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