00001 // -*- C++ -*- // // Package: MuonSegmentMatcher // Class: MuonSegmentMatcher // /**\class MuonSegmentMatcher MuonSegmentMatcher.cc 00002 // Description: <one line class summary> 00003 // Implementation: 00004 // <Notes on implementation> 00005 //*/ 00006 // 00007 // Original Author: Alan Tua 00008 // Created: Wed Jul 9 21:40:17 CEST 2008 00009 // $Id: MuonSegmentMatcher.cc,v 1.4 2009/01/06 13:12:53 ptraczyk Exp $ 00010 // 00011 // 00012 00013 #include "RecoMuon/TrackingTools/interface/MuonSegmentMatcher.h" 00014 00015 00016 // user include files 00017 #include "FWCore/Framework/interface/EDAnalyzer.h" 00018 #include "FWCore/Framework/interface/Event.h" 00019 00020 #include "DataFormats/MuonReco/interface/Muon.h" 00021 #include "RecoMuon/TrackingTools/interface/MuonServiceProxy.h" 00022 00023 #include "DataFormats/TrackingRecHit/interface/TrackingRecHit.h" 00024 #include "DataFormats/TrackReco/interface/Track.h" 00025 00026 #include "DataFormats/MuonDetId/interface/MuonSubdetId.h" 00027 #include "DataFormats/MuonDetId/interface/CSCDetId.h" 00028 00029 #include "DataFormats/DTRecHit/interface/DTRecSegment4DCollection.h" 00030 #include "DataFormats/DTRecHit/interface/DTRecHitCollection.h" 00031 #include "DataFormats/CSCRecHit/interface/CSCSegmentCollection.h" 00032 #include "DataFormats/CSCRecHit/interface/CSCRecHit2DCollection.h" 00033 #include "DataFormats/GeometryVector/interface/GlobalPoint.h" 00034 00035 // system include files 00036 #include <vector> 00037 #include <iostream> 00038 00039 class MuonServiceProxy; 00040 class MuonSegmentMatcher; 00041 00042 using namespace std; 00043 00044 // constructors and destructor 00045 00046 MuonSegmentMatcher::MuonSegmentMatcher(const edm::ParameterSet& matchParameters, MuonServiceProxy* service) 00047 : 00048 theService(service), 00049 DTSegmentTags_(matchParameters.getUntrackedParameter<edm::InputTag>("DTsegments")), 00050 CSCSegmentTags_(matchParameters.getUntrackedParameter<edm::InputTag>("CSCsegments")), 00051 dtTightMatch(matchParameters.getParameter<bool>("TightMatchDT")), 00052 cscTightMatch(matchParameters.getParameter<bool>("TightMatchCSC")) 00053 { 00054 } 00055 00056 MuonSegmentMatcher::~MuonSegmentMatcher() 00057 { 00058 } 00059 00060 // member functions 00061 vector<const DTRecSegment4D*> MuonSegmentMatcher::matchDT(const reco::Track &muon, const edm::Event& event) 00062 { 00063 using namespace edm; 00064 00065 edm::Handle<DTRecSegment4DCollection> dtRecHits; 00066 event.getByLabel(DTSegmentTags_, dtRecHits); 00067 00068 vector<const DTRecSegment4D*> pointerTo4DSegments; 00069 00070 TrackingRecHitRefVector dtHits; 00071 00072 // Loop and select DT recHits 00073 for(trackingRecHit_iterator hit = muon.recHitsBegin(); hit != muon.recHitsEnd(); ++hit) { 00074 if ( (*hit)->geographicalId().det() != DetId::Muon ) continue; 00075 if ( (*hit)->geographicalId().subdetId() != MuonSubdetId::DT ) continue; 00076 if (!(*hit)->isValid()) continue; 00077 dtHits.push_back(*hit); 00078 } 00079 00080 double PhiCutParameter=0.001; 00081 double ZCutParameter=0.001; 00082 double matchRatioZ=0; 00083 double matchRatioPhi=0; 00084 00085 for (DTRecSegment4DCollection::const_iterator rechit = dtRecHits->begin(); rechit!=dtRecHits->end();++rechit) { 00086 00087 double nhitsPhi = 0; 00088 double nhitsZ = 0; 00089 00090 if(rechit->hasZed()) { 00091 double countMuonDTHits = 0; 00092 double countAgreeingHits=0; 00093 00094 const DTRecSegment2D* segmZ; 00095 segmZ = dynamic_cast<const DTRecSegment2D*>(rechit->zSegment()); 00096 nhitsZ = segmZ->recHits().size(); 00097 00098 const vector<DTRecHit1D> hits1d = segmZ->specificRecHits(); 00099 DTChamberId chamberSegIdT((segmZ->geographicalId()).rawId()); 00100 00101 // Loop over muon recHits 00102 for(trackingRecHit_iterator hit = dtHits.begin(); hit != dtHits.end(); ++hit) { 00103 00104 DetId idT = (*hit)->geographicalId(); 00105 DTChamberId dtDetIdHitT(idT.rawId()); 00106 DTSuperLayerId dtDetLayerIdHitT(idT.rawId()); 00107 00108 LocalPoint pointLocal = (*hit)->localPosition(); 00109 00110 if ((chamberSegIdT==dtDetIdHitT) && (dtDetLayerIdHitT.superlayer()==2)) countMuonDTHits++; 00111 00112 for (vector<DTRecHit1D>::const_iterator hiti=hits1d.begin(); hiti!=hits1d.end(); hiti++) { 00113 00114 // Pick the one in the same DT Layer as the 1D hit 00115 if(!(hiti->geographicalId().rawId()==idT.rawId())) continue; 00116 00117 // and compare the local positions 00118 LocalPoint segLocal = hiti->localPosition(); 00119 // cout << "Zed Segment Point = "<<pointLocal<<" Muon Point = "<<segLocal<<" " << endl; 00120 if ((fabs(pointLocal.x()-segLocal.x())<ZCutParameter) && 00121 (fabs(pointLocal.y()-segLocal.y())<ZCutParameter)) 00122 countAgreeingHits++; 00123 } //End Segment Hit Iteration 00124 } //End Muon Hit Iteration 00125 00126 matchRatioZ = countAgreeingHits/countMuonDTHits; 00127 } //End HasZed Check 00128 00129 if(rechit->hasPhi()) { 00130 double countMuonDTHits = 0; 00131 double countAgreeingHits=0; 00132 00133 //PREPARE PARAMETERS FOR SEGMENT DETECTOR GEOMETRY 00134 const DTRecSegment2D* segmPhi; 00135 segmPhi = dynamic_cast<const DTRecSegment2D*>(rechit->phiSegment()); 00136 nhitsPhi = segmPhi->recHits().size(); 00137 00138 const vector<DTRecHit1D> hits1d = segmPhi->specificRecHits(); 00139 DTChamberId chamberSegIdT((segmPhi->geographicalId()).rawId()); 00140 00141 // Loop over muon recHits 00142 for(trackingRecHit_iterator hit = dtHits.begin(); hit != dtHits.end(); ++hit) { 00143 00144 DetId idT = (*hit)->geographicalId(); 00145 DTChamberId dtDetIdHitT(idT.rawId()); 00146 DTSuperLayerId dtDetLayerIdHitT(idT.rawId()); 00147 00148 LocalPoint pointLocal = (*hit)->localPosition(); //Localposition is in DTLayer https://cmslxr.fnal.gov/lxr/source/DataFormats/DTRecHit/interface/DTRecHit1D.h 00149 00150 if ((chamberSegIdT==dtDetIdHitT)&&((dtDetLayerIdHitT.superlayer()==1)||(dtDetLayerIdHitT.superlayer()==3))) 00151 countMuonDTHits++; 00152 00153 for (vector<DTRecHit1D>::const_iterator hiti=hits1d.begin(); hiti!=hits1d.end(); hiti++) { 00154 00155 // Pick the one in the same DT Layer as the 1D hit 00156 if(!(hiti->geographicalId().rawId()==idT.rawId())) continue; 00157 00158 // and compare the local positions 00159 LocalPoint segLocal = hiti->localPosition(); 00160 // cout << " Phi Segment Point = "<<pointLocal<<" Muon Point = "<<segLocal<<" " << endl; 00161 00162 if ((fabs(pointLocal.x()-segLocal.x())<PhiCutParameter) && 00163 (fabs(pointLocal.y()-segLocal.y())<PhiCutParameter)) 00164 countAgreeingHits++; 00165 } // End Segment Hit Iteration 00166 } // End Muon Hit Iteration 00167 00168 matchRatioPhi = countAgreeingHits/countMuonDTHits; 00169 } // End HasPhi Check 00170 00171 // DTChamberId chamberSegId2((rechit->geographicalId()).rawId()); 00172 00173 if (dtTightMatch && nhitsPhi && nhitsZ) { 00174 if((matchRatioPhi>0.9)&&(matchRatioZ>0.9)) { 00175 // cout<<"Making a tight match in Chamber "<<chamberSegId2<<endl; 00176 pointerTo4DSegments.push_back(&(*rechit)); 00177 } 00178 } else { 00179 if((matchRatioPhi>0.9 && nhitsPhi)||(matchRatioZ>0.9 && nhitsZ)) { 00180 // cout<<"Making a loose match in Chamber "<<chamberSegId2<<endl; 00181 pointerTo4DSegments.push_back(&(*rechit)); 00182 } 00183 } 00184 00185 } //End DT Segment Iteration 00186 00187 return pointerTo4DSegments; 00188 } 00189 00190 00191 00192 vector<const CSCSegment*> MuonSegmentMatcher::matchCSC(const reco::Track& muon, const edm::Event& event) 00193 { 00194 00195 using namespace edm; 00196 00197 edm::Handle<CSCSegmentCollection> allSegmentsCSC; 00198 event.getByLabel(CSCSegmentTags_, allSegmentsCSC); 00199 00200 vector<const CSCSegment*> pointerToCSCSegments; 00201 00202 double matchRatioCSC=0; 00203 int numCSC = 0; 00204 double CSCXCut = 0.001; 00205 double CSCYCut = 0.001; 00206 double countMuonCSCHits = 0; 00207 00208 for(CSCSegmentCollection::const_iterator segmentCSC = allSegmentsCSC->begin(); segmentCSC != allSegmentsCSC->end(); segmentCSC++) { 00209 double CSCcountAgreeingHits=0; 00210 00211 numCSC++; 00212 double CSCnhits = segmentCSC->recHits().size(); 00213 const vector<CSCRecHit2D>& CSCRechits2D = segmentCSC->specificRecHits(); 00214 countMuonCSCHits = 0; 00215 CSCDetId myChamber((*segmentCSC).geographicalId().rawId()); 00216 00217 for(trackingRecHit_iterator hitC = muon.recHitsBegin(); hitC != muon.recHitsEnd(); ++hitC) { 00218 if ( (*hitC)->geographicalId().det() != DetId::Muon ) continue; 00219 if ( (*hitC)->geographicalId().subdetId() != MuonSubdetId::CSC ) continue; 00220 if(!(*hitC)->isValid()) continue; 00221 00222 //DETECTOR CONSTRUCTION 00223 DetId id = (*hitC)->geographicalId(); 00224 CSCDetId cscDetIdHit(id.rawId()); 00225 00226 if(!(cscDetIdHit.ring()==myChamber.ring())) continue; 00227 if(!(cscDetIdHit.station()==myChamber.station())) continue; 00228 if(!(cscDetIdHit.endcap()==myChamber.endcap())) continue; 00229 if(!(cscDetIdHit.chamber()==myChamber.chamber())) continue; 00230 00231 countMuonCSCHits++; 00232 00233 LocalPoint positionLocalCSC = (*hitC)->localPosition(); 00234 00235 for (vector<CSCRecHit2D>::const_iterator hiti=CSCRechits2D.begin(); hiti!=CSCRechits2D.end(); hiti++) { 00236 00237 CSCDetId cscDetId((hiti->geographicalId()).rawId()); 00238 00239 if ((*hitC)->geographicalId().rawId()!=(hiti->geographicalId()).rawId()) continue; 00240 00241 LocalPoint segLocalCSC = hiti->localPosition(); 00242 // cout<<"Layer Id (MuonHit) = "<<cscDetIdHit<<" Muon Local Position (det frame) "<<positionLocalCSC <<endl; 00243 // cout<<"Layer Id (CSCHit) = "<<cscDetId<<" Hit Local Position (det frame) "<<segLocalCSC <<endl; 00244 if((fabs(positionLocalCSC.x()-segLocalCSC.x())<CSCXCut) && 00245 (fabs(positionLocalCSC.y()-segLocalCSC.y())<CSCYCut)) { 00246 CSCcountAgreeingHits++; 00247 // cout << " Matched." << endl; 00248 } 00249 }//End 2D rechit iteration 00250 }//End muon hit iteration 00251 00252 matchRatioCSC = CSCcountAgreeingHits/countMuonCSCHits; 00253 00254 // cout<<"Matching Ratio(CSC): "<<matchRatioCSC<<" and Num of Hits in Muon "<<countMuonCSCHits<<endl; 00255 // cout<<"Num of Hits in Segment = "<<CSCnhits<<endl; 00256 00257 if ((matchRatioCSC>0.9) && ((countMuonCSCHits>1) || !cscTightMatch)) pointerToCSCSegments.push_back(&(*segmentCSC)); 00258 00259 } //End CSC Segment Iteration 00260 00261 return pointerToCSCSegments; 00262 00263 } 00264 00265 //define this as a plug-in 00266 //DEFINE_FWK_MODULE(MuonSegmentMatcher);