]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ANALYSIS/AliAnalysisTaskESDfilter.cxx
Changes for #81171 Port request- AliReconstruction: ESD friend filtering - Random...
[u/mrichter/AliRoot.git] / ANALYSIS / AliAnalysisTaskESDfilter.cxx
CommitLineData
f7ed7e88 1/**************************************************************************\r
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
3 * *\r
4 * Author: The ALICE Off-line Project. *\r
5 * Contributors are mentioned in the code where appropriate. *\r
6 * *\r
7 * Permission to use, copy, modify and distribute this software and its *\r
8 * documentation strictly for non-commercial purposes is hereby granted *\r
9 * without fee, provided that the above copyright notice appears in all *\r
10 * copies and that both the copyright notice and this permission notice *\r
11 * appear in the supporting documentation. The authors make no claims *\r
12 * about the suitability of this software for any purpose. It is *\r
13 * provided "as is" without express or implied warranty. *\r
14 **************************************************************************/\r
15\r
c7ef1c47 16/* $Id$ */\r
f7ed7e88 17 \r
18#include <TChain.h>\r
19#include <TTree.h>\r
20#include <TList.h>\r
21#include <TArrayI.h>\r
22#include <TRandom.h>\r
23#include <TParticle.h>\r
a5f7aba4 24#include <TFile.h>\r
f7ed7e88 25\r
26#include "AliAnalysisTaskESDfilter.h"\r
27#include "AliAnalysisManager.h"\r
28#include "AliESDEvent.h"\r
f2f21b24 29#include "AliESDRun.h"\r
f7ed7e88 30#include "AliStack.h"\r
31#include "AliAODEvent.h"\r
32#include "AliMCEvent.h"\r
33#include "AliMCEventHandler.h"\r
34#include "AliESDInputHandler.h"\r
35#include "AliAODHandler.h"\r
36#include "AliAODMCParticle.h"\r
37#include "AliAnalysisFilter.h"\r
38#include "AliESDMuonTrack.h"\r
39#include "AliESDVertex.h"\r
fba4b2e2 40#include "AliCentrality.h"\r
ce7adfe9 41#include "AliEventplane.h"\r
f7ed7e88 42#include "AliESDv0.h"\r
43#include "AliESDkink.h"\r
44#include "AliESDcascade.h"\r
45#include "AliESDPmdTrack.h"\r
46#include "AliESDCaloCluster.h"\r
47#include "AliESDCaloCells.h"\r
48#include "AliMultiplicity.h"\r
49#include "AliLog.h"\r
c7ef1c47 50#include "AliCodeTimer.h"\r
a7b1466d 51#include "AliESDtrackCuts.h"\r
c7ef1c47 52#include "AliESDpid.h"\r
53#include "Riostream.h"\r
54\r
f7ed7e88 55ClassImp(AliAnalysisTaskESDfilter)\r
56\r
57////////////////////////////////////////////////////////////////////////\r
58\r
59AliAnalysisTaskESDfilter::AliAnalysisTaskESDfilter():\r
60 AliAnalysisTaskSE(),\r
61 fTrackFilter(0x0),\r
62 fKinkFilter(0x0),\r
63 fV0Filter(0x0),\r
64 fCascadeFilter(0x0),\r
65 fHighPthreshold(0),\r
12588061 66 fPtshape(0x0),\r
26917b0b 67 fEnableFillAOD(kTRUE),\r
c7ef1c47 68fUsedTrack(0x0),\r
69fUsedKink(0x0),\r
70fUsedV0(0x0),\r
71fAODTrackRefs(0x0),\r
72fAODV0VtxRefs(0x0),\r
73fAODV0Refs(0x0),\r
74fMChandler(0x0),\r
75fNumberOfTracks(0),\r
76fNumberOfPositiveTracks(0),\r
77fNumberOfV0s(0),\r
78fNumberOfVertices(0),\r
79fNumberOfCascades(0),\r
80fNumberOfKinks(0),\r
81fOldESDformat(kFALSE),\r
82fPrimaryVertex(0x0),\r
83fTPCOnlyFilterMask(0),\r
84fIsVZEROEnabled(kTRUE),\r
85fAreCascadesEnabled(kTRUE),\r
86fAreV0sEnabled(kTRUE),\r
87fAreKinksEnabled(kTRUE),\r
88fAreTracksEnabled(kTRUE),\r
89fArePmdClustersEnabled(kTRUE),\r
90fAreCaloClustersEnabled(kTRUE),\r
91fAreEMCALCellsEnabled(kTRUE),\r
92fArePHOSCellsEnabled(kTRUE),\r
93fAreTrackletsEnabled(kTRUE),\r
94fESDpid(0x0),\r
95fIsPidOwner(kFALSE),\r
96fTimeZeroType(AliESDpid::kTOF_T0)\r
f7ed7e88 97{\r
98 // Default constructor\r
99}\r
100\r
c7ef1c47 101//______________________________________________________________________________\r
f7ed7e88 102AliAnalysisTaskESDfilter::AliAnalysisTaskESDfilter(const char* name):\r
103 AliAnalysisTaskSE(name),\r
104 fTrackFilter(0x0),\r
105 fKinkFilter(0x0),\r
106 fV0Filter(0x0),\r
107 fCascadeFilter(0x0),\r
108 fHighPthreshold(0),\r
12588061 109 fPtshape(0x0),\r
26917b0b 110 fEnableFillAOD(kTRUE),\r
c7ef1c47 111fUsedTrack(0x0),\r
112fUsedKink(0x0),\r
113fUsedV0(0x0),\r
114fAODTrackRefs(0x0),\r
115fAODV0VtxRefs(0x0),\r
116fAODV0Refs(0x0),\r
117fMChandler(0x0),\r
118fNumberOfTracks(0),\r
119fNumberOfPositiveTracks(0),\r
120fNumberOfV0s(0),\r
121fNumberOfVertices(0),\r
122fNumberOfCascades(0),\r
123fNumberOfKinks(0),\r
124fOldESDformat(kFALSE),\r
125fPrimaryVertex(0x0),\r
126fTPCOnlyFilterMask(0),\r
127fIsVZEROEnabled(kTRUE),\r
128fAreCascadesEnabled(kTRUE),\r
129fAreV0sEnabled(kTRUE),\r
130fAreKinksEnabled(kTRUE),\r
131fAreTracksEnabled(kTRUE),\r
132fArePmdClustersEnabled(kTRUE),\r
133fAreCaloClustersEnabled(kTRUE),\r
134fAreEMCALCellsEnabled(kTRUE),\r
135fArePHOSCellsEnabled(kTRUE),\r
136fAreTrackletsEnabled(kTRUE),\r
137fESDpid(0x0),\r
138fIsPidOwner(kFALSE),\r
139fTimeZeroType(AliESDpid::kTOF_T0)\r
f7ed7e88 140{\r
141 // Constructor\r
142}\r
75a996ab 143AliAnalysisTaskESDfilter::~AliAnalysisTaskESDfilter(){\r
144 if(fIsPidOwner)delete fESDpid;\r
145}\r
c7ef1c47 146//______________________________________________________________________________\r
f7ed7e88 147void AliAnalysisTaskESDfilter::UserCreateOutputObjects()\r
148{\r
cd5eefa7 149 //\r
150 // Create Output Objects conenct filter to outputtree\r
151 // \r
c7ef1c47 152 if(OutputTree())\r
153 {\r
f7ed7e88 154 OutputTree()->GetUserInfo()->Add(fTrackFilter);\r
d3815ca9 155 }\r
c7ef1c47 156 else\r
157 {\r
d3815ca9 158 AliError("No OutputTree() for adding the track filter");\r
159 }\r
f7ed7e88 160}\r
161\r
c7ef1c47 162//______________________________________________________________________________\r
f7ed7e88 163void AliAnalysisTaskESDfilter::Init()\r
164{\r
c7ef1c47 165 // Initialization\r
166 if (fDebug > 1) AliInfo("Init() \n");\r
167 // Call configuration file\r
f7ed7e88 168}\r
169\r
c7ef1c47 170//______________________________________________________________________________\r
171void AliAnalysisTaskESDfilter::PrintTask(Option_t *option, Int_t indent) const\r
172{\r
173 AliInfo("");\r
174 \r
175 AliAnalysisTaskSE::PrintTask(option,indent);\r
176 \r
177 TString spaces(' ',indent+3);\r
178 \r
179 cout << spaces.Data() << Form("Cascades are %s",fAreCascadesEnabled ? "ENABLED":"DISABLED") << endl;\r
180 cout << spaces.Data() << Form("V0s are %s",fAreV0sEnabled ? "ENABLED":"DISABLED") << endl;\r
181 cout << spaces.Data() << Form("Kinks are %s",fAreKinksEnabled ? "ENABLED":"DISABLED") << endl;\r
182 cout << spaces.Data() << Form("Tracks are %s",fAreTracksEnabled ? "ENABLED":"DISABLED") << endl;\r
183 cout << spaces.Data() << Form("PmdClusters are %s",fArePmdClustersEnabled ? "ENABLED":"DISABLED") << endl;\r
184 cout << spaces.Data() << Form("CaloClusters are %s",fAreCaloClustersEnabled ? "ENABLED":"DISABLED") << endl;\r
185 cout << spaces.Data() << Form("EMCAL cells are %s",fAreEMCALCellsEnabled ? "ENABLED":"DISABLED") << endl;\r
186 cout << spaces.Data() << Form("Tracklets are %s",fAreTrackletsEnabled ? "ENABLED":"DISABLED") << endl; \r
187}\r
f7ed7e88 188\r
c7ef1c47 189//______________________________________________________________________________\r
f7ed7e88 190void AliAnalysisTaskESDfilter::UserExec(Option_t */*option*/)\r
191{\r
192// Execute analysis for current event\r
193//\r
194 \r
195 Long64_t ientry = Entry();\r
f2f21b24 196 \r
e26c0b71 197 if (fDebug > 0) {\r
f2f21b24 198 printf("Filter: Analysing event # %5d\n", (Int_t) ientry);\r
199 if (fHighPthreshold == 0) AliInfo("detector PID signals are stored in each track");\r
200 if (!fPtshape) AliInfo("detector PID signals are not stored below the pt threshold");\r
e26c0b71 201 }\r
5b173801 202 // Filters must explicitely enable AOD filling in their UserExec (AG)\r
203 if (!AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()) AliFatal("Cannot run ESD filter without an output event handler");\r
307b6330 204 if(fEnableFillAOD) {\r
205 AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillAOD(kTRUE);\r
206 AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillExtension(kTRUE);\r
207 } \r
f7ed7e88 208 ConvertESDtoAOD();\r
209}\r
210\r
c7ef1c47 211//______________________________________________________________________________\r
212TClonesArray& AliAnalysisTaskESDfilter::Cascades()\r
213{\r
214 return *(AODEvent()->GetCascades());\r
215}\r
f7ed7e88 216\r
c7ef1c47 217//______________________________________________________________________________\r
218TClonesArray& AliAnalysisTaskESDfilter::Tracks()\r
219{\r
220 return *(AODEvent()->GetTracks());\r
221}\r
f7ed7e88 222\r
c7ef1c47 223//______________________________________________________________________________\r
224TClonesArray& AliAnalysisTaskESDfilter::V0s()\r
225{\r
226 return *(AODEvent()->GetV0s());\r
227}\r
f7ed7e88 228\r
c7ef1c47 229//______________________________________________________________________________\r
230TClonesArray& AliAnalysisTaskESDfilter::Vertices()\r
231{\r
232 return *(AODEvent()->GetVertices());\r
233}\r
f7ed7e88 234\r
c7ef1c47 235//______________________________________________________________________________\r
236AliAODHeader* AliAnalysisTaskESDfilter::ConvertHeader(const AliESDEvent& esd)\r
237{\r
238 AliCodeTimerAuto("",0);\r
239 \r
240 AliAODHeader* header = AODEvent()->GetHeader();\r
241 \r
242 header->SetRunNumber(esd.GetRunNumber());\r
243 header->SetOfflineTrigger(fInputHandler->IsEventSelected()); // propagate the decision of the physics selection\r
f7ed7e88 244\r
c7ef1c47 245 TTree* tree = fInputHandler->GetTree();\r
246 if (tree) {\r
247 TFile* file = tree->GetCurrentFile();\r
248 if (file) header->SetESDFileName(file->GetName());\r
249 }\r
250 \r
251 if (fOldESDformat) {\r
252 header->SetBunchCrossNumber(0);\r
253 header->SetOrbitNumber(0);\r
254 header->SetPeriodNumber(0);\r
255 header->SetEventType(0);\r
256 header->SetMuonMagFieldScale(-999.);\r
257 header->SetCentrality(0); \r
ce7adfe9 258 header->SetEventplane(0);\r
c7ef1c47 259 } else {\r
260 header->SetBunchCrossNumber(esd.GetBunchCrossNumber());\r
261 header->SetOrbitNumber(esd.GetOrbitNumber());\r
262 header->SetPeriodNumber(esd.GetPeriodNumber());\r
263 header->SetEventType(esd.GetEventType());\r
f7ed7e88 264 \r
c7ef1c47 265 header->SetEventNumberESDFile(esd.GetHeader()->GetEventNumberInFile());\r
266 if(const_cast<AliESDEvent&>(esd).GetCentrality()){\r
75a996ab 267 header->SetCentrality(const_cast<AliESDEvent&>(esd).GetCentrality());\r
a98c78e5 268 }\r
c7ef1c47 269 else{\r
270 header->SetCentrality(0);\r
a98c78e5 271 }\r
ce7adfe9 272 if(const_cast<AliESDEvent&>(esd).GetEventplane()){\r
273 header->SetEventplane(const_cast<AliESDEvent&>(esd).GetEventplane());\r
274 }\r
275 else{\r
276 header->SetEventplane(0);\r
277 }\r
c7ef1c47 278 }\r
279 \r
280 // Trigger\r
281 header->SetFiredTriggerClasses(esd.GetFiredTriggerClasses());\r
282 header->SetTriggerMask(esd.GetTriggerMask()); \r
283 header->SetTriggerCluster(esd.GetTriggerCluster());\r
284 header->SetL0TriggerInputs(esd.GetHeader()->GetL0TriggerInputs()); \r
285 header->SetL1TriggerInputs(esd.GetHeader()->GetL1TriggerInputs()); \r
286 header->SetL2TriggerInputs(esd.GetHeader()->GetL2TriggerInputs()); \r
287 \r
288 header->SetMagneticField(esd.GetMagneticField());\r
289 header->SetMuonMagFieldScale(esd.GetCurrentDip()/6000.);\r
290 header->SetZDCN1Energy(esd.GetZDCN1Energy());\r
291 header->SetZDCP1Energy(esd.GetZDCP1Energy());\r
292 header->SetZDCN2Energy(esd.GetZDCN2Energy());\r
293 header->SetZDCP2Energy(esd.GetZDCP2Energy());\r
294 header->SetZDCEMEnergy(esd.GetZDCEMEnergy(0),esd.GetZDCEMEnergy(1));\r
295 \r
f7ce7a37 296 // ITS Cluster Multiplicty\r
297 const AliMultiplicity *mult = esd.GetMultiplicity();\r
298 for (Int_t ilay = 0; ilay < 6; ilay++) header->SetITSClusters(ilay, mult->GetNumberOfITSClusters(ilay));\r
c7ef1c47 299 \r
f7ce7a37 300 // TPC only Reference Multiplicty\r
301 AliESDtrackCuts* cuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();\r
302 Int_t refMult = cuts ? (Short_t)cuts->GetReferenceMultiplicity(&esd, kTRUE) : -1;\r
303 header->SetTPConlyRefMultiplicity(refMult);\r
304 \r
305 //\r
c7ef1c47 306 Float_t diamxy[2]={esd.GetDiamondX(),esd.GetDiamondY()};\r
307 Float_t diamcov[3]; \r
308 esd.GetDiamondCovXY(diamcov);\r
309 header->SetDiamond(diamxy,diamcov);\r
310 header->SetDiamondZ(esd.GetDiamondZ(),esd.GetSigma2DiamondZ());\r
311 \r
312 return header;\r
313}\r
a98c78e5 314\r
c7ef1c47 315//______________________________________________________________________________\r
316void AliAnalysisTaskESDfilter::ConvertCascades(const AliESDEvent& esd) \r
317{\r
318 // Convert the cascades part of the ESD.\r
319 // Return the number of cascades\r
320 \r
321 AliCodeTimerAuto("",0);\r
322 \r
323 // Create vertices starting from the most complex objects\r
324 Double_t chi2 = 0.;\r
325 \r
326 const AliESDVertex* vtx = esd.GetPrimaryVertex();\r
327 Double_t pos[3] = { 0. };\r
328 Double_t covVtx[6] = { 0. };\r
329 Double_t momBach[3]={0.};\r
330 Double_t covTr[21]={0.};\r
331 Double_t pid[10]={0.};\r
332 AliAODPid* detpid(0x0);\r
333 AliAODVertex* vV0FromCascade(0x0);\r
334 AliAODv0* aodV0(0x0);\r
335 AliAODcascade* aodCascade(0x0);\r
336 AliAODTrack* aodTrack(0x0);\r
337 Double_t momPos[3]={0.};\r
338 Double_t momNeg[3] = { 0. };\r
339 Double_t momPosAtV0vtx[3]={0.};\r
340 Double_t momNegAtV0vtx[3]={0.};\r
341\r
342 TClonesArray& verticesArray = Vertices();\r
343 TClonesArray& tracksArray = Tracks();\r
344 TClonesArray& cascadesArray = Cascades();\r
345 \r
346 // Cascades (Modified by A.Maire - February 2009)\r
347 for (Int_t nCascade = 0; nCascade < esd.GetNumberOfCascades(); ++nCascade) {\r
6878d51e 348 \r
c7ef1c47 349 // 0- Preparation\r
350 //\r
351 AliESDcascade *esdCascade = esd.GetCascade(nCascade);\r
352 Int_t idxPosFromV0Dghter = esdCascade->GetPindex();\r
353 Int_t idxNegFromV0Dghter = esdCascade->GetNindex();\r
354 Int_t idxBachFromCascade = esdCascade->GetBindex();\r
b13bcc35 355 \r
c7ef1c47 356 AliESDtrack *esdCascadePos = esd.GetTrack( idxPosFromV0Dghter);\r
357 AliESDtrack *esdCascadeNeg = esd.GetTrack( idxNegFromV0Dghter);\r
358 AliESDtrack *esdCascadeBach = esd.GetTrack( idxBachFromCascade);\r
253cab77 359 \r
c7ef1c47 360 // Identification of the V0 within the esdCascade (via both daughter track indices)\r
361 AliESDv0 * currentV0 = 0x0;\r
362 Int_t idxV0FromCascade = -1;\r
363 \r
364 for (Int_t iV0=0; iV0<esd.GetNumberOfV0s(); ++iV0) {\r
365 \r
366 currentV0 = esd.GetV0(iV0);\r
367 Int_t posCurrentV0 = currentV0->GetPindex();\r
368 Int_t negCurrentV0 = currentV0->GetNindex();\r
369 \r
370 if (posCurrentV0==idxPosFromV0Dghter && negCurrentV0==idxNegFromV0Dghter) {\r
371 idxV0FromCascade = iV0;\r
372 break;\r
373 }\r
b13bcc35 374 }\r
f7ed7e88 375 \r
c7ef1c47 376 if(idxV0FromCascade < 0){\r
377 printf("Cascade - no matching for the V0 (index V0 = -1) ! Skip ... \n");\r
378 continue;\r
379 }// a priori, useless check, but safer ... in case of pb with tracks "out of bounds"\r
380 \r
381 AliESDv0 *esdV0FromCascade = esd.GetV0(idxV0FromCascade);\r
382 \r
383 // 1 - Cascade selection \r
384 \r
385 // AliESDVertex *esdPrimVtx = new AliESDVertex(*(esd.GetPrimaryVertex()));\r
386 // TList cascadeObjects;\r
387 // cascadeObjects.AddAt(esdV0FromCascade, 0);\r
388 // cascadeObjects.AddAt(esdCascadePos, 1);\r
389 // cascadeObjects.AddAt(esdCascadeNeg, 2);\r
390 // cascadeObjects.AddAt(esdCascade, 3);\r
391 // cascadeObjects.AddAt(esdCascadeBach, 4);\r
392 // cascadeObjects.AddAt(esdPrimVtx, 5);\r
393 // \r
394 // UInt_t selectCascade = 0;\r
395 // if (fCascadeFilter) {\r
396 // // selectCascade = fCascadeFilter->IsSelected(&cascadeObjects); \r
397 // // FIXME AliESDCascadeCuts to be implemented ...\r
398 // \r
399 // // Here we may encounter a moot point at the V0 level \r
400 // // between the cascade selections and the V0 ones :\r
401 // // the V0 selected along with the cascade (secondary V0) may \r
402 // // usually be removed from the dedicated V0 selections (prim V0) ...\r
403 // // -> To be discussed !\r
404 // \r
405 // // this is a little awkward but otherwise the \r
406 // // list wants to access the pointer (delete it) \r
407 // // again when going out of scope\r
408 // delete cascadeObjects.RemoveAt(5); // esdPrimVtx created via copy construct\r
409 // esdPrimVtx = 0;\r
410 // if (!selectCascade) \r
411 // continue;\r
412 // }\r
413 // else{\r
414 // delete cascadeObjects.RemoveAt(5); // esdPrimVtx created via copy construct\r
415 // esdPrimVtx = 0;\r
416 // }\r
417 \r
418 // 2 - Add the cascade vertex\r
419 \r
420 esdCascade->GetXYZcascade(pos[0], pos[1], pos[2]);\r
421 esdCascade->GetPosCovXi(covVtx);\r
422 chi2 = esdCascade->GetChi2Xi(); \r
423 \r
424 AliAODVertex *vCascade = new(verticesArray[fNumberOfVertices++]) AliAODVertex( pos,\r
425 covVtx,\r
426 chi2, // FIXME = Chi2/NDF will be needed\r
427 fPrimaryVertex,\r
428 nCascade, // id\r
429 AliAODVertex::kCascade);\r
430 fPrimaryVertex->AddDaughter(vCascade);\r
431 \r
432// if (fDebug > 2) {\r
433// printf("---- Cascade / Cascade Vertex (AOD) : \n");\r
434// vCascade->Print();\r
435// }\r
436 \r
437 if(esd.GetTOFHeader() && fIsPidOwner) fESDpid->SetTOFResponse(const_cast<AliESDEvent*>(&esd), AliESDpid::kBest_T0); //in case of AOD production strating form LHC10e without Tender. \r
f7ed7e88 438\r
f7ed7e88 439\r
c7ef1c47 440 // 3 - Add the bachelor track from the cascade\r
441 \r
442 if (!fUsedTrack[idxBachFromCascade]) {\r
443 \r
444 esdCascadeBach->GetPxPyPz(momBach);\r
445 esdCascadeBach->GetXYZ(pos);\r
446 esdCascadeBach->GetCovarianceXYZPxPyPz(covTr);\r
447 esdCascadeBach->GetESDpid(pid);\r
448 \r
449 fUsedTrack[idxBachFromCascade] = kTRUE;\r
f7ed7e88 450 UInt_t selectInfo = 0;\r
9eeae5d5 451 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadeBach);\r
c7ef1c47 452 if (fMChandler) fMChandler->SelectParticle(esdCascadeBach->GetLabel());\r
453 aodTrack = new(tracksArray[fNumberOfTracks++]) AliAODTrack(esdCascadeBach->GetID(),\r
454 esdCascadeBach->GetLabel(), \r
455 momBach, \r
456 kTRUE,\r
457 pos,\r
458 kFALSE, // Why kFALSE for "isDCA" ? FIXME\r
459 covTr, \r
460 (Short_t)esdCascadeBach->GetSign(),\r
461 esdCascadeBach->GetITSClusterMap(), \r
462 pid,\r
463 vCascade,\r
464 kTRUE, // usedForVtxFit = kFALSE ? FIXME\r
465 vtx->UsesTrack(esdCascadeBach->GetID()),\r
466 AliAODTrack::kSecondary,\r
467 selectInfo);\r
f2f21b24 468 aodTrack->SetTPCClusterMap(esdCascadeBach->GetTPCClusterMap());\r
469 aodTrack->SetTPCSharedMap (esdCascadeBach->GetTPCSharedMap());\r
40a060e5 470 aodTrack->SetChi2perNDF(Chi2perNDF(esdCascadeBach));\r
3c01c166 471 aodTrack->SetTPCPointsF(esdCascadeBach->GetTPCNclsF());\r
c7ef1c47 472 fAODTrackRefs->AddAt(aodTrack,idxBachFromCascade);\r
f2f21b24 473 \r
c7ef1c47 474 if (esdCascadeBach->GetSign() > 0) ++fNumberOfPositiveTracks;\r
f7ed7e88 475 aodTrack->ConvertAliPIDtoAODPID();\r
476 aodTrack->SetFlags(esdCascadeBach->GetStatus());\r
1a3579a6 477 SetAODPID(esdCascadeBach,aodTrack,detpid);\r
c7ef1c47 478 }\r
479 else {\r
480 aodTrack = static_cast<AliAODTrack*>( fAODTrackRefs->At(idxBachFromCascade) );\r
481 }\r
482 \r
483 vCascade->AddDaughter(aodTrack);\r
484 \r
485// if (fDebug > 4) {\r
486// printf("---- Cascade / bach dghter : \n");\r
487// aodTrack->Print();\r
488// }\r
489 \r
490 \r
491 // 4 - Add the V0 from the cascade. \r
492 // = V0vtx + both pos and neg daughter tracks + the aodV0 itself\r
493 //\r
494 \r
495 if ( !fUsedV0[idxV0FromCascade] ) {\r
496 // 4.A - if VO structure hasn't been created yet\r
497 \r
498 // 4.A.1 - Create the V0 vertex of the cascade\r
499 \r
500 esdV0FromCascade->GetXYZ(pos[0], pos[1], pos[2]);\r
501 esdV0FromCascade->GetPosCov(covVtx);\r
502 chi2 = esdV0FromCascade->GetChi2V0(); // = chi2/NDF since NDF = 2*2-3 ?\r
f7ed7e88 503 \r
c7ef1c47 504 vV0FromCascade = new(verticesArray[fNumberOfVertices++]) AliAODVertex(pos,\r
505 covVtx,\r
506 chi2,\r
507 vCascade,\r
508 idxV0FromCascade, //id of ESDv0\r
509 AliAODVertex::kV0);\r
510 // Note:\r
511 // one V0 can be used by several cascades.\r
512 // So, one AOD V0 vtx can have several parent vtx.\r
513 // This is not directly allowed by AliAODvertex.\r
514 // Setting the parent vtx (here = param "vCascade") doesn't lead to a crash\r
515 // but to a problem of consistency within AODEvent.\r
516 // -> See below paragraph 4.B, for the proposed treatment of such a case.\r
517 \r
518 // Add the vV0FromCascade to the aodVOVtxRefs\r
519 fAODV0VtxRefs->AddAt(vV0FromCascade,idxV0FromCascade);\r
520 \r
521 \r
522 // 4.A.2 - Add the positive tracks from the V0\r
523 \r
524 esdCascadePos->GetPxPyPz(momPos);\r
525 esdCascadePos->GetXYZ(pos);\r
526 esdCascadePos->GetCovarianceXYZPxPyPz(covTr);\r
527 esdCascadePos->GetESDpid(pid);\r
528 \r
529 \r
530 if (!fUsedTrack[idxPosFromV0Dghter]) {\r
531 fUsedTrack[idxPosFromV0Dghter] = kTRUE;\r
532 \r
533 UInt_t selectInfo = 0;\r
534 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadePos);\r
535 if(fMChandler) fMChandler->SelectParticle(esdCascadePos->GetLabel());\r
536 aodTrack = new(tracksArray[fNumberOfTracks++]) \r
537 AliAODTrack( esdCascadePos->GetID(),\r
538 esdCascadePos->GetLabel(), \r
539 momPos, \r
540 kTRUE,\r
541 pos,\r
542 kFALSE, // Why kFALSE for "isDCA" ? FIXME\r
543 covTr, \r
544 (Short_t)esdCascadePos->GetSign(),\r
545 esdCascadePos->GetITSClusterMap(), \r
546 pid,\r
547 vV0FromCascade,\r
548 kTRUE, // usedForVtxFit = kFALSE ? FIXME\r
549 vtx->UsesTrack(esdCascadePos->GetID()),\r
550 AliAODTrack::kSecondary,\r
551 selectInfo);\r
552 aodTrack->SetTPCClusterMap(esdCascadePos->GetTPCClusterMap());\r
553 aodTrack->SetTPCSharedMap (esdCascadePos->GetTPCSharedMap());\r
554 aodTrack->SetChi2perNDF(Chi2perNDF(esdCascadePos));\r
3c01c166 555 aodTrack->SetTPCPointsF(esdCascadePos->GetTPCNclsF());\r
c7ef1c47 556 fAODTrackRefs->AddAt(aodTrack,idxPosFromV0Dghter);\r
557 \r
558 if (esdCascadePos->GetSign() > 0) ++fNumberOfPositiveTracks;\r
559 aodTrack->ConvertAliPIDtoAODPID();\r
560 aodTrack->SetFlags(esdCascadePos->GetStatus());\r
1a3579a6 561 SetAODPID(esdCascadePos,aodTrack,detpid);\r
c7ef1c47 562 }\r
563 else {\r
564 aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(idxPosFromV0Dghter));\r
565 }\r
566 vV0FromCascade->AddDaughter(aodTrack);\r
567 \r
568 \r
569 // 4.A.3 - Add the negative tracks from the V0\r
570 \r
571 esdCascadeNeg->GetPxPyPz(momNeg);\r
572 esdCascadeNeg->GetXYZ(pos);\r
573 esdCascadeNeg->GetCovarianceXYZPxPyPz(covTr);\r
574 esdCascadeNeg->GetESDpid(pid);\r
575 \r
576 \r
577 if (!fUsedTrack[idxNegFromV0Dghter]) {\r
578 fUsedTrack[idxNegFromV0Dghter] = kTRUE;\r
579 \r
580 UInt_t selectInfo = 0;\r
581 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadeNeg);\r
582 if(fMChandler)fMChandler->SelectParticle(esdCascadeNeg->GetLabel());\r
583 aodTrack = new(tracksArray[fNumberOfTracks++]) AliAODTrack( esdCascadeNeg->GetID(),\r
584 esdCascadeNeg->GetLabel(),\r
585 momNeg,\r
586 kTRUE,\r
587 pos,\r
588 kFALSE, // Why kFALSE for "isDCA" ? FIXME\r
589 covTr, \r
590 (Short_t)esdCascadeNeg->GetSign(),\r
591 esdCascadeNeg->GetITSClusterMap(), \r
592 pid,\r
593 vV0FromCascade,\r
594 kTRUE, // usedForVtxFit = kFALSE ? FIXME\r
595 vtx->UsesTrack(esdCascadeNeg->GetID()),\r
596 AliAODTrack::kSecondary,\r
597 selectInfo);\r
598 aodTrack->SetTPCClusterMap(esdCascadeNeg->GetTPCClusterMap());\r
599 aodTrack->SetTPCSharedMap (esdCascadeNeg->GetTPCSharedMap());\r
600 aodTrack->SetChi2perNDF(Chi2perNDF(esdCascadeNeg));\r
3c01c166 601 aodTrack->SetTPCPointsF(esdCascadeNeg->GetTPCNclsF());\r
c7ef1c47 602 fAODTrackRefs->AddAt(aodTrack,idxNegFromV0Dghter);\r
603 \r
604 if (esdCascadeNeg->GetSign() > 0) ++fNumberOfPositiveTracks;\r
605 aodTrack->ConvertAliPIDtoAODPID();\r
606 aodTrack->SetFlags(esdCascadeNeg->GetStatus());\r
1a3579a6 607 SetAODPID(esdCascadeNeg,aodTrack,detpid);\r
c7ef1c47 608 }\r
609 else {\r
610 aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(idxNegFromV0Dghter));\r
611 }\r
612 \r
613 vV0FromCascade->AddDaughter(aodTrack);\r
614 \r
f7ed7e88 615 \r
c7ef1c47 616 // 4.A.4 - Add the V0 from cascade to the V0 array\r
617 \r
618 Double_t dcaV0Daughters = esdV0FromCascade->GetDcaV0Daughters();\r
619 Double_t dcaV0ToPrimVertex = esdV0FromCascade->GetD( esd.GetPrimaryVertex()->GetX(),\r
620 esd.GetPrimaryVertex()->GetY(),\r
621 esd.GetPrimaryVertex()->GetZ() );\r
622 esdV0FromCascade->GetPPxPyPz( momPosAtV0vtx[0],momPosAtV0vtx[1],momPosAtV0vtx[2] ); \r
623 esdV0FromCascade->GetNPxPyPz( momNegAtV0vtx[0],momNegAtV0vtx[1],momNegAtV0vtx[2] ); \r
624 \r
625 Double_t dcaDaughterToPrimVertex[2] = { 999., 999.}; // ..[0] = DCA in (x,y) for Pos and ..[1] = Neg\r
626 dcaDaughterToPrimVertex[0] = TMath::Abs(esdCascadePos->GetD( esd.GetPrimaryVertex()->GetX(),\r
627 esd.GetPrimaryVertex()->GetY(),\r
628 esd.GetMagneticField()) );\r
629 dcaDaughterToPrimVertex[1] = TMath::Abs(esdCascadeNeg->GetD( esd.GetPrimaryVertex()->GetX(),\r
630 esd.GetPrimaryVertex()->GetY(),\r
631 esd.GetMagneticField()) );\r
632 \r
633 aodV0 = new(V0s()[fNumberOfV0s++]) AliAODv0( vV0FromCascade, \r
634 dcaV0Daughters,\r
635 dcaV0ToPrimVertex, \r
636 momPosAtV0vtx, \r
637 momNegAtV0vtx, \r
638 dcaDaughterToPrimVertex); \r
639 // set the aod v0 on-the-fly status\r
640 aodV0->SetOnFlyStatus(esdV0FromCascade->GetOnFlyStatus());\r
641 \r
642 // Add the aodV0 to the aodVORefs\r
643 fAODV0Refs->AddAt(aodV0,idxV0FromCascade);\r
644 \r
645 fUsedV0[idxV0FromCascade] = kTRUE;\r
646 \r
647 } else { \r
648 // 4.B - if V0 structure already used\r
649 \r
650 // Note :\r
651 // one V0 can be used by several cascades (frequent in PbPb evts) : \r
652 // same V0 which used but attached to different bachelor tracks\r
653 // -> aodVORefs and fAODV0VtxRefs are needed.\r
654 // Goal : avoid a redundancy of the info in "Vertices" and "v0s" clones array.\r
655 \r
656 vV0FromCascade = static_cast<AliAODVertex*>( fAODV0VtxRefs->At(idxV0FromCascade) );\r
657 aodV0 = static_cast<AliAODv0*> ( fAODV0Refs ->At(idxV0FromCascade) );\r
658 \r
659 // - Treatment of the parent for such a "re-used" V0 :\r
660 // Insert the cascade that reuses the V0 vertex in the lineage chain\r
661 // Before : vV0 -> vCascade1 -> vPrimary\r
662 // - Hyp : cascade2 uses the same V0 as cascade1\r
663 // After : vV0 -> vCascade2 -> vCascade1 -> vPrimary\r
664 \r
665 AliAODVertex *vCascadePreviousParent = static_cast<AliAODVertex*> (vV0FromCascade->GetParent());\r
f7ed7e88 666 vV0FromCascade->SetParent(vCascade);\r
667 vCascade ->SetParent(vCascadePreviousParent);\r
c7ef1c47 668 \r
669// if(fDebug > 2) \r
670// printf("---- Cascade / Lineage insertion\n"\r
671// "Parent of V0 vtx = Cascade vtx %p\n"\r
672// "Parent of the cascade vtx = Cascade vtx %p\n"\r
673// "Parent of the parent cascade vtx = Cascade vtx %p\n", \r
674// static_cast<void*> (vV0FromCascade->GetParent()),\r
675// static_cast<void*> (vCascade->GetParent()),\r
676// static_cast<void*> (vCascadePreviousParent->GetParent()) );\r
677 \r
678 }// end if V0 structure already used\r
679 \r
680// if (fDebug > 2) {\r
681// printf("---- Cascade / V0 vertex: \n");\r
682// vV0FromCascade->Print();\r
683// }\r
684// \r
685// if (fDebug > 4) {\r
686// printf("---- Cascade / pos dghter : \n");\r
687// aodTrack->Print();\r
688// printf("---- Cascade / neg dghter : \n");\r
689// aodTrack->Print();\r
690// printf("---- Cascade / aodV0 : \n");\r
691// aodV0->Print();\r
692// }\r
693 \r
694 // In any case (used V0 or not), add the V0 vertex to the cascade one.\r
695 vCascade->AddDaughter(vV0FromCascade); \r
696 \r
f7ed7e88 697 \r
c7ef1c47 698 // 5 - Add the primary track of the cascade (if any)\r
699 \r
700 \r
701 // 6 - Add the cascade to the AOD array of cascades\r
702 \r
703 Double_t dcaBachToPrimVertexXY = TMath::Abs(esdCascadeBach->GetD(esd.GetPrimaryVertex()->GetX(),\r
704 esd.GetPrimaryVertex()->GetY(),\r
705 esd.GetMagneticField()) );\r
706 \r
707 Double_t momBachAtCascadeVtx[3]={0.};\r
f7ed7e88 708\r
c7ef1c47 709 esdCascade->GetBPxPyPz(momBachAtCascadeVtx[0], momBachAtCascadeVtx[1], momBachAtCascadeVtx[2]);\r
710 \r
711 aodCascade = new(cascadesArray[fNumberOfCascades++]) AliAODcascade( vCascade,\r
712 esdCascade->Charge(),\r
713 esdCascade->GetDcaXiDaughters(),\r
714 -999.,\r
715 // DCAXiToPrimVtx -> needs to be calculated ----|\r
716 // doesn't exist at ESD level;\r
717 // See AODcascade::DcaXiToPrimVertex(Double, Double, Double)\r
718 dcaBachToPrimVertexXY,\r
719 momBachAtCascadeVtx,\r
720 *aodV0);\r
721 \r
40771071 722 if (fDebug > 10) {\r
723 printf("---- Cascade / AOD cascade : \n\n");\r
724 aodCascade->PrintXi(fPrimaryVertex->GetX(), fPrimaryVertex->GetY(), fPrimaryVertex->GetZ());\r
725 }\r
c7ef1c47 726 \r
727 } // end of the loop on cascades\r
728 \r
729 Cascades().Expand(fNumberOfCascades);\r
730}\r
f7ed7e88 731\r
c7ef1c47 732//______________________________________________________________________________\r
733void AliAnalysisTaskESDfilter::ConvertV0s(const AliESDEvent& esd)\r
734{\r
735 // Access to the AOD container of V0s\r
736 \r
737 AliCodeTimerAuto("",0);\r
f7ed7e88 738\r
c7ef1c47 739 //\r
740 // V0s\r
741 //\r
742 \r
743 Double_t pos[3] = { 0. }; \r
744 Double_t chi2(0.0);\r
745 Double_t covVtx[6] = { 0. };\r
746 Double_t momPos[3]={0.};\r
747 Double_t covTr[21]={0.};\r
748 Double_t pid[10]={0.};\r
749 AliAODTrack* aodTrack(0x0);\r
750 AliAODPid* detpid(0x0);\r
751 Double_t momNeg[3]={0.};\r
752 Double_t momPosAtV0vtx[3]={0.};\r
753 Double_t momNegAtV0vtx[3]={0.};\r
754 \r
755 for (Int_t nV0 = 0; nV0 < esd.GetNumberOfV0s(); ++nV0) \r
756 {\r
757 if (fUsedV0[nV0]) continue; // skip if already added to the AOD\r
758 \r
759 AliESDv0 *v0 = esd.GetV0(nV0);\r
760 Int_t posFromV0 = v0->GetPindex();\r
761 Int_t negFromV0 = v0->GetNindex();\r
762 \r
763 // V0 selection \r
f7ed7e88 764 //\r
c7ef1c47 765 AliESDVertex *esdVtx = new AliESDVertex(*(esd.GetPrimaryVertex()));\r
766 AliESDtrack *esdV0Pos = esd.GetTrack(posFromV0);\r
767 AliESDtrack *esdV0Neg = esd.GetTrack(negFromV0);\r
768 TList v0objects;\r
769 v0objects.AddAt(v0, 0);\r
770 v0objects.AddAt(esdV0Pos, 1);\r
771 v0objects.AddAt(esdV0Neg, 2);\r
772 v0objects.AddAt(esdVtx, 3);\r
773 UInt_t selectV0 = 0;\r
774 if (fV0Filter) {\r
775 selectV0 = fV0Filter->IsSelected(&v0objects);\r
776 // this is a little awkward but otherwise the \r
777 // list wants to access the pointer (delete it) \r
778 // again when going out of scope\r
779 delete v0objects.RemoveAt(3); // esdVtx created via copy construct\r
780 esdVtx = 0;\r
781 if (!selectV0) \r
782 continue;\r
783 }\r
784 else{\r
785 delete v0objects.RemoveAt(3); // esdVtx created via copy construct\r
786 esdVtx = 0;\r
787 }\r
788 \r
789 v0->GetXYZ(pos[0], pos[1], pos[2]);\r
f7ed7e88 790 \r
c7ef1c47 791 if (!fOldESDformat) {\r
f7ed7e88 792 chi2 = v0->GetChi2V0(); // = chi2/NDF since NDF = 2*2-3\r
793 v0->GetPosCov(covVtx);\r
c7ef1c47 794 } else {\r
f7ed7e88 795 chi2 = -999.;\r
796 for (Int_t i = 0; i < 6; i++) covVtx[i] = 0.;\r
c7ef1c47 797 }\r
798 \r
799 \r
800 AliAODVertex * vV0 = \r
801 new(Vertices()[fNumberOfVertices++]) AliAODVertex(pos,\r
802 covVtx,\r
803 chi2,\r
804 fPrimaryVertex,\r
805 nV0,\r
806 AliAODVertex::kV0);\r
807 fPrimaryVertex->AddDaughter(vV0);\r
808 \r
809 \r
810 // Add the positive tracks from the V0\r
811 \r
812\r
813 esdV0Pos->GetPxPyPz(momPos);\r
814 esdV0Pos->GetXYZ(pos);\r
815 esdV0Pos->GetCovarianceXYZPxPyPz(covTr);\r
816 esdV0Pos->GetESDpid(pid);\r
817 \r
818 const AliESDVertex *vtx = esd.GetPrimaryVertex();\r
819\r
820 if (!fUsedTrack[posFromV0]) {\r
821 fUsedTrack[posFromV0] = kTRUE;\r
f7ed7e88 822 UInt_t selectInfo = 0;\r
9eeae5d5 823 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdV0Pos);\r
c7ef1c47 824 if(fMChandler)fMChandler->SelectParticle(esdV0Pos->GetLabel());\r
825 aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdV0Pos->GetID(),\r
826 esdV0Pos->GetLabel(), \r
827 momPos, \r
828 kTRUE,\r
829 pos,\r
830 kFALSE,\r
831 covTr, \r
832 (Short_t)esdV0Pos->GetSign(),\r
833 esdV0Pos->GetITSClusterMap(), \r
834 pid,\r
835 vV0,\r
836 kTRUE, // check if this is right\r
837 vtx->UsesTrack(esdV0Pos->GetID()),\r
838 AliAODTrack::kSecondary,\r
839 selectInfo);\r
f2f21b24 840 aodTrack->SetTPCClusterMap(esdV0Pos->GetTPCClusterMap());\r
841 aodTrack->SetTPCSharedMap (esdV0Pos->GetTPCSharedMap());\r
40a060e5 842 aodTrack->SetChi2perNDF(Chi2perNDF(esdV0Pos));\r
3c01c166 843 aodTrack->SetTPCPointsF(esdV0Pos->GetTPCNclsF());\r
c7ef1c47 844 fAODTrackRefs->AddAt(aodTrack,posFromV0);\r
f7ed7e88 845 // if (fDebug > 0) printf("-------------------Bo: pos track from original pt %.3f \n",aodTrack->Pt());\r
c7ef1c47 846 if (esdV0Pos->GetSign() > 0) ++fNumberOfPositiveTracks;\r
f7ed7e88 847 aodTrack->ConvertAliPIDtoAODPID();\r
848 aodTrack->SetFlags(esdV0Pos->GetStatus());\r
1a3579a6 849 SetAODPID(esdV0Pos,aodTrack,detpid);\r
c7ef1c47 850 }\r
851 else {\r
852 aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(posFromV0));\r
f7ed7e88 853 // if (fDebug > 0) printf("-------------------Bo pos track from refArray pt %.3f \n",aodTrack->Pt());\r
c7ef1c47 854 }\r
855 vV0->AddDaughter(aodTrack);\r
856 \r
857 // Add the negative tracks from the V0\r
858 \r
859 esdV0Neg->GetPxPyPz(momNeg);\r
860 esdV0Neg->GetXYZ(pos);\r
861 esdV0Neg->GetCovarianceXYZPxPyPz(covTr);\r
862 esdV0Neg->GetESDpid(pid);\r
863 \r
864 if (!fUsedTrack[negFromV0]) {\r
865 fUsedTrack[negFromV0] = kTRUE;\r
f7ed7e88 866 UInt_t selectInfo = 0;\r
9eeae5d5 867 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdV0Neg);\r
c7ef1c47 868 if(fMChandler)fMChandler->SelectParticle(esdV0Neg->GetLabel());\r
869 aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdV0Neg->GetID(),\r
870 esdV0Neg->GetLabel(),\r
871 momNeg,\r
872 kTRUE,\r
873 pos,\r
874 kFALSE,\r
875 covTr, \r
876 (Short_t)esdV0Neg->GetSign(),\r
877 esdV0Neg->GetITSClusterMap(), \r
878 pid,\r
879 vV0,\r
880 kTRUE, // check if this is right\r
881 vtx->UsesTrack(esdV0Neg->GetID()),\r
882 AliAODTrack::kSecondary,\r
883 selectInfo);\r
f2f21b24 884 aodTrack->SetTPCClusterMap(esdV0Neg->GetTPCClusterMap());\r
885 aodTrack->SetTPCSharedMap (esdV0Neg->GetTPCSharedMap());\r
40a060e5 886 aodTrack->SetChi2perNDF(Chi2perNDF(esdV0Neg));\r
3c01c166 887 aodTrack->SetTPCPointsF(esdV0Neg->GetTPCNclsF());\r
f7ed7e88 888 \r
c7ef1c47 889 fAODTrackRefs->AddAt(aodTrack,negFromV0);\r
f7ed7e88 890 // if (fDebug > 0) printf("-------------------Bo: neg track from original pt %.3f \n",aodTrack->Pt());\r
c7ef1c47 891 if (esdV0Neg->GetSign() > 0) ++fNumberOfPositiveTracks;\r
f7ed7e88 892 aodTrack->ConvertAliPIDtoAODPID();\r
893 aodTrack->SetFlags(esdV0Neg->GetStatus());\r
1a3579a6 894 SetAODPID(esdV0Neg,aodTrack,detpid);\r
c7ef1c47 895 }\r
896 else {\r
897 aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(negFromV0));\r
f7ed7e88 898 // if (fDebug > 0) printf("-------------------Bo neg track from refArray pt %.3f \n",aodTrack->Pt());\r
f7ed7e88 899 }\r
c7ef1c47 900 vV0->AddDaughter(aodTrack);\r
901 \r
902 \r
903 // Add the V0 the V0 array as well\r
904 \r
905 Double_t dcaV0Daughters = v0->GetDcaV0Daughters();\r
906 Double_t dcaV0ToPrimVertex = v0->GetD(esd.GetPrimaryVertex()->GetX(),\r
907 esd.GetPrimaryVertex()->GetY(),\r
908 esd.GetPrimaryVertex()->GetZ());\r
909 \r
910 v0->GetPPxPyPz(momPosAtV0vtx[0],momPosAtV0vtx[1],momPosAtV0vtx[2]); \r
911 v0->GetNPxPyPz(momNegAtV0vtx[0],momNegAtV0vtx[1],momNegAtV0vtx[2]); \r
912 \r
913 Double_t dcaDaughterToPrimVertex[2] = { 999., 999.}; // ..[0] = DCA in (x,y) for Pos and ..[1] = Neg\r
914 dcaDaughterToPrimVertex[0] = TMath::Abs(esdV0Pos->GetD( esd.GetPrimaryVertex()->GetX(),\r
915 esd.GetPrimaryVertex()->GetY(),\r
916 esd.GetMagneticField()) );\r
917 dcaDaughterToPrimVertex[1] = TMath::Abs(esdV0Neg->GetD( esd.GetPrimaryVertex()->GetX(),\r
918 esd.GetPrimaryVertex()->GetY(),\r
919 esd.GetMagneticField()) );\r
f7ed7e88 920 \r
c7ef1c47 921 AliAODv0* aodV0 = new(V0s()[fNumberOfV0s++]) AliAODv0(vV0, \r
922 dcaV0Daughters,\r
923 dcaV0ToPrimVertex,\r
924 momPosAtV0vtx,\r
925 momNegAtV0vtx,\r
926 dcaDaughterToPrimVertex);\r
927 \r
928 // set the aod v0 on-the-fly status\r
929 aodV0->SetOnFlyStatus(v0->GetOnFlyStatus());\r
930 }//End of loop on V0s \r
f7ed7e88 931 \r
c7ef1c47 932 V0s().Expand(fNumberOfV0s); \r
933}\r
f7ed7e88 934\r
c7ef1c47 935//______________________________________________________________________________\r
936void AliAnalysisTaskESDfilter::ConvertTPCOnlyTracks(const AliESDEvent& esd)\r
937{\r
938 // Convert TPC only tracks\r
939 \r
940 AliCodeTimerAuto("",0);\r
941 \r
942 // Loop over the tracks and extract and mask out all aod tracks that pass the selections for AODt racks\r
943 for(int it = 0;it < fNumberOfTracks;++it)\r
944 {\r
945 AliAODTrack *tr = (AliAODTrack*)(Tracks().At(it));\r
946 if(!tr)continue;\r
947 UInt_t map = tr->GetFilterMap();\r
948 if(map&fTPCOnlyFilterMask){\r
949 // we only reset the track select ionfo, no deletion...\r
950 tr->SetFilterMap(map&~fTPCOnlyFilterMask);\r
951 }\r
952 }\r
953 // Loop over the ESD trcks and pick out the tracks passing TPC only cuts\r
954 \r
955 \r
956 const AliESDVertex *vtxSPD = esd.GetPrimaryVertexSPD();\r
957 const AliESDVertex *vtx = esd.GetPrimaryVertex();\r
958\r
959 Double_t pos[3] = { 0. }; \r
960 Double_t covTr[21]={0.};\r
961 Double_t pid[10]={0.}; \r
962 Double_t p[3] = { 0. };\r
4b6f9f24 963\r
964 Double_t pDCA[3] = { 0. }; // momentum at DCA\r
965 Double_t rDCA[3] = { 0. }; // position at DCA\r
966 Float_t dDCA[2] = {0.}; // DCA to the vertex d and z\r
967 Float_t cDCA[3] = {0.}; // covariance of impact parameters\r
968\r
969\r
c7ef1c47 970 AliAODTrack* aodTrack(0x0);\r
971 \r
972 for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack) \r
973 {\r
974 AliESDtrack* esdTrack = esd.GetTrack(nTrack); //carefull do not modify it othwise need to work with a copy \r
975 \r
976 UInt_t selectInfo = 0;\r
977 //\r
978 // Track selection\r
979 if (fTrackFilter) {\r
980 selectInfo = fTrackFilter->IsSelected(esdTrack);\r
981 }\r
982 selectInfo &= fTPCOnlyFilterMask;\r
983 if (!selectInfo)continue;\r
f7ed7e88 984 \r
c7ef1c47 985 // create a tpc only tracl\r
986 AliESDtrack *track = AliESDtrackCuts::GetTPCOnlyTrack(const_cast<AliESDEvent*>(&esd),esdTrack->GetID());\r
987 if(!track) continue;\r
988 \r
989 if(track->Pt()>0.)\r
990 {\r
991 // only constrain tracks above threshold\r
992 AliExternalTrackParam exParam;\r
993 // take the B-field from the ESD, no 3D fieldMap available at this point\r
994 Bool_t relate = false;\r
2fb345ad 995 relate = track->RelateToVertexTPC(vtxSPD,esd.GetMagneticField(),kVeryBig,&exParam);\r
c7ef1c47 996 if(!relate){\r
997 delete track;\r
998 continue;\r
a7b1466d 999 }\r
4b6f9f24 1000 // fetch the track parameters at the DCA (unconstraint)\r
1001 if(track->GetTPCInnerParam()){\r
1002 track->GetTPCInnerParam()->GetPxPyPz(pDCA);\r
1003 track->GetTPCInnerParam()->GetXYZ(rDCA);\r
1004 }\r
1005 // get the DCA to the vertex:\r
1006 track->GetImpactParametersTPC(dDCA,cDCA);\r
1007 // set the constraint parameters to the track\r
c7ef1c47 1008 track->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance());\r
1009 }\r
1010 \r
1011 track->GetPxPyPz(p);\r
1012 track->GetXYZ(pos);\r
1013 track->GetCovarianceXYZPxPyPz(covTr);\r
1014 track->GetESDpid(pid);\r
1015 if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());\r
1016 aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack((track->GetID()+1)*-1,\r
1017 track->GetLabel(),\r
1018 p,\r
1019 kTRUE,\r
1020 pos,\r
1021 kFALSE,\r
1022 covTr, \r
1023 (Short_t)track->GetSign(),\r
1024 track->GetITSClusterMap(), \r
1025 pid,\r
1026 fPrimaryVertex,\r
1027 kTRUE, // check if this is right\r
1028 vtx->UsesTrack(track->GetID()),\r
1029 AliAODTrack::kPrimary, \r
1030 selectInfo);\r
1031 aodTrack->SetTPCClusterMap(track->GetTPCClusterMap());\r
1032 aodTrack->SetTPCSharedMap (track->GetTPCSharedMap());\r
2fb345ad 1033 Float_t ndf = track->GetTPCNcls()+1 - 5 ;\r
1034 if(ndf>0){\r
af769c14 1035 aodTrack->SetChi2perNDF(track->GetConstrainedChi2TPC());\r
2fb345ad 1036 }\r
1037 else{\r
1038 aodTrack->SetChi2perNDF(-1);\r
1039 }\r
4b6f9f24 1040\r
1041 // set the DCA values to the AOD track\r
1042 aodTrack->SetPxPyPzAtDCA(pDCA[0],pDCA[1],pDCA[2]);\r
1043 aodTrack->SetXYAtDCA(rDCA[0],rDCA[1]);\r
1044 aodTrack->SetDCA(dDCA[0],dDCA[1]);\r
1045\r
c7ef1c47 1046 aodTrack->SetFlags(track->GetStatus());\r
3c01c166 1047 aodTrack->SetTPCPointsF(track->GetTPCNclsF());\r
1048\r
c7ef1c47 1049 delete track;\r
1050 } // end of loop on tracks\r
1051 \r
1052}\r
1053\r
1054//______________________________________________________________________________\r
1055void AliAnalysisTaskESDfilter::ConvertTracks(const AliESDEvent& esd)\r
1056{\r
1057 // Tracks (primary and orphan)\r
1058\r
1059 AliCodeTimerAuto("",0);\r
1060 \r
1061 AliDebug(1,Form("NUMBER OF ESD TRACKS %5d\n", esd.GetNumberOfTracks()));\r
1062 \r
1063 const AliESDVertex *vtx = esd.GetPrimaryVertex();\r
1064 Double_t p[3] = { 0. };\r
1065 Double_t pos[3] = { 0. };\r
1066 Double_t covTr[21] = { 0. };\r
1067 Double_t pid[10] = { 0. };\r
1068 AliAODTrack* aodTrack(0x0);\r
1069 AliAODPid* detpid(0x0);\r
1070 \r
1071 for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack) \r
1072 {\r
1073 if (fUsedTrack[nTrack]) continue;\r
1074 \r
1075 AliESDtrack *esdTrack = esd.GetTrack(nTrack);\r
1076 UInt_t selectInfo = 0;\r
1077 //\r
1078 // Track selection\r
1079 if (fTrackFilter) {\r
9eeae5d5 1080 selectInfo = fTrackFilter->IsSelected(esdTrack);\r
f7ed7e88 1081 if (!selectInfo && !vtx->UsesTrack(esdTrack->GetID())) continue;\r
f7ed7e88 1082 }\r
c7ef1c47 1083 \r
1084 \r
1085 esdTrack->GetPxPyPz(p);\r
1086 esdTrack->GetXYZ(pos);\r
1087 esdTrack->GetCovarianceXYZPxPyPz(covTr);\r
1088 esdTrack->GetESDpid(pid);\r
1089 if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());\r
1090 fPrimaryVertex->AddDaughter(aodTrack =\r
1091 new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrack->GetID(),\r
1092 esdTrack->GetLabel(),\r
1093 p,\r
1094 kTRUE,\r
1095 pos,\r
1096 kFALSE,\r
1097 covTr, \r
1098 (Short_t)esdTrack->GetSign(),\r
1099 esdTrack->GetITSClusterMap(), \r
1100 pid,\r
1101 fPrimaryVertex,\r
1102 kTRUE, // check if this is right\r
1103 vtx->UsesTrack(esdTrack->GetID()),\r
1104 AliAODTrack::kPrimary, \r
1105 selectInfo)\r
1106 );\r
1107 aodTrack->SetTPCClusterMap(esdTrack->GetTPCClusterMap());\r
1108 aodTrack->SetTPCSharedMap (esdTrack->GetTPCSharedMap());\r
1109 aodTrack->SetChi2perNDF(Chi2perNDF(esdTrack));\r
3c01c166 1110 aodTrack->SetTPCPointsF(esdTrack->GetTPCNclsF());\r
1111\r
c7ef1c47 1112 fAODTrackRefs->AddAt(aodTrack, nTrack);\r
1113 \r
1114 \r
1115 if (esdTrack->GetSign() > 0) ++fNumberOfPositiveTracks;\r
1116 aodTrack->SetFlags(esdTrack->GetStatus());\r
1117 aodTrack->ConvertAliPIDtoAODPID();\r
1a3579a6 1118 SetAODPID(esdTrack,aodTrack,detpid);\r
c7ef1c47 1119 } // end of loop on tracks\r
1120}\r
f7ed7e88 1121\r
c7ef1c47 1122//______________________________________________________________________________\r
1123void AliAnalysisTaskESDfilter::ConvertPmdClusters(const AliESDEvent& esd)\r
1124{\r
1125 AliCodeTimerAuto("",0);\r
1126 Int_t jPmdClusters=0;\r
1127 // Access to the AOD container of PMD clusters\r
1128 TClonesArray &pmdClusters = *(AODEvent()->GetPmdClusters());\r
1129 for (Int_t iPmd = 0; iPmd < esd.GetNumberOfPmdTracks(); ++iPmd) {\r
1130 // file pmd clusters, to be revised!\r
1131 AliESDPmdTrack *pmdTrack = esd.GetPmdTrack(iPmd);\r
1132 Int_t nLabel = 0;\r
1133 Int_t *label = 0x0;\r
1134 Double_t posPmd[3] = { pmdTrack->GetClusterX(), pmdTrack->GetClusterY(), pmdTrack->GetClusterZ()};\r
1135 Double_t pidPmd[13] = { 0.}; // to be revised!\r
1136 // type not set!\r
1137 // assoc cluster not set\r
1138 new(pmdClusters[jPmdClusters++]) AliAODPmdCluster(iPmd, nLabel, label, pmdTrack->GetClusterADC(), posPmd, pidPmd);\r
1139 }\r
1140}\r
f7ed7e88 1141\r
c7ef1c47 1142//______________________________________________________________________________\r
1143void AliAnalysisTaskESDfilter::ConvertCaloClusters(const AliESDEvent& esd)\r
1144{\r
1145 AliCodeTimerAuto("",0);\r
1146 \r
1147 // Access to the AOD container of clusters\r
1148 TClonesArray &caloClusters = *(AODEvent()->GetCaloClusters());\r
1149 Int_t jClusters(0);\r
1150 \r
1151 for (Int_t iClust=0; iClust<esd.GetNumberOfCaloClusters(); ++iClust) {\r
1152 \r
1153 AliESDCaloCluster * cluster = esd.GetCaloCluster(iClust);\r
1154 \r
1155 Int_t id = cluster->GetID();\r
1156 Int_t nLabel = cluster->GetNLabels();\r
1157 Int_t *labels = cluster->GetLabels();\r
1158 if(labels){ \r
c8fe2783 1159 for(int i = 0;i < nLabel;++i){\r
c7ef1c47 1160 if(fMChandler)fMChandler->SelectParticle(labels[i]);\r
c8fe2783 1161 }\r
1162 } \r
c7ef1c47 1163 \r
1164 Float_t energy = cluster->E();\r
1165 Float_t posF[3] = { 0.};\r
1166 cluster->GetPosition(posF);\r
1167 \r
1168 AliAODCaloCluster *caloCluster = new(caloClusters[jClusters++]) AliAODCaloCluster(id,\r
1169 nLabel,\r
1170 labels,\r
1171 energy,\r
1172 posF,\r
1173 NULL,\r
1174 cluster->GetType(),0);\r
1175 \r
1176 caloCluster->SetCaloCluster(cluster->GetDistanceToBadChannel(),\r
1177 cluster->GetDispersion(),\r
1178 cluster->GetM20(), cluster->GetM02(),\r
1179 cluster->GetEmcCpvDistance(), \r
1180 cluster->GetNExMax(),cluster->GetTOF()) ;\r
1181 \r
1182 caloCluster->SetPIDFromESD(cluster->GetPID());\r
1183 caloCluster->SetNCells(cluster->GetNCells());\r
1184 caloCluster->SetCellsAbsId(cluster->GetCellsAbsId());\r
1185 caloCluster->SetCellsAmplitudeFraction(cluster->GetCellsAmplitudeFraction());\r
1186 \r
1187 TArrayI* matchedT = cluster->GetTracksMatched();\r
1188 if (fNumberOfTracks>0 && matchedT && cluster->GetTrackMatchedIndex() >= 0) { \r
1189 for (Int_t im = 0; im < matchedT->GetSize(); im++) {\r
1190 Int_t iESDtrack = matchedT->At(im);;\r
1191 if (fAODTrackRefs->At(iESDtrack) != 0) {\r
1192 caloCluster->AddTrackMatched((AliAODTrack*)fAODTrackRefs->At(iESDtrack));\r
1193 }\r
1194 }\r
1195 }\r
1196 \r
1197 } \r
1198 caloClusters.Expand(jClusters); // resize TObjArray to 'remove' slots for pseudo clusters \r
1199}\r
f7ed7e88 1200\r
c7ef1c47 1201//______________________________________________________________________________\r
1202void AliAnalysisTaskESDfilter::ConvertEMCALCells(const AliESDEvent& esd)\r
1203{\r
1204 AliCodeTimerAuto("",0);\r
1205 // fill EMCAL cell info\r
1206 if (esd.GetEMCALCells()) { // protection against missing ESD information\r
1207 AliESDCaloCells &esdEMcells = *(esd.GetEMCALCells());\r
1208 Int_t nEMcell = esdEMcells.GetNumberOfCells() ;\r
1209 \r
1210 AliAODCaloCells &aodEMcells = *(AODEvent()->GetEMCALCells());\r
1211 aodEMcells.CreateContainer(nEMcell);\r
1212 aodEMcells.SetType(AliAODCaloCells::kEMCALCell);\r
1213 for (Int_t iCell = 0; iCell < nEMcell; iCell++) { \r
1214 aodEMcells.SetCell(iCell,esdEMcells.GetCellNumber(iCell),esdEMcells.GetAmplitude(iCell));\r
1215 }\r
1216 aodEMcells.Sort();\r
1217 }\r
1218}\r
f7ed7e88 1219\r
c7ef1c47 1220//______________________________________________________________________________\r
1221void AliAnalysisTaskESDfilter::ConvertPHOSCells(const AliESDEvent& esd)\r
1222{\r
1223 AliCodeTimerAuto("",0);\r
1224 // fill PHOS cell info\r
1225 if (esd.GetPHOSCells()) { // protection against missing ESD information\r
1226 AliESDCaloCells &esdPHcells = *(esd.GetPHOSCells());\r
1227 Int_t nPHcell = esdPHcells.GetNumberOfCells() ;\r
1228 \r
1229 AliAODCaloCells &aodPHcells = *(AODEvent()->GetPHOSCells());\r
1230 aodPHcells.CreateContainer(nPHcell);\r
1231 aodPHcells.SetType(AliAODCaloCells::kPHOSCell);\r
1232 for (Int_t iCell = 0; iCell < nPHcell; iCell++) { \r
1233 aodPHcells.SetCell(iCell,esdPHcells.GetCellNumber(iCell),esdPHcells.GetAmplitude(iCell));\r
f7ed7e88 1234 }\r
c7ef1c47 1235 aodPHcells.Sort();\r
1236 }\r
1237}\r
f7ed7e88 1238\r
c7ef1c47 1239//______________________________________________________________________________\r
1240void AliAnalysisTaskESDfilter::ConvertTracklets(const AliESDEvent& esd)\r
1241{\r
1242 // tracklets \r
1243 AliCodeTimerAuto("",0);\r
1244\r
1245 AliAODTracklets &SPDTracklets = *(AODEvent()->GetTracklets());\r
1246 const AliMultiplicity *mult = esd.GetMultiplicity();\r
1247 if (mult) {\r
1248 if (mult->GetNumberOfTracklets()>0) {\r
1249 SPDTracklets.CreateContainer(mult->GetNumberOfTracklets());\r
f7ed7e88 1250 \r
c7ef1c47 1251 for (Int_t n=0; n<mult->GetNumberOfTracklets(); n++) {\r
1252 if(fMChandler){\r
1253 fMChandler->SelectParticle(mult->GetLabel(n, 0));\r
1254 fMChandler->SelectParticle(mult->GetLabel(n, 1));\r
1255 }\r
1256 SPDTracklets.SetTracklet(n, mult->GetTheta(n), mult->GetPhi(n), mult->GetDeltaPhi(n), mult->GetLabel(n, 0),mult->GetLabel(n, 1));\r
f7ed7e88 1257 }\r
f7ed7e88 1258 }\r
c7ef1c47 1259 } else {\r
1260 //Printf("ERROR: AliMultiplicity could not be retrieved from ESD");\r
1261 }\r
1262}\r
f7ed7e88 1263\r
c7ef1c47 1264//______________________________________________________________________________\r
1265void AliAnalysisTaskESDfilter::ConvertKinks(const AliESDEvent& esd)\r
1266{\r
1267 AliCodeTimerAuto("",0);\r
1268 \r
1269 // Kinks: it is a big mess the access to the information in the kinks\r
1270 // The loop is on the tracks in order to find the mother and daugther of each kink\r
1271 \r
1272 Double_t covTr[21]={0.};\r
1273 Double_t pid[10]={0.};\r
1274 AliAODPid* detpid(0x0);\r
1275 \r
1276 fNumberOfKinks = esd.GetNumberOfKinks();\r
1277\r
1278 const AliESDVertex* vtx = esd.GetPrimaryVertex();\r
1279 \r
1280 for (Int_t iTrack=0; iTrack<esd.GetNumberOfTracks(); ++iTrack) \r
1281 {\r
1282 AliESDtrack * esdTrack = esd.GetTrack(iTrack);\r
1283 \r
1284 Int_t ikink = esdTrack->GetKinkIndex(0);\r
1285 \r
1286 if (ikink && fNumberOfKinks) {\r
1287 // Negative kink index: mother, positive: daughter\r
1288 \r
1289 // Search for the second track of the kink\r
1290 \r
1291 for (Int_t jTrack = iTrack+1; jTrack<esd.GetNumberOfTracks(); ++jTrack) {\r
1292 \r
1293 AliESDtrack * esdTrack1 = esd.GetTrack(jTrack);\r
1294 \r
1295 Int_t jkink = esdTrack1->GetKinkIndex(0);\r
1296 \r
1297 if ( TMath::Abs(ikink)==TMath::Abs(jkink) ) {\r
1298 \r
1299 // The two tracks are from the same kink\r
1300 \r
1301 if (fUsedKink[TMath::Abs(ikink)-1]) continue; // skip used kinks\r
1302 \r
1303 Int_t imother = -1;\r
1304 Int_t idaughter = -1;\r
1305 \r
1306 if (ikink<0 && jkink>0) {\r
1307 \r
1308 imother = iTrack;\r
1309 idaughter = jTrack;\r
1310 }\r
1311 else if (ikink>0 && jkink<0) {\r
1312 \r
1313 imother = jTrack;\r
1314 idaughter = iTrack;\r
1315 }\r
1316 else {\r
1317 // cerr << "Error: Wrong combination of kink indexes: "\r
1318 // << ikink << " " << jkink << endl;\r
1319 continue;\r
1320 }\r
1321 \r
1322 // Add the mother track if it passed primary track selection cuts\r
1323 \r
1324 AliAODTrack * mother = NULL;\r
1325 \r
1326 UInt_t selectInfo = 0;\r
1327 if (fTrackFilter) {\r
1328 selectInfo = fTrackFilter->IsSelected(esd.GetTrack(imother));\r
1329 if (!selectInfo) continue;\r
1330 }\r
1331 \r
1332 if (!fUsedTrack[imother]) {\r
1333 \r
1334 fUsedTrack[imother] = kTRUE;\r
1335 \r
1336 AliESDtrack *esdTrackM = esd.GetTrack(imother);\r
1337 Double_t p[3] = { 0. };\r
1338 Double_t pos[3] = { 0. };\r
1339 esdTrackM->GetPxPyPz(p);\r
1340 esdTrackM->GetXYZ(pos);\r
1341 esdTrackM->GetCovarianceXYZPxPyPz(covTr);\r
1342 esdTrackM->GetESDpid(pid);\r
1343 if(fMChandler)fMChandler->SelectParticle(esdTrackM->GetLabel());\r
1344 mother = \r
1345 new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrackM->GetID(),\r
1346 esdTrackM->GetLabel(),\r
1347 p,\r
1348 kTRUE,\r
1349 pos,\r
1350 kFALSE,\r
1351 covTr, \r
1352 (Short_t)esdTrackM->GetSign(),\r
1353 esdTrackM->GetITSClusterMap(), \r
1354 pid,\r
1355 fPrimaryVertex,\r
1356 kTRUE, // check if this is right\r
1357 vtx->UsesTrack(esdTrack->GetID()),\r
1358 AliAODTrack::kPrimary,\r
1359 selectInfo);\r
1360 mother->SetTPCClusterMap(esdTrackM->GetTPCClusterMap());\r
1361 mother->SetTPCSharedMap (esdTrackM->GetTPCSharedMap());\r
1362 mother->SetChi2perNDF(Chi2perNDF(esdTrackM));\r
3c01c166 1363 mother->SetTPCPointsF(esdTrackM->GetTPCNclsF());\r
1364\r
c7ef1c47 1365 fAODTrackRefs->AddAt(mother, imother);\r
1366 \r
1367 if (esdTrackM->GetSign() > 0) ++fNumberOfPositiveTracks;\r
1368 mother->SetFlags(esdTrackM->GetStatus());\r
1369 mother->ConvertAliPIDtoAODPID();\r
1370 fPrimaryVertex->AddDaughter(mother);\r
1371 mother->ConvertAliPIDtoAODPID();\r
1a3579a6 1372 SetAODPID(esdTrackM,mother,detpid);\r
c7ef1c47 1373 }\r
1374 else {\r
1375 // cerr << "Error: event " << esd.GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1\r
1376 // << " track " << imother << " has already been used!" << endl;\r
1377 }\r
1378 \r
1379 // Add the kink vertex\r
1380 AliESDkink * kink = esd.GetKink(TMath::Abs(ikink)-1);\r
1381 \r
1382 AliAODVertex * vkink = \r
1383 new(Vertices()[fNumberOfVertices++]) AliAODVertex(kink->GetPosition(),\r
1384 NULL,\r
1385 0.,\r
1386 mother,\r
1387 esdTrack->GetID(), // This is the track ID of the mother's track!\r
1388 AliAODVertex::kKink);\r
1389 // Add the daughter track\r
1390 \r
1391 AliAODTrack * daughter = NULL;\r
1392 \r
1393 if (!fUsedTrack[idaughter]) {\r
1394 \r
1395 fUsedTrack[idaughter] = kTRUE;\r
1396 \r
1397 AliESDtrack *esdTrackD = esd.GetTrack(idaughter);\r
1398 Double_t p[3] = { 0. };\r
1399 Double_t pos[3] = { 0. };\r
1400\r
1401 esdTrackD->GetPxPyPz(p);\r
1402 esdTrackD->GetXYZ(pos);\r
1403 esdTrackD->GetCovarianceXYZPxPyPz(covTr);\r
1404 esdTrackD->GetESDpid(pid);\r
1405 selectInfo = 0;\r
1406 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdTrackD);\r
1407 if(fMChandler)fMChandler->SelectParticle(esdTrackD->GetLabel());\r
1408 daughter = \r
1409 new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrackD->GetID(),\r
1410 esdTrackD->GetLabel(),\r
1411 p,\r
1412 kTRUE,\r
1413 pos,\r
1414 kFALSE,\r
1415 covTr, \r
1416 (Short_t)esdTrackD->GetSign(),\r
1417 esdTrackD->GetITSClusterMap(), \r
1418 pid,\r
1419 vkink,\r
1420 kTRUE, // check if this is right\r
1421 vtx->UsesTrack(esdTrack->GetID()),\r
1422 AliAODTrack::kSecondary,\r
1423 selectInfo);\r
1424 daughter->SetTPCClusterMap(esdTrackD->GetTPCClusterMap());\r
1425 daughter->SetTPCSharedMap (esdTrackD->GetTPCSharedMap());\r
3c01c166 1426 daughter->SetTPCPointsF(esdTrackD->GetTPCNclsF());\r
c7ef1c47 1427 fAODTrackRefs->AddAt(daughter, idaughter);\r
1428 \r
1429 if (esdTrackD->GetSign() > 0) ++fNumberOfPositiveTracks;\r
1430 daughter->SetFlags(esdTrackD->GetStatus());\r
1431 daughter->ConvertAliPIDtoAODPID();\r
1432 vkink->AddDaughter(daughter);\r
1433 daughter->ConvertAliPIDtoAODPID();\r
1a3579a6 1434 SetAODPID(esdTrackD,daughter,detpid);\r
c7ef1c47 1435 }\r
1436 else {\r
1437 // cerr << "Error: event " << esd.GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1\r
1438 // << " track " << idaughter << " has already been used!" << endl;\r
1439 }\r
1440 }\r
1441 }\r
1442 } \r
1443 }\r
1444}\r
1445\r
1446//______________________________________________________________________________\r
1447void AliAnalysisTaskESDfilter::ConvertPrimaryVertices(const AliESDEvent& esd)\r
1448{\r
1449 AliCodeTimerAuto("",0);\r
1450 \r
1451 // Access to the AOD container of vertices\r
1452 fNumberOfVertices = 0;\r
1453 \r
1454 Double_t pos[3] = { 0. };\r
1455 Double_t covVtx[6] = { 0. };\r
f7ed7e88 1456\r
c7ef1c47 1457 // Add primary vertex. The primary tracks will be defined\r
1458 // after the loops on the composite objects (V0, cascades, kinks)\r
1459 const AliESDVertex *vtx = esd.GetPrimaryVertex();\r
1460 \r
1461 vtx->GetXYZ(pos); // position\r
1462 vtx->GetCovMatrix(covVtx); //covariance matrix\r
1463 \r
1464 fPrimaryVertex = new(Vertices()[fNumberOfVertices++])\r
1465 AliAODVertex(pos, covVtx, vtx->GetChi2toNDF(), NULL, -1, AliAODVertex::kPrimary);\r
1466 fPrimaryVertex->SetName(vtx->GetName());\r
1467 fPrimaryVertex->SetTitle(vtx->GetTitle());\r
1468 \r
1469 TString vtitle = vtx->GetTitle();\r
1470 if (!vtitle.Contains("VertexerTracks")) \r
1471 fPrimaryVertex->SetNContributors(vtx->GetNContributors());\r
1472 \r
1473 if (fDebug > 0) fPrimaryVertex->Print(); \r
1474 \r
1475 // Add SPD "main" vertex \r
1476 const AliESDVertex *vtxS = esd.GetPrimaryVertexSPD();\r
1477 vtxS->GetXYZ(pos); // position\r
1478 vtxS->GetCovMatrix(covVtx); //covariance matrix\r
1479 AliAODVertex * mVSPD = new(Vertices()[fNumberOfVertices++])\r
1480 AliAODVertex(pos, covVtx, vtxS->GetChi2toNDF(), NULL, -1, AliAODVertex::kMainSPD);\r
1481 mVSPD->SetName(vtxS->GetName());\r
1482 mVSPD->SetTitle(vtxS->GetTitle());\r
1483 mVSPD->SetNContributors(vtxS->GetNContributors()); \r
1484 \r
1485 // Add SPD pileup vertices\r
1486 for(Int_t iV=0; iV<esd.GetNumberOfPileupVerticesSPD(); ++iV)\r
1487 {\r
1488 const AliESDVertex *vtxP = esd.GetPileupVertexSPD(iV);\r
1489 vtxP->GetXYZ(pos); // position\r
1490 vtxP->GetCovMatrix(covVtx); //covariance matrix\r
1491 AliAODVertex * pVSPD = new(Vertices()[fNumberOfVertices++])\r
1492 AliAODVertex(pos, covVtx, vtxP->GetChi2toNDF(), NULL, -1, AliAODVertex::kPileupSPD);\r
1493 pVSPD->SetName(vtxP->GetName());\r
1494 pVSPD->SetTitle(vtxP->GetTitle());\r
1495 pVSPD->SetNContributors(vtxP->GetNContributors()); \r
1496 }\r
1497 \r
1498 // Add TRK pileup vertices\r
1499 for(Int_t iV=0; iV<esd.GetNumberOfPileupVerticesTracks(); ++iV)\r
1500 {\r
1501 const AliESDVertex *vtxP = esd.GetPileupVertexTracks(iV);\r
1502 vtxP->GetXYZ(pos); // position\r
1503 vtxP->GetCovMatrix(covVtx); //covariance matrix\r
1504 AliAODVertex * pVTRK = new(Vertices()[fNumberOfVertices++])\r
1505 AliAODVertex(pos, covVtx, vtxP->GetChi2toNDF(), NULL, -1, AliAODVertex::kPileupTracks);\r
1506 pVTRK->SetName(vtxP->GetName());\r
1507 pVTRK->SetTitle(vtxP->GetTitle());\r
1508 pVTRK->SetNContributors(vtxP->GetNContributors());\r
1509 }\r
f7ed7e88 1510}\r
1511\r
c7ef1c47 1512//______________________________________________________________________________\r
1513void AliAnalysisTaskESDfilter::ConvertVZERO(const AliESDEvent& esd)\r
1514{\r
1515 // Convert VZERO data\r
1516 AliAODVZERO* vzeroData = AODEvent()->GetVZEROData();\r
1517 *vzeroData = *(esd.GetVZEROData());\r
1518}\r
f7ed7e88 1519\r
c7ef1c47 1520//______________________________________________________________________________\r
1521void AliAnalysisTaskESDfilter::ConvertESDtoAOD() \r
1522{\r
1523 // ESD Filter analysis task executed for each event\r
1524 \r
1525 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());\r
1526 \r
1527 if(!esd)return;\r
1528\r
1529 AliCodeTimerAuto("",0);\r
1530 \r
1531 fOldESDformat = ( esd->GetAliESDOld() != 0x0 );\r
1532 \r
1533 fNumberOfTracks = 0;\r
1534 fNumberOfPositiveTracks = 0;\r
1535 fNumberOfV0s = 0;\r
1536 fNumberOfVertices = 0;\r
1537 fNumberOfCascades = 0;\r
1538 fNumberOfKinks = 0;\r
1539 \r
1540 AliAODHeader* header = ConvertHeader(*esd);\r
1541\r
1542 if ( fIsVZEROEnabled ) ConvertVZERO(*esd);\r
1543 \r
1544 // Fetch Stack for debuggging if available \r
1545 fMChandler=0x0;\r
1546 if(MCEvent())\r
1547 {\r
1548 fMChandler = (AliMCEventHandler*) ((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler()); \r
1549 }\r
1550 \r
1551 // loop over events and fill them\r
1552 // Multiplicity information needed by the header (to be revised!)\r
1553 Int_t nTracks = esd->GetNumberOfTracks();\r
1554 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) esd->GetTrack(iTrack)->SetESDEvent(esd);\r
1555\r
1556 // Update the header\r
1557\r
1558 Int_t nV0s = esd->GetNumberOfV0s();\r
1559 Int_t nCascades = esd->GetNumberOfCascades();\r
1560 Int_t nKinks = esd->GetNumberOfKinks();\r
1561 Int_t nVertices = nV0s + nCascades /*V0 wihtin cascade already counted*/+ nKinks + 1 /* = prim. vtx*/;\r
1562 Int_t nPileSPDVertices=1+esd->GetNumberOfPileupVerticesSPD(); // also SPD main vertex\r
1563 Int_t nPileTrkVertices=esd->GetNumberOfPileupVerticesTracks();\r
1564 nVertices+=nPileSPDVertices;\r
1565 nVertices+=nPileTrkVertices;\r
1566 Int_t nJets = 0;\r
1567 Int_t nCaloClus = esd->GetNumberOfCaloClusters();\r
1568 Int_t nFmdClus = 0;\r
1569 Int_t nPmdClus = esd->GetNumberOfPmdTracks();\r
1570 \r
1571 AliDebug(1,Form(" NV0=%d NCASCADES=%d NKINKS=%d", nV0s, nCascades, nKinks));\r
1572 \r
1573 AODEvent()->ResetStd(nTracks, nVertices, nV0s, nCascades, nJets, nCaloClus, nFmdClus, nPmdClus);\r
1574\r
1575 if (nV0s > 0) \r
1576 {\r
1577 // RefArray to store a mapping between esd V0 number and newly created AOD-Vertex V0\r
1578 fAODV0VtxRefs = new TRefArray(nV0s);\r
1579 // RefArray to store the mapping between esd V0 number and newly created AOD-V0\r
1580 fAODV0Refs = new TRefArray(nV0s); \r
1581 // Array to take into account the V0s already added to the AOD (V0 within cascades)\r
1582 fUsedV0 = new Bool_t[nV0s];\r
1583 for (Int_t iV0=0; iV0<nV0s; ++iV0) fUsedV0[iV0]=kFALSE;\r
1584 }\r
1585 \r
1586 if (nTracks>0) \r
1587 {\r
1588 // RefArray to store the mapping between esd track number and newly created AOD-Track\r
1589 \r
1590 fAODTrackRefs = new TRefArray(nTracks);\r
1591\r
1592 // Array to take into account the tracks already added to the AOD \r
1593 fUsedTrack = new Bool_t[nTracks];\r
1594 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) fUsedTrack[iTrack]=kFALSE;\r
1595 }\r
1596 \r
1597 // Array to take into account the kinks already added to the AOD\r
1598 if (nKinks>0) \r
1599 {\r
1600 fUsedKink = new Bool_t[nKinks];\r
1601 for (Int_t iKink=0; iKink<nKinks; ++iKink) fUsedKink[iKink]=kFALSE;\r
1602 }\r
1603 \r
1604 ConvertPrimaryVertices(*esd);\r
1605\r
1606 //setting best TOF PID\r
39d6c296 1607 AliESDInputHandler* esdH = dynamic_cast<AliESDInputHandler*>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());\r
1608 if (esdH)\r
1609 fESDpid = esdH->GetESDpid();\r
c7ef1c47 1610\r
39d6c296 1611 if (fIsPidOwner && fESDpid){\r
75a996ab 1612 delete fESDpid;\r
1613 fESDpid = 0;\r
1614 }\r
c7ef1c47 1615 if(!fESDpid)\r
1616 { //in case of no Tender attached \r
1617 fESDpid = new AliESDpid;\r
1618 fIsPidOwner = kTRUE;\r
1619 }\r
1620 \r
1621 if(!esd->GetTOFHeader())\r
1622 { //protection in case the pass2 LHC10b,c,d have been processed without tender. \r
1623 Float_t t0spread[10];\r
1624 Float_t intrinsicTOFres=100; //ps ok for LHC10b,c,d pass2!! \r
1625 for (Int_t i=0; i<10; i++) t0spread[i] = (TMath::Sqrt(esd->GetSigma2DiamondZ()))/0.03; //0.03 to convert from cm to ps\r
1626 fESDpid->GetTOFResponse().SetT0resolution(t0spread);\r
1627 fESDpid->GetTOFResponse().SetTimeResolution(intrinsicTOFres);\r
1628 \r
1629 fESDpid->SetTOFResponse(esd, (AliESDpid::EStartTimeType_t)fTimeZeroType); \r
1630 }\r
1631 \r
1632 if(esd->GetTOFHeader() && fIsPidOwner) fESDpid->SetTOFResponse(esd, (AliESDpid::EStartTimeType_t)fTimeZeroType); //in case of AOD production strating form LHC10e without Tender. \r
1633 \r
1634 if ( fAreCascadesEnabled ) ConvertCascades(*esd);\r
1635\r
1636 if ( fAreV0sEnabled ) ConvertV0s(*esd);\r
1637 \r
1638 if ( fAreKinksEnabled ) ConvertKinks(*esd);\r
1639 \r
1640 if ( fAreTracksEnabled ) ConvertTracks(*esd);\r
1641 \r
1642 // Update number of AOD tracks in header at the end of track loop (M.G.)\r
1643 header->SetRefMultiplicity(fNumberOfTracks);\r
1644 header->SetRefMultiplicityPos(fNumberOfPositiveTracks);\r
1645 header->SetRefMultiplicityNeg(fNumberOfTracks - fNumberOfPositiveTracks);\r
1646\r
1647 if ( fTPCOnlyFilterMask ) ConvertTPCOnlyTracks(*esd);\r
1648 \r
1649 if ( fArePmdClustersEnabled ) ConvertPmdClusters(*esd);\r
1650 \r
1651 if ( fAreCaloClustersEnabled ) ConvertCaloClusters(*esd);\r
1652 \r
1653 if ( fAreEMCALCellsEnabled )ConvertEMCALCells(*esd);\r
1654 \r
1655 if ( fArePHOSCellsEnabled )ConvertPHOSCells(*esd);\r
1656 \r
1657 if ( fAreTrackletsEnabled ) ConvertTracklets(*esd);\r
1658 \r
1659 delete fAODTrackRefs; fAODTrackRefs=0x0;\r
1660 delete fAODV0VtxRefs; fAODV0VtxRefs=0x0;\r
1661 delete fAODV0Refs; fAODV0Refs=0x0;\r
1662 \r
1663 delete[] fUsedTrack; fUsedTrack=0x0;\r
1664 delete[] fUsedV0; fUsedV0=0x0;\r
1665 delete[] fUsedKink; fUsedKink=0x0;\r
67456c92 1666\r
1667 if ( fIsPidOwner){\r
1668 delete fESDpid;\r
1669 fESDpid = 0x0;\r
1670 }\r
1671\r
1672\r
c7ef1c47 1673}\r
1674\r
1675\r
1676//______________________________________________________________________________\r
1a3579a6 1677void AliAnalysisTaskESDfilter::SetAODPID(AliESDtrack *esdtrack, AliAODTrack *aodtrack, AliAODPid *detpid)\r
f7ed7e88 1678{\r
1679 //\r
1680 // Setter for the raw PID detector signals\r
1681 //\r
1682\r
1683 // Save PID object for candidate electrons\r
1684 Bool_t pidSave = kFALSE;\r
1685 if (fTrackFilter) {\r
0c6c629b 1686 Bool_t selectInfo = fTrackFilter->IsSelected((char*) "Electrons");\r
f7ed7e88 1687 if (selectInfo) pidSave = kTRUE;\r
1688 }\r
1689\r
1690\r
1691 // Tracks passing pt cut \r
1692 if(esdtrack->Pt()>fHighPthreshold) {\r
1693 pidSave = kTRUE;\r
1694 } else {\r
1695 if(fPtshape){\r
1696 if(esdtrack->Pt()> fPtshape->GetXmin()){\r
1697 Double_t y = fPtshape->Eval(esdtrack->Pt())/fPtshape->Eval(fHighPthreshold);\r
1698 if(gRandom->Rndm(0)<1./y){\r
1699 pidSave = kTRUE;\r
1700 }//end rndm\r
1701 }//end if p < pmin\r
1702 }//end if p function\r
1703 }// end else\r
1704\r
1705 if (pidSave) {\r
6ece2430 1706 if(!aodtrack->GetDetPid()){// prevent memory leak when calling SetAODPID twice for the same track\r
f7ed7e88 1707 detpid = new AliAODPid();\r
1a3579a6 1708 SetDetectorRawSignals(detpid,esdtrack);\r
f7ed7e88 1709 aodtrack->SetDetPID(detpid);\r
6ece2430 1710 }\r
f7ed7e88 1711 }\r
1712}\r
1713\r
c7ef1c47 1714//______________________________________________________________________________\r
1a3579a6 1715void AliAnalysisTaskESDfilter::SetDetectorRawSignals(AliAODPid *aodpid, AliESDtrack *track)\r
f7ed7e88 1716{\r
1717//\r
1718//assignment of the detector signals (AliXXXesdPID inspired)\r
1719//\r
1720 if(!track){\r
1721 AliInfo("no ESD track found. .....exiting");\r
1722 return;\r
1723 }\r
febefa51 1724 // TPC momentum\r
1725 const AliExternalTrackParam *in=track->GetInnerParam();\r
1726 if (in) {\r
1727 aodpid->SetTPCmomentum(in->GetP());\r
1728 }else{\r
1729 aodpid->SetTPCmomentum(-1.);\r
1730 }\r
1731\r
f7ed7e88 1732\r
1733 aodpid->SetITSsignal(track->GetITSsignal());\r
1734 aodpid->SetTPCsignal(track->GetTPCsignal());\r
597b16d8 1735 aodpid->SetTPCsignalN(track->GetTPCsignalN());\r
f7ed7e88 1736\r
febefa51 1737 //n TRD planes = 6\r
f7ed7e88 1738 Int_t nslices = track->GetNumberOfTRDslices()*6;\r
1739 Double_t *trdslices = new Double_t[nslices];\r
1740 for(Int_t iSl =0; iSl < track->GetNumberOfTRDslices(); iSl++) {\r
febefa51 1741 for(Int_t iPl =0; iPl<6; iPl++) trdslices[iPl*track->GetNumberOfTRDslices()+iSl] = track->GetTRDslice(iPl,iSl);\r
1742 }\r
1743 \r
1744//TRD momentum\r
1745 for(Int_t iPl=0;iPl<6;iPl++){\r
1746 Double_t trdmom=track->GetTRDmomentum(iPl);\r
1747 aodpid->SetTRDmomentum(iPl,trdmom);\r
1748 }\r
f7ed7e88 1749\r
1750 aodpid->SetTRDsignal(track->GetNumberOfTRDslices()*6,trdslices);\r
b13bcc35 1751\r
1752 //TOF PID \r
f7ed7e88 1753 Double_t times[AliAODPid::kSPECIES]; track->GetIntegratedTimes(times);\r
1754 aodpid->SetIntegratedTimes(times);\r
1755\r
c7ef1c47 1756 Float_t tzeroTrack = fESDpid->GetTOFResponse().GetStartTime(track->P());\r
1757 aodpid->SetTOFsignal(track->GetTOFsignal()-tzeroTrack);\r
1758 \r
1759 Double_t tofRes[5];\r
1760 for (Int_t iMass=0; iMass<5; iMass++){\r
1761 tofRes[iMass]=(Double_t)fESDpid->GetTOFResponse().GetExpectedSigma(track->P(), times[iMass], AliPID::ParticleMass(iMass));\r
1762 }\r
1763 aodpid->SetTOFpidResolution(tofRes);\r
1764\r
1765 aodpid->SetHMPIDsignal(track->GetHMPIDsignal());\r
f7ed7e88 1766\r
1767 //Extrapolate track to EMCAL surface for AOD-level track-cluster matching\r
1768 Double_t emcpos[3] = {0.,0.,0.};\r
1769 Double_t emcmom[3] = {0.,0.,0.};\r
1770 aodpid->SetEMCALPosition(emcpos);\r
1771 aodpid->SetEMCALMomentum(emcmom);\r
1772\r
1773 AliExternalTrackParam *outerparam = (AliExternalTrackParam*)track->GetOuterParam();\r
1774 if(!outerparam) return;\r
1775\r
1776 //To be replaced by call to AliEMCALGeoUtils when the class becomes available\r
1a3579a6 1777 Bool_t okpos = outerparam->GetXYZ(emcpos);\r
1778 Bool_t okmom = outerparam->GetPxPyPz(emcmom);\r
f7ed7e88 1779 if(!(okpos && okmom)) return;\r
1780\r
1781 aodpid->SetEMCALPosition(emcpos);\r
1782 aodpid->SetEMCALMomentum(emcmom);\r
1783\r
1784}\r
1785\r
40a060e5 1786Double_t AliAnalysisTaskESDfilter::Chi2perNDF(AliESDtrack* track)\r
1787{\r
1788 // Calculate chi2 per ndf for track\r
1789 Int_t nClustersTPC = track->GetTPCNcls();\r
1790\r
1791 if ( nClustersTPC > 5) {\r
1792 return (track->GetTPCchi2()/Float_t(nClustersTPC - 5));\r
1793 } else {\r
1794 return (-1.);\r
1795 }\r
1796 }\r
1797\r
1798\r
c7ef1c47 1799//______________________________________________________________________________\r
f7ed7e88 1800void AliAnalysisTaskESDfilter::Terminate(Option_t */*option*/)\r
1801{\r
1802// Terminate analysis\r
1803//\r
1804 if (fDebug > 1) printf("AnalysisESDfilter: Terminate() \n");\r
1805}\r
1806\r
c7ef1c47 1807//______________________________________________________________________________\r
f7ed7e88 1808void AliAnalysisTaskESDfilter::PrintMCInfo(AliStack *pStack,Int_t label){\r
f690bf48 1809// Print MC info\r
f7ed7e88 1810 if(!pStack)return;\r
1811 label = TMath::Abs(label);\r
1812 TParticle *part = pStack->Particle(label);\r
1813 Printf("########################");\r
1814 Printf("%s:%d %d UniqueID %d PDG %d P %3.3f",(char*)__FILE__,__LINE__,label,part->GetUniqueID(),part->GetPdgCode(),part->P());\r
1815 part->Print();\r
1816 TParticle* mother = part;\r
1817 Int_t imo = part->GetFirstMother();\r
1818 Int_t nprim = pStack->GetNprimary();\r
1819 // while((imo >= nprim) && (mother->GetUniqueID() == 4)) {\r
1820 while((imo >= nprim)) {\r
1821 mother = pStack->Particle(imo);\r
1822 Printf("Mother %s:%d Label %d UniqueID %d PDG %d P %3.3f",(char*)__FILE__,__LINE__,imo,mother->GetUniqueID(),mother->GetPdgCode(),mother->P());\r
1823 mother->Print();\r
1824 imo = mother->GetFirstMother();\r
1825 }\r
1826 Printf("########################");\r
1827}\r