CMS 3D CMS Logo

List of all members | Public Member Functions | Private Attributes
MuonSegmentMatcher Class Reference

#include <MuonSegmentMatcher.h>

Public Member Functions

std::vector< const CSCSegment * > matchCSC (const reco::Track &muon, const edm::Event &event)
 
std::vector< const DTRecSegment4D * > matchDT (const reco::Track &muon, const edm::Event &event)
 perform the matching More...
 
std::vector< const RPCRecHit * > matchRPC (const reco::Track &muon, const edm::Event &event)
 
 MuonSegmentMatcher (const edm::ParameterSet &, edm::ConsumesCollector &iC)
 constructor with Parameter Set and MuonServiceProxy More...
 
virtual ~MuonSegmentMatcher ()
 destructor More...
 

Private Attributes

edm::EDGetTokenT< CSCSegmentCollectionallSegmentsCSCToken
 
edm::InputTag CSCSegmentTags_
 
bool cscTightMatch
 
double dtRadius_
 
edm::EDGetTokenT< DTRecSegment4DCollectiondtRecHitsToken
 
edm::InputTag DTSegmentTags_
 
bool dtTightMatch
 
edm::InputTag RPCHitTags_
 
edm::EDGetTokenT< RPCRecHitCollectionrpcRecHitsToken
 
const edm::EventtheEvent
 
const MuonServiceProxytheService
 
edm::InputTag TKtrackTags_
 
edm::InputTag trackTags_
 

Detailed Description

Definition at line 29 of file MuonSegmentMatcher.h.

Constructor & Destructor Documentation

◆ MuonSegmentMatcher()

MuonSegmentMatcher::MuonSegmentMatcher ( const edm::ParameterSet matchParameters,
edm::ConsumesCollector iC 
)

constructor with Parameter Set and MuonServiceProxy

Definition at line 41 of file MuonSegmentMatcher.cc.

References allSegmentsCSCToken, edm::ConsumesCollector::consumes(), CSCSegmentTags_, dtRecHitsToken, DTSegmentTags_, edm::ParameterSet::existsAs(), edm::ParameterSet::getParameter(), ProducerED_cfi::InputTag, RPCHitTags_, and rpcRecHitsToken.

42  : DTSegmentTags_(matchParameters.getParameter<edm::InputTag>("DTsegments")),
43  CSCSegmentTags_(matchParameters.getParameter<edm::InputTag>("CSCsegments")),
44  dtRadius_(matchParameters.getParameter<double>("DTradius")),
45  dtTightMatch(matchParameters.getParameter<bool>("TightMatchDT")),
46  cscTightMatch(matchParameters.getParameter<bool>("TightMatchCSC")) {
47  if (matchParameters.existsAs<edm::InputTag>("RPChits")) {
48  RPCHitTags_ = matchParameters.getParameter<edm::InputTag>("RPChits");
49  } else
50  RPCHitTags_ = edm::InputTag("hltRpcRecHits");
51 
55 }
EDGetTokenT< ProductType > consumes(edm::InputTag const &tag)
T getParameter(std::string const &) const
Definition: ParameterSet.h:307
edm::EDGetTokenT< RPCRecHitCollection > rpcRecHitsToken
bool existsAs(std::string const &parameterName, bool trackiness=true) const
checks if a parameter exists as a given type
Definition: ParameterSet.h:172
edm::InputTag RPCHitTags_
edm::EDGetTokenT< CSCSegmentCollection > allSegmentsCSCToken
edm::InputTag DTSegmentTags_
edm::InputTag CSCSegmentTags_
edm::EDGetTokenT< DTRecSegment4DCollection > dtRecHitsToken

◆ ~MuonSegmentMatcher()

MuonSegmentMatcher::~MuonSegmentMatcher ( )
virtual

destructor

Definition at line 57 of file MuonSegmentMatcher.cc.

57 {}

Member Function Documentation

◆ matchCSC()

vector< const CSCSegment * > MuonSegmentMatcher::matchCSC ( const reco::Track muon,
const edm::Event event 
)

