A cut on minimal digit energy added; OnFlight mode added

diff --git a/PHOS/AliPHOSClusterizerv1.cxx b/PHOS/AliPHOSClusterizerv1.cxx
index 8dda769c6373df1c003386fefbf6457e0cc53bab..3b2a10d326a050889240c76fe69113cea53554d0 100644 (file)
--- a/PHOS/AliPHOSClusterizerv1.cxx
+++ b/PHOS/AliPHOSClusterizerv1.cxx
@@ -60,6 +60,7 @@
 #include "AliPHOSCpvRecPoint.h"
 #include "AliPHOSDigit.h"
 #include "AliPHOSDigitizer.h"
+#include "AliPHOSCalibrationDB.h"
 
 ClassImp(AliPHOSClusterizerv1)
   
@@ -88,11 +89,7 @@ AliPHOSClusterizerv1::AliPHOSClusterizerv1(const TString alirunFileName, const T
 {
   // dtor
 
- //  delete fPedestals ;
-//   delete fGains ;
 }
-
-
 //____________________________________________________________________________
 const TString AliPHOSClusterizerv1::BranchName() const 
 {  
@@ -103,14 +100,15 @@ 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 ;       
+  }
 }
 
 //____________________________________________________________________________
@@ -129,20 +127,19 @@ void AliPHOSClusterizerv1::Exec(Option_t *option)
     return ;
   }
 
-  AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ; 
+  AliPHOSGetter * gime = AliPHOSGetter::Instance() ; 
+  GetCalibrationParameters() ;
   
   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 = fFirstEvent; ievent <= fLastEvent; ievent++) {
     gime->Event(ievent, "D");
     
-    GetCalibrationParameters() ;
-
     fNumberOfEmcClusters  = fNumberOfCpvClusters  = 0 ;
     
     MakeClusters() ;
@@ -160,7 +157,8 @@ 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");
@@ -270,16 +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() ; 
+  }  
 }
 
 //____________________________________________________________________________
@@ -287,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);
+  }
 }
 
 //____________________________________________________________________________
@@ -313,7 +317,10 @@ void AliPHOSClusterizerv1::InitParameters()
   
   fEmcLocMaxCut            = 0.03 ;
   fCpvLocMaxCut            = 0.03 ;
-  
+
+  fEmcMinE                 = 0.01 ;
+  fCpvMinE                 = 0.0  ;
+
   fW0                      = 4.5 ;
   fW0CPV                   = 4.0 ;
 
@@ -323,6 +330,8 @@ void AliPHOSClusterizerv1::InitParameters()
     
   fRecPointsInRun          = 0 ;
 
+  fWrite                   = kTRUE ;
+
   SetEventRange(0,-1) ;
 }
 
@@ -369,8 +378,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
 {
@@ -422,24 +444,33 @@ 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<AliPHOSCpvRecPoint *>( cpvRecPoints->At(index) )->EvalAll(fW0CPV,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++)
@@ -447,22 +478,26 @@ void AliPHOSClusterizerv1::WriteRecPoints()
   
   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");
+  }
 }
 
 //____________________________________________________________________________
@@ -482,6 +517,10 @@ void AliPHOSClusterizerv1::MakeClusters()
   
   TClonesArray * digits = gime->Digits() ; 
 
+  //Remove digits below threshold
+  CleanDigits(digits) ;
+
+
   TClonesArray * digitsC =  static_cast<TClonesArray*>( digits->Clone() ) ;
   
  

diff --git a/PHOS/AliPHOSClusterizerv1.h b/PHOS/AliPHOSClusterizerv1.h
index e85bd5288540aebe56a5076f237d0493a434f8f6..d6ecb18a79b344cbdbdf2f5d2cd302305562feb2 100644 (file)
--- a/PHOS/AliPHOSClusterizerv1.h
+++ b/PHOS/AliPHOSClusterizerv1.h
@@ -15,7 +15,7 @@
 //*-- Author: Yves Schutz (SUBATECH)
 
 // --- ROOT system ---
