PWG2/FLOW/macros/runFlowAnalysis.C

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