hooks for PMD flow analysis
[u/mrichter/AliRoot.git] / PWG2 / FLOW / macros / runFlowAnalysis.C
CommitLineData
929098e4 1#include <Riostream.h>
2#include <TStopwatch.h>
3#include <TObjArray.h>
4#include <TFile.h>
5#include <TRandom3.h>
6#include <TTimeStamp.h>
7
8#include <AliFlowCommon/AliFlowCommonConstants.h>
9#include <AliFlowCommon/AliFlowLYZConstants.h>
10#include <AliFlowCommon/AliFlowCumuConstants.h>
11#include <AliFlowCommon/AliFlowVector.h>
12#include <AliFlowCommon/AliFlowTrackSimple.h>
7382279b 13#include <AliFlowCommon/AliFlowEvent.h>
929098e4 14#include <AliFlowCommon/AliFlowEventSimple.h>
15#include <AliFlowCommon/AliFlowTrackSimpleCuts.h>
16#include <AliFlowCommon/AliFlowCommonHist.h>
17#include <AliFlowCommon/AliFlowCommonHistResults.h>
18#include <AliFlowCommon/AliFlowLYZHist1.h>
19#include <AliFlowCommon/AliFlowLYZHist2.h>
20#include <AliFlowCommon/AliCumulantsFunctions.h>
21#include <AliFlowCommon/AliFlowLYZEventPlane.h>
22#include <AliFlowCommon/AliFlowAnalysisWithMCEventPlane.h>
23#include <AliFlowCommon/AliFlowAnalysisWithScalarProduct.h>
24#include <AliFlowCommon/AliFlowAnalysisWithLYZEventPlane.h>
25#include <AliFlowCommon/AliFlowAnalysisWithLeeYangZeros.h>
26#include <AliFlowCommon/AliFlowAnalysisWithCumulants.h>
27#include <AliFlowCommon/AliFlowAnalysisWithQCumulants.h>
28#include <AliFlowCommon/AliFlowAnalysisWithFittingQDistribution.h>
29#include <AliFlowCommon/AliFlowAnalysisWithMixedHarmonics.h>
30#include <AliFlowCommon/AliFlowAnalysisWithNestedLoops.h>
334e3256 31
03a02aca 32//--------------------------------------------------------------------------------------
c852719e 33// Run flow analysis on local data with custom FlowEvent maker
03a02aca 34// RUN SETTINGS
334e3256 35//flow analysis method can be: (set to kTRUE or kFALSE)
c741f5d0 36Bool_t SP = kTRUE;
5b40431d 37Bool_t LYZ1SUM = kTRUE;
38Bool_t LYZ1PROD = kTRUE;
c741f5d0 39Bool_t LYZ2SUM = kFALSE;
40Bool_t LYZ2PROD = kFALSE;
5b40431d 41Bool_t LYZEP = kFALSE;
c741f5d0 42Bool_t GFC = kTRUE;
43Bool_t QC = kTRUE;
b5425e51 44Bool_t FQD = kTRUE;
ecac11c2 45Bool_t MH = kTRUE;
46Bool_t NL = kFALSE;
5b40431d 47Bool_t MCEP = kFALSE; //does not work yet 24/12/08
03a02aca 48//--------------------------------------------------------------------------------------
334e3256 49
20e52807 50// Weights
51// Use weights for Q vector
52Bool_t usePhiWeights = kFALSE; //Phi (correction for non-uniform azimuthal acceptance)
53Bool_t usePtWeights = kFALSE; //v'(pt) (differential flow in pt)
54Bool_t useEtaWeights = kFALSE; //v'(eta) (differential flow in eta)
55
03a02aca 56//--------------------------------------------------------------------------------------
57// CUT SETTINGS
334e3256 58//integrated selection
59Double_t ptMaxInt = 10.;
60Double_t ptMinInt = 0.;
61Double_t etaMaxInt = 1.;
62Double_t etaMinInt = -1.;
63Double_t phiMaxInt = 7.5;
64Double_t phiMinInt = 0.;
65Int_t PIDInt = 211;
66
67//differential selection
68Double_t ptMaxDiff = 10.;
69Double_t ptMinDiff = 0.;
70Double_t etaMaxDiff = 1.;
71Double_t etaMinDiff = -1.;
72Double_t phiMaxDiff = 7.5;
73Double_t phiMinDiff = 0.;
74Int_t PIDDiff = 211;
c741f5d0 75/*
03a02aca 76//--------------------------------------------------------------------------------------
33314629 77// FLOW SETTINGS (R.Rietkerk)
78Int_t nLoops=1; // Number of times to use the same particle (nonflow).
c741f5d0 79Double_t xEllipticFlowValue=0.1;// Add Elliptic Flow. Must be in range [0,1].
80Int_t nMultiplicityOfEvent=500; // Set Average Multiplicity.
33314629 81Double_t xSigmaFlow=0.00; // Add Elliptic Flow. Must be in range [0,1].
c741f5d0 82Int_t nSigmaMult=50; // Set Average Multiplicity.
33314629 83//--------------------------------------------------------------------------------------
c741f5d0 84*/
26c4cbb9 85
fdd8c18f 86enum anaModes {mLocal,mLocalSource,mLocalPAR,};
87//mLocal: Analyze data on your computer using aliroot
88//mLocalPAR: Analyze data on your computer using root + PAR files
89//mLocalSource: Analyze data on your computer using root + source files
90
929098e4 91void LoadLibraries(const anaModes mode);
92
26c4cbb9 93Int_t offset = 0;
94
929098e4 95int runFlowAnalysis(const anaModes mode=mLocal, Int_t aRuns = 100, const char*
c741f5d0 96 dir="/data/alice1/kolk/KineOnly3/")
97 // dir="/Users/snelling/alice_data/KineOnly3/")
33314629 98 // dir="/Users/snelling/alice_data/stoomboot/5b/")
334e3256 99{
100 TStopwatch timer;
101 timer.Start();
102
929098e4 103 if (LYZ1SUM && LYZ2SUM) {cout<<"WARNING: you cannot run LYZ1 and LYZ2 at the same time! LYZ2 needs the output from LYZ1."<<endl; exit(1); }
104 if (LYZ1PROD && LYZ2PROD) {cout<<"WARNING: you cannot run LYZ1 and LYZ2 at the same time! LYZ2 needs the output from LYZ1."<<endl; exit(1); }
105 if (LYZ2SUM && LYZEP) {cout<<"WARNING: you cannot run LYZ2 and LYZEP at the same time! LYZEP needs the output from LYZ2."<<endl; exit(1); }
106 if (LYZ1SUM && LYZEP) {cout<<"WARNING: you cannot run LYZ1 and LYZEP at the same time! LYZEP needs the output from LYZ2."<<endl; exit(1); }
334e3256 107
108
109 cout<<endl;
110 cout<<" ---- BEGIN ANALYSIS ---- "<<endl;
111 cout<<endl;
112
fdd8c18f 113 LoadLibraries(mode);
2a413443 114
33314629 115 TRandom3 random3Temp; //init for manual settings (R.Rietkerk)
116 TTimeStamp dt;
117 Int_t sseed = dt.GetNanoSec()/1000;
118 random3Temp.SetSeed(sseed);
119
fdd8c18f 120 if (mode == mLocal || mode == mLocalPAR) {
fdd8c18f 121 }
122 else if (mode == mLocalSource) {
fdd8c18f 123 }
124 else{
125 cout << "No supported running mode selected!" << endl;
126 break;
127 }
2a413443 128
2a413443 129
334e3256 130 //------------------------------------------------------------------------
e085f1a9 131 //cuts:
334e3256 132 AliFlowTrackSimpleCuts* cutsInt = new AliFlowTrackSimpleCuts();
133 cutsInt->SetPtMax(ptMaxInt);
134 cutsInt->SetPtMin(ptMinInt);
135 cutsInt->SetEtaMax(etaMaxInt);
136 cutsInt->SetEtaMin(etaMinInt);
137 cutsInt->SetPhiMax(phiMaxInt);
138 cutsInt->SetPhiMin(phiMinInt);
139 cutsInt->SetPID(PIDInt);
140
141 AliFlowTrackSimpleCuts* cutsDiff = new AliFlowTrackSimpleCuts();
142 cutsDiff->SetPtMax(ptMaxDiff);
143 cutsDiff->SetPtMin(ptMinDiff);
144 cutsDiff->SetEtaMax(etaMaxDiff);
145 cutsDiff->SetEtaMin(etaMinDiff);
146 cutsDiff->SetPhiMax(phiMaxDiff);
147 cutsDiff->SetPhiMin(phiMinDiff);
148 cutsDiff->SetPID(PIDDiff);
149
03a02aca 150 //if the weights are used:
151 TFile *fileWithWeights = NULL;
152 TList *listWithWeights = NULL;
334e3256 153
14db9c04 154 if(usePhiWeights||usePtWeights||useEtaWeights) {
155 fileWithWeights = TFile::Open("weights.root","READ");
156 if(fileWithWeights) {
157 listWithWeights = (TList*)fileWithWeights->Get("weights");
158 }
159 else
160 {cout << " WARNING: the file <weights.root> with weights from the previous run was not found."<<endl;
161 break;
162 }
e085f1a9 163 }
164
165 //flow methods:
ce4a88f5 166 AliFlowAnalysisWithQCumulants *qc = NULL;
167 AliFlowAnalysisWithCumulants *gfc = NULL;
168 AliFlowAnalysisWithFittingQDistribution *fqd = NULL;
169 AliFlowAnalysisWithLeeYangZeros *lyz1sum = NULL;
170 AliFlowAnalysisWithLeeYangZeros *lyz1prod = NULL;
171 AliFlowAnalysisWithLeeYangZeros *lyz2sum = NULL;
172 AliFlowAnalysisWithLeeYangZeros *lyz2prod = NULL;
173 AliFlowAnalysisWithLYZEventPlane *lyzep = NULL;
174 AliFlowAnalysisWithScalarProduct *sp = NULL;
ecac11c2 175 AliFlowAnalysisWithMCEventPlane *mcep = NULL;
176 AliFlowAnalysisWithMixedHarmonics *mh = NULL;
177 AliFlowAnalysisWithNestedLoops *nl = NULL;
deebd72f 178
334e3256 179 //MCEP = monte carlo event plane
deebd72f 180 if (MCEP) {
181 AliFlowAnalysisWithMCEventPlane *mcep = new AliFlowAnalysisWithMCEventPlane();
182 mcep->Init();
183 }
334e3256 184
185 //QC = Q-cumulants
deebd72f 186 if(QC) {
187 AliFlowAnalysisWithQCumulants* qc = new AliFlowAnalysisWithQCumulants();
03a02aca 188 if(listWithWeights) qc->SetWeightsList(listWithWeights);
189 if(usePhiWeights) qc->SetUsePhiWeights(usePhiWeights);
190 if(usePtWeights) qc->SetUsePtWeights(usePtWeights);
191 if(useEtaWeights) qc->SetUseEtaWeights(useEtaWeights);
ecac11c2 192 qc->Init();
deebd72f 193 }
334e3256 194
195 //GFC = Generating Function Cumulants
deebd72f 196 if(GFC) {
197 AliFlowAnalysisWithCumulants* gfc = new AliFlowAnalysisWithCumulants();
e5e75b58 198 if(listWithWeights) gfc->SetWeightsList(listWithWeights);
199 if(usePhiWeights) gfc->SetUsePhiWeights(usePhiWeights);
200 if(usePtWeights) gfc->SetUsePtWeights(usePtWeights);
201 if(useEtaWeights) gfc->SetUseEtaWeights(useEtaWeights);
ecac11c2 202 gfc->Init();
deebd72f 203 }
334e3256 204
205 //FQD = Fitting q-distribution
deebd72f 206 if(FQD) {
ce4a88f5 207 AliFlowAnalysisWithFittingQDistribution* fqd = new AliFlowAnalysisWithFittingQDistribution();
14db9c04 208 if(listWithWeights) fqd->SetWeightsList(listWithWeights);
209 if(usePhiWeights) fqd->SetUsePhiWeights(usePhiWeights);
ecac11c2 210 fqd->Init();
deebd72f 211 }
60ec66fb 212
213 //SP = Scalar Product
214 if(SP) {
215 AliFlowAnalysisWithScalarProduct* sp = new AliFlowAnalysisWithScalarProduct();
ecac11c2 216 if(usePhiWeights) sp->SetUsePhiWeights(usePhiWeights);
60ec66fb 217 sp->Init();
218 }
219
334e3256 220 //LYZ1 = Lee-Yang Zeroes first run
c741f5d0 221 if(LYZ1SUM) {
222 AliFlowAnalysisWithLeeYangZeros* lyz1sum = new AliFlowAnalysisWithLeeYangZeros();
223 lyz1sum->SetFirstRun(kTRUE);
224 lyz1sum->SetUseSum(kTRUE);
225 lyz1sum->Init();
226 }
227 if(LYZ1PROD) {
228 AliFlowAnalysisWithLeeYangZeros* lyz1prod = new AliFlowAnalysisWithLeeYangZeros();
229 lyz1prod->SetFirstRun(kTRUE);
230 lyz1prod->SetUseSum(kFALSE);
231 lyz1prod->Init();
deebd72f 232 }
334e3256 233 //LYZ2 = Lee-Yang Zeroes second run
c741f5d0 234 if(LYZ2SUM) {
235 AliFlowAnalysisWithLeeYangZeros* lyz2sum = new AliFlowAnalysisWithLeeYangZeros();
deebd72f 236 // read the input file from the first run
c741f5d0 237 TString inputFileNameLYZ2SUM = "outputLYZ1SUManalysis.root" ;
238 TFile* inputFileLYZ2SUM = new TFile(inputFileNameLYZ2SUM.Data(),"READ");
239 if(!inputFileLYZ2SUM || inputFileLYZ2SUM->IsZombie()) {
5b40431d 240 cerr << " ERROR: To run LYZ2SUM you need the output file from LYZ1SUM. This file is not there! Please run LYZ1SUM first." << endl ;
ce78b98c 241 break;
deebd72f 242 }
243 else {
c741f5d0 244 TList* inputListLYZ2SUM = (TList*)inputFileLYZ2SUM->Get("cobjLYZ1SUM");
245 if (!inputListLYZ2SUM) {cout<<"SUM Input list is NULL pointer!"<<endl; break;}
deebd72f 246 else {
c741f5d0 247 cout<<"LYZ2SUM input file/list read..."<<endl;
248 lyz2sum->SetFirstRunList(inputListLYZ2SUM);
249 lyz2sum->SetFirstRun(kFALSE);
250 lyz2sum->SetUseSum(kTRUE);
251 lyz2sum->Init();
252 }
253 }
254 }
255 if(LYZ2PROD) {
256 AliFlowAnalysisWithLeeYangZeros* lyz2prod = new AliFlowAnalysisWithLeeYangZeros();
257 // read the input file from the first run
258 TString inputFileNameLYZ2PROD = "outputLYZ1PRODanalysis.root" ;
259 TFile* inputFileLYZ2PROD = new TFile(inputFileNameLYZ2PROD.Data(),"READ");
260 if(!inputFileLYZ2PROD || inputFileLYZ2PROD->IsZombie()) {
5b40431d 261 cerr << " ERROR: To run LYZ2PROD you need the output file from LYZ1PROD. This file is not there! Please run LYZ1PROD first." << endl ;
c741f5d0 262 break;
263 }
264 else {
265 TList* inputListLYZ2PROD = (TList*)inputFileLYZ2PROD->Get("cobjLYZ1PROD");
266 if (!inputListLYZ2PROD) {cout<<"PROD Input list is NULL pointer!"<<endl; break;}
267 else {
268 cout<<"LYZ2PROD input file/list read..."<<endl;
269 lyz2prod->SetFirstRunList(inputListLYZ2PROD);
270 lyz2prod->SetFirstRun(kFALSE);
271 lyz2prod->SetUseSum(kTRUE);
272 lyz2prod->Init();
334e3256 273 }
274 }
deebd72f 275 }
334e3256 276 //LYZEP = Lee-Yang Zeroes event plane
deebd72f 277 if(LYZEP) {
278 AliFlowLYZEventPlane* ep = new AliFlowLYZEventPlane() ;
279 AliFlowAnalysisWithLYZEventPlane* lyzep = new AliFlowAnalysisWithLYZEventPlane();
280 // read the input file from the second lyz run
c741f5d0 281 TString inputFileNameLYZEP = "outputLYZ2SUManalysis.root" ;
deebd72f 282 TFile* inputFileLYZEP = new TFile(inputFileNameLYZEP.Data(),"READ");
283 if(!inputFileLYZEP || inputFileLYZEP->IsZombie()) {
5b40431d 284 cerr << " ERROR: To run LYZEP you need the output file from LYZ2SUM. This file is not there! Please run LYZ2SUM first." << endl ;
ce78b98c 285 break;
286 }
deebd72f 287 else {
c741f5d0 288 TList* inputListLYZEP = (TList*)inputFileLYZEP->Get("cobjLYZ2SUM");
ace38bad 289 if (!inputListLYZEP) {cout<<"Input list is NULL pointer!"<<endl; break;}
deebd72f 290 else {
291 cout<<"LYZEP input file/list read..."<<endl;
292 ep ->SetSecondRunList(inputListLYZEP);
293 lyzep->SetSecondRunList(inputListLYZEP);
294 ep ->Init();
295 lyzep->Init();
334e3256 296 }
297 }
deebd72f 298 }
ecac11c2 299 // MH = Mixed Harmonics:
300 if(MH) {
301 AliFlowAnalysisWithMixedHarmonics* mh = new AliFlowAnalysisWithMixedHarmonics();
302 if(listWithWeights) mh->SetWeightsList(listWithWeights);
303 //if(usePhiWeights) mh->SetUsePhiWeights(usePhiWeights); // to be improved (enabled)
304 //if(usePtWeights) mh->SetUsePtWeights(usePtWeights); // to be improved (enabled)
305 //if(useEtaWeights) mh->SetUseEtaWeights(useEtaWeights); // to be improved (enabled)
306 mh->Init();
307 }
308 // NL = Nested Loops:
309 if(NL) {
310 AliFlowAnalysisWithNestedLoops* nl = new AliFlowAnalysisWithNestedLoops();
311 if(listWithWeights) nl->SetWeightsList(listWithWeights);
312 //if(usePhiWeights) nl->SetUsePhiWeights(usePhiWeights); // to be improved (enabled)
313 //if(usePtWeights) nl->SetUsePtWeights(usePtWeights); // to be improved (enabled)
314 //if(useEtaWeights) nl->SetUseEtaWeights(useEtaWeights); // to be improved (enabled)
315 nl->Init();
316 }
317
334e3256 318 //------------------------------------------------------------------------
319
26c4cbb9 320
b25cc698 321 //standard code to read files in directory
334e3256 322 Int_t fCount = 0;
323 TString execDir(gSystem->pwd());
324 TString targetDir(dir);
325 TSystemDirectory* baseDir = new TSystemDirectory(".", dir);
326 TList* dirList = baseDir->GetListOfFiles();
327 if (!dirList) {
328 cout << endl << "No input files in: " << targetDir.Data() << endl;
329 break;
330 }
331 Int_t nDirs = dirList->GetEntries();
332 cout<<endl;
333 cout<<"Int_t nDirs = "<<nDirs<<endl;
334 gSystem->cd(execDir);
335
336 for(Int_t iDir=0;iDir<nDirs;++iDir)
337 {
338 TSystemFile* presentDir = (TSystemFile*)dirList->At(iDir);
339 if(!presentDir || !presentDir->IsDirectory() ||
340 strcmp(presentDir->GetName(), ".") == 0 ||
341 strcmp(presentDir->GetName(), "..") == 0)
342 {
343 cout << endl;
344 cout << "Directory (" << iDir << "): " << presentDir->GetName() <<
345 " - Skipping ... " << endl;
346 continue ;
347 }
348
349 if(offset > 0) { --offset ; continue ; }
350 if((aRuns > 0) && (fCount >= aRuns)) { break ; }
351
352 TString presentDirName(dir); // aDataDir
353 presentDirName += presentDir->GetName();
354 presentDirName += "/";
355 //cerr<<" presentDirName = "<<presentDirName<<endl;
356
357 TString fileName = presentDirName;
358 fileName += "galice.root";
359 Long_t *id, *size, *flags, *modtime;
360 if(gSystem->GetPathInfo(fileName.Data(),id,size,flags,modtime))
361 {
362 cout << " File : " << fileName << " does NOT exist ! - Skipping ... "
363 << endl;
364 continue;
365 }
deebd72f 366 // cout << endl ; cout << "Directory (" << iDir << "): " << presentDirName << " ... " << endl;
334e3256 367
368 //loop (simulations in the present dir)
369 TSystemDirectory* evtsDir = new TSystemDirectory(".",
370 presentDirName.Data());
371 TList* fileList = evtsDir->GetListOfFiles();
372 Int_t nFiles = fileList->GetEntries();
deebd72f 373 //cout<<" Int_t nFiles = "<<nFiles<<endl;
334e3256 374 gSystem->cd(execDir);
375 for(Int_t iFiles=0; iFiles<nFiles; ++iFiles)
376 {
377 TSystemFile* presentFile = (TSystemFile*) fileList->At(iFiles);
378 TString presentFileName(presentDirName);
379 presentFileName += presentFile->GetName();
380
381 if(!(presentFileName.Contains("Kinematics") &&
382 presentFileName.Contains("root"))) { continue ; }
383
deebd72f 384 //cout << " found: " << presentFileName.Data() << endl;
334e3256 385
386 TFile* kineFile = new TFile(presentFileName.Data(), "READ");
387 // kineFile->ls();
388 Int_t nEvts = kineFile->GetNkeys() ;
deebd72f 389 //cout << " . found: " << nEvts << " KineTree(s) in " << presentFileName.Data() << endl;
334e3256 390 TList* kineEventsList = (TList*)kineFile->GetListOfKeys();
391 TTree* kTree;
392 TIter next(kineEventsList);
393 TKey* key;
33314629 394
395 Double_t xRPAngle;
334e3256 396 //loop over the events
397 while( key=(TKey *)next() )
398 {
399 TDirectory* tDir = (TDirectory*)key->ReadObj();
400 if(!tDir) break;
401
402 TString evtDir(tDir->GetName());
deebd72f 403 //cout << " . . found: " << tDir->GetName() << endl;
334e3256 404
405 kTree = (TTree *)tDir->Get("TreeK");
406 if(!kTree) break;
407
408 Int_t nPart = kTree->GetEntries();
deebd72f 409 //cout << " . . . kTree " << fCount << " has " << nPart << " particles " << endl;
334e3256 410
411 //-----------------------------------------------------------
26c4cbb9 412 //fill and save the flow event
33314629 413
c741f5d0 414 /*
33314629 415 Int_t nNewMultOfEvent = random3Temp.Gaus(nMultiplicityOfEvent,nSigmaMult);
416 cout << "new multiplicity: " << nNewMultOfEvent << endl;
417 Double_t xNewFlowValue = random3Temp.Gaus(xEllipticFlowValue,xSigmaFlow);
4f678629 418 if ( (fCount % 100) == 0) {
419 cout << "new multiplicity: " << nNewMultOfEvent << endl;
420 cout << "new flow value: " << xNewFlowValue << endl;
421 }
c741f5d0 422 */
33314629 423
34b15925 424 AliFlowEventSimple *fEvent = new AliFlowEventSimple(kTree, cutsInt, cutsDiff);
425 //xRPAngle=random3Temp.Uniform(0.0,TMath::TwoPi());
426 //fEvent->SetMCReactionPlaneAngle(xRPAngle);
427 //fEvent->SetV2(xNewFlowValue);
428 //fEvent->CloneTracks(nLoops);
26c4cbb9 429
b25cc698 430 // do flow analysis for various methods
c741f5d0 431 if(MCEP) mcep->Make(fEvent);
432 if(QC) qc->Make(fEvent);
433 if(GFC) gfc->Make(fEvent);
434 if(FQD) fqd->Make(fEvent);
435 if(LYZ1SUM) lyz1sum->Make(fEvent);
436 if(LYZ1PROD)lyz1prod->Make(fEvent);
437 if(LYZ2SUM) lyz2sum->Make(fEvent);
438 if(LYZ2PROD)lyz2prod->Make(fEvent);
439 if(LYZEP) lyzep->Make(fEvent,ep);
ecac11c2 440 if(SP) sp->Make(fEvent);
441 if(MH) mh->Make(fEvent);
442 if(NL) nl->Make(fEvent);
334e3256 443 //-----------------------------------------------------------
334e3256 444 fCount++;
951817d5 445 //cout << "# " << fCount << " events processed" << endl;
334e3256 446 delete kTree;
447 delete fEvent;
448 }
449 delete kineFile ;
450 }
451 delete evtsDir ;
452 }
453
b5425e51 454 //---------------------------------------------------------------------------------------
455 // create a new file which will hold the final results of all methods:
456 TString outputFileName = "AnalysisResults.root";
457 TFile *outputFile = new TFile(outputFileName.Data(),"RECREATE");
458 // create a new file for each method wich will hold list with final results:
ecac11c2 459 const Int_t nMethods = 12;
460 TString method[nMethods] = {"MCEP","SP","GFC","QC","FQD","LYZ1SUM","LYZ1PROD","LYZ2SUM","LYZ2PROD","LYZEP","MH","NL"};
b5425e51 461 TDirectoryFile *dirFileFinal[nMethods] = {NULL};
462 TString fileNameMethod[nMethods];
463 for(Int_t i=0;i<nMethods;i++)
464 {
465 // form a file name for each method:
466 fileNameMethod[i]+="output";
467 fileNameMethod[i]+=method[i].Data();
468 fileNameMethod[i]+="analysis";
469 dirFileFinal[i] = new TDirectoryFile(fileNameMethod[i].Data(),fileNameMethod[i].Data());
470 }
471
472 // calculating and storing the final results of default flow analysis:
473 if(MCEP) {mcep->Finish(); mcep->WriteHistograms(dirFileFinal[0]);}
474 if(SP) {sp->Finish(); sp->WriteHistograms(dirFileFinal[1]);}
475 if(GFC) {gfc->Finish(); gfc->WriteHistograms(dirFileFinal[2]);}
476 if(QC) {qc->Finish(); qc->WriteHistograms(dirFileFinal[3]);}
477 if(FQD) {fqd->Finish(); fqd->WriteHistograms(dirFileFinal[4]);}
478 if(LYZ1SUM) {lyz1sum->Finish(); lyz1sum->WriteHistograms(dirFileFinal[5]);}
479 if(LYZ1PROD){lyz1prod->Finish();lyz1prod->WriteHistograms(dirFileFinal[6]);}
480 if(LYZ2SUM) {lyz2sum->Finish(); lyz2sum->WriteHistograms(dirFileFinal[7]);}
481 if(LYZ2PROD){lyz2prod->Finish();lyz2prod->WriteHistograms(dirFileFinal[8]);}
482 if(LYZEP) {lyzep->Finish(); lyzep->WriteHistograms(dirFileFinal[9]);}
ecac11c2 483 if(MH) {mh->Finish(); mh->WriteHistograms(dirFileFinal[10]);}
484 if(NL) {nl->Finish(); nl->WriteHistograms(dirFileFinal[11]);}
b5425e51 485 //---------------------------------------------------------------------------------------
486
487 outputFile->Close();
488 delete outputFile;
334e3256 489
490 cout << endl;
ace38bad 491 cout << " Fini ... " << endl;
334e3256 492 cout << endl;
493
494 timer.Stop();
495 cout << endl;
496 timer.Print();
497}
21e694dd 498
499void SetupPar(char* pararchivename)
500{
501 //Load par files, create analysis libraries
502 //For testing, if par file already decompressed and modified
503 //classes then do not decompress.
504
505 TString cdir(Form("%s", gSystem->WorkingDirectory() )) ;
506 TString parpar(Form("%s.par", pararchivename)) ;
507 if ( gSystem->AccessPathName(parpar.Data()) ) {
508 gSystem->ChangeDirectory(gSystem->Getenv("ALICE_ROOT")) ;
509 TString processline(Form(".! make %s", parpar.Data())) ;
510 gROOT->ProcessLine(processline.Data()) ;
511 gSystem->ChangeDirectory(cdir) ;
512 processline = Form(".! mv /tmp/%s .", parpar.Data()) ;
513 gROOT->ProcessLine(processline.Data()) ;
514 }
515 if ( gSystem->AccessPathName(pararchivename) ) {
516 TString processline = Form(".! tar xvzf %s",parpar.Data()) ;
517 gROOT->ProcessLine(processline.Data());
518 }
519
520 TString ocwd = gSystem->WorkingDirectory();
521 gSystem->ChangeDirectory(pararchivename);
522
523 // check for BUILD.sh and execute
524 if (!gSystem->AccessPathName("PROOF-INF/BUILD.sh")) {
525 printf("*******************************\n");
526 printf("*** Building PAR archive ***\n");
527 cout<<pararchivename<<endl;
528 printf("*******************************\n");
529
530 if (gSystem->Exec("PROOF-INF/BUILD.sh")) {
531 Error("runProcess","Cannot Build the PAR Archive! - Abort!");
532 return -1;
533 }
534 }
535 // check for SETUP.C and execute
536 if (!gSystem->AccessPathName("PROOF-INF/SETUP.C")) {
537 printf("*******************************\n");
538 printf("*** Setup PAR archive ***\n");
539 cout<<pararchivename<<endl;
540 printf("*******************************\n");
541 gROOT->Macro("PROOF-INF/SETUP.C");
542 }
543
544 gSystem->ChangeDirectory(ocwd.Data());
545 printf("Current dir: %s\n", ocwd.Data());
546}
547
fdd8c18f 548void LoadLibraries(const anaModes mode) {
549
550 //--------------------------------------
551 // Load the needed libraries most of them already loaded by aliroot
552 //--------------------------------------
b5425e51 553 //gSystem->Load("libTree");
5d040cf3 554 gSystem->Load("libGeom");
555 gSystem->Load("libVMC");
556 gSystem->Load("libXMLIO");
557 gSystem->Load("libPhysics");
fdd8c18f 558
559 //----------------------------------------------------------
560 // >>>>>>>>>>> Local mode <<<<<<<<<<<<<<
561 //----------------------------------------------------------
562 if (mode==mLocal) {
563 //--------------------------------------------------------
564 // If you want to use already compiled libraries
565 // in the aliroot distribution
566 //--------------------------------------------------------
567 gSystem->Load("libSTEERBase");
568 gSystem->Load("libESD");
569 gSystem->Load("libAOD");
570 gSystem->Load("libANALYSIS");
571 gSystem->Load("libANALYSISalice");
5d040cf3 572 gSystem->Load("libCORRFW");
573 cerr<<"libCORRFW loaded..."<<endl;
574 gSystem->Load("libPWG2flowCommon");
575 cerr<<"libPWG2flowCommon loaded..."<<endl;
576 gSystem->Load("libPWG2flowTasks");
577 cerr<<"libPWG2flowTasks loaded..."<<endl;
fdd8c18f 578 }
579
580 else if (mode == mLocalPAR) {
581 //--------------------------------------------------------
582 //If you want to use root and par files from aliroot
583 //--------------------------------------------------------
584 //If you want to use root and par files from aliroot
585 //--------------------------------------------------------
586 SetupPar("STEERBase");
587 SetupPar("ESD");
588 SetupPar("AOD");
589 SetupPar("ANALYSIS");
590 SetupPar("ANALYSISalice");
591 SetupPar("PWG2AOD");
592 SetupPar("CORRFW");
593 SetupPar("PWG2flowCommon");
594 cerr<<"PWG2flowCommon.par loaded..."<<endl;
595 SetupPar("PWG2flowTasks");
596 cerr<<"PWG2flowTasks.par loaded..."<<endl;
597 }
598
599 //---------------------------------------------------------
600 // <<<<<<<<<< Source mode >>>>>>>>>>>>
601 //---------------------------------------------------------
602 else if (mode==mLocalSource) {
603
604 // In root inline compile
605
606
607 // Constants
608 gROOT->LoadMacro("AliFlowCommon/AliFlowCommonConstants.cxx+");
609 gROOT->LoadMacro("AliFlowCommon/AliFlowLYZConstants.cxx+");
fdd8c18f 610
611 // Flow event
612 gROOT->LoadMacro("AliFlowCommon/AliFlowVector.cxx+");
5062cb0f 613 gROOT->LoadMacro("AliFlowCommon/AliFlowTrackSimple.cxx+");
fdd8c18f 614 gROOT->LoadMacro("AliFlowCommon/AliFlowTrackSimpleCuts.cxx+");
5062cb0f 615 gROOT->LoadMacro("AliFlowCommon/AliFlowEventSimple.cxx+");
616
fdd8c18f 617 // Output histosgrams
618 gROOT->LoadMacro("AliFlowCommon/AliFlowCommonHist.cxx+");
619 gROOT->LoadMacro("AliFlowCommon/AliFlowCommonHistResults.cxx+");
620 gROOT->LoadMacro("AliFlowCommon/AliFlowLYZHist1.cxx+");
621 gROOT->LoadMacro("AliFlowCommon/AliFlowLYZHist2.cxx+");
622
623 // Functions needed for various methods
624 gROOT->LoadMacro("AliFlowCommon/AliCumulantsFunctions.cxx+");
fdd8c18f 625 gROOT->LoadMacro("AliFlowCommon/AliFlowLYZEventPlane.cxx+");
626
627 // Flow Analysis code for various methods
628 gROOT->LoadMacro("AliFlowCommon/AliFlowAnalysisWithMCEventPlane.cxx+");
629 gROOT->LoadMacro("AliFlowCommon/AliFlowAnalysisWithScalarProduct.cxx+");
630 gROOT->LoadMacro("AliFlowCommon/AliFlowAnalysisWithLYZEventPlane.cxx+");
631 gROOT->LoadMacro("AliFlowCommon/AliFlowAnalysisWithLeeYangZeros.cxx+");
632 gROOT->LoadMacro("AliFlowCommon/AliFlowAnalysisWithCumulants.cxx+");
633 gROOT->LoadMacro("AliFlowCommon/AliFlowAnalysisWithQCumulants.cxx+");
ecac11c2 634 gROOT->LoadMacro("AliFlowCommon/AliFlowAnalysisWithFittingQDistribution.cxx+");
635 gROOT->LoadMacro("AliFlowCommon/AliFlowAnalysisWithMixedHarmonics.cxx+");
636 gROOT->LoadMacro("AliFlowCommon/AliFlowAnalysisWithNestedLoops.cxx+");
fdd8c18f 637
fdd8c18f 638 cout << "finished loading macros!" << endl;
639
640 }
641
642}
643
644