Number of sigma pedestal cut increased to 4
[u/mrichter/AliRoot.git] / STEER / AliAODEvent.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-2007, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
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  **************************************************************************/
15
16 /* $Id$ */
17
18 //-------------------------------------------------------------------------
19 //     AOD base class
20 //     Author: Markus Oldenburg, CERN
21 //-------------------------------------------------------------------------
22
23 #include <TROOT.h>
24 #include <TTree.h>
25 #include <TFolder.h>
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"
34
35 ClassImp(AliAODEvent)
36
37 // definition of std AOD member names
38   const char* AliAODEvent::fAODListName[kAODListN] = {"header",
39                                                       "tracks",
40                                                       "vertices",
41                                                       "v0s",
42                                                       "cascades",
43                                                       "tracklets",
44                                                       "jets",
45                                                       "emcalCells",
46                                                       "phosCells",
47                                                       "caloClusters",
48                                                       "fmdClusters",
49                                                       "pmdClusters",
50                                                       "dimuons",
51                                                       "AliAODVZERO",
52                                                       "AliAODZDC"
53                                                       
54 };
55 //______________________________________________________________________________
56 AliAODEvent::AliAODEvent() :
57   AliVEvent(),
58   fAODObjects(new TList()),
59   fAODFolder(0),
60   fConnected(kFALSE),
61   fHeader(0),
62   fTracks(0),
63   fVertices(0),
64   fV0s(0),
65   fCascades(0),
66   fTracklets(0),
67   fJets(0),
68   fEmcalCells(0),
69   fPhosCells(0),
70   fCaloClusters(0),
71   fFmdClusters(0),
72   fPmdClusters(0),
73   fDimuons(0),
74   fAODVZERO(0),
75   fAODZDC(0)
76 {
77   // default constructor
78 }
79
80 //______________________________________________________________________________
81 AliAODEvent::AliAODEvent(const AliAODEvent& aod):
82   AliVEvent(aod),
83   fAODObjects(new TList()),
84   fAODFolder(0),
85   fConnected(kFALSE),
86   fHeader(new AliAODHeader(*aod.fHeader)),
87   fTracks(new TClonesArray(*aod.fTracks)),
88   fVertices(new TClonesArray(*aod.fVertices)),
89   fV0s(new TClonesArray(*aod.fV0s)),
90   fCascades(new TClonesArray(*aod.fCascades)),
91   fTracklets(new AliAODTracklets(*aod.fTracklets)),
92   fJets(new TClonesArray(*aod.fJets)),
93   fEmcalCells(new AliAODCaloCells(*aod.fEmcalCells)),
94   fPhosCells(new AliAODCaloCells(*aod.fPhosCells)),
95   fCaloClusters(new TClonesArray(*aod.fCaloClusters)),
96   fFmdClusters(new TClonesArray(*aod.fFmdClusters)),
97   fPmdClusters(new TClonesArray(*aod.fPmdClusters)),
98   fDimuons(new TClonesArray(*aod.fDimuons)),
99   fAODVZERO(new AliAODVZERO(*aod.fAODVZERO)),
100   fAODZDC(new AliAODZDC(*aod.fAODZDC))
101 {
102   // Copy constructor
103   AddObject(fHeader);
104   AddObject(fTracks);
105   AddObject(fVertices);
106   AddObject(fV0s);
107   AddObject(fCascades);
108   AddObject(fTracklets);
109   AddObject(fJets);
110   AddObject(fEmcalCells);
111   AddObject(fPhosCells);
112   AddObject(fCaloClusters);
113   AddObject(fFmdClusters);
114   AddObject(fPmdClusters);
115   AddObject(fDimuons);
116   AddObject(fAODVZERO);
117   AddObject(fAODZDC);
118   fConnected = aod.fConnected;
119   GetStdContent();
120   CreateStdFolders();
121 }
122
123 //______________________________________________________________________________
124 AliAODEvent & AliAODEvent::operator=(const AliAODEvent& aod) {
125
126     // Assignment operator
127
128   if(&aod == this) return *this;
129   AliVEvent::operator=(aod);
130
131   // This assumes that the list is already created
132   // and that the virtual void Copy(Tobject&) function
133   // is correctly implemented in the derived class
134   // otherwise only TObject::Copy() will be used
135
136   if((fAODObjects->GetSize()==0)&&(aod.fAODObjects->GetSize()>=kAODListN)){
137     // We cover the case that we do not yet have the 
138     // standard content but the source has it
139     CreateStdContent();
140   }
141   
142   // Here we have the standard content without user additions, but the content is 
143   // not matching the aod source.
144   
145   // Iterate the list of source objects
146   TIter next(aod.GetList());
147   TObject *its = 0;
148   TString name;
149   while ((its = next())) {
150     name = its->GetName();
151     // Check if we have this object type in out list
152     TObject *mine = fAODObjects->FindObject(name);    
153     if(!mine) {
154       // We have to create the same type of object.
155       TClass* pClass=TClass::GetClass(its->ClassName());     
156       if (!pClass) {
157         AliWarning(Form("Can not find class description for entry %s (%s)\n",
158                    its->ClassName(), name.Data()));
159         continue;
160       }
161       mine=(TObject*)pClass->New();
162       if(!mine){
163         // not in this: can be added to list
164         AliWarning(Form("%s:%d Could not find %s for copying \n",
165                    (char*)__FILE__,__LINE__,name.Data()));
166         continue;
167       }  
168       if(mine->InheritsFrom("TNamed"))  {
169         ((TNamed*)mine)->SetName(name);
170       } else if(mine->InheritsFrom("TCollection")){
171         if(mine->InheritsFrom("TClonesArray")) {
172           TClonesArray *itscl = dynamic_cast<TClonesArray*>(its);
173           if (!itscl) {
174             AliWarning(Form("Class description for entry %s (%s) not TClonesArray\n",
175                    its->ClassName(), name.Data()));
176             continue;
177           
178           }
179                dynamic_cast<TClonesArray*>(mine)->SetClass(itscl->GetClass(), itscl->GetSize());
180         }
181         dynamic_cast<TCollection*>(mine)->SetName(name);
182       }
183       AliDebug(1, Form("adding object %s of type %s", mine->GetName(), mine->ClassName()));
184       AddObject(mine);
185     }
186     // Now we have an object of the same type and name, but different content.        
187     if(!its->InheritsFrom("TCollection")){
188       // simple objects (do they have a Copy method that calls operator= ?)
189       its->Copy(*mine);
190     } else if (its->InheritsFrom("TClonesArray")) {
191       // Create or expand the tclonesarray pointers
192       // so we can directly copy to the object
193       TClonesArray *its_tca = (TClonesArray*)its;
194       TClonesArray *mine_tca = (TClonesArray*)mine;
195       // this leaves the capacity of the TClonesArray the same
196       // except for a factor of 2 increase when size > capacity
197       // does not release any memory occupied by the tca
198       Int_t its_entries = its_tca->GetEntriesFast();
199       mine_tca->ExpandCreate(its_entries);
200       for(int i=0; i<its_entries; i++){
201         // copy 
202         TObject *mine_tca_obj = mine_tca->At(i);
203         TObject *its_tca_obj = its_tca->At(i);
204         // no need to delete first
205         // pointers within the class should be handled by Copy()...
206         // Can there be Empty slots?
207         its_tca_obj->Copy(*mine_tca_obj);
208       }
209     } else {
210       AliWarning(Form("%s:%d cannot copy TCollection \n",
211                       (char*)__FILE__,__LINE__));
212     }
213   }  
214   fConnected = aod.fConnected;
215   return *this;
216 }
217
218
219 //______________________________________________________________________________
220 AliAODEvent::~AliAODEvent() 
221 {
222 // destructor
223     if(!fConnected) delete fAODObjects;
224     delete fAODFolder;
225 }
226
227 //______________________________________________________________________________
228 void AliAODEvent::AddObject(TObject* obj) 
229 {
230   // Add an object to the list of objects.
231   // Please be aware that in order to increase performance you should
232   // refrain from using TObjArrays (if possible). Use TClonesArrays, instead.
233   
234   if ( !fAODObjects->FindObject(obj) ) 
235   {
236     fAODObjects->AddLast(obj);
237   }
238 }
239
240 //______________________________________________________________________________
241 void AliAODEvent::RemoveObject(TObject* obj) 
242 {
243   // Removes an object from the list of objects.
244   
245   fAODObjects->Remove(obj);
246 }
247
248 //______________________________________________________________________________
249 TObject *AliAODEvent::FindListObject(const char *objName) const
250 {
251   // Return the pointer to the object with the given name.
252
253   return fAODObjects->FindObject(objName);
254 }
255
256 //______________________________________________________________________________
257 void AliAODEvent::CreateStdContent() 
258 {
259   // create the standard AOD content and set pointers
260
261   // create standard objects and add them to the TList of objects
262   AddObject(new AliAODHeader());
263   AddObject(new TClonesArray("AliAODTrack", 0));
264   AddObject(new TClonesArray("AliAODVertex", 0));
265   AddObject(new TClonesArray("AliAODv0", 0));
266   AddObject(new TClonesArray("AliAODcascade", 0));
267   AddObject(new AliAODTracklets());
268   AddObject(new TClonesArray("AliAODJet", 0));
269   AddObject(new AliAODCaloCells());
270   AddObject(new AliAODCaloCells());
271   AddObject(new TClonesArray("AliAODCaloCluster", 0));
272   AddObject(new TClonesArray("AliAODFmdCluster", 0));
273   AddObject(new TClonesArray("AliAODPmdCluster", 0));
274   AddObject(new TClonesArray("AliAODDimuon", 0));
275   AddObject(new AliAODVZERO());
276   AddObject(new AliAODZDC());
277   // set names
278   SetStdNames();
279
280   // read back pointers
281   GetStdContent();
282   CreateStdFolders();
283   return;
284 }
285
286 void  AliAODEvent::MakeEntriesReferencable()
287 {
288     // Make all entries referencable in a subsequent process
289     //
290     TIter next(fAODObjects);
291     TObject* obj;
292     while ((obj = next()))
293     {
294         if(obj->InheritsFrom("TCollection"))
295             {
296                 AssignIDtoCollection((TCollection*)obj);
297             }
298     }
299 }
300
301 //______________________________________________________________________________
302 void AliAODEvent::SetStdNames()
303 {
304   // introduce the standard naming
305
306   if(fAODObjects->GetEntries()==kAODListN){
307     for(int i = 0;i < fAODObjects->GetEntries();i++){
308       TObject *fObj = fAODObjects->At(i);
309       if(fObj->InheritsFrom("TNamed")){
310         ((TNamed*)fObj)->SetName(fAODListName[i]);
311       }
312       else if(fObj->InheritsFrom("TClonesArray")){
313         ((TClonesArray*)fObj)->SetName(fAODListName[i]);
314       }
315     }
316   }
317   else{
318     printf("%s:%d SetStdNames() Wrong number of Std Entries \n",(char*)__FILE__,__LINE__);
319   }
320
321
322 void AliAODEvent::CreateStdFolders()
323 {
324     // Create the standard folder structure
325   if(fAODFolder)delete fAODFolder;
326     fAODFolder = gROOT->GetRootFolder()->AddFolder("AOD", "AOD");
327     if(fAODObjects->GetEntries()==kAODListN){
328         for(int i = 0;i < fAODObjects->GetEntries();i++){
329             TObject *fObj = fAODObjects->At(i);
330             if(fObj->InheritsFrom("TClonesArray")){
331                 fAODFolder->AddFolder(fAODListName[i], fAODListName[i], (TCollection*) fObj);
332             } else {
333                 fAODFolder->AddFolder(fAODListName[i], fAODListName[i], 0);
334             }
335         }
336     }
337     else{
338         printf("%s:%d CreateStdFolders() Wrong number of Std Entries \n",(char*)__FILE__,__LINE__);
339     }
340
341
342 //______________________________________________________________________________
343 void AliAODEvent::GetStdContent()
344 {
345   // set pointers for standard content
346
347   fHeader        = (AliAODHeader*)fAODObjects->FindObject("header");
348   fTracks        = (TClonesArray*)fAODObjects->FindObject("tracks");
349   fVertices      = (TClonesArray*)fAODObjects->FindObject("vertices");
350   fV0s           = (TClonesArray*)fAODObjects->FindObject("v0s");
351   fCascades      = (TClonesArray*)fAODObjects->FindObject("cascades");
352   fTracklets     = (AliAODTracklets*)fAODObjects->FindObject("tracklets");
353   fJets          = (TClonesArray*)fAODObjects->FindObject("jets");
354   fEmcalCells    = (AliAODCaloCells*)fAODObjects->FindObject("emcalCells");
355   fPhosCells     = (AliAODCaloCells*)fAODObjects->FindObject("phosCells");
356   fCaloClusters  = (TClonesArray*)fAODObjects->FindObject("caloClusters");
357   fFmdClusters   = (TClonesArray*)fAODObjects->FindObject("fmdClusters");
358   fPmdClusters   = (TClonesArray*)fAODObjects->FindObject("pmdClusters");
359   fDimuons       = (TClonesArray*)fAODObjects->FindObject("dimuons");
360   fAODVZERO      = (AliAODVZERO*)fAODObjects->FindObject("AliAODVZERO");
361   fAODZDC        = (AliAODZDC*)fAODObjects->FindObject("AliAODZDC");
362 }
363
364 //______________________________________________________________________________
365 void AliAODEvent::ResetStd(Int_t trkArrSize, 
366                            Int_t vtxArrSize, 
367                            Int_t v0ArrSize,
368                            Int_t cascadeArrSize,
369                            Int_t jetSize, 
370                            Int_t caloClusSize, 
371                            Int_t fmdClusSize, 
372                            Int_t pmdClusSize,
373                            Int_t dimuonArrSize
374                            )
375 {
376   // deletes content of standard arrays and resets size 
377   
378   if (fTracks) {
379     fTracks->Delete();
380     if (trkArrSize > fTracks->GetSize()) 
381       fTracks->Expand(trkArrSize);
382   }
383   if (fVertices) {
384     fVertices->Delete();
385     if (vtxArrSize > fVertices->GetSize()) 
386       fVertices->Expand(vtxArrSize);
387   }
388   if (fV0s) {
389     fV0s->Delete();
390     if (v0ArrSize > fV0s->GetSize()) 
391       fV0s->Expand(v0ArrSize);
392   }
393   if (fCascades) {
394     fCascades->Delete();
395     if (cascadeArrSize > fCascades->GetSize()) 
396       fCascades->Expand(cascadeArrSize);
397   }
398   if (fJets) {
399     fJets->Delete();
400     if (jetSize > fJets->GetSize())
401       fJets->Expand(jetSize);
402   }
403   if (fCaloClusters) {
404     fCaloClusters->Delete();
405     if (caloClusSize > fCaloClusters->GetSize()) 
406       fCaloClusters->Expand(caloClusSize);
407   }
408   if (fFmdClusters) {
409     fFmdClusters->Delete();
410     if (fmdClusSize > fFmdClusters->GetSize()) 
411       fFmdClusters->Expand(fmdClusSize);
412   }
413   if (fPmdClusters) {
414     fPmdClusters->Delete();
415     if (pmdClusSize > fPmdClusters->GetSize()) 
416       fPmdClusters->Expand(pmdClusSize);
417   }
418   if (fDimuons) {
419     fDimuons->Delete();
420     if (dimuonArrSize > fDimuons->GetSize()) 
421       fDimuons->Expand(dimuonArrSize);
422   }
423   if (fTracklets)
424     fTracklets->DeleteContainer();
425   if (fPhosCells)
426     fPhosCells->DeleteContainer();  
427   if (fEmcalCells)
428     fEmcalCells->DeleteContainer();
429 }
430
431 void AliAODEvent::ClearStd()
432 {
433   // clears the standard arrays
434   if (fHeader)
435     fHeader        ->Clear();
436   if (fTracks)
437     fTracks        ->Delete();
438   if (fVertices)
439     fVertices      ->Delete();
440   if (fV0s)
441     fV0s           ->Delete();
442   if (fCascades)
443     fCascades      ->Delete();
444   if (fTracklets)
445     fTracklets     ->DeleteContainer();
446   if (fJets)
447     fJets          ->Delete();
448   if (fEmcalCells)
449     fEmcalCells    ->DeleteContainer();
450   if (fPhosCells)
451     fPhosCells     ->DeleteContainer();
452   if (fCaloClusters)
453     fCaloClusters  ->Delete();
454   if (fFmdClusters)
455     fFmdClusters   ->Clear();
456   if (fPmdClusters)
457     fPmdClusters   ->Clear();
458   if (fDimuons)
459     fDimuons       ->Clear();
460 }
461
462 //_________________________________________________________________
463 Int_t AliAODEvent::GetPHOSClusters(TRefArray *clusters) const
464 {
465   // fills the provided TRefArray with all found phos clusters
466   
467   clusters->Clear();
468   
469   AliAODCaloCluster *cl = 0;
470   Bool_t first = kTRUE;
471   for (Int_t i = 0; i < GetNumberOfCaloClusters() ; i++) {
472     if ( (cl = GetCaloCluster(i)) ) {
473       if (cl->IsPHOS()){
474         if(first) {
475           new (clusters) TRefArray(TProcessID::GetProcessWithUID(cl)); 
476           first=kFALSE;
477         }
478         clusters->Add(cl);
479         //printf("IsPHOS cluster %d, E %2.3f Size: %d \n",i,cl->E(),clusters->GetEntriesFast());
480       }
481     }
482   }
483   return clusters->GetEntriesFast();
484 }
485
486 //_________________________________________________________________
487 Int_t AliAODEvent::GetEMCALClusters(TRefArray *clusters) const
488 {
489   // fills the provided TRefArray with all found emcal clusters
490
491   clusters->Clear();
492   AliAODCaloCluster *cl = 0;
493   Bool_t first = kTRUE;
494   for (Int_t i = 0; i < GetNumberOfCaloClusters(); i++) {
495     if ( (cl = GetCaloCluster(i)) ) {
496       if (cl->IsEMCAL()){
497         if(first) {
498           new (clusters) TRefArray(TProcessID::GetProcessWithUID(cl)); 
499           first=kFALSE;
500         }
501         clusters->Add(cl);
502         //printf("IsEMCal cluster %d, E %2.3f Size: %d \n",i,cl->E(),clusters->GetEntriesFast());
503       }
504     }
505   }
506   return clusters->GetEntriesFast();
507 }
508
509
510 //______________________________________________________________________________
511 Int_t AliAODEvent::GetMuonTracks(TRefArray *muonTracks) const
512 {
513   // fills the provided TRefArray with all found muon tracks
514
515   muonTracks->Clear();
516
517   AliAODTrack *track = 0;
518   for (Int_t iTrack = 0; iTrack < GetNTracks(); iTrack++) {
519     track = GetTrack(iTrack);
520     if (track->IsMuonTrack()) {
521       muonTracks->Add(track);
522     }
523   }
524   
525   return muonTracks->GetEntriesFast();
526 }
527
528
529 //______________________________________________________________________________
530 Int_t AliAODEvent::GetNumberOfMuonTracks() const
531 {
532   // get number of muon tracks
533   Int_t nMuonTracks=0;
534   for (Int_t iTrack = 0; iTrack < GetNTracks(); iTrack++) {
535     if ((GetTrack(iTrack))->IsMuonTrack()) {
536        nMuonTracks++;
537     }
538   }
539   
540   return nMuonTracks;
541 }
542
543 //______________________________________________________________________________
544 void AliAODEvent::ReadFromTree(TTree *tree, Option_t* opt /*= ""*/)
545 {
546     // Connects aod event to tree
547   
548   if(!tree){
549     AliWarning("Zero Pointer to Tree \n");
550     return;
551   }
552     // load the TTree
553   if(!tree->GetTree())tree->LoadTree(0);
554   
555     // Try to find AliAODEvent
556   AliAODEvent *aodEvent = 0;
557   aodEvent = (AliAODEvent*)tree->GetTree()->GetUserInfo()->FindObject("AliAODEvent");
558   if(aodEvent){
559     // This event is connected to the tree by definition, just say so
560     aodEvent->SetConnected();
561       // Check if already connected to tree
562     TList* connectedList = (TList*) (tree->GetUserInfo()->FindObject("AODObjectsConnectedToTree"));
563     if (connectedList && (strcmp(opt, "reconnect"))) {
564         // If connected use the connected list of objects
565         delete fAODObjects;
566         fAODObjects = connectedList;
567         GetStdContent(); 
568         fConnected = kTRUE;
569         return;
570     } 
571       // Connect to tree
572       // prevent a memory leak when reading back the TList
573 //      if (!(strcmp(opt, "reconnect"))) fAODObjects->Delete();
574     
575       // create a new TList from the UserInfo TList... 
576       // copy constructor does not work...
577     //    fAODObjects = (TList*)(aodEvent->GetList()->Clone());
578     fAODObjects = (TList*)aodEvent->GetList();
579     fAODObjects->SetOwner(kTRUE);
580     if(fAODObjects->GetEntries()<kAODListN)
581     {
582       AliWarning(Form("AliAODEvent::ReadFromTree() TList contains less than the standard contents %d < %d"
583                       " That might be fine though (at least for filtered AODs)",fAODObjects->GetEntries(),kAODListN));
584     }
585       //
586       // Let's find out whether we have friends
587     TList* friendL = tree->GetTree()->GetListOfFriends();
588     if (friendL) 
589     {
590       TIter next(friendL);
591       TFriendElement* fe;
592       while ((fe = (TFriendElement*)next())){
593         aodEvent = (AliAODEvent*)(fe->GetTree()->GetUserInfo()->FindObject("AliAODEvent"));
594         if (!aodEvent) {
595           printf("No UserInfo on tree \n");
596         } else {
597           
598           //          TList* objL = (TList*)(aodEvent->GetList()->Clone());
599           TList* objL = (TList*)aodEvent->GetList();
600           printf("Get list of object from tree %d !!\n", objL->GetEntries());
601           TIter nextobject(objL);
602           TObject* obj =  0;
603           while((obj = nextobject()))
604           {
605             printf("Adding object from friend %s !\n", obj->GetName());
606             fAODObjects->Add(obj);
607           } // object "branch" loop
608         } // has userinfo  
609       } // friend loop
610     } // has friends    
611       // set the branch addresses
612     TIter next(fAODObjects);
613     TNamed *el;
614     while((el=(TNamed*)next())){
615       TString bname(el->GetName());
616         // check if branch exists under this Name
617       TBranch *br = tree->GetTree()->GetBranch(bname.Data());
618       if(br){
619         tree->SetBranchAddress(bname.Data(),fAODObjects->GetObjectRef(el));
620       } else {
621         br = tree->GetBranch(Form("%s.",bname.Data()));
622         if(br){
623           tree->SetBranchAddress(Form("%s.",bname.Data()),fAODObjects->GetObjectRef(el));
624         }
625         else{
626           printf("%s %d AliAODEvent::ReadFromTree() No Branch found with Name %s. \n",
627                  (char*)__FILE__,__LINE__,bname.Data());
628         }       
629       }
630     }
631     GetStdContent();
632       // when reading back we are not owner of the list 
633       // must not delete it
634     fAODObjects->SetOwner(kTRUE);
635     fAODObjects->SetName("AODObjectsConnectedToTree");
636       // we are not owner of the list objects 
637       // must not delete it
638     tree->GetUserInfo()->Add(fAODObjects);
639     fConnected = kTRUE;
640   }// no aodEvent
641   else {
642       // we can't get the list from the user data, create standard content
643       // and set it by hand
644     CreateStdContent();
645     TIter next(fAODObjects);
646     TNamed *el;
647     while((el=(TNamed*)next())){
648       TString bname(el->GetName());    
649       tree->SetBranchAddress(bname.Data(),fAODObjects->GetObjectRef(el));
650     }
651     GetStdContent();
652       // when reading back we are not owner of the list 
653       // must not delete it
654     fAODObjects->SetOwner(kTRUE);
655   }
656 }
657   //______________________________________________________________________________
658 Int_t  AliAODEvent::GetNumberOfPileupVerticesSPD() const{
659   // count number of SPD pileup vertices
660   Int_t nVertices=GetNumberOfVertices();
661   Int_t nPileupVertices=0;
662   for(Int_t iVert=0; iVert<nVertices; iVert++){
663     AliAODVertex *v=GetVertex(iVert);
664     if(v->GetType()==AliAODVertex::kPileupSPD) nPileupVertices++;
665   }
666   return nPileupVertices;
667 }
668 //______________________________________________________________________________
669 Int_t  AliAODEvent::GetNumberOfPileupVerticesTracks() const{
670   // count number of track pileup vertices
671   Int_t nVertices=GetNumberOfVertices();
672   Int_t nPileupVertices=0;
673   for(Int_t iVert=0; iVert<nVertices; iVert++){
674     AliAODVertex *v=GetVertex(iVert);
675     if(v->GetType()==AliAODVertex::kPileupTracks) nPileupVertices++;
676   }
677   return nPileupVertices;
678 }
679 //______________________________________________________________________________
680 AliAODVertex* AliAODEvent::GetPrimaryVertexSPD() const{
681   // Get SPD primary vertex
682   Int_t nVertices=GetNumberOfVertices();
683   for(Int_t iVert=0; iVert<nVertices; iVert++){
684     AliAODVertex *v=GetVertex(iVert);
685     if(v->GetType()==AliAODVertex::kMainSPD) return v;
686   }
687   return 0;
688 }
689 //______________________________________________________________________________
690 AliAODVertex* AliAODEvent::GetPileupVertexSPD(Int_t iV) const{
691   // Get pile-up vertex iV
692   Int_t nVertices=GetNumberOfVertices();
693   Int_t counter=0;
694   for(Int_t iVert=0; iVert<nVertices; iVert++){
695     AliAODVertex *v=GetVertex(iVert);
696     if(v->GetType()==AliAODVertex::kPileupSPD){
697       if(counter==iV) return v;
698       ++counter;
699     }
700   }
701   return 0;
702 }
703 //______________________________________________________________________________
704 AliAODVertex* AliAODEvent::GetPileupVertexTracks(Int_t iV) const{
705   // Get pile-up vertex iV
706   Int_t nVertices=GetNumberOfVertices();
707   Int_t counter=0;
708   for(Int_t iVert=0; iVert<nVertices; iVert++){
709     AliAODVertex *v=GetVertex(iVert);
710     if(v->GetType()==AliAODVertex::kPileupTracks){
711       if(counter==iV) return v;
712       ++counter;
713     }
714   }
715   return 0;
716 }
717 //______________________________________________________________________________
718 Bool_t  AliAODEvent::IsPileupFromSPD(Int_t minContributors, 
719                                      Double_t minZdist, 
720                                      Double_t nSigmaZdist, 
721                                      Double_t nSigmaDiamXY, 
722                                      Double_t nSigmaDiamZ) const{
723   //
724   // This function checks if there was a pile up
725   // reconstructed with SPD
726   //
727   AliAODVertex *mainV=GetPrimaryVertexSPD();
728   if(!mainV) return kFALSE;
729   Int_t nc1=mainV->GetNContributors();
730   if(nc1<1) return kFALSE;
731   Int_t nPileVert=GetNumberOfPileupVerticesSPD();
732   if(nPileVert==0) return kFALSE;
733   Int_t nVertices=GetNumberOfVertices();
734   
735   for(Int_t iVert=0; iVert<nVertices; iVert++){
736     AliAODVertex *pv=GetVertex(iVert);
737     if(pv->GetType()!=AliAODVertex::kPileupSPD) continue;
738     Int_t nc2=pv->GetNContributors();
739     if(nc2>=minContributors){
740       Double_t z1=mainV->GetZ();
741       Double_t z2=pv->GetZ();
742       Double_t distZ=TMath::Abs(z2-z1);
743       Double_t distZdiam=TMath::Abs(z2-GetDiamondZ());
744       Double_t cutZdiam=nSigmaDiamZ*TMath::Sqrt(GetSigma2DiamondZ());
745       if(GetSigma2DiamondZ()<0.0001)cutZdiam=99999.; //protection for missing z diamond information
746       if(distZ>minZdist && distZdiam<cutZdiam){
747         Double_t x2=pv->GetX();
748         Double_t y2=pv->GetY();
749         Double_t distXdiam=TMath::Abs(x2-GetDiamondX());
750         Double_t distYdiam=TMath::Abs(y2-GetDiamondY());
751         Double_t cov1[6],cov2[6];       
752         mainV->GetCovarianceMatrix(cov1);
753         pv->GetCovarianceMatrix(cov2);
754         Double_t errxDist=TMath::Sqrt(cov2[0]+GetSigma2DiamondX());
755         Double_t erryDist=TMath::Sqrt(cov2[2]+GetSigma2DiamondY());
756         Double_t errzDist=TMath::Sqrt(cov1[5]+cov2[5]);
757         Double_t cutXdiam=nSigmaDiamXY*errxDist;
758         if(GetSigma2DiamondX()<0.0001)cutXdiam=99999.; //protection for missing diamond information
759         Double_t cutYdiam=nSigmaDiamXY*erryDist;
760         if(GetSigma2DiamondY()<0.0001)cutYdiam=99999.; //protection for missing diamond information
761         if( (distXdiam<cutXdiam) && (distYdiam<cutYdiam) && (distZ>nSigmaZdist*errzDist) ){
762           return kTRUE;
763         }
764       }
765     }
766   }
767   return kFALSE;
768 }
769
770 //______________________________________________________________________________
771 void AliAODEvent::Print(Option_t *) const
772 {
773   // Print the names of the all branches
774   TIter next(fAODObjects);
775   TNamed *el;
776   Printf(">>>>>  AOD  Content <<<<<");    
777   while((el=(TNamed*)next())){
778     Printf(">> %s ",el->GetName());      
779   }
780   Printf(">>>>>                <<<<<");    
781   
782   return;
783 }
784
785 void AliAODEvent::AssignIDtoCollection(const TCollection* col)
786 {
787     // Static method which assigns a ID to each object in a collection
788     // In this way the objects are marked as referenced and written with 
789     // an ID. This has the advantage that TRefs to this objects can be 
790     // written by a subsequent process.
791     TIter next(col);
792     TObject* obj;
793     while ((obj = next()))
794         TProcessID::AssignID(obj);
795 }
796
797 Bool_t AliAODEvent::IsPileupFromSPDInMultBins() const {
798     Int_t nTracklets=GetTracklets()->GetNumberOfTracklets();
799     if(nTracklets<20) return IsPileupFromSPD(3,0.8);
800     else if(nTracklets<50) return IsPileupFromSPD(4,0.8);
801     else return IsPileupFromSPD(5,0.8);
802 }
803