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