00001 #include "RecoMuon/TrackerSeedGenerator/interface/L1MuonSeedsMerger.h" 00002 00003 #include "FWCore/ParameterSet/interface/ParameterSet.h" 00004 #include "DataFormats/TrackReco/interface/Track.h" 00005 00006 L1MuonSeedsMerger::L1MuonSeedsMerger(const edm::ParameterSet& cfg) 00007 { 00008 theDeltaEtaCut = cfg.getParameter<double>("deltaEtaCut"); 00009 theDiffRelPtCut = cfg.getParameter<double>("diffRelPtCut"); 00010 } 00011 00012 void L1MuonSeedsMerger::resolve(std::vector<TrackAndHits>& tracks) const 00013 { 00014 sort(tracks.begin(),tracks.end(), Less()); 00015 typedef std::vector<TrackAndHits>::iterator Tracks_Itr; 00016 Tracks_Itr it1 = tracks.begin(); 00017 while (it1 != tracks.end() ) { 00018 for (Tracks_Itr it2 = it1+1; it1->first && it2<tracks.end(); it2++) { 00019 if (! it2->first) continue; 00020 if ( it2->first->eta() - it1->first->eta() > theDeltaEtaCut) break; 00021 switch ( compare( &(*it1), &(*it2) ) ) { 00022 case killFirst : 00023 delete it1->first; 00024 it1->first = 0; 00025 break; 00026 case killSecond : 00027 delete it2->first; 00028 it2->first = 0; 00029 break; 00030 case mergeTwo : 00031 *it2 = *(merge(&(*it1),&(*it2))); 00032 it1->first = 0; 00033 break; 00034 case goAhead : default: break; 00035 } 00036 } 00037 if (0 == it1->first) tracks.erase(it1); else it1++; 00038 } 00039 } 00040 00041 bool L1MuonSeedsMerger::Less::operator()(const TrackAndHits& a, const TrackAndHits& b) const 00042 { 00043 return (a.first->eta() < b.first->eta()); 00044 } 00045 00046 const L1MuonSeedsMerger::TrackAndHits* 00047 L1MuonSeedsMerger::merge(const TrackAndHits* a, const TrackAndHits* b) const 00048 { 00049 // temporary algorith, takes track with bigger pt, to be reimplemented 00050 if (a->first->pt() > b->first->pt()) { 00051 delete b->first; 00052 return a; 00053 } else { 00054 delete a->first; 00055 return b; 00056 } 00057 } 00058 00059 L1MuonSeedsMerger::Action 00060 L1MuonSeedsMerger::compare(const TrackAndHits* a, const TrackAndHits* b) const 00061 { 00062 int nshared = 0; 00063 const SeedingHitSet & hitsA = a->second; 00064 const SeedingHitSet & hitsB = b->second; 00065 for (unsigned int iHitA=0, nHitsA=hitsA.size(); iHitA < nHitsA; ++iHitA) { 00066 const TrackingRecHit* trhA = hitsA[iHitA]->hit(); 00067 for (unsigned int iHitB=0, nHitsB=hitsB.size(); iHitB < nHitsB; ++iHitB) { 00068 const TrackingRecHit* trhB = hitsB[iHitB]->hit(); 00069 if (trhA==trhB) nshared++; 00070 } 00071 } 00072 00073 if (nshared >= 2) { 00074 if (hitsA.size() >= 3 && hitsB.size() >= 3 ) 00075 return (a->first->chi2() > b->first->chi2()) ? killFirst : killSecond; 00076 else if (hitsB.size() >= 3) 00077 return killFirst; 00078 else 00079 return killSecond; 00080 } 00081 else if (nshared >= 1) { 00082 if (hitsA.size() != hitsB.size()) 00083 return (hitsA.size() < hitsB.size()) ? killFirst : killSecond; 00084 else if ( hitsA.size() >= 3 00085 && a->first->charge()==b->first->charge() 00086 && fabs(a->first->pt()-b->first->pt())/b->first->pt() < theDiffRelPtCut ) 00087 return (a->first->chi2() > b->first->chi2()) ? killFirst : killSecond; 00088 else if ( hitsA.size() == 2 ) 00089 return (a->first->pt() < b->first->pt()) ? killFirst : killSecond; 00090 else 00091 return goAhead; 00092 } 00093 else return goAhead; 00094 } 00095