ATLAS Offline Software
Loading...
Searching...
No Matches
EgammaCalibrationAndSmearingTool.cxx File Reference
#include <AsgServices/AsgServiceConfig.h>
#include <AsgTools/AsgToolConfig.h>
#include <algorithm>
#include <format>
#include <memory>
#include <string>
#include <utility>
#include "AthContainers/ConstAccessor.h"
#include "PATInterfaces/SystematicRegistry.h"
#include "PathResolver/PathResolver.h"
#include "xAODCaloEvent/CaloCluster.h"
#include "xAODEgamma/Egamma.h"
#include "xAODEgamma/EgammaDefs.h"
#include "xAODEgamma/EgammaxAODHelpers.h"
#include "xAODEventInfo/EventInfo.h"
#include "xAODTracking/TrackParticle.h"
#include "xAODTracking/TrackingPrimitives.h"
#include "xAODTracking/Vertex.h"
#include "EgammaAnalysisInterfaces/GlobalEventInfo.h"
#include "AthAnalysisBaseComps/AthAnalysisHelper.h"
#include <cmath>
#include "ElectronPhotonFourMomentumCorrection/EgammaCalibrationAndSmearingTool.h"
#include "ElectronPhotonFourMomentumCorrection/GainTool.h"
#include "ElectronPhotonFourMomentumCorrection/GainUncertainty.h"
#include "ElectronPhotonFourMomentumCorrection/LinearityADC.h"
#include "ElectronPhotonFourMomentumCorrection/egammaEnergyCorrectionTool.h"
#include "egammaLayerRecalibTool/egammaLayerRecalibTool.h"
#include "ElectronPhotonFourMomentumCorrection/systematics_es2024_Run3_v0.def"
#include "ElectronPhotonFourMomentumCorrection/systematics_S12_2022.def"
#include "ElectronPhotonFourMomentumCorrection/systematics_1NPCOR_PLUS_UNCOR.def"
#include "ElectronPhotonFourMomentumCorrection/systematics.def"
#include "ElectronPhotonFourMomentumCorrection/systematics_S12.def"

Go to the source code of this file.

Macros

#define SYSMACRO(name, fullcorrelated, decorrelation, flagup, flagdown)
#define SYSMACRO(name, fullcorrelated, decorrelation, flagup, flagdown)
#define SYSMACRO(name, fullcorrelated, decorrelation, flagup, flagdown)
#define SYSMACRO(name, fullcorrelated, decorrelation, flagup, flagdown)
#define SYSMACRO(name, fullcorrelated, decorrelation, flagup, flagdown)
#define SYSMACRO(name, fullcorrelated, decorrelation, flagup, flagdown)
#define SYSMACRO(name, fullcorrelated, decorrelation, flagup, flagdown)

Functions

namespace ADD_CMODE (CP)

Macro Definition Documentation

◆ SYSMACRO [1/7]

