CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
PFRecoTauMassPlugin.cc
Go to the documentation of this file.
1 /*
2  * =============================================================================
3  * Filename: PFRecoTauMassPlugin.cc
4  *
5  * Description: Set mass of taus reconstructed in 1prong0pi0 decay mode
6  * to charged pion mass
7  *
8  * Created: 27/10/2015 10:30:00
9  *
10  * Authors: Christian Veelken (Tallinn)
11  *
12  * =============================================================================
13  */
14 
16 
19 
20 #include <TMath.h>
21 
22 namespace reco {
23  namespace tau {
24 
26  public:
28  ~PFRecoTauMassPlugin() override;
29  void operator()(PFTau&) const override;
30  void beginEvent() override;
31  void endEvent() override;
32 
33  private:
35  };
36 
39  verbosity_ = cfg.getParameter<int>("verbosity");
40  }
41 
43 
45 
47  if (verbosity_) {
48  std::cout << "<PFRecoTauMassPlugin::operator()>:" << std::endl;
49  std::cout << "tau: Pt = " << tau.pt() << ", eta = " << tau.eta() << ", phi = " << tau.phi()
50  << ", mass = " << tau.mass() << " (decayMode = " << tau.decayMode() << ")" << std::endl;
51  }
52 
53  if (tau.decayMode() == reco::PFTau::kOneProng0PiZero) {
54  double tauEn = tau.energy();
55  const double chargedPionMass = 0.13957; // GeV
56  if (tauEn < chargedPionMass)
57  tauEn = chargedPionMass;
58  double tauP_modified = TMath::Sqrt(tauEn * tauEn - chargedPionMass * chargedPionMass);
59  double tauPx_modified = TMath::Cos(tau.phi()) * TMath::Sin(tau.theta()) * tauP_modified;
60  double tauPy_modified = TMath::Sin(tau.phi()) * TMath::Sin(tau.theta()) * tauP_modified;
61  double tauPz_modified = TMath::Cos(tau.theta()) * tauP_modified;
62  reco::Candidate::LorentzVector tauP4_modified(tauPx_modified, tauPy_modified, tauPz_modified, tauEn);
63  if (verbosity_) {
64  std::cout << "--> setting tauP4: Pt = " << tauP4_modified.pt() << ", eta = " << tauP4_modified.eta()
65  << ", phi = " << tauP4_modified.phi() << ", mass = " << tauP4_modified.mass() << std::endl;
66  }
67  tau.setP4(tauP4_modified);
68  }
69  }
70 
72 
73  } // namespace tau
74 } // namespace reco
75 
77 
void operator()(PFTau &) const override
PFRecoTauMassPlugin(const edm::ParameterSet &, edm::ConsumesCollector &&iC)
math::XYZTLorentzVector LorentzVector
Lorentz vector.
Definition: Candidate.h:36
fixed size matrix
#define DEFINE_EDM_PLUGIN(factory, type, name)
def move(src, dest)
Definition: eostools.py:511