1 /**************************************************************************
2 * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 //-----------------------------------------------------------------
17 // AliAnalysisTaskPerformanceStrange class
18 // This task is for a performance study of V0 identification.
19 // It works with MC info and ESD tree.
20 // Author: H.Ricaud, H.Ricaud@gsi.de
21 //-----------------------------------------------------------------
23 #include <Riostream.h>
36 #include "AliAnalysisManager.h"
38 #include "AliPhysicsSelection.h"
39 #include "AliBackgroundSelection.h"
41 #include "AliESDVertex.h"
42 #include "AliESDEvent.h"
43 #include "AliESDInputHandler.h"
44 #include "AliESDtrack.h"
46 #include "AliESDtrackCuts.h"
47 #include "AliESDpid.h"
48 #include "AliMultiplicity.h"
50 #include "AliAODEvent.h"
51 #include "AliAODVertex.h"
52 #include "AliAODTrack.h"
54 #include "AliAODMCHeader.h"
55 #include "AliAODInputHandler.h"
57 #include "AliAODMCParticle.h"
59 #include "AliMCEventHandler.h"
60 #include "AliMCEvent.h"
62 #include "AliGenEventHeader.h"
66 #include "AliKFVertex.h"
67 #include "AliVertexerTracks.h"
69 #include "AliAnalysisTaskPerformanceStrange.h"
70 #include "AliAnalysisCentralitySelector.h"
73 ClassImp(AliAnalysisTaskPerformanceStrange)
76 //________________________________________________________________________
77 AliAnalysisTaskPerformanceStrange::AliAnalysisTaskPerformanceStrange()
78 : AliAnalysisTaskSE("TaskStrange"), fAnalysisMC(0), fAnalysisType("infoType"), fCollidingSystems(0), fUsePID("infoPID"), fUseCut("infoCut"), fUseOnTheFly(0), fESD(0), fListHist(0), fCentrSelector(0),
80 // MC histograms ---------------------------------------
81 fHistMCPrimaryVertexX(0),
82 fHistMCPrimaryVertexY(0),
83 fHistMCPrimaryVertexZ(0),
85 fHistMCMultiplicityPrimary(0),
86 fHistMCMultiplicityTracks(0),
88 fHistMCtracksProdRadiusK0s(0),
89 fHistMCtracksProdRadiusLambda(0),
90 fHistMCtracksProdRadiusAntiLambda(0),
92 fHistMCtracksDecayRadiusK0s(0),
93 fHistMCtracksDecayRadiusLambda(0),
94 fHistMCtracksDecayRadiusAntiLambda(0),
97 fHistMCPtAllLambda(0),
98 fHistMCPtAllAntiLambda(0),
100 fHistMCProdRadiusK0s(0),
101 fHistMCProdRadiusLambda(0),
102 fHistMCProdRadiusAntiLambda(0),
105 fHistMCRapInPtRangeK0s(0),
107 fHistMCRapInPtRangeLambda(0),
108 fHistMCRapAntiLambda(0),
109 fHistMCRapInPtRangeAntiLambda(0),
111 fHistMCRapInPtRangeXi(0),
113 fHistMCRapInPtRangePhi(0),
118 fHistMCPtLambdaFromSigma(0),
119 fHistMCPtAntiLambdaFromSigma(0),
121 fHistNTimesRecK0s(0),
122 fHistNTimesRecLambda(0),
123 fHistNTimesRecAntiLambda(0),
124 fHistNTimesRecK0sVsPt(0),
125 fHistNTimesRecLambdaVsPt(0),
126 fHistNTimesRecAntiLambdaVsPt(0),
127 // ------------------------------------------------------
129 // Reconstructed particle histograms -------------------
130 fHistNumberEvents(0),
131 fHistTrackPerEvent(0),
132 fHistTrackletPerEvent(0),
133 fHistMCDaughterTrack(0),
135 fHistSPDPrimaryVertexZ(0),
137 fHistPrimaryVertexX(0),
138 fHistPrimaryVertexY(0),
139 fHistPrimaryVertexZ(0),
141 fHistPrimaryVertexResX(0),
142 fHistPrimaryVertexResY(0),
143 fHistPrimaryVertexResZ(0),
145 fHistPrimaryVertexPosXV0events(0),
146 fHistPrimaryVertexPosYV0events(0),
147 fHistPrimaryVertexPosZV0events(0),
151 fHistDcaPosToPrimVertex(0),
152 fHistDcaNegToPrimVertex(0),
153 fHistDcaPosToPrimVertexZoom(0),
154 fHistDcaNegToPrimVertexZoom(0),
156 fHistDecayLengthV0(0),
157 fHistDcaV0Daughters(0),
159 fHistCosPointAngle(0),
160 fHistCosPointAngleZoom(0),
163 fHistV0Multiplicity(0),
165 fHistChi2KFBeforeCutK0s(0),
166 fHistChi2KFBeforeCutLambda(0),
167 fHistChi2KFBeforeCutAntiLambda(0),
168 fHistChi2KFAfterCutK0s(0),
169 fHistChi2KFAfterCutLambda(0),
170 fHistChi2KFAfterCutAntiLambda(0),
174 fHistMassAntiLambda(0),
175 fHistMassVsRadiusK0(0),
176 fHistMassVsRadiusLambda(0),
177 fHistMassVsRadiusAntiLambda(0),
180 fHistPtVsMassLambda(0),
181 fHistArmenterosPodolanski(0),
182 // ------------------------------------------------------
185 // PID histograms --------------------------------------
186 fHistNsigmaPosPionAntiLambda(0),
187 fHistNsigmaNegProtonAntiLambda(0),
188 fHistNsigmaPosProtonLambda(0),
189 fHistNsigmaNegPionLambda(0),
190 fHistNsigmaPosPionK0(0),
191 fHistNsigmaNegPionK0(0),
192 // ------------------------------------------------------
194 // Associated particles ---------------------------------
196 fHistAsMcRapLambda(0),
197 fHistAsMcRapAntiLambda(0),
200 fHistAsMcPtLambda(0),
202 fHistAsMcPtZoomK0(0),
203 fHistAsMcPtZoomLambda(0),
205 fHistAsMcProdRadiusK0(0),
206 fHistAsMcProdRadiusLambda(0),
207 fHistAsMcProdRadiusAntiLambda(0),
209 fHistAsMcProdRadiusXvsYK0s(0),
210 fHistAsMcProdRadiusXvsYLambda(0),
211 fHistAsMcProdRadiusXvsYAntiLambda(0),
214 fHistPidMcMassLambda(0),
215 fHistPidMcMassAntiLambda(0),
217 fHistAsMcMassLambda(0),
218 fHistAsMcMassAntiLambda(0),
220 fHistAsMcPtVsMassK0(0),
221 fHistAsMcPtVsMassLambda(0),
222 fHistAsMcPtVsMassAntiLambda(0),
224 fHistAsMcMassVsRadiusK0(0),
225 fHistAsMcMassVsRadiusLambda(0),
226 fHistAsMcMassVsRadiusAntiLambda(0),
232 fHistAsMcResrVsRadiusK0(0),
233 fHistAsMcReszVsRadiusK0(0),
235 fHistAsMcResxLambda(0),
236 fHistAsMcResyLambda(0),
237 fHistAsMcReszLambda(0),
239 fHistAsMcResrVsRadiusLambda(0),
240 fHistAsMcReszVsRadiusLambda(0),
242 fHistAsMcResxAntiLambda(0),
243 fHistAsMcResyAntiLambda(0),
244 fHistAsMcReszAntiLambda(0),
246 fHistAsMcResrVsRadiusAntiLambda(0),
247 fHistAsMcReszVsRadiusAntiLambda(0),
250 fHistAsMcResPtLambda(0),
251 fHistAsMcResPtAntiLambda(0),
253 fHistAsMcResPtVsRapK0(0),
254 fHistAsMcResPtVsRapLambda(0),
255 fHistAsMcResPtVsRapAntiLambda(0),
256 fHistAsMcResPtVsPtK0(0),
257 fHistAsMcResPtVsPtLambda(0),
258 fHistAsMcResPtVsPtAntiLambda(0),
260 fHistAsMcMotherPdgCodeK0s(0),
261 fHistAsMcMotherPdgCodeLambda(0),
262 fHistAsMcMotherPdgCodeAntiLambda(0),
264 fHistAsMcPtLambdaFromSigma(0),
265 fHistAsMcPtAntiLambdaFromSigma(0),
266 // ------------------------------------------------------
268 // Associated secondary particle histograms -------------
269 fHistAsMcSecondaryPtVsRapK0s(0),
270 fHistAsMcSecondaryPtVsRapLambda(0),
271 fHistAsMcSecondaryPtVsRapAntiLambda(0),
273 fHistAsMcSecondaryProdRadiusK0s(0),
274 fHistAsMcSecondaryProdRadiusLambda(0),
275 fHistAsMcSecondaryProdRadiusAntiLambda(0),
277 fHistAsMcSecondaryProdRadiusXvsYK0s(0),
278 fHistAsMcSecondaryProdRadiusXvsYLambda(0),
279 fHistAsMcSecondaryProdRadiusXvsYAntiLambda(0),
281 fHistAsMcSecondaryMotherPdgCodeK0s(0),
282 fHistAsMcSecondaryMotherPdgCodeLambda(0),
283 fHistAsMcSecondaryMotherPdgCodeAntiLambda(0),
285 fHistAsMcSecondaryPtLambdaFromSigma(0),
286 fHistAsMcSecondaryPtAntiLambdaFromSigma(0)
291 /* fCuts[0]=33; // max allowed chi2
292 fCuts[1]=0.05; // min allowed impact parameter for the 1st daughter
293 fCuts[2]=0.05; // min allowed impact parameter for the 2nd daughter
294 fCuts[3]=0.5; // max allowed DCA between the daughter tracks
295 fCuts[4]=0.00; // max allowed cosine of V0's pointing angle
296 fCuts[5]=0.2; // min radius of the fiducial volume
297 fCuts[6]=100; // max radius of the fiducial volume
301 //________________________________________________________________________
302 AliAnalysisTaskPerformanceStrange::AliAnalysisTaskPerformanceStrange(const char *name)
303 : AliAnalysisTaskSE(name), fAnalysisMC(0), fAnalysisType("infoType"), fCollidingSystems(0), fUsePID("infoPID"), fUseCut("infoCut"), fUseOnTheFly(0), fESD(0), fListHist(),fCentrSelector(0),
305 // MC histograms ---------------------------------------
306 fHistMCPrimaryVertexX(0),
307 fHistMCPrimaryVertexY(0),
308 fHistMCPrimaryVertexZ(0),
310 fHistMCMultiplicityPrimary(0),
311 fHistMCMultiplicityTracks(0),
313 fHistMCtracksProdRadiusK0s(0),
314 fHistMCtracksProdRadiusLambda(0),
315 fHistMCtracksProdRadiusAntiLambda(0),
317 fHistMCtracksDecayRadiusK0s(0),
318 fHistMCtracksDecayRadiusLambda(0),
319 fHistMCtracksDecayRadiusAntiLambda(0),
322 fHistMCPtAllLambda(0),
323 fHistMCPtAllAntiLambda(0),
325 fHistMCProdRadiusK0s(0),
326 fHistMCProdRadiusLambda(0),
327 fHistMCProdRadiusAntiLambda(0),
330 fHistMCRapInPtRangeK0s(0),
332 fHistMCRapInPtRangeLambda(0),
333 fHistMCRapAntiLambda(0),
334 fHistMCRapInPtRangeAntiLambda(0),
336 fHistMCRapInPtRangeXi(0),
338 fHistMCRapInPtRangePhi(0),
343 fHistMCPtLambdaFromSigma(0),
344 fHistMCPtAntiLambdaFromSigma(0),
346 fHistNTimesRecK0s(0),
347 fHistNTimesRecLambda(0),
348 fHistNTimesRecAntiLambda(0),
349 fHistNTimesRecK0sVsPt(0),
350 fHistNTimesRecLambdaVsPt(0),
351 fHistNTimesRecAntiLambdaVsPt(0),
352 // ------------------------------------------------------
354 // Reconstructed particle histograms -------------------
355 fHistNumberEvents(0),
356 fHistTrackPerEvent(0),
357 fHistTrackletPerEvent(0),
358 fHistMCDaughterTrack(0),
360 fHistSPDPrimaryVertexZ(0),
362 fHistPrimaryVertexX(0),
363 fHistPrimaryVertexY(0),
364 fHistPrimaryVertexZ(0),
366 fHistPrimaryVertexResX(0),
367 fHistPrimaryVertexResY(0),
368 fHistPrimaryVertexResZ(0),
370 fHistPrimaryVertexPosXV0events(0),
371 fHistPrimaryVertexPosYV0events(0),
372 fHistPrimaryVertexPosZV0events(0),
376 fHistDcaPosToPrimVertex(0),
377 fHistDcaNegToPrimVertex(0),
378 fHistDcaPosToPrimVertexZoom(0),
379 fHistDcaNegToPrimVertexZoom(0),
381 fHistDecayLengthV0(0),
382 fHistDcaV0Daughters(0),
384 fHistCosPointAngle(0),
385 fHistCosPointAngleZoom(0),
388 fHistV0Multiplicity(0),
390 fHistChi2KFBeforeCutK0s(0),
391 fHistChi2KFBeforeCutLambda(0),
392 fHistChi2KFBeforeCutAntiLambda(0),
393 fHistChi2KFAfterCutK0s(0),
394 fHistChi2KFAfterCutLambda(0),
395 fHistChi2KFAfterCutAntiLambda(0),
399 fHistMassAntiLambda(0),
400 fHistMassVsRadiusK0(0),
401 fHistMassVsRadiusLambda(0),
402 fHistMassVsRadiusAntiLambda(0),
405 fHistPtVsMassLambda(0),
406 fHistArmenterosPodolanski(0),
407 // ------------------------------------------------------
410 // PID histograms --------------------------------------
411 fHistNsigmaPosPionAntiLambda(0),
412 fHistNsigmaNegProtonAntiLambda(0),
413 fHistNsigmaPosProtonLambda(0),
414 fHistNsigmaNegPionLambda(0),
415 fHistNsigmaPosPionK0(0),
416 fHistNsigmaNegPionK0(0),
417 // ------------------------------------------------------
419 // Associated particles ---------------------------------
421 fHistAsMcRapLambda(0),
422 fHistAsMcRapAntiLambda(0),
425 fHistAsMcPtLambda(0),
427 fHistAsMcPtZoomK0(0),
428 fHistAsMcPtZoomLambda(0),
430 fHistAsMcProdRadiusK0(0),
431 fHistAsMcProdRadiusLambda(0),
432 fHistAsMcProdRadiusAntiLambda(0),
434 fHistAsMcProdRadiusXvsYK0s(0),
435 fHistAsMcProdRadiusXvsYLambda(0),
436 fHistAsMcProdRadiusXvsYAntiLambda(0),
439 fHistPidMcMassLambda(0),
440 fHistPidMcMassAntiLambda(0),
442 fHistAsMcMassLambda(0),
443 fHistAsMcMassAntiLambda(0),
445 fHistAsMcPtVsMassK0(0),
446 fHistAsMcPtVsMassLambda(0),
447 fHistAsMcPtVsMassAntiLambda(0),
449 fHistAsMcMassVsRadiusK0(0),
450 fHistAsMcMassVsRadiusLambda(0),
451 fHistAsMcMassVsRadiusAntiLambda(0),
457 fHistAsMcResrVsRadiusK0(0),
458 fHistAsMcReszVsRadiusK0(0),
460 fHistAsMcResxLambda(0),
461 fHistAsMcResyLambda(0),
462 fHistAsMcReszLambda(0),
464 fHistAsMcResrVsRadiusLambda(0),
465 fHistAsMcReszVsRadiusLambda(0),
467 fHistAsMcResxAntiLambda(0),
468 fHistAsMcResyAntiLambda(0),
469 fHistAsMcReszAntiLambda(0),
471 fHistAsMcResrVsRadiusAntiLambda(0),
472 fHistAsMcReszVsRadiusAntiLambda(0),
475 fHistAsMcResPtLambda(0),
476 fHistAsMcResPtAntiLambda(0),
478 fHistAsMcResPtVsRapK0(0),
479 fHistAsMcResPtVsRapLambda(0),
480 fHistAsMcResPtVsRapAntiLambda(0),
481 fHistAsMcResPtVsPtK0(0),
482 fHistAsMcResPtVsPtLambda(0),
483 fHistAsMcResPtVsPtAntiLambda(0),
485 fHistAsMcMotherPdgCodeK0s(0),
486 fHistAsMcMotherPdgCodeLambda(0),
487 fHistAsMcMotherPdgCodeAntiLambda(0),
489 fHistAsMcPtLambdaFromSigma(0),
490 fHistAsMcPtAntiLambdaFromSigma(0),
491 // ------------------------------------------------------
493 // Associated secondary particle histograms -------------
494 fHistAsMcSecondaryPtVsRapK0s(0),
495 fHistAsMcSecondaryPtVsRapLambda(0),
496 fHistAsMcSecondaryPtVsRapAntiLambda(0),
498 fHistAsMcSecondaryProdRadiusK0s(0),
499 fHistAsMcSecondaryProdRadiusLambda(0),
500 fHistAsMcSecondaryProdRadiusAntiLambda(0),
502 fHistAsMcSecondaryProdRadiusXvsYK0s(0),
503 fHistAsMcSecondaryProdRadiusXvsYLambda(0),
504 fHistAsMcSecondaryProdRadiusXvsYAntiLambda(0),
506 fHistAsMcSecondaryMotherPdgCodeK0s(0),
507 fHistAsMcSecondaryMotherPdgCodeLambda(0),
508 fHistAsMcSecondaryMotherPdgCodeAntiLambda(0),
510 fHistAsMcSecondaryPtLambdaFromSigma(0),
511 fHistAsMcSecondaryPtAntiLambdaFromSigma(0)
516 /* fCuts[0]=33; // max allowed chi2
517 fCuts[1]=0.05; // min allowed impact parameter for the 1st daughter
518 fCuts[2]=0.05; // min allowed impact parameter for the 2nd daughter
519 fCuts[3]=0.5; // max allowed DCA between the daughter tracks
520 fCuts[4]=0.00; // max allowed cosine of V0's pointing angle
521 fCuts[5]=0.2; // min radius of the fiducial volume
522 fCuts[6]=100; // max radius of the fiducial volume
524 // Define output slots only here
525 // Output slot #1 writes into a TList container
526 DefineOutput(1, TList::Class());
527 DefineOutput(2, AliAnalysisCentralitySelector::Class());
531 //________________________________________________________________________
532 void AliAnalysisTaskPerformanceStrange::UserCreateOutputObjects()
536 //*******************
537 fListHist = new TList();
538 fListHist->SetOwner();
540 // Bo: tbd: condition before allocation (i.e. if (!fHistMCMultiplicityPrimary){...} for each histo...
546 // Primary Vertex X,Y,Z:
547 fHistMCPrimaryVertexX = new TH1F("h1MCPrimaryVertexX", "MC Primary Vertex Position X;Primary Vertex Position X (cm);Events",100,-0.5,0.5);
548 fListHist->Add(fHistMCPrimaryVertexX);
550 fHistMCPrimaryVertexY = new TH1F("h1MCPrimaryVertexY", "MC Primary Vertex Position Y;Primary Vertex Position Y (cm);Events",100,-0.5,0.5);
551 fListHist->Add(fHistMCPrimaryVertexY);
553 fHistMCPrimaryVertexZ = new TH1F("h1MCPrimaryVertexZ", "MC Primary Vertex Position Z;Primary Vertex Position Z (cm);Events",200,-20,20);
554 fListHist->Add(fHistMCPrimaryVertexZ);
557 fHistMCMultiplicityPrimary = new TH1F("h1MCMultiplicityPrimary", "MC Primary Particles;NPrimary;Count", 201, -0.5, 200.5);
558 fListHist->Add(fHistMCMultiplicityPrimary);
560 fHistMCMultiplicityTracks = new TH1F("h1MCMultiplicityTracks", "MC Tracks;Ntracks;Count", 201, -0.5, 200.5);
561 fListHist->Add(fHistMCMultiplicityTracks);
563 // Production Radius of non-primary particles:
564 fHistMCtracksProdRadiusK0s = new TH2F("h2MCtracksProdRadiusK0s","Non-primary MC K^{0} Production Radius;x (cm); y (cm)",200,-50,50,200,-50,50);
565 fListHist->Add(fHistMCtracksProdRadiusK0s);
567 fHistMCtracksProdRadiusLambda = new TH2F("h2MCtracksProdRadiusLambda","Non-primary MC #Lambda^{0} Production Radius;x (cm); y (cm)",200,-50,50,200,-50,50);
568 fListHist->Add(fHistMCtracksProdRadiusLambda);
570 fHistMCtracksProdRadiusAntiLambda = new TH2F("h2MCtracksProdRadiusAntiLambda","Non-primary MC #bar{#Lambda}^{0} Production Radius;x (cm); y (cm)",200,-50,50,200,-50,50);
571 fListHist->Add(fHistMCtracksProdRadiusAntiLambda);
573 // Decay Radius of non-primary particles:
574 fHistMCtracksDecayRadiusK0s = new TH1F("h1MCtracksDecayRadiusK0s","Non-primary MC K^{0} Decay Radius;r (cm)",101,-1,100);
575 fListHist->Add(fHistMCtracksDecayRadiusK0s);
577 fHistMCtracksDecayRadiusLambda = new TH1F("h1MCtracksDecayRadiusLambda","Non-primary MC #Lambda^{0} Decay Radius;r (cm)",101,-1,100);
578 fListHist->Add(fHistMCtracksDecayRadiusLambda);
580 fHistMCtracksDecayRadiusAntiLambda = new TH1F("h1MCtracksDecayRadiusAntiLambda","Non-primary #bar{#Lambda}^{0} Decay Radius;r (cm)",100,1,101);
581 fListHist->Add(fHistMCtracksDecayRadiusAntiLambda);
583 // Pt Distribution of non-primary particles:
584 fHistMCPtAllK0s = new TH1F("h1MCPtAllK0s", "Non-primary MC K^{0};p_{t} (GeV/c);Counts",240,0,12);
585 fListHist->Add(fHistMCPtAllK0s);
587 fHistMCPtAllLambda = new TH1F("h1MCPtAllLambda", "Non-primary MC #Lambda^{0};p_{t} (GeV/c);Counts",240,0,12);
588 fListHist->Add(fHistMCPtAllLambda);
590 fHistMCPtAllAntiLambda = new TH1F("h1MCPtAllAntiLambda", "Non-primary MC #bar{#Lambda}^{0};p_{t} (GeV/c);Counts",240,0,12);
591 fListHist->Add(fHistMCPtAllAntiLambda);
594 fHistMCProdRadiusK0s = new TH1F("h1MCProdRadiusK0s", "MC K^{0} Production Radius;r (cm);Count", 400, -2, 2);
595 fListHist->Add(fHistMCProdRadiusK0s);
597 fHistMCProdRadiusLambda = new TH1F("h1MCProdRadiusLambda", "MC #Lambda^{0} Production Radius;r (cm);Count", 400, -2, 2);
598 fListHist->Add(fHistMCProdRadiusLambda);
600 fHistMCProdRadiusAntiLambda = new TH1F("h1MCProdRadiusAntiLambda", "MC #bar{#Lambda}^{0} Production Radius;r (cm);Count", 400, -2, 2);
601 fListHist->Add(fHistMCProdRadiusAntiLambda);
604 // Rapidity distribution:
605 fHistMCRapK0s = new TH1F("h1MCRapK0s", "K^{0};y",160,-4,4);
606 fListHist->Add(fHistMCRapK0s);
608 fHistMCRapInPtRangeK0s = new TH1F("h1MCRapInPtRangeK0s", "K^{0};y",160,-4,4);
609 fListHist->Add(fHistMCRapInPtRangeK0s);
611 fHistMCRapLambda = new TH1F("h1MCRapLambda", "#Lambda;y",160,-4,4);
612 fListHist->Add(fHistMCRapLambda);
614 fHistMCRapInPtRangeLambda = new TH1F("h1MCRapInPtRangeLambda", "#Lambda;y",160,-4,4);
615 fListHist->Add(fHistMCRapInPtRangeLambda);
617 fHistMCRapAntiLambda = new TH1F("h1MCRapAntiLambda", "#bar{#Lambda};y",160,-4,4);
618 fListHist->Add(fHistMCRapAntiLambda);
620 fHistMCRapInPtRangeAntiLambda = new TH1F("h1MCRapInPtRangeAntiLambda", "#bar{#Lambda};y",160,-4,4);
621 fListHist->Add(fHistMCRapInPtRangeAntiLambda);
623 fHistMCRapXi = new TH1F("h1MCRapXi", "#Xi;y",160,-4,4);
624 fListHist->Add(fHistMCRapXi);
626 fHistMCRapInPtRangeXi = new TH1F("h1MCRapInPtRangeXi", "#Xi;y",160,-4,4);
627 fListHist->Add(fHistMCRapInPtRangeXi);
629 fHistMCRapPhi = new TH1F("h1MCRapPhi", "#phi;y",160,-4,4);
630 fListHist->Add(fHistMCRapPhi);
632 fHistMCRapInPtRangePhi = new TH1F("h1MCRapInPtRangePhi", "#phi;y",160,-4,4);
633 fListHist->Add(fHistMCRapInPtRangePhi);
637 fHistMCPtK0s = new TH1F("h1MCPtK0s", "K^{0};p_{t} (GeV/c)",240,0,12);
638 fListHist->Add(fHistMCPtK0s);
640 fHistMCPtLambda = new TH1F("h1MCPtLambda", "#Lambda^{0};p_{t} (GeV/c)",240,0,12);
641 fListHist->Add(fHistMCPtLambda);
643 // Pt distribution of Lambda coming from Sigma decay:
644 fHistMCPtLambdaFromSigma = new TH1F("h1MCPtLambdaFromSigma", "#Lambda^{0};p_{t} (GeV/c)",240,0,12);
645 fListHist->Add(fHistMCPtLambdaFromSigma);
647 fHistMCPtAntiLambdaFromSigma = new TH1F("h1MCPtAntiLambdaFromSigma", "#Lambda^{0};p_{t} (GeV/c)",240,0,12);
648 fListHist->Add(fHistMCPtAntiLambdaFromSigma);
650 // Multiple reconstruction studies:
651 fHistNTimesRecK0s = new TH1F("h1NTimesRecK0s","number of times a K0s is reconstructed in -1<y<1;number of times;counts",500,-0.5,4.5);
652 fListHist->Add(fHistNTimesRecK0s);
654 fHistNTimesRecLambda = new TH1F("h1NTimesRecLambda","number of times a Lambda is reconstructed in -1<y<1;number of times;counts",500,-0.5,4.5);
655 fListHist->Add(fHistNTimesRecLambda);
657 fHistNTimesRecAntiLambda = new TH1F("h1NTimesRecAntiLambda","number of times an AntiLambda is reconstructed in -1<y<1;number of times;counts",500,-0.5,4.5);
658 fListHist->Add(fHistNTimesRecAntiLambda);
660 fHistNTimesRecK0sVsPt = new TH2F("h2NTimesRecK0sVsPt","NTimes versus Pt, K^{0} in -1<y<1;p_{t} (GeV/c);number of times",75,0,15,5,-0.5,4.5);
661 fListHist->Add(fHistNTimesRecK0sVsPt);
663 fHistNTimesRecLambdaVsPt = new TH2F("h2NTimesRecLambdaVsPt","NTimes versus Pt, #Lambda^{0} in -1<y<1;p_{t} (GeV/c);number of times",75,0,15,5,-0.5,4.5);
664 fListHist->Add(fHistNTimesRecLambdaVsPt);
666 fHistNTimesRecAntiLambdaVsPt = new TH2F("h2NTimesRecAntiLambdaVsPt","NTimes versus Pt, #bar{#Lambda}^{0} in -1<y<1;p_{t} (GeV/c);number of times",75,0,15,5,-0.5,4.5);
667 fListHist->Add(fHistNTimesRecAntiLambdaVsPt);
669 //***********************************
670 // Reconstructed particles histograms
671 //***********************************
674 fHistNumberEvents = new TH1F("h1NumberEvents", "Number of events; index;Number of Events",10,0,10);
675 fListHist->Add(fHistNumberEvents);
678 fHistTrackPerEvent = new TH1F("h1TrackPerEvent", "Tracks per event;Number of Tracks;Number of Events",20000,0,20000);
679 fListHist->Add(fHistTrackPerEvent);
681 fHistTrackletPerEvent = new TH1F("h1TrackletPerEvent", "Number of tracklets;Number of tracklets per events;Number of events",1000,0,1000);
682 fListHist->Add(fHistTrackletPerEvent);
684 fHistMCDaughterTrack = new TH1F("h1MCDaughterTrack","Distribution of mc id for daughters;id tags;Counts",15,0,15);
685 fListHist->Add(fHistMCDaughterTrack);
688 fHistSPDPrimaryVertexZ = new TH1F("h1SPDPrimaryVertexZ", "SPD Primary Vertex Position Z;Primary Vertex Position Z (cm);Events",200,-20,20);
689 fListHist->Add(fHistSPDPrimaryVertexZ);
691 fHistPrimaryVertexX = new TH1F("h1PrimaryVertexX", "Primary Vertex Position X;Primary Vertex Position X (cm);Events",100,-0.5,0.5);
692 fListHist->Add(fHistPrimaryVertexX);
694 fHistPrimaryVertexY = new TH1F("h1PrimaryVertexY", "Primary Vertex Position Y;Primary Vertex Position Y (cm);Events",100,-0.5,0.5);
695 fListHist->Add(fHistPrimaryVertexY);
697 fHistPrimaryVertexZ = new TH1F("h1PrimaryVertexZ", "Primary Vertex Position Z;Primary Vertex Position Z (cm);Events",200,-20,20);
698 fListHist->Add(fHistPrimaryVertexZ);
701 // Primary vertex resolution:
702 fHistPrimaryVertexResX = new TH1F("h1PrimaryVertexResX", "Primary Vertex Resolution X;Primary Vertex Resolution X (cm);Events",100,-0.25,0.25);
703 fListHist->Add(fHistPrimaryVertexResX);
705 fHistPrimaryVertexResY = new TH1F("h1PrimaryVertexResY", "Primary Vertex Resolution Y;Primary Vertex Resolution Y (cm);Events",100,-0.25,0.25);
706 fListHist->Add(fHistPrimaryVertexResY);
708 fHistPrimaryVertexResZ = new TH1F("h1PrimaryVertexResZ", "Primary Vertex Resolution Z;Primary Vertex Resolution Z (cm);Events",200,-0.25,0.25);
709 fListHist->Add(fHistPrimaryVertexResZ);
712 // Primary Vertex in events with V0 candidates:
713 fHistPrimaryVertexPosXV0events = new TH1F("h1PrimaryVertexPosXV0events", "Primary Vertex Position X;Primary Vertex Position X (cm);Events",100,-0.5,0.5);
714 fListHist->Add(fHistPrimaryVertexPosXV0events);
715 fHistPrimaryVertexPosYV0events = new TH1F("h1PrimaryVertexPosYV0events", "Primary Vertex Position Y;Primary Vertex Position Y (cm);Events",100,-0.5,0.5);
716 fListHist->Add(fHistPrimaryVertexPosYV0events);
717 fHistPrimaryVertexPosZV0events = new TH1F("h1PrimaryVertexPosZV0events", "Primary Vertex Position Z;Primary Vertex Position Z (cm);Events",200,-20.0,20.0);
718 fListHist->Add(fHistPrimaryVertexPosZV0events);
721 fHistDaughterPt = new TH2F("h2DaughterPt", "Daughter Pt;Positive Daughter Pt; Negative Daughter Pt",200,0,2,200,0,2);
722 fListHist->Add(fHistDaughterPt);
725 fHistDcaPosToPrimVertex = new TH2F("h2DcaPosToPrimVertex", "Positive V0 daughter;dca(cm);Status",500,0,5,2,-0.5,1.5);
726 fListHist->Add(fHistDcaPosToPrimVertex);
728 fHistDcaNegToPrimVertex = new TH2F("h2DcaNegToPrimVertex", "Negative V0 daughter;dca(cm);Status",500,0,5,2,-0.5,1.5);
729 fListHist->Add(fHistDcaNegToPrimVertex);
731 fHistDcaPosToPrimVertexZoom = new TH2F("h2DcaPosToPrimVertexZoom", "Positive V0 daughter;dca(cm);Status",100,0,0.1,2,-0.5,1.5);
732 fListHist->Add(fHistDcaPosToPrimVertexZoom);
734 fHistDcaNegToPrimVertexZoom = new TH2F("h2DcaNegToPrimVertexZoom", "Negative V0 daughter;dca(cm);Status",100,0,0.1,2,-0.5,1.5);
735 fListHist->Add(fHistDcaNegToPrimVertexZoom);
737 fHistRadiusV0 = new TH2F("h2RadiusV0", "Radius;Radius(cm);Status",5000,0,500,2,-0.5,1.5);
738 fListHist->Add(fHistRadiusV0);
740 fHistDecayLengthV0 = new TH2F("h2DecayLengthV0", "V0s decay Length;decay length(cm);Status", 240, 0, 120,2,-0.5,1.5);
741 fListHist->Add(fHistDecayLengthV0);
743 fHistDcaV0Daughters = new TH2F("h2DcaV0Daughters", "DCA between daughters;dca(cm);Status", 400, 0, 4,2,-0.5,1.5);
744 fListHist->Add(fHistDcaV0Daughters);
746 fHistChi2 = new TH2F("h2Chi2", "V0s chi2;chi2;Status", 33, 0, 33,2,-0.5,1.5);
747 fListHist->Add(fHistChi2);
749 fHistCosPointAngle = new TH2F("h2CosPointAngle", "Cosine of V0's pointing angle", 100,0,1,2,-0.5,1.5);
750 fListHist->Add(fHistCosPointAngle);
752 fHistCosPointAngleZoom = new TH2F("h2CosPointAngleZoom", "Cosine of V0's pointing angle", 1000,0.9,1,2,-0.5,1.5);
753 fListHist->Add(fHistCosPointAngleZoom);
755 fHistProdRadius = new TH2F("h2ProdRadius", "Production position;x (cm);y (cm)", 100,-50,50,100,-50,50);
756 fListHist->Add(fHistProdRadius);
759 if (!fHistV0Multiplicity) {
760 if (fCollidingSystems)
761 fHistV0Multiplicity = new TH1F("fHistV0Multiplicity", "Multiplicity distribution;Number of Offline V0s;Events", 200, 0, 40000);
763 fHistV0Multiplicity = new TH1F("fHistV0Multiplicity", "Multiplicity distribution;Number of Offline V0s;Events", 10, 0, 10);
764 fListHist->Add(fHistV0Multiplicity);
767 // Kalman Filter Chi2:
768 fHistChi2KFBeforeCutK0s = new TH2F("h1Chi2KFBeforeCutK0s", "K^{0} candidates;#Chi^{2});Counts", 250, 0, 50, 2,-0.5,1.5);
769 fListHist->Add(fHistChi2KFBeforeCutK0s);
770 fHistChi2KFBeforeCutLambda = new TH2F("h1Chi2KFBeforeCutLambda", "#Lambda^{0} candidates;#Chi^{2};Counts", 250, 0, 50, 2,-0.5,1.5);
771 fListHist->Add(fHistChi2KFBeforeCutLambda);
772 fHistChi2KFBeforeCutAntiLambda = new TH2F("h1Chi2KFBeforeCutAntiLambda", "#bar{#Lambda}^{0} candidates;#Chi^{2};Counts", 250, 0, 50, 2,-0.5,1.5);
773 fListHist->Add(fHistChi2KFBeforeCutAntiLambda);
775 fHistChi2KFAfterCutK0s = new TH2F("h1Chi2KFAfterCutK0s", "K^{0} candidates;#Chi^{2});Counts", 250, 0, 50, 2,-0.5,1.5);
776 fListHist->Add(fHistChi2KFAfterCutK0s);
777 fHistChi2KFAfterCutLambda = new TH2F("h1Chi2KFAfterCutLambda", "#Lambda^{0} candidates;#Chi^{2};Counts", 250, 0, 50, 2,-0.5,1.5);
778 fListHist->Add(fHistChi2KFAfterCutLambda);
779 fHistChi2KFAfterCutAntiLambda = new TH2F("h1Chi2KFAfterCutAntiLambda", "#bar{#Lambda}^{0} candidates;#Chi^{2};Counts", 250, 0, 50, 2,-0.5,1.5);
780 fListHist->Add(fHistChi2KFAfterCutAntiLambda);
783 fHistMassK0 = new TH1F("h1MassK0", "K^{0} candidates;M(#pi^{+}#pi^{-}) (GeV/c^{2});Counts", 100, 0.4, 0.6);
784 fListHist->Add(fHistMassK0);
786 fHistMassLambda = new TH1F("h1MassLambda", "#Lambda^{0} candidates;M(p#pi^{-}) (GeV/c^{2});Counts", 75, 1.05, 1.2);
787 fListHist->Add(fHistMassLambda);
789 fHistMassAntiLambda = new TH1F("h1MassAntiLambda", "#bar{#Lambda}^{0} candidates;M(#bar{p}#pi^{+}) (GeV/c^{2});Counts", 75, 1.05, 1.2);
790 fListHist->Add(fHistMassAntiLambda);
792 // Invariant mass vs. radius:
793 const Double_t radius[10] = {0.0,2.5,2.9,3.9,7.6,15.0,23.9,37.8,42.8,100.0};
794 Int_t lNbinRadius = 9;
795 Int_t lNbinInvMassLambda = 300;
797 fHistMassVsRadiusK0 = new TH2F("h2MassVsRadiusK0", "K^{0} candidates;radius (cm);M(#pi^{+}#pi^{-}) (GeV/c^{2})",lNbinRadius,radius, 200, 0.4, 0.6);
798 fListHist->Add(fHistMassVsRadiusK0);
800 fHistMassVsRadiusLambda = new TH2F("h2MassVsRadiusLambda", "#Lambda candidates;radius (cm);M(p#pi^{-}) (GeV/c^{2})",lNbinRadius,radius, 140, 1.06, 1.2);
801 fListHist->Add(fHistMassVsRadiusLambda);
803 fHistMassVsRadiusAntiLambda = new TH2F("h2MassVsRadiusAntiLambda", "#bar{#Lambda} candidates;radius (cm);M(#bar{p}#pi^{+}) (GeV/c^{2})",lNbinRadius,radius, 140, 1.06, 1.2);
804 fListHist->Add(fHistMassVsRadiusAntiLambda);
807 fHistPtVsMassK0 = new TH2F("h2PtVsMassK0","K^{0} candidates;M(#pi^{+}#pi^{-}) (GeV/c^{2});p_{t} (GeV/c)",200, 0.4, 0.6,240,0,12);
808 fListHist->Add(fHistPtVsMassK0);
810 fHistPtVsMassLambda = new TH2F("h2PtVsMassLambda","#Lambda^{0} candidates;M(p#pi^{-}) (GeV/c^{2});p_{t} (GeV/c)",140, 1.06, 1.2,240,0,12);
811 fListHist->Add(fHistPtVsMassLambda);
813 fHistArmenterosPodolanski = new TH2F("h2ArmenterosPodolanski","Armenteros-Podolanski phase space;#alpha;p_{t} arm",100,-1.0,1.0,50,0,0.5);
814 fListHist->Add(fHistArmenterosPodolanski);
818 fHistNsigmaPosPionAntiLambda = new TH1F("h1NsigmaPosPionAntiLambda", "Positive daughter of Antilambda;NsigmaPion;Counts",25,0,5);
819 fListHist->Add(fHistNsigmaPosPionAntiLambda);
821 fHistNsigmaNegProtonAntiLambda = new TH1F("h1NsigmaNegProtonAntiLambda", "Negative daughter of Antilambda;NsigmaProton;Counts",25,0,5);
822 fListHist->Add(fHistNsigmaNegProtonAntiLambda);
824 fHistNsigmaPosProtonLambda = new TH1F("h1NsigmaPosProtonLambda", "Positive daughter of Lambda;NsigmaProton;Counts",25,0,5);
825 fListHist->Add(fHistNsigmaPosProtonLambda);
827 fHistNsigmaNegPionLambda = new TH1F("h1NsigmaNegPionLambda", "Negative daughter of Lambda;NsigmaPion;Counts",25,0,5);
828 fListHist->Add(fHistNsigmaNegPionLambda);
830 fHistNsigmaPosPionK0 = new TH1F("h1NsigmaPosPionK0", "Positive daughter of K0s;NsigmaPion;Counts",25,0,5);
831 fListHist->Add(fHistNsigmaPosPionK0);
833 fHistNsigmaNegPionK0 = new TH1F("h1NsigmaNegPionK0", "Negative daughter of K0s;NsigmaPion;Counts",25,0,5);
834 fListHist->Add(fHistNsigmaNegPionK0);
836 //********************************
837 // Associated particle histograms
838 //********************************
840 // Rapidity distribution:
841 fHistAsMcRapK0 = new TH1F("h1AsMcRapK0", "K^{0} associated;eta;Counts", 60, -1.5, 1.5);
842 fListHist->Add(fHistAsMcRapK0);
844 fHistAsMcRapLambda = new TH1F("h1AsMcRapLambda", "#Lambda^{0} associated;eta;Counts", 60, -1.5, 1.5);
845 fListHist->Add(fHistAsMcRapLambda);
847 fHistAsMcRapAntiLambda = new TH1F("h1AsMcRapAntiLambda", "#bar{#Lambda}^{0} associated;eta;Counts", 60, -1.5, 1.5);
848 fListHist->Add(fHistAsMcRapAntiLambda);
851 fHistAsMcPtK0 = new TH1F("h1AsMcPtK0", "K^{0} associated;p_{t} (GeV/c);Counts", 240,0,12);
852 fListHist->Add(fHistAsMcPtK0);
854 fHistAsMcPtLambda = new TH1F("h1AsMcPtLambda", "#Lambda^{0} associated;p_{t} (GeV/c);Counts", 240,0,12);
855 fListHist->Add(fHistAsMcPtLambda);
857 fHistAsMcPtZoomK0 = new TH1F("h1AsMcPtZoomK0", "K^{0} candidates in -1 <y<1;p_{t} (GeV/c);Counts",20,0,1);
858 fListHist->Add(fHistAsMcPtZoomK0);
860 fHistAsMcPtZoomLambda = new TH1F("h1AsMcPtZoomLambda", "#Lambda^{0} candidates in -1 <y<1;p_{t} (GeV/c);Counts",20,0,1);
861 fListHist->Add(fHistAsMcPtZoomLambda);
863 // Radius distribution:
864 fHistAsMcProdRadiusK0 = new TH1F("h1AsMcProdRadiusK0", "K^{0} associated;r (cm);Counts", 500, 0, 100);
865 fListHist->Add(fHistAsMcProdRadiusK0);
867 fHistAsMcProdRadiusLambda = new TH1F("h1AsMcProdRadiusLambda", "#Lambda^{0} associated;r (cm);Counts", 500, 0, 100);
868 fListHist->Add(fHistAsMcProdRadiusLambda);
870 fHistAsMcProdRadiusAntiLambda = new TH1F("h1AsMcProdRadiusAntiLambda", "#bar{#Lambda}^{0} associated;r (cm);Counts", 500, 0, 100);
871 fListHist->Add(fHistAsMcProdRadiusAntiLambda);
873 // Radius distribution vs. rapidity:
874 fHistAsMcProdRadiusXvsYK0s = new TH2F("h2AsMcProdRadiusXvsYK0s","Associated Secondary K^{0} Production Radius;x (cm); y (cm)",200,-50,50,200,-50,50);
875 fListHist->Add(fHistAsMcProdRadiusXvsYK0s);
877 fHistAsMcProdRadiusXvsYLambda = new TH2F("h2AsMcProdRadiusXvsYLambda","Associated Secondary #Lambda^{0} Production Radius;x (cm); y (cm)",200,-50,50,200,-50,50);
878 fListHist->Add(fHistAsMcProdRadiusXvsYLambda);
880 fHistAsMcProdRadiusXvsYAntiLambda = new TH2F("h2AsMcProdRadiusXvsYAntiLambda","Associated Secondary #bar{#Lambda}^{0} Production Radius;x (cm); y (cm)",200,-50,50,200,-50,50);
881 fListHist->Add(fHistAsMcProdRadiusXvsYAntiLambda);
883 // Invariant mass distribution with PID checked:
884 fHistPidMcMassK0 = new TH1F("h1PidMcMassK0", "K^{0} MC PId checked;M(#pi^{+}#pi^{-}) (GeV/c^{2});Counts", 100, 0.4, 0.6);
885 fListHist->Add(fHistPidMcMassK0);
887 fHistPidMcMassLambda = new TH1F("h1PidMcMassLambda", "#Lambda^{0} MC PId checked;M(p#pi^{-}) (GeV/c^{2});Counts", 75, 1.05, 1.2);
888 fListHist->Add(fHistPidMcMassLambda);
890 fHistPidMcMassAntiLambda = new TH1F("h1PidMcMassAntiLambda", "#bar{#Lambda}^{0} MC PId checked;M(#bar{p}#pi^{+}) (GeV/c^{2});Counts", 75, 1.05, 1.2);
891 fListHist->Add(fHistPidMcMassAntiLambda);
893 // Invariant mass distribution:
894 fHistAsMcMassK0 = new TH1F("h1AsMcMassK0", "K^{0} associated;M(#pi^{+}#pi^{-}) (GeV/c^{2});Counts", 100, 0.4, 0.6);
895 fListHist->Add(fHistAsMcMassK0);
897 fHistAsMcMassLambda = new TH1F("h1AsMcMassLambda", "#Lambda^{0} associated;M(p#pi^{-}) (GeV/c^{2});Counts", 75, 1.05, 1.2);
898 fListHist->Add(fHistAsMcMassLambda);
900 fHistAsMcMassAntiLambda = new TH1F("h1AsMcMassAntiLambda", "#bar{#Lambda}^{0} associated;M(#bar{p}#pi^{+}) (GeV/c^{2});Counts", 75, 1.05, 1.2);
901 fListHist->Add(fHistAsMcMassAntiLambda);
903 // Pt vs. invariant mass:
904 fHistAsMcPtVsMassK0 = new TH2F("h2AsMcPtVsMassK0","K^{0} associated;M(#pi^{+}#pi^{-}) (GeV/c^{2});p_{t} (GeV/c)",200, 0.4, 0.6,240,0,12);
905 fListHist->Add(fHistAsMcPtVsMassK0);
907 fHistAsMcPtVsMassLambda = new TH2F("h2AsMcPtVsMassLambda","#Lambda^{0} associated;M(p#pi^{-}) (GeV/c^{2});p_{t} (GeV/c)",140, 1.06, 1.2,240,0,12);
908 fListHist->Add(fHistAsMcPtVsMassLambda);
910 fHistAsMcPtVsMassAntiLambda = new TH2F("h2AsMcPtVsMassAntiLambda","#bar{#Lambda}^{0} associated;M(#bar{p}#pi^{+}) (GeV/c^{2});p_{t} (GeV/c)",140, 1.06, 1.2,240,0,12);
911 fListHist->Add(fHistAsMcPtVsMassAntiLambda);
914 // Invariant mass vs. radius:
915 fHistAsMcMassVsRadiusK0 = new TH2F("h2AsMcMassVsRadiusK0", "K^{0} associated;radius (cm);M(#pi^{+}#pi^{-}) (GeV/c^{2})",lNbinRadius,radius, 500, 0.47, 0.52);
916 fListHist->Add(fHistAsMcMassVsRadiusK0);
918 fHistAsMcMassVsRadiusLambda = new TH2F("h2AsMcMassVsRadiusLambda", "#Lambda associated;radius (cm);M(p#pi^{-}) (GeV/c^{2})",lNbinRadius,radius, lNbinInvMassLambda, 1.10, 1.13);
919 fListHist->Add(fHistAsMcMassVsRadiusLambda);
921 fHistAsMcMassVsRadiusAntiLambda = new TH2F("h2AsMcMassVsRadiusAntiLambda", "#bar{#Lambda} associated;radius (cm);M(#bar{p}#pi^{+}) (GeV/c^{2})",lNbinRadius,radius,lNbinInvMassLambda , 1.10, 1.13);
922 fListHist->Add(fHistAsMcMassVsRadiusAntiLambda);
924 // Position resolution for K0s:
925 fHistAsMcResxK0 = new TH1F("h1AsMcResxK0", "K^{0} associated;#Delta x (cm);Counts", 50, -0.25, 0.25);
926 fListHist->Add(fHistAsMcResxK0);
927 fHistAsMcResyK0 = new TH1F("h1AsMcResyK0", "K^{0} associated;#Delta y (cm);Counts", 50, -0.25, 0.25);
928 fListHist->Add(fHistAsMcResyK0);
929 fHistAsMcReszK0 = new TH1F("h1AsMcReszK0", "K^{0} associated;#Delta z (cm);Counts", 50, -0.25, 0.25);
930 fListHist->Add(fHistAsMcReszK0);
932 // Position resolution vs. radius for K0s:
933 fHistAsMcResrVsRadiusK0 = new TH2F("h2AsMcResrVsRadiusK0", "K^{0} associated;Radius (cm);#Delta r (cm)",200,0.0,50., 50, -0.25, 0.25);
934 fListHist->Add(fHistAsMcResrVsRadiusK0);
935 fHistAsMcReszVsRadiusK0 = new TH2F("h2AsMcReszVsRadiusK0", "K^{0} associated;Radius (cm);#Delta z (cm)",200,0.0,50.0, 50, -0.25, 0.25);
936 fListHist->Add(fHistAsMcReszVsRadiusK0);
938 // Position resolution for Lambda:
939 fHistAsMcResxLambda = new TH1F("h1AsMcResxLambda", "#Lambda^{0} associated;#Delta x (cm);Counts", 50, -0.25, 0.25);
940 fListHist->Add(fHistAsMcResxLambda);
941 fHistAsMcResyLambda = new TH1F("h1AsMcResyLambda", "#Lambda^{0} associated;#Delta y (cm);Counts", 50, -0.25, 0.25);
942 fListHist->Add(fHistAsMcResyLambda);
943 fHistAsMcReszLambda = new TH1F("h1AsMcReszLambda", "#Lambda^{0} associated;#Delta z (cm);Counts", 50, -0.25, 0.25);
944 fListHist->Add(fHistAsMcReszLambda);
946 // Position resolution vs. radius for Lambda:
947 fHistAsMcResrVsRadiusLambda = new TH2F("h2AsMcResrVsRadiusLambda", "#Lambda^{0} associated;Radius (cm);#Delta r (cm)",200,0.0,50.0, 50, -0.25, 0.25);
948 fListHist->Add(fHistAsMcResrVsRadiusLambda);
949 fHistAsMcReszVsRadiusLambda = new TH2F("h2AsMcReszVsRadiusLambda", "#Lambda^{0} associated;Radius (cm);#Delta z (cm)",200,0.0,50.0, 50, -0.25, 0.25);
950 fListHist->Add(fHistAsMcReszVsRadiusLambda);
952 // Position resolution for anti-Lambda:
953 fHistAsMcResxAntiLambda = new TH1F("h1AsMcResxAntiLambda", "#bar{#Lambda}^{0} associated;#Delta x (cm);Counts", 50, -0.25, 0.25);
954 fListHist->Add(fHistAsMcResxAntiLambda);
955 fHistAsMcResyAntiLambda = new TH1F("h1AsMcResyAntiLambda", "#bar{#Lambda}^{0} associated;#Delta y (cm);Counts", 50, -0.25, 0.25);
956 fListHist->Add(fHistAsMcResyAntiLambda);
957 fHistAsMcReszAntiLambda = new TH1F("h1AsMcReszAntiLambda", "#bar{#Lambda}^{0} associated;#Delta z (cm);Counts", 50, -0.25, 0.25);
958 fListHist->Add(fHistAsMcReszAntiLambda);
960 // Position resolution vs. radius for anti-Lambda:
961 fHistAsMcResrVsRadiusAntiLambda = new TH2F("h2AsMcResrVsRadiusAntiLambda", "#bar{#Lambda}^{0} associated;Radius (cm);#Delta r (cm)",200,0.0,50.0, 50, -0.25, 0.25);
962 fListHist->Add(fHistAsMcResrVsRadiusAntiLambda);
963 fHistAsMcReszVsRadiusAntiLambda = new TH2F("h2AsMcReszVsRadiusAntiLambda", "#bar{#Lambda}^{0} associated;Radius (cm);#Delta z (cm)",200,0.0,50.0, 50, -0.25, 0.25);
964 fListHist->Add(fHistAsMcReszVsRadiusAntiLambda);
967 fHistAsMcResPtK0 = new TH1F("h1AsMcResPtK0","Pt Resolution K^{0};#Delta Pt;Counts",200,-1,1);
968 fListHist->Add(fHistAsMcResPtK0);
969 fHistAsMcResPtLambda = new TH1F("h1AsMcResPtLambda","Pt Resolution #Lambda^{0};#Delta Pt;Counts",200,-1,1);
970 fListHist->Add(fHistAsMcResPtLambda);
971 fHistAsMcResPtAntiLambda = new TH1F("h1AsMcResPtAntiLambda","Pt Resolution #bar{#Lambda}^{0};#Delta Pt;Counts",200,-1,1);
972 fListHist->Add(fHistAsMcResPtAntiLambda);
974 // Pt Resolution vs. rapidity:
975 fHistAsMcResPtVsRapK0 = new TH2F("h2AsMcResPtVsRapK0","Pt Resolution K^{0};#Delta Pt;Rap",200,-1,1,20,-1,1);
976 fListHist->Add(fHistAsMcResPtVsRapK0);
977 fHistAsMcResPtVsRapLambda = new TH2F("h2AsMcResPtVsRapLambda","Pt Resolution #Lambda^{0};#Delta Pt;Rap",200,-1,1,20,-1,1);
978 fListHist->Add(fHistAsMcResPtVsRapLambda);
979 fHistAsMcResPtVsRapAntiLambda = new TH2F("h2AsMcResPtVsRapAntiLambda","Pt Resolution #bar{#Lambda}^{0};#Delta Pt;Rap",200,-1,1,20,-1,1);
980 fListHist->Add(fHistAsMcResPtVsRapAntiLambda);
982 // Pt Resolution vs. Pt:
983 fHistAsMcResPtVsPtK0 = new TH2F("h2AsMcResPtVsPtK0","Pt Resolution K^{0};#Delta Pt;Pt",600,-0.15,0.15,240,0,12);
984 fListHist->Add(fHistAsMcResPtVsPtK0);
985 fHistAsMcResPtVsPtLambda = new TH2F("h2AsMcResPtVsPtLambda","Pt Resolution #Lambda^{0};#Delta Pt;Pt",600,-0.15,0.15,240,0,12);
986 fListHist->Add(fHistAsMcResPtVsPtLambda);
987 fHistAsMcResPtVsPtAntiLambda = new TH2F("h2AsMcResPtVsPtAntiLambda","Pt Resolution #bar{#Lambda}^{0};#Delta Pt;Pt",300,-0.15,0.15,240,0,12);
988 fListHist->Add(fHistAsMcResPtVsPtAntiLambda);
990 // Pdg code of mother particle:
991 fHistAsMcMotherPdgCodeK0s = new TH1F("h1AsMcMotherPdgCodeK0s","Mother of Associated K^{0};mother;counts",11,0,11);
992 fListHist->Add(fHistAsMcMotherPdgCodeK0s);
993 fHistAsMcMotherPdgCodeLambda = new TH1F("h1AsMcMotherPdgCodeLambda","Mother of Associated #Lambda^{0};mother;counts",11,0,11);
994 fListHist->Add(fHistAsMcMotherPdgCodeLambda);
995 fHistAsMcMotherPdgCodeAntiLambda = new TH1F("h1AsMcMotherPdgCodeAntiLambda","Mother of Associated #bar{#Lambda}^{0};mother;counts",11,0,11);
996 fListHist->Add(fHistAsMcMotherPdgCodeAntiLambda);
998 // Pt distribution of Lambda <- Sigma decay
999 fHistAsMcPtLambdaFromSigma = new TH1F("h1AsMcPtLambdaFromSigma","#Lambda}^{0} associated from Sigma;p_{t} (GeV/c);Count",240,0,12);
1000 fListHist->Add(fHistAsMcPtLambdaFromSigma);
1001 fHistAsMcPtAntiLambdaFromSigma = new TH1F("h1AsMcPtAntiLambdaFromSigma","#bar{#Lambda}^{0} associated from Sigma;p_{t} (GeV/c);Count",240,0,12);
1002 fListHist->Add(fHistAsMcPtAntiLambdaFromSigma);
1004 //*******************************************
1005 // Associated secondary particles histograms
1006 //*******************************************
1008 // Pt vs. rapidity distribution:
1009 fHistAsMcSecondaryPtVsRapK0s = new TH2F("h2AsMcSecondaryPtVsRapK0s", "K^{0} associated secondary;p_{t} (GeV/c);rapidity",240,0,12,30,-1.5,1.5);
1010 fListHist->Add(fHistAsMcSecondaryPtVsRapK0s);
1011 fHistAsMcSecondaryPtVsRapLambda = new TH2F("h2AsMcSecondaryPtVsRapLambda", "#Lambda^{0} associated secondary;p_{t} (GeV/c);rapidity",240,0,12,30,-1.5,1.5);
1012 fListHist->Add(fHistAsMcSecondaryPtVsRapLambda);
1014 fHistAsMcSecondaryPtVsRapAntiLambda = new TH2F("h2AsMcSecondaryPtVsRapAntiLambda", "#bar{#Lambda}^{0} associated secondary;p_{t} (GeV/c);rapidity",240,0,12,30,-1.5,1.5);
1015 fListHist->Add(fHistAsMcSecondaryPtVsRapAntiLambda);
1017 // Production radius
1018 fHistAsMcSecondaryProdRadiusK0s = new TH1F("h1AsMcSecondaryProdRadiusK0s", "K^{0} Production Radius;r (cm);Count", 170, -2, 15);
1019 fListHist->Add(fHistAsMcSecondaryProdRadiusK0s);
1021 fHistAsMcSecondaryProdRadiusLambda = new TH1F("h1AsMcSecondaryProdRadiusLambda", "#Lambda^{0} Production Radius;r (cm);Count", 170, -2, 15);
1022 fListHist->Add(fHistAsMcSecondaryProdRadiusLambda);
1024 fHistAsMcSecondaryProdRadiusAntiLambda = new TH1F("h1AsMcSecondaryProdRadiusAntiLambda", "#bar{#Lambda}^{0} Production Radius;r (cm);Count", 170, -2, 15);
1025 fListHist->Add(fHistAsMcSecondaryProdRadiusAntiLambda);
1027 // Production radius vs. rapidity:
1028 fHistAsMcSecondaryProdRadiusXvsYK0s = new TH2F("h2AsMcSecondaryProdRadiusXvsYK0s","Associated Secondary K^{0} Production Radius;x (cm); y (cm)",200,-20,20,200,-20,20);
1029 fListHist->Add(fHistAsMcSecondaryProdRadiusXvsYK0s);
1030 fHistAsMcSecondaryProdRadiusXvsYLambda = new TH2F("h2AsMcSecondaryProdRadiusXvsYLambda","Associated Secondary #Lambda^{0} Production Radius;x (cm); y (cm)",200,-20,20,200,-20,20);
1031 fListHist->Add(fHistAsMcSecondaryProdRadiusXvsYLambda);
1032 fHistAsMcSecondaryProdRadiusXvsYAntiLambda = new TH2F("h2AsMcSecondaryProdRadiusXvsYAntiLambda","Associated Secondary #bar{#Lambda}^{0} Production Radius;x (cm); y (cm)",200,-20,20,200,-20,20);
1033 fListHist->Add(fHistAsMcSecondaryProdRadiusXvsYAntiLambda);
1035 // Pdg code of mother particle for secondary V0s:
1036 fHistAsMcSecondaryMotherPdgCodeK0s = new TH1F("h1AsMcSecondaryMotherPdgCodeK0s","Mother of Associated Secondary K^{0};mother;counts",11,0,11);
1037 fListHist->Add(fHistAsMcSecondaryMotherPdgCodeK0s);
1038 fHistAsMcSecondaryMotherPdgCodeLambda = new TH1F("h1AsMcSecondaryMotherPdgCodeLambda","Mother of Associated Secondary #Lambda^{0};mother;counts",11,0,11);
1039 fListHist->Add(fHistAsMcSecondaryMotherPdgCodeLambda);
1040 fHistAsMcSecondaryMotherPdgCodeAntiLambda = new TH1F("h1AsMcSecondaryMotherPdgCodeAntiLambda","Mother of Associated Secondary #bar{#Lambda}^{0};mother;counts",11,0,11);
1041 fListHist->Add(fHistAsMcSecondaryMotherPdgCodeAntiLambda);
1043 // Pt distribution of secondary Lambda <- Sigma decay:
1044 fHistAsMcSecondaryPtLambdaFromSigma = new TH1F("h1AsMcSecondaryPtLambdaFromSigma","#Lambda}^{0} associated from Sigma;p_{t} (GeV/c);Count",240,0,12);
1045 fListHist->Add(fHistAsMcSecondaryPtLambdaFromSigma);
1046 fHistAsMcSecondaryPtAntiLambdaFromSigma = new TH1F("h1AsMcSecondaryPtAntiLambdaFromSigma","#bar{#Lambda}^{0} associated from Sigma;p_{t} (GeV/c);Count",240,0,12);
1047 fListHist->Add(fHistAsMcSecondaryPtAntiLambdaFromSigma);
1049 PostData(1, fListHist);
1050 PostData(2, fCentrSelector);
1053 //________________________________________________________________________
1054 void AliAnalysisTaskPerformanceStrange::UserExec(Option_t *)
1057 // Called for each event
1059 AliStack* stack = NULL;
1060 TClonesArray *mcArray = NULL;
1061 TArrayF mcPrimaryVtx;
1063 fESD=(AliESDEvent *)InputEvent();
1066 Printf("ERROR: fESD not available");
1070 // FIXME: levent not used. Can I remove it?
1071 AliVEvent* lEvent = InputEvent();
1074 Printf("ERROR: Event not available");
1078 fHistNumberEvents->Fill(0.5);
1080 //******************
1081 // Trigger Selection ! Warning Works only for ESD, add protection in case of AOD loop
1082 //******************
1085 (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected()
1087 if (!isSelected) return;
1089 // Centrality selection
1090 static AliESDtrackCuts * trackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(); // FIXME: make it a data member
1091 Bool_t isCentralitySelected = fCentrSelector->IsCentralityBinSelected(fESD,trackCuts);
1092 if(!isCentralitySelected) return;
1093 // FIXME: add to hist number events another entry for centrality.
1095 // Done by the AliPhysicsSelection Task ! Only the selected events are passed to this task
1097 fHistNumberEvents->Fill(1.5); // FIXME: use enum here
1100 //*************************
1101 //End track multiplicity
1102 //*************************
1105 // Remove Events with no tracks
1106 //if (!(fESD->GetNumberOfTracks())) return;
1108 fHistNumberEvents->Fill(2.5);
1109 fHistTrackPerEvent->Fill(fESD->GetNumberOfTracks());
1111 //*************************************
1113 //*************************************
1114 // FIXME: Create a cut object, to be configured in the steering macro and to be streamed in the output to reference those cuts
1116 Double_t lCutRap = 0.75;
1118 // Cut AliKF Chi2 for Reconstructed particles
1119 Double_t cutChi2KF = 1E3;
1122 Double_t lLimitPPID = 0.7;
1123 Float_t cutNSigmaLowP = 1E3;
1124 Float_t cutNSigmaHighP = 1E3;
1125 if (fUsePID.Contains("withPID")) {
1126 cutNSigmaLowP = 5.0;
1127 cutNSigmaHighP = 3.0;
1131 // Cut Daughters pt (GeV/c):
1132 Double_t cutMinPtDaughter = 0.160;
1134 // Cut primary vertex:
1135 Double_t cutPrimVertex = 10.0;
1137 // Min number of TPC clusters:
1138 Int_t nbMinTPCclusters = 80;
1140 //*******************
1142 //*******************
1143 // FIXME: OADB or momber TFormula?
1144 Double_t fAlephParameters[5] = {0,0,0,0,0,};
1146 fAlephParameters[0] = 0.0283086;
1147 fAlephParameters[1] = 2.63394e+01;
1148 fAlephParameters[2] = 5.04114e-11;
1149 fAlephParameters[3] = 2.12543e+00;
1150 fAlephParameters[4] = 4.88663e+00;
1153 //*******************
1155 //*******************
1156 // FIXME:: move this two branches directly in the loops below
1158 if(fAnalysisType == "ESD") {
1159 AliMCEventHandler* eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
1160 if (!eventHandler) {
1161 Printf("ERROR: Could not retrieve MC event handler");
1164 AliMCEvent* mcEvent = eventHandler->MCEvent();
1166 Printf("ERROR: Could not retrieve MC event");
1169 stack = mcEvent->Stack();
1171 Printf("ERROR: Could not retrieve stack");
1175 AliGenEventHeader* mcHeader=mcEvent->GenEventHeader();
1176 if(!mcHeader) return;
1177 mcHeader->PrimaryVertex(mcPrimaryVtx);
1181 else if(fAnalysisType == "AOD") {
1183 // load MC particles
1184 mcArray = (TClonesArray*)fESD->GetList()->FindObject(AliAODMCParticle::StdBranchName());
1186 Printf("strange analysis::UserExec: MC particles branch not found!\n");
1191 AliAODMCHeader *mcHeader =
1192 (AliAODMCHeader*)fESD->GetList()->FindObject(AliAODMCHeader::StdBranchName());
1194 Printf("strange analysis::UserExec: MC header branch not found!\n");
1199 // PID parameters for MC simulations:
1200 // FIXME: set above, with the others
1201 fAlephParameters[0] = 2.15898e+00/50.;
1202 fAlephParameters[1] = 1.75295e+01;
1203 fAlephParameters[2] = 3.40030e-09;
1204 fAlephParameters[3] = 1.96178e+00;
1205 fAlephParameters[4] = 3.91720e+00;
1209 //**********************************************
1211 //**********************************************
1213 Double_t lmcPrimVtxR = 0;
1215 Int_t lNbMCPrimary = 0;
1216 Int_t lNbMCPart = 0;
1218 Int_t lPdgcodeCurrentPart = 0;
1219 Double_t lRapCurrentPart = 0;
1220 Double_t lPtCurrentPart = 0;
1222 Int_t lComeFromSigma = 0;
1225 // Production Radius
1226 Double_t mcPosX = 0.0, mcPosY = 0.0, mcPosZ = 0.0;
1227 Double_t mcPosR = 0.0;
1230 Double_t mcDecayPosX = 0, mcDecayPosY = 0, mcDecayPosR = 0;
1232 // current mc particle 's mother
1233 Int_t iCurrentMother = 0, lPdgCurrentMother = 0;
1234 Bool_t lCurrentMotherIsPrimary;
1236 // current mc particles 's daughter:
1237 Int_t lPdgCurrentDaughter0 = 0, lPdgCurrentDaughter1 = 0;
1239 // variables for multiple reconstruction studies:
1240 Int_t id0 = 0, id1 = 0;
1241 //Int_t lLabelTrackN = 0, lLabelTrackP = 0;
1242 //Int_t lPartNMother = 0, lPartPMother = 0;
1243 //Int_t lPartPMotherPDGcode = 0;
1244 Int_t lNtimesReconstructedK0s = 0, lNtimesReconstructedLambda = 0, lNtimesReconstructedAntiLambda = 0;
1245 // Int_t lNtimesReconstructedK0sMI = 0, lNtimesReconstructedLambdaMI = 0, lNtimesReconstructedAntiLambdaMI = 0;
1247 //****************************
1248 // Start loop over MC particles
1252 fHistMCPrimaryVertexX->Fill(mcPrimaryVtx.At(0));
1253 fHistMCPrimaryVertexY->Fill(mcPrimaryVtx.At(1));
1254 fHistMCPrimaryVertexZ->Fill(mcPrimaryVtx.At(2));
1256 lmcPrimVtxR = TMath::Sqrt(mcPrimaryVtx.At(0)*mcPrimaryVtx.At(0)+mcPrimaryVtx.At(1)*mcPrimaryVtx.At(1));
1258 // FIXME: move these loops to other functions, for better readibility?
1259 if(fAnalysisType == "ESD") {
1261 lNbMCPrimary = stack->GetNprimary(); // FIXME: This does not correspond to our definition of primaries, but maybe it is ok for strange particles
1262 lNbMCPart = stack->GetNtrack();
1264 fHistMCMultiplicityPrimary->Fill(lNbMCPrimary);
1265 fHistMCMultiplicityTracks->Fill(lNbMCPart);
1268 for (Int_t iMc = 0; iMc < (stack->GetNtrack()); iMc++) {
1269 TParticle *p0 = stack->Particle(iMc);
1271 //Printf("ERROR: particle with label %d not found in stack (mc loop)", iMc);
1274 lPdgcodeCurrentPart = p0->GetPdgCode();
1276 // Keep only K0s, Lambda and AntiLambda, Xi and Phi:
1277 if ( (lPdgcodeCurrentPart != 310 ) && (lPdgcodeCurrentPart != 3122 ) && (lPdgcodeCurrentPart != -3122 ) && (lPdgcodeCurrentPart != 3312 ) && (lPdgcodeCurrentPart != -3312) && (lPdgcodeCurrentPart != -333) ) continue;
1279 lRapCurrentPart = MyRapidity(p0->Energy(),p0->Pz());
1280 //lEtaCurrentPart = p0->Eta();
1281 lPtCurrentPart = p0->Pt();
1283 iCurrentMother = p0->GetFirstMother();
1285 // lPdgCurrentMother = stack->Particle(iCurrentMother)->GetPdgCode();
1286 if (iCurrentMother == -1){lPdgCurrentMother=0; } else {lPdgCurrentMother = stack->Particle(iCurrentMother)->GetPdgCode();}
1291 mcPosR = TMath::Sqrt(mcPosX*mcPosX+mcPosY*mcPosY);
1293 id0 = p0->GetDaughter(0);
1294 id1 = p0->GetDaughter(1);
1297 if ( id0 <= lNbMCPart && id0 > 0 && id1 <= lNbMCPart && id1 > 0) {
1298 TParticle *pDaughter0 = stack->Particle(id0);
1299 TParticle *pDaughter1 = stack->Particle(id1);
1300 lPdgCurrentDaughter0 = pDaughter0->GetPdgCode();
1301 lPdgCurrentDaughter1 = pDaughter1->GetPdgCode();
1303 mcDecayPosX = pDaughter0->Vx();
1304 mcDecayPosY = pDaughter0->Vy();
1305 mcDecayPosR = TMath::Sqrt(mcDecayPosX*mcDecayPosX+mcDecayPosY*mcDecayPosY);
1308 //FIXME: shouldn't this be a fatal?
1309 //Printf("ERROR: particle with label %d and/or %d not found in stack (mc loop)", id0,id1);
1312 // FIXME using array of histos and conversion PDGCode -> enum would make this much easier to read
1313 // We could also have a function FillMcHistos (pos, radius, rap...) which we call from both the AOD and ESD loops
1314 if (lPdgcodeCurrentPart==310) {
1315 fHistMCtracksProdRadiusK0s->Fill(mcPosX,mcPosY);
1316 fHistMCtracksDecayRadiusK0s->Fill(mcDecayPosR);
1317 if (TMath::Abs(lRapCurrentPart) < lCutRap) fHistMCPtAllK0s->Fill(lPtCurrentPart);
1319 else if (lPdgcodeCurrentPart==3122) {
1320 fHistMCtracksProdRadiusLambda->Fill(mcPosX,mcPosY);
1321 fHistMCtracksDecayRadiusLambda->Fill(mcDecayPosR);
1322 if (TMath::Abs(lRapCurrentPart) < lCutRap) fHistMCPtAllLambda->Fill(lPtCurrentPart);
1324 else if (lPdgcodeCurrentPart==-3122) {
1325 fHistMCtracksProdRadiusAntiLambda->Fill(mcPosX,mcPosY);
1326 fHistMCtracksDecayRadiusAntiLambda->Fill(mcDecayPosR);
1327 if (TMath::Abs(lRapCurrentPart) < lCutRap) fHistMCPtAllAntiLambda->Fill(lPtCurrentPart);
1330 // FIXME: not sure if I understand this: is it correct? (definition of primaries)
1331 if ( ( ( TMath::Abs(lPdgCurrentMother) == 3212) ||
1332 ( TMath::Abs(lPdgCurrentMother) == 3224) ||
1333 ( TMath::Abs(lPdgCurrentMother) == 3214) ||
1334 ( TMath::Abs(lPdgCurrentMother) == 3114) )
1335 && ( iCurrentMother <= lNbMCPrimary )
1336 ) lComeFromSigma = 1;
1337 else lComeFromSigma = 0;
1339 //*********************************************
1340 // Now keep only primary particles
1341 if ( ( iMc > lNbMCPrimary ) && (!lComeFromSigma) ) continue;
1343 //********************************************
1344 //check if V0 is reconstructed several times
1346 lNtimesReconstructedK0s = 0; lNtimesReconstructedLambda = 0; lNtimesReconstructedAntiLambda = 0;
1348 //for (Int_t jV0 = 0; jV0 < fESD->GetNumberOfV0s(); jV0++) {
1350 //lLabelTrackN = 0; lLabelTrackP = 0;
1351 //lPartNMother = 0; lPartPMother = 0;
1353 //AliESDv0 *vertexESD = ((AliESDEvent*)fESD)->GetV0(jV0);
1354 //if (!vertexESD) continue;
1356 //AliESDtrack *trackNESD = ((AliESDEvent*)fESD)->GetTrack(TMath::Abs(vertexESD->GetNindex()));
1357 //lLabelTrackN = (UInt_t)TMath::Abs(trackNESD->GetLabel());
1358 //if (lLabelTrackN!=id0 && lLabelTrackN!=id1) continue;
1360 //AliESDtrack *trackPESD = ((AliESDEvent*)fESD)->GetTrack(TMath::Abs(vertexESD->GetPindex()));
1361 //lLabelTrackP = (UInt_t)TMath::Abs(trackPESD->GetLabel());
1362 //if (lLabelTrackP!=id0 && lLabelTrackP!=id1) continue;
1364 //TParticle *lPartNESD = stack->Particle(lLabelTrackN);
1365 //TParticle *lPartPESD = stack->Particle(lLabelTrackP);
1366 //lPartNMother = lPartNESD->GetFirstMother();
1367 //lPartPMother = lPartPESD->GetFirstMother();
1369 //lPartPMotherPDGcode = stack->Particle(lPartPMother)->GetPdgCode();
1371 //switch (vertexESD->GetOnFlyStatus()){
1374 //if ( (lPartPMother==lPartNMother) && (lPartPMotherPDGcode==310) ) lNtimesReconstructedK0s++;
1375 //else if ( (lPartPMother==lPartNMother) && (lPartPMotherPDGcode==3122) ) lNtimesReconstructedLambda++;
1376 //else if ( (lPartPMother==lPartNMother) && (lPartPMotherPDGcode==-3122) ) lNtimesReconstructedAntiLambda++;
1380 //if ( (lPartPMother==lPartNMother) && (lPartPMotherPDGcode==310) ) lNtimesReconstructedK0sMI++;
1381 //else if ( (lPartPMother==lPartNMother) && (lPartPMotherPDGcode==3122) ) lNtimesReconstructedLambdaMI++;
1382 //else if ( (lPartPMother==lPartNMother) && (lPartPMotherPDGcode==-3122) ) lNtimesReconstructedAntiLambdaMI++;
1386 //} // end loop over reconstructed V0s inside MC loop
1388 // FIXME: same comemtn for array of histos
1390 if (lPdgcodeCurrentPart==310) {
1391 fHistMCRapK0s->Fill(lRapCurrentPart);
1392 if (lPtCurrentPart < 0.2 && lPtCurrentPart < 3.0)
1393 fHistMCRapInPtRangeK0s->Fill(lRapCurrentPart);
1396 if (lPdgcodeCurrentPart==3122) {
1397 fHistMCRapLambda->Fill(lRapCurrentPart);
1398 if (lPtCurrentPart < 0.6 && lPtCurrentPart < 3.5)
1399 fHistMCRapInPtRangeLambda->Fill(lRapCurrentPart);
1402 if (lPdgcodeCurrentPart==-3122) {
1403 fHistMCRapAntiLambda->Fill(lRapCurrentPart);
1404 if (lPtCurrentPart < 0.6 && lPtCurrentPart < 3.5)
1405 fHistMCRapInPtRangeAntiLambda->Fill(lRapCurrentPart);
1408 if (lPdgcodeCurrentPart==3312 || lPdgcodeCurrentPart==-3312) {
1409 fHistMCRapXi->Fill(lRapCurrentPart);
1410 if (lPtCurrentPart < 0.6 && lPtCurrentPart < 3.0)
1411 fHistMCRapInPtRangeXi->Fill(lRapCurrentPart);
1414 if (lPdgcodeCurrentPart==333) {
1415 fHistMCRapPhi->Fill(lRapCurrentPart);
1416 if (lPtCurrentPart < 0.7 && lPtCurrentPart < 3.0)
1417 fHistMCRapInPtRangePhi->Fill(lRapCurrentPart);
1421 if (TMath::Abs(lRapCurrentPart) > lCutRap) continue;
1423 if (lPdgcodeCurrentPart==310) {
1424 fHistMCProdRadiusK0s->Fill(mcPosR);
1426 fHistMCPtK0s->Fill(lPtCurrentPart);
1428 fHistNTimesRecK0s->Fill(lNtimesReconstructedK0s);
1429 fHistNTimesRecK0sVsPt->Fill(lPtCurrentPart,lNtimesReconstructedK0s);
1432 if (lPdgcodeCurrentPart==3122) {
1433 fHistMCProdRadiusLambda->Fill(mcPosR);
1435 fHistMCPtLambda->Fill(lPtCurrentPart);
1438 fHistNTimesRecLambda->Fill(lNtimesReconstructedLambda);
1439 fHistNTimesRecLambdaVsPt->Fill(lPtCurrentPart,lNtimesReconstructedLambda);
1440 if (lComeFromSigma) fHistMCPtLambdaFromSigma->Fill(lPtCurrentPart);
1441 //printf("found Lambda MC pT=%e\n",lPtCurrentPart);
1442 //printf("found Lambda MC Plabel=%d PPDGcode=%d Nlabel=%d NPDGcode=%d\n\n",id0,lPdgCurrentDaughter0,id1,lPdgCurrentDaughter1);
1446 } // end loop ESD MC
1448 } // end ESD condition
1450 // FIXME: I skipped the AOD loop
1451 else if(fAnalysisType == "AOD") {
1452 lNbMCPart = mcArray->GetEntriesFast();
1455 fHistMCMultiplicityTracks->Fill(lNbMCPart);
1457 for (Int_t iMc = 0; iMc < lNbMCPart; iMc++) {
1459 // Primary vertex TO DO !!
1462 AliAODMCParticle *mcAODPart = (AliAODMCParticle*)mcArray->At(iMc);
1464 //Printf("Strange analysis task (mc loop): particle with label %d not found", iMc);
1467 lPdgcodeCurrentPart = mcAODPart->GetPdgCode();
1468 if (mcAODPart->IsPhysicalPrimary()) {lNbMCPrimary = lNbMCPrimary +1;}
1470 // Keep only K0s, Lambda and AntiLambda:
1471 if ( (lPdgcodeCurrentPart != 310 ) && (lPdgcodeCurrentPart != 3122 ) && (lPdgcodeCurrentPart != -3122 ) ) continue;
1473 //lEtaCurrentPart = mcAODPart->Eta();
1474 lRapCurrentPart = mcAODPart->Y();
1475 lPtCurrentPart = mcAODPart->Pt();
1476 iCurrentMother = mcAODPart->GetMother();
1477 lPdgCurrentMother = ((AliAODMCParticle*)mcArray->At(iCurrentMother))->GetPdgCode();
1478 lCurrentMotherIsPrimary = ((AliAODMCParticle*)mcArray->At(iCurrentMother))->IsPhysicalPrimary();
1480 mcPosX = mcAODPart->Xv();
1481 mcPosY = mcAODPart->Yv();
1482 mcPosZ = mcAODPart->Zv();
1483 mcPosR = TMath::Sqrt(mcPosX*mcPosX+mcPosY*mcPosY);
1485 id0 = mcAODPart->GetDaughter(0);
1486 id1 = mcAODPart->GetDaughter(1);
1488 // Decay Radius and Production Radius
1489 if ( id0 <= lNbMCPart && id0 > 0 && id1 <= lNbMCPart && id1 > 0) {
1490 AliAODMCParticle *mcAODDaughter1 = (AliAODMCParticle*)mcArray->At(id1);
1492 //Printf("Strange analysis task (mc loop): daughter not found");
1495 mcDecayPosX = mcAODDaughter1->Xv();
1496 mcDecayPosY = mcAODDaughter1->Yv();
1497 mcDecayPosR = TMath::Sqrt(mcDecayPosX*mcDecayPosX+mcDecayPosY*mcDecayPosY);
1500 //Printf("ERROR: particle with label %d and/or %d not found in stack (mc loop)", id0,id1);
1504 if (lPdgcodeCurrentPart==310) {
1505 fHistMCtracksProdRadiusK0s->Fill(mcPosX,mcPosY);
1506 fHistMCtracksDecayRadiusK0s->Fill(mcDecayPosR);
1507 if (TMath::Abs(lRapCurrentPart) < lCutRap) fHistMCPtAllK0s->Fill(lPtCurrentPart);
1509 else if (lPdgcodeCurrentPart==3122) {
1510 fHistMCtracksProdRadiusLambda->Fill(mcPosX,mcPosY);
1511 fHistMCtracksDecayRadiusLambda->Fill(mcDecayPosR);
1512 if (TMath::Abs(lRapCurrentPart) < lCutRap) fHistMCPtAllLambda->Fill(lPtCurrentPart);
1514 else if (lPdgcodeCurrentPart==-3122) {
1515 fHistMCtracksProdRadiusAntiLambda->Fill(mcPosX,mcPosY);
1516 fHistMCtracksDecayRadiusAntiLambda->Fill(mcDecayPosR);
1517 if (TMath::Abs(lRapCurrentPart) < lCutRap) fHistMCPtAllAntiLambda->Fill(lPtCurrentPart);
1520 if ( ( ( TMath::Abs(lPdgCurrentMother) == 3212) ||
1521 ( TMath::Abs(lPdgCurrentMother) == 3224) ||
1522 ( TMath::Abs(lPdgCurrentMother) == 3214) ||
1523 ( TMath::Abs(lPdgCurrentMother) == 3114) )
1524 && (lCurrentMotherIsPrimary)
1525 ) lComeFromSigma = 1;
1526 else lComeFromSigma = 0;
1528 //*********************************************
1529 // Now keep only primary particles
1531 // FIX IT !!!! iMC is not defined !!!! FIX IT also in ESD/AOD loop !!
1532 if ( ( iMc > lNbMCPrimary ) && (!lComeFromSigma) ) continue;
1534 //********************************************
1535 // check if V0 is reconstructed several times
1537 //lNtimesReconstructedK0s = 0; lNtimesReconstructedLambda = 0; lNtimesReconstructedAntiLambda = 0;
1538 //lNtimesReconstructedK0sMI = 0; lNtimesReconstructedLambdaMI = 0; lNtimesReconstructedAntiLambdaMI = 0;
1540 //for (Int_t jV0 = 0; jV0 < fESD->GetNumberOfV0s(); jV0++) {
1542 //lLabelTrackN = 0; lLabelTrackP = 0;
1543 //lPartNMother = 0; lPartPMother = 0;
1545 //AliAODv0 *vertexAOD= ((AliAODEvent*)fESD)->GetV0(jV0);
1546 //if (!vertexAOD) continue;
1547 //printf("enter!!");
1548 //AliVParticle *trackP = ((AliVEvent*)fESD)->GetTrack(vertexAOD->GetPosID());
1549 //if (!trackP) continue;
1550 //lLabelTrackP = TMath::Abs(trackP->GetLabel());
1551 //if (lLabelTrackP!=id0 && lLabelTrackP!=id1) continue;
1553 //AliVParticle *trackN = ((AliVEvent*)fESD)->GetTrack(vertexAOD->GetNegID());
1554 //if (!trackN) continue;
1555 //lLabelTrackN = TMath::Abs(trackN->GetLabel());
1556 //if (lLabelTrackN!=id0 && lLabelTrackN!=id1) continue;
1558 //AliAODMCParticle *lPartNAOD = (AliAODMCParticle*)mcArray->At(lLabelTrackN);
1559 //if (!lPartNAOD) continue;
1560 //AliAODMCParticle *lPartPAOD = (AliAODMCParticle*)mcArray->At(lLabelTrackP);
1561 //if (!lPartPAOD) continue;
1563 //lPartNMother = lPartNAOD->GetMother();
1564 //lPartPMother = lPartPAOD->GetMother();
1566 //lPartPMotherPDGcode = ((AliAODMCParticle*)mcArray->At(lPartPMother))->GetPdgCode();
1568 //switch (vertexAOD->GetOnFlyStatus()){
1571 //if ( (lPartPMother==lPartNMother) && (lPartPMotherPDGcode==310) ) lNtimesReconstructedK0s++;
1572 //else if ( (lPartPMother==lPartNMother) && (lPartPMotherPDGcode==3122) ) lNtimesReconstructedLambda++;
1573 //else if ( (lPartPMother==lPartNMother) && (lPartPMotherPDGcode==-3122) ) lNtimesReconstructedAntiLambda++;
1577 //if ( (lPartPMother==lPartNMother) && (lPartPMotherPDGcode==310) ) lNtimesReconstructedK0sMI++;
1578 //else if ( (lPartPMother==lPartNMother) && (lPartPMotherPDGcode==3122) ) lNtimesReconstructedLambdaMI++;
1579 //else if ( (lPartPMother==lPartNMother) && (lPartPMotherPDGcode==-3122) ) lNtimesReconstructedAntiLambdaMI++;
1583 //} // end loop over reconstructed V0s inside MC loop
1585 if (TMath::Abs(lRapCurrentPart) > lCutRap) continue;
1587 if (lPdgcodeCurrentPart==310) {
1588 fHistMCProdRadiusK0s->Fill(mcPosR);
1589 fHistMCPtK0s->Fill(lPtCurrentPart);
1590 fHistNTimesRecK0s->Fill(lNtimesReconstructedK0s);
1591 fHistNTimesRecK0sVsPt->Fill(lPtCurrentPart,lNtimesReconstructedK0s);
1593 else if (lPdgcodeCurrentPart==3122) {
1594 fHistMCProdRadiusLambda->Fill(mcPosR);
1595 fHistMCPtLambda->Fill(lPtCurrentPart);
1596 fHistNTimesRecLambda->Fill(lNtimesReconstructedLambda);
1597 fHistNTimesRecLambdaVsPt->Fill(lPtCurrentPart,lNtimesReconstructedLambda);
1598 if (lComeFromSigma) fHistMCPtLambdaFromSigma->Fill(lPtCurrentPart);
1600 else if (lPdgcodeCurrentPart==-3122) {
1601 fHistMCProdRadiusAntiLambda->Fill(mcPosR);
1602 fHistNTimesRecAntiLambda->Fill(lNtimesReconstructedAntiLambda);
1603 fHistNTimesRecAntiLambdaVsPt->Fill(lPtCurrentPart,lNtimesReconstructedAntiLambda);
1604 if (lComeFromSigma) fHistMCPtAntiLambdaFromSigma->Fill(lPtCurrentPart);
1607 } // end loop over AODMC particles
1608 fHistMCMultiplicityPrimary->Fill(lNbMCPrimary);
1610 } // end AOD condition
1612 } // End Loop over MC condition
1617 //************************************
1619 //************************************
1621 Double_t lMagneticField = 999;
1624 Int_t nv0sTot= 0, nv0s = 0;
1627 Double_t lV0Position[3];
1629 Double_t lDcaPosToPrimVertex = 0;
1630 Double_t lDcaNegToPrimVertex = 0;
1631 Double_t lDcaV0Daughters = 0;
1632 Double_t lV0cosPointAngle = 0;
1633 Double_t lChi2V0 = 0;
1634 Double_t lV0DecayLength = 0;
1635 Double_t lV0Radius = 0;
1636 Double_t lDcaV0ToPrimVertex = 0;
1638 Int_t lOnFlyStatus = 0;
1639 //Float_t tdcaPosToPrimVertexXYZ[2], tdcaNegToPrimVertexXYZ[2]; // ..[0] = Impact parameter in XY plane and ..[1] = Impact parameter in Z
1640 //Double_t tdcaDaughterToPrimVertex[2]; // ..[0] = Pos and ..[1] = Neg
1644 Double_t lInvMassK0s = 0, lInvMassLambda = 0, lInvMassAntiLambda = 0;
1645 Double_t lPtK0s = 0, lPtLambda = 0, lPtAntiLambda = 0;
1646 Double_t lRapK0s = 0, lRapLambda = 0, lRapAntiLambda = 0;
1647 Double_t lEtaK0s = 0, lEtaLambda = 0, lEtaAntiLambda = 0;
1648 Double_t lAlphaV0 = 0, lPtArmV0 = 0;
1650 Double_t lPzK0s = 0, lPzLambda = 0, lPzAntiLambda = 0;
1653 Double_t lV0Eta = 999;
1655 // to study Associated V0s:
1656 Int_t lIndexTrackPos = 0, lIndexTrackNeg = 0;
1657 UInt_t lLabelTrackPos = 0, lLabelTrackNeg = 0;
1658 Int_t lCheckPIdK0Short = 0, lCheckMcK0Short = 0;
1659 Int_t lCheckPIdLambda = 0, lCheckMcLambda = 0;
1660 Int_t lCheckPIdAntiLambda = 0, lCheckMcAntiLambda = 0;
1661 Int_t lCheckSecondaryK0s = 0, lCheckSecondaryLambda = 0, lCheckSecondaryAntiLambda = 0;
1662 Int_t lCheckGamma = 0;
1663 Double_t mcPosMotherX = 0, mcPosMotherY = 0, mcPosMotherZ = 0;
1664 Double_t mcPosMotherR = 0;
1665 Double_t mcMotherPt = 0;
1667 Int_t lIndexPosMother = 0;
1668 Int_t lIndexNegMother = 0;
1669 Int_t lIndexMotherOfMother = 0;
1670 Int_t lPDGCodePosDaughter = 0;
1671 Int_t lPDGCodeNegDaughter = 0;
1672 Int_t lPdgcodeMother = 0;
1673 Int_t lPdgcodeMotherOfMother = 0;
1675 // Reconstructed position
1676 Double_t rcPosXK0s = 0, rcPosYK0s = 0, rcPosZK0s = 0;
1677 Double_t rcPosRK0s = 0;
1678 Double_t rcPosXLambda = 0, rcPosYLambda = 0, rcPosZLambda = 0;
1679 Double_t rcPosRLambda = 0;
1680 Double_t rcPosXAntiLambda = 0, rcPosYAntiLambda = 0, rcPosZAntiLambda = 0;
1681 Double_t rcPosRAntiLambda = 0;
1684 Double_t deltaPtK0s = 0, deltaPtLambda = 0, deltaPtAntiLambda = 0;
1687 AliESDtrack *myTrackPos = NULL;
1688 AliESDtrack *myTrackNeg = NULL;
1689 AliVParticle *lVPartPos = NULL;
1690 AliVParticle *lVPartNeg = NULL;
1692 // Daughters' momentum:
1693 Double_t lMomPos[3] = {999,999,999};
1694 Double_t lMomNeg[3] = {999,999,999};
1695 Double_t lPtPos = 999, lPtNeg = 999;
1696 Double_t lPPos = 999, lPNeg = 999;
1698 // Inner Wall parameters:
1699 Double_t lMomInnerWallPos =999, lMomInnerWallNeg = 999;
1701 // AliKF Chi2 and Armenteros variables
1702 Double_t lChi2KFK0s = 0, lChi2KFLambda = 0, lChi2KFAntiLambda = 0;
1703 Double_t lAlphaV0K0s = 0, lAlphaV0Lambda = 0, lAlphaV0AntiLambda = 0;
1704 Double_t lPtArmV0K0s = 0, lPtArmV0Lambda = 0, lPtArmV0AntiLambda = 0;
1705 Double_t lQlPos = 0, lQlNeg = 0;
1709 Float_t nSigmaPosPion = 0;
1710 Float_t nSigmaNegPion = 0;
1712 Float_t nSigmaPosProton = 0;
1713 Float_t nSigmaNegProton = 0;
1715 Int_t lCheckPIDK0sPosDaughter = 0, lCheckPIDK0sNegDaughter = 0;
1716 Int_t lCheckPIDLambdaPosDaughter = 0, lCheckPIDLambdaNegDaughter = 0;
1717 Int_t lCheckPIDAntiLambdaPosDaughter = 0, lCheckPIDAntiLambdaNegDaughter = 0;
1721 //***********************
1722 // Primary Vertex cuts &
1723 // Magnetic field and Quality tracks cuts
1725 Double_t lPrimaryVtxPosition[3];
1726 Double_t lPrimaryVtxCov[6];
1727 Double_t lPrimaryVtxChi2 = 999;
1728 Double_t lResPrimaryVtxX = 999;
1729 Double_t lResPrimaryVtxY = 999;
1730 Double_t lResPrimaryVtxZ = 999;
1732 AliAODVertex *myPrimaryVertex = NULL;
1733 //const AliVVertex *mySPDPrimaryVertex = NULL;
1735 AliESDtrackCuts *myTracksCuts = NULL;
1737 const AliMultiplicity *myMultiplicty = ((AliESDEvent*)fESD)->GetMultiplicity();
1739 if(fAnalysisType == "ESD") {
1741 // Best Primary Vertex:
1742 const AliESDVertex *myBestPrimaryVertex = ((AliESDEvent*)fESD)->GetPrimaryVertex();
1743 myBestPrimaryVertex = ((AliESDEvent*)fESD)->GetPrimaryVertex();
1744 if (!myBestPrimaryVertex) return;
1745 if (!myBestPrimaryVertex->GetStatus()) return;
1746 fHistNumberEvents->Fill(3.5);
1747 myBestPrimaryVertex->GetXYZ(lPrimaryVtxPosition);
1748 myBestPrimaryVertex->GetCovMatrix(lPrimaryVtxCov);
1749 if ( ( TMath::Abs(lPrimaryVtxPosition[2]) ) > cutPrimVertex) return ;
1750 fHistNumberEvents->Fill(4.5);
1751 lPrimaryVtxChi2 = myBestPrimaryVertex->GetChi2toNDF();
1752 lResPrimaryVtxX = myBestPrimaryVertex->GetXRes();
1753 lResPrimaryVtxY = myBestPrimaryVertex->GetYRes();
1754 lResPrimaryVtxZ = myBestPrimaryVertex->GetZRes();
1756 // remove TPC-only primary vertex : retain only events with tracking + SPD vertex
1757 const AliESDVertex *mySPDPrimaryVertex = ((AliESDEvent*)fESD)->GetPrimaryVertexSPD();
1758 if (!mySPDPrimaryVertex) return;
1759 fHistSPDPrimaryVertexZ->Fill(mySPDPrimaryVertex->GetZ());
1760 const AliESDVertex *myPrimaryVertexTracking = ((AliESDEvent*)fESD)->GetPrimaryVertexTracks();
1761 if (!myPrimaryVertexTracking) return;
1762 if (!mySPDPrimaryVertex->GetStatus() && !myPrimaryVertexTracking->GetStatus() ) return;
1763 fHistNumberEvents->Fill(5.5);
1766 myPrimaryVertex = new AliAODVertex(lPrimaryVtxPosition, lPrimaryVtxCov, lPrimaryVtxChi2, NULL, -1, AliAODVertex::kPrimary);
1767 if (!myPrimaryVertex) return;
1770 // Number of Tracklets:
1771 //const AliMultiplicity *myMultiplicty = ((AliESDEvent*)fESD)->GetMultiplicity();
1772 //if (myMultiplicty->GetNumberOfTracklets() < 10) return;
1773 fHistTrackletPerEvent->Fill(myMultiplicty->GetNumberOfTracklets());
1775 lMagneticField = ((AliESDEvent*)fESD)->GetMagneticField();
1777 myTracksCuts = new AliESDtrackCuts();
1778 // require TPC refit
1779 myTracksCuts->SetRequireTPCRefit(kTRUE);
1780 // minimum number of clusters in TPC
1781 myTracksCuts->SetMinNClustersTPC(nbMinTPCclusters);
1785 else if(fAnalysisType == "AOD") {
1786 printf("enter AOD!!");
1787 myPrimaryVertex = ((AliAODEvent*)fESD)->GetPrimaryVertex();
1788 if (!myPrimaryVertex) return;
1790 lPrimaryVtxPosition[0] = myPrimaryVertex->GetX();
1791 lPrimaryVtxPosition[1] = myPrimaryVertex->GetY();
1792 lPrimaryVtxPosition[2] = myPrimaryVertex->GetZ();
1794 // Cut on SPD vertex and fill histo Nevents: FIX it !
1796 // Tracks cuts FIX IT !
1799 lMagneticField = 999;
1804 fHistPrimaryVertexX->Fill(lPrimaryVtxPosition[0]);
1805 fHistPrimaryVertexY->Fill(lPrimaryVtxPosition[1]);
1806 fHistPrimaryVertexZ->Fill(lPrimaryVtxPosition[2]);
1807 //Double_t lrcPrimVtxR = TMath::Sqrt(lPrimaryVtxPosition[0]*lPrimaryVtxPosition[0]+lPrimaryVtxPosition[0]*lPrimaryVtxPosition[0]);
1809 fHistPrimaryVertexResX->Fill(lResPrimaryVtxX);
1810 fHistPrimaryVertexResY->Fill(lResPrimaryVtxY);
1811 fHistPrimaryVertexResZ->Fill(lResPrimaryVtxZ);
1814 //***********************
1815 // AliKF Primary Vertex
1817 AliKFVertex primaryVtxKF( *myPrimaryVertex );
1818 AliKFParticle::SetField(lMagneticField);
1821 //************************************
1824 AliESDpid *fESDpid = new AliESDpid(); // FIXME delete
1825 fESDpid->GetTPCResponse().SetBetheBlochParameters(fAlephParameters[0],fAlephParameters[1],fAlephParameters[2],fAlephParameters[3],fAlephParameters[4]);
1829 //***Rerun the V0 finder
1831 // fESD->ResetV0s();
1832 // AliV0vertexer v0Vertexer;
1833 // v0Vertexer.SetCuts(fCuts);
1834 // v0Vertexer.Tracks2V0vertices(fESD);
1836 //*************************
1839 nv0sTot = fESD->GetNumberOfV0s();
1840 if (!nv0sTot) fHistNumberEvents->Fill(6.5);
1842 for (Int_t iV0 = 0; iV0 < nv0sTot; iV0++) {
1845 AliKFParticle* negPiKF = NULL;
1846 AliKFParticle* posPiKF = NULL;
1847 AliKFParticle* posPKF = NULL;
1848 AliKFParticle* negAPKF = NULL;
1851 lIndexPosMother = 0; lIndexNegMother = 0; lIndexMotherOfMother = 0;
1852 lCheckPIdK0Short = 0; lCheckMcK0Short = 0; lCheckSecondaryK0s = 0;
1853 lCheckPIdLambda = 0; lCheckMcLambda = 0; lCheckSecondaryLambda = 0;
1854 lCheckPIdAntiLambda = 0; lCheckMcAntiLambda = 0; lCheckSecondaryAntiLambda = 0;
1855 lComeFromSigma = -1;
1858 if(fAnalysisType == "ESD") {
1861 AliESDv0 *v0 = ((AliESDEvent*)fESD)->GetV0(iV0);
1865 fHistPrimaryVertexPosXV0events->Fill(lPrimaryVtxPosition[0]);
1866 fHistPrimaryVertexPosYV0events->Fill(lPrimaryVtxPosition[1]);
1867 fHistPrimaryVertexPosZV0events->Fill(lPrimaryVtxPosition[2]);
1870 lIndexTrackPos = TMath::Abs(v0->GetPindex());
1871 lIndexTrackNeg = TMath::Abs(v0->GetNindex());
1872 AliESDtrack *myTrackPosTest = ((AliESDEvent*)fESD)->GetTrack(lIndexTrackPos);
1873 AliESDtrack *myTrackNegTest = ((AliESDEvent*)fESD)->GetTrack(lIndexTrackNeg);
1874 if (!myTrackPosTest || !myTrackNegTest) {
1875 // FIXME: shouldn't this be fatal?
1876 Printf("strange analysis::UserExec:: Error:Could not retreive one of the daughter track\n");
1880 if ( myTrackPosTest->GetSign() == myTrackNegTest->GetSign()){
1881 // FIXME: how can this happen?
1885 // FIXME: are the GetParamN/GetParamP reliable? If so, why do you need to check the sign below?
1886 // VO's main characteristics to check the reconstruction cuts
1887 lOnFlyStatus = v0->GetOnFlyStatus();
1888 lChi2V0 = v0->GetChi2V0();
1889 lDcaV0Daughters = v0->GetDcaV0Daughters();
1890 lDcaV0ToPrimVertex = v0->GetD(lPrimaryVtxPosition[0],lPrimaryVtxPosition[1],lPrimaryVtxPosition[2]);
1891 lV0cosPointAngle = v0->GetV0CosineOfPointingAngle(lPrimaryVtxPosition[0],lPrimaryVtxPosition[1], lPrimaryVtxPosition[2]);
1893 v0->GetXYZ(lV0Position[0], lV0Position[1], lV0Position[2]);
1895 lV0Radius = TMath::Sqrt(lV0Position[0]*lV0Position[0]+lV0Position[1]*lV0Position[1]);
1896 lV0DecayLength = TMath::Sqrt(TMath::Power(lV0Position[0] - lPrimaryVtxPosition[0],2) +
1897 TMath::Power(lV0Position[1] - lPrimaryVtxPosition[1],2) +
1898 TMath::Power(lV0Position[2] - lPrimaryVtxPosition[2],2 ));
1901 if( myTrackPosTest->GetSign() ==1){
1903 myTrackPos = ((AliESDEvent*)fESD)->GetTrack(lIndexTrackPos);
1904 myTrackNeg = ((AliESDEvent*)fESD)->GetTrack(lIndexTrackNeg);
1906 // Daughters' momentum;
1907 v0->GetPPxPyPz(lMomPos[0],lMomPos[1],lMomPos[2]);
1908 v0->GetNPxPyPz(lMomNeg[0],lMomNeg[1],lMomNeg[2]);
1910 if (negPiKF) delete negPiKF; negPiKF=NULL;
1911 if (posPiKF) delete posPiKF; posPiKF=NULL;
1912 if (posPKF) delete posPKF; posPKF=NULL;
1913 if (negAPKF) delete negAPKF; negAPKF=NULL;
1915 negPiKF = new AliKFParticle( *(v0->GetParamN()) ,-211);
1916 posPiKF = new AliKFParticle( *(v0->GetParamP()) ,211);
1917 posPKF = new AliKFParticle( *(v0->GetParamP()) ,2212);
1918 negAPKF = new AliKFParticle( *(v0->GetParamN()) ,-2212);
1922 if( myTrackPosTest->GetSign() ==-1){
1924 myTrackPos = ((AliESDEvent*)fESD)->GetTrack(lIndexTrackNeg);
1925 myTrackNeg = ((AliESDEvent*)fESD)->GetTrack(lIndexTrackPos);
1927 // Daughters' momentum;
1928 v0->GetPPxPyPz(lMomNeg[0],lMomNeg[1],lMomNeg[2]);
1929 v0->GetNPxPyPz(lMomPos[0],lMomPos[1],lMomPos[2]);
1931 if (negPiKF) delete negPiKF; negPiKF=NULL;
1932 if (posPiKF) delete posPiKF; posPiKF=NULL;
1933 if (posPKF) delete posPKF; posPKF=NULL;
1934 if (negAPKF) delete negAPKF; negAPKF=NULL;
1936 negPiKF = new AliKFParticle( *(v0->GetParamP()) ,-211);
1937 posPiKF = new AliKFParticle( *(v0->GetParamN()) ,211);
1938 posPKF = new AliKFParticle( *(v0->GetParamN()) ,2212);
1939 negAPKF = new AliKFParticle( *(v0->GetParamP()) ,-2212);
1943 lLabelTrackPos = (UInt_t)TMath::Abs(myTrackPos->GetLabel());
1944 lLabelTrackNeg = (UInt_t)TMath::Abs(myTrackNeg->GetLabel());
1946 // Daughters Pt and P:
1947 lPtPos = TMath::Sqrt(lMomPos[0]*lMomPos[0] + lMomPos[1]*lMomPos[1]);
1948 lPtNeg = TMath::Sqrt(lMomNeg[0]*lMomNeg[0] + lMomNeg[1]*lMomNeg[1]);
1950 lPPos = TMath::Sqrt(lMomPos[0]*lMomPos[0] + lMomPos[1]*lMomPos[1] + lMomPos[2]*lMomPos[2]);
1951 lPNeg = TMath::Sqrt(lMomNeg[0]*lMomNeg[0] + lMomNeg[1]*lMomNeg[1] + lMomNeg[2]*lMomNeg[2]);
1953 // Inner Wall parameter (used for pid):
1954 const AliExternalTrackParam *myInnerWallTrackPos = myTrackPos->GetInnerParam();
1955 if(myInnerWallTrackPos) lMomInnerWallPos = myInnerWallTrackPos->GetP();
1956 const AliExternalTrackParam *myInnerWallTrackNeg = myTrackNeg->GetInnerParam();
1957 if(myInnerWallTrackNeg) lMomInnerWallNeg = myInnerWallTrackNeg->GetP();
1959 // DCA between daughter and Primary Vertex:
1960 if (myTrackPos) lDcaPosToPrimVertex = TMath::Abs(myTrackPos->GetD(lPrimaryVtxPosition[0],lPrimaryVtxPosition[1],lMagneticField) );
1962 if (myTrackNeg) lDcaNegToPrimVertex = TMath::Abs(myTrackNeg->GetD(lPrimaryVtxPosition[0],lPrimaryVtxPosition[1],lMagneticField) );
1964 // Quality tracks cuts:
1965 if ( !(myTracksCuts->IsSelected(myTrackPos)) || !(myTracksCuts->IsSelected(myTrackNeg)) ) {
1966 if (negPiKF) delete negPiKF; negPiKF=NULL;
1967 if (posPiKF) delete posPiKF; posPiKF=NULL;
1968 if (posPKF) delete posPKF; posPKF=NULL;
1969 if (negAPKF) delete negAPKF; negAPKF=NULL;
1972 // Armenteros variables:
1973 lAlphaV0 = v0->AlphaV0();
1974 lPtArmV0 = v0->PtArmV0();
1980 if (fUsePID.Contains("withPID")) {
1981 nSigmaPosPion = TMath::Abs(fESDpid->NumberOfSigmasTPC(myTrackPos,AliPID::kPion));
1983 nSigmaNegPion = TMath::Abs(fESDpid->NumberOfSigmasTPC(myTrackNeg,AliPID::kPion));
1985 nSigmaPosProton = TMath::Abs(fESDpid->NumberOfSigmasTPC(myTrackPos,AliPID::kProton));
1987 nSigmaNegProton = TMath::Abs(fESDpid->NumberOfSigmasTPC(myTrackNeg,AliPID::kProton));
1990 nSigmaPosPion = 0; nSigmaNegPion =0; nSigmaPosProton = 0; nSigmaNegProton= 0;
1993 // Monte-Carlo particle associated to reconstructed particles:
1995 //if (lLabelTrackPos < 0 || lLabelTrackNeg < 0) continue;
1996 TParticle *lMCESDPartPos = stack->Particle(lLabelTrackPos);
1997 if(!lMCESDPartPos) {
1998 Printf("no MC particle for positive and/or negative daughter\n");
2000 if (negPiKF) delete negPiKF; negPiKF=NULL;
2001 if (posPiKF) delete posPiKF; posPiKF=NULL;
2002 if (posPKF) delete posPKF; posPKF=NULL;
2003 if (negAPKF) delete negAPKF; negAPKF=NULL;
2007 TParticle *lMCESDPartNeg = stack->Particle(lLabelTrackNeg);
2008 if (!lMCESDPartNeg) {
2009 if (negPiKF) delete negPiKF; negPiKF=NULL;
2010 if (posPiKF) delete posPiKF; posPiKF=NULL;
2011 if (posPKF) delete posPKF; posPKF=NULL;
2012 if (negAPKF) delete negAPKF; negAPKF=NULL;
2015 lPDGCodePosDaughter = lMCESDPartPos->GetPdgCode();
2016 lPDGCodeNegDaughter = lMCESDPartNeg->GetPdgCode();
2017 lIndexPosMother = lMCESDPartPos->GetFirstMother();
2018 lIndexNegMother = lMCESDPartNeg->GetFirstMother();
2020 if (lIndexPosMother == -1) {
2024 lIndexMotherOfMother = 0;
2038 TParticle *lMCESDMother = stack->Particle(lIndexPosMother);
2039 if (!lMCESDMother) {
2040 if (negPiKF) delete negPiKF; negPiKF=NULL;
2041 if (posPiKF) delete posPiKF; posPiKF=NULL;
2042 if (posPKF) delete posPKF; posPKF=NULL;
2043 if (negAPKF) delete negAPKF; negAPKF=NULL;
2046 lPdgcodeMother = lMCESDMother->GetPdgCode();
2047 lIndexMotherOfMother = lMCESDMother->GetFirstMother();
2048 if (lIndexMotherOfMother ==-1) lPdgcodeMotherOfMother = 0;
2050 TParticle *lMCESDMotherOfMother = stack->Particle(lIndexMotherOfMother);
2051 if (!lMCESDMotherOfMother) {
2052 if (negPiKF) delete negPiKF; negPiKF=NULL;
2053 if (posPiKF) delete posPiKF; posPiKF=NULL;
2054 if (posPKF) delete posPKF; posPKF=NULL;
2055 if (negAPKF) delete negAPKF; negAPKF=NULL;
2058 lPdgcodeMotherOfMother = lMCESDMotherOfMother->GetPdgCode();
2061 mcPosX = lMCESDPartPos->Vx();
2062 mcPosY = lMCESDPartPos->Vy();
2063 mcPosZ = lMCESDPartPos->Vz();
2064 mcPosR = TMath::Sqrt(mcPosX*mcPosX+mcPosY*mcPosY);
2065 mcPosMotherX = lMCESDMother->Vx();
2066 mcPosMotherY = lMCESDMother->Vy();
2067 mcPosMotherZ = lMCESDMother->Vz();
2068 mcPosMotherR = TMath::Sqrt(mcPosMotherX*mcPosMotherX+mcPosMotherY*mcPosMotherY);
2070 mcMotherPt = lMCESDMother->Pt();
2073 } // end ESD condition
2077 else if(fAnalysisType == "AOD") {
2079 AliAODv0 *myAODv0 = ((AliAODEvent*)fESD)->GetV0(iV0);
2080 if (!myAODv0) continue;
2083 fHistPrimaryVertexPosXV0events->Fill(lPrimaryVtxPosition[0]);
2084 fHistPrimaryVertexPosYV0events->Fill(lPrimaryVtxPosition[1]);
2085 fHistPrimaryVertexPosZV0events->Fill(lPrimaryVtxPosition[2]);
2089 if(lOnFlyStatus == fUseOnTheFly) nv0s++;
2092 lIndexTrackPos = TMath::Abs(myAODv0->GetPosID());
2093 lIndexTrackNeg = TMath::Abs(myAODv0->GetNegID());
2095 AliVParticle *lVPartPosTest = ((AliVEvent*)fESD)->GetTrack(lIndexTrackPos);
2096 AliVParticle *lVPartNegTest = ((AliVEvent*)fESD)->GetTrack(lIndexTrackNeg);
2097 //AliAODTrack *lVPartPos = ((AliAODEvent*)fESD)->GetTrack(lIndexTrackPos);
2098 //AliAODTrack *lVPartNeg = ((AliAODEvent*)fESD)->GetTrack(lIndexTrackNeg);
2100 if (!lVPartPosTest ||(!lVPartNegTest )) {
2101 Printf("strange analysis::UserExec:: Could not retreive one of the daughter track\n");
2108 // TPC refit condition (done during reconstruction for Offline but not for On-the-fly)
2109 //if( !(lVPartPosTest->GetStatus() & AliAODTrack::kTPCrefit)) continue;
2110 //if( !(lVPartNegTest->GetStatus() & AliAODTrack::kTPCrefit)) continue;
2113 lDcaPosToPrimVertex = myAODv0->DcaPosToPrimVertex();
2114 lDcaNegToPrimVertex = myAODv0->DcaNegToPrimVertex();
2115 lOnFlyStatus = myAODv0->GetOnFlyStatus();
2116 lChi2V0 = myAODv0->Chi2V0();
2117 lDcaV0Daughters = myAODv0->DcaV0Daughters();
2118 lDcaV0ToPrimVertex = myAODv0->DcaV0ToPrimVertex();
2119 lV0DecayLength = myAODv0->DecayLengthV0(lPrimaryVtxPosition);
2120 lV0cosPointAngle = myAODv0->CosPointingAngle(lPrimaryVtxPosition);
2121 lV0Radius = myAODv0->RadiusV0();
2123 if( lVPartPosTest->Charge() ==1){
2125 lVPartPos = ((AliVEvent*)fESD)->GetTrack(lIndexTrackPos);
2126 lVPartNeg = ((AliVEvent*)fESD)->GetTrack(lIndexTrackNeg);
2129 if (negPiKF) delete negPiKF; negPiKF=NULL;
2130 if (posPiKF) delete posPiKF; posPiKF=NULL;
2131 if (posPKF) delete posPKF; posPKF=NULL;
2132 if (negAPKF) delete negAPKF; negAPKF=NULL;
2134 //negPiKF = new AliKFParticle( *(myAODv0->GetParamN()) ,-211);
2135 //posPiKF = new AliKFParticle( *(myAODv0->GetParamP()) ,211);
2136 //posPKF = new AliKFParticle( *(myAODv0->GetParamP()) ,2212);
2137 //negAPKF = new AliKFParticle( *(myAODv0->GetParamN()) ,-2212);
2146 if( lVPartPosTest->Charge() ==-1){
2148 lVPartPos = ((AliVEvent*)fESD)->GetTrack(lIndexTrackNeg);
2149 lVPartNeg = ((AliVEvent*)fESD)->GetTrack(lIndexTrackPos);
2151 if (negPiKF) delete negPiKF; negPiKF=NULL;
2152 if (posPiKF) delete posPiKF; posPiKF=NULL;
2153 if (posPKF) delete posPKF; posPKF=NULL;
2154 if (negAPKF) delete negAPKF; negAPKF=NULL;
2156 //negPiKF = new AliKFParticle( *(myAODv0->GetParamP()) ,-211);
2157 //posPiKF = new AliKFParticle( *(myAODv0->GetParamN()) ,211);
2158 //posPKF = new AliKFParticle( *(myAODv0->GetParamN()) ,2212);
2159 //negAPKF = new AliKFParticle( *(myAODv0->GetParamP()) ,-2212);
2166 lLabelTrackPos = TMath::Abs(lVPartPos->GetLabel());
2167 lLabelTrackNeg = TMath::Abs(lVPartNeg->GetLabel());
2169 // Armenteros variables:
2170 lAlphaV0 = myAODv0->AlphaV0();
2171 lPtArmV0 = myAODv0->PtArmV0();
2174 lV0Eta = myAODv0->PseudoRapV0();
2176 // PID not accessible with AOD !
2177 nSigmaPosPion = 0; nSigmaNegPion =0; nSigmaPosProton = 0; nSigmaNegProton= 0;
2180 // Monte-Carlo particle associated to reconstructed particles:
2182 AliAODMCParticle *lMCAODPartPos = (AliAODMCParticle*)mcArray->At(lLabelTrackPos);
2183 if (!lMCAODPartPos) {
2184 if (negPiKF) delete negPiKF; negPiKF=NULL;
2185 if (posPiKF) delete posPiKF; posPiKF=NULL;
2186 if (posPKF) delete posPKF; posPKF=NULL;
2187 if (negAPKF) delete negAPKF; negAPKF=NULL;
2190 AliAODMCParticle *lMCAODPartNeg = (AliAODMCParticle*)mcArray->At(lLabelTrackNeg);
2192 // Printf("strange analysis::UserExec:no MC particle for negative daughter\n");
2194 if (negPiKF) delete negPiKF; negPiKF=NULL;
2195 if (posPiKF) delete posPiKF; posPiKF=NULL;
2196 if (posPKF) delete posPKF; posPKF=NULL;
2197 if (negAPKF) delete negAPKF; negAPKF=NULL;
2201 lPDGCodePosDaughter = lMCAODPartPos->GetPdgCode();
2202 lPDGCodeNegDaughter = lMCAODPartNeg->GetPdgCode();
2203 lIndexPosMother = lMCAODPartPos->GetMother();
2204 lIndexNegMother = lMCAODPartNeg->GetMother();
2206 AliAODMCParticle *lMCAODMother = (AliAODMCParticle*)mcArray->At(lIndexPosMother);
2207 lPdgcodeMother = lMCAODMother->GetPdgCode();
2208 lIndexMotherOfMother = lMCAODMother->GetMother();
2209 if (lIndexMotherOfMother ==-1) lPdgcodeMotherOfMother = 0;
2211 lPdgcodeMotherOfMother = ((AliAODMCParticle*)mcArray->At(lIndexMotherOfMother))->GetPdgCode();
2214 mcPosX = lMCAODPartPos->Xv();
2215 mcPosY = lMCAODPartPos->Yv();
2216 mcPosZ = lMCAODPartPos->Zv();
2217 mcPosR = TMath::Sqrt(mcPosX*mcPosX+mcPosY*mcPosY);
2218 mcPosMotherX = lMCAODMother->Xv();
2219 mcPosMotherY = lMCAODMother->Yv();
2220 mcPosMotherZ = lMCAODMother->Zv();
2221 mcPosMotherR = TMath::Sqrt(mcPosMotherX*mcPosMotherX+mcPosMotherY*mcPosMotherY);
2222 mcMotherPt = lMCAODMother->Pt();
2225 } // end AOD condition
2229 if(lOnFlyStatus == fUseOnTheFly) nv0s++;
2231 // Daughter momentum cut: ! FIX it in case of AOD !
2232 if ( (lPtPos < cutMinPtDaughter ) ||
2233 (lPtNeg < cutMinPtDaughter )
2235 if (negPiKF) delete negPiKF; negPiKF=NULL;
2236 if (posPiKF) delete posPiKF; posPiKF=NULL;
2237 if (posPKF) delete posPKF; posPKF=NULL;
2238 if (negAPKF) delete negAPKF; negAPKF=NULL;
2242 AliKFParticle v0K0sKF;
2243 v0K0sKF+=(*negPiKF);
2244 v0K0sKF+=(*posPiKF);
2245 v0K0sKF.SetProductionVertex(primaryVtxKF);
2247 AliKFParticle v0LambdaKF;
2248 v0LambdaKF+=(*negPiKF);
2249 v0LambdaKF+=(*posPKF);
2250 v0LambdaKF.SetProductionVertex(primaryVtxKF);
2252 AliKFParticle v0AntiLambdaKF;
2253 v0AntiLambdaKF+=(*posPiKF);
2254 v0AntiLambdaKF+=(*negAPKF);
2255 v0AntiLambdaKF.SetProductionVertex(primaryVtxKF);
2258 lInvMassK0s = v0K0sKF.GetMass();
2259 lInvMassLambda = v0LambdaKF.GetMass();
2260 lInvMassAntiLambda = v0AntiLambdaKF.GetMass();
2263 lRapK0s = 0.5*TMath::Log((v0K0sKF.GetE()+v0K0sKF.GetPz())/(v0K0sKF.GetE()-v0K0sKF.GetPz()+1.e-13));
2264 lRapLambda = 0.5*TMath::Log((v0LambdaKF.GetE()+v0LambdaKF.GetPz())/(v0LambdaKF.GetE()-v0LambdaKF.GetPz()+1.e-13));
2265 lRapAntiLambda = 0.5*TMath::Log((v0AntiLambdaKF.GetE()+v0AntiLambdaKF.GetPz())/(v0AntiLambdaKF.GetE()-v0AntiLambdaKF.GetPz()+1.e-13));
2268 lEtaK0s = v0K0sKF.GetEta();
2269 lEtaLambda = v0LambdaKF.GetEta();
2270 lEtaAntiLambda = v0AntiLambdaKF.GetEta();
2273 lPzK0s = v0K0sKF.GetPz();
2274 lPzLambda = v0LambdaKF.GetPz();
2275 lPzAntiLambda = v0AntiLambdaKF.GetPz();
2278 lPtK0s = v0K0sKF.GetPt();
2279 lPtLambda = v0LambdaKF.GetPt();
2280 lPtAntiLambda = v0AntiLambdaKF.GetPt();
2283 if (negPiKF) delete negPiKF; negPiKF=NULL;
2284 if (posPiKF) delete posPiKF; posPiKF=NULL;
2285 if (posPKF) delete posPKF; posPKF=NULL;
2286 if (negAPKF) delete negAPKF; negAPKF=NULL;
2287 // FIXME: should you really continue here and below?
2291 if (negPiKF) delete negPiKF; negPiKF=NULL;
2292 if (posPiKF) delete posPiKF; posPiKF=NULL;
2293 if (posPKF) delete posPKF; posPKF=NULL;
2294 if (negAPKF) delete negAPKF; negAPKF=NULL;
2298 deltaPtK0s = (lPtK0s - mcMotherPt)/mcMotherPt;
2299 deltaPtLambda = (lPtLambda - mcMotherPt)/mcMotherPt;
2300 deltaPtAntiLambda = (lPtAntiLambda - mcMotherPt)/mcMotherPt;
2303 lChi2KFK0s = v0K0sKF.GetChi2();
2304 lChi2KFLambda = v0LambdaKF.GetChi2();
2305 lChi2KFAntiLambda = v0AntiLambdaKF.GetChi2();
2307 // Reconstructed Position
2308 rcPosXK0s = v0K0sKF.GetX();
2309 rcPosYK0s = v0K0sKF.GetY();
2310 rcPosZK0s = v0K0sKF.GetZ();
2311 rcPosRK0s = TMath::Sqrt(rcPosXK0s*rcPosXK0s+rcPosYK0s*rcPosYK0s);
2313 rcPosXLambda = v0LambdaKF.GetX();
2314 rcPosYLambda = v0LambdaKF.GetY();
2315 rcPosZLambda = v0LambdaKF.GetZ();
2316 rcPosRLambda = TMath::Sqrt(rcPosXLambda*rcPosXLambda+rcPosYLambda*rcPosYLambda);
2318 rcPosXAntiLambda = v0AntiLambdaKF.GetX();
2319 rcPosYAntiLambda = v0AntiLambdaKF.GetY();
2320 rcPosZAntiLambda = v0AntiLambdaKF.GetZ();
2321 rcPosRAntiLambda = TMath::Sqrt(rcPosXAntiLambda*rcPosXAntiLambda+rcPosYAntiLambda*rcPosYAntiLambda);
2323 TVector3 momPos(lMomPos[0],lMomPos[1],lMomPos[2]);
2324 TVector3 momNeg(lMomNeg[0],lMomNeg[1],lMomNeg[2]);
2325 TVector3 momTotK0s(v0K0sKF.GetPx(),v0K0sKF.GetPy(),v0K0sKF.GetPz());
2326 TVector3 momTotLambda(v0LambdaKF.GetPx(),v0LambdaKF.GetPy(),v0LambdaKF.GetPz());
2327 TVector3 momTotAntiLambda(v0AntiLambdaKF.GetPx(),v0AntiLambdaKF.GetPy(),v0AntiLambdaKF.GetPz());
2329 lQlPos = momPos.Dot(momTotK0s)/momTotK0s.Mag();
2330 lQlNeg = momNeg.Dot(momTotK0s)/momTotK0s.Mag();
2331 lAlphaV0K0s = 1.-2./(1.+lQlPos/lQlNeg);
2332 lQlPos = momPos.Dot(momTotLambda)/momTotLambda.Mag();
2333 lQlNeg = momNeg.Dot(momTotLambda)/momTotLambda.Mag();
2334 lAlphaV0Lambda = 1.-2./(1.+lQlPos/lQlNeg);
2335 lQlPos = momPos.Dot(momTotAntiLambda)/momTotAntiLambda.Mag();
2336 lQlNeg = momNeg.Dot(momTotAntiLambda)/momTotAntiLambda.Mag();
2337 lAlphaV0AntiLambda = 1.-2./(1.+lQlPos/lQlNeg);
2339 lPtArmV0K0s = momPos.Perp(momTotK0s);
2340 lPtArmV0Lambda = momPos.Perp(momTotLambda);
2341 lPtArmV0AntiLambda = momPos.Perp(momTotAntiLambda);
2343 // Look for associated particles:
2345 if( (lIndexPosMother==-1) || (lIndexNegMother==-1) ) {
2346 fHistMCDaughterTrack->Fill(1);
2349 else if( ( (lPDGCodePosDaughter==+211) && (lPDGCodeNegDaughter==-211) )
2351 lCheckPIdK0Short = 1;
2352 fHistMCDaughterTrack->Fill(3);
2353 if ( (lIndexPosMother==lIndexNegMother) &&
2354 (lPdgcodeMother==310) ) {
2355 if (lIndexPosMother <= lNbMCPrimary) lCheckMcK0Short = 1;
2356 else lCheckSecondaryK0s = 1;
2359 else if( ( (lPDGCodePosDaughter==+2212) && (lPDGCodeNegDaughter==-211) )
2361 lCheckPIdLambda = 1;
2362 fHistMCDaughterTrack->Fill(5);
2363 if ( (lIndexPosMother==lIndexNegMother) &&
2364 (lPdgcodeMother==3122) ){
2365 if ( ( TMath::Abs(lPdgcodeMotherOfMother) == 3212) ||
2366 ( TMath::Abs(lPdgcodeMotherOfMother) == 3224) ||
2367 ( TMath::Abs(lPdgcodeMotherOfMother) == 3214) ||
2368 ( TMath::Abs(lPdgcodeMotherOfMother) == 3114)
2369 ) lComeFromSigma = 1;
2370 else lComeFromSigma = 0;
2371 if ( (lIndexPosMother <= lNbMCPrimary) ||
2372 ( ( lIndexPosMother > lNbMCPrimary) && (lComeFromSigma) )
2373 ) lCheckMcLambda = 1;
2374 else lCheckSecondaryLambda = 1;
2377 else if( ( (lPDGCodePosDaughter==211) && (lPDGCodeNegDaughter==-2212) )
2379 lCheckPIdAntiLambda = 1;
2380 fHistMCDaughterTrack->Fill(7);
2381 if ( (lIndexPosMother==lIndexNegMother) &&
2382 (lPdgcodeMother==-3122) ) {
2383 if ( ( TMath::Abs(lPdgcodeMotherOfMother) == 3212) ||
2384 ( TMath::Abs(lPdgcodeMotherOfMother) == 3224) ||
2385 ( TMath::Abs(lPdgcodeMotherOfMother) == 3214) ||
2386 ( TMath::Abs(lPdgcodeMotherOfMother) == 3114)
2387 ) lComeFromSigma = 1;
2388 else lComeFromSigma = 0;
2389 if ( (lIndexPosMother <= lNbMCPrimary) ||
2390 ( ( lIndexPosMother > lNbMCPrimary) && (lComeFromSigma) )
2391 ) lCheckMcAntiLambda = 1;
2392 else lCheckSecondaryAntiLambda = 1;
2397 else if ( (lPDGCodePosDaughter==11) &&
2398 (lPDGCodeNegDaughter==-11) &&
2399 (lPdgcodeMother==22 ) )
2401 } // end "look for associated particles
2405 /* if (fUseCut.Contains("yes")) {
2406 if ( (lDcaPosToPrimVertex < 0.036 ) ||
2407 (lDcaNegToPrimVertex < 0.036 ) ||
2408 (lDcaV0Daughters > 0.5 ) ||
2409 (lV0cosPointAngle < 0.999 )
2416 if ( (lDcaV0Daughters > 0.3 ) ||
2417 (lV0cosPointAngle < 0.998 )
2422 lCheckPIDK0sPosDaughter = 0, lCheckPIDK0sNegDaughter = 0;
2423 lCheckPIDLambdaPosDaughter = 0, lCheckPIDLambdaNegDaughter = 0;
2424 lCheckPIDAntiLambdaPosDaughter = 0, lCheckPIDAntiLambdaNegDaughter = 0;
2426 if (lMomInnerWallPos < lLimitPPID) {
2427 if (nSigmaPosPion < cutNSigmaLowP) {
2428 lCheckPIDK0sPosDaughter = 1;
2429 lCheckPIDAntiLambdaPosDaughter = 1;
2431 if (nSigmaPosProton < cutNSigmaLowP) lCheckPIDLambdaPosDaughter = 1;
2434 else if (lMomInnerWallPos > lLimitPPID) {
2435 if (nSigmaPosPion < cutNSigmaHighP) {
2436 lCheckPIDK0sPosDaughter = 1;
2437 lCheckPIDAntiLambdaPosDaughter = 1;
2439 if (nSigmaPosProton < cutNSigmaHighP) lCheckPIDLambdaPosDaughter = 1;
2442 if (lMomInnerWallNeg < lLimitPPID) {
2443 if (nSigmaNegPion < cutNSigmaLowP) {
2444 lCheckPIDK0sNegDaughter = 1;
2445 lCheckPIDLambdaNegDaughter = 1;
2447 if (nSigmaNegProton < cutNSigmaLowP) lCheckPIDAntiLambdaNegDaughter = 1;
2450 else if (lMomInnerWallNeg > lLimitPPID) {
2451 if (nSigmaNegPion < cutNSigmaHighP) {
2452 lCheckPIDK0sNegDaughter = 1;
2453 lCheckPIDLambdaNegDaughter = 1;
2455 if (nSigmaNegProton < cutNSigmaHighP) lCheckPIDAntiLambdaNegDaughter = 1;
2459 //*****************************
2460 // filling histograms
2462 fHistDcaPosToPrimVertex->Fill(lDcaPosToPrimVertex,lOnFlyStatus);
2463 fHistDcaNegToPrimVertex->Fill(lDcaNegToPrimVertex,lOnFlyStatus);
2464 fHistDcaPosToPrimVertexZoom->Fill(lDcaPosToPrimVertex,lOnFlyStatus);
2465 fHistDcaNegToPrimVertexZoom->Fill(lDcaNegToPrimVertex,lOnFlyStatus);
2466 fHistRadiusV0->Fill(lV0Radius,lOnFlyStatus);
2467 fHistDecayLengthV0->Fill(lV0DecayLength,lOnFlyStatus);
2468 fHistDcaV0Daughters->Fill(lDcaV0Daughters,lOnFlyStatus);
2469 fHistChi2->Fill(lChi2V0,lOnFlyStatus);
2470 fHistCosPointAngle->Fill(lV0cosPointAngle,lOnFlyStatus);
2471 if (lV0cosPointAngle >= 0.9) fHistCosPointAngleZoom->Fill(lV0cosPointAngle,lOnFlyStatus);
2472 fHistChi2KFBeforeCutK0s->Fill(lChi2KFK0s,lOnFlyStatus);
2473 fHistChi2KFBeforeCutLambda->Fill(lChi2KFLambda,lOnFlyStatus);
2474 fHistChi2KFBeforeCutAntiLambda->Fill(lChi2KFAntiLambda,lOnFlyStatus);
2477 // Histo versus Rap and armenteros plot
2478 if (lOnFlyStatus == fUseOnTheFly){
2479 if (lCheckMcK0Short) fHistAsMcRapK0->Fill(lRapK0s);
2480 if (lCheckMcLambda) fHistAsMcRapLambda->Fill(lRapLambda);
2481 if (lCheckMcAntiLambda) fHistAsMcRapLambda->Fill(lRapAntiLambda);
2482 // fHistArmenterosPodolanski->Fill(lAlphaV0,lPtArmV0);
2483 // fHistDaughterPt->Fill(lPtPos,lPtNeg);
2486 // FIXME: associated histos, what are they used for?
2488 // K0s associated histograms in |rap| < lCutRap:
2491 if ( lCheckPIDK0sPosDaughter && lCheckPIDK0sNegDaughter
2492 && (lChi2KFK0s < cutChi2KF) && (TMath::Abs(lPzK0s/lPtK0s)<0.7) ) {
2495 fHistChi2KFAfterCutK0s->Fill(lChi2KFK0s,lOnFlyStatus);
2497 if (TMath::Abs(lRapK0s) < lCutRap) {
2499 fHistNsigmaPosPionK0->Fill(nSigmaPosPion);
2500 fHistNsigmaNegPionK0->Fill(nSigmaNegPion);
2502 if (lOnFlyStatus == fUseOnTheFly){
2503 fHistMassK0->Fill(lInvMassK0s);
2504 fHistMassVsRadiusK0->Fill(rcPosRK0s,lInvMassK0s);
2505 fHistPtVsMassK0->Fill(lInvMassK0s,lPtK0s);
2508 // fHistMultVsPtVsMassK0->Fill(multiplicity ,lInvMassK0s,lPtK0s);
2509 if(lCheckPIdK0Short) fHistPidMcMassK0->Fill(lInvMassK0s);
2510 if(lCheckMcK0Short) {
2511 fHistAsMcMassK0->Fill(lInvMassK0s);
2512 fHistAsMcPtK0->Fill(lPtK0s);
2515 fHistAsMcPtVsMassK0->Fill(lInvMassK0s,lPtK0s);
2516 if (lPtK0s <= 1) fHistAsMcPtZoomK0->Fill(lPtK0s);
2517 fHistAsMcMassVsRadiusK0->Fill(rcPosRK0s,lInvMassK0s);
2518 fHistAsMcResxK0->Fill(rcPosXK0s-mcPosX);
2519 fHistAsMcResyK0->Fill(rcPosYK0s-mcPosY);
2520 fHistAsMcReszK0->Fill(rcPosZK0s-mcPosZ);
2521 fHistAsMcResrVsRadiusK0->Fill(rcPosRK0s,rcPosRK0s-mcPosR);
2522 fHistAsMcReszVsRadiusK0->Fill(rcPosZK0s,rcPosZK0s-mcPosZ);
2523 fHistAsMcProdRadiusK0->Fill(mcPosMotherR);
2524 fHistAsMcProdRadiusXvsYK0s->Fill(mcPosMotherX,mcPosMotherY);
2525 fHistAsMcResPtK0->Fill(deltaPtK0s);
2526 fHistAsMcResPtVsRapK0->Fill(deltaPtK0s,lRapK0s);
2527 fHistAsMcResPtVsPtK0->Fill(deltaPtK0s,lPtK0s);
2529 else if (lCheckSecondaryK0s) {
2530 fHistAsMcSecondaryPtVsRapK0s->Fill(lPtK0s,lRapK0s);
2531 fHistAsMcSecondaryProdRadiusK0s->Fill(mcPosMotherR);
2532 fHistAsMcSecondaryProdRadiusXvsYK0s->Fill(mcPosMotherX,mcPosMotherY);
2533 switch (lPdgcodeMotherOfMother) {
2534 case 130 : fHistAsMcSecondaryMotherPdgCodeK0s->Fill(0.5);break; // K0L
2535 case 321 : fHistAsMcSecondaryMotherPdgCodeK0s->Fill(1.5);break; // K+
2536 case -321 : fHistAsMcSecondaryMotherPdgCodeK0s->Fill(2.5);break; // K-
2537 case -3122 : fHistAsMcSecondaryMotherPdgCodeK0s->Fill(3.5);break; //AntiLambda
2538 default : fHistAsMcSecondaryMotherPdgCodeK0s->Fill(6.5);break;
2542 } // end rapidity condition
2543 } // end nsigma condition
2546 // Associated Lambda histograms in |rap| < lCutRap
2549 if ( lCheckPIDLambdaPosDaughter && lCheckPIDLambdaNegDaughter
2550 && (lChi2KFLambda < cutChi2KF) && (TMath::Abs(lPzLambda/lPtLambda)<0.7) ) {
2554 fHistChi2KFAfterCutLambda->Fill(lChi2KFLambda,lOnFlyStatus);
2556 if (TMath::Abs(lRapLambda) < lCutRap) {
2558 fHistNsigmaPosProtonLambda->Fill(nSigmaPosProton);
2559 fHistNsigmaNegPionLambda->Fill(nSigmaNegPion);
2560 if (lOnFlyStatus == fUseOnTheFly){
2561 fHistMassLambda->Fill(lInvMassLambda);
2562 fHistMassVsRadiusLambda->Fill(rcPosRLambda,lInvMassLambda);
2563 fHistPtVsMassLambda->Fill(lInvMassLambda,lPtLambda);
2565 // fHistMultVsPtVsMassLambda->Fill(multiplicity ,lInvMassLambda,lPtLambda);
2566 if(lCheckPIdLambda) fHistPidMcMassLambda->Fill(lInvMassLambda);
2568 if(lCheckMcLambda) {
2569 fHistAsMcMassLambda->Fill(lInvMassLambda);
2570 fHistAsMcPtLambda->Fill(lPtLambda);
2573 fHistAsMcPtVsMassLambda->Fill(lInvMassLambda,lPtLambda);
2574 if (lPtLambda <= 1) fHistAsMcPtZoomLambda->Fill(lPtLambda);
2575 fHistAsMcMassVsRadiusLambda->Fill(rcPosRLambda,lInvMassLambda);
2576 fHistAsMcResxLambda->Fill(rcPosXLambda-mcPosX);
2577 fHistAsMcResyLambda->Fill(rcPosYLambda-mcPosY);
2578 fHistAsMcReszLambda->Fill(rcPosZLambda-mcPosZ);
2579 fHistAsMcResrVsRadiusLambda->Fill(rcPosRLambda,rcPosRLambda-mcPosR);
2580 fHistAsMcReszVsRadiusLambda->Fill(rcPosZLambda,rcPosZLambda-mcPosZ);
2581 fHistAsMcProdRadiusLambda->Fill(mcPosMotherR);
2582 fHistAsMcProdRadiusXvsYLambda->Fill(mcPosMotherX,mcPosMotherY);
2583 fHistAsMcResPtLambda->Fill(deltaPtLambda);
2584 fHistAsMcResPtVsRapLambda->Fill(deltaPtLambda,lRapLambda);
2585 fHistAsMcResPtVsPtLambda->Fill(deltaPtLambda,lPtLambda);
2586 if (lComeFromSigma) fHistAsMcPtLambdaFromSigma->Fill(lPtLambda);
2587 switch (lPdgcodeMotherOfMother) {
2588 case 3222 : fHistAsMcMotherPdgCodeLambda->Fill(0.5); break; // Sigma +
2589 case 3212 : fHistAsMcMotherPdgCodeLambda->Fill(1.5); break; // Sigma 0
2590 case 3112 : fHistAsMcMotherPdgCodeLambda->Fill(2.5); break;// Sigma -
2591 case 3224 : fHistAsMcMotherPdgCodeLambda->Fill(3.5); break;// Sigma(1385) +
2592 case 3214 : fHistAsMcMotherPdgCodeLambda->Fill(4.5); break;// Sigma(1385) 0
2593 case 3114 : fHistAsMcMotherPdgCodeLambda->Fill(5.5); break;// Sigma(1385) -
2594 case 3322 : fHistAsMcMotherPdgCodeLambda->Fill(6.5); break; // Xi 0
2595 case 3312 : fHistAsMcMotherPdgCodeLambda->Fill(7.5); break; // Xi -
2596 case 3334 : fHistAsMcMotherPdgCodeLambda->Fill(8.5); break; // Omega
2597 case -1 : fHistAsMcMotherPdgCodeLambda->Fill(9.5); break;
2598 default : fHistAsMcMotherPdgCodeLambda->Fill(10.5);break;
2602 //printf("found Lambda RC dcaPos=%e dcaNeg=%e dcaDau=%e cosP=%e pT=%e mass=%e\n",lDcaPosToPrimVertex ,lDcaNegToPrimVertex ,lDcaV0Daughters,lV0cosPointAngle,lPtLambda,lInvMassLambda);
2603 //printf("found Lambda RC Pindex=%d Nindex=%d Plabel=%d Nlabel=%d\n\n",lIndexTrackPos,lIndexTrackNeg,lLabelTrackPos,lLabelTrackNeg);
2607 else if (lCheckSecondaryLambda) {
2608 fHistAsMcSecondaryPtVsRapLambda->Fill(lPtLambda,lRapLambda);
2609 fHistAsMcSecondaryProdRadiusLambda->Fill(mcPosMotherR);
2610 fHistAsMcSecondaryProdRadiusXvsYLambda->Fill(mcPosMotherX,mcPosMotherY);
2611 if (lComeFromSigma) fHistAsMcSecondaryPtLambdaFromSigma->Fill(lPtLambda);
2612 printf(" lPdgcodeMotherOfMother= %d",lPdgcodeMotherOfMother);
2613 switch (lPdgcodeMotherOfMother) {
2614 case 3222 : fHistAsMcSecondaryMotherPdgCodeLambda->Fill(0.5); break;// Sigma +
2615 case 3212 : fHistAsMcSecondaryMotherPdgCodeLambda->Fill(1.5); break;// Sigma 0
2616 case 3112 : fHistAsMcSecondaryMotherPdgCodeLambda->Fill(2.5); break;// Sigma -
2617 case 3224 : fHistAsMcSecondaryMotherPdgCodeLambda->Fill(3.5); break;// Sigma(1385) +
2618 case 3214 : fHistAsMcSecondaryMotherPdgCodeLambda->Fill(4.5); break;// Sigma(1385) 0
2619 case 3114 : fHistAsMcSecondaryMotherPdgCodeLambda->Fill(5.5); break;// Sigma(1385) -
2620 case 3322 : fHistAsMcSecondaryMotherPdgCodeLambda->Fill(6.5); break; // Xi 0
2621 case 3312 : fHistAsMcSecondaryMotherPdgCodeLambda->Fill(7.5); break; // Xi -
2622 case 3334 : fHistAsMcSecondaryMotherPdgCodeLambda->Fill(8.5); break; // Omega
2623 case -1 : fHistAsMcSecondaryMotherPdgCodeLambda->Fill(9.5); break;
2624 default : fHistAsMcSecondaryMotherPdgCodeLambda->Fill(10.5);break;
2628 } // end rapidity condition
2629 } //end nsigma condition - lambda
2632 if (negPiKF) delete negPiKF; negPiKF= NULL;
2633 if (posPiKF) delete posPiKF; posPiKF= NULL;
2634 if (posPKF) delete posPKF; posPKF = NULL;
2635 if (negAPKF) delete negAPKF; negAPKF= NULL;
2639 fHistV0Multiplicity->Fill(nv0s);
2640 if (fAnalysisType == "ESD") { if(myPrimaryVertex) delete myPrimaryVertex; }
2642 if(myTracksCuts) delete myTracksCuts;
2645 PostData(1, fListHist);
2646 PostData(2, fCentrSelector);
2649 //________________________________________________________________________
2650 void AliAnalysisTaskPerformanceStrange::Terminate(Option_t *)
2652 // Draw result to the screen
2653 // Called once at the end of the query
2655 TList *cRetrievedList = 0x0;
2656 cRetrievedList = (TList*)GetOutputData(1);
2658 if(!cRetrievedList){
2659 AliWarning("ERROR - AliAnalysisTaskPerformanceStrange: output data container list not available\n"); return;
2663 fHistV0Multiplicity = dynamic_cast<TH1F*> ( cRetrievedList->FindObject("fHistV0Multiplicity"));
2664 if (!fHistV0Multiplicity) {
2665 Printf("ERROR: fHistV0Multiplicity not available");
2669 fHistV0MultiplicityMI = dynamic_cast<TH1F*> ( cRetrievedList->FindObject("fHistV0MultiplicityMI"));
2670 if (!fHistV0MultiplicityMI) {
2671 Printf("ERROR: fHistV0MultiplicityMI not available");
2675 TCanvas *canPerformanceStrange = new TCanvas("AliAnalysisTaskCheckV0","Multiplicity",10,10,510,510);
2676 canPerformanceStrange->Divide(2,1);
2677 if (fHistV0Multiplicity->GetMaximum() > 0.) canPerformanceStrange->cd(1)->SetLogy();
2678 fHistV0Multiplicity->SetMarkerStyle(25);
2679 fHistV0Multiplicity->DrawCopy("E");
2680 if (fHistV0MultiplicityMI->GetMaximum() > 0.) canPerformanceStrange->cd(2)->SetLogy();
2681 fHistV0MultiplicityMI->SetMarkerStyle(24);
2682 fHistV0MultiplicityMI->DrawCopy("E");
2688 //----------------------------------------------------------------------------
2690 Double_t AliAnalysisTaskPerformanceStrange::MyRapidity(Double_t rE, Double_t rPz) const
2692 // Local calculation for rapidity
2693 return 0.5*TMath::Log((rE+rPz)/(rE-rPz+1.e-13));
2695 //----------------------------------------------------------------------------