ATLAS Offline Software
TauCalibrateLC.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2025 CERN for the benefit of the ATLAS collaboration
3 */
4 
6 
9 
10 #include "TFile.h"
11 #include "TF1.h"
12 #include "TH1D.h"
13 #include <utility>
14 
15 #define GeV 1000
16 
17 /********************************************************************/
18 TauCalibrateLC::TauCalibrateLC(const std::string& name) :
20 }
21 
22 /********************************************************************/
24 }
25 
26 /********************************************************************/
28 
29  ATH_CHECK( m_vertexInputContainer.initialize() );
30 
31  std::string fullPath = find_file(m_calibrationFile);
32  ATH_MSG_INFO("Using calibration file: " << fullPath);
33 
34  std::unique_ptr<TFile> file(TFile::Open(fullPath.c_str(), "READ"));
35  if (!file) {
36  ATH_MSG_FATAL("Failed to open " << fullPath);
37  return StatusCode::FAILURE;
38  }
39 
40  // get the histogram defining eta binning
41  std::string key = "etaBinning";
42  TH1* histo = dynamic_cast<TH1*>(file->Get(key.c_str()));
43  if (histo) {
44  histo->SetDirectory(nullptr);
45  m_etaBinHist = std::unique_ptr<TH1>(histo);
46  }
47  else {
48  ATH_MSG_FATAL("Failed to get an object with key " << key);
49  return StatusCode::FAILURE;
50  }
51 
52  //retrieve number of eta bins from file
53  m_nEtaBins = m_etaBinHist->GetNbinsX();
54 
55  TString tmpSlopKey[s_nProngBins] = {"slopeNPV1P", "slopeNPV3P"};
56  TString tmpFuncBase[s_nProngBins] = {"OneP_Eta_", "MultiP_Eta_"};
57 
60 
62  m_calibFunc.resize(s_nProngBins);
63  for (int i = 0; i < s_nProngBins; ++i) {
64  m_calibFunc[i].resize(m_nEtaBins);
65  }
66 
67  for (int i = 0; i < s_nProngBins; i++) {
68  histo = dynamic_cast<TH1*>(file->Get(tmpSlopKey[i])); // get pile-up slope histograms
69  if (histo) {
70  histo->SetDirectory(nullptr);
71  m_slopeNPVHist[i] = std::unique_ptr<TH1>(histo);
72  }
73  else {
74  ATH_MSG_FATAL("Failed to get an object with key " << tmpSlopKey[i]);
75  return StatusCode::FAILURE;
76  }
77 
78  for (int j = 0; j < m_nEtaBins; j++) {
79  TString key = tmpFuncBase[i];
80  key += j;
81  TF1* fcn = dynamic_cast<TF1*>(file->Get(key));
82  if (fcn) {
83  m_calibFunc[i][j] = std::unique_ptr<TF1>(fcn);
84  }
85  else {
86  ATH_MSG_FATAL("Failed to get an object with key " << key);
87  return StatusCode::FAILURE;
88  }
89  }
90  }
91  m_averageNPV = m_slopeNPVHist[0]->GetBinContent(0); // underflow is the average number of reconstructed primary vertices
92  ATH_MSG_DEBUG("averageNPV = " << m_averageNPV);
93 
94  file->Close();
95 
96  return StatusCode::SUCCESS;
97 }
98 
99 /********************************************************************/
101 {
102  // get IntermediateAxis or DetectorAxis momentum
104 
105  // skip calibration if negative pt
106  if ( tau_p4.Pt() <= 0.) {
107  ATH_MSG_DEBUG("tau energy at LC scale is " << tau_p4.Pt()/GeV << "--> set energy=0.001");
108  tau.setP4(0.001, tau_p4.Eta(), tau_p4.Phi(), tau.m());
109  tau.setP4(xAOD::TauJetParameters::TauEnergyScale, tau.pt(), tau.eta(), tau.phi(), tau.m());
110  return StatusCode::SUCCESS;
111  }
112 
113  // energy calibration depends on number of tracks - 1p or Mp
114  int prongBin = 1; //Mp
115  if (tau.nTracks() <= 1) prongBin = 0; //1p
116 
117  double absEta = std::abs( tau_p4.Eta() );
118  int etaBin = std::as_const(*m_etaBinHist).GetXaxis()->FindBin(absEta) - 1;
119 
120  if (etaBin>=m_nEtaBins) etaBin = m_nEtaBins-1; // correction from last bin should be applied on all taus outside stored eta range
121 
122  int nVertex = 0;
123 
124  // Obtain pileup
126  if (!vertexInHandle.isValid()) {
127  ATH_MSG_ERROR ("Could not retrieve HiveDataObj with key " << vertexInHandle.key());
128  return StatusCode::FAILURE;
129  }
130  const xAOD::VertexContainer * vxContainer = vertexInHandle.cptr();
131  for (const auto *const vertex : *vxContainer) {
132  if (vertex->vertexType() == xAOD::VxType::PileUp) {
133  ++nVertex;
134  }
135  }
136  ATH_MSG_DEBUG("calculated nVertex " << nVertex );
137 
138  double calibConst = 1.;
139 
140  double slopeNPV = m_slopeNPVHist[prongBin]->GetBinContent(etaBin + 1);
141  double offset = slopeNPV * (nVertex - m_averageNPV);
142 
143  // pt (energy) response parameterized as a function of pileup-corrected pt (energy)
144  double energyLC = tau_p4.Pt()/GeV;
145 
146  if (energyLC - offset <= 0.) {
147  ATH_MSG_DEBUG("after pile-up correction energy would be = " << energyLC - offset << " --> setting offset=0 now!");
148  offset = 0.;
149  }
150 
151  // apply offset correction
152  double energyPileupCorr = energyLC - offset;
153 
154  if (energyPileupCorr > 0. and energyPileupCorr < 10000.) // from 0 to 10 TeV
155  {
156  calibConst = m_calibFunc[prongBin][etaBin]->Eval(energyPileupCorr);
157 
158  if (calibConst <= 0.) {
159  ATH_MSG_DEBUG("calibration constant = " << calibConst);
160  ATH_MSG_DEBUG("prongBin = " << prongBin);
161  ATH_MSG_DEBUG("etaBin = " << etaBin);
162  ATH_MSG_DEBUG("energyPileupCorr = " << energyPileupCorr);
163  ATH_MSG_DEBUG("energyLC = " << energyLC);
164  calibConst = 1.;
165  }
166  }
167 
168  double energyFinal = energyPileupCorr / calibConst;
169 
170  // tau.m() is 0 by convention in tauRecTools
171  // while mDetectorAxis and mIntermediateAxis are not forced to 0, mTauEnergyScale is
172  tau.setP4( energyFinal * GeV, tau_p4.Eta(), tau_p4.Phi(), tau.m());
173  tau.setP4(xAOD::TauJetParameters::TauEnergyScale, tau.pt(), tau.eta(), tau.phi(), tau.m());
174  ATH_MSG_DEBUG("Energy at LC scale = " << energyLC << " pile-up offset " << offset << " calib. const. = " << calibConst << " final energy = " << energyFinal);
175 
176  return StatusCode::SUCCESS;
177 }
ATH_MSG_FATAL
#define ATH_MSG_FATAL(x)
Definition: AthMsgStreamMacros.h:34
TauCalibrateLC.h
xAOD::TauJetParameters::IntermediateAxis
@ IntermediateAxis
Definition: TauDefs.h:338
ATH_MSG_INFO
#define ATH_MSG_INFO(x)
Definition: AthMsgStreamMacros.h:31
SG::ReadHandle::cptr
const_pointer_type cptr()
Dereference the pointer.
xAOD::TauJet_v3::eta
virtual double eta() const
The pseudorapidity ( ) of the particle.
SG::ReadHandle
Definition: StoreGate/StoreGate/ReadHandle.h:67
xAOD::TauJet_v3::m
virtual double m() const
The invariant mass of the particle.
CscCalibQuery.fullPath
string fullPath
Definition: CscCalibQuery.py:359
TauCalibrateLC::m_doVertexCorrection
Gaudi::Property< bool > m_doVertexCorrection
Definition: TauCalibrateLC.h:49
xAOD::TauJet_v3::nTracks
size_t nTracks(TauJetParameters::TauTrackFlag flag=TauJetParameters::TauTrackFlag::classifiedCharged) const
Definition: TauJet_v3.cxx:488
TauCalibrateLC::s_nProngBins
static const int s_nProngBins
Definition: TauCalibrateLC.h:39
TauRecToolBase
The base class for all tau tools.
Definition: TauRecToolBase.h:21
xAOD::TauJet_v3::setP4
void setP4(double pt, double eta, double phi, double m)
Set methods for IParticle values.
Definition: TauJet_v3.cxx:171
TauCalibrateLC::TauCalibrateLC
TauCalibrateLC(const std::string &name="TauCalibrateLC")
Definition: TauCalibrateLC.cxx:18
TauCalibrateLC::m_calibrationFile
Gaudi::Property< std::string > m_calibrationFile
Definition: TauCalibrateLC.h:48
xAOD::TauJet_v3::pt
virtual double pt() const
The transverse momentum ( ) of the particle.
xAOD::TauJet_v3::phi
virtual double phi() const
The azimuthal angle ( ) of the particle.
xAOD::etaBin
setSAddress setEtaMS setDirPhiMS setDirZMS setBarrelRadius setEndcapAlpha setEndcapRadius setInterceptInner setEtaMap etaBin
Definition: L2StandAloneMuon_v1.cxx:149
TauCalibrateLC::m_etaBinHist
std::unique_ptr< TH1 > m_etaBinHist
Definition: TauCalibrateLC.h:43
ATH_MSG_ERROR
#define ATH_MSG_ERROR(x)
Definition: AthMsgStreamMacros.h:33
TauCalibrateLC::initialize
virtual StatusCode initialize() override
Tool initializer.
Definition: TauCalibrateLC.cxx:27
lumiFormat.i
int i
Definition: lumiFormat.py:85
TauCalibrateLC::m_vertexInputContainer
SG::ReadHandleKey< xAOD::VertexContainer > m_vertexInputContainer
Definition: TauCalibrateLC.h:51
EL::StatusCode
::StatusCode StatusCode
StatusCode definition for legacy code.
Definition: PhysicsAnalysis/D3PDTools/EventLoop/EventLoop/StatusCode.h:22
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x)
Definition: AthMsgStreamMacros.h:29
xAOD::TauJet_v3
Class describing a tau jet.
Definition: TauJet_v3.h:41
TauCalibrateLC::m_nEtaBins
int m_nEtaBins
Definition: TauCalibrateLC.h:45
file
TFile * file
Definition: tile_monitor.h:29
TauCalibrateLC::execute
virtual StatusCode execute(xAOD::TauJet &tau) const override
Execute - called for each tau candidate.
Definition: TauCalibrateLC.cxx:100
ATH_CHECK
#define ATH_CHECK
Definition: AthCheckMacros.h:40
DataVector
Derived DataVector<T>.
Definition: DataVector.h:794
SG::ReadHandle::isValid
virtual bool isValid() override final
Can the handle be successfully dereferenced?
xAOD::TauJetParameters::TauEnergyScale
@ TauEnergyScale
Definition: TauDefs.h:339
TauCalibrateLC::m_averageNPV
double m_averageNPV
Definition: TauCalibrateLC.h:46
xAOD::VxType::PileUp
@ PileUp
Pile-up vertex.
Definition: TrackingPrimitives.h:574
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:240
ReadHandle.h
Handle class for reading from StoreGate.
TauRecToolBase::find_file
std::string find_file(const std::string &fname) const
Definition: TauRecToolBase.cxx:19
TauCalibrateLC::m_slopeNPVHist
std::vector< std::unique_ptr< TH1 > > m_slopeNPVHist
Definition: TauCalibrateLC.h:42
SG::VarHandleBase::key
virtual const std::string & key() const override final
Return the StoreGate ID for the referenced object.
Definition: AthToolSupport/AsgDataHandles/Root/VarHandleBase.cxx:64
TauCalibrateLC::~TauCalibrateLC
~TauCalibrateLC()
Definition: TauCalibrateLC.cxx:23
Trk::vertex
@ vertex
Definition: MeasurementType.h:21
xAOD::TauJet_v3::p4
virtual FourMom_t p4() const
The full 4-momentum of the particle.
Definition: TauJet_v3.cxx:96
GeV
#define GeV
Definition: TauCalibrateLC.cxx:15
ReadDecorHandle.h
Handle class for reading a decoration on an object.
convertTimingResiduals.offset
offset
Definition: convertTimingResiduals.py:71
TauGNNUtils::Variables::absEta
bool absEta(const xAOD::TauJet &tau, double &out)
Definition: TauGNNUtils.cxx:245
xAOD::TauJetParameters::DetectorAxis
@ DetectorAxis
Definition: TauDefs.h:337
plotBeamSpotCompare.histo
histo
Definition: plotBeamSpotCompare.py:414
TauCalibrateLC::m_calibFunc
std::vector< std::vector< std::unique_ptr< TF1 > > > m_calibFunc
Definition: TauCalibrateLC.h:41
mapkey::key
key
Definition: TElectronEfficiencyCorrectionTool.cxx:37