27 {
28
29
30
31ElectronChargeEfficiencyCorrectionTool::
32 ElectronChargeEfficiencyCorrectionTool(const std::string& name)
34 , m_dataTypeOverwrite(-1)
35 , m_eventInfoCollectionName("EventInfo")
36 , m_SF_SS()
37 , m_SF_OS()
38 , m_RunNumbers()
39 , m_useRandomRunNumber(true)
40 , m_defaultRandomRunNumber(999999)
41 , m_filename("")
42 , m_workingPoint("")
43 , m_eta_lowlimit(0.0)
44 , m_eta_uplimit(0.0)
45 , m_pt_lowlimit(0.0)
46 , m_pt_uplimit(0.0)
47 , m_gevmev(0.0)
48 , m_filtered_sys_sets()
49 , m_mySysConf()
50 , m_affectingSys()
51 , m_appliedSystematics(nullptr)
52 , m_sf_decoration_name("chargeIDEffiSF")
53 , m_sfDec(nullptr)
54{
55
56 declareProperty("CorrectionFileName",
57 m_filename,
58 "Name of the file with charge flipping rates");
59 declareProperty(
60 "WorkingPoint", m_workingPoint, "Name of working point folder in the file");
61 declareProperty("ScaleFactorDecorationName", m_sf_decoration_name);
62 declareProperty("ForceDataType",
63 m_dataTypeOverwrite,
64 "Force the DataType of the electron to specified value (to "
65 "circumvent problem of incorrect DataType for forward "
66 "electrons in some old releases)");
67 declareProperty("EventInfoCollectionName",
68 m_eventInfoCollectionName,
69 "The EventInfo Collection Name");
70 declareProperty("UseRandomRunNumber", m_useRandomRunNumber);
71 declareProperty("DefaultRandomRunNumber", m_defaultRandomRunNumber);
72}
73
74
75
76
77ElectronChargeEfficiencyCorrectionTool::
78 ~ElectronChargeEfficiencyCorrectionTool()
79{
80 if (m_sfDec)
81 delete m_sfDec;
82}
83
84
85
86
88ElectronChargeEfficiencyCorrectionTool::initialize()
89{
91
92
93
94
95
96 if (m_sfDec)
97 delete m_sfDec;
99
100
102 if (m_filename.empty()) {
103 ATH_MSG_ERROR(
" PathResolver was not able to find the file ... aborting");
104 return StatusCode::FAILURE;
105 }
106
107
108 TFile*
rootFile = TFile::Open(rootfilename.c_str());
109
110
111 if (rootFile == nullptr) {
113 << " ...... aborting");
114 return StatusCode::FAILURE;
115 }
116
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145 m_SF_SS.clear();
146 m_SF_OS.clear();
147 TList* keyListfolder =
rootFile->GetListOfKeys();
148 std::vector<std::string>
names;
149 std::set<std::string> set_systematics;
150
151 names.reserve(keyListfolder->GetEntries());
152 for (
int j = 0;
j < keyListfolder->GetEntries();
j++) {
153 names.emplace_back((keyListfolder->At(j)->GetName()));
154 }
156
157 for (const std::string & name : names) {
158
160 if (
name.contains(
"SFCentral_")) {
162
163 bool isOS = false;
164 if (
name.contains(
"_OS")) {
165 isOS = true;
167 }
168 if (isOS) {
169 std::string histid = (
name);
170 histid.erase(0, 10);
171 histid.erase(histid.size() - 3, 3);
173
174 if (histid.contains("RunNumber")) {
176 std::string runlow = histid;
177 runlow.erase(histid.find("RunNumber"), 9);
178 runlow.erase(runlow.find('_'), runlow.size());
179 m_RunNumbers.push_back(
180 static_cast<unsigned int>(
atoi(runlow.c_str())));
181 std::string runhigh = histid;
182 runhigh.erase(histid.find("RunNumber"), 9);
183 runhigh.erase(0, runhigh.find('_') + 1);
184 m_RunNumbers.push_back(
185 static_cast<unsigned int>(
atoi(runhigh.c_str())));
186 }
188 m_SF_OS[histid].push_back(
static_cast<TH2*
>(
rootFile->Get(
name.c_str())));
189 } else {
190 std::string histid = (
name);
191 histid.erase(0, 10);
192 histid.erase(histid.size() - 3, 3);
194 m_SF_SS[histid].push_back(
static_cast<TH2*
>(
rootFile->Get(
name.c_str())));
195 }
196 }
197
199 if (
name.contains(
"STAT_")) {
201 bool isOS = false;
202 if (
name.contains(
"_OS")) {
203 isOS = true;
205 }
206 if (isOS) {
207 std::string histid = (
name);
208 histid.erase(0, 5);
209 histid.erase(histid.size() - 3, 3);
212 m_SF_OS[histid].push_back(
static_cast<TH2*
>(
rootFile->Get(
name.c_str())));
213 } else {
214 std::string histid = (
name);
216 histid.erase(0, 5);
217 histid.erase(histid.size() - 3, 3);
219 m_SF_SS[histid].push_back(
static_cast<TH2*
>(
rootFile->Get(
name.c_str())));
220 }
221
222 }
223
225 if (
name.contains(
"SYST")) {
227 bool isOS = false;
228 if (
name.contains(
"_OS")) {
229 isOS = true;
231 }
232 if (isOS) {
233 std::string histid = (
name);
234 histid.erase(0, 4);
235 histid.erase(histid.size() - 3, 3);
236
237 std::string sysname = histid;
238 sysname.erase(sysname.find('_'), sysname.size());
239 set_systematics.insert(std::move(sysname));
240
241 histid.erase(0, histid.find('_') + 1);
244 m_SF_OS[histid].push_back(
static_cast<TH2*
>(
rootFile->Get(
name.c_str())));
245 } else {
246 std::string histid = (
name);
247 histid.erase(0, 4);
248 histid.erase(histid.size() - 3, 3);
249 histid.erase(0, histid.find('_') + 1);
251 m_SF_SS[histid].push_back(
static_cast<TH2*
>(
rootFile->Get(
name.c_str())));
252 }
253
254 }
255 }
256
258
259 if (m_SF_OS.empty() || m_SF_SS.empty() || m_SF_SS.size() != m_SF_OS.size()) {
261 "OS/SS SF vectors not filled or of different size. -- Problem with "
262 "files. -- Report to <hn-atlas-EGammaWG@cern.ch>");
263 return StatusCode::FAILURE;
264 }
265
266 m_systematics.insert(m_systematics.end(), set_systematics.begin(), set_systematics.end());
267
268 std::sort(m_RunNumbers.begin(), m_RunNumbers.end());
270
271
274 std::map<std::string, std::vector<TH2*>>::iterator
it = m_SF_OS.begin();
275
276
277 m_eta_lowlimit = (*it).second.at(0)->GetYaxis()->GetXmin();
278 m_eta_uplimit = (*it).second.at(0)->GetYaxis()->GetXmax();
279 ATH_MSG_VERBOSE(
"|eta| limits " << m_eta_lowlimit <<
", " << m_eta_uplimit);
280
281 m_pt_lowlimit = (*it).second.at(0)->GetXaxis()->GetXmin();
282 m_pt_uplimit = (*it).second.at(0)->GetXaxis()->GetXmax();
283 ATH_MSG_VERBOSE(
"pt limits " << m_pt_lowlimit <<
", " << m_pt_uplimit);
284
286 m_gevmev = 0.001;
287
288
289 m_affectingSys = affectingSystematics();
290
291
292 if (registerSystematics() != StatusCode::SUCCESS) {
293 ATH_MSG_ERROR(
"(registerSystematics() != CP::SystematicCode::Ok)");
294 return StatusCode::FAILURE;
295 }
296
297 return StatusCode::SUCCESS;
298}
299
300
301
302
303
304
305
306
308ElectronChargeEfficiencyCorrectionTool::getEfficiencyScaleFactor(
310 double& sf) const
311{
312
313
315
316
317
318 bool goodEle = false;
319 CP::CorrectionCode goodEle_result =
323 ATH_MSG_DEBUG(
"This is the check of goodeleCC in getscalefactor. Scale "
324 "factor set to -999");
325 return goodEle_result;
326 }
327
328 if (!goodEle) {
329
332 }
333
334
335 int reco_ele_charge = ele.
charge();
336 const double ele_pt = ele.
pt() * m_gevmev;
338
339
340 int truth_ele_charge = 9999;
341 CP::CorrectionCode charge_result =
345 ATH_MSG_VERBOSE(
"This is check of geteletruthchargeCC in getscalefactor. "
346 "Scale factor set to -9999");
347 return charge_result;
348 }
349
350 if (truth_ele_charge == 0) {
353 }
354
356 << "; Truth charge = " << truth_ele_charge);
357
358
359 float retVal(0.0);
360
362
363
364 std::string cutRunNumber = "all";
365
366 if (!m_RunNumbers.empty()) {
367 unsigned int runnumber = m_defaultRandomRunNumber;
369 << m_useRandomRunNumber);
370 if (m_useRandomRunNumber) {
372 evtStore()->retrieve<
const xAOD::EventInfo>(m_eventInfoCollectionName);
373 if (!eventInfo) {
377 }
379 "RandomRunNumber");
380 if (!randomrunnumber.isAvailable(*eventInfo)) {
383 "Pileup tool not run before using ElectronEfficiencyTool! SFs do not "
384 "reflect PU distribution in data");
386 }
387 runnumber = randomrunnumber(*(eventInfo));
388 }
389 ATH_MSG_DEBUG(
"Number of RunNumbers in file: " << m_RunNumbers.size());
390 for (std::size_t
r = 0;
r < m_RunNumbers.size();
r++) {
392 }
394
395 bool isInRunNumberRange = false;
396 for ( std::size_t
r=0;
r<m_RunNumbers.size()-1;
r+=2 ){
397
398
399 if (
runnumber >= (
unsigned int)m_RunNumbers.at(
r) &&
400 runnumber <= (
unsigned int)m_RunNumbers.at(
r+1) ) {
401 cutRunNumber.clear();
402 cutRunNumber =
403 Form(
"RunNumber%d_%d", m_RunNumbers.at(
r), m_RunNumbers.at(
r + 1));
404 ATH_MSG_DEBUG(
"Random run number lies in range " << m_RunNumbers.at(
r) <<
" " << m_RunNumbers.at(
r+1));
405 isInRunNumberRange = true;
406 }
407 }
408
410 (
runnumber > m_RunNumbers.at(m_RunNumbers.size() - 1))) {
414 }
415
416 if ( !isInRunNumberRange ) {
418 }
419 }
420
421
422 if ( ele_eta < m_eta_lowlimit || ele_eta > m_eta_uplimit ) {
423
424 ATH_MSG_DEBUG(
"Got an electron outside of the range of eta validity " << ele_eta);
426 }
427
428 if ( ele_pt < m_pt_lowlimit ) {
429
430 ATH_MSG_DEBUG(
"Got an electron outside of the range of pt validity: pt lower than lower limit");
432 }
433
434
435 const std::vector<TH2*>& SShistograms = m_SF_SS.at(cutRunNumber.c_str());
436 const std::vector<TH2*>& OShistograms = m_SF_OS.at(cutRunNumber.c_str());
437
438
439 bool isOS = false;
440
441 if (truth_ele_charge * reco_ele_charge > 0)
442 isOS = true;
443
444 if (isOS) {
445 retVal = this->getChargeFlipRate(ele_eta, ele_pt, OShistograms.at(0), sf);
446 if (retVal != 0) {
449 }
450 } else {
452 retVal = this->getChargeFlipRate(ele_eta, ele_pt, SShistograms.at(0), sf);
453 if (retVal != 0) {
456 }
457 }
458
461
462
463
464 double val_stat;
465
467 if (isOS) {
468 retVal =
469 this->getChargeFlipRate(ele_eta, ele_pt, OShistograms.at(1), val_stat);
470 if (retVal != 0) {
473 }
474 } else {
476 retVal =
477 this->getChargeFlipRate(ele_eta, ele_pt, SShistograms.at(1), val_stat);
478 if (retVal != 0) {
481 }
482 }
483
484 std::vector<float> systs;
485 double val_sys{ 0.0 };
487 for (
unsigned int s = 2;
s < OShistograms.size();
s++) {
488 if (isOS) {
489 retVal =
490 this->getChargeFlipRate(ele_eta, ele_pt, OShistograms.at(s), val_sys);
491 if (retVal != 0) {
492 val_sys = -9999.0;
494 }
495 } else {
497 retVal =
498 this->getChargeFlipRate(ele_eta, ele_pt, SShistograms.at(s), val_sys);
499 if (retVal != 0) {
500 val_sys = -9999.0;
502 }
503 }
504 systs.push_back(static_cast<float>(val_sys));
505 }
506
508
509 if (m_mySysConf.empty()) {
511 } else if (*(m_mySysConf.begin()) ==
512 CP::SystematicVariation("EL_CHARGEID_STAT", 1)) {
513 sf = (
sf + (val_stat));
515 } else if (*(m_mySysConf.begin()) ==
516 CP::SystematicVariation("EL_CHARGEID_STAT", -1)) {
517 sf = (
sf - (val_stat));
519 } else {
520
521 for (
unsigned int i = 0;
i < m_systematics.size();
i++) {
522 if (*(m_mySysConf.begin()) ==
523 CP::SystematicVariation(
524 Form("EL_CHARGEID_SYS%s", m_systematics.at(i).c_str()), 1)) {
525 sf = (
sf + (val_sys));
527 }
528
529 if (*(m_mySysConf.begin()) ==
530 CP::SystematicVariation(
531 Form("EL_CHARGEID_SYS%s", m_systematics.at(i).c_str()), -1)) {
532 sf = (
sf - (val_sys));
534 }
535 }
536
537 }
538
540}
541
542
543
544
545
547ElectronChargeEfficiencyCorrectionTool::applyEfficiencyScaleFactor(
549{
551 "In "
552 "ElectronChargeEfficiencyCorrectionTool::applyEfficiencyScaleFactor("
553 "const xAOD::IParticle& part) const");
555 CP::CorrectionCode
result = this->getEfficiencyScaleFactor(part, sf);
556
557 (*m_sfDec)(
part) =
static_cast<float>(sf);
559}
560
561
562float
563ElectronChargeEfficiencyCorrectionTool::getChargeFlipRate(
565 double pt,
566 TH2* hrates,
567 double& flipRate) const
568{
570 << " TH2, double&)");
571
572 if (pt > m_pt_uplimit)
573 pt = m_pt_uplimit * 0.999;
574
575 int bin2D = hrates->FindBin(pt,
eta);
576 flipRate = hrates->GetBinContent(bin2D);
577
579 << hrates->GetName());
580
581 return 0;
582}
583
585
586bool
587ElectronChargeEfficiencyCorrectionTool::isAffectedBySystematic(
589{
590
591 CP::SystematicSet
sys = affectingSystematics();
592 return sys.find(systematic) !=
sys.end();
593}
594
596
597
599ElectronChargeEfficiencyCorrectionTool::affectingSystematics() const
600{
602 result.insert(CP::SystematicVariation(
"EL_CHARGEID_STAT", 1));
603 result.insert(CP::SystematicVariation(
"EL_CHARGEID_STAT", -1));
604
605 for (
unsigned int i = 0;
i < m_systematics.size();
i++) {
606 result.insert(CP::SystematicVariation(
607 Form("EL_CHARGEID_SYS%s", m_systematics.at(i).c_str()), 1));
608 result.insert(CP::SystematicVariation(
609 Form("EL_CHARGEID_SYS%s", m_systematics.at(i).c_str()), -1));
610 }
612}
613
615
617ElectronChargeEfficiencyCorrectionTool::recommendedSystematics() const
618{
619
620 return affectingSystematics();
621}
622
624
626ElectronChargeEfficiencyCorrectionTool::applySystematicVariation(
628{
629
631 systConfig, m_affectingSys, m_mySysConf)) {
633 "Unsupported combination of systematics passed to the tool! ");
634 return StatusCode::FAILURE;
635 }
636
637 return StatusCode::SUCCESS;
638}
640
641StatusCode ElectronChargeEfficiencyCorrectionTool::registerSystematics() {
643
644 if (
registry.registerSystematics(*
this) != StatusCode::SUCCESS) {
645 ATH_MSG_ERROR(
"Failed to add systematic to list of recommended systematics.");
646 return StatusCode::FAILURE;
647 }
648
649 return StatusCode::SUCCESS;
650}
651
653}
Scalar eta() const
pseudorapidity method
#define ATH_MSG_DEBUG(x,...)
#define ATH_MSG_ERROR(x,...)
#define ATH_MSG_WARNING(x,...)
#define ATH_MSG_VERBOSE(x,...)
std::string PathResolverFindCalibFile(const std::string &logical_file_name)
Return value from object correction CP tools.
@ Error
Some error happened during the object correction.
@ OutOfValidityRange
Input object is out of validity range.
@ Ok
The correction was done successfully.
static SystematicRegistry & getInstance()
Get the singleton instance of the registry for the curren thread.
Class to wrap a set of SystematicVariations.
static StatusCode filterForAffectingSystematics(const SystematicSet &systConfig, const SystematicSet &affectingSystematics, SystematicSet &filteredSystematics)
description: filter the systematics for the affected systematics returns: success guarantee: strong f...
float etaBE(const unsigned layer) const
Get the eta in one layer of the EM Calo.
virtual double pt() const override final
The transverse momentum ( ) of the particle.
const xAOD::CaloCluster * caloCluster(size_t index=0) const
Pointer to the xAOD::CaloCluster/s that define the electron candidate.
float charge() const
Obtain the charge of the object.
static std::vector< uint32_t > runnumber
int atoi(std::string_view str)
Helper functions to unpack numbers decoded in string into integers and doubles The strings are requir...
::StatusCode StatusCode
StatusCode definition for legacy code.
CP::CorrectionCode getEleTruthCharge(const xAOD::Electron &ele, int &truthcharge)
Get the charge of the original electron.
CP::CorrectionCode isGoodEle(const xAOD::Electron &ele, bool &goodEle)
Return true if it's good ele for charge flip measurements.
SG::Accessor< T, ALLOC > Accessor
SG::Decorator< T, ALLOC > Decorator
Helper class to provide type-safe access to aux data, specialized for JaggedVecElt.
float j(const xAOD::IParticle &, const xAOD::TrackMeasurementValidation &hit, const Eigen::Matrix3d &jab_inv)
void sort(typename DataModel_detail::iterator< DVL > beg, typename DataModel_detail::iterator< DVL > end)
Specialization of sort for DataVector/List.
EventInfo_v1 EventInfo
Definition of the latest event info version.
Electron_v1 Electron
Definition of the current "egamma version".