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