]> git.uio.no Git - u/mrichter/AliRoot.git/blob - ANALYSIS/AliAnalysisTaskESDfilter.cxx
Coding rule violations corrected
[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   Double_t p[3] = { 0. };\r
999 \r
1000   Double_t pDCA[3] = { 0. }; // momentum at DCA\r
1001   Double_t rDCA[3] = { 0. }; // position at DCA\r
1002   Float_t  dDCA[2] = {0.};    // DCA to the vertex d and z\r
1003   Float_t  cDCA[3] = {0.};    // covariance of impact parameters\r
1004 \r
1005 \r
1006   AliAODTrack* aodTrack(0x0);\r
1007   \r
1008   // account for change in pT after the constraint\r
1009   Float_t ptMax = 1E10;\r
1010   Float_t ptMin = 0;\r
1011   for(int i = 0;i<32;i++){\r
1012     if(fTPCConstrainedFilterMask&(1<<i)){\r
1013       AliESDtrackCuts*cuts = (AliESDtrackCuts*)fTrackFilter->GetCuts()->At(i);\r
1014       Float_t tmp1= 0,tmp2 = 0;\r
1015       cuts->GetPtRange(tmp1,tmp2);\r
1016       if(tmp1>ptMin)ptMin=tmp1;\r
1017       if(tmp2<ptMax)ptMax=tmp2;\r
1018     }\r
1019   } \r
1020 \r
1021   for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack) \r
1022   {\r
1023     AliESDtrack* esdTrack = esd.GetTrack(nTrack); //carefull do not modify it othwise  need to work with a copy \r
1024     \r
1025     UInt_t selectInfo = 0;\r
1026     Bool_t isHybridITSTPC = false;\r
1027     //\r
1028     // Track selection\r
1029     if (fTrackFilter) {\r
1030       selectInfo = fTrackFilter->IsSelected(esdTrack);\r
1031     }\r
1032 \r
1033     if(!(selectInfo&fHybridFilterMaskTPCCG)){\r
1034       // not already selected tracks, use second part of hybrid tracks\r
1035       isHybridITSTPC = true;\r
1036       // too save space one could only store these...\r
1037     }\r
1038 \r
1039     selectInfo &= fTPCConstrainedFilterMask;\r
1040     if (!selectInfo)continue;\r
1041     if (fWriteHybridTPCCOnly&&!isHybridITSTPC)continue; // write only complementary tracks\r
1042     // create a tpc only tracl\r
1043     AliESDtrack *track = AliESDtrackCuts::GetTPCOnlyTrack(const_cast<AliESDEvent*>(&esd),esdTrack->GetID());\r
1044     if(!track) continue;\r
1045     \r
1046     if(track->Pt()>0.)\r
1047     {\r
1048       // only constrain tracks above threshold\r
1049       AliExternalTrackParam exParam;\r
1050       // take the B-field from the ESD, no 3D fieldMap available at this point\r
1051       Bool_t relate = false;\r
1052       relate = track->RelateToVertexTPC(vtxSPD,esd.GetMagneticField(),kVeryBig,&exParam);\r
1053       if(!relate){\r
1054         delete track;\r
1055         continue;\r
1056       }\r
1057       // fetch the track parameters at the DCA (unconstraint)\r
1058       if(track->GetTPCInnerParam()){\r
1059         track->GetTPCInnerParam()->GetPxPyPz(pDCA);\r
1060         track->GetTPCInnerParam()->GetXYZ(rDCA);\r
1061       }\r
1062       // get the DCA to the vertex:\r
1063       track->GetImpactParametersTPC(dDCA,cDCA);\r
1064       // set the constrained parameters to the track\r
1065       track->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance());\r
1066     }\r
1067     \r
1068     track->GetPxPyPz(p);\r
1069 \r
1070     Float_t pT = track->Pt();\r
1071     if(pT<ptMin||pT>ptMax){\r
1072       delete track;\r
1073       continue;\r
1074     }\r
1075 \r
1076     // \r
1077 \r
1078 \r
1079     track->GetXYZ(pos);\r
1080     track->GetCovarianceXYZPxPyPz(covTr);\r
1081     track->GetESDpid(pid);\r
1082 \r
1083     if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());\r
1084     aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack((track->GetID()+1)*-1,\r
1085                                                             track->GetLabel(),\r
1086                                                             p,\r
1087                                                             kTRUE,\r
1088                                                             pos,\r
1089                                                             kFALSE,\r
1090                                                             covTr, \r
1091                                                             (Short_t)track->GetSign(),\r
1092                                                             track->GetITSClusterMap(), \r
1093                                                             pid,\r
1094                                                             fPrimaryVertex,\r
1095                                                             kTRUE, // check if this is right\r
1096                                                             vtx->UsesTrack(track->GetID()),\r
1097                                                             AliAODTrack::kPrimary, \r
1098                                                             selectInfo);\r
1099     aodTrack->SetIsHybridTPCConstrainedGlobal(isHybridITSTPC);    \r
1100     aodTrack->SetTPCClusterMap(track->GetTPCClusterMap());\r
1101     aodTrack->SetTPCSharedMap (track->GetTPCSharedMap());\r
1102     aodTrack->SetIsTPCConstrained(kTRUE);    \r
1103     Float_t ndf = track->GetTPCNcls()+1 - 5 ;\r
1104     if(ndf>0){\r
1105       aodTrack->SetChi2perNDF(track->GetConstrainedChi2TPC());\r
1106     }\r
1107     else{\r
1108       aodTrack->SetChi2perNDF(-1);\r
1109     }\r
1110 \r
1111     // set the DCA values to the AOD track\r
1112     aodTrack->SetPxPyPzAtDCA(pDCA[0],pDCA[1],pDCA[2]);\r
1113     aodTrack->SetXYAtDCA(rDCA[0],rDCA[1]);\r
1114     aodTrack->SetDCA(dDCA[0],dDCA[1]);\r
1115 \r
1116     aodTrack->SetFlags(track->GetStatus());\r
1117     aodTrack->SetTPCPointsF(track->GetTPCNclsF());\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   const AliESDVertex *vtx = esd.GetPrimaryVertex();\r
1169 \r
1170   // account for change in pT after the constraint\r
1171   Float_t ptMax = 1E10;\r
1172   Float_t ptMin = 0;\r
1173   for(int i = 0;i<32;i++){\r
1174     if(fGlobalConstrainedFilterMask&(1<<i)){\r
1175       AliESDtrackCuts*cuts = (AliESDtrackCuts*)fTrackFilter->GetCuts()->At(i);\r
1176       Float_t tmp1= 0,tmp2 = 0;\r
1177       cuts->GetPtRange(tmp1,tmp2);\r
1178       if(tmp1>ptMin)ptMin=tmp1;\r
1179       if(tmp2<ptMax)ptMax=tmp2;\r
1180     }\r
1181   } \r
1182 \r
1183 \r
1184 \r
1185  for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack) \r
1186   {\r
1187     AliESDtrack* esdTrack = esd.GetTrack(nTrack); //carefull do not modify it othwise  need to work with a copy \r
1188     const AliExternalTrackParam * exParamGC = esdTrack->GetConstrainedParam();\r
1189     if(!exParamGC)continue;\r
1190 \r
1191     UInt_t selectInfo = 0;\r
1192     Bool_t isHybridGC = false;\r
1193 \r
1194     //\r
1195     // Track selection\r
1196     if (fTrackFilter) {\r
1197       selectInfo = fTrackFilter->IsSelected(esdTrack);\r
1198     }\r
1199 \r
1200 \r
1201     if(!(selectInfo&fHybridFilterMaskGCG))isHybridGC = true;\r
1202     if (fWriteHybridGCOnly&&!isHybridGC)continue; // write only complementary tracks\r
1203 \r
1204     selectInfo &= fGlobalConstrainedFilterMask;\r
1205     if (!selectInfo)continue;\r
1206     // fetch the track parameters at the DCA (unconstrained)\r
1207     esdTrack->GetPxPyPz(pDCA);\r
1208     esdTrack->GetXYZ(rDCA);\r
1209     // get the DCA to the vertex:\r
1210     esdTrack->GetImpactParameters(dDCA,cDCA);\r
1211 \r
1212     esdTrack->GetConstrainedPxPyPz(p);\r
1213 \r
1214 \r
1215     Float_t pT = exParamGC->Pt();\r
1216     if(pT<ptMin||pT>ptMax){\r
1217       continue;\r
1218     }\r
1219 \r
1220 \r
1221     esdTrack->GetConstrainedXYZ(pos);\r
1222     exParamGC->GetCovarianceXYZPxPyPz(covTr);\r
1223     esdTrack->GetESDpid(pid);\r
1224     if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());\r
1225     aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack((esdTrack->GetID()+1)*-1,\r
1226                                                             esdTrack->GetLabel(),\r
1227                                                             p,\r
1228                                                             kTRUE,\r
1229                                                             pos,\r
1230                                                             kFALSE,\r
1231                                                             covTr, \r
1232                                                             (Short_t)esdTrack->GetSign(),\r
1233                                                             esdTrack->GetITSClusterMap(), \r
1234                                                             pid,\r
1235                                                             fPrimaryVertex,\r
1236                                                             kTRUE, // check if this is right\r
1237                                                             vtx->UsesTrack(esdTrack->GetID()),\r
1238                                                             AliAODTrack::kPrimary, \r
1239                                                             selectInfo);\r
1240     aodTrack->SetIsHybridGlobalConstrainedGlobal(isHybridGC);    \r
1241     aodTrack->SetIsGlobalConstrained(kTRUE);    \r
1242     aodTrack->SetTPCClusterMap(esdTrack->GetTPCClusterMap());\r
1243     aodTrack->SetTPCSharedMap (esdTrack->GetTPCSharedMap());\r
1244     Float_t ndf = esdTrack->GetTPCNcls()+1 - 5 ;\r
1245     if(ndf>0){\r
1246       aodTrack->SetChi2perNDF(esdTrack->GetConstrainedChi2TPC());\r
1247     }\r
1248     else{\r
1249       aodTrack->SetChi2perNDF(-1);\r
1250     }\r
1251 \r
1252     // set the DCA values to the AOD track\r
1253     aodTrack->SetPxPyPzAtDCA(pDCA[0],pDCA[1],pDCA[2]);\r
1254     aodTrack->SetXYAtDCA(rDCA[0],rDCA[1]);\r
1255     aodTrack->SetDCA(dDCA[0],dDCA[1]);\r
1256 \r
1257     aodTrack->SetFlags(esdTrack->GetStatus());\r
1258     aodTrack->SetTPCPointsF(esdTrack->GetTPCNclsF());\r
1259   } // end of loop on tracks\r
1260   \r
1261 }\r
1262 \r
1263 \r
1264 //______________________________________________________________________________\r
1265 void AliAnalysisTaskESDfilter::ConvertTracks(const AliESDEvent& esd)\r
1266 {\r
1267   // Tracks (primary and orphan)\r
1268 \r
1269   AliCodeTimerAuto("",0);\r
1270   \r
1271   AliDebug(1,Form("NUMBER OF ESD TRACKS %5d\n", esd.GetNumberOfTracks()));\r
1272   \r
1273   const AliESDVertex *vtx = esd.GetPrimaryVertex();\r
1274   Double_t p[3] = { 0. };\r
1275   Double_t pos[3] = { 0. };\r
1276   Double_t covTr[21] = { 0. };\r
1277   Double_t pid[10] = { 0. };\r
1278   AliAODTrack* aodTrack(0x0);\r
1279   AliAODPid* detpid(0x0);\r
1280   \r
1281   for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack) \r
1282   {\r
1283     if (fUsedTrack[nTrack]) continue;\r
1284     \r
1285     AliESDtrack *esdTrack = esd.GetTrack(nTrack);\r
1286     UInt_t selectInfo = 0;\r
1287     //\r
1288     // Track selection\r
1289     if (fTrackFilter) {\r
1290             selectInfo = fTrackFilter->IsSelected(esdTrack);\r
1291             if (!selectInfo && !vtx->UsesTrack(esdTrack->GetID())) continue;\r
1292     }\r
1293     \r
1294     \r
1295     esdTrack->GetPxPyPz(p);\r
1296     esdTrack->GetXYZ(pos);\r
1297     esdTrack->GetCovarianceXYZPxPyPz(covTr);\r
1298     esdTrack->GetESDpid(pid);\r
1299     if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());\r
1300     fPrimaryVertex->AddDaughter(aodTrack =\r
1301                          new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrack->GetID(),\r
1302                                                             esdTrack->GetLabel(),\r
1303                                                             p,\r
1304                                                             kTRUE,\r
1305                                                             pos,\r
1306                                                             kFALSE,\r
1307                                                             covTr, \r
1308                                                             (Short_t)esdTrack->GetSign(),\r
1309                                                             esdTrack->GetITSClusterMap(), \r
1310                                                             pid,\r
1311                                                             fPrimaryVertex,\r
1312                                                             kTRUE, // check if this is right\r
1313                                                             vtx->UsesTrack(esdTrack->GetID()),\r
1314                                                             AliAODTrack::kPrimary, \r
1315                                                             selectInfo)\r
1316                          );\r
1317     aodTrack->SetTPCClusterMap(esdTrack->GetTPCClusterMap());\r
1318     aodTrack->SetTPCSharedMap (esdTrack->GetTPCSharedMap());\r
1319     aodTrack->SetChi2perNDF(Chi2perNDF(esdTrack));\r
1320     aodTrack->SetTPCPointsF(esdTrack->GetTPCNclsF());\r
1321     if(esdTrack->IsEMCAL()) aodTrack->SetEMCALcluster(esdTrack->GetEMCALcluster());\r
1322     if(esdTrack->IsPHOS())  aodTrack->SetPHOScluster(esdTrack->GetPHOScluster());\r
1323 \r
1324     fAODTrackRefs->AddAt(aodTrack, nTrack);\r
1325     \r
1326     \r
1327     if (esdTrack->GetSign() > 0) ++fNumberOfPositiveTracks;\r
1328     aodTrack->SetFlags(esdTrack->GetStatus());\r
1329     aodTrack->ConvertAliPIDtoAODPID();\r
1330     SetAODPID(esdTrack,aodTrack,detpid);\r
1331   } // end of loop on tracks\r
1332 }\r
1333 \r
1334 //______________________________________________________________________________\r
1335 void AliAnalysisTaskESDfilter::ConvertPmdClusters(const AliESDEvent& esd)\r
1336 {\r
1337 // Convert PMD Clusters \r
1338   AliCodeTimerAuto("",0);\r
1339   Int_t jPmdClusters=0;\r
1340   // Access to the AOD container of PMD clusters\r
1341   TClonesArray &pmdClusters = *(AODEvent()->GetPmdClusters());\r
1342   for (Int_t iPmd = 0; iPmd < esd.GetNumberOfPmdTracks(); ++iPmd) {\r
1343     // file pmd clusters, to be revised!\r
1344     AliESDPmdTrack *pmdTrack = esd.GetPmdTrack(iPmd);\r
1345     Int_t nLabel = 0;\r
1346     Int_t *label = 0x0;\r
1347     Double_t posPmd[3] = { pmdTrack->GetClusterX(), pmdTrack->GetClusterY(), pmdTrack->GetClusterZ()};\r
1348     Double_t pidPmd[13] = { 0.}; // to be revised!\r
1349     // type not set!\r
1350     // assoc cluster not set\r
1351     new(pmdClusters[jPmdClusters++]) AliAODPmdCluster(iPmd, nLabel, label, pmdTrack->GetClusterADC(), posPmd, pidPmd);\r
1352   }\r
1353 }\r
1354 \r
1355 //______________________________________________________________________________\r
1356 void AliAnalysisTaskESDfilter::ConvertCaloClusters(const AliESDEvent& esd)\r
1357 {\r
1358 // Convert Calorimeter Clusters\r
1359   AliCodeTimerAuto("",0);\r
1360   \r
1361   // Access to the AOD container of clusters\r
1362   TClonesArray &caloClusters = *(AODEvent()->GetCaloClusters());\r
1363   Int_t jClusters(0);\r
1364   \r
1365   for (Int_t iClust=0; iClust<esd.GetNumberOfCaloClusters(); ++iClust) {\r
1366     \r
1367     AliESDCaloCluster * cluster = esd.GetCaloCluster(iClust);\r
1368     \r
1369     Int_t  id        = cluster->GetID();\r
1370     Int_t  nLabel    = cluster->GetNLabels();\r
1371     Int_t *labels    = cluster->GetLabels();\r
1372     if(labels){ \r
1373                   for(int i = 0;i < nLabel;++i){\r
1374                           if(fMChandler)fMChandler->SelectParticle(labels[i]);\r
1375                   }\r
1376           }             \r
1377     \r
1378     Float_t energy = cluster->E();\r
1379     Float_t posF[3] = { 0.};\r
1380     cluster->GetPosition(posF);\r
1381     \r
1382     AliAODCaloCluster *caloCluster = new(caloClusters[jClusters++]) AliAODCaloCluster(id,\r
1383                                                                                       nLabel,\r
1384                                                                                       labels,\r
1385                                                                                       energy,\r
1386                                                                                       posF,\r
1387                                                                                       NULL,\r
1388                                                                                       cluster->GetType(),0);\r
1389     \r
1390     caloCluster->SetCaloCluster(cluster->GetDistanceToBadChannel(),\r
1391                                 cluster->GetDispersion(),\r
1392                                 cluster->GetM20(), cluster->GetM02(),\r
1393                                 cluster->GetEmcCpvDistance(),  \r
1394                                 cluster->GetNExMax(),cluster->GetTOF()) ;\r
1395     \r
1396     caloCluster->SetPIDFromESD(cluster->GetPID());\r
1397     caloCluster->SetNCells(cluster->GetNCells());\r
1398     caloCluster->SetCellsAbsId(cluster->GetCellsAbsId());\r
1399     caloCluster->SetCellsAmplitudeFraction(cluster->GetCellsAmplitudeFraction());\r
1400     \r
1401     TArrayI* matchedT =         cluster->GetTracksMatched();\r
1402     if (fNumberOfTracks>0 && matchedT && cluster->GetTrackMatchedIndex() >= 0) {        \r
1403       for (Int_t im = 0; im < matchedT->GetSize(); im++) {\r
1404         Int_t iESDtrack = matchedT->At(im);;\r
1405         if (fAODTrackRefs->At(iESDtrack) != 0) {\r
1406           caloCluster->AddTrackMatched((AliAODTrack*)fAODTrackRefs->At(iESDtrack));\r
1407         }\r
1408       }\r
1409     }\r
1410     \r
1411   } \r
1412   caloClusters.Expand(jClusters); // resize TObjArray to 'remove' slots for pseudo clusters      \r
1413 }\r
1414 \r
1415 //______________________________________________________________________________\r
1416 void AliAnalysisTaskESDfilter::ConvertEMCALCells(const AliESDEvent& esd)\r
1417 {\r
1418 // Convert EMCAL Cells\r
1419   AliCodeTimerAuto("",0);\r
1420   // fill EMCAL cell info\r
1421   if (esd.GetEMCALCells()) { // protection against missing ESD information\r
1422     AliESDCaloCells &esdEMcells = *(esd.GetEMCALCells());\r
1423     Int_t nEMcell = esdEMcells.GetNumberOfCells() ;\r
1424     \r
1425     AliAODCaloCells &aodEMcells = *(AODEvent()->GetEMCALCells());\r
1426     aodEMcells.CreateContainer(nEMcell);\r
1427     aodEMcells.SetType(AliAODCaloCells::kEMCALCell);\r
1428     for (Int_t iCell = 0; iCell < nEMcell; iCell++) {      \r
1429       aodEMcells.SetCell(iCell,esdEMcells.GetCellNumber(iCell),esdEMcells.GetAmplitude(iCell));\r
1430     }\r
1431     aodEMcells.Sort();\r
1432   }\r
1433 }\r
1434 \r
1435 //______________________________________________________________________________\r
1436 void AliAnalysisTaskESDfilter::ConvertPHOSCells(const AliESDEvent& esd)\r
1437 {\r
1438 // Convert PHOS Cells\r
1439   AliCodeTimerAuto("",0);\r
1440   // fill PHOS cell info\r
1441   if (esd.GetPHOSCells()) { // protection against missing ESD information\r
1442     AliESDCaloCells &esdPHcells = *(esd.GetPHOSCells());\r
1443     Int_t nPHcell = esdPHcells.GetNumberOfCells() ;\r
1444     \r
1445     AliAODCaloCells &aodPHcells = *(AODEvent()->GetPHOSCells());\r
1446     aodPHcells.CreateContainer(nPHcell);\r
1447     aodPHcells.SetType(AliAODCaloCells::kPHOSCell);\r
1448     for (Int_t iCell = 0; iCell < nPHcell; iCell++) {      \r
1449       aodPHcells.SetCell(iCell,esdPHcells.GetCellNumber(iCell),esdPHcells.GetAmplitude(iCell));\r
1450     }\r
1451     aodPHcells.Sort();\r
1452   }\r
1453 }\r
1454 \r
1455 //______________________________________________________________________________\r
1456 void AliAnalysisTaskESDfilter::ConvertTracklets(const AliESDEvent& esd)\r
1457 {\r
1458   // tracklets    \r
1459   AliCodeTimerAuto("",0);\r
1460 \r
1461   AliAODTracklets &SPDTracklets = *(AODEvent()->GetTracklets());\r
1462   const AliMultiplicity *mult = esd.GetMultiplicity();\r
1463   if (mult) {\r
1464     if (mult->GetNumberOfTracklets()>0) {\r
1465       SPDTracklets.CreateContainer(mult->GetNumberOfTracklets());\r
1466       \r
1467       for (Int_t n=0; n<mult->GetNumberOfTracklets(); n++) {\r
1468         if(fMChandler){\r
1469           fMChandler->SelectParticle(mult->GetLabel(n, 0));\r
1470           fMChandler->SelectParticle(mult->GetLabel(n, 1));\r
1471         }\r
1472         SPDTracklets.SetTracklet(n, mult->GetTheta(n), mult->GetPhi(n), mult->GetDeltaPhi(n), mult->GetLabel(n, 0),mult->GetLabel(n, 1));\r
1473       }\r
1474     }\r
1475   } else {\r
1476     //Printf("ERROR: AliMultiplicity could not be retrieved from ESD");\r
1477   }\r
1478 }\r
1479 \r
1480 //______________________________________________________________________________\r
1481 void AliAnalysisTaskESDfilter::ConvertKinks(const AliESDEvent& esd)\r
1482 {\r
1483   AliCodeTimerAuto("",0);\r
1484   \r
1485   // Kinks: it is a big mess the access to the information in the kinks\r
1486   // The loop is on the tracks in order to find the mother and daugther of each kink\r
1487   \r
1488   Double_t covTr[21]={0.};\r
1489   Double_t pid[10]={0.};\r
1490   AliAODPid* detpid(0x0);\r
1491   \r
1492   fNumberOfKinks = esd.GetNumberOfKinks();\r
1493 \r
1494   const AliESDVertex* vtx = esd.GetPrimaryVertex();\r
1495   \r
1496   for (Int_t iTrack=0; iTrack<esd.GetNumberOfTracks(); ++iTrack) \r
1497   {\r
1498     AliESDtrack * esdTrack = esd.GetTrack(iTrack);\r
1499     \r
1500     Int_t ikink = esdTrack->GetKinkIndex(0);\r
1501     \r
1502     if (ikink && fNumberOfKinks) {\r
1503             // Negative kink index: mother, positive: daughter\r
1504             \r
1505             // Search for the second track of the kink\r
1506             \r
1507             for (Int_t jTrack = iTrack+1; jTrack<esd.GetNumberOfTracks(); ++jTrack) {\r
1508         \r
1509         AliESDtrack * esdTrack1 = esd.GetTrack(jTrack);\r
1510         \r
1511         Int_t jkink = esdTrack1->GetKinkIndex(0);\r
1512         \r
1513         if ( TMath::Abs(ikink)==TMath::Abs(jkink) ) {\r
1514           \r
1515           // The two tracks are from the same kink\r
1516           \r
1517           if (fUsedKink[TMath::Abs(ikink)-1]) continue; // skip used kinks\r
1518           \r
1519           Int_t imother = -1;\r
1520           Int_t idaughter = -1;\r
1521           \r
1522           if (ikink<0 && jkink>0) {\r
1523             \r
1524             imother = iTrack;\r
1525             idaughter = jTrack;\r
1526           }\r
1527           else if (ikink>0 && jkink<0) {\r
1528             \r
1529             imother = jTrack;\r
1530             idaughter = iTrack;\r
1531           }\r
1532           else {\r
1533             //                  cerr << "Error: Wrong combination of kink indexes: "\r
1534             //                       << ikink << " " << jkink << endl;\r
1535             continue;\r
1536           }\r
1537           \r
1538           // Add the mother track if it passed primary track selection cuts\r
1539           \r
1540           AliAODTrack * mother = NULL;\r
1541           \r
1542           UInt_t selectInfo = 0;\r
1543           if (fTrackFilter) {\r
1544             selectInfo = fTrackFilter->IsSelected(esd.GetTrack(imother));\r
1545             if (!selectInfo) continue;\r
1546           }\r
1547           \r
1548           if (!fUsedTrack[imother]) {\r
1549             \r
1550             fUsedTrack[imother] = kTRUE;\r
1551             \r
1552             AliESDtrack *esdTrackM = esd.GetTrack(imother);\r
1553             Double_t p[3] = { 0. };\r
1554             Double_t pos[3] = { 0. };\r
1555             esdTrackM->GetPxPyPz(p);\r
1556             esdTrackM->GetXYZ(pos);\r
1557             esdTrackM->GetCovarianceXYZPxPyPz(covTr);\r
1558             esdTrackM->GetESDpid(pid);\r
1559             if(fMChandler)fMChandler->SelectParticle(esdTrackM->GetLabel());\r
1560             mother = \r
1561             new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrackM->GetID(),\r
1562                                                esdTrackM->GetLabel(),\r
1563                                                p,\r
1564                                                kTRUE,\r
1565                                                pos,\r
1566                                                kFALSE,\r
1567                                                covTr, \r
1568                                                (Short_t)esdTrackM->GetSign(),\r
1569                                                esdTrackM->GetITSClusterMap(), \r
1570                                                pid,\r
1571                                                fPrimaryVertex,\r
1572                                                kTRUE, // check if this is right\r
1573                                                vtx->UsesTrack(esdTrack->GetID()),\r
1574                                                AliAODTrack::kPrimary,\r
1575                                                selectInfo);\r
1576             mother->SetTPCClusterMap(esdTrackM->GetTPCClusterMap());\r
1577             mother->SetTPCSharedMap (esdTrackM->GetTPCSharedMap());\r
1578             mother->SetChi2perNDF(Chi2perNDF(esdTrackM));\r
1579             mother->SetTPCPointsF(esdTrackM->GetTPCNclsF());\r
1580 \r
1581             fAODTrackRefs->AddAt(mother, imother);\r
1582             \r
1583             if (esdTrackM->GetSign() > 0) ++fNumberOfPositiveTracks;\r
1584             mother->SetFlags(esdTrackM->GetStatus());\r
1585             mother->ConvertAliPIDtoAODPID();\r
1586             fPrimaryVertex->AddDaughter(mother);\r
1587             mother->ConvertAliPIDtoAODPID();\r
1588             SetAODPID(esdTrackM,mother,detpid);\r
1589           }\r
1590           else {\r
1591             //                  cerr << "Error: event " << esd.GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1\r
1592             //                       << " track " << imother << " has already been used!" << endl;\r
1593           }\r
1594           \r
1595           // Add the kink vertex\r
1596           AliESDkink * kink = esd.GetKink(TMath::Abs(ikink)-1);\r
1597           \r
1598           AliAODVertex * vkink = \r
1599           new(Vertices()[fNumberOfVertices++]) AliAODVertex(kink->GetPosition(),\r
1600                                                   NULL,\r
1601                                                   0.,\r
1602                                                   mother,\r
1603                                                   esdTrack->GetID(),  // This is the track ID of the mother's track!\r
1604                                                   AliAODVertex::kKink);\r
1605           // Add the daughter track\r
1606           \r
1607           AliAODTrack * daughter = NULL;\r
1608           \r
1609           if (!fUsedTrack[idaughter]) {\r
1610             \r
1611             fUsedTrack[idaughter] = kTRUE;\r
1612             \r
1613             AliESDtrack *esdTrackD = esd.GetTrack(idaughter);\r
1614             Double_t p[3] = { 0. };\r
1615             Double_t pos[3] = { 0. };\r
1616 \r
1617             esdTrackD->GetPxPyPz(p);\r
1618             esdTrackD->GetXYZ(pos);\r
1619             esdTrackD->GetCovarianceXYZPxPyPz(covTr);\r
1620             esdTrackD->GetESDpid(pid);\r
1621             selectInfo = 0;\r
1622             if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdTrackD);\r
1623             if(fMChandler)fMChandler->SelectParticle(esdTrackD->GetLabel());\r
1624             daughter = \r
1625             new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrackD->GetID(),\r
1626                                                esdTrackD->GetLabel(),\r
1627                                                p,\r
1628                                                kTRUE,\r
1629                                                pos,\r
1630                                                kFALSE,\r
1631                                                covTr, \r
1632                                                (Short_t)esdTrackD->GetSign(),\r
1633                                                esdTrackD->GetITSClusterMap(), \r
1634                                                pid,\r
1635                                                vkink,\r
1636                                                kTRUE, // check if this is right\r
1637                                                vtx->UsesTrack(esdTrack->GetID()),\r
1638                                                AliAODTrack::kSecondary,\r
1639                                                selectInfo);\r
1640             daughter->SetTPCClusterMap(esdTrackD->GetTPCClusterMap());\r
1641             daughter->SetTPCSharedMap (esdTrackD->GetTPCSharedMap());\r
1642             daughter->SetTPCPointsF(esdTrackD->GetTPCNclsF());\r
1643             fAODTrackRefs->AddAt(daughter, idaughter);\r
1644             \r
1645             if (esdTrackD->GetSign() > 0) ++fNumberOfPositiveTracks;\r
1646             daughter->SetFlags(esdTrackD->GetStatus());\r
1647             daughter->ConvertAliPIDtoAODPID();\r
1648             vkink->AddDaughter(daughter);\r
1649             daughter->ConvertAliPIDtoAODPID();\r
1650             SetAODPID(esdTrackD,daughter,detpid);\r
1651           }\r
1652           else {\r
1653             //                  cerr << "Error: event " << esd.GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1\r
1654             //                       << " track " << idaughter << " has already been used!" << endl;\r
1655           }\r
1656         }\r
1657             }\r
1658     }      \r
1659   }\r
1660 }\r
1661 \r
1662 //______________________________________________________________________________\r
1663 void AliAnalysisTaskESDfilter::ConvertPrimaryVertices(const AliESDEvent& esd)\r
1664 {\r
1665   AliCodeTimerAuto("",0);\r
1666   \r
1667   // Access to the AOD container of vertices\r
1668   fNumberOfVertices = 0;\r
1669   \r
1670   Double_t pos[3] = { 0. };\r
1671   Double_t covVtx[6] = { 0. };\r
1672 \r
1673   // Add primary vertex. The primary tracks will be defined\r
1674   // after the loops on the composite objects (V0, cascades, kinks)\r
1675   const AliESDVertex *vtx = esd.GetPrimaryVertex();\r
1676   \r
1677   vtx->GetXYZ(pos); // position\r
1678   vtx->GetCovMatrix(covVtx); //covariance matrix\r
1679   \r
1680   fPrimaryVertex = new(Vertices()[fNumberOfVertices++])\r
1681   AliAODVertex(pos, covVtx, vtx->GetChi2toNDF(), NULL, -1, AliAODVertex::kPrimary);\r
1682   fPrimaryVertex->SetName(vtx->GetName());\r
1683   fPrimaryVertex->SetTitle(vtx->GetTitle());\r
1684   \r
1685   TString vtitle = vtx->GetTitle();\r
1686   if (!vtitle.Contains("VertexerTracks")) \r
1687     fPrimaryVertex->SetNContributors(vtx->GetNContributors());\r
1688   \r
1689   if (fDebug > 0) fPrimaryVertex->Print();  \r
1690   \r
1691   // Add SPD "main" vertex \r
1692   const AliESDVertex *vtxS = esd.GetPrimaryVertexSPD();\r
1693   vtxS->GetXYZ(pos); // position\r
1694   vtxS->GetCovMatrix(covVtx); //covariance matrix\r
1695   AliAODVertex * mVSPD = new(Vertices()[fNumberOfVertices++])\r
1696   AliAODVertex(pos, covVtx, vtxS->GetChi2toNDF(), NULL, -1, AliAODVertex::kMainSPD);\r
1697   mVSPD->SetName(vtxS->GetName());\r
1698   mVSPD->SetTitle(vtxS->GetTitle());\r
1699   mVSPD->SetNContributors(vtxS->GetNContributors()); \r
1700   \r
1701   // Add SPD pileup vertices\r
1702   for(Int_t iV=0; iV<esd.GetNumberOfPileupVerticesSPD(); ++iV)\r
1703   {\r
1704     const AliESDVertex *vtxP = esd.GetPileupVertexSPD(iV);\r
1705     vtxP->GetXYZ(pos); // position\r
1706     vtxP->GetCovMatrix(covVtx); //covariance matrix\r
1707     AliAODVertex * pVSPD =  new(Vertices()[fNumberOfVertices++])\r
1708     AliAODVertex(pos, covVtx, vtxP->GetChi2toNDF(), NULL, -1, AliAODVertex::kPileupSPD);\r
1709     pVSPD->SetName(vtxP->GetName());\r
1710     pVSPD->SetTitle(vtxP->GetTitle());\r
1711     pVSPD->SetNContributors(vtxP->GetNContributors()); \r
1712     pVSPD->SetBC(vtxP->GetBC());\r
1713   }\r
1714   \r
1715   // Add TRK pileup vertices\r
1716   for(Int_t iV=0; iV<esd.GetNumberOfPileupVerticesTracks(); ++iV)\r
1717   {\r
1718     const AliESDVertex *vtxP = esd.GetPileupVertexTracks(iV);\r
1719     vtxP->GetXYZ(pos); // position\r
1720     vtxP->GetCovMatrix(covVtx); //covariance matrix\r
1721     AliAODVertex * pVTRK = new(Vertices()[fNumberOfVertices++])\r
1722     AliAODVertex(pos, covVtx, vtxP->GetChi2toNDF(), NULL, -1, AliAODVertex::kPileupTracks);\r
1723     pVTRK->SetName(vtxP->GetName());\r
1724     pVTRK->SetTitle(vtxP->GetTitle());\r
1725     pVTRK->SetNContributors(vtxP->GetNContributors());\r
1726     pVTRK->SetBC(vtxP->GetBC());\r
1727   }\r
1728 }\r
1729 \r
1730 //______________________________________________________________________________\r
1731 void AliAnalysisTaskESDfilter::ConvertVZERO(const AliESDEvent& esd)\r
1732 {\r
1733   // Convert VZERO data\r
1734   AliAODVZERO* vzeroData = AODEvent()->GetVZEROData();\r
1735   *vzeroData = *(esd.GetVZEROData());\r
1736 }\r
1737 \r
1738 //______________________________________________________________________________\r
1739 void AliAnalysisTaskESDfilter::ConvertTZERO(const AliESDEvent& esd)\r
1740 {\r
1741   // Convert TZERO data\r
1742   const AliESDTZERO* esdTzero = esd.GetESDTZERO(); \r
1743   AliAODTZERO* aodTzero = AODEvent()->GetTZEROData();\r
1744 \r
1745   for (Int_t icase=0; icase<3; icase++){ \r
1746     aodTzero->SetT0TOF(    icase, esdTzero->GetT0TOF(icase));\r
1747     aodTzero->SetT0TOFbest(icase, esdTzero->GetT0TOFbest(icase)); \r
1748   }\r
1749   aodTzero->SetBackgroundFlag(esdTzero->GetBackgroundFlag());\r
1750   aodTzero->SetPileupFlag(esdTzero->GetPileupFlag());\r
1751   aodTzero->SetSatelliteFlag(esdTzero->GetSatellite()); \r
1752 \r
1753 }\r
1754 \r
1755 \r
1756 //______________________________________________________________________________\r
1757 void AliAnalysisTaskESDfilter::ConvertZDC(const AliESDEvent& esd)\r
1758 {\r
1759   // Convert ZDC data\r
1760   AliESDZDC* esdZDC = esd.GetZDCData();\r
1761   \r
1762   const Double_t zem1Energy = esdZDC->GetZEM1Energy();\r
1763   const Double_t zem2Energy = esdZDC->GetZEM2Energy();\r
1764    \r
1765   const Double_t *towZNC = esdZDC->GetZNCTowerEnergy();\r
1766   const Double_t *towZPC = esdZDC->GetZPCTowerEnergy();\r
1767   const Double_t *towZNA = esdZDC->GetZNATowerEnergy();\r
1768   const Double_t *towZPA = esdZDC->GetZPATowerEnergy();\r
1769   const Double_t *towZNCLG = esdZDC->GetZNCTowerEnergyLR();\r
1770   const Double_t *towZNALG = esdZDC->GetZNATowerEnergyLR();\r
1771   \r
1772   AliAODZDC* zdcAOD = AODEvent()->GetZDCData();\r
1773 \r
1774   zdcAOD->SetZEM1Energy(zem1Energy);\r
1775   zdcAOD->SetZEM2Energy(zem2Energy);\r
1776   zdcAOD->SetZNCTowers(towZNC, towZNCLG);\r
1777   zdcAOD->SetZNATowers(towZNA, towZNALG);\r
1778   zdcAOD->SetZPCTowers(towZPC);\r
1779   zdcAOD->SetZPATowers(towZPA);\r
1780   \r
1781   zdcAOD->SetZDCParticipants(esdZDC->GetZDCParticipants(), esdZDC->GetZDCPartSideA(), esdZDC->GetZDCPartSideC());\r
1782   zdcAOD->SetZDCImpactParameter(esdZDC->GetImpactParameter(), esdZDC->GetImpactParamSideA(), \r
1783         esdZDC->GetImpactParamSideC());\r
1784   zdcAOD->SetZDCTDCSum(esdZDC->GetZNTDCSum(0)); \r
1785   zdcAOD->SetZDCTDCDiff(esdZDC->GetZNTDCDiff(0));       \r
1786 \r
1787 }\r
1788 \r
1789 //______________________________________________________________________________\r
1790 void AliAnalysisTaskESDfilter::ConvertESDtoAOD() \r
1791 {\r
1792   // ESD Filter analysis task executed for each event\r
1793   \r
1794   AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());\r
1795   \r
1796   if(!esd)return;\r
1797 \r
1798   AliCodeTimerAuto("",0);\r
1799   \r
1800   fOldESDformat = ( esd->GetAliESDOld() != 0x0 );\r
1801   \r
1802   fNumberOfTracks = 0;\r
1803   fNumberOfPositiveTracks = 0;\r
1804   fNumberOfV0s = 0;\r
1805   fNumberOfVertices = 0;\r
1806   fNumberOfCascades = 0;\r
1807   fNumberOfKinks = 0;\r
1808     \r
1809   AliAODHeader* header = ConvertHeader(*esd);\r
1810 \r
1811   if ( fIsVZEROEnabled ) ConvertVZERO(*esd);\r
1812   if ( fIsTZEROEnabled ) ConvertTZERO(*esd);\r
1813   \r
1814   // Fetch Stack for debuggging if available \r
1815   fMChandler=0x0;\r
1816   if(MCEvent())\r
1817   {\r
1818     fMChandler = (AliMCEventHandler*) ((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler()); \r
1819   }\r
1820   \r
1821   // loop over events and fill them\r
1822   // Multiplicity information needed by the header (to be revised!)\r
1823   Int_t nTracks    = esd->GetNumberOfTracks();\r
1824   for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) esd->GetTrack(iTrack)->SetESDEvent(esd);\r
1825 \r
1826   // Update the header\r
1827 \r
1828   Int_t nV0s      = esd->GetNumberOfV0s();\r
1829   Int_t nCascades = esd->GetNumberOfCascades();\r
1830   Int_t nKinks    = esd->GetNumberOfKinks();\r
1831   Int_t nVertices = nV0s + nCascades /*V0 wihtin cascade already counted*/+ nKinks + 1 /* = prim. vtx*/;\r
1832   Int_t nPileSPDVertices=1+esd->GetNumberOfPileupVerticesSPD(); // also SPD main vertex\r
1833   Int_t nPileTrkVertices=esd->GetNumberOfPileupVerticesTracks();\r
1834   nVertices+=nPileSPDVertices;\r
1835   nVertices+=nPileTrkVertices;\r
1836   Int_t nJets     = 0;\r
1837   Int_t nCaloClus = esd->GetNumberOfCaloClusters();\r
1838   Int_t nFmdClus  = 0;\r
1839   Int_t nPmdClus  = esd->GetNumberOfPmdTracks();\r
1840     \r
1841   AliDebug(1,Form("   NV0=%d  NCASCADES=%d  NKINKS=%d", nV0s, nCascades, nKinks));\r
1842        \r
1843   AODEvent()->ResetStd(nTracks, nVertices, nV0s, nCascades, nJets, nCaloClus, nFmdClus, nPmdClus);\r
1844 \r
1845   if (nV0s > 0) \r
1846   {\r
1847     // RefArray to store a mapping between esd V0 number and newly created AOD-Vertex V0\r
1848     fAODV0VtxRefs = new TRefArray(nV0s);\r
1849     // RefArray to store the mapping between esd V0 number and newly created AOD-V0\r
1850     fAODV0Refs = new TRefArray(nV0s); \r
1851     // Array to take into account the V0s already added to the AOD (V0 within cascades)\r
1852     fUsedV0 = new Bool_t[nV0s];\r
1853     for (Int_t iV0=0; iV0<nV0s; ++iV0) fUsedV0[iV0]=kFALSE;\r
1854   }\r
1855   \r
1856   if (nTracks>0) \r
1857   {\r
1858     // RefArray to store the mapping between esd track number and newly created AOD-Track\r
1859     \r
1860     fAODTrackRefs = new TRefArray(nTracks);\r
1861 \r
1862     // Array to take into account the tracks already added to the AOD    \r
1863     fUsedTrack = new Bool_t[nTracks];\r
1864     for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) fUsedTrack[iTrack]=kFALSE;\r
1865   }\r
1866   \r
1867   // Array to take into account the kinks already added to the AOD\r
1868   if (nKinks>0) \r
1869   {\r
1870     fUsedKink = new Bool_t[nKinks];\r
1871     for (Int_t iKink=0; iKink<nKinks; ++iKink) fUsedKink[iKink]=kFALSE;\r
1872   }\r
1873     \r
1874   ConvertPrimaryVertices(*esd);\r
1875 \r
1876   //setting best TOF PID\r
1877   AliESDInputHandler* esdH = dynamic_cast<AliESDInputHandler*>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());\r
1878   if (esdH)\r
1879       fESDpid = esdH->GetESDpid();\r
1880 \r
1881   if (fIsPidOwner && fESDpid){\r
1882     delete fESDpid;\r
1883     fESDpid = 0;\r
1884   }\r
1885   if(!fESDpid)\r
1886   { //in case of no Tender attached \r
1887     fESDpid = new AliESDpid;\r
1888     fIsPidOwner = kTRUE;\r
1889   }\r
1890   \r
1891   if(!esd->GetTOFHeader())\r
1892   { //protection in case the pass2 LHC10b,c,d have been processed without tender. \r
1893     Float_t t0spread[10];\r
1894     Float_t intrinsicTOFres=100; //ps ok for LHC10b,c,d pass2!! \r
1895     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
1896     fESDpid->GetTOFResponse().SetT0resolution(t0spread);\r
1897     fESDpid->GetTOFResponse().SetTimeResolution(intrinsicTOFres);\r
1898     \r
1899     fESDpid->SetTOFResponse(esd, (AliESDpid::EStartTimeType_t)fTimeZeroType);    \r
1900   }\r
1901   \r
1902   if(esd->GetTOFHeader() && fIsPidOwner) fESDpid->SetTOFResponse(esd, (AliESDpid::EStartTimeType_t)fTimeZeroType); //in case of AOD production strating form LHC10e without Tender. \r
1903   \r
1904   if ( fAreCascadesEnabled ) ConvertCascades(*esd);\r
1905 \r
1906   if ( fAreV0sEnabled ) ConvertV0s(*esd);\r
1907   \r
1908   if ( fAreKinksEnabled ) ConvertKinks(*esd);\r
1909   \r
1910   if ( fAreTracksEnabled ) ConvertTracks(*esd);\r
1911   \r
1912   // Update number of AOD tracks in header at the end of track loop (M.G.)\r
1913   header->SetRefMultiplicity(fNumberOfTracks);\r
1914   header->SetRefMultiplicityPos(fNumberOfPositiveTracks);\r
1915   header->SetRefMultiplicityNeg(fNumberOfTracks - fNumberOfPositiveTracks);\r
1916 \r
1917   if ( fTPCConstrainedFilterMask ) ConvertTPCOnlyTracks(*esd);\r
1918   if( fGlobalConstrainedFilterMask) ConvertGlobalConstrainedTracks(*esd);  \r
1919 \r
1920   if ( fArePmdClustersEnabled ) ConvertPmdClusters(*esd);\r
1921   \r
1922   if ( fAreCaloClustersEnabled ) ConvertCaloClusters(*esd);\r
1923   \r
1924   if ( fAreEMCALCellsEnabled )ConvertEMCALCells(*esd);\r
1925   \r
1926   if ( fArePHOSCellsEnabled )ConvertPHOSCells(*esd);\r
1927   \r
1928   if ( fAreTrackletsEnabled ) ConvertTracklets(*esd);\r
1929   \r
1930   delete fAODTrackRefs; fAODTrackRefs=0x0;\r
1931   delete fAODV0VtxRefs; fAODV0VtxRefs=0x0;\r
1932   delete fAODV0Refs; fAODV0Refs=0x0;\r
1933   \r
1934   delete[] fUsedTrack; fUsedTrack=0x0;\r
1935   delete[] fUsedV0; fUsedV0=0x0;\r
1936   delete[] fUsedKink; fUsedKink=0x0;\r
1937 \r
1938   if ( fIsPidOwner){\r
1939     delete fESDpid;\r
1940     fESDpid = 0x0;\r
1941   }\r
1942 \r
1943 \r
1944 }\r
1945 \r
1946 \r
1947 //______________________________________________________________________________\r
1948 void AliAnalysisTaskESDfilter::SetAODPID(AliESDtrack *esdtrack, AliAODTrack *aodtrack, AliAODPid *detpid)\r
1949 {\r
1950   //\r
1951   // Setter for the raw PID detector signals\r
1952   //\r
1953 \r
1954   // Save PID object for candidate electrons\r
1955     Bool_t pidSave = kFALSE;\r
1956     if (fTrackFilter) {\r
1957         Bool_t selectInfo = fTrackFilter->IsSelected((char*) "Electrons");\r
1958         if (selectInfo)  pidSave = kTRUE;\r
1959     }\r
1960 \r
1961 \r
1962     // Tracks passing pt cut \r
1963     if(esdtrack->Pt()>fHighPthreshold) {\r
1964         pidSave = kTRUE;\r
1965     } else {\r
1966         if(fPtshape){\r
1967             if(esdtrack->Pt()> fPtshape->GetXmin()){\r
1968                 Double_t y = fPtshape->Eval(esdtrack->Pt())/fPtshape->Eval(fHighPthreshold);\r
1969                 if(gRandom->Rndm(0)<1./y){\r
1970                     pidSave = kTRUE;\r
1971                 }//end rndm\r
1972             }//end if p < pmin\r
1973         }//end if p function\r
1974     }// end else\r
1975 \r
1976     if (pidSave) {\r
1977       if(!aodtrack->GetDetPid()){// prevent memory leak when calling SetAODPID twice for the same track\r
1978         detpid = new AliAODPid();\r
1979         SetDetectorRawSignals(detpid,esdtrack);\r
1980         aodtrack->SetDetPID(detpid);\r
1981       }\r
1982     }\r
1983 }\r
1984 \r
1985 //______________________________________________________________________________\r
1986 void AliAnalysisTaskESDfilter::SetDetectorRawSignals(AliAODPid *aodpid, AliESDtrack *track)\r
1987 {\r
1988 //\r
1989 //assignment of the detector signals (AliXXXesdPID inspired)\r
1990 //\r
1991  if(!track){\r
1992  AliInfo("no ESD track found. .....exiting");\r
1993  return;\r
1994  }\r
1995  // TPC momentum\r
1996  const AliExternalTrackParam *in=track->GetInnerParam();\r
1997  if (in) {\r
1998    aodpid->SetTPCmomentum(in->GetP());\r
1999  }else{\r
2000    aodpid->SetTPCmomentum(-1.);\r
2001  }\r
2002 \r
2003 \r
2004  aodpid->SetITSsignal(track->GetITSsignal());\r
2005  Double_t itsdedx[4]; // dE/dx samples for individual ITS layers\r
2006  track->GetITSdEdxSamples(itsdedx);\r
2007  aodpid->SetITSdEdxSamples(itsdedx);\r
2008 \r
2009  aodpid->SetTPCsignal(track->GetTPCsignal());\r
2010  aodpid->SetTPCsignalN(track->GetTPCsignalN());\r
2011 \r
2012  //n TRD planes = 6\r
2013  Int_t nslices = track->GetNumberOfTRDslices()*6;\r
2014  TArrayD trdslices(nslices);\r
2015  for(Int_t iSl =0; iSl < track->GetNumberOfTRDslices(); iSl++) {\r
2016    for(Int_t iPl =0; iPl<6; iPl++) trdslices[iPl*track->GetNumberOfTRDslices()+iSl] = track->GetTRDslice(iPl,iSl);\r
2017  }\r
2018  \r
2019 //TRD momentum\r
2020  for(Int_t iPl=0;iPl<6;iPl++){\r
2021    Double_t trdmom=track->GetTRDmomentum(iPl);\r
2022    aodpid->SetTRDmomentum(iPl,trdmom);\r
2023  }\r
2024 \r
2025  aodpid->SetTRDsignal(track->GetNumberOfTRDslices()*6,trdslices.GetArray());\r
2026 \r
2027  //TRD clusters and tracklets\r
2028  aodpid->SetTRDncls(track->GetTRDncls());\r
2029  aodpid->SetTRDntrackletsPID(track->GetTRDntrackletsPID());\r
2030  \r
2031  //TOF PID  \r
2032  Double_t times[AliAODPid::kSPECIES]; track->GetIntegratedTimes(times);\r
2033  aodpid->SetIntegratedTimes(times);\r
2034 \r
2035   Float_t tzeroTrack = fESDpid->GetTOFResponse().GetStartTime(track->P());\r
2036   aodpid->SetTOFsignal(track->GetTOFsignal()-tzeroTrack);\r
2037   \r
2038   Double_t tofRes[5];\r
2039   for (Int_t iMass=0; iMass<5; iMass++){\r
2040     tofRes[iMass]=(Double_t)fESDpid->GetTOFResponse().GetExpectedSigma(track->P(), times[iMass], AliPID::ParticleMass(iMass));\r
2041   }\r
2042   aodpid->SetTOFpidResolution(tofRes);\r
2043 \r
2044   aodpid->SetHMPIDsignal(track->GetHMPIDsignal());\r
2045 \r
2046  //Extrapolate track to EMCAL surface for AOD-level track-cluster matching\r
2047  Double_t emcpos[3] = {0.,0.,0.};\r
2048  Double_t emcmom[3] = {0.,0.,0.};\r
2049  aodpid->SetEMCALPosition(emcpos);\r
2050  aodpid->SetEMCALMomentum(emcmom);\r
2051 \r
2052  Double_t hmpPid[5] = {0};\r
2053  track->GetHMPIDpid(hmpPid);\r
2054  aodpid->SetHMPIDprobs(hmpPid);\r
2055 \r
2056  AliExternalTrackParam *outerparam = (AliExternalTrackParam*)track->GetOuterParam();\r
2057  if(!outerparam) return;\r
2058 \r
2059  //To be replaced by call to AliEMCALGeoUtils when the class becomes available\r
2060  Bool_t okpos = outerparam->GetXYZ(emcpos);\r
2061  Bool_t okmom = outerparam->GetPxPyPz(emcmom);\r
2062  if(!(okpos && okmom)) return;\r
2063 \r
2064  aodpid->SetEMCALPosition(emcpos);\r
2065  aodpid->SetEMCALMomentum(emcmom);\r
2066 \r
2067 }\r
2068 \r
2069 Double_t  AliAnalysisTaskESDfilter::Chi2perNDF(AliESDtrack* track)\r
2070 {\r
2071     // Calculate chi2 per ndf for track\r
2072     Int_t  nClustersTPC = track->GetTPCNcls();\r
2073 \r
2074     if ( nClustersTPC > 5) {\r
2075        return (track->GetTPCchi2()/Float_t(nClustersTPC - 5));\r
2076     } else {\r
2077        return (-1.);\r
2078     }\r
2079  }\r
2080 \r
2081 \r
2082 //______________________________________________________________________________\r
2083 void AliAnalysisTaskESDfilter::Terminate(Option_t */*option*/)\r
2084 {\r
2085 // Terminate analysis\r
2086 //\r
2087     if (fDebug > 1) printf("AnalysisESDfilter: Terminate() \n");\r
2088 }\r
2089 \r
2090 //______________________________________________________________________________\r
2091 void  AliAnalysisTaskESDfilter::PrintMCInfo(AliStack *pStack,Int_t label){\r
2092 // Print MC info\r
2093   if(!pStack)return;\r
2094   label = TMath::Abs(label);\r
2095   TParticle *part = pStack->Particle(label);\r
2096   Printf("########################");\r
2097   Printf("%s:%d %d UniqueID %d PDG %d P %3.3f",(char*)__FILE__,__LINE__,label,part->GetUniqueID(),part->GetPdgCode(),part->P());\r
2098   part->Print();\r
2099   TParticle* mother = part;\r
2100   Int_t imo = part->GetFirstMother();\r
2101   Int_t nprim = pStack->GetNprimary();\r
2102   //  while((imo >= nprim) && (mother->GetUniqueID() == 4)) {\r
2103   while((imo >= nprim)) {\r
2104     mother =  pStack->Particle(imo);\r
2105     Printf("Mother %s:%d Label %d UniqueID %d PDG %d P %3.3f",(char*)__FILE__,__LINE__,imo,mother->GetUniqueID(),mother->GetPdgCode(),mother->P());\r
2106     mother->Print();\r
2107     imo =  mother->GetFirstMother();\r
2108   }\r
2109   Printf("########################");\r
2110 }\r