PHOS/AliPHOSEsdCluster.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 //_________________________________________________________________________
  17 //  AliESDCaloCluster extension for PHOS to recalculate cluster
  18 //  parameters in case of recalibration.
  19 //*--
  20 //*-- Author: Dmitri Peressounko (RRC KI)
  21
  22
  23 // --- ROOT system ---
  24 #include "TVector3.h"
  25 #include "TMath.h"
  26
  27 // --- Standard library ---
  28
  29 // --- AliRoot header files ---
  30 #include "AliLog.h"
  31 #include "AliPHOSGeometry.h"
  32 #include "AliPHOSPIDv1.h"
  33 #include "AliPHOSEsdCluster.h"
  34 #include "AliPHOSCalibData.h"
  35 #include "AliESDCaloCells.h"
  36
  37 ClassImp(AliPHOSEsdCluster)
  38
  39 //____________________________________________________________________________
  40 AliPHOSEsdCluster::AliPHOSEsdCluster() :
  41   AliESDCaloCluster(),fRecalibrated(0)
  42 {
  43   // ctor
  44 }
  45 //____________________________________________________________________________
  46 AliPHOSEsdCluster::AliPHOSEsdCluster(const AliESDCaloCluster & clu) :
  47   AliESDCaloCluster(clu),fRecalibrated(0)
  48 {
  49   // cpy ctor
  50 }
  51
  52 //____________________________________________________________________________
  53 AliPHOSEsdCluster::~AliPHOSEsdCluster()
  54 {
  55   // dtor
  56 }
  57 //____________________________________________________________________________
  58 void AliPHOSEsdCluster::Recalibrate(AliPHOSCalibData * calibData,AliESDCaloCells *phsCells){
  59   //If not done yet, apply recalibration coefficients to energies list
  60   //NOTE that after recalibration fCellsAmpFraction contains not FRACTION but FULL energy
  61
  62   if(fRecalibrated)
  63    return ;
  64
  65   if(!calibData)
  66     return ;
  67
  68   AliPHOSGeometry * phosgeom =  AliPHOSGeometry::GetInstance() ;
  69   if(!phosgeom)
  70     AliFatal("AliPHOSGeometry was not contructed\n") ;
  71
  72   for(Int_t i=0; i<fNCells; i++){
  73     Int_t relId[4];
  74     phosgeom->AbsToRelNumbering(fCellsAbsId[i],relId) ;
  75     Int_t   module = relId[0];
  76     Int_t   column = relId[3];
  77     Int_t   row    = relId[2];
  78     Double_t energy = phsCells->GetCellAmplitude(fCellsAbsId[i]) ;
  79     fCellsAmpFraction[i]*=energy*calibData->GetADCchannelEmc(module,column,row);
  80   }
  81   fRecalibrated=kTRUE;
  82 }
  83 //____________________________________________________________________________
  84 void  AliPHOSEsdCluster::EvalAll(Float_t logWeight, TVector3 &vtx){
  85     //If recalibrated - recalculate all cluster parameters
  86   if(!fRecalibrated)
  87     return ;
  88
  89   EvalEnergy() ; //Energy should be evaluated first
  90   EvalCoord(logWeight, vtx) ;
  91
  92 }
  93 //____________________________________________________________________________
  94 void AliPHOSEsdCluster::EvalEnergy(){
  95   if(!fRecalibrated) // no need to recalibrate
  96     return ;
  97
  98   fEnergy=0. ;
  99   for(Int_t iDigit=0; iDigit<fNCells; iDigit++) {
 100     fEnergy+=fCellsAmpFraction[iDigit] ;
 101   }
 102   //Correct for nonlinearity later
 103
 104 }
 105 //____________________________________________________________________________
 106 void AliPHOSEsdCluster::EnergyCorrection(AliPHOSPIDv1 * pid){
 107   //apply nonlinearity correction same as in AliPHOSPIDv1.
 108   fEnergy = pid->GetCalibratedEnergy(fEnergy) ;
 109 }
 110 //____________________________________________________________________________
 111 void AliPHOSEsdCluster::EvalPID(AliPHOSPIDv1 * /*pid*/){
 112   //re-evaluate identification parameters
 113 //  pid->CalculatePID(fEnergy,fDispersion,fEmcCpvDistance,tof,fPID) ;
 114 //  pid->CalculatePID(fEnergy,fDispersion,fM20,fM02,fEmcCpvDistance,tof,fPID) ;
 115 }
 116 //____________________________________________________________________________
 117 void AliPHOSEsdCluster::EvalCoord(Float_t logWeight, TVector3 &vtx)
 118 {
 119   // Calculates new center of gravity in the local PHOS-module coordinates
 120   // and tranfers into global ALICE coordinates
 121   // Calculates Dispersion and main axis
 122   if(!fRecalibrated) // no need to recalibrate
 123     return ;
 124
 125   Float_t wtot = 0. ;
 126   Int_t relid[4] ;
 127   Int_t phosMod=0 ;
 128   Float_t xMean = 0. ;
 129   Float_t zMean = 0. ;
 130
 131   AliPHOSGeometry * phosgeom =  AliPHOSGeometry::GetInstance() ;
 132   if(!phosgeom)
 133     AliFatal("AliPHOSGeometry was not contructed\n") ;
 134
 135   for(Int_t iDigit=0; iDigit<fNCells; iDigit++) {
 136     Float_t xi ;
 137     Float_t zi ;
 138     phosgeom->AbsToRelNumbering(fCellsAbsId[iDigit], relid) ;
 139     phosgeom->RelPosInModule(relid, xi, zi);
 140     phosMod=relid[0] ;
 141     Double_t ei=fCellsAmpFraction[iDigit] ;
 142     if (fEnergy>0 && ei>0) {
 143       Float_t w = TMath::Max( 0., logWeight + TMath::Log(ei/fEnergy) ) ;
 144       xMean+= xi * w ;
 145       zMean+= zi * w ;
 146       wtot += w ;
 147     }
 148     else
 149       AliError(Form("Wrong energy %f and/or amplitude %f\n", ei, fEnergy));
 150   }
 151   if (wtot>0) {
 152     xMean /= wtot ;
 153     zMean /= wtot ;
 154   }
 155   else
 156     AliError(Form("Wrong weight %f\n", wtot));
 157
 158
 159 // Calculates the dispersion and second momenta
 160   Double_t d=0. ;
 161   Double_t dxx  = 0.;
 162   Double_t dzz  = 0.;
 163   Double_t dxz  = 0.;
 164   for(Int_t iDigit=0; iDigit < fNCells; iDigit++) {
 165     Float_t xi ;
 166     Float_t zi ;
 167     phosgeom->AbsToRelNumbering(fCellsAbsId[iDigit], relid) ;
 168     phosgeom->RelPosInModule(relid, xi, zi);
 169     Double_t ei=fCellsAmpFraction[iDigit] ;
 170     if (fEnergy>0 && ei>0) {
 171       Float_t w = TMath::Max( 0., logWeight + TMath::Log(ei/fEnergy) ) ;
 172       d += w*((xi-xMean)*(xi-xMean) + (zi-zMean)*(zi-zMean) ) ;
 173       dxx  += w * xi * xi ;
 174       dzz  += w * zi * zi ;
 175       dxz  += w * xi * zi ;
 176    }
 177     else
 178       AliError(Form("Wrong energy %f and/or amplitude %f\n", ei, fEnergy));
 179   }
 180
 181   if (wtot>0) {
 182     d /= wtot ;
 183     dxx /= wtot ;
 184     dzz /= wtot ;
 185     dxz /= wtot ;
 186     dxx -= xMean * xMean ;
 187     dzz -= zMean * zMean ;
 188     dxz -= xMean * zMean ;
 189     fM02 =  0.5 * (dxx + dzz) + TMath::Sqrt( 0.25 * (dxx - dzz) * (dxx - dzz) + dxz * dxz )  ;
 190     fM20 =  0.5 * (dxx + dzz) - TMath::Sqrt( 0.25 * (dxx - dzz) * (dxx - dzz) + dxz * dxz )  ;
 191   }
 192   else{
 193     AliError(Form("Wrong weight %f\n", wtot));
 194     d=0. ;
 195     fM20=0. ;
 196     fM02=0. ;
 197   }
 198
 199   if (d>=0)
 200     fDispersion = TMath::Sqrt(d) ;
 201   else
 202     fDispersion = 0 ;
 203
 204
 205   // Correction for the depth of the shower starting point (TDR p 127)
 206   Float_t para = 0.925 ;
 207   Float_t parb = 6.52 ;
 208
 209   TVector3 vInc ;
 210   phosgeom->GetIncidentVector(vtx,phosMod,xMean,zMean,vInc) ;
 211
 212   Float_t depthx = 0.;
 213   Float_t depthz = 0.;
 214   if (fEnergy>0&&vInc.Y()!=0.) {
 215     depthx = ( para * TMath::Log(fEnergy) + parb ) * vInc.X()/TMath::Abs(vInc.Y()) ;
 216     depthz = ( para * TMath::Log(fEnergy) + parb ) * vInc.Z()/TMath::Abs(vInc.Y()) ;
 217   }
 218   else
 219     AliError(Form("Wrong amplitude %f\n", fEnergy));
 220
 221   xMean-= depthx  ;
 222   zMean-= depthz  ;
 223
 224   //Go to the global system
 225   TVector3 gps ;
 226   phosgeom->Local2Global(phosMod, xMean, zMean, gps) ;
 227   fGlobalPos[0]=gps[0] ;
 228   fGlobalPos[1]=gps[1] ;
 229   fGlobalPos[2]=gps[2] ;
 230 }