diff --git a/INSTALL b/INSTALL index 35f9cbaf7fa3ac1fe9affbd6c6f43e1e69b91a2e..2235f0b7d40d65a52d7522f02be5d73ed5a06aa4 100644 --- a/INSTALL +++ b/INSTALL @@ -29,7 +29,7 @@ I) REQUIREMENTS II) WORKING CONFIGURATIONS The NPTool package has been mainly tested (and developped) with the two following configurations: - 1) Mac OS X (10.5.7) + G4 (4.9.2) + CLHEP (2.0.4.2) + ROOT (5.22/00 and 5.24/00) + GSL () + 1) Mac OS X (10.5.7) + G4 (4.9.2 and 4.9.3.b1) + CLHEP (2.0.4.2) + ROOT (5.22/00 and 5.24/00) + GSL () 2) Linux Fedora (kernel 2.6.29) + G4 (4.9.1p3) + CLHEP (2.0.4.2) + ROOT (5.22/00) + GSL (1.12) Please, report any working or non-working configuration. diff --git a/NPAnalysis/Gaspard/RunToTreat.txt b/NPAnalysis/Gaspard/RunToTreat.txt index e1e642cec7edf653af2ca43f8dde76c1ccd3c306..5e6a3389381c7b3d88643f0ffc3311d201f72d07 100644 --- a/NPAnalysis/Gaspard/RunToTreat.txt +++ b/NPAnalysis/Gaspard/RunToTreat.txt @@ -1,5 +1,5 @@ TTreeName SimulatedTree RootFileName -% ../../Outputs/Simulation/mySimul.root - ../../Outputs/Simulation/sn132dp_10MeVA_T1_B0_E0.root + ../../Outputs/Simulation/mySimul.root +% ../../Outputs/Simulation/sn132dp_10MeVA_T1_B0_E0.root diff --git a/NPDocumentation/Gaspard.tex b/NPDocumentation/Gaspard.tex index 5f0ccce02498d0d86f4d34fb6b3e633ce8596262..e486f4bc743e576454702cc44f370f8a31e6e8ae 100755 --- a/NPDocumentation/Gaspard.tex +++ b/NPDocumentation/Gaspard.tex @@ -169,10 +169,10 @@ contains three classes: \subsection{Adding a new detector shape to Gaspard} A special class (GaspardTrackerDummyShape) has been created to show how to define a new module for the Gaspard Tracker. This class describes a -simple 5x5 cm2 square telescope made of three layers of silicon which +simple 5x5 cm$^2$ square telescope made of three layers of silicon which has been used for some preliminary studies of the tracker. So, when considering adding a new module to the Gaspard Tracker, please do not use -this class but create you own instead. However, for the explanations the +this class but create your own instead. However, for the explanations the GaspardTrackerDummyShape case will be considered. When adding a new detector you need to follow several steps: diff --git a/NPSimulation/include/Target.hh b/NPSimulation/include/Target.hh index 5fcee069b67b69476d59c52eb7ca5af51fde514e..92da2be70c450cd297b0229bee4eb68050db72f3 100644 --- a/NPSimulation/include/Target.hh +++ b/NPSimulation/include/Target.hh @@ -79,12 +79,13 @@ public: public: - G4double GetTargetThickness() {return m_TargetThickness;} - 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;} + 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;} private: diff --git a/NPSimulation/include/VEventGenerator.hh b/NPSimulation/include/VEventGenerator.hh index 1364870777c13833f1c8cb741f6a99207d3fc135..f5d040a10d404b7acedbf9fd0c8fb2fdcf7b2da4 100644 --- a/NPSimulation/include/VEventGenerator.hh +++ b/NPSimulation/include/VEventGenerator.hh @@ -11,7 +11,7 @@ * Original Author: Adrien MATTA contact address: matta@ipno.in2p3.fr * * * * Creation Date : January 2009 * - * Last update : * + * Last update : 03/11/2009 * *---------------------------------------------------------------------------* * Decription: * * All event generator added in the project should derive from this virtual * @@ -23,6 +23,8 @@ * * *---------------------------------------------------------------------------* * Comment: * + * + 03/11/09: Adding EffectiveThiknessBeforeInteraction in the * + * CalculateBeamInteraction() method (N. de Sereville) * * * *****************************************************************************/ // C++ header @@ -64,7 +66,8 @@ public: double MeanPosY, double SigmaPosY, double MeanPosPhi, double SigmaPosPhi, Target* target, G4ThreeVector &InterCoord, double &AngleEmittanceTheta, double &AngleEmittancePhi, - double &AngleIncidentTheta, double &AngleIncidentPhi); + double &AngleIncidentTheta, double &AngleIncidentPhi, + double &EffectiveTargetThicknessBeforeInteraction); // Used to simulate beam emmitance effect void RandomGaussian2D(double MeanX, double MeanY, double SigmaX, double SigmaY, double &X, double &Y, double NumberOfSigma = 10000); diff --git a/NPSimulation/src/EventGeneratorBeam.cc b/NPSimulation/src/EventGeneratorBeam.cc index 664ccf1d592f38eeb07ecfcb7233131ea65a443b..c0418d7222318e0ebe0ad748b7bca500a01e806d 100644 --- a/NPSimulation/src/EventGeneratorBeam.cc +++ b/NPSimulation/src/EventGeneratorBeam.cc @@ -180,11 +180,13 @@ void EventGeneratorBeam::GenerateEvent(G4Event* anEvent, G4ParticleGun* particle G4ThreeVector InterCoord; G4double Beam_thetaX = 0, Beam_phiY = 0; G4double Beam_theta = 0, Beam_phi = 0; + G4double EffectiveThickness = 0; CalculateBeamInteraction(0, m_SigmaX, 0, m_SigmaThetaX, 0, m_SigmaY, 0, m_SigmaPhiY, m_Target, InterCoord, Beam_thetaX, Beam_phiY, - Beam_theta, Beam_phi); + Beam_theta, Beam_phi, + EffectiveThickness); // write vertex position to ROOT file G4double x0 = InterCoord.x(); diff --git a/NPSimulation/src/EventGeneratorTransfert.cc b/NPSimulation/src/EventGeneratorTransfert.cc index 87bba8e1144d8b126c41c8e9071a2feb4d795a18..2f4d422519761efdd20a1b0516112803aa16dddb 100644 --- a/NPSimulation/src/EventGeneratorTransfert.cc +++ b/NPSimulation/src/EventGeneratorTransfert.cc @@ -30,6 +30,7 @@ // G4 headers #include "G4ParticleTable.hh" +#include "G4EmCalculator.hh" #include "G4ParticleGun.hh" #include "G4RotationMatrix.hh" @@ -320,18 +321,25 @@ void EventGeneratorTransfert::GenerateEvent(G4Event* anEvent , G4ParticleGun* pa ////////////////////////////////////////////////// //////Define the kind of particle to shoot//////// ////////////////////////////////////////////////// + // Light G4int LightZ = m_Reaction->GetNucleus3()->GetZ() ; G4int LightA = m_Reaction->GetNucleus3()->GetA() ; G4ParticleDefinition* LightName = G4ParticleTable::GetParticleTable()->GetIon(LightZ, LightA, 0.); + // Recoil G4int HeavyZ = m_Reaction->GetNucleus4()->GetZ() ; G4int HeavyA = m_Reaction->GetNucleus4()->GetA() ; G4ParticleDefinition* HeavyName = G4ParticleTable::GetParticleTable()->GetIon(HeavyZ, HeavyA, m_Reaction->GetExcitation()*MeV); + // Beam + G4int BeamZ = m_Reaction->GetNucleus1()->GetZ(); + G4int BeamA = m_Reaction->GetNucleus1()->GetA(); + G4ParticleDefinition* BeamName = G4ParticleTable::GetParticleTable()->GetIon(BeamZ, BeamA, 0); + /////////////////////////////////////////////////////////////////////// ///// Calculate the incident beam direction as well as the vertex ///// ///// of interaction in target ///// @@ -339,11 +347,13 @@ void EventGeneratorTransfert::GenerateEvent(G4Event* anEvent , G4ParticleGun* pa G4ThreeVector InterCoord; G4double Beam_thetaX = 0, Beam_phiY = 0; G4double Beam_theta = 0, Beam_phi = 0; + G4double EffectiveThickness = 0; CalculateBeamInteraction(0, m_SigmaX, 0, m_SigmaThetaX, 0, m_SigmaY, 0, m_SigmaPhiY, m_Target, InterCoord, Beam_thetaX, Beam_phiY, - Beam_theta, Beam_phi); + Beam_theta, Beam_phi, + EffectiveThickness); // write vertex position to ROOT file G4double x0 = InterCoord.x(); @@ -380,9 +390,19 @@ void EventGeneratorTransfert::GenerateEvent(G4Event* anEvent , G4ParticleGun* pa ///// Angles for emitted particles following Cross Section ////// ///// Angles are in the beam frame ////// ///////////////////////////////////////////////////////////////// - // Beam incident energy + // Beam nominal incident energy before interaction with the target G4double NominalBeamEnergy = m_BeamEnergy; G4double IncidentBeamEnergy = RandGauss::shoot(NominalBeamEnergy, m_BeamEnergySpread / 2.35); + // Slowing down the beam to the interaction layer in the target + // Number of Layers + const G4int NbLayers = 50; + G4EmCalculator emCalculator; + for (G4int i = 0; i < NbLayers; i++) { +// G4double dedx = emCalculator.GetDEDX(IncidentBeamEnergy, BeamName, m_Target->GetTargetMaterial()); + G4double dedx = emCalculator.ComputeTotalDEDX(IncidentBeamEnergy, BeamName, m_Target->GetTargetMaterial()); + G4double de = dedx * EffectiveThickness / NbLayers; + IncidentBeamEnergy -= de; + } m_Reaction->SetBeamEnergy(IncidentBeamEnergy); m_InitConditions->SetICIncidentEnergy(IncidentBeamEnergy / MeV); // Angles diff --git a/NPSimulation/src/EventGeneratorTransfertToResonance.cc b/NPSimulation/src/EventGeneratorTransfertToResonance.cc index eeb2e8afc8d16a893e2314af056ad8e1e19943dc..790eb932209a63fd01b80b775f41bb0b2c633931 100644 --- a/NPSimulation/src/EventGeneratorTransfertToResonance.cc +++ b/NPSimulation/src/EventGeneratorTransfertToResonance.cc @@ -397,6 +397,7 @@ void EventGeneratorTransfertToResonance::GenerateEvent(G4Event* anEvent , G4Part G4ThreeVector InterCoord; G4double Beam_thetaX = 0, Beam_phiY = 0; G4double Beam_theta = 0, Beam_phi = 0; + G4double EffectiveThickness = 0; if(m_SigmaX==0 && m_SigmaY==0 && m_SigmaThetaX == 0 && m_SigmaPhiY==0) InterCoord = G4ThreeVector(0,0,0); @@ -405,7 +406,8 @@ void EventGeneratorTransfertToResonance::GenerateEvent(G4Event* anEvent , G4Part 0, m_SigmaY, 0, m_SigmaPhiY, m_Target, InterCoord, Beam_thetaX, Beam_phiY, - Beam_theta, Beam_phi); + Beam_theta, Beam_phi, + EffectiveThickness); else InterCoord = G4ThreeVector(0,0,0); diff --git a/NPSimulation/src/VEventGenerator.cc b/NPSimulation/src/VEventGenerator.cc index 6ab94ebb3f4373f950ed733414510ea43825fa56..4dc2aac3db8508daa908d46eeeafbca5bf0a5c86 100644 --- a/NPSimulation/src/VEventGenerator.cc +++ b/NPSimulation/src/VEventGenerator.cc @@ -9,7 +9,7 @@ * Original Author: Adrien MATTA contact address: matta@ipno.in2p3.fr * * * * Creation Date : January 2009 * - * Last update : * + * Last update : 03/11/2009 * *---------------------------------------------------------------------------* * Decription: * * All event generator added in the project should derive from this virtual * @@ -21,6 +21,8 @@ * * *---------------------------------------------------------------------------* * Comment: * + * + 03/11/09: Adding EffectiveThiknessBeforeInteraction in the * + * CalculateBeamInteraction() method (N. de Sereville) * * * *****************************************************************************/ #include "VEventGenerator.hh" @@ -28,8 +30,8 @@ // C++ headers #include "cmath" -// ROOT headers -//#include"TMath.h" +// G4 headers +#include "G4UnitsTable.hh" // CLHEP headers #include "Randomize.hh" @@ -51,7 +53,8 @@ void VEventGenerator::CalculateBeamInteraction(double MeanPosX, double SigmaPosX double MeanPosY, double SigmaPosY, double MeanPosPhi, double SigmaPosPhi, Target* target, G4ThreeVector &InterCoord, double &AngleEmittanceTheta, double &AngleEmittancePhi, - double &AngleIncidentTheta, double &AngleIncidentPhi) + double &AngleIncidentTheta, double &AngleIncidentPhi, + double &EffectiveTargetThicknessBeforeInteraction) { // target parameters double TargetThickness = target->GetTargetThickness(); @@ -71,6 +74,7 @@ void VEventGenerator::CalculateBeamInteraction(double MeanPosX, double SigmaPosX RandomGaussian2D(MeanPosY, MeanPosPhi, SigmaPosY, SigmaPosPhi, y0, AngleEmittancePhi); } // in case target is tilted, correct the z-position of interaction + // x is the vertical axis dz = x0 * tan(TargetAngle); } else { // if no target radius is given consider a point-like target @@ -78,18 +82,6 @@ void VEventGenerator::CalculateBeamInteraction(double MeanPosX, double SigmaPosX RandomGaussian2D(0, 0, 0, SigmaPosPhi, y0, AngleEmittancePhi); } - // correct for the target angle wrt the beam axis - // this simple correction is only valid if the beam is parallel to the beam axis - // should be improved in a next version - TargetThickness /= cos(TargetAngle); - z0 = dz + (-TargetThickness / 2 + RandFlat::shoot() * TargetThickness); - - // Move to the target position - x0 += target->GetTargetX(); - y0 += target->GetTargetY(); - z0 += target->GetTargetZ(); - InterCoord = G4ThreeVector(x0, y0, z0); - // Calculate incident angle in spherical coordinate, passing by the direction vector dir double Xdir = sin(AngleEmittanceTheta); double Ydir = sin(AngleEmittancePhi); @@ -99,6 +91,29 @@ void VEventGenerator::CalculateBeamInteraction(double MeanPosX, double SigmaPosX AngleIncidentPhi = atan2(Ydir, Xdir) * rad; if (AngleIncidentPhi < 0) AngleIncidentPhi += 2*pi; if (AngleIncidentTheta < 1e-6) AngleIncidentPhi = 0; + + // Calculation of effective target thickness and z-position of interaction + // when the target is tilted wrt the beam axis + // * exact if target is perpendicular to the beam axis (target not tilted) + // in any case of beam emittance + // * exact if target is tilted wrt the beam axis and the beam is parallel + // (no beam divergence) + // * wrong if target is tilted and for general case of beam emittance + // (should be fixed in next release). For small beam divergence this + // should not be such a big problem + TargetThickness /= cos(TargetAngle); + double uniform = RandFlat::shoot(); + z0 = dz + (-TargetThickness / 2 + uniform * TargetThickness); + + // Calculate the effective thickness before interaction in target + // This is useful to slow down the beam + EffectiveTargetThicknessBeforeInteraction = TargetThickness * uniform / cos(AngleIncidentTheta); + + // Move to the target position + x0 += target->GetTargetX(); + y0 += target->GetTargetY(); + z0 += target->GetTargetZ(); + InterCoord = G4ThreeVector(x0, y0, z0); } diff --git a/TODO b/TODO index 5bb0b745cec1dd08f7bfde55842c8dfb82cd118f..a7a00e5e14c17d7a742c0cec06020773c39b4bb5 100644 --- a/TODO +++ b/TODO @@ -22,22 +22,21 @@ TODO for NPTool: ---------------- - + Build a dummy detector as a simple example of how NPTool is working (Adrien) + + Retrieve energyloss tables from G4 to be used for slowing the beam + and in the analysis. + Target work (Nicolas) - + Make homogeneous treatment of materials definition wrt other detectors - (Adrien) + Take correctly into account the emittance when the target is tilted. + IsInsideTarget method for debugging purposes + Add a dedicated class to deal with materials (see example from G4 tutorial) (still under debate) + Split physic list and give the possibility to choose which package to use + Add support for messengers + + Scorers: CopyNumber v.s. DetectorNumber? (to be checked) + Build NPTool with the Autotool/Automake + Add documentation (Adrien + Nicolas) TODO for NPAnalysis: -------------------- - + Add calibration manager class from Julien (Adrien) + Use only one class for the analysis per detector instead of two currently (ROOT feature in private member definition ( //!))