3 #include <CLHEP/Random/RandGaussQ.h>
4 #include <CLHEP/Random/RandFlat.h>
5 #include <CLHEP/Random/Random.h>
29 std::cout <<
"[ElectronEnergyCalibrator] Initialization in DATA mode" << std::endl;
36 <<
"[ElectronEnergyCalibrator] Cannot open the file "
37 <<
pathData_ <<
"\n It is not found, missed or corrupted" ;
42 std::cout <<
"[ElectronEnergyCalibrator] File "
43 <<
pathData_ <<
" succesfully opened" << std::endl;
47 std::vector<std::string> selements;
78 std::cout <<
"[ElectronEnergyCalibrator] File closed" << std::endl;
90 <<
"[ElectronEnergyCalibrator] Cannot open the file "<<
pathLinData_ <<
"\n It is not found, missed or corrupted" ;
95 std::cout<<
"[ElectronEnergyCalibrator] File with Linearity Corrections "<<
pathLinData_<<
" succesfully opened"<<std::endl;
99 std::vector<std::string> selements;
103 while ( !finlin.eof() )
129 std::cout<<
"[ElectronEnergyCalibrator] File closed"<<std::endl;
137 std::cout <<
"[ElectronEnergyCalibrator] Initialization in MC mode" << std::endl;
147 while ( ( std::string::npos != pos ) || ( std::string::npos != lastpos ) )
149 elements.push_back(fullstr.substr(lastpos, pos-lastpos));
150 lastpos = fullstr.find_first_not_of(delimiter, pos);
151 pos = fullstr.find_first_of(delimiter, lastpos);
157 std::istringstream stm;
170 bool isEB = electron.
isEB();
172 double r9 = electron.
getR9();
179 std::cout <<
"[ElectronEnergyCalibrator] Using regression energy for calibration" << std::endl;
187 std::cout <<
"[ElectronEnergyCalibrator] Using scale corrections for new regression" << std::endl;
195 std::cout <<
"[ElectronEnergyCalibrator] Using standard ecal energy for calibration" << std::endl;
206 <<
"XXXXXXX requires the RandomNumberGeneratorService\n"
207 "which is not present in the configuration file. You must add the service\n"
208 "in the configuration file or remove the modules that require it.";
273 if ( run_ <= 203002 )
275 if ( isEB && fabs(eta) < 1 && r9 < 0.94 ) dsigMC = 0.0103;
276 if ( isEB && fabs(eta) < 1 && r9 >= 0.94 ) dsigMC = 0.0090;
277 if ( isEB && fabs(eta) >= 1 && r9 < 0.94 ) dsigMC = 0.0190;
278 if ( isEB && fabs(eta) >= 1 && r9 >= 0.94 ) dsigMC = 0.0156;
279 if ( !isEB && fabs(eta) < 2 && r9 < 0.94 ) dsigMC = 0.0269;
280 if ( !isEB && fabs(eta) < 2 && r9 >= 0.94 ) dsigMC = 0.0287;
281 if ( !isEB && fabs(eta) >= 2 && r9 < 0.94 ) dsigMC = 0.0364;
282 if ( !isEB && fabs(eta) >= 2 && r9 >= 0.94 ) dsigMC = 0.0321;
285 if ( isEB && fabs(eta) < 1 && r9 < 0.94 ) dsigMC = 0.0109;
286 if ( isEB && fabs(eta) < 1 && r9 >= 0.94 ) dsigMC = 0.0099;
287 if ( isEB && fabs(eta) >= 1 && r9 < 0.94 ) dsigMC = 0.0182;
288 if ( isEB && fabs(eta) >= 1 && r9 >= 0.94 ) dsigMC = 0.0200;
289 if ( !isEB && fabs(eta) < 2 && r9 < 0.94 ) dsigMC = 0.0282;
290 if ( !isEB && fabs(eta) < 2 && r9 >= 0.94 ) dsigMC = 0.0309;
291 if ( !isEB && fabs(eta) >= 2 && r9 < 0.94 ) dsigMC = 0.0386;
292 if ( !isEB && fabs(eta) >= 2 && r9 >= 0.94 ) dsigMC = 0.0359;
296 CLHEP::RandFlat flatRandom(rng->
getEngine(streamID));
297 double rn = flatRandom.fire();
300 if ( isEB && fabs(eta) < 1 && r9 < 0.94 ) dsigMC = 0.0109;
301 if ( isEB && fabs(eta) < 1 && r9 >= 0.94 ) dsigMC = 0.0099;
302 if ( isEB && fabs(eta) >= 1 && r9 < 0.94 ) dsigMC = 0.0182;
303 if ( isEB && fabs(eta) >= 1 && r9 >= 0.94 ) dsigMC = 0.0200;
304 if ( !isEB && fabs(eta) < 2 && r9 < 0.94 ) dsigMC = 0.0282;
305 if ( !isEB && fabs(eta) < 2 && r9 >= 0.94 ) dsigMC = 0.0309;
306 if ( !isEB && fabs(eta) >= 2 && r9 < 0.94 ) dsigMC = 0.0386;
307 if ( !isEB && fabs(eta) >= 2 && r9 >= 0.94 ) dsigMC = 0.0359;
310 if ( isEB && fabs(eta) < 1 && r9 < 0.94 ) dsigMC = 0.0103;
311 if ( isEB && fabs(eta) < 1 && r9 >= 0.94 ) dsigMC = 0.0090;
312 if ( isEB && fabs(eta) >= 1 && r9 < 0.94 ) dsigMC = 0.0190;
313 if ( isEB && fabs(eta) >= 1 && r9 >= 0.94 ) dsigMC = 0.0156;
314 if ( !isEB && fabs(eta) < 2 && r9 < 0.94 ) dsigMC = 0.0269;
315 if ( !isEB && fabs(eta) < 2 && r9 >= 0.94 ) dsigMC = 0.0287;
316 if ( !isEB && fabs(eta) >= 2 && r9 < 0.94 ) dsigMC = 0.0364;
317 if ( !isEB && fabs(eta) >= 2 && r9 >= 0.94 ) dsigMC = 0.0321;
321 if ( isEB && fabs(eta) < 1 && r9 < 0.94 ) dsigMC = 0.0103;
322 if ( isEB && fabs(eta) < 1 && r9 >= 0.94 ) dsigMC = 0.0090;
323 if ( isEB && fabs(eta) >= 1 && r9 < 0.94 ) dsigMC = 0.0190;
324 if ( isEB && fabs(eta) >= 1 && r9 >= 0.94 ) dsigMC = 0.0156;
325 if ( !isEB && fabs(eta) < 2 && r9 < 0.94 ) dsigMC = 0.0269;
326 if ( !isEB && fabs(eta) < 2 && r9 >= 0.94 ) dsigMC = 0.0287;
327 if ( !isEB && fabs(eta) >= 2 && r9 < 0.94 ) dsigMC = 0.0364;
328 if ( !isEB && fabs(eta) >= 2 && r9 >= 0.94 ) dsigMC = 0.0321;
332 if ( isEB && fabs(eta) < 1 && r9 < 0.94 ) dsigMC = 0.0109;
333 if ( isEB && fabs(eta) < 1 && r9 >= 0.94 ) dsigMC = 0.0099;
334 if ( isEB && fabs(eta) >= 1 && r9 < 0.94 ) dsigMC = 0.0182;
335 if ( isEB && fabs(eta) >= 1 && r9 >= 0.94 ) dsigMC = 0.0200;
336 if ( !isEB && fabs(eta) < 2 && r9 < 0.94 ) dsigMC = 0.0282;
337 if ( !isEB && fabs(eta) < 2 && r9 >= 0.94 ) dsigMC = 0.0309;
338 if ( !isEB && fabs(eta) >= 2 && r9 < 0.94 ) dsigMC = 0.0386;
339 if ( !isEB && fabs(eta) >= 2 && r9 >= 0.94 ) dsigMC = 0.0359;
349 if ( isEB && fabs(eta) < 1 && r9 < 0.94 ) dsigMC = 0.0094;
350 if ( isEB && fabs(eta) < 1 && r9 >= 0.94 ) dsigMC = 0.0092;
351 if ( isEB && fabs(eta) >= 1 && r9 < 0.94 ) dsigMC = 0.0182;
352 if ( isEB && fabs(eta) >= 1 && r9 >= 0.94 ) dsigMC = 0.0139;
353 if ( !isEB && fabs(eta) < 2 && r9 < 0.94 ) dsigMC = 0.0220;
354 if ( !isEB && fabs(eta) < 2 && r9 >= 0.94 ) dsigMC = 0.0229;
355 if ( !isEB && fabs(eta) >= 2 && r9 < 0.94 ) dsigMC = 0.0290;
356 if ( !isEB && fabs(eta) >= 2 && r9 >= 0.94 ) dsigMC = 0.0234;
364 if ( isEB && fabs(eta) < 1 && r9 < 0.94 ) dsigMC = 0.01;
365 if ( isEB && fabs(eta) < 1 && r9 >= 0.94 ) dsigMC = 0.0099;
366 if ( isEB && fabs(eta) >= 1 && r9 < 0.94 ) dsigMC = 0.0217;
367 if ( isEB && fabs(eta) >= 1 && r9 >= 0.94 ) dsigMC = 0.0157;
368 if ( !isEB && fabs(eta) < 2 && r9 < 0.94 ) dsigMC = 0.0326;
369 if ( !isEB && fabs(eta) < 2 && r9 >= 0.94 ) dsigMC = 0.0330;
370 if ( !isEB && fabs(eta) >= 2 && r9 < 0.94 ) dsigMC = 0.0331;
371 if ( !isEB && fabs(eta) >= 2 && r9 >= 0.94 ) dsigMC = 0.0378;
374 if ( isEB && fabs(eta) < 1 && r9 < 0.94 ) dsigMC = 0.0096;
375 if ( isEB && fabs(eta) < 1 && r9 >= 0.94 ) dsigMC = 0.0074;
376 if ( isEB && fabs(eta) >= 1 && r9 < 0.94 ) dsigMC = 0.0196;
377 if ( isEB && fabs(eta) >= 1 && r9 >= 0.94 ) dsigMC = 0.0141;
378 if ( !isEB && fabs(eta) < 2 && r9 < 0.94 ) dsigMC = 0.0279;
379 if ( !isEB && fabs(eta) < 2 && r9 >= 0.94 ) dsigMC = 0.0268;
380 if ( !isEB && fabs(eta) >= 2 && r9 < 0.94 ) dsigMC = 0.0301;
381 if ( !isEB && fabs(eta) >= 2 && r9 >= 0.94 ) dsigMC = 0.0293;
384 if ( isEB && fabs(eta) < 1 && r9 < 0.94 ) dsigMC = 0.0119;
385 if ( isEB && fabs(eta) < 1 && r9 >= 0.94 ) dsigMC = 0.0107;
386 if ( isEB && fabs(eta) >= 1 && r9 < 0.94 ) dsigMC = 0.0240;
387 if ( isEB && fabs(eta) >= 1 && r9 >= 0.94 ) dsigMC = 0.0149;
388 if ( !isEB && fabs(eta) < 2 && r9 < 0.94 ) dsigMC = 0.0330;
389 if ( !isEB && fabs(eta) < 2 && r9 >= 0.94 ) dsigMC = 0.0375;
390 if ( !isEB && fabs(eta) >= 2 && r9 < 0.94 ) dsigMC = 0.0602;
391 if ( !isEB && fabs(eta) >= 2 && r9 >= 0.94 ) dsigMC = 0.0607;
392 }
else if (
dataset_ ==
"Summer12_DR53X_HCP2012" ||
dataset_ ==
"Moriond2013" )
394 if ( isEB && fabs(eta) < 1 && r9 < 0.94 ) dsigMC = 0.0099;
395 if ( isEB && fabs(eta) < 1 && r9 >= 0.94 ) dsigMC = 0.0103;
396 if ( isEB && fabs(eta) >= 1 && r9 < 0.94 ) dsigMC = 0.0219;
397 if ( isEB && fabs(eta) >= 1 && r9 >= 0.94 ) dsigMC = 0.0158;
398 if ( !isEB && fabs(eta) < 2 && r9 < 0.94 ) dsigMC = 0.0222;
399 if ( !isEB && fabs(eta) < 2 && r9 >= 0.94 ) dsigMC = 0.0298;
400 if ( !isEB && fabs(eta) >= 2 && r9 < 0.94 ) dsigMC = 0.0318;
401 if ( !isEB && fabs(eta) >= 2 && r9 >= 0.94 ) dsigMC = 0.0302;
408 CLHEP::RandGaussQ gaussDistribution(rng->
getEngine(streamID), 1.,dsigMC);
409 corrMC = gaussDistribution.fire();
416 std::cout <<
"[ElectronEnergyCalibrator] "
417 <<
"======================= SYNCRONIZATION MODE! ======================="
437 std::cout <<
"[ElectronEnergyCalibrator] initial energy "
442 std::cout <<
"[ElectronEnergyCalibrator] initial energy error "
455 bool isEB = electron.
isEB();
459 double pt = p * fabs(
sin(theta));
461 double linscale = 0.;
471 if (classification<2)
480 if (classification<2)
490 if (classification<2)
500 double newP = p/(1.+linscale);
503 std::cout <<
"[ElectronEnergyCalibrator] Applying a linearity correction of " << 1./(1.+linscale) <<
" to " << pt <<
" GeV in pt" << std::endl;
508 std::cout <<
"[ElectronEnergyCalibrator] calibrated transverse momentum " << pt <<
" GeV recalibrated for linearity to momentum " << electron.
getCombinedMomentum()*fabs(
sin(theta)) <<
" GeV" << std::endl;
correctionValues corrValArray[100]
double getSCEnergy() const
void setNewEnergyError(double newEnergyError)
double getRegEnergy() const
void setCombinedMomentum(double combinedMomentum)
void setNewEnergy(double newEnergy)
Sin< T >::type sin(const T &t)
Geom::Theta< T > theta() const
void calibrate(SimpleElectron &electron, edm::StreamID const &)
virtual CLHEP::HepRandomEngine & getEngine(StreamID const &)=0
Use this engine in event methods.
double stringToDouble(const std::string &str)
double getCombinedMomentum() const
unsigned int getRunNumber() const
void correctLinearity(SimpleElectron &electron)
linearityCorrectionValues linCorrValArray[100]
double getRegEnergyError() const
bool applyLinearityCorrection_
double getSCEnergyError() const
void splitString(const std::string &fullstr, std::vector< std::string > &elements, const std::string &delimiter)