]>
Commit | Line | Data |
---|---|---|
0c5b726e | 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 | // AliAODCaloCluster extension for EMCAL to recalculate cluster | |
18 | // parameters in case of recalibration. | |
19 | // Copy-paste from methods in AliEMCALRecPoint. | |
20 | //*-- | |
21 | //*-- Author: Dmitri Peressounko (RRC KI) for PHOS | |
22 | //*-- Adapted for EMCAL: Gustavo Conesa (INFN-LNF) | |
23 | ||
24 | // --- ROOT system --- | |
25 | #include "TVector3.h" | |
26 | #include "TMath.h" | |
27 | ||
28 | // --- Standard library --- | |
29 | ||
30 | // --- AliRoot header files --- | |
31 | #include "AliLog.h" | |
32 | #include "AliEMCALGeometry.h" | |
33 | #include "AliEMCALPID.h" | |
34 | #include "AliEMCALAodCluster.h" | |
35 | #include "AliEMCALCalibData.h" | |
36 | #include "AliAODCaloCells.h" | |
37 | ||
38 | ClassImp(AliEMCALAodCluster) | |
39 | ||
40 | //____________________________________________________________________________ | |
41 | AliEMCALAodCluster::AliEMCALAodCluster() : | |
42 | AliAODCaloCluster(),fRecalibrated(0) | |
43 | { | |
44 | // ctor | |
45 | } | |
46 | //____________________________________________________________________________ | |
47 | AliEMCALAodCluster::AliEMCALAodCluster(const AliAODCaloCluster & clu) : | |
48 | AliAODCaloCluster(clu),fRecalibrated(0) | |
49 | { | |
50 | // cpy ctor | |
51 | } | |
52 | ||
53 | //____________________________________________________________________________ | |
54 | AliEMCALAodCluster::~AliEMCALAodCluster() | |
55 | { | |
56 | // dtor | |
57 | } | |
58 | //____________________________________________________________________________ | |
59 | void AliEMCALAodCluster::Recalibrate(AliEMCALCalibData * calibData, AliAODCaloCells *emcCells, TString emcalGeoName){ | |
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 | AliEMCALGeometry * emcalgeo = AliEMCALGeometry::GetInstance(emcalGeoName) ; | |
7e1d9a9b | 70 | if(emcalgeo){ |
0c5b726e | 71 | |
72 | Double32_t * cellsAmpFraction = GetCellsAmplitudeFraction(); | |
73 | Int_t iSupMod = -1; | |
74 | Int_t iTower = -1; | |
75 | Int_t iIphi = -1; | |
76 | Int_t iIeta = -1; | |
77 | Int_t iphi = -1; | |
78 | Int_t ieta = -1; | |
79 | ||
80 | for(Int_t i=0; i < GetNCells(); i++){ | |
81 | ||
82 | //Get from the absid the supermodule, tower and eta/phi numbers | |
83 | emcalgeo->GetCellIndex(GetCellAbsId(i),iSupMod,iTower,iIphi,iIeta); | |
84 | //Gives SuperModule and Tower numbers | |
85 | emcalgeo->GetCellPhiEtaIndexInSModule(iSupMod,iTower, | |
86 | iIphi, iIeta,iphi,ieta); | |
87 | ||
88 | Double_t energy = emcCells->GetCellAmplitude(GetCellAbsId(i)) ; | |
60ea57cf | 89 | //AliDebug(2,Form("Recalibrate: cell %f, calib %f, fraction %f\n",energy,calibData->GetADCchannel(iSupMod,ieta,iphi),cellsAmpFraction[i])); |
90 | if(cellsAmpFraction[i]< 1e-4) cellsAmpFraction[i] = 1; //Unfolding off, took all the cell energy | |
0c5b726e | 91 | cellsAmpFraction[i]*=energy*calibData->GetADCchannel(iSupMod,ieta,iphi); |
92 | } | |
93 | ||
94 | SetCellsAmplitudeFraction(cellsAmpFraction); | |
95 | fRecalibrated=kTRUE; | |
7e1d9a9b | 96 | |
97 | } | |
98 | else AliFatal("AliEMCALGeometry was not constructed\n") ; | |
99 | ||
0c5b726e | 100 | } |
101 | //____________________________________________________________________________ | |
102 | void AliEMCALAodCluster::EvalAll(Float_t logWeight, TString geoname){ | |
103 | //If recalibrated - recalculate all cluster parameters | |
104 | if(!fRecalibrated) | |
105 | return ; | |
106 | //printf("EvalAll e org %f\n",E()); | |
107 | EvalEnergy() ; //Energy should be evaluated first | |
108 | //printf("EvalAll e2 %f\n",E()); | |
109 | EvalPositionAndShowerShape(logWeight, geoname) ; | |
110 | //printf("EvalAll e3 %f\n",E()); | |
111 | EvalPID() ; //Should be evaluated after energy and shower shape recalculation | |
112 | //printf("EvalAll e4 %f\n",E()); | |
113 | } | |
114 | //____________________________________________________________________________ | |
115 | void AliEMCALAodCluster::EvalEnergy(){ | |
116 | //Evaluate energy | |
117 | if(!fRecalibrated) // no need to recalibrate | |
118 | return ; | |
119 | ||
120 | Float_t energy=0. ; | |
121 | for(Int_t iDigit=0; iDigit < GetNCells(); iDigit++) { | |
122 | energy+=GetCellAmplitudeFraction(iDigit) ; | |
123 | } | |
124 | //printf("EvalEnergy: e %f\n", energy); | |
125 | SetE(energy); | |
126 | ||
127 | ||
128 | } | |
129 | ////____________________________________________________________________________ | |
130 | //void AliEMCALAodCluster::EnergyCorrection(AliEMCALPID * pid){ | |
131 | // //apply nonlinearity correction same as in AliEMCALPID. | |
132 | // SetE(pid->GetCalibratedEnergy(E())) ; | |
133 | //} | |
134 | ||
135 | //____________________________________________________________________________ | |
136 | void AliEMCALAodCluster::EvalPID(){ | |
137 | ||
138 | //re-evaluate identification parameters | |
139 | // pid->CalculatePID(E(),GetDispersion(),GetEmcCpvDistance(),GetTOF(),GetPID()) ; | |
140 | // pid->CalculatePID(E(),GetDispersion(),GetM20(),GetM02(),GetEmcCpvDistance(),GetTOF(),GetPID()) ; | |
141 | ||
142 | //With bayesian | |
143 | AliEMCALPID *pid = new AliEMCALPID(kFALSE); | |
144 | pid->SetLowFluxParam(); // Need to be fixed | |
145 | Float_t pidlist[AliPID::kSPECIESN+1]; | |
146 | for(Int_t i = 0; i < AliPID::kSPECIESN+1; i++) pidlist[i] = pid->GetPIDFinal(i); | |
147 | SetPIDFromESD(pidlist); | |
148 | ||
149 | } | |
150 | ||
151 | //____________________________________________________________________________ | |
152 | void AliEMCALAodCluster::EvalPositionAndShowerShape(Float_t logWeight, TString emcalGeoName) | |
153 | { | |
154 | // Calculates new center of gravity in the local EMCAL-module coordinates | |
155 | // and tranfers into global ALICE coordinates | |
156 | // Calculates Dispersion and main axis | |
157 | if(!fRecalibrated) // no need to recalibrate | |
158 | return ; | |
159 | ||
160 | Int_t nstat = 0; | |
161 | Float_t wtot = 0. ; | |
162 | ||
0c5b726e | 163 | Int_t iSupMod = -1; |
164 | Int_t iTower = -1; | |
165 | Int_t iIphi = -1; | |
166 | Int_t iIeta = -1; | |
167 | Int_t iphi = -1; | |
168 | Int_t ieta = -1; | |
169 | Double_t etai = -1.; | |
170 | Double_t phii = -1.; | |
171 | ||
172 | Double_t clXYZ[3] ={0.,0.,0.}; | |
173 | Double_t xyzi[3] ={0.,0.,0.}; | |
174 | ||
175 | Double_t d = 0.; | |
176 | Double_t dxx = 0.; | |
177 | Double_t dzz = 0.; | |
178 | Double_t dxz = 0.; | |
179 | Double_t xmean = 0.; | |
180 | Double_t zmean = 0.; | |
181 | ||
182 | AliEMCALGeometry * emcalgeo = AliEMCALGeometry::GetInstance(emcalGeoName) ; | |
7e1d9a9b | 183 | if(emcalgeo){ |
0c5b726e | 184 | |
60ea57cf | 185 | Double_t dist = TmaxInCm(E(),0); |
0c5b726e | 186 | for(Int_t iDigit=0; iDigit < GetNCells(); iDigit++) { |
7cfcebd3 | 187 | |
0c5b726e | 188 | //Get from the absid the supermodule, tower and eta/phi numbers |
189 | emcalgeo->GetCellIndex(GetCellAbsId(iDigit),iSupMod,iTower,iIphi,iIeta); | |
7cfcebd3 | 190 | emcalgeo->RelPosCellInSModule(GetCellAbsId(iDigit), dist, xyzi[0], xyzi[1], xyzi[2]); |
0c5b726e | 191 | emcalgeo->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi,iIeta, iphi,ieta); |
192 | ||
193 | Double_t ei = GetCellAmplitudeFraction(iDigit) ; | |
194 | if (E() > 0 && ei > 0) { | |
195 | Float_t w = ei; | |
196 | if(logWeight > 0) w = TMath::Max( 0., logWeight + TMath::Log(ei/E()) ) ; | |
197 | etai=(Double_t)ieta; | |
198 | phii=(Double_t)iphi; | |
199 | if(w > 0.0) { | |
200 | wtot += w ; | |
201 | nstat++; | |
202 | for(Int_t i = 0; i < 3; i++ ) clXYZ[i] += (w*xyzi[i]); | |
203 | ||
204 | //Shower shape | |
205 | dxx += w * etai * etai ; | |
206 | xmean+= w * etai ; | |
207 | dzz += w * phii * phii ; | |
208 | zmean+= w * phii ; | |
209 | dxz += w * etai * phii ; | |
210 | } | |
211 | } | |
212 | else | |
213 | AliError(Form("Wrong energy %f and/or amplitude %f\n", ei, E())); | |
214 | } | |
215 | ||
216 | //Normalize to the weight | |
217 | if (wtot > 0) { | |
218 | for(Int_t i=0; i<3; i++ ) clXYZ[i] /= wtot; | |
219 | xmean /= wtot ; | |
220 | zmean /= wtot ; | |
221 | } | |
222 | else | |
223 | AliError(Form("Wrong weight %f\n", wtot)); | |
224 | ||
225 | //Put cluster position in the global system | |
226 | TVector3 gpos ; | |
227 | emcalgeo->GetGlobal(clXYZ, gpos, iSupMod); | |
228 | ||
c8fe2783 | 229 | SetPositionAt(gpos[0],0) ; |
230 | SetPositionAt(gpos[1],1) ; | |
231 | SetPositionAt(gpos[2],2) ; | |
0c5b726e | 232 | |
233 | //Calculate dispersion | |
234 | for(Int_t iDigit=0; iDigit < GetNCells(); iDigit++) { | |
235 | //Get from the absid the supermodule, tower and eta/phi numbers | |
236 | emcalgeo->GetCellIndex(GetCellAbsId(iDigit),iSupMod,iTower,iIphi,iIeta); | |
237 | emcalgeo->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi,iIeta, iphi,ieta); | |
238 | ||
239 | Double_t ei=GetCellAmplitudeFraction(iDigit) ; | |
240 | if (E() > 0 && ei > 0) { | |
241 | Float_t w = ei; | |
242 | if(logWeight > 0) w = TMath::Max( 0., logWeight + TMath::Log(ei/E()) ) ; | |
243 | etai=(Double_t)ieta; | |
244 | phii=(Double_t)iphi; | |
245 | if(w > 0.0) d += w*((etai-xmean)*(etai-xmean)+(phii-zmean)*(phii-zmean)); | |
246 | } | |
247 | else | |
248 | AliError(Form("Wrong energy %f and/or amplitude %f\n", ei, E())); | |
249 | } | |
250 | ||
251 | //Normalize to the weigth and set shower shape parameters | |
252 | if (wtot > 0 && nstat > 1) { | |
253 | d /= wtot ; | |
254 | dxx /= wtot ; | |
255 | dzz /= wtot ; | |
256 | dxz /= wtot ; | |
257 | dxx -= xmean * xmean ; | |
258 | dzz -= zmean * zmean ; | |
259 | dxz -= xmean * zmean ; | |
260 | SetM02(0.5 * (dxx + dzz) + TMath::Sqrt( 0.25 * (dxx - dzz) * (dxx - dzz) + dxz * dxz )) ; | |
261 | SetM20(0.5 * (dxx + dzz) - TMath::Sqrt( 0.25 * (dxx - dzz) * (dxx - dzz) + dxz * dxz )); | |
262 | } | |
263 | else{ | |
264 | d=0. ; | |
265 | SetM20(0.) ; | |
266 | SetM02(0.) ; | |
267 | } | |
268 | ||
269 | if (d>=0) | |
270 | SetDispersion(TMath::Sqrt(d)) ; | |
271 | else | |
272 | SetDispersion(0) ; | |
7e1d9a9b | 273 | } |
274 | else { | |
275 | AliFatal("AliEMCALGeometry was not constructed\n") ; | |
276 | } | |
277 | ||
0c5b726e | 278 | |
279 | } | |
280 | ||
281 | //_____________________________________________________________________ | |
282 | Double_t AliEMCALAodCluster::TmaxInCm(const Double_t e , const Int_t key) const | |
283 | { | |
284 | // e energy in GeV) | |
285 | // key = 0(gamma, default) | |
286 | // != 0(electron) | |
287 | static Double_t ca = 4.82; // shower max parameter - first guess; ca=TMath::Log(1000./8.07) | |
288 | static Double_t x0 = 1.23; // radiation lenght (cm) | |
289 | static Double_t tmax = 0.; // position of electromagnetic shower max in cm | |
290 | ||
291 | tmax = 0.0; | |
292 | if(e>0.1) { | |
293 | tmax = TMath::Log(e) + ca; | |
294 | if (key==0) tmax += 0.5; | |
295 | else tmax -= 0.5; | |
296 | tmax *= x0; // convert to cm | |
297 | } | |
298 | return tmax; | |
299 | } | |
300 |