]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - EMCAL/AliEMCALPIDv1.cxx
Small changes
[u/mrichter/AliRoot.git] / EMCAL / AliEMCALPIDv1.cxx
index 7184d7b1c78e9643fc385b004be30e06e41026ca..848c3086521ea102037670f39ba8c16bac82206b 100644 (file)
 #include "TFolder.h"
 #include "TSystem.h"
 #include "TBenchmark.h"
-
 #include "TSystem.h"
-
-// --- Standard library ---
-
-#include <Riostream.h>
-
-// --- AliRoot header files ---
-
-#include "AliRun.h"
+  
+  // --- Standard library ---
+  
+  // --- AliRoot header files ---
+  
 #include "AliGenerator.h"
 #include "AliEMCAL.h"
 #include "AliEMCALPIDv1.h"
 #include "AliEMCALRecParticle.h"
 #include "AliEMCALGeometry.h"
 #include "AliEMCALGetter.h"
-
-ClassImp( AliEMCALPIDv1) 
+  
+  ClassImp( AliEMCALPIDv1) 
 
 //____________________________________________________________________________
 AliEMCALPIDv1::AliEMCALPIDv1():AliEMCALPID()
 { 
   // default ctor
+  
   InitParameters() ; 
   fDefaultInit = kTRUE ; 
+  
+}
 
+//____________________________________________________________________________
+AliEMCALPIDv1::AliEMCALPIDv1(const AliEMCALPIDv1 & pid ):AliEMCALPID(pid)
+{ 
+  // ctor
+  InitParameters() ; 
+  Init() ;
+  
 }
 
 //____________________________________________________________________________
-AliEMCALPIDv1::AliEMCALPIDv1(const char * headerFile,const char * name, const Bool_t toSplit)
-:AliEMCALPID(headerFile, name,toSplit)
+AliEMCALPIDv1::AliEMCALPIDv1(const TString alirunFileName, const TString eventFolderName):AliEMCALPID(alirunFileName, eventFolderName)
 { 
   //ctor with the indication on where to look for the track segments
  
   InitParameters() ; 
-
   Init() ;
   fDefaultInit = kFALSE ; 
 
@@ -78,18 +81,13 @@ AliEMCALPIDv1::AliEMCALPIDv1(const char * headerFile,const char * name, const Bo
 AliEMCALPIDv1::~AliEMCALPIDv1()
 { 
   // dtor
-
-  if (!fDefaultInit) {  
-    fSplitFile = 0 ; 
-  }
 }
 
 //____________________________________________________________________________
 const TString AliEMCALPIDv1::BranchName() const 
 {  
-  TString branchName(GetName() ) ;
-  branchName.Remove(branchName.Index(Version())-1) ;
-  return branchName ;
+
+  return GetName() ;
 }
  
 //____________________________________________________________________________
@@ -98,42 +96,11 @@ void AliEMCALPIDv1::Init()
   // Make all memory allocations that are not possible in default constructor
   // Add the PID task to the list of EMCAL tasks
 
-  if ( strcmp(GetTitle(), "") == 0 )
-    SetTitle("galice.root") ;
-
-  TString branchname(GetName()) ;
-  branchname.Remove(branchname.Index(Version())-1) ;    
-  AliEMCALGetter * gime = AliEMCALGetter::GetInstance(GetTitle(),branchname.Data(),fToSplit ) ; 
 
-  //  gime->SetRecParticlesTitle(BranchName()) ;
-  if ( gime == 0 ) {
-    Error("Init", "Could not obtain the Getter object !" ) ;  
-    return ;
-  } 
+  AliEMCALGetter * gime = AliEMCALGetter::Instance(GetTitle(), fEventFolderName.Data()) ; 
 
-  fSplitFile = 0 ;
-  if(fToSplit){
-    //First - extract full path if necessary
-    TString fileName(GetTitle()) ;
-    Ssiz_t islash = fileName.Last('/') ;
-    if(islash<fileName.Length())
-      fileName.Remove(islash+1,fileName.Length()) ;
-    else
-      fileName="" ;
-    fileName+="EMCAL.RecData." ;
-    if((strcmp(branchname.Data(),"Default")!=0)&&(strcmp(branchname.Data(),"")!=0)){
-      fileName+=branchname.Data() ;
-      fileName+="." ;
-    }
-    fileName+="root" ;
-    fSplitFile = static_cast<TFile*>(gROOT->GetFile(fileName.Data()));   
-    if(!fSplitFile)
-      fSplitFile =  TFile::Open(fileName.Data(),"update") ;
-  }
-  
-  gime->PostPID(this) ;
-  gime->PostRecParticles(branchname) ; 
-  
+  if ( !gime->PID() ) 
+    gime->PostPID(this) ;
 }
 
 //____________________________________________________________________________
@@ -143,13 +110,6 @@ void AliEMCALPIDv1::InitParameters()
   fRecParticlesInRun = 0 ; 
   fNEvent            = 0 ;            
   fRecParticlesInRun = 0 ;
-  TString pidName( GetName()) ;
-  if (pidName.IsNull() ) 
-    pidName = "Default" ; 
-  pidName.Append(":") ; 
-  pidName.Append(Version()) ; 
-  SetName(pidName) ;
-  fPi0Analysis = kFALSE ;
 }
 
 //____________________________________________________________________________
