]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TUHKMgen/TUHKMgen.cxx
Xrootd verbosity switched on for AliAnalysisManager debug level >2
[u/mrichter/AliRoot.git] / TUHKMgen / TUHKMgen.cxx
index 967f0a9f93144ce35f25f1e19effa3f4d68e92f8..cf81d4a32a891f50e0847c5dedeb62332e80b79b 100755 (executable)
@@ -15,7 +15,7 @@
 //                                                                         //
 ////////////////////////////////////////////////////////////////////////////
 
-#ifndef TUHKMgen_H
+#ifndef TUHKMGEN_H
 #include "TUHKMgen.h"
 #endif
 #include "TObjArray.h"
 #include <iostream>
 using namespace std;
 
+//class TObjArray;
+
 ClassImp(TUHKMgen)
 
 TUHKMgen::TUHKMgen() : 
   TGenerator("UHKM","UHKM"),
   fInitialState(0x0),
-  fUseCharmParticles(kFALSE),
-  fMinWidth(0.0),
-  fMaxWidth(10.0),
-  fMinMass(0.0001),
-  fMaxMass(10.0)
+  fAllocator(),
+  fSecondariesList(),
+  fNPprim(0),
+  fNPsec(0),
+  fHydjetParams(),
+  fStableFlagged(0)
 {
-  cout << "TUHKMgen::TUHKMgen() IN" << endl;
+  // default constructor setting reasonable defaults for initial parameters (central Pb+Pb at 5.5 TeV)
+
+  //  ParticleAllocator fAllocator;
+  //  List_t fSecondariesList;
 
   // Set reasonable default values for LHC
   
@@ -122,20 +128,21 @@ TUHKMgen::TUHKMgen() :
     fStableFlagStatus[i] = kFALSE;
   }
   fStableFlagged = 0;
