1 /**************************************************************************
\r
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
\r
4 * Author: The ALICE Off-line Project. *
\r
5 * Contributors are mentioned in the code where appropriate. *
\r
7 * Permission to use, copy, modify and distribute this software and its *
\r
8 * documentation strictly for non-commercial purposes is hereby granted *
\r
9 * without fee, provided that the above copyright notice appears in all *
\r
10 * copies and that both the copyright notice and this permission notice *
\r
11 * appear in the supporting documentation. The authors make no claims *
\r
12 * about the suitability of this software for any purpose. It is *
\r
13 * provided "as is" without express or implied warranty. *
\r
14 **************************************************************************/
\r
21 #include <TArrayI.h>
\r
22 #include <TRandom.h>
\r
23 #include <TParticle.h>
\r
26 #include "AliAnalysisTaskESDfilter.h"
\r
27 #include "AliAnalysisManager.h"
\r
28 #include "AliESDEvent.h"
\r
29 #include "AliESDRun.h"
\r
30 #include "AliStack.h"
\r
31 #include "AliAODEvent.h"
\r
32 #include "AliMCEvent.h"
\r
33 #include "AliMCEventHandler.h"
\r
34 #include "AliESDInputHandler.h"
\r
35 #include "AliAODHandler.h"
\r
36 #include "AliAODMCParticle.h"
\r
37 #include "AliAnalysisFilter.h"
\r
38 #include "AliESDMuonTrack.h"
\r
39 #include "AliESDVertex.h"
\r
40 #include "AliCentrality.h"
\r
41 #include "AliEventplane.h"
\r
42 #include "AliESDv0.h"
\r
43 #include "AliESDkink.h"
\r
44 #include "AliESDcascade.h"
\r
45 #include "AliESDPmdTrack.h"
\r
46 #include "AliESDCaloCluster.h"
\r
47 #include "AliESDCaloCells.h"
\r
48 #include "AliMultiplicity.h"
\r
50 #include "AliCodeTimer.h"
\r
51 #include "AliESDtrackCuts.h"
\r
52 #include "AliESDpid.h"
\r
53 #include "Riostream.h"
\r
55 ClassImp(AliAnalysisTaskESDfilter)
\r
57 ////////////////////////////////////////////////////////////////////////
\r
59 AliAnalysisTaskESDfilter::AliAnalysisTaskESDfilter():
\r
60 AliAnalysisTaskSE(),
\r
64 fCascadeFilter(0x0),
\r
67 fEnableFillAOD(kTRUE),
\r
76 fNumberOfPositiveTracks(0),
\r
78 fNumberOfVertices(0),
\r
79 fNumberOfCascades(0),
\r
81 fOldESDformat(kFALSE),
\r
82 fPrimaryVertex(0x0),
\r
83 fTPCConstrainedFilterMask(0),
\r
84 fHybridFilterMaskTPCCG(0),
\r
85 fWriteHybridTPCCOnly(kFALSE),
\r
86 fGlobalConstrainedFilterMask(0),
\r
87 fHybridFilterMaskGCG(0),
\r
88 fWriteHybridGCOnly(kFALSE),
\r
89 fIsVZEROEnabled(kTRUE),
\r
90 fIsTZEROEnabled(kTRUE),
\r
91 fIsZDCEnabled(kTRUE),
\r
92 fAreCascadesEnabled(kTRUE),
\r
93 fAreV0sEnabled(kTRUE),
\r
94 fAreKinksEnabled(kTRUE),
\r
95 fAreTracksEnabled(kTRUE),
\r
96 fArePmdClustersEnabled(kTRUE),
\r
97 fAreCaloClustersEnabled(kTRUE),
\r
98 fAreEMCALCellsEnabled(kTRUE),
\r
99 fArePHOSCellsEnabled(kTRUE),
\r
100 fAreTrackletsEnabled(kTRUE),
\r
102 fIsPidOwner(kFALSE),
\r
103 fTimeZeroType(AliESDpid::kTOF_T0),
\r
104 fTPCaloneTrackCuts(0)
\r
106 // Default constructor
\r
109 //______________________________________________________________________________
\r
110 AliAnalysisTaskESDfilter::AliAnalysisTaskESDfilter(const char* name):
\r
111 AliAnalysisTaskSE(name),
\r
115 fCascadeFilter(0x0),
\r
116 fHighPthreshold(0),
\r
118 fEnableFillAOD(kTRUE),
\r
122 fAODTrackRefs(0x0),
\r
123 fAODV0VtxRefs(0x0),
\r
126 fNumberOfTracks(0),
\r
127 fNumberOfPositiveTracks(0),
\r
129 fNumberOfVertices(0),
\r
130 fNumberOfCascades(0),
\r
132 fOldESDformat(kFALSE),
\r
133 fPrimaryVertex(0x0),
\r
134 fTPCConstrainedFilterMask(0),
\r
135 fHybridFilterMaskTPCCG(0),
\r
136 fWriteHybridTPCCOnly(kFALSE),
\r
137 fGlobalConstrainedFilterMask(0),
\r
138 fHybridFilterMaskGCG(0),
\r
139 fWriteHybridGCOnly(kFALSE),
\r
140 fIsVZEROEnabled(kTRUE),
\r
141 fIsTZEROEnabled(kTRUE),
\r
142 fIsZDCEnabled(kTRUE),
\r
143 fAreCascadesEnabled(kTRUE),
\r
144 fAreV0sEnabled(kTRUE),
\r
145 fAreKinksEnabled(kTRUE),
\r
146 fAreTracksEnabled(kTRUE),
\r
147 fArePmdClustersEnabled(kTRUE),
\r
148 fAreCaloClustersEnabled(kTRUE),
\r
149 fAreEMCALCellsEnabled(kTRUE),
\r
150 fArePHOSCellsEnabled(kTRUE),
\r
151 fAreTrackletsEnabled(kTRUE),
\r
153 fIsPidOwner(kFALSE),
\r
154 fTimeZeroType(AliESDpid::kTOF_T0),
\r
155 fTPCaloneTrackCuts(0)
\r
159 AliAnalysisTaskESDfilter::~AliAnalysisTaskESDfilter(){
\r
160 if(fIsPidOwner) delete fESDpid;
\r
162 //______________________________________________________________________________
\r
163 void AliAnalysisTaskESDfilter::UserCreateOutputObjects()
\r
166 // Create Output Objects conenct filter to outputtree
\r
170 OutputTree()->GetUserInfo()->Add(fTrackFilter);
\r
174 AliError("No OutputTree() for adding the track filter");
\r
176 fTPCaloneTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
\r
179 //______________________________________________________________________________
\r
180 void AliAnalysisTaskESDfilter::Init()
\r
183 if (fDebug > 1) AliInfo("Init() \n");
\r
184 // Call configuration file
\r
187 //______________________________________________________________________________
\r
188 void AliAnalysisTaskESDfilter::PrintTask(Option_t *option, Int_t indent) const
\r
190 // Print selection task information
\r
193 AliAnalysisTaskSE::PrintTask(option,indent);
\r
195 TString spaces(' ',indent+3);
\r
197 cout << spaces.Data() << Form("Cascades are %s",fAreCascadesEnabled ? "ENABLED":"DISABLED") << endl;
\r
198 cout << spaces.Data() << Form("V0s are %s",fAreV0sEnabled ? "ENABLED":"DISABLED") << endl;
\r
199 cout << spaces.Data() << Form("Kinks are %s",fAreKinksEnabled ? "ENABLED":"DISABLED") << endl;
\r
200 cout << spaces.Data() << Form("Tracks are %s",fAreTracksEnabled ? "ENABLED":"DISABLED") << endl;
\r
201 cout << spaces.Data() << Form("PmdClusters are %s",fArePmdClustersEnabled ? "ENABLED":"DISABLED") << endl;
\r
202 cout << spaces.Data() << Form("CaloClusters are %s",fAreCaloClustersEnabled ? "ENABLED":"DISABLED") << endl;
\r
203 cout << spaces.Data() << Form("EMCAL cells are %s",fAreEMCALCellsEnabled ? "ENABLED":"DISABLED") << endl;
\r
204 cout << spaces.Data() << Form("Tracklets are %s",fAreTrackletsEnabled ? "ENABLED":"DISABLED") << endl;
\r
207 //______________________________________________________________________________
\r
208 void AliAnalysisTaskESDfilter::UserExec(Option_t */*option*/)
\r
210 // Execute analysis for current event
\r
213 Long64_t ientry = Entry();
\r
216 printf("Filter: Analysing event # %5d\n", (Int_t) ientry);
\r
217 if (fHighPthreshold == 0) AliInfo("detector PID signals are stored in each track");
\r
218 if (!fPtshape) AliInfo("detector PID signals are not stored below the pt threshold");
\r
220 // Filters must explicitely enable AOD filling in their UserExec (AG)
\r
221 if (!AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()) AliFatal("Cannot run ESD filter without an output event handler");
\r
222 if(fEnableFillAOD) {
\r
223 AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillAOD(kTRUE);
\r
224 AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillExtension(kTRUE);
\r
229 //______________________________________________________________________________
\r
230 TClonesArray& AliAnalysisTaskESDfilter::Cascades()
\r
232 return *(AODEvent()->GetCascades());
\r
235 //______________________________________________________________________________
\r
236 TClonesArray& AliAnalysisTaskESDfilter::Tracks()
\r
238 return *(AODEvent()->GetTracks());
\r
241 //______________________________________________________________________________
\r
242 TClonesArray& AliAnalysisTaskESDfilter::V0s()
\r
244 return *(AODEvent()->GetV0s());
\r
247 //______________________________________________________________________________
\r
248 TClonesArray& AliAnalysisTaskESDfilter::Vertices()
\r
250 return *(AODEvent()->GetVertices());
\r
253 //______________________________________________________________________________
\r
254 AliAODHeader* AliAnalysisTaskESDfilter::ConvertHeader(const AliESDEvent& esd)
\r
256 // Convert header information
\r
258 AliCodeTimerAuto("",0);
\r
260 AliAODHeader* header = AODEvent()->GetHeader();
\r
262 header->SetRunNumber(esd.GetRunNumber());
\r
263 header->SetOfflineTrigger(fInputHandler->IsEventSelected()); // propagate the decision of the physics selection
\r
265 TTree* tree = fInputHandler->GetTree();
\r
267 TFile* file = tree->GetCurrentFile();
\r
268 if (file) header->SetESDFileName(file->GetName());
\r
271 if (fOldESDformat) {
\r
272 header->SetBunchCrossNumber(0);
\r
273 header->SetOrbitNumber(0);
\r
274 header->SetPeriodNumber(0);
\r
275 header->SetEventType(0);
\r
276 header->SetMuonMagFieldScale(-999.);
\r
277 header->SetCentrality(0);
\r
278 header->SetEventplane(0);
\r
280 header->SetBunchCrossNumber(esd.GetBunchCrossNumber());
\r
281 header->SetOrbitNumber(esd.GetOrbitNumber());
\r
282 header->SetPeriodNumber(esd.GetPeriodNumber());
\r
283 header->SetEventType(esd.GetEventType());
\r
285 header->SetEventNumberESDFile(esd.GetHeader()->GetEventNumberInFile());
\r
286 if(const_cast<AliESDEvent&>(esd).GetCentrality()){
\r
287 header->SetCentrality(const_cast<AliESDEvent&>(esd).GetCentrality());
\r
290 header->SetCentrality(0);
\r
292 if(const_cast<AliESDEvent&>(esd).GetEventplane()){
\r
293 header->SetEventplane(const_cast<AliESDEvent&>(esd).GetEventplane());
\r
296 header->SetEventplane(0);
\r
301 header->SetFiredTriggerClasses(esd.GetFiredTriggerClasses());
\r
302 header->SetTriggerMask(esd.GetTriggerMask());
\r
303 header->SetTriggerCluster(esd.GetTriggerCluster());
\r
304 header->SetL0TriggerInputs(esd.GetHeader()->GetL0TriggerInputs());
\r
305 header->SetL1TriggerInputs(esd.GetHeader()->GetL1TriggerInputs());
\r
306 header->SetL2TriggerInputs(esd.GetHeader()->GetL2TriggerInputs());
\r
308 header->SetMagneticField(esd.GetMagneticField());
\r
309 header->SetMuonMagFieldScale(esd.GetCurrentDip()/6000.);
\r
310 header->SetZDCN1Energy(esd.GetZDCN1Energy());
\r
311 header->SetZDCP1Energy(esd.GetZDCP1Energy());
\r
312 header->SetZDCN2Energy(esd.GetZDCN2Energy());
\r
313 header->SetZDCP2Energy(esd.GetZDCP2Energy());
\r
314 header->SetZDCEMEnergy(esd.GetZDCEMEnergy(0),esd.GetZDCEMEnergy(1));
\r
316 // ITS Cluster Multiplicty
\r
317 const AliMultiplicity *mult = esd.GetMultiplicity();
\r
318 for (Int_t ilay = 0; ilay < 6; ilay++) header->SetITSClusters(ilay, mult->GetNumberOfITSClusters(ilay));
\r
320 // TPC only Reference Multiplicty
\r
321 Int_t refMult = fTPCaloneTrackCuts ? (Short_t)fTPCaloneTrackCuts->GetReferenceMultiplicity(&esd, kTRUE) : -1;
\r
322 header->SetTPConlyRefMultiplicity(refMult);
\r
325 Float_t diamxy[2]={esd.GetDiamondX(),esd.GetDiamondY()};
\r
326 Float_t diamcov[3];
\r
327 esd.GetDiamondCovXY(diamcov);
\r
328 header->SetDiamond(diamxy,diamcov);
\r
329 header->SetDiamondZ(esd.GetDiamondZ(),esd.GetSigma2DiamondZ());
\r
331 // VZERO channel equalization factors for event-plane reconstruction
\r
332 header->SetVZEROEqFactors(esd.GetVZEROEqFactors());
\r
337 //______________________________________________________________________________
\r
338 void AliAnalysisTaskESDfilter::ConvertCascades(const AliESDEvent& esd)
\r
340 // Convert the cascades part of the ESD.
\r
341 // Return the number of cascades
\r
343 AliCodeTimerAuto("",0);
\r
345 // Create vertices starting from the most complex objects
\r
346 Double_t chi2 = 0.;
\r
348 const AliESDVertex* vtx = esd.GetPrimaryVertex();
\r
349 Double_t pos[3] = { 0. };
\r
350 Double_t covVtx[6] = { 0. };
\r
351 Double_t momBach[3]={0.};
\r
352 Double_t covTr[21]={0.};
\r
353 Double_t pid[10]={0.};
\r
354 AliAODPid* detpid(0x0);
\r
355 AliAODVertex* vV0FromCascade(0x0);
\r
356 AliAODv0* aodV0(0x0);
\r
357 AliAODcascade* aodCascade(0x0);
\r
358 AliAODTrack* aodTrack(0x0);
\r
359 Double_t momPos[3]={0.};
\r
360 Double_t momNeg[3] = { 0. };
\r
361 Double_t momPosAtV0vtx[3]={0.};
\r
362 Double_t momNegAtV0vtx[3]={0.};
\r
364 TClonesArray& verticesArray = Vertices();
\r
365 TClonesArray& tracksArray = Tracks();
\r
366 TClonesArray& cascadesArray = Cascades();
\r
368 // Cascades (Modified by A.Maire - February 2009)
\r
369 for (Int_t nCascade = 0; nCascade < esd.GetNumberOfCascades(); ++nCascade) {
\r
373 AliESDcascade *esdCascade = esd.GetCascade(nCascade);
\r
374 Int_t idxPosFromV0Dghter = esdCascade->GetPindex();
\r
375 Int_t idxNegFromV0Dghter = esdCascade->GetNindex();
\r
376 Int_t idxBachFromCascade = esdCascade->GetBindex();
\r
378 AliESDtrack *esdCascadePos = esd.GetTrack( idxPosFromV0Dghter);
\r
379 AliESDtrack *esdCascadeNeg = esd.GetTrack( idxNegFromV0Dghter);
\r
380 AliESDtrack *esdCascadeBach = esd.GetTrack( idxBachFromCascade);
\r
382 // Identification of the V0 within the esdCascade (via both daughter track indices)
\r
383 AliESDv0 * currentV0 = 0x0;
\r
384 Int_t idxV0FromCascade = -1;
\r
386 for (Int_t iV0=0; iV0<esd.GetNumberOfV0s(); ++iV0) {
\r
388 currentV0 = esd.GetV0(iV0);
\r
389 Int_t posCurrentV0 = currentV0->GetPindex();
\r
390 Int_t negCurrentV0 = currentV0->GetNindex();
\r
392 if (posCurrentV0==idxPosFromV0Dghter && negCurrentV0==idxNegFromV0Dghter) {
\r
393 idxV0FromCascade = iV0;
\r
398 if(idxV0FromCascade < 0){
\r
399 printf("Cascade - no matching for the V0 (index V0 = -1) ! Skip ... \n");
\r
401 }// a priori, useless check, but safer ... in case of pb with tracks "out of bounds"
\r
403 AliESDv0 *esdV0FromCascade = esd.GetV0(idxV0FromCascade);
\r
405 // 1 - Cascade selection
\r
407 // AliESDVertex *esdPrimVtx = new AliESDVertex(*(esd.GetPrimaryVertex()));
\r
408 // TList cascadeObjects;
\r
409 // cascadeObjects.AddAt(esdV0FromCascade, 0);
\r
410 // cascadeObjects.AddAt(esdCascadePos, 1);
\r
411 // cascadeObjects.AddAt(esdCascadeNeg, 2);
\r
412 // cascadeObjects.AddAt(esdCascade, 3);
\r
413 // cascadeObjects.AddAt(esdCascadeBach, 4);
\r
414 // cascadeObjects.AddAt(esdPrimVtx, 5);
\r
416 // UInt_t selectCascade = 0;
\r
417 // if (fCascadeFilter) {
\r
418 // // selectCascade = fCascadeFilter->IsSelected(&cascadeObjects);
\r
419 // // FIXME AliESDCascadeCuts to be implemented ...
\r
421 // // Here we may encounter a moot point at the V0 level
\r
422 // // between the cascade selections and the V0 ones :
\r
423 // // the V0 selected along with the cascade (secondary V0) may
\r
424 // // usually be removed from the dedicated V0 selections (prim V0) ...
\r
425 // // -> To be discussed !
\r
427 // // this is a little awkward but otherwise the
\r
428 // // list wants to access the pointer (delete it)
\r
429 // // again when going out of scope
\r
430 // delete cascadeObjects.RemoveAt(5); // esdPrimVtx created via copy construct
\r
432 // if (!selectCascade)
\r
436 // delete cascadeObjects.RemoveAt(5); // esdPrimVtx created via copy construct
\r
440 // 2 - Add the cascade vertex
\r
442 esdCascade->GetXYZcascade(pos[0], pos[1], pos[2]);
\r
443 esdCascade->GetPosCovXi(covVtx);
\r
444 chi2 = esdCascade->GetChi2Xi();
\r
446 AliAODVertex *vCascade = new(verticesArray[fNumberOfVertices++]) AliAODVertex( pos,
\r
448 chi2, // FIXME = Chi2/NDF will be needed
\r
451 AliAODVertex::kCascade);
\r
452 fPrimaryVertex->AddDaughter(vCascade);
\r
454 // if (fDebug > 2) {
\r
455 // printf("---- Cascade / Cascade Vertex (AOD) : \n");
\r
456 // vCascade->Print();
\r
459 if(esd.GetTOFHeader() && fIsPidOwner) fESDpid->SetTOFResponse(const_cast<AliESDEvent*>(&esd), (AliESDpid::EStartTimeType_t)fTimeZeroType); //in case of AOD production strating form LHC10e without Tender.
\r
462 // 3 - Add the bachelor track from the cascade
\r
464 if (!fUsedTrack[idxBachFromCascade]) {
\r
466 esdCascadeBach->GetPxPyPz(momBach);
\r
467 esdCascadeBach->GetXYZ(pos);
\r
468 esdCascadeBach->GetCovarianceXYZPxPyPz(covTr);
\r
469 esdCascadeBach->GetESDpid(pid);
\r
471 fUsedTrack[idxBachFromCascade] = kTRUE;
\r
472 UInt_t selectInfo = 0;
\r
473 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadeBach);
\r
474 if (fMChandler) fMChandler->SelectParticle(esdCascadeBach->GetLabel());
\r
475 aodTrack = new(tracksArray[fNumberOfTracks++]) AliAODTrack(esdCascadeBach->GetID(),
\r
476 esdCascadeBach->GetLabel(),
\r
480 kFALSE, // Why kFALSE for "isDCA" ? FIXME
\r
482 (Short_t)esdCascadeBach->GetSign(),
\r
483 esdCascadeBach->GetITSClusterMap(),
\r
486 kTRUE, // usedForVtxFit = kFALSE ? FIXME
\r
487 vtx->UsesTrack(esdCascadeBach->GetID()),
\r
488 AliAODTrack::kSecondary,
\r
490 aodTrack->SetTPCFitMap(esdCascadeBach->GetTPCFitMap());
\r
491 aodTrack->SetTPCClusterMap(esdCascadeBach->GetTPCClusterMap());
\r
492 aodTrack->SetTPCSharedMap (esdCascadeBach->GetTPCSharedMap());
\r
493 aodTrack->SetChi2perNDF(Chi2perNDF(esdCascadeBach));
\r
494 aodTrack->SetTPCPointsF(esdCascadeBach->GetTPCNclsF());
\r
495 fAODTrackRefs->AddAt(aodTrack,idxBachFromCascade);
\r
497 if (esdCascadeBach->GetSign() > 0) ++fNumberOfPositiveTracks;
\r
498 aodTrack->ConvertAliPIDtoAODPID();
\r
499 aodTrack->SetFlags(esdCascadeBach->GetStatus());
\r
500 SetAODPID(esdCascadeBach,aodTrack,detpid);
\r
503 aodTrack = static_cast<AliAODTrack*>( fAODTrackRefs->At(idxBachFromCascade) );
\r
506 vCascade->AddDaughter(aodTrack);
\r
508 // if (fDebug > 4) {
\r
509 // printf("---- Cascade / bach dghter : \n");
\r
510 // aodTrack->Print();
\r
514 // 4 - Add the V0 from the cascade.
\r
515 // = V0vtx + both pos and neg daughter tracks + the aodV0 itself
\r
518 if ( !fUsedV0[idxV0FromCascade] ) {
\r
519 // 4.A - if VO structure hasn't been created yet
\r
521 // 4.A.1 - Create the V0 vertex of the cascade
\r
523 esdV0FromCascade->GetXYZ(pos[0], pos[1], pos[2]);
\r
524 esdV0FromCascade->GetPosCov(covVtx);
\r
525 chi2 = esdV0FromCascade->GetChi2V0(); // = chi2/NDF since NDF = 2*2-3 ?
\r
527 vV0FromCascade = new(verticesArray[fNumberOfVertices++]) AliAODVertex(pos,
\r
531 idxV0FromCascade, //id of ESDv0
\r
532 AliAODVertex::kV0);
\r
534 // one V0 can be used by several cascades.
\r
535 // So, one AOD V0 vtx can have several parent vtx.
\r
536 // This is not directly allowed by AliAODvertex.
\r
537 // Setting the parent vtx (here = param "vCascade") doesn't lead to a crash
\r
538 // but to a problem of consistency within AODEvent.
\r
539 // -> See below paragraph 4.B, for the proposed treatment of such a case.
\r
541 // Add the vV0FromCascade to the aodVOVtxRefs
\r
542 fAODV0VtxRefs->AddAt(vV0FromCascade,idxV0FromCascade);
\r
545 // 4.A.2 - Add the positive tracks from the V0
\r
547 esdCascadePos->GetPxPyPz(momPos);
\r
548 esdCascadePos->GetXYZ(pos);
\r
549 esdCascadePos->GetCovarianceXYZPxPyPz(covTr);
\r
550 esdCascadePos->GetESDpid(pid);
\r
553 if (!fUsedTrack[idxPosFromV0Dghter]) {
\r
554 fUsedTrack[idxPosFromV0Dghter] = kTRUE;
\r
556 UInt_t selectInfo = 0;
\r
557 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadePos);
\r
558 if(fMChandler) fMChandler->SelectParticle(esdCascadePos->GetLabel());
\r
559 aodTrack = new(tracksArray[fNumberOfTracks++])
\r
560 AliAODTrack( esdCascadePos->GetID(),
\r
561 esdCascadePos->GetLabel(),
\r
565 kFALSE, // Why kFALSE for "isDCA" ? FIXME
\r
567 (Short_t)esdCascadePos->GetSign(),
\r
568 esdCascadePos->GetITSClusterMap(),
\r
571 kTRUE, // usedForVtxFit = kFALSE ? FIXME
\r
572 vtx->UsesTrack(esdCascadePos->GetID()),
\r
573 AliAODTrack::kSecondary,
\r
575 aodTrack->SetTPCFitMap(esdCascadePos->GetTPCFitMap());
\r
576 aodTrack->SetTPCClusterMap(esdCascadePos->GetTPCClusterMap());
\r
577 aodTrack->SetTPCSharedMap (esdCascadePos->GetTPCSharedMap());
\r
578 aodTrack->SetChi2perNDF(Chi2perNDF(esdCascadePos));
\r
579 aodTrack->SetTPCPointsF(esdCascadePos->GetTPCNclsF());
\r
580 fAODTrackRefs->AddAt(aodTrack,idxPosFromV0Dghter);
\r
582 if (esdCascadePos->GetSign() > 0) ++fNumberOfPositiveTracks;
\r
583 aodTrack->ConvertAliPIDtoAODPID();
\r
584 aodTrack->SetFlags(esdCascadePos->GetStatus());
\r
585 SetAODPID(esdCascadePos,aodTrack,detpid);
\r
588 aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(idxPosFromV0Dghter));
\r
590 vV0FromCascade->AddDaughter(aodTrack);
\r
593 // 4.A.3 - Add the negative tracks from the V0
\r
595 esdCascadeNeg->GetPxPyPz(momNeg);
\r
596 esdCascadeNeg->GetXYZ(pos);
\r
597 esdCascadeNeg->GetCovarianceXYZPxPyPz(covTr);
\r
598 esdCascadeNeg->GetESDpid(pid);
\r
601 if (!fUsedTrack[idxNegFromV0Dghter]) {
\r
602 fUsedTrack[idxNegFromV0Dghter] = kTRUE;
\r
604 UInt_t selectInfo = 0;
\r
605 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdCascadeNeg);
\r
606 if(fMChandler)fMChandler->SelectParticle(esdCascadeNeg->GetLabel());
\r
607 aodTrack = new(tracksArray[fNumberOfTracks++]) AliAODTrack( esdCascadeNeg->GetID(),
\r
608 esdCascadeNeg->GetLabel(),
\r
612 kFALSE, // Why kFALSE for "isDCA" ? FIXME
\r
614 (Short_t)esdCascadeNeg->GetSign(),
\r
615 esdCascadeNeg->GetITSClusterMap(),
\r
618 kTRUE, // usedForVtxFit = kFALSE ? FIXME
\r
619 vtx->UsesTrack(esdCascadeNeg->GetID()),
\r
620 AliAODTrack::kSecondary,
\r
622 aodTrack->SetTPCFitMap(esdCascadeNeg->GetTPCFitMap());
\r
623 aodTrack->SetTPCClusterMap(esdCascadeNeg->GetTPCClusterMap());
\r
624 aodTrack->SetTPCSharedMap (esdCascadeNeg->GetTPCSharedMap());
\r
625 aodTrack->SetChi2perNDF(Chi2perNDF(esdCascadeNeg));
\r
626 aodTrack->SetTPCPointsF(esdCascadeNeg->GetTPCNclsF());
\r
627 fAODTrackRefs->AddAt(aodTrack,idxNegFromV0Dghter);
\r
629 if (esdCascadeNeg->GetSign() > 0) ++fNumberOfPositiveTracks;
\r
630 aodTrack->ConvertAliPIDtoAODPID();
\r
631 aodTrack->SetFlags(esdCascadeNeg->GetStatus());
\r
632 SetAODPID(esdCascadeNeg,aodTrack,detpid);
\r
635 aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(idxNegFromV0Dghter));
\r
638 vV0FromCascade->AddDaughter(aodTrack);
\r
641 // 4.A.4 - Add the V0 from cascade to the V0 array
\r
643 Double_t dcaV0Daughters = esdV0FromCascade->GetDcaV0Daughters();
\r
644 Double_t dcaV0ToPrimVertex = esdV0FromCascade->GetD( esd.GetPrimaryVertex()->GetX(),
\r
645 esd.GetPrimaryVertex()->GetY(),
\r
646 esd.GetPrimaryVertex()->GetZ() );
\r
647 esdV0FromCascade->GetPPxPyPz( momPosAtV0vtx[0],momPosAtV0vtx[1],momPosAtV0vtx[2] );
\r
648 esdV0FromCascade->GetNPxPyPz( momNegAtV0vtx[0],momNegAtV0vtx[1],momNegAtV0vtx[2] );
\r
650 Double_t dcaDaughterToPrimVertex[2] = { 999., 999.}; // ..[0] = DCA in (x,y) for Pos and ..[1] = Neg
\r
651 dcaDaughterToPrimVertex[0] = TMath::Abs(esdCascadePos->GetD( esd.GetPrimaryVertex()->GetX(),
\r
652 esd.GetPrimaryVertex()->GetY(),
\r
653 esd.GetMagneticField()) );
\r
654 dcaDaughterToPrimVertex[1] = TMath::Abs(esdCascadeNeg->GetD( esd.GetPrimaryVertex()->GetX(),
\r
655 esd.GetPrimaryVertex()->GetY(),
\r
656 esd.GetMagneticField()) );
\r
658 aodV0 = new(V0s()[fNumberOfV0s++]) AliAODv0( vV0FromCascade,
\r
660 dcaV0ToPrimVertex,
\r
663 dcaDaughterToPrimVertex);
\r
664 // set the aod v0 on-the-fly status
\r
665 aodV0->SetOnFlyStatus(esdV0FromCascade->GetOnFlyStatus());
\r
667 // Add the aodV0 to the aodVORefs
\r
668 fAODV0Refs->AddAt(aodV0,idxV0FromCascade);
\r
670 fUsedV0[idxV0FromCascade] = kTRUE;
\r
673 // 4.B - if V0 structure already used
\r
676 // one V0 can be used by several cascades (frequent in PbPb evts) :
\r
677 // same V0 which used but attached to different bachelor tracks
\r
678 // -> aodVORefs and fAODV0VtxRefs are needed.
\r
679 // Goal : avoid a redundancy of the info in "Vertices" and "v0s" clones array.
\r
681 vV0FromCascade = static_cast<AliAODVertex*>( fAODV0VtxRefs->At(idxV0FromCascade) );
\r
682 aodV0 = static_cast<AliAODv0*> ( fAODV0Refs ->At(idxV0FromCascade) );
\r
684 // - Treatment of the parent for such a "re-used" V0 :
\r
685 // Insert the cascade that reuses the V0 vertex in the lineage chain
\r
686 // Before : vV0 -> vCascade1 -> vPrimary
\r
687 // - Hyp : cascade2 uses the same V0 as cascade1
\r
688 // After : vV0 -> vCascade2 -> vCascade1 -> vPrimary
\r
690 AliAODVertex *vCascadePreviousParent = static_cast<AliAODVertex*> (vV0FromCascade->GetParent());
\r
691 vV0FromCascade->SetParent(vCascade);
\r
692 vCascade ->SetParent(vCascadePreviousParent);
\r
695 // printf("---- Cascade / Lineage insertion\n"
\r
696 // "Parent of V0 vtx = Cascade vtx %p\n"
\r
697 // "Parent of the cascade vtx = Cascade vtx %p\n"
\r
698 // "Parent of the parent cascade vtx = Cascade vtx %p\n",
\r
699 // static_cast<void*> (vV0FromCascade->GetParent()),
\r
700 // static_cast<void*> (vCascade->GetParent()),
\r
701 // static_cast<void*> (vCascadePreviousParent->GetParent()) );
\r
703 }// end if V0 structure already used
\r
705 // if (fDebug > 2) {
\r
706 // printf("---- Cascade / V0 vertex: \n");
\r
707 // vV0FromCascade->Print();
\r
710 // if (fDebug > 4) {
\r
711 // printf("---- Cascade / pos dghter : \n");
\r
712 // aodTrack->Print();
\r
713 // printf("---- Cascade / neg dghter : \n");
\r
714 // aodTrack->Print();
\r
715 // printf("---- Cascade / aodV0 : \n");
\r
719 // In any case (used V0 or not), add the V0 vertex to the cascade one.
\r
720 vCascade->AddDaughter(vV0FromCascade);
\r
723 // 5 - Add the primary track of the cascade (if any)
\r
726 // 6 - Add the cascade to the AOD array of cascades
\r
728 Double_t dcaBachToPrimVertexXY = TMath::Abs(esdCascadeBach->GetD(esd.GetPrimaryVertex()->GetX(),
\r
729 esd.GetPrimaryVertex()->GetY(),
\r
730 esd.GetMagneticField()) );
\r
732 Double_t momBachAtCascadeVtx[3]={0.};
\r
734 esdCascade->GetBPxPyPz(momBachAtCascadeVtx[0], momBachAtCascadeVtx[1], momBachAtCascadeVtx[2]);
\r
736 aodCascade = new(cascadesArray[fNumberOfCascades++]) AliAODcascade( vCascade,
\r
737 esdCascade->Charge(),
\r
738 esdCascade->GetDcaXiDaughters(),
\r
740 // DCAXiToPrimVtx -> needs to be calculated ----|
\r
741 // doesn't exist at ESD level;
\r
742 // See AODcascade::DcaXiToPrimVertex(Double, Double, Double)
\r
743 dcaBachToPrimVertexXY,
\r
744 momBachAtCascadeVtx,
\r
748 printf("---- Cascade / AOD cascade : \n\n");
\r
749 aodCascade->PrintXi(fPrimaryVertex->GetX(), fPrimaryVertex->GetY(), fPrimaryVertex->GetZ());
\r
752 } // end of the loop on cascades
\r
754 Cascades().Expand(fNumberOfCascades);
\r
757 //______________________________________________________________________________
\r
758 void AliAnalysisTaskESDfilter::ConvertV0s(const AliESDEvent& esd)
\r
760 // Access to the AOD container of V0s
\r
762 AliCodeTimerAuto("",0);
\r
768 Double_t pos[3] = { 0. };
\r
769 Double_t chi2(0.0);
\r
770 Double_t covVtx[6] = { 0. };
\r
771 Double_t momPos[3]={0.};
\r
772 Double_t covTr[21]={0.};
\r
773 Double_t pid[10]={0.};
\r
774 AliAODTrack* aodTrack(0x0);
\r
775 AliAODPid* detpid(0x0);
\r
776 Double_t momNeg[3]={0.};
\r
777 Double_t momPosAtV0vtx[3]={0.};
\r
778 Double_t momNegAtV0vtx[3]={0.};
\r
780 for (Int_t nV0 = 0; nV0 < esd.GetNumberOfV0s(); ++nV0)
\r
782 if (fUsedV0[nV0]) continue; // skip if already added to the AOD
\r
784 AliESDv0 *v0 = esd.GetV0(nV0);
\r
785 Int_t posFromV0 = v0->GetPindex();
\r
786 Int_t negFromV0 = v0->GetNindex();
\r
790 AliESDVertex *esdVtx = new AliESDVertex(*(esd.GetPrimaryVertex()));
\r
791 AliESDtrack *esdV0Pos = esd.GetTrack(posFromV0);
\r
792 AliESDtrack *esdV0Neg = esd.GetTrack(negFromV0);
\r
794 v0objects.AddAt(v0, 0);
\r
795 v0objects.AddAt(esdV0Pos, 1);
\r
796 v0objects.AddAt(esdV0Neg, 2);
\r
797 v0objects.AddAt(esdVtx, 3);
\r
798 UInt_t selectV0 = 0;
\r
800 selectV0 = fV0Filter->IsSelected(&v0objects);
\r
801 // this is a little awkward but otherwise the
\r
802 // list wants to access the pointer (delete it)
\r
803 // again when going out of scope
\r
804 delete v0objects.RemoveAt(3); // esdVtx created via copy construct
\r
810 delete v0objects.RemoveAt(3); // esdVtx created via copy construct
\r
814 v0->GetXYZ(pos[0], pos[1], pos[2]);
\r
816 if (!fOldESDformat) {
\r
817 chi2 = v0->GetChi2V0(); // = chi2/NDF since NDF = 2*2-3
\r
818 v0->GetPosCov(covVtx);
\r
821 for (Int_t i = 0; i < 6; i++) covVtx[i] = 0.;
\r
825 AliAODVertex * vV0 =
\r
826 new(Vertices()[fNumberOfVertices++]) AliAODVertex(pos,
\r
831 AliAODVertex::kV0);
\r
832 fPrimaryVertex->AddDaughter(vV0);
\r
835 // Add the positive tracks from the V0
\r
838 esdV0Pos->GetPxPyPz(momPos);
\r
839 esdV0Pos->GetXYZ(pos);
\r
840 esdV0Pos->GetCovarianceXYZPxPyPz(covTr);
\r
841 esdV0Pos->GetESDpid(pid);
\r
843 const AliESDVertex *vtx = esd.GetPrimaryVertex();
\r
845 if (!fUsedTrack[posFromV0]) {
\r
846 fUsedTrack[posFromV0] = kTRUE;
\r
847 UInt_t selectInfo = 0;
\r
848 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdV0Pos);
\r
849 if(fMChandler)fMChandler->SelectParticle(esdV0Pos->GetLabel());
\r
850 aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdV0Pos->GetID(),
\r
851 esdV0Pos->GetLabel(),
\r
857 (Short_t)esdV0Pos->GetSign(),
\r
858 esdV0Pos->GetITSClusterMap(),
\r
861 kTRUE, // check if this is right
\r
862 vtx->UsesTrack(esdV0Pos->GetID()),
\r
863 AliAODTrack::kSecondary,
\r
865 aodTrack->SetTPCFitMap(esdV0Pos->GetTPCFitMap());
\r
866 aodTrack->SetTPCClusterMap(esdV0Pos->GetTPCClusterMap());
\r
867 aodTrack->SetTPCSharedMap (esdV0Pos->GetTPCSharedMap());
\r
868 aodTrack->SetChi2perNDF(Chi2perNDF(esdV0Pos));
\r
869 aodTrack->SetTPCPointsF(esdV0Pos->GetTPCNclsF());
\r
870 fAODTrackRefs->AddAt(aodTrack,posFromV0);
\r
871 // if (fDebug > 0) printf("-------------------Bo: pos track from original pt %.3f \n",aodTrack->Pt());
\r
872 if (esdV0Pos->GetSign() > 0) ++fNumberOfPositiveTracks;
\r
873 aodTrack->ConvertAliPIDtoAODPID();
\r
874 aodTrack->SetFlags(esdV0Pos->GetStatus());
\r
875 SetAODPID(esdV0Pos,aodTrack,detpid);
\r
878 aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(posFromV0));
\r
879 // if (fDebug > 0) printf("-------------------Bo pos track from refArray pt %.3f \n",aodTrack->Pt());
\r
881 vV0->AddDaughter(aodTrack);
\r
883 // Add the negative tracks from the V0
\r
885 esdV0Neg->GetPxPyPz(momNeg);
\r
886 esdV0Neg->GetXYZ(pos);
\r
887 esdV0Neg->GetCovarianceXYZPxPyPz(covTr);
\r
888 esdV0Neg->GetESDpid(pid);
\r
890 if (!fUsedTrack[negFromV0]) {
\r
891 fUsedTrack[negFromV0] = kTRUE;
\r
892 UInt_t selectInfo = 0;
\r
893 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdV0Neg);
\r
894 if(fMChandler)fMChandler->SelectParticle(esdV0Neg->GetLabel());
\r
895 aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdV0Neg->GetID(),
\r
896 esdV0Neg->GetLabel(),
\r
902 (Short_t)esdV0Neg->GetSign(),
\r
903 esdV0Neg->GetITSClusterMap(),
\r
906 kTRUE, // check if this is right
\r
907 vtx->UsesTrack(esdV0Neg->GetID()),
\r
908 AliAODTrack::kSecondary,
\r
910 aodTrack->SetTPCFitMap(esdV0Neg->GetTPCFitMap());
\r
911 aodTrack->SetTPCClusterMap(esdV0Neg->GetTPCClusterMap());
\r
912 aodTrack->SetTPCSharedMap (esdV0Neg->GetTPCSharedMap());
\r
913 aodTrack->SetChi2perNDF(Chi2perNDF(esdV0Neg));
\r
914 aodTrack->SetTPCPointsF(esdV0Neg->GetTPCNclsF());
\r
916 fAODTrackRefs->AddAt(aodTrack,negFromV0);
\r
917 // if (fDebug > 0) printf("-------------------Bo: neg track from original pt %.3f \n",aodTrack->Pt());
\r
918 if (esdV0Neg->GetSign() > 0) ++fNumberOfPositiveTracks;
\r
919 aodTrack->ConvertAliPIDtoAODPID();
\r
920 aodTrack->SetFlags(esdV0Neg->GetStatus());
\r
921 SetAODPID(esdV0Neg,aodTrack,detpid);
\r
924 aodTrack = static_cast<AliAODTrack*>(fAODTrackRefs->At(negFromV0));
\r
925 // if (fDebug > 0) printf("-------------------Bo neg track from refArray pt %.3f \n",aodTrack->Pt());
\r
927 vV0->AddDaughter(aodTrack);
\r
930 // Add the V0 the V0 array as well
\r
932 Double_t dcaV0Daughters = v0->GetDcaV0Daughters();
\r
933 Double_t dcaV0ToPrimVertex = v0->GetD(esd.GetPrimaryVertex()->GetX(),
\r
934 esd.GetPrimaryVertex()->GetY(),
\r
935 esd.GetPrimaryVertex()->GetZ());
\r
937 v0->GetPPxPyPz(momPosAtV0vtx[0],momPosAtV0vtx[1],momPosAtV0vtx[2]);
\r
938 v0->GetNPxPyPz(momNegAtV0vtx[0],momNegAtV0vtx[1],momNegAtV0vtx[2]);
\r
940 Double_t dcaDaughterToPrimVertex[2] = { 999., 999.}; // ..[0] = DCA in (x,y) for Pos and ..[1] = Neg
\r
941 dcaDaughterToPrimVertex[0] = TMath::Abs(esdV0Pos->GetD( esd.GetPrimaryVertex()->GetX(),
\r
942 esd.GetPrimaryVertex()->GetY(),
\r
943 esd.GetMagneticField()) );
\r
944 dcaDaughterToPrimVertex[1] = TMath::Abs(esdV0Neg->GetD( esd.GetPrimaryVertex()->GetX(),
\r
945 esd.GetPrimaryVertex()->GetY(),
\r
946 esd.GetMagneticField()) );
\r
948 AliAODv0* aodV0 = new(V0s()[fNumberOfV0s++]) AliAODv0(vV0,
\r
953 dcaDaughterToPrimVertex);
\r
955 // set the aod v0 on-the-fly status
\r
956 aodV0->SetOnFlyStatus(v0->GetOnFlyStatus());
\r
957 }//End of loop on V0s
\r
959 V0s().Expand(fNumberOfV0s);
\r
962 //______________________________________________________________________________
\r
963 void AliAnalysisTaskESDfilter::ConvertTPCOnlyTracks(const AliESDEvent& esd)
\r
965 // Convert TPC only tracks
\r
966 // Here we have wo hybrid appraoch to remove fakes
\r
967 // ******* ITSTPC ********
\r
968 // Uses a cut on the ITS properties to select global tracks
\r
969 // which are than marked as HybdridITSTPC for the remainder
\r
970 // the TPC only tracks are flagged as HybridITSTPConly.
\r
971 // Note, in order not to get fakes back in the TPC cuts, one needs
\r
972 // two "ITS" cuts one tight (1) (to throw out fakes) and one lose (2) (to NOT flag the trakcs in the TPC only)
\r
973 // using cut number (3)
\r
974 // so fHybridFilterMask == (1)|(2) fTPCFilterMask = (3), Usercode needs to slect with mask = (1)|(3) and track->IsHybridITSTPC()
\r
975 // ******* TPC ********
\r
976 // Here, only TPC tracks are flagged that pass the tight ITS cuts and tracks that pass the TPC cuts and NOT the loose ITS cuts
\r
977 // the ITS cuts neeed to be added to the filter as extra cuts, since here the selections info is reset in the global and put to the TPC only track
\r
979 AliCodeTimerAuto("",0);
\r
981 // Loop over the tracks and extract and mask out all aod tracks that pass the selections for AODt racks
\r
982 for(int it = 0;it < fNumberOfTracks;++it)
\r
984 AliAODTrack *tr = (AliAODTrack*)(Tracks().At(it));
\r
986 UInt_t map = tr->GetFilterMap();
\r
987 if(map&fTPCConstrainedFilterMask){
\r
988 // we only reset the track select ionfo, no deletion...
\r
989 tr->SetFilterMap(map&~fTPCConstrainedFilterMask);
\r
991 if(map&fHybridFilterMaskTPCCG){
\r
992 // this is one part of the hybrid tracks
\r
993 // the others not passing the selection will be TPC only selected below
\r
994 tr->SetIsHybridTPCConstrainedGlobal(kTRUE);
\r
997 // Loop over the ESD trcks and pick out the tracks passing TPC only cuts
\r
1000 const AliESDVertex *vtxSPD = esd.GetPrimaryVertexSPD();
\r
1001 const AliESDVertex *vtx = esd.GetPrimaryVertex();
\r
1003 Double_t pos[3] = { 0. };
\r
1004 Double_t covTr[21]={0.};
\r
1005 Double_t pid[10]={0.};
\r
1008 Double_t p[3] = { 0. };
\r
1010 Double_t pDCA[3] = { 0. }; // momentum at DCA
\r
1011 Double_t rDCA[3] = { 0. }; // position at DCA
\r
1012 Float_t dDCA[2] = {0.}; // DCA to the vertex d and z
\r
1013 Float_t cDCA[3] = {0.}; // covariance of impact parameters
\r
1016 AliAODTrack* aodTrack(0x0);
\r
1017 // AliAODPid* detpid(0x0);
\r
1019 // account for change in pT after the constraint
\r
1020 Float_t ptMax = 1E10;
\r
1021 Float_t ptMin = 0;
\r
1022 for(int i = 0;i<32;i++){
\r
1023 if(fTPCConstrainedFilterMask&(1<<i)){
\r
1024 AliESDtrackCuts*cuts = (AliESDtrackCuts*)fTrackFilter->GetCuts()->At(i);
\r
1025 Float_t tmp1= 0,tmp2 = 0;
\r
1026 cuts->GetPtRange(tmp1,tmp2);
\r
1027 if(tmp1>ptMin)ptMin=tmp1;
\r
1028 if(tmp2<ptMax)ptMax=tmp2;
\r
1032 for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack)
\r
1034 AliESDtrack* esdTrack = esd.GetTrack(nTrack); //carefull do not modify it othwise need to work with a copy
\r
1036 UInt_t selectInfo = 0;
\r
1037 Bool_t isHybridITSTPC = false;
\r
1039 // Track selection
\r
1040 if (fTrackFilter) {
\r
1041 selectInfo = fTrackFilter->IsSelected(esdTrack);
\r
1044 if(!(selectInfo&fHybridFilterMaskTPCCG)){
\r
1045 // not already selected tracks, use second part of hybrid tracks
\r
1046 isHybridITSTPC = true;
\r
1047 // too save space one could only store these...
\r
1050 selectInfo &= fTPCConstrainedFilterMask;
\r
1051 if (!selectInfo)continue;
\r
1052 if (fWriteHybridTPCCOnly&&!isHybridITSTPC)continue; // write only complementary tracks
\r
1053 // create a tpc only tracl
\r
1054 AliESDtrack *track = AliESDtrackCuts::GetTPCOnlyTrack(const_cast<AliESDEvent*>(&esd),esdTrack->GetID());
\r
1055 if(!track) continue;
\r
1057 if(track->Pt()>0.)
\r
1059 // only constrain tracks above threshold
\r
1060 AliExternalTrackParam exParam;
\r
1061 // take the B-field from the ESD, no 3D fieldMap available at this point
\r
1062 Bool_t relate = false;
\r
1063 relate = track->RelateToVertexTPC(vtxSPD,esd.GetMagneticField(),kVeryBig,&exParam);
\r
1068 // fetch the track parameters at the DCA (unconstraint)
\r
1069 if(track->GetTPCInnerParam()){
\r
1070 track->GetTPCInnerParam()->GetPxPyPz(pDCA);
\r
1071 track->GetTPCInnerParam()->GetXYZ(rDCA);
\r
1073 // get the DCA to the vertex:
\r
1074 track->GetImpactParametersTPC(dDCA,cDCA);
\r
1075 // set the constrained parameters to the track
\r
1076 track->Set(exParam.GetX(),exParam.GetAlpha(),exParam.GetParameter(),exParam.GetCovariance());
\r
1079 track->GetPxPyPz(p);
\r
1081 Float_t pT = track->Pt();
\r
1082 if(pT<ptMin||pT>ptMax){
\r
1090 track->GetXYZ(pos);
\r
1091 track->GetCovarianceXYZPxPyPz(covTr);
\r
1092 esdTrack->GetESDpid(pid);// original PID
\r
1094 if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());
\r
1095 aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack((track->GetID()+1)*-1,
\r
1096 track->GetLabel(),
\r
1102 (Short_t)track->GetSign(),
\r
1103 track->GetITSClusterMap(),
\r
1106 kTRUE, // check if this is right
\r
1107 vtx->UsesTrack(track->GetID()),
\r
1108 AliAODTrack::kPrimary,
\r
1110 aodTrack->SetIsHybridTPCConstrainedGlobal(isHybridITSTPC);
\r
1111 aodTrack->SetTPCFitMap(track->GetTPCFitMap());
\r
1112 aodTrack->SetTPCClusterMap(track->GetTPCClusterMap());
\r
1113 aodTrack->SetTPCSharedMap (track->GetTPCSharedMap());
\r
1114 aodTrack->SetIsTPCConstrained(kTRUE);
\r
1115 aodTrack->SetChi2perNDF(Chi2perNDF(esdTrack)); // original track
\r
1116 // set the DCA values to the AOD track
\r
1117 aodTrack->SetPxPyPzAtDCA(pDCA[0],pDCA[1],pDCA[2]);
\r
1118 aodTrack->SetXYAtDCA(rDCA[0],rDCA[1]);
\r
1119 aodTrack->SetDCA(dDCA[0],dDCA[1]);
\r
1121 aodTrack->SetFlags(track->GetStatus());
\r
1122 aodTrack->SetTPCPointsF(track->GetTPCNclsF());
\r
1124 // do not duplicate PID information
\r
1125 // aodTrack->ConvertAliPIDtoAODPID();
\r
1126 // SetAODPID(esdTrack,aodTrack,detpid);
\r
1129 } // end of loop on tracks
\r
1134 void AliAnalysisTaskESDfilter::ConvertGlobalConstrainedTracks(const AliESDEvent& esd)
\r
1137 // Here we have the option to store the complement from global constraint information
\r
1138 // to tracks passing tight cuts (1) in order not to get fakes back in, one needs
\r
1139 // two sets of cuts one tight (1) (to throw out fakes) and one lose (2) (fakes/bad tracks would pass (2) but not (1))
\r
1140 // using cut number (3) selects the tracks that complement (1) e.g. tracks witout ITS refit or cluster requirement
\r
1143 AliCodeTimerAuto("",0);
\r
1145 // Loop over the tracks and extract and mask out all aod tracks that pass the selections for AODt racks
\r
1146 for(int it = 0;it < fNumberOfTracks;++it)
\r
1148 AliAODTrack *tr = (AliAODTrack*)(Tracks().At(it));
\r
1150 UInt_t map = tr->GetFilterMap();
\r
1151 if(map&fGlobalConstrainedFilterMask){
\r
1152 // we only reset the track select info, no deletion...
\r
1153 // mask reset mask in case track is already taken
\r
1154 tr->SetFilterMap(map&~fGlobalConstrainedFilterMask);
\r
1156 if(map&fHybridFilterMaskGCG){
\r
1157 // this is one part of the hybrid tracks
\r
1158 // the others not passing the selection will be the ones selected below
\r
1159 tr->SetIsHybridGlobalConstrainedGlobal(kTRUE);
\r
1162 // Loop over the ESD trcks and pick out the tracks passing the GlobalConstraint cuts
\r
1165 Double_t pos[3] = { 0. };
\r
1166 Double_t covTr[21]={0.};
\r
1167 Double_t pid[10]={0.};
\r
1168 Double_t p[3] = { 0. };
\r
1170 Double_t pDCA[3] = { 0. }; // momentum at DCA
\r
1171 Double_t rDCA[3] = { 0. }; // position at DCA
\r
1172 Float_t dDCA[2] = {0.}; // DCA to the vertex d and z
\r
1173 Float_t cDCA[3] = {0.}; // covariance of impact parameters
\r
1176 AliAODTrack* aodTrack(0x0);
\r
1177 AliAODPid* detpid(0x0);
\r
1178 const AliESDVertex *vtx = esd.GetPrimaryVertex();
\r
1180 // account for change in pT after the constraint
\r
1181 Float_t ptMax = 1E10;
\r
1182 Float_t ptMin = 0;
\r
1183 for(int i = 0;i<32;i++){
\r
1184 if(fGlobalConstrainedFilterMask&(1<<i)){
\r
1185 AliESDtrackCuts*cuts = (AliESDtrackCuts*)fTrackFilter->GetCuts()->At(i);
\r
1186 Float_t tmp1= 0,tmp2 = 0;
\r
1187 cuts->GetPtRange(tmp1,tmp2);
\r
1188 if(tmp1>ptMin)ptMin=tmp1;
\r
1189 if(tmp2<ptMax)ptMax=tmp2;
\r
1195 for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack)
\r
1197 AliESDtrack* esdTrack = esd.GetTrack(nTrack); //carefull do not modify it othwise need to work with a copy
\r
1198 const AliExternalTrackParam * exParamGC = esdTrack->GetConstrainedParam();
\r
1199 if(!exParamGC)continue;
\r
1201 UInt_t selectInfo = 0;
\r
1202 Bool_t isHybridGC = false;
\r
1205 // Track selection
\r
1206 if (fTrackFilter) {
\r
1207 selectInfo = fTrackFilter->IsSelected(esdTrack);
\r
1211 if(!(selectInfo&fHybridFilterMaskGCG))isHybridGC = true;
\r
1212 if (fWriteHybridGCOnly&&!isHybridGC)continue; // write only complementary tracks
\r
1214 selectInfo &= fGlobalConstrainedFilterMask;
\r
1215 if (!selectInfo)continue;
\r
1216 // fetch the track parameters at the DCA (unconstrained)
\r
1217 esdTrack->GetPxPyPz(pDCA);
\r
1218 esdTrack->GetXYZ(rDCA);
\r
1219 // get the DCA to the vertex:
\r
1220 esdTrack->GetImpactParameters(dDCA,cDCA);
\r
1222 esdTrack->GetConstrainedPxPyPz(p);
\r
1225 Float_t pT = exParamGC->Pt();
\r
1226 if(pT<ptMin||pT>ptMax){
\r
1231 esdTrack->GetConstrainedXYZ(pos);
\r
1232 exParamGC->GetCovarianceXYZPxPyPz(covTr);
\r
1233 esdTrack->GetESDpid(pid);
\r
1234 if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());
\r
1235 aodTrack = new(Tracks()[fNumberOfTracks++]) AliAODTrack((esdTrack->GetID()+1)*-1,
\r
1236 esdTrack->GetLabel(),
\r
1242 (Short_t)esdTrack->GetSign(),
\r
1243 esdTrack->GetITSClusterMap(),
\r
1246 kTRUE, // check if this is right
\r
1247 vtx->UsesTrack(esdTrack->GetID()),
\r
1248 AliAODTrack::kPrimary,
\r
1250 aodTrack->SetIsHybridGlobalConstrainedGlobal(isHybridGC);
\r
1251 aodTrack->SetIsGlobalConstrained(kTRUE);
\r
1252 aodTrack->SetTPCFitMap(esdTrack->GetTPCFitMap());
\r
1253 aodTrack->SetTPCClusterMap(esdTrack->GetTPCClusterMap());
\r
1254 aodTrack->SetTPCSharedMap (esdTrack->GetTPCSharedMap());
\r
1255 aodTrack->SetChi2perNDF(Chi2perNDF(esdTrack));
\r
1258 // set the DCA values to the AOD track
\r
1259 aodTrack->SetPxPyPzAtDCA(pDCA[0],pDCA[1],pDCA[2]);
\r
1260 aodTrack->SetXYAtDCA(rDCA[0],rDCA[1]);
\r
1261 aodTrack->SetDCA(dDCA[0],dDCA[1]);
\r
1263 aodTrack->SetFlags(esdTrack->GetStatus());
\r
1264 aodTrack->SetTPCPointsF(esdTrack->GetTPCNclsF());
\r
1267 // only copy AOD information for hybrid, no duplicate information
\r
1268 aodTrack->ConvertAliPIDtoAODPID();
\r
1269 SetAODPID(esdTrack,aodTrack,detpid);
\r
1271 } // end of loop on tracks
\r
1276 //______________________________________________________________________________
\r
1277 void AliAnalysisTaskESDfilter::ConvertTracks(const AliESDEvent& esd)
\r
1279 // Tracks (primary and orphan)
\r
1281 AliCodeTimerAuto("",0);
\r
1283 AliDebug(1,Form("NUMBER OF ESD TRACKS %5d\n", esd.GetNumberOfTracks()));
\r
1285 const AliESDVertex *vtx = esd.GetPrimaryVertex();
\r
1286 Double_t p[3] = { 0. };
\r
1287 Double_t pos[3] = { 0. };
\r
1288 Double_t covTr[21] = { 0. };
\r
1289 Double_t pid[10] = { 0. };
\r
1290 AliAODTrack* aodTrack(0x0);
\r
1291 AliAODPid* detpid(0x0);
\r
1293 for (Int_t nTrack = 0; nTrack < esd.GetNumberOfTracks(); ++nTrack)
\r
1295 if (fUsedTrack[nTrack]) continue;
\r
1297 AliESDtrack *esdTrack = esd.GetTrack(nTrack);
\r
1298 UInt_t selectInfo = 0;
\r
1300 // Track selection
\r
1301 if (fTrackFilter) {
\r
1302 selectInfo = fTrackFilter->IsSelected(esdTrack);
\r
1303 if (!selectInfo && !vtx->UsesTrack(esdTrack->GetID())) continue;
\r
1307 esdTrack->GetPxPyPz(p);
\r
1308 esdTrack->GetXYZ(pos);
\r
1309 esdTrack->GetCovarianceXYZPxPyPz(covTr);
\r
1310 esdTrack->GetESDpid(pid);
\r
1311 if(fMChandler)fMChandler->SelectParticle(esdTrack->GetLabel());
\r
1312 fPrimaryVertex->AddDaughter(aodTrack =
\r
1313 new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrack->GetID(),
\r
1314 esdTrack->GetLabel(),
\r
1320 (Short_t)esdTrack->GetSign(),
\r
1321 esdTrack->GetITSClusterMap(),
\r
1324 kTRUE, // check if this is right
\r
1325 vtx->UsesTrack(esdTrack->GetID()),
\r
1326 AliAODTrack::kPrimary,
\r
1329 aodTrack->SetTPCFitMap(esdTrack->GetTPCFitMap());
\r
1330 aodTrack->SetTPCClusterMap(esdTrack->GetTPCClusterMap());
\r
1331 aodTrack->SetTPCSharedMap (esdTrack->GetTPCSharedMap());
\r
1332 aodTrack->SetChi2perNDF(Chi2perNDF(esdTrack));
\r
1333 aodTrack->SetTPCPointsF(esdTrack->GetTPCNclsF());
\r
1334 if(esdTrack->IsEMCAL()) aodTrack->SetEMCALcluster(esdTrack->GetEMCALcluster());
\r
1335 if(esdTrack->IsPHOS()) aodTrack->SetPHOScluster(esdTrack->GetPHOScluster());
\r
1337 fAODTrackRefs->AddAt(aodTrack, nTrack);
\r
1340 if (esdTrack->GetSign() > 0) ++fNumberOfPositiveTracks;
\r
1341 aodTrack->SetFlags(esdTrack->GetStatus());
\r
1342 aodTrack->ConvertAliPIDtoAODPID();
\r
1343 SetAODPID(esdTrack,aodTrack,detpid);
\r
1344 } // end of loop on tracks
\r
1347 //______________________________________________________________________________
\r
1348 void AliAnalysisTaskESDfilter::ConvertPmdClusters(const AliESDEvent& esd)
\r
1350 // Convert PMD Clusters
\r
1351 AliCodeTimerAuto("",0);
\r
1352 Int_t jPmdClusters=0;
\r
1353 // Access to the AOD container of PMD clusters
\r
1354 TClonesArray &pmdClusters = *(AODEvent()->GetPmdClusters());
\r
1355 for (Int_t iPmd = 0; iPmd < esd.GetNumberOfPmdTracks(); ++iPmd) {
\r
1356 // file pmd clusters, to be revised!
\r
1357 AliESDPmdTrack *pmdTrack = esd.GetPmdTrack(iPmd);
\r
1359 Int_t *label = 0x0;
\r
1360 Double_t posPmd[3] = { pmdTrack->GetClusterX(), pmdTrack->GetClusterY(), pmdTrack->GetClusterZ()};
\r
1361 Double_t pidPmd[13] = { 0.}; // to be revised!
\r
1363 // assoc cluster not set
\r
1364 new(pmdClusters[jPmdClusters++]) AliAODPmdCluster(iPmd, nLabel, label, pmdTrack->GetClusterADC(), posPmd, pidPmd);
\r
1368 //______________________________________________________________________________
\r
1369 void AliAnalysisTaskESDfilter::ConvertCaloClusters(const AliESDEvent& esd)
\r
1371 // Convert Calorimeter Clusters
\r
1372 AliCodeTimerAuto("",0);
\r
1374 // Access to the AOD container of clusters
\r
1375 TClonesArray &caloClusters = *(AODEvent()->GetCaloClusters());
\r
1376 Int_t jClusters(0);
\r
1378 for (Int_t iClust=0; iClust<esd.GetNumberOfCaloClusters(); ++iClust) {
\r
1380 AliESDCaloCluster * cluster = esd.GetCaloCluster(iClust);
\r
1382 Int_t id = cluster->GetID();
\r
1383 Int_t nLabel = cluster->GetNLabels();
\r
1384 Int_t *labels = cluster->GetLabels();
\r
1386 for(int i = 0;i < nLabel;++i){
\r
1387 if(fMChandler)fMChandler->SelectParticle(labels[i]);
\r
1391 Float_t energy = cluster->E();
\r
1392 Float_t posF[3] = { 0.};
\r
1393 cluster->GetPosition(posF);
\r
1395 AliAODCaloCluster *caloCluster = new(caloClusters[jClusters++]) AliAODCaloCluster(id,
\r
1401 cluster->GetType(),0);
\r
1403 caloCluster->SetCaloCluster(cluster->GetDistanceToBadChannel(),
\r
1404 cluster->GetDispersion(),
\r
1405 cluster->GetM20(), cluster->GetM02(),
\r
1406 cluster->GetEmcCpvDistance(),
\r
1407 cluster->GetNExMax(),cluster->GetTOF()) ;
\r
1409 caloCluster->SetPIDFromESD(cluster->GetPID());
\r
1410 caloCluster->SetNCells(cluster->GetNCells());
\r
1411 caloCluster->SetCellsAbsId(cluster->GetCellsAbsId());
\r
1412 caloCluster->SetCellsAmplitudeFraction(cluster->GetCellsAmplitudeFraction());
\r
1414 caloCluster->SetTrackDistance(cluster->GetTrackDx(), cluster->GetTrackDz());
\r
1416 TArrayI* matchedT = cluster->GetTracksMatched();
\r
1417 if (fNumberOfTracks>0 && matchedT && cluster->GetTrackMatchedIndex() >= 0) {
\r
1418 for (Int_t im = 0; im < matchedT->GetSize(); im++) {
\r
1419 Int_t iESDtrack = matchedT->At(im);;
\r
1420 if (fAODTrackRefs->At(iESDtrack) != 0) {
\r
1421 caloCluster->AddTrackMatched((AliAODTrack*)fAODTrackRefs->At(iESDtrack));
\r
1427 caloClusters.Expand(jClusters); // resize TObjArray to 'remove' slots for pseudo clusters
\r
1430 //______________________________________________________________________________
\r
1431 void AliAnalysisTaskESDfilter::ConvertEMCALCells(const AliESDEvent& esd)
\r
1433 // Convert EMCAL Cells
\r
1434 AliCodeTimerAuto("",0);
\r
1435 // fill EMCAL cell info
\r
1436 if (esd.GetEMCALCells()) { // protection against missing ESD information
\r
1437 AliESDCaloCells &esdEMcells = *(esd.GetEMCALCells());
\r
1438 Int_t nEMcell = esdEMcells.GetNumberOfCells() ;
\r
1440 AliAODCaloCells &aodEMcells = *(AODEvent()->GetEMCALCells());
\r
1441 aodEMcells.CreateContainer(nEMcell);
\r
1442 aodEMcells.SetType(AliAODCaloCells::kEMCALCell);
\r
1443 for (Int_t iCell = 0; iCell < nEMcell; iCell++) {
\r
1444 aodEMcells.SetCell(iCell,esdEMcells.GetCellNumber(iCell),esdEMcells.GetAmplitude(iCell));
\r
1446 aodEMcells.Sort();
\r
1450 //______________________________________________________________________________
\r
1451 void AliAnalysisTaskESDfilter::ConvertPHOSCells(const AliESDEvent& esd)
\r
1453 // Convert PHOS Cells
\r
1454 AliCodeTimerAuto("",0);
\r
1455 // fill PHOS cell info
\r
1456 if (esd.GetPHOSCells()) { // protection against missing ESD information
\r
1457 AliESDCaloCells &esdPHcells = *(esd.GetPHOSCells());
\r
1458 Int_t nPHcell = esdPHcells.GetNumberOfCells() ;
\r
1460 AliAODCaloCells &aodPHcells = *(AODEvent()->GetPHOSCells());
\r
1461 aodPHcells.CreateContainer(nPHcell);
\r
1462 aodPHcells.SetType(AliAODCaloCells::kPHOSCell);
\r
1463 for (Int_t iCell = 0; iCell < nPHcell; iCell++) {
\r
1464 aodPHcells.SetCell(iCell,esdPHcells.GetCellNumber(iCell),esdPHcells.GetAmplitude(iCell));
\r
1466 aodPHcells.Sort();
\r
1470 //______________________________________________________________________________
\r
1471 void AliAnalysisTaskESDfilter::ConvertTracklets(const AliESDEvent& esd)
\r
1474 AliCodeTimerAuto("",0);
\r
1476 AliAODTracklets &SPDTracklets = *(AODEvent()->GetTracklets());
\r
1477 const AliMultiplicity *mult = esd.GetMultiplicity();
\r
1479 if (mult->GetNumberOfTracklets()>0) {
\r
1480 SPDTracklets.CreateContainer(mult->GetNumberOfTracklets());
\r
1482 for (Int_t n=0; n<mult->GetNumberOfTracklets(); n++) {
\r
1484 fMChandler->SelectParticle(mult->GetLabel(n, 0));
\r
1485 fMChandler->SelectParticle(mult->GetLabel(n, 1));
\r
1487 SPDTracklets.SetTracklet(n, mult->GetTheta(n), mult->GetPhi(n), mult->GetDeltaPhi(n), mult->GetLabel(n, 0),mult->GetLabel(n, 1));
\r
1491 //Printf("ERROR: AliMultiplicity could not be retrieved from ESD");
\r
1495 //______________________________________________________________________________
\r
1496 void AliAnalysisTaskESDfilter::ConvertKinks(const AliESDEvent& esd)
\r
1498 AliCodeTimerAuto("",0);
\r
1500 // Kinks: it is a big mess the access to the information in the kinks
\r
1501 // The loop is on the tracks in order to find the mother and daugther of each kink
\r
1503 Double_t covTr[21]={0.};
\r
1504 Double_t pid[10]={0.};
\r
1505 AliAODPid* detpid(0x0);
\r
1507 fNumberOfKinks = esd.GetNumberOfKinks();
\r
1509 const AliESDVertex* vtx = esd.GetPrimaryVertex();
\r
1511 for (Int_t iTrack=0; iTrack<esd.GetNumberOfTracks(); ++iTrack)
\r
1513 AliESDtrack * esdTrack = esd.GetTrack(iTrack);
\r
1515 Int_t ikink = esdTrack->GetKinkIndex(0);
\r
1517 if (ikink && fNumberOfKinks) {
\r
1518 // Negative kink index: mother, positive: daughter
\r
1520 // Search for the second track of the kink
\r
1522 for (Int_t jTrack = iTrack+1; jTrack<esd.GetNumberOfTracks(); ++jTrack) {
\r
1524 AliESDtrack * esdTrack1 = esd.GetTrack(jTrack);
\r
1526 Int_t jkink = esdTrack1->GetKinkIndex(0);
\r
1528 if ( TMath::Abs(ikink)==TMath::Abs(jkink) ) {
\r
1530 // The two tracks are from the same kink
\r
1532 if (fUsedKink[TMath::Abs(ikink)-1]) continue; // skip used kinks
\r
1534 Int_t imother = -1;
\r
1535 Int_t idaughter = -1;
\r
1537 if (ikink<0 && jkink>0) {
\r
1540 idaughter = jTrack;
\r
1542 else if (ikink>0 && jkink<0) {
\r
1545 idaughter = iTrack;
\r
1548 // cerr << "Error: Wrong combination of kink indexes: "
\r
1549 // << ikink << " " << jkink << endl;
\r
1553 // Add the mother track if it passed primary track selection cuts
\r
1555 AliAODTrack * mother = NULL;
\r
1557 UInt_t selectInfo = 0;
\r
1558 if (fTrackFilter) {
\r
1559 selectInfo = fTrackFilter->IsSelected(esd.GetTrack(imother));
\r
1560 if (!selectInfo) continue;
\r
1563 if (!fUsedTrack[imother]) {
\r
1565 fUsedTrack[imother] = kTRUE;
\r
1567 AliESDtrack *esdTrackM = esd.GetTrack(imother);
\r
1568 Double_t p[3] = { 0. };
\r
1569 Double_t pos[3] = { 0. };
\r
1570 esdTrackM->GetPxPyPz(p);
\r
1571 esdTrackM->GetXYZ(pos);
\r
1572 esdTrackM->GetCovarianceXYZPxPyPz(covTr);
\r
1573 esdTrackM->GetESDpid(pid);
\r
1574 if(fMChandler)fMChandler->SelectParticle(esdTrackM->GetLabel());
\r
1576 new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrackM->GetID(),
\r
1577 esdTrackM->GetLabel(),
\r
1583 (Short_t)esdTrackM->GetSign(),
\r
1584 esdTrackM->GetITSClusterMap(),
\r
1587 kTRUE, // check if this is right
\r
1588 vtx->UsesTrack(esdTrack->GetID()),
\r
1589 AliAODTrack::kPrimary,
\r
1591 mother->SetTPCFitMap(esdTrackM->GetTPCFitMap());
\r
1592 mother->SetTPCClusterMap(esdTrackM->GetTPCClusterMap());
\r
1593 mother->SetTPCSharedMap (esdTrackM->GetTPCSharedMap());
\r
1594 mother->SetChi2perNDF(Chi2perNDF(esdTrackM));
\r
1595 mother->SetTPCPointsF(esdTrackM->GetTPCNclsF());
\r
1597 fAODTrackRefs->AddAt(mother, imother);
\r
1599 if (esdTrackM->GetSign() > 0) ++fNumberOfPositiveTracks;
\r
1600 mother->SetFlags(esdTrackM->GetStatus());
\r
1601 mother->ConvertAliPIDtoAODPID();
\r
1602 fPrimaryVertex->AddDaughter(mother);
\r
1603 mother->ConvertAliPIDtoAODPID();
\r
1604 SetAODPID(esdTrackM,mother,detpid);
\r
1607 // cerr << "Error: event " << esd.GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1
\r
1608 // << " track " << imother << " has already been used!" << endl;
\r
1611 // Add the kink vertex
\r
1612 AliESDkink * kink = esd.GetKink(TMath::Abs(ikink)-1);
\r
1614 AliAODVertex * vkink =
\r
1615 new(Vertices()[fNumberOfVertices++]) AliAODVertex(kink->GetPosition(),
\r
1619 esdTrack->GetID(), // This is the track ID of the mother's track!
\r
1620 AliAODVertex::kKink);
\r
1621 // Add the daughter track
\r
1623 AliAODTrack * daughter = NULL;
\r
1625 if (!fUsedTrack[idaughter]) {
\r
1627 fUsedTrack[idaughter] = kTRUE;
\r
1629 AliESDtrack *esdTrackD = esd.GetTrack(idaughter);
\r
1630 Double_t p[3] = { 0. };
\r
1631 Double_t pos[3] = { 0. };
\r
1633 esdTrackD->GetPxPyPz(p);
\r
1634 esdTrackD->GetXYZ(pos);
\r
1635 esdTrackD->GetCovarianceXYZPxPyPz(covTr);
\r
1636 esdTrackD->GetESDpid(pid);
\r
1638 if (fTrackFilter) selectInfo = fTrackFilter->IsSelected(esdTrackD);
\r
1639 if(fMChandler)fMChandler->SelectParticle(esdTrackD->GetLabel());
\r
1641 new(Tracks()[fNumberOfTracks++]) AliAODTrack(esdTrackD->GetID(),
\r
1642 esdTrackD->GetLabel(),
\r
1648 (Short_t)esdTrackD->GetSign(),
\r
1649 esdTrackD->GetITSClusterMap(),
\r
1652 kTRUE, // check if this is right
\r
1653 vtx->UsesTrack(esdTrack->GetID()),
\r
1654 AliAODTrack::kSecondary,
\r
1656 daughter->SetTPCFitMap(esdTrackD->GetTPCFitMap());
\r
1657 daughter->SetTPCClusterMap(esdTrackD->GetTPCClusterMap());
\r
1658 daughter->SetTPCSharedMap (esdTrackD->GetTPCSharedMap());
\r
1659 daughter->SetTPCPointsF(esdTrackD->GetTPCNclsF());
\r
1660 fAODTrackRefs->AddAt(daughter, idaughter);
\r
1662 if (esdTrackD->GetSign() > 0) ++fNumberOfPositiveTracks;
\r
1663 daughter->SetFlags(esdTrackD->GetStatus());
\r
1664 daughter->ConvertAliPIDtoAODPID();
\r
1665 vkink->AddDaughter(daughter);
\r
1666 daughter->ConvertAliPIDtoAODPID();
\r
1667 SetAODPID(esdTrackD,daughter,detpid);
\r
1670 // cerr << "Error: event " << esd.GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1
\r
1671 // << " track " << idaughter << " has already been used!" << endl;
\r
1679 //______________________________________________________________________________
\r
1680 void AliAnalysisTaskESDfilter::ConvertPrimaryVertices(const AliESDEvent& esd)
\r
1682 AliCodeTimerAuto("",0);
\r
1684 // Access to the AOD container of vertices
\r
1685 fNumberOfVertices = 0;
\r
1687 Double_t pos[3] = { 0. };
\r
1688 Double_t covVtx[6] = { 0. };
\r
1690 // Add primary vertex. The primary tracks will be defined
\r
1691 // after the loops on the composite objects (V0, cascades, kinks)
\r
1692 const AliESDVertex *vtx = esd.GetPrimaryVertex();
\r
1694 vtx->GetXYZ(pos); // position
\r
1695 vtx->GetCovMatrix(covVtx); //covariance matrix
\r
1697 fPrimaryVertex = new(Vertices()[fNumberOfVertices++])
\r
1698 AliAODVertex(pos, covVtx, vtx->GetChi2toNDF(), NULL, -1, AliAODVertex::kPrimary);
\r
1699 fPrimaryVertex->SetName(vtx->GetName());
\r
1700 fPrimaryVertex->SetTitle(vtx->GetTitle());
\r
1702 TString vtitle = vtx->GetTitle();
\r
1703 if (!vtitle.Contains("VertexerTracks"))
\r
1704 fPrimaryVertex->SetNContributors(vtx->GetNContributors());
\r
1706 if (fDebug > 0) fPrimaryVertex->Print();
\r
1708 // Add SPD "main" vertex
\r
1709 const AliESDVertex *vtxS = esd.GetPrimaryVertexSPD();
\r
1710 vtxS->GetXYZ(pos); // position
\r
1711 vtxS->GetCovMatrix(covVtx); //covariance matrix
\r
1712 AliAODVertex * mVSPD = new(Vertices()[fNumberOfVertices++])
\r
1713 AliAODVertex(pos, covVtx, vtxS->GetChi2toNDF(), NULL, -1, AliAODVertex::kMainSPD);
\r
1714 mVSPD->SetName(vtxS->GetName());
\r
1715 mVSPD->SetTitle(vtxS->GetTitle());
\r
1716 mVSPD->SetNContributors(vtxS->GetNContributors());
\r
1718 // Add SPD pileup vertices
\r
1719 for(Int_t iV=0; iV<esd.GetNumberOfPileupVerticesSPD(); ++iV)
\r
1721 const AliESDVertex *vtxP = esd.GetPileupVertexSPD(iV);
\r
1722 vtxP->GetXYZ(pos); // position
\r
1723 vtxP->GetCovMatrix(covVtx); //covariance matrix
\r
1724 AliAODVertex * pVSPD = new(Vertices()[fNumberOfVertices++])
\r
1725 AliAODVertex(pos, covVtx, vtxP->GetChi2toNDF(), NULL, -1, AliAODVertex::kPileupSPD);
\r
1726 pVSPD->SetName(vtxP->GetName());
\r
1727 pVSPD->SetTitle(vtxP->GetTitle());
\r
1728 pVSPD->SetNContributors(vtxP->GetNContributors());
\r
1729 pVSPD->SetBC(vtxP->GetBC());
\r
1732 // Add TRK pileup vertices
\r
1733 for(Int_t iV=0; iV<esd.GetNumberOfPileupVerticesTracks(); ++iV)
\r
1735 const AliESDVertex *vtxP = esd.GetPileupVertexTracks(iV);
\r
1736 vtxP->GetXYZ(pos); // position
\r
1737 vtxP->GetCovMatrix(covVtx); //covariance matrix
\r
1738 AliAODVertex * pVTRK = new(Vertices()[fNumberOfVertices++])
\r
1739 AliAODVertex(pos, covVtx, vtxP->GetChi2toNDF(), NULL, -1, AliAODVertex::kPileupTracks);
\r
1740 pVTRK->SetName(vtxP->GetName());
\r
1741 pVTRK->SetTitle(vtxP->GetTitle());
\r
1742 pVTRK->SetNContributors(vtxP->GetNContributors());
\r
1743 pVTRK->SetBC(vtxP->GetBC());
\r
1747 //______________________________________________________________________________
\r
1748 void AliAnalysisTaskESDfilter::ConvertVZERO(const AliESDEvent& esd)
\r
1750 // Convert VZERO data
\r
1751 AliAODVZERO* vzeroData = AODEvent()->GetVZEROData();
\r
1752 *vzeroData = *(esd.GetVZEROData());
\r
1755 //______________________________________________________________________________
\r
1756 void AliAnalysisTaskESDfilter::ConvertTZERO(const AliESDEvent& esd)
\r
1758 // Convert TZERO data
\r
1759 const AliESDTZERO* esdTzero = esd.GetESDTZERO();
\r
1760 AliAODTZERO* aodTzero = AODEvent()->GetTZEROData();
\r
1762 for (Int_t icase=0; icase<3; icase++){
\r
1763 aodTzero->SetT0TOF( icase, esdTzero->GetT0TOF(icase));
\r
1764 aodTzero->SetT0TOFbest(icase, esdTzero->GetT0TOFbest(icase));
\r
1766 aodTzero->SetBackgroundFlag(esdTzero->GetBackgroundFlag());
\r
1767 aodTzero->SetPileupFlag(esdTzero->GetPileupFlag());
\r
1768 aodTzero->SetSatelliteFlag(esdTzero->GetSatellite());
\r
1773 //______________________________________________________________________________
\r
1774 void AliAnalysisTaskESDfilter::ConvertZDC(const AliESDEvent& esd)
\r
1776 // Convert ZDC data
\r
1777 AliESDZDC* esdZDC = esd.GetZDCData();
\r
1779 const Double_t zem1Energy = esdZDC->GetZEM1Energy();
\r
1780 const Double_t zem2Energy = esdZDC->GetZEM2Energy();
\r
1782 const Double_t *towZNC = esdZDC->GetZNCTowerEnergy();
\r
1783 const Double_t *towZPC = esdZDC->GetZPCTowerEnergy();
\r
1784 const Double_t *towZNA = esdZDC->GetZNATowerEnergy();
\r
1785 const Double_t *towZPA = esdZDC->GetZPATowerEnergy();
\r
1786 const Double_t *towZNCLG = esdZDC->GetZNCTowerEnergyLR();
\r
1787 const Double_t *towZNALG = esdZDC->GetZNATowerEnergyLR();
\r
1789 AliAODZDC* zdcAOD = AODEvent()->GetZDCData();
\r
1791 zdcAOD->SetZEM1Energy(zem1Energy);
\r
1792 zdcAOD->SetZEM2Energy(zem2Energy);
\r
1793 zdcAOD->SetZNCTowers(towZNC, towZNCLG);
\r
1794 zdcAOD->SetZNATowers(towZNA, towZNALG);
\r
1795 zdcAOD->SetZPCTowers(towZPC);
\r
1796 zdcAOD->SetZPATowers(towZPA);
\r
1798 zdcAOD->SetZDCParticipants(esdZDC->GetZDCParticipants(), esdZDC->GetZDCPartSideA(), esdZDC->GetZDCPartSideC());
\r
1799 zdcAOD->SetZDCImpactParameter(esdZDC->GetImpactParameter(), esdZDC->GetImpactParamSideA(),
\r
1800 esdZDC->GetImpactParamSideC());
\r
1801 zdcAOD->SetZDCTDCSum(esdZDC->GetZNTDCSum(0));
\r
1802 zdcAOD->SetZDCTDCDiff(esdZDC->GetZNTDCDiff(0));
\r
1806 //______________________________________________________________________________
\r
1807 void AliAnalysisTaskESDfilter::ConvertESDtoAOD()
\r
1809 // ESD Filter analysis task executed for each event
\r
1811 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
\r
1815 AliCodeTimerAuto("",0);
\r
1817 fOldESDformat = ( esd->GetAliESDOld() != 0x0 );
\r
1819 fNumberOfTracks = 0;
\r
1820 fNumberOfPositiveTracks = 0;
\r
1822 fNumberOfVertices = 0;
\r
1823 fNumberOfCascades = 0;
\r
1824 fNumberOfKinks = 0;
\r
1826 AliAODHeader* header = ConvertHeader(*esd);
\r
1828 if ( fIsVZEROEnabled ) ConvertVZERO(*esd);
\r
1829 if ( fIsTZEROEnabled ) ConvertTZERO(*esd);
\r
1831 // Fetch Stack for debuggging if available
\r
1835 fMChandler = (AliMCEventHandler*) ((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler());
\r
1838 // loop over events and fill them
\r
1839 // Multiplicity information needed by the header (to be revised!)
\r
1840 Int_t nTracks = esd->GetNumberOfTracks();
\r
1841 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) esd->GetTrack(iTrack)->SetESDEvent(esd);
\r
1843 // Update the header
\r
1845 Int_t nV0s = esd->GetNumberOfV0s();
\r
1846 Int_t nCascades = esd->GetNumberOfCascades();
\r
1847 Int_t nKinks = esd->GetNumberOfKinks();
\r
1848 Int_t nVertices = nV0s + nCascades /*V0 wihtin cascade already counted*/+ nKinks + 1 /* = prim. vtx*/;
\r
1849 Int_t nPileSPDVertices=1+esd->GetNumberOfPileupVerticesSPD(); // also SPD main vertex
\r
1850 Int_t nPileTrkVertices=esd->GetNumberOfPileupVerticesTracks();
\r
1851 nVertices+=nPileSPDVertices;
\r
1852 nVertices+=nPileTrkVertices;
\r
1854 Int_t nCaloClus = esd->GetNumberOfCaloClusters();
\r
1855 Int_t nFmdClus = 0;
\r
1856 Int_t nPmdClus = esd->GetNumberOfPmdTracks();
\r
1858 AliDebug(1,Form(" NV0=%d NCASCADES=%d NKINKS=%d", nV0s, nCascades, nKinks));
\r
1860 AODEvent()->ResetStd(nTracks, nVertices, nV0s, nCascades, nJets, nCaloClus, nFmdClus, nPmdClus);
\r
1864 // RefArray to store a mapping between esd V0 number and newly created AOD-Vertex V0
\r
1865 fAODV0VtxRefs = new TRefArray(nV0s);
\r
1866 // RefArray to store the mapping between esd V0 number and newly created AOD-V0
\r
1867 fAODV0Refs = new TRefArray(nV0s);
\r
1868 // Array to take into account the V0s already added to the AOD (V0 within cascades)
\r
1869 fUsedV0 = new Bool_t[nV0s];
\r
1870 for (Int_t iV0=0; iV0<nV0s; ++iV0) fUsedV0[iV0]=kFALSE;
\r
1875 // RefArray to store the mapping between esd track number and newly created AOD-Track
\r
1877 fAODTrackRefs = new TRefArray(nTracks);
\r
1879 // Array to take into account the tracks already added to the AOD
\r
1880 fUsedTrack = new Bool_t[nTracks];
\r
1881 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) fUsedTrack[iTrack]=kFALSE;
\r
1884 // Array to take into account the kinks already added to the AOD
\r
1887 fUsedKink = new Bool_t[nKinks];
\r
1888 for (Int_t iKink=0; iKink<nKinks; ++iKink) fUsedKink[iKink]=kFALSE;
\r
1891 ConvertPrimaryVertices(*esd);
\r
1893 //setting best TOF PID
\r
1894 AliESDInputHandler* esdH = dynamic_cast<AliESDInputHandler*>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
\r
1896 fESDpid = esdH->GetESDpid();
\r
1898 if (fIsPidOwner && fESDpid){
\r
1903 { //in case of no Tender attached
\r
1904 fESDpid = new AliESDpid;
\r
1905 fIsPidOwner = kTRUE;
\r
1908 if(!esd->GetTOFHeader())
\r
1909 { //protection in case the pass2 LHC10b,c,d have been processed without tender.
\r
1910 Float_t t0spread[10];
\r
1911 Float_t intrinsicTOFres=100; //ps ok for LHC10b,c,d pass2!!
\r
1912 for (Int_t i=0; i<10; i++) t0spread[i] = (TMath::Sqrt(esd->GetSigma2DiamondZ()))/0.03; //0.03 to convert from cm to ps
\r
1913 fESDpid->GetTOFResponse().SetT0resolution(t0spread);
\r
1914 fESDpid->GetTOFResponse().SetTimeResolution(intrinsicTOFres);
\r
1916 fESDpid->SetTOFResponse(esd, (AliESDpid::EStartTimeType_t)fTimeZeroType);
\r
1919 if(esd->GetTOFHeader() && fIsPidOwner) fESDpid->SetTOFResponse(esd, (AliESDpid::EStartTimeType_t)fTimeZeroType); //in case of AOD production strating form LHC10e without Tender.
\r
1921 if ( fAreCascadesEnabled ) ConvertCascades(*esd);
\r
1923 if ( fAreV0sEnabled ) ConvertV0s(*esd);
\r
1925 if ( fAreKinksEnabled ) ConvertKinks(*esd);
\r
1927 if ( fAreTracksEnabled ) ConvertTracks(*esd);
\r
1929 // Update number of AOD tracks in header at the end of track loop (M.G.)
\r
1930 header->SetRefMultiplicity(fNumberOfTracks);
\r
1931 header->SetRefMultiplicityPos(fNumberOfPositiveTracks);
\r
1932 header->SetRefMultiplicityNeg(fNumberOfTracks - fNumberOfPositiveTracks);
\r
1934 if ( fTPCConstrainedFilterMask ) ConvertTPCOnlyTracks(*esd);
\r
1935 if( fGlobalConstrainedFilterMask) ConvertGlobalConstrainedTracks(*esd);
\r
1937 if ( fArePmdClustersEnabled ) ConvertPmdClusters(*esd);
\r
1939 if ( fAreCaloClustersEnabled ) ConvertCaloClusters(*esd);
\r
1941 if ( fAreEMCALCellsEnabled )ConvertEMCALCells(*esd);
\r
1943 if ( fArePHOSCellsEnabled )ConvertPHOSCells(*esd);
\r
1945 if ( fAreTrackletsEnabled ) ConvertTracklets(*esd);
\r
1947 delete fAODTrackRefs; fAODTrackRefs=0x0;
\r
1948 delete fAODV0VtxRefs; fAODV0VtxRefs=0x0;
\r
1949 delete fAODV0Refs; fAODV0Refs=0x0;
\r
1951 delete[] fUsedTrack; fUsedTrack=0x0;
\r
1952 delete[] fUsedV0; fUsedV0=0x0;
\r
1953 delete[] fUsedKink; fUsedKink=0x0;
\r
1955 if ( fIsPidOwner){
\r
1964 //______________________________________________________________________________
\r
1965 void AliAnalysisTaskESDfilter::SetAODPID(AliESDtrack *esdtrack, AliAODTrack *aodtrack, AliAODPid *detpid)
\r
1968 // Setter for the raw PID detector signals
\r
1971 // Save PID object for candidate electrons
\r
1972 Bool_t pidSave = kFALSE;
\r
1973 if (fTrackFilter) {
\r
1974 Bool_t selectInfo = fTrackFilter->IsSelected((char*) "Electrons");
\r
1975 if (selectInfo) pidSave = kTRUE;
\r
1979 // Tracks passing pt cut
\r
1980 if(esdtrack->Pt()>fHighPthreshold) {
\r
1984 if(esdtrack->Pt()> fPtshape->GetXmin()){
\r
1985 Double_t y = fPtshape->Eval(esdtrack->Pt())/fPtshape->Eval(fHighPthreshold);
\r
1986 if(gRandom->Rndm(0)<1./y){
\r
1989 }//end if p < pmin
\r
1990 }//end if p function
\r
1994 if(!aodtrack->GetDetPid()){// prevent memory leak when calling SetAODPID twice for the same track
\r
1995 detpid = new AliAODPid();
\r
1996 SetDetectorRawSignals(detpid,esdtrack);
\r
1997 aodtrack->SetDetPID(detpid);
\r
2002 //______________________________________________________________________________
\r
2003 void AliAnalysisTaskESDfilter::SetDetectorRawSignals(AliAODPid *aodpid, AliESDtrack *track)
\r
2006 //assignment of the detector signals (AliXXXesdPID inspired)
\r
2009 AliInfo("no ESD track found. .....exiting");
\r
2013 const AliExternalTrackParam *in=track->GetInnerParam();
\r
2015 aodpid->SetTPCmomentum(in->GetP());
\r
2017 aodpid->SetTPCmomentum(-1.);
\r
2021 aodpid->SetITSsignal(track->GetITSsignal());
\r
2022 Double_t itsdedx[4]; // dE/dx samples for individual ITS layers
\r
2023 track->GetITSdEdxSamples(itsdedx);
\r
2024 aodpid->SetITSdEdxSamples(itsdedx);
\r
2026 aodpid->SetTPCsignal(track->GetTPCsignal());
\r
2027 aodpid->SetTPCsignalN(track->GetTPCsignalN());
\r
2029 //n TRD planes = 6
\r
2030 Int_t nslices = track->GetNumberOfTRDslices()*6;
\r
2031 TArrayD trdslices(nslices);
\r
2032 for(Int_t iSl =0; iSl < track->GetNumberOfTRDslices(); iSl++) {
\r
2033 for(Int_t iPl =0; iPl<6; iPl++) trdslices[iPl*track->GetNumberOfTRDslices()+iSl] = track->GetTRDslice(iPl,iSl);
\r
2037 for(Int_t iPl=0;iPl<6;iPl++){
\r
2038 Double_t trdmom=track->GetTRDmomentum(iPl);
\r
2039 aodpid->SetTRDmomentum(iPl,trdmom);
\r
2042 aodpid->SetTRDsignal(track->GetNumberOfTRDslices()*6,trdslices.GetArray());
\r
2044 //TRD clusters and tracklets
\r
2045 aodpid->SetTRDncls(track->GetTRDncls());
\r
2046 aodpid->SetTRDntrackletsPID(track->GetTRDntrackletsPID());
\r
2049 Double_t times[AliAODPid::kSPECIES]; track->GetIntegratedTimes(times);
\r
2050 aodpid->SetIntegratedTimes(times);
\r
2052 Float_t tzeroTrack = fESDpid->GetTOFResponse().GetStartTime(track->P());
\r
2053 aodpid->SetTOFsignal(track->GetTOFsignal()-tzeroTrack);
\r
2055 Double_t tofRes[5];
\r
2056 for (Int_t iMass=0; iMass<5; iMass++){
\r
2057 tofRes[iMass]=(Double_t)fESDpid->GetTOFResponse().GetExpectedSigma(track->P(), times[iMass], AliPID::ParticleMass(iMass));
\r
2059 aodpid->SetTOFpidResolution(tofRes);
\r
2061 aodpid->SetHMPIDsignal(track->GetHMPIDsignal());
\r
2063 //Extrapolate track to EMCAL surface for AOD-level track-cluster matching
\r
2064 Double_t emcpos[3] = {0.,0.,0.};
\r
2065 Double_t emcmom[3] = {0.,0.,0.};
\r
2066 aodpid->SetEMCALPosition(emcpos);
\r
2067 aodpid->SetEMCALMomentum(emcmom);
\r
2069 Double_t hmpPid[5] = {0};
\r
2070 track->GetHMPIDpid(hmpPid);
\r
2071 aodpid->SetHMPIDprobs(hmpPid);
\r
2073 AliExternalTrackParam *outerparam = (AliExternalTrackParam*)track->GetOuterParam();
\r
2074 if(!outerparam) return;
\r
2076 //To be replaced by call to AliEMCALGeoUtils when the class becomes available
\r
2077 Bool_t okpos = outerparam->GetXYZ(emcpos);
\r
2078 Bool_t okmom = outerparam->GetPxPyPz(emcmom);
\r
2079 if(!(okpos && okmom)) return;
\r
2081 aodpid->SetEMCALPosition(emcpos);
\r
2082 aodpid->SetEMCALMomentum(emcmom);
\r
2086 Double_t AliAnalysisTaskESDfilter::Chi2perNDF(AliESDtrack* track)
\r
2088 // Calculate chi2 per ndf for track
\r
2089 Int_t nClustersTPC = track->GetTPCNcls();
\r
2091 if ( nClustersTPC > 5) {
\r
2092 return (track->GetTPCchi2()/Float_t(nClustersTPC - 5));
\r
2099 //______________________________________________________________________________
\r
2100 void AliAnalysisTaskESDfilter::Terminate(Option_t */*option*/)
\r
2102 // Terminate analysis
\r
2104 if (fDebug > 1) printf("AnalysisESDfilter: Terminate() \n");
\r
2107 //______________________________________________________________________________
\r
2108 void AliAnalysisTaskESDfilter::PrintMCInfo(AliStack *pStack,Int_t label){
\r
2110 if(!pStack)return;
\r
2111 label = TMath::Abs(label);
\r
2112 TParticle *part = pStack->Particle(label);
\r
2113 Printf("########################");
\r
2114 Printf("%s:%d %d UniqueID %d PDG %d P %3.3f",(char*)__FILE__,__LINE__,label,part->GetUniqueID(),part->GetPdgCode(),part->P());
\r
2116 TParticle* mother = part;
\r
2117 Int_t imo = part->GetFirstMother();
\r
2118 Int_t nprim = pStack->GetNprimary();
\r
2119 // while((imo >= nprim) && (mother->GetUniqueID() == 4)) {
\r
2120 while((imo >= nprim)) {
\r
2121 mother = pStack->Particle(imo);
\r
2122 Printf("Mother %s:%d Label %d UniqueID %d PDG %d P %3.3f",(char*)__FILE__,__LINE__,imo,mother->GetUniqueID(),mother->GetPdgCode(),mother->P());
\r
2124 imo = mother->GetFirstMother();
\r
2126 Printf("########################");
\r