ATLAS Offline Software
ShowerDepthTool.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 // Local include(s).
7 
8 // Project include(s).
11 
12 // ROOT include(s).
13 #include <TFile.h>
14 #include <TH1.h>
15 
16 // System include(s).
17 #include <cmath>
18 #include <string>
19 
20 namespace CP {
21 
22 // Calibration file / histogram name(s).
23 static const char* const CONFIG_FILE_NAME =
24  "ElectronIsolationSelection/v1/CaloDeltaRZ.root";
25 static const char* const DATA_HISTO_NAME = "hData";
26 static const char* const MC_HISTO_NAME = "hMC";
27 
28 // Make the messaging functions available.
29 ANA_MSG_SOURCE(ShowerDepthToolMessaging, "CP::ShowerDepthTool");
30 using namespace ShowerDepthToolMessaging;
32 
34 
36  const std::string filename = PathResolverFindCalibFile(CONFIG_FILE_NAME);
37 
38  m_hData = getHistoFromFile(filename.c_str(), DATA_HISTO_NAME);
39  m_hMC = getHistoFromFile(filename.c_str(), MC_HISTO_NAME);
40 
41  return (m_hData && m_hMC);
42 }
43 
45 float ShowerDepthTool::getCorrectedShowerDepthEM1(float etas1, float phi,
46  bool isData) const {
47  return getShowerDepthEM1(etas1) - getRZCorrection(etas1, phi, isData);
48 }
49 
51 float ShowerDepthTool::getCorrectedShowerDepthEM2(float etas2, float phi,
52  bool isData) const {
53  return getShowerDepthEM2(etas2) - getRZCorrection(etas2, phi, isData);
54 }
55 
62  float radius, aetas1 = std::abs(etas1);
63  if (aetas1 < 0.8) {
64  radius = (1558.859292 - 4.990838 * aetas1 - 21.144279 * aetas1 * aetas1);
65  } else if (aetas1 < 1.5) {
66  radius = (1522.775373 + 27.970192 * aetas1 - 21.104108 * aetas1 * aetas1);
67  } else {
68  // endcap so Z
69  radius = 3790.671754;
70  }
71  if (etas1 < 0. and aetas1 > 1.5) {
72  // negative endcap
73  return -radius;
74  }
75  return radius;
76 }
77 
84  float radius, aetas2 = std::abs(etas2);
85  if (aetas2 < 1.425) { // Barrel
86  radius = (1698.990944 - 49.431767 * aetas2 - 24.504976 * aetas2 * aetas2);
87  } else if (aetas2 < 1.5) { // EME2 in tool
88  radius = (8027.574119 - 2717.653528 * aetas2);
89  } else {
90  // endcap so Z
91  radius = (3473.473909 + 453.941515 * aetas2 - 119.101945 * aetas2 * aetas2);
92  }
93  if (etas2 < 0. and aetas2 > 1.5) {
94  // negative endcap
95  return -radius;
96  }
97  return radius;
98 }
99 
101  float phi, bool isData,
102  int sampling) const {
103  std::pair<float, float> RZ = getCorrectedRZ(eta, phi, isData, sampling);
104  return getEtaDirection(zvertex, RZ.first, RZ.second);
105 }
106 
107 std::pair<float, float> ShowerDepthTool::getRZ(float eta, int sampling) {
108  if ((sampling != 1 && sampling != 2) || (std::abs(eta) > 10)) {
109  ANA_MSG_LVL_SERIOUS(MSG::WARNING,
110  "Invalid sampling, eta: " << sampling << ", " << eta);
111  return std::make_pair(0., 0.);
112  }
113  float depth =
114  (sampling == 1 ? getShowerDepthEM1(eta) : getShowerDepthEM2(eta));
115  if (std::abs(eta) < 1.5)
116  return std::make_pair(depth, depth * std::sinh(eta));
117  return std::make_pair(depth / std::sinh(eta), depth);
118 }
119 
120 std::optional<float> ShowerDepthTool::getCaloPointingEta(float etas1,
121  float etas2, float phi,
122  bool isData) const {
123  std::pair<float, float> RZ1 = getCorrectedRZ(etas1, phi, isData, 1);
124  std::pair<float, float> RZ2 = getCorrectedRZ(etas2, phi, isData, 2);
125 
126  // Sanity check
127  constexpr float epsilon = 1e-6;
128  if (std::abs(RZ2.first - RZ1.first) < epsilon){
129  return std::nullopt;
130  }
131 
132  return {std::asinh((RZ2.second - RZ1.second) / (RZ2.first - RZ1.first))};
133 }
134 
135 std::pair<float, float> ShowerDepthTool::getCorrectedRZ(float eta, float phi,
136  bool isData,
137  int sampling) const {
138  if ((sampling != 1 && sampling != 2) || (std::abs(eta) > 10)) {
139  ANA_MSG_LVL_SERIOUS(MSG::WARNING,
140  "Invalid sampling, eta: " << sampling << ", " << eta);
141  return std::make_pair(0., 0.);
142  }
143  float depth = (sampling == 1 ? getCorrectedShowerDepthEM1(eta, phi, isData)
144  : getCorrectedShowerDepthEM2(eta, phi, isData));
145  if (std::abs(eta) < 1.5)
146  return std::make_pair(depth, depth * std::sinh(eta));
147  return std::make_pair(depth / std::sinh(eta), depth);
148 }
149 
151 float ShowerDepthTool::getRZCorrection(float eta, float phi,
152  bool isData) const {
153  // Get the correct histogram.
154  const TH1* histo = (isData ? m_hData.get() : m_hMC.get());
155  if (!histo) {
156  return 0;
157  }
158  // Make sure that we can perform the interpolation in both eta and phi.
159  // Note that std::numeric_limits<float>::epsilon() is just not large enough
160  // for the following. :-(
161  static constexpr float epsilon = 1e-6f;
162  const Int_t etaBin = histo->GetXaxis()->FindFixBin(eta);
163  if (etaBin < 1) {
164  const float etaOld = eta;
165  eta = histo->GetXaxis()->GetBinLowEdge(1) + epsilon;
166  ANA_MSG_LVL_SERIOUS(MSG::WARNING,
167  "Using eta " << eta << " instead of " << etaOld);
168  } else if (etaBin > histo->GetNbinsX()) {
169  const float etaOld = eta;
170  eta = histo->GetXaxis()->GetBinUpEdge(histo->GetNbinsX()) - epsilon;
171  ANA_MSG_LVL_SERIOUS(MSG::WARNING,
172  "Using eta " << eta << " instead of " << etaOld);
173  }
174  const Int_t phiBin = histo->GetYaxis()->FindFixBin(phi);
175  if (phiBin < 1) {
176  const float phiOld = phi;
177  phi = histo->GetYaxis()->GetBinLowEdge(1) + epsilon;
178  ANA_MSG_LVL_SERIOUS(MSG::WARNING,
179  "Using phi " << phi << " instead of " << phiOld);
180  } else if (phiBin > histo->GetNbinsY()) {
181  const float phiOld = phi;
182  phi = histo->GetYaxis()->GetBinUpEdge(histo->GetNbinsY()) - epsilon;
183  ANA_MSG_LVL_SERIOUS(MSG::WARNING,
184  "Using phi " << phi << " instead of " << phiOld);
185  }
186  // Get the correction as an interpolation.
187  return histo->Interpolate(eta, phi);
188 }
189 
190 float ShowerDepthTool::getEtaDirection(float zvertex, float R, float z) {
191  return std::asinh((z - zvertex) / R);
192 }
193 
194 std::unique_ptr<TH1> ShowerDepthTool::getHistoFromFile(const char* fileName,
195  const char* histoName) {
196  std::unique_ptr<TFile> f(TFile::Open(fileName, "READ"));
197  if (!f) {
198  ANA_MSG_LVL_SERIOUS(MSG::WARNING,
199  "Could not open file: \"" << fileName << "\"");
200  return {};
201  }
202  TH1* h = dynamic_cast<TH1*>(f->Get(histoName));
203  if (!h) {
204  ANA_MSG_LVL_SERIOUS(MSG::WARNING, "Could not get histogram: \""
205  << histoName << "\" from file: \""
206  << fileName << "\"");
207  return {};
208  }
209  // The file we be deleted so use SetDirectory
210  h->SetDirectory(nullptr);
211  return std::unique_ptr<TH1>(h);
212 }
213 
214 } // namespace CP
AllowedVariables::e
e
Definition: AsgElectronSelectorTool.cxx:37
CP::ShowerDepthTool::getShowerDepthEM1
static float getShowerDepthEM1(float etas1)
Shower depth (in mm) vs.
Definition: ShowerDepthTool.cxx:61
egammaParameters::depth
@ depth
pointing depth of the shower as calculated in egammaqgcld
Definition: egammaParamDefs.h:276
CP::ShowerDepthTool::getCaloPointingEta
std::optional< float > getCaloPointingEta(float etas1, float etas2, float phi, bool isData=true) const
Eta direction from samplings 1 and 2 (pointing)
Definition: ShowerDepthTool.cxx:120
ShowerDepthTool.h
CP::ShowerDepthTool::getRZCorrection
float getRZCorrection(float eta, float phi, bool isData=true) const
Return the calorimeter displacement in R(Z) for barrel (endcap)
Definition: ShowerDepthTool.cxx:151
egammaParameters::zvertex
@ zvertex
pointing z at vertex reconstructed from the cluster
Definition: egammaParamDefs.h:270
CP::ShowerDepthTool::getCorrectedShowerDepthEM1
float getCorrectedShowerDepthEM1(float etas1, float phi, bool isData=true) const
Shower depth (in mm) on EM1 vs.
Definition: ShowerDepthTool.cxx:45
CP
Select isolated Photons, Electrons and Muons.
Definition: Control/xAODRootAccess/xAODRootAccess/TEvent.h:48
CP::ShowerDepthTool::getCorrectedRZ
std::pair< float, float > getCorrectedRZ(float eta, float phi, bool isData=true, int sampling=1) const
Return the shower depth in R,Z considering misalignments.
Definition: ShowerDepthTool.cxx:135
xAOD::phi
setEt phi
Definition: TrigEMCluster_v1.cxx:29
xAOD::etaBin
setSAddress setEtaMS setDirPhiMS setDirZMS setBarrelRadius setEndcapAlpha setEndcapRadius setInterceptInner setEtaMap etaBin
Definition: L2StandAloneMuon_v1.cxx:148
FortranAlgorithmOptions.fileName
fileName
Definition: FortranAlgorithmOptions.py:13
z
#define z
h
CP::ShowerDepthTool::getRZ
static std::pair< float, float > getRZ(float eta, int sampling)
Shower depth in R,Z for the given sampling.
Definition: ShowerDepthTool.cxx:107
MessageCheck.h
macros for messaging and checking status codes
CP::ShowerDepthTool::getCorrectedShowerDepthEM2
float getCorrectedShowerDepthEM2(float etas2, float phi, bool isData=true) const
Shower depth (in mm) on EM2 vs.
Definition: ShowerDepthTool.cxx:51
hist_file_dump.f
f
Definition: hist_file_dump.py:135
python.Dumpers.asinh
def asinh(x)
helper methods ---------------------------------------------------------—
Definition: Dumpers.py:89
CP::ShowerDepthTool::ShowerDepthTool
ShowerDepthTool()
CP::ShowerDepthTool::~ShowerDepthTool
~ShowerDepthTool()
PathResolver.h
xAOD::phiBin
setSAddress setEtaMS setDirPhiMS setDirZMS setBarrelRadius setEndcapAlpha setEndcapRadius setPhiMap phiBin
Definition: L2StandAloneMuon_v2.cxx:144
ANA_MSG_LVL_SERIOUS
#define ANA_MSG_LVL_SERIOUS(lvl, xmsg)
Macro used to print "serious" messages.
Definition: Control/AthToolSupport/AsgMessaging/AsgMessaging/MessageCheck.h:270
ParticleGun_SamplingFraction.radius
radius
Definition: ParticleGun_SamplingFraction.py:96
CP::ShowerDepthTool::getShowerDepthEM2
static float getShowerDepthEM2(float etas2)
Shower depth (in mm) vs.
Definition: ShowerDepthTool.cxx:83
PathResolverFindCalibFile
std::string PathResolverFindCalibFile(const std::string &logical_file_name)
Definition: PathResolver.cxx:431
CP::ShowerDepthTool::getCorrectedEtaDirection
float getCorrectedEtaDirection(float zvertex, float eta, float phi, bool isData=true, int sampling=1) const
Eta direction from zvertex to the shower in the given sampling.
Definition: ShowerDepthTool.cxx:100
CP::ANA_MSG_SOURCE
ANA_MSG_SOURCE(ShowerDepthToolMessaging, "CP::ShowerDepthTool")
CP::ShowerDepthTool::getEtaDirection
static float getEtaDirection(float zvertex, float R, float z)
Definition: ShowerDepthTool.cxx:190
CaloCellTimeCorrFiller.filename
filename
Definition: CaloCellTimeCorrFiller.py:24
plotBeamSpotCompare.histo
histo
Definition: plotBeamSpotCompare.py:415
CP::ShowerDepthTool::initialize
bool initialize()
Function initialising the tool.
Definition: ShowerDepthTool.cxx:35
CP::ShowerDepthTool::getHistoFromFile
static std::unique_ptr< TH1 > getHistoFromFile(const char *fileName, const char *histoName)
Return TH1* from file given fileName, histoName.
Definition: ShowerDepthTool.cxx:194