]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - STEER/AliMCEventHandler.cxx
Optimisation
[u/mrichter/AliRoot.git] / STEER / AliMCEventHandler.cxx
index 6fb6528077c83f96ba1fa953d3d5a1ad5882001c..866b5c6750dc0ca06cbd08899fd53ec3bace754f 100644 (file)
@@ -1,5 +1,5 @@
-/**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+/************************************************************************** 
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved  *
  *                                                                        *
  * Author: The ALICE Off-line Project.                                    *
  * Contributors are mentioned in the code where appropriate.              *
@@ -26,6 +26,9 @@
 
 
 #include "AliMCEventHandler.h"
+#include "AliMCEvent.h"
+#include "AliMCParticle.h"
+#include "AliPDG.h"
 #include "AliTrackReference.h"
 #include "AliHeader.h"
 #include "AliStack.h"
 
 #include <TTree.h>
 #include <TFile.h>
+#include <TList.h>
 #include <TParticle.h>
+#include <TString.h>
 #include <TClonesArray.h>
 #include <TDirectoryFile.h>
-#include <TArrow.h>
-#include <TMarker.h>
-#include <TH2F.h>
-
 
 ClassImp(AliMCEventHandler)
 
 AliMCEventHandler::AliMCEventHandler() :
     AliVEventHandler(),
+    fMCEvent(new AliMCEvent()),
     fFileE(0),
     fFileK(0),
     fFileTR(0),
-    fTmpFileTR(0),
     fTreeE(0),
     fTreeK(0),
     fTreeTR(0),
-    fTmpTreeTR(0),
-    fStack(0),
-    fHeader(0),
-    fTrackReferences(0),
+    fDirK(0),
+    fDirTR(0),
+    fParticleSelected(0),
+    fLabelMap(0),
     fNEvent(-1),
     fEvent(-1),
-    fNprimaries(-1),
-    fNparticles(-1),
-    fPathName("./"),
+    fPathName(new TString("./")),
     fExtension(""),
     fFileNumber(0),
-    fEventsPerFile(0)
+    fEventsPerFile(0),
+    fReadTR(kTRUE),
+    fInitOk(kFALSE),
+    fSubsidiaryHandlers(0),
+    fEventsInContainer(0),
+    fPreReadMode(kNoPreRead)
 {
-    // Default constructor
+  //
+  // Default constructor
+  //
+  // Be sure to add all particles to the PDG database
+  AliPDG::AddParticlesToPdgDataBase();
 }
 
 AliMCEventHandler::AliMCEventHandler(const char* name, const char* title) :
     AliVEventHandler(name, title),
+    fMCEvent(new AliMCEvent()),
     fFileE(0),
     fFileK(0),
     fFileTR(0),
-    fTmpFileTR(0),
     fTreeE(0),
     fTreeK(0),
     fTreeTR(0),
-    fTmpTreeTR(0),
-    fStack(0),
-    fHeader(new AliHeader()),
-    fTrackReferences(new TClonesArray("AliTrackReference", 200)),
+    fDirK(0),
+    fDirTR(0),
+    fParticleSelected(0),
+    fLabelMap(0),
     fNEvent(-1),
     fEvent(-1),
-    fNprimaries(-1),
-    fNparticles(-1),
-    fPathName("./"),
+    fPathName(new TString("./")),
     fExtension(""),
     fFileNumber(0),
-    fEventsPerFile(0)
+    fEventsPerFile(0),
+    fReadTR(kTRUE),
+    fInitOk(kFALSE),
+    fSubsidiaryHandlers(0),
+    fEventsInContainer(0),
+    fPreReadMode(kNoPreRead)
 {
-    // Constructor
+  //
+  // Constructor
+  //
+  // Be sure to add all particles to the PDG database
+  AliPDG::AddParticlesToPdgDataBase();
 }
 AliMCEventHandler::~AliMCEventHandler()
 { 
     // Destructor
+  delete fPathName;
+    delete fMCEvent;
     delete fFileE;
     delete fFileK;
     delete fFileTR;
 }
 
-Bool_t AliMCEventHandler::InitIO(Option_t* /*opt*/) 
+Bool_t AliMCEventHandler::Init(Option_t* opt)
 { 
     // Initialize input
     //
-    fFileE = TFile::Open(Form("%sgalice.root", fPathName));
-    if (!fFileE) AliFatal(Form("AliMCEventHandler:galice.root not found in directory %s ! \n", fPathName));
-
+    if (!(strcmp(opt, "proof")) || !(strcmp(opt, "local"))) return kTRUE;
+    //
+    fFileE = TFile::Open(Form("%sgalice.root", fPathName->Data()));
+    if (!fFileE) {
+       AliError(Form("AliMCEventHandler:galice.root not found in directory %s ! \n", fPathName->Data()));
+       fInitOk = kFALSE;
+       return kFALSE;
+    }
+    
+    //
+    // Tree E
     fFileE->GetObject("TE", fTreeE);
-    fTreeE->SetBranchAddress("Header", &fHeader);
+    // Connect Tree E to the MCEvent
+    fMCEvent->ConnectTreeE(fTreeE);
     fNEvent = fTreeE->GetEntries();
     //
     // Tree K
-    fFileK = TFile::Open(Form("%sKinematics%s.root", fPathName, fExtension));
-    if (!fFileK) AliFatal(Form("AliMCEventHandler:Kinematics.root not found in directory %s ! \n", fPathName));
+    fFileK = TFile::Open(Form("%sKinematics%s.root", fPathName->Data(), fExtension));
+    if (!fFileK) {
+       AliError(Form("AliMCEventHandler:Kinematics.root not found in directory %s ! \n", fPathName->Data()));
+       fInitOk = kFALSE;
+       return kTRUE;
+    }
+    
     fEventsPerFile = fFileK->GetNkeys() - fFileK->GetNProcessIDs();
     //
     // Tree TR
-    fFileTR = TFile::Open(Form("%sTrackRefs%s.root", fPathName, fExtension));
-    if (!fFileTR) AliWarning(Form("AliMCEventHandler:TrackRefs.root not found in directory %s ! \n", fPathName));
+    if (fReadTR) {
+       fFileTR = TFile::Open(Form("%sTrackRefs%s.root", fPathName->Data(), fExtension));
+       if (!fFileTR) {
+           AliError(Form("AliMCEventHandler:TrackRefs.root not found in directory %s ! \n", fPathName->Data()));
+           fInitOk = kFALSE;
+           return kTRUE;
+       }
+    }
     //
     // Reset the event number
     fEvent      = -1;
     fFileNumber =  0;
-    
-    AliInfo(Form("AliMCEventHandler:Number of events in this directory %5d \n", fNEvent));
+    AliInfo(Form("Number of events in this directory %5d \n", fNEvent));
+    fInitOk = kTRUE;
+
+
+    if (fSubsidiaryHandlers) {
+       TIter next(fSubsidiaryHandlers);
+       AliMCEventHandler *handler;
+       while((handler = (AliMCEventHandler*)next())) {
+           handler->Init(opt);
+           handler->SetNumberOfEventsInContainer(fNEvent);
+       }
+    }
+
     return kTRUE;
-    
 }
 
 Bool_t AliMCEventHandler::GetEvent(Int_t iev)