Definition at line 236 of file MuonSegmentMatcher.cc.

References allSegmentsCSCToken, CSCDetId::chamber(), MuonSubdetId::CSC, cscTightMatch, CSCDetId::endcap(), DetId::Muon, nano_mu_digi_cff::rawId, DetId::rawId(), CSCDetId::ring(), CSCDetId::station(), PV3DBase< T, PVType, FrameType >::x(), and PV3DBase< T, PVType, FrameType >::y().

Referenced by CSCHaloAlgo::Calculate(), and CSCTimingExtractor::fillTiming().

236  {
237  using namespace edm;
238 
239  edm::Handle<CSCSegmentCollection> allSegmentsCSC;
240  event.getByToken(allSegmentsCSCToken, allSegmentsCSC);
241 
242  vector<const CSCSegment*> pointerToCSCSegments;
243 
244  double matchRatioCSC = 0;
245  double CSCXCut = 0.001;
246  double CSCYCut = 0.001;
247  double countMuonCSCHits = 0;
248 
249  for (CSCSegmentCollection::const_iterator segmentCSC = allSegmentsCSC->begin(); segmentCSC != allSegmentsCSC->end();
250  segmentCSC++) {
251  double CSCcountAgreeingHits = 0;
252 
253  if (!segmentCSC->isValid())
254  continue;
255 
256  const vector<CSCRecHit2D>& CSCRechits2D = segmentCSC->specificRecHits();
257  countMuonCSCHits = 0;
258  CSCDetId myChamber((*segmentCSC).geographicalId().rawId());
259 
260  bool segments = false;
261 
262  for (auto const& hitC : muon.recHits()) {
263  if (!hitC->isValid())
264  continue;
265  if (hitC->geographicalId().det() != DetId::Muon)
266  continue;
267  if (hitC->geographicalId().subdetId() != MuonSubdetId::CSC)
268  continue;
269  if (!hitC->isValid())
270  continue;
271  if (hitC->recHits().size() > 1)
272  segments = true;
273 
274  //DETECTOR CONSTRUCTION
275  DetId id = hitC->geographicalId();
276  CSCDetId cscDetIdHit(id.rawId());
277 
278  if (segments) {
279  if (!(myChamber.rawId() == cscDetIdHit.rawId()))
280  continue;
281 
282  // and compare the local positions
283  LocalPoint positionLocalCSC = hitC->localPosition();
284  LocalPoint segLocalCSC = segmentCSC->localPosition();
285  if ((fabs(positionLocalCSC.x() - segLocalCSC.x()) < CSCXCut) &&
286  (fabs(positionLocalCSC.y() - segLocalCSC.y()) < CSCYCut))
287  pointerToCSCSegments.push_back(&(*segmentCSC));
288  continue;
289  }
290 
291  if (!(cscDetIdHit.ring() == myChamber.ring()))
292  continue;
293  if (!(cscDetIdHit.station() == myChamber.station()))
294  continue;
295  if (!(cscDetIdHit.endcap() == myChamber.endcap()))
296  continue;
297  if (!(cscDetIdHit.chamber() == myChamber.chamber()))
298  continue;
299 
300  countMuonCSCHits++;
301 
302  LocalPoint positionLocalCSC = hitC->localPosition();
303 
304  for (vector<CSCRecHit2D>::const_iterator hiti = CSCRechits2D.begin(); hiti != CSCRechits2D.end(); hiti++) {
305  if (!hiti->isValid())
306  continue;
307  CSCDetId cscDetId((hiti->geographicalId()).rawId());
308 
309  if (hitC->geographicalId().rawId() != (hiti->geographicalId()).rawId())
310  continue;
311 
312  LocalPoint segLocalCSC = hiti->localPosition();
313  // cout<<"Layer Id (MuonHit) = "<<cscDetIdHit<<" Muon Local Position (det frame) "<<positionLocalCSC <<endl;
314  // cout<<"Layer Id (CSCHit) = "<<cscDetId<<" Hit Local Position (det frame) "<<segLocalCSC <<endl;
315  if ((fabs(positionLocalCSC.x() - segLocalCSC.x()) < CSCXCut) &&
316  (fabs(positionLocalCSC.y() - segLocalCSC.y()) < CSCYCut)) {
317  CSCcountAgreeingHits++;
318  // cout << " Matched." << endl;
319  }
320  } //End 2D rechit iteration
321  } //End muon hit iteration
322 
323  matchRatioCSC = countMuonCSCHits == 0 ? 0 : CSCcountAgreeingHits / countMuonCSCHits;
324 
325  if ((matchRatioCSC > 0.9) && ((countMuonCSCHits > 1) || !cscTightMatch))
326  pointerToCSCSegments.push_back(&(*segmentCSC));
327 
328  } //End CSC Segment Iteration
329 
330  return pointerToCSCSegments;
331 }
C::const_iterator const_iterator
constant access iterator type
Definition: RangeMap.h:43
T x() const
Definition: PV3DBase.h:59
T y() const
Definition: PV3DBase.h:60
edm::EDGetTokenT< CSCSegmentCollection > allSegmentsCSCToken
Definition: DetId.h:17
HLT enums.
static constexpr int CSC
Definition: MuonSubdetId.h:12

