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