]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PHOS/AliPHOSClusterizerv1.cxx
Coding convention
[u/mrichter/AliRoot.git] / PHOS / AliPHOSClusterizerv1.cxx
index 75a4420c623620169d9b020eeca8c837671b8d5e..014cd6e345c4075310d6271b7f1df87dc69a761c 100644 (file)
 
 /* $Id$ */
 
-/* $Log:
-   1 October 2000. Yuri Kharlov:
-     AreNeighbours()
-     PPSD upper layer is considered if number of layers>1
-
-   18 October 2000. Yuri Kharlov:
-     AliPHOSClusterizerv1()
-     CPV clusterizing parameters added
-
-     MakeClusters()
-     After first PPSD digit remove EMC digits only once
-*/
 //*-- Author: Yves Schutz (SUBATECH)  & Dmitri Peressounko (SUBATECH & Kurchatov Institute)
 //////////////////////////////////////////////////////////////////////////////
 //  Clusterization class. Performs clusterization (collects neighbouring active cells) and 
 //  This TTask is normally called from Reconstructioner, but can as well be used in 
 //  standalone mode.
 // Use Case:
-//  root [0] AliPHOSClusterizerv1 * cl = new AliPHOSClusterizerv1("galice.root")  
+//  root [0] AliPHOSClusterizerv1 * cl = new AliPHOSClusterizerv1("galice.root", "recpointsname", "digitsname")  
 //  Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
-//               //reads gAlice from header file "..."                      
+//               // reads gAlice from header file "galice.root", uses digits stored in the branch names "digitsname" (default = "Default")
+//               // and saves recpoints in branch named "recpointsname" (default = "digitsname")                       
 //  root [1] cl->ExecuteTask()  
 //               //finds RecPoints in all events stored in galice.root
 //  root [2] cl->SetDigitsBranch("digits2") 
 
 // --- ROOT system ---
 
-#include "TROOT.h" 
-#include "TFile.h" 
-#include "TFolder.h" 
 #include "TMath.h" 
 #include "TMinuit.h"
 #include "TTree.h" 
-#include "TSystem.h" 
 #include "TBenchmark.h"
 
 // --- Standard library ---
 
-#include <iostream.h>
-#include <iomanip.h>
-
 // --- AliRoot header files ---
-
+#include "AliPHOSGetter.h"
+#include "AliPHOSGeometry.h" 
 #include "AliPHOSClusterizerv1.h"
+#include "AliPHOSEmcRecPoint.h"
 #include "AliPHOSCpvRecPoint.h"
 #include "AliPHOSDigit.h"
 #include "AliPHOSDigitizer.h"
-#include "AliPHOSEmcRecPoint.h"
-#include "AliPHOS.h"
-#include "AliPHOSGetter.h"
-#include "AliRun.h"
 
 ClassImp(AliPHOSClusterizerv1)
   
@@ -88,165 +68,119 @@ ClassImp(AliPHOSClusterizerv1)
 {
   // default ctor (to be used mainly by Streamer)
   
-  fNumberOfCpvClusters     = 0 ; 
-  fNumberOfEmcClusters     = 0 ; 
-  
-  fCpvClusteringThreshold  = 0.0;
-  fEmcClusteringThreshold  = 0.2;   
-  
-  fEmcLocMaxCut            = 0.03 ;
-  fCpvLocMaxCut            = 0.03 ;
-  
-  fW0                      = 4.5 ;
-  fW0CPV                   = 4.0 ;
-
-  fEmcTimeGate             = 1.e-8 ; 
-
-  fToUnfold                = kTRUE ;
-
-  fHeaderFileName          = "" ; 
-  fDigitsBranchTitle       = "" ;
-  fRecPointsInRun          = 0 ;   
+  InitParameters() ; 
+  fDefaultInit = kTRUE ; 
 }
 
 //____________________________________________________________________________
-AliPHOSClusterizerv1::AliPHOSClusterizerv1(const char* headerFile,const char* name)
-:AliPHOSClusterizer(headerFile, name)
+AliPHOSClusterizerv1::AliPHOSClusterizerv1(const TString alirunFileName, const TString eventFolderName)
+:AliPHOSClusterizer(alirunFileName, eventFolderName)
 {
   // ctor with the indication of the file where header Tree and digits Tree are stored
   
-
-  fNumberOfCpvClusters     = 0 ; 
-  fNumberOfEmcClusters     = 0 ; 
-  
-  fCpvClusteringThreshold  = 0.0;
-  fEmcClusteringThreshold  = 0.2;   
-  
-  fEmcLocMaxCut            = 0.03 ;
-  fCpvLocMaxCut            = 0.03 ;
-  
-  fW0                      = 4.5 ;
-  fW0CPV                   = 4.0 ;
-
-  fEmcTimeGate             = 1.e-8 ; 
-  
-  fToUnfold                = kTRUE ;
-  
-  fHeaderFileName          = GetTitle() ; 
-  fDigitsBranchTitle       = GetName() ;
-  
-  TString clusterizerName( GetName()) ; 
-  clusterizerName.Append(":") ; 
-  clusterizerName.Append(Version()) ; 
-  SetName(clusterizerName) ;
-  fRecPointsInRun          = 0 ; 
-
+  InitParameters() ;
   Init() ;
-
+  fDefaultInit = kFALSE ; 
 }
+
 //____________________________________________________________________________
   AliPHOSClusterizerv1::~AliPHOSClusterizerv1()
 {
+  // dtor
+
+ //  delete fPedestals ;
+//   delete fGains ;
+}
+
+
+//____________________________________________________________________________
+const TString AliPHOSClusterizerv1::BranchName() const 
+{  
+  return GetName();
 }
 //____________________________________________________________________________
 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 //Digitize as CPV
+    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
-
-  if( strcmp(GetName(), "")== 0 ) 
-    Init() ;
+  // 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"))
-    Print("") ; 
+  if(strstr(option,"print")) {
+    Print() ; 
+    return ;
+  }
 
-  gAlice->GetEvent(0) ;
+  AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ; 
   
- //check, if the branch with name of this" already exits?
-  TObjArray * lob = (TObjArray*)gAlice->TreeR()->GetListOfBranches() ;
-  TIter next(lob) ; 
-  TBranch * branch = 0 ;  
-  Bool_t phosemcfound = kFALSE, phoscpvfound = kFALSE, clusterizerfound = kFALSE ; 
-
-  TString branchname = GetName() ;
-  branchname.Remove(branchname.Index(Version())-1) ;
+  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;
   