◆ matchDT()

vector< const DTRecSegment4D * > MuonSegmentMatcher::matchDT ( const reco::Track muon,
const edm::Event event 
)

perform the matching

Definition at line 60 of file MuonSegmentMatcher.cc.

References MuonSubdetId::DT, dtRadius_, dtRecHitsToken, dtTightMatch, TrackingRecHit::geographicalId(), DetId::Muon, nano_mu_digi_cff::rawId, DetId::rawId(), FastTrackerRecHitMaskProducer_cfi::recHits, DTRecSegment2D::recHits(), DTRecSegment2D::specificRecHits(), PV3DBase< T, PVType, FrameType >::x(), and PV3DBase< T, PVType, FrameType >::y().

Referenced by DTTimingExtractor::fillTiming().

60  {
61  using namespace edm;
62 
64  event.getByToken(dtRecHitsToken, dtRecHits);
65 
66  vector<const DTRecSegment4D*> pointerTo4DSegments;
67 
68  std::vector<TrackingRecHit const*> dtHits;
69 
70  bool segments = false;
71 
72  // Loop and select DT recHits
73  for (auto const& hit : muon.recHits()) {
74  if (!hit->isValid())
75  continue;
76  if (hit->geographicalId().det() != DetId::Muon)
77  continue;
78  if (hit->geographicalId().subdetId() != MuonSubdetId::DT)
79  continue;
80  if (!hit->recHits().empty())
81  if ((*hit->recHits().begin())->recHits().size() > 1)
82  segments = true;
83  dtHits.push_back(hit);
84  }
85 
86  // cout << "Muon DT hits found: " << dtHits.size() << " segments " << segments << endl;
87 
88  double PhiCutParameter = dtRadius_;
89  double ZCutParameter = dtRadius_;
90  double matchRatioZ = 0;
91  double matchRatioPhi = 0;
92 
93  for (DTRecSegment4DCollection::const_iterator rechit = dtRecHits->begin(); rechit != dtRecHits->end(); ++rechit) {
94  LocalPoint pointLocal = rechit->localPosition();
95 
96  if (segments) {
97  // Loop over muon recHits
98  for (auto hit = dtHits.begin(); hit != dtHits.end(); ++hit) {
99  // Pick the one in the same DT Chamber as the muon
100  DetId idT = (*hit)->geographicalId();
101  if (!(rechit->geographicalId().rawId() == idT.rawId()))
102  continue;
103 
104  // and compare the local positions
105  LocalPoint segLocal = (*hit)->localPosition();
106  if ((fabs(pointLocal.x() - segLocal.x()) < ZCutParameter) &&
107  (fabs(pointLocal.y() - segLocal.y()) < ZCutParameter))
108  pointerTo4DSegments.push_back(&(*rechit));
109  }
110  continue;
111  }
112 
113  double nhitsPhi = 0;
114  double nhitsZ = 0;
115 
116  if (rechit->hasZed()) {
117  double countMuonDTHits = 0;
118  double countAgreeingHits = 0;
119 
120  const DTRecSegment2D* segmZ;
121  segmZ = dynamic_cast<const DTRecSegment2D*>(rechit->zSegment());
122  nhitsZ = segmZ->recHits().size();
123 
124  const vector<DTRecHit1D> hits1d = segmZ->specificRecHits();
125  DTChamberId chamberSegIdT((segmZ->geographicalId()).rawId());
126 
127  // Loop over muon recHits
128  for (auto hit = dtHits.begin(); hit != dtHits.end(); ++hit) {
129  if (!(*hit)->isValid())
130  continue;
131 
132  DetId idT = (*hit)->geographicalId();
133  DTChamberId dtDetIdHitT(idT.rawId());
134  DTSuperLayerId dtDetLayerIdHitT(idT.rawId());
135 
136  LocalPoint pointLocal = (*hit)->localPosition();
137 
138  if ((chamberSegIdT == dtDetIdHitT) && (dtDetLayerIdHitT.superlayer() == 2))
139  countMuonDTHits++;
140 
141  for (vector<DTRecHit1D>::const_iterator hiti = hits1d.begin(); hiti != hits1d.end(); hiti++) {
142  if (!hiti->isValid())
143  continue;
144 
145  // Pick the one in the same DT Layer as the 1D hit
146  if (!(hiti->geographicalId().rawId() == idT.rawId()))
147  continue;
148 
149  // and compare the local positions
150  LocalPoint segLocal = hiti->localPosition();
151  // cout << "Zed Segment Point = "<<pointLocal<<" Muon Point = "<<segLocal<<" Dist: "
152  // << (fabs(pointLocal.x()-segLocal.x()))+(fabs(pointLocal.y()-segLocal.y()))<< endl;
153  if ((fabs(pointLocal.x() - segLocal.x()) < ZCutParameter) &&
154  (fabs(pointLocal.y() - segLocal.y()) < ZCutParameter))
155  countAgreeingHits++;
156  } //End Segment Hit Iteration
157  } //End Muon Hit Iteration
158 
159  matchRatioZ = countMuonDTHits == 0 ? 0 : countAgreeingHits / countMuonDTHits;
160  if (nhitsZ)
161  if (countAgreeingHits / nhitsZ > matchRatioZ)
162  matchRatioZ = countAgreeingHits / nhitsZ;
163  } //End HasZed Check
164 
165  if (rechit->hasPhi()) {
166  double countMuonDTHits = 0;
167  double countAgreeingHits = 0;
168 
169  //PREPARE PARAMETERS FOR SEGMENT DETECTOR GEOMETRY
170  const DTRecSegment2D* segmPhi;
171  segmPhi = dynamic_cast<const DTRecSegment2D*>(rechit->phiSegment());
172  nhitsPhi = segmPhi->recHits().size();
173 
174  const vector<DTRecHit1D> hits1d = segmPhi->specificRecHits();
175  DTChamberId chamberSegIdT((segmPhi->geographicalId()).rawId());
176 
177  // Loop over muon recHits
178  for (auto hit = dtHits.begin(); hit != dtHits.end(); ++hit) {
179  if (!(*hit)->isValid())
180  continue;
181 
182  DetId idT = (*hit)->geographicalId();
183  DTChamberId dtDetIdHitT(idT.rawId());
184  DTSuperLayerId dtDetLayerIdHitT(idT.rawId());
185 
186  LocalPoint pointLocal =
187  (*hit)
188  ->localPosition(); //Localposition is in DTLayer http://cmslxr.fnal.gov/lxr/source/DataFormats/DTRecHit/interface/DTRecHit1D.h
189 
190  if ((chamberSegIdT == dtDetIdHitT) &&
191  ((dtDetLayerIdHitT.superlayer() == 1) || (dtDetLayerIdHitT.superlayer() == 3)))
192  countMuonDTHits++;
193 
194  for (vector<DTRecHit1D>::const_iterator hiti = hits1d.begin(); hiti != hits1d.end(); hiti++) {
195  if (!hiti->isValid())
196  continue;
197 
198  // Pick the one in the same DT Layer as the 1D hit
199  if (!(hiti->geographicalId().rawId() == idT.rawId()))
200  continue;
201 
202  // and compare the local positions
203  LocalPoint segLocal = hiti->localPosition();
204  // cout << " Phi Segment Point = "<<pointLocal<<" Muon Point = "<<segLocal<<" Dist: "
205  // << (fabs(pointLocal.x()-segLocal.x()))+(fabs(pointLocal.y()-segLocal.y()))<< endl;
206 
207  if ((fabs(pointLocal.x() - segLocal.x()) < PhiCutParameter) &&
208  (fabs(pointLocal.y() - segLocal.y()) < PhiCutParameter))
209  countAgreeingHits++;
210  } // End Segment Hit Iteration
211  } // End Muon Hit Iteration
212 
213  matchRatioPhi = countMuonDTHits != 0 ? countAgreeingHits / countMuonDTHits : 0;
214  if (nhitsPhi)
215  if (countAgreeingHits / nhitsPhi > matchRatioPhi)
216  matchRatioPhi = countAgreeingHits / nhitsPhi;
217  } // End HasPhi Check
218  // DTChamberId chamberSegId2((rechit->geographicalId()).rawId());
219  if (dtTightMatch && nhitsPhi && nhitsZ) {
220  if ((matchRatioPhi > 0.9) && (matchRatioZ > 0.9)) {
221  // cout<<"Making a tight match in Chamber "<<chamberSegId2<<endl;
222  pointerTo4DSegments.push_back(&(*rechit));
223  }
224  } else {
225  if ((matchRatioPhi > 0.9 && nhitsPhi) || (matchRatioZ > 0.9 && nhitsZ)) {
226  // cout<<"Making a loose match in Chamber "<<chamberSegId2<<endl;
227  pointerTo4DSegments.push_back(&(*rechit));
228  }
229  }
230 
231  } //End DT Segment Iteration
232 
233  return pointerTo4DSegments;
234 }
C::const_iterator const_iterator
constant access iterator type
Definition: RangeMap.h:43
T x() const
Definition: PV3DBase.h:59
T y() const
Definition: PV3DBase.h:60
std::vector< const TrackingRecHit * > recHits() const override
Access to component RecHits (if any)
Definition: DetId.h:17
DetId geographicalId() const
std::vector< DTRecHit1D > specificRecHits() const
Access to specific components.
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:57
HLT enums.
edm::EDGetTokenT< DTRecSegment4DCollection > dtRecHitsToken
static constexpr int DT
Definition: MuonSubdetId.h:11

