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