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