[u/mrichter/AliRoot.git] / STEER / AliMCEventHandler.cxx

/************************************************************************** 
 * 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.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/

/* $Id$ */
//---------------------------------------------------------------------------------
//                          Class AliMCEventHandler
// This class gives access to MC truth during the analysis.
// Monte Carlo truth is containe in the kinematics tree (produced particles) and 
// the tree of reference hits.
//      
// Origin: Andreas Morsch, CERN, andreas.morsch@cern.ch 
//---------------------------------------------------------------------------------


#include "AliMCEventHandler.h"
#include "AliMCEvent.h"
#include "AliMCParticle.h"
#include "AliPDG.h"
#include "AliTrackReference.h"
#include "AliHeader.h"
#include "AliStack.h"
#include "AliLog.h"

#include <TTree.h>
#include <TFile.h>
#include <TList.h>
#include <TParticle.h>
#include <TString.h>
#include <TClonesArray.h>
#include <TDirectoryFile.h>

ClassImp(AliMCEventHandler)

AliMCEventHandler::AliMCEventHandler() :
    AliVEventHandler(),
    fMCEvent(new AliMCEvent()),
    fFileE(0),
    fFileK(0),
    fFileTR(0),
    fTreeE(0),
    fTreeK(0),
    fTreeTR(0),
    fDirK(0),
    fDirTR(0),
    fParticleSelected(0),
    fLabelMap(0),
    fNEvent(-1),
    fEvent(-1),
    fPathName(new TString("./")),
    fExtension(""),
    fFileNumber(0),
    fEventsPerFile(0),
    fReadTR(kTRUE),
    fInitOk(kFALSE),
    fSubsidiaryHandlers(0),
    fEventsInContainer(0),
    fPreReadMode(kNoPreRead)
{
  //
  // 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),
    fTreeE(0),
    fTreeK(0),
    fTreeTR(0),
    fDirK(0),
    fDirTR(0),
    fParticleSelected(0),
    fLabelMap(0),
    fNEvent(-1),
    fEvent(-1),
    fPathName(new TString("./")),
    fExtension(""),
    fFileNumber(0),
    fEventsPerFile(0),
    fReadTR(kTRUE),
    fInitOk(kFALSE),
    fSubsidiaryHandlers(0),
    fEventsInContainer(0),
    fPreReadMode(kNoPreRead)
{
  //
  // Constructor
  //
  // Be sure to add all particles to the PDG database
  AliPDG::AddParticlesToPdgDataBase();
}
AliMCEventHandler::~AliMCEventHandler()
{ 
    // Destructor
    delete fMCEvent;
    delete fFileE;
    delete fFileK;
    delete fFileTR;
}

Bool_t AliMCEventHandler::Init(Option_t* opt)
{ 
    // Initialize input
    //
    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);
    // Connect Tree E to the MCEvent
    fMCEvent->ConnectTreeE(fTreeE);
    fNEvent = fTreeE->GetEntries();
    //
    // Tree K
    fFileK = TFile::Open(Form("%sKinematics%s.root", fPathName->Data(), fExtension));
    if (!fFileK) {
	AliError(Form("AliMCEventHandler:Kinematics.root not found in directory %s ! \n", fPathName));
	fInitOk = kFALSE;
	return kTRUE;
    }
    
    fEventsPerFile = fFileK->GetNkeys() - fFileK->GetNProcessIDs();
    //
    // Tree TR
    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("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)
{
    // Load the event number iev
    //
    // Calculate the file number
  if (!fInitOk) return kFALSE;
    
  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);
  // Tree K
  fFileK->GetObject(folder, fDirK);
  if (!fDirK) {
    AliWarning(Form("AliMCEventHandler: Event #%5d - Cannot get kinematics\n", iev));
    return kFALSE;
  }
    
  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
    if (i > 0) {
	fExtension = Form("%d", i);
    } else {
	fExtension = "";
    }
    
    
    delete fFileK;
    fFileK = TFile::Open(Form("%sKinematics%s.root", fPathName->Data(), fExtension));
    if (!fFileK) {
	AliError(Form("AliMCEventHandler:Kinematics%s.root not found in directory %s ! \n", fExtension, fPathName->Data()));
	fInitOk = kFALSE;
	return 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;
	}
    }
    
    fInitOk = kTRUE;

    return kTRUE;
}

