CMS 3D CMS Logo

List of all members | Public Member Functions | Protected Member Functions | Private Attributes
PFMEtSignInterfaceBase Class Reference

#include <PFMEtSignInterfaceBase.h>

Public Member Functions

template<typename T >
metsig::SigInputObj compResolution (const T *particle) const
 
template<typename T >
std::vector< metsig::SigInputObjcompResolution (const std::list< T * > &particles) const
 
reco::METCovMatrix operator() (const std::vector< metsig::SigInputObj > &) const
 
template<typename T >
reco::METCovMatrix operator() (const std::list< T * > &particles) const
 
 PFMEtSignInterfaceBase (const edm::ParameterSet &)
 
 ~PFMEtSignInterfaceBase ()
 

Protected Member Functions

template<typename T >
void addPFMEtSignObjects (std::vector< metsig::SigInputObj > &metSignObjects, const std::list< T * > &particles) const
 

Private Attributes

TFile * inputFile_
 
TH2 * lut_
 
metsig::SignAlgoResolutionspfMEtResolution_
 
int verbosity_
 

Detailed Description

Auxiliary class interfacing the TauAnalysis software to RecoMET/METAlgorithms/interface/significanceAlgo.h for computing (PF)MEt significance (see CMS AN-10/400 for description of the (PF)MEt significance computation)

Author
Christian Veelken, LLR

Definition at line 41 of file PFMEtSignInterfaceBase.h.

Constructor & Destructor Documentation

PFMEtSignInterfaceBase::PFMEtSignInterfaceBase ( const edm::ParameterSet cfg)

Definition at line 15 of file PFMEtSignInterfaceBase.cc.

References Exception, edm::ParameterSet::exists(), edm::ParameterSet::getParameter(), inputFile_, writeAntiElectronDiscrMVA_cfi::inputFileName, edm::FileInPath::Local, lut_, pfMEtResolution_, AlCaHLTBitMon_QueryRunRegistry::string, and verbosity_.

16  : pfMEtResolution_(nullptr),
17  inputFile_(nullptr),
18  lut_(nullptr)
19 {
21 
22  if ( cfg.exists("addJERcorr") ) {
23  edm::ParameterSet cfgJERcorr = cfg.getParameter<edm::ParameterSet>("addJERcorr");
24  edm::FileInPath inputFileName = cfgJERcorr.getParameter<edm::FileInPath>("inputFileName");
25  std::string lutName = cfgJERcorr.getParameter<std::string>("lutName");
26  if ( inputFileName.location()!=edm::FileInPath::Local )
27  throw cms::Exception("PFMEtSignInterfaceBase")
28  << " Failed to find File = " << inputFileName << " !!\n";
29 
30  inputFile_ = new TFile(inputFileName.fullPath().data());
31  lut_ = dynamic_cast<TH2*>(inputFile_->Get(lutName.data()));
32  if ( !lut_ )
33  throw cms::Exception("PFMEtSignInterfaceBase")
34  << " Failed to load LUT = " << lutName.data() << " from file = " << inputFileName.fullPath().data() << " !!\n";
35  }
36 
37  verbosity_ = cfg.exists("verbosity") ?
38  cfg.getParameter<int>("verbosity") : 0;
39 }
T getParameter(std::string const &) const
metsig::SignAlgoResolutions * pfMEtResolution_
bool exists(std::string const &parameterName) const
checks if a parameter exists
PFMEtSignInterfaceBase::~PFMEtSignInterfaceBase ( )

Definition at line 41 of file PFMEtSignInterfaceBase.cc.

References inputFile_, lut_, and pfMEtResolution_.

42 {
43  delete pfMEtResolution_;
44  delete inputFile_;
45  delete lut_;
46 }
metsig::SignAlgoResolutions * pfMEtResolution_

Member Function Documentation

template<typename T >
void PFMEtSignInterfaceBase::addPFMEtSignObjects ( std::vector< metsig::SigInputObj > &  metSignObjects,
const std::list< T * > &  particles 
) const
inlineprotected

Definition at line 178 of file PFMEtSignInterfaceBase.h.

References compResolution().

Referenced by compResolution().

180  {
181  for ( typename std::list<T*>::const_iterator particle = particles.begin();
182  particle != particles.end(); ++particle ) {
183  metSignObjects.push_back(this->compResolution(*particle));
184  }
185  }
metsig::SigInputObj compResolution(const T *particle) const
template<typename T >
metsig::SigInputObj PFMEtSignInterfaceBase::compResolution ( const T particle) const
inline

