]> git.uio.no Git - u/mrichter/AliRoot.git/blob - STEER/AliAODHandler.cxx
Make it find libraries independently of location
[u/mrichter/AliRoot.git] / STEER / AliAODHandler.cxx
1
2 /**************************************************************************
3  * Copyright(c) 1998-2007, ALICE Experiment at CERN, All rights reserved. *
4  *                                                                        *
5  * Author: The ALICE Off-line Project.                                    *
6  * Contributors are mentioned in the code where appropriate.              *
7  *                                                                        *
8  * Permission to use, copy, modify and distribute this software and its   *
9  * documentation strictly for non-commercial purposes is hereby granted   *
10  * without fee, provided that the above copyright notice appears in all   *
11  * copies and that both the copyright notice and this permission notice   *
12  * appear in the supporting documentation. The authors make no claims     *
13  * about the suitability of this software for any purpose. It is          *
14  * provided "as is" without express or implied warranty.                  *
15  **************************************************************************/
16
17 /* $Id$ */
18
19 //-------------------------------------------------------------------------
20 //     Implementation of the Virtual Event Handler Interface for AOD
21 //     Author: Andreas Morsch, CERN
22 //-------------------------------------------------------------------------
23
24
25 #include <TTree.h>
26 #include <TFile.h>
27 #include <TString.h>
28 #include <TList.h>
29 #include <TROOT.h>
30
31 #include "AliLog.h"
32 #include "AliAODHandler.h"
33 #include "AliAODEvent.h"
34 #include "AliAODTracklets.h"
35 #include "AliStack.h"
36 #include "AliAODMCParticle.h"
37 #include "AliAODMCHeader.h"
38 #include "AliMCEventHandler.h"
39 #include "AliMCEvent.h"
40 #include "AliGenEventHeader.h"
41 #include "AliGenHijingEventHeader.h"
42 #include "AliGenDPMjetEventHeader.h"
43 #include "AliGenPythiaEventHeader.h"
44 #include "AliGenCocktailEventHeader.h"
45
46
47
48 ClassImp(AliAODHandler)
49
50 //______________________________________________________________________________
51 AliAODHandler::AliAODHandler() :
52     AliVEventHandler(),
53     fIsStandard(kTRUE),
54     fFillAOD(kTRUE),
55     fFillAODRun(kTRUE),
56     fNeedsHeaderReplication(kFALSE),
57     fNeedsTracksBranchReplication(kFALSE),
58     fNeedsVerticesBranchReplication(kFALSE),
59     fNeedsV0sBranchReplication(kFALSE),
60     fNeedsCascadesBranchReplication(kFALSE),
61     fNeedsTrackletsBranchReplication(kFALSE),
62     fNeedsPMDClustersBranchReplication(kFALSE),
63     fNeedsJetsBranchReplication(kFALSE),
64     fNeedsFMDClustersBranchReplication(kFALSE),
65     fNeedsCaloClustersBranchReplication(kFALSE),
66     fNeedsMCParticlesBranchReplication(kFALSE),
67     fNeedsDimuonsBranchReplication(kFALSE),
68     fAODIsReplicated(kFALSE),
69     fAODEvent(NULL),
70     fMCEventH(NULL),
71     fTreeA(NULL),
72     fFileA(NULL),
73     fFileName(""),
74     fExtensions(NULL),
75     fFilters(NULL)
76 {
77   // default constructor
78 }
79
80 //______________________________________________________________________________
81 AliAODHandler::AliAODHandler(const char* name, const char* title):
82     AliVEventHandler(name, title),
83     fIsStandard(kTRUE),
84     fFillAOD(kTRUE),
85     fFillAODRun(kTRUE),
86     fNeedsHeaderReplication(kFALSE),
87     fNeedsTracksBranchReplication(kFALSE),
88     fNeedsVerticesBranchReplication(kFALSE),
89     fNeedsV0sBranchReplication(kFALSE),
90     fNeedsCascadesBranchReplication(kFALSE),
91     fNeedsTrackletsBranchReplication(kFALSE),
92     fNeedsPMDClustersBranchReplication(kFALSE),
93     fNeedsJetsBranchReplication(kFALSE),
94     fNeedsFMDClustersBranchReplication(kFALSE),
95     fNeedsCaloClustersBranchReplication(kFALSE),
96     fNeedsMCParticlesBranchReplication(kFALSE),
97     fNeedsDimuonsBranchReplication(kFALSE),
98     fAODIsReplicated(kFALSE),
99     fAODEvent(NULL),
100     fMCEventH(NULL),
101     fTreeA(NULL),
102     fFileA(NULL),
103     fFileName(""),
104     fExtensions(NULL),
105     fFilters(NULL)
106 {
107 // Normal constructor.
108 }
109
110 //______________________________________________________________________________
111 AliAODHandler::~AliAODHandler() 
112 {
113  // Destructor.
114   if (fAODEvent) delete fAODEvent;
115   if(fFileA){
116     // is already handled in TerminateIO
117     fFileA->Close();
118     delete fFileA;
119     fTreeA = 0;
120   }
121   if (fTreeA) delete fTreeA;
122   if (fExtensions) {fExtensions->Delete(); delete fExtensions;}
123   if (fFilters)    {fFilters->Delete();    delete fFilters;}
124 }
125
126 //______________________________________________________________________________
127 Bool_t AliAODHandler::Init(Option_t* opt)
128 {
129   // Initialize IO
130   //
131   // Create the AODevent object
132     
133   Bool_t createStdAOD = fIsStandard || fFillAOD;
134   if(!fAODEvent && createStdAOD){
135     fAODEvent = new AliAODEvent();
136     if (fIsStandard) 
137       fAODEvent->CreateStdContent();
138   }
139   //
140   // File opening according to execution mode
141   TString option(opt);
142   option.ToLower();
143   if (createStdAOD) {
144     TDirectory *owd = gDirectory;
145     if (option.Contains("proof")) {
146       // proof
147       // Merging via files. Need to access analysis manager via interpreter.
148       gROOT->ProcessLine(Form("AliAnalysisDataContainer *c_common_out = AliAnalysisManager::GetAnalysisManager()->GetCommonOutputContainer();"));
149       gROOT->ProcessLine(Form("AliAnalysisManager::GetAnalysisManager()->OpenProofFile(c_common_out, \"RECREATE\");"));
150       fFileA = gFile;
151     } else {
152       // local and grid
153       fFileA = new TFile(fFileName.Data(), "RECREATE");
154     }
155     CreateTree(1);
156     owd->cd();
157   }  
158   if (fExtensions) {
159      TIter next(fExtensions);
160      AliAODExtension *ext;
161      while ((ext=(AliAODExtension*)next())) ext->Init(option);
162   }   
163   if (fFilters) {
164      TIter nextf(fFilters);
165      AliAODExtension *filteredAOD;
166      while ((filteredAOD=(AliAODExtension*)nextf())) {
167         filteredAOD->SetEvent(fAODEvent);
168         filteredAOD->Init(option);
169      }   
170   }   
171   return kTRUE;
172 }
173
174 //______________________________________________________________________________
175 void AliAODHandler::StoreMCParticles(){
176
177   // 
178   // Remap the labels from ESD stack and store
179   // the AODMCParticles, makes only sense if we have
180   // the mcparticles branch
181   // has to be done here since we cannot know in advance 
182   // which particles are needed (e.g. by the tracks etc.)
183   //
184   // Particles have been selected by AliMCEventhanlder->SelectParticle()
185   // To use the MCEventhandler here we need to set it from the outside
186   // can vanish when Handler go to the ANALYSISalice library
187   //
188   // The Branch booking for mcParticles and mcHeader has to happen 
189   // in an external task for now since the AODHandler does not have access
190   // the AnalysisManager. For the same reason the pointer t o the MCEventH
191   // has to passed to the AOD Handler by this task 
192   // (doing this in the steering macro would not work on PROOF)
193
194   if (!fAODEvent) return;
195   TClonesArray *mcarray = (TClonesArray*)fAODEvent->FindListObject(AliAODMCParticle::StdBranchName()); 
196   if(!mcarray)return;
197
198   AliAODMCHeader *mcHeader = (AliAODMCHeader*)fAODEvent->FindListObject(AliAODMCHeader::StdBranchName()); 
199   if(!mcHeader)return;
200
201   // Get the MC Infos.. Handler needs to be set before 
202   // while adding the branch
203   // This needs to be done, not to depend on the AnalysisManager
204
205   if(!fMCEventH)return;
206   if(!fMCEventH->MCEvent())return;
207   AliStack *pStack = fMCEventH->MCEvent()->Stack();
208   if(!pStack)return;
209
210   fMCEventH->CreateLabelMap();
211
212   //
213   // Get the Event Header 
214   // 
215
216   AliHeader* header = fMCEventH->MCEvent()->Header();
217    // get the MC vertex
218   AliGenEventHeader* genHeader = 0;
219   if (header) genHeader = header->GenEventHeader();
220   if (genHeader) {
221       TArrayF vtxMC(3);
222       genHeader->PrimaryVertex(vtxMC);
223       mcHeader->SetVertex(vtxMC[0],vtxMC[1],vtxMC[2]);
224
225       // we search the MCEventHeaders first 
226       // Two cases, cocktail or not...
227       AliGenCocktailEventHeader* genCocktailHeader = dynamic_cast<AliGenCocktailEventHeader*>(genHeader);
228       if(genCocktailHeader){
229           // we have a coktail header
230           mcHeader->AddGeneratorName(genHeader->GetName());
231           // Loop from the back so that the first one sets the process type
232           TList* headerList = genCocktailHeader->GetHeaders();
233           for(int i = headerList->GetEntries()-1;i>=0;--i){
234               AliGenEventHeader *headerEntry = dynamic_cast<AliGenEventHeader*>(headerList->At(i));
235               SetMCHeaderInfo(mcHeader,headerEntry);
236           }
237       }
238       else{
239           // No Cocktail just take the first one
240           SetMCHeaderInfo(mcHeader,genHeader);
241       }
242   }
243   
244
245
246
247
248   // Store the AliAODParticlesMC
249   AliMCEvent* mcEvent = fMCEventH->MCEvent();
250   
251   Int_t np    = mcEvent->GetNumberOfTracks();
252   Int_t nprim = mcEvent->GetNumberOfPrimaries();
253
254
255   Int_t j = 0;
256   TClonesArray& l = *mcarray;
257
258   for(int i = 0; i < np; ++i){
259       if(fMCEventH->IsParticleSelected(i)){
260           Int_t flag = 0;
261           AliMCParticle* mcpart =  (AliMCParticle*) mcEvent->GetTrack(i);
262           if(i<nprim)flag |= AliAODMCParticle::kPrimary;
263           
264           if(mcEvent->IsPhysicalPrimary(i))flag |= AliAODMCParticle::kPhysicalPrim;
265           
266           if(fMCEventH->GetNewLabel(i)!=j){
267               AliError(Form("MISMATCH New label %d j: %d",fMCEventH->GetNewLabel(i),j));
268           }
269
270           AliAODMCParticle mcpart_tmp(mcpart,i,flag);
271           
272           mcpart_tmp.SetStatus(mcpart->Particle()->GetStatusCode());
273           // 
274           Int_t d0 =  mcpart_tmp.GetDaughter(0);
275           Int_t d1 =  mcpart_tmp.GetDaughter(1);
276           Int_t m =   mcpart_tmp.GetMother();
277           
278           // other than for the track labels, negative values mean
279           // no daughter/mother so preserve it
280           
281           if(d0<0 && d1<0){
282               // no first daughter -> no second daughter
283               // nothing to be done
284               // second condition not needed just for sanity check at the end
285               mcpart_tmp.SetDaughter(0,d0);
286               mcpart_tmp.SetDaughter(1,d1);
287           } else if(d1 < 0 && d0 >= 0) {
288               // Only one daughter
289               // second condition not needed just for sanity check at the end
290               if(fMCEventH->IsParticleSelected(d0)){
291                   mcpart_tmp.SetDaughter(0,fMCEventH->GetNewLabel(d0));
292               } else {
293                   mcpart_tmp.SetDaughter(0,-1);
294               }
295               mcpart_tmp.SetDaughter(1,d1);
296           }
297           else if (d0 > 0 && d1 > 0 ){
298               // we have two or more daughters loop on the stack to see if they are
299               // selected
300               Int_t d0_tmp = -1;
301               Int_t d1_tmp = -1;
302               for(int id = d0; id<=d1;++id){
303                   if(fMCEventH->IsParticleSelected(id)){
304                       if(d0_tmp==-1){
305                           // first time
306                           d0_tmp = fMCEventH->GetNewLabel(id);
307                           d1_tmp = d0_tmp; // this is to have the same schema as on the stack i.e. with one daugther d0 and d1 are the same 
308                       }
309                       else d1_tmp = fMCEventH->GetNewLabel(id);
310                   }
311               }
312               mcpart_tmp.SetDaughter(0,d0_tmp);
313               mcpart_tmp.SetDaughter(1,d1_tmp);
314           } else {
315               AliError(Form("Unxpected indices %d %d",d0,d1));
316           }
317           
318           if(m<0){
319               mcpart_tmp.SetMother(m);
320           } else {
321               if(fMCEventH->IsParticleSelected(m))mcpart_tmp.SetMother(fMCEventH->GetNewLabel(m));
322               else AliError(Form("PROBLEM Mother not selected %d \n", m));
323           }
324
325           new (l[j++]) AliAODMCParticle(mcpart_tmp);
326           
327       }
328   }
329   AliInfo(Form("AliAODHandler::StoreMCParticles: Selected %d (Primaries %d / total %d) after validation \n",
330                j,nprim,np));
331   
332   // Set the labels in the AOD output...
333   // Remapping
334
335   // AODTracks
336   TClonesArray* tracks = fAODEvent->GetTracks();
337   if(tracks){
338     for(int it = 0; it < fAODEvent->GetNTracks();++it){
339       AliAODTrack *track = fAODEvent->GetTrack(it);
340       
341       Int_t sign = 1;
342       Int_t label = track->GetLabel();
343       if(label<0){ // preserve the sign for later usage
344         label *= -1;
345         sign  = -1;
346       }
347
348       if (label >= AliMCEvent::BgLabelOffset()) label =  mcEvent->BgLabelToIndex(label);
349       if(label > np || track->GetLabel() == 0){
350         AliWarning(Form("Wrong ESD track label %5d (%5d)",track->GetLabel(), label));
351       }
352       if(fMCEventH->GetNewLabel(label) == 0){
353         AliWarning(Form("New label not found for %5d (%5d)",track->GetLabel(), label));
354       }
355       track->SetLabel(sign*fMCEventH->GetNewLabel(label));
356     }
357   }
358   
359   // AOD calo cluster
360   TClonesArray *clusters = fAODEvent->GetCaloClusters();
361   if(clusters){
362     for (Int_t iClust = 0;iClust < fAODEvent->GetNumberOfCaloClusters(); ++iClust) {
363       AliAODCaloCluster * cluster = fAODEvent->GetCaloCluster(iClust);
364       UInt_t nLabel    = cluster->GetNLabels();
365       // Ugly but do not want to fragment memory by creating 
366       // new Int_t (nLabel)
367       Int_t* labels    = const_cast<Int_t*>(cluster->GetLabels());
368       if (labels){
369         for(UInt_t i = 0;i < nLabel;++i){
370           labels[i] = fMCEventH->GetNewLabel(cluster->GetLabelAt(i));
371         }
372       }
373       //      cluster->SetLabels(labels,nLabel);
374     }// iClust
375   }// clusters
376
377   // AOD tracklets
378   AliAODTracklets *tracklets = fAODEvent->GetTracklets();
379   if(tracklets){
380     for(int it = 0;it < tracklets->GetNumberOfTracklets();++it){
381       int label0 = tracklets->GetLabel(it,0);
382       int label1 = tracklets->GetLabel(it,1);
383       if(label0>=0)label0 = fMCEventH->GetNewLabel(label0);      
384       if(label1>=0)label1 = fMCEventH->GetNewLabel(label1);
385       tracklets->SetLabel(it,0,label0);
386       tracklets->SetLabel(it,1,label1);
387     }
388   }
389
390 }
391
392 //______________________________________________________________________________
393 Bool_t AliAODHandler::FinishEvent()
394 {
395   // Fill data structures
396   if(fFillAOD && fFillAODRun && fAODEvent){
397       fAODEvent->MakeEntriesReferencable();
398       fTreeA->BranchRef();
399       FillTree();
400   }
401
402   if (fFillAOD && fFillAODRun) {      
403       if (fExtensions) {
404           TIter next(fExtensions);
405           AliAODExtension *ext;
406           while ((ext=(AliAODExtension*)next())) ext->FinishEvent();
407       }
408       if (fFilters) {   
409           TIter nextf(fFilters);
410           AliAODExtension *ext;
411           while ((ext=(AliAODExtension*)nextf())) {
412               ext->FinishEvent();
413           }  
414       }       
415   }  
416   if (fIsStandard) fAODEvent->ResetStd();
417   if (fAODEvent) {
418     TClonesArray *mcarray = (TClonesArray*)fAODEvent->FindListObject(AliAODMCParticle::StdBranchName()); 
419     if(mcarray) mcarray->Delete();
420     
421     AliAODMCHeader *mcHeader = (AliAODMCHeader*)fAODEvent->FindListObject(AliAODMCHeader::StdBranchName()); 
422     if(mcHeader) mcHeader->Reset();
423   }
424   // Reset AOD replication flag
425   fAODIsReplicated = kFALSE;
426   return kTRUE;
427 }
428
429 //______________________________________________________________________________
430 Bool_t AliAODHandler::Terminate()
431 {
432   // Terminate 
433   AddAODtoTreeUserInfo();
434   if (fExtensions) {
435     TIter next(fExtensions);
436     AliAODExtension *ext;
437     while ((ext=(AliAODExtension*)next())) ext->GetTree()->GetUserInfo()->Add(ext->GetAOD());
438   }  
439   return kTRUE;
440 }
441
442 //______________________________________________________________________________
443 Bool_t AliAODHandler::TerminateIO()
444 {
445   // Terminate IO
446   if (fFileA) {
447     fFileA->Write();
448     fFileA->Close();
449     delete fFileA;
450     fFileA = 0;
451     // When closing the file, the tree is also deleted.
452     fTreeA = 0;
453   }
454   if (fExtensions) {
455     TIter next(fExtensions);
456     AliAODExtension *ext;
457     while ((ext=(AliAODExtension*)next())) ext->TerminateIO();
458   }  
459   if (fFilters) {
460     TIter nextf(fFilters);
461     AliAODExtension *ext;
462     while ((ext=(AliAODExtension*)nextf())) ext->TerminateIO();
463   }  
464   return kTRUE;
465 }
466
467 //______________________________________________________________________________
468 void AliAODHandler::CreateTree(Int_t flag)
469 {
470     // Creates the AOD Tree
471     fTreeA = new TTree("aodTree", "AliAOD tree");
472     fTreeA->Branch(fAODEvent->GetList());
473     if (flag == 0) fTreeA->SetDirectory(0);
474 }
475
476 //______________________________________________________________________________
477 void AliAODHandler::FillTree()
478 {
479  
480     // Fill the AOD Tree
481     fTreeA->Fill();
482 }
483
484 //______________________________________________________________________________
485 void AliAODHandler::AddAODtoTreeUserInfo()
486 {
487   // Add aod event to tree user info
488   if (fTreeA) fTreeA->GetUserInfo()->Add(fAODEvent);
489   // Now the tree owns our fAODEvent...
490   fAODEvent = 0;
491 }
492
493 //______________________________________________________________________________
494 void AliAODHandler::AddBranch(const char* cname, void* addobj, const char* filename)
495 {
496     // Add a new branch to the aod. Added optional filename parameter if the
497     // branch should be written to a separate file.
498     if (strlen(filename)) {
499        AliAODExtension *ext = AddExtension(filename);
500        ext->AddBranch(cname, addobj);
501        return;
502     }
503     TDirectory *owd = gDirectory;
504     if (fFileA) {
505       fFileA->cd();
506     }
507     char** apointer = (char**) addobj;
508     TObject* obj = (TObject*) *apointer;
509
510     fAODEvent->AddObject(obj);
511  
512     const Int_t kSplitlevel = 99; // default value in TTree::Branch()
513     const Int_t kBufsize = 32000; // default value in TTree::Branch()
514
515     if (!fTreeA->FindBranch(obj->GetName())) {
516       // Do the same as if we book via 
517       // TTree::Branch(TCollection*)
518       
519       fTreeA->Bronch(obj->GetName(), cname, fAODEvent->GetList()->GetObjectRef(obj),
520                      kBufsize, kSplitlevel - 1);
521     }
522     owd->cd();
523 }
524
525 //______________________________________________________________________________
526 AliAODExtension *AliAODHandler::AddExtension(const char *filename, const char *title)
527 {
528 // Add an AOD extension with some branches in a different file.
529    TString fname(filename);
530    if (!fname.EndsWith(".root")) fname += ".root";
531    if (!fExtensions) {
532       fExtensions = new TObjArray();
533       fExtensions->SetOwner();
534    }   
535    AliAODExtension *ext = (AliAODExtension*)fExtensions->FindObject(fname);
536    if (!ext) {
537       ext = new AliAODExtension(fname, title);
538       fExtensions->Add(ext);
539    }   
540    return ext;
541 }
542
543 //______________________________________________________________________________
544 AliAODExtension *AliAODHandler::GetExtension(const char *filename) const
545 {
546 // Getter for AOD extensions via file name.
547    if (!fExtensions) return NULL;
548    return (AliAODExtension*)fExtensions->FindObject(filename);
549 }   
550
551 //______________________________________________________________________________
552 AliAODExtension *AliAODHandler::AddFilteredAOD(const char *filename, const char *filtername)
553 {
554 // Add an AOD extension that can write only AOD events that pass a user filter.
555    if (!fFilters) {
556       fFilters = new TObjArray();
557       fFilters->SetOwner();
558    } 
559    AliAODExtension *filter = (AliAODExtension*)fFilters->FindObject(filename);
560    if (!filter) {
561       filter = new AliAODExtension(filename, filtername, kTRUE);
562       fFilters->Add(filter);
563    }
564    return filter;
565 }      
566
567 //______________________________________________________________________________
568 AliAODExtension *AliAODHandler::GetFilteredAOD(const char *filename) const
569 {
570 // Getter for AOD filters via file name.
571    if (!fFilters) return NULL;
572    return (AliAODExtension*)fFilters->FindObject(filename);
573 }   
574    
575 //______________________________________________________________________________
576 void AliAODHandler::SetOutputFileName(const char* fname)
577 {
578 // Set file name.
579    fFileName = fname;
580 }
581
582 //______________________________________________________________________________
583 const char *AliAODHandler::GetOutputFileName()
584 {
585 // Get file name.
586    return fFileName.Data();
587 }
588
589 //______________________________________________________________________________
590 const char *AliAODHandler::GetExtraOutputs() const
591 {
592 // Get extra outputs as a string separated by commas.
593    static TString eoutputs;
594    eoutputs = "";
595    TObject *obj;
596    if (fExtensions) {
597       TIter next1(fExtensions);
598       while ((obj=next1())) {
599          if (!eoutputs.IsNull()) eoutputs += ",";
600          eoutputs += obj->GetName();
601       }
602    }
603    if (fFilters) {
604       TIter next2(fFilters);
605       while ((obj=next2())) {
606          if (!eoutputs.IsNull()) eoutputs += ",";
607          eoutputs += obj->GetName();
608       }
609    }
610    return eoutputs.Data();
611 }
612
613 //______________________________________________________________________________
614 void  AliAODHandler::SetMCHeaderInfo(AliAODMCHeader *mcHeader,AliGenEventHeader *genHeader){
615
616
617   // Utility function to cover different cases for the AliGenEventHeader
618   // Needed since different ProcessType and ImpactParamter are not 
619   // in the base class...
620   // We don't encode process types for event cocktails yet
621   // could be done e.g. by adding offsets depnding on the generator
622
623   mcHeader->AddGeneratorName(genHeader->GetName());
624   if(!genHeader)return;
625   AliGenPythiaEventHeader *pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
626     if (pythiaGenHeader) {
627       mcHeader->SetEventType(pythiaGenHeader->ProcessType());
628       mcHeader->SetPtHard(pythiaGenHeader->GetPtHard());
629       return;
630     }
631     
632     AliGenDPMjetEventHeader* dpmJetGenHeader = dynamic_cast<AliGenDPMjetEventHeader*>(genHeader);
633
634   if (dpmJetGenHeader){
635     mcHeader->SetEventType(dpmJetGenHeader->ProcessType());
636     return;
637   } 
638
639   AliGenHijingEventHeader* hijingGenHeader = dynamic_cast<AliGenHijingEventHeader*>(genHeader);
640   if(hijingGenHeader){
641     mcHeader->SetImpactParameter(hijingGenHeader->ImpactParameter());
642     return;
643   }
644   
645   AliWarning(Form("MC Eventheader not known: %s",genHeader->GetName()));
646
647 }
648
649 ClassImp(AliAODExtension)
650
651 //-------------------------------------------------------------------------
652 //     Support class for AOD extensions. This is created by the user analysis
653 //     that requires a separate file for some AOD branches. The name of the 
654 //     AliAODExtension object is the file name where the AOD branches will be
655 //     stored.
656 //-------------------------------------------------------------------------
657
658 //______________________________________________________________________________
659 AliAODExtension::AliAODExtension(const char* name, const char* title, Bool_t isfilter)
660                 :TNamed(name,title), 
661                  fAODEvent(0), 
662                  fTreeE(0), 
663                  fFileE(0), 
664                  fNtotal(0), 
665                  fNpassed(0),
666                  fSelected(kFALSE)
667 {
668 // Constructor.
669   if (isfilter) {
670     TObject::SetBit(kFilteredAOD);
671     printf("####### Added AOD filter %s\n", name);
672   } else printf("####### Added AOD extension %s\n", name);
673 }   
674
675 //______________________________________________________________________________
676 AliAODExtension::~AliAODExtension()
677 {
678 // Destructor.
679   if(fFileE){
680     // is already handled in TerminateIO
681     fFileE->Close();
682     delete fFileE;
683     fTreeE = 0;
684     fAODEvent = 0;
685   }
686   if (fTreeE) delete fTreeE;
687 }
688
689 //______________________________________________________________________________
690 void AliAODExtension::AddBranch(const char* cname, void* addobj)
691 {
692     // Add a new branch to the aod 
693     if (IsFilteredAOD()) {
694        Error("AddBranch", "Not allowed to add branched to filtered AOD's.");
695        return;
696     }   
697     if (!fAODEvent) {
698        char type[20];
699        gROOT->ProcessLine(Form("TString s_tmp; AliAnalysisManager::GetAnalysisManager()->GetAnalysisTypeString(s_tmp); sprintf((char*)%p, \"%%s\", s_tmp.Data());", type));
700        Init(type);
701     }
702     TDirectory *owd = gDirectory;
703     if (fFileE) {
704       fFileE->cd();
705     }
706     char** apointer = (char**) addobj;
707     TObject* obj = (TObject*) *apointer;
708
709     fAODEvent->AddObject(obj);
710  
711     const Int_t kSplitlevel = 99; // default value in TTree::Branch()
712     const Int_t kBufsize = 32000; // default value in TTree::Branch()
713
714     if (!fTreeE->FindBranch(obj->GetName())) {
715       // Do the same as if we book via 
716       // TTree::Branch(TCollection*)
717       fTreeE->Bronch(obj->GetName(), cname, fAODEvent->GetList()->GetObjectRef(obj),
718                      kBufsize, kSplitlevel - 1);
719     }
720     owd->cd();
721 }
722
723 //______________________________________________________________________________
724 Bool_t AliAODExtension::FinishEvent()
725 {
726 // Fill current event.
727   fNtotal++;
728   if (!IsFilteredAOD()) {
729     fAODEvent->MakeEntriesReferencable();
730     fTreeE->Fill();
731     return kTRUE;
732   }  
733   // Filtered AOD. Fill only if event is selected.
734   if (!fSelected) return kTRUE;
735   fNpassed++;
736   fTreeE->Fill();
737   fSelected = kFALSE; // so that next event will not be selected unless demanded
738   return kTRUE;
739 }  
740
741 //______________________________________________________________________________
742 Bool_t AliAODExtension::Init(Option_t *option)
743 {
744 // Initialize IO.
745   if(!fAODEvent) fAODEvent = new AliAODEvent();
746   TDirectory *owd = gDirectory;
747   TString opt(option);
748   opt.ToLower();
749   if (opt.Contains("proof")) {
750     // proof
751     // Merging via files. Need to access analysis manager via interpreter.
752     gROOT->ProcessLine(Form("AliAnalysisManager::GetAnalysisManager()->OpenProofFile(\"%s\", \"RECREATE\");", fName.Data()));
753     fFileE = gFile;
754   } else {
755     fFileE = new TFile(GetName(), "RECREATE");
756   }  
757   fTreeE = new TTree("aodTree", "AliAOD tree");
758   fTreeE->Branch(fAODEvent->GetList());
759   fTreeE->BranchRef();
760   owd->cd();
761   return kTRUE;
762 }
763
764 //______________________________________________________________________________
765 void AliAODExtension::SetEvent(AliAODEvent *event)
766 {
767 // Connects to an external event
768    if (!IsFilteredAOD()) {
769       Error("SetEvent", "Not allowed to set external event for filtered AOD's");   
770       return;
771    }
772    // Use the copy constructor or assignment operator to synchronize with external event.
773 //   AliAODEvent &other = *event;
774 //   if (!fAODEvent)     fAODEvent = new AliAODEvent(other);
775 //   else if (fSelected) *fAODEvent = other;
776    fAODEvent = event;
777 }
778    
779 //______________________________________________________________________________
780 Bool_t AliAODExtension::TerminateIO()
781 {
782   // Terminate IO
783   if (TObject::TestBit(kFilteredAOD))
784     printf("AOD Filter %s: events processed: %d   passed: %d\n", GetName(), fNtotal, fNpassed);
785   else
786     printf("AOD extension %s: events processed: %d\n", GetName(), fNtotal);
787   if (fFileE) {
788     fFileE->Write();
789     fFileE->Close();
790     delete fFileE;
791     fFileE = 0;
792     fTreeE = 0;
793     fAODEvent = 0;
794   }
795   return kTRUE;
796 }