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