CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
MuonMesh.cc
Go to the documentation of this file.
1 #include "RecoMuon/MuonIdentification/interface/MuonMesh.h"
2 #include "DataFormats/CSCRecHit/interface/CSCSegmentCollection.h"
3 #include "DataFormats/DTRecHit/interface/DTRecSegment4DCollection.h"
4 #include "Geometry/CSCGeometry/interface/CSCGeometry.h"
5 
6 #include <utility>
7 #include <algorithm>
8 
9 MuonMesh::MuonMesh(const edm::ParameterSet& parm)
10  : doME1a(parm.getParameter<bool>("ME1a")),
11  doOverlaps(parm.getParameter<bool>("Overlap")),
12  doClustering(parm.getParameter<bool>("Clustering")),
13  OverlapDPhi(parm.getParameter<double>("OverlapDPhi")),
14  OverlapDTheta(parm.getParameter<double>("OverlapDTheta")),
15  ClusterDPhi(parm.getParameter<double>("ClusterDPhi")),
16  ClusterDTheta(parm.getParameter<double>("ClusterDTheta")) {}
17 
18 void MuonMesh::fillMesh(std::vector<reco::Muon>* inputMuons) {
19  for (std::vector<reco::Muon>::iterator muonIter1 = inputMuons->begin(); muonIter1 != inputMuons->end(); ++muonIter1) {
20  if (muonIter1->isTrackerMuon()) {
21  mesh_[&*muonIter1]; // create this entry if it's a tracker muon
22  for (std::vector<reco::Muon>::iterator muonIter2 = inputMuons->begin(); muonIter2 != inputMuons->end();
23  ++muonIter2) {
24  if (muonIter2->isTrackerMuon()) {
25  if (muonIter2 != muonIter1) {
26  // now fill all the edges for muon1 based on overlaps with muon2
27  for (std::vector<reco::MuonChamberMatch>::iterator chamberIter1 = muonIter1->matches().begin();
28  chamberIter1 != muonIter1->matches().end();
29  ++chamberIter1) {
30  for (std::vector<reco::MuonSegmentMatch>::iterator segmentIter1 = chamberIter1->segmentMatches.begin();
31  segmentIter1 != chamberIter1->segmentMatches.end();
32  ++segmentIter1) {
33  for (std::vector<reco::MuonChamberMatch>::iterator chamberIter2 = muonIter2->matches().begin();
34  chamberIter2 != muonIter2->matches().end();
35  ++chamberIter2) {
36  for (std::vector<reco::MuonSegmentMatch>::iterator segmentIter2 =
37  chamberIter2->segmentMatches.begin();
38  segmentIter2 != chamberIter2->segmentMatches.end();
39  ++segmentIter2) {
40  //if(segmentIter1->mask & 0x1e0000 && segmentIter2->mask &0x1e0000) {
41 
42  bool addsegment(false);
43 
44  if (segmentIter1->cscSegmentRef.isNonnull() && segmentIter2->cscSegmentRef.isNonnull()) {
45  if (doME1a && isDuplicateOf(segmentIter1->cscSegmentRef, segmentIter2->cscSegmentRef) &&
46  CSCDetId(chamberIter1->id).ring() == 4 && CSCDetId(chamberIter2->id).ring() == 4 &&
47  chamberIter1->id == chamberIter2->id) {
48  addsegment = true;
49  //std::cout << "\tME1/a sharing detected." << std::endl;
50  }
51 
52  if (doOverlaps &&
53  isDuplicateOf(std::make_pair(CSCDetId(chamberIter1->id), segmentIter1->cscSegmentRef),
54  std::make_pair(CSCDetId(chamberIter2->id), segmentIter2->cscSegmentRef))) {
55  addsegment = true;
56  //std::cout << "\tChamber Overlap sharing detected." << std::endl;
57  }
58 
59  if (doClustering &&
60  isClusteredWith(std::make_pair(CSCDetId(chamberIter1->id), segmentIter1->cscSegmentRef),
61  std::make_pair(CSCDetId(chamberIter2->id), segmentIter2->cscSegmentRef))) {
62  addsegment = true;
63  //std::cout << "\tCluster sharing detected." << std::endl;
64  }
65  //std::cout << std::endl;
66  } // has valid csc segment ref
67 
68  if (addsegment) { // add segment if clusters/overlaps/replicant and doesn't already exist
69 
70  if (find(mesh_[&*muonIter1].begin(),
71  mesh_[&*muonIter1].end(),
72  std::make_pair(&*muonIter2, std::make_pair(&*chamberIter2, &*segmentIter2))) ==
73  mesh_[&*muonIter1].end()) {
74  mesh_[&*muonIter1].push_back(
75  std::make_pair(&*muonIter2, std::make_pair(&*chamberIter2, &*segmentIter2)));
76  } // find
77  } // add segment?
78  //} // both segments won arbitration
79  } // segmentIter 2
80  } // chamberIter2
81  } //segmentIter1
82  } // chamberIter1
83 
84  } // if different muon
85  } // is tracker muon
86  } // muonIter2
87 
88  } // is tracker muon
89  } // muonIter1
90 
91  // special cases
92 
93  // one muon: mark all segments belonging to a muon as cleaned as there are no other muons to fight with
94  if (mesh_.size() == 1) {
95  for (std::vector<reco::MuonChamberMatch>::iterator chamberIter1 = mesh_.begin()->first->matches().begin();
96  chamberIter1 != mesh_.begin()->first->matches().end();
97  ++chamberIter1) {
98  for (std::vector<reco::MuonSegmentMatch>::iterator segmentIter1 = chamberIter1->segmentMatches.begin();
99  segmentIter1 != chamberIter1->segmentMatches.end();
100  ++segmentIter1) {
101  segmentIter1->setMask(reco::MuonSegmentMatch::BelongsToTrackByCleaning);
102  } // segmentIter1
103  } // chamberIter1
104  } // if only one tracker muon set winner bit boosted arbitration
105 
106  // segments that are not shared amongst muons and the have won all segment arbitration flags need to be promoted
107  // also promote DT segments
108  if (mesh_.size() > 1) {
109  for (MeshType::iterator i = mesh_.begin(); i != mesh_.end(); ++i) {
110  for (std::vector<reco::MuonChamberMatch>::iterator chamberIter1 = i->first->matches().begin();
111  chamberIter1 != i->first->matches().end();
112  ++chamberIter1) {
113  for (std::vector<reco::MuonSegmentMatch>::iterator segmentIter1 = chamberIter1->segmentMatches.begin();
114  segmentIter1 != chamberIter1->segmentMatches.end();
115  ++segmentIter1) {
116  bool shared(false);
117 
118  for (AssociationType::iterator j = i->second.begin(); j != i->second.end(); ++j) {
119  if (segmentIter1->cscSegmentRef.isNonnull() && j->second.second->cscSegmentRef.isNonnull()) {
120  if (chamberIter1->id.subdetId() == MuonSubdetId::CSC &&
121  j->second.first->id.subdetId() == MuonSubdetId::CSC) {
122  CSCDetId segIterId(chamberIter1->id), shareId(j->second.first->id);
123 
124  if (doOverlaps && isDuplicateOf(std::make_pair(segIterId, segmentIter1->cscSegmentRef),
125  std::make_pair(shareId, j->second.second->cscSegmentRef)))
126  shared = true;
127 
128  if (doME1a && isDuplicateOf(segmentIter1->cscSegmentRef, j->second.second->cscSegmentRef) &&
129  segIterId.ring() == 4 && shareId.ring() == 4 && segIterId == segIterId)
130  shared = true;
131 
132  if (doClustering &&
133  isClusteredWith(std::make_pair(CSCDetId(chamberIter1->id), segmentIter1->cscSegmentRef),
134  std::make_pair(CSCDetId(j->second.first->id), j->second.second->cscSegmentRef)))
135  shared = true;
136  } // in CSCs?
137  } // cscSegmentRef non null?
138  } // j
139 
140  // Promote segments which have won all arbitration and are not shared or are DT segments
141  if (((segmentIter1->mask & 0x1e0000) == 0x1e0000 && !shared) ||
142  (chamberIter1->id.subdetId() == MuonSubdetId::DT && (segmentIter1->mask & 0x1e0000)))
143  segmentIter1->setMask(reco::MuonSegmentMatch::BelongsToTrackByCleaning);
144 
145  } // segmentIter1
146  } // chamberIter1
147  } // i
148  } // if non-trivial case
149 }
150 
151 void MuonMesh::pruneMesh() {
152  for (MeshType::iterator i = mesh_.begin(); i != mesh_.end(); ++i) {
153  for (AssociationType::iterator j = i->second.begin(); j != i->second.end(); ++j) {
154  for (std::vector<reco::MuonChamberMatch>::iterator chamberIter1 = i->first->matches().begin();
155  chamberIter1 != i->first->matches().end();
156  ++chamberIter1) {
157  for (std::vector<reco::MuonSegmentMatch>::iterator segmentIter1 = chamberIter1->segmentMatches.begin();
158  segmentIter1 != chamberIter1->segmentMatches.end();
159  ++segmentIter1) {
160  if (j->second.second->cscSegmentRef.isNonnull() && segmentIter1->cscSegmentRef.isNonnull()) {
161  //UNUSED: bool me1a(false), overlap(false), cluster(false);
162 
163  // remove physical overlap duplicates first
164  if (doOverlaps &&
165  isDuplicateOf(std::make_pair(CSCDetId(chamberIter1->id), segmentIter1->cscSegmentRef),
166  std::make_pair(CSCDetId(j->second.first->id), j->second.second->cscSegmentRef))) {
167  if (i->first->numberOfMatches((reco::Muon::ArbitrationType)0x1e0000) >
168  j->first->numberOfMatches((reco::Muon::ArbitrationType)0x1e0000)) {
169  segmentIter1->setMask(reco::MuonSegmentMatch::BelongsToTrackByOvlClean);
170  segmentIter1->setMask(reco::MuonSegmentMatch::BelongsToTrackByCleaning);
171 
172  //UNUSED: overlap = true;
173  } else if (i->first->numberOfMatches((reco::Muon::ArbitrationType)0x1e0000) ==
174  j->first->numberOfMatches(
175  (reco::Muon::ArbitrationType)0x1e0000)) { // muon with more matched stations wins
176 
177  if ((segmentIter1->mask & 0x1e0000) >
178  (j->second.second->mask & 0x1e0000)) { // segment with better match wins
179 
180  segmentIter1->setMask(reco::MuonSegmentMatch::BelongsToTrackByOvlClean);
181  segmentIter1->setMask(reco::MuonSegmentMatch::BelongsToTrackByCleaning);
182 
183  //UNUSED: overlap = true;
184  } else { // ??
185  // leave this available for later
186  }
187  } // overlap duplicate resolution
188  } // is overlap duplicate
189 
190  // do ME1/a arbitration second since the tie breaker depends on other stations
191  // Unlike the other cleanings this one removes the bits from segments associated to tracks which
192  // fail cleaning. (Instead of setting bits for the segments which win.)
193  if (doME1a && isDuplicateOf(segmentIter1->cscSegmentRef, j->second.second->cscSegmentRef) &&
194  CSCDetId(chamberIter1->id).ring() == 4 && CSCDetId(j->second.first->id).ring() == 4 &&
195  chamberIter1->id == j->second.first->id) {
196  if (j->first->numberOfMatches((reco::Muon::ArbitrationType)0x1e0000) <
197  i->first->numberOfMatches((reco::Muon::ArbitrationType)0x1e0000)) {
198  for (AssociationType::iterator AsscIter1 = i->second.begin(); AsscIter1 != i->second.end();
199  ++AsscIter1) {
200  if (AsscIter1->second.second->cscSegmentRef.isNonnull())
201  if (j->first == AsscIter1->first && j->second.first == AsscIter1->second.first &&
202  isDuplicateOf(segmentIter1->cscSegmentRef, AsscIter1->second.second->cscSegmentRef)) {
203  AsscIter1->second.second->mask &= ~reco::MuonSegmentMatch::BelongsToTrackByME1aClean;
204  }
205  }
206 
207  //UNUSED: me1a = true;
208  } else if (j->first->numberOfMatches((reco::Muon::ArbitrationType)0x1e0000) ==
209  i->first->numberOfMatches(
210  (reco::Muon::ArbitrationType)0x1e0000)) { // muon with best arbitration wins
211 
212  bool bestArb(true);
213 
214  for (AssociationType::iterator AsscIter1 = i->second.begin(); AsscIter1 != i->second.end();
215  ++AsscIter1) {
216  if (AsscIter1->second.second->cscSegmentRef.isNonnull())
217  if (j->first == AsscIter1->first && j->second.first == AsscIter1->second.first &&
218  isDuplicateOf(segmentIter1->cscSegmentRef, AsscIter1->second.second->cscSegmentRef) &&
219  (segmentIter1->mask & 0x1e0000) < (AsscIter1->second.second->mask & 0x1e0000))
220  bestArb = false;
221  }
222 
223  if (bestArb) {
224  for (AssociationType::iterator AsscIter1 = i->second.begin(); AsscIter1 != i->second.end();
225  ++AsscIter1) {
226  if (AsscIter1->second.second->cscSegmentRef.isNonnull())
227  if (j->first == AsscIter1->first && j->second.first == AsscIter1->second.first &&
228  isDuplicateOf(segmentIter1->cscSegmentRef, AsscIter1->second.second->cscSegmentRef)) {
229  AsscIter1->second.second->mask &= ~reco::MuonSegmentMatch::BelongsToTrackByME1aClean;
230  }
231  }
232  }
233  //UNUSED me1a = true;
234 
235  } // ME1/a duplicate resolution
236  } // is ME1/aduplicate?
237 
238  if (doClustering &&
239  isClusteredWith(std::make_pair(CSCDetId(chamberIter1->id), segmentIter1->cscSegmentRef),
240  std::make_pair(CSCDetId(j->second.first->id), j->second.second->cscSegmentRef))) {
241  if (i->first->numberOfMatches((reco::Muon::ArbitrationType)0x1e0000) >
242  j->first->numberOfMatches((reco::Muon::ArbitrationType)0x1e0000)) {
243  segmentIter1->setMask(reco::MuonSegmentMatch::BelongsToTrackByClusClean);
244  segmentIter1->setMask(reco::MuonSegmentMatch::BelongsToTrackByCleaning);
245 
246  //UNUSED: cluster = true;
247  } else if (i->first->numberOfMatches((reco::Muon::ArbitrationType)0x1e0000) <
248  j->first->numberOfMatches((reco::Muon::ArbitrationType)0x1e0000)) {
249  j->second.second->setMask(reco::MuonSegmentMatch::BelongsToTrackByClusClean);
250  j->second.second->setMask(reco::MuonSegmentMatch::BelongsToTrackByCleaning);
251 
252  //UNUSED: cluster = true;
253  } else { // muon with more matched stations wins
254 
255  if ((segmentIter1->mask & 0x1e0000) >
256  (j->second.second->mask & 0x1e0000)) { // segment with better match wins
257 
258  segmentIter1->setMask(reco::MuonSegmentMatch::BelongsToTrackByClusClean);
259  segmentIter1->setMask(reco::MuonSegmentMatch::BelongsToTrackByCleaning);
260 
261  //UNUSED: cluster = true;
262  } else if ((segmentIter1->mask & 0x1e0000) < (j->second.second->mask & 0x1e0000)) { //
263 
264  j->second.second->setMask(reco::MuonSegmentMatch::BelongsToTrackByClusClean);
265  j->second.second->setMask(reco::MuonSegmentMatch::BelongsToTrackByCleaning);
266 
267  //UNUSED: cluster = true;
268  } else {
269  }
270  } // cluster sharing resolution
271 
272  } // is clustered with?
273 
274  } // csc ref nonnull
275  } // segmentIter1
276  } // chamberIter1
277  } // j, associated segments iterator
278  } // i, map iterator
279 
280  // final step: make sure everything that's won a cleaning flag has the "BelongsToTrackByCleaning" flag
281 
282  for (MeshType::iterator i = mesh_.begin(); i != mesh_.end(); ++i) {
283  for (std::vector<reco::MuonChamberMatch>::iterator chamberIter1 = i->first->matches().begin();
284  chamberIter1 != i->first->matches().end();
285  ++chamberIter1) {
286  for (std::vector<reco::MuonSegmentMatch>::iterator segmentIter1 = chamberIter1->segmentMatches.begin();
287  segmentIter1 != chamberIter1->segmentMatches.end();
288  ++segmentIter1) {
289  // set cleaning bit if initial no cleaning bit but there are cleaning algorithm bits set.
290  if (!segmentIter1->isMask(reco::MuonSegmentMatch::BelongsToTrackByCleaning) && segmentIter1->isMask(0xe00000))
291  segmentIter1->setMask(reco::MuonSegmentMatch::BelongsToTrackByCleaning);
292  } // segmentIter1
293  } // chamberIter1
294  } // i
295 }
296 
297 bool MuonMesh::isDuplicateOf(const CSCSegmentRef& lhs, const CSCSegmentRef& rhs)
298  const // this isDuplicateOf() deals with duplicate segments in ME1/a
299 {
300  bool result(false);
301 
302  if (!lhs->isME11a_duplicate())
303  return result;
304 
305  std::vector<CSCSegment> lhs_duplicates = lhs->duplicateSegments();
306 
307  if (fabs(lhs->localPosition().x() - rhs->localPosition().x()) < 1E-3 &&
308  fabs(lhs->localPosition().y() - rhs->localPosition().y()) < 1E-3 &&
309  fabs(lhs->localDirection().x() / lhs->localDirection().z() -
310  rhs->localDirection().x() / rhs->localDirection().z()) < 1E-3 &&
311  fabs(lhs->localDirection().y() / lhs->localDirection().z() -
312  rhs->localDirection().y() / rhs->localDirection().z()) < 1E-3 &&
313  fabs(lhs->localPositionError().xx() - rhs->localPositionError().xx()) < 1E-3 &&
314  fabs(lhs->localPositionError().yy() - rhs->localPositionError().yy()) < 1E-3 &&
315  fabs(lhs->localDirectionError().xx() - rhs->localDirectionError().xx()) < 1E-3 &&
316  fabs(lhs->localDirectionError().yy() - rhs->localDirectionError().yy()) < 1E-3)
317  result = true;
318 
319  for (std::vector<CSCSegment>::const_iterator segIter1 = lhs_duplicates.begin(); segIter1 != lhs_duplicates.end();
320  ++segIter1) { // loop over lhs duplicates
321 
322  if (fabs(segIter1->localPosition().x() - rhs->localPosition().x()) < 1E-3 &&
323  fabs(segIter1->localPosition().y() - rhs->localPosition().y()) < 1E-3 &&
324  fabs(segIter1->localDirection().x() / segIter1->localDirection().z() -
325  rhs->localDirection().x() / rhs->localDirection().z()) < 1E-3 &&
326  fabs(segIter1->localDirection().y() / segIter1->localDirection().z() -
327  rhs->localDirection().y() / rhs->localDirection().z()) < 1E-3 &&
328  fabs(segIter1->localPositionError().xx() - rhs->localPositionError().xx()) < 1E-3 &&
329  fabs(segIter1->localPositionError().yy() - rhs->localPositionError().yy()) < 1E-3 &&
330  fabs(segIter1->localDirectionError().xx() - rhs->localDirectionError().xx()) < 1E-3 &&
331  fabs(segIter1->localDirectionError().yy() - rhs->localDirectionError().yy()) < 1E-3)
332  result = true;
333  /*
334  if(fabs(segIter1->localPosition().x() - rhs->localPosition().x() ) < 2*sqrt(segIter1->localPositionError().xx()) &&
335  fabs(segIter1->localPosition().y() - rhs->localPosition().y() ) < 2*sqrt(segIter1->localPositionError().yy()) &&
336  fabs(segIter1->localDirection().x()/segIter1->localDirection().z() - rhs->localDirection().x()/rhs->localDirection().z() )
337  < 2*std::sqrt(std::max(segIter1->localDirectionError().yy(),rhs->localDirectionError().xx())) &&
338  fabs(segIter1->localDirection().y()/segIter1->localDirection().z() - rhs->localDirection().y()/rhs->localDirection().z() )
339  < 2*std::sqrt(std::max(segIter1->localDirectionError().yy(),rhs->localDirectionError().yy())))
340  result = true;
341  */
342 
343  } // loop over duplicates
344 
345  return result;
346 }
347 
348 bool MuonMesh::isDuplicateOf(const std::pair<CSCDetId, CSCSegmentRef>& rhs,
349  const std::pair<CSCDetId, CSCSegmentRef>& lhs)
350  const // this isDuplicateOf() deals with "overlapping chambers" duplicates
351 {
352  bool result(false);
353 
354  // try to implement the simple case first just back-to-back segments without treatment of ME1/a ganging
355  // ME1a should be a simple extension of this
356 
357  if (rhs.first.endcap() == lhs.first.endcap() && rhs.first.station() == lhs.first.station() &&
358  rhs.first.ring() == lhs.first.ring()) { // if same endcap,station,ring (minimal requirement for ovl candidate)
359  /*
360  std::cout << "Chamber 1: " << rhs.first << std::endl
361  << "Chamber 2: " << lhs.first << std::endl;
362 
363  std::cout << "Same endcap,ring,station." << std::endl;
364  */
365  //create neighboring chamber labels, treat ring as (Z mod 36 or 18) + 1 number line: left, right defined accordingly.
366  unsigned modulus = ((rhs.first.ring() != 1 || rhs.first.station() == 1) ? 36 : 18);
367  int left_neighbor = (((rhs.first.chamber() - 1 + modulus) % modulus == 0)
368  ? modulus
369  : (rhs.first.chamber() - 1 + modulus) % modulus); // + modulus to ensure positivity
370  int right_neighbor = (((rhs.first.chamber() + 1) % modulus == 0) ? modulus : (rhs.first.chamber() + 1) % modulus);
371 
372  if (lhs.first.chamber() == left_neighbor ||
373  lhs.first.chamber() == right_neighbor) { // if this is a neighboring chamber then it can be an overlap
374 
375  std::vector<CSCSegment> thesegments;
376  thesegments.push_back(*(lhs.second));
377  /*
378  if(lhs.second->isME11a_duplicate())
379  thesegments.insert(thesegments.begin(),
380  lhs.second->duplicateSegments().begin(),
381  lhs.second->duplicateSegments().end());
382  */
383 
384  //std::cout << "lhs is in neighoring chamber of rhs." << std::endl;
385 
386  // rhs local direction info
387  /*
388  double rhs_dydz = geometry_->chamber(rhs.first)->toGlobal(rhs.second->localDirection()).y()/
389  geometry_->chamber(rhs.first)->toGlobal(rhs.second->localDirection()).z();
390  double rhs_dxdz = geometry_->chamber(rhs.first)->toGlobal(rhs.second->localDirection()).x()/
391  geometry_->chamber(rhs.first)->toGlobal(rhs.second->localDirection()).z();
392  double rhs_dydz_err = rhs.second->localDirectionError().yy();
393  double rhs_dxdz_err = rhs.second->localDirectionError().xx();
394  */
395 
396  //rhs global position info
397  double rhs_phi = geometry_->chamber(rhs.first)->toGlobal(rhs.second->localPosition()).phi();
398  double rhs_theta = geometry_->chamber(rhs.first)->toGlobal(rhs.second->localPosition()).theta();
399  double rhs_z = geometry_->chamber(rhs.first)->toGlobal(rhs.second->localPosition()).z();
400 
401  for (std::vector<CSCSegment>::const_iterator ilhs = thesegments.begin(); ilhs != thesegments.end(); ++ilhs) {
402  // lhs local direction info
403  /*
404  double lhs_dydz = geometry_->chamber(lhs.first)->toGlobal(ilhs->localDirection()).y()/
405  geometry_->chamber(lhs.first)->toGlobal(ilhs->localDirection()).z();
406  double lhs_dxdz = geometry_->chamber(lhs.first)->toGlobal(ilhs->localDirection()).x()/
407  geometry_->chamber(lhs.first)->toGlobal(ilhs->localDirection()).z();
408  double lhs_dydz_err = ilhs->localDirectionError().yy();
409  double lhs_dxdz_err = ilhs->localDirectionError().xx();
410  */
411 
412  //lhs global position info
413  double lhs_phi = geometry_->chamber(lhs.first)->toGlobal(ilhs->localPosition()).phi();
414  double lhs_theta = geometry_->chamber(lhs.first)->toGlobal(ilhs->localPosition()).theta();
415  double lhs_z = geometry_->chamber(lhs.first)->toGlobal(ilhs->localPosition()).z();
416  /*
417  std::cout << "RHS Segment Parameters:" << std::endl
418  << "\t RHS Phi : " << rhs_phi << std::endl
419  << "\t RHS Theta : " << rhs_theta << std::endl
420  << "\t RHS dx/dz : " << rhs_dxdz << " +- " << rhs_dxdz_err << std::endl
421  << "\t RHS dy/dz : " << rhs_dydz << " +- " << rhs_dydz_err << std::endl;
422 
423  std::cout << "LHS Segment Parameters:" << std::endl
424  << "\t LHS Phi : " << lhs_phi << std::endl
425  << "\t LHS Theta : " << lhs_theta << std::endl
426  << "\t LHS dx/dz : " << lhs_dxdz << " +- " << lhs_dxdz_err << std::endl
427  << "\t LHS dy/dz : " << lhs_dydz << " +- " << lhs_dydz_err << std::endl;
428  */
429 
430  double phidiff =
431  (fabs(rhs_phi - lhs_phi) > 2 * M_PI ? fabs(rhs_phi - lhs_phi) - 2 * M_PI : fabs(rhs_phi - lhs_phi));
432 
433  if (phidiff < OverlapDPhi && fabs(rhs_theta - lhs_theta) < OverlapDTheta && fabs(rhs_z) < fabs(lhs_z) &&
434  rhs_z * lhs_z > 0) // phi overlap region is 3.5 degrees and rhs is infront of lhs
435  result = true;
436  } // loop over duplicate segments
437  } // neighboring chamber
438  } // same endcap,station,ring
439 
440  return result;
441 }
442 
443 bool MuonMesh::isClusteredWith(const std::pair<CSCDetId, CSCSegmentRef>& lhs,
444  const std::pair<CSCDetId, CSCSegmentRef>& rhs) const {
445  bool result(false);
446 
447  // try to implement the simple case first just back-to-back segments without treatment of ME1/a ganging
448  // ME1a should be a simple extension of this
449 
450  //std::cout << lhs.first << ' ' << rhs.first << std::endl;
451 
452  if (rhs.first.endcap() == lhs.first.endcap() && lhs.first.station() < rhs.first.station()) {
453  std::vector<CSCSegment> thesegments;
454  thesegments.push_back(*(lhs.second));
455  /*
456  if(lhs.second->isME11a_duplicate())
457  thesegments.insert(thesegments.begin(),
458  lhs.second->duplicateSegments().begin(),
459  lhs.second->duplicateSegments().end());
460  */
461  //std::cout << "lhs is in neighoring chamber of rhs." << std::endl;
462 
463  // rhs local direction info
464  /*
465  double rhs_dydz = geometry_->chamber(rhs.first)->toGlobal(rhs.second->localDirection()).y()/
466  geometry_->chamber(rhs.first)->toGlobal(rhs.second->localDirection()).z();
467  double rhs_dxdz = geometry_->chamber(rhs.first)->toGlobal(rhs.second->localDirection()).x()/
468  geometry_->chamber(rhs.first)->toGlobal(rhs.second->localDirection()).z();
469  double rhs_dydz_err = rhs.second->localDirectionError().yy();
470  double rhs_dxdz_err = rhs.second->localDirectionError().xx();
471  */
472 
473  //rhs global position info
474  double rhs_phi = geometry_->chamber(rhs.first)->toGlobal(rhs.second->localPosition()).phi();
475  double rhs_theta = geometry_->chamber(rhs.first)->toGlobal(rhs.second->localPosition()).theta();
476 
477  for (std::vector<CSCSegment>::const_iterator ilhs = thesegments.begin(); ilhs != thesegments.end(); ++ilhs) {
478  // lhs local direction info
479  /*
480  double lhs_dydz = geometry_->chamber(lhs.first)->toGlobal(ilhs->localDirection()).y()/
481  geometry_->chamber(lhs.first)->toGlobal(ilhs->localDirection()).z();
482  double lhs_dxdz = geometry_->chamber(lhs.first)->toGlobal(ilhs->localDirection()).x()/
483  geometry_->chamber(lhs.first)->toGlobal(ilhs->localDirection()).z();
484  double lhs_dydz_err = ilhs->localDirectionError().yy();
485  double lhs_dxdz_err = ilhs->localDirectionError().xx();
486  */
487 
488  //lhs global position info
489  double lhs_phi = geometry_->chamber(lhs.first)->toGlobal(ilhs->localPosition()).phi();
490  double lhs_theta = geometry_->chamber(lhs.first)->toGlobal(ilhs->localPosition()).theta();
491  /*
492  std::cout << "RHS Segment Parameters:" << std::endl
493  << "\t RHS Phi : " << rhs_phi << std::endl
494  << "\t RHS Theta : " << rhs_theta << std::endl
495  << "\t RHS dx/dz : " << rhs_dxdz << " +- " << rhs_dxdz_err << std::endl
496  << "\t RHS dy/dz : " << rhs_dydz << " +- " << rhs_dydz_err << std::endl;
497 
498  std::cout << "LHS Segment Parameters:" << std::endl
499  << "\t LHS Phi : " << lhs_phi << std::endl
500  << "\t LHS Theta : " << lhs_theta << std::endl
501  << "\t LHS dx/dz : " << lhs_dxdz << " +- " << lhs_dxdz_err << std::endl
502  << "\t LHS dy/dz : " << lhs_dydz << " +- " << lhs_dydz_err << std::endl;
503  */
504 
505  double phidiff =
506  (fabs(rhs_phi - lhs_phi) > 2 * M_PI ? fabs(rhs_phi - lhs_phi) - 2 * M_PI : fabs(rhs_phi - lhs_phi));
507 
508  if (phidiff < ClusterDPhi && fabs(rhs_theta - lhs_theta) < ClusterDTheta) // phi overlap region is 37 degrees
509  result = true;
510  } // loop over duplicate segments
511  } // same endcap,station,ring
512 
513  return result;
514 }
mps_fire.i
i
Definition: mps_fire.py:428
CSCDetId
CSCDetId
Definition: RecoLocalMuon_RECO.doi:36
MuonMesh::ClusterDTheta
const double ClusterDTheta
Definition: MuonMesh.h:63
MuonMesh::geometry_
const CSCGeometry * geometry_
Definition: MuonMesh.h:58
GeomDet::toGlobal
GlobalPoint toGlobal(const Local2DPoint &lp) const
Conversion to the global R.F. from the R.F. of the GeomDet.
Definition: GeomDet.h:49
MuonMesh::doClustering
const bool doClustering
Definition: MuonMesh.h:61
theta
Geom::Theta< T > theta() const
Definition: Basic3DVectorLD.h:150
MuonMesh::doOverlaps
const bool doOverlaps
Definition: MuonMesh.h:61
hitfit::phidiff
double phidiff(double phi)
Normalized difference in azimuthal angles to a range between .
Definition: fourvec.cc:220
MuonMesh::ClusterDPhi
const double ClusterDPhi
Definition: MuonMesh.h:63
dt_dqm_sourceclient_common_cff.reco
tuple reco
Definition: dt_dqm_sourceclient_common_cff.py:111
spr::find
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:19
MuonMesh::pruneMesh
void pruneMesh()
Definition: MuonMesh.cc:151
MuonMesh::isClusteredWith
bool isClusteredWith(const std::pair< CSCDetId, CSCSegmentRef > &lhs, const std::pair< CSCDetId, CSCSegmentRef > &rhs) const
Definition: MuonMesh.cc:443
detailsBasic3DVector::z
float float float z
Definition: extBasic3DVector.h:14
MuonMesh::MuonMesh
MuonMesh(const edm::ParameterSet &)
Definition: MuonMesh.cc:9
mps_fire.result
tuple result
Definition: mps_fire.py:311
reco::Muon::ArbitrationType
ArbitrationType
define arbitration schemes
Definition: Muon.h:187
MuonMesh::fillMesh
void fillMesh(std::vector< reco::Muon > *)
Definition: MuonMesh.cc:18
reco::MuonSegmentMatch::BelongsToTrackByCleaning
static const unsigned int BelongsToTrackByCleaning
Definition: MuonSegmentMatch.h:32
M_PI
#define M_PI
Definition: BXVectorInputProducer.cc:50
CSCDetId::ring
int ring() const
Definition: CSCDetId.h:68
reco::MuonSegmentMatch::BelongsToTrackByOvlClean
static const unsigned int BelongsToTrackByOvlClean
Definition: MuonSegmentMatch.h:29
MuonMesh::mesh_
MeshType mesh_
Definition: MuonMesh.h:55
MuonMesh::isDuplicateOf
bool isDuplicateOf(const CSCSegmentRef &lhs, const CSCSegmentRef &rhs) const
Definition: MuonMesh.cc:297
CSCGeometry::chamber
const CSCChamber * chamber(CSCDetId id) const
Return the chamber corresponding to given DetId.
Definition: CSCGeometry.cc:100
MuonMesh::doME1a
const bool doME1a
Definition: MuonMesh.h:61
dataset.end
string end
Definition: dataset.py:937
MuonSubdetId::DT
static constexpr int DT
Definition: MuonSubdetId.h:11
reco::MuonSegmentMatch::BelongsToTrackByClusClean
static const unsigned int BelongsToTrackByClusClean
Definition: MuonSegmentMatch.h:30
MuonSubdetId::CSC
static constexpr int CSC
Definition: MuonSubdetId.h:12
Generated for CMSSW Reference Manual by doxygen 1.8.5