-  cout << "TUHKMgen::TUHKMgen() OUT" << endl;
+  //  cout << "TUHKMgen::TUHKMgen() OUT" << endl;
 }
 
 //______________________________________________________________________________
 TUHKMgen::~TUHKMgen()
 {
-  cout << "TUHKMgen::~TUHKMgen() IN" << endl;
+  // destructor, deletes the InitialStateHydjet object
   if(fInitialState)
     delete fInitialState;
-  cout << "TUHKMgen::~TUHKMgen() OUT" << endl;
 }
 
 void TUHKMgen::SetAllParametersRHIC()
 {
+  // Set reasonable input parameters for 0-5% central Au+Au collisions
+  // at 200 GeV at RHIC
   SetEcms(200.0);                // RHIC top energy
   SetAw(197);                    // Au+Au
   SetIfb(1);
@@ -171,6 +178,7 @@ void TUHKMgen::SetAllParametersRHIC()
 
 void TUHKMgen::SetAllParametersLHC()
 {
+  // Set reasonable input parameters for 0-5% Pb+Pb collisions at 5.5 TeV at LHC
   SetEcms(5500.0);               // LHC
   SetAw(207);                    // Pb+Pb
   SetIfb(1);
@@ -204,16 +212,16 @@ void TUHKMgen::SetAllParametersLHC()
   SetPyquenIanglu(0);            // small gluon angular distribution
 }
 
-TObjArray* TUHKMgen::ImportParticles(Option_t *option)
+TObjArray* TUHKMgen::ImportParticles(const Option_t *)
 {
-  cout << "TObjArray* TUHKMgen::ImportParticles(Option_t) IN" << endl;
+  // Function overloading the TGenerator::ImportParticles() member function.
+  // The particles from the local particle list (fSecondariesList) are
+  // forwarded to the TGenerator::fParticles 
   fParticles->Clear();
-  cout << "TObjArray* TUHKMgen::ImportParticles(Option_t): UHKM stack contains " << fNPsec << " particles." << endl;
   Int_t nump = 0;
   LPIT_t it,e;
-  Double_t Fm2Gev = 1.;//0.197327; //UHKM has coordinates in fm
    
-  for(it = fSecondariesList.begin(), e = fSecondariesList.end(); it != e; ++it) {
+  for(it = fSecondariesList.begin(), e = fSecondariesList.end(); it != e; it++) {
     TVector3 pos(it->Pos().Vect());
     TVector3 mom(it->Mom().Vect());
     Float_t m1 = it->TableMass();
@@ -221,17 +229,14 @@ TObjArray* TUHKMgen::ImportParticles(Option_t *option)
     Int_t im2 = -1;
     Int_t id1 = -1;
     Int_t id2 = -1;
-    Int_t id3 = -1;
 
-    Int_t nd = it->GetNDaughters();
-    
     Int_t type = it->GetType();  // 0-hydro, 1-jets
 
     if (im1> -1) {
-      TParticle *mother = (TParticle*) (fParticles->UncheckedAt(im1));    
+      TParticle *mother = (TParticle*) (fParticles->UncheckedAt(im1+1));          
       mother->SetLastDaughter(nump);
       if (mother->GetFirstDaughter()==-1)
-       mother->SetFirstDaughter(nump);
+       mother->SetFirstDaughter(nump+1);
     }
 
     nump++;
@@ -239,119 +244,113 @@ TObjArray* TUHKMgen::ImportParticles(Option_t *option)
     TParticle* p = new TParticle(it->Encoding(), type,                         //pdg,stat
                                 im1, im2, id1, id2,                                            //m1,m2,d1,d2
                                 mom[0], mom[1], mom[2], TMath::Sqrt(mom.Mag2() + m1 * m1),     //px,py,pz,e
-                                pos[0]*1.e-13*Fm2Gev, pos[1]*1.e-13*Fm2Gev, pos[2]*1.e-13*Fm2Gev
-                                it->T()*1.e-13*Fm2Gev/3e10);                           //x,y,z,t
+                                pos[0]*1.e-13, pos[1]*1.e-13, pos[2]*1.e-13
+                                it->T()*1.e-13/3e10);                          //x,y,z,t
   
      p->SetUniqueID(nump);
      fParticles->Add(p);
   } //end for
  
-  fAllocator.FreeList(fSourceList);
   fAllocator.FreeList(fSecondariesList);
-  cout << "TObjArray* TUHKMgen::ImportParticles(Option_t) OUT" << endl;
+
   return fParticles;
 }
 
-Int_t TUHKMgen::ImportParticles(TClonesArray *particles,Option_t* option)
+Int_t TUHKMgen::ImportParticles(TClonesArray *particles, const Option_t* option)
 {
-  //z-start import
-  cout << "Int_t TUHKMgen::ImportParticles(TClonesArray*, Option_t) IN" << endl;
+  // Function overloading the TGenerator::ImportParticles() member function.
+  // The particles from the local particle list (fSecondariesList) are
+  // forwarded to the TGenerator::fParticles
+  option = option;   // just to avoid the warning
+  
+
   if(particles==0) return 0;
   TClonesArray &particlesR=*particles;
   particlesR.Clear();
-
  
   Int_t numprim,numsec;  numprim=numsec=0;
   Int_t nump = 0;
   LPIT_t it,e;
-  Double_t Fm2Gev = 1.; //UHKM has coordinates in fm //0.197327;
-
-  cout << "TUHKMgen::ImportParticles() option(All or Sec) = " << option << endl;
-  for(it = fSecondariesList.begin(), e = fSecondariesList.end(); it != e; ++it) {
-    //      cout << "TUHKMgen::ImportParticles() import particle pdg(" << it->Encoding() << ")" << endl;
-    TVector3 pos(it->Pos().Vect());
-    TVector3 mom(it->Mom().Vect());
-    Float_t m1 = it->TableMass();
-    Int_t im1 = it->GetMother();
-    Int_t nd = 0;
+  
+  for(it = fSecondariesList.begin(), e = fSecondariesList.end(); it != e; it++) {
+    TVector3 pos(it->Pos().Vect());  
+    TVector3 mom(it->Mom().Vect());  
+    Float_t m1 = it->TableMass();  
+    Int_t im1 = it->GetMother();  
     Int_t im2 = -1;
     Int_t id1 = -1;
     Int_t id2 = -1;
-    Int_t id3 = -1;
-    
-    nd = it->GetNDaughters();
-
-    Int_t type = it->GetType();  // 0-hydro, 1-jets
     
+    Int_t type = it->GetType();  
+      
     if (im1> -1) {
-      TParticle *mother = (TParticle*) (particlesR.UncheckedAt(im1+1));           
+     // particle not a primary -> set the daughter indexes for the mother particle"<< endl;
+      TParticle *mother = (TParticle*) (particlesR.UncheckedAt(im1));
       mother->SetLastDaughter(nump);
-      if (mother->GetFirstDaughter()==-1)
-       mother->SetFirstDaughter(nump);
+      if(mother->GetFirstDaughter()==-1)
+       mother->SetFirstDaughter(nump);
     }
       
-    nump++;
-
+    
     new (particlesR[nump]) TParticle(it->Encoding(), type,                                             //pdg,stat
                                     im1, im2, id1, id2,                                                //m1,m2,d1,d2
                                     mom[0], mom[1], mom[2], TMath::Sqrt(mom.Mag2() + m1 * m1), //px,py,pz,e
-                                    pos[0]*1.e-13*Fm2Gev, pos[1]*1.e-13*Fm2Gev, pos[2]*1.e-13*Fm2Gev, 
-                                    it->T()*1.e-13*Fm2Gev/3e10);                       //x,y,z,t
-
-
+                                    pos[0]*1.e-13, pos[1]*1.e-13, pos[2]*1.e-13, 
+                                    it->T()*1.e-13/3e10);                      //x,y,z,t
+    
     particlesR[nump]->SetUniqueID(nump);
+    nump++;
     numsec++;
   }//end for
   
-  fAllocator.FreeList(fSourceList);
-  fAllocator.FreeList(fSecondariesList);
+  fSecondariesList.clear();
   printf("Scan and add prim %d sec %d and all %d particles\n",
-        numprim,numsec,nump);
-  cout << "Int_t TUHKMgen::ImportParticles(TClonesArray*, Option_t) OUT" << endl;
+        numprim,numsec,nump);
   return nump;
 }
 
 //______________________________________________________________________________
 void TUHKMgen::Initialize()
 {
-  cout << "TUHKMgen::Initialize() IN" << endl;
-  fInitialState = new InitialStateHydjet();
-  SetAllParameters();
-  fInitialState->LoadPDGInfo();
+  // Function overloading the TGenerator::Initialize() member function.
+  // The Monte-Carlo model is initialized (input parameters are transmitted, 
+  // particle list and decay channels are loaded, average multiplicities are calculated, etc.)
+  fInitialState = new InitialStateHydjet();  
+  SetAllParameters();  
+  fInitialState->LoadPDGInfo();  
   // set the stable flags
   for(Int_t i=0; i<fStableFlagged; i++) 
-    fInitialState->SetPDGParticleStable(fStableFlagPDG[i], fStableFlagStatus[i]);
+    fInitialState->SetPDGParticleStable(fStableFlagPDG[i], fStableFlagStatus[i]); 
 
-  if(!fInitialState->MultIni())
-    Error("TUHKMgen::Initialize()", "Bad status return from MultIni(). Check it out!! \n");
-  cout << "TUHKMgen::Initialize() OUT" << endl;
+  if(!fInitialState->MultIni()) 
+    Error("TUHKMgen::Initialize()", "Bad status return from MultIni(). Check it out!! \n");  //
 }
 
