]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGJE/StrangenessInJets/AliAnalysisTaskV0sInJets.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Merge branch 'master' of https://git.cern.ch/reps/AliRoot
[u/mrichter/AliRoot.git] / PWGJE / StrangenessInJets / AliAnalysisTaskV0sInJets.cxx
CommitLineData
257f0eee 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
7729cc34 16// task for analysis of V0s (K0S, (anti-)Lambda) in charged jets
17// Author: Vit Kucera (vit.kucera@cern.ch)
18
19#include "TChain.h"
20#include "TTree.h"
21#include "TH1D.h"
22#include "TH2D.h"
23#include "THnSparse.h"
24#include "TCanvas.h"
25
26#include "AliAnalysisTask.h"
27#include "AliAnalysisManager.h"
28
29#include "AliESDEvent.h"
30#include "AliAODEvent.h"
31#include "AliAODTrack.h"
32#include <TDatabasePDG.h>
33#include <TPDGCode.h>
34#include "AliPIDResponse.h"
35#include "AliInputEventHandler.h"
36#include "AliAODMCHeader.h"
37#include "AliAODMCParticle.h"
38#include "TClonesArray.h"
39//#include "AliEventInfoObject.cxx"
40//#include "AliV0Object.cxx"
41//#include "AliJetObject.cxx"
42#include "TRandom3.h"
43
44#include "AliAnalysisTaskV0sInJets.h"
45
46ClassImp(AliAnalysisTaskV0sInJets)
47
48// upper edges of centrality bins
9706bde3 49//const Int_t AliAnalysisTaskV0sInJets::fgkiCentBinRanges[AliAnalysisTaskV0sInJets::fgkiNBinsCent] = {10, 30, 50, 80}; // Alice Zimmermann
7729cc34 50//const Int_t AliAnalysisTaskV0sInJets::fgkiCentBinRanges[AliAnalysisTaskV0sInJets::fgkiNBinsCent] = {10, 20, 40, 60, 80}; // Vit Kucera, initial binning
51//const Int_t AliAnalysisTaskV0sInJets::fgkiCentBinRanges[AliAnalysisTaskV0sInJets::fgkiNBinsCent] = {5, 10, 20, 40, 60, 80}; // Iouri Belikov, LF analysis
9706bde3 52const Int_t AliAnalysisTaskV0sInJets::fgkiCentBinRanges[AliAnalysisTaskV0sInJets::fgkiNBinsCent] = {10}; // only central
7729cc34 53
54// axis: pT of V0
9706bde3 55const Double_t AliAnalysisTaskV0sInJets::fgkdBinsPtV0[2] = {0, 12};
7729cc34 56const Int_t AliAnalysisTaskV0sInJets::fgkiNBinsPtV0 = sizeof(AliAnalysisTaskV0sInJets::fgkdBinsPtV0)/sizeof((AliAnalysisTaskV0sInJets::fgkdBinsPtV0)[0])-1;
57const Int_t AliAnalysisTaskV0sInJets::fgkiNBinsPtV0Init = int(((AliAnalysisTaskV0sInJets::fgkdBinsPtV0)[AliAnalysisTaskV0sInJets::fgkiNBinsPtV0]-(AliAnalysisTaskV0sInJets::fgkdBinsPtV0)[0])/0.1); // bin width 0.1 GeV/c
58// axis: pT of jets
59const Double_t AliAnalysisTaskV0sInJets::fgkdBinsPtJet[2] = {0, 100};
60const Int_t AliAnalysisTaskV0sInJets::fgkiNBinsPtJet = sizeof(AliAnalysisTaskV0sInJets::fgkdBinsPtJet)/sizeof(AliAnalysisTaskV0sInJets::fgkdBinsPtJet[0])-1;
61const Int_t AliAnalysisTaskV0sInJets::fgkiNBinsPtJetInit = int(((AliAnalysisTaskV0sInJets::fgkdBinsPtJet)[AliAnalysisTaskV0sInJets::fgkiNBinsPtJet]-(AliAnalysisTaskV0sInJets::fgkdBinsPtJet)[0])/5.); // bin width 5 GeV/c
62// axis: K0S invariant mass
9706bde3 63const Int_t AliAnalysisTaskV0sInJets::fgkiNBinsMassK0s = 300;
64const Double_t AliAnalysisTaskV0sInJets::fgkdMassK0sMin = 0.35;
65const Double_t AliAnalysisTaskV0sInJets::fgkdMassK0sMax = 0.65;
7729cc34 66// axis: Lambda invariant mass
67const Int_t AliAnalysisTaskV0sInJets::fgkiNBinsMassLambda = 200;
87e23d3e 68const Double_t AliAnalysisTaskV0sInJets::fgkdMassLambdaMin = 1.05;
69const Double_t AliAnalysisTaskV0sInJets::fgkdMassLambdaMax = 1.25;
7729cc34 70
71
72// Default constructor
73AliAnalysisTaskV0sInJets::AliAnalysisTaskV0sInJets():
74 AliAnalysisTaskSE(),
75 fAODIn(0),
76 fAODOut(0),
77 fOutputListStd(0),
78 fOutputListQA(0),
79 fOutputListCuts(0),
80 fOutputListMC(0),
81// ftreeOut(0),
82
83 fiAODAnalysis(1),
84
85 fdCutDCAToPrimVtxMin(0.1),
86 fdCutDCADaughtersMax(1.),
87 fdCutNSigmadEdxMax(3),
88 fdCutCPAMin(0.998),
89 fdCutNTauMax(5),
90
91 fsJetBranchName(0),
9706bde3 92 fsJetBgBranchName(0),
87e23d3e 93 fdCutPtJetMin(0),
94 fdCutPtTrackMin(5),
95 fdRadiusJet(0.4),
7729cc34 96 fbJetSelection(0),
97 fbMCAnalysis(0),
98// fbTreeOutput(0),
99 fRandom(0),
100
87e23d3e 101 fdCutVertexZ(10),
102 fdCutVertexR2(1),
103 fdCutCentLow(0),
104 fdCutCentHigh(80),
7729cc34 105
106/*
107 fBranchV0Rec(0),
108 fBranchV0Gen(0),
109 fBranchJet(0),
110 fEventInfo(0),
111*/
87e23d3e 112 fdCentrality(0),
7729cc34 113 fh1EventCounterCut(0),
114 fh1EventCent(0),
115 fh1EventCent2(0),
169867f8 116 fh1EventCent2Jets(0),
117 fh1EventCent2NoJets(0),
7729cc34 118 fh2EventCentTracks(0),
119 fh1V0CandPerEvent(0),
120 fh1NRndConeCent(0),
9706bde3 121 fh1NMedConeCent(0),
7729cc34 122 fh1AreaExcluded(0),
123
124 fh2CCK0s(0),
125 fh2CCLambda(0),
126 fh3CCMassCorrelBoth(0),
127 fh3CCMassCorrelKNotL(0),
128 fh3CCMassCorrelLNotK(0)
129{
130 for (Int_t i =0; i < fgkiNQAIndeces; i++)
131 {
132 fh1QAV0Status[i] = 0;
133 fh1QAV0TPCRefit[i] = 0;
134 fh1QAV0TPCRows[i] = 0;
135 fh1QAV0TPCFindable[i] = 0;
136 fh1QAV0TPCRowsFind[i] = 0;
137 fh1QAV0Eta[i] = 0;
138 fh2QAV0EtaRows[i] = 0;
139 fh2QAV0PtRows[i] = 0;
140 fh2QAV0PhiRows[i] = 0;
141 fh2QAV0NClRows[i] = 0;
142 fh2QAV0EtaNCl[i] = 0;
143
144 fh2QAV0EtaPtK0sPeak[i] = 0;
145 fh2QAV0EtaEtaK0s[i] = 0;
146 fh2QAV0PhiPhiK0s[i] = 0;
147 fh1QAV0RapK0s[i] = 0;
148 fh2QAV0PtPtK0sPeak[i] = 0;
149 fh2ArmPodK0s[i] = 0;
150
151 fh2QAV0EtaPtLambdaPeak[i] = 0;
152 fh2QAV0EtaEtaLambda[i] = 0;
153 fh2QAV0PhiPhiLambda[i] = 0;
154 fh1QAV0RapLambda[i] = 0;
155 fh2QAV0PtPtLambdaPeak[i] = 0;
156 fh2ArmPodLambda[i] = 0;
157
158 fh2QAV0EtaPtALambdaPeak[i] = 0;
159 fh2QAV0EtaEtaALambda[i] = 0;
160 fh2QAV0PhiPhiALambda[i] = 0;
161 fh1QAV0RapALambda[i] = 0;
162 fh2QAV0PtPtALambdaPeak[i] = 0;
163 fh2ArmPodALambda[i] = 0;
164
165 fh1QAV0Pt[i] = 0;
166 fh1QAV0Charge[i] = 0;
167 fh1QAV0DCAVtx[i] = 0;
168 fh1QAV0DCAV0[i] = 0;
169 fh1QAV0Cos[i] = 0;
170 fh1QAV0R[i] = 0;
171 fh1QACTau2D[i] = 0;
172 fh1QACTau3D[i] = 0;
173
174 fh2ArmPod[i] = 0;
175
176 fh2CutTPCRowsK0s[i] = 0;
177 fh2CutTPCRowsLambda[i] = 0;
178 fh2CutPtPosK0s[i] = 0;
179 fh2CutPtNegK0s[i] = 0;
180 fh2CutPtPosLambda[i] = 0;
181 fh2CutPtNegLambda[i] = 0;
182 fh2CutDCAVtx[i] = 0;
183 fh2CutDCAV0[i] = 0;
184 fh2CutCos[i] = 0;
185 fh2CutR[i] = 0;
186 fh2CutEtaK0s[i] = 0;
187 fh2CutEtaLambda[i] = 0;
188 fh2CutRapK0s[i] = 0;
189 fh2CutRapLambda[i] = 0;
190 fh2CutCTauK0s[i] = 0;
191 fh2CutCTauLambda[i] = 0;
192 fh2CutPIDPosK0s[i] = 0;
193 fh2CutPIDNegK0s[i] = 0;
194 fh2CutPIDPosLambda[i] = 0;
195 fh2CutPIDNegLambda[i] = 0;
196
197 fh2Tau3DVs2D[i] = 0;
198 }
199 for (Int_t i = 0; i<fgkiNCategV0; i++)
200 {
201 fh1V0InvMassK0sAll[i] = 0;
202 fh1V0InvMassLambdaAll[i] = 0;
203 fh1V0InvMassALambdaAll[i] = 0;
204 }
205 for (Int_t i = 0; i < fgkiNBinsCent; i++)
206 {
207 fh1EventCounterCutCent[i] = 0;
208 fh1V0CounterCentK0s[i] = 0;
209 fh1V0CounterCentLambda[i] = 0;
210 fh1V0CounterCentALambda[i] = 0;
211 fh1V0CandPerEventCentK0s[i] = 0;
212 fh1V0CandPerEventCentLambda[i] = 0;
213 fh1V0CandPerEventCentALambda[i] = 0;
214 fh1V0InvMassK0sCent[i] = 0;
215 fh1V0InvMassLambdaCent[i] = 0;
216 fh1V0InvMassALambdaCent[i] = 0;
217 fh1V0K0sPtMCGen[i] = 0;
218 fh2V0K0sPtMassMCRec[i] = 0;
219 fh1V0K0sPtMCRecFalse[i] = 0;
220 fh2V0K0sEtaPtMCGen[i] = 0;
221 fh3V0K0sEtaPtMassMCRec[i] = 0;
222 fh2V0K0sInJetPtMCGen[i] = 0;
223 fh3V0K0sInJetPtMassMCRec[i] = 0;
224 fh3V0K0sInJetEtaPtMCGen[i] = 0;
225 fh4V0K0sInJetEtaPtMassMCRec[i] = 0;
226 fh2V0K0sMCResolMPt[i] = 0;
227 fh2V0K0sMCPtGenPtRec[i] = 0;
228 fh1V0LambdaPtMCGen[i] = 0;
229 fh2V0LambdaPtMassMCRec[i] = 0;
230 fh1V0LambdaPtMCRecFalse[i] = 0;
231 fh2V0LambdaEtaPtMCGen[i] = 0;
232 fh3V0LambdaEtaPtMassMCRec[i] = 0;
233 fh2V0LambdaInJetPtMCGen[i] = 0;
234 fh3V0LambdaInJetPtMassMCRec[i] = 0;
235 fh3V0LambdaInJetEtaPtMCGen[i] = 0;
236 fh4V0LambdaInJetEtaPtMassMCRec[i] = 0;
237 fh2V0LambdaMCResolMPt[i] = 0;
238 fh2V0LambdaMCPtGenPtRec[i] = 0;
239 fhnV0LambdaInclMCFD[i] = 0;
240 fhnV0LambdaInJetsMCFD[i] = 0;
241 fhnV0LambdaBulkMCFD[i] = 0;
242 fh1V0XiPtMCGen[i] = 0;
243 fh1V0ALambdaPt[i] = 0;
244 fh1V0ALambdaPtMCGen[i] = 0;
245 fh1V0ALambdaPtMCRec[i] = 0;
246 fh2V0ALambdaPtMassMCRec[i] = 0;
247 fh1V0ALambdaPtMCRecFalse[i] = 0;
248 fh2V0ALambdaEtaPtMCGen[i] = 0;
249 fh3V0ALambdaEtaPtMassMCRec[i] = 0;
250 fh2V0ALambdaInJetPtMCGen[i] = 0;
251 fh2V0ALambdaInJetPtMCRec[i] = 0;
252 fh3V0ALambdaInJetPtMassMCRec[i] = 0;
253 fh3V0ALambdaInJetEtaPtMCGen[i] = 0;
254 fh4V0ALambdaInJetEtaPtMassMCRec[i] = 0;
255 fh2V0ALambdaMCResolMPt[i] = 0;
256 fh2V0ALambdaMCPtGenPtRec[i] = 0;
257 fhnV0ALambdaInclMCFD[i] = 0;
258 fhnV0ALambdaInJetsMCFD[i] = 0;
259 fhnV0ALambdaBulkMCFD[i] = 0;
260 fh1V0AXiPtMCGen[i] = 0;
261
262 // eta daughters
263// fhnV0K0sInclDaughterEtaPtPtMCGen[i] = 0;
264 fhnV0K0sInclDaughterEtaPtPtMCRec[i] = 0;
265// fhnV0K0sInJetsDaughterEtaPtPtMCGen[i] = 0;
266 fhnV0K0sInJetsDaughterEtaPtPtMCRec[i] = 0;
267// fhnV0LambdaInclDaughterEtaPtPtMCGen[i] = 0;
268 fhnV0LambdaInclDaughterEtaPtPtMCRec[i] = 0;
269// fhnV0LambdaInJetsDaughterEtaPtPtMCGen[i] = 0;
270 fhnV0LambdaInJetsDaughterEtaPtPtMCRec[i] = 0;
271// fhnV0ALambdaInclDaughterEtaPtPtMCGen[i] = 0;
272 fhnV0ALambdaInclDaughterEtaPtPtMCRec[i] = 0;
273// fhnV0ALambdaInJetsDaughterEtaPtPtMCGen[i] = 0;
274 fhnV0ALambdaInJetsDaughterEtaPtPtMCRec[i] = 0;
275
276 // Inclusive
277 fhnV0InclusiveK0s[i] = 0;
278 fhnV0InclusiveLambda[i] = 0;
279 fhnV0InclusiveALambda[i] = 0;
280 // Cones
281 fhnV0InJetK0s[i] = 0;
282 fhnV0InPerpK0s[i] = 0;
283 fhnV0InRndK0s[i] = 0;
9706bde3 284 fhnV0InMedK0s[i] = 0;
7729cc34 285 fhnV0OutJetK0s[i] = 0;
286 fhnV0NoJetK0s[i] = 0;
287 fhnV0InJetLambda[i] = 0;
288 fhnV0InPerpLambda[i] = 0;
289 fhnV0InRndLambda[i] = 0;
9706bde3 290 fhnV0InMedLambda[i] = 0;
7729cc34 291 fhnV0OutJetLambda[i] = 0;
292 fhnV0NoJetLambda[i] = 0;
293 fhnV0InJetALambda[i] = 0;
294 fhnV0InPerpALambda[i] = 0;
295 fhnV0InRndALambda[i] = 0;
9706bde3 296 fhnV0InMedALambda[i] = 0;
7729cc34 297 fhnV0OutJetALambda[i] = 0;
298 fhnV0NoJetALambda[i] = 0;
299
300 fh2V0PtJetAngleK0s[i] = 0;
301 fh2V0PtJetAngleLambda[i] = 0;
302 fh2V0PtJetAngleALambda[i] = 0;
303 fh1DCAInK0s[i] = 0;
304 fh1DCAInLambda[i] = 0;
305 fh1DCAInALambda[i] = 0;
306 fh1DCAOutK0s[i] = 0;
307 fh1DCAOutLambda[i] = 0;
308 fh1DCAOutALambda[i] = 0;
309 fh1DeltaZK0s[i] = 0;
310 fh1DeltaZLambda[i] = 0;
311 fh1DeltaZALambda[i] = 0;
312
313 fh1PtJet[i] = 0;
314 fh1EtaJet[i] = 0;
315 fh2EtaPtJet[i] = 0;
316 fh1PhiJet[i] = 0;
317 fh1NJetPerEvent[i] = 0;
318 fh2EtaPhiRndCone[i] = 0;
9706bde3 319 fh2EtaPhiMedCone[i] = 0;
7729cc34 320
321 fh1VtxZ[i] = 0;
322 fh2VtxXY[i] = 0;
323 }
324}
325
326// Constructor
327AliAnalysisTaskV0sInJets::AliAnalysisTaskV0sInJets(const char* name):
328 AliAnalysisTaskSE(name),
329 fAODIn(0),
330 fAODOut(0),
331 fOutputListStd(0),
332 fOutputListQA(0),
333 fOutputListCuts(0),
334 fOutputListMC(0),
335// ftreeOut(0),
336
337 fiAODAnalysis(1),
338
339 fdCutDCAToPrimVtxMin(0.1),
340 fdCutDCADaughtersMax(1.),
341 fdCutNSigmadEdxMax(3),
342 fdCutCPAMin(0.998),
343 fdCutNTauMax(5),
344
345 fsJetBranchName(0),
9706bde3 346 fsJetBgBranchName(0),
87e23d3e 347 fdCutPtJetMin(0),
348 fdCutPtTrackMin(5),
349 fdRadiusJet(0.4),
7729cc34 350 fbJetSelection(0),
351 fbMCAnalysis(0),
352// fbTreeOutput(0),
353 fRandom(0),
354
87e23d3e 355 fdCutVertexZ(10),
356 fdCutVertexR2(1),
357 fdCutCentLow(0),
358 fdCutCentHigh(80),
7729cc34 359/*
360 fBranchV0Rec(0),
361 fBranchV0Gen(0),
362 fBranchJet(0),
363 fEventInfo(0),
364*/
87e23d3e 365 fdCentrality(0),
7729cc34 366 fh1EventCounterCut(0),
367 fh1EventCent(0),
368 fh1EventCent2(0),
169867f8 369 fh1EventCent2Jets(0),
370 fh1EventCent2NoJets(0),
7729cc34 371 fh2EventCentTracks(0),
372 fh1V0CandPerEvent(0),
373 fh1NRndConeCent(0),
9706bde3 374 fh1NMedConeCent(0),
7729cc34 375 fh1AreaExcluded(0),
376
377 fh2CCK0s(0),
378 fh2CCLambda(0),
379 fh3CCMassCorrelBoth(0),
380 fh3CCMassCorrelKNotL(0),
381 fh3CCMassCorrelLNotK(0)
382{
383 for (Int_t i =0; i < fgkiNQAIndeces; i++)
384 {
385 fh1QAV0Status[i] = 0;
386 fh1QAV0TPCRefit[i] = 0;
387 fh1QAV0TPCRows[i] = 0;
388 fh1QAV0TPCFindable[i] = 0;
389 fh1QAV0TPCRowsFind[i] = 0;
390 fh1QAV0Eta[i] = 0;
391 fh2QAV0EtaRows[i] = 0;
392 fh2QAV0PtRows[i] = 0;
393 fh2QAV0PhiRows[i] = 0;
394 fh2QAV0NClRows[i] = 0;
395 fh2QAV0EtaNCl[i] = 0;
396
397 fh2QAV0EtaPtK0sPeak[i] = 0;
398 fh2QAV0EtaEtaK0s[i] = 0;
399 fh2QAV0PhiPhiK0s[i] = 0;
400 fh1QAV0RapK0s[i] = 0;
401 fh2QAV0PtPtK0sPeak[i] = 0;
402 fh2ArmPodK0s[i] = 0;
403
404 fh2QAV0EtaPtLambdaPeak[i] = 0;
405 fh2QAV0EtaEtaLambda[i] = 0;
406 fh2QAV0PhiPhiLambda[i] = 0;
407 fh1QAV0RapLambda[i] = 0;
408 fh2QAV0PtPtLambdaPeak[i] = 0;
409 fh2ArmPodLambda[i] = 0;
410
411 fh2QAV0EtaPtALambdaPeak[i] = 0;
412 fh2QAV0EtaEtaALambda[i] = 0;
413 fh2QAV0PhiPhiALambda[i] = 0;
414 fh1QAV0RapALambda[i] = 0;
415 fh2QAV0PtPtALambdaPeak[i] = 0;
416 fh2ArmPodALambda[i] = 0;
417
418 fh1QAV0Pt[i] = 0;
419 fh1QAV0Charge[i] = 0;
420 fh1QAV0DCAVtx[i] = 0;
421 fh1QAV0DCAV0[i] = 0;
422 fh1QAV0Cos[i] = 0;
423 fh1QAV0R[i] = 0;
424 fh1QACTau2D[i] = 0;
425 fh1QACTau3D[i] = 0;
426
427 fh2ArmPod[i] = 0;
428
429 fh2CutTPCRowsK0s[i] = 0;
430 fh2CutTPCRowsLambda[i] = 0;
431 fh2CutPtPosK0s[i] = 0;
432 fh2CutPtNegK0s[i] = 0;
433 fh2CutPtPosLambda[i] = 0;
434 fh2CutPtNegLambda[i] = 0;
435 fh2CutDCAVtx[i] = 0;
436 fh2CutDCAV0[i] = 0;
437 fh2CutCos[i] = 0;
438 fh2CutR[i] = 0;
439 fh2CutEtaK0s[i] = 0;
440 fh2CutEtaLambda[i] = 0;
441 fh2CutRapK0s[i] = 0;
442 fh2CutRapLambda[i] = 0;
443 fh2CutCTauK0s[i] = 0;
444 fh2CutCTauLambda[i] = 0;
445 fh2CutPIDPosK0s[i] = 0;
446 fh2CutPIDNegK0s[i] = 0;
447 fh2CutPIDPosLambda[i] = 0;
448 fh2CutPIDNegLambda[i] = 0;
449
450 fh2Tau3DVs2D[i] = 0;
451 }
452 for (Int_t i = 0; i<fgkiNCategV0; i++)
453 {
454 fh1V0InvMassK0sAll[i] = 0;
455 fh1V0InvMassLambdaAll[i] = 0;
456 fh1V0InvMassALambdaAll[i] = 0;
457 }
458 for (Int_t i = 0; i < fgkiNBinsCent; i++)
459 {
460 fh1EventCounterCutCent[i] = 0;
461 fh1V0CounterCentK0s[i] = 0;
462 fh1V0CounterCentLambda[i] = 0;
463 fh1V0CounterCentALambda[i] = 0;
464 fh1V0CandPerEventCentK0s[i] = 0;
465 fh1V0CandPerEventCentLambda[i] = 0;
466 fh1V0CandPerEventCentALambda[i] = 0;
467 fh1V0InvMassK0sCent[i] = 0;
468 fh1V0InvMassLambdaCent[i] = 0;
469 fh1V0InvMassALambdaCent[i] = 0;
470 fh1V0K0sPtMCGen[i] = 0;
471 fh2V0K0sPtMassMCRec[i] = 0;
472 fh1V0K0sPtMCRecFalse[i] = 0;
473 fh2V0K0sEtaPtMCGen[i] = 0;
474 fh3V0K0sEtaPtMassMCRec[i] = 0;
475 fh2V0K0sInJetPtMCGen[i] = 0;
476 fh3V0K0sInJetPtMassMCRec[i] = 0;
477 fh3V0K0sInJetEtaPtMCGen[i] = 0;
478 fh4V0K0sInJetEtaPtMassMCRec[i] = 0;
479 fh2V0K0sMCResolMPt[i] = 0;
480 fh2V0K0sMCPtGenPtRec[i] = 0;
481 fh1V0LambdaPtMCGen[i] = 0;
482 fh2V0LambdaPtMassMCRec[i] = 0;
483 fh1V0LambdaPtMCRecFalse[i] = 0;
484 fh2V0LambdaEtaPtMCGen[i] = 0;
485 fh3V0LambdaEtaPtMassMCRec[i] = 0;
486 fh2V0LambdaInJetPtMCGen[i] = 0;
487 fh3V0LambdaInJetPtMassMCRec[i] = 0;
488 fh3V0LambdaInJetEtaPtMCGen[i] = 0;
489 fh4V0LambdaInJetEtaPtMassMCRec[i] = 0;
490 fh2V0LambdaMCResolMPt[i] = 0;
491 fh2V0LambdaMCPtGenPtRec[i] = 0;
492 fhnV0LambdaInclMCFD[i] = 0;
493 fhnV0LambdaInJetsMCFD[i] = 0;
494 fhnV0LambdaBulkMCFD[i] = 0;
495 fh1V0XiPtMCGen[i] = 0;
496 fh1V0ALambdaPt[i] = 0;
497 fh1V0ALambdaPtMCGen[i] = 0;
498 fh1V0ALambdaPtMCRec[i] = 0;
499 fh2V0ALambdaPtMassMCRec[i] = 0;
500 fh1V0ALambdaPtMCRecFalse[i] = 0;
501 fh2V0ALambdaEtaPtMCGen[i] = 0;
502 fh3V0ALambdaEtaPtMassMCRec[i] = 0;
503 fh2V0ALambdaInJetPtMCGen[i] = 0;
504 fh2V0ALambdaInJetPtMCRec[i] = 0;
505 fh3V0ALambdaInJetPtMassMCRec[i] = 0;
506 fh3V0ALambdaInJetEtaPtMCGen[i] = 0;
507 fh4V0ALambdaInJetEtaPtMassMCRec[i] = 0;
508 fh2V0ALambdaMCResolMPt[i] = 0;
509 fh2V0ALambdaMCPtGenPtRec[i] = 0;
510 fhnV0ALambdaInclMCFD[i] = 0;
511 fhnV0ALambdaInJetsMCFD[i] = 0;
512 fhnV0ALambdaBulkMCFD[i] = 0;
513 fh1V0AXiPtMCGen[i] = 0;
514
515 // eta daughters
516// fhnV0K0sInclDaughterEtaPtPtMCGen[i] = 0;
517 fhnV0K0sInclDaughterEtaPtPtMCRec[i] = 0;
518// fhnV0K0sInJetsDaughterEtaPtPtMCGen[i] = 0;
519 fhnV0K0sInJetsDaughterEtaPtPtMCRec[i] = 0;
520// fhnV0LambdaInclDaughterEtaPtPtMCGen[i] = 0;
521 fhnV0LambdaInclDaughterEtaPtPtMCRec[i] = 0;
522// fhnV0LambdaInJetsDaughterEtaPtPtMCGen[i] = 0;
523 fhnV0LambdaInJetsDaughterEtaPtPtMCRec[i] = 0;
524// fhnV0ALambdaInclDaughterEtaPtPtMCGen[i] = 0;
525 fhnV0ALambdaInclDaughterEtaPtPtMCRec[i] = 0;
526// fhnV0ALambdaInJetsDaughterEtaPtPtMCGen[i] = 0;
527 fhnV0ALambdaInJetsDaughterEtaPtPtMCRec[i] = 0;
528
529 // Inclusive
530 fhnV0InclusiveK0s[i] = 0;
531 fhnV0InclusiveLambda[i] = 0;
532 fhnV0InclusiveALambda[i] = 0;
533 // Cones
534 fhnV0InJetK0s[i] = 0;
535 fhnV0InPerpK0s[i] = 0;
536 fhnV0InRndK0s[i] = 0;
9706bde3 537 fhnV0InMedK0s[i] = 0;
7729cc34 538 fhnV0OutJetK0s[i] = 0;
539 fhnV0NoJetK0s[i] = 0;
540 fhnV0InJetLambda[i] = 0;
541 fhnV0InPerpLambda[i] = 0;
542 fhnV0InRndLambda[i] = 0;
9706bde3 543 fhnV0InMedLambda[i] = 0;
7729cc34 544 fhnV0OutJetLambda[i] = 0;
545 fhnV0NoJetLambda[i] = 0;
546 fhnV0InJetALambda[i] = 0;
547 fhnV0InPerpALambda[i] = 0;
548 fhnV0InRndALambda[i] = 0;
9706bde3 549 fhnV0InMedALambda[i] = 0;
7729cc34 550 fhnV0OutJetALambda[i] = 0;
551 fhnV0NoJetALambda[i] = 0;
552
553 fh2V0PtJetAngleK0s[i] = 0;
554 fh2V0PtJetAngleLambda[i] = 0;
555 fh2V0PtJetAngleALambda[i] = 0;
556 fh1DCAInK0s[i] = 0;
557 fh1DCAInLambda[i] = 0;
558 fh1DCAInALambda[i] = 0;
559 fh1DCAOutK0s[i] = 0;
560 fh1DCAOutLambda[i] = 0;
561 fh1DCAOutALambda[i] = 0;
562 fh1DeltaZK0s[i] = 0;
563 fh1DeltaZLambda[i] = 0;
564 fh1DeltaZALambda[i] = 0;
565
566 fh1PtJet[i] = 0;
567 fh1EtaJet[i] = 0;
568 fh2EtaPtJet[i] = 0;
569 fh1PhiJet[i] = 0;
570 fh1NJetPerEvent[i] = 0;
571 fh2EtaPhiRndCone[i] = 0;
9706bde3 572 fh2EtaPhiMedCone[i] = 0;
7729cc34 573
574 fh1VtxZ[i] = 0;
575 fh2VtxXY[i] = 0;
576 }
577 // Define input and output slots here
578 // Input slot #0 works with a TChain
579 DefineInput(0, TChain::Class());
580 // Output slot #0 id reserved by the base class for AOD
581 // Output slot #1 writes into a TList container
582 DefineOutput(1, TList::Class());
583 DefineOutput(2, TList::Class());
584 DefineOutput(3, TList::Class());
585 DefineOutput(4, TList::Class());
586 DefineOutput(5, TTree::Class());
587}
588
589AliAnalysisTaskV0sInJets::~AliAnalysisTaskV0sInJets()
590{
9706bde3 591 /*
592 if (fBranchV0Rec)
593 fBranchV0Rec->Delete();
594 delete fBranchV0Rec;
595 fBranchV0Rec = 0;
596 if (fBranchV0Gen)
597 fBranchV0Gen->Delete();
598 delete fBranchV0Gen;
599 fBranchV0Gen = 0;
600 if (fBranchJet)
601 fBranchJet->Delete();
602 delete fBranchJet;
603 fBranchJet = 0;
604 if (fEventInfo)
605 fEventInfo->Delete();
606 delete fEventInfo;
607 fEventInfo = 0;
608 */
7729cc34 609 delete fRandom;
610 fRandom = 0;
611}
612
613void AliAnalysisTaskV0sInJets::UserCreateOutputObjects()
614{
615 // Create histograms
616 // Called once
617
618 fRandom = new TRandom3(0);
619
9706bde3 620 /*
621 if (!fBranchV0Rec && fbTreeOutput)
622 {
623 // fBranchV0Rec = new TClonesArray("AliAODv0",0);
624 fBranchV0Rec = new TClonesArray("AliV0Object",0);
625 fBranchV0Rec->SetName("branch_V0Rec");
626 }
627 if (!fBranchV0Gen && fbTreeOutput)
628 {
629 fBranchV0Gen = new TClonesArray("AliAODMCParticle",0);
630 fBranchV0Gen->SetName("branch_V0Gen");
631 }
632 if (!fBranchJet && fbTreeOutput)
633 {
634 // fBranchJet = new TClonesArray("AliAODJet",0);
635 fBranchJet = new TClonesArray("AliJetObject",0);
636 fBranchJet->SetName("branch_Jet");
637 }
638 if (!fEventInfo && fbTreeOutput)
639 {
640 fEventInfo = new AliEventInfoObject();
641 fEventInfo->SetName("eventInfo");
642 }
643 Int_t dSizeBuffer = 32000; // default 32000
644 if (fbTreeOutput)
7729cc34 645 {
9706bde3 646 ftreeOut = new TTree("treeV0","Tree V0");
647 ftreeOut->Branch("branch_V0Rec",&fBranchV0Rec,dSizeBuffer,2);
648 ftreeOut->Branch("branch_V0Gen",&fBranchV0Gen,dSizeBuffer,2);
649 ftreeOut->Branch("branch_Jet",&fBranchJet,dSizeBuffer,2);
650 ftreeOut->Branch("eventInfo",&fEventInfo,dSizeBuffer,2);
7729cc34 651 }
9706bde3 652 */
7729cc34 653
654 fOutputListStd = new TList();
655 fOutputListStd->SetOwner();
656 fOutputListQA = new TList();
657 fOutputListQA->SetOwner();
658 fOutputListCuts = new TList();
659 fOutputListCuts->SetOwner();
660 fOutputListMC = new TList();
661 fOutputListMC->SetOwner();
662
663 // event categories
664 const Int_t iNCategEvent = 6;
665 TString categEvent[iNCategEvent] = {"coll. candid.","AOD OK","vtx & cent","with V0","with jets","jet selection"};
666 // labels for stages of V0 selection
667 TString categV0[fgkiNCategV0] = {"all"/*0*/,"mass range"/*1*/,"rec. method"/*2*/,"tracks TPC"/*3*/,"track pt"/*4*/,"DCA prim v"/*5*/,"DCA daughters"/*6*/,"CPA"/*7*/,"volume"/*8*/,"track #it{#eta}"/*9*/,"V0 #it{y} & #it{#eta}"/*10*/,"lifetime"/*11*/,"PID"/*12*/,"Arm.-Pod."/*13*/,"inclusive"/*14*/,"in jet event"/*15*/,"in jet"/*16*/};
668
669 fh1EventCounterCut = new TH1D("fh1EventCounterCut","Number of events after filtering;selection filter;counts",iNCategEvent,0,iNCategEvent);
670 for (Int_t i = 0; i < iNCategEvent; i++)
671 fh1EventCounterCut->GetXaxis()->SetBinLabel(i+1,categEvent[i].Data());
672 fh1EventCent2 = new TH1D("fh1EventCent2","Number of events vs centrality;centrality;counts",100,0,100);
169867f8 673 fh1EventCent2Jets = new TH1D("fh1EventCent2Jets","Number of sel.-jet events vs centrality;centrality;counts",100,0,100);
674 fh1EventCent2NoJets = new TH1D("fh1EventCent2NoJets","Number of no-jet events vs centrality;centrality;counts",100,0,100);
7729cc34 675 fh2EventCentTracks = new TH2D("fh2EventCentTracks","Number of tracks vs centrality;centrality;tracks;counts",100,0,100,150,0,15e3);
676 fh1EventCent = new TH1D("fh1EventCent","Number of events in centrality bins;centrality;counts",fgkiNBinsCent,0,fgkiNBinsCent);
677 for (Int_t i = 0; i < fgkiNBinsCent; i++)
678 fh1EventCent->GetXaxis()->SetBinLabel(i+1,GetCentBinLabel(i).Data());
679 fh1NRndConeCent = new TH1D("fh1NRndConeCent","Number of rnd. cones in centrality bins;centrality;counts",fgkiNBinsCent,0,fgkiNBinsCent);
680 for (Int_t i = 0; i < fgkiNBinsCent; i++)
681 fh1NRndConeCent->GetXaxis()->SetBinLabel(i+1,GetCentBinLabel(i).Data());
9706bde3 682 fh1NMedConeCent = new TH1D("fh1NMedConeCent","Number of med.-cl. cones in centrality bins;centrality;counts",fgkiNBinsCent,0,fgkiNBinsCent);
683 for (Int_t i = 0; i < fgkiNBinsCent; i++)
684 fh1NMedConeCent->GetXaxis()->SetBinLabel(i+1,GetCentBinLabel(i).Data());
7729cc34 685 fh1AreaExcluded = new TH1D("fh1AreaExcluded","Area of excluded cones in centrality bins;centrality;area",fgkiNBinsCent,0,fgkiNBinsCent);
686 for (Int_t i = 0; i < fgkiNBinsCent; i++)
687 fh1AreaExcluded->GetXaxis()->SetBinLabel(i+1,GetCentBinLabel(i).Data());
688 fOutputListStd->Add(fh1EventCounterCut);
689 fOutputListStd->Add(fh1EventCent);
690 fOutputListStd->Add(fh1EventCent2);
169867f8 691 fOutputListStd->Add(fh1EventCent2Jets);
692 fOutputListStd->Add(fh1EventCent2NoJets);
7729cc34 693 fOutputListStd->Add(fh1NRndConeCent);
9706bde3 694 fOutputListStd->Add(fh1NMedConeCent);
7729cc34 695 fOutputListStd->Add(fh1AreaExcluded);
696 fOutputListStd->Add(fh2EventCentTracks);
697
698 fh1V0CandPerEvent = new TH1D("fh1V0CandPerEvent","Number of all V0 candidates per event;candidates;events",1000,0,1000);
699 fOutputListStd->Add(fh1V0CandPerEvent);
700
701 for (Int_t i = 0; i < fgkiNBinsCent; i++)
702 {
703 fh1EventCounterCutCent[i] = new TH1D(Form("fh1EventCounterCutCent_%d",i),Form("Number of events after filtering, cent %s;selection filter;counts",GetCentBinLabel(i).Data()),iNCategEvent,0,iNCategEvent);
704 for (Int_t j = 0; j < iNCategEvent; j++)
705 fh1EventCounterCutCent[i]->GetXaxis()->SetBinLabel(j+1,categEvent[j].Data());
706 fh1V0CandPerEventCentK0s[i] = new TH1D(Form("fh1V0CandPerEventCentK0s_%d",i),Form("Number of selected K0s candidates per event, cent %s;candidates;events",GetCentBinLabel(i).Data()),100,0,100);
707 fh1V0CandPerEventCentLambda[i] = new TH1D(Form("fh1V0CandPerEventCentLambda_%d",i),Form("Number of selected Lambda candidates per event, cent %s;candidates;events",GetCentBinLabel(i).Data()),100,0,100);
708 fh1V0CandPerEventCentALambda[i] = new TH1D(Form("fh1V0CandPerEventCentALambda_%d",i),Form("Number of selected ALambda candidates per event, cent %s;candidates;events",GetCentBinLabel(i).Data()),100,0,100);
709 fh1V0CounterCentK0s[i] = new TH1D(Form("fh1V0CounterCentK0s_%d",i),Form("Number of K0s candidates after cuts, cent %s;cut;counts",GetCentBinLabel(i).Data()),fgkiNCategV0,0,fgkiNCategV0);
710 fh1V0CounterCentLambda[i] = new TH1D(Form("fh1V0CounterCentLambda_%d",i),Form("Number of Lambda candidates after cuts, cent %s;cut;counts",GetCentBinLabel(i).Data()),fgkiNCategV0,0,fgkiNCategV0);
711 fh1V0CounterCentALambda[i] = new TH1D(Form("fh1V0CounterCentALambda_%d",i),Form("Number of ALambda candidates after cuts, cent %s;cut;counts",GetCentBinLabel(i).Data()),fgkiNCategV0,0,fgkiNCategV0);
712 for (Int_t j = 0; j < fgkiNCategV0; j++)
713 {
714 fh1V0CounterCentK0s[i]->GetXaxis()->SetBinLabel(j+1,categV0[j].Data());
715 fh1V0CounterCentLambda[i]->GetXaxis()->SetBinLabel(j+1,categV0[j].Data());
716 fh1V0CounterCentALambda[i]->GetXaxis()->SetBinLabel(j+1,categV0[j].Data());
717 }
718 fOutputListStd->Add(fh1EventCounterCutCent[i]);
719 fOutputListStd->Add(fh1V0CandPerEventCentK0s[i]);
720 fOutputListStd->Add(fh1V0CandPerEventCentLambda[i]);
721 fOutputListStd->Add(fh1V0CandPerEventCentALambda[i]);
722 fOutputListStd->Add(fh1V0CounterCentK0s[i]);
723 fOutputListStd->Add(fh1V0CounterCentLambda[i]);
724 fOutputListStd->Add(fh1V0CounterCentALambda[i]);
725 }
726 // pt binning for V0 and jets
727 Int_t iNBinsPtV0 = fgkiNBinsPtV0Init;
87e23d3e 728 Double_t dPtV0Min = fgkdBinsPtV0[0];
729 Double_t dPtV0Max = fgkdBinsPtV0[fgkiNBinsPtV0];
7729cc34 730 Int_t iNJetPtBins = fgkiNBinsPtJetInit;
87e23d3e 731 Double_t dJetPtMin = fgkdBinsPtJet[0];
732 Double_t dJetPtMax = fgkdBinsPtJet[fgkiNBinsPtJet];
7729cc34 733
87e23d3e 734 fh2CCK0s = new TH2D("fh2CCK0s","K0s candidates in Lambda peak",fgkiNBinsMassK0s,fgkdMassK0sMin,fgkdMassK0sMax,iNBinsPtV0,dPtV0Min,dPtV0Max);
735 fh2CCLambda = new TH2D("fh2CCLambda","Lambda candidates in K0s peak",fgkiNBinsMassLambda,fgkdMassLambdaMin,fgkdMassLambdaMax,iNBinsPtV0,dPtV0Min,dPtV0Max);
7729cc34 736 Int_t binsCorrel[3] = {fgkiNBinsMassK0s, fgkiNBinsMassLambda, iNBinsPtV0};
87e23d3e 737 Double_t xminCorrel[3] = {fgkdMassK0sMin, fgkdMassLambdaMin, dPtV0Min};
738 Double_t xmaxCorrel[3] = {fgkdMassK0sMax, fgkdMassLambdaMax, dPtV0Max};
7729cc34 739// Int_t binsCorrel[3] = {200, 200, iNBinsPtV0};
87e23d3e 740// Double_t xminCorrel[3] = {0, 0, dPtV0Min};
741// Double_t xmaxCorrel[3] = {2, 2, dPtV0Max};
7729cc34 742 fh3CCMassCorrelBoth = new THnSparseD("fh3CCMassCorrelBoth","Mass correlation: K0S && Lambda;m K0S;m Lambda;pT",3,binsCorrel,xminCorrel,xmaxCorrel);
743 fh3CCMassCorrelKNotL = new THnSparseD("fh3CCMassCorrelKNotL","Mass correlation: K0S, not Lambda;m K0S;m Lambda;pT",3,binsCorrel,xminCorrel,xmaxCorrel);
744 fh3CCMassCorrelLNotK = new THnSparseD("fh3CCMassCorrelLNotK","Mass correlation: Lambda, not K0S;m K0S;m Lambda;pT",3,binsCorrel,xminCorrel,xmaxCorrel);
745 fOutputListQA->Add(fh2CCK0s);
746 fOutputListQA->Add(fh2CCLambda);
747 fOutputListQA->Add(fh3CCMassCorrelBoth);
748 fOutputListQA->Add(fh3CCMassCorrelKNotL);
749 fOutputListQA->Add(fh3CCMassCorrelLNotK);
750
87e23d3e 751 Double_t dStepEtaV0 = 0.025;
752 Double_t dRangeEtaV0Max = 0.8;
753 const Int_t iNBinsEtaV0 = 2*Int_t(dRangeEtaV0Max/dStepEtaV0);
7729cc34 754 // inclusive
755 const Int_t iNDimIncl = 3;
756 Int_t binsKIncl[iNDimIncl] = {fgkiNBinsMassK0s, iNBinsPtV0, iNBinsEtaV0};
87e23d3e 757 Double_t xminKIncl[iNDimIncl] = {fgkdMassK0sMin, dPtV0Min, -dRangeEtaV0Max};
758 Double_t xmaxKIncl[iNDimIncl] = {fgkdMassK0sMax, dPtV0Max, dRangeEtaV0Max};
7729cc34 759 Int_t binsLIncl[iNDimIncl] = {fgkiNBinsMassLambda, iNBinsPtV0, iNBinsEtaV0};
87e23d3e 760 Double_t xminLIncl[iNDimIncl] = {fgkdMassLambdaMin, dPtV0Min, -dRangeEtaV0Max};
761 Double_t xmaxLIncl[iNDimIncl] = {fgkdMassLambdaMax, dPtV0Max, dRangeEtaV0Max};
7729cc34 762 // binning in jets
763 const Int_t iNDimInJC = 4;
764 Int_t binsKInJC[iNDimInJC] = {fgkiNBinsMassK0s, iNBinsPtV0, iNBinsEtaV0, iNJetPtBins};
87e23d3e 765 Double_t xminKInJC[iNDimInJC] = {fgkdMassK0sMin, dPtV0Min, -dRangeEtaV0Max, dJetPtMin};
766 Double_t xmaxKInJC[iNDimInJC] = {fgkdMassK0sMax, dPtV0Max, dRangeEtaV0Max, dJetPtMax};
7729cc34 767 Int_t binsLInJC[iNDimInJC] = {fgkiNBinsMassLambda, iNBinsPtV0, iNBinsEtaV0, iNJetPtBins};
87e23d3e 768 Double_t xminLInJC[iNDimInJC] = {fgkdMassLambdaMin, dPtV0Min, -dRangeEtaV0Max, dJetPtMin};
769 Double_t xmaxLInJC[iNDimInJC] = {fgkdMassLambdaMax, dPtV0Max, dRangeEtaV0Max, dJetPtMax};
7729cc34 770
771 // binning eff inclusive vs eta-pT
87e23d3e 772 Double_t dStepDeltaEta = 0.1;
773 Double_t dRangeDeltaEtaMax = 0.5;
774 const Int_t iNBinsDeltaEta = 2*Int_t(dRangeDeltaEtaMax/dStepDeltaEta);
7729cc34 775 Int_t binsEtaK[3] = {fgkiNBinsMassK0s, iNBinsPtV0, iNBinsEtaV0};
87e23d3e 776 Double_t xminEtaK[3] = {fgkdMassK0sMin, dPtV0Min, -dRangeEtaV0Max};
777 Double_t xmaxEtaK[3] = {fgkdMassK0sMax, dPtV0Max, dRangeEtaV0Max};
7729cc34 778 Int_t binsEtaL[3] = {fgkiNBinsMassLambda, iNBinsPtV0, iNBinsEtaV0};
87e23d3e 779 Double_t xminEtaL[3] = {fgkdMassLambdaMin, dPtV0Min, -dRangeEtaV0Max};
780 Double_t xmaxEtaL[3] = {fgkdMassLambdaMax, dPtV0Max, dRangeEtaV0Max};
7729cc34 781 // binning eff in jets vs eta-pT
782 // associated
783 Int_t binsEtaKInRec[5] = {fgkiNBinsMassK0s, iNBinsPtV0, iNBinsEtaV0, iNJetPtBins, iNBinsDeltaEta};
87e23d3e 784 Double_t xminEtaKInRec[5] = {fgkdMassK0sMin, dPtV0Min, -dRangeEtaV0Max, dJetPtMin, -dRangeDeltaEtaMax};
785 Double_t xmaxEtaKInRec[5] = {fgkdMassK0sMax, dPtV0Max, dRangeEtaV0Max, dJetPtMax, dRangeDeltaEtaMax};
7729cc34 786 Int_t binsEtaLInRec[5] = {fgkiNBinsMassLambda, iNBinsPtV0, iNBinsEtaV0, iNJetPtBins, iNBinsDeltaEta};
87e23d3e 787 Double_t xminEtaLInRec[5] = {fgkdMassLambdaMin, dPtV0Min, -dRangeEtaV0Max, dJetPtMin, -dRangeDeltaEtaMax};
788 Double_t xmaxEtaLInRec[5] = {fgkdMassLambdaMax, dPtV0Max, dRangeEtaV0Max, dJetPtMax, dRangeDeltaEtaMax};
7729cc34 789 // generated
790 Int_t binsEtaInGen[4] = {iNBinsPtV0, iNBinsEtaV0, iNJetPtBins, iNBinsDeltaEta};
87e23d3e 791 Double_t xminEtaInGen[4] = {dPtV0Min, -dRangeEtaV0Max, dJetPtMin, -dRangeDeltaEtaMax};
792 Double_t xmaxEtaInGen[4] = {dPtV0Max, dRangeEtaV0Max, dJetPtMax, dRangeDeltaEtaMax};
7729cc34 793 // daughter eta: charge-etaD-ptD-etaV0-ptV0-ptJet
794 const Int_t iNDimEtaD = 6;
795 Int_t binsEtaDaughter[iNDimEtaD] = {2, 20, iNBinsPtV0, iNBinsEtaV0, iNBinsPtV0, iNJetPtBins};
87e23d3e 796 Double_t xminEtaDaughter[iNDimEtaD] = {0, -1, dPtV0Min, -dRangeEtaV0Max, dPtV0Min, dJetPtMin};
797 Double_t xmaxEtaDaughter[iNDimEtaD] = {2, 1, dPtV0Max, dRangeEtaV0Max, dPtV0Max, dJetPtMax};
7729cc34 798
799 for (Int_t i = 0; i < fgkiNBinsCent; i++)
800 {
87e23d3e 801 fh1V0InvMassK0sCent[i] = new TH1D(Form("fh1V0InvMassK0sCent_%d",i),Form("K0s: V0 invariant mass, cent %s;#it{m}_{inv} (GeV/#it{c}^{2});counts",GetCentBinLabel(i).Data()),fgkiNBinsMassK0s,fgkdMassK0sMin,fgkdMassK0sMax);
802 fh1V0InvMassLambdaCent[i] = new TH1D(Form("fh1V0InvMassLambdaCent_%d",i),Form("Lambda: V0 invariant mass, cent %s;#it{m}_{inv} (GeV/#it{c}^{2});counts",GetCentBinLabel(i).Data()),fgkiNBinsMassLambda,fgkdMassLambdaMin,fgkdMassLambdaMax);
803 fh1V0InvMassALambdaCent[i] = new TH1D(Form("fh1V0InvMassALambdaCent_%d",i),Form("ALambda: V0 invariant mass, cent %s;#it{m}_{inv} (GeV/#it{c}^{2});counts",GetCentBinLabel(i).Data()),fgkiNBinsMassLambda,fgkdMassLambdaMin,fgkdMassLambdaMax);
7729cc34 804 fOutputListStd->Add(fh1V0InvMassK0sCent[i]);
805 fOutputListStd->Add(fh1V0InvMassLambdaCent[i]);
806 fOutputListStd->Add(fh1V0InvMassALambdaCent[i]);
807 // Inclusive
808 fhnV0InclusiveK0s[i] = new THnSparseD(Form("fhnV0InclusiveK0s_C%d",i), "K0s: V0 invariant mass vs pt;#it{m}_{inv} (GeV/#it{c}^{2});pt (GeV/#it{c});counts",iNDimIncl,binsKIncl,xminKIncl,xmaxKIncl);
809 fhnV0InclusiveLambda[i] = new THnSparseD(Form("fhnV0InclusiveLambda_C%d",i), "Lambda: V0 invariant mass vs pt;#it{m}_{inv} (GeV/#it{c}^{2});pt (GeV/#it{c});counts",iNDimIncl,binsLIncl,xminLIncl,xmaxLIncl);
810 fhnV0InclusiveALambda[i] = new THnSparseD(Form("fhnV0InclusiveALambda_C%d",i), "ALambda: V0 invariant mass vs pt;#it{m}_{inv} (GeV/#it{c}^{2});pt (GeV/#it{c});counts",iNDimIncl,binsLIncl,xminLIncl,xmaxLIncl);
811 fOutputListStd->Add(fhnV0InclusiveK0s[i]);
812 fOutputListStd->Add(fhnV0InclusiveLambda[i]);
813 fOutputListStd->Add(fhnV0InclusiveALambda[i]);
814 // In cones
815 fhnV0InJetK0s[i] = new THnSparseD(Form("fhnV0InJetK0s_%d",i),Form("K0s: Mass vs Pt in jets, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});#it{p}_{T}^{V0} (GeV/#it{c});#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimInJC,binsKInJC,xminKInJC,xmaxKInJC);
816 fOutputListStd->Add(fhnV0InJetK0s[i]);
817 fhnV0InPerpK0s[i] = new THnSparseD(Form("fhnV0InPerpK0s_%d",i),Form("K0s: Mass vs Pt in perp. cones, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});#it{p}_{T}^{V0} (GeV/#it{c});#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimInJC,binsKInJC,xminKInJC,xmaxKInJC);
818 fOutputListStd->Add(fhnV0InPerpK0s[i]);
819 fhnV0InRndK0s[i] = new THnSparseD(Form("fhnV0InRndK0s_%d",i),Form("K0s: Mass vs Pt in rnd. cones, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});#it{p}_{T}^{V0} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimIncl,binsKIncl,xminKIncl,xmaxKIncl);
820 fOutputListStd->Add(fhnV0InRndK0s[i]);
9706bde3 821 fhnV0InMedK0s[i] = new THnSparseD(Form("fhnV0InMedK0s_%d",i),Form("K0s: Mass vs Pt in med.-cl. cones, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});#it{p}_{T}^{V0} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimIncl,binsKIncl,xminKIncl,xmaxKIncl);
822 fOutputListStd->Add(fhnV0InMedK0s[i]);
7729cc34 823 fhnV0OutJetK0s[i] = new THnSparseD(Form("fhnV0OutJetK0s_%d",i),Form("K0s: Pt outside jets, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});#it{p}_{T}^{V0} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimIncl,binsKIncl,xminKIncl,xmaxKIncl);
824 fOutputListStd->Add(fhnV0OutJetK0s[i]);
825 fhnV0NoJetK0s[i] = new THnSparseD(Form("fhnV0NoJetK0s_%d",i),Form("K0s: Pt in jet-less events, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});#it{p}_{T}^{V0} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimIncl,binsKIncl,xminKIncl,xmaxKIncl);
826 fOutputListStd->Add(fhnV0NoJetK0s[i]);
827 fhnV0InJetLambda[i] = new THnSparseD(Form("fhnV0InJetLambda_%d",i),Form("Lambda: Mass vs Pt in jets, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});#it{p}_{T}^{V0} (GeV/#it{c});#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimInJC,binsLInJC,xminLInJC,xmaxLInJC);
828 fOutputListStd->Add(fhnV0InJetLambda[i]);
829 fhnV0InPerpLambda[i] = new THnSparseD(Form("fhnV0InPerpLambda_%d",i),Form("Lambda: Mass vs Pt in perp. cones, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});#it{p}_{T}^{V0} (GeV/#it{c});#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimInJC,binsLInJC,xminLInJC,xmaxLInJC);
830 fOutputListStd->Add(fhnV0InPerpLambda[i]);
831 fhnV0InRndLambda[i] = new THnSparseD(Form("fhnV0InRndLambda_%d",i),Form("Lambda: Mass vs Pt in rnd. cones, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});#it{p}_{T}^{V0} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimIncl,binsLIncl,xminLIncl,xmaxLIncl);
832 fOutputListStd->Add(fhnV0InRndLambda[i]);
9706bde3 833 fhnV0InMedLambda[i] = new THnSparseD(Form("fhnV0InMedLambda_%d",i),Form("Lambda: Mass vs Pt in med.-cl. cones, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});#it{p}_{T}^{V0} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimIncl,binsLIncl,xminLIncl,xmaxLIncl);
834 fOutputListStd->Add(fhnV0InMedLambda[i]);
7729cc34 835 fhnV0OutJetLambda[i] = new THnSparseD(Form("fhnV0OutJetLambda_%d",i),Form("Lambda: Pt outside jets, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});#it{p}_{T}^{V0} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimIncl,binsLIncl,xminLIncl,xmaxLIncl);
836 fOutputListStd->Add(fhnV0OutJetLambda[i]);
837 fhnV0NoJetLambda[i] = new THnSparseD(Form("fhnV0NoJetLambda_%d",i),Form("Lambda: Pt in jet-less events, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});#it{p}_{T}^{V0} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimIncl,binsLIncl,xminLIncl,xmaxLIncl);
838 fOutputListStd->Add(fhnV0NoJetLambda[i]);
839 fhnV0InJetALambda[i] = new THnSparseD(Form("fhnV0InJetALambda_%d",i),Form("ALambda: Mass vs Pt in jets, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});#it{p}_{T}^{V0} (GeV/#it{c});#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimInJC,binsLInJC,xminLInJC,xmaxLInJC);
840 fOutputListStd->Add(fhnV0InJetALambda[i]);
841 fhnV0InPerpALambda[i] = new THnSparseD(Form("fhnV0InPerpALambda_%d",i),Form("ALambda: Mass vs Pt in perp. cones, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});#it{p}_{T}^{V0} (GeV/#it{c});#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimInJC,binsLInJC,xminLInJC,xmaxLInJC);
842 fOutputListStd->Add(fhnV0InPerpALambda[i]);
843 fhnV0InRndALambda[i] = new THnSparseD(Form("fhnV0InRndALambda_%d",i),Form("ALambda: Mass vs Pt in rnd. cones, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});#it{p}_{T}^{V0} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimIncl,binsLIncl,xminLIncl,xmaxLIncl);
844 fOutputListStd->Add(fhnV0InRndALambda[i]);
9706bde3 845 fhnV0InMedALambda[i] = new THnSparseD(Form("fhnV0InMedALambda_%d",i),Form("ALambda: Mass vs Pt in med.-cl. cones, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});#it{p}_{T}^{V0} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimIncl,binsLIncl,xminLIncl,xmaxLIncl);
846 fOutputListStd->Add(fhnV0InMedALambda[i]);
7729cc34 847 fhnV0OutJetALambda[i] = new THnSparseD(Form("fhnV0OutJetALambda_%d",i),Form("ALambda: Pt outside jets, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});#it{p}_{T}^{V0} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimIncl,binsLIncl,xminLIncl,xmaxLIncl);
848 fOutputListStd->Add(fhnV0OutJetALambda[i]);
849 fhnV0NoJetALambda[i] = new THnSparseD(Form("fhnV0NoJetALambda_%d",i),Form("ALambda: Pt in jet-less events, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});#it{p}_{T}^{V0} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimIncl,binsLIncl,xminLIncl,xmaxLIncl);
850 fOutputListStd->Add(fhnV0NoJetALambda[i]);
851
87e23d3e 852 fh2V0PtJetAngleK0s[i] = new TH2D(Form("fh2V0PtJetAngleK0s_%d",i),Form("K0s: #it{p}_{T}^{jet} vs angle V0-jet, cent: %s;#it{p}_{T}^{jet};#it{#alpha}",GetCentBinLabel(i).Data()),iNJetPtBins,dJetPtMin,dJetPtMax,100,0,fdRadiusJet+0.1);
7729cc34 853 fOutputListStd->Add(fh2V0PtJetAngleK0s[i]);
87e23d3e 854 fh2V0PtJetAngleLambda[i] = new TH2D(Form("fh2V0PtJetAngleLambda_%d",i),Form("Lambda: #it{p}_{T}^{jet} vs angle V0-jet, cent: %s;#it{p}_{T}^{jet};#it{#alpha}",GetCentBinLabel(i).Data()),iNJetPtBins,dJetPtMin,dJetPtMax,100,0,fdRadiusJet+0.1);
7729cc34 855 fOutputListStd->Add(fh2V0PtJetAngleLambda[i]);
87e23d3e 856 fh2V0PtJetAngleALambda[i] = new TH2D(Form("fh2V0PtJetAngleALambda_%d",i),Form("ALambda: #it{p}_{T}^{jet} vs angle V0-jet, cent: %s;#it{p}_{T}^{jet};#it{#alpha}",GetCentBinLabel(i).Data()),iNJetPtBins,dJetPtMin,dJetPtMax,100,0,fdRadiusJet+0.1);
7729cc34 857 fOutputListStd->Add(fh2V0PtJetAngleALambda[i]);
858
859 fh1DCAInK0s[i] = new TH1D(Form("fh1DCAInK0s_%d",i),Form("K0s in jets: DCA daughters, cent %s;DCA (#sigma)",GetCentBinLabel(i).Data()),50,0,1);
860 fOutputListQA->Add(fh1DCAInK0s[i]);
861 fh1DCAInLambda[i] = new TH1D(Form("fh1DCAInLambda_%d",i),Form("Lambda in jets: DCA daughters, cent %s;DCA (#sigma)",GetCentBinLabel(i).Data()),50,0,1);
862 fOutputListQA->Add(fh1DCAInLambda[i]);
863 fh1DCAInALambda[i] = new TH1D(Form("fh1DCAInALambda_%d",i),Form("ALambda in jets: DCA daughters, cent %s;DCA (#sigma)",GetCentBinLabel(i).Data()),50,0,1);
864 fOutputListQA->Add(fh1DCAInALambda[i]);
865
866 fh1DCAOutK0s[i] = new TH1D(Form("fh1DCAOutK0s_%d",i),Form("K0s outside jets: DCA daughters, cent %s;DCA (#sigma)",GetCentBinLabel(i).Data()),50,0,1);
867 fOutputListQA->Add(fh1DCAOutK0s[i]);
868 fh1DCAOutLambda[i] = new TH1D(Form("fh1DCAOutLambda_%d",i),Form("Lambda outside jets: DCA daughters, cent %s;DCA (#sigma)",GetCentBinLabel(i).Data()),50,0,1);
869 fOutputListQA->Add(fh1DCAOutLambda[i]);
870 fh1DCAOutALambda[i] = new TH1D(Form("fh1DCAOutALambda_%d",i),Form("ALambda outside jets: DCA daughters, cent %s;DCA (#sigma)",GetCentBinLabel(i).Data()),50,0,1);
871 fOutputListQA->Add(fh1DCAOutALambda[i]);
872
873 fh1DeltaZK0s[i] = new TH1D(Form("fh1DeltaZK0s_%d",i),Form("K0s: #Delta#it{z} vertices, cent %s;#it{z} (cm)",GetCentBinLabel(i).Data()),50,-10,10);
874 fOutputListQA->Add(fh1DeltaZK0s[i]);
875 fh1DeltaZLambda[i] = new TH1D(Form("fh1DeltaZLambda_%d",i),Form("Lambda: #Delta#it{z} vertices, cent %s;#it{z} (cm)",GetCentBinLabel(i).Data()),50,-10,10);
876 fOutputListQA->Add(fh1DeltaZLambda[i]);
877 fh1DeltaZALambda[i] = new TH1D(Form("fh1DeltaZALambda_%d",i),Form("ALambda: #Delta#it{z} vertices, cent %s;#it{z} (cm)",GetCentBinLabel(i).Data()),50,-10,10);
878 fOutputListQA->Add(fh1DeltaZALambda[i]);
879
880 // jet histograms
87e23d3e 881 fh1PtJet[i] = new TH1D(Form("fh1PtJet_%d",i),Form("Jet pt spectrum, cent: %s;#it{p}_{T} jet (GeV/#it{c})",GetCentBinLabel(i).Data()),iNJetPtBins,dJetPtMin,dJetPtMax);
7729cc34 882 fOutputListStd->Add(fh1PtJet[i]);
883 fh1EtaJet[i] = new TH1D(Form("fh1EtaJet_%d",i),Form("Jet eta spectrum, cent: %s;#it{#eta} jet",GetCentBinLabel(i).Data()),80,-1.,1.);
884 fOutputListStd->Add(fh1EtaJet[i]);
87e23d3e 885 fh2EtaPtJet[i] = new TH2D(Form("fh2EtaPtJet_%d",i),Form("Jet eta vs pT spectrum, cent: %s;#it{#eta} jet;#it{p}_{T} jet (GeV/#it{c})",GetCentBinLabel(i).Data()),80,-1.,1.,iNJetPtBins,dJetPtMin,dJetPtMax);
7729cc34 886 fOutputListStd->Add(fh2EtaPtJet[i]);
887 fh2EtaPhiRndCone[i] = new TH2D(Form("fh2EtaPhiRndCone_%d",i),Form("Rnd. cones: eta vs phi, cent: %s;#it{#eta} cone;#it{#phi} cone",GetCentBinLabel(i).Data()),80,-1.,1.,100,0.,TMath::TwoPi());
888 fOutputListStd->Add(fh2EtaPhiRndCone[i]);
9706bde3 889 fh2EtaPhiMedCone[i] = new TH2D(Form("fh2EtaPhiMedCone_%d",i),Form("Med.-cl. cones: eta vs phi, cent: %s;#it{#eta} cone;#it{#phi} cone",GetCentBinLabel(i).Data()),80,-1.,1.,100,0.,TMath::TwoPi());
890 fOutputListStd->Add(fh2EtaPhiMedCone[i]);
7729cc34 891 fh1PhiJet[i] = new TH1D(Form("fh1PhiJet_%d",i),Form("Jet phi spectrum, cent: %s;#it{#phi} jet",GetCentBinLabel(i).Data()),100,0.,TMath::TwoPi());
892 fOutputListStd->Add(fh1PhiJet[i]);
893 fh1NJetPerEvent[i] = new TH1D(Form("fh1NJetPerEvent_%d",i),Form("Number of selected jets per event, cent: %s;# jets;# events",GetCentBinLabel(i).Data()),100,0.,100.);
894 fOutputListStd->Add(fh1NJetPerEvent[i]);
895 // event histograms
87e23d3e 896 fh1VtxZ[i] = new TH1D(Form("fh1VtxZ_%d",i),Form("#it{z} coordinate of the primary vertex, cent: %s;#it{z} (cm)",GetCentBinLabel(i).Data()),150,-1.5*fdCutVertexZ,1.5*fdCutVertexZ);
7729cc34 897 fOutputListQA->Add(fh1VtxZ[i]);
898 fh2VtxXY[i] = new TH2D(Form("fh2VtxXY_%d",i),Form("#it{xy} coordinate of the primary vertex, cent: %s;#it{x} (cm);#it{y} (cm)",GetCentBinLabel(i).Data()),200,-0.2,0.2,500,-0.5,0.5);
899 fOutputListQA->Add(fh2VtxXY[i]);
900 // fOutputListStd->Add([i]);
901 if (fbMCAnalysis)
902 {
903 // inclusive pt
87e23d3e 904 fh1V0K0sPtMCGen[i] = new TH1D(Form("fh1V0K0sPtMCGen_%d",i),Form("MC K0s generated: pt spectrum, cent: %s;MC #it{p}_{T} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNBinsPtV0,dPtV0Min,dPtV0Max);
7729cc34 905 fOutputListMC->Add(fh1V0K0sPtMCGen[i]);
87e23d3e 906 fh2V0K0sPtMassMCRec[i] = new TH2D(Form("fh2V0K0sPtMassMCRec_%d",i),Form("MC K0s associated: pt-m spectrum, cent: %s;MC #it{p}_{T} (GeV/#it{c});#it{m}_{inv} (GeV/#it{c}^{2})",GetCentBinLabel(i).Data()),iNBinsPtV0,dPtV0Min,dPtV0Max,fgkiNBinsMassK0s,fgkdMassK0sMin,fgkdMassK0sMax);
7729cc34 907 fOutputListMC->Add(fh2V0K0sPtMassMCRec[i]);
87e23d3e 908 fh1V0K0sPtMCRecFalse[i] = new TH1D(Form("fh1V0K0sPtMCRecFalse_%d",i),Form("MC K0s false: pt spectrum, cent: %s;MC #it{p}_{T} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNBinsPtV0,dPtV0Min,dPtV0Max);
7729cc34 909 fOutputListMC->Add(fh1V0K0sPtMCRecFalse[i]);
910 // inclusive pt-eta
87e23d3e 911 fh2V0K0sEtaPtMCGen[i] = new TH2D(Form("fh2V0K0sEtaPtMCGen_%d",i),Form("MC K0s generated: pt-eta spectrum, cent: %s;MC #it{p}_{T} (GeV/#it{c});#eta",GetCentBinLabel(i).Data()),iNBinsPtV0,dPtV0Min,dPtV0Max,iNBinsEtaV0,-dRangeEtaV0Max,dRangeEtaV0Max);
7729cc34 912 fOutputListMC->Add(fh2V0K0sEtaPtMCGen[i]);
913 fh3V0K0sEtaPtMassMCRec[i] = new THnSparseD(Form("fh3V0K0sEtaPtMassMCRec_%d",i),Form("MC K0s associated: m-pt-eta spectrum, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});MC #it{p}_{T} (GeV/#it{c});#eta",GetCentBinLabel(i).Data()),3,binsEtaK,xminEtaK,xmaxEtaK);
914 fOutputListMC->Add(fh3V0K0sEtaPtMassMCRec[i]);
915 // in jet pt
87e23d3e 916 fh2V0K0sInJetPtMCGen[i] = new TH2D(Form("fh2V0K0sInJetPtMCGen_%d",i),Form("MC K0s in jet generated: pt-ptJet spectrum, cent: %s;MC #it{p}_{T} (GeV/#it{c});#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNBinsPtV0,dPtV0Min,dPtV0Max,iNJetPtBins,dJetPtMin,dJetPtMax);
7729cc34 917 fOutputListMC->Add(fh2V0K0sInJetPtMCGen[i]);
918 fh3V0K0sInJetPtMassMCRec[i] = new THnSparseD(Form("fh3V0K0sInJetPtMassMCRec_%d",i),Form("MC K0s in jet associated: m-pt-ptJet spectrum, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});MC #it{p}_{T} (GeV/#it{c});#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimInJC,binsKInJC,xminKInJC,xmaxKInJC);
919 fOutputListMC->Add(fh3V0K0sInJetPtMassMCRec[i]);
920 // in jet pt-eta
921 fh3V0K0sInJetEtaPtMCGen[i] = new THnSparseD(Form("fh3V0K0sInJetEtaPtMCGen_%d",i),Form("MC K0s generated: pt-eta-ptJet spectrum, cent: %s;MC #it{p}_{T} (GeV/#it{c});#eta;#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),4,binsEtaInGen,xminEtaInGen,xmaxEtaInGen);
922 fOutputListMC->Add(fh3V0K0sInJetEtaPtMCGen[i]);
923 fh4V0K0sInJetEtaPtMassMCRec[i] = new THnSparseD(Form("fh4V0K0sInJetEtaPtMassMCRec_%d",i),Form("MC K0s associated: m-pt-eta-ptJet spectrum, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});MC #it{p}_{T} (GeV/#it{c});#eta;#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),5,binsEtaKInRec,xminEtaKInRec,xmaxEtaKInRec);
924 fOutputListMC->Add(fh4V0K0sInJetEtaPtMassMCRec[i]);
925
87e23d3e 926 fh2V0K0sMCResolMPt[i] = new TH2D(Form("fh2V0K0sMCResolMPt_%d",i),Form("MC K0s associated: #Delta#it{m} vs pt, cent %s;#Delta#it{m} = #it{m}_{inv} - #it{m}_{true} (GeV/#it{c}^{2});#it{p}_{T}^{rec} (GeV/#it{c})",GetCentBinLabel(i).Data()),100,-0.02,0.02,iNBinsPtV0,dPtV0Min,dPtV0Max);
7729cc34 927 fOutputListMC->Add(fh2V0K0sMCResolMPt[i]);
87e23d3e 928 fh2V0K0sMCPtGenPtRec[i] = new TH2D(Form("fh2V0K0sMCPtGenPtRec_%d",i),Form("MC K0s associated: pt gen vs pt rec, cent %s;#it{p}_{T}^{gen} (GeV/#it{c});#it{p}_{T}^{rec} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNBinsPtV0,dPtV0Min,dPtV0Max,iNBinsPtV0,dPtV0Min,dPtV0Max);
7729cc34 929 fOutputListMC->Add(fh2V0K0sMCPtGenPtRec[i]);
930
931 // inclusive pt
87e23d3e 932 fh1V0LambdaPtMCGen[i] = new TH1D(Form("fh1V0LambdaPtMCGen_%d",i),Form("MC Lambda generated: pt spectrum, cent: %s;MC #it{p}_{T} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNBinsPtV0,dPtV0Min,dPtV0Max);
7729cc34 933 fOutputListMC->Add(fh1V0LambdaPtMCGen[i]);
87e23d3e 934 fh2V0LambdaPtMassMCRec[i] = new TH2D(Form("fh2V0LambdaPtMassMCRec_%d",i),Form("MC Lambda associated: pt-m spectrum, cent: %s;MC #it{p}_{T} (GeV/#it{c});#it{m}_{inv} (GeV/#it{c}^{2})",GetCentBinLabel(i).Data()),iNBinsPtV0,dPtV0Min,dPtV0Max,fgkiNBinsMassLambda,fgkdMassLambdaMin,fgkdMassLambdaMax);
7729cc34 935 fOutputListMC->Add(fh2V0LambdaPtMassMCRec[i]);
87e23d3e 936 fh1V0LambdaPtMCRecFalse[i] = new TH1D(Form("fh1V0LambdaPtMCRecFalse_%d",i),Form("MC Lambda false: pt spectrum, cent: %s;MC #it{p}_{T} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNBinsPtV0,dPtV0Min,dPtV0Max);
7729cc34 937 fOutputListMC->Add(fh1V0LambdaPtMCRecFalse[i]);
938 // inclusive pt-eta
87e23d3e 939 fh2V0LambdaEtaPtMCGen[i] = new TH2D(Form("fh2V0LambdaEtaPtMCGen_%d",i),Form("MC Lambda generated: pt-eta spectrum, cent: %s;MC #it{p}_{T} (GeV/#it{c});#eta",GetCentBinLabel(i).Data()),iNBinsPtV0,dPtV0Min,dPtV0Max,iNBinsEtaV0,-dRangeEtaV0Max,dRangeEtaV0Max);
7729cc34 940 fOutputListMC->Add(fh2V0LambdaEtaPtMCGen[i]);
941 fh3V0LambdaEtaPtMassMCRec[i] = new THnSparseD(Form("fh3V0LambdaEtaPtMassMCRec_%d",i),Form("MC Lambda associated: m-pt-eta spectrum, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});MC #it{p}_{T} (GeV/#it{c});#eta",GetCentBinLabel(i).Data()),3,binsEtaL,xminEtaL,xmaxEtaL);
942 fOutputListMC->Add(fh3V0LambdaEtaPtMassMCRec[i]);
943 // in jet pt
87e23d3e 944 fh2V0LambdaInJetPtMCGen[i] = new TH2D(Form("fh2V0LambdaInJetPtMCGen_%d",i),Form("MC Lambda in jet generated: pt-ptJet spectrum, cent: %s;MC #it{p}_{T} (GeV/#it{c});#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNBinsPtV0,dPtV0Min,dPtV0Max,iNJetPtBins,dJetPtMin,dJetPtMax);
7729cc34 945 fOutputListMC->Add(fh2V0LambdaInJetPtMCGen[i]);
946 fh3V0LambdaInJetPtMassMCRec[i] = new THnSparseD(Form("fh3V0LambdaInJetPtMassMCRec_%d",i),Form("MC Lambda in jet associated: m-pt-ptJet spectrum, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});MC #it{p}_{T} (GeV/#it{c});#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimInJC,binsLInJC,xminLInJC,xmaxLInJC);
947 fOutputListMC->Add(fh3V0LambdaInJetPtMassMCRec[i]);
948 // in jet pt-eta
949 fh3V0LambdaInJetEtaPtMCGen[i] = new THnSparseD(Form("fh3V0LambdaInJetEtaPtMCGen_%d",i),Form("MC Lambda generated: pt-eta-ptJet spectrum, cent: %s;MC #it{p}_{T} (GeV/#it{c});#eta;#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),4,binsEtaInGen,xminEtaInGen,xmaxEtaInGen);
950 fOutputListMC->Add(fh3V0LambdaInJetEtaPtMCGen[i]);
951 fh4V0LambdaInJetEtaPtMassMCRec[i] = new THnSparseD(Form("fh4V0LambdaInJetEtaPtMassMCRec_%d",i),Form("MC Lambda associated: m-pt-eta-ptJet spectrum, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});MC #it{p}_{T} (GeV/#it{c});#eta;#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),5,binsEtaLInRec,xminEtaLInRec,xmaxEtaLInRec);
952 fOutputListMC->Add(fh4V0LambdaInJetEtaPtMassMCRec[i]);
953
87e23d3e 954 fh2V0LambdaMCResolMPt[i] = new TH2D(Form("fh2V0LambdaMCResolMPt_%d",i),Form("MC Lambda associated: #Delta#it{m} vs pt, cent %s;#Delta#it{m} = #it{m}_{inv} - #it{m}_{true} (GeV/#it{c}^{2});#it{p}_{T}^{rec} (GeV/#it{c})",GetCentBinLabel(i).Data()),100,-0.02,0.02,iNBinsPtV0,dPtV0Min,dPtV0Max);
7729cc34 955 fOutputListMC->Add(fh2V0LambdaMCResolMPt[i]);
87e23d3e 956 fh2V0LambdaMCPtGenPtRec[i] = new TH2D(Form("fh2V0LambdaMCPtGenPtRec_%d",i),Form("MC Lambda associated: pt gen vs pt rec, cent %s;#it{p}_{T}^{gen} (GeV/#it{c});#it{p}_{T}^{rec} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNBinsPtV0,dPtV0Min,dPtV0Max,iNBinsPtV0,dPtV0Min,dPtV0Max);
7729cc34 957 fOutputListMC->Add(fh2V0LambdaMCPtGenPtRec[i]);
958
959 // inclusive pt
87e23d3e 960 fh1V0ALambdaPtMCGen[i] = new TH1D(Form("fh1V0ALambdaPtMCGen_%d",i),Form("MC ALambda generated: pt spectrum, cent: %s;MC #it{p}_{T} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNBinsPtV0,dPtV0Min,dPtV0Max);
7729cc34 961 fOutputListMC->Add(fh1V0ALambdaPtMCGen[i]);
87e23d3e 962 fh2V0ALambdaPtMassMCRec[i] = new TH2D(Form("fh2V0ALambdaPtMassMCRec_%d",i),Form("MC ALambda associated: pt-m spectrum, cent: %s;MC #it{p}_{T} (GeV/#it{c});#it{m}_{inv} (GeV/#it{c}^{2})",GetCentBinLabel(i).Data()),iNBinsPtV0,dPtV0Min,dPtV0Max,fgkiNBinsMassLambda,fgkdMassLambdaMin,fgkdMassLambdaMax);
7729cc34 963 fOutputListMC->Add(fh2V0ALambdaPtMassMCRec[i]);
87e23d3e 964 fh1V0ALambdaPtMCRecFalse[i] = new TH1D(Form("fh1V0ALambdaPtMCRecFalse_%d",i),Form("MC ALambda false: pt spectrum, cent: %s;MC #it{p}_{T} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNBinsPtV0,dPtV0Min,dPtV0Max);
7729cc34 965 fOutputListMC->Add(fh1V0ALambdaPtMCRecFalse[i]);
966 // inclusive pt-eta
87e23d3e 967 fh2V0ALambdaEtaPtMCGen[i] = new TH2D(Form("fh2V0ALambdaEtaPtMCGen_%d",i),Form("MC ALambda generated: pt-eta spectrum, cent: %s;MC #it{p}_{T} (GeV/#it{c});#eta",GetCentBinLabel(i).Data()),iNBinsPtV0,dPtV0Min,dPtV0Max,iNBinsEtaV0,-dRangeEtaV0Max,dRangeEtaV0Max);
7729cc34 968 fOutputListMC->Add(fh2V0ALambdaEtaPtMCGen[i]);
969 fh3V0ALambdaEtaPtMassMCRec[i] = new THnSparseD(Form("fh3V0ALambdaEtaPtMassMCRec_%d",i),Form("MC ALambda associated: m-pt-eta spectrum, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});MC #it{p}_{T} (GeV/#it{c});#eta",GetCentBinLabel(i).Data()),3,binsEtaL,xminEtaL,xmaxEtaL);
970 fOutputListMC->Add(fh3V0ALambdaEtaPtMassMCRec[i]);
971 // in jet pt
87e23d3e 972 fh2V0ALambdaInJetPtMCGen[i] = new TH2D(Form("fh2V0ALambdaInJetPtMCGen_%d",i),Form("MC ALambda in jet generated: pt-ptJet spectrum, cent: %s;MC #it{p}_{T} (GeV/#it{c});#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNBinsPtV0,dPtV0Min,dPtV0Max,iNJetPtBins,dJetPtMin,dJetPtMax);
7729cc34 973 fOutputListMC->Add(fh2V0ALambdaInJetPtMCGen[i]);
974 fh3V0ALambdaInJetPtMassMCRec[i] = new THnSparseD(Form("fh3V0ALambdaInJetPtMassMCRec_%d",i),Form("MC ALambda in jet associated: m-pt-ptJet spectrum, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});MC #it{p}_{T} (GeV/#it{c});#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimInJC,binsLInJC,xminLInJC,xmaxLInJC);
975 fOutputListMC->Add(fh3V0ALambdaInJetPtMassMCRec[i]);
976 // in jet pt-eta
977 fh3V0ALambdaInJetEtaPtMCGen[i] = new THnSparseD(Form("fh3V0ALambdaInJetEtaPtMCGen_%d",i),Form("MC ALambda generated: pt-eta-ptJet spectrum, cent: %s;MC #it{p}_{T} (GeV/#it{c});#eta;#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),4,binsEtaInGen,xminEtaInGen,xmaxEtaInGen);
978 fOutputListMC->Add(fh3V0ALambdaInJetEtaPtMCGen[i]);
979 fh4V0ALambdaInJetEtaPtMassMCRec[i] = new THnSparseD(Form("fh4V0ALambdaInJetEtaPtMassMCRec_%d",i),Form("MC ALambda associated: m-pt-eta-ptJet spectrum, cent: %s;#it{m}_{inv} (GeV/#it{c}^{2});MC #it{p}_{T} (GeV/#it{c});#eta;#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),5,binsEtaLInRec,xminEtaLInRec,xmaxEtaLInRec);
980 fOutputListMC->Add(fh4V0ALambdaInJetEtaPtMassMCRec[i]);
981
87e23d3e 982 fh2V0ALambdaMCResolMPt[i] = new TH2D(Form("fh2V0ALambdaMCResolMPt_%d",i),Form("MC ALambda associated: #Delta#it{m} vs pt, cent %s;#Delta#it{m} = #it{m}_{inv} - #it{m}_{true} (GeV/#it{c}^{2});#it{p}_{T}^{rec} (GeV/#it{c})",GetCentBinLabel(i).Data()),100,-0.02,0.02,iNBinsPtV0,dPtV0Min,dPtV0Max);
7729cc34 983 fOutputListMC->Add(fh2V0ALambdaMCResolMPt[i]);
87e23d3e 984 fh2V0ALambdaMCPtGenPtRec[i] = new TH2D(Form("fh2V0ALambdaMCPtGenPtRec_%d",i),Form("MC ALambda associated: pt gen vs pt rec, cent %s;#it{p}_{T}^{gen} (GeV/#it{c});#it{p}_{T}^{rec} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNBinsPtV0,dPtV0Min,dPtV0Max,iNBinsPtV0,dPtV0Min,dPtV0Max);
7729cc34 985 fOutputListMC->Add(fh2V0ALambdaMCPtGenPtRec[i]);
986
987 Int_t iNBinsPtXi = 80;
988 Double_t dPtXiMin = 0;
989 Double_t dPtXiMax = 8;
990 const Int_t iNDimFD = 3;
991 Int_t binsFD[iNDimFD] = {iNBinsPtV0, iNBinsPtXi, iNJetPtBins};
87e23d3e 992 Double_t xminFD[iNDimFD] = {dPtV0Min, dPtXiMin, dJetPtMin};
993 Double_t xmaxFD[iNDimFD] = {dPtV0Max, dPtXiMax, dJetPtMax};
7729cc34 994 fhnV0LambdaInclMCFD[i] = new THnSparseD(Form("fhnV0LambdaInclMCFD_%d",i),Form("MC Lambda associated, inclusive, from Xi: pt-pt, cent %s;#it{p}_{T}^{#Lambda,gen.} (GeV/#it{c});#it{p}_{T}^{#Xi,gen.} (GeV/#it{c});#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimFD,binsFD,xminFD,xmaxFD);
995 fOutputListMC->Add(fhnV0LambdaInclMCFD[i]);
996 fhnV0LambdaInJetsMCFD[i] = new THnSparseD(Form("fhnV0LambdaInJetsMCFD_%d",i),Form("MC Lambda associated, in JC, from Xi: pt-pt-ptJet, cent %s;#it{p}_{T}^{#Lambda,gen.} (GeV/#it{c});#it{p}_{T}^{#Xi,gen.} (GeV/#it{c});#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimFD,binsFD,xminFD,xmaxFD);
997 fOutputListMC->Add(fhnV0LambdaInJetsMCFD[i]);
998 fhnV0LambdaBulkMCFD[i] = new THnSparseD(Form("fhnV0LambdaBulkMCFD_%d",i),Form("MC Lambda associated, in no jet events, from Xi: pt-pt, cent %s;#it{p}_{T}^{#Lambda,gen.} (GeV/#it{c});#it{p}_{T}^{#Xi,gen.} (GeV/#it{c});#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimFD,binsFD,xminFD,xmaxFD);
999 fOutputListMC->Add(fhnV0LambdaBulkMCFD[i]);
1000 fh1V0XiPtMCGen[i] = new TH1D(Form("fh1V0XiPtMCGen_%d",i),Form("MC Xi^{-} generated: Pt spectrum, cent %s;#it{p}_{T}^{#Xi^{-},gen.} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNBinsPtXi,dPtXiMin,dPtXiMax);
1001 fOutputListMC->Add(fh1V0XiPtMCGen[i]);
1002 fhnV0ALambdaInclMCFD[i] = new THnSparseD(Form("fhnV0ALambdaInclMCFD_%d",i),Form("MC ALambda associated, from AXi: pt-pt, cent %s;#it{p}_{T}^{A#Lambda,gen.} (GeV/#it{c});#it{p}_{T}^{A#Xi,gen.} (GeV/#it{c});#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimFD,binsFD,xminFD,xmaxFD);
1003 fOutputListMC->Add(fhnV0ALambdaInclMCFD[i]);
1004 fhnV0ALambdaInJetsMCFD[i] = new THnSparseD(Form("fhnV0ALambdaInJetsMCFD_%d",i),Form("MC ALambda associated, in JC, from AXi: pt-pt-ptJet, cent %s;#it{p}_{T}^{A#Lambda,gen.} (GeV/#it{c});#it{p}_{T}^{A#Xi,gen.} (GeV/#it{c});#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimFD,binsFD,xminFD,xmaxFD);
1005 fOutputListMC->Add(fhnV0ALambdaInJetsMCFD[i]);
1006 fhnV0ALambdaBulkMCFD[i] = new THnSparseD(Form("fhnV0ALambdaBulkMCFD_%d",i),Form("MC ALambda associated, in no jet events, from AXi: pt-pt-ptJet, cent %s;#it{p}_{T}^{A#Lambda,gen.} (GeV/#it{c});#it{p}_{T}^{A#Xi,gen.} (GeV/#it{c});#it{p}_{T}^{jet} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNDimFD,binsFD,xminFD,xmaxFD);
1007 fOutputListMC->Add(fhnV0ALambdaBulkMCFD[i]);
1008 fh1V0AXiPtMCGen[i] = new TH1D(Form("fh1V0AXiPtMCGen_%d",i),Form("MC AXi^{-} generated: Pt spectrum, cent %s;#it{p}_{T}^{A#Xi^{-},gen.} (GeV/#it{c})",GetCentBinLabel(i).Data()),iNBinsPtXi,dPtXiMin,dPtXiMax);
1009 fOutputListMC->Add(fh1V0AXiPtMCGen[i]);
1010
1011 // daughter eta
1012// fhnV0K0sInclDaughterEtaPtPtMCGen[i] = new THnSparseD(Form("fhnV0K0sInclDaughterEtaPtPtMCGen_%d",i),Form("MC K0S, inclusive, gen., daughters: charge-etaD-ptD-etaV0-ptV0-ptJet, cent: %s;charge;eta daughter;pT daughter;eta V0;pT V0;pT jet",GetCentBinLabel(i).Data()),iNDimEtaD,binsEtaDaughter,xminEtaDaughter,xmaxEtaDaughter);
1013 fhnV0K0sInclDaughterEtaPtPtMCRec[i] = new THnSparseD(Form("fhnV0K0sInclDaughterEtaPtPtMCRec_%d",i),Form("MC K0S, inclusive, assoc., daughters: charge-etaD-ptD-etaV0-ptV0-ptJet, cent: %s;charge;eta daughter;pT daughter;eta V0;pT V0;pT jet",GetCentBinLabel(i).Data()),iNDimEtaD,binsEtaDaughter,xminEtaDaughter,xmaxEtaDaughter);
1014// fhnV0K0sInJetsDaughterEtaPtPtMCGen[i] = new THnSparseD(Form("fhnV0K0sInJetsDaughterEtaPtPtMCGen_%d",i),Form("MC K0S, in JC, gen., daughters: charge-etaD-ptD-etaV0-ptV0-ptJet, cent: %s;charge;eta daughter;pT daughter;eta V0;pT V0;pT jet",GetCentBinLabel(i).Data()),iNDimEtaD,binsEtaDaughter,xminEtaDaughter,xmaxEtaDaughter);
1015 fhnV0K0sInJetsDaughterEtaPtPtMCRec[i] = new THnSparseD(Form("fhnV0K0sInJetsDaughterEtaPtPtMCRec_%d",i),Form("MC K0S, in JC, assoc., daughters: charge-etaD-ptD-etaV0-ptV0-ptJet, cent: %s;charge;eta daughter;pT daughter;eta V0;pT V0;pT jet",GetCentBinLabel(i).Data()),iNDimEtaD,binsEtaDaughter,xminEtaDaughter,xmaxEtaDaughter);
1016// fhnV0LambdaInclDaughterEtaPtPtMCGen[i] = new THnSparseD(Form("fhnV0LambdaInclDaughterEtaPtPtMCGen_%d",i),Form("MC Lambda, inclusive, gen., daughters: charge-etaD-ptD-etaV0-ptV0-ptJet, cent: %s;charge;eta daughter;pT daughter;eta V0;pT V0;pT jet",GetCentBinLabel(i).Data()),iNDimEtaD,binsEtaDaughter,xminEtaDaughter,xmaxEtaDaughter);
1017 fhnV0LambdaInclDaughterEtaPtPtMCRec[i] = new THnSparseD(Form("fhnV0LambdaInclDaughterEtaPtPtMCRec_%d",i),Form("MC Lambda, inclusive, assoc., daughters: charge-etaD-ptD-etaV0-ptV0-ptJet, cent: %s;charge;eta daughter;pT daughter;eta V0;pT V0;pT jet",GetCentBinLabel(i).Data()),iNDimEtaD,binsEtaDaughter,xminEtaDaughter,xmaxEtaDaughter);
1018// fhnV0LambdaInJetsDaughterEtaPtPtMCGen[i] = new THnSparseD(Form("fhnV0LambdaInJetsDaughterEtaPtPtMCGen_%d",i),Form("MC Lambda, in JC, gen., daughters: charge-etaD-ptD-etaV0-ptV0-ptJet, cent: %s;charge;eta daughter;pT daughter;eta V0;pT V0;pT jet",GetCentBinLabel(i).Data()),iNDimEtaD,binsEtaDaughter,xminEtaDaughter,xmaxEtaDaughter);
1019 fhnV0LambdaInJetsDaughterEtaPtPtMCRec[i] = new THnSparseD(Form("fhnV0LambdaInJetsDaughterEtaPtPtMCRec_%d",i),Form("MC Lambda, in JC, assoc., daughters: charge-etaD-ptD-etaV0-ptV0-ptJet, cent: %s;charge;eta daughter;pT daughter;eta V0;pT V0;pT jet",GetCentBinLabel(i).Data()),iNDimEtaD,binsEtaDaughter,xminEtaDaughter,xmaxEtaDaughter);
1020// fhnV0ALambdaInclDaughterEtaPtPtMCGen[i] = new THnSparseD(Form("fhnV0ALambdaInclDaughterEtaPtPtMCGen_%d",i),Form("MC ALambda, inclusive, gen., daughters: charge-etaD-ptD-etaV0-ptV0-ptJet, cent: %s;charge;eta daughter;pT daughter;eta V0;pT V0;pT jet",GetCentBinLabel(i).Data()),iNDimEtaD,binsEtaDaughter,xminEtaDaughter,xmaxEtaDaughter);
1021 fhnV0ALambdaInclDaughterEtaPtPtMCRec[i] = new THnSparseD(Form("fhnV0ALambdaInclDaughterEtaPtPtMCRec_%d",i),Form("MC ALambda, inclusive, assoc., daughters: charge-etaD-ptD-etaV0-ptV0-ptJet, cent: %s;charge;eta daughter;pT daughter;eta V0;pT V0;pT jet",GetCentBinLabel(i).Data()),iNDimEtaD,binsEtaDaughter,xminEtaDaughter,xmaxEtaDaughter);
1022// fhnV0ALambdaInJetsDaughterEtaPtPtMCGen[i] = new THnSparseD(Form("fhnV0ALambdaInJetsDaughterEtaPtPtMCGen_%d",i),Form("MC ALambda, in JC, gen., daughters: charge-etaD-ptD-etaV0-ptV0-ptJet, cent: %s;charge;eta daughter;pT daughter;eta V0;pT V0;pT jet",GetCentBinLabel(i).Data()),iNDimEtaD,binsEtaDaughter,xminEtaDaughter,xmaxEtaDaughter);
1023 fhnV0ALambdaInJetsDaughterEtaPtPtMCRec[i] = new THnSparseD(Form("fhnV0ALambdaInJetsDaughterEtaPtPtMCRec_%d",i),Form("MC ALambda, in JC, assoc., daughters: charge-etaD-ptD-etaV0-ptV0-ptJet, cent: %s;charge;eta daughter;pT daughter;eta V0;pT V0;pT jet",GetCentBinLabel(i).Data()),iNDimEtaD,binsEtaDaughter,xminEtaDaughter,xmaxEtaDaughter);
1024
1025// fOutputListMC->Add(fhnV0K0sInclDaughterEtaPtPtMCGen[i]);
1026 fOutputListMC->Add(fhnV0K0sInclDaughterEtaPtPtMCRec[i]);
1027// fOutputListMC->Add(fhnV0K0sInJetsDaughterEtaPtPtMCGen[i]);
1028 fOutputListMC->Add(fhnV0K0sInJetsDaughterEtaPtPtMCRec[i]);
1029// fOutputListMC->Add(fhnV0LambdaInclDaughterEtaPtPtMCGen[i]);
1030 fOutputListMC->Add(fhnV0LambdaInclDaughterEtaPtPtMCRec[i]);
1031// fOutputListMC->Add(fhnV0LambdaInJetsDaughterEtaPtPtMCGen[i]);
1032 fOutputListMC->Add(fhnV0LambdaInJetsDaughterEtaPtPtMCRec[i]);
1033// fOutputListMC->Add(fhnV0ALambdaInclDaughterEtaPtPtMCGen[i]);
1034 fOutputListMC->Add(fhnV0ALambdaInclDaughterEtaPtPtMCRec[i]);
1035// fOutputListMC->Add(fhnV0ALambdaInJetsDaughterEtaPtPtMCGen[i]);
1036 fOutputListMC->Add(fhnV0ALambdaInJetsDaughterEtaPtPtMCRec[i]);
1037 }
1038 }
1039
1040 // QA Histograms
1041 for (Int_t i = 0; i < fgkiNQAIndeces; i++)
1042 {
1043// [i] = new TH1D(Form("%d",i),";;Counts",,,);
1044 fh1QAV0Status[i] = new TH1D(Form("fh1QAV0Status_%d",i),"QA: V0 status",2,0,2);
1045 fh1QAV0TPCRefit[i] = new TH1D(Form("fh1QAV0TPCRefit_%d",i),"QA: TPC refit",2,0,2);
1046 fh1QAV0TPCRows[i] = new TH1D(Form("fh1QAV0TPCRows_%d",i),"QA: TPC Rows",160,0,160);
1047 fh1QAV0TPCFindable[i] = new TH1D(Form("fh1QAV0TPCFindable_%d",i),"QA: TPC Findable",160,0,160);
1048 fh1QAV0TPCRowsFind[i] = new TH1D(Form("fh1QAV0TPCRowsFind_%d",i),"QA: TPC Rows/Findable",100,0,2);
1049 fh1QAV0Eta[i] = new TH1D(Form("fh1QAV0Eta_%d",i),"QA: Daughter Eta",200,-2,2);
1050 fh2QAV0EtaRows[i] = new TH2D(Form("fh2QAV0EtaRows_%d",i),"QA: Daughter Eta vs TPC rows;#eta;TPC rows",200,-2,2,160,0,160);
1051 fh2QAV0PtRows[i] = new TH2D(Form("fh2QAV0PtRows_%d",i),"QA: Daughter Pt vs TPC rows;pt;TPC rows",100,0,10,160,0,160);
1052 fh2QAV0PhiRows[i] = new TH2D(Form("fh2QAV0PhiRows_%d",i),"QA: Daughter Phi vs TPC rows;#phi;TPC rows",100,0,TMath::TwoPi(),160,0,160);
1053 fh2QAV0NClRows[i] = new TH2D(Form("fh2QAV0NClRows_%d",i),"QA: Daughter NCl vs TPC rows;findable clusters;TPC rows",100,0,160,160,0,160);
1054 fh2QAV0EtaNCl[i] = new TH2D(Form("fh2QAV0EtaNCl_%d",i),"QA: Daughter Eta vs NCl;#eta;findable clusters",200,-2,2,160,0,160);
1055
87e23d3e 1056 fh2QAV0EtaPtK0sPeak[i] = new TH2D(Form("fh2QAV0EtaPtK0sPeak_%d",i),"QA: K0s: Daughter Eta vs V0 pt, peak;track eta;V0 pt",200,-2,2,iNBinsPtV0,dPtV0Min,dPtV0Max);
7729cc34 1057 fh2QAV0EtaEtaK0s[i] = new TH2D(Form("fh2QAV0EtaEtaK0s_%d",i),"QA: K0s: Eta vs Eta Daughter",200,-2,2,200,-2,2);
1058 fh2QAV0PhiPhiK0s[i] = new TH2D(Form("fh2QAV0PhiPhiK0s_%d",i),"QA: K0s: Phi vs Phi Daughter",200,0,TMath::TwoPi(),200,0,TMath::TwoPi());
1059 fh1QAV0RapK0s[i] = new TH1D(Form("fh1QAV0RapK0s_%d",i),"QA: K0s: V0 Rapidity",200,-2,2);
1060 fh2QAV0PtPtK0sPeak[i] = new TH2D(Form("fh2QAV0PtPtK0sPeak_%d",i),"QA: K0s: Daughter Pt vs Pt;neg pt;pos pt",100,0,5,100,0,5);
1061
87e23d3e 1062 fh2QAV0EtaPtLambdaPeak[i] = new TH2D(Form("fh2QAV0EtaPtLambdaPeak_%d",i),"QA: Lambda: Daughter Eta vs V0 pt, peak;track eta;V0 pt",200,-2,2,iNBinsPtV0,dPtV0Min,dPtV0Max);
7729cc34 1063 fh2QAV0EtaEtaLambda[i] = new TH2D(Form("fh2QAV0EtaEtaLambda_%d",i),"QA: Lambda: Eta vs Eta Daughter",200,-2,2,200,-2,2);
1064 fh2QAV0PhiPhiLambda[i] = new TH2D(Form("fh2QAV0PhiPhiLambda_%d",i),"QA: Lambda: Phi vs Phi Daughter",200,0,TMath::TwoPi(),200,0,TMath::TwoPi());
1065 fh1QAV0RapLambda[i] = new TH1D(Form("fh1QAV0RapLambda_%d",i),"QA: Lambda: V0 Rapidity",200,-2,2);
1066 fh2QAV0PtPtLambdaPeak[i] = new TH2D(Form("fh2QAV0PtPtLambdaPeak_%d",i),"QA: Lambda: Daughter Pt vs Pt;neg pt;pos pt",100,0,5,100,0,5);
1067
1068 fh1QAV0Pt[i] = new TH1D(Form("fh1QAV0Pt_%d",i),"QA: Daughter Pt",100,0,5);
1069 fh1QAV0Charge[i] = new TH1D(Form("fh1QAV0Charge_%d",i),"QA: V0 Charge",3,-1,2);
1070 fh1QAV0DCAVtx[i] = new TH1D(Form("fh1QAV0DCAVtx_%d",i),"QA: DCA daughters to primary vertex",100,0,10);
1071 fh1QAV0DCAV0[i] = new TH1D(Form("fh1QAV0DCAV0_%d",i),"QA: DCA daughters",100,0,2);
1072 fh1QAV0Cos[i] = new TH1D(Form("fh1QAV0Cos_%d",i),"QA: CPA",10000,0.9,1);
1073 fh1QAV0R[i] = new TH1D(Form("fh1QAV0R_%d",i),"QA: R",1500,0,150);
1074 fh1QACTau2D[i] = new TH1D(Form("fh1QACTau2D_%d",i),"QA: K0s: c#tau 2D;mR/pt#tau",100,0,10);
1075 fh1QACTau3D[i] = new TH1D(Form("fh1QACTau3D_%d",i),"QA: K0s: c#tau 3D;mL/p#tau",100,0,10);
1076
1077 fh2ArmPod[i] = new TH2D(Form("fh2ArmPod_%d",i),"Armenteros-Podolanski;#alpha;#it{p}_{T}^{Arm}",100,-1.,1.,50,0.,0.25);
1078 fh2ArmPodK0s[i] = new TH2D(Form("fh2ArmPodK0s_%d",i),"K0s: Armenteros-Podolanski;#alpha;#it{p}_{T}^{Arm}",100,-1.,1.,50,0.,0.25);
1079 fh2ArmPodLambda[i] = new TH2D(Form("fh2ArmPodLambda_%d",i),"Lambda: Armenteros-Podolanski;#alpha;#it{p}_{T}^{Arm}",100,-1.,1.,50,0.,0.25);
1080 fh2ArmPodALambda[i] = new TH2D(Form("fh2ArmPodALambda_%d",i),"ALambda: Armenteros-Podolanski;#alpha;#it{p}_{T}^{Arm}",100,-1.,1.,50,0.,0.25);
1081
87e23d3e 1082 fh2CutTPCRowsK0s[i] = new TH2D(Form("fh2CutTPCRowsK0s_%d",i),"Cuts: K0s: TPC Rows vs mass;#it{m}_{inv} (GeV/#it{c}^{2});TPC rows",fgkiNBinsMassK0s,fgkdMassK0sMin,fgkdMassK0sMax,160,0,160);
1083 fh2CutTPCRowsLambda[i] = new TH2D(Form("fh2CutTPCRowsLambda_%d",i),"Cuts: Lambda: TPC Rows vs mass;#it{m}_{inv} (GeV/#it{c}^{2});TPC rows",fgkiNBinsMassLambda,fgkdMassLambdaMin,fgkdMassLambdaMax,160,0,160);
1084 fh2CutPtPosK0s[i] = new TH2D(Form("fh2CutPtPosK0s_%d",i),"Cuts: K0s: Pt pos;#it{m}_{inv} (GeV/#it{c}^{2});pt pos",fgkiNBinsMassK0s,fgkdMassK0sMin,fgkdMassK0sMax,100,0,5);
1085 fh2CutPtNegK0s[i] = new TH2D(Form("fh2CutPtNegK0s_%d",i),"Cuts: K0s: Pt neg;#it{m}_{inv} (GeV/#it{c}^{2});pt neg",fgkiNBinsMassK0s,fgkdMassK0sMin,fgkdMassK0sMax,100,0,5);
1086 fh2CutPtPosLambda[i] = new TH2D(Form("fh2CutPtPosLambda_%d",i),"Cuts: Lambda: Pt pos;#it{m}_{inv} (GeV/#it{c}^{2});pt pos",fgkiNBinsMassLambda,fgkdMassLambdaMin,fgkdMassLambdaMax,100,0,5);
1087 fh2CutPtNegLambda[i] = new TH2D(Form("fh2CutPtNegLambda_%d",i),"Cuts: Lambda: Pt neg;#it{m}_{inv} (GeV/#it{c}^{2});pt neg",fgkiNBinsMassLambda,fgkdMassLambdaMin,fgkdMassLambdaMax,100,0,5);
1088 fh2CutDCAVtx[i] = new TH2D(Form("fh2CutDCAVtx_%d",i),"Cuts: DCA daughters to prim. vtx.;#it{m}_{inv} (GeV/#it{c}^{2});DCA daughter to prim. vtx. (cm)",fgkiNBinsMassK0s,fgkdMassK0sMin,fgkdMassK0sMax,100,0,10);
1089 fh2CutDCAV0[i] = new TH2D(Form("fh2CutDCAV0_%d",i),"Cuts: DCA daughters;#it{m}_{inv} (GeV/#it{c}^{2});DCA daughters / #sigma_{TPC}",fgkiNBinsMassK0s,fgkdMassK0sMin,fgkdMassK0sMax,100,0,2);
1090 fh2CutCos[i] = new TH2D(Form("fh2CutCos_%d",i),"Cuts: CPA;#it{m}_{inv} (GeV/#it{c}^{2});CPA",fgkiNBinsMassK0s,fgkdMassK0sMin,fgkdMassK0sMax,10000,0.9,1);
1091 fh2CutR[i] = new TH2D(Form("fh2CutR_%d",i),"Cuts: R;#it{m}_{inv} (GeV/#it{c}^{2});R (cm)",fgkiNBinsMassK0s,fgkdMassK0sMin,fgkdMassK0sMax,1500,0,150);
1092 fh2CutEtaK0s[i] = new TH2D(Form("fh2CutEtaK0s_%d",i),"Cuts: K0s: Eta;#it{m}_{inv} (GeV/#it{c}^{2});#eta",fgkiNBinsMassK0s,fgkdMassK0sMin,fgkdMassK0sMax,200,-2,2);
1093 fh2CutEtaLambda[i] = new TH2D(Form("fh2CutEtaLambda_%d",i),"Cuts: Lambda: Eta;#it{m}_{inv} (GeV/#it{c}^{2});#eta",fgkiNBinsMassLambda,fgkdMassLambdaMin,fgkdMassLambdaMax,200,-2,2);
1094 fh2CutRapK0s[i] = new TH2D(Form("fh2CutRapK0s_%d",i),"Cuts: K0s: Rapidity;#it{m}_{inv} (GeV/#it{c}^{2});y",fgkiNBinsMassK0s,fgkdMassK0sMin,fgkdMassK0sMax,200,-2,2);
1095 fh2CutRapLambda[i] = new TH2D(Form("fh2CutRapLambda_%d",i),"Cuts: Lambda: Rapidity;#it{m}_{inv} (GeV/#it{c}^{2});y",fgkiNBinsMassLambda,fgkdMassLambdaMin,fgkdMassLambdaMax,200,-2,2);
1096 fh2CutCTauK0s[i] = new TH2D(Form("fh2CutCTauK0s_%d",i),"Cuts: K0s: #it{c#tau};#it{m}_{inv} (GeV/#it{c}^{2});#it{mL/p#tau}",fgkiNBinsMassK0s,fgkdMassK0sMin,fgkdMassK0sMax,100,0,10);
1097 fh2CutCTauLambda[i] = new TH2D(Form("fh2CutCTauLambda_%d",i),"Cuts: Lambda: #it{c#tau};#it{m}_{inv} (GeV/#it{c}^{2});#it{mL/p#tau}",fgkiNBinsMassLambda,fgkdMassLambdaMin,fgkdMassLambdaMax,100,0,10);
1098 fh2CutPIDPosK0s[i] = new TH2D(Form("fh2CutPIDPosK0s_%d",i),"Cuts: K0s: PID pos;#it{m}_{inv} (GeV/#it{c}^{2});##sigma_{d#it{E}/d#it{x}}",fgkiNBinsMassK0s,fgkdMassK0sMin,fgkdMassK0sMax,100,0,10);
1099 fh2CutPIDNegK0s[i] = new TH2D(Form("fh2CutPIDNegK0s_%d",i),"Cuts: K0s: PID neg;#it{m}_{inv} (GeV/#it{c}^{2});##sigma_{d#it{E}/d#it{x}}",fgkiNBinsMassK0s,fgkdMassK0sMin,fgkdMassK0sMax,100,0,10);
1100 fh2CutPIDPosLambda[i] = new TH2D(Form("fh2CutPIDPosLambda_%d",i),"Cuts: Lambda: PID pos;#it{m}_{inv} (GeV/#it{c}^{2});##sigma_{d#it{E}/d#it{x}}",fgkiNBinsMassLambda,fgkdMassLambdaMin,fgkdMassLambdaMax,100,0,10);
1101 fh2CutPIDNegLambda[i] = new TH2D(Form("fh2CutPIDNegLambda_%d",i),"Cuts: Lambda: PID neg;#it{m}_{inv} (GeV/#it{c}^{2});##sigma_{d#it{E}/d#it{x}}",fgkiNBinsMassLambda,fgkdMassLambdaMin,fgkdMassLambdaMax,100,0,10);
7729cc34 1102 fh2Tau3DVs2D[i] = new TH2D(Form("fh2Tau3DVs2D_%d",i),"Decay 3D vs 2D;pt;3D/2D",100,0,10,200,0.5,1.5);
1103
1104 fOutputListQA->Add(fh1QAV0Status[i]);
1105 fOutputListQA->Add(fh1QAV0TPCRefit[i]);
1106 fOutputListQA->Add(fh1QAV0TPCRows[i]);
1107 fOutputListQA->Add(fh1QAV0TPCFindable[i]);
1108 fOutputListQA->Add(fh1QAV0TPCRowsFind[i]);
1109 fOutputListQA->Add(fh1QAV0Eta[i]);
1110 fOutputListQA->Add(fh2QAV0EtaRows[i]);
1111 fOutputListQA->Add(fh2QAV0PtRows[i]);
1112 fOutputListQA->Add(fh2QAV0PhiRows[i]);
1113 fOutputListQA->Add(fh2QAV0NClRows[i]);
1114 fOutputListQA->Add(fh2QAV0EtaNCl[i]);
1115
1116 fOutputListQA->Add(fh2QAV0EtaPtK0sPeak[i]);
1117 fOutputListQA->Add(fh2QAV0EtaEtaK0s[i]);
1118 fOutputListQA->Add(fh2QAV0PhiPhiK0s[i]);
1119 fOutputListQA->Add(fh1QAV0RapK0s[i]);
1120 fOutputListQA->Add(fh2QAV0PtPtK0sPeak[i]);
1121
1122 fOutputListQA->Add(fh2QAV0EtaPtLambdaPeak[i]);
1123 fOutputListQA->Add(fh2QAV0EtaEtaLambda[i]);
1124 fOutputListQA->Add(fh2QAV0PhiPhiLambda[i]);
1125 fOutputListQA->Add(fh1QAV0RapLambda[i]);
1126 fOutputListQA->Add(fh2QAV0PtPtLambdaPeak[i]);
1127
1128 fOutputListQA->Add(fh1QAV0Pt[i]);
1129 fOutputListQA->Add(fh1QAV0Charge[i]);
1130 fOutputListQA->Add(fh1QAV0DCAVtx[i]);
1131 fOutputListQA->Add(fh1QAV0DCAV0[i]);
1132 fOutputListQA->Add(fh1QAV0Cos[i]);
1133 fOutputListQA->Add(fh1QAV0R[i]);
1134 fOutputListQA->Add(fh1QACTau2D[i]);
1135 fOutputListQA->Add(fh1QACTau3D[i]);
1136
1137 fOutputListQA->Add(fh2ArmPod[i]);
1138 fOutputListQA->Add(fh2ArmPodK0s[i]);
1139 fOutputListQA->Add(fh2ArmPodLambda[i]);
1140 fOutputListQA->Add(fh2ArmPodALambda[i]);
1141
1142 fOutputListCuts->Add(fh2CutTPCRowsK0s[i]);
1143 fOutputListCuts->Add(fh2CutTPCRowsLambda[i]);
1144 fOutputListCuts->Add(fh2CutPtPosK0s[i]);
1145 fOutputListCuts->Add(fh2CutPtNegK0s[i]);
1146 fOutputListCuts->Add(fh2CutPtPosLambda[i]);
1147 fOutputListCuts->Add(fh2CutPtNegLambda[i]);
1148 fOutputListCuts->Add(fh2CutDCAVtx[i]);
1149 fOutputListCuts->Add(fh2CutDCAV0[i]);
1150 fOutputListCuts->Add(fh2CutCos[i]);
1151 fOutputListCuts->Add(fh2CutR[i]);
1152 fOutputListCuts->Add(fh2CutEtaK0s[i]);
1153 fOutputListCuts->Add(fh2CutEtaLambda[i]);
1154 fOutputListCuts->Add(fh2CutRapK0s[i]);
1155 fOutputListCuts->Add(fh2CutRapLambda[i]);
1156 fOutputListCuts->Add(fh2CutCTauK0s[i]);
1157 fOutputListCuts->Add(fh2CutCTauLambda[i]);
1158 fOutputListCuts->Add(fh2CutPIDPosK0s[i]);
1159 fOutputListCuts->Add(fh2CutPIDNegK0s[i]);
1160 fOutputListCuts->Add(fh2CutPIDPosLambda[i]);
1161 fOutputListCuts->Add(fh2CutPIDNegLambda[i]);
1162 fOutputListCuts->Add(fh2Tau3DVs2D[i]);
1163 }
1164
1165 for (Int_t i = 0; i < fgkiNCategV0; i++)
1166 {
87e23d3e 1167 fh1V0InvMassK0sAll[i] = new TH1D(Form("fh1V0InvMassK0sAll_%d",i), Form("K0s: V0 invariant mass, %s;#it{m}_{inv} (GeV/#it{c}^{2});counts",categV0[i].Data()),fgkiNBinsMassK0s,fgkdMassK0sMin,fgkdMassK0sMax);
1168 fh1V0InvMassLambdaAll[i] = new TH1D(Form("fh1V0InvMassLambdaAll_%d",i), Form("Lambda: V0 invariant mass, %s;#it{m}_{inv} (GeV/#it{c}^{2});counts",categV0[i].Data()),fgkiNBinsMassLambda,fgkdMassLambdaMin,fgkdMassLambdaMax);
1169 fh1V0InvMassALambdaAll[i] = new TH1D(Form("fh1V0InvMassALambdaAll_%d",i), Form("ALambda: V0 invariant mass, %s;#it{m}_{inv} (GeV/#it{c}^{2});counts",categV0[i].Data()),fgkiNBinsMassLambda,fgkdMassLambdaMin,fgkdMassLambdaMax);
7729cc34 1170 fOutputListStd->Add(fh1V0InvMassK0sAll[i]);
1171 fOutputListStd->Add(fh1V0InvMassLambdaAll[i]);
1172 fOutputListStd->Add(fh1V0InvMassALambdaAll[i]);
1173 }
1174
1175 for (Int_t i = 0; i < fOutputListStd->GetEntries(); ++i)
1176 {
1177 TH1* h1 = dynamic_cast<TH1*>(fOutputListStd->At(i));
1178 if (h1)
1179 {
1180 h1->Sumw2();
1181 continue;
1182 }
1183 THnSparse* hn = dynamic_cast<THnSparse*>(fOutputListStd->At(i));
1184 if(hn) hn->Sumw2();
1185 }
1186 for (Int_t i = 0; i < fOutputListQA->GetEntries(); ++i)
1187 {
1188 TH1* h1 = dynamic_cast<TH1*>(fOutputListQA->At(i));
1189 if (h1)
1190 {
1191 h1->Sumw2();
1192 continue;
1193 }
1194 THnSparse* hn = dynamic_cast<THnSparse*>(fOutputListQA->At(i));
1195 if(hn) hn->Sumw2();
1196 }
1197 for (Int_t i = 0; i < fOutputListCuts->GetEntries(); ++i)
1198 {
1199 TH1* h1 = dynamic_cast<TH1*>(fOutputListCuts->At(i));
1200 if (h1)
1201 {
1202 h1->Sumw2();
1203 continue;
1204 }
1205 THnSparse* hn = dynamic_cast<THnSparse*>(fOutputListCuts->At(i));
1206 if(hn) hn->Sumw2();
1207 }
1208 for (Int_t i = 0; i < fOutputListMC->GetEntries(); ++i)
1209 {
1210 TH1* h1 = dynamic_cast<TH1*>(fOutputListMC->At(i));
1211 if (h1)
1212 {
1213 h1->Sumw2();
1214 continue;
1215 }
1216 THnSparse* hn = dynamic_cast<THnSparse*>(fOutputListMC->At(i));
1217 if(hn) hn->Sumw2();
1218 }
1219
1220 PostData(1,fOutputListStd);
1221 PostData(2,fOutputListQA);
1222 PostData(3,fOutputListCuts);
1223 PostData(4,fOutputListMC);
1224// if (fbTreeOutput)
1225// PostData(5,ftreeOut);
1226}
1227
1228void AliAnalysisTaskV0sInJets::UserExec(Option_t *)
1229{
1230 // Main loop, called for each event
1231 if(fDebug>5) printf("TaskV0sInJets: UserExec: Start\n");
9706bde3 1232 /*
1233 // reset branches for each event
1234 if (fBranchV0Rec)
1235 fBranchV0Rec->Clear();
1236 if (fBranchV0Gen)
1237 fBranchV0Gen->Clear();
1238 if (fBranchJet)
1239 fBranchJet->Clear();
1240 if (fEventInfo)
1241 fEventInfo->Reset();
1242 */
7729cc34 1243 if (!fiAODAnalysis)
1244 return;
1245
1246 if(fDebug>2) printf("TaskV0sInJets: AOD analysis\n");
1247 fh1EventCounterCut->Fill(0); // all available selected events (collision candidates)
1248
1249 if(fDebug>5) printf("TaskV0sInJets: UserExec: Loading AOD\n");
1250 fAODIn = dynamic_cast<AliAODEvent*>(InputEvent()); // input AOD
1251 fAODOut = AODEvent(); // output AOD
1252 if (!fAODOut)
1253 {
87e23d3e 1254 if(fDebug>0) printf("TaskV0sInJets: No output AOD found\n");
7729cc34 1255 return;
1256 }
1257 if (!fAODIn)
1258 {
87e23d3e 1259 if(fDebug>0) printf("TaskV0sInJets: No input AOD found\n");
7729cc34 1260 return;
1261 }
1262 if(fDebug>5) printf("TaskV0sInJets: UserExec: Loading AOD OK\n");
1263
1264 TClonesArray* arrayMC = 0; // array particles in the MC event
1265 AliAODMCHeader* headerMC = 0; // MC header
1266 Int_t iNTracksMC = 0; // number of MC tracks
1267 Double_t dPrimVtxMCX=0., dPrimVtxMCY=0., dPrimVtxMCZ=0.; // position of the MC primary vertex
1268
1269 // Simulation info
1270 if (fbMCAnalysis)
1271 {
1272 arrayMC = (TClonesArray*)fAODIn->FindListObject(AliAODMCParticle::StdBranchName());
1273 if (!arrayMC)
1274 {
87e23d3e 1275 if(fDebug>0) printf("TaskV0sInJets: No MC array found\n");
7729cc34 1276 return;
1277 }
1278 if(fDebug>5) printf("TaskV0sInJets: MC array found\n");
1279 iNTracksMC = arrayMC->GetEntriesFast();
1280 if(fDebug>5) printf("TaskV0sInJets: There are %d MC tracks in this event\n",iNTracksMC);
1281// if (!iNTracksMC)
1282// return;
1283 headerMC = (AliAODMCHeader*)fAODIn->FindListObject(AliAODMCHeader::StdBranchName());
1284 if (!headerMC)
1285 {
87e23d3e 1286 if(fDebug>0) printf("TaskV0sInJets: No MC header found\n");
7729cc34 1287 return;
1288 }
1289 // get position of the MC primary vertex
1290 dPrimVtxMCX=headerMC->GetVtxX();
1291 dPrimVtxMCY=headerMC->GetVtxY();
1292 dPrimVtxMCZ=headerMC->GetVtxZ();
1293 }
1294
1295 // PID Response Task object
1296 AliAnalysisManager* mgr = AliAnalysisManager::GetAnalysisManager();
1297 AliInputEventHandler* inputHandler = (AliInputEventHandler*)mgr->GetInputEventHandler();
1298 AliPIDResponse* fPIDResponse = inputHandler->GetPIDResponse();
1299 if (!fPIDResponse)
1300 {
87e23d3e 1301 if(fDebug>0) printf("TaskV0sInJets: No PID response object found\n");
7729cc34 1302 return;
1303 }
1304
1305 // AOD files are OK
1306 fh1EventCounterCut->Fill(1);
1307
1308 // Event selection
87e23d3e 1309 if (!IsSelectedForJets(fAODIn,fdCutVertexZ,fdCutVertexR2,fdCutCentLow,fdCutCentHigh,1,0.1)) // cut on |delta z| in 2011 data between SPD vertex and nominal primary vertex
1310// if (!IsSelectedForJets(fAODIn,fdCutVertexZ,fdCutVertexR2,fdCutCentLow,fdCutCentHigh)) // no need for cutting in 2010 data
7729cc34 1311 {
1312 if(fDebug>5) printf("TaskV0sInJets: Event rejected\n");
1313 return;
1314 }
1315
87e23d3e 1316// fdCentrality = fAODIn->GetHeader()->GetCentrality(); // event centrality
1317 fdCentrality = fAODIn->GetHeader()->GetCentralityP()->GetCentralityPercentile("V0M"); // event centrality
1318 Int_t iCentIndex = GetCentralityBinIndex(fdCentrality); // get index of centrality bin
9706bde3 1319 if (iCentIndex<0)
1320 {
1321 if(fDebug>5) printf("TaskV0sInJets: Event is out of histogram range\n");
1322 return;
1323 }
7729cc34 1324 fh1EventCounterCut->Fill(2); // selected events (vertex, centrality)
1325 fh1EventCounterCutCent[iCentIndex]->Fill(2);
1326
1327 UInt_t iNTracks = fAODIn->GetNumberOfTracks(); // get number of tracks in event
1328 if(fDebug>5) printf("TaskV0sInJets: There are %d tracks in this event\n",iNTracks);
1329// if (!iNTracks)
1330// return;
1331
1332 Int_t iNV0s = fAODIn->GetNumberOfV0s(); // get the number of V0 candidates
1333 if (!iNV0s)
1334 {
87e23d3e 1335 if(fDebug>2) printf("TaskV0sInJets: No V0s found in event\n");
7729cc34 1336// return;
1337 }
1338
1339 /*===== Event is OK for the analysis =====*/
1340 fh1EventCent->Fill(iCentIndex);
87e23d3e 1341 fh1EventCent2->Fill(fdCentrality);
1342 fh2EventCentTracks->Fill(fdCentrality,iNTracks);
7729cc34 1343
1344// if (fbTreeOutput)
1345// fEventInfo->SetAODEvent(fAODIn);
1346// printf("V0sInJets: EventInfo: Centrality: %f\n",fEventInfo->GetCentrality());
1347
1348 if (iNV0s)
1349 {
1350 fh1EventCounterCut->Fill(3); // events with V0s
1351 fh1EventCounterCutCent[iCentIndex]->Fill(3);
1352 }
1353
1354// Int_t iNV0SelV0Rec = 0;
1355// Int_t iNV0SelV0Gen = 0;
1356
7729cc34 1357 AliAODv0* v0 = 0; // pointer to V0 candidates
1358// AliV0Object* objectV0 = 0;
1359 TVector3 vecV0Momentum; // 3D vector of V0 momentum
87e23d3e 1360 Double_t dMassV0K0s = 0; // invariant mass of the K0s candidate
1361 Double_t dMassV0Lambda = 0; // invariant mass of the Lambda candidate
1362 Double_t dMassV0ALambda = 0; // invariant mass of the Lambda candidate
7729cc34 1363 Int_t iNV0CandTot = 0; // counter of all V0 candidates at the beginning
1364 Int_t iNV0CandK0s = 0; // counter of K0s candidates at the end
1365 Int_t iNV0CandLambda = 0; // counter of Lambda candidates at the end
1366 Int_t iNV0CandALambda = 0; // counter of Lambda candidates at the end
1367
1368 Bool_t bUseOldCuts = 0; // old reconstruction cuts
1369 Bool_t bUseAliceCuts = 0; // cuts used by Alice Zimmermann
1370 Bool_t bUseIouriCuts = 0; // cuts used by Iouri
87e23d3e 1371 Bool_t bPrintCuts = 0; // print out which cuts are applied
7729cc34 1372 Bool_t bPrintJetSelection = 0; // print out which jets are selected
1373
1374 // Values of V0 reconstruction cuts:
1375 // Daughter tracks
1376 Int_t iRefit = AliAODTrack::kTPCrefit; // TPC refit for daughter tracks
87e23d3e 1377 Double_t dDCAToPrimVtxMin = fdCutDCAToPrimVtxMin; // 0.1; // [cm] min DCA of daughters to the prim vtx
1378 Double_t dDCADaughtersMax = fdCutDCADaughtersMax; // 1.; // [sigma of TPC tracking] max DCA between daughters
1379 Double_t dEtaDaughterMax = 0.8; // max |pseudorapidity| of daughter tracks
1380 Double_t dNSigmadEdxMax = fdCutNSigmadEdxMax;// 3.; // [sigma dE/dx] max difference between measured and expected signal of dE/dx in the TPC
1381 Double_t dPtProtonPIDMax = 1.; // [GeV/c] maxium pT of proton for applying PID cut
7729cc34 1382 // V0 candidate
1383 Bool_t bOnFly = 0; // on-the-fly (yes) or offline (no) reconstructed
87e23d3e 1384 Double_t dCPAMin = fdCutCPAMin;// 0.998; // min cosine of the pointing angle
1385 Double_t dRadiusDecayMin = 5.; // [cm] min radial distance of the decay vertex
1386 Double_t dRadiusDecayMax = 100.; // [cm] max radial distance of the decay vertex
1387 Double_t dEtaMax = 0.7; // max |pseudorapidity| of V0
1388 Double_t dNTauMax = fdCutNTauMax; // 5.0; // [tau] max proper lifetime in multiples of the mean lifetime
7729cc34 1389
1390 // Old cuts Start
87e23d3e 1391 Double_t dNCrossedRowsTPCMin = 70.; // min number of crossed TPC rows (turned off)
1392// Double_t dCrossedRowsOverFindMin = 0.8; // min ratio crossed rows / findable clusters (turned off)
1393// Double_t dCrossedRowsOverFindMax = 1e3; // max ratio crossed rows / findable clusters (turned off)
1394 Double_t dPtDaughterMin = 0.150; // [GeV/c] min transverse momentum of daughter tracks (turned off)
1395 Double_t dRapMax = 0.75; // max |rapidity| of V0 (turned off)
7729cc34 1396 // Old cuts End
1397
1398 // Other cuts
87e23d3e 1399 Double_t dNSigmaMassMax = 3.; // [sigma m] max difference between candidate mass and real particle mass (used only for mass peak method of signal extraction)
1400 Double_t dDistPrimaryMax = 0.01; // [cm] max distance of production point to the primary vertex (criterion for choice of MC particles considered as primary)
7729cc34 1401
1402 // Selection of active cuts
1403 Bool_t bCutEtaDaughter = 1; // daughter pseudorapidity
1404 Bool_t bCutRapV0 = 0; // V0 rapidity
1405 Bool_t bCutEtaV0 = 1; // V0 pseudorapidity
1406 Bool_t bCutTau = 1; // V0 lifetime
1407 Bool_t bCutPid = 1; // PID (TPC dE/dx)
1408 Bool_t bCutArmPod = 1; // Armenteros-Podolanski for K0S
1409// Bool_t bCutCross = 0; // cross contamination
1410
1411 if (bUseOldCuts)
1412 {
1413 bCutRapV0 = 1;
87e23d3e 1414 dEtaMax = 0.75;
1415 dNTauMax = 3.0;
7729cc34 1416 }
1417 else if (bUseAliceCuts)
1418 {
1419// bOnFly = 1;
87e23d3e 1420 dEtaMax = 0.75;
1421 dNTauMax = 5.0;
7729cc34 1422 }
1423 else if (bUseIouriCuts)
1424 {
1425 bCutRapV0 = 1;
1426 bCutEtaV0 = 0;
87e23d3e 1427 dNTauMax = 3.0;
1428 dRapMax = 0.5;
7729cc34 1429 }
1430
87e23d3e 1431 Double_t dCTauK0s = 2.6844; // [cm] c tau of K0S
1432 Double_t dCTauLambda = 7.89; // [cm] c tau of Lambda
7729cc34 1433
1434 // Load PDG values of particle masses
87e23d3e 1435 Double_t dMassPDGK0s = TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass();
1436 Double_t dMassPDGLambda = TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass();
7729cc34 1437
1438 // PDG codes of used particles
1439 Int_t iPdgCodePion = 211;
1440 Int_t iPdgCodeProton = 2212;
1441 Int_t iPdgCodeK0s = 310;
1442 Int_t iPdgCodeLambda = 3122;
1443
87e23d3e 1444 // Jet selection: fdCutPtJetMin, fdCutPtTrackMin
1445 Double_t dJetEtaWindow = dEtaMax-fdRadiusJet; // max jet |pseudorapidity|, to make sure that V0s can appear in the entire jet area
1446 Double_t dCutJetAreaMin = 0.6*TMath::Pi()*fdRadiusJet*fdRadiusJet; // minimum jet area
1447 Double_t dRadiusExcludeCone = 2*fdRadiusJet; // radius of cones around jets excluded for V0 outside jets
7729cc34 1448 Bool_t bLeadingJetOnly = 0;
1449
1450 if (bUseAliceCuts)
1451 {
87e23d3e 1452 fdCutPtJetMin = 5;
1453 fdCutPtTrackMin = 5;
1454 dCutJetAreaMin = 0;
7729cc34 1455 bLeadingJetOnly = 0;
1456 }
1457
1458// Int_t iNJetAll = 0; // number of reconstructed jets in fBranchJet
1459// iNJetAll = fBranchJet->GetEntriesFast(); // number of reconstructed jets
1460 TClonesArray* jetArray = 0; // object where the input jets are stored
9706bde3 1461 TClonesArray* jetArrayBg = 0; // object where the kt clusters are stored
7729cc34 1462 Int_t iNJet = 0; // number of reconstructed jets in the input
1463 TClonesArray* jetArraySel = new TClonesArray("AliAODJet",0); // object where the selected jets are copied
1464 Int_t iNJetSel = 0; // number of selected reconstructed jets
1465// iNJetSel = jetArraySel->GetEntriesFast(); // number of selected reconstructed jets
1466 TClonesArray* jetArrayPerp = new TClonesArray("AliAODJet",0); // object where the perp. cones are stored
1467 Int_t iNJetPerp = 0; // number of perpendicular cones
1468
1469 AliAODJet* jet = 0; // pointer to a jet
1470// AliJetObject* objectJet = 0;
1471 AliAODJet* jetPerp = 0; // pointer to a perp. cone
1472 AliAODJet* jetRnd = 0; // pointer to a rand. cone
9706bde3 1473 AliAODJet* jetMed = 0; // pointer to a median cluster
7729cc34 1474 TVector3 vecJetMomentum; // 3D vector of jet momentum
1475// TVector3 vecPerpConeMomentum; // 3D vector of perpendicular cone momentum
1476// TVector3 vecRndConeMomentum; // 3D vector of random cone momentum
1477 Bool_t bJetEventGood = kTRUE; // indicator of good jet events
1478
1479// printf("iNJetAll, iNJetSel: %d %d\n",iNJetAll,iNJetSel);
1480
1481 if (fbJetSelection) // analysis of V0s in jets is switched on
1482 {
1483 jetArray = (TClonesArray*)(fAODOut->FindListObject(fsJetBranchName.Data())); // find object with jets in the output AOD
1484 if (!jetArray)
1485 {
87e23d3e 1486 if(fDebug>0) printf("TaskV0sInJets: No array of name: %s\n",fsJetBranchName.Data());
7729cc34 1487 bJetEventGood = kFALSE;
1488 }
1489 if (bJetEventGood)
1490 iNJet = jetArray->GetEntriesFast();
1491 if (bJetEventGood && !iNJet) // check whether there are some jets
1492 {
87e23d3e 1493 if(fDebug>2) printf("TaskV0sInJets: No jets in array\n");
7729cc34 1494 bJetEventGood = kFALSE;
1495 }
9706bde3 1496 if (bJetEventGood)
1497 {
1498// printf("TaskV0sInJets: Loading bg array of name: %s\n",fsJetBgBranchName.Data());
1499 jetArrayBg = (TClonesArray*)(fAODOut->FindListObject(fsJetBgBranchName.Data())); // find object with jets in the output AOD
1500 if (!jetArrayBg)
1501 {
1502 if(fDebug>0) printf("TaskV0sInJets: No bg array of name: %s\n",fsJetBgBranchName.Data());
1503// bJetEventGood = kFALSE;
1504 }
1505 }
7729cc34 1506 }
1507 else // no in-jet analysis
1508 bJetEventGood = kFALSE;
1509
1510 // select good jets and copy them to another array
1511 if (bJetEventGood)
1512 {
1513 if (bLeadingJetOnly)
1514 iNJet = 1; // only leading jets
87e23d3e 1515 if(fDebug>5) printf("TaskV0sInJets: Jet selection for %d jets\n",iNJet);
7729cc34 1516 for (Int_t iJet = 0; iJet<iNJet; iJet++)
1517 {
1518 AliAODJet* jetSel = (AliAODJet*)jetArray->At(iJet); // load a jet in the list
1519 if (!jetSel)
1520 {
87e23d3e 1521 if(fDebug>0) printf("TaskV0sInJets: Cannot load jet %d\n",iJet);
7729cc34 1522 continue;
1523 }
1524 if (bPrintJetSelection)
87e23d3e 1525 if(fDebug>7) printf("jet: i = %d, pT = %f, eta = %f, phi = %f, pt lead tr = %f ",iJet,jetSel->Pt(),jetSel->Eta(),jetSel->Phi(),jetSel->GetPtLeading());
1526// printf("TaskV0sInJets: Checking pt > %.2f for jet %d with pt %.2f\n",fdCutPtJetMin,iJet,jetSel->Pt());
1527 if (jetSel->Pt() < fdCutPtJetMin) // selection of high-pt jets
7729cc34 1528 {
1529 if (bPrintJetSelection)
87e23d3e 1530 if(fDebug>7) printf("rejected (pt)\n");
7729cc34 1531 continue;
1532 }
87e23d3e 1533// printf("TaskV0sInJets: Checking |eta| < %.2f for jet %d with |eta| %.2f\n",dJetEtaWindow,iJet,TMath::Abs(jetSel->Eta()));
1534 if (TMath::Abs(jetSel->Eta()) > dJetEtaWindow) // selection of jets in the chosen pseudorapidity range
7729cc34 1535 {
1536 if (bPrintJetSelection)
87e23d3e 1537 if(fDebug>7) printf("rejected (eta)\n");
7729cc34 1538 continue;
1539 }
1540 if (!bUseOldCuts)
1541 {
87e23d3e 1542 if (jetSel->EffectiveAreaCharged() < dCutJetAreaMin)
7729cc34 1543 continue;
1544 }
1545 Int_t iNTracksInJet = 0;
87e23d3e 1546 Double_t dPtLeadTrack = 0; // pt of the leading track
7729cc34 1547// Int_t iLeadTrack = 0;
1548 iNTracksInJet = jetSel->GetRefTracks()->GetEntriesFast(); // number od tracks that constitute the jet
1549// printf("TaskV0sInJets: Searching for leading track from %d tracks in jet %d\n",iNTracksInJet,iJet);
87e23d3e 1550 if (fdCutPtTrackMin > 0) // a positive min leading track pt is set
7729cc34 1551 {
1552 for (Int_t j = 0; j < iNTracksInJet; j++) // find the track with the highest pt
1553 {
1554 AliAODTrack* track = (AliAODTrack*)jetSel->GetTrack(j); // is this the leading track?
1555 if (!track)
1556 continue;
1557// printf("TaskV0sInJets: %d: %.2f\n",j,track->Pt());
87e23d3e 1558 if (track->Pt() > dPtLeadTrack)
7729cc34 1559 {
87e23d3e 1560 dPtLeadTrack = track->Pt();
7729cc34 1561// iLeadTrack = j;
1562 }
1563 }
87e23d3e 1564// printf("Leading track pT: my: %f, ali: %f\n",dPtLeadTrack,jetSel->GetPtLeading());
1565// printf("TaskV0sInJets: Checking leading track pt > %.2f for pt %.2f of track %d in jet %d\n",fdCutPtTrackMin,dPtLeadTrack,iLeadTrack,iJet);
1566 if (dPtLeadTrack < fdCutPtTrackMin) // selection of high-pt jet-track events
7729cc34 1567 {
1568 if (bPrintJetSelection)
87e23d3e 1569 if(fDebug>7) printf("rejected (track pt)\n");
7729cc34 1570 continue;
1571 }
1572 }
1573 if (bPrintJetSelection)
87e23d3e 1574 if(fDebug>7) printf("accepted\n");
7729cc34 1575 if(fDebug>5) printf("TaskV0sInJets: Jet %d with pt %.2f passed selection\n",iJet,jetSel->Pt());
9706bde3 1576 /*
1577 if (fbTreeOutput)
1578 {
1579 // new ((*fBranchJet)[iNJetAll++]) AliAODJet(*((AliAODJet*)jetSel));
1580 objectJet = new ((*fBranchJet)[iNJetAll++]) AliJetObject(jetSel); // copy selected jet to the array
1581 // objectJet->SetPtLeadingTrack(dPtLeadTrack);
1582 objectJet->SetRadius(fdRadiusJet);
1583 objectJet = 0;
1584 }
1585 */
7729cc34 1586 TLorentzVector vecPerpPlus(*(jetSel->MomentumVector()));
1587 vecPerpPlus.RotateZ(TMath::Pi()/2.); // rotate vector by 90 deg around z
1588 TLorentzVector vecPerpMinus(*(jetSel->MomentumVector()));
1589 vecPerpMinus.RotateZ(-TMath::Pi()/2.); // rotate vector by -90 deg around z
1590// AliAODJet jetTmp = AliAODJet(vecPerp);
1591 if(fDebug>5) printf("TaskV0sInJets: Adding perp. cones number %d, %d\n",iNJetPerp,iNJetPerp+1);
1592// printf("TaskV0sInJets: Adding perp. cone number %d: pT %f, phi %f, eta %f, pT %f, phi %f, eta %f\n",iNJetSel,vecPerp.Pt(),vecPerp.Phi(),vecPerp.Eta(),jetTmp.Pt(),jetTmp.Phi(),jetTmp.Eta());
1593 new ((*jetArrayPerp)[iNJetPerp++]) AliAODJet(vecPerpPlus); // write perp. cone to the array
1594 new ((*jetArrayPerp)[iNJetPerp++]) AliAODJet(vecPerpMinus); // write perp. cone to the array
1595 if(fDebug>5) printf("TaskV0sInJets: Adding jet number %d\n",iNJetSel);
1596// printf("TaskV0sInJets: Adding jet number %d: pT %f, phi %f, eta %f\n",iNJetSel,jetSel->Pt(),jetSel->Phi(),jetSel->Eta());
1597 new ((*jetArraySel)[iNJetSel++]) AliAODJet(*((AliAODJet*)jetSel)); // copy selected jet to the array
1598 }
1599 if(fDebug>5) printf("TaskV0sInJets: Added jets: %d\n",iNJetSel);
1600 iNJetSel = jetArraySel->GetEntriesFast();
87e23d3e 1601 if(fDebug>2) printf("TaskV0sInJets: Selected jets in array: %d\n",iNJetSel);
7729cc34 1602 fh1NJetPerEvent[iCentIndex]->Fill(iNJetSel);
1603 // fill jet spectra
1604 for (Int_t iJet = 0; iJet<iNJetSel; iJet++)
1605 {
1606 jet = (AliAODJet*)jetArraySel->At(iJet); // load a jet in the list
1607 fh1PtJet[iCentIndex]->Fill(jet->Pt()); // pt spectrum of selected jets
1608 fh1EtaJet[iCentIndex]->Fill(jet->Eta()); // eta spectrum of selected jets
1609 fh2EtaPtJet[iCentIndex]->Fill(jet->Eta(),jet->Pt()); // eta-pT spectrum of selected jets
1610 fh1PhiJet[iCentIndex]->Fill(jet->Phi()); // phi spectrum of selected jets
1611 Double_t dAreaExcluded = TMath::Pi()*dRadiusExcludeCone*dRadiusExcludeCone; // area of the cone
87e23d3e 1612 dAreaExcluded -= AreaCircSegment(dRadiusExcludeCone,dEtaMax-jet->Eta()); // positive eta overhang
1613 dAreaExcluded -= AreaCircSegment(dRadiusExcludeCone,dEtaMax+jet->Eta()); // negative eta overhang
7729cc34 1614 fh1AreaExcluded->Fill(iCentIndex,dAreaExcluded);
1615 }
1616 jet = 0;
1617 }
1618
1619 if (bJetEventGood) // there should be some reconstructed jets
1620 {
1621 fh1EventCounterCut->Fill(4); // events with jet(s)
1622 fh1EventCounterCutCent[iCentIndex]->Fill(4); // events with jet(s)
1623 if (iNJetSel)
1624 {
1625 fh1EventCounterCut->Fill(5); // events with selected jets
1626 fh1EventCounterCutCent[iCentIndex]->Fill(5);
1627 }
1628 }
169867f8 1629 if (iNJetSel)
1630 fh1EventCent2Jets->Fill(fdCentrality);
1631 else
1632 fh1EventCent2NoJets->Fill(fdCentrality);
7729cc34 1633
1634 if (iNJetSel)
1635 {
87e23d3e 1636 jetRnd = GetRandomCone(jetArraySel,dJetEtaWindow,2*fdRadiusJet);
7729cc34 1637 if (jetRnd)
1638 {
1639 fh1NRndConeCent->Fill(iCentIndex);
1640 fh2EtaPhiRndCone[iCentIndex]->Fill(jetRnd->Eta(),jetRnd->Phi());
1641 }
9706bde3 1642 jetMed = GetMedianCluster(jetArrayBg,dJetEtaWindow);
1643 if (jetMed)
1644 {
1645 fh1NMedConeCent->Fill(iCentIndex);
1646 fh2EtaPhiMedCone[iCentIndex]->Fill(jetMed->Eta(),jetMed->Phi());
1647 }
7729cc34 1648 }
1649
1650 // Loading primary vertex info
1651 AliAODVertex* primVtx = fAODIn->GetPrimaryVertex(); // get the primary vertex
1652 Double_t dPrimVtxPos[3]; // primary vertex position {x,y,z}
1653 primVtx->GetXYZ(dPrimVtxPos);
1654 fh1VtxZ[iCentIndex]->Fill(dPrimVtxPos[2]);
1655 fh2VtxXY[iCentIndex]->Fill(dPrimVtxPos[0],dPrimVtxPos[1]);
1656
1657 /*===== Start of loop over V0 candidates =====*/
87e23d3e 1658 if(fDebug>2) printf("TaskV0sInJets: Start of V0 loop\n");
7729cc34 1659 for (Int_t iV0 = 0; iV0 < iNV0s; iV0++)
1660 {
1661 v0 = fAODIn->GetV0(iV0); // get next candidate from the list in AOD
1662 if (!v0)
1663 continue;
1664
1665 iNV0CandTot++;
1666
1667 // Initialization of status indicators
1668 Bool_t bIsCandidateK0s = kTRUE; // candidate for K0s
1669 Bool_t bIsCandidateLambda = kTRUE; // candidate for Lambda
1670 Bool_t bIsCandidateALambda = kTRUE; // candidate for Lambda
1671 Bool_t bIsInPeakK0s = kFALSE; // candidate within the K0s mass peak
1672 Bool_t bIsInPeakLambda = kFALSE; // candidate within the Lambda mass peak
1673 Bool_t bIsInConeJet = kFALSE; // candidate within the jet cones
1674 Bool_t bIsInConePerp = kFALSE; // candidate within the perpendicular cone
1675 Bool_t bIsInConeRnd = kFALSE; // candidate within the random cone
9706bde3 1676 Bool_t bIsInConeMed = kFALSE; // candidate within the median-cluster cone
7729cc34 1677 Bool_t bIsOutsideCones = kFALSE; // candidate outside excluded cones
1678
1679 // Invariant mass calculation
87e23d3e 1680 dMassV0K0s = v0->MassK0Short();
1681 dMassV0Lambda = v0->MassLambda();
1682 dMassV0ALambda = v0->MassAntiLambda();
7729cc34 1683
1684 Int_t iCutIndex = 0; // indicator of current selection step
1685 // 0
1686 // All V0 candidates
87e23d3e 1687 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, bIsCandidateK0s, bIsCandidateLambda, bIsCandidateALambda, iCutIndex, iCentIndex);
7729cc34 1688 iCutIndex++;
1689
1690 // Skip candidates outside the histogram range
87e23d3e 1691 if ( (dMassV0K0s < fgkdMassK0sMin) || (dMassV0K0s >= fgkdMassK0sMax) )
7729cc34 1692 bIsCandidateK0s = kFALSE;
87e23d3e 1693 if ( (dMassV0Lambda < fgkdMassLambdaMin) || (dMassV0Lambda >= fgkdMassLambdaMax) )
7729cc34 1694 bIsCandidateLambda = kFALSE;
87e23d3e 1695 if ( (dMassV0ALambda < fgkdMassLambdaMin) || (dMassV0ALambda >= fgkdMassLambdaMax) )
7729cc34 1696 bIsCandidateALambda = kFALSE;
1697 if (!bIsCandidateK0s && !bIsCandidateLambda && !bIsCandidateALambda)
1698 continue;
1699
87e23d3e 1700 Double_t dPtV0 = TMath::Sqrt(v0->Pt2V0()); // transverse momentum of V0
7729cc34 1701 vecV0Momentum = TVector3(v0->Px(),v0->Py(),v0->Pz()); // set the vector of V0 momentum
1702
1703 // Sigma of the mass peak window
87e23d3e 1704 Double_t dMassPeakWindowK0s = dNSigmaMassMax*MassPeakSigmaOld(dPtV0,0);
1705 Double_t dMassPeakWindowLambda = dNSigmaMassMax*MassPeakSigmaOld(dPtV0,1);
1706// Double_t dMassPeakWindowK0s = dNSigmaMassMax*MassPeakSigma(iCentIndex,dPtV0,0);
1707// Double_t dMassPeakWindowLambda = dNSigmaMassMax*MassPeakSigma(iCentIndex,dPtV0,1);
7729cc34 1708
1709 // Invariant mass peak selection
87e23d3e 1710 if (TMath::Abs(dMassV0K0s-dMassPDGK0s) < dMassPeakWindowK0s)
7729cc34 1711 bIsInPeakK0s = kTRUE;
87e23d3e 1712 if (TMath::Abs(dMassV0Lambda-dMassPDGLambda) < dMassPeakWindowLambda)
7729cc34 1713 bIsInPeakLambda = kTRUE;
1714
1715 // Retrieving all relevant properties of the V0 candidate
1716 Bool_t bOnFlyStatus = v0->GetOnFlyStatus(); // online (on fly) reconstructed vs offline reconstructed
1717 const AliAODTrack* trackPos = (AliAODTrack*)v0->GetDaughter(0); // positive daughter track
1718 const AliAODTrack* trackNeg = (AliAODTrack*)v0->GetDaughter(1); // negative daughter track
87e23d3e 1719 Double_t dPtDaughterPos = trackPos->Pt(); // transverse momentum of a daughter track
1720 Double_t dPtDaughterNeg = trackNeg->Pt();
1721 Double_t dNRowsPos = trackPos->GetTPCClusterInfo(2,1); // crossed TPC pad rows of a daughter track
1722 Double_t dNRowsNeg = trackNeg->GetTPCClusterInfo(2,1);
1723 Double_t dDCAToPrimVtxPos = TMath::Abs(v0->DcaPosToPrimVertex()); // dca of a daughter to the primary vertex
1724 Double_t dDCAToPrimVtxNeg = TMath::Abs(v0->DcaNegToPrimVertex());
1725 Double_t dDCADaughters = v0->DcaV0Daughters(); // dca between daughters
1726 Double_t dCPA = v0->CosPointingAngle(primVtx); // cosine of the pointing angle
7729cc34 1727 Double_t dSecVtxPos[3]; // V0 vertex position {x,y,z}
1728// Double_t dSecVtxPos[3] = {v0->DecayVertexV0X(),v0->DecayVertexV0Y(),v0->DecayVertexV0Z()}; // V0 vertex position
1729 v0->GetSecondaryVtx(dSecVtxPos);
87e23d3e 1730 Double_t dRadiusDecay = TMath::Sqrt(dSecVtxPos[0]*dSecVtxPos[0] + dSecVtxPos[1]*dSecVtxPos[1]); // distance of the V0 vertex from the z-axis
1731 Double_t dEtaDaughterNeg = trackNeg->Eta(); // = v0->EtaProng(1), pseudorapidity of a daughter track
1732 Double_t dEtaDaughterPos = trackPos->Eta(); // = v0->EtaProng(0)
1733 Double_t dRapK0s = v0->RapK0Short(); // rapidity calculated for K0s assumption
1734 Double_t dRapLambda = v0->RapLambda(); // rapidity calculated for Lambda assumption
1735 Double_t dEtaV0 = v0->Eta(); // V0 pseudorapidity
1736// Double_t dPhiV0 = v0->Phi(); // V0 pseudorapidity
7729cc34 1737 Double_t dDecayPath[3];
1738 for (Int_t iPos = 0; iPos < 3; iPos++)
1739 dDecayPath[iPos] = dSecVtxPos[iPos]-dPrimVtxPos[iPos]; // vector of the V0 path
87e23d3e 1740 Double_t dDecLen = TMath::Sqrt(dDecayPath[0]*dDecayPath[0]+dDecayPath[1]*dDecayPath[1]+dDecayPath[2]*dDecayPath[2]); // path length L
1741 Double_t dDecLen2D = TMath::Sqrt(dDecayPath[0]*dDecayPath[0]+dDecayPath[1]*dDecayPath[1]); // transverse path length R
1742 Double_t dLOverP = dDecLen/v0->P(); // L/p
1743 Double_t dROverPt = dDecLen2D/dPtV0; // R/pT
1744 Double_t dMLOverPK0s = dMassPDGK0s*dLOverP; // m*L/p = c*(proper lifetime)
1745// Double_t dMLOverPLambda = dMassPDGLambda*dLOverP; // m*L/p
1746 Double_t dMROverPtK0s = dMassPDGK0s*dROverPt; // m*R/pT
1747 Double_t dMROverPtLambda = dMassPDGLambda*dROverPt; // m*R/pT
1748 Double_t dNSigmaPosPion = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(trackPos,AliPID::kPion)); // difference between measured and expected signal of the dE/dx in the TPC
1749 Double_t dNSigmaPosProton = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(trackPos,AliPID::kProton));
1750 Double_t dNSigmaNegPion = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(trackNeg,AliPID::kPion));
1751 Double_t dNSigmaNegProton = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(trackNeg,AliPID::kProton));
1752 Double_t dAlpha = v0->AlphaV0(); // Armenteros-Podolanski alpha
1753 Double_t dPtArm = v0->PtArmV0(); // Armenteros-Podolanski pT
7729cc34 1754 AliAODVertex* prodVtxDaughterPos = (AliAODVertex*)(trackPos->GetProdVertex()); // production vertex of the positive daughter track
1755 Char_t cTypeVtxProdPos = prodVtxDaughterPos->GetType(); // type of the production vertex
1756 AliAODVertex* prodVtxDaughterNeg = (AliAODVertex*)(trackNeg->GetProdVertex()); // production vertex of the negative daughter track
1757 Char_t cTypeVtxProdNeg = prodVtxDaughterNeg->GetType(); // type of the production vertex
1758
87e23d3e 1759 fh2Tau3DVs2D[0]->Fill(dPtV0,dLOverP/dROverPt);
7729cc34 1760
1761 // QA histograms before cuts
1762 FillQAHistogramV0(primVtx,v0,0,bIsCandidateK0s,bIsCandidateLambda,bIsInPeakK0s,bIsInPeakLambda);
1763 // Cut vs mass histograms before cuts
1764 if (bIsCandidateK0s)
1765 {
87e23d3e 1766 fh2CutTPCRowsK0s[0]->Fill(dMassV0K0s,dNRowsPos);
1767 fh2CutTPCRowsK0s[0]->Fill(dMassV0K0s,dNRowsNeg);
1768 fh2CutPtPosK0s[0]->Fill(dMassV0K0s,dPtDaughterPos);
1769 fh2CutPtNegK0s[0]->Fill(dMassV0K0s,dPtDaughterNeg);
1770 fh2CutDCAVtx[0]->Fill(dMassV0K0s,dDCAToPrimVtxPos);
1771 fh2CutDCAVtx[0]->Fill(dMassV0K0s,dDCAToPrimVtxNeg);
1772 fh2CutDCAV0[0]->Fill(dMassV0K0s,dDCADaughters);
1773 fh2CutCos[0]->Fill(dMassV0K0s,dCPA);
1774 fh2CutR[0]->Fill(dMassV0K0s,dRadiusDecay);
1775 fh2CutEtaK0s[0]->Fill(dMassV0K0s,dEtaDaughterPos);
1776 fh2CutEtaK0s[0]->Fill(dMassV0K0s,dEtaDaughterNeg);
1777 fh2CutRapK0s[0]->Fill(dMassV0K0s,dRapK0s);
1778 fh2CutCTauK0s[0]->Fill(dMassV0K0s,dMROverPtK0s/dCTauK0s);
1779 fh2CutPIDPosK0s[0]->Fill(dMassV0K0s,dNSigmaPosPion);
1780 fh2CutPIDNegK0s[0]->Fill(dMassV0K0s,dNSigmaNegPion);
7729cc34 1781 }
1782 if (bIsCandidateLambda)
1783 {
87e23d3e 1784 fh2CutTPCRowsLambda[0]->Fill(dMassV0Lambda,dNRowsPos);
1785 fh2CutTPCRowsLambda[0]->Fill(dMassV0Lambda,dNRowsNeg);
1786 fh2CutPtPosLambda[0]->Fill(dMassV0Lambda,dPtDaughterPos);
1787 fh2CutPtNegLambda[0]->Fill(dMassV0Lambda,dPtDaughterNeg);
1788 fh2CutEtaLambda[0]->Fill(dMassV0Lambda,dEtaDaughterPos);
1789 fh2CutEtaLambda[0]->Fill(dMassV0Lambda,dEtaDaughterNeg);
1790 fh2CutRapLambda[0]->Fill(dMassV0Lambda,dRapLambda);
1791 fh2CutCTauLambda[0]->Fill(dMassV0Lambda,dMROverPtLambda/dCTauLambda);
1792 fh2CutPIDPosLambda[0]->Fill(dMassV0Lambda,dNSigmaPosProton);
1793 fh2CutPIDNegLambda[0]->Fill(dMassV0Lambda,dNSigmaNegPion);
7729cc34 1794 }
1795
1796 /*===== Start of reconstruction cutting =====*/
1797
1798 // 1
1799 // All V0 candidates
87e23d3e 1800 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, bIsCandidateK0s, bIsCandidateLambda, bIsCandidateALambda, iCutIndex, iCentIndex);
7729cc34 1801 iCutIndex++;
1802
1803 /* Start of global cuts */
1804 // 2
1805 // Reconstruction method
1806 if (bPrintCuts) printf("Rec: Applying cut: Reconstruction method: on-the-fly? %s\n",(bOnFly ? "yes" : "no"));
1807 if (bOnFlyStatus!=bOnFly)
1808 continue;
87e23d3e 1809 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, bIsCandidateK0s, bIsCandidateLambda, bIsCandidateALambda, iCutIndex, iCentIndex);
7729cc34 1810 iCutIndex++;
1811
1812 // 3
1813 // Tracks TPC OK
1814 if (bPrintCuts) printf("Rec: Applying cut: Correct charge of daughters\n");
1815 if ( !trackNeg || !trackPos )
1816 continue;
1817 if (trackNeg->Charge() == trackPos->Charge()) // daughters have different charge?
1818 continue;
1819 if (trackNeg->Charge() != -1) // daughters have expected charge?
1820 continue;
1821 if (trackPos->Charge() != 1) // daughters have expected charge?
1822 continue;
1823
1824 if (bPrintCuts) printf("Rec: Applying cut: TPC refit: %d\n",iRefit);
1825 if (!trackNeg->IsOn(iRefit)) // TPC refit is ON?
1826 continue;
1827 if (bPrintCuts) printf("Rec: Applying cut: Type of production vertex of daughter: Not %d\n",AliAODVertex::kKink);
1828 if(cTypeVtxProdNeg == AliAODVertex::kKink) // kink daughter rejection
1829 continue;
1830 // Old cuts Start
1831 if (bUseOldCuts)
1832 {
87e23d3e 1833 if (bPrintCuts) printf("Rec: Applying cut: Number of TPC rows: > %f\n",dNCrossedRowsTPCMin);
1834 if (dNRowsNeg < dNCrossedRowsTPCMin) // Crossed TPC padrows
7729cc34 1835 continue;
1836// Int_t findable = trackNeg->GetTPCNclsF(); // Findable clusters
1837// if (findable <= 0)
1838// continue;
87e23d3e 1839// if (dNRowsNeg/findable < dCrossedRowsOverFindMin)
7729cc34 1840// continue;
87e23d3e 1841// if (dNRowsNeg/findable > dCrossedRowsOverFindMax)
7729cc34 1842// continue;
1843 }
1844 // Old cuts End
1845
1846 if (!trackPos->IsOn(iRefit))
1847 continue;
1848 if(cTypeVtxProdPos == AliAODVertex::kKink) // kink daughter rejection
1849 continue;
1850 // Old cuts Start
1851 if (bUseOldCuts)
1852 {
87e23d3e 1853 if (dNRowsPos < dNCrossedRowsTPCMin)
7729cc34 1854 continue;
1855// findable = trackPos->GetTPCNclsF();
1856// if (findable <= 0)
1857// continue;
87e23d3e 1858// if (dNRowsPos/findable < dCrossedRowsOverFindMin)
7729cc34 1859// continue;
87e23d3e 1860// if (dNRowsPos/findable > dCrossedRowsOverFindMax)
7729cc34 1861// continue;
1862 }
1863 // Old cuts End
1864
87e23d3e 1865 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, bIsCandidateK0s, bIsCandidateLambda, bIsCandidateALambda, iCutIndex, iCentIndex);
7729cc34 1866 iCutIndex++;
1867
1868 // 4
1869 // Daughters: transverse momentum cut
1870 if (bUseOldCuts)
1871 {
87e23d3e 1872 if (bPrintCuts) printf("Rec: Applying cut: Daughter pT: > %f\n",dPtDaughterMin);
1873 if ( ( dPtDaughterNeg < dPtDaughterMin ) || ( dPtDaughterPos < dPtDaughterMin ) )
7729cc34 1874 continue;
87e23d3e 1875 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, bIsCandidateK0s, bIsCandidateLambda, bIsCandidateALambda, iCutIndex, iCentIndex);
7729cc34 1876 }
1877 iCutIndex++;
1878
1879 // 5
1880 // Daughters: Impact parameter of daughters to prim vtx
87e23d3e 1881 if (bPrintCuts) printf("Rec: Applying cut: Daughter DCA to prim vtx: > %f\n",dDCAToPrimVtxMin);
1882 if ( ( dDCAToPrimVtxNeg < dDCAToPrimVtxMin ) || ( dDCAToPrimVtxPos < dDCAToPrimVtxMin ) )
7729cc34 1883 continue;
87e23d3e 1884 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, bIsCandidateK0s, bIsCandidateLambda, bIsCandidateALambda, iCutIndex, iCentIndex);
7729cc34 1885 iCutIndex++;
1886
1887 // 6
1888 // Daughters: DCA
87e23d3e 1889 if (bPrintCuts) printf("Rec: Applying cut: DCA between daughters: < %f\n",dDCADaughtersMax);
1890 if (dDCADaughters > dDCADaughtersMax)
7729cc34 1891 continue;
87e23d3e 1892 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, bIsCandidateK0s, bIsCandidateLambda, bIsCandidateALambda, iCutIndex, iCentIndex);
7729cc34 1893 iCutIndex++;
1894
1895 // 7
1896 // V0: Cosine of the pointing angle
87e23d3e 1897 if (bPrintCuts) printf("Rec: Applying cut: CPA: > %f\n",dCPAMin);
1898 if (dCPA < dCPAMin)
7729cc34 1899 continue;
87e23d3e 1900 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, bIsCandidateK0s, bIsCandidateLambda, bIsCandidateALambda, iCutIndex, iCentIndex);
7729cc34 1901 iCutIndex++;
1902
1903 // 8
1904 // V0: Fiducial volume
87e23d3e 1905 if (bPrintCuts) printf("Rec: Applying cut: Decay radius: > %f, < %f\n",dRadiusDecayMin,dRadiusDecayMax);
1906 if ( (dRadiusDecay < dRadiusDecayMin) || (dRadiusDecay > dRadiusDecayMax) )
7729cc34 1907 continue;
87e23d3e 1908 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, bIsCandidateK0s, bIsCandidateLambda, bIsCandidateALambda, iCutIndex, iCentIndex);
7729cc34 1909 iCutIndex++;
1910
1911 // 9
1912 // Daughters: pseudorapidity cut
1913 if (bCutEtaDaughter)
1914 {
87e23d3e 1915 if (bPrintCuts) printf("Rec: Applying cut: Daughter |eta|: < %f\n",dEtaDaughterMax);
1916 if ( (TMath::Abs(dEtaDaughterNeg) > dEtaDaughterMax) || (TMath::Abs(dEtaDaughterPos) > dEtaDaughterMax) )
7729cc34 1917 continue;
87e23d3e 1918 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, bIsCandidateK0s, bIsCandidateLambda, bIsCandidateALambda, iCutIndex, iCentIndex);
7729cc34 1919 }
1920 iCutIndex++;
1921 /* End of global cuts */
1922
1923 /* Start of particle-dependent cuts */
1924 // 10
1925 // V0: rapidity cut & pseudorapidity cut
1926 if (bCutRapV0)
1927 {
87e23d3e 1928 if (bPrintCuts) printf("Rec: Applying cut: V0 |y|: < %f\n",dRapMax);
1929 if (TMath::Abs(dRapK0s) > dRapMax)
7729cc34 1930 bIsCandidateK0s = kFALSE;
87e23d3e 1931 if (TMath::Abs(dRapLambda) > dRapMax)
7729cc34 1932 {
1933 bIsCandidateLambda = kFALSE;
1934 bIsCandidateALambda = kFALSE;
1935 }
1936 }
1937 if (bCutEtaV0)
1938 {
87e23d3e 1939 if (bPrintCuts) printf("Rec: Applying cut: V0 |eta|: < %f\n",dEtaMax);
1940 if (TMath::Abs(dEtaV0) > dEtaMax)
7729cc34 1941 {
1942 bIsCandidateK0s = kFALSE;
1943 bIsCandidateLambda = kFALSE;
1944 bIsCandidateALambda = kFALSE;
1945 }
87e23d3e 1946 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, bIsCandidateK0s, bIsCandidateLambda, bIsCandidateALambda, iCutIndex, iCentIndex);
7729cc34 1947 }
1948 iCutIndex++;
1949
1950 // 11
1951 // Lifetime cut
1952 if (bCutTau)
1953 {
87e23d3e 1954 if (bPrintCuts) printf("Rec: Applying cut: Proper lifetime: < %f\n",dNTauMax);
1955 if (dMROverPtK0s > dNTauMax*dCTauK0s)
7729cc34 1956 bIsCandidateK0s = kFALSE;
87e23d3e 1957 if (dMROverPtLambda > dNTauMax*dCTauLambda)
7729cc34 1958 {
1959 bIsCandidateLambda = kFALSE;
1960 bIsCandidateALambda = kFALSE;
1961 }
87e23d3e 1962 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, bIsCandidateK0s, bIsCandidateLambda, bIsCandidateALambda, iCutIndex, iCentIndex);
7729cc34 1963 }
1964 iCutIndex++;
1965
1966 // 12
1967 // Daughter PID
1968 if (bCutPid)
1969 {
1970 if (bUseOldCuts)
1971 {
87e23d3e 1972 if (bPrintCuts) printf("Rec: Applying cut: Delta dE/dx (both daughters): < %f\n",dNSigmadEdxMax);
1973 if (dNSigmaPosPion > dNSigmadEdxMax || dNSigmaNegPion > dNSigmadEdxMax) // pi+, pi-
7729cc34 1974 bIsCandidateK0s = kFALSE;
87e23d3e 1975 if (dNSigmaPosProton > dNSigmadEdxMax || dNSigmaNegPion > dNSigmadEdxMax) // p+, pi-
7729cc34 1976 bIsCandidateLambda = kFALSE;
87e23d3e 1977 if (dNSigmaNegProton > dNSigmadEdxMax || dNSigmaPosPion > dNSigmadEdxMax) // p-, pi+
7729cc34 1978 bIsCandidateALambda = kFALSE;
1979 }
1980 else
1981 {
87e23d3e 1982 if (bPrintCuts) printf("Rec: Applying cut: Delta dE/dx (proton below %f GeV/c): < %f\n",dPtProtonPIDMax,dNSigmadEdxMax);
1983 if ( (dPtDaughterPos < dPtProtonPIDMax) && (dNSigmaPosProton > dNSigmadEdxMax) ) // p+
7729cc34 1984 bIsCandidateLambda = kFALSE;
87e23d3e 1985 if ( (dPtDaughterNeg < dPtProtonPIDMax) && (dNSigmaNegProton > dNSigmadEdxMax) ) // p-
7729cc34 1986 bIsCandidateALambda = kFALSE;
1987 }
87e23d3e 1988 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, bIsCandidateK0s, bIsCandidateLambda, bIsCandidateALambda, iCutIndex, iCentIndex);
7729cc34 1989 }
1990 iCutIndex++;
1991
87e23d3e 1992 Double_t valueCorrel[3] = {dMassV0K0s,dMassV0Lambda,dPtV0};
7729cc34 1993 if (bIsCandidateK0s && bIsCandidateLambda)
1994 fh3CCMassCorrelBoth->Fill(valueCorrel); // correlation of mass distribution of candidates selected as both K0s and Lambda
1995 if (bIsCandidateK0s && !bIsCandidateLambda)
1996 fh3CCMassCorrelKNotL->Fill(valueCorrel); // correlation of mass distribution of candidates selected as K0s and not Lambda
1997 if (!bIsCandidateK0s && bIsCandidateLambda)
1998 fh3CCMassCorrelLNotK->Fill(valueCorrel); // correlation of mass distribution of candidates selected as not K0s and Lambda
1999
2000 // 13
2001 // Armenteros-Podolanski cut
2002 if (bCutArmPod)
2003 {
2004 if (bPrintCuts) printf("Rec: Applying cut: Armenteros-Podolanski (K0S): pT > %f * |alpha|\n",0.2);
87e23d3e 2005 if(dPtArm < TMath::Abs(0.2*dAlpha))
7729cc34 2006 bIsCandidateK0s = kFALSE;
87e23d3e 2007 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, bIsCandidateK0s, bIsCandidateLambda, bIsCandidateALambda, iCutIndex, iCentIndex);
7729cc34 2008 }
2009 iCutIndex++;
2010
2011 // Cross contamination
2012 if (bIsInPeakK0s)
2013 {
2014 if (bIsCandidateLambda) // Lambda candidates in K0s peak, excluded from Lambda candidates by CC cut
87e23d3e 2015 fh2CCLambda->Fill(dMassV0Lambda,dPtV0);
7729cc34 2016 }
2017 if (bIsInPeakLambda)
2018 {
2019 if (bIsCandidateK0s) // K0s candidates in Lambda peak, excluded from K0s candidates by CC cut
87e23d3e 2020 fh2CCK0s->Fill(dMassV0K0s,dPtV0);
7729cc34 2021 }
2022// if (bCutCross)
2023// {
2024// if (bIsInPeakK0s)
2025// bIsCandidateLambda = kFALSE;
2026// if (bIsInPeakLambda)
2027// bIsCandidateK0s = kFALSE;
87e23d3e 2028// FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, bIsCandidateK0s, bIsCandidateLambda, bIsCandidateALambda, iCutIndex, iCentIndex);
7729cc34 2029// }
2030// iCutIndex++;
2031
2032 /* End of particle-dependent cuts */
2033
2034 /*===== End of reconstruction cutting =====*/
2035
2036 if (!bIsCandidateK0s && !bIsCandidateLambda && !bIsCandidateALambda)
2037 continue;
2038
9706bde3 2039 /*
2040 if(fDebug>5) printf("TaskV0sInJets: Adding selected V0 to branch\n");
2041 // Add selected candidates to the output tree branch
2042 if ((bIsCandidateK0s || bIsCandidateLambda || bIsCandidateALambda) && fbTreeOutput)
2043 {
2044 objectV0 = new ((*fBranchV0Rec)[iNV0SelV0Rec++]) AliV0Object(v0,primVtx);
2045 // new ((*fBranchV0Rec)[iNV0SelV0Rec++]) AliAODv0(*((AliAODv0*)v0));
2046 objectV0->SetIsCandidateK0S(bIsCandidateK0s);
2047 objectV0->SetIsCandidateLambda(bIsCandidateLambda);
2048 objectV0->SetIsCandidateALambda(bIsCandidateALambda);
2049 objectV0->SetNSigmaPosProton(dNSigmaPosProton);
2050 objectV0->SetNSigmaNegProton(dNSigmaNegProton);
2051 }
2052 */
7729cc34 2053
2054 // Selection of V0s in jet cones, perpendicular cones, random cones, outside cones
2055 if (bJetEventGood && iNJetSel && (bIsCandidateK0s || bIsCandidateLambda || bIsCandidateALambda))
2056 {
2057 // Selection of V0s in jet cones
2058 if(fDebug>5) printf("TaskV0sInJets: Searching for V0 %d %d in %d jet cones\n",bIsCandidateK0s,bIsCandidateLambda,iNJetSel);
2059 for (Int_t iJet = 0; iJet<iNJetSel; iJet++)
2060 {
2061 jet = (AliAODJet*)jetArraySel->At(iJet); // load a jet in the list
2062 vecJetMomentum = TVector3(jet->Px(),jet->Py(),jet->Pz()); // set the vector of jet momentum
2063 if(fDebug>5) printf("TaskV0sInJets: Checking if V0 %d %d in jet cone %d\n",bIsCandidateK0s,bIsCandidateLambda,iJet);
87e23d3e 2064 if (IsParticleInCone(v0,jet,fdRadiusJet)) // If good jet in event, find out whether V0 is in that jet
7729cc34 2065 {
2066 if(fDebug>5) printf("TaskV0sInJets: V0 %d %d found in jet cone %d\n",bIsCandidateK0s,bIsCandidateLambda,iJet);
2067 bIsInConeJet = kTRUE;
2068 break;
2069 }
2070 }
2071 // Selection of V0s in perp. cones
2072 if(fDebug>5) printf("TaskV0sInJets: Searching for V0 %d %d in %d perp. cones\n",bIsCandidateK0s,bIsCandidateLambda,iNJetSel);
2073 for (Int_t iJet = 0; iJet<iNJetPerp; iJet++)
2074 {
2075 jetPerp = (AliAODJet*)jetArrayPerp->At(iJet); // load a jet in the list
2076 if(fDebug>5) printf("TaskV0sInJets: Checking if V0 %d %d in perp. cone %d\n",bIsCandidateK0s,bIsCandidateLambda,iJet);
87e23d3e 2077 if (IsParticleInCone(v0,jetPerp,fdRadiusJet)) // V0 in perp. cone
7729cc34 2078 {
2079 if(fDebug>5) printf("TaskV0sInJets: V0 %d %d found in perp. cone %d\n",bIsCandidateK0s,bIsCandidateLambda,iJet);
2080 bIsInConePerp = kTRUE;
2081 break;
2082 }
2083 }
2084 // Selection of V0s in random cones
2085 if (jetRnd)
2086 {
2087 if(fDebug>5) printf("TaskV0sInJets: Searching for V0 %d %d in the rnd. cone\n",bIsCandidateK0s,bIsCandidateLambda);
87e23d3e 2088 if (IsParticleInCone(v0,jetRnd,fdRadiusJet)) // V0 in rnd. cone?
7729cc34 2089 {
2090 if(fDebug>5) printf("TaskV0sInJets: V0 %d %d found in the rnd. cone\n",bIsCandidateK0s,bIsCandidateLambda);
2091 bIsInConeRnd = kTRUE;
2092 }
2093 }
9706bde3 2094 // Selection of V0s in median-cluster cones
2095 if (jetMed)
2096 {
2097 if(fDebug>5) printf("TaskV0sInJets: Searching for V0 %d %d in the med. cone\n",bIsCandidateK0s,bIsCandidateLambda);
f8afd72f 2098 if (IsParticleInCone(v0,jetMed,fdRadiusJet)) // V0 in med. cone?
9706bde3 2099 {
2100 if(fDebug>5) printf("TaskV0sInJets: V0 %d %d found in the med. cone\n",bIsCandidateK0s,bIsCandidateLambda);
2101 bIsInConeMed = kTRUE;
2102 }
2103 }
7729cc34 2104 // Selection of V0s outside jet cones
2105 if(fDebug>5) printf("TaskV0sInJets: Searching for V0 %d %d outside jet cones\n",bIsCandidateK0s,bIsCandidateLambda);
2106 if (!OverlapWithJets(jetArraySel,v0,dRadiusExcludeCone)) // V0 oustide jet cones
2107 {
2108 if(fDebug>5) printf("TaskV0sInJets: V0 %d %d found outside jet cones\n",bIsCandidateK0s,bIsCandidateLambda);
2109 bIsOutsideCones = kTRUE;
2110 }
2111 }
2112
2113 // QA histograms after cuts
2114 FillQAHistogramV0(primVtx,v0,1,bIsCandidateK0s,bIsCandidateLambda,bIsInPeakK0s,bIsInPeakLambda);
2115 // Cut vs mass histograms after cuts
2116 if (bIsCandidateK0s)
2117 {
87e23d3e 2118 fh2CutTPCRowsK0s[1]->Fill(dMassV0K0s,dNRowsPos);
2119 fh2CutTPCRowsK0s[1]->Fill(dMassV0K0s,dNRowsNeg);
2120 fh2CutPtPosK0s[1]->Fill(dMassV0K0s,dPtDaughterPos);
2121 fh2CutPtNegK0s[1]->Fill(dMassV0K0s,dPtDaughterNeg);
2122 fh2CutDCAVtx[1]->Fill(dMassV0K0s,dDCAToPrimVtxPos);
2123 fh2CutDCAVtx[1]->Fill(dMassV0K0s,dDCAToPrimVtxNeg);
2124 fh2CutDCAV0[1]->Fill(dMassV0K0s,dDCADaughters);
2125 fh2CutCos[1]->Fill(dMassV0K0s,dCPA);
2126 fh2CutR[1]->Fill(dMassV0K0s,dRadiusDecay);
2127 fh2CutEtaK0s[1]->Fill(dMassV0K0s,dEtaDaughterPos);
2128 fh2CutEtaK0s[1]->Fill(dMassV0K0s,dEtaDaughterNeg);
2129 fh2CutRapK0s[1]->Fill(dMassV0K0s,dRapK0s);
2130 fh2CutCTauK0s[1]->Fill(dMassV0K0s,dMROverPtK0s/dCTauK0s);
2131 fh2CutPIDPosK0s[1]->Fill(dMassV0K0s,dNSigmaPosPion);
2132 fh2CutPIDNegK0s[1]->Fill(dMassV0K0s,dNSigmaNegPion);
7729cc34 2133 fh1DeltaZK0s[iCentIndex]->Fill(dDecayPath[2]);
2134 }
2135 if (bIsCandidateLambda)
2136 {
87e23d3e 2137 fh2CutTPCRowsLambda[1]->Fill(dMassV0Lambda,dNRowsPos);
2138 fh2CutTPCRowsLambda[1]->Fill(dMassV0Lambda,dNRowsNeg);
2139 fh2CutPtPosLambda[1]->Fill(dMassV0Lambda,dPtDaughterPos);
2140 fh2CutPtNegLambda[1]->Fill(dMassV0Lambda,dPtDaughterNeg);
2141 fh2CutEtaLambda[1]->Fill(dMassV0Lambda,dEtaDaughterPos);
2142 fh2CutEtaLambda[1]->Fill(dMassV0Lambda,dEtaDaughterNeg);
2143 fh2CutRapLambda[1]->Fill(dMassV0Lambda,dRapLambda);
2144 fh2CutCTauLambda[1]->Fill(dMassV0Lambda,dMROverPtLambda/dCTauLambda);
2145 fh2CutPIDPosLambda[1]->Fill(dMassV0Lambda,dNSigmaPosProton);
2146 fh2CutPIDNegLambda[1]->Fill(dMassV0Lambda,dNSigmaNegPion);
7729cc34 2147 fh1DeltaZLambda[iCentIndex]->Fill(dDecayPath[2]);
2148 }
2149
2150 /*===== Start of filling V0 spectra =====*/
2151
87e23d3e 2152 Double_t dAngle = TMath::Pi(); // angle between V0 momentum and jet momentum
7729cc34 2153 if (bIsInConeJet)
2154 {
87e23d3e 2155 dAngle = vecV0Momentum.Angle(vecJetMomentum);
7729cc34 2156 }
2157
2158 // iCutIndex = 14
2159 if (bIsCandidateK0s)
2160 {
2161 // 14 K0s candidates after cuts
87e23d3e 2162// printf("K0S: i = %d, m = %f, pT = %f, eta = %f, phi = %f\n",iV0,dMassV0K0s,dPtV0,dEtaV0,dPhiV0);
2163 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, bIsCandidateK0s, kFALSE, kFALSE, iCutIndex, iCentIndex);
2164 Double_t valueKIncl[3] = {dMassV0K0s,dPtV0,dEtaV0};
7729cc34 2165 fhnV0InclusiveK0s[iCentIndex]->Fill(valueKIncl);
87e23d3e 2166 fh1V0InvMassK0sCent[iCentIndex]->Fill(dMassV0K0s);
7729cc34 2167
87e23d3e 2168 fh1QACTau2D[1]->Fill(dMROverPtK0s/dCTauK0s);
2169 fh1QACTau3D[1]->Fill(dMLOverPK0s/dCTauK0s);
2170 fh2Tau3DVs2D[1]->Fill(dPtV0,dLOverP/dROverPt);
7729cc34 2171
2172 if (iNJetSel)
2173 {
2174 // 15 K0s in jet events
87e23d3e 2175 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, bIsCandidateK0s, kFALSE, kFALSE, iCutIndex+1, iCentIndex);
7729cc34 2176 }
2177 if (bIsInConeJet)
2178 {
2179 // 16 K0s in jets
87e23d3e 2180 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, bIsCandidateK0s, kFALSE, kFALSE, iCutIndex+2, iCentIndex);
2181 Double_t valueKInJC[4] = {dMassV0K0s,dPtV0,dEtaV0,jet->Pt()};
7729cc34 2182 fhnV0InJetK0s[iCentIndex]->Fill(valueKInJC);
87e23d3e 2183 fh2V0PtJetAngleK0s[iCentIndex]->Fill(jet->Pt(),dAngle);
7729cc34 2184 }
2185 if (bIsOutsideCones)
2186 {
87e23d3e 2187 Double_t valueKOutJC[3] = {dMassV0K0s,dPtV0,dEtaV0};
7729cc34 2188 fhnV0OutJetK0s[iCentIndex]->Fill(valueKOutJC);
2189 }
2190 if (bIsInConePerp)
2191 {
87e23d3e 2192 Double_t valueKInPC[4] = {dMassV0K0s,dPtV0,dEtaV0,jetPerp->Pt()};
7729cc34 2193 fhnV0InPerpK0s[iCentIndex]->Fill(valueKInPC);
2194 }
2195 if (bIsInConeRnd)
2196 {
87e23d3e 2197 Double_t valueKInRnd[3] = {dMassV0K0s,dPtV0,dEtaV0};
7729cc34 2198 fhnV0InRndK0s[iCentIndex]->Fill(valueKInRnd);
2199 }
9706bde3 2200 if (bIsInConeMed)
2201 {
2202 Double_t valueKInMed[3] = {dMassV0K0s,dPtV0,dEtaV0};
2203 fhnV0InMedK0s[iCentIndex]->Fill(valueKInMed);
2204 }
7729cc34 2205 if (!iNJetSel)
2206 {
87e23d3e 2207 Double_t valueKNoJet[3] = {dMassV0K0s,dPtV0,dEtaV0};
7729cc34 2208 fhnV0NoJetK0s[iCentIndex]->Fill(valueKNoJet);
2209 }
2210 iNV0CandK0s++;
2211 }
2212 if (bIsCandidateLambda)
2213 {
2214 // 14 Lambda candidates after cuts
87e23d3e 2215// printf("La: i = %d, m = %f, pT = %f, eta = %f, phi = %f\n",iV0,dMassV0Lambda,dPtV0,dEtaV0,dPhiV0);
2216 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, kFALSE, bIsCandidateLambda, kFALSE, iCutIndex, iCentIndex);
2217 Double_t valueLIncl[3] = {dMassV0Lambda,dPtV0,dEtaV0};
7729cc34 2218 fhnV0InclusiveLambda[iCentIndex]->Fill(valueLIncl);
87e23d3e 2219 fh1V0InvMassLambdaCent[iCentIndex]->Fill(dMassV0Lambda);
7729cc34 2220 if (iNJetSel)
2221 {
2222 // 15 Lambda in jet events
87e23d3e 2223 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, kFALSE, bIsCandidateLambda, kFALSE, iCutIndex+1, iCentIndex);
7729cc34 2224 }
2225 if (bIsInConeJet)
2226 {
2227 // 16 Lambda in jets
87e23d3e 2228 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, kFALSE, bIsCandidateLambda, kFALSE, iCutIndex+2, iCentIndex);
2229 Double_t valueLInJC[4] = {dMassV0Lambda,dPtV0,dEtaV0,jet->Pt()};
7729cc34 2230 fhnV0InJetLambda[iCentIndex]->Fill(valueLInJC);
87e23d3e 2231 fh2V0PtJetAngleLambda[iCentIndex]->Fill(jet->Pt(),dAngle);
7729cc34 2232 }
2233 if (bIsOutsideCones)
2234 {
87e23d3e 2235 Double_t valueLOutJet[3] = {dMassV0Lambda,dPtV0,dEtaV0};
7729cc34 2236 fhnV0OutJetLambda[iCentIndex]->Fill(valueLOutJet);
2237 }
2238 if (bIsInConePerp)
2239 {
87e23d3e 2240 Double_t valueLInPC[4] = {dMassV0Lambda,dPtV0,dEtaV0,jetPerp->Pt()};
7729cc34 2241 fhnV0InPerpLambda[iCentIndex]->Fill(valueLInPC);
2242 }
2243 if (bIsInConeRnd)
2244 {
87e23d3e 2245 Double_t valueLInRnd[3] = {dMassV0Lambda,dPtV0,dEtaV0};
7729cc34 2246 fhnV0InRndLambda[iCentIndex]->Fill(valueLInRnd);
2247 }
9706bde3 2248 if (bIsInConeMed)
2249 {
2250 Double_t valueLInMed[3] = {dMassV0Lambda,dPtV0,dEtaV0};
2251 fhnV0InMedLambda[iCentIndex]->Fill(valueLInMed);
2252 }
7729cc34 2253 if (!iNJetSel)
2254 {
87e23d3e 2255 Double_t valueLNoJet[3] = {dMassV0Lambda,dPtV0,dEtaV0};
7729cc34 2256 fhnV0NoJetLambda[iCentIndex]->Fill(valueLNoJet);
2257 }
2258 iNV0CandLambda++;
2259 }
2260 if (bIsCandidateALambda)
2261 {
2262 // 14 ALambda candidates after cuts
87e23d3e 2263// printf("AL: i = %d, m = %f, pT = %f, eta = %f, phi = %f\n",iV0,dMassV0ALambda,dPtV0,dEtaV0,dPhiV0);
2264 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, kFALSE, kFALSE, bIsCandidateALambda, iCutIndex, iCentIndex);
2265 Double_t valueALIncl[3] = {dMassV0ALambda,dPtV0,dEtaV0};
7729cc34 2266 fhnV0InclusiveALambda[iCentIndex]->Fill(valueALIncl);
87e23d3e 2267 fh1V0InvMassALambdaCent[iCentIndex]->Fill(dMassV0ALambda);
7729cc34 2268 if (iNJetSel)
2269 {
2270 // 15 ALambda in jet events
87e23d3e 2271 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, kFALSE, kFALSE, bIsCandidateALambda, iCutIndex+1, iCentIndex);
7729cc34 2272 }
2273 if (bIsInConeJet)
2274 {
2275 // 16 ALambda in jets
87e23d3e 2276 FillCandidates(dMassV0K0s, dMassV0Lambda, dMassV0ALambda, kFALSE, kFALSE, bIsCandidateALambda, iCutIndex+2, iCentIndex);
2277 Double_t valueLInJC[4] = {dMassV0ALambda,dPtV0,dEtaV0,jet->Pt()};
7729cc34 2278 fhnV0InJetALambda[iCentIndex]->Fill(valueLInJC);
87e23d3e 2279 fh2V0PtJetAngleALambda[iCentIndex]->Fill(jet->Pt(),dAngle);
7729cc34 2280 }
2281 if (bIsOutsideCones)
2282 {
87e23d3e 2283 Double_t valueALOutJet[3] = {dMassV0ALambda,dPtV0,dEtaV0};
7729cc34 2284 fhnV0OutJetALambda[iCentIndex]->Fill(valueALOutJet);
2285 }
2286 if (bIsInConePerp)
2287 {
87e23d3e 2288 Double_t valueLInPC[4] = {dMassV0ALambda,dPtV0,dEtaV0,jetPerp->Pt()};
7729cc34 2289 fhnV0InPerpALambda[iCentIndex]->Fill(valueLInPC);
2290 }
2291 if (bIsInConeRnd)
2292 {
87e23d3e 2293 Double_t valueALInRnd[3] = {dMassV0ALambda,dPtV0,dEtaV0};
7729cc34 2294 fhnV0InRndALambda[iCentIndex]->Fill(valueALInRnd);
2295 }
9706bde3 2296 if (bIsInConeMed)
2297 {
2298 Double_t valueALInMed[3] = {dMassV0ALambda,dPtV0,dEtaV0};
2299 fhnV0InMedALambda[iCentIndex]->Fill(valueALInMed);
2300 }
7729cc34 2301 if (!iNJetSel)
2302 {
87e23d3e 2303 Double_t valueALNoJet[3] = {dMassV0ALambda,dPtV0,dEtaV0};
7729cc34 2304 fhnV0NoJetALambda[iCentIndex]->Fill(valueALNoJet);
2305 }
2306 iNV0CandALambda++;
2307 }
2308 /*===== End of filling V0 spectra =====*/
2309
2310
2311 /*===== Association of reconstructed V0 candidates with MC particles =====*/
2312 if (fbMCAnalysis)
2313 {
2314 // Associate selected candidates only
2315// if ( !(bIsCandidateK0s && bIsInPeakK0s) && !(bIsCandidateLambda && bIsInPeakLambda) ) // signal candidates
2316 if ( !(bIsCandidateK0s) && !(bIsCandidateLambda) && !(bIsCandidateALambda) ) // chosen candidates with any mass
2317 continue;
2318
2319 // Get MC labels of reconstructed daughter tracks
2320 Int_t iLabelPos = TMath::Abs(trackPos->GetLabel());
2321 Int_t iLabelNeg = TMath::Abs(trackNeg->GetLabel());
2322
2323 // Make sure MC daughters are in the array range
2324 if ( (iLabelNeg<0) || (iLabelNeg>=iNTracksMC) || (iLabelPos<0) || (iLabelPos>=iNTracksMC) )
2325 continue;
2326
2327 // Get MC particles corresponding to reconstructed daughter tracks
2328 AliAODMCParticle* particleMCDaughterNeg = (AliAODMCParticle*)arrayMC->At(iLabelNeg);
2329 AliAODMCParticle* particleMCDaughterPos = (AliAODMCParticle*)arrayMC->At(iLabelPos);
2330 if (!particleMCDaughterNeg || !particleMCDaughterPos)
2331 continue;
2332
2333 // Make sure MC daughter particles are not physical primary
2334 if ( (particleMCDaughterNeg->IsPhysicalPrimary()) || (particleMCDaughterPos->IsPhysicalPrimary()) )
2335 continue;
2336
2337 // Get identities of MC daughter particles
2338 Int_t iPdgCodeDaughterPos = particleMCDaughterPos->GetPdgCode();
2339 Int_t iPdgCodeDaughterNeg = particleMCDaughterNeg->GetPdgCode();
2340
2341 // Get index of the mother particle for each MC daughter particle
2342 Int_t iIndexMotherPos = particleMCDaughterPos->GetMother();
2343 Int_t iIndexMotherNeg = particleMCDaughterNeg->GetMother();
2344
2345 if ( (iIndexMotherNeg<0) || (iIndexMotherNeg>=iNTracksMC) || (iIndexMotherPos<0) || (iIndexMotherPos>=iNTracksMC) )
2346 continue;
2347
2348 // Check whether MC daughter particles have the same mother
2349 if (iIndexMotherNeg != iIndexMotherPos)
2350 continue;
2351
2352 // Get the MC mother particle of both MC daughter particles
2353 AliAODMCParticle* particleMCMother = (AliAODMCParticle*)arrayMC->At(iIndexMotherPos);
2354 if (!particleMCMother)
2355 continue;
2356
2357 // Get identity of the MC mother particle
2358 Int_t iPdgCodeMother = particleMCMother->GetPdgCode();
2359
2360 // Skip not interesting particles
2361 if ( (iPdgCodeMother != iPdgCodeK0s) && (TMath::Abs(iPdgCodeMother) != iPdgCodeLambda) )
2362 continue;
2363
2364 // Check identity of the MC mother particle and the decay channel
2365 // Is MC mother particle K0S?
2366 Bool_t bV0MCIsK0s = ( (iPdgCodeMother==iPdgCodeK0s) && (iPdgCodeDaughterPos==+iPdgCodePion) && (iPdgCodeDaughterNeg==-iPdgCodePion) );
2367 // Is MC mother particle Lambda?
2368 Bool_t bV0MCIsLambda = ( (iPdgCodeMother==+iPdgCodeLambda) && (iPdgCodeDaughterPos==+iPdgCodeProton) && (iPdgCodeDaughterNeg==-iPdgCodePion) );
87e23d3e 2369 // Is MC mother particle anti-Lambda?
7729cc34 2370 Bool_t bV0MCIsALambda = ( (iPdgCodeMother==-iPdgCodeLambda) && (iPdgCodeDaughterPos==+iPdgCodePion) && (iPdgCodeDaughterNeg==-iPdgCodeProton) );
2371
87e23d3e 2372 Double_t dPtV0Gen = particleMCMother->Pt();
2373// Double_t dRapV0MC = particleMCMother->Y();
2374 Double_t dEtaV0Gen = particleMCMother->Eta();
2375// Double_t dPhiV0Gen = particleMCMother->Phi();
7729cc34 2376
2377 // Is MC mother particle physical primary? Attention!! Definition of IsPhysicalPrimary may change!!
2378// Bool_t bV0MCIsPrimary = particleMCMother->IsPhysicalPrimary();
2379 // Get the MC mother particle of the MC mother particle
2380 Int_t iIndexMotherOfMother = particleMCMother->GetMother();
2381 AliAODMCParticle* particleMCMotherOfMother = 0;
2382 if (iIndexMotherOfMother >= 0)
2383 particleMCMotherOfMother = (AliAODMCParticle*)arrayMC->At(iIndexMotherOfMother);
2384 // Get identity of the MC mother particle of the MC mother particle if it exists
2385 Int_t iPdgCodeMotherOfMother = 0;
2386 if (particleMCMotherOfMother)
2387 iPdgCodeMotherOfMother = particleMCMotherOfMother->GetPdgCode();
2388 // Check if the MC mother particle of the MC mother particle is a physical primary Sigma (3212 - Sigma0, 3224 - Sigma*+, 3214 - Sigma*0, 3114 - Sigma*-)
2389// Bool_t bV0MCComesFromSigma = kFALSE; // Is MC mother particle daughter of a Sigma?
2390// if ( (particleMCMotherOfMother && particleMCMotherOfMother->IsPhysicalPrimary()) && ( (TMath::Abs(iPdgCodeMotherOfMother)==3212) || (TMath::Abs(iPdgCodeMotherOfMother)==3224) || (TMath::Abs(iPdgCodeMotherOfMother)==3214) || (TMath::Abs(iPdgCodeMotherOfMother)==3114) ) )
2391// bV0MCComesFromSigma = kTRUE;
2392 // Should MC mother particle be considered as primary when it is Lambda?
2393// Bool_t bV0MCIsPrimaryLambda = (bV0MCIsPrimary || bV0MCComesFromSigma);
2394 // Check if the MC mother particle of the MC mother particle is a Xi (3322 - Xi0, 3312 - Xi-)
2395 Bool_t bV0MCComesFromXi = ( (particleMCMotherOfMother) && ( (iPdgCodeMotherOfMother==3322) || (iPdgCodeMotherOfMother==3312) ) ); // Is MC mother particle daughter of a Xi?
2396 Bool_t bV0MCComesFromAXi = ( (particleMCMotherOfMother) && ( (iPdgCodeMotherOfMother==-3322) || (iPdgCodeMotherOfMother==-3312) ) ); // Is MC mother particle daughter of a anti-Xi?
2397
2398 // Get the distance between production point of the MC mother particle and the primary vertex
2399 Double_t dx = dPrimVtxMCX-particleMCMother->Xv();
2400 Double_t dy = dPrimVtxMCY-particleMCMother->Yv();
2401 Double_t dz = dPrimVtxMCZ-particleMCMother->Zv();
87e23d3e 2402 Double_t dDistPrimary = TMath::Sqrt(dx*dx + dy*dy + dz*dz);
2403 Bool_t bV0MCIsPrimaryDist = (dDistPrimary < dDistPrimaryMax); // Is close enough to be considered primary-like?
7729cc34 2404
2405 // phi, eta resolution for V0-reconstruction
87e23d3e 2406// Double_t dResolutionV0Eta = particleMCMother->Eta()-v0->Eta();
2407// Double_t dResolutionV0Phi = particleMCMother->Phi()-v0->Phi();
7729cc34 2408
9706bde3 2409 /*
2410 if (fbTreeOutput)
2411 {
2412 objectV0->SetPtTrue(dPtV0Gen);
2413 objectV0->SetEtaTrue(dEtaV0Gen);
2414 objectV0->SetPhiTrue(dPhiV0Gen);
2415 objectV0->SetPDGCode(iPdgCodeMother);
2416 objectV0->SetPDGCodeMother(iPdgCodeMotherOfMother);
2417 }
2418 */
7729cc34 2419
2420 // K0s
2421// if (bIsCandidateK0s && bIsInPeakK0s) // selected candidates in peak
2422 if (bIsCandidateK0s) // selected candidates with any mass
2423 {
2424// if (bV0MCIsK0s && bV0MCIsPrimary) // well reconstructed candidates
2425 if (bV0MCIsK0s && bV0MCIsPrimaryDist) // well reconstructed candidates
2426 {
2427// if (fbTreeOutput)
2428// objectV0->SetOrigin(1);
87e23d3e 2429 fh2V0K0sPtMassMCRec[iCentIndex]->Fill(dPtV0Gen,dMassV0K0s);
2430 Double_t valueEtaK[3] = {dMassV0K0s,dPtV0Gen,dEtaV0Gen};
7729cc34 2431 fh3V0K0sEtaPtMassMCRec[iCentIndex]->Fill(valueEtaK);
2432
87e23d3e 2433 Double_t valueEtaDKNeg[6] = {0,particleMCDaughterNeg->Eta(),particleMCDaughterNeg->Pt(),dEtaV0Gen,dPtV0Gen,0};
7729cc34 2434 fhnV0K0sInclDaughterEtaPtPtMCRec[iCentIndex]->Fill(valueEtaDKNeg);
87e23d3e 2435 Double_t valueEtaDKPos[6] = {1,particleMCDaughterPos->Eta(),particleMCDaughterPos->Pt(),dEtaV0Gen,dPtV0Gen,0};
7729cc34 2436 fhnV0K0sInclDaughterEtaPtPtMCRec[iCentIndex]->Fill(valueEtaDKPos);
2437
87e23d3e 2438 fh2V0K0sMCResolMPt[iCentIndex]->Fill(dMassV0K0s-dMassPDGK0s,dPtV0);
2439 fh2V0K0sMCPtGenPtRec[iCentIndex]->Fill(dPtV0Gen,dPtV0);
7729cc34 2440 if (bIsInConeJet) // true V0 associated to a candidate in jet
2441 {
87e23d3e 2442 Double_t valueKInJCMC[4] = {dMassV0K0s,dPtV0Gen,dEtaV0Gen,jet->Pt()};
7729cc34 2443 fh3V0K0sInJetPtMassMCRec[iCentIndex]->Fill(valueKInJCMC);
87e23d3e 2444 Double_t valueEtaKIn[5] = {dMassV0K0s,dPtV0Gen,dEtaV0Gen,jet->Pt(),dEtaV0Gen-jet->Eta()};
7729cc34 2445 fh4V0K0sInJetEtaPtMassMCRec[iCentIndex]->Fill(valueEtaKIn);
2446
87e23d3e 2447 Double_t valueEtaDKJCNeg[6] = {0,particleMCDaughterNeg->Eta(),particleMCDaughterNeg->Pt(),dEtaV0Gen,dPtV0Gen,jet->Pt()};
7729cc34 2448 fhnV0K0sInJetsDaughterEtaPtPtMCRec[iCentIndex]->Fill(valueEtaDKJCNeg);
87e23d3e 2449 Double_t valueEtaDKJCPos[6] = {1,particleMCDaughterPos->Eta(),particleMCDaughterPos->Pt(),dEtaV0Gen,dPtV0Gen,jet->Pt()};
7729cc34 2450 fhnV0K0sInJetsDaughterEtaPtPtMCRec[iCentIndex]->Fill(valueEtaDKJCPos);
2451 }
2452 }
2453 if (bV0MCIsK0s && !bV0MCIsPrimaryDist) // not primary K0s
2454 {
2455// if (fbTreeOutput)
2456// objectV0->SetOrigin(-1);
87e23d3e 2457 fh1V0K0sPtMCRecFalse[iCentIndex]->Fill(dPtV0Gen);
7729cc34 2458 }
2459 }
2460 // Lambda
2461// if (bIsCandidateLambda && bIsInPeakLambda) // selected candidates in peak
2462 if (bIsCandidateLambda) // selected candidates with any mass
2463 {
2464// if (bV0MCIsLambda && bV0MCIsPrimaryLambda) // well reconstructed candidates
2465 if (bV0MCIsLambda && bV0MCIsPrimaryDist) // well reconstructed candidates
2466 {
2467// if (fbTreeOutput)
2468// objectV0->SetOrigin(1);
87e23d3e 2469 fh2V0LambdaPtMassMCRec[iCentIndex]->Fill(dPtV0Gen,dMassV0Lambda);
2470 Double_t valueEtaL[3] = {dMassV0Lambda,dPtV0Gen,dEtaV0Gen};
7729cc34 2471 fh3V0LambdaEtaPtMassMCRec[iCentIndex]->Fill(valueEtaL);
2472
87e23d3e 2473 Double_t valueEtaDLNeg[6] = {0,particleMCDaughterNeg->Eta(),particleMCDaughterNeg->Pt(),dEtaV0Gen,dPtV0Gen,0};
7729cc34 2474 fhnV0LambdaInclDaughterEtaPtPtMCRec[iCentIndex]->Fill(valueEtaDLNeg);
87e23d3e 2475 Double_t valueEtaDLPos[6] = {1,particleMCDaughterPos->Eta(),particleMCDaughterPos->Pt(),dEtaV0Gen,dPtV0Gen,0};
7729cc34 2476 fhnV0LambdaInclDaughterEtaPtPtMCRec[iCentIndex]->Fill(valueEtaDLPos);
2477
87e23d3e 2478 fh2V0LambdaMCResolMPt[iCentIndex]->Fill(dMassV0Lambda-dMassPDGLambda,dPtV0);
2479 fh2V0LambdaMCPtGenPtRec[iCentIndex]->Fill(dPtV0Gen,dPtV0);
7729cc34 2480 if (bIsInConeJet) // true V0 associated to a reconstructed candidate in jet
2481 {
87e23d3e 2482 Double_t valueLInJCMC[4] = {dMassV0Lambda,dPtV0Gen,dEtaV0Gen,jet->Pt()};
7729cc34 2483 fh3V0LambdaInJetPtMassMCRec[iCentIndex]->Fill(valueLInJCMC);
87e23d3e 2484 Double_t valueEtaLIn[5] = {dMassV0Lambda,dPtV0Gen,dEtaV0Gen,jet->Pt(),dEtaV0Gen-jet->Eta()};
7729cc34 2485 fh4V0LambdaInJetEtaPtMassMCRec[iCentIndex]->Fill(valueEtaLIn);
2486
87e23d3e 2487 Double_t valueEtaDLJCNeg[6] = {0,particleMCDaughterNeg->Eta(),particleMCDaughterNeg->Pt(),dEtaV0Gen,dPtV0Gen,jet->Pt()};
7729cc34 2488 fhnV0LambdaInJetsDaughterEtaPtPtMCRec[iCentIndex]->Fill(valueEtaDLJCNeg);
87e23d3e 2489 Double_t valueEtaDLJCPos[6] = {1,particleMCDaughterPos->Eta(),particleMCDaughterPos->Pt(),dEtaV0Gen,dPtV0Gen,jet->Pt()};
7729cc34 2490 fhnV0LambdaInJetsDaughterEtaPtPtMCRec[iCentIndex]->Fill(valueEtaDLJCPos);
2491 }
2492 }
2493 // Fill the feed-down histograms
2494 if (bV0MCIsLambda && bV0MCComesFromXi)
2495 {
2496// if (fbTreeOutput)
2497// objectV0->SetOrigin(2);
87e23d3e 2498 Double_t valueFDLIncl[3] = {dPtV0Gen,particleMCMotherOfMother->Pt(),0.};
7729cc34 2499 fhnV0LambdaInclMCFD[iCentIndex]->Fill(valueFDLIncl);
2500 if (bIsInConeRnd)
2501 {
2502 fhnV0LambdaBulkMCFD[iCentIndex]->Fill(valueFDLIncl);
2503 }
2504 if (bIsInConeJet)
2505 {
87e23d3e 2506 Double_t valueFDLInJets[3] = {dPtV0Gen,particleMCMotherOfMother->Pt(),jet->Pt()};
7729cc34 2507 fhnV0LambdaInJetsMCFD[iCentIndex]->Fill(valueFDLInJets);
2508 }
2509 }
2510 if (bV0MCIsLambda && !bV0MCIsPrimaryDist && !bV0MCComesFromXi) // not primary Lambda
2511 {
2512// if (fbTreeOutput)
2513// objectV0->SetOrigin(-1);
87e23d3e 2514 fh1V0LambdaPtMCRecFalse[iCentIndex]->Fill(dPtV0Gen);
7729cc34 2515 }
2516 }
2517 // anti-Lambda
2518// if (bIsCandidateALambda && bIsInPeakALambda) // selected candidates in peak
2519 if (bIsCandidateALambda) // selected candidates with any mass
2520 {
2521// if (bV0MCIsALambda && bV0MCIsPrimaryALambda) // well reconstructed candidates
2522 if (bV0MCIsALambda && bV0MCIsPrimaryDist) // well reconstructed candidates
2523 {
2524// if (fbTreeOutput)
2525// objectV0->SetOrigin(1);
87e23d3e 2526 fh2V0ALambdaPtMassMCRec[iCentIndex]->Fill(dPtV0Gen,dMassV0ALambda);
2527 Double_t valueEtaAL[3] = {dMassV0ALambda,dPtV0Gen,dEtaV0Gen};
7729cc34 2528 fh3V0ALambdaEtaPtMassMCRec[iCentIndex]->Fill(valueEtaAL);
2529
87e23d3e 2530 Double_t valueEtaDALNeg[6] = {0,particleMCDaughterNeg->Eta(),particleMCDaughterNeg->Pt(),dEtaV0Gen,dPtV0Gen,0};
7729cc34 2531 fhnV0ALambdaInclDaughterEtaPtPtMCRec[iCentIndex]->Fill(valueEtaDALNeg);
87e23d3e 2532 Double_t valueEtaDALPos[6] = {1,particleMCDaughterPos->Eta(),particleMCDaughterPos->Pt(),dEtaV0Gen,dPtV0Gen,0};
7729cc34 2533 fhnV0ALambdaInclDaughterEtaPtPtMCRec[iCentIndex]->Fill(valueEtaDALPos);
2534
87e23d3e 2535 fh2V0ALambdaMCResolMPt[iCentIndex]->Fill(dMassV0ALambda-dMassPDGLambda,dPtV0);
2536 fh2V0ALambdaMCPtGenPtRec[iCentIndex]->Fill(dPtV0Gen,dPtV0);
7729cc34 2537 if (bIsInConeJet) // true V0 associated to a reconstructed candidate in jet
2538 {
87e23d3e 2539 Double_t valueALInJCMC[4] = {dMassV0ALambda,dPtV0Gen,dEtaV0Gen,jet->Pt()};
7729cc34 2540 fh3V0ALambdaInJetPtMassMCRec[iCentIndex]->Fill(valueALInJCMC);
87e23d3e 2541 Double_t valueEtaALIn[5] = {dMassV0ALambda,dPtV0Gen,dEtaV0Gen,jet->Pt(),dEtaV0Gen-jet->Eta()};
7729cc34 2542 fh4V0ALambdaInJetEtaPtMassMCRec[iCentIndex]->Fill(valueEtaALIn);
2543
87e23d3e 2544 Double_t valueEtaDALJCNeg[6] = {0,particleMCDaughterNeg->Eta(),particleMCDaughterNeg->Pt(),dEtaV0Gen,dPtV0Gen,jet->Pt()};
7729cc34 2545 fhnV0ALambdaInJetsDaughterEtaPtPtMCRec[iCentIndex]->Fill(valueEtaDALJCNeg);
87e23d3e 2546 Double_t valueEtaDALJCPos[6] = {1,particleMCDaughterPos->Eta(),particleMCDaughterPos->Pt(),dEtaV0Gen,dPtV0Gen,jet->Pt()};
7729cc34 2547 fhnV0ALambdaInJetsDaughterEtaPtPtMCRec[iCentIndex]->Fill(valueEtaDALJCPos);
2548 }
2549 }
2550 // Fill the feed-down histograms
2551 if (bV0MCIsALambda && bV0MCComesFromAXi)
2552 {
2553// if (fbTreeOutput)
2554// objectV0->SetOrigin(2);
87e23d3e 2555 Double_t valueFDALIncl[3] = {dPtV0Gen,particleMCMotherOfMother->Pt(),0.};
7729cc34 2556 fhnV0ALambdaInclMCFD[iCentIndex]->Fill(valueFDALIncl);
2557 if (bIsInConeRnd)
2558 {
2559 fhnV0ALambdaBulkMCFD[iCentIndex]->Fill(valueFDALIncl);
2560 }
2561 if (bIsInConeJet)
2562 {
87e23d3e 2563 Double_t valueFDALInJets[3] = {dPtV0Gen,particleMCMotherOfMother->Pt(),jet->Pt()};
7729cc34 2564 fhnV0ALambdaInJetsMCFD[iCentIndex]->Fill(valueFDALInJets);
2565 }
2566 }
87e23d3e 2567 if (bV0MCIsALambda && !bV0MCIsPrimaryDist && !bV0MCComesFromAXi) // not primary anti-Lambda
7729cc34 2568 {
2569// if (fbTreeOutput)
2570// objectV0->SetOrigin(-1);
87e23d3e 2571 fh1V0ALambdaPtMCRecFalse[iCentIndex]->Fill(dPtV0Gen);
7729cc34 2572 }
2573 }
2574 }
2575 /*===== End Association of reconstructed V0 candidates with MC particles =====*/
2576 }
2577 /*===== End of V0 loop =====*/
2578
2579 fh1V0CandPerEvent->Fill(iNV0CandTot);
2580 fh1V0CandPerEventCentK0s[iCentIndex]->Fill(iNV0CandK0s);
2581 fh1V0CandPerEventCentLambda[iCentIndex]->Fill(iNV0CandLambda);
2582 fh1V0CandPerEventCentALambda[iCentIndex]->Fill(iNV0CandALambda);
2583
87e23d3e 2584 if(fDebug>2) printf("TaskV0sInJets: End of V0 loop\n");
7729cc34 2585
2586 // Spectra of generated particles
2587 if (fbMCAnalysis)
2588 {
2589 for (Int_t iPartMC = 0; iPartMC < iNTracksMC; iPartMC++)
2590 {
2591 // Get MC particle
2592 AliAODMCParticle* particleMC = (AliAODMCParticle*)arrayMC->At(iPartMC);
2593 if(!particleMC)
2594 continue;
2595
2596 // Get identity of MC particle
2597 Int_t iPdgCodeParticleMC = particleMC->GetPdgCode();
2598 // Fill Xi spectrum (3322 - Xi0, 3312 - Xi-)
2599// if ( (iPdgCodeParticleMC==3322) || (iPdgCodeParticleMC==3312) )
2600 if ( (iPdgCodeParticleMC==3312) && (TMath::Abs(particleMC->Y())<0.5) )
2601 {
2602// if (fbTreeOutput)
2603// new ((*fBranchV0Gen)[iNV0SelV0Gen++]) AliAODMCParticle(*((AliAODMCParticle*)particleMC));
2604 fh1V0XiPtMCGen[iCentIndex]->Fill(particleMC->Pt());
2605 }
2606 if ( (iPdgCodeParticleMC==-3312) && (TMath::Abs(particleMC->Y())<0.5) )
2607 {
2608// if (fbTreeOutput)
2609// new ((*fBranchV0Gen)[iNV0SelV0Gen++]) AliAODMCParticle(*((AliAODMCParticle*)particleMC));
2610 fh1V0AXiPtMCGen[iCentIndex]->Fill(particleMC->Pt());
2611 }
2612 // Skip not interesting particles
2613 if ( (iPdgCodeParticleMC != iPdgCodeK0s) && (TMath::Abs(iPdgCodeParticleMC) != iPdgCodeLambda) )
2614 continue;
2615
2616 // Check identity of the MC V0 particle
2617 // Is MC V0 particle K0S?
2618 Bool_t bV0MCIsK0s = (iPdgCodeParticleMC==iPdgCodeK0s);
2619 // Is MC V0 particle Lambda?
2620 Bool_t bV0MCIsLambda = (iPdgCodeParticleMC==+iPdgCodeLambda);
2621 // Is MC V0 particle anti-Lambda?
2622 Bool_t bV0MCIsALambda = (iPdgCodeParticleMC==-iPdgCodeLambda);
2623
87e23d3e 2624 Double_t dPtV0Gen = particleMC->Pt();
2625 Double_t dRapV0Gen = particleMC->Y();
2626 Double_t dEtaV0Gen = particleMC->Eta();
7729cc34 2627
2628 // V0 rapidity cut
2629 if (bCutRapV0)
2630 {
87e23d3e 2631 if (bPrintCuts) printf("Gen: Applying cut: V0 |y|: < %f\n",dRapMax);
2632 if ( (TMath::Abs(dRapV0Gen) > dRapMax) )
7729cc34 2633 continue;
2634 }
2635 // V0 pseudorapidity cut
2636 if (bCutEtaV0)
2637 {
87e23d3e 2638 if (bPrintCuts) printf("Gen: Applying cut: V0 |eta|: < %f\n",dEtaMax);
2639 if ( (TMath::Abs(dEtaV0Gen) > dEtaMax) )
7729cc34 2640 continue;
2641 }
2642 /*
2643 // Is MC V0 particle physical primary? Attention!! Definition of IsPhysicalPrimary may change!!
2644 Bool_t bV0MCIsPrimary = particleMC->IsPhysicalPrimary();
2645
2646 // Get the MC mother particle of the MC V0 particle
2647 Int_t iIndexMotherOfMother = particleMC->GetMother();
2648 AliAODMCParticle* particleMCMotherOfMother = 0;
2649 if (iIndexMotherOfMother >= 0)
2650 particleMCMotherOfMother = (AliAODMCParticle*)arrayMC->At(iIndexMotherOfMother);
2651 // Get identity of the MC mother particle of the MC V0 particle if it exists
2652 Int_t iPdgCodeMotherOfMother = 0;
2653 if (particleMCMotherOfMother)
2654 iPdgCodeMotherOfMother = particleMCMotherOfMother->GetPdgCode();
2655 // Check if the MC mother particle is a physical primary Sigma
2656 Bool_t bV0MCComesFromSigma = kFALSE;
2657 if ((particleMCMotherOfMother && particleMCMotherOfMother->IsPhysicalPrimary()) && (TMath::Abs(iPdgCodeMotherOfMother)==3212) || (TMath::Abs(iPdgCodeMotherOfMother)==3224) || (TMath::Abs(iPdgCodeMotherOfMother)==3214) || (TMath::Abs(iPdgCodeMotherOfMother)==3114) )
2658 bV0MCComesFromSigma = kTRUE;
2659 // Should the MC V0 particle be considered as primary when it is Lambda?
2660 Bool_t bV0MCIsPrimaryLambda = (bV0MCIsPrimary || bV0MCComesFromSigma);
2661 */
2662 // Reject non primary particles
2663// if (!bV0MCIsPrimaryLambda)
2664// continue;
2665
2666 // Get the distance between the production point of the MC V0 particle and the primary vertex
2667 Double_t dx = dPrimVtxMCX-particleMC->Xv();
2668 Double_t dy = dPrimVtxMCY-particleMC->Yv();
2669 Double_t dz = dPrimVtxMCZ-particleMC->Zv();
87e23d3e 2670 Double_t dDistPrimary = TMath::Sqrt(dx*dx + dy*dy + dz*dz);
2671 Bool_t bV0MCIsPrimaryDist = (dDistPrimary < dDistPrimaryMax); // Is close enough to be considered primary-like?
7729cc34 2672
2673 // Check whether the MC V0 particle is in a MC jet
2674 AliAODJet* jetMC = 0;
2675 Bool_t bIsMCV0InJet = kFALSE;
2676 if (iNJetSel)
2677 {
2678 if(fDebug>5) printf("TaskV0sInJets: Searching for gen V0 in %d MC jets\n",iNJetSel);
2679 for (Int_t iJet = 0; iJet<iNJetSel; iJet++)
2680 {
2681 jetMC = (AliAODJet*)jetArraySel->At(iJet); // load a jet in the list
2682 if(fDebug>5) printf("TaskV0sInJets: Checking if gen V0 in MC jet %d\n",iJet);
87e23d3e 2683 if (IsParticleInCone(particleMC,jetMC,fdRadiusJet)) // If good jet in event, find out whether V0 is in that jet
7729cc34 2684 {
2685 if(fDebug>5) printf("TaskV0sInJets: gen V0 found in MC jet %d\n",iJet);
2686 bIsMCV0InJet = kTRUE;
2687 break;
2688 }
2689 }
2690 }
2691
2692 // Select only primary-like MC V0 particles
2693 // K0s
2694// if (bV0MCIsK0s && bV0MCIsPrimary) // well reconstructed candidates
2695 if (bV0MCIsK0s && bV0MCIsPrimaryDist) // well reconstructed candidates
2696 {
2697// if (fbTreeOutput)
2698// new ((*fBranchV0Gen)[iNV0SelV0Gen++]) AliAODMCParticle(*((AliAODMCParticle*)particleMC));
87e23d3e 2699 fh1V0K0sPtMCGen[iCentIndex]->Fill(dPtV0Gen);
2700 fh2V0K0sEtaPtMCGen[iCentIndex]->Fill(dPtV0Gen,dEtaV0Gen);
7729cc34 2701 if (bIsMCV0InJet)
2702 {
87e23d3e 2703 fh2V0K0sInJetPtMCGen[iCentIndex]->Fill(dPtV0Gen,jetMC->Pt());
2704 Double_t valueEtaKInGen[4] = {dPtV0Gen,dEtaV0Gen,jetMC->Pt(),dEtaV0Gen-jetMC->Eta()};
7729cc34 2705 fh3V0K0sInJetEtaPtMCGen[iCentIndex]->Fill(valueEtaKInGen);
2706 }
2707 }
2708 // Lambda
2709// if (bV0MCIsLambda && bV0MCIsPrimaryLambda) // well reconstructed candidates
2710 if (bV0MCIsLambda && bV0MCIsPrimaryDist) // well reconstructed candidates
2711 {
2712// if (fbTreeOutput)
2713// new ((*fBranchV0Gen)[iNV0SelV0Gen++]) AliAODMCParticle(*((AliAODMCParticle*)particleMC));
87e23d3e 2714 fh1V0LambdaPtMCGen[iCentIndex]->Fill(dPtV0Gen);
2715 fh2V0LambdaEtaPtMCGen[iCentIndex]->Fill(dPtV0Gen,dEtaV0Gen);
7729cc34 2716 if (bIsMCV0InJet)
2717 {
87e23d3e 2718 fh2V0LambdaInJetPtMCGen[iCentIndex]->Fill(dPtV0Gen,jetMC->Pt());
2719 Double_t valueEtaLInGen[4] = {dPtV0Gen,dEtaV0Gen,jetMC->Pt(),dEtaV0Gen-jetMC->Eta()};
7729cc34 2720 fh3V0LambdaInJetEtaPtMCGen[iCentIndex]->Fill(valueEtaLInGen);
2721 }
2722 }
2723 // anti-Lambda
2724// if (bV0MCIsALambda && bV0MCIsPrimaryALambda) // well reconstructed candidates
2725 if (bV0MCIsALambda && bV0MCIsPrimaryDist) // well reconstructed candidates
2726 {
2727// if (fbTreeOutput)
2728// new ((*fBranchV0Gen)[iNV0SelV0Gen++]) AliAODMCParticle(*((AliAODMCParticle*)particleMC));
87e23d3e 2729 fh1V0ALambdaPtMCGen[iCentIndex]->Fill(dPtV0Gen);
2730 fh2V0ALambdaEtaPtMCGen[iCentIndex]->Fill(dPtV0Gen,dEtaV0Gen);
7729cc34 2731 if (bIsMCV0InJet)
2732 {
87e23d3e 2733 fh2V0ALambdaInJetPtMCGen[iCentIndex]->Fill(dPtV0Gen,jetMC->Pt());
2734 Double_t valueEtaALInGen[4] = {dPtV0Gen,dEtaV0Gen,jetMC->Pt(),dEtaV0Gen-jetMC->Eta()};
7729cc34 2735 fh3V0ALambdaInJetEtaPtMCGen[iCentIndex]->Fill(valueEtaALInGen);
2736 }
2737 }
2738 }
2739 }
2740
2741// if (fbTreeOutput)
2742// ftreeOut->Fill();
2743
2744 jetArraySel->Delete();
2745 delete jetArraySel;
2746 jetArrayPerp->Delete();
2747 delete jetArrayPerp;
2748 if (jetRnd)
2749 delete jetRnd;
2750 jetRnd = 0;
2751
2752 PostData(1,fOutputListStd);
2753 PostData(2,fOutputListQA);
2754 PostData(3,fOutputListCuts);
2755 PostData(4,fOutputListMC);
2756// if (fbTreeOutput)
2757// PostData(5,ftreeOut);
2758// if(fDebug>5) printf("TaskV0sInJets: UserExec: End\n");
2759}
2760
2761void AliAnalysisTaskV0sInJets::FillQAHistogramV0(AliAODVertex* vtx, const AliAODv0* vZero, Int_t iIndexHisto, Bool_t IsCandK0s, Bool_t IsCandLambda, Bool_t IsInPeakK0s, Bool_t IsInPeakLambda)
2762{
2763 if (!IsCandK0s && !IsCandLambda)
2764 return;
2765
87e23d3e 2766// Double_t dMassK0s = vZero->MassK0Short();
2767// Double_t dMassLambda = vZero->MassLambda();
7729cc34 2768
2769 fh1QAV0Status[iIndexHisto]->Fill(vZero->GetOnFlyStatus());
2770
2771 AliAODTrack* trackNeg=(AliAODTrack*)vZero->GetDaughter(1); // negative track
2772 AliAODTrack* trackPos=(AliAODTrack*)vZero->GetDaughter(0); // positive track
2773
2774 Short_t fTotalCharge = 0;
2775 for (Int_t i = 0; i < 2; i++)
2776 {
2777 AliAODTrack* track = (AliAODTrack*)vZero->GetDaughter(i); // track
2778 // Tracks TPC OK
2779 fh1QAV0TPCRefit[iIndexHisto]->Fill(track->IsOn(AliAODTrack::kTPCrefit));
2780 Double_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
2781 fh1QAV0TPCRows[iIndexHisto]->Fill(nCrossedRowsTPC);
2782 Int_t findable = track->GetTPCNclsF();
2783 fh1QAV0TPCFindable[iIndexHisto]->Fill(findable);
2784 if (findable != 0)
2785 {
2786 fh1QAV0TPCRowsFind[iIndexHisto]->Fill(nCrossedRowsTPC/findable);
2787 }
2788 // Daughters: pseudo-rapidity cut
2789 fh1QAV0Eta[iIndexHisto]->Fill(track->Eta());
2790 if ( (nCrossedRowsTPC > (160./(250.-85.)*(255.*TMath::Abs(tan(track->Theta()))-85.))+20.) && (track->Eta() < 0) && (track->Pt() > 0.15) )
2791// if (IsCandK0s)
2792 {
2793 fh2QAV0EtaRows[iIndexHisto]->Fill(track->Eta(),nCrossedRowsTPC);
2794 fh2QAV0PtRows[iIndexHisto]->Fill(track->Pt(),nCrossedRowsTPC);
2795 fh2QAV0PhiRows[iIndexHisto]->Fill(track->Phi(),nCrossedRowsTPC);
2796 fh2QAV0NClRows[iIndexHisto]->Fill(findable,nCrossedRowsTPC);
2797 fh2QAV0EtaNCl[iIndexHisto]->Fill(track->Eta(),findable);
2798 }
2799
2800 // Daughters: transverse momentum cut
2801 fh1QAV0Pt[iIndexHisto]->Fill(track->Pt());
2802 fTotalCharge+=track->Charge();
2803 }
2804 fh1QAV0Charge[iIndexHisto]->Fill(fTotalCharge);
2805
2806 // Daughters: Impact parameter of daughters to prim vtx
2807 fh1QAV0DCAVtx[iIndexHisto]->Fill(TMath::Abs(vZero->DcaNegToPrimVertex()));
2808 fh1QAV0DCAVtx[iIndexHisto]->Fill(TMath::Abs(vZero->DcaPosToPrimVertex()));
87e23d3e 2809// fh2CutDCAVtx[iIndexHisto]->Fill(dMassK0s,TMath::Abs(vZero->DcaNegToPrimVertex()));
7729cc34 2810
2811 // Daughters: DCA
2812 fh1QAV0DCAV0[iIndexHisto]->Fill(vZero->DcaV0Daughters());
87e23d3e 2813// fh2CutDCAV0[iIndexHisto]->Fill(dMassK0s,vZero->DcaV0Daughters());
7729cc34 2814
2815 // V0: Cosine of the pointing angle
2816 fh1QAV0Cos[iIndexHisto]->Fill(vZero->CosPointingAngle(vtx));
87e23d3e 2817// fh2CutCos[iIndexHisto]->Fill(dMassK0s,vZero->CosPointingAngle(vtx));
7729cc34 2818
2819 // V0: Fiducial volume
2820 Double_t xyz[3];
2821 vZero->GetSecondaryVtx(xyz);
2822 Double_t r2=xyz[0]*xyz[0] + xyz[1]*xyz[1];
2823 fh1QAV0R[iIndexHisto]->Fill(TMath::Sqrt(r2));
2824
87e23d3e 2825 Double_t dAlpha = vZero->AlphaV0();
2826 Double_t dPtArm = vZero->PtArmV0();
7729cc34 2827
2828 if (IsCandK0s)
2829 {
2830 if (IsInPeakK0s)
2831 {
2832// fh2QAV0EtaPtK0sPeak[iIndexHisto]->Fill(trackNeg->Eta(),vZero->Pt());
2833// fh2QAV0EtaPtK0sPeak[iIndexHisto]->Fill(trackPos->Eta(),vZero->Pt());
2834 fh2QAV0EtaPtK0sPeak[iIndexHisto]->Fill(vZero->Eta(),vZero->Pt());
2835 fh2QAV0PtPtK0sPeak[iIndexHisto]->Fill(trackNeg->Pt(),trackPos->Pt());
87e23d3e 2836 fh2ArmPodK0s[iIndexHisto]->Fill(dAlpha,dPtArm);
7729cc34 2837 }
2838 fh2QAV0EtaEtaK0s[iIndexHisto]->Fill(trackNeg->Eta(),trackPos->Eta());
2839 fh2QAV0PhiPhiK0s[iIndexHisto]->Fill(trackNeg->Phi(),trackPos->Phi());
2840 fh1QAV0RapK0s[iIndexHisto]->Fill(vZero->RapK0Short());
2841 }
2842
2843 if (IsCandLambda)
2844 {
2845 if (IsInPeakLambda)
2846 {
2847// fh2QAV0EtaPtLambdaPeak[iIndexHisto]->Fill(trackNeg->Eta(),vZero->Pt());
2848// fh2QAV0EtaPtLambdaPeak[iIndexHisto]->Fill(trackPos->Eta(),vZero->Pt());
2849 fh2QAV0EtaPtLambdaPeak[iIndexHisto]->Fill(vZero->Eta(),vZero->Pt());
2850 fh2QAV0PtPtLambdaPeak[iIndexHisto]->Fill(trackNeg->Pt(),trackPos->Pt());
87e23d3e 2851 fh2ArmPodLambda[iIndexHisto]->Fill(dAlpha,dPtArm);
7729cc34 2852 }
2853 fh2QAV0EtaEtaLambda[iIndexHisto]->Fill(trackNeg->Eta(),trackPos->Eta());
2854 fh2QAV0PhiPhiLambda[iIndexHisto]->Fill(trackNeg->Phi(),trackPos->Phi());
2855 fh1QAV0RapLambda[iIndexHisto]->Fill(vZero->RapLambda());
2856 }
2857
87e23d3e 2858 fh2ArmPod[iIndexHisto]->Fill(dAlpha,dPtArm);
7729cc34 2859
2860}
2861
2862void AliAnalysisTaskV0sInJets::FillCandidates(Double_t mK, Double_t mL, Double_t mAL, Bool_t isK, Bool_t isL, Bool_t isAL, Int_t iCut/*cut index*/, Int_t iCent/*cent index*/)
2863{
2864 if (isK)
2865 {
2866 fh1V0CounterCentK0s[iCent]->Fill(iCut);
2867 fh1V0InvMassK0sAll[iCut]->Fill(mK);
2868 }
2869 if (isL)
2870 {
2871 fh1V0CounterCentLambda[iCent]->Fill(iCut);
2872 fh1V0InvMassLambdaAll[iCut]->Fill(mL);
2873 }
2874 if (isAL)
2875 {
2876 fh1V0CounterCentALambda[iCent]->Fill(iCut);
2877 fh1V0InvMassALambdaAll[iCut]->Fill(mAL);
2878 }
2879}
2880
2881Bool_t AliAnalysisTaskV0sInJets::IsParticleInCone(const AliVParticle* part1, const AliVParticle* part2, Double_t dRMax) const
2882{
2883// decides whether a particle is inside a jet cone
2884 if (!part1 || !part2)
2885 return kFALSE;
2886
2887 TVector3 vecMom2(part2->Px(),part2->Py(),part2->Pz());
2888 TVector3 vecMom1(part1->Px(),part1->Py(),part1->Pz());
2889 Double_t dR = vecMom2.DeltaR(vecMom1); // = sqrt(dEta*dEta+dPhi*dPhi)
2890 if(dR<dRMax) // momentum vectors of part1 and part2 are closer than dRMax
2891 return kTRUE;
2892 return kFALSE;
2893}
2894
2895Bool_t AliAnalysisTaskV0sInJets::OverlapWithJets(const TClonesArray* array, const AliVParticle* part, Double_t dDistance) const
2896{
2897// decides whether a cone overlaps with other jets
2898 if (!part)
2899 {
87e23d3e 2900 if(fDebug>0) printf("AliAnalysisTaskV0sInJets::OverlapWithJets: Error: No part\n");
7729cc34 2901 return kFALSE;
2902 }
2903 if (!array)
2904 {
87e23d3e 2905 if(fDebug>0) printf("AliAnalysisTaskV0sInJets::OverlapWithJets: Error: No array\n");
7729cc34 2906 return kFALSE;
2907 }
2908 Int_t iNJets = array->GetEntriesFast();
2909 if (iNJets<=0)
2910 {
87e23d3e 2911 if(fDebug>2) printf("AliAnalysisTaskV0sInJets::OverlapWithJets: Warning: No jets\n");
7729cc34 2912 return kFALSE;
2913 }
2914 AliVParticle* jet = 0;
2915 for (Int_t iJet=0; iJet<iNJets; iJet++)
2916 {
2917 jet = (AliVParticle*)array->At(iJet);
2918 if (!jet)
2919 {
87e23d3e 2920 if(fDebug>0) printf("AliAnalysisTaskV0sInJets::OverlapWithJets: Error: Failed to load jet %d/%d\n",iJet,iNJets);
7729cc34 2921 continue;
2922 }
2923 if (IsParticleInCone(part,jet,dDistance))
2924 return kTRUE;
2925 }
2926 return kFALSE;
2927}
2928
9706bde3 2929AliAODJet* AliAnalysisTaskV0sInJets::GetRandomCone(const TClonesArray* array, Double_t dEtaConeMax, Double_t dDistance) const
7729cc34 2930{
2931// generate a random cone which does not overlap with selected jets
2932// printf("Generating random cone...\n");
2933 TLorentzVector vecCone;
2934 AliAODJet* part = 0;
2935 Double_t dEta, dPhi;
2936 Int_t iNTrialsMax = 10;
2937 Bool_t bStatus = kFALSE;
2938 for (Int_t iTry=0; iTry<iNTrialsMax; iTry++)
2939 {
2940// printf("Try %d\n",iTry);
9706bde3 2941 dEta = dEtaConeMax*(2*fRandom->Rndm()-1.); // random eta in [-dEtaConeMax,+dEtaConeMax]
7729cc34 2942 dPhi = TMath::TwoPi()*fRandom->Rndm(); // random phi in [0,2*Pi]
2943 vecCone.SetPtEtaPhiM(1.,dEta,dPhi,0.);
2944 part = new AliAODJet(vecCone);
2945 if (!OverlapWithJets(array,part,dDistance))
2946 {
2947 bStatus = kTRUE;
2948// printf("Success\n");
2949 break;
2950 }
2951 else
2952 delete part;
2953 }
2954 if (!bStatus)
2955 part = 0;
2956 return part;
2957}
2958
9706bde3 2959AliAODJet* AliAnalysisTaskV0sInJets::GetMedianCluster(const TClonesArray* array, Double_t dEtaConeMax) const
2960{
f8afd72f 2961// sort kt clusters by pT/area and return the middle one, based on code in AliAnalysisTaskJetChem
9706bde3 2962 if (!array)
2963 {
2964 if(fDebug>0) printf("AliAnalysisTaskV0sInJets::GetMedianCluster: Error: No array\n");
2965 return NULL;
2966 }
2967 Int_t iNCl = array->GetEntriesFast();
f8afd72f 2968// Int_t iNClE = array->GetEntries();
9706bde3 2969 if (iNCl<3) // need at least 3 clusters (skipping 2 highest)
2970 {
2971 if(fDebug>2) printf("AliAnalysisTaskV0sInJets::GetMedianCluster: Warning: Too little clusters\n");
2972 return NULL;
2973 }
f8afd72f 2974// printf("AliAnalysisTaskV0sInJets::GetMedianCluster: EntriesFast: %d, Entries: %d\n",iNCl,iNClE);
9706bde3 2975
2976 // get list of densities
2977 Double_t* dBgDensity = new Double_t[iNCl];
2978 Int_t* iIndexList = new Int_t[iNCl];
2979 for (Int_t ij=0; ij<iNCl; ij++)
2980 {
2981 AliAODJet* clusterBg = (AliAODJet*)(array->At(ij));
2982 if (!clusterBg)
2983 {
f8afd72f 2984// printf("AliAnalysisTaskV0sInJets::GetMedianCluster: cluster %d/%d not ok\n",ij,iNCl);
9706bde3 2985 delete[] dBgDensity;
2986 delete[] iIndexList;
2987 return NULL;
2988 }
2989 Double_t dPtBg = clusterBg->Pt();
2990 Double_t dAreaBg = clusterBg->EffectiveAreaCharged();
2991
2992 Double_t dDensityBg = 0;
2993 if(dAreaBg>0)
2994 dDensityBg = dPtBg/dAreaBg;
2995 dBgDensity[ij] = dDensityBg;
2996 iIndexList[ij] = ij;
2997 }
f8afd72f 2998 // sort iIndexList by dBgDensity in descending order
9706bde3 2999 TMath::Sort(iNCl, dBgDensity, iIndexList);
3000
3001 // get median cluster with median density
3002 AliAODJet* clusterMed = 0;
3003 Int_t iIndexMed = 0;
3004 if (TMath::Odd(iNCl)) // odd number of clusters
3005 {
3006 iIndexMed = iIndexList[(Int_t) (0.5*(iNCl+1))]; // = (n - skip + 1)/2 + 1, skip = 2
3007 }
3008 else // even number: picking randomly one of the two closest to median
3009 {
3010 Int_t iIndexMed1 = iIndexList[(Int_t) (0.5*iNCl)]; // = (n - skip)/2 + 1, skip = 2
f8afd72f 3011 Int_t iIndexMed2 = iIndexList[(Int_t) (0.5*iNCl+1)]; // = (n - skip)/2 + 1 + 1, skip = 2
3012 iIndexMed = ( (fRandom->Rndm()>0.5) ? iIndexMed1 : iIndexMed2 ); // select one randomly to avoid adding areas
9706bde3 3013 }
f8afd72f 3014// printf("AliAnalysisTaskV0sInJets::GetMedianCluster: getting median cluster %d/%d ok\n",iIndexMed,iNCl);
9706bde3 3015 clusterMed = (AliAODJet*)(array->At(iIndexMed));
3016
3017 delete[] dBgDensity;
3018 delete[] iIndexList;
3019
f8afd72f 3020// printf("AliAnalysisTaskV0sInJets::GetMedianCluster: checking eta cut %g\n",dEtaConeMax);
3021// printf("AliAnalysisTaskV0sInJets::GetMedianCluster: checking eta cut |%g| < %g?\n",clusterMed->Eta(),dEtaConeMax);
9706bde3 3022 if (TMath::Abs(clusterMed->Eta())>dEtaConeMax)
3023 return NULL;
f8afd72f 3024// printf("AliAnalysisTaskV0sInJets::GetMedianCluster: cluster %d/%d passed\n",iIndexMed,iNCl);
9706bde3 3025 return clusterMed;
3026}
3027
7729cc34 3028Double_t AliAnalysisTaskV0sInJets::AreaCircSegment(Double_t dRadius, Double_t dDistance) const
3029{
3030// calculate area of a circular segment defined by the circle radius and the (oriented) distance between the secant line and the circle centre
3031 Double_t dEpsilon = 1e-2;
3032 Double_t dR = dRadius;
3033 Double_t dD = dDistance;
3034 if (TMath::Abs(dR)<dEpsilon)
3035 {
87e23d3e 3036 if(fDebug>0) printf("AliAnalysisTaskV0sInJets::AreaCircSegment: Error: Too small radius: %f < %f\n",dR,dEpsilon);
7729cc34 3037 return 0.;
3038 }
3039 if (dD>=dR)
3040 return 0.;
3041 if (dD<=-dR)
3042 return TMath::Pi()*dR*dR;
3043 return dR*dR*TMath::ACos(dD/dR)-dD*TMath::Sqrt(dR*dR-dD*dD);
3044}
3045
87e23d3e 3046Bool_t AliAnalysisTaskV0sInJets::IsSelectedForJets(AliAODEvent* fAOD, Double_t dVtxZCut, Double_t dVtxR2Cut, Double_t dCentCutLo, Double_t dCentCutUp, Bool_t bCutDeltaZ, Double_t dDeltaZMax)
7729cc34 3047{
3048// event selection
3049 AliAODVertex* vertex = fAOD->GetPrimaryVertex();
3050 if (!vertex)
3051 return kFALSE;
3052 if (vertex->GetNContributors() < 3)
3053 return kFALSE;
3054 TString vtxTitle(vertex->GetTitle());
3055 if (vtxTitle.Contains("TPCVertex"))
3056 return kFALSE;
3057 Double_t zVertex = vertex->GetZ();
87e23d3e 3058 if (TMath::Abs(zVertex) > dVtxZCut)
7729cc34 3059 return kFALSE;
3060 if (bCutDeltaZ)
3061 {
3062 AliAODVertex* vertexSPD = fAOD->GetPrimaryVertexSPD();
3063 if (!vertexSPD)
3064 {
3065// printf("IsSelectedForJets: Error: No SPD vertex\n");
3066 return kFALSE;
3067 }
3068 Double_t zVertexSPD = vertexSPD->GetZ();
87e23d3e 3069 if (TMath::Abs(zVertex-zVertexSPD) > dDeltaZMax)
7729cc34 3070 {
3071// printf("IsSelectedForJets: Rejecting event due to delta z = %f - %f = %f\n",zVertex,zVertexSPD,zVertex-zVertexSPD);
3072 return kFALSE;
3073 }
3074// printf("IsSelectedForJets: Event OK: %f - %f = %f\n",zVertex,zVertexSPD,zVertex-zVertexSPD);
3075 }
3076 Double_t xVertex = vertex->GetX();
3077 Double_t yVertex = vertex->GetY();
3078 Double_t radiusSq = yVertex*yVertex+xVertex*xVertex;
87e23d3e 3079 if (radiusSq > dVtxR2Cut)
7729cc34 3080 return kFALSE;
3081 Double_t centrality;
3082// centrality = fAOD->GetHeader()->GetCentrality();
3083 centrality = fAOD->GetHeader()->GetCentralityP()->GetCentralityPercentile("V0M");
3084 if (centrality < 0)
3085 return kFALSE;
87e23d3e 3086 if( (dCentCutUp < 0) || (dCentCutLo < 0) || (dCentCutUp > 100) || (dCentCutLo > 100) || (dCentCutLo > dCentCutUp) )
7729cc34 3087 return kFALSE;
87e23d3e 3088 if ( (centrality < dCentCutLo) || (centrality > dCentCutUp) )
7729cc34 3089 return kFALSE;
3090 return kTRUE;
3091}
3092
3093Int_t AliAnalysisTaskV0sInJets::GetCentralityBinIndex(Double_t centrality)
3094{
3095// returns index of the centrality bin corresponding to the provided value of centrality
3096 if (centrality < 0 || centrality > fgkiCentBinRanges[fgkiNBinsCent-1])
3097 return -1;
3098 for (Int_t i = 0; i < fgkiNBinsCent; i++)
3099 {
3100 if (centrality <= fgkiCentBinRanges[i])
3101 return i;
3102 }
3103 return -1;
3104}
3105
3106Int_t AliAnalysisTaskV0sInJets::GetCentralityBinEdge(Int_t index)
3107{
3108// returns the upper edge of the centrality bin corresponding to the provided value of index
3109 if (index < 0 || index >= fgkiNBinsCent)
3110 return -1;
3111 return fgkiCentBinRanges[index];
3112}
3113
3114TString AliAnalysisTaskV0sInJets::GetCentBinLabel(Int_t index)
3115{
3116// get string with centrality range for given bin
3117 TString lowerEdge = ( (index == 0) ? "0" : Form("%d",GetCentralityBinEdge(index-1)));
3118 TString upperEdge = Form("%d",GetCentralityBinEdge(index));
3119 return Form("%s-%s %%",lowerEdge.Data(),upperEdge.Data());
3120}
3121
3122Double_t AliAnalysisTaskV0sInJets::MassPeakSigmaOld(Double_t pt, Int_t particle)
3123{
3124// estimation of the sigma of the invariant-mass peak as a function of pT and particle type
3125 switch (particle)
3126 {
3127 case 0: // K0S
3128 return 0.0044 + 0.0004*(pt - 1.);
3129 break;
3130 case 1: // Lambda
3131 return 0.0023 + 0.00034*(pt - 1.);
3132 break;
3133 default:
3134 return 0;
3135 break;
3136 }
3137}