diff --git a/Inputs/DetectorConfiguration/Riken_65mm.detector b/Inputs/DetectorConfiguration/Riken_65mm.detector
index 0855d67ee4c0fd72fe87854cf6d1fbdfc055f31f..5dd3b2d0cff7b7899e5d9f287d35f3aa619fe0dc 100644
--- a/Inputs/DetectorConfiguration/Riken_65mm.detector
+++ b/Inputs/DetectorConfiguration/Riken_65mm.detector
@@ -5,10 +5,10 @@
%Material name according to the target library
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-GeneralTarget
+%GeneralTarget
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Target
- THICKNESS= 10000
+ THICKNESS= 1
RADIUS= 45
MATERIAL= CD2
ANGLE= 45
diff --git a/Inputs/EventGenerator/10He.reaction b/Inputs/EventGenerator/10He.reaction
index 2c1407b0f4d2339fef7a3ead00cc6e449d203b9a..1987df73649f40bdf632faad34251b0c2a9d64b3 100644
--- a/Inputs/EventGenerator/10He.reaction
+++ b/Inputs/EventGenerator/10He.reaction
@@ -7,7 +7,7 @@ TransfertToResonance
Target= 2H
Light= 3He
Heavy= 10He
- ExcitationEnergy= 2.0
+ ExcitationEnergy= 0.0
BeamEnergy= 550
BeamEnergySpread= 0
SigmaThetaX= 0.6921330164
@@ -19,7 +19,7 @@ TransfertToResonance
ResonanceDecayA= 8
CrossSectionPath= 11Li(d,3He)10He.txt
ShootLight= 1
- ShootHeavy= 1
+ ShootHeavy= 0
ShootDecayProduct= 0
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
diff --git a/NPAnalysis/10He_Riken/Analysis b/NPAnalysis/10He_Riken/Analysis
index ef3df396e03189c2958aeb881d7d327c569c9ee1..75f4cb005797b9560b71bcab202eac2dbc49b5fb 100755
Binary files a/NPAnalysis/10He_Riken/Analysis and b/NPAnalysis/10He_Riken/Analysis differ
diff --git a/NPAnalysis/10He_Riken/include/ObjectManager.hh b/NPAnalysis/10He_Riken/include/ObjectManager.hh
index c61e480db4743573a1b0290c7757247286083df2..ab0d7e99fb64cf0ed8893087375509fb786af239 100644
--- a/NPAnalysis/10He_Riken/include/ObjectManager.hh
+++ b/NPAnalysis/10He_Riken/include/ObjectManager.hh
@@ -108,40 +108,40 @@ namespace ENERGYLOSS
// 3He Energy Loss
EnergyLoss He3TargetWind = EnergyLoss ( "3He_Mylar.txt" ,
- 1000 ,
- 3 ,
+ 10000 ,
+ 1 ,
3 );
- EnergyLoss He3TargetGaz = EnergyLoss ( "3He_D2solid_1b_26K.txt" ,
- 1000 ,
- 3 ,
+ EnergyLoss He3TargetGaz = EnergyLoss ( "3He_D2gaz_1b_26K.txt" ,
+ 10000 ,
+ 1 ,
3 );
EnergyLoss He3StripAl = EnergyLoss ( "3He_Al.txt" ,
- 1000 ,
- 3 ,
+ 10000 ,
+ 1 ,
3 );
EnergyLoss He3StripSi = EnergyLoss ( "3He_Si.txt" ,
- 1000 ,
- 3 ,
+ 10000 ,
+ 1 ,
3 );
// proton Energy Loss
EnergyLoss protonTargetWind = EnergyLoss ( "proton_Mylar.txt" ,
1000 ,
- 3 ,
+ 1 ,
1 );
EnergyLoss protonTargetGaz = EnergyLoss ( "proton_D2gaz_1b_26K.txt" ,
1000 ,
- 3 ,
+ 1 ,
1 );
EnergyLoss protonStripAl = EnergyLoss ( "proton_Al.txt" ,
100 ,
- 3 ,
+ 1 ,
1 );
diff --git a/NPAnalysis/10He_Riken/src/Analysis.cc b/NPAnalysis/10He_Riken/src/Analysis.cc
index 09c9e496b761fb4e4ee10544d5a6787c2396ac46..4015e62a980c68cb2fa625e272695163ccd5744b 100644
--- a/NPAnalysis/10He_Riken/src/Analysis.cc
+++ b/NPAnalysis/10He_Riken/src/Analysis.cc
@@ -48,6 +48,7 @@ int main(int argc,char** argv)
//E lab et Theta lab
RootOutput::getInstance()->GetTree()->Branch("ThetaCM",ThetaCM,"ThetaCM[2]/D") ;
RootOutput::getInstance()->GetTree()->Branch("ELab",ELab,"ELab[2]/D") ;
+ RootOutput::getInstance()->GetTree()->Branch("ThetaLab",ThetaLab,"ThetaLab[2]/D") ;
RootOutput::getInstance()->GetTree()->Branch("X",X,"X/D[2]") ;
RootOutput::getInstance()->GetTree()->Branch("Y",Y,"Y/D[2]") ;
@@ -110,7 +111,7 @@ int main(int argc,char** argv)
-/* // Target (from initial condition)
+ // Target (from initial condition)
XTarget = Init->GetICPositionX(0);
YTarget = Init->GetICPositionY(0);
// XTarget = RandomEngine.Gaus(Init->GetICPositionX(0),1);
@@ -120,9 +121,11 @@ int main(int argc,char** argv)
////
// Must 2 And ThinSi //
- for(int hit = 0; hit < M2 -> GetEventMultiplicity() ; hit ++)
+ //for(int hit = 0; hit < M2 -> GetEventMultiplicity() ; hit ++)
+ for(int hit = 0; hit < M2 -> Si_E.size() ; hit ++)
{
- // ELab[hit] = M2 -> GetEnergyDeposit(hit);
+ ELab[hit] = -1 ; ThetaLab[hit] = -1;
+ // Get Hit Direction
TVector3 HitDirection = M2 -> GetPositionOfInteraction(hit) - TVector3(XTarget,YTarget,0);
// Angle between beam and particle
@@ -132,67 +135,77 @@ int main(int argc,char** argv)
// ANgle between particule and Must2 Si surface
double ThetaMM2Surface = ThetaCalculation ( HitDirection , M2 -> GetTelescopeNormal(hit) );
- if(M2 -> GetPositionOfInteraction(hit).Z()>0)
+ if(M2 -> GetPositionOfInteraction(hit).Z() > 0)
{
- if(ELab[hit]>-1000 && ThinSi_E>0 )
- if(M2Physics.CsI.size)
- {
- ELab[hit]= He3StripAl.EvaluateInitialEnergy( ELab[hit] , // Energy of the detected particle
- 2*0.4*micrometer , // Target Thickness at 0 degree
- ThetaMM2Surface );
-
-
- ELab[hit]= He3StripSi.EvaluateInitialEnergy( ELab[hit] , // Energy of the detected particle
- 20*micrometer , // Target Thickness at 0 degree
- ThetaMM2Surface );
+ if( M2 -> CsI_E[hit] == 0 && M2 -> Si_E[hit] > 0)
+ {
+ ELab[hit] = M2 -> Si_E[hit] ;
+ ELab[hit]= He3StripAl.EvaluateInitialEnergy( ELab[hit] , // Energy of the detected particle
+ 2*0.4*micrometer , // Target Thickness at 0 degree
+ ThetaMM2Surface );
- ELab[hit]= He3StripAl.EvaluateInitialEnergy( ELab[hit] , // Energy of the detected particle
- 0.4*micrometer , // Target Thickness at 0 degree
- ThetaMM2Surface );
+
+// ELab[hit]= He3StripSi.EvaluateInitialEnergy( ELab[hit] , // Energy of the detected particle
+// 20*micrometer , // Target Thickness at 0 degree
+// ThetaMM2Surface );
+
+ if(ThinSi -> Energy.size() > 0)ELab[hit] += ThinSi-> Energy[hit];
+
+ ELab[hit]= He3StripAl.EvaluateInitialEnergy( ELab[hit] , // Energy of the detected particle
+ 0.4*micrometer , // Target Thickness at 0 degree
+ ThetaMM2Surface );
- ELab[hit]= He3TargetWind.EvaluateInitialEnergy( ELab[hit] , // Energy of the detected particle
- 15*micrometer , // Target Thickness at 0 degree
- ThetaN );
+ ELab[hit]= He3TargetWind.EvaluateInitialEnergy( ELab[hit] , // Energy of the detected particle
+ 15*micrometer , // Target Thickness at 0 degree
+ ThetaN );
- ELab[hit]= He3TargetGaz.EvaluateInitialEnergy( ELab[hit] , // Energy of the detected particle
- 1.5*mm , // Target Thickness at 0 degree
- ThetaN );
+ ELab[hit]= He3TargetGaz.EvaluateInitialEnergy( ELab[hit] , // Energy of the detected particle
+ 1.5*mm , // Target Thickness at 0 degree
+ ThetaN );
- ThetaCM[hit] = myReaction -> EnergyLabToThetaCM( ELab[hit] ) /deg ;
- ExcitationEnergy[hit] = myReaction -> ReconstructRelativistic( ELab[hit] , ThetaLab[hit] ) ;
- X[hit] = HitDirection . X();
+ ThetaCM[hit] = He10Reaction -> EnergyLabToThetaCM( ELab[hit] ) /deg ;
+ ExcitationEnergy[hit] = He10Reaction -> ReconstructRelativistic( ELab[hit] , ThetaLab[hit] ) ;
X[hit] = HitDirection . X();
Y[hit] = HitDirection . Y();
ThetaLab[hit] = ThetaLab[hit] / deg ;
}
-
- else if(ELab[hit]>-1000 )
+
+ else if(M2 ->CsI_E[hit] > 0 && M2 -> Si_E[hit] > 0)
{
- if(ELab[hit]>21.66)//CsI are inside a Mylar foil, plus rear alu strip
- {
- ELab[hit]= He3TargetWind.EvaluateInitialEnergy( ELab[hit] , // Energy of the detected particle
- 3*micrometer , // Target Thickness at 0 degree
- ThetaMM2Surface );
- ELab[hit]= He3StripAl.EvaluateInitialEnergy( ELab[hit] , // Energy of the detected particle
- 0.4*micrometer , // Target Thickness at 0 degree
- ThetaMM2Surface );
- }
- ELab[hit]= He3StripAl.EvaluateInitialEnergy( ELab[hit] , // Energy of the detected particle
- 0.4*micrometer , // Target Thickness at 0 degree
- ThetaMM2Surface );
+ ELab[hit]= M2 ->CsI_E[hit] ;
- ELab[hit]= He3TargetWind.EvaluateInitialEnergy( ELab[hit] , // Energy of the detected particle
- 15*micrometer , // Target Thickness at 0 degree
- ThetaN );
+ ELab[hit]= He3TargetWind.EvaluateInitialEnergy( ELab[hit] , // Energy of the detected particle
+ 3*micrometer , // Target Thickness at 0 degree
+ ThetaMM2Surface );
+
+ ELab[hit]= He3StripAl.EvaluateInitialEnergy( ELab[hit] , // Energy of the detected particle
+ 0.4*micrometer , // Target Thickness at 0 degree
+ ThetaMM2Surface );
+ ELab[hit]+= M2 ->Si_E[hit];
- ELab[hit]= He3TargetGaz.EvaluateInitialEnergy( ELab[hit] , // Energy of the detected particle
- 1.5*mm , // Target Thickness at 0 degree
- ThetaN );
+ ELab[hit]= He3StripAl.EvaluateInitialEnergy( ELab[hit] , // Energy of the detected particle
+ 0.4*micrometer , // Target Thickness at 0 degree
+ ThetaMM2Surface );
+
+ if(ThinSi -> Energy.size() > 0)ELab[hit] += ThinSi-> Energy[hit];
+
+
+ ELab[hit]= He3StripAl.EvaluateInitialEnergy( ELab[hit] , // Energy of the detected particle
+ 0.4*micrometer , // Target Thickness at 0 degree
+ ThetaMM2Surface );
+
+ ELab[hit]= He3TargetWind.EvaluateInitialEnergy( ELab[hit] , // Energy of the detected particle
+ 15*micrometer , // Target Thickness at 0 degree
+ ThetaN );
+
+ ELab[hit]= He3TargetGaz.EvaluateInitialEnergy( ELab[hit] , // Energy of the detected particle
+ 1.5*mm , // Target Thickness at 0 degree
+ ThetaN );
- ThetaCM[hit]= myReaction -> EnergyLabToThetaCM( ELab[hit] , 1 ) /deg ;
- ExcitationEnergy[hit] = myReaction -> ReconstructRelativistic( ELab[hit], ThetaLab[hit] ) ;
+ ThetaCM[hit]= He10Reaction -> EnergyLabToThetaCM( ELab[hit] ) /deg ;
+ ExcitationEnergy[hit] = He10Reaction -> ReconstructRelativistic( ELab[hit], ThetaLab[hit] ) ;
X[hit] = HitDirection . X();
Y[hit] = HitDirection . Y();
ThetaLab[hit] = ThetaLab[hit] / deg ;
@@ -202,7 +215,7 @@ int main(int argc,char** argv)
}
- else if(M2 -> GetPositionOfInteraction(hit).Z()<0)
+ /*else if(M2 -> GetPositionOfInteraction(hit).Z()<0)
{
if(ELab[hit]>-1000 )
@@ -234,19 +247,11 @@ int main(int argc,char** argv)
else {ExcitationEnergy[hit]=-100 ; X[hit] = -100 ; Y[hit] = -100 ;ThetaLab[hit]=-100; ThetaCM[hit]=-100 ;}
- }
+ } */
}
- // Plastic
- for(int yy = 0 ; yy < Plastic->GetEnergySize() ; yy++)
- {
- if(Plastic->GetPlasticNumber(yy)==1) EPl1[yy]=Plastic->GetEnergy(yy);
- else if(Plastic->GetPlasticNumber(yy)==2) EPl2[yy]=Plastic->GetEnergy(yy);
-
- }*/
-
RootOutput::getInstance()->GetTree()->Fill() ;
}
diff --git a/NPDocumentation/NPTool_UserGuide.tex b/NPDocumentation/NPTool_UserGuide.tex
index 08af1c4665b826214ce58f6595d4c77294cfb3ec..1d9e94b71907bf95457147a170af10c43b55b7ae 100755
--- a/NPDocumentation/NPTool_UserGuide.tex
+++ b/NPDocumentation/NPTool_UserGuide.tex
@@ -92,9 +92,34 @@ It also help to understand what happen during analysis.
NPLib is widely build on top of ROOT,GSL and CLHEP (by order of importance), but not on Geant4.
\section{the directory layout}
-
+ NPLib is made of several folder, one for each detector basicly and some for other shared object such as reaction and nucleus class. A common makefile is present in the NPLib directory calling the individual make file of every folder. In this folder you will also find a README file and the liblist file, containing a small macro used in NPA and NPS makefile to return the list of avaible librairies. Our phylosophy is to have two class for each detector, inheritted fron the TObject ROOT class so theyr can be put in a TTree. The first one is called the Data class, named following the TDetectorNameData, where T is there to remind the TObject dependance. This Data class hold the raw parameter given by the detector during an experiment, understood before any treatment. NPSimulation will output those kind of raw data, the only difference between a simulated and an experiement data object is the one comming from NPS are already calibrated. The second class is a Physics class, named TDetectorNamePhysics. This class hold the results of the treated Data object, ie after application of threshold, selection of good event, application of calibration,... and so one, depending on your detector. Note that the Physics class also derived from the VDetector class, and therefore can be used in NPA via the DetectorManager class. This system allow your annalysis programm to be more flexible and lisible. Let's have a look to the Directory layout:
+
+
+ \begin{itemize}
+ \item[DummyDetector] Template base for making your own Data and Physics class
+ \item[AnnularS1] Containing the Data object for the Annular detector "S1" by micron
+ \item[GASPARD] Data and Physics object associate to the Gamma and Segmented Particule ARray Detector
+ \item[MUST2] Data and Physics object associate to the MUr a STrip 2 detector
+ \item[Plastic] Data and Physics object associate to the Plastic scintillator detector
+ \item[SSSD] Data and Physics object associate to the Single Sided Striped Detector
+ \item[InitialConditions] Used in NPS to output the initial physics condition of the generated event (Beam interaction position and energy, angle and energy of particule,...)
+ \item[InteractionCoordinates] Used in NPS to output the true physical point of interaction in the differents detector
+ \item[CalibrationManager] This class allowed to manage different file and token and associate them to detector
+ \item[IORoot] Hold the ROOTInput and ROOTOutput classes that allow an easier used of TTRee in NPS and NPA
+ \item[Tools] A few usefull classes such as Nuclear Reaction kinematics, Nucleus object and Energy Loss
+ \item[VDetector] Virtual Detector class from wich inherrited all Physics class
+ \item[include] All header are copy in this folder during compilation
+ \item[lib] All libraries are copy in this folder during compilation
+ \end{itemize}
+
\section{Step by step from TDUMMYDetectorData class}
\subsection{step 1}
+ First, go to your consol and open the NPLib folder, there type \begin{verbatim}make distclean\end{verbatim}
+ This command will remove all the compilation generated file and only the source will remain.
+ Copy the DummyDetector folder and rename it, let's say, WonderdullDetector.
+ Ebter the directory and rename the TDummyDetectorData.h and .cxx with your detector name, you should have now only three file: Makefile, TWonderfullDetectorData.h and TWonderfullDetectorData.cxx
+ If any other file remain, just remove it.
+ No open the .h and .cxx file with your favorite editor (gedit or vim for instance). Replace the the DUMMYDetector by your detector name. Just have look in the .h file, the class is very basic, a couple of value are stored in stl vector (allowing smaller tree and increasing speed access). A few methods allowed to set or get those value and to know the size of the different vector.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\chapter[NPSimulation]{NPSimulation}
diff --git a/NPLib/MUST2/TMust2Physics.cxx b/NPLib/MUST2/TMust2Physics.cxx
index 70b6fcf02de0c1652595e3302fd5f1e19a151a8c..4cd918223b2ae842225aa0c02b6217b5b48634cd 100644
--- a/NPLib/MUST2/TMust2Physics.cxx
+++ b/NPLib/MUST2/TMust2Physics.cxx
@@ -57,15 +57,23 @@ void TMust2Physics::BuildSimplePhysicalEvent()
void TMust2Physics::BuildPhysicalEvent()
{
+ bool check_SILI = false ;
+ bool check_CSI = false ;
+
+
if( CheckEvent() == 1 )
{
vector< TVector2 > couple = Match_X_Y() ;
for(unsigned int i = 0 ; i < couple.size() ; i++)
{
+ check_SILI = false ;
+ check_CSI = false ;
+
int N = EventData->GetMMStripXEDetectorNbr(couple[i].X()) ;
- int X = EventData->GetMMStripXEStripNbr(couple[i].X()) ;
- int Y = EventData->GetMMStripXEStripNbr(couple[i].Y()) ;
+
+ int X = EventData->GetMMStripXEStripNbr(couple[i].X()) ;
+ int Y = EventData->GetMMStripYEStripNbr(couple[i].Y()) ;
double Si_X_E = fSi_X_E(EventData , couple[i].X()) ;
double Si_Y_E = fSi_Y_E(EventData , couple[i].Y()) ;
@@ -95,6 +103,8 @@ void TMust2Physics::BuildPhysicalEvent()
SiLi_N.push_back(EventData->GetMMSiLiEPadNbr(j)) ;
SiLi_E.push_back(fSiLi_E(EventData , j)) ;
SiLi_T.push_back(fSiLi_T(EventData , j)) ;
+ check_SILI = true ;
+
}
}
}
@@ -112,10 +122,27 @@ void TMust2Physics::BuildPhysicalEvent()
CsI_N.push_back(EventData->GetMMCsIECristalNbr(j)) ;
CsI_E.push_back(fCsI_E(EventData , j)) ;
CsI_T.push_back(fCsI_T(EventData , j)) ;
+ check_CSI = true ;
}
}
}
}
+
+
+ if(!check_SILI)
+ {
+ SiLi_N.push_back(0) ;
+ SiLi_E.push_back(0) ;
+ SiLi_T.push_back(0) ;
+ }
+
+ if(!check_CSI)
+ {
+ CsI_N.push_back(0) ;
+ CsI_E.push_back(0) ;
+ CsI_T.push_back(0) ;
+ }
+
}
}
@@ -872,8 +899,8 @@ void TMust2Physics::AddTelescope( double theta ,
TVector3 TMust2Physics::GetPositionOfInteraction(int i)
{
TVector3 Position = TVector3 ( GetStripPositionX( TelescopeNumber[i] , Si_X[i] , Si_Y[i] ) ,
- GetStripPositionY( TelescopeNumber[i] , Si_X[i] , Si_Y[i] ) ,
- GetStripPositionZ( TelescopeNumber[i] , Si_X[i] , Si_Y[i] ) ) ;
+ GetStripPositionY( TelescopeNumber[i] , Si_X[i] , Si_Y[i] ) ,
+ GetStripPositionZ( TelescopeNumber[i] , Si_X[i] , Si_Y[i] ) ) ;
return(Position) ;
@@ -881,21 +908,21 @@ TVector3 TMust2Physics::GetPositionOfInteraction(int i)
TVector3 TMust2Physics::GetTelescopeNormal( int i)
{
- TVector3 U = TVector3 ( GetStripPositionX( TelescopeNumber[i] , 128 , 1 ) ,
- GetStripPositionY( TelescopeNumber[i] , 128 , 1 ) ,
- GetStripPositionZ( TelescopeNumber[i] , 128 , 1 ) )
+ TVector3 U = TVector3 ( GetStripPositionX( TelescopeNumber[i] , 128 , 1 ) ,
+ GetStripPositionY( TelescopeNumber[i] , 128 , 1 ) ,
+ GetStripPositionZ( TelescopeNumber[i] , 128 , 1 ) )
- - TVector3 ( GetStripPositionX( TelescopeNumber[i] , 1 , 1 ) ,
- GetStripPositionY( TelescopeNumber[i] , 1 , 1 ) ,
- GetStripPositionZ( TelescopeNumber[i] , 1 , 1 ) );
+ - TVector3 ( GetStripPositionX( TelescopeNumber[i] , 1 , 1 ) ,
+ GetStripPositionY( TelescopeNumber[i] , 1 , 1 ) ,
+ GetStripPositionZ( TelescopeNumber[i] , 1 , 1 ) );
- TVector3 V = TVector3 ( GetStripPositionX( TelescopeNumber[i] , 128 , 128 ) ,
- GetStripPositionY( TelescopeNumber[i] , 128 , 128 ) ,
- GetStripPositionZ( TelescopeNumber[i] , 128 , 128 ) )
+ TVector3 V = TVector3 ( GetStripPositionX( TelescopeNumber[i] , 128 , 128 ) ,
+ GetStripPositionY( TelescopeNumber[i] , 128 , 128 ) ,
+ GetStripPositionZ( TelescopeNumber[i] , 128 , 128 ) )
- - TVector3 ( GetStripPositionX( TelescopeNumber[i] , 128 , 1 ) ,
- GetStripPositionY( TelescopeNumber[i] , 128 , 1 ) ,
- GetStripPositionZ( TelescopeNumber[i] , 128 , 1 ) );
+ - TVector3 ( GetStripPositionX( TelescopeNumber[i] , 128 , 1 ) ,
+ GetStripPositionY( TelescopeNumber[i] , 128 , 1 ) ,
+ GetStripPositionZ( TelescopeNumber[i] , 128 , 1 ) );
TVector3 Normal = U.Cross(V);
diff --git a/NPSimulation/include/EventGeneratorBeam.hh b/NPSimulation/include/EventGeneratorBeam.hh
index 04e2cd7420f7a51e69ead73f7807f676ad139153..915aa13cc2bca265e69c9675f3c338685c747228 100644
--- a/NPSimulation/include/EventGeneratorBeam.hh
+++ b/NPSimulation/include/EventGeneratorBeam.hh
@@ -50,7 +50,7 @@ public: // Inherit from VEventGenerator Class
void GenerateEvent(G4Event*, G4ParticleGun*);
void InitializeRootOutput();
- void SetTarget(Target* Target) {m_Target = Target;}
+ void SetTarget(Target* Target) {if(Target!=0)m_Target = Target;}
private: // TTree to store initial value of beam and reaction
diff --git a/NPSimulation/include/EventGeneratorTransfert.hh b/NPSimulation/include/EventGeneratorTransfert.hh
index dde7b3f2d34e6e7a8094b13bdb4f6f70d20eb480..cfe4b7be666978167e85d39c8e600733867bdf28 100644
--- a/NPSimulation/include/EventGeneratorTransfert.hh
+++ b/NPSimulation/include/EventGeneratorTransfert.hh
@@ -69,7 +69,7 @@ public: // Constructors and Destructors
public: // Inherit from VEventGenerator class
void ReadConfiguration(string);
void GenerateEvent(G4Event*, G4ParticleGun*);
- void SetTarget(Target* Target) {m_Target = Target;}
+ void SetTarget(Target* Target) {if(Target!=0)m_Target = Target;};
private: // Particle Shoot Option
diff --git a/NPSimulation/include/EventGeneratorTransfertToResonance.hh b/NPSimulation/include/EventGeneratorTransfertToResonance.hh
index 63492f180ae23b47924d13040dcaf7232baa1cb7..eba9c9ef8bf74cf360523791310771d9f0c7f76b 100644
--- a/NPSimulation/include/EventGeneratorTransfertToResonance.hh
+++ b/NPSimulation/include/EventGeneratorTransfertToResonance.hh
@@ -74,7 +74,7 @@ class EventGeneratorTransfertToResonance : public VEventGenerator
public: // Inherit from VEventGenerator class
void ReadConfiguration(string) ;
void GenerateEvent(G4Event*, G4ParticleGun*) ;
- void SetTarget(Target* Target) {m_Target = Target;}
+ void SetTarget(Target* Target) {if(Target!=0)m_Target = Target;};
private: // Particle Shoot Option
bool m_ShootLight ;
diff --git a/NPSimulation/include/Target.hh b/NPSimulation/include/Target.hh
index 3b24daf07f67a18046ac3a7e5128282d22159996..170f7394cb4769e9fc416af1668266d5d88841ba 100644
--- a/NPSimulation/include/Target.hh
+++ b/NPSimulation/include/Target.hh
@@ -43,87 +43,89 @@
using namespace std;
-
-
class Target : public VDetector
{
-public:
- Target();
-
-
-public:
- // Read stream at Configfile to pick-up parameters of detector (Position,...)
- // Called in DetecorConstruction::ReadDetextorConfiguration Method
- void ReadConfiguration(string Path);
-
- // Construct detector and inialise sensitive part.
- // Called After DetecorConstruction::AddDetector Method
- void ConstructDetector(G4LogicalVolume* world);
-
- // Add Detector branch to the EventTree.
- // Called After DetecorConstruction::AddDetector Method
- void InitializeRootOutput();
-
- // Read sensitive part and fill the Root tree.
- // Called at in the EventAction::EndOfEventAvtion
- void ReadSensitive(const G4Event* event);
-
-public:
- // method for debug purpose (still to be implemented)
- // This method should check if the results of the beam interaction within the target
- // (interaction coordinates) are well located inside the target volume
- bool IsInsideTarget() {return false;};
-
- // Used to calculate the incident beam direction (taking into account
- // the emittance) and the vertex of interaction in target
- // Also compute the energy lost by the beam in the target before interacting
- void CalculateBeamInteraction(double MeanPosX, double SigmaPosX, double MeanPosTheta, double SigmaPosTheta,
- double MeanPosY, double SigmaPosY, double MeanPosPhi, double SigmaPosPhi,
- double IncidentBeamEnergy,
- G4ParticleDefinition* BeamName,
- G4ThreeVector &InterCoord, double &AngleEmittanceTheta, double &AngleEmittancePhi,
- double &AngleIncidentTheta, double &AngleIncidentPhi,
- double &FinalBeamEnergy);
-
- // Used to simulate beam emmitance effect
- void RandomGaussian2D(double MeanX, double MeanY, double SigmaX, double SigmaY, double &X, double &Y, double NumberOfSigma = 10000);
-
-public:
- G4Material* GetMaterialFromLibrary(G4String MaterialName, G4double Temperature = 0, G4double Pressure = 0);
-
-
-public:
- G4double GetTargetThickness() {return m_TargetThickness;}
- G4Material* GetTargetMaterial() {return m_TargetMaterial;}
- G4double GetTargetRadius() {return m_TargetRadius;}
- G4double GetTargetAngle() {return m_TargetAngle;}
- G4double GetTargetX() {return m_TargetX;}
- G4double GetTargetY() {return m_TargetY;}
- G4double GetTargetZ() {return m_TargetZ;}
- G4int GetTargetNbLayers() {return m_TargetNbLayers;}
-
-
-private:
- // Target type : true = normal ; false = cryo
- bool m_TargetType;
-
- // Standard parameter
- G4double m_TargetThickness;
- G4double m_TargetRadius;
- G4double m_TargetAngle;
- G4Material* m_TargetMaterial;
- G4int m_TargetNbLayers;
-
- // For Cryo Target
- G4double m_TargetTemperature;
- G4double m_TargetPressure;
- G4double m_WindowsThickness;
- G4Material* m_WindowsMaterial;
-
- // Positioning
- G4double m_TargetX;
- G4double m_TargetY;
- G4double m_TargetZ;
+ public:
+ Target();
+ ~Target(){};
+
+
+ public:
+ // Read stream at Configfile to pick-up parameters of detector (Position,...)
+ // Called in DetecorConstruction::ReadDetextorConfiguration Method
+ void ReadConfiguration(string Path);
+
+ // Construct detector and inialise sensitive part.
+ // Called After DetecorConstruction::AddDetector Method
+ void ConstructDetector(G4LogicalVolume* world);
+
+ // Add Detector branch to the EventTree.
+ // Called After DetecorConstruction::AddDetector Method
+ void InitializeRootOutput();
+
+ // Read sensitive part and fill the Root tree.
+ // Called at in the EventAction::EndOfEventAvtion
+ void ReadSensitive(const G4Event* event);
+
+ public:
+ // method for debug purpose (still to be implemented)
+ // This method should check if the results of the beam interaction within the target
+ // (interaction coordinates) are well located inside the target volume
+ bool IsInsideTarget() {return false;};
+
+ // Used to calculate the incident beam direction (taking into account
+ // the emittance) and the vertex of interaction in target
+ // Also compute the energy lost by the beam in the target before interacting
+ void CalculateBeamInteraction(double MeanPosX, double SigmaPosX, double MeanPosTheta, double SigmaPosTheta,
+ double MeanPosY, double SigmaPosY, double MeanPosPhi, double SigmaPosPhi,
+ double IncidentBeamEnergy,
+ G4ParticleDefinition* BeamName,
+ G4ThreeVector &InterCoord, double &AngleEmittanceTheta, double &AngleEmittancePhi,
+ double &AngleIncidentTheta, double &AngleIncidentPhi,
+ double &FinalBeamEnergy);
+
+ // Used to simulate beam emmitance effect
+ void RandomGaussian2D(double MeanX, double MeanY, double SigmaX, double SigmaY, double &X, double &Y, double NumberOfSigma = 10000);
+
+ public:
+ // Return Material from the Target Material Library
+ G4Material* GetMaterialFromLibrary(G4String MaterialName, G4double Temperature = 0, G4double Pressure = 0);
+
+ // Generate a DEDX file table using the material used in the target
+ void WriteDEDXTable(G4ParticleDefinition* Particle,G4double Emin,G4double Emax);
+
+ public:
+ G4double GetTargetThickness() {return m_TargetThickness;}
+ G4Material* GetTargetMaterial() {return m_TargetMaterial;}
+ G4double GetTargetRadius() {return m_TargetRadius;}
+ G4double GetTargetAngle() {return m_TargetAngle;}
+ G4double GetTargetX() {return m_TargetX;}
+ G4double GetTargetY() {return m_TargetY;}
+ G4double GetTargetZ() {return m_TargetZ;}
+ G4int GetTargetNbLayers() {return m_TargetNbLayers;}
+
+
+ private:
+ // Target type : true = normal ; false = cryo
+ bool m_TargetType;
+
+ // Standard parameter
+ G4double m_TargetThickness;
+ G4double m_TargetRadius;
+ G4double m_TargetAngle;
+ G4Material* m_TargetMaterial;
+ G4int m_TargetNbLayers;
+
+ // For Cryo Target
+ G4double m_TargetTemperature;
+ G4double m_TargetPressure;
+ G4double m_WindowsThickness;
+ G4Material* m_WindowsMaterial;
+
+ // Positioning
+ G4double m_TargetX;
+ G4double m_TargetY;
+ G4double m_TargetZ;
};
#endif
diff --git a/NPSimulation/src/EventGeneratorTransfertToResonance.cc b/NPSimulation/src/EventGeneratorTransfertToResonance.cc
index c12b8b13f146f4a1ad652d507c082046b76b565b..0990ad3abcc8a33b57bfbec613b2b0d6075652fb 100644
--- a/NPSimulation/src/EventGeneratorTransfertToResonance.cc
+++ b/NPSimulation/src/EventGeneratorTransfertToResonance.cc
@@ -57,7 +57,7 @@ EventGeneratorTransfertToResonance::EventGeneratorTransfertToResonance()
//------------- Default Constructor -------------
m_InitConditions = new TInitialConditions() ;
m_Reaction = new Reaction() ;
- m_Target = 0;
+ m_Target = new Target();
m_SigmaX = 0 ;
m_SigmaY = 0 ;
m_SigmaThetaX = 0 ;
diff --git a/NPSimulation/src/Target.cc b/NPSimulation/src/Target.cc
index af9f4998b6a4a1c12f0619d30232b73439311d7c..2d9e53ed4e5360ee8d5f4c9379338019a2adde9c 100644
--- a/NPSimulation/src/Target.cc
+++ b/NPSimulation/src/Target.cc
@@ -28,6 +28,7 @@
// C++ header
#include <fstream>
#include <limits>
+
// G4 geometry header
#include "G4Tubs.hh"
@@ -42,6 +43,8 @@
#include "G4VisAttributes.hh"
#include "G4Colour.hh"
#include "G4EmCalculator.hh"
+#include "G4ParticleDefinition.hh"
+#include "G4ParticleTable.hh"
#include "Randomize.hh"
using namespace CLHEP ;
// NPTool header
@@ -480,6 +483,9 @@ void Target::ConstructDetector(G4LogicalVolume* world)
}
}
+
+
+WriteDEDXTable(G4ParticleTable::GetParticleTable()->GetIon(2, 3, 0.) ,0,300);
}
// Add Detector branch to the EventTree.
@@ -553,20 +559,28 @@ void Target::CalculateBeamInteraction( double MeanPosX, double SigmaPosX, double
y0 += m_TargetY;
z0 += m_TargetZ;
InterCoord = G4ThreeVector(x0, y0, z0);
-
- G4EmCalculator emCalculator;
+
if(m_TargetType)
{
- for (G4int i = 0; i < m_TargetNbLayers; i++)
+ G4EmCalculator emCalculator;
+ if(m_TargetThickness!=0)
{
- G4double dedx = emCalculator.ComputeTotalDEDX(IncidentBeamEnergy, BeamName, m_TargetMaterial);
- G4double de = dedx * EffectiveTargetThicknessBeforeInteraction / m_TargetNbLayers;
- IncidentBeamEnergy -= de;
+ for (G4int i = 0; i < m_TargetNbLayers; i++)
+ {
+ G4double dedx = emCalculator.ComputeTotalDEDX(IncidentBeamEnergy, BeamName, m_TargetMaterial);
+ G4double de = dedx * EffectiveTargetThicknessBeforeInteraction / m_TargetNbLayers;
+ IncidentBeamEnergy -= de;
+ }
+
}
}
+
+
else
- { // Windows
+ { G4EmCalculator emCalculator;
+ // Windows
+ if(m_WindowsThickness!=0)
for (G4int i = 0; i < m_TargetNbLayers; i++)
{
G4double dedx = emCalculator.ComputeTotalDEDX(IncidentBeamEnergy, BeamName, m_WindowsMaterial);
@@ -575,6 +589,7 @@ void Target::CalculateBeamInteraction( double MeanPosX, double SigmaPosX, double
}
// Target
+ if(m_TargetThickness!=0)
for (G4int i = 0; i < m_TargetNbLayers; i++)
{
G4double dedx = emCalculator.ComputeTotalDEDX(IncidentBeamEnergy, BeamName, m_TargetMaterial);
@@ -607,3 +622,46 @@ void Target::RandomGaussian2D(double MeanX, double MeanY, double SigmaX, double
Y = RandGauss::shoot(MeanY, SigmaY);
}
}
+
+// Generate a DEDX file table using the material used in the target
+void Target::WriteDEDXTable(G4ParticleDefinition* Particle ,G4double Emin,G4double Emax)
+ {
+ // Opening hte output file
+ G4String GlobalPath = getenv("NPTOOL");
+ G4String Path = GlobalPath + "/Inputs/EnergyLoss/" + Particle->GetParticleName() + "_" + m_TargetMaterial->GetName() + ".G4table";
+
+ ofstream File ;
+ File.open(Path) ;
+
+ if(!File) return ;
+
+ File << "Table from Geant4 generate using NPSimulation" << endl
+ << "Particle: " << Particle->GetParticleName() << "\tMaterial: " << m_TargetMaterial->GetName() << endl ;
+
+ G4EmCalculator emCalculator;
+
+ for (G4double E=Emin*MeV; E < Emax*MeV; E+=(Emax-Emin)*MeV/10000.)
+ {
+ G4double dedx = emCalculator.ComputeTotalDEDX(E, Particle, m_TargetMaterial);
+ File << E/MeV << "\t" << dedx/(MeV/um) << endl ;
+ }
+ File.close();
+
+ if(!m_TargetType)
+ {
+ G4String Path = GlobalPath + "/Inputs/EnergyLoss/" + Particle->GetParticleName() + "_" + m_WindowsMaterial->GetName() + ".G4table";
+ File.open(Path) ;
+ if(!File) return ;
+ File << "Table from Geant4 generate using NPSimulation" << endl
+ << "Particle: " << Particle->GetParticleName() << "\tMaterial: " << m_TargetMaterial->GetName() << endl ;
+
+ for (G4double E=Emin*MeV; E < Emax*MeV; E+=(Emax-Emin)*MeV/10000.)
+ {
+ G4double dedx = emCalculator.ComputeTotalDEDX(E, Particle, m_WindowsMaterial);
+ File << E/MeV << "\t" << dedx/(MeV/um) << endl ;
+ }
+ }
+
+
+ File.close();
+ }