-
Notifications
You must be signed in to change notification settings - Fork 4.2k
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
Adding chsMET and trackMET to miniAOD #20655
Changes from 9 commits
046cf84
ce18926
79b2107
8aaf704
123a086
07f5ac9
b226afb
159dd89
d81d212
d139b77
28f06fe
aab6709
bd91dae
7fb017a
File filter
Filter by extension
Conversations
Jump to
Diff view
Diff view
There are no files selected for viewing
Original file line number | Diff line number | Diff line change |
---|---|---|
|
@@ -57,12 +57,12 @@ namespace pat { | |
/// constructor for corrected METs (keeping specific informations from src MET) | ||
MET(const reco::MET & corMET, const MET& srcMET ); | ||
/// destructor | ||
virtual ~MET(); | ||
~MET() override; | ||
|
||
MET& operator=(MET const&); | ||
|
||
/// required reimplementation of the Candidate's clone method | ||
virtual MET * clone() const { return new MET(*this); } | ||
MET * clone() const override { return new MET(*this); } | ||
|
||
// ---- methods for generated MET link ---- | ||
/// return the associated GenMET | ||
|
@@ -76,7 +76,6 @@ namespace pat { | |
// get the MET significance | ||
double metSignificance() const; | ||
|
||
|
||
// ---- methods for uncorrected MET ---- | ||
// Methods not yet defined | ||
//float uncorrectedPt() const; | ||
|
@@ -159,12 +158,12 @@ namespace pat { | |
enum METCorrectionLevel { | ||
Raw=0, Type1=1, Type01=2, TypeXY=3, Type1XY=4, Type01XY=5, | ||
Type1Smear=6, Type01Smear=7, Type1SmearXY=8, | ||
Type01SmearXY=9, RawCalo=10, METCorrectionLevelSize=11 | ||
Type01SmearXY=9, RawCalo=10, RawChs=11, RawTrk=12, METCorrectionLevelSize=13 | ||
}; | ||
enum METCorrectionType { | ||
None=0, T1=1, T0=2, TXY=3, TXYForRaw=4, | ||
TXYForT01=5, TXYForT1Smear=6, TXYForT01Smear=7, | ||
Smear=8, Calo=9, METCorrectionTypeSize=10 | ||
Smear=8, Calo=9, Chs=10, Trk=11, METCorrectionTypeSize=12 | ||
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Please indent as for the other lines |
||
}; | ||
|
||
struct Vector2 { | ||
|
Original file line number | Diff line number | Diff line change |
---|---|---|
|
@@ -172,6 +172,56 @@ def miniAOD_customizeCommon(process): | |
del process.slimmedMETsNoHF.caloMET | ||
# ================== NoHF pfMET | ||
|
||
# ================== CHSMET | ||
process.CHSCands = cms.EDFilter("CandPtrSelector", | ||
src=cms.InputTag("packedPFCandidates"), | ||
cut=cms.string("fromPV(0) > 0") | ||
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. following up on today's discussion and last week mail thread, are we sure we want to use this definition for TkMET? how about using this definition There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. given that we do not have a better definition right now (studies started), we keep it as is, and we can systematically update when we decide. This definition is hand in hand correlated with what jets and is used in the JEC computation. therefore the met has to be built from the pf candidates identical to what jets are reconstructed with. |
||
) | ||
task.add(process.CHSCands) | ||
|
||
process.pfMetCHS = cms.EDProducer("PFMETProducer", | ||
src = cms.InputTag("CHSCands"), | ||
alias = cms.string('pfMet'), | ||
globalThreshold = cms.double(0.0), | ||
calculateSignificance = cms.bool(False), | ||
) | ||
task.add(process.pfMetCHS) | ||
|
||
addMETCollection(process, | ||
labelName = "patChsMet", | ||
metSource = "pfMetCHS" | ||
) | ||
|
||
process.patChsMet.computeMETSignificance = cms.bool(False) | ||
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. hi @zdemirag - please follow a convention here. Eg CHS or Chs (seems CHS is the established one) (and similarly for TrkMet) There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. The differences between patXXXMet and MetXXX is the pat vs reco modules. therefore the notation is self consistent. There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Its not XXXMet vs MetXXX - rather its "XXX" being consistent that I'm requesting.It looks like CHS is currently universal in CMSSW python: [dlange@lxplus098 CMSSW_9_4_0_pre2]$ cmsglimpse CHS | grep python | wc |
||
|
||
# ================== CHSMET | ||
|
||
# ================== TrkMET | ||
process.TrkCands = cms.EDFilter("CandPtrSelector", | ||
src=cms.InputTag("packedPFCandidates"), | ||
cut=cms.string("fromPV(0) > 0 && charge()!=0") | ||
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. @zdemirag I meant to comment on this one (tkMET) not on the one above (chsMET), for this I do not see any reason to be consistent with CHS, we should rather have something that has the features expected by the users (i.e. stability against PU) hence we could use a tighter PV association such as There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. ... or perhaps even There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. ah thats fine then. I thought you were referring to CHS. Give me (tops) one-two days, I have preliminary studies that can back up the decision we make. but in general, I see no problem updating this portion of the code in this PR and will update. There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. i'll move the definition to ==> "charge()!=0 && pvAssociationQuality()>=4 && vertexRef().key()==0" , it is because >=4 has slightly lower energy removal for barrel jets (~5%) at all pts in a non-pu sample. There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. I guess mini was updated after nano was already implemented |
||
) | ||
task.add(process.TrkCands) | ||
|
||
process.pfMetTrk = cms.EDProducer("PFMETProducer", | ||
src = cms.InputTag("TrkCands"), | ||
alias = cms.string('pfMet'), | ||
globalThreshold = cms.double(0.0), | ||
calculateSignificance = cms.bool(False), | ||
) | ||
|
||
task.add(process.pfMetTrk) | ||
|
||
addMETCollection(process, | ||
labelName = "patTrkMet", | ||
metSource = "pfMetTrk" | ||
) | ||
|
||
process.patTrkMet.computeMETSignificance = cms.bool(False) | ||
|
||
# ================== TrkMET | ||
|
||
|
||
## PU JetID | ||
process.load("RecoJets.JetProducers.PileupJetID_cfi") | ||
task.add(process.pileUpJetIDTask) | ||
|
Original file line number | Diff line number | Diff line change |
---|---|---|
|
@@ -330,8 +330,8 @@ def toolCode(self, process): | |
postfix) | ||
|
||
#pre-preparation to run over miniAOD | ||
if onMiniAOD: | ||
self.miniAODConfigurationPre(process, patMetModuleSequence, postfix) | ||
if onMiniAOD: | ||
self.miniAODConfigurationPre(process, patMetModuleSequence, pfCandCollection, postfix) | ||
|
||
#default MET production | ||
self.produceMET(process, metType,patMetModuleSequence, postfix) | ||
|
@@ -1509,7 +1509,7 @@ def ak4JetReclustering(self,process, pfCandCollection, patMetModuleSequence, pos | |
|
||
|
||
|
||
def miniAODConfigurationPre(self, process, patMetModuleSequence, postfix): | ||
def miniAODConfigurationPre(self, process, patMetModuleSequence, pfCandCollection, postfix): | ||
|
||
#extractor for caloMET === temporary for the beginning of the data taking | ||
self.extractMET(process,"rawCalo",patMetModuleSequence,postfix) | ||
|
@@ -1520,6 +1520,56 @@ def miniAODConfigurationPre(self, process, patMetModuleSequence, postfix): | |
) | ||
getattr(process,"patCaloMet").addGenMET = False | ||
|
||
##adding the necessary chs and track met configuration | ||
task = getPatAlgosToolsTask(process) | ||
|
||
pfChs = cms.EDFilter("CandPtrSelector", src = pfCandCollection, cut = cms.string("fromPV(0)>0")) | ||
addToProcessAndTask("pfChs", pfChs, process, task) | ||
pfMetChs = cms.EDProducer("PFMETProducer", | ||
src = cms.InputTag('pfChs'), | ||
alias = cms.string('pfMet'), | ||
globalThreshold = cms.double(0.0), | ||
calculateSignificance = cms.bool(False), | ||
) | ||
|
||
addToProcessAndTask("pfMetChs", pfMetChs, process, task) | ||
|
||
addMETCollection(process, | ||
labelName = "patChsMet", | ||
metSource = "pfMetChs" | ||
) | ||
|
||
getattr(process,"patChsMet").computeMETSignificance = cms.bool(False) | ||
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. here and below - these getattr seem unnecessary. Eg, process.patChsMet does the job and is more readable. There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. this was done so we don't waste time running the significance code since all we need is the central value of the met. There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. hi @zdemirag - i think you answered a different question. I'm purely suggesting a more clear syntax not a functionality change Eg process.patChsMet.computeMETSignificant = .... There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. ah sorry, will do it right away. |
||
getattr(process,"patChsMet").addGenMET = False | ||
|
||
patMetModuleSequence += getattr(process, "pfChs") | ||
patMetModuleSequence += getattr(process, "pfMetChs") | ||
patMetModuleSequence += getattr(process, "patChsMet") | ||
|
||
pfTrk = cms.EDFilter("CandPtrSelector", src = pfCandCollection, cut = cms.string("fromPV(0) > 0 && charge()!=0")) | ||
addToProcessAndTask("pfTrk", pfTrk, process, task) | ||
pfMetTrk = cms.EDProducer("PFMETProducer", | ||
src = cms.InputTag('pfTrk'), | ||
alias = cms.string('pfMet'), | ||
globalThreshold = cms.double(0.0), | ||
calculateSignificance = cms.bool(False), | ||
) | ||
|
||
addToProcessAndTask("pfMetTrk", pfMetTrk, process, task) | ||
|
||
addMETCollection(process, | ||
labelName = "patTrkMet", | ||
metSource = "pfMetTrk" | ||
) | ||
|
||
getattr(process,"patTrkMet").computeMETSignificance = cms.bool(False) | ||
getattr(process,"patTrkMet").addGenMET = False | ||
|
||
patMetModuleSequence += getattr(process, "pfTrk") | ||
patMetModuleSequence += getattr(process, "pfMetTrk") | ||
patMetModuleSequence += getattr(process, "patTrkMet") | ||
|
||
|
||
def miniAODConfigurationPost(self, process, postfix): | ||
|
||
if self._parameters["metType"].value == "PF": | ||
|
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Please indent as for the other lines