Correct handling of seeds for MC on the fly trains
authorjgrosseo <jgrosseo@cern.ch>
Wed, 3 Dec 2014 16:43:42 +0000 (17:43 +0100)
committerjgrosseo <jgrosseo@cern.ch>
Wed, 3 Dec 2014 16:43:42 +0000 (17:43 +0100)
Adding AliGenerator::SetSeed which AliFatal by default
It has to be overwritten in the deriving classes

Implemented in PYTHIA6, 8, Hijing, AMPT
Change of signature in TUHKMgen which already had SetSeed function

16 files changed:
PYTHIA6/AliGenPythia.cxx
PYTHIA6/AliGenPythia.h
PYTHIA6/AliGenPythiaPlus.cxx
PYTHIA6/AliGenPythiaPlus.h
PYTHIA6/AliPythiaBase.cxx
PYTHIA6/AliPythiaBase.h
PYTHIA8/AliPythia8.cxx
PYTHIA8/AliPythia8.h
STEER/STEER/AliGenerator.cxx
STEER/STEER/AliGenerator.h
STEER/STEER/AliMCGenHandler.cxx
TAmpt/AliGenAmpt.cxx
TAmpt/AliGenAmpt.h
THijing/AliGenHijing.cxx
THijing/AliGenHijing.h
TUHKMgen/AliGenUHKM.h

index 5ad46d2..9ea21ed 100644 (file)
@@ -2672,3 +2672,7 @@ Bool_t AliGenPythia::TriggerOnSelectedParticles(Int_t np)
   return ok;
 }
 
+void AliGenPythia::SetSeed(UInt_t seed)
+{
+  GetRandom()->SetSeed(seed);
+}
index 01d05d1..9e69b9a 100644 (file)
@@ -45,6 +45,7 @@ class AliGenPythia : public AliGenMC
     // Select process type
     virtual void    SetProcess(Process_t proc = kPyCharm) {fProcess = proc;}
     virtual void    SetTune(Int_t itune) {fItune = itune;}
+    virtual void    SetSeed(UInt_t seed);
 
     // Select structure function
     virtual void    SetStrucFunc(StrucFunc_t func =  kCTEQ5L) {fStrucFunc = func;}
index 76b6fce..038666c 100644 (file)
@@ -508,6 +508,11 @@ void AliGenPythiaPlus::Init()
 //    }
 }
 
+void AliGenPythiaPlus::SetSeed(UInt_t seed)
+{
+  fPythia->SetSeed(seed);
+}
+
 void AliGenPythiaPlus::Generate()
 {
 // Generate one event
@@ -1486,4 +1491,3 @@ void AliGenPythiaPlus::Streamer(TBuffer &R__b)
 #endif
 
 
-
index bd13f8e..e4d3a30 100644 (file)
@@ -39,6 +39,8 @@ class AliGenPythiaPlus : public AliGenMC
     virtual ~AliGenPythiaPlus();
     virtual void    Generate();
     virtual void    Init();
+    virtual void    SetSeed(UInt_t seed);
+    
     // Range of events to be printed
     virtual void    SetEventListRange(Int_t eventFirst=-1, Int_t eventLast=-1);
     // Select process type
index 898d195..41af0f3 100644 (file)
@@ -17,6 +17,7 @@
 /* $Id$ */
 
 #include "AliPythiaBase.h"
+#include <cstdlib>
 
 ClassImp(AliPythiaBase)
 
@@ -25,5 +26,10 @@ AliPythiaBase::AliPythiaBase()
     // Default constructor
 }
 
+void AliPythiaBase::SetSeed(UInt_t)
+{
+  Printf("FATAL -> SetSeed not implemented");
+  ::abort();
+}
 
 
index b4a28b6..ef2a308 100644 (file)
@@ -26,6 +26,7 @@ class AliPythiaBase : public AliRndm
     virtual Int_t CheckedLuComp(Int_t /*kf*/) {return -1;}   
     // Pythia initialisation for selected processes
     virtual void  ProcInit (Process_t /*process*/, Float_t /*energy*/, StrucFunc_t /*strucfunc*/, Int_t /* tune */) {;}
+    virtual void  SetSeed(UInt_t seed);
     virtual void  GenerateEvent() {;}
     virtual void  GenerateMIEvent() {;}
     virtual Int_t GetNumberOfParticles() {return -1;};
index 45030e1..9f642a3 100644 (file)
@@ -543,6 +543,16 @@ void AliPythia8::ProcInit(Process_t process, Float_t energy, StrucFunc_t strucfu
     Initialize(2212, 2212, fEcms);
 }
 
