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