1 /**************************************************************************
2 * Copyright(c) 1998-2007, 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 //-------------------------------------------------------------------------
20 // Author: Markus Oldenburg, CERN
21 //-------------------------------------------------------------------------
26 #include <TFriendElement.h>
27 #include <TProcessID.h>
28 #include <TCollection.h>
29 #include "Riostream.h"
30 #include "AliAODEvent.h"
31 #include "AliAODHeader.h"
32 #include "AliAODTrack.h"
33 #include "AliAODDimuon.h"
37 // definition of std AOD member names
38 const char* AliAODEvent::fAODListName[kAODListN] = {"header",
53 //______________________________________________________________________________
54 AliAODEvent::AliAODEvent() :
56 fAODObjects(new TList()),
73 // default constructor
76 //______________________________________________________________________________
77 AliAODEvent::AliAODEvent(const AliAODEvent& aod):
79 fAODObjects(new TList()),
82 fHeader(new AliAODHeader(*aod.fHeader)),
83 fTracks(new TClonesArray(*aod.fTracks)),
84 fVertices(new TClonesArray(*aod.fVertices)),
85 fV0s(new TClonesArray(*aod.fV0s)),
86 fCascades(new TClonesArray(*aod.fCascades)),
87 fTracklets(new AliAODTracklets(*aod.fTracklets)),
88 fJets(new TClonesArray(*aod.fJets)),
89 fEmcalCells(new AliAODCaloCells(*aod.fEmcalCells)),
90 fPhosCells(new AliAODCaloCells(*aod.fPhosCells)),
91 fCaloClusters(new TClonesArray(*aod.fCaloClusters)),
92 fFmdClusters(new TClonesArray(*aod.fFmdClusters)),
93 fPmdClusters(new TClonesArray(*aod.fPmdClusters)),
94 fDimuons(new TClonesArray(*aod.fDimuons))
101 AddObject(fCascades);
102 AddObject(fTracklets);
104 AddObject(fEmcalCells);
105 AddObject(fPhosCells);
106 AddObject(fCaloClusters);
107 AddObject(fFmdClusters);
108 AddObject(fPmdClusters);
110 fConnected = aod.fConnected;
115 //______________________________________________________________________________
116 AliAODEvent & AliAODEvent::operator=(const AliAODEvent& aod) {
118 // Assignment operator
120 if(&aod == this) return *this;
121 AliVEvent::operator=(aod);
123 // This assumes that the list is already created
124 // and that the virtual void Copy(Tobject&) function
125 // is correctly implemented in the derived class
126 // otherwise only TObject::Copy() will be used
128 if((fAODObjects->GetSize()==0)&&(aod.fAODObjects->GetSize()>=kAODListN)){
129 // We cover the case that we do not yet have the
130 // standard content but the source has it
134 // Here we have the standard content without user additions, but the content is
135 // not matching the aod source.
137 // Iterate the list of source objects
138 TIter next(aod.GetList());
141 while ((its = next())) {
142 name = its->GetName();
143 // Check if we have this object type in out list
144 TObject *mine = fAODObjects->FindObject(name);
146 // We have to create the same type of object.
147 TClass* pClass=TClass::GetClass(its->ClassName());
149 AliWarning(Form("Can not find class description for entry %s (%s)\n",
150 its->ClassName(), name.Data()));
153 mine=(TObject*)pClass->New();
155 // not in this: can be added to list
156 AliWarning(Form("%s:%d Could not find %s for copying \n",
157 (char*)__FILE__,__LINE__,name.Data()));
160 if(mine->InheritsFrom("TNamed")) {
161 ((TNamed*)mine)->SetName(name);
162 } else if(mine->InheritsFrom("TCollection")){
163 if(mine->InheritsFrom("TClonesArray")) {
164 TClonesArray *itscl = dynamic_cast<TClonesArray*>(its);
166 AliWarning(Form("Class description for entry %s (%s) not TClonesArray\n",
167 its->ClassName(), name.Data()));
171 dynamic_cast<TClonesArray*>(mine)->SetClass(itscl->GetClass(), itscl->GetSize());
173 dynamic_cast<TCollection*>(mine)->SetName(name);
175 AliDebug(1, Form("adding object %s of type %s", mine->GetName(), mine->ClassName()));
178 // Now we have an object of the same type and name, but different content.
179 if(!its->InheritsFrom("TCollection")){
180 // simple objects (do they have a Copy method that calls operator= ?)
182 } else if (its->InheritsFrom("TClonesArray")) {
183 // Create or expand the tclonesarray pointers
184 // so we can directly copy to the object
185 TClonesArray *its_tca = (TClonesArray*)its;
186 TClonesArray *mine_tca = (TClonesArray*)mine;
187 // this leaves the capacity of the TClonesArray the same
188 // except for a factor of 2 increase when size > capacity
189 // does not release any memory occupied by the tca
190 Int_t its_entries = its_tca->GetEntriesFast();
191 mine_tca->ExpandCreate(its_entries);
192 for(int i=0; i<its_entries; i++){
194 TObject *mine_tca_obj = mine_tca->At(i);
195 TObject *its_tca_obj = its_tca->At(i);
196 // no need to delete first
197 // pointers within the class should be handled by Copy()...
198 // Can there be Empty slots?
199 its_tca_obj->Copy(*mine_tca_obj);
202 AliWarning(Form("%s:%d cannot copy TCollection \n",
203 (char*)__FILE__,__LINE__));
206 fConnected = aod.fConnected;
211 //______________________________________________________________________________
212 AliAODEvent::~AliAODEvent()
215 if(fAODObjects&&!fConnected)
225 //______________________________________________________________________________
226 void AliAODEvent::AddObject(TObject* obj)
228 // Add an object to the list of objects.
229 // Please be aware that in order to increase performance you should
230 // refrain from using TObjArrays (if possible). Use TClonesArrays, instead.
232 fAODObjects->AddLast(obj);
235 //______________________________________________________________________________
236 void AliAODEvent::RemoveObject(TObject* obj)
238 // Removes an object from the list of objects.
240 fAODObjects->Remove(obj);
243 //______________________________________________________________________________
244 TObject *AliAODEvent::FindListObject(const char *objName) const
246 // Return the pointer to the object with the given name.
248 return fAODObjects->FindObject(objName);
251 //______________________________________________________________________________
252 void AliAODEvent::CreateStdContent()
254 // create the standard AOD content and set pointers
256 // create standard objects and add them to the TList of objects
257 AddObject(new AliAODHeader());
258 AddObject(new TClonesArray("AliAODTrack", 0));
259 AddObject(new TClonesArray("AliAODVertex", 0));
260 AddObject(new TClonesArray("AliAODv0", 0));
261 AddObject(new TClonesArray("AliAODcascade", 0));
262 AddObject(new AliAODTracklets());
263 AddObject(new TClonesArray("AliAODJet", 0));
264 AddObject(new AliAODCaloCells());
265 AddObject(new AliAODCaloCells());
266 AddObject(new TClonesArray("AliAODCaloCluster", 0));
267 AddObject(new TClonesArray("AliAODFmdCluster", 0));
268 AddObject(new TClonesArray("AliAODPmdCluster", 0));
269 AddObject(new TClonesArray("AliAODDimuon", 0));
273 // read back pointers
279 void AliAODEvent::MakeEntriesReferencable()
281 // Make all entries referencable in a subsequent process
283 TIter next(fAODObjects);
285 while ((obj = next()))
287 if(obj->InheritsFrom("TCollection"))
289 AssignIDtoCollection((TCollection*)obj);
294 //______________________________________________________________________________
295 void AliAODEvent::SetStdNames()
297 // introduce the standard naming
299 if(fAODObjects->GetEntries()==kAODListN){
300 for(int i = 0;i < fAODObjects->GetEntries();i++){
301 TObject *fObj = fAODObjects->At(i);
302 if(fObj->InheritsFrom("TNamed")){
303 ((TNamed*)fObj)->SetName(fAODListName[i]);
305 else if(fObj->InheritsFrom("TClonesArray")){
306 ((TClonesArray*)fObj)->SetName(fAODListName[i]);
311 printf("%s:%d SetStdNames() Wrong number of Std Entries \n",(char*)__FILE__,__LINE__);
315 void AliAODEvent::CreateStdFolders()
317 // Create the standard folder structure
318 if(fAODFolder)delete fAODFolder;
319 fAODFolder = gROOT->GetRootFolder()->AddFolder("AOD", "AOD");
320 if(fAODObjects->GetEntries()==kAODListN){
321 for(int i = 0;i < fAODObjects->GetEntries();i++){
322 TObject *fObj = fAODObjects->At(i);
323 if(fObj->InheritsFrom("TClonesArray")){
324 fAODFolder->AddFolder(fAODListName[i], fAODListName[i], (TCollection*) fObj);
326 fAODFolder->AddFolder(fAODListName[i], fAODListName[i], 0);
331 printf("%s:%d CreateStdFolders() Wrong number of Std Entries \n",(char*)__FILE__,__LINE__);
335 //______________________________________________________________________________
336 void AliAODEvent::GetStdContent()
338 // set pointers for standard content
340 fHeader = (AliAODHeader*)fAODObjects->FindObject("header");
341 fTracks = (TClonesArray*)fAODObjects->FindObject("tracks");
342 fVertices = (TClonesArray*)fAODObjects->FindObject("vertices");
343 fV0s = (TClonesArray*)fAODObjects->FindObject("v0s");
344 fCascades = (TClonesArray*)fAODObjects->FindObject("cascades");
345 fTracklets = (AliAODTracklets*)fAODObjects->FindObject("tracklets");
346 fJets = (TClonesArray*)fAODObjects->FindObject("jets");
347 fEmcalCells = (AliAODCaloCells*)fAODObjects->FindObject("emcalCells");
348 fPhosCells = (AliAODCaloCells*)fAODObjects->FindObject("phosCells");
349 fCaloClusters = (TClonesArray*)fAODObjects->FindObject("caloClusters");
350 fFmdClusters = (TClonesArray*)fAODObjects->FindObject("fmdClusters");
351 fPmdClusters = (TClonesArray*)fAODObjects->FindObject("pmdClusters");
352 fDimuons = (TClonesArray*)fAODObjects->FindObject("dimuons");
355 //______________________________________________________________________________
356 void AliAODEvent::ResetStd(Int_t trkArrSize,
359 Int_t cascadeArrSize,
367 // deletes content of standard arrays and resets size
370 if (trkArrSize > fTracks->GetSize())
371 fTracks->Expand(trkArrSize);
374 if (vtxArrSize > fVertices->GetSize())
375 fVertices->Expand(vtxArrSize);
378 if (v0ArrSize > fV0s->GetSize())
379 fV0s->Expand(v0ArrSize);
382 if (cascadeArrSize > fCascades->GetSize())
383 fCascades->Expand(cascadeArrSize);
386 if (jetSize > fJets->GetSize())
387 fJets->Expand(jetSize);
389 fCaloClusters->Delete();
390 if (caloClusSize > fCaloClusters->GetSize())
391 fCaloClusters->Expand(caloClusSize);
393 fFmdClusters->Delete();
394 if (fmdClusSize > fFmdClusters->GetSize())
395 fFmdClusters->Expand(fmdClusSize);
397 fPmdClusters->Delete();
398 if (pmdClusSize > fPmdClusters->GetSize())
399 fPmdClusters->Expand(pmdClusSize);
402 if (dimuonArrSize > fDimuons->GetSize())
403 fDimuons->Expand(dimuonArrSize);
405 // Reset the tracklets
406 fTracklets->DeleteContainer();
407 fPhosCells->DeleteContainer();
408 fEmcalCells->DeleteContainer();
412 void AliAODEvent::ClearStd()
414 // clears the standard arrays
415 fHeader ->RemoveQTheta();
417 fVertices ->Delete();
419 fCascades ->Delete();
420 fTracklets ->DeleteContainer();
422 fEmcalCells ->DeleteContainer();
423 fPhosCells ->DeleteContainer();
424 fCaloClusters ->Delete();
425 fFmdClusters ->Clear();
426 fPmdClusters ->Clear();
430 //_________________________________________________________________
431 Int_t AliAODEvent::GetPHOSClusters(TRefArray *clusters) const
433 // fills the provided TRefArray with all found phos clusters
437 AliAODCaloCluster *cl = 0;
438 Bool_t first = kTRUE;
439 for (Int_t i = 0; i < GetNumberOfCaloClusters() ; i++) {
440 if ( (cl = GetCaloCluster(i)) ) {
443 new (clusters) TRefArray(TProcessID::GetProcessWithUID(cl));
447 //printf("IsPHOS cluster %d, E %2.3f Size: %d \n",i,cl->E(),clusters->GetEntriesFast());
451 return clusters->GetEntriesFast();
454 //_________________________________________________________________
455 Int_t AliAODEvent::GetEMCALClusters(TRefArray *clusters) const
457 // fills the provided TRefArray with all found emcal clusters
460 AliAODCaloCluster *cl = 0;
461 Bool_t first = kTRUE;
462 for (Int_t i = 0; i < GetNumberOfCaloClusters(); i++) {
463 if ( (cl = GetCaloCluster(i)) ) {
466 new (clusters) TRefArray(TProcessID::GetProcessWithUID(cl));
470 //printf("IsEMCal cluster %d, E %2.3f Size: %d \n",i,cl->E(),clusters->GetEntriesFast());
474 return clusters->GetEntriesFast();
478 //______________________________________________________________________________
479 Int_t AliAODEvent::GetMuonTracks(TRefArray *muonTracks) const
481 // fills the provided TRefArray with all found muon tracks
485 AliAODTrack *track = 0;
486 for (Int_t iTrack = 0; iTrack < GetNTracks(); iTrack++) {
487 track = GetTrack(iTrack);
488 if (track->IsMuonTrack()) {
489 muonTracks->Add(track);
493 return muonTracks->GetEntriesFast();
497 //______________________________________________________________________________
498 Int_t AliAODEvent::GetNumberOfMuonTracks() const
500 // get number of muon tracks
502 for (Int_t iTrack = 0; iTrack < GetNTracks(); iTrack++) {
503 if ((GetTrack(iTrack))->IsMuonTrack()) {
511 void AliAODEvent::ReadFromTree(TTree *tree, Option_t* opt /*= ""*/)
513 // Connects aod event to tree
516 Printf("%s %d AliAODEvent::ReadFromTree() Zero Pointer to Tree \n",(char*)__FILE__,__LINE__);
520 if(!tree->GetTree())tree->LoadTree(0);
522 // Try to find AliAODEvent
523 AliAODEvent *aodEvent = 0;
524 aodEvent = (AliAODEvent*)tree->GetTree()->GetUserInfo()->FindObject("AliAODEvent");
526 // Check if already connected to tree
527 TList* connectedList = (TList*) (tree->GetUserInfo()->FindObject("AODObjectsConnectedToTree"));
528 if (connectedList && (strcmp(opt, "reconnect"))) {
529 // If connected use the connected list of objects
530 fAODObjects->Delete();
531 fAODObjects = connectedList;
537 // prevent a memory leak when reading back the TList
538 // if (!(strcmp(opt, "reconnect"))) fAODObjects->Delete();
540 // create a new TList from the UserInfo TList...
541 // copy constructor does not work...
542 fAODObjects = (TList*)(aodEvent->GetList()->Clone());
543 fAODObjects->SetOwner(kTRUE);
544 if(fAODObjects->GetEntries()<kAODListN){
545 printf("%s %d AliAODEvent::ReadFromTree() TList contains less than the standard contents %d < %d \n",
546 (char*)__FILE__,__LINE__,fAODObjects->GetEntries(),kAODListN);
549 // Let's find out whether we have friends
550 TList* friendL = tree->GetTree()->GetListOfFriends();
555 while ((fe = (TFriendElement*)next())){
556 aodEvent = (AliAODEvent*)(fe->GetTree()->GetUserInfo()->FindObject("AliAODEvent"));
558 printf("No UserInfo on tree \n");
561 TList* objL = (TList*)(aodEvent->GetList()->Clone());
562 printf("Get list of object from tree %d !!\n", objL->GetEntries());
563 TIter nextobject(objL);
565 while((obj = nextobject()))
567 printf("Adding object from friend %s !\n", obj->GetName());
568 fAODObjects->Add(obj);
569 } // object "branch" loop
575 // set the branch addresses
576 TIter next(fAODObjects);
578 while((el=(TNamed*)next())){
579 TString bname(el->GetName());
580 // check if branch exists under this Name
581 TBranch *br = tree->GetTree()->GetBranch(bname.Data());
583 tree->SetBranchAddress(bname.Data(),fAODObjects->GetObjectRef(el));
585 br = tree->GetBranch(Form("%s.",bname.Data()));
587 tree->SetBranchAddress(Form("%s.",bname.Data()),fAODObjects->GetObjectRef(el));
590 printf("%s %d AliAODEvent::ReadFromTree() No Branch found with Name %s. \n",
591 (char*)__FILE__,__LINE__,bname.Data());
596 // when reading back we are not owner of the list
597 // must not delete it
598 fAODObjects->SetOwner(kTRUE);
599 fAODObjects->SetName("AODObjectsConnectedToTree");
600 // we are not owner of the list objects
601 // must not delete it
602 tree->GetUserInfo()->Add(fAODObjects);
606 // we can't get the list from the user data, create standard content
607 // and set it by hand
609 TIter next(fAODObjects);
611 while((el=(TNamed*)next())){
612 TString bname(el->GetName());
613 tree->SetBranchAddress(bname.Data(),fAODObjects->GetObjectRef(el));
616 // when reading back we are not owner of the list
617 // must not delete it
618 fAODObjects->SetOwner(kTRUE);
621 //______________________________________________________________________________
622 Int_t AliAODEvent::GetNumberOfPileupVerticesSPD() const{
623 // count number of SPD pileup vertices
624 Int_t nVertices=GetNumberOfVertices();
625 Int_t nPileupVertices=0;
626 for(Int_t iVert=0; iVert<nVertices; iVert++){
627 AliAODVertex *v=GetVertex(iVert);
628 if(v->GetType()==AliAODVertex::kPileupSPD) nPileupVertices++;
630 return nPileupVertices;
632 //______________________________________________________________________________
633 Int_t AliAODEvent::GetNumberOfPileupVerticesTracks() const{
634 // count number of track pileup vertices
635 Int_t nVertices=GetNumberOfVertices();
636 Int_t nPileupVertices=0;
637 for(Int_t iVert=0; iVert<nVertices; iVert++){
638 AliAODVertex *v=GetVertex(iVert);
639 if(v->GetType()==AliAODVertex::kPileupTracks) nPileupVertices++;
641 return nPileupVertices;
643 //______________________________________________________________________________
644 AliAODVertex* AliAODEvent::GetPrimaryVertexSPD() const{
646 Int_t nVertices=GetNumberOfVertices();
647 for(Int_t iVert=0; iVert<nVertices; iVert++){
648 AliAODVertex *v=GetVertex(iVert);
649 if(v->GetType()==AliAODVertex::kMainSPD) return v;
653 //______________________________________________________________________________
654 AliAODVertex* AliAODEvent::GetPileupVertexSPD(Int_t iV) const{
656 Int_t nVertices=GetNumberOfVertices();
658 for(Int_t iVert=0; iVert<nVertices; iVert++){
659 AliAODVertex *v=GetVertex(iVert);
660 if(v->GetType()==AliAODVertex::kPileupSPD){
661 if(counter==iV) return v;
667 //______________________________________________________________________________
668 AliAODVertex* AliAODEvent::GetPileupVertexTracks(Int_t iV) const{
670 Int_t nVertices=GetNumberOfVertices();
672 for(Int_t iVert=0; iVert<nVertices; iVert++){
673 AliAODVertex *v=GetVertex(iVert);
674 if(v->GetType()==AliAODVertex::kPileupTracks){
675 if(counter==iV) return v;
681 //______________________________________________________________________________
682 Bool_t AliAODEvent::IsPileupFromSPD(Int_t minContributors,
684 Double_t nSigmaZdist,
685 Double_t nSigmaDiamXY,
686 Double_t nSigmaDiamZ) const{
688 // This function checks if there was a pile up
689 // reconstructed with SPD
691 AliAODVertex *mainV=GetPrimaryVertexSPD();
692 if(!mainV) return kFALSE;
693 Int_t nc1=mainV->GetNContributors();
694 if(nc1<1) return kFALSE;
695 Int_t nPileVert=GetNumberOfPileupVerticesSPD();
696 if(nPileVert==0) return kFALSE;
697 Int_t nVertices=GetNumberOfVertices();
699 for(Int_t iVert=0; iVert<nVertices; iVert++){
700 AliAODVertex *pv=GetVertex(iVert);
701 if(pv->GetType()!=AliAODVertex::kPileupSPD) continue;
702 Int_t nc2=pv->GetNContributors();
703 if(nc2>=minContributors){
704 Double_t z1=mainV->GetZ();
705 Double_t z2=pv->GetZ();
706 Double_t distZ=TMath::Abs(z2-z1);
707 Double_t distZdiam=TMath::Abs(z2-GetDiamondZ());
708 Double_t cutZdiam=nSigmaDiamZ*TMath::Sqrt(GetSigma2DiamondZ());
709 if(GetSigma2DiamondZ()<0.0001)cutZdiam=99999.; //protection for missing z diamond information
710 if(distZ>minZdist && distZdiam<cutZdiam){
711 Double_t x2=pv->GetX();
712 Double_t y2=pv->GetY();
713 Double_t distXdiam=TMath::Abs(x2-GetDiamondX());
714 Double_t distYdiam=TMath::Abs(y2-GetDiamondY());
715 Double_t cov1[6],cov2[6];
716 mainV->GetCovarianceMatrix(cov1);
717 pv->GetCovarianceMatrix(cov2);
718 Double_t errxDist=TMath::Sqrt(cov2[0]+GetSigma2DiamondX());
719 Double_t erryDist=TMath::Sqrt(cov2[2]+GetSigma2DiamondY());
720 Double_t errzDist=TMath::Sqrt(cov1[5]+cov2[5]);
721 Double_t cutXdiam=nSigmaDiamXY*errxDist;
722 if(GetSigma2DiamondX()<0.0001)cutXdiam=99999.; //protection for missing diamond information
723 Double_t cutYdiam=nSigmaDiamXY*erryDist;
724 if(GetSigma2DiamondY()<0.0001)cutYdiam=99999.; //protection for missing diamond information
725 if( (distXdiam<cutXdiam) && (distYdiam<cutYdiam) && (distZ>nSigmaZdist*errzDist) ){
734 //______________________________________________________________________________
735 void AliAODEvent::Print(Option_t *) const
737 // Print the names of the all branches
738 TIter next(fAODObjects);
740 Printf(">>>>> AOD Content <<<<<");
741 while((el=(TNamed*)next())){
742 Printf(">> %s ",el->GetName());
744 Printf(">>>>> <<<<<");
749 void AliAODEvent::AssignIDtoCollection(TCollection* col)
751 // Static method which assigns a ID to each object in a collection
752 // In this way the objects are marked as referenced and written with
753 // an ID. This has the advantage that TRefs to this objects can be
754 // written by a subsequent process.
757 while ((obj = next()))
758 TProcessID::AssignID(obj);
761 Bool_t AliAODEvent::IsPileupFromSPDInMultBins() const {
762 Int_t nTracklets=GetTracklets()->GetNumberOfTracklets();
763 if(nTracklets<20) return IsPileupFromSPD(3,0.8);
764 else if(nTracklets<50) return IsPileupFromSPD(4,0.8);
765 else return IsPileupFromSPD(5,0.8);