diff --git a/Inputs/EventGenerator/3He.source b/Inputs/EventGenerator/3He.source index d85a5bce8286662497a458a9fad9ce991b45c15f..b80e9a4a1fec2107cdf38b2507847f66ca9022fa 100644 --- a/Inputs/EventGenerator/3He.source +++ b/Inputs/EventGenerator/3He.source @@ -5,9 +5,9 @@ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Isotropic EnergyLow= 0 - EnergyHigh= 520 + EnergyHigh= 710 HalfOpenAngleMin= 0 - HalfOpenAngleMax= 90 + HalfOpenAngleMax= 10 x0= 0 y0= 0 z0= 0 diff --git a/Inputs/EventGenerator/alpha.source b/Inputs/EventGenerator/alpha.source index efbad0a97245a4b50c5111e894a5886a989f3aa8..9c9bc63fd260253abb84404ccf4e4b31a487ef90 100644 --- a/Inputs/EventGenerator/alpha.source +++ b/Inputs/EventGenerator/alpha.source @@ -5,9 +5,9 @@ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Isotropic EnergyLow= 0 - EnergyHigh= 585 + EnergyHigh= 795 HalfOpenAngleMin= 0 - HalfOpenAngleMax= 90 + HalfOpenAngleMax= 10 x0= 0 y0= 0 z0= 0 diff --git a/Inputs/EventGenerator/deuton.source b/Inputs/EventGenerator/deuton.source index 2730c05f1604e224d664a63851a9296419a93228..7a864b0fb3018a654583b78a9f4ab1e2d1e00a19 100644 --- a/Inputs/EventGenerator/deuton.source +++ b/Inputs/EventGenerator/deuton.source @@ -5,9 +5,9 @@ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Isotropic EnergyLow= 0 - EnergyHigh= 250 + EnergyHigh= 265 HalfOpenAngleMin= 0 - HalfOpenAngleMax= 90 + HalfOpenAngleMax= 10 x0= 0 y0= 0 z0= 0 diff --git a/Inputs/EventGenerator/proton.source b/Inputs/EventGenerator/proton.source index 9bb49875926c71e424ea7975009273cdd77649b2..7e08b60bc93daf9b29800306fbd9f142bad79b03 100644 --- a/Inputs/EventGenerator/proton.source +++ b/Inputs/EventGenerator/proton.source @@ -5,9 +5,9 @@ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Isotropic EnergyLow= 0 - EnergyHigh= 190 + EnergyHigh= 195 HalfOpenAngleMin= 0 - HalfOpenAngleMax= 90 + HalfOpenAngleMax= 10 x0= 0 y0= 0 z0= 0 diff --git a/Inputs/EventGenerator/triton.source b/Inputs/EventGenerator/triton.source index a3430584ab09a241a519a1340fb6b5a5cc1fc32b..f67b1d5fb4d199a7907c33275ed56469da4ba73c 100644 --- a/Inputs/EventGenerator/triton.source +++ b/Inputs/EventGenerator/triton.source @@ -5,9 +5,9 @@ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Isotropic EnergyLow= 0 - EnergyHigh= 280 + EnergyHigh= 315 HalfOpenAngleMin= 0 - HalfOpenAngleMax= 90 + HalfOpenAngleMax= 10 x0= 0 y0= 0 z0= 0 diff --git a/NPLib/Physics/NPNucleus.cxx b/NPLib/Physics/NPNucleus.cxx index f2fb75b9869b6a08874e55905387926f40d2d7df..84af50d3690b86f265d902af5e4c51339dc02383 100644 --- a/NPLib/Physics/NPNucleus.cxx +++ b/NPLib/Physics/NPNucleus.cxx @@ -348,6 +348,37 @@ double Nucleus::GetEnergyCM(double EnergyLab, double ThetaLab, double PhiLab, do return EnergyCM; } +//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... +double Nucleus::GetThetaCM(double EnergyLab, double ThetaLab, double PhiLab, double relativisticboost) +{ + SetKineticEnergy(EnergyLab); + double EnergyCM; + double ThetaCM; + double ImpulsionLab; + double ImpulsionLabX; + double ImpulsionLabY; + double ImpulsionLabZ; + TVector3 VImpulsionLAB; + TLorentzVector LVEnergyImpulsionLAB; + TLorentzVector LVEnergyImpulsionCM; + + ImpulsionLab = sqrt(EnergyLab*EnergyLab + 2*EnergyLab*Mass()); + ImpulsionLabX = ImpulsionLab*sin(ThetaLab)*cos(PhiLab); + ImpulsionLabY = ImpulsionLab*sin(ThetaLab)*sin(PhiLab); + ImpulsionLabZ = ImpulsionLab*cos(ThetaLab); + + VImpulsionLAB = TVector3(ImpulsionLabX, ImpulsionLabY, ImpulsionLabZ); + LVEnergyImpulsionLAB = TLorentzVector(VImpulsionLAB,Mass()+EnergyLab); + + LVEnergyImpulsionCM = LVEnergyImpulsionLAB; + LVEnergyImpulsionCM.Boost(0,0,-relativisticboost); + + EnergyCM = LVEnergyImpulsionCM.Energy() - Mass(); + ThetaCM = LVEnergyImpulsionCM.Theta(); + + return ThetaCM; +} + diff --git a/NPLib/Physics/NPNucleus.h b/NPLib/Physics/NPNucleus.h index 015bda45d34f4b1459cff7a7f8f9aedd15c0c9fd..f4f8678739d8d92cc41c14f41d7f5fe413d0aefa 100644 --- a/NPLib/Physics/NPNucleus.h +++ b/NPLib/Physics/NPNucleus.h @@ -123,6 +123,7 @@ namespace NPL { BetaToVelocity();} void SetBeta(double beta) {fBeta = beta; BetaToGamma(); BetaToEnergy(); EnergyToTof(); EnergyToBrho();BetaToVelocity();} double GetEnergyCM(double EnergyLab, double ThetaLab, double PhiLab, double relativisticboost); + double GetThetaCM(double EnergyLab, double ThetaLab, double PhiLab, double relativisticboost); // Nuclear mass in MeV double Mass() const {return (fAtomicWeight*amu_c2 + fMassExcess/1000. - fCharge*electron_mass_c2);} diff --git a/NPSimulation/Core/PhysicsList.hh b/NPSimulation/Core/PhysicsList.hh index f263cf6388e805e6de195f65269d970f86a6cc84..8e0cd5f81df29e27452eeda07c32ecd73621a560 100644 --- a/NPSimulation/Core/PhysicsList.hh +++ b/NPSimulation/Core/PhysicsList.hh @@ -85,7 +85,7 @@ class PhysicsList: public G4VModularPhysicsList{ G4VPhysicsConstructor* emPhysicsList; G4VPhysicsConstructor* decay_List; G4VPhysicsConstructor* radioactiveDecay_List; - + private: // Physics option std::string m_EmList; double m_IonBinaryCascadePhysics; diff --git a/NPSimulation/Detectors/Hira/Hira.cc b/NPSimulation/Detectors/Hira/Hira.cc index b5d08b82082deb49c898d8508985072429301ef8..c49fde8854be64d313ac843e118156563df2a5ea 100644 --- a/NPSimulation/Detectors/Hira/Hira.cc +++ b/NPSimulation/Detectors/Hira/Hira.cc @@ -68,15 +68,16 @@ Hira::Hira(){ InitializeMaterial(); m_EventHira = new THiraData(); - // Dark Grey - m_SiliconVisAtt = new G4VisAttributes(G4Colour(0.3, 0.3, 0.3)) ; - // Green - m_CsIVisAtt = new G4VisAttributes(G4Colour(0.2, 0.5, 0.2)) ; + // Silicon + m_SiliconVisAtt = new G4VisAttributes(G4Colour(0.839216, 0.839216, 0.839216)) ; + m_SiliconVisAtt2 = new G4VisAttributes(G4Colour(0.839216, 0.839216, 0.839216)) ; + // CsI Color + m_CsIVisAtt = new G4VisAttributes(G4Colour(0.529412, 0.807843, 0.980392, 0.9)) ; //m_CsIVisAtt->SetForceWireframe(true); - m_LogicThinSi = 0; - m_LogicThickSi = 0; - m_LogicCsICrystal = 0; - m_LogicCluster = 0; + m_LogicThinSi = 0; + m_LogicThickSi = 0; + m_LogicCsICrystal = 0; + m_LogicCluster = 0; } //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... Hira::~Hira(){ @@ -523,7 +524,7 @@ void Hira::VolumeMaker(G4int DetectorNumber, m_LogicThickSi->SetSensitiveDetector(m_ThickSiStripScorer); // Visualisation of ThickSi - m_LogicThickSi->SetVisAttributes(m_SiliconVisAtt) ; + m_LogicThickSi->SetVisAttributes(m_SiliconVisAtt2) ; } /////////////////////////////////////////////////// @@ -563,8 +564,8 @@ void Hira::VolumeMaker(G4int DetectorNumber, // Sub Mother Volume G4Trd* solidCluster = new G4Trd("SolidCluster", 0.5*ClusterFaceFront,0.5*ClusterFaceBack,0.5*ClusterFaceFront,0.5*ClusterFaceBack, 0.5*CsIThickness); m_LogicCluster = new G4LogicalVolume(solidCluster, m_MaterialVacuum, "LogicSolidCluster", 0, 0, 0); - //m_LogicCluster->SetVisAttributes(G4VisAttributes::Invisible); - G4VisAttributes* TempVisAtt = new G4VisAttributes(G4Colour(0.6, 0.6, 0.3)) ; + m_LogicCluster->SetVisAttributes(G4VisAttributes::Invisible); + G4VisAttributes* TempVisAtt = new G4VisAttributes(G4Colour(0.415686, 0.352941, 0.803922, 0.1)) ; TempVisAtt->SetForceWireframe(true); m_LogicCluster->SetVisAttributes(TempVisAtt); diff --git a/NPSimulation/Detectors/Hira/Hira.hh b/NPSimulation/Detectors/Hira/Hira.hh index 56edbb193f420a2e31eeed05eed376615b6a72c6..b198cbb29fc1c2fa3f4dd8fc7fb9b4631549d338 100644 --- a/NPSimulation/Detectors/Hira/Hira.hh +++ b/NPSimulation/Detectors/Hira/Hira.hh @@ -222,8 +222,9 @@ private: private:/// Visualisation Attribute: - G4VisAttributes* m_SiliconVisAtt ; - G4VisAttributes* m_CsIVisAtt; + G4VisAttributes* m_SiliconVisAtt ; + G4VisAttributes* m_SiliconVisAtt2 ; + G4VisAttributes* m_CsIVisAtt; public: static NPS::VDetector* Construct(); diff --git a/NPSimulation/PhysicsListOption.txt b/NPSimulation/PhysicsListOption.txt index 6e6fc7c040fc7fe15a6a81a2c6400d52d2ec1711..42fbf0b5f36097997a801301737c23b68c197e38 100644 --- a/NPSimulation/PhysicsListOption.txt +++ b/NPSimulation/PhysicsListOption.txt @@ -1,10 +1,10 @@ EmPhysicsList Option4 DefaultCutOff 1 -IonBinaryCascadePhysics 0 +IonBinaryCascadePhysics 1 EmExtraPhysics 0 HadronElasticPhysics 0 StoppingPhysics 0 OpticalPhysics 0 -HadronPhysicsQGSP_BIC_HP 0 -Decay 0 +HadronPhysicsQGSP_BIC_HP 1 +Decay 1 diff --git a/NPSimulation/ressources/macro/vis.mac b/NPSimulation/ressources/macro/vis.mac index bb1c0b15b2d80e764a118d4db75e5a2ac406d0c8..9ce524983b522db1fef7359b3003dde98915af13 100644 --- a/NPSimulation/ressources/macro/vis.mac +++ b/NPSimulation/ressources/macro/vis.mac @@ -47,7 +47,8 @@ /vis/viewer/set/autoRefresh true /vis/verbose 0 -/vis/viewer/set/background black +#/vis/viewer/set/background black +/vis/viewer/set/background white # print Option #/vis/ogl/set/printMode vectored diff --git a/Projects/Hira2/Analysis.cxx b/Projects/Hira2/Analysis.cxx deleted file mode 100644 index 60d0f38879e42fea579f06eed13510ef4407db9b..0000000000000000000000000000000000000000 --- a/Projects/Hira2/Analysis.cxx +++ /dev/null @@ -1,244 +0,0 @@ - /***************************************************************************** - * Copyright (C) 2009-2014 this file is part of the NPTool Project * - * * - * For the licensing terms see $NPTOOL/Licence/NPTool_Licence * - * For the list of contributors see $NPTOOL/Licence/Contributors * - *****************************************************************************/ - -/***************************************************************************** - * Original Author: Adrien MATTA contact address: a.matta@surrey.ac.uk * - * * - * Creation Date : march 2025 * - * Last update : * - *---------------------------------------------------------------------------* - * Decription: * - * Class describing the property of an Analysis object * - * * - *---------------------------------------------------------------------------* - * Comment: * - * * - * * - *****************************************************************************/ -#include<iostream> -using namespace std; -#include"Analysis.h" -#include"NPAnalysisFactory.h" -#include"NPDetectorManager.h" -#include"NPOptionManager.h" -#include"RootOutput.h" -#include"RootInput.h" -//////////////////////////////////////////////////////////////////////////////// -Analysis::Analysis(){ -} -//////////////////////////////////////////////////////////////////////////////// -Analysis::~Analysis(){ -} - -//////////////////////////////////////////////////////////////////////////////// -void Analysis::Init(){ - Hira = (THiraPhysics*) m_DetectorManager->GetDetector("HIRAArray"); - InitialConditions=new TInitialConditions(); - InitOutputBranch(); - InitInputBranch(); - Rand = TRandom3(); - TransferReaction= new NPL::Reaction("8B(p,4He)5Be@400"); - //TransferReaction->ReadConfigurationFile(NPOptionManager::getInstance()->GetReactionFile()); - - E_ThinSi = 0; - E_ThickSi = 0; - E_CsI = 0; - ELab = 0; - ThetaLab = 0; - PhiLab = 0; - ThetaLab_simu = 0; - ThetaCM = 0; - ExcitationEnergy = 0; - X = 0; - Y = 0; - Z = 0; - TelescopeNumber = 0; - EnergyThreshold = 1; - - event = 0; - good_event=0; - ProtonEnergy = 50; - - TargetThickness = m_DetectorManager->GetTargetThickness()*micrometer; - // Energy loss table: the G4Table are generated by the simulation - Proton_CsI = EnergyLoss("proton_CsI.G4table","G4Table",100 ); - Triton_CD2 = EnergyLoss("triton_CD2.G4table","G4Table",100 ); - Triton_CH2 = EnergyLoss("triton_CD2.G4table","G4Table",100 ); - Deuteron_CH2 = EnergyLoss("deuteron_CH2.G4table","G4Table",100 ); - He3_CD2 = EnergyLoss("He3_CD2.G4table","G4Table",100 ); - He4_CH2 = EnergyLoss("alpha_CH2.G4table","G4Table",100 ); - - f_proton = new TF1("f_proton","1 - TMath::Exp(-(x-0.0538203)/28.3125)",0,10); - -} - -//////////////////////////////////////////////////////////////////////////////// -void Analysis::TreatEvent(){ - // Reinitiate calculated variable - ReInitValue(); - - - //double BeamEnergy = Rand.Gaus(Initial->GetIncidentInitialKineticEnergy(),4.5); - double BeamEnergy = InitialConditions->GetIncidentInitialKineticEnergy(); - double Thickness; - - - //cout << Thickness << endl; - //cout << eval << endl; - - //TransferReaction->SetBeamEnergy(BeamEnergy); - //////////////////////////// LOOP on Hira Hit ////////////////// - if(Hira -> ThickSi_E.size() == 1){ - for(unsigned int countHira = 0 ; countHira < Hira->ThickSi_E.size() ; countHira++){ - event += 1; - TelescopeNumber = Hira->TelescopeNumber[countHira]; - - TargetThickness = m_DetectorManager->GetTargetThickness(); - - X = Hira->GetPositionOfInteraction(0).X(); - Y = Hira->GetPositionOfInteraction(0).Y(); - Z = Hira->GetPositionOfInteraction(0).Z(); - - TVector3 PositionOnHira = TVector3(X,Y,Z); - TVector3 ZUnit = TVector3(0,0,1); - - X_target = InitialConditions->GetIncidentPositionX();// + Rand.Gaus(0,2.); - Y_target = InitialConditions->GetIncidentPositionY();// + Rand.Gaus(0,2.); - double Z_target = InitialConditions->GetIncidentPositionZ(); - - //TVector3 PositionOnTarget = TVector3(X_target,Y_target,Z_target); - TVector3 PositionOnTarget = TVector3(0,0,0); - TVector3 HitDirection = PositionOnHira - PositionOnTarget; - TVector3 HitDirectionUnit = HitDirection.Unit(); - - /*TVector3 BeamDirection = InitialConditions->GetBeamDirection(); - double XBeam = BeamDirection.X(); - double YBeam = BeamDirection.Y(); - double ZBeam = BeamDirection.Z();*/ - - //ThetaLab = BeamDirection.Angle(HitDirection); - ThetaLab = ZUnit.Angle(HitDirection); - ThetaLab = ThetaLab/deg; - - PhiLab = HitDirection.Phi()/deg; - //if(PhiLab>0) PhiLab = 180-PhiLab; - //if(PhiLab<0) PhiLab = -PhiLab-180; - - E_ThickSi = Hira->ThickSi_E[countHira]; - E_ThinSi = Hira->ThinSi_E[countHira]; - //ELab = E_ThinSi + E_ThickSi; - if(Hira->CsI_E.size() == 1){ - for(int countCsI =0; countCsI<Hira->CsI_E.size(); countCsI++){ - //Try to simulate the nuclear reaction loss - Thickness = Proton_CsI.EvaluateMaterialThickness(0*MeV, Hira->CsI_E[countCsI]*MeV, 200*millimeter, 0.1*millimeter); - //cout << Thickness << endl; - double eval = f_proton->Eval(Thickness/10); - double Random_value = Rand.Uniform(0,1); - - //E_CsI = Hira->CsI_E[countCsI]; - if(Random_value>eval)E_CsI = Hira->CsI_E[countCsI]; - else E_CsI = Rand.Uniform(0,Hira->CsI_E[countCsI]); - - //ELab += E_CsI; - if(E_CsI>EnergyThreshold) ELab = E_ThinSi + E_ThickSi + E_CsI; - } - } - else ELab = E_ThinSi + E_ThickSi; - - if(ELab>0){ - //ELab = Deuteron_CH2.EvaluateInitialEnergy(ELab,(double)TargetThickness/2*micrometer,0); - //ELab = Triton_CH2.EvaluateInitialEnergy(ELab,(double)TargetThickness/2*micrometer,0); - ELab = He4_CH2.EvaluateInitialEnergy(ELab,(double)TargetThickness/2*micrometer,0); - //cout << TargetThickness << endl; - // ********************** Angle in the CM frame ***************************** - TransferReaction -> SetNuclei3(ELab, ThetaLab*deg); - ThetaCM = TransferReaction->GetThetaCM()/deg; - ExcitationEnergy = TransferReaction->GetExcitation4(); - } - else{ - ThetaCM = -1000; - ExcitationEnergy = -1000; - } - }//end loop Hira - }//end if Hira - - -} - -//////////////////////////////////////////////////////////////////////////////// -void Analysis::End(){ - /*cout << "Incident Energy = " << ProtonEnergy << " MeV" << endl; - cout << "event = " << event << endl; - cout << "Good event = " << good_event << endl; - cout << "Efficiency = " << good_event/event*100 << "%" << endl;*/ -} - -//////////////////////////////////////////////////////////////////////////////// -void Analysis::InitOutputBranch() { - RootOutput::getInstance()->GetTree()->Branch( "ThetaLab" , &ThetaLab , "ThetaLab/D" ); - RootOutput::getInstance()->GetTree()->Branch( "PhiLab" , &PhiLab , "PhiLab/D" ) ; - RootOutput::getInstance()->GetTree()->Branch("ThetaCM", &ThetaCM,"ThetaCM/D") ; - RootOutput::getInstance()->GetTree()->Branch( "E_ThinSi" , &E_ThinSi , "E_ThinSi/D" ) ; - RootOutput::getInstance()->GetTree()->Branch( "E_ThickSi" , &E_ThickSi , "E_ThickSi/D" ) ; - RootOutput::getInstance()->GetTree()->Branch( "E_CsI" , &E_CsI , "E_CsI/D" ) ; - RootOutput::getInstance()->GetTree()->Branch( "ELab" , &ELab , "ELab/D" ) ; - RootOutput::getInstance()->GetTree()->Branch("ExcitationEnergy", &ExcitationEnergy,"ExcitationEnergy/D") ; - RootOutput::getInstance()->GetTree()->Branch( "X" , &X , "X/D" ) ; - RootOutput::getInstance()->GetTree()->Branch( "Y" , &Y , "Y/D" ) ; - RootOutput::getInstance()->GetTree()->Branch( "Z" , &Z , "Z/D" ) ; - RootOutput::getInstance()->GetTree()->Branch( "TelescopeNumber" , &TelescopeNumber , "TelescopeNumber/D" ) ; - RootOutput::getInstance()->GetTree()->Branch( "X_target" , &X_target , "X_target/D" ) ; - RootOutput::getInstance()->GetTree()->Branch( "Y_target" , &Y_target , "Y_target/D" ) ; - //RootOutput::getInstance()->GetTree()-> Branch("InteractionCoordinates","TInteractionCoordinates",&InteractionCoordinates); - //RootOutput::getInstance()->GetTree()->Branch("InitialConditions","TInitialConditions",&InitialConditions); -} - -//////////////////////////////////////////////////////////////////////////////// -void Analysis::InitInputBranch(){ - RootInput:: getInstance()->GetChain()->SetBranchStatus("InitialConditions",true ); - RootInput:: getInstance()->GetChain()->SetBranchStatus("fIC_*",true ); - RootInput:: getInstance()->GetChain()->SetBranchAddress("InitialConditions",&InitialConditions); -} - -//////////////////////////////////////////////////////////////////////////////// -void Analysis::ReInitValue(){ - E_ThinSi = -1000; - E_ThickSi = -1000; - E_CsI = -1000; - ELab = -1000; - ThetaLab = -1000; - PhiLab = -1000; - ThetaCM = -1000; - ExcitationEnergy = -1000; - X = -1000; - Y = -1000; - Z = -1000; - TelescopeNumber = -1; - X_target = -1000; - Y_target = -1000; -} - -//////////////////////////////////////////////////////////////////////////////// -// Construct Method to be pass to the DetectorFactory // -//////////////////////////////////////////////////////////////////////////////// -NPL::VAnalysis* Analysis::Construct(){ - return (NPL::VAnalysis*) new Analysis(); -} - -//////////////////////////////////////////////////////////////////////////////// -// Registering the construct method to the factory // -//////////////////////////////////////////////////////////////////////////////// -extern "C"{ - class proxy{ - public: - proxy(){ - NPL::AnalysisFactory::getInstance()->SetConstructor(Analysis::Construct); - } - }; - - proxy p; -} \ No newline at end of file diff --git a/Projects/Hira2/Analysis.h b/Projects/Hira2/Analysis.h deleted file mode 100644 index 3263899d6cd1e9737e3fc5031d931141b5828859..0000000000000000000000000000000000000000 --- a/Projects/Hira2/Analysis.h +++ /dev/null @@ -1,87 +0,0 @@ -#ifndef Analysis_h -#define Analysis_h -/***************************************************************************** - * Copyright (C) 2009-2014 this file is part of the NPTool Project * - * * - * For the licensing terms see $NPTOOL/Licence/NPTool_Licence * - * For the list of contributors see $NPTOOL/Licence/Contributors * - *****************************************************************************/ - -/***************************************************************************** - * Original Author: Adrien MATTA contact address: a.matta@surrey.ac.uk * - * * - * Creation Date : march 2025 * - * Last update : * - *---------------------------------------------------------------------------* - * Decription: * - * Class describing the property of an Analysis object * - * * - *---------------------------------------------------------------------------* - * Comment: * - * * - * * - *****************************************************************************/ -#include"NPVAnalysis.h" -#include"THiraPhysics.h" -#include "TInitialConditions.h" -#include "TInteractionCoordinates.h" -#include "NPEnergyLoss.h" -#include "NPReaction.h" -#include "TRandom3.h" -#include "TF1.h" -class Analysis: public NPL::VAnalysis{ -public: - Analysis(); - ~Analysis(); - -public: - void Init(); - void TreatEvent(); - void End(); - void InitOutputBranch(); - void InitInputBranch(); - void ReInitValue(); - static NPL::VAnalysis* Construct(); - - double event; - double good_event; - double ProtonEnergy; - -private: - double TargetThickness; - double ExcitationEnergy; - double ELab; - double E_ThinSi; - double E_ThickSi; - double E_CsI; - double PhiLab; - double ThetaLab; - double ThetaLab_simu; - double ThetaCM; - double X,Y,Z; - double TelescopeNumber; - double EnergyThreshold; - double X_target; - double Y_target; - - - - NPL::Reaction* TransferReaction; - - // intermediate variable - TRandom3 Rand; - - - TF1* f_proton; - NPL::EnergyLoss Proton_CsI ; - NPL::EnergyLoss Triton_CD2 ; - NPL::EnergyLoss Triton_CH2 ; - NPL::EnergyLoss Deuteron_CH2 ; - NPL::EnergyLoss He3_CD2; - NPL::EnergyLoss He4_CH2; - - THiraPhysics* Hira; - TInitialConditions* InitialConditions; - TInteractionCoordinates* InteractionCoordinates; -}; -#endif \ No newline at end of file diff --git a/Projects/Hira2/CMakeLists.txt b/Projects/Hira2/CMakeLists.txt deleted file mode 100644 index 6806d778cdd5240538c93278f2ac5dcce3429083..0000000000000000000000000000000000000000 --- a/Projects/Hira2/CMakeLists.txt +++ /dev/null @@ -1,2 +0,0 @@ -cmake_minimum_required (VERSION 2.8) -include("../../NPLib/ressources/CMake/NPAnalysis.cmake") diff --git a/Projects/Hira2/RunToTreat.txt b/Projects/Hira2/RunToTreat.txt deleted file mode 100755 index f91784f3c00f836b5e922f7e5282c9b18daae412..0000000000000000000000000000000000000000 --- a/Projects/Hira2/RunToTreat.txt +++ /dev/null @@ -1,6 +0,0 @@ -TTreeName - SimulatedTree -RootFileName - ../../Outputs/Simulation/dumb.root - %../../Outputs/Simulation/12Be_pt_iso.root - %../../Outputs/Simulation/12Be_pt_1st.root diff --git a/Projects/Lassa/Analysis.cxx b/Projects/Lassa/Analysis.cxx index 32ab354ce799b826dc95a9b8660f6e42613fdd5d..7b40f5a52b8986cefb9fa650fdf3f429ebf5caf3 100644 --- a/Projects/Lassa/Analysis.cxx +++ b/Projects/Lassa/Analysis.cxx @@ -82,7 +82,13 @@ void Analysis::TreatEvent(){ double ParticleEnergy = InitialConditions->GetKineticEnergy(0); double EDelta = 2.0; - if(Lassa->ThickSi_E.size()>0) InitialEnergy = ParticleEnergy; + InitialEnergy = ParticleEnergy; + //if(Lassa->ThickSi_E.size()>0) InitialEnergy = ParticleEnergy; + double phi_in = acos(InitialConditions->GetMomentumDirectionX(0)/sin(InitialConditions->GetThetaCM(0)*deg)); + if(InitialEnergy>0){ + ECM_initial = proton->GetEnergyCM(InitialEnergy, InitialConditions->GetThetaCM(0)*deg, phi_in, BetaCM); + } + else ECM_initial = -100; ///////////////////////////LOOP on Lassa Hit////////////////////////////////// if(Lassa->ThickSi_E.size() == 1){ detectedEvents++; @@ -142,14 +148,10 @@ void Analysis::TreatEvent(){ } - if(fabs(ParticleEnergy-ELab)<EDelta){ - peakEvents++; - } - else{ + if(fabs(ParticleEnergy-ELab)>EDelta){ ELab = -100; - //E_CsI = -100; } - + if(fabs(ParticleEnergy-ELab_nucl)>EDelta) ELab_nucl = -100; if(ELab>0){ @@ -157,11 +159,6 @@ void Analysis::TreatEvent(){ } else ECM = -100; - if(InitialEnergy>0){ - ECM_initial = proton->GetEnergyCM(InitialEnergy, ThetaLab, PhiLab, BetaCM); - } - else ECM = -100; - ThetaLab = ThetaLab/deg; PhiLab = PhiLab/deg; } @@ -169,7 +166,7 @@ void Analysis::TreatEvent(){ //////////////////////////////////////////////////////////////////////////////// void Analysis::End(){ - geomEff = 100*((double)detectedEvents)/((double)totalEvents); + /*geomEff = 100*((double)detectedEvents)/((double)totalEvents); peakEff = 100*((double)peakEvents)/((double)detectedEvents); @@ -179,7 +176,7 @@ void Analysis::End(){ cout << "PeakEvents: " << peakEvents << endl; cout << "Geometric Efficiency: " << geomEff << endl; - cout << "Peak Efficiency: " << peakEff << endl; + cout << "Peak Efficiency: " << peakEff << endl;*/ } diff --git a/Projects/Lassa/RunToTreat.txt b/Projects/Lassa/RunToTreat.txt index 65eaf703b21cf251646179abb994e90ba46099c5..82cfe37dc41941c8e5c2033e7b67ec8e194c5112 100755 --- a/Projects/Lassa/RunToTreat.txt +++ b/Projects/Lassa/RunToTreat.txt @@ -1,7 +1,7 @@ TTreeName SimulatedTree RootFileName - ../../Outputs/Simulation/lassa_proton_pbuu.root + ../../Outputs/Simulation/lassa_proton.root %../../Outputs/Simulation/lassa_pid_d.root %../../Outputs/Simulation/lassa_pid_t.root