]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PHOS/AliPHOSClusterizerv1.cxx
Implement the new logging system (AliLog)
[u/mrichter/AliRoot.git] / PHOS / AliPHOSClusterizerv1.cxx
index 44b067e3d3132681b1b9f9f1eccb1579d773f57a..677633a55eebae604106543de48ab827f3999841 100644 (file)
@@ -53,6 +53,7 @@
 // --- Standard library ---
 
 // --- AliRoot header files ---
+#include "AliLog.h"
 #include "AliPHOSGetter.h"
 #include "AliPHOSGeometry.h" 
 #include "AliPHOSClusterizerv1.h"
@@ -60,6 +61,7 @@
 #include "AliPHOSCpvRecPoint.h"
 #include "AliPHOSDigit.h"
 #include "AliPHOSDigitizer.h"
+#include "AliPHOSCalibrationDB.h"
 
 ClassImp(AliPHOSClusterizerv1)
   
@@ -88,11 +90,7 @@ AliPHOSClusterizerv1::AliPHOSClusterizerv1(const TString alirunFileName, const T
 {
   // dtor
 
- //  delete fPedestals ;
-//   delete fGains ;
 }
-
-
 //____________________________________________________________________________
 const TString AliPHOSClusterizerv1::BranchName() const 
 {  
@@ -103,38 +101,46 @@ const TString AliPHOSClusterizerv1::BranchName() const
 Float_t  AliPHOSClusterizerv1::Calibrate(Int_t amp, Int_t absId) const
 {  
   //To be replaced later by the method, reading individual parameters from the database
-  //   if(fCalibrationDB)
-  //     return fCalibrationDB->Calibrate(amp,absId) ;
-  //   else{ //simulation
-  if(absId <= fEmcCrystals) //calibrate as EMC 
-    return fADCpedestalEmc + amp*fADCchanelEmc ;        
-  else //calibrate as CPV
-    return fADCpedestalCpv+ amp*fADCchanelCpv ;       
-  //  }
+  if(fCalibrationDB){
+    return fCalibrationDB->Calibrate(amp,absId) ;
+  }
+  else{ //simulation
+    if(absId <= fEmcCrystals) //calibrate as EMC 
+      return fADCpedestalEmc + amp*fADCchanelEmc ;        
+    else //calibrate as CPV
+      return fADCpedestalCpv+ amp*fADCchanelCpv ;       
+  }
 }
 
 //____________________________________________________________________________
-void AliPHOSClusterizerv1::Exec(Option_t * option)
+void AliPHOSClusterizerv1::Exec(Option_t *option)
 {
-  // Steering method
+  // Steering method to perform clusterization for events
+  // in the range from fFirstEvent to fLastEvent.
+  // This range is optionally set by SetEventRange().
+  // if fLastEvent=-1 (by default), then process events until the end.
 
   if(strstr(option,"tim"))
     gBenchmark->Start("PHOSClusterizer"); 
   
-  if(strstr(option,"print"))
+  if(strstr(option,"print")) {
     Print() ; 
+    return ;
+  }
 
   AliPHOSGetter * gime = AliPHOSGetter::Instance() ; 
+  GetCalibrationParameters() ;
   
-  Int_t nevents = gime->MaxEvent() ;
+  if (fLastEvent == -1) 
+    fLastEvent = gime->MaxEvent() - 1 ;
+  else 
+    fLastEvent = TMath::Min(fFirstEvent, gime->MaxEvent()); // one event at the time 
+  Int_t nEvents   = fLastEvent - fFirstEvent + 1;
+
   Int_t ievent ;
-  
-  for(ievent = 0; ievent < nevents; ievent++)
-   {
+  for (ievent = fFirstEvent; ievent <= fLastEvent; ievent++) {
     gime->Event(ievent, "D");
-
-    GetCalibrationParameters() ;
-
+    
     fNumberOfEmcClusters  = fNumberOfCpvClusters  = 0 ;
     
     MakeClusters() ;
@@ -152,13 +158,14 @@ void AliPHOSClusterizerv1::Exec(Option_t * option)
     fRecPointsInRun += gime->CpvRecPoints()->GetEntriesFast() ;  
     }
   
-  Unload();
+  if(fWrite) //do not unload in "on flight" mode
+    Unload();
   
   if(strstr(option,"tim")){
     gBenchmark->Stop("PHOSClusterizer");
-    Info("Exec", "  took %f seconds for Clusterizing %f seconds per event \n",
+    AliInfo(Form("took %f seconds for Clusterizing %f seconds per event \n",
         gBenchmark->GetCpuTime("PHOSClusterizer"), 
-        gBenchmark->GetCpuTime("PHOSClusterizer")/nevents ) ; 
+        gBenchmark->GetCpuTime("PHOSClusterizer")/nEvents )) ; 
   } 
 }
 
@@ -261,18 +268,19 @@ void AliPHOSClusterizerv1::GetCalibrationParameters()
 {
 
   AliPHOSGetter * gime = AliPHOSGetter::Instance();
-
-  if ( !gime->Digitizer() ) 
-    gime->LoadDigitizer();
-  AliPHOSDigitizer * dig = gime->Digitizer(); 
-  
-  fADCchanelEmc   = dig->GetEMCchannel() ;
-  fADCpedestalEmc = dig->GetEMCpedestal();
-  
-  fADCchanelCpv   = dig->GetCPVchannel() ;
-  fADCpedestalCpv = dig->GetCPVpedestal() ; 
-  
-  //   fCalibrationDB = gime->CalibrationDB();
+  if(gime->IsRawDigits()){
+    fCalibrationDB = gime->CalibrationDB();    
+  }
+  else{
+    if ( !gime->Digitizer() ) 
+      gime->LoadDigitizer();
+    AliPHOSDigitizer * dig = gime->Digitizer(); 
+    fADCchanelEmc   = dig->GetEMCchannel() ;
+    fADCpedestalEmc = dig->GetEMCpedestal();
+    
+    fADCchanelCpv   = dig->GetCPVchannel() ;
+    fADCpedestalCpv = dig->GetCPVpedestal() ; 
+  }  
 }
 
 //____________________________________________________________________________
@@ -280,18 +288,21 @@ void AliPHOSClusterizerv1::Init()
 {
   // Make all memory allocations which can not be done in default constructor.
   // Attach the Clusterizer task to the list of PHOS tasks
-   
-  AliPHOSGetter* gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName.Data());
+  AliPHOSGetter* gime = AliPHOSGetter::Instance() ;
+  if(!gime)
+    gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName.Data());
 
   AliPHOSGeometry * geom = gime->PHOSGeometry();
-  
+
   fEmcCrystals = geom->GetNModules() *  geom->GetNCristalsInModule() ;
 
   if(!gMinuit) 
     gMinuit = new TMinuit(100);
 
-  if ( !gime->Clusterizer() ) 
+  if ( !gime->Clusterizer() ) {
     gime->PostClusterizer(this);
+  }
 }
 
 //____________________________________________________________________________
@@ -306,7 +317,10 @@ void AliPHOSClusterizerv1::InitParameters()
   
   fEmcLocMaxCut            = 0.03 ;
   fCpvLocMaxCut            = 0.03 ;
-  
+
+  fEmcMinE                 = 0.01 ;
+  fCpvMinE                 = 0.0  ;
+
   fW0                      = 4.5 ;
   fW0CPV                   = 4.0 ;
 
@@ -316,6 +330,11 @@ void AliPHOSClusterizerv1::InitParameters()
     
   fRecPointsInRun          = 0 ;
 
+  fWrite                   = kTRUE ;
+
+  fCalibrationDB           = 0 ;
+
+  SetEventRange(0,-1) ;
 }
 
 //____________________________________________________________________________