@@ -157,10 +117,7 @@ void AliEMCALPIDv1::InitParameters()
 void  AliEMCALPIDv1::Exec(Option_t * option) 
 {
   //Steering method
-  
-  if( strcmp(GetName(), "")== 0 ) 
-    Init() ;
-  
+
   if(strstr(option,"tim"))
     gBenchmark->Start("EMCALPID");
   
@@ -168,34 +125,32 @@ void  AliEMCALPIDv1::Exec(Option_t * option)
     Print("") ; 
     return ; 
   }
-  AliEMCALGetter * gime = AliEMCALGetter::GetInstance() ; 
-  if(gime->BranchExists("RecParticles") )
-    return ;
-  Int_t nevents = gime->MaxEvent() ;       //(Int_t) gAlice->TreeE()->GetEntries() ;
+  AliEMCALGetter * gime = AliEMCALGetter::Instance() ; 
+
+  Int_t nevents = gime->MaxEvent() ;      
   Int_t ievent ;
 
 
   for(ievent = 0; ievent < nevents; ievent++){
-    gime->Event(ievent,"R") ;
-    MakeRecParticles() ;
-    
-    WriteRecParticles(ievent);
-    
-    if(strstr(option,"deb"))
-      PrintRecParticles(option) ;
-
-    //increment the total number of rec particles per run 
-    fRecParticlesInRun += gime->RecParticles(BranchName())->GetEntriesFast() ; 
-
+    gime->Event(ievent,"TR") ;
+    if(gime->TrackSegments() && //Skip events, where no track segments made
+       gime->TrackSegments()->GetEntriesFast()) {   
+      MakeRecParticles() ;
+      WriteRecParticles();
+      if(strstr(option,"deb"))
+       PrintRecParticles(option) ;
+      //increment the total number of rec particles per run 
+      fRecParticlesInRun += gime->RecParticles()->GetEntriesFast() ; 
+    }
   }
-  
   if(strstr(option,"tim")){
     gBenchmark->Stop("EMCALPID");
     Info("Exec", "took %f seconds for PID %f seconds per event", 
         gBenchmark->GetCpuTime("EMCALPID"),  
         gBenchmark->GetCpuTime("EMCALPID")/nevents) ;
   } 
+
+  Unload();
 }
 
 //____________________________________________________________________________
