ATLAS Offline Software
Loading...
Searching...
No Matches
MuonCalibIntScaleSmearTool.cxx
Go to the documentation of this file.
1/*
2 Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration
3*/
4
5// Framework include(s):
6#include "AsgDataHandles/ReadHandle.h"
7#include "PathResolver/PathResolver.h"
8
9// Local include(s):
10#include <cmath>
11#include "TRandom3.h"
12
13#include "MuonMomentumCorrections/MuonCalibIntScaleSmearTool.h"
14#include "MuonMomentumCorrections/EnumDef.h"
15#include <MuonMomentumCorrections/CalibInitializer.h>
16
17namespace CP
18{
19
22
23 StatusCode MuonCalibIntScaleSmearTool::initialize()
24 {
25 // Greet the user:
26 ATH_MSG_INFO("Initializing MuonCalibIntScaleSmearTool");
27
28 // Get the m_eventinfo container
29 ATH_CHECK(m_eventInfo.initialize());
30
31
32 // Load the constants
33 for(const auto & year: MCP::dataYearList)
34 {
35 for(const auto & param: m_paramList)
36 {
37 m_IDparams[year][param] = MCP::CalibInitializer::createScaleResCorrMap(year, MCP::TrackType::ID, m_release, param);
38 m_MEparams[year][param] = MCP::CalibInitializer::createScaleResCorrMap(year, MCP::TrackType::ME, m_release, param);
39 m_CBparams[year][param] = MCP::CalibInitializer::createScaleResCorrMap(year, MCP::TrackType::CB, m_release, param);
40 }
41
42 m_IDExpectedResparams[year] = MCP::CalibInitializer::createExpectedResMap(year, MCP::TrackType::ID, m_release);
43 m_MEExpectedResparams[year] = MCP::CalibInitializer::createExpectedResMap(year, MCP::TrackType::ME, m_release);
44 }
45
46
47 m_currentParameters = nullptr;
48 // Init the systematics
49 m_Parameters.initialize(affectingSystematics(), [this](const SystematicSet &systConfig, ParameterSetScaleSmear &param)
50 { return calcSystematicVariation(systConfig, param); });
51 if (!applySystematicVariation(SystematicSet()))
52 {
53 ATH_MSG_ERROR("Unable to run with no systematic");
54 return StatusCode::FAILURE;
55 }
56 SystematicRegistry &registry = SystematicRegistry::getInstance();
57 if (registry.registerSystematics(*this) != StatusCode::SUCCESS)
58 {
59 ATH_MSG_ERROR("Unkown systematic list");
60 return StatusCode::FAILURE;
61 }
62 // Return gracefully:
63 return StatusCode::SUCCESS;
64 }
65
66 CorrectionCode MuonCalibIntScaleSmearTool::applyCorrection(MCP::MuonObj& mu) const
67 {
68 // Do nothing, if it is data
69 if(mu.CB.isData) return CorrectionCode::Ok;
70
71 // Get the correction constants for each year
72 auto IDcorrConstants = getConstants(mu.ID, m_IDparams.at(mu.ID.year), m_currentParameters->m_IDparams);
73 auto MEcorrConstants = getConstants(mu.ME, m_MEparams.at(mu.ME.year), m_currentParameters->m_MEparams);
74 auto CBcorrConstants = getConstants(mu.CB, m_CBparams.at(mu.CB.year), m_currentParameters->m_CBparams);
75
76 double corrIDpT = getCorrectedPt(mu, mu.ID, IDcorrConstants);
77 double corrMEpT = getCorrectedPt(mu, mu.ME, MEcorrConstants);
78 double corrCBpT = getCorrectedPt(mu, mu.CB, CBcorrConstants);
79
80 //Smearing can be positive and negative since it depends on Gaussians centered in 0 with sigma 1.
81 //However, the resolution correction factor depends on 1/(1+smearing)
82 //and therefore if the smearing is less than -1 the pT turns out to have a negative value.
83 //Physically it means that the resolution has caused the pT to be reconstructed with a flip of the charge.
84 //Therefore in addition to correcting the pT, a flip of the charge is also done through the calib_charge variable.
85
86 // First compute the calibrated CB pT under the recombination scheme with uncalibrated ID & ME pT
87 double corrCBpTWithIDME = getCorrectedCBPtWithIDMSComb(mu, IDcorrConstants, MEcorrConstants);
88
89 // Write the pT into the object (if negative pT, multiply it by -1)
90 //< 0.1 because sometimes the pT is -0.
91 mu.ID.calib_pt = corrIDpT * ((corrIDpT < -0.1) ? -1 : 1);
92 mu.ME.calib_pt = corrMEpT * ((corrMEpT < -0.1) ? -1 : 1);
93 mu.CB.calib_pt = corrCBpT * ((corrCBpT < -0.1) ? -1 : 1);
94 //charge calibration (flip of the charge if needed)
95 mu.ID.calib_charge = mu.ID.uncalib_charge * ((corrIDpT < -0.1) ? -1 : 1);
96 mu.ME.calib_charge = mu.ME.uncalib_charge * ((corrMEpT < -0.1) ? -1 : 1);
97 mu.CB.calib_charge = mu.CB.uncalib_charge * ((corrCBpT < -0.1) ? -1 : 1);
98
99 if(!m_doDirectCBCalib) {
100 mu.CB.calib_pt = corrCBpTWithIDME * ((corrCBpTWithIDME < -0.1) ? -1 : 1);
101 mu.CB.calib_charge = mu.CB.uncalib_charge * ((corrCBpTWithIDME < -0.1) ? -1 : 1);
102 }
103
104 // Return gracefully:
105 return CorrectionCode::Ok;
106 }
107
108 std::map<MCP::ScaleSmearParam, double> MuonCalibIntScaleSmearTool::getConstants(const MCP::TrackCalibObj& trk, const std::map<MCP::ScaleSmearParam, ScaleSmearCorrConstMap>& constants, const std::map<MCP::ScaleSmearParam, double>& direction) const
109 {
110
111 std::map<MCP::ScaleSmearParam, double> calibConstants;
112
113 // Extra the constants from container into a simple map
114 for(const auto& param: m_paramList)
115 {
116 const auto & contantList = constants.at(param);
117
118 double val = contantList.at(MCP::ScaleResCorrection::Nominal)->getCalibConstant(trk);
119 if(direction.at(param) == MCP::SystVariation::Up)
120 {
121 val += contantList.at(MCP::ScaleResCorrection::SystErr__1up)->getCalibConstant(trk);
122 }
123 else if (direction.at(param) == MCP::SystVariation::Down)
124 {
125 val -= contantList.at(MCP::ScaleResCorrection::SystErr__1down)->getCalibConstant(trk);
126 }
127
128 // Remove unphysical smearing
129 if(param == MCP::ScaleSmearParam::r0 || param == MCP::ScaleSmearParam::r1 || param == MCP::ScaleSmearParam::r2)
130 {
131 val = std::max(0.0, val);
132 }
133
134 calibConstants[param] = val;
135 }
136
137 return calibConstants;
138 }
139
140
141 double MuonCalibIntScaleSmearTool::getCorrectedPt(const MCP::MuonObj& mu, const MCP::TrackCalibObj& trk, const std::map<MCP::ScaleSmearParam, double>& calibConstant) const
142 {
143 if(trk.calib_pt == 0) return 0;
144
145 // For debugging:: Todo add checking if in verbose mode
146 ATH_MSG_VERBOSE("MuonType: "<<MCP::toString(trk.type));
147 for(const auto& var: calibConstant) ATH_MSG_VERBOSE("var: "<<MCP::toString(var.first)<<" = "<<var.second);
148
149 double pT = trk.calib_pt;
150
151 // Calculate the denominator for the scale correction
152 double smearCorr = getSmearCorr(mu, trk, calibConstant.at(MCP::ScaleSmearParam::r0), calibConstant.at(MCP::ScaleSmearParam::r1), calibConstant.at(MCP::ScaleSmearParam::r2));
153 ATH_MSG_VERBOSE("smearCorr: "<<smearCorr);
154
155 // apply the full correction
156 double corrpT = (pT + pT * calibConstant.at(MCP::ScaleSmearParam::s1) + calibConstant.at(MCP::ScaleSmearParam::s0))/(1 + smearCorr);
157
158 return corrpT;
159 }
160
161 double MuonCalibIntScaleSmearTool::getSmearCorr(const MCP::MuonObj& mu, const MCP::TrackCalibObj& trk, double r0, double r1, double r2) const
162 {
163 if(trk.calib_pt == 0) return 1;
164
165 double pT = trk.calib_pt;
166 if(trk.type == MCP::TrackType::ID)
167 {
168 // ID applies a tan(theta) correction for r2 for high eta muons
169 double additional_weight = 1.;
170 if (std::abs(trk.eta) > 2) additional_weight = sinh(trk.eta);
171
172 return r1 * mu.rnd_g3 + r2 * mu.rnd_g4 * pT * additional_weight;
173 }
174 return r0 * mu.rnd_g0 / pT + r1 * mu.rnd_g1 + r2 * mu.rnd_g2 * pT;
175 }
176
177
178 double MuonCalibIntScaleSmearTool::getCorrectedCBPtWithIDMSComb(const MCP::MuonObj& mu, const std::map<MCP::ScaleSmearParam, double>& calibIDConstant, const std::map<MCP::ScaleSmearParam, double>& calibMEConstant) const
179 {
180 // calculate the relative of ID and ME for the corrections
181 double weightID = 0.5;
182 double weightME = 0.5;
183
184 double deltaCBME = mu.CB.calib_pt - mu.ME.calib_pt;
185 double deltaCBID = mu.CB.calib_pt - mu.ID.calib_pt;
186
187 if (mu.ME.calib_pt == 0)
188 {
189 weightID = 1.0;
190 weightME = 0.0;
191 }
192 else if (mu.ID.calib_pt == 0)
193 {
194 weightID = 0.0;
195 weightME = 1.0;
196 }
197 else if (mu.CB.calib_pt != 0)
198 {
199 if (std::abs(deltaCBME) > 0 || std::abs(deltaCBID) > 0)
200 {
201 double R = 1, Rplus = 1;
202 if (std::abs(deltaCBME) == std::abs(deltaCBID))
203 {
204 // do nothing
205 }
206 else if (std::abs(deltaCBME) != 0 &&
207 std::abs(deltaCBME) > std::abs(deltaCBID)) {
208 R = (-deltaCBID) / deltaCBME; /* R~wMS/wID */
209 Rplus = 1 + R;
210 if (Rplus != 0 && R > 0) {
211 weightID = 1 / Rplus;
212 weightME = R / Rplus;
213 }
214 }
215 else if (std::abs(deltaCBID) != 0 &&
216 std::abs(deltaCBME) < std::abs(deltaCBID)) {
217 R = (-deltaCBME) / (deltaCBID); /* R~wID/wMS */
218 Rplus = 1 + R;
219 if (Rplus != 0 && R > 0) {
220 weightID = R / Rplus;
221 weightME = 1 / Rplus;
222 }
223 }
224 }
225 }
226 // If no correction was found
227 if(weightID == 0.5 && weightME == 0.5)
228 {
229 if(mu.expectedPercentResME<std::numeric_limits<double>::epsilon() ||
230 mu.CB.calib_pt<std::numeric_limits<double>::epsilon() || mu.expectedPercentResID<std::numeric_limits<double>::epsilon() ){
231 ATH_MSG_VERBOSE("Potential FPE caught! Continue with averaging ID+MS");
232 ATH_MSG_VERBOSE("mu.expectedPercentResME = "<<mu.expectedPercentResME);
233 ATH_MSG_VERBOSE("mu.expectedPercentResID = "<<mu.expectedPercentResID);
234 ATH_MSG_VERBOSE("mu.CB.calib_pt = "<<mu.CB.calib_pt);
235 } else {
236 double wME = mu.ME.calib_pt / mu.CB.calib_pt / std::pow(mu.expectedPercentResME, 2);
237 double wID = mu.ID.calib_pt / mu.CB.calib_pt / std::pow(mu.expectedPercentResID, 2);
238 weightID = wID / (wME + wID);
239 weightME = wME / (wME + wID);
240 }
241 }
242
243
244 // Calculate the denominator for the scale correction
245 double smearIDCorr = getSmearCorr(mu, mu.ID, calibIDConstant.at(MCP::ScaleSmearParam::r0), calibIDConstant.at(MCP::ScaleSmearParam::r1), calibIDConstant.at(MCP::ScaleSmearParam::r2));
246 double smearMECorr = getSmearCorr(mu, mu.ME, calibMEConstant.at(MCP::ScaleSmearParam::r0), calibMEConstant.at(MCP::ScaleSmearParam::r1), calibMEConstant.at(MCP::ScaleSmearParam::r2));
247 double smearCorr = weightID * smearIDCorr + weightME * smearMECorr;
248 double scaleCB = 0;
249
250 // apply the full correction
251 if(weightID == 0) scaleCB = (calibMEConstant.at(MCP::ScaleSmearParam::s1) + calibMEConstant.at(MCP::ScaleSmearParam::s0) / mu.ME.calib_pt) * weightME;
252 else if(weightME == 0) scaleCB = calibIDConstant.at(MCP::ScaleSmearParam::s1) * weightID ;
253 else scaleCB = (calibIDConstant.at(MCP::ScaleSmearParam::s1) * weightID + (calibMEConstant.at(MCP::ScaleSmearParam::s1) + calibMEConstant.at(MCP::ScaleSmearParam::s0) / mu.ME.calib_pt) * weightME);
254
255 ATH_MSG_VERBOSE("mu.ID.calib_pt: "<<mu.ID.calib_pt);
256 ATH_MSG_VERBOSE("mu.ME.calib_pt: "<<mu.ME.calib_pt);
257 ATH_MSG_VERBOSE("mu.CB.calib_pt: "<<mu.CB.calib_pt);
258 ATH_MSG_VERBOSE("mu.expectedPercentResME: "<<mu.expectedPercentResME);
259 ATH_MSG_VERBOSE("mu.expectedPercentResID: "<<mu.expectedPercentResID);
260
261 ATH_MSG_VERBOSE("weightID: "<<weightID);
262 ATH_MSG_VERBOSE("weightME: "<<weightME);
263 ATH_MSG_VERBOSE("smearIDCorr: "<<smearIDCorr);
264 ATH_MSG_VERBOSE("smearMECorr: "<<smearMECorr);
265 ATH_MSG_VERBOSE("smearCorr: "<<smearCorr);
266 ATH_MSG_VERBOSE("scaleCB: "<<scaleCB);
267
268 double pT = mu.CB.calib_pt;
269 double corrpT = (pT + pT * scaleCB)/(1 + smearCorr);
270
271 return corrpT;
272 }
273
274
275 bool MuonCalibIntScaleSmearTool::isAffectedBySystematic(const SystematicVariation &systematic) const
276 {
277 SystematicSet sys = affectingSystematics();
278 return sys.find(systematic) != sys.end();
279 }
280
281 SystematicSet MuonCalibIntScaleSmearTool::affectingSystematics() const
282 {
283 SystematicSet result;
285 // Resolution systematics
287 if (m_doDirectCBCalib || m_sysScheme == "AllSys") {
288 // CB systematics
289 result.insert(SystematicVariation("MUON_CB", 1));
290 result.insert(SystematicVariation("MUON_CB", -1));
291 }
292 if (!m_doDirectCBCalib || m_sysScheme == "AllSys") {
293 // ID systematics
294 result.insert(SystematicVariation("MUON_ID", 1));
295 result.insert(SystematicVariation("MUON_ID", -1));
296
297 // MS systematics
298 result.insert(SystematicVariation("MUON_MS", 1));
299 result.insert(SystematicVariation("MUON_MS", -1));
300 }
301
305 if (m_sysScheme == "Corr_Scale") {
306 result.insert(SystematicVariation("MUON_SCALE", 1));
307 result.insert(SystematicVariation("MUON_SCALE", -1));
308 }
309 else if (m_sysScheme == "Decorr_Scale" || m_sysScheme == "AllSys") {
310 // Either doing direct calib of CB or asking for all the sys
311 if (m_doDirectCBCalib || m_sysScheme == "AllSys") {
312 result.insert(SystematicVariation("MUON_SCALE_CB", 1));
313 result.insert(SystematicVariation("MUON_SCALE_CB", -1));
314
315 result.insert(SystematicVariation("MUON_SCALE_CB_ELOSS", 1));
316 result.insert(SystematicVariation("MUON_SCALE_CB_ELOSS", -1));
317 }
318
319 // Either not doing direct calib of CB or asking for all the sys
320 if (!m_doDirectCBCalib || m_sysScheme == "AllSys") {
321 result.insert(SystematicVariation("MUON_SCALE_ID", 1));
322 result.insert(SystematicVariation("MUON_SCALE_ID", -1));
323
324 result.insert(SystematicVariation("MUON_SCALE_MS", 1));
325 result.insert(SystematicVariation("MUON_SCALE_MS", -1));
326
327 result.insert(SystematicVariation("MUON_SCALE_MS_ELOSS", 1));
328 result.insert(SystematicVariation("MUON_SCALE_MS_ELOSS", -1));
329 }
330 }
331
332
333
334
335 return result;
336 }
337
338 SystematicSet MuonCalibIntScaleSmearTool::recommendedSystematics() const { return affectingSystematics(); }
339
340 StatusCode MuonCalibIntScaleSmearTool::calcSystematicVariation(const SystematicSet &systConfig, ParameterSetScaleSmear &param) const
341 {
342 param.m_IDparams[MCP::ScaleSmearParam::s0] = MCP::SystVariation::Default;
343 param.m_IDparams[MCP::ScaleSmearParam::s1] = MCP::SystVariation::Default;
344 param.m_IDparams[MCP::ScaleSmearParam::r0] = MCP::SystVariation::Default;
345 param.m_IDparams[MCP::ScaleSmearParam::r1] = MCP::SystVariation::Default;
346 param.m_IDparams[MCP::ScaleSmearParam::r2] = MCP::SystVariation::Default;
347
348 param.m_MEparams[MCP::ScaleSmearParam::s0] = MCP::SystVariation::Default;
349 param.m_MEparams[MCP::ScaleSmearParam::s1] = MCP::SystVariation::Default;
350 param.m_MEparams[MCP::ScaleSmearParam::r0] = MCP::SystVariation::Default;
351 param.m_MEparams[MCP::ScaleSmearParam::r1] = MCP::SystVariation::Default;
352 param.m_MEparams[MCP::ScaleSmearParam::r2] = MCP::SystVariation::Default;
353
354 param.m_CBparams[MCP::ScaleSmearParam::s0] = MCP::SystVariation::Default;
355 param.m_CBparams[MCP::ScaleSmearParam::s1] = MCP::SystVariation::Default;
356 param.m_CBparams[MCP::ScaleSmearParam::r0] = MCP::SystVariation::Default;
357 param.m_CBparams[MCP::ScaleSmearParam::r1] = MCP::SystVariation::Default;
358 param.m_CBparams[MCP::ScaleSmearParam::r2] = MCP::SystVariation::Default;
359
360
361 // ID systematics
362 SystematicVariation syst = systConfig.getSystematicByBaseName("MUON_ID");
363
364 if (syst == SystematicVariation("MUON_ID", 1)) {
365 param.m_IDparams[MCP::ScaleSmearParam::r0] = MCP::SystVariation::Up;
366 param.m_IDparams[MCP::ScaleSmearParam::r1] = MCP::SystVariation::Up;
367 param.m_IDparams[MCP::ScaleSmearParam::r2] = MCP::SystVariation::Up;
368 } else if (syst == SystematicVariation("MUON_ID", -1)) {
369 param.m_IDparams[MCP::ScaleSmearParam::r0] = MCP::SystVariation::Down;
370 param.m_IDparams[MCP::ScaleSmearParam::r1] = MCP::SystVariation::Down;
371 param.m_IDparams[MCP::ScaleSmearParam::r2] = MCP::SystVariation::Down;
372 } else if (!syst.empty()) return StatusCode::FAILURE;
373
374 // MS systematics
375 syst = systConfig.getSystematicByBaseName("MUON_MS");
376
377 if (syst == SystematicVariation("MUON_MS", 1)) {
378 param.m_MEparams[MCP::ScaleSmearParam::r0] = MCP::SystVariation::Up;
379 param.m_MEparams[MCP::ScaleSmearParam::r1] = MCP::SystVariation::Up;
380 param.m_MEparams[MCP::ScaleSmearParam::r2] = MCP::SystVariation::Up;
381 } else if (syst == SystematicVariation("MUON_MS", -1)) {
382 param.m_MEparams[MCP::ScaleSmearParam::r0] = MCP::SystVariation::Down;
383 param.m_MEparams[MCP::ScaleSmearParam::r1] = MCP::SystVariation::Down;
384 param.m_MEparams[MCP::ScaleSmearParam::r2] = MCP::SystVariation::Down;
385 } else if (!syst.empty()) return StatusCode::FAILURE;
386
387 // CB systematics
388 syst = systConfig.getSystematicByBaseName("MUON_CB");
389
390 if (syst == SystematicVariation("MUON_CB", 1)) {
391 param.m_CBparams[MCP::ScaleSmearParam::r0] = MCP::SystVariation::Up;
392 param.m_CBparams[MCP::ScaleSmearParam::r1] = MCP::SystVariation::Up;
393 param.m_CBparams[MCP::ScaleSmearParam::r2] = MCP::SystVariation::Up;
394 } else if (syst == SystematicVariation("MUON_CB", -1)) {
395 param.m_CBparams[MCP::ScaleSmearParam::r0] = MCP::SystVariation::Down;
396 param.m_CBparams[MCP::ScaleSmearParam::r1] = MCP::SystVariation::Down;
397 param.m_CBparams[MCP::ScaleSmearParam::r2] = MCP::SystVariation::Down;
398 } else if (!syst.empty()) return StatusCode::FAILURE;
399
400 // Scale systematics
401 syst = systConfig.getSystematicByBaseName("MUON_SCALE");
402
403 if (syst == SystematicVariation("MUON_SCALE", 1)) {
404 param.m_IDparams[MCP::ScaleSmearParam::s0] = MCP::SystVariation::Down;
405 param.m_IDparams[MCP::ScaleSmearParam::s1] = MCP::SystVariation::Down;
406
407 param.m_MEparams[MCP::ScaleSmearParam::s0] = MCP::SystVariation::Down;
408 param.m_MEparams[MCP::ScaleSmearParam::s1] = MCP::SystVariation::Down;
409
410 param.m_CBparams[MCP::ScaleSmearParam::s0] = MCP::SystVariation::Down;
411 param.m_CBparams[MCP::ScaleSmearParam::s1] = MCP::SystVariation::Down;
412
413 } else if (syst == SystematicVariation("MUON_SCALE", -1)) {
414 param.m_IDparams[MCP::ScaleSmearParam::s0] = MCP::SystVariation::Up;
415 param.m_IDparams[MCP::ScaleSmearParam::s1] = MCP::SystVariation::Up;
416
417 param.m_MEparams[MCP::ScaleSmearParam::s0] = MCP::SystVariation::Up;
418 param.m_MEparams[MCP::ScaleSmearParam::s1] = MCP::SystVariation::Up;
419
420 param.m_CBparams[MCP::ScaleSmearParam::s0] = MCP::SystVariation::Up;
421 param.m_CBparams[MCP::ScaleSmearParam::s1] = MCP::SystVariation::Up;
422 } else if (!syst.empty()) return StatusCode::FAILURE;
423
424 // Split scale ID/MS/EGloss
425 syst = systConfig.getSystematicByBaseName("MUON_SCALE_ID");
426
427 if (syst == SystematicVariation("MUON_SCALE_ID", 1)) {
428 param.m_IDparams[MCP::ScaleSmearParam::s0] = MCP::SystVariation::Down;
429 param.m_IDparams[MCP::ScaleSmearParam::s1] = MCP::SystVariation::Down;
430 } else if (syst == SystematicVariation("MUON_SCALE_ID", -1)) {
431 param.m_IDparams[MCP::ScaleSmearParam::s0] = MCP::SystVariation::Up;
432 param.m_IDparams[MCP::ScaleSmearParam::s1] = MCP::SystVariation::Up;
433 } else if (!syst.empty()) return StatusCode::FAILURE;
434
435 syst = systConfig.getSystematicByBaseName("MUON_SCALE_MS");
436
437 if (syst == SystematicVariation("MUON_SCALE_MS", 1)) {
438 param.m_MEparams[MCP::ScaleSmearParam::s1] = MCP::SystVariation::Down;
439 } else if (syst == SystematicVariation("MUON_SCALE_MS", -1)) {
440 param.m_MEparams[MCP::ScaleSmearParam::s1] = MCP::SystVariation::Up;
441 } else if (!syst.empty()) return StatusCode::FAILURE;
442
443 syst = systConfig.getSystematicByBaseName("MUON_SCALE_MS_ELOSS");
444
445 if (syst == SystematicVariation("MUON_SCALE_MS_ELOSS", 1)) {
446 param.m_MEparams[MCP::ScaleSmearParam::s0] = MCP::SystVariation::Down;
447 } else if (syst == SystematicVariation("MUON_SCALE_MS_ELOSS", -1)) {
448 param.m_MEparams[MCP::ScaleSmearParam::s0] = MCP::SystVariation::Up;
449 } else if (!syst.empty()) return StatusCode::FAILURE;
450
451 syst = systConfig.getSystematicByBaseName("MUON_SCALE_CB");
452
453 if (syst == SystematicVariation("MUON_SCALE_CB", 1)) {
454 param.m_CBparams[MCP::ScaleSmearParam::s1] = MCP::SystVariation::Down;
455 } else if (syst == SystematicVariation("MUON_SCALE_CB", -1)) {
456 param.m_CBparams[MCP::ScaleSmearParam::s1] = MCP::SystVariation::Up;
457 } else if (!syst.empty()) return StatusCode::FAILURE;
458
459 syst = systConfig.getSystematicByBaseName("MUON_SCALE_CB_ELOSS");
460
461 if (syst == SystematicVariation("MUON_SCALE_CB_ELOSS", 1)) {
462 param.m_CBparams[MCP::ScaleSmearParam::s0] = MCP::SystVariation::Down;
463 } else if (syst == SystematicVariation("MUON_SCALE_CB_ELOSS", -1)) {
464 param.m_CBparams[MCP::ScaleSmearParam::s0] = MCP::SystVariation::Up;
465 } else if (!syst.empty()) return StatusCode::FAILURE;
466
467 return StatusCode::SUCCESS;
468 }
469
470 StatusCode MuonCalibIntScaleSmearTool::applySystematicVariation(const SystematicSet &systConfig)
471 {
472 return m_Parameters.get(systConfig, m_currentParameters);
473 }
474
475
476 double MuonCalibIntScaleSmearTool::getExpectedResolution(const int &DetType, double pT, double eta, double phi, MCP::DataYear year, bool addMCCorrectionSmearing) const
477 {
478
479 const auto & IDcorrConstants = m_IDparams.at(year);
480 const auto & MEcorrConstants = m_MEparams.at(year);
481 const auto & IDExpectedResConstants = m_IDExpectedResparams.at(year);
482 const auto & MEExpectedResConstants = m_MEExpectedResparams.at(year);
483
484 MCP::TrackType type = MCP::TrackType::ID;
485 if(DetType == MCP::DetectorType::MS) type = MCP::TrackType::ME;
486
487 auto trk = MCP::TrackCalibObj(type, pT * GeVtoMeV, eta, phi, year, false);
488
489
490 double expRes = 0.;
491 if (DetType == MCP::DetectorType::MS)
492 {
493 if (pT == 0) return 1e12;
494
495 double p0 = 0;
496 double p1 = 0;
497 double p2 = 0;
498
499 if(!addMCCorrectionSmearing)
500 {
501 p0 = MEExpectedResConstants.at(MCP::ExpectedResParam::r0)->getCalibConstant(trk);
502 p1 = MEExpectedResConstants.at(MCP::ExpectedResParam::r1)->getCalibConstant(trk);
503 p2 = MEExpectedResConstants.at(MCP::ExpectedResParam::r2)->getCalibConstant(trk);
504 }
505 else
506 {
507 double expectedP0 = MEExpectedResConstants.at(MCP::ExpectedResParam::r0)->getCalibConstant(trk);
508 double expectedP1 = MEExpectedResConstants.at(MCP::ExpectedResParam::r1)->getCalibConstant(trk);
509 double expectedP2 = MEExpectedResConstants.at(MCP::ExpectedResParam::r2)->getCalibConstant(trk);
510
511 double r0 = MEcorrConstants.at(MCP::ScaleSmearParam::r0).at(MCP::ScaleResCorrection::Nominal)->getCalibConstant(trk);
512 double r1 = MEcorrConstants.at(MCP::ScaleSmearParam::r1).at(MCP::ScaleResCorrection::Nominal)->getCalibConstant(trk);
513 double r2 = MEcorrConstants.at(MCP::ScaleSmearParam::r2).at(MCP::ScaleResCorrection::Nominal)->getCalibConstant(trk);
514
515 p0 = std::sqrt(std::pow(expectedP0, 2) + std::pow(r0, 2));
516 p1 = std::sqrt(std::pow(expectedP1, 2) + std::pow(r1, 2));
517 p2 = std::sqrt(std::pow(expectedP2, 2) + std::pow(r2, 2));
518 }
519
520 expRes = std::sqrt(std::pow(p0 / pT, 2) + std::pow(p1, 2) + std::pow(p2 * pT, 2));
521
522 }
523 else if (DetType == MCP::DetectorType::ID)
524 {
525 if (pT == 0) return 1e12;
526
527 double p1 = 0;
528 double p2 = 0;
529
530 if(!addMCCorrectionSmearing)
531 {
532 p1 = IDExpectedResConstants.at(MCP::ExpectedResParam::r1)->getCalibConstant(trk);
533 p2 = IDExpectedResConstants.at(MCP::ExpectedResParam::r2)->getCalibConstant(trk);
534 double p2Tan = IDExpectedResConstants.at(MCP::ExpectedResParam::r2tan2)->getCalibConstant(trk);
535 if(p2Tan) p2 = p2Tan;
536 }
537 else
538 {
539 double expectedP1 = IDExpectedResConstants.at(MCP::ExpectedResParam::r1)->getCalibConstant(trk);
540 double expectedP2 = IDExpectedResConstants.at(MCP::ExpectedResParam::r2)->getCalibConstant(trk);
541 double p2Tan = IDExpectedResConstants.at(MCP::ExpectedResParam::r2tan2)->getCalibConstant(trk);
542 if(p2Tan) expectedP2 = p2Tan;
543
544 double r1 = IDcorrConstants.at(MCP::ScaleSmearParam::r1).at(MCP::ScaleResCorrection::Nominal)->getCalibConstant(trk);
545 double r2 = IDcorrConstants.at(MCP::ScaleSmearParam::r2).at(MCP::ScaleResCorrection::Nominal)->getCalibConstant(trk);
546
547 p1 = std::sqrt(std::pow(expectedP1, 2) + std::pow(r1, 2));
548 p2 = std::sqrt(std::pow(expectedP2, 2) + std::pow(r2, 2));
549
550 if(p2Tan) p2 = p2 * std::sinh(eta) * std::sinh(eta);
551 }
552
553 expRes = std::sqrt(std::pow(p1, 2) + std::pow(p2 * pT, 2));
554 return expRes;
555 }
556 else
557 {
558 ATH_MSG_ERROR("wrong DetType in input " << DetType);
559 return 0.;
560 }
561
562 return expRes;
563 }
564
565
566} // namespace CP
eta
Scalar eta() const
pseudorapidity method
Definition AmgMatrixBasePlugin.h:83
ATH_CHECK
#define ATH_CHECK
Evaluate an expression and check for errors.
Definition AthCheckMacros.h:40
ATH_MSG_ERROR
#define ATH_MSG_ERROR(x)
Definition AthMsgStreamMacros.h:33
ATH_MSG_INFO
#define ATH_MSG_INFO(x)
Definition AthMsgStreamMacros.h:31
ATH_MSG_VERBOSE
#define ATH_MSG_VERBOSE(x)
Definition AthMsgStreamMacros.h:28
ReadHandle.h
Handle class for reading from StoreGate.
CP::CorrectionCode
Return value from object correction CP tools.
Definition CorrectionCode.h:31
CP::CorrectionCode::Ok
@ Ok
The correction was done successfully.
Definition CorrectionCode.h:38
CP::IMuonCalibIntTool::applyCorrection
virtual CorrectionCode applyCorrection(MCP::MuonObj &mu) const =0
Declare the interface that the class provides.
CP::MuonCalibIntScaleSmearTool::m_MEExpectedResparams
std::map< MCP::DataYear, std::map< MCP::ExpectedResParam, std::shared_ptr< MCP::CalibContainer > > > m_MEExpectedResparams
Definition MuonCalibIntScaleSmearTool.h:96
CP::MuonCalibIntScaleSmearTool::m_release
Gaudi::Property< std::string > m_release
Definition MuonCalibIntScaleSmearTool.h:64
CP::MuonCalibIntScaleSmearTool::initialize
virtual StatusCode initialize() override
Dummy implementation of the initialisation function.
Definition MuonCalibIntScaleSmearTool.cxx:23
CP::MuonCalibIntScaleSmearTool::affectingSystematics
virtual SystematicSet affectingSystematics() const override
the list of all systematics this tool can be affected by
Definition MuonCalibIntScaleSmearTool.cxx:281
CP::MuonCalibIntScaleSmearTool::m_paramList
static constexpr std::array< MCP::ScaleSmearParam, 5 > m_paramList
Definition MuonCalibIntScaleSmearTool.h:98
CP::MuonCalibIntScaleSmearTool::m_IDExpectedResparams
std::map< MCP::DataYear, std::map< MCP::ExpectedResParam, std::shared_ptr< MCP::CalibContainer > > > m_IDExpectedResparams
Definition MuonCalibIntScaleSmearTool.h:95
CP::MuonCalibIntScaleSmearTool::recommendedSystematics
virtual SystematicSet recommendedSystematics() const override
the list of all systematics this tool recommends to use
Definition MuonCalibIntScaleSmearTool.cxx:338
CP::MuonCalibIntScaleSmearTool::m_CBparams
std::map< MCP::DataYear, std::map< MCP::ScaleSmearParam, ScaleSmearCorrConstMap > > m_CBparams
Definition MuonCalibIntScaleSmearTool.h:91
CP::MuonCalibIntScaleSmearTool::isAffectedBySystematic
virtual bool isAffectedBySystematic(const SystematicVariation &systematic) const override
Declare the interface that this class provides.
Definition MuonCalibIntScaleSmearTool.cxx:275
CP::MuonCalibIntScaleSmearTool::m_sysScheme
Gaudi::Property< std::string > m_sysScheme
Definition MuonCalibIntScaleSmearTool.h:67
CP::MuonCalibIntScaleSmearTool::getCorrectedPt
double getCorrectedPt(const MCP::MuonObj &mu, const MCP::TrackCalibObj &trk, const std::map< MCP::ScaleSmearParam, double > &calibConstant) const
Definition MuonCalibIntScaleSmearTool.cxx:141
CP::MuonCalibIntScaleSmearTool::m_Parameters
SystematicsCache< ParameterSetScaleSmear > m_Parameters
Definition MuonCalibIntScaleSmearTool.h:85
CP::MuonCalibIntScaleSmearTool::getConstants
std::map< MCP::ScaleSmearParam, double > getConstants(const MCP::TrackCalibObj &trk, const std::map< MCP::ScaleSmearParam, ScaleSmearCorrConstMap > &constants, const std::map< MCP::ScaleSmearParam, double > &direction) const
Definition MuonCalibIntScaleSmearTool.cxx:108
CP::MuonCalibIntScaleSmearTool::getExpectedResolution
virtual ASG_TOOL_CLASS3(MuonCalibIntScaleSmearTool, CP::IMuonCalibIntScaleSmearTool, CP::ISystematicsTool, CP::IReentrantSystematicsTool) public double getExpectedResolution(const int &DetType, double pT, double eta, double phi, MCP::DataYear year, bool isData) const override
Declare the interface that the class provides.
Definition MuonCalibIntScaleSmearTool.cxx:476
CP::MuonCalibIntScaleSmearTool::m_currentParameters
const ParameterSetScaleSmear * m_currentParameters
Definition MuonCalibIntScaleSmearTool.h:86
CP::MuonCalibIntScaleSmearTool::m_eventInfo
SG::ReadHandleKey< xAOD::EventInfo > m_eventInfo
Definition MuonCalibIntScaleSmearTool.h:61
CP::MuonCalibIntScaleSmearTool::m_MEparams
std::map< MCP::DataYear, std::map< MCP::ScaleSmearParam, ScaleSmearCorrConstMap > > m_MEparams
Definition MuonCalibIntScaleSmearTool.h:90
CP::MuonCalibIntScaleSmearTool::MuonCalibIntScaleSmearTool
MuonCalibIntScaleSmearTool(const std::string &name)
Definition MuonCalibIntScaleSmearTool.cxx:20
CP::MuonCalibIntScaleSmearTool::applySystematicVariation
virtual StatusCode applySystematicVariation(const SystematicSet &systConfig) override
effects: configure this tool for the given list of systematic variations.
Definition MuonCalibIntScaleSmearTool.cxx:470
CP::MuonCalibIntScaleSmearTool::calcSystematicVariation
StatusCode calcSystematicVariation(const SystematicSet &systConfig, ParameterSetScaleSmear &param) const
Definition MuonCalibIntScaleSmearTool.cxx:340
CP::MuonCalibIntScaleSmearTool::getSmearCorr
double getSmearCorr(const MCP::MuonObj &mu, const MCP::TrackCalibObj &trk, double r0, double r1, double r2) const
Definition MuonCalibIntScaleSmearTool.cxx:161
CP::MuonCalibIntScaleSmearTool::getCorrectedCBPtWithIDMSComb
double getCorrectedCBPtWithIDMSComb(const MCP::MuonObj &mu, const std::map< MCP::ScaleSmearParam, double > &calibIDConstant, const std::map< MCP::ScaleSmearParam, double > &calibMEConstant) const
Definition MuonCalibIntScaleSmearTool.cxx:178
CP::MuonCalibIntScaleSmearTool::m_IDparams
std::map< MCP::DataYear, std::map< MCP::ScaleSmearParam, ScaleSmearCorrConstMap > > m_IDparams
Definition MuonCalibIntScaleSmearTool.h:89
CP::SystematicRegistry
This module implements the central registry for handling systematic uncertainties with CP tools.
Definition SystematicRegistry.h:25
CP::SystematicRegistry::getInstance
static SystematicRegistry & getInstance()
Get the singleton instance of the registry for the curren thread.
Definition SystematicRegistry.cxx:25
CP::SystematicRegistry::registerSystematics
StatusCode registerSystematics(const IReentrantSystematicsTool &tool)
effects: register all the systematics from the tool
Definition SystematicRegistry.cxx:45
CP::SystematicSet
Class to wrap a set of SystematicVariations.
Definition SystematicSet.h:31
CP::SystematicSet::getSystematicByBaseName
SystematicVariation getSystematicByBaseName(const std::string &basename) const
description: get the first systematic matching basename
Definition SystematicSet.cxx:170
CP::SystematicVariation::empty
bool empty() const
returns: whether this is an empty systematic, i.e.
Definition SystematicVariation.cxx:294
CorrectionCode::Ok
@ Ok
The correction was done successfully.
Definition CorrectionCode.h:38
asg::AsgTool::AsgTool
AsgTool(const std::string &name)
Constructor specifying the tool instance's name.
Definition AsgTool.cxx:58
CP
Select isolated Photons, Electrons and Muons.
Definition Control/xAODRootAccess/xAODRootAccess/TEvent.h:27
MCP::CalibInitializer::createExpectedResMap
std::map< ExpectedResParam, std::shared_ptr< CalibContainer > > createExpectedResMap(DataYear dataYear, TrackType type, const std::string &recommendationPath)
Definition CalibInitializer.cxx:74
MCP::CalibInitializer::createScaleResCorrMap
std::map< ScaleResCorrection, std::shared_ptr< CalibContainer > > createScaleResCorrMap(DataYear dataYear, TrackType type, const std::string &recommendationPath, ScaleSmearParam param)
Definition CalibInitializer.cxx:41
MCP::toString
std::string toString(TrackType trkType)
Definition EnumDef.h:39
MCP::DataYear
DataYear
Definition EnumDef.h:28
MCP::dataYearList
static constexpr std::array< MCP::DataYear, 7 > dataYearList
Definition EnumDef.h:34
MCP::TrackType
TrackType
Definition EnumDef.h:13
type
MCP::TrackCalibObj
Basic object to cache all relevant information from the track.
Definition MuonObj.h:74
MCP::TrackCalibObj::eta
const double eta
Value of the track-eta.
Definition MuonObj.h:157
MCP::TrackCalibObj::calib_pt
double calib_pt
Smeared track pt.
Definition MuonObj.h:155
MCP::TrackCalibObj::type
const TrackType type
Flag telling the code whether this is CB/ME/ID.
Definition MuonObj.h:149
result
Generated on for ATLAS Offline Software by doxygen 1.14.0