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 
24 #include <cmath>
25 
27 public:
28  explicit L1MuonMatcherAlgo(const edm::ParameterSet &iConfig);
30 
33  void init(const edm::EventSetup &iSetup);
34 
38 
42 
46 
50  const PropagateToMuon &propagatorToMuon() const { return prop_; }
51 
55  bool match(const reco::Track &tk,
56  const l1extra::L1MuonParticle &l1,
57  float &deltaR,
58  float &deltaPhi,
59  TrajectoryStateOnSurface &propagated) const {
60  propagated = extrapolate(tk);
61  return propagated.isValid() ? match(propagated, l1, deltaR, deltaPhi) : false;
62  }
63 
67  bool match(const reco::Candidate &c,
68  const l1extra::L1MuonParticle &l1,
69  float &deltaR,
70  float &deltaPhi,
71  TrajectoryStateOnSurface &propagated) const {
72  propagated = extrapolate(c);
73  return propagated.isValid() ? match(propagated, l1, deltaR, deltaPhi) : false;
74  }
75 
79  bool match(TrajectoryStateOnSurface &propagated,
80  const l1extra::L1MuonParticle &l1,
81  float &deltaR,
82  float &deltaPhi) const;
83 
88  int match(const reco::Track &tk,
89  const std::vector<l1extra::L1MuonParticle> &l1,
90  float &deltaR,
91  float &deltaPhi,
92  TrajectoryStateOnSurface &propagated) const {
93  propagated = extrapolate(tk);
94  return propagated.isValid() ? match(propagated, l1, deltaR, deltaPhi) : -1;
95  }
96 
101  int match(const reco::Candidate &c,
102  const std::vector<l1extra::L1MuonParticle> &l1,
103  float &deltaR,
104  float &deltaPhi,
105  TrajectoryStateOnSurface &propagated) const {
106  propagated = extrapolate(c);
107  return propagated.isValid() ? match(propagated, l1, deltaR, deltaPhi) : -1;
108  }
109 
114  int match(TrajectoryStateOnSurface &propagated,
115  const std::vector<l1extra::L1MuonParticle> &l1,
116  float &deltaR,
117  float &deltaPhi) const;
118 
124  template <typename Collection, typename Selector>
125  int matchGeneric(const reco::Track &tk,
126  const Collection &l1,
127  const Selector &sel,
128  float &deltaR,
129  float &deltaPhi,
130  TrajectoryStateOnSurface &propagated) const {
131  propagated = extrapolate(tk);
132  return propagated.isValid() ? matchGeneric(propagated, l1, sel, deltaR, deltaPhi) : -1;
133  }
134 
140  template <typename Collection, typename Selector>
142  const Collection &l1,
143  const Selector &sel,
144  float &deltaR,
145  float &deltaPhi,
146  TrajectoryStateOnSurface &propagated) const {
147  propagated = extrapolate(c);
148  return propagated.isValid() ? matchGeneric(propagated, l1, sel, deltaR, deltaPhi) : -1;
149  }
150 
156  template <typename Collection, typename Selector>
157  int matchGeneric(TrajectoryStateOnSurface &propagated,
158  const Collection &l1,
159  const Selector &sel,
160  float &deltaR,
161  float &deltaPhi) const;
162 
163 private:
165 
168  L1Selector preselectionCut_;
169 
172 
175 };
176 
177 template <typename Collection, typename Selector>
179  const Collection &l1s,
180  const Selector &sel,
181  float &deltaR,
182  float &deltaPhi) const {
183  typedef typename Collection::value_type obj;
184  int match = -1;
185  double minDeltaPhi = deltaPhi_;
186  double minDeltaR2 = deltaR2_;
187  GlobalPoint pos = propagated.globalPosition();
188  for (int i = 0, n = l1s.size(); i < n; ++i) {
189  const obj &l1 = l1s[i];
190  if (sel(l1)) {
191  double thisDeltaPhi = ::deltaPhi(double(pos.phi()), l1.phi());
192  double thisDeltaR2 = ::deltaR2(double(pos.eta()), double(pos.phi()), l1.eta(), l1.phi());
193  if ((fabs(thisDeltaPhi) < deltaPhi_) && (thisDeltaR2 < deltaR2_)) { // check both
194  if (sortByDeltaPhi_ ? (fabs(thisDeltaPhi) < fabs(minDeltaPhi)) : (thisDeltaR2 < minDeltaR2)) { // sort on one
195  match = i;
196  deltaR = std::sqrt(thisDeltaR2);
197  deltaPhi = thisDeltaPhi;
198  if (sortByDeltaPhi_)
199  minDeltaPhi = thisDeltaPhi;
200  else
201  minDeltaR2 = thisDeltaR2;
202  }
203  }
204  }
205  }
206  return match;
207 }
208 
209 #endif
L1Selector preselectionCut_
Preselection cut to apply to L1 candidates before matching.
Definition: L1Scalers.h:17
TrajectoryStateOnSurface extrapolate(const FreeTrajectoryState &state) const
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
PropagateToMuon prop_
TrajectoryStateOnSurface extrapolate(const reco::Track &tk) const
Geom::Phi< T > phi() const
Definition: PV3DBase.h:69
GlobalPoint globalPosition() const
bool sortByDeltaPhi_
Sort by deltaPhi instead of deltaR.
const PropagateToMuon & propagatorToMuon() const
Return the propagator to second muon station (in case it&#39;s needed)
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
Functor that operates on <T>
Definition: Selector.h:23
PropagateToMuon & propagatorToMuon()
Return the propagator to second muon station (in case it&#39;s needed)
bool match(const reco::Candidate &c, const l1extra::L1MuonParticle &l1, float &deltaR, float &deltaPhi, TrajectoryStateOnSurface &propagated) const
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
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
StringCutObjectSelector< l1extra::L1MuonParticle > L1Selector
L1MuonMatcherAlgo(const edm::ParameterSet &iConfig)
Matcher of reconstructed objects to L1 Muons.
TrajectoryStateOnSurface extrapolate(const reco::Candidate &tk) const
double deltaR2_
Matching cuts.