few extra options for momentum dependence v_n
[u/mrichter/AliRoot.git] / PWG2 / FLOW / AliFlowCommon / AliFlowEventSimpleMakerOnTheFly.cxx
CommitLineData
93ff27bd 1/*************************************************************************
2* Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
3* *
4* Author: The ALICE Off-line Project. *
5* Contributors are mentioned in the code where appropriate. *
6* *
7* Permission to use, copy, modify and distribute this software and its *
8* documentation strictly for non-commercial purposes is hereby granted *
9* without fee, provided that the above copyright notice appears in all *
10* copies and that both the copyright notice and this permission notice *
11* appear in the supporting documentation. The authors make no claims *
12* about the suitability of this software for any purpose. It is *
13* provided "as is" without express or implied warranty. *
14**************************************************************************/
15
16/**********************************
17 * create an event and perform *
18 * flow analysis 'on the fly' *
19 * *
20 * authors: Raimond Snellings *
21 * (snelling@nikhef.nl) *
22 * Ante Bilandzic *
23 * (anteb@nikhef.nl) *
24 *********************************/
25
26#include "Riostream.h"
27#include "TMath.h"
28#include "TF1.h"
29#include "TRandom3.h"
4d603348 30//#include "TUnuran.h"
93ff27bd 31
32#include "AliFlowEventSimpleMakerOnTheFly.h"
33#include "AliFlowEventSimple.h"
34#include "AliFlowTrackSimple.h"
35
36ClassImp(AliFlowEventSimpleMakerOnTheFly)
37
38
39//========================================================================
40
41
5b8dccde 42AliFlowEventSimpleMakerOnTheFly::AliFlowEventSimpleMakerOnTheFly(UInt_t iseed):
43 fMultiplicityOfRP(0),
44 fMultiplicitySpreadOfRP(0.),
4d603348 45 fUseConstantHarmonics(kFALSE),
e8c7e66c 46 fV1RP(0.),
47 fV1SpreadRP(0.),
48 fV2RP(0.),
49 fV2SpreadRP(0.),
4d603348 50 fV4RP(0.),
51 fV4SpreadRP(0.),
52 fV2RPMax(0.),
53 fPtCutOff(0.),
e8c7e66c 54 fPtSpectra(NULL),
55 fPhiDistribution(NULL),
5b8dccde 56 fMyTRandom3(NULL),
aaa62cf0 57 fCount(0),
58 fNoOfLoops(1)
93ff27bd 59 {
60 // constructor
4d603348 61 fMyTRandom3 = new TRandom3(iseed);
62 gRandom->SetSeed(fMyTRandom3->Integer(65539));
63
64 //fMyUnuran = new TUnuran();
65 //fMyUnuran->SetSeed(iseed);
93ff27bd 66 }
67
68
69//========================================================================
70
71
72AliFlowEventSimpleMakerOnTheFly::~AliFlowEventSimpleMakerOnTheFly()
73{
74 // destructor
4d603348 75 if (fPtSpectra) delete fPtSpectra;
76 if (fPhiDistribution) delete fPhiDistribution;
77 if (fMyTRandom3) delete fMyTRandom3;
93ff27bd 78}
79
80
81//========================================================================
82
83
4d603348 84void AliFlowEventSimpleMakerOnTheFly::Init()
85{
86 // define the pt spectra and phi distribution
87
88 // pt spectra:
89 Double_t dPtMin = 0.; // to be improved (move this to the body of contstructor?)
90 Double_t dPtMax = 10.; // to be improved (move this to the body of contstructor?)
91
92 fPtSpectra = new TF1("fPtSpectra","[0]*x*TMath::Exp(-x*x)",dPtMin,dPtMax);
93 fPtSpectra->SetParName(0,"Multiplicity of RPs");
94
95 // phi distribution:
96 Double_t dPhiMin = 0.; // to be improved (move this to the body of contstructor?)
97 Double_t dPhiMax = TMath::TwoPi(); // to be improved (move this to the body of contstructor?)
98
99 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);
100 fPhiDistribution->SetParName(0,"directed flow");
101 fPhiDistribution->SetParName(1,"elliptic flow");
102 fPhiDistribution->SetParName(2,"Reaction Plane");
103 fPhiDistribution->SetParName(3,"harmonic 4"); // to be improved (name)
104}
105
106//========================================================================
107
93ff27bd 108AliFlowEventSimple* AliFlowEventSimpleMakerOnTheFly::CreateEventOnTheFly()
109{
110 // method to create event on the fly
111
112 AliFlowEventSimple* pEvent = new AliFlowEventSimple(fMultiplicityOfRP);
4d603348 113
e8c7e66c 114 // sampling the multiplicity:
4d603348 115 Int_t iNewMultiplicityOfRP = fMultiplicityOfRP;
116 if(fMultiplicitySpreadOfRP>0.0) iNewMultiplicityOfRP = (Int_t)fMyTRandom3->Gaus(fMultiplicityOfRP,fMultiplicitySpreadOfRP);
117 fPtSpectra->SetParameter(0,iNewMultiplicityOfRP);
93ff27bd 118
4d603348 119 // sampling the reaction plane
120 Double_t dMCReactionPlaneAngle = fMyTRandom3->Uniform(0.,TMath::TwoPi());
121 fPhiDistribution->SetParameter(2,dMCReactionPlaneAngle);
e8c7e66c 122
123 // sampling the V1:
4d603348 124 Double_t dNewV1RP=fV1RP;
125 if(fV1SpreadRP>0.0) {dNewV1RP = fMyTRandom3->Gaus(fV1RP,fV1SpreadRP);}
126 fPhiDistribution->SetParameter(0,dNewV1RP);
93ff27bd 127
e8c7e66c 128 // sampling the V2:
4d603348 129 if(fUseConstantHarmonics)
130 {
131 Double_t dNewV2RP = fV2RP;
132 if(fV2SpreadRP>0.0) dNewV2RP = fMyTRandom3->Gaus(fV2RP,fV2SpreadRP);
133 fPhiDistribution->SetParameter(1,dNewV2RP);
134 }
135
136 // sampling the V4:
137 Double_t dNewV4RP = fV4RP;
138 if(fV4SpreadRP>0.0) dNewV4RP = fMyTRandom3->Gaus(fV4RP,fV4SpreadRP);
139 fPhiDistribution->SetParameter(3,dNewV4RP);
e8c7e66c 140
93ff27bd 141 // eta:
142 Double_t dEtaMin = -1.; // to be improved
143 Double_t dEtaMax = 1.; // to be improved
e8c7e66c 144
93ff27bd 145 Int_t iGoodTracks = 0;
146 Int_t iSelParticlesRP = 0;
147 Int_t iSelParticlesPOI = 0;
4d603348 148 Double_t dTmpPt = 0.;
149 Double_t dTmpEta = 0.;
150 Double_t dTmpV2 = 0.;
151 Double_t dTmpPhi = 0.;
152 for(Int_t i=0;i<iNewMultiplicityOfRP;i++) {
153 dTmpPt = fPtSpectra->GetRandom();
154 dTmpEta = fMyTRandom3->Uniform(dEtaMin,dEtaMax);
155 // to be improved:
156 if(!fUseConstantHarmonics) {
157 if(dTmpPt >= fPtCutOff) {
158 dTmpV2 = fV2RPMax;
159 } else {
160 dTmpV2 = fV2RPMax*(dTmpPt/fPtCutOff);
161 }
162 fPhiDistribution->SetParameter(1,dTmpV2);
163 }
164 dTmpPhi = fPhiDistribution->GetRandom();
165 // add the track to the event
aaa62cf0 166 for(Int_t d=0;d<fNoOfLoops;d++) {
167 AliFlowTrackSimple* pTrack = new AliFlowTrackSimple();
4d603348 168 pTrack->SetPt(dTmpPt);
169 pTrack->SetEta(dTmpEta);
170 pTrack->SetPhi(dTmpPhi);
aaa62cf0 171 pTrack->SetForRPSelection(kTRUE);
172 iSelParticlesRP++;
173 pTrack->SetForPOISelection(kTRUE);
174 iSelParticlesPOI++;
175 pEvent->TrackCollection()->Add(pTrack);
176 iGoodTracks++;
177 }
93ff27bd 178 }
4d603348 179
180 // update the event quantities
93ff27bd 181 pEvent->SetEventNSelTracksRP(iSelParticlesRP);
182 pEvent->SetNumberOfTracks(iGoodTracks);//tracks used either for RP or for POI selection
4d603348 183 pEvent->SetMCReactionPlaneAngle(dMCReactionPlaneAngle);
93ff27bd 184
4d603348 185 if (!dMCReactionPlaneAngle == 0) cout<<" MC Reaction Plane Angle = "<< dMCReactionPlaneAngle << endl;
5b8dccde 186 else cout<<" MC Reaction Plane Angle = unknown "<< endl;
187
188 cout<<" iGoodTracks = "<< iGoodTracks << endl;
189 cout<<" # of RP selected tracks = "<<iSelParticlesRP<<endl;
190 cout<<" # of POI selected tracks = "<<iSelParticlesPOI<<endl;
191 cout << "# " << ++fCount << " events processed" << endl;
4d603348 192
28ece29b 193 return pEvent;
93ff27bd 194
195} // end of CreateEventOnTheFly()
196
197
198