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.

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

Add DeepMET into NanoAOD for UltraLegacy in 10_6_X (backport of 30291 and 30926) #31120

Merged
merged 9 commits into from
Aug 26, 2020
Merged
27 changes: 27 additions & 0 deletions PhysicsTools/NanoAOD/python/met_cff.py
Original file line number Diff line number Diff line change
Expand Up @@ -106,6 +106,32 @@
),
)

deepMetResolutionTuneTable = cms.EDProducer("SimpleCandidateFlatTableProducer",
# current deepMets are saved in slimmedMETs in MiniAOD,
# in the same way as chsMet/TkMET
src = metTable.src,
name = cms.string("DeepMETResolutionTune"),
doc = cms.string("Deep MET trained with resolution tune"),
singleton = cms.bool(True), # there's always exactly one MET per event
extension = cms.bool(False), # this is the main table for the MET
variables = cms.PSet(#NOTA BENE: we don't copy PTVars here!
pt = Var("corPt('RawDeepResolutionTune')", float, doc="DeepMET ResolutionTune pt",precision=-1),
phi = Var("corPhi('RawDeepResolutionTune')", float, doc="DeepmET ResolutionTune phi",precision=12),
),
)

deepMetResponseTuneTable = cms.EDProducer("SimpleCandidateFlatTableProducer",
src = metTable.src,
name = cms.string("DeepMETResponseTune"),
doc = cms.string("Deep MET trained with extra response tune"),
singleton = cms.bool(True), # there's always exactly one MET per event
extension = cms.bool(False), # this is the main table for the MET
variables = cms.PSet(#NOTA BENE: we don't copy PTVars here!
pt = Var("corPt('RawDeepResponseTune')", float, doc="DeepMET ResponseTune pt",precision=-1),
phi = Var("corPhi('RawDeepResponseTune')", float, doc="DeepMET ResponseTune phi",precision=12),
),
)

metFixEE2017Table = metTable.clone()
metFixEE2017Table.src = cms.InputTag("slimmedMETsFixEE2017")
metFixEE2017Table.name = cms.string("METFixEE2017")
Expand All @@ -127,6 +153,7 @@


