CMS 3D CMS Logo

L1MuonMatcherAlgo.h
Go to the documentation of this file.
1 #ifndef MuonAnalysis_MuonAssociators_interface_L1MuonMatcherAlgo_h
2 #define MuonAnalysis_MuonAssociators_interface_L1MuonMatcherAlgo_h
3 //
4 //
5 
15 #include <cmath>
16 #include <type_traits>
28 
29 class L1MuonMatcherAlgo {
30  public:
31  explicit L1MuonMatcherAlgo(const edm::ParameterSet & iConfig) ;
33 
35  void init(const edm::EventSetup &iSetup) ;
36 
39 
42 
45 
48 
52  const PropagateToMuon & propagatorToMuon() const { return prop_; }
53 
56  bool match(const reco::Track &tk, const l1extra::L1MuonParticle &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const {
57  propagated = extrapolate(tk);
58  return propagated.isValid() ? match(propagated, l1, deltaR, deltaPhi) : false;
59  }
60 
63  bool match(const reco::Candidate &c, const l1extra::L1MuonParticle &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const {
64  propagated = extrapolate(c);
65  return propagated.isValid() ? match(propagated, l1, deltaR, deltaPhi) : false;
66  }
67 
70  bool match(const SimTrack &tk, const edm::SimVertexContainer &vtxs, const l1extra::L1MuonParticle &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const {
71  propagated = extrapolate(tk, vtxs);
72  return propagated.isValid() ? match(propagated, l1, deltaR, deltaPhi) : false;
73  }
74 
77  bool match(TrajectoryStateOnSurface & propagated, const l1extra::L1MuonParticle &l1, float &deltaR, float &deltaPhi) const ;
78 
79 
80 
81  // Methods to match with vectors of legacy L1 muon trigger objects
82 
86  int match(const reco::Track &tk, const std::vector<l1extra::L1MuonParticle> &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const {
87  propagated = extrapolate(tk);
88  return propagated.isValid() ? match(propagated, l1, deltaR, deltaPhi) : -1;
89  }
90 
94  int match(const reco::Candidate &c, const std::vector<l1extra::L1MuonParticle> &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const {
95  propagated = extrapolate(c);
96  return propagated.isValid() ? match(propagated, l1, deltaR, deltaPhi) : -1;
97  }
98 
102  int match(const SimTrack &tk, const edm::SimVertexContainer &vtxs, const std::vector<l1extra::L1MuonParticle> &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const {
103  propagated = extrapolate(tk, vtxs);
104  return propagated.isValid() ? match(propagated, l1, deltaR, deltaPhi) : -1;
105  }
106 
110  int match(TrajectoryStateOnSurface &propagated, const std::vector<l1extra::L1MuonParticle> &l1, float &deltaR, float &deltaPhi) const ;
111 
112 
113 
114 
115  // Methods to match with vectors of stage2 L1 muon trigger objects
116 
120  int match(const reco::Track &tk, const std::vector<l1t::Muon> &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const {
121  propagated = extrapolate(tk);
122  return propagated.isValid() ? match(propagated, l1, deltaR, deltaPhi) : -1;
123  }
124 
128  int match(const reco::Candidate &c, const std::vector<l1t::Muon> &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const {
129  propagated = extrapolate(c);
130  return propagated.isValid() ? match(propagated, l1, deltaR, deltaPhi) : -1;
131  }
132 
136  int match(const SimTrack &tk, const edm::SimVertexContainer &vtxs, const std::vector<l1t::Muon> &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const {
137  propagated = extrapolate(tk, vtxs);
138  return propagated.isValid() ? match(propagated, l1, deltaR, deltaPhi) : -1;
139  }
140 
144  int match(TrajectoryStateOnSurface &propagated, const std::vector<l1t::Muon> &l1, float &deltaR, float &deltaPhi) const ;
145 
146 
147 
148  // Generic matching
149 
154  template<typename Collection, typename Selector>
155  int matchGeneric(const reco::Track &tk, const Collection &l1, const Selector &sel, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const {
156  propagated = extrapolate(tk);
157  return propagated.isValid() ? matchGeneric(propagated, l1, sel, deltaR, deltaPhi) : -1;
158  }
159 
164  template<typename Collection, typename Selector>
165  int matchGeneric(const reco::Candidate &c, const Collection &l1, const Selector &sel, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const {
166  propagated = extrapolate(c);
167  return propagated.isValid() ? matchGeneric(propagated, l1, sel, deltaR, deltaPhi) : -1;
168  }
169 
174  template<typename Collection, typename Selector>
175  int matchGeneric(TrajectoryStateOnSurface &propagated, const Collection &l1, const Selector &sel, float &deltaR, float &deltaPhi) const ;
176 
178  void setL1PhiOffset(double l1PhiOffset) { l1PhiOffset_ = l1PhiOffset; }
179 
180  private:
181 
182  template<class T> int genericQuality(const T & cand) const;
183 
185 
187 
190  L1Selector preselectionCut_;
191 
194 
198 
200  double l1PhiOffset_;
201 
202 };
203 
204 template<>
205 inline int L1MuonMatcherAlgo::genericQuality<l1extra::L1MuonParticle>(const l1extra::L1MuonParticle & cand) const
206 {
207  return cand.gmtMuonCand().rank();
208 }
209 
210 template<>
211 inline int L1MuonMatcherAlgo::genericQuality<l1t::Muon>(const l1t::Muon & cand) const
212 {
213  return cand.hwQual();
214 }
215 
216 template<class T>
217 inline int L1MuonMatcherAlgo::genericQuality(const T & cand) const
218 {
219  return 0;
220 }
221 
222 
223 
224 template<typename Collection, typename Selector>
225 int
226 L1MuonMatcherAlgo::matchGeneric(TrajectoryStateOnSurface &propagated, const Collection &l1s, const Selector &sel, float &deltaR, float &deltaPhi) const {
227  typedef typename Collection::value_type obj;
228 
229  int match = -1;
230  double minDeltaPhi = deltaPhi_;
231  double minDeltaEta = deltaEta_;
232  double minDeltaR2 = deltaR2_;
233  double maxPt = -1.0;
234  int maxQual = 0;
235  GlobalPoint pos = propagated.globalPosition();
236  for (int i = 0, n = l1s.size(); i < n; ++i) {
237  const obj &l1 = l1s[i];
238  if (sel(l1)) {
239  double thisDeltaPhi = ::deltaPhi(double(pos.phi()), l1.phi()+l1PhiOffset_);
240  double thisDeltaEta = pos.eta() - l1.eta();
241  double thisDeltaR2 = ::deltaR2(double(pos.eta()), double(pos.phi()), l1.eta(), l1.phi()+l1PhiOffset_);
242  int thisQual = genericQuality<obj>(l1);
243  double thisPt = l1.pt();
244  if ((fabs(thisDeltaPhi) < deltaPhi_) && (fabs(thisDeltaEta) < deltaEta_) && (thisDeltaR2 < deltaR2_)) { // check all
245  bool betterMatch = (match == -1);
246  switch (sortBy_) {
247  case SortByDeltaR: betterMatch = (thisDeltaR2 < minDeltaR2); break;
248  case SortByDeltaEta: betterMatch = (fabs(thisDeltaEta) < fabs(minDeltaEta)); break;
249  case SortByDeltaPhi: betterMatch = (fabs(thisDeltaPhi) < fabs(minDeltaPhi)); break;
250  case SortByPt: betterMatch = (thisPt > maxPt); break;
251  // Quality is an int, adding sorting by pT in case of identical qualities
252  case SortByQual: betterMatch = (thisQual > maxQual || ((thisQual == maxQual) && (thisPt > maxPt))); break;
253  }
254  if (betterMatch) { // sort on one
255  match = i;
256  deltaR = std::sqrt(thisDeltaR2);
257  deltaPhi = thisDeltaPhi;
258  minDeltaR2 = thisDeltaR2;
259  minDeltaEta = thisDeltaEta;
260  minDeltaPhi = thisDeltaPhi;
261  maxQual = thisQual;
262  maxPt = thisPt;
263  }
264  }
265  }
266  }
267  return match;
268 }
269 
270 #endif
L1Selector preselectionCut_
Preselection cut to apply to L1 candidates before matching.
Definition: L1Scalers.h:17
TrajectoryStateOnSurface extrapolate(const FreeTrajectoryState &state) const
Extrapolate a FreeTrajectoryState to the muon station 2, return an invalid TSOS if it fails...
void init(const edm::EventSetup &iSetup)
TrajectoryStateOnSurface extrapolate(const reco::Track &tk) const
int matchGeneric(const reco::Candidate &c, const Collection &l1, const Selector &sel, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const
int match(const SimTrack &tk, const edm::SimVertexContainer &vtxs, const std::vector< l1extra::L1MuonParticle > &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const
StringCutObjectSelector< reco::Candidate, true > L1Selector
PropagateToMuon prop_
TrajectoryStateOnSurface extrapolate(const reco::Track &tk) const
Extrapolate reco::Track to the muon station 2, return an invalid TSOS if it fails.
Geom::Phi< T > phi() const
Definition: PV3DBase.h:69
TrajectoryStateOnSurface extrapolate(const SimTrack &tk, const edm::SimVertexContainer &vtx) const
Extrapolate a SimTrack to the muon station 2, return an invalid TSOS if it fails. Requires SimVertice...
GlobalPoint globalPosition() const
const PropagateToMuon & propagatorToMuon() const
Return the propagator to second muon station (in case it&#39;s needed)
bool match(const SimTrack &tk, const edm::SimVertexContainer &vtxs, const l1extra::L1MuonParticle &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const
int match(const SimTrack &tk, const edm::SimVertexContainer &vtxs, const std::vector< l1t::Muon > &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const
int match(const reco::Candidate &c, const std::vector< l1t::Muon > &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const
Propagate an object (usually a track) to the second muon station. Support for other muon stations wil...
int match(const reco::Track &tk, const std::vector< l1extra::L1MuonParticle > &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const
int match(const reco::Candidate &c, const std::vector< l1extra::L1MuonParticle > &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const
T sqrt(T t)
Definition: SSEVec.h:18
int genericQuality(const T &cand) const
Functor that operates on <T>
Definition: Selector.h:23
PropagateToMuon & propagatorToMuon()
Return the propagator to second muon station (in case it&#39;s needed)
SortBy
Sort by deltaPhi or deltaEta instead of deltaR.
Definition: Muon.h:21
bool match(const reco::Candidate &c, const l1extra::L1MuonParticle &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const
void setL1PhiOffset(double l1PhiOffset)
Add this offset to the L1 phi before doing the match, to correct for different scales in L1 vs offlin...
constexpr auto deltaR2(const T1 &t1, const T2 &t2) -> decltype(t1.eta())
Definition: deltaR.h:16
bool match(const reco::Track &tk, const l1extra::L1MuonParticle &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const
std::vector< SimVertex > SimVertexContainer
int matchGeneric(const reco::Track &tk, const Collection &l1, const Selector &sel, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const
T eta() const
Definition: PV3DBase.h:76
L1MuonMatcherAlgo(const edm::ParameterSet &iConfig)
int match(const reco::Track &tk, const std::vector< l1t::Muon > &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const
Matcher of reconstructed objects to L1 Muons.
TrajectoryStateOnSurface extrapolate(const reco::Candidate &tk) const
Extrapolate reco::Candidate to the muon station 2, return an invalid TSOS if it fails.
double l1PhiOffset_
offset to be added to the L1 phi before the match
long double T
double deltaR2_
Matching cuts.