Definition at line 49 of file PFMEtSignInterfaceBase.h.

References funct::abs(), metsig::ET, PVValHelper::eta, metsig::SignAlgoResolutions::eval(), metsig::SignAlgoResolutions::evalPF(), metsig::SignAlgoResolutions::evalPFJet(), Exception, edm::Ref< C, T, F >::get(), metsig::SigInputObj::get_energy(), metsig::SigInputObj::get_phi(), metsig::SigInputObj::get_sigma_e(), metsig::SigInputObj::get_sigma_tan(), metsig::SigInputObj::get_type(), reco::Jet::getJetConstituents(), edm::Ref< C, T, F >::isAvailable(), edm::Ref< C, T, F >::isNonnull(), pat::Jet::isPFJet(), metsig::jet, reco::PFTau::jetRef(), lut_, reco::LeafCandidate::p4(), objects.autophobj::particleType, pat::Tau::pfJetRef(), pfMEtResolution_, pat::Jet::pfSpecific(), metsig::PFtype2, metsig::PFtype3, phi, metsig::PHI, reco::TrackBase::phiError(), EnergyCorrector::pt, reco::TrackBase::ptError(), metsig::SigInputObj::set(), AlCaHLTBitMon_QueryRunRegistry::string, reco::Muon::track(), pat::Muon::track(), reco::LeafCandidate::vertex(), x, and y.

Referenced by addPFMEtSignObjects(), operator()(), NoPileUpPFMEtProducer::produce(), and NoPileUpPFMEtDataProducer::produce().