@@ -361,8 +380,21 @@ Int_t AliPHOSClusterizerv1::AreNeighbours(AliPHOSDigit * d1, AliPHOSDigit * d2)c
 
   return rv ; 
 }
+//____________________________________________________________________________
+void AliPHOSClusterizerv1::CleanDigits(TClonesArray * digits){
+  for(Int_t i=0; i<digits->GetEntriesFast(); i++){
+    AliPHOSDigit * digit = static_cast<AliPHOSDigit*>(digits->At(i)) ;
+    Float_t cut = IsInEmc(digit) ? fEmcMinE : fCpvMinE ;
+    if(Calibrate(digit->GetAmp(),digit->GetId()) < cut)
+      digits->RemoveAt(i) ;
+  }
+  digits->Compress() ;
+  for (Int_t i = 0 ; i < digits->GetEntriesFast() ; i++) { 
+    AliPHOSDigit *digit = static_cast<AliPHOSDigit*>( digits->At(i) ) ; 
+    digit->SetIndexInList(i) ;     
+  }
 
-
+}
 //____________________________________________________________________________
 Bool_t AliPHOSClusterizerv1::IsInEmc(AliPHOSDigit * digit) const
 {
@@ -414,47 +446,60 @@ void AliPHOSClusterizerv1::WriteRecPoints()
   TObjArray * emcRecPoints = gime->EmcRecPoints() ; 
   TObjArray * cpvRecPoints = gime->CpvRecPoints() ; 
   TClonesArray * digits = gime->Digits() ; 
-  
-  TTree * treeR = gime->TreeR();
-
   Int_t index ;
-  //Evaluate position, dispersion and other RecPoint properties...
-  for(index = 0; index < emcRecPoints->GetEntries(); index++)
-    dynamic_cast<AliPHOSEmcRecPoint *>( emcRecPoints->At(index) )->EvalAll(fW0,digits) ;
-  
+  //Evaluate position, dispersion and other RecPoint properties..
+  Int_t nEmc = emcRecPoints->GetEntriesFast();
+  for(index = 0; index < nEmc; index++){
+    AliPHOSEmcRecPoint * rp = dynamic_cast<AliPHOSEmcRecPoint *>( emcRecPoints->At(index) );
+    rp->Purify(fEmcMinE) ;
+    if(rp->GetMultiplicity()>0.) //If this RP is not empty
+      rp->EvalAll(fW0,digits) ;
+    else{
+      emcRecPoints->RemoveAt(index) ;
+      delete rp ;
+    }
+  }
+  emcRecPoints->Compress() ;
   emcRecPoints->Sort() ;
-  for(index = 0; index < emcRecPoints->GetEntries(); index++)
+  //  emcRecPoints->Expand(emcRecPoints->GetEntriesFast()) ;
+  for(index = 0; index < emcRecPoints->GetEntries(); index++){
     dynamic_cast<AliPHOSEmcRecPoint *>( emcRecPoints->At(index) )->SetIndexInList(index) ;
-  
-  emcRecPoints->Expand(emcRecPoints->GetEntriesFast()) ; 
+  }
   
   //Now the same for CPV
-  for(index = 0; index < cpvRecPoints->GetEntries(); index++)
-    dynamic_cast<AliPHOSRecPoint *>( cpvRecPoints->At(index) )->EvalAll(digits)  ;
-  
+  for(index = 0; index < cpvRecPoints->GetEntries(); index++){
+    AliPHOSCpvRecPoint * rp = dynamic_cast<AliPHOSCpvRecPoint *>( cpvRecPoints->At(index) );
+    rp->EvalAll(fW0CPV,digits) ;
+  }
   cpvRecPoints->Sort() ;
   
   for(index = 0; index < cpvRecPoints->GetEntries(); index++)
-    dynamic_cast<AliPHOSRecPoint *>( cpvRecPoints->At(index) )->SetIndexInList(index) ;
+    dynamic_cast<AliPHOSCpvRecPoint *>( cpvRecPoints->At(index) )->SetIndexInList(index) ;
   
   cpvRecPoints->Expand(cpvRecPoints->GetEntriesFast()) ;
   
-  Int_t bufferSize = 32000 ;    
-  Int_t splitlevel = 0 ;
-  //First EMC
-  TBranch * emcBranch = treeR->Branch("PHOSEmcRP","TObjArray",&emcRecPoints,bufferSize,splitlevel);
-  emcBranch->SetTitle(BranchName());
+  if(fWrite){ //We write TreeR
+    TTree * treeR = gime->TreeR();
     
-  //Now CPV branch
-  TBranch * cpvBranch = treeR->Branch("PHOSCpvRP","TObjArray",&cpvRecPoints,bufferSize,splitlevel);
-  cpvBranch->SetTitle(BranchName());
+    Int_t bufferSize = 32000 ;    
+    Int_t splitlevel = 0 ;
     
-  emcBranch        ->Fill() ;
-  cpvBranch        ->Fill() ;
-  
-  gime->WriteRecPoints("OVERWRITE");
-  gime->WriteClusterizer("OVERWRITE");
+    //First EMC
+    TBranch * emcBranch = treeR->Branch("PHOSEmcRP","TObjArray",&emcRecPoints,bufferSize,splitlevel);
+    emcBranch->SetTitle(BranchName());
+    
+    //Now CPV branch
+    TBranch * cpvBranch = treeR->Branch("PHOSCpvRP","TObjArray",&cpvRecPoints,bufferSize,splitlevel);
+    cpvBranch->SetTitle(BranchName());
+    
+    emcBranch        ->Fill() ;
+    cpvBranch        ->Fill() ;
+    
+    gime->WriteRecPoints("OVERWRITE");
+    gime->WriteClusterizer("OVERWRITE");
+  }
 }
 
 //____________________________________________________________________________
