1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 //-----------------------------------------------------------------
19 // Implementation of the AliESDEvent class
20 // This is the class to deal with during the physics analysis of data.
21 // It also ensures the backward compatibility with the old ESD format.
23 AliESDEvent *ev= new AliESDEvent();
24 ev->ReadFromTree(esdTree);
26 for (Int_t i=0; i<nev; i++) {
28 if(ev->GetAliESDOld())ev->CopyFromOldESD();
30 // The AliESDInputHandler does this automatically for you
32 // Origin: Christian Klein-Boesing, CERN, Christian.Klein-Boesing@cern.ch
33 //-----------------------------------------------------------------
36 #include "TRefArray.h"
39 #include "AliESDEvent.h"
40 #include "AliESDfriend.h"
41 #include "AliESDVZERO.h"
42 #include "AliESDFMD.h"
44 #include "AliESDMuonTrack.h"
45 #include "AliESDPmdTrack.h"
46 #include "AliESDTrdTrack.h"
47 #include "AliESDVertex.h"
48 #include "AliESDcascade.h"
49 #include "AliESDPmdTrack.h"
50 #include "AliESDTrdTrack.h"
51 #include "AliESDVertex.h"
52 #include "AliESDcascade.h"
53 #include "AliESDkink.h"
54 #include "AliESDtrack.h"
55 #include "AliESDHLTtrack.h"
56 #include "AliESDCaloCluster.h"
57 #include "AliESDCaloCells.h"
59 #include "AliESDFMD.h"
60 #include "AliESDVZERO.h"
61 #include "AliMultiplicity.h"
62 #include "AliRawDataErrorLog.h"
64 #include "AliESDACORDE.h"
69 // here we define the names, some classes are no TNamed, therefore the classnames
71 const char* AliESDEvent::fgkESDListName[kESDListN] = {"AliESDRun",
95 "AliRawDataErrorLogs",
98 //______________________________________________________________________________
99 AliESDEvent::AliESDEvent():
101 fESDObjects(new TList()),
115 fSPDPileupVertices(0),
116 fTrkPileupVertices(0),
125 fEMCALCells(0), fPHOSCells(0),
132 fFirstEMCALCluster(-1),
134 fFirstPHOSCluster(-1)
137 //______________________________________________________________________________
138 AliESDEvent::AliESDEvent(const AliESDEvent& esd):
140 fESDObjects(new TList()),
141 fESDRun(new AliESDRun(*esd.fESDRun)),
142 fHeader(new AliESDHeader(*esd.fHeader)),
143 fESDZDC(new AliESDZDC(*esd.fESDZDC)),
144 fESDFMD(new AliESDFMD(*esd.fESDFMD)),
145 fESDVZERO(new AliESDVZERO(*esd.fESDVZERO)),
146 fESDTZERO(new AliESDTZERO(*esd.fESDTZERO)),
147 fTPCVertex(new AliESDVertex(*esd.fTPCVertex)),
148 fSPDVertex(new AliESDVertex(*esd.fSPDVertex)),
149 fPrimaryVertex(new AliESDVertex(*esd.fPrimaryVertex)),
150 fSPDMult(new AliMultiplicity(*esd.fSPDMult)),
151 fPHOSTrigger(new AliESDCaloTrigger(*esd.fPHOSTrigger)),
152 fEMCALTrigger(new AliESDCaloTrigger(*esd.fEMCALTrigger)),
153 fESDACORDE(new AliESDACORDE(*esd.fESDACORDE)),
154 fSPDPileupVertices(new TClonesArray(*esd.fSPDPileupVertices)),
155 fTrkPileupVertices(new TClonesArray(*esd.fTrkPileupVertices)),
156 fTracks(new TClonesArray(*esd.fTracks)),
157 fMuonTracks(new TClonesArray(*esd.fMuonTracks)),
158 fPmdTracks(new TClonesArray(*esd.fPmdTracks)),
159 fTrdTracks(new TClonesArray(*esd.fTrdTracks)),
160 fV0s(new TClonesArray(*esd.fV0s)),
161 fCascades(new TClonesArray(*esd.fCascades)),
162 fKinks(new TClonesArray(*esd.fKinks)),
163 fCaloClusters(new TClonesArray(*esd.fCaloClusters)),
164 fEMCALCells(new AliESDCaloCells(*esd.fEMCALCells)),
165 fPHOSCells(new AliESDCaloCells(*esd.fPHOSCells)),
166 fErrorLogs(new TClonesArray(*esd.fErrorLogs)),
167 fESDOld(new AliESD(*esd.fESDOld)),
168 fESDFriendOld(new AliESDfriend(*esd.fESDFriendOld)),
169 fConnected(esd.fConnected),
170 fUseOwnList(esd.fUseOwnList),
171 fEMCALClusters(esd.fEMCALClusters),
172 fFirstEMCALCluster(esd.fFirstEMCALCluster),
173 fPHOSClusters(esd.fPHOSClusters),
174 fFirstPHOSCluster(esd.fFirstPHOSCluster)
177 // CKB init in the constructor list and only add here ...
182 AddObject(fESDVZERO);
183 AddObject(fESDTZERO);
184 AddObject(fTPCVertex);
185 AddObject(fSPDVertex);
186 AddObject(fPrimaryVertex);
188 AddObject(fPHOSTrigger);
189 AddObject(fEMCALTrigger);
190 AddObject(fSPDPileupVertices);
191 AddObject(fTrkPileupVertices);
193 AddObject(fMuonTracks);
194 AddObject(fPmdTracks);
195 AddObject(fTrdTracks);
197 AddObject(fCascades);
199 AddObject(fCaloClusters);
200 AddObject(fEMCALCells);
201 AddObject(fPHOSCells);
202 AddObject(fErrorLogs);
203 AddObject(fESDACORDE);
209 //______________________________________________________________________________
210 AliESDEvent & AliESDEvent::operator=(const AliESDEvent& source) {
212 // Assignment operator
214 if(&source == this) return *this;
215 AliVEvent::operator=(source);
217 // This assumes that the list is already created
218 // and that the virtual void Copy(Tobject&) function
219 // is correctly implemented in the derived class
220 // otherwise only TObject::Copy() will be used
224 if((fESDObjects->GetSize()==0)&&(source.fESDObjects->GetSize()>=kESDListN)){
225 // We cover the case that we do not yet have the
226 // standard content but the source has it
230 TIter next(source.GetList());
233 while ((its = next())) {
234 name.Form("%s", its->GetName());
235 TObject *mine = fESDObjects->FindObject(name.Data());
237 TClass* pClass=TClass::GetClass(its->ClassName());
239 AliWarning(Form("Can not find class description for entry %s (%s)\n",
240 its->ClassName(), name.Data()));
244 mine=(TObject*)pClass->New();
246 // not in this: can be added to list
247 AliWarning(Form("%s:%d Could not find %s for copying \n",
248 (char*)__FILE__,__LINE__,name.Data()));
251 if(mine->InheritsFrom("TNamed")){
252 ((TNamed*)mine)->SetName(name);
254 else if(mine->InheritsFrom("TCollection")){
255 if(mine->InheritsFrom("TClonesArray"))
256 dynamic_cast<TClonesArray*>(mine)->SetClass(dynamic_cast<TClonesArray*>(its)->GetClass());
257 dynamic_cast<TCollection*>(mine)->SetName(name);
259 AliDebug(1, Form("adding object %s of type %s", mine->GetName(), mine->ClassName()));
263 if(!its->InheritsFrom("TCollection")){
267 else if(its->InheritsFrom("TClonesArray")){
268 // Create or expand the tclonesarray pointers
269 // so we can directly copy to the object
270 TClonesArray *its_tca = (TClonesArray*)its;
271 TClonesArray *mine_tca = (TClonesArray*)mine;
273 // this leaves the capacity of the TClonesArray the same
274 // except for a factor of 2 increase when size > capacity
275 // does not release any memory occupied by the tca
276 mine_tca->ExpandCreate(its_tca->GetEntriesFast());
277 for(int i = 0;i < its_tca->GetEntriesFast();++i){
279 TObject *mine_tca_obj = mine_tca->At(i);
280 TObject *its_tca_obj = its_tca->At(i);
281 // no need to delete first
282 // pointers within the class should be handled by Copy()...
283 // Can there be Empty slots?
284 its_tca_obj->Copy(*mine_tca_obj);
288 AliWarning(Form("%s:%d cannot copy TCollection \n",
289 (char*)__FILE__,__LINE__));
293 fConnected = source.fConnected;
294 fUseOwnList = source.fUseOwnList;
295 fEMCALClusters = source.fEMCALClusters;
296 fFirstEMCALCluster = source.fFirstEMCALCluster;
297 fPHOSClusters = source.fPHOSClusters;
298 fFirstPHOSCluster = source.fFirstPHOSCluster;
306 //______________________________________________________________________________
307 AliESDEvent::~AliESDEvent()
310 // Standard destructor
313 // everthing on the list gets deleted automatically
316 if(fESDObjects&&!fConnected)
325 void AliESDEvent::Copy(TObject &obj) const {
327 // interface to TOBject::Copy
328 // Copies the content of this into obj!
329 // bascially obj = *this
331 if(this==&obj)return;
332 AliESDEvent *robj = dynamic_cast<AliESDEvent*>(&obj);
333 if(!robj)return; // not an AliESEvent
338 //______________________________________________________________________________
339 void AliESDEvent::Reset()
343 // Std content + Non std content
345 // Reset the standard contents
348 // reset for the old data without AliESDEvent...
349 if(fESDOld)fESDOld->Reset();
351 fESDFriendOld->~AliESDfriend();
352 new (fESDFriendOld) AliESDfriend();
356 if(fESDObjects->GetSize()>kESDListN){
357 // we have non std content
358 // this also covers esdfriends
359 for(int i = kESDListN;i < fESDObjects->GetSize();++i){
360 TObject *pObject = fESDObjects->At(i);
362 if(pObject->InheritsFrom(TClonesArray::Class())){
363 ((TClonesArray*)pObject)->Delete();
365 else if(!pObject->InheritsFrom(TCollection::Class())){
366 ResetWithPlacementNew(pObject);
369 AliWarning(Form("No reset for %s (%s)\n",
370 pObject->ClassName()));
377 Bool_t AliESDEvent::ResetWithPlacementNew(TObject *pObject){
378 Long_t dtoronly = TObject::GetDtorOnly();
379 TClass *pClass = TClass::GetClass(pObject->ClassName());
380 TObject::SetDtorOnly(pObject);
382 // Recreate with placement new
383 pClass->New(pObject);
384 // Restore the state.
385 TObject::SetDtorOnly((void*)dtoronly);
389 void AliESDEvent::ResetStdContent()
391 // Reset the standard contents
392 if(fESDRun) fESDRun->Reset();
393 if(fHeader) fHeader->Reset();
394 if(fESDZDC) fESDZDC->Reset();
399 // reset by callin d'to /c'tor keep the pointer
400 fESDVZERO->~AliESDVZERO();
401 new (fESDVZERO) AliESDVZERO();
404 fESDACORDE->~AliESDACORDE();
405 new (fESDACORDE) AliESDACORDE();
407 if(fESDTZERO) fESDTZERO->Reset();
408 // CKB no clear/reset implemented
410 fTPCVertex->~AliESDVertex();
411 new (fTPCVertex) AliESDVertex();
412 fTPCVertex->SetName(fgkESDListName[kTPCVertex]);
415 fSPDVertex->~AliESDVertex();
416 new (fSPDVertex) AliESDVertex();
417 fSPDVertex->SetName(fgkESDListName[kSPDVertex]);
420 fPrimaryVertex->~AliESDVertex();
421 new (fPrimaryVertex) AliESDVertex();
422 fPrimaryVertex->SetName(fgkESDListName[kPrimaryVertex]);
425 fSPDMult->~AliMultiplicity();
426 new (fSPDMult) AliMultiplicity();
428 if(fPHOSTrigger)fPHOSTrigger->Reset();
429 if(fEMCALTrigger)fEMCALTrigger->Reset();
430 if(fSPDPileupVertices)fSPDPileupVertices->Delete();
431 if(fTrkPileupVertices)fTrkPileupVertices->Delete();
432 if(fTracks)fTracks->Delete();
433 if(fMuonTracks)fMuonTracks->Delete();
434 if(fPmdTracks)fPmdTracks->Delete();
435 if(fTrdTracks)fTrdTracks->Delete();
436 if(fV0s)fV0s->Delete();
437 if(fCascades)fCascades->Delete();
438 if(fKinks)fKinks->Delete();
439 if(fCaloClusters)fCaloClusters->Delete();
440 if(fPHOSCells)fPHOSCells->DeleteContainer();
441 if(fEMCALCells)fEMCALCells->DeleteContainer();
442 if(fErrorLogs) fErrorLogs->Delete();
444 // don't reset fconnected fConnected and the list
447 fFirstEMCALCluster=-1;
449 fFirstPHOSCluster=-1;
453 Int_t AliESDEvent::AddV0(const AliESDv0 *v) {
457 TClonesArray &fv = *fV0s;
458 Int_t idx=fV0s->GetEntriesFast();
459 new(fv[idx]) AliESDv0(*v);
463 //______________________________________________________________________________
464 void AliESDEvent::Print(Option_t *) const
467 // Print header information of the event
469 printf("ESD run information\n");
470 printf("Event # in file %d Bunch crossing # %d Orbit # %d Period # %d Run # %d Trigger %lld Magnetic field %f \n",
471 GetEventNumberInFile(),
472 GetBunchCrossNumber(),
477 GetMagneticField() );
478 printf("Vertex: (%.4f +- %.4f, %.4f +- %.4f, %.4f +- %.4f) cm\n",
479 fPrimaryVertex->GetXv(), fPrimaryVertex->GetXRes(),
480 fPrimaryVertex->GetYv(), fPrimaryVertex->GetYRes(),
481 fPrimaryVertex->GetZv(), fPrimaryVertex->GetZRes());
482 printf("Mean vertex in RUN: X=%.4f Y=%.4f cm\n",
483 GetDiamondX(),GetDiamondY());
484 printf("SPD Multiplicity. Number of tracklets %d \n",
485 fSPDMult->GetNumberOfTracklets());
486 printf("Number of pileup primary vertices reconstructed with SPD %d\n",
487 GetNumberOfPileupVerticesSPD());
488 printf("Number of pileup primary vertices reconstructed using the tracks %d\n",
489 GetNumberOfPileupVerticesTracks());
490 printf("Number of tracks: \n");
491 printf(" charged %d\n", GetNumberOfTracks());
492 printf(" muon %d\n", GetNumberOfMuonTracks());
493 printf(" pmd %d\n", GetNumberOfPmdTracks());
494 printf(" trd %d\n", GetNumberOfTrdTracks());
495 printf(" v0 %d\n", GetNumberOfV0s());
496 printf(" cascades %d\n", GetNumberOfCascades());
497 printf(" kinks %d\n", GetNumberOfKinks());
498 if(fPHOSCells)printf(" PHOSCells %d\n", fPHOSCells->GetNumberOfCells());
499 else printf(" PHOSCells not in the Event\n");
500 if(fEMCALCells)printf(" EMCALCells %d\n", fEMCALCells->GetNumberOfCells());
501 else printf(" EMCALCells not in the Event\n");
502 printf(" CaloClusters %d\n", GetNumberOfCaloClusters());
503 printf(" phos %d\n", GetNumberOfPHOSClusters());
504 printf(" emcal %d\n", GetNumberOfEMCALClusters());
505 printf(" FMD %s\n", (fESDFMD ? "yes" : "no"));
506 printf(" VZERO %s\n", (fESDVZERO ? "yes" : "no"));
511 void AliESDEvent::SetESDfriend(const AliESDfriend *ev) const {
513 // Attaches the complementary info to the ESD
517 // to be sure that we set the tracks also
518 // in case of old esds
519 // if(fESDOld)CopyFromOldESD();
521 Int_t ntrk=ev->GetNumberOfTracks();
523 for (Int_t i=0; i<ntrk; i++) {
524 const AliESDfriendTrack *f=ev->GetTrack(i);
525 GetTrack(i)->SetFriendTrack(f);
529 Bool_t AliESDEvent::RemoveKink(Int_t rm) const {
530 // ---------------------------------------------------------
531 // Remove a kink candidate and references to it from ESD,
532 // if this candidate does not come from a reconstructed decay
533 // Not yet implemented...
534 // ---------------------------------------------------------
535 Int_t last=GetNumberOfKinks()-1;
536 if ((rm<0)||(rm>last)) return kFALSE;
541 Bool_t AliESDEvent::RemoveV0(Int_t rm) const {
542 // ---------------------------------------------------------
543 // Remove a V0 candidate and references to it from ESD,
544 // if this candidate does not come from a reconstructed decay
545 // ---------------------------------------------------------
546 Int_t last=GetNumberOfV0s()-1;
547 if ((rm<0)||(rm>last)) return kFALSE;
549 AliESDv0 *v0=GetV0(rm);
550 Int_t idxP=v0->GetPindex(), idxN=v0->GetNindex();
553 Int_t lastIdxP=v0->GetPindex(), lastIdxN=v0->GetNindex();
557 // Check if this V0 comes from a reconstructed decay
558 Int_t ncs=GetNumberOfCascades();
559 for (Int_t n=0; n<ncs; n++) {
560 AliESDcascade *cs=GetCascade(n);
562 Int_t csIdxP=cs->GetPindex();
563 Int_t csIdxN=cs->GetNindex();
566 if (idxN==csIdxN) return kFALSE;
568 if (csIdxP==lastIdxP)
569 if (csIdxN==lastIdxN) used++;
572 //Replace the removed V0 with the last V0
573 TClonesArray &a=*fV0s;
574 delete a.RemoveAt(rm);
576 if (rm==last) return kTRUE;
578 //v0 is pointing to the last V0 candidate...
579 new (a[rm]) AliESDv0(*v0);
580 delete a.RemoveAt(last);
582 if (!used) return kTRUE;
585 // Remap the indices of the daughters of reconstructed decays
586 for (Int_t n=0; n<ncs; n++) {
587 AliESDcascade *cs=GetCascade(n);
590 Int_t csIdxP=cs->GetPindex();
591 Int_t csIdxN=cs->GetNindex();
593 if (csIdxP==lastIdxP)
594 if (csIdxN==lastIdxN) {
595 cs->AliESDv0::SetIndex(1,idxP);
596 cs->AliESDv0::SetIndex(0,idxN);
598 if (!used) return kTRUE;
605 Bool_t AliESDEvent::RemoveTrack(Int_t rm) const {
606 // ---------------------------------------------------------
607 // Remove a track and references to it from ESD,
608 // if this track does not come from a reconstructed decay
609 // ---------------------------------------------------------
610 Int_t last=GetNumberOfTracks()-1;
611 if ((rm<0)||(rm>last)) return kFALSE;
615 // Check if this track comes from the reconstructed primary vertices
616 if (fTPCVertex && fTPCVertex->GetStatus()) {
617 UShort_t *primIdx=fTPCVertex->GetIndices();
618 Int_t n=fTPCVertex->GetNIndices();
620 Int_t idx=Int_t(primIdx[n]);
621 if (rm==idx) return kFALSE;
622 if (idx==last) used++;
625 if (fPrimaryVertex && fPrimaryVertex->GetStatus()) {
626 UShort_t *primIdx=fPrimaryVertex->GetIndices();
627 Int_t n=fPrimaryVertex->GetNIndices();
629 Int_t idx=Int_t(primIdx[n]);
630 if (rm==idx) return kFALSE;
631 if (idx==last) used++;
635 // Check if this track comes from a reconstructed decay
636 Int_t nv0=GetNumberOfV0s();
637 for (Int_t n=0; n<nv0; n++) {
638 AliESDv0 *v0=GetV0(n);
640 Int_t idx=v0->GetNindex();
641 if (rm==idx) return kFALSE;
642 if (idx==last) used++;
645 if (rm==idx) return kFALSE;
646 if (idx==last) used++;
649 Int_t ncs=GetNumberOfCascades();
650 for (Int_t n=0; n<ncs; n++) {
651 AliESDcascade *cs=GetCascade(n);
653 Int_t idx=cs->GetIndex();
654 if (rm==idx) return kFALSE;
655 if (idx==last) used++;
659 if (rm==idx) return kFALSE;
660 if (idx==last) used++;
663 if (rm==idx) return kFALSE;
664 if (idx==last) used++;
667 Int_t nkn=GetNumberOfKinks();
668 for (Int_t n=0; n<nkn; n++) {
669 AliESDkink *kn=GetKink(n);
671 Int_t idx=kn->GetIndex(0);
672 if (rm==idx) return kFALSE;
673 if (idx==last) used++;
676 if (rm==idx) return kFALSE;
677 if (idx==last) used++;
680 // Check if this track is associated with a CaloCluster
681 Int_t ncl=GetNumberOfCaloClusters();
682 for (Int_t n=0; n<ncl; n++) {
683 AliESDCaloCluster *cluster=GetCaloCluster(n);
684 TArrayI *arr=cluster->GetTracksMatched();
685 Int_t s=arr->GetSize();
687 Int_t idx=arr->At(s);
688 if (rm==idx) return kFALSE;
689 if (idx==last) used++;
695 //Replace the removed track with the last track
696 TClonesArray &a=*fTracks;
697 delete a.RemoveAt(rm);
699 if (rm==last) return kTRUE;
701 AliESDtrack *t=GetTrack(last);
703 new (a[rm]) AliESDtrack(*t);
704 delete a.RemoveAt(last);
707 if (!used) return kTRUE;
710 // Remap the indices of the tracks used for the primary vertex reconstruction
711 if (fTPCVertex && fTPCVertex->GetStatus()) {
712 UShort_t *primIdx=fTPCVertex->GetIndices();
713 Int_t n=fTPCVertex->GetNIndices();
715 Int_t idx=Int_t(primIdx[n]);
717 primIdx[n]=Short_t(rm);
719 if (!used) return kTRUE;
723 if (fPrimaryVertex && fPrimaryVertex->GetStatus()) {
724 UShort_t *primIdx=fPrimaryVertex->GetIndices();
725 Int_t n=fPrimaryVertex->GetNIndices();
727 Int_t idx=Int_t(primIdx[n]);
729 primIdx[n]=Short_t(rm);
731 if (!used) return kTRUE;
736 // Remap the indices of the daughters of reconstructed decays
737 for (Int_t n=0; n<nv0; n++) {
738 AliESDv0 *v0=GetV0(n);
739 if (v0->GetIndex(0)==last) {
742 if (!used) return kTRUE;
744 if (v0->GetIndex(1)==last) {
747 if (!used) return kTRUE;
751 for (Int_t n=0; n<ncs; n++) {
752 AliESDcascade *cs=GetCascade(n);
753 if (cs->GetIndex()==last) {
756 if (!used) return kTRUE;
759 if (v0->GetIndex(0)==last) {
762 if (!used) return kTRUE;
764 if (v0->GetIndex(1)==last) {
767 if (!used) return kTRUE;
771 for (Int_t n=0; n<nkn; n++) {
772 AliESDkink *kn=GetKink(n);
773 if (kn->GetIndex(0)==last) {
776 if (!used) return kTRUE;
778 if (kn->GetIndex(1)==last) {
781 if (!used) return kTRUE;
785 // Remap the indices of the tracks accosicated with CaloClusters
786 for (Int_t n=0; n<ncl; n++) {
787 AliESDCaloCluster *cluster=GetCaloCluster(n);
788 TArrayI *arr=cluster->GetTracksMatched();
789 Int_t s=arr->GetSize();
791 Int_t idx=arr->At(s);
795 if (!used) return kTRUE;
804 Bool_t AliESDEvent::Clean(Float_t *cleanPars) {
806 // Remove the data which are not needed for the physics analysis.
808 // 1) Cleaning the V0 candidates
809 // ---------------------------
810 // If the cosine of the V0 pointing angle "csp" and
811 // the DCA between the daughter tracks "dca" does not satisfy
814 // csp > cleanPars[1] + dca/cleanPars[0]*(1.- cleanPars[1])
816 // an attempt to remove this V0 candidate from ESD is made.
818 // The V0 candidate gets removed if it does not belong to any
819 // recosntructed cascade decay
821 // 12.11.2007, optimal values: cleanPars[0]=0.5, cleanPars[1]=0.999
823 // 2) Cleaning the tracks
824 // ----------------------
825 // If track's transverse parameter is larger than cleanPars[2]
827 // track's longitudinal parameter is larger than cleanPars[3]
828 // an attempt to remove this track from ESD is made.
830 // The track gets removed if it does not come
831 // from a reconstructed decay
835 Float_t dcaMax=cleanPars[0];
836 Float_t cspMin=cleanPars[1];
838 Int_t nV0s=GetNumberOfV0s();
839 for (Int_t i=nV0s-1; i>=0; i--) {
840 AliESDv0 *v0=GetV0(i);
842 Float_t dca=v0->GetDcaV0Daughters();
843 Float_t csp=v0->GetV0CosineOfPointingAngle();
844 Float_t cspcut=cspMin + dca/dcaMax*(1.-cspMin);
845 if (csp > cspcut) continue;
846 if (RemoveV0(i)) rc=kTRUE;
850 Float_t dmax=cleanPars[2], zmax=cleanPars[3];
852 const AliESDVertex *vertex=GetPrimaryVertexSPD();
853 Bool_t vtxOK=vertex->GetStatus();
855 Int_t nTracks=GetNumberOfTracks();
856 for (Int_t i=nTracks-1; i>=0; i--) {
857 AliESDtrack *track=GetTrack(i);
858 Float_t xy,z; track->GetImpactParameters(xy,z);
859 if ((TMath::Abs(xy) > dmax) || (vtxOK && (TMath::Abs(z) > zmax))) {
860 if (RemoveTrack(i)) rc=kTRUE;
867 void AliESDEvent::AddPileupVertexSPD(const AliESDVertex *vtx)
869 // Add a pileup primary vertex reconstructed with SPD
870 TClonesArray &ftr = *fSPDPileupVertices;
871 new(ftr[fSPDPileupVertices->GetEntriesFast()]) AliESDVertex(*vtx);
874 void AliESDEvent::AddPileupVertexTracks(const AliESDVertex *vtx)
876 // Add a pileup primary vertex reconstructed with SPD
877 TClonesArray &ftr = *fTrkPileupVertices;
878 new(ftr[fTrkPileupVertices->GetEntriesFast()]) AliESDVertex(*vtx);
881 Int_t AliESDEvent::AddTrack(const AliESDtrack *t)
884 TClonesArray &ftr = *fTracks;
885 AliESDtrack * track = new(ftr[fTracks->GetEntriesFast()])AliESDtrack(*t);
886 track->SetID(fTracks->GetEntriesFast()-1);
887 return track->GetID();
890 void AliESDEvent::AddMuonTrack(const AliESDMuonTrack *t)
892 TClonesArray &fmu = *fMuonTracks;
893 new(fmu[fMuonTracks->GetEntriesFast()]) AliESDMuonTrack(*t);
896 void AliESDEvent::AddPmdTrack(const AliESDPmdTrack *t)
898 TClonesArray &fpmd = *fPmdTracks;
899 new(fpmd[fPmdTracks->GetEntriesFast()]) AliESDPmdTrack(*t);
902 void AliESDEvent::AddTrdTrack(const AliESDTrdTrack *t)
904 TClonesArray &ftrd = *fTrdTracks;
905 new(ftrd[fTrdTracks->GetEntriesFast()]) AliESDTrdTrack(*t);
911 Int_t AliESDEvent::AddKink(const AliESDkink *c)
914 TClonesArray &fk = *fKinks;
915 AliESDkink * kink = new(fk[fKinks->GetEntriesFast()]) AliESDkink(*c);
916 kink->SetID(fKinks->GetEntriesFast()); // CKB different from the other imps..
917 return fKinks->GetEntriesFast()-1;
921 void AliESDEvent::AddCascade(const AliESDcascade *c)
923 TClonesArray &fc = *fCascades;
924 new(fc[fCascades->GetEntriesFast()]) AliESDcascade(*c);
928 Int_t AliESDEvent::AddCaloCluster(const AliESDCaloCluster *c)
931 TClonesArray &fc = *fCaloClusters;
932 AliESDCaloCluster *clus = new(fc[fCaloClusters->GetEntriesFast()]) AliESDCaloCluster(*c);
933 clus->SetID(fCaloClusters->GetEntriesFast()-1);
934 return fCaloClusters->GetEntriesFast()-1;
938 void AliESDEvent::AddRawDataErrorLog(const AliRawDataErrorLog *log) const {
939 TClonesArray &errlogs = *fErrorLogs;
940 new(errlogs[errlogs.GetEntriesFast()]) AliRawDataErrorLog(*log);
943 void AliESDEvent::SetPrimaryVertexTPC(const AliESDVertex *vertex)
945 // Set the TPC vertex
946 // use already allocated space
948 *fTPCVertex = *vertex;
949 fTPCVertex->SetName(fgkESDListName[kTPCVertex]);
953 void AliESDEvent::SetPrimaryVertexSPD(const AliESDVertex *vertex)
955 // Set the SPD vertex
956 // use already allocated space
958 *fSPDVertex = *vertex;
959 fSPDVertex->SetName(fgkESDListName[kSPDVertex]);
963 void AliESDEvent::SetPrimaryVertexTracks(const AliESDVertex *vertex)
965 // Set the primary vertex reconstructed using he ESD tracks.
966 // use already allocated space
968 *fPrimaryVertex = *vertex;
969 fPrimaryVertex->SetName(fgkESDListName[kPrimaryVertex]);
973 const AliESDVertex * AliESDEvent::GetPrimaryVertex() const
976 // Get the "best" available reconstructed primary vertex.
979 if (fPrimaryVertex->GetStatus()) return fPrimaryVertex;
982 if (fSPDVertex->GetStatus()) return fSPDVertex;
984 if(fTPCVertex) return fTPCVertex;
986 AliWarning("No primary vertex available. Returning the \"default\"...");
990 void AliESDEvent::SetMultiplicity(const AliMultiplicity *mul)
992 // Set the SPD Multiplicity
999 void AliESDEvent::SetFMDData(AliESDFMD * obj)
1001 // use already allocated space
1007 void AliESDEvent::SetVZEROData(AliESDVZERO * obj)
1009 // use already allocated space
1014 void AliESDEvent::SetACORDEData(AliESDACORDE * obj)
1021 void AliESDEvent::GetESDfriend(AliESDfriend *ev) const
1024 // Extracts the complementary info from the ESD
1028 Int_t ntrk=GetNumberOfTracks();
1030 for (Int_t i=0; i<ntrk; i++) {
1031 AliESDtrack *t=GetTrack(i);
1032 const AliESDfriendTrack *f=t->GetFriendTrack();
1035 t->ReleaseESDfriendTrack();// Not to have two copies of "friendTrack"
1039 AliESDfriend *fr = (AliESDfriend*)(const_cast<AliESDEvent*>(this)->FindListObject("AliESDfriend"));
1040 if (fr) ev->SetVZEROfriend(fr->GetVZEROfriend());
1043 void AliESDEvent::AddObject(TObject* obj)
1045 // Add an object to the list of object.
1046 // Please be aware that in order to increase performance you should
1047 // refrain from using TObjArrays (if possible). Use TClonesArrays, instead.
1048 fESDObjects->SetOwner(kTRUE);
1049 fESDObjects->AddLast(obj);
1053 void AliESDEvent::GetStdContent()
1055 // set pointers for standard content
1056 // get by name much safer and not a big overhead since not called very often
1058 fESDRun = (AliESDRun*)fESDObjects->FindObject(fgkESDListName[kESDRun]);
1059 fHeader = (AliESDHeader*)fESDObjects->FindObject(fgkESDListName[kHeader]);
1060 fESDZDC = (AliESDZDC*)fESDObjects->FindObject(fgkESDListName[kESDZDC]);
1061 fESDFMD = (AliESDFMD*)fESDObjects->FindObject(fgkESDListName[kESDFMD]);
1062 fESDVZERO = (AliESDVZERO*)fESDObjects->FindObject(fgkESDListName[kESDVZERO]);
1063 fESDTZERO = (AliESDTZERO*)fESDObjects->FindObject(fgkESDListName[kESDTZERO]);
1064 fTPCVertex = (AliESDVertex*)fESDObjects->FindObject(fgkESDListName[kTPCVertex]);
1065 fSPDVertex = (AliESDVertex*)fESDObjects->FindObject(fgkESDListName[kSPDVertex]);
1066 fPrimaryVertex = (AliESDVertex*)fESDObjects->FindObject(fgkESDListName[kPrimaryVertex]);
1067 fSPDMult = (AliMultiplicity*)fESDObjects->FindObject(fgkESDListName[kSPDMult]);
1068 fPHOSTrigger = (AliESDCaloTrigger*)fESDObjects->FindObject(fgkESDListName[kPHOSTrigger]);
1069 fEMCALTrigger = (AliESDCaloTrigger*)fESDObjects->FindObject(fgkESDListName[kEMCALTrigger]);
1070 fSPDPileupVertices = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kSPDPileupVertices]);
1071 fTrkPileupVertices = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kTrkPileupVertices]);
1072 fTracks = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kTracks]);
1073 fMuonTracks = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kMuonTracks]);
1074 fPmdTracks = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kPmdTracks]);
1075 fTrdTracks = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kTrdTracks]);
1076 fV0s = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kV0s]);
1077 fCascades = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kCascades]);
1078 fKinks = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kKinks]);
1079 fCaloClusters = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kCaloClusters]);
1080 fEMCALCells = (AliESDCaloCells*)fESDObjects->FindObject(fgkESDListName[kEMCALCells]);
1081 fPHOSCells = (AliESDCaloCells*)fESDObjects->FindObject(fgkESDListName[kPHOSCells]);
1082 fErrorLogs = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kErrorLogs]);
1083 fESDACORDE = (AliESDACORDE*)fESDObjects->FindObject(fgkESDListName[kESDACORDE]);
1087 void AliESDEvent::SetStdNames(){
1088 // Set the names of the standard contents
1090 if(fESDObjects->GetEntries()>=kESDListN){
1091 for(int i = 0;i < fESDObjects->GetEntries() && i<kESDListN;i++){
1092 TObject *fObj = fESDObjects->At(i);
1093 if(fObj->InheritsFrom("TNamed")){
1094 ((TNamed*)fObj)->SetName(fgkESDListName[i]);
1096 else if(fObj->InheritsFrom("TClonesArray")){
1097 ((TClonesArray*)fObj)->SetName(fgkESDListName[i]);
1102 AliWarning("Std Entries missing");
1107 void AliESDEvent::CreateStdContent(Bool_t bUseThisList){
1108 fUseOwnList = bUseThisList;
1112 void AliESDEvent::CreateStdContent()
1114 // create the standard AOD content and set pointers
1116 // create standard objects and add them to the TList of objects
1117 AddObject(new AliESDRun());
1118 AddObject(new AliESDHeader());
1119 AddObject(new AliESDZDC());
1120 AddObject(new AliESDFMD());
1121 AddObject(new AliESDVZERO());
1122 AddObject(new AliESDTZERO());
1123 AddObject(new AliESDVertex());
1124 AddObject(new AliESDVertex());
1125 AddObject(new AliESDVertex());
1126 AddObject(new AliMultiplicity());
1127 AddObject(new AliESDCaloTrigger());
1128 AddObject(new AliESDCaloTrigger());
1129 AddObject(new TClonesArray("AliESDVertex",0));
1130 AddObject(new TClonesArray("AliESDVertex",0));
1131 AddObject(new TClonesArray("AliESDtrack",0));
1132 AddObject(new TClonesArray("AliESDMuonTrack",0));
1133 AddObject(new TClonesArray("AliESDPmdTrack",0));
1134 AddObject(new TClonesArray("AliESDTrdTrack",0));
1135 AddObject(new TClonesArray("AliESDv0",0));
1136 AddObject(new TClonesArray("AliESDcascade",0));
1137 AddObject(new TClonesArray("AliESDkink",0));
1138 AddObject(new TClonesArray("AliESDCaloCluster",0));
1139 AddObject(new AliESDCaloCells());
1140 AddObject(new AliESDCaloCells());
1141 AddObject(new TClonesArray("AliRawDataErrorLog",0));
1142 AddObject(new AliESDACORDE());
1144 // check the order of the indices against enum...
1148 // read back pointers
1152 TObject* AliESDEvent::FindListObject(const char *name){
1154 // Find object with name "name" in the list of branches
1157 return fESDObjects->FindObject(name);
1162 Int_t AliESDEvent::GetPHOSClusters(TRefArray *clusters) const
1164 // fills the provided TRefArray with all found phos clusters
1168 AliESDCaloCluster *cl = 0;
1169 for (Int_t i = 0; i < GetNumberOfCaloClusters(); i++) {
1171 if ( (cl = GetCaloCluster(i)) ) {
1174 AliDebug(1,Form("IsPHOS cluster %d Size: %d \n",i,clusters->GetEntriesFast()));
1178 return clusters->GetEntriesFast();
1181 Int_t AliESDEvent::GetEMCALClusters(TRefArray *clusters) const
1183 // fills the provided TRefArray with all found emcal clusters
1187 AliESDCaloCluster *cl = 0;
1188 for (Int_t i = 0; i < GetNumberOfCaloClusters(); i++) {
1190 if ( (cl = GetCaloCluster(i)) ) {
1193 AliDebug(1,Form("IsEMCAL cluster %d Size: %d \n",i,clusters->GetEntriesFast()));
1197 return clusters->GetEntriesFast();
1200 void AliESDEvent::WriteToTree(TTree* tree) const {
1201 // Book the branches as in TTree::Branch(TCollection*)
1202 // but add a "." at the end of top level branches which are
1203 // not a TClonesArray
1207 TIter next(fESDObjects);
1208 const Int_t kSplitlevel = 99; // default value in TTree::Branch()
1209 const Int_t kBufsize = 32000; // default value in TTree::Branch()
1212 while ((obj = next())) {
1213 branchname.Form("%s", obj->GetName());
1214 if(branchname.CompareTo("AliESDfriend")==0)branchname = "ESDfriend.";
1215 if ((kSplitlevel > 1) && !obj->InheritsFrom(TClonesArray::Class())) {
1216 if(!branchname.EndsWith("."))branchname += ".";
1218 if (!tree->FindBranch(branchname)) {
1219 tree->Bronch(branchname, obj->ClassName(), fESDObjects->GetObjectRef(obj),
1220 kBufsize, kSplitlevel - 1);
1226 void AliESDEvent::ReadFromTree(TTree *tree, Option_t* /*opt*/){
1228 // Connect the ESDEvent to a tree
1231 AliWarning("AliESDEvent::ReadFromTree() Zero Pointer to Tree \n");
1235 if(!tree->GetTree())tree->LoadTree(0);
1237 // if we find the "ESD" branch on the tree we do have the old structure
1238 if(tree->GetBranch("ESD")) {
1239 char ** address = (char **)(tree->GetBranch("ESD")->GetAddress());
1240 // do we have the friend branch
1241 TBranch * esdFB = tree->GetBranch("ESDfriend.");
1242 char ** addressF = 0;
1243 if(esdFB)addressF = (char **)(esdFB->GetAddress());
1245 AliInfo("AliESDEvent::ReadFromTree() Reading old Tree");
1246 tree->SetBranchAddress("ESD", &fESDOld);
1248 tree->SetBranchAddress("ESDfriend.",&fESDFriendOld);
1251 AliInfo("AliESDEvent::ReadFromTree() Reading old Tree");
1252 AliInfo("Branch already connected. Using existing branch address.");
1253 fESDOld = (AliESD*) (*address);
1254 // addressF can still be 0, since branch needs to switched on
1255 if(addressF)fESDFriendOld = (AliESDfriend*) (*addressF);
1258 // have already connected the old ESD structure... ?
1259 // reuse also the pointer of the AlliESDEvent
1260 // otherwise create new ones
1261 TList* connectedList = (TList*) (tree->GetUserInfo()->FindObject("ESDObjectsConnectedToTree"));
1264 // If connected use the connected list of objects
1265 if(fESDObjects!= connectedList){
1266 // protect when called twice
1267 fESDObjects->Delete();
1268 fESDObjects = connectedList;
1273 // The pointer to the friend changes when called twice via InitIO
1274 // since AliESDEvent is deleted
1275 TObject* oldf = FindListObject("AliESDfriend");
1278 newf = (TObject*)*addressF;
1280 if(newf!=0&&oldf!=newf){
1281 // remove the old reference
1282 // Should we also delete it? Or is this handled in TTree I/O
1283 // since it is created by the first SetBranchAddress
1284 fESDObjects->Remove(oldf);
1286 fESDObjects->Add(newf);
1293 CreateStdContent(); // create for copy
1294 // if we have the esdfriend add it, so we always can access it via the userinfo
1295 if(fESDFriendOld)AddObject(fESDFriendOld);
1296 // we are not owner of the list objects
1297 // must not delete it
1298 fESDObjects->SetOwner(kFALSE);
1299 fESDObjects->SetName("ESDObjectsConnectedToTree");
1300 tree->GetUserInfo()->Add(fESDObjects);
1308 // Try to find AliESDEvent
1309 AliESDEvent *esdEvent = 0;
1310 esdEvent = (AliESDEvent*)tree->GetTree()->GetUserInfo()->FindObject("AliESDEvent");
1312 // Check if already connected to tree
1313 TList* connectedList = (TList*) (tree->GetUserInfo()->FindObject("ESDObjectsConnectedToTree"));
1314 if (connectedList) {
1315 // If connected use the connected list if objects
1316 fESDObjects->Delete();
1317 fESDObjects = connectedList;
1324 // prevent a memory leak when reading back the TList
1329 // create a new TList from the UserInfo TList...
1330 // copy constructor does not work...
1331 fESDObjects = (TList*)(esdEvent->GetList()->Clone());
1332 fESDObjects->SetOwner(kFALSE);
1334 else if ( fESDObjects->GetEntries()==0){
1335 // at least create the std content if we want to read to our list
1340 // we only need new things in the list if we do no already have it..
1341 // TODO just add new entries
1343 if(fESDObjects->GetEntries()<kESDListN){
1344 AliWarning(Form("AliESDEvent::ReadFromTree() TList contains less than the standard contents %d < %d \n",
1345 fESDObjects->GetEntries(),kESDListN));
1347 // set the branch addresses
1348 TIter next(fESDObjects);
1350 while((el=(TNamed*)next())){
1351 TString bname(el->GetName());
1352 if(bname.CompareTo("AliESDfriend")==0)
1354 // AliESDfriend does not have a name ...
1355 tree->SetBranchAddress("ESDfriend.",fESDObjects->GetObjectRef(el));
1358 // check if branch exists under this Name
1359 TBranch *br = tree->GetBranch(bname.Data());
1361 tree->SetBranchAddress(bname.Data(),fESDObjects->GetObjectRef(el));
1364 br = tree->GetBranch(Form("%s.",bname.Data()));
1366 tree->SetBranchAddress(Form("%s.",bname.Data()),fESDObjects->GetObjectRef(el));
1369 AliWarning(Form("AliESDEvent::ReadFromTree() No Branch found with Name %s or %s.",bname.Data(),bname.Data()));
1376 // when reading back we are not owner of the list
1377 // must not delete it
1378 fESDObjects->SetOwner(kFALSE);
1379 fESDObjects->SetName("ESDObjectsConnectedToTree");
1380 // we are not owner of the list objects
1381 // must not delete it
1382 tree->GetUserInfo()->Add(fESDObjects);
1386 // we can't get the list from the user data, create standard content
1387 // and set it by hand (no ESDfriend at the moment
1389 TIter next(fESDObjects);
1391 while((el=(TNamed*)next())){
1392 TString bname(el->GetName());
1393 TBranch *br = tree->GetBranch(bname.Data());
1395 tree->SetBranchAddress(bname.Data(),fESDObjects->GetObjectRef(el));
1398 br = tree->GetBranch(Form("%s.",bname.Data()));
1400 tree->SetBranchAddress(Form("%s.",bname.Data()),fESDObjects->GetObjectRef(el));
1405 // when reading back we are not owner of the list
1406 // must not delete it
1407 fESDObjects->SetOwner(kFALSE);
1412 void AliESDEvent::CopyFromOldESD()
1414 // Method which copies over everthing from the old esd structure to the
1419 SetRunNumber(fESDOld->GetRunNumber());
1420 SetPeriodNumber(fESDOld->GetPeriodNumber());
1421 SetMagneticField(fESDOld->GetMagneticField());
1423 // leave out diamond ...
1424 // SetDiamond(const AliESDVertex *vertex) { fESDRun->SetDiamond(vertex);}
1427 SetTriggerMask(fESDOld->GetTriggerMask());
1428 SetOrbitNumber(fESDOld->GetOrbitNumber());
1429 SetTimeStamp(fESDOld->GetTimeStamp());
1430 SetEventType(fESDOld->GetEventType());
1431 SetEventNumberInFile(fESDOld->GetEventNumberInFile());
1432 SetBunchCrossNumber(fESDOld->GetBunchCrossNumber());
1433 SetTriggerCluster(fESDOld->GetTriggerCluster());
1437 SetZDC(fESDOld->GetZDCN1Energy(),
1438 fESDOld->GetZDCP1Energy(),
1439 fESDOld->GetZDCEMEnergy(),
1441 fESDOld->GetZDCN2Energy(),
1442 fESDOld->GetZDCP2Energy(),
1443 fESDOld->GetZDCParticipants(),
1448 if(fESDOld->GetFMDData())SetFMDData(fESDOld->GetFMDData());
1452 SetT0zVertex(fESDOld->GetT0zVertex());
1453 SetT0(fESDOld->GetT0());
1457 if (fESDOld->GetVZEROData()) SetVZEROData(fESDOld->GetVZEROData());
1459 if(fESDOld->GetVertex())SetPrimaryVertexSPD(fESDOld->GetVertex());
1461 if(fESDOld->GetPrimaryVertex())SetPrimaryVertexTracks(fESDOld->GetPrimaryVertex());
1463 if(fESDOld->GetMultiplicity())SetMultiplicity(fESDOld->GetMultiplicity());
1465 for(int i = 0;i<fESDOld->GetNumberOfTracks();i++){
1466 AddTrack(fESDOld->GetTrack(i));
1469 for(int i = 0;i<fESDOld->GetNumberOfMuonTracks();i++){
1470 AddMuonTrack(fESDOld->GetMuonTrack(i));
1473 for(int i = 0;i<fESDOld->GetNumberOfPmdTracks();i++){
1474 AddPmdTrack(fESDOld->GetPmdTrack(i));
1477 for(int i = 0;i<fESDOld->GetNumberOfTrdTracks();i++){
1478 AddTrdTrack(fESDOld->GetTrdTrack(i));
1481 for(int i = 0;i<fESDOld->GetNumberOfV0s();i++){
1482 AddV0(fESDOld->GetV0(i));
1485 for(int i = 0;i<fESDOld->GetNumberOfCascades();i++){
1486 AddCascade(fESDOld->GetCascade(i));
1489 for(int i = 0;i<fESDOld->GetNumberOfKinks();i++){
1490 AddKink(fESDOld->GetKink(i));
1494 for(int i = 0;i<fESDOld->GetNumberOfCaloClusters();i++){
1495 AddCaloCluster(fESDOld->GetCaloCluster(i));