-
Notifications
You must be signed in to change notification settings - Fork 9
Expand file tree
/
Copy pathelectron.h
More file actions
158 lines (146 loc) · 7.31 KB
/
electron.h
File metadata and controls
158 lines (146 loc) · 7.31 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
#pragma once
#include "correction.h"
#include "corrections.h"
namespace correction {
class EleCorrProvider : public CorrectionsBase<EleCorrProvider> {
public:
enum class UncSource : int {
Central = -1,
EleID = 0,
EleES = 1,
Ele_dEsigma = 2,
};
static std::string getESScaleStr(UncScale scale) {
static const std::map<UncScale, std::string> scale_names = {
{UncScale::Down, "scaledown"},
{UncScale::Up, "scaleup"},
};
return scale_names.at(scale);
}
static std::string getIDScaleStr(UncScale scale) {
static const std::map<UncScale, std::string> scale_names = {
{UncScale::Down, "sfdown"},
{UncScale::Central, "sf"},
{UncScale::Up, "sfup"},
};
return scale_names.at(scale);
}
static bool sourceApplies(UncSource source) {
if (source == UncSource::EleID)
return true;
if (source == UncSource::EleES)
return true;
return false;
}
EleCorrProvider(const std::string& EleIDFile,
const std::string& EleESFile,
const std::string& EleIDFile_key,
const std::string& EleESFile_key)
: corrections_(CorrectionSet::from_file(EleIDFile)),
correctionsES_(CorrectionSet::from_file(EleESFile)),
EleIDSF_(corrections_->at(EleIDFile_key)),
EleES_(correctionsES_->at(EleESFile_key)) {}
float getID_SF(const LorentzVectorM& Electron_p4,
std::string working_point,
std::string period,
UncSource source,
UncScale scale) const {
const UncScale jet_scale = sourceApplies(source) ? scale : UncScale::Central;
float value = 1.0;
if (period.starts_with("2023")) {
value = safeEvaluate(EleIDSF_,
period,
getIDScaleStr(jet_scale),
working_point,
Electron_p4.eta(),
Electron_p4.pt(),
Electron_p4.phi());
} else {
value = safeEvaluate(
EleIDSF_, period, getIDScaleStr(jet_scale), working_point, Electron_p4.eta(), Electron_p4.pt());
}
return value;
}
// https://gitlab.cern.ch/cms-analysis-corrections/EGM/examples/-/blob/latest/egmScaleAndSmearingExample.py?ref_type=heads
RVecLV getESEtDep_data(const RVecLV& Electron_p4,
const RVecI& Electron_genMatch,
const RVecUC& Electron_seedGain,
const RVecF& Electron_SCeta,
unsigned int run,
const RVecUC& Electron_r9,
UncSource source,
UncScale scale) const {
RVecLV final_p4 = Electron_p4;
for (size_t n = 0; n < Electron_p4.size(); ++n) {
const GenLeptonMatch genMatch = static_cast<GenLeptonMatch>(Electron_genMatch.at(n));
if (scale != UncScale::Central &&
(genMatch == GenLeptonMatch::Electron || genMatch == GenLeptonMatch::TauElectron)) {
double sf = Electron_p4[n].pt() < 15 ? 0. : correctionsES_->compound().at("Scale")->evaluate({"scale", // or SmearAndSyst ??? and smear?
static_cast<double>(run),
Electron_SCeta[n],
static_cast<double>(Electron_r9.at(n)),
Electron_p4[n].pt(),
Electron_seedGain.at(n)});
final_p4[n] *= 1 + static_cast<int>(scale) * sf;
}
}
return final_p4;
}
RVecLV getESEtDep_MC(const RVecLV& Electron_p4,
const RVecI& Electron_genMatch,
const RVecUC& Electron_seedGain,
const RVecF& Electron_SCeta,
unsigned int run,
const RVecUC& Electron_r9,
UncSource source,
UncScale scale) const {
TRandom3 rng(0);
RVecLV final_p4 = Electron_p4;
for (size_t n = 0; n < Electron_p4.size(); ++n) {
const GenLeptonMatch genMatch = static_cast<GenLeptonMatch>(Electron_genMatch.at(n));
if (scale != UncScale::Central &&
(genMatch == GenLeptonMatch::Electron || genMatch == GenLeptonMatch::TauElectron)) {
double smear = Electron_p4[n].pt() < 15 ? 0. : EleES_->evaluate({"smear", // or SmearAndSyst ??? and smear?
// static_cast<double>(run),
Electron_p4[n].pt(),
static_cast<double>(Electron_r9.at(n)),
Electron_SCeta[n]});
const double random_number = rng.Gaus(0.0, 1.0);
const double sf = 1.0 + static_cast<int>(scale) * smear * random_number;
final_p4[n] = LorentzVectorM(Electron_p4[n].pt() * sf, Electron_p4[n].eta(),Electron_p4[n].phi(),Electron_p4[n].M());
// const double energyErr_corr = std::sqrt(
// std::pow(mc_electrons_energyErr[i], 2) +
// std::pow(mc_electrons_energy[i] * smear, 2)
// ) * smearing;
}
}
return final_p4;
}
RVecLV getES(const RVecLV& Electron_p4,
const RVecI& Electron_genMatch,
const RVecUC& Electron_seedGain,
int run,
const RVecUC& Electron_r9,
UncSource source,
UncScale scale) const {
RVecLV final_p4 = Electron_p4;
for (size_t n = 0; n < Electron_p4.size(); ++n) {
const GenLeptonMatch genMatch = static_cast<GenLeptonMatch>(Electron_genMatch.at(n));
if (scale != UncScale::Central &&
(genMatch == GenLeptonMatch::Electron || genMatch == GenLeptonMatch::TauElectron)) {
double sf = EleES_->evaluate({"total_uncertainty",
static_cast<int>(Electron_seedGain.at(n)),
static_cast<double>(run),
Electron_p4[n].eta(),
static_cast<double>(Electron_r9.at(n)),
Electron_p4[n].pt()});
final_p4[n] *= 1 + static_cast<int>(scale) * sf;
}
}
return final_p4;
}
private:
std::unique_ptr<CorrectionSet> corrections_, correctionsES_;
Correction::Ref EleIDSF_, EleES_;
};
} //namespace correction