Skip to content

Commit

Permalink
Ad stage2 L1 trigger matching and possibility to sort matchings by L1…
Browse files Browse the repository at this point in the history 
… quality.
  • Loading branch information
Carlo Battilana committed Sep 28, 2016
1 parent d79b60d commit 94c279c
Show file tree
Hide file tree
Showing 7 changed files with 268 additions and 59 deletions.
92 changes: 82 additions & 10 deletions MuonAnalysis/MuonAssociators/interface/L1MuonMatcherAlgo.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -13,11 +13,13 @@


#include <cmath>
#include <type_traits>
#include "DataFormats/Math/interface/deltaR.h"
#include "DataFormats/Math/interface/deltaPhi.h"
#include "DataFormats/Candidate/interface/Candidate.h"
#include "DataFormats/TrackReco/interface/Track.h"
#include "DataFormats/L1Trigger/interface/L1MuonParticle.h"
#include "DataFormats/L1Trigger/interface/Muon.h"
#include "FWCore/Framework/interface/EventSetup.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"
#include "CommonTools/Utils/interface/StringCutObjectSelector.h"
Expand Down Expand Up @@ -70,11 +72,14 @@ class L1MuonMatcherAlgo {
return propagated.isValid() ? match(propagated, l1, deltaR, deltaPhi) : false;
}


/// Try to match one track to one L1. Return true if succeeded (and update deltaR, deltaPhi accordingly)
/// The preselection cut on L1, if specified in the config, is applied before the match
bool match(TrajectoryStateOnSurface & propagated, const l1extra::L1MuonParticle &l1, float &deltaR, float &deltaPhi) const ;



// Methods to match with vectors of legacy L1 muon trigger objects

/// Find the best match to L1, and return its index in the vector (and update deltaR, deltaPhi and propagated TSOS accordingly)
/// Returns -1 if the match fails
/// The preselection cut on L1, if specified in the config, is applied before the match
Expand All @@ -99,13 +104,49 @@ class L1MuonMatcherAlgo {
return propagated.isValid() ? match(propagated, l1, deltaR, deltaPhi) : -1;
}


/// Find the best match to L1, and return its index in the vector (and update deltaR, deltaPhi accordingly)
/// Returns -1 if the match fails
/// The preselection cut on L1, if specified in the config, is applied before the match
int match(TrajectoryStateOnSurface &propagated, const std::vector<l1extra::L1MuonParticle> &l1, float &deltaR, float &deltaPhi) const ;




// Methods to match with vectors of stage2 L1 muon trigger objects