#define SYSMACRO ( name,
fullcorrelated,
decorrelation,
flagup,
flagdown )
Value:
if (bool(fullcorrelated)) { \
m_syst_description[CP::SystematicVariation(#name, +1)] = \
SysInfo{always, flagup}; \
m_syst_description[CP::SystematicVariation(#name, -1)] = \
SysInfo{always, flagdown}; \
} else { \
int i = 0; \
for (const auto& p : AbsEtaCaloPredicatesFactory(decorrelation)) { \
m_syst_description[CP::SystematicVariation( \
#name "__ETABIN" + std::to_string(i), +1)] = SysInfo{p, flagup}; \
m_syst_description[CP::SystematicVariation( \
#name "__ETABIN" + std::to_string(i), -1)] = SysInfo{p, flagdown}; \
i += 1; \
} \
}

◆ SYSMACRO [2/7]

#define SYSMACRO ( name,
fullcorrelated,
decorrelation,
flagup,
flagdown )
Value:
if (bool(fullcorrelated)) { \
m_syst_description[CP::SystematicVariation(#name, +1)] = \
SysInfo{always, flagup}; \
m_syst_description[CP::SystematicVariation(#name, -1)] = \
SysInfo{always, flagdown}; \
} else { \
int i = 0; \
for (const auto& p : AbsEtaCaloPredicatesFactory(decorrelation)) { \
m_syst_description[CP::SystematicVariation( \
#name "__ETABIN" + std::to_string(i), +1)] = SysInfo{p, flagup}; \
m_syst_description[CP::SystematicVariation( \
#name "__ETABIN" + std::to_string(i), -1)] = SysInfo{p, flagdown}; \
i += 1; \
} \
}

◆ SYSMACRO [3/7]

#define SYSMACRO ( name,
fullcorrelated,
decorrelation,
flagup,
flagdown )
Value:
if (bool(fullcorrelated)) { \
m_syst_description[CP::SystematicVariation(#name, +1)] = \
SysInfo{always, flagup}; \
m_syst_description[CP::SystematicVariation(#name, -1)] = \
SysInfo{always, flagdown}; \
} else { \
int i = 0; \
for (const auto& p : AbsEtaCaloPredicatesFactory(decorrelation)) { \
m_syst_description[CP::SystematicVariation( \
#name "__ETABIN" + std::to_string(i), +1)] = SysInfo{p, flagup}; \
m_syst_description[CP::SystematicVariation( \
#name "__ETABIN" + std::to_string(i), -1)] = SysInfo{p, flagdown}; \
i += 1; \
} \
}

◆ SYSMACRO [4/7]

#define SYSMACRO ( name,
fullcorrelated,
decorrelation,
flagup,
flagdown )
Value:
if (bool(fullcorrelated)) { \
m_syst_description[CP::SystematicVariation(#name, +1)] = \
SysInfo{always, flagup}; \
m_syst_description[CP::SystematicVariation(#name, -1)] = \
SysInfo{always, flagdown}; \
} else { \
int i = 0; \
for (const auto& p : AbsEtaCaloPredicatesFactory(decorrelation)) { \
m_syst_description[CP::SystematicVariation( \
#name "__ETABIN" + std::to_string(i), +1)] = SysInfo{p, flagup}; \
m_syst_description[CP::SystematicVariation( \
#name "__ETABIN" + std::to_string(i), -1)] = SysInfo{p, flagdown}; \
i += 1; \
} \
}

◆ SYSMACRO [5/7]

#define SYSMACRO ( name,
fullcorrelated,
decorrelation,
flagup,
flagdown )
Value:
m_syst_description[CP::SystematicVariation(#name, +1)] = \
SysInfo{always, flagup}; \
m_syst_description[CP::SystematicVariation(#name, -1)] = \
SysInfo{always, flagdown};

◆ SYSMACRO [6/7]

#define SYSMACRO ( name,
fullcorrelated,
decorrelation,
flagup,
flagdown )
Value:
m_syst_description[CP::SystematicVariation(#name, +1)] = \
SysInfo{always, flagup}; \
m_syst_description[CP::SystematicVariation(#name, -1)] = \
SysInfo{always, flagdown};

◆ SYSMACRO [7/7]

#define SYSMACRO ( name,
fullcorrelated,
decorrelation,
flagup,
flagdown )
Value:
m_syst_description[CP::SystematicVariation(#name, +1)] = \
SysInfo{always, flagup}; \
m_syst_description[CP::SystematicVariation(#name, -1)] = \
SysInfo{always, flagdown};

Function Documentation

◆ ADD_CMODE()

namespace ADD_CMODE ( CP )

Definition at line 1 of file EgammaCalibrationAndSmearingTool.cxx.

41 {
42
43const double GeV = 1000.;
44
45std::unique_ptr<egGain::GainTool> gainToolFactory(egEnergyCorr::ESModel model) {
46 switch (model) {
59 const std::string gain_filename1 = PathResolverFindCalibFile(
60 "ElectronPhotonFourMomentumCorrection/v8/FunctionsTO.root");
61 const std::string gain_filename2 = PathResolverFindCalibFile(
62 "ElectronPhotonFourMomentumCorrection/v8/FunctionsG_all.root");
63 return std::make_unique<egGain::GainTool>(gain_filename1, gain_filename2);
64 }
80 return nullptr;
81 default:
82 return nullptr;
83 }
84}
85
86std::string egammaMVAToolFolder(egEnergyCorr::ESModel model) {
87 std::string folder;
88 switch (model) {
92 folder = "egammaMVACalib/v1";
93 break;
95 folder = "egammaMVACalib/v1";
96 break;
101 folder = "egammaMVACalib/offline/v3";
102 break;
105 folder = "egammaMVACalib/offline/v3_E4crack_bis";
106 break;
114 folder = "egammaMVACalib/offline/v4.0";
115 break;
124 folder = "egammaMVACalib/offline/v7";
125 break;
128 folder = "egammaMVACalib/offline/v9";
129 break;
131 folder = "egammaMVACalib/offline/v10";
132 break;
133 default:
134 folder = "";
135 }
136
137 return folder;
138}
139
140std::unique_ptr<egammaLayerRecalibTool> egammaLayerRecalibToolFactory(
141 egEnergyCorr::ESModel model, int enableSacc) {
142 std::string tune = "";
143 switch (model) {
147 tune = "2011_alt_with_layer2";
148 break;
161 tune = "2012_alt_with_layer2";
162 break;
165 tune = "es2017_20.7_improved";
166 break;
168 tune = "es2017_20.7_final";
169 break;
176 tune = "es2017_21.0_v0";
177 break;
179 tune = "es2018_21.0_v0";
180 break;
182 tune = "es2022_22.0_Precision";
183 break;
185 tune = "es2022_22.0_Precision_v1";
186 break;
188 tune = "es2024_run3_extrapolate_v0";
189 break;
190 default:
191 return nullptr;
192 }
193 return std::make_unique<egammaLayerRecalibTool>(tune, enableSacc);
194}
195
196bool use_intermodule_correction(egEnergyCorr::ESModel model) {
197 switch (model) {
204 return false;
233 return true;
234 case egEnergyCorr::UNDEFINED: // TODO: find better logic
235 return false;
236 }
237 assert(false);
238 return false;
239}
240
241bool use_phi_uniform_correction(egEnergyCorr::ESModel model) {
242 return use_intermodule_correction(model); // they are equal
243}
244
245bool is_after_run1(egEnergyCorr::ESModel model) {
246 switch (model) {
256 return false;
282 return true;
283 case egEnergyCorr::UNDEFINED: // TODO: find better logic
284 return false;
285 }
286 assert(false);
287 return false;
288}
289
290EgammaCalibrationAndSmearingTool::EgammaCalibrationAndSmearingTool(
291 const std::string& name)
292 : asg::AsgMetadataTool(name),
293 m_TESModel(egEnergyCorr::UNDEFINED),
294 m_TResolutionType(egEnergyCorr::Resolution::SigmaEff90),
295 m_use_mapping_correction(false),
296 m_currentScaleVariation_MC(egEnergyCorr::Scale::None),
297 m_currentScaleVariation_data(egEnergyCorr::Scale::Nominal),
298 m_currentResolutionVariation_MC(egEnergyCorr::Resolution::Nominal),
299 m_currentResolutionVariation_data(egEnergyCorr::Resolution::None),
300 m_set_seed_function([](const EgammaCalibrationAndSmearingTool& tool,
302 columnar::EventInfoId ei) {
303 const Accessors& acc = *tool.m_accessors;
304 // avoid 0 as result, see
305 // https://root.cern.ch/root/html/TRandom3.html#TRandom3:SetSeed
306 auto cluster = acc.caloClusterAcc(egamma)[0].value();
307 return 1 + static_cast<RandomNumber>(
308 std::abs(acc.clusterPhiAcc(cluster)) * 1E6 +
309 std::abs(acc.clusterEtaAcc(cluster)) * 1E3 +
310 acc.eventNumberAcc(ei));
311 }),
312 m_accessors(std::make_unique<Accessors>(*this)) {
313
314 declareProperty("ESModel", m_ESModel = "");
315 declareProperty("decorrelationModel", m_decorrelation_model_name = "");
316 declareProperty("decorrelationModelScale",
317 m_decorrelation_model_scale_name = "");
318 declareProperty("decorrelationModelResolution",
319 m_decorrelation_model_resolution_name = "");
320 declareProperty("ResolutionType", m_ResolutionType = "SigmaEff90");
321 declareProperty("varSF", m_varSF = 1.0);
322 declareProperty("doScaleCorrection", m_doScaleCorrection = AUTO);
323 declareProperty("doSmearing", m_doSmearing = AUTO);
324 declareProperty("useLayerCorrection", m_useLayerCorrection = AUTO);
325 declareProperty("usePSCorrection", m_usePSCorrection = AUTO);
326 declareProperty("useS12Correction", m_useS12Correction = AUTO);
327 declareProperty("useSaccCorrection", m_useSaccCorrection = AUTO);
328 declareProperty("useIntermoduleCorrection",
329 m_useIntermoduleCorrection = AUTO);
330 declareProperty("usePhiUniformCorrection", m_usePhiUniformCorrection = AUTO);
331 declareProperty("useCaloDistPhiUnifCorrection",
332 m_useCaloDistPhiUnifCorrection = AUTO);
333 declareProperty("useGainCorrection", m_useGainCorrection = AUTO);
334 declareProperty("useGainInterpolation", m_useGainInterpolation = AUTO);
335 declareProperty("doADCLinearityCorrection",
336 m_doADCLinearityCorrection = AUTO);
337 declareProperty("doLeakageCorrection", m_doLeakageCorrection = AUTO);
338 declareProperty("MVAfolder", m_MVAfolder = "");
339 declareProperty("layerRecalibrationTune", m_layer_recalibration_tune = "");
340 declareProperty("useEPCombination", m_use_ep_combination = false);
341 declareProperty("useMVACalibration", m_use_mva_calibration = AUTO);
342 declareProperty("use_full_statistical_error",
343 m_use_full_statistical_error = false);
344 declareProperty("use_temp_correction201215",
345 m_use_temp_correction201215 = AUTO);
346 declareProperty("use_uA2MeV_2015_first2weeks_correction",
347 m_use_uA2MeV_2015_first2weeks_correction = AUTO);
348 declareProperty("randomRunNumber", m_user_random_run_number = 0);
349 // this is the user input, it is never changed by the tool. The tool uses
350 // m_simulation.
351 declareProperty("useFastSim", m_useFastSim = -1,
352 "This should be explicitly set by the user depending on the "
353 "data type (int)0=full sim, (int)1=fast sim");
354 declareProperty(
355 "useAFII", m_use_AFII = -1,
356 "This is now deprecated. Kept for explicit error message for now");
357 declareProperty("decorateEmva", m_decorateEmva = false, "whether to decorate the eMVA value");
358}
359
360EgammaCalibrationAndSmearingTool::~EgammaCalibrationAndSmearingTool() {
361 ATH_MSG_DEBUG("destructor");
362 delete m_layer_recalibration_tool;
363 delete m_gain_tool;
364}
365
366StatusCode EgammaCalibrationAndSmearingTool::initialize() {
367 ATH_MSG_INFO("Initialization");
368
369 if (m_ESModel == "es2015XX") {
370 ATH_MSG_ERROR("es2015XX is deprecated. Use es2015PRE");
371 }
372
373 if (m_ESModel == "es2010") {
374 m_TESModel = egEnergyCorr::es2010;
375 } // legacy
376 else if (m_ESModel == "es2011c") {
377 m_TESModel = egEnergyCorr::es2011c;
378 } // mc11c : faulty G4; old geometry
379 else if (m_ESModel == "es2011d") {
380 m_TESModel = egEnergyCorr::es2011d;
381 } // mc11d : corrected G4; new geometry == final Run1 scheme
382 else if (m_ESModel == "es2012a") {
383 m_TESModel = egEnergyCorr::es2012a;
384 } // mc12a : "crude" G4 fix; old geometry
385 else if (m_ESModel == "es2012c") {
386 m_TESModel = egEnergyCorr::es2012c;
387 } // mc12c : corrected G4; new geometry == final Run1 scheme
388 else if (m_ESModel == "es2012XX") {
389 m_TESModel = egEnergyCorr::es2012XX;
390 } else if (m_ESModel == "es2015PRE") {
391 m_TESModel = egEnergyCorr::es2015PRE;
392 } else if (m_ESModel == "es2015PRE_res_improved") {
394 } else if (m_ESModel == "es2015cPRE") {
395 m_TESModel = egEnergyCorr::es2015cPRE;
396 } else if (m_ESModel == "es2015cPRE_res_improved") {
398 } else if (m_ESModel == "es2015c_summer") {
399 m_TESModel = egEnergyCorr::es2015c_summer;
400 } else if (m_ESModel == "es2016PRE") {
401 m_TESModel = egEnergyCorr::es2016PRE;
402 } else if (m_ESModel == "es2016data_mc15c") {
403 m_TESModel = egEnergyCorr::es2017;
404 } else if (m_ESModel == "es2016data_mc15c_summer") {
405 m_TESModel = egEnergyCorr::es2017_summer;
406 } else if (m_ESModel == "es2016data_mc15c_summer_improved") {
408 } else if (m_ESModel == "es2016data_mc15c_final") {
410 } else if (m_ESModel == "es2015_5TeV") {
411 m_TESModel = egEnergyCorr::es2015_5TeV;
412 } else if (m_ESModel == "es2017_R21_PRE") {
413 m_TESModel = egEnergyCorr::es2017_R21_PRE;
414 } else if (m_ESModel == "es2017_R21_v0") {
415 m_TESModel = egEnergyCorr::es2017_R21_v0;
416 } else if (m_ESModel == "es2017_R21_v1") {
417 m_TESModel = egEnergyCorr::es2017_R21_v1;
418 } else if (m_ESModel == "es2017_R21_ofc0_v1") {
420 } else if (m_ESModel == "es2018_R21_v0") {
421 m_TESModel = egEnergyCorr::es2018_R21_v0;
422 } else if (m_ESModel == "es2018_R21_v1") {
423 m_TESModel = egEnergyCorr::es2018_R21_v1;
424 } else if (m_ESModel == "es2022_R22_PRE") {
425 m_TESModel = egEnergyCorr::es2022_R22_PRE;
426 } else if (m_ESModel == "es2023_R22_Run2_v0") {
428 } else if (m_ESModel == "es2023_R22_Run2_v1") {
430 } else if (m_ESModel == "es2024_Run3_ofc0_v0") {
432 } else if (m_ESModel == "es2024_Run3_v0") {
433 m_TESModel = egEnergyCorr::es2024_Run3_v0;
434 } else if (m_ESModel == "es2025_Run3_GNN_v0") {
436 } else if (m_ESModel.empty()) {
437 ATH_MSG_ERROR("you must set ESModel property");
438 return StatusCode::FAILURE;
439 } else {
440 ATH_MSG_ERROR("Cannot understand model " << m_ESModel);
441 return StatusCode::FAILURE;
442 }
443
444 if (m_ResolutionType == "Gaussian") {
445 m_TResolutionType = egEnergyCorr::Resolution::Gaussian;
446 } else if (m_ResolutionType == "SigmaEff80") {
447 m_TResolutionType = egEnergyCorr::Resolution::SigmaEff80;
448 } else if (m_ResolutionType == "SigmaEff90") {
449 m_TResolutionType = egEnergyCorr::Resolution::SigmaEff90;
450 } else {
451 ATH_MSG_ERROR("Cannot understand resolution " << m_ResolutionType);
452 return StatusCode::FAILURE;
453 }
454
455 if (m_use_AFII != -1) {
457 "Property useAFII is deprecated. It is now replaced with useFastSim, "
458 "which should be explicitly configured");
459 return StatusCode::FAILURE;
460 }
461
462 if (m_useFastSim == 1) {
463 m_simulation = PATCore::ParticleDataType::Fast;
464 } else if (m_useFastSim == 0) {
465 m_simulation = PATCore::ParticleDataType::Full;
466 } else {
467 ATH_MSG_ERROR("Property useFastSim should be explicitly configured");
468 return StatusCode::FAILURE;
469 }
470
471 if ( (m_TESModel == egEnergyCorr::es2022_R22_PRE || m_TESModel == egEnergyCorr::es2024_Run3_ofc0_v0) &&
472 m_simulation == PATCore::ParticleDataType::Fast) {
474 "Sample is FastSim but no AF3 calibration is supported with "
475 "MC23 pre-recommendations (es2022_R22_PRE and es2024_Run3_ofc0_v0). "
476 "Please swtich to Run3 consolidated recommendations (es2024_Run3_v0), "
477 "or get in touch with the EGamma CP group in case you are using this");
478 return StatusCode::FAILURE;
479 }
480
481 // configure decorrelation model, translate string property to internal class
482 // enum
483 /* S R SR
484 0. 0 0 0 WARNING Full, Full (this is the default without configuration)
485 1. 0 0 1 SR
486 2. 0 1 0 FATAL
487 3. 0 1 1 WARNING SR then R
488 4. 1 0 0 FATAL
489 5. 1 0 1 WARNING SR then S
490 6. 1 1 0 S, R
491 7. 1 1 1 FATAL
492 */
493 if (m_decorrelation_model_name.empty() and
494 m_decorrelation_model_scale_name.empty() and
495 m_decorrelation_model_resolution_name.empty()) {
496 // case 0
497 ATH_MSG_WARNING("no decorrelation model specified, assuming full model");
498 m_decorrelation_model_scale = ScaleDecorrelation::FULL;
499 m_decorrelation_model_resolution = ResolutionDecorrelation::FULL;
500 m_decorrelation_model_name = "FULL_v1";
501 } else if (not m_decorrelation_model_name.empty() and
502 not m_decorrelation_model_scale_name.empty() and
503 not m_decorrelation_model_resolution_name.empty()) {
504 // case 7
505 ATH_MSG_FATAL("too many flags for the decorrelation model");
506 return StatusCode::FAILURE;
507 } else {
508 // set scale decorrelation model
509 if (not m_decorrelation_model_scale_name.empty()) { // case 4, 5, 6, (7)
510 if (not m_decorrelation_model_name.empty()) {
512 "flag decorrelation model ignored for scale decorrelation model");
513 } // case 5
514 if (m_decorrelation_model_scale_name == "1NP_v1")
515 m_decorrelation_model_scale = ScaleDecorrelation::ONENP;
516 else if (m_decorrelation_model_scale_name == "FULL_ETACORRELATED_v1")
517 m_decorrelation_model_scale = ScaleDecorrelation::FULL_ETA_CORRELATED;
518 else if (m_decorrelation_model_scale_name == "1NPCOR_PLUS_UNCOR")
519 m_decorrelation_model_scale = ScaleDecorrelation::ONENP_PLUS_UNCONR;
520 else if (m_decorrelation_model_scale_name == "FULL_v1")
521 m_decorrelation_model_scale = ScaleDecorrelation::FULL;
522 else {
523 ATH_MSG_FATAL("cannot understand the scale decorrelation model '"
524 << m_decorrelation_model_scale_name << "'(typo?)");
525 return StatusCode::FAILURE;
526 }
527 } else if (not m_decorrelation_model_name.empty()) { // case 1, 3
528 if (m_decorrelation_model_name == "1NP_v1")
529 m_decorrelation_model_scale = ScaleDecorrelation::ONENP;
530 else if (m_decorrelation_model_name == "FULL_ETACORRELATED_v1")
531 m_decorrelation_model_scale = ScaleDecorrelation::FULL_ETA_CORRELATED;
532 else if (m_decorrelation_model_name == "1NPCOR_PLUS_UNCOR")
533 m_decorrelation_model_scale = ScaleDecorrelation::ONENP_PLUS_UNCONR;
534 else if (m_decorrelation_model_name == "FULL_v1")
535 m_decorrelation_model_scale = ScaleDecorrelation::FULL;
536 else {
537 ATH_MSG_FATAL("cannot understand the decorrelation model '"
538 << m_decorrelation_model_name << "'(typo?)");
539 return StatusCode::FAILURE;
540 }
541 } else { // case 2, (7)
543 "not information how to initialize the scale decorrelation model");
544 return StatusCode::FAILURE;
545 }
546
547 // set resolution decorralation model
548 if (not m_decorrelation_model_resolution_name
549 .empty()) { // case 2, 3, 6, (7)
550 if (not m_decorrelation_model_name.empty()) {
552 "flag decorrelation model ignored for resolution decorrelation "
553 "model");
554 } // case 3
555 if (m_decorrelation_model_resolution_name == "1NP_v1")
556 m_decorrelation_model_resolution = ResolutionDecorrelation::ONENP;
557 else if (m_decorrelation_model_resolution_name == "FULL_v1")
558 m_decorrelation_model_resolution = ResolutionDecorrelation::FULL;
559 else {
560 ATH_MSG_FATAL("cannot understand the resolution decorrelation model '"
561 << m_decorrelation_model_resolution_name << "'(typo?)");
562 return StatusCode::FAILURE;
563 }
564 } else if (not m_decorrelation_model_name.empty()) { // case 1, 5
565 if (m_decorrelation_model_name == "1NP_v1")
566 m_decorrelation_model_resolution = ResolutionDecorrelation::ONENP;
567 else if (m_decorrelation_model_name == "FULL_ETACORRELATED_v1")
568 m_decorrelation_model_resolution = ResolutionDecorrelation::FULL;
569 else if (m_decorrelation_model_name == "1NPCOR_PLUS_UNCOR")
570 m_decorrelation_model_resolution = ResolutionDecorrelation::ONENP;
571 else if (m_decorrelation_model_name == "FULL_v1")
572 m_decorrelation_model_resolution = ResolutionDecorrelation::FULL;
573 else {
574 ATH_MSG_FATAL("cannot understand the decorrelation model '"
575 << m_decorrelation_model_name << "'(typo?)");
576 return StatusCode::FAILURE;
577 }
578 }
579 }
580
581 // create correction tool
582 ATH_MSG_DEBUG("creating internal correction tool");
583 m_rootTool = std::make_unique<AtlasRoot::egammaEnergyCorrectionTool>();
584 if (!m_rootTool) {
585 ATH_MSG_ERROR("Cannot initialize underlying tool");
586 return StatusCode::FAILURE;
587 }
588 m_rootTool->setESModel(m_TESModel);
589
590 if ( (m_TESModel == egEnergyCorr::es2023_R22_Run2_v0 ||
591 m_TESModel == egEnergyCorr::es2023_R22_Run2_v1 ||
592 m_TESModel == egEnergyCorr::es2024_Run3_v0) &&
593 (m_useGainInterpolation == AUTO || m_useGainInterpolation == 1)) {
595 "Using linear interpolation in the gain tool (uncertainties only)");
596 m_useGainInterpolation = 1;
597 m_rootTool->setApplyL2GainInterpolation();
598 }
599 m_rootTool->msg().setLevel(this->msg().level());
600 m_rootTool->initialize();
601
602 // configure MVA calibration
603 if (m_use_mva_calibration != 0) {
604 ATH_MSG_DEBUG("creating MVA calibration tool (if needed)");
605 if (m_MVAfolder.empty()) { // automatically configure MVA tool
606 m_MVAfolder = egammaMVAToolFolder(m_TESModel);
607 }
608
609 if (not m_MVAfolder.empty()) {
610
611 // electron MVA tool
612 asg::AsgToolConfig config_mva_electron(
613 "egammaMVACalibTool/tool_mva_electron");
614 config_mva_electron.setPropertyFromString("folder", m_MVAfolder);
615 ATH_CHECK(config_mva_electron.setProperty("use_layer_corrected", true));
616 ATH_CHECK(config_mva_electron.setProperty(
617 "ParticleType", xAOD::EgammaParameters::electron));
618
619 // unconverted photon MVA tool
620 asg::AsgToolConfig config_mva_unconverted(
621 "egammaMVACalibTool/tool_mva_unconverted");
622 config_mva_unconverted.setPropertyFromString("folder", m_MVAfolder);
623 ATH_CHECK(
624 config_mva_unconverted.setProperty("use_layer_corrected", true));
625 ATH_CHECK(config_mva_unconverted.setProperty(
627 ATH_CHECK(config_mva_unconverted.setProperty("OutputLevel",
628 this->msg().level()));
629
630 // converted photon MVA tool
631 asg::AsgToolConfig config_mva_converted(
632 "egammaMVACalibTool/tool_mva_converted");
633 config_mva_converted.setPropertyFromString("folder", m_MVAfolder);
634 ATH_CHECK(config_mva_converted.setProperty("use_layer_corrected", true));
635 ATH_CHECK(config_mva_converted.setProperty(
637 ATH_CHECK(config_mva_converted.setProperty("OutputLevel",
638 this->msg().level()));
639
640 // initialize the ServiceHandler egammaMVASvc
641 // make the name unique
642 std::ostringstream mva_service_name;
643 size_t dot = this->name().find('.');
644 size_t underscore = this->name().rfind('_', dot);
645 std::string collection = "";
646 if (dot != std::string::npos &&
647 underscore != std::string::npos &&
648 underscore < dot)
649 {
650 collection = this->name().substr(underscore + 1, dot - underscore - 1);
651 mva_service_name << "egammaMVASvc/service_mva_egamma_id"
652 << collection;;
653 }
654 else
655 {
656 mva_service_name << "egammaMVASvc/service_mva_egamma_id"
657 << (void const*)this;
658 }
659 asg::AsgServiceConfig config_mva_service(mva_service_name.str());
660 // each tool loads its BDTs at initialize; egammaMVASvc disables missing ones
661 if (m_doElectronMVA) {
662 ATH_CHECK(config_mva_service.addPrivateTool("ElectronTool",
663 std::move(config_mva_electron)));
664 }
665 if (m_doPhotonMVA) {
666 ATH_CHECK(config_mva_service.addPrivateTool("UnconvertedPhotonTool",
667 std::move(config_mva_unconverted)));
668 ATH_CHECK(config_mva_service.addPrivateTool("ConvertedPhotonTool",
669 std::move(config_mva_converted)));
670 }
671 // fwd electron MVA tool
672 if (m_doFwdCalib) {
673 asg::AsgToolConfig config_mva_fwdelectron(
674 "egammaMVACalibTool/tool_mva_fwdelectron");
675 config_mva_fwdelectron.setPropertyFromString("folder", m_MVAfolder);
676 ATH_CHECK(config_mva_fwdelectron.setProperty(
678 ATH_CHECK(config_mva_fwdelectron.setProperty("ShiftType", 0));
679 ATH_CHECK(config_mva_fwdelectron.setProperty("OutputLevel", this->msg().level()));
680 ATH_CHECK(config_mva_service.addPrivateTool("FwdElectronTool",
681 std::move(config_mva_fwdelectron)));
682 }
683 config_mva_service.setPropertyFromString("folder", m_MVAfolder);
684 ATH_CHECK(
685 config_mva_service.setProperty("OutputLevel", this->msg().level()));
686 ATH_CHECK(config_mva_service.makeService(m_MVACalibSvc));
687 } else if (m_TESModel == egEnergyCorr::es2025_Run3_GNN_v0) {
689 "WIP: testing GNN based calibration for Run3,"
690 "requring decorated GNN energy (gnn_energy) from input,"
691 "which should already have the layer calibration applied");
692 } else {
693 m_use_mva_calibration = false;
694 }
695 }
696
697 // configure layer recalibration tool
698 // For now: layer recalibration not applied to PRE release 21 (using run 1
699 // based calibration applied at reco level)
700 // for following R21 recommendations, need to apply the run2/run1 layer
701 // calibration ratio
702 if (m_ESModel == "es2017_R21_PRE") {
703 ATH_MSG_INFO("Layer recalibration already applied at cell level");
704 m_useLayerCorrection = false;
705 } else if (!m_useLayerCorrection) {
706 ATH_MSG_INFO("Layer corrections disabled!");
707 } else {
708 ATH_MSG_DEBUG("initializing layer recalibration tool (if needed)");
709 if (m_layer_recalibration_tune
710 .empty()) { // automatically configure layer recalibration tool
711 m_layer_recalibration_tool =
712 egammaLayerRecalibToolFactory(m_TESModel, m_useSaccCorrection)
713 .release();
714 if (!m_layer_recalibration_tool) {
715 ATH_MSG_INFO("not using layer recalibration");
716 }
717 } else {
718 m_layer_recalibration_tool =
719 new egammaLayerRecalibTool(m_layer_recalibration_tune, m_useSaccCorrection);
720 }
721 if (m_layer_recalibration_tool) {
722 m_layer_recalibration_tool->msg().setLevel(this->msg().level());
723 m_layer_recalibration_tool->fixForMissingCells(m_fixForMissingCells);
724 if (!m_usePSCorrection) {
725 ATH_MSG_INFO("PS corrections disabled!");
726 m_layer_recalibration_tool->disable_PSCorrections();
727 }
728 if (!m_useS12Correction) {
729 ATH_MSG_INFO("S12 corrections disabled!");
730 m_layer_recalibration_tool->disable_S12Corrections();
731 }
732 if (!m_useSaccCorrection) {
733 ATH_MSG_INFO("Sacc corrections disabled!");
734 m_layer_recalibration_tool->disable_SaccCorrections();
735 }
736 }
737 }
738
739 if (m_use_temp_correction201215 != AUTO)
740 m_rootTool->use_temp_correction201215(m_use_temp_correction201215);
741 if (m_use_uA2MeV_2015_first2weeks_correction != AUTO)
742 m_rootTool->use_uA2MeV_2015_first2weeks_correction(
743 m_use_uA2MeV_2015_first2weeks_correction);
744 if (not m_use_full_statistical_error and
745 m_decorrelation_model_scale == ScaleDecorrelation::FULL) {
746 m_rootTool->useStatErrorScaling(true);
747 }
748
749 if (m_use_ep_combination) {
750 ATH_MSG_ERROR("ep combination not supported yet");
751 throw std::runtime_error("ep combination not supported yet");
752 }
753
754 if (m_useIntermoduleCorrection == AUTO) {
755 m_useIntermoduleCorrection = use_intermodule_correction(m_TESModel);
756 }
757 if (m_usePhiUniformCorrection == AUTO) {
758 m_usePhiUniformCorrection = use_phi_uniform_correction(m_TESModel);
759 }
760 m_use_mapping_correction = not is_after_run1(m_TESModel);
761
762 if (m_useGainCorrection == AUTO) {
763 if (m_TESModel == egEnergyCorr::es2023_R22_Run2_v0 ||
764 m_TESModel == egEnergyCorr::es2023_R22_Run2_v1 ||
765 m_TESModel == egEnergyCorr::es2024_Run3_v0) {
766 m_useGainCorrection = 0;
767 }
768 else {
769 ATH_MSG_DEBUG("initializing gain tool");
770 m_gain_tool = gainToolFactory(m_TESModel).release();
771 m_useGainCorrection = bool(m_gain_tool);
772 }
773 }
774 else if (m_useGainCorrection == 1) {
775 if (m_TESModel != egEnergyCorr::es2023_R22_Run2_v0 &&
776 m_TESModel != egEnergyCorr::es2023_R22_Run2_v1 &&
777 m_TESModel != egEnergyCorr::es2024_Run3_v0) {
778 m_useGainCorrection = 0;
780 "cannot instantiate gain tool for this model (you can only disable "
781 "the gain tool, but not enable it)");
782 }
783 else {
785 "initializing gain tool for run2 final precision recommendations");
787 "Gain corrections required but Zee scales are derived without Gain, "
788 "will cause inconsistency!");
789 std::string gain_tool_run_2_filename = PathResolverFindCalibFile(
790 "ElectronPhotonFourMomentumCorrection/v29/"
791 "gain_uncertainty_specialRun.root");
792 m_gain_tool_run2 = std::make_unique<egGain::GainUncertainty>(
793 gain_tool_run_2_filename, false, "GainCorrection",
794 m_useGainInterpolation);
795 m_gain_tool_run2->msg().setLevel(this->msg().level());
796 }
797 }
798
799 if (m_TESModel == egEnergyCorr::es2023_R22_Run2_v0 ||
800 m_TESModel == egEnergyCorr::es2023_R22_Run2_v1 ||
801 m_TESModel == egEnergyCorr::es2024_Run3_v0 ||
802 m_TESModel == egEnergyCorr::es2025_Run3_GNN_v0) {
803 // ADC non linearity correction
804 if (m_doADCLinearityCorrection == AUTO || m_doADCLinearityCorrection == 1) {
805 m_doADCLinearityCorrection = 1;
806 std::string adcLinearityCorr_filename = PathResolverFindCalibFile(
807 "ElectronPhotonFourMomentumCorrection/v25/linearity_ADC.root");
808 m_ADCLinearity_tool = std::make_shared<LinearityADC>(adcLinearityCorr_filename);
809 m_ADCLinearity_tool->msg().setLevel(this->msg().level());
810 m_rootTool->setADCTool(m_ADCLinearity_tool);
811 } else {
813 m_ESModel + " recommendations use ADC corrections for scale "
814 "derivation. Disabling the ADCLinearity flag will create "
815 "inconsistency!");
816 }
817
818 if (m_doLeakageCorrection == AUTO || m_doLeakageCorrection == 1) {
819 m_doLeakageCorrection = 1;
820 m_rootTool->setApplyLeakageCorrection(true);
821 }
822
823 // Calo distortion phi unif correction
824 if (m_useCaloDistPhiUnifCorrection == AUTO ||
825 m_useCaloDistPhiUnifCorrection == 1) {
826 m_useCaloDistPhiUnifCorrection = 1;
827 std::string phiUnifCorrfileName = PathResolverFindCalibFile(
828 "ElectronPhotonFourMomentumCorrection/v33/"
829 "egammaEnergyCorrectionData.root");
830 std::unique_ptr<TFile> fCorr(
831 TFile::Open(phiUnifCorrfileName.c_str(), "READ"));
832 m_caloDistPhiUnifCorr.reset(
833 dynamic_cast<TH2*>(fCorr->Get("CaloDistortionPhiUniformityCorrection/"
834 "es2023_R22_Run2_v0/h2DcorrPhiUnif")));
835 m_caloDistPhiUnifCorr->SetDirectory(nullptr);
836 } else {
838 m_ESModel + " recommendations use CaloDistPhiUnif for scale "
839 "derivation. Disabling the CaloDistPhiUnif flag will create "
840 "inconsistency!");
841 }
842 }
843
844 // No scale correction for release 21 ==> obsolete
845 /*if (m_ESModel == "es2017_R21_PRE"){
846 m_doScaleCorrection = 0;
847 }
848 */
849
850 ATH_MSG_INFO("ESModel: " << m_ESModel);
851 ATH_MSG_INFO("ResolutionType: " << m_ResolutionType);
852 ATH_MSG_INFO("decorrelation Model: " << m_decorrelation_model_name);
853 ATH_MSG_DEBUG("layer correction = " << m_useLayerCorrection);
854 ATH_MSG_DEBUG("PS correction = " << m_usePSCorrection);
855 ATH_MSG_DEBUG("S12 correction = " << m_useS12Correction);
856 ATH_MSG_DEBUG("Sacc correction = " << m_useSaccCorrection);
857 ATH_MSG_DEBUG("intermodule correction = " << m_useIntermoduleCorrection);
858 ATH_MSG_DEBUG("phi uniformity correction = " << m_usePhiUniformCorrection);
859 ATH_MSG_DEBUG("distorted calo phi uniformity correction = "
860 << m_useCaloDistPhiUnifCorrection);
861 ATH_MSG_DEBUG("gain correction = " << m_useGainCorrection);
862 ATH_MSG_DEBUG("ADC non-linearity correction = " << m_doADCLinearityCorrection);
863 ATH_MSG_DEBUG("leakage correction for photons = " << m_doLeakageCorrection);
864 ATH_MSG_DEBUG("smearing = " << m_doSmearing);
865 ATH_MSG_DEBUG("insitu scales = " << m_doScaleCorrection);
866 ATH_MSG_DEBUG("ep combination = " << m_use_ep_combination);
867 ATH_MSG_DEBUG("use MVA calibration = " << m_use_mva_calibration);
869 "use temperature correction 2015 = " << m_use_temp_correction201215);
870 ATH_MSG_DEBUG("use uA2MeV correction 2015 1/2 week = "
871 << m_use_uA2MeV_2015_first2weeks_correction);
872
873 setupSystematics();
874
875 applySystematicVariation(CP::SystematicSet())
876 .ignore(); // this set the flags for the internal tool without
877 // systematics
879 if (registry.registerSystematics(*this) != StatusCode::SUCCESS)
880 return StatusCode::FAILURE;
881
882 // For columnar, it is important only to set this accessor if it is
883 // needed, as it creates a hard data dependency on the column, which
884 // is not present in older PHYSLITE files, causing the tool to fail.
885 // An alternative would be to mark the column with `isOptional`, but
886 // since it is always required for the GNN calibration (and never
887 // otherwise), declaring it only for the GNN calibration seemed
888 // cleaner.
889 if (m_TESModel == egEnergyCorr::es2025_Run3_GNN_v0) {
890 resetAccessor (m_accessors->gnn_energy_Acc, *this, "TransformerEnergy");
891 }
892 if (m_onlyElectrons.value() && m_onlyPhotons.value()) {
893 ATH_MSG_ERROR("Cannot select both onlyElectrons and onlyPhotons");
894 return StatusCode::FAILURE;
895 }
896 if (m_onlyElectrons.value()) {
897 resetElectron (m_accessors->momAcc, *m_accessors);
898 if (m_TESModel == egEnergyCorr::es2011c) {
899 resetAccessor (m_accessors->electronTrackAcc, *this, "trackParticleLinks");
900 }
901 }
902 if (m_onlyPhotons.value()) {
903 resetPhoton (m_accessors->momAcc, *m_accessors);
904 resetAccessor (m_accessors->photonVertexAcc, *this, "vertexLinks");
905 }
906 if (m_decorateEmva)
907 resetAccessor (m_accessors->decEmva, *this, "E_mva_only");
908
909 ANA_CHECK (initializeColumns ());
910
911 return StatusCode::SUCCESS;
912}
913
914
915PATCore::ParticleType::Type EgammaCalibrationAndSmearingTool::xAOD2ptype(columnar::EgammaId particle) const
916{
917 const Accessors& acc = *m_accessors;
918
919 // this is departing from the logic below, as we are now requiring the
920 // user to specify at configuration time whether we run on electrons
921 // or photons. this is necessary to configure the columns we need
922 // correctly.
923 if (m_onlyElectrons.value())
924 {
926 }
927 if (m_onlyPhotons.value())
928 {
929 if (acc.photonVertexAcc(particle).size() > 0)
930 {
932 }
933 else
934 {
936 }
937 }
938
939 // this is the old logic and should not be visited in columnar mode
940 // (disabled by turning on onlyElectrons or onlyPhotons)
942 //no ForwardElectron ptype: consider them as Electron
944 else if (xAOD::EgammaHelpers::isPhoton(&particle.getXAODObject())) {
947 }
948 else {
949 ATH_MSG_ERROR("particle is not electron of photon");
950 throw std::runtime_error("particle is not electron or photon");
951 }
952 return ptype;
953}
954
955double EgammaCalibrationAndSmearingTool::getResolution(
956 const xAOD::Egamma& particle, bool withCT) const {
957 const auto ptype = xAOD2ptype(particle);
958 const auto cl_etaCalo =
959 xAOD::get_eta_calo(*particle.caloCluster(), particle.author());
960
961 return m_rootTool->resolution(particle.e(), particle.caloCluster()->eta(),
962 cl_etaCalo, ptype, withCT,
963 false); // TODO: always for full simulation
964}
965
966double EgammaCalibrationAndSmearingTool::resolution(
967 double energy, double cl_eta, double cl_etaCalo,
968 PATCore::ParticleType::Type ptype, bool withCT) const {
969 return m_rootTool->resolution(energy, cl_eta, cl_etaCalo, ptype, withCT,
970 false);
971}
972
973CP::CorrectionCode EgammaCalibrationAndSmearingTool::applyCorrection(
974 xAOD::Egamma& input) const {
975 // Retrieve the event information:
976 const xAOD::EventInfo* event_info = nullptr;
977 if (evtStore()->retrieve(event_info, "EventInfo").isFailure()) {
978 ATH_MSG_ERROR("No EventInfo object could be retrieved");
980 }
981 return applyCorrection(input, *event_info);
982}
983
984CP::CorrectionCode EgammaCalibrationAndSmearingTool::correctedCopy(
985 const xAOD::Electron& input, xAOD::Electron*& output) const {
986 // A sanity check:
987 if (output)
989 "Non-null pointer received. "
990 "There's a possible memory leak!");
991
992 output = new xAOD::Electron();
993 output->makePrivateStore(input);
994 return applyCorrection(*output);
995}
996
997CP::CorrectionCode EgammaCalibrationAndSmearingTool::correctedCopy(
998 const xAOD::Photon& input, xAOD::Photon*& output) const {
999 // A sanity check:
1000 if (output)
1002 "Non-null pointer received. "
1003 "There's a possible memory leak!");
1004
1005 output = new xAOD::Photon();
1006 output->makePrivateStore(input);
1007 return applyCorrection(*output);
1008}
1009
1010double EgammaCalibrationAndSmearingTool::getEnergy(
1011 const xAOD::Photon& input) const {
1012 xAOD::Photon* new_particle = nullptr;
1013 ANA_CHECK_THROW(correctedCopy(input, new_particle));
1014 const double e = new_particle->e();
1015 delete new_particle;
1016 return e;
1017}
1018
1019double EgammaCalibrationAndSmearingTool::getEnergy(
1020 const xAOD::Electron& input) const {
1021 xAOD::Electron* new_particle = nullptr;
1022 ANA_CHECK_THROW(correctedCopy(input, new_particle));
1023 const double e = new_particle->e();
1024 delete new_particle;
1025 return e;
1026}
1027
1028CP::CorrectionCode EgammaCalibrationAndSmearingTool::applyCorrection(
1029 columnar::MutableEgammaId input, columnar::EventInfoId event_info) const {
1030 const Accessors& acc = *m_accessors;
1031
1032 // only used in simulation (for the smearing)
1033 RandomNumber seed = m_set_seed_function(*this, input, event_info);
1034
1035 columnar::ClusterId inputCluster = acc.caloClusterAcc (input)[0].value();
1036
1037 if (m_layer_recalibration_tool && acc.authorAcc (input) !=
1039 ATH_MSG_DEBUG("applying energy recalibration before E0|E1|E2|E3 = "
1040 << acc.energyBEAcc (inputCluster, 0) << "|"
1041 << acc.energyBEAcc (inputCluster, 1) << "|"
1042 << acc.energyBEAcc (inputCluster, 2) << "|"
1043 << acc.energyBEAcc (inputCluster, 3));
1044 // for now just go back to the xAOD object to access the subtool
1045 const CP::CorrectionCode status_layer_recalibration = m_layer_recalibration_tool->applyCorrection(input.getXAODObject(), event_info.getXAODObject());
1046 if (status_layer_recalibration == CP::CorrectionCode::Error) { return CP::CorrectionCode::Error; }
1047 ATH_MSG_DEBUG("eta|phi = " << acc.etaAcc (input) << "|" << acc.phiAcc (input));
1048 if (status_layer_recalibration == CP::CorrectionCode::Ok) {
1049 ATH_MSG_DEBUG("decoration E0|E1|E2|E3 = "
1050 << acc.Es0Acc(inputCluster) << "|"
1051 << acc.Es1Acc(inputCluster) << "|"
1052 << acc.Es2Acc(inputCluster) << "|"
1053 << acc.Es3Acc(inputCluster) << "|");
1054 if (acc.Es2Acc(inputCluster) == 0 and acc.Es1Acc(inputCluster) == 0 and
1055 acc.Es3Acc(inputCluster) == 0 and acc.Es0Acc(inputCluster) == 0 and
1056 (std::abs(acc.etaAcc (input)) < 1.37 or (std::abs(acc.etaAcc (input)) > 1.55 and std::abs(acc.etaAcc (input)) < 2.47)))
1057 {
1058 ATH_MSG_WARNING("all layer energies are zero");
1059 }
1060 }
1061 }
1062
1063 double energy = acc.momAcc.e (input);
1064 // apply MVA calibration
1065 if (!m_MVACalibSvc.empty()) {
1067 if (acc.authorAcc (input) ==
1069 const xAOD::VertexContainer* pVtxCont = nullptr;
1070 if (evtStore()->retrieve(pVtxCont, m_pVtxKey).isFailure()) {
1071 ATH_MSG_ERROR("No primary vertex container " << m_pVtxKey << " could be retrieved");
1073 }
1074 unsigned int npv(0);
1075 for (const auto *vtx : *pVtxCont) {
1076 if (vtx->vertexType() == xAOD::VxType::PriVtx ||
1077 vtx->vertexType() == xAOD::VxType::PileUp) { ++npv; }
1078 }
1079 gei.nPV = npv;
1080 gei.acmu = acc.actIntPerXingAcc(event_info);
1081 ATH_MSG_DEBUG("Retrieved nPV = " << gei.nPV << " and mu = " << gei.acmu);
1082 }
1083 if (acc.authorAcc (input) !=
1085 if (m_MVACalibSvc->getEnergy(inputCluster.getXAODObject(), input.getXAODObject(), energy, gei)
1086 .isFailure()) {
1087 ATH_MSG_ERROR("Failure in MVACalib service");
1089 }
1090 ATH_MSG_DEBUG("energy after MVA calibration = " << std::format("{:.2f}", energy));
1091 }
1092 }
1093 else if (m_TESModel == egEnergyCorr::es2025_Run3_GNN_v0) {
1094 // for now assume the input has already the layer calibration applied
1095 // and just take the decorated gnn_energy
1096 if (!acc.gnn_energy_Acc.isAvailable(input)) {
1097 ATH_MSG_ERROR("GNN energy requested but decoration TransformerEnergy not found");
1099 }
1100 energy = acc.gnn_energy_Acc(input);
1101 ATH_MSG_DEBUG("energy after GNN calibration = " << std::format("{:.2f}", energy));
1102
1103 }
1104 if (m_decorateEmva)
1105 {
1106 acc.decEmva(input) = energy;
1107 }
1108
1109 // For the time being, it is just the MVA calib
1110 if (acc.authorAcc (input) ==
1112 setPt(input, energy);
1114 }
1115
1116 if (m_TESModel == egEnergyCorr::es2011c) {
1117 // Crack calibation correction for es2011c (calibration hits calibration)
1118 const auto ptype = xAOD2ptype(input);
1119 const double etaden =
1121 ? static_cast<const xAOD::Electron&>(input.getXAODObject()).trackParticle()->eta()
1122 : acc.clusterEtaAcc(inputCluster);
1123 energy *= m_rootTool->applyMCCalibration(acc.clusterEtaAcc(inputCluster),
1124 energy / cosh(etaden), ptype);
1125 ATH_MSG_DEBUG("energy after crack calibration es2011c = "
1126 << std::format("{:.2f}", energy));
1127 }
1128
1129 /*
1130 * Here we check for each event the kind of data DATA vs FullSim
1131 * The m_simulation flavour has already been configured
1132 */
1134 (acc.eventTypeAcc(event_info,xAOD::EventInfo::IS_SIMULATION))
1135 ? m_simulation
1137
1138 unsigned int runNumber_for_tool = 0;
1139
1140 // apply uniformity corrections to data
1141 if (dataType == PATCore::ParticleDataType::Data) {
1142 // Get run number
1143 runNumber_for_tool = acc.runNumberAcc(event_info);
1144 // Get etaCalo, phiCalo
1145 const auto cl_eta = acc.clusterEtaAcc(inputCluster);
1146 double etaCalo = 0, phiCalo = 0;
1147 if (m_ADCLinearity_tool || m_gain_tool_run2 || m_usePhiUniformCorrection) {
1148 etaCalo = acc.etaCaloAcc(inputCluster, acc.authorAcc(input), false);
1149 if (m_usePhiUniformCorrection) {
1150 phiCalo =
1151 acc.phiCaloAcc(inputCluster, acc.authorAcc(input), false);
1152 }
1153 }
1154
1155 // Intermodule
1156 if (m_useIntermoduleCorrection) {
1157 energy =
1158 intermodule_correction(energy, acc.clusterPhiAcc(inputCluster), cl_eta);
1159 ATH_MSG_DEBUG("energy after intermodule correction = "
1160 << std::format("{:.2f}", energy));
1161 }
1162
1163 // Calo distortion
1164 if ( (m_TESModel == egEnergyCorr::es2023_R22_Run2_v0 ||
1165 m_TESModel == egEnergyCorr::es2023_R22_Run2_v1 ||
1166 m_TESModel == egEnergyCorr::es2024_Run3_v0 ||
1167 m_TESModel == egEnergyCorr::es2025_Run3_GNN_v0) &&
1168 m_useCaloDistPhiUnifCorrection) {
1169 double etaC = acc.clusterEtaAcc(inputCluster);
1170 double phiC = acc.clusterPhiAcc(inputCluster);
1171 int ieta = m_caloDistPhiUnifCorr->GetXaxis()->FindBin(etaC);
1172 ieta = ieta == 0 ? 1
1173 : (ieta > m_caloDistPhiUnifCorr->GetNbinsX()
1174 ? m_caloDistPhiUnifCorr->GetNbinsX()
1175 : ieta);
1176 int iphi = m_caloDistPhiUnifCorr->GetYaxis()->FindBin(phiC);
1177 iphi = iphi == 0 ? 1
1178 : (iphi > m_caloDistPhiUnifCorr->GetNbinsY()
1179 ? m_caloDistPhiUnifCorr->GetNbinsY()
1180 : iphi);
1181 energy *= m_caloDistPhiUnifCorr->GetBinContent(ieta, iphi);
1183 "energy after phi uniformity correction (for calo distortion) = "
1184 << std::format("{:.2f}", energy));
1185 }
1186
1187 // Phi
1188 if (m_usePhiUniformCorrection) {
1189 energy *= correction_phi_unif(etaCalo, phiCalo);
1190 ATH_MSG_DEBUG("energy after uniformity correction = "
1191 << std::format("{:.2f}", energy));
1192 }
1193
1194 // ADC
1195 if (m_ADCLinearity_tool) {
1196 double et = energy / std::cosh(cl_eta);
1197 double corr =
1198 m_ADCLinearity_tool->getCorr(etaCalo, et, xAOD2ptype(input));
1199 energy *= corr;
1200 ATH_MSG_DEBUG("energy after ADC linearity correction = "
1201 << std::format("{:.2f}", energy));
1202 }
1203
1204 // Gain
1205 if (m_gain_tool) {
1206 const auto es2 = acc.Es2Acc.isAvailable(inputCluster)
1207 ? acc.Es2Acc(inputCluster)
1208 : acc.energyBEAcc (inputCluster,2);
1209 if (!(std::abs(cl_eta) < 1.52 and std::abs(cl_eta) > 1.37) and
1210 std::abs(cl_eta) < 2.4)
1211 energy = m_gain_tool->CorrectionGainTool(
1212 cl_eta, energy / GeV, es2 / GeV,
1213 xAOD2ptype(input)); // cl_eta ok, TODO: check corrected E2
1214 } else if (m_gain_tool_run2) {
1215 double et = energy / std::cosh(cl_eta);
1216 double corr = m_gain_tool_run2->getUncertainty(etaCalo, et,
1217 xAOD2ptype(input), true);
1218 energy /= (1 + corr);
1219 }
1220 ATH_MSG_DEBUG("energy after gain correction = " << std::format("{:.2f}", energy));
1221 } else {
1222 if (m_user_random_run_number == 0) {
1223 if (acc.randomrunnumber_getter.isAvailable(event_info)) {
1224 runNumber_for_tool = acc.randomrunnumber_getter(event_info);
1225 } else {
1227 "Pileup tool not run before using "
1228 "ElectronPhotonFourMomentumCorrection! Assuming it is 2016. If you "
1229 "want to force a specific period set the property randomRunNumber "
1230 "of the tool, e.g. in the job option: "
1231 "tool.randomRunNumber = 123456 or "
1232 "tool.randomRunNumber = "
1233 "EgammaCalibrationAndSmearingToolRunNumbersExample.run_2016");
1235 }
1236 } else {
1237 runNumber_for_tool = m_user_random_run_number;
1238 }
1239 }
1240
1241 const double eraw = ((acc.Es0Acc.isAvailable(inputCluster)
1242 ? acc.Es0Acc(inputCluster)
1243 : acc.energyBEAcc(inputCluster,0)) +
1244 (acc.Es1Acc.isAvailable(inputCluster)
1245 ? acc.Es1Acc(inputCluster)
1246 : acc.energyBEAcc(inputCluster,1)) +
1247 (acc.Es2Acc.isAvailable(inputCluster)
1248 ? acc.Es2Acc(inputCluster)
1249 : acc.energyBEAcc(inputCluster,2)) +
1250 (acc.Es3Acc.isAvailable(inputCluster)
1251 ? acc.Es3Acc(inputCluster)
1252 : acc.energyBEAcc(inputCluster,3)));
1253
1254
1255 if (dataType == PATCore::ParticleDataType::Fast)
1256 ATH_MSG_DEBUG("is fast");
1257 else if (dataType == PATCore::ParticleDataType::Full)
1258 ATH_MSG_DEBUG("is full");
1259 else if (dataType == PATCore::ParticleDataType::Data)
1260 ATH_MSG_DEBUG("is data");
1261
1262 // apply scale factors or systematics
1263 energy = m_rootTool->getCorrectedEnergy(
1264 runNumber_for_tool, dataType, xAOD2ptype(input),
1265 inputCluster(acc.clusterEtaAcc),
1266 inputCluster(acc.clusterEtaBEAcc,2),
1267 acc.etaCaloAcc(inputCluster, acc.authorAcc (input), false), energy,
1268 acc.Es2Acc.isAvailable(inputCluster)
1269 ? acc.Es2Acc(inputCluster)
1270 : inputCluster(acc.energyBEAcc,2),
1271 eraw, seed, oldtool_scale_flag_this_event(input, event_info),
1272 oldtool_resolution_flag_this_event(input, event_info), m_TResolutionType,
1273 m_varSF);
1274
1275 ATH_MSG_DEBUG("energy after scale/systematic correction = " << std::format("{:.2f}", energy));
1276
1277 // TODO: this check should be done before systematics variations
1278 setPt(input, energy);
1280}
1281
1282void EgammaCalibrationAndSmearingTool::setPt(columnar::MutableEgammaId input, double energy) const {
1283 const double new_energy2 = energy * energy;
1284 const auto ptype = xAOD2ptype(input);
1286 const double m2 = m * m;
1287 const double p2 = new_energy2 > m2 ? new_energy2 - m2 : 0.;
1288 m_accessors->ptOutDec (input) = sqrt(p2) / cosh(m_accessors->etaAcc (input));
1289 ATH_MSG_DEBUG("after setting pt, energy = " << m_accessors->momAcc.e (input));
1290}
1291
1292double EgammaCalibrationAndSmearingTool::getEnergy(
1293 xAOD::Egamma* p, const xAOD::EventInfo* event_info) {
1294 ANA_CHECK_THROW(applyCorrection(*p, *event_info));
1295 ATH_MSG_DEBUG("returning " << p->e());
1296 return p->e();
1297}
1298
1300EgammaCalibrationAndSmearingTool::oldtool_scale_flag_this_event(
1301 columnar::EgammaId p, columnar::EventInfoId event_info) const {
1302 const Accessors& acc = *m_accessors;
1303 if (!acc.eventTypeAcc (event_info, xAOD::EventInfo::IS_SIMULATION))
1304 return m_currentScaleVariation_data;
1305 if (m_currentScalePredicate(*this,p))
1306 return m_currentScaleVariation_MC;
1307 else
1309}
1310
1312EgammaCalibrationAndSmearingTool::oldtool_resolution_flag_this_event(
1313 columnar::EgammaId, columnar::EventInfoId event_info) const {
1314 const Accessors& acc = *m_accessors;
1315 return acc.eventTypeAcc (event_info, xAOD::EventInfo::IS_SIMULATION)
1316 ? m_currentResolutionVariation_MC
1317 : m_currentResolutionVariation_data;
1318}
1319
1320double EgammaCalibrationAndSmearingTool::getElectronMomentum(
1321 const xAOD::Electron* el, const xAOD::EventInfo* event_info) {
1324 ? m_simulation
1326
1327 const xAOD::TrackParticle* eTrack = el->trackParticle();
1328
1329 // track momentum and eta
1330 const float el_tracketa = eTrack->eta();
1331 const float el_trackmomentum = eTrack->pt() * cosh(el->eta());
1332
1333 return m_rootTool->getCorrectedMomentum(
1334 dataType, PATCore::ParticleType::Electron, el_trackmomentum, el_tracketa,
1335 oldtool_scale_flag_this_event(*el, *event_info), m_varSF);
1336}
1337
1338bool EgammaCalibrationAndSmearingTool::isAffectedBySystematic(
1339 const CP::SystematicVariation& systematic) const {
1340 CP::SystematicSet sys = affectingSystematics();
1341 return sys.find(systematic) != sys.end();
1342}
1343
1344CP::SystematicSet EgammaCalibrationAndSmearingTool::affectingSystematics()
1345 const {
1346 CP::SystematicSet affecting_systematics;
1347 for (const auto& it : m_syst_description) {
1348 affecting_systematics.insert(it.first);
1349 }
1350 for (const auto& it : m_syst_description_resolution) {
1351 affecting_systematics.insert(it.first);
1352 }
1353
1354 return affecting_systematics;
1355}
1356
1357void EgammaCalibrationAndSmearingTool::setupSystematics() {
1358 const EgammaPredicate always = [](const EgammaCalibrationAndSmearingTool&, columnar::EgammaId) { return true; };
1359
1360 // Try to simplify a bit for the ones that are fully correlate in eta,
1361 // whatever the model and that are not included in the macros including
1362 // - ADC non linearity
1363 // - L2Gain
1364 // - Leakage
1365 // - Conversion related
1366 // - TopoCluster threshold
1367 // - AF2
1368 // - PS_BARREL_B12
1369 // - S12EXTRALASTETABINRUN2
1370 // - ZEESTAT
1371 // - Run3 pre OFC + EXTRA
1372 if (m_decorrelation_model_scale == ScaleDecorrelation::FULL_ETA_CORRELATED ||
1373 m_decorrelation_model_scale == ScaleDecorrelation::FULL) {
1374 // Electron leakage, ADCLin, convReco only in final run2 recommendations
1375 if (m_TESModel == egEnergyCorr::es2023_R22_Run2_v0 ||
1376 m_TESModel == egEnergyCorr::es2023_R22_Run2_v1 ||
1377 m_TESModel == egEnergyCorr::es2024_Run3_v0) {
1378 // systematic related to ADC non linearity correction. Before 2022, there
1379 // was not correction, nor related systematic
1380 if (m_doADCLinearityCorrection) {
1381 m_syst_description[CP::SystematicVariation("EG_SCALE_ADCLIN", +1)] =
1382 SysInfo{always, egEnergyCorr::Scale::ADCLinUp};
1383 m_syst_description[CP::SystematicVariation("EG_SCALE_ADCLIN", -1)] =
1384 SysInfo{always, egEnergyCorr::Scale::ADCLinDown};
1385 }
1386 // Gain splitted uncertainty
1387 m_syst_description[CP::SystematicVariation("EG_SCALE_L2MEDIUMGAIN", +1)] =
1388 SysInfo{always, egEnergyCorr::Scale::L2MediumGainUp};
1389 m_syst_description[CP::SystematicVariation("EG_SCALE_L2MEDIUMGAIN", -1)] =
1391 m_syst_description[CP::SystematicVariation("EG_SCALE_L2LOWGAIN", +1)] =
1392 SysInfo{always, egEnergyCorr::Scale::L2LowGainUp};
1393 m_syst_description[CP::SystematicVariation("EG_SCALE_L2LOWGAIN", -1)] =
1394 SysInfo{always, egEnergyCorr::Scale::L2LowGainDown};
1395
1396 // Electron leakage
1397 m_syst_description[CP::SystematicVariation("EG_SCALE_LEAKAGEELEC", +1)] =
1398 SysInfo{always, egEnergyCorr::Scale::LeakageElecUp};
1399 m_syst_description[CP::SystematicVariation("EG_SCALE_LEAKAGEELEC", -1)] =
1400 SysInfo{always, egEnergyCorr::Scale::LeakageElecDown};
1401
1402 // Conversion related
1403 m_syst_description[CP::SystematicVariation("PH_SCALE_CONVRECO", +1)] =
1404 SysInfo{always, egEnergyCorr::Scale::ConvRecoUp};
1405 m_syst_description[CP::SystematicVariation("PH_SCALE_CONVRECO", -1)] =
1406 SysInfo{always, egEnergyCorr::Scale::ConvRecoDown};
1407 }
1408 // The equivalent of convReco (convefficiency and convfakerate) for other
1409 // models
1410 else {
1411 m_syst_description[CP::SystematicVariation("PH_SCALE_CONVEFFICIENCY",
1412 +1)] =
1414 m_syst_description[CP::SystematicVariation("PH_SCALE_CONVEFFICIENCY",
1415 -1)] =
1417 m_syst_description[CP::SystematicVariation("PH_SCALE_CONVFAKERATE", +1)] =
1418 SysInfo{always, egEnergyCorr::Scale::ConvFakeRateUp};
1419 m_syst_description[CP::SystematicVariation("PH_SCALE_CONVFAKERATE", -1)] =
1421 }
1422
1423 // additional systematics for R22 OFC and MC21 pre and bulk
1424 if (m_TESModel == egEnergyCorr::es2022_R22_PRE) {
1425 m_syst_description[CP::SystematicVariation("EG_SCALE_OFC", +1)] =
1426 SysInfo{always, egEnergyCorr::Scale::OFCUp};
1427 m_syst_description[CP::SystematicVariation("EG_SCALE_OFC", -1)] =
1428 SysInfo{always, egEnergyCorr::Scale::OFCDown};
1429
1430 m_syst_description[CP::SystematicVariation("EG_SCALE_EXTRARUN3PRE", +1)] =
1431 SysInfo{always, egEnergyCorr::Scale::EXTRARUN3PREUp};
1432 m_syst_description[CP::SystematicVariation("EG_SCALE_EXTRARUN3PRE", -1)] =
1434 }
1435
1436 if (m_TESModel == egEnergyCorr::es2017_R21_v0 or
1437 m_TESModel == egEnergyCorr::es2017_R21_v1 or
1438 m_TESModel == egEnergyCorr::es2017_R21_ofc0_v1 or
1439 m_TESModel == egEnergyCorr::es2024_Run3_ofc0_v0 or
1440 m_TESModel == egEnergyCorr::es2018_R21_v0 or
1441 m_TESModel == egEnergyCorr::es2018_R21_v1 or
1442 m_TESModel == egEnergyCorr::es2022_R22_PRE or
1443 m_TESModel == egEnergyCorr::es2023_R22_Run2_v0 or
1444 m_TESModel == egEnergyCorr::es2023_R22_Run2_v1 or
1445 m_TESModel == egEnergyCorr::es2024_Run3_v0) {
1446
1447 // topo clustr threshold systematics aded to release 21 recommendations
1448 m_syst_description[CP::SystematicVariation("EG_SCALE_TOPOCLUSTER_THRES",
1449 +1)] =
1451 m_syst_description[CP::SystematicVariation("EG_SCALE_TOPOCLUSTER_THRES",
1452 -1)] =
1454
1455 // AF3 for run3 models: es2022_R22_PRE and esmodel >= es2023_R22_Run2_v1
1456 // although we prevent AF for es2022_R22_PRE and es2024_Run3_ofc0_v0
1457 // AF3 are still technically added in the tool
1458 // but normally uncertainty will be 0
1459 if (m_TESModel >= egEnergyCorr::es2022_R22_PRE) {
1460 m_syst_description[CP::SystematicVariation("EG_SCALE_AF3", +1)] =
1461 SysInfo{always, egEnergyCorr::Scale::afUp};
1462 m_syst_description[CP::SystematicVariation("EG_SCALE_AF3", -1)] =
1463 SysInfo{always, egEnergyCorr::Scale::afDown};
1464 }
1465 else {
1466 // and extra AF2 systematics for release 21 recommendations - Moriond 2018
1467 // - pending proper AF2 to FullSim correction with release 21
1468 m_syst_description[CP::SystematicVariation("EG_SCALE_AF2", +1)] =
1469 SysInfo{always, egEnergyCorr::Scale::afUp};
1470 m_syst_description[CP::SystematicVariation("EG_SCALE_AF2", -1)] =
1471 SysInfo{always, egEnergyCorr::Scale::afDown};
1472 }
1473 }
1474
1475 // PS correlated barrel uncertainty
1476 if (m_TESModel == egEnergyCorr::es2017_summer_final or
1477 m_TESModel == egEnergyCorr::es2017_R21_v0 or
1478 m_TESModel == egEnergyCorr::es2017_R21_v1 or
1479 m_TESModel == egEnergyCorr::es2017_R21_ofc0_v1 or
1480 m_TESModel == egEnergyCorr::es2024_Run3_ofc0_v0 or
1481 m_TESModel == egEnergyCorr::es2018_R21_v0 or
1482 m_TESModel == egEnergyCorr::es2018_R21_v1 or
1483 m_TESModel == egEnergyCorr::es2022_R22_PRE) {
1484 m_syst_description[CP::SystematicVariation("EG_SCALE_PS_BARREL_B12",
1485 +1)] =
1486 SysInfo{always, egEnergyCorr::Scale::PSb12Up};
1487 m_syst_description[CP::SystematicVariation("EG_SCALE_PS_BARREL_B12",
1488 -1)] =
1489 SysInfo{always, egEnergyCorr::Scale::PSb12Down};
1490 }
1491
1492 // additional systematic for S12 last eta bin run2
1493 if (m_TESModel == egEnergyCorr::es2017 or
1494 m_TESModel == egEnergyCorr::es2017_summer or
1496 m_TESModel == egEnergyCorr::es2015_5TeV) {
1497 m_syst_description[CP::SystematicVariation(
1498 "EG_SCALE_S12EXTRALASTETABINRUN2", +1)] =
1500 m_syst_description[CP::SystematicVariation(
1501 "EG_SCALE_S12EXTRALASTETABINRUN2", -1)] =
1503 }
1504
1505 // Zee stat, if for FULL we do not ask for m_use_full_statistical_error
1506 if (m_decorrelation_model_scale ==
1507 ScaleDecorrelation::FULL_ETA_CORRELATED or
1508 !m_use_full_statistical_error) {
1509 // return 1 variation only, fully correlated in eta, equal to the correct
1510 // value but scaled by sqrt(number of bins) the scaling is done by the old
1511 // tool
1512 m_syst_description[CP::SystematicVariation("EG_SCALE_ZEESTAT", +1)] =
1513 SysInfo{always, egEnergyCorr::Scale::ZeeStatUp};
1514 m_syst_description[CP::SystematicVariation("EG_SCALE_ZEESTAT", -1)] =
1515 SysInfo{always, egEnergyCorr::Scale::ZeeStatDown};
1516 }
1517 }
1518 if (m_decorrelation_model_scale == ScaleDecorrelation::ONENP) {
1519 // TODO: independet implementation of ALL UP looping on all the variations
1520 m_syst_description[CP::SystematicVariation("EG_SCALE_ALL", +1)] =
1521 SysInfo{always, egEnergyCorr::Scale::AllUp};
1522 m_syst_description[CP::SystematicVariation("EG_SCALE_ALL", -1)] =
1523 SysInfo{always, egEnergyCorr::Scale::AllDown};
1524
1525 // to be consistent with other schemes, we add
1526 // extra AF systematics in addition to the 1NP
1527 if (m_TESModel == egEnergyCorr::es2017_R21_v0 ||
1528 m_TESModel == egEnergyCorr::es2017_R21_v1 ||
1529 m_TESModel == egEnergyCorr::es2017_R21_ofc0_v1 ||
1530 m_TESModel == egEnergyCorr::es2024_Run3_ofc0_v0 ||
1531 m_TESModel == egEnergyCorr::es2018_R21_v0 ||
1532 m_TESModel == egEnergyCorr::es2018_R21_v1) {
1533 m_syst_description[CP::SystematicVariation("EG_SCALE_AF2", +1)] =
1534 SysInfo{always, egEnergyCorr::Scale::afUp};
1535 m_syst_description[CP::SystematicVariation("EG_SCALE_AF2", -1)] =
1536 SysInfo{always, egEnergyCorr::Scale::afDown};
1537 }
1538 else if (m_TESModel >= egEnergyCorr::es2022_R22_PRE) {
1539 m_syst_description[CP::SystematicVariation("EG_SCALE_AF3", +1)] =
1540 SysInfo{always, egEnergyCorr::Scale::afUp};
1541 m_syst_description[CP::SystematicVariation("EG_SCALE_AF3", -1)] =
1542 SysInfo{always, egEnergyCorr::Scale::afDown};
1543 }
1544 }
1545 else if (m_decorrelation_model_scale ==
1546 ScaleDecorrelation::FULL_ETA_CORRELATED) {
1547// all the physical effects separately, considered as fully correlated in eta
1548 if (m_TESModel == egEnergyCorr::es2024_Run3_v0) {
1549#define SYSMACRO(name, fullcorrelated, decorrelation, flagup, flagdown) \
1550 m_syst_description[CP::SystematicVariation(#name, +1)] = \
1551 SysInfo{always, flagup}; \
1552 m_syst_description[CP::SystematicVariation(#name, -1)] = \
1553 SysInfo{always, flagdown};
1554#include "ElectronPhotonFourMomentumCorrection/systematics_es2024_Run3_v0.def"
1555#undef SYSMACRO
1556 }
1557 else {
1558// common systematics for all the esmodels
1559#define SYSMACRO(name, fullcorrelated, decorrelation, flagup, flagdown) \
1560 m_syst_description[CP::SystematicVariation(#name, +1)] = \
1561 SysInfo{always, flagup}; \
1562 m_syst_description[CP::SystematicVariation(#name, -1)] = \
1563 SysInfo{always, flagdown};
1564#include "ElectronPhotonFourMomentumCorrection/systematics_S12_2022.def"
1565#undef SYSMACRO
1566 }
1567
1568 if (m_TESModel != egEnergyCorr::es2023_R22_Run2_v0 &&
1569 m_TESModel != egEnergyCorr::es2023_R22_Run2_v1 &&
1570 m_TESModel != egEnergyCorr::es2024_Run3_v0) {
1571 m_syst_description[CP::SystematicVariation("EG_SCALE_LARCALIB", +1)] =
1572 SysInfo{always, egEnergyCorr::Scale::LArCalibUp};
1573 m_syst_description[CP::SystematicVariation("EG_SCALE_LARCALIB", -1)] =
1574 SysInfo{always, egEnergyCorr::Scale::LArCalibDown};
1575 m_syst_description[CP::SystematicVariation("EG_SCALE_L2GAIN", +1)] =
1576 SysInfo{always, egEnergyCorr::Scale::L2GainUp};
1577 m_syst_description[CP::SystematicVariation("EG_SCALE_L2GAIN", -1)] =
1578 SysInfo{always, egEnergyCorr::Scale::L2GainDown};
1579 }
1580
1581 // additional systematics for S12 run2
1582 if (m_TESModel == egEnergyCorr::es2015PRE_res_improved or
1583 m_TESModel == egEnergyCorr::es2015PRE or
1584 m_TESModel == egEnergyCorr::es2015cPRE or
1585 m_TESModel == egEnergyCorr::es2015c_summer or
1586 m_TESModel == egEnergyCorr::es2016PRE or
1587 m_TESModel == egEnergyCorr::es2017 or
1588 m_TESModel == egEnergyCorr::es2017_summer or
1590 m_TESModel == egEnergyCorr::es2015_5TeV) {
1591 m_syst_description[CP::SystematicVariation(
1592 "EG_SCALE_LARCALIB_EXTRA2015PRE", +1)] =
1594 m_syst_description[CP::SystematicVariation(
1595 "EG_SCALE_LARCALIB_EXTRA2015PRE", -1)] =
1597 }
1598
1599 // additional systematics for temperature run1->run2
1600 if (m_TESModel == egEnergyCorr::es2015PRE_res_improved or
1601 m_TESModel == egEnergyCorr::es2015PRE or
1602 m_TESModel == egEnergyCorr::es2015cPRE or
1603 m_TESModel == egEnergyCorr::es2015c_summer or
1604 m_TESModel == egEnergyCorr::es2016PRE) {
1605 m_syst_description[CP::SystematicVariation(
1606 "EG_SCALE_LARTEMPERATURE_EXTRA2015PRE", +1)] =
1608 m_syst_description[CP::SystematicVariation(
1609 "EG_SCALE_LARTEMPERATURE_EXTRA2015PRE", -1)] =
1611 }
1612
1613 // additional systematic for temperature 2015->2016
1614 if (m_TESModel == egEnergyCorr::es2016PRE) {
1615 m_syst_description[CP::SystematicVariation(
1616 "EG_SCALE_LARTEMPERATURE_EXTRA2016PRE", +1)] =
1618 m_syst_description[CP::SystematicVariation(
1619 "EG_SCALE_LARTEMPERATURE_EXTRA2016PRE", -1)] =
1621 }
1622
1623 // additional systematic for PP0 region
1624 if (m_TESModel == egEnergyCorr::es2017 or
1625 m_TESModel == egEnergyCorr::es2017_summer or
1627 m_TESModel == egEnergyCorr::es2017_summer_final or
1628 m_TESModel == egEnergyCorr::es2015_5TeV or
1629 m_TESModel == egEnergyCorr::es2017_R21_v0 or
1630 m_TESModel == egEnergyCorr::es2017_R21_v1 or
1631 m_TESModel == egEnergyCorr::es2017_R21_ofc0_v1 or
1632 m_TESModel == egEnergyCorr::es2024_Run3_ofc0_v0 or
1633 m_TESModel == egEnergyCorr::es2018_R21_v0 or
1634 m_TESModel == egEnergyCorr::es2018_R21_v1 or
1635 m_TESModel == egEnergyCorr::es2022_R22_PRE or
1636 m_TESModel == egEnergyCorr::es2023_R22_Run2_v0 or
1637 m_TESModel == egEnergyCorr::es2023_R22_Run2_v1 or
1638 m_TESModel == egEnergyCorr::es2024_Run3_v0) {
1639 m_syst_description[CP::SystematicVariation("EG_SCALE_MATPP0", +1)] =
1640 SysInfo{always, egEnergyCorr::Scale::MatPP0Up};
1641 m_syst_description[CP::SystematicVariation("EG_SCALE_MATPP0", -1)] =
1642 SysInfo{always, egEnergyCorr::Scale::MatPP0Down};
1643 }
1644
1645 // systematic related to wtots1
1646 if (m_TESModel == egEnergyCorr::es2017 or
1647 m_TESModel == egEnergyCorr::es2017_summer or
1649 m_TESModel == egEnergyCorr::es2017_summer_final or
1650 m_TESModel == egEnergyCorr::es2015_5TeV or
1651 m_TESModel == egEnergyCorr::es2017_R21_v0 or
1652 m_TESModel == egEnergyCorr::es2017_R21_v1 or
1653 m_TESModel == egEnergyCorr::es2017_R21_ofc0_v1 or
1654 m_TESModel == egEnergyCorr::es2024_Run3_ofc0_v0 or
1655 m_TESModel == egEnergyCorr::es2018_R21_v0 or
1656 m_TESModel == egEnergyCorr::es2018_R21_v1 or
1657 m_TESModel == egEnergyCorr::es2022_R22_PRE or
1658 m_TESModel == egEnergyCorr::es2023_R22_Run2_v0 or
1659 m_TESModel == egEnergyCorr::es2023_R22_Run2_v1 or
1660 m_TESModel == egEnergyCorr::es2024_Run3_v0) {
1661 m_syst_description[CP::SystematicVariation("EG_SCALE_WTOTS1", +1)] =
1662 SysInfo{always, egEnergyCorr::Scale::Wtots1Up};
1663 m_syst_description[CP::SystematicVariation("EG_SCALE_WTOTS1", -1)] =
1664 SysInfo{always, egEnergyCorr::Scale::Wtots1Down};
1665 }
1666
1667 // systematic for the scintillators
1668 if (m_TESModel == egEnergyCorr::es2015cPRE or
1669 m_TESModel == egEnergyCorr::es2015c_summer or
1670 m_TESModel == egEnergyCorr::es2016PRE or
1671 m_TESModel == egEnergyCorr::es2017 or
1672 m_TESModel == egEnergyCorr::es2017_summer or
1674 m_TESModel == egEnergyCorr::es2017_summer_final or
1675 m_TESModel == egEnergyCorr::es2015_5TeV or
1676 m_TESModel == egEnergyCorr::es2017_R21_v0 or
1677 m_TESModel == egEnergyCorr::es2017_R21_v1 or
1678 m_TESModel == egEnergyCorr::es2017_R21_ofc0_v1 or
1679 m_TESModel == egEnergyCorr::es2024_Run3_ofc0_v0 or
1680 m_TESModel == egEnergyCorr::es2018_R21_v0 or
1681 m_TESModel == egEnergyCorr::es2018_R21_v1 or
1682 m_TESModel == egEnergyCorr::es2022_R22_PRE or
1683 m_TESModel == egEnergyCorr::es2023_R22_Run2_v0 or
1684 m_TESModel == egEnergyCorr::es2023_R22_Run2_v1 or
1685 m_TESModel == egEnergyCorr::es2024_Run3_v0) {
1686 // scintillator systematics
1687 m_syst_description[CP::SystematicVariation("EG_SCALE_E4SCINTILLATOR",
1688 +1)] =
1690 m_syst_description[CP::SystematicVariation("EG_SCALE_E4SCINTILLATOR",
1691 -1)] =
1693 }
1694
1695 } else if (m_decorrelation_model_scale ==
1696 ScaleDecorrelation::ONENP_PLUS_UNCONR) {
1697// qsum of all variations correlated 8/13 TeV + uncorrelated (additional
1698// systematics for 2015PRE or 2016) all the physical effects separately,
1699// considered as fully correlated in eta
1700// TODO: fix for es2017
1701#define SYSMACRO(name, fullcorrelated, decorrelation, flagup, flagdown) \
1702 m_syst_description[CP::SystematicVariation(#name, +1)] = \
1703 SysInfo{always, flagup}; \
1704 m_syst_description[CP::SystematicVariation(#name, -1)] = \
1705 SysInfo{always, flagdown};
1706#include "ElectronPhotonFourMomentumCorrection/systematics_1NPCOR_PLUS_UNCOR.def"
1707#undef SYSMACRO
1708
1709 // additional systematic for S12 last eta bin run2 - not needed anymore for
1710 // last 20.7 model since it is part of bin per bin E1/E2 uncertainty in root
1711 // file
1712 if (m_TESModel == egEnergyCorr::es2017 or
1713 m_TESModel == egEnergyCorr::es2017_summer or
1715 m_syst_description[CP::SystematicVariation(
1716 "EG_SCALE_S12EXTRALASTETABINRUN2", +1)] =
1718 m_syst_description[CP::SystematicVariation(
1719 "EG_SCALE_S12EXTRALASTETABINRUN2", -1)] =
1721 }
1722
1723 } else if (m_decorrelation_model_scale == ScaleDecorrelation::FULL) {
1724 using pairvector = std::vector<std::pair<double, double>>;
1725 const pairvector decorrelation_bins_BE = {{0., 1.45}, {1.52, 2.5}};
1726 const std::vector<double> decorrelation_edges_TWELVE = {
1727 0., 0.2, 0.4, 0.6, 0.8, 1.0, 1.2, 1.4, 1.6, 1.8, 2.0, 2.2, 2.4};
1728 std::vector<double> decorrelation_edges_MODULE = {
1729 0., 0.2, 0.4, 0.6, 0.8, 1.0, 1.2, 1.37, 1.52, 1.8};
1730 const std::vector<double> decorrelation_edges_MATERIAL = {0.0, 1.1, 1.5,
1731 2.1, 2.5};
1732 std::vector<double> decorrelation_edges_S12_EXTRARUN3 = {
1733 0., 0.8, 1.5, 2.5};
1734
1735 std::vector<double> decorrelation_edges_S12;
1736 // for es2018_R21_v1 : 4 eta bins for muon E1/E2 uncertainty correlation
1737 if (m_TESModel == egEnergyCorr::es2018_R21_v1) {
1738 decorrelation_edges_S12.resize(5);
1739 decorrelation_edges_S12 = {0., 1.35, 1.5, 2.4, 2.5};
1740 } else if (m_TESModel == egEnergyCorr::es2023_R22_Run2_v0 or
1741 m_TESModel == egEnergyCorr::es2023_R22_Run2_v1 or
1742 m_TESModel == egEnergyCorr::es2024_Run3_v0) {
1743 decorrelation_edges_S12.resize(8);
1744 decorrelation_edges_S12 = {0., 0.6, 1.0, 1.35, 1.5, 1.8, 2.4, 2.5};
1745 //
1746 // PS scale from muons, so "crack" is a bit different
1747 decorrelation_edges_MODULE[7] = 1.4;
1748 decorrelation_edges_MODULE[8] = 1.5;
1749 }
1750 // for previous run 2 muon calibration with 20.7, 5 eta bins for E1/E2
1751 // uncertainty correlation
1752 else {
1753 decorrelation_edges_S12.resize(6);
1754 decorrelation_edges_S12 = {0., 0.6, 1.4, 1.5, 2.4, 2.5};
1755 }
1756
1757 if (m_TESModel != egEnergyCorr::es2017_summer_final and
1758 m_TESModel != egEnergyCorr::es2017_R21_v0 and
1759 m_TESModel != egEnergyCorr::es2017_R21_v1 and
1760 m_TESModel != egEnergyCorr::es2017_R21_ofc0_v1 and
1761 m_TESModel != egEnergyCorr::es2024_Run3_ofc0_v0 and
1762 m_TESModel != egEnergyCorr::es2018_R21_v0 and
1763 m_TESModel != egEnergyCorr::es2018_R21_v1 and
1764 m_TESModel != egEnergyCorr::es2022_R22_PRE and
1765 m_TESModel != egEnergyCorr::es2023_R22_Run2_v0 and
1766 m_TESModel != egEnergyCorr::es2023_R22_Run2_v1 and
1767 m_TESModel != egEnergyCorr::es2024_Run3_v0) {
1768#define SYSMACRO(name, fullcorrelated, decorrelation, flagup, flagdown) \
1769 if (bool(fullcorrelated)) { \
1770 m_syst_description[CP::SystematicVariation(#name, +1)] = \
1771 SysInfo{always, flagup}; \
1772 m_syst_description[CP::SystematicVariation(#name, -1)] = \
1773 SysInfo{always, flagdown}; \
1774 } else { \
1775 int i = 0; \
1776 for (const auto& p : AbsEtaCaloPredicatesFactory(decorrelation)) { \
1777 m_syst_description[CP::SystematicVariation( \
1778 #name "__ETABIN" + std::to_string(i), +1)] = SysInfo{p, flagup}; \
1779 m_syst_description[CP::SystematicVariation( \
1780 #name "__ETABIN" + std::to_string(i), -1)] = SysInfo{p, flagdown}; \
1781 i += 1; \
1782 } \
1783 }
1784#include "ElectronPhotonFourMomentumCorrection/systematics.def"
1785#undef SYSMACRO
1786 } else if (m_TESModel == egEnergyCorr::es2023_R22_Run2_v0 or
1787 m_TESModel == egEnergyCorr::es2023_R22_Run2_v1) {
1788#define SYSMACRO(name, fullcorrelated, decorrelation, flagup, flagdown) \
1789 if (bool(fullcorrelated)) { \
1790 m_syst_description[CP::SystematicVariation(#name, +1)] = \
1791 SysInfo{always, flagup}; \
1792 m_syst_description[CP::SystematicVariation(#name, -1)] = \
1793 SysInfo{always, flagdown}; \
1794 } else { \
1795 int i = 0; \
1796 for (const auto& p : AbsEtaCaloPredicatesFactory(decorrelation)) { \
1797 m_syst_description[CP::SystematicVariation( \
1798 #name "__ETABIN" + std::to_string(i), +1)] = SysInfo{p, flagup}; \
1799 m_syst_description[CP::SystematicVariation( \
1800 #name "__ETABIN" + std::to_string(i), -1)] = SysInfo{p, flagdown}; \
1801 i += 1; \
1802 } \
1803 }
1804#include "ElectronPhotonFourMomentumCorrection/systematics_S12_2022.def"
1805#undef SYSMACRO
1806 } else if (m_TESModel == egEnergyCorr::es2024_Run3_v0) {
1807#define SYSMACRO(name, fullcorrelated, decorrelation, flagup, flagdown) \
1808 if (bool(fullcorrelated)) { \
1809 m_syst_description[CP::SystematicVariation(#name, +1)] = \
1810 SysInfo{always, flagup}; \
1811 m_syst_description[CP::SystematicVariation(#name, -1)] = \
1812 SysInfo{always, flagdown}; \
1813 } else { \
1814 int i = 0; \
1815 for (const auto& p : AbsEtaCaloPredicatesFactory(decorrelation)) { \
1816 m_syst_description[CP::SystematicVariation( \
1817 #name "__ETABIN" + std::to_string(i), +1)] = SysInfo{p, flagup}; \
1818 m_syst_description[CP::SystematicVariation( \
1819 #name "__ETABIN" + std::to_string(i), -1)] = SysInfo{p, flagdown}; \
1820 i += 1; \
1821 } \
1822 }
1823#include "ElectronPhotonFourMomentumCorrection/systematics_es2024_Run3_v0.def"
1824#undef SYSMACRO
1825 } else {
1826#define SYSMACRO(name, fullcorrelated, decorrelation, flagup, flagdown) \
1827 if (bool(fullcorrelated)) { \
1828 m_syst_description[CP::SystematicVariation(#name, +1)] = \
1829 SysInfo{always, flagup}; \
1830 m_syst_description[CP::SystematicVariation(#name, -1)] = \
1831 SysInfo{always, flagdown}; \
1832 } else { \
1833 int i = 0; \
1834 for (const auto& p : AbsEtaCaloPredicatesFactory(decorrelation)) { \
1835 m_syst_description[CP::SystematicVariation( \
1836 #name "__ETABIN" + std::to_string(i), +1)] = SysInfo{p, flagup}; \
1837 m_syst_description[CP::SystematicVariation( \
1838 #name "__ETABIN" + std::to_string(i), -1)] = SysInfo{p, flagdown}; \
1839 i += 1; \
1840 } \
1841 }
1842#include "ElectronPhotonFourMomentumCorrection/systematics_S12.def"
1843#undef SYSMACRO
1844 } // else
1845
1846 if (m_use_full_statistical_error) {
1847 // statistical error, decorrelate in *all* the bins
1848 int i = 0;
1849 const TAxis& axis_statistical_error(m_rootTool->get_ZeeStat_eta_axis());
1850 for (int ibin = 1; ibin <= axis_statistical_error.GetNbins(); ++ibin) {
1851 auto p = EtaCaloPredicateFactory(
1852 axis_statistical_error.GetBinLowEdge(ibin),
1853 axis_statistical_error.GetBinLowEdge(ibin + 1));
1854 m_syst_description[CP::SystematicVariation(
1855 "EG_SCALE_ZEESTAT__ETABIN" + std::to_string(i), +1)] =
1857 m_syst_description[CP::SystematicVariation(
1858 "EG_SCALE_ZEESTAT__ETABIN" + std::to_string(i), -1)] =
1860 ++i;
1861 }
1862 }
1863
1864 // additional systematics for S12 run2
1865 if (m_TESModel == egEnergyCorr::es2015PRE_res_improved or
1866 m_TESModel == egEnergyCorr::es2015PRE or
1867 m_TESModel == egEnergyCorr::es2015cPRE or
1868 m_TESModel == egEnergyCorr::es2015c_summer or
1869 m_TESModel == egEnergyCorr::es2016PRE or
1870 m_TESModel == egEnergyCorr::es2017 or
1871 m_TESModel == egEnergyCorr::es2017_summer or
1873 m_TESModel == egEnergyCorr::es2015_5TeV) {
1874 m_syst_description[CP::SystematicVariation(
1875 "EG_SCALE_LARCALIB_EXTRA2015PRE__ETABIN0", +1)] =
1876 SysInfo{AbsEtaCaloPredicateFactory({0, 1.45}),
1878 m_syst_description[CP::SystematicVariation(
1879 "EG_SCALE_LARCALIB_EXTRA2015PRE__ETABIN0", -1)] =
1880 SysInfo{AbsEtaCaloPredicateFactory({0, 1.45}),
1882 m_syst_description[CP::SystematicVariation(
1883 "EG_SCALE_LARCALIB_EXTRA2015PRE__ETABIN1", +1)] =
1884 SysInfo{AbsEtaCaloPredicateFactory({1.45, 2.47}),
1886 m_syst_description[CP::SystematicVariation(
1887 "EG_SCALE_LARCALIB_EXTRA2015PRE__ETABIN1", -1)] =
1888 SysInfo{AbsEtaCaloPredicateFactory({1.45, 2.47}),
1890 m_syst_description[CP::SystematicVariation(
1891 "EG_SCALE_LARCALIB_EXTRA2015PRE__ETABIN2", +1)] =
1892 SysInfo{AbsEtaCaloPredicateFactory({2.47, 3.2}),
1894 m_syst_description[CP::SystematicVariation(
1895 "EG_SCALE_LARCALIB_EXTRA2015PRE__ETABIN2", -1)] =
1896 SysInfo{AbsEtaCaloPredicateFactory({2.47, 3.2}),
1898 }
1899
1900 // additional systematics for temperature run1->run2
1901 if (m_TESModel == egEnergyCorr::es2015PRE_res_improved or
1902 m_TESModel == egEnergyCorr::es2015PRE or
1903 m_TESModel == egEnergyCorr::es2015cPRE or
1904 m_TESModel == egEnergyCorr::es2015c_summer or
1905 m_TESModel == egEnergyCorr::es2016PRE) {
1906 m_syst_description[CP::SystematicVariation(
1907 "EG_SCALE_LARTEMPERATURE_EXTRA2015PRE__ETABIN0", +1)] =
1908 SysInfo{AbsEtaCaloPredicateFactory(decorrelation_bins_BE[0]),
1910 m_syst_description[CP::SystematicVariation(
1911 "EG_SCALE_LARTEMPERATURE_EXTRA2015PRE__ETABIN0", -1)] =
1912 SysInfo{AbsEtaCaloPredicateFactory(decorrelation_bins_BE[0]),
1914 m_syst_description[CP::SystematicVariation(
1915 "EG_SCALE_LARTEMPERATURE_EXTRA2015PRE__ETABIN1", +1)] =
1916 SysInfo{AbsEtaCaloPredicateFactory(decorrelation_bins_BE[1]),
1918 m_syst_description[CP::SystematicVariation(
1919 "EG_SCALE_LARTEMPERATURE_EXTRA2015PRE__ETABIN1", -1)] =
1920 SysInfo{AbsEtaCaloPredicateFactory(decorrelation_bins_BE[1]),
1922 }
1923
1924 // additional systematic for temperature 2015->2016
1925 if (m_TESModel == egEnergyCorr::es2016PRE) {
1926 m_syst_description[CP::SystematicVariation(
1927 "EG_SCALE_LARTEMPERATURE_EXTRA2016PRE__ETABIN0", +1)] =
1928 SysInfo{AbsEtaCaloPredicateFactory(decorrelation_bins_BE[0]),
1930 m_syst_description[CP::SystematicVariation(
1931 "EG_SCALE_LARTEMPERATURE_EXTRA2016PRE__ETABIN1", +1)] =
1932 SysInfo{AbsEtaCaloPredicateFactory(decorrelation_bins_BE[1]),
1934 m_syst_description[CP::SystematicVariation(
1935 "EG_SCALE_LARTEMPERATURE_EXTRA2016PRE__ETABIN0", -1)] =
1936 SysInfo{AbsEtaCaloPredicateFactory(decorrelation_bins_BE[0]),
1938 m_syst_description[CP::SystematicVariation(
1939 "EG_SCALE_LARTEMPERATURE_EXTRA2016PRE__ETABIN1", -1)] =
1940 SysInfo{AbsEtaCaloPredicateFactory(decorrelation_bins_BE[1]),
1942 }
1943
1944 // additional systematic for PP0 region
1945 if (m_TESModel == egEnergyCorr::es2017 or
1946 m_TESModel == egEnergyCorr::es2017_summer or
1948 m_TESModel == egEnergyCorr::es2017_summer_final or
1949 m_TESModel == egEnergyCorr::es2015_5TeV or
1950 m_TESModel == egEnergyCorr::es2017_R21_v0 or
1951 m_TESModel == egEnergyCorr::es2017_R21_v1 or
1952 m_TESModel == egEnergyCorr::es2017_R21_ofc0_v1 or
1953 m_TESModel == egEnergyCorr::es2024_Run3_ofc0_v0 or
1954 m_TESModel == egEnergyCorr::es2018_R21_v0 or
1955 m_TESModel == egEnergyCorr::es2018_R21_v1 or
1956 m_TESModel == egEnergyCorr::es2022_R22_PRE or
1957 m_TESModel == egEnergyCorr::es2023_R22_Run2_v0 or
1958 m_TESModel == egEnergyCorr::es2023_R22_Run2_v1 or
1959 m_TESModel == egEnergyCorr::es2024_Run3_v0) {
1960 m_syst_description[CP::SystematicVariation("EG_SCALE_MATPP0__ETABIN0",
1961 +1)] = SysInfo{
1962 AbsEtaCaloPredicateFactory(0, 1.5), egEnergyCorr::Scale::MatPP0Up};
1963 m_syst_description[CP::SystematicVariation("EG_SCALE_MATPP0__ETABIN1",
1964 +1)] = SysInfo{
1965 AbsEtaCaloPredicateFactory(1.5, 2.5), egEnergyCorr::Scale::MatPP0Up};
1966 m_syst_description[CP::SystematicVariation("EG_SCALE_MATPP0__ETABIN0",
1967 -1)] = SysInfo{
1968 AbsEtaCaloPredicateFactory(0, 1.5), egEnergyCorr::Scale::MatPP0Down};
1969 m_syst_description[CP::SystematicVariation("EG_SCALE_MATPP0__ETABIN1",
1970 -1)] =
1971 SysInfo{AbsEtaCaloPredicateFactory(1.5, 2.5),
1973 }
1974
1975 // systematic related to wtots1
1976 if (m_TESModel == egEnergyCorr::es2017 or
1977 m_TESModel == egEnergyCorr::es2017_summer or
1979 m_TESModel == egEnergyCorr::es2017_summer_final or
1980 m_TESModel == egEnergyCorr::es2015_5TeV or
1981 m_TESModel == egEnergyCorr::es2017_R21_v0 or
1982 m_TESModel == egEnergyCorr::es2017_R21_v1 or
1983 m_TESModel == egEnergyCorr::es2017_R21_ofc0_v1 or
1984 m_TESModel == egEnergyCorr::es2024_Run3_ofc0_v0 or
1985 m_TESModel == egEnergyCorr::es2018_R21_v0 or
1986 m_TESModel == egEnergyCorr::es2018_R21_v1 or
1987 m_TESModel == egEnergyCorr::es2022_R22_PRE) {
1988 m_syst_description[CP::SystematicVariation("EG_SCALE_WTOTS1", +1)] =
1989 SysInfo{always, egEnergyCorr::Scale::Wtots1Up};
1990 m_syst_description[CP::SystematicVariation("EG_SCALE_WTOTS1", -1)] =
1991 SysInfo{always, egEnergyCorr::Scale::Wtots1Down};
1992 }
1993
1994 // systematic related to wtots1, decorrelate eta bin [1.52,1.82] from rest
1995 if (m_TESModel == egEnergyCorr::es2023_R22_Run2_v0 or
1996 m_TESModel == egEnergyCorr::es2023_R22_Run2_v1 or
1997 m_TESModel == egEnergyCorr::es2024_Run3_v0) {
1998 m_syst_description[CP::SystematicVariation("EG_SCALE_WTOTS1__ETABIN0",
1999 +1)] =
2000 SysInfo{DoubleOrAbsEtaCaloPredicate(0, 1.52, 1.82, 2.47),
2002 m_syst_description[CP::SystematicVariation("EG_SCALE_WTOTS1__ETABIN0",
2003 -1)] =
2004 SysInfo{DoubleOrAbsEtaCaloPredicate(0, 1.52, 1.82, 2.47),
2006 m_syst_description[CP::SystematicVariation("EG_SCALE_WTOTS1__ETABIN1",
2007 +1)] =
2008 SysInfo{AbsEtaCaloPredicateFactory({1.52, 1.82}),
2010 m_syst_description[CP::SystematicVariation("EG_SCALE_WTOTS1__ETABIN1",
2011 -1)] =
2012 SysInfo{AbsEtaCaloPredicateFactory({1.52, 1.82}),
2014 }
2015
2016 // systematic for the scintillators
2017 if (m_TESModel == egEnergyCorr::es2015cPRE or
2018 m_TESModel == egEnergyCorr::es2015c_summer or
2019 m_TESModel == egEnergyCorr::es2016PRE or
2020 m_TESModel == egEnergyCorr::es2017 or
2021 m_TESModel == egEnergyCorr::es2017_summer or
2023 m_TESModel == egEnergyCorr::es2017_summer_final or
2024 m_TESModel == egEnergyCorr::es2015_5TeV or
2025 m_TESModel == egEnergyCorr::es2017_R21_v0 or
2026 m_TESModel == egEnergyCorr::es2017_R21_v1 or
2027 m_TESModel == egEnergyCorr::es2017_R21_ofc0_v1 or
2028 m_TESModel == egEnergyCorr::es2024_Run3_ofc0_v0 or
2029 m_TESModel == egEnergyCorr::es2018_R21_v0 or
2030 m_TESModel == egEnergyCorr::es2018_R21_v1 or
2031 m_TESModel == egEnergyCorr::es2022_R22_PRE or
2032 m_TESModel == egEnergyCorr::es2023_R22_Run2_v0 or
2033 m_TESModel == egEnergyCorr::es2023_R22_Run2_v1 or
2034 m_TESModel == egEnergyCorr::es2024_Run3_v0) {
2035 m_syst_description[CP::SystematicVariation(
2036 "EG_SCALE_E4SCINTILLATOR__ETABIN0", +1)] =
2037 SysInfo{AbsEtaCaloPredicateFactory(1.4, 1.46),
2039 m_syst_description[CP::SystematicVariation(
2040 "EG_SCALE_E4SCINTILLATOR__ETABIN1", +1)] =
2041 SysInfo{AbsEtaCaloPredicateFactory(1.46, 1.52),
2043 m_syst_description[CP::SystematicVariation(
2044 "EG_SCALE_E4SCINTILLATOR__ETABIN2", +1)] =
2045 SysInfo{AbsEtaCaloPredicateFactory(1.52, 1.6),
2047 m_syst_description[CP::SystematicVariation(
2048 "EG_SCALE_E4SCINTILLATOR__ETABIN0", -1)] =
2049 SysInfo{AbsEtaCaloPredicateFactory(1.4, 1.46),
2051 m_syst_description[CP::SystematicVariation(
2052 "EG_SCALE_E4SCINTILLATOR__ETABIN1", -1)] =
2053 SysInfo{AbsEtaCaloPredicateFactory(1.46, 1.52),
2055 m_syst_description[CP::SystematicVariation(
2056 "EG_SCALE_E4SCINTILLATOR__ETABIN2", -1)] =
2057 SysInfo{AbsEtaCaloPredicateFactory(1.52, 1.6),
2059 if (m_TESModel == egEnergyCorr::es2024_Run3_v0) { // extended E4 in Run 3
2060 m_syst_description[CP::SystematicVariation(
2061 "EG_SCALE_E4SCINTILLATOR__ETABIN3", +1)] =
2062 SysInfo{AbsEtaCaloPredicateFactory(1.6, 1.72),
2064 m_syst_description[CP::SystematicVariation(
2065 "EG_SCALE_E4SCINTILLATOR__ETABIN3", -1)] =
2066 SysInfo{AbsEtaCaloPredicateFactory(1.6, 1.72),
2068 }
2069 }
2070 } else {
2071 ATH_MSG_FATAL("scale decorrelation model invalid");
2072 }
2073
2074 // resolution systematics
2075 if (m_decorrelation_model_resolution == ResolutionDecorrelation::ONENP) {
2076 // ALL will not include AF2/AF3 systematic
2077 // individual AF NP is always provided
2078 // linghua.guo@cern.ch 2025-04-23
2079 m_syst_description_resolution[CP::SystematicVariation(
2080 "EG_RESOLUTION_ALL", +1)] = egEnergyCorr::Resolution::AllUp;
2081 m_syst_description_resolution[CP::SystematicVariation(
2082 "EG_RESOLUTION_ALL", -1)] = egEnergyCorr::Resolution::AllDown;
2083 } else if (m_decorrelation_model_resolution ==
2084 ResolutionDecorrelation::FULL) {
2085 m_syst_description_resolution[CP::SystematicVariation(
2086 "EG_RESOLUTION_ZSMEARING", +1)] = egEnergyCorr::Resolution::ZSmearingUp;
2087 m_syst_description_resolution[CP::SystematicVariation(
2088 "EG_RESOLUTION_ZSMEARING", -1)] =
2090 m_syst_description_resolution[CP::SystematicVariation(
2091 "EG_RESOLUTION_SAMPLINGTERM", +1)] =
2093 m_syst_description_resolution[CP::SystematicVariation(
2094 "EG_RESOLUTION_SAMPLINGTERM", -1)] =
2096 m_syst_description_resolution[CP::SystematicVariation(
2097 "EG_RESOLUTION_MATERIALID", +1)] =
2099 m_syst_description_resolution[CP::SystematicVariation(
2100 "EG_RESOLUTION_MATERIALID", -1)] =
2102 m_syst_description_resolution[CP::SystematicVariation(
2103 "EG_RESOLUTION_MATERIALCALO", +1)] =
2105 m_syst_description_resolution[CP::SystematicVariation(
2106 "EG_RESOLUTION_MATERIALCALO", -1)] =
2108 m_syst_description_resolution[CP::SystematicVariation(
2109 "EG_RESOLUTION_MATERIALGAP", +1)] =
2111 m_syst_description_resolution[CP::SystematicVariation(
2112 "EG_RESOLUTION_MATERIALGAP", -1)] =
2114 m_syst_description_resolution[CP::SystematicVariation(
2115 "EG_RESOLUTION_MATERIALCRYO", +1)] =
2117 m_syst_description_resolution[CP::SystematicVariation(
2118 "EG_RESOLUTION_MATERIALCRYO", -1)] =
2120 m_syst_description_resolution[CP::SystematicVariation(
2121 "EG_RESOLUTION_PILEUP", +1)] = egEnergyCorr::Resolution::PileUpUp;
2122 m_syst_description_resolution[CP::SystematicVariation(
2123 "EG_RESOLUTION_PILEUP", -1)] = egEnergyCorr::Resolution::PileUpDown;
2124 if (m_TESModel == egEnergyCorr::es2017 or
2125 m_TESModel == egEnergyCorr::es2017_summer or
2127 m_TESModel == egEnergyCorr::es2017_summer_final or
2128 m_TESModel == egEnergyCorr::es2015_5TeV or
2129 m_TESModel == egEnergyCorr::es2017_R21_v0 or
2130 m_TESModel == egEnergyCorr::es2017_R21_v1 or
2131 m_TESModel == egEnergyCorr::es2017_R21_ofc0_v1 or
2132 m_TESModel == egEnergyCorr::es2024_Run3_ofc0_v0 or
2133 m_TESModel == egEnergyCorr::es2018_R21_v0 or
2134 m_TESModel == egEnergyCorr::es2018_R21_v1 or
2135 m_TESModel == egEnergyCorr::es2022_R22_PRE or
2136 m_TESModel == egEnergyCorr::es2023_R22_Run2_v0 or
2137 m_TESModel == egEnergyCorr::es2023_R22_Run2_v1 or
2138 m_TESModel == egEnergyCorr::es2024_Run3_v0) {
2139 m_syst_description_resolution[CP::SystematicVariation(
2140 "EG_RESOLUTION_MATERIALIBL", +1)] =
2142 m_syst_description_resolution[CP::SystematicVariation(
2143 "EG_RESOLUTION_MATERIALIBL", -1)] =
2145 m_syst_description_resolution[CP::SystematicVariation(
2146 "EG_RESOLUTION_MATERIALPP0", +1)] =
2148 m_syst_description_resolution[CP::SystematicVariation(
2149 "EG_RESOLUTION_MATERIALPP0", -1)] =
2151
2152 if (m_TESModel == egEnergyCorr::es2022_R22_PRE) { // exta sys. for Run-3
2153 // pre-recommendations
2154 m_syst_description_resolution[CP::SystematicVariation(
2155 "EG_RESOLUTION_OFC", +1)] = egEnergyCorr::Resolution::OFCUp;
2156 m_syst_description_resolution[CP::SystematicVariation(
2157 "EG_RESOLUTION_OFC", -1)] = egEnergyCorr::Resolution::OFCDown;
2158 }
2159 }
2160 } else {
2161 ATH_MSG_FATAL("resolution decorrelation model invalid");
2162 }
2163
2164 // Always use individual AF2/AF3 systematics for resolution
2165 if (m_TESModel == egEnergyCorr::es2017_R21_v1 ||
2166 m_TESModel == egEnergyCorr::es2017_R21_ofc0_v1 ||
2167 m_TESModel == egEnergyCorr::es2018_R21_v0 ||
2168 m_TESModel == egEnergyCorr::es2018_R21_v1) {
2169 m_syst_description_resolution[CP::SystematicVariation(
2170 "EG_RESOLUTION_AF2", +1)] = egEnergyCorr::Resolution::afUp;
2171 m_syst_description_resolution[CP::SystematicVariation(
2172 "EG_RESOLUTION_AF2", -1)] = egEnergyCorr::Resolution::afDown;
2173 }
2174 else if (m_TESModel >= egEnergyCorr::es2022_R22_PRE){
2175 m_syst_description_resolution[CP::SystematicVariation(
2176 "EG_RESOLUTION_AF3", +1)] = egEnergyCorr::Resolution::afUp;
2177 m_syst_description_resolution[CP::SystematicVariation(
2178 "EG_RESOLUTION_AF3", -1)] = egEnergyCorr::Resolution::afDown;
2179 }
2180
2181 // ep combination systematics
2182 if (m_use_ep_combination) {
2183 m_syst_description[CP::SystematicVariation("EL_SCALE_MOMENTUM", +1)] =
2184 SysInfo{always, egEnergyCorr::Scale::MomentumUp};
2185 m_syst_description[CP::SystematicVariation("EL_SCALE_MOMENTUM", -1)] =
2186 SysInfo{always, egEnergyCorr::Scale::MomentumDown};
2187 }
2188}
2189
2190CP::SystematicSet EgammaCalibrationAndSmearingTool::recommendedSystematics()
2191 const {
2192 return affectingSystematics();
2193}
2194
2195StatusCode EgammaCalibrationAndSmearingTool::applySystematicVariation(
2196 const CP::SystematicSet& systConfig) {
2197
2198 // set the nominal one (no systematics)
2199 m_currentScaleVariation_MC = egEnergyCorr::Scale::None;
2200 m_currentScaleVariation_data = m_doScaleCorrection
2203 m_currentResolutionVariation_MC = m_doSmearing
2206 m_currentResolutionVariation_data = egEnergyCorr::Resolution::None;
2207 m_currentScalePredicate = [](const EgammaCalibrationAndSmearingTool&, columnar::EgammaId) { return true; };
2208
2209 if (systConfig.empty())
2210 return StatusCode::SUCCESS;
2211
2212 // the following code allows only ONE systematic variation at a time (1 for
2213 // scale, 1 for resolution)
2214
2215 bool first_scale = true;
2216 bool first_resolution = true;
2217 for (const auto& it : systConfig) {
2218 const auto found_scale = m_syst_description.find(it);
2219 if (found_scale != m_syst_description.end()) {
2220 if (not first_scale) {
2221 ATH_MSG_ERROR("multiple scale variations not supported");
2222 throw std::runtime_error("multiple scale variations not supported");
2223 }
2224 first_scale = false;
2225 m_currentScaleVariation_MC = found_scale->second.effect;
2226 m_currentScalePredicate = found_scale->second.predicate;
2227 }
2228
2229 const auto found_resolution = m_syst_description_resolution.find(it);
2230 if (found_resolution != m_syst_description_resolution.end()) {
2231 if (not first_resolution) {
2232 ATH_MSG_ERROR("multiple resolution variations not supported");
2233 throw std::runtime_error(
2234 "multiple resolution variations not supported");
2235 }
2236 first_resolution = false;
2237 m_currentResolutionVariation_MC = found_resolution->second;
2238 }
2239 }
2240
2241 return StatusCode::SUCCESS;
2242}
2243
2244double EgammaCalibrationAndSmearingTool::intermodule_correction(
2245 double Ecl, double phi, double eta) const {
2246
2247 // Intermodule Widening Correction: E_corr = E / (a' - b' * ((1 / (1 +
2248 // exp((phi_mod - 2 * pi / 32) * c))) * (1 / (1 + exp((phi_mod - 2 * pi / 32)
2249 // * (d)))))) (phi_min, phi_max) : [a' = a / a, b' = b / a, c, d]
2250
2251 double Ecl_corr = 0.;
2252 int DivInt = 0;
2253 double pi = M_PI;
2254
2255 if (m_TESModel == egEnergyCorr::es2017_summer_improved ||
2256 m_TESModel == egEnergyCorr::es2017_summer_final ||
2257 m_TESModel == egEnergyCorr::es2017_R21_v0 ||
2258 m_TESModel == egEnergyCorr::es2017_R21_v1 ||
2259 m_TESModel == egEnergyCorr::es2017_R21_ofc0_v1 ||
2260 m_TESModel == egEnergyCorr::es2024_Run3_ofc0_v0 ||
2261 m_TESModel == egEnergyCorr::es2018_R21_v0 ||
2262 m_TESModel == egEnergyCorr::es2018_R21_v1 ||
2263 m_TESModel == egEnergyCorr::es2022_R22_PRE ||
2264 m_TESModel == egEnergyCorr::es2023_R22_Run2_v0 ||
2265 m_TESModel == egEnergyCorr::es2023_R22_Run2_v1 ||
2266 m_TESModel == egEnergyCorr::es2024_Run3_v0 ||
2267 m_TESModel == egEnergyCorr::es2025_Run3_GNN_v0) {
2268
2269 double phi_mod = 0;
2270 if (phi < 0)
2271 phi_mod = fmod(phi, 2 * pi / 16.) + pi / 8.;
2272 else
2273 phi_mod = fmod(phi, 2 * pi / 16.);
2274
2275 // The correction concerns only the barrel
2276 if (std::abs(eta) <= 1.37) {
2277
2278 if (phi < (-7 * pi / 8) && phi > (-1 * pi))
2279 Ecl_corr =
2280 Ecl /
2281 (1 - 0.1086 *
2282 ((1 / (1 + exp((phi_mod - 2 * pi / 32.) * 175.2759))) *
2283 (1 / (1 + exp((phi_mod - 2 * pi / 32.) * (-189.3612))))));
2284 if (phi < (-6 * pi / 8) && phi > (-7 * pi / 8))
2285 Ecl_corr =
2286 Ecl /
2287 (1 - 0.0596 *
2288 ((1 / (1 + exp((phi_mod - 2 * pi / 32.) * 170.8305))) *
2289 (1 / (1 + exp((phi_mod - 2 * pi / 32.) * (-233.3782))))));
2290 if (phi < (-5 * pi / 8) && phi > (-6 * pi / 8))
2291 Ecl_corr =
2292 Ecl /
2293 (1 - 0.0596 *
2294 ((1 / (1 + exp((phi_mod - 2 * pi / 32.) * 147.1451))) *
2295 (1 / (1 + exp((phi_mod - 2 * pi / 32.) * (-139.3386))))));
2296 if (phi < (-4 * pi / 8) && phi > (-5 * pi / 8))
2297 Ecl_corr =
2298 Ecl /
2299 (1 - 0.0583 *
2300 ((1 / (1 + exp((phi_mod - 2 * pi / 32.) * 168.4644))) *
2301 (1 / (1 + exp((phi_mod - 2 * pi / 32.) * (-246.2897))))));
2302 if (phi < (-3 * pi / 8) && phi > (-4 * pi / 8))
2303 Ecl_corr =
2304 Ecl /
2305 (1 - 0.0530 *
2306 ((1 / (1 + exp((phi_mod - 2 * pi / 32.) * 177.6703))) *
2307 (1 / (1 + exp((phi_mod - 2 * pi / 32.) * (-198.3227))))));
2308 if (phi < (-2 * pi / 8) && phi > (-3 * pi / 8))
2309 Ecl_corr =
2310 Ecl /
2311 (1 - 0.0672 *
2312 ((1 / (1 + exp((phi_mod - 2 * pi / 32.) * 145.0693))) *
2313 (1 / (1 + exp((phi_mod - 2 * pi / 32.) * (-242.1771))))));
2314 if (phi < (-1 * pi / 8) && phi > (-2 * pi / 8))
2315 Ecl_corr =
2316 Ecl /
2317 (1 - 0.0871 *
2318 ((1 / (1 + exp((phi_mod - 2 * pi / 32.) * 132.3303))) *
2319 (1 / (1 + exp((phi_mod - 2 * pi / 32.) * (-166.1833))))));
2320 if (phi < (0 * pi / 8) && phi > (-1 * pi / 8))
2321 Ecl_corr =
2322 Ecl /
2323 (1 - 0.0948 *
2324 ((1 / (1 + exp((phi_mod - 2 * pi / 32.) * 127.6780))) *
2325 (1 / (1 + exp((phi_mod - 2 * pi / 32.) * (-150.0700))))));
2326 if (phi < (1 * pi / 8) && phi > (0 * pi / 8))
2327 Ecl_corr =
2328 Ecl /
2329 (1 - 0.1166 *
2330 ((1 / (1 + exp((phi_mod - 2 * pi / 32.) * 172.0679))) *
2331 (1 / (1 + exp((phi_mod - 2 * pi / 32.) * (-235.3293))))));
2332 if (phi < (2 * pi / 8) && phi > (1 * pi / 8))
2333 Ecl_corr =
2334 Ecl /
2335 (1 - 0.1172 *
2336 ((1 / (1 + exp((phi_mod - 2 * pi / 32.) * 190.3524))) *
2337 (1 / (1 + exp((phi_mod - 2 * pi / 32.) * (-198.9400))))));
2338 if (phi < (3 * pi / 8) && phi > (2 * pi / 8))
2339 Ecl_corr =
2340 Ecl /
2341 (1 - 0.1292 *
2342 ((1 / (1 + exp((phi_mod - 2 * pi / 32.) * 158.0540))) *
2343 (1 / (1 + exp((phi_mod - 2 * pi / 32.) * (-165.3893))))));
2344 if (phi < (4 * pi / 8) && phi > (3 * pi / 8))
2345 Ecl_corr =
2346 Ecl /
2347 (1 - 0.1557 *
2348 ((1 / (1 + exp((phi_mod - 2 * pi / 32.) * 162.2793))) *
2349 (1 / (1 + exp((phi_mod - 2 * pi / 32.) * (-133.5131))))));
2350 if (phi < (5 * pi / 8) && phi > (4 * pi / 8))
2351 Ecl_corr =
2352 Ecl /
2353 (1 - 0.1659 *
2354 ((1 / (1 + exp((phi_mod - 2 * pi / 32.) * 180.5270))) *
2355 (1 / (1 + exp((phi_mod - 2 * pi / 32.) * (-168.5074))))));
2356 if (phi < (6 * pi / 8) && phi > (5 * pi / 8))
2357 Ecl_corr =
2358 Ecl /
2359 (1 - 0.1123 *
2360 ((1 / (1 + exp((phi_mod - 2 * pi / 32.) * 128.2277))) *
2361 (1 / (1 + exp((phi_mod - 2 * pi / 32.) * (-154.4455))))));
2362 if (phi < (7 * pi / 8) && phi > (6 * pi / 8))
2363 Ecl_corr =
2364 Ecl /
2365 (1 - 0.1394 *
2366 ((1 / (1 + exp((phi_mod - 2 * pi / 32.) * 192.1216))) *
2367 (1 / (1 + exp((phi_mod - 2 * pi / 32.) * (-198.0727))))));
2368 if (phi < (8 * pi / 8) && phi > (7 * pi / 8))
2369 Ecl_corr =
2370 Ecl /
2371 (1 - 0.1001 *
2372 ((1 / (1 + exp((phi_mod - 2 * pi / 32.) * 199.1735))) *
2373 (1 / (1 + exp((phi_mod - 2 * pi / 32.) * (-176.4056))))));
2374 }
2375
2376 // No correction for the EC
2377 else {
2378 Ecl_corr = Ecl;
2379 }
2380
2381 }
2382
2383 else {
2384
2385 // Definitions of module folding into four quarters (top, left, bottom and
2386 // right)
2387
2388 DivInt = (int)(phi / ((2 * pi) / 16.));
2389 double phi_mod = phi - DivInt * (2 * pi / 16.);
2390
2391 // Centring on the intermodule --> phi_mod will now be in [0,0.4]
2392 if (phi_mod < 0)
2393 phi_mod += pi / 8.;
2394
2395 // The correction concerns only the barrel
2396 if (std::abs(eta) <= 1.4) {
2397
2398 // Top quarter
2399 if (phi < (3 * pi) / 4. && phi >= pi / 4.) {
2400 Ecl_corr =
2401 Ecl / (1 - 0.131 * ((1 / (1 + exp((phi_mod - 0.2) * 199.08))) *
2402 (1 / (1 + exp((phi_mod - 0.2) * (-130.36))))));
2403 }
2404
2405 // Right quarter
2406 if (phi < pi / 4. && phi >= -pi / 4.) {
2407 Ecl_corr =
2408 Ecl / (1 - 0.0879 * ((1 / (1 + exp((phi_mod - 0.2) * 221.01))) *
2409 (1 / (1 + exp((phi_mod - 0.2) * (-149.51))))));
2410 }
2411 // Bottom quarter
2412 if (phi < -pi / 4. && phi >= (-3 * pi) / 4.) {
2413 Ecl_corr =
2414 Ecl / (1 - 0.0605 * ((1 / (1 + exp((phi_mod - 0.2) * 281.37))) *
2415 (1 / (1 + exp((phi_mod - 0.2) * (-170.29))))));
2416 }
2417 // Left quarter
2418 if ((phi < (-3 * pi) / 4.) || (phi >= (3 * pi) / 4.)) {
2419 Ecl_corr =
2420 Ecl / (1 - 0.102 * ((1 / (1 + exp((phi_mod - 0.2) * 235.37))) *
2421 (1 / (1 + exp((phi_mod - 0.2) * (-219.04))))));
2422 }
2423 }
2424
2425 // No correction for the EC
2426 else {
2427 Ecl_corr = Ecl;
2428 }
2429 }
2430
2431 return Ecl_corr;
2432}
2433
2434double EgammaCalibrationAndSmearingTool::correction_phi_unif(double eta,
2435 double phi) const {
2436 constexpr double PI = M_PI;
2437 double Fcorr = 1.0;
2438
2439 if (m_use_mapping_correction) {
2440 // wrong mapping HV -> sectors in run1
2441 if (eta < -0.4 && eta > -0.6) {
2442 if (phi < (14 * PI / 32.) && phi > (13 * PI / 32.)) {
2443 Fcorr += 0.035;
2444 } else if (phi < (13 * PI / 32.) && phi > (12 * PI / 32.)) {
2445 Fcorr -= 0.035;
2446 }
2447 }
2448 }
2449
2450 if (m_TESModel == egEnergyCorr::es2017_summer_improved ||
2451 m_TESModel == egEnergyCorr::es2017_summer_final ||
2452 m_TESModel == egEnergyCorr::es2017_R21_v0 ||
2453 m_TESModel == egEnergyCorr::es2017_R21_v1 or
2454 m_TESModel == egEnergyCorr::es2017_R21_ofc0_v1 ||
2455 m_TESModel == egEnergyCorr::es2024_Run3_ofc0_v0 ||
2456 m_TESModel == egEnergyCorr::es2018_R21_v0 ||
2457 m_TESModel == egEnergyCorr::es2018_R21_v1 ||
2458 m_TESModel == egEnergyCorr::es2022_R22_PRE ||
2459 m_TESModel == egEnergyCorr::es2023_R22_Run2_v0 ||
2460 m_TESModel == egEnergyCorr::es2023_R22_Run2_v1 ||
2461 m_TESModel == egEnergyCorr::es2024_Run3_v0 ||
2462 m_TESModel == egEnergyCorr::es2025_Run3_GNN_v0) {
2463
2464 if (eta < 0.2 && eta > 0.) {
2465 if (phi < (-7 * 2 * PI / 32.) && phi > (-8 * 2 * PI / 32.)) {
2466 Fcorr = 1.016314;
2467 }
2468 }
2469
2470 else if (eta < 0.6 && eta > 0.4) {
2471 if (phi < 0 && phi > (-2 * PI / 32.)) {
2472 Fcorr = 1.041591;
2473 } else if (phi < (-4 * 2 * PI / 32.) && phi > (-5 * 2 * PI / 32.)) {
2474 Fcorr = 1.067346;
2475 }
2476 }
2477
2478 else if (eta < 0.8 && eta > 0.6) {
2479 if (phi < (7 * 2 * PI / 32.) && phi > (6 * 2 * PI / 32.)) {
2480 Fcorr = 1.027980;
2481 }
2482 }
2483
2484 else if (eta < 1.4 && eta > 1.2) {
2485 if (phi < (-9 * 2 * PI / 32.) && phi > (-10 * 2 * PI / 32.)) {
2486 Fcorr = 1.020299;
2487 } else if (phi < (-11 * 2 * PI / 32.) && phi > (-12 * 2 * PI / 32.)) {
2488 Fcorr = 1.051426;
2489 }
2490 }
2491
2492 else if (eta < 2.3 && eta > 2.1) {
2493 if (phi < (-12 * 2 * PI / 32.) && phi > (-13 * 2 * PI / 32.)) {
2494 Fcorr = 1.071695;
2495 }
2496 }
2497
2498 else if (eta < 0. && eta > -0.2) {
2499 if (phi < (-12 * 2 * PI / 32.) && phi > (-13 * 2 * PI / 32.)) {
2500 Fcorr = 1.008227;
2501 } else if (phi < (-8 * 2 * PI / 32.) && phi > (-9 * 2 * PI / 32.)) {
2502 Fcorr = 1.013929;
2503 }
2504 }
2505
2506 else if (eta < -0.2 && eta > -0.4) {
2507 if (phi < (-9 * 2 * PI / 32.) && phi > (-10 * 2 * PI / 32.)) {
2508 Fcorr = 1.015749;
2509 }
2510 }
2511
2512 else if (eta < -1.2 && eta > -1.4) {
2513 if (phi < (-6 * 2 * PI / 32.) && phi > (-7 * 2 * PI / 32.)) {
2514 Fcorr = 1.064954;
2515 }
2516 }
2517
2518 else if (eta < -1.6 && eta > -1.8) {
2519 if (phi < (9 * 2 * PI / 32.) && phi > (8 * 2 * PI / 32.)) {
2520 Fcorr = 1.027448;
2521 }
2522 }
2523
2524 else if (eta < -2.3 && eta > -2.5) {
2525 if (phi < (-8 * 2 * PI / 32.) && phi > (-9 * 2 * PI / 32.)) {
2526 Fcorr = 1.025882;
2527 } else if (phi < (5 * 2 * PI / 32.) && phi > (4 * 2 * PI / 32.)) {
2528 Fcorr = 1.036616;
2529 } else if (phi < (9 * 2 * PI / 32.) && phi > (8 * 2 * PI / 32.)) {
2530 Fcorr = 1.053838;
2531 } else if (phi < (10 * 2 * PI / 32.) && phi > (9 * 2 * PI / 32.)) {
2532 Fcorr = 1.026856;
2533 } else if (phi < (11 * 2 * PI / 32.) && phi > (10 * 2 * PI / 32.)) {
2534 Fcorr = 0.994382;
2535 }
2536 }
2537
2538 } // es2017_summer_improved end
2539
2540 else {
2541 if (eta < 0.6 && eta > 0.4) {
2542 if (phi < 0 && phi > (-2 * PI / 32.)) {
2543 Fcorr = 1.028;
2544 } else if (phi < (-4 * 2 * PI / 32.) && phi > (-5 * 2 * PI / 32.)) {
2545 Fcorr = 1.044;
2546 }
2547 }
2548
2549 else if (eta < 0.8 && eta > 0.6) {
2550 if (phi < (7 * 2 * PI / 32.) && phi > (6 * 2 * PI / 32.)) {
2551 Fcorr = 1.022;
2552 }
2553 }
2554
2555 else if (eta < 1.4 && eta > 1.2) {
2556 if (phi < (-11 * 2 * PI / 32.) && phi > (-12 * 2 * PI / 32.)) {
2557 Fcorr = 1.038;
2558 }
2559 }
2560
2561 else if (eta < 2.0 && eta > 1.9) {
2562 if (phi < (10 * 2 * PI / 32.) && phi > (9 * 2 * PI / 32.)) {
2563 Fcorr = 1.029;
2564 }
2565 }
2566
2567 else if (eta < -1.2 && eta > -1.4) {
2568 if (phi < (-4 * 2 * PI / 32.) && phi > (-5 * 2 * PI / 32.)) {
2569 Fcorr = 1.048;
2570 } else if (phi < (-6 * 2 * PI / 32.) && phi > (-7 * 2 * PI / 32.)) {
2571 Fcorr = 1.048;
2572 }
2573 }
2574
2575 else if (eta < -1.6 && eta > -1.8) {
2576 if (phi < (9 * 2 * PI / 32.) && phi > (8 * 2 * PI / 32.)) {
2577 Fcorr = 1.024;
2578 }
2579 }
2580
2581 else if (eta < -2.3 && eta > -2.5) {
2582 if (phi < (-8 * 2 * PI / 32.) && phi > (-9 * 2 * PI / 32.)) {
2583 Fcorr = 1.037;
2584 } else if (phi < (5 * 2 * PI / 32.) && phi > (4 * 2 * PI / 32.)) {
2585 Fcorr = 1.031;
2586 } else if (phi < (9 * 2 * PI / 32.) && phi > (8 * 2 * PI / 32.)) {
2587 Fcorr = 1.040;
2588 } else if (phi < (10 * 2 * PI / 32.) && phi > (9 * 2 * PI / 32.)) {
2589 Fcorr = 1.030;
2590 } else if (phi < (11 * 2 * PI / 32.) && phi > (10 * 2 * PI / 32.)) {
2591 Fcorr = 1.020;
2592 }
2593 }
2594 }
2595
2596 return Fcorr;
2597}
2598
2599void EgammaCalibrationAndSmearingTool ::
2600callSingleEvent (columnar::MutableEgammaRange egammas, columnar::EventInfoId event) const
2601{
2602 for (auto egamma : egammas) {
2603 if (applyCorrection (egamma, event) != CP::CorrectionCode::Ok)
2604 throw std::runtime_error ("EgammaCalibrationAndSmearingTool::callEvents: apply failed");
2605 }
2606}
2607
2608void EgammaCalibrationAndSmearingTool ::
2609callEvents (columnar::EventContextRange events) const
2610{
2611 const Accessors& acc = *m_accessors;
2612 for (auto event : events) {
2613 auto eventInfo = acc.m_eventHandle(event);
2614 callSingleEvent (acc.m_egammaHandle(event), eventInfo);
2615 }
2616}
2617
2618} // namespace ADD_CMODE(CP)
#define M_PI
Scalar eta() const
pseudorapidity method
Scalar phi() const
phi method
#define ATH_CHECK
Evaluate an expression and check for errors.
#define ATH_MSG_DEBUG(x,...)
#define ATH_MSG_ERROR(x,...)
#define ATH_MSG_WARNING(x,...)
#define ATH_MSG_INFO(x,...)
#define ATH_MSG_FATAL(x,...)
#define ANA_CHECK(EXP)
check whether the given expression was successful
#define ANA_CHECK_THROW(EXP)
check whether the given expression was successful, throwing an exception on failure
std::string PathResolverFindCalibFile(const std::string &logical_file_name)
@ None
#define pi
static const Attributes_t empty
Return value from object correction CP tools.
@ Error
Some error happened during the object correction.
@ Ok
The correction was done successfully.
This module implements the central registry for handling systematic uncertainties with CP tools.
static SystematicRegistry & getInstance()
Get the singleton instance of the registry for the curren thread.
Class to wrap a set of SystematicVariations.
bool empty() const
returns: whether the set is empty
void insert(const SystematicVariation &systematic)
description: insert a systematic into the set
an object that can create a AsgService
an object that can create a AsgTool
virtual double e() const override
The total energy of the particle.
Definition Egamma_v1.cxx:86
bool eventType(EventType type) const
Check for one particular bitmask value.
@ IS_SIMULATION
true: simulation, false: data
virtual double e() const override final
The total energy of the particle.
Definition Photon_v1.cxx:46
virtual double pt() const override final
The transverse momentum ( ) of the particle.
virtual double eta() const override final
The pseudorapidity ( ) of the particle.
void Scale(TH1 *h, double d=1)
::StatusCode StatusCode
StatusCode definition for legacy code.
constexpr double electronMassInMeV
the mass of the electron (in MeV)
constexpr double photonMassInMeV
various mass-less particles
constexpr ParticleHypothesis particle[PARTICLEHYPOTHESES]
the array of masses
ObjectRange< MutableEgammaDef > MutableEgammaRange
Definition EgammaDef.h:52
ObjectId< MutableEgammaDef > MutableEgammaId
Definition EgammaDef.h:53
void resetElectron(MomentumAccessors< CI, CM > &accessors, ColumnarTool< CM > &columnarTool)
ObjectRange< EventContextDef > EventContextRange
void resetPhoton(MomentumAccessors< CI, CM > &accessors, ColumnarTool< CM > &columnarTool)
ObjectId< ClusterDef > ClusterId
Definition ClusterDef.h:23
ObjectId< EgammaDef > EgammaId
Definition EgammaDef.h:47
ObjectId< EventInfoDef > EventInfoId
void resetAccessor(AccessorTemplate< CI, CT, CAM, CM > &accessor, ColumnarTool< CM > &columnBase, const std::string &name, ColumnAccessorOptions &&options={})
reset a column accessor to point to a new column
Definition dot.py:1
output
Definition merge.py:16
bool isConvertedPhoton(const xAOD::Egamma *eg, bool excludeTRT=false)
is the object a converted photon
bool isElectron(const xAOD::Egamma *eg)
is the object an electron (not Fwd)
bool isPhoton(const xAOD::Egamma *eg)
is the object a photon
const uint16_t AuthorFwdElectron
Electron reconstructed by the Forward cluster-based algorithm.
Definition EgammaDefs.h:30
@ PileUp
Pile-up vertex.
@ PriVtx
Primary vertex.
EventInfo_v1 EventInfo
Definition of the latest event info version.
TrackParticle_v1 TrackParticle
Reference the current persistent version:
VertexContainer_v1 VertexContainer
Definition of the current "Vertex container version".
Egamma_v1 Egamma
Definition of the current "egamma version".
Definition Egamma.h:17
float get_eta_calo(const xAOD::CaloCluster &cluster, int author, bool do_throw=false)
Photon_v1 Photon
Definition of the current "egamma version".
Definition Photon.h:17
setBGCode setTAP setLVL2ErrorBits bool
Electron_v1 Electron
Definition of the current "egamma version".
Definition Electron.h:17
A structure holding some global event information.
Extra patterns decribing particle interation process.
MsgStream & msg
Definition testRead.cxx:32
const float PI