50  {
51  double pt = particle->pt();
52  double eta = particle->eta();
53  double phi = particle->phi();
54 
55  if ( dynamic_cast<const reco::GsfElectron*>(particle) != 0 ||
56  dynamic_cast<const pat::Electron*>(particle) != 0 ) {
57  std::string particleType = "electron";
58  // WARNING: SignAlgoResolutions::PFtype2 needs to be kept in sync with reco::PFCandidate::e !!
59  double dpt = pfMEtResolution_->eval(metsig::PFtype2, metsig::ET, pt, phi, eta);
60  double dphi = pfMEtResolution_->eval(metsig::PFtype2, metsig::PHI, pt, phi, eta);
61  //std::cout << "electron: pt = " << pt << ", eta = " << eta << ", phi = " << phi
62  // << " --> dpt = " << dpt << ", dphi = " << dphi << std::endl;
63  return metsig::SigInputObj(particleType, pt, phi, dpt, dphi);
64  } else if ( dynamic_cast<const reco::Photon*>(particle) != 0 ||
65  dynamic_cast<const pat::Photon*>(particle) != 0 ) {
66  // CV: assume resolutions for photons to be identical to electron resolutions
67  std::string particleType = "electron";
68  // WARNING: SignAlgoResolutions::PFtype2 needs to be kept in sync with reco::PFCandidate::e !!
69  double dpt = pfMEtResolution_->eval(metsig::PFtype2, metsig::ET, pt, phi, eta);
70  double dphi = pfMEtResolution_->eval(metsig::PFtype2, metsig::PHI, pt, phi, eta);
71  //std::cout << "photon: pt = " << pt << ", eta = " << eta << ", phi = " << phi
72  // << " --> dpt = " << dpt << ", dphi = " << dphi << std::endl;
73  return metsig::SigInputObj(particleType, pt, phi, dpt, dphi);
74  } else if ( dynamic_cast<const reco::Muon*>(particle) != 0 ||
75  dynamic_cast<const pat::Muon*>(particle) != 0 ) {
76  std::string particleType = "muon";
77  double dpt, dphi;
78  const reco::Track* muonTrack = nullptr;
79  if ( dynamic_cast<const pat::Muon*>(particle) != 0 ) {
80  const pat::Muon* muon = dynamic_cast<const pat::Muon*>(particle);
81  if ( muon->track().isNonnull() && muon->track().isAvailable() ) muonTrack = muon->track().get();
82  } else if ( dynamic_cast<const reco::Muon*>(particle) != 0 ) {
83  const reco::Muon* muon = dynamic_cast<const reco::Muon*>(particle);
84  if ( muon->track().isNonnull() && muon->track().isAvailable() ) muonTrack = muon->track().get();
85  } else assert(0);
86  if ( muonTrack ) {
87  dpt = muonTrack->ptError();
88  dphi = pt*muonTrack->phiError(); // CV: pt*dphi is indeed correct
89  } else {
90  // WARNING: SignAlgoResolutions::PFtype3 needs to be kept in sync with reco::PFCandidate::mu !!
91  dpt = pfMEtResolution_->eval(metsig::PFtype3, metsig::ET, pt, phi, eta);
92  dphi = pfMEtResolution_->eval(metsig::PFtype3, metsig::PHI, pt, phi, eta);
93  }
94  //std::cout << "muon: pt = " << pt << ", eta = " << eta << ", phi = " << phi
95  // << " --> dpt = " << dpt << ", dphi = " << dphi << std::endl;
96  return metsig::SigInputObj(particleType, pt, phi, dpt, dphi);
97  } else if ( dynamic_cast<const reco::PFTau*>(particle) != 0 ||
98  dynamic_cast<const pat::Tau*>(particle) != 0 ) {
99  // CV: use PFJet resolutions for PFTaus for now...
100  // (until PFTau specific resolutions are available)
101  if ( dynamic_cast<const pat::Tau*>(particle) != 0 ) {
102  const pat::Tau* pfTau = dynamic_cast<const pat::Tau*>(particle);
103  //std::cout << "tau: pt = " << pt << ", eta = " << eta << ", phi = " << phi << std::endl;
104  return pfMEtResolution_->evalPFJet(pfTau->pfJetRef().get());
105  } else if ( dynamic_cast<const reco::PFTau*>(particle) != 0 ) {
106  const reco::PFTau* pfTau = dynamic_cast<const reco::PFTau*>(particle);
107  //std::cout << "tau: pt = " << pt << ", eta = " << eta << ", phi = " << phi << std::endl;
108  return pfMEtResolution_->evalPFJet(pfTau->jetRef().get());
109  } else assert(0);
110  } else if ( dynamic_cast<const reco::PFJet*>(particle) != 0 ||
111  dynamic_cast<const pat::Jet*>(particle) != 0 ) {
112  metsig::SigInputObj pfJetResolution;
113  if ( dynamic_cast<const reco::PFJet*>(particle) != 0 ) {
114  const reco::PFJet* pfJet = dynamic_cast<const reco::PFJet*>(particle);
115  pfJetResolution = pfMEtResolution_->evalPFJet(pfJet);
116  } else if ( dynamic_cast<const pat::Jet*>(particle) != 0 ) {
117  const pat::Jet* jet = dynamic_cast<const pat::Jet*>(particle);
118  if ( jet->isPFJet() ) {
119  reco::PFJet pfJet(jet->p4(), jet->vertex(), jet->pfSpecific(), jet->getJetConstituents());
120  pfJetResolution = pfMEtResolution_->evalPFJet(&pfJet);
121  } else throw cms::Exception("addPFMEtSignObjects")
122  << "PAT jet not of PF-type !!\n";
123  } else assert(0);
124  //std::cout << "pfJet: pt = " << pt << ", eta = " << eta << ", phi = " << phi << std::endl;
125  // CV: apply additional jet energy resolution corrections
126  // not included in (PF)MEt significance algorithm yet
127  // (cf. CMS AN-11/400 vs. CMS AN-11/330)
128  if ( lut_ && pt > 10. ) {
129  double x = std::abs(eta);
130  double y = pt;
131  if ( x > lut_->GetXaxis()->GetXmin() && x < lut_->GetXaxis()->GetXmax() &&
132  y > lut_->GetYaxis()->GetXmin() && y < lut_->GetYaxis()->GetXmax() ) {
133  int binIndex = lut_->FindBin(x, y);
134  double addJERcorrFactor = lut_->GetBinContent(binIndex);
135  //std::cout << " addJERcorrFactor = " << addJERcorrFactor << std::endl;
136  pfJetResolution.set(pfJetResolution.get_type(),
137  pfJetResolution.get_energy(),
138  pfJetResolution.get_phi(),
139  addJERcorrFactor*pfJetResolution.get_sigma_e(),
140  pfJetResolution.get_sigma_tan());
141  }
142  }
143  return pfJetResolution;
144  } else if ( dynamic_cast<const reco::PFCandidate*>(particle) != 0 ) {
145  const reco::PFCandidate* pfCandidate = dynamic_cast<const reco::PFCandidate*>(particle);
146  //std::cout << "pfCandidate: pt = " << pt << ", eta = " << eta << ", phi = " << phi << std::endl;
147  return pfMEtResolution_->evalPF(pfCandidate);
148  } else throw cms::Exception("addPFMEtSignObjects")
149  << "Invalid type of particle:"
150  << " valid types = { reco::GsfElectron/pat::Electron, reco::Photon/pat::Photon, reco::Muon/pat::Muon, reco::PFTau/pat::Tau,"
151  << " reco::PFJet/pat::Jet, reco::PFCandidate } !!\n";
152  }
bool isAvailable() const
Definition: Ref.h:577
double get_sigma_tan() const
Definition: SigInputObj.h:46
metsig::SignAlgoResolutions * pfMEtResolution_
bool isNonnull() const
Checks for non-null.
Definition: Ref.h:253
const PFJetRef & jetRef() const
Definition: PFTau.cc:58
reco::TrackRef track() const
reference to Track reconstructed in the tracker only (reimplemented from reco::Muon) ...
double get_phi() const
Definition: SigInputObj.h:44
const PFSpecific & pfSpecific() const
retrieve the pf specific part of the jet
Definition: Jet.h:267
virtual const Point & vertex() const
vertex position (overwritten by PF...)
metsig::SigInputObj evalPF(const reco::PFCandidate *candidate) const
virtual TrackRef track() const
reference to a Track
Definition: Muon.h:49
virtual Constituents getJetConstituents() const
list of constituents
double eval(const resolutionType &type, const resolutionFunc &func, const double &et, const double &phi, const double &eta, const double &p) const
Jets made from PFObjects.
Definition: PFJet.h:21
void set(const std::string &m_type, const double &m_energy, const double &m_phi, const double &m_sigma_e, const double &m_sigma_tan)
Definition: SigInputObj.h:48
double ptError() const
error on Pt (set to 1000 TeV if charge==0 for safety)
Definition: TrackBase.h:763
double phiError() const
error on phi
Definition: TrackBase.h:790
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
T const * get() const
Returns C++ pointer to the item.
Definition: Ref.h:245
Analysis-level tau class.
Definition: Tau.h:55
bool isPFJet() const
check to see if the jet is a reco::PFJet
Definition: Jet.h:253
const reco::PFJetRef & pfJetRef() const
Definition: Tau.h:163
double get_energy() const
Definition: SigInputObj.h:43
Analysis-level calorimeter jet class.
Definition: Jet.h:78
metsig::SigInputObj evalPFJet(const reco::PFJet *jet) const
Particle reconstructed by the particle flow algorithm.
Definition: PFCandidate.h:39
virtual const LorentzVector & p4() const final
four-momentum Lorentz vector
Definition: LeafCandidate.h:99
std::string get_type() const
Definition: SigInputObj.h:42
Analysis-level muon class.
Definition: Muon.h:49
double get_sigma_e() const
Definition: SigInputObj.h:45
template<typename T >
std::vector<metsig::SigInputObj> PFMEtSignInterfaceBase::compResolution ( const std::list< T * > &  particles) const
inline

