]> 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
Bug fixed on decoding flags
[u/mrichter/AliRoot.git] / ANALYSIS / AliAnalysisTaskESDfilter.cxx
CommitLineData
96b85d73 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: AliAnalysisTaskESDfilter.cxx 24535 2008-03-16 22:43:30Z fca $ */
17
18#include <TChain.h>
19#include <TFile.h>
e6500713 20#include <TArrayI.h>
d7bdc804 21#include <TRandom.h>
96b85d73 22
23#include "AliAnalysisTaskESDfilter.h"
24#include "AliAnalysisManager.h"
25#include "AliESDEvent.h"
26#include "AliAODEvent.h"
27#include "AliESDInputHandler.h"
28#include "AliAODHandler.h"
29#include "AliAnalysisFilter.h"
96b85d73 30#include "AliESDMuonTrack.h"
31#include "AliESDVertex.h"
32#include "AliESDv0.h"
33#include "AliESDkink.h"
34#include "AliESDcascade.h"
35#include "AliESDPmdTrack.h"
36#include "AliESDCaloCluster.h"
37#include "AliESDCaloCells.h"
38#include "AliMultiplicity.h"
39#include "AliLog.h"
40
41ClassImp(AliAnalysisTaskESDfilter)
42
43////////////////////////////////////////////////////////////////////////
44
45AliAnalysisTaskESDfilter::AliAnalysisTaskESDfilter():
46 AliAnalysisTaskSE(),
47 fTrackFilter(0x0),
48 fKinkFilter(0x0),
d7bdc804 49 fV0Filter(0x0),
50 fHighPthreshold(0),
51 fPtshape(0x0)
96b85d73 52{
53 // Default constructor
54}
55
56AliAnalysisTaskESDfilter::AliAnalysisTaskESDfilter(const char* name):
57 AliAnalysisTaskSE(name),
58 fTrackFilter(0x0),
59 fKinkFilter(0x0),
d7bdc804 60 fV0Filter(0x0),
61 fHighPthreshold(0),
62 fPtshape(0x0)
96b85d73 63{
64 // Constructor
65}
66
67void AliAnalysisTaskESDfilter::UserCreateOutputObjects()
68{
69// Create the output container
70 OutputTree()->GetUserInfo()->Add(fTrackFilter);
71}
72
73void AliAnalysisTaskESDfilter::Init()
74{
75 // Initialization
76 if (fDebug > 1) AliInfo("Init() \n");
77 // Call configuration file
78}
79
80
81void AliAnalysisTaskESDfilter::UserExec(Option_t */*option*/)
82{
83// Execute analysis for current event
84//
85
86 Long64_t ientry = Entry();
2827d046 87 if (fDebug > 0) printf("Filter: Analysing event # %5d\n", (Int_t) ientry);
d7bdc804 88 if (fHighPthreshold == 0) AliInfo("detector PID signals are stored in each track");
89 if (!fPtshape) AliInfo("detector PID signals are not stored below the pt threshold");
96b85d73 90 ConvertESDtoAOD();
91}
92
93void AliAnalysisTaskESDfilter::ConvertESDtoAOD() {
94 // ESD Filter analysis task executed for each event
95 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
96 AliESD* old = esd->GetAliESDOld();
97
98 // set arrays and pointers
99 Float_t posF[3];
100 Double_t pos[3];
101 Double_t p[3];
102 Double_t p_pos[3];
103 Double_t p_neg[3];
8d69965b 104 Double_t p_pos_atv0[3];
105 Double_t p_neg_atv0[3];
96b85d73 106 Double_t covVtx[6];
107 Double_t covTr[21];
108 Double_t pid[10];
109
110 for (Int_t i = 0; i < 6; i++) covVtx[i] = 0.;
111 for (Int_t i = 0; i < 21; i++) covTr [i] = 0.;
112
113
114 // loop over events and fill them
115
116 // Multiplicity information needed by the header (to be revised!)
117 Int_t nTracks = esd->GetNumberOfTracks();
8d69965b 118 // if (fDebug > 0) printf("-------------------Bo: Number of ESD tracks %d \n",nTracks);
119
96b85d73 120 Int_t nPosTracks = 0;
3e492e99 121// for (Int_t iTrack = 0; iTrack < nTracks; ++iTrack)
122// if (esd->GetTrack(iTrack)->GetSign()> 0) nPosTracks++;
96b85d73 123
124 // Update the header
125
126 AliAODHeader* header = AODEvent()->GetHeader();
127 header->SetRunNumber(esd->GetRunNumber());
128 if (old) {
129 header->SetBunchCrossNumber(0);
130 header->SetOrbitNumber(0);
131 header->SetPeriodNumber(0);
132 header->SetEventType(0);
133 header->SetMuonMagFieldScale(-999.); // FIXME
134 header->SetCentrality(-999.); // FIXME
135 } else {
136 header->SetBunchCrossNumber(esd->GetBunchCrossNumber());
137 header->SetOrbitNumber(esd->GetOrbitNumber());
138 header->SetPeriodNumber(esd->GetPeriodNumber());
139 header->SetEventType(esd->GetEventType());
140 header->SetMuonMagFieldScale(-999.); // FIXME
141 header->SetCentrality(-999.); // FIXME
142 }
143
144 header->SetTriggerMask(esd->GetTriggerMask());
145 header->SetTriggerCluster(esd->GetTriggerCluster());
146 header->SetMagneticField(esd->GetMagneticField());
147 header->SetZDCN1Energy(esd->GetZDCN1Energy());
148 header->SetZDCP1Energy(esd->GetZDCP1Energy());
149 header->SetZDCN2Energy(esd->GetZDCN2Energy());
150 header->SetZDCP2Energy(esd->GetZDCP2Energy());
151 header->SetZDCEMEnergy(esd->GetZDCEMEnergy(0),esd->GetZDCEMEnergy(1));
96b85d73 152//
153//
154 Int_t nV0s = esd->GetNumberOfV0s();
155 Int_t nCascades = esd->GetNumberOfCascades();
156 Int_t nKinks = esd->GetNumberOfKinks();
157 Int_t nVertices = nV0s + 2*nCascades /*could lead to two vertices, one V0 and the Xi */+ nKinks + 1 /* = prim. vtx*/;
158 Int_t nJets = 0;
159 Int_t nCaloClus = esd->GetNumberOfCaloClusters();
160 Int_t nFmdClus = 0;
161 Int_t nPmdClus = esd->GetNumberOfPmdTracks();
162
2827d046 163 if (fDebug > 0)
164 printf(" NV0=%d NCASCADES=%d NKINKS=%d\n", nV0s, nCascades, nKinks);
96b85d73 165
166 AODEvent()->ResetStd(nTracks, nVertices, nV0s+nCascades, nJets, nCaloClus, nFmdClus, nPmdClus);
167
168 AliAODTrack *aodTrack = 0x0;
d7bdc804 169 AliAODPid *detpid = 0x0;
170 Double_t timezero = 0; //TO BE FIXED
7dd90a3c 171 AliAODv0 *aodV0 = 0x0;
172
8d69965b 173 // RefArray to take into account the tracks associated to V0s
174 TRefArray *v0DaughterTracks = NULL;
175 if (nTracks>0) v0DaughterTracks = new TRefArray(nTracks);
176
96b85d73 177 // Array to take into account the tracks already added to the AOD
178 Bool_t * usedTrack = NULL;
179 if (nTracks>0) {
180 usedTrack = new Bool_t[nTracks];
181 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) usedTrack[iTrack]=kFALSE;
182 }
183 // Array to take into account the V0s already added to the AOD
184 Bool_t * usedV0 = NULL;
185 if (nV0s>0) {
186 usedV0 = new Bool_t[nV0s];
187 for (Int_t iV0=0; iV0<nV0s; ++iV0) usedV0[iV0]=kFALSE;
188 }
189 // Array to take into account the kinks already added to the AOD
190 Bool_t * usedKink = NULL;
191 if (nKinks>0) {
192 usedKink = new Bool_t[nKinks];
193 for (Int_t iKink=0; iKink<nKinks; ++iKink) usedKink[iKink]=kFALSE;
194 }
195
196 // Access to the AOD container of vertices
197 TClonesArray &vertices = *(AODEvent()->GetVertices());
198 Int_t jVertices=0;
199
200 // Access to the AOD container of tracks
201 TClonesArray &tracks = *(AODEvent()->GetTracks());
202 Int_t jTracks=0;
203
204 // Access to the AOD container of V0s
205 TClonesArray &V0s = *(AODEvent()->GetV0s());
206 Int_t jV0s=0;
207
208 // Add primary vertex. The primary tracks will be defined
209 // after the loops on the composite objects (V0, cascades, kinks)
210 const AliESDVertex *vtx = esd->GetPrimaryVertex();
211
212 vtx->GetXYZ(pos); // position
213 vtx->GetCovMatrix(covVtx); //covariance matrix
214
215 AliAODVertex * primary = new(vertices[jVertices++])
c06d243f 216 AliAODVertex(pos, covVtx, vtx->GetChi2toNDF(), NULL, -1, AliAODVertex::kPrimary);
2827d046 217 if (fDebug > 0) primary->Print();
96b85d73 218
219 // Create vertices starting from the most complex objects
220 Double_t chi2 = 0.;
221
222 // Cascades
223 for (Int_t nCascade = 0; nCascade < nCascades; ++nCascade) {
224 AliESDcascade *cascade = esd->GetCascade(nCascade);
225
226 cascade->GetXYZ(pos[0], pos[1], pos[2]);
227
228 if (!old) {
229 chi2 = cascade->GetChi2Xi(); // = chi2/NDF since NDF = 2*2-3
230 cascade->GetPosCovXi(covVtx);
231 } else {
232 chi2 = -999.;
233 for (Int_t i = 0; i < 6; i++) covVtx[i] = 0.;
234 }
235 // Add the cascade vertex
236 AliAODVertex * vcascade = new(vertices[jVertices++]) AliAODVertex(pos,
237 covVtx,
238 chi2,
239 primary,
240 nCascade,
241 AliAODVertex::kCascade);
242
243 primary->AddDaughter(vcascade);
244
245 // Add the V0 from the cascade. The ESD class have to be optimized...
246 // Now we have to search for the corresponding Vo in the list of V0s
247 // using the indeces of the positive and negative tracks
248
249 Int_t posFromV0 = cascade->GetPindex();
250 Int_t negFromV0 = cascade->GetNindex();
251
252
253 AliESDv0 * v0 = 0x0;
254 Int_t indV0 = -1;
255
256 for (Int_t iV0=0; iV0<nV0s; ++iV0) {
257
258 v0 = esd->GetV0(iV0);
259 Int_t posV0 = v0->GetPindex();
260 Int_t negV0 = v0->GetNindex();
261
262 if (posV0==posFromV0 && negV0==negFromV0) {
263 indV0 = iV0;
264 break;
265 }
266 }
267
268 AliAODVertex * vV0FromCascade = 0x0;
269
270 if (indV0>-1 && !usedV0[indV0] ) {
271
272 // the V0 exists in the array of V0s and is not used
273
274 usedV0[indV0] = kTRUE;
275
276 v0->GetXYZ(pos[0], pos[1], pos[2]);
277 if (!old) {
278 chi2 = v0->GetChi2V0(); // = chi2/NDF since NDF = 2*2-3
279 v0->GetPosCov(covVtx);
280 } else {
281 chi2 = -999.;
282 for (Int_t i = 0; i < 6; i++) covVtx[i] = 0.;
283 }
284
285 vV0FromCascade = new(vertices[jVertices++]) AliAODVertex(pos,
286 covVtx,
287 chi2,
288 vcascade,
289 indV0,
290 AliAODVertex::kV0);
291 } else {
292
293 // the V0 doesn't exist in the array of V0s or was used
294// cerr << "Error: event " << esd->GetEventNumberInFile() << " cascade " << nCascade
295// << " The V0 " << indV0
296// << " doesn't exist in the array of V0s or was used!" << endl;
297
298 cascade->GetXYZ(pos[0], pos[1], pos[2]);
299
300 if (!old) {
301 chi2 = v0->GetChi2V0();
302 cascade->GetPosCov(covVtx);
303 } else {
304 chi2 = -999.;
305 for (Int_t i = 0; i < 6; i++) covVtx[i] = 0.;
306 }
307
308 vV0FromCascade = new(vertices[jVertices++]) AliAODVertex(pos,
309 covVtx,
310 chi2, // = chi2/NDF since NDF = 2*2-3 (AM)
311 vcascade,
312 indV0,
313 AliAODVertex::kV0);
314 vcascade->AddDaughter(vV0FromCascade);
315 }
316
317 // Add the positive tracks from the V0
318
319 if (posFromV0>-1 && !usedTrack[posFromV0]) {
320
321 usedTrack[posFromV0] = kTRUE;
322
323 AliESDtrack *esdTrack = esd->GetTrack(posFromV0);
324 esdTrack->GetPxPyPz(p_pos);
325 esdTrack->GetXYZ(pos);
326 esdTrack->GetCovarianceXYZPxPyPz(covTr);
327 esdTrack->GetESDpid(pid);
328
329 vV0FromCascade->AddDaughter(aodTrack =
330 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
331 esdTrack->GetLabel(),
332 p_pos,
333 kTRUE,
334 pos,
335 kFALSE,
336 covTr,
337 (Short_t)esdTrack->GetSign(),
338 esdTrack->GetITSClusterMap(),
339 pid,
340 vV0FromCascade,
341 kTRUE, // check if this is right
342 kFALSE, // check if this is right
343 AliAODTrack::kSecondary)
344 );
3e492e99 345 if (esdTrack->GetSign() > 0) nPosTracks++;
96b85d73 346 aodTrack->ConvertAliPIDtoAODPID();
4fcd837c 347 aodTrack->SetFlags(esdTrack->GetStatus());
d7bdc804 348 SetAODPID(esdTrack,aodTrack,detpid,timezero);
96b85d73 349 }
350 else {
351// cerr << "Error: event " << esd->GetEventNumberInFile() << " cascade " << nCascade
352// << " track " << posFromV0 << " has already been used!" << endl;
353 }
354
355 // Add the negative tracks from the V0
356
357 if (negFromV0>-1 && !usedTrack[negFromV0]) {
358
359 usedTrack[negFromV0] = kTRUE;
360
361 AliESDtrack *esdTrack = esd->GetTrack(negFromV0);
362 esdTrack->GetPxPyPz(p_neg);
363 esdTrack->GetXYZ(pos);
364 esdTrack->GetCovarianceXYZPxPyPz(covTr);
365 esdTrack->GetESDpid(pid);
366
367 vV0FromCascade->AddDaughter(aodTrack =
368 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
369 esdTrack->GetLabel(),
370 p_neg,
371 kTRUE,
372 pos,
373 kFALSE,
374 covTr,
375 (Short_t)esdTrack->GetSign(),
376 esdTrack->GetITSClusterMap(),
377 pid,
378 vV0FromCascade,
379 kTRUE, // check if this is right
380 kFALSE, // check if this is right
381 AliAODTrack::kSecondary)
382 );
3e492e99 383 if (esdTrack->GetSign() > 0) nPosTracks++;
96b85d73 384 aodTrack->ConvertAliPIDtoAODPID();
4fcd837c 385 aodTrack->SetFlags(esdTrack->GetStatus());
d7bdc804 386 SetAODPID(esdTrack,aodTrack,detpid,timezero);
96b85d73 387 }
388 else {
389// cerr << "Error: event " << esd->GetEventNumberInFile() << " cascade " << nCascade
390// << " track " << negFromV0 << " has already been used!" << endl;
391 }
392
393 // add it to the V0 array as well
394 Double_t d0[2] = { -999., -99.};
395 // counting is probably wrong
396 new(V0s[jV0s++]) AliAODv0(vV0FromCascade, -999., -99., p_pos, p_neg, d0); // to be refined
397
398 // Add the bachelor track from the cascade
399
400 Int_t bachelor = cascade->GetBindex();
401
402 if(bachelor>-1 && !usedTrack[bachelor]) {
403
404 usedTrack[bachelor] = kTRUE;
405
406 AliESDtrack *esdTrack = esd->GetTrack(bachelor);
407 esdTrack->GetPxPyPz(p);
408 esdTrack->GetXYZ(pos);
409 esdTrack->GetCovarianceXYZPxPyPz(covTr);
410 esdTrack->GetESDpid(pid);
411
412 vcascade->AddDaughter(aodTrack =
413 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
414 esdTrack->GetLabel(),
415 p,
416 kTRUE,
417 pos,
418 kFALSE,
419 covTr,
420 (Short_t)esdTrack->GetSign(),
421 esdTrack->GetITSClusterMap(),
422 pid,
423 vcascade,
424 kTRUE, // check if this is right
425 kFALSE, // check if this is right
426 AliAODTrack::kSecondary)
427 );
3e492e99 428 if (esdTrack->GetSign() > 0) nPosTracks++;
96b85d73 429 aodTrack->ConvertAliPIDtoAODPID();
4fcd837c 430 aodTrack->SetFlags(esdTrack->GetStatus());
d7bdc804 431 SetAODPID(esdTrack,aodTrack,detpid,timezero);
96b85d73 432 }
433 else {
434// cerr << "Error: event " << esd->GetEventNumberInFile() << " cascade " << nCascade
435// << " track " << bachelor << " has already been used!" << endl;
436 }
437
438 // Add the primary track of the cascade (if any)
439
440 } // end of the loop on cascades
441
442 // V0s
443
444 for (Int_t nV0 = 0; nV0 < nV0s; ++nV0) {
445
7dd90a3c 446 if (usedV0[nV0]) continue; // skip if already added to the AOD
96b85d73 447
448 AliESDv0 *v0 = esd->GetV0(nV0);
449
450 v0->GetXYZ(pos[0], pos[1], pos[2]);
451
452 if (!old) {
453 chi2 = v0->GetChi2V0(); // = chi2/NDF since NDF = 2*2-3
454 v0->GetPosCov(covVtx);
455 } else {
456 chi2 = -999.;
457 for (Int_t i = 0; i < 6; i++) covVtx[i] = 0.;
458 }
459
460
461 AliAODVertex * vV0 =
462 new(vertices[jVertices++]) AliAODVertex(pos,
463 covVtx,
464 chi2,
465 primary,
466 nV0,
467 AliAODVertex::kV0);
468 primary->AddDaughter(vV0);
469
470 Int_t posFromV0 = v0->GetPindex();
471 Int_t negFromV0 = v0->GetNindex();
7dd90a3c 472
473 Float_t dcaPosToPrimVertexXYZ[2] = { 999., 999.}; // ..[0] = in XY plane and ..[1] = in Z
474 Float_t dcaNegToPrimVertexXYZ[2] = { 999., 999.}; // ..[0] = in XY plane and ..[1] = in Z
475 Double_t dcaDaughterToPrimVertex[2] = { 999., 999.}; // ..[0] = Pos and ..[1] = Neg
96b85d73 476
7dd90a3c 477 Double_t dcaV0Daughters = v0->GetDcaV0Daughters();
478 Double_t dcaV0ToPrimVertex = v0->GetD();
479
8d69965b 480 v0->GetPPxPyPz(p_pos_atv0[0],p_pos_atv0[1],p_pos_atv0[2]);
481 v0->GetNPxPyPz(p_neg_atv0[0],p_neg_atv0[1],p_neg_atv0[2]);
7dd90a3c 482
96b85d73 483 // Add the positive tracks from the V0
484
8d69965b 485 if (posFromV0>-1){
486 AliESDtrack *esdTrack = esd->GetTrack(posFromV0);
487 esdTrack->GetPxPyPz(p_pos);
488 esdTrack->GetXYZ(pos);
489 esdTrack->GetCovarianceXYZPxPyPz(covTr);
490 esdTrack->GetESDpid(pid);
491 esdTrack->GetImpactParameters(dcaPosToPrimVertexXYZ[0],dcaPosToPrimVertexXYZ[1]);
492 if (!usedTrack[posFromV0]) {
96b85d73 493 usedTrack[posFromV0] = kTRUE;
8d69965b 494 aodTrack = new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
495 esdTrack->GetLabel(),
496 p_pos,
497 kTRUE,
498 pos,
499 kFALSE,
500 covTr,
501 (Short_t)esdTrack->GetSign(),
502 esdTrack->GetITSClusterMap(),
503 pid,
504 vV0,
505 kTRUE, // check if this is right
506 kFALSE, // check if this is right
507 AliAODTrack::kSecondary);
508 v0DaughterTracks->AddAt(aodTrack,posFromV0);
509 // if (fDebug > 0) printf("-------------------Bo: pos track from original pt %.3f \n",aodTrack->Pt());
510 if (esdTrack->GetSign() > 0) nPosTracks++;
96b85d73 511 aodTrack->ConvertAliPIDtoAODPID();
4fcd837c 512 aodTrack->SetFlags(esdTrack->GetStatus());
d7bdc804 513 SetAODPID(esdTrack,aodTrack,detpid,timezero);
8d69965b 514 }
515 else {
516 aodTrack = dynamic_cast<AliAODTrack*>(v0DaughterTracks->At(posFromV0));
517 // if (fDebug > 0) printf("-------------------Bo pos track from refArray pt %.3f \n",aodTrack->Pt());
518 }
519 vV0->AddDaughter(aodTrack);
96b85d73 520 }
521
522 // Add the negative tracks from the V0
523
8d69965b 524 if (negFromV0>-1){
525 AliESDtrack *esdTrack = esd->GetTrack(negFromV0);
526 esdTrack->GetPxPyPz(p_neg);
527 esdTrack->GetXYZ(pos);
528 esdTrack->GetCovarianceXYZPxPyPz(covTr);
529 esdTrack->GetESDpid(pid);
530 esdTrack->GetImpactParameters(dcaNegToPrimVertexXYZ[0],dcaNegToPrimVertexXYZ[1]);
531
532 if (!usedTrack[negFromV0]) {
96b85d73 533 usedTrack[negFromV0] = kTRUE;
8d69965b 534 aodTrack = new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
535 esdTrack->GetLabel(),
536 p_neg,
537 kTRUE,
538 pos,
539 kFALSE,
540 covTr,
541 (Short_t)esdTrack->GetSign(),
542 esdTrack->GetITSClusterMap(),
543 pid,
544 vV0,
545 kTRUE, // check if this is right
546 kFALSE, // check if this is right
547 AliAODTrack::kSecondary);
548
549 v0DaughterTracks->AddAt(aodTrack,negFromV0);
550 // if (fDebug > 0) printf("-------------------Bo: neg track from original pt %.3f \n",aodTrack->Pt());
551 if (esdTrack->GetSign() > 0) nPosTracks++;
96b85d73 552 aodTrack->ConvertAliPIDtoAODPID();
4fcd837c 553 aodTrack->SetFlags(esdTrack->GetStatus());
d7bdc804 554 SetAODPID(esdTrack,aodTrack,detpid,timezero);
8d69965b 555 }
556 else {
557 aodTrack = dynamic_cast<AliAODTrack*>(v0DaughterTracks->At(negFromV0));
558 // if (fDebug > 0) printf("-------------------Bo neg track from refArray pt %.3f \n",aodTrack->Pt());
559 }
560 vV0->AddDaughter(aodTrack);
96b85d73 561 }
7dd90a3c 562 dcaDaughterToPrimVertex[0] =
8d69965b 563 TMath::Sqrt(dcaPosToPrimVertexXYZ[0]*dcaPosToPrimVertexXYZ[0]
564 +dcaPosToPrimVertexXYZ[1]*dcaPosToPrimVertexXYZ[1]);
7dd90a3c 565 dcaDaughterToPrimVertex[1] =
8d69965b 566 TMath::Sqrt(dcaNegToPrimVertexXYZ[0]*dcaNegToPrimVertexXYZ[0]
567 +dcaNegToPrimVertexXYZ[1]*dcaNegToPrimVertexXYZ[1]);
96b85d73 568 // add it to the V0 array as well
7dd90a3c 569 aodV0 = new(V0s[jV0s++])
8d69965b 570 AliAODv0(vV0, dcaV0Daughters, dcaV0ToPrimVertex, p_pos_atv0, p_neg_atv0, dcaDaughterToPrimVertex); // to be refined
7dd90a3c 571 // set the aod v0 on-the-fly status
572 aodV0->SetOnFlyStatus(v0->GetOnFlyStatus());
8d69965b 573// if (fDebug > 0)
574// printf("-------------------Bo: In the v0 loop: onFlyStatus=%d,dcaV0Daughters=%.3f dcaV0ToPrimVertex=%.3f dcaDaughterToPrimVertex=[%.3f,%.3f]\n",
575// aodV0->GetOnFlyStatus(), aodV0->DcaV0Daughters(), aodV0->DcaV0ToPrimVertex(),
576// aodV0->DcaPosToPrimVertex(),aodV0->DcaNegToPrimVertex());
96b85d73 577 }
578 V0s.Expand(jV0s);
579 // end of the loop on V0s
580
581 // Kinks: it is a big mess the access to the information in the kinks
582 // The loop is on the tracks in order to find the mother and daugther of each kink
583
584
585 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) {
586
587 AliESDtrack * esdTrack = esd->GetTrack(iTrack);
588
589 Int_t ikink = esdTrack->GetKinkIndex(0);
590
591 if (ikink && nKinks) {
592 // Negative kink index: mother, positive: daughter
593
594 // Search for the second track of the kink
595
596 for (Int_t jTrack = iTrack+1; jTrack<nTracks; ++jTrack) {
597
598 AliESDtrack * esdTrack1 = esd->GetTrack(jTrack);
599
600 Int_t jkink = esdTrack1->GetKinkIndex(0);
601
602 if ( TMath::Abs(ikink)==TMath::Abs(jkink) ) {
603
604 // The two tracks are from the same kink
605
606 if (usedKink[TMath::Abs(ikink)-1]) continue; // skip used kinks
607
608 Int_t imother = -1;
609 Int_t idaughter = -1;
610
611 if (ikink<0 && jkink>0) {
612
613 imother = iTrack;
614 idaughter = jTrack;
615 }
616 else if (ikink>0 && jkink<0) {
617
618 imother = jTrack;
619 idaughter = iTrack;
620 }
621 else {
622// cerr << "Error: Wrong combination of kink indexes: "
623// << ikink << " " << jkink << endl;
624 continue;
625 }
626
627 // Add the mother track
628
629 AliAODTrack * mother = NULL;
630
631 if (!usedTrack[imother]) {
632
633 usedTrack[imother] = kTRUE;
634
03a5cc9f 635 AliESDtrack *esdTrackM = esd->GetTrack(imother);
636 esdTrackM->GetPxPyPz(p);
637 esdTrackM->GetXYZ(pos);
638 esdTrackM->GetCovarianceXYZPxPyPz(covTr);
639 esdTrackM->GetESDpid(pid);
96b85d73 640
641 mother =
03a5cc9f 642 new(tracks[jTracks++]) AliAODTrack(esdTrackM->GetID(),
643 esdTrackM->GetLabel(),
96b85d73 644 p,
645 kTRUE,
646 pos,
647 kFALSE,
648 covTr,
03a5cc9f 649 (Short_t)esdTrackM->GetSign(),
650 esdTrackM->GetITSClusterMap(),
96b85d73 651 pid,
652 primary,
653 kTRUE, // check if this is right
654 kTRUE, // check if this is right
655 AliAODTrack::kPrimary);
03a5cc9f 656 if (esdTrackM->GetSign() > 0) nPosTracks++;
657 mother->SetFlags(esdTrackM->GetStatus());
635ef56c 658 mother->ConvertAliPIDtoAODPID();
96b85d73 659 primary->AddDaughter(mother);
660 mother->ConvertAliPIDtoAODPID();
d7bdc804 661 SetAODPID(esdTrackM,mother,detpid,timezero);
96b85d73 662 }
663 else {
664// cerr << "Error: event " << esd->GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1
665// << " track " << imother << " has already been used!" << endl;
666 }
667
668 // Add the kink vertex
669 AliESDkink * kink = esd->GetKink(TMath::Abs(ikink)-1);
670
671 AliAODVertex * vkink =
672 new(vertices[jVertices++]) AliAODVertex(kink->GetPosition(),
673 NULL,
674 0.,
675 mother,
676 esdTrack->GetID(), // This is the track ID of the mother's track!
677 AliAODVertex::kKink);
678 // Add the daughter track
679
680 AliAODTrack * daughter = NULL;
681
682 if (!usedTrack[idaughter]) {
683
684 usedTrack[idaughter] = kTRUE;
685
03a5cc9f 686 AliESDtrack *esdTrackD = esd->GetTrack(idaughter);
687 esdTrackD->GetPxPyPz(p);
688 esdTrackD->GetXYZ(pos);
689 esdTrackD->GetCovarianceXYZPxPyPz(covTr);
690 esdTrackD->GetESDpid(pid);
96b85d73 691
692 daughter =
03a5cc9f 693 new(tracks[jTracks++]) AliAODTrack(esdTrackD->GetID(),
694 esdTrackD->GetLabel(),
96b85d73 695 p,
696 kTRUE,
697 pos,
698 kFALSE,
699 covTr,
03a5cc9f 700 (Short_t)esdTrackD->GetSign(),
701 esdTrackD->GetITSClusterMap(),
96b85d73 702 pid,
703 vkink,
704 kTRUE, // check if this is right
705 kTRUE, // check if this is right
3123180c 706 AliAODTrack::kSecondary);
03a5cc9f 707 if (esdTrackD->GetSign() > 0) nPosTracks++;
708 daughter->SetFlags(esdTrackD->GetStatus());
635ef56c 709 daughter->ConvertAliPIDtoAODPID();
96b85d73 710 vkink->AddDaughter(daughter);
711 daughter->ConvertAliPIDtoAODPID();
d7bdc804 712 SetAODPID(esdTrackD,daughter,detpid,timezero);
96b85d73 713 }
714 else {
715// cerr << "Error: event " << esd->GetEventNumberInFile() << " kink " << TMath::Abs(ikink)-1
716// << " track " << idaughter << " has already been used!" << endl;
717 }
718 }
719 }
720 }
721 }
722
723
724 // Tracks (primary and orphan)
725
3e492e99 726 if (fDebug > 0) printf("NUMBER OF ESD TRACKS %5d\n", nTracks);
96b85d73 727
728 for (Int_t nTrack = 0; nTrack < nTracks; ++nTrack) {
729
730
731 if (usedTrack[nTrack]) continue;
732
733 AliESDtrack *esdTrack = esd->GetTrack(nTrack);
734 UInt_t selectInfo = 0;
735 //
736 // Track selection
737 if (fTrackFilter) {
738 selectInfo = fTrackFilter->IsSelected(esdTrack);
739 if (!selectInfo) continue;
740 }
741
742 //
743 esdTrack->GetPxPyPz(p);
744 esdTrack->GetXYZ(pos);
745 esdTrack->GetCovarianceXYZPxPyPz(covTr);
746 esdTrack->GetESDpid(pid);
747
748 Float_t impactXY, impactZ;
749
750 esdTrack->GetImpactParameters(impactXY,impactZ);
751
752 if (impactXY<3) {
753 // track inside the beam pipe
754
755 primary->AddDaughter(aodTrack =
756 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
757 esdTrack->GetLabel(),
758 p,
759 kTRUE,
760 pos,
761 kFALSE,
762 covTr,
763 (Short_t)esdTrack->GetSign(),
764 esdTrack->GetITSClusterMap(),
765 pid,
766 primary,
767 kTRUE, // check if this is right
768 kTRUE, // check if this is right
769 AliAODTrack::kPrimary,
770 selectInfo)
771 );
3e492e99 772 if (esdTrack->GetSign() > 0) nPosTracks++;
4fcd837c 773 aodTrack->SetFlags(esdTrack->GetStatus());
96b85d73 774 aodTrack->ConvertAliPIDtoAODPID();
d7bdc804 775 SetAODPID(esdTrack,aodTrack,detpid,timezero);
96b85d73 776 }
777 else {
778 // outside the beam pipe: orphan track
779 // Don't write them anymore!
780 continue;
781 /*
782 aodTrack =
783 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
784 esdTrack->GetLabel(),
785 p,
786 kTRUE,
787 pos,
788 kFALSE,
789 covTr,
790 (Short_t)esdTrack->GetSign(),
791 esdTrack->GetITSClusterMap(),
792 pid,
793 NULL,
794 kFALSE, // check if this is right
795 kFALSE, // check if this is right
796 AliAODTrack::kOrphan,
797 selectInfo);
3e492e99 798 if (esdTrack->GetSign() > 0) nPosTracks++;
4fcd837c 799 aodTrack->SetFlags(esdTrack->GetStatus());
96b85d73 800 aodTrack->ConvertAliPIDtoAODPID();
801 */
802 }
803 } // end of loop on tracks
804
3e492e99 805 // Update number of AOD tracks in header at the end of track loop (M.G.)
806 header->SetRefMultiplicity(jTracks);
807 header->SetRefMultiplicityPos(nPosTracks);
808 header->SetRefMultiplicityNeg(jTracks - nPosTracks);
809 if (fDebug > 0)
810 printf(" NAODTRACKS=%d NPOS=%d NNEG=%d\n", jTracks, nPosTracks, jTracks - nPosTracks);
811 // Do not shrink the array of tracks - other filters may add to it (M.G)
812// tracks.Expand(jTracks); // remove 'empty slots' due to unwritten tracks
96b85d73 813
814 // Access to the AOD container of PMD clusters
815 TClonesArray &pmdClusters = *(AODEvent()->GetPmdClusters());
816 Int_t jPmdClusters=0;
817
818 for (Int_t iPmd = 0; iPmd < nPmdClus; ++iPmd) {
819 // file pmd clusters, to be revised!
820 AliESDPmdTrack *pmdTrack = esd->GetPmdTrack(iPmd);
821 Int_t nLabel = 0;
822 Int_t *label = 0x0;
03a5cc9f 823 Double_t posPmd[3] = { pmdTrack->GetClusterX(), pmdTrack->GetClusterY(), pmdTrack->GetClusterZ()};
824 Double_t pidPmd[9] = { 0., 0., 0., 0., 0., 0., 0., 0., 0. }; // to be revised!
96b85d73 825 // type not set!
826 // assoc cluster not set
03a5cc9f 827 new(pmdClusters[jPmdClusters++]) AliAODPmdCluster(iPmd, nLabel, label, pmdTrack->GetClusterADC(), posPmd, pidPmd);
96b85d73 828 }
829
830 // Access to the AOD container of clusters
831 TClonesArray &caloClusters = *(AODEvent()->GetCaloClusters());
832 Int_t jClusters=0;
833
834 for (Int_t iClust=0; iClust<nCaloClus; ++iClust) {
835
836 AliESDCaloCluster * cluster = esd->GetCaloCluster(iClust);
837
e6500713 838 Int_t id = cluster->GetID();
839 Int_t nLabel = cluster->GetNLabels();
840 TArrayI* labels = cluster->GetLabels();
841 Int_t *label = 0;
842 if (labels) label = (cluster->GetLabels())->GetArray();
843
96b85d73 844 Float_t energy = cluster->E();
845 cluster->GetPosition(posF);
e6500713 846 Char_t ttype = AliAODCluster::kUndef;
96b85d73 847
848 if (cluster->GetClusterType() == AliESDCaloCluster::kPHOSCluster) {
849 ttype=AliAODCluster::kPHOSNeutral;
850 }
851 else if (cluster->GetClusterType() == AliESDCaloCluster::kEMCALClusterv1) {
852 ttype = AliAODCluster::kEMCALClusterv1;
853 }
854
855
856 AliAODCaloCluster *caloCluster = new(caloClusters[jClusters++]) AliAODCaloCluster(id,
857 nLabel,
858 label,
859 energy,
860 pos,
861 NULL,
862 ttype);
863
e6500713 864 caloCluster->SetCaloCluster(cluster->GetDistanceToBadChannel(),
865 cluster->GetClusterDisp(),
78902954 866 cluster->GetM20(), cluster->GetM02(),
867 cluster->GetEmcCpvDistance(),
868 cluster->GetNExMax(),cluster->GetTOF()) ;
e6500713 869
3cbe9a4d 870 caloCluster->SetPIDFromESD(cluster->GetPid());
e6500713 871 caloCluster->SetNCells(cluster->GetNCells());
872 caloCluster->SetCellsAbsId(cluster->GetCellsAbsId());
873 caloCluster->SetCellsAmplitudeFraction(cluster->GetCellsAmplitudeFraction());
78902954 874
875 TArrayI* matchedT = cluster->GetTracksMatched();
876 if (matchedT) {
877 for (Int_t im = 0; im < matchedT->GetSize(); im++) {
878 caloCluster->AddTrackMatched((esd->GetTrack(im)));
879 }
880 }
881
96b85d73 882 }
883 caloClusters.Expand(jClusters); // resize TObjArray to 'remove' slots for pseudo clusters
884 // end of loop on calo clusters
885
886 // fill EMCAL cell info
887 if (esd->GetEMCALCells()) { // protection against missing ESD information
888 AliESDCaloCells &esdEMcells = *(esd->GetEMCALCells());
889 Int_t nEMcell = esdEMcells.GetNumberOfCells() ;
890
891 AliAODCaloCells &aodEMcells = *(AODEvent()->GetEMCALCells());
892 aodEMcells.CreateContainer(nEMcell);
893 aodEMcells.SetType(AliAODCaloCells::kEMCAL);
894 for (Int_t iCell = 0; iCell < nEMcell; iCell++) {
895 aodEMcells.SetCell(iCell,esdEMcells.GetCellNumber(iCell),esdEMcells.GetAmplitude(iCell));
896 }
897 aodEMcells.Sort();
898 }
899
900 // fill PHOS cell info
901 if (esd->GetPHOSCells()) { // protection against missing ESD information
902 AliESDCaloCells &esdPHcells = *(esd->GetPHOSCells());
903 Int_t nPHcell = esdPHcells.GetNumberOfCells() ;
904
905 AliAODCaloCells &aodPHcells = *(AODEvent()->GetPHOSCells());
906 aodPHcells.CreateContainer(nPHcell);
907 aodPHcells.SetType(AliAODCaloCells::kPHOS);
908 for (Int_t iCell = 0; iCell < nPHcell; iCell++) {
909 aodPHcells.SetCell(iCell,esdPHcells.GetCellNumber(iCell),esdPHcells.GetAmplitude(iCell));
910 }
911 aodPHcells.Sort();
912 }
913
914 // tracklets
915 AliAODTracklets &SPDTracklets = *(AODEvent()->GetTracklets());
916 const AliMultiplicity *mult = esd->GetMultiplicity();
917 if (mult) {
918 if (mult->GetNumberOfTracklets()>0) {
919 SPDTracklets.CreateContainer(mult->GetNumberOfTracklets());
920
921 for (Int_t n=0; n<mult->GetNumberOfTracklets(); n++) {
0939e22a 922 SPDTracklets.SetTracklet(n, mult->GetTheta(n), mult->GetPhi(n), mult->GetDeltaPhi(n), mult->GetLabel(n, 0), mult->GetLabel(n, 1));
96b85d73 923 }
924 }
925 } else {
926 //Printf("ERROR: AliMultiplicity could not be retrieved from ESD");
927 }
928
929 delete [] usedTrack;
930 delete [] usedV0;
931 delete [] usedKink;
8d69965b 932 delete v0DaughterTracks;
96b85d73 933
934 return;
935}
936
d7bdc804 937void AliAnalysisTaskESDfilter::SetAODPID(AliESDtrack *esdtrack, AliAODTrack *aodtrack, AliAODPid *detpid, Double_t timezero)
938{
939 //
940 // Setter for the raw PID detector signals
941 //
942
943 if(esdtrack->Pt()>fHighPthreshold) {
944 detpid = new AliAODPid();
945 detpid->SetDetectorRawSignals(esdtrack,timezero);
946 aodtrack->SetDetPID(detpid);
947 } else {
948 if(fPtshape){
949 if(esdtrack->Pt()> fPtshape->GetXmin()){
950 Double_t y = fPtshape->Eval(esdtrack->Pt())/fPtshape->Eval(fHighPthreshold);
951 if(gRandom->Rndm(0)<1./y){
952 detpid = new AliAODPid();
953 detpid->SetDetectorRawSignals(esdtrack,timezero);
954 aodtrack->SetDetPID(detpid);
955 }//end rndm
956 }//end if p < pmin
957 }//end if p function
958 }// end else
959}
960
961
96b85d73 962void AliAnalysisTaskESDfilter::Terminate(Option_t */*option*/)
963{
964// Terminate analysis
965//
966 if (fDebug > 1) printf("AnalysisESDfilter: Terminate() \n");
967}
968