From 7acf600847d7c3706c34409560cce379580d954e Mon Sep 17 00:00:00 2001 From: schutz Date: Wed, 4 Apr 2001 13:21:56 +0000 Subject: [PATCH] new design: derived from TTask --- PHOS/AliPHOSPID.cxx | 18 +- PHOS/AliPHOSPID.h | 43 +- PHOS/AliPHOSPIDv1.cxx | 669 ++++++++++++++++++++++++----- PHOS/AliPHOSPIDv1.h | 79 +++- PHOS/AliPHOSReconstructioner.cxx | 694 +++++++++++++++++++------------ PHOS/AliPHOSReconstructioner.h | 80 ++-- PHOS/AliPHOSSDigitizer.cxx | 307 ++++++++++---- PHOS/AliPHOSSDigitizer.h | 36 +- 8 files changed, 1385 insertions(+), 541 deletions(-) diff --git a/PHOS/AliPHOSPID.cxx b/PHOS/AliPHOSPID.cxx index a36d5f1091a..92ea9aca051 100644 --- a/PHOS/AliPHOSPID.cxx +++ b/PHOS/AliPHOSPID.cxx @@ -23,7 +23,7 @@ // the code checker // -//*-- Author: Yves Schutz (SUBATECH) +//*-- Author: Yves Schutz (SUBATECH) & Dmitri Peressounko // --- ROOT system --- @@ -31,7 +31,7 @@ // --- Standard library --- - +#include // --- AliRoot header files --- @@ -41,12 +41,22 @@ ClassImp(AliPHOSPID) //____________________________________________________________________________ -AliPHOSPID::AliPHOSPID() + AliPHOSPID::AliPHOSPID():TTask() { // ctor - fGeom = AliPHOSGeometry::GetInstance() ; } +//____________________________________________________________________________ +AliPHOSPID::AliPHOSPID(const char* header,const char * branch ):TTask() +{ + // ctor + cout << "AliPHOSPID: This constructor should be overwritten! "<< endl ; + abort() ; + +} + + + //____________________________________________________________________________ AliPHOSPID::~AliPHOSPID() diff --git a/PHOS/AliPHOSPID.h b/PHOS/AliPHOSPID.h index 5d0d82ae1a6..9e82f64529b 100644 --- a/PHOS/AliPHOSPID.h +++ b/PHOS/AliPHOSPID.h @@ -13,38 +13,47 @@ // --- ROOT system --- -#include "TObject.h" -#include "TClonesArray.h" +#include "TTask.h" +class TFormula ; +class TClonesArray ; // --- Standard library --- // --- AliRoot header files --- -#include "AliPHOSTrackSegment.h" -#include "AliPHOSRecParticle.h" -#include "AliPHOSGeometry.h" +class AliPHOSGeometry ; +class AliPHOSClusterizer ; +class AliPHOSTrackSegmentMaker ; - - -class AliPHOSPID : public TObject { +class AliPHOSPID : public TTask { public: AliPHOSPID() ; // ctor + AliPHOSPID(const char* headerFile,const char * tsBranch = 0) ; virtual ~AliPHOSPID() ; // dtor - virtual void MakeParticles(AliPHOSTrackSegment::TrackSegmentsList * trsl, - AliPHOSRecParticle::RecParticlesList * rpl) {} ; - virtual void Print(const char *){} ; - virtual void SetShowerProfileCuts(Float_t, Float_t, Float_t, Float_t) {} ; - virtual void SetDispersionCutOff(Float_t ) {} ; - virtual void SetRelativeDistanceCut(Float_t ) {}; + virtual void Exec(Option_t * option) = 0 ; + virtual char * GetRecParticlesBranch()const = 0 ; + virtual char * GetTrackSegmentsBranch()const = 0 ; + virtual void Init()= 0 ; + + virtual void Print(Option_t * option) const = 0 ; + virtual void PlotDispersionCuts()const = 0; + virtual Bool_t ReadTrackSegments()= 0 ; + + virtual void SetIdentificationMethod(char * option = "CPV DISP" ) = 0 ; + + virtual void SetShowerProfileCut(char * formula) = 0 ; + virtual void SetDispersionCut(Float_t cut) = 0 ; + virtual void SetCpvtoEmcDistanceCut(Float_t cut ) = 0; - protected: - - AliPHOSGeometry * fGeom ; // pointer to PHOS geometry + virtual void SetTrackSegmentsBranch(const char* title) = 0 ; + virtual void SetRecParticlesBranch (const char* title) = 0 ; + virtual void WriteRecParticles()= 0 ; +protected: ClassDef(AliPHOSPID,1) // Particle Identifier algorithm (base class) diff --git a/PHOS/AliPHOSPIDv1.cxx b/PHOS/AliPHOSPIDv1.cxx index 5ba64fb2a6a..e8909e142d6 100644 --- a/PHOS/AliPHOSPIDv1.cxx +++ b/PHOS/AliPHOSPIDv1.cxx @@ -17,166 +17,633 @@ //_________________________________________________________________________ // Implementation version v1 of the PHOS particle identifier -// Identification is based on information from PPSD and EMC -// Why should I put meaningless comments -// just to satisfy -// the code checker -// -//*-- Author: Yves Schutz (SUBATECH) & Gines Martinez (SUBATECH) +// Particle identification based on the +// - CPV information, +// - Preshower information (in MIX or GPS2 geometries) +// - shower width. +// CPV or Preshower cluster should be clother in PHOS plane than fCpvEmcDistance (in cm). +// This variable can be set by method SetCpvtoEmcDistanceCut(Float_t cut) +// +// One can set desirable ID method by the function SetIdentificationMethod(option). +// Now the following options can be used together or separately : +// - "disp": use dispersion cut on shower width +// (width can be set by method SetDispersionCut(Float_t cut) +// - "ell" : use cut on the axis of the ellipse, drown around shower +// (this cut can be changed by SetShowerProfileCut(char* formula), +// where formula - any function of two variables f(lambda[0],lambda[1]). +// Shower is considered as EM if f() > 0 ) +// One can visualize current cuts calling method PlotDispersionCuts(). +// +// Below we present usercase: +// root [0] AliPHOSPIDv1 * p1 = new AliPHOSPIDv1("galice.root") +// Warning in : object already instantiated +// root [1] p1->SetIdentificationMethod("disp ellipse") +// root [2] p1->ExecuteTask() +// root [3] AliPHOSPIDv1 * p2 = new AliPHOSPIDv1("galice1.root","ts1") +// Warning in : object already instantiated +// // reading headers from file galice1.root and TrackSegments +// // with title "ts1" +// root [4] p2->SetRecParticlesBranch("rp1") +// // set file name for the branch RecParticles +// root [5] p2->ExecuteTask("deb all time") +// // available options +// // "deb" - prints # of reconstructed particles +// // "deb all" - prints # and list of RecParticles +// // "time" - prints benchmarking results +// +//*-- Author: Yves Schutz (SUBATECH) & Gines Martinez (SUBATECH) & +// Dmitri Peressounko (SUBATECH & Kurchatov Institute) +// Complitely redesined by Dmitri Peressounko, March 2001 // --- ROOT system --- - +#include "TROOT.h" +#include "TTree.h" +#include "TFile.h" +#include "TF2.h" +#include "TFormula.h" +#include "TCanvas.h" +#include "TFolder.h" +#include "TSystem.h" +#include "TBenchmark.h" // --- Standard library --- #include +#include // --- AliRoot header files --- +#include "AliRun.h" +#include "AliGenerator.h" +#include "AliPHOS.h" #include "AliPHOSPIDv1.h" +#include "AliPHOSClusterizerv1.h" #include "AliPHOSTrackSegment.h" +#include "AliPHOSTrackSegmentMakerv1.h" #include "AliPHOSRecParticle.h" -#include "AliPHOSIndexToObject.h" ClassImp( AliPHOSPIDv1) //____________________________________________________________________________ AliPHOSPIDv1::AliPHOSPIDv1():AliPHOSPID() { - fCutOnDispersion = 2.0; - fCutOnRelativeDistance = 3.0 ; + fIsInitialized = kFALSE ; +} + +//____________________________________________________________________________ +AliPHOSPIDv1::AliPHOSPIDv1(const char * headeFile,const char * tsBranchTitle):AliPHOSPID() +{ + + fHeaderFileName = headeFile ; + + fTSTitle = tsBranchTitle ; + + SetName("AliPHOSPID") ; + SetTitle("version1") ; + + TFile * file = (TFile*) gROOT->GetFile(fHeaderFileName.Data() ) ; + + if(file == 0){ + file = new TFile(fHeaderFileName.Data(),"update") ; + gAlice = (AliRun *) file->Get("gAlice") ; + } + + AliPHOS * phos = (AliPHOS *) gAlice->GetDetector("PHOS") ; + fGeom = AliPHOSGeometry::GetInstance(phos->GetGeometry()->GetName(),phos->GetGeometry()->GetTitle() ); + + fTrackSegments = new TClonesArray("AliPHOSTrackSegment",1) ; + fTSMaker = 0 ; + fEmcRecPoints = new TObjArray(1) ; + fCpvRecPoints = new TObjArray(1) ; + fClusterizer = 0 ; + fRecParticles = new TClonesArray("AliPHOSRecParticle",100) ; + + fFormula = new TFormula("LambdaCuts","(x>1)*(x<3)*(y>0)*(yGetRootFolder()->FindObject("Tasks") ; + roottasks->Add(this) ; + + fDispersion = 2.0; + fCpvEmcDistance = 3.0 ; + fIsInitialized = kTRUE ; + +} +//____________________________________________________________________________ +AliPHOSPIDv1::~AliPHOSPIDv1() +{ + +} +//____________________________________________________________________________ +void AliPHOSPIDv1::Init() +{ + if(!fIsInitialized){ + if(fHeaderFileName.IsNull()) + fHeaderFileName = "galice.root" ; + + TFile * file = (TFile*) gROOT->GetFile(fHeaderFileName.Data() ) ; + + if(file == 0){ + file = new TFile(fHeaderFileName.Data(),"update") ; + gAlice = (AliRun *) file->Get("gAlice") ; + } + + AliPHOS * phos = (AliPHOS *) gAlice->GetDetector("PHOS") ; + fGeom = AliPHOSGeometry::GetInstance(phos->GetGeometry()->GetName(),phos->GetGeometry()->GetTitle() ); + + fTrackSegments = new TClonesArray("AliPHOSTrackSegment",1) ; + fTSMaker = new AliPHOSTrackSegmentMakerv1() ; + fEmcRecPoints = new TObjArray(1) ; + fCpvRecPoints = new TObjArray(1) ; + fClusterizer = new AliPHOSClusterizerv1() ; + fRecParticles = new TClonesArray("AliPHOSRecParticle",100) ; + + fFormula = new TFormula("LambdaCuts","(x>1)*(x<2.5)*(y>0)*(yGetRootFolder()->FindObject("Tasks") ; + roottasks->Add(this) ; + + fDispersion = 2.0; + fCpvEmcDistance = 3.0 ; + + fIsInitialized = kTRUE ; + } + + } +//____________________________________________________________________________ +Bool_t AliPHOSPIDv1::ReadTrackSegments() +{ + //Fist read Track Segment Branch and extract RecPointsBranch from fTSMaker + + fTrackSegments->Clear() ; + fEmcRecPoints->Clear() ; + fCpvRecPoints->Clear() ; + fRecParticles->Clear() ; + + gAlice->GetEvent(fNEvent) ; + + TTree * treeR = gAlice->TreeR() ; + + if(treeR==0){ + char treeName[20]; + sprintf(treeName,"TreeR%d",fNEvent); + cout << "Error in AliPHOSClusterizerv1 : no "<GetListOfBranches() ; + Int_t ibranch; + Bool_t tsMakerNotFound = kTRUE ; + Bool_t tsNotFound = kTRUE ; + + for(ibranch = 0;(ibranch GetEntries())&&(tsMakerNotFound||tsNotFound);ibranch++){ + if(tsMakerNotFound){ + tsMakerBranch=(TBranch *) branches->At(ibranch) ; + if( fTSTitle.CompareTo(tsMakerBranch->GetTitle())==0 ) + if( strcmp(tsMakerBranch->GetName(),"AliPHOSTrackSegmentMaker") == 0) + tsMakerNotFound = kFALSE ; + } + if(tsNotFound){ + tsBranch=(TBranch *) branches->At(ibranch) ; + if( fTSTitle.CompareTo(tsBranch->GetTitle())==0 ) + if( strcmp(tsBranch->GetName(),"PHOSTS") == 0) + tsNotFound = kFALSE ; + } + } + + if(tsMakerNotFound ||tsNotFound ){ + cout << "Can't find Branch with TrackSegmentMaker and TrackSegments " ; + cout << "Do nothing" <SetAddress(&fTSMaker) ; + tsBranch->SetAddress(&fTrackSegments) ; + + treeR->GetEvent(0) ; + + fRecPointsTitle = fTSMaker->GetRecPointsBranch() ; + //reading now recponts branches + TBranch * emcBranch = 0; + TBranch * cpvBranch = 0; + TBranch * cluBranch = 0; + + Bool_t emcNotFound = kTRUE ; + Bool_t cpvNotFound = kTRUE ; + Bool_t cluNotFound = kTRUE ; + + for(ibranch = 0;(ibranch GetEntries())&&(emcNotFound||cpvNotFound||cluNotFound);ibranch++){ + if(emcNotFound){ + emcBranch=(TBranch *) branches->At(ibranch) ; + if( fRecPointsTitle.CompareTo(emcBranch->GetTitle())==0 ) + if( strcmp(emcBranch->GetName(),"PHOSEmcRP") == 0) + emcNotFound = kFALSE ; + } + if(cpvNotFound){ + cpvBranch=(TBranch *) branches->At(ibranch) ; + if( fRecPointsTitle.CompareTo(cpvBranch->GetTitle())==0 ) + if( strcmp(cpvBranch->GetName(),"PHOSCpvRP") == 0) + cpvNotFound = kFALSE ; + } + if(cluNotFound){ + cluBranch=(TBranch *) branches->At(ibranch) ; + if( fRecPointsTitle.CompareTo(cluBranch->GetTitle())==0 ) + if( strcmp(cluBranch->GetName(),"AliPHOSClusterizer") == 0) + cluNotFound = kFALSE ; + } + } + + if(emcNotFound ||cpvNotFound ||cluNotFound ){ + cout << "Can't find Branch with RecPoints or AliPHOSClusterizer " ; + cout << "Do nothing" <SetAddress(&fEmcRecPoints) ; + cpvBranch->SetAddress(&fCpvRecPoints) ; + cluBranch->SetAddress(&fClusterizer) ; + + treeR->GetEvent(0) ; + return kTRUE ; + + + +} //____________________________________________________________________________ -Float_t AliPHOSPIDv1::GetDistanceInPHOSPlane(AliPHOSEmcRecPoint * emcclu,AliPHOSPpsdRecPoint * PpsdClu, Bool_t &toofar, Option_t * Axis) +Float_t AliPHOSPIDv1::GetDistance(AliPHOSEmcRecPoint * emc,AliPHOSRecPoint * cpv, Option_t * Axis)const { // Calculates the distance between the EMC RecPoint and the PPSD RecPoint - Float_t r = 1.e+10; TVector3 vecEmc ; - TVector3 vecPpsd ; - - emcclu->GetLocalPosition(vecEmc) ; - PpsdClu->GetLocalPosition(vecPpsd) ; - if(emcclu->GetPHOSMod() == PpsdClu->GetPHOSMod()) - { - - // Correct to difference in CPV and EMC position due to different distance to center. - // we assume, that particle moves from center - Float_t dCPV = fGeom->GetIPtoOuterCoverDistance(); - Float_t dEMC = fGeom->GetIPtoCrystalSurface() ; - dEMC = dEMC / dCPV ; - vecPpsd = dEMC * vecPpsd - vecEmc ; - r = vecPpsd.Mag() ; - if (Axis == "X") r = vecPpsd.X(); - if (Axis == "Y") r = vecPpsd.Y(); - if (Axis == "Z") r = vecPpsd.Z(); - if (Axis == "R") r = vecPpsd.Mag(); - - } - else - { - toofar = kTRUE ; - } - return r ; + TVector3 vecCpv ; + + emc->GetLocalPosition(vecEmc) ; + cpv->GetLocalPosition(vecCpv) ; + if(emc->GetPHOSMod() == cpv->GetPHOSMod()){ + + // Correct to difference in CPV and EMC position due to different distance to center. + // we assume, that particle moves from center + Float_t dCPV = fGeom->GetIPtoOuterCoverDistance(); + Float_t dEMC = fGeom->GetIPtoCrystalSurface() ; + dEMC = dEMC / dCPV ; + vecCpv = dEMC * vecCpv - vecEmc ; + if (Axis == "X") return vecCpv.X(); + if (Axis == "Y") return vecCpv.Y(); + if (Axis == "Z") return vecCpv.Z(); + if (Axis == "R") return vecCpv.Mag(); + } + + return 100000000 ; } //____________________________________________________________________________ -void AliPHOSPIDv1::MakeParticles(AliPHOSTrackSegment::TrackSegmentsList * trsl, - AliPHOSRecParticle::RecParticlesList * rpl) +void AliPHOSPIDv1::Exec(Option_t * option) { + if(!fIsInitialized) + Init() ; + + if(strstr(option,"tim")) + gBenchmark->Start("PHOSPID"); + + + Int_t nEvents = (Int_t) gAlice->TreeE()->GetEntries() ; + + for(fNEvent = 0 ;fNEvent Stop("PHOSPID"); + cout << "AliPHOSPID:" << endl ; + cout << " took " << gBenchmark->GetCpuTime("PHOSPID") << " seconds for PID " + << gBenchmark->GetCpuTime("PHOSPID")/nEvents << " seconds per event " << endl ; + cout << endl ; + } + +} +//____________________________________________________________________________ +void AliPHOSPIDv1::MakeRecParticles(){ + // Makes a RecParticle out of a TrackSegment - TIter next(trsl) ; - AliPHOSTrackSegment * tracksegment ; + TIter next(fTrackSegments) ; + AliPHOSTrackSegment * ts ; Int_t index = 0 ; AliPHOSRecParticle * rp ; - Bool_t tDistance; - Int_t type ; - Int_t showerprofile; // 0 narrow and 1 wide - Int_t cpvdetector ; // 1 hit and 0 no hit - Int_t pcdetector ; // 1 hit and 0 no hit - - while ( (tracksegment = (AliPHOSTrackSegment *)next()) ) { - new( (*rpl)[index] ) AliPHOSRecParticle(tracksegment) ; - rp = (AliPHOSRecParticle *)rpl->At(index) ; - AliPHOSEmcRecPoint * recp = tracksegment->GetEmcRecPoint() ; - AliPHOSPpsdRecPoint * rpcpv = tracksegment->GetPpsdUpRecPoint() ; - AliPHOSPpsdRecPoint * rppc = tracksegment->GetPpsdLowRecPoint() ; + + Bool_t ellips = fIDOptions.Contains("ell",TString::kIgnoreCase ) ; + Bool_t disp = fIDOptions.Contains("dis",TString::kIgnoreCase ) ; + + while ( (ts = (AliPHOSTrackSegment *)next()) ) { - // Float_t * lambda = new Float_t[2]; - // recp->GetElipsAxis(lambda) ; + new( (*fRecParticles)[index] ) AliPHOSRecParticle() ; + rp = (AliPHOSRecParticle *)fRecParticles->At(index) ; + rp->SetTraskSegment(index) ; + + AliPHOSEmcRecPoint * emc = 0 ; + if(ts->GetEmcIndex()>=0) + emc = (AliPHOSEmcRecPoint *) fEmcRecPoints->At(ts->GetEmcIndex()) ; - // // Looking at the lateral development of the shower - // if ( ( lambda[0] > fLambda1m && lambda[0] < fLambda1M ) && // shower profile cut - // ( lambda[1] > fLambda2m && lambda[1] < fLambda2M ) ) - // // Float_t R ; - // //R=(lambda[0]-1.386)*(lambda[0]-1.386)+1.707*1.707*(lambda[1]-1.008)*(lambda[1]-1.008) ; - // //if(R<0.35*0.35) - - Float_t dispersion; - dispersion = recp->GetDispersion(); - if (dispersion < fCutOnDispersion) - showerprofile = 0 ; // NARROW PROFILE - else - showerprofile = 1 ;// WIDE PROFILE + AliPHOSRecPoint * cpv = 0 ; + if(ts->GetCpvIndex()>=0) + cpv = (AliPHOSRecPoint *) fCpvRecPoints->At(ts->GetCpvIndex()) ; + AliPHOSRecPoint * ppsd = 0 ; + if(ts->GetPpsdIndex()>=0) + ppsd= (AliPHOSRecPoint *) fCpvRecPoints->At(ts->GetPpsdIndex()) ; + + //set momentum and energy first + Float_t e = emc->GetEnergy() ; + TVector3 dir = GetMomentumDirection(emc,cpv,ppsd) ; + dir.SetMag(e) ; + + rp->SetMomentum(dir.X(),dir.Y(),dir.Z(),e) ; + rp->SetCalcMass(0); + + //now set type (reconstructed) of the particle + Int_t showerprofile = 0; // 0 narrow and 1 wide - // Looking at the photon conversion detector - if( tracksegment->GetPpsdLowRecPoint() == 0 ) - pcdetector = 0 ; // No hit - else{ - if (GetDistanceInPHOSPlane(recp, rppc, tDistance, "R") < fCutOnRelativeDistance) - pcdetector = 1 ; // hit - else - pcdetector = 0 ; + if(ellips){ + Float_t lambda[2] ; + emc->GetElipsAxis(lambda) ; + if(fFormula->Eval(lambda[0],lambda[1]) <= 0 ) + showerprofile = 1 ; // not narrow } + if(disp) + if(emc->GetDispersion() > fDispersion ) + showerprofile = 1 ; // not narrow + + // Looking at the photon conversion detector - if( tracksegment->GetPpsdUpRecPoint() == 0 ) - cpvdetector = 0 ; // No hit - else{ - if (GetDistanceInPHOSPlane(recp, rpcpv, tDistance, "R")< fCutOnRelativeDistance) - cpvdetector = 1 ; // Hit - else - cpvdetector = 0 ; - } - - type = showerprofile + 2 * pcdetector + 4 * cpvdetector ; + Int_t pcdetector= 0 ; //1 hit and 0 no hit + if(ppsd) + if(GetDistance(emc, ppsd, "R") < fCpvEmcDistance) + pcdetector = 1 ; + + // Looking at the CPV detector + Int_t cpvdetector= 0 ; //1 hit and 0 no hit + if(cpv) + if(GetDistance(emc, cpv, "R") < fCpvEmcDistance) + cpvdetector = 1 ; + + Int_t type = showerprofile + 2 * pcdetector + 4 * cpvdetector ; rp->SetType(type) ; index++ ; } - + } //____________________________________________________________________________ -void AliPHOSPIDv1:: Print(const char * opt) +void AliPHOSPIDv1:: Print(Option_t * option) const { // Print the parameters used for the particle type identification - cout << "AliPHOSPIDv1 : cuts for the particle idendification based on the shower profile " << endl - << fLambda1m << " < value1 < " << fLambda1M << endl - << fLambda2m << " < value2 < " << fLambda2M << endl ; + cout << "AliPHOSPIDv1 : cuts for the particle idendification based on the shower profile " << endl ; + + cout << "Eliptic cuts function " << endl ; + cout << " " << fFormula->GetTitle() << endl ; } //____________________________________________________________________________ -void AliPHOSPIDv1::SetShowerProfileCuts(Float_t l1m, Float_t l1M, Float_t l2m, Float_t l2M) +void AliPHOSPIDv1::SetShowerProfileCut(char * formula){ + //set shape of the cut on the axis of ellipce, drown around shouer + //shower considered "narrow" if Formula(lambda[0],lambda[1]) > 0. + if(fFormula) + delete fFormula; + fFormula = new TF2("Lambda Cut",formula,0,3,0,3) ; +} +//____________________________________________________________________________ +void AliPHOSPIDv1::WriteRecParticles() { - // Modifies the parameters used for the particle type identification - fLambda1m = l1m ; - fLambda1M = l1M ; - fLambda2m = l2m ; - fLambda2M = l2M ; -} + //check, if these branches already exist + TBranch * pidBranch = 0; + TBranch * rpBranch = 0; + + TObjArray * branches = gAlice->TreeR()->GetListOfBranches() ; + Int_t ibranch; + Bool_t pidNotFound = kTRUE ; + Bool_t rpNotFound = kTRUE ; + + for(ibranch = 0;(ibranch GetEntries())&& pidNotFound && rpNotFound;ibranch++){ + if(pidNotFound){ + pidBranch=(TBranch *) branches->At(ibranch) ; + if( (strcmp(pidBranch->GetName(),"PHOSPID") == 0) && + (fRecparticlesTitle.CompareTo(pidBranch->GetTitle()) == 0) ) + pidNotFound = kFALSE ; + } + if(rpNotFound){ + rpBranch=(TBranch *) branches->At(ibranch) ; + if( (strcmp(rpBranch->GetName(),"PHOSRP") == 0) && + (fRecparticlesTitle.CompareTo(rpBranch->GetTitle())==0 )) + rpNotFound = kFALSE ; + } + } + + if(!pidNotFound || !rpNotFound) { + cout << "AliPHOSPIDv1 error: " << endl ; + cout << " Branch PHOSRP and PHOSPID with title '"<Getenv("CONFIG_SPLIT_FILE")!=0){ //generating file name + filename = new char[strlen(gAlice->GetBaseFile())+20] ; + sprintf(filename,"%s/PHOS.Reco.root",gAlice->GetBaseFile()) ; + } + + TDirectory *cwd = gDirectory; + + //First rp + Int_t bufferSize = 32000 ; + rpBranch = gAlice->TreeR()->Branch("PHOSRP",&fRecParticles,bufferSize); + rpBranch->SetTitle(fRecparticlesTitle.Data()); + if (filename) { + rpBranch->SetFile(filename); + TIter next( rpBranch->GetListOfBranches()); + while ((rpBranch=(TBranch*)next())) { + rpBranch->SetFile(filename); + } + cwd->cd(); + } + + //second, pid + Int_t splitlevel = 0 ; + AliPHOSPIDv1 * pid = this ; + pidBranch = gAlice->TreeR()->Branch("AliPHOSPID","AliPHOSPIDv1",&pid,bufferSize,splitlevel); + pidBranch->SetTitle(fRecparticlesTitle.Data()); + if (filename) { + pidBranch->SetFile(filename); + TIter next( pidBranch->GetListOfBranches()); + while ((pidBranch=(TBranch*)next())) { + pidBranch->SetFile(filename); + } + cwd->cd(); + } + + gAlice->TreeR()->Fill() ; + gAlice->TreeR()->Write(0,kOverwrite) ; + +} //____________________________________________________________________________ -void AliPHOSPIDv1::SetRelativeDistanceCut(Float_t CutOnRelativeDistance) +void AliPHOSPIDv1::PlotDispersionCuts()const { - // Modifies the parameters used for the particle type identification + TCanvas* lambdas = new TCanvas("lambdas","Cuts on the elipse axise",200,10,700,500); + + if(fIDOptions.Contains("ell",TString::kIgnoreCase ) ){ + TF2 * ell = new TF2("Elliptic Cuts",fFormula->GetName(),0,3,0,3) ; + ell->SetMinimum(0.0000001) ; + ell->SetMaximum(0.001) ; + ell->SetLineStyle(1) ; + ell->SetLineWidth(2) ; + ell->Draw() ; + } + + if( fIDOptions.Contains("dis",TString::kIgnoreCase ) ){ + TF2 * dsp = new TF2("dispersion","(ySetParameter(0,fDispersion) ; + dsp->SetMinimum(0.0000001) ; + dsp->SetMaximum(0.001) ; + dsp->SetLineStyle(1) ; + dsp->SetLineColor(2) ; + dsp->SetLineWidth(2) ; + dsp->SetNpx(200) ; + dsp->SetNpy(200) ; + if(fIDOptions.Contains("ell",TString::kIgnoreCase ) ) + dsp->Draw("same") ; + else + dsp->Draw() ; + } + lambdas->Update(); +} - fCutOnRelativeDistance = CutOnRelativeDistance ; +//____________________________________________________________________________ +TVector3 AliPHOSPIDv1::GetMomentumDirection(AliPHOSEmcRecPoint * emc, AliPHOSRecPoint * cpv,AliPHOSRecPoint * ppsd)const +{ + // Calculates the momentum direction: + // 1. if only a EMC RecPoint, direction is given by IP and this RecPoint + // 2. if a EMC RecPoint and one PPSD RecPoint, direction is given by the line through the 2 recpoints + // 3. if a EMC RecPoint and two PPSD RecPoints, dirrection is given by the average line through + // the 2 pairs of recpoints + // However because of the poor position resolution of PPSD the direction is always taken as if we were + // in case 1. + + TVector3 dir(0,0,0) ; + + TVector3 emcglobalpos ; + TMatrix dummy ; + + emc->GetGlobalPosition(emcglobalpos, dummy) ; + + // The following commeneted code becomes valid once the PPSD provides + // a reasonable position resolution, at least as good as EMC ! + // TVector3 ppsdlglobalpos ; + // TVector3 ppsduglobalpos ; + // if( fPpsdLowRecPoint ){ // certainly a photon that has concerted + // fPpsdLowRecPoint->GetGlobalPosition(ppsdlglobalpos, mdummy) ; + // dir = emcglobalpos - ppsdlglobalpos ; + // if( fPpsdUpRecPoint ){ // not looks like a charged + // fPpsdUpRecPoint->GetGlobalPosition(ppsduglobalpos, mdummy) ; + // dir = ( dir + emcglobalpos - ppsduglobalpos ) * 0.5 ; + // } + // } + // else { // looks like a neutral + // dir = emcglobalpos ; + // } + + dir = emcglobalpos ; + dir.SetZ( -dir.Z() ) ; // why ? + dir.SetMag(1.) ; + + //account correction to the position of IP + Float_t xo,yo,zo ; //Coordinates of the origin + gAlice->Generator()->GetOrigin(xo,yo,zo) ; + TVector3 origin(xo,yo,zo); + dir = dir - origin ; + + return dir ; +} +//____________________________________________________________________________ +void AliPHOSPIDv1::PrintRecParticles(Option_t * option){ + + cout << "AliPHOSPIDv1: " << endl ; + cout << " found " << fRecParticles->GetEntriesFast() << " RecParticles " << endl ; + + if(strstr(option,"all")) { // printing found TS + + cout << " PARTICLE " + << " Index " << endl ; + // << " X " + // << " Y " + // << " Z " + // << " # of primaries " + // << " Primaries list " << endl; + + Int_t index ; + for (index = 0 ; index < fRecParticles->GetEntries() ; index++) { + AliPHOSRecParticle * rp = (AliPHOSRecParticle * ) fRecParticles->At(index) ; + + Text_t particle[11]; + switch(rp->GetType()) { + case AliPHOSFastRecParticle::kNEUTRALEM: + strcpy( particle, "NEUTRAL_EM"); + break; + case AliPHOSFastRecParticle::kNEUTRALHA: + strcpy(particle, "NEUTRAL_HA"); + break; + case AliPHOSFastRecParticle::kGAMMA: + strcpy(particle, "GAMMA"); + break ; + case AliPHOSFastRecParticle::kGAMMAHA: + strcpy(particle, "GAMMA_H"); + break ; + case AliPHOSFastRecParticle::kABSURDEM: + strcpy(particle, "ABSURD_EM") ; + break ; + case AliPHOSFastRecParticle::kABSURDHA: + strcpy(particle, "ABSURD_HA") ; + break ; + case AliPHOSFastRecParticle::kELECTRON: + strcpy(particle, "ELECTRON") ; + break ; + case AliPHOSFastRecParticle::kCHARGEDHA: + strcpy(particle, "CHARGED_HA") ; + break ; + } + + // Int_t * primaries; + // Int_t nprimaries; + // primaries = rp->GetPrimaries(nprimaries); + + cout << setw(15) << particle << " " + << setw(3) << rp->GetIndexInList() << " " ; + // << setw(4) << nprimaries << " "; + // for (Int_t iprimary=0; iprimaryExecuteTask() +// +// // One can specify the title for each branch +// root [2] r->SetBranchFileName("RecPoints","RecPoints1") ; +// // By default branches are stored in galice.root (in non-split mode) +// // or PHOS.SDigits.root, PHOS.Digits.root etc. +// +// // 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" // --- Standard library --- - -#include +#include // --- AliRoot header files --- - +#include "AliRun.h" #include "AliPHOSReconstructioner.h" -#include "AliPHOSClusterizer.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(AliPHOSClusterizer * Clusterizer, - AliPHOSTrackSegmentMaker * Tracker, - AliPHOSPID * Pid) + AliPHOSReconstructioner::AliPHOSReconstructioner():TTask("AliPHOSReconstructioner","") { // ctor - - fClusterizer = Clusterizer ; - fTrackSegmentMaker = Tracker ; - fPID = Pid ; - fDebugReconstruction = kFALSE ; -} + fDigitizer = 0 ; + fClusterizer = 0 ; + fTSMaker = 0 ; + fPID = 0 ; + fSDigitizer = 0 ; + fIsInitialized = kFALSE ; -//____________________________________________________________________________ - void AliPHOSReconstructioner::Init(AliPHOSClusterizer * Clusterizer, - AliPHOSTrackSegmentMaker * Tracker, - AliPHOSPID * Pid) -{ - // Initialisation - - fClusterizer = Clusterizer ; - fTrackSegmentMaker = Tracker ; - fPID = Pid ; - fDebugReconstruction = kFALSE ; } //____________________________________________________________________________ - void AliPHOSReconstructioner::Make(DigitsList * dl, - AliPHOSRecPoint::RecPointsList * emccl, - AliPHOSRecPoint::RecPointsList * ppsdl, - AliPHOSTrackSegment::TrackSegmentsList * trsl, - AliPHOSRecParticle::RecParticlesList * rpl) +AliPHOSReconstructioner::AliPHOSReconstructioner(const char* headerFile):TTask("AliPHOSReconstructioner","") { - // Launches the Reconstruction process in the sequence: Make the reconstructed poins (clusterize) - // Make the track segments - // Make the reconstructed particles - Int_t index ; - if (fDebugReconstruction) - cout << "\n\nDebugReconstruction>>> " << "Start making reconstructed points (clusterizing!!)" << endl; + // ctor + + fHeaderFileName = headerFile ; + + fSDigitsBranch="" ; + fSDigitizer = new AliPHOSSDigitizer(fHeaderFileName.Data(),fSDigitsBranch.Data()) ; + Add(fSDigitizer) ; + + fDigitsBranch="" ; + fDigitizer = new AliPHOSDigitizer(fHeaderFileName.Data(),fDigitsBranch.Data()) ; + Add(fDigitizer) ; + + + fRecPointBranch="" ; + fClusterizer = new AliPHOSClusterizerv1(fHeaderFileName.Data(),fRecPointBranch.Data()) ; + Add(fClusterizer) ; - fClusterizer->MakeClusters(dl, emccl, ppsdl); - if (fDebugReconstruction){ - cout << "DebugReconstruction>>> " << "AliPHOSReconstructioner: Digit list entries is " << dl->GetEntries() << endl ; - cout << "AliPHOSReconstructioner: Emc list entries is " << emccl->GetEntries() << endl ; - cout << "AliPHOSReconstructioner: Ppsd list entries is " << ppsdl->GetEntries() << endl ; + fTSBranch="" ; + fTSMaker = new AliPHOSTrackSegmentMakerv1(fHeaderFileName.Data(),fTSBranch.Data()) ; + Add(fTSMaker) ; + + + fRecPartBranch="" ; + fPID = new AliPHOSPIDv1(fHeaderFileName.Data(),fRecPartBranch.Data()) ; + Add(fPID) ; + + fIsInitialized = kTRUE ; + +} +//____________________________________________________________________________ +void AliPHOSReconstructioner::Exec(Option_t *option){ + //chesk, if the names of branches, which should be made conicide with already + //existing + if(!fIsInitialized) + Init() ; + + gAlice->GetEvent(0) ; + + if(fSDigitizer->IsActive()&& gAlice->TreeS()){ //Will produce SDigits + + TBranch * sdigitsBranch = 0; + TBranch * sdigitizerBranch = 0; + + TObjArray * branches = gAlice->TreeS()->GetListOfBranches() ; + Int_t ibranch; + Bool_t phosNotFound = kTRUE ; + Bool_t sdigitizerNotFound = kTRUE ; + + for(ibranch = 0;ibranch GetEntries();ibranch++){ + if(phosNotFound){ + sdigitsBranch=(TBranch *) branches->At(ibranch) ; + if(( strcmp("PHOS",sdigitsBranch->GetName())==0 ) && + (fSDigitsBranch.CompareTo(sdigitsBranch->GetTitle())== 0 )) + phosNotFound = kFALSE ; + } + if(sdigitizerNotFound){ + sdigitizerBranch = (TBranch *) branches->At(ibranch) ; + if(( strcmp(sdigitizerBranch->GetName(),"AliPHOSSDigitizer") == 0) && + (fSDigitsBranch.CompareTo(sdigitizerBranch->GetTitle())== 0 ) ) + sdigitizerNotFound = kFALSE ; + } + } + + if(!(sdigitizerNotFound && phosNotFound)){ + cout << "AliPHOSReconstructioner error: "<< endl ; + cout << " Branches ''PHOS'' or ''AliPHOSSDigitizer'' with title ``" << fSDigitsBranch.Data() << "''" << endl ; + cout << " already exist in TreeS. ROOT does not allow updating/overwriting." << endl ; + cout << " Specify another title for branches or use ''StartFrom()'' method" << endl ; + + //mark all tasks as inactive + TIter next(fTasks); + TTask *task; + while((task=(TTask*)next())) + task->SetActive(kFALSE) ; + + return ; + } } - // Digit Debuging - if (fDebugReconstruction) { - cout << ">>>>>>>>>>>>>>>>>>>>>> DebugReconstruction <<<<<<<<<<<<<<<<<<<<<<<<<<" << endl ; - cout << "DebugReconstruction>>> Digit list entries is " << dl->GetEntries() << endl ; - AliPHOSDigit * digit; - Bool_t calorimeter ; - Float_t factor; - cout << "DebugReconstruction>>> Vol Id " << - " Ene(MeV, KeV) " << - " Index " << - " Nprim " << - " Primaries list " << endl; - for (index = 0 ; index < dl->GetEntries() ; index++) { - digit = (AliPHOSDigit * ) dl->At(index) ; - calorimeter = fClusterizer->IsInEmc(digit); - if (calorimeter) factor =1000. ; else factor=1000000.; - cout << "DebugReconstruction>>> " << - setw(8) << digit->GetId() << " " << - setw(3) << (Int_t) calorimeter << - setw(10) << factor*fClusterizer->Calibrate(digit->GetAmp()) << " " << - setw(6) << digit->GetIndexInList() << " " << - setw(5) << digit->GetNprimary() <<" "; - for (Int_t iprimary=0; iprimaryGetNprimary(); iprimary++) - cout << setw(5) << digit->GetPrimary(iprimary+1) << " "; - cout << endl; + if(fDigitizer->IsActive() && gAlice->TreeD()){ //Will produce Digits + TBranch * digitsBranch = 0; + TBranch * digitizerBranch = 0; + + TObjArray * branches = gAlice->TreeD()->GetListOfBranches() ; + Int_t ibranch; + Bool_t phosNotFound = kTRUE ; + Bool_t digitizerNotFound = kTRUE ; + + for(ibranch = 0;ibranch GetEntries();ibranch++){ + if(phosNotFound){ + digitsBranch=(TBranch *) branches->At(ibranch) ; + if(( strcmp("PHOS",digitsBranch->GetName())==0 ) && + (fDigitsBranch.CompareTo(digitsBranch->GetTitle())== 0 )) + phosNotFound = kFALSE ; + } + if(digitizerNotFound){ + digitizerBranch = (TBranch *) branches->At(ibranch) ; + if(( strcmp(digitizerBranch->GetName(),"AliPHOSDigitizer") == 0) && + (fDigitsBranch.CompareTo(digitizerBranch->GetTitle())== 0 ) ) + digitizerNotFound = kFALSE ; + } } + if(!(digitizerNotFound && phosNotFound)){ + cout << "AliPHOSReconstructioner error: "<< endl ; + cout << " Branches ''PHOS'' or ''AliPHOSDigitizer'' with title ``" << fDigitsBranch.Data() << "''" << endl ; + cout << " already exist in TreeD. ROOT does not allow updating/overwriting." << endl ; + cout << " Specify another title for branches or use ''StartFrom()'' method" << endl ; + + //mark all tasks as inactive + TIter next(fTasks); + TTask *task; + while((task=(TTask*)next())) + task->SetActive(kFALSE) ; + + return ; + } } - // Making Clusters - if (fDebugReconstruction) cout << "DebugReconstruction>>> Start making reconstructed points (clusterizing)" << endl; + if(fClusterizer->IsActive() && gAlice->TreeR()){ //Will produce RecPoints + TBranch * emcBranch = 0; + TBranch * cpvBranch = 0; + TBranch * clusterizerBranch = 0; + + TObjArray * branches = gAlice->TreeR()->GetListOfBranches() ; + Int_t ibranch; + Bool_t emcNotFound = kTRUE ; + Bool_t cpvNotFound = kTRUE ; + Bool_t clusterizerNotFound = kTRUE ; + + for(ibranch = 0;ibranch GetEntries();ibranch++){ + + if(emcNotFound){ + emcBranch=(TBranch *) branches->At(ibranch) ; + if(fRecPointBranch.CompareTo(emcBranch->GetTitle())==0 ) + if( strcmp(emcBranch->GetName(),"PHOSEmcRP") == 0) + emcNotFound = kFALSE ; + } + if(cpvNotFound){ + cpvBranch=(TBranch *) branches->At(ibranch) ; + if(fRecPointBranch.CompareTo(cpvBranch->GetTitle())==0 ) + if( strcmp(cpvBranch->GetName(),"PHOSCpvRP") == 0) + cpvNotFound = kFALSE ; + } + if(clusterizerNotFound){ + clusterizerBranch = (TBranch *) branches->At(ibranch) ; + if( fRecPointBranch.CompareTo(clusterizerBranch->GetTitle()) == 0) + if( strcmp(clusterizerBranch->GetName(),"AliPHOSClusterizer") == 0) + clusterizerNotFound = kFALSE ; + } + } - // mark the position of the RecPoints in the array - AliPHOSEmcRecPoint * emcrp ; - for (index = 0 ; index < emccl->GetEntries() ; index++) { - emcrp = (AliPHOSEmcRecPoint * )emccl->At(index) ; - emcrp->SetIndexInList(index) ; - } - AliPHOSPpsdRecPoint * ppsdrp ; - for (index = 0 ; index < ppsdl->GetEntries() ; index++) { - ppsdrp = (AliPHOSPpsdRecPoint * )ppsdl->At(index) ; - ppsdrp->SetIndexInList(index) ; + if(!(clusterizerNotFound && emcNotFound && cpvNotFound)){ + cout << "AliPHOSReconstructioner error: "<< endl ; + cout << " Branches ''PHOSEmcRP'', ''PHOSCpvRP'' or ''AliPHOSClusterizer'' with title ``" + << fRecPointBranch.Data() << "''" << endl ; + cout << " already exist in TreeR. ROOT does not allow updating/overwriting." << endl ; + cout << " Specify another title for branches or use ''StartFrom()'' method" << endl ; + + //mark all tasks as inactive + TIter next(fTasks); + TTask *task; + while((task=(TTask*)next())) + task->SetActive(kFALSE) ; + return ; + } } - if (fDebugReconstruction) { - cout << "DebugReconstruction>>> Cluster emc list entries is " << emccl->GetEntries() << endl ; - AliPHOSEmcRecPoint * recpoint; - cout << "DebugReconstruction>>> Module " << - "Ene(MeV) " << - "Index " << - "Multi " << - " X " << - " Y " << - " Z " << - " Lambda 1 " << - " Lambda 2 " << - "MaxEnergy(MeV) " << - "Nprim " << - " Primaries list " << endl; - for (index = 0 ; index < emccl->GetEntries() ; index++) { - recpoint = (AliPHOSEmcRecPoint * )emccl->At(index) ; - TVector3 locpos; recpoint->GetLocalPosition(locpos); - Float_t lambda[2]; recpoint->GetElipsAxis(lambda); - Int_t * primaries; - Int_t nprimaries; - primaries = recpoint->GetPrimaries(nprimaries); - cout << "DebugReconstruction>>> " << - setw(2) <GetPHOSMod() << " " << - setw(9) << 1000.*recpoint->GetEnergy() << " " << - setw(6) << recpoint->GetIndexInList() << " " << - setw(5) << recpoint->GetMultiplicity() <<" " << - setw(8) << locpos.X() <<" " << - setw(8) << locpos.Y() <<" " << - setw(8) << locpos.Z() << " " << - setw(10) << lambda[0] << " " << - setw(10) << lambda[1] << " " << - setw(9) << 1000*recpoint->GetMaximalEnergy() << " " << - setw(9) << nprimaries << " "; - for (Int_t iprimary=0; iprimaryIsActive() && gAlice->TreeR()){ //Produce TrackSegments + + TBranch * tsMakerBranch = 0; + TBranch * tsBranch = 0; + + TObjArray * branches = gAlice->TreeR()->GetListOfBranches() ; + Int_t ibranch; + Bool_t tsMakerNotFound = kTRUE ; + Bool_t tsNotFound = kTRUE ; + + for(ibranch = 0;(ibranch GetEntries())&&(tsMakerNotFound||tsNotFound);ibranch++){ + if(tsMakerNotFound){ + tsMakerBranch=(TBranch *) branches->At(ibranch) ; + if( fTSBranch.CompareTo(tsMakerBranch->GetTitle())==0 ) + if( strcmp(tsMakerBranch->GetName(),"AliPHOSTrackSegmentMaker") == 0) + tsMakerNotFound = kFALSE ; + } + if(tsNotFound){ + tsBranch=(TBranch *) branches->At(ibranch) ; + if( fTSBranch.CompareTo(tsBranch->GetTitle())==0 ) + if( strcmp(tsBranch->GetName(),"PHOSTS") == 0) + tsNotFound = kFALSE ; + } } + + if(!(tsMakerNotFound &&tsNotFound) ){ + cout << "AliPHOSReconstructioner error: "<< endl ; + cout << " Branches ''PHOSTS'' or ''AliPHOSTrackSegmentMaker'' with title ``" + << fTSBranch.Data() << "''" << endl ; + cout << " already exist in TreeR. ROOT does not allow updating/overwriting." << endl ; + cout << " Specify another title for branches or use ''StartFrom()'' method" << endl ; + + //mark all tasks as inactive + TIter next(fTasks); + TTask *task; + while((task=(TTask*)next())) + task->SetActive(kFALSE) ; + return ; - cout << "DebugReconstruction>>> Cluster ppsd list entries is " << ppsdl->GetEntries() << endl ; - AliPHOSPpsdRecPoint * ppsdrecpoint; - Text_t detector[4]; - cout << "DebugReconstruction>>> Module " << - "Det " << - "Ene(KeV) " << - "Index " << - "Multi " << - " X " << - " Y " << - " Z " << - "Nprim " << - " Primaries list " << endl; - for (index = 0 ; index < ppsdl->GetEntries() ; index++) { - ppsdrecpoint = (AliPHOSPpsdRecPoint * ) ppsdl->At(index) ; - TVector3 locpos; ppsdrecpoint->GetLocalPosition(locpos); - Int_t * primaries; - Int_t nprimaries; - if (ppsdrecpoint->GetUp()) - strcpy(detector, "CPV"); - else - strcpy(detector, "PC "); - primaries = ppsdrecpoint->GetPrimaries(nprimaries); - cout << "DebugReconstruction>>> " << - setw(4) << ppsdrecpoint->GetPHOSMod() << " " << - setw(4) << detector << " " << - setw(9) << 1000000.*ppsdrecpoint->GetEnergy() << " " << - setw(6) << ppsdrecpoint->GetIndexInList() << " " << - setw(5) << ppsdrecpoint->GetMultiplicity() <<" " << - setw(8) << locpos.X() <<" " << - setw(8) << locpos.Y() <<" " << - setw(8) << locpos.Z() << " " << - setw(9) << nprimaries << " "; - for (Int_t iprimary=0; iprimary>>> Start making track segments(unfolding+tracksegments)" << endl; - fTrackSegmentMaker->MakeTrackSegments(dl, emccl, ppsdl, trsl) ; - - // mark the position of the TrackSegments in the array - AliPHOSTrackSegment * trs ; - for (index = 0 ; index < trsl->GetEntries() ; index++) { - trs = (AliPHOSTrackSegment * )trsl->At(index) ; - trs->SetIndexInList(index) ; + } - if (fDebugReconstruction){ - cout << "DebugReconstruction>>> Track segment list entries is " << trsl->GetEntries() << endl ; - cout << "DebugReconstruction>>> Module " << - "Ene(KeV) " << - "Index " << - " X " << - " Y " << - " Z " << - " rX " << - " rY " << - " rZ " << - "Nprim " << - " Primaries list " << endl; + + if(fPID->IsActive() && gAlice->TreeR()){ //Produce RecParticles + TBranch * pidBranch = 0; + TBranch * rpBranch = 0; - for (index = 0 ; index < trsl->GetEntries() ; index++) { - trs = (AliPHOSTrackSegment * )trsl->At(index) ; - TVector3 locpos; trs->GetPosition(locpos); - Int_t * primaries; - Int_t nprimaries; - primaries = trs->GetPrimariesEmc(nprimaries); - cout << "DebugReconstruction>>> " << - setw(4) << trs->GetPHOSMod() << " " << - setw(9) << 1000.*trs->GetEnergy() << " " << - setw(3) << trs->GetIndexInList() << " " << - setw(9) << locpos.X() <<" " << - setw(9) << locpos.Y() <<" " << - setw(9) << locpos.Z() << " " << - setw(10) << (trs->GetMomentumDirection()).X() << " " << - setw(10) << (trs->GetMomentumDirection()).Y() << " " << - setw(10) << (trs->GetMomentumDirection()).Z() << " " << - setw(4) << nprimaries << " "; - for (Int_t iprimary=0; iprimaryTreeR()->GetListOfBranches() ; + Int_t ibranch; + Bool_t pidNotFound = kTRUE ; + Bool_t rpNotFound = kTRUE ; + + for(ibranch = 0;(ibranch GetEntries()) && pidNotFound && rpNotFound ;ibranch++){ + if(pidNotFound){ + pidBranch=(TBranch *) branches->At(ibranch) ; + if( (strcmp(fRecPartBranch,pidBranch->GetTitle())==0 ) && + (strcmp(pidBranch->GetName(),"AliPHOSPID") == 0) ) + pidNotFound = kFALSE ; + } + if(rpNotFound){ + rpBranch=(TBranch *) branches->At(ibranch) ; + if( (strcmp(fRecPartBranch,rpBranch->GetTitle())==0 ) && + (strcmp(rpBranch->GetName(),"PHOSRP") == 0) ) + rpNotFound = kFALSE ; + } + } + + if(!pidNotFound || !rpNotFound ){ + cout << "AliPHOSReconstructioner error: "<< endl ; + cout << " Branches ''PHOSRP'' or ''AliPHOSPID'' with title ``" + << fRecPartBranch.Data() << "''" << endl ; + cout << " already exist in TreeR. ROOT does not allow updating/overwriting." << endl ; + cout << " Specify another title for branches." << endl ; + + //mark all tasks as inactive + TIter next(fTasks); + TTask *task; + while((task=(TTask*)next())) + task->SetActive(kFALSE) ; + return ; } } - if (fDebugReconstruction) cout << "DebugReconstruction>>>> Start making reconstructed particles" << endl; +} +//____________________________________________________________________________ + void AliPHOSReconstructioner::Init() +{ + //initiase Reconstructioner if necessary: we can not do this in default constructor + + if(!fIsInitialized){ + // Initialisation + + fSDigitsBranch="" ; + fSDigitizer = new AliPHOSSDigitizer(fHeaderFileName.Data(),fSDigitsBranch.Data()) ; + Add(fSDigitizer) ; + + fDigitsBranch="" ; + fDigitizer = new AliPHOSDigitizer(fHeaderFileName.Data(),fDigitsBranch.Data()) ; + Add(fDigitizer) ; + + fRecPointBranch="" ; + fClusterizer = new AliPHOSClusterizerv1(fHeaderFileName.Data(),fRecPointBranch.Data()) ; + Add(fClusterizer) ; + + fTSBranch="" ; + fTSMaker = new AliPHOSTrackSegmentMakerv1(fHeaderFileName.Data(),fTSBranch.Data()) ; + Add(fTSMaker) ; + + + fRecPartBranch="" ; + fPID = new AliPHOSPIDv1(fHeaderFileName.Data(),fRecPartBranch.Data()) ; + Add(fPID) ; + + fIsInitialized = kTRUE ; + } +} +//____________________________________________________________________________ +AliPHOSReconstructioner::~AliPHOSReconstructioner() +{ - if (fPID) { - fPID->MakeParticles(trsl, rpl) ; + if(fSDigitizer) + delete fSDigitizer ; - // mark the position of the RecParticles in the array - AliPHOSRecParticle * rp ; - for (index = 0 ; index < rpl->GetEntries() ; index++) { - rp = (AliPHOSRecParticle * )rpl->At(index) ; - rp->SetIndexInList(index) ; - } - //Debugger of RecParticles - if (fDebugReconstruction){ - cout << "DebugReconstruction>>> Reconstructed particle list entries is " << rpl->GetEntries() << endl ; - cout << "DebugReconstruction>>> Module " << - " PARTICLE " << - "Ene(KeV) " << - "Index " << - " X " << - " Y " << - " Z " << - "Nprim " << - " Primaries list " << endl; - for (index = 0 ; index < rpl->GetEntries() ; index++) { - rp = (AliPHOSRecParticle * ) rpl->At(index) ; - TVector3 locpos; (rp->GetPHOSTrackSegment())->GetPosition(locpos); - Int_t * primaries; - Int_t nprimaries; - Text_t particle[11]; - primaries = (rp->GetPHOSTrackSegment())->GetPrimariesEmc(nprimaries); - switch(rp->GetType()) { - case AliPHOSFastRecParticle::kNEUTRALEM: - strcpy( particle, "NEUTRAL_EM"); - break; - case AliPHOSFastRecParticle::kNEUTRALHA: - strcpy(particle, "NEUTRAL_HA"); - break; - case AliPHOSFastRecParticle::kGAMMA: - strcpy(particle, "GAMMA"); - break ; - case AliPHOSFastRecParticle::kGAMMAHA: - strcpy(particle, "GAMMA_H"); - break ; - case AliPHOSFastRecParticle::kABSURDEM: - strcpy(particle, "ABSURD_EM") ; - break ; - case AliPHOSFastRecParticle::kABSURDHA: - strcpy(particle, "ABSURD_HA") ; - break ; - case AliPHOSFastRecParticle::kELECTRON: - strcpy(particle, "ELECTRON") ; - break ; - case AliPHOSFastRecParticle::kCHARGEDHA: - strcpy(particle, "CHARGED_HA") ; - break ; - } - - cout << "DebugReconstruction>>> " << - setw(4) << (rp->GetPHOSTrackSegment())->GetPHOSMod() << " " << - setw(15) << particle << " " << - setw(9) << 1000.*(rp->GetPHOSTrackSegment())->GetEnergy() << " " << - setw(3) << rp->GetIndexInList() << " " << - setw(9) << locpos.X() <<" " << - setw(9) << locpos.Y() <<" " << - setw(9) << locpos.Z() << " " << - setw(4) << nprimaries << " "; - for (Int_t iprimary=0; iprimarySetSDigitsBranch(title) ; + fDigitizer->SetSDigitsBranch(title) ; + fSDigitsBranch = title ; + return ; } + + if(strcmp(branch,"Digits") == 0){ + fDigitizer->SetDigitsBranch(title) ; + fClusterizer->SetDigitsBranch(title) ; + fDigitsBranch = title ; + return ; + } + + if(strcmp(branch,"RecPoints") == 0){ + fClusterizer->SetRecPointsBranch(title) ; + fTSMaker->SetRecPointsBranch(title) ; + fRecPointBranch = title ; + return ; + } + + if(strcmp(branch,"TrackSegments") == 0){ + fTSMaker->SetTrackSegmentsBranch(title) ; + fPID->SetTrackSegmentsBranch(title) ; + fTSBranch = title ; + return ; + } + + if(strcmp(branch,"RecParticles") == 0){ + fPID->SetRecParticlesBranch(title) ; + fRecPartBranch = title ; + 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 " <GetName() << endl ; + } +} + +//____________________________________________________________________________ +void AliPHOSReconstructioner::Print(Option_t * option)const { + + cout << "-----------------AliPHOSReconstructioner---------------" << endl ; + cout << " Reconstruction of the header file " <IsActive()){ + cout << " (+) " << fSDigitizer->GetName() << " to branch : " << fSDigitsBranch.Data() << endl ; + cout << endl ; + } + if(fDigitizer->IsActive()){ + cout << " (+) " << fDigitizer->GetName() << " to branch : " << fDigitsBranch.Data() << endl ; + cout << endl ; + } + + if(fClusterizer->IsActive()){ + cout << " (+) " <GetName() << " to branch : " <IsActive()){ + cout << " (+) " << fTSMaker->GetName() << " to branch : " << fTSBranch.Data() << endl ; + cout << endl ; + } + + + if(fPID->IsActive()){ + cout << " (+) " << fPID->GetName() << " to branch : " < // --- AliRoot header files --- - +#include "AliRun.h" #include "AliPHOSDigit.h" #include "AliPHOSHit.h" -#include "AliPHOSv1.h" #include "AliPHOSSDigitizer.h" -#include "TROOT.h" -#include "TFolder.h" ClassImp(AliPHOSSDigitizer) @@ -49,58 +59,88 @@ ClassImp(AliPHOSSDigitizer) fA = 0; fB = 10000000. ; fPrimThreshold = 0.01 ; - fNevents = 0 ; // Number of events to digitize, 0 means all evens in current file - // add Task to //root/Tasks folder - TTask * roottasks = (TTask*)gROOT->GetRootFolder()->FindObject("Tasks") ; - roottasks->Add(this) ; + fNevents = 0 ; + fSDigits = 0 ; + fHits = 0 ; + fIsInitialized = kFALSE ; } + //____________________________________________________________________________ -AliPHOSSDigitizer::AliPHOSSDigitizer(char* HeaderFile, char *SDigitsFile):TTask("AliPHOSSDigitizer","") +AliPHOSSDigitizer::AliPHOSSDigitizer(const char* headerFile, const char *sDigitsTitle):TTask("AliPHOSSDigitizer","") { // ctor fA = 0; fB = 10000000.; fPrimThreshold = 0.01 ; - fNevents = 0 ; // Number of events to digitize, 0 means all events in current file - fSDigitsFile = SDigitsFile ; - fHeadersFile = HeaderFile ; + fNevents = 0 ; + fSDigitsTitle = sDigitsTitle ; + fHeadersFile = headerFile ; + fSDigits = new TClonesArray("AliPHOSDigit",1000); + fHits = new TClonesArray("AliPHOSHit",1000); + + TFile * file = (TFile*) gROOT->GetFile(fHeadersFile.Data() ) ; + + //File was not opened yet + if(file == 0){ + file = new TFile(fHeadersFile.Data(),"update") ; + gAlice = (AliRun *) file->Get("gAlice") ; + } + //add Task to //root/Tasks folder TTask * roottasks = (TTask*)gROOT->GetRootFolder()->FindObject("Tasks") ; roottasks->Add(this) ; - + + fIsInitialized = kTRUE ; } //____________________________________________________________________________ - AliPHOSSDigitizer::~AliPHOSSDigitizer() +AliPHOSSDigitizer::~AliPHOSSDigitizer() { // dtor + if(fSDigits) + delete fSDigits ; + if(fHits) + delete fHits ; } +//____________________________________________________________________________ +void AliPHOSSDigitizer::Init(){ + //Initialization can not be done in the default constructor + + if(!fIsInitialized){ + if(fHeadersFile.IsNull()) + fHeadersFile="galice.root" ; + TFile * file = (TFile*) gROOT->GetFile(fHeadersFile.Data() ) ; + + //if file was not opened yet, read gAlice + if(file == 0){ + file = new TFile(fHeadersFile.Data(),"update") ; + gAlice = (AliRun *) file->Get("gAlice") ; + } + + fHits = new TClonesArray("AliPHOSHit",1000); + fSDigits = new TClonesArray("AliPHOSDigit",1000); + + // add Task to //root/Tasks folder + TTask * roottasks = (TTask*)gROOT->GetRootFolder()->FindObject("Tasks") ; + roottasks->Add(this) ; + + fIsInitialized = kTRUE ; + } +} //____________________________________________________________________________ void AliPHOSSDigitizer::Exec(Option_t *option) { //Collects all hits in the same active volume into digit - - TFile * file = 0; - - // if(gAlice->TreeE()==0) //If gAlice not yet read/constructed - if(fHeadersFile.IsNull()) - file = new TFile("galice.root", "update") ; - else - file = new TFile(fHeadersFile.Data(),"update") ; - - gAlice = (AliRun *) file->Get("gAlice") ; - - - TClonesArray * sdigits = new TClonesArray("AliPHOSDigit",1000) ; + if(!fIsInitialized) + Init() ; - AliPHOS * PHOS = (AliPHOS *) gAlice->GetDetector("PHOS") ; + if(strstr(option,"tim")) + gBenchmark->Start("PHOSSDigitizer"); - - if(fNevents == 0) - fNevents = (Int_t) gAlice->TreeE()->GetEntries() ; + fNevents = (Int_t) gAlice->TreeE()->GetEntries() ; Int_t ievent ; for(ievent = 0; ievent < fNevents; ievent++){ @@ -112,42 +152,22 @@ void AliPHOSSDigitizer::Exec(Option_t *option) { return ; } - if(gAlice->TreeS() == 0) - gAlice->MakeTree("S") ; - - TClonesArray * hits = PHOS->Hits() ; + //set address of the hits + TBranch * branch = gAlice->TreeH()->GetBranch("PHOS"); + if (branch) + branch->SetAddress(&fHits); + else{ + cout << "ERROR in AliPHOSSDigitizer: "<< endl ; + cout << " no branch PHOS in TreeH"<< endl ; + cout << " do nothing " << endl ; + return ; + } - sdigits->Clear(); + fSDigits->Clear(); Int_t nSdigits = 0 ; - //Make branches - char branchname[20]; - sprintf(branchname,"%s",PHOS->GetName()); - - Int_t bufferSize = 16000 ; - char * file =0; - if(!fSDigitsFile.IsNull()) - file = (char*) fSDigitsFile.Data() ; //ievent ; - else - if(gSystem->Getenv("CONFIG_SPLIT_FILE")){ //generating file name - file = new char[30] ; - // sprintf(file,"PHOS.SDigits%d.root",ievent) ; - sprintf(file,"PHOS.SDigits.root") ; - } - else - file = 0 ; - - gAlice->MakeBranchInTree(gAlice->TreeS(),branchname,&sdigits,bufferSize,file); - - - Int_t splitlevel = 0 ; - sprintf(branchname,"AliPHOSSDigitizer"); - AliPHOSSDigitizer * sd = this ; - gAlice->MakeBranchInTree(gAlice->TreeS(),branchname,"AliPHOSSDigitizer",&sd, bufferSize, splitlevel,file); - - //Now made SDigits from hits, for PHOS it is the same - + //Now made SDigits from hits, for PHOS it is the same, so just copy Int_t itrack ; for (itrack=0; itrackGetNtrack(); itrack++){ @@ -156,8 +176,8 @@ void AliPHOSSDigitizer::Exec(Option_t *option) { gAlice->TreeH()->GetEvent(itrack); Int_t i; - for ( i = 0 ; i < hits->GetEntries() ; i++ ) { - AliPHOSHit * hit = (AliPHOSHit*)hits->At(i) ; + for ( i = 0 ; i < fHits->GetEntries() ; i++ ) { + AliPHOSHit * hit = (AliPHOSHit*)fHits->At(i) ; AliPHOSDigit * newdigit ; // Assign primary number only if contribution is significant @@ -166,7 +186,7 @@ void AliPHOSSDigitizer::Exec(Option_t *option) { else newdigit = new AliPHOSDigit( -1 , hit->GetId(), Digitize( hit->GetEnergy() ) ) ; - new((*sdigits)[nSdigits]) AliPHOSDigit(* newdigit) ; + new((*fSDigits)[nSdigits]) AliPHOSDigit(* newdigit) ; nSdigits++ ; delete newdigit ; @@ -174,44 +194,123 @@ void AliPHOSSDigitizer::Exec(Option_t *option) { } // loop over tracks - sdigits->Sort() ; + fSDigits->Sort() ; - nSdigits = sdigits->GetEntries() ; - sdigits->Expand(nSdigits) ; + nSdigits = fSDigits->GetEntriesFast() ; + fSDigits->Expand(nSdigits) ; Int_t i ; for (i = 0 ; i < nSdigits ; i++) { - AliPHOSDigit * digit = (AliPHOSDigit *) sdigits->At(i) ; + AliPHOSDigit * digit = (AliPHOSDigit *) fSDigits->At(i) ; digit->SetIndexInList(i) ; } + + if(gAlice->TreeS() == 0) + gAlice->MakeTree("S") ; + + //check, if this branch already exits? + TBranch * sdigitsBranch = 0; + TBranch * sdigitizerBranch = 0; + + TObjArray * branches = gAlice->TreeS()->GetListOfBranches() ; + Int_t ibranch; + Bool_t phosNotFound = kTRUE ; + Bool_t sdigitizerNotFound = kTRUE ; + + for(ibranch = 0;ibranch GetEntries();ibranch++){ + + if(phosNotFound){ + sdigitsBranch=(TBranch *) branches->At(ibranch) ; + if( (strcmp("PHOS",sdigitsBranch->GetName())==0 ) && + (fSDigitsTitle.CompareTo(sdigitsBranch->GetTitle()) == 0) ) + phosNotFound = kFALSE ; + } + if(sdigitizerNotFound){ + sdigitizerBranch = (TBranch *) branches->At(ibranch) ; + if( (strcmp(sdigitizerBranch->GetName(),"AliPHOSSDigitizer") == 0)&& + (fSDigitsTitle.CompareTo(sdigitizerBranch->GetTitle()) == 0) ) + sdigitizerNotFound = kFALSE ; + } + } + + if(!(sdigitizerNotFound && phosNotFound)){ + cout << "AliPHOSSdigitizer error:" << endl ; + cout << "Can not overwrite existing branches: do not write" << endl ; + return ; + } + //Make (if necessary) branches + char * file =0; + if(gSystem->Getenv("CONFIG_SPLIT_FILE")){ //generating file name + file = new char[strlen(gAlice->GetBaseFile())+20] ; + sprintf(file,"%s/PHOS.SDigits.root",gAlice->GetBaseFile()) ; + } + + TDirectory *cwd = gDirectory; + + //First list of sdigits + Int_t bufferSize = 32000 ; + sdigitsBranch = gAlice->TreeS()->Branch("PHOS",&fSDigits,bufferSize); + sdigitsBranch->SetTitle(fSDigitsTitle.Data()); + if (file) { + sdigitsBranch->SetFile(file); + TIter next( sdigitsBranch->GetListOfBranches()); + while ((sdigitsBranch=(TBranch*)next())) { + sdigitsBranch->SetFile(file); + } + cwd->cd(); + } + + //second - SDigitizer + Int_t splitlevel = 0 ; + AliPHOSSDigitizer * sd = this ; + sdigitizerBranch = gAlice->TreeS()->Branch("AliPHOSSDigitizer","AliPHOSSDigitizer", + &sd,bufferSize,splitlevel); + sdigitizerBranch->SetTitle(fSDigitsTitle.Data()); + if (file) { + sdigitizerBranch->SetFile(file); + TIter next( sdigitizerBranch->GetListOfBranches()); + while ((sdigitizerBranch=(TBranch*)next())) { + sdigitizerBranch->SetFile(file); + } + cwd->cd(); + delete file; + } + gAlice->TreeS()->Fill() ; gAlice->TreeS()->Write(0,TObject::kOverwrite) ; + + if(strstr(option,"deb")) + PrintSDigits(option) ; + } - - delete sdigits ; - if(file) - file->Close() ; - + + if(strstr(option,"tim")){ + gBenchmark->Stop("PHOSSDigitizer"); + cout << "AliPHOSSDigitizer:" << endl ; + cout << " took " << gBenchmark->GetCpuTime("PHOSSDigitizer") << " seconds for SDigitizing " + << gBenchmark->GetCpuTime("PHOSSDigitizer")/fNevents << " seconds per event " << endl ; + cout << endl ; + } + + } //__________________________________________________________________ -void AliPHOSSDigitizer::SetSDigitsFile(char * file ){ - if(!fSDigitsFile.IsNull()) - cout << "Changing SDigits file from " <<(char *)fSDigitsFile.Data() << " to " << file << endl ; - fSDigitsFile=file ; +void AliPHOSSDigitizer::SetSDigitsBranch(const char * title ){ + //Seting title to branch SDigits + if(!fSDigitsTitle.IsNull()) + cout << "AliPHOSSdigitizer: changing SDigits file from " <GetEntriesFast() << endl ; + cout << endl ; + + if(strstr(option,"all")){// print all digits + + //loop over digits + AliPHOSDigit * digit; + cout << "SDigit Id " << " Amplitude " << " Index " << " Nprim " << " Primaries list " << endl; + Int_t index ; + for (index = 0 ; index < fSDigits->GetEntries() ; index++) { + digit = (AliPHOSDigit * ) fSDigits->At(index) ; + cout << setw(8) << digit->GetId() << " " << setw(3) << digit->GetAmp() << " " + << setw(6) << digit->GetIndexInList() << " " + << setw(5) << digit->GetNprimary() <<" "; + + Int_t iprimary; + for (iprimary=0; iprimaryGetNprimary(); iprimary++) + cout << setw(5) << digit->GetPrimary(iprimary+1) << " "; + cout << endl; + } + + } +} diff --git a/PHOS/AliPHOSSDigitizer.h b/PHOS/AliPHOSSDigitizer.h index c45259170e3..c3f4ae07a05 100644 --- a/PHOS/AliPHOSSDigitizer.h +++ b/PHOS/AliPHOSSDigitizer.h @@ -22,30 +22,40 @@ class AliPHOSSDigitizer: public TTask { public: AliPHOSSDigitizer() ; // ctor - AliPHOSSDigitizer(char* HeaderFile,char *SdigitsFile = 0) ; - + AliPHOSSDigitizer(const char* HeaderFile,const char *SdigitsTitle = 0) ; virtual ~AliPHOSSDigitizer() ; // dtor - Float_t Calibrate(Int_t amp){return (amp - fA)/fB ; } - Int_t Digitize(Float_t Energy){ return (Int_t ) ( fA + Energy*fB); } - Float_t GetPedestalParameter(){return fA;} - Float_t GetCalibrationParameter(){return fB;} - char *GetSDigitsFile()const{return (char*) fSDigitsFile.Data();} + Float_t Calibrate(Int_t amp)const {return (amp - fA)/fB ; } + Int_t Digitize(Float_t Energy)const { return (Int_t ) ( fA + Energy*fB); } + virtual void Exec(Option_t *option); - void SetNEvents(Int_t Nevents){fNevents = Nevents;} - void SetPedestalParameter(Float_t A){fA = A ;} - void SetSlopeParameter(Float_t B){fB = B ;} - void SetSDigitsFile(char * file ) ; + + Float_t GetPedestalParameter()const {return fA;} + Float_t GetCalibrationParameter()const{return fB;} + char * GetSDigitsBranch()const{return (char*) fSDigitsTitle.Data();} + virtual void Print(Option_t* option) const ; - Bool_t operator == (const AliPHOSSDigitizer & sd) const ; + + void SetPedestalParameter(Float_t A){fA = A ;} + void SetSlopeParameter(Float_t B){fB = B ;} + void SetSDigitsBranch(const char * title ) ; + + Bool_t operator == (const AliPHOSSDigitizer & sd) const ; + +private: + void Init() ; + void PrintSDigits(Option_t * option) ; private: Float_t fA ; //Pedestal parameter Float_t fB ; //Slope Digitizition parameters Int_t fNevents ; // Number of events to digitize Float_t fPrimThreshold ; // To store primari if Elos > threshold - TString fSDigitsFile ; //output file + TString fSDigitsTitle ; // title of SDigits branch TString fHeadersFile ; //input file + Bool_t fIsInitialized ; + TClonesArray * fSDigits ; //! list of SDigits + TClonesArray * fHits ; //! ClassDef(AliPHOSSDigitizer,1) // description -- 2.31.1