added fluctuations in multiplicity in fixed range
authorsnelling <snelling@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 8 May 2009 13:27:48 +0000 (13:27 +0000)
committersnelling <snelling@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 8 May 2009 13:27:48 +0000 (13:27 +0000)
PWG2/FLOW/AliFlowCommon/AliFlowEventSimpleMakerOnTheFly.cxx
PWG2/FLOW/AliFlowCommon/AliFlowEventSimpleMakerOnTheFly.h
PWG2/FLOW/macros/runFlowAnalysisOnTheFly.C

index a5c3469..d0490b8 100644 (file)
@@ -39,8 +39,11 @@ ClassImp(AliFlowEventSimpleMakerOnTheFly)
 
 
 AliFlowEventSimpleMakerOnTheFly::AliFlowEventSimpleMakerOnTheFly(UInt_t iseed):
+  fMultDistrOfRPsIsGauss(kFALSE),
   fMultiplicityOfRP(0),
   fMultiplicitySpreadOfRP(0.),
+  fMinMultOfRP(0),
+  fMaxMultOfRP(0),  
   fTemperatureOfRP(0.),  
   fUseConstantHarmonics(kFALSE),
   fV1RP(0.), 
@@ -55,12 +58,12 @@ AliFlowEventSimpleMakerOnTheFly::AliFlowEventSimpleMakerOnTheFly(UInt_t iseed):
   fPhiDistribution(NULL),
   fMyTRandom3(NULL),
   fCount(0),
-  fNoOfLoops(1)
- {
+  fNoOfLoops(1) 
+{
   // constructor
   fMyTRandom3 = new TRandom3(iseed);   
   gRandom->SetSeed(fMyTRandom3->Integer(65539));
- }
+}
 
 
 //========================================================================
@@ -72,7 +75,7 @@ AliFlowEventSimpleMakerOnTheFly::~AliFlowEventSimpleMakerOnTheFly()
  if (fPtSpectra) delete fPtSpectra;
  if (fPhiDistribution) delete fPhiDistribution;
  if (fMyTRandom3) delete fMyTRandom3;
-}
+}      
 
 
 //========================================================================
@@ -99,6 +102,7 @@ void AliFlowEventSimpleMakerOnTheFly::Init()
  fPhiDistribution->SetParName(1,"elliptic flow"); 
  fPhiDistribution->SetParName(2,"Reaction Plane");
  fPhiDistribution->SetParName(3,"harmonic 4"); // to be improved (name)
 }
 
 //========================================================================
@@ -111,9 +115,20 @@ AliFlowEventSimple* AliFlowEventSimpleMakerOnTheFly::CreateEventOnTheFly()
    
   // sampling the multiplicity:
   Int_t iNewMultiplicityOfRP = fMultiplicityOfRP;
-  if(fMultiplicitySpreadOfRP>0.0) iNewMultiplicityOfRP = (Int_t)fMyTRandom3->Gaus(fMultiplicityOfRP,fMultiplicitySpreadOfRP);
-  fPtSpectra->SetParameter(0,iNewMultiplicityOfRP);
   
+  if(fMultDistrOfRPsIsGauss) {
+    if (fMultiplicitySpreadOfRP>0.0) iNewMultiplicityOfRP = (Int_t)fMyTRandom3->Gaus(fMultiplicityOfRP,fMultiplicitySpreadOfRP);
+    fPtSpectra->SetParameter(0,iNewMultiplicityOfRP);
+  } else {
+    if (fMinMultOfRP != fMaxMultOfRP) {
+      iNewMultiplicityOfRP = (Int_t)fMyTRandom3->Uniform(fMinMultOfRP,fMaxMultOfRP);
+      fPtSpectra->SetParameter(0,iNewMultiplicityOfRP);
+    }
+    else {
+      fPtSpectra->SetParameter(0,fMinMultOfRP);
+    }
+  }
+    
   // set the 'temperature' of RPs
   fPtSpectra->SetParameter(1,fTemperatureOfRP);  
   
