PWG2/FLOW/macros/runAliAnalysisTaskFlow.C

   1 /////////////////////////////////////////////////////////////////////////////////
   2 //
   3 // HOW TO USE THIS MACRO:
   4 //
   5 // With this macro several flow analysis can be run.
   6 // SP    = Scalar Product                (for PbPb or pp)
   7 // LYZ1  = Lee Yang Zeroes first run     (for PbPb)
   8 // LYZ2  = Lee Yang Zeroes second run    (for PbPb)
   9 // LYZEP = Lee Yang Zeroes Event Plane   (for PbPb)
  10 // GFC   = Generating Function Cumulants (for PbPb)
  11 // QC    = Q-cumulants                   (for PbPb)
  12 // FQD   = Fitting q-distribution        (for PbPb)
  13 // MCEP  = Flow calculated from the real MC event plane (for PbPb only!)
  14 //
  15 // The LYZ analysis should be done in the following order;
  16 // LYZ1 -> LYZ2 -> LYZEP,
  17 // because LYZ2 depends on the outputfile of LYZ1 and LYZEP on the outputfile
  18 // of LYZ2.
  19 //
  20 // The MCEP method is a reference method.
  21 // It can only be run when MC information (kinematics.root & galice.root file) is available
  22 // in which the reaction plane is stored.
  23 //
  24 // One can run on ESD, AOD or MC.
  25 // Additional options are ESDMC0, ESDMC1. In these options the ESD and MC
  26 // information is combined. Tracks are selected in the ESD, the PID information
  27 // is taken from the MC (perfect PID). For ESDMC0 the track kinematics is taken
  28 // from the ESD and for ESDMC1 it is taken from the MC information.
  29 //
  30 ///////////////////////////////////////////////////////////////////////////////////
  31
  32
  33 // SETTING THE ANALYSIS
  34
  35 // Flow analysis methods can be: (set to kTRUE or kFALSE)
  36 Bool_t SP    = kTRUE;
  37 Bool_t LYZ1  = kTRUE;
  38 Bool_t LYZ2  = kFALSE;
  39 Bool_t LYZEP = kFALSE;
  40 Bool_t GFC   = kTRUE;
  41 Bool_t QC    = kTRUE;
  42 Bool_t FQD   = kTRUE;
  43 Bool_t MCEP  = kTRUE;
  44
  45
  46 // Type of analysis can be:
  47 // ESD, AOD, MC, ESDMC0, ESDMC1
  48 const TString type = "ESD";
  49
  50 // Boolean to fill/not fill the QA histograms
  51 Bool_t QA = kFALSE;
  52
  53
  54 // WEIGHTS SETTINGS:
  55 // to use weigths for the Q vector
  56 Bool_t usePhiWeights = kFALSE; //Phi
  57 Bool_t usePtWeights  = kFALSE; //v'(pt)
  58 Bool_t useEtaWeights = kFALSE; //v'(eta)
  59 Bool_t useWeights = usePhiWeights||usePtWeights||usePtWeights;
  60
  61 // SETTING THE CUTS
  62
  63 // For integrated flow
  64 const Double_t ptmin1 = 0.0;
  65 const Double_t ptmax1 = 10.0;
  66 const Double_t ymin1  = -1.;
  67 const Double_t ymax1  = 1.;
  68 const Int_t mintrackrefsTPC1 = 2;
  69 const Int_t mintrackrefsITS1 = 3;
  70 const Int_t charge1 = 1;
  71 Bool_t UsePIDIntegratedFlow = kFALSE;
  72 const Int_t PDG1 = 211;
  73 const Int_t minclustersTPC1 = 50;
  74 const Int_t maxnsigmatovertex1 = 3;
  75
  76 // For differential flow
  77 const Double_t ptmin2 = 0.0;
  78 const Double_t ptmax2 = 10.0;
  79 const Double_t ymin2  = -1.;
  80 const Double_t ymax2  = 1.;
  81 const Int_t mintrackrefsTPC2 = 2;
  82 const Int_t mintrackrefsITS2 = 3;
  83 const Int_t charge2 = 1;
  84 Bool_t UsePIDDifferentialFlow = kFALSE;
  85 const Int_t PDG2 = 211;
  86 const Int_t minclustersTPC2 = 50;
  87 const Int_t maxnsigmatovertex2 = 3;
  88
  89 //void runAliAnalysisTaskFlow(Int_t nRuns = 50, const Char_t* dataDir="/data/alice2/kolk/Therminator_midcentral", Int_t offset = 0)
  90 void runAliAnalysisTaskFlow(Int_t nRuns = -1, const Char_t* dataDir="/Users/snelling/alice_data/Therminator_midcentral", Int_t offset = 0)
  91
  92 {
  93   TStopwatch timer;
  94   timer.Start();
  95
  96   if (LYZ1 && LYZ2) {cout<<"WARNING: you cannot run LYZ1 and LYZ2 at the same time! LYZ2 needs the output from LYZ1."<<endl; exit(); }
  97   if (LYZ2 && LYZEP) {cout<<"WARNING: you cannot run LYZ2 and LYZEP at the same time! LYZEP needs the output from LYZ2."<<endl; exit(); }
  98   if (LYZ1 && LYZEP) {cout<<"WARNING: you cannot run LYZ1 and LYZEP at the same time! LYZEP needs the output from LYZ2."<<endl; exit(); }
  99
 100   // Include path (to find the .h files when compiling)
 101   gSystem->AddIncludePath("-I$ROOTSYS/include") ;
 102
 103   // Load needed root libraries
 104   gSystem->Load("libTree.so");
 105   gSystem->Load("libGeom.so");
 106   gSystem->Load("libVMC.so");
 107   gSystem->Load("libPhysics.so");
 108
 109   // For using root and par files
 110   SetupPar("PWG2flowCommon");
 111   cerr<<"PWG2flowCommon.par loaded..."<<endl;
 112
 113   // Alice specific stuff
 114   SetupPar("STEERBase");
 115   SetupPar("ESD");
 116   SetupPar("AOD");
 117
 118   SetupPar("ANALYSIS");
 119   SetupPar("ANALYSISalice");
 120   SetupPar("PWG2AOD");
 121
 122   SetupPar("CORRFW");
 123
 124   SetupPar("PWG2flowTasks");
 125   cerr<<"PWG2flowTasks.par loaded..."<<endl;
 126
 127   // For using aliroot
 128   /*
 129   //  gSystem->AddIncludePath("-I$ALICE_ROOT/include") ;
 130
 131   gSystem->Load("libESD.so");
 132   cerr<<"libESD loaded..."<<endl;
 133   gSystem->Load("libANALYSIS.so");
 134   cerr<<"libANALYSIS.so loaded..."<<endl;
 135   gSystem->Load("libANALYSISRL.so");
 136   cerr<<"libANALYSISRL.so loaded..."<<endl;
 137   gSystem->Load("libCORRFW.so");
 138   cerr<<"libCORRFW.so loaded..."<<endl;
 139   gSystem->Load("libPWG2flowCommon.so");
 140   cerr<<"libPWG2flowCommon.so loaded..."<<endl;
 141   gSystem->Load("libPWG2flowTasks.so");
 142   cerr<<"libPWG2flowTasks.so loaded..."<<endl;
 143   */
 144
 145   // Create the TChain. CreateESDChain() is defined in CreateESDChain.C
 146   if (type!="AOD") { TChain* chain = CreateESDChain(dataDir, nRuns, offset);}
 147   //  cout<<"chain ("<<chain<<")"<<endl; }
 148   else { TChain* chain = CreateAODChain(dataDir, nRuns, offset);}
 149   //  cout<<"chain ("<<chain<<")"<<endl; }
 150
 151   //____________________________________________//
 152   // Create cuts using correction framework
 153
 154   // Set TList for the QA histograms
 155   if (QA) {
 156     if (SP){
 157       TList* qaIntSP = new TList();
 158       TList* qaDiffSP = new TList(); }
 159     if (LYZ1) {
 160       TList* qaIntLYZ1 = new TList();
 161       TList* qaDiffLYZ1 = new TList(); }
 162     if (LYZ2) {
 163       TList* qaIntLYZ2 = new TList();
 164       TList* qaDiffLYZ2 = new TList(); }
 165     if (LYZEP) {
 166       TList* qaIntLYZEP = new TList();
 167       TList* qaDiffLYZEP = new TList(); }
 168     if (GFC) {
 169       TList* qaIntGFC = new TList();
 170       TList* qaDiffGFC = new TList(); }
 171     if (QC) {
 172       TList* qaIntQC = new TList();
 173       TList* qaDiffQC = new TList(); }
 174     if (FQD) {
 175       TList* qaIntFQD = new TList();
 176       TList* qaDiffFQD = new TList(); }
 177     if (MCEP) {
 178       TList* qaIntMCEP = new TList();
 179       TList* qaDiffMCEP = new TList(); }
 180   }
 181
 182   //############# cuts on MC
 183   AliCFTrackKineCuts* mcKineCuts1 = new AliCFTrackKineCuts("mcKineCuts1","MC-level kinematic cuts");
 184   mcKineCuts1->SetPtRange(ptmin1,ptmax1);
 185   mcKineCuts1->SetRapidityRange(ymin1,ymax1);
 186   mcKineCuts1->SetChargeMC(charge1);
 187   if (QA) {
 188     if (SP)   { mcKineCuts1->SetQAOn(qaIntSP); }
 189     if (LYZ1) { mcKineCuts1->SetQAOn(qaIntLYZ1); }
 190     if (LYZ2) { mcKineCuts1->SetQAOn(qaIntLYZ2); }
 191     if (LYZEP){ mcKineCuts1->SetQAOn(qaIntLYZEP); }
 192     if (GFC)  { mcKineCuts1->SetQAOn(qaIntGFC); }
 193     if (QC)   { mcKineCuts1->SetQAOn(qaIntQC); }
 194     if (FQD)  { mcKineCuts1->SetQAOn(qaIntFQD); }
 195     if (MCEP) { mcKineCuts1->SetQAOn(qaIntMCEP); }
 196   }
 197
 198   AliCFTrackKineCuts* mcKineCuts2 = new AliCFTrackKineCuts("mcKineCuts2","MC-level kinematic cuts");
 199   mcKineCuts2->SetPtRange(ptmin2,ptmax2);
 200   mcKineCuts2->SetRapidityRange(ymin2,ymax2);
 201   mcKineCuts2->SetChargeMC(charge2);
 202   if (QA) {
 203     if (SP)   { mcKineCuts2->SetQAOn(qaDiffSP); }
 204     if (LYZ1) { mcKineCuts2->SetQAOn(qaDiffLYZ1); }
 205     if (LYZ2) { mcKineCuts2->SetQAOn(qaDiffLYZ2); }
 206     if (LYZEP){ mcKineCuts2->SetQAOn(qaDiffLYZEP); }
 207     if (GFC)  { mcKineCuts2->SetQAOn(qaDiffGFC); }
 208     if (QC)   { mcKineCuts2->SetQAOn(qaDiffQC); }
 209     if (FQD)  { mcKineCuts2->SetQAOn(qaDiffFQD); }
 210     if (MCEP) { mcKineCuts2->SetQAOn(qaDiffMCEP); }
 211   }
 212
 213   AliCFParticleGenCuts* mcGenCuts1 = new AliCFParticleGenCuts("mcGenCuts1","MC particle generation cuts for integrated flow");
 214   mcGenCuts1->SetRequireIsPrimary();
 215   if (UsePIDIntegratedFlow) {mcGenCuts1->SetRequirePdgCode(PDG1);}
 216   if (QA) {
 217     if (SP)   { mcGenCuts1->SetQAOn(qaIntSP); }
 218     if (LYZ1) { mcGenCuts1->SetQAOn(qaIntLYZ1); }
 219     if (LYZ2) { mcGenCuts1->SetQAOn(qaIntLYZ2); }
 220     if (LYZEP){ mcGenCuts1->SetQAOn(qaIntLYZEP); }
 221     if (GFC)  { mcGenCuts1->SetQAOn(qaIntGFC); }
 222     if (QC)   { mcGenCuts1->SetQAOn(qaIntQC); }
 223     if (FQD)  { mcGenCuts1->SetQAOn(qaIntFQD); }
 224     if (MCEP) { mcGenCuts1->SetQAOn(qaIntMCEP); }
 225   }
 226
 227   AliCFParticleGenCuts* mcGenCuts2 = new AliCFParticleGenCuts("mcGenCuts2","MC particle generation cuts for differential flow");
 228   mcGenCuts2->SetRequireIsPrimary();
 229   if (UsePIDDifferentialFlow) {mcGenCuts2->SetRequirePdgCode(PDG2);}
 230   if (QA) {
 231     if (SP)   { mcGenCuts2->SetQAOn(qaDiffSP); }
 232     if (LYZ1) { mcGenCuts2->SetQAOn(qaDiffLYZ1); }
 233     if (LYZ2) { mcGenCuts2->SetQAOn(qaDiffLYZ2); }
 234     if (LYZEP){ mcGenCuts2->SetQAOn(qaDiffLYZEP); }
 235     if (GFC)  { mcGenCuts2->SetQAOn(qaDiffGFC); }
 236     if (QC)   { mcGenCuts2->SetQAOn(qaDiffQC); }
 237     if (FQD)  { mcGenCuts2->SetQAOn(qaDiffFQD); }
 238     if (MCEP) { mcGenCuts2->SetQAOn(qaDiffMCEP); }
 239   }
 240
 241   //############# Acceptance Cuts
 242   AliCFAcceptanceCuts *mcAccCuts1 = new AliCFAcceptanceCuts("mcAccCuts1","MC acceptance cuts");
 243   mcAccCuts1->SetMinNHitITS(mintrackrefsITS1);
 244   mcAccCuts1->SetMinNHitTPC(mintrackrefsTPC1);
 245   if (QA) {
 246     if (SP)   { mcAccCuts1->SetQAOn(qaIntSP); }
 247     if (LYZ1) { mcAccCuts1->SetQAOn(qaIntLYZ1); }
 248     if (LYZ2) { mcAccCuts1->SetQAOn(qaIntLYZ2); }
 249     if (LYZEP){ mcAccCuts1->SetQAOn(qaIntLYZEP); }
 250     if (GFC)  { mcAccCuts1->SetQAOn(qaIntGFC); }
 251     if (QC)   { mcAccCuts1->SetQAOn(qaIntQC); }
 252     if (FQD)  { mcAccCuts1->SetQAOn(qaIntFQD); }
 253     if (MCEP) { mcAccCuts1->SetQAOn(qaIntMCEP); }
 254   }
 255
 256   AliCFAcceptanceCuts *mcAccCuts2 = new AliCFAcceptanceCuts("mcAccCuts2","MC acceptance cuts");
 257   mcAccCuts2->SetMinNHitITS(mintrackrefsITS2);
 258   mcAccCuts2->SetMinNHitTPC(mintrackrefsTPC2);
 259   if (QA) {
 260     if (SP)   { mcAccCuts2->SetQAOn(qaDiffSP); }
 261     if (LYZ1) { mcAccCuts2->SetQAOn(qaDiffLYZ1); }
 262     if (LYZ2) { mcAccCuts2->SetQAOn(qaDiffLYZ2); }
 263     if (LYZEP){ mcAccCuts2->SetQAOn(qaDiffLYZEP); }
 264     if (GFC)  { mcAccCuts2->SetQAOn(qaDiffGFC); }
 265     if (QC)   { mcAccCuts2->SetQAOn(qaDiffQC); }
 266     if (FQD)  { mcAccCuts2->SetQAOn(qaDiffFQD); }
 267     if (MCEP) { mcAccCuts2->SetQAOn(qaDiffMCEP); }
 268   }
 269
 270   //############# Rec-Level kinematic cuts
 271   AliCFTrackKineCuts *recKineCuts1 = new AliCFTrackKineCuts("recKineCuts1","rec-level kine cuts");
 272   recKineCuts1->SetPtRange(ptmin1,ptmax1);
 273   recKineCuts1->SetRapidityRange(ymin1,ymax1);
 274   recKineCuts1->SetChargeRec(charge1);
 275   if (QA) {
 276     if (SP)   { recKineCuts1->SetQAOn(qaIntSP); }
 277     if (LYZ1) { recKineCuts1->SetQAOn(qaIntLYZ1); }
 278     if (LYZ2) { recKineCuts1->SetQAOn(qaIntLYZ2); }
 279     if (LYZEP){ recKineCuts1->SetQAOn(qaIntLYZEP); }
 280     if (GFC)  { recKineCuts1->SetQAOn(qaIntGFC); }
 281     if (QC)   { recKineCuts1->SetQAOn(qaIntQC); }
 282     if (FQD)  { recKineCuts1->SetQAOn(qaIntFQD); }
 283     if (MCEP) { recKineCuts1->SetQAOn(qaIntMCEP); }
 284   }
 285
 286   AliCFTrackKineCuts *recKineCuts2 = new AliCFTrackKineCuts("recKineCuts2","rec-level kine cuts");
 287   recKineCuts2->SetPtRange(ptmin2,ptmax2);
 288   recKineCuts2->SetRapidityRange(ymin2,ymax2);
 289   recKineCuts2->SetChargeRec(charge2);
 290   if (QA) {
 291     if (SP)   { recKineCuts2->SetQAOn(qaDiffSP); }
 292     if (LYZ1) { recKineCuts2->SetQAOn(qaDiffLYZ1); }
 293     if (LYZ2) { recKineCuts2->SetQAOn(qaDiffLYZ2); }
 294     if (LYZEP){ recKineCuts2->SetQAOn(qaDiffLYZEP); }
 295     if (GFC)  { recKineCuts2->SetQAOn(qaDiffGFC); }
 296     if (QC)   { recKineCuts2->SetQAOn(qaDiffQC); }
 297     if (FQD)  { recKineCuts2->SetQAOn(qaDiffFQD); }
 298     if (MCEP) { recKineCuts2->SetQAOn(qaDiffMCEP); }
 299   }
 300
 301   AliCFTrackQualityCuts *recQualityCuts1 = new AliCFTrackQualityCuts("recQualityCuts1","rec-level quality cuts");
 302   recQualityCuts1->SetMinNClusterTPC(minclustersTPC1);
 303   recQualityCuts1->SetStatus(AliESDtrack::kITSrefit);
 304   if (QA) {
 305     if (SP)   { recQualityCuts1->SetQAOn(qaIntSP); }
 306     if (LYZ1) { recQualityCuts1->SetQAOn(qaIntLYZ1); }
 307     if (LYZ2) { recQualityCuts1->SetQAOn(qaIntLYZ2); }
 308     if (LYZEP){ recQualityCuts1->SetQAOn(qaIntLYZEP); }
 309     if (GFC)  { recQualityCuts1->SetQAOn(qaIntGFC); }
 310     if (QC)   { recQualityCuts1->SetQAOn(qaIntQC); }
 311     if (FQD)  { recQualityCuts1->SetQAOn(qaIntFQD); }
 312     if (MCEP) { recQualityCuts1->SetQAOn(qaIntMCEP); }
 313   }
 314
 315   AliCFTrackQualityCuts *recQualityCuts2 = new AliCFTrackQualityCuts("recQualityCuts2","rec-level quality cuts");
 316   recQualityCuts2->SetMinNClusterTPC(minclustersTPC2);
 317   recQualityCuts2->SetStatus(AliESDtrack::kITSrefit);
 318   if (QA) {
 319     if (SP)   { recQualityCuts2->SetQAOn(qaDiffSP); }
 320     if (LYZ1) { recQualityCuts2->SetQAOn(qaDiffLYZ1); }
 321     if (LYZ2) { recQualityCuts2->SetQAOn(qaDiffLYZ2); }
 322     if (LYZEP){ recQualityCuts2->SetQAOn(qaDiffLYZEP); }
 323     if (GFC)  { recQualityCuts2->SetQAOn(qaDiffGFC); }
 324     if (QC)   { recQualityCuts2->SetQAOn(qaDiffQC); }
 325     if (FQD)  { recQualityCuts2->SetQAOn(qaDiffFQD); }
 326     if (MCEP) { recQualityCuts2->SetQAOn(qaDiffMCEP); }
 327   }
 328
 329   AliCFTrackIsPrimaryCuts *recIsPrimaryCuts1 = new AliCFTrackIsPrimaryCuts("recIsPrimaryCuts1","rec-level isPrimary cuts");
 330   recIsPrimaryCuts1->SetMaxNSigmaToVertex(maxnsigmatovertex1);
 331   if (QA) {
 332     if (SP)   { recIsPrimaryCuts1->SetQAOn(qaIntSP); }
 333     if (LYZ1) { recIsPrimaryCuts1->SetQAOn(qaIntLYZ1); }
 334     if (LYZ2) { recIsPrimaryCuts1->SetQAOn(qaIntLYZ2); }
 335     if (LYZEP){ recIsPrimaryCuts1->SetQAOn(qaIntLYZEP); }
 336     if (GFC)  { recIsPrimaryCuts1->SetQAOn(qaIntGFC); }
 337     if (QC)   { recIsPrimaryCuts1->SetQAOn(qaIntQC); }
 338     if (FQD)  { recIsPrimaryCuts1->SetQAOn(qaIntFQD); }
 339     if (MCEP) { recIsPrimaryCuts1->SetQAOn(qaIntMCEP); }
 340   }
 341
 342   AliCFTrackIsPrimaryCuts *recIsPrimaryCuts2 = new AliCFTrackIsPrimaryCuts("recIsPrimaryCuts2","rec-level isPrimary cuts");
 343   recIsPrimaryCuts2->SetMaxNSigmaToVertex(maxnsigmatovertex2);
 344   if (QA) {
 345     if (SP)   { recIsPrimaryCuts2->SetQAOn(qaDiffSP); }
 346     if (LYZ1) { recIsPrimaryCuts2->SetQAOn(qaDiffLYZ1); }
 347     if (LYZ2) { recIsPrimaryCuts2->SetQAOn(qaDiffLYZ2); }
 348     if (LYZEP){ recIsPrimaryCuts2->SetQAOn(qaDiffLYZEP); }
 349     if (GFC)  { recIsPrimaryCuts2->SetQAOn(qaDiffGFC); }
 350     if (QC)   { recIsPrimaryCuts2->SetQAOn(qaDiffQC); }
 351     if (FQD)  { recIsPrimaryCuts2->SetQAOn(qaDiffFQD); }
 352     if (MCEP) { recIsPrimaryCuts2->SetQAOn(qaDiffMCEP); }
 353   }
 354
 355
 356
 357   int n_species = AliPID::kSPECIES ;
 358   Double_t* prior = new Double_t[n_species];
 359
 360   prior[0] = 0.0244519 ;
 361   prior[1] = 0.0143988 ;
 362   prior[2] = 0.805747  ;
 363   prior[3] = 0.0928785 ;
 364   prior[4] = 0.0625243 ;
 365
 366   AliCFTrackCutPid* cutPID1 = NULL;
 367   if(UsePIDIntegratedFlow) {
 368     AliCFTrackCutPid* cutPID1 = new AliCFTrackCutPid("cutPID1","ESD_PID for integrated flow") ;
 369     cutPID1->SetPriors(prior);
 370     cutPID1->SetProbabilityCut(0.0);
 371     cutPID1->SetDetectors("TPC TOF");
 372     switch(TMath::Abs(PDG1)) {
 373     case 11   : cutPID1->SetParticleType(AliPID::kElectron, kTRUE); break;
 374     case 13   : cutPID1->SetParticleType(AliPID::kMuon    , kTRUE); break;
 375     case 211  : cutPID1->SetParticleType(AliPID::kPion    , kTRUE); break;
 376     case 321  : cutPID1->SetParticleType(AliPID::kKaon    , kTRUE); break;
 377     case 2212 : cutPID1->SetParticleType(AliPID::kProton  , kTRUE); break;
 378     default   : printf("UNDEFINED PID\n"); break;
 379     }
 380     if (QA) {
 381       if (SP) {cutPID1->SetQAOn(qaIntSP);}
 382       if (LYZ1) {cutPID1->SetQAOn(qaIntLYZ1);}
 383       if (LYZ2) {cutPID1->SetQAOn(qaIntLYZ2);}
 384       if (LYZEP){cutPID1->SetQAOn(qaIntLYZEP);}
 385       if (GFC) {cutPID1->SetQAOn(qaIntGFC);}
 386       if (QC) {cutPID1->SetQAOn(qaIntQC);}
 387       if (FQD) {cutPID1->SetQAOn(qaIntFQD);}
 388       if (MCEP) {cutPID1->SetQAOn(qaIntMCEP);}
 389     }
 390   }
 391
 392   AliCFTrackCutPid* cutPID2 = NULL;
 393   if (UsePIDDifferentialFlow) {
 394     AliCFTrackCutPid* cutPID2 = new AliCFTrackCutPid("cutPID2","ESD_PID for differential flow") ;
 395     cutPID2->SetPriors(prior);
 396     cutPID2->SetProbabilityCut(0.0);
 397     cutPID2->SetDetectors("TPC TOF");
 398     switch(TMath::Abs(PDG2)) {
 399     case 11   : cutPID2->SetParticleType(AliPID::kElectron, kTRUE); break;
 400     case 13   : cutPID2->SetParticleType(AliPID::kMuon    , kTRUE); break;
 401     case 211  : cutPID2->SetParticleType(AliPID::kPion    , kTRUE); break;
 402     case 321  : cutPID2->SetParticleType(AliPID::kKaon    , kTRUE); break;
 403     case 2212 : cutPID2->SetParticleType(AliPID::kProton  , kTRUE); break;
 404     default   : printf("UNDEFINED PID\n"); break;
 405     }
 406     if (QA) {
 407       if (SP) {cutPID2->SetQAOn(qaIntSP);}
 408       if (LYZ1) {cutPID2->SetQAOn(qaIntLYZ1);}
 409       if (LYZ2) {cutPID2->SetQAOn(qaIntLYZ2);}
 410       if (LYZEP){cutPID2->SetQAOn(qaIntLYZEP);}
 411       if (GFC) {cutPID2->SetQAOn(qaIntGFC);}
 412       if (QC) {cutPID2->SetQAOn(qaIntQC);}
 413       if (FQD) {cutPID2->SetQAOn(qaIntFQD);}
 414       if (MCEP) {cutPID2->SetQAOn(qaIntMCEP);}
 415     }
 416   }
 417
 418   printf("CREATE MC KINE CUTS\n");
 419   TObjArray* mcList1 = new TObjArray(0);
 420   mcList1->AddLast(mcKineCuts1);
 421   mcList1->AddLast(mcGenCuts1);
 422
 423   TObjArray* mcList2 = new TObjArray(0);
 424   mcList2->AddLast(mcKineCuts2);
 425   mcList2->AddLast(mcGenCuts2);
 426
 427   printf("CREATE ACCEPTANCE CUTS\n");
 428   TObjArray* accList1 = new TObjArray(0) ;
 429   accList1->AddLast(mcAccCuts1);
 430
 431   TObjArray* accList2 = new TObjArray(0) ;
 432   accList2->AddLast(mcAccCuts2);
 433
 434   printf("CREATE RECONSTRUCTION CUTS\n");
 435   TObjArray* recList1 = new TObjArray(0) ;
 436   recList1->AddLast(recKineCuts1);
 437   recList1->AddLast(recQualityCuts1);
 438   recList1->AddLast(recIsPrimaryCuts1);
 439
 440   TObjArray* recList2 = new TObjArray(0) ;
 441   recList2->AddLast(recKineCuts2);
 442   recList2->AddLast(recQualityCuts2);
 443   recList2->AddLast(recIsPrimaryCuts2);
 444
 445   printf("CREATE PID CUTS\n");
 446   TObjArray* fPIDCutList1 = new TObjArray(0) ;
 447   if(UsePIDIntegratedFlow) {fPIDCutList1->AddLast(cutPID1);}
 448
 449   TObjArray* fPIDCutList2 = new TObjArray(0) ;
 450   if (UsePIDDifferentialFlow)  {fPIDCutList2->AddLast(cutPID2);}
 451
 452   printf("CREATE INTERFACE AND CUTS FOR INTEGRATED FLOW\n");
 453   AliCFManager* cfmgr1 = new AliCFManager();
 454   //  cfmgr1->SetNStepParticle(4); // needed for trunk not in 4.16
 455   cfmgr1->SetParticleCutsList(AliCFManager::kPartGenCuts,mcList1);       //on MC
 456   cfmgr1->SetParticleCutsList(AliCFManager::kPartAccCuts,accList1);      //on MC
 457   cfmgr1->SetParticleCutsList(AliCFManager::kPartRecCuts,recList1);      //on ESD
 458   cfmgr1->SetParticleCutsList(AliCFManager::kPartSelCuts,fPIDCutList1);  //on ESD
 459
 460   printf("CREATE INTERFACE AND CUTS FOR DIFFERENTIAL FLOW\n");
 461   AliCFManager* cfmgr2 = new AliCFManager();
 462   //  cfmgr2->SetNStepParticle(4); // needed for trunk not in 4.16
 463   cfmgr2->SetParticleCutsList(AliCFManager::kPartGenCuts,mcList2);
 464   cfmgr2->SetParticleCutsList(AliCFManager::kPartAccCuts,accList2);
 465   cfmgr2->SetParticleCutsList(AliCFManager::kPartRecCuts,recList2);
 466   cfmgr2->SetParticleCutsList(AliCFManager::kPartSelCuts,fPIDCutList2);
 467
 468   //weights:
 469   TFile *weightsFile = NULL;
 470   TList *weightsList = NULL;
 471
 472   if(useWeights)
 473   {
 474    //open the file with the weights:
 475    weightsFile = TFile::Open("weights.root","READ");
 476    if(weightsFile) {
 477      //access the list which holds the histos with weigths:
 478      weightsList = (TList*)weightsFile->Get("weights");
 479    }
 480    else {
 481      cout<<" WARNING: the file <weights.root> with weights from the previous run was not available."<<endl;
 482      break;
 483    }
 484   }//end of if(useWeights)
 485
 486   if (LYZ2){
 487     // read the input file from the first run
 488     TString inputFileNameLYZ2 = "outputLYZ1analysis" ;
 489     inputFileNameLYZ2 += type;
 490     inputFileNameLYZ2 += ".root";
 491     cout<<"The input file is "<<inputFileNameLYZ2.Data()<<endl;
 492     TFile* fInputFileLYZ2 = new TFile(inputFileNameLYZ2.Data(),"READ");
 493     if(!fInputFileLYZ2 || fInputFileLYZ2->IsZombie()) {
 494       cerr << " ERROR: NO First Run file... " << endl ; }
 495     else {
 496       TList* fInputListLYZ2 = (TList*)fInputFileLYZ2->Get("cobjLYZ1");
 497       if (!fInputListLYZ2) {cout<<"list is NULL pointer!"<<endl;}
 498     }
 499     cout<<"LYZ2 input file/list read..."<<endl;
 500   }
 501
 502   if (LYZEP) {
 503     // read the input file from the second LYZ run
 504     TString inputFileNameLYZEP = "outputLYZ2analysis" ;
 505     inputFileNameLYZEP += type;
 506     inputFileNameLYZEP += ".root";
 507     cout<<"The input file is "<<inputFileNameLYZEP.Data()<<endl;
 508     TFile* fInputFileLYZEP = new TFile(inputFileNameLYZEP.Data(),"READ");
 509     if(!fInputFileLYZEP || fInputFileLYZEP->IsZombie()) { cerr << " ERROR: NO First Run file... " << endl ; }
 510     else {
 511       TList* fInputListLYZEP = (TList*)fInputFileLYZEP->Get("cobjLYZ2");
 512       if (!fInputListLYZEP) {cout<<"list is NULL pointer!"<<endl;}
 513     }
 514     cout<<"LYZEP input file/list read..."<<endl;
 515   }
 516
 517
 518   //____________________________________________//
 519   // Make the analysis manager
 520   AliAnalysisManager *mgr = new AliAnalysisManager("FlowAnalysisManager");
 521
 522   if (type == "ESD"){
 523     AliVEventHandler* esdH = new AliESDInputHandler;
 524     mgr->SetInputEventHandler(esdH);
 525     if (MCEP) {
 526       AliMCEventHandler *mc = new AliMCEventHandler();
 527       mgr->SetMCtruthEventHandler(mc);}  }
 528
 529   if (type == "AOD"){
 530     AliVEventHandler* aodH = new AliAODInputHandler;
 531     mgr->SetInputEventHandler(aodH);
 532     if (MCEP) {
 533       AliMCEventHandler *mc = new AliMCEventHandler();
 534       mgr->SetMCtruthEventHandler(mc);}  }
 535
 536   if (type == "MC" || type == "ESDMC0" || type == "ESDMC1"){
 537     AliVEventHandler* esdH = new AliESDInputHandler;
 538     mgr->SetInputEventHandler(esdH);
 539
 540     AliMCEventHandler *mc = new AliMCEventHandler();
 541     mgr->SetMCtruthEventHandler(mc); }
 542
 543   //____________________________________________//
 544   // Task
 545
 546   if (SP){
 547     if (QA) { AliAnalysisTaskScalarProduct *taskSP = new AliAnalysisTaskScalarProduct("TaskScalarProduct",kTRUE); }
 548     else { AliAnalysisTaskScalarProduct *taskSP = new AliAnalysisTaskScalarProduct("TaskScalarProduct",kFALSE); }
 549     taskSP->SetAnalysisType(type);
 550     taskSP->SetCFManager1(cfmgr1);
 551     taskSP->SetCFManager2(cfmgr2);
 552     if (QA) {
 553       taskSP->SetQAList1(qaIntSP);
 554       taskSP->SetQAList2(qaDiffSP); }
 555     mgr->AddTask(taskSP);
 556   }
 557   if (LYZ1){
 558     if (QA) { AliAnalysisTaskLeeYangZeros *taskLYZ1 = new AliAnalysisTaskLeeYangZeros("TaskLeeYangZeros",kTRUE,kTRUE);}
 559     else { AliAnalysisTaskLeeYangZeros *taskLYZ1 = new AliAnalysisTaskLeeYangZeros("TaskLeeYangZeros",kTRUE,kFALSE);}
 560     taskLYZ1->SetAnalysisType(type);
 561     taskLYZ1->SetFirstRunLYZ(kTRUE);
 562     taskLYZ1->SetUseSumLYZ(kTRUE);
 563     taskLYZ1->SetCFManager1(cfmgr1);
 564     taskLYZ1->SetCFManager2(cfmgr2);
 565     if (QA) {
 566       taskLYZ1->SetQAList1(qaIntLYZ1);
 567       taskLYZ1->SetQAList2(qaDiffLYZ1);}
 568     mgr->AddTask(taskLYZ1);
 569   }
 570   if (LYZ2){
 571     if (QA) { AliAnalysisTaskLeeYangZeros *taskLYZ2 = new AliAnalysisTaskLeeYangZeros("TaskLeeYangZeros",kFALSE,kTRUE);}
 572     else { AliAnalysisTaskLeeYangZeros *taskLYZ2 = new AliAnalysisTaskLeeYangZeros("TaskLeeYangZeros",kFALSE,kFALSE); }
 573     taskLYZ2->SetAnalysisType(type);
 574     taskLYZ2->SetFirstRunLYZ(kFALSE);
 575     taskLYZ2->SetUseSumLYZ(kTRUE);
 576     taskLYZ2->SetCFManager1(cfmgr1);
 577     taskLYZ2->SetCFManager2(cfmgr2);
 578     if (QA) {
 579       taskLYZ2->SetQAList1(qaIntLYZ2);
 580       taskLYZ2->SetQAList2(qaDiffLYZ2); }
 581     mgr->AddTask(taskLYZ2);
 582   }
 583   if (LYZEP){
 584     if (QA) { AliAnalysisTaskLYZEventPlane *taskLYZEP = new AliAnalysisTaskLYZEventPlane("TaskLYZEventPlane",kTRUE); }
 585     else { AliAnalysisTaskLYZEventPlane *taskLYZEP = new AliAnalysisTaskLYZEventPlane("TaskLYZEventPlane",kFALSE); }
 586     taskLYZEP->SetAnalysisType(type);
 587     taskLYZEP->SetCFManager1(cfmgr1);
 588     taskLYZEP->SetCFManager2(cfmgr2);
 589     if (QA) {
 590       taskLYZEP->SetQAList1(qaIntLYZEP);
 591       taskLYZEP->SetQAList2(qaDiffLYZEP); }
 592     mgr->AddTask(taskLYZEP);
 593   }
 594   if (GFC){
 595     if (QA) { AliAnalysisTaskCumulants *taskGFC = new AliAnalysisTaskCumulants("TaskCumulants",kTRUE);}
 596     else { AliAnalysisTaskCumulants *taskGFC = new AliAnalysisTaskCumulants("TaskCumulants",kFALSE);}
 597     taskGFC->SetAnalysisType(type);
 598     taskGFC->SetCFManager1(cfmgr1);
 599     taskGFC->SetCFManager2(cfmgr2);
 600     if (QA) {
 601       taskGFC->SetQAList1(qaIntGFC);
 602       taskGFC->SetQAList2(qaDiffGFC); }
 603     mgr->AddTask(taskGFC);
 604   }
 605   if (QC){
 606     if (QA) { AliAnalysisTaskQCumulants *taskQC = new AliAnalysisTaskQCumulants("TaskQCumulants",kTRUE,useWeights);}
 607     else { AliAnalysisTaskQCumulants *taskQC = new AliAnalysisTaskQCumulants("TaskQCumulants",kFALSE,useWeights);}
 608     taskQC->SetAnalysisType(type);
 609     taskQC->SetUsePhiWeights(usePhiWeights);
 610     taskQC->SetUsePtWeights(usePtWeights);
 611     taskQC->SetUseEtaWeights(useEtaWeights);
 612     taskQC->SetCFManager1(cfmgr1);
 613     taskQC->SetCFManager2(cfmgr2);
 614     if (QA) {
 615       taskQC->SetQAList1(qaIntQC);
 616       taskQC->SetQAList2(qaDiffQC); }
 617     mgr->AddTask(taskQC);
 618   }
 619   if (FQD){
 620     if (QA) { AliAnalysisTaskFittingQDistribution *taskFQD = new AliAnalysisTaskFittingQDistribution("TaskFittingQDistribution",kTRUE);}
 621     else { AliAnalysisTaskFittingQDistribution *taskFQD = new AliAnalysisTaskFittingQDistribution("TaskFittingQDistribution",kFALSE);}
 622     taskFQD->SetAnalysisType(type);
 623     taskFQD->SetCFManager1(cfmgr1);
 624     taskFQD->SetCFManager2(cfmgr2);
 625     if (QA) {
 626       taskFQD->SetQAList1(qaIntFQD);
 627       taskFQD->SetQAList2(qaDiffFQD); }
 628     mgr->AddTask(taskFQD);
 629   }
 630   if (MCEP){
 631     if (QA) { AliAnalysisTaskMCEventPlane *taskMCEP = new AliAnalysisTaskMCEventPlane("TaskMCEventPlane",kTRUE);}
 632     else { AliAnalysisTaskMCEventPlane *taskMCEP = new AliAnalysisTaskMCEventPlane("TaskMCEventPlane",kFALSE);}
 633     taskMCEP->SetAnalysisType(type);
 634     taskMCEP->SetCFManager1(cfmgr1);
 635     taskMCEP->SetCFManager2(cfmgr2);
 636     if (QA) {
 637       taskMCEP->SetQAList1(qaIntMCEP);
 638       taskMCEP->SetQAList2(qaDiffMCEP); }
 639     mgr->AddTask(taskMCEP);
 640   }
 641
 642
 643   // Create containers for input/output
 644   AliAnalysisDataContainer *cinput1 =
 645     mgr->CreateContainer("cchain1",TChain::Class(),AliAnalysisManager::kInputContainer);
 646
 647   //AliAnalysisDataContainer *cinputWeights = NULL;
 648   if(useWeights)
 649   {
 650    cinputWeights = mgr->CreateContainer("cobjWeights",TList::Class(),AliAnalysisManager::kInputContainer);
 651   }
 652
 653   if (LYZ2){
 654     AliAnalysisDataContainer *cinputLYZ2 =
 655       mgr->CreateContainer("cobjLYZ2in",TList::Class(),AliAnalysisManager::kInputContainer); }
 656   if (LYZEP){
 657     AliAnalysisDataContainer *cinputLYZEP =
 658       mgr->CreateContainer("cobjLYZEPin",TList::Class(),AliAnalysisManager::kInputContainer); }
 659
 660   if(SP) {
 661     TString outputSP = "outputSPanalysis";
 662     outputSP+= type;
 663     outputSP+= ".root";
 664     AliAnalysisDataContainer *coutputSP =
 665       mgr->CreateContainer("cobjSP", TList::Class(),AliAnalysisManager::kOutputContainer,outputSP);
 666   }
 667
 668   if(LYZ1) {
 669     TString outputLYZ1 = "outputLYZ1analysis";
 670     outputLYZ1+= type;
 671     outputLYZ1+= ".root";
 672     AliAnalysisDataContainer *coutputLYZ1 =
 673       mgr->CreateContainer("cobjLYZ1", TList::Class(),AliAnalysisManager::kOutputContainer,outputLYZ1);
 674   }
 675
 676   if(LYZ2) {
 677     TString outputLYZ2 = "outputLYZ2analysis";
 678     outputLYZ2+= type;
 679     outputLYZ2+= ".root";
 680     AliAnalysisDataContainer *coutputLYZ2 =
 681       mgr->CreateContainer("cobjLYZ2", TList::Class(),AliAnalysisManager::kOutputContainer,outputLYZ2);
 682   }
 683
 684   if(LYZEP) {
 685     TString outputLYZEP = "outputLYZEPanalysis";
 686     outputLYZEP+= type;
 687     outputLYZEP+= ".root";
 688     AliAnalysisDataContainer *coutputLYZEP =
 689       mgr->CreateContainer("cobjLYZEP", TList::Class(),AliAnalysisManager::kOutputContainer,outputLYZEP);
 690   }
 691
 692   if(GFC) {
 693     TString outputGFC = "outputGFCanalysis";
 694     outputGFC+= type;
 695     outputGFC+= ".root";
 696     AliAnalysisDataContainer *coutputGFC =
 697       mgr->CreateContainer("cobjGFC", TList::Class(),AliAnalysisManager::kOutputContainer,outputGFC);
 698   }
 699
 700   if(QC) {
 701     TString outputQC = "outputQCanalysis";
 702     outputQC+= type;
 703     outputQC+= ".root";
 704     AliAnalysisDataContainer *coutputQC =
 705       mgr->CreateContainer("cobjQC", TList::Class(),AliAnalysisManager::kOutputContainer,outputQC);
 706   }
 707
 708   if(FQD) {
 709     TString outputFQD = "outputFQDanalysis";
 710     outputFQD+= type;
 711     outputFQD+= ".root";
 712     AliAnalysisDataContainer *coutputFQD =
 713       mgr->CreateContainer("cobjFQD", TList::Class(),AliAnalysisManager::kOutputContainer,outputFQD);
 714   }
 715
 716   if(MCEP) {
 717     TString outputMCEP = "outputMCEPanalysis";
 718     outputMCEP+= type;
 719     outputMCEP+= ".root";
 720     AliAnalysisDataContainer *coutputMCEP =
 721       mgr->CreateContainer("cobjMCEP", TList::Class(),AliAnalysisManager::kOutputContainer,outputMCEP);
 722   }
 723
 724   if (QA) {
 725     if(SP) {
 726       TString qaNameIntSP = "QAforInt_SP_";
 727       qaNameIntSP += type;
 728       qaNameIntSP += ".root";
 729       AliAnalysisDataContainer *coutputQA1SP =
 730         mgr->CreateContainer("QAintSP", TList::Class(),AliAnalysisManager::kOutputContainer,qaNameIntSP);
 731
 732       TString qaNameDiffSP = "QAforDiff_SP_";
 733       qaNameDiffSP += type;
 734       qaNameDiffSP += ".root";
 735       AliAnalysisDataContainer *coutputQA2SP =
 736         mgr->CreateContainer("QAdiffSP", TList::Class(),AliAnalysisManager::kOutputContainer,qaNameDiffSP);
 737     }
 738     if(LYZ1) {
 739       TString qaNameIntLYZ1 = "QAforInt_LYZ1_";
 740       qaNameIntLYZ1 += type;
 741       qaNameIntLYZ1 += ".root";
 742       AliAnalysisDataContainer *coutputQA1LYZ1 =
 743         mgr->CreateContainer("QAintLYZ1", TList::Class(),AliAnalysisManager::kOutputContainer,qaNameIntLYZ1);
 744
 745       TString qaNameDiffLYZ1 = "QAforDiff_LYZ1_";
 746       qaNameDiffLYZ1 += type;
 747       qaNameDiffLYZ1 += ".root";
 748       AliAnalysisDataContainer *coutputQA2LYZ1 =
 749         mgr->CreateContainer("QAdiffLYZ1", TList::Class(),AliAnalysisManager::kOutputContainer,qaNameDiffLYZ1);
 750     }
 751     if(LYZ2) {
 752       TString qaNameIntLYZ2 = "QAforInt_LYZ2_";
 753       qaNameIntLYZ2 += type;
 754       qaNameIntLYZ2 += ".root";
 755       AliAnalysisDataContainer *coutputQA1LYZ2 =
 756         mgr->CreateContainer("QAintLYZ2", TList::Class(),AliAnalysisManager::kOutputContainer,qaNameIntLYZ2);
 757
 758       TString qaNameDiffLYZ2 = "QAforDiff_LYZ2_";
 759       qaNameDiffLYZ2 += type;
 760       qaNameDiffLYZ2 += ".root";
 761       AliAnalysisDataContainer *coutputQA2LYZ2 =
 762         mgr->CreateContainer("QAdiffLYZ2", TList::Class(),AliAnalysisManager::kOutputContainer,qaNameDiffLYZ2);
 763     }
 764     if(LYZEP) {
 765       TString qaNameIntLYZEP = "QAforInt_LYZEP_";
 766       qaNameIntLYZEP += type;
 767       qaNameIntLYZEP += ".root";
 768       AliAnalysisDataContainer *coutputQA1LYZEP =
 769         mgr->CreateContainer("QAintLYZEP", TList::Class(),AliAnalysisManager::kOutputContainer,qaNameIntLYZEP);
 770
 771       TString qaNameDiffLYZEP = "QAforDiff_LYZEP_";
 772       qaNameDiffLYZEP += type;
 773       qaNameDiffLYZEP += ".root";
 774       AliAnalysisDataContainer *coutputQA2LYZEP =
 775         mgr->CreateContainer("QAdiffLYZEP", TList::Class(),AliAnalysisManager::kOutputContainer,qaNameDiffLYZEP);
 776     }
 777     if(GFC) {
 778       TString qaNameIntGFC = "QAforInt_GFC_";
 779       qaNameIntGFC += type;
 780       qaNameIntGFC += ".root";
 781       AliAnalysisDataContainer *coutputQA1GFC =
 782         mgr->CreateContainer("QAintGFC", TList::Class(),AliAnalysisManager::kOutputContainer,qaNameIntGFC);
 783
 784       TString qaNameDiffGFC = "QAforDiff_GFC_";
 785       qaNameDiffGFC += type;
 786       qaNameDiffGFC += ".root";
 787       AliAnalysisDataContainer *coutputQA2GFC =
 788         mgr->CreateContainer("QAdiffGFC", TList::Class(),AliAnalysisManager::kOutputContainer,qaNameDiffGFC);
 789     }
 790     if(QC) {
 791       TString qaNameIntQC = "QAforInt_QC_";
 792       qaNameIntQC += type;
 793       qaNameIntQC += ".root";
 794       AliAnalysisDataContainer *coutputQA1QC =
 795         mgr->CreateContainer("QAintQC", TList::Class(),AliAnalysisManager::kOutputContainer,qaNameIntQC);
 796
 797       TString qaNameDiffQC = "QAforDiff_QC_";
 798       qaNameDiffQC += type;
 799       qaNameDiffQC += ".root";
 800       AliAnalysisDataContainer *coutputQA2QC =
 801         mgr->CreateContainer("QAdiffQC", TList::Class(),AliAnalysisManager::kOutputContainer,qaNameDiffQC);
 802     }
 803     if(FQD) {
 804       TString qaNameIntFQD = "QAforInt_FQD_";
 805       qaNameIntFQD += type;
 806       qaNameIntFQD += ".root";
 807       AliAnalysisDataContainer *coutputQA1FQD =
 808         mgr->CreateContainer("QAintFQD", TList::Class(),AliAnalysisManager::kOutputContainer,qaNameIntFQD);
 809
 810       TString qaNameDiffFQD = "QAforDiff_FQD_";
 811       qaNameDiffFQD += type;
 812       qaNameDiffFQD += ".root";
 813       AliAnalysisDataContainer *coutputQA2FQD =
 814         mgr->CreateContainer("QAdiffFQD", TList::Class(),AliAnalysisManager::kOutputContainer,qaNameDiffFQD);
 815     }
 816     if(MCEP) {
 817       TString qaNameIntMCEP = "QAforInt_MCEP_";
 818       qaNameIntMCEP += type;
 819       qaNameIntMCEP += ".root";
 820       AliAnalysisDataContainer *coutputQA1MCEP =
 821         mgr->CreateContainer("QAintMCEP", TList::Class(),AliAnalysisManager::kOutputContainer,qaNameIntMCEP);
 822
 823       TString qaNameDiffMCEP = "QAforDiff_MCEP_";
 824       qaNameDiffMCEP += type;
 825       qaNameDiffMCEP += ".root";
 826       AliAnalysisDataContainer *coutputQA2MCEP =
 827         mgr->CreateContainer("QAdiffMCEP", TList::Class(),AliAnalysisManager::kOutputContainer,qaNameDiffMCEP);
 828     }
 829   }
 830
 831   //____________________________________________//
 832
 833   if (SP)   {
 834     mgr->ConnectInput(taskSP,0,cinput1);
 835     mgr->ConnectOutput(taskSP,0,coutputSP);
 836     if (QA) { mgr->ConnectOutput(taskSP,1,coutputQA1SP);
 837     mgr->ConnectOutput(taskSP,2,coutputQA2SP); }
 838   }
 839   if (LYZ1) {
 840     mgr->ConnectInput(taskLYZ1,0,cinput1);
 841     mgr->ConnectOutput(taskLYZ1,0,coutputLYZ1);
 842     if (QA) { mgr->ConnectOutput(taskLYZ1,1,coutputQA1LYZ1);
 843     mgr->ConnectOutput(taskLYZ1,2,coutputQA2LYZ1); }
 844   }
 845   if (LYZ2) {
 846     mgr->ConnectInput(taskLYZ2,0,cinput1);
 847     mgr->ConnectInput(taskLYZ2,1,cinputLYZ2);
 848     mgr->ConnectOutput(taskLYZ2,0,coutputLYZ2);
 849     if (QA) { mgr->ConnectOutput(taskLYZ2,1,coutputQA1LYZ2);
 850     mgr->ConnectOutput(taskLYZ2,2,coutputQA2LYZ2); }
 851     cinputLYZ2->SetData(fInputListLYZ2);
 852   }
 853   if (LYZEP) {
 854     mgr->ConnectInput(taskLYZEP,0,cinput1);
 855     mgr->ConnectInput(taskLYZEP,1,cinputLYZEP);
 856     mgr->ConnectOutput(taskLYZEP,0,coutputLYZEP);
 857     if (QA) { mgr->ConnectOutput(taskLYZEP,1,coutputQA1LYZEP);
 858     mgr->ConnectOutput(taskLYZEP,2,coutputQA2LYZEP); }
 859     cinputLYZEP->SetData(fInputListLYZEP);
 860   }
 861   if (GFC)   {
 862     mgr->ConnectInput(taskGFC,0,cinput1);
 863     mgr->ConnectOutput(taskGFC,0,coutputGFC);
 864     if (QA) { mgr->ConnectOutput(taskGFC,1,coutputQA1GFC);
 865     mgr->ConnectOutput(taskGFC,2,coutputQA2GFC); }
 866   }
 867   if (QC)   {
 868     mgr->ConnectInput(taskQC,0,cinput1);
 869     mgr->ConnectOutput(taskQC,0,coutputQC);
 870     if (QA) { mgr->ConnectOutput(taskQC,1,coutputQA1QC);
 871     mgr->ConnectOutput(taskQC,2,coutputQA2QC); }
 872     if(useWeights)
 873     {
 874      mgr->ConnectInput(taskQC,1,cinputWeights);
 875      cinputWeights->SetData(weightsList);
 876     }
 877   }
 878   if (FQD)   {
 879     mgr->ConnectInput(taskFQD,0,cinput1);
 880     mgr->ConnectOutput(taskFQD,0,coutputFQD);
 881     if (QA) { mgr->ConnectOutput(taskFQD,1,coutputQA1FQD);
 882     mgr->ConnectOutput(taskFQD,2,coutputQA2FQD); }
 883   }
 884   if (MCEP)  {
 885     mgr->ConnectInput(taskMCEP,0,cinput1);
 886     mgr->ConnectOutput(taskMCEP,0,coutputMCEP);
 887     if (QA) { mgr->ConnectOutput(taskMCEP,1,coutputQA1MCEP);
 888     mgr->ConnectOutput(taskMCEP,2,coutputQA2MCEP); }
 889   }
 890
 891   if (!mgr->InitAnalysis()) return;
 892   mgr->PrintStatus();
 893   mgr->StartAnalysis("local",chain);
 894
 895   timer.Stop();
 896   timer.Print();
 897 }
 898
 899
 900 // Helper macros for creating chains
 901 // from: CreateESDChain.C,v 1.10 jgrosseo Exp
 902
 903 TChain* CreateESDChain(const char* aDataDir, Int_t aRuns, Int_t offset)
 904 {
 905   // creates chain of files in a given directory or file containing a list.
 906   // In case of directory the structure is expected as:
 907   // <aDataDir>/<dir0>/AliESDs.root
 908   // <aDataDir>/<dir1>/AliESDs.root
 909   // ...
 910
 911   if (!aDataDir)
 912     return 0;
 913
 914   Long_t id, size, flags, modtime;
 915   if (gSystem->GetPathInfo(aDataDir, &id, &size, &flags, &modtime))
 916     {
 917       printf("%s not found.\n", aDataDir);
 918       return 0;
 919     }
 920
 921   TChain* chain = new TChain("esdTree");
 922   TChain* chaingAlice = 0;
 923
 924   if (flags & 2)
 925     {
 926       TString execDir(gSystem->pwd());
 927       TSystemDirectory* baseDir = new TSystemDirectory(".", aDataDir);
 928       TList* dirList            = baseDir->GetListOfFiles();
 929       Int_t nDirs               = dirList->GetEntries();
 930       gSystem->cd(execDir);
 931
 932       Int_t count = 0;
 933
 934       for (Int_t iDir=0; iDir<nDirs; ++iDir)
 935         {
 936           TSystemFile* presentDir = (TSystemFile*) dirList->At(iDir);
 937           if (!presentDir || !presentDir->IsDirectory() || strcmp(presentDir->GetName(), ".") == 0 || strcmp(presentDir->GetName(), "..") == 0)
 938             continue;
 939
 940           if (offset > 0)
 941             {
 942               --offset;
 943               continue;
 944             }
 945
 946           if (count++ == aRuns)
 947             break;
 948
 949           TString presentDirName(aDataDir);
 950           presentDirName += "/";
 951           presentDirName += presentDir->GetName();
 952           chain->Add(presentDirName + "/AliESDs.root/esdTree");
 953           //  cerr<<presentDirName<<endl;
 954         }
 955
 956     }
 957   else
 958     {
 959       // Open the input stream
 960       ifstream in;
 961       in.open(aDataDir);
 962
 963       Int_t count = 0;
 964
 965       // Read the input list of files and add them to the chain
 966       TString esdfile;
 967       while(in.good()) {
 968         in >> esdfile;
 969         if (!esdfile.Contains("root")) continue; // protection
 970
 971         if (offset > 0)
 972           {
 973             --offset;
 974             continue;
 975           }
 976
 977         if (count++ == aRuns)
 978           break;
 979
 980         // add esd file
 981         chain->Add(esdfile);
 982       }
 983
 984       in.close();
 985     }
 986
 987   return chain;
 988 }
 989
 990 // Helper macros for creating chains
 991 // from: CreateESDChain.C,v 1.10 jgrosseo Exp
 992
 993 TChain* CreateAODChain(const char* aDataDir, Int_t aRuns, Int_t offset)
 994 {
 995   // creates chain of files in a given directory or file containing a list.
 996   // In case of directory the structure is expected as:
 997   // <aDataDir>/<dir0>/AliAOD.root
 998   // <aDataDir>/<dir1>/AliAOD.root
 999   // ...
1000
1001   if (!aDataDir)
1002     return 0;
1003
1004   Long_t id, size, flags, modtime;
1005   if (gSystem->GetPathInfo(aDataDir, &id, &size, &flags, &modtime))
1006     {
1007       printf("%s not found.\n", aDataDir);
1008       return 0;
1009     }
1010
1011   TChain* chain = new TChain("aodTree");
1012   TChain* chaingAlice = 0;
1013
1014   if (flags & 2)
1015     {
1016       TString execDir(gSystem->pwd());
1017       TSystemDirectory* baseDir = new TSystemDirectory(".", aDataDir);
1018       TList* dirList            = baseDir->GetListOfFiles();
1019       Int_t nDirs               = dirList->GetEntries();
1020       gSystem->cd(execDir);
1021
1022       Int_t count = 0;
1023
1024       for (Int_t iDir=0; iDir<nDirs; ++iDir)
1025         {
1026           TSystemFile* presentDir = (TSystemFile*) dirList->At(iDir);
1027           if (!presentDir || !presentDir->IsDirectory() || strcmp(presentDir->GetName(), ".") == 0 || strcmp(presentDir->GetName(), "..") == 0)
1028             continue;
1029
1030           if (offset > 0)
1031             {
1032               --offset;
1033               continue;
1034             }
1035
1036           if (count++ == aRuns)
1037             break;
1038
1039           TString presentDirName(aDataDir);
1040           presentDirName += "/";
1041           presentDirName += presentDir->GetName();
1042           chain->Add(presentDirName + "/AliAOD.root/aodTree");
1043           // cerr<<presentDirName<<endl;
1044         }
1045
1046     }
1047   else
1048     {
1049       // Open the input stream
1050       ifstream in;
1051       in.open(aDataDir);
1052
1053       Int_t count = 0;
1054
1055       // Read the input list of files and add them to the chain
1056       TString aodfile;
1057       while(in.good()) {
1058         in >> aodfile;
1059         if (!aodfile.Contains("root")) continue; // protection
1060
1061         if (offset > 0)
1062           {
1063             --offset;
1064             continue;
1065           }
1066
1067         if (count++ == aRuns)
1068           break;
1069
1070         // add aod file
1071         chain->Add(aodfile);
1072       }
1073
1074       in.close();
1075     }
1076
1077   return chain;
1078 }
1079
1080 void SetupPar(char* pararchivename)
1081 {
1082   //Load par files, create analysis libraries
1083   //For testing, if par file already decompressed and modified
1084   //classes then do not decompress.
1085
1086   TString cdir(Form("%s", gSystem->WorkingDirectory() )) ;
1087   TString parpar(Form("%s.par", pararchivename)) ;
1088   if ( gSystem->AccessPathName(parpar.Data()) ) {
1089     gSystem->ChangeDirectory(gSystem->Getenv("ALICE_ROOT")) ;
1090     TString processline(Form(".! make %s", parpar.Data())) ;
1091     gROOT->ProcessLine(processline.Data()) ;
1092     gSystem->ChangeDirectory(cdir) ;
1093     processline = Form(".! mv /tmp/%s .", parpar.Data()) ;
1094     gROOT->ProcessLine(processline.Data()) ;
1095   }
1096   if ( gSystem->AccessPathName(pararchivename) ) {
1097     TString processline = Form(".! tar xvzf %s",parpar.Data()) ;
1098     gROOT->ProcessLine(processline.Data());
1099   }
1100
1101   TString ocwd = gSystem->WorkingDirectory();
1102   gSystem->ChangeDirectory(pararchivename);
1103
1104   // check for BUILD.sh and execute
1105   if (!gSystem->AccessPathName("PROOF-INF/BUILD.sh")) {
1106     printf("*******************************\n");
1107     printf("*** Building PAR archive    ***\n");
1108     cout<<pararchivename<<endl;
1109     printf("*******************************\n");
1110
1111     if (gSystem->Exec("PROOF-INF/BUILD.sh")) {
1112       Error("runProcess","Cannot Build the PAR Archive! - Abort!");
1113       return -1;
1114     }
1115   }
1116   // check for SETUP.C and execute
1117   if (!gSystem->AccessPathName("PROOF-INF/SETUP.C")) {
1118     printf("*******************************\n");
1119     printf("*** Setup PAR archive       ***\n");
1120     cout<<pararchivename<<endl;
1121     printf("*******************************\n");
1122     gROOT->Macro("PROOF-INF/SETUP.C");
1123   }
1124
1125   gSystem->ChangeDirectory(ocwd.Data());
1126   printf("Current dir: %s\n", ocwd.Data());
1127 }
1128