-  while ( (branch = (TBranch*)next()) && (!phosemcfound || !phoscpvfound || !clusterizerfound) ) {
-    if ( (strcmp(branch->GetName(), "PHOSEmcRP")==0) && (strcmp(branch->GetTitle(), branchname.Data())==0) ) 
-      phosemcfound = kTRUE ;
-    
-    else if ( (strcmp(branch->GetName(), "PHOSCpvRP")==0) && (strcmp(branch->GetTitle(), branchname.Data())==0) ) 
-      phoscpvfound = kTRUE ;
-   
-    else if ((strcmp(branch->GetName(), "AliPHOSClusterizer")==0) && (strcmp(branch->GetTitle(), GetName())==0) ) 
-      clusterizerfound = kTRUE ; 
-  }
-
-  if ( phoscpvfound || phosemcfound || clusterizerfound ) {
-    cerr << "WARNING: AliPHOSClusterizer::Exec -> Emc(Cpv)RecPoints and/or Clusterizer branch with name " 
-        << branchname.Data() << " already exits" << endl ;
-    return ; 
-  }       
-
-  AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ;
-  Int_t nevents = (Int_t) gAlice->TreeE()->GetEntries() ;
   Int_t ievent ;
+  for (ievent = fFirstEvent; ievent <= fLastEvent; ievent++) {
+    gime->Event(ievent, "D");
 
-  for(ievent = 0; ievent < nevents; ievent++){
-
-    if(ievent == 0) GetCalibrationParameters() ;
+    GetCalibrationParameters() ;
 
-    fNumberOfEmcClusters  = 0 ;
-    fNumberOfCpvClusters  = 0 ;
-   
-    gime->Event(ievent,"D") ;
-    
-    //    if(!ReadDigits(ievent))   continue;  //reads digits for event ievent
+    fNumberOfEmcClusters  = fNumberOfCpvClusters  = 0 ;
     
     MakeClusters() ;
-    
-    if(fToUnfold)             MakeUnfolding() ;
 
-    WriteRecPoints(ievent) ;
+    if(fToUnfold)             
+      MakeUnfolding() ;
+
+    WriteRecPoints();
 
     if(strstr(option,"deb"))  
       PrintRecPoints(option) ;
 
-    //increment the total number of digits per run 
+    //increment the total number of recpoints per run 
     fRecPointsInRun += gime->EmcRecPoints()->GetEntriesFast() ;  
     fRecPointsInRun += gime->CpvRecPoints()->GetEntriesFast() ;  
- }
+    }
+  
+  Unload();
   
   if(strstr(option,"tim")){
     gBenchmark->Stop("PHOSClusterizer");
-    cout << "AliPHOSClusterizer:" << endl ;
-    cout << "  took " << gBenchmark->GetCpuTime("PHOSClusterizer") << " seconds for Clusterizing " 
-        <<  gBenchmark->GetCpuTime("PHOSClusterizer")/nevents << " seconds per event " << endl ;
-    cout << endl ;
-  }
-  
+    Info("Exec", "  took %f seconds for Clusterizing %f seconds per event \n",
+        gBenchmark->GetCpuTime("PHOSClusterizer"), 
+        gBenchmark->GetCpuTime("PHOSClusterizer")/nEvents ) ; 
+  } 
 }
 
 //____________________________________________________________________________