index ca0f767..888120e 100644 (file)
@@ -39,12 +39,21 @@ class AliFlowEventSimpleMakerOnTheFly {
   //                        *****************************
   //................................................................................................
   // setters and getters for global parameters:
+  void SetMultDistrOfRPsIsGauss(Bool_t const mdorig) {this->fMultDistrOfRPsIsGauss = mdorig;};
+  Bool_t GetMultDistrOfRPsIsGauss() const {return this->fMultDistrOfRPsIsGauss;};
+  
   void SetMultiplicityOfRP(Int_t multRP) {this->fMultiplicityOfRP = multRP;}
   Int_t GetMultiplicityOfRP() const {return this->fMultiplicityOfRP;} 
   
   void SetMultiplicitySpreadOfRP(Double_t multSpreadRP) {this->fMultiplicitySpreadOfRP = multSpreadRP;}
   Double_t GetMultiplicitySpreadOfRP() const {return this->fMultiplicitySpreadOfRP;} 
   
+  void SetMinMultOfRP(Int_t minmr) {this->fMinMultOfRP = minmr;}
+  Int_t GetMinMultOfRP() const {return this->fMinMultOfRP;} 
+  
+  void SetMaxMultOfRP(Int_t maxmr) {this->fMaxMultOfRP = maxmr;}
+  Int_t GetMaxMultOfRP() const {return this->fMaxMultOfRP;} 
+  
   void SetTemperatureOfRP(Double_t temperatureRP) {this->fTemperatureOfRP = temperatureRP;}
   Double_t GetTemperatureOfRP() const {return this->fTemperatureOfRP;} 
   
@@ -88,8 +97,14 @@ class AliFlowEventSimpleMakerOnTheFly {
   
   //................................................................................................
   // global parameters:
-  Int_t     fMultiplicityOfRP;       // multiplicity of RPs
-  Double_t  fMultiplicitySpreadOfRP; // multiplicity spread of RPs 
+  Bool_t    fMultDistrOfRPsIsGauss;  // 1.) if kTRUE  = multiplicitiy of RPs is sampled e-b-e from Gaussian distribution with
+                                     //                 mean = fMultiplicityOfRP and spread = fMultiplicitySpreadOfRP
+                                     // 2.) if kFALSE = multiplicitiy of RPs is sampled e-b-e uniformly from 
+                                     //                 interval [fMinMultOfRP,fMaxMultOfRP]
+  Int_t     fMultiplicityOfRP;       // mean multiplicity of RPs (if sampled from Gaussian)
+  Double_t  fMultiplicitySpreadOfRP; // multiplicity spread of RPs (if sampled from Gaussian)
+  Int_t     fMinMultOfRP;            // minimal multiplicity of RPs (if sampled uniformly)
+  Int_t     fMaxMultOfRP;            // maximum multiplicity of RPs (if sampled uniformly)
   Double_t  fTemperatureOfRP;        // "temperature" of RPs in GeV/c (increase this parameter to get more high pt RPs) 
   Bool_t    fUseConstantHarmonics;   // harmonics V1, V2, V4... are constant (kTRUE) or functions of pt and eta (kFALSE)     
   // constant harmonics: 
@@ -113,6 +128,7 @@ class AliFlowEventSimpleMakerOnTheFly {
   TRandom3* fMyTRandom3; // our TRandom3 generator
   Int_t     fCount;      // count number of events 
   Int_t     fNoOfLoops;  // number of times to use the same particle (nonflow)
+  
 
   ClassDef(AliFlowEventSimpleMakerOnTheFly,0) // macro for rootcint
 };
index 2fdb009..4856f5b 100644 (file)
@@ -13,7 +13,7 @@ Bool_t LYZEP = kFALSE;
 Bool_t GFC   = kTRUE;
 Bool_t QC    = kTRUE;
 Bool_t FQD   = kTRUE;
-Bool_t MCEP  = kTRUE; 
+Bool_t MCEP  = kTRUE;
 //--------------------------------------------------------------------------------------
 
 // Weights 
@@ -27,12 +27,22 @@ Bool_t bSameSeed = kFALSE; // use always the same seed for random generators.
                            // usage od same seed (kTRUE) is relevant in two cases:
                            // 1.) If you want to use LYZ method to calcualte differential flow;
                            // 2.) If you want to use phi weights for GFC, QC and FQD
+                           
 Bool_t bConstantHarmonics = kTRUE; // harmonics V1, V2, V4... are constant (kTRUE) or functions of pt and eta (kFALSE)
 
 Int_t iLoops = 1; // number of times to use each track (to simulate nonflow)
 
-Int_t iMultiplicityOfRP = 500; // multiplicity of RPs
-Double_t dMultiplicitySpreadOfRP = 0; // multiplicity spread of RPs
+Bool_t bMultDistrOfRPsIsGauss = kTRUE; // 1.) if kTRUE  = multiplicitiy of RPs is sampled e-b-e from Gaussian distribution with
+                                        //                 mean = iMultiplicityOfRP and spread = dMultiplicitySpreadOfRP
+                                        // 2.) if kFALSE = multiplicitiy of RPs is sampled e-b-e uniformly from 
+                                        //                 interval [iMinMultOfRP,iMaxMultOfRP]
+                                        // 3.) for a fixed multiplicity use Gaussian with zero spread or use uniform with iMinMult=iMaxMult
+                                         
+Int_t iMultiplicityOfRP = 500;        // mean multiplicity of RPs (if sampled from Gaussian)
+Double_t dMultiplicitySpreadOfRP = 0; // multiplicity spread of RPs (if sampled from Gaussian)
+Int_t iMinMultOfRP = 400;             // minimal multiplicity of RPs (if sampled uniformly)
+Int_t iMaxMultOfRP = 600;             // maximal multiplicity of RPs (if sampled uniformly)
+
 Double_t dTemperatureOfRP = 0.44; // 'temperature' of RPs in GeV/c (increase this parameter to get more high pt RPs) 
 
 //......................................................................................  
@@ -58,7 +68,7 @@ enum anaModes {mLocal,mLocalSource,mLocalPAR};
 // mLocalPAR: Analyze data on your computer using root + PAR files
 // mLocalSource: Analyze data on your computer using root + source files
                                           
-int runFlowAnalysisOnTheFly(Int_t mode=mLocal, Int_t nEvts=1000)
+int runFlowAnalysisOnTheFly(Int_t mode=mLocal, Int_t nEvts=440)
 {
  TStopwatch timer;
  timer.Start();
@@ -223,8 +233,19 @@ int runFlowAnalysisOnTheFly(Int_t mode=mLocal, Int_t nEvts=1000)
   
  // set the global event parameters: 
  eventMakerOnTheFly->SetNoOfLoops(iLoops);
- eventMakerOnTheFly->SetMultiplicityOfRP(iMultiplicityOfRP);
- eventMakerOnTheFly->SetMultiplicitySpreadOfRP(dMultiplicitySpreadOfRP);
+ if(bMultDistrOfRPsIsGauss)
+ {
+  eventMakerOnTheFly->SetMultDistrOfRPsIsGauss(bMultDistrOfRPsIsGauss);
+  eventMakerOnTheFly->SetMultiplicityOfRP(iMultiplicityOfRP);
+  eventMakerOnTheFly->SetMultiplicitySpreadOfRP(dMultiplicitySpreadOfRP);
+ } else
+   {
+    eventMakerOnTheFly->SetMultDistrOfRPsIsGauss(bMultDistrOfRPsIsGauss);
+    eventMakerOnTheFly->SetMinMultOfRP(iMinMultOfRP);
+    eventMakerOnTheFly->SetMaxMultOfRP(iMaxMultOfRP); 
+   }
+  
  eventMakerOnTheFly->SetTemperatureOfRP(dTemperatureOfRP);
 
  eventMakerOnTheFly->SetV1RP(dV1RP);