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 <TInterpreter.h>
41 #include "AliESDEvent.h"
42 #include "AliESDfriend.h"
43 #include "AliESDVZERO.h"
44 #include "AliESDFMD.h"
46 #include "AliESDMuonTrack.h"
47 #include "AliESDPmdTrack.h"
48 #include "AliESDTrdTrack.h"
49 #include "AliESDVertex.h"
50 #include "AliESDcascade.h"
51 #include "AliESDPmdTrack.h"
52 #include "AliESDTrdTrigger.h"
53 #include "AliESDTrdTrack.h"
54 #include "AliESDTrdTracklet.h"
55 #include "AliESDVertex.h"
56 #include "AliVertexerTracks.h"
57 #include "AliESDcascade.h"
58 #include "AliESDkink.h"
59 #include "AliESDtrack.h"
60 #include "AliESDHLTtrack.h"
61 #include "AliESDCaloCluster.h"
62 #include "AliESDCaloCells.h"
64 #include "AliESDFMD.h"
65 #include "AliESDVZERO.h"
66 #include "AliMultiplicity.h"
67 #include "AliRawDataErrorLog.h"
69 #include "AliESDACORDE.h"
70 #include "AliESDHLTDecision.h"
71 #include "AliCentrality.h"
72 #include "AliESDCosmicTrack.h"
74 #include "AliESDMFT.h"
76 #include "AliEventplane.h"
83 // here we define the names, some classes are no TNamed, therefore the classnames
85 const char* AliESDEvent::fgkESDListName[kESDListN] = {"AliESDRun",
111 "AliRawDataErrorLogs",
120 //______________________________________________________________________________
121 AliESDEvent::AliESDEvent():
123 fESDObjects(new TList()),
138 fSPDPileupVertices(0),
139 fTrkPileupVertices(0),
149 fEMCALCells(0), fPHOSCells(0),
159 fDetectorStatus(0xFFFFFFFF)
165 //______________________________________________________________________________
166 AliESDEvent::AliESDEvent(const AliESDEvent& esd):
168 fESDObjects(new TList()),
169 fESDRun(new AliESDRun(*esd.fESDRun)),
170 fHeader(new AliESDHeader(*esd.fHeader)),
171 fESDZDC(new AliESDZDC(*esd.fESDZDC)),
172 fESDFMD(new AliESDFMD(*esd.fESDFMD)),
173 fESDVZERO(new AliESDVZERO(*esd.fESDVZERO)),
174 fESDTZERO(new AliESDTZERO(*esd.fESDTZERO)),
175 fTPCVertex(new AliESDVertex(*esd.fTPCVertex)),
176 fSPDVertex(new AliESDVertex(*esd.fSPDVertex)),
177 fPrimaryVertex(new AliESDVertex(*esd.fPrimaryVertex)),
178 fSPDMult(new AliMultiplicity(*esd.fSPDMult)),
179 fPHOSTrigger(new AliESDCaloTrigger(*esd.fPHOSTrigger)),
180 fEMCALTrigger(new AliESDCaloTrigger(*esd.fEMCALTrigger)),
181 fESDACORDE(new AliESDACORDE(*esd.fESDACORDE)),
182 fTrdTrigger(new AliESDTrdTrigger(*esd.fTrdTrigger)),
183 fSPDPileupVertices(new TClonesArray(*esd.fSPDPileupVertices)),
184 fTrkPileupVertices(new TClonesArray(*esd.fTrkPileupVertices)),
185 fTracks(new TClonesArray(*esd.fTracks)),
186 fMuonTracks(new TClonesArray(*esd.fMuonTracks)),
187 fPmdTracks(new TClonesArray(*esd.fPmdTracks)),
188 fTrdTracks(new TClonesArray(*esd.fTrdTracks)),
189 fTrdTracklets(new TClonesArray(*esd.fTrdTracklets)),
190 fV0s(new TClonesArray(*esd.fV0s)),
191 fCascades(new TClonesArray(*esd.fCascades)),
192 fKinks(new TClonesArray(*esd.fKinks)),
193 fCaloClusters(new TClonesArray(*esd.fCaloClusters)),
194 fEMCALCells(new AliESDCaloCells(*esd.fEMCALCells)),
195 fPHOSCells(new AliESDCaloCells(*esd.fPHOSCells)),
196 fCosmicTracks(new TClonesArray(*esd.fCosmicTracks)),
197 fErrorLogs(new TClonesArray(*esd.fErrorLogs)),
198 fESDOld(esd.fESDOld ? new AliESD(*esd.fESDOld) : 0),
199 fESDFriendOld(esd.fESDFriendOld ? new AliESDfriend(*esd.fESDFriendOld) : 0),
200 fConnected(esd.fConnected),
201 fUseOwnList(esd.fUseOwnList),
202 fTOFHeader(new AliTOFHeader(*esd.fTOFHeader)),
203 fCentrality(new AliCentrality(*esd.fCentrality)),
204 fEventplane(new AliEventplane(*esd.fEventplane)),
205 fDetectorStatus(esd.fDetectorStatus)
207 // , fESDMFT(new AliESDMFT(*esd.fESDMFT))
212 printf("copying ESD event...\n"); // AU
213 // CKB init in the constructor list and only add here ...
218 AddObject(fESDVZERO);
219 AddObject(fESDTZERO);
220 AddObject(fTPCVertex);
221 AddObject(fSPDVertex);
222 AddObject(fPrimaryVertex);
224 AddObject(fPHOSTrigger);
225 AddObject(fEMCALTrigger);
226 AddObject(fTrdTrigger);
227 AddObject(fSPDPileupVertices);
228 AddObject(fTrkPileupVertices);
230 AddObject(fMuonTracks);
231 AddObject(fPmdTracks);
232 AddObject(fTrdTracks);
233 AddObject(fTrdTracklets);
235 AddObject(fCascades);
237 AddObject(fCaloClusters);
238 AddObject(fEMCALCells);
239 AddObject(fPHOSCells);
240 AddObject(fCosmicTracks);
241 AddObject(fErrorLogs);
242 AddObject(fESDACORDE);
243 AddObject(fTOFHeader);
245 // AddObject(fESDMFT);
251 //______________________________________________________________________________
252 AliESDEvent & AliESDEvent::operator=(const AliESDEvent& source) {
254 // Assignment operator
256 if(&source == this) return *this;
257 AliVEvent::operator=(source);
259 // This assumes that the list is already created
260 // and that the virtual void Copy(Tobject&) function
261 // is correctly implemented in the derived class
262 // otherwise only TObject::Copy() will be used
266 if((fESDObjects->GetSize()==0)&&(source.fESDObjects->GetSize()>=kESDListN)){
267 // We cover the case that we do not yet have the
268 // standard content but the source has it
272 TIter next(source.GetList());
275 while ((its = next())) {
276 name.Form("%s", its->GetName());
277 TObject *mine = fESDObjects->FindObject(name.Data());
279 TClass* pClass=TClass::GetClass(its->ClassName());
281 AliWarning(Form("Can not find class description for entry %s (%s)\n",
282 its->ClassName(), name.Data()));
286 mine=(TObject*)pClass->New();
288 // not in this: can be added to list
289 AliWarning(Form("%s:%d Could not find %s for copying \n",
290 (char*)__FILE__,__LINE__,name.Data()));
293 if(mine->InheritsFrom("TNamed")){
294 ((TNamed*)mine)->SetName(name);
296 else if(mine->InheritsFrom("TCollection")){
297 if(mine->InheritsFrom("TClonesArray")) {
298 TClonesArray* tcits = dynamic_cast<TClonesArray*>(its);
300 dynamic_cast<TClonesArray*>(mine)->SetClass(tcits->GetClass());
302 dynamic_cast<TCollection*>(mine)->SetName(name);
304 AliDebug(1, Form("adding object %s of type %s", mine->GetName(), mine->ClassName()));
308 if(!its->InheritsFrom("TCollection")){
312 else if(its->InheritsFrom("TClonesArray")){
313 // Create or expand the tclonesarray pointers
314 // so we can directly copy to the object
315 TClonesArray *itstca = (TClonesArray*)its;
316 TClonesArray *minetca = (TClonesArray*)mine;
318 // this leaves the capacity of the TClonesArray the same
319 // except for a factor of 2 increase when size > capacity
320 // does not release any memory occupied by the tca
321 minetca->ExpandCreate(itstca->GetEntriesFast());
322 for(int i = 0;i < itstca->GetEntriesFast();++i){
324 TObject *minetcaobj = minetca->At(i);
325 TObject *itstcaobj = itstca->At(i);
326 // no need to delete first
327 // pointers within the class should be handled by Copy()...
328 // Can there be Empty slots?
329 itstcaobj->Copy(*minetcaobj);
333 AliWarning(Form("%s:%d cannot copy TCollection \n",
334 (char*)__FILE__,__LINE__));
338 fCentrality = source.fCentrality;
339 fEventplane = source.fEventplane;
341 fConnected = source.fConnected;
342 fUseOwnList = source.fUseOwnList;
344 fDetectorStatus = source.fDetectorStatus;
350 //______________________________________________________________________________
351 AliESDEvent::~AliESDEvent()
354 // Standard destructor
357 // everthing on the list gets deleted automatically
360 if(fESDObjects&&!fConnected)
365 if (fCentrality) delete fCentrality;
366 if (fEventplane) delete fEventplane;
370 void AliESDEvent::Copy(TObject &obj) const {
372 // interface to TOBject::Copy
373 // Copies the content of this into obj!
374 // bascially obj = *this
376 if(this==&obj)return;
377 AliESDEvent *robj = dynamic_cast<AliESDEvent*>(&obj);
378 if(!robj)return; // not an AliESEvent
383 //______________________________________________________________________________
384 void AliESDEvent::Reset()
388 // Std content + Non std content
390 // Reset the standard contents
393 // reset for the old data without AliESDEvent...
394 if(fESDOld)fESDOld->Reset();
396 fESDFriendOld->~AliESDfriend();
397 new (fESDFriendOld) AliESDfriend();
401 if(fESDObjects->GetSize()>kESDListN){
402 // we have non std content
403 // this also covers esdfriends
404 for(int i = kESDListN;i < fESDObjects->GetSize();++i){
405 TObject *pObject = fESDObjects->At(i);
407 if(pObject->InheritsFrom(TClonesArray::Class())){
408 ((TClonesArray*)pObject)->Delete();
410 else if(!pObject->InheritsFrom(TCollection::Class())){
411 TClass *pClass = TClass::GetClass(pObject->ClassName());
412 if (pClass && pClass->GetListOfMethods()->FindObject("Clear")) {
413 AliDebug(1, Form("Clear for object %s class %s", pObject->GetName(), pObject->ClassName()));
417 AliDebug(1, Form("ResetWithPlacementNew for object %s class %s", pObject->GetName(), pObject->ClassName()));
418 ResetWithPlacementNew(pObject);
422 AliWarning(Form("No reset for %s \n",
423 pObject->ClassName()));
430 Bool_t AliESDEvent::ResetWithPlacementNew(TObject *pObject){
432 // funtion to reset using the already allocated space
434 Long_t dtoronly = TObject::GetDtorOnly();
435 TClass *pClass = TClass::GetClass(pObject->ClassName());
436 TObject::SetDtorOnly(pObject);
438 // Recreate with placement new
439 pClass->New(pObject);
440 // Restore the state.
441 TObject::SetDtorOnly((void*)dtoronly);
445 void AliESDEvent::ResetStdContent()
447 // Reset the standard contents
448 if(fESDRun) fESDRun->Reset();
449 if(fHeader) fHeader->Reset();
450 if(fCentrality) fCentrality->Reset();
451 if(fEventplane) fEventplane->Reset();
452 if(fESDZDC) fESDZDC->Reset();
457 // reset by callin d'to /c'tor keep the pointer
458 fESDVZERO->~AliESDVZERO();
459 new (fESDVZERO) AliESDVZERO();
462 fESDACORDE->~AliESDACORDE();
463 new (fESDACORDE) AliESDACORDE();
465 if(fESDTZERO) fESDTZERO->Reset();
466 // CKB no clear/reset implemented
468 fTPCVertex->~AliESDVertex();
469 new (fTPCVertex) AliESDVertex();
470 fTPCVertex->SetName(fgkESDListName[kTPCVertex]);
473 fSPDVertex->~AliESDVertex();
474 new (fSPDVertex) AliESDVertex();
475 fSPDVertex->SetName(fgkESDListName[kSPDVertex]);
478 fPrimaryVertex->~AliESDVertex();
479 new (fPrimaryVertex) AliESDVertex();
480 fPrimaryVertex->SetName(fgkESDListName[kPrimaryVertex]);
483 fSPDMult->~AliMultiplicity();
484 new (fSPDMult) AliMultiplicity();
487 fTOFHeader->~AliTOFHeader();
488 new (fTOFHeader) AliTOFHeader();
489 //fTOFHeader->SetName(fgkESDListName[kTOFHeader]);
492 fTrdTrigger->~AliESDTrdTrigger();
493 new (fTrdTrigger) AliESDTrdTrigger();
497 // fESDMFT->~AliESDMFT();
498 // new (fESDMFT) AliESDMFT();
502 if(fPHOSTrigger)fPHOSTrigger->DeAllocate();
503 if(fEMCALTrigger)fEMCALTrigger->DeAllocate();
504 if(fSPDPileupVertices)fSPDPileupVertices->Delete();
505 if(fTrkPileupVertices)fTrkPileupVertices->Delete();
506 if(fTracks)fTracks->Clear("C"); //Delete();
507 if(fMuonTracks)fMuonTracks->Clear("C"); //Delete();
508 if(fPmdTracks)fPmdTracks->Delete();
509 if(fTrdTracks)fTrdTracks->Delete();
510 if(fTrdTracklets)fTrdTracklets->Delete();
511 if(fV0s)fV0s->Delete();
512 if(fCascades)fCascades->Delete();
513 if(fKinks)fKinks->Delete();
514 if(fCaloClusters)fCaloClusters->Delete();
515 if(fPHOSCells)fPHOSCells->DeleteContainer();
516 if(fEMCALCells)fEMCALCells->DeleteContainer();
517 if(fCosmicTracks)fCosmicTracks->Delete();
518 if(fErrorLogs) fErrorLogs->Delete();
520 // don't reset fconnected fConnected and the list
525 Int_t AliESDEvent::AddV0(const AliESDv0 *v) {
529 TClonesArray &fv = *fV0s;
530 Int_t idx=fV0s->GetEntriesFast();
531 new(fv[idx]) AliESDv0(*v);
535 //______________________________________________________________________________
536 void AliESDEvent::Print(Option_t *) const
539 // Print header information of the event
541 printf("ESD run information\n");
542 printf("Event # in file %d Bunch crossing # %d Orbit # %d Period # %d Run # %d Trigger %lld Magnetic field %f \n",
543 GetEventNumberInFile(),
544 GetBunchCrossNumber(),
549 GetMagneticField() );
551 printf("Vertex: (%.4f +- %.4f, %.4f +- %.4f, %.4f +- %.4f) cm\n",
552 fPrimaryVertex->GetXv(), fPrimaryVertex->GetXRes(),
553 fPrimaryVertex->GetYv(), fPrimaryVertex->GetYRes(),
554 fPrimaryVertex->GetZv(), fPrimaryVertex->GetZRes());
555 printf("Mean vertex in RUN: X=%.4f Y=%.4f Z=%.4f cm\n",
556 GetDiamondX(),GetDiamondY(),GetDiamondZ());
558 printf("SPD Multiplicity. Number of tracklets %d \n",
559 fSPDMult->GetNumberOfTracklets());
560 printf("Number of pileup primary vertices reconstructed with SPD %d\n",
561 GetNumberOfPileupVerticesSPD());
562 printf("Number of pileup primary vertices reconstructed using the tracks %d\n",
563 GetNumberOfPileupVerticesTracks());
564 printf("Number of tracks: \n");
565 printf(" charged %d\n", GetNumberOfTracks());
566 printf(" muon %d\n", GetNumberOfMuonTracks());
567 printf(" pmd %d\n", GetNumberOfPmdTracks());
568 printf(" trd %d\n", GetNumberOfTrdTracks());
569 printf(" trd trkl %d\n", GetNumberOfTrdTracklets());
570 printf(" v0 %d\n", GetNumberOfV0s());
571 printf(" cascades %d\n", GetNumberOfCascades());
572 printf(" kinks %d\n", GetNumberOfKinks());
573 if(fPHOSCells)printf(" PHOSCells %d\n", fPHOSCells->GetNumberOfCells());
574 else printf(" PHOSCells not in the Event\n");
575 if(fEMCALCells)printf(" EMCALCells %d\n", fEMCALCells->GetNumberOfCells());
576 else printf(" EMCALCells not in the Event\n");
577 printf(" CaloClusters %d\n", GetNumberOfCaloClusters());
578 printf(" FMD %s\n", (fESDFMD ? "yes" : "no"));
579 printf(" VZERO %s\n", (fESDVZERO ? "yes" : "no"));
580 if (fCosmicTracks) printf(" Cosmics %d\n", GetNumberOfCosmicTracks());
582 //printf(" MFT %s\n", (fESDMFT ? "yes" : "no"));
586 TObject* pHLTDecision=GetHLTTriggerDecision();
587 printf("HLT trigger decision: %s\n", pHLTDecision?pHLTDecision->GetOption():"not available");
588 if (pHLTDecision) pHLTDecision->Print("compact");
593 void AliESDEvent::SetESDfriend(const AliESDfriend *ev) const {
595 // Attaches the complementary info to the ESD
599 // to be sure that we set the tracks also
600 // in case of old esds
601 // if(fESDOld)CopyFromOldESD();
603 Int_t ntrk=ev->GetNumberOfTracks();
605 for (Int_t i=0; i<ntrk; i++) {
606 const AliESDfriendTrack *f=ev->GetTrack(i);
607 if (!f) {AliFatal(Form("NULL pointer for ESD track %d",i));}
608 GetTrack(i)->SetFriendTrack(f);
612 Bool_t AliESDEvent::RemoveKink(Int_t rm) const {
613 // ---------------------------------------------------------
614 // Remove a kink candidate and references to it from ESD,
615 // if this candidate does not come from a reconstructed decay
616 // Not yet implemented...
617 // ---------------------------------------------------------
618 Int_t last=GetNumberOfKinks()-1;
619 if ((rm<0)||(rm>last)) return kFALSE;
624 Bool_t AliESDEvent::RemoveV0(Int_t rm) const {
625 // ---------------------------------------------------------
626 // Remove a V0 candidate and references to it from ESD,
627 // if this candidate does not come from a reconstructed decay
628 // ---------------------------------------------------------
629 Int_t last=GetNumberOfV0s()-1;
630 if ((rm<0)||(rm>last)) return kFALSE;
632 AliESDv0 *v0=GetV0(rm);
633 Int_t idxP=v0->GetPindex(), idxN=v0->GetNindex();
636 Int_t lastIdxP=v0->GetPindex(), lastIdxN=v0->GetNindex();
640 // Check if this V0 comes from a reconstructed decay
641 Int_t ncs=GetNumberOfCascades();
642 for (Int_t n=0; n<ncs; n++) {
643 AliESDcascade *cs=GetCascade(n);
645 Int_t csIdxP=cs->GetPindex();
646 Int_t csIdxN=cs->GetNindex();
649 if (idxN==csIdxN) return kFALSE;
651 if (csIdxP==lastIdxP)
652 if (csIdxN==lastIdxN) used++;
655 //Replace the removed V0 with the last V0
656 TClonesArray &a=*fV0s;
657 delete a.RemoveAt(rm);
659 if (rm==last) return kTRUE;
661 //v0 is pointing to the last V0 candidate...
662 new (a[rm]) AliESDv0(*v0);
663 delete a.RemoveAt(last);
665 if (!used) return kTRUE;
668 // Remap the indices of the daughters of reconstructed decays
669 for (Int_t n=0; n<ncs; n++) {
670 AliESDcascade *cs=GetCascade(n);
673 Int_t csIdxP=cs->GetPindex();
674 Int_t csIdxN=cs->GetNindex();
676 if (csIdxP==lastIdxP)
677 if (csIdxN==lastIdxN) {
678 cs->AliESDv0::SetIndex(1,idxP);
679 cs->AliESDv0::SetIndex(0,idxN);
681 if (!used) return kTRUE;
688 Bool_t AliESDEvent::RemoveTrack(Int_t rm) const {
689 // ---------------------------------------------------------
690 // Remove a track and references to it from ESD,
691 // if this track does not come from a reconstructed decay
692 // ---------------------------------------------------------
693 Int_t last=GetNumberOfTracks()-1;
694 if ((rm<0)||(rm>last)) return kFALSE;
698 // Check if this track comes from the reconstructed primary vertices
699 if (fTPCVertex && fTPCVertex->GetStatus()) {
700 UShort_t *primIdx=fTPCVertex->GetIndices();
701 Int_t n=fTPCVertex->GetNIndices();
703 Int_t idx=Int_t(primIdx[n]);
704 if (rm==idx) return kFALSE;
705 if (idx==last) used++;
708 if (fPrimaryVertex && fPrimaryVertex->GetStatus()) {
709 UShort_t *primIdx=fPrimaryVertex->GetIndices();
710 Int_t n=fPrimaryVertex->GetNIndices();
712 Int_t idx=Int_t(primIdx[n]);
713 if (rm==idx) return kFALSE;
714 if (idx==last) used++;
718 // Check if this track comes from a reconstructed decay
719 Int_t nv0=GetNumberOfV0s();
720 for (Int_t n=0; n<nv0; n++) {
721 AliESDv0 *v0=GetV0(n);
723 Int_t idx=v0->GetNindex();
724 if (rm==idx) return kFALSE;
725 if (idx==last) used++;
728 if (rm==idx) return kFALSE;
729 if (idx==last) used++;
732 Int_t ncs=GetNumberOfCascades();
733 for (Int_t n=0; n<ncs; n++) {
734 AliESDcascade *cs=GetCascade(n);
736 Int_t idx=cs->GetIndex();
737 if (rm==idx) return kFALSE;
738 if (idx==last) used++;
742 if (rm==idx) return kFALSE;
743 if (idx==last) used++;
746 if (rm==idx) return kFALSE;
747 if (idx==last) used++;
750 Int_t nkn=GetNumberOfKinks();
751 for (Int_t n=0; n<nkn; n++) {
752 AliESDkink *kn=GetKink(n);
754 Int_t idx=kn->GetIndex(0);
755 if (rm==idx) return kFALSE;
756 if (idx==last) used++;
759 if (rm==idx) return kFALSE;
760 if (idx==last) used++;
763 // Check if this track is associated with a CaloCluster
764 Int_t ncl=GetNumberOfCaloClusters();
765 for (Int_t n=0; n<ncl; n++) {
766 AliESDCaloCluster *cluster=GetCaloCluster(n);
767 TArrayI *arr=cluster->GetTracksMatched();
768 Int_t s=arr->GetSize();
770 Int_t idx=arr->At(s);
771 if (rm==idx) return kFALSE;
772 if (idx==last) used++;
778 //Replace the removed track with the last track
779 TClonesArray &a=*fTracks;
780 delete a.RemoveAt(rm);
782 if (rm==last) return kTRUE;
784 AliESDtrack *t=GetTrack(last);
785 if (!t) {AliFatal(Form("NULL pointer for ESD track %d",last));}
787 new (a[rm]) AliESDtrack(*t);
788 delete a.RemoveAt(last);
791 if (!used) return kTRUE;
794 // Remap the indices of the tracks used for the primary vertex reconstruction
795 if (fTPCVertex && fTPCVertex->GetStatus()) {
796 UShort_t *primIdx=fTPCVertex->GetIndices();
797 Int_t n=fTPCVertex->GetNIndices();
799 Int_t idx=Int_t(primIdx[n]);
801 primIdx[n]=Short_t(rm);
803 if (!used) return kTRUE;
807 if (fPrimaryVertex && fPrimaryVertex->GetStatus()) {
808 UShort_t *primIdx=fPrimaryVertex->GetIndices();
809 Int_t n=fPrimaryVertex->GetNIndices();
811 Int_t idx=Int_t(primIdx[n]);
813 primIdx[n]=Short_t(rm);
815 if (!used) return kTRUE;
820 // Remap the indices of the daughters of reconstructed decays
821 for (Int_t n=0; n<nv0; n++) {
822 AliESDv0 *v0=GetV0(n);
823 if (v0->GetIndex(0)==last) {
826 if (!used) return kTRUE;
828 if (v0->GetIndex(1)==last) {
831 if (!used) return kTRUE;
835 for (Int_t n=0; n<ncs; n++) {
836 AliESDcascade *cs=GetCascade(n);
837 if (cs->GetIndex()==last) {
840 if (!used) return kTRUE;
843 if (v0->GetIndex(0)==last) {
846 if (!used) return kTRUE;
848 if (v0->GetIndex(1)==last) {
851 if (!used) return kTRUE;
855 for (Int_t n=0; n<nkn; n++) {
856 AliESDkink *kn=GetKink(n);
857 if (kn->GetIndex(0)==last) {
860 if (!used) return kTRUE;
862 if (kn->GetIndex(1)==last) {
865 if (!used) return kTRUE;
869 // Remap the indices of the tracks accosicated with CaloClusters
870 for (Int_t n=0; n<ncl; n++) {
871 AliESDCaloCluster *cluster=GetCaloCluster(n);
872 TArrayI *arr=cluster->GetTracksMatched();
873 Int_t s=arr->GetSize();
875 Int_t idx=arr->At(s);
879 if (!used) return kTRUE;
888 Bool_t AliESDEvent::Clean(Float_t *cleanPars) {
890 // Remove the data which are not needed for the physics analysis.
892 // 1) Cleaning the V0 candidates
893 // ---------------------------
894 // If the cosine of the V0 pointing angle "csp" and
895 // the DCA between the daughter tracks "dca" does not satisfy
898 // csp > cleanPars[1] + dca/cleanPars[0]*(1.- cleanPars[1])
900 // an attempt to remove this V0 candidate from ESD is made.
902 // The V0 candidate gets removed if it does not belong to any
903 // recosntructed cascade decay
905 // 12.11.2007, optimal values: cleanPars[0]=0.5, cleanPars[1]=0.999
907 // 2) Cleaning the tracks
908 // ----------------------
909 // If track's transverse parameter is larger than cleanPars[2]
911 // track's longitudinal parameter is larger than cleanPars[3]
912 // an attempt to remove this track from ESD is made.
914 // The track gets removed if it does not come
915 // from a reconstructed decay
919 Float_t dcaMax=cleanPars[0];
920 Float_t cspMin=cleanPars[1];
922 Int_t nV0s=GetNumberOfV0s();
923 for (Int_t i=nV0s-1; i>=0; i--) {
924 AliESDv0 *v0=GetV0(i);
926 Float_t dca=v0->GetDcaV0Daughters();
927 Float_t csp=v0->GetV0CosineOfPointingAngle();
928 Float_t cspcut=cspMin + dca/dcaMax*(1.-cspMin);
929 if (csp > cspcut) continue;
930 if (RemoveV0(i)) rc=kTRUE;
934 Float_t dmax=cleanPars[2], zmax=cleanPars[3];
936 const AliESDVertex *vertex=GetPrimaryVertexSPD();
937 Bool_t vtxOK=vertex->GetStatus();
939 Int_t nTracks=GetNumberOfTracks();
940 for (Int_t i=nTracks-1; i>=0; i--) {
941 AliESDtrack *track=GetTrack(i);
942 if (!track) {AliFatal(Form("NULL pointer for ESD track %d",i));}
943 Float_t xy,z; track->GetImpactParameters(xy,z);
944 if ((TMath::Abs(xy) > dmax) || (vtxOK && (TMath::Abs(z) > zmax))) {
945 if (RemoveTrack(i)) rc=kTRUE;
952 Char_t AliESDEvent::AddPileupVertexSPD(const AliESDVertex *vtx)
954 // Add a pileup primary vertex reconstructed with SPD
955 TClonesArray &ftr = *fSPDPileupVertices;
956 Char_t n=Char_t(ftr.GetEntriesFast());
957 AliESDVertex *vertex = new(ftr[n]) AliESDVertex(*vtx);
962 Char_t AliESDEvent::AddPileupVertexTracks(const AliESDVertex *vtx)
964 // Add a pileup primary vertex reconstructed with SPD
965 TClonesArray &ftr = *fTrkPileupVertices;
966 Char_t n=Char_t(ftr.GetEntriesFast());
967 AliESDVertex *vertex = new(ftr[n]) AliESDVertex(*vtx);
972 Int_t AliESDEvent::AddTrack(const AliESDtrack *t)
975 TClonesArray &ftr = *fTracks;
976 AliESDtrack * track = new(ftr[fTracks->GetEntriesFast()])AliESDtrack(*t);
977 track->SetID(fTracks->GetEntriesFast()-1);
978 return track->GetID();
981 AliESDtrack* AliESDEvent::NewTrack()
984 TClonesArray &ftr = *fTracks;
985 AliESDtrack * track = new(ftr[fTracks->GetEntriesFast()])AliESDtrack();
986 track->SetID(fTracks->GetEntriesFast()-1);
990 void AliESDEvent::AddMuonTrack(const AliESDMuonTrack *t)
992 TClonesArray &fmu = *fMuonTracks;
993 new(fmu[fMuonTracks->GetEntriesFast()]) AliESDMuonTrack(*t);
996 void AliESDEvent::AddPmdTrack(const AliESDPmdTrack *t)
998 TClonesArray &fpmd = *fPmdTracks;
999 new(fpmd[fPmdTracks->GetEntriesFast()]) AliESDPmdTrack(*t);
1002 void AliESDEvent::SetTrdTrigger(const AliESDTrdTrigger *t)
1007 void AliESDEvent::AddTrdTrack(const AliESDTrdTrack *t)
1009 TClonesArray &ftrd = *fTrdTracks;
1010 new(ftrd[fTrdTracks->GetEntriesFast()]) AliESDTrdTrack(*t);
1013 void AliESDEvent::AddTrdTracklet(const AliESDTrdTracklet *trkl)
1015 new ((*fTrdTracklets)[fTrdTracklets->GetEntriesFast()]) AliESDTrdTracklet(*trkl);
1018 Int_t AliESDEvent::AddKink(const AliESDkink *c)
1021 TClonesArray &fk = *fKinks;
1022 AliESDkink * kink = new(fk[fKinks->GetEntriesFast()]) AliESDkink(*c);
1023 kink->SetID(fKinks->GetEntriesFast()); // CKB different from the other imps..
1024 return fKinks->GetEntriesFast()-1;
1028 void AliESDEvent::AddCascade(const AliESDcascade *c)
1030 TClonesArray &fc = *fCascades;
1031 new(fc[fCascades->GetEntriesFast()]) AliESDcascade(*c);
1034 void AliESDEvent::AddCosmicTrack(const AliESDCosmicTrack *t)
1036 TClonesArray &ft = *fCosmicTracks;
1037 new(ft[fCosmicTracks->GetEntriesFast()]) AliESDCosmicTrack(*t);
1041 Int_t AliESDEvent::AddCaloCluster(const AliESDCaloCluster *c)
1044 TClonesArray &fc = *fCaloClusters;
1045 AliESDCaloCluster *clus = new(fc[fCaloClusters->GetEntriesFast()]) AliESDCaloCluster(*c);
1046 clus->SetID(fCaloClusters->GetEntriesFast()-1);
1047 return fCaloClusters->GetEntriesFast()-1;
1051 void AliESDEvent::AddRawDataErrorLog(const AliRawDataErrorLog *log) const {
1052 TClonesArray &errlogs = *fErrorLogs;
1053 new(errlogs[errlogs.GetEntriesFast()]) AliRawDataErrorLog(*log);
1056 void AliESDEvent::SetZDCData(const AliESDZDC * obj)
1058 // use already allocated space
1063 void AliESDEvent::SetPrimaryVertexTPC(const AliESDVertex *vertex)
1065 // Set the TPC vertex
1066 // use already allocated space
1068 *fTPCVertex = *vertex;
1069 fTPCVertex->SetName(fgkESDListName[kTPCVertex]);
1073 void AliESDEvent::SetPrimaryVertexSPD(const AliESDVertex *vertex)
1075 // Set the SPD vertex
1076 // use already allocated space
1078 *fSPDVertex = *vertex;
1079 fSPDVertex->SetName(fgkESDListName[kSPDVertex]);
1083 void AliESDEvent::SetPrimaryVertexTracks(const AliESDVertex *vertex)
1085 // Set the primary vertex reconstructed using he ESD tracks.
1086 // use already allocated space
1088 *fPrimaryVertex = *vertex;
1089 fPrimaryVertex->SetName(fgkESDListName[kPrimaryVertex]);
1093 const AliESDVertex * AliESDEvent::GetPrimaryVertex() const
1096 // Get the "best" available reconstructed primary vertex.
1099 if (fPrimaryVertex->GetStatus()) return fPrimaryVertex;
1102 if (fSPDVertex->GetStatus()) return fSPDVertex;
1104 if(fTPCVertex) return fTPCVertex;
1106 AliWarning("No primary vertex available. Returning the \"default\"...");
1110 AliESDVertex * AliESDEvent::PrimaryVertexTracksUnconstrained() const
1113 // Removes diamond constraint from fPrimaryVertex (reconstructed with tracks)
1114 // Returns a AliESDVertex which has to be deleted by the user
1116 if(!fPrimaryVertex) {
1117 AliWarning("No primary vertex from tracks available.");
1120 if(!fPrimaryVertex->GetStatus()) {
1121 AliWarning("No primary vertex from tracks available.");
1125 AliVertexerTracks vertexer(GetMagneticField());
1126 Float_t diamondxyz[3]={(Float_t)GetDiamondX(),(Float_t)GetDiamondY(),0.};
1127 Float_t diamondcovxy[3]; GetDiamondCovXY(diamondcovxy);
1128 Float_t diamondcov[6]={diamondcovxy[0],diamondcovxy[1],diamondcovxy[2],0.,0.,7.};
1129 AliESDVertex *vertex =
1130 (AliESDVertex*)vertexer.RemoveConstraintFromVertex(fPrimaryVertex,diamondxyz,diamondcov);
1135 void AliESDEvent::SetMultiplicity(const AliMultiplicity *mul)
1137 // Set the SPD Multiplicity
1144 void AliESDEvent::SetFMDData(AliESDFMD * obj)
1146 // use already allocated space
1152 void AliESDEvent::SetVZEROData(const AliESDVZERO * obj)
1154 // use already allocated space
1159 void AliESDEvent::SetTZEROData(const AliESDTZERO * obj)
1161 // use already allocated space
1167 //void AliESDEvent::SetMFTData(AliESDMFT * obj)
1174 void AliESDEvent::SetACORDEData(AliESDACORDE * obj)
1181 void AliESDEvent::GetESDfriend(AliESDfriend *ev) const
1184 // Extracts the complementary info from the ESD
1188 Int_t ntrk=GetNumberOfTracks();
1190 for (Int_t i=0; i<ntrk; i++) {
1191 AliESDtrack *t=GetTrack(i);
1192 if (!t) {AliFatal(Form("NULL pointer for ESD track %d",i));}
1193 const AliESDfriendTrack *f=t->GetFriendTrack();
1196 t->ReleaseESDfriendTrack();// Not to have two copies of "friendTrack"
1200 AliESDfriend *fr = (AliESDfriend*)(const_cast<AliESDEvent*>(this)->FindListObject("AliESDfriend"));
1201 if (fr) ev->SetVZEROfriend(fr->GetVZEROfriend());
1204 void AliESDEvent::AddObject(TObject* obj)
1206 // Add an object to the list of object.
1207 // Please be aware that in order to increase performance you should
1208 // refrain from using TObjArrays (if possible). Use TClonesArrays, instead.
1209 fESDObjects->SetOwner(kTRUE);
1210 fESDObjects->AddLast(obj);
1214 void AliESDEvent::GetStdContent()
1216 // set pointers for standard content
1217 // get by name much safer and not a big overhead since not called very often
1219 fESDRun = (AliESDRun*)fESDObjects->FindObject(fgkESDListName[kESDRun]);
1220 fHeader = (AliESDHeader*)fESDObjects->FindObject(fgkESDListName[kHeader]);
1221 fESDZDC = (AliESDZDC*)fESDObjects->FindObject(fgkESDListName[kESDZDC]);
1222 fESDFMD = (AliESDFMD*)fESDObjects->FindObject(fgkESDListName[kESDFMD]);
1223 fESDVZERO = (AliESDVZERO*)fESDObjects->FindObject(fgkESDListName[kESDVZERO]);
1224 fESDTZERO = (AliESDTZERO*)fESDObjects->FindObject(fgkESDListName[kESDTZERO]);
1225 fTPCVertex = (AliESDVertex*)fESDObjects->FindObject(fgkESDListName[kTPCVertex]);
1226 fSPDVertex = (AliESDVertex*)fESDObjects->FindObject(fgkESDListName[kSPDVertex]);
1227 fPrimaryVertex = (AliESDVertex*)fESDObjects->FindObject(fgkESDListName[kPrimaryVertex]);
1228 fSPDMult = (AliMultiplicity*)fESDObjects->FindObject(fgkESDListName[kSPDMult]);
1229 fPHOSTrigger = (AliESDCaloTrigger*)fESDObjects->FindObject(fgkESDListName[kPHOSTrigger]);
1230 fEMCALTrigger = (AliESDCaloTrigger*)fESDObjects->FindObject(fgkESDListName[kEMCALTrigger]);
1231 fSPDPileupVertices = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kSPDPileupVertices]);
1232 fTrkPileupVertices = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kTrkPileupVertices]);
1233 fTracks = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kTracks]);
1234 fMuonTracks = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kMuonTracks]);
1235 fPmdTracks = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kPmdTracks]);
1236 fTrdTrigger = (AliESDTrdTrigger*)fESDObjects->FindObject(fgkESDListName[kTrdTrigger]);
1237 fTrdTracks = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kTrdTracks]);
1238 fTrdTracklets = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kTrdTracklets]);
1239 fV0s = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kV0s]);
1240 fCascades = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kCascades]);
1241 fKinks = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kKinks]);
1242 fCaloClusters = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kCaloClusters]);
1243 fEMCALCells = (AliESDCaloCells*)fESDObjects->FindObject(fgkESDListName[kEMCALCells]);
1244 fPHOSCells = (AliESDCaloCells*)fESDObjects->FindObject(fgkESDListName[kPHOSCells]);
1245 fErrorLogs = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kErrorLogs]);
1246 fESDACORDE = (AliESDACORDE*)fESDObjects->FindObject(fgkESDListName[kESDACORDE]);
1247 fTOFHeader = (AliTOFHeader*)fESDObjects->FindObject(fgkESDListName[kTOFHeader]);
1248 fCosmicTracks = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kCosmicTracks]);
1250 // fESDMFT = (AliESDMFT*)fESDObjects->FindObject(fgkESDListName[kESDMFT]);
1254 void AliESDEvent::SetStdNames(){
1255 // Set the names of the standard contents
1257 if(fESDObjects->GetEntries()>=kESDListN){
1258 for(int i = 0;i < fESDObjects->GetEntries() && i<kESDListN;i++){
1259 TObject *fObj = fESDObjects->At(i);
1260 if(fObj->InheritsFrom("TNamed")){
1261 ((TNamed*)fObj)->SetName(fgkESDListName[i]);
1263 else if(fObj->InheritsFrom("TClonesArray")){
1264 ((TClonesArray*)fObj)->SetName(fgkESDListName[i]);
1269 AliWarning("Std Entries missing");
1274 void AliESDEvent::CreateStdContent(Bool_t bUseThisList){
1275 fUseOwnList = bUseThisList;
1279 void AliESDEvent::CreateStdContent()
1281 // create the standard AOD content and set pointers
1283 // create standard objects and add them to the TList of objects
1284 AddObject(new AliESDRun());
1285 AddObject(new AliESDHeader());
1286 AddObject(new AliESDZDC());
1287 AddObject(new AliESDFMD());
1288 AddObject(new AliESDVZERO());
1289 AddObject(new AliESDTZERO());
1290 AddObject(new AliESDVertex());
1291 AddObject(new AliESDVertex());
1292 AddObject(new AliESDVertex());
1293 AddObject(new AliMultiplicity());
1294 AddObject(new AliESDCaloTrigger());
1295 AddObject(new AliESDCaloTrigger());
1296 AddObject(new TClonesArray("AliESDVertex",0));
1297 AddObject(new TClonesArray("AliESDVertex",0));
1298 AddObject(new TClonesArray("AliESDtrack",0));
1299 AddObject(new TClonesArray("AliESDMuonTrack",0));
1300 AddObject(new TClonesArray("AliESDPmdTrack",0));
1301 AddObject(new AliESDTrdTrigger());
1302 AddObject(new TClonesArray("AliESDTrdTrack",0));
1303 AddObject(new TClonesArray("AliESDTrdTracklet",0));
1304 AddObject(new TClonesArray("AliESDv0",0));
1305 AddObject(new TClonesArray("AliESDcascade",0));
1306 AddObject(new TClonesArray("AliESDkink",0));
1307 AddObject(new TClonesArray("AliESDCaloCluster",0));
1308 AddObject(new AliESDCaloCells());
1309 AddObject(new AliESDCaloCells());
1310 AddObject(new TClonesArray("AliRawDataErrorLog",0));
1311 AddObject(new AliESDACORDE());
1312 AddObject(new AliTOFHeader());
1313 AddObject(new TClonesArray("AliESDCosmicTrack",0));
1315 //AddObject(new AliESDMFT());
1318 // check the order of the indices against enum...
1322 // read back pointers
1326 TObject* AliESDEvent::FindListObject(const char *name) const {
1328 // Find object with name "name" in the list of branches
1331 return fESDObjects->FindObject(name);
1336 Int_t AliESDEvent::GetPHOSClusters(TRefArray *clusters) const
1338 // fills the provided TRefArray with all found phos clusters
1342 AliESDCaloCluster *cl = 0;
1343 for (Int_t i = 0; i < GetNumberOfCaloClusters(); i++) {
1345 if ( (cl = GetCaloCluster(i)) ) {
1348 AliDebug(1,Form("IsPHOS cluster %d Size: %d \n",i,clusters->GetEntriesFast()));
1352 return clusters->GetEntriesFast();
1355 Int_t AliESDEvent::GetEMCALClusters(TRefArray *clusters) const
1357 // fills the provided TRefArray with all found emcal clusters
1361 AliESDCaloCluster *cl = 0;
1362 for (Int_t i = 0; i < GetNumberOfCaloClusters(); i++) {
1364 if ( (cl = GetCaloCluster(i)) ) {
1367 AliDebug(1,Form("IsEMCAL cluster %d Size: %d \n",i,clusters->GetEntriesFast()));
1371 return clusters->GetEntriesFast();
1374 void AliESDEvent::WriteToTree(TTree* tree) const {
1375 // Book the branches as in TTree::Branch(TCollection*)
1376 // but add a "." at the end of top level branches which are
1377 // not a TClonesArray
1381 TIter next(fESDObjects);
1382 const Int_t kSplitlevel = 99; // default value in TTree::Branch()
1383 const Int_t kBufsize = 32000; // default value in TTree::Branch()
1386 while ((obj = next())) {
1387 branchname.Form("%s", obj->GetName());
1388 if(branchname.CompareTo("AliESDfriend")==0)branchname = "ESDfriend.";
1389 if ((kSplitlevel > 1) && !obj->InheritsFrom(TClonesArray::Class())) {
1390 if(!branchname.EndsWith("."))branchname += ".";
1392 if (!tree->FindBranch(branchname)) {
1393 // For the custom streamer to be called splitlevel
1394 // has to be negative, only needed for HLT
1395 Int_t splitLevel = (TString(obj->ClassName()) == "AliHLTGlobalTriggerDecision") ? -1 : kSplitlevel - 1;
1396 tree->Bronch(branchname, obj->ClassName(), fESDObjects->GetObjectRef(obj),kBufsize, splitLevel);
1402 void AliESDEvent::ReadFromTree(TTree *tree, Option_t* opt){
1404 // Connect the ESDEvent to a tree
1407 AliWarning("AliESDEvent::ReadFromTree() Zero Pointer to Tree \n");
1411 if(!tree->GetTree())tree->LoadTree(0);
1413 // if we find the "ESD" branch on the tree we do have the old structure
1414 if(tree->GetBranch("ESD")) {
1415 char ** address = (char **)(tree->GetBranch("ESD")->GetAddress());
1416 // do we have the friend branch
1417 TBranch * esdFB = tree->GetBranch("ESDfriend.");
1418 char ** addressF = 0;
1419 if(esdFB)addressF = (char **)(esdFB->GetAddress());
1421 AliInfo("AliESDEvent::ReadFromTree() Reading old Tree");
1422 tree->SetBranchAddress("ESD", &fESDOld);
1424 tree->SetBranchAddress("ESDfriend.",&fESDFriendOld);
1427 AliInfo("AliESDEvent::ReadFromTree() Reading old Tree");
1428 AliInfo("Branch already connected. Using existing branch address.");
1429 fESDOld = (AliESD*) (*address);
1430 // addressF can still be 0, since branch needs to switched on
1431 if(addressF)fESDFriendOld = (AliESDfriend*) (*addressF);
1434 // have already connected the old ESD structure... ?
1435 // reuse also the pointer of the AlliESDEvent
1436 // otherwise create new ones
1437 TList* connectedList = (TList*) (tree->GetUserInfo()->FindObject("ESDObjectsConnectedToTree"));
1440 // If connected use the connected list of objects
1441 if(fESDObjects!= connectedList){
1442 // protect when called twice
1443 fESDObjects->Delete();
1444 fESDObjects = connectedList;
1449 // The pointer to the friend changes when called twice via InitIO
1450 // since AliESDEvent is deleted
1451 TObject* oldf = FindListObject("AliESDfriend");
1454 newf = (TObject*)*addressF;
1456 if(newf!=0&&oldf!=newf){
1457 // remove the old reference
1458 // Should we also delete it? Or is this handled in TTree I/O
1459 // since it is created by the first SetBranchAddress
1460 fESDObjects->Remove(oldf);
1462 fESDObjects->Add(newf);
1469 CreateStdContent(); // create for copy
1470 // if we have the esdfriend add it, so we always can access it via the userinfo
1471 if(fESDFriendOld)AddObject(fESDFriendOld);
1472 // we are not owner of the list objects
1473 // must not delete it
1474 fESDObjects->SetOwner(kTRUE);
1475 fESDObjects->SetName("ESDObjectsConnectedToTree");
1476 tree->GetUserInfo()->Add(fESDObjects);
1484 // Try to find AliESDEvent
1485 AliESDEvent *esdEvent = 0;
1486 esdEvent = (AliESDEvent*)tree->GetTree()->GetUserInfo()->FindObject("AliESDEvent");
1488 // Check if already connected to tree
1490 TList* connectedList = (TList*) (tree->GetUserInfo()->FindObject("ESDObjectsConnectedToTree"));
1493 if (connectedList && (strcmp(opt, "reconnect"))) {
1494 // If connected use the connected list if objects
1495 fESDObjects->Delete();
1496 fESDObjects = connectedList;
1503 // prevent a memory leak when reading back the TList
1504 // if (!(strcmp(opt, "reconnect"))) fESDObjects->Delete();
1507 // create a new TList from the UserInfo TList...
1508 // copy constructor does not work...
1509 fESDObjects = (TList*)(esdEvent->GetList()->Clone());
1510 fESDObjects->SetOwner(kTRUE);
1512 else if ( fESDObjects->GetEntries()==0){
1513 // at least create the std content if we want to read to our list
1518 // we only need new things in the list if we do no already have it..
1519 // TODO just add new entries
1521 if(fESDObjects->GetEntries()<kESDListN){
1522 AliWarning(Form("AliESDEvent::ReadFromTree() TList contains less than the standard contents %d < %d \n",
1523 fESDObjects->GetEntries(),kESDListN));
1525 // set the branch addresses
1526 TIter next(fESDObjects);
1528 while((el=(TNamed*)next())){
1529 TString bname(el->GetName());
1530 if(bname.CompareTo("AliESDfriend")==0)
1532 // AliESDfriend does not have a name ...
1533 TBranch *br = tree->GetBranch("ESDfriend.");
1534 if (br) tree->SetBranchAddress("ESDfriend.",fESDObjects->GetObjectRef(el));
1537 // check if branch exists under this Name
1538 TBranch *br = tree->GetBranch(bname.Data());
1540 tree->SetBranchAddress(bname.Data(),fESDObjects->GetObjectRef(el));
1543 br = tree->GetBranch(Form("%s.",bname.Data()));
1545 tree->SetBranchAddress(Form("%s.",bname.Data()),fESDObjects->GetObjectRef(el));
1548 AliWarning(Form("AliESDEvent::ReadFromTree() No Branch found with Name %s or %s.",bname.Data(),bname.Data()));
1555 // when reading back we are not owner of the list
1556 // must not delete it
1557 fESDObjects->SetOwner(kTRUE);
1558 fESDObjects->SetName("ESDObjectsConnectedToTree");
1559 // we are not owner of the list objects
1560 // must not delete it
1561 tree->GetUserInfo()->Add(fESDObjects);
1562 tree->GetUserInfo()->SetOwner(kFALSE);
1566 // we can't get the list from the user data, create standard content
1567 // and set it by hand (no ESDfriend at the moment
1569 TIter next(fESDObjects);
1571 while((el=(TNamed*)next())){
1572 TString bname(el->GetName());
1573 TBranch *br = tree->GetBranch(bname.Data());
1575 tree->SetBranchAddress(bname.Data(),fESDObjects->GetObjectRef(el));
1578 br = tree->GetBranch(Form("%s.",bname.Data()));
1580 tree->SetBranchAddress(Form("%s.",bname.Data()),fESDObjects->GetObjectRef(el));
1585 // when reading back we are not owner of the list
1586 // must not delete it
1587 fESDObjects->SetOwner(kTRUE);
1592 void AliESDEvent::CopyFromOldESD()
1594 // Method which copies over everthing from the old esd structure to the
1599 SetRunNumber(fESDOld->GetRunNumber());
1600 SetPeriodNumber(fESDOld->GetPeriodNumber());
1601 SetMagneticField(fESDOld->GetMagneticField());
1603 // leave out diamond ...
1604 // SetDiamond(const AliESDVertex *vertex) { fESDRun->SetDiamond(vertex);}
1607 SetTriggerMask(fESDOld->GetTriggerMask());
1608 SetOrbitNumber(fESDOld->GetOrbitNumber());
1609 SetTimeStamp(fESDOld->GetTimeStamp());
1610 SetEventType(fESDOld->GetEventType());
1611 SetEventNumberInFile(fESDOld->GetEventNumberInFile());
1612 SetBunchCrossNumber(fESDOld->GetBunchCrossNumber());
1613 SetTriggerCluster(fESDOld->GetTriggerCluster());
1617 SetZDC(fESDOld->GetZDCN1Energy(),
1618 fESDOld->GetZDCP1Energy(),
1619 fESDOld->GetZDCEMEnergy(),
1621 fESDOld->GetZDCN2Energy(),
1622 fESDOld->GetZDCP2Energy(),
1623 fESDOld->GetZDCParticipants(),
1633 if(fESDOld->GetFMDData())SetFMDData(fESDOld->GetFMDData());
1637 SetT0zVertex(fESDOld->GetT0zVertex());
1638 SetT0(fESDOld->GetT0());
1642 if (fESDOld->GetVZEROData()) SetVZEROData(fESDOld->GetVZEROData());
1644 if(fESDOld->GetVertex())SetPrimaryVertexSPD(fESDOld->GetVertex());
1646 if(fESDOld->GetPrimaryVertex())SetPrimaryVertexTracks(fESDOld->GetPrimaryVertex());
1648 if(fESDOld->GetMultiplicity())SetMultiplicity(fESDOld->GetMultiplicity());
1650 for(int i = 0;i<fESDOld->GetNumberOfTracks();i++){
1651 AddTrack(fESDOld->GetTrack(i));
1654 for(int i = 0;i<fESDOld->GetNumberOfMuonTracks();i++){
1655 AddMuonTrack(fESDOld->GetMuonTrack(i));
1658 for(int i = 0;i<fESDOld->GetNumberOfPmdTracks();i++){
1659 AddPmdTrack(fESDOld->GetPmdTrack(i));
1662 for(int i = 0;i<fESDOld->GetNumberOfTrdTracks();i++){
1663 AddTrdTrack(fESDOld->GetTrdTrack(i));
1666 for(int i = 0;i<fESDOld->GetNumberOfV0s();i++){
1667 AddV0(fESDOld->GetV0(i));
1670 for(int i = 0;i<fESDOld->GetNumberOfCascades();i++){
1671 AddCascade(fESDOld->GetCascade(i));
1674 for(int i = 0;i<fESDOld->GetNumberOfKinks();i++){
1675 AddKink(fESDOld->GetKink(i));
1679 for(int i = 0;i<fESDOld->GetNumberOfCaloClusters();i++){
1680 AddCaloCluster(fESDOld->GetCaloCluster(i));
1685 // if (fESDOld->GetMFTData()) SetMFTData(fESDOld->GetMFTData());
1691 Bool_t AliESDEvent::IsEventSelected(const char *trigExpr) const
1693 // Check if the event satisfies the trigger
1694 // selection expression trigExpr.
1695 // trigExpr can be any logical expression
1696 // of the trigger classes defined in AliESDRun
1697 // In case of wrong syntax return kTRUE.
1699 TString expr(trigExpr);
1700 if (expr.IsNull()) return kTRUE;
1702 ULong64_t mask = GetTriggerMask();
1703 for(Int_t itrig = 0; itrig < AliESDRun::kNTriggerClasses; itrig++) {
1704 if (mask & (1ull << itrig)) {
1705 expr.ReplaceAll(GetESDRun()->GetTriggerClass(itrig),"1");
1708 expr.ReplaceAll(GetESDRun()->GetTriggerClass(itrig),"0");
1713 if ((gROOT->ProcessLineFast(expr.Data(),&error) == 0) &&
1714 (error == TInterpreter::kNoError)) {
1722 TObject* AliESDEvent::GetHLTTriggerDecision() const
1724 // get the HLT trigger decission object
1726 // cast away const'nes because the FindListObject method
1728 AliESDEvent* pNonConst=const_cast<AliESDEvent*>(this);
1729 return pNonConst->FindListObject("HLTGlobalTrigger");
1732 TString AliESDEvent::GetHLTTriggerDescription() const
1734 // get the HLT trigger decission description
1735 TString description;
1736 TObject* pDecision=GetHLTTriggerDecision();
1738 description=pDecision->GetTitle();
1744 Bool_t AliESDEvent::IsHLTTriggerFired(const char* name) const
1746 // get the HLT trigger decission description
1747 TObject* pDecision=GetHLTTriggerDecision();
1748 if (!pDecision) return kFALSE;
1750 Option_t* option=pDecision->GetOption();
1751 if (option==NULL || *option!='1') return kFALSE;
1754 TString description=GetHLTTriggerDescription();
1755 Int_t index=description.Index(name);
1756 if (index<0) return kFALSE;
1757 index+=strlen(name);
1758 if (index>=description.Length()) return kFALSE;
1759 if (description[index]!=0 && description[index]!=' ') return kFALSE;
1764 //______________________________________________________________________________
1765 Bool_t AliESDEvent::IsPileupFromSPD(Int_t minContributors,
1767 Double_t nSigmaZdist,
1768 Double_t nSigmaDiamXY,
1769 Double_t nSigmaDiamZ) const{
1771 // This function checks if there was a pile up
1772 // reconstructed with SPD
1774 Int_t nc1=fSPDVertex->GetNContributors();
1775 if(nc1<1) return kFALSE;
1776 Int_t nPileVert=GetNumberOfPileupVerticesSPD();
1777 if(nPileVert==0) return kFALSE;
1779 for(Int_t i=0; i<nPileVert;i++){
1780 const AliESDVertex* pv=GetPileupVertexSPD(i);
1781 Int_t nc2=pv->GetNContributors();
1782 if(nc2>=minContributors){
1783 Double_t z1=fSPDVertex->GetZ();
1784 Double_t z2=pv->GetZ();
1785 Double_t distZ=TMath::Abs(z2-z1);
1786 Double_t distZdiam=TMath::Abs(z2-GetDiamondZ());
1787 Double_t cutZdiam=nSigmaDiamZ*TMath::Sqrt(GetSigma2DiamondZ());
1788 if(GetSigma2DiamondZ()<0.0001)cutZdiam=99999.; //protection for missing z diamond information
1789 if(distZ>minZdist && distZdiam<cutZdiam){
1790 Double_t x2=pv->GetX();
1791 Double_t y2=pv->GetY();
1792 Double_t distXdiam=TMath::Abs(x2-GetDiamondX());
1793 Double_t distYdiam=TMath::Abs(y2-GetDiamondY());
1794 Double_t cov1[6],cov2[6];
1795 fSPDVertex->GetCovarianceMatrix(cov1);
1796 pv->GetCovarianceMatrix(cov2);
1797 Double_t errxDist=TMath::Sqrt(cov2[0]+GetSigma2DiamondX());
1798 Double_t erryDist=TMath::Sqrt(cov2[2]+GetSigma2DiamondY());
1799 Double_t errzDist=TMath::Sqrt(cov1[5]+cov2[5]);
1800 Double_t cutXdiam=nSigmaDiamXY*errxDist;
1801 if(GetSigma2DiamondX()<0.0001)cutXdiam=99999.; //protection for missing diamond information
1802 Double_t cutYdiam=nSigmaDiamXY*erryDist;
1803 if(GetSigma2DiamondY()<0.0001)cutYdiam=99999.; //protection for missing diamond information
1804 if( (distXdiam<cutXdiam) && (distYdiam<cutYdiam) && (distZ>nSigmaZdist*errzDist) ){
1813 //______________________________________________________________________________
1814 void AliESDEvent::EstimateMultiplicity(Int_t &tracklets, Int_t &trITSTPC, Int_t &trITSSApure, Double_t eta, Bool_t useDCAFlag,Bool_t useV0Flag) const
1817 // calculates 3 estimators for the multiplicity in the -eta:eta range
1818 // tracklets : using SPD tracklets only
1819 // trITSTPC : using TPC/ITS + complementary ITS SA tracks + tracklets from clusters not used by tracks
1820 // trITSSApure : using ITS standalone tracks + tracklets from clusters not used by tracks
1821 // if useDCAFlag is true: account for the ESDtrack flag marking the tracks with large DCA
1822 // if useV0Flag is true: account for the ESDtrack flag marking conversion and K0's V0s
1823 tracklets = trITSSApure = trITSTPC = 0;
1824 int ntr = fSPDMult ? fSPDMult->GetNumberOfTracklets() : 0;
1827 for (int itr=ntr;itr--;) {
1828 if (TMath::Abs(fSPDMult->GetEta(itr))>eta) continue;
1830 if (fSPDMult->FreeClustersTracklet(itr,0)) trITSTPC++; // not used in ITS/TPC or ITS_SA track
1831 if (fSPDMult->FreeClustersTracklet(itr,1)) trITSSApure++; // not used in ITS_SA_Pure track
1834 // count real tracks
1835 ntr = GetNumberOfTracks();
1836 for (int itr=ntr;itr--;) {
1837 AliESDtrack *t = GetTrack(itr);
1838 if (!t) {AliFatal(Form("NULL pointer for ESD track %d",itr));}
1839 if (TMath::Abs(t->Eta())>eta) continue;
1840 if (!t->IsOn(AliESDtrack::kITSin)) continue;
1841 if (useDCAFlag && t->IsOn(AliESDtrack::kMultSec)) continue;
1842 if (useV0Flag && t->IsOn(AliESDtrack::kMultInV0)) continue;
1843 if (t->IsOn(AliESDtrack::kITSpureSA)) trITSSApure++;
1849 Bool_t AliESDEvent::IsPileupFromSPDInMultBins() const {
1850 Int_t nTracklets=GetMultiplicity()->GetNumberOfTracklets();
1851 if(nTracklets<20) return IsPileupFromSPD(3,0.8);
1852 else if(nTracklets<50) return IsPileupFromSPD(4,0.8);
1853 else return IsPileupFromSPD(5,0.8);
1856 void AliESDEvent::SetTOFHeader(const AliTOFHeader *header)
1859 // Set the TOF event_time
1863 *fTOFHeader=*header;
1864 //fTOFHeader->SetName(fgkESDListName[kTOFHeader]);
1867 // for analysis of reconstructed events
1868 // when this information is not avaliable
1869 fTOFHeader = new AliTOFHeader(*header);
1870 //AddObject(fTOFHeader);
1875 AliCentrality* AliESDEvent::GetCentrality()
1877 if (!fCentrality) fCentrality = new AliCentrality();
1881 AliEventplane* AliESDEvent::GetEventplane()
1883 if (!fEventplane) fEventplane = new AliEventplane();
1887 Float_t AliESDEvent::GetVZEROEqMultiplicity(Int_t i) const
1889 // Get VZERO Multiplicity for channel i
1890 // Themethod uses the equalization factors
1891 // stored in the ESD-run object in order to
1892 // get equal multiplicities within a VZERO rins (1/8 of VZERO)
1893 if (!fESDVZERO || !fESDRun) return -1;
1896 Float_t factorSum = 0;
1897 for(Int_t j = 8*ring; j < (8*ring+8); ++j) {
1898 factorSum += fESDRun->GetVZEROEqFactors(j);
1900 Float_t factor = fESDRun->GetVZEROEqFactors(i)*8./factorSum;
1902 return (fESDVZERO->GetMultiplicity(i)/factor);