AOD CaloCluster in PHOS.

diff --git a/PHOS/AliPHOSAodCluster.cxx b/PHOS/AliPHOSAodCluster.cxx
new file mode 100755 (executable)
index 0000000..9bd0fc5
--- /dev/null
+++ b/PHOS/AliPHOSAodCluster.cxx
@@ -0,0 +1,236 @@
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+//_________________________________________________________________________
+//  AliAODCaloCluster extension for PHOS to recalculate cluster 
+//  parameters in case of recalibration.
+//*--
+//*-- Author: Dmitri Peressounko (RRC KI)
+
+
+// --- ROOT system ---
+#include "TVector3.h"
+#include "TMath.h"
+
+// --- Standard library ---
+
+// --- AliRoot header files ---
+#include "AliLog.h" 
+#include "AliPHOSGeometry.h" 
+#include "AliPHOSPIDv1.h" 
+#include "AliPHOSAodCluster.h" 
+#include "AliPHOSCalibData.h"
+#include "AliAODCaloCells.h"
+
+ClassImp(AliPHOSAodCluster)
+
+//____________________________________________________________________________
+AliPHOSAodCluster::AliPHOSAodCluster() : 
+  AliAODCaloCluster(),fRecalibrated(0)
+{
+  // ctor
+}
+//____________________________________________________________________________
+AliPHOSAodCluster::AliPHOSAodCluster(const AliAODCaloCluster & clu) : 
+  AliAODCaloCluster(clu),fRecalibrated(0)
+{
+  // cpy ctor
+}
+
+//____________________________________________________________________________
+AliPHOSAodCluster::~AliPHOSAodCluster()
+{
+  // dtor
+}
+//____________________________________________________________________________
+void AliPHOSAodCluster::Recalibrate(AliPHOSCalibData * calibData,AliAODCaloCells *phsCells){
+  //If not done yet, apply recalibration coefficients to energies list
+  //NOTE that after recalibration fCellsAmpFraction contains not FRACTION but FULL energy 
+  
+  if(fRecalibrated)
+   return ;
+  
+  if(!calibData)
+    return ;
+
+  AliPHOSGeometry * phosgeom =  AliPHOSGeometry::GetInstance() ;
+  if(!phosgeom)
+    AliFatal("AliPHOSGeometry was not contructed\n") ;
+       
+  Double32_t * cellsAmpFraction = GetCellsAmplitudeFraction(); 
+       
+  for(Int_t i=0; i < GetNCells(); i++){
+    Int_t relId[4];
+    phosgeom->AbsToRelNumbering(GetCellAbsId(i),relId) ;
+    Int_t   module = relId[0];
+    Int_t   column = relId[3];
+    Int_t   row    = relId[2];
+    Double_t energy = phsCells->GetCellAmplitude(GetCellAbsId(i)) ;
+    cellsAmpFraction[i]*=energy*calibData->GetADCchannelEmc(module,column,row);
+  }
+       
+  SetCellsAmplitudeFraction(cellsAmpFraction);
+  fRecalibrated=kTRUE; 
+}
+//____________________________________________________________________________
+void  AliPHOSAodCluster::EvalAll(Float_t logWeight, TVector3 &vtx){
+    //If recalibrated - recalculate all cluster parameters
+  if(!fRecalibrated)
+    return ;
+
+  EvalEnergy() ; //Energy should be evaluated first
+  EvalCoord(logWeight, vtx) ;
+  
+}
+//____________________________________________________________________________
+void AliPHOSAodCluster::EvalEnergy(){
+  if(!fRecalibrated) // no need to recalibrate
+    return ;
+    
+  Float_t energy=0. ;
+  for(Int_t iDigit=0; iDigit < GetNCells(); iDigit++) {
+    energy+=GetCellAmplitudeFraction(iDigit) ;
+  }
+  
+  SetE(energy);
+  //Correct for nonlinearity later
+   
+}
+//____________________________________________________________________________
+void AliPHOSAodCluster::EnergyCorrection(AliPHOSPIDv1 * pid){
+  //apply nonlinearity correction same as in AliPHOSPIDv1.
+  SetE(pid->GetCalibratedEnergy(E())) ;
+}
+//____________________________________________________________________________
+void AliPHOSAodCluster::EvalPID(AliPHOSPIDv1 * /*pid*/){           
+  //re-evaluate identification parameters
+//  pid->CalculatePID(E(),GetDispersion(),GetEmcCpvDistance(),GetTOF(),GetPID()) ;  
+//  pid->CalculatePID(E(),GetDispersion(),GetM20(),GetM02(),GetEmcCpvDistance(),GetTOF(),GetPID()) ;
+}
+//____________________________________________________________________________
+void AliPHOSAodCluster::EvalCoord(Float_t logWeight, TVector3 &vtx)
+{
+  // Calculates new center of gravity in the local PHOS-module coordinates 
+  // and tranfers into global ALICE coordinates
+  // Calculates Dispersion and main axis
+  if(!fRecalibrated) // no need to recalibrate
+    return ;
+ 
+  Float_t wtot = 0. ;
+  Int_t relid[4] ;
+  Int_t phosMod=0 ;
+  Float_t xMean = 0. ;
+  Float_t zMean = 0. ;
+
+  AliPHOSGeometry * phosgeom =  AliPHOSGeometry::GetInstance() ;
+  if(!phosgeom)
+    AliFatal("AliPHOSGeometry was not contructed\n") ;
+
+  for(Int_t iDigit=0; iDigit < GetNCells(); iDigit++) {
+    Float_t xi ;
+    Float_t zi ;
+    phosgeom->AbsToRelNumbering(GetCellAbsId(iDigit), relid) ;
+    phosgeom->RelPosInModule(relid, xi, zi);
+    phosMod=relid[0] ;
+    Double_t ei=GetCellAmplitudeFraction(iDigit) ;
+    if (E() > 0 && ei > 0) {
+      Float_t w = TMath::Max( 0., logWeight + TMath::Log(ei/E()) ) ;
+      xMean+= xi * w ;
+      zMean+= zi * w ;
+      wtot += w ;
+    }
+    else
+      AliError(Form("Wrong energy %f and/or amplitude %f\n", ei, E()));
+  }
+  if (wtot>0) {
+    xMean /= wtot ;
+    zMean /= wtot ;
+  }
+  else
+    AliError(Form("Wrong weight %f\n", wtot));
+
+
+// Calculates the dispersion and second momenta
+  Double_t d=0. ;
+  Double_t dxx  = 0.;
+  Double_t dzz  = 0.;
+  Double_t dxz  = 0.;
+  for(Int_t iDigit=0; iDigit < GetNCells(); iDigit++) {
+    Float_t xi ;
+    Float_t zi ;
+    phosgeom->AbsToRelNumbering(GetCellAbsId(iDigit), relid) ;
+    phosgeom->RelPosInModule(relid, xi, zi);
+    Double_t ei=GetCellAmplitudeFraction(iDigit) ;
+    if (E() > 0 && ei > 0) {
+      Float_t w = TMath::Max( 0., logWeight + TMath::Log(ei/E()) ) ;
+      d += w*((xi-xMean)*(xi-xMean) + (zi-zMean)*(zi-zMean) ) ; 
+      dxx  += w * xi * xi ;
+      dzz  += w * zi * zi ;
+      dxz  += w * xi * zi ; 
+   }
+    else
+      AliError(Form("Wrong energy %f and/or amplitude %f\n", ei, E()));
+  }
+  
+  if (wtot>0) {
+    d /= wtot ;
+    dxx /= wtot ;
+    dzz /= wtot ;
+    dxz /= wtot ;
+    dxx -= xMean * xMean ;
+    dzz -= zMean * zMean ;
+    dxz -= xMean * zMean ;
+    SetM02(0.5 * (dxx + dzz) + TMath::Sqrt( 0.25 * (dxx - dzz) * (dxx - dzz) + dxz * dxz )) ;
+    SetM20(0.5 * (dxx + dzz) - TMath::Sqrt( 0.25 * (dxx - dzz) * (dxx - dzz) + dxz * dxz ));
+  }
+  else{
+    AliError(Form("Wrong weight %f\n", wtot));
+    d=0. ;
+    SetM20(0.) ;
+    SetM02(0.) ;
+  }
+
+  if (d>=0)
+    SetDispersion(TMath::Sqrt(d)) ;
+  else    
+    SetDispersion(0) ;
+
+
+  // Correction for the depth of the shower starting point (TDR p 127)  
+  Float_t para = 0.925 ; 
+  Float_t parb = 6.52 ; 
+
+  TVector3 vInc ;
+  phosgeom->GetIncidentVector(vtx,phosMod,xMean,zMean,vInc) ;
+
+  Float_t depthx = 0.; 
+  Float_t depthz = 0.;
+  if (E()>0&&vInc.Y()!=0.) {
+    depthx = ( para * TMath::Log(E()) + parb ) * vInc.X()/TMath::Abs(vInc.Y()) ;
+    depthz = ( para * TMath::Log(E()) + parb ) * vInc.Z()/TMath::Abs(vInc.Y()) ;
+  }
+  else 
+    AliError(Form("Wrong amplitude %f\n", E()));
+
+  xMean-= depthx  ;
+  zMean-= depthz  ;
+
+  //Go to the global system
+  TVector3 gps ;
+  phosgeom->Local2Global(phosMod, xMean, zMean, gps) ;
+  SetPosition(0, gps[0]) ;
+  SetPosition(1, gps[1]) ;  
+  SetPosition(2, gps[2]) ;
+}

