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