2 // This macro add an analysis task for computing efficiency.
3 // It will have as output an AliCFContainer with several steps:
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
11 Bool_t AddRsnAnalysisTaskEffPhi
15 const char *path = "$(ALICE_ROOT)/PWG2/RESONANCES/macros/train"
18 // ==================================================================================================================
19 // == OPTIONS =======================================================================================================
20 // ==================================================================================================================
22 // Instead or getting confused with plenty of arguments in the macro (with default values),
23 // we use a unique string of options with a set of conventional strings to set up the job:
24 // -- "MC"/"DATA" --> what kind of sample
25 // -- "ITS"/"TPC" --> what tracks to use (ITS standalone and/or TPC+ITS)
26 // -- "xxxPID" --> add the PID cuts for the detector xxx.
28 // In this point, these options are converted into boolean variables.
32 opt.ReplaceAll(" ", "");
34 Bool_t addITS = opt.Contains("ITS");
35 Bool_t addTPC = opt.Contains("TPC");
36 Bool_t useITS = opt.Contains("ITSPID");
37 Bool_t useTPC = opt.Contains("TPCPID");
38 Bool_t useTOF = opt.Contains("TOFPID");
40 // correct options when needed
41 if (!addITS) useITS = kFALSE;
42 if (!addTPC) useTPC = useTOF = kFALSE;
44 // ==================================================================================================================
45 // == DEFINITIONS ===================================================================================================
46 // ==================================================================================================================
48 // We put here the definitions of all objects which are needed in the following, in order to have then
49 // a more readable code in the points where these objects are added to the analysis manager.
53 AliRsnPairDef *pairPhi = new AliRsnPairDef(AliPID::kKaon, '+', AliPID::kKaon, '-', 333, 1.019455);
55 // ==================================================================================================================
56 // == CUTS AND AXES =================================================================================================
57 // ==================================================================================================================
60 // Track quality for ITS standalone:
61 // this cut is used to select tracks of good quality, irrespective of the PID.
62 // When adding status flags, the second argument tells if each considered flag
63 // must be active or not in the track status, since the ITS-SA tracks need that
64 // some of them are OFF (e.g.: kTPCin)
66 AliRsnCutTrackQuality *cutQualityITS = new AliRsnCutTrackQuality("cutQualityITS");
67 cutQualityITS->AddStatusFlag(AliESDtrack::kITSin , kTRUE);
68 cutQualityITS->AddStatusFlag(AliESDtrack::kTPCin , kFALSE);
69 cutQualityITS->AddStatusFlag(AliESDtrack::kITSrefit , kTRUE);
70 cutQualityITS->AddStatusFlag(AliESDtrack::kTPCrefit , kFALSE);
71 cutQualityITS->AddStatusFlag(AliESDtrack::kITSpureSA, kFALSE);
72 cutQualityITS->AddStatusFlag(AliESDtrack::kITSpid , kTRUE);
73 cutQualityITS->SetPtRange(0.15, 1E+20);
74 cutQualityITS->SetEtaRange(-0.8, 0.8);
75 cutQualityITS->SetDCARPtFormula("0.0595+0.0182/pt^1.55");
76 cutQualityITS->SetDCAZmax(2.0);
77 cutQualityITS->SetSPDminNClusters(1);
78 cutQualityITS->SetITSminNClusters(4);
79 cutQualityITS->SetITSmaxChi2(2.0);
80 cutQualityITS->SetTPCminNClusters(0);
81 cutQualityITS->SetTPCmaxChi2(1E+10);
82 cutQualityITS->SetRejectKinkDaughters();
85 // Track quality for TPC+ITS:
86 // works exactly like the one above, but has settings for selecting TPC+ITS tracks
87 // in this case, the flags required are all necessary, so here the procedure is simpler
89 AliRsnCutTrackQuality *cutQualityTPC = new AliRsnCutTrackQuality("cutQualityTPC");
90 cutQualityTPC->AddStatusFlag(AliESDtrack::kTPCin , kTRUE);
91 cutQualityTPC->AddStatusFlag(AliESDtrack::kTPCrefit, kTRUE);
92 cutQualityTPC->AddStatusFlag(AliESDtrack::kITSrefit, kTRUE);
93 cutQualityTPC->SetPtRange(0.15, 1E+20);
94 cutQualityTPC->SetEtaRange(-0.8, 0.8);
95 cutQualityTPC->SetDCARPtFormula("0.0182+0.0350/pt^1.01");
96 cutQualityTPC->SetDCAZmax(2.0);
97 cutQualityTPC->SetSPDminNClusters(1);
98 cutQualityTPC->SetITSminNClusters(0);
99 cutQualityTPC->SetITSmaxChi2(1E+20);
100 cutQualityTPC->SetTPCminNClusters(70);
101 cutQualityTPC->SetTPCmaxChi2(4.0);
102 cutQualityTPC->SetRejectKinkDaughters();
106 // In this implementation, it is a 3sigma cut around the Bethe-Bloch value.
109 // --> The initialization of the BB is different between data and MC.
110 // --> The cut is the same for all momenta.
112 AliRsnCutPIDITS *cutPIDITSkaon = new AliRsnCutPIDITS("cutPIDITSkaon", AliPID::kKaon, -3.0, 3.0);
113 cutPIDITSkaon->SetMC(kTRUE);
114 cutPIDITSkaon->SetMomentumRange(0.0, 1E20);
118 // In this implementation, there are two instances:
119 // - below 350 MeV --> 5 sigma cut
120 // - above 350 MeV --> 3 sigma cut
123 // --> The initialization of the BB is different between data and MC.
125 AliRsnCutPIDTPC *cutPIDTPCkaonLow = new AliRsnCutPIDTPC("cutPIDTPCkaonLow" , AliPID::kKaon, -5.0, 5.0);
126 AliRsnCutPIDTPC *cutPIDTPCkaonHigh = new AliRsnCutPIDTPC("cutPIDTPCkaonHigh", AliPID::kKaon, -3.0, 3.0);
128 // assign the momentum range and tell to reject tracks outside it
129 cutPIDTPCkaonLow ->SetMomentumRange(0.00, 0.35);
130 cutPIDTPCkaonHigh->SetMomentumRange(0.35, 1E20);
131 cutPIDTPCkaonLow ->SetRejectOutside(kTRUE);
132 cutPIDTPCkaonHigh->SetRejectOutside(kTRUE);
134 // BB parameterization depends on data sample (MC, data)
135 // the momentum range is passed and tracks outside it are rejected
137 bbPar[0] = 2.15898 / 50.0;
138 bbPar[1] = 1.75295E1;
139 bbPar[2] = 3.40030E-9;
142 cutPIDTPCkaonLow ->SetBBParam(bbPar);
143 cutPIDTPCkaonHigh->SetBBParam(bbPar);
147 // In this implementation it is a 3sigma cout aroung expected kaon time.
150 // --> It is important to choose if this cut must reject tracks not matched in TOF.
151 // Usually, if TPC pid is used, we can accept them, otherwise we must reject.
152 // (here we assume TPC is used)
154 AliRsnCutPIDTOF *cutPIDTOFkaon = new AliRsnCutPIDTOF("cutPIDTOFkaon", AliPID::kKaon, -3.0, 3.0);
155 cutPIDTOFkaon->SetRejectUnmatched(!useTPC);
159 // Only thing to consider is that it needs a support object to define mass
161 AliRsnCutValue *cutRapidity = new AliRsnCutValue("cutY", AliRsnValue::kPairY, -0.5, 0.5);
162 cutRapidity->GetValueObj()->SetSupportObject(pairPhi);
168 // --> multiplicity has variable bins, defined by array below
171 Double_t mult[] = { 0., 1., 2., 3., 4., 5., 6., 7., 8., 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19.,
172 20., 21., 22., 23., 24., 25., 30., 35., 40., 50., 60., 70., 80., 90., 100., 120., 140., 160., 180., 200., 500.};
173 Int_t nmult = sizeof(mult) / sizeof(mult[0]);
175 AliRsnValue *axisIM = new AliRsnValue("IM" , AliRsnValue::kPairInvMass , 0.9, 1.4, 0.001);
176 AliRsnValue *axisRes = new AliRsnValue("Res" , AliRsnValue::kPairInvMassRes, -0.5, 0.5, 0.001);
177 AliRsnValue *axisPt = new AliRsnValue("PT" , AliRsnValue::kPairPt , 0.0, 5.0, 0.1 );
178 AliRsnValue *axisY = new AliRsnValue("Y" , AliRsnValue::kPairY , -1.1, 1.1, 0.1 );
179 AliRsnValue *axisMultESD = new AliRsnValue("MESD", AliRsnValue::kEventMultESDCuts, nmult, mult);
180 AliRsnValue *axisMultSPD = new AliRsnValue("MSPD", AliRsnValue::kEventMultSPD , nmult, mult);
181 AliRsnValue *axisMultMC = new AliRsnValue("MMC" , AliRsnValue::kEventMultMC , nmult, mult);
183 // ==================================================================================================================
184 // == PRELIMINARY OPERATIONS ========================================================================================
185 // ==================================================================================================================
187 // retrieve analysis manager
188 AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
191 AliRsnAnalysisTaskEffPair *task = new AliRsnAnalysisTaskEffPair(Form("RsnEff_%s", options));
192 task->SelectCollisionCandidates();
194 // add pair definition, to choose the checked resonance
195 task->AddDef(pairPhi);
197 // add the output histogram axis
198 task->AddAxis(axisIM);
199 task->AddAxis(axisRes);
200 task->AddAxis(axisPt);
201 task->AddAxis(axisY);
202 task->AddAxis(axisMultSPD);
203 task->AddAxis(axisMultMC);
205 // ==================================================================================================================
206 // == EVENT CUTS ====================================================================================================
207 // ==================================================================================================================
209 gROOT->LoadMacro(Form("%s/AddRsnEventComputations.C", path));
210 AddRsnEventComputations(kTRUE, evtopts);
212 // ==================================================================================================================
213 // == STEPS =========================================================================================================
214 // ==================================================================================================================
216 Char_t qualityITS[255], qualityTPC[255];
217 Char_t pidITS[255], pidTPC[255], pidTOF[255];
218 Char_t schemeITS[255], schemeTPC[255], scheme[255];
220 sprintf(qualityITS, "%s" , cutQualityITS->GetName());
221 sprintf(qualityTPC, "%s" , cutQualityTPC->GetName());
222 sprintf(pidITS , "%s" , cutPIDITSkaon->GetName());
223 sprintf(pidTPC , "(%s|%s)", cutPIDTPCkaonHigh->GetName(), cutPIDTPCkaonLow->GetName());
224 sprintf(pidTOF , "%s" , cutPIDTOFkaon->GetName());
225 sprintf(schemeITS , "");
226 sprintf(schemeTPC , "");
227 sprintf(scheme , "");
230 // *** STEP 0 - All resonances which decay in the specified pair
232 // This step does not need any kind of definition, since
233 // its requirement is automatically checked during execution,
234 // but to avoid segfaults, it is needed to initialize a cut manager.
236 AliRsnCutManager *mgr_step0 = new AliRsnCutManager("mc_step0", "");
239 // *** STEP 1 - Track quality
241 // All resonances whose daughters were reconstructed
242 // and pass quality track cuts will enter this step
244 AliRsnCutManager *mgr_step1 = new AliRsnCutManager("rec_step1", "");
245 AliRsnCutSet *set_step1 = mgr_step1->GetCommonDaughterCuts();
247 if (addTPC && addITS) {
248 set_step1->AddCut(cutQualityTPC);
249 set_step1->AddCut(cutQualityITS);
250 set_step1->SetCutScheme(Form("%s|%s", qualityTPC, qualityITS));
252 set_step1->AddCut(cutQualityTPC);
253 set_step1->SetCutScheme(qualityTPC);
255 set_step1->AddCut(cutQualityITS);
256 set_step1->SetCutScheme(qualityITS);
258 ::Error("Need to ad at least one between ITS and TPC tracks");
261 ::Info("AddRsnAnalysisTaskEffPhi", "Cut scheme for step #1: %s", set_step1->GetCutScheme().Data());
266 // Add all TPC cuts, according to options
268 AliRsnCutManager *mgr_step2 = new AliRsnCutManager("esd_step2", "");
269 AliRsnCutSet *set_step2 = mgr_step2->GetCommonDaughterCuts();
271 if (addITS && useITS) {
272 sprintf(schemeITS, "%s & %s", qualityITS, pidITS);
273 set_step2->AddCut(cutPIDITSkaon);
276 if (useTPC && useTOF) {
277 set_step2->AddCut(cutPIDTPCkaonLow);
278 set_step2->AddCut(cutPIDTPCkaonHigh);
279 set_step2->AddCut(cutPIDTOFkaon);
280 sprintf(schemeTPC, "%s & %s", pidTPC, pidTOF);
282 set_step2->AddCut(cutPIDTPCkaonLow);
283 set_step2->AddCut(cutPIDTPCkaonHigh);
284 sprintf(schemeTPC, "%s & %s", pidTPC);
286 set_step2->AddCut(cutPIDTOFkaon);
287 sprintf(schemeTPC, "%s & %s", pidTOF);
291 // final scheme depends on what of the above were added
292 // in case both ITS-SA and TPC tracks are added, the scheme
293 // is the OR of the cuts for the first and those for the second
294 // category of tracks
295 if (strlen(schemeITS) > 0 && strlen(schemeTPC) > 0) {
296 sprintf(scheme, "(%s) | (%s)", schemeITS, schemeTPC);
297 } else if (strlen(schemeITS) > 0) {
298 sprintf(scheme, "%s", schemeITS);
299 } else if (strlen(schemeTPC) > 0) {
300 sprintf(scheme, "%s", schemeTPC);
302 ::Error("Scheme is empty!");
306 set_step2->SetCutScheme(scheme);
307 ::Info("AddRsnAnalysisTaskEffPhi", "Cut scheme for step #1: %s", set_step2->GetCutScheme().Data());
309 // add all steps to the task:
310 // - first step computed on MC
311 // - all other steps computed on reconstruction
312 task->AddStepMC(mgr_step0);
313 task->AddStepRec(mgr_step1);
314 task->AddStepRec(mgr_step2);
316 // add the task to the manager and connect to input
318 mgr->ConnectInput(task, 0, mgr->GetCommonInputContainer());
320 // create paths for the output in the common file
321 TString infoname(task->GetName());
322 TString histname(task->GetName());
323 infoname.Append("_Info");
324 histname.Append("_Hist");
325 AliAnalysisDataContainer *outputInfo = mgr->CreateContainer(infoname.Data(), TList::Class(), AliAnalysisManager::kOutputContainer, mgr->GetCommonFileName());
326 AliAnalysisDataContainer *outputHist = mgr->CreateContainer(histname.Data(), TList::Class(), AliAnalysisManager::kOutputContainer, mgr->GetCommonFileName());
327 mgr->ConnectOutput(task, 1, outputInfo);
328 mgr->ConnectOutput(task, 2, outputHist);