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 **************************************************************************/
22 #include <TParameter.h>
24 #include <TParticle.h>
27 #include "AliAnalysisTaskESDfilter.h"
28 #include "AliAnalysisManager.h"
29 #include "AliESDEvent.h"
30 #include "AliESDRun.h"
32 #include "AliAODEvent.h"
33 #include "AliMCEvent.h"
34 #include "AliMCEventHandler.h"
35 #include "AliESDInputHandler.h"
36 #include "AliAODHandler.h"
37 #include "AliAODMCParticle.h"
38 #include "AliAnalysisFilter.h"
39 #include "AliESDMuonTrack.h"
40 #include "AliESDVertex.h"
41 #include "AliCentrality.h"
42 #include "AliEventplane.h"
44 #include "AliESDkink.h"
45 #include "AliESDcascade.h"
46 #include "AliESDPmdTrack.h"
47 #include "AliESDCaloCluster.h"
48 #include "AliESDCaloCells.h"
49 #include "AliMultiplicity.h"
51 #include "AliCodeTimer.h"
52 #include "AliESDtrackCuts.h"
53 #include "AliESDpid.h"
54 #include "AliAODHMPIDrings.h"
55 #include "AliV0vertexer.h"
56 #include "AliCascadeVertexer.h"
57 #include "Riostream.h"
58 #include "AliExternalTrackParam.h"
59 #include "AliTrackerBase.h"
61 #include "AliTPCdEdxInfo.h"
65 ClassImp(AliAnalysisTaskESDfilter)
67 ////////////////////////////////////////////////////////////////////////
69 AliAnalysisTaskESDfilter::AliAnalysisTaskESDfilter():
77 fEnableFillAOD(kTRUE),
86 fNumberOfPositiveTracks(0),
91 fOldESDformat(kFALSE),
93 fTPCConstrainedFilterMask(0),
94 fHybridFilterMaskTPCCG(0),
95 fWriteHybridTPCCOnly(kFALSE),
96 fGlobalConstrainedFilterMask(0),
97 fHybridFilterMaskGCG(0),
98 fWriteHybridGCOnly(kFALSE),
99 fIsVZEROEnabled(kTRUE),
100 fIsTZEROEnabled(kTRUE),
101 fIsZDCEnabled(kTRUE),
102 fIsHMPIDEnabled(kTRUE),
103 fIsV0CascadeRecoEnabled(kFALSE),
104 fAreCascadesEnabled(kTRUE),
105 fAreV0sEnabled(kTRUE),
106 fAreKinksEnabled(kTRUE),
107 fAreTracksEnabled(kTRUE),
108 fArePmdClustersEnabled(kTRUE),
109 fAreCaloClustersEnabled(kTRUE),
110 fAreEMCALCellsEnabled(kTRUE),
111 fArePHOSCellsEnabled(kTRUE),
112 fAreEMCALTriggerEnabled(kTRUE),
113 fArePHOSTriggerEnabled(kTRUE),
114 fAreTrackletsEnabled(kTRUE),
117 fTPCaloneTrackCuts(0),
118 fDoPropagateTrackToEMCal(kTRUE)
120 // Default constructor
121 fV0Cuts[0] = 33. ; // max allowed chi2
122 fV0Cuts[1] = 0.1 ; // min allowed impact parameter for the 1st daughter
123 fV0Cuts[2] = 0.1 ; // min allowed impact parameter for the 2nd daughter
124 fV0Cuts[3] = 1. ; // max allowed DCA between the daughter tracks
125 fV0Cuts[4] = .998; // min allowed cosine of V0's pointing angle
126 fV0Cuts[5] = 0.9 ; // min radius of the fiducial volume
127 fV0Cuts[6] = 100. ; // max radius of the fiducial volume
129 fCascadeCuts[0] = 33. ; // max allowed chi2 (same as PDC07)
130 fCascadeCuts[1] = 0.05 ; // min allowed V0 impact parameter
131 fCascadeCuts[2] = 0.008; // "window" around the Lambda mass
132 fCascadeCuts[3] = 0.03 ; // min allowed bachelor's impact parameter
133 fCascadeCuts[4] = 0.3 ; // max allowed DCA between the V0 and the bachelor
134 fCascadeCuts[5] = 0.999; // min allowed cosine of the cascade pointing angle
135 fCascadeCuts[6] = 0.9 ; // min radius of the fiducial volume
136 fCascadeCuts[7] = 100. ; // max radius of the fiducial volume
139 //______________________________________________________________________________
140 AliAnalysisTaskESDfilter::AliAnalysisTaskESDfilter(const char* name):
141 AliAnalysisTaskSE(name),
148 fEnableFillAOD(kTRUE),
157 fNumberOfPositiveTracks(0),
159 fNumberOfVertices(0),
160 fNumberOfCascades(0),
162 fOldESDformat(kFALSE),
164 fTPCConstrainedFilterMask(0),
165 fHybridFilterMaskTPCCG(0),
166 fWriteHybridTPCCOnly(kFALSE),
167 fGlobalConstrainedFilterMask(0),
168 fHybridFilterMaskGCG(0),
169 fWriteHybridGCOnly(kFALSE),
170 fIsVZEROEnabled(kTRUE),
171 fIsTZEROEnabled(kTRUE),
172 fIsZDCEnabled(kTRUE),
173 fIsHMPIDEnabled(kTRUE),
174 fIsV0CascadeRecoEnabled(kFALSE),
175 fAreCascadesEnabled(kTRUE),
176 fAreV0sEnabled(kTRUE),
177 fAreKinksEnabled(kTRUE),
178 fAreTracksEnabled(kTRUE),
179 fArePmdClustersEnabled(kTRUE),
180 fAreCaloClustersEnabled(kTRUE),
181 fAreEMCALCellsEnabled(kTRUE),
182 fArePHOSCellsEnabled(kTRUE),
183 fAreEMCALTriggerEnabled(kTRUE),
184 fArePHOSTriggerEnabled(kTRUE),
185 fAreTrackletsEnabled(kTRUE),
188 fTPCaloneTrackCuts(0),
189 fDoPropagateTrackToEMCal(kTRUE)
193 fV0Cuts[0] = 33. ; // max allowed chi2
194 fV0Cuts[1] = 0.1 ; // min allowed impact parameter for the 1st daughter
195 fV0Cuts[2] = 0.1 ; // min allowed impact parameter for the 2nd daughter
196 fV0Cuts[3] = 1. ; // max allowed DCA between the daughter tracks
197 fV0Cuts[4] = .998; // min allowed cosine of V0's pointing angle
198 fV0Cuts[5] = 0.9 ; // min radius of the fiducial volume
199 fV0Cuts[6] = 100. ; // max radius of the fiducial volume
201 fCascadeCuts[0] = 33. ; // max allowed chi2 (same as PDC07)
202 fCascadeCuts[1] = 0.05 ; // min allowed V0 impact parameter
203 fCascadeCuts[2] = 0.008; // "window" around the Lambda mass
204 fCascadeCuts[3] = 0.03 ; // min allowed bachelor's impact parameter
205 fCascadeCuts[4] = 0.3 ; // max allowed DCA between the V0 and the bachelor
206 fCascadeCuts[5] = 0.999; // min allowed cosine of the cascade pointing angle
207 fCascadeCuts[6] = 0.9 ; // min radius of the fiducial volume
208 fCascadeCuts[7] = 100. ; // max radius of the fiducial volume
213 AliAnalysisTaskESDfilter::~AliAnalysisTaskESDfilter(){
214 if(fIsPidOwner) delete fESDpid;
216 //______________________________________________________________________________
217 void AliAnalysisTaskESDfilter::UserCreateOutputObjects()
220 // Create Output Objects conenct filter to outputtree
224 OutputTree()->GetUserInfo()->Add(fTrackFilter);
228 AliError("No OutputTree() for adding the track filter");
230 fTPCaloneTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
233 //______________________________________________________________________________
234 void AliAnalysisTaskESDfilter::Init()
237 if (fDebug > 1) AliInfo("Init() \n");
238 // Call configuration file
241 //______________________________________________________________________________
242 Bool_t AliAnalysisTaskESDfilter::Notify()
245 AddMetadataToUserInfo();
249 //______________________________________________________________________________
250 Bool_t AliAnalysisTaskESDfilter::AddMetadataToUserInfo()
252 // Copy metadata to AOD user info.
253 static Bool_t copyFirst = kFALSE;
255 AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
257 AliError("AliAnalysisTaskESDfilter::AddMetadataToUserInfo() : No analysis manager !");
260 TTree *esdTree = mgr->GetTree()->GetTree();
261 if (!esdTree) return kFALSE;
262 TNamed *alirootVersion = (TNamed*)esdTree->GetUserInfo()->FindObject("alirootVersion");
263 if (!alirootVersion) return kFALSE;
264 AliAODHandler *aodHandler = dynamic_cast<AliAODHandler*>(mgr->GetOutputEventHandler());
265 if (!aodHandler) return kFALSE;
266 TTree *aodTree = aodHandler->GetTree();
267 if (!aodTree) return kFALSE;
268 aodTree->GetUserInfo()->Add(new TNamed(*alirootVersion));
274 //______________________________________________________________________________
275 void AliAnalysisTaskESDfilter::PrintTask(Option_t *option, Int_t indent) const
277 // Print selection task information
280 AliAnalysisTaskSE::PrintTask(option,indent);
282 TString spaces(' ',indent+3);
284 cout << spaces.Data() << Form("Cascades are %s",fAreCascadesEnabled ? "ENABLED":"DISABLED") << endl;
285 cout << spaces.Data() << Form("V0s are %s",fAreV0sEnabled ? "ENABLED":"DISABLED") << endl;
286 cout << spaces.Data() << Form("Kinks are %s",fAreKinksEnabled ? "ENABLED":"DISABLED") << endl;
287 cout << spaces.Data() << Form("Tracks are %s",fAreTracksEnabled ? "ENABLED":"DISABLED") << endl;
288 cout << spaces.Data() << Form("PmdClusters are %s",fArePmdClustersEnabled ? "ENABLED":"DISABLED") << endl;
289 cout << spaces.Data() << Form("CaloClusters are %s",fAreCaloClustersEnabled ? "ENABLED":"DISABLED") << endl;
290 cout << spaces.Data() << Form("EMCAL cells are %s",fAreEMCALCellsEnabled ? "ENABLED":"DISABLED") << endl;
291 cout << spaces.Data() << Form("EMCAL triggers are %s",fAreEMCALTriggerEnabled ? "ENABLED":"DISABLED") << endl;
292 cout << spaces.Data() << Form("PHOS triggers are %s",fArePHOSTriggerEnabled ? "ENABLED":"DISABLED") << endl;
293 cout << spaces.Data() << Form("Tracklets are %s",fAreTrackletsEnabled ? "ENABLED":"DISABLED") << endl;
294 cout << spaces.Data() << Form("PropagateTrackToEMCal is %s", fDoPropagateTrackToEMCal ? "ENABLED":"DISABLED") << endl;
297 //______________________________________________________________________________
298 void AliAnalysisTaskESDfilter::UserExec(Option_t */*option*/)
300 // Execute analysis for current event
303 Long64_t ientry = Entry();
306 printf("Filter: Analysing event # %5d\n", (Int_t) ientry);
307 if (fHighPthreshold == 0) AliInfo("detector PID signals are stored in each track");
308 if (!fPtshape) AliInfo("detector PID signals are not stored below the pt threshold");
310 // Filters must explicitely enable AOD filling in their UserExec (AG)
311 if (!AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()) AliFatal("Cannot run ESD filter without an output event handler");
313 AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillAOD(kTRUE);
314 AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillExtension(kTRUE);
319 //______________________________________________________________________________
320 TClonesArray& AliAnalysisTaskESDfilter::Cascades()
322 return *(AODEvent()->GetCascades());
325 //______________________________________________________________________________
326 TClonesArray& AliAnalysisTaskESDfilter::Tracks()
328 return *(AODEvent()->GetTracks());
331 //______________________________________________________________________________
332 TClonesArray& AliAnalysisTaskESDfilter::V0s()
334 return *(AODEvent()->GetV0s());
337 //______________________________________________________________________________
338 TClonesArray& AliAnalysisTaskESDfilter::Vertices()
340 return *(AODEvent()->GetVertices());
343 //______________________________________________________________________________
344 AliAODHeader* AliAnalysisTaskESDfilter::ConvertHeader(const AliESDEvent& esd)
346 // Convert header information
348 AliCodeTimerAuto("",0);
350 AliAODHeader* header = AODEvent()->GetHeader();
352 header->SetRunNumber(esd.GetRunNumber());
353 header->SetOfflineTrigger(fInputHandler->IsEventSelected()); // propagate the decision of the physics selection
355 TTree* tree = fInputHandler->GetTree();
357 TFile* file = tree->GetCurrentFile();
358 if (file) header->SetESDFileName(file->GetName());
362 header->SetBunchCrossNumber(0);
363 header->SetOrbitNumber(0);
364 header->SetPeriodNumber(0);
365 header->SetEventType(0);
366 header->SetMuonMagFieldScale(-999.);
367 header->SetCentrality(0);
368 header->SetEventplane(0);
370 header->SetBunchCrossNumber(esd.GetBunchCrossNumber());
371 header->SetOrbitNumber(esd.GetOrbitNumber());
372 header->SetPeriodNumber(esd.GetPeriodNumber());
373 header->SetEventType(esd.GetEventType());
375 header->SetEventNumberESDFile(esd.GetHeader()->GetEventNumberInFile());
376 if(const_cast<AliESDEvent&>(esd).GetCentrality()){
377 header->SetCentrality(const_cast<AliESDEvent&>(esd).GetCentrality());
380 header->SetCentrality(0);
382 if(const_cast<AliESDEvent&>(esd).GetEventplane()){
383 header->SetEventplane(const_cast<AliESDEvent&>(esd).GetEventplane());
386 header->SetEventplane(0);
391 header->SetFiredTriggerClasses(esd.GetFiredTriggerClasses());
392 header->SetTriggerMask(esd.GetTriggerMask());
393 header->SetTriggerCluster(esd.GetTriggerCluster());
394 header->SetL0TriggerInputs(esd.GetHeader()->GetL0TriggerInputs());
395 header->SetL1TriggerInputs(esd.GetHeader()->GetL1TriggerInputs());
396 header->SetL2TriggerInputs(esd.GetHeader()->GetL2TriggerInputs());
398 header->SetMagneticField(esd.GetMagneticField());
399 header->SetMuonMagFieldScale(esd.GetCurrentDip()/6000.);
400 header->SetZDCN1Energy(esd.GetZDCN1Energy());
401 header->SetZDCP1Energy(esd.GetZDCP1Energy());
402 header->SetZDCN2Energy(esd.GetZDCN2Energy());
403 header->SetZDCP2Energy(esd.GetZDCP2Energy());
404 header->SetZDCEMEnergy(esd.GetZDCEMEnergy(0),esd.GetZDCEMEnergy(1));
406 // ITS Cluster Multiplicty
407 const AliMultiplicity *mult = esd.GetMultiplicity();
408 for (Int_t ilay = 0; ilay < 6; ilay++) header->SetITSClusters(ilay, mult->GetNumberOfITSClusters(ilay));
410 // TPC only Reference Multiplicty
411 Int_t refMult = fTPCaloneTrackCuts ? (Short_t)fTPCaloneTrackCuts->GetReferenceMultiplicity(&esd, kTRUE) : -1;
412 header->SetTPConlyRefMultiplicity(refMult);
415 Float_t diamxy[2]={esd.GetDiamondX(),esd.GetDiamondY()};
417 esd.GetDiamondCovXY(diamcov);
418 header->SetDiamond(diamxy,diamcov);
419 header->SetDiamondZ(esd.GetDiamondZ(),esd.GetSigma2DiamondZ());
421 // VZERO channel equalization factors for event-plane reconstruction
422 header->SetVZEROEqFactors(esd.GetVZEROEqFactors());
424 // T0 Resolution information
425 const AliESDRun* esdRun = esd.GetESDRun();
426 for (Int_t i=0;i<AliESDRun::kT0spreadSize;i++) header->SetT0spread(i,esdRun->GetT0spread(i));
431 //______________________________________________________________________________
432 void AliAnalysisTaskESDfilter::ConvertCascades(const AliESDEvent& esd)
435 // Convert the cascades part of the ESD.
436 // Return the number of cascades
438 AliCodeTimerAuto("",0);
440 // Create vertices starting from the most complex objects
443 const AliESDVertex* vtx = esd.GetPrimaryVertex();
444 Double_t pos[3] = { 0. };
445 Double_t covVtx[6] = { 0. };
446 Double_t momBach[3]={0.};
447 Double_t covTr[21]={0.};
448 Double_t pid[10]={0.};
449 AliAODPid* detpid(0x0);
450 AliAODVertex* vV0FromCascade(0x0);
451 AliAODv0* aodV0(0x0);
452 AliAODcascade* aodCascade(0x0);
453 AliAODTrack* aodTrack(0x0);
454 Double_t momPos[3]={0.};
455 Double_t momNeg[3] = { 0. };
456 Double_t momPosAtV0vtx[3]={0.};
457 Double_t momNegAtV0vtx[3]={0.};
459 TClonesArray& verticesArray = Vertices();
460 TClonesArray& tracksArray = Tracks();
461 TClonesArray& cascadesArray = Cascades();
463 // Cascades (Modified by A.Maire - February 2009)
464 for (Int_t nCascade = 0; nCascade < esd.GetNumberOfCascades(); ++nCascade) {
468 AliESDcascade *esdCascade = esd.GetCascade(nCascade);
469 Int_t idxPosFromV0Dghter = esdCascade->GetPindex();
470 Int_t idxNegFromV0Dghter = esdCascade->GetNindex();
471 Int_t idxBachFromCascade = esdCascade->GetBindex();
473 AliESDtrack *esdCascadePos = esd.GetTrack( idxPosFromV0Dghter);
474 AliESDtrack *esdCascadeNeg = esd.GetTrack( idxNegFromV0Dghter);
475 AliESDtrack *esdCascadeBach = esd.GetTrack( idxBachFromCascade);
477 // Identification of the V0 within the esdCascade (via both daughter track indices)
478 AliESDv0 * currentV0 = 0x0;
479 Int_t idxV0FromCascade = -1;
481 for (Int_t iV0=0; iV0<esd.GetNumberOfV0s(); ++iV0) {
483 currentV0 = esd.GetV0(iV0);
484 Int_t posCurrentV0 = currentV0->GetPindex();
485 Int_t negCurrentV0 = currentV0->GetNindex();
487 if (posCurrentV0==idxPosFromV0Dghter && negCurrentV0==idxNegFromV0Dghter) {
488 idxV0FromCascade = iV0;
493 if(idxV0FromCascade < 0){
494 printf("Cascade - no matching for the V0 (index V0 = -1) ! Skip ... \n");
496 }// a priori, useless check, but safer ... in case of pb with tracks "out of bounds"
498 AliESDv0 *esdV0FromCascade = esd.GetV0(idxV0FromCascade);
500 // 1 - Cascade selection
502 // AliESDVertex *esdPrimVtx = new AliESDVertex(*(esd.GetPrimaryVertex()));
503 // TList cascadeObjects;
504 // cascadeObjects.AddAt(esdV0FromCascade, 0);
505 // cascadeObjects.AddAt(esdCascadePos, 1);
506 // cascadeObjects.AddAt(esdCascadeNeg, 2);
507 // cascadeObjects.AddAt(esdCascade, 3);
508 // cascadeObjects.AddAt(esdCascadeBach, 4);
509 // cascadeObjects.AddAt(esdPrimVtx, 5);
511 // UInt_t selectCascade = 0;
512 // if (fCascadeFilter) {
513 // // selectCascade = fCascadeFilter->IsSelected(&cascadeObjects);
514 // // FIXME AliESDCascadeCuts to be implemented ...
516 // // Here we may encounter a moot point at the V0 level
517 // // between the cascade selections and the V0 ones :
518 // // the V0 selected along with the cascade (secondary V0) may
519 // // usually be removed from the dedicated V0 selections (prim V0) ...
520 // // -> To be discussed !
522 // // this is a little awkward but otherwise the
523 // // list wants to access the pointer (delete it)
524 // // again when going out of scope
525 // delete cascadeObjects.RemoveAt(5); // esdPrimVtx created via copy construct
527 // if (!selectCascade)
531 // delete cascadeObjects.RemoveAt(5); // esdPrimVtx created via copy construct
535 // 2 - Add the cascade vertex
537 esdCascade->GetXYZcascade(pos[0], pos[1], pos[2]);
538 esdCascade->GetPosCovXi(covVtx);
539 chi2 = esdCascade->GetChi2Xi();
541 AliAODVertex *vCascade = new(verticesArray[fNumberOfVertices++]) AliAODVertex( pos,
543 chi2, // FIXME = Chi2/NDF will be needed
546 AliAODVertex::kCascade);
547 fPrimaryVertex->AddDaughter(vCascade);
550 // printf("---- Cascade / Cascade Vertex (AOD) : \n");
551 // vCascade->Print();
554 // if(esd.GetTOFHeader() && fIsPidOwner) fESDpid->SetTOFResponse(const_cast<AliESDEvent*>(&esd), (AliESDpid::EStartTimeType_t)fTimeZeroType); //in case of AOD production starting form LHC10e without Tender.
557 // 3 - Add the bachelor track from the cascade
559 if (!fUsedTrack[idxBachFromCascade]) {
561 esdCascadeBach->GetPxPyPz(momBach);
562 esdCascadeBach->GetXYZ(pos);
563 esdCascadeBach->GetCovarianceXYZPxPyPz(covTr);
564 esdCascadeBach->GetESDpid(pid);
566 fUsedTrack[idxBachFromCascade] = kTRUE;
567 UInt_t selectInfo = 0;
568 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadeBach);
569 if (fMChandler) fMChandler->SelectParticle(esdCascadeBach->GetLabel());
570 aodTrack = new(tracksArray[fNumberOfTracks++]) AliAODTrack(esdCascadeBach->GetID(),
571 esdCascadeBach->GetLabel(),
575 kFALSE, // Why kFALSE for "isDCA" ? FIXME
577 (Short_t)esdCascadeBach->GetSign(),
578 esdCascadeBach->GetITSClusterMap(),
581 kTRUE, // usedForVtxFit = kFALSE ? FIXME
582 vtx->UsesTrack(esdCascadeBach->GetID()),
583 AliAODTrack::kSecondary,
585 aodTrack->SetTPCFitMap(esdCascadeBach->GetTPCFitMap());
586 aodTrack->SetTPCClusterMap(esdCascadeBach->GetTPCClusterMap());
587 aodTrack->SetTPCSharedMap (esdCascadeBach->GetTPCSharedMap());
588 aodTrack->SetChi2perNDF(Chi2perNDF(esdCascadeBach));
589 aodTrack->SetTPCPointsF(esdCascadeBach->GetTPCNclsF());
590 aodTrack->SetTPCNCrossedRows(UShort_t(esdCascadeBach->GetTPCCrossedRows()));
591 fAODTrackRefs->AddAt(aodTrack,idxBachFromCascade);
593 if (esdCascadeBach->GetSign() > 0) ++fNumberOfPositiveTracks;
594 aodTrack->ConvertAliPIDtoAODPID();
595 aodTrack->SetFlags(esdCascadeBach->GetStatus());
596 SetAODPID(esdCascadeBach,aodTrack,detpid);
599 aodTrack = static_cast<AliAODTrack*>( fAODTrackRefs->At(idxBachFromCascade) );
602 vCascade->AddDaughter(aodTrack);
605 // printf("---- Cascade / bach dghter : \n");
606 // aodTrack->Print();
610 // 4 - Add the V0 from the cascade.
611 // = V0vtx + both pos and neg daughter tracks + the aodV0 itself
614 if ( !fUsedV0[idxV0FromCascade] ) {
615 // 4.A - if VO structure hasn't been created yet
617 // 4.A.1 - Create the V0 vertex of the cascade
619 esdV0FromCascade->GetXYZ(pos[0], pos[1], pos[2]);
620 esdV0FromCascade->GetPosCov(covVtx);
621 chi2 = esdV0FromCascade->GetChi2V0(); // = chi2/NDF since NDF = 2*2-3 ?
623 vV0FromCascade = new(verticesArray[fNumberOfVertices++]) AliAODVertex(pos,
627 idxV0FromCascade, //id of ESDv0
630 // one V0 can be used by several cascades.
631 // So, one AOD V0 vtx can have several parent vtx.
632 // This is not directly allowed by AliAODvertex.
633 // Setting the parent vtx (here = param "vCascade") doesn't lead to a crash
634 // but to a problem of consistency within AODEvent.
635 // -> See below paragraph 4.B, for the proposed treatment of such a case.
637 // Add the vV0FromCascade to the aodVOVtxRefs
638 fAODV0VtxRefs->AddAt(vV0FromCascade,idxV0FromCascade);
641 // 4.A.2 - Add the positive tracks from the V0
643 esdCascadePos->GetPxPyPz(momPos);
644 esdCascadePos->GetXYZ(pos);
645 esdCascadePos->GetCovarianceXYZPxPyPz(covTr);
646 esdCascadePos->GetESDpid(pid);
649 if (!fUsedTrack[idxPosFromV0Dghter]) {
650 fUsedTrack[idxPosFromV0Dghter] = kTRUE;
652 UInt_t selectInfo = 0;
653 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadePos);
654 if(fMChandler) fMChandler->SelectParticle(esdCascadePos->GetLabel());
655 aodTrack = new(tracksArray[fNumberOfTracks++])
656 AliAODTrack( esdCascadePos->GetID(),
657 esdCascadePos->GetLabel(),
661 kFALSE, // Why kFALSE for "isDCA" ? FIXME
663 (Short_t)esdCascadePos->GetSign(),
664 esdCascadePos->GetITSClusterMap(),
667 kTRUE, // usedForVtxFit = kFALSE ? FIXME
668 vtx->UsesTrack(esdCascadePos->GetID()),
669 AliAODTrack::kSecondary,
671 aodTrack->SetTPCFitMap(esdCascadePos->GetTPCFitMap());
672 aodTrack->SetTPCClusterMap(esdCascadePos->GetTPCClusterMap());
673 aodTrack->SetTPCSharedMap (esdCascadePos->GetTPCSharedMap());
674 aodTrack->SetChi2perNDF(Chi2perNDF(esdCascadePos));
675 aodTrack->SetTPCPointsF(esdCascadePos->GetTPCNclsF());
676 aodTrack->SetTPCNCrossedRows(UShort_t(esdCascadePos->GetTPCCrossedRows()));
677 fAODTrackRefs->AddAt(aodTrack,idxPosFromV0Dghter);
679 if (esdCascadePos->GetSign() > 0) ++fNumberOfPositiveTracks;
680 aodTrack->ConvertAliPIDtoAODPID();
681 aodTrack->SetFlags(esdCascadePos->GetStatus());
682 SetAODPID(esdCascadePos,aodTrack,detpid);
685 aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(idxPosFromV0Dghter));
687 vV0FromCascade->AddDaughter(aodTrack);
690 // 4.A.3 - Add the negative tracks from the V0
692 esdCascadeNeg->GetPxPyPz(momNeg);
693 esdCascadeNeg->GetXYZ(pos);
694 esdCascadeNeg->GetCovarianceXYZPxPyPz(covTr);
695 esdCascadeNeg->GetESDpid(pid);
698 if (!fUsedTrack[idxNegFromV0Dghter]) {
699 fUsedTrack[idxNegFromV0Dghter] = kTRUE;
701 UInt_t selectInfo = 0;
702 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadeNeg);
703 if(fMChandler)fMChandler->SelectParticle(esdCascadeNeg->GetLabel());
704 aodTrack = new(tracksArray[fNumberOfTracks++]) AliAODTrack( esdCascadeNeg->GetID(),
705 esdCascadeNeg->GetLabel(),
709 kFALSE, // Why kFALSE for "isDCA" ? FIXME
711 (Short_t)esdCascadeNeg->GetSign(),
712 esdCascadeNeg->GetITSClusterMap(),
715 kTRUE, // usedForVtxFit = kFALSE ? FIXME
716 vtx->UsesTrack(esdCascadeNeg->GetID()),
717 AliAODTrack::kSecondary,
719 aodTrack->SetTPCFitMap(esdCascadeNeg->GetTPCFitMap());
720 aodTrack->SetTPCClusterMap(esdCascadeNeg->GetTPCClusterMap());
721 aodTrack->SetTPCSharedMap (esdCascadeNeg->GetTPCSharedMap());
722 aodTrack->SetChi2perNDF(Chi2perNDF(esdCascadeNeg));
723 aodTrack->SetTPCPointsF(esdCascadeNeg->GetTPCNclsF());
724 aodTrack->SetTPCNCrossedRows(UShort_t(esdCascadeNeg->GetTPCCrossedRows()));
725 fAODTrackRefs->AddAt(aodTrack,idxNegFromV0Dghter);
727 if (esdCascadeNeg->GetSign() > 0) ++fNumberOfPositiveTracks;
728 aodTrack->ConvertAliPIDtoAODPID();
729 aodTrack->SetFlags(esdCascadeNeg->GetStatus());
730 SetAODPID(esdCascadeNeg,aodTrack,detpid);
733 aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(idxNegFromV0Dghter));
736 vV0FromCascade->AddDaughter(aodTrack);
739 // 4.A.4 - Add the V0 from cascade to the V0 array
741 Double_t dcaV0Daughters = esdV0FromCascade->GetDcaV0Daughters();
742 Double_t dcaV0ToPrimVertex = esdV0FromCascade->GetD( esd.GetPrimaryVertex()->GetX(),
743 esd.GetPrimaryVertex()->GetY(),
744 esd.GetPrimaryVertex()->GetZ() );
745 esdV0FromCascade->GetPPxPyPz( momPosAtV0vtx[0],momPosAtV0vtx[1],momPosAtV0vtx[2] );
746 esdV0FromCascade->GetNPxPyPz( momNegAtV0vtx[0],momNegAtV0vtx[1],momNegAtV0vtx[2] );
748 Double_t dcaDaughterToPrimVertex[2] = { 999., 999.}; // ..[0] = DCA in (x,y) for Pos and ..[1] = Neg
749 dcaDaughterToPrimVertex[0] = TMath::Abs(esdCascadePos->GetD( esd.GetPrimaryVertex()->GetX(),
750 esd.GetPrimaryVertex()->GetY(),
751 esd.GetMagneticField()) );
752 dcaDaughterToPrimVertex[1] = TMath::Abs(esdCascadeNeg->GetD( esd.GetPrimaryVertex()->GetX(),
753 esd.GetPrimaryVertex()->GetY(),
754 esd.GetMagneticField()) );
756 aodV0 = new(V0s()[fNumberOfV0s++]) AliAODv0( vV0FromCascade,
761 dcaDaughterToPrimVertex);
762 // set the aod v0 on-the-fly status
763 aodV0->SetOnFlyStatus(esdV0FromCascade->GetOnFlyStatus());
765 // Add the aodV0 to the aodVORefs
766 fAODV0Refs->AddAt(aodV0,idxV0FromCascade);
768 fUsedV0[idxV0FromCascade] = kTRUE;
771 // 4.B - if V0 structure already used
774 // one V0 can be used by several cascades (frequent in PbPb evts) :
775 // same V0 which used but attached to different bachelor tracks
776 // -> aodVORefs and fAODV0VtxRefs are needed.
777 // Goal : avoid a redundancy of the info in "Vertices" and "v0s" clones array.
779 vV0FromCascade = static_cast<AliAODVertex*>( fAODV0VtxRefs->At(idxV0FromCascade) );
780 aodV0 = static_cast<AliAODv0*> ( fAODV0Refs ->At(idxV0FromCascade) );
782 // - Treatment of the parent for such a "re-used" V0 :
783 // Insert the cascade that reuses the V0 vertex in the lineage chain
784 // Before : vV0 -> vCascade1 -> vPrimary
785 // - Hyp : cascade2 uses the same V0 as cascade1
786 // After : vV0 -> vCascade2 -> vCascade1 -> vPrimary
788 AliAODVertex *vCascadePreviousParent = static_cast<AliAODVertex*> (vV0FromCascade->GetParent());
789 vV0FromCascade->SetParent(vCascade);
790 vCascade ->SetParent(vCascadePreviousParent);
793 // printf("---- Cascade / Lineage insertion\n"
794 // "Parent of V0 vtx = Cascade vtx %p\n"
795 // "Parent of the cascade vtx = Cascade vtx %p\n"
796 // "Parent of the parent cascade vtx = Cascade vtx %p\n",
797 // static_cast<void*> (vV0FromCascade->GetParent()),
798 // static_cast<void*> (vCascade->GetParent()),
799 // static_cast<void*> (vCascadePreviousParent->GetParent()) );
801 }// end if V0 structure already used
804 // printf("---- Cascade / V0 vertex: \n");
805 // vV0FromCascade->Print();
809 // printf("---- Cascade / pos dghter : \n");
810 // aodTrack->Print();
811 // printf("---- Cascade / neg dghter : \n");
812 // aodTrack->Print();
813 // printf("---- Cascade / aodV0 : \n");
817 // In any case (used V0 or not), add the V0 vertex to the cascade one.
818 vCascade->AddDaughter(vV0FromCascade);
821 // 5 - Add the primary track of the cascade (if any)
824 // 6 - Add the cascade to the AOD array of cascades
826 Double_t dcaBachToPrimVertexXY = TMath::Abs(esdCascadeBach->GetD(esd.GetPrimaryVertex()->GetX(),
827 esd.GetPrimaryVertex()->GetY(),
828 esd.GetMagneticField()) );
830 Double_t momBachAtCascadeVtx[3]={0.};
832 esdCascade->GetBPxPyPz(momBachAtCascadeVtx[0], momBachAtCascadeVtx[1], momBachAtCascadeVtx[2]);
834 aodCascade = new(cascadesArray[fNumberOfCascades++]) AliAODcascade( vCascade,
835 esdCascade->Charge(),
836 esdCascade->GetDcaXiDaughters(),
838 // DCAXiToPrimVtx -> needs to be calculated ----|
839 // doesn't exist at ESD level;
840 // See AODcascade::DcaXiToPrimVertex(Double, Double, Double)
841 dcaBachToPrimVertexXY,
846 printf("---- Cascade / AOD cascade : \n\n");
847 aodCascade->PrintXi(fPrimaryVertex->GetX(), fPrimaryVertex->GetY(), fPrimaryVertex->GetZ());
850 } // end of the loop on cascades
852 Cascades().Expand(fNumberOfCascades);
855 //______________________________________________________________________________
856 void AliAnalysisTaskESDfilter::ConvertV0s(const AliESDEvent& esd)
858 // Access to the AOD container of V0s
860 AliCodeTimerAuto("",0);
866 Double_t pos[3] = { 0. };
868 Double_t covVtx[6] = { 0. };
869 Double_t momPos[3]={0.};
870 Double_t covTr[21]={0.};
871 Double_t pid[10]={0.};
872 AliAODTrack* aodTrack(0x0);
873 AliAODPid* detpid(0x0);
874 Double_t momNeg[3]={0.};
875 Double_t momPosAtV0vtx[3]={0.};
876 Double_t momNegAtV0vtx[3]={0.};
878 for (Int_t nV0 = 0; nV0 < esd.GetNumberOfV0s(); ++nV0)
880 if (fUsedV0[nV0]) continue; // skip if already added to the AOD
882 AliESDv0 *v0 = esd.GetV0(nV0);
883 Int_t posFromV0 = v0->GetPindex();
884 Int_t negFromV0 = v0->GetNindex();
888 AliESDVertex *esdVtx = new AliESDVertex(*(esd.GetPrimaryVertex()));
889 AliESDtrack *esdV0Pos = esd.GetTrack(posFromV0);
890 AliESDtrack *esdV0Neg = esd.GetTrack(negFromV0);
892 v0objects.AddAt(v0, 0);
893 v0objects.AddAt(esdV0Pos, 1);
894 v0objects.AddAt(esdV0Neg, 2);
895 v0objects.AddAt(esdVtx, 3);
898 selectV0 = fV0Filter->IsSelected(&v0objects);
899 // this is a little awkward but otherwise the
900 // list wants to access the pointer (delete it)
901 // again when going out of scope
902 delete v0objects.RemoveAt(3); // esdVtx created via copy construct
908 delete v0objects.RemoveAt(3); // esdVtx created via copy construct
912 v0->GetXYZ(pos[0], pos[1], pos[2]);
914 if (!fOldESDformat) {
915 chi2 = v0->GetChi2V0(); // = chi2/NDF since NDF = 2*2-3
916 v0->GetPosCov(covVtx);
919 for (Int_t i = 0; i < 6; i++) covVtx[i] = 0.;
924 new(Vertices()[fNumberOfVertices++]) AliAODVertex(pos,
930 fPrimaryVertex->AddDaughter(vV0);
933 // Add the positive tracks from the V0
936 esdV0Pos->GetPxPyPz(momPos);
937 esdV0Pos->GetXYZ(pos);
938 esdV0Pos->GetCovarianceXYZPxPyPz(covTr);
939 esdV0Pos->GetESDpid(pid);
941 const AliESDVertex *vtx = esd.GetPrimaryVertex();
943 if (!fUsedTrack[posFromV0]) {
944 fUsedTrack[posFromV0] = kTRUE;
945 UInt_t selectInfo = 0;
946 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdV0Pos);
947 if(fMChandler)fMChandler->SelectParticle(esdV0Pos->GetLabel());
948 aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdV0Pos->GetID(),
949 esdV0Pos->GetLabel(),
955 (Short_t)esdV0Pos->GetSign(),
956 esdV0Pos->GetITSClusterMap(),
959 kTRUE, // check if this is right
960 vtx->UsesTrack(esdV0Pos->GetID()),
961 AliAODTrack::kSecondary,
963 aodTrack->SetTPCFitMap(esdV0Pos->GetTPCFitMap());
964 aodTrack->SetTPCClusterMap(esdV0Pos->GetTPCClusterMap());
965 aodTrack->SetTPCSharedMap (esdV0Pos->GetTPCSharedMap());
966 aodTrack->SetChi2perNDF(Chi2perNDF(esdV0Pos));
967 aodTrack->SetTPCPointsF(esdV0Pos->GetTPCNclsF());
968 aodTrack->SetTPCNCrossedRows(UShort_t(esdV0Pos->GetTPCCrossedRows()));
969 fAODTrackRefs->AddAt(aodTrack,posFromV0);
970 // if (fDebug > 0) printf("-------------------Bo: pos track from original pt %.3f \n",aodTrack->Pt());
971 if (esdV0Pos->GetSign() > 0) ++fNumberOfPositiveTracks;
972 aodTrack->ConvertAliPIDtoAODPID();
973 aodTrack->SetFlags(esdV0Pos->GetStatus());
974 SetAODPID(esdV0Pos,aodTrack,detpid);
977 aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(posFromV0));
978 // if (fDebug > 0) printf("-------------------Bo pos track from refArray pt %.3f \n",aodTrack->Pt());
980 vV0->AddDaughter(aodTrack);
982 // Add the negative tracks from the V0
984 esdV0Neg->GetPxPyPz(momNeg);
985 esdV0Neg->GetXYZ(pos);
986 esdV0Neg->GetCovarianceXYZPxPyPz(covTr);
987 esdV0Neg->GetESDpid(pid);
989 if (!fUsedTrack[negFromV0]) {
990 fUsedTrack[negFromV0] = kTRUE;
991 UInt_t selectInfo = 0;
992 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdV0Neg);
993 if(fMChandler)fMChandler->SelectParticle(esdV0Neg->GetLabel());
994 aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdV0Neg->GetID(),
995 esdV0Neg->GetLabel(),
1001 (Short_t)esdV0Neg->GetSign(),
1002 esdV0Neg->GetITSClusterMap(),
1005 kTRUE, // check if this is right
1006 vtx->UsesTrack(esdV0Neg->GetID()),
1007 AliAODTrack::kSecondary,
1009 aodTrack->SetTPCFitMap(esdV0Neg->GetTPCFitMap());
1010 aodTrack->SetTPCClusterMap(esdV0Neg->GetTPCClusterMap());
1011 aodTrack->SetTPCSharedMap (esdV0Neg->GetTPCSharedMap());
1012 aodTrack->SetChi2perNDF(Chi2perNDF(esdV0Neg));
1013 aodTrack->SetTPCPointsF(esdV0Neg->GetTPCNclsF());
1014 aodTrack->SetTPCNCrossedRows(UShort_t(esdV0Neg->GetTPCCrossedRows()));
1016 fAODTrackRefs->AddAt(aodTrack,negFromV0);
1017 // if (fDebug > 0) printf("-------------------Bo: neg track from original pt %.3f \n",aodTrack->Pt());
1018 if (esdV0Neg->GetSign() > 0) ++fNumberOfPositiveTracks;
1019 aodTrack->ConvertAliPIDtoAODPID();
1020 aodTrack->SetFlags(esdV0Neg->GetStatus());
1021 SetAODPID(esdV0Neg,aodTrack,detpid);
1024 aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(negFromV0));
1025 // if (fDebug > 0) printf("-------------------Bo neg track from refArray pt %.3f \n",aodTrack->Pt());
1027 vV0->AddDaughter(aodTrack);
1030 // Add the V0 the V0 array as well
1032 Double_t dcaV0Daughters = v0->GetDcaV0Daughters();
1033 Double_t dcaV0ToPrimVertex = v0->GetD(esd.GetPrimaryVertex()->GetX(),
1034 esd.GetPrimaryVertex()->GetY(),
1035 esd.GetPrimaryVertex()->GetZ());
1037 v0->GetPPxPyPz(momPosAtV0vtx[0],momPosAtV0vtx[1],momPosAtV0vtx[2]);
1038 v0->GetNPxPyPz(momNegAtV0vtx[0],momNegAtV0vtx[1],momNegAtV0vtx[2]);
1040 Double_t dcaDaughterToPrimVertex[2] = { 999., 999.}; // ..[0] = DCA in (x,y) for Pos and ..[1] = Neg
1041 dcaDaughterToPrimVertex[0] = TMath::Abs(esdV0Pos->GetD( esd.GetPrimaryVertex()->GetX(),
1042 esd.GetPrimaryVertex()->GetY(),
1043 esd.GetMagneticField()) );
1044 dcaDaughterToPrimVertex[1] = TMath::Abs(esdV0Neg->GetD( esd.GetPrimaryVertex()->GetX(),
1045 esd.GetPrimaryVertex()->GetY(),
1046 esd.GetMagneticField()) );
1048 AliAODv0* aodV0 = new(V0s()[fNumberOfV0s++]) AliAODv0(vV0,
1053 dcaDaughterToPrimVertex);
1055 // set the aod v0 on-the-fly status
1056 aodV0->SetOnFlyStatus(v0->GetOnFlyStatus());
1057 }//End of loop on V0s
1059 V0s().Expand(fNumberOfV0s);
1062 //______________________________________________________________________________
1063 void AliAnalysisTaskESDfilter::ConvertTPCOnlyTracks(const AliESDEvent& esd)
1065 // Convert TPC only tracks
1066 // Here we have wo hybrid appraoch to remove fakes
1067 // ******* ITSTPC ********
1068 // Uses a cut on the ITS properties to select global tracks
1069 // which are than marked as HybdridITSTPC for the remainder
1070 // the TPC only tracks are flagged as HybridITSTPConly.
1071 // Note, in order not to get fakes back in the TPC cuts, one needs
1072 // two "ITS" cuts one tight (1) (to throw out fakes) and one lose (2) (to NOT flag the trakcs in the TPC only)
1073 // using cut number (3)
1074 // so fHybridFilterMask == (1)|(2) fTPCFilterMask = (3), Usercode needs to slect with mask = (1)|(3) and track->IsHybridITSTPC()
1075 // ******* TPC ********
1076 // 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
1077 // 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
1079 AliCodeTimerAuto("",0);
1081 // Loop over the tracks and extract and mask out all aod tracks that pass the selections for AODt racks
1082 for(int it = 0;it < fNumberOfTracks;++it)
1084 AliAODTrack *tr = (AliAODTrack*)(Tracks().At(it));
1086 UInt_t map = tr->GetFilterMap();
1087 if(map&fTPCConstrainedFilterMask){
1088 // we only reset the track select ionfo, no deletion...
1089 tr->SetFilterMap(map&~fTPCConstrainedFilterMask);
1091 if(map&fHybridFilterMaskTPCCG){
1092 // this is one part of the hybrid tracks
1093 // the others not passing the selection will be TPC only selected below
1094 tr->SetIsHybridTPCConstrainedGlobal(kTRUE);
1097 // Loop over the ESD trcks and pick out the tracks passing TPC only cuts
1100 const AliESDVertex *vtxSPD = esd.GetPrimaryVertexSPD();
1101 const AliESDVertex *vtx = esd.GetPrimaryVertex();
1103 Double_t pos[3] = { 0. };
1104 Double_t covTr[21]={0.};
1105 Double_t pid[10]={0.};
1108 Double_t p[3] = { 0. };
1110 Double_t pDCA[3] = { 0. }; // momentum at DCA
1111 Double_t rDCA[3] = { 0. }; // position at DCA
1112 Float_t dDCA[2] = {0.}; // DCA to the vertex d and z
1113 Float_t cDCA[3] = {0.}; // covariance of impact parameters
1116 AliAODTrack* aodTrack(0x0);
1117 // AliAODPid* detpid(0x0);
1119 // account for change in pT after the constraint
1120 Float_t ptMax = 1E10;
1122 for(int i = 0;i<32;i++){
1123 if(fTPCConstrainedFilterMask&(1<<i)){
1124 AliESDtrackCuts*cuts = (AliESDtrackCuts*)fTrackFilter->GetCuts()->At(i);
1125 Float_t tmp1= 0,tmp2 = 0;
1126 cuts->GetPtRange(tmp1,tmp2);
1127 if(tmp1>ptMin)ptMin=tmp1;
1128 if(tmp2<ptMax)ptMax=tmp2;
1132 for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack)
1134 AliESDtrack* esdTrack = esd.GetTrack(nTrack); //carefull do not modify it othwise need to work with a copy
1136 UInt_t selectInfo = 0;
1137 Bool_t isHybridITSTPC = false;
1141 selectInfo = fTrackFilter->IsSelected(esdTrack);
1144 if(!(selectInfo&fHybridFilterMaskTPCCG)){
1145 // not already selected tracks, use second part of hybrid tracks
1146 isHybridITSTPC = true;
1147 // too save space one could only store these...
1150 selectInfo &= fTPCConstrainedFilterMask;
1151 if (!selectInfo)continue;
1152 if (fWriteHybridTPCCOnly&&!isHybridITSTPC)continue; // write only complementary tracks
1153 // create a tpc only tracl
1154 AliESDtrack *track = AliESDtrackCuts::GetTPCOnlyTrack(const_cast<AliESDEvent*>(&esd),esdTrack->GetID());
1155 if(!track) continue;
1159 // only constrain tracks above threshold
1160 AliExternalTrackParam exParam;
1161 // take the B-field from the ESD, no 3D fieldMap available at this point
1162 Bool_t relate = false;
1163 relate = track->RelateToVertexTPC(vtxSPD,esd.GetMagneticField(),kVeryBig,&exParam);
1168 // fetch the track parameters at the DCA (unconstraint)
1169 if(track->GetTPCInnerParam()){
1170 track->GetTPCInnerParam()->GetPxPyPz(pDCA);
1171 track->GetTPCInnerParam()->GetXYZ(rDCA);
1173 // get the DCA to the vertex:
1174 track->GetImpactParametersTPC(dDCA,cDCA);
1175 // set the constrained parameters to the track
1176 track->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance());
1179 track->GetPxPyPz(p);
1181 Float_t pT = track->Pt();
1182 if(pT<ptMin||pT>ptMax){
1191 track->GetCovarianceXYZPxPyPz(covTr);
1192 esdTrack->GetESDpid(pid);// original PID
1194 if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());
1195 aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack((track->GetID()+1)*-1,
1202 (Short_t)track->GetSign(),
1203 track->GetITSClusterMap(),
1206 kTRUE, // check if this is right
1207 vtx->UsesTrack(track->GetID()),
1208 AliAODTrack::kPrimary,
1210 aodTrack->SetIsHybridTPCConstrainedGlobal(isHybridITSTPC);
1211 aodTrack->SetTPCFitMap(track->GetTPCFitMap());
1212 aodTrack->SetTPCClusterMap(track->GetTPCClusterMap());
1213 aodTrack->SetTPCSharedMap (track->GetTPCSharedMap());
1214 aodTrack->SetIsTPCConstrained(kTRUE);
1215 aodTrack->SetChi2perNDF(Chi2perNDF(esdTrack)); // original track
1216 // set the DCA values to the AOD track
1217 aodTrack->SetPxPyPzAtDCA(pDCA[0],pDCA[1],pDCA[2]);
1218 aodTrack->SetXYAtDCA(rDCA[0],rDCA[1]);
1219 aodTrack->SetDCA(dDCA[0],dDCA[1]);
1221 aodTrack->SetFlags(track->GetStatus());
1222 aodTrack->SetTPCPointsF(track->GetTPCNclsF());
1223 aodTrack->SetTPCNCrossedRows(UShort_t(track->GetTPCCrossedRows()));
1225 //Perform progagation of tracks if needed
1226 if(fDoPropagateTrackToEMCal) PropagateTrackToEMCal(esdTrack);
1227 aodTrack->SetTrackPhiEtaOnEMCal(esdTrack->GetTrackPhiOnEMCal(),esdTrack->GetTrackEtaOnEMCal());
1229 // do not duplicate PID information
1230 // aodTrack->ConvertAliPIDtoAODPID();
1231 // SetAODPID(esdTrack,aodTrack,detpid);
1234 } // end of loop on tracks
1238 //______________________________________________________________________________
1239 void AliAnalysisTaskESDfilter::ConvertGlobalConstrainedTracks(const AliESDEvent& esd)
1242 // Here we have the option to store the complement from global constraint information
1243 // to tracks passing tight cuts (1) in order not to get fakes back in, one needs
1244 // two sets of cuts one tight (1) (to throw out fakes) and one lose (2) (fakes/bad tracks would pass (2) but not (1))
1245 // using cut number (3) selects the tracks that complement (1) e.g. tracks witout ITS refit or cluster requirement
1248 AliCodeTimerAuto("",0);
1250 // Loop over the tracks and extract and mask out all aod tracks that pass the selections for AODt racks
1251 for(int it = 0;it < fNumberOfTracks;++it)
1253 AliAODTrack *tr = (AliAODTrack*)(Tracks().At(it));
1255 UInt_t map = tr->GetFilterMap();
1256 if(map&fGlobalConstrainedFilterMask){
1257 // we only reset the track select info, no deletion...
1258 // mask reset mask in case track is already taken
1259 tr->SetFilterMap(map&~fGlobalConstrainedFilterMask);
1261 if(map&fHybridFilterMaskGCG){
1262 // this is one part of the hybrid tracks
1263 // the others not passing the selection will be the ones selected below
1264 tr->SetIsHybridGlobalConstrainedGlobal(kTRUE);
1267 // Loop over the ESD trcks and pick out the tracks passing the GlobalConstraint cuts
1270 Double_t pos[3] = { 0. };
1271 Double_t covTr[21]={0.};
1272 Double_t pid[10]={0.};
1273 Double_t p[3] = { 0. };
1275 Double_t pDCA[3] = { 0. }; // momentum at DCA
1276 Double_t rDCA[3] = { 0. }; // position at DCA
1277 Float_t dDCA[2] = {0.}; // DCA to the vertex d and z
1278 Float_t cDCA[3] = {0.}; // covariance of impact parameters
1281 AliAODTrack* aodTrack(0x0);
1282 AliAODPid* detpid(0x0);
1283 const AliESDVertex *vtx = esd.GetPrimaryVertex();
1285 // account for change in pT after the constraint
1286 Float_t ptMax = 1E10;
1288 for(int i = 0;i<32;i++){
1289 if(fGlobalConstrainedFilterMask&(1<<i)){
1290 AliESDtrackCuts*cuts = (AliESDtrackCuts*)fTrackFilter->GetCuts()->At(i);
1291 Float_t tmp1= 0,tmp2 = 0;
1292 cuts->GetPtRange(tmp1,tmp2);
1293 if(tmp1>ptMin)ptMin=tmp1;
1294 if(tmp2<ptMax)ptMax=tmp2;
1300 for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack)
1302 AliESDtrack* esdTrack = esd.GetTrack(nTrack); //carefull do not modify it othwise need to work with a copy
1303 const AliExternalTrackParam * exParamGC = esdTrack->GetConstrainedParam();
1304 if(!exParamGC)continue;
1306 UInt_t selectInfo = 0;
1307 Bool_t isHybridGC = false;
1312 selectInfo = fTrackFilter->IsSelected(esdTrack);
1316 if(!(selectInfo&fHybridFilterMaskGCG))isHybridGC = true;
1317 if (fWriteHybridGCOnly&&!isHybridGC)continue; // write only complementary tracks
1319 selectInfo &= fGlobalConstrainedFilterMask;
1320 if (!selectInfo)continue;
1321 // fetch the track parameters at the DCA (unconstrained)
1322 esdTrack->GetPxPyPz(pDCA);
1323 esdTrack->GetXYZ(rDCA);
1324 // get the DCA to the vertex:
1325 esdTrack->GetImpactParameters(dDCA,cDCA);
1327 if (!esdTrack->GetConstrainedPxPyPz(p)) continue;
1330 Float_t pT = exParamGC->Pt();
1331 if(pT<ptMin||pT>ptMax){
1336 esdTrack->GetConstrainedXYZ(pos);
1337 exParamGC->GetCovarianceXYZPxPyPz(covTr);
1338 esdTrack->GetESDpid(pid);
1339 if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());
1340 aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack((esdTrack->GetID()+1)*-1,
1341 esdTrack->GetLabel(),
1347 (Short_t)esdTrack->GetSign(),
1348 esdTrack->GetITSClusterMap(),
1351 kTRUE, // check if this is right
1352 vtx->UsesTrack(esdTrack->GetID()),
1353 AliAODTrack::kPrimary,
1355 aodTrack->SetIsHybridGlobalConstrainedGlobal(isHybridGC);
1356 aodTrack->SetIsGlobalConstrained(kTRUE);
1357 aodTrack->SetTPCFitMap(esdTrack->GetTPCFitMap());
1358 aodTrack->SetTPCClusterMap(esdTrack->GetTPCClusterMap());
1359 aodTrack->SetTPCSharedMap (esdTrack->GetTPCSharedMap());
1360 aodTrack->SetChi2perNDF(Chi2perNDF(esdTrack));
1363 // set the DCA values to the AOD track
1364 aodTrack->SetPxPyPzAtDCA(pDCA[0],pDCA[1],pDCA[2]);
1365 aodTrack->SetXYAtDCA(rDCA[0],rDCA[1]);
1366 aodTrack->SetDCA(dDCA[0],dDCA[1]);
1368 aodTrack->SetFlags(esdTrack->GetStatus());
1369 aodTrack->SetTPCPointsF(esdTrack->GetTPCNclsF());
1370 aodTrack->SetTPCNCrossedRows(UShort_t(esdTrack->GetTPCCrossedRows()));
1373 // only copy AOD information for hybrid, no duplicate information
1374 aodTrack->ConvertAliPIDtoAODPID();
1375 SetAODPID(esdTrack,aodTrack,detpid);
1378 //Perform progagation of tracks if needed
1379 if(fDoPropagateTrackToEMCal) PropagateTrackToEMCal(esdTrack);
1380 aodTrack->SetTrackPhiEtaOnEMCal(esdTrack->GetTrackPhiOnEMCal(),esdTrack->GetTrackEtaOnEMCal());
1381 } // end of loop on tracks
1386 //______________________________________________________________________________
1387 void AliAnalysisTaskESDfilter::ConvertTracks(const AliESDEvent& esd)
1389 // Tracks (primary and orphan)
1391 AliCodeTimerAuto("",0);
1393 AliDebug(1,Form("NUMBER OF ESD TRACKS %5d\n", esd.GetNumberOfTracks()));
1395 const AliESDVertex *vtx = esd.GetPrimaryVertex();
1396 Double_t p[3] = { 0. };
1397 Double_t pos[3] = { 0. };
1398 Double_t covTr[21] = { 0. };
1399 Double_t pid[10] = { 0. };
1400 AliAODTrack* aodTrack(0x0);
1401 AliAODPid* detpid(0x0);
1403 for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack)
1405 if (fUsedTrack[nTrack]) continue;
1407 AliESDtrack *esdTrack = esd.GetTrack(nTrack);
1408 UInt_t selectInfo = 0;
1412 selectInfo = fTrackFilter->IsSelected(esdTrack);
1413 if (!selectInfo && !vtx->UsesTrack(esdTrack->GetID())) continue;
1417 esdTrack->GetPxPyPz(p);
1418 esdTrack->GetXYZ(pos);
1419 esdTrack->GetCovarianceXYZPxPyPz(covTr);
1420 esdTrack->GetESDpid(pid);
1421 if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());
1422 fPrimaryVertex->AddDaughter(aodTrack =
1423 new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrack->GetID(),
1424 esdTrack->GetLabel(),
1430 (Short_t)esdTrack->GetSign(),
1431 esdTrack->GetITSClusterMap(),
1434 kTRUE, // check if this is right
1435 vtx->UsesTrack(esdTrack->GetID()),
1436 AliAODTrack::kPrimary,
1439 aodTrack->SetTPCFitMap(esdTrack->GetTPCFitMap());
1440 aodTrack->SetTPCClusterMap(esdTrack->GetTPCClusterMap());
1441 aodTrack->SetTPCSharedMap (esdTrack->GetTPCSharedMap());
1442 aodTrack->SetChi2perNDF(Chi2perNDF(esdTrack));
1443 aodTrack->SetTPCPointsF(esdTrack->GetTPCNclsF());
1444 aodTrack->SetTPCNCrossedRows(UShort_t(esdTrack->GetTPCCrossedRows()));
1445 if(esdTrack->IsEMCAL()) aodTrack->SetEMCALcluster(esdTrack->GetEMCALcluster());
1446 if(esdTrack->IsPHOS()) aodTrack->SetPHOScluster(esdTrack->GetPHOScluster());
1448 //Perform progagation of tracks if needed
1449 if(fDoPropagateTrackToEMCal) PropagateTrackToEMCal(esdTrack);
1450 aodTrack->SetTrackPhiEtaOnEMCal(esdTrack->GetTrackPhiOnEMCal(),esdTrack->GetTrackEtaOnEMCal());
1452 fAODTrackRefs->AddAt(aodTrack, nTrack);
1455 if (esdTrack->GetSign() > 0) ++fNumberOfPositiveTracks;
1456 aodTrack->SetFlags(esdTrack->GetStatus());
1457 aodTrack->ConvertAliPIDtoAODPID();
1458 SetAODPID(esdTrack,aodTrack,detpid);
1459 } // end of loop on tracks
1462 //______________________________________________________________________________
1463 void AliAnalysisTaskESDfilter::PropagateTrackToEMCal(AliESDtrack *esdTrack)
1465 Double_t trkPos[3] = {0.,0.,0.};
1466 Double_t EMCalEta=-999, EMCalPhi=-999;
1467 Double_t trkphi = esdTrack->Phi()*TMath::RadToDeg();
1468 if(TMath::Abs(esdTrack->Eta())<0.9 && trkphi > 10 && trkphi < 250 )
1470 AliExternalTrackParam *trkParam = const_cast<AliExternalTrackParam*>(esdTrack->GetInnerParam());
1473 AliExternalTrackParam trkParamTmp(*trkParam);
1474 if(AliTrackerBase::PropagateTrackToBxByBz(&trkParamTmp, 430, esdTrack->GetMass(), 20, kTRUE, 0.8, -1))
1476 trkParamTmp.GetXYZ(trkPos);
1477 TVector3 trkPosVec(trkPos[0],trkPos[1],trkPos[2]);
1478 EMCalEta = trkPosVec.Eta();
1479 EMCalPhi = trkPosVec.Phi();
1480 if(EMCalPhi<0) EMCalPhi += 2*TMath::Pi();
1481 esdTrack->SetTrackPhiEtaOnEMCal(EMCalPhi,EMCalEta);
1487 //______________________________________________________________________________
1488 void AliAnalysisTaskESDfilter::ConvertPmdClusters(const AliESDEvent& esd)
1490 // Convert PMD Clusters
1491 AliCodeTimerAuto("",0);
1492 Int_t jPmdClusters=0;
1493 // Access to the AOD container of PMD clusters
1494 TClonesArray &pmdClusters = *(AODEvent()->GetPmdClusters());
1495 for (Int_t iPmd = 0; iPmd < esd.GetNumberOfPmdTracks(); ++iPmd) {
1496 // file pmd clusters, to be revised!
1497 AliESDPmdTrack *pmdTrack = esd.GetPmdTrack(iPmd);
1500 Double_t posPmd[3] = { pmdTrack->GetClusterX(), pmdTrack->GetClusterY(), pmdTrack->GetClusterZ()};
1501 Double_t pidPmd[13] = { 0.}; // to be revised!
1503 // assoc cluster not set
1504 new(pmdClusters[jPmdClusters++]) AliAODPmdCluster(iPmd, nLabel, label, pmdTrack->GetClusterADC(), posPmd, pidPmd);
1509 //______________________________________________________________________________
1510 void AliAnalysisTaskESDfilter::ConvertCaloClusters(const AliESDEvent& esd)
1512 // Convert Calorimeter Clusters
1513 AliCodeTimerAuto("",0);
1515 // Access to the AOD container of clusters
1516 TClonesArray &caloClusters = *(AODEvent()->GetCaloClusters());
1519 for (Int_t iClust=0; iClust<esd.GetNumberOfCaloClusters(); ++iClust) {
1521 AliESDCaloCluster * cluster = esd.GetCaloCluster(iClust);
1523 Int_t id = cluster->GetID();
1524 Int_t nLabel = cluster->GetNLabels();
1525 Int_t *labels = cluster->GetLabels();
1527 for(int i = 0;i < nLabel;++i){
1528 if(fMChandler)fMChandler->SelectParticle(labels[i]);
1532 Float_t energy = cluster->E();
1533 Float_t posF[3] = { 0.};
1534 cluster->GetPosition(posF);
1536 AliAODCaloCluster *caloCluster = new(caloClusters[jClusters++]) AliAODCaloCluster(id,
1542 cluster->GetType(),0);
1544 caloCluster->SetCaloCluster(cluster->GetDistanceToBadChannel(),
1545 cluster->GetDispersion(),
1546 cluster->GetM20(), cluster->GetM02(),
1547 cluster->GetEmcCpvDistance(),
1548 cluster->GetNExMax(),cluster->GetTOF()) ;
1550 caloCluster->SetPIDFromESD(cluster->GetPID());
1551 caloCluster->SetNCells(cluster->GetNCells());
1552 caloCluster->SetCellsAbsId(cluster->GetCellsAbsId());
1553 caloCluster->SetCellsAmplitudeFraction(cluster->GetCellsAmplitudeFraction());
1555 caloCluster->SetTrackDistance(cluster->GetTrackDx(), cluster->GetTrackDz());
1557 Int_t nMatchCount = 0;
1558 TArrayI* matchedT = cluster->GetTracksMatched();
1559 if (fNumberOfTracks>0 && matchedT && cluster->GetTrackMatchedIndex() >= 0) {
1560 for (Int_t im = 0; im < matchedT->GetSize(); im++) {
1561 Int_t iESDtrack = matchedT->At(im);;
1562 if (fAODTrackRefs->At(iESDtrack) != 0) {
1563 caloCluster->AddTrackMatched((AliAODTrack*)fAODTrackRefs->At(iESDtrack));
1569 caloCluster->SetTrackDistance(-999,-999);
1572 caloClusters.Expand(jClusters); // resize TObjArray to 'remove' slots for pseudo clusters
1575 //______________________________________________________________________________
1576 void AliAnalysisTaskESDfilter::ConvertCaloTrigger(TString calo, const AliESDEvent& esd)
1578 AliCodeTimerAuto("",0);
1582 AliAODCaloTrigger &aodTrigger = *(AODEvent()->GetCaloTrigger(calo));
1583 AliESDCaloTrigger &esdTrigger = *(esd.GetCaloTrigger(calo));
1585 aodTrigger.Allocate(esdTrigger.GetEntries());
1591 while (esdTrigger.Next()) {
1592 esdTrigger.GetPosition(tmod,tabsId);
1593 esdTrigger.GetAmplitude(a);
1594 aodTrigger.Add(tmod,tabsId,a,0.,(Int_t*)NULL,0,0,0);
1600 AliAODHandler *aodHandler = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
1604 TTree *aodTree = aodHandler->GetTree();
1608 Int_t *type = esd.GetCaloTriggerType();
1610 for (Int_t i = 0; i < 15; i++)
1612 aodTree->GetUserInfo()->Add(new TParameter<int>(Form("EMCALCaloTrigger%d",i), type[i]));
1617 AliAODCaloTrigger &aodTrigger = *(AODEvent()->GetCaloTrigger(calo));
1619 AliESDCaloTrigger &esdTrigger = *(esd.GetCaloTrigger(calo));
1621 aodTrigger.Allocate(esdTrigger.GetEntries());
1624 while (esdTrigger.Next())
1626 Int_t px, py, ts, nTimes, times[10], b;
1629 esdTrigger.GetPosition(px, py);
1631 esdTrigger.GetAmplitude(a);
1632 esdTrigger.GetTime(t);
1634 esdTrigger.GetL0Times(times);
1635 esdTrigger.GetNL0Times(nTimes);
1637 esdTrigger.GetL1TimeSum(ts);
1639 esdTrigger.GetTriggerBits(b);
1641 aodTrigger.Add(px, py, a, t, times, nTimes, ts, b);
1644 for (int i = 0; i < 4; i++) aodTrigger.SetL1Threshold(i, esdTrigger.GetL1Threshold(i));
1648 esdTrigger.GetL1V0(0),
1649 esdTrigger.GetL1V0(1)
1652 aodTrigger.SetL1V0(v0);
1653 aodTrigger.SetL1FrameMask(esdTrigger.GetL1FrameMask());
1656 //______________________________________________________________________________
1657 void AliAnalysisTaskESDfilter::ConvertEMCALCells(const AliESDEvent& esd)
1659 // Convert EMCAL Cells
1660 AliCodeTimerAuto("",0);
1661 // fill EMCAL cell info
1662 if (esd.GetEMCALCells()) { // protection against missing ESD information
1663 AliESDCaloCells &esdEMcells = *(esd.GetEMCALCells());
1664 Int_t nEMcell = esdEMcells.GetNumberOfCells() ;
1666 AliAODCaloCells &aodEMcells = *(AODEvent()->GetEMCALCells());
1667 aodEMcells.CreateContainer(nEMcell);
1668 aodEMcells.SetType(AliAODCaloCells::kEMCALCell);
1669 for (Int_t iCell = 0; iCell < nEMcell; iCell++) {
1670 aodEMcells.SetCell(iCell,esdEMcells.GetCellNumber(iCell),esdEMcells.GetAmplitude(iCell),
1671 esdEMcells.GetTime(iCell), esdEMcells.GetMCLabel(iCell), esdEMcells.GetEFraction(iCell));
1677 //______________________________________________________________________________
1678 void AliAnalysisTaskESDfilter::ConvertPHOSCells(const AliESDEvent& esd)
1680 // Convert PHOS Cells
1681 AliCodeTimerAuto("",0);
1682 // fill PHOS cell info
1683 if (esd.GetPHOSCells()) { // protection against missing ESD information
1684 AliESDCaloCells &esdPHcells = *(esd.GetPHOSCells());
1685 Int_t nPHcell = esdPHcells.GetNumberOfCells() ;
1687 AliAODCaloCells &aodPHcells = *(AODEvent()->GetPHOSCells());
1688 aodPHcells.CreateContainer(nPHcell);
1689 aodPHcells.SetType(AliAODCaloCells::kPHOSCell);
1690 for (Int_t iCell = 0; iCell < nPHcell; iCell++) {
1691 aodPHcells.SetCell(iCell,esdPHcells.GetCellNumber(iCell),esdPHcells.GetAmplitude(iCell),
1692 esdPHcells.GetTime(iCell), esdPHcells.GetMCLabel(iCell), esdPHcells.GetEFraction(iCell));
1698 //______________________________________________________________________________
1699 void AliAnalysisTaskESDfilter::ConvertTracklets(const AliESDEvent& esd)
1702 AliCodeTimerAuto("",0);
1704 AliAODTracklets &SPDTracklets = *(AODEvent()->GetTracklets());
1705 const AliMultiplicity *mult = esd.GetMultiplicity();
1707 if (mult->GetNumberOfTracklets()>0) {
1708 SPDTracklets.CreateContainer(mult->GetNumberOfTracklets());
1710 for (Int_t n=0; n<mult->GetNumberOfTracklets(); n++) {
1712 fMChandler->SelectParticle(mult->GetLabel(n, 0));
1713 fMChandler->SelectParticle(mult->GetLabel(n, 1));
1715 SPDTracklets.SetTracklet(n, mult->GetTheta(n), mult->GetPhi(n), mult->GetDeltaPhi(n), mult->GetLabel(n, 0),mult->GetLabel(n, 1));
1719 //Printf("ERROR: AliMultiplicity could not be retrieved from ESD");
1723 //______________________________________________________________________________
1724 void AliAnalysisTaskESDfilter::ConvertKinks(const AliESDEvent& esd)
1726 AliCodeTimerAuto("",0);
1728 // Kinks: it is a big mess the access to the information in the kinks
1729 // The loop is on the tracks in order to find the mother and daugther of each kink
1731 Double_t covTr[21]={0.};
1732 Double_t pid[10]={0.};
1733 AliAODPid* detpid(0x0);
1735 fNumberOfKinks = esd.GetNumberOfKinks();
1737 const AliESDVertex* vtx = esd.GetPrimaryVertex();
1739 for (Int_t iTrack=0; iTrack<esd.GetNumberOfTracks(); ++iTrack)
1741 AliESDtrack * esdTrack = esd.GetTrack(iTrack);
1743 Int_t ikink = esdTrack->GetKinkIndex(0);
1745 if (ikink && fNumberOfKinks) {
1746 // Negative kink index: mother, positive: daughter
1748 // Search for the second track of the kink
1750 for (Int_t jTrack = iTrack+1; jTrack<esd.GetNumberOfTracks(); ++jTrack) {
1752 AliESDtrack * esdTrack1 = esd.GetTrack(jTrack);
1754 Int_t jkink = esdTrack1->GetKinkIndex(0);
1756 if ( TMath::Abs(ikink)==TMath::Abs(jkink) ) {
1758 // The two tracks are from the same kink
1760 if (fUsedKink[TMath::Abs(ikink)-1]) continue; // skip used kinks
1763 Int_t idaughter = -1;
1765 if (ikink<0 && jkink>0) {
1770 else if (ikink>0 && jkink<0) {
1776 // cerr << "Error: Wrong combination of kink indexes: "
1777 // << ikink << " " << jkink << endl;
1781 // Add the mother track if it passed primary track selection cuts
1783 AliAODTrack * mother = NULL;
1785 UInt_t selectInfo = 0;
1787 selectInfo = fTrackFilter->IsSelected(esd.GetTrack(imother));
1788 if (!selectInfo) continue;
1791 if (!fUsedTrack[imother]) {
1793 fUsedTrack[imother] = kTRUE;
1795 AliESDtrack *esdTrackM = esd.GetTrack(imother);
1796 Double_t p[3] = { 0. };
1797 Double_t pos[3] = { 0. };
1798 esdTrackM->GetPxPyPz(p);
1799 esdTrackM->GetXYZ(pos);
1800 esdTrackM->GetCovarianceXYZPxPyPz(covTr);
1801 esdTrackM->GetESDpid(pid);
1802 if(fMChandler)fMChandler->SelectParticle(esdTrackM->GetLabel());
1804 new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrackM->GetID(),
1805 esdTrackM->GetLabel(),
1811 (Short_t)esdTrackM->GetSign(),
1812 esdTrackM->GetITSClusterMap(),
1815 kTRUE, // check if this is right
1816 vtx->UsesTrack(esdTrack->GetID()),
1817 AliAODTrack::kPrimary,
1819 mother->SetTPCFitMap(esdTrackM->GetTPCFitMap());
1820 mother->SetTPCClusterMap(esdTrackM->GetTPCClusterMap());
1821 mother->SetTPCSharedMap (esdTrackM->GetTPCSharedMap());
1822 mother->SetChi2perNDF(Chi2perNDF(esdTrackM));
1823 mother->SetTPCPointsF(esdTrackM->GetTPCNclsF());
1824 mother->SetTPCNCrossedRows(UShort_t(esdTrackM->GetTPCCrossedRows()));
1826 fAODTrackRefs->AddAt(mother, imother);
1828 if (esdTrackM->GetSign() > 0) ++fNumberOfPositiveTracks;
1829 mother->SetFlags(esdTrackM->GetStatus());
1830 mother->ConvertAliPIDtoAODPID();
1831 fPrimaryVertex->AddDaughter(mother);
1832 mother->ConvertAliPIDtoAODPID();
1833 SetAODPID(esdTrackM,mother,detpid);
1836 // cerr << "Error: event " << esd.GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1
1837 // << " track " << imother << " has already been used!" << endl;
1840 // Add the kink vertex
1841 AliESDkink * kink = esd.GetKink(TMath::Abs(ikink)-1);
1843 AliAODVertex * vkink =
1844 new(Vertices()[fNumberOfVertices++]) AliAODVertex(kink->GetPosition(),
1848 esdTrack->GetID(), // This is the track ID of the mother's track!
1849 AliAODVertex::kKink);
1850 // Add the daughter track
1852 AliAODTrack * daughter = NULL;
1854 if (!fUsedTrack[idaughter]) {
1856 fUsedTrack[idaughter] = kTRUE;
1858 AliESDtrack *esdTrackD = esd.GetTrack(idaughter);
1859 Double_t p[3] = { 0. };
1860 Double_t pos[3] = { 0. };
1862 esdTrackD->GetPxPyPz(p);
1863 esdTrackD->GetXYZ(pos);
1864 esdTrackD->GetCovarianceXYZPxPyPz(covTr);
1865 esdTrackD->GetESDpid(pid);
1867 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdTrackD);
1868 if(fMChandler)fMChandler->SelectParticle(esdTrackD->GetLabel());
1870 new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrackD->GetID(),
1871 esdTrackD->GetLabel(),
1877 (Short_t)esdTrackD->GetSign(),
1878 esdTrackD->GetITSClusterMap(),
1881 kTRUE, // check if this is right
1882 vtx->UsesTrack(esdTrack->GetID()),
1883 AliAODTrack::kSecondary,
1885 daughter->SetTPCFitMap(esdTrackD->GetTPCFitMap());
1886 daughter->SetTPCClusterMap(esdTrackD->GetTPCClusterMap());
1887 daughter->SetTPCSharedMap (esdTrackD->GetTPCSharedMap());
1888 daughter->SetTPCPointsF(esdTrackD->GetTPCNclsF());
1889 daughter->SetTPCNCrossedRows(UShort_t(esdTrackD->GetTPCCrossedRows()));
1890 fAODTrackRefs->AddAt(daughter, idaughter);
1892 if (esdTrackD->GetSign() > 0) ++fNumberOfPositiveTracks;
1893 daughter->SetFlags(esdTrackD->GetStatus());
1894 daughter->ConvertAliPIDtoAODPID();
1895 vkink->AddDaughter(daughter);
1896 daughter->ConvertAliPIDtoAODPID();
1897 SetAODPID(esdTrackD,daughter,detpid);
1900 // cerr << "Error: event " << esd.GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1
1901 // << " track " << idaughter << " has already been used!" << endl;
1909 //______________________________________________________________________________
1910 void AliAnalysisTaskESDfilter::ConvertPrimaryVertices(const AliESDEvent& esd)
1912 AliCodeTimerAuto("",0);
1914 // Access to the AOD container of vertices
1915 fNumberOfVertices = 0;
1917 Double_t pos[3] = { 0. };
1918 Double_t covVtx[6] = { 0. };
1920 // Add primary vertex. The primary tracks will be defined
1921 // after the loops on the composite objects (V0, cascades, kinks)
1922 const AliESDVertex *vtx = esd.GetPrimaryVertex();
1924 vtx->GetXYZ(pos); // position
1925 vtx->GetCovMatrix(covVtx); //covariance matrix
1927 fPrimaryVertex = new(Vertices()[fNumberOfVertices++])
1928 AliAODVertex(pos, covVtx, vtx->GetChi2toNDF(), NULL, -1, AliAODVertex::kPrimary);
1929 fPrimaryVertex->SetName(vtx->GetName());
1930 fPrimaryVertex->SetTitle(vtx->GetTitle());
1931 fPrimaryVertex->SetBC(vtx->GetBC());
1933 TString vtitle = vtx->GetTitle();
1934 if (!vtitle.Contains("VertexerTracks"))
1935 fPrimaryVertex->SetNContributors(vtx->GetNContributors());
1937 if (fDebug > 0) fPrimaryVertex->Print();
1939 // Add SPD "main" vertex
1940 const AliESDVertex *vtxS = esd.GetPrimaryVertexSPD();
1941 vtxS->GetXYZ(pos); // position
1942 vtxS->GetCovMatrix(covVtx); //covariance matrix
1943 AliAODVertex * mVSPD = new(Vertices()[fNumberOfVertices++])
1944 AliAODVertex(pos, covVtx, vtxS->GetChi2toNDF(), NULL, -1, AliAODVertex::kMainSPD);
1945 mVSPD->SetName(vtxS->GetName());
1946 mVSPD->SetTitle(vtxS->GetTitle());
1947 mVSPD->SetNContributors(vtxS->GetNContributors());
1949 // Add SPD pileup vertices
1950 for(Int_t iV=0; iV<esd.GetNumberOfPileupVerticesSPD(); ++iV)
1952 const AliESDVertex *vtxP = esd.GetPileupVertexSPD(iV);
1953 vtxP->GetXYZ(pos); // position
1954 vtxP->GetCovMatrix(covVtx); //covariance matrix
1955 AliAODVertex * pVSPD = new(Vertices()[fNumberOfVertices++])
1956 AliAODVertex(pos, covVtx, vtxP->GetChi2toNDF(), NULL, -1, AliAODVertex::kPileupSPD);
1957 pVSPD->SetName(vtxP->GetName());
1958 pVSPD->SetTitle(vtxP->GetTitle());
1959 pVSPD->SetNContributors(vtxP->GetNContributors());
1960 pVSPD->SetBC(vtxP->GetBC());
1963 // Add TRK pileup vertices
1964 for(Int_t iV=0; iV<esd.GetNumberOfPileupVerticesTracks(); ++iV)
1966 const AliESDVertex *vtxP = esd.GetPileupVertexTracks(iV);
1967 vtxP->GetXYZ(pos); // position
1968 vtxP->GetCovMatrix(covVtx); //covariance matrix
1969 AliAODVertex * pVTRK = new(Vertices()[fNumberOfVertices++])
1970 AliAODVertex(pos, covVtx, vtxP->GetChi2toNDF(), NULL, -1, AliAODVertex::kPileupTracks);
1971 pVTRK->SetName(vtxP->GetName());
1972 pVTRK->SetTitle(vtxP->GetTitle());
1973 pVTRK->SetNContributors(vtxP->GetNContributors());
1974 pVTRK->SetBC(vtxP->GetBC());
1978 //______________________________________________________________________________
1979 void AliAnalysisTaskESDfilter::ConvertVZERO(const AliESDEvent& esd)
1981 // Convert VZERO data
1982 AliAODVZERO* vzeroData = AODEvent()->GetVZEROData();
1983 *vzeroData = *(esd.GetVZEROData());
1986 //______________________________________________________________________________
1987 void AliAnalysisTaskESDfilter::ConvertTZERO(const AliESDEvent& esd)
1989 // Convert TZERO data
1990 const AliESDTZERO* esdTzero = esd.GetESDTZERO();
1991 AliAODTZERO* aodTzero = AODEvent()->GetTZEROData();
1993 for (Int_t icase=0; icase<3; icase++){
1994 aodTzero->SetT0TOF( icase, esdTzero->GetT0TOF(icase));
1995 aodTzero->SetT0TOFbest(icase, esdTzero->GetT0TOFbest(icase));
1997 aodTzero->SetBackgroundFlag(esdTzero->GetBackgroundFlag());
1998 aodTzero->SetPileupFlag(esdTzero->GetPileupFlag());
1999 aodTzero->SetSatelliteFlag(esdTzero->GetSatellite());
2001 Float_t rawTime[24];
2002 for(Int_t ipmt=0; ipmt<24; ipmt++)
2003 rawTime[ipmt] = esdTzero->GetTimeFull(ipmt,0);
2005 Int_t idxOfFirstPmtA = -1, idxOfFirstPmtC = -1;
2006 Float_t timeOfFirstPmtA = 9999, timeOfFirstPmtC = 9999;
2007 for(int ipmt=0; ipmt<12; ipmt++){
2008 if( rawTime[ipmt] > -200 && rawTime[ipmt] < timeOfFirstPmtC && rawTime[ipmt]!=0){
2009 timeOfFirstPmtC = rawTime[ipmt];
2010 idxOfFirstPmtC = ipmt;
2013 for(int ipmt=12; ipmt<24; ipmt++){
2014 if( rawTime[ipmt] > -200 && rawTime[ipmt] < timeOfFirstPmtA && rawTime[ipmt]!=0 ){
2015 timeOfFirstPmtA = rawTime[ipmt];
2016 idxOfFirstPmtA = ipmt;
2020 if(idxOfFirstPmtA != -1 && idxOfFirstPmtC != -1){
2021 //speed of light in cm/ns TMath::C()*1e-7
2022 Float_t vertexraw = TMath::C()*1e-7 * (rawTime[idxOfFirstPmtA] - rawTime[idxOfFirstPmtC])/2;
2023 aodTzero->SetT0VertexRaw( vertexraw );
2025 aodTzero->SetT0VertexRaw(99999);
2028 aodTzero->SetT0zVertex(esdTzero->GetT0zVertex());
2032 //______________________________________________________________________________
2033 void AliAnalysisTaskESDfilter::ConvertZDC(const AliESDEvent& esd)
2036 AliESDZDC* esdZDC = esd.GetZDCData();
2038 const Double_t zem1Energy = esdZDC->GetZEM1Energy();
2039 const Double_t zem2Energy = esdZDC->GetZEM2Energy();
2041 const Double_t *towZNC = esdZDC->GetZNCTowerEnergy();
2042 const Double_t *towZPC = esdZDC->GetZPCTowerEnergy();
2043 const Double_t *towZNA = esdZDC->GetZNATowerEnergy();
2044 const Double_t *towZPA = esdZDC->GetZPATowerEnergy();
2045 const Double_t *towZNCLG = esdZDC->GetZNCTowerEnergyLR();
2046 const Double_t *towZNALG = esdZDC->GetZNATowerEnergyLR();
2048 AliAODZDC* zdcAOD = AODEvent()->GetZDCData();
2050 zdcAOD->SetZEM1Energy(zem1Energy);
2051 zdcAOD->SetZEM2Energy(zem2Energy);
2052 zdcAOD->SetZNCTowers(towZNC, towZNCLG);
2053 zdcAOD->SetZNATowers(towZNA, towZNALG);
2054 zdcAOD->SetZPCTowers(towZPC);
2055 zdcAOD->SetZPATowers(towZPA);
2057 zdcAOD->SetZDCParticipants(esdZDC->GetZDCParticipants(), esdZDC->GetZDCPartSideA(), esdZDC->GetZDCPartSideC());
2058 zdcAOD->SetZDCImpactParameter(esdZDC->GetImpactParameter(), esdZDC->GetImpactParamSideA(),
2059 esdZDC->GetImpactParamSideC());
2060 zdcAOD->SetZDCTDCSum(esdZDC->GetZNTDCSum(0));
2061 zdcAOD->SetZDCTDCDiff(esdZDC->GetZNTDCDiff(0));
2065 //_______________________________________________________________________________________________________________________________________
2066 Int_t AliAnalysisTaskESDfilter::ConvertHMPID(const AliESDEvent& esd) // clm
2069 // Convtert ESD HMPID info to AOD and return the number of good tracks with HMPID signal.
2070 // We need to return an int since there is no signal counter in the ESD.
2073 AliCodeTimerAuto("",0);
2075 Int_t cntHmpidGoodTracks = 0;
2084 Float_t thetaTrk = 0;
2087 Double_t hmpPid[5]={0};
2088 Double_t hmpMom[3]={0};
2090 TClonesArray &hmpidRings = *(AODEvent()->GetHMPIDrings());
2092 for (Int_t iTrack=0; iTrack<esd.GetNumberOfTracks(); ++iTrack)
2094 if(! esd.GetTrack(iTrack) ) continue;
2096 if(esd.GetTrack(iTrack)->GetHMPIDsignal() > -20 ) { //
2098 (esd.GetTrack(iTrack))->GetHMPIDmip(xMip, yMip, qMip, nphMip); // Get MIP properties
2099 (esd.GetTrack(iTrack))->GetHMPIDtrk(xTrk,yTrk,thetaTrk,phiTrk);
2100 (esd.GetTrack(iTrack))->GetHMPIDpid(hmpPid);
2101 if((esd.GetTrack(iTrack))->GetOuterHmpParam()) (esd.GetTrack(iTrack))->GetOuterHmpPxPyPz(hmpMom);
2103 if(esd.GetTrack(iTrack)->GetHMPIDsignal() == 0 && thetaTrk == 0 && qMip == 0 && nphMip ==0 ) continue; //
2105 new(hmpidRings[cntHmpidGoodTracks++]) AliAODHMPIDrings(
2106 (esd.GetTrack(iTrack))->GetID(), // Unique track id to attach the ring to
2107 1000000*nphMip+qMip, // MIP charge and number of photons
2108 (esd.GetTrack(iTrack))->GetHMPIDcluIdx(), // 1000000*chamber id + cluster idx of the assigned MIP cluster
2109 thetaTrk, // track inclination angle theta
2110 phiTrk, // track inclination angle phi
2111 (esd.GetTrack(iTrack))->GetHMPIDsignal(), // Cherenkov angle
2112 (esd.GetTrack(iTrack))->GetHMPIDoccupancy(), // Occupancy claculated for the given chamber
2113 (esd.GetTrack(iTrack))->GetHMPIDchi2(), // Ring resolution squared
2114 xTrk, // Track x coordinate (LORS)
2115 yTrk, // Track y coordinate (LORS)
2116 xMip, // MIP x coordinate (LORS)
2117 yMip, // MIP y coordinate (LORS)
2118 hmpPid, // PID probablities from ESD, remove later once it is in CombinedPid
2119 hmpMom // Track momentum in HMPID at ring reconstruction
2122 // Printf(Form("+++++++++ yes/no: %d %lf %lf %lf %lf %lf %lf ",(esd.GetTrack(iTrack))->IsHMPID(),thetaTrk, (esd.GetTrack(iTrack))->GetHMPIDchi2(),xTrk, yTrk , xMip, yMip));
2125 }// HMPID signal > -20
2126 }//___esd track loop
2128 return cntHmpidGoodTracks;
2131 //______________________________________________________________________________
2132 void AliAnalysisTaskESDfilter::ConvertESDtoAOD()
2134 // ESD Filter analysis task executed for each event
2136 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
2140 AliCodeTimerAuto("",0);
2142 fOldESDformat = ( esd->GetAliESDOld() != 0x0 );
2144 // Reconstruct cascades and V0 here
2145 if (fIsV0CascadeRecoEnabled) {
2146 esd->ResetCascades();
2149 AliV0vertexer lV0vtxer;
2150 AliCascadeVertexer lCascVtxer;
2152 lV0vtxer.SetCuts(fV0Cuts);
2153 lCascVtxer.SetCuts(fCascadeCuts);
2156 lV0vtxer.Tracks2V0vertices(esd);
2157 lCascVtxer.V0sTracks2CascadeVertices(esd);
2161 fNumberOfTracks = 0;
2162 fNumberOfPositiveTracks = 0;
2164 fNumberOfVertices = 0;
2165 fNumberOfCascades = 0;
2168 AliAODHeader* header = ConvertHeader(*esd);
2170 if ( fIsVZEROEnabled ) ConvertVZERO(*esd);
2171 if ( fIsTZEROEnabled ) ConvertTZERO(*esd);
2173 // Fetch Stack for debuggging if available
2177 fMChandler = (AliMCEventHandler*) ((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler());
2180 // loop over events and fill them
2181 // Multiplicity information needed by the header (to be revised!)
2182 Int_t nTracks = esd->GetNumberOfTracks();
2183 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) esd->GetTrack(iTrack)->SetESDEvent(esd);
2185 // Update the header
2187 Int_t nV0s = esd->GetNumberOfV0s();
2188 Int_t nCascades = esd->GetNumberOfCascades();
2189 Int_t nKinks = esd->GetNumberOfKinks();
2190 Int_t nVertices = nV0s + nCascades /*V0 wihtin cascade already counted*/+ nKinks + 1 /* = prim. vtx*/;
2191 Int_t nPileSPDVertices=1+esd->GetNumberOfPileupVerticesSPD(); // also SPD main vertex
2192 Int_t nPileTrkVertices=esd->GetNumberOfPileupVerticesTracks();
2193 nVertices+=nPileSPDVertices;
2194 nVertices+=nPileTrkVertices;
2196 Int_t nCaloClus = esd->GetNumberOfCaloClusters();
2198 Int_t nPmdClus = esd->GetNumberOfPmdTracks();
2199 Int_t nHmpidRings = 0;
2201 AliDebug(1,Form(" NV0=%d NCASCADES=%d NKINKS=%d", nV0s, nCascades, nKinks));
2203 AODEvent()->ResetStd(nTracks, nVertices, nV0s, nCascades, nJets, nCaloClus, nFmdClus, nPmdClus,nHmpidRings);
2207 // RefArray to store a mapping between esd V0 number and newly created AOD-Vertex V0
2208 fAODV0VtxRefs = new TRefArray(nV0s);
2209 // RefArray to store the mapping between esd V0 number and newly created AOD-V0
2210 fAODV0Refs = new TRefArray(nV0s);
2211 // Array to take into account the V0s already added to the AOD (V0 within cascades)
2212 fUsedV0 = new Bool_t[nV0s];
2213 for (Int_t iV0=0; iV0<nV0s; ++iV0) fUsedV0[iV0]=kFALSE;
2218 // RefArray to store the mapping between esd track number and newly created AOD-Track
2220 fAODTrackRefs = new TRefArray(nTracks);
2222 // Array to take into account the tracks already added to the AOD
2223 fUsedTrack = new Bool_t[nTracks];
2224 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) fUsedTrack[iTrack]=kFALSE;
2227 // Array to take into account the kinks already added to the AOD
2230 fUsedKink = new Bool_t[nKinks];
2231 for (Int_t iKink=0; iKink<nKinks; ++iKink) fUsedKink[iKink]=kFALSE;
2234 ConvertPrimaryVertices(*esd);
2236 //setting best TOF PID
2237 AliESDInputHandler* esdH = dynamic_cast<AliESDInputHandler*>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
2239 fESDpid = esdH->GetESDpid();
2241 if (fIsPidOwner && fESDpid){
2246 { //in case of no Tender attached
2247 fESDpid = new AliESDpid;
2248 fIsPidOwner = kTRUE;
2251 if(!esd->GetTOFHeader())
2252 { //protection in case the pass2 LHC10b,c,d have been processed without tender.
2253 Float_t t0spread[10];
2254 Float_t intrinsicTOFres=100; //ps ok for LHC10b,c,d pass2!!
2255 for (Int_t i=0; i<10; i++) t0spread[i] = (TMath::Sqrt(esd->GetSigma2DiamondZ()))/0.03; //0.03 to convert from cm to ps
2256 fESDpid->GetTOFResponse().SetT0resolution(t0spread);
2257 fESDpid->GetTOFResponse().SetTimeResolution(intrinsicTOFres);
2258 // fESDpid->SetTOFResponse(esd, (AliESDpid::EStartTimeType_t)fTimeZeroType);
2259 AliTOFHeader tmpTOFHeader(0,t0spread[0],0,NULL,NULL,NULL,intrinsicTOFres,t0spread[0]);
2260 AODEvent()->SetTOFHeader(&tmpTOFHeader); // write dummy TOF header in AOD
2262 AODEvent()->SetTOFHeader(esd->GetTOFHeader()); // write TOF header in AOD
2265 // if(esd->GetTOFHeader() && fIsPidOwner) fESDpid->SetTOFResponse(esd, (AliESDpid::EStartTimeType_t)fTimeZeroType); //in case of AOD production strating form LHC10e without Tender.
2267 if ( fAreCascadesEnabled ) ConvertCascades(*esd);
2269 if ( fAreV0sEnabled ) ConvertV0s(*esd);
2271 if ( fAreKinksEnabled ) ConvertKinks(*esd);
2273 if ( fAreTracksEnabled ) ConvertTracks(*esd);
2275 // Update number of AOD tracks in header at the end of track loop (M.G.)
2276 header->SetRefMultiplicity(fNumberOfTracks);
2277 header->SetRefMultiplicityPos(fNumberOfPositiveTracks);
2278 header->SetRefMultiplicityNeg(fNumberOfTracks - fNumberOfPositiveTracks);
2280 if ( fTPCConstrainedFilterMask ) ConvertTPCOnlyTracks(*esd);
2281 if( fGlobalConstrainedFilterMask) ConvertGlobalConstrainedTracks(*esd);
2283 if ( fArePmdClustersEnabled ) ConvertPmdClusters(*esd);
2285 if ( fAreCaloClustersEnabled ) ConvertCaloClusters(*esd);
2287 if ( fAreEMCALCellsEnabled )ConvertEMCALCells(*esd);
2289 if ( fArePHOSCellsEnabled )ConvertPHOSCells(*esd);
2291 if ( fAreEMCALTriggerEnabled )ConvertCaloTrigger(TString("EMCAL"), *esd);
2293 if ( fArePHOSTriggerEnabled )ConvertCaloTrigger(TString("PHOS"), *esd);
2295 if ( fAreTrackletsEnabled ) ConvertTracklets(*esd);
2296 if ( fIsZDCEnabled ) ConvertZDC(*esd);
2298 if(fIsHMPIDEnabled) nHmpidRings = ConvertHMPID(*esd);
2300 delete fAODTrackRefs; fAODTrackRefs=0x0;
2301 delete fAODV0VtxRefs; fAODV0VtxRefs=0x0;
2302 delete fAODV0Refs; fAODV0Refs=0x0;
2304 delete[] fUsedTrack; fUsedTrack=0x0;
2305 delete[] fUsedV0; fUsedV0=0x0;
2306 delete[] fUsedKink; fUsedKink=0x0;
2317 //______________________________________________________________________________
2318 void AliAnalysisTaskESDfilter::SetAODPID(AliESDtrack *esdtrack, AliAODTrack *aodtrack, AliAODPid *detpid)
2321 // Setter for the raw PID detector signals
2324 // Save PID object for candidate electrons
2325 Bool_t pidSave = kFALSE;
2327 Bool_t selectInfo = fTrackFilter->IsSelected((char*) "Electrons");
2328 if (selectInfo) pidSave = kTRUE;
2332 // Tracks passing pt cut
2333 if(esdtrack->Pt()>fHighPthreshold) {
2337 if(esdtrack->Pt()> fPtshape->GetXmin()){
2338 Double_t y = fPtshape->Eval(esdtrack->Pt())/fPtshape->Eval(fHighPthreshold);
2339 if(gRandom->Rndm(0)<1./y){
2343 }//end if p function
2347 if(!aodtrack->GetDetPid()){// prevent memory leak when calling SetAODPID twice for the same track
2348 detpid = new AliAODPid();
2349 SetDetectorRawSignals(detpid,esdtrack);
2350 aodtrack->SetDetPID(detpid);
2355 //______________________________________________________________________________
2356 void AliAnalysisTaskESDfilter::SetDetectorRawSignals(AliAODPid *aodpid, AliESDtrack *track)
2359 //assignment of the detector signals (AliXXXesdPID inspired)
2362 AliInfo("no ESD track found. .....exiting");
2366 const AliExternalTrackParam *in=track->GetInnerParam();
2368 aodpid->SetTPCmomentum(in->GetP());
2370 aodpid->SetTPCmomentum(-1.);
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::kSPECIES]; track->GetIntegratedTimes(times);
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);
2420 // aodpid->SetHMPIDsignal(0); // set to zero for compression but it will be removed later
2424 Double_t AliAnalysisTaskESDfilter::Chi2perNDF(AliESDtrack* track)
2426 // Calculate chi2 per ndf for track
2427 Int_t nClustersTPC = track->GetTPCNcls();
2429 if ( nClustersTPC > 5) {
2430 return (track->GetTPCchi2()/Float_t(nClustersTPC - 5));
2437 //______________________________________________________________________________
2438 void AliAnalysisTaskESDfilter::Terminate(Option_t */*option*/)
2440 // Terminate analysis
2442 if (fDebug > 1) printf("AnalysisESDfilter: Terminate() \n");
2445 //______________________________________________________________________________
2446 void AliAnalysisTaskESDfilter::PrintMCInfo(AliStack *pStack,Int_t label){
2449 label = TMath::Abs(label);
2450 TParticle *part = pStack->Particle(label);
2451 Printf("########################");
2452 Printf("%s:%d %d UniqueID %d PDG %d P %3.3f",(char*)__FILE__,__LINE__,label,part->GetUniqueID(),part->GetPdgCode(),part->P());
2454 TParticle* mother = part;
2455 Int_t imo = part->GetFirstMother();
2456 Int_t nprim = pStack->GetNprimary();
2457 // while((imo >= nprim) && (mother->GetUniqueID() == 4)) {
2458 while((imo >= nprim)) {
2459 mother = pStack->Particle(imo);
2460 Printf("Mother %s:%d Label %d UniqueID %d PDG %d P %3.3f",(char*)__FILE__,__LINE__,imo,mother->GetUniqueID(),mother->GetPdgCode(),mother->P());
2462 imo = mother->GetFirstMother();
2464 Printf("########################");
2467 //______________________________________________________