-void TUHKMgen::Print(
+void TUHKMgen::Print(const Option_t*) const
 {
   cout << "TUHKMgen::Print() method not implemented yet!!" << endl;
 }
 
 void TUHKMgen::GenerateEvent()
 {
-  cout << "TUHKMgen::GenerateEvent() IN" << endl;
-  // Generate the initial particles
-  fInitialState->Initialize(fSourceList, fAllocator);
-  cout << "TUHKMgen::fInitialState->Initialize" << endl;
+  // Member function overloading the TGenerator::GenerateEvent()
+  // The HYDJET++ model is run and the particle lists (fSourceList and fSecondariesList) are filled
+  
+  fInitialState->Initialize(fSecondariesList, fAllocator);  
  
-  if(fSourceList.empty()) 
-    Error("TUHKMgen::GenerateEvent()", "Source particle list empty after fireball initialization!! \n");
+  if(fSecondariesList.empty())Error("TUHKMgen::GenerateEvent()", "Source particle list empty after fireball initialization!! \n");  //
 
   // Run the decays
-  cout << "after Initilize: get_time, weak_decay_limit" <<fInitialState->GetTime()<<fInitialState->GetWeakDecayLimit()<< endl; 
-  if(fInitialState->GetTime() >= 0.)
-    fInitialState->Evolve(fSourceList, fSecondariesList, fAllocator, fInitialState->GetWeakDecayLimit());
-  cout << "TUHKMgen::GenerateEvent() OUT" << endl;
+
+  if(fInitialState->RunDecays())  
+    fInitialState->Evolve(fSecondariesList, fAllocator, fInitialState->GetWeakDecayLimit());  
+  
 }
 
 void TUHKMgen::SetAllParameters() {
-  cout << "TUHKMgen::SetAllParameters() IN" << endl;
+  // forward all model input parameters to the InitialStateHydjet object
+  // which will handle all the Monte-Carlo simulation using HYDJET++ model
 
   fInitialState->fParams.fSqrtS = fHydjetParams.fSqrtS;
   fInitialState->fParams.fAw = fHydjetParams.fAw;
@@ -375,8 +374,8 @@ void TUHKMgen::SetAllParameters() {
   fInitialState->fParams.fDelta = fHydjetParams.fDelta;
   fInitialState->fParams.fEpsilon = fHydjetParams.fEpsilon;
 
- // fInitialState->fParams.fTime = fHydjetParams.fTime;
   fInitialState->fParams.fWeakDecay = fHydjetParams.fWeakDecay;
+  fInitialState->fParams.fDecay = 1;                   // run decays
   fInitialState->fParams.fCorrS = fHydjetParams.fCorrS;
   fInitialState->fParams.fEtaType = fHydjetParams.fEtaType;
  
@@ -391,10 +390,10 @@ void TUHKMgen::SetAllParameters() {
 
   fInitialState->SetPDGParticleFilename(fParticleFilename);
   fInitialState->SetPDGDecayFilename(fDecayFilename);
-  fInitialState->SetUseCharmParticles(fUseCharmParticles);
-  fInitialState->SetWidthRange(fMinWidth, fMaxWidth);
-  fInitialState->SetMassRange(fMinMass, fMaxMass);
-  //for(Int_t i=0; i<fStableFlagged; i++) fInitialState->SetPDGParticleStable(fStableFlagPDG[i], fStableFlagStatus[i]);
-  cout << "TUHKMgen::SetAllParameters() OUT" << endl;
+  //  fInitialState->SetUseCharmParticles(fUseCharmParticles);
+  //  fInitialState->SetWidthRange(fMinWidth, fMaxWidth);
+  //  fInitialState->SetMassRange(fMinMass, fMaxMass);
+  //  for(Int_t i=0; i<fStableFlagged; i++) 
+  //    fInitialState->SetPDGParticleStable(fStableFlagPDG[i], fStableFlagStatus[i]);
+  //  cout << "TUHKMgen::SetAllParameters() OUT" << endl;
 }
-