hooks for PMD flow analysis
[u/mrichter/AliRoot.git] / PWG2 / FLOW / macros / runFlowAnalysisOnTheFly.C
CommitLineData
93ff27bd 1#include "TStopwatch.h"
2#include "TObjArray"
3#include "Riostream.h"
4#include "TFile.h"
5
6//--------------------------------------------------------------------------------------
7// RUN SETTINGS
8// flow analysis method can be: (set to kTRUE or kFALSE)
06868ea2 9Bool_t MCEP = kTRUE;
c741f5d0 10Bool_t SP = kTRUE;
06868ea2 11Bool_t GFC = kTRUE;
12Bool_t QC = kTRUE;
13Bool_t FQD = kTRUE;
14Bool_t LYZ1SUM = kTRUE;
15Bool_t LYZ1PROD = kTRUE;
c741f5d0 16Bool_t LYZ2SUM = kFALSE;
17Bool_t LYZ2PROD = kFALSE;
5b40431d 18Bool_t LYZEP = kFALSE;
65f19276 19Bool_t MH = kFALSE; // mixed harmonics
83bc3e95 20Bool_t NL = kFALSE; // nested loops
d66d46f7 21Bool_t MCEP_AH = kFALSE; // MCEP in another harmonic
93ff27bd 22//--------------------------------------------------------------------------------------
23
7fb984f6 24// Weights
25// use weights for Q-vector:
26Bool_t usePhiWeights = kFALSE; // phi weights (correction for non-uniform azimuthal acceptance)
27Bool_t usePtWeights = kFALSE; // pt weights
28Bool_t useEtaWeights = kFALSE; // eta weights
29
7a01f4a7 30// Define the range for eta subevents
31Double_t minA = -0.9;
32Double_t maxA = -0.01;
33Double_t minB = 0.01;
34Double_t maxB = 0.9;
35
e1911c19 36// Define simple cuts for RP selection:
37Double_t ptMinRP = 0.0; // in GeV
38Double_t ptMaxRP = 10.0; // in GeV
d66d46f7 39Double_t etaMinRP = -1.;
40Double_t etaMaxRP = 1.;
e1911c19 41Double_t phiMinRP = 0.0; // in degrees
42Double_t phiMaxRP = 360.0; // in degrees
43Int_t pidRP = -1; // to be improved (supported eventually)
44
45// Define simple cuts for POI selection:
46Double_t ptMinPOI = 0.0; // in GeV
47Double_t ptMaxPOI = 10.0; // in GeV
d66d46f7 48Double_t etaMinPOI = -1.;
49Double_t etaMaxPOI = 1.;
e1911c19 50Double_t phiMinPOI = 0.0; // in degrees
51Double_t phiMaxPOI = 360.0; // in degrees
52Int_t pidPOI = -1; // to be improved (supported eventually)
53
7fb984f6 54// Parameters for the simulation of events 'on the fly':
7a01f4a7 55//===SEED========================================================
56// use always the same seed for random generators.
57// usage of same seed (kTRUE) is relevant in two cases:
58// 1.) If you want to use LYZ method to calcualte differential flow;
59// 2.) If you want to use phi weights for GFC, QC and FQD
2a54d329 60Bool_t bSameSeed = kFALSE;
4d603348 61
7a01f4a7 62//===NONFLOW=============================================
63// number of times to use each track (to simulate nonflow)
64Int_t iLoops = 1;
524798bf 65Double_t phiRange = 0.0; // If the original track with azimuthal angle phi is splitted,
66 // the azimuthal angle of splitted track is sampled uniformly
67 // from (phi-phiRange, phi+phiRange). If phiRange = 0.0, the
68 // azimuthal angle of splitted track is the same as azimuthal
69 // angle of original track. (Remark: phiRange's unit is degree.)
2a54d329 70Double_t ptRange = 0.0; // If the original track with momentum pt is splitted,
71 // the momentum of splitted track is sampled uniformly
72 // from (pt-ptRange, pt+ptRange). If ptRange = 0.0, the
73 // momentum of splitted track is the same as momentum
524798bf 74 // of original track. (Remark: ptRange's unit is GeV.)
2a54d329 75Double_t etaRange = 0.0; // If the original track with pseudorapidity eta is splitted,
524798bf 76 // the pseudorapidity of splitted track is sampled uniformly
77 // from (eta-etaRange, eta+etaRange). If etaRange = 0.0, the
78 // pseudorapidity of splitted track will be the same as the
79 // pseudorapidity of original track.
e1911c19 80// in addition one can simulate nonflow only in a certain detector's sector
81// ranging from nonflowSectorMin to nonflowSectorMax. Outside this sector the
82// tracks will not be splitted. Use the following two settings only if iLoops>1:
83Double_t nonflowSectorMin = 0.0; // detector's sector in which tracks will be splitted starts at this angle (in degrees)
1f89f68d 84Double_t nonflowSectorMax = 360.0; // detector's sector in which tracks will be splitted ends at this angle (in degrees)
85// to study nonflow fluctuations event-by-event one can create the jet-like structure in the following way:
86Bool_t bCreateJets = kFALSE; // set kTRUE if you want to study nonflow fluctuations via jet creation
87Double_t jetProbability = 0.2; // probability than event will have a jet
88Double_t jetTracksFraction = 0.1; // percentage of total tracks in an event belonging to the jet
89Double_t jetCone = 20; // angular size of jet cone (in degrees)
93ff27bd 90
598d292e 91//===GLAUBER MODEL================================================================
92Bool_t bUseGlauberModel = kFALSE; // 1.) if kTRUE = multiplicity and constant flow harmonics are
93 // determined event-by-event from Glauber model
94 // (pt and/or eta dependence of flow harmonics not supported;
95 // settings for multiplicity and flow harmonics bellow are irrelevant)
96 // 2.) if kFALSE = multiplicity and flow harmonics are determined
97 // independently event-by-event with bellow settings
98 // (pt and/or eta dependence of flow harmonics supported)
99
7a01f4a7 100//===FLOW HARMONICS===============================================================
bdb5c432 101// harmonics V1, V2, V4... are constants or functions of pt and eta:
102Bool_t bPtDependentHarmonicV1 = kFALSE;
103Bool_t bEtaDependentHarmonicV1 = kFALSE;
104Bool_t bPtDependentHarmonicV2 = kFALSE;
105Bool_t bEtaDependentHarmonicV2 = kFALSE;
106Bool_t bPtDependentHarmonicV4 = kFALSE;
107Bool_t bEtaDependentHarmonicV4 = kFALSE;
108// 1.) if you use constant harmonics (bPtDependentHarmonic* = kFALSE, bEtaDependentHarmonic* = kFALSE)
109// you can additionally select if V2 will be sampled e-b-e from Gaussian or from uniform distribution
110// (constant harmonics V1 and V4 will be always sampled from Gaussian):
524798bf 111Bool_t bConstantV2IsSampledFromGauss = kTRUE;
112 // 1a.) if kTRUE = elliptic flow of RPs is sampled e-b-e from Gaussian distribution with
113 // mean = dV2RP and spread = dV2SpreadRP (analogously for V1 and V4).
114 // 1b.) if kFALSE = elliptic flow of RPs is sampled e-b-e uniformly from
115 // interval [dMinV2RP,dMaxV2RP] (not supported for V1 and V4).
116 // 1c.) for a fixed elliptic flow e-b-e use Gaussian with zero spread or use uniform with dMinV2RP=dMaxV2RP
117 // V2:
118 Double_t dV2RP = 0.05; // elliptic flow of RPs (if sampled from Gaussian)
119 Double_t dV2SpreadRP = 0.0; // elliptic flow spread of RPs (if sampled from Gaussian)
d66d46f7 120 Double_t dMinV2RP = 0.04; // minimal elliptic flow of RPs (if sampled uniformly)
121 Double_t dMaxV2RP = 0.06; // maximal elliptic flow of RPs (if sampled uniformly)
524798bf 122 // V1:
123 Double_t dV1RP = 0.0; // directed flow of RPs
124 Double_t dV1SpreadRP = 0.0; // directed flow spread of RPs
65f19276 125 // V3:
126 Double_t dV3RP = 0.0; // triangular flow of RPs
127 Double_t dV3SpreadRP = 0.0; // triangular flow spread of RP
524798bf 128 // V4:
129 Double_t dV4RP = 0.0; // harmonic V4 of RPs (to be improved: name needed)
130 Double_t dV4SpreadRP = 0.0; // harmonic V4's spread of RPs (to be improved: name needed)
bdb5c432 131// 2.) if you use (pt,eta) dependent harmonic V1 (bPtDependentHarmonicV1 = kTRUE or bEtaDependentHarmonicV1 = kTRUE):
132// 2a.) V1(pt) is linear up to pt = dV1PtCutOff and for pt > dV1PtCutOff it is constant, V1(pt) = dV1vsPtEtaMax:
133// 2b.) V1(eta) is determined from formula V1(eta) = -eta
134 Double_t dV1vsPtEtaMax = 0.10; // maximum value of V1(pt,eta) (for pt >= dV1PtCutOff and at eta = 0)
135 Double_t dV1PtCutOff = 2.0; // up to this pt value V1(pt) is growing linearly versus pt
136// 3.) if you use (pt,eta) dependent harmonic V2 (bPtDependentHarmonicV2 = kTRUE or bEtaDependentHarmonicV2 = kTRUE):
137// 2a.) V2(pt) is linear up to pt = dV2PtCutOff and for pt > dV2PtCutOff it is constant, V2(pt) = dV2vsPtEtaMax:
524798bf 138// 2b.) V2(eta) is Gaussian centered at midrapidity (eta=0) with V2(0) = dV2vsPtEtaMax and spread = dV2vsEtaSpread:
bdb5c432 139 Double_t dV2vsPtEtaMax = 0.20; // maximum value of V2(pt,eta) (for pt >= dV2PtCutOff and at eta = 0)
140 Double_t dV2PtCutOff = 2.0; // up to this pt value V2(pt) is growing linearly versus pt
524798bf 141 Double_t dV2vsEtaSpread = 0.75; // V2(eta) is Gaussian centered at midrapidity (eta=0) with spread = dV2vsEtaSpread
bdb5c432 142// 4.) if you use (pt,eta) dependent harmonic V4 (bPtDependentHarmonicV4 = kTRUE or bEtaDependentHarmonicV4 = kTRUE):
143// V4(pt,eta) is determined as V4(pt,eta) = pow(V2(pt,eta),2.)
7a01f4a7 144
145//===MULTIPLICITY===============================================================
146// 1.) if kTRUE = multiplicitiy of RPs is sampled e-b-e from Gaussian distribution with
147// mean = iMultiplicityOfRP and spread = dMultiplicitySpreadOfRP
148// 2.) if kFALSE = multiplicitiy of RPs is sampled e-b-e uniformly from
149// interval [iMinMultOfRP,iMaxMultOfRP]
150// 3.) for a fixed multiplicity use Gaussian with zero spread or use uniform with iMinMult=iMaxMult
151Bool_t bMultDistrOfRPsIsGauss = kTRUE;
152
153Int_t iMultiplicityOfRP = 500; // mean multiplicity of RPs (if sampled from Gaussian)
154Double_t dMultiplicitySpreadOfRP = 0; // multiplicity spread of RPs (if sampled from Gaussian)
598d292e 155Int_t iMinMultOfRP = 400; // minimal multiplicity of RPs(if sampled uniformly)
d66d46f7 156Int_t iMaxMultOfRP = 600; // maximal multiplicity of RPs (if sampled uniformly)
7a01f4a7 157
158//===DETECTOR ACCEPTANCE===============================================================
159
160// 1.) if kTRUE = detectors has uniform azimuthal acceptance
161// 2.) if kFALSE = you will simulate detector with non-uniform acceptance in one or two sectors.
162// For each of two sectors you specify phi_min, phi_max and probability p. Then all particles
163// going in direction phi_min < phi < phi_max will be taken with probability p. If p = 0, that
164// sector is blocked. Set bellow phimin1, phimax1, p1 for the first sector and phimin2, phimax2, p2
165// for the second sector. If you set phimin2 = phimax2 = p2 = 0, only first non-uniform sector is
166// simulated.
167Bool_t uniformAcceptance = kTRUE;
168
e23cbd39 169// settings for non-uniform acceptance:
170// Remark: set the angles in degrees from interval [0,360] and probability from interval [0,1]
171
172// 1st non-uniform sector:
173Double_t phimin1 = 60; // first non-uniform sector starts at this azimuth
174Double_t phimax1 = 120; // first non-uniform sector ends at this azimuth
e1911c19 175Double_t p1 = 0.5; // e.g. if p1 = 0 all particles emitted in phimin1 < phi < phimax1 are blocked
e23cbd39 176 // e.g. if p1 = 0.5 half of the particles emitted in phimin1 < phi < phimax1 are blocked
177
178// 2nd non-uniform sector (Remark: if you do NOT want to simulate this sector, set phimin2 = phimax2 = p2 = 0):
d66d46f7 179Double_t phimin2 = 0.0; // second non-uniform sector starts at this azimuth (make sure phimin2 > phimax1 !!!!)
180Double_t phimax2 = 0.0; // second non-uniform sector ends at this azimuth
524798bf 181Double_t p2 = 0.0;
7a01f4a7 182
524798bf 183
184//===MODIFYING Pt SPECTRA===============================================================
185Double_t dTemperatureOfRP = 0.44; // 'temperature' of RPs in GeV/c (increase this parameter to get more high pt RPs)
e23cbd39 186
4d603348 187enum anaModes {mLocal,mLocalSource,mLocalPAR};
93ff27bd 188// mLocal: Analyze data on your computer using aliroot
189// mLocalPAR: Analyze data on your computer using root + PAR files
190// mLocalSource: Analyze data on your computer using root + source files
191
524798bf 192int runFlowAnalysisOnTheFly(Int_t nEvts=1000, Int_t mode=mLocal)
93ff27bd 193{
bdb5c432 194 CheckUserSettings();
93ff27bd 195 TStopwatch timer;
196 timer.Start();
197
93ff27bd 198 cout<<endl;
199 cout<<endl;
200 cout<<" ---- ARE YOU READY TO FLY ? ---- "<<endl;
201 cout<<endl;
202
203 cout<<endl;
204 cout<<" ---- BEGIN FLOW ANALYSIS 'ON THE FLY' ---- "<<endl;
205 cout<<endl;
206 cout<<endl;
207
93ff27bd 208 LoadLibraries(mode);
209
4d603348 210 // Initialize the seed for random generator
211 UInt_t sseed = 0;
212
213 if(bSameSeed)
214 {
215 sseed = 44; // the default constant value for seed for random generators
216 }
5b8dccde 217
4d603348 218 if(!bSameSeed)
219 {
220 TTimeStamp dt;
221 sseed = dt.GetNanoSec()/1000;
222 }
223
e23cbd39 224 cout<<endl;
225 cout<<"Seed for the random generators is "<<sseed<<endl;
226 cout<<endl;
227
7fb984f6 228 //---------------------------------------------------------------------------------------
229 // If the weights are used:
230 TFile *fileWithWeights = NULL;
231 TList *listWithWeights = NULL;
232
233 if(usePhiWeights||usePtWeights||useEtaWeights) {
234 fileWithWeights = TFile::Open("weights.root","READ");
235 if(fileWithWeights) {
236 listWithWeights = (TList*)fileWithWeights->Get("weights");
237 }
238 else
239 {cout << " WARNING: the file <weights.root> with weights from the previous run was not found."<<endl;
240 break;
241 }
242 }
243
93ff27bd 244 //---------------------------------------------------------------------------------------
245 // Initialize the flowevent maker
5b8dccde 246 AliFlowEventSimpleMakerOnTheFly* eventMakerOnTheFly = new AliFlowEventSimpleMakerOnTheFly(sseed);
4d603348 247 eventMakerOnTheFly->Init();
7a01f4a7 248 //set the range for eta subevents
249 eventMakerOnTheFly->SetSubeventEtaRange(minA,maxA,minB,maxB);
93ff27bd 250
93ff27bd 251 //---------------------------------------------------------------------------------------
ce4a88f5 252 // Initialize all the flow methods for default analysis:
253 AliFlowAnalysisWithQCumulants *qc = NULL;
254 AliFlowAnalysisWithCumulants *gfc = NULL;
255 AliFlowAnalysisWithFittingQDistribution *fqd = NULL;
256 AliFlowAnalysisWithLeeYangZeros *lyz1sum = NULL;
257 AliFlowAnalysisWithLeeYangZeros *lyz1prod = NULL;
258 AliFlowAnalysisWithLeeYangZeros *lyz2sum = NULL;
259 AliFlowAnalysisWithLeeYangZeros *lyz2prod = NULL;
260 AliFlowAnalysisWithLYZEventPlane *lyzep = NULL;
261 AliFlowAnalysisWithScalarProduct *sp = NULL;
d66d46f7 262 AliFlowAnalysisWithMixedHarmonics *mh = NULL;
83bc3e95 263 AliFlowAnalysisWithNestedLoops *nl = NULL;
ce4a88f5 264 AliFlowAnalysisWithMCEventPlane *mcep = NULL;
b80d67e1 265 AliFlowAnalysisWithMCEventPlane *mcep_ah = NULL;
93ff27bd 266
267 // MCEP = monte carlo event plane
268 if (MCEP) {
269 AliFlowAnalysisWithMCEventPlane *mcep = new AliFlowAnalysisWithMCEventPlane();
bdb5c432 270 //mcep->SetHarmonic(2); // default is v2
b3a29c75 271 //mcep->SetFlowOfResonances(kTRUE);
93ff27bd 272 mcep->Init();
273 }
274
d66d46f7 275 // MCEP = monte carlo event plane in another harmonic:
276 if(MCEP_AH)
277 {
278 AliFlowAnalysisWithMCEventPlane *mcep_ah = new AliFlowAnalysisWithMCEventPlane();
279 mcep_ah->SetHarmonic(1);
280 mcep_ah->Init();
281 }
282
283 // Mixed harmonics:
284 if(MH)
285 {
286 AliFlowAnalysisWithMixedHarmonics *mh = new AliFlowAnalysisWithMixedHarmonics();
7e570884 287 mh->SetHarmonic(1); // integer n in expression cos[n(2phi1-phi2-phi3)] = v2n*vn^2
83bc3e95 288 mh->SetMinMultiplicity(100);
289 mh->SetNoOfMultipicityBins(5);
290 mh->SetMultipicityBinWidth(200);
d66d46f7 291 mh->Init();
292 }
293
83bc3e95 294 // NL = nested loops:
295 if(NL) {
296 AliFlowAnalysisWithNestedLoops *nl = new AliFlowAnalysisWithNestedLoops();
297 nl->Init();
298 }
299
ce4a88f5 300 // QC = Q-cumulants
93ff27bd 301 if(QC) {
302 AliFlowAnalysisWithQCumulants* qc = new AliFlowAnalysisWithQCumulants();
7fb984f6 303 if(listWithWeights) qc->SetWeightsList(listWithWeights);
304 if(usePhiWeights) qc->SetUsePhiWeights(usePhiWeights);
305 if(usePtWeights) qc->SetUsePtWeights(usePtWeights);
306 if(useEtaWeights) qc->SetUseEtaWeights(useEtaWeights);
afa4af05 307 // qc->SetHarmonic(2); // default is v2
dd442cd2 308 qc->SetApplyCorrectionForNUA(kFALSE);
309 qc->SetFillMultipleControlHistograms(kFALSE);
a5b7efd0 310 // qc->SetCalculate2DFlow(kFALSE); // default
311 // qc->SetMultiplicityWeight("combinations"); // default
312 // qc->SetMultiplicityWeight("unit");
313 // qc->SetMultiplicityWeight("multiplicity");
dd442cd2 314 qc->SetCalculateCumulantsVsM(kFALSE);
3688e459 315 qc->SetnBinsMult(10000);
316 qc->SetMinMult(0);
317 qc->SetMaxMult(10000);
ce4a88f5 318 qc->Init();
93ff27bd 319 }
320
321 // GFC = Generating Function Cumulants
322 if(GFC) {
323 AliFlowAnalysisWithCumulants* gfc = new AliFlowAnalysisWithCumulants();
7fb984f6 324 if(listWithWeights) gfc->SetWeightsList(listWithWeights);
325 if(usePhiWeights) gfc->SetUsePhiWeights(usePhiWeights);
326 if(usePtWeights) gfc->SetUsePtWeights(usePtWeights);
327 if(useEtaWeights) gfc->SetUseEtaWeights(useEtaWeights);
3688e459 328 // calculation vs multiplicity:
329 gfc->SetCalculateVsMultiplicity(kFALSE);
330 gfc->SetnBinsMult(10000);
331 gfc->SetMinMult(0);
332 gfc->SetMaxMult(10000);
333 // tuning of interpolating parameters:
334 gfc->SetTuneParameters(kFALSE);
335 Double_t r0[10] = {1.8,1.9,2.0,2.1,2.2,2.3,2.4,2.5,2.6,2.7}; // up to 10 values allowed
336 for(Int_t r=0;r<10;r++) {gfc->SetTuningR0(r0[r],r);}
ce4a88f5 337 gfc->Init();
93ff27bd 338 }
339
340 // FQD = Fitting q-distribution
341 if(FQD) {
ce4a88f5 342 AliFlowAnalysisWithFittingQDistribution* fqd = new AliFlowAnalysisWithFittingQDistribution();
7fb984f6 343 if(listWithWeights) fqd->SetWeightsList(listWithWeights);
344 if(usePhiWeights) fqd->SetUsePhiWeights(usePhiWeights);
ce4a88f5 345 fqd->Init();
93ff27bd 346 }
347
348 // SP = Scalar Product
349 if(SP) {
350 AliFlowAnalysisWithScalarProduct* sp = new AliFlowAnalysisWithScalarProduct();
7a01f4a7 351 if(listWithWeights) sp->SetWeightsList(listWithWeights);
352 sp->SetUsePhiWeights(usePhiWeights);
93ff27bd 353 sp->Init();
354 }
355
356 // LYZ1 = Lee-Yang Zeroes first run
c741f5d0 357 if(LYZ1SUM) {
358 AliFlowAnalysisWithLeeYangZeros* lyz1sum = new AliFlowAnalysisWithLeeYangZeros();
359 lyz1sum->SetFirstRun(kTRUE);
360 lyz1sum->SetUseSum(kTRUE);
361 lyz1sum->Init();
362 }
363 if(LYZ1PROD) {
364 AliFlowAnalysisWithLeeYangZeros* lyz1prod = new AliFlowAnalysisWithLeeYangZeros();
365 lyz1prod->SetFirstRun(kTRUE);
366 lyz1prod->SetUseSum(kFALSE);
367 lyz1prod->Init();
93ff27bd 368 }
93ff27bd 369 // LYZ2 = Lee-Yang Zeroes second run
c741f5d0 370 if(LYZ2SUM) {
371 AliFlowAnalysisWithLeeYangZeros* lyz2sum = new AliFlowAnalysisWithLeeYangZeros();
93ff27bd 372 // read the input file from the first run
c741f5d0 373 TString inputFileNameLYZ2SUM = "outputLYZ1SUManalysis.root" ;
374 TFile* inputFileLYZ2SUM = new TFile(inputFileNameLYZ2SUM.Data(),"READ");
375 if(!inputFileLYZ2SUM || inputFileLYZ2SUM->IsZombie()) {
5b40431d 376 cerr << " ERROR: To run LYZ2SUM you need the output file from LYZ1SUM. This file is not there! Please run LYZ1SUM first." << endl ;
93ff27bd 377 break;
378 }
379 else {
c741f5d0 380 TList* inputListLYZ2SUM = (TList*)inputFileLYZ2SUM->Get("cobjLYZ1SUM");
381 if (!inputListLYZ2SUM) {cout<<"Input list LYZ2SUM is NULL pointer!"<<endl; break;}
93ff27bd 382 else {
c741f5d0 383 cout<<"LYZ2SUM input file/list read..."<<endl;
384 lyz2sum->SetFirstRunList(inputListLYZ2SUM);
385 lyz2sum->SetFirstRun(kFALSE);
386 lyz2sum->SetUseSum(kTRUE);
387 lyz2sum->Init();
93ff27bd 388 }
389 }
390 }
c741f5d0 391 if(LYZ2PROD) {
392 AliFlowAnalysisWithLeeYangZeros* lyz2prod = new AliFlowAnalysisWithLeeYangZeros();
393 // read the input file from the first run
394 TString inputFileNameLYZ2PROD = "outputLYZ1PRODanalysis.root" ;
395 TFile* inputFileLYZ2PROD = new TFile(inputFileNameLYZ2PROD.Data(),"READ");
396 if(!inputFileLYZ2PROD || inputFileLYZ2PROD->IsZombie()) {
5b40431d 397 cerr << " ERROR: To run LYZ2PROD you need the output file from LYZ1PROD. This file is not there! Please run LYZ1PROD first." << endl ;
c741f5d0 398 break;
399 }
400 else {
401 TList* inputListLYZ2PROD = (TList*)inputFileLYZ2PROD->Get("cobjLYZ1PROD");
402 if (!inputListLYZ2PROD) {cout<<"Input list LYZ2PROD is NULL pointer!"<<endl; break;}
403 else {
404 cout<<"LYZ2PROD input file/list read..."<<endl;
405 lyz2prod->SetFirstRunList(inputListLYZ2PROD);
406 lyz2prod->SetFirstRun(kFALSE);
407 lyz2prod->SetUseSum(kFALSE);
408 lyz2prod->Init();
409 }
410 }
411 }
412
93ff27bd 413 // LYZEP = Lee-Yang Zeroes event plane
414 if(LYZEP) {
415 AliFlowLYZEventPlane* ep = new AliFlowLYZEventPlane() ;
416 AliFlowAnalysisWithLYZEventPlane* lyzep = new AliFlowAnalysisWithLYZEventPlane();
417 // read the input file from the second lyz run
c741f5d0 418 TString inputFileNameLYZEP = "outputLYZ2SUManalysis.root" ;
93ff27bd 419 TFile* inputFileLYZEP = new TFile(inputFileNameLYZEP.Data(),"READ");
420 if(!inputFileLYZEP || inputFileLYZEP->IsZombie()) {
5b40431d 421 cerr << " ERROR: To run LYZEP you need the output file from LYZ2SUM. This file is not there! Please run LYZ2SUM first." << endl ;
93ff27bd 422 break;
423 }
424 else {
c741f5d0 425 TList* inputListLYZEP = (TList*)inputFileLYZEP->Get("cobjLYZ2SUM");
426 if (!inputListLYZEP) {cout<<"Input list LYZEP is NULL pointer!"<<endl; break;}
93ff27bd 427 else {
428 cout<<"LYZEP input file/list read..."<<endl;
429 ep ->SetSecondRunList(inputListLYZEP);
430 lyzep->SetSecondRunList(inputListLYZEP);
431 ep ->Init();
432 lyzep->Init();
433 }
434 }
435 }
436 //---------------------------------------------------------------------------------------
ce4a88f5 437
524798bf 438 // set the global event parameters:
aaa62cf0 439 eventMakerOnTheFly->SetNoOfLoops(iLoops);
524798bf 440 eventMakerOnTheFly->SetPhiRange(phiRange*TMath::Pi()/180.);
2a54d329 441 eventMakerOnTheFly->SetPtRange(ptRange);
442 eventMakerOnTheFly->SetEtaRange(etaRange);
e1911c19 443 eventMakerOnTheFly->SetNonflowSectorMin(nonflowSectorMin*TMath::Pi()/180.);
444 eventMakerOnTheFly->SetNonflowSectorMax(nonflowSectorMax*TMath::Pi()/180.);
598d292e 445 eventMakerOnTheFly->SetUseGlauberModel(bUseGlauberModel);
1f89f68d 446 eventMakerOnTheFly->SetCreateJets(bCreateJets);
447 eventMakerOnTheFly->SetJetProbability(jetProbability);
448 eventMakerOnTheFly->SetJetTracksFraction(jetTracksFraction);
449 eventMakerOnTheFly->SetJetCone(jetCone);
598d292e 450
451 if(!bUseGlauberModel)
678f8252 452 {
598d292e 453 if(bMultDistrOfRPsIsGauss)
454 {
455 eventMakerOnTheFly->SetMultDistrOfRPsIsGauss(bMultDistrOfRPsIsGauss);
456 eventMakerOnTheFly->SetMultiplicityOfRP(iMultiplicityOfRP);
457 eventMakerOnTheFly->SetMultiplicitySpreadOfRP(dMultiplicitySpreadOfRP);
458 } else
459 {
460 eventMakerOnTheFly->SetMultDistrOfRPsIsGauss(bMultDistrOfRPsIsGauss);
461 eventMakerOnTheFly->SetMinMultOfRP(iMinMultOfRP);
462 eventMakerOnTheFly->SetMaxMultOfRP(iMaxMultOfRP);
463 }
464 } // end of if(!bUseGlauberModel)
678f8252 465
cf90787f 466 eventMakerOnTheFly->SetTemperatureOfRP(dTemperatureOfRP);
bdb5c432 467 eventMakerOnTheFly->SetPtDependentHarmonicV1(bPtDependentHarmonicV1);
468 eventMakerOnTheFly->SetEtaDependentHarmonicV1(bEtaDependentHarmonicV1);
469 eventMakerOnTheFly->SetPtDependentHarmonicV2(bPtDependentHarmonicV2);
470 eventMakerOnTheFly->SetEtaDependentHarmonicV2(bEtaDependentHarmonicV2);
471 eventMakerOnTheFly->SetPtDependentHarmonicV4(bPtDependentHarmonicV4);
472 eventMakerOnTheFly->SetEtaDependentHarmonicV4(bEtaDependentHarmonicV4);
473 // V1:
474 if(!(bPtDependentHarmonicV1 || bEtaDependentHarmonicV1))
524798bf 475 {
598d292e 476 if(!bUseGlauberModel)
477 {
478 eventMakerOnTheFly->SetV1RP(dV1RP);
479 eventMakerOnTheFly->SetV1SpreadRP(dV1SpreadRP);
480 }
bdb5c432 481 } else // (pt,eta) dependent V1
482 {
483 eventMakerOnTheFly->SetV1vsPtEtaMax(dV1vsPtEtaMax);
484 if(bPtDependentHarmonicV1)
485 {
486 eventMakerOnTheFly->SetV1PtCutOff(dV1PtCutOff);
487 }
488 }
489 // V2:
490 if(!(bPtDependentHarmonicV2 || bEtaDependentHarmonicV2)) // constant V2
4d603348 491 {
598d292e 492 if(!bUseGlauberModel)
55381065 493 {
598d292e 494 eventMakerOnTheFly->SetConstantV2IsSampledFromGauss(bConstantV2IsSampledFromGauss);
495 if(bConstantV2IsSampledFromGauss) // Gauss
496 {
497 eventMakerOnTheFly->SetV2RP(dV2RP);
498 eventMakerOnTheFly->SetV2SpreadRP(dV2SpreadRP);
499 } else // uniform
500 {
501 eventMakerOnTheFly->SetMinV2RP(dMinV2RP);
502 eventMakerOnTheFly->SetMaxV2RP(dMaxV2RP);
503 }
504 } // end of if(!bUseGlauberModel)
524798bf 505 } else // (pt,eta) dependent V2
506 {
507 eventMakerOnTheFly->SetV2vsPtEtaMax(dV2vsPtEtaMax);
bdb5c432 508 if(bPtDependentHarmonicV2)
524798bf 509 {
bdb5c432 510 eventMakerOnTheFly->SetV2PtCutOff(dV2PtCutOff);
524798bf 511 }
bdb5c432 512 if(bEtaDependentHarmonicV2)
524798bf 513 {
514 eventMakerOnTheFly->SetV2vsEtaSpread(dV2vsEtaSpread);
515 }
bdb5c432 516 }
65f19276 517 // V3:
518 if(!bUseGlauberModel)
519 {
520 eventMakerOnTheFly->SetV3RP(dV3RP);
521 eventMakerOnTheFly->SetV3SpreadRP(dV3SpreadRP);
522 }
bdb5c432 523 // V4:
524 if(!(bPtDependentHarmonicV4 || bEtaDependentHarmonicV4))
525 {
598d292e 526 if(!bUseGlauberModel)
527 {
528 eventMakerOnTheFly->SetV4RP(dV4RP);
529 eventMakerOnTheFly->SetV4SpreadRP(dV4SpreadRP);
530 }
bdb5c432 531 }
e23cbd39 532
533 // non-uniform acceptance:
534 if(!uniformAcceptance)
535 {
536 eventMakerOnTheFly->SetFirstSectorPhiMin(phimin1);
537 eventMakerOnTheFly->SetFirstSectorPhiMax(phimax1);
538 eventMakerOnTheFly->SetFirstSectorProbability(p1);
539 eventMakerOnTheFly->SetSecondSectorPhiMin(phimin2);
540 eventMakerOnTheFly->SetSecondSectorPhiMax(phimax2);
541 eventMakerOnTheFly->SetSecondSectorProbability(p2);
542 }
e1911c19 543
544 // simple cuts for RPs and POIs:
545 AliFlowTrackSimpleCuts *cutsRP = new AliFlowTrackSimpleCuts();
546 cutsRP->SetPtMax(ptMaxRP);
547 cutsRP->SetPtMin(ptMinRP);
548 cutsRP->SetEtaMax(etaMaxRP);
549 cutsRP->SetEtaMin(etaMinRP);
550 cutsRP->SetPhiMax(phiMaxRP*TMath::Pi()/180.);
551 cutsRP->SetPhiMin(phiMinRP*TMath::Pi()/180.);
552 cutsRP->SetPID(pidRP);
553
554 AliFlowTrackSimpleCuts *cutsPOI = new AliFlowTrackSimpleCuts();
555 cutsPOI->SetPtMax(ptMaxPOI);
556 cutsPOI->SetPtMin(ptMinPOI);
557 cutsPOI->SetEtaMax(etaMaxPOI);
558 cutsPOI->SetEtaMin(etaMinPOI);
559 cutsPOI->SetPhiMax(phiMaxPOI*TMath::Pi()/180.);
560 cutsPOI->SetPhiMin(phiMinPOI*TMath::Pi()/180.);
561 cutsPOI->SetPID(pidPOI);
562
93ff27bd 563 //---------------------------------------------------------------------------------------
564 // create and analyze events 'on the fly':
565
93ff27bd 566 for(Int_t i=0;i<nEvts;i++) {
567 // creating the event with above settings:
e1911c19 568 AliFlowEventSimple *event = eventMakerOnTheFly->CreateEventOnTheFly(cutsRP,cutsPOI);
93ff27bd 569
570 // analyzing the created event 'on the fly':
93ff27bd 571 // do flow analysis for various methods:
c741f5d0 572 if(MCEP) mcep->Make(event);
573 if(QC) qc->Make(event);
574 if(GFC) gfc->Make(event);
575 if(FQD) fqd->Make(event);
576 if(LYZ1SUM) lyz1sum->Make(event);
577 if(LYZ1PROD)lyz1prod->Make(event);
578 if(LYZ2SUM) lyz2sum->Make(event);
579 if(LYZ2PROD)lyz2prod->Make(event);
580 if(LYZEP) lyzep->Make(event,ep);
581 if(SP) sp->Make(event);
d66d46f7 582 if(MH) mh->Make(event);
83bc3e95 583 if(NL) nl->Make(event);
d66d46f7 584 if(MCEP_AH) mcep_ah->Make(event);
ce4a88f5 585
93ff27bd 586 delete event;
587 } // end of for(Int_t i=0;i<nEvts;i++)
588 //---------------------------------------------------------------------------------------
589
93ff27bd 590 //---------------------------------------------------------------------------------------
ad87ae62 591 // create a new file which will hold the final results of all methods:
592 TString outputFileName = "AnalysisResults.root";
593 TFile *outputFile = new TFile(outputFileName.Data(),"RECREATE");
594 // create a new file for each method wich will hold list with final results:
83bc3e95 595 const Int_t nMethods = 13;
596 TString method[nMethods] = {"MCEP","SP","GFC","QC","FQD","LYZ1SUM","LYZ1PROD","LYZ2SUM","LYZ2PROD","LYZEP","MH","NL","MCEP_AH"};
ad87ae62 597 TDirectoryFile *dirFileFinal[nMethods] = {NULL};
598 TString fileName[nMethods];
599 for(Int_t i=0;i<nMethods;i++)
600 {
601 // form a file name for each method:
602 fileName[i]+="output";
603 fileName[i]+=method[i].Data();
604 fileName[i]+="analysis";
605 dirFileFinal[i] = new TDirectoryFile(fileName[i].Data(),fileName[i].Data());
606 }
607
ce4a88f5 608 // calculating and storing the final results of default flow analysis:
ad87ae62 609 if(MCEP) {mcep->Finish(); mcep->WriteHistograms(dirFileFinal[0]);}
610 if(SP) {sp->Finish(); sp->WriteHistograms(dirFileFinal[1]);}
611 if(GFC) {gfc->Finish(); gfc->WriteHistograms(dirFileFinal[2]);}
612 if(QC) {qc->Finish(); qc->WriteHistograms(dirFileFinal[3]);}
613 if(FQD) {fqd->Finish(); fqd->WriteHistograms(dirFileFinal[4]);}
614 if(LYZ1SUM) {lyz1sum->Finish(); lyz1sum->WriteHistograms(dirFileFinal[5]);}
615 if(LYZ1PROD){lyz1prod->Finish();lyz1prod->WriteHistograms(dirFileFinal[6]);}
616 if(LYZ2SUM) {lyz2sum->Finish(); lyz2sum->WriteHistograms(dirFileFinal[7]);}
617 if(LYZ2PROD){lyz2prod->Finish();lyz2prod->WriteHistograms(dirFileFinal[8]);}
618 if(LYZEP) {lyzep->Finish(); lyzep->WriteHistograms(dirFileFinal[9]);}
d66d46f7 619 if(MH) {mh->Finish(); mh->WriteHistograms(dirFileFinal[10]);}
83bc3e95 620 if(NL) {nl->Finish(); nl->WriteHistograms(dirFileFinal[11]);}
621 if(MCEP_AH) {mcep_ah->Finish(); mcep_ah->WriteHistograms(dirFileFinal[12]);}
93ff27bd 622 //---------------------------------------------------------------------------------------
623
ad87ae62 624 outputFile->Close();
625 delete outputFile;
93ff27bd 626
627 cout<<endl;
628 cout<<endl;
629 cout<<" ---- LANDED SUCCESSFULLY ---- "<<endl;
630 cout<<endl;
631
632 timer.Stop();
633 cout << endl;
634 timer.Print();
635}
636
637void SetupPar(char* pararchivename)
638{
639 //Load par files, create analysis libraries
640 //For testing, if par file already decompressed and modified
641 //classes then do not decompress.
642
643 TString cdir(Form("%s", gSystem->WorkingDirectory() )) ;
644 TString parpar(Form("%s.par", pararchivename)) ;
645 if ( gSystem->AccessPathName(parpar.Data()) ) {
646 gSystem->ChangeDirectory(gSystem->Getenv("ALICE_ROOT")) ;
647 TString processline(Form(".! make %s", parpar.Data())) ;
648 gROOT->ProcessLine(processline.Data()) ;
649 gSystem->ChangeDirectory(cdir) ;
650 processline = Form(".! mv /tmp/%s .", parpar.Data()) ;
651 gROOT->ProcessLine(processline.Data()) ;
652 }
653 if ( gSystem->AccessPathName(pararchivename) ) {
654 TString processline = Form(".! tar xvzf %s",parpar.Data()) ;
655 gROOT->ProcessLine(processline.Data());
656 }
657
658 TString ocwd = gSystem->WorkingDirectory();
659 gSystem->ChangeDirectory(pararchivename);
660
661 // check for BUILD.sh and execute
662 if (!gSystem->AccessPathName("PROOF-INF/BUILD.sh")) {
663 printf("*******************************\n");
664 printf("*** Building PAR archive ***\n");
665 cout<<pararchivename<<endl;
666 printf("*******************************\n");
667
668 if (gSystem->Exec("PROOF-INF/BUILD.sh")) {
669 Error("runProcess","Cannot Build the PAR Archive! - Abort!");
670 return -1;
671 }
672 }
673 // check for SETUP.C and execute
674 if (!gSystem->AccessPathName("PROOF-INF/SETUP.C")) {
675 printf("*******************************\n");
676 printf("*** Setup PAR archive ***\n");
677 cout<<pararchivename<<endl;
678 printf("*******************************\n");
679 gROOT->Macro("PROOF-INF/SETUP.C");
680 }
681
682 gSystem->ChangeDirectory(ocwd.Data());
683 printf("Current dir: %s\n", ocwd.Data());
684}
685
bdb5c432 686void CheckUserSettings()
687{
688 // Check if user settings make sense before taking off:
689
690 if (LYZ1SUM && LYZ2SUM) {cout<<"WARNING: you cannot run LYZ1 and LYZ2 at the same time! LYZ2 needs the output from LYZ1. "<<endl; exit(); }
691 if (LYZ1PROD && LYZ2PROD) {cout<<"WARNING: you cannot run LYZ1 and LYZ2 at the same time! LYZ2 needs the output from LYZ1. "<<endl; exit(); }
692 if (LYZ2SUM && LYZEP) {cout<<"WARNING: you cannot run LYZ2 and LYZEP at the same time! LYZEP needs the output from LYZ2."<<endl; exit(); }
693 if (LYZ1SUM && LYZEP) {cout<<"WARNING: you cannot run LYZ1 and LYZEP at the same time! LYZEP needs the output from LYZ2."<<endl; exit(); }
694
695 if(!uniformAcceptance && phimin1 > phimax1)
696 {
697 cout<<"WARNING: you must have phimin1 < phimax1 !!!!"<<endl;
698 break;
699 }
700
701 if (!uniformAcceptance && !((phimin2 == 0.) && (phimax2 == 0.) && (p2 == 0.)) && (phimin2 < phimax1 || phimin2 > phimax2))
702 {
703 cout<<"WARNING: you must have phimin2 > phimax1 and phimin2 < phimax2 !!!!"<<endl;
704 break;
705 }
706
707 if((phimin1 < 0 || phimin1 > 360) || (phimax1 < 0 || phimax1 > 360) ||
708 (phimin2 < 0 || phimin2 > 360) || (phimax2 < 0 || phimax2 > 360) )
709 {
710 cout<<"WARNING: you must take azimuthal angles from interval [0,360] !!!!"<<endl;
711 break;
712 }
713
714 if((p1 < 0 || p1 > 1) || (p2 < 0 || p2 > 1))
715 {
716 cout<<"WARNING: you must take p1 and p2 from interval [0,1] !!!!"<<endl;
717 break;
718 }
719
598d292e 720 if(bUseGlauberModel)
721 {
722 if(bPtDependentHarmonicV1||bEtaDependentHarmonicV1||
723 bPtDependentHarmonicV2||bEtaDependentHarmonicV2||
724 bPtDependentHarmonicV4||bEtaDependentHarmonicV4)
725 {
726 cout<<"WARNING: When using Glauber model pt and/or eta dependence of flow harmonics is not supported !!!!"<<endl;
727 cout<<" Set all booleans bPtDependentHarmonic* to kFALSE in the macro."<<endl;
728 exit(0);
729 }
730 }
731
bdb5c432 732 if(bPtDependentHarmonicV4 && !bPtDependentHarmonicV2))
733 {
734 cout<<"WARNING: V4(pt,eta) is determined as pow(V2(pt,eta),2.) !!!!"<<endl;
735 cout<<" You must also set bPtDependentHarmonicV2 = kTRUE in the macro."<<endl;
736 exit(0);
737 }
738
739 if(bEtaDependentHarmonicV4 && !bEtaDependentHarmonicV2))
740 {
741 cout<<"WARNING: V4(pt,eta) is determined as pow(V2(pt,eta),2.) !!!!"<<endl;
742 cout<<" You must also set bEtaDependentHarmonicV2 = kTRUE in the macro."<<endl;
743 exit(0);
744 }
745
746 if(bPtDependentHarmonicV4 && (bEtaDependentHarmonicV4 != bEtaDependentHarmonicV2))
747 {
748 cout<<"WARNING: V4(pt,eta) is determined as pow(V2(pt,eta),2.) !!!!"<<endl;
749 cout<<" You must also have bEtaDependentHarmonicV2 = bEtaDependentHarmonicV4 in the macro."<<endl;
750 exit(0);
751 }
752
753 if(bEtaDependentHarmonicV4 && (bPtDependentHarmonicV4 != bPtDependentHarmonicV2))
754 {
755 cout<<"WARNING: V4(pt,eta) is determined as pow(V2(pt,eta),2.) !!!!"<<endl;
756 cout<<" You must also have bPtDependentHarmonicV2 = bPtDependentHarmonicV4 in the macro."<<endl;
757 exit(0);
758 }
1f89f68d 759
760 if(!uniformAcceptance && bCreateJets)
761 {
762 cout<<"WARNING: Jet creation is supported only for uniform acceptance !!!!"<<endl;
763 exit(0);
764 }
598d292e 765
bdb5c432 766} // end of void CheckUserSettings()
767
93ff27bd 768void LoadLibraries(const anaModes mode) {
769
770 //--------------------------------------
771 // Load the needed libraries most of them already loaded by aliroot
772 //--------------------------------------
ad87ae62 773 //gSystem->Load("libTree");
5d040cf3 774 gSystem->Load("libGeom");
775 gSystem->Load("libVMC");
776 gSystem->Load("libXMLIO");
777 gSystem->Load("libPhysics");
93ff27bd 778
779 //----------------------------------------------------------
780 // >>>>>>>>>>> Local mode <<<<<<<<<<<<<<
781 //----------------------------------------------------------
782 if (mode==mLocal) {
783 //--------------------------------------------------------
784 // If you want to use already compiled libraries
785 // in the aliroot distribution
786 //--------------------------------------------------------
787 gSystem->Load("libSTEERBase");
788 gSystem->Load("libESD");
789 gSystem->Load("libAOD");
790 gSystem->Load("libANALYSIS");
791 gSystem->Load("libANALYSISalice");
5d040cf3 792 gSystem->Load("libCORRFW");
793 cerr<<"libCORRFW loaded..."<<endl;
794 gSystem->Load("libPWG2flowCommon");
795 cerr<<"libPWG2flowCommon loaded..."<<endl;
796 gSystem->Load("libPWG2flowTasks");
797 cerr<<"libPWG2flowTasks loaded..."<<endl;
93ff27bd 798 }
799
800 else if (mode == mLocalPAR) {
801 //--------------------------------------------------------
802 //If you want to use root and par files from aliroot
803 //--------------------------------------------------------
804 //If you want to use root and par files from aliroot
805 //--------------------------------------------------------
806 SetupPar("STEERBase");
807 SetupPar("ESD");
808 SetupPar("AOD");
809 SetupPar("ANALYSIS");
810 SetupPar("ANALYSISalice");
811 SetupPar("PWG2AOD");
812 SetupPar("CORRFW");
813 SetupPar("PWG2flowCommon");
814 cerr<<"PWG2flowCommon.par loaded..."<<endl;
815 SetupPar("PWG2flowTasks");
816 cerr<<"PWG2flowTasks.par loaded..."<<endl;
817 }
818
819 //---------------------------------------------------------
820 // <<<<<<<<<< Source mode >>>>>>>>>>>>
821 //---------------------------------------------------------
822 else if (mode==mLocalSource) {
823
824 // In root inline compile
825
93ff27bd 826 // Constants
827 gROOT->LoadMacro("AliFlowCommon/AliFlowCommonConstants.cxx+");
828 gROOT->LoadMacro("AliFlowCommon/AliFlowLYZConstants.cxx+");
93ff27bd 829
830 // Flow event
831 gROOT->LoadMacro("AliFlowCommon/AliFlowVector.cxx+");
5062cb0f 832 gROOT->LoadMacro("AliFlowCommon/AliFlowTrackSimple.cxx+");
93ff27bd 833 gROOT->LoadMacro("AliFlowCommon/AliFlowTrackSimpleCuts.cxx+");
5062cb0f 834 gROOT->LoadMacro("AliFlowCommon/AliFlowEventSimple.cxx+");
93ff27bd 835
836 // Output histosgrams
837 gROOT->LoadMacro("AliFlowCommon/AliFlowCommonHist.cxx+");
838 gROOT->LoadMacro("AliFlowCommon/AliFlowCommonHistResults.cxx+");
839 gROOT->LoadMacro("AliFlowCommon/AliFlowLYZHist1.cxx+");
840 gROOT->LoadMacro("AliFlowCommon/AliFlowLYZHist2.cxx+");
841
842 // Functions needed for various methods
843 gROOT->LoadMacro("AliFlowCommon/AliCumulantsFunctions.cxx+");
93ff27bd 844 gROOT->LoadMacro("AliFlowCommon/AliFlowLYZEventPlane.cxx+");
845
846 // Flow Analysis code for various methods
847 gROOT->LoadMacro("AliFlowCommon/AliFlowAnalysisWithMCEventPlane.cxx+");
848 gROOT->LoadMacro("AliFlowCommon/AliFlowAnalysisWithScalarProduct.cxx+");
849 gROOT->LoadMacro("AliFlowCommon/AliFlowAnalysisWithLYZEventPlane.cxx+");
850 gROOT->LoadMacro("AliFlowCommon/AliFlowAnalysisWithLeeYangZeros.cxx+");
851 gROOT->LoadMacro("AliFlowCommon/AliFlowAnalysisWithCumulants.cxx+");
852 gROOT->LoadMacro("AliFlowCommon/AliFlowAnalysisWithQCumulants.cxx+");
ce4a88f5 853 gROOT->LoadMacro("AliFlowCommon/AliFlowAnalysisWithFittingQDistribution.cxx+");
ecac11c2 854 gROOT->LoadMacro("AliFlowCommon/AliFlowAnalysisWithMixedHarmonics.cxx+");
855 gROOT->LoadMacro("AliFlowCommon/AliFlowAnalysisWithNestedLoops.cxx+");
93ff27bd 856
a47dc7e4 857 // Class to fill the FlowEvent on the fly (generate Monte Carlo events)
858 gROOT->LoadMacro("AliFlowCommon/AliFlowEventSimpleMakerOnTheFly.cxx+");
93ff27bd 859
860 cout << "finished loading macros!" << endl;
861
862 }
863
864}
865
866