First commit of experimental task to test VZERO percentiles from AliCentrality for...
[u/mrichter/AliRoot.git] / PWGLF / STRANGENESS / Cascades / AliAnalysisTaskStrangenessVsMultiplicity.cxx
CommitLineData
af0d4791 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// 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:
21//
22// 1) (under development) it should act as an auxiliary task and provide a
23// calibrated estimator
24//
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:
28//
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)
37// --- Run Number
38//
39// +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
40
41class TTree;
42class TParticle;
43class TVector3;
44
45//class AliMCEventHandler;
46//class AliMCEvent;
47//class AliStack;
48
49class AliESDVertex;
50class AliAODVertex;
51class AliESDv0;
52class AliAODv0;
53
54#include <Riostream.h>
55#include "TList.h"
56#include "TH1.h"
57#include "TH2.h"
58#include "TH3.h"
59#include "TFile.h"
60#include "THnSparse.h"
61#include "TVector3.h"
62#include "TCanvas.h"
63#include "TMath.h"
64#include "TLegend.h"
65//#include "AliLog.h"
66
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"
78#include "AliStack.h"
79#include "AliCentrality.h"
80
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"
91
92using std::cout;
93using std::endl;
94
95ClassImp(AliAnalysisTaskStrangenessVsMultiplicity)
96
97AliAnalysisTaskStrangenessVsMultiplicity::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
103 fAmplitude_V0A (0),
104 fAmplitude_V0C (0),
105 fAmplitude_V0AEq (0),
106 fAmplitude_V0CEq (0),
107 fCentrality_V0A(0),
108 fCentrality_V0C(0),
109 fCentrality_V0M(0),
110 fCentrality_V0AEq(0),
111 fCentrality_V0CEq(0),
112 fCentrality_V0MEq(0),
113 fRefMultEta5(0),
114 fRefMultEta8(0),
115 fRunNumber(0),
116 //---> Variables for fTreeCascade
117 fTreeCascVarCharge(0),
118 fTreeCascVarMassAsXi(0),
119 fTreeCascVarMassAsOmega(0),
120 fTreeCascVarPt(0),
121 fTreeCascVarRapXi(0),
122 fTreeCascVarRapOmega(0),
123 fTreeCascVarNegEta(0),
124 fTreeCascVarPosEta(0),
125 fTreeCascVarBachEta(0),
126 fTreeCascVarDCACascDaughters(0),
127 fTreeCascVarDCABachToPrimVtx(0),
128 fTreeCascVarDCAV0Daughters(0),
129 fTreeCascVarDCAV0ToPrimVtx(0),
130 fTreeCascVarDCAPosToPrimVtx(0),
131 fTreeCascVarDCANegToPrimVtx(0),
132 fTreeCascVarCascCosPointingAngle(0),
133 fTreeCascVarCascRadius(0),
134 fTreeCascVarV0Mass(0),
135 fTreeCascVarV0CosPointingAngle(0),
136 fTreeCascVarV0CosPointingAngleSpecial(0),
137 fTreeCascVarV0Radius(0),
138 fTreeCascVarLeastNbrClusters(0),
139 fTreeCascVarDistOverTotMom(0),
140 fTreeCascVarNegNSigmaPion(0),
141 fTreeCascVarNegNSigmaProton(0),
142 fTreeCascVarPosNSigmaPion(0),
143 fTreeCascVarPosNSigmaProton(0),
144 fTreeCascVarBachNSigmaPion(0),
145 fTreeCascVarBachNSigmaKaon(0),
146 fTreeCascVarCentV0A(0),
147 fTreeCascVarCentV0C(0),
148 fTreeCascVarCentV0M(0),
149 fTreeCascVarCentV0AEq(0),
150 fTreeCascVarCentV0CEq(0),
151 fTreeCascVarCentV0MEq(0),
152 fTreeCascVarAmpV0A(0),
153 fTreeCascVarAmpV0C(0),
154 fTreeCascVarAmpV0AEq(0),
155 fTreeCascVarAmpV0CEq(0),
156 fTreeCascVarRefMultEta8(0),
157 fTreeCascVarRefMultEta5(0),
158 //Histos
159 fHistEventCounter(0)
160//------------------------------------------------
161// Tree Variables
162{
163
164}
165
166AliAnalysisTaskStrangenessVsMultiplicity::AliAnalysisTaskStrangenessVsMultiplicity(const char *name)
167 : AliAnalysisTaskSE(name), fListHist(0), fTreeEvent(0), fTreeV0(0), fTreeCascade(0), fPIDResponse(0), fESDtrackCuts(0),
168 fkSaveV0Tree ( kFALSE ),
169 fkSaveCascadeTree ( kTRUE ),
170 fkRunVertexers ( kTRUE ),
171 //---> Variables for fTreeEvent
172 fAmplitude_V0A (0),
173 fAmplitude_V0C (0),
174 fAmplitude_V0AEq (0),
175 fAmplitude_V0CEq (0),
176 fCentrality_V0A(0),
177 fCentrality_V0C(0),
178 fCentrality_V0M(0),
179 fCentrality_V0AEq(0),
180 fCentrality_V0CEq(0),
181 fCentrality_V0MEq(0),
182 fRefMultEta5(0),
183 fRefMultEta8(0),
184 fRunNumber(0),
185 //---> Variables for fTreeCascade
186 fTreeCascVarCharge(0),
187 fTreeCascVarMassAsXi(0),
188 fTreeCascVarMassAsOmega(0),
189 fTreeCascVarPt(0),
190 fTreeCascVarRapXi(0),
191 fTreeCascVarRapOmega(0),
192 fTreeCascVarNegEta(0),
193 fTreeCascVarPosEta(0),
194 fTreeCascVarBachEta(0),
195 fTreeCascVarDCACascDaughters(0),
196 fTreeCascVarDCABachToPrimVtx(0),
197 fTreeCascVarDCAV0Daughters(0),
198 fTreeCascVarDCAV0ToPrimVtx(0),
199 fTreeCascVarDCAPosToPrimVtx(0),
200 fTreeCascVarDCANegToPrimVtx(0),
201 fTreeCascVarCascCosPointingAngle(0),
202 fTreeCascVarCascRadius(0),
203 fTreeCascVarV0Mass(0),
204 fTreeCascVarV0CosPointingAngle(0),
205 fTreeCascVarV0CosPointingAngleSpecial(0),
206 fTreeCascVarV0Radius(0),
207 fTreeCascVarLeastNbrClusters(0),
208 fTreeCascVarDistOverTotMom(0),
209 fTreeCascVarNegNSigmaPion(0),
210 fTreeCascVarNegNSigmaProton(0),
211 fTreeCascVarPosNSigmaPion(0),
212 fTreeCascVarPosNSigmaProton(0),
213 fTreeCascVarBachNSigmaPion(0),
214 fTreeCascVarBachNSigmaKaon(0),
215 fTreeCascVarCentV0A(0),
216 fTreeCascVarCentV0C(0),
217 fTreeCascVarCentV0M(0),
218 fTreeCascVarCentV0AEq(0),
219 fTreeCascVarCentV0CEq(0),
220 fTreeCascVarCentV0MEq(0),
221 fTreeCascVarAmpV0A(0),
222 fTreeCascVarAmpV0C(0),
223 fTreeCascVarAmpV0AEq(0),
224 fTreeCascVarAmpV0CEq(0),
225 fTreeCascVarRefMultEta8(0),
226 fTreeCascVarRefMultEta5(0),
227 //Histos
228 fHistEventCounter(0)
229{
230
231 //Re-vertex: Will only apply for cascade candidates
232
233 fV0VertexerSels[0] = 33. ; // max allowed chi2
234 fV0VertexerSels[1] = 0.02; // min allowed impact parameter for the 1st daughter (LHC09a4 : 0.05)
235 fV0VertexerSels[2] = 0.02; // min allowed impact parameter for the 2nd daughter (LHC09a4 : 0.05)
236 fV0VertexerSels[4] = 0.95; // min allowed cosine of V0's pointing angle (LHC09a4 : 0.99)
237 fV0VertexerSels[5] = 1.0 ; // min radius of the fiducial volume (LHC09a4 : 0.2)
238 fV0VertexerSels[6] = 200. ; // max radius of the fiducial volume (LHC09a4 : 100.0)
239
240 fCascadeVertexerSels[0] = 33. ; // max allowed chi2 (same as PDC07)
241 fCascadeVertexerSels[1] = 0.05 ; // min allowed V0 impact parameter (PDC07 : 0.05 / LHC09a4 : 0.025 )
242 fCascadeVertexerSels[2] = 0.010; // "window" around the Lambda mass (PDC07 : 0.008 / LHC09a4 : 0.010 )
243 fCascadeVertexerSels[3] = 0.03 ; // min allowed bachelor's impact parameter (PDC07 : 0.035 / LHC09a4 : 0.025 )
244 fCascadeVertexerSels[4] = 2.0 ; // max allowed DCA between the V0 and the bachelor (PDC07 : 0.1 / LHC09a4 : 0.2 )
245 fCascadeVertexerSels[5] = 0.95 ; // min allowed cosine of the cascade pointing angle (PDC07 : 0.9985 / LHC09a4 : 0.998 )
246 fCascadeVertexerSels[6] = 0.4 ; // min radius of the fiducial volume (PDC07 : 0.9 / LHC09a4 : 0.2 )
247 fCascadeVertexerSels[7] = 100. ; // max radius of the fiducial volume (PDC07 : 100 / LHC09a4 : 100 )
248
249
250 DefineOutput(1, TList::Class()); // Event Counter Histo
251 DefineOutput(2, TTree::Class()); // Event Tree
252 DefineOutput(3, TTree::Class()); // V0 Tree
253 DefineOutput(4, TTree::Class()); // Cascade Tree
254}
255
256
257AliAnalysisTaskStrangenessVsMultiplicity::~AliAnalysisTaskStrangenessVsMultiplicity()
258{
259//------------------------------------------------
260// DESTRUCTOR
261//------------------------------------------------
262
263 if (fListHist){
264 delete fListHist;
265 fListHist = 0x0;
266 }
267 if (fTreeEvent){
268 delete fTreeEvent;
269 fTreeEvent = 0x0;
270 }
271 if (fTreeV0){
272 delete fTreeV0;
273 fTreeV0 = 0x0;
274 }
275 if (fTreeCascade){
276 delete fTreeCascade;
277 fTreeCascade = 0x0;
278 }
279}
280
281//________________________________________________________________________
282void AliAnalysisTaskStrangenessVsMultiplicity::UserCreateOutputObjects()
283{
284
285 OpenFile(2);
286 // Called once
287
288//------------------------------------------------
289
290 fTreeEvent = new TTree("fTreeEvent","Event");
291
292//------------------------------------------------
293// fTree Branch definitions - Event by Event info
294//------------------------------------------------
295
296//-----------BASIC-INFO---------------------------
297
298 //--- VZERO Data (Integrated)
299 fTreeEvent->Branch("fAmplitude_V0A",&fAmplitude_V0A,"fAmplitude_V0A/F");
300 fTreeEvent->Branch("fAmplitude_V0C",&fAmplitude_V0C,"fAmplitude_V0C/F");
301 fTreeEvent->Branch("fAmplitude_V0AEq",&fAmplitude_V0AEq,"fAmplitude_V0AEq/F");
302 fTreeEvent->Branch("fAmplitude_V0CEq",&fAmplitude_V0CEq,"fAmplitude_V0CEq/F");
303
304 //Info from AliCentrality (not necessarily 'centrality' per se)
305 fTreeEvent->Branch("fCentrality_V0A",&fCentrality_V0A,"fCentrality_V0A/F");
306 fTreeEvent->Branch("fCentrality_V0C",&fCentrality_V0C,"fCentrality_V0C/F");
307 fTreeEvent->Branch("fCentrality_V0M",&fCentrality_V0A,"fCentrality_V0M/F");
308 fTreeEvent->Branch("fCentrality_V0AEq",&fCentrality_V0AEq,"fCentrality_V0AEq/F");
309 fTreeEvent->Branch("fCentrality_V0CEq",&fCentrality_V0CEq,"fCentrality_V0CEq/F");
310 fTreeEvent->Branch("fCentrality_V0MEq",&fCentrality_V0AEq,"fCentrality_V0MEq/F");
311
312 //Official GetReferenceMultiplicity
313 fTreeEvent->Branch("fRefMultEta5",&fRefMultEta5,"fRefMultEta5/I");
314 fTreeEvent->Branch("fRefMultEta8",&fRefMultEta8,"fRefMultEta8/I");
315
316 //Run Number
317 fTreeEvent->Branch("fRunNumber", &fRunNumber, "fRunNumber/I");
318
319 //Create Basic V0 Output Tree
320 fTreeV0 = new TTree( "fTreeV0", "V0 Candidates");
321
322 //Create Cascade output tree
323 fTreeCascade = new TTree("fTreeCascade","CascadeCandidates");
324
325//------------------------------------------------
326// fTreeCascade Branch definitions - Cascade Tree
327//------------------------------------------------
328
329//-----------BASIC-INFO---------------------------
330 fTreeCascade->Branch("fTreeCascVarCharge",&fTreeCascVarCharge,"fTreeCascVarCharge/I");
331 fTreeCascade->Branch("fTreeCascVarMassAsXi",&fTreeCascVarMassAsXi,"fTreeCascVarMassAsXi/F");
332 fTreeCascade->Branch("fTreeCascVarMassAsOmega",&fTreeCascVarMassAsOmega,"fTreeCascVarMassAsOmega/F");
333 fTreeCascade->Branch("fTreeCascVarPt",&fTreeCascVarPt,"fTreeCascVarPt/F");
334 fTreeCascade->Branch("fTreeCascVarRapXi",&fTreeCascVarRapXi,"fTreeCascVarRapXi/F");
335 fTreeCascade->Branch("fTreeCascVarRapOmega",&fTreeCascVarRapOmega,"fTreeCascVarRapOmega/F");
336 fTreeCascade->Branch("fTreeCascVarNegEta",&fTreeCascVarNegEta,"fTreeCascVarNegEta/F");
337 fTreeCascade->Branch("fTreeCascVarPosEta",&fTreeCascVarPosEta,"fTreeCascVarPosEta/F");
338 fTreeCascade->Branch("fTreeCascVarBachEta",&fTreeCascVarBachEta,"fTreeCascVarBachEta/F");
339//-----------INFO-FOR-CUTS------------------------
340 fTreeCascade->Branch("fTreeCascVarDCACascDaughters",&fTreeCascVarDCACascDaughters,"fTreeCascVarDCACascDaughters/F");
341 fTreeCascade->Branch("fTreeCascVarDCABachToPrimVtx",&fTreeCascVarDCABachToPrimVtx,"fTreeCascVarDCABachToPrimVtx/F");
342 fTreeCascade->Branch("fTreeCascVarDCAV0Daughters",&fTreeCascVarDCAV0Daughters,"fTreeCascVarDCAV0Daughters/F");
343 fTreeCascade->Branch("fTreeCascVarDCAV0ToPrimVtx",&fTreeCascVarDCAV0ToPrimVtx,"fTreeCascVarDCAV0ToPrimVtx/F");
344 fTreeCascade->Branch("fTreeCascVarDCAPosToPrimVtx",&fTreeCascVarDCAPosToPrimVtx,"fTreeCascVarDCAPosToPrimVtx/F");
345 fTreeCascade->Branch("fTreeCascVarDCANegToPrimVtx",&fTreeCascVarDCANegToPrimVtx,"fTreeCascVarDCANegToPrimVtx/F");
346 fTreeCascade->Branch("fTreeCascVarCascCosPointingAngle",&fTreeCascVarCascCosPointingAngle,"fTreeCascVarCascCosPointingAngle/F");
347 fTreeCascade->Branch("fTreeCascVarCascRadius",&fTreeCascVarCascRadius,"fTreeCascVarCascRadius/F");
348 fTreeCascade->Branch("fTreeCascVarV0Mass",&fTreeCascVarV0Mass,"fTreeCascVarV0Mass/F");
349 fTreeCascade->Branch("fTreeCascVarV0CosPointingAngle",&fTreeCascVarV0CosPointingAngle,"fTreeCascVarV0CosPointingAngle/F");
350 fTreeCascade->Branch("fTreeCascVarV0CosPointingAngleSpecial",&fTreeCascVarV0CosPointingAngleSpecial,"fTreeCascVarV0CosPointingAngleSpecial/F");
351 fTreeCascade->Branch("fTreeCascVarV0Radius",&fTreeCascVarV0Radius,"fTreeCascVarV0Radius/F");
352 fTreeCascade->Branch("fTreeCascVarLeastNbrClusters",&fTreeCascVarLeastNbrClusters,"fTreeCascVarLeastNbrClusters/I");
353//-----------MULTIPLICITY-INFO--------------------
354 fTreeCascade->Branch("fTreeCascVarCentV0A",&fTreeCascVarCentV0A,"fTreeCascVarCentV0A/F");
355 fTreeCascade->Branch("fTreeCascVarCentV0C",&fTreeCascVarCentV0C,"fTreeCascVarCentV0C/F");
356 fTreeCascade->Branch("fTreeCascVarCentV0M",&fTreeCascVarCentV0M,"fTreeCascVarCentV0M/F");
357 fTreeCascade->Branch("fTreeCascVarCentV0AEq",&fTreeCascVarCentV0AEq,"fTreeCascVarCentV0AEq/F");
358 fTreeCascade->Branch("fTreeCascVarCentV0CEq",&fTreeCascVarCentV0CEq,"fTreeCascVarCentV0CEq/F");
359 fTreeCascade->Branch("fTreeCascVarCentV0MEq",&fTreeCascVarCentV0MEq,"fTreeCascVarCentV0MEq/F");
360 fTreeCascade->Branch("fTreeCascVarAmpV0A",&fTreeCascVarAmpV0A,"fTreeCascVarAmpV0A/F");
361 fTreeCascade->Branch("fTreeCascVarAmpV0C",&fTreeCascVarAmpV0C,"fTreeCascVarAmpV0C/F");
362 fTreeCascade->Branch("fTreeCascVarAmpV0AEq",&fTreeCascVarAmpV0AEq,"fTreeCascVarAmpV0AEq/F");
363 fTreeCascade->Branch("fTreeCascVarAmpV0CEq",&fTreeCascVarAmpV0CEq,"fTreeCascVarAmpV0CEq/F");
364 fTreeCascade->Branch("fTreeCascVarRefMultEta8",&fTreeCascVarRefMultEta8,"fTreeCascVarRefMultEta8/I");
365 fTreeCascade->Branch("fTreeCascVarRefMultEta5",&fTreeCascVarRefMultEta5,"fTreeCascVarRefMultEta5/I");
366//-----------DECAY-LENGTH-INFO--------------------
367 fTreeCascade->Branch("fTreeCascVarDistOverTotMom",&fTreeCascVarDistOverTotMom,"fTreeCascVarDistOverTotMom/F");
368//------------------------------------------------
369 fTreeCascade->Branch("fTreeCascVarNegNSigmaPion",&fTreeCascVarNegNSigmaPion,"fTreeCascVarNegNSigmaPion/F");
370 fTreeCascade->Branch("fTreeCascVarNegNSigmaProton",&fTreeCascVarNegNSigmaProton,"fTreeCascVarNegNSigmaProton/F");
371 fTreeCascade->Branch("fTreeCascVarPosNSigmaPion",&fTreeCascVarPosNSigmaPion,"fTreeCascVarPosNSigmaPion/F");
372 fTreeCascade->Branch("fTreeCascVarPosNSigmaProton",&fTreeCascVarPosNSigmaProton,"fTreeCascVarPosNSigmaProton/F");
373 fTreeCascade->Branch("fTreeCascVarBachNSigmaPion",&fTreeCascVarBachNSigmaPion,"fTreeCascVarBachNSigmaPion/F");
374 fTreeCascade->Branch("fTreeCascVarBachNSigmaKaon",&fTreeCascVarBachNSigmaKaon,"fTreeCascVarBachNSigmaKaon/F");
375
376//------------------------------------------------
377// Particle Identification Setup
378//------------------------------------------------
379
380 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
381 AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
382 fPIDResponse = inputHandler->GetPIDResponse();
383
384 // Multiplicity
385 if(! fESDtrackCuts ){
386 fESDtrackCuts = new AliESDtrackCuts();
387 }
388
389//------------------------------------------------
390// V0 Multiplicity Histograms
391//------------------------------------------------
392
393 // Create histograms
394 OpenFile(1);
395 fListHist = new TList();
396 fListHist->SetOwner(); // See http://root.cern.ch/root/html/TCollection.html#TCollection:SetOwner
397
398 if(! fHistEventCounter ) {
399 //Histogram Output: Event-by-Event
400 fHistEventCounter = new TH1D( "fHistEventCounter", ";Evt. Sel. Step;Count",5,0,5);
401 fHistEventCounter->GetXaxis()->SetBinLabel(1, "Processed");
402 fHistEventCounter->GetXaxis()->SetBinLabel(2, "Phys-Sel");
403 fHistEventCounter->GetXaxis()->SetBinLabel(3, "Has Vtx");
404 fHistEventCounter->GetXaxis()->SetBinLabel(4, "Vtx |z|<10cm");
405 fHistEventCounter->GetXaxis()->SetBinLabel(5, "Isn't Pileup");
406 fListHist->Add(fHistEventCounter);
407 }
408
409 //List of Histograms: Normal
410 PostData(1, fListHist);
411
412 //TTree Object: Saved to base directory. Should cache to disk while saving.
413 //(Important to avoid excessive memory usage, particularly when merging)
414 PostData(2, fTreeEvent);
415 PostData(3, fTreeV0);
416 PostData(4, fTreeCascade);
417
418}// end UserCreateOutputObjects
419
420
421//________________________________________________________________________
422void AliAnalysisTaskStrangenessVsMultiplicity::UserExec(Option_t *)
423{
424 // Main loop
425 // Called for each event
426
427 AliESDEvent *lESDevent = 0x0;
428
429 // Connect to the InputEvent
430 // After these lines, we should have an ESD/AOD event + the number of V0s in it.
431
432 // Appropriate for ESD analysis!
433
434 lESDevent = dynamic_cast<AliESDEvent*>( InputEvent() );
435 if (!lESDevent) {
436 AliWarning("ERROR: lESDevent not available \n");
437 return;
438 }
439
440 //Get VZERO Information for multiplicity later
441 AliVVZERO* esdV0 = lESDevent->GetVZEROData();
442 if (!esdV0) {
443 AliError("AliVVZERO not available");
444 return;
445 }
446
447 fRunNumber = lESDevent->GetRunNumber();
448
449 Double_t lMagneticField = -10;
450 lMagneticField = lESDevent->GetMagneticField( );
451
452//------------------------------------------------
453// Variable Definition
454//------------------------------------------------
455
456//------------------------------------------------
457// Physics Selection
458//------------------------------------------------
459
460 fHistEventCounter->Fill(0.5);
461
462 UInt_t maskIsSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
463 Bool_t isSelected = 0;
464 isSelected = (maskIsSelected & AliVEvent::kMB) == AliVEvent::kMB;
465
466 //Standard Min-Bias Selection
467 if ( ! isSelected ) {
468 PostData(1, fListHist);
469 PostData(2, fTreeEvent);
470 PostData(3, fTreeV0);
471 PostData(4, fTreeCascade);
472 return;
473 }
474
475 fHistEventCounter->Fill(1.5);
476
477 //------------------------------------------------
478 // Primary Vertex Requirements Section:
479 // ---> pp: has vertex, |z|<10cm
480 //------------------------------------------------
481
482 //classical Proton-proton like selection
483 const AliESDVertex *lPrimaryBestESDVtx = lESDevent->GetPrimaryVertex();
484 const AliESDVertex *lPrimaryTrackingESDVtx = lESDevent->GetPrimaryVertexTracks();
485 const AliESDVertex *lPrimarySPDVtx = lESDevent->GetPrimaryVertexSPD();
486
487 Double_t lBestPrimaryVtxPos[3] = {-100.0, -100.0, -100.0};
488 lPrimaryBestESDVtx->GetXYZ( lBestPrimaryVtxPos );
489
490 //Only accept if Tracking or SPD vertex is fine
491 if (!lPrimarySPDVtx->GetStatus() && !lPrimaryTrackingESDVtx->GetStatus() ){
492 AliWarning("Pb / No SPD prim. vertex nor prim. Tracking vertex ... return !");
493 PostData(1, fListHist);
494 PostData(2, fTreeEvent);
495 PostData(3, fTreeV0);
496 PostData(4, fTreeCascade);
497 return;
498 }
499
500 //Has SPD or Tracking Vertex
501 fHistEventCounter -> Fill(2.5);
502
503 //Always do Primary Vertex Selection
504 if(TMath::Abs(lBestPrimaryVtxPos[2]) > 10.0) {
505 AliWarning("Pb / | Z position of Best Prim Vtx | > 10.0 cm ... return !");
506 PostData(1, fListHist);
507 PostData(2, fTreeEvent);
508 PostData(3, fTreeV0);
509 PostData(4, fTreeCascade);
510 return;
511 }
512
513 //Fill Event selected counter
514 fHistEventCounter -> Fill(3.5);
515
516 //------------------------------------------------
517 // Check if this isn't pileup
518 //------------------------------------------------
519
520 if(lESDevent->IsPileupFromSPDInMultBins() ){
521 // minContributors=3, minZdist=0.8, nSigmaZdist=3., nSigmaDiamXY=2., nSigmaDiamZ=5.
522 //-> see http://alisoft.cern.ch/viewvc/trunk/STEER/AliESDEvent.h?root=AliRoot&r1=41914&r2=42199&pathrev=42199
523 AliWarning("Pb / Event tagged as pile-up by SPD... return !");
524 PostData(1, fListHist);
525 PostData(2, fTreeEvent);
526 PostData(3, fTreeV0);
527 PostData(4, fTreeCascade);
528 return;
529 }
530 //Fill Event isn't pileup counter
531 fHistEventCounter -> Fill(4.5);
532
533//------------------------------------------------
534// Multiplicity Information Acquistion
535//------------------------------------------------
536
537 //Standard GetReferenceMultiplicity Estimator (0.5 and 0.8)
538 fRefMultEta5 = fESDtrackCuts->GetReferenceMultiplicity(lESDevent, AliESDtrackCuts::kTrackletsITSTPC,0.5);
539 fRefMultEta8 = fESDtrackCuts->GetReferenceMultiplicity(lESDevent, AliESDtrackCuts::kTrackletsITSTPC,0.8);
540
541 // VZERO PART
542 Float_t multV0A = 0; // multiplicity from V0 reco side A
543 Float_t multV0C = 0; // multiplicity from V0 reco side C
544 Float_t multV0AEq = 0; // multiplicity from V0 reco side A
545 Float_t multV0CEq = 0; // multiplicity from V0 reco side C
546 Float_t multV0ACorr = 0; // multiplicity from V0 reco side A
547 Float_t multV0CCorr = 0; // multiplicity from V0 reco side C
548
549 //Non-Equalized Signal: copy of multV0ACorr and multV0CCorr from AliCentralitySelectionTask
550 //Getters for uncorrected multiplicity
551 multV0A=esdV0->GetMTotV0A();
552 multV0C=esdV0->GetMTotV0C();
553
554 //Get Z vertex position of SPD vertex (why not Tracking if available?)
555 Float_t zvtx = lPrimarySPDVtx->GetZ();
556
557 //Acquire Corrected multV0A
558 multV0ACorr = AliESDUtils::GetCorrV0A(multV0A,zvtx);
559 multV0CCorr = AliESDUtils::GetCorrV0C(multV0C,zvtx);
560
561 //Copy to Event Tree for extra information
562 fAmplitude_V0A = multV0ACorr;
563 fAmplitude_V0C = multV0CCorr;
564
565 // Equalized signals // From AliCentralitySelectionTask
566 for(Int_t iCh = 4; iCh < 7; ++iCh) {
567 Double_t mult = lESDevent->GetVZEROEqMultiplicity(iCh);
568 multV0AEq += mult;
569 }
570 for(Int_t iCh = 0; iCh < 3; ++iCh) {
571 Double_t mult = lESDevent->GetVZEROEqMultiplicity(iCh);
572 multV0CEq += mult;
573 }
574 fAmplitude_V0AEq = multV0AEq;
575 fAmplitude_V0CEq = multV0CEq;
576
577 fCentrality_V0A = -100;
578 fCentrality_V0C = -100;
579 fCentrality_V0M = -100;
580 fCentrality_V0AEq = -100;
581 fCentrality_V0CEq = -100;
582 fCentrality_V0MEq = -100;
583
584 //AliCentrality... Check if working?
585 AliCentrality* centrality;
586 centrality = lESDevent->GetCentrality();
587 if ( !(centrality->GetQuality()>1) ){
588 fCentrality_V0A = centrality->GetCentralityPercentile( "V0A" );
589 fCentrality_V0C = centrality->GetCentralityPercentile( "V0C" );
590 fCentrality_V0M = centrality->GetCentralityPercentile( "V0M" );
591 fCentrality_V0AEq = centrality->GetCentralityPercentile( "V0AEq" );
592 fCentrality_V0CEq = centrality->GetCentralityPercentile( "V0CEq" );
593 fCentrality_V0MEq = centrality->GetCentralityPercentile( "V0MEq" );
594 }
595
596 //Event-level fill
597 fTreeEvent->Fill() ;
598
599//------------------------------------------------
600
601//------------------------------------------------
602// Rerun cascade vertexer!
603//------------------------------------------------
604
605 if( fkRunVertexers ){
606 lESDevent->ResetCascades();
607 lESDevent->ResetV0s();
608
609 AliV0vertexer lV0vtxer;
610 AliCascadeVertexer lCascVtxer;
611
612 lV0vtxer.SetDefaultCuts(fV0VertexerSels);
613 lCascVtxer.SetDefaultCuts(fCascadeVertexerSels);
614
615 lV0vtxer.Tracks2V0vertices(lESDevent);
616 lCascVtxer.V0sTracks2CascadeVertices(lESDevent);
617 }
618
619//------------------------------------------------
620// MAIN CASCADE LOOP STARTS HERE
621//------------------------------------------------
622// Code Credit: Antonin Maire (thanks^100)
623// ---> This is an adaptation
624
625 Long_t ncascades = 0;
626 ncascades = lESDevent->GetNumberOfCascades();
627
628 for (Int_t iXi = 0; iXi < ncascades; iXi++){
629 //------------------------------------------------
630 // Initializations
631 //------------------------------------------------
632 //Double_t lTrkgPrimaryVtxRadius3D = -500.0;
633 //Double_t lBestPrimaryVtxRadius3D = -500.0;
634
635 // - 1st part of initialisation : variables needed to store AliESDCascade data members
636 Double_t lEffMassXi = 0. ;
637 //Double_t lChi2Xi = -1. ;
638 Double_t lDcaXiDaughters = -1. ;
639 Double_t lXiCosineOfPointingAngle = -1. ;
640 Double_t lPosXi[3] = { -1000.0, -1000.0, -1000.0 };
641 Double_t lXiRadius = -1000. ;
642
643 // - 2nd part of initialisation : Nbr of clusters within TPC for the 3 daughter cascade tracks
644 Int_t lPosTPCClusters = -1; // For ESD only ...//FIXME : wait for availability in AOD
645 Int_t lNegTPCClusters = -1; // For ESD only ...
646 Int_t lBachTPCClusters = -1; // For ESD only ...
647
648 // - 3rd part of initialisation : about V0 part in cascades
649 Double_t lInvMassLambdaAsCascDghter = 0.;
650 //Double_t lV0Chi2Xi = -1. ;
651 Double_t lDcaV0DaughtersXi = -1.;
652
653 Double_t lDcaBachToPrimVertexXi = -1., lDcaV0ToPrimVertexXi = -1.;
654 Double_t lDcaPosToPrimVertexXi = -1.;
655 Double_t lDcaNegToPrimVertexXi = -1.;
656 Double_t lV0CosineOfPointingAngleXi = -1. ;
657 Double_t lV0CosineOfPointingAngleXiSpecial = -1. ;
658 Double_t lPosV0Xi[3] = { -1000. , -1000., -1000. }; // Position of VO coming from cascade
659 Double_t lV0RadiusXi = -1000.0;
660 Double_t lV0quality = 0.;
661
662 // - 4th part of initialisation : Effective masses
663 Double_t lInvMassXiMinus = 0.;
664 Double_t lInvMassXiPlus = 0.;
665 Double_t lInvMassOmegaMinus = 0.;
666 Double_t lInvMassOmegaPlus = 0.;
667
668 // - 6th part of initialisation : extra info for QA
669 Double_t lXiMomX = 0. , lXiMomY = 0., lXiMomZ = 0.;
670 Double_t lXiTransvMom = 0. ;
671 //Double_t lXiTransvMomMC= 0. ;
672 Double_t lXiTotMom = 0. ;
673
674 Double_t lBachMomX = 0., lBachMomY = 0., lBachMomZ = 0.;
675 //Double_t lBachTransvMom = 0.;
676 //Double_t lBachTotMom = 0.;
677
678 fTreeCascVarNegNSigmaPion = -100;
679 fTreeCascVarNegNSigmaProton = -100;
680 fTreeCascVarPosNSigmaPion = -100;
681 fTreeCascVarPosNSigmaProton = -100;
682 fTreeCascVarBachNSigmaPion = -100;
683 fTreeCascVarBachNSigmaKaon = -100;
684
685 Short_t lChargeXi = -2;
686 //Double_t lV0toXiCosineOfPointingAngle = 0. ;
687
688 Double_t lRapXi = -20.0, lRapOmega = -20.0; // lEta = -20.0, lTheta = 360., lPhi = 720. ;
689 //Double_t lAlphaXi = -200., lPtArmXi = -200.0;
690
691 // -------------------------------------
692 // II.ESD - Calculation Part dedicated to Xi vertices (ESD)
693
694 AliESDcascade *xi = lESDevent->GetCascade(iXi);
695 if (!xi) continue;
696
697 // - II.Step 2 : Assigning the necessary variables for specific AliESDcascade data members (ESD)
698 //-------------
699 lV0quality = 0.;
700 xi->ChangeMassHypothesis(lV0quality , 3312); // default working hypothesis : cascade = Xi- decay
701
702 lEffMassXi = xi->GetEffMassXi();
703 //lChi2Xi = xi->GetChi2Xi();
704 lDcaXiDaughters = xi->GetDcaXiDaughters();
705 lXiCosineOfPointingAngle = xi->GetCascadeCosineOfPointingAngle( lBestPrimaryVtxPos[0],
706 lBestPrimaryVtxPos[1],
707 lBestPrimaryVtxPos[2] );
708 // Take care : the best available vertex should be used (like in AliCascadeVertexer)
709
710 xi->GetXYZcascade( lPosXi[0], lPosXi[1], lPosXi[2] );
711 lXiRadius = TMath::Sqrt( lPosXi[0]*lPosXi[0] + lPosXi[1]*lPosXi[1] );
712
713 // - II.Step 3 : around the tracks : Bach + V0 (ESD)
714 // ~ Necessary variables for ESDcascade data members coming from the ESDv0 part (inheritance)
715 //-------------
716
717 UInt_t lIdxPosXi = (UInt_t) TMath::Abs( xi->GetPindex() );
718 UInt_t lIdxNegXi = (UInt_t) TMath::Abs( xi->GetNindex() );
719 UInt_t lBachIdx = (UInt_t) TMath::Abs( xi->GetBindex() );
720 // Care track label can be negative in MC production (linked with the track quality)
721 // However = normally, not the case for track index ...
722
723 // FIXME : rejection of a double use of a daughter track (nothing but just a crosscheck of what is done in the cascade vertexer)
724 if(lBachIdx == lIdxNegXi) {
725 AliWarning("Pb / Idx(Bach. track) = Idx(Neg. track) ... continue!"); continue;
726 }
727 if(lBachIdx == lIdxPosXi) {
728 AliWarning("Pb / Idx(Bach. track) = Idx(Pos. track) ... continue!"); continue;
729 }
730
731 AliESDtrack *pTrackXi = lESDevent->GetTrack( lIdxPosXi );
732 AliESDtrack *nTrackXi = lESDevent->GetTrack( lIdxNegXi );
733 AliESDtrack *bachTrackXi = lESDevent->GetTrack( lBachIdx );
734
735 if (!pTrackXi || !nTrackXi || !bachTrackXi ) {
736 AliWarning("ERROR: Could not retrieve one of the 3 ESD daughter tracks of the cascade ...");
737 continue;
738 }
739
740 fTreeCascVarPosEta = pTrackXi->Eta();
741 fTreeCascVarNegEta = nTrackXi->Eta();
742 fTreeCascVarBachEta = bachTrackXi->Eta();
743
744 Double_t lBMom[3], lNMom[3], lPMom[3];
745 xi->GetBPxPyPz( lBMom[0], lBMom[1], lBMom[2] );
746 xi->GetPPxPyPz( lPMom[0], lPMom[1], lPMom[2] );
747 xi->GetNPxPyPz( lNMom[0], lNMom[1], lNMom[2] );
748
749 //fTreeCascVarBachTransMom = TMath::Sqrt( lBMom[0]*lBMom[0] + lBMom[1]*lBMom[1] );
750 //fTreeCascVarPosTransMom = TMath::Sqrt( lPMom[0]*lPMom[0] + lPMom[1]*lPMom[1] );
751 //fTreeCascVarNegTransMom = TMath::Sqrt( lNMom[0]*lNMom[0] + lNMom[1]*lNMom[1] );
752
753 //------------------------------------------------
754 // TPC dEdx information
755 //------------------------------------------------
756 fTreeCascVarNegNSigmaPion = fPIDResponse->NumberOfSigmasTPC( nTrackXi, AliPID::kPion );
757 fTreeCascVarNegNSigmaProton = fPIDResponse->NumberOfSigmasTPC( nTrackXi, AliPID::kProton );
758 fTreeCascVarPosNSigmaPion = fPIDResponse->NumberOfSigmasTPC( pTrackXi, AliPID::kPion );
759 fTreeCascVarPosNSigmaProton = fPIDResponse->NumberOfSigmasTPC( pTrackXi, AliPID::kProton );
760 fTreeCascVarBachNSigmaPion = fPIDResponse->NumberOfSigmasTPC( bachTrackXi, AliPID::kPion );
761 fTreeCascVarBachNSigmaKaon = fPIDResponse->NumberOfSigmasTPC( bachTrackXi, AliPID::kKaon );
762
763 //------------------------------------------------
764 // TPC Number of clusters info
765 // --- modified to save the smallest number
766 // --- of TPC clusters for the 3 tracks
767 //------------------------------------------------
768
769 lPosTPCClusters = pTrackXi->GetTPCNcls();
770 lNegTPCClusters = nTrackXi->GetTPCNcls();
771 lBachTPCClusters = bachTrackXi->GetTPCNcls();
772
773 // 1 - Poor quality related to TPCrefit
774 ULong_t pStatus = pTrackXi->GetStatus();
775 ULong_t nStatus = nTrackXi->GetStatus();
776 ULong_t bachStatus = bachTrackXi->GetStatus();
777
778 //fTreeCascVarkITSRefitBachelor = kTRUE;
779 //fTreeCascVarkITSRefitNegative = kTRUE;
780 //fTreeCascVarkITSRefitPositive = kTRUE;
781
782 if ((pStatus&AliESDtrack::kTPCrefit) == 0) { AliWarning("Pb / V0 Pos. track has no TPCrefit ... continue!"); continue; }
783 if ((nStatus&AliESDtrack::kTPCrefit) == 0) { AliWarning("Pb / V0 Neg. track has no TPCrefit ... continue!"); continue; }
784 if ((bachStatus&AliESDtrack::kTPCrefit) == 0) { AliWarning("Pb / Bach. track has no TPCrefit ... continue!"); continue; }
785
786 //Extra Debug Information: booleans for ITS refit
787 //if ((pStatus&AliESDtrack::kITSrefit) == 0) { fTreeCascVarkITSRefitPositive = kFALSE; }
788 //if ((nStatus&AliESDtrack::kITSrefit) == 0) { fTreeCascVarkITSRefitNegative = kFALSE; }
789 //if ((bachStatus&AliESDtrack::kITSrefit) == 0) { fTreeCascVarkITSRefitBachelor = kFALSE; }
790
791 // 2 - Poor quality related to TPC clusters: lowest cut of 70 clusters
792 if(lPosTPCClusters < 70) { AliWarning("Pb / V0 Pos. track has less than 70 TPC clusters ... continue!"); continue; }
793 if(lNegTPCClusters < 70) { AliWarning("Pb / V0 Neg. track has less than 70 TPC clusters ... continue!"); continue; }
794 if(lBachTPCClusters < 70) { AliWarning("Pb / Bach. track has less than 70 TPC clusters ... continue!"); continue; }
795 Int_t leastnumberofclusters = 1000;
796 if( lPosTPCClusters < leastnumberofclusters ) leastnumberofclusters = lPosTPCClusters;
797 if( lNegTPCClusters < leastnumberofclusters ) leastnumberofclusters = lNegTPCClusters;
798 if( lBachTPCClusters < leastnumberofclusters ) leastnumberofclusters = lBachTPCClusters;
799
800 lInvMassLambdaAsCascDghter = xi->GetEffMass();
801 // This value shouldn't change, whatever the working hyp. is : Xi-, Xi+, Omega-, Omega+
802 lDcaV0DaughtersXi = xi->GetDcaV0Daughters();
803 //lV0Chi2Xi = xi->GetChi2V0();
804
805 lV0CosineOfPointingAngleXi = xi->GetV0CosineOfPointingAngle( lBestPrimaryVtxPos[0],
806 lBestPrimaryVtxPos[1],
807 lBestPrimaryVtxPos[2] );
808 //Modification: V0 CosPA wrt to Cascade decay vertex
809 lV0CosineOfPointingAngleXiSpecial = xi->GetV0CosineOfPointingAngle( lPosXi[0],
810 lPosXi[1],
811 lPosXi[2] );
812
813 lDcaV0ToPrimVertexXi = xi->GetD( lBestPrimaryVtxPos[0],
814 lBestPrimaryVtxPos[1],
815 lBestPrimaryVtxPos[2] );
816
817 lDcaBachToPrimVertexXi = TMath::Abs( bachTrackXi->GetD( lBestPrimaryVtxPos[0],
818 lBestPrimaryVtxPos[1],
819 lMagneticField ) );
820 // Note : AliExternalTrackParam::GetD returns an algebraic value ...
821
822 xi->GetXYZ( lPosV0Xi[0], lPosV0Xi[1], lPosV0Xi[2] );
823 lV0RadiusXi = TMath::Sqrt( lPosV0Xi[0]*lPosV0Xi[0] + lPosV0Xi[1]*lPosV0Xi[1] );
824
825 lDcaPosToPrimVertexXi = TMath::Abs( pTrackXi ->GetD( lBestPrimaryVtxPos[0],
826 lBestPrimaryVtxPos[1],
827 lMagneticField ) );
828
829 lDcaNegToPrimVertexXi = TMath::Abs( nTrackXi ->GetD( lBestPrimaryVtxPos[0],
830 lBestPrimaryVtxPos[1],
831 lMagneticField ) );
832
833 // - II.Step 4 : around effective masses (ESD)
834 // ~ change mass hypotheses to cover all the possibilities : Xi-/+, Omega -/+
835
836 if( bachTrackXi->Charge() < 0 ) {
837 lV0quality = 0.;
838 xi->ChangeMassHypothesis(lV0quality , 3312);
839 // Calculate the effective mass of the Xi- candidate.
840 // pdg code 3312 = Xi-
841 lInvMassXiMinus = xi->GetEffMassXi();
842
843 lV0quality = 0.;
844 xi->ChangeMassHypothesis(lV0quality , 3334);
845 // Calculate the effective mass of the Xi- candidate.
846 // pdg code 3334 = Omega-
847 lInvMassOmegaMinus = xi->GetEffMassXi();
848
849 lV0quality = 0.;
850 xi->ChangeMassHypothesis(lV0quality , 3312); // Back to default hyp.
851 }// end if negative bachelor
852
853
854 if( bachTrackXi->Charge() > 0 ){
855 lV0quality = 0.;
856 xi->ChangeMassHypothesis(lV0quality , -3312);
857 // Calculate the effective mass of the Xi+ candidate.
858 // pdg code -3312 = Xi+
859 lInvMassXiPlus = xi->GetEffMassXi();
860
861 lV0quality = 0.;
862 xi->ChangeMassHypothesis(lV0quality , -3334);
863 // Calculate the effective mass of the Xi+ candidate.
864 // pdg code -3334 = Omega+
865 lInvMassOmegaPlus = xi->GetEffMassXi();
866
867 lV0quality = 0.;
868 xi->ChangeMassHypothesis(lV0quality , -3312); // Back to "default" hyp.
869 }// end if positive bachelor
870 // - II.Step 6 : extra info for QA (ESD)
871 // miscellaneous pieces of info that may help regarding data quality assessment.
872 //-------------
873
874 xi->GetPxPyPz( lXiMomX, lXiMomY, lXiMomZ );
875 lXiTransvMom = TMath::Sqrt( lXiMomX*lXiMomX + lXiMomY*lXiMomY );
876 lXiTotMom = TMath::Sqrt( lXiMomX*lXiMomX + lXiMomY*lXiMomY + lXiMomZ*lXiMomZ );
877
878 xi->GetBPxPyPz( lBachMomX, lBachMomY, lBachMomZ );
879 //lBachTransvMom = TMath::Sqrt( lBachMomX*lBachMomX + lBachMomY*lBachMomY );
880 //lBachTotMom = TMath::Sqrt( lBachMomX*lBachMomX + lBachMomY*lBachMomY + lBachMomZ*lBachMomZ );
881
882 lChargeXi = xi->Charge();
883
884 //lV0toXiCosineOfPointingAngle = xi->GetV0CosineOfPointingAngle( lPosXi[0], lPosXi[1], lPosXi[2] );
885
886 lRapXi = xi->RapXi();
887 lRapOmega = xi->RapOmega();
888 //lEta = xi->Eta();
889 //lTheta = xi->Theta() *180.0/TMath::Pi();
890 //lPhi = xi->Phi() *180.0/TMath::Pi();
891 //lAlphaXi = xi->AlphaXi();
892 //lPtArmXi = xi->PtArmXi();
893
894 //------------------------------------------------
895 // Set Variables for adding to tree
896 //------------------------------------------------
897
898 fTreeCascVarCharge = lChargeXi;
899 if(lInvMassXiMinus!=0) fTreeCascVarMassAsXi = lInvMassXiMinus;
900 if(lInvMassXiPlus!=0) fTreeCascVarMassAsXi = lInvMassXiPlus;
901 if(lInvMassOmegaMinus!=0) fTreeCascVarMassAsOmega = lInvMassOmegaMinus;
902 if(lInvMassOmegaPlus!=0) fTreeCascVarMassAsOmega = lInvMassOmegaPlus;
903 fTreeCascVarPt = lXiTransvMom;
904 fTreeCascVarRapXi = lRapXi ;
905 fTreeCascVarRapOmega = lRapOmega ;
906 fTreeCascVarDCACascDaughters = lDcaXiDaughters;
907 fTreeCascVarDCABachToPrimVtx = lDcaBachToPrimVertexXi;
908 fTreeCascVarDCAV0Daughters = lDcaV0DaughtersXi;
909 fTreeCascVarDCAV0ToPrimVtx = lDcaV0ToPrimVertexXi;
910 fTreeCascVarDCAPosToPrimVtx = lDcaPosToPrimVertexXi;
911 fTreeCascVarDCANegToPrimVtx = lDcaNegToPrimVertexXi;
912 fTreeCascVarCascCosPointingAngle = lXiCosineOfPointingAngle;
913 fTreeCascVarCascRadius = lXiRadius;
914 fTreeCascVarV0Mass = lInvMassLambdaAsCascDghter;
915 fTreeCascVarV0CosPointingAngle = lV0CosineOfPointingAngleXi;
916 fTreeCascVarV0CosPointingAngleSpecial = lV0CosineOfPointingAngleXiSpecial;
917 fTreeCascVarV0Radius = lV0RadiusXi;
918 fTreeCascVarLeastNbrClusters = leastnumberofclusters;
919
920 //Copy Multiplicity information
921 fTreeCascVarCentV0A = fCentrality_V0A;
922 fTreeCascVarCentV0C = fCentrality_V0C;
923 fTreeCascVarCentV0M = fCentrality_V0M;
924 fTreeCascVarCentV0AEq = fCentrality_V0AEq;
925 fTreeCascVarCentV0CEq = fCentrality_V0CEq;
926 fTreeCascVarCentV0MEq = fCentrality_V0MEq;
927 fTreeCascVarAmpV0A = fAmplitude_V0A;
928 fTreeCascVarAmpV0C = fAmplitude_V0C;
929 fTreeCascVarAmpV0AEq = fAmplitude_V0AEq;
930 fTreeCascVarAmpV0CEq = fAmplitude_V0CEq;
931 fTreeCascVarRefMultEta8 = fRefMultEta8;
932 fTreeCascVarRefMultEta5 = fRefMultEta5;
933
934 fTreeCascVarDistOverTotMom = TMath::Sqrt(
935 TMath::Power( lPosXi[0] - lBestPrimaryVtxPos[0] , 2) +
936 TMath::Power( lPosXi[1] - lBestPrimaryVtxPos[1] , 2) +
937 TMath::Power( lPosXi[2] - lBestPrimaryVtxPos[2] , 2)
938 );
939 fTreeCascVarDistOverTotMom /= (lXiTotMom+1e-13);
940
941//All vars not specified here: specified elsewhere!
942
943//------------------------------------------------
944// Fill Tree!
945//------------------------------------------------
946
947// The conditional is meant to decrease excessive
948// memory usage! Be careful when loosening the
949// cut!
950
951 //Xi Mass window: 150MeV wide
952 //Omega mass window: 150MeV wide
953
954 if( fkSaveCascadeTree && ( (fTreeCascVarMassAsXi<1.32+0.075&&fTreeCascVarMassAsXi>1.32-0.075) ||
955 (fTreeCascVarMassAsOmega<1.68+0.075&&fTreeCascVarMassAsOmega>1.68-0.075) ) ){
956 fTreeCascade->Fill();
957 }
958
959//------------------------------------------------
960// Fill tree over.
961//------------------------------------------------
962
963 }// end of the Cascade loop (ESD or AOD)
964
965 // Post output data.
966 PostData(1, fListHist);
967 PostData(2, fTreeEvent);
968 PostData(3, fTreeV0);
969 PostData(4, fTreeCascade);
970}
971
972//________________________________________________________________________
973void AliAnalysisTaskStrangenessVsMultiplicity::Terminate(Option_t *)
974{
975 // Draw result to the screen
976 // Called once at the end of the query
977
978 TList *cRetrievedList = 0x0;
979 cRetrievedList = (TList*)GetOutputData(1);
980 if(!cRetrievedList){
981 Printf("ERROR - AliAnalysisTaskStrangenessVsMultiplicity : ouput data container list not available\n");
982 return;
983 }
984
985 fHistEventCounter = dynamic_cast<TH1D*> ( cRetrievedList->FindObject("fHistEventCounter") );
986 if (!fHistEventCounter) {
987 Printf("ERROR - AliAnalysisTaskStrangenessVsMultiplicity : fHistEventCounter not available");
988 return;
989 }
990
991 TCanvas *canCheck = new TCanvas("AliAnalysisTaskStrangenessVsMultiplicity","V0 Multiplicity",10,10,510,510);
992 canCheck->cd(1)->SetLogy();
993
994 fHistEventCounter->SetMarkerStyle(22);
995 fHistEventCounter->DrawCopy("E");
996}
997
998//----------------------------------------------------------------------------
999
1000Double_t AliAnalysisTaskStrangenessVsMultiplicity::MyRapidity(Double_t rE, Double_t rPz) const
1001{
1002 // Local calculation for rapidity
1003 Double_t ReturnValue = -100;
1004 if( (rE-rPz+1.e-13) != 0 && (rE+rPz) != 0 ){
1005 ReturnValue = 0.5*TMath::Log((rE+rPz)/(rE-rPz+1.e-13));
1006 }
1007 return ReturnValue;
1008}