PWGJE/EMCALJetTasks/macros/AddTaskPtEMCalTriggerV1.C

   1 #if !defined (__CINT__) || defined (__MAKECINT__)
   2 #include "AliAnalysisManager.h"
   3 #include "AliAnalysisTaskPtEMCalTriggerV1.h"
   4 #include "AliESDtrackCuts.h"
   5 #include "AliJetContainer.h"
   6 #include <TList.h>
   7 #include <TString.h>
   8 #include <cstring>
   9 #endif
  10
  11 void AddClusterComponent(EMCalTriggerPtAnalysis::AliEMCalTriggerTaskGroup *group, bool isMC);
  12 void AddTrackComponent(EMCalTriggerPtAnalysis::AliEMCalTriggerTaskGroup *group, AliESDtrackCuts *trackcuts, bool isMC, bool isSwapEta);
  13 void AddMCParticleComponent(EMCalTriggerPtAnalysis::AliEMCalTriggerTaskGroup *group);
  14 void AddEventCounterComponent(EMCalTriggerPtAnalysis::AliEMCalTriggerTaskGroup *group, bool isMC);
  15 void AddMCJetComponent(EMCalTriggerPtAnalysis::AliEMCalTriggerTaskGroup *group, double minJetPt);
  16 void AddRecJetComponent(EMCalTriggerPtAnalysis::AliEMCalTriggerTaskGroup *group, AliESDtrackCuts *trackcuts, double minJetPt, bool isMC, bool isSwapEta);
  17 void CreateJetPtBinning(EMCalTriggerPtAnalysis::AliEMCalTriggerBinningComponent *binning);
  18 AliESDtrackCuts *CreateDefaultTrackCuts();
  19 AliESDtrackCuts *CreateHybridTrackCuts();
  20 AliESDtrackCuts *TrackCutsFactory(const char *trackCutsName);
  21
  22 AliAnalysisTask* AddTaskPtEMCalTriggerV1(
  23     bool isMC,
  24     bool usePythiaHard,
  25     const char *period ="LHC13d",
  26     const char *ntrackContainer = "",
  27     const char *nclusterContainer = "",
  28     const char *njetcontainerData = "",
  29     const char *njetcontainerMC = "",
  30     const char *ntriggerContainer = "",
  31     double jetradius = 0.5,
  32     const char *ntrackcuts = "standard"
  33 )
  34 {
  35   //AliLog::SetClassDebugLevel("EMCalTriggerPtAnalysis::AliAnalysisTaskPtEMCalTrigger", 2);
  36   AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
  37
  38   if (!mgr) {
  39     ::Error("AddTaskPtEMCalTrigger", "No analysis manager to connect to.");
  40     return NULL;
  41   }
  42
  43   if (!mgr->GetInputEventHandler()) {
  44     ::Error("AddTaskPtEMCalTrigger", "This task requires an input event handler");
  45     return NULL;
  46   }
  47
  48   bool isSwapEta = TString(period).CompareTo("LHC13f") ? kFALSE : kTRUE;
  49   EMCalTriggerPtAnalysis::AliAnalysisTaskPtEMCalTriggerV1 *pttriggertask = new EMCalTriggerPtAnalysis::AliAnalysisTaskPtEMCalTriggerV1("ptemcaltriggertask");
  50   //pttriggertask->SelectCollisionCandidates(AliVEvent::kINT7 | AliVEvent::kEMC7);                          // Select both INT7 or EMC7 triggered events
  51   pttriggertask->SelectCollisionCandidates(AliVEvent::kAny);
  52   if(isMC) pttriggertask->SetSwapThresholds();
  53
  54   mgr->AddTask(pttriggertask);
  55   if(usePythiaHard){
  56     pttriggertask->SetIsPythia(kTRUE);
  57   }
  58
  59   // Add components
  60   if(strlen(ntriggerContainer)){
  61     pttriggertask->SetCaloTriggerPatchInfoName(ntriggerContainer);
  62     EMCalTriggerPtAnalysis::AliEMCalTriggerTaskGroup *noselect = new EMCalTriggerPtAnalysis::AliEMCalTriggerTaskGroup("noselect");
  63     noselect->AddAnalysisComponent(new EMCalTriggerPtAnalysis::AliEMCalTriggerPatchAnalysisComponent("patchanalysis"));
  64     pttriggertask->AddAnalysisGroup(noselect);
  65   } else {
  66     pttriggertask->SetCaloTriggerPatchInfoName("");
  67   }
  68
  69   double jetpt[4] = {40., 60., 80., 100.};
  70   EMCalTriggerPtAnalysis::AliEMCalTriggerTaskGroup *defaultselect = new EMCalTriggerPtAnalysis::AliEMCalTriggerTaskGroup("defaultselect");
  71   defaultselect->SetEventSelection(new EMCalTriggerPtAnalysis::AliEMCalTriggerEventSelection());
  72   EMCalTriggerPtAnalysis::AliEMCalTriggerKineCuts *kineCuts = new EMCalTriggerPtAnalysis::AliEMCalTriggerKineCuts();
  73   kineCuts->SetPtRange(2., 100.);
  74   defaultselect->SetKineCuts(kineCuts);
  75   AddEventCounterComponent(defaultselect, isMC);
  76   if(isMC){
  77     AddMCParticleComponent(defaultselect);
  78     AddMCJetComponent(defaultselect, 20.);
  79     /*
  80     for(int ijpt = 0; ijpt < 4; ijpt++)
  81       AddMCJetComponent(defaultselect, jetpt[ijpt]);
  82     */
  83   }
  84   AddClusterComponent(defaultselect, isMC);
  85   AddTrackComponent(defaultselect, TrackCutsFactory(ntrackcuts), isMC, isSwapEta);
  86   AddRecJetComponent(defaultselect, TrackCutsFactory(ntrackcuts), 20., isMC, isSwapEta);
  87   /*
  88    * for(int ijpt = 0; ijpt < 4; ijpt++)
  89        AddRecJetComponent(defaultselect, TrackCutsFactory(ntrackcuts), jetpt[ijpt], isMC, isSwapEta);
  90    */
  91
  92   pttriggertask->AddAnalysisGroup(defaultselect);
  93
  94   // Add containers
  95   Bool_t isAOD = mgr->GetInputEventHandler()->IsA() == AliAODInputHandler::Class();
  96   AliParticleContainer *trackContainer = pttriggertask->AddParticleContainer(ntrackContainer);
  97   //trackContainer->SetClassName("AliVTrack");
  98   AliClusterContainer *clusterContainer = pttriggertask->AddClusterContainer(nclusterContainer);
  99   AliParticleContainer *mcpartcont = isMC ? pttriggertask->AddParticleContainer("MCParticlesSelected") : NULL;
 100
 101
 102   // Handle Jet Containers
 103   if(strlen(njetcontainerData)){
 104     AliJetContainer *jetcontainerData = pttriggertask->AddJetContainer(njetcontainerData, "EMCAL", jetradius);
 105     pttriggertask->SetDataJetContainerName("PtTriggerTaskJetsData");
 106     jetcontainerData->ConnectParticleContainer(trackContainer);
 107     jetcontainerData->SetName("PtTriggerTaskJetsData");
 108     jetcontainerData->SetJetPtCut(20.);
 109     printf("jet container added for Data\n");
 110   }
 111   if(isMC && strlen(njetcontainerMC)){
 112     AliJetContainer *jetcontainerMC = pttriggertask->AddJetContainer(njetcontainerMC, "EMCAL", jetradius);
 113     pttriggertask->SetMCJetContainerName("PtTriggerTaskJetsMC");
 114     jetcontainerMC->ConnectParticleContainer(mcpartcont);
 115     jetcontainerMC->SetName("PtTriggerTaskJetsMC");
 116     jetcontainerMC->SetJetPtCut(20.);
 117     printf("Jet container added for MC");
 118   }
 119
 120   TString containerName = mgr->GetCommonFileName();
 121   containerName += ":PtEMCalTriggerTask" + TString(ntrackcuts);
 122   printf("container name: %s\n", containerName.Data());
 123
 124   AliAnalysisDataContainer *cinput = mgr->GetCommonInputContainer();
 125   AliAnalysisDataContainer *coutput = mgr->CreateContainer("results", TList::Class(),    AliAnalysisManager::kOutputContainer, containerName.Data());
 126
 127   //Connect input/output
 128   mgr->ConnectInput(pttriggertask, 0, cinput);
 129   mgr->ConnectOutput(pttriggertask, 1, coutput);
 130
 131   return pttriggertask;
 132 }
 133
 134 void AddClusterComponent(EMCalTriggerPtAnalysis::AliEMCalTriggerTaskGroup *group, bool isMC){
 135   EMCalTriggerPtAnalysis::AliEMCalTriggerClusterAnalysisComponent *clusteranalysis = new EMCalTriggerPtAnalysis::AliEMCalTriggerClusterAnalysisComponent("clusterAnalysis");
 136   clusteranalysis->SetEnergyRange(2., 100.);
 137   if(isMC) clusteranalysis->SetUsePatches();
 138   group->AddAnalysisComponent(clusteranalysis);
 139 }
 140
 141 void AddTrackComponent(EMCalTriggerPtAnalysis::AliEMCalTriggerTaskGroup *group, AliESDtrackCuts * trackcuts, bool isMC, bool isSwapEta){
 142   EMCalTriggerPtAnalysis::AliEMCalTriggerRecTrackAnalysisComponent *trackanalysis = new EMCalTriggerPtAnalysis::AliEMCalTriggerRecTrackAnalysisComponent("trackAnalysisStandard");
 143   group->AddAnalysisComponent(trackanalysis);
 144   // Create charged hadrons pPb standard track cuts
 145   trackanalysis->SetTrackSelection(new EMCalTriggerPtAnalysis::AliEMCalPtTaskTrackSelectionESD(trackcuts));
 146
 147   if(isMC){
 148     trackanalysis->SetRequestMCtrueTracks();
 149     trackanalysis->SetUsePatches();
 150   }
 151   if(isSwapEta) trackanalysis->SetSwapEta();
 152 }
 153
 154 void AddEventCounterComponent(EMCalTriggerPtAnalysis::AliEMCalTriggerTaskGroup *group, Bool_t isMC){
 155   EMCalTriggerPtAnalysis::AliEMCalTriggerEventCounterAnalysisComponent * evcount = new EMCalTriggerPtAnalysis::AliEMCalTriggerEventCounterAnalysisComponent("eventCounter");
 156   evcount->SetUsePatches(isMC);
 157   group->AddAnalysisComponent(evcount);
 158 }
 159
 160 void AddMCParticleComponent(EMCalTriggerPtAnalysis::AliEMCalTriggerTaskGroup *group){
 161   group->AddAnalysisComponent(new EMCalTriggerPtAnalysis::AliEMCalTriggerMCParticleAnalysisComponent("partana"));
 162 }
 163
 164 void AddMCJetComponent(EMCalTriggerPtAnalysis::AliEMCalTriggerTaskGroup *group, double minJetPt){
 165   EMCalTriggerPtAnalysis::AliEMCalTriggerMCJetAnalysisComponent *jetana = new EMCalTriggerPtAnalysis::AliEMCalTriggerMCJetAnalysisComponent(Form("MCJetAna%f", minJetPt));
 166   jetana->SetMinimumJetPt(minJetPt);
 167   jetana->SetUsePatches();
 168   group->AddAnalysisComponent(jetana);
 169 }
 170
 171 void AddRecJetComponent(EMCalTriggerPtAnalysis::AliEMCalTriggerTaskGroup *group, AliESDtrackCuts *trackcuts, double minJetPt, bool isMC, bool isSwapEta){
 172   EMCalTriggerPtAnalysis::AliEMCalTriggerRecJetAnalysisComponent *jetana = new EMCalTriggerPtAnalysis::AliEMCalTriggerRecJetAnalysisComponent(Form("RecJetAna%f", minJetPt));
 173   jetana->SetMinimumJetPt(minJetPt);
 174   jetana->SetUsePatches();
 175   jetana->SetSingleTrackCuts(new EMCalTriggerPtAnalysis::AliEMCalPtTaskTrackSelectionESD(trackcuts));
 176   group->AddAnalysisComponent(jetana);
 177 }
 178
 179 void CreateJetPtBinning(EMCalTriggerPtAnalysis::AliEMCalTriggerBinningComponent *binning){
 180   // Linear binning in steps of 10 GeV/c up to 200 GeV/c
 181   TArrayD binlimits(21);
 182   for(int i = 0; i < 21; i++) binlimits[i] = 10.*i;
 183   binning->SetBinning("jetpt", binlimits);
 184 }
 185
 186 AliESDtrackCuts *CreateDefaultTrackCuts(){
 187   AliESDtrackCuts *standardTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2011(true, 1);
 188   standardTrackCuts->SetName("Standard Track cuts");
 189   standardTrackCuts->SetMinNCrossedRowsTPC(120);
 190   standardTrackCuts->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01");
 191   return standardTrackCuts;
 192 }
 193
 194 AliESDtrackCuts *CreateHybridTrackCuts(){
 195   AliESDtrackCuts* hybridTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2011(kFALSE);
 196   hybridTrackCuts->SetName("Global Hybrid tracks, loose DCA");
 197   hybridTrackCuts->SetMaxDCAToVertexXY(2.4);
 198   hybridTrackCuts->SetMaxDCAToVertexZ(3.2);
 199   hybridTrackCuts->SetDCAToVertex2D(kTRUE);
 200   hybridTrackCuts->SetMaxChi2TPCConstrainedGlobal(36);
 201   hybridTrackCuts->SetMaxFractionSharedTPCClusters(0.4);
 202   return hybridTrackCuts;
 203 }
 204
 205 AliESDtrackCuts* TrackCutsFactory(const char* trackCutsName) {
 206   if(!strcmp(trackCutsName, "standard")) return CreateDefaultTrackCuts();
 207   else if(!strcmp(trackCutsName, "hybrid")) return CreateHybridTrackCuts();
 208   return NULL;
 209 }