3b1144aa0e338b80b1dd3767bb1b07d89f03e3fe
[u/mrichter/AliRoot.git] / STEER / ESD / AliESDEvent.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, 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 //           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.
22 /*
23    AliESDEvent *ev= new AliESDEvent();
24    ev->ReadFromTree(esdTree);
25    ...
26     for (Int_t i=0; i<nev; i++) {
27       esdTree->GetEntry(i);
28       if(ev->GetAliESDOld())ev->CopyFromOldESD();
29 */
30 //   The AliESDInputHandler does this automatically for you
31 //
32 // Origin: Christian Klein-Boesing, CERN, Christian.Klein-Boesing@cern.ch
33 //-----------------------------------------------------------------
34
35 #include "TList.h"
36 #include "TRefArray.h"
37 #include <TNamed.h>
38 #include <TROOT.h>
39 #include <TInterpreter.h>
40
41 #include "AliESDEvent.h"
42 #include "AliESDfriend.h"
43 #include "AliESDVZERO.h"
44 #include "AliESDFMD.h"
45 #include "AliESD.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"
63 #include "AliESDv0.h"
64 #include "AliESDFMD.h"
65 #include "AliESDVZERO.h"
66 #include "AliMultiplicity.h"
67 #include "AliRawDataErrorLog.h"
68 #include "AliLog.h"
69 #include "AliESDACORDE.h"
70 #include "AliESDHLTDecision.h"
71 #include "AliCentrality.h"
72 #ifdef MFT_UPGRADE
73 #include "AliESDMFT.h"
74 #endif
75 #include "AliEventplane.h"
76
77
78 ClassImp(AliESDEvent)
79
80
81
82 // here we define the names, some classes are no TNamed, therefore the classnames 
83 // are the Names
84   const char* AliESDEvent::fgkESDListName[kESDListN] = {"AliESDRun",
85                                                         "AliESDHeader",
86                                                         "AliESDZDC",
87                                                         "AliESDFMD",
88                                                         "AliESDVZERO",
89                                                         "AliESDTZERO",
90                                                         "TPCVertex",
91                                                         "SPDVertex",
92                                                         "PrimaryVertex",
93                                                         "AliMultiplicity",
94                                                         "PHOSTrigger",
95                                                         "EMCALTrigger",
96                                                         "SPDPileupVertices",
97                                                         "TrkPileupVertices",
98                                                         "Tracks",
99                                                         "MuonTracks",
100                                                         "PmdTracks",
101                                                         "AliESDTrdTrigger",
102                                                         "TrdTracks",
103                                                         "TrdTracklets",
104                                                         "V0s",
105                                                         "Cascades",
106                                                         "Kinks",
107                                                         "CaloClusters",
108                                                         "EMCALCells",
109                                                         "PHOSCells",
110                                                         "AliRawDataErrorLogs",
111                                                         "AliESDACORDE",
112                                                         "AliTOFHeader"
113                                 #ifdef MFT_UPGRADE
114 //                              , "AliESDMFT"
115                                                         #endif
116   };
117
118 //______________________________________________________________________________
119 AliESDEvent::AliESDEvent():
120   AliVEvent(),
121   fESDObjects(new TList()),
122   fESDRun(0),
123   fHeader(0),
124   fESDZDC(0),
125   fESDFMD(0),
126   fESDVZERO(0),
127   fESDTZERO(0),
128   fTPCVertex(0),
129   fSPDVertex(0),
130   fPrimaryVertex(0),
131   fSPDMult(0),
132   fPHOSTrigger(0),
133   fEMCALTrigger(0),
134   fESDACORDE(0),
135   fTrdTrigger(0),
136   fSPDPileupVertices(0),
137   fTrkPileupVertices(0),
138   fTracks(0),
139   fMuonTracks(0),
140   fPmdTracks(0),
141   fTrdTracks(0),
142   fTrdTracklets(0),
143   fV0s(0),  
144   fCascades(0),
145   fKinks(0),
146   fCaloClusters(0),
147   fEMCALCells(0), fPHOSCells(0),
148   fErrorLogs(0),
149   fESDOld(0),
150   fESDFriendOld(0),
151   fConnected(kFALSE),
152   fUseOwnList(kFALSE),
153   fTOFHeader(0),
154   fCentrality(0),
155   fEventplane(0),
156   fDetectorStatus(0xFFFFFFFF)
157   #ifdef MFT_UPGRADE
158 //  , fESDMFT(0)
159   #endif
160 {
161 }
162 //______________________________________________________________________________
163 AliESDEvent::AliESDEvent(const AliESDEvent& esd):
164   AliVEvent(esd),
165   fESDObjects(new TList()),
166   fESDRun(new AliESDRun(*esd.fESDRun)),
167   fHeader(new AliESDHeader(*esd.fHeader)),
168   fESDZDC(new AliESDZDC(*esd.fESDZDC)),
169   fESDFMD(new AliESDFMD(*esd.fESDFMD)),
170   fESDVZERO(new AliESDVZERO(*esd.fESDVZERO)),
171   fESDTZERO(new AliESDTZERO(*esd.fESDTZERO)),
172   fTPCVertex(new AliESDVertex(*esd.fTPCVertex)),
173   fSPDVertex(new AliESDVertex(*esd.fSPDVertex)),
174   fPrimaryVertex(new AliESDVertex(*esd.fPrimaryVertex)),
175   fSPDMult(new AliMultiplicity(*esd.fSPDMult)),
176   fPHOSTrigger(new AliESDCaloTrigger(*esd.fPHOSTrigger)),
177   fEMCALTrigger(new AliESDCaloTrigger(*esd.fEMCALTrigger)),
178   fESDACORDE(new AliESDACORDE(*esd.fESDACORDE)),
179   fTrdTrigger(new AliESDTrdTrigger(*esd.fTrdTrigger)),
180   fSPDPileupVertices(new TClonesArray(*esd.fSPDPileupVertices)),
181   fTrkPileupVertices(new TClonesArray(*esd.fTrkPileupVertices)),
182   fTracks(new TClonesArray(*esd.fTracks)),
183   fMuonTracks(new TClonesArray(*esd.fMuonTracks)),
184   fPmdTracks(new TClonesArray(*esd.fPmdTracks)),
185   fTrdTracks(new TClonesArray(*esd.fTrdTracks)),
186   fTrdTracklets(new TClonesArray(*esd.fTrdTracklets)),
187   fV0s(new TClonesArray(*esd.fV0s)),  
188   fCascades(new TClonesArray(*esd.fCascades)),
189   fKinks(new TClonesArray(*esd.fKinks)),
190   fCaloClusters(new TClonesArray(*esd.fCaloClusters)),
191   fEMCALCells(new AliESDCaloCells(*esd.fEMCALCells)),
192   fPHOSCells(new AliESDCaloCells(*esd.fPHOSCells)),
193   fErrorLogs(new TClonesArray(*esd.fErrorLogs)),
194   fESDOld(esd.fESDOld ? new AliESD(*esd.fESDOld) : 0),
195   fESDFriendOld(esd.fESDFriendOld ? new AliESDfriend(*esd.fESDFriendOld) : 0),
196   fConnected(esd.fConnected),
197   fUseOwnList(esd.fUseOwnList),
198   fTOFHeader(new AliTOFHeader(*esd.fTOFHeader)),
199   fCentrality(new AliCentrality(*esd.fCentrality)),
200   fEventplane(new AliEventplane(*esd.fEventplane)),
201   fDetectorStatus(esd.fDetectorStatus)
202   #ifdef MFT_UPGRADE
203 //  , fESDMFT(new AliESDMFT(*esd.fESDMFT))
204   #endif
205
206
207 {
208   printf("copying ESD event...\n");   // AU
209   // CKB init in the constructor list and only add here ...
210   AddObject(fESDRun);
211   AddObject(fHeader);
212   AddObject(fESDZDC);
213   AddObject(fESDFMD);
214   AddObject(fESDVZERO);
215   AddObject(fESDTZERO);
216   AddObject(fTPCVertex);
217   AddObject(fSPDVertex);
218   AddObject(fPrimaryVertex);
219   AddObject(fSPDMult);
220   AddObject(fPHOSTrigger);
221   AddObject(fEMCALTrigger);
222   AddObject(fTrdTrigger);
223   AddObject(fSPDPileupVertices);
224   AddObject(fTrkPileupVertices);
225   AddObject(fTracks);
226   AddObject(fMuonTracks);
227   AddObject(fPmdTracks);
228   AddObject(fTrdTracks);
229   AddObject(fTrdTracklets);
230   AddObject(fV0s);
231   AddObject(fCascades);
232   AddObject(fKinks);
233   AddObject(fCaloClusters);
234   AddObject(fEMCALCells);
235   AddObject(fPHOSCells);
236   AddObject(fErrorLogs);
237   AddObject(fESDACORDE);
238   AddObject(fTOFHeader);
239   #ifdef MFT_UPGRADE
240 //  AddObject(fESDMFT);
241   #endif
242   GetStdContent();
243
244 }
245
246 //______________________________________________________________________________
247 AliESDEvent & AliESDEvent::operator=(const AliESDEvent& source) {
248
249   // Assignment operator
250
251   if(&source == this) return *this;
252   AliVEvent::operator=(source);
253
254   // This assumes that the list is already created
255   // and that the virtual void Copy(Tobject&) function
256   // is correctly implemented in the derived class
257   // otherwise only TObject::Copy() will be used
258
259
260
261   if((fESDObjects->GetSize()==0)&&(source.fESDObjects->GetSize()>=kESDListN)){
262     // We cover the case that we do not yet have the 
263     // standard content but the source has it
264     CreateStdContent();
265   }
266
267   TIter next(source.GetList());
268   TObject *its = 0;
269   TString name;
270   while ((its = next())) {
271     name.Form("%s", its->GetName());
272     TObject *mine = fESDObjects->FindObject(name.Data());
273     if(!mine){
274       TClass* pClass=TClass::GetClass(its->ClassName());
275       if (!pClass) {
276         AliWarning(Form("Can not find class description for entry %s (%s)\n",
277                         its->ClassName(), name.Data()));
278         continue;
279       }
280
281       mine=(TObject*)pClass->New();
282       if(!mine){
283       // not in this: can be added to list
284         AliWarning(Form("%s:%d Could not find %s for copying \n",
285                         (char*)__FILE__,__LINE__,name.Data()));
286         continue;
287       }  
288       if(mine->InheritsFrom("TNamed")){
289         ((TNamed*)mine)->SetName(name);
290       }
291       else if(mine->InheritsFrom("TCollection")){
292         if(mine->InheritsFrom("TClonesArray")) {
293           TClonesArray* tcits = dynamic_cast<TClonesArray*>(its);
294           if (tcits)
295             dynamic_cast<TClonesArray*>(mine)->SetClass(tcits->GetClass());
296         }
297         dynamic_cast<TCollection*>(mine)->SetName(name);
298       }
299       AliDebug(1, Form("adding object %s of type %s", mine->GetName(), mine->ClassName()));
300       AddObject(mine);
301     }  
302    
303     if(!its->InheritsFrom("TCollection")){
304       // simple objects
305       its->Copy(*mine);
306     }
307     else if(its->InheritsFrom("TClonesArray")){
308       // Create or expand the tclonesarray pointers
309       // so we can directly copy to the object
310       TClonesArray *its_tca = (TClonesArray*)its;
311       TClonesArray *mine_tca = (TClonesArray*)mine;
312
313       // this leaves the capacity of the TClonesArray the same
314       // except for a factor of 2 increase when size > capacity
315       // does not release any memory occupied by the tca
316       mine_tca->ExpandCreate(its_tca->GetEntriesFast());
317       for(int i = 0;i < its_tca->GetEntriesFast();++i){
318         // copy 
319         TObject *mine_tca_obj = mine_tca->At(i);
320         TObject *its_tca_obj = its_tca->At(i);
321         // no need to delete first
322         // pointers within the class should be handled by Copy()...
323         // Can there be Empty slots?
324         its_tca_obj->Copy(*mine_tca_obj);
325       }
326     }
327     else{
328       AliWarning(Form("%s:%d cannot copy TCollection \n",
329                       (char*)__FILE__,__LINE__));
330     }
331   }
332
333   fCentrality = source.fCentrality;
334   fEventplane = source.fEventplane;
335
336   fConnected  = source.fConnected;
337   fUseOwnList = source.fUseOwnList;
338   
339   fDetectorStatus = source.fDetectorStatus;
340
341   return *this;
342 }
343
344
345 //______________________________________________________________________________
346 AliESDEvent::~AliESDEvent()
347 {
348   //
349   // Standard destructor
350   //
351
352   // everthing on the list gets deleted automatically
353
354   
355   if(fESDObjects&&!fConnected)
356     {
357       delete fESDObjects;
358       fESDObjects = 0;
359     }
360   if (fCentrality) delete fCentrality;
361   if (fEventplane) delete fEventplane;
362   
363 }
364
365 void AliESDEvent::Copy(TObject &obj) const {
366
367   // interface to TOBject::Copy
368   // Copies the content of this into obj!
369   // bascially obj = *this
370
371   if(this==&obj)return;
372   AliESDEvent *robj = dynamic_cast<AliESDEvent*>(&obj);
373   if(!robj)return; // not an AliESEvent
374   *robj = *this;
375   return;
376 }
377
378 //______________________________________________________________________________
379 void AliESDEvent::Reset()
380 {
381
382   // Handle the cases
383   // Std content + Non std content
384
385   // Reset the standard contents
386   ResetStdContent(); 
387
388   //  reset for the old data without AliESDEvent...
389   if(fESDOld)fESDOld->Reset();
390   if(fESDFriendOld){
391     fESDFriendOld->~AliESDfriend();
392     new (fESDFriendOld) AliESDfriend();
393   }
394   // 
395
396   if(fESDObjects->GetSize()>kESDListN){
397     // we have non std content
398     // this also covers esdfriends
399     for(int i = kESDListN;i < fESDObjects->GetSize();++i){
400       TObject *pObject = fESDObjects->At(i);
401       // TClonesArrays
402       if(pObject->InheritsFrom(TClonesArray::Class())){
403         ((TClonesArray*)pObject)->Delete();
404       }
405       else if(!pObject->InheritsFrom(TCollection::Class())){
406         TClass *pClass = TClass::GetClass(pObject->ClassName());
407         if (pClass && pClass->GetListOfMethods()->FindObject("Clear")) {
408           AliDebug(1, Form("Clear for object %s class %s", pObject->GetName(), pObject->ClassName()));
409           pObject->Clear();
410         }
411         else {
412           AliDebug(1, Form("ResetWithPlacementNew for object %s class %s", pObject->GetName(), pObject->ClassName()));
413           ResetWithPlacementNew(pObject);
414         }
415       }
416       else{
417         AliWarning(Form("No reset for %s \n",
418                         pObject->ClassName()));
419       }
420     }
421   }
422
423 }
424
425 Bool_t AliESDEvent::ResetWithPlacementNew(TObject *pObject){
426   Long_t dtoronly = TObject::GetDtorOnly();
427   TClass *pClass = TClass::GetClass(pObject->ClassName()); 
428   TObject::SetDtorOnly(pObject);
429   delete pObject;
430   // Recreate with placement new
431   pClass->New(pObject);
432   // Restore the state.
433   TObject::SetDtorOnly((void*)dtoronly);
434   return kTRUE;
435 }
436
437 void AliESDEvent::ResetStdContent()
438 {
439   // Reset the standard contents
440   if(fESDRun) fESDRun->Reset();
441   if(fHeader) fHeader->Reset();
442   if(fCentrality) fCentrality->Reset();
443   if(fEventplane) fEventplane->Reset();
444   if(fESDZDC) fESDZDC->Reset();
445   if(fESDFMD) {
446     fESDFMD->Clear();
447   }
448   if(fESDVZERO){
449     // reset by callin d'to /c'tor keep the pointer
450     fESDVZERO->~AliESDVZERO();
451     new (fESDVZERO) AliESDVZERO();
452   }  
453   if(fESDACORDE){
454     fESDACORDE->~AliESDACORDE();
455     new (fESDACORDE) AliESDACORDE();    
456   } 
457   if(fESDTZERO) fESDTZERO->Reset(); 
458   // CKB no clear/reset implemented
459   if(fTPCVertex){
460     fTPCVertex->~AliESDVertex();
461     new (fTPCVertex) AliESDVertex();
462     fTPCVertex->SetName(fgkESDListName[kTPCVertex]);
463   }
464   if(fSPDVertex){
465     fSPDVertex->~AliESDVertex();
466     new (fSPDVertex) AliESDVertex();
467     fSPDVertex->SetName(fgkESDListName[kSPDVertex]);
468   }
469   if(fPrimaryVertex){
470     fPrimaryVertex->~AliESDVertex();
471     new (fPrimaryVertex) AliESDVertex();
472     fPrimaryVertex->SetName(fgkESDListName[kPrimaryVertex]);
473   }
474   if(fSPDMult){
475     fSPDMult->~AliMultiplicity();
476     new (fSPDMult) AliMultiplicity();
477   }
478   if(fTOFHeader){
479     fTOFHeader->~AliTOFHeader();
480     new (fTOFHeader) AliTOFHeader();
481     //fTOFHeader->SetName(fgkESDListName[kTOFHeader]);
482   }
483   if (fTrdTrigger) {
484     fTrdTrigger->~AliESDTrdTrigger();
485     new (fTrdTrigger) AliESDTrdTrigger();
486   }
487   #ifdef MFT_UPGRADE
488   //if(fESDMFT){
489 //      fESDMFT->~AliESDMFT();
490 //      new (fESDMFT) AliESDMFT();
491  // }  
492   #endif
493         
494   if(fPHOSTrigger)fPHOSTrigger->DeAllocate(); 
495   if(fEMCALTrigger)fEMCALTrigger->DeAllocate(); 
496   if(fSPDPileupVertices)fSPDPileupVertices->Delete();
497   if(fTrkPileupVertices)fTrkPileupVertices->Delete();
498   if(fTracks)fTracks->Delete();
499   if(fMuonTracks)fMuonTracks->Delete();
500   if(fPmdTracks)fPmdTracks->Delete();
501   if(fTrdTracks)fTrdTracks->Delete();
502   if(fTrdTracklets)fTrdTracklets->Delete();
503   if(fV0s)fV0s->Delete();
504   if(fCascades)fCascades->Delete();
505   if(fKinks)fKinks->Delete();
506   if(fCaloClusters)fCaloClusters->Delete();
507   if(fPHOSCells)fPHOSCells->DeleteContainer();
508   if(fEMCALCells)fEMCALCells->DeleteContainer();
509   if(fErrorLogs) fErrorLogs->Delete();
510
511   // don't reset fconnected fConnected and the list
512
513 }
514
515
516 Int_t AliESDEvent::AddV0(const AliESDv0 *v) {
517   //
518   // Add V0
519   //
520   TClonesArray &fv = *fV0s;
521   Int_t idx=fV0s->GetEntriesFast();
522   new(fv[idx]) AliESDv0(*v);
523   return idx;
524 }  
525
526 //______________________________________________________________________________
527 void AliESDEvent::Print(Option_t *) const 
528 {
529   //
530   // Print header information of the event
531   //
532   printf("ESD run information\n");
533   printf("Event # in file %d Bunch crossing # %d Orbit # %d Period # %d Run # %d Trigger %lld Magnetic field %f \n",
534          GetEventNumberInFile(),
535          GetBunchCrossNumber(),
536          GetOrbitNumber(),
537          GetPeriodNumber(),
538          GetRunNumber(),
539          GetTriggerMask(),
540          GetMagneticField() );
541   if (fPrimaryVertex)
542     printf("Vertex: (%.4f +- %.4f, %.4f +- %.4f, %.4f +- %.4f) cm\n",
543            fPrimaryVertex->GetXv(), fPrimaryVertex->GetXRes(),
544            fPrimaryVertex->GetYv(), fPrimaryVertex->GetYRes(),
545            fPrimaryVertex->GetZv(), fPrimaryVertex->GetZRes());
546   printf("Mean vertex in RUN: X=%.4f Y=%.4f Z=%.4f cm\n",
547          GetDiamondX(),GetDiamondY(),GetDiamondZ());
548   if(fSPDMult)
549     printf("SPD Multiplicity. Number of tracklets %d \n",
550            fSPDMult->GetNumberOfTracklets());
551   printf("Number of pileup primary vertices reconstructed with SPD %d\n", 
552          GetNumberOfPileupVerticesSPD());
553   printf("Number of pileup primary vertices reconstructed using the tracks %d\n",
554          GetNumberOfPileupVerticesTracks());
555   printf("Number of tracks: \n");
556   printf("                 charged   %d\n", GetNumberOfTracks());
557   printf("                 muon      %d\n", GetNumberOfMuonTracks());
558   printf("                 pmd       %d\n", GetNumberOfPmdTracks());
559   printf("                 trd       %d\n", GetNumberOfTrdTracks());
560   printf("                 trd trkl  %d\n", GetNumberOfTrdTracklets());
561   printf("                 v0        %d\n", GetNumberOfV0s());
562   printf("                 cascades  %d\n", GetNumberOfCascades());
563   printf("                 kinks     %d\n", GetNumberOfKinks());
564   if(fPHOSCells)printf("                 PHOSCells %d\n", fPHOSCells->GetNumberOfCells());
565   else printf("                 PHOSCells not in the Event\n");
566   if(fEMCALCells)printf("                 EMCALCells %d\n", fEMCALCells->GetNumberOfCells());
567   else printf("                 EMCALCells not in the Event\n");
568   printf("                 CaloClusters %d\n", GetNumberOfCaloClusters());
569   printf("                 FMD       %s\n", (fESDFMD ? "yes" : "no"));
570   printf("                 VZERO     %s\n", (fESDVZERO ? "yes" : "no"));
571   #ifdef MFT_UPGRADE
572   //printf("                 MFT     %s\n", (fESDMFT ? "yes" : "no"));
573   #endif
574         
575         
576   TObject* pHLTDecision=GetHLTTriggerDecision();
577   printf("HLT trigger decision: %s\n", pHLTDecision?pHLTDecision->GetOption():"not available");
578   if (pHLTDecision) pHLTDecision->Print("compact");
579
580   return;
581 }
582
583 void AliESDEvent::SetESDfriend(const AliESDfriend *ev) const {
584   //
585   // Attaches the complementary info to the ESD
586   //
587   if (!ev) return;
588
589   // to be sure that we set the tracks also
590   // in case of old esds 
591   // if(fESDOld)CopyFromOldESD();
592
593   Int_t ntrk=ev->GetNumberOfTracks();
594  
595   for (Int_t i=0; i<ntrk; i++) {
596     const AliESDfriendTrack *f=ev->GetTrack(i);
597     GetTrack(i)->SetFriendTrack(f);
598   }
599 }
600
601 Bool_t  AliESDEvent::RemoveKink(Int_t rm) const {
602   // ---------------------------------------------------------
603   // Remove a kink candidate and references to it from ESD,
604   // if this candidate does not come from a reconstructed decay
605   // Not yet implemented...
606   // ---------------------------------------------------------
607   Int_t last=GetNumberOfKinks()-1;
608   if ((rm<0)||(rm>last)) return kFALSE;
609
610   return kTRUE;
611 }
612
613 Bool_t  AliESDEvent::RemoveV0(Int_t rm) const {
614   // ---------------------------------------------------------
615   // Remove a V0 candidate and references to it from ESD,
616   // if this candidate does not come from a reconstructed decay
617   // ---------------------------------------------------------
618   Int_t last=GetNumberOfV0s()-1;
619   if ((rm<0)||(rm>last)) return kFALSE;
620
621   AliESDv0 *v0=GetV0(rm);
622   Int_t idxP=v0->GetPindex(), idxN=v0->GetNindex();
623
624   v0=GetV0(last);
625   Int_t lastIdxP=v0->GetPindex(), lastIdxN=v0->GetNindex();
626
627   Int_t used=0;
628
629   // Check if this V0 comes from a reconstructed decay
630   Int_t ncs=GetNumberOfCascades();
631   for (Int_t n=0; n<ncs; n++) {
632     AliESDcascade *cs=GetCascade(n);
633
634     Int_t csIdxP=cs->GetPindex();
635     Int_t csIdxN=cs->GetNindex();
636
637     if (idxP==csIdxP)
638        if (idxN==csIdxN) return kFALSE;
639
640     if (csIdxP==lastIdxP)
641        if (csIdxN==lastIdxN) used++;
642   }
643
644   //Replace the removed V0 with the last V0 
645   TClonesArray &a=*fV0s;
646   delete a.RemoveAt(rm);
647
648   if (rm==last) return kTRUE;
649
650   //v0 is pointing to the last V0 candidate... 
651   new (a[rm]) AliESDv0(*v0);
652   delete a.RemoveAt(last);
653
654   if (!used) return kTRUE;
655   
656
657   // Remap the indices of the daughters of reconstructed decays
658   for (Int_t n=0; n<ncs; n++) {
659     AliESDcascade *cs=GetCascade(n);
660
661
662     Int_t csIdxP=cs->GetPindex();
663     Int_t csIdxN=cs->GetNindex();
664
665     if (csIdxP==lastIdxP)
666       if (csIdxN==lastIdxN) {
667          cs->AliESDv0::SetIndex(1,idxP);
668          cs->AliESDv0::SetIndex(0,idxN);
669          used--;
670          if (!used) return kTRUE;
671       }
672   }
673
674   return kTRUE;
675 }
676
677 Bool_t  AliESDEvent::RemoveTrack(Int_t rm) const {
678   // ---------------------------------------------------------
679   // Remove a track and references to it from ESD,
680   // if this track does not come from a reconstructed decay
681   // ---------------------------------------------------------
682   Int_t last=GetNumberOfTracks()-1;
683   if ((rm<0)||(rm>last)) return kFALSE;
684
685   Int_t used=0;
686
687   // Check if this track comes from the reconstructed primary vertices
688   if (fTPCVertex && fTPCVertex->GetStatus()) {
689      UShort_t *primIdx=fTPCVertex->GetIndices();
690      Int_t n=fTPCVertex->GetNIndices();
691      while (n--) {
692        Int_t idx=Int_t(primIdx[n]);
693        if (rm==idx) return kFALSE;
694        if (idx==last) used++; 
695      }
696   }
697   if (fPrimaryVertex && fPrimaryVertex->GetStatus()) {
698      UShort_t *primIdx=fPrimaryVertex->GetIndices();
699      Int_t n=fPrimaryVertex->GetNIndices();
700      while (n--) {
701        Int_t idx=Int_t(primIdx[n]);
702        if (rm==idx) return kFALSE;
703        if (idx==last) used++; 
704      }
705   }
706   
707   // Check if this track comes from a reconstructed decay
708   Int_t nv0=GetNumberOfV0s();
709   for (Int_t n=0; n<nv0; n++) {
710     AliESDv0 *v0=GetV0(n);
711
712     Int_t idx=v0->GetNindex();
713     if (rm==idx) return kFALSE;
714     if (idx==last) used++;
715
716     idx=v0->GetPindex();
717     if (rm==idx) return kFALSE;
718     if (idx==last) used++;
719   }
720
721   Int_t ncs=GetNumberOfCascades();
722   for (Int_t n=0; n<ncs; n++) {
723     AliESDcascade *cs=GetCascade(n);
724
725     Int_t idx=cs->GetIndex();
726     if (rm==idx) return kFALSE;
727     if (idx==last) used++;
728
729     AliESDv0 *v0=cs;
730     idx=v0->GetNindex();
731     if (rm==idx) return kFALSE;
732     if (idx==last) used++;
733
734     idx=v0->GetPindex();
735     if (rm==idx) return kFALSE;
736     if (idx==last) used++;
737   }
738
739   Int_t nkn=GetNumberOfKinks();
740   for (Int_t n=0; n<nkn; n++) {
741     AliESDkink *kn=GetKink(n);
742
743     Int_t idx=kn->GetIndex(0);
744     if (rm==idx) return kFALSE;
745     if (idx==last) used++;
746
747     idx=kn->GetIndex(1);
748     if (rm==idx) return kFALSE;
749     if (idx==last) used++;
750   }
751
752   // Check if this track is associated with a CaloCluster
753   Int_t ncl=GetNumberOfCaloClusters();
754   for (Int_t n=0; n<ncl; n++) {
755     AliESDCaloCluster *cluster=GetCaloCluster(n);
756     TArrayI *arr=cluster->GetTracksMatched();
757     Int_t s=arr->GetSize();
758     while (s--) {
759       Int_t idx=arr->At(s);
760       if (rm==idx) return kFALSE;
761       if (idx==last) used++;     
762     }
763   }
764
765
766
767   //Replace the removed track with the last track 
768   TClonesArray &a=*fTracks;
769   delete a.RemoveAt(rm);
770
771   if (rm==last) return kTRUE;
772
773   AliESDtrack *t=GetTrack(last);
774   t->SetID(rm);
775   new (a[rm]) AliESDtrack(*t);
776   delete a.RemoveAt(last);
777
778
779   if (!used) return kTRUE;
780   
781
782   // Remap the indices of the tracks used for the primary vertex reconstruction
783   if (fTPCVertex && fTPCVertex->GetStatus()) {
784      UShort_t *primIdx=fTPCVertex->GetIndices();
785      Int_t n=fTPCVertex->GetNIndices();
786      while (n--) {
787        Int_t idx=Int_t(primIdx[n]);
788        if (idx==last) {
789           primIdx[n]=Short_t(rm); 
790           used--;
791           if (!used) return kTRUE;
792        }
793      }
794   }  
795   if (fPrimaryVertex && fPrimaryVertex->GetStatus()) {
796      UShort_t *primIdx=fPrimaryVertex->GetIndices();
797      Int_t n=fPrimaryVertex->GetNIndices();
798      while (n--) {
799        Int_t idx=Int_t(primIdx[n]);
800        if (idx==last) {
801           primIdx[n]=Short_t(rm); 
802           used--;
803           if (!used) return kTRUE;
804        }
805      }
806   }  
807
808   // Remap the indices of the daughters of reconstructed decays
809   for (Int_t n=0; n<nv0; n++) {
810     AliESDv0 *v0=GetV0(n);
811     if (v0->GetIndex(0)==last) {
812        v0->SetIndex(0,rm);
813        used--;
814        if (!used) return kTRUE;
815     }
816     if (v0->GetIndex(1)==last) {
817        v0->SetIndex(1,rm);
818        used--;
819        if (!used) return kTRUE;
820     }
821   }
822
823   for (Int_t n=0; n<ncs; n++) {
824     AliESDcascade *cs=GetCascade(n);
825     if (cs->GetIndex()==last) {
826        cs->SetIndex(rm);
827        used--;
828        if (!used) return kTRUE;
829     }
830     AliESDv0 *v0=cs;
831     if (v0->GetIndex(0)==last) {
832        v0->SetIndex(0,rm);
833        used--;
834        if (!used) return kTRUE;
835     }
836     if (v0->GetIndex(1)==last) {
837        v0->SetIndex(1,rm);
838        used--;
839        if (!used) return kTRUE;
840     }
841   }
842
843   for (Int_t n=0; n<nkn; n++) {
844     AliESDkink *kn=GetKink(n);
845     if (kn->GetIndex(0)==last) {
846        kn->SetIndex(rm,0);
847        used--;
848        if (!used) return kTRUE;
849     }
850     if (kn->GetIndex(1)==last) {
851        kn->SetIndex(rm,1);
852        used--;
853        if (!used) return kTRUE;
854     }
855   }
856
857   // Remap the indices of the tracks accosicated with CaloClusters
858   for (Int_t n=0; n<ncl; n++) {
859     AliESDCaloCluster *cluster=GetCaloCluster(n);
860     TArrayI *arr=cluster->GetTracksMatched();
861     Int_t s=arr->GetSize();
862     while (s--) {
863       Int_t idx=arr->At(s);
864       if (idx==last) {
865          arr->AddAt(rm,s);
866          used--; 
867          if (!used) return kTRUE;
868       }
869     }
870   }
871
872   return kTRUE;
873 }
874
875
876 Bool_t AliESDEvent::Clean(Float_t *cleanPars) {
877   //
878   // Remove the data which are not needed for the physics analysis.
879   //
880   // 1) Cleaning the V0 candidates
881   //    ---------------------------
882   //    If the cosine of the V0 pointing angle "csp" and 
883   //    the DCA between the daughter tracks "dca" does not satisfy 
884   //    the conditions 
885   //
886   //     csp > cleanPars[1] + dca/cleanPars[0]*(1.- cleanPars[1])
887   //
888   //    an attempt to remove this V0 candidate from ESD is made.
889   //
890   //    The V0 candidate gets removed if it does not belong to any 
891   //    recosntructed cascade decay
892   //
893   //    12.11.2007, optimal values: cleanPars[0]=0.5, cleanPars[1]=0.999
894   //
895   // 2) Cleaning the tracks
896   //    ----------------------
897   //    If track's transverse parameter is larger than cleanPars[2]
898   //                       OR
899   //    track's longitudinal parameter is larger than cleanPars[3]
900   //    an attempt to remove this track from ESD is made.
901   //
902   //    The track gets removed if it does not come 
903   //    from a reconstructed decay
904   //
905   Bool_t rc=kFALSE;
906
907   Float_t dcaMax=cleanPars[0];
908   Float_t cspMin=cleanPars[1];
909
910   Int_t nV0s=GetNumberOfV0s();
911   for (Int_t i=nV0s-1; i>=0; i--) {
912     AliESDv0 *v0=GetV0(i);
913
914     Float_t dca=v0->GetDcaV0Daughters();
915     Float_t csp=v0->GetV0CosineOfPointingAngle();
916     Float_t cspcut=cspMin + dca/dcaMax*(1.-cspMin);
917     if (csp > cspcut) continue;
918     if (RemoveV0(i)) rc=kTRUE;
919   }
920
921
922   Float_t dmax=cleanPars[2], zmax=cleanPars[3];
923
924   const AliESDVertex *vertex=GetPrimaryVertexSPD();
925   Bool_t vtxOK=vertex->GetStatus();
926   
927   Int_t nTracks=GetNumberOfTracks();
928   for (Int_t i=nTracks-1; i>=0; i--) {
929     AliESDtrack *track=GetTrack(i);
930     Float_t xy,z; track->GetImpactParameters(xy,z);
931     if ((TMath::Abs(xy) > dmax) || (vtxOK && (TMath::Abs(z) > zmax))) {
932       if (RemoveTrack(i)) rc=kTRUE;
933     }
934   }
935
936   return rc;
937 }
938
939 Char_t  AliESDEvent::AddPileupVertexSPD(const AliESDVertex *vtx) 
940 {
941     // Add a pileup primary vertex reconstructed with SPD
942     TClonesArray &ftr = *fSPDPileupVertices;
943     Char_t n=Char_t(ftr.GetEntriesFast());
944     AliESDVertex *vertex = new(ftr[n]) AliESDVertex(*vtx);
945     vertex->SetID(n);
946     return n;
947 }
948
949 Char_t  AliESDEvent::AddPileupVertexTracks(const AliESDVertex *vtx) 
950 {
951     // Add a pileup primary vertex reconstructed with SPD
952     TClonesArray &ftr = *fTrkPileupVertices;
953     Char_t n=Char_t(ftr.GetEntriesFast());
954     AliESDVertex *vertex = new(ftr[n]) AliESDVertex(*vtx);
955     vertex->SetID(n);
956     return n;
957 }
958
959 Int_t  AliESDEvent::AddTrack(const AliESDtrack *t) 
960 {
961     // Add track
962     TClonesArray &ftr = *fTracks;
963     AliESDtrack * track = new(ftr[fTracks->GetEntriesFast()])AliESDtrack(*t);
964     track->SetID(fTracks->GetEntriesFast()-1);
965     return  track->GetID();    
966 }
967
968 AliESDtrack*  AliESDEvent::NewTrack() 
969 {
970     // Add a new track
971     TClonesArray &ftr = *fTracks;
972     AliESDtrack * track = new(ftr[fTracks->GetEntriesFast()])AliESDtrack();
973     track->SetID(fTracks->GetEntriesFast()-1);
974     return  track;
975 }
976
977  void AliESDEvent::AddMuonTrack(const AliESDMuonTrack *t) 
978 {
979     TClonesArray &fmu = *fMuonTracks;
980     new(fmu[fMuonTracks->GetEntriesFast()]) AliESDMuonTrack(*t);
981 }
982
983 void AliESDEvent::AddPmdTrack(const AliESDPmdTrack *t) 
984 {
985   TClonesArray &fpmd = *fPmdTracks;
986   new(fpmd[fPmdTracks->GetEntriesFast()]) AliESDPmdTrack(*t);
987 }
988
989 void AliESDEvent::SetTrdTrigger(const AliESDTrdTrigger *t)
990 {
991   *fTrdTrigger = *t;
992 }
993
994 void AliESDEvent::AddTrdTrack(const AliESDTrdTrack *t) 
995 {
996   TClonesArray &ftrd = *fTrdTracks;
997   new(ftrd[fTrdTracks->GetEntriesFast()]) AliESDTrdTrack(*t);
998 }
999
1000 void AliESDEvent::AddTrdTracklet(const AliESDTrdTracklet *trkl)
1001 {
1002   new ((*fTrdTracklets)[fTrdTracklets->GetEntriesFast()]) AliESDTrdTracklet(*trkl);
1003 }
1004
1005 Int_t AliESDEvent::AddKink(const AliESDkink *c) 
1006 {
1007     // Add kink
1008     TClonesArray &fk = *fKinks;
1009     AliESDkink * kink = new(fk[fKinks->GetEntriesFast()]) AliESDkink(*c);
1010     kink->SetID(fKinks->GetEntriesFast()); // CKB different from the other imps..
1011     return fKinks->GetEntriesFast()-1;
1012 }
1013
1014
1015 void AliESDEvent::AddCascade(const AliESDcascade *c) 
1016 {
1017   TClonesArray &fc = *fCascades;
1018   new(fc[fCascades->GetEntriesFast()]) AliESDcascade(*c);
1019 }
1020
1021
1022 Int_t AliESDEvent::AddCaloCluster(const AliESDCaloCluster *c) 
1023 {
1024     // Add calocluster
1025     TClonesArray &fc = *fCaloClusters;
1026     AliESDCaloCluster *clus = new(fc[fCaloClusters->GetEntriesFast()]) AliESDCaloCluster(*c);
1027     clus->SetID(fCaloClusters->GetEntriesFast()-1);
1028     return fCaloClusters->GetEntriesFast()-1;
1029   }
1030
1031
1032 void  AliESDEvent::AddRawDataErrorLog(const AliRawDataErrorLog *log) const {
1033   TClonesArray &errlogs = *fErrorLogs;
1034   new(errlogs[errlogs.GetEntriesFast()])  AliRawDataErrorLog(*log);
1035 }
1036
1037 void AliESDEvent::SetZDCData(AliESDZDC * obj)
1038
1039   // use already allocated space
1040   if(fESDZDC)
1041     *fESDZDC = *obj;
1042 }
1043
1044 void  AliESDEvent::SetPrimaryVertexTPC(const AliESDVertex *vertex) 
1045 {
1046   // Set the TPC vertex
1047   // use already allocated space
1048   if(fTPCVertex){
1049     *fTPCVertex = *vertex;
1050     fTPCVertex->SetName(fgkESDListName[kTPCVertex]);
1051   }
1052 }
1053
1054 void  AliESDEvent::SetPrimaryVertexSPD(const AliESDVertex *vertex) 
1055 {
1056   // Set the SPD vertex
1057   // use already allocated space
1058   if(fSPDVertex){
1059     *fSPDVertex = *vertex;
1060     fSPDVertex->SetName(fgkESDListName[kSPDVertex]);
1061   }
1062 }
1063
1064 void  AliESDEvent::SetPrimaryVertexTracks(const AliESDVertex *vertex) 
1065 {
1066   // Set the primary vertex reconstructed using he ESD tracks.
1067   // use already allocated space
1068   if(fPrimaryVertex){
1069     *fPrimaryVertex = *vertex;
1070     fPrimaryVertex->SetName(fgkESDListName[kPrimaryVertex]);
1071   }
1072 }
1073
1074 const AliESDVertex * AliESDEvent::GetPrimaryVertex() const 
1075 {
1076   //
1077   // Get the "best" available reconstructed primary vertex.
1078   //
1079   if(fPrimaryVertex){
1080     if (fPrimaryVertex->GetStatus()) return fPrimaryVertex;
1081   }
1082   if(fSPDVertex){
1083     if (fSPDVertex->GetStatus()) return fSPDVertex;
1084   }
1085   if(fTPCVertex) return fTPCVertex;
1086   
1087   AliWarning("No primary vertex available. Returning the \"default\"...");
1088   return fSPDVertex;
1089 }
1090
1091 AliESDVertex * AliESDEvent::PrimaryVertexTracksUnconstrained() const 
1092 {
1093   //
1094   // Removes diamond constraint from fPrimaryVertex (reconstructed with tracks)
1095   // Returns a AliESDVertex which has to be deleted by the user
1096   //
1097   if(!fPrimaryVertex) {
1098     AliWarning("No primary vertex from tracks available.");
1099     return 0;
1100   }
1101   if(!fPrimaryVertex->GetStatus()) {
1102     AliWarning("No primary vertex from tracks available.");
1103     return 0;
1104   }
1105
1106   AliVertexerTracks vertexer(GetMagneticField());
1107   Float_t diamondxyz[3]={(Float_t)GetDiamondX(),(Float_t)GetDiamondY(),0.};
1108   Float_t diamondcovxy[3]; GetDiamondCovXY(diamondcovxy);
1109   Float_t diamondcov[6]={diamondcovxy[0],diamondcovxy[1],diamondcovxy[2],0.,0.,7.};
1110   AliESDVertex *vertex = 
1111     (AliESDVertex*)vertexer.RemoveConstraintFromVertex(fPrimaryVertex,diamondxyz,diamondcov);
1112
1113   return vertex;
1114 }
1115
1116 void AliESDEvent::SetMultiplicity(const AliMultiplicity *mul) 
1117 {
1118   // Set the SPD Multiplicity
1119   if(fSPDMult){
1120     *fSPDMult = *mul;
1121   }
1122 }
1123
1124
1125 void AliESDEvent::SetFMDData(AliESDFMD * obj) 
1126
1127   // use already allocated space
1128   if(fESDFMD){
1129     *fESDFMD = *obj;
1130   }
1131 }
1132
1133 void AliESDEvent::SetVZEROData(AliESDVZERO * obj)
1134
1135   // use already allocated space
1136   if(fESDVZERO)
1137     *fESDVZERO = *obj;
1138 }
1139
1140 void AliESDEvent::SetTZEROData(AliESDTZERO * obj)
1141
1142   // use already allocated space
1143   if(fESDTZERO)
1144     *fESDTZERO = *obj;
1145 }
1146
1147 #ifdef MFT_UPGRADE
1148 //void AliESDEvent::SetMFTData(AliESDMFT * obj)
1149 //{ 
1150 //  if(fESDMFT)
1151 //      *fESDMFT = *obj;
1152 //}
1153 #endif
1154
1155 void AliESDEvent::SetACORDEData(AliESDACORDE * obj)
1156 {
1157   if(fESDACORDE)
1158     *fESDACORDE = *obj;
1159 }
1160
1161
1162 void AliESDEvent::GetESDfriend(AliESDfriend *ev) const 
1163 {
1164   //
1165   // Extracts the complementary info from the ESD
1166   //
1167   if (!ev) return;
1168
1169   Int_t ntrk=GetNumberOfTracks();
1170
1171   for (Int_t i=0; i<ntrk; i++) {
1172     AliESDtrack *t=GetTrack(i);
1173     const AliESDfriendTrack *f=t->GetFriendTrack();
1174     ev->AddTrack(f);
1175
1176     t->ReleaseESDfriendTrack();// Not to have two copies of "friendTrack"
1177
1178   }
1179
1180   AliESDfriend *fr = (AliESDfriend*)(const_cast<AliESDEvent*>(this)->FindListObject("AliESDfriend"));
1181   if (fr) ev->SetVZEROfriend(fr->GetVZEROfriend());
1182 }
1183
1184 void AliESDEvent::AddObject(TObject* obj) 
1185 {
1186   // Add an object to the list of object.
1187   // Please be aware that in order to increase performance you should
1188   // refrain from using TObjArrays (if possible). Use TClonesArrays, instead.
1189   fESDObjects->SetOwner(kTRUE);
1190   fESDObjects->AddLast(obj);
1191 }
1192
1193
1194 void AliESDEvent::GetStdContent() 
1195 {
1196   // set pointers for standard content
1197   // get by name much safer and not a big overhead since not called very often
1198  
1199   fESDRun = (AliESDRun*)fESDObjects->FindObject(fgkESDListName[kESDRun]);
1200   fHeader = (AliESDHeader*)fESDObjects->FindObject(fgkESDListName[kHeader]);
1201   fESDZDC = (AliESDZDC*)fESDObjects->FindObject(fgkESDListName[kESDZDC]);
1202   fESDFMD = (AliESDFMD*)fESDObjects->FindObject(fgkESDListName[kESDFMD]);
1203   fESDVZERO = (AliESDVZERO*)fESDObjects->FindObject(fgkESDListName[kESDVZERO]);
1204   fESDTZERO = (AliESDTZERO*)fESDObjects->FindObject(fgkESDListName[kESDTZERO]);
1205   fTPCVertex = (AliESDVertex*)fESDObjects->FindObject(fgkESDListName[kTPCVertex]);
1206   fSPDVertex = (AliESDVertex*)fESDObjects->FindObject(fgkESDListName[kSPDVertex]);
1207   fPrimaryVertex = (AliESDVertex*)fESDObjects->FindObject(fgkESDListName[kPrimaryVertex]);
1208   fSPDMult =       (AliMultiplicity*)fESDObjects->FindObject(fgkESDListName[kSPDMult]);
1209   fPHOSTrigger = (AliESDCaloTrigger*)fESDObjects->FindObject(fgkESDListName[kPHOSTrigger]);
1210   fEMCALTrigger = (AliESDCaloTrigger*)fESDObjects->FindObject(fgkESDListName[kEMCALTrigger]);
1211   fSPDPileupVertices = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kSPDPileupVertices]);
1212   fTrkPileupVertices = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kTrkPileupVertices]);
1213   fTracks = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kTracks]);
1214   fMuonTracks = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kMuonTracks]);
1215   fPmdTracks = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kPmdTracks]);
1216   fTrdTrigger = (AliESDTrdTrigger*)fESDObjects->FindObject(fgkESDListName[kTrdTrigger]);
1217   fTrdTracks = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kTrdTracks]);
1218   fTrdTracklets = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kTrdTracklets]);
1219   fV0s = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kV0s]);
1220   fCascades = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kCascades]);
1221   fKinks = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kKinks]);
1222   fCaloClusters = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kCaloClusters]);
1223   fEMCALCells = (AliESDCaloCells*)fESDObjects->FindObject(fgkESDListName[kEMCALCells]);
1224   fPHOSCells = (AliESDCaloCells*)fESDObjects->FindObject(fgkESDListName[kPHOSCells]);
1225   fErrorLogs = (TClonesArray*)fESDObjects->FindObject(fgkESDListName[kErrorLogs]);
1226   fESDACORDE = (AliESDACORDE*)fESDObjects->FindObject(fgkESDListName[kESDACORDE]);
1227   fTOFHeader = (AliTOFHeader*)fESDObjects->FindObject(fgkESDListName[kTOFHeader]);
1228   #ifdef MFT_UPGRADE
1229  // fESDMFT = (AliESDMFT*)fESDObjects->FindObject(fgkESDListName[kESDMFT]);
1230   #endif
1231 }
1232
1233 void AliESDEvent::SetStdNames(){
1234   // Set the names of the standard contents
1235   // 
1236   if(fESDObjects->GetEntries()>=kESDListN){
1237     for(int i = 0;i < fESDObjects->GetEntries() && i<kESDListN;i++){
1238       TObject *fObj = fESDObjects->At(i);
1239       if(fObj->InheritsFrom("TNamed")){
1240         ((TNamed*)fObj)->SetName(fgkESDListName[i]);
1241       }
1242       else if(fObj->InheritsFrom("TClonesArray")){
1243         ((TClonesArray*)fObj)->SetName(fgkESDListName[i]);
1244       }
1245     }
1246   }
1247   else{
1248      AliWarning("Std Entries missing");
1249   }
1250
1251
1252
1253 void AliESDEvent::CreateStdContent(Bool_t bUseThisList){
1254   fUseOwnList = bUseThisList;
1255   CreateStdContent();
1256 }
1257
1258 void AliESDEvent::CreateStdContent() 
1259 {
1260   // create the standard AOD content and set pointers
1261
1262   // create standard objects and add them to the TList of objects
1263   AddObject(new AliESDRun());
1264   AddObject(new AliESDHeader());
1265   AddObject(new AliESDZDC());
1266   AddObject(new AliESDFMD());
1267   AddObject(new AliESDVZERO());
1268   AddObject(new AliESDTZERO());
1269   AddObject(new AliESDVertex());
1270   AddObject(new AliESDVertex());
1271   AddObject(new AliESDVertex());
1272   AddObject(new AliMultiplicity());
1273   AddObject(new AliESDCaloTrigger());
1274   AddObject(new AliESDCaloTrigger());
1275   AddObject(new TClonesArray("AliESDVertex",0));
1276   AddObject(new TClonesArray("AliESDVertex",0));
1277   AddObject(new TClonesArray("AliESDtrack",0));
1278   AddObject(new TClonesArray("AliESDMuonTrack",0));
1279   AddObject(new TClonesArray("AliESDPmdTrack",0));
1280   AddObject(new AliESDTrdTrigger());
1281   AddObject(new TClonesArray("AliESDTrdTrack",0));
1282   AddObject(new TClonesArray("AliESDTrdTracklet",0));
1283   AddObject(new TClonesArray("AliESDv0",0));
1284   AddObject(new TClonesArray("AliESDcascade",0));
1285   AddObject(new TClonesArray("AliESDkink",0));
1286   AddObject(new TClonesArray("AliESDCaloCluster",0));
1287   AddObject(new AliESDCaloCells());
1288   AddObject(new AliESDCaloCells());
1289   AddObject(new TClonesArray("AliRawDataErrorLog",0));
1290   AddObject(new AliESDACORDE()); 
1291   AddObject(new AliTOFHeader());
1292   #ifdef MFT_UPGRADE
1293   //AddObject(new AliESDMFT());
1294   #endif
1295         
1296   // check the order of the indices against enum...
1297
1298   // set names
1299   SetStdNames();
1300   // read back pointers
1301   GetStdContent();
1302 }
1303
1304 TObject* AliESDEvent::FindListObject(const char *name) const {
1305 //
1306 // Find object with name "name" in the list of branches
1307 //
1308   if(fESDObjects){
1309     return fESDObjects->FindObject(name);
1310   }
1311   return 0;
1312
1313
1314 Int_t AliESDEvent::GetPHOSClusters(TRefArray *clusters) const
1315 {
1316   // fills the provided TRefArray with all found phos clusters
1317   
1318   clusters->Clear();
1319   
1320   AliESDCaloCluster *cl = 0;
1321   for (Int_t i = 0; i < GetNumberOfCaloClusters(); i++) {
1322     
1323     if ( (cl = GetCaloCluster(i)) ) {
1324       if (cl->IsPHOS()){
1325         clusters->Add(cl);
1326         AliDebug(1,Form("IsPHOS cluster %d Size: %d \n",i,clusters->GetEntriesFast()));
1327       }
1328     }
1329   }
1330   return clusters->GetEntriesFast();
1331 }
1332
1333 Int_t AliESDEvent::GetEMCALClusters(TRefArray *clusters) const
1334 {
1335   // fills the provided TRefArray with all found emcal clusters
1336
1337   clusters->Clear();
1338
1339   AliESDCaloCluster *cl = 0;
1340   for (Int_t i = 0; i < GetNumberOfCaloClusters(); i++) {
1341
1342     if ( (cl = GetCaloCluster(i)) ) {
1343       if (cl->IsEMCAL()){
1344         clusters->Add(cl);
1345         AliDebug(1,Form("IsEMCAL cluster %d Size: %d \n",i,clusters->GetEntriesFast()));
1346       }
1347     }
1348   }
1349   return clusters->GetEntriesFast();
1350 }
1351
1352 void AliESDEvent::WriteToTree(TTree* tree) const {
1353   // Book the branches as in TTree::Branch(TCollection*)
1354   // but add a "." at the end of top level branches which are
1355   // not a TClonesArray
1356
1357
1358   TString branchname;
1359   TIter next(fESDObjects);
1360   const Int_t kSplitlevel = 99; // default value in TTree::Branch()
1361   const Int_t kBufsize = 32000; // default value in TTree::Branch()
1362   TObject *obj = 0;
1363
1364   while ((obj = next())) {
1365     branchname.Form("%s", obj->GetName());
1366     if(branchname.CompareTo("AliESDfriend")==0)branchname = "ESDfriend.";
1367     if ((kSplitlevel > 1) &&  !obj->InheritsFrom(TClonesArray::Class())) {
1368       if(!branchname.EndsWith("."))branchname += ".";
1369     }
1370     if (!tree->FindBranch(branchname)) {
1371       // For the custom streamer to be called splitlevel
1372       // has to be negative, only needed for HLT
1373       Int_t splitLevel = (TString(obj->ClassName()) == "AliHLTGlobalTriggerDecision") ? -1 : kSplitlevel - 1;
1374       tree->Bronch(branchname, obj->ClassName(), fESDObjects->GetObjectRef(obj),kBufsize, splitLevel);
1375     }
1376   }
1377 }
1378
1379
1380 void AliESDEvent::ReadFromTree(TTree *tree, Option_t* opt){
1381 //
1382 // Connect the ESDEvent to a tree
1383 //
1384   if(!tree){
1385     AliWarning("AliESDEvent::ReadFromTree() Zero Pointer to Tree \n");
1386     return;
1387   }
1388   // load the TTree
1389   if(!tree->GetTree())tree->LoadTree(0);
1390
1391   // if we find the "ESD" branch on the tree we do have the old structure
1392   if(tree->GetBranch("ESD")) {
1393     char ** address  = (char **)(tree->GetBranch("ESD")->GetAddress());
1394     // do we have the friend branch
1395     TBranch * esdFB = tree->GetBranch("ESDfriend.");
1396     char ** addressF = 0;
1397     if(esdFB)addressF = (char **)(esdFB->GetAddress());
1398     if (!address) {
1399       AliInfo("AliESDEvent::ReadFromTree() Reading old Tree");
1400       tree->SetBranchAddress("ESD",       &fESDOld);
1401       if(esdFB){
1402         tree->SetBranchAddress("ESDfriend.",&fESDFriendOld);
1403       }
1404     } else {
1405       AliInfo("AliESDEvent::ReadFromTree() Reading old Tree");
1406       AliInfo("Branch already connected. Using existing branch address.");
1407       fESDOld       = (AliESD*)       (*address);
1408       // addressF can still be 0, since branch needs to switched on
1409       if(addressF)fESDFriendOld = (AliESDfriend*) (*addressF);
1410     }
1411                                        
1412     //  have already connected the old ESD structure... ?
1413     // reuse also the pointer of the AlliESDEvent
1414     // otherwise create new ones
1415     TList* connectedList = (TList*) (tree->GetUserInfo()->FindObject("ESDObjectsConnectedToTree"));
1416   
1417     if(connectedList){
1418       // If connected use the connected list of objects
1419       if(fESDObjects!= connectedList){
1420         // protect when called twice 
1421         fESDObjects->Delete();
1422         fESDObjects = connectedList;
1423       }
1424       GetStdContent(); 
1425
1426       
1427       // The pointer to the friend changes when called twice via InitIO
1428       // since AliESDEvent is deleted
1429       TObject* oldf = FindListObject("AliESDfriend");
1430       TObject* newf = 0;
1431       if(addressF){
1432         newf = (TObject*)*addressF;
1433       }
1434       if(newf!=0&&oldf!=newf){
1435         // remove the old reference
1436         // Should we also delete it? Or is this handled in TTree I/O
1437         // since it is created by the first SetBranchAddress
1438         fESDObjects->Remove(oldf);
1439         // add the new one 
1440         fESDObjects->Add(newf);
1441       }
1442       
1443       fConnected = true;
1444       return;
1445     }
1446     // else...    
1447     CreateStdContent(); // create for copy
1448     // if we have the esdfriend add it, so we always can access it via the userinfo
1449     if(fESDFriendOld)AddObject(fESDFriendOld);
1450     // we are not owner of the list objects 
1451     // must not delete it
1452     fESDObjects->SetOwner(kTRUE);
1453     fESDObjects->SetName("ESDObjectsConnectedToTree");
1454     tree->GetUserInfo()->Add(fESDObjects);
1455     fConnected = true;
1456     return;
1457   }
1458   
1459
1460     delete fESDOld;
1461     fESDOld = 0;
1462   // Try to find AliESDEvent
1463   AliESDEvent *esdEvent = 0;
1464   esdEvent = (AliESDEvent*)tree->GetTree()->GetUserInfo()->FindObject("AliESDEvent");
1465   if(esdEvent){   
1466       // Check if already connected to tree
1467     esdEvent->Reset();
1468     TList* connectedList = (TList*) (tree->GetUserInfo()->FindObject("ESDObjectsConnectedToTree"));
1469
1470     
1471     if (connectedList && (strcmp(opt, "reconnect"))) {
1472       // If connected use the connected list if objects
1473       fESDObjects->Delete();
1474       fESDObjects = connectedList;
1475       GetStdContent(); 
1476       fConnected = true;
1477       return;
1478     }
1479
1480     // Connect to tree
1481     // prevent a memory leak when reading back the TList
1482     // if (!(strcmp(opt, "reconnect"))) fESDObjects->Delete();
1483     
1484     if(!fUseOwnList){
1485       // create a new TList from the UserInfo TList... 
1486       // copy constructor does not work...
1487       fESDObjects = (TList*)(esdEvent->GetList()->Clone());
1488       fESDObjects->SetOwner(kTRUE);
1489     }
1490     else if ( fESDObjects->GetEntries()==0){
1491       // at least create the std content if we want to read to our list
1492       CreateStdContent(); 
1493     }
1494
1495     // in principle
1496     // we only need new things in the list if we do no already have it..
1497     // TODO just add new entries
1498
1499     if(fESDObjects->GetEntries()<kESDListN){
1500       AliWarning(Form("AliESDEvent::ReadFromTree() TList contains less than the standard contents %d < %d \n",
1501                       fESDObjects->GetEntries(),kESDListN));
1502     }
1503     // set the branch addresses
1504     TIter next(fESDObjects);
1505     TNamed *el;
1506     while((el=(TNamed*)next())){
1507       TString bname(el->GetName());
1508       if(bname.CompareTo("AliESDfriend")==0)
1509         {
1510           // AliESDfriend does not have a name ...
1511             TBranch *br = tree->GetBranch("ESDfriend.");
1512             if (br) tree->SetBranchAddress("ESDfriend.",fESDObjects->GetObjectRef(el));
1513         }
1514       else{
1515         // check if branch exists under this Name
1516         TBranch *br = tree->GetBranch(bname.Data());
1517         if(br){
1518           tree->SetBranchAddress(bname.Data(),fESDObjects->GetObjectRef(el));
1519         }
1520         else{
1521           br = tree->GetBranch(Form("%s.",bname.Data()));
1522           if(br){
1523             tree->SetBranchAddress(Form("%s.",bname.Data()),fESDObjects->GetObjectRef(el));
1524           }
1525           else{
1526             AliWarning(Form("AliESDEvent::ReadFromTree() No Branch found with Name %s or %s.",bname.Data(),bname.Data()));
1527           }
1528
1529         }
1530       }
1531     }
1532     GetStdContent();
1533     // when reading back we are not owner of the list 
1534     // must not delete it
1535     fESDObjects->SetOwner(kTRUE);
1536     fESDObjects->SetName("ESDObjectsConnectedToTree");
1537     // we are not owner of the list objects 
1538     // must not delete it
1539     tree->GetUserInfo()->Add(fESDObjects);
1540     tree->GetUserInfo()->SetOwner(kFALSE);
1541     fConnected = true;
1542   }// no esdEvent -->
1543   else {
1544     // we can't get the list from the user data, create standard content
1545     // and set it by hand (no ESDfriend at the moment
1546     CreateStdContent();
1547     TIter next(fESDObjects);
1548     TNamed *el;
1549     while((el=(TNamed*)next())){
1550       TString bname(el->GetName());    
1551       TBranch *br = tree->GetBranch(bname.Data());
1552       if(br){
1553         tree->SetBranchAddress(bname.Data(),fESDObjects->GetObjectRef(el));
1554       }
1555       else{
1556         br = tree->GetBranch(Form("%s.",bname.Data()));
1557         if(br){
1558           tree->SetBranchAddress(Form("%s.",bname.Data()),fESDObjects->GetObjectRef(el));
1559         }
1560       }
1561     }
1562     GetStdContent();
1563     // when reading back we are not owner of the list 
1564     // must not delete it
1565     fESDObjects->SetOwner(kTRUE);
1566   }
1567 }
1568
1569
1570 void AliESDEvent::CopyFromOldESD()
1571 {
1572   // Method which copies over everthing from the old esd structure to the 
1573   // new  
1574   if(fESDOld){
1575     ResetStdContent();
1576      // Run
1577     SetRunNumber(fESDOld->GetRunNumber());
1578     SetPeriodNumber(fESDOld->GetPeriodNumber());
1579     SetMagneticField(fESDOld->GetMagneticField());
1580   
1581     // leave out diamond ...
1582     // SetDiamond(const AliESDVertex *vertex) { fESDRun->SetDiamond(vertex);}
1583
1584     // header
1585     SetTriggerMask(fESDOld->GetTriggerMask());
1586     SetOrbitNumber(fESDOld->GetOrbitNumber());
1587     SetTimeStamp(fESDOld->GetTimeStamp());
1588     SetEventType(fESDOld->GetEventType());
1589     SetEventNumberInFile(fESDOld->GetEventNumberInFile());
1590     SetBunchCrossNumber(fESDOld->GetBunchCrossNumber());
1591     SetTriggerCluster(fESDOld->GetTriggerCluster());
1592
1593     // ZDC
1594
1595     SetZDC(fESDOld->GetZDCN1Energy(),
1596            fESDOld->GetZDCP1Energy(),
1597            fESDOld->GetZDCEMEnergy(),
1598            0,
1599            fESDOld->GetZDCN2Energy(),
1600            fESDOld->GetZDCP2Energy(),
1601            fESDOld->GetZDCParticipants(),
1602            0,
1603            0,
1604            0,
1605            0,
1606            0,
1607            0);
1608
1609     // FMD
1610     
1611     if(fESDOld->GetFMDData())SetFMDData(fESDOld->GetFMDData());
1612
1613     // T0
1614
1615     SetT0zVertex(fESDOld->GetT0zVertex());
1616     SetT0(fESDOld->GetT0());
1617     //  leave amps out
1618
1619     // VZERO
1620     if (fESDOld->GetVZEROData()) SetVZEROData(fESDOld->GetVZEROData());
1621
1622     if(fESDOld->GetVertex())SetPrimaryVertexSPD(fESDOld->GetVertex());
1623
1624     if(fESDOld->GetPrimaryVertex())SetPrimaryVertexTracks(fESDOld->GetPrimaryVertex());
1625
1626     if(fESDOld->GetMultiplicity())SetMultiplicity(fESDOld->GetMultiplicity());
1627
1628     for(int i = 0;i<fESDOld->GetNumberOfTracks();i++){
1629       AddTrack(fESDOld->GetTrack(i));
1630     }
1631
1632     for(int i = 0;i<fESDOld->GetNumberOfMuonTracks();i++){
1633       AddMuonTrack(fESDOld->GetMuonTrack(i));
1634     }
1635
1636     for(int i = 0;i<fESDOld->GetNumberOfPmdTracks();i++){
1637       AddPmdTrack(fESDOld->GetPmdTrack(i));
1638     }
1639
1640     for(int i = 0;i<fESDOld->GetNumberOfTrdTracks();i++){
1641       AddTrdTrack(fESDOld->GetTrdTrack(i));
1642     }
1643
1644     for(int i = 0;i<fESDOld->GetNumberOfV0s();i++){
1645       AddV0(fESDOld->GetV0(i));
1646     }
1647
1648     for(int i = 0;i<fESDOld->GetNumberOfCascades();i++){
1649       AddCascade(fESDOld->GetCascade(i));
1650     }
1651
1652     for(int i = 0;i<fESDOld->GetNumberOfKinks();i++){
1653       AddKink(fESDOld->GetKink(i));
1654     }
1655
1656
1657     for(int i = 0;i<fESDOld->GetNumberOfCaloClusters();i++){
1658       AddCaloCluster(fESDOld->GetCaloCluster(i));
1659     }
1660           
1661         #ifdef MFT_UPGRADE  
1662         // MFT
1663 //      if (fESDOld->GetMFTData()) SetMFTData(fESDOld->GetMFTData());
1664     #endif
1665
1666   }// if fesdold
1667 }
1668
1669 Bool_t AliESDEvent::IsEventSelected(const char *trigExpr) const
1670 {
1671   // Check if the event satisfies the trigger
1672   // selection expression trigExpr.
1673   // trigExpr can be any logical expression
1674   // of the trigger classes defined in AliESDRun
1675   // In case of wrong syntax return kTRUE.
1676
1677   TString expr(trigExpr);
1678   if (expr.IsNull()) return kTRUE;
1679
1680   ULong64_t mask = GetTriggerMask();
1681   for(Int_t itrig = 0; itrig < AliESDRun::kNTriggerClasses; itrig++) {
1682     if (mask & (1ull << itrig)) {
1683       expr.ReplaceAll(GetESDRun()->GetTriggerClass(itrig),"1");
1684     }
1685     else {
1686       expr.ReplaceAll(GetESDRun()->GetTriggerClass(itrig),"0");
1687     }
1688   }
1689
1690   Int_t error;
1691   if ((gROOT->ProcessLineFast(expr.Data(),&error) == 0) &&
1692       (error == TInterpreter::kNoError)) {
1693     return kFALSE;
1694   }
1695
1696   return kTRUE;
1697
1698 }
1699
1700 TObject*  AliESDEvent::GetHLTTriggerDecision() const
1701 {
1702   // get the HLT trigger decission object
1703
1704   // cast away const'nes because the FindListObject method
1705   // is not const
1706   AliESDEvent* pNonConst=const_cast<AliESDEvent*>(this);
1707   return pNonConst->FindListObject("HLTGlobalTrigger");
1708 }
1709
1710 TString   AliESDEvent::GetHLTTriggerDescription() const
1711 {
1712   // get the HLT trigger decission description
1713   TString description;
1714   TObject* pDecision=GetHLTTriggerDecision();
1715   if (pDecision) {
1716     description=pDecision->GetTitle();
1717   }
1718
1719   return description;
1720 }
1721
1722 Bool_t    AliESDEvent::IsHLTTriggerFired(const char* name) const
1723 {
1724   // get the HLT trigger decission description
1725   TObject* pDecision=GetHLTTriggerDecision();
1726   if (!pDecision) return kFALSE;
1727
1728   Option_t* option=pDecision->GetOption();
1729   if (option==NULL || *option!='1') return kFALSE;
1730
1731   if (name) {
1732     TString description=GetHLTTriggerDescription();
1733     Int_t index=description.Index(name);
1734     if (index<0) return kFALSE;
1735     index+=strlen(name);
1736     if (index>=description.Length()) return kFALSE;
1737     if (description[index]!=0 && description[index]!=' ') return kFALSE;
1738   }
1739   return kTRUE;
1740 }
1741
1742 //______________________________________________________________________________
1743 Bool_t  AliESDEvent::IsPileupFromSPD(Int_t minContributors, 
1744                                      Double_t minZdist, 
1745                                      Double_t nSigmaZdist, 
1746                                      Double_t nSigmaDiamXY, 
1747                                      Double_t nSigmaDiamZ) const{
1748   //
1749   // This function checks if there was a pile up
1750   // reconstructed with SPD
1751   //
1752   Int_t nc1=fSPDVertex->GetNContributors();
1753   if(nc1<1) return kFALSE;
1754   Int_t nPileVert=GetNumberOfPileupVerticesSPD();
1755   if(nPileVert==0) return kFALSE;
1756   
1757   for(Int_t i=0; i<nPileVert;i++){
1758     const AliESDVertex* pv=GetPileupVertexSPD(i);
1759     Int_t nc2=pv->GetNContributors();
1760     if(nc2>=minContributors){
1761       Double_t z1=fSPDVertex->GetZ();
1762       Double_t z2=pv->GetZ();
1763       Double_t distZ=TMath::Abs(z2-z1);
1764       Double_t distZdiam=TMath::Abs(z2-GetDiamondZ());
1765       Double_t cutZdiam=nSigmaDiamZ*TMath::Sqrt(GetSigma2DiamondZ());
1766       if(GetSigma2DiamondZ()<0.0001)cutZdiam=99999.; //protection for missing z diamond information
1767       if(distZ>minZdist && distZdiam<cutZdiam){
1768         Double_t x2=pv->GetX();
1769         Double_t y2=pv->GetY();
1770         Double_t distXdiam=TMath::Abs(x2-GetDiamondX());
1771         Double_t distYdiam=TMath::Abs(y2-GetDiamondY());
1772         Double_t cov1[6],cov2[6];       
1773         fSPDVertex->GetCovarianceMatrix(cov1);
1774         pv->GetCovarianceMatrix(cov2);
1775         Double_t errxDist=TMath::Sqrt(cov2[0]+GetSigma2DiamondX());
1776         Double_t erryDist=TMath::Sqrt(cov2[2]+GetSigma2DiamondY());
1777         Double_t errzDist=TMath::Sqrt(cov1[5]+cov2[5]);
1778         Double_t cutXdiam=nSigmaDiamXY*errxDist;
1779         if(GetSigma2DiamondX()<0.0001)cutXdiam=99999.; //protection for missing diamond information
1780         Double_t cutYdiam=nSigmaDiamXY*erryDist;
1781         if(GetSigma2DiamondY()<0.0001)cutYdiam=99999.; //protection for missing diamond information
1782         if( (distXdiam<cutXdiam) && (distYdiam<cutYdiam) && (distZ>nSigmaZdist*errzDist) ){
1783           return kTRUE;
1784         }
1785       }
1786     }
1787   }
1788   return kFALSE;
1789 }
1790
1791 //______________________________________________________________________________
1792 void AliESDEvent::EstimateMultiplicity(Int_t &tracklets, Int_t &trITSTPC, Int_t &trITSSApure, Double_t eta, Bool_t useDCAFlag,Bool_t useV0Flag) const
1793 {
1794   //
1795   // calculates 3 estimators for the multiplicity in the -eta:eta range
1796   // tracklets   : using SPD tracklets only
1797   // trITSTPC    : using TPC/ITS + complementary ITS SA tracks + tracklets from clusters not used by tracks
1798   // trITSSApure : using ITS standalone tracks + tracklets from clusters not used by tracks
1799   // if useDCAFlag is true: account for the ESDtrack flag marking the tracks with large DCA
1800   // if useV0Flag  is true: account for the ESDtrack flag marking conversion and K0's V0s
1801   tracklets = trITSSApure = trITSTPC = 0;
1802   int ntr = fSPDMult ? fSPDMult->GetNumberOfTracklets() : 0;
1803   //
1804   // count tracklets
1805   for (int itr=ntr;itr--;) { 
1806     if (TMath::Abs(fSPDMult->GetEta(itr))>eta) continue;
1807     tracklets++;
1808     if (fSPDMult->FreeClustersTracklet(itr,0)) trITSTPC++;    // not used in ITS/TPC or ITS_SA track
1809     if (fSPDMult->FreeClustersTracklet(itr,1)) trITSSApure++; // not used in ITS_SA_Pure track
1810   }
1811   //
1812   // count real tracks
1813   ntr = GetNumberOfTracks();
1814   for (int itr=ntr;itr--;) {
1815     AliESDtrack *t = GetTrack(itr);
1816     if (TMath::Abs(t->Eta())>eta) continue;
1817     if (!t->IsOn(AliESDtrack::kITSin)) continue;
1818     if (useDCAFlag && t->IsOn(AliESDtrack::kMultSec))  continue;
1819     if (useV0Flag  && t->IsOn(AliESDtrack::kMultInV0)) continue;    
1820     if (t->IsOn(AliESDtrack::kITSpureSA)) trITSSApure++;
1821     else                                  trITSTPC++;
1822   }
1823   //
1824 }
1825
1826 Bool_t AliESDEvent::IsPileupFromSPDInMultBins() const {
1827     Int_t nTracklets=GetMultiplicity()->GetNumberOfTracklets();
1828     if(nTracklets<20) return IsPileupFromSPD(3,0.8);
1829     else if(nTracklets<50) return IsPileupFromSPD(4,0.8);
1830     else return IsPileupFromSPD(5,0.8);
1831 }
1832
1833 void  AliESDEvent::SetTOFHeader(const AliTOFHeader *header)
1834 {
1835   //
1836   // Set the TOF event_time
1837   //
1838
1839   if (fTOFHeader) {
1840     *fTOFHeader=*header;
1841     //fTOFHeader->SetName(fgkESDListName[kTOFHeader]);
1842   }
1843   else {
1844     // for analysis of reconstructed events
1845     // when this information is not avaliable
1846     fTOFHeader = new AliTOFHeader(*header);
1847     //AddObject(fTOFHeader);
1848   }
1849
1850 }
1851
1852 AliCentrality* AliESDEvent::GetCentrality()
1853 {
1854     if (!fCentrality) fCentrality = new AliCentrality();
1855     return  fCentrality;
1856 }
1857
1858 AliEventplane* AliESDEvent::GetEventplane()
1859 {
1860     if (!fEventplane) fEventplane = new AliEventplane();
1861     return  fEventplane;
1862 }
1863
1864 Float_t AliESDEvent::GetVZEROEqMultiplicity(Int_t i) const
1865 {
1866   // Get VZERO Multiplicity for channel i
1867   // Themethod uses the equalization factors
1868   // stored in the ESD-run object in order to
1869   // get equal multiplicities within a VZERO rins (1/8 of VZERO)
1870   if (!fESDVZERO || !fESDRun) return -1;
1871
1872   Int_t ring = i/8;
1873   Float_t factorSum = 0;
1874   for(Int_t j = 8*ring; j < (8*ring+8); ++j) {
1875     factorSum += fESDRun->GetVZEROEqFactors(j);
1876   }
1877   Float_t factor = fESDRun->GetVZEROEqFactors(i)*8./factorSum;
1878
1879   return (fESDVZERO->GetMultiplicity(i)/factor);
1880 }