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