@@ -134,60 +181,60 @@ Bool_t AliMCEventHandler::GetEvent(Int_t iev)
     // Load the event number iev
     //
     // Calculate the file number
-    Int_t inew  = iev/fEventsPerFile;
-    if (inew != fFileNumber) {
-       fFileNumber = inew;
-       if (!OpenFile(fFileNumber)){
-           return kFALSE;
-       }
-    }
-    // Folder name
-    char folder[20];
-    sprintf(folder, "Event%d", iev);
-
-    // TreeE
-    fTreeE->GetEntry(iev);
-    fStack = fHeader->Stack();
-    // Tree K
-    TDirectoryFile* dirK  = 0;
-    fFileK->GetObject(folder, dirK);
-    if (!dirK) {
-       AliWarning(Form("AliMCEventHandler: Event #%5d not found\n", iev));
-       return kFALSE;
-    }
-    dirK->GetObject("TreeK", fTreeK);
-    fStack->ConnectTree(fTreeK);
-    fStack->GetEvent();
+  if (!fInitOk) return kFALSE;
     
-    //Tree TR 
-
-    if (fFileTR) {
-       TDirectoryFile* dirTR = 0;
-       fFileTR->GetObject(folder, dirTR);
-       dirTR->GetObject("TreeTR", fTreeTR);
-       if (fTreeTR->GetBranch("AliRun")) {
-           // This is an old format with one branch per detector not in synch with TreeK
-           ReorderAndExpandTreeTR();
-           delete dirTR;
-       } else {
-           // New format 
-           fTreeTR->SetBranchAddress("TrackReferences", &fTrackReferences);
-       }
+  Int_t inew  = iev / fEventsPerFile;
+  if (inew != fFileNumber) {
+    fFileNumber = inew;
+    if (!OpenFile(fFileNumber)){
+      return kFALSE;
     }
-
-    //
-    fNparticles = fStack->GetNtrack();
-    fNprimaries = fStack->GetNprimary();
-    AliInfo(Form("AliMCEventHandler: Event#: %5d Number of particles: %5d (all) %5d (primaries)\n", 
-                fEvent, fNparticles, fNprimaries));
+  }
+  // Folder name
+  char folder[20];
+  snprintf(folder, 20, "Event%d", iev);
+  // TreeE
+  fTreeE->GetEntry(iev);
+  // Tree K
+  fFileK->GetObject(folder, fDirK);
+  if (!fDirK) {
+    AliWarning(Form("AliMCEventHandler: Event #%5d - Cannot get kinematics\n", iev));
+    return kFALSE;
+  }
     
-    return kTRUE;
+  fDirK ->GetObject("TreeK", fTreeK);
+  if (!fTreeK) {
+    AliError(Form("AliMCEventHandler: Event #%5d - Cannot get TreeK\n",iev));
+    return kFALSE;
+  }  
+  // Connect TreeK to MCEvent
+  fMCEvent->ConnectTreeK(fTreeK);
+  //Tree TR 
+  if (fFileTR) {
+    // Check which format has been read
+    fFileTR->GetObject(folder, fDirTR);
+    if (!fDirTR) {
+      AliError(Form("AliMCEventHandler: Event #%5d - Cannot get track references\n",iev));
+      return kFALSE;
+    }  
+     
+    fDirTR->GetObject("TreeTR", fTreeTR);
+    //
+    if (!fTreeTR) {
+      AliError(Form("AliMCEventHandler: Event #%5d - Cannot get TreeTR\n",iev));
+      return kFALSE;
+    }  
+    // Connect TR to MCEvent
+    fMCEvent->ConnectTreeTR(fTreeTR);
+  }
+
+  //
+  return kTRUE;
 }
 
 Bool_t AliMCEventHandler::OpenFile(Int_t i)
 {
     // Open file i
-    Bool_t ok = kTRUE;
     if (i > 0) {
        fExtension = Form("%d", i);
     } else {
@@ -196,129 +243,263 @@ Bool_t AliMCEventHandler::OpenFile(Int_t i)
     
     
     delete fFileK;
-    fFileK = TFile::Open(Form("%sKinematics%s.root", fPathName, fExtension));
+    fFileK = TFile::Open(Form("%sKinematics%s.root", fPathName->Data(), fExtension));
     if (!fFileK) {
-       AliFatal(Form("AliMCEventHandler:Kinematics%s.root not found in directory %s ! \n", fExtension, fPathName));
-       ok = kFALSE;
+       AliError(Form("AliMCEventHandler:Kinematics%s.root not found in directory %s ! \n", fExtension, fPathName->Data()));
+       fInitOk = kFALSE;
+       return kFALSE;
     }
     
-    delete fFileTR;
-    fFileTR = TFile::Open(Form("%sTrackRefs%s.root", fPathName, fExtension));
-    if (!fFileTR) {
-       AliWarning(Form("AliMCEventHandler:TrackRefs%s.root not found in directory %s ! \n", fExtension, fPathName));
-       ok = kFALSE;
+    if (fReadTR) {
+       delete fFileTR;
+       fFileTR = TFile::Open(Form("%sTrackRefs%s.root", fPathName->Data(), fExtension));
+       if (!fFileTR) {
+           AliWarning(Form("AliMCEventHandler:TrackRefs%s.root not found in directory %s ! \n", fExtension, fPathName->Data()));
+           fInitOk = kFALSE;
+           return kFALSE;
+       }
     }
     
-    return ok;
+    fInitOk = kTRUE;
+
+    return kTRUE;
 }
 
-Bool_t AliMCEventHandler::BeginEvent()
+Bool_t AliMCEventHandler::BeginEvent(Long64_t entry)
 { 
+    // Begin event
+    fParticleSelected.Delete();
+    fLabelMap.Delete();
     // Read the next event
-    fEvent++;
-    if (fEvent >= fNEvent) {
-       AliWarning(Form("AliMCEventHandler: Event number out of range %5d\n", fEvent));
+
+    if (fEventsInContainer != 0) {
+       entry = (Long64_t) ( entry * Float_t(fNEvent) / Float_t (fEventsInContainer));
+    }
+
+
+    if (entry == -1) {
+       fEvent++;
+       entry = fEvent;
+    } else {
+       fEvent = entry;
+    }
+
+    if (entry >= fNEvent) {
+       AliWarning(Form("AliMCEventHandler: Event number out of range %5lld %5d\n", entry, fNEvent));
        return kFALSE;
     }
-    return GetEvent(fEvent);
+    
+    Bool_t result = GetEvent(entry);
+
+    if (fSubsidiaryHandlers) {
+       TIter next(fSubsidiaryHandlers);
+       AliMCEventHandler *handler;
+       while((handler = (AliMCEventHandler*)next())) {
+           handler->BeginEvent(entry);
+       }
+       next.Reset();
+       while((handler = (AliMCEventHandler*)next())) {
+           fMCEvent->AddSubsidiaryEvent(handler->MCEvent());
+       }
+       fMCEvent->InitEvent();
+    }
+    
+    if (fPreReadMode == kLmPreRead) {
+       fMCEvent->PreReadAll();
+    }
+
+    return result;
+    
+}
+
+void AliMCEventHandler::SelectParticle(Int_t i){
+  // taking the absolute values here, need to take care 
+  // of negative daughter and mother
+  // IDs when setting!
+    if (TMath::Abs(i) >= AliMCEvent::BgLabelOffset()) i =  fMCEvent->BgLabelToIndex(TMath::Abs(i));
+    if(!IsParticleSelected(TMath::Abs(i)))fParticleSelected.Add(TMath::Abs(i),1);
+}
+
+Bool_t AliMCEventHandler::IsParticleSelected(Int_t i)  {
+  // taking the absolute values here, need to take 
+  // care with negative daughter and mother
+  // IDs when setting!
+  return (fParticleSelected.GetValue(TMath::Abs(i))==1);
+}
+
+
+void AliMCEventHandler::CreateLabelMap(){
+
+  //
+  // this should be called once all selections where done 
+  //
+
+  fLabelMap.Delete();
+  if(!fMCEvent){
+    fParticleSelected.Delete();
+    return;
+  }
+
+  VerifySelectedParticles();
+
+  Int_t iNew = 0;
+  for(int i = 0;i < fMCEvent->GetNumberOfTracks();++i){
+    if(IsParticleSelected(i)){
+      fLabelMap.Add(i,iNew);
+      iNew++;
+    }
+  }
+}
+
+Int_t AliMCEventHandler::GetNewLabel(Int_t i) {
+  // Gets the label from the new created Map
+  // Call CreatLabelMap before
+  // otherwise only 0 returned
+  return fLabelMap.GetValue(TMath::Abs(i));
+}
+
+void  AliMCEventHandler::VerifySelectedParticles(){
+
+  //  
+  // Make sure that each particle has at least it's predecessors
+  // selected so that we have the complete ancestry tree
+  // Private, should be only called by CreateLabelMap
+
+  if(!fMCEvent){
+      fParticleSelected.Delete();
+      return;
+  }
+
+  Int_t nprim = fMCEvent->GetNumberOfPrimaries();
+
+  for(int i = 0;i < fMCEvent->GetNumberOfTracks(); ++i){
+      if(i < nprim){
+         SelectParticle(i);// take all primaries
+         continue;
+      }
+
+      if(!IsParticleSelected(i))continue;
+
+      AliMCParticle* mcpart = (AliMCParticle*) fMCEvent->GetTrack(i);
+      Int_t imo = mcpart->GetMother();
+      while((imo >= nprim)&&!IsParticleSelected(imo)){
+         // Mother not yet selected
+         SelectParticle(imo);
+         mcpart = (AliMCParticle*) fMCEvent->GetTrack(imo);
+         imo = mcpart->GetMother();
+      }
+    // after last step we may have an unselected primary
+    // mother
+    if(imo>=0){
+      if(!IsParticleSelected(imo))
+       SelectParticle(imo);
+    } 
+  }// loop over all tracks
 }
 
 Int_t AliMCEventHandler::GetParticleAndTR(Int_t i, TParticle*& particle, TClonesArray*& trefs)
 {
     // Retrieve entry i
-    if (i < 0 || i >= fNparticles) {
-       AliWarning(Form("AliMCEventHandler::GetEntry: Index out of range"));
-       particle = 0;
-       trefs    = 0;
-       return (-1);
-    }
-    particle = fStack->Particle(i);
-    if (fFileTR) {
-       fTreeTR->GetEntry(fStack->TreeKEntry(i));
-       trefs    = fTrackReferences;
-       return trefs->GetEntries();
+    if (!fInitOk) {
+       return 0;
     } else {
-       trefs = 0;
-       return -1;
+       return (fMCEvent->GetParticleAndTR(i, particle, trefs));
     }
 }
 
 void AliMCEventHandler::DrawCheck(Int_t i, Int_t search)
 {
     // Retrieve entry i and draw momentum vector and hits
-    if (!fFileTR) {
-       AliWarning("No Track Reference information available");
-       return;
-    } 
-
-    if (i > -1 && i < fNparticles) {
-       fTreeTR->GetEntry(fStack->TreeKEntry(i));
-    } else {
-       AliWarning("AliMCEventHandler::GetEntry: Index out of range");
-    }
-    
-    Int_t nh = fTrackReferences->GetEntries();
-    
-    
-    if (search) {
-       while(nh <= search && i < fNparticles - 1) {
-           i++;
-           fTreeTR->GetEntry(fStack->TreeKEntry(i));
-           nh =  fTrackReferences->GetEntries();
-       }
-       printf("Found Hits at %5d\n", i);
-    }
-    TParticle* particle = fStack->Particle(i);
-
-    TH2F*    h = new TH2F("", "", 100, -500, 500, 100, -500, 500);
-    Float_t x0 = particle->Vx();
-    Float_t y0 = particle->Vy();
-
-    Float_t x1 = particle->Vx() + particle->Px() * 50.;
-    Float_t y1 = particle->Vy() + particle->Py() * 50.;
-    
-    TArrow*  a = new TArrow(x0, y0, x1, y1, 0.01);
-    h->Draw();
-    a->SetLineColor(2);
-    
-    a->Draw();
-    
-    for (Int_t ih = 0; ih < nh; ih++) {
-       AliTrackReference* ref = (AliTrackReference*) fTrackReferences->At(ih);
-       TMarker* m = new TMarker(ref->X(), ref->Y(), 20);
-       m->Draw();
-       m->SetMarkerSize(0.4);
-       
-    }
+    fMCEvent->DrawCheck(i, search);
 }
 
 Bool_t AliMCEventHandler::Notify(const char *path)
 {
-    // Notify about directory change
-    // The directory is taken from the 'path' argument
-    // Reconnect trees
+  // Notify about directory change
+  // The directory is taken from the 'path' argument
+  // Reconnect trees
+    TString fileName(path);
+    if(fileName.Contains("AliESDs.root")){
+       fileName.ReplaceAll("AliESDs.root", "");
+    }
+    else if(fileName.Contains("AliAOD.root")){
+       fileName.ReplaceAll("AliAOD.root", "");
+    }
+    else if(fileName.Contains("galice.root")){
+       // for running with galice and kinematics alone...
+       fileName.ReplaceAll("galice.root", "");
+    }
+    else if (fileName.BeginsWith("root:")) {
+      fileName.Append("?ZIP=");
+    }
 
-    printf("AliMCEventHandler::Notify() file: %s\n", path);
-    fPathName = Form("%s",  path);
+    *fPathName = fileName;
+    AliInfo(Form("Path: -%s-\n", fPathName->Data()));
+    
     ResetIO();
     InitIO("");
+
+// Handle subsidiary handlers
+    if (fSubsidiaryHandlers) {
+       TIter next(fSubsidiaryHandlers);
+       AliMCEventHandler *handler;
+       while((handler = (AliMCEventHandler*) next())) {
+           TString* spath = handler->GetInputPath();
+           if (spath->Contains("merged")) {
+               if (! fPathName->IsNull()) {
+                   handler->Notify(Form("%s/../.", fPathName->Data()));
+               } else {
+                   handler->Notify("../");
+               }
+           }
+       }
+    }
+    
     return kTRUE;
 }
-    
+
 void AliMCEventHandler::ResetIO()
 {
-    // Reset files
-    if (fFileE)  delete fFileE;
-    if (fFileK)  delete fFileK;
-    if (fFileTR) delete fFileTR;
+//  Clear header and stack
+    
+    if (fInitOk) fMCEvent->Clean();
+    
+// Delete Tree E
+    delete fTreeE; fTreeE = 0;
+     
+// Reset files
+    if (fFileE)  {delete fFileE;  fFileE  = 0;}
+    if (fFileK)  {delete fFileK;  fFileK  = 0;}
+    if (fFileTR) {delete fFileTR; fFileTR = 0;}
+    fExtension="";
+    fInitOk = kFALSE;
+
+    if (fSubsidiaryHandlers) {
+       TIter next(fSubsidiaryHandlers);
+       AliMCEventHandler *handler;
+       while((handler = (AliMCEventHandler*)next())) {
+           handler->ResetIO();
+       }
+    }
+
 }
 
                            
 Bool_t AliMCEventHandler::FinishEvent()
 {
-    // Reset the stack 
-    Stack()->Reset();
-    
+    // Clean-up after each event
+    delete fDirTR;  fDirTR = 0;
+    delete fDirK;   fDirK  = 0;    
+    if (fInitOk) fMCEvent->FinishEvent();
+
+    if (fSubsidiaryHandlers) {
+       TIter next(fSubsidiaryHandlers);
+       AliMCEventHandler *handler;
+       while((handler = (AliMCEventHandler*)next())) {
+           handler->FinishEvent();
+       }
+    }
+
     return kTRUE;
 }
 
@@ -334,198 +515,18 @@ Bool_t AliMCEventHandler::TerminateIO()
     return kTRUE;
 }
     
-void AliMCEventHandler::ReorderAndExpandTreeTR()
-{
-//
-//  Reorder and expand the track reference tree in order to match the kinematics tree.
-//  Copy the information from different branches into one
-//
-//  TreeTR
-    if (fTmpTreeTR) {
-       fTmpTreeTR->Delete("all");
-    }
-    
-    if (fTmpFileTR) {
-       fTmpFileTR->Close();
-       delete fTmpFileTR;
-    }
 
-    fTmpFileTR = new TFile("TrackRefsTmp.root", "recreate");
-    fTmpTreeTR = new TTree("TreeTR", "Track References");
-    if (!fTrackReferences)  fTrackReferences = new TClonesArray("AliTrackReference", 100);
-    fTmpTreeTR->Branch("TrackReferences", "TClonesArray", &fTrackReferences, 4000);
-
-//
-    fTreeTR->SetBranchStatus("*",      0);
-    fTreeTR->SetBranchStatus("AliRun", 1);
-    fTreeTR->SetBranchStatus("ITS",    1);
-    fTreeTR->SetBranchStatus("TPC",    1);
-    fTreeTR->SetBranchStatus("TRD",    1);
-    fTreeTR->SetBranchStatus("TOF",    1);
-    fTreeTR->SetBranchStatus("FRAME",  1);
-    fTreeTR->SetBranchStatus("MUON",   1);
-
-    TClonesArray* trefs[7];
-    for (Int_t i = 0; i < 7; i++) trefs[i] = 0;
-    if (fTreeTR){
-       // make branch for central track references
-       if (fTreeTR->GetBranch("AliRun")) fTreeTR->SetBranchAddress("AliRun", &trefs[0]);
-       if (fTreeTR->GetBranch("ITS"))    fTreeTR->SetBranchAddress("ITS",    &trefs[1]);
-       if (fTreeTR->GetBranch("TPC"))    fTreeTR->SetBranchAddress("TPC",    &trefs[2]);
-       if (fTreeTR->GetBranch("TRD"))    fTreeTR->SetBranchAddress("TRD",    &trefs[3]);
-       if (fTreeTR->GetBranch("TOF"))    fTreeTR->SetBranchAddress("TOF",    &trefs[4]);
-       if (fTreeTR->GetBranch("FRAME"))  fTreeTR->SetBranchAddress("FRAME",  &trefs[5]);
-       if (fTreeTR->GetBranch("MUON"))   fTreeTR->SetBranchAddress("MUON",   &trefs[6]);
-    }
+void AliMCEventHandler::SetInputPath(const char* fname)
+{
+    // Set the input path name
+    delete fPathName;
+    fPathName = new TString(fname);
+}
 
-    Int_t np = fStack->GetNprimary();
-    Int_t nt = fTreeTR->GetEntries();
-    //
-    // Loop over tracks and find the secondaries with the help of the kine tree
-    Int_t ifills = 0;
-    Int_t it     = 0;
-    Int_t itlast = 0;
-    
-    for (Int_t ip = np - 1; ip > -1; ip--) {
-       TParticle *part = fStack->Particle(ip);
-//     printf("Particle %5d %5d %5d %5d %5d %5d \n", 
-//            ip, part->GetPdgCode(), part->GetFirstMother(), part->GetFirstDaughter(), 
-//            part->GetLastDaughter(), part->TestBit(kTransportBit));
-
-       // Determine range of secondaries produced by this primary during transport     
-       Int_t dau1  = part->GetFirstDaughter();
-       if (dau1 < np) continue;  // This particle has no secondaries produced during transport
-       Int_t dau2  = -1;
-       if (dau1 > -1) {
-           Int_t inext = ip - 1;
-           while (dau2 < 0) {
-               if (inext >= 0) {
-                   part = fStack->Particle(inext);
-                   dau2 =  part->GetFirstDaughter();
-                   if (dau2 == -1 || dau2 < np) {
-                       dau2 = -1;
-                   } else {
-                       dau2--;
-                   }
-               } else {
-                   dau2 = fStack->GetNtrack() - 1;
-               }
-               inext--;
-           } // find upper bound
-       }  // dau2 < 0
-       
-
-//     printf("Check (1) %5d %5d %5d %5d %5d \n", ip, np, it, dau1, dau2);
-//
-// Loop over reference hits and find secondary label
-// First the tricky part: find the entry in treeTR than contains the hits or
-// make sure that no hits exist.
-//
-       Bool_t hasHits   = kFALSE;
-       Bool_t isOutside = kFALSE;
-
-       it = itlast;
-       while (!hasHits && !isOutside && it < nt) {
-           fTreeTR->GetEntry(it++);
-           for (Int_t ib = 0; ib < 7; ib++) {
-               if (!trefs[ib]) continue;
-               Int_t nh = trefs[ib]->GetEntries();
-               for (Int_t ih = 0; ih < nh; ih++) {
-                   AliTrackReference* tr = (AliTrackReference*) trefs[ib]->At(ih);
-                   Int_t label = tr->Label();
-                   if (label >= dau1 && label <= dau2) {
-                       hasHits = kTRUE;
-                       itlast = it - 1;
-                       break;
-                   }
-                   if (label > dau2 || label < ip) {
-                       isOutside = kTRUE;
-                       itlast = it - 1;
-                       break;
-                   }
-               } // hits
-               if (hasHits || isOutside) break;
-           } // branches
-       } // entries
-
-       if (!hasHits) {
-           // Write empty entries
-           for (Int_t id = dau1; (id <= dau2); id++) {
-               fTmpTreeTR->Fill();
-               ifills++;
-           } 
-       } else {
-           // Collect all hits
-           fTreeTR->GetEntry(itlast);
-           for (Int_t id = dau1; (id <= dau2) && (dau1 > -1); id++) {
-               for (Int_t ib = 0; ib < 7; ib++) {
-                   if (!trefs[ib]) continue;
-                   Int_t nh = trefs[ib]->GetEntries();
-                   for (Int_t ih = 0; ih < nh; ih++) {
-                       AliTrackReference* tr = (AliTrackReference*) trefs[ib]->At(ih);
-                       Int_t label = tr->Label();
-                       // Skip primaries
-                       if (label == ip) continue;
-                       if (label > dau2 || label < dau1) 
-                           printf("AliMCEventHandler::Track Reference Label out of range !: %5d %5d %5d %5d \n", 
-                                  itlast, label, dau1, dau2);
-                       if (label == id) {
-                           // secondary found
-                           tr->SetDetectorId(ib-1);
-                           Int_t nref =  fTrackReferences->GetEntriesFast();
-                           TClonesArray &lref = *fTrackReferences;
-                           new(lref[nref]) AliTrackReference(*tr);
-                       }
-                   } // hits
-               } // branches
-               fTmpTreeTR->Fill();
-               fTrackReferences->Clear();
-               ifills++;
-           } // daughters
-       } // has hits
-    } // tracks
-    //
-    // Now loop again and write the primaries
-    //
-    it = nt - 1;
-    for (Int_t ip = 0; ip < np; ip++) {
-       Int_t labmax = -1;
-       while (labmax < ip && it > -1) {
-           fTreeTR->GetEntry(it--);
-           for (Int_t ib = 0; ib < 7; ib++) {
-               if (!trefs[ib]) continue;
-               Int_t nh = trefs[ib]->GetEntries();
-               // 
-               // Loop over reference hits and find primary labels
-               for (Int_t ih = 0; ih < nh; ih++) {
-                   AliTrackReference* tr = (AliTrackReference*)  trefs[ib]->At(ih);
-                   Int_t label = tr->Label();
-                   if (label < np && label > labmax) {
-                       labmax = label;
-                   }
-                   
-                   if (label == ip) {
-                       tr->SetDetectorId(ib-1);
-                       Int_t nref = fTrackReferences->GetEntriesFast();
-                       TClonesArray &lref = *fTrackReferences;
-                       new(lref[nref]) AliTrackReference(*tr);
-                   }
-               } // hits
-           } // branches
-       } // entries
-       it++;
-       fTmpTreeTR->Fill();
-       fTrackReferences->Clear();
-       ifills++;
-    } // tracks
-    // Check
-    delete fTreeTR;
-    for (Int_t ib = 0; ib < 7; ib++) {
-       if (trefs[ib]) delete trefs[ib];
-    }
+void AliMCEventHandler::AddSubsidiaryHandler(AliMCEventHandler* handler)
+{
+    // Add a subsidiary handler. For example for background events
 
-    if (ifills != fStack->GetNtrack()) 
-       printf("AliMCEventHandler:Number of entries in TreeTR (%5d) unequal to TreeK (%5d) \n", 
-              ifills, fStack->GetNtrack());
-    fTreeTR = fTmpTreeTR;
+    if (!fSubsidiaryHandlers) fSubsidiaryHandlers = new TList();
+    fSubsidiaryHandlers->Add(handler);
 }