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