Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

JER MET uncertainties and METSignifiance re-estimation #16435

Merged
merged 3 commits into from Nov 17, 2016
Merged

JER MET uncertainties and METSignifiance re-estimation #16435

Show file tree
Hide file tree
Changes from 2 commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
b6bb6a0
enable update of MET sgnificance when correcting the MET
Oct 31, 2016
9d3044a
enable the JER uncertainty variation when MET is not smeared
Aug 26, 2016
27b0a91
cleaning
Nov 7, 2016
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
4 changes: 2 additions & 2 deletions PhysicsTools/PatAlgos/plugins/PATMETSlimmer.cc
View file
Open in desktop
Expand Up @@ -117,8 +117,8 @@ pat::PATMETSlimmer::OneMETShift::OneMETShift(pat::MET::METUncertainty shift_, pa
case pat::MET::NoShift : snprintf(buff, 1023, baseTagStr.c_str(), ""); break;
case pat::MET::JetEnUp : snprintf(buff, 1023, baseTagStr.c_str(), "JetEnUp"); break;
case pat::MET::JetEnDown: snprintf(buff, 1023, baseTagStr.c_str(), "JetEnDown"); break;
case pat::MET::JetResUp : snprintf(buff, 1023, baseTagStr.c_str(), isSmeared?"JetResUp":""); break;
case pat::MET::JetResDown: snprintf(buff, 1023, baseTagStr.c_str(), isSmeared?"JetResDown":""); break;
case pat::MET::JetResUp : snprintf(buff, 1023, baseTagStr.c_str(), "JetResUp"); break;
case pat::MET::JetResDown: snprintf(buff, 1023, baseTagStr.c_str(), "JetResDown"); break;
case pat::MET::MuonEnUp : snprintf(buff, 1023, baseTagStr.c_str(), "MuonEnUp"); break;
case pat::MET::MuonEnDown: snprintf(buff, 1023, baseTagStr.c_str(), "MuonEnDown"); break;
case pat::MET::ElectronEnUp : snprintf(buff, 1023, baseTagStr.c_str(), "ElectronEnUp"); break;
Expand Down
17 changes: 14 additions & 3 deletions PhysicsTools/PatUtils/interface/SmearedJetProducerT.h
View file
Open in desktop
Expand Up @@ -131,14 +131,23 @@ class SmearedJetProducerT : public edm::stream::EDProducer<> {
m_genJetMatcher = std::make_shared<pat::GenJetMatcher>(cfg, consumesCollector());

std::int32_t variation = cfg.getParameter<std::int32_t>("variation");
_nomVar=1;
if (variation == 0)
m_systematic_variation = Variation::NOMINAL;
else if (variation == 1)
m_systematic_variation = Variation::UP;
else if (variation == -1)
m_systematic_variation = Variation::DOWN;
else if (variation == 101) {
m_systematic_variation = Variation::NOMINAL;
_nomVar=1;
}
else if (variation == -101) {
m_systematic_variation = Variation::NOMINAL;
_nomVar=-1;
}
else
throw edm::Exception(edm::errors::ConfigFileReadError, "Invalid value for 'variation' parameter. Only -1, 0 or 1 are supported.");
throw edm::Exception(edm::errors::ConfigFileReadError, "Invalid value for 'variation' parameter. Only -1, 0, 1 or 101, -101 are supported.");
}

produces<JetCollection>();
Expand Down Expand Up @@ -234,7 +243,7 @@ class SmearedJetProducerT : public edm::stream::EDProducer<> {
}

double dPt = jet.pt() - genJet->pt();
smearFactor = 1 + (jer_sf - 1.) * dPt / jet.pt();
smearFactor = 1 + _nomVar*(jer_sf - 1.) * dPt / jet.pt();
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'm probably missing something.
What is the purpose of this "101" shift?
For +/-101 we get the same NOMINAL jet_sf value and you will smear by the full +/- |jet_sf -1|
while the systematics on the JER can be smaller or larger

Looking at
https://twiki.cern.ch/twiki/bin/view/CMS/JetResolution#JER_Scaling_factors_and_Uncertai
for eta 0 to 0.5 JER SF is 1.122 +-0.026 or "jer_sf - 1 = 0.122 +-0.026".

So, smearing by +-0.122 in this case is not meaningful, but maybe there is some other use case outside the nominal definition of JER.
Please clarify.

Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@slava77 Here the total semaring factor is basically smear by 1-0.122 or 1+0.122. This is the same logic as the previous code, but the scale factor is used instead of the uncertainties associated to the scale factor


} else if (jer_sf > 1) {
/*
Expand All @@ -247,7 +256,7 @@ class SmearedJetProducerT : public edm::stream::EDProducer<> {
}

std::normal_distribution<> d(0, sigma);
smearFactor = 1. + d(m_random_generator);
smearFactor = 1. + _nomVar*d(m_random_generator);
} else if (m_debug) {
std::cout << "Impossible to smear this jet" << std::endl;
}
Expand Down Expand Up @@ -298,5 +307,7 @@ class SmearedJetProducerT : public edm::stream::EDProducer<> {
std::mt19937 m_random_generator;

GreaterByPt<T> jetPtComparator;

int _nomVar;
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

leading underscore is frowned upon.

Please change to m_nomVar to be more in style with all other member data

};
#endif
10 changes: 9 additions & 1 deletion PhysicsTools/PatUtils/plugins/CorrectedPATMETProducer.cc
View file
Open in desktop
Expand Up @@ -15,6 +15,7 @@
#include "DataFormats/METReco/interface/CorrMETData.h"

#include "JetMETCorrections/Type1MET/interface/AddCorrectionsToGenericMET.h"
#include "RecoMET/METAlgorithms/interface/METSignificance.h"

#include <vector>

Expand Down Expand Up @@ -49,7 +50,6 @@ class CorrectedPATMETProducer : public edm::stream::EDProducer<>

edm::EDGetTokenT<METCollection> token_;


void produce(edm::Event& evt, const edm::EventSetup& es) override
{
edm::Handle<METCollection> srcMETCollection;
Expand All @@ -62,6 +62,14 @@ class CorrectedPATMETProducer : public edm::stream::EDProducer<>
pat::MET outMET(corrMET, srcMET);

auto product = std::make_unique<METCollection>();

reco::METCovMatrix cov=srcMET.getSignificanceMatrix();
if( !(cov(0,0)==0 && cov(0,1)==0 && cov(1,0)==0 && cov(1,1)==0) ) {
outMET.setSignificanceMatrix(cov);
double metSig=metsig::METSignificance::getSignificance(cov, outMET);
outMET.setMETSignificance(metSig);
}

product->push_back(outMET);
evt.put(std::move(product));
}
Expand Down
8 changes: 7 additions & 1 deletion PhysicsTools/PatUtils/python/tools/runMETCorrectionsAndUncertainties.py
View file
Open in desktop
Expand Up @@ -966,7 +966,10 @@ def getVariations(self, process, metModName, identifier,preId, objectCollection,
smear=False
if "Smear" in metModName:
smear=True

else:
smear=True
varyByNsigmas=101

shiftedCollModules['Up'] = self.createShiftedJetResModule(process, smear, objectCollection, +1.*varyByNsigmas,
"Up", postfix)
shiftedCollModules['Down'] = self.createShiftedJetResModule(process, smear, objectCollection, -1.*varyByNsigmas,
Expand Down Expand Up @@ -1464,6 +1467,9 @@ def miniAODConfigurationPost(self, process, postfix):
if hasattr(process, "patPFMetTxyCorr"+postfix):
getattr(process, "patPFMetTxyCorr"+postfix).vertexCollection = cms.InputTag("offlineSlimmedPrimaryVertices")

if self._parameters['computeUncertainties'].value:
getattr(process, "shiftedPatJetResDown"+postfix).genJets = cms.InputTag("slimmedGenJets")
getattr(process, "shiftedPatJetResUp"+postfix).genJets = cms.InputTag("slimmedGenJets")


def miniAODConfiguration(self, process, pfCandCollection, jetCollection,
Expand Down
2 changes: 1 addition & 1 deletion RecoMET/METAlgorithms/interface/METSignificance.h
View file
Open in desktop
Expand Up @@ -46,7 +46,7 @@ namespace metsig {
JME::JetResolutionScaleFactor & resSFObj,
bool isRealData);

double getSignificance(const reco::METCovMatrix& cov, const reco::MET& met ) const;
static double getSignificance(const reco::METCovMatrix& cov, const reco::MET& met );

private:
bool cleanJet(const reco::Jet& jet,
Expand Down
20 changes: 11 additions & 9 deletions RecoMET/METAlgorithms/src/METSignificance.cc
View file
Open in desktop
Expand Up @@ -20,16 +20,18 @@ Description: [one line class summary]

metsig::METSignificance::METSignificance(const edm::ParameterSet& iConfig) {

edm::ParameterSet cfgParams = iConfig.getParameter<edm::ParameterSet>("parameters");
// if(iConfig.exists("parameters")) {
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

is this commented out code needed?
Please remove or add comments inline in the code why the commented out block is relevant

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

looks like this commented out "if" is also the reason for the unnecessary re-indentation of the code in metsig::METSignificance::METSignificance

Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

not needed anymore, forgot to remove it during the development period, going to be removed

edm::ParameterSet cfgParams = iConfig.getParameter<edm::ParameterSet>("parameters");

double dRmatch = cfgParams.getParameter<double>("dRMatch");
dR2match_ = dRmatch*dRmatch;
double dRmatch = cfgParams.getParameter<double>("dRMatch");
dR2match_ = dRmatch*dRmatch;

jetThreshold_ = cfgParams.getParameter<double>("jetThreshold");
jetEtas_ = cfgParams.getParameter<std::vector<double> >("jeta");
jetParams_ = cfgParams.getParameter<std::vector<double> >("jpar");
pjetParams_ = cfgParams.getParameter<std::vector<double> >("pjpar");
// }

jetThreshold_ = cfgParams.getParameter<double>("jetThreshold");
jetEtas_ = cfgParams.getParameter<std::vector<double> >("jeta");
jetParams_ = cfgParams.getParameter<std::vector<double> >("jpar");
pjetParams_ = cfgParams.getParameter<std::vector<double> >("pjpar");

}

metsig::METSignificance::~METSignificance() {
Expand Down Expand Up @@ -172,7 +174,7 @@ metsig::METSignificance::getCovariance(const edm::View<reco::Jet>& jets,
}

double
metsig::METSignificance::getSignificance(const reco::METCovMatrix& cov, const reco::MET& met) const {
metsig::METSignificance::getSignificance(const reco::METCovMatrix& cov, const reco::MET& met) {

// covariance matrix determinant
double det = cov(0,0)*cov(1,1) - cov(0,1)*cov(1,0);
Expand Down