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