◆ matchRPC()

vector< const RPCRecHit * > MuonSegmentMatcher::matchRPC ( const reco::Track muon,
const edm::Event event 
)

Definition at line 333 of file MuonSegmentMatcher.cc.

References muonTagProbeFilters_cff::matched, DetId::Muon, nano_mu_digi_cff::rawId, MuonSubdetId::RPC, dtTriggerPhase2PrimitiveDigis_cfi::rpcRecHits, rpcRecHitsToken, and PV3DBase< T, PVType, FrameType >::x().

333  {
334  using namespace edm;
335 
337  event.getByToken(rpcRecHitsToken, rpcRecHits);
338 
339  vector<const RPCRecHit*> pointerToRPCRecHits;
340  double RPCCut = 0.001;
341 
342  for (RPCRecHitCollection::const_iterator hitRPC = rpcRecHits->begin(); hitRPC != rpcRecHits->end(); hitRPC++) {
343  if (!hitRPC->isValid())
344  continue;
345 
346  RPCDetId myChamber((*hitRPC).geographicalId().rawId());
347  LocalPoint posLocalRPC = hitRPC->localPosition();
348  bool matched = false;
349 
350  for (auto const& hitC : muon.recHits()) {
351  if (!hitC->isValid())
352  continue;
353  if (hitC->geographicalId().det() != DetId::Muon)
354  continue;
355  if (hitC->geographicalId().subdetId() != MuonSubdetId::RPC)
356  continue;
357  if (!hitC->isValid())
358  continue;
359 
360  //DETECTOR CONSTRUCTION
361  DetId id = hitC->geographicalId();
362  RPCDetId rpcDetIdHit(id.rawId());
363 
364  if (rpcDetIdHit != myChamber)
365  continue;
366  LocalPoint posLocalMuon = hitC->localPosition();
367 
368  // cout<<"Layer Id (MuonHit) = "<<rpcDetIdHit<<" Muon Local Position (det frame) "<<posLocalMuon <<endl;
369  // cout<<"Layer Id (RPCHit) = "<<myChamber<<" Hit Local Position (det frame) "<<posLocalRPC <<endl;
370  if ((fabs(posLocalMuon.x() - posLocalRPC.x()) < RPCCut)) {
371  matched = true;
372  break;
373  }
374  }
375 
376  if (matched)
377  pointerToRPCRecHits.push_back(&(*hitRPC));
378  }
379 
380  return pointerToRPCRecHits;
381 }
edm::EDGetTokenT< RPCRecHitCollection > rpcRecHitsToken
C::const_iterator const_iterator
constant access iterator type
Definition: RangeMap.h:43
T x() const
Definition: PV3DBase.h:59
Definition: DetId.h:17
static constexpr int RPC
Definition: MuonSubdetId.h:13
HLT enums.

