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