1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
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.
23 // Origin: Andreas Morsch, CERN, andreas.morsch@cern.ch
24 //---------------------------------------------------------------------------------
28 #include "AliMCEventHandler.h"
29 #include "AliMCEvent.h"
30 #include "AliMCParticle.h"
32 #include "AliTrackReference.h"
33 #include "AliHeader.h"
40 #include <TParticle.h>
42 #include <TClonesArray.h>
43 #include <TDirectoryFile.h>
45 ClassImp(AliMCEventHandler)
47 AliMCEventHandler::AliMCEventHandler() :
49 fMCEvent(new AliMCEvent()),
62 fPathName(new TString("./")),
68 fSubsidiaryHandlers(0),
69 fEventsInContainer(0),
70 fPreReadMode(kNoPreRead)
73 // Default constructor
75 // Be sure to add all particles to the PDG database
76 AliPDG::AddParticlesToPdgDataBase();
79 AliMCEventHandler::AliMCEventHandler(const char* name, const char* title) :
80 AliVEventHandler(name, title),
81 fMCEvent(new AliMCEvent()),
94 fPathName(new TString("./")),
100 fSubsidiaryHandlers(0),
101 fEventsInContainer(0),
102 fPreReadMode(kNoPreRead)
107 // Be sure to add all particles to the PDG database
108 AliPDG::AddParticlesToPdgDataBase();
110 AliMCEventHandler::~AliMCEventHandler()
120 Bool_t AliMCEventHandler::Init(Option_t* opt)
124 if (!(strcmp(opt, "proof")) || !(strcmp(opt, "local"))) return kTRUE;
126 fFileE = TFile::Open(Form("%sgalice.root", fPathName->Data()));
128 AliError(Form("AliMCEventHandler:galice.root not found in directory %s ! \n", fPathName->Data()));
135 fFileE->GetObject("TE", fTreeE);
136 // Connect Tree E to the MCEvent
137 fMCEvent->ConnectTreeE(fTreeE);
138 fNEvent = fTreeE->GetEntries();
141 fFileK = TFile::Open(Form("%sKinematics%s.root", fPathName->Data(), fExtension));
143 AliError(Form("AliMCEventHandler:Kinematics.root not found in directory %s ! \n", fPathName->Data()));
148 fEventsPerFile = fFileK->GetNkeys() - fFileK->GetNProcessIDs();
152 fFileTR = TFile::Open(Form("%sTrackRefs%s.root", fPathName->Data(), fExtension));
154 AliError(Form("AliMCEventHandler:TrackRefs.root not found in directory %s ! \n", fPathName->Data()));
160 // Reset the event number
163 AliInfo(Form("Number of events in this directory %5d \n", fNEvent));
167 if (fSubsidiaryHandlers) {
168 TIter next(fSubsidiaryHandlers);
169 AliMCEventHandler *handler;
170 while((handler = (AliMCEventHandler*)next())) {
172 handler->SetNumberOfEventsInContainer(fNEvent);
179 Bool_t AliMCEventHandler::GetEvent(Int_t iev)
181 // Load the event number iev
183 // Calculate the file number
184 if (!fInitOk) return kFALSE;
186 Int_t inew = iev / fEventsPerFile;
187 if (inew != fFileNumber) {
189 if (!OpenFile(fFileNumber)){
195 snprintf(folder, 20, "Event%d", iev);
197 fTreeE->GetEntry(iev);
199 fFileK->GetObject(folder, fDirK);
201 AliWarning(Form("AliMCEventHandler: Event #%5d - Cannot get kinematics\n", iev));
205 fDirK ->GetObject("TreeK", fTreeK);
207 AliError(Form("AliMCEventHandler: Event #%5d - Cannot get TreeK\n",iev));
210 // Connect TreeK to MCEvent
211 fMCEvent->ConnectTreeK(fTreeK);
214 // Check which format has been read
215 fFileTR->GetObject(folder, fDirTR);
217 AliError(Form("AliMCEventHandler: Event #%5d - Cannot get track references\n",iev));
221 fDirTR->GetObject("TreeTR", fTreeTR);
224 AliError(Form("AliMCEventHandler: Event #%5d - Cannot get TreeTR\n",iev));
227 // Connect TR to MCEvent
228 fMCEvent->ConnectTreeTR(fTreeTR);
235 Bool_t AliMCEventHandler::OpenFile(Int_t i)
239 fExtension = Form("%d", i);
246 fFileK = TFile::Open(Form("%sKinematics%s.root", fPathName->Data(), fExtension));
248 AliError(Form("AliMCEventHandler:Kinematics%s.root not found in directory %s ! \n", fExtension, fPathName->Data()));
255 fFileTR = TFile::Open(Form("%sTrackRefs%s.root", fPathName->Data(), fExtension));
257 AliWarning(Form("AliMCEventHandler:TrackRefs%s.root not found in directory %s ! \n", fExtension, fPathName->Data()));
268 Bool_t AliMCEventHandler::BeginEvent(Long64_t entry)
271 fParticleSelected.Delete();
273 // Read the next event
275 if (fEventsInContainer != 0) {
276 entry = (Long64_t) ( entry * Float_t(fNEvent) / Float_t (fEventsInContainer));
287 if (entry >= fNEvent) {
288 AliWarning(Form("AliMCEventHandler: Event number out of range %5lld %5d\n", entry, fNEvent));
292 Bool_t result = GetEvent(entry);
294 if (fSubsidiaryHandlers) {
295 TIter next(fSubsidiaryHandlers);
296 AliMCEventHandler *handler;
297 while((handler = (AliMCEventHandler*)next())) {
298 handler->BeginEvent(entry);
301 while((handler = (AliMCEventHandler*)next())) {
302 fMCEvent->AddSubsidiaryEvent(handler->MCEvent());
304 fMCEvent->InitEvent();
307 if (fPreReadMode == kLmPreRead) {
308 fMCEvent->PreReadAll();
315 void AliMCEventHandler::SelectParticle(Int_t i){
316 // taking the absolute values here, need to take care
317 // of negative daughter and mother
319 if (TMath::Abs(i) >= AliMCEvent::BgLabelOffset()) i = fMCEvent->BgLabelToIndex(TMath::Abs(i));
320 if(!IsParticleSelected(TMath::Abs(i)))fParticleSelected.Add(TMath::Abs(i),1);
323 Bool_t AliMCEventHandler::IsParticleSelected(Int_t i) {
324 // taking the absolute values here, need to take
325 // care with negative daughter and mother
327 return (fParticleSelected.GetValue(TMath::Abs(i))==1);
331 void AliMCEventHandler::CreateLabelMap(){
334 // this should be called once all selections where done
339 fParticleSelected.Delete();
343 VerifySelectedParticles();
346 for(int i = 0;i < fMCEvent->GetNumberOfTracks();++i){
347 if(IsParticleSelected(i)){
348 fLabelMap.Add(i,iNew);
354 Int_t AliMCEventHandler::GetNewLabel(Int_t i) {
355 // Gets the label from the new created Map
356 // Call CreatLabelMap before
357 // otherwise only 0 returned
358 return fLabelMap.GetValue(TMath::Abs(i));
361 void AliMCEventHandler::VerifySelectedParticles(){
364 // Make sure that each particle has at least it's predecessors
365 // selected so that we have the complete ancestry tree
366 // Private, should be only called by CreateLabelMap
369 fParticleSelected.Delete();
373 Int_t nprim = fMCEvent->GetNumberOfPrimaries();
375 for(int i = 0;i < fMCEvent->GetNumberOfTracks(); ++i){
377 SelectParticle(i);// take all primaries
381 if(!IsParticleSelected(i))continue;
383 AliMCParticle* mcpart = (AliMCParticle*) fMCEvent->GetTrack(i);
384 Int_t imo = mcpart->GetMother();
385 while((imo >= nprim)&&!IsParticleSelected(imo)){
386 // Mother not yet selected
388 mcpart = (AliMCParticle*) fMCEvent->GetTrack(imo);
389 imo = mcpart->GetMother();
391 // after last step we may have an unselected primary
394 if(!IsParticleSelected(imo))
397 }// loop over all tracks
400 Int_t AliMCEventHandler::GetParticleAndTR(Int_t i, TParticle*& particle, TClonesArray*& trefs)
406 return (fMCEvent->GetParticleAndTR(i, particle, trefs));
410 void AliMCEventHandler::DrawCheck(Int_t i, Int_t search)
412 // Retrieve entry i and draw momentum vector and hits
413 fMCEvent->DrawCheck(i, search);
416 Bool_t AliMCEventHandler::Notify(const char *path)
418 // Notify about directory change
419 // The directory is taken from the 'path' argument
421 TString fileName(path);
422 if(fileName.Contains("AliESDs.root")){
423 fileName.ReplaceAll("AliESDs.root", "");
425 else if(fileName.Contains("AliESDs_wSDD.root")){
426 fileName.ReplaceAll("AliESDs_wSDD.root", "");
428 else if(fileName.Contains("AliAOD.root")){
429 fileName.ReplaceAll("AliAOD.root", "");
431 else if(fileName.Contains("galice.root")){
432 // for running with galice and kinematics alone...
433 fileName.ReplaceAll("galice.root", "");
435 else if (fileName.BeginsWith("root:")) {
436 fileName.Append("?ZIP=");
439 *fPathName = fileName;
440 AliInfo(Form("Path: -%s-\n", fPathName->Data()));
445 // Handle subsidiary handlers
446 if (fSubsidiaryHandlers) {
447 TIter next(fSubsidiaryHandlers);
448 AliMCEventHandler *handler;
449 while((handler = (AliMCEventHandler*) next())) {
450 TString* spath = handler->GetInputPath();
451 if (spath->Contains("merged")) {
452 if (! fPathName->IsNull()) {
453 handler->Notify(Form("%s/../.", fPathName->Data()));
455 handler->Notify("../");
464 void AliMCEventHandler::ResetIO()
466 // Clear header and stack
468 if (fInitOk) fMCEvent->Clean();
471 delete fTreeE; fTreeE = 0;
474 if (fFileE) {delete fFileE; fFileE = 0;}
475 if (fFileK) {delete fFileK; fFileK = 0;}
476 if (fFileTR) {delete fFileTR; fFileTR = 0;}
480 if (fSubsidiaryHandlers) {
481 TIter next(fSubsidiaryHandlers);
482 AliMCEventHandler *handler;
483 while((handler = (AliMCEventHandler*)next())) {
491 Bool_t AliMCEventHandler::FinishEvent()
493 // Clean-up after each event
494 delete fDirTR; fDirTR = 0;
495 delete fDirK; fDirK = 0;
496 if (fInitOk) fMCEvent->FinishEvent();
498 if (fSubsidiaryHandlers) {
499 TIter next(fSubsidiaryHandlers);
500 AliMCEventHandler *handler;
501 while((handler = (AliMCEventHandler*)next())) {
502 handler->FinishEvent();
509 Bool_t AliMCEventHandler::Terminate()
515 Bool_t AliMCEventHandler::TerminateIO()
522 void AliMCEventHandler::SetInputPath(const char* fname)
524 // Set the input path name
526 fPathName = new TString(fname);
529 void AliMCEventHandler::AddSubsidiaryHandler(AliMCEventHandler* handler)
531 // Add a subsidiary handler. For example for background events
533 if (!fSubsidiaryHandlers) fSubsidiaryHandlers = new TList();
534 fSubsidiaryHandlers->Add(handler);