31 static const double muon_barrel_endcap_boundary = 1.05;
45 m_appliedSystematics(nullptr),
48 m_efficiencyMapReplicaArray(),
49 m_muonquality(
"Medium"),
50 m_calibration_version(
"240717_mc23ad"),
54 m_experimental(false),
57 m_replicaTriggerList(),
60 m_ReplicaRandomSeed(12345) {
73 declareProperty(
"AllowZeroSF",
m_allowZeroSF,
"If a trigger is not available will return 0 instead of throwing an error. More difficult to spot configuration issues. Use at own risk");
83 if (
year == 2015)
fileName =
"muontrigger_sf_2015_mc20a_v1.root";
84 else if (
year == 2016)
fileName =
"muontrigger_sf_2016_mc20a_v1.root";
85 else if (
year == 2017)
fileName =
"muontrigger_sf_2017_mc20d_v1.root";
86 else if (
year == 2018)
fileName =
"muontrigger_sf_2018_mc20e_v2.root";
87 else if (
year == 2022)
fileName =
"muontrigger_sf_2022_mc23a_v02.root";
88 else if (
year == 2023)
fileName =
"muontrigger_sf_2023_mc23d_v02.root";
91 return StatusCode::SUCCESS;
106 ATH_MSG_INFO(
"Note: setting up with user specified input file location " <<
m_custom_dir <<
" - this is not encouraged!");
112 if (
file ==
nullptr || !
file->IsOpen()) {
114 return StatusCode::FAILURE;
119 static const std::vector<std::string>
type {
"data",
"mc" };
120 static const std::vector<std::string> region {
"barrel",
"endcap" };
121 static const std::vector<std::string> systematic {
"nominal",
"stat_up",
"stat_down",
"syst_up",
"syst_down" };
126 if (qualityDirectory ==
nullptr) {
127 ATH_MSG_FATAL(
"MuonTriggerScaleFactors::initialize cannot find directory with selected quality");
128 return StatusCode::FAILURE;
131 TIter nextPeriod(qualityDirectory->GetListOfKeys());
132 while ((periodKey = (TKey*) nextPeriod())) {
133 if (not periodKey->IsFolder())
continue;
134 TDirectory* periodDirectory = qualityDirectory->GetDirectory(periodKey->GetName());
135 std::string
periodName = std::string(periodKey->GetName());
139 TIter nextTrigger(periodDirectory->GetListOfKeys());
140 while ((triggerKey = (TKey*) nextTrigger())) {
141 if (not triggerKey->IsFolder())
continue;
142 TDirectory* triggerDirectory = periodDirectory->GetDirectory(triggerKey->GetName());
143 std::string triggerName = std::string(triggerKey->GetName());
144 if(!std::set<std::string>{
"HLT_mu26_ivarmedium",
"HLT_mu50",
"HLT_mu26_ivarmedium_OR_HLT_mu50"}.count(triggerName) &&
m_binning ==
"coarse"){
145 ATH_MSG_DEBUG(
"Coarse binning not supported for di-muon trigger legs at the moment");
148 for (
const auto& iregion : region) {
149 bool isBarrel = iregion.find(
"barrel") != std::string::npos;
150 for (
const auto& itype :
type) {
151 bool isData = itype.find(
"data") != std::string::npos;
152 std::string
histname = (
"_MuonTrigEff_" +
periodName +
"_" + triggerName +
"_" + quality +
"_" +
"_EtaPhi_" +
m_binning +
"_" + iregion +
"_" + itype);
153 for (
const auto& isys : systematic) {
154 if (itype.find(
"data") != std::string::npos && isys.find(
"syst") != std::string::npos)
continue;
155 std::string
path =
"eff_etaphi_" +
m_binning +
"_" + iregion +
"_" + itype +
"_" + isys;
156 TH2*
hist =
dynamic_cast<TH2*
>(triggerDirectory->Get(
path.c_str()));
159 ATH_MSG_FATAL(
"MuonTriggerScaleFactors::initialize " <<
path <<
" not found under trigger " << triggerName <<
" and period " <<
periodName <<
" for year: " <<
year);
162 hist->SetDirectory(0);
167 return StatusCode::FAILURE;
169 m_efficiencyMap.insert(std::pair<EffiHistoIdent, TH1_Ptr>(HistoId, std::shared_ptr < TH1 > (
hist)));
179 TH1_Ptr tmp_h2 =
TH1_Ptr(
dynamic_cast<TH2F*
>(Nominal_H->Clone(Form(
"tmp_h2_%s", Nominal_H->GetName()))));
180 const int xbins = tmp_h2->GetNbinsX(),
ybins = tmp_h2->GetNbinsY();
181 for (
int x_i = 0; x_i <=
xbins; ++x_i) {
182 for (
int y_i = 0; y_i <=
ybins; ++y_i) {
183 double statErr = std::abs(tmp_h2->GetBinContent(x_i, y_i) - StatUp_H->GetBinContent(x_i, y_i));
184 tmp_h2->SetBinError(x_i, y_i, statErr);
196 return StatusCode::SUCCESS;
213 return StatusCode::FAILURE;
218 return StatusCode::FAILURE;
225 constexpr
auto years_to_run = std::to_array<int>({2015, 2016, 2017, 2018, 2022, 2023});
226 for (
const int &
year: years_to_run) {
229 return StatusCode::SUCCESS;
234 ATH_MSG_ERROR(
"MuonTriggerScaleFactors::getTriggerScaleFactor This is an experimental function. If you really know what you are doing set UseExperimental property.");
238 if (trigger.empty()) {
239 ATH_MSG_ERROR(
"MuonTriggerScaleFactors::getTriggerScaleFactor Trigger must have value.");
245 if (trigger ==
"HLT_mu8noL1")
246 ATH_MSG_WARNING(
"What you are trying to do is not correct. For di-muon triggers you should get the efficiency with getTriggerEfficiency and compute the SF by yourself.");
247 else if (trigger.find(
"HLT_2mu10") != std::string::npos || trigger.find(
"HLT_2mu14") != std::string::npos)
248 ATH_MSG_WARNING(
"Di-muon trigger scale factors for single reco muons are not supported!");
255 if (trigger.empty()) {
256 ATH_MSG_ERROR(
"MuonTriggerScaleFactors::getTriggerScaleFactor Trigger must have value.");
262 if (trigger ==
"HLT_mu8noL1") {
263 ATH_MSG_WARNING(
"What you are trying to do is not correct. For di-muon triggers you should get the efficiency with getTriggerEfficiency and compute the SF by yourself.");
265 else if (trigger.find(
"HLT_2mu10") != std::string::npos || trigger.find(
"HLT_2mu14") != std::string::npos) {
283 std::size_t
pos = sysBaseName.find(
"MCTOY");
284 if (
pos == std::string::npos)
return -1;
285 return atoi(sysBaseName.substr(
pos + 5,
pos + 8).c_str());
289 if (trigger.empty()) {
290 ATH_MSG_ERROR(
"MuonTriggerScaleFactors::getTriggerEfficiency Trigger must have value.");
310 std::string systype =
"";
312 systype =
"syst_down";
316 systype =
"stat_down";
330 if (configuration.
replicaIndex != -1) systype =
"replicas";
347 std::vector<TH1_Ptr> replica_v;
348 const int xbins =
h->GetNbinsX(),
ybins =
h->GetNbinsY();
350 for (
int t = 0;
t < nrep; ++
t) {
351 TH2* replica =
dynamic_cast<TH2*
>(
h->Clone(Form(
"rep%d_%s",
t,
h->GetName())));
353 for (
int x_i = 0; x_i <=
xbins; ++x_i) {
354 for (
int y_i = 0; y_i <=
ybins; ++y_i) {
355 replica->SetBinContent(x_i, y_i, Rndm.Gaus(
h->GetBinContent(x_i, y_i),
h->GetBinError(x_i, y_i)));
358 replica_v.push_back(
TH1_Ptr(replica));
365 return H1.get() !=
nullptr;
370 ATH_MSG_ERROR(
"MuonTriggerScaleFactors::getTriggerScaleFactor This is an experimental function. If you really know what you are doing set UseExperimental property.");
380 ATH_MSG_ERROR(
"Could not find efficiency map for muon with eta: " <<
mu_eta <<
" and phi: " <<
mu_phi <<
". Something is inconsistent. Please check your settings for year, mc and trigger." );
384 double mu_phi_corr =
mu_phi;
385 if (mu_phi_corr < eff_h2->GetYaxis()->GetXmin()) mu_phi_corr += 2.0 *
M_PI;
386 if (mu_phi_corr > eff_h2->GetYaxis()->GetXmax()) mu_phi_corr -= 2.0 *
M_PI;
387 return eff_h2->FindFixBin(
mu_eta, mu_phi_corr);
393 return std::hash<std::string>()(
histName);
428 ATH_MSG_WARNING(
"Could not find what you are looking for in the efficiency map. The trigger you are looking for, year and mc are not consistent, or the trigger is unavailable in this data period. Returning efficiency = 0.");
434 ATH_MSG_ERROR(
"Could not find what you are looking for in the efficiency map. The trigger you are looking for, year and mc are not consistent, or the trigger is unavailable in this data period. Please check how you set up the tool.");
439 if (configuration.
replicaIndex >= (
int) cit->second.size()) {
440 ATH_MSG_ERROR(
"MuonTriggerScaleFactors::getMuonEfficiency ; index for replicated histograms is out of range.");
447 if (cit.get() ==
nullptr) {
449 ATH_MSG_WARNING(
"Could not find what you are looking for in the efficiency map. The trigger you are looking for, year and mc are not consistent, or the trigger is unavailable in this data period. Returning efficiency = 0.");
453 ATH_MSG_ERROR(
"Could not find what you are looking for in the efficiency map. The trigger you are looking for, year and mc are not consistent, or the trigger is unavailable in this data period. Please check how you set up the tool.");
460 double mu_phi_corr =
mu_phi;
461 if (mu_phi_corr < eff_h2->GetYaxis()->GetXmin()) mu_phi_corr += 2.0 *
M_PI;
462 if (mu_phi_corr > eff_h2->GetYaxis()->GetXmax()) mu_phi_corr -= 2.0 *
M_PI;
464 const int bin = eff_h2->FindFixBin(
mu_eta, mu_phi_corr);
477 if (mucont.
size() != 2) {
478 ATH_MSG_FATAL(
"MuonTriggerScaleFactors::GetTriggerSF;Currently dimuon trigger chains only implemented for events with exactly 2 muons.");
482 Double_t eff_data = 0;
485 std::string data_err =
"";
486 std::string mc_err =
"";
489 data_err =
"nominal";
492 data_err =
"nominal";
493 mc_err =
"syst_down";
495 data_err =
"stat_down";
498 data_err =
"stat_up";
501 data_err =
"nominal";
512 if (configuration.
replicaIndex != -1) data_err =
"replicas";
515 configuration.
isData =
true;
519 configuration.
isData =
false;
524 double event_SF = 1.;
526 if (std::abs(1. - eff_mc) > 0.0001) {
527 event_SF = eff_data / eff_mc;
530 TriggerSF = event_SF;
539 double rate_not_fired_data = 1.;
540 double rate_not_fired_mc = 1.;
542 for (
const auto mu : mucont) {
544 double eff_data = 0., eff_mc = 0.;
553 std::string muon_trigger_name = trigger;
554 std::string data_err =
"";
555 std::string mc_err =
"";
564 data_err =
"nominal";
567 data_err =
"nominal";
568 mc_err =
"syst_down";
570 data_err =
"stat_down";
573 data_err =
"stat_up";
576 data_err =
"nominal";
587 if (configuration.
replicaIndex != -1) data_err =
"replicas";
590 configuration.
isData =
true;
593 configuration.
isData =
false;
598 rate_not_fired_data *= (1. - eff_data);
599 rate_not_fired_mc *= (1. - eff_mc);
602 double event_SF = 1.;
603 if (1 - rate_not_fired_data == 0) event_SF = 0;
604 if ((mucont.size()) and (std::abs(1. - rate_not_fired_mc) > 0.0001)) {
606 event_SF = (1. - rate_not_fired_data) / (1. - rate_not_fired_mc);
608 TriggerSF = event_SF;
619 double eff_data = 0., eff_mc = 0.;
628 std::string muon_trigger_name = trigger;
629 std::string data_err =
"";
630 std::string mc_err =
"";
639 data_err =
"nominal";
642 data_err =
"nominal";
643 mc_err =
"syst_down";
645 data_err =
"stat_down";
648 data_err =
"stat_up";
651 data_err =
"nominal";
657 if (configuration.
replicaIndex != -1) data_err =
"replicas";
660 configuration.
isData =
true;
664 configuration.
isData =
false;
671 if (std::abs(eff_mc) > 0.0001)
672 TriggerSF = eff_data / eff_mc;
702 if (trigger.find(
"2mu10") != std::string::npos)
return "HLT_mu10";
703 if (trigger.find(
"2mu14") != std::string::npos)
return "HLT_mu14";
704 throw std::runtime_error(
"Unknown dimuon trigger");
708 std::size_t
index = trigger.find(
"HLT_mu");
709 if (
index != std::string::npos) {
710 std::string rawNumber = trigger.substr(
index + 6);
711 if (!rawNumber.empty() && isdigit(rawNumber[0])) {
712 std::stringstream(rawNumber) >>
threshold;
718 ATH_MSG_ERROR(
"MuonTriggerScaleFactors::getThreshold Could not extract threshold for trigger " << trigger);
756 else if (
year == 2016) {
769 else if (
year == 2017) {
780 else if (
year == 2018) {
798 else if (
year == 2022) {
803 else if (
year == 2023) {
808 ATH_MSG_FATAL(
"RunNumber: " <<
runNumber <<
" not known! Will stop the code to prevent using wrong SFs.");
809 throw std::invalid_argument{
""};
815 if (
info.operator->()==
nullptr) {
817 throw std::invalid_argument{
""};
820 ATH_MSG_DEBUG(
"The current event is a data event. Return runNumber instead.");
821 return info->runNumber();
825 ATH_MSG_FATAL(
"Failed to find the RandomRunNumber decoration. Please call the apply() method from the PileupReweightingTool beforehand in order to get period dependent SFs");
826 throw std::invalid_argument{
""};
827 }
else if (acc_rnd(*
info) == 0) {
828 ATH_MSG_FATAL(
"Pile up tool has given runNumber 0. Exiting the code.");
829 throw std::invalid_argument{
""};
833 return acc_rnd(*
info);
838 TDirectory* tempDir = 0;
840 while (not tempDir) {
843 if (gROOT->GetDirectory((
dirname.str()).c_str())) {
847 tempDir = gROOT->mkdir((
dirname.str()).c_str());
849 ATH_MSG_ERROR(
"getTemporaryDirectory::Temporary directory could not be created");
859 if (!systematic.
empty()) {
861 return sys.find(systematic) !=
sys.end();
889 if (
registry.registerSystematics(*
this) != StatusCode::SUCCESS) {
890 ATH_MSG_ERROR(
"Failed to add systematic to list of recommended systematics.");
891 return StatusCode::FAILURE;
893 return StatusCode::SUCCESS;
912 ATH_MSG_ERROR(
"Unsupported combination of systematics passed to the tool!");
913 return StatusCode::FAILURE;
917 itr =
m_systFilter.insert(std::make_pair(systConfig, filteredSys)).first;
930 return StatusCode::FAILURE;
934 return StatusCode::SUCCESS;