CMS 3D CMS Logo

CMSSW_4_4_3_patch1/src/RecoMuon/TrackingTools/src/MuonSegmentMatcher.cc

Go to the documentation of this file.
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.10 2011/01/25 22:48:19 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.getParameter<edm::InputTag>("DTsegments")),
00050   CSCSegmentTags_(matchParameters.getParameter<edm::InputTag>("CSCsegments")),
00051   dtRadius_(matchParameters.getParameter<double>("DTradius")),
00052   dtTightMatch(matchParameters.getParameter<bool>("TightMatchDT")),
00053   cscTightMatch(matchParameters.getParameter<bool>("TightMatchCSC"))
00054 {
00055 }
00056 
00057 MuonSegmentMatcher::~MuonSegmentMatcher()
00058 {
00059 }
00060 
00061 // member functions
00062 vector<const DTRecSegment4D*> MuonSegmentMatcher::matchDT(const reco::Track &muon, const edm::Event& event)
00063 {
00064   using namespace edm;
00065 
00066   edm::Handle<DTRecSegment4DCollection> dtRecHits;
00067   event.getByLabel(DTSegmentTags_, dtRecHits);  
00068   
00069   vector<const DTRecSegment4D*> pointerTo4DSegments;
00070 
00071   TrackingRecHitRefVector dtHits;
00072 
00073   // Loop and select DT recHits
00074   for(trackingRecHit_iterator hit = muon.recHitsBegin(); hit != muon.recHitsEnd(); ++hit) {
00075     if ( !(*hit)->isValid()) continue; 
00076     if ( (*hit)->geographicalId().det() != DetId::Muon ) continue; 
00077     if ( (*hit)->geographicalId().subdetId() != MuonSubdetId::DT ) continue; 
00078     if (!(*hit)->isValid()) continue; 
00079     dtHits.push_back(*hit);
00080   }
00081   
00082   double PhiCutParameter=dtRadius_;
00083   double ZCutParameter=dtRadius_;
00084   double matchRatioZ=0;
00085   double matchRatioPhi=0;
00086 
00087   for (DTRecSegment4DCollection::const_iterator rechit = dtRecHits->begin(); rechit!=dtRecHits->end();++rechit) {
00088   
00089     if ( !rechit->isValid()) continue; 
00090 
00091     double nhitsPhi = 0;
00092     double nhitsZ = 0;
00093           
00094     if(rechit->hasZed()) {
00095       double countMuonDTHits = 0;
00096       double countAgreeingHits=0;
00097 
00098       const DTRecSegment2D* segmZ;
00099       segmZ = dynamic_cast<const DTRecSegment2D*>(rechit->zSegment());
00100       nhitsZ = segmZ->recHits().size(); 
00101                 
00102       const vector<DTRecHit1D> hits1d = segmZ->specificRecHits();
00103       DTChamberId chamberSegIdT((segmZ->geographicalId()).rawId());
00104                 
00105       // Loop over muon recHits
00106       for(trackingRecHit_iterator hit = dtHits.begin(); hit != dtHits.end(); ++hit) {
00107 
00108         if ( !(*hit)->isValid()) continue; 
00109         
00110         DetId idT = (*hit)->geographicalId();
00111         DTChamberId dtDetIdHitT(idT.rawId());
00112         DTSuperLayerId dtDetLayerIdHitT(idT.rawId());
00113 
00114         LocalPoint  pointLocal = (*hit)->localPosition();
00115                                         
00116         if ((chamberSegIdT==dtDetIdHitT) && (dtDetLayerIdHitT.superlayer()==2)) countMuonDTHits++;
00117 
00118         for (vector<DTRecHit1D>::const_iterator hiti=hits1d.begin(); hiti!=hits1d.end(); hiti++) {
00119 
00120           if ( !hiti->isValid()) continue; 
00121 
00122           // Pick the one in the same DT Layer as the 1D hit
00123           if(!(hiti->geographicalId().rawId()==idT.rawId())) continue; 
00124 
00125           // and compare the local positions
00126           LocalPoint segLocal = hiti->localPosition();
00127 //        cout << "Zed Segment Point = "<<pointLocal<<"    Muon Point = "<<segLocal<<"  Dist:  "
00128 //             << (fabs(pointLocal.x()-segLocal.x()))+(fabs(pointLocal.y()-segLocal.y()))<< endl;
00129           if ((fabs(pointLocal.x()-segLocal.x())<ZCutParameter) && 
00130               (fabs(pointLocal.y()-segLocal.y())<ZCutParameter)) 
00131             countAgreeingHits++;
00132         } //End Segment Hit Iteration
00133       } //End Muon Hit Iteration
00134                 
00135       matchRatioZ = countMuonDTHits == 0 ? 0 : countAgreeingHits/countMuonDTHits;
00136     } //End HasZed Check
00137                         
00138     if(rechit->hasPhi()) {
00139       double countMuonDTHits = 0;
00140       double countAgreeingHits=0;
00141 
00142       //PREPARE PARAMETERS FOR SEGMENT DETECTOR GEOMETRY
00143       const DTRecSegment2D* segmPhi;
00144       segmPhi = dynamic_cast<const DTRecSegment2D*>(rechit->phiSegment());
00145       nhitsPhi = segmPhi->recHits().size();
00146                 
00147       const vector<DTRecHit1D> hits1d = segmPhi->specificRecHits();
00148       DTChamberId chamberSegIdT((segmPhi->geographicalId()).rawId());
00149                 
00150       // Loop over muon recHits
00151       for(trackingRecHit_iterator hit = dtHits.begin(); hit != dtHits.end(); ++hit) {
00152 
00153         if ( !(*hit)->isValid()) continue; 
00154 
00155         DetId idT = (*hit)->geographicalId();
00156         DTChamberId dtDetIdHitT(idT.rawId());
00157         DTSuperLayerId dtDetLayerIdHitT(idT.rawId());
00158 
00159         LocalPoint pointLocal = (*hit)->localPosition(); //Localposition is in DTLayer https://cmslxr.fnal.gov/lxr/source/DataFormats/DTRecHit/interface/DTRecHit1D.h
00160 
00161         if ((chamberSegIdT==dtDetIdHitT)&&((dtDetLayerIdHitT.superlayer()==1)||(dtDetLayerIdHitT.superlayer()==3))) 
00162           countMuonDTHits++;
00163 
00164         for (vector<DTRecHit1D>::const_iterator hiti=hits1d.begin(); hiti!=hits1d.end(); hiti++) {
00165 
00166           if ( !hiti->isValid()) continue; 
00167 
00168           // Pick the one in the same DT Layer as the 1D hit
00169           if(!(hiti->geographicalId().rawId()==idT.rawId())) continue; 
00170                                          
00171           // and compare the local positions
00172           LocalPoint segLocal = hiti->localPosition();
00173 //        cout << "     Phi Segment Point = "<<pointLocal<<"    Muon Point = "<<segLocal<<"  Dist:   " 
00174 //             << (fabs(pointLocal.x()-segLocal.x()))+(fabs(pointLocal.y()-segLocal.y()))<< endl;
00175 
00176           if ((fabs(pointLocal.x()-segLocal.x())<PhiCutParameter) && 
00177               (fabs(pointLocal.y()-segLocal.y())<PhiCutParameter))
00178             countAgreeingHits++; 
00179         } // End Segment Hit Iteration
00180       } // End Muon Hit Iteration
00181 
00182       matchRatioPhi = countMuonDTHits != 0 ? countAgreeingHits/countMuonDTHits : 0;
00183     } // End HasPhi Check
00184 
00185 //    DTChamberId chamberSegId2((rechit->geographicalId()).rawId());
00186     
00187     if (dtTightMatch && nhitsPhi && nhitsZ) {
00188       if((matchRatioPhi>0.9)&&(matchRatioZ>0.9)) {
00189 //      cout<<"Making a tight match in Chamber "<<chamberSegId2<<endl;
00190         pointerTo4DSegments.push_back(&(*rechit));
00191       }
00192     } else {
00193       if((matchRatioPhi>0.9 && nhitsPhi)||(matchRatioZ>0.9 && nhitsZ)) {
00194 //      cout<<"Making a loose match in Chamber "<<chamberSegId2<<endl;
00195         pointerTo4DSegments.push_back(&(*rechit));
00196       }
00197     }
00198     
00199   } //End DT Segment Iteration
00200 
00201   return pointerTo4DSegments;
00202 }
00203 
00204 
00205 
00206 vector<const CSCSegment*> MuonSegmentMatcher::matchCSC(const reco::Track& muon, const edm::Event& event)
00207 {
00208 
00209   using namespace edm;
00210 
00211   edm::Handle<CSCSegmentCollection> allSegmentsCSC;
00212   event.getByLabel(CSCSegmentTags_, allSegmentsCSC);
00213 
00214   vector<const CSCSegment*> pointerToCSCSegments;
00215 
00216   double matchRatioCSC=0;
00217   int numCSC = 0;
00218   double CSCXCut = 0.001;
00219   double CSCYCut = 0.001;
00220   double countMuonCSCHits = 0;
00221 
00222   for(CSCSegmentCollection::const_iterator segmentCSC = allSegmentsCSC->begin(); segmentCSC != allSegmentsCSC->end(); segmentCSC++) {
00223     double CSCcountAgreeingHits=0;
00224 
00225     if ( !segmentCSC->isValid()) continue; 
00226 
00227     numCSC++;
00228     const vector<CSCRecHit2D>& CSCRechits2D = segmentCSC->specificRecHits();
00229     countMuonCSCHits = 0;
00230     CSCDetId myChamber((*segmentCSC).geographicalId().rawId());
00231 
00232     for(trackingRecHit_iterator hitC = muon.recHitsBegin(); hitC != muon.recHitsEnd(); ++hitC) {
00233       if (!(*hitC)->isValid()) continue; 
00234       if ( (*hitC)->geographicalId().det() != DetId::Muon ) continue; 
00235       if ( (*hitC)->geographicalId().subdetId() != MuonSubdetId::CSC ) continue;
00236       if (!(*hitC)->isValid()) continue;
00237 
00238       //DETECTOR CONSTRUCTION
00239       DetId id = (*hitC)->geographicalId();
00240       CSCDetId cscDetIdHit(id.rawId());
00241 
00242       if(!(cscDetIdHit.ring()==myChamber.ring())) continue;
00243       if(!(cscDetIdHit.station()==myChamber.station())) continue;
00244       if(!(cscDetIdHit.endcap()==myChamber.endcap())) continue;
00245       if(!(cscDetIdHit.chamber()==myChamber.chamber())) continue;
00246 
00247       countMuonCSCHits++;
00248 
00249       LocalPoint positionLocalCSC = (*hitC)->localPosition();
00250         
00251       for (vector<CSCRecHit2D>::const_iterator hiti=CSCRechits2D.begin(); hiti!=CSCRechits2D.end(); hiti++) {
00252 
00253         if ( !hiti->isValid()) continue; 
00254         CSCDetId cscDetId((hiti->geographicalId()).rawId());
00255                 
00256         if ((*hitC)->geographicalId().rawId()!=(hiti->geographicalId()).rawId()) continue;
00257 
00258         LocalPoint segLocalCSC = hiti->localPosition();
00259         //              cout<<"Layer Id (MuonHit) =  "<<cscDetIdHit<<" Muon Local Position (det frame) "<<positionLocalCSC <<endl;
00260         //              cout<<"Layer Id  (CSCHit) =  "<<cscDetId<<"  Hit Local Position (det frame) "<<segLocalCSC <<endl;
00261         if((fabs(positionLocalCSC.x()-segLocalCSC.x())<CSCXCut) && 
00262            (fabs(positionLocalCSC.y()-segLocalCSC.y())<CSCYCut)) {
00263           CSCcountAgreeingHits++;
00264           //              cout << "   Matched." << endl;
00265         }  
00266       }//End 2D rechit iteration
00267     }//End muon hit iteration
00268     
00269     matchRatioCSC = countMuonCSCHits == 0 ? 0 : CSCcountAgreeingHits/countMuonCSCHits;
00270                 
00271     if ((matchRatioCSC>0.9) && ((countMuonCSCHits>1) || !cscTightMatch)) pointerToCSCSegments.push_back(&(*segmentCSC));
00272 
00273   } //End CSC Segment Iteration 
00274 
00275   return pointerToCSCSegments;
00276 
00277 }
00278 
00279 //define this as a plug-in
00280 //DEFINE_FWK_MODULE(MuonSegmentMatcher);