@@ -474,6 +519,10 @@ void AliPHOSClusterizerv1::MakeClusters()
   
   TClonesArray * digits = gime->Digits() ; 
 
+  //Remove digits below threshold
+  CleanDigits(digits) ;
+
+
   TClonesArray * digitsC =  static_cast<TClonesArray*>( digits->Clone() ) ;
   
  
@@ -913,7 +962,7 @@ void AliPHOSClusterizerv1::Print()const
   else 
     message = " AliPHOSClusterizerv1 not initialized " ;
   
-  Info("Print", message.Data(),  
+  AliInfo(Form("%s, %s %s %s %s %s %s %s %s %s %s", message.Data(),  
        taskName.Data(), 
        GetTitle(),
        taskName.Data(), 
@@ -923,7 +972,7 @@ void AliPHOSClusterizerv1::Print()const
        fW0, 
        fCpvClusteringThreshold, 
        fCpvLocMaxCut, 
-       fW0CPV ) ; 
+       fW0CPV )) ; 
 }
 
 
@@ -937,9 +986,10 @@ void AliPHOSClusterizerv1::PrintRecPoints(Option_t * option)
   TObjArray * emcRecPoints = gime->EmcRecPoints() ; 
   TObjArray * cpvRecPoints = gime->CpvRecPoints() ; 
 
