Go to the documentation of this file.00001 #include "RecoEgamma/EgammaPhotonAlgos/interface/EnergyUncertaintyPhotonSpecific.h"
00002 #include "TMath.h"
00003
00004 #include <iostream>
00005
00006
00007 EnergyUncertaintyPhotonSpecific::EnergyUncertaintyPhotonSpecific( const edm::ParameterSet& config ) {
00008
00009
00010
00011 }
00012
00013
00014 EnergyUncertaintyPhotonSpecific::~EnergyUncertaintyPhotonSpecific() {
00015
00016 }
00017
00018
00019
00020 void EnergyUncertaintyPhotonSpecific::init ( const edm::EventSetup& theEventSetup ) {
00021
00022 }
00023
00024
00025 double EnergyUncertaintyPhotonSpecific::computePhotonEnergyUncertainty_lowR9(double eta, double brem, double energy){
00026
00027
00028
00029 double et = energy/cosh(eta);
00030
00031 const int nBinsEta=6;
00032 const double EtaBins[nBinsEta+1] = {0.0, 0.7, 1.15, 1.44, 1.56, 2.0, 2.5};
00033
00034 const int nBinsBrem=6;
00035 const double BremBins [nBinsBrem+1]= {0.8, 1.0, 2.0, 3.0, 4.0, 5.0, 10.0};
00036
00037
00038 float par0[nBinsEta][nBinsBrem];
00039 float par1[nBinsEta][nBinsBrem];
00040 float par2[nBinsEta][nBinsBrem];
00041 float par3[nBinsEta][nBinsBrem];
00042
00043 par0[0][0]=0.0232291;
00044 par1[0][0]=0;
00045 par2[0][0]=0;
00046 par3[0][0]=0;
00047
00048 par0[0][1]=0.00703187;
00049 par1[0][1]=0.646644;
00050 par2[0][1]=-7.4698;
00051 par3[0][1]=5.53373e-08;
00052
00053 par0[0][2]=0.00692465;
00054 par1[0][2]=0.292698;
00055 par2[0][2]=4.16907;
00056 par3[0][2]=5.61149e-06;
00057
00058 par0[0][3]=0.00855993;
00059 par1[0][3]=0.280843;
00060 par2[0][3]=4.25527;
00061 par3[0][3]=9.6404e-07;
00062
00063 par0[0][4]=0.00795058;
00064 par1[0][4]=0.370007;
00065 par2[0][4]=3.03429;
00066 par3[0][4]=4.43986e-07;
00067
00068 par0[0][5]=0.0107494;
00069 par1[0][5]=0.276159;
00070 par2[0][5]=4.44532;
00071 par3[0][5]=2.58822e-06;
00072
00073 par0[1][0]=0.0614866;
00074 par1[1][0]=0;
00075 par2[1][0]=0;
00076 par3[1][0]=0;
00077
00078 par0[1][1]=0.00894211;
00079 par1[1][1]=0.466937;
00080 par2[1][1]=3.33434;
00081 par3[1][1]=0.000114835;
00082
00083 par0[1][2]=0.0102959;
00084 par1[1][2]=0.313568;
00085 par2[1][2]=6.34301;
00086 par3[1][2]=2.86726e-07;
00087
00088 par0[1][3]=0.0128934;
00089 par1[1][3]=0.302943;
00090 par2[1][3]=6.35598;
00091 par3[1][3]=0.00190694;
00092
00093 par0[1][4]=0.0130199;
00094 par1[1][4]=0.505135;
00095 par2[1][4]=2.52964;
00096 par3[1][4]=0.120204;
00097
00098 par0[1][5]=0.0180839;
00099 par1[1][5]=0.382134;
00100 par2[1][5]=5.3388;
00101 par3[1][5]=3.59921e-07;
00102
00103 par0[2][0]=0.0291343;
00104 par1[2][0]=0;
00105 par2[2][0]=0;
00106 par3[2][0]=0;
00107
00108 par0[2][1]=0.00876269;
00109 par1[2][1]=0.375159;
00110 par2[2][1]=7.11411;
00111 par3[2][1]=0.0438575;
00112
00113 par0[2][2]=0.0120863;
00114 par1[2][2]=0.397635;
00115 par2[2][2]=5.97451;
00116 par3[2][2]=0.0469782;
00117
00118 par0[2][3]=0.0112655;
00119 par1[2][3]=0.856565;
00120 par2[2][3]=-5.76122;
00121 par3[2][3]=4.99993;
00122
00123 par0[2][4]=0.0168267;
00124 par1[2][4]=0.636468;
00125 par2[2][4]=-1.54548;
00126 par3[2][4]=4.99992;
00127
00128 par0[2][5]=0.0168059;
00129 par1[2][5]=1.09268;
00130 par2[2][5]=-0.547554;
00131 par3[2][5]=0.0952985;
00132
00133 par0[3][0]=0.158403;
00134 par1[3][0]=0;
00135 par2[3][0]=0;
00136 par3[3][0]=0;
00137
00138 par0[3][1]=0.0717431;
00139 par1[3][1]=1.66981;
00140 par2[3][1]=6.86275;
00141 par3[3][1]=0.00543544;
00142
00143 par0[3][2]=0.0385666;
00144 par1[3][2]=3.6319;
00145 par2[3][2]=-3.76633;
00146 par3[3][2]=6.56718e-05;
00147
00148 par0[3][3]=0.0142631;
00149 par1[3][3]=8.85991;
00150 par2[3][3]=-32.6073;
00151 par3[3][3]=0.00119538;
00152
00153 par0[3][4]=0.0421638;
00154 par1[3][4]=3.1289;
00155 par2[3][4]=-6.58653;
00156 par3[3][4]=1.10125e-05;
00157
00158 par0[3][5]=0.046331;
00159 par1[3][5]=1.29951;
00160 par2[3][5]=1.76117;
00161 par3[3][5]=0.00204206;
00162
00163 par0[4][0]=0.0483944;
00164 par1[4][0]=0;
00165 par2[4][0]=0;
00166 par3[4][0]=0;
00167
00168 par0[4][1]=0.0168516;
00169 par1[4][1]=1.19617;
00170 par2[4][1]=-6.78666;
00171 par3[4][1]=4.98192;
00172
00173 par0[4][2]=0.0243039;
00174 par1[4][2]=0.994626;
00175 par2[4][2]=-4.26073;
00176 par3[4][2]=4.99984;
00177
00178 par0[4][3]=0.031795;
00179 par1[4][3]=0.875925;
00180 par2[4][3]=1.43183;
00181 par3[4][3]=0.0920944;
00182
00183 par0[4][4]=0.0414953;
00184 par1[4][4]=0.654605;
00185 par2[4][4]=4.45367;
00186 par3[4][4]=0.030385;
00187
00188 par0[4][5]=0.058031;
00189 par1[4][5]=0.292915;
00190 par2[4][5]=8.48307;
00191 par3[4][5]=0.0134321;
00192
00193 par0[5][0]=0.107158;
00194 par1[5][0]=0;
00195 par2[5][0]=0;
00196 par3[5][0]=0;
00197
00198 par0[5][1]=0.021685;
00199 par1[5][1]=0.574207;
00200 par2[5][1]=-0.566981;
00201 par3[5][1]=0.0120609;
00202
00203 par0[5][2]=0.0196619;
00204 par1[5][2]=0.940217;
00205 par2[5][2]=-6.05845;
00206 par3[5][2]=0.000193818;
00207
00208 par0[5][3]=0.0324734;
00209 par1[5][3]=0.574766;
00210 par2[5][3]=-5.23571;
00211 par3[5][3]=4.9419;
00212
00213 par0[5][4]=0.0414953;
00214 par1[5][4]=0.654605;
00215 par2[5][4]=4.45367;
00216 par3[5][4]=0.030385;
00217
00218 par0[5][5]=0.058031;
00219 par1[5][5]=0.292915;
00220 par2[5][5]=8.48307;
00221 par3[5][5]=0.0134321;
00222
00223
00224 int iEtaSl = -1;
00225 for (int iEta = 0; iEta < nBinsEta; ++iEta){
00226 if ( EtaBins[iEta] <= TMath::Abs(eta) && TMath::Abs(eta) <EtaBins[iEta+1] ){
00227 iEtaSl = iEta;
00228 }
00229 }
00230
00231
00232 int iBremSl = -1;
00233 for (int iBrem = 0; iBrem < nBinsBrem; ++iBrem){
00234 if ( BremBins[iBrem] <= brem && brem <BremBins[iBrem+1] ){
00235 iBremSl = iBrem;
00236 }
00237 }
00238
00239 if (TMath::Abs(eta)>2.5) iEtaSl = nBinsEta-1;
00240 if (brem<BremBins[0]) iBremSl = 0;
00241 if (brem>BremBins[nBinsBrem-1]) iBremSl = nBinsBrem-1;
00242
00243 float uncertainty = 0;
00244 if (et<5) uncertainty = par0[iEtaSl][iBremSl] + par1[iEtaSl][iBremSl]/(5-par2[iEtaSl][iBremSl]) + par3[iEtaSl][iBremSl]/((5-par2[iEtaSl][iBremSl])*(5-par2[iEtaSl][iBremSl]));
00245 if (et>200) uncertainty = par0[iEtaSl][iBremSl] + par1[iEtaSl][iBremSl]/(200-par2[iEtaSl][iBremSl]) + par3[iEtaSl][iBremSl]/((200-par2[iEtaSl][iBremSl])*(200-par2[iEtaSl][iBremSl]));
00246
00247 if (et>5 && et<200) uncertainty = par0[iEtaSl][iBremSl] + par1[iEtaSl][iBremSl]/(et-par2[iEtaSl][iBremSl]) + par3[iEtaSl][iBremSl]/((et-par2[iEtaSl][iBremSl])*(et-par2[iEtaSl][iBremSl]));
00248
00249 return (uncertainty*energy);
00250
00251 }
00252
00253 double EnergyUncertaintyPhotonSpecific::computePhotonEnergyUncertainty_highR9(double eta, double brem, double energy){
00254
00255
00256
00257 double et = energy/cosh(eta);
00258
00259 const int nBinsEta=6;
00260 const double EtaBins[nBinsEta+1] = {0.0, 0.7, 1.15, 1.44, 1.56, 2.0, 2.5};
00261
00262 const int nBinsBrem=6;
00263 const double BremBins [nBinsBrem+1]= {0.8, 1.0, 2.0, 3.0, 4.0, 5.0, 10.0};
00264
00265 float par0[nBinsEta][nBinsBrem];
00266 float par1[nBinsEta][nBinsBrem];
00267 float par2[nBinsEta][nBinsBrem];
00268 float par3[nBinsEta][nBinsBrem];
00269
00270 par0[0][0]=0.00806753;
00271 par1[0][0]=0.143754;
00272 par2[0][0]=-0.00368104;
00273 par3[0][0]=0.219829;
00274
00275 par0[0][1]=0.00899298;
00276 par1[0][1]=0.10159;
00277 par2[0][1]=4.70884;
00278 par3[0][1]=9.07419e-08;
00279
00280 par0[0][2]=0;
00281 par1[0][2]=0;
00282 par2[0][2]=0;
00283 par3[0][2]=0;
00284
00285 par0[0][3]=0;
00286 par1[0][3]=0;
00287 par2[0][3]=0;
00288 par3[0][3]=0;
00289
00290 par0[0][4]=0;
00291 par1[0][4]=0;
00292 par2[0][4]=0;
00293 par3[0][4]=0;
00294
00295 par0[0][5]=0;
00296 par1[0][5]=0;
00297 par2[0][5]=0;
00298 par3[0][5]=0;
00299
00300 par0[1][0]=0.00880649;
00301 par1[1][0]=0.0716169;
00302 par2[1][0]=5.23856;
00303 par3[1][0]=0.00632907;
00304
00305 par0[1][1]=0.00972275;
00306 par1[1][1]=0.0752675;
00307 par2[1][1]=3.35623;
00308 par3[1][1]=2.49397e-07;
00309
00310 par0[1][2]=0;
00311 par1[1][2]=0;
00312 par2[1][2]=0;
00313 par3[1][2]=0;
00314
00315 par0[1][3]=0;
00316 par1[1][3]=0;
00317 par2[1][3]=0;
00318 par3[1][3]=0;
00319
00320 par0[1][4]=0;
00321 par1[1][4]=0;
00322 par2[1][4]=0;
00323 par3[1][4]=0;
00324
00325 par0[1][5]=0;
00326 par1[1][5]=0;
00327 par2[1][5]=0;
00328 par3[1][5]=0;
00329
00330 par0[2][0]=0.0101474;
00331 par1[2][0]=-0.332171;
00332 par2[2][0]=-31.8456;
00333 par3[2][0]=22.543;
00334
00335 par0[2][1]=0.0109109;
00336 par1[2][1]=0.0425903;
00337 par2[2][1]=6.52561;
00338 par3[2][1]=2.18593e-08;
00339
00340 par0[2][2]=0;
00341 par1[2][2]=0;
00342 par2[2][2]=0;
00343 par3[2][2]=0;
00344
00345 par0[2][3]=0;
00346 par1[2][3]=0;
00347 par2[2][3]=0;
00348 par3[2][3]=0;
00349
00350 par0[2][4]=0;
00351 par1[2][4]=0;
00352 par2[2][4]=0;
00353 par3[2][4]=0;
00354
00355 par0[2][5]=0;
00356 par1[2][5]=0;
00357 par2[2][5]=0;
00358 par3[2][5]=0;
00359
00360 par0[3][0]=0.00343003;
00361 par1[3][0]=11.5791;
00362 par2[3][0]=-112.084;
00363 par3[3][0]=-863.968;
00364
00365 par0[3][1]=0.0372159;
00366 par1[3][1]=1.44028;
00367 par2[3][1]=-40;
00368 par3[3][1]=0.00102639;
00369
00370 par0[3][2]=0;
00371 par1[3][2]=0;
00372 par2[3][2]=0;
00373 par3[3][2]=0;
00374
00375 par0[3][3]=0;
00376 par1[3][3]=0;
00377 par2[3][3]=0;
00378 par3[3][3]=0;
00379
00380 par0[3][4]=0;
00381 par1[3][4]=0;
00382 par2[3][4]=0;
00383 par3[3][4]=0;
00384
00385 par0[3][5]=0;
00386 par1[3][5]=0;
00387 par2[3][5]=0;
00388 par3[3][5]=0;
00389
00390 par0[4][0]=0.0192411;
00391 par1[4][0]=0.0511006;
00392 par2[4][0]=7.56304;
00393 par3[4][0]=0.00331583;
00394
00395 par0[4][1]=0.0195124;
00396 par1[4][1]=0.104321;
00397 par2[4][1]=5.71476;
00398 par3[4][1]=6.12472e-06;
00399
00400 par0[4][2]=0;
00401 par1[4][2]=0;
00402 par2[4][2]=0;
00403 par3[4][2]=0;
00404
00405 par0[4][3]=0;
00406 par1[4][3]=0;
00407 par2[4][3]=0;
00408 par3[4][3]=0;
00409
00410 par0[4][4]=0;
00411 par1[4][4]=0;
00412 par2[4][4]=0;
00413 par3[4][4]=0;
00414
00415 par0[4][5]=0;
00416 par1[4][5]=0;
00417 par2[4][5]=0;
00418 par3[4][5]=0;
00419
00420 par0[5][0]=0.0203644;
00421 par1[5][0]=-0.050789;
00422 par2[5][0]=-7.96854;
00423 par3[5][0]=4.71223;
00424
00425 par0[5][1]=0.0198718;
00426 par1[5][1]=0.106859;
00427 par2[5][1]=3.54235;
00428 par3[5][1]=6.89631e-06;
00429
00430 par0[5][2]=0;
00431 par1[5][2]=0;
00432 par2[5][2]=0;
00433 par3[5][2]=0;
00434
00435 par0[5][3]=0;
00436 par1[5][3]=0;
00437 par2[5][3]=0;
00438 par3[5][3]=0;
00439
00440 par0[5][4]=0;
00441 par1[5][4]=0;
00442 par2[5][4]=0;
00443 par3[5][4]=0;
00444
00445 par0[5][5]=0;
00446 par1[5][5]=0;
00447 par2[5][5]=0;
00448 par3[5][5]=0;
00449
00450
00451 int iEtaSl = -1;
00452 for (int iEta = 0; iEta < nBinsEta; ++iEta){
00453 if ( EtaBins[iEta] <= TMath::Abs(eta) && TMath::Abs(eta) <EtaBins[iEta+1] ){
00454 iEtaSl = iEta;
00455 }
00456 }
00457
00458
00459 int iBremSl = -1;
00460 for (int iBrem = 0; iBrem < nBinsBrem; ++iBrem){
00461 if ( BremBins[iBrem] <= brem && brem <BremBins[iBrem+1] ){
00462 iBremSl = iBrem;
00463 }
00464 }
00465
00466 if (TMath::Abs(eta)>2.5) iEtaSl = nBinsEta-1;
00467 if (brem<BremBins[0]) iBremSl = 0;
00468 if (brem>BremBins[nBinsBrem-1]) iBremSl = nBinsBrem-1;
00469 if (brem>2) iBremSl = 1;
00470
00471 float uncertainty = 0;
00472 if (et<5) uncertainty = par0[iEtaSl][iBremSl] + par1[iEtaSl][iBremSl]/(5-par2[iEtaSl][iBremSl]) + par3[iEtaSl][iBremSl]/((5-par2[iEtaSl][iBremSl])*(5-par2[iEtaSl][iBremSl]));
00473 if (et>200) uncertainty = par0[iEtaSl][iBremSl] + par1[iEtaSl][iBremSl]/(200-par2[iEtaSl][iBremSl]) + par3[iEtaSl][iBremSl]/((200-par2[iEtaSl][iBremSl])*(200-par2[iEtaSl][iBremSl]));
00474
00475 if (et>5 && et<200) uncertainty = par0[iEtaSl][iBremSl] + par1[iEtaSl][iBremSl]/(et-par2[iEtaSl][iBremSl]) + par3[iEtaSl][iBremSl]/((et-par2[iEtaSl][iBremSl])*(et-par2[iEtaSl][iBremSl]));
00476
00477 return (uncertainty*energy);
00478
00479 }