PWGGA/CaloTrackCorrelations/macros/ana.C

   1 /* $Id:  $ */
   2 //--------------------------------------------------
   3 // Example macro to do analysis with the
   4 // analysis classes in CaloTrackCorrelations
   5 // Can be executed with Root and AliRoot
   6 //
   7 // Pay attention to the options and definitions
   8 // set in the lines below
   9 //
  10 //  Author : Gustavo Conesa Balbastre (INFN-LNF)
  11 //
  12 //-------------------------------------------------
  13 enum anaModes {mLocal=0, mPROOF=1, mPlugin=2, mGRID=3};
  14 //mLocal    = 0: Analyze locally files in your computer
  15 //mPROOF    = 1: Analyze files on GRID with Plugin
  16 //mPlugin   = 2: Analyze files on GRID with Plugin
  17 //mGRID     = 3: Analyze files on GRID, jobs launched from aliensh
  18
  19 //---------------------------------------------------------------------------
  20 // Settings to read locally several files, only for "mLocal" mode
  21 // The different values are default, they can be set with environmental
  22 // variables: INDIR, PATTERN, NFILES, respectivelly
  23
  24 char * kInDir   = "/user/data/files/";
  25 char * kPattern = ""; // Data are in files kInDir/kPattern+i
  26 Int_t  kFile    = 2;
  27
  28 //---------------------------------------------------------------------------
  29 // Dataset for proof analysis, mode=mPROOF
  30 // char * kDataset = "/alice/vernet/PbPb_LHC10h_ESD";
  31
  32 char *  kDatasetPROOF     = "/alice/vernet/LHC11b_149646";
  33 Int_t   kDatasetNMaxFiles = 20;
  34 TString ccin2p3UserName   = "arbor" ;
  35 TString alienUserName     = "narbor" ;
  36
  37 //---------------------------------------------------------------------------
  38 // Collection file for grid analysis
  39
  40 char * kXML = "collection.xml";
  41
  42 //---------------------------------------------------------------------------
  43 //Scale histograms from file. Change to kTRUE when xsection file exists
  44 //Put name of file containing xsection
  45 //Put number of events per ESD file
  46 //This is an specific case for normalization of Pythia files.
  47 const char * kXSFileName = "pyxsec.root";
  48
  49 // container of xs if xs in file pyxsec_hist.root
  50 TArrayF* xsArr;
  51 TArrayI* trArr;
  52
  53 //---------------------------------------------------------------------------
  54
  55 //Set some default values, but used values are set in the code!
  56
  57 Bool_t  kMC        = kFALSE; //With real data kMC = kFALSE
  58 TString kInputData = "ESD"; //ESD, AOD, MC, deltaAOD
  59 Int_t   kYear      = 2011;
  60 TString kCollision = "pp";
  61 Bool_t  outAOD     = kFALSE; //Some tasks doesnt need it.
  62 TString kTreeName;
  63 TString kPass      = "";
  64 char    kTrigger[1024];
  65 Int_t   kRun       = 0;
  66
  67 //________________________
  68 void ana(Int_t mode=mGRID)
  69 {
  70   // Main
  71
  72   //--------------------------------------------------------------------
  73   // Load analysis libraries
  74   // Look at the method below,
  75   // change whatever you need for your analysis case
  76   // ------------------------------------------------------------------
  77
  78   LoadLibraries(mode) ;
  79   //gSystem->ListLibraries();
  80
  81   //-------------------------------------------------------------------------------------------------
  82   //Create chain from ESD and from cross sections files, look below for options.
  83   //-------------------------------------------------------------------------------------------------
  84
  85   // Set kInputData and kTreeName looking to the kINDIR
  86
  87   CheckInputData(mode);
  88
  89   // Check global analysis settings
  90
  91   CheckEnvironmentVariables();
  92
  93   printf("*********************************************\n");
  94   printf("*** Input data < %s >, pass %s, tree < %s >, MC?  < %d > ***\n",kInputData.Data(),kPass.Data(),kTreeName.Data(),kMC);
  95   printf("*********************************************\n");
  96
  97   TChain * chain   = new TChain(kTreeName) ;
  98   TChain * chainxs = new TChain("Xsection") ;
  99   CreateChain(mode, chain, chainxs);
 100
 101   Double_t scale  = -1;
 102   printf("===== kMC %d, chainxs %p\n",kMC,chainxs);
 103
 104   if(kMC)
 105   {
 106     //Get the cross section
 107     Double_t xsection = 0;
 108     Float_t  ntrials  = 0;
 109     Int_t    nfiles =  0;
 110
 111     Bool_t ok = GetAverageXsection(chainxs, xsection, ntrials, nfiles);
 112
 113     printf("n xs files %d",nfiles);
 114
 115     if(ok)
 116     {
 117       Int_t  nEventsPerFile = chain->GetEntries() / nfiles;
 118
 119       Double_t trials = ntrials / nEventsPerFile ;
 120
 121       scale = xsection / trials;
 122
 123       printf("Get Cross section : nfiles  %d, nevents %d, nevents per file %d \n",nfiles, chain->GetEntries(),nEventsPerFile);
 124       printf("                    ntrials %d, trials %2.2f, xs %2.2e, scale factor %2.2e\n", ntrials,trials,xsection,scale);
 125
 126       if(chainxs->GetEntries()!=chain->GetEntries()) printf("CAREFUL: Number of files in data chain %d, in cross section chain %d \n",
 127                                                             chainxs->GetEntries(),chain->GetEntries());
 128     } // ok
 129
 130     // comment out this line in case the simulation did not have the cross section files produced in the directory
 131     if( scale <= 0  || !ok)
 132     { printf( "STOP, cross section not available! nfiles %d \n", chainxs->GetEntries() ) ; return ; }
 133
 134   }
 135
 136
 137   printf("*********************************************\n");
 138   printf("number of entries # %lld, skipped %d\n", chain->GetEntries()) ;
 139   printf("*********************************************\n");
 140
 141   if(!chain)
 142   {
 143     printf("STOP, no chain available\n");
 144     return;
 145   }
 146
 147   AliLog::SetGlobalLogLevel(AliLog::kError);//Minimum prints on screen
 148
 149   //------------------------------------------
 150   //  Alien handler part
 151   //------------------------------------------
 152   AliAnalysisGrid *alienHandler=0x0;
 153   if(mode==mPlugin)
 154   {
 155     // Create and configure the alien handler plugin
 156     gROOT->LoadMacro("CreateAlienHandler.C");
 157     alienHandler = CreateAlienHandler();
 158     if (!alienHandler) return;
 159   }
 160
 161   //--------------------------------------
 162   // Make the analysis manager
 163   //-------------------------------------
 164   AliAnalysisManager *mgr  = new AliAnalysisManager("Manager", "Manager");
 165   //AliAnalysisManager::SetUseProgressBar(kTRUE);
 166   //mgr->SetSkipTerminate(kTRUE);
 167   //mgr->SetNSysInfo(1);
 168
 169   if(mode==mPlugin)
 170   {
 171     // Connect plugin to the analysis manager
 172     mgr->SetGridHandler(alienHandler);
 173   }
 174
 175   // MC handler
 176   if((kMC || kInputData == "MC") && !kInputData.Contains("AOD"))
 177   {
 178     AliMCEventHandler* mcHandler = new AliMCEventHandler();
 179     mcHandler->SetReadTR(kFALSE);//Do not search TrackRef file
 180     mgr->SetMCtruthEventHandler(mcHandler);
 181     if( kInputData == "MC")
 182     {
 183       cout<<"MC INPUT EVENT HANDLER"<<endl;
 184       mgr->SetInputEventHandler(NULL);
 185     }
 186   }
 187
 188
 189   // AOD output handler
 190   if(kInputData!="deltaAOD" && outAOD)
 191   {
 192     cout<<"Init output handler"<<endl;
 193     AliAODHandler* aodoutHandler   = new AliAODHandler();
 194     aodoutHandler->SetOutputFileName("aod.root");
 195     ////aodoutHandler->SetCreateNonStandardAOD();
 196     mgr->SetOutputEventHandler(aodoutHandler);
 197   }
 198
 199   //input
 200
 201   if(kInputData == "ESD")
 202   {
 203     // ESD handler
 204     AliESDInputHandler *esdHandler = new AliESDInputHandler();
 205     esdHandler->SetReadFriends(kFALSE);
 206     mgr->SetInputEventHandler(esdHandler);
 207     cout<<"ESD handler "<<mgr->GetInputEventHandler()<<endl;
 208   }
 209   else if(kInputData.Contains("AOD"))
 210   {
 211     // AOD handler
 212     AliAODInputHandler *aodHandler = new AliAODInputHandler();
 213     mgr->SetInputEventHandler(aodHandler);
 214     if(kInputData == "deltaAOD") aodHandler->AddFriend("deltaAODCaloTrackCorr.root");
 215     cout<<"AOD handler "<<mgr->GetInputEventHandler()<<endl;
 216   }
 217
 218   //mgr->RegisterExternalFile("deltaAODCaloTrackCorr.root");
 219   //mgr->SetDebugLevel(1); // For debugging, do not uncomment if you want no messages.
 220
 221   TString outputFile = AliAnalysisManager::GetCommonFileName();
 222
 223   //-------------------------------------------------------------------------
 224   // Define task, put here any other task that you want to use.
 225   //-------------------------------------------------------------------------
 226
 227   // Physics selection
 228   if(kInputData=="ESD" && !kMC)
 229   {
 230     gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPhysicsSelection.C");
 231     AliPhysicsSelectionTask* physSelTask = AddTaskPhysicsSelection(kMC);
 232   }
 233
 234   // Centrality
 235   if(kCollision=="PbPb" && kInputData=="ESD")
 236   {
 237     gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskCentrality.C");
 238     AliCentralitySelectionTask *taskCentrality = AddTaskCentrality();
 239   }
 240
 241   if(kCollision=="PbPb")
 242   {
 243     gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskVZEROEPSelection.C");
 244     AliVZEROEPSelectionTask  * EPV0 = AddTaskVZEROEPSelection();
 245
 246     gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskEventplane.C");
 247     AliEPSelectionTask * EP = AddTaskEventplane();
 248   }
 249
 250   // Simple event counting tasks
 251
 252   gROOT->LoadMacro("$ALICE_ROOT/PWGGA/CaloTrackCorrelations/macros/AddTaskCounter.C");
 253
 254   AliAnalysisTask* count    = AddTaskCounter("",kMC);   // All, fill histo with cross section and trials if kMC is true
 255   AliAnalysisTask* countmb  = AddTaskCounter("MB"); // Min Bias
 256   AliAnalysisTask* countany = AddTaskCounter("Any");
 257   AliAnalysisTask* countint = AddTaskCounter("AnyINT");// Min Bias
 258
 259   if(!kMC)
 260   {
 261     AliAnalysisTaskCounter* countemg = AddTaskCounter("EMCEGA");
 262     AliAnalysisTaskCounter* countemj = AddTaskCounter("EMCEJE");
 263     if(kCollision=="PbPb")
 264     {
 265       AliAnalysisTaskCounter* countcen = AddTaskCounter("Central");
 266       AliAnalysisTaskCounter* countsce = AddTaskCounter("SemiCentral");
 267       AliAnalysisTaskCounter* countssce= AddTaskCounter("SemiOrCentral");
 268       AliAnalysisTaskCounter* countphP = AddTaskCounter("PHOSPb");
 269     }
 270     else
 271     {
 272       AliAnalysisTaskCounter* countem1 = AddTaskCounter("EMC1"); // Trig Th > 1.5 GeV approx
 273       AliAnalysisTaskCounter* countem7 = AddTaskCounter("EMC7"); // Trig Th > 4-5 GeV
 274       AliAnalysisTaskCounter* countphp = AddTaskCounter("PHOS");
 275     }
 276   }
 277   // -----------------
 278   // Photon conversion
 279   // -----------------
 280 /*
 281   if(kInputData=="ESD"){
 282     printf("* Configure photon conversion analysis in macro \n");
 283     TString arguments = "-run-on-train -use-own-xyz  -force-aod -mc-off ";
 284     gROOT->LoadMacro("$ALICE_ROOT/PWGGA/GammaConversion/macros/ConfigGammaConversion.C");
 285     AliAnalysisTaskGammaConversion * taskGammaConversion =
 286     ConfigGammaConversion(arguments,mgr->GetCommonInputContainer());
 287     taskGammaConversion->SelectCollisionCandidates();
 288
 289     // Gamma Conversion AOD to AODPWG4Particle
 290     AliAnalysisTaskGCPartToPWG4Part * taskGCToPC = new AliAnalysisTaskGCPartToPWG4Part("GCPartToPWG4Part");
 291     taskGCToPC->SetGammaCutId("90035620401003321022000000090");
 292     mgr->AddTask(taskGCToPC);
 293     mgr->ConnectInput  (taskGCToPC, 0, mgr->GetCommonInputContainer() );
 294     mgr->ConnectOutput (taskGCToPC, 0, mgr->GetCommonOutputContainer());
 295   }
 296 */
 297
 298   Bool_t kPrint   = kFALSE;
 299   Bool_t deltaAOD = kFALSE;
 300   gROOT->LoadMacro("AddTaskCaloTrackCorr.C");   // $ALICE_ROOT/PWGGA/CaloTrackCorrelations/macros
 301   gROOT->LoadMacro("$ALICE_ROOT/PWGGA/EMCALTasks/macros/AddTaskEMCALClusterize.C"); // $ALICE_ROOT/PWGGA/EMCALTasks/macros
 302
 303   //gROOT->LoadMacro("$ALICE_ROOT/PWGGA/CaloTrackCorrelations/macros/QA/AddTaskCalorimeterQA.C");
 304   //AliAnalysisTaskCaloTrackCorrelation * qatask = AddTaskCalorimeterQA(kInputData,kYear,kPrint,kMC);
 305
 306   //gROOT->LoadMacro("$ALICE_ROOT/PWGGA/EMCALTasks/macros/AddTaskEMCALTriggerQA.C");
 307   //AliAnalysisTaskEMCALTriggerQA * qatrigtask = AddTaskEMCALTriggerQA();
 308
 309   // Calibration, bad map ...
 310
 311   Bool_t calibEE = kTRUE; // It is set automatically, but here we force to use it or not in any case
 312   Bool_t calibTT = kTRUE; // It is set automatically, but here we force to use it or not in any case
 313   if(kRun < 122195 || (kRun > 126437 && kRun < 136851) || kMC) calibTT=kFALSE ; // Recalibration parameters not available for LHC10a,b,c,e,f,g
 314   Bool_t badMap  = kTRUE; // It is set automatically, but here we force to use it or not in any case
 315
 316   if(kCollision=="pp")
 317   {
 318     printf("==================================== \n");
 319     printf("CONFIGURE ANALYSIS FOR PP COLLISIONS \n");
 320     printf("==================================== \n");
 321
 322     Bool_t  mixing    = kTRUE; // mixing in reader and hadron correlation, and pi0
 323     Bool_t  clTM      = kTRUE;
 324     Bool_t  reTM      = kFALSE; // Recalculate matches if not already done in clusterizer
 325     Bool_t  anTM      = kTRUE;  // Remove matched
 326     Bool_t  exo       = kTRUE;  // Remove exotic cells
 327     Bool_t  clnonlin  = kTRUE;  // Apply non linearity (clusterization)
 328     Bool_t  annonlin  = kFALSE; // Apply non linearity (analysis)
 329     Int_t   minEcell  = 50;     // 50  MeV (10 MeV used in reconstruction)
 330     Int_t   minEseed  = 100;    // 100 MeV
 331     Int_t   dTime     = 0;      // default, 250 ns
 332     Int_t   wTime     = 30;     // default 425 < T < 825 ns, careful if time calibration is on
 333     Int_t   unfMinE   = 15;     // Remove cells with less than 15 MeV from cluster after unfolding
 334     Int_t   unfFrac   = 1;      // Remove cells with less than 1% of cluster energy after unfolding
 335
 336     //Trigger
 337     TString clTrigger   = "";
 338     TString anTrigger   = "EMC7";
 339
 340     if(kMC)
 341     {
 342       clTrigger   = "";
 343       anTrigger   = "";
 344       dTime       = 0;
 345       wTime       = 0;
 346     }
 347
 348     Bool_t  selectEvents = kFALSE; // Select events depending on V0, pile-up and vertex quality
 349     Bool_t  qa     = kTRUE; // Do besides calorimeter QA analysis
 350     Bool_t  hadron = kTRUE; // Do besides charged track correlations analysis
 351
 352     //Analysis with clusterizer V1
 353
 354     TString arrayNameV1 = "";
 355
 356     AliAnalysisTaskEMCALClusterize * clv1 = AddTaskEMCALClusterize(arrayNameV1,outAOD,kMC,exo,"V1",clTrigger, clTM,
 357                                                                    minEcell,minEseed,dTime,wTime,unfMinE,unfFrac,
 358                                                                    calibEE,badMap,calibTT,clnonlin);
 359
 360     printf("Name of clusterizer1 array: %s\n",arrayNameV1.Data());
 361
 362     if(!kMC)
 363     {
 364
 365       AliAnalysisTaskCaloTrackCorrelation *anav1trig   = AddTaskCaloTrackCorr(kInputData, "EMCAL", kMC, selectEvents, exo, annonlin, outputFile.Data(),
 366                                                                               kYear,kCollision,anTrigger,arrayNameV1,mixing,reTM,anTM,
 367                                                                               -1,-1, qa, hadron,calibEE,badMap,calibTT,deltaAOD,kPrint,scale,kRun);
 368     }
 369
 370     AliAnalysisTaskCaloTrackCorrelation *anav1mb     = AddTaskCaloTrackCorr(kInputData, "EMCAL", kMC, selectEvents, exo, annonlin, outputFile.Data(),
 371                                                                             kYear,kCollision,"AnyINT",arrayNameV1,mixing,reTM,anTM,
 372                                                                             -1,-1, qa, hadron,calibEE,badMap,calibTT,deltaAOD,kPrint,scale,kRun);
 373
 374
 375
 376     //Analysis with clusterizer V2
 377     TString arrayNameV2 = "";
 378     AliAnalysisTaskEMCALClusterize * clv2 = AddTaskEMCALClusterize(arrayNameV2,outAOD,kMC,exo,"V2",clTrigger, clTM,
 379                                                                    minEcell,minEseed,dTime,wTime,
 380                                                                    calibEE,badMap,calibTT,clnonlin);
 381
 382     printf("Name of clusterizer2 array: %s\n",arrayNameV2.Data());
 383
 384     hadron = kFALSE;
 385     if(!kMC)
 386     {
 387
 388
 389       AliAnalysisTaskCaloTrackCorrelation *anav2tr   = AddTaskCaloTrackCorr(kInputData, "EMCAL", kMC, selectEvents, exo, annonlin, outputFile.Data(),
 390                                                                             kYear,kCollision,anTrigger,arrayNameV2,mixing,reTM,anTM,
 391                                                                             -1,-1,qa,hadron,calibEE,badMap,calibTT,deltaAOD,kPrint,scale,kRun);
 392     }
 393
 394
 395     AliAnalysisTaskCaloTrackCorrelation *anav2mb     = AddTaskCaloTrackCorr(kInputData, "EMCAL", kMC, selectEvents, exo, annonlin, outputFile.Data(),
 396                                                                             kYear,kCollision,"AnyINT",arrayNameV2,mixing,reTM,anTM,
 397                                                                             -1,-1, qa, hadron,calibEE,badMap,calibTT,deltaAOD,kPrint,scale,kRun);
 398   }
 399
 400   if(kCollision=="PbPb")
 401   {
 402     printf("====================================== \n");
 403     printf("CONFIGURE ANALYSIS FOR PbPb COLLISIONS \n");
 404     printf("====================================== \n");
 405     Bool_t  mixing    = kTRUE;
 406     Bool_t  clTM      = kTRUE;
 407     Bool_t  reTM      = kFALSE; // Recalculate matches if not already done in clusterizer
 408     Bool_t  anTM      = kTRUE;  // Remove matched
 409     Bool_t  exo       = kTRUE;  // Remove exotic cells
 410     Bool_t  clnonlin  = kTRUE;  // Apply non linearity (clusterization)
 411     Bool_t  annonlin  = kFALSE; // Apply non linearity (analysis)
 412     Int_t   minEcell  = 150;    // 50  MeV (10 MeV used in reconstruction)
 413     Int_t   minEseed  = 300;    // 100 MeV
 414     Int_t   dTime     = 0;      // default, 250 ns
 415     Int_t   wTime     = 0;      // default 425 < T < 825 ns
 416     Int_t   unfMinE   = 15;     // Remove cells with less than 15 MeV from cluster after unfolding
 417     Int_t   unfFrac   = 1;      // Remove cells with less than 1% of cluster energy after unfolding
 418
 419     // Trigger
 420     TString clTrigger = "";
 421     TString anTrigger = "EMCGA";
 422     if(kMC)
 423     {
 424       clTrigger = "";
 425       anTrigger = "";
 426       dTime       = 0;
 427       wTime       = 0;
 428     }
 429
 430     Bool_t  selectEvents = kFALSE; // Select events depending on V0, pile-up and vertex quality
 431     Bool_t  qa     = kTRUE; // Do besides calorimeter QA analysis
 432     Bool_t  hadron = kTRUE; // Do besides charged track correlations analysis
 433
 434     //Analysis with clusterizer V1
 435
 436     TString arrayNameV1 = "";
 437     AliAnalysisTaskEMCALClusterize * clv1 = AddTaskEMCALClusterize(arrayNameV1,outAOD,kMC,exo,"V1",clTrigger, clTM,
 438                                                                    minEcell,minEseed,dTime,wTime,unfMinE,unfFrac,
 439                                                                    calibEE,badMap,calibTT,clnonlin);
 440
 441     printf("Name of clusterizer1 array: %s\n",arrayNameV1.Data());
 442
 443     if(!kMC)
 444     {
 445       AliAnalysisTaskCaloTrackCorrelation *anav1c   = AddTaskCaloTrackCorr(kInputData, "EMCAL", kMC, selectEvents, exo, annonlin, outputFile.Data(),
 446                                                                            kYear,kCollision,anTrigger,arrayNameV1,mixing,reTM,anTM,
 447                                                                            0,20,qa,hadron,calibEE,badMap,calibTT,deltaAOD,kPrint,scale,kRun);
 448       AliAnalysisTaskCaloTrackCorrelation *anav1m   = AddTaskCaloTrackCorr(kInputData, "EMCAL", kMC, selectEvents, exo, annonlin, outputFile.Data(),
 449                                                                            kYear,kCollision,anTrigger,arrayNameV1,mixing,reTM,anTM,
 450                                                                            20,40,qa,hadron,calibEE,badMap,calibTT,deltaAOD,kPrint,scale,kRun);
 451       AliAnalysisTaskCaloTrackCorrelation *anav1p   = AddTaskCaloTrackCorr(kInputData, "EMCAL", kMC,  selectEvents, exo, annonlin, outputFile.Data(),
 452                                                                            kYear,kCollision,anTrigger,arrayNameV1,mixing,reTM,anTM,
 453                                                                            60,80,qa,hadron,calibEE,badMap,calibTT,deltaAOD,kPrint,scale,kRun);
 454     }
 455
 456     //Analysis with clusterizer V2
 457
 458     TString arrayNameV2 = "";
 459     AliAnalysisTaskEMCALClusterize * clv2 = AddTaskEMCALClusterize(arrayNameV2,outAOD,kMC,exo,"V2",clTrigger, clTM,
 460                                                                    minEcell,minEseed,dTime,wTime,unfMinE,unfFrac,
 461                                                                    calibEE,badMap,calibTT,clnonlin);
 462
 463     printf("Name of clusterizer2 array: %s\n",arrayNameV2.Data());
 464
 465     hadron = kFALSE;
 466
 467     if(!kMC)
 468     {
 469
 470       AliAnalysisTaskCaloTrackCorrelation *anav2cT   = AddTaskCaloTrackCorr(kInputData, "EMCAL", kMC, selectEvents, exo, annonlin, outputFile.Data(),
 471                                                                             kYear,kCollision,anTrigger,arrayNameV2,mixing,reTM,anTM,
 472                                                                             0,20,qa,hadron,calibEE,badMap,calibTT,deltaAOD,kPrint,scale,kRun);
 473       AliAnalysisTaskCaloTrackCorrelation *anav2mT   = AddTaskCaloTrackCorr(kInputData, "EMCAL", kMC, selectEvents, exo, annonlin, outputFile.Data(),
 474                                                                             kYear,kCollision,anTrigger,arrayNameV2,mixing,reTM,anTM,
 475                                                                             20,40,qa,hadron,calibEE,badMap,calibTT,deltaAOD,kPrint,scale,kRun);
 476       AliAnalysisTaskCaloTrackCorrelation *anav2pT   = AddTaskCaloTrackCorr(kInputData, "EMCAL", kMC, selectEvents, exo, annonlin, outputFile.Data(),
 477                                                                             kYear,kCollision,anTrigger,arrayNameV2,mixing,reTM,anTM,
 478                                                                             60,80,qa,hadron,calibEE,badMap,calibTT,deltaAOD,kPrint,scale,kRun);
 479     }
 480
 481     AliAnalysisTaskCaloTrackCorrelation *anav2cMB   = AddTaskCaloTrackCorr(kInputData, "EMCAL", kMC, selectEvents, exo, annonlin, outputFile.Data(),
 482                                                                            kYear,kCollision,"AnyINT",arrayNameV2,mixing,reTM,anTM,
 483                                                                            0,20,qa,hadron,calibEE,badMap,calibTT,deltaAOD,kPrint,scale,kRun);
 484     AliAnalysisTaskCaloTrackCorrelation *anav2mMB   = AddTaskCaloTrackCorr(kInputData, "EMCAL", kMC, selectEvents, exo, annonlin, outputFile.Data(),
 485                                                                            kYear,kCollision,"AnyINT",arrayNameV2,mixing,reTM,anTM,
 486                                                                            20,40,qa,hadron,calibEE,badMap,calibTT,deltaAOD,kPrint,scale,kRun);
 487     AliAnalysisTaskCaloTrackCorrelation *anav2pMB   = AddTaskCaloTrackCorr(kInputData, "EMCAL", kMC, selectEvents, exo, annonlin, outputFile.Data(),
 488                                                                            kYear,kCollision,"AnyINT",arrayNameV2,mixing,reTM,anTM,
 489                                                                            60,80,qa,hadron,calibEE,badMap,calibTT,deltaAOD,kPrint,scale,kRun);
 490
 491
 492   }
 493
 494
 495
 496   //-----------------------
 497   // Run the analysis
 498   //-----------------------
 499   mgr->InitAnalysis();
 500   mgr->PrintStatus();
 501
 502   if      (mode == mPlugin) mgr->StartAnalysis("grid");
 503   else if (mode == mPROOF ) mgr->StartAnalysis("proof",chain);
 504   else                      mgr->StartAnalysis("local",chain);
 505
 506   cout <<" Analysis ended sucessfully "<< endl ;
 507
 508 }
 509
 510 //_____________________________
 511 void  LoadLibraries(Int_t mode)
 512 {
 513
 514   if (mode == mPROOF) {
 515     //TProof::Mgr("ccalpmaster")->SetROOTVersion("ALICE_v5-27-06b");
 516     gROOT->LoadMacro("/afs/in2p3.fr/group/alice/laf/EnableAliRootForLAF.C");
 517     TProof* proof = EnableAliRootForLAF("ccaplmaster",nPROOFWorkers.Data(),ccin2p3UserName.Data(),alienUserName.Data(),"",kFALSE,kTRUE,kTRUE,"OADB:ANALYSIS:ANALYSISalice:AOD:ESD:CORRFW:STEERBase:EMCALUtils:PHOSUtils:PWGCaloTrackCorrBase:PWGGACaloTrackCorrelations:PWGEMCAL:PWGGAEMCALTasks");
 518
 519     //  TProof* proof = TProof::Open("ccaplmaster",Form("workers=%s",nPROOFWorkers.Data()));
 520
 521     //     //proof->ClearPackages();
 522     //     proof->UploadPackage("STEERBase");
 523     //     proof->UploadPackage("ESD");
 524     //     proof->UploadPackage("AOD");
 525     //     proof->UploadPackage("ANALYSIS");
 526     //     proof->UploadPackage("OADB");
 527     //     proof->UploadPackage("ANALYSISalice");
 528     //     proof->UploadPackage("CORRFW");
 529     //     //proof->UploadPackage("JETAN");
 530     //     proof->UploadPackage("PHOSUtils");
 531     //     proof->UploadPackage("EMCALUtils");
 532     //     proof->UploadPackage("PWGCaloTrackCorrBase");
 533     //     proof->UploadPackage("PWGGACaloTrackCorrelations");
 534     //     proof->UploadPackage("PWGEMCAL");
 535     //     proof->UploadPackage("PWGGAEMCALTasks");
 536
 537     //     proof->EnablePackage("STEERBase");
 538     //     proof->EnablePackage("ESD");
 539     //     proof->EnablePackage("AOD");
 540     //     proof->EnablePackage("ANALYSIS");
 541     //     proof->EnablePackage("OADB");
 542     //     proof->EnablePackage("ANALYSISalice");
 543     //     proof->EnablePackage("CORRFW");
 544     //     //proof->EnablePackage("JETAN");
 545     //     proof->EnablePackage("PHOSUtils");
 546     //     proof->EnablePackage("EMCALUtils");
 547     //     proof->EnablePackage("PWGCaloTrackCorrBase");
 548     //     proof->EnablePackage("PWGGACaloTrackCorrelations");
 549     //     proof->EnablePackage("PWGEMCAL");
 550     //     proof->EnablePackage("PWGGAEMCALTasks");
 551     return;
 552   }
 553
 554   //--------------------------------------
 555   // Load the needed libraries most of them already loaded by aliroot
 556   //--------------------------------------
 557   gSystem->Load("libTree");
 558   gSystem->Load("libGeom");
 559   gSystem->Load("libVMC");
 560   gSystem->Load("libXMLIO");
 561   gSystem->Load("libMatrix");
 562   gSystem->Load("libPhysics");
 563   gSystem->Load("libMinuit"); // Root + libraries to if reclusterization is done
 564
 565   gSystem->Load("libSTEERBase");
 566   gSystem->Load("libGui"); // Root + libraries to if reclusterization is done
 567   gSystem->Load("libCDB"); // Root + libraries to if reclusterization is done
 568   gSystem->Load("libESD"); // Root + libraries to if reclusterization is done
 569   gSystem->Load("libAOD");
 570   gSystem->Load("libRAWDatabase"); // Root + libraries to if reclusterization is done
 571   gSystem->Load("libProof");
 572   gSystem->Load("libOADB");
 573   gSystem->Load("libANALYSIS");
 574   gSystem->Load("libSTEER"); // Root + libraries to if reclusterization is done
 575
 576   gSystem->Load("libRAWDatarec"); // Root + libraries to if reclusterization is done
 577   gSystem->Load("libRAWDatasim"); // Root + libraries to if reclusterization is done
 578   gSystem->Load("libVZERObase");  // Root + libraries to if reclusterization is done
 579   gSystem->Load("libVZEROrec");   // Root + libraries to if reclusterization is done
 580
 581   gSystem->Load("libPHOSUtils");
 582
 583   gSystem->Load("libEMCALUtils");
 584   //SetupPar("EMCALUtils");
 585   gSystem->Load("libEMCALraw");  // Root + libraries to if reclusterization is done
 586   gSystem->Load("libEMCALbase"); // Root + libraries to if reclusterization is done
 587   gSystem->Load("libEMCALsim");  // Root + libraries to if reclusterization is done
 588   gSystem->Load("libEMCALrec");  // Root + libraries to if reclusterization is done
 589   //SetupPar("EMCALraw");
 590   //SetupPar("EMCALbase");
 591   //SetupPar("EMCALsim");
 592   //SetupPar("EMCALrec");
 593
 594   gSystem->Load("libANALYSISalice");
 595   gSystem->Load("libESDfilter");
 596
 597   gSystem->Load("libTender");
 598   gSystem->Load("libTenderSupplies");
 599
 600   gSystem->Load("libCORRFW");
 601   gSystem->Load("libPWGTools");
 602
 603   gSystem->Load("libPWGEMCAL");
 604   gSystem->Load("libPWGGAEMCALTasks");
 605   //SetupPar("PWGEMCAL");
 606   //SetupPar("PWGGAEMCALTasks");
 607
 608   gSystem->Load("libPWGCaloTrackCorrBase");
 609   gSystem->Load("libPWGGACaloTrackCorrelations");
 610   //SetupPar("PWGCaloTrackCorrBase");
 611   //SetupPar("PWGGACaloTrackCorrelations");
 612
 613   //gSystem->Load("libJETAN");
 614   //gSystem->Load("FASTJETAN");
 615   //gSystem->Load("PWGJE");
 616
 617   //gSystem->Load("libCORRFW");
 618   //gSystem->Load("libPWGGAGammaConv");
 619   //SetupPar("PWGGAGammaConv");
 620
 621   // needed for plugin?
 622   gSystem->AddIncludePath("-I$ALICE_ROOT");
 623   gSystem->AddIncludePath("-I./");
 624
 625 }
 626
 627 //_________________________________
 628 void SetupPar(char* pararchivename)
 629 {
 630   //Load par files, create analysis libraries
 631   //For testing, if par file already decompressed and modified
 632   //classes then do not decompress.
 633
 634   TString cdir(Form("%s", gSystem->WorkingDirectory() )) ;
 635   TString parpar(Form("%s.par", pararchivename)) ;
 636
 637   if ( gSystem->AccessPathName(pararchivename) ) {
 638     TString processline = Form(".! tar xvzf %s",parpar.Data()) ;
 639     gROOT->ProcessLine(processline.Data());
 640   }
 641
 642   TString ocwd = gSystem->WorkingDirectory();
 643   gSystem->ChangeDirectory(pararchivename);
 644
 645   // check for BUILD.sh and execute
 646   if (!gSystem->AccessPathName("PROOF-INF/BUILD.sh")) {
 647     printf("*******************************\n");
 648     printf("*** Building PAR archive    ***\n");
 649     cout<<pararchivename<<endl;
 650     printf("*******************************\n");
 651
 652     if (gSystem->Exec("PROOF-INF/BUILD.sh")) {
 653       Error("runProcess","Cannot Build the PAR Archive! - Abort!");
 654       return -1;
 655     }
 656   }
 657   // check for SETUP.C and execute
 658   if (!gSystem->AccessPathName("PROOF-INF/SETUP.C")) {
 659     printf("*******************************\n");
 660     printf("*** Setup PAR archive       ***\n");
 661     cout<<pararchivename<<endl;
 662     printf("*******************************\n");
 663     gROOT->Macro("PROOF-INF/SETUP.C");
 664   }
 665
 666   gSystem->ChangeDirectory(ocwd.Data());
 667   printf("Current dir: %s\n", ocwd.Data());
 668 }
 669
 670 //______________________________________
 671 void CheckInputData(const anaModes mode)
 672 {
 673   //Sets input data and tree
 674
 675   TString ocwd = gSystem->WorkingDirectory();
 676
 677   //---------------------------------------
 678   //Local files analysis
 679   //---------------------------------------
 680   if(mode == mLocal){
 681     //If you want to add several ESD files sitting in a common directory INDIR
 682     //Specify as environmental variables the directory (INDIR), the number of files
 683     //to analyze (NFILES) and the pattern name of the directories with files (PATTERN)
 684
 685     if(gSystem->Getenv("INDIR"))
 686       kInDir = gSystem->Getenv("INDIR") ;
 687     else cout<<"INDIR not set, use default: "<<kInDir<<endl;
 688
 689     TString sindir(kInDir);
 690     if     (sindir.Contains("pass1")) kPass = "pass1";
 691     else if(sindir.Contains("pass2")) kPass = "pass2";
 692     else if(sindir.Contains("pass3")) kPass = "pass3";
 693
 694     if(gSystem->Getenv("PATTERN"))
 695       kPattern = gSystem->Getenv("PATTERN") ;
 696     else  cout<<"PATTERN not set, use default: "<<kPattern<<endl;
 697
 698     cout<<"INDIR   : "<<kInDir<<endl;
 699     cout<<"NFILES  : "<<kFile<<endl;
 700
 701     char fileE[120] ;
 702     char fileA[120] ;
 703     char fileG[120] ;
 704     char fileEm[120] ;
 705     for (Int_t event = 0 ; event < kFile ; event++) {
 706       sprintf(fileE,  "%s/%s%d/AliESDs.root",    kInDir,kPattern,event) ;
 707       sprintf(fileA,  "%s/%s%d/AliAOD.root",     kInDir,kPattern,event) ;
 708       sprintf(fileG,  "%s/%s%d/galice.root",     kInDir,kPattern,event) ;
 709       sprintf(fileEm, "%s/%s%d/embededAOD.root", kInDir,kPattern,event) ;
 710
 711       TFile * fESD = TFile::Open(fileE) ;
 712       TFile * fAOD = TFile::Open(fileA) ;
 713
 714       //Check if file exists and add it, if not skip it
 715       if (fESD)
 716       {
 717         kTreeName  = "esdTree";
 718         kInputData = "ESD";
 719         TFile * fG = TFile::Open(fileG);
 720         if(fG) { kMC = kTRUE; fG->Close();}
 721         else     kMC = kFALSE;
 722
 723         // Get run number
 724         TTree* esdTree = (TTree*)fESD->Get("esdTree");
 725         AliESDEvent* esd = new AliESDEvent();
 726         esd->ReadFromTree(esdTree);
 727         esdTree->GetEvent(0);
 728         kRun = esd->GetRunNumber();
 729
 730         return;
 731       }
 732       else if(fAOD)
 733       {
 734         kTreeName  = "aodTree";
 735         kInputData = "AOD";
 736         if(((TTree*) fAOD->Get("aodTree"))->GetBranch("mcparticles")) kMC=kTRUE;
 737         else kMC = kFALSE;
 738
 739         // Get run number
 740         TTree* aodTree = (TTree*)fAOD->Get("aodTree");
 741         AliAODEvent* aod = new AliAODEvent();
 742         aod->ReadFromTree(aodTree);
 743         aodTree->GetEvent(0);
 744         kRun = aod->GetRunNumber();
 745         return;
 746       }
 747       else if(TFile::Open(fileEm))
 748       {
 749         kTreeName  = "aodTree";
 750         kInputData = "AOD";
 751         kMC        = kTRUE;
 752
 753         return;
 754       }
 755       else if(TFile::Open(fileG))
 756       {
 757         kTreeName  = "TE";
 758         kInputData = "MC";
 759         kMC        = kTRUE;
 760         return;
 761       }
 762     }
 763
 764     if(fESD) fESD->Close();
 765     if(fAOD) fAOD->Close();
 766
 767   }// local files analysis
 768
 769   //------------------------------
 770   //GRID xml files
 771   //-----------------------------
 772   else if(mode == mGRID){
 773     //Get colection file. It is specified by the environmental
 774     //variable XML
 775
 776     if(gSystem->Getenv("XML") )
 777       kXML = gSystem->Getenv("XML");
 778     else
 779       sprintf(kXML, "collection.xml") ;
 780
 781     if (!TFile::Open(kXML)) {
 782       printf("No collection file with name -- %s -- was found\n",kXML);
 783       return ;
 784     }
 785     else cout<<"XML file "<<kXML<<endl;
 786
 787     //Load necessary libraries and connect to the GRID
 788     gSystem->Load("libNetx") ;
 789     gSystem->Load("libRAliEn");
 790     TGrid::Connect("alien://") ;
 791
 792     //Feed Grid with collection file
 793     TGridCollection * collection = (TGridCollection*) TAlienCollection::Open(kXML);
 794     if (! collection) {
 795       AliError(Form("%s not found", kXML)) ;
 796       return kFALSE ;
 797     }
 798     TGridResult* result = collection->GetGridResult("",0 ,0);
 799
 800     for (Int_t index = 0; index < result->GetEntries(); index++) {
 801       TString alienURL = result->GetKey(index, "turl") ;
 802       cout << "================== " << alienURL << endl ;
 803
 804       if     (alienURL.Contains("pass1")) kPass = "pass1";
 805       else if(alienURL.Contains("pass2")) kPass = "pass2";
 806       else if(alienURL.Contains("pass3")) kPass = "pass3";
 807
 808       kRun = AliAnalysisManager::GetRunFromAlienPath(alienURL.Data());
 809       printf("Run number from alien path = %d\n",kRun);
 810
 811       TFile * fAOD = 0 ;
 812       //Check if file exists and add it, if not skip it
 813       if (alienURL.Contains("AliESDs.root"))
 814       {
 815         kTreeName  = "esdTree";
 816         kInputData = "ESD";
 817         alienURL.ReplaceAll("AliESDs.root","galice.root");
 818         if(TFile::Open(alienURL)) kMC=kTRUE;
 819         else kMC = kFALSE;
 820         return;
 821       }
 822       else if(alienURL.Contains("AliAOD.root"))
 823       {
 824         kTreeName  = "aodTree";
 825         kInputData = "AOD";
 826         fAOD = TFile::Open(alienURL);
 827         if(((TTree*) fAOD->Get("aodTree"))->GetBranch("mcparticles")) kMC=kTRUE;
 828         else kMC = kFALSE;
 829         return;
 830       }
 831       else if(alienURL.Contains("embededAOD.root"))
 832       {
 833         kTreeName  = "aodTree";
 834         kInputData = "AOD";
 835         kMC=kTRUE;
 836         return;
 837       }
 838       else if(alienURL.Contains("galice.root"))
 839       {
 840         kTreeName  = "TE";
 841         kInputData = "MC";
 842         kMC=kTRUE;
 843         return;
 844       }
 845     }
 846   }// xml analysis
 847   //------------------------------
 848   //PROOF files
 849   //-----------------------------
 850   else if(mode == mPROOF){
 851
 852     TFileCollection* coll  = gProof->GetDataSet(kDatasetPROOF)->GetStagedSubset();
 853
 854     TIter iter(coll->GetList());
 855
 856     TFileInfo* fileInfo = 0;
 857     while ((fileInfo = dynamic_cast<TFileInfo*> (iter())))
 858     {
 859       if (fileInfo->GetFirstUrl()) {
 860         TString ProofURL = fileInfo->GetFirstUrl()->GetUrl();
 861         cout << "================== " << ProofURL << endl ;
 862
 863         if     (ProofURL.Contains("pass1")) kPass = "pass1";
 864         else if(ProofURL.Contains("pass2")) kPass = "pass2";
 865         else if(ProofURL.Contains("pass3")) kPass = "pass3";
 866
 867         kRun = AliAnalysisManager::GetRunFromAlienPath(ProofURL.Data());
 868         printf("Run number from alien path = %d\n",kRun);
 869
 870         TFile * fAOD = 0 ;
 871         //Check if file exists and add it, if not skip it
 872         if (ProofURL.Contains("AliESDs.root"))
 873         {
 874           kTreeName  = "esdTree";
 875           kInputData = "ESD";
 876           alienURL.ReplaceAll("AliESDs.root","galice.root");
 877           if(TFile::Open(ProofURL)) kMC=kTRUE;
 878           else kMC = kFALSE;
 879
 880           return;
 881         }
 882         else if(ProofURL.Contains("AliAOD.root"))
 883         {
 884           kTreeName  = "aodTree";
 885           kInputData = "AOD";
 886           fAOD = TFile::Open(ProofURL);
 887           if(((TTree*) fAOD->Get("aodTree"))->GetBranch("mcparticles")) kMC=kTRUE;
 888           else kMC = kFALSE;
 889           return;
 890         }
 891         else if(ProofURL.Contains("embededAOD.root"))
 892         {
 893           kTreeName  = "aodTree";
 894           kInputData = "AOD";
 895           kMC=kTRUE;
 896           return;
 897         }
 898         else if(ProofURL.Contains("galice.root"))
 899         {
 900           kTreeName  = "TE";
 901           kInputData = "MC";
 902           kMC=kTRUE;
 903           return;
 904         }
 905       }
 906     }
 907   }// proof analysis
 908
 909   gSystem->ChangeDirectory(ocwd.Data());
 910
 911 }
 912
 913 //_____________________________________________________________________
 914 void CreateChain(const anaModes mode, TChain * chain, TChain * chainxs)
 915 {
 916   //Fills chain with data
 917   TString ocwd = gSystem->WorkingDirectory();
 918
 919   if(kInputData == "AOD")
 920   {
 921     xsArr = new TArrayF;
 922     trArr = new TArrayI;
 923   }
 924
 925   //---------------------------------------
 926   // Local files analysis
 927   //---------------------------------------
 928   if(mode == mLocal){
 929     //If you want to add several ESD files sitting in a common directory INDIR
 930     //Specify as environmental variables the directory (INDIR), the number of files
 931     //to analyze (NFILES) and the pattern name of the directories with files (PATTERN)
 932
 933     if(gSystem->Getenv("INDIR"))
 934       kInDir = gSystem->Getenv("INDIR") ;
 935     else cout<<"INDIR not set, use default: "<<kInDir<<endl;
 936
 937     if(gSystem->Getenv("PATTERN"))
 938       kPattern = gSystem->Getenv("PATTERN") ;
 939     else  cout<<"PATTERN not set, use default: "<<kPattern<<endl;
 940
 941     if(gSystem->Getenv("NFILES"))
 942       kFile = atoi(gSystem->Getenv("NFILES")) ;
 943     else cout<<"NFILES not set, use default: "<<kFile<<endl;
 944
 945     //Check if env variables are set and are correct
 946     if ( kInDir  && kFile) {
 947       printf("Get %d files from directory %s\n",kFile,kInDir);
 948       if ( ! gSystem->cd(kInDir) ) {//check if ESDs directory exist
 949         printf("%s does not exist\n", kInDir) ;
 950         return ;
 951       }
 952
 953       //if(gSystem->Getenv("XSFILE"))
 954       //kXSFileName = gSystem->Getenv("XSFILE") ;
 955       //else cout<<" XS file name not set, use default: "<<kXSFileName<<endl;
 956       char * kGener = gSystem->Getenv("GENER");
 957       if(kGener) {
 958         cout<<"GENER "<<kGener<<endl;
 959         if     (!strcmp(kGener,"PYTHIA")) kXSFileName = "pyxsec.root";
 960         else if(!strcmp(kGener,"HERWIG")) kXSFileName = "hexsec.root";
 961         else cout<<" UNKNOWN GENER, use default: "<<kXSFileName<<endl;
 962       }
 963       else cout<<" GENER not set, use default xs file name: "<<kXSFileName<<endl;
 964
 965       if(kInputData=="AOD")
 966       {
 967         kXSFileName = "pyxsec_hists.root";
 968         xsArr->Set(kFile);
 969         trArr->Set(kFile);
 970       }
 971
 972       cout<<"INDIR   : "<<kInDir     <<endl;
 973       cout<<"NFILES  : "<<kFile      <<endl;
 974       cout<<"PATTERN : "<<kPattern   <<endl;
 975       cout<<"XSFILE  : "<<kXSFileName<<endl;
 976
 977       TString datafile="";
 978       if     (kInputData == "ESD")        datafile = "AliESDs.root" ;
 979       else if(kInputData.Contains("AOD")) datafile = "AliAOD.root"  ;
 980       else if(kInputData == "MC")         datafile = "galice.root"  ;
 981
 982       //Loop on ESD/AOD/MC files, add them to chain
 983       Int_t event =0;
 984       Int_t skipped=0 ;
 985       char file[120] ;
 986       char filexs[120] ;
 987
 988       for (event = 0 ; event < kFile ; event++) {
 989         sprintf(file,   "%s/%s%d/%s", kInDir,kPattern,event,datafile.Data()) ;
 990         sprintf(filexs, "%s/%s%d/%s", kInDir,kPattern,event,kXSFileName) ;
 991         TFile * fData = 0 ;
 992         //Check if file exists and add it, if not skip it
 993         if ( fData = TFile::Open(file)) {
 994           if ( fData->Get(kTreeName) ) {
 995             printf("++++ Adding %s\n", file) ;
 996             chain->AddFile(file);
 997
 998             if(kInputData != "AOD")
 999             {
1000               chainxs->Add(filexs) ;
1001             }
1002             else
1003             {
1004               TFile*  fxsec = TFile::Open(filexs);
1005               if(fxsec)
1006               {
1007                 TKey* key = (TKey*)fxsec->GetListOfKeys()->At(0);
1008                 if(!key)
1009                 {
1010                   fxsec->Close();
1011                   printf("No key!");
1012                   continue;
1013                 }
1014
1015                 TList *list = dynamic_cast<TList*>(key->ReadObj());
1016                 if(!list)
1017                 {
1018                   fxsec->Close();
1019                   printf("No list!");
1020                   continue;
1021                 }
1022
1023                 Float_t xsection = ((TProfile*)list->FindObject("h1Xsec"))  ->GetBinContent(1);
1024                 Int_t   ntrials  = ((TH1F*)    list->FindObject("h1Trials"))->GetBinContent(1);
1025                 fxsec->Close();
1026
1027                 xsArr->SetAt(xsection,event);
1028                 trArr->SetAt(ntrials,event);
1029
1030                 printf("recovered xs %f, ntrials %d, event %d\n",xsection,ntrials, event);
1031                 //chainxs->Add(tree);
1032                 //fileTMP->Close();
1033
1034
1035               } // fxsec exists
1036             } // xs in AODs
1037
1038           }
1039         }
1040         else {
1041           printf("---- Skipping %s\n", file) ;
1042           skipped++ ;
1043         }
1044       }
1045     }
1046     else {
1047       TString input = "AliESDs.root" ;
1048       cout<<">>>>>> No list added, take a single file <<<<<<<<< "<<input<<endl;
1049       chain->AddFile(input);
1050     }
1051
1052   }// local files analysis
1053
1054   //------------------------------
1055   // GRID xml files
1056   //------------------------------
1057   else if(mode == mGRID){
1058     //Get colection file. It is specified by the environmental
1059     //variable XML
1060
1061     //Feed Grid with collection file
1062     TGridCollection * collection = (TGridCollection*) TAlienCollection::Open(kXML);
1063     if (! collection) {
1064       AliError(Form("%s not found", kXML)) ;
1065       return kFALSE ;
1066     }
1067
1068     TGridResult* result = collection->GetGridResult("",0 ,0);
1069
1070     // Makes the ESD chain
1071     printf("*** Getting the Chain       ***\n");
1072     for (Int_t index = 0; index < result->GetEntries(); index++) {
1073       TString alienURL = result->GetKey(index, "turl") ;
1074       cout << "================== " << alienURL << endl ;
1075       chain->Add(alienURL) ;
1076
1077       if(kInputData != "AOD")
1078       {
1079         alienURL.ReplaceAll("AliESDs.root",kXSFileName);
1080         alienURL.ReplaceAll("AliAOD.root",kXSFileName);
1081         chainxs->Add(alienURL) ;
1082       }
1083       else
1084       {
1085         alienURL.ReplaceAll("AliESDs.root","pyxsec_hists.root");
1086         alienURL.ReplaceAll("AliAOD.root", "pyxsec_hists.root");
1087         TFile*  fxsec = TFile::Open(alienURL);
1088         if(fxsec)
1089         {
1090           TKey* key = (TKey*)fxsec->GetListOfKeys()->At(0);
1091           if(!key)
1092           {
1093             fxsec->Close();
1094             printf("No key!");
1095             continue;
1096           }
1097
1098           TList *list = dynamic_cast<TList*>(key->ReadObj());
1099           if(!list)
1100           {
1101             fxsec->Close();
1102             printf("No list!");
1103             continue;
1104           }
1105
1106           Float_t xsection = ((TProfile*)list->FindObject("h1Xsec"))  ->GetBinContent(1);
1107           Int_t   ntrials  = ((TH1F*)    list->FindObject("h1Trials"))->GetBinContent(1);
1108           fxsec->Close();
1109
1110           xsArr->SetAt(xsection,event);
1111           trArr->SetAt(ntrials,event);
1112
1113           printf("recovered xs %f, ntrials %d, event %d\n",xsection,ntrials, event);
1114
1115         } // fxsec exists
1116       } // xs in AODs
1117
1118     }
1119   }// xml analysis
1120
1121   //------------------------------
1122   // PROOF
1123   //------------------------------
1124   else if (mode == mPROOF) {
1125
1126     TFileCollection* ds= gProof->GetDataSet(kDatasetPROOF)->GetStagedSubset();
1127
1128     gROOT->LoadMacro("/afs/in2p3.fr/group/alice/laf/dataset_management/CreateChainFromDataSet.C");
1129     chain = CreateChainFromDataSet(ds, kTreeName , kDatasetNMaxFiles);
1130     printf("chain has %d entries\n",chain->GetEntries());
1131   }
1132
1133   gSystem->ChangeDirectory(ocwd.Data());
1134
1135 }
1136
1137 //______________________________
1138 void CheckEnvironmentVariables()
1139 {
1140
1141   sprintf(kTrigger,"");
1142
1143   Bool_t bRecalibrate = kFALSE;
1144   Bool_t bBadChannel = kFALSE;
1145
1146   for (int i=0; i< gApplication->Argc();i++){
1147
1148 #ifdef VERBOSEARGS
1149
1150     printf("Arg  %d:  %s\n",i,gApplication->Argv(i));
1151
1152 #endif
1153
1154     TString sRun = "";
1155
1156     if (!(strcmp(gApplication->Argv(i),"--trigger")))
1157       sprintf(trigger,gApplication->Argv(i+1));
1158
1159     if (!(strcmp(gApplication->Argv(i),"--recalibrate")))
1160       bRecalibrate = atoi(gApplication->Argv(i+1));
1161
1162     if (!(strcmp(gApplication->Argv(i),"--badchannel")))
1163       bBadChannel = atoi(gApplication->Argv(i+1));
1164
1165     if (!(strcmp(gApplication->Argv(i),"--run"))){
1166       sRun = gApplication->Argv(i+1);
1167       if(sRun.Contains("LHC10")) {
1168         kYear = 2010;
1169       }
1170       else {
1171         if(kRun <=0){
1172           kRun = atoi(gApplication->Argv(i+1));
1173         }
1174         else printf("** Run number already set  to %d, do not set to %d\n",kRun,atoi(gApplication->Argv(i+1)));
1175       }//numeric run
1176     }//--run available
1177
1178   }// args loop
1179
1180   if(!sRun.Contains("LHC10")){
1181     if     ( kRun < 140000)
1182     {
1183       kYear = 2010;
1184       if( kRun >= 136851 ) kCollision = "PbPb";
1185     }
1186     else if( kRun < 170600)
1187     {
1188       kYear = 2011;
1189       if( kRun >= 166500 ) kCollision = "PbPb";
1190     }
1191     else
1192     {
1193       kYear = 2012;
1194
1195     }
1196   }
1197
1198   if(kMC) sprintf(kTrigger,"");
1199
1200   printf("*********************************************\n");
1201   //printf("*** Settings trigger %s, recalibrate %d, remove bad channels %d, year %d, collision %s, run %d ***\n",
1202   //       kTrigger,bRecalibrate,bBadChannel, kYear,kCollision.Data(), kRun);
1203   printf("*** Settings year %d, collision %s, run %d ***\n",kYear,kCollision.Data(), kRun);
1204   printf("*********************************************\n");
1205
1206 }
1207
1208 //______________________________________________________________________________
1209 Bool_t GetAverageXsection(TTree * tree, Double_t & xs, Float_t & ntr, Int_t & n)
1210 {
1211   // Read the PYTHIA statistics from the file pyxsec.root created by
1212   // the function WriteXsection():
1213   // integrated cross section (xsection) and
1214   // the  number of Pyevent() calls (ntrials)
1215   // and calculate the weight per one event xsection/ntrials
1216   // The spectrum calculated by a user should be
1217   // multiplied by this weight, something like this:
1218   // TH1F *userSpectrum ... // book and fill the spectrum
1219   // userSpectrum->Scale(weight)
1220   //
1221   // Yuri Kharlov 19 June 2007
1222   // Gustavo Conesa 15 April 2008
1223   // Add recovery of xs from pyxsec_hists.root file 15/jan/2015
1224
1225   Double_t xsection = 0 ;
1226   UInt_t    ntrials = 0 ;
1227   Int_t      nfiles = 0 ;
1228
1229   xs  = 0;
1230   ntr = 0;
1231   n   = 0;
1232   if( kInputData != "AOD" &&  tree))
1233   {
1234     nfiles =  tree->GetEntries()  ;
1235
1236     tree->SetBranchAddress("xsection",&xsection);
1237     tree->SetBranchAddress("ntrials" ,&ntrials );
1238     for(Int_t i = 0; i < nfiles; i++)
1239     {
1240       tree->GetEntry(i);
1241       if(xsection > 0)
1242       {
1243         xs  += xsection ;
1244         ntr += ntrials ;
1245         n++;
1246       }
1247       cout << "xsection " <<xsection<<" ntrials "<<ntrials<<endl;
1248     } // loop
1249   }
1250   else if( kInputData == "AOD" && xsArr))
1251   {
1252     nfiles = xsArr->GetSize();
1253
1254     for(Int_t i = 0; i < nfiles; i++)
1255     {
1256       if(xsArr->GetAt(i) > 0)
1257       {
1258         xs  += xsArr->GetAt(i) ;
1259         ntr += trArr->GetAt(i) ;
1260         n++;
1261       }
1262       cout << "xsection " <<xsArr->GetAt(i)<<" ntrials "<<trArr->GetAt(i)<<endl;
1263     } // loop
1264   }
1265   else return kFALSE;
1266
1267   xs =   xs /  n;
1268   ntr =  ntr / n;
1269   cout << "-----------------------------------------------------------------"<<endl;
1270   cout << "Average of "<< n <<" files: xsection " <<xs<<" ntrials "<<ntr<<endl;
1271   cout << "-----------------------------------------------------------------"<<endl;
1272
1273   return kTRUE;
1274
1275 }
1276
1277