23 constexpr
float epsilon = 1
e-6;
25 const char*
const CONFIG_FILE_NAME =
26 "ElectronIsolationSelection/v1/CaloDeltaRZ.root";
27 const char*
const DATA_HISTO_NAME =
"hData";
28 const char*
const MC_HISTO_NAME =
"hMC";
36 using namespace ShowerDepthUtilMessaging;
41 m_hData = getHistoFromFile(
filename.c_str(), DATA_HISTO_NAME);
42 m_hMC = getHistoFromFile(
filename.c_str(), MC_HISTO_NAME);
43 if(!m_hData || !m_hMC){
44 throw std::runtime_error(
"Could not create instance of ShowerDepthUtil");
51 return getShowerDepthEM1(etas1) - getRZCorrection(etas1,
phi, isData);
57 return getShowerDepthEM2(etas2) - getRZCorrection(etas2,
phi, isData);
66 float radius, aetas1 = std::abs(etas1);
68 radius = (1558.859292 - 4.990838 * aetas1 - 21.144279 * aetas1 * aetas1);
69 }
else if (aetas1 < 1.5) {
70 radius = (1522.775373 + 27.970192 * aetas1 - 21.104108 * aetas1 * aetas1);
75 if (etas1 < 0. and aetas1 > 1.5) {
88 float radius, aetas2 = std::abs(etas2);
90 radius = (1698.990944 - 49.431767 * aetas2 - 24.504976 * aetas2 * aetas2);
91 }
else if (aetas2 < 1.5) {
92 radius = (8027.574119 - 2717.653528 * aetas2);
95 radius = (3473.473909 + 453.941515 * aetas2 - 119.101945 * aetas2 * aetas2);
97 if (etas2 < 0. and aetas2 > 1.5) {
105 float phi,
bool isData,
106 int sampling)
const {
107 std::pair<float, float> RZ = getCorrectedRZ(eta,
phi, isData, sampling);
108 return getEtaDirection(
zvertex, RZ.first, RZ.second);
112 if ((sampling != 1 && sampling != 2) || (std::abs(eta) > 10)) {
114 "Invalid sampling, eta: " << sampling <<
", " << eta);
115 return std::make_pair(0., 0.);
118 (sampling == 1 ? getShowerDepthEM1(eta) : getShowerDepthEM2(eta));
119 if (std::abs(eta) < 1.5)
120 return std::make_pair(
depth,
depth * std::sinh(eta));
121 return std::make_pair(
depth / std::sinh(eta),
depth);
125 float etas2,
float phi,
127 std::pair<float, float> RZ1 = getCorrectedRZ(etas1,
phi, isData, 1);
128 std::pair<float, float> RZ2 = getCorrectedRZ(etas2,
phi, isData, 2);
130 if (std::abs(RZ2.first - RZ1.first) < epsilon) {
134 return {
std::asinh((RZ2.second - RZ1.second) / (RZ2.first - RZ1.first))};
139 int sampling)
const {
140 if ((sampling != 1 && sampling != 2) || (std::abs(eta) > 10)) {
142 "Invalid sampling, eta: " << sampling <<
", " << eta);
143 return std::make_pair(0., 0.);
145 float depth = (sampling == 1 ? getCorrectedShowerDepthEM1(eta,
phi, isData)
146 : getCorrectedShowerDepthEM2(eta,
phi, isData));
147 if (std::abs(eta) < 1.5)
148 return std::make_pair(
depth,
depth * std::sinh(eta));
149 return std::make_pair(
depth / std::sinh(eta),
depth);
156 const TH1*
histo = (isData ? m_hData.get() : m_hMC.get());
161 const Int_t
etaBin =
histo->GetXaxis()->FindFixBin(eta);
163 const float etaOld = eta;
164 eta =
histo->GetXaxis()->GetBinLowEdge(1) + epsilon;
166 "Using eta " << eta <<
" instead of " << etaOld);
168 const float etaOld = eta;
169 eta =
histo->GetXaxis()->GetBinUpEdge(
histo->GetNbinsX()) - epsilon;
171 "Using eta " << eta <<
" instead of " << etaOld);
175 const float phiOld =
phi;
176 phi =
histo->GetYaxis()->GetBinLowEdge(1) + epsilon;
178 "Using phi " <<
phi <<
" instead of " << phiOld);
180 const float phiOld =
phi;
181 phi =
histo->GetYaxis()->GetBinUpEdge(
histo->GetNbinsY()) - epsilon;
183 "Using phi " <<
phi <<
" instead of " << phiOld);
186 return histo->Interpolate(eta,
phi);
194 const char* histoName) {
195 std::unique_ptr<TFile>
f(TFile::Open(
fileName,
"READ"));
198 "Could not open file: \"" <<
fileName <<
"\"");
201 TH1*
h =
dynamic_cast<TH1*
>(
f->Get(histoName));
204 << histoName <<
"\" from file: \""
209 h->SetDirectory(
nullptr);
210 return std::unique_ptr<TH1>(
h);