made abstract functions in AliVVvertex non-abstract, becasue AliFlatESDVertex does...
[u/mrichter/AliRoot.git] / 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 "AliPHOSReconstructor.h"
34 #include "AliPHOSEsdCluster.h" 
35 #include "AliPHOSCalibData.h"
36 #include "AliESDCaloCells.h"
37
38 ClassImp(AliPHOSEsdCluster)
39
40 //____________________________________________________________________________
41 AliPHOSEsdCluster::AliPHOSEsdCluster() : 
42   AliESDCaloCluster(),fRecalibrated(0)
43 {
44   // ctor
45 }
46 //____________________________________________________________________________
47 AliPHOSEsdCluster::AliPHOSEsdCluster(const AliESDCaloCluster & clu) : 
48   AliESDCaloCluster(clu),fRecalibrated(0)
49 {
50   // cpy ctor
51 }
52
53 //____________________________________________________________________________
54 AliPHOSEsdCluster::~AliPHOSEsdCluster()
55 {
56   // dtor
57 }
58 //____________________________________________________________________________
59 void AliPHOSEsdCluster::Recalibrate(AliPHOSCalibData * calibData,AliESDCaloCells *phsCells){
60   //If not done yet, apply recalibration coefficients to energies list
61   //NOTE that after recalibration fCellsAmpFraction contains not FRACTION but FULL energy 
62   
63   if(fRecalibrated)
64     return ;
65   
66   if(!calibData)
67     return ;
68   
69   AliPHOSGeometry * phosgeom =  AliPHOSGeometry::GetInstance() ;
70   if(!phosgeom)
71     AliFatal("AliPHOSGeometry was not contructed\n") ;
72   
73   for(Int_t i=0; i<fNCells; i++){
74     Int_t relId[4];
75     phosgeom->AbsToRelNumbering(fCellsAbsId[i],relId) ;
76     Int_t   module = relId[0];
77     Int_t   column = relId[3];
78     Int_t   row    = relId[2];
79     Double_t energy = phsCells->GetCellAmplitude(fCellsAbsId[i]) ;
80     fCellsAmpFraction[i]*=energy*calibData->GetADCchannelEmc(module,column,row);
81   }
82   
83   fRecalibrated=kTRUE; 
84 }
85 //____________________________________________________________________________
86 void  AliPHOSEsdCluster::EvalAll(Float_t logWeight, TVector3 &vtx){
87     //If recalibrated - recalculate all cluster parameters
88   if(!fRecalibrated)
89     return ;
90
91   EvalEnergy() ; //Energy should be evaluated first
92   EvalCoord(logWeight, vtx) ;
93   
94 }
95 //____________________________________________________________________________
96 void AliPHOSEsdCluster::EvalEnergy(){
97   if(!fRecalibrated) // no need to recalibrate
98     return ;
99     
100   fEnergy=0. ;
101   for(Int_t iDigit=0; iDigit<fNCells; iDigit++) {
102     fEnergy+=fCellsAmpFraction[iDigit] ;
103   }
104   //Correct for nonlinearity later   
105 }
106 //____________________________________________________________________________
107 void AliPHOSEsdCluster::EnergyCorrection(){
108   //apply nonlinearity correction 
109   fEnergy = AliPHOSReconstructor::CorrectNonlinearity(fEnergy) ;
110 }
111 //____________________________________________________________________________
112 void AliPHOSEsdCluster::EvalPID(AliPHOSPIDv1 * /*pid*/){           
113   //re-evaluate identification parameters
114 //  pid->CalculatePID(fEnergy,fDispersion,fEmcCpvDistance,tof,fPID) ;  
115 //  pid->CalculatePID(fEnergy,fDispersion,fM20,fM02,fEmcCpvDistance,tof,fPID) ;
116 }
117 //____________________________________________________________________________
118 void AliPHOSEsdCluster::EvalCoord(Float_t logWeight, TVector3 &vtx)
119 {
120   // Calculates new center of gravity in the local PHOS-module coordinates 
121   // and tranfers into global ALICE coordinates
122   // Calculates Dispersion and main axis
123   if(!fRecalibrated) // no need to recalibrate
124     return ;
125  
126   Float_t wtot = 0. ;
127   Int_t relid[4] ;
128   Int_t phosMod=0 ;
129   Float_t xMean = 0. ;
130   Float_t zMean = 0. ;
131
132   AliPHOSGeometry * phosgeom =  AliPHOSGeometry::GetInstance() ;
133   if(!phosgeom)
134     AliFatal("AliPHOSGeometry was not contructed\n") ;
135
136   for(Int_t iDigit=0; iDigit<fNCells; iDigit++) {
137     Float_t xi ;
138     Float_t zi ;
139     phosgeom->AbsToRelNumbering(fCellsAbsId[iDigit], relid) ;
140     phosgeom->RelPosInModule(relid, xi, zi);
141     phosMod=relid[0] ;
142     Double_t ei=fCellsAmpFraction[iDigit] ;
143     if (fEnergy>0 && ei>0) {
144       Float_t w = TMath::Max( 0., logWeight + TMath::Log(ei/fEnergy) ) ;
145       xMean+= xi * w ;
146       zMean+= zi * w ;
147       wtot += w ;
148     }
149     else
150       AliError(Form("Wrong energy %f and/or amplitude %f\n", ei, fEnergy));
151   }
152   if (wtot>0) {
153     xMean /= wtot ;
154     zMean /= wtot ;
155   }
156   else
157     AliError(Form("Wrong weight %f\n", wtot));
158
159
160 // Calculates the dispersion and second momenta
161   Double_t d=0. ;
162   Double_t dxx  = 0.;
163   Double_t dzz  = 0.;
164   Double_t dxz  = 0.;
165   for(Int_t iDigit=0; iDigit < fNCells; iDigit++) {
166     Float_t xi ;
167     Float_t zi ;
168     phosgeom->AbsToRelNumbering(fCellsAbsId[iDigit], relid) ;
169     phosgeom->RelPosInModule(relid, xi, zi);
170     Double_t ei=fCellsAmpFraction[iDigit] ;
171     if (fEnergy>0 && ei>0) {
172       Float_t w = TMath::Max( 0., logWeight + TMath::Log(ei/fEnergy) ) ;
173       d += w*((xi-xMean)*(xi-xMean) + (zi-zMean)*(zi-zMean) ) ; 
174       dxx  += w * xi * xi ;
175       dzz  += w * zi * zi ;
176       dxz  += w * xi * zi ; 
177    }
178     else
179       AliError(Form("Wrong energy %f and/or amplitude %f\n", ei, fEnergy));
180   }
181   
182   if (wtot>0) {
183     d /= wtot ;
184     dxx /= wtot ;
185     dzz /= wtot ;
186     dxz /= wtot ;
187     dxx -= xMean * xMean ;
188     dzz -= zMean * zMean ;
189     dxz -= xMean * zMean ;
190     fM02 =  0.5 * (dxx + dzz) + TMath::Sqrt( 0.25 * (dxx - dzz) * (dxx - dzz) + dxz * dxz )  ;
191     fM20 =  0.5 * (dxx + dzz) - TMath::Sqrt( 0.25 * (dxx - dzz) * (dxx - dzz) + dxz * dxz )  ;
192   }
193   else{
194     AliError(Form("Wrong weight %f\n", wtot));
195     d=0. ;
196     fM20=0. ;
197     fM02=0. ;
198   }
199
200   if (d>=0)
201     fDispersion = TMath::Sqrt(d) ;
202   else    
203     fDispersion = 0 ;
204
205
206   // Correction for the depth of the shower starting point (TDR p 127)  
207   Float_t para = 0.925 ; 
208   Float_t parb = 6.52 ; 
209
210   TVector3 vInc ;
211   phosgeom->GetIncidentVector(vtx,phosMod,xMean,zMean,vInc) ;
212
213   Float_t depthx = 0.; 
214   Float_t depthz = 0.;
215   if (fEnergy>0&&vInc.Y()!=0.) {
216     depthx = ( para * TMath::Log(fEnergy) + parb ) * vInc.X()/TMath::Abs(vInc.Y()) ;
217     depthz = ( para * TMath::Log(fEnergy) + parb ) * vInc.Z()/TMath::Abs(vInc.Y()) ;
218   }
219   else 
220     AliError(Form("Wrong amplitude %f\n", fEnergy));
221
222   xMean-= depthx  ;
223   zMean-= depthz  ;
224
225   //Go to the global system
226   TVector3 gps ;
227   phosgeom->Local2Global(phosMod, xMean, zMean, gps) ;
228   fGlobalPos[0]=gps[0] ;
229   fGlobalPos[1]=gps[1] ;  
230   fGlobalPos[2]=gps[2] ;
231 }