From 4d4d02028932c77469c2bf4a01bf6f5bb17fb922 Mon Sep 17 00:00:00 2001
From: matta <matta@npt>
Date: Thu, 31 Jan 2013 12:56:51 +0000
Subject: [PATCH] * NPBeam now inherit from NPNucleus

---
 NPLib/Physics/NPBeam.cxx                      |   7 +-
 NPLib/Physics/NPBeam.h                        |   8 +-
 NPLib/Physics/NPNucleus.cxx                   |  59 +++---
 NPLib/Physics/NPNucleus.h                     | 177 +++++++++---------
 NPLib/Physics/NPReaction.cxx                  |  22 +--
 NPLib/Physics/NPReaction.h                    |   9 +-
 NPSimulation/src/EventGeneratorBeam.cc        |   2 +-
 .../src/EventGeneratorParticleDecay.cc        |   2 +-
 8 files changed, 140 insertions(+), 146 deletions(-)

diff --git a/NPLib/Physics/NPBeam.cxx b/NPLib/Physics/NPBeam.cxx
index 23df5d87c..a7850f136 100644
--- a/NPLib/Physics/NPBeam.cxx
+++ b/NPLib/Physics/NPBeam.cxx
@@ -46,7 +46,6 @@ using namespace NPL;
 
 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 
 Beam::Beam(){
-  fBeamNucleus = new Nucleus();
   fEnergy = 0;
   fSigmaEnergy = -1 ;
   fMeanX = 0 ;
@@ -78,7 +77,6 @@ Beam::Beam(){
 
 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
 Beam::~Beam(){
-  delete fBeamNucleus ;
   delete fEnergyHist  ;
   delete fXThetaXHist ;
   delete fYPhiYHist   ;
@@ -137,9 +135,8 @@ void Beam::ReadConfigurationFile(string Path){
       else if (DataBuffer == "Particle=") {
         check_BeamName = true ;
         ReactionFile >> DataBuffer;
-        delete fBeamNucleus;
-        fBeamNucleus = new Nucleus(DataBuffer);
-        if(fVerboseLevel==1) cout << "Beam Particle: " << fBeamNucleus->GetName() << endl;
+        SetUp(DataBuffer);
+        if(fVerboseLevel==1) cout << "Beam Particle: " << GetName() << endl;
       }
       
       else if (DataBuffer == "Energy=") {
diff --git a/NPLib/Physics/NPBeam.h b/NPLib/Physics/NPBeam.h
index db6b991dc..5b3f85fdd 100755
--- a/NPLib/Physics/NPBeam.h
+++ b/NPLib/Physics/NPBeam.h
@@ -41,7 +41,7 @@ using namespace NPL;
 
 namespace NPL{
   
-  class Beam{
+  class Beam:public NPL::Nucleus{
     
   public:  // Constructors and Destructors
     Beam();
@@ -54,7 +54,7 @@ namespace NPL{
     int fVerboseLevel;
   
   private:
-    Nucleus* fBeamNucleus;
+    //Nucleus* fBeamNucleus;
     double fEnergy;
     double fSigmaEnergy;
     double fMeanX;
@@ -74,7 +74,7 @@ namespace NPL{
   public:
     // Getters and Setters
     // Set
-    void SetBeamNucleus (Nucleus* BeamNucleus)  {delete fBeamNucleus ; fBeamNucleus = new Nucleus(BeamNucleus->GetZ(),BeamNucleus->GetA());}
+    // void SetBeamNucleus (Nucleus* BeamNucleus)  {delete fBeamNucleus ; fBeamNucleus = new Nucleus(BeamNucleus->GetZ(),BeamNucleus->GetA());}
     void SetEnergy      (double Energy)         {fEnergy=Energy;}
     void SetSigmaEnergy (double SigmaEnergy)    {fSigmaEnergy=SigmaEnergy;}
     void SetMeanX       (double MeanX)          {fMeanX=MeanX;}
@@ -91,7 +91,7 @@ namespace NPL{
     void SetVerboseLevel(int verbose)           {fVerboseLevel = verbose;}
 
     // Get
-    Nucleus*  GetNucleus     () const {return fBeamNucleus;}
+    // Nucleus*  GetNucleus     () const {return fBeamNucleus;}
     double    GetEnergy      () const {return fEnergy;}
     double    GetSigmaEnergy () const {return fSigmaEnergy;}
     double    GetMeanX       () const {return fMeanX;}
diff --git a/NPLib/Physics/NPNucleus.cxx b/NPLib/Physics/NPNucleus.cxx
index bc4d608ab..abd853aa4 100644
--- a/NPLib/Physics/NPNucleus.cxx
+++ b/NPLib/Physics/NPNucleus.cxx
@@ -56,34 +56,39 @@ Nucleus::Nucleus()
 }
 
 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 
-Nucleus::Nucleus(string isotope)
-{
-   //----------- Constructor Using nubtab03.asc by name----------
-   // open the file to read and check if it is open
-   
-   const char* Isotope = isotope.c_str();
-   
-   ifstream inFile;
-   string Path = getenv("NPTOOL") ;
-   string FileName = Path + "/NPLib/Physics/nubtab03.asc";
-   inFile.open(FileName.c_str());
-   
-   // reading the file
-   string line, s_name;
-   size_t space;
-   if (inFile.is_open()) {
-      while (!inFile.eof()) {
-         getline(inFile,line);
-
-    s_name = line.substr(11,7);
-    space = s_name.find_first_of(" "); 
-    s_name.resize(space);
+Nucleus::Nucleus(string isotope){
+  SetUp(isotope);
+}
 
-    if (s_name.find(Isotope) != string::npos && s_name.length() == isotope.length()) break;
-      }
-      Extract(line.data());
-   }
-   else cout << "Unable to open file nuclear data base file " << FileName << endl;
+void Nucleus::SetUp(string isotope){
+  //----------- Constructor Using nubtab03.asc by name----------
+  // open the file to read and check if it is open
+  
+  const char* Isotope = isotope.c_str();
+  
+  ifstream inFile;
+  string Path = getenv("NPTOOL") ;
+  string FileName = Path + "/NPLib/Physics/nubtab03.asc";
+  inFile.open(FileName.c_str());
+  
+  // reading the file
+  string line, s_name;
+  size_t space;
+  if (inFile.is_open()) {
+    while (!inFile.eof()) {
+      getline(inFile,line);
+      
+      s_name = line.substr(11,7);
+      space = s_name.find_first_of(" ");
+      s_name.resize(space);
+      
+      if (s_name.find(Isotope) != string::npos && s_name.length() == isotope.length()) break;
+    }
+    Extract(line.data());
+  }
+  else cout << "Unable to open file nuclear data base file " << FileName << endl;
+  
+  
 }
 
 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 
diff --git a/NPLib/Physics/NPNucleus.h b/NPLib/Physics/NPNucleus.h
index 650ab6cfb..7cbee59a9 100644
--- a/NPLib/Physics/NPNucleus.h
+++ b/NPLib/Physics/NPNucleus.h
@@ -27,93 +27,92 @@
 #include "TLorentzVector.h"
 using namespace std;
 
-namespace NPL 
-   {
-      class Nucleus {
-
-      public:
-         Nucleus();
-         Nucleus(string isotope);
-         Nucleus(int Z, int A);
-         ~Nucleus();
-
-       private :
-		  //intrinsic properties
-		  const char* fName;         // Nucleus name
-		  string	    fNucleusName;
-		  int         fCharge;      // Nucleus charge
-		  int         fAtomicWeight;   // Nucleus atomic weight
-		  double      fMassExcess;   // Nucleus mass excess in keV
-		  const char* fSpinParity;   // Nucleus spin and parity
-		  double      fSpin;         // Nucleus spin
-		  const char* fParity;      // Nucleus parity
-		  //kinematic properties
-		  double fKineticEnergy;
-		  double fBeta;
-		  double fGamma;
-		  double fBrho;
-		  double fTimeOfFlight;
-		  double fVelocity;
-		  TLorentzVector fEnergyImpulsion;
-		  
-	  public:
-		  void  EnergyToBrho();
-		  void  EnergyToTof();
-		  void	BetaToVelocity();
-		  void  BrhoToEnergy();
-		  void  BrhoToTof()    {BrhoToEnergy(); EnergyToTof();}
-		  void	TofToEnergy();
-		  void	TofToBrho()    {TofToEnergy(); EnergyToBrho();}
-		  void	EnergyToBeta();
-		  void	BetaToEnergy();
-		  void	BetaToGamma();
-		  double DopplerCorrection(double EnergyLabGamma, double ThetaLabGamma);
-		  
-		  
-       protected :
-		  void Extract(const char* line);
-
-       public :
-		  void				GetNucleusName();
-		  string			GetName()			const				{return fNucleusName;}
-		  int				GetZ()				const				{return fCharge;}
-		  int				GetA()				const				{return fAtomicWeight;}
-		  double			GetMassExcess()		const				{return fMassExcess;}
-		  const char*		GetSpinParity()		const				{return fSpinParity;}
-		  double			GetSpin()			const				{return fSpin;}
-		  const char*		GetParity()			const				{return fParity;}
-		  double			GetEnergy()			const				{return fKineticEnergy;}
-		  double			GetBrho()			const				{return fBrho;}
-		  double			GetTimeOfFlight()	const				{return fTimeOfFlight;}
-		  double			GetBeta()			const				{return fBeta;}
-		  double			GetGamma()			const				{return fGamma;}
-		  double			GetVelocity()		const				{return fVelocity;}
-		  TLorentzVector	GetEnergyImpulsion() const				{return fEnergyImpulsion;}
-		  void				SetName(const char* name)				{fName = name;}
-		  void				SetZ(int charge)						{fCharge = charge;}
-		  void				SetA(int mass)							{fAtomicWeight = mass;}
-		  void				SetMassExcess(double massexcess)		{fMassExcess = massexcess;}
-		  void				SetSpinParity(const char* spinparity)	{fSpinParity = spinparity;}
-		  void				SetSpin(double spin)					{fSpin = spin;}
-		  void				SetParity(const char* parity)			{fParity = parity;}
-		  void				SetKineticEnergy(double energy)			{fKineticEnergy = energy; EnergyToBrho(); EnergyToTof(); EnergyToBeta(); BetaToGamma();BetaToVelocity();}
-		  void				SetBrho(double brho)					{fBrho = brho; BrhoToEnergy(); BrhoToTof(); EnergyToBeta(); BetaToGamma();BetaToVelocity();}
-		  void				SetTimeOfFlight(double tof)				{fTimeOfFlight = tof; TofToEnergy(); TofToBrho(); EnergyToBeta(); BetaToGamma();BetaToVelocity();}
-		  void				SetEnergyImpulsion(TLorentzVector P) 	{fEnergyImpulsion = P; 
-																	fKineticEnergy = P.E() - Mass();
-																	EnergyToBrho();
-																	EnergyToTof();
-																	EnergyToBeta();
-																	BetaToGamma();
-																	BetaToVelocity();}	
-		  void				SetBeta(double beta)					{fBeta = beta; BetaToGamma(); BetaToEnergy(); EnergyToTof(); EnergyToBrho();BetaToVelocity();}
-													
-			  
-										
-		  
-		  // Nuclear mass in MeV
-         double      Mass() const {return (fAtomicWeight*uma + fMassExcess/1000.);}
-         void        Print() const   ;
-      };
-   }
+namespace NPL {
+  class Nucleus {
+    
+  public:
+    Nucleus();
+    Nucleus(string isotope);
+    Nucleus(int Z, int A);
+    ~Nucleus();
+  
+  public:
+    void SetUp(string isotope);
+    
+  private :
+    //intrinsic properties
+    const char* fName;         // Nucleus name
+    string	    fNucleusName;
+    int         fCharge;      // Nucleus charge
+    int         fAtomicWeight;   // Nucleus atomic weight
+    double      fMassExcess;   // Nucleus mass excess in keV
+    const char* fSpinParity;   // Nucleus spin and parity
+    double      fSpin;         // Nucleus spin
+    const char* fParity;      // Nucleus parity
+                              //kinematic properties
+    double fKineticEnergy;
+    double fBeta;
+    double fGamma;
+    double fBrho;
+    double fTimeOfFlight;
+    double fVelocity;
+    TLorentzVector fEnergyImpulsion;
+    
+  public:
+    void  EnergyToBrho();
+    void  EnergyToTof();
+    void	BetaToVelocity();
+    void  BrhoToEnergy();
+    void  BrhoToTof()    {BrhoToEnergy(); EnergyToTof();}
+    void	TofToEnergy();
+    void	TofToBrho()    {TofToEnergy(); EnergyToBrho();}
+    void	EnergyToBeta();
+    void	BetaToEnergy();
+    void	BetaToGamma();
+    double DopplerCorrection(double EnergyLabGamma, double ThetaLabGamma);
+    
+    
+  protected :
+    void Extract(const char* line);
+    
+    public :
+    void				GetNucleusName();
+    string			GetName()			const				{return fNucleusName;}
+    int				GetZ()				const				{return fCharge;}
+    int				GetA()				const				{return fAtomicWeight;}
+    double			GetMassExcess()		const				{return fMassExcess;}
+    const char*		GetSpinParity()		const				{return fSpinParity;}
+    double			GetSpin()			const				{return fSpin;}
+    const char*		GetParity()			const				{return fParity;}
+    double			GetEnergy()			const				{return fKineticEnergy;}
+    double			GetBrho()			const				{return fBrho;}
+    double			GetTimeOfFlight()	const				{return fTimeOfFlight;}
+    double			GetBeta()			const				{return fBeta;}
+    double			GetGamma()			const				{return fGamma;}
+    double			GetVelocity()		const				{return fVelocity;}
+    TLorentzVector	GetEnergyImpulsion() const				{return fEnergyImpulsion;}
+    void				SetName(const char* name)				{fName = name;}
+    void				SetZ(int charge)						{fCharge = charge;}
+    void				SetA(int mass)							{fAtomicWeight = mass;}
+    void				SetMassExcess(double massexcess)		{fMassExcess = massexcess;}
+    void				SetSpinParity(const char* spinparity)	{fSpinParity = spinparity;}
+    void				SetSpin(double spin)					{fSpin = spin;}
+    void				SetParity(const char* parity)			{fParity = parity;}
+    void				SetKineticEnergy(double energy)			{fKineticEnergy = energy; EnergyToBrho(); EnergyToTof(); EnergyToBeta(); BetaToGamma();BetaToVelocity();}
+    void				SetBrho(double brho)					{fBrho = brho; BrhoToEnergy(); BrhoToTof(); EnergyToBeta(); BetaToGamma();BetaToVelocity();}
+    void				SetTimeOfFlight(double tof)				{fTimeOfFlight = tof; TofToEnergy(); TofToBrho(); EnergyToBeta(); BetaToGamma();BetaToVelocity();}
+    void				SetEnergyImpulsion(TLorentzVector P) 	{fEnergyImpulsion = P;
+      fKineticEnergy = P.E() - Mass();
+      EnergyToBrho();
+      EnergyToTof();
+      EnergyToBeta();
+      BetaToGamma();
+      BetaToVelocity();}
+    void				SetBeta(double beta)					{fBeta = beta; BetaToGamma(); BetaToEnergy(); EnergyToTof(); EnergyToBrho();BetaToVelocity();}
+    
+    // Nuclear mass in MeV
+    double      Mass() const {return (fAtomicWeight*uma + fMassExcess/1000.);}
+    void        Print() const   ;
+  };
+}
 #endif
diff --git a/NPLib/Physics/NPReaction.cxx b/NPLib/Physics/NPReaction.cxx
index a16674eb8..ea85665e3 100644
--- a/NPLib/Physics/NPReaction.cxx
+++ b/NPLib/Physics/NPReaction.cxx
@@ -42,7 +42,6 @@
 #include <vector>
 
 #include "NPReaction.h"
-#include "NPBeam.h"
 #include "NPOptionManager.h"
 #include "NPFunction.h"
 
@@ -50,15 +49,11 @@
 #include "CLHEP/Units/GlobalSystemOfUnits.h"
 #include "CLHEP/Units/PhysicalConstants.h"
 
-
-using namespace NPL;
-
-
 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
 Reaction::Reaction(){
   //------------- Default Constructor -------------
   
-  fNuclei1              = new Nucleus();
+  fNuclei1              = new Beam();
   fNuclei2              = new Nucleus();
   fNuclei3              = new Nucleus();
   fNuclei4              = new Nucleus();
@@ -258,9 +253,11 @@ void Reaction::ReadConfigurationFile(string Path){
       else if (DataBuffer=="Beam=") {
         check_Beam = true ;
         ReactionFile >> DataBuffer;
-        fNuclei1         = new Nucleus(DataBuffer);
-        Beam = DataBuffer;
-        if(fVerboseLevel==1) cout << "Beam " << fNuclei1->GetName() << endl;
+        // Pick up the beam energy from the Beam event generator
+        fNuclei1->SetVerboseLevel(0);
+        fNuclei1->ReadConfigurationFile(Path);
+        fBeamEnergy= fNuclei1->GetEnergy();
+        if(fVerboseLevel==1) cout << "Beam " << fNuclei1->GetName() << " @ " << fBeamEnergy << " MeV" << endl;
       }
       
       else if (DataBuffer=="Target=") {
@@ -352,13 +349,6 @@ void Reaction::ReadConfigurationFile(string Path){
         ReadingStatus = false;
     }
   }
-  
-  // Pick up the beam energy from the Beam event generator
-  NPL::Beam* localBeam= new NPL::Beam();
-  localBeam->SetVerboseLevel(0);
-  localBeam->ReadConfigurationFile(Path);
-  fBeamEnergy= localBeam->GetEnergy();
-  delete localBeam;
 
   // Modifiy the CS to shoot only within ]HalfOpenAngleMin,HalfOpenAngleMax[
   SetCSAngle(CSHalfOpenAngleMin,CSHalfOpenAngleMax);
diff --git a/NPLib/Physics/NPReaction.h b/NPLib/Physics/NPReaction.h
index ee58d18b7..d39be1835 100644
--- a/NPLib/Physics/NPReaction.h
+++ b/NPLib/Physics/NPReaction.h
@@ -34,12 +34,15 @@
  *                                                                           *
  *                                                                           *
  *****************************************************************************/
-  // C++ header
+// C++ header
 #include <string>
 
+// NPL
 #include "NPNucleus.h"
+#include "NPBeam.h"
+using namespace NPL;
 
-  // ROOT header
+// ROOT header
 #include "TLorentzVector.h"
 #include "TLorentzRotation.h"
 #include "TVector3.h"
@@ -68,7 +71,7 @@ namespace NPL{
     bool fshoot4;
     
   private:
-    Nucleus *fNuclei1;                 // Beam
+    Beam    *fNuclei1;                 // Beam
     Nucleus *fNuclei2;                 // Target
     Nucleus *fNuclei3;                 // Light ejectile
     Nucleus *fNuclei4;                 // Heavy ejectile
diff --git a/NPSimulation/src/EventGeneratorBeam.cc b/NPSimulation/src/EventGeneratorBeam.cc
index 4a8b3f36d..39be201db 100644
--- a/NPSimulation/src/EventGeneratorBeam.cc
+++ b/NPSimulation/src/EventGeneratorBeam.cc
@@ -78,7 +78,7 @@ void EventGeneratorBeam::GenerateEvent(G4Event* anEvent){
   //--------------write the DeDx Table -------------------
   if( anEvent->GetEventID()==0){
     // Define the particle to be shoot
-    m_particle = G4ParticleTable::GetParticleTable()->GetIon(m_Beam->GetNucleus()->GetZ(), m_Beam->GetNucleus()->GetA() , 0.);
+    m_particle = G4ParticleTable::GetParticleTable()->GetIon(m_Beam->GetZ(), m_Beam->GetA() , 0.);
 
     if(m_Target!=0 ){
       m_Target->WriteDEDXTable(m_particle ,0, m_Beam->GetEnergy()+4*m_Beam->GetSigmaEnergy());
diff --git a/NPSimulation/src/EventGeneratorParticleDecay.cc b/NPSimulation/src/EventGeneratorParticleDecay.cc
index 5d937dcc0..e50b250d1 100644
--- a/NPSimulation/src/EventGeneratorParticleDecay.cc
+++ b/NPSimulation/src/EventGeneratorParticleDecay.cc
@@ -388,7 +388,7 @@ void EventGeneratorParticleDecay::SetDecay(vector<string> DaughterName, vector<b
   
   m_DifferentialCrossSection = CSPath;
   if(CSPath!="TGenPhaseSpace") {
-    if(m_CrossSectionPath!="_void_"){
+    if(m_CrossSectionPath!="_void_")
       m_CrossSectionHist = Read1DProfile(m_CrossSectionPath,m_CrossSectionName);
       
       else{
-- 
GitLab