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),
119 fEMCalSurfaceDistance(440)
121 // Default constructor
122 fV0Cuts[0] = 33. ; // max allowed chi2
123 fV0Cuts[1] = 0.1 ; // min allowed impact parameter for the 1st daughter
124 fV0Cuts[2] = 0.1 ; // min allowed impact parameter for the 2nd daughter
125 fV0Cuts[3] = 1. ; // max allowed DCA between the daughter tracks
126 fV0Cuts[4] = .998; // min allowed cosine of V0's pointing angle
127 fV0Cuts[5] = 0.9 ; // min radius of the fiducial volume
128 fV0Cuts[6] = 100. ; // max radius of the fiducial volume
130 fCascadeCuts[0] = 33. ; // max allowed chi2 (same as PDC07)
131 fCascadeCuts[1] = 0.05 ; // min allowed V0 impact parameter
132 fCascadeCuts[2] = 0.008; // "window" around the Lambda mass
133 fCascadeCuts[3] = 0.03 ; // min allowed bachelor's impact parameter
134 fCascadeCuts[4] = 0.3 ; // max allowed DCA between the V0 and the bachelor
135 fCascadeCuts[5] = 0.999; // min allowed cosine of the cascade pointing angle
136 fCascadeCuts[6] = 0.9 ; // min radius of the fiducial volume
137 fCascadeCuts[7] = 100. ; // max radius of the fiducial volume
140 //______________________________________________________________________________
141 AliAnalysisTaskESDfilter::AliAnalysisTaskESDfilter(const char* name):
142 AliAnalysisTaskSE(name),
149 fEnableFillAOD(kTRUE),
158 fNumberOfPositiveTracks(0),
160 fNumberOfVertices(0),
161 fNumberOfCascades(0),
163 fOldESDformat(kFALSE),
165 fTPCConstrainedFilterMask(0),
166 fHybridFilterMaskTPCCG(0),
167 fWriteHybridTPCCOnly(kFALSE),
168 fGlobalConstrainedFilterMask(0),
169 fHybridFilterMaskGCG(0),
170 fWriteHybridGCOnly(kFALSE),
171 fIsVZEROEnabled(kTRUE),
172 fIsTZEROEnabled(kTRUE),
173 fIsZDCEnabled(kTRUE),
174 fIsHMPIDEnabled(kTRUE),
175 fIsV0CascadeRecoEnabled(kFALSE),
176 fAreCascadesEnabled(kTRUE),
177 fAreV0sEnabled(kTRUE),
178 fAreKinksEnabled(kTRUE),
179 fAreTracksEnabled(kTRUE),
180 fArePmdClustersEnabled(kTRUE),
181 fAreCaloClustersEnabled(kTRUE),
182 fAreEMCALCellsEnabled(kTRUE),
183 fArePHOSCellsEnabled(kTRUE),
184 fAreEMCALTriggerEnabled(kTRUE),
185 fArePHOSTriggerEnabled(kTRUE),
186 fAreTrackletsEnabled(kTRUE),
189 fTPCaloneTrackCuts(0),
190 fDoPropagateTrackToEMCal(kTRUE),
191 fEMCalSurfaceDistance(440)
195 fV0Cuts[0] = 33. ; // max allowed chi2
196 fV0Cuts[1] = 0.1 ; // min allowed impact parameter for the 1st daughter
197 fV0Cuts[2] = 0.1 ; // min allowed impact parameter for the 2nd daughter
198 fV0Cuts[3] = 1. ; // max allowed DCA between the daughter tracks
199 fV0Cuts[4] = .998; // min allowed cosine of V0's pointing angle
200 fV0Cuts[5] = 0.9 ; // min radius of the fiducial volume
201 fV0Cuts[6] = 100. ; // max radius of the fiducial volume
203 fCascadeCuts[0] = 33. ; // max allowed chi2 (same as PDC07)
204 fCascadeCuts[1] = 0.05 ; // min allowed V0 impact parameter
205 fCascadeCuts[2] = 0.008; // "window" around the Lambda mass
206 fCascadeCuts[3] = 0.03 ; // min allowed bachelor's impact parameter
207 fCascadeCuts[4] = 0.3 ; // max allowed DCA between the V0 and the bachelor
208 fCascadeCuts[5] = 0.999; // min allowed cosine of the cascade pointing angle
209 fCascadeCuts[6] = 0.9 ; // min radius of the fiducial volume
210 fCascadeCuts[7] = 100. ; // max radius of the fiducial volume
215 AliAnalysisTaskESDfilter::~AliAnalysisTaskESDfilter(){
216 if(fIsPidOwner) delete fESDpid;
218 //______________________________________________________________________________
219 void AliAnalysisTaskESDfilter::UserCreateOutputObjects()
222 // Create Output Objects conenct filter to outputtree
226 OutputTree()->GetUserInfo()->Add(fTrackFilter);
230 AliError("No OutputTree() for adding the track filter");
232 fTPCaloneTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
235 //______________________________________________________________________________
236 void AliAnalysisTaskESDfilter::Init()
239 if (fDebug > 1) AliInfo("Init() \n");
240 // Call configuration file
243 //______________________________________________________________________________
244 Bool_t AliAnalysisTaskESDfilter::Notify()
247 AddMetadataToUserInfo();
251 //______________________________________________________________________________
252 Bool_t AliAnalysisTaskESDfilter::AddMetadataToUserInfo()
254 // Copy metadata to AOD user info.
255 static Bool_t copyFirst = kFALSE;
257 AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
259 AliError("AliAnalysisTaskESDfilter::AddMetadataToUserInfo() : No analysis manager !");
262 TTree *esdTree = mgr->GetTree()->GetTree();
263 if (!esdTree) return kFALSE;
264 TNamed *alirootVersion = (TNamed*)esdTree->GetUserInfo()->FindObject("alirootVersion");
265 if (!alirootVersion) return kFALSE;
266 AliAODHandler *aodHandler = dynamic_cast<AliAODHandler*>(mgr->GetOutputEventHandler());
267 if (!aodHandler) return kFALSE;
268 TTree *aodTree = aodHandler->GetTree();
269 if (!aodTree) return kFALSE;
270 aodTree->GetUserInfo()->Add(new TNamed(*alirootVersion));
276 //______________________________________________________________________________
277 void AliAnalysisTaskESDfilter::PrintTask(Option_t *option, Int_t indent) const
279 // Print selection task information
282 AliAnalysisTaskSE::PrintTask(option,indent);
284 TString spaces(' ',indent+3);
286 cout << spaces.Data() << Form("Cascades are %s",fAreCascadesEnabled ? "ENABLED":"DISABLED") << endl;
287 cout << spaces.Data() << Form("V0s are %s",fAreV0sEnabled ? "ENABLED":"DISABLED") << endl;
288 cout << spaces.Data() << Form("Kinks are %s",fAreKinksEnabled ? "ENABLED":"DISABLED") << endl;
289 cout << spaces.Data() << Form("Tracks are %s",fAreTracksEnabled ? "ENABLED":"DISABLED") << endl;
290 cout << spaces.Data() << Form("PmdClusters are %s",fArePmdClustersEnabled ? "ENABLED":"DISABLED") << endl;
291 cout << spaces.Data() << Form("CaloClusters are %s",fAreCaloClustersEnabled ? "ENABLED":"DISABLED") << endl;
292 cout << spaces.Data() << Form("EMCAL cells are %s",fAreEMCALCellsEnabled ? "ENABLED":"DISABLED") << endl;
293 cout << spaces.Data() << Form("EMCAL triggers are %s",fAreEMCALTriggerEnabled ? "ENABLED":"DISABLED") << endl;
294 cout << spaces.Data() << Form("PHOS triggers are %s",fArePHOSTriggerEnabled ? "ENABLED":"DISABLED") << endl;
295 cout << spaces.Data() << Form("Tracklets are %s",fAreTrackletsEnabled ? "ENABLED":"DISABLED") << endl;
296 cout << spaces.Data() << Form("PropagateTrackToEMCal is %s", fDoPropagateTrackToEMCal ? "ENABLED":"DISABLED") << endl;
299 //______________________________________________________________________________
300 void AliAnalysisTaskESDfilter::UserExec(Option_t */*option*/)
302 // Execute analysis for current event
305 Long64_t ientry = Entry();
308 printf("Filter: Analysing event # %5d\n", (Int_t) ientry);
309 if (fHighPthreshold == 0) AliInfo("detector PID signals are stored in each track");
310 if (!fPtshape) AliInfo("detector PID signals are not stored below the pt threshold");
312 // Filters must explicitely enable AOD filling in their UserExec (AG)
313 if (!AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()) AliFatal("Cannot run ESD filter without an output event handler");
315 AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillAOD(kTRUE);
316 AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillExtension(kTRUE);
321 //______________________________________________________________________________
322 TClonesArray& AliAnalysisTaskESDfilter::Cascades()
324 return *(AODEvent()->GetCascades());
327 //______________________________________________________________________________
328 TClonesArray& AliAnalysisTaskESDfilter::Tracks()
330 return *(AODEvent()->GetTracks());
333 //______________________________________________________________________________
334 TClonesArray& AliAnalysisTaskESDfilter::V0s()
336 return *(AODEvent()->GetV0s());
339 //______________________________________________________________________________
340 TClonesArray& AliAnalysisTaskESDfilter::Vertices()
342 return *(AODEvent()->GetVertices());
345 //______________________________________________________________________________
346 AliAODHeader* AliAnalysisTaskESDfilter::ConvertHeader(const AliESDEvent& esd)
348 // Convert header information
350 AliCodeTimerAuto("",0);
352 AliAODHeader* header = AODEvent()->GetHeader();
354 header->SetRunNumber(esd.GetRunNumber());
355 header->SetOfflineTrigger(fInputHandler->IsEventSelected()); // propagate the decision of the physics selection
356 header->SetNumberOfESDTracks(esd.GetNumberOfTracks());
358 TTree* tree = fInputHandler->GetTree();
360 TFile* file = tree->GetCurrentFile();
361 if (file) header->SetESDFileName(file->GetName());
365 header->SetBunchCrossNumber(0);
366 header->SetOrbitNumber(0);
367 header->SetPeriodNumber(0);
368 header->SetEventType(0);
369 header->SetMuonMagFieldScale(-999.);
370 header->SetCentrality(0);
371 header->SetEventplane(0);
373 header->SetBunchCrossNumber(esd.GetBunchCrossNumber());
374 header->SetOrbitNumber(esd.GetOrbitNumber());
375 header->SetPeriodNumber(esd.GetPeriodNumber());
376 header->SetEventType(esd.GetEventType());
378 header->SetEventNumberESDFile(esd.GetHeader()->GetEventNumberInFile());
379 if(const_cast<AliESDEvent&>(esd).GetCentrality()){
380 header->SetCentrality(const_cast<AliESDEvent&>(esd).GetCentrality());
383 header->SetCentrality(0);
385 if(const_cast<AliESDEvent&>(esd).GetEventplane()){
386 header->SetEventplane(const_cast<AliESDEvent&>(esd).GetEventplane());
389 header->SetEventplane(0);
394 header->SetFiredTriggerClasses(esd.GetFiredTriggerClasses());
395 header->SetTriggerMask(esd.GetTriggerMask());
396 header->SetTriggerCluster(esd.GetTriggerCluster());
397 header->SetL0TriggerInputs(esd.GetHeader()->GetL0TriggerInputs());
398 header->SetL1TriggerInputs(esd.GetHeader()->GetL1TriggerInputs());
399 header->SetL2TriggerInputs(esd.GetHeader()->GetL2TriggerInputs());
401 header->SetMagneticField(esd.GetMagneticField());
402 header->SetMuonMagFieldScale(esd.GetCurrentDip()/6000.);
403 header->SetZDCN1Energy(esd.GetZDCN1Energy());
404 header->SetZDCP1Energy(esd.GetZDCP1Energy());
405 header->SetZDCN2Energy(esd.GetZDCN2Energy());
406 header->SetZDCP2Energy(esd.GetZDCP2Energy());
407 header->SetZDCEMEnergy(esd.GetZDCEMEnergy(0),esd.GetZDCEMEnergy(1));
409 header->SetIRInt2InteractionMap(esd.GetHeader()->GetIRInt2InteractionMap());
410 header->SetIRInt1InteractionMap(esd.GetHeader()->GetIRInt1InteractionMap());
412 // ITS Cluster Multiplicty
413 const AliMultiplicity *mult = esd.GetMultiplicity();
414 for (Int_t ilay = 0; ilay < 6; ilay++) header->SetITSClusters(ilay, mult->GetNumberOfITSClusters(ilay));
416 // TPC only Reference Multiplicty
417 Int_t refMult = fTPCaloneTrackCuts ? (Short_t)fTPCaloneTrackCuts->GetReferenceMultiplicity(&esd, kTRUE) : -1;
418 header->SetTPConlyRefMultiplicity(refMult);
420 AliESDtrackCuts::MultEstTrackType estType = esd.GetPrimaryVertexTracks()->GetStatus() ? AliESDtrackCuts::kTrackletsITSTPC : AliESDtrackCuts::kTracklets;
421 header->SetRefMultiplicityComb05(AliESDtrackCuts::GetReferenceMultiplicity(&esd,estType,0.5));
422 header->SetRefMultiplicityComb08(AliESDtrackCuts::GetReferenceMultiplicity(&esd,estType,0.8));
424 Float_t diamxy[2]={esd.GetDiamondX(),esd.GetDiamondY()};
426 esd.GetDiamondCovXY(diamcov);
427 header->SetDiamond(diamxy,diamcov);
428 header->SetDiamondZ(esd.GetDiamondZ(),esd.GetSigma2DiamondZ());
430 // VZERO channel equalization factors for event-plane reconstruction
431 header->SetVZEROEqFactors(esd.GetVZEROEqFactors());
433 // T0 Resolution information
434 const AliESDRun* esdRun = esd.GetESDRun();
435 for (Int_t i=0;i<AliESDRun::kT0spreadSize;i++) header->SetT0spread(i,esdRun->GetT0spread(i));
440 //______________________________________________________________________________
441 void AliAnalysisTaskESDfilter::ConvertCascades(const AliESDEvent& esd)
444 // Convert the cascades part of the ESD.
445 // Return the number of cascades
447 AliCodeTimerAuto("",0);
449 // Create vertices starting from the most complex objects
452 const AliESDVertex* vtx = esd.GetPrimaryVertex();
453 Double_t pos[3] = { 0. };
454 Double_t covVtx[6] = { 0. };
455 Double_t momBach[3]={0.};
456 Double_t covTr[21]={0.};
457 Double_t pid[10]={0.};
458 AliAODPid* detpid(0x0);
459 AliAODVertex* vV0FromCascade(0x0);
460 AliAODv0* aodV0(0x0);
461 AliAODcascade* aodCascade(0x0);
462 AliAODTrack* aodTrack(0x0);
463 Double_t momPos[3]={0.};
464 Double_t momNeg[3] = { 0. };
465 Double_t momPosAtV0vtx[3]={0.};
466 Double_t momNegAtV0vtx[3]={0.};
467 Int_t tofLabel[3] = {0};
468 TClonesArray& verticesArray = Vertices();
469 TClonesArray& tracksArray = Tracks();
470 TClonesArray& cascadesArray = Cascades();
472 // Cascades (Modified by A.Maire - February 2009)
473 for (Int_t nCascade = 0; nCascade < esd.GetNumberOfCascades(); ++nCascade) {
477 AliESDcascade *esdCascade = esd.GetCascade(nCascade);
478 Int_t idxPosFromV0Dghter = esdCascade->GetPindex();
479 Int_t idxNegFromV0Dghter = esdCascade->GetNindex();
480 Int_t idxBachFromCascade = esdCascade->GetBindex();
482 AliESDtrack *esdCascadePos = esd.GetTrack( idxPosFromV0Dghter);
483 AliESDtrack *esdCascadeNeg = esd.GetTrack( idxNegFromV0Dghter);
484 AliESDtrack *esdCascadeBach = esd.GetTrack( idxBachFromCascade);
486 // Identification of the V0 within the esdCascade (via both daughter track indices)
487 AliESDv0 * currentV0 = 0x0;
488 Int_t idxV0FromCascade = -1;
490 for (Int_t iV0=0; iV0<esd.GetNumberOfV0s(); ++iV0) {
492 currentV0 = esd.GetV0(iV0);
493 Int_t posCurrentV0 = currentV0->GetPindex();
494 Int_t negCurrentV0 = currentV0->GetNindex();
496 if (posCurrentV0==idxPosFromV0Dghter && negCurrentV0==idxNegFromV0Dghter) {
497 idxV0FromCascade = iV0;
502 if(idxV0FromCascade < 0){
503 printf("Cascade - no matching for the V0 (index V0 = -1) ! Skip ... \n");
505 }// a priori, useless check, but safer ... in case of pb with tracks "out of bounds"
507 AliESDv0 *esdV0FromCascade = esd.GetV0(idxV0FromCascade);
509 // 1 - Cascade selection
511 // AliESDVertex *esdPrimVtx = new AliESDVertex(*(esd.GetPrimaryVertex()));
512 // TList cascadeObjects;
513 // cascadeObjects.AddAt(esdV0FromCascade, 0);
514 // cascadeObjects.AddAt(esdCascadePos, 1);
515 // cascadeObjects.AddAt(esdCascadeNeg, 2);
516 // cascadeObjects.AddAt(esdCascade, 3);
517 // cascadeObjects.AddAt(esdCascadeBach, 4);
518 // cascadeObjects.AddAt(esdPrimVtx, 5);
520 // UInt_t selectCascade = 0;
521 // if (fCascadeFilter) {
522 // // selectCascade = fCascadeFilter->IsSelected(&cascadeObjects);
523 // // FIXME AliESDCascadeCuts to be implemented ...
525 // // Here we may encounter a moot point at the V0 level
526 // // between the cascade selections and the V0 ones :
527 // // the V0 selected along with the cascade (secondary V0) may
528 // // usually be removed from the dedicated V0 selections (prim V0) ...
529 // // -> To be discussed !
531 // // this is a little awkward but otherwise the
532 // // list wants to access the pointer (delete it)
533 // // again when going out of scope
534 // delete cascadeObjects.RemoveAt(5); // esdPrimVtx created via copy construct
536 // if (!selectCascade)
540 // delete cascadeObjects.RemoveAt(5); // esdPrimVtx created via copy construct
544 // 2 - Add the cascade vertex
546 esdCascade->GetXYZcascade(pos[0], pos[1], pos[2]);
547 esdCascade->GetPosCovXi(covVtx);
548 chi2 = esdCascade->GetChi2Xi();
550 AliAODVertex *vCascade = new(verticesArray[fNumberOfVertices++]) AliAODVertex( pos,
552 chi2, // FIXME = Chi2/NDF will be needed
555 AliAODVertex::kCascade);
556 fPrimaryVertex->AddDaughter(vCascade);
559 // printf("---- Cascade / Cascade Vertex (AOD) : \n");
560 // vCascade->Print();
563 // if(esd.GetTOFHeader() && fIsPidOwner) fESDpid->SetTOFResponse(const_cast<AliESDEvent*>(&esd), (AliESDpid::EStartTimeType_t)fTimeZeroType); //in case of AOD production starting form LHC10e without Tender.
566 // 3 - Add the bachelor track from the cascade
568 if (!fUsedTrack[idxBachFromCascade]) {
570 esdCascadeBach->GetPxPyPz(momBach);
571 esdCascadeBach->GetXYZ(pos);
572 esdCascadeBach->GetCovarianceXYZPxPyPz(covTr);
573 esdCascadeBach->GetESDpid(pid);
574 esdCascadeBach->GetTOFLabel(tofLabel);
576 fUsedTrack[idxBachFromCascade] = kTRUE;
577 UInt_t selectInfo = 0;
578 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadeBach);
579 if (fMChandler) fMChandler->SelectParticle(esdCascadeBach->GetLabel());
580 aodTrack = new(tracksArray[fNumberOfTracks++]) AliAODTrack(esdCascadeBach->GetID(),
581 esdCascadeBach->GetLabel(),
585 kFALSE, // Why kFALSE for "isDCA" ? FIXME
587 (Short_t)esdCascadeBach->GetSign(),
588 esdCascadeBach->GetITSClusterMap(),
591 kTRUE, // usedForVtxFit = kFALSE ? FIXME
592 vtx->UsesTrack(esdCascadeBach->GetID()),
593 AliAODTrack::kSecondary,
595 aodTrack->SetTPCFitMap(esdCascadeBach->GetTPCFitMap());
596 aodTrack->SetTPCClusterMap(esdCascadeBach->GetTPCClusterMap());
597 aodTrack->SetTPCSharedMap (esdCascadeBach->GetTPCSharedMap());
598 aodTrack->SetChi2perNDF(Chi2perNDF(esdCascadeBach));
599 aodTrack->SetTPCPointsF(esdCascadeBach->GetTPCNclsF());
600 aodTrack->SetTPCNCrossedRows(UShort_t(esdCascadeBach->GetTPCCrossedRows()));
601 aodTrack->SetIntegratedLength(esdCascadeBach->GetIntegratedLength());
602 aodTrack->SetTOFLabel(tofLabel);
603 fAODTrackRefs->AddAt(aodTrack,idxBachFromCascade);
605 if (esdCascadeBach->GetSign() > 0) ++fNumberOfPositiveTracks;
606 aodTrack->ConvertAliPIDtoAODPID();
607 aodTrack->SetFlags(esdCascadeBach->GetStatus());
608 SetAODPID(esdCascadeBach,aodTrack,detpid);
611 aodTrack = static_cast<AliAODTrack*>( fAODTrackRefs->At(idxBachFromCascade) );
614 vCascade->AddDaughter(aodTrack);
617 // printf("---- Cascade / bach dghter : \n");
618 // aodTrack->Print();
622 // 4 - Add the V0 from the cascade.
623 // = V0vtx + both pos and neg daughter tracks + the aodV0 itself
626 if ( !fUsedV0[idxV0FromCascade] ) {
627 // 4.A - if VO structure hasn't been created yet
629 // 4.A.1 - Create the V0 vertex of the cascade
631 esdV0FromCascade->GetXYZ(pos[0], pos[1], pos[2]);
632 esdV0FromCascade->GetPosCov(covVtx);
633 chi2 = esdV0FromCascade->GetChi2V0(); // = chi2/NDF since NDF = 2*2-3 ?
635 vV0FromCascade = new(verticesArray[fNumberOfVertices++]) AliAODVertex(pos,
639 idxV0FromCascade, //id of ESDv0
642 // one V0 can be used by several cascades.
643 // So, one AOD V0 vtx can have several parent vtx.
644 // This is not directly allowed by AliAODvertex.
645 // Setting the parent vtx (here = param "vCascade") doesn't lead to a crash
646 // but to a problem of consistency within AODEvent.
647 // -> See below paragraph 4.B, for the proposed treatment of such a case.
649 // Add the vV0FromCascade to the aodVOVtxRefs
650 fAODV0VtxRefs->AddAt(vV0FromCascade,idxV0FromCascade);
653 // 4.A.2 - Add the positive tracks from the V0
655 esdCascadePos->GetPxPyPz(momPos);
656 esdCascadePos->GetXYZ(pos);
657 esdCascadePos->GetCovarianceXYZPxPyPz(covTr);
658 esdCascadePos->GetESDpid(pid);
659 esdCascadePos->GetTOFLabel(tofLabel);
661 if (!fUsedTrack[idxPosFromV0Dghter]) {
662 fUsedTrack[idxPosFromV0Dghter] = kTRUE;
664 UInt_t selectInfo = 0;
665 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadePos);
666 if(fMChandler) fMChandler->SelectParticle(esdCascadePos->GetLabel());
667 aodTrack = new(tracksArray[fNumberOfTracks++])
668 AliAODTrack( esdCascadePos->GetID(),
669 esdCascadePos->GetLabel(),
673 kFALSE, // Why kFALSE for "isDCA" ? FIXME
675 (Short_t)esdCascadePos->GetSign(),
676 esdCascadePos->GetITSClusterMap(),
679 kTRUE, // usedForVtxFit = kFALSE ? FIXME
680 vtx->UsesTrack(esdCascadePos->GetID()),
681 AliAODTrack::kSecondary,
683 aodTrack->SetTPCFitMap(esdCascadePos->GetTPCFitMap());
684 aodTrack->SetTPCClusterMap(esdCascadePos->GetTPCClusterMap());
685 aodTrack->SetTPCSharedMap (esdCascadePos->GetTPCSharedMap());
686 aodTrack->SetChi2perNDF(Chi2perNDF(esdCascadePos));
687 aodTrack->SetTPCPointsF(esdCascadePos->GetTPCNclsF());
688 aodTrack->SetTPCNCrossedRows(UShort_t(esdCascadePos->GetTPCCrossedRows()));
689 aodTrack->SetIntegratedLength(esdCascadePos->GetIntegratedLength());
690 aodTrack->SetTOFLabel(tofLabel);
691 fAODTrackRefs->AddAt(aodTrack,idxPosFromV0Dghter);
693 if (esdCascadePos->GetSign() > 0) ++fNumberOfPositiveTracks;
694 aodTrack->ConvertAliPIDtoAODPID();
695 aodTrack->SetFlags(esdCascadePos->GetStatus());
696 SetAODPID(esdCascadePos,aodTrack,detpid);
699 aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(idxPosFromV0Dghter));
701 vV0FromCascade->AddDaughter(aodTrack);
704 // 4.A.3 - Add the negative tracks from the V0
706 esdCascadeNeg->GetPxPyPz(momNeg);
707 esdCascadeNeg->GetXYZ(pos);
708 esdCascadeNeg->GetCovarianceXYZPxPyPz(covTr);
709 esdCascadeNeg->GetESDpid(pid);
710 esdCascadeNeg->GetTOFLabel(tofLabel);
713 if (!fUsedTrack[idxNegFromV0Dghter]) {
714 fUsedTrack[idxNegFromV0Dghter] = kTRUE;
716 UInt_t selectInfo = 0;
717 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadeNeg);
718 if(fMChandler)fMChandler->SelectParticle(esdCascadeNeg->GetLabel());
719 aodTrack = new(tracksArray[fNumberOfTracks++]) AliAODTrack( esdCascadeNeg->GetID(),
720 esdCascadeNeg->GetLabel(),
724 kFALSE, // Why kFALSE for "isDCA" ? FIXME
726 (Short_t)esdCascadeNeg->GetSign(),
727 esdCascadeNeg->GetITSClusterMap(),
730 kTRUE, // usedForVtxFit = kFALSE ? FIXME
731 vtx->UsesTrack(esdCascadeNeg->GetID()),
732 AliAODTrack::kSecondary,
734 aodTrack->SetTPCFitMap(esdCascadeNeg->GetTPCFitMap());
735 aodTrack->SetTPCClusterMap(esdCascadeNeg->GetTPCClusterMap());
736 aodTrack->SetTPCSharedMap (esdCascadeNeg->GetTPCSharedMap());
737 aodTrack->SetChi2perNDF(Chi2perNDF(esdCascadeNeg));
738 aodTrack->SetTPCPointsF(esdCascadeNeg->GetTPCNclsF());
739 aodTrack->SetTPCNCrossedRows(UShort_t(esdCascadeNeg->GetTPCCrossedRows()));
740 aodTrack->SetIntegratedLength(esdCascadeNeg->GetIntegratedLength());
741 aodTrack->SetTOFLabel(tofLabel);
742 fAODTrackRefs->AddAt(aodTrack,idxNegFromV0Dghter);
744 if (esdCascadeNeg->GetSign() > 0) ++fNumberOfPositiveTracks;
745 aodTrack->ConvertAliPIDtoAODPID();
746 aodTrack->SetFlags(esdCascadeNeg->GetStatus());
747 SetAODPID(esdCascadeNeg,aodTrack,detpid);
750 aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(idxNegFromV0Dghter));
753 vV0FromCascade->AddDaughter(aodTrack);
756 // 4.A.4 - Add the V0 from cascade to the V0 array
758 Double_t dcaV0Daughters = esdV0FromCascade->GetDcaV0Daughters();
759 Double_t dcaV0ToPrimVertex = esdV0FromCascade->GetD( esd.GetPrimaryVertex()->GetX(),
760 esd.GetPrimaryVertex()->GetY(),
761 esd.GetPrimaryVertex()->GetZ() );
762 esdV0FromCascade->GetPPxPyPz( momPosAtV0vtx[0],momPosAtV0vtx[1],momPosAtV0vtx[2] );
763 esdV0FromCascade->GetNPxPyPz( momNegAtV0vtx[0],momNegAtV0vtx[1],momNegAtV0vtx[2] );
765 Double_t dcaDaughterToPrimVertex[2] = { 999., 999.}; // ..[0] = DCA in (x,y) for Pos and ..[1] = Neg
766 dcaDaughterToPrimVertex[0] = TMath::Abs(esdCascadePos->GetD( esd.GetPrimaryVertex()->GetX(),
767 esd.GetPrimaryVertex()->GetY(),
768 esd.GetMagneticField()) );
769 dcaDaughterToPrimVertex[1] = TMath::Abs(esdCascadeNeg->GetD( esd.GetPrimaryVertex()->GetX(),
770 esd.GetPrimaryVertex()->GetY(),
771 esd.GetMagneticField()) );
773 aodV0 = new(V0s()[fNumberOfV0s++]) AliAODv0( vV0FromCascade,
778 dcaDaughterToPrimVertex);
779 // set the aod v0 on-the-fly status
780 aodV0->SetOnFlyStatus(esdV0FromCascade->GetOnFlyStatus());
782 // Add the aodV0 to the aodVORefs
783 fAODV0Refs->AddAt(aodV0,idxV0FromCascade);
785 fUsedV0[idxV0FromCascade] = kTRUE;
788 // 4.B - if V0 structure already used
791 // one V0 can be used by several cascades (frequent in PbPb evts) :
792 // same V0 which used but attached to different bachelor tracks
793 // -> aodVORefs and fAODV0VtxRefs are needed.
794 // Goal : avoid a redundancy of the info in "Vertices" and "v0s" clones array.
796 vV0FromCascade = static_cast<AliAODVertex*>( fAODV0VtxRefs->At(idxV0FromCascade) );
797 aodV0 = static_cast<AliAODv0*> ( fAODV0Refs ->At(idxV0FromCascade) );
799 // - Treatment of the parent for such a "re-used" V0 :
800 // Insert the cascade that reuses the V0 vertex in the lineage chain
801 // Before : vV0 -> vCascade1 -> vPrimary
802 // - Hyp : cascade2 uses the same V0 as cascade1
803 // After : vV0 -> vCascade2 -> vCascade1 -> vPrimary
805 AliAODVertex *vCascadePreviousParent = static_cast<AliAODVertex*> (vV0FromCascade->GetParent());
806 vV0FromCascade->SetParent(vCascade);
807 vCascade ->SetParent(vCascadePreviousParent);
810 // printf("---- Cascade / Lineage insertion\n"
811 // "Parent of V0 vtx = Cascade vtx %p\n"
812 // "Parent of the cascade vtx = Cascade vtx %p\n"
813 // "Parent of the parent cascade vtx = Cascade vtx %p\n",
814 // static_cast<void*> (vV0FromCascade->GetParent()),
815 // static_cast<void*> (vCascade->GetParent()),
816 // static_cast<void*> (vCascadePreviousParent->GetParent()) );
818 }// end if V0 structure already used
821 // printf("---- Cascade / V0 vertex: \n");
822 // vV0FromCascade->Print();
826 // printf("---- Cascade / pos dghter : \n");
827 // aodTrack->Print();
828 // printf("---- Cascade / neg dghter : \n");
829 // aodTrack->Print();
830 // printf("---- Cascade / aodV0 : \n");
834 // In any case (used V0 or not), add the V0 vertex to the cascade one.
835 vCascade->AddDaughter(vV0FromCascade);
838 // 5 - Add the primary track of the cascade (if any)
841 // 6 - Add the cascade to the AOD array of cascades
843 Double_t dcaBachToPrimVertexXY = TMath::Abs(esdCascadeBach->GetD(esd.GetPrimaryVertex()->GetX(),
844 esd.GetPrimaryVertex()->GetY(),
845 esd.GetMagneticField()) );
847 Double_t momBachAtCascadeVtx[3]={0.};
849 esdCascade->GetBPxPyPz(momBachAtCascadeVtx[0], momBachAtCascadeVtx[1], momBachAtCascadeVtx[2]);
851 aodCascade = new(cascadesArray[fNumberOfCascades++]) AliAODcascade( vCascade,
852 esdCascade->Charge(),
853 esdCascade->GetDcaXiDaughters(),
855 // DCAXiToPrimVtx -> needs to be calculated ----|
856 // doesn't exist at ESD level;
857 // See AODcascade::DcaXiToPrimVertex(Double, Double, Double)
858 dcaBachToPrimVertexXY,
863 printf("---- Cascade / AOD cascade : \n\n");
864 aodCascade->PrintXi(fPrimaryVertex->GetX(), fPrimaryVertex->GetY(), fPrimaryVertex->GetZ());
867 } // end of the loop on cascades
869 Cascades().Expand(fNumberOfCascades);
872 //______________________________________________________________________________
873 void AliAnalysisTaskESDfilter::ConvertV0s(const AliESDEvent& esd)
875 // Access to the AOD container of V0s
877 AliCodeTimerAuto("",0);
883 Double_t pos[3] = { 0. };
885 Double_t covVtx[6] = { 0. };
886 Double_t momPos[3]={0.};
887 Double_t covTr[21]={0.};
888 Double_t pid[10]={0.};
889 AliAODTrack* aodTrack(0x0);
890 AliAODPid* detpid(0x0);
891 Double_t momNeg[3]={0.};
892 Double_t momPosAtV0vtx[3]={0.};
893 Double_t momNegAtV0vtx[3]={0.};
894 Int_t tofLabel[3] = {0};
895 for (Int_t nV0 = 0; nV0 < esd.GetNumberOfV0s(); ++nV0)
897 if (fUsedV0[nV0]) continue; // skip if already added to the AOD
899 AliESDv0 *v0 = esd.GetV0(nV0);
900 Int_t posFromV0 = v0->GetPindex();
901 Int_t negFromV0 = v0->GetNindex();
905 AliESDVertex *esdVtx = new AliESDVertex(*(esd.GetPrimaryVertex()));
906 AliESDtrack *esdV0Pos = esd.GetTrack(posFromV0);
907 AliESDtrack *esdV0Neg = esd.GetTrack(negFromV0);
909 v0objects.AddAt(v0, 0);
910 v0objects.AddAt(esdV0Pos, 1);
911 v0objects.AddAt(esdV0Neg, 2);
912 v0objects.AddAt(esdVtx, 3);
915 selectV0 = fV0Filter->IsSelected(&v0objects);
916 // this is a little awkward but otherwise the
917 // list wants to access the pointer (delete it)
918 // again when going out of scope
919 delete v0objects.RemoveAt(3); // esdVtx created via copy construct
925 delete v0objects.RemoveAt(3); // esdVtx created via copy construct
929 v0->GetXYZ(pos[0], pos[1], pos[2]);
931 if (!fOldESDformat) {
932 chi2 = v0->GetChi2V0(); // = chi2/NDF since NDF = 2*2-3
933 v0->GetPosCov(covVtx);
936 for (Int_t i = 0; i < 6; i++) covVtx[i] = 0.;
941 new(Vertices()[fNumberOfVertices++]) AliAODVertex(pos,
947 fPrimaryVertex->AddDaughter(vV0);
950 // Add the positive tracks from the V0
953 esdV0Pos->GetPxPyPz(momPos);
954 esdV0Pos->GetXYZ(pos);
955 esdV0Pos->GetCovarianceXYZPxPyPz(covTr);
956 esdV0Pos->GetESDpid(pid);
957 esdV0Pos->GetTOFLabel(tofLabel);
959 const AliESDVertex *vtx = esd.GetPrimaryVertex();
961 if (!fUsedTrack[posFromV0]) {
962 fUsedTrack[posFromV0] = kTRUE;
963 UInt_t selectInfo = 0;
964 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdV0Pos);
965 if(fMChandler)fMChandler->SelectParticle(esdV0Pos->GetLabel());
966 aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdV0Pos->GetID(),
967 esdV0Pos->GetLabel(),
973 (Short_t)esdV0Pos->GetSign(),
974 esdV0Pos->GetITSClusterMap(),
977 kTRUE, // check if this is right
978 vtx->UsesTrack(esdV0Pos->GetID()),
979 AliAODTrack::kSecondary,
981 aodTrack->SetTPCFitMap(esdV0Pos->GetTPCFitMap());
982 aodTrack->SetTPCClusterMap(esdV0Pos->GetTPCClusterMap());
983 aodTrack->SetTPCSharedMap (esdV0Pos->GetTPCSharedMap());
984 aodTrack->SetChi2perNDF(Chi2perNDF(esdV0Pos));
985 aodTrack->SetTPCPointsF(esdV0Pos->GetTPCNclsF());
986 aodTrack->SetTPCNCrossedRows(UShort_t(esdV0Pos->GetTPCCrossedRows()));
987 aodTrack->SetIntegratedLength(esdV0Pos->GetIntegratedLength());
988 aodTrack->SetTOFLabel(tofLabel);
989 fAODTrackRefs->AddAt(aodTrack,posFromV0);
990 // if (fDebug > 0) printf("-------------------Bo: pos track from original pt %.3f \n",aodTrack->Pt());
991 if (esdV0Pos->GetSign() > 0) ++fNumberOfPositiveTracks;
992 aodTrack->ConvertAliPIDtoAODPID();
993 aodTrack->SetFlags(esdV0Pos->GetStatus());
994 SetAODPID(esdV0Pos,aodTrack,detpid);
997 aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(posFromV0));
998 // if (fDebug > 0) printf("-------------------Bo pos track from refArray pt %.3f \n",aodTrack->Pt());
1000 vV0->AddDaughter(aodTrack);
1002 // Add the negative tracks from the V0
1004 esdV0Neg->GetPxPyPz(momNeg);
1005 esdV0Neg->GetXYZ(pos);
1006 esdV0Neg->GetCovarianceXYZPxPyPz(covTr);
1007 esdV0Neg->GetESDpid(pid);
1008 esdV0Neg->GetTOFLabel(tofLabel);
1010 if (!fUsedTrack[negFromV0]) {
1011 fUsedTrack[negFromV0] = kTRUE;
1012 UInt_t selectInfo = 0;
1013 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdV0Neg);
1014 if(fMChandler)fMChandler->SelectParticle(esdV0Neg->GetLabel());
1015 aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdV0Neg->GetID(),
1016 esdV0Neg->GetLabel(),
1022 (Short_t)esdV0Neg->GetSign(),
1023 esdV0Neg->GetITSClusterMap(),
1026 kTRUE, // check if this is right
1027 vtx->UsesTrack(esdV0Neg->GetID()),
1028 AliAODTrack::kSecondary,
1030 aodTrack->SetTPCFitMap(esdV0Neg->GetTPCFitMap());
1031 aodTrack->SetTPCClusterMap(esdV0Neg->GetTPCClusterMap());
1032 aodTrack->SetTPCSharedMap (esdV0Neg->GetTPCSharedMap());
1033 aodTrack->SetChi2perNDF(Chi2perNDF(esdV0Neg));
1034 aodTrack->SetTPCPointsF(esdV0Neg->GetTPCNclsF());
1035 aodTrack->SetTPCNCrossedRows(UShort_t(esdV0Neg->GetTPCCrossedRows()));
1036 aodTrack->SetIntegratedLength(esdV0Neg->GetIntegratedLength());
1037 aodTrack->SetTOFLabel(tofLabel);
1038 fAODTrackRefs->AddAt(aodTrack,negFromV0);
1039 // if (fDebug > 0) printf("-------------------Bo: neg track from original pt %.3f \n",aodTrack->Pt());
1040 if (esdV0Neg->GetSign() > 0) ++fNumberOfPositiveTracks;
1041 aodTrack->ConvertAliPIDtoAODPID();
1042 aodTrack->SetFlags(esdV0Neg->GetStatus());
1043 SetAODPID(esdV0Neg,aodTrack,detpid);
1046 aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(negFromV0));
1047 // if (fDebug > 0) printf("-------------------Bo neg track from refArray pt %.3f \n",aodTrack->Pt());
1049 vV0->AddDaughter(aodTrack);
1052 // Add the V0 the V0 array as well
1054 Double_t dcaV0Daughters = v0->GetDcaV0Daughters();
1055 Double_t dcaV0ToPrimVertex = v0->GetD(esd.GetPrimaryVertex()->GetX(),
1056 esd.GetPrimaryVertex()->GetY(),
1057 esd.GetPrimaryVertex()->GetZ());
1059 v0->GetPPxPyPz(momPosAtV0vtx[0],momPosAtV0vtx[1],momPosAtV0vtx[2]);
1060 v0->GetNPxPyPz(momNegAtV0vtx[0],momNegAtV0vtx[1],momNegAtV0vtx[2]);
1062 Double_t dcaDaughterToPrimVertex[2] = { 999., 999.}; // ..[0] = DCA in (x,y) for Pos and ..[1] = Neg
1063 dcaDaughterToPrimVertex[0] = TMath::Abs(esdV0Pos->GetD( esd.GetPrimaryVertex()->GetX(),
1064 esd.GetPrimaryVertex()->GetY(),
1065 esd.GetMagneticField()) );
1066 dcaDaughterToPrimVertex[1] = TMath::Abs(esdV0Neg->GetD( esd.GetPrimaryVertex()->GetX(),
1067 esd.GetPrimaryVertex()->GetY(),
1068 esd.GetMagneticField()) );
1070 AliAODv0* aodV0 = new(V0s()[fNumberOfV0s++]) AliAODv0(vV0,
1075 dcaDaughterToPrimVertex);
1077 // set the aod v0 on-the-fly status
1078 aodV0->SetOnFlyStatus(v0->GetOnFlyStatus());
1079 }//End of loop on V0s
1081 V0s().Expand(fNumberOfV0s);
1084 //______________________________________________________________________________
1085 void AliAnalysisTaskESDfilter::ConvertTPCOnlyTracks(const AliESDEvent& esd)
1087 // Convert TPC only tracks
1088 // Here we have wo hybrid appraoch to remove fakes
1089 // ******* ITSTPC ********
1090 // Uses a cut on the ITS properties to select global tracks
1091 // which are than marked as HybdridITSTPC for the remainder
1092 // the TPC only tracks are flagged as HybridITSTPConly.
1093 // Note, in order not to get fakes back in the TPC cuts, one needs
1094 // two "ITS" cuts one tight (1) (to throw out fakes) and one lose (2) (to NOT flag the trakcs in the TPC only)
1095 // using cut number (3)
1096 // so fHybridFilterMask == (1)|(2) fTPCFilterMask = (3), Usercode needs to slect with mask = (1)|(3) and track->IsHybridITSTPC()
1097 // ******* TPC ********
1098 // 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
1099 // 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
1101 AliCodeTimerAuto("",0);
1103 // Loop over the tracks and extract and mask out all aod tracks that pass the selections for AODt racks
1104 for(int it = 0;it < fNumberOfTracks;++it)
1106 AliAODTrack *tr = (AliAODTrack*)(Tracks().At(it));
1108 UInt_t map = tr->GetFilterMap();
1109 if(map&fTPCConstrainedFilterMask){
1110 // we only reset the track select ionfo, no deletion...
1111 tr->SetFilterMap(map&~fTPCConstrainedFilterMask);
1113 if(map&fHybridFilterMaskTPCCG){
1114 // this is one part of the hybrid tracks
1115 // the others not passing the selection will be TPC only selected below
1116 tr->SetIsHybridTPCConstrainedGlobal(kTRUE);
1119 // Loop over the ESD trcks and pick out the tracks passing TPC only cuts
1122 const AliESDVertex *vtxSPD = esd.GetPrimaryVertexSPD();
1123 const AliESDVertex *vtx = esd.GetPrimaryVertex();
1125 Double_t pos[3] = { 0. };
1126 Double_t covTr[21]={0.};
1127 Double_t pid[10]={0.};
1130 Double_t p[3] = { 0. };
1132 Double_t pDCA[3] = { 0. }; // momentum at DCA
1133 Double_t rDCA[3] = { 0. }; // position at DCA
1134 Float_t dDCA[2] = {0.}; // DCA to the vertex d and z
1135 Float_t cDCA[3] = {0.}; // covariance of impact parameters
1136 Int_t tofLabel[3] = {0};
1138 AliAODTrack* aodTrack(0x0);
1139 // AliAODPid* detpid(0x0);
1141 // account for change in pT after the constraint
1142 Float_t ptMax = 1E10;
1144 for(int i = 0;i<32;i++){
1145 if(fTPCConstrainedFilterMask&(1<<i)){
1146 AliESDtrackCuts*cuts = (AliESDtrackCuts*)fTrackFilter->GetCuts()->At(i);
1147 Float_t tmp1= 0,tmp2 = 0;
1148 cuts->GetPtRange(tmp1,tmp2);
1149 if(tmp1>ptMin)ptMin=tmp1;
1150 if(tmp2<ptMax)ptMax=tmp2;
1154 for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack)
1156 AliESDtrack* esdTrack = esd.GetTrack(nTrack); //carefull do not modify it othwise need to work with a copy
1158 UInt_t selectInfo = 0;
1159 Bool_t isHybridITSTPC = false;
1163 selectInfo = fTrackFilter->IsSelected(esdTrack);
1166 if(!(selectInfo&fHybridFilterMaskTPCCG)){
1167 // not already selected tracks, use second part of hybrid tracks
1168 isHybridITSTPC = true;
1169 // too save space one could only store these...
1172 selectInfo &= fTPCConstrainedFilterMask;
1173 if (!selectInfo)continue;
1174 if (fWriteHybridTPCCOnly&&!isHybridITSTPC)continue; // write only complementary tracks
1175 // create a tpc only tracl
1176 AliESDtrack *track = AliESDtrackCuts::GetTPCOnlyTrack(const_cast<AliESDEvent*>(&esd),esdTrack->GetID());
1177 if(!track) continue;
1181 // only constrain tracks above threshold
1182 AliExternalTrackParam exParam;
1183 // take the B-field from the ESD, no 3D fieldMap available at this point
1184 Bool_t relate = false;
1185 relate = track->RelateToVertexTPC(vtxSPD,esd.GetMagneticField(),kVeryBig,&exParam);
1190 // fetch the track parameters at the DCA (unconstraint)
1191 if(track->GetTPCInnerParam()){
1192 track->GetTPCInnerParam()->GetPxPyPz(pDCA);
1193 track->GetTPCInnerParam()->GetXYZ(rDCA);
1195 // get the DCA to the vertex:
1196 track->GetImpactParametersTPC(dDCA,cDCA);
1197 // set the constrained parameters to the track
1198 track->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance());
1201 track->GetPxPyPz(p);
1203 Float_t pT = track->Pt();
1204 if(pT<ptMin||pT>ptMax){
1213 track->GetCovarianceXYZPxPyPz(covTr);
1214 esdTrack->GetESDpid(pid);// original PID
1215 esdTrack->GetTOFLabel(tofLabel);
1216 if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());
1217 aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack((track->GetID()+1)*-1,
1224 (Short_t)track->GetSign(),
1225 track->GetITSClusterMap(),
1228 kTRUE, // check if this is right
1229 vtx->UsesTrack(track->GetID()),
1230 AliAODTrack::kPrimary,
1232 aodTrack->SetIsHybridTPCConstrainedGlobal(isHybridITSTPC);
1233 aodTrack->SetTPCFitMap(track->GetTPCFitMap());
1234 aodTrack->SetTPCClusterMap(track->GetTPCClusterMap());
1235 aodTrack->SetTPCSharedMap (track->GetTPCSharedMap());
1236 aodTrack->SetIsTPCConstrained(kTRUE);
1237 aodTrack->SetChi2perNDF(Chi2perNDF(esdTrack)); // original track
1238 // set the DCA values to the AOD track
1239 aodTrack->SetPxPyPzAtDCA(pDCA[0],pDCA[1],pDCA[2]);
1240 aodTrack->SetXYAtDCA(rDCA[0],rDCA[1]);
1241 aodTrack->SetDCA(dDCA[0],dDCA[1]);
1243 aodTrack->SetFlags(track->GetStatus());
1244 aodTrack->SetTPCPointsF(track->GetTPCNclsF());
1245 aodTrack->SetTPCNCrossedRows(UShort_t(track->GetTPCCrossedRows()));
1246 aodTrack->SetIntegratedLength(track->GetIntegratedLength());
1247 aodTrack->SetTOFLabel(tofLabel);
1248 //Perform progagation of tracks if needed
1249 if(fDoPropagateTrackToEMCal) PropagateTrackToEMCal(esdTrack);
1250 aodTrack->SetTrackPhiEtaPtOnEMCal(esdTrack->GetTrackPhiOnEMCal(),esdTrack->GetTrackEtaOnEMCal(),esdTrack->GetTrackPtOnEMCal());
1252 // do not duplicate PID information
1253 // aodTrack->ConvertAliPIDtoAODPID();
1254 // SetAODPID(esdTrack,aodTrack,detpid);
1257 } // end of loop on tracks
1261 //______________________________________________________________________________
1262 void AliAnalysisTaskESDfilter::ConvertGlobalConstrainedTracks(const AliESDEvent& esd)
1265 // Here we have the option to store the complement from global constraint information
1266 // to tracks passing tight cuts (1) in order not to get fakes back in, one needs
1267 // two sets of cuts one tight (1) (to throw out fakes) and one lose (2) (fakes/bad tracks would pass (2) but not (1))
1268 // using cut number (3) selects the tracks that complement (1) e.g. tracks witout ITS refit or cluster requirement
1271 AliCodeTimerAuto("",0);
1273 // Loop over the tracks and extract and mask out all aod tracks that pass the selections for AODt racks
1274 for(int it = 0;it < fNumberOfTracks;++it)
1276 AliAODTrack *tr = (AliAODTrack*)(Tracks().At(it));
1278 UInt_t map = tr->GetFilterMap();
1279 if(map&fGlobalConstrainedFilterMask){
1280 // we only reset the track select info, no deletion...
1281 // mask reset mask in case track is already taken
1282 tr->SetFilterMap(map&~fGlobalConstrainedFilterMask);
1284 if(map&fHybridFilterMaskGCG){
1285 // this is one part of the hybrid tracks
1286 // the others not passing the selection will be the ones selected below
1287 tr->SetIsHybridGlobalConstrainedGlobal(kTRUE);
1290 // Loop over the ESD trcks and pick out the tracks passing the GlobalConstraint cuts
1293 Double_t pos[3] = { 0. };
1294 Double_t covTr[21]={0.};
1295 Double_t pid[10]={0.};
1296 Double_t p[3] = { 0. };
1298 Double_t pDCA[3] = { 0. }; // momentum at DCA
1299 Double_t rDCA[3] = { 0. }; // position at DCA
1300 Float_t dDCA[2] = {0.}; // DCA to the vertex d and z
1301 Float_t cDCA[3] = {0.}; // covariance of impact parameters
1302 Int_t tofLabel[3] = {0};
1305 AliAODTrack* aodTrack(0x0);
1306 AliAODPid* detpid(0x0);
1307 const AliESDVertex *vtx = esd.GetPrimaryVertex();
1309 // account for change in pT after the constraint
1310 Float_t ptMax = 1E10;
1312 for(int i = 0;i<32;i++){
1313 if(fGlobalConstrainedFilterMask&(1<<i)){
1314 AliESDtrackCuts*cuts = (AliESDtrackCuts*)fTrackFilter->GetCuts()->At(i);
1315 Float_t tmp1= 0,tmp2 = 0;
1316 cuts->GetPtRange(tmp1,tmp2);
1317 if(tmp1>ptMin)ptMin=tmp1;
1318 if(tmp2<ptMax)ptMax=tmp2;
1324 for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack)
1326 AliESDtrack* esdTrack = esd.GetTrack(nTrack); //carefull do not modify it othwise need to work with a copy
1327 const AliExternalTrackParam * exParamGC = esdTrack->GetConstrainedParam();
1328 if(!exParamGC)continue;
1330 UInt_t selectInfo = 0;
1331 Bool_t isHybridGC = false;
1336 selectInfo = fTrackFilter->IsSelected(esdTrack);
1340 if(!(selectInfo&fHybridFilterMaskGCG))isHybridGC = true;
1341 if (fWriteHybridGCOnly&&!isHybridGC)continue; // write only complementary tracks
1343 selectInfo &= fGlobalConstrainedFilterMask;
1344 if (!selectInfo)continue;
1345 // fetch the track parameters at the DCA (unconstrained)
1346 esdTrack->GetPxPyPz(pDCA);
1347 esdTrack->GetXYZ(rDCA);
1348 // get the DCA to the vertex:
1349 esdTrack->GetImpactParameters(dDCA,cDCA);
1351 if (!esdTrack->GetConstrainedPxPyPz(p)) continue;
1354 Float_t pT = exParamGC->Pt();
1355 if(pT<ptMin||pT>ptMax){
1360 esdTrack->GetConstrainedXYZ(pos);
1361 exParamGC->GetCovarianceXYZPxPyPz(covTr);
1362 esdTrack->GetESDpid(pid);
1363 esdTrack->GetTOFLabel(tofLabel);
1364 if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());
1365 aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack((esdTrack->GetID()+1)*-1,
1366 esdTrack->GetLabel(),
1372 (Short_t)esdTrack->GetSign(),
1373 esdTrack->GetITSClusterMap(),
1376 kTRUE, // check if this is right
1377 vtx->UsesTrack(esdTrack->GetID()),
1378 AliAODTrack::kPrimary,
1380 aodTrack->SetIsHybridGlobalConstrainedGlobal(isHybridGC);
1381 aodTrack->SetIsGlobalConstrained(kTRUE);
1382 aodTrack->SetTPCFitMap(esdTrack->GetTPCFitMap());
1383 aodTrack->SetTPCClusterMap(esdTrack->GetTPCClusterMap());
1384 aodTrack->SetTPCSharedMap (esdTrack->GetTPCSharedMap());
1385 aodTrack->SetChi2perNDF(Chi2perNDF(esdTrack));
1388 // set the DCA values to the AOD track
1389 aodTrack->SetPxPyPzAtDCA(pDCA[0],pDCA[1],pDCA[2]);
1390 aodTrack->SetXYAtDCA(rDCA[0],rDCA[1]);
1391 aodTrack->SetDCA(dDCA[0],dDCA[1]);
1393 aodTrack->SetFlags(esdTrack->GetStatus());
1394 aodTrack->SetTPCPointsF(esdTrack->GetTPCNclsF());
1395 aodTrack->SetTPCNCrossedRows(UShort_t(esdTrack->GetTPCCrossedRows()));
1396 aodTrack->SetIntegratedLength(esdTrack->GetIntegratedLength());
1397 aodTrack->SetTOFLabel(tofLabel);
1399 // only copy AOD information for hybrid, no duplicate information
1400 aodTrack->ConvertAliPIDtoAODPID();
1401 SetAODPID(esdTrack,aodTrack,detpid);
1404 //Perform progagation of tracks if needed
1405 if(fDoPropagateTrackToEMCal) PropagateTrackToEMCal(esdTrack);
1406 aodTrack->SetTrackPhiEtaPtOnEMCal(esdTrack->GetTrackPhiOnEMCal(),esdTrack->GetTrackEtaOnEMCal(),esdTrack->GetTrackPtOnEMCal());
1407 } // end of loop on tracks
1412 //______________________________________________________________________________
1413 void AliAnalysisTaskESDfilter::ConvertTracks(const AliESDEvent& esd)
1415 // Tracks (primary and orphan)
1417 AliCodeTimerAuto("",0);
1419 AliDebug(1,Form("NUMBER OF ESD TRACKS %5d\n", esd.GetNumberOfTracks()));
1421 const AliESDVertex *vtx = esd.GetPrimaryVertex();
1422 Double_t p[3] = { 0. };
1423 Double_t pos[3] = { 0. };
1424 Double_t covTr[21] = { 0. };
1425 Double_t pid[10] = { 0. };
1426 Int_t tofLabel[3] = {0};
1427 AliAODTrack* aodTrack(0x0);
1428 AliAODPid* detpid(0x0);
1430 for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack)
1432 if (fUsedTrack[nTrack]) continue;
1434 AliESDtrack *esdTrack = esd.GetTrack(nTrack);
1435 UInt_t selectInfo = 0;
1439 selectInfo = fTrackFilter->IsSelected(esdTrack);
1440 if (!selectInfo && !vtx->UsesTrack(esdTrack->GetID())) continue;
1444 esdTrack->GetPxPyPz(p);
1445 esdTrack->GetXYZ(pos);
1446 esdTrack->GetCovarianceXYZPxPyPz(covTr);
1447 esdTrack->GetESDpid(pid);
1448 esdTrack->GetTOFLabel(tofLabel);
1449 if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());
1450 fPrimaryVertex->AddDaughter(aodTrack =
1451 new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrack->GetID(),
1452 esdTrack->GetLabel(),
1458 (Short_t)esdTrack->GetSign(),
1459 esdTrack->GetITSClusterMap(),
1462 kTRUE, // check if this is right
1463 vtx->UsesTrack(esdTrack->GetID()),
1464 AliAODTrack::kPrimary,
1467 aodTrack->SetTPCFitMap(esdTrack->GetTPCFitMap());
1468 aodTrack->SetTPCClusterMap(esdTrack->GetTPCClusterMap());
1469 aodTrack->SetTPCSharedMap (esdTrack->GetTPCSharedMap());
1470 aodTrack->SetChi2perNDF(Chi2perNDF(esdTrack));
1471 aodTrack->SetTPCPointsF(esdTrack->GetTPCNclsF());
1472 aodTrack->SetTPCNCrossedRows(UShort_t(esdTrack->GetTPCCrossedRows()));
1473 aodTrack->SetIntegratedLength(esdTrack->GetIntegratedLength());
1474 aodTrack->SetTOFLabel(tofLabel);
1475 if(esdTrack->IsEMCAL()) aodTrack->SetEMCALcluster(esdTrack->GetEMCALcluster());
1476 if(esdTrack->IsPHOS()) aodTrack->SetPHOScluster(esdTrack->GetPHOScluster());
1478 //Perform progagation of tracks if needed
1479 if(fDoPropagateTrackToEMCal) PropagateTrackToEMCal(esdTrack);
1480 aodTrack->SetTrackPhiEtaPtOnEMCal(esdTrack->GetTrackPhiOnEMCal(),esdTrack->GetTrackEtaOnEMCal(),esdTrack->GetTrackPtOnEMCal());
1482 fAODTrackRefs->AddAt(aodTrack, nTrack);
1485 if (esdTrack->GetSign() > 0) ++fNumberOfPositiveTracks;
1486 aodTrack->SetFlags(esdTrack->GetStatus());
1487 aodTrack->ConvertAliPIDtoAODPID();
1488 SetAODPID(esdTrack,aodTrack,detpid);
1489 } // end of loop on tracks
1492 //______________________________________________________________________________
1493 void AliAnalysisTaskESDfilter::PropagateTrackToEMCal(AliESDtrack *esdTrack)
1495 Double_t trkPos[3] = {0.,0.,0.};
1496 Double_t EMCalEta=-999, EMCalPhi=-999, EMCalPt=-999;
1497 Double_t trkphi = esdTrack->Phi()*TMath::RadToDeg();
1498 if(TMath::Abs(esdTrack->Eta())<0.9 && trkphi > 10 && trkphi < 250 )
1500 AliExternalTrackParam *trkParam = const_cast<AliExternalTrackParam*>(esdTrack->GetInnerParam());
1503 AliExternalTrackParam trkParamTmp(*trkParam);
1504 if(AliTrackerBase::PropagateTrackToBxByBz(&trkParamTmp, fEMCalSurfaceDistance, esdTrack->GetMass(), 20, kTRUE, 0.8, -1))
1506 trkParamTmp.GetXYZ(trkPos);
1507 TVector3 trkPosVec(trkPos[0],trkPos[1],trkPos[2]);
1508 EMCalEta = trkPosVec.Eta();
1509 EMCalPhi = trkPosVec.Phi();
1510 EMCalPt = trkParamTmp.Pt();
1511 if(EMCalPhi<0) EMCalPhi += 2*TMath::Pi();
1512 esdTrack->SetTrackPhiEtaPtOnEMCal(EMCalPhi,EMCalEta,EMCalPt);
1518 //______________________________________________________________________________
1519 void AliAnalysisTaskESDfilter::ConvertPmdClusters(const AliESDEvent& esd)
1521 // Convert PMD Clusters
1522 AliCodeTimerAuto("",0);
1523 Int_t jPmdClusters=0;
1524 // Access to the AOD container of PMD clusters
1525 TClonesArray &pmdClusters = *(AODEvent()->GetPmdClusters());
1526 for (Int_t iPmd = 0; iPmd < esd.GetNumberOfPmdTracks(); ++iPmd) {
1527 // file pmd clusters, to be revised!
1528 AliESDPmdTrack *pmdTrack = esd.GetPmdTrack(iPmd);
1531 Double_t posPmd[3] = { pmdTrack->GetClusterX(), pmdTrack->GetClusterY(), pmdTrack->GetClusterZ()};
1532 Double_t pidPmd[13] = { 0.}; // to be revised!
1534 // assoc cluster not set
1535 new(pmdClusters[jPmdClusters++]) AliAODPmdCluster(iPmd, nLabel, label, pmdTrack->GetClusterADC(), posPmd, pidPmd);
1540 //______________________________________________________________________________
1541 void AliAnalysisTaskESDfilter::ConvertCaloClusters(const AliESDEvent& esd)
1543 // Convert Calorimeter Clusters
1544 AliCodeTimerAuto("",0);
1546 // Access to the AOD container of clusters
1547 TClonesArray &caloClusters = *(AODEvent()->GetCaloClusters());
1550 for (Int_t iClust=0; iClust<esd.GetNumberOfCaloClusters(); ++iClust) {
1552 AliESDCaloCluster * cluster = esd.GetCaloCluster(iClust);
1554 Int_t id = cluster->GetID();
1555 Int_t nLabel = cluster->GetNLabels();
1556 Int_t *labels = cluster->GetLabels();
1558 for(int i = 0;i < nLabel;++i){
1559 if(fMChandler)fMChandler->SelectParticle(labels[i]);
1563 Float_t energy = cluster->E();
1564 Float_t posF[3] = { 0.};
1565 cluster->GetPosition(posF);
1567 AliAODCaloCluster *caloCluster = new(caloClusters[jClusters++]) AliAODCaloCluster(id,
1573 cluster->GetType(),0);
1575 caloCluster->SetCaloCluster(cluster->GetDistanceToBadChannel(),
1576 cluster->GetDispersion(),
1577 cluster->GetM20(), cluster->GetM02(),
1578 cluster->GetEmcCpvDistance(),
1579 cluster->GetNExMax(),cluster->GetTOF()) ;
1581 caloCluster->SetPIDFromESD(cluster->GetPID());
1582 caloCluster->SetNCells(cluster->GetNCells());
1583 caloCluster->SetCellsAbsId(cluster->GetCellsAbsId());
1584 caloCluster->SetCellsAmplitudeFraction(cluster->GetCellsAmplitudeFraction());
1586 caloCluster->SetTrackDistance(cluster->GetTrackDx(), cluster->GetTrackDz());
1588 Int_t nMatchCount = 0;
1589 TArrayI* matchedT = cluster->GetTracksMatched();
1590 if (fNumberOfTracks>0 && matchedT && cluster->GetTrackMatchedIndex() >= 0) {
1591 for (Int_t im = 0; im < matchedT->GetSize(); im++) {
1592 Int_t iESDtrack = matchedT->At(im);;
1593 if (fAODTrackRefs->At(iESDtrack) != 0) {
1594 caloCluster->AddTrackMatched((AliAODTrack*)fAODTrackRefs->At(iESDtrack));
1600 caloCluster->SetTrackDistance(-999,-999);
1603 caloClusters.Expand(jClusters); // resize TObjArray to 'remove' slots for pseudo clusters
1606 //______________________________________________________________________________
1607 void AliAnalysisTaskESDfilter::ConvertCaloTrigger(TString calo, const AliESDEvent& esd)
1609 AliCodeTimerAuto("",0);
1613 AliAODCaloTrigger &aodTrigger = *(AODEvent()->GetCaloTrigger(calo));
1614 AliESDCaloTrigger &esdTrigger = *(esd.GetCaloTrigger(calo));
1616 aodTrigger.Allocate(esdTrigger.GetEntries());
1622 while (esdTrigger.Next()) {
1623 esdTrigger.GetPosition(tmod,tabsId);
1624 esdTrigger.GetAmplitude(a);
1625 aodTrigger.Add(tmod,tabsId,a,0.,(Int_t*)NULL,0,0,0);
1631 AliAODHandler *aodHandler = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
1635 TTree *aodTree = aodHandler->GetTree();
1639 Int_t *type = esd.GetCaloTriggerType();
1641 for (Int_t i = 0; i < 15; i++)
1643 aodTree->GetUserInfo()->Add(new TParameter<int>(Form("EMCALCaloTrigger%d",i), type[i]));
1648 AliAODCaloTrigger &aodTrigger = *(AODEvent()->GetCaloTrigger(calo));
1650 AliESDCaloTrigger &esdTrigger = *(esd.GetCaloTrigger(calo));
1652 aodTrigger.Allocate(esdTrigger.GetEntries());
1655 while (esdTrigger.Next())
1657 Int_t px, py, ts, nTimes, times[10], b;
1660 esdTrigger.GetPosition(px, py);
1662 esdTrigger.GetAmplitude(a);
1663 esdTrigger.GetTime(t);
1665 esdTrigger.GetL0Times(times);
1666 esdTrigger.GetNL0Times(nTimes);
1668 esdTrigger.GetL1TimeSum(ts);
1670 esdTrigger.GetTriggerBits(b);
1672 aodTrigger.Add(px, py, a, t, times, nTimes, ts, b);
1675 for (int i = 0; i < 4; i++) aodTrigger.SetL1Threshold(i, esdTrigger.GetL1Threshold(i));
1679 esdTrigger.GetL1V0(0),
1680 esdTrigger.GetL1V0(1)
1683 aodTrigger.SetL1V0(v0);
1684 aodTrigger.SetL1FrameMask(esdTrigger.GetL1FrameMask());
1687 //______________________________________________________________________________
1688 void AliAnalysisTaskESDfilter::ConvertEMCALCells(const AliESDEvent& esd)
1690 // Convert EMCAL Cells
1691 AliCodeTimerAuto("",0);
1692 // fill EMCAL cell info
1693 if (esd.GetEMCALCells()) { // protection against missing ESD information
1694 AliESDCaloCells &esdEMcells = *(esd.GetEMCALCells());
1695 Int_t nEMcell = esdEMcells.GetNumberOfCells() ;
1697 AliAODCaloCells &aodEMcells = *(AODEvent()->GetEMCALCells());
1698 aodEMcells.CreateContainer(nEMcell);
1699 aodEMcells.SetType(AliAODCaloCells::kEMCALCell);
1700 for (Int_t iCell = 0; iCell < nEMcell; iCell++) {
1701 aodEMcells.SetCell(iCell,esdEMcells.GetCellNumber(iCell),esdEMcells.GetAmplitude(iCell),
1702 esdEMcells.GetTime(iCell), esdEMcells.GetMCLabel(iCell), esdEMcells.GetEFraction(iCell));
1708 //______________________________________________________________________________
1709 void AliAnalysisTaskESDfilter::ConvertPHOSCells(const AliESDEvent& esd)
1711 // Convert PHOS Cells
1712 AliCodeTimerAuto("",0);
1713 // fill PHOS cell info
1714 if (esd.GetPHOSCells()) { // protection against missing ESD information
1715 AliESDCaloCells &esdPHcells = *(esd.GetPHOSCells());
1716 Int_t nPHcell = esdPHcells.GetNumberOfCells() ;
1718 AliAODCaloCells &aodPHcells = *(AODEvent()->GetPHOSCells());
1719 aodPHcells.CreateContainer(nPHcell);
1720 aodPHcells.SetType(AliAODCaloCells::kPHOSCell);
1721 for (Int_t iCell = 0; iCell < nPHcell; iCell++) {
1722 aodPHcells.SetCell(iCell,esdPHcells.GetCellNumber(iCell),esdPHcells.GetAmplitude(iCell),
1723 esdPHcells.GetTime(iCell), esdPHcells.GetMCLabel(iCell), esdPHcells.GetEFraction(iCell));
1729 //______________________________________________________________________________
1730 void AliAnalysisTaskESDfilter::ConvertTracklets(const AliESDEvent& esd)
1733 AliCodeTimerAuto("",0);
1735 AliAODTracklets &SPDTracklets = *(AODEvent()->GetTracklets());
1736 const AliMultiplicity *mult = esd.GetMultiplicity();
1738 if (mult->GetNumberOfTracklets()>0) {
1739 SPDTracklets.CreateContainer(mult->GetNumberOfTracklets());
1741 for (Int_t n=0; n<mult->GetNumberOfTracklets(); n++) {
1743 fMChandler->SelectParticle(mult->GetLabel(n, 0));
1744 fMChandler->SelectParticle(mult->GetLabel(n, 1));
1746 SPDTracklets.SetTracklet(n, mult->GetTheta(n), mult->GetPhi(n), mult->GetDeltaPhi(n), mult->GetLabel(n, 0),mult->GetLabel(n, 1));
1750 //Printf("ERROR: AliMultiplicity could not be retrieved from ESD");
1754 //______________________________________________________________________________
1755 void AliAnalysisTaskESDfilter::ConvertKinks(const AliESDEvent& esd)
1757 AliCodeTimerAuto("",0);
1759 // Kinks: it is a big mess the access to the information in the kinks
1760 // The loop is on the tracks in order to find the mother and daugther of each kink
1762 Double_t covTr[21]={0.};
1763 Double_t pid[10]={0.};
1764 AliAODPid* detpid(0x0);
1765 Int_t tofLabel[3] = {0};
1767 fNumberOfKinks = esd.GetNumberOfKinks();
1769 const AliESDVertex* vtx = esd.GetPrimaryVertex();
1771 for (Int_t iTrack=0; iTrack<esd.GetNumberOfTracks(); ++iTrack)
1773 AliESDtrack * esdTrack = esd.GetTrack(iTrack);
1775 Int_t ikink = esdTrack->GetKinkIndex(0);
1777 if (ikink && fNumberOfKinks) {
1778 // Negative kink index: mother, positive: daughter
1780 // Search for the second track of the kink
1782 for (Int_t jTrack = iTrack+1; jTrack<esd.GetNumberOfTracks(); ++jTrack) {
1784 AliESDtrack * esdTrack1 = esd.GetTrack(jTrack);
1786 Int_t jkink = esdTrack1->GetKinkIndex(0);
1788 if ( TMath::Abs(ikink)==TMath::Abs(jkink) ) {
1790 // The two tracks are from the same kink
1792 if (fUsedKink[TMath::Abs(ikink)-1]) continue; // skip used kinks
1795 Int_t idaughter = -1;
1797 if (ikink<0 && jkink>0) {
1802 else if (ikink>0 && jkink<0) {
1808 // cerr << "Error: Wrong combination of kink indexes: "
1809 // << ikink << " " << jkink << endl;
1813 // Add the mother track if it passed primary track selection cuts
1815 AliAODTrack * mother = NULL;
1817 UInt_t selectInfo = 0;
1819 selectInfo = fTrackFilter->IsSelected(esd.GetTrack(imother));
1820 if (!selectInfo) continue;
1823 if (!fUsedTrack[imother]) {
1825 fUsedTrack[imother] = kTRUE;
1827 AliESDtrack *esdTrackM = esd.GetTrack(imother);
1828 Double_t p[3] = { 0. };
1829 Double_t pos[3] = { 0. };
1830 esdTrackM->GetPxPyPz(p);
1831 esdTrackM->GetXYZ(pos);
1832 esdTrackM->GetCovarianceXYZPxPyPz(covTr);
1833 esdTrackM->GetESDpid(pid);
1834 esdTrackM->GetTOFLabel(tofLabel);
1835 if(fMChandler)fMChandler->SelectParticle(esdTrackM->GetLabel());
1837 new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrackM->GetID(),
1838 esdTrackM->GetLabel(),
1844 (Short_t)esdTrackM->GetSign(),
1845 esdTrackM->GetITSClusterMap(),
1848 kTRUE, // check if this is right
1849 vtx->UsesTrack(esdTrack->GetID()),
1850 AliAODTrack::kPrimary,
1852 mother->SetTPCFitMap(esdTrackM->GetTPCFitMap());
1853 mother->SetTPCClusterMap(esdTrackM->GetTPCClusterMap());
1854 mother->SetTPCSharedMap (esdTrackM->GetTPCSharedMap());
1855 mother->SetChi2perNDF(Chi2perNDF(esdTrackM));
1856 mother->SetTPCPointsF(esdTrackM->GetTPCNclsF());
1857 mother->SetTPCNCrossedRows(UShort_t(esdTrackM->GetTPCCrossedRows()));
1858 mother->SetIntegratedLength(esdTrackM->GetIntegratedLength());
1859 mother->SetTOFLabel(tofLabel);
1860 fAODTrackRefs->AddAt(mother, imother);
1862 if (esdTrackM->GetSign() > 0) ++fNumberOfPositiveTracks;
1863 mother->SetFlags(esdTrackM->GetStatus());
1864 mother->ConvertAliPIDtoAODPID();
1865 fPrimaryVertex->AddDaughter(mother);
1866 mother->ConvertAliPIDtoAODPID();
1867 SetAODPID(esdTrackM,mother,detpid);
1870 // cerr << "Error: event " << esd.GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1
1871 // << " track " << imother << " has already been used!" << endl;
1874 // Add the kink vertex
1875 AliESDkink * kink = esd.GetKink(TMath::Abs(ikink)-1);
1877 AliAODVertex * vkink =
1878 new(Vertices()[fNumberOfVertices++]) AliAODVertex(kink->GetPosition(),
1882 esdTrack->GetID(), // This is the track ID of the mother's track!
1883 AliAODVertex::kKink);
1884 // Add the daughter track
1886 AliAODTrack * daughter = NULL;
1888 if (!fUsedTrack[idaughter]) {
1890 fUsedTrack[idaughter] = kTRUE;
1892 AliESDtrack *esdTrackD = esd.GetTrack(idaughter);
1893 Double_t p[3] = { 0. };
1894 Double_t pos[3] = { 0. };
1896 esdTrackD->GetPxPyPz(p);
1897 esdTrackD->GetXYZ(pos);
1898 esdTrackD->GetCovarianceXYZPxPyPz(covTr);
1899 esdTrackD->GetESDpid(pid);
1900 esdTrackD->GetTOFLabel(tofLabel);
1902 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdTrackD);
1903 if(fMChandler)fMChandler->SelectParticle(esdTrackD->GetLabel());
1905 new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrackD->GetID(),
1906 esdTrackD->GetLabel(),
1912 (Short_t)esdTrackD->GetSign(),
1913 esdTrackD->GetITSClusterMap(),
1916 kTRUE, // check if this is right
1917 vtx->UsesTrack(esdTrack->GetID()),
1918 AliAODTrack::kSecondary,
1920 daughter->SetTPCFitMap(esdTrackD->GetTPCFitMap());
1921 daughter->SetTPCClusterMap(esdTrackD->GetTPCClusterMap());
1922 daughter->SetTPCSharedMap (esdTrackD->GetTPCSharedMap());
1923 daughter->SetTPCPointsF(esdTrackD->GetTPCNclsF());
1924 daughter->SetTPCNCrossedRows(UShort_t(esdTrackD->GetTPCCrossedRows()));
1925 daughter->SetIntegratedLength(esdTrackD->GetIntegratedLength());
1926 daughter->SetTOFLabel(tofLabel);
1927 fAODTrackRefs->AddAt(daughter, idaughter);
1929 if (esdTrackD->GetSign() > 0) ++fNumberOfPositiveTracks;
1930 daughter->SetFlags(esdTrackD->GetStatus());
1931 daughter->ConvertAliPIDtoAODPID();
1932 vkink->AddDaughter(daughter);
1933 daughter->ConvertAliPIDtoAODPID();
1934 SetAODPID(esdTrackD,daughter,detpid);
1937 // cerr << "Error: event " << esd.GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1
1938 // << " track " << idaughter << " has already been used!" << endl;
1946 //______________________________________________________________________________
1947 void AliAnalysisTaskESDfilter::ConvertPrimaryVertices(const AliESDEvent& esd)
1949 AliCodeTimerAuto("",0);
1951 // Access to the AOD container of vertices
1952 fNumberOfVertices = 0;
1954 Double_t pos[3] = { 0. };
1955 Double_t covVtx[6] = { 0. };
1957 // Add primary vertex. The primary tracks will be defined
1958 // after the loops on the composite objects (V0, cascades, kinks)
1959 const AliESDVertex *vtx = esd.GetPrimaryVertex();
1961 vtx->GetXYZ(pos); // position
1962 vtx->GetCovMatrix(covVtx); //covariance matrix
1964 fPrimaryVertex = new(Vertices()[fNumberOfVertices++])
1965 AliAODVertex(pos, covVtx, vtx->GetChi2toNDF(), NULL, -1, AliAODVertex::kPrimary);
1966 fPrimaryVertex->SetName(vtx->GetName());
1967 fPrimaryVertex->SetTitle(vtx->GetTitle());
1968 fPrimaryVertex->SetBC(vtx->GetBC());
1970 TString vtitle = vtx->GetTitle();
1971 if (!vtitle.Contains("VertexerTracks"))
1972 fPrimaryVertex->SetNContributors(vtx->GetNContributors());
1974 if (fDebug > 0) fPrimaryVertex->Print();
1976 // Add SPD "main" vertex
1977 const AliESDVertex *vtxS = esd.GetPrimaryVertexSPD();
1978 vtxS->GetXYZ(pos); // position
1979 vtxS->GetCovMatrix(covVtx); //covariance matrix
1980 AliAODVertex * mVSPD = new(Vertices()[fNumberOfVertices++])
1981 AliAODVertex(pos, covVtx, vtxS->GetChi2toNDF(), NULL, -1, AliAODVertex::kMainSPD);
1982 mVSPD->SetName(vtxS->GetName());
1983 mVSPD->SetTitle(vtxS->GetTitle());
1984 mVSPD->SetNContributors(vtxS->GetNContributors());
1986 // Add SPD pileup vertices
1987 for(Int_t iV=0; iV<esd.GetNumberOfPileupVerticesSPD(); ++iV)
1989 const AliESDVertex *vtxP = esd.GetPileupVertexSPD(iV);
1990 vtxP->GetXYZ(pos); // position
1991 vtxP->GetCovMatrix(covVtx); //covariance matrix
1992 AliAODVertex * pVSPD = new(Vertices()[fNumberOfVertices++])
1993 AliAODVertex(pos, covVtx, vtxP->GetChi2toNDF(), NULL, -1, AliAODVertex::kPileupSPD);
1994 pVSPD->SetName(vtxP->GetName());
1995 pVSPD->SetTitle(vtxP->GetTitle());
1996 pVSPD->SetNContributors(vtxP->GetNContributors());
1997 pVSPD->SetBC(vtxP->GetBC());
2000 // Add TRK pileup vertices
2001 for(Int_t iV=0; iV<esd.GetNumberOfPileupVerticesTracks(); ++iV)
2003 const AliESDVertex *vtxP = esd.GetPileupVertexTracks(iV);
2004 vtxP->GetXYZ(pos); // position
2005 vtxP->GetCovMatrix(covVtx); //covariance matrix
2006 AliAODVertex * pVTRK = new(Vertices()[fNumberOfVertices++])
2007 AliAODVertex(pos, covVtx, vtxP->GetChi2toNDF(), NULL, -1, AliAODVertex::kPileupTracks);
2008 pVTRK->SetName(vtxP->GetName());
2009 pVTRK->SetTitle(vtxP->GetTitle());
2010 pVTRK->SetNContributors(vtxP->GetNContributors());
2011 pVTRK->SetBC(vtxP->GetBC());
2014 // Add TPC "main" vertex
2015 const AliESDVertex *vtxT = esd.GetPrimaryVertexTPC();
2016 vtxT->GetXYZ(pos); // position
2017 vtxT->GetCovMatrix(covVtx); //covariance matrix
2018 AliAODVertex * mVTPC = new(Vertices()[fNumberOfVertices++])
2019 AliAODVertex(pos, covVtx, vtxT->GetChi2toNDF(), NULL, -1, AliAODVertex::kMainTPC);
2020 mVTPC->SetName(vtxT->GetName());
2021 mVTPC->SetTitle(vtxT->GetTitle());
2022 mVTPC->SetNContributors(vtxT->GetNContributors());
2027 //______________________________________________________________________________
2028 void AliAnalysisTaskESDfilter::ConvertVZERO(const AliESDEvent& esd)
2030 // Convert VZERO data
2031 AliAODVZERO* vzeroData = AODEvent()->GetVZEROData();
2032 *vzeroData = *(esd.GetVZEROData());
2035 //______________________________________________________________________________
2036 void AliAnalysisTaskESDfilter::ConvertTZERO(const AliESDEvent& esd)
2038 // Convert TZERO data
2039 const AliESDTZERO* esdTzero = esd.GetESDTZERO();
2040 AliAODTZERO* aodTzero = AODEvent()->GetTZEROData();
2042 for (Int_t icase=0; icase<3; icase++){
2043 aodTzero->SetT0TOF( icase, esdTzero->GetT0TOF(icase));
2044 aodTzero->SetT0TOFbest(icase, esdTzero->GetT0TOFbest(icase));
2046 aodTzero->SetBackgroundFlag(esdTzero->GetBackgroundFlag());
2047 aodTzero->SetPileupFlag(esdTzero->GetPileupFlag());
2048 aodTzero->SetSatelliteFlag(esdTzero->GetSatellite());
2050 Float_t rawTime[24];
2051 for(Int_t ipmt=0; ipmt<24; ipmt++)
2052 rawTime[ipmt] = esdTzero->GetTimeFull(ipmt,0);
2054 Int_t idxOfFirstPmtA = -1, idxOfFirstPmtC = -1;
2055 Float_t timeOfFirstPmtA = 9999, timeOfFirstPmtC = 9999;
2056 for(int ipmt=0; ipmt<12; ipmt++){
2057 if( rawTime[ipmt] > -200 && rawTime[ipmt] < timeOfFirstPmtC && rawTime[ipmt]!=0){
2058 timeOfFirstPmtC = rawTime[ipmt];
2059 idxOfFirstPmtC = ipmt;
2062 for(int ipmt=12; ipmt<24; ipmt++){
2063 if( rawTime[ipmt] > -200 && rawTime[ipmt] < timeOfFirstPmtA && rawTime[ipmt]!=0 ){
2064 timeOfFirstPmtA = rawTime[ipmt];
2065 idxOfFirstPmtA = ipmt;
2069 if(idxOfFirstPmtA != -1 && idxOfFirstPmtC != -1){
2070 //speed of light in cm/ns TMath::C()*1e-7
2071 Float_t vertexraw = TMath::C()*1e-7 * (rawTime[idxOfFirstPmtA] - rawTime[idxOfFirstPmtC])/2;
2072 aodTzero->SetT0VertexRaw( vertexraw );
2074 aodTzero->SetT0VertexRaw(99999);
2077 aodTzero->SetT0zVertex(esdTzero->GetT0zVertex());
2080 const Double32_t *amp=esdTzero->GetT0amplitude();
2081 for(int ipmt=0; ipmt<24; ipmt++)
2082 aodTzero->SetAmp(ipmt, amp[ipmt]);
2083 aodTzero->SetAmp(24,esdTzero->GetMultC() );
2084 aodTzero->SetAmp(25,esdTzero->GetMultA() );
2089 //______________________________________________________________________________
2090 void AliAnalysisTaskESDfilter::ConvertZDC(const AliESDEvent& esd)
2093 AliESDZDC* esdZDC = esd.GetZDCData();
2095 const Double_t zem1Energy = esdZDC->GetZEM1Energy();
2096 const Double_t zem2Energy = esdZDC->GetZEM2Energy();
2098 const Double_t *towZNC = esdZDC->GetZNCTowerEnergy();
2099 const Double_t *towZPC = esdZDC->GetZPCTowerEnergy();
2100 const Double_t *towZNA = esdZDC->GetZNATowerEnergy();
2101 const Double_t *towZPA = esdZDC->GetZPATowerEnergy();
2102 const Double_t *towZNCLG = esdZDC->GetZNCTowerEnergyLR();
2103 const Double_t *towZNALG = esdZDC->GetZNATowerEnergyLR();
2105 AliAODZDC* zdcAOD = AODEvent()->GetZDCData();
2107 zdcAOD->SetZEM1Energy(zem1Energy);
2108 zdcAOD->SetZEM2Energy(zem2Energy);
2109 zdcAOD->SetZNCTowers(towZNC, towZNCLG);
2110 zdcAOD->SetZNATowers(towZNA, towZNALG);
2111 zdcAOD->SetZPCTowers(towZPC);
2112 zdcAOD->SetZPATowers(towZPA);
2114 zdcAOD->SetZDCParticipants(esdZDC->GetZDCParticipants(), esdZDC->GetZDCPartSideA(), esdZDC->GetZDCPartSideC());
2115 zdcAOD->SetZDCImpactParameter(esdZDC->GetImpactParameter(), esdZDC->GetImpactParamSideA(),
2116 esdZDC->GetImpactParamSideC());
2117 zdcAOD->SetZDCTDCSum(esdZDC->GetZNTDCSum(0));
2118 zdcAOD->SetZDCTDCDiff(esdZDC->GetZNTDCDiff(0));
2119 if(esdZDC->IsZNChit()) zdcAOD->SetZNCTDC(esdZDC->GetZDCTDCCorrected(10,0));
2120 if(esdZDC->IsZNAhit()) zdcAOD->SetZNATDC(esdZDC->GetZDCTDCCorrected(12,0));
2123 //_______________________________________________________________________________________________________________________________________
2124 Int_t AliAnalysisTaskESDfilter::ConvertHMPID(const AliESDEvent& esd) // clm
2127 // Convtert ESD HMPID info to AOD and return the number of good tracks with HMPID signal.
2128 // We need to return an int since there is no signal counter in the ESD.
2131 AliCodeTimerAuto("",0);
2133 Int_t cntHmpidGoodTracks = 0;
2142 Float_t thetaTrk = 0;
2145 Double_t hmpPid[5]={0};
2146 Double_t hmpMom[3]={0};
2148 TClonesArray &hmpidRings = *(AODEvent()->GetHMPIDrings());
2150 for (Int_t iTrack=0; iTrack<esd.GetNumberOfTracks(); ++iTrack)
2152 if(! esd.GetTrack(iTrack) ) continue;
2154 if(esd.GetTrack(iTrack)->GetHMPIDsignal() > -20 ) { //
2156 (esd.GetTrack(iTrack))->GetHMPIDmip(xMip, yMip, qMip, nphMip); // Get MIP properties
2157 (esd.GetTrack(iTrack))->GetHMPIDtrk(xTrk,yTrk,thetaTrk,phiTrk);
2158 (esd.GetTrack(iTrack))->GetHMPIDpid(hmpPid);
2159 if((esd.GetTrack(iTrack))->GetOuterHmpParam()) (esd.GetTrack(iTrack))->GetOuterHmpPxPyPz(hmpMom);
2161 if(esd.GetTrack(iTrack)->GetHMPIDsignal() == 0 && thetaTrk == 0 && qMip == 0 && nphMip ==0 ) continue; //
2163 new(hmpidRings[cntHmpidGoodTracks++]) AliAODHMPIDrings(
2164 (esd.GetTrack(iTrack))->GetID(), // Unique track id to attach the ring to
2165 1000000*nphMip+qMip, // MIP charge and number of photons
2166 (esd.GetTrack(iTrack))->GetHMPIDcluIdx(), // 1000000*chamber id + cluster idx of the assigned MIP cluster
2167 thetaTrk, // track inclination angle theta
2168 phiTrk, // track inclination angle phi
2169 (esd.GetTrack(iTrack))->GetHMPIDsignal(), // Cherenkov angle
2170 (esd.GetTrack(iTrack))->GetHMPIDoccupancy(), // Occupancy claculated for the given chamber
2171 (esd.GetTrack(iTrack))->GetHMPIDchi2(), // Ring resolution squared
2172 xTrk, // Track x coordinate (LORS)
2173 yTrk, // Track y coordinate (LORS)
2174 xMip, // MIP x coordinate (LORS)
2175 yMip, // MIP y coordinate (LORS)
2176 hmpPid, // PID probablities from ESD, remove later once it is in CombinedPid
2177 hmpMom // Track momentum in HMPID at ring reconstruction
2180 // Printf(Form("+++++++++ yes/no: %d %lf %lf %lf %lf %lf %lf ",(esd.GetTrack(iTrack))->IsHMPID(),thetaTrk, (esd.GetTrack(iTrack))->GetHMPIDchi2(),xTrk, yTrk , xMip, yMip));
2183 }// HMPID signal > -20
2184 }//___esd track loop
2186 return cntHmpidGoodTracks;
2189 //______________________________________________________________________________
2190 void AliAnalysisTaskESDfilter::ConvertESDtoAOD()
2192 // ESD Filter analysis task executed for each event
2194 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
2198 AliCodeTimerAuto("",0);
2200 fOldESDformat = ( esd->GetAliESDOld() != 0x0 );
2202 // Reconstruct cascades and V0 here
2203 if (fIsV0CascadeRecoEnabled) {
2204 esd->ResetCascades();
2207 AliV0vertexer lV0vtxer;
2208 AliCascadeVertexer lCascVtxer;
2210 lV0vtxer.SetCuts(fV0Cuts);
2211 lCascVtxer.SetCuts(fCascadeCuts);
2214 lV0vtxer.Tracks2V0vertices(esd);
2215 lCascVtxer.V0sTracks2CascadeVertices(esd);
2219 fNumberOfTracks = 0;
2220 fNumberOfPositiveTracks = 0;
2222 fNumberOfVertices = 0;
2223 fNumberOfCascades = 0;
2226 AliAODHeader* header = ConvertHeader(*esd);
2228 if ( fIsVZEROEnabled ) ConvertVZERO(*esd);
2229 if ( fIsTZEROEnabled ) ConvertTZERO(*esd);
2231 // Fetch Stack for debuggging if available
2235 fMChandler = (AliMCEventHandler*) ((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler());
2238 // loop over events and fill them
2239 // Multiplicity information needed by the header (to be revised!)
2240 Int_t nTracks = esd->GetNumberOfTracks();
2241 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) esd->GetTrack(iTrack)->SetESDEvent(esd);
2243 // Update the header
2245 Int_t nV0s = esd->GetNumberOfV0s();
2246 Int_t nCascades = esd->GetNumberOfCascades();
2247 Int_t nKinks = esd->GetNumberOfKinks();
2248 Int_t nVertices = nV0s + nCascades /*V0 wihtin cascade already counted*/+ nKinks + 1 /* = prim. vtx*/;
2249 Int_t nPileSPDVertices=1+esd->GetNumberOfPileupVerticesSPD(); // also SPD main vertex
2250 Int_t nPileTrkVertices=esd->GetNumberOfPileupVerticesTracks();
2251 nVertices+=nPileSPDVertices;
2252 nVertices+=nPileTrkVertices;
2254 Int_t nCaloClus = esd->GetNumberOfCaloClusters();
2256 Int_t nPmdClus = esd->GetNumberOfPmdTracks();
2257 Int_t nHmpidRings = 0;
2259 AliDebug(1,Form(" NV0=%d NCASCADES=%d NKINKS=%d", nV0s, nCascades, nKinks));
2261 AODEvent()->ResetStd(nTracks, nVertices, nV0s, nCascades, nJets, nCaloClus, nFmdClus, nPmdClus,nHmpidRings);
2265 // RefArray to store a mapping between esd V0 number and newly created AOD-Vertex V0
2266 fAODV0VtxRefs = new TRefArray(nV0s);
2267 // RefArray to store the mapping between esd V0 number and newly created AOD-V0
2268 fAODV0Refs = new TRefArray(nV0s);
2269 // Array to take into account the V0s already added to the AOD (V0 within cascades)
2270 fUsedV0 = new Bool_t[nV0s];
2271 for (Int_t iV0=0; iV0<nV0s; ++iV0) fUsedV0[iV0]=kFALSE;
2276 // RefArray to store the mapping between esd track number and newly created AOD-Track
2278 fAODTrackRefs = new TRefArray(nTracks);
2280 // Array to take into account the tracks already added to the AOD
2281 fUsedTrack = new Bool_t[nTracks];
2282 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) fUsedTrack[iTrack]=kFALSE;
2285 // Array to take into account the kinks already added to the AOD
2288 fUsedKink = new Bool_t[nKinks];
2289 for (Int_t iKink=0; iKink<nKinks; ++iKink) fUsedKink[iKink]=kFALSE;
2292 ConvertPrimaryVertices(*esd);
2294 //setting best TOF PID
2295 AliESDInputHandler* esdH = dynamic_cast<AliESDInputHandler*>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
2297 fESDpid = esdH->GetESDpid();
2299 if (fIsPidOwner && fESDpid){
2304 { //in case of no Tender attached
2305 fESDpid = new AliESDpid;
2306 fIsPidOwner = kTRUE;
2309 if(!esd->GetTOFHeader())
2310 { //protection in case the pass2 LHC10b,c,d have been processed without tender.
2311 Float_t t0spread[10];
2312 Float_t intrinsicTOFres=100; //ps ok for LHC10b,c,d pass2!!
2313 for (Int_t i=0; i<10; i++) t0spread[i] = (TMath::Sqrt(esd->GetSigma2DiamondZ()))/0.03; //0.03 to convert from cm to ps
2314 fESDpid->GetTOFResponse().SetT0resolution(t0spread);
2315 fESDpid->GetTOFResponse().SetTimeResolution(intrinsicTOFres);
2316 // fESDpid->SetTOFResponse(esd, (AliESDpid::EStartTimeType_t)fTimeZeroType);
2317 AliTOFHeader tmpTOFHeader(0,t0spread[0],0,NULL,NULL,NULL,intrinsicTOFres,t0spread[0]);
2318 AODEvent()->SetTOFHeader(&tmpTOFHeader); // write dummy TOF header in AOD
2320 AODEvent()->SetTOFHeader(esd->GetTOFHeader()); // write TOF header in AOD
2323 // if(esd->GetTOFHeader() && fIsPidOwner) fESDpid->SetTOFResponse(esd, (AliESDpid::EStartTimeType_t)fTimeZeroType); //in case of AOD production strating form LHC10e without Tender.
2325 if ( fAreCascadesEnabled ) ConvertCascades(*esd);
2327 if ( fAreV0sEnabled ) ConvertV0s(*esd);
2329 if ( fAreKinksEnabled ) ConvertKinks(*esd);
2331 if ( fAreTracksEnabled ) ConvertTracks(*esd);
2333 // Update number of AOD tracks in header at the end of track loop (M.G.)
2334 header->SetRefMultiplicity(fNumberOfTracks);
2335 header->SetRefMultiplicityPos(fNumberOfPositiveTracks);
2336 header->SetRefMultiplicityNeg(fNumberOfTracks - fNumberOfPositiveTracks);
2338 if ( fTPCConstrainedFilterMask ) ConvertTPCOnlyTracks(*esd);
2339 if( fGlobalConstrainedFilterMask) ConvertGlobalConstrainedTracks(*esd);
2341 if ( fArePmdClustersEnabled ) ConvertPmdClusters(*esd);
2343 if ( fAreCaloClustersEnabled ) ConvertCaloClusters(*esd);
2345 if ( fAreEMCALCellsEnabled )ConvertEMCALCells(*esd);
2347 if ( fArePHOSCellsEnabled )ConvertPHOSCells(*esd);
2349 if ( fAreEMCALTriggerEnabled )ConvertCaloTrigger(TString("EMCAL"), *esd);
2351 if ( fArePHOSTriggerEnabled )ConvertCaloTrigger(TString("PHOS"), *esd);
2353 if ( fAreTrackletsEnabled ) ConvertTracklets(*esd);
2354 if ( fIsZDCEnabled ) ConvertZDC(*esd);
2356 if(fIsHMPIDEnabled) nHmpidRings = ConvertHMPID(*esd);
2358 delete fAODTrackRefs; fAODTrackRefs=0x0;
2359 delete fAODV0VtxRefs; fAODV0VtxRefs=0x0;
2360 delete fAODV0Refs; fAODV0Refs=0x0;
2362 delete[] fUsedTrack; fUsedTrack=0x0;
2363 delete[] fUsedV0; fUsedV0=0x0;
2364 delete[] fUsedKink; fUsedKink=0x0;
2375 //______________________________________________________________________________
2376 void AliAnalysisTaskESDfilter::SetAODPID(AliESDtrack *esdtrack, AliAODTrack *aodtrack, AliAODPid *detpid)
2379 // Setter for the raw PID detector signals
2382 // Save PID object for candidate electrons
2383 Bool_t pidSave = kFALSE;
2385 Bool_t selectInfo = fTrackFilter->IsSelected((char*) "Electrons");
2386 if (selectInfo) pidSave = kTRUE;
2390 // Tracks passing pt cut
2391 if(esdtrack->Pt()>fHighPthreshold) {
2395 if(esdtrack->Pt()> fPtshape->GetXmin()){
2396 Double_t y = fPtshape->Eval(esdtrack->Pt())/fPtshape->Eval(fHighPthreshold);
2397 if(gRandom->Rndm(0)<1./y){
2401 }//end if p function
2405 if(!aodtrack->GetDetPid()){// prevent memory leak when calling SetAODPID twice for the same track
2406 detpid = new AliAODPid();
2407 SetDetectorRawSignals(detpid,esdtrack);
2408 aodtrack->SetDetPID(detpid);
2413 //______________________________________________________________________________
2414 void AliAnalysisTaskESDfilter::SetDetectorRawSignals(AliAODPid *aodpid, AliESDtrack *track)
2417 //assignment of the detector signals (AliXXXesdPID inspired)
2420 AliInfo("no ESD track found. .....exiting");
2424 aodpid->SetTPCmomentum(track->GetTPCmomentum());
2425 aodpid->SetTPCTgl(track->GetTPCTgl());
2426 // const AliExternalTrackParam *in=track->GetInnerParam();
2428 // aodpid->SetTPCmomentum(in->GetP());
2430 // aodpid->SetTPCmomentum(-1.);
2434 aodpid->SetITSsignal(track->GetITSsignal());
2435 Double_t itsdedx[4]; // dE/dx samples for individual ITS layers
2436 track->GetITSdEdxSamples(itsdedx);
2437 aodpid->SetITSdEdxSamples(itsdedx);
2439 aodpid->SetTPCsignal(track->GetTPCsignal());
2440 aodpid->SetTPCsignalN(track->GetTPCsignalN());
2441 if(track->GetTPCdEdxInfo()) aodpid->SetTPCdEdxInfo(track->GetTPCdEdxInfo());
2444 Int_t nslices = track->GetNumberOfTRDslices()*6;
2445 TArrayD trdslices(nslices);
2446 for(Int_t iSl =0; iSl < track->GetNumberOfTRDslices(); iSl++) {
2447 for(Int_t iPl =0; iPl<6; iPl++) trdslices[iPl*track->GetNumberOfTRDslices()+iSl] = track->GetTRDslice(iPl,iSl);
2451 for(Int_t iPl=0;iPl<6;iPl++){
2452 Double_t trdmom=track->GetTRDmomentum(iPl);
2453 aodpid->SetTRDmomentum(iPl,trdmom);
2456 aodpid->SetTRDslices(track->GetNumberOfTRDslices()*6,trdslices.GetArray());
2457 aodpid->SetTRDsignal(track->GetTRDsignal());
2459 //TRD clusters and tracklets
2460 aodpid->SetTRDncls(track->GetTRDncls());
2461 aodpid->SetTRDntrackletsPID(track->GetTRDntrackletsPID());
2463 aodpid->SetTRDChi2(track->GetTRDchi2());
2466 Double_t times[AliPID::kSPECIES]; track->GetIntegratedTimes(times);
2467 aodpid->SetIntegratedTimes(times);
2469 // Float_t tzeroTrack = fESDpid->GetTOFResponse().GetStartTime(track->P());
2470 // aodpid->SetTOFsignal(track->GetTOFsignal()-tzeroTrack);
2471 aodpid->SetTOFsignal(track->GetTOFsignal());
2474 for (Int_t iMass=0; iMass<5; iMass++){
2475 // tofRes[iMass]=(Double_t)fESDpid->GetTOFResponse().GetExpectedSigma(track->P(), times[iMass], AliPID::ParticleMass(iMass));
2476 tofRes[iMass]=0; //backward compatibility
2478 aodpid->SetTOFpidResolution(tofRes);
2480 // aodpid->SetHMPIDsignal(0); // set to zero for compression but it will be removed later
2484 Double_t AliAnalysisTaskESDfilter::Chi2perNDF(AliESDtrack* track)
2486 // Calculate chi2 per ndf for track
2487 Int_t nClustersTPC = track->GetTPCNcls();
2489 if ( nClustersTPC > 5) {
2490 return (track->GetTPCchi2()/Float_t(nClustersTPC - 5));
2497 //______________________________________________________________________________
2498 void AliAnalysisTaskESDfilter::Terminate(Option_t */*option*/)
2500 // Terminate analysis
2502 if (fDebug > 1) printf("AnalysisESDfilter: Terminate() \n");
2505 //______________________________________________________________________________
2506 void AliAnalysisTaskESDfilter::PrintMCInfo(AliStack *pStack,Int_t label){
2509 label = TMath::Abs(label);
2510 TParticle *part = pStack->Particle(label);
2511 Printf("########################");
2512 Printf("%s:%d %d UniqueID %d PDG %d P %3.3f",(char*)__FILE__,__LINE__,label,part->GetUniqueID(),part->GetPdgCode(),part->P());
2514 TParticle* mother = part;
2515 Int_t imo = part->GetFirstMother();
2516 Int_t nprim = pStack->GetNprimary();
2517 // while((imo >= nprim) && (mother->GetUniqueID() == 4)) {
2518 while((imo >= nprim)) {
2519 mother = pStack->Particle(imo);
2520 Printf("Mother %s:%d Label %d UniqueID %d PDG %d P %3.3f",(char*)__FILE__,__LINE__,imo,mother->GetUniqueID(),mother->GetPdgCode(),mother->P());
2522 imo = mother->GetFirstMother();
2524 Printf("########################");
2527 //______________________________________________________