Member Data Documentation

◆ allSegmentsCSCToken

edm::EDGetTokenT<CSCSegmentCollection> MuonSegmentMatcher::allSegmentsCSCToken
private

Definition at line 56 of file MuonSegmentMatcher.h.

Referenced by matchCSC(), and MuonSegmentMatcher().

◆ CSCSegmentTags_

edm::InputTag MuonSegmentMatcher::CSCSegmentTags_
private

Definition at line 52 of file MuonSegmentMatcher.h.

Referenced by MuonSegmentMatcher().

◆ cscTightMatch

bool MuonSegmentMatcher::cscTightMatch
private

Definition at line 62 of file MuonSegmentMatcher.h.

Referenced by matchCSC().

◆ dtRadius_

double MuonSegmentMatcher::dtRadius_
private

Definition at line 59 of file MuonSegmentMatcher.h.

Referenced by matchDT().

◆ dtRecHitsToken

edm::EDGetTokenT<DTRecSegment4DCollection> MuonSegmentMatcher::dtRecHitsToken
private

Definition at line 55 of file MuonSegmentMatcher.h.

Referenced by matchDT(), and MuonSegmentMatcher().

◆ DTSegmentTags_

edm::InputTag MuonSegmentMatcher::DTSegmentTags_
private

Definition at line 51 of file MuonSegmentMatcher.h.

