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