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