Definition at line 155 of file PFMEtSignInterfaceBase.h.

References addPFMEtSignObjects(), LogDebug, and operator()().

156  {
157  LogDebug("compResolution")
158  << " particles: entries = " << particles.size() << std::endl;
159 
160  std::vector<metsig::SigInputObj> pfMEtSignObjects;
161  addPFMEtSignObjects(pfMEtSignObjects, particles);
162 
163  return pfMEtSignObjects;
164  }
#define LogDebug(id)
void addPFMEtSignObjects(std::vector< metsig::SigInputObj > &metSignObjects, const std::list< T * > &particles) const
reco::METCovMatrix PFMEtSignInterfaceBase::operator() ( const std::vector< metsig::SigInputObj > &  pfMEtSignObjects) const

Definition at line 48 of file PFMEtSignInterfaceBase.cc.

References funct::abs(), metsig::significanceAlgo::addObjects(), defaultPFMEtResolutionX, defaultPFMEtResolutionY, epsilon, and metsig::significanceAlgo::getSignifMatrix().

Referenced by compResolution(), and operator()().

49 {
50  // if ( this->verbosity_ ) {
51  // std::cout << "<PFMEtSignInterfaceBase::operator()>:" << std::endl;
52  // std::cout << " pfMEtSignObjects: entries = " << pfMEtSignObjects.size() << std::endl;
53  // double dpt2Sum = 0.;
54  // for ( std::vector<metsig::SigInputObj>::const_iterator pfMEtSignObject = pfMEtSignObjects.begin();
55  // pfMEtSignObject != pfMEtSignObjects.end(); ++pfMEtSignObject ) {
56  // std::cout << pfMEtSignObject->get_type() << ": pt = " << pfMEtSignObject->get_energy() << ","
57  // << " phi = " << pfMEtSignObject->get_phi() << " --> dpt = " << pfMEtSignObject->get_sigma_e() << std::endl;
58  // dpt2Sum += pfMEtSignObject->get_sigma_e();
59  // }
60  // std::cout << "--> sqrt(sum(dpt^2)) = " << sqrt(dpt2Sum) << std::endl;
61  // }
62 
63  reco::METCovMatrix pfMEtCov;
64  if ( pfMEtSignObjects.size() >= 2 ) {
65  metsig::significanceAlgo pfMEtSignAlgorithm;
66  pfMEtSignAlgorithm.addObjects(pfMEtSignObjects);
67  pfMEtCov = pfMEtSignAlgorithm.getSignifMatrix();
68 
69  double det=0;
70  pfMEtCov.Det(det);
71 
72  // if ( this->verbosity_ && std::abs(det) > epsilon ) {
73  // //keep TMatrixD as it is much easier to find
74  // //eigenvectors and values than with SMatrix;
75  // //not used anyway, except for debugging
76  // TMatrixD tmpMatrix(2,2);
77  // tmpMatrix(0,0) = pfMEtCov(0,0);
78  // tmpMatrix(0,1) = pfMEtCov(0,1);
79  // tmpMatrix(1,0) = pfMEtCov(1,0);
80  // tmpMatrix(1,1) = pfMEtCov(1,1);
81 
82  // TVectorD eigenValues(2);
83  // TMatrixD eigenVectors = tmpMatrix.EigenVectors(eigenValues);
84  // // CV: eigenvectors are stored in columns
85  // // and are sorted such that the one corresponding to the highest eigenvalue is in the **first** column
86  // for ( unsigned iEigenVector = 0; iEigenVector < 2; ++iEigenVector ) {
87  // std::cout << "eigenVector #" << iEigenVector << " (eigenValue = " << eigenValues(iEigenVector) << "):"
88  // << " x = " << eigenVectors(0, iEigenVector) << ", y = " << eigenVectors(1, iEigenVector) << std::endl;
89  // }
90  // }
91 
92  //--- substitute (PF)MEt resolution matrix by default values
93  // in case resolution matrix cannot be inverted
94  if (std::abs(det) < epsilon ) {
95  edm::LogWarning("PFMEtSignInterfaceBase::operator()")
96  << "Inversion of PFMEt covariance matrix failed, det = " << det
97  << " --> replacing covariance matrix by resolution defaults !!";
99  pfMEtCov(0,1) = 0.;
100  pfMEtCov(1,0) = 0.;
102  }
103  } else {
104  pfMEtCov(0,0) = 0.;
105  pfMEtCov(0,1) = 0.;
106  pfMEtCov(1,0) = 0.;
107  pfMEtCov(1,1) = 0.;
108  }
109 
110  return pfMEtCov;
111 }
const void addObjects(const std::vector< metsig::SigInputObj > &EventVec)
const double defaultPFMEtResolutionX
ROOT::Math::SMatrix< double, 2 > METCovMatrix
Definition: MET.h:40
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
const double defaultPFMEtResolutionY
reco::METCovMatrix getSignifMatrix() const
const double epsilon
template<typename T >
reco::METCovMatrix PFMEtSignInterfaceBase::operator() ( const std::list< T * > &  particles) const
inline

Definition at line 169 of file PFMEtSignInterfaceBase.h.

References compResolution(), and operator()().

170  {
171  std::vector<metsig::SigInputObj> pfMEtSignObjects = compResolution(particles);
172  return this->operator()(pfMEtSignObjects);
173  }
reco::METCovMatrix operator()(const std::vector< metsig::SigInputObj > &) const
metsig::SigInputObj compResolution(const T *particle) const

Member Data Documentation

TFile* PFMEtSignInterfaceBase::inputFile_
private

Definition at line 194 of file PFMEtSignInterfaceBase.h.

Referenced by PFMEtSignInterfaceBase(), and ~PFMEtSignInterfaceBase().

TH2* PFMEtSignInterfaceBase::lut_
private
metsig::SignAlgoResolutions* PFMEtSignInterfaceBase::pfMEtResolution_
private
int PFMEtSignInterfaceBase::verbosity_
private

Definition at line 197 of file PFMEtSignInterfaceBase.h.

Referenced by PFMEtSignInterfaceBase().