Bool_t AliMCEventHandler::BeginEvent(Long64_t entry)
{ 
    fParticleSelected.Delete();
    fLabelMap.Delete();
    // Read the next event

    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 %5d %5d\n", entry, fNEvent));
	return kFALSE;
    }
    
    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 (!fInitOk) {
	return 0;
    } else {
	return (fMCEvent->GetParticleAndTR(i, particle, trefs));
    }
}

void AliMCEventHandler::DrawCheck(Int_t i, Int_t search)
{
    // Retrieve entry i and draw momentum vector and hits
    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
    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=");
    }

    *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()
{
//  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()
{
    // 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;
}

Bool_t AliMCEventHandler::Terminate()
{ 
    // Dummy 
    return kTRUE;
}

Bool_t AliMCEventHandler::TerminateIO()
{ 
    // Dummy
    return kTRUE;
}
    

void AliMCEventHandler::SetInputPath(const char* fname)
{
    // Set the input path name
    delete fPathName;
    fPathName = new TString(fname);
}

void AliMCEventHandler::AddSubsidiaryHandler(AliMCEventHandler* handler)
{
    // Add a subsidiary handler. For example for background events

    if (!fSubsidiaryHandlers) fSubsidiaryHandlers = new TList();
    fSubsidiaryHandlers->Add(handler);
}
Commit	Line	Data
	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	// Class AliMCEventHandler
	19	// This class gives access to MC truth during the analysis.
	20	// Monte Carlo truth is containe in the kinematics tree (produced particles) and
	21	// the tree of reference hits.
	22	//
	23	// Origin: Andreas Morsch, CERN, andreas.morsch@cern.ch
	24	//---------------------------------------------------------------------------------
	25
	26
	27
	28	#include "AliMCEventHandler.h"
	29	#include "AliMCEvent.h"
	30	#include "AliMCParticle.h"
	31	#include "AliPDG.h"
	32	#include "AliTrackReference.h"
	33	#include "AliHeader.h"
	34	#include "AliStack.h"
	35	#include "AliLog.h"
	36
	37	#include <TTree.h>
	38	#include <TFile.h>
	39	#include <TList.h>
	40	#include <TParticle.h>
	41	#include <TString.h>
	42	#include <TClonesArray.h>
	43	#include <TDirectoryFile.h>
	44
	45	ClassImp(AliMCEventHandler)
	46
	47	AliMCEventHandler::AliMCEventHandler() :
	48	AliVEventHandler(),
	49	fMCEvent(new AliMCEvent()),
	50	fFileE(0),
	51	fFileK(0),
	52	fFileTR(0),
	53	fTreeE(0),
	54	fTreeK(0),
	55	fTreeTR(0),
	56	fDirK(0),
	57	fDirTR(0),
	58	fParticleSelected(0),
	59	fLabelMap(0),
	60	fNEvent(-1),
	61	fEvent(-1),
	62	fPathName(new TString("./")),
	63	fExtension(""),
	64	fFileNumber(0),
	65	fEventsPerFile(0),
	66	fReadTR(kTRUE),
	67	fInitOk(kFALSE),
	68	fSubsidiaryHandlers(0),
	69	fEventsInContainer(0),
	70	fPreReadMode(kNoPreRead)
	71	{
	72	//
	73	// Default constructor
	74	//
	75	// Be sure to add all particles to the PDG database
	76	AliPDG::AddParticlesToPdgDataBase();
	77	}
	78
	79	AliMCEventHandler::AliMCEventHandler(const char* name, const char* title) :
	80	AliVEventHandler(name, title),
	81	fMCEvent(new AliMCEvent()),
	82	fFileE(0),
	83	fFileK(0),
	84	fFileTR(0),
	85	fTreeE(0),
	86	fTreeK(0),
	87	fTreeTR(0),
	88	fDirK(0),
	89	fDirTR(0),
	90	fParticleSelected(0),
	91	fLabelMap(0),
	92	fNEvent(-1),
	93	fEvent(-1),
	94	fPathName(new TString("./")),
	95	fExtension(""),
	96	fFileNumber(0),
	97	fEventsPerFile(0),
	98	fReadTR(kTRUE),
	99	fInitOk(kFALSE),
	100	fSubsidiaryHandlers(0),
	101	fEventsInContainer(0),
	102	fPreReadMode(kNoPreRead)
	103	{
	104	//
	105	// Constructor
	106	//
	107	// Be sure to add all particles to the PDG database
	108	AliPDG::AddParticlesToPdgDataBase();
	109	}
	110	AliMCEventHandler::~AliMCEventHandler()
	111	{
	112	// Destructor
	113	delete fMCEvent;
	114	delete fFileE;
	115	delete fFileK;
	116	delete fFileTR;
	117	}
	118
	119	Bool_t AliMCEventHandler::Init(Option_t* opt)
	120	{
	121	// Initialize input
	122	//
	123	if (!(strcmp(opt, "proof")) \|\| !(strcmp(opt, "local"))) return kTRUE;
	124	//
	125	fFileE = TFile::Open(Form("%sgalice.root", fPathName->Data()));
	126	if (!fFileE) {
	127	AliError(Form("AliMCEventHandler:galice.root not found in directory %s ! \n", fPathName->Data()));
	128	fInitOk = kFALSE;
	129	return kFALSE;
	130	}
	131
	132	//
	133	// Tree E
	134	fFileE->GetObject("TE", fTreeE);
	135	// Connect Tree E to the MCEvent
	136	fMCEvent->ConnectTreeE(fTreeE);
	137	fNEvent = fTreeE->GetEntries();
	138	//
	139	// Tree K
	140	fFileK = TFile::Open(Form("%sKinematics%s.root", fPathName->Data(), fExtension));
	141	if (!fFileK) {
	142	AliError(Form("AliMCEventHandler:Kinematics.root not found in directory %s ! \n", fPathName));
	143	fInitOk = kFALSE;
	144	return kTRUE;
	145	}
	146
	147	fEventsPerFile = fFileK->GetNkeys() - fFileK->GetNProcessIDs();
	148	//
	149	// Tree TR
	150	if (fReadTR) {
	151	fFileTR = TFile::Open(Form("%sTrackRefs%s.root", fPathName->Data(), fExtension));
	152	if (!fFileTR) {
	153	AliError(Form("AliMCEventHandler:TrackRefs.root not found in directory %s ! \n", fPathName->Data()));
	154	fInitOk = kFALSE;
	155	return kTRUE;
	156	}
	157	}
	158	//
	159	// Reset the event number
	160	fEvent = -1;
	161	fFileNumber = 0;
	162	AliInfo(Form("Number of events in this directory %5d \n", fNEvent));
	163	fInitOk = kTRUE;
	164
	165
	166	if (fSubsidiaryHandlers) {
	167	TIter next(fSubsidiaryHandlers);
	168	AliMCEventHandler *handler;
	169	while((handler = (AliMCEventHandler*)next())) {
	170	handler->Init(opt);
	171	handler->SetNumberOfEventsInContainer(fNEvent);
	172	}
	173	}
	174
	175	return kTRUE;
	176	}
	177
	178	Bool_t AliMCEventHandler::GetEvent(Int_t iev)
	179	{
	180	// Load the event number iev
	181	//
	182	// Calculate the file number
	183	if (!fInitOk) return kFALSE;
	184
	185	Int_t inew = iev / fEventsPerFile;
	186	if (inew != fFileNumber) {
	187	fFileNumber = inew;
	188	if (!OpenFile(fFileNumber)){
	189	return kFALSE;
	190	}
	191	}
	192	// Folder name
	193	char folder[20];
	194	sprintf(folder, "Event%d", iev);
	195	// TreeE
	196	fTreeE->GetEntry(iev);
	197	// Tree K
	198	fFileK->GetObject(folder, fDirK);
	199	if (!fDirK) {
	200	AliWarning(Form("AliMCEventHandler: Event #%5d - Cannot get kinematics\n", iev));
	201	return kFALSE;
	202	}
	203
	204	fDirK ->GetObject("TreeK", fTreeK);
	205	if (!fTreeK) {
	206	AliError(Form("AliMCEventHandler: Event #%5d - Cannot get TreeK\n",iev));
	207	return kFALSE;
	208	}
	209	// Connect TreeK to MCEvent
	210	fMCEvent->ConnectTreeK(fTreeK);
	211	//Tree TR
	212	if (fFileTR) {
	213	// Check which format has been read
	214	fFileTR->GetObject(folder, fDirTR);
	215	if (!fDirTR) {
	216	AliError(Form("AliMCEventHandler: Event #%5d - Cannot get track references\n",iev));
	217	return kFALSE;
	218	}
	219
	220	fDirTR->GetObject("TreeTR", fTreeTR);
	221	//
	222	if (!fTreeTR) {
	223	AliError(Form("AliMCEventHandler: Event #%5d - Cannot get TreeTR\n",iev));
	224	return kFALSE;
	225	}
	226	// Connect TR to MCEvent
	227	fMCEvent->ConnectTreeTR(fTreeTR);
	228	}
	229
	230	//
	231	return kTRUE;
	232	}
	233
	234	Bool_t AliMCEventHandler::OpenFile(Int_t i)
	235	{
	236	// Open file i
	237	if (i > 0) {
	238	fExtension = Form("%d", i);
	239	} else {
	240	fExtension = "";
	241	}
	242
	243
	244	delete fFileK;
	245	fFileK = TFile::Open(Form("%sKinematics%s.root", fPathName->Data(), fExtension));
	246	if (!fFileK) {
	247	AliError(Form("AliMCEventHandler:Kinematics%s.root not found in directory %s ! \n", fExtension, fPathName->Data()));
	248	fInitOk = kFALSE;
	249	return kFALSE;
	250	}
	251
	252	if (fReadTR) {
	253	delete fFileTR;
	254	fFileTR = TFile::Open(Form("%sTrackRefs%s.root", fPathName->Data(), fExtension));
	255	if (!fFileTR) {
	256	AliWarning(Form("AliMCEventHandler:TrackRefs%s.root not found in directory %s ! \n", fExtension, fPathName->Data()));
	257	fInitOk = kFALSE;
	258	return kFALSE;
	259	}
	260	}
	261
	262	fInitOk = kTRUE;
	263
	264	return kTRUE;
	265	}
	266
	267	Bool_t AliMCEventHandler::BeginEvent(Long64_t entry)
	268	{
	269	fParticleSelected.Delete();
	270	fLabelMap.Delete();
	271	// Read the next event
	272
	273	if (fEventsInContainer != 0) {
	274	entry = (Long64_t) ( entry * Float_t(fNEvent) / Float_t (fEventsInContainer));
	275	}
	276
	277
	278	if (entry == -1) {
	279	fEvent++;
	280	entry = fEvent;
	281	} else {
	282	fEvent = entry;
	283	}
	284
	285	if (entry >= fNEvent) {
	286	AliWarning(Form("AliMCEventHandler: Event number out of range %5d %5d\n", entry, fNEvent));
	287	return kFALSE;
	288	}
	289
	290	Bool_t result = GetEvent(entry);
	291
	292	if (fSubsidiaryHandlers) {
	293	TIter next(fSubsidiaryHandlers);
	294	AliMCEventHandler *handler;
	295	while((handler = (AliMCEventHandler*)next())) {
	296	handler->BeginEvent(entry);
	297	}
	298	next.Reset();
	299	while((handler = (AliMCEventHandler*)next())) {
	300	fMCEvent->AddSubsidiaryEvent(handler->MCEvent());
	301	}
	302	fMCEvent->InitEvent();
	303	}
	304
	305	if (fPreReadMode == kLmPreRead) {
	306	fMCEvent->PreReadAll();
	307	}
	308
	309	return result;
	310
	311	}
	312
	313	void AliMCEventHandler::SelectParticle(Int_t i){
	314	// taking the absolute values here, need to take care
	315	// of negative daughter and mother
	316	// IDs when setting!
	317	if (TMath::Abs(i) >= AliMCEvent::BgLabelOffset()) i = fMCEvent->BgLabelToIndex(TMath::Abs(i));
	318	if(!IsParticleSelected(TMath::Abs(i)))fParticleSelected.Add(TMath::Abs(i),1);
	319	}
	320
	321	Bool_t AliMCEventHandler::IsParticleSelected(Int_t i) {
	322	// taking the absolute values here, need to take
	323	// care with negative daughter and mother
	324	// IDs when setting!
	325	return (fParticleSelected.GetValue(TMath::Abs(i))==1);
	326	}
	327
	328
	329	void AliMCEventHandler::CreateLabelMap(){
	330
	331	//
	332	// this should be called once all selections where done
	333	//
	334
	335	fLabelMap.Delete();
	336	if(!fMCEvent){
	337	fParticleSelected.Delete();
	338	return;
	339	}
	340
	341	VerifySelectedParticles();
	342
	343	Int_t iNew = 0;
	344	for(int i = 0;i < fMCEvent->GetNumberOfTracks();++i){
	345	if(IsParticleSelected(i)){
	346	fLabelMap.Add(i,iNew);
	347	iNew++;
	348	}
	349	}
	350	}
	351
	352	Int_t AliMCEventHandler::GetNewLabel(Int_t i) {
	353	// Gets the label from the new created Map
	354	// Call CreatLabelMap before
	355	// otherwise only 0 returned
	356	return fLabelMap.GetValue(TMath::Abs(i));
	357	}
	358
	359	void AliMCEventHandler::VerifySelectedParticles(){
	360
	361	//
	362	// Make sure that each particle has at least it's predecessors
	363	// selected so that we have the complete ancestry tree
	364	// Private, should be only called by CreateLabelMap
	365
	366	if(!fMCEvent){
	367	fParticleSelected.Delete();
	368	return;
	369	}
	370
	371	Int_t nprim = fMCEvent->GetNumberOfPrimaries();
	372
	373	for(int i = 0;i < fMCEvent->GetNumberOfTracks(); ++i){
	374	if(i < nprim){
	375	SelectParticle(i);// take all primaries
	376	continue;
	377	}
	378
	379	if(!IsParticleSelected(i))continue;
	380
	381	AliMCParticle* mcpart = (AliMCParticle*) fMCEvent->GetTrack(i);
	382	Int_t imo = mcpart->GetMother();
	383	while((imo >= nprim)&&!IsParticleSelected(imo)){
	384	// Mother not yet selected
	385	SelectParticle(imo);
	386	mcpart = (AliMCParticle*) fMCEvent->GetTrack(imo);
	387	imo = mcpart->GetMother();
	388	}
	389	// after last step we may have an unselected primary
	390	// mother
	391	if(imo>=0){
	392	if(!IsParticleSelected(imo))
	393	SelectParticle(imo);
	394	}
	395	}// loop over all tracks
	396	}
	397
	398	Int_t AliMCEventHandler::GetParticleAndTR(Int_t i, TParticle& particle, TClonesArray& trefs)
	399	{
	400	// Retrieve entry i
	401	if (!fInitOk) {
	402	return 0;
	403	} else {
	404	return (fMCEvent->GetParticleAndTR(i, particle, trefs));
	405	}
	406	}
	407
	408	void AliMCEventHandler::DrawCheck(Int_t i, Int_t search)
	409	{
	410	// Retrieve entry i and draw momentum vector and hits
	411	fMCEvent->DrawCheck(i, search);
	412	}
	413
	414	Bool_t AliMCEventHandler::Notify(const char *path)
	415	{
	416	// Notify about directory change
	417	// The directory is taken from the 'path' argument
	418	// Reconnect trees
	419	TString fileName(path);
	420	if(fileName.Contains("AliESDs.root")){
	421	fileName.ReplaceAll("AliESDs.root", "");
	422	}
	423	else if(fileName.Contains("AliAOD.root")){
	424	fileName.ReplaceAll("AliAOD.root", "");
	425	}
	426	else if(fileName.Contains("galice.root")){
	427	// for running with galice and kinematics alone...
	428	fileName.ReplaceAll("galice.root", "");
	429	}
	430	else if (fileName.BeginsWith("root:")) {
	431	fileName.Append("?ZIP=");
	432	}
	433
	434	*fPathName = fileName;
	435	AliInfo(Form("Path: -%s-\n", fPathName->Data()));
	436
	437	ResetIO();
	438	InitIO("");
	439
	440	// Handle subsidiary handlers
	441	if (fSubsidiaryHandlers) {
	442	TIter next(fSubsidiaryHandlers);
	443	AliMCEventHandler *handler;
	444	while((handler = (AliMCEventHandler*) next())) {
	445	TString* spath = handler->GetInputPath();
	446	if (spath->Contains("merged")) {
	447	if (! fPathName->IsNull()) {
	448	handler->Notify(Form("%s/../.", fPathName->Data()));
	449	} else {
	450	handler->Notify("../");
	451	}
	452	}
	453	}
	454	}
	455
	456	return kTRUE;
	457	}
	458
	459	void AliMCEventHandler::ResetIO()
	460	{
	461	// Clear header and stack
	462
	463	if (fInitOk) fMCEvent->Clean();
	464
	465	// Delete Tree E
	466	delete fTreeE; fTreeE = 0;
	467
	468	// Reset files
	469	if (fFileE) {delete fFileE; fFileE = 0;}
	470	if (fFileK) {delete fFileK; fFileK = 0;}
	471	if (fFileTR) {delete fFileTR; fFileTR = 0;}
	472	fExtension="";
	473	fInitOk = kFALSE;
	474
	475	if (fSubsidiaryHandlers) {
	476	TIter next(fSubsidiaryHandlers);
	477	AliMCEventHandler *handler;
	478	while((handler = (AliMCEventHandler*)next())) {
	479	handler->ResetIO();
	480	}
	481	}
	482
	483	}
	484
	485
	486	Bool_t AliMCEventHandler::FinishEvent()
	487	{
	488	// Clean-up after each event
	489	delete fDirTR; fDirTR = 0;
	490	delete fDirK; fDirK = 0;
	491	if (fInitOk) fMCEvent->FinishEvent();
	492
	493	if (fSubsidiaryHandlers) {
	494	TIter next(fSubsidiaryHandlers);
	495	AliMCEventHandler *handler;
	496	while((handler = (AliMCEventHandler*)next())) {
	497	handler->FinishEvent();
	498	}
	499	}
	500
	501	return kTRUE;
	502	}
	503
	504	Bool_t AliMCEventHandler::Terminate()
	505	{
	506	// Dummy
	507	return kTRUE;
	508	}
	509
	510	Bool_t AliMCEventHandler::TerminateIO()
	511	{
	512	// Dummy
	513	return kTRUE;
	514	}
	515
	516
	517	void AliMCEventHandler::SetInputPath(const char* fname)
	518	{
	519	// Set the input path name
	520	delete fPathName;
	521	fPathName = new TString(fname);
	522	}
	523
	524	void AliMCEventHandler::AddSubsidiaryHandler(AliMCEventHandler* handler)
	525	{
	526	// Add a subsidiary handler. For example for background events
	527
	528	if (!fSubsidiaryHandlers) fSubsidiaryHandlers = new TList();
	529	fSubsidiaryHandlers->Add(handler);
	530	}