1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 #include <Riostream.h>
23 #include <TParameter.h>
25 #include <TParticle.h>
29 #include "AliAnalysisTaskESDfilter.h"
30 #include "AliAnalysisManager.h"
31 #include "AliESDEvent.h"
32 #include "AliESDRun.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"
46 #include "AliESDkink.h"
47 #include "AliESDcascade.h"
48 #include "AliESDPmdTrack.h"
49 #include "AliESDCaloCluster.h"
50 #include "AliESDCaloCells.h"
51 #include "AliMultiplicity.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"
63 #include "AliESDTrdTrack.h"
64 #include "AliESDTrdTracklet.h"
65 #include "AliAODTrdTrack.h"
66 #include "AliAODTrdTracklet.h"
67 #include "AliEMCALRecoUtils.h"
68 #include "AliESDUtils.h"
72 ClassImp(AliAnalysisTaskESDfilter)
74 ////////////////////////////////////////////////////////////////////////
76 AliAnalysisTaskESDfilter::AliAnalysisTaskESDfilter():
84 fEnableFillAOD(kTRUE),
94 fNumberOfPositiveTracks(0),
99 fOldESDformat(kFALSE),
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),
126 fTPCaloneTrackCuts(0),
127 fDoPropagateTrackToEMCal(kTRUE),
128 fEMCalSurfaceDistance(440),
129 fRefitVertexTracks(-1),
130 fRefitVertexTracksNCuts(0),
131 fRefitVertexTracksCuts(0)
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
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
152 //______________________________________________________________________________
153 AliAnalysisTaskESDfilter::AliAnalysisTaskESDfilter(const char* name):
154 AliAnalysisTaskSE(name),
161 fEnableFillAOD(kTRUE),
171 fNumberOfPositiveTracks(0),
173 fNumberOfVertices(0),
174 fNumberOfCascades(0),
176 fOldESDformat(kFALSE),
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),
203 fTPCaloneTrackCuts(0),
204 fDoPropagateTrackToEMCal(kTRUE),
205 fEMCalSurfaceDistance(440),
206 fRefitVertexTracks(-1),
207 fRefitVertexTracksNCuts(0),
208 fRefitVertexTracksCuts(0)
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
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
230 AliAnalysisTaskESDfilter::~AliAnalysisTaskESDfilter()
232 if(fIsPidOwner) delete fESDpid;
233 delete[] fRefitVertexTracksCuts;
236 //______________________________________________________________________________
237 void AliAnalysisTaskESDfilter::UserCreateOutputObjects()
240 // Create Output Objects conenct filter to outputtree
244 OutputTree()->GetUserInfo()->Add(fTrackFilter);
248 AliError("No OutputTree() for adding the track filter");
250 fTPCaloneTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
253 //______________________________________________________________________________
254 void AliAnalysisTaskESDfilter::Init()
257 if (fDebug > 1) AliInfo("Init() \n");
260 //______________________________________________________________________________
261 Bool_t AliAnalysisTaskESDfilter::Notify()
264 AddMetadataToUserInfo();
268 //______________________________________________________________________________
269 Bool_t AliAnalysisTaskESDfilter::AddMetadataToUserInfo()
271 // Copy metadata to AOD user info.
272 static Bool_t copyFirst = kFALSE;
274 AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
276 AliError("AliAnalysisTaskESDfilter::AddMetadataToUserInfo() : No analysis manager !");
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));
293 //______________________________________________________________________________
294 void AliAnalysisTaskESDfilter::PrintTask(Option_t *option, Int_t indent) const
296 // Print selection task information
299 AliAnalysisTaskSE::PrintTask(option,indent);
301 TString spaces(' ',indent+3);
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;
318 //______________________________________________________________________________
319 void AliAnalysisTaskESDfilter::UserExec(Option_t */*option*/)
321 // Execute analysis for current event
323 Long64_t ientry = Entry();
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");
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");
333 AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillAOD(kTRUE);
334 AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillExtension(kTRUE);
339 //______________________________________________________________________________
340 TClonesArray& AliAnalysisTaskESDfilter::Cascades()
342 return *(AODEvent()->GetCascades());
345 //______________________________________________________________________________
346 TClonesArray& AliAnalysisTaskESDfilter::Tracks()
348 return *(AODEvent()->GetTracks());
351 //______________________________________________________________________________
352 TClonesArray& AliAnalysisTaskESDfilter::V0s()
354 return *(AODEvent()->GetV0s());
357 //______________________________________________________________________________
358 TClonesArray& AliAnalysisTaskESDfilter::Vertices()
360 return *(AODEvent()->GetVertices());
363 //______________________________________________________________________________
364 AliAODHeader* AliAnalysisTaskESDfilter::ConvertHeader(const AliESDEvent& esd)
366 // Convert header information
368 AliCodeTimerAuto("",0);
370 AliAODHeader* header = dynamic_cast<AliAODHeader*>(AODEvent()->GetHeader());
371 if(!header) AliFatal("Not a standard AOD");
373 header->SetRunNumber(esd.GetRunNumber());
374 header->SetOfflineTrigger(fInputHandler->IsEventSelected()); // propagate the decision of the physics selection
375 header->SetNumberOfESDTracks(esd.GetNumberOfTracks());
377 TTree* tree = fInputHandler->GetTree();
379 TFile* file = tree->GetCurrentFile();
380 if (file) header->SetESDFileName(file->GetName());
384 header->SetBunchCrossNumber(0);
385 header->SetOrbitNumber(0);
386 header->SetPeriodNumber(0);
387 header->SetEventType(0);
388 header->SetMuonMagFieldScale(-999.);
389 header->SetCentrality(0);
390 header->SetEventplane(0);
392 header->SetBunchCrossNumber(esd.GetBunchCrossNumber());
393 header->SetOrbitNumber(esd.GetOrbitNumber());
394 header->SetPeriodNumber(esd.GetPeriodNumber());
395 header->SetEventType(esd.GetEventType());
397 header->SetEventNumberESDFile(esd.GetHeader()->GetEventNumberInFile());
398 if(const_cast<AliESDEvent&>(esd).GetCentrality()){
399 header->SetCentrality(const_cast<AliESDEvent&>(esd).GetCentrality());
402 header->SetCentrality(0);
404 if(const_cast<AliESDEvent&>(esd).GetEventplane()){
405 header->SetEventplane(const_cast<AliESDEvent&>(esd).GetEventplane());
408 header->SetEventplane(0);
413 header->SetFiredTriggerClasses(esd.GetFiredTriggerClasses());
414 header->SetTriggerMask(esd.GetTriggerMask());
415 header->SetTriggerCluster(esd.GetTriggerCluster());
416 header->SetL0TriggerInputs(esd.GetHeader()->GetL0TriggerInputs());
417 header->SetL1TriggerInputs(esd.GetHeader()->GetL1TriggerInputs());
418 header->SetL2TriggerInputs(esd.GetHeader()->GetL2TriggerInputs());
420 header->SetMagneticField(esd.GetMagneticField());
421 header->SetMuonMagFieldScale(esd.GetCurrentDip()/6000.);
422 header->SetZDCN1Energy(esd.GetZDCN1Energy());
423 header->SetZDCP1Energy(esd.GetZDCP1Energy());
424 header->SetZDCN2Energy(esd.GetZDCN2Energy());
425 header->SetZDCP2Energy(esd.GetZDCP2Energy());
426 header->SetZDCEMEnergy(esd.GetZDCEMEnergy(0),esd.GetZDCEMEnergy(1));
428 header->SetIRInt2InteractionMap(esd.GetHeader()->GetIRInt2InteractionMap());
429 header->SetIRInt1InteractionMap(esd.GetHeader()->GetIRInt1InteractionMap());
431 // ITS Cluster Multiplicty
432 const AliMultiplicity *mult = esd.GetMultiplicity();
433 for (Int_t ilay = 0; ilay < 6; ilay++) header->SetITSClusters(ilay, mult->GetNumberOfITSClusters(ilay));
435 // TPC only Reference Multiplicty
436 Int_t refMult = fTPCaloneTrackCuts ? (Short_t)fTPCaloneTrackCuts->GetReferenceMultiplicity(&esd, kTRUE) : -1;
437 header->SetTPConlyRefMultiplicity(refMult);
439 AliESDtrackCuts::MultEstTrackType estType = esd.GetPrimaryVertexTracks()->GetStatus() ? AliESDtrackCuts::kTrackletsITSTPC : AliESDtrackCuts::kTracklets;
440 header->SetRefMultiplicityComb05(AliESDtrackCuts::GetReferenceMultiplicity(&esd,estType,0.5));
441 header->SetRefMultiplicityComb08(AliESDtrackCuts::GetReferenceMultiplicity(&esd,estType,0.8));
443 Float_t diamxy[2]={(Float_t)esd.GetDiamondX(),(Float_t)esd.GetDiamondY()};
445 esd.GetDiamondCovXY(diamcov);
446 header->SetDiamond(diamxy,diamcov);
447 header->SetDiamondZ(esd.GetDiamondZ(),esd.GetSigma2DiamondZ());
449 // VZERO channel equalization factors for event-plane reconstruction
450 header->SetVZEROEqFactors(esd.GetVZEROEqFactors());
452 // T0 Resolution information
453 const AliESDRun* esdRun = esd.GetESDRun();
454 for (Int_t i=0;i<AliESDRun::kT0spreadSize;i++) header->SetT0spread(i,esdRun->GetT0spread(i));
459 //______________________________________________________________________________
460 void AliAnalysisTaskESDfilter::ConvertCascades(const AliESDEvent& esd)
462 // Convert the cascades part of the ESD.
463 // Return the number of cascades
465 AliCodeTimerAuto("",0);
467 // Create vertices starting from the most complex objects
470 const AliESDVertex* vtx = esd.GetPrimaryVertex();
471 Double_t pos[3] = { 0. };
472 Double_t covVtx[6] = { 0. };
473 Double_t momBach[3]={0.};
474 Double_t covTr[21]={0.};
475 // Double_t pid[10]={0.};
476 AliAODPid* detpid(0x0);
477 AliAODVertex* vV0FromCascade(0x0);
478 AliAODv0* aodV0(0x0);
479 AliAODcascade* aodCascade(0x0);
480 AliAODTrack* aodTrack(0x0);
481 Double_t momPos[3]={0.};
482 Double_t momNeg[3] = { 0. };
483 Double_t momPosAtV0vtx[3]={0.};
484 Double_t momNegAtV0vtx[3]={0.};
485 Int_t tofLabel[3] = {0};
486 TClonesArray& verticesArray = Vertices();
487 TClonesArray& tracksArray = Tracks();
488 TClonesArray& cascadesArray = Cascades();
490 // Cascades (Modified by A.Maire - February 2009)
491 for (Int_t nCascade = 0; nCascade < esd.GetNumberOfCascades(); ++nCascade) {
495 AliESDcascade *esdCascade = esd.GetCascade(nCascade);
496 Int_t idxPosFromV0Dghter = esdCascade->GetPindex();
497 Int_t idxNegFromV0Dghter = esdCascade->GetNindex();
498 Int_t idxBachFromCascade = esdCascade->GetBindex();
500 AliESDtrack *esdCascadePos = esd.GetTrack( idxPosFromV0Dghter);
501 AliESDtrack *esdCascadeNeg = esd.GetTrack( idxNegFromV0Dghter);
502 AliESDtrack *esdCascadeBach = esd.GetTrack( idxBachFromCascade);
504 // Identification of the V0 within the esdCascade (via both daughter track indices)
505 AliESDv0 * currentV0 = 0x0;
506 Int_t idxV0FromCascade = -1;
508 for (Int_t iV0=0; iV0<esd.GetNumberOfV0s(); ++iV0) {
510 currentV0 = esd.GetV0(iV0);
511 Int_t posCurrentV0 = currentV0->GetPindex();
512 Int_t negCurrentV0 = currentV0->GetNindex();
514 if (posCurrentV0==idxPosFromV0Dghter && negCurrentV0==idxNegFromV0Dghter) {
515 idxV0FromCascade = iV0;
520 if(idxV0FromCascade < 0){
521 printf("Cascade - no matching for the V0 (index V0 = -1) ! Skip ... \n");
523 }// a priori, useless check, but safer ... in case of pb with tracks "out of bounds"
525 AliESDv0 *esdV0FromCascade = esd.GetV0(idxV0FromCascade);
527 // 1 - Cascade selection
529 // AliESDVertex *esdPrimVtx = new AliESDVertex(*(esd.GetPrimaryVertex()));
530 // TList cascadeObjects;
531 // cascadeObjects.AddAt(esdV0FromCascade, 0);
532 // cascadeObjects.AddAt(esdCascadePos, 1);
533 // cascadeObjects.AddAt(esdCascadeNeg, 2);
534 // cascadeObjects.AddAt(esdCascade, 3);
535 // cascadeObjects.AddAt(esdCascadeBach, 4);
536 // cascadeObjects.AddAt(esdPrimVtx, 5);
538 // UInt_t selectCascade = 0;
539 // if (fCascadeFilter) {
540 // // selectCascade = fCascadeFilter->IsSelected(&cascadeObjects);
541 // // FIXME AliESDCascadeCuts to be implemented ...
543 // // Here we may encounter a moot point at the V0 level
544 // // between the cascade selections and the V0 ones :
545 // // the V0 selected along with the cascade (secondary V0) may
546 // // usually be removed from the dedicated V0 selections (prim V0) ...
547 // // -> To be discussed !
549 // // this is a little awkward but otherwise the
550 // // list wants to access the pointer (delete it)
551 // // again when going out of scope
552 // delete cascadeObjects.RemoveAt(5); // esdPrimVtx created via copy construct
554 // if (!selectCascade)
558 // delete cascadeObjects.RemoveAt(5); // esdPrimVtx created via copy construct
562 // 2 - Add the cascade vertex
564 esdCascade->GetXYZcascade(pos[0], pos[1], pos[2]);
565 esdCascade->GetPosCovXi(covVtx);
566 chi2 = esdCascade->GetChi2Xi();
568 AliAODVertex *vCascade = new(verticesArray[fNumberOfVertices++]) AliAODVertex(pos,
570 chi2, // FIXME = Chi2/NDF will be needed
573 AliAODVertex::kCascade);
574 fPrimaryVertex->AddDaughter(vCascade);
576 // 3 - Add the bachelor track from the cascade
578 if (!fUsedTrack[idxBachFromCascade]) {
580 esdCascadeBach->GetPxPyPz(momBach);
581 esdCascadeBach->GetXYZ(pos);
582 esdCascadeBach->GetCovarianceXYZPxPyPz(covTr);
583 // esdCascadeBach->GetESDpid(pid);
584 esdCascadeBach->GetTOFLabel(tofLabel);
586 fUsedTrack[idxBachFromCascade] = kTRUE;
587 UInt_t selectInfo = 0;
588 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadeBach);
589 if (fMChandler) fMChandler->SelectParticle(esdCascadeBach->GetLabel());
590 aodTrack = new(tracksArray[fNumberOfTracks++]) AliAODTrack(esdCascadeBach->GetID(),
591 esdCascadeBach->GetLabel(),
595 kFALSE, // Why kFALSE for "isDCA" ? FIXME
597 (Short_t)esdCascadeBach->GetSign(),
598 esdCascadeBach->GetITSClusterMap(),
601 kTRUE, // usedForVtxFit = kFALSE ? FIXME
602 vtx->UsesTrack(esdCascadeBach->GetID()),
603 AliAODTrack::kFromDecayVtx,
605 aodTrack->SetPIDForTracking(esdCascadeBach->GetPIDForTracking());
606 aodTrack->SetTPCFitMap(esdCascadeBach->GetTPCFitMap());
607 aodTrack->SetTPCClusterMap(esdCascadeBach->GetTPCClusterMap());
608 aodTrack->SetTPCSharedMap (esdCascadeBach->GetTPCSharedMap());
609 aodTrack->SetChi2perNDF(Chi2perNDF(esdCascadeBach));
610 aodTrack->SetTPCPointsF(esdCascadeBach->GetTPCNclsF());
611 aodTrack->SetTPCNCrossedRows(UShort_t(esdCascadeBach->GetTPCCrossedRows()));
612 aodTrack->SetIntegratedLength(esdCascadeBach->GetIntegratedLength());
613 aodTrack->SetTOFLabel(tofLabel);
614 CopyCaloProps(esdCascadeBach,aodTrack);
615 fAODTrackRefs->AddAt(aodTrack,idxBachFromCascade);
617 if (esdCascadeBach->GetSign() > 0) ++fNumberOfPositiveTracks;
618 aodTrack->ConvertAliPIDtoAODPID();
619 aodTrack->SetFlags(esdCascadeBach->GetStatus());
620 SetAODPID(esdCascadeBach,aodTrack,detpid);
623 aodTrack = static_cast<AliAODTrack*>( fAODTrackRefs->At(idxBachFromCascade) );
626 vCascade->AddDaughter(aodTrack);
628 // 4 - Add the V0 from the cascade.
629 // = V0vtx + both pos and neg daughter tracks + the aodV0 itself
632 if ( !fUsedV0[idxV0FromCascade] ) {
633 // 4.A - if VO structure hasn't been created yet
635 // 4.A.1 - Create the V0 vertex of the cascade
636 esdV0FromCascade->GetXYZ(pos[0], pos[1], pos[2]);
637 esdV0FromCascade->GetPosCov(covVtx);
638 chi2 = esdV0FromCascade->GetChi2V0(); // = chi2/NDF since NDF = 2*2-3 ?
640 vV0FromCascade = new(verticesArray[fNumberOfVertices++]) AliAODVertex(pos,
644 idxV0FromCascade, //id of ESDv0
647 // one V0 can be used by several cascades.
648 // So, one AOD V0 vtx can have several parent vtx.
649 // This is not directly allowed by AliAODvertex.
650 // Setting the parent vtx (here = param "vCascade") doesn't lead to a crash
651 // but to a problem of consistency within AODEvent.
652 // -> See below paragraph 4.B, for the proposed treatment of such a case.
654 // Add the vV0FromCascade to the aodVOVtxRefs
655 fAODV0VtxRefs->AddAt(vV0FromCascade,idxV0FromCascade);
657 // 4.A.2 - Add the positive tracks from the V0
659 esdCascadePos->GetPxPyPz(momPos);
660 esdCascadePos->GetXYZ(pos);
661 esdCascadePos->GetCovarianceXYZPxPyPz(covTr);
662 // esdCascadePos->GetESDpid(pid);
663 esdCascadePos->GetTOFLabel(tofLabel);
665 if (!fUsedTrack[idxPosFromV0Dghter]) {
666 fUsedTrack[idxPosFromV0Dghter] = kTRUE;
668 UInt_t selectInfo = 0;
669 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadePos);
670 if(fMChandler) fMChandler->SelectParticle(esdCascadePos->GetLabel());
671 aodTrack = new(tracksArray[fNumberOfTracks++]) AliAODTrack(esdCascadePos->GetID(),
672 esdCascadePos->GetLabel(),
676 kFALSE, // Why kFALSE for "isDCA" ? FIXME
678 (Short_t)esdCascadePos->GetSign(),
679 esdCascadePos->GetITSClusterMap(),
682 kTRUE, // usedForVtxFit = kFALSE ? FIXME
683 vtx->UsesTrack(esdCascadePos->GetID()),
684 AliAODTrack::kFromDecayVtx,
686 aodTrack->SetPIDForTracking(esdCascadePos->GetPIDForTracking());
687 aodTrack->SetTPCFitMap(esdCascadePos->GetTPCFitMap());
688 aodTrack->SetTPCClusterMap(esdCascadePos->GetTPCClusterMap());
689 aodTrack->SetTPCSharedMap (esdCascadePos->GetTPCSharedMap());
690 aodTrack->SetChi2perNDF(Chi2perNDF(esdCascadePos));
691 aodTrack->SetTPCPointsF(esdCascadePos->GetTPCNclsF());
692 aodTrack->SetTPCNCrossedRows(UShort_t(esdCascadePos->GetTPCCrossedRows()));
693 aodTrack->SetIntegratedLength(esdCascadePos->GetIntegratedLength());
694 aodTrack->SetTOFLabel(tofLabel);
695 CopyCaloProps(esdCascadePos,aodTrack);
696 fAODTrackRefs->AddAt(aodTrack,idxPosFromV0Dghter);
698 if (esdCascadePos->GetSign() > 0) ++fNumberOfPositiveTracks;
699 aodTrack->ConvertAliPIDtoAODPID();
700 aodTrack->SetFlags(esdCascadePos->GetStatus());
701 SetAODPID(esdCascadePos,aodTrack,detpid);
704 aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(idxPosFromV0Dghter));
706 vV0FromCascade->AddDaughter(aodTrack);
708 // 4.A.3 - Add the negative tracks from the V0
710 esdCascadeNeg->GetPxPyPz(momNeg);
711 esdCascadeNeg->GetXYZ(pos);
712 esdCascadeNeg->GetCovarianceXYZPxPyPz(covTr);
713 // esdCascadeNeg->GetESDpid(pid);
714 esdCascadeNeg->GetTOFLabel(tofLabel);
716 if (!fUsedTrack[idxNegFromV0Dghter]) {
717 fUsedTrack[idxNegFromV0Dghter] = kTRUE;
719 UInt_t selectInfo = 0;
720 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadeNeg);
722 fMChandler->SelectParticle(esdCascadeNeg->GetLabel());
723 aodTrack = new(tracksArray[fNumberOfTracks++]) AliAODTrack(esdCascadeNeg->GetID(),
724 esdCascadeNeg->GetLabel(),
728 kFALSE, // Why kFALSE for "isDCA" ? FIXME
730 (Short_t)esdCascadeNeg->GetSign(),
731 esdCascadeNeg->GetITSClusterMap(),
734 kTRUE, // usedForVtxFit = kFALSE ? FIXME
735 vtx->UsesTrack(esdCascadeNeg->GetID()),
736 AliAODTrack::kFromDecayVtx,
738 aodTrack->SetPIDForTracking(esdCascadeNeg->GetPIDForTracking());
739 aodTrack->SetTPCFitMap(esdCascadeNeg->GetTPCFitMap());
740 aodTrack->SetTPCClusterMap(esdCascadeNeg->GetTPCClusterMap());
741 aodTrack->SetTPCSharedMap (esdCascadeNeg->GetTPCSharedMap());
742 aodTrack->SetChi2perNDF(Chi2perNDF(esdCascadeNeg));
743 aodTrack->SetTPCPointsF(esdCascadeNeg->GetTPCNclsF());
744 aodTrack->SetTPCNCrossedRows(UShort_t(esdCascadeNeg->GetTPCCrossedRows()));
745 aodTrack->SetIntegratedLength(esdCascadeNeg->GetIntegratedLength());
746 aodTrack->SetTOFLabel(tofLabel);
747 CopyCaloProps(esdCascadeNeg,aodTrack);
748 fAODTrackRefs->AddAt(aodTrack,idxNegFromV0Dghter);
750 if (esdCascadeNeg->GetSign() > 0) ++fNumberOfPositiveTracks;
751 aodTrack->ConvertAliPIDtoAODPID();
752 aodTrack->SetFlags(esdCascadeNeg->GetStatus());
753 SetAODPID(esdCascadeNeg,aodTrack,detpid);
756 aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(idxNegFromV0Dghter));
759 vV0FromCascade->AddDaughter(aodTrack);
761 // 4.A.4 - Add the V0 from cascade to the V0 array
763 Double_t dcaV0Daughters = esdV0FromCascade->GetDcaV0Daughters();
764 Double_t dcaV0ToPrimVertex = esdV0FromCascade->GetD(esd.GetPrimaryVertex()->GetX(),
765 esd.GetPrimaryVertex()->GetY(),
766 esd.GetPrimaryVertex()->GetZ() );
767 esdV0FromCascade->GetPPxPyPz( momPosAtV0vtx[0],momPosAtV0vtx[1],momPosAtV0vtx[2] );
768 esdV0FromCascade->GetNPxPyPz( momNegAtV0vtx[0],momNegAtV0vtx[1],momNegAtV0vtx[2] );
770 Double_t dcaDaughterToPrimVertex[2] = { 999., 999.}; // ..[0] = DCA in (x,y) for Pos and ..[1] = Neg
771 dcaDaughterToPrimVertex[0] = TMath::Abs(esdCascadePos->GetD(esd.GetPrimaryVertex()->GetX(),
772 esd.GetPrimaryVertex()->GetY(),
773 esd.GetMagneticField()) );
774 dcaDaughterToPrimVertex[1] = TMath::Abs(esdCascadeNeg->GetD(esd.GetPrimaryVertex()->GetX(),
775 esd.GetPrimaryVertex()->GetY(),
776 esd.GetMagneticField()) );
778 aodV0 = new(V0s()[fNumberOfV0s++]) AliAODv0(vV0FromCascade,
783 dcaDaughterToPrimVertex);
784 // set the aod v0 on-the-fly status
785 aodV0->SetOnFlyStatus(esdV0FromCascade->GetOnFlyStatus());
787 // Add the aodV0 to the aodVORefs
788 fAODV0Refs->AddAt(aodV0,idxV0FromCascade);
790 fUsedV0[idxV0FromCascade] = kTRUE;
793 // 4.B - if V0 structure already used
796 // one V0 can be used by several cascades (frequent in PbPb evts) :
797 // same V0 which used but attached to different bachelor tracks
798 // -> aodVORefs and fAODV0VtxRefs are needed.
799 // Goal : avoid a redundancy of the info in "Vertices" and "v0s" clones array.
801 vV0FromCascade = static_cast<AliAODVertex*>(fAODV0VtxRefs->At(idxV0FromCascade));
802 aodV0 = static_cast<AliAODv0*> (fAODV0Refs ->At(idxV0FromCascade));
804 // - Treatment of the parent for such a "re-used" V0 :
805 // Insert the cascade that reuses the V0 vertex in the lineage chain
806 // Before : vV0 -> vCascade1 -> vPrimary
807 // - Hyp : cascade2 uses the same V0 as cascade1
808 // After : vV0 -> vCascade2 -> vCascade1 -> vPrimary
810 AliAODVertex *vCascadePreviousParent = static_cast<AliAODVertex*> (vV0FromCascade->GetParent());
811 vV0FromCascade->SetParent(vCascade);
812 vCascade ->SetParent(vCascadePreviousParent);
814 }// end if V0 structure already used
816 // In any case (used V0 or not), add the V0 vertex to the cascade one.
817 vCascade->AddDaughter(vV0FromCascade);
819 // 5 - Add the primary track of the cascade (if any)
821 // 6 - Add the cascade to the AOD array of cascades
822 Double_t dcaBachToPrimVertexXY = TMath::Abs(esdCascadeBach->GetD(esd.GetPrimaryVertex()->GetX(),
823 esd.GetPrimaryVertex()->GetY(),
824 esd.GetMagneticField()) );
826 Double_t momBachAtCascadeVtx[3]={0.};
828 esdCascade->GetBPxPyPz(momBachAtCascadeVtx[0], momBachAtCascadeVtx[1], momBachAtCascadeVtx[2]);
830 aodCascade = new(cascadesArray[fNumberOfCascades++]) AliAODcascade(vCascade,
831 esdCascade->Charge(),
832 esdCascade->GetDcaXiDaughters(),
834 // DCAXiToPrimVtx -> needs to be calculated ----|
835 // doesn't exist at ESD level;
836 // See AODcascade::DcaXiToPrimVertex(Double, Double, Double)
837 dcaBachToPrimVertexXY,
841 printf("---- Cascade / AOD cascade : \n\n");
842 aodCascade->PrintXi(fPrimaryVertex->GetX(), fPrimaryVertex->GetY(), fPrimaryVertex->GetZ());
845 } // end of the loop on cascades
847 Cascades().Expand(fNumberOfCascades);
850 //______________________________________________________________________________
851 void AliAnalysisTaskESDfilter::ConvertV0s(const AliESDEvent& esd)
853 // Access to the AOD container of V0s
855 AliCodeTimerAuto("",0);
860 Double_t pos[3] = { 0. };
862 Double_t covVtx[6] = { 0. };
863 Double_t momPos[3]={0.};
864 Double_t covTr[21]={0.};
865 // Double_t pid[10]={0.};
866 AliAODTrack* aodTrack(0x0);
867 AliAODPid* detpid(0x0);
868 Double_t momNeg[3]={0.};
869 Double_t momPosAtV0vtx[3]={0.};
870 Double_t momNegAtV0vtx[3]={0.};
871 Int_t tofLabel[3] = {0};
872 for (Int_t nV0 = 0; nV0 < esd.GetNumberOfV0s(); ++nV0) {
873 if (fUsedV0[nV0]) continue; // skip if already added to the AOD
875 AliESDv0 *v0 = esd.GetV0(nV0);
876 Int_t posFromV0 = v0->GetPindex();
877 Int_t negFromV0 = v0->GetNindex();
881 AliESDVertex *esdVtx = new AliESDVertex(*(esd.GetPrimaryVertex()));
882 AliESDtrack *esdV0Pos = esd.GetTrack(posFromV0);
883 AliESDtrack *esdV0Neg = esd.GetTrack(negFromV0);
885 v0objects.AddAt(v0, 0);
886 v0objects.AddAt(esdV0Pos, 1);
887 v0objects.AddAt(esdV0Neg, 2);
888 v0objects.AddAt(esdVtx, 3);
891 selectV0 = fV0Filter->IsSelected(&v0objects);
892 // this is a little awkward but otherwise the
893 // list wants to access the pointer (delete it)
894 // again when going out of scope
895 delete v0objects.RemoveAt(3); // esdVtx created via copy construct
900 delete v0objects.RemoveAt(3); // esdVtx created via copy construct
904 v0->GetXYZ(pos[0], pos[1], pos[2]);
906 if (!fOldESDformat) {
907 chi2 = v0->GetChi2V0(); // = chi2/NDF since NDF = 2*2-3
908 v0->GetPosCov(covVtx);
911 for (Int_t i = 0; i < 6; i++) covVtx[i] = 0.;
916 new(Vertices()[fNumberOfVertices++]) AliAODVertex(pos,
922 fPrimaryVertex->AddDaughter(vV0);
925 // Add the positive tracks from the V0
927 esdV0Pos->GetPxPyPz(momPos);
928 esdV0Pos->GetXYZ(pos);
929 esdV0Pos->GetCovarianceXYZPxPyPz(covTr);
930 // esdV0Pos->GetESDpid(pid);
931 esdV0Pos->GetTOFLabel(tofLabel);
933 const AliESDVertex *vtx = esd.GetPrimaryVertex();
935 if (!fUsedTrack[posFromV0]) {
936 fUsedTrack[posFromV0] = kTRUE;
937 UInt_t selectInfo = 0;
938 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdV0Pos);
939 if(fMChandler)fMChandler->SelectParticle(esdV0Pos->GetLabel());
940 aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdV0Pos->GetID(),
941 esdV0Pos->GetLabel(),
947 (Short_t)esdV0Pos->GetSign(),
948 esdV0Pos->GetITSClusterMap(),
951 kTRUE, // check if this is right
952 vtx->UsesTrack(esdV0Pos->GetID()),
953 AliAODTrack::kFromDecayVtx,
955 aodTrack->SetPIDForTracking(esdV0Pos->GetPIDForTracking());
956 aodTrack->SetTPCFitMap(esdV0Pos->GetTPCFitMap());
957 aodTrack->SetTPCClusterMap(esdV0Pos->GetTPCClusterMap());
958 aodTrack->SetTPCSharedMap (esdV0Pos->GetTPCSharedMap());
959 aodTrack->SetChi2perNDF(Chi2perNDF(esdV0Pos));
960 aodTrack->SetTPCPointsF(esdV0Pos->GetTPCNclsF());
961 aodTrack->SetTPCNCrossedRows(UShort_t(esdV0Pos->GetTPCCrossedRows()));
962 aodTrack->SetIntegratedLength(esdV0Pos->GetIntegratedLength());
963 aodTrack->SetTOFLabel(tofLabel);
964 CopyCaloProps(esdV0Pos,aodTrack);
965 fAODTrackRefs->AddAt(aodTrack,posFromV0);
966 if (esdV0Pos->GetSign() > 0) ++fNumberOfPositiveTracks;
967 aodTrack->ConvertAliPIDtoAODPID();
968 aodTrack->SetFlags(esdV0Pos->GetStatus());
969 SetAODPID(esdV0Pos,aodTrack,detpid);
971 aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(posFromV0));
973 vV0->AddDaughter(aodTrack);
975 // Add the negative tracks from the V0
976 esdV0Neg->GetPxPyPz(momNeg);
977 esdV0Neg->GetXYZ(pos);
978 esdV0Neg->GetCovarianceXYZPxPyPz(covTr);
979 // esdV0Neg->GetESDpid(pid);
980 esdV0Neg->GetTOFLabel(tofLabel);
982 if (!fUsedTrack[negFromV0]) {
983 fUsedTrack[negFromV0] = kTRUE;
984 UInt_t selectInfo = 0;
985 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdV0Neg);
986 if(fMChandler)fMChandler->SelectParticle(esdV0Neg->GetLabel());
987 aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdV0Neg->GetID(),
988 esdV0Neg->GetLabel(),
994 (Short_t)esdV0Neg->GetSign(),
995 esdV0Neg->GetITSClusterMap(),
998 kTRUE, // check if this is right
999 vtx->UsesTrack(esdV0Neg->GetID()),
1000 AliAODTrack::kFromDecayVtx,
1002 aodTrack->SetPIDForTracking(esdV0Neg->GetPIDForTracking());
1003 aodTrack->SetTPCFitMap(esdV0Neg->GetTPCFitMap());
1004 aodTrack->SetTPCClusterMap(esdV0Neg->GetTPCClusterMap());
1005 aodTrack->SetTPCSharedMap (esdV0Neg->GetTPCSharedMap());
1006 aodTrack->SetChi2perNDF(Chi2perNDF(esdV0Neg));
1007 aodTrack->SetTPCPointsF(esdV0Neg->GetTPCNclsF());
1008 aodTrack->SetTPCNCrossedRows(UShort_t(esdV0Neg->GetTPCCrossedRows()));
1009 aodTrack->SetIntegratedLength(esdV0Neg->GetIntegratedLength());
1010 aodTrack->SetTOFLabel(tofLabel);
1011 CopyCaloProps(esdV0Neg,aodTrack);
1012 fAODTrackRefs->AddAt(aodTrack,negFromV0);
1013 if (esdV0Neg->GetSign() > 0) ++fNumberOfPositiveTracks;
1014 aodTrack->ConvertAliPIDtoAODPID();
1015 aodTrack->SetFlags(esdV0Neg->GetStatus());
1016 SetAODPID(esdV0Neg,aodTrack,detpid);
1018 aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(negFromV0));
1020 vV0->AddDaughter(aodTrack);
1022 // Add the V0 the V0 array as well
1023 Double_t dcaV0Daughters = v0->GetDcaV0Daughters();
1024 Double_t dcaV0ToPrimVertex = v0->GetD(esd.GetPrimaryVertex()->GetX(),
1025 esd.GetPrimaryVertex()->GetY(),
1026 esd.GetPrimaryVertex()->GetZ());
1028 v0->GetPPxPyPz(momPosAtV0vtx[0],momPosAtV0vtx[1],momPosAtV0vtx[2]);
1029 v0->GetNPxPyPz(momNegAtV0vtx[0],momNegAtV0vtx[1],momNegAtV0vtx[2]);
1031 Double_t dcaDaughterToPrimVertex[2] = { 999., 999.}; // ..[0] = DCA in (x,y) for Pos and ..[1] = Neg
1032 dcaDaughterToPrimVertex[0] = TMath::Abs(esdV0Pos->GetD( esd.GetPrimaryVertex()->GetX(),
1033 esd.GetPrimaryVertex()->GetY(),
1034 esd.GetMagneticField()) );
1035 dcaDaughterToPrimVertex[1] = TMath::Abs(esdV0Neg->GetD( esd.GetPrimaryVertex()->GetX(),
1036 esd.GetPrimaryVertex()->GetY(),
1037 esd.GetMagneticField()) );
1039 AliAODv0* aodV0 = new(V0s()[fNumberOfV0s++]) AliAODv0(vV0,
1044 dcaDaughterToPrimVertex);
1046 // set the aod v0 on-the-fly status
1047 aodV0->SetOnFlyStatus(v0->GetOnFlyStatus());
1048 }//End of loop on V0s
1050 V0s().Expand(fNumberOfV0s);
1053 //______________________________________________________________________________
1054 void AliAnalysisTaskESDfilter::ConvertTPCOnlyTracks(const AliESDEvent& esd)
1056 // Convert TPC only tracks
1057 // Here we have wo hybrid appraoch to remove fakes
1058 // ******* ITSTPC ********
1059 // Uses a cut on the ITS properties to select global tracks
1060 // which are than marked as HybdridITSTPC for the remainder
1061 // the TPC only tracks are flagged as HybridITSTPConly.
1062 // Note, in order not to get fakes back in the TPC cuts, one needs
1063 // two "ITS" cuts one tight (1) (to throw out fakes) and one lose (2) (to NOT flag the trakcs in the TPC only)
1064 // using cut number (3)
1065 // so fHybridFilterMask == (1)|(2) fTPCFilterMask = (3), Usercode needs to slect with mask = (1)|(3) and track->IsHybridITSTPC()
1066 // ******* TPC ********
1067 // 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
1068 // 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
1070 AliCodeTimerAuto("",0);
1072 // Loop over the tracks and extract and mask out all aod tracks that pass the selections for AODt racks
1073 for(int it = 0;it < fNumberOfTracks;++it)
1075 AliAODTrack *tr = (AliAODTrack*)(Tracks().At(it));
1077 UInt_t map = tr->GetFilterMap();
1078 if(map&fTPCConstrainedFilterMask){
1079 // we only reset the track select ionfo, no deletion...
1080 tr->SetFilterMap(map&~fTPCConstrainedFilterMask);
1082 if(map&fHybridFilterMaskTPCCG){
1083 // this is one part of the hybrid tracks
1084 // the others not passing the selection will be TPC only selected below
1085 tr->SetIsHybridTPCConstrainedGlobal(kTRUE);
1089 // Loop over the ESD trcks and pick out the tracks passing TPC only cuts
1090 const AliESDVertex *vtxSPD = esd.GetPrimaryVertexSPD();
1091 const AliESDVertex *vtx = esd.GetPrimaryVertex();
1093 Double_t pos[3] = { 0. };
1094 Double_t covTr[21]={0.};
1095 // Double_t pid[10]={0.};
1096 Double_t p[3] = { 0. };
1097 Double_t pDCA[3] = { 0. }; // momentum at DCA
1098 Double_t rDCA[3] = { 0. }; // position at DCA
1099 Float_t dDCA[2] = {0.}; // DCA to the vertex d and z
1100 Float_t cDCA[3] = {0.}; // covariance of impact parameters
1101 Int_t tofLabel[3] = {0};
1103 AliAODTrack* aodTrack(0x0);
1105 // account for change in pT after the constraint
1106 Float_t ptMax = 1E10;
1108 for(int i = 0;i<32;i++){
1109 if(fTPCConstrainedFilterMask&(1<<i)){
1110 AliESDtrackCuts*cuts = (AliESDtrackCuts*)fTrackFilter->GetCuts()->At(i);
1111 Float_t tmp1= 0,tmp2 = 0;
1112 cuts->GetPtRange(tmp1,tmp2);
1113 if(tmp1>ptMin)ptMin=tmp1;
1114 if(tmp2<ptMax)ptMax=tmp2;
1118 for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack)
1120 AliESDtrack* esdTrack = esd.GetTrack(nTrack); //carefull do not modify it othwise need to work with a copy
1122 UInt_t selectInfo = 0;
1123 Bool_t isHybridITSTPC = false;
1127 selectInfo = fTrackFilter->IsSelected(esdTrack);
1130 if(!(selectInfo&fHybridFilterMaskTPCCG)){
1131 // not already selected tracks, use second part of hybrid tracks
1132 isHybridITSTPC = true;
1133 // too save space one could only store these...
1136 selectInfo &= fTPCConstrainedFilterMask;
1137 if (!selectInfo) continue;
1138 if (fWriteHybridTPCCOnly&&!isHybridITSTPC) continue; // write only complementary tracks
1139 // create a tpc only tracl
1140 AliESDtrack *track = AliESDtrackCuts::GetTPCOnlyTrack(const_cast<AliESDEvent*>(&esd),esdTrack->GetID());
1141 if (!track) continue;
1143 if (track->Pt()>0.) {
1144 // only constrain tracks above threshold
1145 AliExternalTrackParam exParam;
1146 // take the B-field from the ESD, no 3D fieldMap available at this point
1147 Bool_t relate = false;
1148 relate = track->RelateToVertexTPC(vtxSPD,esd.GetMagneticField(),kVeryBig,&exParam);
1153 // fetch the track parameters at the DCA (unconstraint)
1154 if(track->GetTPCInnerParam()){
1155 track->GetTPCInnerParam()->GetPxPyPz(pDCA);
1156 track->GetTPCInnerParam()->GetXYZ(rDCA);
1158 // get the DCA to the vertex:
1159 track->GetImpactParametersTPC(dDCA,cDCA);
1160 // set the constrained parameters to the track
1161 track->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance());
1164 track->GetPxPyPz(p);
1166 Float_t pT = track->Pt();
1167 if(pT<ptMin||pT>ptMax){
1173 track->GetCovarianceXYZPxPyPz(covTr);
1174 // esdTrack->GetESDpid(pid);// original PID
1175 esdTrack->GetTOFLabel(tofLabel);
1176 if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());
1177 fUsedTrackCopy[nTrack] |= selectInfo;
1178 aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack((track->GetID()+1)*-1,
1185 (Short_t)track->GetSign(),
1186 track->GetITSClusterMap(),
1189 kTRUE, // check if this is right
1190 vtx->UsesTrack(track->GetID()),
1191 AliAODTrack::kPrimary,
1193 aodTrack->SetPIDForTracking(track->GetPIDForTracking());
1194 aodTrack->SetIsHybridTPCConstrainedGlobal(isHybridITSTPC);
1195 aodTrack->SetTPCFitMap(track->GetTPCFitMap());
1196 aodTrack->SetTPCClusterMap(track->GetTPCClusterMap());
1197 aodTrack->SetTPCSharedMap (track->GetTPCSharedMap());
1198 aodTrack->SetIsTPCConstrained(kTRUE);
1199 aodTrack->SetChi2perNDF(Chi2perNDF(esdTrack)); // original track
1200 // set the DCA values to the AOD track
1201 aodTrack->SetPxPyPzAtDCA(pDCA[0],pDCA[1],pDCA[2]);
1202 aodTrack->SetXYAtDCA(rDCA[0],rDCA[1]);
1203 aodTrack->SetDCA(dDCA[0],dDCA[1]);
1204 aodTrack->SetFlags(track->GetStatus());
1205 aodTrack->SetTPCPointsF(track->GetTPCNclsF());
1206 aodTrack->SetTPCNCrossedRows(UShort_t(track->GetTPCCrossedRows()));
1207 aodTrack->SetIntegratedLength(track->GetIntegratedLength());
1208 aodTrack->SetTOFLabel(tofLabel);
1209 CopyCaloProps(track,aodTrack);
1210 // do not duplicate PID information
1211 // aodTrack->ConvertAliPIDtoAODPID();
1212 // SetAODPID(esdTrack,aodTrack,detpid);
1214 } // end of loop on tracks
1217 //______________________________________________________________________________
1218 void AliAnalysisTaskESDfilter::ConvertGlobalConstrainedTracks(const AliESDEvent& esd)
1220 // Here we have the option to store the complement from global constraint information
1221 // to tracks passing tight cuts (1) in order not to get fakes back in, one needs
1222 // two sets of cuts one tight (1) (to throw out fakes) and one lose (2) (fakes/bad tracks would pass (2) but not (1))
1223 // using cut number (3) selects the tracks that complement (1) e.g. tracks witout ITS refit or cluster requirement
1225 AliCodeTimerAuto("",0);
1227 // Loop over the tracks and extract and mask out all aod tracks that pass the selections for AODt racks
1228 for(int it = 0;it < fNumberOfTracks;++it)
1230 AliAODTrack *tr = (AliAODTrack*)(Tracks().At(it));
1232 UInt_t map = tr->GetFilterMap();
1233 if(map&fGlobalConstrainedFilterMask){
1234 // we only reset the track select info, no deletion...
1235 // mask reset mask in case track is already taken
1236 tr->SetFilterMap(map&~fGlobalConstrainedFilterMask);
1238 if(map&fHybridFilterMaskGCG){
1239 // this is one part of the hybrid tracks
1240 // the others not passing the selection will be the ones selected below
1241 tr->SetIsHybridGlobalConstrainedGlobal(kTRUE);
1244 // Loop over the ESD trcks and pick out the tracks passing the GlobalConstraint cuts
1246 Double_t pos[3] = { 0. };
1247 Double_t covTr[21]={0.};
1248 // Double_t pid[10]={0.};
1249 Double_t p[3] = { 0. };
1251 Double_t pDCA[3] = { 0. }; // momentum at DCA
1252 Double_t rDCA[3] = { 0. }; // position at DCA
1253 Float_t dDCA[2] = {0.}; // DCA to the vertex d and z
1254 Float_t cDCA[3] = {0.}; // covariance of impact parameters
1255 Int_t tofLabel[3] = {0};
1257 AliAODTrack* aodTrack(0x0);
1258 AliAODPid* detpid(0x0);
1259 const AliESDVertex *vtx = esd.GetPrimaryVertex();
1261 // account for change in pT after the constraint
1262 Float_t ptMax = 1E10;
1264 for(int i = 0;i<32;i++){
1265 if(fGlobalConstrainedFilterMask&(1<<i)){
1266 AliESDtrackCuts*cuts = (AliESDtrackCuts*)fTrackFilter->GetCuts()->At(i);
1267 Float_t tmp1= 0,tmp2 = 0;
1268 cuts->GetPtRange(tmp1,tmp2);
1269 if(tmp1>ptMin)ptMin=tmp1;
1270 if(tmp2<ptMax)ptMax=tmp2;
1274 for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack)
1276 AliESDtrack* esdTrack = esd.GetTrack(nTrack); //carefull do not modify it othwise need to work with a copy
1277 const AliExternalTrackParam * exParamGC = esdTrack->GetConstrainedParam();
1278 if(!exParamGC)continue;
1280 UInt_t selectInfo = 0;
1281 Bool_t isHybridGC = false;
1286 selectInfo = fTrackFilter->IsSelected(esdTrack);
1289 if (!(selectInfo&fHybridFilterMaskGCG)) isHybridGC = true;
1290 if (fWriteHybridGCOnly&&!isHybridGC) continue; // write only complementary tracks
1292 selectInfo &= fGlobalConstrainedFilterMask;
1293 if (!selectInfo) continue;
1294 // fetch the track parameters at the DCA (unconstrained)
1295 esdTrack->GetPxPyPz(pDCA);
1296 esdTrack->GetXYZ(rDCA);
1297 // get the DCA to the vertex:
1298 esdTrack->GetImpactParameters(dDCA,cDCA);
1299 if (!esdTrack->GetConstrainedPxPyPz(p)) continue;
1301 Float_t pT = exParamGC->Pt();
1302 if(pT<ptMin||pT>ptMax){
1306 esdTrack->GetConstrainedXYZ(pos);
1307 exParamGC->GetCovarianceXYZPxPyPz(covTr);
1308 // esdTrack->GetESDpid(pid);
1309 esdTrack->GetTOFLabel(tofLabel);
1310 if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());
1311 fUsedTrackCopy[nTrack] |= selectInfo;
1312 aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack((esdTrack->GetID()+1)*-1,
1313 esdTrack->GetLabel(),
1319 (Short_t)esdTrack->GetSign(),
1320 esdTrack->GetITSClusterMap(),
1323 kTRUE, // check if this is right
1324 vtx->UsesTrack(esdTrack->GetID()),
1325 AliAODTrack::kPrimary,
1327 aodTrack->SetPIDForTracking(esdTrack->GetPIDForTracking());
1328 aodTrack->SetIsHybridGlobalConstrainedGlobal(isHybridGC);
1329 aodTrack->SetIsGlobalConstrained(kTRUE);
1330 aodTrack->SetTPCFitMap(esdTrack->GetTPCFitMap());
1331 aodTrack->SetTPCClusterMap(esdTrack->GetTPCClusterMap());
1332 aodTrack->SetTPCSharedMap (esdTrack->GetTPCSharedMap());
1333 aodTrack->SetChi2perNDF(Chi2perNDF(esdTrack));
1335 // set the DCA values to the AOD track
1336 aodTrack->SetPxPyPzAtDCA(pDCA[0],pDCA[1],pDCA[2]);
1337 aodTrack->SetXYAtDCA(rDCA[0],rDCA[1]);
1338 aodTrack->SetDCA(dDCA[0],dDCA[1]);
1339 aodTrack->SetFlags(esdTrack->GetStatus());
1340 aodTrack->SetTPCPointsF(esdTrack->GetTPCNclsF());
1341 aodTrack->SetTPCNCrossedRows(UShort_t(esdTrack->GetTPCCrossedRows()));
1342 aodTrack->SetIntegratedLength(esdTrack->GetIntegratedLength());
1343 aodTrack->SetTOFLabel(tofLabel);
1344 CopyCaloProps(esdTrack,aodTrack);
1346 // only copy AOD information for hybrid, no duplicate information
1347 aodTrack->ConvertAliPIDtoAODPID();
1348 SetAODPID(esdTrack,aodTrack,detpid);
1350 } // end of loop on tracks
1353 //______________________________________________________________________________
1354 void AliAnalysisTaskESDfilter::ConvertTracks(const AliESDEvent& esd)
1356 // Tracks (primary and orphan)
1358 AliCodeTimerAuto("",0);
1360 AliDebug(1,Form("NUMBER OF ESD TRACKS %5d\n", esd.GetNumberOfTracks()));
1362 const AliESDVertex *vtx = esd.GetPrimaryVertex();
1363 Double_t p[3] = { 0. };
1364 Double_t pos[3] = { 0. };
1365 Double_t covTr[21] = { 0. };
1366 // Double_t pid[10] = { 0. };
1367 Int_t tofLabel[3] = {0};
1368 AliAODTrack* aodTrack(0x0);
1369 AliAODPid* detpid(0x0);
1371 for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack)
1373 if (fUsedTrack[nTrack]) continue;
1375 AliESDtrack *esdTrack = esd.GetTrack(nTrack);
1376 UInt_t selectInfo = 0;
1380 selectInfo = fTrackFilter->IsSelected(esdTrack);
1381 if (!selectInfo && !vtx->UsesTrack(esdTrack->GetID())) continue;
1384 esdTrack->GetPxPyPz(p);
1385 esdTrack->GetXYZ(pos);
1386 esdTrack->GetCovarianceXYZPxPyPz(covTr);
1387 // esdTrack->GetESDpid(pid);
1388 esdTrack->GetTOFLabel(tofLabel);
1389 if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());
1390 fUsedTrack[nTrack] = kTRUE;
1391 aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrack->GetID(),
1392 esdTrack->GetLabel(),
1398 (Short_t)esdTrack->GetSign(),
1399 esdTrack->GetITSClusterMap(),
1402 kTRUE, // check if this is right
1403 vtx->UsesTrack(esdTrack->GetID()),
1404 AliAODTrack::kPrimary,
1406 fPrimaryVertex->AddDaughter(aodTrack);
1407 aodTrack->SetPIDForTracking(esdTrack->GetPIDForTracking());
1408 aodTrack->SetTPCFitMap(esdTrack->GetTPCFitMap());
1409 aodTrack->SetTPCClusterMap(esdTrack->GetTPCClusterMap());
1410 aodTrack->SetTPCSharedMap (esdTrack->GetTPCSharedMap());
1411 aodTrack->SetChi2perNDF(Chi2perNDF(esdTrack));
1412 aodTrack->SetTPCPointsF(esdTrack->GetTPCNclsF());
1413 aodTrack->SetTPCNCrossedRows(UShort_t(esdTrack->GetTPCCrossedRows()));
1414 aodTrack->SetIntegratedLength(esdTrack->GetIntegratedLength());
1415 aodTrack->SetTOFLabel(tofLabel);
1416 CopyCaloProps(esdTrack,aodTrack);
1417 fAODTrackRefs->AddAt(aodTrack, nTrack);
1418 if (esdTrack->GetSign() > 0) ++fNumberOfPositiveTracks;
1419 aodTrack->SetFlags(esdTrack->GetStatus());
1420 aodTrack->ConvertAliPIDtoAODPID();
1421 SetAODPID(esdTrack,aodTrack,detpid);
1422 } // end of loop on tracks
1425 //______________________________________________________________________________
1426 void AliAnalysisTaskESDfilter::ConvertPmdClusters(const AliESDEvent& esd)
1428 // Convert PMD Clusters
1429 AliCodeTimerAuto("",0);
1430 Int_t jPmdClusters=0;
1431 // Access to the AOD container of PMD clusters
1432 TClonesArray &pmdClusters = *(AODEvent()->GetPmdClusters());
1433 for (Int_t iPmd = 0; iPmd < esd.GetNumberOfPmdTracks(); ++iPmd) {
1434 // file pmd clusters, to be revised!
1435 AliESDPmdTrack *pmdTrack = esd.GetPmdTrack(iPmd);
1438 Double_t posPmd[3] = { pmdTrack->GetClusterX(), pmdTrack->GetClusterY(), pmdTrack->GetClusterZ()};
1439 Double_t pidPmd[13] = { 0.}; // to be revised!
1441 // assoc cluster not set
1442 new(pmdClusters[jPmdClusters++]) AliAODPmdCluster(iPmd, nLabel, label, pmdTrack->GetClusterADC(), posPmd, pidPmd);
1446 //______________________________________________________________________________
1447 void AliAnalysisTaskESDfilter::ConvertCaloClusters(const AliESDEvent& esd)
1449 // Convert Calorimeter Clusters
1450 AliCodeTimerAuto("",0);
1452 // Access to the AOD container of clusters
1453 TClonesArray &caloClusters = *(AODEvent()->GetCaloClusters());
1456 for (Int_t iClust=0; iClust<esd.GetNumberOfCaloClusters(); ++iClust) {
1457 AliESDCaloCluster * cluster = esd.GetCaloCluster(iClust);
1459 Int_t id = cluster->GetID();
1460 Int_t nLabel = cluster->GetNLabels();
1461 Int_t *labels = cluster->GetLabels();
1463 for(int i = 0;i < nLabel;++i) {
1464 if(fMChandler)fMChandler->SelectParticle(labels[i]);
1468 Float_t energy = cluster->E();
1469 Float_t posF[3] = { 0.};
1470 cluster->GetPosition(posF);
1472 AliAODCaloCluster *caloCluster = new(caloClusters[jClusters++]) AliAODCaloCluster(id,
1478 cluster->GetType(),0);
1480 caloCluster->SetCaloCluster(cluster->GetDistanceToBadChannel(),
1481 cluster->GetDispersion(),
1482 cluster->GetM20(), cluster->GetM02(),
1483 cluster->GetEmcCpvDistance(),
1484 cluster->GetNExMax(),cluster->GetTOF()) ;
1485 caloCluster->SetPIDFromESD(cluster->GetPID());
1486 caloCluster->SetNCells(cluster->GetNCells());
1487 caloCluster->SetCellsAbsId(cluster->GetCellsAbsId());
1488 caloCluster->SetCellsAmplitudeFraction(cluster->GetCellsAmplitudeFraction());
1489 caloCluster->SetTrackDistance(cluster->GetTrackDx(), cluster->GetTrackDz());
1491 Int_t nMatchCount = 0;
1492 TArrayI* matchedT = cluster->GetTracksMatched();
1493 if (fNumberOfTracks>0 && matchedT && cluster->GetTrackMatchedIndex() >= 0) {
1494 for (Int_t im = 0; im < matchedT->GetSize(); im++) {
1495 Int_t iESDtrack = matchedT->At(im);;
1496 if (fAODTrackRefs->At(iESDtrack) != 0) {
1497 caloCluster->AddTrackMatched((AliAODTrack*)fAODTrackRefs->At(iESDtrack));
1503 caloCluster->SetTrackDistance(-999,-999);
1506 caloClusters.Expand(jClusters); // resize TObjArray to 'remove' slots for pseudo clusters
1509 //______________________________________________________________________________
1510 void AliAnalysisTaskESDfilter::ConvertCaloTrigger(TString calo, const AliESDEvent& esd)
1512 AliCodeTimerAuto("",0);
1514 if (calo == "PHOS") {
1515 AliAODCaloTrigger &aodTrigger = *(AODEvent()->GetCaloTrigger(calo));
1516 AliESDCaloTrigger &esdTrigger = *(esd.GetCaloTrigger(calo));
1518 aodTrigger.Allocate(esdTrigger.GetEntries());
1523 while (esdTrigger.Next()) {
1524 esdTrigger.GetPosition(tmod,tabsId);
1525 esdTrigger.GetAmplitude(a);
1526 aodTrigger.Add(tmod,tabsId,a,0.,(Int_t*)NULL,0,0,0);
1531 AliAODHandler *aodHandler = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
1533 TTree *aodTree = aodHandler->GetTree();
1535 Int_t *type = esd.GetCaloTriggerType();
1536 for (Int_t i = 0; i < 15; i++) {
1537 aodTree->GetUserInfo()->Add(new TParameter<int>(Form("EMCALCaloTrigger%d",i), type[i]));
1542 AliAODCaloTrigger &aodTrigger = *(AODEvent()->GetCaloTrigger(calo));
1543 AliESDCaloTrigger &esdTrigger = *(esd.GetCaloTrigger(calo));
1544 aodTrigger.Allocate(esdTrigger.GetEntries());
1547 while (esdTrigger.Next()) {
1548 Int_t px, py, ts, nTimes, times[10], b;
1550 esdTrigger.GetPosition(px, py);
1551 esdTrigger.GetAmplitude(a);
1552 esdTrigger.GetTime(t);
1553 esdTrigger.GetL0Times(times);
1554 esdTrigger.GetNL0Times(nTimes);
1555 esdTrigger.GetL1TimeSum(ts);
1556 esdTrigger.GetTriggerBits(b);
1557 aodTrigger.Add(px, py, a, t, times, nTimes, ts, b);
1560 for (int i = 0; i < 4; i++) aodTrigger.SetL1Threshold(i, esdTrigger.GetL1Threshold(i));
1562 esdTrigger.GetL1V0(0),
1563 esdTrigger.GetL1V0(1)
1565 aodTrigger.SetL1V0(v0);
1566 aodTrigger.SetL1FrameMask(esdTrigger.GetL1FrameMask());
1569 //______________________________________________________________________________
1570 void AliAnalysisTaskESDfilter::ConvertEMCALCells(const AliESDEvent& esd)
1572 // Convert EMCAL Cells
1573 AliCodeTimerAuto("",0);
1575 // fill EMCAL cell info
1576 if (esd.GetEMCALCells()) { // protection against missing ESD information
1577 AliESDCaloCells &esdEMcells = *(esd.GetEMCALCells());
1578 Int_t nEMcell = esdEMcells.GetNumberOfCells() ;
1579 AliAODCaloCells &aodEMcells = *(AODEvent()->GetEMCALCells());
1580 aodEMcells.CreateContainer(nEMcell);
1581 aodEMcells.SetType(AliAODCaloCells::kEMCALCell);
1582 for (Int_t iCell = 0; iCell < nEMcell; iCell++) {
1583 aodEMcells.SetCell(iCell,esdEMcells.GetCellNumber(iCell),esdEMcells.GetAmplitude(iCell),
1584 esdEMcells.GetTime(iCell), esdEMcells.GetMCLabel(iCell), esdEMcells.GetEFraction(iCell),
1585 esdEMcells.GetHighGain(iCell) );
1591 //______________________________________________________________________________
1592 void AliAnalysisTaskESDfilter::ConvertPHOSCells(const AliESDEvent& esd)
1594 // Convert PHOS Cells
1595 AliCodeTimerAuto("",0);
1597 // fill PHOS cell info
1598 if (esd.GetPHOSCells()) { // protection against missing ESD information
1599 AliESDCaloCells &esdPHcells = *(esd.GetPHOSCells());
1600 Int_t nPHcell = esdPHcells.GetNumberOfCells() ;
1602 AliAODCaloCells &aodPHcells = *(AODEvent()->GetPHOSCells());
1603 aodPHcells.CreateContainer(nPHcell);
1604 aodPHcells.SetType(AliAODCaloCells::kPHOSCell);
1605 for (Int_t iCell = 0; iCell < nPHcell; iCell++) {
1606 aodPHcells.SetCell(iCell,esdPHcells.GetCellNumber(iCell),esdPHcells.GetAmplitude(iCell),
1607 esdPHcells.GetTime(iCell), esdPHcells.GetMCLabel(iCell), esdPHcells.GetEFraction(iCell),
1608 esdPHcells.GetHighGain(iCell) );
1614 //______________________________________________________________________________
1615 void AliAnalysisTaskESDfilter::ConvertTracklets(const AliESDEvent& esd)
1618 AliCodeTimerAuto("",0);
1620 AliAODTracklets &SPDTracklets = *(AODEvent()->GetTracklets());
1621 const AliMultiplicity *mult = esd.GetMultiplicity();
1623 if (mult->GetNumberOfTracklets()>0) {
1624 SPDTracklets.CreateContainer(mult->GetNumberOfTracklets());
1625 SPDTracklets.SetScaleDThetaBySin2T(mult->GetScaleDThetaBySin2T());
1626 for (Int_t n=0; n<mult->GetNumberOfTracklets(); n++) {
1628 fMChandler->SelectParticle(mult->GetLabel(n, 0));
1629 fMChandler->SelectParticle(mult->GetLabel(n, 1));
1631 SPDTracklets.SetTracklet(n, mult->GetTheta(n), mult->GetPhi(n), mult->GetDeltaPhi(n), mult->GetLabel(n, 0),mult->GetLabel(n, 1));
1635 //Printf("ERROR: AliMultiplicity could not be retrieved from ESD");
1639 //______________________________________________________________________________
1640 void AliAnalysisTaskESDfilter::ConvertKinks(const AliESDEvent& esd)
1642 AliCodeTimerAuto("",0);
1644 // Kinks: it is a big mess the access to the information in the kinks
1645 // The loop is on the tracks in order to find the mother and daugther of each kink
1647 Double_t covTr[21]={0.};
1648 // Double_t pid[10]={0.};
1649 AliAODPid* detpid(0x0);
1650 Int_t tofLabel[3] = {0};
1652 fNumberOfKinks = esd.GetNumberOfKinks();
1654 const AliESDVertex* vtx = esd.GetPrimaryVertex();
1656 for (Int_t iTrack=0; iTrack<esd.GetNumberOfTracks(); ++iTrack)
1658 AliESDtrack * esdTrack = esd.GetTrack(iTrack);
1660 Int_t ikink = esdTrack->GetKinkIndex(0);
1662 if (ikink && fNumberOfKinks) {
1663 // Negative kink index: mother, positive: daughter
1664 // Search for the second track of the kink
1666 for (Int_t jTrack = iTrack+1; jTrack<esd.GetNumberOfTracks(); ++jTrack) {
1667 AliESDtrack * esdTrack1 = esd.GetTrack(jTrack);
1668 Int_t jkink = esdTrack1->GetKinkIndex(0);
1670 if ( TMath::Abs(ikink)==TMath::Abs(jkink) ) {
1671 // The two tracks are from the same kink
1672 if (fUsedKink[TMath::Abs(ikink)-1]) continue; // skip used kinks
1675 Int_t idaughter = -1;
1677 if (ikink<0 && jkink>0) {
1680 } else if (ikink>0 && jkink<0) {
1684 //cerr << "Error: Wrong combination of kink indexes: "
1685 // << ikink << " " << jkink << endl;
1689 // Add the mother track if it passed primary track selection cuts
1690 AliAODTrack * mother = NULL;
1692 UInt_t selectInfo = 0;
1694 selectInfo = fTrackFilter->IsSelected(esd.GetTrack(imother));
1695 if (!selectInfo) continue;
1698 if (!fUsedTrack[imother]) {
1699 fUsedTrack[imother] = kTRUE;
1700 AliESDtrack *esdTrackM = esd.GetTrack(imother);
1701 Double_t p[3] = { 0. };
1702 Double_t pos[3] = { 0. };
1703 esdTrackM->GetPxPyPz(p);
1704 esdTrackM->GetXYZ(pos);
1705 esdTrackM->GetCovarianceXYZPxPyPz(covTr);
1706 // esdTrackM->GetESDpid(pid);
1707 esdTrackM->GetTOFLabel(tofLabel);
1708 if(fMChandler)fMChandler->SelectParticle(esdTrackM->GetLabel());
1709 mother = new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrackM->GetID(),
1710 esdTrackM->GetLabel(),
1716 (Short_t)esdTrackM->GetSign(),
1717 esdTrackM->GetITSClusterMap(),
1720 kTRUE, // check if this is right
1721 vtx->UsesTrack(esdTrack->GetID()),
1722 AliAODTrack::kPrimary,
1724 mother->SetPIDForTracking(esdTrackM->GetPIDForTracking());
1725 mother->SetTPCFitMap(esdTrackM->GetTPCFitMap());
1726 mother->SetTPCClusterMap(esdTrackM->GetTPCClusterMap());
1727 mother->SetTPCSharedMap (esdTrackM->GetTPCSharedMap());
1728 mother->SetChi2perNDF(Chi2perNDF(esdTrackM));
1729 mother->SetTPCPointsF(esdTrackM->GetTPCNclsF());
1730 mother->SetTPCNCrossedRows(UShort_t(esdTrackM->GetTPCCrossedRows()));
1731 mother->SetIntegratedLength(esdTrackM->GetIntegratedLength());
1732 mother->SetTOFLabel(tofLabel);
1733 CopyCaloProps(esdTrackM,mother);
1734 fAODTrackRefs->AddAt(mother, imother);
1735 if (esdTrackM->GetSign() > 0) ++fNumberOfPositiveTracks;
1736 mother->SetFlags(esdTrackM->GetStatus());
1737 mother->ConvertAliPIDtoAODPID();
1738 fPrimaryVertex->AddDaughter(mother);
1739 mother->ConvertAliPIDtoAODPID();
1740 SetAODPID(esdTrackM,mother,detpid);
1743 //cerr << "Error: event " << esd.GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1
1744 // << " track " << imother << " has already been used!" << endl;
1747 // Add the kink vertex
1748 AliESDkink * kink = esd.GetKink(TMath::Abs(ikink)-1);
1750 AliAODVertex * vkink = new(Vertices()[fNumberOfVertices++]) AliAODVertex(kink->GetPosition(),
1754 esdTrack->GetID(), // ID of mother's track!
1755 AliAODVertex::kKink);
1756 // Add the daughter track
1757 AliAODTrack * daughter = NULL;
1758 if (!fUsedTrack[idaughter]) {
1759 fUsedTrack[idaughter] = kTRUE;
1760 AliESDtrack *esdTrackD = esd.GetTrack(idaughter);
1761 Double_t p[3] = { 0. };
1762 Double_t pos[3] = { 0. };
1763 esdTrackD->GetPxPyPz(p);
1764 esdTrackD->GetXYZ(pos);
1765 esdTrackD->GetCovarianceXYZPxPyPz(covTr);
1766 // esdTrackD->GetESDpid(pid);
1767 esdTrackD->GetTOFLabel(tofLabel);
1769 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdTrackD);
1770 if(fMChandler)fMChandler->SelectParticle(esdTrackD->GetLabel());
1771 daughter = new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrackD->GetID(),
1772 esdTrackD->GetLabel(),
1778 (Short_t)esdTrackD->GetSign(),
1779 esdTrackD->GetITSClusterMap(),
1782 kTRUE, // check if this is right
1783 vtx->UsesTrack(esdTrack->GetID()),
1784 AliAODTrack::kFromDecayVtx,
1786 daughter->SetPIDForTracking(esdTrackD->GetPIDForTracking());
1787 daughter->SetTPCFitMap(esdTrackD->GetTPCFitMap());
1788 daughter->SetTPCClusterMap(esdTrackD->GetTPCClusterMap());
1789 daughter->SetTPCSharedMap (esdTrackD->GetTPCSharedMap());
1790 daughter->SetTPCPointsF(esdTrackD->GetTPCNclsF());
1791 daughter->SetTPCNCrossedRows(UShort_t(esdTrackD->GetTPCCrossedRows()));
1792 daughter->SetIntegratedLength(esdTrackD->GetIntegratedLength());
1793 daughter->SetTOFLabel(tofLabel);
1794 CopyCaloProps(esdTrackD,daughter);
1795 fAODTrackRefs->AddAt(daughter, idaughter);
1796 if (esdTrackD->GetSign() > 0) ++fNumberOfPositiveTracks;
1797 daughter->SetFlags(esdTrackD->GetStatus());
1798 daughter->ConvertAliPIDtoAODPID();
1799 vkink->AddDaughter(daughter);
1800 daughter->ConvertAliPIDtoAODPID();
1801 SetAODPID(esdTrackD,daughter,detpid);
1803 //cerr << "Error: event " << esd.GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1
1804 // << " track " << idaughter << " has already been used!" << endl;
1812 //______________________________________________________________________________
1813 void AliAnalysisTaskESDfilter::ConvertPrimaryVertices(const AliESDEvent& esd)
1815 AliCodeTimerAuto("",0);
1817 // Access to the AOD container of vertices
1818 fNumberOfVertices = 0;
1820 Double_t pos[3] = { 0. };
1821 Double_t covVtx[6] = { 0. };
1823 // Add primary vertex. The primary tracks will be defined
1824 // after the loops on the composite objects (V0, cascades, kinks)
1825 const AliESDVertex *vtx = esd.GetPrimaryVertex();
1827 vtx->GetXYZ(pos); // position
1828 vtx->GetCovMatrix(covVtx); //covariance matrix
1830 fPrimaryVertex = new(Vertices()[fNumberOfVertices++])
1831 AliAODVertex(pos, covVtx, vtx->GetChi2toNDF(), NULL, -1, AliAODVertex::kPrimary);
1832 fPrimaryVertex->SetName(vtx->GetName());
1833 fPrimaryVertex->SetTitle(vtx->GetTitle());
1834 fPrimaryVertex->SetBC(vtx->GetBC());
1836 TString vtitle = vtx->GetTitle();
1837 if (!vtitle.Contains("VertexerTracks"))
1838 fPrimaryVertex->SetNContributors(vtx->GetNContributors());
1840 if (fDebug > 0) fPrimaryVertex->Print();
1842 // Add SPD "main" vertex
1843 const AliESDVertex *vtxS = esd.GetPrimaryVertexSPD();
1844 vtxS->GetXYZ(pos); // position
1845 vtxS->GetCovMatrix(covVtx); //covariance matrix
1846 AliAODVertex * mVSPD = new(Vertices()[fNumberOfVertices++])
1847 AliAODVertex(pos, covVtx, vtxS->GetChi2toNDF(), NULL, -1, AliAODVertex::kMainSPD);
1848 mVSPD->SetName(vtxS->GetName());
1849 mVSPD->SetTitle(vtxS->GetTitle());
1850 mVSPD->SetNContributors(vtxS->GetNContributors());
1852 // Add SPD pileup vertices
1853 for(Int_t iV=0; iV<esd.GetNumberOfPileupVerticesSPD(); ++iV) {
1854 const AliESDVertex *vtxP = esd.GetPileupVertexSPD(iV);
1855 vtxP->GetXYZ(pos); // position
1856 vtxP->GetCovMatrix(covVtx); //covariance matrix
1857 AliAODVertex * pVSPD = new(Vertices()[fNumberOfVertices++])
1858 AliAODVertex(pos, covVtx, vtxP->GetChi2toNDF(), NULL, -1, AliAODVertex::kPileupSPD);
1859 pVSPD->SetName(vtxP->GetName());
1860 pVSPD->SetTitle(vtxP->GetTitle());
1861 pVSPD->SetNContributors(vtxP->GetNContributors());
1862 pVSPD->SetBC(vtxP->GetBC());
1865 // Add TRK pileup vertices
1866 for(Int_t iV=0; iV<esd.GetNumberOfPileupVerticesTracks(); ++iV) {
1867 const AliESDVertex *vtxP = esd.GetPileupVertexTracks(iV);
1868 vtxP->GetXYZ(pos); // position
1869 vtxP->GetCovMatrix(covVtx); //covariance matrix
1870 AliAODVertex * pVTRK = new(Vertices()[fNumberOfVertices++])
1871 AliAODVertex(pos, covVtx, vtxP->GetChi2toNDF(), NULL, -1, AliAODVertex::kPileupTracks);
1872 pVTRK->SetName(vtxP->GetName());
1873 pVTRK->SetTitle(vtxP->GetTitle());
1874 pVTRK->SetNContributors(vtxP->GetNContributors());
1875 pVTRK->SetBC(vtxP->GetBC());
1878 // Add TPC "main" vertex
1879 const AliESDVertex *vtxT = esd.GetPrimaryVertexTPC();
1880 vtxT->GetXYZ(pos); // position
1881 vtxT->GetCovMatrix(covVtx); //covariance matrix
1882 AliAODVertex * mVTPC = new(Vertices()[fNumberOfVertices++])
1883 AliAODVertex(pos, covVtx, vtxT->GetChi2toNDF(), NULL, -1, AliAODVertex::kMainTPC);
1884 mVTPC->SetName(vtxT->GetName());
1885 mVTPC->SetTitle(vtxT->GetTitle());
1886 mVTPC->SetNContributors(vtxT->GetNContributors());
1889 //______________________________________________________________________________
1890 void AliAnalysisTaskESDfilter::ConvertVZERO(const AliESDEvent& esd)
1892 // Convert VZERO data
1893 AliAODVZERO* vzeroData = AODEvent()->GetVZEROData();
1894 *vzeroData = *(esd.GetVZEROData());
1897 //______________________________________________________________________________
1898 void AliAnalysisTaskESDfilter::ConvertTZERO(const AliESDEvent& esd)
1900 // Convert TZERO data
1901 const AliESDTZERO* esdTzero = esd.GetESDTZERO();
1902 AliAODTZERO* aodTzero = AODEvent()->GetTZEROData();
1904 for (Int_t icase=0; icase<3; icase++){
1905 aodTzero->SetT0TOF( icase, esdTzero->GetT0TOF(icase));
1906 aodTzero->SetT0TOFbest(icase, esdTzero->GetT0TOFbest(icase));
1908 aodTzero->SetBackgroundFlag(esdTzero->GetBackgroundFlag());
1909 aodTzero->SetPileupFlag(esdTzero->GetPileupFlag());
1910 aodTzero->SetSatelliteFlag(esdTzero->GetSatellite());
1912 Float_t rawTime[24];
1913 for(Int_t ipmt=0; ipmt<24; ipmt++)
1914 rawTime[ipmt] = esdTzero->GetTimeFull(ipmt,0);
1916 Int_t idxOfFirstPmtA = -1, idxOfFirstPmtC = -1;
1917 Float_t timeOfFirstPmtA = 9999, timeOfFirstPmtC = 9999;
1918 for(int ipmt=0; ipmt<12; ipmt++){
1919 if( rawTime[ipmt] > -200 && rawTime[ipmt] < timeOfFirstPmtC && rawTime[ipmt]!=0){
1920 timeOfFirstPmtC = rawTime[ipmt];
1921 idxOfFirstPmtC = ipmt;
1924 for(int ipmt=12; ipmt<24; ipmt++){
1925 if( rawTime[ipmt] > -200 && rawTime[ipmt] < timeOfFirstPmtA && rawTime[ipmt]!=0 ){
1926 timeOfFirstPmtA = rawTime[ipmt];
1927 idxOfFirstPmtA = ipmt;
1931 if(idxOfFirstPmtA != -1 && idxOfFirstPmtC != -1){
1932 //speed of light in cm/ns TMath::C()*1e-7
1933 Float_t vertexraw = TMath::C()*1e-7 * (rawTime[idxOfFirstPmtA] - rawTime[idxOfFirstPmtC])/2;
1934 aodTzero->SetT0VertexRaw( vertexraw );
1936 aodTzero->SetT0VertexRaw(99999);
1939 aodTzero->SetT0zVertex(esdTzero->GetT0zVertex());
1942 const Double32_t *amp=esdTzero->GetT0amplitude();
1943 for(int ipmt=0; ipmt<24; ipmt++)
1944 aodTzero->SetAmp(ipmt, amp[ipmt]);
1945 aodTzero->SetAmp(24,esdTzero->GetMultC() );
1946 aodTzero->SetAmp(25,esdTzero->GetMultA() );
1949 //______________________________________________________________________________
1950 void AliAnalysisTaskESDfilter::ConvertZDC(const AliESDEvent& esd)
1953 AliESDZDC* esdZDC = esd.GetZDCData();
1955 const Double_t zem1Energy = esdZDC->GetZEM1Energy();
1956 const Double_t zem2Energy = esdZDC->GetZEM2Energy();
1958 const Double_t *towZNC = esdZDC->GetZNCTowerEnergy();
1959 const Double_t *towZPC = esdZDC->GetZPCTowerEnergy();
1960 const Double_t *towZNA = esdZDC->GetZNATowerEnergy();
1961 const Double_t *towZPA = esdZDC->GetZPATowerEnergy();
1962 const Double_t *towZNCLG = esdZDC->GetZNCTowerEnergyLR();
1963 const Double_t *towZNALG = esdZDC->GetZNATowerEnergyLR();
1965 AliAODZDC* zdcAOD = AODEvent()->GetZDCData();
1967 zdcAOD->SetZEM1Energy(zem1Energy);
1968 zdcAOD->SetZEM2Energy(zem2Energy);
1969 zdcAOD->SetZNCTowers(towZNC, towZNCLG);
1970 zdcAOD->SetZNATowers(towZNA, towZNALG);
1971 zdcAOD->SetZPCTowers(towZPC);
1972 zdcAOD->SetZPATowers(towZPA);
1974 zdcAOD->SetZDCParticipants(esdZDC->GetZDCParticipants(), esdZDC->GetZDCPartSideA(), esdZDC->GetZDCPartSideC());
1975 zdcAOD->SetZDCImpactParameter(esdZDC->GetImpactParameter(), esdZDC->GetImpactParamSideA(), esdZDC->GetImpactParamSideC());
1976 zdcAOD->SetZDCTDCSum(esdZDC->GetZNTDCSum(0));
1977 zdcAOD->SetZDCTDCDiff(esdZDC->GetZNTDCDiff(0));
1978 if(esdZDC->IsZNChit()) zdcAOD->SetZNCTDC(esdZDC->GetZDCTDCCorrected(10,0));
1979 if(esdZDC->IsZNAhit()) zdcAOD->SetZNATDC(esdZDC->GetZDCTDCCorrected(12,0));
1982 //_____________________________________________________________________________
1983 Int_t AliAnalysisTaskESDfilter::ConvertHMPID(const AliESDEvent& esd) // clm
1986 // Convtert ESD HMPID info to AOD and return the number of good tracks with HMPID signal.
1987 // We need to return an int since there is no signal counter in the ESD.
1990 AliCodeTimerAuto("",0);
1992 Int_t cntHmpidGoodTracks = 0;
2001 Float_t thetaTrk = 0;
2004 Double_t hmpPid[5]={0};
2005 Double_t hmpMom[3]={0};
2007 TClonesArray &hmpidRings = *(AODEvent()->GetHMPIDrings());
2009 for (Int_t iTrack=0; iTrack<esd.GetNumberOfTracks(); ++iTrack) {
2010 if(! esd.GetTrack(iTrack) ) continue;
2012 if(esd.GetTrack(iTrack)->GetHMPIDsignal() > -20 ) { //
2013 (esd.GetTrack(iTrack))->GetHMPIDmip(xMip, yMip, qMip, nphMip); // Get MIP properties
2014 (esd.GetTrack(iTrack))->GetHMPIDtrk(xTrk,yTrk,thetaTrk,phiTrk);
2015 (esd.GetTrack(iTrack))->GetHMPIDpid(hmpPid);
2016 if((esd.GetTrack(iTrack))->GetOuterHmpParam()) (esd.GetTrack(iTrack))->GetOuterHmpPxPyPz(hmpMom);
2018 if(esd.GetTrack(iTrack)->GetHMPIDsignal() == 0 && thetaTrk == 0 && qMip == 0 && nphMip ==0 ) continue; //
2020 new(hmpidRings[cntHmpidGoodTracks++]) AliAODHMPIDrings((esd.GetTrack(iTrack))->GetID(), // Unique track id to attach the ring to
2021 1000000*nphMip+qMip, // MIP charge and number of photons
2022 (esd.GetTrack(iTrack))->GetHMPIDcluIdx(), // 1000000*chamber id + cluster idx of assigned MIP cluster
2023 thetaTrk, // track inclination angle theta
2024 phiTrk, // track inclination angle phi
2025 (esd.GetTrack(iTrack))->GetHMPIDsignal(), // Cherenkov angle
2026 (esd.GetTrack(iTrack))->GetHMPIDoccupancy(), // Occupancy claculated for the given chamber
2027 (esd.GetTrack(iTrack))->GetHMPIDchi2(), // Ring resolution squared
2028 xTrk, // Track x coordinate (LORS)
2029 yTrk, // Track y coordinate (LORS)
2030 xMip, // MIP x coordinate (LORS)
2031 yMip, // MIP y coordinate (LORS)
2032 hmpPid, // PID probablities from ESD, remove once it is CombinedPid
2033 hmpMom // Track momentum in HMPID at ring reconstruction
2038 return cntHmpidGoodTracks;
2041 void AliAnalysisTaskESDfilter::ConvertTRD(const AliESDEvent& esd)
2043 // fill TRD on-line tracks with assiocated tracklets
2044 // as used for the TRD level-1 triggers
2046 const Int_t nTrdTracks = esd.GetNumberOfTrdTracks();
2047 const Int_t nLayers = 6;
2049 for (Int_t iTrdTrack = 0; iTrdTrack < nTrdTracks; ++iTrdTrack) {
2050 // copy information from ESD track to AOD track
2051 const AliESDTrdTrack *esdTrdTrk = esd.GetTrdTrack(iTrdTrack);
2052 AliAODTrdTrack &aodTrdTrk = AODEvent()->AddTrdTrack(esdTrdTrk);
2054 // copy the contributing tracklets
2055 for (Int_t iTracklet = 0; iTracklet < nLayers; ++iTracklet) {
2056 if (const AliESDTrdTracklet *esdTrdTrkl = esdTrdTrk->GetTracklet(iTracklet))
2057 aodTrdTrk.AddTracklet(*esdTrdTrkl, iTracklet);
2060 // add the reference to the matched global track
2061 AliAODTrack *aodTrkMatch = 0x0;
2062 AliESDtrack *esdTrkMatch = (AliESDtrack*) esdTrdTrk->GetTrackMatch();
2064 Int_t idx = esdTrkMatch->GetID();
2067 AliError("track has a matched track that was not found");
2068 else if (esdTrkMatch != esd.GetTrack(idx))
2069 AliError("wrong track found for ESD track index");
2071 UInt_t selectInfo = fTrackFilter ? fTrackFilter->IsSelected(esdTrkMatch) : 0;
2073 if (fUsedTrack[idx]) {
2074 aodTrkMatch = (AliAODTrack*) (*fAODTrackRefs)[idx];
2075 AliDebug(2, Form("event %lld: existing track (idx %i, pt = %f) matched to TRD track (idx %i, pt = %f), cut flags: 0x%08x",
2076 Entry(), idx, esdTrkMatch->Pt(), iTrdTrack, esdTrdTrk->Pt(),
2079 if (selectInfo & fUsedTrackCopy[idx]) {
2080 // mask filter bits already used in track copies
2081 selectInfo &= ~fUsedTrackCopy[idx];
2082 AliWarning(Form("event %lld: copied 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));
2086 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",
2087 Entry(), idx, esdTrkMatch->Pt(), iTrdTrack, esdTrdTrk->Pt(),
2088 fTrackFilter->IsSelected(esdTrkMatch), selectInfo));
2090 Double_t mom[3]={0.};
2091 Double_t pos[3]={0.};
2092 Double_t covTr[21]={0.};
2093 // Double_t pid[10]={0.};
2095 esdTrkMatch->GetPxPyPz(mom);
2096 esdTrkMatch->GetXYZ(pos);
2097 esdTrkMatch->GetCovarianceXYZPxPyPz(covTr);
2098 // esdTrkMatch->GetESDpid(pid);
2099 const AliESDVertex* vtx = esd.GetPrimaryVertex();
2101 fUsedTrack[idx] = kTRUE;
2102 if(fMChandler) fMChandler->SelectParticle(esdTrkMatch->GetLabel());
2104 aodTrkMatch = new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrkMatch->GetID(),
2105 esdTrkMatch->GetLabel(),
2111 (Short_t)esdTrkMatch->GetSign(),
2112 esdTrkMatch->GetITSClusterMap(),
2116 vtx->UsesTrack(esdTrkMatch->GetID()),
2117 AliAODTrack::kUndef,
2119 aodTrkMatch->SetPIDForTracking(esdTrkMatch->GetPIDForTracking());
2120 aodTrkMatch->SetTPCFitMap(esdTrkMatch->GetTPCFitMap());
2121 aodTrkMatch->SetTPCClusterMap(esdTrkMatch->GetTPCClusterMap());
2122 aodTrkMatch->SetTPCSharedMap (esdTrkMatch->GetTPCSharedMap());
2123 aodTrkMatch->SetChi2perNDF(Chi2perNDF(esdTrkMatch));
2124 aodTrkMatch->SetTPCPointsF(esdTrkMatch->GetTPCNclsF());
2125 aodTrkMatch->SetTPCNCrossedRows(UShort_t(esdTrkMatch->GetTPCCrossedRows()));
2126 aodTrkMatch->SetIntegratedLength(esdTrkMatch->GetIntegratedLength());
2127 CopyCaloProps(esdTrkMatch,aodTrkMatch);
2128 fAODTrackRefs->AddAt(aodTrkMatch,idx);
2129 if (esdTrkMatch->GetSign() > 0) ++fNumberOfPositiveTracks;
2130 aodTrkMatch->ConvertAliPIDtoAODPID();
2131 aodTrkMatch->SetFlags(esdTrkMatch->GetStatus());
2135 aodTrdTrk.SetTrackMatchReference(aodTrkMatch);
2139 //______________________________________________________________________________
2140 void AliAnalysisTaskESDfilter::ConvertESDtoAOD()
2142 // ESD Filter analysis task executed for each event
2144 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
2148 AliCodeTimerAuto("",0);
2150 if (fRefitVertexTracks>=0) AliESDUtils::RefitESDVertexTracks(esd,fRefitVertexTracks,
2151 fRefitVertexTracksNCuts ? fRefitVertexTracksCuts:0);
2153 fOldESDformat = ( esd->GetAliESDOld() != 0x0 );
2155 // Reconstruct cascades and V0 here
2156 if (fIsV0CascadeRecoEnabled) {
2157 esd->ResetCascades();
2160 AliV0vertexer lV0vtxer;
2161 AliCascadeVertexer lCascVtxer;
2163 lV0vtxer.SetCuts(fV0Cuts);
2164 lCascVtxer.SetCuts(fCascadeCuts);
2167 lV0vtxer.Tracks2V0vertices(esd);
2168 lCascVtxer.V0sTracks2CascadeVertices(esd);
2171 // Perform progagation of tracks if needed
2172 if (fDoPropagateTrackToEMCal) {
2173 const Int_t ntrack = esd->GetNumberOfTracks();
2174 for (Int_t i=0;i<ntrack;++i) {
2175 AliESDtrack *t = esd->GetTrack(i);
2176 AliEMCALRecoUtils::ExtrapolateTrackToEMCalSurface(t,fEMCalSurfaceDistance);
2180 fNumberOfTracks = 0;
2181 fNumberOfPositiveTracks = 0;
2183 fNumberOfVertices = 0;
2184 fNumberOfCascades = 0;
2187 AliAODHeader* header = ConvertHeader(*esd);
2189 if ( fIsVZEROEnabled ) ConvertVZERO(*esd);
2190 if ( fIsTZEROEnabled ) ConvertTZERO(*esd);
2192 // Fetch Stack for debuggging if available
2195 fMChandler = (AliMCEventHandler*) ((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler());
2198 // loop over events and fill them
2199 // Multiplicity information needed by the header (to be revised!)
2200 Int_t nTracks = esd->GetNumberOfTracks();
2202 // The line below should not be needed anymore (tracks already connected)
2203 // for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) esd->GetTrack(iTrack)->SetESDEvent(esd);
2205 // Update the header
2206 Int_t nV0s = esd->GetNumberOfV0s();
2207 Int_t nCascades = esd->GetNumberOfCascades();
2208 Int_t nKinks = esd->GetNumberOfKinks();
2209 Int_t nVertices = nV0s + nCascades /*V0 wihtin cascade already counted*/+ nKinks + 1 /* = prim. vtx*/;
2210 Int_t nPileSPDVertices=1+esd->GetNumberOfPileupVerticesSPD(); // also SPD main vertex
2211 Int_t nPileTrkVertices=esd->GetNumberOfPileupVerticesTracks();
2212 nVertices+=nPileSPDVertices;
2213 nVertices+=nPileTrkVertices;
2215 Int_t nCaloClus = esd->GetNumberOfCaloClusters();
2217 Int_t nPmdClus = esd->GetNumberOfPmdTracks();
2218 Int_t nHmpidRings = 0;
2220 AliDebug(1,Form(" NV0=%d NCASCADES=%d NKINKS=%d", nV0s, nCascades, nKinks));
2222 AODEvent()->ResetStd(nTracks, nVertices, nV0s, nCascades, nJets, nCaloClus, nFmdClus, nPmdClus, nHmpidRings);
2225 // RefArray to store a mapping between esd V0 number and newly created AOD-Vertex V0
2226 fAODV0VtxRefs = new TRefArray(nV0s);
2227 // RefArray to store the mapping between esd V0 number and newly created AOD-V0
2228 fAODV0Refs = new TRefArray(nV0s);
2229 // Array to take into account the V0s already added to the AOD (V0 within cascades)
2230 fUsedV0 = new Bool_t[nV0s];
2231 for (Int_t iV0=0; iV0<nV0s; ++iV0) fUsedV0[iV0]=kFALSE;
2235 // RefArray to store the mapping between esd track number and newly created AOD-Track
2237 fAODTrackRefs = new TRefArray(nTracks);
2239 // Array to take into account the tracks already added to the AOD
2240 fUsedTrack = new Bool_t[nTracks];
2241 fUsedTrackCopy = new UInt_t[nTracks];
2242 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) {
2243 fUsedTrack[iTrack]=kFALSE;
2244 fUsedTrackCopy[iTrack] = 0;
2248 // Array to take into account the kinks already added to the AOD
2250 fUsedKink = new Bool_t[nKinks];
2251 for (Int_t iKink=0; iKink<nKinks; ++iKink) fUsedKink[iKink]=kFALSE;
2254 ConvertPrimaryVertices(*esd);
2256 //setting best TOF PID
2257 AliESDInputHandler* esdH = dynamic_cast<AliESDInputHandler*>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
2259 fESDpid = esdH->GetESDpid();
2261 if (fIsPidOwner && fESDpid) {
2265 if (!fESDpid) { //in case of no Tender attached
2266 fESDpid = new AliESDpid;
2267 fIsPidOwner = kTRUE;
2270 if (!esd->GetTOFHeader()) { //protection in case the pass2 LHC10b,c,d have been processed without tender.
2271 Float_t t0spread[10];
2272 Float_t intrinsicTOFres=100; //ps ok for LHC10b,c,d pass2!!
2273 for (Int_t i=0; i<10; i++) t0spread[i] = (TMath::Sqrt(esd->GetSigma2DiamondZ()))/0.03; //0.03 to convert from cm to ps
2274 fESDpid->GetTOFResponse().SetT0resolution(t0spread);
2275 fESDpid->GetTOFResponse().SetTimeResolution(intrinsicTOFres);
2276 // fESDpid->SetTOFResponse(esd, (AliESDpid::EStartTimeType_t)fTimeZeroType);
2277 AliTOFHeader tmpTOFHeader(0,t0spread[0],0,NULL,NULL,NULL,intrinsicTOFres,t0spread[0]);
2278 AODEvent()->SetTOFHeader(&tmpTOFHeader); // write dummy TOF header in AOD
2280 AODEvent()->SetTOFHeader(esd->GetTOFHeader()); // write TOF header in AOD
2283 // In case of AOD production strating form LHC10e without Tender.
2284 //if(esd->GetTOFHeader() && fIsPidOwner) fESDpid->SetTOFResponse(esd, (AliESDpid::EStartTimeType_t)fTimeZeroType);
2286 if (fAreCascadesEnabled) ConvertCascades(*esd);
2287 if (fAreV0sEnabled) ConvertV0s(*esd);
2288 if (fAreKinksEnabled) ConvertKinks(*esd);
2289 if (fAreTracksEnabled) ConvertTracks(*esd);
2291 // Update number of AOD tracks in header at the end of track loop (M.G.)
2292 header->SetRefMultiplicity(fNumberOfTracks);
2293 header->SetRefMultiplicityPos(fNumberOfPositiveTracks);
2294 header->SetRefMultiplicityNeg(fNumberOfTracks - fNumberOfPositiveTracks);
2296 if (fTPCConstrainedFilterMask) ConvertTPCOnlyTracks(*esd);
2297 if (fGlobalConstrainedFilterMask) ConvertGlobalConstrainedTracks(*esd);
2298 if (fArePmdClustersEnabled) ConvertPmdClusters(*esd);
2299 if (fAreCaloClustersEnabled) ConvertCaloClusters(*esd);
2300 if (fAreEMCALCellsEnabled) ConvertEMCALCells(*esd);
2301 if (fArePHOSCellsEnabled) ConvertPHOSCells(*esd);
2302 if (fAreEMCALTriggerEnabled) ConvertCaloTrigger(TString("EMCAL"), *esd);
2303 if (fArePHOSTriggerEnabled) ConvertCaloTrigger(TString("PHOS"), *esd);
2304 if (fAreTrackletsEnabled) ConvertTracklets(*esd);
2305 if (fIsZDCEnabled) ConvertZDC(*esd);
2306 if (fIsHMPIDEnabled) nHmpidRings = ConvertHMPID(*esd);
2307 if (fIsTRDEnabled) ConvertTRD(*esd);
2309 delete fAODTrackRefs; fAODTrackRefs=0x0;
2310 delete fAODV0VtxRefs; fAODV0VtxRefs=0x0;
2311 delete fAODV0Refs; fAODV0Refs=0x0;
2312 delete[] fUsedTrack; fUsedTrack=0x0;
2313 delete[] fUsedTrackCopy; fUsedTrackCopy=0x0;
2314 delete[] fUsedV0; fUsedV0=0x0;
2315 delete[] fUsedKink; fUsedKink=0x0;
2321 AODEvent()->ConnectTracks();
2324 //______________________________________________________________________________
2325 void AliAnalysisTaskESDfilter::SetAODPID(AliESDtrack *esdtrack, AliAODTrack *aodtrack, AliAODPid *detpid)
2328 // Setter for the raw PID detector signals
2331 // Save PID object for candidate electrons
2332 Bool_t pidSave = kFALSE;
2334 Bool_t selectInfo = fTrackFilter->IsSelected((char*) "Electrons");
2335 if (selectInfo) pidSave = kTRUE;
2338 // Tracks passing pt cut
2339 if(esdtrack->Pt()>fHighPthreshold) {
2343 if(esdtrack->Pt()> fPtshape->GetXmin()){
2344 Double_t y = fPtshape->Eval(esdtrack->Pt())/fPtshape->Eval(fHighPthreshold);
2345 if(gRandom->Rndm(0)<1./y){
2349 }//end if p function
2353 if(!aodtrack->GetDetPid()){// prevent memory leak when calling SetAODPID twice for the same track
2354 detpid = new AliAODPid();
2355 SetDetectorRawSignals(detpid,esdtrack);
2356 aodtrack->SetDetPID(detpid);
2361 //______________________________________________________________________________
2362 void AliAnalysisTaskESDfilter::SetDetectorRawSignals(AliAODPid *aodpid, AliESDtrack *track)
2364 // Assignment of the detector signals (AliXXXesdPID inspired)
2367 AliInfo("no ESD track found. .....exiting");
2372 aodpid->SetTPCmomentum(track->GetTPCmomentum());
2373 aodpid->SetTPCTgl(track->GetTPCTgl());
2374 aodpid->SetITSsignal(track->GetITSsignal());
2375 Double_t itsdedx[4]; // dE/dx samples for individual ITS layers
2376 track->GetITSdEdxSamples(itsdedx);
2377 aodpid->SetITSdEdxSamples(itsdedx);
2379 aodpid->SetTPCsignal(track->GetTPCsignal());
2380 aodpid->SetTPCsignalN(track->GetTPCsignalN());
2381 if (track->GetTPCdEdxInfo()) aodpid->SetTPCdEdxInfo(track->GetTPCdEdxInfo());
2384 Int_t nslices = track->GetNumberOfTRDslices()*6;
2385 TArrayD trdslices(nslices);
2386 for(Int_t iSl =0; iSl < track->GetNumberOfTRDslices(); iSl++) {
2387 for(Int_t iPl =0; iPl<6; iPl++) trdslices[iPl*track->GetNumberOfTRDslices()+iSl] = track->GetTRDslice(iPl,iSl);
2391 for(Int_t iPl=0;iPl<6;iPl++){
2392 Double_t trdmom=track->GetTRDmomentum(iPl);
2393 aodpid->SetTRDmomentum(iPl,trdmom);
2396 aodpid->SetTRDslices(track->GetNumberOfTRDslices()*6,trdslices.GetArray());
2397 aodpid->SetTRDsignal(track->GetTRDsignal());
2399 //TRD clusters and tracklets
2400 aodpid->SetTRDncls(track->GetTRDncls());
2401 aodpid->SetTRDntrackletsPID(track->GetTRDntrackletsPID());
2403 aodpid->SetTRDChi2(track->GetTRDchi2());
2406 Double_t times[AliPID::kSPECIESC]; track->GetIntegratedTimes(times,AliPID::kSPECIESC);
2407 aodpid->SetIntegratedTimes(times);
2409 // Float_t tzeroTrack = fESDpid->GetTOFResponse().GetStartTime(track->P());
2410 // aodpid->SetTOFsignal(track->GetTOFsignal()-tzeroTrack);
2411 aodpid->SetTOFsignal(track->GetTOFsignal());
2414 for (Int_t iMass=0; iMass<5; iMass++){
2415 // tofRes[iMass]=(Double_t)fESDpid->GetTOFResponse().GetExpectedSigma(track->P(), times[iMass], AliPID::ParticleMass(iMass));
2416 tofRes[iMass]=0; //backward compatibility
2418 aodpid->SetTOFpidResolution(tofRes);
2419 //aodpid->SetHMPIDsignal(0); // set to zero for compression but it will be removed later
2422 Double_t AliAnalysisTaskESDfilter::Chi2perNDF(AliESDtrack* track)
2424 // Calculate chi2 per ndf for track
2426 Int_t nClustersTPC = track->GetTPCNcls();
2427 if ( nClustersTPC > 5) {
2428 return (track->GetTPCchi2()/Float_t(nClustersTPC - 5));
2434 //______________________________________________________________________________
2435 void AliAnalysisTaskESDfilter::Terminate(Option_t */*option*/)
2437 // Terminate analysis
2439 if (fDebug > 1) printf("AnalysisESDfilter: Terminate() \n");
2442 //______________________________________________________________________________
2443 void AliAnalysisTaskESDfilter::PrintMCInfo(AliStack *pStack,Int_t label)
2446 if (!pStack) return;
2447 label = TMath::Abs(label);
2448 TParticle *part = pStack->Particle(label);
2449 Printf("########################");
2450 Printf("%s:%d %d UniqueID %d PDG %d P %3.3f",(char*)__FILE__,__LINE__,label,part->GetUniqueID(),part->GetPdgCode(),part->P());
2452 TParticle* mother = part;
2453 Int_t imo = part->GetFirstMother();
2454 Int_t nprim = pStack->GetNprimary();
2455 // while((imo >= nprim) && (mother->GetUniqueID() == 4)) {
2456 while((imo >= nprim)) {
2457 mother = pStack->Particle(imo);
2458 Printf("Mother %s:%d Label %d UniqueID %d PDG %d P %3.3f",(char*)__FILE__,__LINE__,imo,mother->GetUniqueID(),mother->GetPdgCode(),mother->P());
2460 imo = mother->GetFirstMother();
2462 Printf("########################");
2465 //______________________________________________________________________________
2466 void AliAnalysisTaskESDfilter::CopyCaloProps(AliESDtrack *tre, AliAODTrack *tra)
2468 // Copy calo properties from ESD track to AOD track
2469 tra->SetTrackPhiEtaPtOnEMCal(tre->GetTrackPhiOnEMCal(),tre->GetTrackEtaOnEMCal(),tre->GetTrackPtOnEMCal());
2470 if (tre->IsEMCAL()) tra->SetEMCALcluster(tre->GetEMCALcluster());
2471 if (tre->IsPHOS()) tra->SetPHOScluster(tre->GetPHOScluster());
2474 //______________________________________________________________________________
2475 void AliAnalysisTaskESDfilter::SetRefitVertexTracks(Int_t algo, Double_t* cuts)
2477 // request vertexTrack reprocessing from ESDtracks
2478 // if algo>=0 and cuts==0 then algo is interpreted as the algorithm ID to be run with default cuts
2479 // otherwise it is number of cuts to digest
2480 fRefitVertexTracks = algo;
2481 if (algo>0 && cuts) {
2482 fRefitVertexTracksCuts = new Double_t[fRefitVertexTracks];
2483 for (int i=fRefitVertexTracks;i--;) fRefitVertexTracksCuts[i] = cuts[i];
2484 fRefitVertexTracksNCuts = fRefitVertexTracks;