PWGGA/CaloTrackCorrelations/macros/AddTaskPi0EMCALPbPb.C

   1
   2 AliAnalysisTaskCaloTrackCorrelation *AddTaskPi0EMCALPbPb(
   3                                                          const TString data          = "AOD",
   4                                                          const TString calorimeter   = "EMCAL",
   5                                                          const Bool_t  printSettings = kFALSE,
   6                                                          const Bool_t  simulation    = kFALSE,
   7                                                          const Bool_t  outputAOD     = kFALSE,
   8                                                          const TString outputfile    = "",
   9                                                          const Int_t   year          = 2011,
  10                                                          const TString col           = "PbPb",
  11                                                          const TString trig          = "",
  12                                                          const TString clustersArray = ""
  13                                                         )
  14 {
  15   // Creates a PartCorr task, configures it and adds it to the analysis manager.
  16
  17   // Get the pointer to the existing analysis manager via the static access method.
  18
  19   AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
  20   if (!mgr) {
  21     ::Error("AddTask", "No analysis manager to connect to.");
  22     return NULL;
  23   }
  24
  25   // Check the analysis type using the event handlers connected to the analysis manager.
  26
  27   if (!mgr->GetInputEventHandler()) {
  28     ::Error("AddTask", "This task requires an input event handler");
  29     return NULL;
  30   }
  31   Bool_t kInputDataType = "AOD";
  32   if(!data.Contains("delta"))
  33     kInputDataType = mgr->GetInputEventHandler()->GetDataType(); // can be "ESD" or "AOD"
  34
  35   Bool_t kUseKinematics = kFALSE;
  36   if(simulation) {
  37     kUseKinematics = (mgr->GetMCtruthEventHandler())?kTRUE:kFALSE;
  38     if (!kUseKinematics && data=="AOD" && kInputDataType != "ESD") kUseKinematics = kTRUE; //AOD primary should be available ...
  39   }
  40
  41   //cout<<"********* ACCESS KINE? "<<kUseKinematics<<endl;
  42
  43   // #### Configure analysis ####
  44
  45   AliAnaPartCorrMaker * maker = new AliAnaPartCorrMaker();
  46
  47   // General frame setting and configuration
  48   maker->SetReader   (ConfigureReader(data,kInputDataType,calorimeter,clustersArray,col,outputAOD,kUseKinematics,printSettings)   );
  49   maker->SetCaloUtils(ConfigureCaloUtils(year,col,clustersArray, simulation, printSettings));
  50
  51   // Analysis tasks setting and configuration
  52   Int_t n = 0;//Analysis number, order is important
  53   maker->AddAnalysis(ConfigurePhotonAnalysis(data,calorimeter,year,col,clustersArray,simulation,"",printSettings), n++); // Photon cluster selection
  54   maker->AddAnalysis(ConfigurePi0Analysis   (data,calorimeter,year,col,clustersArray,simulation,"",printSettings), n++); // Pi0 invariant mass analysis
  55
  56   maker->SetAnaDebug(-1)  ;
  57   maker->SwitchOnHistogramsMaker()  ;
  58   if(data.Contains("delta")) maker->SwitchOffAODsMaker() ;
  59   else                       maker->SwitchOnAODsMaker()  ;
  60
  61   if(printSettings) maker->Print("");
  62
  63   //printf("<< End Configuration of %d analysis for calorimeter %s >>\n",n, calorimeter.Data());
  64
  65   // Create task
  66
  67   AliAnalysisTaskCaloTrackCorrelation * task = new AliAnalysisTaskCaloTrackCorrelation (Form("PartCorr%s_Trig%s_Cl%s",calorimeter.Data(),trig.Data(),clustersArray.Data()));
  68   task->SetConfigFileName(""); //Don't configure the analysis via configuration file.
  69   //task->SetDebugLevel(-1);
  70   task->SetBranches("ESD:AliESDRun.,AliESDHeader"); //just a trick to get Constantin's analysis to work
  71   task->SetAnalysisMaker(maker);
  72   mgr->AddTask(task);
  73
  74   //Create containers
  75
  76   if(outputfile.Length()==0)outputfile = AliAnalysisManager::GetCommonFileName();
  77   TString name(Form("%s_Trig%s_Cl%s",calorimeter.Data(),trig.Data(),clustersArray.Data()));
  78   AliAnalysisDataContainer *cout_pc   = mgr->CreateContainer(name.Data(), TList::Class(),
  79                                                              AliAnalysisManager::kOutputContainer,
  80                                                              Form("%s:%s",outputfile.Data(),name.Data()));
  81
  82   AliAnalysisDataContainer *cout_cuts = mgr->CreateContainer(Form("%sCuts",name.Data()), TList::Class(),
  83                                                              AliAnalysisManager::kParamContainer,
  84                                                              Form("%s:%sCuts",outputfile.Data(),name.Data()));
  85
  86   // Create ONLY the output containers for the data produced by the task.
  87   // Get and connect other common input/output containers via the manager as below
  88   //==============================================================================
  89   mgr->ConnectInput  (task, 0, mgr->GetCommonInputContainer());
  90   // AOD output slot will be used in a different way in future
  91   if(!data.Contains("delta")   && outputAOD) mgr->ConnectOutput (task, 0, mgr->GetCommonOutputContainer());
  92   mgr->ConnectOutput (task, 1, cout_pc);
  93   mgr->ConnectOutput (task, 2, cout_cuts);
  94
  95   return task;
  96 }
  97
  98 //____________________________________
  99 AliCaloTrackReader * ConfigureReader(TString kData,TString kInputDataType, TString kCalorimeter,
 100                                      TString kClusterArray, TString kCollisions,
 101                                      Bool_t kOutputAOD, Bool_t kUseKinematics, Bool_t kPrint)
 102 {
 103
 104   AliCaloTrackReader * reader = 0;
 105   if     (kData.Contains("AOD"))   reader = new AliCaloTrackAODReader();
 106   else if(kData=="ESD")            reader = new AliCaloTrackESDReader();
 107   else if(kData=="MC" &&
 108           kInputDataType == "ESD") reader = new AliCaloTrackMCReader();
 109
 110   reader->SetDebug(-1);//10 for lots of messages
 111
 112   //Delta AOD?
 113   //reader->SetDeltaAODFileName("");
 114   if(kOutputAOD) reader->SwitchOnWriteDeltaAOD()  ;
 115
 116   // MC settings
 117   if(kUseKinematics){
 118     if(kInputDataType == "ESD"){
 119       reader->SwitchOnStack();
 120       reader->SwitchOffAODMCParticles();
 121     }
 122     else if(kInputDataType == "AOD"){
 123       reader->SwitchOffStack();
 124       reader->SwitchOnAODMCParticles();
 125     }
 126   }
 127
 128   //------------------------
 129   // Detector input filling
 130   //------------------------
 131
 132   //Min cluster/track E
 133   reader->SetEMCALEMin(0.5);
 134   reader->SetEMCALEMax(1000);
 135   reader->SetPHOSEMin(0.3);
 136   reader->SetPHOSEMax(1000);
 137   reader->SetCTSPtMin(0.1);
 138   reader->SetCTSPtMax(1000);
 139
 140   reader->SwitchOffFiducialCut();
 141
 142   // Tracks
 143   reader->SwitchOffCTS();
 144   //gROOT->LoadMacro("$ALICE_ROOT/PWG4/macros/CreateTrackCutsPWG4.C");
 145   //AliESDtrackCuts * esdTrackCuts = CreateTrackCutsPWG4(10011004);
 146   //reader->SetTrackCuts(esdTrackCuts);
 147
 148   // Calorimeter
 149
 150   reader->SetEMCALClusterListName(kClusterArray);
 151   if(kClusterArray == "") {
 152     //printf("**************** Normal analysis **************** \n");
 153     reader->SwitchOnClusterRecalculation(); // Bad map removal
 154   }
 155   else {
 156     printf("**************** Input for analysis is Clusterizer %s **************** \n", kClusterArray.Data());
 157     reader->SwitchOffClusterRecalculation();
 158   }
 159
 160   if(kCalorimeter == "EMCAL") {
 161     reader->SwitchOnEMCALCells();
 162     reader->SwitchOnEMCAL();
 163   }
 164   if(kCalorimeter == "PHOS") {
 165     reader->SwitchOnPHOSCells();
 166     reader->SwitchOnPHOS();
 167   }
 168
 169   // for case data="deltaAOD", no need to fill the EMCAL/PHOS cluster lists
 170   if(kData.Contains("delta")){
 171     reader->SwitchOffEMCAL();
 172     reader->SwitchOffPHOS();
 173     reader->SwitchOffEMCALCells();
 174     reader->SwitchOffPHOSCells();
 175   }
 176
 177   //-----------------
 178   // Event selection
 179   //-----------------
 180
 181   //if(!kUseKinematics) reader->SetFiredTriggerClassName("CEMC7EGA-B-NOPF-CENTNOTRD"); // L1 Gamma
 182
 183   if     (kCollisions=="pp"  ) {
 184     reader->SwitchOffPileUpEventRejection();   // remove pileup by default
 185     reader->SwitchOffV0ANDSelection() ;        // and besides v0 AND
 186     reader->SwitchOffPrimaryVertexSelection(); // and besides primary vertex
 187     reader->SetZvertexCut(50.);                // Open cut
 188   }
 189   else if(kCollisions=="PbPb") {
 190     reader->SwitchOffPileUpEventRejection();   // remove pileup by default
 191     reader->SwitchOffV0ANDSelection() ;        // and besides v0 AND
 192     reader->SwitchOffPrimaryVertexSelection(); // and besides primary vertex
 193     reader->SetZvertexCut(10.);                // Centrality defined in this range.
 194
 195     // Centrality
 196     reader->SetCentralityClass("V0M");
 197     reader->SetCentralityOpt(10);  // 10 centrality bins
 198     reader->SetCentralityBin(-1,-1); // Accept all events, if not select range
 199
 200     // Event plane (only used in AliAnaPi0 for the moment)
 201     reader->SetEventPlaneMethod("Q");
 202   }
 203
 204   if(kPrint) reader->Print("");
 205
 206   return reader;
 207
 208 }
 209
 210 //_______________________________________
 211 AliCalorimeterUtils* ConfigureCaloUtils(Int_t kYears, TString kCollisions, TString kClusterArray,
 212                                         Bool_t kSimulation, Bool_t kPrint)
 213 {
 214
 215   AliCalorimeterUtils *cu = new AliCalorimeterUtils;
 216   cu->SetDebug(-1);
 217
 218   if(kYears==2010) cu->SetEMCALGeometryName("EMCAL_FIRSTYEARV1");
 219   else             cu->SetEMCALGeometryName("EMCAL_COMPLETEV1");
 220
 221   // Remove clusters close to borders, at least max energy cell is 1 cell away
 222   cu->SetNumberOfCellsFromEMCALBorder(1);
 223   cu->SetNumberOfCellsFromPHOSBorder(2);
 224
 225   //if(kClusterArray == "")
 226   //  cu->SwitchOffRecalculateClusterTrackMatching(); // Done in clusterization
 227   //else
 228   //  cu->SwitchOnRecalculateClusterTrackMatching();
 229
 230   //EMCAL only settings
 231   AliEMCALRecoUtils * recou = cu->GetEMCALRecoUtils();
 232
 233
 234   if(kCollisions=="pp"  ) { // Do only for pp for the moment
 235     cu->SwitchOnCorrectClusterLinearity();
 236     if(!kSimulation) recou->SetNonLinearityFunction(AliEMCALRecoUtils::kBeamTestCorrected);
 237     else             recou->SetNonLinearityFunction(AliEMCALRecoUtils::kPi0MC);
 238   }
 239
 240   recou->SwitchOnRejectExoticCell();
 241   if(kClusterArray == "") recou->SwitchOnRejectExoticCluster();
 242   else                    recou->SwitchOffRejectExoticCluster();
 243
 244   if(kCalorimeter=="PHOS")
 245   {
 246     if (kYears <  2014) cu->SetNumberOfSuperModulesUsed(3);
 247     else                cu->SetNumberOfSuperModulesUsed(4);
 248   }
 249   else
 250   {
 251     if      (kYears == 2010) cu->SetNumberOfSuperModulesUsed(4); //EMCAL first year
 252     else if (kYears <  2014) cu->SetNumberOfSuperModulesUsed(10);
 253     else                     cu->SetNumberOfSuperModulesUsed(20);
 254   }
 255
 256
 257   if(kPrint) cu->Print("");
 258
 259   return cu;
 260
 261 }
 262
 263 //_____________________________________
 264 AliAnaPhoton* ConfigurePhotonAnalysis(TString kData, TString kCalorimeter, Int_t kYears, TString kCollisions, TString kClusterArray,
 265                                       Bool_t kSimulation, TString kTrig = "", Bool_t kPrint)
 266 {
 267
 268   AliAnaPhoton *anaphoton = new AliAnaPhoton();
 269   anaphoton->SetDebug(-1); //10 for lots of messages
 270
 271   // cluster selection cuts
 272
 273   anaphoton->SwitchOffFiducialCut();
 274
 275   anaphoton->SetCalorimeter(kCalorimeter);
 276
 277   if(kCalorimeter == "PHOS"){
 278     anaphoton->SetNCellCut(2);// At least 2 cells
 279     anaphoton->SetMinPt(0.3);
 280     anaphoton->SetMinDistanceToBadChannel(2, 4, 5);
 281   }
 282   else {//EMCAL
 283     anaphoton->SetNCellCut(1);// At least 2 cells
 284     anaphoton->SetMinPt(0.5); // avoid mip peak at E = 260 MeV
 285     anaphoton->SetMaxPt(1000);
 286     anaphoton->SetTimeCut(-1000,1000); // open cut, usual time window of [425-825] ns if time recalibration is off
 287     anaphoton->SetMinDistanceToBadChannel(1, 2, 3); // For filtered AODs, new releases.
 288   }
 289
 290   anaphoton->SwitchOnTrackMatchRejection() ;
 291
 292   //PID cuts (shower shape)
 293   anaphoton->SwitchOnCaloPID(); // do PID selection, unless specified in GetCaloPID, selection not based on bayesian
 294   AliCaloPID* caloPID = anaphoton->GetCaloPID();
 295   //Not used in bayesian
 296
 297   //EMCAL
 298   caloPID->SetEMCALLambda0CutMax(0.30);
 299   caloPID->SetEMCALLambda0CutMin(0.10);
 300
 301   caloPID->SetEMCALDEtaCut(0.025);
 302   caloPID->SetEMCALDPhiCut(0.05);
 303   // In case of official AODs when dX and dZ was not stored, open the cuts
 304   // and rely on having a match recorded. In case of reclusterization, try.
 305   if(kData=="AOD" && kClusterArray==""){
 306     caloPID->SetEMCALDEtaCut(2000);
 307     caloPID->SetEMCALDPhiCut(2000);
 308   }
 309
 310   //PHOS
 311   caloPID->SetPHOSDispersionCut(2.5);
 312   caloPID->SetPHOSRCut(2.);
 313   if(kData=="AOD") caloPID->SetPHOSRCut(2000.); // Open cut since dX, dZ not stored
 314
 315   if(kCalorimeter=="PHOS"){
 316     caloPID->SetHistoDEtaRangeAndNBins(-200, 200, 200); // dZ
 317     caloPID->SetHistoDPhiRangeAndNBins(-200, 200, 200); // dX
 318   }
 319
 320   //caloPID->SetTOFCut(10000000); // Not used, only to set PID bits
 321
 322   anaphoton->SwitchOffFillShowerShapeHistograms();  // Filled before photon shower shape selection
 323
 324   // Input / output delta AOD settings
 325
 326   if(!kData.Contains("delta")) {
 327     anaphoton->SetOutputAODName(Form("Photon%s_Trig%s_Cl%s",kCalorimeter.Data(), kTrig.Data(),kClusterArray.Data()));
 328     anaphoton->SetOutputAODClassName("AliAODPWG4ParticleCorrelation");
 329     //anaphoton->SetOutputAODClassName("AliAODPWG4Particle"); // use if no correlation done
 330   }
 331   else anaphoton->SetInputAODName(Form("Photon%s_Trig%s_Cl%s",kCalorimeter.Data(), kTrig.Data(),kClusterArray.Data()));
 332
 333   //Set Histograms name tag, bins and ranges
 334
 335   anaphoton->AddToHistogramsName("AnaPhoton_");
 336   SetHistoRangeAndNBins(anaphoton); // see method below
 337
 338   // Number of particle type MC histograms
 339   anaphoton->FillNOriginHistograms(8);
 340   anaphoton->FillNPrimaryHistograms(4);
 341
 342   if(kPrint) anaphoton->Print("");
 343
 344   return anaphoton;
 345
 346 }
 347
 348
 349
 350 //_______________________________
 351 AliAnaPi0* ConfigurePi0Analysis(TString kData, TString kCalorimeter, Int_t kYears, TString kCollisions, TString kClusterArray,
 352                                 Bool_t kSimulation, TString kTrig = "", Bool_t kPrint)
 353 {
 354
 355   AliAnaPi0 *anapi0 = new AliAnaPi0();
 356
 357   anapi0->SetDebug(-1);//10 for lots of messages
 358   if(kPrint) anapi0->Print("");
 359
 360   // Input delta AOD settings
 361   anapi0->SetInputAODName(Form("Photon%s_Trig%s_Cl%s",kCalorimeter.Data(), kTrig.Data(),kClusterArray.Data()));
 362
 363   // Calorimeter settings
 364   anapi0->SetCalorimeter(kCalorimeter);
 365
 366   //settings for pp collision mixing
 367   anapi0->SwitchOnOwnMix(); //Off when mixing done with general mixing frame
 368
 369   // Cuts
 370   if(kCalorimeter=="EMCAL") anapi0->SetPairTimeCut(70);
 371
 372   if     (kCollisions=="pp"  ) {
 373     printf("****** Configure Pi0 for pp analysis\n");
 374     anapi0->SetNCentrBin(1);
 375     anapi0->SetNZvertBin(10);
 376     anapi0->SetNRPBin(1);
 377     anapi0->SetNMaxEvMix(100);
 378   }
 379   else if(kCollisions=="PbPb") {
 380     printf("****** Configure Pi0 for PbPb analysis\n");
 381     anapi0->SetNCentrBin(5);
 382     anapi0->SetNZvertBin(3);
 383     anapi0->SetNRPBin(1);
 384     anapi0->SetNMaxEvMix(5);
 385   }
 386
 387   anapi0->SwitchOffMultipleCutAnalysis();
 388   anapi0->SwitchOffSMCombinations();
 389
 390   //Set Histograms name tag, bins and ranges
 391
 392   anapi0->AddToHistogramsName("AnaPi0_");
 393   SetHistoRangeAndNBins(anapi0, kCalorimeter, kYears, kCollisions, kSimulation); // see method below
 394
 395   return anapi0;
 396
 397 }
 398
 399 //________________________________________________________
 400 void SetHistoRangeAndNBins (AliAnaPartCorrBaseClass* ana, TString kCalorimeter,
 401                             Int_t kYears, TString kCollisions, Bool_t kSimulation)
 402 {
 403   // Set common bins for all analysis and MC histograms filling
 404
 405   if(kSimulation) ana->SwitchOnDataMC() ;//Access MC stack and fill more histograms, AOD MC not implemented yet.
 406   else            ana->SwitchOffDataMC() ;
 407
 408   ana->SetHistoPtRangeAndNBins(0, 100, 250) ; // Energy and pt histograms
 409
 410   if(kCalorimeter=="EMCAL"){
 411     if(kYears==2010){
 412       ana->SetHistoPhiRangeAndNBins(78*TMath::DegToRad(), 122*TMath::DegToRad(), 78) ;
 413       ana->SetHistoXRangeAndNBins(-230,90,120); // QA
 414       ana->SetHistoYRangeAndNBins(370,450,40);  // QA
 415     }
 416     else {
 417       ana->SetHistoPhiRangeAndNBins(78*TMath::DegToRad(), 182*TMath::DegToRad(), 108) ;
 418       ana->SetHistoXRangeAndNBins(-600,90,200); // QA
 419       ana->SetHistoYRangeAndNBins(100,450,100); // QA
 420     }
 421
 422     ana->SetHistoEtaRangeAndNBins(-0.72, 0.72, 144) ;
 423   }
 424   else{
 425     ana->SetHistoPhiRangeAndNBins(260*TMath::DegToRad(), 320*TMath::DegToRad(), 60) ;
 426     ana->SetHistoEtaRangeAndNBins(-0.13, 0.13, 130) ;
 427
 428   }
 429
 430   ana->SetHistoShowerShapeRangeAndNBins(0, 3, 300);
 431
 432   // Invariant mass analysis
 433   ana->SetHistoMassRangeAndNBins(0., 1., 200) ;
 434   ana->SetHistoAsymmetryRangeAndNBins(0., 1. , 100) ;
 435
 436   // check if time calibration is on
 437   ana->SetHistoTimeRangeAndNBins(-1000.,1000,1000);
 438   ana->SetHistoDiffTimeRangeAndNBins(-200, 200, 800);
 439
 440   // QA, electron, charged
 441   ana->SetHistoPOverERangeAndNBins(0,10.,100);
 442   ana->SetHistodEdxRangeAndNBins(0.,200.,200);
 443
 444   // QA
 445   ana->SetHistoFinePtRangeAndNBins(0, 10, 200) ; // bining for fhAmpId
 446   ana->SetHistodRRangeAndNBins(0.,TMath::Pi(),150);
 447   ana->SetHistoRatioRangeAndNBins(0.,2.,100);
 448   ana->SetHistoVertexDistRangeAndNBins(0.,500.,500);
 449   ana->SetHistoNClusterCellRangeAndNBins(0,500,500);
 450   ana->SetHistoZRangeAndNBins(-400,400,200);
 451   ana->SetHistoRRangeAndNBins(400,450,25);
 452   ana->SetHistoV0SignalRangeAndNBins(0,5000,500);
 453   ana->SetHistoV0MultiplicityRangeAndNBins(0,5000,500);
 454   ana->SetHistoTrackMultiplicityRangeAndNBins(0,5000,500);
 455
 456 }
 457