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