-
+class TClonesArray ;
 // --- Standard library ---
 
 // --- AliRoot header files ---
@@ -62,6 +62,8 @@ public:
 
   void Print()const ;
 
+  void SetEmcMinE(Float_t e){fEmcMinE = e ;}
+  void SetCpvMinE(Float_t e){fCpvMinE = e ;}
   virtual void SetEmcClusteringThreshold(Float_t cluth)  { fEmcClusteringThreshold = cluth ; }
   virtual void SetEmcLocalMaxCut(Float_t cut)            { fEmcLocMaxCut = cut ; }
   virtual void SetEmcLogWeight(Float_t w)                { fW0 = w ; }
@@ -69,11 +71,13 @@ public:
   virtual void SetCpvClusteringThreshold(Float_t cluth)  { fCpvClusteringThreshold = cluth ; }
   virtual void SetCpvLocalMaxCut(Float_t cut)            { fCpvLocMaxCut = cut ; }
   virtual void SetCpvLogWeight(Float_t w)                { fW0CPV = w ; }
-  virtual void SetUnfolding(Bool_t toUnfold = kTRUE )    { fToUnfold = toUnfold ;}  
+  virtual void SetUnfolding(Bool_t toUnfold = kTRUE )    { fToUnfold = toUnfold ;}
+  //Switch to "on flyght" mode, without writing to TreeR and file  
+  void SetWriting(const Bool_t toWrite = kFALSE){fWrite = toWrite;} 
   static Double_t ShowerShape(Double_t r) ; // Shape of EM shower used in unfolding; 
                                             //class member function (not object member function)
   static void UnfoldingChiSquare(Int_t & nPar, Double_t * Grad, Double_t & fret, Double_t * x, Int_t iflag)  ;
-                                            // Chi^2 of the fit. Should be static to be passes to MINUIT
+                                            // Chi^2 of the fit. Should be static to be passed to MINUIT
   void Unload() ; 
   virtual const char * Version() const { return "clu-v1"; }  
 
@@ -83,7 +87,7 @@ protected:
   virtual void   MakeClusters( ) ;            
   virtual Bool_t IsInEmc (AliPHOSDigit * digit)const ;     // Tells if id digit is in EMC
   virtual Bool_t IsInCpv (AliPHOSDigit * digit)const ;     // Tells if id digit is in CPV
-
+  void           CleanDigits(TClonesArray * digits) ;
   
 private:
 
@@ -107,6 +111,7 @@ private:
   Int_t   fEmcCrystals ;             // number of EMC cristalls in PHOS
 
   Bool_t  fToUnfold ;                // To perform unfolding 
+  Bool_t  fWrite ;                   // Write RecPoints to TreeR  
 
   Int_t   fNumberOfEmcClusters ;     // number of EMC clusters found 
   Int_t   fNumberOfCpvClusters ;     // number of CPV clusters found
@@ -118,8 +123,10 @@ private:
   Float_t fADCchanelCpv ;           // width of one ADC channel in CPV 'popugais'
   Float_t fADCpedestalCpv ;         // 
 
-  Float_t fEmcClusteringThreshold ;  // minimum energy to include a EMC digit in a cluster
-  Float_t fCpvClusteringThreshold ;  // minimum energy to include a CPV digit in a cluster
+  Float_t fEmcClusteringThreshold ;  // minimum energy to start EMC cluster
+  Float_t fCpvClusteringThreshold ;  // minimum energy to start CPV cluster
+  Float_t fEmcMinE ;                 // minimum energy of digit to be included into cluster
+  Float_t fCpvMinE ;                 // minimum energy of digit to be included into cluster
   Float_t fEmcLocMaxCut ;            // minimum energy difference to distinguish local maxima in a cluster
   Float_t fW0 ;                      // logarithmic weight for the cluster center of gravity calculation
   Float_t fCpvLocMaxCut ;            // minimum energy difference to distinguish local maxima in a CPV cluster