CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
L1MuonSeedsMerger.cc
Go to the documentation of this file.
2 
5 
7 {
8  theDeltaEtaCut = cfg.getParameter<double>("deltaEtaCut");
9  theDiffRelPtCut = cfg.getParameter<double>("diffRelPtCut");
10 }
11 
12 void L1MuonSeedsMerger::resolve(std::vector<TrackAndHits>& tracks) const
13 {
14  sort(tracks.begin(),tracks.end(), Less());
15  typedef std::vector<TrackAndHits>::iterator Tracks_Itr;
16  Tracks_Itr it1 = tracks.begin();
17  while (it1 != tracks.end() ) {
18  for (Tracks_Itr it2 = it1+1; it1->first && it2<tracks.end(); it2++) {
19  if (! it2->first) continue;
20  if ( it2->first->eta() - it1->first->eta() > theDeltaEtaCut) break;
21  switch ( compare( &(*it1), &(*it2) ) ) {
22  case killFirst :
23  delete it1->first;
24  it1->first = 0;
25  break;
26  case killSecond :
27  delete it2->first;
28  it2->first = 0;
29  break;
30  case mergeTwo :
31  *it2 = *(merge(&(*it1),&(*it2)));
32  it1->first = 0;
33  break;
34  case goAhead : default: break;
35  }
36  }
37  if (0 == it1->first) tracks.erase(it1); else it1++;
38  }
39 }
40 
42 {
43  return (a.first->eta() < b.first->eta());
44 }
45 
48 {
49 // temporary algorith, takes track with bigger pt, to be reimplemented
50  if (a->first->pt() > b->first->pt()) {
51  delete b->first;
52  return a;
53  } else {
54  delete a->first;
55  return b;
56  }
57 }
58 
61 {
62  int nshared = 0;
63  const SeedingHitSet & hitsA = a->second;
64  const SeedingHitSet & hitsB = b->second;
65  for (unsigned int iHitA=0, nHitsA=hitsA.size(); iHitA < nHitsA; ++iHitA) {
66  const TrackingRecHit* trhA = hitsA[iHitA]->hit();
67  for (unsigned int iHitB=0, nHitsB=hitsB.size(); iHitB < nHitsB; ++iHitB) {
68  const TrackingRecHit* trhB = hitsB[iHitB]->hit();
69  if (trhA==trhB) nshared++;
70  }
71  }
72 
73  if (nshared >= 2) {
74  if (hitsA.size() >= 3 && hitsB.size() >= 3 )
75  return (a->first->chi2() > b->first->chi2()) ? killFirst : killSecond;
76  else if (hitsB.size() >= 3)
77  return killFirst;
78  else
79  return killSecond;
80  }
81  else if (nshared >= 1) {
82  if (hitsA.size() != hitsB.size())
83  return (hitsA.size() < hitsB.size()) ? killFirst : killSecond;
84  else if ( hitsA.size() >= 3
85  && a->first->charge()==b->first->charge()
86  && fabs(a->first->pt()-b->first->pt())/b->first->pt() < theDiffRelPtCut )
87  return (a->first->chi2() > b->first->chi2()) ? killFirst : killSecond;
88  else if ( hitsA.size() == 2 )
89  return (a->first->pt() < b->first->pt()) ? killFirst : killSecond;
90  else
91  return goAhead;
92  }
93  else return goAhead;
94 }
95 
T getParameter(std::string const &) const
tuple cfg
Definition: looper.py:259
std::pair< const reco::Track *, SeedingHitSet > TrackAndHits
bool operator()(const TrackAndHits &, const TrackAndHits &) const
Action compare(const TrackAndHits *, const TrackAndHits *) const
virtual void resolve(TracksAndHits &) const
L1MuonSeedsMerger(const edm::ParameterSet &cfg)
tuple tracks
Definition: testEve_cfg.py:39
double b
Definition: hdecay.h:120
const TrackAndHits * merge(const TrackAndHits *, const TrackAndHits *) const
double a
Definition: hdecay.h:121
unsigned int size() const
Definition: SeedingHitSet.h:44