metTables = cms.Sequence( metTable + rawMetTable + caloMetTable + puppiMetTable + rawPuppiMetTable+ tkMetTable + chsMetTable)
deepMetTables = cms.Sequence( deepMetResolutionTuneTable + deepMetResponseTuneTable )
_withFixEE2017_sequence = cms.Sequence(metTables.copy() + metFixEE2017Table)
for modifier in run2_nanoAOD_94XMiniAODv1, run2_nanoAOD_94XMiniAODv2:
modifier.toReplaceWith(metTables,_withFixEE2017_sequence) # only in old miniAOD, the new ones will come from the UL rereco
Expand Down
14 changes: 14 additions & 0 deletions PhysicsTools/NanoAOD/python/nanoDQM_cfi.py
Original file line number Diff line number Diff line change
Expand Up @@ -32,6 +32,20 @@
Plot1D('rawPt', 'rawPt', 20, 5, 25, "pt()*jecFactor('Uncorrected')"),
)
),
DeepMETResolutionTune = cms.PSet(
sels = cms.PSet(),
plots = cms.VPSet(
Plot1D('phi', 'phi', 20, -3.14159, 3.14159, 'Deep MET Resolution Tune phi'),
Plot1D('pt', 'pt', 20, 0, 400, 'Deep MET Response Tune pt'),
)
),
DeepMETResponseTune = cms.PSet(
sels = cms.PSet(),
plots = cms.VPSet(
Plot1D('phi', 'phi', 20, -3.14159, 3.14159, 'Deep MET Response Tune phi'),
Plot1D('pt', 'pt', 20, 0, 400, 'Deep MET Response Tune pt'),
)
),
Electron = cms.PSet(
sels = cms.PSet(
Good = cms.string('pt > 15 && abs(dxy) < 0.2 && abs(dz) < 0.5 && cutBased >= 3 && miniPFRelIso_all < 0.4')
Expand Down
32 changes: 31 additions & 1 deletion PhysicsTools/NanoAOD/python/nano_cff.py
Original file line number Diff line number Diff line change
Expand Up @@ -179,10 +179,40 @@ def nanoAOD_addDeepInfo(process,addDeepBTag,addDeepFlavour):
process.updatedJets.jetSource="selectedUpdatedPatJetsWithDeepInfo"
return process

def nanoAOD_addDeepMET(process, addDeepMETProducer, ResponseTune_Graph):
if addDeepMETProducer:
# produce DeepMET on the fly if it is not in MiniAOD
print("add DeepMET Producers")
process.load('RecoMET.METPUSubtraction.deepMETProducer_cfi')
process.deepMETsResolutionTune = process.deepMETProducer.clone()
process.deepMETsResponseTune = process.deepMETProducer.clone()
#process.deepMETsResponseTune.graph_path = 'RecoMET/METPUSubtraction/data/deepmet/deepmet_resp_v1_2018.pb'
process.deepMETsResponseTune.graph_path = ResponseTune_Graph.value()
process.metTables += process.deepMetTables
return process

from PhysicsTools.PatUtils.tools.runMETCorrectionsAndUncertainties import runMetCorAndUncFromMiniAOD
#from PhysicsTools.PatAlgos.slimming.puppiForMET_cff import makePuppiesFromMiniAOD
def nanoAOD_recalibrateMETs(process,isData):
runMetCorAndUncFromMiniAOD(process,isData=isData)
# add DeepMETs
nanoAOD_DeepMET_switch = cms.PSet(
nanoAOD_addDeepMET_switch = cms.untracked.bool(True), # decide if DeeMET should be included in Nano
nanoAOD_produceDeepMET_switch = cms.untracked.bool(False), # decide if DeepMET should be computed on the fly
ResponseTune_Graph = cms.untracked.string('RecoMET/METPUSubtraction/data/deepmet/deepmet_resp_v1_2018.pb')
)
# compute DeepMETs for the eras before UL-ReminiAOD
(~run2_miniAOD_devel).toModify(nanoAOD_DeepMET_switch, nanoAOD_produceDeepMET_switch = cms.untracked.bool(True))
for modifier in run2_miniAOD_80XLegacy, run2_nanoAOD_94X2016:
modifier.toModify(nanoAOD_DeepMET_switch, ResponseTune_Graph=cms.untracked.string("RecoMET/METPUSubtraction/data/deepmet/deepmet_resp_v1_2016.pb"))
if nanoAOD_DeepMET_switch.nanoAOD_addDeepMET_switch:
process = nanoAOD_addDeepMET(process,
addDeepMETProducer=nanoAOD_DeepMET_switch.nanoAOD_produceDeepMET_switch,
ResponseTune_Graph=nanoAOD_DeepMET_switch.ResponseTune_Graph)

# if included in Nano, and not computed in the fly, then it should be extracted from minAOD
extractDeepMETs = nanoAOD_DeepMET_switch.nanoAOD_addDeepMET_switch and not nanoAOD_DeepMET_switch.nanoAOD_produceDeepMET_switch

runMetCorAndUncFromMiniAOD(process,isData=isData, extractDeepMETs=extractDeepMETs)
process.nanoSequenceCommon.insert(process.nanoSequenceCommon.index(process.jetSequence),cms.Sequence(process.fullPatMetSequence))
process.basicJetsForMetForT1METNano = process.basicJetsForMet.clone(
src = process.updatedJetsWithUserData.src,
Expand Down
67 changes: 37 additions & 30 deletions PhysicsTools/PatAlgos/plugins/RecoMETExtractor.cc
Original file line number Diff line number Diff line change
Expand Up @@ -8,57 +8,64 @@
using namespace pat;

RecoMETExtractor::RecoMETExtractor(const edm::ParameterSet& iConfig) {

edm::InputTag metIT = iConfig.getParameter<edm::InputTag>("metSource");
metSrcToken_ = consumes<pat::METCollection>(metIT);

std::string corLevel = iConfig.getParameter<std::string>("correctionLevel");

//all possible met flavors
if(corLevel=="raw") { corLevel_=pat::MET::Raw;}
else if(corLevel=="type1") { corLevel_=pat::MET::Type1;}
else if(corLevel=="type01") { corLevel_=pat::MET::Type01;}
else if(corLevel=="typeXY") { corLevel_=pat::MET::TypeXY;}
else if(corLevel=="type1XY") { corLevel_=pat::MET::Type1XY;}
else if(corLevel=="type01XY") { corLevel_=pat::MET::Type01XY;}
else if(corLevel=="type1Smear") { corLevel_=pat::MET::Type1Smear;}
else if(corLevel=="type01Smear") { corLevel_=pat::MET::Type01Smear;}
else if(corLevel=="type1SmearXY") { corLevel_=pat::MET::Type1SmearXY;}
else if(corLevel=="type01SmearXY") { corLevel_=pat::MET::Type01SmearXY;}
else if(corLevel=="rawCalo") { corLevel_=pat::MET::RawCalo;}
else {
if (corLevel == "raw") {
corLevel_ = pat::MET::Raw;
} else if (corLevel == "type1") {
corLevel_ = pat::MET::Type1;
} else if (corLevel == "type01") {
corLevel_ = pat::MET::Type01;
} else if (corLevel == "typeXY") {
corLevel_ = pat::MET::TypeXY;
} else if (corLevel == "type1XY") {
corLevel_ = pat::MET::Type1XY;
} else if (corLevel == "type01XY") {
corLevel_ = pat::MET::Type01XY;
} else if (corLevel == "type1Smear") {
corLevel_ = pat::MET::Type1Smear;
} else if (corLevel == "type01Smear") {
corLevel_ = pat::MET::Type01Smear;
} else if (corLevel == "type1SmearXY") {
corLevel_ = pat::MET::Type1SmearXY;
} else if (corLevel == "type01SmearXY") {
corLevel_ = pat::MET::Type01SmearXY;
} else if (corLevel == "rawCalo") {
corLevel_ = pat::MET::RawCalo;
} else if (corLevel == "rawDeepResponseTune") {
corLevel_ = pat::MET::RawDeepResponseTune;
} else if (corLevel == "rawDeepResolutionTune") {
corLevel_ = pat::MET::RawDeepResolutionTune;
} else {
//throw exception

}

// produces vector of recoMet
produces<std::vector<reco::MET> >();
}

RecoMETExtractor::~RecoMETExtractor() {}

RecoMETExtractor::~RecoMETExtractor() {

}
void RecoMETExtractor::produce(edm::StreamID streamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const {
edm::Handle<std::vector<pat::MET> > src;
iEvent.getByToken(metSrcToken_, src);
if (src->empty())
edm::LogError("RecoMETExtractor::produce") << "input reco MET collection is empty";

std::vector<reco::MET>* metCol = new std::vector<reco::MET>();

void RecoMETExtractor::produce(edm::StreamID streamID, edm::Event & iEvent,
const edm::EventSetup & iSetup) const {
reco::MET met(src->front().corSumEt(corLevel_), src->front().corP4(corLevel_), src->front().vertex());

edm::Handle<std::vector<pat::MET> > src;
iEvent.getByToken(metSrcToken_, src);
if(src->empty()) edm::LogError("RecoMETExtractor::produce") << "input reco MET collection is empty" ;
metCol->push_back(met);

std::vector<reco::MET> *metCol = new std::vector<reco::MET>();

reco::MET met(src->front().corSumEt(corLevel_), src->front().corP4(corLevel_), src->front().vertex() );

metCol->push_back( met );

std::unique_ptr<std::vector<reco::MET> > recoMETs(metCol);
iEvent.put(std::move(recoMETs));
}


#include "FWCore/Framework/interface/MakerMacros.h"

DEFINE_FWK_MODULE(RecoMETExtractor);
Original file line number Diff line number Diff line change
Expand Up @@ -91,6 +91,8 @@ def __init__(self):
"Exclude jets and PF candidates with EE noise characteristics (fix for 2017 run)", Type=bool)
self.addParameter(self._defaultParameters,'fixEE2017Params', {'userawPt': True, 'ptThreshold': 50.0, 'minEtaThreshold': 2.65, 'maxEtaThreshold': 3.139},
"Parameters dict for fixEE2017: userawPt, ptThreshold, minEtaThreshold, maxEtaThreshold", Type=dict)
self.addParameter(self._defaultParameters, 'extractDeepMETs', False,
"Extract DeepMETs from miniAOD, instead of recomputing them.", Type=bool)

#private parameters
self.addParameter(self._defaultParameters, 'Puppi', False,
Expand Down Expand Up @@ -136,6 +138,7 @@ def __call__(self, process,
onMiniAOD =None,
fixEE2017 =None,
fixEE2017Params =None,
extractDeepMETs =None,
postfix =None):
electronCollection = self.initializeInputTag(electronCollection, 'electronCollection')
photonCollection = self.initializeInputTag(photonCollection, 'photonCollection')
Expand Down Expand Up @@ -209,6 +212,8 @@ def __call__(self, process,
fixEE2017 = self._defaultParameters['fixEE2017'].value
if fixEE2017Params is None :
fixEE2017Params = self._defaultParameters['fixEE2017Params'].value
if extractDeepMETs is None :
extractDeepMETs = self._defaultParameters['extractDeepMETs'].value

self.setParameter('metType',metType),
self.setParameter('correctionLevel',correctionLevel),
Expand Down Expand Up @@ -242,6 +247,7 @@ def __call__(self, process,
self.setParameter('postfix',postfix),
self.setParameter('fixEE2017',fixEE2017),
self.setParameter('fixEE2017Params',fixEE2017Params),
self.setParameter('extractDeepMETs',extractDeepMETs),

#if mva/puppi MET, autoswitch to std jets
if metType == "MVA" or metType == "Puppi":
Expand Down Expand Up @@ -310,6 +316,7 @@ def toolCode(self, process):
postfix = self._parameters['postfix'].value
fixEE2017 = self._parameters['fixEE2017'].value
fixEE2017Params = self._parameters['fixEE2017Params'].value
extractDeepMETs = self._parameters['extractDeepMETs'].value

#prepare jet configuration
jetUncInfos = { "jCorrPayload":jetFlavor, "jCorLabelUpToL3":jetCorLabelUpToL3,
Expand Down Expand Up @@ -469,6 +476,11 @@ def toolCode(self, process):
#adding the slimmed MET
if hasattr(process, "patCaloMet"):
fullPatMetSequence +=getattr(process, "patCaloMet")
# include deepMETsResolutionTune and deepMETsResponseTune into fullPatMetSequence
if hasattr(process, "deepMETsResolutionTune"):
fullPatMetSequence +=getattr(process, "deepMETsResolutionTune")
if hasattr(process, "deepMETsResponseTune"):
fullPatMetSequence += getattr(process, "deepMETsResponseTune")
if hasattr(process, "slimmedMETs"+postfix):
fullPatMetSequence +=getattr(process, "slimmedMETs"+postfix)

Expand Down Expand Up @@ -1683,6 +1695,26 @@ def miniAODConfigurationPre(self, process, patMetModuleSequence, pfCandCollectio
)
getattr(process,"patCaloMet").addGenMET = False

# extract DeepMETs (ResponseTune and ResolutionTune) from miniAOD
if self._parameters["extractDeepMETs"].value:
self.extractMET(process, "rawDeepResponseTune", patMetModuleSequence, postfix)
deepMetResponseTuneName = "metrawDeepResponseTune" if hasattr(process, "metrawDeepResponseTune") else "metrawDeepResponseTune"+postfix
addMETCollection(process,
labelName = "deepMETsResponseTune",
metSource = deepMetResponseTuneName
)
getattr(process, "deepMETsResponseTune").addGenMET = False
getattr(process, "deepMETsResponseTune").computeMETSignificance = cms.bool(False)

self.extractMET(process, "rawDeepResolutionTune", patMetModuleSequence, postfix)
deepMetResolutionTuneName = "metrawDeepResolutionTune" if hasattr(process, "metrawDeepResolutionTune") else "metrawDeepResolutionTune"+postfix
addMETCollection(process,
labelName = "deepMETsResolutionTune",
metSource = deepMetResolutionTuneName
)
getattr(process, "deepMETsResolutionTune").addGenMET = False
getattr(process, "deepMETsResolutionTune").computeMETSignificance = cms.bool(False)

##adding the necessary chs and track met configuration
task = getPatAlgosToolsTask(process)

Expand Down Expand Up @@ -1783,6 +1815,10 @@ def miniAODConfiguration(self, process, pfCandCollection, jetCollection,
getattr(process,"slimmedMETs"+postfix).runningOnMiniAOD = True
getattr(process,"slimmedMETs"+postfix).t01Variation = cms.InputTag("slimmedMETs" if not self._parameters["Puppi"].value else "slimmedMETsPuppi",processName=cms.InputTag.skipCurrentProcess())

if hasattr(process, "deepMETsResolutionTune") and hasattr(process, "deepMETsResponseTune"):
# process includes producing/extracting deepMETsResolutionTune and deepMETsResponseTune
# add them to the slimmedMETs
getattr(process,"slimmedMETs"+postfix).addDeepMETs = True

#smearing and type0 variations not yet supported in reprocessing
#del getattr(process,"slimmedMETs"+postfix).t1SmearedVarsAndUncs
Expand Down Expand Up @@ -2054,6 +2090,7 @@ def runMetCorAndUncFromMiniAOD(process, metType="PF",
computeMETSignificance=True,
fixEE2017=False,
fixEE2017Params=None,
extractDeepMETs=False,
postfix=""):

runMETCorrectionsAndUncertainties = RunMETCorrectionsAndUncertainties()
Expand Down Expand Up @@ -2087,6 +2124,7 @@ def runMetCorAndUncFromMiniAOD(process, metType="PF",
postfix=postfix,
fixEE2017=fixEE2017,
fixEE2017Params=fixEE2017Params,
extractDeepMETs=extractDeepMETs,
)

#MET T1+Txy / Smear
Expand Down Expand Up @@ -2118,6 +2156,7 @@ def runMetCorAndUncFromMiniAOD(process, metType="PF",
postfix=postfix,
fixEE2017=fixEE2017,
fixEE2017Params=fixEE2017Params,
extractDeepMETs=extractDeepMETs,
)
#MET T1+Smear + uncertainties
runMETCorrectionsAndUncertainties(process, metType=metType,
Expand Down Expand Up @@ -2148,4 +2187,5 @@ def runMetCorAndUncFromMiniAOD(process, metType="PF",
postfix=postfix,
fixEE2017=fixEE2017,
fixEE2017Params=fixEE2017Params,
extractDeepMETs=extractDeepMETs,
)
2 changes: 1 addition & 1 deletion RecoMET/METPUSubtraction/plugins/DeepMETProducer.cc
Original file line number Diff line number Diff line change
Expand Up @@ -108,7 +108,7 @@ void DeepMETProducer::produce(edm::Event& event, const edm::EventSetup& setup) {
input_cat2_.tensor<float, 3>()(0, i_pf, 0) = pf.fromPV();

++i_pf;
if (i_pf > max_n_pf_) {
if (i_pf == max_n_pf_) {
break; // output a warning?
}
}
Expand Down