/**************************************************************************
* 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$ */
//_________________________________________________________________________
// Algorithm class for the reconstruction:
//
//
//*--
//*-- Author: Gines Martinez & Yves Schutz (SUBATECH)
//*-- Complitely redisigned by Dmitri Peressounko (SUBATECH & RRC KI) March 2001
/////////////////////////////////////////////////////////////////////////////////////
// Supervising class for reconstruction
// use case:
//
// root [0] AliPHOSReconstructioner * r = new AliPHOSReconstructioner("galice.root")
// // Set the header file
// root [1] r->ExecuteTask()
//
// // One can specify, to what file each branch should be written
// root [2] r->SetBranchFileName("RecPoints","PHOS.RecPoints1.root") ;
// // By default branches are stored in galice.root (in non-split mode)
// // or PHOS.SDigits.root, PHOS.Digits.root etc.
// // Note, that if you set already existing names of branch files,
// // these branches will be OVERWRITTEN!
//
// // One can specify the starting point of the reconstruction
// root [3] r->StartFrom("AliPHOSClusterizer")
// // means that SDigits and Digits will not be regenerated, only RecPoints,
// // TS and RecParticles
//
// // And finally one can call ExecuteTask() with the following options
// root [4] r->ExecuteTask("debug all timing")
// // deb - prints the numbers of produced SDigits, Digits etc.
// // deb all - prints in addition list of made SDigits, digits etc.
// // timing - prints benchmarking results
// --- ROOT system ---
#include "TClonesArray.h"
#include "TROOT.h"
#include "TTree.h"
#include "TSystem.h"
// --- Standard library ---
#include <iostream.h>
// --- AliRoot header files ---
#include "AliRun.h"
#include "AliPHOSReconstructioner.h"
#include "AliPHOSClusterizerv1.h"
#include "AliPHOSDigitizer.h"
#include "AliPHOSSDigitizer.h"
#include "AliPHOSTrackSegmentMakerv1.h"
#include "AliPHOSPIDv1.h"
#include "AliPHOSFastRecParticle.h"
#include "AliPHOSCpvRecPoint.h"
ClassImp(AliPHOSReconstructioner)
//____________________________________________________________________________
AliPHOSReconstructioner::AliPHOSReconstructioner():TTask("AliPHOSReconstructioner","")
{
// ctor
fDigitizer = 0 ;
fClusterizer = 0 ;
fTSMaker = 0 ;
fPID = 0 ;
fSDigitizer = 0 ;
fIsInitialized = kFALSE ;
}
//____________________________________________________________________________
AliPHOSReconstructioner::AliPHOSReconstructioner(const char* headerFile):TTask("AliPHOSReconstructioner","")
{
// ctor
fHeaderFileName = headerFile ;
if(gSystem->Getenv("CONFIG_SPLIT_FILE")) //generating file name
fSDigitsBranch="PHOS.SDigits.root" ;
else
fSDigitsBranch="" ;
fSDigitizer = new AliPHOSSDigitizer(fHeaderFileName.Data(),fSDigitsBranch.Data()) ;
Add(fSDigitizer) ;
if(gSystem->Getenv("CONFIG_SPLIT_FILE")) //generating file name
fDigitsBranch="PHOS.Digits.root" ;
else
fDigitsBranch="" ;
fDigitizer = new AliPHOSDigitizer(fHeaderFileName.Data(),fDigitsBranch.Data()) ;
Add(fDigitizer) ;
if(gSystem->Getenv("CONFIG_SPLIT_FILE")) //generating file name
fRecPointBranch="PHOS.Reco.root" ;
else
fRecPointBranch="" ;
fClusterizer = new AliPHOSClusterizerv1(fHeaderFileName.Data(),fRecPointBranch.Data()) ;
Add(fClusterizer) ;
if(gSystem->Getenv("CONFIG_SPLIT_FILE")) //generating file name
fTSBranch="PHOS.Reco.root" ;
else
fTSBranch="" ;
fTSMaker = new AliPHOSTrackSegmentMakerv1(fHeaderFileName.Data(),fTSBranch.Data()) ;
Add(fTSMaker) ;
if(gSystem->Getenv("CONFIG_SPLIT_FILE")) //generating file name
fRecPartBranch="PHOS.Reco.root" ;
else
fRecPartBranch="" ;
fPID = new AliPHOSPIDv1(fHeaderFileName.Data(),fRecPartBranch.Data()) ;
Add(fPID) ;
fIsInitialized = kTRUE ;
}
//____________________________________________________________________________
void AliPHOSReconstructioner::Init()
{
//initiase Reconstructioner if necessary: we can not do this in default constructor
if(!fIsInitialized){
// Initialisation
if(gSystem->Getenv("CONFIG_SPLIT_FILE")) //generating file name
fSDigitsBranch="PHOS.SDigits.root" ;
else
fSDigitsBranch="" ;
fSDigitizer = new AliPHOSSDigitizer(fHeaderFileName.Data(),fSDigitsBranch.Data()) ;
Add(fSDigitizer) ;
if(gSystem->Getenv("CONFIG_SPLIT_FILE")) //generating file name
fDigitsBranch="PHOS.Digits.root" ;
else
fDigitsBranch="" ;
fDigitizer = new AliPHOSDigitizer(fHeaderFileName.Data(),fDigitsBranch.Data()) ;
Add(fDigitizer) ;
if(gSystem->Getenv("CONFIG_SPLIT_FILE")) //generating file name
fRecPointBranch="PHOS.Reco.root" ;
else
fRecPointBranch="" ;
fClusterizer = new AliPHOSClusterizerv1(fHeaderFileName.Data(),fRecPointBranch.Data()) ;
Add(fClusterizer) ;
if(gSystem->Getenv("CONFIG_SPLIT_FILE")) //generating file name
fTSBranch="PHOS.Reco.root" ;
else
fTSBranch="" ;
fTSMaker = new AliPHOSTrackSegmentMakerv1(fHeaderFileName.Data(),fTSBranch.Data()) ;
Add(fTSMaker) ;
if(gSystem->Getenv("CONFIG_SPLIT_FILE")) //generating file name
fRecPartBranch="PHOS.Reco.root" ;
else
fRecPartBranch="" ;
fPID = new AliPHOSPIDv1(fHeaderFileName.Data(),fRecPartBranch.Data()) ;
Add(fPID) ;
fIsInitialized = kTRUE ;
}
}
//____________________________________________________________________________
AliPHOSReconstructioner::~AliPHOSReconstructioner()
{
if(fSDigitizer)
delete fSDigitizer ;
if(fDigitizer)
delete fDigitizer ;
if(fClusterizer)
delete fClusterizer ;
if(fTSMaker)
delete fTSMaker ;
if(fPID)
delete fPID ;
}
//____________________________________________________________________________
void AliPHOSReconstructioner::SetBranchFileName(const char* branch, const char * fileName){
//Diverge correcpoinding branch to the file "fileName"
if(strcmp(branch,"SDigits") == 0){
fSDigitizer->SetSDigitsBranch(fileName) ;
fDigitizer->SetSDigitsBranch(fileName) ;
fSDigitsBranch = fileName ;
return ;
}
if(strcmp(branch,"Digits") == 0){
fDigitizer->SetDigitsBranch(fileName) ;
fClusterizer->SetDigitsBranch(fileName) ;
fDigitsBranch = fileName ;
return ;
}
if(strcmp(branch,"RecPoints") == 0){
fClusterizer->SetRecPointsBranch(fileName) ;
fTSMaker->SetRecPointsBranch(fileName) ;
fRecPointBranch = fileName ;
return ;
}
if(strcmp(branch,"TrackSegments") == 0){
fTSMaker->SetTrackSegmentsBranch(fileName) ;
fPID->SetTrackSegmentsBranch(fileName) ;
fTSBranch = fileName ;
return ;
}
if(strcmp(branch,"RecParticles") == 0){
fPID->SetRecParticlesBranch(fileName) ;
fRecPartBranch = fileName ;
return ;
}
cout << "There is no branch " << branch << "!"<< endl ;
cout << "Available branches `SDigits', `Digits', `RecPoints', `TrackSegments' and `RecParticles' " << endl ;
}
//____________________________________________________________________________
void AliPHOSReconstructioner::StartFrom(Option_t * module){
//in the next ExecuteTask() reconstruction starts from the module "module"
if(!fIsInitialized)
Init() ;
TIter next(fTasks);
TTask *task;
Bool_t active = kFALSE ;
while((task=(TTask*)next())){
if (strcmp(module,task->GetName())==0)
active = kTRUE;
task->SetActive(active) ;
}
if(!active){
cout << "There is no task " <<module<< endl ;
cout << "Available tasks are: " << endl ;
next.Reset() ;
while((task=(TTask*)next()))
cout<<" " << task->GetName() << endl ;
}
}
//____________________________________________________________________________
void AliPHOSReconstructioner::Print(Option_t * option)const {
cout << "-----------------AliPHOSReconstructioner---------------" << endl ;
cout << " Reconstruction of the header file " <<fHeaderFileName.Data() << endl ;
cout << " with the following modules: " << endl ;
if(fSDigitizer->IsActive())
cout << " (+) " << fSDigitizer->GetName() << " to branch : " << fSDigitsBranch.Data() << endl ;
else
cout << " (-) " << fSDigitizer->GetName() << " to branch : " << fSDigitsBranch.Data() << endl ;
cout << endl ;
if(fDigitizer->IsActive())
cout << " (+) " << fDigitizer->GetName() << " to branch : " << fDigitsBranch.Data() << endl ;
else
cout << " (-) " << fDigitizer->GetName() << " to branch : " << fDigitsBranch.Data() << endl ;
cout << endl ;
if(fClusterizer->IsActive())
cout << " (+) " <<fClusterizer->GetName() << " to branch : " <<fRecPointBranch.Data() << endl ;
else
cout << " (-) " <<fClusterizer->GetName() << " to branch : " <<fRecPointBranch.Data() << endl ;
cout << endl ;
if(fTSMaker->IsActive())
cout << " (+) " << fTSMaker->GetName() << " to branch : " << fTSBranch.Data() << endl ;
else
cout << " (-) " << fTSMaker->GetName() << " to branch : " << fTSBranch.Data() << endl ;
cout << endl ;
if(fPID)
cout << " (+) " << fPID->GetName() << " to branch : " <<fRecPartBranch.Data() << endl ;
else
cout << " (-) " << fPID->GetName() << " to branch : " <<fRecPartBranch.Data() << endl ;
cout << endl ;
}
ROOT page - Class index - Top of the page
This page has been automatically generated. If you have any comments or suggestions about the page layout send a mail to ROOT support, or contact the developers with any questions or problems regarding ROOT.