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