-Bool_t AliPHOSClusterizerv1::FindFit(AliPHOSEmcRecPoint * emcRP, int * maxAt, Float_t * maxAtEnergy,
+Bool_t AliPHOSClusterizerv1::FindFit(AliPHOSEmcRecPoint * emcRP, AliPHOSDigit ** maxAt, Float_t * maxAtEnergy,
                                    Int_t nPar, Float_t * fitparameters) const
 { 
   // Calls TMinuit to fit the energy distribution of a cluster with several maxima 
   // The initial values for fitting procedure are set equal to the positions of local maxima.
   // Cluster will be fitted as a superposition of nPar/3 electromagnetic showers
 
-  AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ; 
-  TClonesArray * digits = gime->Digits() ; 
+  
+  AliPHOSGetter * gime = AliPHOSGetter::Instance();
+  TClonesArray * digits = gime->Digits(); 
   
 
   gMinuit->mncler();                     // Reset Minuit's list of paramters
@@ -272,7 +206,7 @@ Bool_t AliPHOSClusterizerv1::FindFit(AliPHOSEmcRecPoint * emcRP, int * maxAt, Fl
   const AliPHOSGeometry * geom = gime->PHOSGeometry() ; 
 
   for(iDigit = 0; iDigit < nDigits; iDigit++){
-    digit = (AliPHOSDigit *) maxAt[iDigit]; 
+    digit = maxAt[iDigit]; 
 
     Int_t relid[4] ;
     Float_t x = 0.;
@@ -285,19 +219,19 @@ Bool_t AliPHOSClusterizerv1::FindFit(AliPHOSEmcRecPoint * emcRP, int * maxAt, Fl
     gMinuit->mnparm(index, "x",  x, 0.1, 0, 0, ierflg) ;
     index++ ;   
     if(ierflg != 0){ 
-      cout << "PHOS Unfolding>  Unable to set initial value for fit procedure : x = " << x << endl ;
+      Warning("FindFit", "PHOS Unfolding unable to set initial value for fit procedure : x = %f\n", x ) ;
       return kFALSE;
     }
     gMinuit->mnparm(index, "z",  z, 0.1, 0, 0, ierflg) ;
     index++ ;   
     if(ierflg != 0){
-      cout << "PHOS Unfolding>  Unable to set initial value for fit procedure : z = " << z << endl ;
+       Warning("FindFit", "PHOS Unfolding unable to set initial value for fit procedure : z =%f\n", z ) ;
       return kFALSE;
     }
     gMinuit->mnparm(index, "Energy",  energy , 0.05*energy, 0., 4.*energy, ierflg) ;
     index++ ;   
     if(ierflg != 0){
-      cout << "PHOS Unfolding>  Unable to set initial value for fit procedure : energy = " << energy << endl ;      
+      Warning("FindFit", "PHOS Unfolding unable to set initial value for fit procedure : energy = %f\n", energy ) ;      
       return kFALSE;
     }
   }
@@ -315,7 +249,7 @@ Bool_t AliPHOSClusterizerv1::FindFit(AliPHOSEmcRecPoint * emcRP, int * maxAt, Fl
   gMinuit->mnexcm("MIGRAD", &p0, 0, ierflg) ;    // minimize 
 
   if(ierflg == 4){  // Minimum not found   
-    cout << "PHOS Unfolding>  Fit not converged, cluster abandoned "<< endl ;      
+    Warning("FindFit", "PHOS Unfolding fit not converged, cluster abandoned\n" );      
     return kFALSE ;
   }            
   for(index = 0; index < nPar; index++){
@@ -333,48 +267,63 @@ Bool_t AliPHOSClusterizerv1::FindFit(AliPHOSEmcRecPoint * emcRP, int * maxAt, Fl
 //____________________________________________________________________________
 void AliPHOSClusterizerv1::GetCalibrationParameters() 
 {
-  AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ;
-  TString branchname = GetName() ;
-  branchname.Remove(branchname.Index(Version())-1) ;
-  AliPHOSDigitizer * dig = gime->Digitizer(branchname) ;
+
+  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();
 }
+
 //____________________________________________________________________________
 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
-  
-  if ( strcmp(GetTitle(), "") == 0 )
-    SetTitle("galice.root") ;
-
-  TString branchname = GetName() ;
-  branchname.Remove(branchname.Index(Version())-1) ;
+   
+  AliPHOSGetter* gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName.Data());
 
-  AliPHOSGetter * gime = AliPHOSGetter::GetInstance(GetTitle(), branchname) ; 
-  if ( gime == 0 ) {
-    cerr << "ERROR: AliPHOSClusterizerv1::Init -> Could not obtain the Getter object !" << endl ; 
-    return ;
-  } 
-    
-  const AliPHOSGeometry * geom = gime->PHOSGeometry() ;
+  AliPHOSGeometry * geom = gime->PHOSGeometry();
+  
   fEmcCrystals = geom->GetNModules() *  geom->GetNCristalsInModule() ;
 
   if(!gMinuit) 
-    gMinuit = new TMinuit(100) ;
+    gMinuit = new TMinuit(100);
 
-  gime->PostClusterizer(this) ;
-  // create a folder on the white board //YSAlice/WhiteBoard/RecPoints/PHOS/recpointsName
-  gime->PostRecPoints(branchname ) ;
+  if ( !gime->Clusterizer() ) 
+    gime->PostClusterizer(this);
+}
 
-  gime->PostDigits(branchname) ;
-  gime->PostDigitizer(branchname) ;
+//____________________________________________________________________________
+void AliPHOSClusterizerv1::InitParameters()
+{
+
+  fNumberOfCpvClusters     = 0 ; 
+  fNumberOfEmcClusters     = 0 ; 
   
+  fCpvClusteringThreshold  = 0.0;
+  fEmcClusteringThreshold  = 0.2;   
+  
+  fEmcLocMaxCut            = 0.03 ;
+  fCpvLocMaxCut            = 0.03 ;
+  
+  fW0                      = 4.5 ;
+  fW0CPV                   = 4.0 ;
+
+  fEmcTimeGate             = 1.e-8 ; 
+  
+  fToUnfold                = kTRUE ;
+    
+  fRecPointsInRun          = 0 ;
+
+  SetEventRange(0,-1) ;
 }
 
 //____________________________________________________________________________
@@ -387,7 +336,7 @@ Int_t AliPHOSClusterizerv1::AreNeighbours(AliPHOSDigit * d1, AliPHOSDigit * d2)c
   // The order of d1 and d2 is important: first (d1) should be a digit already in a cluster 
   //                                      which is compared to a digit (d2)  not yet in a cluster  
 
-  const AliPHOSGeometry * geom = AliPHOSGetter::GetInstance()->PHOSGeometry() ;
+  AliPHOSGeometry * geom = AliPHOSGetter::Instance()->PHOSGeometry() ;
 
   Int_t rv = 0 ; 
 
@@ -407,7 +356,7 @@ Int_t AliPHOSClusterizerv1::AreNeighbours(AliPHOSDigit * d1, AliPHOSDigit * d2)c
     }
     else {
       if((relid2[2] > relid1[2]) && (relid2[3] > relid1[3]+1)) 
-       rv = 2; //  Difference in row numbers is too large to look further 
+        rv = 2; //  Difference in row numbers is too large to look further 
     }
 
   } 
@@ -428,7 +377,7 @@ Bool_t AliPHOSClusterizerv1::IsInEmc(AliPHOSDigit * digit) const
   // Tells if (true) or not (false) the digit is in a PHOS-EMC module
  
   Bool_t rv = kFALSE ; 
-  const AliPHOSGeometry * geom = AliPHOSGetter::GetInstance()->PHOSGeometry() ;
+  AliPHOSGeometry * geom = AliPHOSGetter::Instance()->PHOSGeometry() ;
 
   Int_t nEMC = geom->GetNModules()*geom->GetNPhi()*geom->GetNZ();  
 
@@ -443,7 +392,8 @@ Bool_t AliPHOSClusterizerv1::IsInCpv(AliPHOSDigit * digit) const
   // Tells if (true) or not (false) the digit is in a PHOS-CPV module
  
   Bool_t rv = kFALSE ; 
-  const AliPHOSGeometry * geom = AliPHOSGetter::GetInstance()->PHOSGeometry() ;
+  
+  AliPHOSGeometry * geom = AliPHOSGetter::Instance()->PHOSGeometry() ;
 
   Int_t nEMC = geom->GetNModules()*geom->GetNPhi()*geom->GetNZ();
 
@@ -453,104 +403,66 @@ Bool_t AliPHOSClusterizerv1::IsInCpv(AliPHOSDigit * digit) const
 }
 
 //____________________________________________________________________________
-void AliPHOSClusterizerv1::WriteRecPoints(Int_t event)
+void AliPHOSClusterizerv1::Unload() 
+{
+  AliPHOSGetter * gime = AliPHOSGetter::Instance() ; 
+  gime->PhosLoader()->UnloadDigits() ; 
+  gime->PhosLoader()->UnloadRecPoints() ; 
+}
+
+//____________________________________________________________________________
+void AliPHOSClusterizerv1::WriteRecPoints()
 {
 
   // Creates new branches with given title
   // fills and writes into TreeR.
 
-  TString branchName(GetName() ) ;
-  branchName.Remove(branchName.Index(Version())-1) ;
+  AliPHOSGetter * gime = AliPHOSGetter::Instance();
 
-  AliPHOSGetter *gime = AliPHOSGetter::GetInstance() ; 
-  TObjArray * emcRecPoints = gime->EmcRecPoints(branchName) ; 
-  TObjArray * cpvRecPoints = gime->CpvRecPoints(branchName) ; 
-  TClonesArray * digits = gime->Digits(branchName) ; 
+  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++)
-    ((AliPHOSEmcRecPoint *)emcRecPoints->At(index))->EvalAll(fW0,digits) ;
-
+    dynamic_cast<AliPHOSEmcRecPoint *>( emcRecPoints->At(index) )->EvalAll(fW0,digits) ;
+  
   emcRecPoints->Sort() ;
-
   for(index = 0; index < emcRecPoints->GetEntries(); index++)
-    ((AliPHOSEmcRecPoint *)emcRecPoints->At(index))->SetIndexInList(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++)
-    ((AliPHOSRecPoint *)cpvRecPoints->At(index))->EvalAll(fW0CPV,digits)  ;
-
+    dynamic_cast<AliPHOSCpvRecPoint *>( cpvRecPoints->At(index) )->EvalAll(fW0CPV,digits) ;
+  
   cpvRecPoints->Sort() ;
-
-  for(index = 0; index < cpvRecPoints->GetEntries(); index++)
-    ((AliPHOSRecPoint *)cpvRecPoints->At(index))->SetIndexInList(index) ;
-
-  cpvRecPoints->Expand(cpvRecPoints->GetEntriesFast()) ;
   
-  //Make branches in TreeR for RecPoints and Clusterizer
-  char * filename = 0;
-  if(gSystem->Getenv("CONFIG_SPLIT_FILE")!=0){   //generating file name
-    filename = new char[strlen(gAlice->GetBaseFile())+20] ;
-    sprintf(filename,"%s/PHOS.Reco.root",gAlice->GetBaseFile()) ;
-  }
+  for(index = 0; index < cpvRecPoints->GetEntries(); index++)
+    dynamic_cast<AliPHOSCpvRecPoint *>( cpvRecPoints->At(index) )->SetIndexInList(index) ;
   
-  //Make new branches
-  TDirectory *cwd = gDirectory;
+  cpvRecPoints->Expand(cpvRecPoints->GetEntriesFast()) ;
   
   Int_t bufferSize = 32000 ;    
   Int_t splitlevel = 0 ;
-
   //First EMC
-  TBranch * emcBranch = gAlice->TreeR()->Branch("PHOSEmcRP","TObjArray",&emcRecPoints,bufferSize,splitlevel);
-  emcBranch->SetTitle(branchName);
-  if (filename) {
-    emcBranch->SetFile(filename);
-    TIter next( emcBranch->GetListOfBranches());
-    TBranch * sb ;
-    while ((sb=(TBranch*)next())) {
-      sb->SetFile(filename);
-    }   
-    
-    cwd->cd();
-  }
+  TBranch * emcBranch = treeR->Branch("PHOSEmcRP","TObjArray",&emcRecPoints,bufferSize,splitlevel);
+  emcBranch->SetTitle(BranchName());
     
   //Now CPV branch
-  TBranch * cpvBranch = gAlice->TreeR()->Branch("PHOSCpvRP","TObjArray",&cpvRecPoints,bufferSize,splitlevel);
-  cpvBranch->SetTitle(branchName);
-  if (filename) {
-    cpvBranch->SetFile(filename);
-    TIter next( cpvBranch->GetListOfBranches());
-    TBranch * sb;
-    while ((sb=(TBranch*)next())) {
-      sb->SetFile(filename);
-    }   
-    cwd->cd();
-  } 
+  TBranch * cpvBranch = treeR->Branch("PHOSCpvRP","TObjArray",&cpvRecPoints,bufferSize,splitlevel);
+  cpvBranch->SetTitle(BranchName());
     
-  //And Finally  clusterizer branch
-  AliPHOSClusterizerv1 * cl = (AliPHOSClusterizerv1*)gime->Clusterizer(branchName) ;
-  TBranch * clusterizerBranch = gAlice->TreeR()->Branch("AliPHOSClusterizer","AliPHOSClusterizerv1",
-                                             &cl,bufferSize,splitlevel);
-  clusterizerBranch->SetTitle(branchName);
-  if (filename) {
-    clusterizerBranch->SetFile(filename);
-    TIter next( clusterizerBranch->GetListOfBranches());
-    TBranch * sb ;
-    while ((sb=(TBranch*)next())) {
-      sb->SetFile(filename);
-    }   
-    cwd->cd();
-  }
   emcBranch        ->Fill() ;
   cpvBranch        ->Fill() ;
-  clusterizerBranch->Fill() ;
-
-  gAlice->TreeR()->Write(0,kOverwrite) ;  
   
+  gime->WriteRecPoints("OVERWRITE");
+  gime->WriteClusterizer("OVERWRITE");
 }
 
 //____________________________________________________________________________
@@ -558,18 +470,19 @@ void AliPHOSClusterizerv1::MakeClusters()
 {
   // Steering method to construct the clusters stored in a list of Reconstructed Points
   // A cluster is defined as a list of neighbour digits
-  TString branchName(GetName()) ; 
-  branchName.Remove(branchName.Index(Version())-1) ; 
 
-  AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ; 
+  
+  AliPHOSGetter * gime = AliPHOSGetter::Instance();
+
+  TObjArray * emcRecPoints = gime->EmcRecPoints() ; 
+  TObjArray * cpvRecPoints = gime->CpvRecPoints() ; 
 
-  TObjArray * emcRecPoints = gime->EmcRecPoints(branchName) ; 
-  TObjArray * cpvRecPoints = gime->CpvRecPoints(branchName) ; 
   emcRecPoints->Delete() ;
   cpvRecPoints->Delete() ;
   
-  TClonesArray * digits = gime->Digits(branchName) ; 
-  TClonesArray * digitsC =  (TClonesArray*)digits->Clone() ;
+  TClonesArray * digits = gime->Digits() ; 
+
+  TClonesArray * digitsC =  static_cast<TClonesArray*>( digits->Clone() ) ;
   
  
   // Clusterization starts  
@@ -578,7 +491,9 @@ void AliPHOSClusterizerv1::MakeClusters()
   AliPHOSDigit * digit ; 
   Bool_t notremoved = kTRUE ;
 
-  while ( (digit = (AliPHOSDigit *)nextdigit()) ) { // scan over the list of digitsC
+  while ( (digit = dynamic_cast<AliPHOSDigit *>( nextdigit()) ) ) { // scan over the list of digitsC
+    
+
     AliPHOSRecPoint * clu = 0 ; 
 
     TArrayI clusterdigitslist(1500) ;   
@@ -586,50 +501,49 @@ void AliPHOSClusterizerv1::MakeClusters()
 
     if (( IsInEmc (digit) && Calibrate(digit->GetAmp(),digit->GetId()) > fEmcClusteringThreshold  ) || 
         ( IsInCpv (digit) && Calibrate(digit->GetAmp(),digit->GetId()) > fCpvClusteringThreshold  ) ) {
-         
       Int_t iDigitInCluster = 0 ; 
       
       if  ( IsInEmc(digit) ) {   
-       // start a new EMC RecPoint
-       if(fNumberOfEmcClusters >= emcRecPoints->GetSize()) 
-         emcRecPoints->Expand(2*fNumberOfEmcClusters+1) ;
-       
-       emcRecPoints->AddAt(new  AliPHOSEmcRecPoint(), fNumberOfEmcClusters) ;
-       clu = (AliPHOSEmcRecPoint *) emcRecPoints->At(fNumberOfEmcClusters) ; 
-       fNumberOfEmcClusters++ ; 
-       clu->AddDigit(*digit, Calibrate(digit->GetAmp(),digit->GetId())) ; 
-       clusterdigitslist[iDigitInCluster] = digit->GetIndexInList() ;  
-       iDigitInCluster++ ; 
-       digitsC->Remove(digit) ; 
+        // start a new EMC RecPoint
+        if(fNumberOfEmcClusters >= emcRecPoints->GetSize()) 
+          emcRecPoints->Expand(2*fNumberOfEmcClusters+1) ;
+          
+        emcRecPoints->AddAt(new  AliPHOSEmcRecPoint(""), fNumberOfEmcClusters) ;
+        clu = dynamic_cast<AliPHOSEmcRecPoint *>( emcRecPoints->At(fNumberOfEmcClusters) ) ; 
+          fNumberOfEmcClusters++ ; 
+        clu->AddDigit(*digit, Calibrate(digit->GetAmp(),digit->GetId())) ; 
+        clusterdigitslist[iDigitInCluster] = digit->GetIndexInList() ;        
+        iDigitInCluster++ ; 
+        digitsC->Remove(digit) ; 
 
       } else { 
-       
-       // start a new CPV cluster
-       if(fNumberOfCpvClusters >= cpvRecPoints->GetSize()) 
-         cpvRecPoints->Expand(2*fNumberOfCpvClusters+1);
-
-       cpvRecPoints->AddAt(new AliPHOSCpvRecPoint(), fNumberOfCpvClusters) ;
-
-       clu =  (AliPHOSCpvRecPoint *) cpvRecPoints->At(fNumberOfCpvClusters)  ;  
-       fNumberOfCpvClusters++ ; 
-       clu->AddDigit(*digit, Calibrate(digit->GetAmp(),digit->GetId()) ) ;     
-       clusterdigitslist[iDigitInCluster] = digit->GetIndexInList()  ; 
-       iDigitInCluster++ ; 
-       digitsC->Remove(digit) ; 
+        
+        // start a new CPV cluster
+        if(fNumberOfCpvClusters >= cpvRecPoints->GetSize()) 
+          cpvRecPoints->Expand(2*fNumberOfCpvClusters+1);
+
+        cpvRecPoints->AddAt(new AliPHOSCpvRecPoint(""), fNumberOfCpvClusters) ;
+
+        clu =  dynamic_cast<AliPHOSCpvRecPoint *>( cpvRecPoints->At(fNumberOfCpvClusters) ) ;  
+        fNumberOfCpvClusters++ ; 
+        clu->AddDigit(*digit, Calibrate(digit->GetAmp(),digit->GetId()) ) ;        
+        clusterdigitslist[iDigitInCluster] = digit->GetIndexInList()  ;        
+        iDigitInCluster++ ; 
+        digitsC->Remove(digit) ; 
         nextdigit.Reset() ;
-       
-       // Here we remove remaining EMC digits, which cannot make a cluster
-       
+        
+        // Here we remove remaining EMC digits, which cannot make a cluster
+        
         if( notremoved ) { 
-         while( ( digit = (AliPHOSDigit *)nextdigit() ) ) {
+          while( ( digit = dynamic_cast<AliPHOSDigit *>( nextdigit() ) ) ) {
             if( IsInEmc(digit) ) 
-             digitsC->Remove(digit) ;
+              digitsC->Remove(digit) ;
             else 
-             break ;
-         }
-         notremoved = kFALSE ;
-       }
-       
+              break ;
+          }
+          notremoved = kFALSE ;
+        }
+        
       } // else        
       
       nextdigit.Reset() ;
@@ -637,28 +551,28 @@ void AliPHOSClusterizerv1::MakeClusters()
       AliPHOSDigit * digitN ; 
       index = 0 ;
       while (index < iDigitInCluster){ // scan over digits already in cluster 
-       digit =  (AliPHOSDigit*)digits->At(clusterdigitslist[index])  ;      
-       index++ ; 
-        while ( (digitN = (AliPHOSDigit *)nextdigit()) ) { // scan over the reduced list of digits 
-         Int_t ineb = AreNeighbours(digit, digitN);       // call (digit,digitN) in THAT oder !!!!!
+        digit =  dynamic_cast<AliPHOSDigit*>( digits->At(clusterdigitslist[index]) )  ;      
+        index++ ; 
+        while ( (digitN = dynamic_cast<AliPHOSDigit *>( nextdigit() ) ) ) { // scan over the reduced list of digits 
+          Int_t ineb = AreNeighbours(digit, digitN);       // call (digit,digitN) in THAT oder !!!!!
           switch (ineb ) {
           case 0 :   // not a neighbour
-           break ;
-         case 1 :   // are neighbours 
-           clu->AddDigit(*digitN, Calibrate( digitN->GetAmp(), digitN->GetId() ) ) ;
-           clusterdigitslist[iDigitInCluster] = digitN->GetIndexInList() ; 
-           iDigitInCluster++ ; 
-           digitsC->Remove(digitN) ;
-           break ;
+            break ;
+          case 1 :   // are neighbours 
+            clu->AddDigit(*digitN, Calibrate( digitN->GetAmp(), digitN->GetId() ) ) ;
+            clusterdigitslist[iDigitInCluster] = digitN->GetIndexInList() ; 
+            iDigitInCluster++ ; 
+            digitsC->Remove(digitN) ;
+            break ;
           case 2 :   // too far from each other
-           goto endofloop;   
-         } // switch
-         
-       } // while digitN
-       
+            goto endofloop;   
+          } // switch
+          
+        } // while digitN
+        
       endofloop: ;
-       nextdigit.Reset() ; 
-       
+        nextdigit.Reset() ; 
+        
       } // loop over cluster     
 
     } // energy theshold  
@@ -676,9 +590,10 @@ void AliPHOSClusterizerv1::MakeUnfolding()
   // Unfolds clusters using the shape of an ElectroMagnetic shower
   // Performs unfolding of all EMC/CPV clusters
 
-  AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ; 
-  
+  AliPHOSGetter * gime = AliPHOSGetter::Instance();
+
   const AliPHOSGeometry * geom = gime->PHOSGeometry() ;
+
   TObjArray * emcRecPoints = gime->EmcRecPoints() ; 
   TObjArray * cpvRecPoints = gime->CpvRecPoints() ; 
   TClonesArray * digits = gime->Digits() ; 
@@ -692,22 +607,22 @@ void AliPHOSClusterizerv1::MakeUnfolding()
     Int_t index ;   
     for(index = 0 ; index < numberofNotUnfolded ; index++){
       
-      AliPHOSEmcRecPoint * emcRecPoint = (AliPHOSEmcRecPoint *) emcRecPoints->At(index) ;
+      AliPHOSEmcRecPoint * emcRecPoint = dynamic_cast<AliPHOSEmcRecPoint *>( emcRecPoints->At(index) ) ;
       if(emcRecPoint->GetPHOSMod()> nModulesToUnfold)
-       break ;
+        break ;
       
       Int_t nMultipl = emcRecPoint->GetMultiplicity() ; 
-      Int_t * maxAt = new Int_t[nMultipl] ;
+      AliPHOSDigit ** maxAt = new AliPHOSDigit*[nMultipl] ;
       Float_t * maxAtEnergy = new Float_t[nMultipl] ;
       Int_t nMax = emcRecPoint->GetNumberOfLocalMax(maxAt, maxAtEnergy,fEmcLocMaxCut,digits) ;
       
       if( nMax > 1 ) {     // if cluster is very flat (no pronounced maximum) then nMax = 0       
-       UnfoldCluster(emcRecPoint, nMax, maxAt, maxAtEnergy) ;
-       emcRecPoints->Remove(emcRecPoint); 
-       emcRecPoints->Compress() ;
-       index-- ;
-       fNumberOfEmcClusters -- ;
-       numberofNotUnfolded-- ;
+        UnfoldCluster(emcRecPoint, nMax, maxAt, maxAtEnergy) ;
+        emcRecPoints->Remove(emcRecPoint); 
+        emcRecPoints->Compress() ;
+        index-- ;
+        fNumberOfEmcClusters -- ;
+        numberofNotUnfolded-- ;
       }
       
       delete[] maxAt ; 
@@ -726,25 +641,25 @@ void AliPHOSClusterizerv1::MakeUnfolding()
     Int_t index ;   
     for(index = 0 ; index < numberofCpvNotUnfolded ; index++){
       
-      AliPHOSRecPoint * recPoint = (AliPHOSRecPoint *) cpvRecPoints->At(index) ;
+      AliPHOSRecPoint * recPoint = dynamic_cast<AliPHOSRecPoint *>( cpvRecPoints->At(index) ) ;
 
       if(recPoint->GetPHOSMod()> nModulesToUnfold)
-       break ;
+        break ;
       
-      AliPHOSEmcRecPoint * emcRecPoint = (AliPHOSEmcRecPoint*) recPoint ; 
+      AliPHOSEmcRecPoint * emcRecPoint = dynamic_cast<AliPHOSEmcRecPoint*>(recPoint) ; 
       
       Int_t nMultipl = emcRecPoint->GetMultiplicity() ; 
-      Int_t * maxAt = new Int_t[nMultipl] ;
+      AliPHOSDigit ** maxAt = new AliPHOSDigit*[nMultipl] ;
       Float_t * maxAtEnergy = new Float_t[nMultipl] ;
       Int_t nMax = emcRecPoint->GetNumberOfLocalMax(maxAt, maxAtEnergy,fCpvLocMaxCut,digits) ;
       
       if( nMax > 1 ) {     // if cluster is very flat (no pronounced maximum) then nMax = 0       
-       UnfoldCluster(emcRecPoint, nMax, maxAt, maxAtEnergy) ;
-       cpvRecPoints->Remove(emcRecPoint); 
-       cpvRecPoints->Compress() ;
-       index-- ;
-       numberofCpvNotUnfolded-- ;
-       fNumberOfCpvClusters-- ;
+        UnfoldCluster(emcRecPoint, nMax, maxAt, maxAtEnergy) ;
+        cpvRecPoints->Remove(emcRecPoint); 
+        cpvRecPoints->Compress() ;
+        index-- ;
+        numberofCpvNotUnfolded-- ;
+        fNumberOfCpvClusters-- ;
       }
       
       delete[] maxAt ; 
@@ -769,15 +684,17 @@ Double_t  AliPHOSClusterizerv1::ShowerShape(Double_t r)
 
 //____________________________________________________________________________
 void  AliPHOSClusterizerv1::UnfoldCluster(AliPHOSEmcRecPoint * iniEmc, 
-                                                Int_t nMax, 
-                                                int * maxAt, 
-                                                Float_t * maxAtEnergy)
+                                          Int_t nMax, 
+                                          AliPHOSDigit ** maxAt, 
+                                          Float_t * maxAtEnergy)
 { 
   // Performs the unfolding of a cluster with nMax overlapping showers 
 
-  AliPHOSGetter * gime = AliPHOSGetter::GetInstance() ; 
+  AliPHOSGetter * gime = AliPHOSGetter::Instance();
+
   const AliPHOSGeometry * geom = gime->PHOSGeometry() ;
-  TClonesArray * digits = gime->Digits() ; 
+
+  const TClonesArray * digits = gime->Digits() ; 
   TObjArray * emcRecPoints = gime->EmcRecPoints() ; 
   TObjArray * cpvRecPoints = gime->CpvRecPoints() ; 
 
@@ -806,7 +723,7 @@ void  AliPHOSClusterizerv1::UnfoldCluster(AliPHOSEmcRecPoint * iniEmc,
   Int_t iparam ;
   Int_t iDigit ;
   for(iDigit = 0 ; iDigit < nDigits ; iDigit ++){
-    digit = (AliPHOSDigit*) digits->At(emcDigits[iDigit] ) ;   
+    digit = dynamic_cast<AliPHOSDigit*>( digits->At(emcDigits[iDigit] ) ) ;   
     geom->AbsToRelNumbering(digit->GetId(), relid) ;
     geom->RelPosInModule(relid, xDigit, zDigit) ;
     efit[iDigit] = 0;
@@ -843,24 +760,24 @@ void  AliPHOSClusterizerv1::UnfoldCluster(AliPHOSEmcRecPoint * iniEmc,
     if(iniEmc->IsEmc()){ //create new entries in fEmcRecPoints...
       
       if(fNumberOfEmcClusters >= emcRecPoints->GetSize())
-       emcRecPoints->Expand(2*fNumberOfEmcClusters) ;
+        emcRecPoints->Expand(2*fNumberOfEmcClusters) ;
       
-      (*emcRecPoints)[fNumberOfEmcClusters] = new AliPHOSEmcRecPoint() ;
-      emcRP = (AliPHOSEmcRecPoint *) emcRecPoints->At(fNumberOfEmcClusters);
+      (*emcRecPoints)[fNumberOfEmcClusters] = new AliPHOSEmcRecPoint("") ;
+      emcRP = dynamic_cast<AliPHOSEmcRecPoint *>( emcRecPoints->At(fNumberOfEmcClusters) ) ;
       fNumberOfEmcClusters++ ;
     }
     else{//create new entries in fCpvRecPoints
       if(fNumberOfCpvClusters >= cpvRecPoints->GetSize())
-       cpvRecPoints->Expand(2*fNumberOfCpvClusters) ;
+        cpvRecPoints->Expand(2*fNumberOfCpvClusters) ;
       
-      (*cpvRecPoints)[fNumberOfCpvClusters] = new AliPHOSCpvRecPoint() ;
-      emcRP = (AliPHOSEmcRecPoint *) cpvRecPoints->At(fNumberOfCpvClusters);
+      (*cpvRecPoints)[fNumberOfCpvClusters] = new AliPHOSCpvRecPoint("") ;
+      emcRP = dynamic_cast<AliPHOSEmcRecPoint *>( cpvRecPoints->At(fNumberOfCpvClusters) ) ;
       fNumberOfCpvClusters++ ;
     }
     
     Float_t eDigit ;
     for(iDigit = 0 ; iDigit < nDigits ; iDigit ++){
-      digit = (AliPHOSDigit*) digits->At( emcDigits[iDigit] ) ; 
+      digit = dynamic_cast<AliPHOSDigit*>( digits->At( emcDigits[iDigit] ) ) ; 
       geom->AbsToRelNumbering(digit->GetId(), relid) ;
       geom->RelPosInModule(relid, xDigit, zDigit) ;
       distance = (xDigit - xpar) * (xDigit - xpar) + (zDigit - zpar) * (zDigit - zpar)  ;
@@ -868,7 +785,7 @@ void  AliPHOSClusterizerv1::UnfoldCluster(AliPHOSEmcRecPoint * iniEmc,
       ratio = epar * ShowerShape(distance) / efit[iDigit] ; 
       eDigit = emcEnergies[iDigit] * ratio ;
       emcRP->AddDigit( *digit, eDigit ) ;
-    }  
+    }        
   }
  
   delete[] fitparameters ; 
@@ -882,23 +799,22 @@ void AliPHOSClusterizerv1::UnfoldingChiSquare(Int_t & nPar, Double_t * Grad, Dou
   // Calculates the Chi square for the cluster unfolding minimization
   // Number of parameters, Gradient, Chi squared, parameters, what to do
 
+  TList * toMinuit = dynamic_cast<TList*>( gMinuit->GetObjectFit() ) ;
 
-  TList * toMinuit = (TList*) gMinuit->GetObjectFit() ;
-
-  AliPHOSEmcRecPoint * emcRP = (AliPHOSEmcRecPoint*) toMinuit->At(0)  ;
-  TClonesArray * digits = (TClonesArray*)toMinuit->At(1)  ;
+  AliPHOSEmcRecPoint * emcRP = dynamic_cast<AliPHOSEmcRecPoint*>( toMinuit->At(0) )  ;
+  TClonesArray * digits = dynamic_cast<TClonesArray*>( toMinuit->At(1) )  ;
 
 
   
-  //  AliPHOSEmcRecPoint * emcRP = (AliPHOSEmcRecPoint *) gMinuit->GetObjectFit() ; // EmcRecPoint to fit
+  //  AliPHOSEmcRecPoint * emcRP = dynamic_cast<AliPHOSEmcRecPoint *>( gMinuit->GetObjectFit() ) ; // EmcRecPoint to fit
 
   Int_t * emcDigits     = emcRP->GetDigitsList() ;
 
   Int_t nOdigits = emcRP->GetDigitsMultiplicity() ; 
 
   Float_t * emcEnergies = emcRP->GetEnergiesList() ;
-
-  const AliPHOSGeometry * geom = AliPHOSGetter::GetInstance()->PHOSGeometry() ; 
+  
+  const AliPHOSGeometry * geom = AliPHOSGetter::Instance()->PHOSGeometry() ; 
   fret = 0. ;     
   Int_t iparam ;
 
@@ -913,7 +829,7 @@ void AliPHOSClusterizerv1::UnfoldingChiSquare(Int_t & nPar, Double_t * Grad, Dou
 
   for( iDigit = 0 ; iDigit < nOdigits ; iDigit++) {
 
-    digit = (AliPHOSDigit*) digits->At( emcDigits[iDigit] ) 
+    digit = dynamic_cast<AliPHOSDigit*>( digits->At( emcDigits[iDigit] ) )
 
     Int_t relid[4] ;
     Float_t xDigit ;
@@ -927,33 +843,33 @@ void AliPHOSClusterizerv1::UnfoldingChiSquare(Int_t & nPar, Double_t * Grad, Dou
        Int_t iParam = 0 ;
        efit = 0 ;
        while(iParam < nPar ){
-        Double_t distance = (xDigit - x[iParam]) * (xDigit - x[iParam]) ;
-        iParam++ ; 
-        distance += (zDigit - x[iParam]) * (zDigit - x[iParam]) ; 
-        distance = TMath::Sqrt( distance ) ; 
-        iParam++ ;      
-        efit += x[iParam] * ShowerShape(distance) ;
-        iParam++ ;
+         Double_t distance = (xDigit - x[iParam]) * (xDigit - x[iParam]) ;
+         iParam++ ; 
+         distance += (zDigit - x[iParam]) * (zDigit - x[iParam]) ; 
+         distance = TMath::Sqrt( distance ) ; 
+         iParam++ ;          
+         efit += x[iParam] * ShowerShape(distance) ;
+         iParam++ ;
        }
        Double_t sum = 2. * (efit - emcEnergies[iDigit]) / emcEnergies[iDigit] ; // Here we assume, that sigma = sqrt(E) 
        iParam = 0 ;
        while(iParam < nPar ){
-        Double_t xpar = x[iParam] ;
-        Double_t zpar = x[iParam+1] ;
-        Double_t epar = x[iParam+2] ;
-        Double_t dr = TMath::Sqrt( (xDigit - xpar) * (xDigit - xpar) + (zDigit - zpar) * (zDigit - zpar) );
-        Double_t shape = sum * ShowerShape(dr) ;
-        Double_t r4 = dr*dr*dr*dr ;
-        Double_t r295 = TMath::Power(dr,2.95) ;
-        Double_t deriv =-4. * dr*dr * ( 2.32 / ( (2.32 + 0.26 * r4) * (2.32 + 0.26 * r4) ) +
-                                        0.0316 * (1. + 0.0171 * r295) / ( ( 1. + 0.0652 * r295) * (1. + 0.0652 * r295) ) ) ;
-        
-        Grad[iParam] += epar * shape * deriv * (xpar - xDigit) ;  // Derivative over x    
-        iParam++ ; 
-        Grad[iParam] += epar * shape * deriv * (zpar - zDigit) ;  // Derivative over z         
-        iParam++ ; 
-        Grad[iParam] += shape ;                                  // Derivative over energy             
-        iParam++ ; 
+         Double_t xpar = x[iParam] ;
+         Double_t zpar = x[iParam+1] ;
+         Double_t epar = x[iParam+2] ;
+         Double_t dr = TMath::Sqrt( (xDigit - xpar) * (xDigit - xpar) + (zDigit - zpar) * (zDigit - zpar) );
+         Double_t shape = sum * ShowerShape(dr) ;
+         Double_t r4 = dr*dr*dr*dr ;
+         Double_t r295 = TMath::Power(dr,2.95) ;
+         Double_t deriv =-4. * dr*dr * ( 2.32 / ( (2.32 + 0.26 * r4) * (2.32 + 0.26 * r4) ) +
+                                         0.0316 * (1. + 0.0171 * r295) / ( ( 1. + 0.0652 * r295) * (1. + 0.0652 * r295) ) ) ;
+         
+         Grad[iParam] += epar * shape * deriv * (xpar - xDigit) ;  // Derivative over x    
+         iParam++ ; 
+         Grad[iParam] += epar * shape * deriv * (zpar - zDigit) ;  // Derivative over z         
+         iParam++ ; 
+         Grad[iParam] += shape ;                                  // Derivative over energy             
+         iParam++ ; 
        }
      }
      efit = 0;
@@ -976,58 +892,70 @@ void AliPHOSClusterizerv1::UnfoldingChiSquare(Int_t & nPar, Double_t * Grad, Dou
 }
 
 //____________________________________________________________________________
-void AliPHOSClusterizerv1::Print(Option_t * option)const
+void AliPHOSClusterizerv1::Print()const
 {
   // Print clusterizer parameters
 
-  if( strcmp(GetName(), "") !=0 ){
-    
+  TString message ; 
+  TString taskName(GetName()) ; 
+  taskName.ReplaceAll(Version(), "") ;
+  
+  if( strcmp(GetName(), "") !=0 ) {  
     // Print parameters
-    TString taskName(GetName()) ; 
-    taskName.ReplaceAll(Version(), "") ;
-
-    cout << "---------------"<< taskName.Data() << " " << GetTitle()<< "-----------" << endl 
-        << "Clusterizing digits from the file: " << fHeaderFileName.Data() << endl 
-        << "                           Branch: " << fDigitsBranchTitle.Data() << endl 
-        << endl 
-        << "                       EMC Clustering threshold = " << fEmcClusteringThreshold << endl
-        << "                       EMC Local Maximum cut    = " << fEmcLocMaxCut << endl
-        << "                       EMC Logarothmic weight   = " << fW0 << endl
-        << endl
-        << "                       CPV Clustering threshold = " << fCpvClusteringThreshold << endl
-        << "                       CPV Local Maximum cut    = " << fCpvLocMaxCut << endl
-       << "                       CPV Logarothmic weight   = " << fW0CPV << endl
-        << endl ;
+    message  = "\n--------------- %s %s -----------\n" ; 
+    message += "Clusterizing digits from the file: %s\n" ;
+    message += "                           Branch: %s\n" ; 
+    message += "                       EMC Clustering threshold = %f\n" ; 
+    message += "                       EMC Local Maximum cut    = %f\n" ; 
+    message += "                       EMC Logarothmic weight   = %f\n" ;
+    message += "                       CPV Clustering threshold = %f\n" ;
+    message += "                       CPV Local Maximum cut    = %f\n" ;
+    message += "                       CPV Logarothmic weight   = %f\n" ;
     if(fToUnfold)
-      cout << " Unfolding on " << endl ;
+      message += " Unfolding on\n" ;
     else
-      cout << " Unfolding off " << endl ;
+      message += " Unfolding off\n" ;
     
-    cout << "------------------------------------------------------------------" <<endl ;
+    message += "------------------------------------------------------------------" ;
   }
-  else
-    cout << " AliPHOSClusterizerv1 not initialized " << endl ;
+  else 
+    message = " AliPHOSClusterizerv1 not initialized " ;
+  
+  Info("Print", message.Data(),  
+       taskName.Data(), 
+       GetTitle(),
+       taskName.Data(), 
+       GetName(), 
+       fEmcClusteringThreshold, 
+       fEmcLocMaxCut, 
+       fW0, 
+       fCpvClusteringThreshold, 
+       fCpvLocMaxCut, 
+       fW0CPV ) ; 
 }
+
+
 //____________________________________________________________________________
 void AliPHOSClusterizerv1::PrintRecPoints(Option_t * option)
 {
   // Prints list of RecPoints produced at the current pass of AliPHOSClusterizer
 
-  TObjArray * emcRecPoints = AliPHOSGetter::GetInstance()->EmcRecPoints() ; 
-  TObjArray * cpvRecPoints = AliPHOSGetter::GetInstance()->CpvRecPoints() ; 
-
-  cout << "AliPHOSClusterizerv1: : event "<<gAlice->GetEvNumber() << endl ;
-  cout << "       Found "<< emcRecPoints->GetEntriesFast() << " EMC Rec Points and " 
-          << cpvRecPoints->GetEntriesFast() << " CPV RecPoints" << endl ;
+  AliPHOSGetter * gime = AliPHOSGetter::Instance();
+  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())  ; 
   fRecPointsInRun +=  emcRecPoints->GetEntriesFast() ; 
   fRecPointsInRun +=  cpvRecPoints->GetEntriesFast() ; 
-
+  
+  
   if(strstr(option,"all")) {
-    cout << "EMC clusters " << endl ;
-    cout << " Index  Ene(MeV)   Multi  Module     X      Y      Z    Lambda 1   Lambda 2  # of prim  Primaries list "      <<  endl;      
-    
+    printf("\n EMC clusters \n") ;
+    printf("Index    Ene(MeV) Multi Module     X    Y   Z    Lambdas_1  Lambda_2  # of prim  Primaries list\n") ;      
     Int_t index ;
     for (index = 0 ; index < emcRecPoints->GetEntries() ; index++) {
       AliPHOSEmcRecPoint * rp = (AliPHOSEmcRecPoint * )emcRecPoints->At(index) ; 
@@ -1038,57 +966,29 @@ void AliPHOSClusterizerv1::PrintRecPoints(Option_t * option)
       Int_t * primaries; 
       Int_t nprimaries;
       primaries = rp->GetPrimaries(nprimaries);
-
-      cout << setw(4) << rp->GetIndexInList() << "   " 
-          << setw(7) << setprecision(3) << rp->GetEnergy() << "           " 
-          << setw(3) <<         rp->GetMultiplicity() << "  " 
-          << setw(1) <<              rp->GetPHOSMod() << "  " 
-          << setw(6) << setprecision(2) << locpos.X() << "  " 
-          << setw(6) << setprecision(2) << locpos.Y() << "  " 
-          << setw(6) << setprecision(2) << locpos.Z() << "  "
-          << setw(4) << setprecision(2) << lambda[0]  << "  "
-          << setw(4) << setprecision(2) << lambda[1]  << "  "
-          << setw(2) << nprimaries << "  " ;
-     
-      for (Int_t iprimary=0; iprimary<nprimaries; iprimary++)
-       cout << setw(4) <<   primaries[iprimary] << "  "  ;
-      cout << endl ;    
+      printf("\n%6d  %8.2f  %3d     %2d     %4.1f   %4.1f   %4.1f   %4f  %4f    %2d     : ", 
+             rp->GetIndexInList(), rp->GetEnergy(), rp->GetMultiplicity(), rp->GetPHOSMod(), 
+             locpos.X(), locpos.Y(), locpos.Z(), lambda[0], lambda[1], nprimaries) ; 
+      
+      for (Int_t iprimary=0; iprimary<nprimaries; iprimary++) {
+       printf("%d ", primaries[iprimary] ) ; 
+      }
+      printf("\n") ;
     }
 
     //Now plot CPV recPoints
-    cout << "EMC clusters " << endl ;
-    cout << "  Index    " 
-        << "  Multi    "
-        << "  Module   "  
-        << "    X      "
-        << "    Y      "
-        << "    Z      "
-        << " # of prim "
-        << " Primaries list "      <<  endl;      
-    
+    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) ; 
-      cout << setw(6) << rp->GetIndexInList() << "     ";
-      cout << setw(6) << rp->GetPHOSMod()     << "        CPV     ";
+      AliPHOSCpvRecPoint * rp = (AliPHOSCpvRecPoint * )cpvRecPoints->At(index) ; 
       
       TVector3  locpos;  
       rp->GetLocalPosition(locpos);
-      cout << setw(6) <<  locpos.X()          << "     ";
-      cout << setw(6) <<  locpos.Y()          << "     ";
-      cout << setw(6) <<  locpos.Z()          << "     ";
       
-      Int_t * primaries; 
-      Int_t nprimaries ; 
-      primaries = rp->GetPrimaries(nprimaries);
-      cout << setw(6) <<    nprimaries         << "     ";
-
-      for (Int_t iprimary=0; iprimary<nprimaries; iprimary++)
-       cout << setw(4)  <<  primaries[iprimary] << " ";
-      cout << endl;     
+      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()) ; 
     }
-
-
-    cout << "-----------------------------------------------------------------------"<<endl ;
   }
 }