]> git.uio.no Git - u/mrichter/AliRoot.git/blame - STEER/CreateAODfromESD.C
Changes to take into account the recent modifications in libESD.pkg
[u/mrichter/AliRoot.git] / STEER / CreateAODfromESD.C
CommitLineData
e985a35a 1#if !defined(__CINT__) || defined(__MAKECINT__)
2
cf22b3fc 3#include <Riostream.h>
4#include <TFile.h>
5#include <TTree.h>
6#include <TMath.h>
df9db588 7
8#include "AliAODEvent.h"
cf22b3fc 9#include "AliAODVertex.h"
10#include "AliAODTrack.h"
a9255000 11#include "AliAODCluster.h"
cf22b3fc 12
df9db588 13#include "AliESD.h"
14#include "AliESDtrack.h"
15#include "AliESDVertex.h"
16#include "AliESDv0.h"
cf22b3fc 17#include "AliESDCascade.h"
18#include "AliESDCaloCluster.h"
df9db588 19
e985a35a 20#endif
21
cf22b3fc 22void CreateAODfromESD(const char *inFileName = "AliESDs.root",
23 const char *outFileName = "AliAOD.root") {
df9db588 24
25 // create an AliAOD object
26 AliAODEvent *aod = new AliAODEvent();
27 aod->CreateStdContent();
28
29 // open the file
30 TFile *outFile = TFile::Open(outFileName, "RECREATE");
31
32 // create the tree
33 TTree *aodTree = new TTree("AOD", "AliAOD tree");
34 aodTree->Branch(aod->GetList());
35
36 // connect to ESD
37 TFile *inFile = TFile::Open(inFileName, "READ");
38 TTree *t = (TTree*) inFile->Get("esdTree");
39 TBranch *b = t->GetBranch("ESD");
40 AliESD *esd = 0;
41 b->SetAddress(&esd);
42
43 Int_t nEvents = b->GetEntries();
44
45 // loop over events and fill them
46 for (Int_t iEvent = 0; iEvent < nEvents; ++iEvent) {
47 b->GetEntry(iEvent);
48
49 // Multiplicity information needed by the header (to be revised!)
50 Int_t nTracks = esd->GetNumberOfTracks();
51 Int_t nPosTracks = 0;
52 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack)
53 if (esd->GetTrack(iTrack)->GetSign()> 0) nPosTracks++;
54
55 // create the header
56 aod->AddHeader(new AliAODHeader(esd ->GetEventNumber(),
57 esd->GetRunNumber(),
58 nTracks,
59 nPosTracks,
60 nTracks-nPosTracks,
61 esd->GetMagneticField(),
62 -999., // centrality; to be filled, still
63 esd->GetTriggerMask(),
64 esd->GetTriggerCluster(),
65 esd->GetEventType()));
66
67 Int_t nV0s = esd->GetNumberOfV0s();
68 Int_t nCascades = esd->GetNumberOfCascades();
69 Int_t nKinks = esd->GetNumberOfKinks();
70 Int_t nVertices = nV0s + nCascades + nKinks;
71
72 aod->ResetStd(nTracks, nVertices);
73
74
75 // Array to take into account the tracks already added to the AOD
76 Bool_t * usedTrack = NULL;
77 if (nTracks>0) {
78 usedTrack = new Bool_t[nTracks];
79 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) usedTrack[iTrack]=kFALSE;
80 }
81 // Array to take into account the V0s already added to the AOD
82 Bool_t * usedV0 = NULL;
83 if (nV0s>0) {
84 usedV0 = new Bool_t[nV0s];
85 for (Int_t iV0=0; iV0<nV0s; ++iV0) usedV0[iV0]=kFALSE;
86 }
cf22b3fc 87 // Array to take into account the kinks already added to the AOD
88 Bool_t * usedKink = NULL;
89 if (nKinks>0) {
90 usedKink = new Bool_t[nKinks];
91 for (Int_t iKink=0; iKink<nKinks; ++iKink) usedKink[iKink]=kFALSE;
92 }
df9db588 93
94 // Access to the AOD container of vertices
95 TClonesArray &vertices = *(aod->GetVertices());
96 Int_t jVertices=0;
97
98 // Access to the AOD container of tracks
99 TClonesArray &tracks = *(aod->GetTracks());
100 Int_t jTracks=0;
101
102 // Add primary vertex. The primary tracks will be defined
103 // after the loops on the composite objects (V0, cascades, kinks)
104 const AliESDVertex *vtx = esd->GetPrimaryVertex();
105
106 Double_t pos[3];
107 vtx->GetXYZ(pos); // position
108 Double_t cov[6];
109 vtx->GetCovMatrix(cov); //covariance matrix
110
111 AliAODVertex * primary = new(vertices[jVertices++])
112 AliAODVertex(pos, cov, vtx->GetChi2(), NULL, AliAODVertex::kPrimary);
113
114 // Create vertices starting from the most complex objects
115
116 // Cascades
117 for (Int_t nCascade = 0; nCascade < nCascades; ++nCascade) {
118 AliESDcascade *cascade = esd->GetCascade(nCascade);
119
120 Double_t posXi[3];
121 cascade->GetXYZ(posXi[0], posXi[1], posXi[2]);
122 Double_t covXi[6];
123 cascade->GetPosCovXi(covXi);
124
125 // Add the cascade vertex
126 AliAODVertex * vcascade = new(vertices[jVertices++]) AliAODVertex(posXi,
127 covXi,
128 cascade->GetChi2Xi(),
129 primary,
130 AliAODVertex::kCascade);
131
132 primary->AddDaughter(vcascade);
133
134 // Add the V0 from the cascade. The ESD class have to be optimized...
135 // Now we have to search for the corresponding Vo in the list of V0s
136 // using the indeces of the positive and negative tracks
137
138 Int_t posFromV0 = cascade->GetPindex();
139 Int_t negFromV0 = cascade->GetNindex();
140
141
142 AliESDv0 * v0 = 0x0;
143 Int_t indV0 = -1;
144
145 for (Int_t iV0=0; iV0<nV0s; ++iV0) {
146
147 v0 = esd->GetV0(iV0);
148 Int_t pos = v0->GetPindex();
149 Int_t neg = v0->GetNindex();
150
151 if (pos==posFromV0 && neg==negFromV0) {
152 indV0 = iV0;
153 break;
154 }
155 }
156
157 AliAODVertex * vV0FromCascade = 0x0;
158
cf22b3fc 159 if (indV0>-1 && !usedV0[indV0] ) {
df9db588 160
cf22b3fc 161 // the V0 exists in the array of V0s and is not used
162
df9db588 163 usedV0[indV0] = kTRUE;
164
165 Double_t posV0[3];
166 v0->GetXYZ(posV0[0], posV0[1], posV0[2]);
167 Double_t covV0[6];
168 v0->GetPosCov(covV0);
cf22b3fc 169
df9db588 170 vV0FromCascade = new(vertices[jVertices++]) AliAODVertex(posV0,
171 covV0,
172 v0->GetChi2V0(),
173 vcascade,
174 AliAODVertex::kV0);
175 } else {
176
cf22b3fc 177 // the V0 doesn't exist in the array of V0s or was used
178 cerr << "Error: event " << iEvent << " cascade " << nCascade
179 << " The V0 " << indV0
180 << " doesn't exist in the array of V0s or was used!" << endl;
181
df9db588 182 Double_t posV0[3];
183 cascade->GetXYZ(posV0[0], posV0[1], posV0[2]);
184 Double_t covV0[6];
185 cascade->GetPosCov(covV0);
186
187 vV0FromCascade = new(vertices[jVertices++]) AliAODVertex(posV0,
188 covV0,
189 v0->GetChi2V0(),
190 vcascade,
191 AliAODVertex::kV0);
192 vcascade->AddDaughter(vV0FromCascade);
193 }
194
195 // Add the positive tracks from the V0
196
cf22b3fc 197 if (! usedTrack[posFromV0]) {
198
199 usedTrack[posFromV0] = kTRUE;
df9db588 200
201 AliESDtrack *esdTrack = esd->GetTrack(posFromV0);
cf22b3fc 202
df9db588 203 Double_t p[3];
204 esdTrack->GetPxPyPz(p);
205
206 Double_t x[3];
207 esdTrack->GetXYZ(x);
208
209 Double_t cov[21];
210 esdTrack->GetCovarianceXYZPxPyPz(cov);
cf22b3fc 211
df9db588 212 Double_t pid[10];
213 esdTrack->GetESDpid(pid);
214
215 vV0FromCascade->AddDaughter(
cf22b3fc 216 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
df9db588 217 esdTrack->GetLabel(),
218 p, kTRUE,
219 x,
220 kFALSE,
221 cov,
222 (Short_t)esdTrack->GetSign(),
223 esdTrack->GetITSClusterMap(),
224 pid,
225 vV0FromCascade,
226 AliAODTrack::kSecondary)
227 );
228 }
cf22b3fc 229 else {
230 cerr << "Error: event " << iEvent << " cascade " << nCascade
231 << " track " << posFromV0 << " has already been used!" << endl;
232 }
df9db588 233
234 // Add the negative tracks from the V0
235
cf22b3fc 236 if (!usedTrack[negFromV0]) {
237
238 usedTrack[negFromV0] = kTRUE;
239
df9db588 240 AliESDtrack *esdTrack = esd->GetTrack(negFromV0);
cf22b3fc 241
df9db588 242 Double_t p[3];
243 esdTrack->GetPxPyPz(p);
244
245 Double_t x[3];
246 esdTrack->GetXYZ(x);
247
248 Double_t cov[21];
249 esdTrack->GetCovarianceXYZPxPyPz(cov);
250
251 Double_t pid[10];
252 esdTrack->GetESDpid(pid);
cf22b3fc 253
df9db588 254 vV0FromCascade->AddDaughter(
255 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
256 esdTrack->GetLabel(),
257 p,
258 kTRUE,
259 x,
260 kFALSE,
261 cov,
262 (Short_t)esdTrack->GetSign(),
263 esdTrack->GetITSClusterMap(),
264 pid,
265 vV0FromCascade,
266 AliAODTrack::kSecondary)
267 );
268 }
cf22b3fc 269 else {
270 cerr << "Error: event " << iEvent << " cascade " << nCascade
271 << " track " << negFromV0 << " has already been used!" << endl;
272 }
df9db588 273
274 // Add the bachelor track from the cascade
275
276 Int_t bachelor = cascade->GetBindex();
df9db588 277
cf22b3fc 278 if(!usedTrack[bachelor]) {
279
280 usedTrack[bachelor] = kTRUE;
281
282 AliESDtrack *esdTrack = esd->GetTrack(bachelor);
283
df9db588 284 Double_t p[3];
285 esdTrack->GetPxPyPz(p);
286
287 Double_t x[3];
288 esdTrack->GetXYZ(x);
289
290 Double_t cov[21];
291 esdTrack->GetCovarianceXYZPxPyPz(cov);
292
293 Double_t pid[10];
294 esdTrack->GetESDpid(pid);
295
296 vcascade->AddDaughter(
297 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
298 esdTrack->GetLabel(),
299 p,
300 kTRUE,
301 x,
302 kFALSE,
303 cov,
304 (Short_t)esdTrack->GetSign(),
305 esdTrack->GetITSClusterMap(),
306 pid,
307 vcascade,
308 AliAODTrack::kSecondary)
309 );
310 }
cf22b3fc 311 else {
312 cerr << "Error: event " << iEvent << " cascade " << nCascade
313 << " track " << bachelor << " has already been used!" << endl;
314 }
df9db588 315
316 // Add the primary track of the cascade (if any)
317
cf22b3fc 318 } // end of the loop on cascades
df9db588 319
320 // V0s
321
322 for (Int_t nV0 = 0; nV0 < nV0s; ++nV0) {
323
324 if (usedV0[nV0]) continue; // skip if aready added to the AOD
325
326 AliESDv0 *v0 = esd->GetV0(nV0);
327
328 Double_t posV0[3];
329 v0->GetXYZ(posV0[0], posV0[1], posV0[2]);
330 Double_t covV0[6];
331 v0->GetPosCov(covV0);
332
333 AliAODVertex * vV0 =
334 new(vertices[jVertices++]) AliAODVertex(posV0,
335 covV0,
336 v0->GetChi2V0(),
337 primary,
338 AliAODVertex::kV0);
339 primary->AddDaughter(vV0);
340
341 Int_t posFromV0 = v0->GetPindex();
342 Int_t negFromV0 = v0->GetNindex();
343
344 // Add the positive tracks from the V0
345
cf22b3fc 346 if (!usedTrack[posFromV0]) {
347
348 usedTrack[posFromV0] = kTRUE;
df9db588 349
350 AliESDtrack *esdTrack = esd->GetTrack(posFromV0);
351
352 Double_t p[3];
353 esdTrack->GetPxPyPz(p);
354
355 Double_t x[3];
356 esdTrack->GetXYZ(x);
357
358 Double_t cov[21];
359 esdTrack->GetCovarianceXYZPxPyPz(cov);
cf22b3fc 360
df9db588 361 Double_t pid[10];
362 esdTrack->GetESDpid(pid);
363
364 vV0->AddDaughter(
365 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
366 esdTrack->GetLabel(),
367 p, kTRUE,
368 x,
369 kFALSE,
370 cov,
371 (Short_t)esdTrack->GetSign(),
372 esdTrack->GetITSClusterMap(),
373 pid,
374 vV0,
375 AliAODTrack::kSecondary)
376 );
377 }
cf22b3fc 378 else {
379 cerr << "Error: event " << iEvent << " V0 " << nV0
380 << " track " << posFromV0 << " has already been used!" << endl;
381 }
df9db588 382
383 // Add the negative tracks from the V0
384
cf22b3fc 385 if (!usedTrack[negFromV0]) {
386
387 usedTrack[negFromV0] = kTRUE;
df9db588 388
389 AliESDtrack *esdTrack = esd->GetTrack(negFromV0);
390
391 Double_t p[3];
392 esdTrack->GetPxPyPz(p);
393
394 Double_t x[3];
395 esdTrack->GetXYZ(x);
396
397 Double_t cov[21];
398 esdTrack->GetCovarianceXYZPxPyPz(cov);
399
400 Double_t pid[10];
401 esdTrack->GetESDpid(pid);
402
403 vV0->AddDaughter(
404 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
405 esdTrack->GetLabel(),
406 p,
407 kTRUE,
408 x,
409 kFALSE,
410 cov,
411 (Short_t)esdTrack->GetSign(),
412 esdTrack->GetITSClusterMap(),
413 pid,
414 vV0,
415 AliAODTrack::kSecondary)
416 );
417 }
cf22b3fc 418 else {
419 cerr << "Error: event " << iEvent << " V0 " << nV0
420 << " track " << negFromV0 << " has already been used!" << endl;
421 }
df9db588 422
cf22b3fc 423 } // end of the loop on V0s
df9db588 424
cf22b3fc 425 // Kinks: it is a big mess the access to the information in the kinks
426 // The loop is on the tracks in order to find the mother and daugther of each kink
427
428
429 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) {
430
431
432 AliESDtrack * track = esd->GetTrack(iTrack);
433
434
435 Int_t ikink = track->GetKinkIndex(0);
436
437 if (ikink) {
438 // Negative kink index: mother, positive: daughter
439
440 // Search for the second track of the kink
441
442 for (Int_t jTrack = iTrack+1; jTrack<nTracks; ++jTrack) {
443
444 AliESDtrack * track1 = esd->GetTrack(jTrack);
445
446 Int_t jkink = track1->GetKinkIndex(0);
447
448 if ( TMath::Abs(ikink)==TMath::Abs(jkink) ) {
449
450 // The two tracks are from the same kink
451
452 if (usedKink[TMath::Abs(ikink)-1]) continue; // skip used kinks
453
454 Int_t imother = -1;
455 Int_t idaughter = -1;
456
457 if (ikink<0 && jkink>0) {
458
459 imother = iTrack;
460 idaughter = jTrack;
461 }
462 else if (ikink>0 && jkink<0) {
463
464 imother = jTrack;
465 idaughter = iTrack;
466 }
467 else {
468 cerr << "Error: Wrong combination of kink indexes: "
469 << ikink << " " << jkink << endl;
470 continue;
471 }
472
473 // Add the mother track
474
475 AliAODTrack * mother = NULL;
476
477 if (!usedTrack[imother]) {
478
479 usedTrack[imother] = kTRUE;
480
481 AliESDtrack *esdTrack = esd->GetTrack(imother);
482
483 Double_t p[3];
484 esdTrack->GetPxPyPz(p);
485
486 Double_t x[3];
487 esdTrack->GetXYZ(x);
488
489 Double_t cov[21];
490 esdTrack->GetCovarianceXYZPxPyPz(cov);
491
492 Double_t pid[10];
493 esdTrack->GetESDpid(pid);
494
495
496 mother =
497 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
498 esdTrack->GetLabel(),
499 p,
500 kTRUE,
501 x,
502 kFALSE,
503 cov,
504 (Short_t)esdTrack->GetSign(),
505 esdTrack->GetITSClusterMap(),
506 pid,
507 primary,
508 AliAODTrack::kPrimary);
509 primary->AddDaughter(mother);
510 }
511 else {
512 cerr << "Error: event " << iEvent << " kink " << TMath::Abs(ikink)-1
513 << " track " << imother << " has already been used!" << endl;
514 }
515
516 // Add the kink vertex
517 AliESDkink * kink = esd->GetKink(TMath::Abs(ikink)-1);
518
519 AliAODVertex * vkink =
520 new(vertices[jVertices++]) AliAODVertex(kink->GetPosition(),
521 NULL,
522 0.,
523 mother,
524 AliAODVertex::kKink);
525 // Add the daughter track
526
527 AliAODTrack * daughter = NULL;
528
529 if (!usedTrack[idaughter]) {
530
531 usedTrack[idaughter] = kTRUE;
532
533 AliESDtrack *esdTrack = esd->GetTrack(idaughter);
534
535 Double_t p[3];
536 esdTrack->GetPxPyPz(p);
537
538 Double_t x[3];
539 esdTrack->GetXYZ(x);
540
541 Double_t cov[21];
542 esdTrack->GetCovarianceXYZPxPyPz(cov);
543
544 Double_t pid[10];
545 esdTrack->GetESDpid(pid);
546
547
548 daughter =
549 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
550 esdTrack->GetLabel(),
551 p,
552 kTRUE,
553 x,
554 kFALSE,
555 cov,
556 (Short_t)esdTrack->GetSign(),
557 esdTrack->GetITSClusterMap(),
558 pid,
559 vkink,
560 AliAODTrack::kPrimary);
561 vkink->AddDaughter(daughter);
562 }
563 else {
564 cerr << "Error: event " << iEvent << " kink " << TMath::Abs(ikink)-1
565 << " track " << idaughter << " has already been used!" << endl;
566 }
567
568
569 }
570 }
571
572 }
df9db588 573
df9db588 574 }
cf22b3fc 575
df9db588 576
cf22b3fc 577 // Tracks (primary and orphan)
df9db588 578
579 for (Int_t nTrack = 0; nTrack < nTracks; ++nTrack) {
580
581
582 if (usedTrack[nTrack]) continue;
583
584 AliESDtrack *esdTrack = esd->GetTrack(nTrack);
585
586 Double_t p[3];
587 esdTrack->GetPxPyPz(p);
588
589 Double_t x[3];
590 esdTrack->GetXYZ(x);
591
592 Double_t cov[21];
593 esdTrack->GetCovarianceXYZPxPyPz(cov);
594
595 Double_t pid[10];
596 esdTrack->GetESDpid(pid);
cf22b3fc 597
598 Float_t impactXY, impactZ;
599
600 esdTrack->GetImpactParameters(impactXY,impactZ);
601
602 if (impactXY<3) {
603 // track inside the beam pipe
df9db588 604
cf22b3fc 605 primary->AddDaughter(
606 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
df9db588 607 esdTrack->GetLabel(),
608 p,
609 kTRUE,
610 x,
611 kFALSE,
612 cov,
613 (Short_t)esdTrack->GetSign(),
614 esdTrack->GetITSClusterMap(),
615 pid,
616 primary,
617 AliAODTrack::kPrimary)
cf22b3fc 618 );
619 }
620 else {
621 // outside the beam pipe: orphan track
622 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
623 esdTrack->GetLabel(),
624 p,
625 kTRUE,
626 x,
627 kFALSE,
628 cov,
629 (Short_t)esdTrack->GetSign(),
630 esdTrack->GetITSClusterMap(),
631 pid,
632 NULL,
633 AliAODTrack::kOrphan);
634 }
635 } // end of loop on tracks
636
637
638 // Access to the AOD container of vertices
a9255000 639 TClonesArray &clusters = *(aod->GetClusters());
cf22b3fc 640 Int_t jClusters=0;
641
642 // Clusters
643 Int_t nClusters = esd->GetNumberOfCaloClusters();
644
645 for (Int_t iClust=0; iClust<nClusters; ++iClust) {
646
647 AliESDCaloCluster * cluster = esd->GetCaloCluster(iClust);
648
649 Int_t id = cluster->GetID();
650 Int_t label = -1;
651 Float_t energy = cluster->GetClusterEnergy();
652 Float_t x[3];
653 cluster->GetGlobalPosition(x);
654 Float_t * covMatrix = NULL;
655 Float_t * pid = NULL;
656 AliAODVertex *prodVertex = primary;
657 AliAODTrack *primTrack = NULL;
a9255000 658 Char_t ttype=AliAODCluster::kUndef;
cf22b3fc 659
a9255000 660 if (cluster->IsPHOS()) ttype=AliAODCluster::kPHOSNeutral;
cf22b3fc 661 else if (cluster->IsEMCAL()) {
662
663 if (cluster->GetClusterType() == AliESDCaloCluster::kPseudoCluster)
a9255000 664 ttype = AliAODCluster::kEMCALPseudoCluster;
cf22b3fc 665 else
a9255000 666 ttype = AliAODCluster::kEMCALClusterv1;
cf22b3fc 667
668 }
df9db588 669
a9255000 670 new(clusters[jClusters++]) AliAODCluster(id,
cf22b3fc 671 label,
672 energy,
673 x,
674 covMatrix,
675 pid,
676 prodVertex,
677 primTrack,
678 ttype);
679
680 } // end of loop on clusters
681
682
683 delete [] usedTrack;
684 delete [] usedV0;
685 delete [] usedKink;
686
df9db588 687 // fill the tree for this event
688 aodTree->Fill();
cf22b3fc 689 } // end of event loop
df9db588 690
691 aodTree->GetUserInfo()->Add(aod);
692
693 // close ESD file
694 inFile->Close();
695
696 // write the tree to the specified file
697 outFile = aodTree->GetCurrentFile();
698 outFile->cd();
699 aodTree->Write();
700 outFile->Close();
701
702}