CMS 3D CMS Logo
PFRecoTauEnergyAlgorithmPlugin.cc
Go to the documentation of this file.
1 /*
2  * =============================================================================
3  * Filename: RecoTauEnergyAlgorithmPlugin.cc
4  *
5  * Description: Determine best estimate for tau energy
6  * for tau candidates reconstructed in different decay modes
7  *
8  * Created: 04/09/2013 11:40:00
9  *
10  * Authors: Christian Veelken (LLR)
11  *
12  * =============================================================================
13  */
14 
15 #include "FWCore/MessageLogger/interface/MessageLogger.h"
16 
17 #include "RecoTauTag/RecoTau/interface/RecoTauBuilderPlugins.h"
18 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
19 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidateFwd.h"
20 #include "DataFormats/PatCandidates/interface/PackedCandidate.h"
21 #include "DataFormats/JetReco/interface/PFJet.h"
22 #include "DataFormats/TrackReco/interface/Track.h"
23 #include "DataFormats/Common/interface/Ptr.h"
24 #include "DataFormats/Common/interface/RefToPtr.h"
25 #include "DataFormats/TauReco/interface/PFRecoTauChargedHadron.h"
26 #include "DataFormats/TauReco/interface/PFRecoTauChargedHadronFwd.h"
27 #include "RecoTauTag/RecoTau/interface/pfRecoTauChargedHadronAuxFunctions.h"
28 
29 #include "DataFormats/Math/interface/deltaR.h"
30 #include "FWCore/Utilities/interface/isFinite.h"
31 
32 #include <vector>
33 #include <cmath>
34 
35 namespace reco {
36  namespace tau {
37 
38  class PFRecoTauEnergyAlgorithmPlugin : public RecoTauModifierPlugin {
39  public:
40  explicit PFRecoTauEnergyAlgorithmPlugin(const edm::ParameterSet&, edm::ConsumesCollector&& iC);
41  ~PFRecoTauEnergyAlgorithmPlugin() override;
42  void operator()(reco::PFTau&) const override;
43  void beginEvent() override;
44  void endEvent() override;
45 
46  private:
47  double dRaddNeutralHadron_;
48  double minNeutralHadronEt_;
49  double dRaddPhoton_;
50  double minGammaEt_;
51 
52  int verbosity_;
53  };
54 
55  PFRecoTauEnergyAlgorithmPlugin::PFRecoTauEnergyAlgorithmPlugin(const edm::ParameterSet& cfg,
56  edm::ConsumesCollector&& iC)
57  : RecoTauModifierPlugin(cfg, std::move(iC)),
58  dRaddNeutralHadron_(cfg.getParameter<double>("dRaddNeutralHadron")),
59  minNeutralHadronEt_(cfg.getParameter<double>("minNeutralHadronEt")),
60  dRaddPhoton_(cfg.getParameter<double>("dRaddPhoton")),
61  minGammaEt_(cfg.getParameter<double>("minGammaEt")),
62  verbosity_(cfg.getParameter<int>("verbosity")) {}
63 
64  PFRecoTauEnergyAlgorithmPlugin::~PFRecoTauEnergyAlgorithmPlugin() {}
65 
66  void PFRecoTauEnergyAlgorithmPlugin::beginEvent() {}
67 
68  namespace {
69  double getTrackPerr2(const reco::Track& track) {
70  double trackPerr = track.p() * (track.ptError() / track.pt());
71  return trackPerr * trackPerr;
72  }
73 
74  void updateTauP4(reco::PFTau& tau, double sf, const reco::Candidate::LorentzVector& addP4) {
75  // preserve tau candidate mass when adding extra neutral energy
76  double tauPx_modified = tau.px() + sf * addP4.px();
77  double tauPy_modified = tau.py() + sf * addP4.py();
78  double tauPz_modified = tau.pz() + sf * addP4.pz();
79  double tauMass = tau.mass();
80  double tauEn_modified = sqrt(tauPx_modified * tauPx_modified + tauPy_modified * tauPy_modified +
81  tauPz_modified * tauPz_modified + tauMass * tauMass);
82  reco::Candidate::LorentzVector tauP4_modified(tauPx_modified, tauPy_modified, tauPz_modified, tauEn_modified);
83  tau.setP4(tauP4_modified);
84  }
85 
86  void killTau(reco::PFTau& tau) {
87  reco::Candidate::LorentzVector tauP4_modified(0., 0., 0., 0.);
88  tau.setP4(tauP4_modified);
89  tau.setStatus(-1);
90  }
91  } // namespace
92 
93  template <class Base, class Der>
94  bool isPtrEqual(const edm::Ptr<Base>& b, const edm::Ptr<Der>& d) {
95  return edm::Ptr<Der>(b) == d;
96  }
97 
98  template <class Base>
99  bool isPtrEqual(const edm::Ptr<Base>& b, const edm::Ptr<Base>& d) {
100  return b == d;
101  }
102 
103  void PFRecoTauEnergyAlgorithmPlugin::operator()(reco::PFTau& tau) const {
104  if (verbosity_) {
105  std::cout << "<PFRecoTauEnergyAlgorithmPlugin::operator()>:" << std::endl;
106  std::cout << "tau: Pt = " << tau.pt() << ", eta = " << tau.eta() << ", phi = " << tau.phi()
107  << " (En = " << tau.energy() << ", decayMode = " << tau.decayMode() << ")" << std::endl;
108  }
109 
110  // Add high Pt PFNeutralHadrons and PFGammas that are not "used" by tau decay mode object
111  std::vector<reco::CandidatePtr> addNeutrals;
112  reco::Candidate::LorentzVector addNeutralsSumP4;
113  const auto& jetConstituents = tau.jetRef()->daughterPtrVector();
114  for (const auto& jetConstituent : jetConstituents) {
115  int jetConstituentPdgId = std::abs(jetConstituent->pdgId());
116  if (!((jetConstituentPdgId == 130 && jetConstituent->et() > minNeutralHadronEt_) ||
117  (jetConstituentPdgId == 22 && jetConstituent->et() > minGammaEt_)))
118  continue;
119 
120  bool isSignalCand = false;
121  const auto& signalCands = tau.signalCands();
122  for (const auto& signalCand : signalCands) {
123  if (isPtrEqual(jetConstituent, signalCand))
124  isSignalCand = true;
125  }
126  if (isSignalCand)
127  continue;
128 
129  double dR = deltaR(jetConstituent->p4(), tau.p4());
130  double dRadd = -1.;
131  if (jetConstituentPdgId == 130)
132  dRadd = dRaddNeutralHadron_;
133  else if (jetConstituentPdgId == 22)
134  dRadd = dRaddPhoton_;
135  if (dR < dRadd) {
136  addNeutrals.push_back(jetConstituent);
137  addNeutralsSumP4 += jetConstituent->p4();
138  }
139  }
140  if (verbosity_) {
141  std::cout << "addNeutralsSumP4: En = " << addNeutralsSumP4.energy() << std::endl;
142  }
143 
144  unsigned numNonPFCandTracks = 0;
145  double nonPFCandTracksSumP = 0.;
146  double nonPFCandTracksSumPerr2 = 0.;
147  const std::vector<PFRecoTauChargedHadron>& chargedHadrons = tau.signalTauChargedHadronCandidatesRestricted();
148  for (std::vector<PFRecoTauChargedHadron>::const_iterator chargedHadron = chargedHadrons.begin();
149  chargedHadron != chargedHadrons.end();
150  ++chargedHadron) {
151  if (chargedHadron->algoIs(PFRecoTauChargedHadron::kTrack)) {
152  ++numNonPFCandTracks;
153  const reco::Track* chargedHadronTrack = getTrackFromChargedHadron(*chargedHadron);
154  if (chargedHadronTrack != nullptr) {
155  nonPFCandTracksSumP += chargedHadronTrack->p();
156  nonPFCandTracksSumPerr2 += getTrackPerr2(*chargedHadronTrack);
157  }
158  }
159  }
160  if (verbosity_) {
161  std::cout << "nonPFCandTracksSumP = " << nonPFCandTracksSumP << " +/- " << sqrt(nonPFCandTracksSumPerr2)
162  << " (numNonPFCandTracks = " << numNonPFCandTracks << ")" << std::endl;
163  }
164 
165  if (numNonPFCandTracks == 0) {
166  // This is the easy case:
167  // All tau energy is taken from PFCandidates reconstructed by PFlow algorithm
168  // and there is no issue with double-counting of energy.
169  if (verbosity_) {
170  std::cout << "easy case: all tracks are associated to PFCandidates --> leaving tau momentum untouched."
171  << std::endl;
172  }
173  updateTauP4(tau, 1., addNeutralsSumP4);
174  return;
175  } else {
176  // This is the difficult case:
177  // The tau energy needs to be computed for an arbitrary mix of charged and neutral PFCandidates plus reco::Tracks.
178  // We need to make sure not to double-count energy deposited by reco::Track in ECAL and/or HCAL as neutral PFCandidates.
179 
180  // Check if we have enough energy in collection of PFNeutralHadrons and PFGammas that are not "used" by tau decay mode object
181  // to balance track momenta:
182  if (nonPFCandTracksSumP < addNeutralsSumP4.energy()) {
183  double scaleFactor = 1. - nonPFCandTracksSumP / addNeutralsSumP4.energy();
184  if (!(scaleFactor >= 0. && scaleFactor <= 1.)) {
185  edm::LogWarning("PFRecoTauEnergyAlgorithmPlugin::operator()")
186  << "Failed to compute tau energy --> killing tau candidate !!" << std::endl;
187  killTau(tau);
188  return;
189  }
190  if (verbosity_) {
191  std::cout << "case (2): addNeutralsSumEn > nonPFCandTracksSumP --> adjusting tau momentum." << std::endl;
192  }
193  updateTauP4(tau, scaleFactor, addNeutralsSumP4);
194  return;
195  }
196 
197  // Determine which neutral PFCandidates are close to PFChargedHadrons
198  // and have been merged into ChargedHadrons
199  std::vector<reco::CandidatePtr> mergedNeutrals;
200  reco::Candidate::LorentzVector mergedNeutralsSumP4;
201  for (std::vector<PFRecoTauChargedHadron>::const_iterator chargedHadron = chargedHadrons.begin();
202  chargedHadron != chargedHadrons.end();
203  ++chargedHadron) {
204  if (chargedHadron->algoIs(PFRecoTauChargedHadron::kTrack)) {
205  const std::vector<reco::CandidatePtr>& neutralPFCands = chargedHadron->getNeutralPFCandidates();
206  for (std::vector<reco::CandidatePtr>::const_iterator neutralPFCand = neutralPFCands.begin();
207  neutralPFCand != neutralPFCands.end();
208  ++neutralPFCand) {
209  mergedNeutrals.push_back(*neutralPFCand);
210  mergedNeutralsSumP4 += (*neutralPFCand)->p4();
211  }
212  }
213  }
214  if (verbosity_) {
215  std::cout << "mergedNeutralsSumP4: En = " << mergedNeutralsSumP4.energy() << std::endl;
216  }
217 
218  // Check if track momenta are balanced by sum of PFNeutralHadrons and PFGammas that are not "used" by tau decay mode object
219  // plus neutral PFCandidates close to PFChargedHadrons:
220  if (nonPFCandTracksSumP < (addNeutralsSumP4.energy() + mergedNeutralsSumP4.energy())) {
221  double scaleFactor = ((addNeutralsSumP4.energy() + mergedNeutralsSumP4.energy()) - nonPFCandTracksSumP) /
222  mergedNeutralsSumP4.energy();
223  if (!(scaleFactor >= 0. && scaleFactor <= 1.)) {
224  edm::LogWarning("PFRecoTauEnergyAlgorithmPlugin::operator()")
225  << "Failed to compute tau energy --> killing tau candidate !!" << std::endl;
226  killTau(tau);
227  return;
228  }
229  reco::Candidate::LorentzVector diffP4;
230  size_t numChargedHadrons = chargedHadrons.size();
231  for (size_t iChargedHadron = 0; iChargedHadron < numChargedHadrons; ++iChargedHadron) {
232  const PFRecoTauChargedHadron& chargedHadron = chargedHadrons[iChargedHadron];
233  if (!chargedHadron.getNeutralPFCandidates().empty()) {
234  PFRecoTauChargedHadron chargedHadron_modified = chargedHadron;
235  setChargedHadronP4(chargedHadron_modified, scaleFactor);
236  tau.signalTauChargedHadronCandidatesRestricted()[iChargedHadron] = chargedHadron_modified;
237  diffP4 += (chargedHadron.p4() - chargedHadron_modified.p4());
238  }
239  }
240  if (verbosity_) {
241  std::cout << "case (3): (addNeutralsSumEn + mergedNeutralsSumEn) > nonPFCandTracksSumP --> adjusting tau "
242  "momentum."
243  << std::endl;
244  }
245  updateTauP4(tau, -1., diffP4);
246  return;
247  }
248 
249  // Determine energy sum of all PFNeutralHadrons interpreted as ChargedHadrons with missing track
250  unsigned numChargedHadronNeutrals = 0;
251  std::vector<reco::CandidatePtr> chargedHadronNeutrals;
252  reco::Candidate::LorentzVector chargedHadronNeutralsSumP4;
253  for (std::vector<PFRecoTauChargedHadron>::const_iterator chargedHadron = chargedHadrons.begin();
254  chargedHadron != chargedHadrons.end();
255  ++chargedHadron) {
256  if (chargedHadron->algoIs(PFRecoTauChargedHadron::kPFNeutralHadron)) {
257  ++numChargedHadronNeutrals;
258  chargedHadronNeutrals.push_back(chargedHadron->getChargedPFCandidate());
259  chargedHadronNeutralsSumP4 += chargedHadron->getChargedPFCandidate()->p4();
260  }
261  }
262  if (verbosity_) {
263  std::cout << "chargedHadronNeutralsSumP4: En = " << chargedHadronNeutralsSumP4.energy()
264  << " (numChargedHadronNeutrals = " << numChargedHadronNeutrals << ")" << std::endl;
265  }
266 
267  // Check if sum of PFNeutralHadrons and PFGammas that are not "used" by tau decay mode object
268  // plus neutral PFCandidates close to PFChargedHadrons plus PFNeutralHadrons interpreted as ChargedHadrons with missing track balances track momenta
269  if (nonPFCandTracksSumP <
270  (addNeutralsSumP4.energy() + mergedNeutralsSumP4.energy() + chargedHadronNeutralsSumP4.energy())) {
271  double scaleFactor =
272  ((addNeutralsSumP4.energy() + mergedNeutralsSumP4.energy() + chargedHadronNeutralsSumP4.energy()) -
273  nonPFCandTracksSumP) /
274  chargedHadronNeutralsSumP4.energy();
275  if (!(scaleFactor >= 0. && scaleFactor <= 1.)) {
276  edm::LogWarning("PFRecoTauEnergyAlgorithmPlugin::operator()")
277  << "Failed to compute tau energy --> killing tau candidate !!" << std::endl;
278  killTau(tau);
279  return;
280  }
281  reco::Candidate::LorentzVector diffP4;
282  size_t numChargedHadrons = chargedHadrons.size();
283  for (size_t iChargedHadron = 0; iChargedHadron < numChargedHadrons; ++iChargedHadron) {
284  const PFRecoTauChargedHadron& chargedHadron = chargedHadrons[iChargedHadron];
285  if (chargedHadron.algoIs(PFRecoTauChargedHadron::kPFNeutralHadron)) {
286  PFRecoTauChargedHadron chargedHadron_modified = chargedHadron;
287  chargedHadron_modified.neutralPFCandidates_.clear();
288  const CandidatePtr& chargedPFCand = chargedHadron.getChargedPFCandidate();
289  double chargedHadronPx_modified = scaleFactor * chargedPFCand->px();
290  double chargedHadronPy_modified = scaleFactor * chargedPFCand->py();
291  double chargedHadronPz_modified = scaleFactor * chargedPFCand->pz();
292  reco::Candidate::LorentzVector chargedHadronP4_modified = compChargedHadronP4fromPxPyPz(
293  chargedHadronPx_modified, chargedHadronPy_modified, chargedHadronPz_modified);
294  chargedHadron_modified.setP4(chargedHadronP4_modified);
295  tau.signalTauChargedHadronCandidatesRestricted()[iChargedHadron] = chargedHadron_modified;
296  diffP4 += (chargedHadron.p4() - chargedHadron_modified.p4());
297  }
298  }
299  if (verbosity_) {
300  std::cout << "case (4): (addNeutralsSumEn + mergedNeutralsSumEn + chargedHadronNeutralsSumEn) > "
301  "nonPFCandTracksSumP --> adjusting momenta of tau and chargedHadrons."
302  << std::endl;
303  }
304  updateTauP4(tau, -1., diffP4);
305  return;
306  } else {
307  double allTracksSumP = 0.;
308  double allTracksSumPerr2 = 0.;
309  const std::vector<PFRecoTauChargedHadron> chargedHadrons = tau.signalTauChargedHadronCandidatesRestricted();
310  for (std::vector<PFRecoTauChargedHadron>::const_iterator chargedHadron = chargedHadrons.begin();
311  chargedHadron != chargedHadrons.end();
312  ++chargedHadron) {
313  if (chargedHadron->algoIs(PFRecoTauChargedHadron::kChargedPFCandidate) ||
314  chargedHadron->algoIs(PFRecoTauChargedHadron::kTrack)) {
315  const reco::Track* chargedHadronTrack = getTrackFromChargedHadron(*chargedHadron);
316  if (chargedHadronTrack != nullptr) {
317  allTracksSumP += chargedHadronTrack->p();
318  allTracksSumPerr2 += getTrackPerr2(*chargedHadronTrack);
319  } else {
320  edm::LogInfo("PFRecoTauEnergyAlgorithmPlugin::operator()")
321  << "PFRecoTauChargedHadron has no associated reco::Track !!";
322  if (verbosity_) {
323  chargedHadron->print();
324  }
325  }
326  }
327  }
328  if (verbosity_) {
329  std::cout << "allTracksSumP = " << allTracksSumP << " +/- " << sqrt(allTracksSumPerr2) << std::endl;
330  }
331  double allNeutralsSumEn = 0.;
332  const auto& signalCands = tau.signalCands();
333  for (const auto& signalCand : signalCands) {
334  if (verbosity_) {
335  std::cout << "Candidate #" << signalCand.id() << ":" << signalCand.key() << ":"
336  << " Pt = " << (signalCand)->pt() << ", eta = " << (signalCand)->eta()
337  << ", phi = " << (signalCand)->phi() << std::endl;
338  }
339  const PFCandidate* pfCand = dynamic_cast<const PFCandidate*>(&*signalCand);
340  if (pfCand) {
341  if (verbosity_) {
342  std::cout << "calorimeter energy:"
343  << " ECAL = " << (pfCand)->ecalEnergy() << ","
344  << " HCAL = " << (pfCand)->hcalEnergy() << ","
345  << " HO = " << (pfCand)->hoEnergy() << std::endl;
346  }
347  // TauReco@MiniAOD: This info is not yet available in miniAOD.
348  if (edm::isFinite(pfCand->ecalEnergy()))
349  allNeutralsSumEn += pfCand->ecalEnergy();
350  if (edm::isFinite(pfCand->hcalEnergy()))
351  allNeutralsSumEn += pfCand->hcalEnergy();
352  if (edm::isFinite(pfCand->hoEnergy()))
353  allNeutralsSumEn += pfCand->hoEnergy();
354  }
355  }
356  allNeutralsSumEn += addNeutralsSumP4.energy();
357  if (allNeutralsSumEn < 0.)
358  allNeutralsSumEn = 0.;
359  if (verbosity_) {
360  std::cout << "allNeutralsSumEn = " << allNeutralsSumEn << std::endl;
361  }
362  if (allNeutralsSumEn > allTracksSumP) {
363  // Adjust momenta of neutral PFCandidates merged into ChargedHadrons
364  size_t numChargedHadrons = chargedHadrons.size();
365  for (size_t iChargedHadron = 0; iChargedHadron < numChargedHadrons; ++iChargedHadron) {
366  const PFRecoTauChargedHadron& chargedHadron = chargedHadrons[iChargedHadron];
367  if (chargedHadron.algoIs(PFRecoTauChargedHadron::kChargedPFCandidate)) {
368  PFRecoTauChargedHadron chargedHadron_modified = chargedHadron;
369  chargedHadron_modified.neutralPFCandidates_.clear();
370  chargedHadron_modified.setP4(chargedHadron.getChargedPFCandidate()->p4());
371  if (verbosity_) {
372  std::cout << "chargedHadron #" << iChargedHadron << ": changing En = " << chargedHadron.energy()
373  << " to " << chargedHadron_modified.energy() << std::endl;
374  }
375  tau.signalTauChargedHadronCandidatesRestricted()[iChargedHadron] = chargedHadron_modified;
376  } else if (chargedHadron.algoIs(PFRecoTauChargedHadron::kTrack)) {
377  PFRecoTauChargedHadron chargedHadron_modified = chargedHadron;
378  chargedHadron_modified.neutralPFCandidates_.clear();
379  reco::Candidate::LorentzVector chargedHadronP4_modified(0., 0., 0., 0.);
380  const reco::Track* chTrack = getTrackFromChargedHadron(chargedHadron);
381  if (chTrack != nullptr) {
382  double chargedHadronPx_modified = chTrack->px();
383  double chargedHadronPy_modified = chTrack->py();
384  double chargedHadronPz_modified = chTrack->pz();
385  chargedHadronP4_modified = compChargedHadronP4fromPxPyPz(
386  chargedHadronPx_modified, chargedHadronPy_modified, chargedHadronPz_modified);
387  } else {
388  edm::LogWarning("PFRecoTauEnergyAlgorithmPlugin::operator()")
389  << "PFRecoTauChargedHadron has no associated reco::Track !!" << std::endl;
390  if (verbosity_) {
391  chargedHadron.print();
392  }
393  }
394  chargedHadron_modified.setP4(chargedHadronP4_modified);
395  if (verbosity_) {
396  std::cout << "chargedHadron #" << iChargedHadron << ": changing En = " << chargedHadron.energy()
397  << " to " << chargedHadron_modified.energy() << std::endl;
398  }
399  tau.signalTauChargedHadronCandidatesRestricted()[iChargedHadron] = chargedHadron_modified;
400  }
401  }
402  double scaleFactor = allNeutralsSumEn / tau.energy();
403  if (verbosity_) {
404  std::cout << "case (5): allNeutralsSumEn > allTracksSumP --> adjusting momenta of tau and chargedHadrons."
405  << std::endl;
406  }
407  updateTauP4(tau, scaleFactor - 1., tau.p4());
408  return;
409  } else {
410  if (numChargedHadronNeutrals == 0 && tau.signalPiZeroCandidates().empty()) {
411  // Adjust momenta of ChargedHadrons build from reco::Tracks to match sum of energy deposits in ECAL + HCAL + HO
412  size_t numChargedHadrons = chargedHadrons.size();
413  for (size_t iChargedHadron = 0; iChargedHadron < numChargedHadrons; ++iChargedHadron) {
414  const PFRecoTauChargedHadron& chargedHadron = chargedHadrons[iChargedHadron];
415  if (chargedHadron.algoIs(PFRecoTauChargedHadron::kChargedPFCandidate) ||
416  chargedHadron.algoIs(PFRecoTauChargedHadron::kTrack)) {
417  PFRecoTauChargedHadron chargedHadron_modified = chargedHadron;
418  chargedHadron_modified.neutralPFCandidates_.clear();
419  reco::Candidate::LorentzVector chargedHadronP4_modified(0., 0., 0., 0.);
420  const reco::Track* chargedHadronTrack = getTrackFromChargedHadron(chargedHadron);
421  if (chargedHadronTrack != nullptr) {
422  double trackP = chargedHadronTrack->p();
423  double trackPerr2 = getTrackPerr2(*chargedHadronTrack);
424  if (verbosity_) {
425  std::cout << "trackP = " << trackP << " +/- " << sqrt(trackPerr2) << std::endl;
426  }
427  // CV: adjust track momenta such that difference beeen (measuredTrackP - adjustedTrackP)/sigmaMeasuredTrackP is minimal
428  // (expression derived using Mathematica)
429  double trackP_modified = (trackP * (allTracksSumPerr2 - trackPerr2) +
430  trackPerr2 * (allNeutralsSumEn - (allTracksSumP - trackP))) /
431  allTracksSumPerr2;
432  // CV: trackP_modified may actually become negative in case sum of energy deposits in ECAL + HCAL + HO is small
433  // and one of the tracks has a significantly larger momentum uncertainty than the other tracks.
434  // In this case set track momentum to small positive value.
435  if (trackP_modified < 1.e-1)
436  trackP_modified = 1.e-1;
437  if (verbosity_) {
438  std::cout << "trackP (modified) = " << trackP_modified << std::endl;
439  }
440  double scaleFactor = trackP_modified / trackP;
441  if (!(scaleFactor >= 0. && scaleFactor <= 1.)) {
442  edm::LogWarning("PFRecoTauEnergyAlgorithmPlugin::operator()")
443  << "Failed to compute tau energy --> killing tau candidate !!" << std::endl;
444  killTau(tau);
445  return;
446  }
447  double chargedHadronPx_modified = scaleFactor * chargedHadronTrack->px();
448  double chargedHadronPy_modified = scaleFactor * chargedHadronTrack->py();
449  double chargedHadronPz_modified = scaleFactor * chargedHadronTrack->pz();
450  chargedHadronP4_modified = compChargedHadronP4fromPxPyPz(
451  chargedHadronPx_modified, chargedHadronPy_modified, chargedHadronPz_modified);
452  } else {
453  edm::LogInfo("PFRecoTauEnergyAlgorithmPlugin::operator()")
454  << "PFRecoTauChargedHadron has no associated reco::Track !!";
455  if (verbosity_) {
456  chargedHadron.print();
457  }
458  }
459  chargedHadron_modified.setP4(chargedHadronP4_modified);
460  if (verbosity_) {
461  std::cout << "chargedHadron #" << iChargedHadron << ": changing En = " << chargedHadron.energy()
462  << " to " << chargedHadron_modified.energy() << std::endl;
463  }
464  tau.signalTauChargedHadronCandidatesRestricted()[iChargedHadron] = chargedHadron_modified;
465  }
466  }
467  double scaleFactor = allNeutralsSumEn / tau.energy();
468  if (verbosity_) {
469  std::cout
470  << "case (6): allNeutralsSumEn < allTracksSumP --> adjusting momenta of tau and chargedHadrons."
471  << std::endl;
472  }
473  updateTauP4(tau, scaleFactor - 1., tau.p4());
474  return;
475  } else {
476  // Interpretation of PFNeutralHadrons as ChargedHadrons with missing track and/or reconstruction of extra PiZeros
477  // is not compatible with the fact that sum of reco::Track momenta exceeds sum of energy deposits in ECAL + HCAL + HO:
478  // kill tau candidate (by setting its four-vector to zero)
479  if (verbosity_) {
480  std::cout << "case (7): allNeutralsSumEn < allTracksSumP not compatible with tau decay mode hypothesis "
481  "--> killing tau candidate."
482  << std::endl;
483  }
484  killTau(tau);
485  return;
486  }
487  }
488  }
489  }
490 
491  // CV: You should never come here.
492  if (verbosity_) {
493  std::cout << "undefined case: you should never come here !!" << std::endl;
494  }
495  assert(0);
496  }
497 
498  void PFRecoTauEnergyAlgorithmPlugin::endEvent() {}
499 
500  } // namespace tau
501 } // namespace reco
502 
503 #include "FWCore/Framework/interface/MakerMacros.h"
504 
505 DEFINE_EDM_PLUGIN(RecoTauModifierPluginFactory,
506  reco::tau::PFRecoTauEnergyAlgorithmPlugin,
507  "PFRecoTauEnergyAlgorithmPlugin");
reco::tau::PFRecoTauEnergyAlgorithmPlugin::endEvent
void endEvent() override
Definition: PFRecoTauEnergyAlgorithmPlugin.cc:498
MessageLogger.h
reco::tau::PFRecoTauEnergyAlgorithmPlugin::beginEvent
void beginEvent() override
Definition: PFRecoTauEnergyAlgorithmPlugin.cc:66
metsig::tau
Definition: SignAlgoResolutions.h:49
PFCandidate.h
DiDispStaMuonMonitor_cfi.pt
pt
Definition: DiDispStaMuonMonitor_cfi.py:39
reco::TrackBase::p
double p() const
momentum vector magnitude
Definition: TrackBase.h:605
gather_cfg.cout
cout
Definition: gather_cfg.py:144
PFJet.h
edm::LogInfo
Definition: MessageLogger.h:254
reco::tau::isPtrEqual
bool isPtrEqual(const edm::Ptr< Base > &b, const edm::Ptr< Der > &d)
Definition: PFRecoTauEnergyAlgorithmPlugin.cc:94
cms::cuda::assert
assert(be >=bs)
reco::PFTau
Definition: PFTau.h:36
reco::TrackBase::px
double px() const
x coordinate of momentum vector
Definition: TrackBase.h:611
reco::tau::getTrackFromChargedHadron
const reco::Track * getTrackFromChargedHadron(const reco::PFRecoTauChargedHadron &chargedHadron)
Definition: pfRecoTauChargedHadronAuxFunctions.cc:12
reco
fixed size matrix
Definition: AlignmentAlgorithmBase.h:45
reco::tau::PFRecoTauEnergyAlgorithmPlugin::dRaddNeutralHadron_
double dRaddNeutralHadron_
Definition: PFRecoTauEnergyAlgorithmPlugin.cc:47
reco::tau::PFRecoTauEnergyAlgorithmPlugin::~PFRecoTauEnergyAlgorithmPlugin
~PFRecoTauEnergyAlgorithmPlugin() override
Definition: PFRecoTauEnergyAlgorithmPlugin.cc:64
HPSPFTauProducerPuppi_cfi.chargedHadron
chargedHadron
Definition: HPSPFTauProducerPuppi_cfi.py:7
reco::PFRecoTauChargedHadron
Definition: PFRecoTauChargedHadron.h:23
deltaR.h
MakerMacros.h
reco::tau::setChargedHadronP4
void setChargedHadronP4(reco::PFRecoTauChargedHadron &chargedHadron, double scaleFactor_neutralPFCands=1.0)
Definition: pfRecoTauChargedHadronAuxFunctions.cc:31
Track.h
reco::TrackBase::py
double py() const
y coordinate of momentum vector
Definition: TrackBase.h:614
PVValHelper::eta
Definition: PVValidationHelpers.h:69
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
reco::Track
Definition: Track.h:27
reco::tau::PFRecoTauEnergyAlgorithmPlugin::PFRecoTauEnergyAlgorithmPlugin
PFRecoTauEnergyAlgorithmPlugin(const edm::ParameterSet &, edm::ConsumesCollector &&iC)
Definition: PFRecoTauEnergyAlgorithmPlugin.cc:55
DEFINE_EDM_PLUGIN
#define DEFINE_EDM_PLUGIN(factory, type, name)
Definition: PluginFactory.h:124
b
double b
Definition: hdecay.h:118
edm::LogWarning
Definition: MessageLogger.h:141
edm::ParameterSet
Definition: ParameterSet.h:36
reco::tau::compChargedHadronP4fromPxPyPz
reco::Candidate::LorentzVector compChargedHadronP4fromPxPyPz(double, double, double)
Definition: pfRecoTauChargedHadronAuxFunctions.cc:80
edmplugin::PluginFactory
Definition: PluginFactory.h:34
PackedCandidate.h
reco::tau::PFRecoTauEnergyAlgorithmPlugin::minGammaEt_
double minGammaEt_
Definition: PFRecoTauEnergyAlgorithmPlugin.cc:50
reco::PFRecoTauChargedHadron::kChargedPFCandidate
Definition: PFRecoTauChargedHadron.h:30
Ptr.h
createfilelist.int
int
Definition: createfilelist.py:10
RefToPtr.h
reco::LeafCandidate::p4
const LorentzVector & p4() const final
four-momentum Lorentz vector
Definition: LeafCandidate.h:114
reco::tau::PFRecoTauEnergyAlgorithmPlugin::operator()
void operator()(reco::PFTau &) const override
Definition: PFRecoTauEnergyAlgorithmPlugin.cc:103
reco::tau::PFRecoTauEnergyAlgorithmPlugin
Definition: PFRecoTauEnergyAlgorithmPlugin.cc:38
reco::PFCandidate::ecalEnergy
double ecalEnergy() const
return corrected Ecal energy
Definition: PFCandidate.h:220
edm::Ptr
Definition: AssociationVector.h:31
looper.cfg
cfg
Definition: looper.py:297
RecoTauBuilderPlugins.h
reco::PFRecoTauChargedHadron::neutralPFCandidates_
std::vector< CandidatePtr > neutralPFCandidates_
Definition: PFRecoTauChargedHadron.h:88
reco::LeafCandidate::setP4
void setP4(const LorentzVector &p4) final
set 4-momentum
Definition: LeafCandidate.h:158
eostools.move
def move(src, dest)
Definition: eostools.py:511
std
Definition: JetResolutionObject.h:76
reco::PFRecoTauChargedHadron::kPFNeutralHadron
Definition: PFRecoTauChargedHadron.h:32
isFinite.h
reco::tau::RecoTauModifierPlugin
Definition: RecoTauBuilderPlugins.h:104
reco::PFCandidate::hcalEnergy
double hcalEnergy() const
return corrected Hcal energy
Definition: PFCandidate.h:232
reco::tau::PFRecoTauEnergyAlgorithmPlugin::dRaddPhoton_
double dRaddPhoton_
Definition: PFRecoTauEnergyAlgorithmPlugin.cc:49
multiplicitycorr_cfi.scaleFactor
scaleFactor
Definition: multiplicitycorr_cfi.py:7
reco::deltaR
constexpr auto deltaR(const T1 &t1, const T2 &t2) -> decltype(t1.eta())
Definition: deltaR.h:30
reco::LeafCandidate::energy
double energy() const final
energy
Definition: LeafCandidate.h:125
pfRecoTauChargedHadronAuxFunctions.h
PFRecoTauChargedHadron.h
reco::PFCandidate
Particle reconstructed by the particle flow algorithm.
Definition: PFCandidate.h:40
HLT_2018_cff.track
track
Definition: HLT_2018_cff.py:10352
ztail.d
d
Definition: ztail.py:151
reco::TrackBase::pz
double pz() const
z coordinate of momentum vector
Definition: TrackBase.h:617
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
edm::isFinite
constexpr bool isFinite(T x)
HGC3DClusterGenMatchSelector_cfi.dR
dR
Definition: HGC3DClusterGenMatchSelector_cfi.py:7
reco::tau::PFRecoTauEnergyAlgorithmPlugin::minNeutralHadronEt_
double minNeutralHadronEt_
Definition: PFRecoTauEnergyAlgorithmPlugin.cc:48
reco::Candidate::LorentzVector
math::XYZTLorentzVector LorentzVector
Lorentz vector.
Definition: Candidate.h:36
reco::tau::PFRecoTauEnergyAlgorithmPlugin::verbosity_
int verbosity_
Definition: PFRecoTauEnergyAlgorithmPlugin.cc:52
edm::ConsumesCollector
Definition: ConsumesCollector.h:39
PFCandidateFwd.h
PFRecoTauChargedHadronFwd.h
reco::PFCandidate::hoEnergy
double hoEnergy() const
return corrected Hcal energy
Definition: PFCandidate.h:244
MillePedeFileConverter_cfg.e
e
Definition: MillePedeFileConverter_cfg.py:37
reco::PFRecoTauChargedHadron::kTrack
Definition: PFRecoTauChargedHadron.h:31
Generated for CMSSW Reference Manual by doxygen 1.8.16