AliFlowEventSimpleMakerOnTheFly::AliFlowEventSimpleMakerOnTheFly(UInt_t iseed):
+ fMultDistrOfRPsIsGauss(kFALSE),
fMultiplicityOfRP(0),
fMultiplicitySpreadOfRP(0.),
- fPtFormula(NULL),
- fPhiFormula(NULL),
+ fMinMultOfRP(0),
+ fMaxMultOfRP(0),
+ fTemperatureOfRP(0.),
+ fUseConstantHarmonics(kFALSE),
+ fV1RP(0.),
+ fV1SpreadRP(0.),
+ fV2DistrOfRPsIsGauss(kFALSE),
+ fV2RP(0.),
+ fV2SpreadRP(0.),
+ fMinV2RP(0.),
+ fMaxV2RP(0.),
+ fV4RP(0.),
+ fV4SpreadRP(0.),
+ fV2RPMax(0.),
+ fPtCutOff(0.),
+ fPhiMin1(0.),
+ fPhiMax1(0.),
+ fProbability1(0.),
+ fPhiMin2(0.),
+ fPhiMax2(0.),
+ fProbability2(0.),
+ fPtSpectra(NULL),
+ fPhiDistribution(NULL),
fMyTRandom3(NULL),
- fCount(0)
- {
+ fCount(0),
+ fNoOfLoops(1),
+ fEtaMinA(-1.0),
+ fEtaMaxA(-0.01),
+ fEtaMinB(0.01),
+ fEtaMaxB(1.0)
+{
// constructor
- fMyTRandom3 = new TRandom3(iseed);
- }
+ fMyTRandom3 = new TRandom3(iseed);
+ gRandom->SetSeed(fMyTRandom3->Integer(65539));
+}
//========================================================================
AliFlowEventSimpleMakerOnTheFly::~AliFlowEventSimpleMakerOnTheFly()
{
// destructor
- if (fPtFormula) delete fPtFormula;
- if (fPhiFormula) delete fPhiFormula;
- if (fMyTRandom3) delete fMyTRandom3;
-}
+ if (fPtSpectra) delete fPtSpectra;
+ if (fPhiDistribution) delete fPhiDistribution;
+ if (fMyTRandom3) delete fMyTRandom3;
+}
//========================================================================
+void AliFlowEventSimpleMakerOnTheFly::Init()
+{
+ // define the pt spectra and phi distribution
+
+ // pt spectra of pions (Boltzman):
+ Double_t dPtMin = 0.; // to be improved (move this to the body of contstructor?)
+ Double_t dPtMax = 10.; // to be improved (move this to the body of contstructor?)
+
+ fPtSpectra = new TF1("fPtSpectra","[0]*x*TMath::Exp(-pow(0.13957*0.13957+x*x,0.5)/[1])",dPtMin,dPtMax);
+ fPtSpectra->SetParName(0,"Multiplicity of RPs");
+ fPtSpectra->SetParName(1,"Temperature of RPs");
+
+ // phi distribution:
+ Double_t dPhiMin = 0.; // to be improved (move this to the body of contstructor?)
+ Double_t dPhiMax = TMath::TwoPi(); // to be improved (move this to the body of contstructor?)
+
+ fPhiDistribution = new TF1("fPhiDistribution","1+2.*[0]*TMath::Cos(x-[2])+2.*[1]*TMath::Cos(2*(x-[2]))+2.*[3]*TMath::Cos(4*(x-[2]))",dPhiMin,dPhiMax);
+ fPhiDistribution->SetParName(0,"directed flow");
+ fPhiDistribution->SetParName(1,"elliptic flow");
+ fPhiDistribution->SetParName(2,"Reaction Plane");
+ fPhiDistribution->SetParName(3,"harmonic 4"); // to be improved (name)
+
+}
+
+//========================================================================
+
AliFlowEventSimple* AliFlowEventSimpleMakerOnTheFly::CreateEventOnTheFly()
{
// method to create event on the fly
AliFlowEventSimple* pEvent = new AliFlowEventSimple(fMultiplicityOfRP);
+
+ // sampling the multiplicity:
+ Int_t iNewMultiplicityOfRP = fMultiplicityOfRP;
- //reaction plane
- Double_t fMCReactionPlaneAngle = fMyTRandom3->Uniform(0.,TMath::TwoPi());
-
- // pt:
- Double_t dPtMin = 0.; // to be improved
- Double_t dPtMax = 10.; // to be improved
-
- fPtFormula = new TF1("PtFormula","[0]*x*TMath::Exp(-x*x)",dPtMin,dPtMax);
-
- fPtFormula->SetParName(0,"Multiplicity of RPs");
- fPtFormula->SetParameter(0,fMultiplicityOfRP);
+ 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);
+ // sampling the reaction plane
+ Double_t dMCReactionPlaneAngle = fMyTRandom3->Uniform(0.,TMath::TwoPi());
+ fPhiDistribution->SetParameter(2,dMCReactionPlaneAngle);
- // phi:
- Double_t dPhiMin = 0.; // to be improved
- Double_t dPhiMax = TMath::TwoPi(); // to be improved
-
- fPhiFormula = new TF1("phiDistribution","(1)*(1+2.*[0]*TMath::Cos(2*x))",dPhiMin,dPhiMax);
-
- Double_t dV2 = 0.044; // to be improved
-
- fPhiFormula->SetParName(0,"elliptic flow");
- fPhiFormula->SetParameter(0,dV2);
+ // sampling the V1:
+ Double_t dNewV1RP=fV1RP;
+ if(fV1SpreadRP>0.0) {dNewV1RP = fMyTRandom3->Gaus(fV1RP,fV1SpreadRP);}
+ fPhiDistribution->SetParameter(0,dNewV1RP);
+ // sampling the V2:
+ if(fUseConstantHarmonics)
+ {
+ Double_t dNewV2RP = fV2RP;
+ if(fV2DistrOfRPsIsGauss)
+ {
+ if(fV2SpreadRP>0.0) dNewV2RP = fMyTRandom3->Gaus(fV2RP,fV2SpreadRP);
+ fPhiDistribution->SetParameter(1,dNewV2RP);
+ } else
+ {
+ if(fMinV2RP != fMaxV2RP)
+ {
+ dNewV2RP = fMyTRandom3->Uniform(fMinV2RP,fMaxV2RP);
+ fPhiDistribution->SetParameter(1,dNewV2RP);
+ } else
+ {
+ dNewV2RP = fMinV2RP;
+ fPhiDistribution->SetParameter(1,dNewV2RP);
+ }
+ }
+ }
+ // sampling the V4:
+ Double_t dNewV4RP = fV4RP;
+ if(fV4SpreadRP>0.0) dNewV4RP = fMyTRandom3->Gaus(fV4RP,fV4SpreadRP);
+ fPhiDistribution->SetParameter(3,dNewV4RP);
+
// eta:
Double_t dEtaMin = -1.; // to be improved
Double_t dEtaMax = 1.; // to be improved
-
-
+
Int_t iGoodTracks = 0;
Int_t iSelParticlesRP = 0;
Int_t iSelParticlesPOI = 0;
- for(Int_t i=0;i<fMultiplicityOfRP;i++) {
- AliFlowTrackSimple* pTrack = new AliFlowTrackSimple();
- pTrack->SetPt(fPtFormula->GetRandom());
- pTrack->SetEta(fMyTRandom3->Uniform(dEtaMin,dEtaMax));
- pTrack->SetPhi(fPhiFormula->GetRandom()+fMCReactionPlaneAngle);
- pTrack->SetForRPSelection(kTRUE);
- iSelParticlesRP++;
- pTrack->SetForPOISelection(kTRUE);
- iSelParticlesPOI++;
-
- pEvent->TrackCollection()->Add(pTrack);
- iGoodTracks++;
+ Double_t dTmpPt = 0.;
+ Double_t dTmpEta = 0.;
+ Double_t dTmpV2 = 0.;
+ Double_t dTmpPhi = 0.;
+ Bool_t bUniformAcceptance = kTRUE;
+ Double_t Pi = TMath::Pi();
+ if(!((fPhiMin1==0.) && (fPhiMax1==0.) && (fPhiMin2==0.) && (fPhiMax2==0.)))
+ {
+ bUniformAcceptance = kFALSE;
}
-
+ for(Int_t i=0;i<iNewMultiplicityOfRP;i++) {
+ dTmpPt = fPtSpectra->GetRandom();
+ dTmpEta = fMyTRandom3->Uniform(dEtaMin,dEtaMax);
+ // to be improved:
+ if(!fUseConstantHarmonics) {
+ if(dTmpPt >= fPtCutOff) {
+ dTmpV2 = fV2RPMax;
+ } else {
+ dTmpV2 = fV2RPMax*(dTmpPt/fPtCutOff);
+ }
+ fPhiDistribution->SetParameter(1,dTmpV2);
+ }
+ dTmpPhi = fPhiDistribution->GetRandom();
+ // add the track to the event
+ for(Int_t d=0;d<fNoOfLoops;d++) {
+ AliFlowTrackSimple* pTrack = new AliFlowTrackSimple();
+
+ // uniform acceptance:
+ if(bUniformAcceptance)
+ {
+ pTrack->SetPt(dTmpPt);
+ pTrack->SetEta(dTmpEta);
+ pTrack->SetPhi(dTmpPhi);
+ pTrack->SetForRPSelection(kTRUE);
+ iSelParticlesRP++;
+ // assign particles to subevents
+ if (pTrack->Eta()>=fEtaMinA && pTrack->Eta()<=fEtaMaxA) {
+ pTrack->SetForSubevent(0);
+ }
+ if (pTrack->Eta()>=fEtaMinB && pTrack->Eta()<=fEtaMaxB) {
+ pTrack->SetForSubevent(1);
+ }
+ pTrack->SetForPOISelection(kTRUE);
+ iSelParticlesPOI++;
+ pEvent->TrackCollection()->Add(pTrack);
+ iGoodTracks++;
+ }
+ // non-uniform acceptance, 1st sector:
+ else if ((dTmpPhi > fPhiMin1*Pi/180) && (dTmpPhi < fPhiMax1*Pi/180))
+ {
+ if(fMyTRandom3->Uniform(0,1) > 1 - fProbability1)
+ {
+ pTrack->SetPt(dTmpPt);
+ pTrack->SetEta(dTmpEta);
+ pTrack->SetPhi(dTmpPhi);
+ pTrack->SetForRPSelection(kTRUE);
+ iSelParticlesRP++;
+ // assign particles to subevents
+ if (pTrack->Eta()>=fEtaMinA && pTrack->Eta()<=fEtaMaxA) {
+ pTrack->SetForSubevent(0);
+ }
+ if (pTrack->Eta()>=fEtaMinB && pTrack->Eta()<=fEtaMaxB) {
+ pTrack->SetForSubevent(1);
+ }
+ pTrack->SetForPOISelection(kTRUE);
+ iSelParticlesPOI++;
+ pEvent->TrackCollection()->Add(pTrack);
+ iGoodTracks++;
+ }
+ }
+ // non-uniform acceptance, 2nd sector:
+ else if ((dTmpPhi > fPhiMin2*Pi/180) && (dTmpPhi < fPhiMax2*Pi/180))
+ {
+ if(fMyTRandom3->Uniform(0,1) > 1 - fProbability2)
+ {
+ pTrack->SetPt(dTmpPt);
+ pTrack->SetEta(dTmpEta);
+ pTrack->SetPhi(dTmpPhi);
+ pTrack->SetForRPSelection(kTRUE);
+ iSelParticlesRP++;
+ // assign particles to subevents
+ if (pTrack->Eta()>=fEtaMinA && pTrack->Eta()<=fEtaMaxA) {
+ pTrack->SetForSubevent(0);
+ }
+ if (pTrack->Eta()>=fEtaMinB && pTrack->Eta()<=fEtaMaxB) {
+ pTrack->SetForSubevent(1);
+ }
+
+ pTrack->SetForPOISelection(kTRUE);
+ iSelParticlesPOI++;
+ pEvent->TrackCollection()->Add(pTrack);
+ iGoodTracks++;
+ }
+ } else
+ {
+ pTrack->SetPt(dTmpPt);
+ pTrack->SetEta(dTmpEta);
+ pTrack->SetPhi(dTmpPhi);
+ pTrack->SetForRPSelection(kTRUE);
+ iSelParticlesRP++;
+ // assign particles to subevents
+ if (pTrack->Eta()>=fEtaMinA && pTrack->Eta()<=fEtaMaxA) {
+ pTrack->SetForSubevent(0);
+ }
+ if (pTrack->Eta()>=fEtaMinB && pTrack->Eta()<=fEtaMaxB) {
+ pTrack->SetForSubevent(1);
+ }
+ pTrack->SetForPOISelection(kTRUE);
+ iSelParticlesPOI++;
+ pEvent->TrackCollection()->Add(pTrack);
+ iGoodTracks++;
+ }
+ }
+ }
+
+ // update the event quantities
pEvent->SetEventNSelTracksRP(iSelParticlesRP);
pEvent->SetNumberOfTracks(iGoodTracks);//tracks used either for RP or for POI selection
- pEvent->SetMCReactionPlaneAngle(fMCReactionPlaneAngle);
-
- if (!fMCReactionPlaneAngle == 0) cout<<" MC Reaction Plane Angle = "<< fMCReactionPlaneAngle << endl;
- else cout<<" MC Reaction Plane Angle = unknown "<< endl;
+ pEvent->SetMCReactionPlaneAngle(dMCReactionPlaneAngle);
- cout<<" iGoodTracks = "<< iGoodTracks << endl;
- cout<<" # of RP selected tracks = "<<iSelParticlesRP<<endl;
- cout<<" # of POI selected tracks = "<<iSelParticlesPOI<<endl;
- cout << "# " << ++fCount << " events processed" << endl;
+ if ( (++fCount % 100) == 0) {
+ if (!dMCReactionPlaneAngle == 0) cout<<" MC Reaction Plane Angle = "<< dMCReactionPlaneAngle << endl;
+ else cout<<" MC Reaction Plane Angle = unknown "<< endl;
+ cout<<" iGoodTracks = "<< iGoodTracks << endl;
+ cout<<" # of RP selected tracks = "<<iSelParticlesRP<<endl;
+ cout<<" # of POI selected tracks = "<<iSelParticlesPOI<<endl;
+ cout << "# " << fCount << " events processed" << endl;
+ }
- return pEvent;
+ return pEvent;
} // end of CreateEventOnTheFly()