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"
36 // definition of std AOD member names
37 const char* AliAODEvent::fAODListName[kAODListN] = {"header",
49 //______________________________________________________________________________
50 AliAODEvent::AliAODEvent() :
52 fAODObjects(new TList()),
67 // default constructor
70 //______________________________________________________________________________
71 AliAODEvent::AliAODEvent(const AliAODEvent& aod):
73 fAODObjects(new TList()),
74 fAODFolder(new TFolder()),
76 fHeader(new AliAODHeader(*aod.fHeader)),
77 fTracks(new TClonesArray(*aod.fTracks)),
78 fVertices(new TClonesArray(*aod.fVertices)),
79 fV0s(new TClonesArray(*aod.fV0s)),
80 fTracklets(new AliAODTracklets(*aod.fTracklets)),
81 fJets(new TClonesArray(*aod.fJets)),
82 fEmcalCells(new AliAODCaloCells(*aod.fEmcalCells)),
83 fPhosCells(new AliAODCaloCells(*aod.fPhosCells)),
84 fCaloClusters(new TClonesArray(*aod.fCaloClusters)),
85 fFmdClusters(new TClonesArray(*aod.fFmdClusters)),
86 fPmdClusters(new TClonesArray(*aod.fPmdClusters))
93 AddObject(fTracklets);
95 AddObject(fEmcalCells);
96 AddObject(fPhosCells);
97 AddObject(fCaloClusters);
98 AddObject(fFmdClusters);
99 AddObject(fPmdClusters);
100 fConnected = aod.fConnected;
104 //______________________________________________________________________________
105 AliAODEvent & AliAODEvent::operator=(const AliAODEvent& aod) {
107 // Assignment operator
109 if(&aod == this) return *this;
110 AliVEvent::operator=(aod);
112 fAODObjects = new TList();
113 fAODFolder = new TFolder();
114 fConnected = aod.fConnected;
115 fHeader = new AliAODHeader(*aod.fHeader);
116 fTracks = new TClonesArray(*aod.fTracks);
117 fVertices = new TClonesArray(*aod.fVertices);
118 fV0s = new TClonesArray(*aod.fV0s);
119 fTracklets = new AliAODTracklets(*aod.fTracklets);
120 fJets = new TClonesArray(*aod.fJets);
121 fEmcalCells = new AliAODCaloCells(*aod.fEmcalCells);
122 fPhosCells = new AliAODCaloCells(*aod.fPhosCells);
123 fCaloClusters = new TClonesArray(*aod.fCaloClusters);
124 fFmdClusters = new TClonesArray(*aod.fFmdClusters);
125 fPmdClusters = new TClonesArray(*aod.fPmdClusters);
127 fAODObjects = new TList();
131 AddObject(fVertices);
133 AddObject(fTracklets);
135 AddObject(fEmcalCells);
136 AddObject(fPhosCells);
137 AddObject(fCaloClusters);
138 AddObject(fFmdClusters);
139 AddObject(fPmdClusters);
145 //______________________________________________________________________________
146 AliAODEvent::~AliAODEvent()
149 if(fAODObjects&&!fConnected)
158 //______________________________________________________________________________
159 void AliAODEvent::AddObject(TObject* obj)
161 // Add an object to the list of objects.
162 // Please be aware that in order to increase performance you should
163 // refrain from using TObjArrays (if possible). Use TClonesArrays, instead.
165 fAODObjects->AddLast(obj);
168 //______________________________________________________________________________
169 void AliAODEvent::RemoveObject(TObject* obj)
171 // Removes an object from the list of objects.
173 fAODObjects->Remove(obj);
176 //______________________________________________________________________________
177 TObject *AliAODEvent::FindListObject(const char *objName)
179 // Return the pointer to the object with the given name.
181 return fAODObjects->FindObject(objName);
184 //______________________________________________________________________________
185 void AliAODEvent::CreateStdContent()
187 // create the standard AOD content and set pointers
189 // create standard objects and add them to the TList of objects
190 AddObject(new AliAODHeader());
191 AddObject(new TClonesArray("AliAODTrack", 0));
192 AddObject(new TClonesArray("AliAODVertex", 0));
193 AddObject(new TClonesArray("AliAODv0", 0));
194 AddObject(new AliAODTracklets());
195 AddObject(new TClonesArray("AliAODJet", 0));
196 AddObject(new AliAODCaloCells());
197 AddObject(new AliAODCaloCells());
198 AddObject(new TClonesArray("AliAODCaloCluster", 0));
199 AddObject(new TClonesArray("AliAODFmdCluster", 0));
200 AddObject(new TClonesArray("AliAODPmdCluster", 0));
204 // read back pointers
210 void AliAODEvent::MakeEntriesReferencable()
212 // Make all entries referencable in a subsequent process
214 TIter next(fAODObjects);
216 while ((obj = next()))
218 if(obj->InheritsFrom("TCollection"))
220 AssignIDtoCollection((TCollection*)obj);
225 //______________________________________________________________________________
226 void AliAODEvent::SetStdNames()
228 // introduce the standard naming
230 if(fAODObjects->GetEntries()==kAODListN){
231 for(int i = 0;i < fAODObjects->GetEntries();i++){
232 TObject *fObj = fAODObjects->At(i);
233 if(fObj->InheritsFrom("TNamed")){
234 ((TNamed*)fObj)->SetName(fAODListName[i]);
236 else if(fObj->InheritsFrom("TClonesArray")){
237 ((TClonesArray*)fObj)->SetName(fAODListName[i]);
242 printf("%s:%d SetStdNames() Wrong number of Std Entries \n",(char*)__FILE__,__LINE__);
246 void AliAODEvent::CreateStdFolders()
248 // Create the standard folder structure
249 fAODFolder = gROOT->GetRootFolder()->AddFolder("AOD", "AOD");
250 if(fAODObjects->GetEntries()==kAODListN){
251 for(int i = 0;i < fAODObjects->GetEntries();i++){
252 TObject *fObj = fAODObjects->At(i);
253 if(fObj->InheritsFrom("TClonesArray")){
254 fAODFolder->AddFolder(fAODListName[i], fAODListName[i], (TCollection*) fObj);
256 fAODFolder->AddFolder(fAODListName[i], fAODListName[i], 0);
261 printf("%s:%d CreateStdFolders() Wrong number of Std Entries \n",(char*)__FILE__,__LINE__);
265 //______________________________________________________________________________
266 void AliAODEvent::GetStdContent()
268 // set pointers for standard content
270 fHeader = (AliAODHeader*)fAODObjects->FindObject("header");
271 fTracks = (TClonesArray*)fAODObjects->FindObject("tracks");
272 fVertices = (TClonesArray*)fAODObjects->FindObject("vertices");
273 fV0s = (TClonesArray*)fAODObjects->FindObject("v0s");
274 fTracklets = (AliAODTracklets*)fAODObjects->FindObject("tracklets");
275 fJets = (TClonesArray*)fAODObjects->FindObject("jets");
276 fEmcalCells = (AliAODCaloCells*)fAODObjects->FindObject("emcalCells");
277 fPhosCells = (AliAODCaloCells*)fAODObjects->FindObject("phosCells");
278 fCaloClusters = (TClonesArray*)fAODObjects->FindObject("caloClusters");
279 fFmdClusters = (TClonesArray*)fAODObjects->FindObject("fmdClusters");
280 fPmdClusters = (TClonesArray*)fAODObjects->FindObject("pmdClusters");
283 //______________________________________________________________________________
284 void AliAODEvent::ResetStd(Int_t trkArrSize,
292 // deletes content of standard arrays and resets size
294 if (trkArrSize > fTracks->GetSize())
295 fTracks->Expand(trkArrSize);
298 if (vtxArrSize > fVertices->GetSize())
299 fVertices->Expand(vtxArrSize);
302 if (v0ArrSize > fV0s->GetSize())
303 fV0s->Expand(v0ArrSize);
306 if (jetSize > fJets->GetSize())
307 fJets->Expand(jetSize);
309 fCaloClusters->Delete();
310 if (caloClusSize > fCaloClusters->GetSize())
311 fCaloClusters->Expand(caloClusSize);
313 fFmdClusters->Delete();
314 if (fmdClusSize > fFmdClusters->GetSize())
315 fFmdClusters->Expand(fmdClusSize);
317 fPmdClusters->Delete();
318 if (pmdClusSize > fPmdClusters->GetSize())
319 fPmdClusters->Expand(pmdClusSize);
321 // Reset the tracklets
322 fTracklets->DeleteContainer();
323 fPhosCells->DeleteContainer();
324 fEmcalCells->DeleteContainer();
328 void AliAODEvent::ClearStd()
330 // clears the standard arrays
331 fHeader ->RemoveQTheta();
333 fVertices ->Delete();
335 fTracklets ->DeleteContainer();
337 fEmcalCells ->DeleteContainer();
338 fPhosCells ->DeleteContainer();
339 fCaloClusters ->Delete();
340 fFmdClusters ->Clear();
341 fPmdClusters ->Clear();
344 //_________________________________________________________________
345 Int_t AliAODEvent::GetPHOSClusters(TRefArray *clusters) const
347 // fills the provided TRefArray with all found phos clusters
351 AliAODCaloCluster *cl = 0;
352 Bool_t first = kTRUE;
353 for (Int_t i = 0; i < GetNCaloClusters() ; i++) {
354 if ( (cl = GetCaloCluster(i)) ) {
355 if (cl->IsPHOSCluster()){
357 new (clusters) TRefArray(TProcessID::GetProcessWithUID(cl));
361 //printf("IsPHOS cluster %d, E %2.3f Size: %d \n",i,cl->E(),clusters->GetEntriesFast());
365 return clusters->GetEntriesFast();
368 //_________________________________________________________________
369 Int_t AliAODEvent::GetEMCALClusters(TRefArray *clusters) const
371 // fills the provided TRefArray with all found emcal clusters
374 cout<<"AOD event 1: nclus "<<GetNCaloClusters()<<endl;
375 AliAODCaloCluster *cl = 0;
376 Bool_t first = kTRUE;
377 for (Int_t i = 0; i < GetNCaloClusters(); i++) {
378 if ( (cl = GetCaloCluster(i)) ) {
379 if (cl->IsEMCALCluster()){
381 new (clusters) TRefArray(TProcessID::GetProcessWithUID(cl));
385 //printf("IsEMCal cluster %d, E %2.3f Size: %d \n",i,cl->E(),clusters->GetEntriesFast());
389 return clusters->GetEntriesFast();
393 //______________________________________________________________________________
394 Int_t AliAODEvent::GetMuonTracks(TRefArray *muonTracks) const
396 // fills the provided TRefArray with all found muon tracks
400 AliAODTrack *track = 0;
401 for (Int_t iTrack = 0; iTrack < GetNTracks(); iTrack++) {
402 if ((track = GetTrack(iTrack))->IsMuonTrack()) {
403 muonTracks->Add(track);
407 return muonTracks->GetEntriesFast();
411 void AliAODEvent::ReadFromTree(TTree *tree, Option_t* opt /*= ""*/)
413 // Connects aod event to tree
416 Printf("%s %d AliAODEvent::ReadFromTree() Zero Pointer to Tree \n",(char*)__FILE__,__LINE__);
420 if(!tree->GetTree())tree->LoadTree(0);
422 // Try to find AliAODEvent
423 AliAODEvent *aodEvent = 0;
424 aodEvent = (AliAODEvent*)tree->GetTree()->GetUserInfo()->FindObject("AliAODEvent");
426 // Check if already connected to tree
427 TList* connectedList = (TList*) (tree->GetUserInfo()->FindObject("AODObjectsConnectedToTree"));
428 if (connectedList && (strcmp(opt, "reconnect"))) {
429 // If connected use the connected list of objects
430 printf("Delete and reconnect \n");
432 fAODObjects->Delete();
433 fAODObjects = connectedList;
439 // prevent a memory leak when reading back the TList
440 if (!(strcmp(opt, "reconnect"))) fAODObjects->Delete();
443 // create a new TList from the UserInfo TList...
444 // copy constructor does not work...
445 fAODObjects = (TList*)(aodEvent->GetList()->Clone());
446 fAODObjects->SetOwner(kFALSE);
447 if(fAODObjects->GetEntries()<kAODListN){
448 printf("%s %d AliAODEvent::ReadFromTree() TList contains less than the standard contents %d < %d \n",
449 (char*)__FILE__,__LINE__,fAODObjects->GetEntries(),kAODListN);
452 // Let's find out whether we have friends
453 TList* friendL = tree->GetTree()->GetListOfFriends();
458 while ((fe = (TFriendElement*)next())){
459 aodEvent = (AliAODEvent*)(fe->GetTree()->GetUserInfo()->FindObject("AliAODEvent"));
461 printf("No UserInfo on tree \n");
464 TList* objL = (TList*)(aodEvent->GetList()->Clone());
465 printf("Get list of object from tree %d !!\n", objL->GetEntries());
466 TIter nextobject(objL);
468 while((obj = nextobject()))
470 printf("Adding object from friend %s !\n", obj->GetName());
471 fAODObjects->Add(obj);
472 } // object "branch" loop
478 // set the branch addresses
479 TIter next(fAODObjects);
481 while((el=(TNamed*)next())){
482 TString bname(el->GetName());
483 // check if branch exists under this Name
484 TBranch *br = tree->GetTree()->GetBranch(bname.Data());
486 tree->SetBranchAddress(bname.Data(),fAODObjects->GetObjectRef(el));
488 br = tree->GetBranch(Form("%s.",bname.Data()));
490 tree->SetBranchAddress(Form("%s.",bname.Data()),fAODObjects->GetObjectRef(el));
493 printf("%s %d AliAODEvent::ReadFromTree() No Branch found with Name %s. \n",
494 (char*)__FILE__,__LINE__,bname.Data());
499 // when reading back we are not owner of the list
500 // must not delete it
501 fAODObjects->SetOwner(kFALSE);
502 fAODObjects->SetName("AODObjectsConnectedToTree");
503 // we are not owner of the list objects
504 // must not delete it
505 tree->GetUserInfo()->Add(fAODObjects);
509 // we can't get the list from the user data, create standard content
510 // and set it by hand
512 TIter next(fAODObjects);
514 while((el=(TNamed*)next())){
515 TString bname(el->GetName());
516 tree->SetBranchAddress(bname.Data(),fAODObjects->GetObjectRef(el));
519 // when reading back we are not owner of the list
520 // must not delete it
521 fAODObjects->SetOwner(kFALSE);
525 //______________________________________________________________________________
526 void AliAODEvent::Print(Option_t *) const
528 // Something meaningful should be implemented here.
533 void AliAODEvent::AssignIDtoCollection(TCollection* col)
535 // Static method which assigns a ID to each object in a collection
536 // In this way the objects are marked as referenced and written with
537 // an ID. This has the advantage that TRefs to this objects can be
538 // written by a subsequent process.
541 while ((obj = next()))
542 TProcessID::AssignID(obj);