Referenced by MuonSegmentMatcher().

◆ dtTightMatch

bool MuonSegmentMatcher::dtTightMatch
private

Definition at line 61 of file MuonSegmentMatcher.h.

Referenced by matchDT().

◆ RPCHitTags_

edm::InputTag MuonSegmentMatcher::RPCHitTags_
private

Definition at line 53 of file MuonSegmentMatcher.h.

Referenced by MuonSegmentMatcher().

◆ rpcRecHitsToken

edm::EDGetTokenT<RPCRecHitCollection> MuonSegmentMatcher::rpcRecHitsToken
private

Definition at line 57 of file MuonSegmentMatcher.h.

Referenced by matchRPC(), and MuonSegmentMatcher().

◆ theEvent

const edm::Event* MuonSegmentMatcher::theEvent
private

Definition at line 47 of file MuonSegmentMatcher.h.

◆ theService

const MuonServiceProxy* MuonSegmentMatcher::theService
private

Definition at line 46 of file MuonSegmentMatcher.h.

◆ TKtrackTags_

edm::InputTag MuonSegmentMatcher::TKtrackTags_
private

Definition at line 49 of file MuonSegmentMatcher.h.

◆ trackTags_

edm::InputTag MuonSegmentMatcher::trackTags_
private

Definition at line 50 of file MuonSegmentMatcher.h.