@@ -203,12 +158,12 @@ void  AliEMCALPIDv1::MakeRecParticles(){
 
   // Makes a RecParticle out of a TrackSegment
   
-  AliEMCALGetter * gime = AliEMCALGetter::GetInstance() ; 
-  TObjArray * aECRecPoints = gime->ECALRecPoints() ; 
-  TObjArray * aPRRecPoints = gime->PRERecPoints() ; 
-  TObjArray * aHCRecPoints = gime->HCALRecPoints() ; 
+  AliEMCALGetter * gime = AliEMCALGetter::Instance() ; 
+  TObjArray * aECARecPoints = gime->ECARecPoints() ; 
+  TObjArray * aPRERecPoints = gime->PRERecPoints() ; 
+  TObjArray * aHCARecPoints = gime->HCARecPoints() ; 
   TClonesArray * trackSegments = gime->TrackSegments() ; 
-  if ( !aECRecPoints || !aPRRecPoints || !aHCRecPoints || !trackSegments ) {
+  if ( !aECARecPoints || !aPRERecPoints || !aHCARecPoints || !trackSegments ) {
     Fatal("MakeRecParticles", "RecPoints or TrackSegments not found !") ;  
   }
   TClonesArray * recParticles  = gime->RecParticles() ; 
@@ -225,30 +180,30 @@ void  AliEMCALPIDv1::MakeRecParticles(){
     rp->SetTrackSegment(index) ;
     rp->SetIndexInList(index) ;
        
-    AliEMCALTowerRecPoint * ecal = 0 ;
-    if(ts->GetECIndex()>=0)
-      ecal = dynamic_cast<AliEMCALTowerRecPoint *>(aECRecPoints->At(ts->GetECIndex())) ;
+    AliEMCALTowerRecPoint * eca = 0 ;
+    if(ts->GetECAIndex()>=0)
+      eca = dynamic_cast<AliEMCALTowerRecPoint *>(aECARecPoints->At(ts->GetECAIndex())) ;
     
-    AliEMCALTowerRecPoint    * pre = 0 ;
+    AliEMCALTowerRecPoint * pre = 0 ;
     if(ts->GetPREIndex()>=0)
-      pre = dynamic_cast<AliEMCALTowerRecPoint *>(aPRRecPoints->At(ts->GetPREIndex())) ;
+      pre = dynamic_cast<AliEMCALTowerRecPoint *>(aPRERecPoints->At(ts->GetPREIndex())) ;
     
-    AliEMCALTowerRecPoint    * hcal = 0 ;
-    if(ts->GetHCIndex()>=0)
-      hcal = dynamic_cast<AliEMCALTowerRecPoint *>(aHCRecPoints->At(ts->GetHCIndex())) ;
+    AliEMCALTowerRecPoint * hca = 0 ;
+    if(ts->GetHCAIndex()>=0)
+      hca = dynamic_cast<AliEMCALTowerRecPoint *>(aHCARecPoints->At(ts->GetHCAIndex())) ;
 
     // Now set type (reconstructed) of the particle
 
     // Choose the cluster energy range
     
-    if (!ecal) {
-      Fatal("MakeRecParticles", "-> emcal(%d) = %d", ts->GetECIndex(), ecal ) ;
+    if (!eca) {
+      Fatal("MakeRecParticles", "-> emcal(%d) = %d", ts->GetECAIndex(), eca ) ;
     }
 
-    Float_t    e = ecal->GetEnergy() ;   
+    Float_t    e = eca->GetEnergy() ;   
     
     Float_t  lambda[2] ;
-    ecal->GetElipsAxis(lambda) ;
+    eca->GetElipsAxis(lambda) ;
     
     if((lambda[0]>0.01) && (lambda[1]>0.01)){
       // Looking PCA. Define and calculate the data (X),
@@ -260,16 +215,16 @@ void  AliEMCALPIDv1::MakeRecParticles(){
       if((lambda[0]+lambda[1])!=0) 
        Spher=fabs(lambda[0]-lambda[1])/(lambda[0]+lambda[1]); 
       
-      Emaxdtotal=ecal->GetMaximalEnergy()/ecal->GetEnergy(); 
+      Emaxdtotal=eca->GetMaximalEnergy()/eca->GetEnergy(); 
     }
     
     //    Float_t time = ecal->GetTime() ;
       
     //Set momentum, energy and other parameters 
-    Float_t  encal = e;
+    Float_t  enca = e;
     TVector3 dir(0., 0., 0.) ; 
-    dir.SetMag(encal) ;
-    rp->SetMomentum(dir.X(),dir.Y(),dir.Z(),encal) ;
+    dir.SetMag(enca) ;
+    rp->SetMomentum(dir.X(),dir.Y(),dir.Z(),enca) ;
     rp->SetCalcMass(0);
     rp->Name(); //If photon sets the particle pdg name to gamma
     rp->SetProductionVertex(0,0,0,0);
@@ -284,69 +239,14 @@ void  AliEMCALPIDv1::MakeRecParticles(){
 }
 
 //____________________________________________________________________________
-void  AliEMCALPIDv1:: Print()
+void  AliEMCALPIDv1:: Print(Option_t * /*option*/) const
 {
   // Print the parameters used for the particle type identification
 
-  TString message ; 
-    message  = "\n=============== AliEMCALPID1 ================\n" ;
-    message += "Making PID\n";
-    message += "    Pricipal analysis file from 0.5 to 100 %s\n" ; 
-    message += "    Name of parameters file     %s\n" ;
-    message += "    Matrix of Parameters: 9x4\n" ;
-    message += "        RCPV 2x3 rows x and z, columns function cut parameters\n" ;
-    message += "        TOF  1x3 [High Eff-Low Pur,Medium Eff-Pur, Low Eff-High Pur]\n" ;
-    message += "        PCA  5x4 [5 ellipse parametres and 4 parametres to calculate them: A/Sqrt(E) + B* E + C * E^2 + D]\n" ;
-    message += "        Energy Calibration  1x3 [3 parametres to calibrate energy: A + B* E + C * E^2]\n" ;
-    Info("Print", message.Data(), fFileName.Data(), fFileNamePar.Data() ) ; 
-}
-
-//____________________________________________________________________________
-void  AliEMCALPIDv1::WriteRecParticles(Int_t event)
-{
-  AliEMCALGetter *gime = AliEMCALGetter::GetInstance() ; 
-
-  TClonesArray * recParticles = gime->RecParticles() ; 
-  recParticles->Expand(recParticles->GetEntriesFast() ) ;
-  TTree * treeR ;
-
-  if(fToSplit){
-    if(!fSplitFile)
-      return ;
-    fSplitFile->cd() ;
-    char name[10] ;
-    sprintf(name,"%s%d", "TreeR",event) ;
-    treeR = dynamic_cast<TTree*>(fSplitFile->Get(name)); 
-  }
-  else{
-    treeR = gAlice->TreeR();
-  }
-  
-  if(!treeR){
-    gAlice->MakeTree("R", fSplitFile);
-    treeR = gAlice->TreeR() ;
-  }
-  
-  //First rp
-  Int_t bufferSize = 32000 ;    
-  TBranch * rpBranch = treeR->Branch("EMCALRP",&recParticles,bufferSize);
-  rpBranch->SetTitle(BranchName());
-
-  
-  //second, pid
-  Int_t splitlevel = 0 ; 
-  AliEMCALPIDv1 * pid = this ;
-  TBranch * pidBranch = treeR->Branch("AliEMCALPID","AliEMCALPIDv1",&pid,bufferSize,splitlevel);
-  pidBranch->SetTitle(BranchName());
-  
-  rpBranch->Fill() ;
-  pidBranch->Fill() ; 
-  
-  treeR->AutoSave() ; //Write(0,kOverwrite) ;  
-  if(gAlice->TreeR()!=treeR){
-    treeR->Delete();
-  }
+    Info("Print", "=============== AliEMCALPID1 ================") ;
+    printf("Making PID\n");
+    printf("    Pricipal analysis file from 0.5 to 100 %s\n", fFileName.Data() ) ; 
+    printf("    Name of parameters file     %s\n", fFileNamePar.Data() )  ;
 }
 
 //____________________________________________________________________________
@@ -354,9 +254,9 @@ void AliEMCALPIDv1::PrintRecParticles(Option_t * option)
 {
   // Print table of reconstructed particles
 
-  AliEMCALGetter *gime = AliEMCALGetter::GetInstance() ; 
+  AliEMCALGetter *gime = AliEMCALGetter::Instance() ; 
 
-  TClonesArray * recParticles = gime->RecParticles(BranchName()) ; 
+  TClonesArray * recParticles = gime->RecParticles() ; 
 
   TString message ; 
   message  = "\nevent " ;
@@ -382,5 +282,36 @@ void AliEMCALPIDv1::PrintRecParticles(Option_t * option)
   Info("Print", message.Data() ) ; 
 }
 
+//____________________________________________________________________________
+void AliEMCALPIDv1::Unload() 
+{
+  AliEMCALGetter * gime = AliEMCALGetter::Instance() ;  
+  gime->EmcalLoader()->UnloadRecPoints() ;
+  gime->EmcalLoader()->UnloadTracks() ;
+  gime->EmcalLoader()->UnloadRecParticles() ;
+}
+
+//____________________________________________________________________________
+void  AliEMCALPIDv1::WriteRecParticles()
+{
+  AliEMCALGetter *gime = AliEMCALGetter::Instance() ; 
+
+  TClonesArray * recParticles = gime->RecParticles() ; 
+  recParticles->Expand(recParticles->GetEntriesFast() ) ;
+  TTree * treeP =  gime->TreeP() ;
 
+  
+  
+  //First rp
+  Int_t bufferSize = 32000 ;    
+  TBranch * rpBranch = treeP->Branch("EMCALRP",&recParticles,bufferSize);
+  rpBranch->SetTitle(BranchName());
+
+  rpBranch->Fill() ;
+  gime->WriteRecParticles("OVERWRITE");
+  gime->WritePID("OVERWRITE");
+
+}