]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGLF/STRANGENESS/LambdaK0/AliAnalysisTaskExtractPerformanceV0.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
analysis code and AliMuonForwardTrackPair updated
[u/mrichter/AliRoot.git] / PWGLF / STRANGENESS / LambdaK0 / AliAnalysisTaskExtractPerformanceV0.cxx
CommitLineData
bcb6ffdb 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16// +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
17//
18// Modified version of AliAnalysisTaskCheckCascade.cxx.
19// This is a 'hybrid' output version, in that it uses a classic TTree
20// ROOT object to store the candidates, plus a couple of histograms filled on
21// a per-event basis for storing variables too numerous to put in a tree.
22//
23// --- Adapted to look for lambdas as well, using code from
24// AliAnalysisTaskCheckPerformanceStrange.cxx
25//
26// --- Algorithm Description
27// 1. Loop over primaries in stack to acquire generated charged Xi
28// 2. Loop over stack to find V0s, fill TH3Fs "PrimRawPt"s for Efficiency
29// 3. Perform Physics Selection
30// 4. Perform Primary Vertex |z|<10cm selection
31// 5. Perform Primary Vertex NoTPCOnly vertexing selection (>0 contrib.)
32// 6. Perform Pileup Rejection
33// 7. Analysis Loops:
34// 7a. Fill TH3Fs "PrimAnalysisPt" for control purposes only
35// 7b. Fill TTree object with V0 information, candidates
36//
37// Please Report Any Bugs!
38//
39// --- David Dobrigkeit Chinellato
40// (david.chinellato@gmail.com)
41//
42// +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
43
44class TTree;
45class TParticle;
46class TVector3;
47
48//class AliMCEventHandler;
49//class AliMCEvent;
50//class AliStack;
51
52class AliESDVertex;
53class AliAODVertex;
54class AliESDv0;
55class AliAODv0;
56
57#include <Riostream.h>
58#include "TList.h"
59#include "TH1.h"
60#include "TH2.h"
61#include "TH3.h"
62#include "TFile.h"
63#include "THnSparse.h"
64#include "TVector3.h"
65#include "TCanvas.h"
66#include "TMath.h"
67#include "TLegend.h"
68#include "AliLog.h"
69
70#include "AliESDEvent.h"
71#include "AliAODEvent.h"
72#include "AliV0vertexer.h"
73#include "AliCascadeVertexer.h"
74#include "AliESDpid.h"
75#include "AliESDtrack.h"
76#include "AliESDtrackCuts.h"
77#include "AliInputEventHandler.h"
78#include "AliAnalysisManager.h"
79#include "AliMCEventHandler.h"
80#include "AliMCEvent.h"
81#include "AliStack.h"
82
83#include "AliCFContainer.h"
84#include "AliMultiplicity.h"
85#include "AliAODMCParticle.h"
86#include "AliESDcascade.h"
87#include "AliAODcascade.h"
88#include "AliESDUtils.h"
7e81bb32 89#include "AliGenEventHeader.h"
bcb6ffdb 90
91#include "AliAnalysisTaskExtractPerformanceV0.h"
92
c64cb1f6 93using std::cout;
94using std::endl;
95
bcb6ffdb 96ClassImp(AliAnalysisTaskExtractPerformanceV0)
97
98AliAnalysisTaskExtractPerformanceV0::AliAnalysisTaskExtractPerformanceV0()
76029adc 99 : AliAnalysisTaskSE(), fListHistV0(0), fTree(0), fPIDResponse(0), fESDtrackCuts(0),
548bd092 100 fkIsNuclear ( kFALSE ),
101 fkLowEnergyPP ( kFALSE ),
102 fkUseOnTheFly ( kFALSE ),
bcb6ffdb 103
104//------------------------------------------------
105// HISTOGRAMS
106// --- Filled on an Event-by-event basis
107//------------------------------------------------
108 fHistV0MultiplicityBeforeTrigSel(0),
109 fHistV0MultiplicityForTrigEvt(0),
110 fHistV0MultiplicityForSelEvt(0),
111 fHistV0MultiplicityForSelEvtNoTPCOnly(0),
112 fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup(0),
ac159816 113 fHistMultiplicityBeforeTrigSel(0),
114 fHistMultiplicityForTrigEvt(0),
bcb6ffdb 115 fHistMultiplicity(0),
116 fHistMultiplicityNoTPCOnly(0),
117 fHistMultiplicityNoTPCOnlyNoPileup(0),
118
119//------------------------------------------------
120// PARTICLE HISTOGRAMS
121// --- Filled on a Particle-by-Particle basis
122//------------------------------------------------
123 f3dHistPrimAnalysisPtVsYVsMultLambda(0),
124 f3dHistPrimAnalysisPtVsYVsMultAntiLambda(0),
125 f3dHistPrimAnalysisPtVsYVsMultK0Short(0),
126 f3dHistPrimRawPtVsYVsMultLambda(0),
127 f3dHistPrimRawPtVsYVsMultAntiLambda(0),
128 f3dHistPrimRawPtVsYVsMultK0Short(0),
7e81bb32 129 f3dHistPrimCloseToPVPtVsYVsMultLambda(0),
130 f3dHistPrimCloseToPVPtVsYVsMultAntiLambda(0),
131 f3dHistPrimCloseToPVPtVsYVsMultK0Short(0),
bcb6ffdb 132 f3dHistPrimRawPtVsYVsDecayLengthLambda(0),
133 f3dHistPrimRawPtVsYVsDecayLengthAntiLambda(0),
134 f3dHistPrimRawPtVsYVsDecayLengthK0Short(0),
135 f3dHistGenPtVsYVsMultXiMinus(0),
136 f3dHistGenPtVsYVsMultXiPlus(0),
c228a5af 137 f3dHistGenPtVsYVsMultOmegaMinus(0),
138 f3dHistGenPtVsYVsMultOmegaPlus(0),
139 f3dHistGenSelectedPtVsYVsMultXiMinus(0),
140 f3dHistGenSelectedPtVsYVsMultXiPlus(0),
141 f3dHistGenSelectedPtVsYVsMultOmegaMinus(0),
142 f3dHistGenSelectedPtVsYVsMultOmegaPlus(0),
bcb6ffdb 143 fHistPVx(0),
144 fHistPVy(0),
145 fHistPVz(0),
146 fHistPVxAnalysis(0),
147 fHistPVyAnalysis(0),
148 fHistPVzAnalysis(0),
149 fHistPVxAnalysisHasHighPtLambda(0),
150 fHistPVyAnalysisHasHighPtLambda(0),
4f8ed367 151 fHistPVzAnalysisHasHighPtLambda(0),
152 fHistSwappedV0Counter(0)
bcb6ffdb 153{
154 // Dummy Constructor
155}
156
157AliAnalysisTaskExtractPerformanceV0::AliAnalysisTaskExtractPerformanceV0(const char *name)
76029adc 158 : AliAnalysisTaskSE(name), fListHistV0(0), fTree(0), fPIDResponse(0), fESDtrackCuts(0),
548bd092 159 fkIsNuclear ( kFALSE ),
160 fkLowEnergyPP ( kFALSE ),
161 fkUseOnTheFly ( kFALSE ),
bcb6ffdb 162
163//------------------------------------------------
164// HISTOGRAMS
165// --- Filled on an Event-by-event basis
166//------------------------------------------------
167 fHistV0MultiplicityBeforeTrigSel(0),
168 fHistV0MultiplicityForTrigEvt(0),
169 fHistV0MultiplicityForSelEvt(0),
170 fHistV0MultiplicityForSelEvtNoTPCOnly(0),
171 fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup(0),
ac159816 172 fHistMultiplicityBeforeTrigSel(0),
173 fHistMultiplicityForTrigEvt(0),
bcb6ffdb 174 fHistMultiplicity(0),
175 fHistMultiplicityNoTPCOnly(0),
176 fHistMultiplicityNoTPCOnlyNoPileup(0),
177
178
179//------------------------------------------------
180// PARTICLE HISTOGRAMS
181// --- Filled on a Particle-by-Particle basis
182//------------------------------------------------
183 f3dHistPrimAnalysisPtVsYVsMultLambda(0),
184 f3dHistPrimAnalysisPtVsYVsMultAntiLambda(0),
185 f3dHistPrimAnalysisPtVsYVsMultK0Short(0),
186 f3dHistPrimRawPtVsYVsMultLambda(0),
187 f3dHistPrimRawPtVsYVsMultAntiLambda(0),
188 f3dHistPrimRawPtVsYVsMultK0Short(0),
7e81bb32 189 f3dHistPrimCloseToPVPtVsYVsMultLambda(0),
190 f3dHistPrimCloseToPVPtVsYVsMultAntiLambda(0),
191 f3dHistPrimCloseToPVPtVsYVsMultK0Short(0),
bcb6ffdb 192 f3dHistPrimRawPtVsYVsDecayLengthLambda(0),
193 f3dHistPrimRawPtVsYVsDecayLengthAntiLambda(0),
194 f3dHistPrimRawPtVsYVsDecayLengthK0Short(0),
195 f3dHistGenPtVsYVsMultXiMinus(0),
196 f3dHistGenPtVsYVsMultXiPlus(0),
c228a5af 197 f3dHistGenPtVsYVsMultOmegaMinus(0),
198 f3dHistGenPtVsYVsMultOmegaPlus(0),
199 f3dHistGenSelectedPtVsYVsMultXiMinus(0),
200 f3dHistGenSelectedPtVsYVsMultXiPlus(0),
201 f3dHistGenSelectedPtVsYVsMultOmegaMinus(0),
202 f3dHistGenSelectedPtVsYVsMultOmegaPlus(0),
bcb6ffdb 203 fHistPVx(0),
204 fHistPVy(0),
205 fHistPVz(0),
206 fHistPVxAnalysis(0),
207 fHistPVyAnalysis(0),
208 fHistPVzAnalysis(0),
209 fHistPVxAnalysisHasHighPtLambda(0),
210 fHistPVyAnalysisHasHighPtLambda(0),
4f8ed367 211 fHistPVzAnalysisHasHighPtLambda(0),
212 fHistSwappedV0Counter(0)
bcb6ffdb 213{
214 // Constructor
215 // Output slot #0 writes into a TList container (Cascade)
216 DefineOutput(1, TList::Class());
217 DefineOutput(2, TTree::Class());
218}
219
220
221AliAnalysisTaskExtractPerformanceV0::~AliAnalysisTaskExtractPerformanceV0()
222{
223//------------------------------------------------
224// DESTRUCTOR
225//------------------------------------------------
226
227 if (fListHistV0){
228 delete fListHistV0;
229 fListHistV0 = 0x0;
230 }
231 if (fTree){
232 delete fTree;
233 fTree = 0x0;
234 }
76029adc 235 //cleanup esd track cuts object too...
236 if (fESDtrackCuts){
237 delete fESDtrackCuts;
238 fESDtrackCuts = 0x0;
239 }
bcb6ffdb 240}
241
242//________________________________________________________________________
243void AliAnalysisTaskExtractPerformanceV0::UserCreateOutputObjects()
244{
bcb6ffdb 245
eff0b6fc 246 OpenFile(2);
bcb6ffdb 247 // Called once
c228a5af 248
249//------------------------------------------------
250
eff0b6fc 251 fTree = new TTree("fTree","V0Candidates");
bcb6ffdb 252
253//------------------------------------------------
c228a5af 254// fTree Branch definitions - V0 Tree
bcb6ffdb 255//------------------------------------------------
256
257//-----------BASIC-INFO---------------------------
258/* 1*/ fTree->Branch("fTreeVariablePrimaryStatus",&fTreeVariablePrimaryStatus,"fTreeVariablePrimaryStatus/I");
259/* 1*/ fTree->Branch("fTreeVariablePrimaryStatusMother",&fTreeVariablePrimaryStatusMother,"fTreeVariablePrimaryStatusMother/I");
260/* 2*/ fTree->Branch("fTreeVariableChi2V0",&fTreeVariableChi2V0,"Chi2V0/F");
261/* 3*/ fTree->Branch("fTreeVariableDcaV0Daughters",&fTreeVariableDcaV0Daughters,"fTreeVariableDcaV0Daughters/F");
262/* 4*/ fTree->Branch("fTreeVariableDcaPosToPrimVertex",&fTreeVariableDcaPosToPrimVertex,"fTreeVariableDcaPosToPrimVertex/F");
263/* 5*/ fTree->Branch("fTreeVariableDcaNegToPrimVertex",&fTreeVariableDcaNegToPrimVertex,"fTreeVariableDcaNegToPrimVertex/F");
264/* 6*/ fTree->Branch("fTreeVariableV0Radius",&fTreeVariableV0Radius,"fTreeVariableV0Radius/F");
265/* 7*/ fTree->Branch("fTreeVariablePt",&fTreeVariablePt,"fTreeVariablePt/F");
266/* 7*/ fTree->Branch("fTreeVariablePtMC",&fTreeVariablePtMC,"fTreeVariablePtMC/F");
267/* 8*/ fTree->Branch("fTreeVariableRapK0Short",&fTreeVariableRapK0Short,"fTreeVariableRapK0Short/F");
268/* 9*/ fTree->Branch("fTreeVariableRapLambda",&fTreeVariableRapLambda,"fTreeVariableRapLambda/F");
269/*10*/ fTree->Branch("fTreeVariableRapMC",&fTreeVariableRapMC,"fTreeVariableRapMC/F");
270/*11*/ fTree->Branch("fTreeVariableInvMassK0s",&fTreeVariableInvMassK0s,"fTreeVariableInvMassK0s/F");
271/*12*/ fTree->Branch("fTreeVariableInvMassLambda",&fTreeVariableInvMassLambda,"fTreeVariableInvMassLambda/F");
272/*13*/ fTree->Branch("fTreeVariableInvMassAntiLambda",&fTreeVariableInvMassAntiLambda,"fTreeVariableInvMassAntiLambda/F");
273/*14*/ fTree->Branch("fTreeVariableAlphaV0",&fTreeVariableAlphaV0,"fTreeVariableAlphaV0/F");
274/*15*/ fTree->Branch("fTreeVariablePtArmV0",&fTreeVariablePtArmV0,"fTreeVariablePtArmV0/F");
275/*16*/ fTree->Branch("fTreeVariableNegTransvMomentum",&fTreeVariableNegTransvMomentum,"fTreeVariableNegTransvMomentum/F");
276/*17*/ fTree->Branch("fTreeVariablePosTransvMomentum",&fTreeVariablePosTransvMomentum,"fTreeVariablePosTransvMomentum/F");
277/*18*/ fTree->Branch("fTreeVariableNegTransvMomentumMC",&fTreeVariableNegTransvMomentumMC,"fTreeVariableNegTransvMomentumMC/F");
278/*19*/ fTree->Branch("fTreeVariablePosTransvMomentumMC",&fTreeVariablePosTransvMomentumMC,"fTreeVariablePosTransvMomentumMC/F");
279/*20*/ fTree->Branch("fTreeVariableLeastNbrCrossedRows",&fTreeVariableLeastNbrCrossedRows,"fTreeVariableLeastNbrCrossedRows/I");
280/*21*/ fTree->Branch("fTreeVariableLeastRatioCrossedRowsOverFindable",&fTreeVariableLeastRatioCrossedRowsOverFindable,"fTreeVariableLeastRatioCrossedRowsOverFindable/F");
281/*22*/ fTree->Branch("fTreeVariablePID",&fTreeVariablePID,"fTreeVariablePID/I");
282/*23*/ fTree->Branch("fTreeVariablePIDPositive",&fTreeVariablePIDPositive,"fTreeVariablePIDPositive/I");
283/*24*/ fTree->Branch("fTreeVariablePIDNegative",&fTreeVariablePIDNegative,"fTreeVariablePIDNegative/I");
284/*25*/ fTree->Branch("fTreeVariablePIDMother",&fTreeVariablePIDMother,"fTreeVariablePIDMother/I");
285/*26*/ fTree->Branch("fTreeVariablePtXiMother",&fTreeVariablePtMother,"fTreeVariablePtMother/F");
286/*27*/ fTree->Branch("fTreeVariableV0CosineOfPointingAngle",&fTreeVariableV0CosineOfPointingAngle,"fTreeVariableV0CosineOfPointingAngle/F");
287//-----------MULTIPLICITY-INFO--------------------
288/*28*/ fTree->Branch("fTreeVariableMultiplicity",&fTreeVariableMultiplicity,"fTreeVariableMultiplicity/I");
289//------------------------------------------------
290/*29*/ fTree->Branch("fTreeVariableDistOverTotMom",&fTreeVariableDistOverTotMom,"fTreeVariableDistOverTotMom/F");
291/*30*/ fTree->Branch("fTreeVariableNSigmasPosProton",&fTreeVariableNSigmasPosProton,"fTreeVariableNSigmasPosProton/F");
292/*31*/ fTree->Branch("fTreeVariableNSigmasPosPion",&fTreeVariableNSigmasPosPion,"fTreeVariableNSigmasPosPion/F");
293/*32*/ fTree->Branch("fTreeVariableNSigmasNegProton",&fTreeVariableNSigmasNegProton,"fTreeVariableNSigmasNegProton/F");
294/*33*/ fTree->Branch("fTreeVariableNSigmasNegPion",&fTreeVariableNSigmasNegPion,"fTreeVariableNSigmasNegPion/F");
295//------------------------------------------------
296/*34*/ fTree->Branch("fTreeVariableNegEta",&fTreeVariableNegEta,"fTreeVariableNegEta/F");
297/*35*/ fTree->Branch("fTreeVariablePosEta",&fTreeVariablePosEta,"fTreeVariablePosEta/F");
298/*36*/ fTree->Branch("fTreeVariableV0CreationRadius",&fTreeVariableV0CreationRadius,"fTreeVariableV0CreationRadius/F");
299/*37*/ fTree->Branch("fTreeVariableIndexStatus",&fTreeVariableIndexStatus,"fTreeVariableIndexStatus/I");
300/*38*/ fTree->Branch("fTreeVariableIndexStatusMother",&fTreeVariableIndexStatusMother,"fTreeVariableIndexStatusMother/I");
301
0eebce93 302/*39*/ fTree->Branch("fTreeVariableRunNumber",&fTreeVariableRunNumber,"fTreeVariableRunNumber/I");
303/*40*/ fTree->Branch("fTreeVariableEventNumber",&fTreeVariableEventNumber,"fTreeVariableEventNumber/l");
304
bcb6ffdb 305//------------------------------------------------
306// Particle Identification Setup
307//------------------------------------------------
308
309 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
310 AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
311 fPIDResponse = inputHandler->GetPIDResponse();
312
76029adc 313 // Multiplicity
314
315 if(! fESDtrackCuts ){
316 fESDtrackCuts = new AliESDtrackCuts();
317 }
318
bcb6ffdb 319//------------------------------------------------
320// V0 Multiplicity Histograms
321//------------------------------------------------
322
eff0b6fc 323 // Create histograms
324 OpenFile(1);
325 fListHistV0 = new TList();
326 fListHistV0->SetOwner(); // See http://root.cern.ch/root/html/TCollection.html#TCollection:SetOwner
327
328
bcb6ffdb 329 if(! fHistV0MultiplicityBeforeTrigSel) {
330 fHistV0MultiplicityBeforeTrigSel = new TH1F("fHistV0MultiplicityBeforeTrigSel",
331 "V0s per event (before Trig. Sel.);Nbr of V0s/Evt;Events",
332 25, 0, 25);
333 fListHistV0->Add(fHistV0MultiplicityBeforeTrigSel);
334 }
335
336 if(! fHistV0MultiplicityForTrigEvt) {
337 fHistV0MultiplicityForTrigEvt = new TH1F("fHistV0MultiplicityForTrigEvt",
338 "V0s per event (for triggered evt);Nbr of V0s/Evt;Events",
339 25, 0, 25);
340 fListHistV0->Add(fHistV0MultiplicityForTrigEvt);
341 }
342
343 if(! fHistV0MultiplicityForSelEvt) {
344 fHistV0MultiplicityForSelEvt = new TH1F("fHistV0MultiplicityForSelEvt",
345 "V0s per event;Nbr of V0s/Evt;Events",
346 25, 0, 25);
347 fListHistV0->Add(fHistV0MultiplicityForSelEvt);
348 }
349
350 if(! fHistV0MultiplicityForSelEvtNoTPCOnly) {
351 fHistV0MultiplicityForSelEvtNoTPCOnly = new TH1F("fHistV0MultiplicityForSelEvtNoTPCOnly",
352 "V0s per event;Nbr of V0s/Evt;Events",
353 25, 0, 25);
354 fListHistV0->Add(fHistV0MultiplicityForSelEvtNoTPCOnly);
355 }
356 if(! fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup) {
357 fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup = new TH1F("fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup",
358 "V0s per event;Nbr of V0s/Evt;Events",
359 25, 0, 25);
360 fListHistV0->Add(fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup);
361 }
362
363//------------------------------------------------
364// Track Multiplicity Histograms
365//------------------------------------------------
366
367 if(! fHistMultiplicityBeforeTrigSel) {
368 fHistMultiplicityBeforeTrigSel = new TH1F("fHistMultiplicityBeforeTrigSel",
369 "Tracks per event;Nbr of Tracks;Events",
370 200, 0, 200);
371 fListHistV0->Add(fHistMultiplicityBeforeTrigSel);
372 }
373 if(! fHistMultiplicityForTrigEvt) {
374 fHistMultiplicityForTrigEvt = new TH1F("fHistMultiplicityForTrigEvt",
375 "Tracks per event;Nbr of Tracks;Events",
376 200, 0, 200);
377 fListHistV0->Add(fHistMultiplicityForTrigEvt);
378 }
379 if(! fHistMultiplicity) {
380 fHistMultiplicity = new TH1F("fHistMultiplicity",
381 "Tracks per event;Nbr of Tracks;Events",
382 200, 0, 200);
383 fListHistV0->Add(fHistMultiplicity);
384 }
385 if(! fHistMultiplicityNoTPCOnly) {
386 fHistMultiplicityNoTPCOnly = new TH1F("fHistMultiplicityNoTPCOnly",
387 "Tracks per event;Nbr of Tracks;Events",
388 200, 0, 200);
389 fListHistV0->Add(fHistMultiplicityNoTPCOnly);
390 }
391 if(! fHistMultiplicityNoTPCOnlyNoPileup) {
392 fHistMultiplicityNoTPCOnlyNoPileup = new TH1F("fHistMultiplicityNoTPCOnlyNoPileup",
393 "Tracks per event;Nbr of Tracks;Events",
394 200, 0, 200);
395 fListHistV0->Add(fHistMultiplicityNoTPCOnlyNoPileup);
396 }
397
398//------------------------------------------------
399// Generated Particle Histograms
400//------------------------------------------------
401
402 Int_t lCustomNBins = 200;
403 Double_t lCustomPtUpperLimit = 20;
404 Int_t lCustomNBinsMultiplicity = 100;
405
406//----------------------------------
407// Raw Generated (Pre-physics-selection)
408//----------------------------------
409
410//--- 3D Histo (Pt, Y, Multiplicity)
411
412 if(! f3dHistPrimRawPtVsYVsMultLambda) {
413 f3dHistPrimRawPtVsYVsMultLambda = new TH3F( "f3dHistPrimRawPtVsYVsMultLambda", "Pt_{lambda} Vs Y_{#Lambda} Vs Multiplicity; Pt_{lambda} (GeV/c); Y_{#Lambda} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
414 fListHistV0->Add(f3dHistPrimRawPtVsYVsMultLambda);
415 }
416 if(! f3dHistPrimRawPtVsYVsMultAntiLambda) {
417 f3dHistPrimRawPtVsYVsMultAntiLambda = new TH3F( "f3dHistPrimRawPtVsYVsMultAntiLambda", "Pt_{antilambda} Vs Y_{#Lambda} Vs Multiplicity; Pt_{antilambda} (GeV/c); Y_{#Lambda} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
418 fListHistV0->Add(f3dHistPrimRawPtVsYVsMultAntiLambda);
419 }
420 if(! f3dHistPrimRawPtVsYVsMultK0Short) {
421 f3dHistPrimRawPtVsYVsMultK0Short = new TH3F( "f3dHistPrimRawPtVsYVsMultK0Short", "Pt_{K0S} Vs Y_{K0S} Vs Multiplicity; Pt_{K0S} (GeV/c); Y_{K0S} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
422 fListHistV0->Add(f3dHistPrimRawPtVsYVsMultK0Short);
423 }
424
7e81bb32 425//--- 3D Histo (Pt, Y, Multiplicity), close to PV criterion
426
427 if(! f3dHistPrimCloseToPVPtVsYVsMultLambda) {
428 f3dHistPrimCloseToPVPtVsYVsMultLambda = new TH3F( "f3dHistPrimCloseToPVPtVsYVsMultLambda", "Pt_{lambda} Vs Y_{#Lambda} Vs Multiplicity; Pt_{lambda} (GeV/c); Y_{#Lambda} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
429 fListHistV0->Add(f3dHistPrimCloseToPVPtVsYVsMultLambda);
430 }
431 if(! f3dHistPrimCloseToPVPtVsYVsMultAntiLambda) {
432 f3dHistPrimCloseToPVPtVsYVsMultAntiLambda = new TH3F( "f3dHistPrimCloseToPVPtVsYVsMultAntiLambda", "Pt_{antilambda} Vs Y_{#Lambda} Vs Multiplicity; Pt_{antilambda} (GeV/c); Y_{#Lambda} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
433 fListHistV0->Add(f3dHistPrimCloseToPVPtVsYVsMultAntiLambda);
434 }
435 if(! f3dHistPrimCloseToPVPtVsYVsMultK0Short) {
436 f3dHistPrimCloseToPVPtVsYVsMultK0Short = new TH3F( "f3dHistPrimCloseToPVPtVsYVsMultK0Short", "Pt_{K0S} Vs Y_{K0S} Vs Multiplicity; Pt_{K0S} (GeV/c); Y_{K0S} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
437 fListHistV0->Add(f3dHistPrimCloseToPVPtVsYVsMultK0Short);
438 }
439
440
bcb6ffdb 441//--- 3D Histo (Pt, Y, Proper Decay Length)
442
443 if(! f3dHistPrimRawPtVsYVsDecayLengthLambda) {
444 f3dHistPrimRawPtVsYVsDecayLengthLambda = new TH3F( "f3dHistPrimRawPtVsYVsDecayLengthLambda", "Pt_{lambda} Vs Y_{#Lambda} Vs DecayLength; Pt_{lambda} (GeV/c); Y_{#Lambda} ; DecayLength", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,200,0,50);
445 fListHistV0->Add(f3dHistPrimRawPtVsYVsDecayLengthLambda);
446 }
447 if(! f3dHistPrimRawPtVsYVsDecayLengthAntiLambda) {
448 f3dHistPrimRawPtVsYVsDecayLengthAntiLambda = new TH3F( "f3dHistPrimRawPtVsYVsDecayLengthAntiLambda", "Pt_{antilambda} Vs Y_{#Lambda} Vs DecayLength; Pt_{antilambda} (GeV/c); Y_{#Lambda} ; DecayLength", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,200,0,50);
449 fListHistV0->Add(f3dHistPrimRawPtVsYVsDecayLengthAntiLambda);
450 }
451 if(! f3dHistPrimRawPtVsYVsDecayLengthK0Short) {
452 f3dHistPrimRawPtVsYVsDecayLengthK0Short = new TH3F( "f3dHistPrimRawPtVsYVsDecayLengthK0Short", "Pt_{K0S} Vs Y_{K0S} Vs DecayLength; Pt_{K0S} (GeV/c); Y_{K0S} ; DecayLength", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,200,0,50);
453 fListHistV0->Add(f3dHistPrimRawPtVsYVsDecayLengthK0Short);
454 }
455
c228a5af 456//--------------------------------------------------------------------------------------
457//--- 3D Histo (Pt, Y, Multiplicity) for generated XiMinus/Plus, all generated
bcb6ffdb 458
459 if(! f3dHistGenPtVsYVsMultXiMinus) {
460 f3dHistGenPtVsYVsMultXiMinus = new TH3F( "f3dHistGenPtVsYVsMultXiMinus", "Pt_{#Xi} Vs Y_{#Xi} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Xi} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
461 fListHistV0->Add(f3dHistGenPtVsYVsMultXiMinus);
462 }
463 if(! f3dHistGenPtVsYVsMultXiPlus) {
464 f3dHistGenPtVsYVsMultXiPlus = new TH3F( "f3dHistGenPtVsYVsMultXiPlus", "Pt_{#Xi} Vs Y_{#Xi} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Xi} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
465 fListHistV0->Add(f3dHistGenPtVsYVsMultXiPlus);
466 }
c228a5af 467//--- 3D Histo (Pt, Y, Multiplicity) for generated OmegaMinus/Plus
468
469 if(! f3dHistGenPtVsYVsMultOmegaMinus) {
470 f3dHistGenPtVsYVsMultOmegaMinus = new TH3F( "f3dHistGenPtVsYVsMultOmegaMinus", "Pt_{#Omega} Vs Y_{#Omega} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Omega} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
471 fListHistV0->Add(f3dHistGenPtVsYVsMultOmegaMinus);
472 }
473 if(! f3dHistGenPtVsYVsMultOmegaPlus) {
474 f3dHistGenPtVsYVsMultOmegaPlus = new TH3F( "f3dHistGenPtVsYVsMultOmegaPlus", "Pt_{#Omega} Vs Y_{#Omega} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Omega} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
475 fListHistV0->Add(f3dHistGenPtVsYVsMultOmegaPlus);
476 }
477
478//--------------------------------------------------------------------------------------
479//--- 3D Histo (Pt, Y, Multiplicity) for generated XiMinus/Plus, at selected analysis evts
480
481 if(! f3dHistGenSelectedPtVsYVsMultXiMinus) {
482 f3dHistGenSelectedPtVsYVsMultXiMinus = new TH3F( "f3dHistGenSelectedPtVsYVsMultXiMinus", "Pt_{#Xi} Vs Y_{#Xi} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Xi} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
483 fListHistV0->Add(f3dHistGenSelectedPtVsYVsMultXiMinus);
484 }
485 if(! f3dHistGenSelectedPtVsYVsMultXiPlus) {
486 f3dHistGenSelectedPtVsYVsMultXiPlus = new TH3F( "f3dHistGenSelectedPtVsYVsMultXiPlus", "Pt_{#Xi} Vs Y_{#Xi} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Xi} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
487 fListHistV0->Add(f3dHistGenSelectedPtVsYVsMultXiPlus);
488 }
489//--- 3D Histo (Pt, Y, Multiplicity) for generated OmegaMinus/Plus
490
491 if(! f3dHistGenSelectedPtVsYVsMultOmegaMinus) {
492 f3dHistGenSelectedPtVsYVsMultOmegaMinus = new TH3F( "f3dHistGenSelectedPtVsYVsMultOmegaMinus", "Pt_{#Omega} Vs Y_{#Omega} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Omega} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
493 fListHistV0->Add(f3dHistGenSelectedPtVsYVsMultOmegaMinus);
494 }
495 if(! f3dHistGenSelectedPtVsYVsMultOmegaPlus) {
496 f3dHistGenSelectedPtVsYVsMultOmegaPlus = new TH3F( "f3dHistGenSelectedPtVsYVsMultOmegaPlus", "Pt_{#Omega} Vs Y_{#Omega} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Omega} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
497 fListHistV0->Add(f3dHistGenSelectedPtVsYVsMultOmegaPlus);
498 }
499
bcb6ffdb 500
501//----------------------------------
502// Histos at analysis level
503//----------------------------------
504
505 if(! f3dHistPrimAnalysisPtVsYVsMultLambda) {
506 f3dHistPrimAnalysisPtVsYVsMultLambda = new TH3F( "f3dHistPrimAnalysisPtVsYVsMultLambda", "Pt_{lambda} Vs Y_{#Lambda} Vs Multiplicity; Pt_{lambda} (GeV/c); Y_{#Lambda} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
507 fListHistV0->Add(f3dHistPrimAnalysisPtVsYVsMultLambda);
508 }
509 if(! f3dHistPrimAnalysisPtVsYVsMultAntiLambda) {
510 f3dHistPrimAnalysisPtVsYVsMultAntiLambda = new TH3F( "f3dHistPrimAnalysisPtVsYVsMultAntiLambda", "Pt_{antilambda} Vs Y_{#Lambda} Vs Multiplicity; Pt_{antilambda} (GeV/c); Y_{#Lambda} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
511 fListHistV0->Add(f3dHistPrimAnalysisPtVsYVsMultAntiLambda);
512 }
513 if(! f3dHistPrimAnalysisPtVsYVsMultK0Short) {
514 f3dHistPrimAnalysisPtVsYVsMultK0Short = new TH3F( "f3dHistPrimAnalysisPtVsYVsMultK0Short", "Pt_{K0S} Vs Y_{K0S} Vs Multiplicity; Pt_{K0S} (GeV/c); Y_{K0S} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
515 fListHistV0->Add(f3dHistPrimAnalysisPtVsYVsMultK0Short);
516 }
517
518//----------------------------------
519// Primary Vertex Position Histos
520//----------------------------------
521
522 if(! fHistPVx) {
523 fHistPVx = new TH1F("fHistPVx",
524 "PV x position;Nbr of Evts;x",
525 2000, -0.5, 0.5);
526 fListHistV0->Add(fHistPVx);
527 }
528 if(! fHistPVy) {
529 fHistPVy = new TH1F("fHistPVy",
530 "PV y position;Nbr of Evts;y",
531 2000, -0.5, 0.5);
532 fListHistV0->Add(fHistPVy);
533 }
534 if(! fHistPVz) {
535 fHistPVz = new TH1F("fHistPVz",
536 "PV z position;Nbr of Evts;z",
537 400, -20, 20);
538 fListHistV0->Add(fHistPVz);
539 }
540
541 if(! fHistPVxAnalysis) {
542 fHistPVxAnalysis = new TH1F("fHistPVxAnalysis",
543 "PV x position;Nbr of Evts;x",
544 2000, -0.5, 0.5);
545 fListHistV0->Add(fHistPVxAnalysis);
546 }
547 if(! fHistPVyAnalysis) {
548 fHistPVyAnalysis = new TH1F("fHistPVyAnalysis",
549 "PV y position;Nbr of Evts;y",
550 2000, -0.5, 0.5);
551 fListHistV0->Add(fHistPVyAnalysis);
552 }
553 if(! fHistPVzAnalysis) {
554 fHistPVzAnalysis = new TH1F("fHistPVzAnalysis",
555 "PV z position;Nbr of Evts;z",
556 400, -20, 20);
557 fListHistV0->Add(fHistPVzAnalysis);
558 }
559
560 if(! fHistPVxAnalysisHasHighPtLambda) {
561 fHistPVxAnalysisHasHighPtLambda = new TH1F("fHistPVxAnalysisHasHighPtLambda",
562 "PV x position;Nbr of Evts;x",
563 2000, -0.5, 0.5);
564 fListHistV0->Add(fHistPVxAnalysisHasHighPtLambda);
565 }
566 if(! fHistPVyAnalysisHasHighPtLambda) {
567 fHistPVyAnalysisHasHighPtLambda = new TH1F("fHistPVyAnalysisHasHighPtLambda",
568 "PV y position;Nbr of Evts;y",
569 2000, -0.5, 0.5);
570 fListHistV0->Add(fHistPVyAnalysisHasHighPtLambda);
571 }
572 if(! fHistPVzAnalysisHasHighPtLambda) {
573 fHistPVzAnalysisHasHighPtLambda = new TH1F("fHistPVzAnalysisHasHighPtLambda",
574 "PV z position;Nbr of Evts;z",
575 400, -20, 20);
576 fListHistV0->Add(fHistPVzAnalysisHasHighPtLambda);
577 }
4f8ed367 578 if(! fHistSwappedV0Counter) {
579 fHistSwappedV0Counter = new TH1F("fHistSwappedV0Counter",
580 "Swap or not histo;Swapped (1) or not (0); count",
581 2, 0, 2);
582 fListHistV0->Add(fHistSwappedV0Counter);
583 }
584
bcb6ffdb 585 //List of Histograms: Normal
586 PostData(1, fListHistV0);
587
588 //TTree Object: Saved to base directory. Should cache to disk while saving.
589 //(Important to avoid excessive memory usage, particularly when merging)
590 PostData(2, fTree);
591
592}// end UserCreateOutputObjects
593
594
595//________________________________________________________________________
596void AliAnalysisTaskExtractPerformanceV0::UserExec(Option_t *)
597{
598 // Main loop
599 // Called for each event
600
601 AliESDEvent *lESDevent = 0x0;
602 AliMCEvent *lMCevent = 0x0;
603 AliStack *lMCstack = 0x0;
604
605 Int_t lNumberOfV0s = -1;
606 Double_t lTrkgPrimaryVtxPos[3] = {-100.0, -100.0, -100.0};
607 Double_t lBestPrimaryVtxPos[3] = {-100.0, -100.0, -100.0};
608 Double_t lMagneticField = -10.;
609
610 // Connect to the InputEvent
611 // After these lines, we should have an ESD/AOD event + the number of V0s in it.
612
613 // Appropriate for ESD analysis!
614
615 lESDevent = dynamic_cast<AliESDEvent*>( InputEvent() );
616 if (!lESDevent) {
617 AliWarning("ERROR: lESDevent not available \n");
618 return;
619 }
620
0eebce93 621 fTreeVariableRunNumber = lESDevent->GetRunNumber();
622 fTreeVariableEventNumber =
623 ( ( ((ULong64_t)lESDevent->GetPeriodNumber() ) << 36 ) |
624 ( ((ULong64_t)lESDevent->GetOrbitNumber () ) << 12 ) |
625 ((ULong64_t)lESDevent->GetBunchCrossNumber() ) );
626
bcb6ffdb 627 lMCevent = MCEvent();
628 if (!lMCevent) {
629 Printf("ERROR: Could not retrieve MC event \n");
630 cout << "Name of the file with pb :" << fInputHandler->GetTree()->GetCurrentFile()->GetName() << endl;
631 return;
632 }
633
634 lMCstack = lMCevent->Stack();
635 if (!lMCstack) {
636 Printf("ERROR: Could not retrieve MC stack \n");
637 cout << "Name of the file with pb :" << fInputHandler->GetTree()->GetCurrentFile()->GetName() << endl;
638 return;
639 }
7e81bb32 640 TArrayF mcPrimaryVtx;
641 AliGenEventHeader* mcHeader=lMCevent->GenEventHeader();
642 if(!mcHeader) return;
643 mcHeader->PrimaryVertex(mcPrimaryVtx);
bcb6ffdb 644
645//------------------------------------------------
646// Multiplicity Information Acquistion
647//------------------------------------------------
648
4f8ed367 649 //REVISED multiplicity estimator after 'multiplicity day' (2011)
650 Int_t lMultiplicity = -100;
651
ad0d3645 652 //testing purposes
76029adc 653 if(fkIsNuclear == kFALSE) lMultiplicity = fESDtrackCuts->GetReferenceMultiplicity(lESDevent, AliESDtrackCuts::kTrackletsITSTPC,0.5);
bcb6ffdb 654
548bd092 655 //---> If this is a nuclear collision, then go nuclear on "multiplicity" variable...
656 //---> Warning: Experimental
657 if(fkIsNuclear == kTRUE){
658 AliCentrality* centrality;
659 centrality = lESDevent->GetCentrality();
660 lMultiplicity = ( ( Int_t ) ( centrality->GetCentralityPercentile( "V0M" ) ) );
661 if (centrality->GetQuality()>1) {
662 PostData(1, fListHistV0);
663 PostData(2, fTree);
664 return;
665 }
666 }
667
668 //Set variable for filling tree afterwards!
669 //---> pp case......: GetReferenceMultiplicity
670 //---> Pb-Pb case...: Centrality by V0M
671 fTreeVariableMultiplicity = lMultiplicity;
672
bcb6ffdb 673 fHistV0MultiplicityBeforeTrigSel->Fill ( lESDevent->GetNumberOfV0s() );
674 fHistMultiplicityBeforeTrigSel->Fill ( lMultiplicity );
675
676//------------------------------------------------
677// MC Information Acquistion
678//------------------------------------------------
679
680 Int_t iNumberOfPrimaries = -1;
681 iNumberOfPrimaries = lMCstack->GetNprimary();
682 if(iNumberOfPrimaries < 1) return;
683 Bool_t lHasHighPtLambda = kFALSE;
684
685//------------------------------------------------
686// Variable Definition
687//------------------------------------------------
688
689 Int_t lNbMCPrimary = 0;
690
691 Int_t lPdgcodeCurrentPart = 0;
692 Double_t lRapCurrentPart = 0;
693 Double_t lPtCurrentPart = 0;
694
695 //Int_t lComeFromSigma = 0;
696
697 // current mc particle 's mother
698 //Int_t iCurrentMother = 0;
699 lNbMCPrimary = lMCstack->GetNprimary();
700
701//------------------------------------------------
702// Pre-Physics Selection
703//------------------------------------------------
704
705//----- Loop on primary Xi, Omega --------------------------------------------------------------
706 for (Int_t iCurrentLabelStack = 0; iCurrentLabelStack < lNbMCPrimary; iCurrentLabelStack++)
707 {// This is the begining of the loop on primaries
708
709 TParticle* lCurrentParticlePrimary = 0x0;
710 lCurrentParticlePrimary = lMCstack->Particle( iCurrentLabelStack );
711 if(!lCurrentParticlePrimary){
712 Printf("Cascade loop %d - MC TParticle pointer to current stack particle = 0x0 ! Skip ...\n", iCurrentLabelStack );
713 continue;
714 }
c228a5af 715 if ( TMath::Abs(lCurrentParticlePrimary->GetPdgCode()) == 3312 || TMath::Abs(lCurrentParticlePrimary->GetPdgCode()) == 3334 ) {
ad0d3645 716 Double_t lRapXiMCPrimary = -100;
717 if( (lCurrentParticlePrimary->Energy() - lCurrentParticlePrimary->Pz() +1.e-13) != 0 ) {
718 if ( (lCurrentParticlePrimary->Energy() + lCurrentParticlePrimary->Pz()) / (lCurrentParticlePrimary->Energy() - lCurrentParticlePrimary->Pz() +1.e-13) !=0 ){
719 lRapXiMCPrimary = 0.5*TMath::Log( (lCurrentParticlePrimary->Energy() + lCurrentParticlePrimary->Pz()) / (lCurrentParticlePrimary->Energy() - lCurrentParticlePrimary->Pz() +1.e-13) );
720 }
721 }
bcb6ffdb 722
723 //=================================================================================
c228a5af 724 // Xi Histograms
bcb6ffdb 725 if( lCurrentParticlePrimary->GetPdgCode() == 3312 ){
726 lPtCurrentPart = lCurrentParticlePrimary->Pt();
727 f3dHistGenPtVsYVsMultXiMinus->Fill(lPtCurrentPart, lRapXiMCPrimary, lMultiplicity);
728 }
729 if( lCurrentParticlePrimary->GetPdgCode() == -3312 ){
730 lPtCurrentPart = lCurrentParticlePrimary->Pt();
731 f3dHistGenPtVsYVsMultXiPlus->Fill(lPtCurrentPart, lRapXiMCPrimary, lMultiplicity);
732 }
c228a5af 733 // Omega Histograms
734 if( lCurrentParticlePrimary->GetPdgCode() == 3334 ){
735 lPtCurrentPart = lCurrentParticlePrimary->Pt();
736 f3dHistGenPtVsYVsMultOmegaMinus->Fill(lPtCurrentPart, lRapXiMCPrimary, lMultiplicity);
737 }
738 if( lCurrentParticlePrimary->GetPdgCode() == -3334 ){
739 lPtCurrentPart = lCurrentParticlePrimary->Pt();
740 f3dHistGenPtVsYVsMultOmegaPlus->Fill(lPtCurrentPart, lRapXiMCPrimary, lMultiplicity);
741 }
bcb6ffdb 742 }
743 }
744//----- End Loop on primary Xi, Omega ----------------------------------------------------------
745
746//----- Loop on Lambda, K0Short ----------------------------------------------------------------
747 for (Int_t iCurrentLabelStack = 0; iCurrentLabelStack < (lMCstack->GetNtrack()); iCurrentLabelStack++)
748 {// This is the begining of the loop on tracks
749
750 TParticle* lCurrentParticleForLambdaCheck = 0x0;
751 lCurrentParticleForLambdaCheck = lMCstack->Particle( iCurrentLabelStack );
752 if(!lCurrentParticleForLambdaCheck){
753 Printf("V0s loop %d - MC TParticle pointer to current stack particle = 0x0 ! Skip ...\n", iCurrentLabelStack );
754 continue;
755 }
756
757 //=================================================================================
758 //Single-Strange checks
759 // Keep only K0s, Lambda and AntiLambda:
760 lPdgcodeCurrentPart = lCurrentParticleForLambdaCheck->GetPdgCode();
761
762 if ( (lCurrentParticleForLambdaCheck->GetPdgCode() == 310 ) ||
763 (lCurrentParticleForLambdaCheck->GetPdgCode() == 3122 ) ||
764 (lCurrentParticleForLambdaCheck->GetPdgCode() == -3122 ) )
765 {
766 lRapCurrentPart = MyRapidity(lCurrentParticleForLambdaCheck->Energy(),lCurrentParticleForLambdaCheck->Pz());
767 lPtCurrentPart = lCurrentParticleForLambdaCheck->Pt();
768
7e81bb32 769 //Use Close to PV for filling CloseToPV histograms!
770 Double_t dx, dy, dz;
771
772 dx = ( (mcPrimaryVtx.At(0)) - (lCurrentParticleForLambdaCheck->Vx()) );
773 dy = ( (mcPrimaryVtx.At(1)) - (lCurrentParticleForLambdaCheck->Vy()) );
774 dz = ( (mcPrimaryVtx.At(2)) - (lCurrentParticleForLambdaCheck->Vz()) );
775 Double_t lDistToPV = TMath::Sqrt(dx*dx + dy*dy + dz*dz);
776 if( lDistToPV <= 0.001){
777 if( lPdgcodeCurrentPart == 3122 ){
778 f3dHistPrimCloseToPVPtVsYVsMultLambda->Fill(lPtCurrentPart, lRapCurrentPart, lMultiplicity);
779 }
780 if( lPdgcodeCurrentPart == -3122 ){
781 f3dHistPrimCloseToPVPtVsYVsMultAntiLambda->Fill(lPtCurrentPart, lRapCurrentPart, lMultiplicity);
782 }
783 if( lPdgcodeCurrentPart == 310 ){
784 f3dHistPrimCloseToPVPtVsYVsMultK0Short->Fill(lPtCurrentPart, lRapCurrentPart, lMultiplicity);
785 }
786 }
787
bcb6ffdb 788 //Use Physical Primaries only for filling PrimRaw Histograms!
789 if ( lMCstack->IsPhysicalPrimary(iCurrentLabelStack)!=kTRUE ) continue;
790
791 if( lPdgcodeCurrentPart == 3122 ){
792 f3dHistPrimRawPtVsYVsMultLambda->Fill(lPtCurrentPart, lRapCurrentPart, lMultiplicity);
793 if( TMath::Abs( lCurrentParticleForLambdaCheck->Eta() )<1.2 && lPtCurrentPart>2 ){
794 lHasHighPtLambda = kTRUE; //Keep track of events with Lambda within |eta|<1.2 and pt>2
795 }
796 }
797 if( lPdgcodeCurrentPart == -3122 ){
798 f3dHistPrimRawPtVsYVsMultAntiLambda->Fill(lPtCurrentPart, lRapCurrentPart, lMultiplicity);
799 }
800 if( lPdgcodeCurrentPart == 310 ){
801 f3dHistPrimRawPtVsYVsMultK0Short->Fill(lPtCurrentPart, lRapCurrentPart, lMultiplicity);
802 }
803 //Decay Length Acquisition=====================================================
804 Double_t decaylength = -1;
805 Double_t lV0Mass = -1;
806
807 if( !(lCurrentParticleForLambdaCheck->GetDaughter(0) < 0) ) {
808 TParticle* lDght0ofV0 = lMCstack->Particle( lCurrentParticleForLambdaCheck->GetDaughter(0) ); //get first daughter
809 if(lDght0ofV0){ // skip if not defined.
810 decaylength = TMath::Sqrt(
811 TMath::Power( lCurrentParticleForLambdaCheck->Vx() - lDght0ofV0->Vx() , 2) +
812 TMath::Power( lCurrentParticleForLambdaCheck->Vy() - lDght0ofV0->Vy() , 2) +
813 TMath::Power( lCurrentParticleForLambdaCheck->Vz() - lDght0ofV0->Vz() , 2)
814 );
815 //Need to correct for relativitity! Involves multiplying by mass and dividing by momentum.
816 if(TMath::Abs( lPdgcodeCurrentPart ) == 3122 ) { lV0Mass = 1.115683; }
817 if(TMath::Abs( lPdgcodeCurrentPart ) == 310 ) { lV0Mass = 0.497614; }
ad0d3645 818 if( lCurrentParticleForLambdaCheck->P() + 1e-10 != 0 ) decaylength = ( lV0Mass * decaylength ) / ( lCurrentParticleForLambdaCheck->P() + 1e-10 );
819 if( lCurrentParticleForLambdaCheck->P() + 1e-10 == 0 ) decaylength = 1e+5;
bcb6ffdb 820 }
821 }
822 if( lPdgcodeCurrentPart == 3122) f3dHistPrimRawPtVsYVsDecayLengthLambda ->Fill( lPtCurrentPart, lRapCurrentPart , decaylength );
823 if( lPdgcodeCurrentPart == -3122) f3dHistPrimRawPtVsYVsDecayLengthAntiLambda ->Fill( lPtCurrentPart, lRapCurrentPart , decaylength );
824 if( lPdgcodeCurrentPart == 310) f3dHistPrimRawPtVsYVsDecayLengthK0Short ->Fill( lPtCurrentPart, lRapCurrentPart , decaylength );
825 }
826 }//End of loop on tracks
827//----- End Loop on Lambda, K0Short ------------------------------------------------------------
828
829// ---> Set Variables to Zero again
830// ---> Variable Definition
831
832 lPdgcodeCurrentPart = 0;
833 lRapCurrentPart = 0;
834 lPtCurrentPart = 0;
835
836//------------------------------------------------
837// Physics Selection
838//------------------------------------------------
839
840 UInt_t maskIsSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
841 Bool_t isSelected = 0;
842 isSelected = (maskIsSelected & AliVEvent::kMB) == AliVEvent::kMB;
548bd092 843
844 //pp at 2.76TeV: special case, ignore FastOnly
845 if ( (fkLowEnergyPP == kTRUE) && (maskIsSelected& AliVEvent::kFastOnly) ){
846 PostData(1, fListHistV0);
847 PostData(2, fTree);
848 return;
849 }
850 //Standard Min-Bias Selection
bcb6ffdb 851 if ( ! isSelected ) {
852 PostData(1, fListHistV0);
853 PostData(2, fTree);
854 return;
855 }
856
857//------------------------------------------------
858// After Trigger Selection
859//------------------------------------------------
860
861 lNumberOfV0s = lESDevent->GetNumberOfV0s();
862
863 //Set variable for filling tree afterwards!
bcb6ffdb 864 fHistV0MultiplicityForTrigEvt->Fill(lNumberOfV0s);
865 fHistMultiplicityForTrigEvt->Fill ( lMultiplicity );
866
867//------------------------------------------------
868// Getting: Primary Vertex + MagField Info
869//------------------------------------------------
870
871 const AliESDVertex *lPrimaryTrackingESDVtx = lESDevent->GetPrimaryVertexTracks();
872 // get the vtx stored in ESD found with tracks
873 lPrimaryTrackingESDVtx->GetXYZ( lTrkgPrimaryVtxPos );
874
875 const AliESDVertex *lPrimaryBestESDVtx = lESDevent->GetPrimaryVertex();
876 // get the best primary vertex available for the event
877 // As done in AliCascadeVertexer, we keep the one which is the best one available.
878 // between : Tracking vertex > SPD vertex > TPC vertex > default SPD vertex
879 // This one will be used for next calculations (DCA essentially)
880 lPrimaryBestESDVtx->GetXYZ( lBestPrimaryVtxPos );
881
882 Double_t lPrimaryVtxPosition[3];
883 const AliVVertex *primaryVtx = lESDevent->GetPrimaryVertex();
884 lPrimaryVtxPosition[0] = primaryVtx->GetX();
885 lPrimaryVtxPosition[1] = primaryVtx->GetY();
886 lPrimaryVtxPosition[2] = primaryVtx->GetZ();
887 fHistPVx->Fill( lPrimaryVtxPosition[0] );
888 fHistPVy->Fill( lPrimaryVtxPosition[1] );
889 fHistPVz->Fill( lPrimaryVtxPosition[2] );
890
891//------------------------------------------------
892// Primary Vertex Z position: SKIP
893//------------------------------------------------
894
895 if(TMath::Abs(lBestPrimaryVtxPos[2]) > 10.0 ) {
896 AliWarning("Pb / | Z position of Best Prim Vtx | > 10.0 cm ... return !");
897 PostData(1, fListHistV0);
898 PostData(2, fTree);
899 return;
900 }
901
902 lMagneticField = lESDevent->GetMagneticField( );
903 fHistV0MultiplicityForSelEvt ->Fill( lNumberOfV0s );
904 fHistMultiplicity->Fill(lMultiplicity);
905
906//------------------------------------------------
907// SKIP: Events with well-established PVtx
908//------------------------------------------------
909
910 const AliESDVertex *lPrimaryTrackingESDVtxCheck = lESDevent->GetPrimaryVertexTracks();
911 const AliESDVertex *lPrimarySPDVtx = lESDevent->GetPrimaryVertexSPD();
912 if (!lPrimarySPDVtx->GetStatus() && !lPrimaryTrackingESDVtxCheck->GetStatus() ){
913 AliWarning("Pb / No SPD prim. vertex nor prim. Tracking vertex ... return !");
914 PostData(1, fListHistV0);
915 PostData(2, fTree);
916 return;
917 }
918 fHistV0MultiplicityForSelEvtNoTPCOnly ->Fill( lNumberOfV0s );
919 fHistMultiplicityNoTPCOnly->Fill(lMultiplicity);
920
921//------------------------------------------------
922// Pileup Rejection Studies
923//------------------------------------------------
924
925 // FIXME : quality selection regarding pile-up rejection
4f8ed367 926 if(lESDevent->IsPileupFromSPD() && !fkIsNuclear ){// minContributors=3, minZdist=0.8, nSigmaZdist=3., nSigmaDiamXY=2., nSigmaDiamZ=5. -> see http://alisoft.cern.ch/viewvc/trunk/STEER/AliESDEvent.h?root=AliRoot&r1=41914&r2=42199&pathrev=42199
bcb6ffdb 927 AliWarning("Pb / This is tagged as Pileup from SPD... return !");
928 PostData(1, fListHistV0);
929 PostData(2, fTree);
930 return;
931 }
932 fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup ->Fill( lNumberOfV0s );
933 fHistMultiplicityNoTPCOnlyNoPileup->Fill(lMultiplicity);
934
935 //Do control histograms without the IsFromVertexerZ events, but consider them in analysis...
936 if( ! (lESDevent->GetPrimaryVertex()->IsFromVertexerZ() ) ){
937 fHistPVxAnalysis->Fill( lPrimaryVtxPosition[0] );
938 fHistPVyAnalysis->Fill( lPrimaryVtxPosition[1] );
939 fHistPVzAnalysis->Fill( lPrimaryVtxPosition[2] );
940 if ( lHasHighPtLambda == kTRUE ){
941 fHistPVxAnalysisHasHighPtLambda->Fill( lPrimaryVtxPosition[0] );
942 fHistPVyAnalysisHasHighPtLambda->Fill( lPrimaryVtxPosition[1] );
943 fHistPVzAnalysisHasHighPtLambda->Fill( lPrimaryVtxPosition[2] );
944 }
945 }
946
947//------------------------------------------------
948// stack loop starts here
949//------------------------------------------------
950
951//---> Loop over ALL PARTICLES
952
953 for (Int_t iMc = 0; iMc < (lMCstack->GetNtrack()); iMc++) {
954 TParticle *p0 = lMCstack->Particle(iMc);
955 if (!p0) {
956 //Printf("ERROR: particle with label %d not found in lMCstack (mc loop)", iMc);
957 continue;
958 }
959 lPdgcodeCurrentPart = p0->GetPdgCode();
960
961 // Keep only K0s, Lambda and AntiLambda:
962 if ( (lPdgcodeCurrentPart != 310 ) && (lPdgcodeCurrentPart != 3122 ) && (lPdgcodeCurrentPart != -3122 ) ) continue;
963
964 lRapCurrentPart = MyRapidity(p0->Energy(),p0->Pz());
965 lPtCurrentPart = p0->Pt();
966
967 //Use Physical Primaries only for filling PrimRaw Histograms!
968 if ( lMCstack->IsPhysicalPrimary(iMc)!=kTRUE ) continue;
969
970 if( lPdgcodeCurrentPart == 3122 ){
971 f3dHistPrimAnalysisPtVsYVsMultLambda->Fill(lPtCurrentPart, lRapCurrentPart, lMultiplicity);
972 }
973 if( lPdgcodeCurrentPart == -3122 ){
974 f3dHistPrimAnalysisPtVsYVsMultAntiLambda->Fill(lPtCurrentPart, lRapCurrentPart, lMultiplicity);
975 }
976 if( lPdgcodeCurrentPart == 310 ){
977 f3dHistPrimAnalysisPtVsYVsMultK0Short->Fill(lPtCurrentPart, lRapCurrentPart, lMultiplicity);
978 }
979 }
980
c228a5af 981//----- Loop on primary Xi, Omega --------------------------------------------------------------
982 for (Int_t iCurrentLabelStack = 0; iCurrentLabelStack < lNbMCPrimary; iCurrentLabelStack++)
983 {// This is the begining of the loop on primaries
984
985 TParticle* lCurrentParticlePrimary = 0x0;
986 lCurrentParticlePrimary = lMCstack->Particle( iCurrentLabelStack );
987 if(!lCurrentParticlePrimary){
988 Printf("Cascade loop %d - MC TParticle pointer to current stack particle = 0x0 ! Skip ...\n", iCurrentLabelStack );
989 continue;
990 }
991 if ( TMath::Abs(lCurrentParticlePrimary->GetPdgCode()) == 3312 || TMath::Abs(lCurrentParticlePrimary->GetPdgCode()) == 3334 ) {
992 Double_t lRapXiMCPrimary = -100;
993 if( (lCurrentParticlePrimary->Energy() - lCurrentParticlePrimary->Pz() +1.e-13) != 0 ) {
994 if ( (lCurrentParticlePrimary->Energy() + lCurrentParticlePrimary->Pz()) / (lCurrentParticlePrimary->Energy() - lCurrentParticlePrimary->Pz() +1.e-13) !=0 ){
995 lRapXiMCPrimary = 0.5*TMath::Log( (lCurrentParticlePrimary->Energy() + lCurrentParticlePrimary->Pz()) / (lCurrentParticlePrimary->Energy() - lCurrentParticlePrimary->Pz() +1.e-13) );
996 }
997 }
998
999 //=================================================================================
1000 // Xi Histograms
1001 if( lCurrentParticlePrimary->GetPdgCode() == 3312 ){
1002 lPtCurrentPart = lCurrentParticlePrimary->Pt();
1003 f3dHistGenSelectedPtVsYVsMultXiMinus->Fill(lPtCurrentPart, lRapXiMCPrimary, lMultiplicity);
1004 }
1005 if( lCurrentParticlePrimary->GetPdgCode() == -3312 ){
1006 lPtCurrentPart = lCurrentParticlePrimary->Pt();
1007 f3dHistGenSelectedPtVsYVsMultXiPlus->Fill(lPtCurrentPart, lRapXiMCPrimary, lMultiplicity);
1008 }
1009 // Omega Histograms
1010 if( lCurrentParticlePrimary->GetPdgCode() == 3334 ){
1011 lPtCurrentPart = lCurrentParticlePrimary->Pt();
1012 f3dHistGenSelectedPtVsYVsMultOmegaMinus->Fill(lPtCurrentPart, lRapXiMCPrimary, lMultiplicity);
1013 }
1014 if( lCurrentParticlePrimary->GetPdgCode() == -3334 ){
1015 lPtCurrentPart = lCurrentParticlePrimary->Pt();
1016 f3dHistGenSelectedPtVsYVsMultOmegaPlus->Fill(lPtCurrentPart, lRapXiMCPrimary, lMultiplicity);
1017 }
1018 }
1019 }
1020//----- End Loop on primary Xi, Omega ----------------------------------------------------------
1021
bcb6ffdb 1022//------------------------------------------------
1023// MAIN LAMBDA LOOP STARTS HERE
1024//------------------------------------------------
1025
1026 //Variable definition
1027 Int_t lOnFlyStatus = 0;
1028 Double_t lChi2V0 = 0;
1029 Double_t lDcaV0Daughters = 0, lDcaV0ToPrimVertex = 0;
1030 Double_t lDcaPosToPrimVertex = 0, lDcaNegToPrimVertex = 0;
1031 Double_t lV0CosineOfPointingAngle = 0;
1032 Double_t lV0Radius = 0, lPt = 0;
1033 Double_t lRapK0Short = 0, lRapLambda = 0;
1034 Double_t lInvMassK0s = 0, lInvMassLambda = 0, lInvMassAntiLambda = 0;
1035 Double_t lAlphaV0 = 0, lPtArmV0 = 0;
1036 Double_t fMinV0Pt = 0;
1037 Double_t fMaxV0Pt = 100;
1038
1039 Int_t nv0s = 0;
1040 nv0s = lESDevent->GetNumberOfV0s();
1041
1042 for (Int_t iV0 = 0; iV0 < nv0s; iV0++)
1043 {// This is the begining of the V0 loop
1044 AliESDv0 *v0 = ((AliESDEvent*)lESDevent)->GetV0(iV0);
1045 if (!v0) continue;
1046
4f8ed367 1047 //---> Fix On-the-Fly candidates
1048 if( v0->GetParamN()->Charge() > 0 && v0->GetParamP()->Charge() < 0 ){
1049 fHistSwappedV0Counter -> Fill( 1 );
1050 }else{
1051 fHistSwappedV0Counter -> Fill( 0 );
1052 }
1053 if ( fkUseOnTheFly ) CheckChargeV0(v0);
1054
1055
bcb6ffdb 1056 Double_t tV0mom[3];
1057 v0->GetPxPyPz( tV0mom[0],tV0mom[1],tV0mom[2] );
1058 Double_t lV0TotalMomentum = TMath::Sqrt(
1059 tV0mom[0]*tV0mom[0]+tV0mom[1]*tV0mom[1]+tV0mom[2]*tV0mom[2] );
1060
1061 Double_t tDecayVertexV0[3]; v0->GetXYZ(tDecayVertexV0[0],tDecayVertexV0[1],tDecayVertexV0[2]);
1062 lV0Radius = TMath::Sqrt(tDecayVertexV0[0]*tDecayVertexV0[0]+tDecayVertexV0[1]*tDecayVertexV0[1]);
1063 lPt = v0->Pt();
1064 lRapK0Short = v0->RapK0Short();
1065 lRapLambda = v0->RapLambda();
1066 if ((lPt<fMinV0Pt)||(fMaxV0Pt<lPt)) continue;
1067
1068 UInt_t lKeyPos = (UInt_t)TMath::Abs(v0->GetPindex());
1069 UInt_t lKeyNeg = (UInt_t)TMath::Abs(v0->GetNindex());
1070
1071 Double_t lMomPos[3]; v0->GetPPxPyPz(lMomPos[0],lMomPos[1],lMomPos[2]);
1072 Double_t lMomNeg[3]; v0->GetNPxPyPz(lMomNeg[0],lMomNeg[1],lMomNeg[2]);
1073
1074 AliESDtrack *pTrack=((AliESDEvent*)lESDevent)->GetTrack(lKeyPos);
1075 AliESDtrack *nTrack=((AliESDEvent*)lESDevent)->GetTrack(lKeyNeg);
1076 if (!pTrack || !nTrack) {
1077 Printf("ERROR: Could not retreive one of the daughter track");
1078 continue;
1079 }
1080
1081 fTreeVariableNegEta = nTrack->Eta();
1082 fTreeVariablePosEta = pTrack->Eta();
1083
1084 // Filter like-sign V0 (next: add counter and distribution)
1085 if ( pTrack->GetSign() == nTrack->GetSign()){
1086 continue;
1087 }
1088
1089 //________________________________________________________________________
1090 // Track quality cuts
1091 Float_t lPosTrackCrossedRows = pTrack->GetTPCClusterInfo(2,1);
1092 Float_t lNegTrackCrossedRows = nTrack->GetTPCClusterInfo(2,1);
4f8ed367 1093 fTreeVariableLeastNbrCrossedRows = (Int_t) lPosTrackCrossedRows;
bcb6ffdb 1094 if( lNegTrackCrossedRows < fTreeVariableLeastNbrCrossedRows )
4f8ed367 1095 fTreeVariableLeastNbrCrossedRows = (Int_t) lNegTrackCrossedRows;
bcb6ffdb 1096
1097 // TPC refit condition (done during reconstruction for Offline but not for On-the-fly)
1098 if( !(pTrack->GetStatus() & AliESDtrack::kTPCrefit)) continue;
1099 if( !(nTrack->GetStatus() & AliESDtrack::kTPCrefit)) continue;
1100
1101 if ( ( ( pTrack->GetTPCClusterInfo(2,1) ) < 70 ) || ( ( nTrack->GetTPCClusterInfo(2,1) ) < 70 ) ) continue;
1102
1103 //GetKinkIndex condition
1104 if( pTrack->GetKinkIndex(0)>0 || nTrack->GetKinkIndex(0)>0 ) continue;
1105
1106 //Findable clusters > 0 condition
1107 if( pTrack->GetTPCNclsF()<=0 || nTrack->GetTPCNclsF()<=0 ) continue;
1108
1109 //Compute ratio Crossed Rows / Findable clusters
1110 //Note: above test avoids division by zero!
ad0d3645 1111 Float_t lPosTrackCrossedRowsOverFindable = -1;
1112 Float_t lNegTrackCrossedRowsOverFindable = -1;
1113 if ( ((double)(pTrack->GetTPCNclsF()) ) != 0 ) lPosTrackCrossedRowsOverFindable = lPosTrackCrossedRows / ((double)(pTrack->GetTPCNclsF()));
1114 if ( ((double)(nTrack->GetTPCNclsF()) ) != 0 ) lNegTrackCrossedRowsOverFindable = lNegTrackCrossedRows / ((double)(nTrack->GetTPCNclsF()));
bcb6ffdb 1115
1116 fTreeVariableLeastRatioCrossedRowsOverFindable = lPosTrackCrossedRowsOverFindable;
1117 if( lNegTrackCrossedRowsOverFindable < fTreeVariableLeastRatioCrossedRowsOverFindable )
1118 fTreeVariableLeastRatioCrossedRowsOverFindable = lNegTrackCrossedRowsOverFindable;
1119
1120 //Lowest Cut Level for Ratio Crossed Rows / Findable = 0.8, set here
1121 if ( fTreeVariableLeastRatioCrossedRowsOverFindable < 0.8 ) continue;
1122
1123 //End track Quality Cuts
1124 //________________________________________________________________________
1125
1126 lDcaPosToPrimVertex = TMath::Abs(pTrack->GetD(lPrimaryVtxPosition[0],
1127 lPrimaryVtxPosition[1],
1128 lMagneticField) );
1129
1130 lDcaNegToPrimVertex = TMath::Abs(nTrack->GetD(lPrimaryVtxPosition[0],
1131 lPrimaryVtxPosition[1],
1132 lMagneticField) );
1133
1134 lOnFlyStatus = v0->GetOnFlyStatus();
1135 lChi2V0 = v0->GetChi2V0();
1136 lDcaV0Daughters = v0->GetDcaV0Daughters();
1137 lDcaV0ToPrimVertex = v0->GetD(lPrimaryVtxPosition[0],lPrimaryVtxPosition[1],lPrimaryVtxPosition[2]);
1138 lV0CosineOfPointingAngle = v0->GetV0CosineOfPointingAngle(lPrimaryVtxPosition[0],lPrimaryVtxPosition[1],lPrimaryVtxPosition[2]);
1139 fTreeVariableV0CosineOfPointingAngle=lV0CosineOfPointingAngle;
1140
1141 // Getting invariant mass infos directly from ESD
1142 v0->ChangeMassHypothesis(310);
1143 lInvMassK0s = v0->GetEffMass();
1144 v0->ChangeMassHypothesis(3122);
1145 lInvMassLambda = v0->GetEffMass();
1146 v0->ChangeMassHypothesis(-3122);
1147 lInvMassAntiLambda = v0->GetEffMass();
1148 lAlphaV0 = v0->AlphaV0();
1149 lPtArmV0 = v0->PtArmV0();
1150
1151 //fTreeVariableOnFlyStatus = lOnFlyStatus;
1152 //fHistV0OnFlyStatus->Fill(lOnFlyStatus);
1153
1154//===============================================
1155// Monte Carlo Association starts here
1156//===============================================
1157
1158 //---> Set Everything to "I don't know" before starting
1159
1160 fTreeVariablePIDPositive = 0;
1161 fTreeVariablePIDNegative = 0;
1162
1163 fTreeVariableIndexStatus = 0;
1164 fTreeVariableIndexStatusMother = 0;
1165
1166 fTreeVariablePtMother = -1;
1167 fTreeVariablePtMC = -1;
1168 fTreeVariableRapMC = -100;
1169
1170 fTreeVariablePID = -1;
1171 fTreeVariablePIDMother = -1;
1172
1173 fTreeVariablePrimaryStatus = 0;
1174 fTreeVariablePrimaryStatusMother = 0;
1175 fTreeVariableV0CreationRadius = -1;
1176
1177 Int_t lblPosV0Dghter = (Int_t) TMath::Abs( pTrack->GetLabel() );
1178 Int_t lblNegV0Dghter = (Int_t) TMath::Abs( nTrack->GetLabel() );
1179
1180 TParticle* mcPosV0Dghter = lMCstack->Particle( lblPosV0Dghter );
1181 TParticle* mcNegV0Dghter = lMCstack->Particle( lblNegV0Dghter );
1182
1183 fTreeVariablePosTransvMomentumMC = mcPosV0Dghter->Pt();
1184 fTreeVariableNegTransvMomentumMC = mcNegV0Dghter->Pt();
1185
1186 Int_t lPIDPositive = mcPosV0Dghter -> GetPdgCode();
1187 Int_t lPIDNegative = mcNegV0Dghter -> GetPdgCode();
1188
1189 fTreeVariablePIDPositive = lPIDPositive;
1190 fTreeVariablePIDNegative = lPIDNegative;
1191
1192 Int_t lblMotherPosV0Dghter = mcPosV0Dghter->GetFirstMother() ;
1193 Int_t lblMotherNegV0Dghter = mcNegV0Dghter->GetFirstMother();
1194
1195 if( lblMotherPosV0Dghter == lblMotherNegV0Dghter && lblMotherPosV0Dghter > -1 ){
1196 //either label is fine, they're equal at this stage
1197 TParticle* pThisV0 = lMCstack->Particle( lblMotherPosV0Dghter );
1198 //Set tree variables
1199 fTreeVariablePID = pThisV0->GetPdgCode(); //PDG Code
1200 fTreeVariablePtMC = pThisV0->Pt(); //Perfect Pt
1201 //Only Interested if it's a Lambda, AntiLambda or K0s
1202 //Avoid the Junction Bug! PYTHIA has particles with Px=Py=Pz=E=0 occasionally,
1203 //having particle code 88 (unrecognized by PDG), for documentation purposes.
1204 //Even ROOT's TParticle::Y() is not prepared to deal with that exception!
1205 //Note that TParticle::Pt() is immune (that would just return 0)...
1206 //Though granted that that should be extremely rare in this precise condition...
1207 if( TMath::Abs(fTreeVariablePID) == 3122 || fTreeVariablePID==310 ){
1208 fTreeVariableRapMC = pThisV0->Y(); //Perfect Y
1209 }
7e81bb32 1210 fTreeVariableV0CreationRadius = TMath::Sqrt(
1211 TMath::Power( ( (mcPrimaryVtx.At(0)) - (pThisV0->Vx()) ) , 2) +
1212 TMath::Power( ( (mcPrimaryVtx.At(1)) - (pThisV0->Vy()) ) , 2) +
1213 TMath::Power( ( (mcPrimaryVtx.At(2)) - (pThisV0->Vz()) ) , 2)
1214 );
bcb6ffdb 1215 if( lblMotherPosV0Dghter < lNbMCPrimary ) fTreeVariableIndexStatus = 1; //looks primary
1216 if( lblMotherPosV0Dghter >= lNbMCPrimary ) fTreeVariableIndexStatus = 2; //looks secondary
1217 if( lMCstack->IsPhysicalPrimary (lblMotherPosV0Dghter) ) fTreeVariablePrimaryStatus = 1; //Is Primary!
1218 if( lMCstack->IsSecondaryFromWeakDecay(lblMotherPosV0Dghter) ) fTreeVariablePrimaryStatus = 2; //Weak Decay!
1219 if( lMCstack->IsSecondaryFromMaterial (lblMotherPosV0Dghter) ) fTreeVariablePrimaryStatus = 3; //Material Int!
1220
1221 //Now we try to acquire the V0 parent particle, if possible
1222 Int_t lblThisV0Parent = pThisV0->GetFirstMother();
1223 if ( lblThisV0Parent > -1 ){ //if it has a parent, get it and store specs
1224 TParticle* pThisV0Parent = lMCstack->Particle( lblThisV0Parent );
1225 fTreeVariablePIDMother = pThisV0Parent->GetPdgCode(); //V0 Mother PDG
1226 fTreeVariablePtMother = pThisV0Parent->Pt(); //V0 Mother Pt
1227 //Primary Status for the V0 Mother particle
1228 if( lblThisV0Parent < lNbMCPrimary ) fTreeVariableIndexStatusMother = 1; //looks primary
1229 if( lblThisV0Parent >= lNbMCPrimary ) fTreeVariableIndexStatusMother = 2; //looks secondary
1230 if( lMCstack->IsPhysicalPrimary (lblThisV0Parent) ) fTreeVariablePrimaryStatusMother = 1; //Is Primary!
1231 if( lMCstack->IsSecondaryFromWeakDecay(lblThisV0Parent) ) fTreeVariablePrimaryStatusMother = 2; //Weak Decay!
1232 if( lMCstack->IsSecondaryFromMaterial (lblThisV0Parent) ) fTreeVariablePrimaryStatusMother = 3; //Material Int!
1233 }
1234 }
1235
1236 fTreeVariablePt = v0->Pt();
1237 fTreeVariableChi2V0 = lChi2V0;
1238 fTreeVariableDcaV0ToPrimVertex = lDcaV0ToPrimVertex;
1239 fTreeVariableDcaV0Daughters = lDcaV0Daughters;
1240 fTreeVariableV0CosineOfPointingAngle = lV0CosineOfPointingAngle;
1241 fTreeVariableV0Radius = lV0Radius;
1242 fTreeVariableDcaPosToPrimVertex = lDcaPosToPrimVertex;
1243 fTreeVariableDcaNegToPrimVertex = lDcaNegToPrimVertex;
1244 fTreeVariableInvMassK0s = lInvMassK0s;
1245 fTreeVariableInvMassLambda = lInvMassLambda;
1246 fTreeVariableInvMassAntiLambda = lInvMassAntiLambda;
1247 fTreeVariableRapK0Short = lRapK0Short;
1248
1249 fTreeVariableRapLambda = lRapLambda;
1250 fTreeVariableAlphaV0 = lAlphaV0;
1251 fTreeVariablePtArmV0 = lPtArmV0;
1252
1253 //Official means of acquiring N-sigmas
1254 fTreeVariableNSigmasPosProton = fPIDResponse->NumberOfSigmasTPC( pTrack, AliPID::kProton );
1255 fTreeVariableNSigmasPosPion = fPIDResponse->NumberOfSigmasTPC( pTrack, AliPID::kPion );
1256 fTreeVariableNSigmasNegProton = fPIDResponse->NumberOfSigmasTPC( nTrack, AliPID::kProton );
1257 fTreeVariableNSigmasNegPion = fPIDResponse->NumberOfSigmasTPC( nTrack, AliPID::kPion );
1258
1259//tDecayVertexV0[0],tDecayVertexV0[1],tDecayVertexV0[2]
ab9d4ee1 1260 Double_t lDistanceTravelled = TMath::Sqrt(
bcb6ffdb 1261 TMath::Power( tDecayVertexV0[0] - lBestPrimaryVtxPos[0] , 2) +
1262 TMath::Power( tDecayVertexV0[1] - lBestPrimaryVtxPos[1] , 2) +
1263 TMath::Power( tDecayVertexV0[2] - lBestPrimaryVtxPos[2] , 2)
1264 );
ad0d3645 1265 fTreeVariableDistOverTotMom = 1e+5;
ab9d4ee1 1266 if( lV0TotalMomentum + 1e-10 != 0 ) fTreeVariableDistOverTotMom = lDistanceTravelled / (lV0TotalMomentum + 1e-10); //avoid division by zero, to be sure
bcb6ffdb 1267
1268 Double_t lMomentumPosTemp[3];
1269 pTrack->GetPxPyPz(lMomentumPosTemp);
1270 Double_t lPtPosTemporary = sqrt(pow(lMomentumPosTemp[0],2) + pow(lMomentumPosTemp[1],2));
1271
1272 Double_t lMomentumNegTemp[3];
1273 nTrack->GetPxPyPz(lMomentumNegTemp);
1274 Double_t lPtNegTemporary = sqrt(pow(lMomentumNegTemp[0],2) + pow(lMomentumNegTemp[1],2));
1275
1276 fTreeVariablePosTransvMomentum = lPtPosTemporary;
1277 fTreeVariableNegTransvMomentum = lPtNegTemporary;
1278
1279
1280//------------------------------------------------
1281// Fill Tree!
1282//------------------------------------------------
1283
1284 // The conditionals are meant to decrease excessive
1285 // memory usage!
1286
1287 //Modified version: Keep only OnFlyStatus == 0
1288 //Keep only if included in a parametric InvMass Region 20 sigmas away from peak
1289
1290 //First Selection: Reject OnFly
4f8ed367 1291 if( (lOnFlyStatus == 0 && fkUseOnTheFly == kFALSE) || (lOnFlyStatus != 0 && fkUseOnTheFly == kTRUE ) ){
bcb6ffdb 1292 //Second Selection: rough 20-sigma band, parametric.
1293 //K0Short: Enough to parametrize peak broadening with linear function.
1294 Double_t lUpperLimitK0Short = (5.63707e-01) + (1.14979e-02)*fTreeVariablePt;
1295 Double_t lLowerLimitK0Short = (4.30006e-01) - (1.10029e-02)*fTreeVariablePt;
1296 //Lambda: Linear (for higher pt) plus exponential (for low-pt broadening)
1297 //[0]+[1]*x+[2]*TMath::Exp(-[3]*x)
1298 Double_t lUpperLimitLambda = (1.13688e+00) + (5.27838e-03)*fTreeVariablePt + (8.42220e-02)*TMath::Exp(-(3.80595e+00)*fTreeVariablePt);
1299 Double_t lLowerLimitLambda = (1.09501e+00) - (5.23272e-03)*fTreeVariablePt - (7.52690e-02)*TMath::Exp(-(3.46339e+00)*fTreeVariablePt);
1300 //Do Selection
1301 if( (fTreeVariableInvMassLambda < lUpperLimitLambda && fTreeVariableInvMassLambda > lLowerLimitLambda ) ||
1302 (fTreeVariableInvMassAntiLambda < lUpperLimitLambda && fTreeVariableInvMassAntiLambda > lLowerLimitLambda ) ||
1303 (fTreeVariableInvMassK0s < lUpperLimitK0Short && fTreeVariableInvMassK0s > lLowerLimitK0Short ) ){
4f8ed367 1304 //Pre-selection in case this is AA...
1305 if( fkIsNuclear == kFALSE ) fTree->Fill();
1306 if( fkIsNuclear == kTRUE){
1307 //If this is a nuclear collision___________________
1308 // ... pre-filter with daughter eta selection only (not TPC)
1309 if ( TMath::Abs(fTreeVariableNegEta)<0.8 && TMath::Abs(fTreeVariablePosEta)<0.8 ) fTree->Fill();
1310 }//end nuclear_____________________________________
bcb6ffdb 1311 }
1312 }
1313
1314//------------------------------------------------
1315// Fill tree over.
1316//------------------------------------------------
1317
1318
1319 }// This is the end of the V0 loop
1320
c228a5af 1321//------------------------------------------------
1322
bcb6ffdb 1323 // Post output data.
1324 PostData(1, fListHistV0);
1325 PostData(2, fTree);
1326}
1327
1328//________________________________________________________________________
1329void AliAnalysisTaskExtractPerformanceV0::Terminate(Option_t *)
1330{
1331 // Draw result to the screen
1332 // Called once at the end of the query
1333
1334 TList *cRetrievedList = 0x0;
1335 cRetrievedList = (TList*)GetOutputData(1);
1336 if(!cRetrievedList){
1337 Printf("ERROR - AliAnalysisTaskExtractV0 : ouput data container list not available\n");
1338 return;
1339 }
1340
1341 fHistV0MultiplicityForTrigEvt = dynamic_cast<TH1F*> ( cRetrievedList->FindObject("fHistV0MultiplicityForTrigEvt") );
1342 if (!fHistV0MultiplicityForTrigEvt) {
1343 Printf("ERROR - AliAnalysisTaskExtractV0 : fHistV0MultiplicityForTrigEvt not available");
1344 return;
1345 }
1346
1347 TCanvas *canCheck = new TCanvas("AliAnalysisTaskExtractV0","V0 Multiplicity",10,10,510,510);
1348 canCheck->cd(1)->SetLogy();
1349
1350 fHistV0MultiplicityForTrigEvt->SetMarkerStyle(22);
1351 fHistV0MultiplicityForTrigEvt->DrawCopy("E");
1352}
1353
1354//----------------------------------------------------------------------------
1355
1356Double_t AliAnalysisTaskExtractPerformanceV0::MyRapidity(Double_t rE, Double_t rPz) const
1357{
1358 // Local calculation for rapidity
ad0d3645 1359 Double_t ReturnValue = -100;
1360 if( (rE-rPz+1.e-13) != 0 && (rE+rPz) != 0 ){
1361 ReturnValue = 0.5*TMath::Log((rE+rPz)/(rE-rPz+1.e-13));
1362 }
1363 return ReturnValue;
bcb6ffdb 1364}
4f8ed367 1365
1366//________________________________________________________________________
1367void AliAnalysisTaskExtractPerformanceV0::CheckChargeV0(AliESDv0 *v0)
1368{
1369 // This function checks charge of negative and positive daughter tracks.
1370 // If incorrectly defined (onfly vertexer), swaps out.
1371 if( v0->GetParamN()->Charge() > 0 && v0->GetParamP()->Charge() < 0 ){
1372 //V0 daughter track swapping is required! Note: everything is swapped here... P->N, N->P
1373 Long_t lCorrectNidx = v0->GetPindex();
1374 Long_t lCorrectPidx = v0->GetNindex();
1375 Double32_t lCorrectNmom[3];
1376 Double32_t lCorrectPmom[3];
1377 v0->GetPPxPyPz( lCorrectNmom[0], lCorrectNmom[1], lCorrectNmom[2] );
1378 v0->GetNPxPyPz( lCorrectPmom[0], lCorrectPmom[1], lCorrectPmom[2] );
1379
1380 AliExternalTrackParam lCorrectParamN(
1381 v0->GetParamP()->GetX() ,
1382 v0->GetParamP()->GetAlpha() ,
1383 v0->GetParamP()->GetParameter() ,
1384 v0->GetParamP()->GetCovariance()
1385 );
1386 AliExternalTrackParam lCorrectParamP(
1387 v0->GetParamN()->GetX() ,
1388 v0->GetParamN()->GetAlpha() ,
1389 v0->GetParamN()->GetParameter() ,
1390 v0->GetParamN()->GetCovariance()
1391 );
1392 lCorrectParamN.SetMostProbablePt( v0->GetParamP()->GetMostProbablePt() );
1393 lCorrectParamP.SetMostProbablePt( v0->GetParamN()->GetMostProbablePt() );
1394
1395 //Get Variables___________________________________________________
1396 Double_t lDcaV0Daughters = v0 -> GetDcaV0Daughters();
1397 Double_t lCosPALocal = v0 -> GetV0CosineOfPointingAngle();
1398 Bool_t lOnFlyStatusLocal = v0 -> GetOnFlyStatus();
1399
1400 //Create Replacement Object_______________________________________
1401 AliESDv0 *v0correct = new AliESDv0(lCorrectParamN,lCorrectNidx,lCorrectParamP,lCorrectPidx);
1402 v0correct->SetDcaV0Daughters ( lDcaV0Daughters );
1403 v0correct->SetV0CosineOfPointingAngle ( lCosPALocal );
1404 v0correct->ChangeMassHypothesis ( kK0Short );
1405 v0correct->SetOnFlyStatus ( lOnFlyStatusLocal );
1406
1407 //Reverse Cluster info..._________________________________________
1408 v0correct->SetClusters( v0->GetClusters( 1 ), v0->GetClusters ( 0 ) );
1409
1410 *v0 = *v0correct;
1411 //Proper cleanup..._______________________________________________
1412 v0correct->Delete();
1413 v0correct = 0x0;
1414
1415 //Just another cross-check and output_____________________________
1416 if( v0->GetParamN()->Charge() > 0 && v0->GetParamP()->Charge() < 0 ) {
1417 AliWarning("Found Swapped Charges, tried to correct but something FAILED!");
1418 }else{
1419 //AliWarning("Found Swapped Charges and fixed.");
1420 }
1421 //________________________________________________________________
1422 }else{
1423 //Don't touch it! ---
1424 //Printf("Ah, nice. Charges are already ordered...");
1425 }
1426 return;
1427}
bcb6ffdb 1428