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

Standard muon selectors #20612

Merged
merged 32 commits into from Oct 8, 2017
Merged

Standard muon selectors #20612

Show file tree
Hide file tree
Changes from 15 commits
Commits
Show all changes
32 commits
Select commit Hold shift + click to select a range
b9c86bc
Added bitmap and needed methods and types for standard muon selectors
drkovalskyi Sep 19, 2017
e317148
Added muon MVA value as well as some jet-related variables widely use…
drkovalskyi Sep 19, 2017
4886664
Added a function to set selector flags for selections based on reco::…
drkovalskyi Sep 19, 2017
8bf0616
Compute standard muon selector flags in RECO
drkovalskyi Sep 19, 2017
ec7e863
Integrated Muon MVA and added calculation of MiniIsolation and MVA wo…
drkovalskyi Sep 19, 2017
1832ba3
Added bitmap and needed methods and types for standard muon selectors
drkovalskyi Sep 19, 2017
ba73ec9
Added muon MVA value as well as some jet-related variables widely use…
drkovalskyi Sep 19, 2017
0c5ab4a
Added a function to set selector flags for selections based on reco::…
drkovalskyi Sep 19, 2017
dd8206b
Compute standard muon selector flags in RECO
drkovalskyi Sep 19, 2017
04cfab4
Integrated Muon MVA and added calculation of MiniIsolation and MVA wo…
drkovalskyi Sep 19, 2017
ba190d9
updated checksum and version
drkovalskyi Sep 19, 2017
e8d20f7
updated checksum
drkovalskyi Sep 21, 2017
bf7f6b5
wrong header file
drkovalskyi Sep 21, 2017
738fd09
protection
drkovalskyi Sep 21, 2017
42f4170
disable muon mva unless we are producing miniaod
drkovalskyi Sep 21, 2017
718dc74
style improvements
drkovalskyi Sep 21, 2017
9dd45cf
removed unnecesary protection
drkovalskyi Sep 21, 2017
8689825
renamed define statement to follow the convention
drkovalskyi Sep 21, 2017
7fbd779
removed commented out code
drkovalskyi Sep 21, 2017
a26e056
use plain float type
drkovalskyi Sep 21, 2017
a8dcb6f
made primary vertex collection configurable
drkovalskyi Sep 22, 2017
59ad7ca
changed the primary vertex collection
drkovalskyi Sep 22, 2017
57e47ad
added TrackerHighPtMuon
drkovalskyi Sep 22, 2017
30998ea
added hip mitigation
drkovalskyi Sep 28, 2017
f5af074
hip mitigation for Run2016BCDEF
drkovalskyi Sep 28, 2017
9cfcfb7
cleanup
drkovalskyi Sep 29, 2017
7902379
cleanup
drkovalskyi Sep 29, 2017
75d3480
extract timing for tracker-only muons
drkovalskyi Sep 30, 2017
f4cb782
cleanup
drkovalskyi Sep 30, 2017
9cc9137
reverted changes
drkovalskyi Oct 1, 2017
79958f0
cleanup
drkovalskyi Oct 1, 2017
0dce48d
resolved conflict
drkovalskyi Oct 6, 2017
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
40 changes: 39 additions & 1 deletion DataFormats/MuonReco/interface/Muon.h
View file
Open in desktop
Expand Up @@ -182,6 +182,43 @@ namespace reco {
enum ArbitrationType { NoArbitration, SegmentArbitration, SegmentAndTrackArbitration, SegmentAndTrackArbitrationCleaned,
RPCHitAndTrackArbitration, GEMSegmentAndTrackArbitration, ME0SegmentAndTrackArbitration };

///
/// ====================== STANDARD SELECTORS ===========================
///
enum Selector {
Copy link
Contributor

Choose a reason for hiding this comment

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

IDSelector may be more descriptive

Copy link
Contributor Author

Choose a reason for hiding this comment

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

It's not ID. It's any selector including isolation.

CutBasedIdLoose = 1UL<< 0,
CutBasedIdMedium = 1UL<< 1,
CutBasedIdMediumPrompt = 1UL<< 2, // IP cuts
CutBasedIdTight = 1UL<< 3,
CutBasedIdHighPt = 1UL<< 4,
PFIsoVeryLoose = 1UL<< 5, // reliso<0.40
PFIsoLoose = 1UL<< 6, // reliso<0.25
PFIsoMedium = 1UL<< 7, // reliso<0.20
PFIsoTight = 1UL<< 8, // reliso<0.15
PFIsoVeryTight = 1UL<< 9, // reliso<0.10
TkIsoLoose = 1UL<<10, // reliso<0.10
TkIsoTight = 1UL<<11, // reliso<0.05
SoftCutBasedId = 1UL<<12,
SoftMvaId = 1UL<<13,
MvaLoose = 1UL<<14,
MvaMedium = 1UL<<15,
MvaTight = 1UL<<16,
MiniIsoLoose = 1UL<<17, // reliso<0.40
MiniIsoMedium = 1UL<<18, // reliso<0.20
MiniIsoTight = 1UL<<19, // reliso<0.10
MiniIsoVeryTight = 1UL<<20 // reliso<0.05
};

bool passed( unsigned int selection ) const { return (selectors_ & selection)==selection; }
Copy link
Contributor

Choose a reason for hiding this comment

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

enums are known in ROOT as well.
It would be much easier to track the meaning of the call to "passed(id)" if it's done using the Selector enum as type.

Copy link
Contributor Author

Choose a reason for hiding this comment

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

Here the choice of type is intentional to be able to do passed(TkIsoLoose|CutBasedIdHighPt), which would check that all flags are set. That's the most typical use case and it's a bit ugly and error prone to do it by hand on the bitmap itself.

unsigned int getSelectionMask() const { return selectors_; }
Copy link
Contributor

Choose a reason for hiding this comment

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

there are no accessors starting with "get" in this class.
Better keep the style and call this e.g. selectionMask()

Copy link
Contributor

Choose a reason for hiding this comment

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

also, there is some naming inconsistency selectors_ vs selectionMask.
better use one to reduce the number of variants to remember in possible future updates

Copy link
Contributor Author

Choose a reason for hiding this comment

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

fixed

void setSelectionMask( unsigned int mask ){ selectors_ = mask; }
void setSelector(Selector selector, bool passed){
if (passed)
selectors_ |= selector;
else
selectors_ &= ~selector;
}

///
/// ====================== USEFUL METHODs ===========================
///
Expand Down Expand Up @@ -288,7 +325,8 @@ namespace reco {
/// get pointers to best segment and corresponding chamber in vector of chambers
std::pair<const MuonChamberMatch*,const MuonSegmentMatch*> pair( const std::vector<const MuonChamberMatch*> &,
ArbitrationType type = SegmentAndTrackArbitration ) const;

/// selector bitmap
unsigned int selectors_;
Copy link
Contributor

Choose a reason for hiding this comment

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

is it safe to commit to 32 bits now, considering that 22 are already taken?
uint64_t could save from some headache in not so distant future.

Copy link
Contributor Author

Choose a reason for hiding this comment

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

32 is good enough. We need to resist a temptation to introduce any minor modification as a separate selector. If that's done I don't think we will reach the limit any time soon. In other words it's good to have some limits to avoid id proliferation.

Copy link
Contributor

Choose a reason for hiding this comment

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

Well, OK.
It's not obvious what happens with phase-2 specific IDs.
It's almost zero effort now to leave more room.

public:
/// get number of segments
int numberOfSegments( int station, int muonSubdetId, ArbitrationType type = SegmentAndTrackArbitration ) const;
Expand Down
2 changes: 2 additions & 0 deletions DataFormats/MuonReco/interface/MuonSelectors.h
View file
Open in desktop
Expand Up @@ -116,5 +116,7 @@ namespace muon {
int sharedSegments( const reco::Muon& muon1, const reco::Muon& muon2,
unsigned int segmentArbitrationMask = reco::MuonSegmentMatch::BestInChamberByDR ) ;

void setCutBasedSelectorFlags(reco::Muon& muon,
const reco::Vertex* vertex=0);
}
#endif
3 changes: 2 additions & 1 deletion DataFormats/MuonReco/src/Muon.cc
View file
Open in desktop
Expand Up @@ -17,6 +17,7 @@ Muon::Muon( Charge q, const LorentzVector & p4, const Point & vtx ) :
type_ = 0;
bestTunePTrackType_ = reco::Muon::None;
bestTrackType_ = reco::Muon::None;
selectors_ = 0;
}

Muon::Muon() {
Expand All @@ -29,7 +30,7 @@ Muon::Muon() {
type_ = 0;
bestTrackType_ = reco::Muon::None;
bestTunePTrackType_ = reco::Muon::None;

selectors_ = 0;
}

bool Muon::overlap( const Candidate & c ) const {
Expand Down
43 changes: 43 additions & 0 deletions DataFormats/MuonReco/src/MuonSelectors.cc
View file
Open in desktop
Expand Up @@ -929,3 +929,46 @@ int muon::sharedSegments( const reco::Muon& mu, const reco::Muon& mu2, unsigned

return ret;
}

void muon::setCutBasedSelectorFlags(reco::Muon& muon,
const reco::Vertex* vertex)
{
// https://twiki.cern.ch/twiki/bin/viewauth/CMS/SWGuideMuonIdRun2
unsigned int selection = muon.getSelectionMask();
// Compute Id and Isolation variables
double chIso = muon.pfIsolationR04().sumChargedHadronPt;
double nIso = muon.pfIsolationR04().sumNeutralHadronEt;
double phoIso = muon.pfIsolationR04().sumPhotonEt;
double puIso = muon.pfIsolationR04().sumPUPt;
double dbCorrectedIsolation = chIso + std::max( nIso + phoIso - .5*puIso, 0. ) ;
double dbCorectedRelIso = dbCorrectedIsolation/muon.pt();
double tkRelIso = muon.isolationR03().sumPt/muon.pt();
Copy link
Contributor

Choose a reason for hiding this comment

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

all original values are floats. Is there a need to go to double precision here?

Copy link
Contributor Author

Choose a reason for hiding this comment

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

Is there any harm? Impressed that you have noticed that, but I don't see any compiler warnings. These variables are used in calculations and it's better to use double precision for any calculations, so sooner or later translation from float to double will happen. See no reason to change that.

Copy link
Contributor

Choose a reason for hiding this comment

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

the harm is in more CPU cycles.
I don't believe that there is any gain in precision in this context.
I guess the argument is that the CPU loss is tiny and you wrote it this way.


// Base selectors
if (muon::isLooseMuon(muon)) selection |= reco::Muon::CutBasedIdLoose;
if (vertex){
if (muon::isTightMuon(muon,*vertex)) selection |= reco::Muon::CutBasedIdTight;
if (muon::isSoftMuon(muon,*vertex)) selection |= reco::Muon::SoftCutBasedId;
if (muon::isHighPtMuon(muon,*vertex)) selection |= reco::Muon::CutBasedIdHighPt;
}
if (muon::isMediumMuon(muon)){
selection |= reco::Muon::CutBasedIdMedium;
if ( vertex and
fabs(muon.muonBestTrack()->dz( vertex->position()))<0.1 and
fabs(muon.muonBestTrack()->dxy(vertex->position()))< 0.02 )
selection |= reco::Muon::CutBasedIdMediumPrompt;
}

// PF isolation
if (dbCorectedRelIso<0.40) selection |= reco::Muon::PFIsoVeryLoose;
if (dbCorectedRelIso<0.25) selection |= reco::Muon::PFIsoLoose;
if (dbCorectedRelIso<0.20) selection |= reco::Muon::PFIsoMedium;
if (dbCorectedRelIso<0.15) selection |= reco::Muon::PFIsoTight;
if (dbCorectedRelIso<0.10) selection |= reco::Muon::PFIsoVeryTight;
Copy link
Contributor

Choose a reason for hiding this comment

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

these magic constants (0.4, 0.25 etc) should be made consts of e.g. reco::Muon:: or in some other place and used consistently instead of entering values explicitly in multiple places and producers

Copy link
Contributor Author

Choose a reason for hiding this comment

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

They are used in one place only and that place sets the flags. Where do you see them used beside that place?

Copy link
Contributor

Choose a reason for hiding this comment

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

0.10 as a "Tight" appears multiple times in this PR, different variants of "Tight" though

Copy link
Contributor Author

Choose a reason for hiding this comment

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

0.10 doesn't mean tight, it's just a number that was selected for different isolations that happens to represent tight, i.e. not a magic number.

Copy link
Contributor

Choose a reason for hiding this comment

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

OK


// Tracker isolation
if (tkRelIso<0.10) selection |= reco::Muon::TkIsoLoose;
if (tkRelIso<0.05) selection |= reco::Muon::TkIsoTight;

muon.setSelectionMask(selection);
}
3 changes: 2 additions & 1 deletion DataFormats/MuonReco/src/classes_def.xml
View file
Open in desktop
@@ -1,5 +1,6 @@
<lcgdict>
<class name="reco::Muon" ClassVersion="18">
<class name="reco::Muon" ClassVersion="19">
<version ClassVersion="19" checksum="3201675353"/>
<version ClassVersion="18" checksum="2257088106"/>
<version ClassVersion="17" checksum="2296096670"/>
<version ClassVersion="16" checksum="1808419014"/>
Expand Down
16 changes: 16 additions & 0 deletions DataFormats/PatCandidates/interface/Muon.h
View file
Open in desktop
Expand Up @@ -265,6 +265,15 @@ namespace pat {
float pfEcalEnergy() const { return pfEcalEnergy_; }
void setPfEcalEnergy(float pfEcalEnergy) { pfEcalEnergy_ = pfEcalEnergy; }

/// near-by jet information
float jetPtRatio() const { return jetPtRatio_; }
float jetPtRel() const { return jetPtRel_; }
void setJetPtRatio(float jetPtRatio){ jetPtRatio_ = jetPtRatio; }
void setJetPtRel(float jetPtRel){ jetPtRel_ = jetPtRel; }

/// Muon MVA
float mvaValue() const { return mvaValue_; }
void setMvaValue(float mva){ mvaValue_ = mva; }
protected:

// ---- for content embedding ----
Expand Down Expand Up @@ -332,6 +341,13 @@ namespace pat {
float puppiNoLeptonsPhotonIso_;

float pfEcalEnergy_;

/// near-by jet information
float jetPtRatio_;
float jetPtRel_;

/// Muon MVA
float mvaValue_;
};


Expand Down
20 changes: 16 additions & 4 deletions DataFormats/PatCandidates/src/Muon.cc
View file
Open in desktop
Expand Up @@ -29,7 +29,10 @@ Muon::Muon() :
normChi2_(0.0),
cachedNumberOfValidHits_(false),
numberOfValidHits_(0),
pfEcalEnergy_(0)
pfEcalEnergy_(0),
jetPtRatio_(0),
jetPtRel_(0),
mvaValue_(0)
{
initImpactParameters();
}
Expand All @@ -53,7 +56,10 @@ Muon::Muon(const reco::Muon & aMuon) :
normChi2_(0.0),
cachedNumberOfValidHits_(false),
numberOfValidHits_(0),
pfEcalEnergy_(0)
pfEcalEnergy_(0),
jetPtRatio_(0),
jetPtRel_(0),
mvaValue_(0)
{
initImpactParameters();
}
Expand All @@ -77,7 +83,10 @@ Muon::Muon(const edm::RefToBase<reco::Muon> & aMuonRef) :
normChi2_(0.0),
cachedNumberOfValidHits_(0),
numberOfValidHits_(0),
pfEcalEnergy_(0)
pfEcalEnergy_(0),
jetPtRatio_(0),
jetPtRel_(0),
mvaValue_(0)
{
initImpactParameters();
}
Expand All @@ -101,7 +110,10 @@ Muon::Muon(const edm::Ptr<reco::Muon> & aMuonRef) :
normChi2_(0.0),
cachedNumberOfValidHits_(0),
numberOfValidHits_(0),
pfEcalEnergy_(0)
pfEcalEnergy_(0),
jetPtRatio_(0),
jetPtRel_(0),
mvaValue_(0)
{
initImpactParameters();
}
Expand Down
3 changes: 2 additions & 1 deletion DataFormats/PatCandidates/src/classes_def_objects.xml
View file
Open in desktop
Expand Up @@ -60,7 +60,8 @@
</ioread>


<class name="pat::Muon" ClassVersion="22">
<class name="pat::Muon" ClassVersion="23">
<version ClassVersion="23" checksum="3471390119"/>
<version ClassVersion="22" checksum="1816533594"/>
<version ClassVersion="21" checksum="1539691612"/>
<version ClassVersion="20" checksum="357097717"/>
Expand Down
46 changes: 46 additions & 0 deletions PhysicsTools/PatAlgos/interface/MuonMvaEstimator.h
View file
Open in desktop
@@ -0,0 +1,46 @@
#ifndef __RecoMuon_MuonIdentification_MuonMvaEstimator__
#define __RecoMuon_MuonIdentification_MuonMvaEstimator__
Copy link
Contributor

Choose a reason for hiding this comment

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

PhysicsTools_PatAlgos not RecoMuon_MuonIdentification

Copy link
Contributor Author

Choose a reason for hiding this comment

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

fixed

#include "DataFormats/MuonReco/interface/Muon.h"
#include "DataFormats/PatCandidates/interface/Muon.h"
#include "DataFormats/BTauReco/interface/JetTag.h"
#include "TMVA/Reader.h"

namespace reco {
class JetCorrector;
}
namespace pat {
class MuonMvaEstimator{
public:
MuonMvaEstimator();
void initialize(std::string weightsfile,
Copy link
Contributor

Choose a reason for hiding this comment

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

initialize(std::string const& weightsfile,

Copy link
Contributor Author

Choose a reason for hiding this comment

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

It's called once in the PATMuonProducer::produce. It's a waste of time to change it.

float dRmax);
void computeMva(const pat::Muon& imuon,
const reco::Vertex& vertex,
const reco::JetTagCollection& bTags,
const reco::JetCorrector* correctorL1=0,
const reco::JetCorrector* correctorL1L2L3Res=0);
float mva() const {return mva_;}
float jetPtRatio() const {return jetPtRatio_;}
float jetPtRel() const {return jetPtRel_;}
private:
TMVA::Reader tmvaReader_;
bool initialized_;
float mva_;
float dRmax_;

/// MVA VAriables
Float_t pt_;
Copy link
Contributor

Choose a reason for hiding this comment

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

Float_t -> float

Copy link
Contributor Author

Choose a reason for hiding this comment

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

The MVA reader chocks if it's not Float_t.

Copy link
Contributor Author

Choose a reason for hiding this comment

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

Seems like I was wrong (it chocked on double, not float). Float seems to be Ok, so I will change it to float.

Float_t eta_;
Float_t jetNDauCharged_;
Float_t miniRelIsoCharged_;
Float_t miniRelIsoNeutral_;
Float_t jetPtRel_;
Float_t jetPtRatio_;
Float_t jetBTagCSV_;
Float_t sip_;
Float_t log_abs_dxyBS_;
Float_t log_abs_dzPV_;
Float_t segmentCompatibility_;
};
}
#endif