diff --git a/PHOS/AliPHOSAodCluster.h b/PHOS/AliPHOSAodCluster.h
new file mode 100755 (executable)
index 0000000..dfe3cb6
--- /dev/null
+++ b/PHOS/AliPHOSAodCluster.h
@@ -0,0 +1,48 @@
+#ifndef ALIPHOSAODCLUSTER_H
+#define ALIPHOSAODCLUSTER_H
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+
+//_________________________________________________________________________
+//  AliAodCaloCluster version for PHOS (used for recalibration)
+//           
+//*-- Author: Dmitri Peressounko (RRC KI)
+
+// --- ROOT system ---
+
+// --- Standard library ---
+
+// --- AliRoot header files ---
+class AliPHOSPIDv1 ;
+class AliPHOSCalibData ;
+class AliAODCaloCells ;
+
+#include "AliAODCaloCluster.h"
+
+class AliPHOSAodCluster : public AliAODCaloCluster  {
+
+public:
+
+  AliPHOSAodCluster() ;
+  AliPHOSAodCluster(const AliAODCaloCluster & clu) ; 
+ 
+  virtual ~AliPHOSAodCluster() ;  
+
+  void  EvalAll(Float_t logWeight, TVector3 &vtx) ; //re-calculate all cluster parameters
+  void  Recalibrate(AliPHOSCalibData * calibData,AliAODCaloCells *phsCells) ; //Apply recalibration to this cluster
+  void  EnergyCorrection(AliPHOSPIDv1 * pid) ;  //Apply non-linearity correction
+  void  EvalPID(AliPHOSPIDv1 * pid) ;           //re-evaluate identification parameters
+
+protected:
+ 
+  void EvalCoord(Float_t logWeight, TVector3 &vtx) ;  //calculate coordinate-related parameters (position, dispersion)
+  void EvalEnergy() ; //re-calculate energy of the cluster
+
+  Bool_t fRecalibrated ;  //Have this cluster been recalibrated
+    
+  ClassDef(AliPHOSAodCluster,1)  // (PHOS AOD cluster)
+
+};
+
+#endif // AliPHOSAODCLUSTER_H