/// Find the best match to stage2 L1, and return its index in the vector (and update deltaR, deltaPhi and propagated TSOS accordingly)
/// Returns -1 if the match fails
/// The preselection cut on stage2 L1, if specified in the config, is applied before the match
int match(const reco::Track &tk, const std::vector<l1t::Muon> &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const {
propagated = extrapolate(tk);
return propagated.isValid() ? match(propagated, l1, deltaR, deltaPhi) : -1;
}

/// Find the best match to stage2 L1, and return its index in the vector (and update deltaR, deltaPhi and propagated TSOS accordingly)
/// Returns -1 if the match fails
/// The preselection cut on stage2 L1, if specified in the config, is applied before the match
int match(const reco::Candidate &c, const std::vector<l1t::Muon> &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const {
propagated = extrapolate(c);
return propagated.isValid() ? match(propagated, l1, deltaR, deltaPhi) : -1;
}

/// Find the best match to stage2 L1, and return its index in the vector (and update deltaR, deltaPhi and propagated TSOS accordingly)
/// Returns -1 if the match fails
/// The preselection cut on stage 2 L1, if specified in the config, is applied before the match
int match(const SimTrack &tk, const edm::SimVertexContainer &vtxs, const std::vector<l1t::Muon> &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const {
propagated = extrapolate(tk, vtxs);
return propagated.isValid() ? match(propagated, l1, deltaR, deltaPhi) : -1;
}

/// Find the best match to stage 2 L1, and return its index in the vector (and update deltaR, deltaPhi accordingly)
/// Returns -1 if the match fails
/// The preselection cut on stage 2 L1, if specified in the config, is applied before the match
int match(TrajectoryStateOnSurface &propagated, const std::vector<l1t::Muon> &l1, float &deltaR, float &deltaPhi) const ;



// Generic matching

/// Find the best match to L1, and return its index in the vector (and update deltaR, deltaPhi and propagated TSOS accordingly)
/// Returns -1 if the match fails
/// Only the objects passing the selector will be allowed for the match.
Expand Down Expand Up @@ -133,52 +174,82 @@ class L1MuonMatcherAlgo {
template<typename Collection, typename Selector>
int matchGeneric(TrajectoryStateOnSurface &propagated, const Collection &l1, const Selector &sel, float &deltaR, float &deltaPhi) const ;


/// Add this offset to the L1 phi before doing the match, to correct for different scales in L1 vs offline
void setL1PhiOffset(double l1PhiOffset) { l1PhiOffset_ = l1PhiOffset; }

private:

template<class T> int genericQuality(const T & cand) const;

PropagateToMuon prop_;

typedef StringCutObjectSelector<l1extra::L1MuonParticle> L1Selector;
bool useStage2L1_;

typedef StringCutObjectSelector<reco::Candidate,true> L1Selector;
/// Preselection cut to apply to L1 candidates before matching
L1Selector preselectionCut_;

/// Matching cuts
double deltaR2_, deltaPhi_, deltaEta_;

/// Sort by deltaPhi or deltaEta instead of deltaR
enum SortBy { SortByDeltaR, SortByDeltaPhi, SortByDeltaEta, SortByPt };
enum SortBy { SortByDeltaR = 0, SortByDeltaPhi, SortByDeltaEta, SortByPt, SortByQual};
SortBy sortBy_;

/// offset to be added to the L1 phi before the match
double l1PhiOffset_;

};

template<>
inline int L1MuonMatcherAlgo::genericQuality<l1extra::L1MuonParticle>(const l1extra::L1MuonParticle & cand) const
{
return cand.gmtMuonCand().rank();
}

template<>
inline int L1MuonMatcherAlgo::genericQuality<l1t::Muon>(const l1t::Muon & cand) const
{
return cand.hwQual();
}

template<class T>
inline int L1MuonMatcherAlgo::genericQuality(const T & cand) const
{
return 0;
}



template<typename Collection, typename Selector>
int
L1MuonMatcherAlgo::matchGeneric(TrajectoryStateOnSurface &propagated, const Collection &l1s, const Selector &sel, float &deltaR, float &deltaPhi) const {
typedef typename Collection::value_type obj;
typedef typename Collection::value_type obj;

int match = -1;
double minDeltaPhi = deltaPhi_;
double minDeltaEta = deltaEta_;
double minDeltaR2 = deltaR2_;
double minPt = -1.0;
double maxPt = -1.0;
int maxQual = 0;
GlobalPoint pos = propagated.globalPosition();
for (int i = 0, n = l1s.size(); i < n; ++i) {
const obj &l1 = l1s[i];
if (sel(l1)) {
double thisDeltaPhi = ::deltaPhi(double(pos.phi()), l1.phi()+l1PhiOffset_);
double thisDeltaEta = pos.eta() - l1.eta();
double thisDeltaR2 = ::deltaR2(double(pos.eta()), double(pos.phi()), l1.eta(), l1.phi()+l1PhiOffset_);
int thisQual = genericQuality<obj>(l1);
double thisPt = l1.pt();
if ((fabs(thisDeltaPhi) < deltaPhi_) && (fabs(thisDeltaEta) < deltaEta_) && (thisDeltaR2 < deltaR2_)) { // check all
bool betterMatch = (match == -1);
switch (sortBy_) {
switch (sortBy_) {
case SortByDeltaR: betterMatch = (thisDeltaR2 < minDeltaR2); break;
case SortByDeltaEta: betterMatch = (fabs(thisDeltaEta) < fabs(minDeltaEta)); break;
case SortByDeltaPhi: betterMatch = (fabs(thisDeltaPhi) < fabs(minDeltaPhi)); break;
case SortByPt: betterMatch = (thisPt > minPt); break;
case SortByPt: betterMatch = (thisPt > maxPt); break;
// Quality is an int, adding sorting by pT in case of identical qualities
case SortByQual: betterMatch = (thisQual > maxQual || ((thisQual == maxQual) && (thisPt > maxPt))); break;
}
if (betterMatch) { // sort on one
match = i;
Expand All @@ -187,7 +258,8 @@ L1MuonMatcherAlgo::matchGeneric(TrajectoryStateOnSurface &propagated, const Coll
minDeltaR2 = thisDeltaR2;
minDeltaEta = thisDeltaEta;
minDeltaPhi = thisDeltaPhi;
minPt = thisPt;
maxQual = thisQual;
maxPt = thisPt;
}
}
}
Expand Down
2 changes: 2 additions & 0 deletions MuonAnalysis/MuonAssociators/interface/PropagateToMuon.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -70,6 +70,8 @@ class PropagateToMuon {
/// for cosmics, some things change: the along-opposite is not in-out, nor the innermost/outermost states are in-out really
bool cosmicPropagation_;

bool useMB2InOverlap_;

// needed services for track propagation
edm::ESHandle<MagneticField> magfield_;
edm::ESHandle<Propagator> propagator_, propagatorAny_, propagatorOpposite_;
Expand Down
128 changes: 107 additions & 21 deletions MuonAnalysis/MuonAssociators/plugins/L1MuonMatcher.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -48,13 +48,21 @@ namespace pat {
/// Tokens for input collections
edm::EDGetTokenT<edm::View<reco::Candidate> > recoToken_;
edm::EDGetTokenT<std::vector<l1extra::L1MuonParticle> > l1Token_;

edm::EDGetTokenT<l1t::MuonBxCollection> l1tToken_;

/// Labels to set as filter names in the output
std::string labelL1_, labelProp_;

/// Write out additional info as ValueMaps
bool writeExtraInfo_;

/// Allow to run both on legacy or stage2 (2016) L1 Muon trigger output
bool useStage2L1_;

/// Skim stage2 BX vector
int firstBX_;
int lastBX_;

/// Store extra information in a ValueMap
template<typename Hand, typename T>
void storeExtraInfo(edm::Event &iEvent,
Expand All @@ -68,11 +76,18 @@ namespace pat {
pat::L1MuonMatcher::L1MuonMatcher(const edm::ParameterSet & iConfig) :
matcher_(iConfig),
recoToken_(consumes<edm::View<reco::Candidate> >(iConfig.getParameter<edm::InputTag>("src"))),
l1Token_(consumes<std::vector<l1extra::L1MuonParticle> >(iConfig.getParameter<edm::InputTag>("matched"))),
labelL1_(iConfig.getParameter<std::string>( "setL1Label")),
labelL1_(iConfig.getParameter<std::string>("setL1Label")),
labelProp_(iConfig.getParameter<std::string>("setPropLabel")),
writeExtraInfo_(iConfig.getParameter<bool>("writeExtraInfo"))
writeExtraInfo_(iConfig.getParameter<bool>("writeExtraInfo")),
useStage2L1_(iConfig.getParameter<bool>("useStage2L1")),
firstBX_(iConfig.getParameter<int>("firstBX")),
lastBX_(iConfig.getParameter<int>("lastBX"))
{
if (useStage2L1_) {
l1tToken_ = consumes<l1t::MuonBxCollection>(iConfig.getParameter<edm::InputTag>("matched"));
} else {
l1Token_ = consumes<std::vector<l1extra::L1MuonParticle> >(iConfig.getParameter<edm::InputTag>("matched"));
}
produces<PATPrimitiveCollection>("l1muons"); // l1 in PAT format
produces<PATPrimitiveCollection>("propagatedReco"); // reco to muon station 2
produces<PATTriggerAssociation>("propagatedReco"); // asso reco to propagated reco
Expand All @@ -82,6 +97,10 @@ pat::L1MuonMatcher::L1MuonMatcher(const edm::ParameterSet & iConfig) :
produces<edm::ValueMap<float> >("deltaPhi");
produces<edm::ValueMap<int> >("quality");
produces<edm::ValueMap<int> >("bx");
if(useStage2L1_) {
produces<edm::ValueMap<int> >("iPhi");
produces<edm::ValueMap<int> >("tfIndex");
}
produces<edm::ValueMap<int> >("isolated");
produces<edm::ValueMap<reco::CandidatePtr> >();
produces<edm::ValueMap<reco::CandidatePtr> >("l1ToReco");
Expand All @@ -95,23 +114,50 @@ pat::L1MuonMatcher::produce(edm::Event & iEvent, const edm::EventSetup & iSetup)

Handle<View<reco::Candidate> > reco;
Handle<vector<l1extra::L1MuonParticle> > l1s;

Handle<l1t::MuonBxCollection> l1tBX;

std::vector<l1t::Muon> l1ts;
std::vector<size_t> bxIdxs;

int minBxIdx = 0;
size_t l1size = 0;

iEvent.getByToken(recoToken_, reco);
iEvent.getByToken(l1Token_, l1s);

if (useStage2L1_) {
iEvent.getByToken(l1tToken_, l1tBX);
l1size = l1tBX->size();

int minBX = max(firstBX_,l1tBX->getFirstBX());
int maxBX = min(lastBX_,l1tBX->getLastBX());

minBxIdx = l1tBX->begin(minBX) - l1tBX->begin();
std::copy(l1tBX->begin(minBX), l1tBX->end(maxBX), std::back_inserter(l1ts));

for (int ibx = l1tBX->getFirstBX(); ibx <= l1tBX->getLastBX(); ++ibx) {
bxIdxs.push_back(l1tBX->end(ibx) - l1tBX->begin());
}
} else {
iEvent.getByToken(l1Token_, l1s);
l1size = l1s->size();
}

unique_ptr<PATPrimitiveCollection> propOut(new PATPrimitiveCollection());
unique_ptr<PATPrimitiveCollection> l1Out(new PATPrimitiveCollection());
std::vector<edm::Ptr<reco::Candidate> > l1rawMatches(reco->size());
vector<int> isSelected(l1s->size(), -1);
std::vector<edm::Ptr<reco::Candidate> > whichRecoMatch(l1s->size());
vector<int> isSelected(l1size, -1);
std::vector<edm::Ptr<reco::Candidate> > whichRecoMatch(l1size);
vector<int> propMatches(reco->size(), -1);
vector<int> fullMatches(reco->size(), -1);
vector<float> deltaRs(reco->size(), 999), deltaPhis(reco->size(), 999);
vector<int> quality(reco->size(), 0), bx(reco->size(), -999), isolated(reco->size(), -999);
vector<int> iPhi(reco->size(), 0), tfIndex(reco->size(), -999);
for (int i = 0, n = reco->size(); i < n; ++i) {
TrajectoryStateOnSurface propagated;
const reco::Candidate &mu = (*reco)[i];
int match = matcher_.match(mu, *l1s, deltaRs[i], deltaPhis[i], propagated);
int match = useStage2L1_ ?
matcher_.match(mu, l1ts, deltaRs[i], deltaPhis[i], propagated) :
matcher_.match(mu, *l1s, deltaRs[i], deltaPhis[i], propagated);
if (propagated.isValid()) {
GlobalPoint pos = propagated.globalPosition();
propMatches[i] = propOut->size();
Expand All @@ -120,21 +166,53 @@ pat::L1MuonMatcher::produce(edm::Event & iEvent, const edm::EventSetup & iSetup)
propOut->back().setCharge(mu.charge());
}
if (match != -1) {
const l1extra::L1MuonParticle & l1 = (*l1s)[match];
whichRecoMatch[match] = reco->ptrAt(i);

if(useStage2L1_) {
match += minBxIdx;
}

whichRecoMatch[match] = reco->ptrAt(i);

int charge = 0;
math::PtEtaPhiMLorentzVector p4;

if (useStage2L1_) {
const l1t::Muon & l1t = (*l1tBX)[match];
charge = l1t.charge();
p4 = l1t.polarP4();
}
else {
const l1extra::L1MuonParticle & l1 = (*l1s)[match];
charge = l1.charge();
p4 = l1.polarP4();
}

if (isSelected[match] == -1) { // copy to output if needed
isSelected[match] = l1Out->size();
l1Out->push_back(PATPrimitive(l1.polarP4()));
l1Out->push_back(PATPrimitive(p4));
l1Out->back().addFilterLabel(labelL1_);
l1Out->back().setCharge(l1.charge());
l1Out->back().setCharge(charge);
}

fullMatches[i] = isSelected[match]; // index in the output collection
const L1MuGMTCand & gmt = l1.gmtMuonCand();
quality[i] = gmt.quality();
bx[i] = gmt.bx();
isolated[i] = gmt.isol();
l1rawMatches[i] = edm::Ptr<reco::Candidate>(l1s, size_t(match));
}

if (useStage2L1_) {
const l1t::Muon & l1t = (*l1tBX)[match];
quality[i] = l1t.hwQual();
bx[i] = l1tBX->getFirstBX() + (std::upper_bound(bxIdxs.begin(),bxIdxs.end(), match) - bxIdxs.begin());
isolated[i] = l1t.hwIso();
l1rawMatches[i] = edm::Ptr<reco::Candidate>(l1tBX, size_t(match));
iPhi[i] = l1t.hwPhi();
tfIndex[i] = l1t.tfMuonIndex();
}
else {
const L1MuGMTCand & gmt = (*l1s)[match].gmtMuonCand();
quality[i] = gmt.quality();
bx[i] = gmt.bx();
isolated[i] = gmt.isol();
l1rawMatches[i] = edm::Ptr<reco::Candidate>(l1s, size_t(match));
}
}
}

OrphanHandle<PATPrimitiveCollection> l1Done = iEvent.put(std::move(l1Out), "l1muons");
Expand All @@ -158,9 +236,17 @@ pat::L1MuonMatcher::produce(edm::Event & iEvent, const edm::EventSetup & iSetup)
storeExtraInfo(iEvent, reco, bx, "bx");
storeExtraInfo(iEvent, reco, isolated, "isolated");
storeExtraInfo(iEvent, reco, quality, "quality");
storeExtraInfo(iEvent, reco, l1rawMatches, "");
storeExtraInfo(iEvent, l1s, whichRecoMatch, "l1ToReco");
storeExtraInfo(iEvent, reco, l1rawMatches, "");
if (useStage2L1_) {
storeExtraInfo(iEvent, l1tBX, whichRecoMatch, "l1ToReco");
storeExtraInfo(iEvent, reco, tfIndex, "tfIndex");
storeExtraInfo(iEvent, reco, iPhi, "iPhi");
}
else {
storeExtraInfo(iEvent, l1s, whichRecoMatch, "l1ToReco");
}
}

}

template<typename Hand, typename T>
Expand Down
Loading

0 comments on commit 94c279c

Please sign in to comment.