1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
16 // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
18 // This task is meant to explore the possibility of using a VZERO amplitude
19 // based multiplicity estimator for proton-proton collisions. For this, two
20 // main operation methods for this task are foreseen:
22 // 1) (under development) it should act as an auxiliary task and provide a
23 // calibrated estimator
25 // 2) "Debug mode" which will also create a ROOT TTree object with event
26 // by event info potentially used for exploration / calibration. This
27 // includes the following info:
29 // --- All VZERO Amplitudes (saved as Float_t)
30 // --- (optional) time for each channel
31 // --- (optional) time width for each channel
32 // --- GetReferenceMultiplicity Estimator, |eta|<0.5
33 // --- GetReferenceMultiplicity Estimator, |eta|<0.8
34 // --- (if MC) True Multiplicity, |eta|<0.5
35 // --- (if MC) True Multiplicity, 2.8 < eta < 5.8 (VZEROA region)
36 // --- (if MC) True Multiplicity, -3.7 < eta <-1.7 (VZEROC region)
39 // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
45 //class AliMCEventHandler;
54 #include <Riostream.h>
60 #include "THnSparse.h"
67 #include "AliESDEvent.h"
68 #include "AliAODEvent.h"
69 #include "AliV0vertexer.h"
70 #include "AliCascadeVertexer.h"
71 #include "AliESDpid.h"
72 #include "AliESDtrack.h"
73 #include "AliESDtrackCuts.h"
74 #include "AliInputEventHandler.h"
75 #include "AliAnalysisManager.h"
76 #include "AliMCEventHandler.h"
77 #include "AliMCEvent.h"
79 #include "AliCentrality.h"
81 #include "AliCFContainer.h"
82 #include "AliMultiplicity.h"
83 #include "AliAODMCParticle.h"
84 #include "AliESDcascade.h"
85 #include "AliAODcascade.h"
86 #include "AliESDUtils.h"
87 #include "AliGenEventHeader.h"
88 #include "AliAnalysisTaskSE.h"
89 #include "AliAnalysisUtils.h"
90 #include "AliAnalysisTaskStrangenessVsMultiplicity.h"
95 ClassImp(AliAnalysisTaskStrangenessVsMultiplicity)
97 AliAnalysisTaskStrangenessVsMultiplicity::AliAnalysisTaskStrangenessVsMultiplicity()
98 : AliAnalysisTaskSE(), fListHist(0), fTreeEvent(0), fTreeV0(0), fTreeCascade(0), fPIDResponse(0), fESDtrackCuts(0),
99 fkSaveV0Tree ( kFALSE ),
100 fkSaveCascadeTree ( kTRUE ),
101 fkRunVertexers ( kTRUE ),
102 //---> Variables for fTreeEvent
105 fAmplitude_V0AEq (0),
106 fAmplitude_V0CEq (0),
110 fCentrality_V0AEq(0),
111 fCentrality_V0CEq(0),
112 fCentrality_V0MEq(0),
116 //---> Variables for fTreeV0
117 fTreeVariableChi2V0(0),
118 fTreeVariableDcaV0Daughters(0),
119 fTreeVariableDcaV0ToPrimVertex(0),
120 fTreeVariableDcaPosToPrimVertex(0),
121 fTreeVariableDcaNegToPrimVertex(0),
122 fTreeVariableV0CosineOfPointingAngle(0),
123 fTreeVariableV0Radius(0),
125 fTreeVariableRapK0Short(0),
126 fTreeVariableRapLambda(0),
127 fTreeVariableInvMassK0s(0),
128 fTreeVariableInvMassLambda(0),
129 fTreeVariableInvMassAntiLambda(0),
130 fTreeVariableAlphaV0(0),
131 fTreeVariablePtArmV0(0),
132 fTreeVariableNegEta(0),
133 fTreeVariablePosEta(0),
135 fTreeVariableNSigmasPosProton(0),
136 fTreeVariableNSigmasPosPion(0),
137 fTreeVariableNSigmasNegProton(0),
138 fTreeVariableNSigmasNegPion(0),
140 fTreeVariableDistOverTotMom(0),
141 fTreeVariableLeastNbrCrossedRows(0),
142 fTreeVariableLeastRatioCrossedRowsOverFindable(0),
144 fTreeVariableCentV0A(0),
145 fTreeVariableCentV0C(0),
146 fTreeVariableCentV0M(0),
147 fTreeVariableCentV0AEq(0),
148 fTreeVariableCentV0CEq(0),
149 fTreeVariableCentV0MEq(0),
150 fTreeVariableAmpV0A(0),
151 fTreeVariableAmpV0C(0),
152 fTreeVariableAmpV0AEq(0),
153 fTreeVariableAmpV0CEq(0),
154 fTreeVariableRefMultEta8(0),
155 fTreeVariableRefMultEta5(0),
156 fTreeVariableRunNumber(0),
157 //---> Variables for fTreeCascade
158 fTreeCascVarCharge(0),
159 fTreeCascVarMassAsXi(0),
160 fTreeCascVarMassAsOmega(0),
162 fTreeCascVarRapXi(0),
163 fTreeCascVarRapOmega(0),
164 fTreeCascVarNegEta(0),
165 fTreeCascVarPosEta(0),
166 fTreeCascVarBachEta(0),
167 fTreeCascVarDCACascDaughters(0),
168 fTreeCascVarDCABachToPrimVtx(0),
169 fTreeCascVarDCAV0Daughters(0),
170 fTreeCascVarDCAV0ToPrimVtx(0),
171 fTreeCascVarDCAPosToPrimVtx(0),
172 fTreeCascVarDCANegToPrimVtx(0),
173 fTreeCascVarCascCosPointingAngle(0),
174 fTreeCascVarCascRadius(0),
175 fTreeCascVarV0Mass(0),
176 fTreeCascVarV0CosPointingAngle(0),
177 fTreeCascVarV0CosPointingAngleSpecial(0),
178 fTreeCascVarV0Radius(0),
179 fTreeCascVarLeastNbrClusters(0),
180 fTreeCascVarDistOverTotMom(0),
181 fTreeCascVarNegNSigmaPion(0),
182 fTreeCascVarNegNSigmaProton(0),
183 fTreeCascVarPosNSigmaPion(0),
184 fTreeCascVarPosNSigmaProton(0),
185 fTreeCascVarBachNSigmaPion(0),
186 fTreeCascVarBachNSigmaKaon(0),
187 fTreeCascVarCentV0A(0),
188 fTreeCascVarCentV0C(0),
189 fTreeCascVarCentV0M(0),
190 fTreeCascVarCentV0AEq(0),
191 fTreeCascVarCentV0CEq(0),
192 fTreeCascVarCentV0MEq(0),
193 fTreeCascVarAmpV0A(0),
194 fTreeCascVarAmpV0C(0),
195 fTreeCascVarAmpV0AEq(0),
196 fTreeCascVarAmpV0CEq(0),
197 fTreeCascVarRefMultEta8(0),
198 fTreeCascVarRefMultEta5(0),
199 fTreeCascVarRunNumber(0),
202 //------------------------------------------------
208 AliAnalysisTaskStrangenessVsMultiplicity::AliAnalysisTaskStrangenessVsMultiplicity(const char *name)
209 : AliAnalysisTaskSE(name), fListHist(0), fTreeEvent(0), fTreeV0(0), fTreeCascade(0), fPIDResponse(0), fESDtrackCuts(0),
210 fkSaveV0Tree ( kFALSE ),
211 fkSaveCascadeTree ( kTRUE ),
212 fkRunVertexers ( kTRUE ),
213 //---> Variables for fTreeEvent
216 fAmplitude_V0AEq (0),
217 fAmplitude_V0CEq (0),
221 fCentrality_V0AEq(0),
222 fCentrality_V0CEq(0),
223 fCentrality_V0MEq(0),
227 //---> Variables for fTreeV0
228 fTreeVariableChi2V0(0),
229 fTreeVariableDcaV0Daughters(0),
230 fTreeVariableDcaV0ToPrimVertex(0),
231 fTreeVariableDcaPosToPrimVertex(0),
232 fTreeVariableDcaNegToPrimVertex(0),
233 fTreeVariableV0CosineOfPointingAngle(0),
234 fTreeVariableV0Radius(0),
236 fTreeVariableRapK0Short(0),
237 fTreeVariableRapLambda(0),
238 fTreeVariableInvMassK0s(0),
239 fTreeVariableInvMassLambda(0),
240 fTreeVariableInvMassAntiLambda(0),
241 fTreeVariableAlphaV0(0),
242 fTreeVariablePtArmV0(0),
243 fTreeVariableNegEta(0),
244 fTreeVariablePosEta(0),
246 fTreeVariableNSigmasPosProton(0),
247 fTreeVariableNSigmasPosPion(0),
248 fTreeVariableNSigmasNegProton(0),
249 fTreeVariableNSigmasNegPion(0),
251 fTreeVariableDistOverTotMom(0),
252 fTreeVariableLeastNbrCrossedRows(0),
253 fTreeVariableLeastRatioCrossedRowsOverFindable(0),
255 fTreeVariableCentV0A(0),
256 fTreeVariableCentV0C(0),
257 fTreeVariableCentV0M(0),
258 fTreeVariableCentV0AEq(0),
259 fTreeVariableCentV0CEq(0),
260 fTreeVariableCentV0MEq(0),
261 fTreeVariableAmpV0A(0),
262 fTreeVariableAmpV0C(0),
263 fTreeVariableAmpV0AEq(0),
264 fTreeVariableAmpV0CEq(0),
265 fTreeVariableRefMultEta8(0),
266 fTreeVariableRefMultEta5(0),
267 fTreeVariableRunNumber(0),
268 //---> Variables for fTreeCascade
269 fTreeCascVarCharge(0),
270 fTreeCascVarMassAsXi(0),
271 fTreeCascVarMassAsOmega(0),
273 fTreeCascVarRapXi(0),
274 fTreeCascVarRapOmega(0),
275 fTreeCascVarNegEta(0),
276 fTreeCascVarPosEta(0),
277 fTreeCascVarBachEta(0),
278 fTreeCascVarDCACascDaughters(0),
279 fTreeCascVarDCABachToPrimVtx(0),
280 fTreeCascVarDCAV0Daughters(0),
281 fTreeCascVarDCAV0ToPrimVtx(0),
282 fTreeCascVarDCAPosToPrimVtx(0),
283 fTreeCascVarDCANegToPrimVtx(0),
284 fTreeCascVarCascCosPointingAngle(0),
285 fTreeCascVarCascRadius(0),
286 fTreeCascVarV0Mass(0),
287 fTreeCascVarV0CosPointingAngle(0),
288 fTreeCascVarV0CosPointingAngleSpecial(0),
289 fTreeCascVarV0Radius(0),
290 fTreeCascVarLeastNbrClusters(0),
291 fTreeCascVarDistOverTotMom(0),
292 fTreeCascVarNegNSigmaPion(0),
293 fTreeCascVarNegNSigmaProton(0),
294 fTreeCascVarPosNSigmaPion(0),
295 fTreeCascVarPosNSigmaProton(0),
296 fTreeCascVarBachNSigmaPion(0),
297 fTreeCascVarBachNSigmaKaon(0),
298 fTreeCascVarCentV0A(0),
299 fTreeCascVarCentV0C(0),
300 fTreeCascVarCentV0M(0),
301 fTreeCascVarCentV0AEq(0),
302 fTreeCascVarCentV0CEq(0),
303 fTreeCascVarCentV0MEq(0),
304 fTreeCascVarAmpV0A(0),
305 fTreeCascVarAmpV0C(0),
306 fTreeCascVarAmpV0AEq(0),
307 fTreeCascVarAmpV0CEq(0),
308 fTreeCascVarRefMultEta8(0),
309 fTreeCascVarRefMultEta5(0),
310 fTreeCascVarRunNumber(0),
315 //Re-vertex: Will only apply for cascade candidates
317 fV0VertexerSels[0] = 33. ; // max allowed chi2
318 fV0VertexerSels[1] = 0.02; // min allowed impact parameter for the 1st daughter (LHC09a4 : 0.05)
319 fV0VertexerSels[2] = 0.02; // min allowed impact parameter for the 2nd daughter (LHC09a4 : 0.05)
320 fV0VertexerSels[4] = 0.95; // min allowed cosine of V0's pointing angle (LHC09a4 : 0.99)
321 fV0VertexerSels[5] = 1.0 ; // min radius of the fiducial volume (LHC09a4 : 0.2)
322 fV0VertexerSels[6] = 200. ; // max radius of the fiducial volume (LHC09a4 : 100.0)
324 fCascadeVertexerSels[0] = 33. ; // max allowed chi2 (same as PDC07)
325 fCascadeVertexerSels[1] = 0.05 ; // min allowed V0 impact parameter (PDC07 : 0.05 / LHC09a4 : 0.025 )
326 fCascadeVertexerSels[2] = 0.010; // "window" around the Lambda mass (PDC07 : 0.008 / LHC09a4 : 0.010 )
327 fCascadeVertexerSels[3] = 0.03 ; // min allowed bachelor's impact parameter (PDC07 : 0.035 / LHC09a4 : 0.025 )
328 fCascadeVertexerSels[4] = 2.0 ; // max allowed DCA between the V0 and the bachelor (PDC07 : 0.1 / LHC09a4 : 0.2 )
329 fCascadeVertexerSels[5] = 0.95 ; // min allowed cosine of the cascade pointing angle (PDC07 : 0.9985 / LHC09a4 : 0.998 )
330 fCascadeVertexerSels[6] = 0.4 ; // min radius of the fiducial volume (PDC07 : 0.9 / LHC09a4 : 0.2 )
331 fCascadeVertexerSels[7] = 100. ; // max radius of the fiducial volume (PDC07 : 100 / LHC09a4 : 100 )
334 DefineOutput(1, TList::Class()); // Event Counter Histo
335 DefineOutput(2, TTree::Class()); // Event Tree
336 DefineOutput(3, TTree::Class()); // V0 Tree
337 DefineOutput(4, TTree::Class()); // Cascade Tree
341 AliAnalysisTaskStrangenessVsMultiplicity::~AliAnalysisTaskStrangenessVsMultiplicity()
343 //------------------------------------------------
345 //------------------------------------------------
365 //________________________________________________________________________
366 void AliAnalysisTaskStrangenessVsMultiplicity::UserCreateOutputObjects()
372 //------------------------------------------------
374 fTreeEvent = new TTree("fTreeEvent","Event");
376 //------------------------------------------------
377 // fTree Branch definitions - Event by Event info
378 //------------------------------------------------
380 //-----------BASIC-INFO---------------------------
382 //--- VZERO Data (Integrated)
383 fTreeEvent->Branch("fAmplitude_V0A",&fAmplitude_V0A,"fAmplitude_V0A/F");
384 fTreeEvent->Branch("fAmplitude_V0C",&fAmplitude_V0C,"fAmplitude_V0C/F");
385 fTreeEvent->Branch("fAmplitude_V0AEq",&fAmplitude_V0AEq,"fAmplitude_V0AEq/F");
386 fTreeEvent->Branch("fAmplitude_V0CEq",&fAmplitude_V0CEq,"fAmplitude_V0CEq/F");
388 //Info from AliCentrality (not necessarily 'centrality' per se)
389 fTreeEvent->Branch("fCentrality_V0A",&fCentrality_V0A,"fCentrality_V0A/F");
390 fTreeEvent->Branch("fCentrality_V0C",&fCentrality_V0C,"fCentrality_V0C/F");
391 fTreeEvent->Branch("fCentrality_V0M",&fCentrality_V0A,"fCentrality_V0M/F");
392 fTreeEvent->Branch("fCentrality_V0AEq",&fCentrality_V0AEq,"fCentrality_V0AEq/F");
393 fTreeEvent->Branch("fCentrality_V0CEq",&fCentrality_V0CEq,"fCentrality_V0CEq/F");
394 fTreeEvent->Branch("fCentrality_V0MEq",&fCentrality_V0AEq,"fCentrality_V0MEq/F");
396 //Official GetReferenceMultiplicity
397 fTreeEvent->Branch("fRefMultEta5",&fRefMultEta5,"fRefMultEta5/I");
398 fTreeEvent->Branch("fRefMultEta8",&fRefMultEta8,"fRefMultEta8/I");
401 fTreeEvent->Branch("fRunNumber", &fRunNumber, "fRunNumber/I");
403 //Create Basic V0 Output Tree
404 fTreeV0 = new TTree( "fTreeV0", "V0 Candidates");
406 //------------------------------------------------
407 // fTreeV0 Branch definitions
408 //------------------------------------------------
410 //-----------BASIC-INFO---------------------------
411 fTreeV0->Branch("fTreeVariableChi2V0",&fTreeVariableChi2V0,"fTreeVariableChi2V0/F");
412 fTreeV0->Branch("fTreeVariableDcaV0Daughters",&fTreeVariableDcaV0Daughters,"fTreeVariableDcaV0Daughters/F");
413 fTreeV0->Branch("fTreeVariableDcaV0ToPrimVertex",&fTreeVariableDcaV0ToPrimVertex,"fTreeVariableDcaV0ToPrimVertex/F");
414 fTreeV0->Branch("fTreeVariableDcaPosToPrimVertex",&fTreeVariableDcaPosToPrimVertex,"fTreeVariableDcaPosToPrimVertex/F");
415 fTreeV0->Branch("fTreeVariableDcaNegToPrimVertex",&fTreeVariableDcaNegToPrimVertex,"fTreeVariableDcaNegToPrimVertex/F");
416 fTreeV0->Branch("fTreeVariableV0Radius",&fTreeVariableV0Radius,"fTreeVariableV0Radius/F");
417 fTreeV0->Branch("fTreeVariablePt",&fTreeVariablePt,"fTreeVariablePt/F");
418 fTreeV0->Branch("fTreeVariableRapK0Short",&fTreeVariableRapK0Short,"fTreeVariableRapK0Short/F");
419 fTreeV0->Branch("fTreeVariableRapLambda",&fTreeVariableRapLambda,"fTreeVariableRapLambda/F");
420 fTreeV0->Branch("fTreeVariableInvMassK0s",&fTreeVariableInvMassK0s,"fTreeVariableInvMassK0s/F");
421 fTreeV0->Branch("fTreeVariableInvMassLambda",&fTreeVariableInvMassLambda,"fTreeVariableInvMassLambda/F");
422 fTreeV0->Branch("fTreeVariableInvMassAntiLambda",&fTreeVariableInvMassAntiLambda,"fTreeVariableInvMassAntiLambda/F");
423 fTreeV0->Branch("fTreeVariableV0CosineOfPointingAngle",&fTreeVariableV0CosineOfPointingAngle,"fTreeVariableV0CosineOfPointingAngle/F");
424 fTreeV0->Branch("fTreeVariableAlphaV0",&fTreeVariableAlphaV0,"fTreeVariableAlphaV0/F");
425 fTreeV0->Branch("fTreeVariablePtArmV0",&fTreeVariablePtArmV0,"fTreeVariablePtArmV0/F");
426 fTreeV0->Branch("fTreeVariableLeastNbrCrossedRows",&fTreeVariableLeastNbrCrossedRows,"fTreeVariableLeastNbrCrossedRows/I");
427 fTreeV0->Branch("fTreeVariableLeastRatioCrossedRowsOverFindable",&fTreeVariableLeastRatioCrossedRowsOverFindable,"fTreeVariableLeastRatioCrossedRowsOverFindable/F");
428 fTreeV0->Branch("fTreeVariableDistOverTotMom",&fTreeVariableDistOverTotMom,"fTreeVariableDistOverTotMom/F");
429 fTreeV0->Branch("fTreeVariableNSigmasPosProton",&fTreeVariableNSigmasPosProton,"fTreeVariableNSigmasPosProton/F");
430 fTreeV0->Branch("fTreeVariableNSigmasPosPion",&fTreeVariableNSigmasPosPion,"fTreeVariableNSigmasPosPion/F");
431 fTreeV0->Branch("fTreeVariableNSigmasNegProton",&fTreeVariableNSigmasNegProton,"fTreeVariableNSigmasNegProton/F");
432 fTreeV0->Branch("fTreeVariableNSigmasNegPion",&fTreeVariableNSigmasNegPion,"fTreeVariableNSigmasNegPion/F");
433 fTreeV0->Branch("fTreeVariableNegEta",&fTreeVariableNegEta,"fTreeVariableNegEta/F");
434 fTreeV0->Branch("fTreeVariablePosEta",&fTreeVariablePosEta,"fTreeVariablePosEta/F");
435 //-----------MULTIPLICITY-INFO--------------------
436 fTreeV0->Branch("fTreeVariableCentV0A",&fTreeVariableCentV0A,"fTreeVariableCentV0A/F");
437 fTreeV0->Branch("fTreeVariableCentV0C",&fTreeVariableCentV0C,"fTreeVariableCentV0C/F");
438 fTreeV0->Branch("fTreeVariableCentV0M",&fTreeVariableCentV0M,"fTreeVariableCentV0M/F");
439 fTreeV0->Branch("fTreeVariableCentV0AEq",&fTreeVariableCentV0AEq,"fTreeVariableCentV0AEq/F");
440 fTreeV0->Branch("fTreeVariableCentV0CEq",&fTreeVariableCentV0CEq,"fTreeVariableCentV0CEq/F");
441 fTreeV0->Branch("fTreeVariableCentV0MEq",&fTreeVariableCentV0MEq,"fTreeVariableCentV0MEq/F");
442 fTreeV0->Branch("fTreeVariableAmpV0A",&fTreeVariableAmpV0A,"fTreeVariableAmpV0A/F");
443 fTreeV0->Branch("fTreeVariableAmpV0C",&fTreeVariableAmpV0C,"fTreeVariableAmpV0C/F");
444 fTreeV0->Branch("fTreeVariableAmpV0AEq",&fTreeVariableAmpV0AEq,"fTreeVariableAmpV0AEq/F");
445 fTreeV0->Branch("fTreeVariableAmpV0CEq",&fTreeVariableAmpV0CEq,"fTreeVariableAmpV0CEq/F");
446 fTreeV0->Branch("fTreeVariableRefMultEta8",&fTreeVariableRefMultEta8,"fTreeVariableRefMultEta8/I");
447 fTreeV0->Branch("fTreeVariableRefMultEta5",&fTreeVariableRefMultEta5,"fTreeVariableRefMultEta5/I");
448 fTreeV0->Branch("fTreeVariableRunNumber",&fTreeVariableRunNumber,"fTreeVariableRunNumber/I");
449 //------------------------------------------------
451 //Create Cascade output tree
452 fTreeCascade = new TTree("fTreeCascade","CascadeCandidates");
454 //------------------------------------------------
455 // fTreeCascade Branch definitions - Cascade Tree
456 //------------------------------------------------
458 //-----------BASIC-INFO---------------------------
459 fTreeCascade->Branch("fTreeCascVarCharge",&fTreeCascVarCharge,"fTreeCascVarCharge/I");
460 fTreeCascade->Branch("fTreeCascVarMassAsXi",&fTreeCascVarMassAsXi,"fTreeCascVarMassAsXi/F");
461 fTreeCascade->Branch("fTreeCascVarMassAsOmega",&fTreeCascVarMassAsOmega,"fTreeCascVarMassAsOmega/F");
462 fTreeCascade->Branch("fTreeCascVarPt",&fTreeCascVarPt,"fTreeCascVarPt/F");
463 fTreeCascade->Branch("fTreeCascVarRapXi",&fTreeCascVarRapXi,"fTreeCascVarRapXi/F");
464 fTreeCascade->Branch("fTreeCascVarRapOmega",&fTreeCascVarRapOmega,"fTreeCascVarRapOmega/F");
465 fTreeCascade->Branch("fTreeCascVarNegEta",&fTreeCascVarNegEta,"fTreeCascVarNegEta/F");
466 fTreeCascade->Branch("fTreeCascVarPosEta",&fTreeCascVarPosEta,"fTreeCascVarPosEta/F");
467 fTreeCascade->Branch("fTreeCascVarBachEta",&fTreeCascVarBachEta,"fTreeCascVarBachEta/F");
468 //-----------INFO-FOR-CUTS------------------------
469 fTreeCascade->Branch("fTreeCascVarDCACascDaughters",&fTreeCascVarDCACascDaughters,"fTreeCascVarDCACascDaughters/F");
470 fTreeCascade->Branch("fTreeCascVarDCABachToPrimVtx",&fTreeCascVarDCABachToPrimVtx,"fTreeCascVarDCABachToPrimVtx/F");
471 fTreeCascade->Branch("fTreeCascVarDCAV0Daughters",&fTreeCascVarDCAV0Daughters,"fTreeCascVarDCAV0Daughters/F");
472 fTreeCascade->Branch("fTreeCascVarDCAV0ToPrimVtx",&fTreeCascVarDCAV0ToPrimVtx,"fTreeCascVarDCAV0ToPrimVtx/F");
473 fTreeCascade->Branch("fTreeCascVarDCAPosToPrimVtx",&fTreeCascVarDCAPosToPrimVtx,"fTreeCascVarDCAPosToPrimVtx/F");
474 fTreeCascade->Branch("fTreeCascVarDCANegToPrimVtx",&fTreeCascVarDCANegToPrimVtx,"fTreeCascVarDCANegToPrimVtx/F");
475 fTreeCascade->Branch("fTreeCascVarCascCosPointingAngle",&fTreeCascVarCascCosPointingAngle,"fTreeCascVarCascCosPointingAngle/F");
476 fTreeCascade->Branch("fTreeCascVarCascRadius",&fTreeCascVarCascRadius,"fTreeCascVarCascRadius/F");
477 fTreeCascade->Branch("fTreeCascVarV0Mass",&fTreeCascVarV0Mass,"fTreeCascVarV0Mass/F");
478 fTreeCascade->Branch("fTreeCascVarV0CosPointingAngle",&fTreeCascVarV0CosPointingAngle,"fTreeCascVarV0CosPointingAngle/F");
479 fTreeCascade->Branch("fTreeCascVarV0CosPointingAngleSpecial",&fTreeCascVarV0CosPointingAngleSpecial,"fTreeCascVarV0CosPointingAngleSpecial/F");
480 fTreeCascade->Branch("fTreeCascVarV0Radius",&fTreeCascVarV0Radius,"fTreeCascVarV0Radius/F");
481 fTreeCascade->Branch("fTreeCascVarLeastNbrClusters",&fTreeCascVarLeastNbrClusters,"fTreeCascVarLeastNbrClusters/I");
482 //-----------MULTIPLICITY-INFO--------------------
483 fTreeCascade->Branch("fTreeCascVarCentV0A",&fTreeCascVarCentV0A,"fTreeCascVarCentV0A/F");
484 fTreeCascade->Branch("fTreeCascVarCentV0C",&fTreeCascVarCentV0C,"fTreeCascVarCentV0C/F");
485 fTreeCascade->Branch("fTreeCascVarCentV0M",&fTreeCascVarCentV0M,"fTreeCascVarCentV0M/F");
486 fTreeCascade->Branch("fTreeCascVarCentV0AEq",&fTreeCascVarCentV0AEq,"fTreeCascVarCentV0AEq/F");
487 fTreeCascade->Branch("fTreeCascVarCentV0CEq",&fTreeCascVarCentV0CEq,"fTreeCascVarCentV0CEq/F");
488 fTreeCascade->Branch("fTreeCascVarCentV0MEq",&fTreeCascVarCentV0MEq,"fTreeCascVarCentV0MEq/F");
489 fTreeCascade->Branch("fTreeCascVarAmpV0A",&fTreeCascVarAmpV0A,"fTreeCascVarAmpV0A/F");
490 fTreeCascade->Branch("fTreeCascVarAmpV0C",&fTreeCascVarAmpV0C,"fTreeCascVarAmpV0C/F");
491 fTreeCascade->Branch("fTreeCascVarAmpV0AEq",&fTreeCascVarAmpV0AEq,"fTreeCascVarAmpV0AEq/F");
492 fTreeCascade->Branch("fTreeCascVarAmpV0CEq",&fTreeCascVarAmpV0CEq,"fTreeCascVarAmpV0CEq/F");
493 fTreeCascade->Branch("fTreeCascVarRefMultEta8",&fTreeCascVarRefMultEta8,"fTreeCascVarRefMultEta8/I");
494 fTreeCascade->Branch("fTreeCascVarRefMultEta5",&fTreeCascVarRefMultEta5,"fTreeCascVarRefMultEta5/I");
495 fTreeCascade->Branch("fTreeCascVarRunNumber",&fTreeCascVarRunNumber,"fTreeCascVarRunNumber/I");
496 //-----------DECAY-LENGTH-INFO--------------------
497 fTreeCascade->Branch("fTreeCascVarDistOverTotMom",&fTreeCascVarDistOverTotMom,"fTreeCascVarDistOverTotMom/F");
498 //------------------------------------------------
499 fTreeCascade->Branch("fTreeCascVarNegNSigmaPion",&fTreeCascVarNegNSigmaPion,"fTreeCascVarNegNSigmaPion/F");
500 fTreeCascade->Branch("fTreeCascVarNegNSigmaProton",&fTreeCascVarNegNSigmaProton,"fTreeCascVarNegNSigmaProton/F");
501 fTreeCascade->Branch("fTreeCascVarPosNSigmaPion",&fTreeCascVarPosNSigmaPion,"fTreeCascVarPosNSigmaPion/F");
502 fTreeCascade->Branch("fTreeCascVarPosNSigmaProton",&fTreeCascVarPosNSigmaProton,"fTreeCascVarPosNSigmaProton/F");
503 fTreeCascade->Branch("fTreeCascVarBachNSigmaPion",&fTreeCascVarBachNSigmaPion,"fTreeCascVarBachNSigmaPion/F");
504 fTreeCascade->Branch("fTreeCascVarBachNSigmaKaon",&fTreeCascVarBachNSigmaKaon,"fTreeCascVarBachNSigmaKaon/F");
506 //------------------------------------------------
507 // Particle Identification Setup
508 //------------------------------------------------
510 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
511 AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
512 fPIDResponse = inputHandler->GetPIDResponse();
515 if(! fESDtrackCuts ){
516 fESDtrackCuts = new AliESDtrackCuts();
519 //------------------------------------------------
520 // V0 Multiplicity Histograms
521 //------------------------------------------------
525 fListHist = new TList();
526 fListHist->SetOwner(); // See http://root.cern.ch/root/html/TCollection.html#TCollection:SetOwner
528 if(! fHistEventCounter ) {
529 //Histogram Output: Event-by-Event
530 fHistEventCounter = new TH1D( "fHistEventCounter", ";Evt. Sel. Step;Count",5,0,5);
531 fHistEventCounter->GetXaxis()->SetBinLabel(1, "Processed");
532 fHistEventCounter->GetXaxis()->SetBinLabel(2, "Phys-Sel");
533 fHistEventCounter->GetXaxis()->SetBinLabel(3, "Has Vtx");
534 fHistEventCounter->GetXaxis()->SetBinLabel(4, "Vtx |z|<10cm");
535 fHistEventCounter->GetXaxis()->SetBinLabel(5, "Isn't Pileup");
536 fListHist->Add(fHistEventCounter);
539 //List of Histograms: Normal
540 PostData(1, fListHist);
542 //TTree Object: Saved to base directory. Should cache to disk while saving.
543 //(Important to avoid excessive memory usage, particularly when merging)
544 PostData(2, fTreeEvent);
545 PostData(3, fTreeV0);
546 PostData(4, fTreeCascade);
548 }// end UserCreateOutputObjects
551 //________________________________________________________________________
552 void AliAnalysisTaskStrangenessVsMultiplicity::UserExec(Option_t *)
555 // Called for each event
557 AliESDEvent *lESDevent = 0x0;
559 // Connect to the InputEvent
560 // After these lines, we should have an ESD/AOD event + the number of V0s in it.
562 // Appropriate for ESD analysis!
564 lESDevent = dynamic_cast<AliESDEvent*>( InputEvent() );
566 AliWarning("ERROR: lESDevent not available \n");
570 //Get VZERO Information for multiplicity later
571 AliVVZERO* esdV0 = lESDevent->GetVZEROData();
573 AliError("AliVVZERO not available");
577 fRunNumber = lESDevent->GetRunNumber();
579 Double_t lMagneticField = -10;
580 lMagneticField = lESDevent->GetMagneticField( );
582 //------------------------------------------------
583 // Variable Definition
584 //------------------------------------------------
586 //------------------------------------------------
588 //------------------------------------------------
590 fHistEventCounter->Fill(0.5);
592 UInt_t maskIsSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
593 Bool_t isSelected = 0;
594 isSelected = (maskIsSelected & AliVEvent::kMB) == AliVEvent::kMB;
596 //Standard Min-Bias Selection
597 if ( ! isSelected ) {
598 PostData(1, fListHist);
599 PostData(2, fTreeEvent);
600 PostData(3, fTreeV0);
601 PostData(4, fTreeCascade);
605 fHistEventCounter->Fill(1.5);
607 //------------------------------------------------
608 // Primary Vertex Requirements Section:
609 // ---> pp: has vertex, |z|<10cm
610 //------------------------------------------------
612 //classical Proton-proton like selection
613 const AliESDVertex *lPrimaryBestESDVtx = lESDevent->GetPrimaryVertex();
614 const AliESDVertex *lPrimaryTrackingESDVtx = lESDevent->GetPrimaryVertexTracks();
615 const AliESDVertex *lPrimarySPDVtx = lESDevent->GetPrimaryVertexSPD();
617 Double_t lBestPrimaryVtxPos[3] = {-100.0, -100.0, -100.0};
618 lPrimaryBestESDVtx->GetXYZ( lBestPrimaryVtxPos );
620 //Only accept if Tracking or SPD vertex is fine
621 if (!lPrimarySPDVtx->GetStatus() && !lPrimaryTrackingESDVtx->GetStatus() ){
622 AliWarning("Pb / No SPD prim. vertex nor prim. Tracking vertex ... return !");
623 PostData(1, fListHist);
624 PostData(2, fTreeEvent);
625 PostData(3, fTreeV0);
626 PostData(4, fTreeCascade);
630 //Has SPD or Tracking Vertex
631 fHistEventCounter -> Fill(2.5);
633 //Always do Primary Vertex Selection
634 if(TMath::Abs(lBestPrimaryVtxPos[2]) > 10.0) {
635 AliWarning("Pb / | Z position of Best Prim Vtx | > 10.0 cm ... return !");
636 PostData(1, fListHist);
637 PostData(2, fTreeEvent);
638 PostData(3, fTreeV0);
639 PostData(4, fTreeCascade);
643 //Fill Event selected counter
644 fHistEventCounter -> Fill(3.5);
646 //------------------------------------------------
647 // Check if this isn't pileup
648 //------------------------------------------------
650 if(lESDevent->IsPileupFromSPDInMultBins() ){
651 // minContributors=3, minZdist=0.8, nSigmaZdist=3., nSigmaDiamXY=2., nSigmaDiamZ=5.
652 //-> see http://alisoft.cern.ch/viewvc/trunk/STEER/AliESDEvent.h?root=AliRoot&r1=41914&r2=42199&pathrev=42199
653 AliWarning("Pb / Event tagged as pile-up by SPD... return !");
654 PostData(1, fListHist);
655 PostData(2, fTreeEvent);
656 PostData(3, fTreeV0);
657 PostData(4, fTreeCascade);
660 //Fill Event isn't pileup counter
661 fHistEventCounter -> Fill(4.5);
663 //------------------------------------------------
664 // Multiplicity Information Acquistion
665 //------------------------------------------------
667 //Standard GetReferenceMultiplicity Estimator (0.5 and 0.8)
668 fRefMultEta5 = fESDtrackCuts->GetReferenceMultiplicity(lESDevent, AliESDtrackCuts::kTrackletsITSTPC,0.5);
669 fRefMultEta8 = fESDtrackCuts->GetReferenceMultiplicity(lESDevent, AliESDtrackCuts::kTrackletsITSTPC,0.8);
672 Float_t multV0A = 0; // multiplicity from V0 reco side A
673 Float_t multV0C = 0; // multiplicity from V0 reco side C
674 Float_t multV0AEq = 0; // multiplicity from V0 reco side A
675 Float_t multV0CEq = 0; // multiplicity from V0 reco side C
676 Float_t multV0ACorr = 0; // multiplicity from V0 reco side A
677 Float_t multV0CCorr = 0; // multiplicity from V0 reco side C
679 //Non-Equalized Signal: copy of multV0ACorr and multV0CCorr from AliCentralitySelectionTask
680 //Getters for uncorrected multiplicity
681 multV0A=esdV0->GetMTotV0A();
682 multV0C=esdV0->GetMTotV0C();
684 //Get Z vertex position of SPD vertex (why not Tracking if available?)
685 Float_t zvtx = lPrimarySPDVtx->GetZ();
687 //Acquire Corrected multV0A
688 multV0ACorr = AliESDUtils::GetCorrV0A(multV0A,zvtx);
689 multV0CCorr = AliESDUtils::GetCorrV0C(multV0C,zvtx);
691 //Copy to Event Tree for extra information
692 fAmplitude_V0A = multV0ACorr;
693 fAmplitude_V0C = multV0CCorr;
695 // Equalized signals // From AliCentralitySelectionTask
696 for(Int_t iCh = 4; iCh < 7; ++iCh) {
697 Double_t mult = lESDevent->GetVZEROEqMultiplicity(iCh);
700 for(Int_t iCh = 0; iCh < 3; ++iCh) {
701 Double_t mult = lESDevent->GetVZEROEqMultiplicity(iCh);
704 fAmplitude_V0AEq = multV0AEq;
705 fAmplitude_V0CEq = multV0CEq;
707 fCentrality_V0A = -100;
708 fCentrality_V0C = -100;
709 fCentrality_V0M = -100;
710 fCentrality_V0AEq = -100;
711 fCentrality_V0CEq = -100;
712 fCentrality_V0MEq = -100;
714 //AliCentrality... Check if working?
715 AliCentrality* centrality;
716 centrality = lESDevent->GetCentrality();
717 if ( !(centrality->GetQuality()>1) ){
718 fCentrality_V0A = centrality->GetCentralityPercentile( "V0A" );
719 fCentrality_V0C = centrality->GetCentralityPercentile( "V0C" );
720 fCentrality_V0M = centrality->GetCentralityPercentile( "V0M" );
721 fCentrality_V0AEq = centrality->GetCentralityPercentile( "V0AEq" );
722 fCentrality_V0CEq = centrality->GetCentralityPercentile( "V0CEq" );
723 fCentrality_V0MEq = centrality->GetCentralityPercentile( "V0MEq" );
729 //------------------------------------------------
731 //------------------------------------------------
732 // Fill V0 Tree as needed
733 //------------------------------------------------
735 //Variable definition
736 Int_t lOnFlyStatus = 0;// nv0sOn = 0, nv0sOff = 0;
737 Double_t lChi2V0 = 0;
738 Double_t lDcaV0Daughters = 0, lDcaV0ToPrimVertex = 0;
739 Double_t lDcaPosToPrimVertex = 0, lDcaNegToPrimVertex = 0;
740 Double_t lV0CosineOfPointingAngle = 0;
741 Double_t lV0Radius = 0, lPt = 0;
742 Double_t lRapK0Short = 0, lRapLambda = 0;
743 Double_t lInvMassK0s = 0, lInvMassLambda = 0, lInvMassAntiLambda = 0;
744 Double_t lAlphaV0 = 0, lPtArmV0 = 0;
746 Double_t fMinV0Pt = 0;
747 Double_t fMaxV0Pt = 100;
750 nv0s = lESDevent->GetNumberOfV0s();
752 for (Int_t iV0 = 0; iV0 < nv0s; iV0++) //extra-crazy test
753 {// This is the begining of the V0 loop
754 AliESDv0 *v0 = ((AliESDEvent*)lESDevent)->GetV0(iV0);
757 Double_t tDecayVertexV0[3]; v0->GetXYZ(tDecayVertexV0[0],tDecayVertexV0[1],tDecayVertexV0[2]);
760 v0->GetPxPyPz( tV0mom[0],tV0mom[1],tV0mom[2] );
761 Double_t lV0TotalMomentum = TMath::Sqrt(
762 tV0mom[0]*tV0mom[0]+tV0mom[1]*tV0mom[1]+tV0mom[2]*tV0mom[2] );
764 lV0Radius = TMath::Sqrt(tDecayVertexV0[0]*tDecayVertexV0[0]+tDecayVertexV0[1]*tDecayVertexV0[1]);
767 lRapK0Short = v0->RapK0Short();
768 lRapLambda = v0->RapLambda();
769 if ((lPt<fMinV0Pt)||(fMaxV0Pt<lPt)) continue;
771 UInt_t lKeyPos = (UInt_t)TMath::Abs(v0->GetPindex());
772 UInt_t lKeyNeg = (UInt_t)TMath::Abs(v0->GetNindex());
774 Double_t lMomPos[3]; v0->GetPPxPyPz(lMomPos[0],lMomPos[1],lMomPos[2]);
775 Double_t lMomNeg[3]; v0->GetNPxPyPz(lMomNeg[0],lMomNeg[1],lMomNeg[2]);
777 AliESDtrack *pTrack=((AliESDEvent*)lESDevent)->GetTrack(lKeyPos);
778 AliESDtrack *nTrack=((AliESDEvent*)lESDevent)->GetTrack(lKeyNeg);
779 if (!pTrack || !nTrack) {
780 Printf("ERROR: Could not retreive one of the daughter track");
784 //Daughter Eta for Eta selection, afterwards
785 fTreeVariableNegEta = nTrack->Eta();
786 fTreeVariablePosEta = pTrack->Eta();
788 // Filter like-sign V0 (next: add counter and distribution)
789 if ( pTrack->GetSign() == nTrack->GetSign()){
793 //________________________________________________________________________
794 // Track quality cuts
795 Float_t lPosTrackCrossedRows = pTrack->GetTPCClusterInfo(2,1);
796 Float_t lNegTrackCrossedRows = nTrack->GetTPCClusterInfo(2,1);
797 fTreeVariableLeastNbrCrossedRows = (Int_t) lPosTrackCrossedRows;
798 if( lNegTrackCrossedRows < fTreeVariableLeastNbrCrossedRows )
799 fTreeVariableLeastNbrCrossedRows = (Int_t) lNegTrackCrossedRows;
801 // TPC refit condition (done during reconstruction for Offline but not for On-the-fly)
802 if( !(pTrack->GetStatus() & AliESDtrack::kTPCrefit)) continue;
803 if( !(nTrack->GetStatus() & AliESDtrack::kTPCrefit)) continue;
806 if ( ( ( pTrack->GetTPCClusterInfo(2,1) ) < 70 ) || ( ( nTrack->GetTPCClusterInfo(2,1) ) < 70 ) ) continue;
808 //GetKinkIndex condition
809 if( pTrack->GetKinkIndex(0)>0 || nTrack->GetKinkIndex(0)>0 ) continue;
811 //Findable clusters > 0 condition
812 if( pTrack->GetTPCNclsF()<=0 || nTrack->GetTPCNclsF()<=0 ) continue;
814 //Compute ratio Crossed Rows / Findable clusters
815 //Note: above test avoids division by zero!
816 Float_t lPosTrackCrossedRowsOverFindable = lPosTrackCrossedRows / ((double)(pTrack->GetTPCNclsF()));
817 Float_t lNegTrackCrossedRowsOverFindable = lNegTrackCrossedRows / ((double)(nTrack->GetTPCNclsF()));
819 fTreeVariableLeastRatioCrossedRowsOverFindable = lPosTrackCrossedRowsOverFindable;
820 if( lNegTrackCrossedRowsOverFindable < fTreeVariableLeastRatioCrossedRowsOverFindable )
821 fTreeVariableLeastRatioCrossedRowsOverFindable = lNegTrackCrossedRowsOverFindable;
823 //Lowest Cut Level for Ratio Crossed Rows / Findable = 0.8, set here
824 if ( fTreeVariableLeastRatioCrossedRowsOverFindable < 0.8 ) continue;
826 //End track Quality Cuts
827 //________________________________________________________________________
829 lDcaPosToPrimVertex = TMath::Abs(pTrack->GetD(lBestPrimaryVtxPos[0],
830 lBestPrimaryVtxPos[1],
833 lDcaNegToPrimVertex = TMath::Abs(nTrack->GetD(lBestPrimaryVtxPos[0],
834 lBestPrimaryVtxPos[1],
837 lOnFlyStatus = v0->GetOnFlyStatus();
838 lChi2V0 = v0->GetChi2V0();
839 lDcaV0Daughters = v0->GetDcaV0Daughters();
840 lDcaV0ToPrimVertex = v0->GetD(lBestPrimaryVtxPos[0],lBestPrimaryVtxPos[1],lBestPrimaryVtxPos[2]);
841 lV0CosineOfPointingAngle = v0->GetV0CosineOfPointingAngle(lBestPrimaryVtxPos[0],lBestPrimaryVtxPos[1],lBestPrimaryVtxPos[2]);
842 fTreeVariableV0CosineOfPointingAngle=lV0CosineOfPointingAngle;
844 // Getting invariant mass infos directly from ESD
845 v0->ChangeMassHypothesis(310);
846 lInvMassK0s = v0->GetEffMass();
847 v0->ChangeMassHypothesis(3122);
848 lInvMassLambda = v0->GetEffMass();
849 v0->ChangeMassHypothesis(-3122);
850 lInvMassAntiLambda = v0->GetEffMass();
851 lAlphaV0 = v0->AlphaV0();
852 lPtArmV0 = v0->PtArmV0();
854 fTreeVariablePt = v0->Pt();
855 fTreeVariableChi2V0 = lChi2V0;
856 fTreeVariableDcaV0ToPrimVertex = lDcaV0ToPrimVertex;
857 fTreeVariableDcaV0Daughters = lDcaV0Daughters;
858 fTreeVariableV0CosineOfPointingAngle = lV0CosineOfPointingAngle;
859 fTreeVariableV0Radius = lV0Radius;
860 fTreeVariableDcaPosToPrimVertex = lDcaPosToPrimVertex;
861 fTreeVariableDcaNegToPrimVertex = lDcaNegToPrimVertex;
862 fTreeVariableInvMassK0s = lInvMassK0s;
863 fTreeVariableInvMassLambda = lInvMassLambda;
864 fTreeVariableInvMassAntiLambda = lInvMassAntiLambda;
865 fTreeVariableRapK0Short = lRapK0Short;
866 fTreeVariableRapLambda = lRapLambda;
867 fTreeVariableAlphaV0 = lAlphaV0;
868 fTreeVariablePtArmV0 = lPtArmV0;
870 //Official means of acquiring N-sigmas
871 fTreeVariableNSigmasPosProton = fPIDResponse->NumberOfSigmasTPC( pTrack, AliPID::kProton );
872 fTreeVariableNSigmasPosPion = fPIDResponse->NumberOfSigmasTPC( pTrack, AliPID::kPion );
873 fTreeVariableNSigmasNegProton = fPIDResponse->NumberOfSigmasTPC( nTrack, AliPID::kProton );
874 fTreeVariableNSigmasNegPion = fPIDResponse->NumberOfSigmasTPC( nTrack, AliPID::kPion );
876 //This requires an Invariant Mass Hypothesis afterwards
877 fTreeVariableDistOverTotMom = TMath::Sqrt(
878 TMath::Power( tDecayVertexV0[0] - lBestPrimaryVtxPos[0] , 2) +
879 TMath::Power( tDecayVertexV0[1] - lBestPrimaryVtxPos[1] , 2) +
880 TMath::Power( tDecayVertexV0[2] - lBestPrimaryVtxPos[2] , 2)
882 fTreeVariableDistOverTotMom /= (lV0TotalMomentum+1e-10); //avoid division by zero, to be sure
884 //Copy Multiplicity information
885 fTreeVariableCentV0A = fCentrality_V0A;
886 fTreeVariableCentV0C = fCentrality_V0C;
887 fTreeVariableCentV0M = fCentrality_V0M;
888 fTreeVariableCentV0AEq = fCentrality_V0AEq;
889 fTreeVariableCentV0CEq = fCentrality_V0CEq;
890 fTreeVariableCentV0MEq = fCentrality_V0MEq;
891 fTreeVariableAmpV0A = fAmplitude_V0A;
892 fTreeVariableAmpV0C = fAmplitude_V0C;
893 fTreeVariableAmpV0AEq = fAmplitude_V0AEq;
894 fTreeVariableAmpV0CEq = fAmplitude_V0CEq;
895 fTreeVariableRefMultEta8 = fRefMultEta8;
896 fTreeVariableRefMultEta5 = fRefMultEta5;
897 fTreeVariableRunNumber = fRunNumber;
899 //------------------------------------------------
901 //------------------------------------------------
903 // The conditionals are meant to decrease excessive
906 //First Selection: Reject OnFly
907 if( lOnFlyStatus == 0 ){
908 //Second Selection: rough 20-sigma band, parametric.
909 //K0Short: Enough to parametrize peak broadening with linear function.
910 Double_t lUpperLimitK0Short = (5.63707e-01) + (1.14979e-02)*fTreeVariablePt;
911 Double_t lLowerLimitK0Short = (4.30006e-01) - (1.10029e-02)*fTreeVariablePt;
912 //Lambda: Linear (for higher pt) plus exponential (for low-pt broadening)
913 //[0]+[1]*x+[2]*TMath::Exp(-[3]*x)
914 Double_t lUpperLimitLambda = (1.13688e+00) + (5.27838e-03)*fTreeVariablePt + (8.42220e-02)*TMath::Exp(-(3.80595e+00)*fTreeVariablePt);
915 Double_t lLowerLimitLambda = (1.09501e+00) - (5.23272e-03)*fTreeVariablePt - (7.52690e-02)*TMath::Exp(-(3.46339e+00)*fTreeVariablePt);
917 if( (fTreeVariableInvMassLambda < lUpperLimitLambda && fTreeVariableInvMassLambda > lLowerLimitLambda ) ||
918 (fTreeVariableInvMassAntiLambda < lUpperLimitLambda && fTreeVariableInvMassAntiLambda > lLowerLimitLambda ) ||
919 (fTreeVariableInvMassK0s < lUpperLimitK0Short && fTreeVariableInvMassK0s > lLowerLimitK0Short ) ){
920 //Pre-selection in case this is AA...
921 if ( TMath::Abs(fTreeVariableNegEta)<0.8 && TMath::Abs(fTreeVariablePosEta)<0.8 && fkSaveV0Tree ) fTreeV0->Fill();
924 }// This is the end of the V0 loop
926 //------------------------------------------------
927 // Fill V0 tree over.
928 //------------------------------------------------
932 //------------------------------------------------
933 // Rerun cascade vertexer!
934 //------------------------------------------------
936 if( fkRunVertexers ){
937 lESDevent->ResetCascades();
938 lESDevent->ResetV0s();
940 AliV0vertexer lV0vtxer;
941 AliCascadeVertexer lCascVtxer;
943 lV0vtxer.SetDefaultCuts(fV0VertexerSels);
944 lCascVtxer.SetDefaultCuts(fCascadeVertexerSels);
946 lV0vtxer.Tracks2V0vertices(lESDevent);
947 lCascVtxer.V0sTracks2CascadeVertices(lESDevent);
950 //------------------------------------------------
951 // MAIN CASCADE LOOP STARTS HERE
952 //------------------------------------------------
953 // Code Credit: Antonin Maire (thanks^100)
954 // ---> This is an adaptation
956 Long_t ncascades = 0;
957 ncascades = lESDevent->GetNumberOfCascades();
959 for (Int_t iXi = 0; iXi < ncascades; iXi++){
960 //------------------------------------------------
962 //------------------------------------------------
963 //Double_t lTrkgPrimaryVtxRadius3D = -500.0;
964 //Double_t lBestPrimaryVtxRadius3D = -500.0;
966 // - 1st part of initialisation : variables needed to store AliESDCascade data members
967 Double_t lEffMassXi = 0. ;
968 //Double_t lChi2Xi = -1. ;
969 Double_t lDcaXiDaughters = -1. ;
970 Double_t lXiCosineOfPointingAngle = -1. ;
971 Double_t lPosXi[3] = { -1000.0, -1000.0, -1000.0 };
972 Double_t lXiRadius = -1000. ;
974 // - 2nd part of initialisation : Nbr of clusters within TPC for the 3 daughter cascade tracks
975 Int_t lPosTPCClusters = -1; // For ESD only ...//FIXME : wait for availability in AOD
976 Int_t lNegTPCClusters = -1; // For ESD only ...
977 Int_t lBachTPCClusters = -1; // For ESD only ...
979 // - 3rd part of initialisation : about V0 part in cascades
980 Double_t lInvMassLambdaAsCascDghter = 0.;
981 //Double_t lV0Chi2Xi = -1. ;
982 Double_t lDcaV0DaughtersXi = -1.;
984 Double_t lDcaBachToPrimVertexXi = -1., lDcaV0ToPrimVertexXi = -1.;
985 Double_t lDcaPosToPrimVertexXi = -1.;
986 Double_t lDcaNegToPrimVertexXi = -1.;
987 Double_t lV0CosineOfPointingAngleXi = -1. ;
988 Double_t lV0CosineOfPointingAngleXiSpecial = -1. ;
989 Double_t lPosV0Xi[3] = { -1000. , -1000., -1000. }; // Position of VO coming from cascade
990 Double_t lV0RadiusXi = -1000.0;
991 Double_t lV0quality = 0.;
993 // - 4th part of initialisation : Effective masses
994 Double_t lInvMassXiMinus = 0.;
995 Double_t lInvMassXiPlus = 0.;
996 Double_t lInvMassOmegaMinus = 0.;
997 Double_t lInvMassOmegaPlus = 0.;
999 // - 6th part of initialisation : extra info for QA
1000 Double_t lXiMomX = 0. , lXiMomY = 0., lXiMomZ = 0.;
1001 Double_t lXiTransvMom = 0. ;
1002 //Double_t lXiTransvMomMC= 0. ;
1003 Double_t lXiTotMom = 0. ;
1005 Double_t lBachMomX = 0., lBachMomY = 0., lBachMomZ = 0.;
1006 //Double_t lBachTransvMom = 0.;
1007 //Double_t lBachTotMom = 0.;
1009 fTreeCascVarNegNSigmaPion = -100;
1010 fTreeCascVarNegNSigmaProton = -100;
1011 fTreeCascVarPosNSigmaPion = -100;
1012 fTreeCascVarPosNSigmaProton = -100;
1013 fTreeCascVarBachNSigmaPion = -100;
1014 fTreeCascVarBachNSigmaKaon = -100;
1016 Short_t lChargeXi = -2;
1017 //Double_t lV0toXiCosineOfPointingAngle = 0. ;
1019 Double_t lRapXi = -20.0, lRapOmega = -20.0; // lEta = -20.0, lTheta = 360., lPhi = 720. ;
1020 //Double_t lAlphaXi = -200., lPtArmXi = -200.0;
1022 // -------------------------------------
1023 // II.ESD - Calculation Part dedicated to Xi vertices (ESD)
1025 AliESDcascade *xi = lESDevent->GetCascade(iXi);
1028 // - II.Step 2 : Assigning the necessary variables for specific AliESDcascade data members (ESD)
1031 xi->ChangeMassHypothesis(lV0quality , 3312); // default working hypothesis : cascade = Xi- decay
1033 lEffMassXi = xi->GetEffMassXi();
1034 //lChi2Xi = xi->GetChi2Xi();
1035 lDcaXiDaughters = xi->GetDcaXiDaughters();
1036 lXiCosineOfPointingAngle = xi->GetCascadeCosineOfPointingAngle( lBestPrimaryVtxPos[0],
1037 lBestPrimaryVtxPos[1],
1038 lBestPrimaryVtxPos[2] );
1039 // Take care : the best available vertex should be used (like in AliCascadeVertexer)
1041 xi->GetXYZcascade( lPosXi[0], lPosXi[1], lPosXi[2] );
1042 lXiRadius = TMath::Sqrt( lPosXi[0]*lPosXi[0] + lPosXi[1]*lPosXi[1] );
1044 // - II.Step 3 : around the tracks : Bach + V0 (ESD)
1045 // ~ Necessary variables for ESDcascade data members coming from the ESDv0 part (inheritance)
1048 UInt_t lIdxPosXi = (UInt_t) TMath::Abs( xi->GetPindex() );
1049 UInt_t lIdxNegXi = (UInt_t) TMath::Abs( xi->GetNindex() );
1050 UInt_t lBachIdx = (UInt_t) TMath::Abs( xi->GetBindex() );
1051 // Care track label can be negative in MC production (linked with the track quality)
1052 // However = normally, not the case for track index ...
1054 // FIXME : rejection of a double use of a daughter track (nothing but just a crosscheck of what is done in the cascade vertexer)
1055 if(lBachIdx == lIdxNegXi) {
1056 AliWarning("Pb / Idx(Bach. track) = Idx(Neg. track) ... continue!"); continue;
1058 if(lBachIdx == lIdxPosXi) {
1059 AliWarning("Pb / Idx(Bach. track) = Idx(Pos. track) ... continue!"); continue;
1062 AliESDtrack *pTrackXi = lESDevent->GetTrack( lIdxPosXi );
1063 AliESDtrack *nTrackXi = lESDevent->GetTrack( lIdxNegXi );
1064 AliESDtrack *bachTrackXi = lESDevent->GetTrack( lBachIdx );
1066 if (!pTrackXi || !nTrackXi || !bachTrackXi ) {
1067 AliWarning("ERROR: Could not retrieve one of the 3 ESD daughter tracks of the cascade ...");
1071 fTreeCascVarPosEta = pTrackXi->Eta();
1072 fTreeCascVarNegEta = nTrackXi->Eta();
1073 fTreeCascVarBachEta = bachTrackXi->Eta();
1075 Double_t lBMom[3], lNMom[3], lPMom[3];
1076 xi->GetBPxPyPz( lBMom[0], lBMom[1], lBMom[2] );
1077 xi->GetPPxPyPz( lPMom[0], lPMom[1], lPMom[2] );
1078 xi->GetNPxPyPz( lNMom[0], lNMom[1], lNMom[2] );
1080 //fTreeCascVarBachTransMom = TMath::Sqrt( lBMom[0]*lBMom[0] + lBMom[1]*lBMom[1] );
1081 //fTreeCascVarPosTransMom = TMath::Sqrt( lPMom[0]*lPMom[0] + lPMom[1]*lPMom[1] );
1082 //fTreeCascVarNegTransMom = TMath::Sqrt( lNMom[0]*lNMom[0] + lNMom[1]*lNMom[1] );
1084 //------------------------------------------------
1085 // TPC dEdx information
1086 //------------------------------------------------
1087 fTreeCascVarNegNSigmaPion = fPIDResponse->NumberOfSigmasTPC( nTrackXi, AliPID::kPion );
1088 fTreeCascVarNegNSigmaProton = fPIDResponse->NumberOfSigmasTPC( nTrackXi, AliPID::kProton );
1089 fTreeCascVarPosNSigmaPion = fPIDResponse->NumberOfSigmasTPC( pTrackXi, AliPID::kPion );
1090 fTreeCascVarPosNSigmaProton = fPIDResponse->NumberOfSigmasTPC( pTrackXi, AliPID::kProton );
1091 fTreeCascVarBachNSigmaPion = fPIDResponse->NumberOfSigmasTPC( bachTrackXi, AliPID::kPion );
1092 fTreeCascVarBachNSigmaKaon = fPIDResponse->NumberOfSigmasTPC( bachTrackXi, AliPID::kKaon );
1094 //------------------------------------------------
1095 // TPC Number of clusters info
1096 // --- modified to save the smallest number
1097 // --- of TPC clusters for the 3 tracks
1098 //------------------------------------------------
1100 lPosTPCClusters = pTrackXi->GetTPCNcls();
1101 lNegTPCClusters = nTrackXi->GetTPCNcls();
1102 lBachTPCClusters = bachTrackXi->GetTPCNcls();
1104 // 1 - Poor quality related to TPCrefit
1105 ULong_t pStatus = pTrackXi->GetStatus();
1106 ULong_t nStatus = nTrackXi->GetStatus();
1107 ULong_t bachStatus = bachTrackXi->GetStatus();
1109 //fTreeCascVarkITSRefitBachelor = kTRUE;
1110 //fTreeCascVarkITSRefitNegative = kTRUE;
1111 //fTreeCascVarkITSRefitPositive = kTRUE;
1113 if ((pStatus&AliESDtrack::kTPCrefit) == 0) { AliWarning("Pb / V0 Pos. track has no TPCrefit ... continue!"); continue; }
1114 if ((nStatus&AliESDtrack::kTPCrefit) == 0) { AliWarning("Pb / V0 Neg. track has no TPCrefit ... continue!"); continue; }
1115 if ((bachStatus&AliESDtrack::kTPCrefit) == 0) { AliWarning("Pb / Bach. track has no TPCrefit ... continue!"); continue; }
1117 // 2 - Poor quality related to TPC clusters: lowest cut of 70 clusters
1118 if(lPosTPCClusters < 70) { AliWarning("Pb / V0 Pos. track has less than 70 TPC clusters ... continue!"); continue; }
1119 if(lNegTPCClusters < 70) { AliWarning("Pb / V0 Neg. track has less than 70 TPC clusters ... continue!"); continue; }
1120 if(lBachTPCClusters < 70) { AliWarning("Pb / Bach. track has less than 70 TPC clusters ... continue!"); continue; }
1121 Int_t leastnumberofclusters = 1000;
1122 if( lPosTPCClusters < leastnumberofclusters ) leastnumberofclusters = lPosTPCClusters;
1123 if( lNegTPCClusters < leastnumberofclusters ) leastnumberofclusters = lNegTPCClusters;
1124 if( lBachTPCClusters < leastnumberofclusters ) leastnumberofclusters = lBachTPCClusters;
1126 lInvMassLambdaAsCascDghter = xi->GetEffMass();
1127 // This value shouldn't change, whatever the working hyp. is : Xi-, Xi+, Omega-, Omega+
1128 lDcaV0DaughtersXi = xi->GetDcaV0Daughters();
1129 //lV0Chi2Xi = xi->GetChi2V0();
1131 lV0CosineOfPointingAngleXi = xi->GetV0CosineOfPointingAngle( lBestPrimaryVtxPos[0],
1132 lBestPrimaryVtxPos[1],
1133 lBestPrimaryVtxPos[2] );
1134 //Modification: V0 CosPA wrt to Cascade decay vertex
1135 lV0CosineOfPointingAngleXiSpecial = xi->GetV0CosineOfPointingAngle( lPosXi[0],
1139 lDcaV0ToPrimVertexXi = xi->GetD( lBestPrimaryVtxPos[0],
1140 lBestPrimaryVtxPos[1],
1141 lBestPrimaryVtxPos[2] );
1143 lDcaBachToPrimVertexXi = TMath::Abs( bachTrackXi->GetD( lBestPrimaryVtxPos[0],
1144 lBestPrimaryVtxPos[1],
1146 // Note : AliExternalTrackParam::GetD returns an algebraic value ...
1148 xi->GetXYZ( lPosV0Xi[0], lPosV0Xi[1], lPosV0Xi[2] );
1149 lV0RadiusXi = TMath::Sqrt( lPosV0Xi[0]*lPosV0Xi[0] + lPosV0Xi[1]*lPosV0Xi[1] );
1151 lDcaPosToPrimVertexXi = TMath::Abs( pTrackXi ->GetD( lBestPrimaryVtxPos[0],
1152 lBestPrimaryVtxPos[1],
1155 lDcaNegToPrimVertexXi = TMath::Abs( nTrackXi ->GetD( lBestPrimaryVtxPos[0],
1156 lBestPrimaryVtxPos[1],
1159 // - II.Step 4 : around effective masses (ESD)
1160 // ~ change mass hypotheses to cover all the possibilities : Xi-/+, Omega -/+
1162 if( bachTrackXi->Charge() < 0 ) {
1164 xi->ChangeMassHypothesis(lV0quality , 3312);
1165 // Calculate the effective mass of the Xi- candidate.
1166 // pdg code 3312 = Xi-
1167 lInvMassXiMinus = xi->GetEffMassXi();
1170 xi->ChangeMassHypothesis(lV0quality , 3334);
1171 // Calculate the effective mass of the Xi- candidate.
1172 // pdg code 3334 = Omega-
1173 lInvMassOmegaMinus = xi->GetEffMassXi();
1176 xi->ChangeMassHypothesis(lV0quality , 3312); // Back to default hyp.
1177 }// end if negative bachelor
1180 if( bachTrackXi->Charge() > 0 ){
1182 xi->ChangeMassHypothesis(lV0quality , -3312);
1183 // Calculate the effective mass of the Xi+ candidate.
1184 // pdg code -3312 = Xi+
1185 lInvMassXiPlus = xi->GetEffMassXi();
1188 xi->ChangeMassHypothesis(lV0quality , -3334);
1189 // Calculate the effective mass of the Xi+ candidate.
1190 // pdg code -3334 = Omega+
1191 lInvMassOmegaPlus = xi->GetEffMassXi();
1194 xi->ChangeMassHypothesis(lV0quality , -3312); // Back to "default" hyp.
1195 }// end if positive bachelor
1196 // - II.Step 6 : extra info for QA (ESD)
1197 // miscellaneous pieces of info that may help regarding data quality assessment.
1200 xi->GetPxPyPz( lXiMomX, lXiMomY, lXiMomZ );
1201 lXiTransvMom = TMath::Sqrt( lXiMomX*lXiMomX + lXiMomY*lXiMomY );
1202 lXiTotMom = TMath::Sqrt( lXiMomX*lXiMomX + lXiMomY*lXiMomY + lXiMomZ*lXiMomZ );
1204 xi->GetBPxPyPz( lBachMomX, lBachMomY, lBachMomZ );
1205 //lBachTransvMom = TMath::Sqrt( lBachMomX*lBachMomX + lBachMomY*lBachMomY );
1206 //lBachTotMom = TMath::Sqrt( lBachMomX*lBachMomX + lBachMomY*lBachMomY + lBachMomZ*lBachMomZ );
1208 lChargeXi = xi->Charge();
1210 //lV0toXiCosineOfPointingAngle = xi->GetV0CosineOfPointingAngle( lPosXi[0], lPosXi[1], lPosXi[2] );
1212 lRapXi = xi->RapXi();
1213 lRapOmega = xi->RapOmega();
1215 //lTheta = xi->Theta() *180.0/TMath::Pi();
1216 //lPhi = xi->Phi() *180.0/TMath::Pi();
1217 //lAlphaXi = xi->AlphaXi();
1218 //lPtArmXi = xi->PtArmXi();
1220 //------------------------------------------------
1221 // Set Variables for adding to tree
1222 //------------------------------------------------
1224 fTreeCascVarCharge = lChargeXi;
1225 if(lInvMassXiMinus!=0) fTreeCascVarMassAsXi = lInvMassXiMinus;
1226 if(lInvMassXiPlus!=0) fTreeCascVarMassAsXi = lInvMassXiPlus;
1227 if(lInvMassOmegaMinus!=0) fTreeCascVarMassAsOmega = lInvMassOmegaMinus;
1228 if(lInvMassOmegaPlus!=0) fTreeCascVarMassAsOmega = lInvMassOmegaPlus;
1229 fTreeCascVarPt = lXiTransvMom;
1230 fTreeCascVarRapXi = lRapXi ;
1231 fTreeCascVarRapOmega = lRapOmega ;
1232 fTreeCascVarDCACascDaughters = lDcaXiDaughters;
1233 fTreeCascVarDCABachToPrimVtx = lDcaBachToPrimVertexXi;
1234 fTreeCascVarDCAV0Daughters = lDcaV0DaughtersXi;
1235 fTreeCascVarDCAV0ToPrimVtx = lDcaV0ToPrimVertexXi;
1236 fTreeCascVarDCAPosToPrimVtx = lDcaPosToPrimVertexXi;
1237 fTreeCascVarDCANegToPrimVtx = lDcaNegToPrimVertexXi;
1238 fTreeCascVarCascCosPointingAngle = lXiCosineOfPointingAngle;
1239 fTreeCascVarCascRadius = lXiRadius;
1240 fTreeCascVarV0Mass = lInvMassLambdaAsCascDghter;
1241 fTreeCascVarV0CosPointingAngle = lV0CosineOfPointingAngleXi;
1242 fTreeCascVarV0CosPointingAngleSpecial = lV0CosineOfPointingAngleXiSpecial;
1243 fTreeCascVarV0Radius = lV0RadiusXi;
1244 fTreeCascVarLeastNbrClusters = leastnumberofclusters;
1246 //Copy Multiplicity information
1247 fTreeCascVarCentV0A = fCentrality_V0A;
1248 fTreeCascVarCentV0C = fCentrality_V0C;
1249 fTreeCascVarCentV0M = fCentrality_V0M;
1250 fTreeCascVarCentV0AEq = fCentrality_V0AEq;
1251 fTreeCascVarCentV0CEq = fCentrality_V0CEq;
1252 fTreeCascVarCentV0MEq = fCentrality_V0MEq;
1253 fTreeCascVarAmpV0A = fAmplitude_V0A;
1254 fTreeCascVarAmpV0C = fAmplitude_V0C;
1255 fTreeCascVarAmpV0AEq = fAmplitude_V0AEq;
1256 fTreeCascVarAmpV0CEq = fAmplitude_V0CEq;
1257 fTreeCascVarRefMultEta8 = fRefMultEta8;
1258 fTreeCascVarRefMultEta5 = fRefMultEta5;
1259 fTreeCascVarRunNumber = fRunNumber;
1261 fTreeCascVarDistOverTotMom = TMath::Sqrt(
1262 TMath::Power( lPosXi[0] - lBestPrimaryVtxPos[0] , 2) +
1263 TMath::Power( lPosXi[1] - lBestPrimaryVtxPos[1] , 2) +
1264 TMath::Power( lPosXi[2] - lBestPrimaryVtxPos[2] , 2)
1266 fTreeCascVarDistOverTotMom /= (lXiTotMom+1e-13);
1268 //All vars not specified here: specified elsewhere!
1270 //------------------------------------------------
1272 //------------------------------------------------
1274 // The conditional is meant to decrease excessive
1275 // memory usage! Be careful when loosening the
1278 //Xi Mass window: 150MeV wide
1279 //Omega mass window: 150MeV wide
1281 if( fkSaveCascadeTree && ( (fTreeCascVarMassAsXi<1.32+0.075&&fTreeCascVarMassAsXi>1.32-0.075) ||
1282 (fTreeCascVarMassAsOmega<1.68+0.075&&fTreeCascVarMassAsOmega>1.68-0.075) ) ){
1283 fTreeCascade->Fill();
1286 //------------------------------------------------
1288 //------------------------------------------------
1290 }// end of the Cascade loop (ESD or AOD)
1292 // Post output data.
1293 PostData(1, fListHist);
1294 PostData(2, fTreeEvent);
1295 PostData(3, fTreeV0);
1296 PostData(4, fTreeCascade);
1299 //________________________________________________________________________
1300 void AliAnalysisTaskStrangenessVsMultiplicity::Terminate(Option_t *)
1302 // Draw result to the screen
1303 // Called once at the end of the query
1305 TList *cRetrievedList = 0x0;
1306 cRetrievedList = (TList*)GetOutputData(1);
1307 if(!cRetrievedList){
1308 Printf("ERROR - AliAnalysisTaskStrangenessVsMultiplicity : ouput data container list not available\n");
1312 fHistEventCounter = dynamic_cast<TH1D*> ( cRetrievedList->FindObject("fHistEventCounter") );
1313 if (!fHistEventCounter) {
1314 Printf("ERROR - AliAnalysisTaskStrangenessVsMultiplicity : fHistEventCounter not available");
1318 TCanvas *canCheck = new TCanvas("AliAnalysisTaskStrangenessVsMultiplicity","V0 Multiplicity",10,10,510,510);
1319 canCheck->cd(1)->SetLogy();
1321 fHistEventCounter->SetMarkerStyle(22);
1322 fHistEventCounter->DrawCopy("E");
1325 //----------------------------------------------------------------------------
1327 Double_t AliAnalysisTaskStrangenessVsMultiplicity::MyRapidity(Double_t rE, Double_t rPz) const
1329 // Local calculation for rapidity
1330 Double_t ReturnValue = -100;
1331 if( (rE-rPz+1.e-13) != 0 && (rE+rPz) != 0 ){
1332 ReturnValue = 0.5*TMath::Log((rE+rPz)/(rE-rPz+1.e-13));