File indexing completed on 2023-03-17 10:58:15
0001
0002
0003
0004
0005
0006
0007
0008 #include "DQMOffline/Muon/interface/SegmentTrackAnalyzer.h"
0009
0010 #include "DataFormats/Common/interface/Handle.h"
0011 #include "DataFormats/MuonDetId/interface/MuonSubdetId.h"
0012 #include "RecoMuon/TransientTrackingRecHit/interface/MuonTransientTrackingRecHit.h"
0013 #include "DataFormats/DTRecHit/interface/DTRecSegment4DCollection.h"
0014 #include "DataFormats/DTRecHit/interface/DTRecHitCollection.h"
0015 #include "DataFormats/SiStripDetId/interface/SiStripDetId.h"
0016 #include "DataFormats/SiPixelDetId/interface/PixelSubdetector.h"
0017 #include "RecoMuon/TrackingTools/interface/MuonServiceProxy.h"
0018
0019 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0020 #include <string>
0021
0022 using namespace std;
0023 using namespace edm;
0024
0025
0026
0027 SegmentTrackAnalyzer::SegmentTrackAnalyzer(const edm::ParameterSet& pSet) {
0028 parameters = pSet;
0029
0030
0031 theMuTrackCollectionLabel_ =
0032 consumes<reco::TrackCollection>(parameters.getParameter<edm::InputTag>("MuTrackCollection"));
0033
0034
0035 const ParameterSet SegmentsTrackAssociatorParameters =
0036 parameters.getParameter<ParameterSet>("SegmentsTrackAssociatorParameters");
0037
0038 edm::ConsumesCollector iC = consumesCollector();
0039 theSegmentsAssociator = new SegmentsTrackAssociator(SegmentsTrackAssociatorParameters, iC);
0040
0041 trackCollection = parameters.getParameter<edm::InputTag>("MuTrackCollection").label() +
0042 parameters.getParameter<edm::InputTag>("MuTrackCollection").instance();
0043
0044 etaBin = parameters.getParameter<int>("etaBin");
0045 etaMin = parameters.getParameter<double>("etaMin");
0046 etaMax = parameters.getParameter<double>("etaMax");
0047 phiBin = parameters.getParameter<int>("phiBin");
0048 phiMin = parameters.getParameter<double>("phiMin");
0049 phiMax = parameters.getParameter<double>("phiMax");
0050 ptBin = parameters.getParameter<int>("ptBin");
0051 ptMin = parameters.getParameter<double>("ptMin");
0052 ptMax = parameters.getParameter<double>("ptMax");
0053 }
0054 void SegmentTrackAnalyzer::bookHistograms(DQMStore::IBooker& ibooker,
0055 const edm::Run& ,
0056 const edm::EventSetup& ) {
0057 ibooker.cd();
0058 ibooker.setCurrentFolder("Muons/SegmentTrackAnalyzer");
0059
0060
0061 hitsNotUsed = ibooker.book1D("HitsNotUsedForGlobalTracking_" + trackCollection,
0062 "recHits not used for GLB [" + trackCollection + "]",
0063 50,
0064 -0.5,
0065 49.5);
0066 hitsNotUsedPercentual = ibooker.book1D("HitsNotUsedForGlobalTrackingDvHitUsed_" + trackCollection,
0067 "(recHits_{notUsedForGLB}) / (recHits_{GLB}) [" + trackCollection + "]",
0068 100,
0069 0,
0070 1.);
0071
0072 TrackSegm = ibooker.book2D("trackSegments_" + trackCollection,
0073 "Number of segments associated to the track [" + trackCollection + "]",
0074 3,
0075 0.5,
0076 3.5,
0077 8,
0078 0,
0079 8);
0080 TrackSegm->setBinLabel(1, "DT+CSC", 1);
0081 TrackSegm->setBinLabel(2, "DT", 1);
0082 TrackSegm->setBinLabel(3, "CSC", 1);
0083
0084 hitStaProvenance = ibooker.book1D("trackHitStaProvenance_" + trackCollection,
0085 "Number of recHits_{STAinTrack} [" + trackCollection + "]",
0086 7,
0087 0.5,
0088 7.5);
0089 hitStaProvenance->setBinLabel(1, "DT");
0090 hitStaProvenance->setBinLabel(2, "CSC");
0091 hitStaProvenance->setBinLabel(3, "RPC");
0092 hitStaProvenance->setBinLabel(4, "DT+CSC");
0093 hitStaProvenance->setBinLabel(5, "DT+RPC");
0094 hitStaProvenance->setBinLabel(6, "CSC+RPC");
0095 hitStaProvenance->setBinLabel(7, "DT+CSC+RPC");
0096
0097 if (trackCollection != "standAloneMuons") {
0098 hitTkrProvenance = ibooker.book1D("trackHitTkrProvenance_" + trackCollection,
0099 "Number of recHits_{TKinTrack} [" + trackCollection + "]",
0100 6,
0101 0.5,
0102 6.5);
0103 hitTkrProvenance->setBinLabel(1, "PixBarrel");
0104 hitTkrProvenance->setBinLabel(2, "PixEndCap");
0105 hitTkrProvenance->setBinLabel(3, "TIB");
0106 hitTkrProvenance->setBinLabel(4, "TID");
0107 hitTkrProvenance->setBinLabel(5, "TOB");
0108 hitTkrProvenance->setBinLabel(6, "TEC");
0109 }
0110
0111 trackHitPercentualVsEta =
0112 ibooker.book2D("trackHitDivtrackSegmHitVsEta_" + trackCollection,
0113 "(recHits_{Track} / recHits_{associatedSegm}) vs #eta [" + trackCollection + "]",
0114 etaBin,
0115 etaMin,
0116 etaMax,
0117 20,
0118 0,
0119 1);
0120 dtTrackHitPercentualVsEta =
0121 ibooker.book2D("dtTrackHitDivtrackSegmHitVsEta_" + trackCollection,
0122 "(recHits_{DTinTrack} / recHits_{associatedSegm}) vs #eta [" + trackCollection + "]",
0123 etaBin,
0124 etaMin,
0125 etaMax,
0126 20,
0127 0,
0128 1);
0129 cscTrackHitPercentualVsEta =
0130 ibooker.book2D("cscTrackHitDivtrackSegmHitVsEta_" + trackCollection,
0131 "(recHits_{CSCinTrack} / recHits_{associatedSegm}) vs #eta [" + trackCollection + "]",
0132 etaBin,
0133 etaMin,
0134 etaMax,
0135 20,
0136 0,
0137 1);
0138
0139 trackHitPercentualVsPhi =
0140 ibooker.book2D("trackHitDivtrackSegmHitVsPhi_" + trackCollection,
0141 "(recHits_{Track} / recHits_{associatedSegm}) vs #phi [" + trackCollection + "]",
0142 phiBin,
0143 phiMin,
0144 phiMax,
0145 20,
0146 0,
0147 1);
0148 trackHitPercentualVsPhi->setAxisTitle("rad", 2);
0149 dtTrackHitPercentualVsPhi =
0150 ibooker.book2D("dtTrackHitDivtrackSegmHitVsPhi_" + trackCollection,
0151 "(recHits_{DTinTrack} / recHits_{associatedSegm}) vs #phi [" + trackCollection + "]",
0152 phiBin,
0153 phiMin,
0154 phiMax,
0155 20,
0156 0,
0157 1);
0158 dtTrackHitPercentualVsPhi->setAxisTitle("rad", 2);
0159 cscTrackHitPercentualVsPhi =
0160 ibooker.book2D("cscTrackHitDivtrackSegmHitVsPhi_" + trackCollection,
0161 "(recHits_{CSCinTrack} / recHits_{associatedSegm}) vs #phi [" + trackCollection + "]",
0162 phiBin,
0163 phiMin,
0164 phiMax,
0165 20,
0166 0,
0167 1);
0168 cscTrackHitPercentualVsPhi->setAxisTitle("rad", 2);
0169
0170 trackHitPercentualVsPt =
0171 ibooker.book2D("trackHitDivtrackSegmHitVsPt_" + trackCollection,
0172 "(recHits_{Track} / recHits_{associatedSegm}) vs 1/p_{t} [" + trackCollection + "]",
0173 ptBin,
0174 ptMin,
0175 ptMax,
0176 20,
0177 0,
0178 1);
0179 trackHitPercentualVsPt->setAxisTitle("GeV", 2);
0180 dtTrackHitPercentualVsPt =
0181 ibooker.book2D("dtTrackHitDivtrackSegmHitVsPt_" + trackCollection,
0182 "(recHits_{DTinTrack} / recHits_{associatedSegm}) vs 1/p_{t} [" + trackCollection + "]",
0183 ptBin,
0184 ptMin,
0185 ptMax,
0186 20,
0187 0,
0188 1);
0189 dtTrackHitPercentualVsPt->setAxisTitle("GeV", 2);
0190 cscTrackHitPercentualVsPt =
0191 ibooker.book2D("cscTrackHitDivtrackSegmHitVsPt_" + trackCollection,
0192 "(recHits_{CSCinTrack} / recHits_{associatedSegm}) vs 1/p_{t} [" + trackCollection + "]",
0193 ptBin,
0194 ptMin,
0195 ptMax,
0196 20,
0197 0,
0198 1);
0199 cscTrackHitPercentualVsPt->setAxisTitle("GeV", 2);
0200 }
0201
0202 void SegmentTrackAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0203 Handle<reco::TrackCollection> glbTracks;
0204 iEvent.getByToken(theMuTrackCollectionLabel_, glbTracks);
0205
0206 for (reco::TrackCollection::const_iterator recoTrack = glbTracks->begin(); recoTrack != glbTracks->end();
0207 ++recoTrack) {
0208 MuonTransientTrackingRecHit::MuonRecHitContainer segments =
0209 theSegmentsAssociator->associate(iEvent, iSetup, *recoTrack);
0210
0211 #ifdef DEBUG
0212 cout << "[SegmentTrackAnalyzer] # of segments associated to the track: " << (segments).size() << endl;
0213 #endif
0214
0215
0216 int hitsFromDt = 0;
0217 int hitsFromCsc = 0;
0218 int hitsFromRpc = 0;
0219 int hitsFromTk = 0;
0220 int hitsFromTrack = 0;
0221 int hitsFromSegmDt = 0;
0222 int hitsFromSegmCsc = 0;
0223
0224 int segmFromDt = 0;
0225 int segmFromCsc = 0;
0226
0227 for (MuonTransientTrackingRecHit::MuonRecHitContainer::const_iterator segment = segments.begin();
0228 segment != segments.end();
0229 segment++) {
0230 DetId id = (*segment)->geographicalId();
0231
0232
0233 if (id.det() == DetId::Muon && id.subdetId() == MuonSubdetId::DT) {
0234 ++segmFromDt;
0235
0236 const DTRecSegment4D* seg4D = dynamic_cast<const DTRecSegment4D*>((*segment)->hit());
0237 if ((*seg4D).hasPhi())
0238 hitsFromSegmDt += (*seg4D).phiSegment()->specificRecHits().size();
0239 if ((*seg4D).hasZed())
0240 hitsFromSegmDt += (*seg4D).zSegment()->specificRecHits().size();
0241 }
0242
0243
0244 if (id.det() == DetId::Muon && id.subdetId() == MuonSubdetId::CSC) {
0245 hitsFromSegmCsc += (*segment)->recHits().size();
0246 segmFromCsc++;
0247 }
0248
0249 #ifdef DEBUG
0250 cout << "[SegmentTrackAnalyzer] # of segments from DT and CSC: " << segmFromDt << ", " << segmFromCsc << endl;
0251 cout << "[SegmentTrackAnalyzer] # of HITS from segments from DT: " << hitsFromSegmDt << endl;
0252 cout << "[SegmentTrackAnalyzer] # of HITS from segments from CSC " << hitsFromSegmCsc << endl;
0253 #endif
0254 }
0255
0256
0257 for (trackingRecHit_iterator recHit = recoTrack->recHitsBegin(); recHit != recoTrack->recHitsEnd(); ++recHit) {
0258 hitsFromTrack++;
0259 DetId id = (*recHit)->geographicalId();
0260
0261 if (id.det() == DetId::Muon && id.subdetId() == MuonSubdetId::DT)
0262 hitsFromDt++;
0263
0264 if (id.det() == DetId::Muon && id.subdetId() == MuonSubdetId::CSC)
0265 hitsFromCsc++;
0266
0267 if (id.det() == DetId::Muon && id.subdetId() == MuonSubdetId::RPC)
0268 hitsFromRpc++;
0269
0270 if (id.det() == DetId::Tracker) {
0271 hitsFromTk++;
0272 if (id.subdetId() == PixelSubdetector::PixelBarrel)
0273 hitTkrProvenance->Fill(1);
0274 if (id.subdetId() == PixelSubdetector::PixelEndcap)
0275 hitTkrProvenance->Fill(2);
0276 if (id.subdetId() == SiStripDetId::TIB)
0277 hitTkrProvenance->Fill(3);
0278 if (id.subdetId() == SiStripDetId::TID)
0279 hitTkrProvenance->Fill(4);
0280 if (id.subdetId() == SiStripDetId::TOB)
0281 hitTkrProvenance->Fill(5);
0282 if (id.subdetId() == SiStripDetId::TEC)
0283 hitTkrProvenance->Fill(6);
0284 }
0285 }
0286
0287
0288 hitsNotUsed->Fill(hitsFromSegmDt + hitsFromSegmCsc + hitsFromRpc + hitsFromTk - hitsFromTrack);
0289 hitsNotUsedPercentual->Fill(double(hitsFromSegmDt + hitsFromSegmCsc + hitsFromRpc + hitsFromTk - hitsFromTrack) /
0290 hitsFromTrack);
0291
0292 if (hitsFromDt != 0 && hitsFromCsc != 0)
0293 TrackSegm->Fill(1, segmFromDt + segmFromCsc);
0294 if (hitsFromDt != 0 && hitsFromCsc == 0)
0295 TrackSegm->Fill(2, segmFromDt);
0296 if (hitsFromDt == 0 && hitsFromCsc != 0)
0297 TrackSegm->Fill(3, segmFromCsc);
0298
0299 if (hitsFromDt != 0 && hitsFromCsc == 0 && hitsFromRpc == 0)
0300 hitStaProvenance->Fill(1);
0301 if (hitsFromCsc != 0 && hitsFromDt == 0 && hitsFromRpc == 0)
0302 hitStaProvenance->Fill(2);
0303 if (hitsFromRpc != 0 && hitsFromDt == 0 && hitsFromCsc == 0)
0304 hitStaProvenance->Fill(3);
0305 if (hitsFromDt != 0 && hitsFromCsc != 0 && hitsFromRpc == 0)
0306 hitStaProvenance->Fill(4);
0307 if (hitsFromDt != 0 && hitsFromRpc != 0 && hitsFromCsc == 0)
0308 hitStaProvenance->Fill(5);
0309 if (hitsFromCsc != 0 && hitsFromRpc != 0 && hitsFromDt == 0)
0310 hitStaProvenance->Fill(6);
0311 if (hitsFromDt != 0 && hitsFromCsc != 0 && hitsFromRpc != 0)
0312 hitStaProvenance->Fill(7);
0313
0314 if (hitsFromSegmDt + hitsFromSegmCsc != 0) {
0315 trackHitPercentualVsEta->Fill(recoTrack->eta(),
0316 double(hitsFromDt + hitsFromCsc) / (hitsFromSegmDt + hitsFromSegmCsc));
0317 trackHitPercentualVsPhi->Fill(recoTrack->phi(),
0318 double(hitsFromDt + hitsFromCsc) / (hitsFromSegmDt + hitsFromSegmCsc));
0319 trackHitPercentualVsPt->Fill(recoTrack->pt(),
0320 double(hitsFromDt + hitsFromCsc) / (hitsFromSegmDt + hitsFromSegmCsc));
0321 }
0322
0323 if (hitsFromSegmDt != 0) {
0324 dtTrackHitPercentualVsEta->Fill(recoTrack->eta(), double(hitsFromDt) / hitsFromSegmDt);
0325 dtTrackHitPercentualVsPhi->Fill(recoTrack->phi(), double(hitsFromDt) / hitsFromSegmDt);
0326 dtTrackHitPercentualVsPt->Fill(recoTrack->pt(), double(hitsFromDt) / hitsFromSegmDt);
0327 }
0328
0329 if (hitsFromSegmCsc != 0) {
0330 cscTrackHitPercentualVsEta->Fill(recoTrack->eta(), double(hitsFromCsc) / hitsFromSegmCsc);
0331 cscTrackHitPercentualVsPhi->Fill(recoTrack->phi(), double(hitsFromCsc) / hitsFromSegmCsc);
0332 cscTrackHitPercentualVsPt->Fill(recoTrack->pt(), double(hitsFromCsc) / hitsFromSegmCsc);
0333 }
0334 }
0335 }