PWG2/RESONANCES/macros/train/AddAnalysisTaskRsnEff.C

   1 //
   2 // This macro add an analysis task for computing efficiency.
   3 // It will have as output an AliCFContainer with several steps:
   4 //
   5 //  0) all resonances in MC which decay in the pair specified
   6 //  1) subset of (0) whose daughters are in acceptance
   7 //  2) subset of (1) whose daughters satisfy quality track cuts (covariance, chi square && nTPCclusters)
   8 //  3) subset of (2) whose daughters satisfy primary track cuts (nsigma to vertex, no kink daughters)
   9 //  4) subset of (3) whose daughters satisty the BB TPC compatibility cut at 3 sigma
  10 //
  11 Bool_t AddAnalysisTaskRsnEff
  12 (
  13   Bool_t      useBB    = kFALSE,
  14   Double_t    sigmaTPC = 0.065,
  15   const char *outFile  = "eff"
  16 )
  17 {
  18   // retrieve analysis manager
  19   AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
  20
  21   // common suffixes
  22   TString suf[2];
  23   suf[0] = "nopid";
  24   suf[1] = "pid";
  25
  26   // create task
  27   AliRsnAnalysisEffSE *task[2];
  28   task[0] = new AliRsnAnalysisEffSE("EffNoPID");
  29   task[1] = new AliRsnAnalysisEffSE("EffPID");
  30
  31   // set prior probabilities for PID
  32   for (Int_t i = 0; i < 2; i++)
  33   {
  34     task[i]->SetPriorProbability(AliPID::kElectron, 0.02);
  35     task[i]->SetPriorProbability(AliPID::kMuon,     0.02);
  36     task[i]->SetPriorProbability(AliPID::kPion,     0.83);
  37     task[i]->SetPriorProbability(AliPID::kKaon,     0.07);
  38     task[i]->SetPriorProbability(AliPID::kProton,   0.06);
  39     task[i]->DumpPriors();
  40   }
  41
  42   // pair definitions:
  43   // phi   --> K+ K-
  44   // kstar --> K+ pi- & K- pi+
  45   AliRsnPairDef *pairPhi    = new AliRsnPairDef('+', AliPID::kKaon, '-', AliPID::kKaon, 333);
  46   AliRsnPairDef *pairKStar1 = new AliRsnPairDef('-', AliPID::kKaon, '+', AliPID::kPion, 313);
  47   AliRsnPairDef *pairKStar2 = new AliRsnPairDef('+', AliPID::kKaon, '-', AliPID::kPion, 313);
  48   for (Int_t i = 0; i < 2; i++)
  49   {
  50     task[i]->AddPairDef(pairPhi);
  51     task[i]->AddPairDef(pairKStar1);
  52     task[i]->AddPairDef(pairKStar2);
  53   }
  54
  55   // axis definition
  56   //  0) transverse momentum
  57   //  1) pseudo-rapidity
  58   //  2) multiplicity (estimated with SPD tracklets - uncorrected)
  59   AliRsnFunctionAxis *axisPt   = new AliRsnFunctionAxis(AliRsnFunctionAxis::kPairPt,       50,  0.0,  10.0);
  60   AliRsnFunctionAxis *axisEta  = new AliRsnFunctionAxis(AliRsnFunctionAxis::kPairEta,      20, -1.5,   1.5);
  61   AliRsnFunctionAxis *axisMult = new AliRsnFunctionAxis(AliRsnFunctionAxis::kEventMult,     8,  0.0, 200.0);
  62   for (Int_t i = 0; i < 2; i++)
  63   {
  64     task[i]->AddAxis(axisMult);
  65     task[i]->AddAxis(axisPt);
  66     task[i]->AddAxis(axisEta);
  67   }
  68
  69   // define cuts for event selection:
  70   // this will determine the filling of bins in the "info" histograms
  71   // and should be computed as additional correction factor in efficiency
  72   AliRsnCutPrimaryVertex *cutVertex   = new AliRsnCutPrimaryVertex("cutVertex", 3);
  73   AliRsnCutSet           *cutSetEvent = new AliRsnCutSet("eventCuts");
  74   cutSetEvent->AddCut(cutVertex);
  75   cutSetEvent->SetCutScheme("cutVertex");
  76   for (Int_t i = 0; i < 2; i++)
  77   {
  78     task[i]->SetEventCuts(cutSetEvent);
  79   }
  80
  81   //
  82   // *** STEP 0 - All resonances which decay in the specified pairs
  83   //
  84   // This step does not need any kind of definition, since
  85   // its requirement is automatically checked during execution,
  86   // but to avoid segfaults, it is better to initialize a cut manager.
  87   //
  88   AliRsnCutMgr *cutMgrMC_step0 = new AliRsnCutMgr("mc_step0", "");
  89
  90   //
  91   // *** STEP 1 - Acceptance
  92   //
  93   // Here we add a cut on the pseudorapidity for both tracks
  94   //
  95   AliRsnCutStd *cutEta = new AliRsnCutStd("cutEta", AliRsnCutStd::kEta, -0.9, 0.9);
  96
  97   AliRsnCutMgr *cutMgrMC_step1    = new AliRsnCutMgr("mc_step1", "");
  98   AliRsnCutSet *cutSetTrack_step1 = new AliRsnCutSet("mc_step1_tracks");
  99
 100   cutSetTrack_step1->AddCut(cutEta);
 101   cutSetTrack_step1->SetCutScheme("cutEta");
 102   cutMgrMC_step1   ->SetCutSet(AliRsnCut::kParticle, cutSetTrack_step1);
 103
 104   //
 105   // *** STEP 2 - Reconstruction & track quality
 106   //
 107   // Use the interface to AliESDtrackCuts
 108   // and set only the cuts we are interested in
 109   AliRsnCutESDPrimary *cutQuality = new AliRsnCutESDPrimary("cutCov");
 110   cutQuality->GetCuts()->SetMaxCovDiagonalElements(2.0, 2.0, 0.5, 0.5, 2.0);
 111   cutQuality->GetCuts()->SetRequireSigmaToVertex(kTRUE);
 112   cutQuality->GetCuts()->SetMaxNsigmaToVertex(10000.0);
 113   cutQuality->GetCuts()->SetRequireTPCRefit(kTRUE);
 114   cutQuality->GetCuts()->SetAcceptKinkDaughters(kTRUE);
 115   cutQuality->GetCuts()->SetMinNClustersTPC(50);
 116   cutQuality->GetCuts()->SetMaxChi2PerClusterTPC(3.5);
 117
 118   AliRsnCutMgr *cutMgrESD_step2   = new AliRsnCutMgr("esd_step2", "");
 119   AliRsnCutSet *cutSetTrack_step2 = new AliRsnCutSet("esd_step2_tracks");
 120
 121   cutSetTrack_step2->AddCut(cutQuality);
 122   cutSetTrack_step2->SetCutScheme("cutQuality");
 123   cutMgrESD_step2  ->SetCutSet(AliRsnCut::kParticle, cutSetTrack_step2);
 124
 125   //
 126   // *** STEP 3 - Primary tracks
 127   //
 128   // Use the interface to AliESDtrackCuts
 129   // and set only the cuts we are interested in
 130   // we also disable the cuts we applied before, for clarity
 131   AliRsnCutESDPrimary *cutESDPrimary = new AliRsnCutESDPrimary("cutESDPrimary");
 132   cutESDPrimary->GetCuts()->SetMaxCovDiagonalElements(100.0, 100.0, 100.0, 100.0, 100.0);
 133   cutESDPrimary->GetCuts()->SetRequireSigmaToVertex(kTRUE);
 134   cutESDPrimary->GetCuts()->SetMaxNsigmaToVertex(3.0);
 135   cutESDPrimary->GetCuts()->SetRequireTPCRefit(kFALSE);
 136   cutESDPrimary->GetCuts()->SetAcceptKinkDaughters(kFALSE);
 137   cutESDPrimary->GetCuts()->SetMinNClustersTPC(0);
 138   cutESDPrimary->GetCuts()->SetMaxChi2PerClusterTPC(100000000.0);
 139
 140   AliRsnCutMgr *cutMgrESD_step3   = new AliRsnCutMgr("esd_step3", "");
 141   AliRsnCutSet *cutSetTrack_step3 = new AliRsnCutSet("esd_step3_tracks");
 142
 143   cutSetTrack_step3->AddCut(cutESDPrimary);
 144   cutSetTrack_step3->SetCutScheme("cutESDPrimary");
 145   cutMgrESD_step3  ->SetCutSet(AliRsnCut::kParticle, cutSetTrack_step3);
 146
 147   //
 148   // *** STEP 4 - Two possibilities (depend on the first macro argument)
 149   //
 150   // option 1 = Bethe-Bloch cut in 3 sigma (the sigma is one argument)
 151   // option 2 = realistic Bayesian PID with all detectors
 152   //
 153   AliRsnCutMgr *cutMgrESD_step4[2];
 154   AliRsnCutSet *cutSetTrack_step4[2];
 155
 156   cutMgrESD_step4[0] = new AliRsnCutMgr("esd_step4_nopid", "");
 157   cutMgrESD_step4[1] = new AliRsnCutMgr("esd_step4_pid", "");
 158   cutSetTrack_step4[0] = new AliRsnCutSet("esd_step4_tracks_nopid");
 159   cutSetTrack_step4[1] = new AliRsnCutSet("esd_step4_tracks_pid");
 160
 161   // Bethe-Bloch with kaon mass hypothesis
 162   AliRsnCutBetheBloch *cutKaonBB = new AliRsnCutBetheBloch("cutKaonBB", 3.0 * sigmaTPC, AliPID::kKaon);
 163   cutKaonBB->SetCalibConstant(0, 0.76176e-1);
 164   cutKaonBB->SetCalibConstant(1, 10.632);
 165   cutKaonBB->SetCalibConstant(2, 0.13279e-4);
 166   cutKaonBB->SetCalibConstant(3, 1.8631);
 167   cutKaonBB->SetCalibConstant(4, 1.9479);
 168
 169   // cuts for realistic PID match
 170   AliRsnCutStd *cutRealisticPID = new AliRsnCutStd("cutKaonPID", AliRsnCutStd::kRealisticPIDMatch, 0);
 171
 172   cutSetTrack_step4[0]->AddCut(cutKaonBB);
 173   cutSetTrack_step4[0]->SetCutScheme("cutKaonBB");
 174
 175   cutSetTrack_step4[1]->AddCut(cutRealisticPID);
 176   cutSetTrack_step4[1]->SetCutScheme("cutRealisticPID");
 177
 178   cutMgrESD_step4[0]->SetCutSet(AliRsnCut::kParticle, cutSetTrack_step4[0]);
 179   cutMgrESD_step4[1]->SetCutSet(AliRsnCut::kParticle, cutSetTrack_step4[1]);
 180
 181   // add all steps to the task
 182   for (Int_t i = 0; i < 2; i++)
 183   {
 184     task[i]->AddStepMC (cutMgrMC_step0);
 185     task[i]->AddStepMC (cutMgrMC_step1);
 186     task[i]->AddStepESD(cutMgrESD_step2);
 187     task[i]->AddStepESD(cutMgrESD_step3);
 188     task[i]->AddStepESD(cutMgrESD_step4[i]);
 189   }
 190
 191   // connect containers and finalize
 192   for (Int_t i = 0; i < 2; i++)
 193   {
 194     mgr->AddTask(task[i]);
 195     mgr->ConnectInput(task[i], 0, mgr->GetCommonInputContainer());
 196
 197     // create paths for the output in the common file
 198     Char_t infoPath[500], effPath[500];
 199     sprintf(infoPath , "%s:PWG2RSNINFO" , AliAnalysisManager::GetCommonFileName());
 200     sprintf(effPath  , "%s:PWG2RSNEFF%s", AliAnalysisManager::GetCommonFileName(), suf[i].Data());
 201
 202     // initialize and connect container for the output
 203     AliAnalysisDataContainer *info = 0x0, *out = 0x0;
 204     info = mgr->CreateContainer(Form("EffInfo_%s", suf[i].Data()), TList::Class(), AliAnalysisManager::kOutputContainer, infoPath);
 205     out  = mgr->CreateContainer(Form("EFF_%s", suf[i].Data()), TList::Class(), AliAnalysisManager::kOutputContainer, effPath);
 206
 207     mgr->ConnectOutput(task[i], 1, info);
 208     mgr->ConnectOutput(task[i], 2, out);
 209   }
 210
 211   return kTRUE;
 212 }