-  printf("\nevent %d \n", gAlice->GetEvNumber())  ;
-  printf("       Found %d EMC RecPoints and %d CPV RecPoints \n",
-         emcRecPoints->GetEntriesFast(),cpvRecPoints->GetEntriesFast())  ; 
+  AliInfo(Form("\nevent %d \n    Found %d EMC RecPoints and %d CPV RecPoints", 
+              gAlice->GetEvNumber(),
+              emcRecPoints->GetEntriesFast(),
+              cpvRecPoints->GetEntriesFast() ))  ;
  
   fRecPointsInRun +=  emcRecPoints->GetEntriesFast() ; 
   fRecPointsInRun +=  cpvRecPoints->GetEntriesFast() ; 
@@ -972,14 +1022,14 @@ void AliPHOSClusterizerv1::PrintRecPoints(Option_t * option)
     printf("\n CPV clusters \n") ;
     printf("Index    Ene(MeV) Module     X     Y    Z  \n") ;      
     for (index = 0 ; index < cpvRecPoints->GetEntries() ; index++) {
-       AliPHOSRecPoint * rp = (AliPHOSRecPoint * )cpvRecPoints->At(index) ; 
-       
-       TVector3  locpos;  
-       rp->GetLocalPosition(locpos);
-       
-       printf("\n%6d  %8.2f  %2d     %4.1f    %4.1f %4.1f \n", 
-               rp->GetIndexInList(), rp->GetEnergy(), rp->GetPHOSMod(), 
-               locpos.X(), locpos.Y(), locpos.Z()) ; 
+      AliPHOSCpvRecPoint * rp = (AliPHOSCpvRecPoint * )cpvRecPoints->At(index) ; 
+      
+      TVector3  locpos;  
+      rp->GetLocalPosition(locpos);
+      
+      printf("\n%6d  %8.2f  %2d     %4.1f    %4.1f %4.1f \n", 
+            rp->GetIndexInList(), rp->GetEnergy(), rp->GetPHOSMod(), 
+            locpos.X(), locpos.Y(), locpos.Z()) ; 
     }
   }
 }