CMS 3D CMS Logo

L1MuonMatcherAlgo.cc
Go to the documentation of this file.
2 
4 
6  : prop_(iConfig),
7  preselectionCut_(iConfig.existsAs<std::string>("preselection") ? iConfig.getParameter<std::string>("preselection")
8  : ""),
9  deltaR2_(std::pow(iConfig.getParameter<double>("maxDeltaR"), 2)),
10  deltaPhi_(iConfig.existsAs<double>("maxDeltaPhi") ? iConfig.getParameter<double>("maxDeltaPhi") : 10),
11  sortByDeltaPhi_(iConfig.existsAs<bool>("sortByDeltaPhi") ? iConfig.getParameter<bool>("sortByDeltaPhi") : false) {
12 }
13 
15 
16 void L1MuonMatcherAlgo::init(const edm::EventSetup &iSetup) { prop_.init(iSetup); }
17 
19  const l1extra::L1MuonParticle &l1,
20  float &deltaR,
21  float &deltaPhi) const {
22  if (preselectionCut_(l1)) {
23  GlobalPoint pos = propagated.globalPosition();
24  double thisDeltaPhi = ::deltaPhi(double(pos.phi()), l1.phi());
25  double thisDeltaR2 = ::deltaR2(double(pos.eta()), double(pos.phi()), l1.eta(), l1.phi());
26  if ((fabs(thisDeltaPhi) < deltaPhi_) && (thisDeltaR2 < deltaR2_)) {
27  deltaR = std::sqrt(thisDeltaR2);
28  deltaPhi = thisDeltaPhi;
29  return true;
30  }
31  }
32  return false;
33 }
34 
36  const std::vector<l1extra::L1MuonParticle> &l1s,
37  float &deltaR,
38  float &deltaPhi) const {
39  return matchGeneric(propagated, l1s, preselectionCut_, deltaR, deltaPhi);
40  /*
41  int match = -1;
42  double minDeltaPhi = deltaPhi_;
43  double minDeltaR2 = deltaR2_;
44  GlobalPoint pos = propagated.globalPosition();
45  for (int i = 0, n = l1s.size(); i < n; ++i) {
46  const l1extra::L1MuonParticle &l1 = l1s[i];
47  if (preselectionCut_(l1)) {
48  double thisDeltaPhi = ::deltaPhi(double(pos.phi()), l1.phi());
49  double thisDeltaR2 = ::deltaR2(double(pos.eta()),
50  double(pos.phi()), l1.eta(), l1.phi()); if ((fabs(thisDeltaPhi) <
51  deltaPhi_) && (thisDeltaR2 < deltaR2_)) { // check both if (sortByDeltaPhi_
52  ? (fabs(thisDeltaPhi) < fabs(minDeltaPhi)) : (thisDeltaR2 < minDeltaR2)) {
53  // sort on one match = i; deltaR = std::sqrt(thisDeltaR2); deltaPhi =
54  thisDeltaPhi; if (sortByDeltaPhi_) minDeltaPhi = thisDeltaPhi; else
55  minDeltaR2 = thisDeltaR2;
56  }
57  }
58  }
59  }
60  return match;
61  */
62 }
L1Selector preselectionCut_
Preselection cut to apply to L1 candidates before matching.
Definition: L1Scalers.h:17
void init(const edm::EventSetup &iSetup)
double eta() const final
momentum pseudorapidity
PropagateToMuon prop_
Geom::Phi< T > phi() const
Definition: PV3DBase.h:69
GlobalPoint globalPosition() const
T sqrt(T t)
Definition: SSEVec.h:18
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
void init(const edm::EventSetup &iSetup)
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)
double phi() const final
momentum azimuthal angle
double deltaR2_
Matching cuts.
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40