Changes for #93916 EMCAL commit attached patch and port to the release
[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));\r
1612     }\r
1613     aodEMcells.Sort();\r
1614   }\r
1615 }\r
1616 \r
1617 //______________________________________________________________________________\r
1618 void AliAnalysisTaskESDfilter::ConvertPHOSCells(const AliESDEvent& esd)\r
1619 {\r
1620 // Convert PHOS Cells\r
1621   AliCodeTimerAuto("",0);\r
1622   // fill PHOS cell info\r
1623   if (esd.GetPHOSCells()) { // protection against missing ESD information\r
1624     AliESDCaloCells &esdPHcells = *(esd.GetPHOSCells());\r
1625     Int_t nPHcell = esdPHcells.GetNumberOfCells() ;\r
1626     \r
1627     AliAODCaloCells &aodPHcells = *(AODEvent()->GetPHOSCells());\r
1628     aodPHcells.CreateContainer(nPHcell);\r
1629     aodPHcells.SetType(AliAODCaloCells::kPHOSCell);\r
1630     for (Int_t iCell = 0; iCell < nPHcell; iCell++) {      \r
1631       aodPHcells.SetCell(iCell,esdPHcells.GetCellNumber(iCell),esdPHcells.GetAmplitude(iCell));\r
1632     }\r
1633     aodPHcells.Sort();\r
1634   }\r
1635 }\r
1636 \r
1637 //______________________________________________________________________________\r
1638 void AliAnalysisTaskESDfilter::ConvertTracklets(const AliESDEvent& esd)\r
1639 {\r
1640   // tracklets    \r
1641   AliCodeTimerAuto("",0);\r
1642 \r
1643   AliAODTracklets &SPDTracklets = *(AODEvent()->GetTracklets());\r
1644   const AliMultiplicity *mult = esd.GetMultiplicity();\r
1645   if (mult) {\r
1646     if (mult->GetNumberOfTracklets()>0) {\r
1647       SPDTracklets.CreateContainer(mult->GetNumberOfTracklets());\r
1648       \r
1649       for (Int_t n=0; n<mult->GetNumberOfTracklets(); n++) {\r
1650         if(fMChandler){\r
1651           fMChandler->SelectParticle(mult->GetLabel(n, 0));\r
1652           fMChandler->SelectParticle(mult->GetLabel(n, 1));\r
1653         }\r
1654         SPDTracklets.SetTracklet(n, mult->GetTheta(n), mult->GetPhi(n), mult->GetDeltaPhi(n), mult->GetLabel(n, 0),mult->GetLabel(n, 1));\r
1655       }\r
1656     }\r
1657   } else {\r
1658     //Printf("ERROR: AliMultiplicity could not be retrieved from ESD");\r
1659   }\r
1660 }\r
1661 \r
1662 //______________________________________________________________________________\r
1663 void AliAnalysisTaskESDfilter::ConvertKinks(const AliESDEvent& esd)\r
1664 {\r
1665   AliCodeTimerAuto("",0);\r
1666   \r
1667   // Kinks: it is a big mess the access to the information in the kinks\r
1668   // The loop is on the tracks in order to find the mother and daugther of each kink\r
1669   \r
1670   Double_t covTr[21]={0.};\r
1671   Double_t pid[10]={0.};\r
1672   AliAODPid* detpid(0x0);\r
1673   \r
1674   fNumberOfKinks = esd.GetNumberOfKinks();\r
1675 \r
1676   const AliESDVertex* vtx = esd.GetPrimaryVertex();\r
1677   \r
1678   for (Int_t iTrack=0; iTrack<esd.GetNumberOfTracks(); ++iTrack) \r
1679   {\r
1680     AliESDtrack * esdTrack = esd.GetTrack(iTrack);\r
1681     \r
1682     Int_t ikink = esdTrack->GetKinkIndex(0);\r
1683     \r
1684     if (ikink && fNumberOfKinks) {\r
1685             // Negative kink index: mother, positive: daughter\r
1686             \r
1687             // Search for the second track of the kink\r
1688             \r
1689             for (Int_t jTrack = iTrack+1; jTrack<esd.GetNumberOfTracks(); ++jTrack) {\r
1690         \r
1691         AliESDtrack * esdTrack1 = esd.GetTrack(jTrack);\r
1692         \r
1693         Int_t jkink = esdTrack1->GetKinkIndex(0);\r
1694         \r
1695         if ( TMath::Abs(ikink)==TMath::Abs(jkink) ) {\r
1696           \r
1697           // The two tracks are from the same kink\r
1698           \r
1699           if (fUsedKink[TMath::Abs(ikink)-1]) continue; // skip used kinks\r
1700           \r
1701           Int_t imother = -1;\r
1702           Int_t idaughter = -1;\r
1703           \r
1704           if (ikink<0 && jkink>0) {\r
1705             \r
1706             imother = iTrack;\r
1707             idaughter = jTrack;\r
1708           }\r
1709           else if (ikink>0 && jkink<0) {\r
1710             \r
1711             imother = jTrack;\r
1712             idaughter = iTrack;\r
1713           }\r
1714           else {\r
1715             //                  cerr << "Error: Wrong combination of kink indexes: "\r
1716             //                       << ikink << " " << jkink << endl;\r
1717             continue;\r
1718           }\r
1719           \r
1720           // Add the mother track if it passed primary track selection cuts\r
1721           \r
1722           AliAODTrack * mother = NULL;\r
1723           \r
1724           UInt_t selectInfo = 0;\r
1725           if (fTrackFilter) {\r
1726             selectInfo = fTrackFilter->IsSelected(esd.GetTrack(imother));\r
1727             if (!selectInfo) continue;\r
1728           }\r
1729           \r
1730           if (!fUsedTrack[imother]) {\r
1731             \r
1732             fUsedTrack[imother] = kTRUE;\r
1733             \r
1734             AliESDtrack *esdTrackM = esd.GetTrack(imother);\r
1735             Double_t p[3] = { 0. };\r
1736             Double_t pos[3] = { 0. };\r
1737             esdTrackM->GetPxPyPz(p);\r
1738             esdTrackM->GetXYZ(pos);\r
1739             esdTrackM->GetCovarianceXYZPxPyPz(covTr);\r
1740             esdTrackM->GetESDpid(pid);\r
1741             if(fMChandler)fMChandler->SelectParticle(esdTrackM->GetLabel());\r
1742             mother = \r
1743             new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrackM->GetID(),\r
1744                                                esdTrackM->GetLabel(),\r
1745                                                p,\r
1746                                                kTRUE,\r
1747                                                pos,\r
1748                                                kFALSE,\r
1749                                                covTr, \r
1750                                                (Short_t)esdTrackM->GetSign(),\r
1751                                                esdTrackM->GetITSClusterMap(), \r
1752                                                pid,\r
1753                                                fPrimaryVertex,\r
1754                                                kTRUE, // check if this is right\r
1755                                                vtx->UsesTrack(esdTrack->GetID()),\r
1756                                                AliAODTrack::kPrimary,\r
1757                                                selectInfo);\r
1758             mother->SetTPCFitMap(esdTrackM->GetTPCFitMap());\r
1759             mother->SetTPCClusterMap(esdTrackM->GetTPCClusterMap());\r
1760             mother->SetTPCSharedMap (esdTrackM->GetTPCSharedMap());\r
1761             mother->SetChi2perNDF(Chi2perNDF(esdTrackM));\r
1762             mother->SetTPCPointsF(esdTrackM->GetTPCNclsF());\r
1763 \r
1764             fAODTrackRefs->AddAt(mother, imother);\r
1765             \r
1766             if (esdTrackM->GetSign() > 0) ++fNumberOfPositiveTracks;\r
1767             mother->SetFlags(esdTrackM->GetStatus());\r
1768             mother->ConvertAliPIDtoAODPID();\r
1769             fPrimaryVertex->AddDaughter(mother);\r
1770             mother->ConvertAliPIDtoAODPID();\r
1771             SetAODPID(esdTrackM,mother,detpid);\r
1772           }\r
1773           else {\r
1774             //                  cerr << "Error: event " << esd.GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1\r
1775             //                       << " track " << imother << " has already been used!" << endl;\r
1776           }\r
1777           \r
1778           // Add the kink vertex\r
1779           AliESDkink * kink = esd.GetKink(TMath::Abs(ikink)-1);\r
1780           \r
1781           AliAODVertex * vkink = \r
1782           new(Vertices()[fNumberOfVertices++]) AliAODVertex(kink->GetPosition(),\r
1783                                                   NULL,\r
1784                                                   0.,\r
1785                                                   mother,\r
1786                                                   esdTrack->GetID(),  // This is the track ID of the mother's track!\r
1787                                                   AliAODVertex::kKink);\r
1788           // Add the daughter track\r
1789           \r
1790           AliAODTrack * daughter = NULL;\r
1791           \r
1792           if (!fUsedTrack[idaughter]) {\r
1793             \r
1794             fUsedTrack[idaughter] = kTRUE;\r
1795             \r
1796             AliESDtrack *esdTrackD = esd.GetTrack(idaughter);\r
1797             Double_t p[3] = { 0. };\r
1798             Double_t pos[3] = { 0. };\r
1799 \r
1800             esdTrackD->GetPxPyPz(p);\r
1801             esdTrackD->GetXYZ(pos);\r
1802             esdTrackD->GetCovarianceXYZPxPyPz(covTr);\r
1803             esdTrackD->GetESDpid(pid);\r
1804             selectInfo = 0;\r
1805             if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdTrackD);\r
1806             if(fMChandler)fMChandler->SelectParticle(esdTrackD->GetLabel());\r
1807             daughter = \r
1808             new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrackD->GetID(),\r
1809                                                esdTrackD->GetLabel(),\r
1810                                                p,\r
1811                                                kTRUE,\r
1812                                                pos,\r
1813                                                kFALSE,\r
1814                                                covTr, \r
1815                                                (Short_t)esdTrackD->GetSign(),\r
1816                                                esdTrackD->GetITSClusterMap(), \r
1817                                                pid,\r
1818                                                vkink,\r
1819                                                kTRUE, // check if this is right\r
1820                                                vtx->UsesTrack(esdTrack->GetID()),\r
1821                                                AliAODTrack::kSecondary,\r
1822                                                selectInfo);\r
1823             daughter->SetTPCFitMap(esdTrackD->GetTPCFitMap());\r
1824             daughter->SetTPCClusterMap(esdTrackD->GetTPCClusterMap());\r
1825             daughter->SetTPCSharedMap (esdTrackD->GetTPCSharedMap());\r
1826             daughter->SetTPCPointsF(esdTrackD->GetTPCNclsF());\r
1827             fAODTrackRefs->AddAt(daughter, idaughter);\r
1828             \r
1829             if (esdTrackD->GetSign() > 0) ++fNumberOfPositiveTracks;\r
1830             daughter->SetFlags(esdTrackD->GetStatus());\r
1831             daughter->ConvertAliPIDtoAODPID();\r
1832             vkink->AddDaughter(daughter);\r
1833             daughter->ConvertAliPIDtoAODPID();\r
1834             SetAODPID(esdTrackD,daughter,detpid);\r
1835           }\r
1836           else {\r
1837             //                  cerr << "Error: event " << esd.GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1\r
1838             //                       << " track " << idaughter << " has already been used!" << endl;\r
1839           }\r
1840         }\r
1841             }\r
1842     }      \r
1843   }\r
1844 }\r
1845 \r
1846 //______________________________________________________________________________\r
1847 void AliAnalysisTaskESDfilter::ConvertPrimaryVertices(const AliESDEvent& esd)\r
1848 {\r
1849   AliCodeTimerAuto("",0);\r
1850   \r
1851   // Access to the AOD container of vertices\r
1852   fNumberOfVertices = 0;\r
1853   \r
1854   Double_t pos[3] = { 0. };\r
1855   Double_t covVtx[6] = { 0. };\r
1856 \r
1857   // Add primary vertex. The primary tracks will be defined\r
1858   // after the loops on the composite objects (V0, cascades, kinks)\r
1859   const AliESDVertex *vtx = esd.GetPrimaryVertex();\r
1860   \r
1861   vtx->GetXYZ(pos); // position\r
1862   vtx->GetCovMatrix(covVtx); //covariance matrix\r
1863   \r
1864   fPrimaryVertex = new(Vertices()[fNumberOfVertices++])\r
1865   AliAODVertex(pos, covVtx, vtx->GetChi2toNDF(), NULL, -1, AliAODVertex::kPrimary);\r
1866   fPrimaryVertex->SetName(vtx->GetName());\r
1867   fPrimaryVertex->SetTitle(vtx->GetTitle());\r
1868   \r
1869   TString vtitle = vtx->GetTitle();\r
1870   if (!vtitle.Contains("VertexerTracks")) \r
1871     fPrimaryVertex->SetNContributors(vtx->GetNContributors());\r
1872   \r
1873   if (fDebug > 0) fPrimaryVertex->Print();  \r
1874   \r
1875   // Add SPD "main" vertex \r
1876   const AliESDVertex *vtxS = esd.GetPrimaryVertexSPD();\r
1877   vtxS->GetXYZ(pos); // position\r
1878   vtxS->GetCovMatrix(covVtx); //covariance matrix\r
1879   AliAODVertex * mVSPD = new(Vertices()[fNumberOfVertices++])\r
1880   AliAODVertex(pos, covVtx, vtxS->GetChi2toNDF(), NULL, -1, AliAODVertex::kMainSPD);\r
1881   mVSPD->SetName(vtxS->GetName());\r
1882   mVSPD->SetTitle(vtxS->GetTitle());\r
1883   mVSPD->SetNContributors(vtxS->GetNContributors()); \r
1884   \r
1885   // Add SPD pileup vertices\r
1886   for(Int_t iV=0; iV<esd.GetNumberOfPileupVerticesSPD(); ++iV)\r
1887   {\r
1888     const AliESDVertex *vtxP = esd.GetPileupVertexSPD(iV);\r
1889     vtxP->GetXYZ(pos); // position\r
1890     vtxP->GetCovMatrix(covVtx); //covariance matrix\r
1891     AliAODVertex * pVSPD =  new(Vertices()[fNumberOfVertices++])\r
1892     AliAODVertex(pos, covVtx, vtxP->GetChi2toNDF(), NULL, -1, AliAODVertex::kPileupSPD);\r
1893     pVSPD->SetName(vtxP->GetName());\r
1894     pVSPD->SetTitle(vtxP->GetTitle());\r
1895     pVSPD->SetNContributors(vtxP->GetNContributors()); \r
1896     pVSPD->SetBC(vtxP->GetBC());\r
1897   }\r
1898   \r
1899   // Add TRK pileup vertices\r
1900   for(Int_t iV=0; iV<esd.GetNumberOfPileupVerticesTracks(); ++iV)\r
1901   {\r
1902     const AliESDVertex *vtxP = esd.GetPileupVertexTracks(iV);\r
1903     vtxP->GetXYZ(pos); // position\r
1904     vtxP->GetCovMatrix(covVtx); //covariance matrix\r
1905     AliAODVertex * pVTRK = new(Vertices()[fNumberOfVertices++])\r
1906     AliAODVertex(pos, covVtx, vtxP->GetChi2toNDF(), NULL, -1, AliAODVertex::kPileupTracks);\r
1907     pVTRK->SetName(vtxP->GetName());\r
1908     pVTRK->SetTitle(vtxP->GetTitle());\r
1909     pVTRK->SetNContributors(vtxP->GetNContributors());\r
1910     pVTRK->SetBC(vtxP->GetBC());\r
1911   }\r
1912 }\r
1913 \r
1914 //______________________________________________________________________________\r
1915 void AliAnalysisTaskESDfilter::ConvertVZERO(const AliESDEvent& esd)\r
1916 {\r
1917   // Convert VZERO data\r
1918   AliAODVZERO* vzeroData = AODEvent()->GetVZEROData();\r
1919   *vzeroData = *(esd.GetVZEROData());\r
1920 }\r
1921 \r
1922 //______________________________________________________________________________\r
1923 void AliAnalysisTaskESDfilter::ConvertTZERO(const AliESDEvent& esd)\r
1924 {\r
1925   // Convert TZERO data\r
1926   const AliESDTZERO* esdTzero = esd.GetESDTZERO(); \r
1927   AliAODTZERO* aodTzero = AODEvent()->GetTZEROData();\r
1928 \r
1929   for (Int_t icase=0; icase<3; icase++){ \r
1930     aodTzero->SetT0TOF(    icase, esdTzero->GetT0TOF(icase));\r
1931     aodTzero->SetT0TOFbest(icase, esdTzero->GetT0TOFbest(icase)); \r
1932   }\r
1933   aodTzero->SetBackgroundFlag(esdTzero->GetBackgroundFlag());\r
1934   aodTzero->SetPileupFlag(esdTzero->GetPileupFlag());\r
1935   aodTzero->SetSatelliteFlag(esdTzero->GetSatellite()); \r
1936 \r
1937   Float_t rawTime[24];\r
1938   for(Int_t ipmt=0; ipmt<24; ipmt++)\r
1939     rawTime[ipmt] = esdTzero->GetTimeFull(ipmt,0);\r
1940    \r
1941   Int_t idxOfFirstPmtA = -1,       idxOfFirstPmtC = -1;\r
1942   Float_t timeOfFirstPmtA = 9999, timeOfFirstPmtC = 9999;\r
1943   for(int ipmt=0;  ipmt<12; ipmt++){\r
1944     if( rawTime[ipmt] > -200 && rawTime[ipmt] < timeOfFirstPmtC && rawTime[ipmt]!=0){\r
1945       timeOfFirstPmtC = rawTime[ipmt];\r
1946       idxOfFirstPmtC  = ipmt;\r
1947     }\r
1948   }\r
1949   for(int ipmt=12; ipmt<24; ipmt++){\r
1950     if( rawTime[ipmt] > -200 && rawTime[ipmt] < timeOfFirstPmtA && rawTime[ipmt]!=0 ){\r
1951       timeOfFirstPmtA = rawTime[ipmt];\r
1952       idxOfFirstPmtA  = ipmt;\r
1953     }\r
1954   }\r
1955 \r
1956   if(idxOfFirstPmtA != -1 && idxOfFirstPmtC != -1){\r
1957     //speed of light in cm/ns   TMath::C()*1e-7 \r
1958     Float_t vertexraw = TMath::C()*1e-7 * (rawTime[idxOfFirstPmtA] - rawTime[idxOfFirstPmtC])/2;\r
1959     aodTzero->SetT0VertexRaw( vertexraw );\r
1960   }else{\r
1961     aodTzero->SetT0VertexRaw(99999);\r
1962   }\r
1963 \r
1964 }\r
1965 \r
1966 \r
1967 //______________________________________________________________________________\r
1968 void AliAnalysisTaskESDfilter::ConvertZDC(const AliESDEvent& esd)\r
1969 {\r
1970   // Convert ZDC data\r
1971   AliESDZDC* esdZDC = esd.GetZDCData();\r
1972   \r
1973   const Double_t zem1Energy = esdZDC->GetZEM1Energy();\r
1974   const Double_t zem2Energy = esdZDC->GetZEM2Energy();\r
1975    \r
1976   const Double_t *towZNC = esdZDC->GetZNCTowerEnergy();\r
1977   const Double_t *towZPC = esdZDC->GetZPCTowerEnergy();\r
1978   const Double_t *towZNA = esdZDC->GetZNATowerEnergy();\r
1979   const Double_t *towZPA = esdZDC->GetZPATowerEnergy();\r
1980   const Double_t *towZNCLG = esdZDC->GetZNCTowerEnergyLR();\r
1981   const Double_t *towZNALG = esdZDC->GetZNATowerEnergyLR();\r
1982   \r
1983   AliAODZDC* zdcAOD = AODEvent()->GetZDCData();\r
1984 \r
1985   zdcAOD->SetZEM1Energy(zem1Energy);\r
1986   zdcAOD->SetZEM2Energy(zem2Energy);\r
1987   zdcAOD->SetZNCTowers(towZNC, towZNCLG);\r
1988   zdcAOD->SetZNATowers(towZNA, towZNALG);\r
1989   zdcAOD->SetZPCTowers(towZPC);\r
1990   zdcAOD->SetZPATowers(towZPA);\r
1991   \r
1992   zdcAOD->SetZDCParticipants(esdZDC->GetZDCParticipants(), esdZDC->GetZDCPartSideA(), esdZDC->GetZDCPartSideC());\r
1993   zdcAOD->SetZDCImpactParameter(esdZDC->GetImpactParameter(), esdZDC->GetImpactParamSideA(), \r
1994         esdZDC->GetImpactParamSideC());\r
1995   zdcAOD->SetZDCTDCSum(esdZDC->GetZNTDCSum(0)); \r
1996   zdcAOD->SetZDCTDCDiff(esdZDC->GetZNTDCDiff(0));       \r
1997 \r
1998 }\r
1999 \r
2000 //______________________________________________________________________________\r
2001 void AliAnalysisTaskESDfilter::ConvertESDtoAOD() \r
2002 {\r
2003   // ESD Filter analysis task executed for each event\r
2004   \r
2005   AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());\r
2006   \r
2007   if(!esd)return;\r
2008 \r
2009   AliCodeTimerAuto("",0);\r
2010   \r
2011   fOldESDformat = ( esd->GetAliESDOld() != 0x0 );\r
2012  \r
2013       // Reconstruct cascades and V0 here\r
2014   if (fIsV0CascadeRecoEnabled) {\r
2015     esd->ResetCascades();\r
2016     esd->ResetV0s();\r
2017 \r
2018     AliV0vertexer lV0vtxer;\r
2019     AliCascadeVertexer lCascVtxer;\r
2020 \r
2021     lV0vtxer.SetCuts(fV0Cuts);\r
2022     lCascVtxer.SetCuts(fCascadeCuts);\r
2023 \r
2024 \r
2025     lV0vtxer.Tracks2V0vertices(esd);\r
2026     lCascVtxer.V0sTracks2CascadeVertices(esd);\r
2027   }\r
2028 \r
2029  \r
2030   fNumberOfTracks = 0;\r
2031   fNumberOfPositiveTracks = 0;\r
2032   fNumberOfV0s = 0;\r
2033   fNumberOfVertices = 0;\r
2034   fNumberOfCascades = 0;\r
2035   fNumberOfKinks = 0;\r
2036     \r
2037   AliAODHeader* header = ConvertHeader(*esd);\r
2038 \r
2039   if ( fIsVZEROEnabled ) ConvertVZERO(*esd);\r
2040   if ( fIsTZEROEnabled ) ConvertTZERO(*esd);\r
2041   \r
2042   // Fetch Stack for debuggging if available \r
2043   fMChandler=0x0;\r
2044   if(MCEvent())\r
2045   {\r
2046     fMChandler = (AliMCEventHandler*) ((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler()); \r
2047   }\r
2048   \r
2049   // loop over events and fill them\r
2050   // Multiplicity information needed by the header (to be revised!)\r
2051   Int_t nTracks    = esd->GetNumberOfTracks();\r
2052   for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) esd->GetTrack(iTrack)->SetESDEvent(esd);\r
2053 \r
2054   // Update the header\r
2055 \r
2056   Int_t nV0s      = esd->GetNumberOfV0s();\r
2057   Int_t nCascades = esd->GetNumberOfCascades();\r
2058   Int_t nKinks    = esd->GetNumberOfKinks();\r
2059   Int_t nVertices = nV0s + nCascades /*V0 wihtin cascade already counted*/+ nKinks + 1 /* = prim. vtx*/;\r
2060   Int_t nPileSPDVertices=1+esd->GetNumberOfPileupVerticesSPD(); // also SPD main vertex\r
2061   Int_t nPileTrkVertices=esd->GetNumberOfPileupVerticesTracks();\r
2062   nVertices+=nPileSPDVertices;\r
2063   nVertices+=nPileTrkVertices;\r
2064   Int_t nJets     = 0;\r
2065   Int_t nCaloClus = esd->GetNumberOfCaloClusters();\r
2066   Int_t nFmdClus  = 0;\r
2067   Int_t nPmdClus  = esd->GetNumberOfPmdTracks();\r
2068     \r
2069   AliDebug(1,Form("   NV0=%d  NCASCADES=%d  NKINKS=%d", nV0s, nCascades, nKinks));\r
2070        \r
2071   AODEvent()->ResetStd(nTracks, nVertices, nV0s, nCascades, nJets, nCaloClus, nFmdClus, nPmdClus);\r
2072 \r
2073   if (nV0s > 0) \r
2074   {\r
2075     // RefArray to store a mapping between esd V0 number and newly created AOD-Vertex V0\r
2076     fAODV0VtxRefs = new TRefArray(nV0s);\r
2077     // RefArray to store the mapping between esd V0 number and newly created AOD-V0\r
2078     fAODV0Refs = new TRefArray(nV0s); \r
2079     // Array to take into account the V0s already added to the AOD (V0 within cascades)\r
2080     fUsedV0 = new Bool_t[nV0s];\r
2081     for (Int_t iV0=0; iV0<nV0s; ++iV0) fUsedV0[iV0]=kFALSE;\r
2082   }\r
2083   \r
2084   if (nTracks>0) \r
2085   {\r
2086     // RefArray to store the mapping between esd track number and newly created AOD-Track\r
2087     \r
2088     fAODTrackRefs = new TRefArray(nTracks);\r
2089 \r
2090     // Array to take into account the tracks already added to the AOD    \r
2091     fUsedTrack = new Bool_t[nTracks];\r
2092     for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) fUsedTrack[iTrack]=kFALSE;\r
2093   }\r
2094   \r
2095   // Array to take into account the kinks already added to the AOD\r
2096   if (nKinks>0) \r
2097   {\r
2098     fUsedKink = new Bool_t[nKinks];\r
2099     for (Int_t iKink=0; iKink<nKinks; ++iKink) fUsedKink[iKink]=kFALSE;\r
2100   }\r
2101     \r
2102   ConvertPrimaryVertices(*esd);\r
2103 \r
2104   //setting best TOF PID\r
2105   AliESDInputHandler* esdH = dynamic_cast<AliESDInputHandler*>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());\r
2106   if (esdH)\r
2107       fESDpid = esdH->GetESDpid();\r
2108 \r
2109   if (fIsPidOwner && fESDpid){\r
2110     delete fESDpid;\r
2111     fESDpid = 0;\r
2112   }\r
2113   if(!fESDpid)\r
2114   { //in case of no Tender attached \r
2115     fESDpid = new AliESDpid;\r
2116     fIsPidOwner = kTRUE;\r
2117   }\r
2118   \r
2119   if(!esd->GetTOFHeader())\r
2120   { //protection in case the pass2 LHC10b,c,d have been processed without tender. \r
2121     Float_t t0spread[10];\r
2122     Float_t intrinsicTOFres=100; //ps ok for LHC10b,c,d pass2!! \r
2123     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
2124     fESDpid->GetTOFResponse().SetT0resolution(t0spread);\r
2125     fESDpid->GetTOFResponse().SetTimeResolution(intrinsicTOFres);\r
2126     fESDpid->SetTOFResponse(esd, (AliESDpid::EStartTimeType_t)fTimeZeroType);    \r
2127     AliTOFHeader tmpTOFHeader(0,t0spread[0],0,NULL,NULL,NULL,intrinsicTOFres,t0spread[0]);   \r
2128     AODEvent()->SetTOFHeader(&tmpTOFHeader);         // write dummy TOF header in AOD\r
2129   } else {\r
2130     AODEvent()->SetTOFHeader(esd->GetTOFHeader());    // write TOF header in AOD\r
2131   }\r
2132   \r
2133   if(esd->GetTOFHeader() && fIsPidOwner) fESDpid->SetTOFResponse(esd, (AliESDpid::EStartTimeType_t)fTimeZeroType); //in case of AOD production strating form LHC10e without Tender. \r
2134   \r
2135   if ( fAreCascadesEnabled ) ConvertCascades(*esd);\r
2136 \r
2137   if ( fAreV0sEnabled ) ConvertV0s(*esd);\r
2138   \r
2139   if ( fAreKinksEnabled ) ConvertKinks(*esd);\r
2140   \r
2141   if ( fAreTracksEnabled ) ConvertTracks(*esd);\r
2142   \r
2143   // Update number of AOD tracks in header at the end of track loop (M.G.)\r
2144   header->SetRefMultiplicity(fNumberOfTracks);\r
2145   header->SetRefMultiplicityPos(fNumberOfPositiveTracks);\r
2146   header->SetRefMultiplicityNeg(fNumberOfTracks - fNumberOfPositiveTracks);\r
2147 \r
2148   if ( fTPCConstrainedFilterMask ) ConvertTPCOnlyTracks(*esd);\r
2149   if( fGlobalConstrainedFilterMask) ConvertGlobalConstrainedTracks(*esd);  \r
2150 \r
2151   if ( fArePmdClustersEnabled ) ConvertPmdClusters(*esd);\r
2152   \r
2153   if ( fAreCaloClustersEnabled ) ConvertCaloClusters(*esd);\r
2154   \r
2155   if ( fAreEMCALCellsEnabled )ConvertEMCALCells(*esd);\r
2156   \r
2157   if ( fArePHOSCellsEnabled )ConvertPHOSCells(*esd);\r
2158         \r
2159         if ( fAreEMCALTriggerEnabled )ConvertCaloTrigger(TString("EMCAL"), *esd);\r
2160 \r
2161         if ( fArePHOSTriggerEnabled )ConvertCaloTrigger(TString("PHOS"), *esd);\r
2162   \r
2163   if ( fAreTrackletsEnabled ) ConvertTracklets(*esd);\r
2164   if ( fIsZDCEnabled ) ConvertZDC(*esd);\r
2165   \r
2166   delete fAODTrackRefs; fAODTrackRefs=0x0;\r
2167   delete fAODV0VtxRefs; fAODV0VtxRefs=0x0;\r
2168   delete fAODV0Refs; fAODV0Refs=0x0;\r
2169   \r
2170   delete[] fUsedTrack; fUsedTrack=0x0;\r
2171   delete[] fUsedV0; fUsedV0=0x0;\r
2172   delete[] fUsedKink; fUsedKink=0x0;\r
2173 \r
2174   if ( fIsPidOwner){\r
2175     delete fESDpid;\r
2176     fESDpid = 0x0;\r
2177   }\r
2178 \r
2179 \r
2180 }\r
2181 \r
2182 \r
2183 //______________________________________________________________________________\r
2184 void AliAnalysisTaskESDfilter::SetAODPID(AliESDtrack *esdtrack, AliAODTrack *aodtrack, AliAODPid *detpid)\r
2185 {\r
2186   //\r
2187   // Setter for the raw PID detector signals\r
2188   //\r
2189 \r
2190   // Save PID object for candidate electrons\r
2191     Bool_t pidSave = kFALSE;\r
2192     if (fTrackFilter) {\r
2193         Bool_t selectInfo = fTrackFilter->IsSelected((char*) "Electrons");\r
2194         if (selectInfo)  pidSave = kTRUE;\r
2195     }\r
2196 \r
2197 \r
2198     // Tracks passing pt cut \r
2199     if(esdtrack->Pt()>fHighPthreshold) {\r
2200         pidSave = kTRUE;\r
2201     } else {\r
2202         if(fPtshape){\r
2203             if(esdtrack->Pt()> fPtshape->GetXmin()){\r
2204                 Double_t y = fPtshape->Eval(esdtrack->Pt())/fPtshape->Eval(fHighPthreshold);\r
2205                 if(gRandom->Rndm(0)<1./y){\r
2206                     pidSave = kTRUE;\r
2207                 }//end rndm\r
2208             }//end if p < pmin\r
2209         }//end if p function\r
2210     }// end else\r
2211 \r
2212     if (pidSave) {\r
2213       if(!aodtrack->GetDetPid()){// prevent memory leak when calling SetAODPID twice for the same track\r
2214         detpid = new AliAODPid();\r
2215         SetDetectorRawSignals(detpid,esdtrack);\r
2216         aodtrack->SetDetPID(detpid);\r
2217       }\r
2218     }\r
2219 }\r
2220 \r
2221 //______________________________________________________________________________\r
2222 void AliAnalysisTaskESDfilter::SetDetectorRawSignals(AliAODPid *aodpid, AliESDtrack *track)\r
2223 {\r
2224 //\r
2225 //assignment of the detector signals (AliXXXesdPID inspired)\r
2226 //\r
2227  if(!track){\r
2228  AliInfo("no ESD track found. .....exiting");\r
2229  return;\r
2230  }\r
2231  // TPC momentum\r
2232  const AliExternalTrackParam *in=track->GetInnerParam();\r
2233  if (in) {\r
2234    aodpid->SetTPCmomentum(in->GetP());\r
2235  }else{\r
2236    aodpid->SetTPCmomentum(-1.);\r
2237  }\r
2238 \r
2239 \r
2240  aodpid->SetITSsignal(track->GetITSsignal());\r
2241  Double_t itsdedx[4]; // dE/dx samples for individual ITS layers\r
2242  track->GetITSdEdxSamples(itsdedx);\r
2243  aodpid->SetITSdEdxSamples(itsdedx);\r
2244 \r
2245  aodpid->SetTPCsignal(track->GetTPCsignal());\r
2246  aodpid->SetTPCsignalN(track->GetTPCsignalN());\r
2247  if(track->GetTPCdEdxInfo()) aodpid->SetTPCdEdxInfo(track->GetTPCdEdxInfo());\r
2248 \r
2249  //n TRD planes = 6\r
2250  Int_t nslices = track->GetNumberOfTRDslices()*6;\r
2251  TArrayD trdslices(nslices);\r
2252  for(Int_t iSl =0; iSl < track->GetNumberOfTRDslices(); iSl++) {\r
2253    for(Int_t iPl =0; iPl<6; iPl++) trdslices[iPl*track->GetNumberOfTRDslices()+iSl] = track->GetTRDslice(iPl,iSl);\r
2254  }\r
2255  \r
2256 //TRD momentum\r
2257  for(Int_t iPl=0;iPl<6;iPl++){\r
2258    Double_t trdmom=track->GetTRDmomentum(iPl);\r
2259    aodpid->SetTRDmomentum(iPl,trdmom);\r
2260  }\r
2261 \r
2262  aodpid->SetTRDsignal(track->GetNumberOfTRDslices()*6,trdslices.GetArray());\r
2263 \r
2264  //TRD clusters and tracklets\r
2265  aodpid->SetTRDncls(track->GetTRDncls());\r
2266  aodpid->SetTRDntrackletsPID(track->GetTRDntrackletsPID());\r
2267  \r
2268  //TOF PID  \r
2269  Double_t times[AliAODPid::kSPECIES]; track->GetIntegratedTimes(times);\r
2270  aodpid->SetIntegratedTimes(times);\r
2271 \r
2272   Float_t tzeroTrack = fESDpid->GetTOFResponse().GetStartTime(track->P());\r
2273   aodpid->SetTOFsignal(track->GetTOFsignal()-tzeroTrack);\r
2274   \r
2275   Double_t tofRes[5];\r
2276   for (Int_t iMass=0; iMass<5; iMass++){\r
2277     tofRes[iMass]=(Double_t)fESDpid->GetTOFResponse().GetExpectedSigma(track->P(), times[iMass], AliPID::ParticleMass(iMass));\r
2278   }\r
2279   aodpid->SetTOFpidResolution(tofRes);\r
2280 \r
2281   aodpid->SetHMPIDsignal(track->GetHMPIDsignal());\r
2282 \r
2283  //Extrapolate track to EMCAL surface for AOD-level track-cluster matching\r
2284  Double_t emcpos[3] = {0.,0.,0.};\r
2285  Double_t emcmom[3] = {0.,0.,0.};\r
2286  aodpid->SetEMCALPosition(emcpos);\r
2287  aodpid->SetEMCALMomentum(emcmom);\r
2288 \r
2289  Double_t hmpPid[5] = {0};\r
2290  track->GetHMPIDpid(hmpPid);\r
2291  aodpid->SetHMPIDprobs(hmpPid);\r
2292 \r
2293  AliExternalTrackParam *outerparam = (AliExternalTrackParam*)track->GetOuterParam();\r
2294  if(!outerparam) return;\r
2295 \r
2296  //To be replaced by call to AliEMCALGeoUtils when the class becomes available\r
2297  Bool_t okpos = outerparam->GetXYZ(emcpos);\r
2298  Bool_t okmom = outerparam->GetPxPyPz(emcmom);\r
2299  if(!(okpos && okmom)) return;\r
2300 \r
2301  aodpid->SetEMCALPosition(emcpos);\r
2302  aodpid->SetEMCALMomentum(emcmom);\r
2303 \r
2304 }\r
2305 \r
2306 Double_t  AliAnalysisTaskESDfilter::Chi2perNDF(AliESDtrack* track)\r
2307 {\r
2308     // Calculate chi2 per ndf for track\r
2309     Int_t  nClustersTPC = track->GetTPCNcls();\r
2310 \r
2311     if ( nClustersTPC > 5) {\r
2312        return (track->GetTPCchi2()/Float_t(nClustersTPC - 5));\r
2313     } else {\r
2314        return (-1.);\r
2315     }\r
2316  }\r
2317 \r
2318 \r
2319 //______________________________________________________________________________\r
2320 void AliAnalysisTaskESDfilter::Terminate(Option_t */*option*/)\r
2321 {\r
2322 // Terminate analysis\r
2323 //\r
2324     if (fDebug > 1) printf("AnalysisESDfilter: Terminate() \n");\r
2325 }\r
2326 \r
2327 //______________________________________________________________________________\r
2328 void  AliAnalysisTaskESDfilter::PrintMCInfo(AliStack *pStack,Int_t label){\r
2329 // Print MC info\r
2330   if(!pStack)return;\r
2331   label = TMath::Abs(label);\r
2332   TParticle *part = pStack->Particle(label);\r
2333   Printf("########################");\r
2334   Printf("%s:%d %d UniqueID %d PDG %d P %3.3f",(char*)__FILE__,__LINE__,label,part->GetUniqueID(),part->GetPdgCode(),part->P());\r
2335   part->Print();\r
2336   TParticle* mother = part;\r
2337   Int_t imo = part->GetFirstMother();\r
2338   Int_t nprim = pStack->GetNprimary();\r
2339   //  while((imo >= nprim) && (mother->GetUniqueID() == 4)) {\r
2340   while((imo >= nprim)) {\r
2341     mother =  pStack->Particle(imo);\r
2342     Printf("Mother %s:%d Label %d UniqueID %d PDG %d P %3.3f",(char*)__FILE__,__LINE__,imo,mother->GetUniqueID(),mother->GetPdgCode(),mother->P());\r
2343     mother->Print();\r
2344     imo =  mother->GetFirstMother();\r
2345   }\r
2346   Printf("########################");\r
2347 }\r
2348 \r
2349 //______________________________________________________\r
2350 \r