+void AliPythia8::SetSeed(UInt_t seed)
+{
+  //
+  // set seed in PYTHIA8
+  // NB. 900000000 is the maximum seed (0 is not allowed)
+  //
+  ReadString("Random:setSeed = on");
+  ReadString(Form("Random:seed = %d", (seed % 900000000) + 1));
+}
+
 void AliPythia8::SetNuclei(Int_t /*a1*/, Int_t /*a2*/)
 {
 // Treat protons as inside nuclei with mass numbers a1 and a2  
index 53ee49c..d7aff47 100644 (file)
@@ -22,6 +22,7 @@ class AliPythia8 :public AliTPythia8, public AliPythiaBase
     virtual Int_t CheckedLuComp(Int_t /*kf*/) {return -1;}
     // Pythia initialisation for selected processes
     virtual void  ProcInit (Process_t process, Float_t energy, StrucFunc_t strucfunc, Int_t tune);
+    virtual void  SetSeed(UInt_t seed);
     virtual void  GenerateEvent();
     virtual void  GenerateMIEvent();
     virtual void  HadronizeEvent();
index 7dfba31..21a8c7b 100644 (file)
@@ -458,3 +458,13 @@ void AliGenerator::FinishRun()
 {
     ;
 }
+
+//_______________________________________________________________________
+void AliGenerator::SetSeed(UInt_t seed)
+{
+  // 
+  // function to set the seed in the random number generator used by this generator
+  //
+  
+  AliFatal("Not implemented!");
+}
index 0bfd1be..dd0a9b7 100644 (file)
@@ -108,6 +108,8 @@ class AliGenerator : public TNamed, public AliRndm
     virtual void    SetTarget(TString tar="", Int_t a = 0, Int_t z = 0)
        {fTarget = tar; fATarget = a; fZTarget = z;}
     
+    virtual void    SetSeed(UInt_t seed);
+    
  protected:
     virtual  void  PushTrack(Int_t done, Int_t parent, Int_t pdg,
                                Float_t *pmom, Float_t *vpos, Float_t *polar,
index fcab1b5..ead4535 100644 (file)
@@ -149,6 +149,7 @@ Bool_t AliMCGenHandler::Init(Option_t* /*opt*/)
 
       Printf("AliMCGenHandler::Init: Using seed: %d", fSeed);
       gRandom->SetSeed(fSeed);
+      fGenerator->SetSeed(fSeed);
     }
 
     AliRunLoader* rl = AliRunLoader::Open("galice.root","FASTRUN","recreate");
index 7b65574..24f61b2 100644 (file)
@@ -234,6 +234,11 @@ void AliGenAmpt::Init()
   fAmpt->SetReactionPlaneAngle(0.0);
 }
 
+void AliGenAmpt::SetSeed(UInt_t seed)
+{
+  AliAmptRndm::GetAmptRandom()->SetSeed(seed);
+}
+
 void AliGenAmpt::Generate()
 {
   // Generate one event
index 2ef9b96..a55618d 100644 (file)
@@ -27,6 +27,7 @@ class AliGenAmpt : public AliGenMC
     virtual TAmpt  *Ampt() { return fAmpt; }
     virtual void    Generate();
     virtual void    Init();
+    virtual void    SetSeed(UInt_t seed);
     virtual void    SetEnergyCMS(Float_t energy=5500) {fEnergyCMS=energy;}
     virtual void    SetReferenceFrame(TString frame="CMS")
        {fFrame=frame;}
index 43cbbd8..e97f778 100644 (file)
@@ -223,6 +223,11 @@ void AliGenHijing::Init()
 //
 }
 
+void AliGenHijing::SetSeed(UInt_t seed)
+{
+  AliHijingRndm::GetHijingRandom()->SetSeed(seed);
+}
+
 void AliGenHijing::Generate()
 {
 // Generate one event
index 0dd4366..8db8161 100644 (file)
@@ -29,6 +29,7 @@ class AliGenHijing : public AliGenMC
     virtual ~AliGenHijing();
     virtual void    Generate();
     virtual void    Init();
+    virtual void    SetSeed(UInt_t seed);
     // set centre of mass energy
     virtual void    SetEnergyCMS(Float_t energy=5500) {fEnergyCMS=energy;}
     virtual void    SetReferenceFrame(TString frame="CMS")
index 401f826..4d9f0f4 100755 (executable)
@@ -48,7 +48,7 @@ class AliGenUHKM : public AliGenMC
   void SetMuQ(Double_t value) {fHydjetParams.fMuI3 = value;}  // Isospin chemical potential [GeV]
   void SetThFrzTemperature(Double_t value) {fHydjetParams.fThFO = value;}  // Temperature for the thermal freezeout [GeV]
   void SetMuPionThermal(Double_t value) {fHydjetParams.fMu_th_pip = value;} // Chemical potential for pi+ at thermal freezeout [GeV]
-  void SetSeed(Int_t value) {fHydjetParams.fSeed = value;} //parameter to set the random nuber seed (=0 the current time is used
+  virtual void SetSeed(UInt_t value) {fHydjetParams.fSeed = value;} //parameter to set the random nuber seed (=0 the current time is used
   //to set the random generator seed, !=0 the value fSeed is
   //used to set the random generator seed and then the state of random
   //number generator in PYTHIA MRPY(1)=fSeed