File indexing completed on 2024-04-06 12:33:10
0001 #include "Validation/RecoMuon/src/RecoDisplacedMuonValidator.h"
0002
0003 #include "DataFormats/MuonReco/interface/MuonSelectors.h"
0004
0005 #include "FWCore/ServiceRegistry/interface/Service.h"
0006 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0007
0008
0009 #include "TrackingTools/TransientTrackingRecHit/interface/TransientTrackingRecHit.h"
0010 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
0011 #include "RecoMuon/TrackingTools/interface/MuonServiceProxy.h"
0012
0013
0014 #include "CommonTools/Utils/interface/StringCutObjectSelector.h"
0015
0016 #include "TMath.h"
0017
0018 using namespace std;
0019 using namespace edm;
0020 using namespace reco;
0021
0022 typedef TrajectoryStateOnSurface TSOS;
0023 typedef FreeTrajectoryState FTS;
0024
0025
0026
0027
0028 struct RecoDisplacedMuonValidator::MuonME {
0029 typedef MonitorElement* MEP;
0030
0031
0032 MEP hSimP_, hSimPt_, hSimEta_, hSimPhi_, hSimDxy_, hSimDz_;
0033
0034 MEP hP_, hPt_, hEta_, hPhi_;
0035 MEP hNSim_, hNMuon_;
0036
0037
0038 MEP hNTrks_, hNTrksEta_, hNTrksPt_;
0039 MEP hMisQPt_, hMisQEta_;
0040
0041
0042 MEP hErrP_, hErrPt_, hErrEta_, hErrPhi_;
0043 MEP hErrPBarrel_, hErrPOverlap_, hErrPEndcap_;
0044 MEP hErrPtBarrel_, hErrPtOverlap_, hErrPtEndcap_;
0045 MEP hErrDxy_, hErrDz_;
0046
0047 MEP hErrP_vs_Eta_, hErrPt_vs_Eta_, hErrQPt_vs_Eta_;
0048 MEP hErrP_vs_P_, hErrPt_vs_Pt_, hErrQPt_vs_Pt_, hErrEta_vs_Eta_;
0049
0050
0051 MEP hErrPt_PF_;
0052 MEP hErrQPt_PF_;
0053 MEP hdPt_vs_Eta_;
0054 MEP hdPt_vs_Pt_;
0055 MEP hPFMomAssCorrectness;
0056 MEP hPt_vs_PFMomAssCorrectness;
0057
0058
0059 MEP hNSimHits_;
0060 MEP hNSimToReco_, hNRecoToSim_;
0061
0062 MEP hNHits_, hNLostHits_, hNTrackerHits_, hNMuonHits_;
0063 MEP hNHits_vs_Pt_, hNHits_vs_Eta_;
0064 MEP hNLostHits_vs_Pt_, hNLostHits_vs_Eta_;
0065 MEP hNTrackerHits_vs_Pt_, hNTrackerHits_vs_Eta_;
0066 MEP hNMuonHits_vs_Pt_, hNMuonHits_vs_Eta_;
0067
0068
0069 MEP hPullPt_, hPullEta_, hPullPhi_, hPullQPt_, hPullDxy_, hPullDz_;
0070 MEP hPullPt_vs_Eta_, hPullPt_vs_Pt_, hPullEta_vs_Eta_, hPullPhi_vs_Eta_, hPullEta_vs_Pt_;
0071
0072
0073 MEP hNDof_, hChi2_, hChi2Norm_, hChi2Prob_;
0074 MEP hNDof_vs_Eta_, hChi2_vs_Eta_, hChi2Norm_vs_Eta_, hChi2Prob_vs_Eta_;
0075
0076 bool doAbsEta_;
0077 bool usePFMuon_;
0078
0079
0080
0081
0082 void bookHistos(DQMStore::IBooker& ibooker, const string& dirName, const HistoDimensions& hDim)
0083
0084 {
0085 ibooker.cd();
0086 ibooker.setCurrentFolder(dirName);
0087
0088 doAbsEta_ = hDim.doAbsEta;
0089 usePFMuon_ = hDim.usePFMuon;
0090
0091
0092 hP_ = ibooker.book1D("P", "p of recoTracks", hDim.nBinP, hDim.minP, hDim.maxP);
0093 hPt_ = ibooker.book1D("Pt", "p_{T} of recoTracks", hDim.nBinPt, hDim.minPt, hDim.maxPt);
0094 hEta_ = ibooker.book1D("Eta", "#eta of recoTracks", hDim.nBinEta, hDim.minEta, hDim.maxEta);
0095 hPhi_ = ibooker.book1D("Phi", "#phi of recoTracks", hDim.nBinPhi, hDim.minPhi, hDim.maxPhi);
0096
0097 hSimP_ = ibooker.book1D("SimP", "p of simTracks", hDim.nBinP, hDim.minP, hDim.maxP);
0098 hSimPt_ = ibooker.book1D("SimPt", "p_{T} of simTracks", hDim.nBinPt, hDim.minPt, hDim.maxPt);
0099 hSimEta_ = ibooker.book1D("SimEta", "#eta of simTracks", hDim.nBinEta, hDim.minEta, hDim.maxEta);
0100 hSimPhi_ = ibooker.book1D("SimPhi", "#phi of simTracks", hDim.nBinPhi, hDim.minPhi, hDim.maxPhi);
0101 hSimDxy_ = ibooker.book1D("SimDxy", "Dxy of simTracks", hDim.nBinDxy, hDim.minDxy, hDim.maxDxy);
0102 hSimDz_ = ibooker.book1D("SimDz", "Dz of simTracks", hDim.nBinDz, hDim.minDz, hDim.maxDz);
0103
0104
0105 hNSim_ = ibooker.book1D("NSim", "Number of particles per event", hDim.nTrks, -0.5, hDim.nTrks + 0.5);
0106 hNMuon_ = ibooker.book1D("NMuon", "Number of muons per event", hDim.nTrks, -0.5, hDim.nTrks + 0.5);
0107
0108
0109 hNTrks_ = ibooker.book1D("NTrks", "Number of reco tracks per event", hDim.nTrks, -0.5, hDim.nTrks + 0.5);
0110 hNTrksEta_ = ibooker.book1D("NTrksEta", "Number of reco tracks vs #eta", hDim.nBinEta, hDim.minEta, hDim.maxEta);
0111 hNTrksPt_ = ibooker.book1D("NTrksPt", "Number of reco tracks vs p_{T}", hDim.nBinPt, hDim.minPt, hDim.maxPt);
0112
0113 hMisQPt_ = ibooker.book1D("MisQPt", "Charge mis-id vs Pt", hDim.nBinPt, hDim.minPt, hDim.maxPt);
0114 hMisQEta_ = ibooker.book1D("MisQEta", "Charge mis-id vs Eta", hDim.nBinEta, hDim.minEta, hDim.maxEta);
0115
0116
0117 hErrP_ = ibooker.book1D("ErrP", "#Delta(p)/p", hDim.nBinErr, hDim.minErrP, hDim.maxErrP);
0118 hErrPBarrel_ = ibooker.book1D("ErrP_barrel", "#Delta(p)/p", hDim.nBinErr, hDim.minErrP, hDim.maxErrP);
0119 hErrPOverlap_ = ibooker.book1D("ErrP_overlap", "#Delta(p)/p", hDim.nBinErr, hDim.minErrP, hDim.maxErrP);
0120 hErrPEndcap_ = ibooker.book1D("ErrP_endcap", "#Delta(p)/p", hDim.nBinErr, hDim.minErrP, hDim.maxErrP);
0121 hErrPt_ = ibooker.book1D("ErrPt", "#Delta(p_{T})/p_{T}", hDim.nBinErr, hDim.minErrPt, hDim.maxErrPt);
0122 hErrPtBarrel_ = ibooker.book1D("ErrPt_barrel", "#Delta(p_{T})/p_{T}", hDim.nBinErr, hDim.minErrPt, hDim.maxErrPt);
0123 hErrPtOverlap_ = ibooker.book1D("ErrPt_overlap", "#Delta(p_{T})/p_{T}", hDim.nBinErr, hDim.minErrPt, hDim.maxErrPt);
0124 hErrPtEndcap_ = ibooker.book1D("ErrPt_endcap", "#Delta(p_{T})/p_{T}", hDim.nBinErr, hDim.minErrPt, hDim.maxErrPt);
0125 hErrEta_ = ibooker.book1D("ErrEta", "#sigma(#eta))", hDim.nBinErr, hDim.minErrEta, hDim.maxErrEta);
0126 hErrPhi_ = ibooker.book1D("ErrPhi", "#sigma(#phi)", hDim.nBinErr, hDim.minErrPhi, hDim.maxErrPhi);
0127 hErrDxy_ = ibooker.book1D("ErrDxy", "#sigma(d_{xy})", hDim.nBinErr, hDim.minErrDxy, hDim.maxErrDxy);
0128 hErrDz_ = ibooker.book1D("ErrDz", "#sigma(d_{z})", hDim.nBinErr, hDim.minErrDz, hDim.maxErrDz);
0129
0130
0131 if (usePFMuon_) {
0132 hErrPt_PF_ = ibooker.book1D("ErrPt_PF", "#Delta(p_{T})|_{PF}/p_{T}", hDim.nBinErr, hDim.minErrPt, hDim.maxErrPt);
0133 hErrQPt_PF_ =
0134 ibooker.book1D("ErrQPt_PF", "#Delta(q/p_{T})|_{PF}/(q/p_{T})", hDim.nBinErr, hDim.minErrQPt, hDim.maxErrQPt);
0135
0136 hPFMomAssCorrectness =
0137 ibooker.book1D("hPFMomAssCorrectness", "Corrected momentum assignement PF/RECO", 2, 0.5, 2.5);
0138 hPt_vs_PFMomAssCorrectness = ibooker.book2D("hPt_vs_PFMomAssCorrectness",
0139 "Corrected momentum assignement PF/RECO",
0140 hDim.nBinPt,
0141 hDim.minPt,
0142 hDim.maxP,
0143 2,
0144 0.5,
0145 2.5);
0146
0147 hdPt_vs_Pt_ = ibooker.book2D("dPt_vs_Pt",
0148 "#Delta(p_{T}) vs p_{T}",
0149 hDim.nBinPt,
0150 hDim.minPt,
0151 hDim.maxPt,
0152 hDim.nBinErr,
0153 hDim.minErrPt,
0154 hDim.maxErrPt);
0155 hdPt_vs_Eta_ = ibooker.book2D("dPt_vs_Eta",
0156 "#Delta(p_{T}) vs #eta",
0157 hDim.nBinEta,
0158 hDim.minEta,
0159 hDim.maxEta,
0160 hDim.nBinErr,
0161 hDim.minErrPt,
0162 hDim.maxErrPt);
0163 }
0164
0165
0166 hErrP_vs_Eta_ = ibooker.book2D("ErrP_vs_Eta",
0167 "#Delta(p)/p vs #eta",
0168 hDim.nBinEta,
0169 hDim.minEta,
0170 hDim.maxEta,
0171 hDim.nBinErr,
0172 hDim.minErrP,
0173 hDim.maxErrP);
0174 hErrPt_vs_Eta_ = ibooker.book2D("ErrPt_vs_Eta",
0175 "#Delta(p_{T})/p_{T} vs #eta",
0176 hDim.nBinEta,
0177 hDim.minEta,
0178 hDim.maxEta,
0179 hDim.nBinErr,
0180 hDim.minErrPt,
0181 hDim.maxErrPt);
0182 hErrQPt_vs_Eta_ = ibooker.book2D("ErrQPt_vs_Eta",
0183 "#Delta(q/p_{T})/(q/p_{T}) vs #eta",
0184 hDim.nBinEta,
0185 hDim.minEta,
0186 hDim.maxEta,
0187 hDim.nBinErr,
0188 hDim.minErrQPt,
0189 hDim.maxErrQPt);
0190 hErrEta_vs_Eta_ = ibooker.book2D("ErrEta_vs_Eta",
0191 "#sigma(#eta) vs #eta",
0192 hDim.nBinEta,
0193 hDim.minEta,
0194 hDim.maxEta,
0195 hDim.nBinErr,
0196 hDim.minErrEta,
0197 hDim.maxErrEta);
0198
0199
0200 hErrP_vs_P_ = ibooker.book2D(
0201 "ErrP_vs_P", "#Delta(p)/p vs p", hDim.nBinP, hDim.minP, hDim.maxP, hDim.nBinErr, hDim.minErrP, hDim.maxErrP);
0202 hErrPt_vs_Pt_ = ibooker.book2D("ErrPt_vs_Pt",
0203 "#Delta(p_{T})/p_{T} vs p_{T}",
0204 hDim.nBinPt,
0205 hDim.minPt,
0206 hDim.maxPt,
0207 hDim.nBinErr,
0208 hDim.minErrPt,
0209 hDim.maxErrPt);
0210 hErrQPt_vs_Pt_ = ibooker.book2D("ErrQPt_vs_Pt",
0211 "#Delta(q/p_{T})/(q/p_{T}) vs p_{T}",
0212 hDim.nBinPt,
0213 hDim.minPt,
0214 hDim.maxPt,
0215 hDim.nBinErr,
0216 hDim.minErrQPt,
0217 hDim.maxErrQPt);
0218
0219
0220 hPullPt_ = ibooker.book1D("PullPt", "Pull(#p_{T})", hDim.nBinPull, -hDim.wPull, hDim.wPull);
0221 hPullEta_ = ibooker.book1D("PullEta", "Pull(#eta)", hDim.nBinPull, -hDim.wPull, hDim.wPull);
0222 hPullPhi_ = ibooker.book1D("PullPhi", "Pull(#phi)", hDim.nBinPull, -hDim.wPull, hDim.wPull);
0223 hPullQPt_ = ibooker.book1D("PullQPt", "Pull(q/p_{T})", hDim.nBinPull, -hDim.wPull, hDim.wPull);
0224 hPullDxy_ = ibooker.book1D("PullDxy", "Pull(D_{xy})", hDim.nBinPull, -hDim.wPull, hDim.wPull);
0225 hPullDz_ = ibooker.book1D("PullDz", "Pull(D_{z})", hDim.nBinPull, -hDim.wPull, hDim.wPull);
0226
0227
0228 hPullPt_vs_Eta_ = ibooker.book2D("PullPt_vs_Eta",
0229 "Pull(p_{T}) vs #eta",
0230 hDim.nBinEta,
0231 hDim.minEta,
0232 hDim.maxEta,
0233 hDim.nBinPull,
0234 -hDim.wPull,
0235 hDim.wPull);
0236 hPullEta_vs_Eta_ = ibooker.book2D("PullEta_vs_Eta",
0237 "Pull(#eta) vs #eta",
0238 hDim.nBinEta,
0239 hDim.minEta,
0240 hDim.maxEta,
0241 hDim.nBinPull,
0242 -hDim.wPull,
0243 hDim.wPull);
0244 hPullPhi_vs_Eta_ = ibooker.book2D("PullPhi_vs_Eta",
0245 "Pull(#phi) vs #eta",
0246 hDim.nBinEta,
0247 hDim.minEta,
0248 hDim.maxEta,
0249 hDim.nBinPull,
0250 -hDim.wPull,
0251 hDim.wPull);
0252
0253
0254 hPullPt_vs_Pt_ = ibooker.book2D("PullPt_vs_Pt",
0255 "Pull(p_{T}) vs p_{T}",
0256 hDim.nBinPt,
0257 hDim.minPt,
0258 hDim.maxPt,
0259 hDim.nBinPull,
0260 -hDim.wPull,
0261 hDim.wPull);
0262 hPullEta_vs_Pt_ = ibooker.book2D("PullEta_vs_Pt",
0263 "Pull(#eta) vs p_{T}",
0264 hDim.nBinPt,
0265 hDim.minPt,
0266 hDim.maxPt,
0267 hDim.nBinPull,
0268 -hDim.wPull,
0269 hDim.wPull);
0270
0271
0272 const int nHits = 100;
0273 hNHits_ = ibooker.book1D("NHits", "Number of hits", nHits + 1, -0.5, nHits + 0.5);
0274 hNHits_vs_Pt_ = ibooker.book2D("NHits_vs_Pt",
0275 "Number of hits vs p_{T}",
0276 hDim.nBinPt,
0277 hDim.minPt,
0278 hDim.maxPt,
0279 nHits / 4 + 1,
0280 -0.25,
0281 nHits + 0.25);
0282 hNHits_vs_Eta_ = ibooker.book2D("NHits_vs_Eta",
0283 "Number of hits vs #eta",
0284 hDim.nBinEta,
0285 hDim.minEta,
0286 hDim.maxEta,
0287 nHits / 4 + 1,
0288 -0.25,
0289 nHits + 0.25);
0290 hNSimHits_ = ibooker.book1D("NSimHits", "Number of simHits", nHits + 1, -0.5, nHits + 0.5);
0291
0292 const int nLostHits = 5;
0293 hNLostHits_ = ibooker.book1D("NLostHits", "Number of Lost hits", nLostHits + 1, -0.5, nLostHits + 0.5);
0294 hNLostHits_vs_Pt_ = ibooker.book2D("NLostHits_vs_Pt",
0295 "Number of lost Hits vs p_{T}",
0296 hDim.nBinPt,
0297 hDim.minPt,
0298 hDim.maxPt,
0299 nLostHits + 1,
0300 -0.5,
0301 nLostHits + 0.5);
0302 hNLostHits_vs_Eta_ = ibooker.book2D("NLostHits_vs_Eta",
0303 "Number of lost Hits vs #eta",
0304 hDim.nBinEta,
0305 hDim.minEta,
0306 hDim.maxEta,
0307 nLostHits + 1,
0308 -0.5,
0309 nLostHits + 0.5);
0310
0311 const int nTrackerHits = 40;
0312 hNTrackerHits_ =
0313 ibooker.book1D("NTrackerHits", "Number of valid tracker hits", nTrackerHits + 1, -0.5, nTrackerHits + 0.5);
0314 hNTrackerHits_vs_Pt_ = ibooker.book2D("NTrackerHits_vs_Pt",
0315 "Number of valid traker hits vs p_{T}",
0316 hDim.nBinPt,
0317 hDim.minPt,
0318 hDim.maxPt,
0319 nTrackerHits / 4 + 1,
0320 -0.25,
0321 nTrackerHits + 0.25);
0322 hNTrackerHits_vs_Eta_ = ibooker.book2D("NTrackerHits_vs_Eta",
0323 "Number of valid tracker hits vs #eta",
0324 hDim.nBinEta,
0325 hDim.minEta,
0326 hDim.maxEta,
0327 nTrackerHits / 4 + 1,
0328 -0.25,
0329 nTrackerHits + 0.25);
0330
0331 const int nMuonHits = 60;
0332 hNMuonHits_ = ibooker.book1D("NMuonHits", "Number of valid muon hits", nMuonHits + 1, -0.5, nMuonHits + 0.5);
0333 hNMuonHits_vs_Pt_ = ibooker.book2D("NMuonHits_vs_Pt",
0334 "Number of valid muon hits vs p_{T}",
0335 hDim.nBinPt,
0336 hDim.minPt,
0337 hDim.maxPt,
0338 nMuonHits / 4 + 1,
0339 -0.25,
0340 nMuonHits + 0.25);
0341 hNMuonHits_vs_Eta_ = ibooker.book2D("NMuonHits_vs_Eta",
0342 "Number of valid muon hits vs #eta",
0343 hDim.nBinEta,
0344 hDim.minEta,
0345 hDim.maxEta,
0346 nMuonHits / 4 + 1,
0347 -0.25,
0348 nMuonHits + 0.25);
0349
0350 hNDof_ = ibooker.book1D("NDof", "Number of DoF", hDim.nDof + 1, -0.5, hDim.nDof + 0.5);
0351 hChi2_ = ibooker.book1D("Chi2", "#Chi^{2}", hDim.nBinErr, 0, 200);
0352 hChi2Norm_ = ibooker.book1D("Chi2Norm", "Normalized #Chi^{2}", hDim.nBinErr, 0, 50);
0353 hChi2Prob_ = ibooker.book1D("Chi2Prob", "Prob(#Chi^{2})", hDim.nBinErr, 0, 1);
0354
0355 hNDof_vs_Eta_ = ibooker.book2D("NDof_vs_Eta",
0356 "Number of DoF vs #eta",
0357 hDim.nBinEta,
0358 hDim.minEta,
0359 hDim.maxEta,
0360 hDim.nDof + 1,
0361 -0.5,
0362 hDim.nDof + 0.5);
0363 hChi2_vs_Eta_ = ibooker.book2D(
0364 "Chi2_vs_Eta", "#Chi^{2} vs #eta", hDim.nBinEta, hDim.minEta, hDim.maxEta, hDim.nBinErr, 0., 200.);
0365 hChi2Norm_vs_Eta_ = ibooker.book2D(
0366 "Chi2Norm_vs_Eta", "Normalized #Chi^{2} vs #eta", hDim.nBinEta, hDim.minEta, hDim.maxEta, hDim.nBinErr, 0., 50.);
0367 hChi2Prob_vs_Eta_ = ibooker.book2D(
0368 "Chi2Prob_vs_Eta", "Prob(#Chi^{2}) vs #eta", hDim.nBinEta, hDim.minEta, hDim.maxEta, hDim.nBinErr, 0., 1.);
0369
0370 hNSimToReco_ =
0371 ibooker.book1D("NSimToReco", "Number of associated reco tracks", hDim.nAssoc + 1, -0.5, hDim.nAssoc + 0.5);
0372 hNRecoToSim_ =
0373 ibooker.book1D("NRecoToSim", "Number of associated sim TP's", hDim.nAssoc + 1, -0.5, hDim.nAssoc + 0.5);
0374 };
0375
0376
0377
0378
0379 void fill(const TrackingParticle* simRef, const Muon* muonRef) {
0380 const double simP = simRef->p();
0381 const double simPt = simRef->pt();
0382 const double simEta = doAbsEta_ ? fabs(simRef->eta()) : simRef->eta();
0383 const double simPhi = simRef->phi();
0384 const double simQ = simRef->charge();
0385 const double simQPt = simQ / simPt;
0386
0387 GlobalPoint simVtx(simRef->vertex().x(), simRef->vertex().y(), simRef->vertex().z());
0388 GlobalVector simMom(simRef->momentum().x(), simRef->momentum().y(), simRef->momentum().z());
0389 const double simDxy = -simVtx.x() * sin(simPhi) + simVtx.y() * cos(simPhi);
0390 const double simDz = simVtx.z() - (simVtx.x() * simMom.x() + simVtx.y() * simMom.y()) * simMom.z() / simMom.perp2();
0391
0392 const double recoQ = muonRef->charge();
0393 if (simQ * recoQ < 0) {
0394 hMisQPt_->Fill(simPt);
0395 hMisQEta_->Fill(simEta);
0396 }
0397
0398 double recoP, recoPt, recoEta, recoPhi, recoQPt;
0399 if (usePFMuon_) {
0400
0401 const double origRecoPt = muonRef->pt();
0402
0403
0404 const double origRecoQPt = recoQ / origRecoPt;
0405 recoP = muonRef->pfP4().P();
0406 recoPt = muonRef->pfP4().Pt();
0407 recoEta = muonRef->pfP4().Eta();
0408 recoPhi = muonRef->pfP4().Phi();
0409 recoQPt = recoQ / recoPt;
0410 hErrPt_PF_->Fill((recoPt - origRecoPt) / origRecoPt);
0411 hErrQPt_PF_->Fill((recoQPt - origRecoQPt) / origRecoQPt);
0412
0413 hdPt_vs_Eta_->Fill(recoEta, recoPt - origRecoPt);
0414 hdPt_vs_Pt_->Fill(recoPt, recoPt - origRecoPt);
0415
0416 int theCorrectPFAss = (fabs(recoPt - simPt) < fabs(origRecoPt - simPt)) ? 1 : 2;
0417 hPFMomAssCorrectness->Fill(theCorrectPFAss);
0418 hPt_vs_PFMomAssCorrectness->Fill(simPt, theCorrectPFAss);
0419 }
0420
0421 else {
0422 recoP = muonRef->p();
0423 recoPt = muonRef->pt();
0424 recoEta = muonRef->eta();
0425 recoPhi = muonRef->phi();
0426 recoQPt = recoQ / recoPt;
0427 }
0428
0429 const double errP = (recoP - simP) / simP;
0430 const double errPt = (recoPt - simPt) / simPt;
0431 const double errEta = (recoEta - simEta) / simEta;
0432 const double errPhi = (recoPhi - simPhi) / simPhi;
0433 const double errQPt = (recoQPt - simQPt) / simQPt;
0434
0435 hP_->Fill(simP);
0436 hPt_->Fill(simPt);
0437 hEta_->Fill(simEta);
0438 hPhi_->Fill(simPhi);
0439
0440 hErrP_->Fill(errP);
0441 hErrPt_->Fill(errPt);
0442 hErrEta_->Fill(errEta);
0443 hErrPhi_->Fill(errPhi);
0444
0445 if (fabs(simEta) > 0. && fabs(simEta) < 0.8) {
0446 hErrPBarrel_->Fill(errP);
0447 hErrPtBarrel_->Fill(errPt);
0448 } else if (fabs(simEta) > 0.8 && fabs(simEta) < 1.2) {
0449 hErrPOverlap_->Fill(errP);
0450 hErrPtOverlap_->Fill(errPt);
0451 } else if (fabs(simEta) > 1.2) {
0452 hErrPEndcap_->Fill(errP);
0453 hErrPtEndcap_->Fill(errPt);
0454 }
0455
0456 hErrP_vs_Eta_->Fill(simEta, errP);
0457 hErrPt_vs_Eta_->Fill(simEta, errPt);
0458 hErrQPt_vs_Eta_->Fill(simEta, errQPt);
0459
0460 hErrP_vs_P_->Fill(simP, errP);
0461 hErrPt_vs_Pt_->Fill(simPt, errPt);
0462 hErrQPt_vs_Pt_->Fill(simQPt, errQPt);
0463
0464 hErrEta_vs_Eta_->Fill(simEta, errEta);
0465
0466
0467 reco::TrackRef recoRef = muonRef->track();
0468 if (recoRef.isNonnull()) {
0469
0470 const int nRecoHits = recoRef->numberOfValidHits();
0471 const int nLostHits = recoRef->numberOfLostHits();
0472
0473 hNHits_->Fill(nRecoHits);
0474 hNHits_vs_Pt_->Fill(simPt, nRecoHits);
0475 hNHits_vs_Eta_->Fill(simEta, nRecoHits);
0476
0477 hNLostHits_->Fill(nLostHits);
0478 hNLostHits_vs_Pt_->Fill(simPt, nLostHits);
0479 hNLostHits_vs_Eta_->Fill(simEta, nLostHits);
0480
0481 const double recoNDof = recoRef->ndof();
0482 const double recoChi2 = recoRef->chi2();
0483 const double recoChi2Norm = recoRef->normalizedChi2();
0484 const double recoChi2Prob = TMath::Prob(recoRef->chi2(), static_cast<int>(recoRef->ndof()));
0485
0486 hNDof_->Fill(recoNDof);
0487 hChi2_->Fill(recoChi2);
0488 hChi2Norm_->Fill(recoChi2Norm);
0489 hChi2Prob_->Fill(recoChi2Prob);
0490
0491 hNDof_vs_Eta_->Fill(simEta, recoNDof);
0492 hChi2_vs_Eta_->Fill(simEta, recoChi2);
0493 hChi2Norm_vs_Eta_->Fill(simEta, recoChi2Norm);
0494 hChi2Prob_vs_Eta_->Fill(simEta, recoChi2Prob);
0495
0496 const double recoDxy = recoRef->dxy();
0497 const double recoDz = recoRef->dz();
0498
0499 const double errDxy = (recoDxy - simDxy) / simDxy;
0500 const double errDz = (recoDz - simDz) / simDz;
0501 hErrDxy_->Fill(errDxy);
0502 hErrDz_->Fill(errDz);
0503
0504 const double pullPt = (recoPt - simPt) / recoRef->ptError();
0505 const double pullQPt = (recoQPt - simQPt) / recoRef->qoverpError();
0506 const double pullEta = (recoEta - simEta) / recoRef->etaError();
0507 const double pullPhi = (recoPhi - simPhi) / recoRef->phiError();
0508 const double pullDxy = (recoDxy - simDxy) / recoRef->dxyError();
0509 const double pullDz = (recoDz - simDz) / recoRef->dzError();
0510
0511 hPullPt_->Fill(pullPt);
0512 hPullEta_->Fill(pullEta);
0513 hPullPhi_->Fill(pullPhi);
0514 hPullQPt_->Fill(pullQPt);
0515 hPullDxy_->Fill(pullDxy);
0516 hPullDz_->Fill(pullDz);
0517
0518 hPullPt_vs_Eta_->Fill(simEta, pullPt);
0519 hPullPt_vs_Pt_->Fill(simPt, pullPt);
0520
0521 hPullEta_vs_Eta_->Fill(simEta, pullEta);
0522 hPullPhi_vs_Eta_->Fill(simEta, pullPhi);
0523
0524 hPullEta_vs_Pt_->Fill(simPt, pullEta);
0525 }
0526 };
0527 };
0528
0529
0530
0531
0532 struct RecoDisplacedMuonValidator::CommonME {
0533 typedef MonitorElement* MEP;
0534
0535
0536 MEP hTrkToGlbDiffNTrackerHits_, hStaToGlbDiffNMuonHits_;
0537 MEP hTrkToGlbDiffNTrackerHitsEta_, hStaToGlbDiffNMuonHitsEta_;
0538 MEP hTrkToGlbDiffNTrackerHitsPt_, hStaToGlbDiffNMuonHitsPt_;
0539
0540
0541 MEP hNInvalidHitsGTHitPattern_, hNInvalidHitsITHitPattern_, hNInvalidHitsOTHitPattern_;
0542 MEP hNDeltaInvalidHitsHitPattern_;
0543
0544
0545 MEP hMuonP_, hMuonPt_, hMuonEta_, hMuonPhi_;
0546
0547 MEP hMuonTrackP_, hMuonTrackPt_, hMuonTrackEta_, hMuonTrackPhi_, hMuonTrackDxy_, hMuonTrackDz_;
0548
0549 MEP hMuonAllP_, hMuonAllPt_, hMuonAllEta_, hMuonAllPhi_;
0550 };
0551
0552
0553
0554
0555 RecoDisplacedMuonValidator::RecoDisplacedMuonValidator(const edm::ParameterSet& pset)
0556 : selector_(pset.getParameter<std::string>("selection")) {
0557
0558 edm::LogVerbatim("RecoDisplacedMuonValidator") << "constructing RecoDisplacedMuonValidator: " << pset.dump();
0559 verbose_ = pset.getUntrackedParameter<unsigned int>("verbose", 0);
0560
0561 outputFileName_ = pset.getUntrackedParameter<string>("outputFileName", "");
0562
0563 wantTightMuon_ = pset.getParameter<bool>("wantTightMuon");
0564 beamspotLabel_ = pset.getParameter<edm::InputTag>("beamSpot");
0565 primvertexLabel_ = pset.getParameter<edm::InputTag>("primaryVertex");
0566 beamspotToken_ = consumes<reco::BeamSpot>(beamspotLabel_);
0567 primvertexToken_ = consumes<reco::VertexCollection>(primvertexLabel_);
0568
0569
0570
0571 hDim.nBinP = pset.getUntrackedParameter<unsigned int>("nBinP");
0572 hDim.minP = pset.getUntrackedParameter<double>("minP");
0573 hDim.maxP = pset.getUntrackedParameter<double>("maxP");
0574
0575 hDim.nBinPt = pset.getUntrackedParameter<unsigned int>("nBinPt");
0576 hDim.minPt = pset.getUntrackedParameter<double>("minPt");
0577 hDim.maxPt = pset.getUntrackedParameter<double>("maxPt");
0578
0579 doAbsEta_ = pset.getUntrackedParameter<bool>("doAbsEta");
0580 hDim.doAbsEta = doAbsEta_;
0581 hDim.nBinEta = pset.getUntrackedParameter<unsigned int>("nBinEta");
0582 hDim.minEta = pset.getUntrackedParameter<double>("minEta");
0583 hDim.maxEta = pset.getUntrackedParameter<double>("maxEta");
0584
0585 hDim.nBinDxy = pset.getUntrackedParameter<unsigned int>("nBinDxy");
0586 hDim.minDxy = pset.getUntrackedParameter<double>("minDxy");
0587 hDim.maxDxy = pset.getUntrackedParameter<double>("maxDxy");
0588
0589 hDim.nBinDz = pset.getUntrackedParameter<unsigned int>("nBinDz");
0590 hDim.minDz = pset.getUntrackedParameter<double>("minDz");
0591 hDim.maxDz = pset.getUntrackedParameter<double>("maxDz");
0592
0593 hDim.nBinPhi = pset.getUntrackedParameter<unsigned int>("nBinPhi");
0594 hDim.minPhi = pset.getUntrackedParameter<double>("minPhi", -TMath::Pi());
0595 hDim.maxPhi = pset.getUntrackedParameter<double>("maxPhi", TMath::Pi());
0596
0597 hDim.nBinErr = pset.getUntrackedParameter<unsigned int>("nBinErr");
0598 hDim.nBinPull = pset.getUntrackedParameter<unsigned int>("nBinPull");
0599
0600 hDim.wPull = pset.getUntrackedParameter<double>("wPull");
0601
0602 hDim.minErrP = pset.getUntrackedParameter<double>("minErrP");
0603 hDim.maxErrP = pset.getUntrackedParameter<double>("maxErrP");
0604
0605 hDim.minErrPt = pset.getUntrackedParameter<double>("minErrPt");
0606 hDim.maxErrPt = pset.getUntrackedParameter<double>("maxErrPt");
0607
0608 hDim.minErrQPt = pset.getUntrackedParameter<double>("minErrQPt");
0609 hDim.maxErrQPt = pset.getUntrackedParameter<double>("maxErrQPt");
0610
0611 hDim.minErrEta = pset.getUntrackedParameter<double>("minErrEta");
0612 hDim.maxErrEta = pset.getUntrackedParameter<double>("maxErrEta");
0613
0614 hDim.minErrPhi = pset.getUntrackedParameter<double>("minErrPhi");
0615 hDim.maxErrPhi = pset.getUntrackedParameter<double>("maxErrPhi");
0616
0617 hDim.minErrDxy = pset.getUntrackedParameter<double>("minErrDxy");
0618 hDim.maxErrDxy = pset.getUntrackedParameter<double>("maxErrDxy");
0619
0620 hDim.minErrDz = pset.getUntrackedParameter<double>("minErrDz");
0621 hDim.maxErrDz = pset.getUntrackedParameter<double>("maxErrDz");
0622
0623 hDim.nTrks = pset.getUntrackedParameter<unsigned int>("nTrks");
0624 hDim.nAssoc = pset.getUntrackedParameter<unsigned int>("nAssoc");
0625 hDim.nDof = pset.getUntrackedParameter<unsigned int>("nDof", 55);
0626
0627
0628 simLabel_ = pset.getParameter<InputTag>("simLabel");
0629 tpRefVector = pset.getParameter<bool>("tpRefVector");
0630 if (tpRefVector)
0631 tpRefVectorToken_ = consumes<TrackingParticleRefVector>(simLabel_);
0632 else
0633 simToken_ = consumes<TrackingParticleCollection>(simLabel_);
0634
0635 muonLabel_ = pset.getParameter<InputTag>("muonLabel");
0636 muonToken_ = consumes<edm::View<reco::Muon> >(muonLabel_);
0637
0638
0639 doAssoc_ = pset.getUntrackedParameter<bool>("doAssoc", true);
0640 muAssocLabel_ = pset.getParameter<InputTag>("muAssocLabel");
0641 if (doAssoc_) {
0642 muAssocToken_ = consumes<reco::MuonToTrackingParticleAssociator>(muAssocLabel_);
0643 }
0644
0645
0646 usePFMuon_ = pset.getUntrackedParameter<bool>("usePFMuon");
0647 hDim.usePFMuon = usePFMuon_;
0648
0649
0650 std::string trackType = pset.getParameter<std::string>("trackType");
0651 if (trackType == "inner")
0652 trackType_ = reco::InnerTk;
0653 else if (trackType == "outer")
0654 trackType_ = reco::OuterTk;
0655 else if (trackType == "global")
0656 trackType_ = reco::GlobalTk;
0657 else if (trackType == "segments")
0658 trackType_ = reco::Segments;
0659 else
0660 throw cms::Exception("Configuration") << "Track type '" << trackType << "' not supported.\n";
0661
0662
0663
0664 ParameterSet tpset = pset.getParameter<ParameterSet>("tpSelector");
0665 tpSelector_ = TrackingParticleSelector(tpset.getParameter<double>("ptMin"),
0666 tpset.getParameter<double>("ptMax"),
0667 tpset.getParameter<double>("minRapidity"),
0668 tpset.getParameter<double>("maxRapidity"),
0669 tpset.getParameter<double>("tip"),
0670 tpset.getParameter<double>("lip"),
0671 tpset.getParameter<int>("minHit"),
0672 tpset.getParameter<bool>("signalOnly"),
0673 tpset.getParameter<bool>("intimeOnly"),
0674 tpset.getParameter<bool>("chargedOnly"),
0675 tpset.getParameter<bool>("stableOnly"),
0676 tpset.getParameter<std::vector<int> >("pdgId"));
0677
0678
0679 ParameterSet serviceParameters = pset.getParameter<ParameterSet>("ServiceParameters");
0680
0681
0682 dbe_ = Service<DQMStore>().operator->();
0683 subsystemname_ = pset.getUntrackedParameter<std::string>("subSystemFolder", "YourSubsystem");
0684
0685 subDir_ = pset.getUntrackedParameter<string>("subDir");
0686 if (subDir_.empty())
0687 subDir_ = "RecoDisplacedMuonV";
0688 if (subDir_[subDir_.size() - 1] == '/')
0689 subDir_.erase(subDir_.size() - 1);
0690 }
0691
0692 void RecoDisplacedMuonValidator::bookHistograms(DQMStore::IBooker& ibooker,
0693 edm::Run const& iRun,
0694 edm::EventSetup const& ) {
0695
0696 ibooker.cd();
0697
0698 ibooker.setCurrentFolder(subDir_);
0699
0700 commonME_ = new CommonME;
0701 muonME_ = new MuonME;
0702
0703
0704 const int nHits = 100;
0705
0706
0707 commonME_->hTrkToGlbDiffNTrackerHits_ = ibooker.book1D("TrkGlbDiffNTrackerHits",
0708 "Difference of number of tracker hits (tkMuon - globalMuon)",
0709 2 * nHits + 1,
0710 -nHits - 0.5,
0711 nHits + 0.5);
0712 commonME_->hStaToGlbDiffNMuonHits_ = ibooker.book1D("StaGlbDiffNMuonHits",
0713 "Difference of number of muon hits (staMuon - globalMuon)",
0714 2 * nHits + 1,
0715 -nHits - 0.5,
0716 nHits + 0.5);
0717
0718 commonME_->hTrkToGlbDiffNTrackerHitsEta_ =
0719 ibooker.book2D("TrkGlbDiffNTrackerHitsEta",
0720 "Difference of number of tracker hits (tkMuon - globalMuon)",
0721 hDim.nBinEta,
0722 hDim.minEta,
0723 hDim.maxEta,
0724 2 * nHits + 1,
0725 -nHits - 0.5,
0726 nHits + 0.5);
0727 commonME_->hStaToGlbDiffNMuonHitsEta_ = ibooker.book2D("StaGlbDiffNMuonHitsEta",
0728 "Difference of number of muon hits (staMuon - globalMuon)",
0729 hDim.nBinEta,
0730 hDim.minEta,
0731 hDim.maxEta,
0732 2 * nHits + 1,
0733 -nHits - 0.5,
0734 nHits + 0.5);
0735
0736 commonME_->hTrkToGlbDiffNTrackerHitsPt_ = ibooker.book2D("TrkGlbDiffNTrackerHitsPt",
0737 "Difference of number of tracker hits (tkMuon - globalMuon)",
0738 hDim.nBinPt,
0739 hDim.minPt,
0740 hDim.maxPt,
0741 2 * nHits + 1,
0742 -nHits - 0.5,
0743 nHits + 0.5);
0744 commonME_->hStaToGlbDiffNMuonHitsPt_ = ibooker.book2D("StaGlbDiffNMuonHitsPt",
0745 "Difference of number of muon hits (staMuon - globalMuon)",
0746 hDim.nBinPt,
0747 hDim.minPt,
0748 hDim.maxPt,
0749 2 * nHits + 1,
0750 -nHits - 0.5,
0751 nHits + 0.5);
0752
0753
0754 commonME_->hNInvalidHitsGTHitPattern_ = ibooker.book1D(
0755 "NInvalidHitsGTHitPattern", "Number of invalid hits on a global track", nHits + 1, -0.5, nHits + 0.5);
0756 commonME_->hNInvalidHitsITHitPattern_ = ibooker.book1D(
0757 "NInvalidHitsITHitPattern", "Number of invalid hits on an inner track", nHits + 1, -0.5, nHits + 0.5);
0758 commonME_->hNInvalidHitsOTHitPattern_ = ibooker.book1D(
0759 "NInvalidHitsOTHitPattern", "Number of invalid hits on an outer track", nHits + 1, -0.5, nHits + 0.5);
0760 commonME_->hNDeltaInvalidHitsHitPattern_ =
0761 ibooker.book1D("hNDeltaInvalidHitsHitPattern",
0762 "The discrepancy for Number of invalid hits on an global track and inner and outer tracks",
0763 2 * nHits + 1,
0764 -nHits - 0.5,
0765 nHits + 0.5);
0766
0767
0768 commonME_->hMuonP_ = ibooker.book1D("PMuon", "p of muon", hDim.nBinP, hDim.minP, hDim.maxP);
0769 commonME_->hMuonPt_ = ibooker.book1D("PtMuon", "p_{T} of muon", hDim.nBinPt, hDim.minPt, hDim.maxPt);
0770 commonME_->hMuonEta_ = ibooker.book1D("EtaMuon", "#eta of muon", hDim.nBinEta, hDim.minEta, hDim.maxEta);
0771 commonME_->hMuonPhi_ = ibooker.book1D("PhiMuon", "#phi of muon", hDim.nBinPhi, hDim.minPhi, hDim.maxPhi);
0772
0773 commonME_->hMuonTrackP_ = ibooker.book1D("PMuonTrack", "p of reco muon track", hDim.nBinP, hDim.minP, hDim.maxP);
0774 commonME_->hMuonTrackPt_ =
0775 ibooker.book1D("PtMuonTrack", "p_{T} of reco muon track", hDim.nBinPt, hDim.minPt, hDim.maxPt);
0776 commonME_->hMuonTrackEta_ =
0777 ibooker.book1D("EtaMuonTrack", "#eta of reco muon track", hDim.nBinEta, hDim.minEta, hDim.maxEta);
0778 commonME_->hMuonTrackPhi_ =
0779 ibooker.book1D("PhiMuonTrack", "#phi of reco muon track", hDim.nBinPhi, hDim.minPhi, hDim.maxPhi);
0780 commonME_->hMuonTrackDxy_ =
0781 ibooker.book1D("DxyMuonTrack", "Dxy of reco muon track", hDim.nBinDxy, hDim.minDxy, hDim.maxDxy);
0782 commonME_->hMuonTrackDz_ =
0783 ibooker.book1D("DzMuonTrack", "Dz of reco muon track", hDim.nBinDz, hDim.minDz, hDim.maxDz);
0784
0785
0786 commonME_->hMuonAllP_ = ibooker.book1D("PMuonAll", "p of muons of all types", hDim.nBinP, hDim.minP, hDim.maxP);
0787 commonME_->hMuonAllPt_ =
0788 ibooker.book1D("PtMuonAll", "p_{T} of muon of all types", hDim.nBinPt, hDim.minPt, hDim.maxPt);
0789 commonME_->hMuonAllEta_ =
0790 ibooker.book1D("EtaMuonAll", "#eta of muon of all types", hDim.nBinEta, hDim.minEta, hDim.maxEta);
0791 commonME_->hMuonAllPhi_ =
0792 ibooker.book1D("PhiMuonAll", "#phi of muon of all types", hDim.nBinPhi, hDim.minPhi, hDim.maxPhi);
0793
0794 muonME_->bookHistos(ibooker, subDir_, hDim);
0795 }
0796
0797
0798
0799
0800 RecoDisplacedMuonValidator::~RecoDisplacedMuonValidator() {}
0801
0802
0803
0804
0805
0806 void RecoDisplacedMuonValidator::dqmBeginRun(const edm::Run&, const EventSetup& eventSetup) {}
0807
0808
0809
0810
0811 void RecoDisplacedMuonValidator::dqmEndRun(edm::Run const&, edm::EventSetup const&) {
0812 if (dbe_ && !outputFileName_.empty())
0813 dbe_->save(outputFileName_);
0814 }
0815
0816
0817
0818
0819 void RecoDisplacedMuonValidator::analyze(const Event& event, const EventSetup& eventSetup) {
0820
0821 reco::Vertex::Point posVtx;
0822 reco::Vertex::Error errVtx;
0823 edm::Handle<reco::VertexCollection> recVtxs;
0824 event.getByToken(primvertexToken_, recVtxs);
0825 unsigned int theIndexOfThePrimaryVertex = 999.;
0826 for (unsigned int ind = 0; ind < recVtxs->size(); ++ind) {
0827 if ((*recVtxs)[ind].isValid() && !((*recVtxs)[ind].isFake())) {
0828 theIndexOfThePrimaryVertex = ind;
0829 break;
0830 }
0831 }
0832 if (theIndexOfThePrimaryVertex < 100) {
0833 posVtx = ((*recVtxs)[theIndexOfThePrimaryVertex]).position();
0834 errVtx = ((*recVtxs)[theIndexOfThePrimaryVertex]).error();
0835 } else {
0836 LogInfo("RecoDisplacedMuonValidator") << "reco::PrimaryVertex not found, use BeamSpot position instead\n";
0837 edm::Handle<reco::BeamSpot> recoBeamSpotHandle;
0838 event.getByToken(beamspotToken_, recoBeamSpotHandle);
0839 reco::BeamSpot bs = *recoBeamSpotHandle;
0840 posVtx = bs.position();
0841 errVtx(0, 0) = bs.BeamWidthX();
0842 errVtx(1, 1) = bs.BeamWidthY();
0843 errVtx(2, 2) = bs.sigmaZ();
0844 }
0845 const reco::Vertex thePrimaryVertex(posVtx, errVtx);
0846
0847
0848 TrackingParticleRefVector TPrefV;
0849 const TrackingParticleRefVector* ptr_TPrefV = nullptr;
0850 Handle<TrackingParticleCollection> simHandle;
0851 Handle<TrackingParticleRefVector> TPCollectionRefVector_H;
0852
0853 if (tpRefVector) {
0854 event.getByToken(tpRefVectorToken_, TPCollectionRefVector_H);
0855 ptr_TPrefV = TPCollectionRefVector_H.product();
0856 TPrefV = *ptr_TPrefV;
0857 } else {
0858 event.getByToken(simToken_, simHandle);
0859 size_t nTP = simHandle->size();
0860 for (size_t i = 0; i < nTP; ++i) {
0861 TPrefV.push_back(TrackingParticleRef(simHandle, i));
0862 }
0863 ptr_TPrefV = &TPrefV;
0864 }
0865
0866
0867 Handle<edm::View<Muon> > muonHandle;
0868 event.getByToken(muonToken_, muonHandle);
0869 View<Muon> muonColl = *(muonHandle.product());
0870
0871 reco::MuonToTrackingParticleAssociator const* assoByHits = nullptr;
0872 if (doAssoc_) {
0873 edm::Handle<reco::MuonToTrackingParticleAssociator> associatorBase;
0874 event.getByToken(muAssocToken_, associatorBase);
0875 assoByHits = associatorBase.product();
0876 }
0877
0878 const size_t nSim = ptr_TPrefV->size();
0879
0880 edm::RefToBaseVector<reco::Muon> Muons;
0881 for (size_t i = 0; i < muonHandle->size(); ++i) {
0882 Muons.push_back(muonHandle->refAt(i));
0883 }
0884
0885 muonME_->hNSim_->Fill(nSim);
0886 muonME_->hNMuon_->Fill(muonColl.size());
0887
0888 reco::MuonToSimCollection muonToSimColl;
0889 reco::SimToMuonCollection simToMuonColl;
0890
0891 if (doAssoc_) {
0892 edm::LogVerbatim("RecoDisplacedMuonValidator")
0893 << "\n >>> MuonToSim association : " << muAssocLabel_ << " <<< \n"
0894 << " muon collection : " << muonLabel_ << " (size = " << muonHandle->size() << ") \n"
0895 << " TrackingParticle collection : " << simLabel_ << " (size = " << nSim << ")";
0896
0897 assoByHits->associateMuons(muonToSimColl, simToMuonColl, Muons, trackType_, TPrefV);
0898 } else {
0899
0900
0901
0902
0903
0904
0905
0906
0907
0908
0909
0910
0911
0912
0913
0914
0915
0916
0917
0918
0919
0920
0921
0922
0923
0924
0925
0926 }
0927
0928 int glbNTrackerHits = 0;
0929 int trkNTrackerHits = 0;
0930 int glbNMuonHits = 0;
0931 int staNMuonHits = 0;
0932 int NTrackerHits = 0;
0933 int NMuonHits = 0;
0934
0935
0936 for (View<Muon>::const_iterator iMuon = muonColl.begin(); iMuon != muonColl.end(); ++iMuon) {
0937 double muonP, muonPt, muonEta, muonPhi;
0938 if (usePFMuon_) {
0939 muonP = iMuon->pfP4().P();
0940 muonPt = iMuon->pfP4().Pt();
0941 muonEta = iMuon->pfP4().Eta();
0942 muonPhi = iMuon->pfP4().Phi();
0943 } else {
0944 muonP = iMuon->p();
0945 muonPt = iMuon->pt();
0946 muonEta = iMuon->eta();
0947 muonPhi = iMuon->phi();
0948 }
0949
0950
0951 commonME_->hMuonAllP_->Fill(muonP);
0952 commonME_->hMuonAllPt_->Fill(muonPt);
0953 commonME_->hMuonAllEta_->Fill(muonEta);
0954 commonME_->hMuonAllPhi_->Fill(muonPhi);
0955
0956 if (!selector_(*iMuon))
0957 continue;
0958 if (wantTightMuon_) {
0959 if (!muon::isTightMuon(*iMuon, thePrimaryVertex))
0960 continue;
0961 }
0962
0963 TrackRef Track = iMuon->track();
0964
0965 if (trackType_ == reco::InnerTk) {
0966 Track = iMuon->track();
0967 trkNTrackerHits = countTrackerHits(*Track);
0968 } else if (trackType_ == reco::OuterTk) {
0969 Track = iMuon->standAloneMuon();
0970 } else if (trackType_ == reco::GlobalTk) {
0971 Track = iMuon->combinedMuon();
0972 glbNTrackerHits = countTrackerHits(*Track);
0973 glbNMuonHits = countMuonHits(*Track);
0974 }
0975
0976 commonME_->hMuonTrackP_->Fill(Track->p());
0977 commonME_->hMuonTrackPt_->Fill(Track->pt());
0978 commonME_->hMuonTrackEta_->Fill(Track->eta());
0979 commonME_->hMuonTrackPhi_->Fill(Track->phi());
0980
0981 commonME_->hMuonTrackDxy_->Fill(Track->dxy());
0982 commonME_->hMuonTrackDz_->Fill(Track->dz());
0983
0984 NTrackerHits = countTrackerHits(*Track);
0985 muonME_->hNTrackerHits_->Fill(NTrackerHits);
0986 muonME_->hNTrackerHits_vs_Pt_->Fill(Track->pt(), NTrackerHits);
0987 muonME_->hNTrackerHits_vs_Eta_->Fill(Track->eta(), NTrackerHits);
0988
0989 NMuonHits = countMuonHits(*Track);
0990 muonME_->hNMuonHits_->Fill(NMuonHits);
0991 muonME_->hNMuonHits_vs_Pt_->Fill(Track->pt(), NMuonHits);
0992 muonME_->hNMuonHits_vs_Eta_->Fill(Track->eta(), NMuonHits);
0993
0994
0995
0996
0997 muonME_->hNTrksEta_->Fill(Track->eta());
0998 muonME_->hNTrksPt_->Fill(Track->pt());
0999
1000 commonME_->hMuonP_->Fill(muonP);
1001 commonME_->hMuonPt_->Fill(muonPt);
1002 commonME_->hMuonEta_->Fill(muonEta);
1003 commonME_->hMuonPhi_->Fill(muonPhi);
1004
1005 if (iMuon->isGlobalMuon()) {
1006 double gtHitPat = iMuon->globalTrack()->hitPattern().numberOfAllHits(HitPattern::TRACK_HITS) -
1007 iMuon->globalTrack()->hitPattern().numberOfValidHits();
1008
1009 double itHitPat = iMuon->innerTrack()->hitPattern().numberOfAllHits(HitPattern::TRACK_HITS) -
1010 iMuon->innerTrack()->hitPattern().numberOfValidHits();
1011
1012 double otHitPat = iMuon->outerTrack()->hitPattern().numberOfAllHits(HitPattern::TRACK_HITS) -
1013 iMuon->outerTrack()->hitPattern().numberOfValidHits();
1014
1015 commonME_->hNInvalidHitsGTHitPattern_->Fill(gtHitPat);
1016 commonME_->hNInvalidHitsITHitPattern_->Fill(itHitPat);
1017 commonME_->hNInvalidHitsOTHitPattern_->Fill(otHitPat);
1018 commonME_->hNDeltaInvalidHitsHitPattern_->Fill(gtHitPat - itHitPat - otHitPat);
1019
1020
1021 if (iMuon->isStandAloneMuon()) {
1022 commonME_->hStaToGlbDiffNMuonHitsEta_->Fill(Track->eta(), staNMuonHits - glbNMuonHits);
1023 commonME_->hStaToGlbDiffNMuonHitsPt_->Fill(Track->pt(), staNMuonHits - glbNMuonHits);
1024 commonME_->hStaToGlbDiffNMuonHits_->Fill(staNMuonHits - glbNMuonHits);
1025 }
1026
1027
1028 if (iMuon->isTrackerMuon()) {
1029 commonME_->hTrkToGlbDiffNTrackerHitsEta_->Fill(Track->eta(), trkNTrackerHits - glbNTrackerHits);
1030 commonME_->hTrkToGlbDiffNTrackerHitsPt_->Fill(Track->pt(), trkNTrackerHits - glbNTrackerHits);
1031 commonME_->hTrkToGlbDiffNTrackerHits_->Fill(trkNTrackerHits - glbNTrackerHits);
1032 }
1033 }
1034
1035 }
1036
1037
1038 for (size_t i = 0; i < nSim; i++) {
1039 TrackingParticleRef simRef = TPrefV[i];
1040 const TrackingParticle* simTP = simRef.get();
1041 if (!tpSelector_(*simTP))
1042 continue;
1043
1044
1045 const double simP = simRef->p();
1046 const double simPt = simRef->pt();
1047 const double simEta = doAbsEta_ ? fabs(simRef->eta()) : simRef->eta();
1048 const double simPhi = simRef->phi();
1049
1050 GlobalPoint simVtx(simRef->vertex().x(), simRef->vertex().y(), simRef->vertex().z());
1051 GlobalVector simMom(simRef->momentum().x(), simRef->momentum().y(), simRef->momentum().z());
1052 const double simDxy = -simVtx.x() * sin(simPhi) + simVtx.y() * cos(simPhi);
1053 const double simDz = simVtx.z() - (simVtx.x() * simMom.x() + simVtx.y() * simMom.y()) * simMom.z() / simMom.perp2();
1054
1055 const unsigned int nSimHits = simRef->numberOfHits();
1056
1057 muonME_->hSimP_->Fill(simP);
1058 muonME_->hSimPt_->Fill(simPt);
1059 muonME_->hSimEta_->Fill(simEta);
1060 muonME_->hSimPhi_->Fill(simPhi);
1061 muonME_->hSimDxy_->Fill(simDxy);
1062 muonME_->hSimDz_->Fill(simDz);
1063 muonME_->hNSimHits_->Fill(nSimHits);
1064
1065
1066 vector<pair<RefToBase<Muon>, double> > MuRefV;
1067 if (simToMuonColl.find(simRef) != simToMuonColl.end()) {
1068 MuRefV = simToMuonColl[simRef];
1069
1070 if (!MuRefV.empty()) {
1071 muonME_->hNSimToReco_->Fill(MuRefV.size());
1072 const Muon* Mu = MuRefV.begin()->first.get();
1073 if (!selector_(*Mu))
1074 continue;
1075 if (wantTightMuon_) {
1076 if (!muon::isTightMuon(*Mu, thePrimaryVertex))
1077 continue;
1078 }
1079
1080 muonME_->fill(&*simTP, Mu);
1081 }
1082 }
1083 }
1084 }
1085
1086 int RecoDisplacedMuonValidator::countMuonHits(const reco::Track& track) const {
1087 TransientTrackingRecHit::ConstRecHitContainer result;
1088
1089 int count = 0;
1090
1091 for (trackingRecHit_iterator hit = track.recHitsBegin(); hit != track.recHitsEnd(); ++hit) {
1092 if ((*hit)->isValid()) {
1093 DetId recoid = (*hit)->geographicalId();
1094 if (recoid.det() == DetId::Muon)
1095 count++;
1096 }
1097 }
1098 return count;
1099 }
1100
1101 int RecoDisplacedMuonValidator::countTrackerHits(const reco::Track& track) const {
1102 TransientTrackingRecHit::ConstRecHitContainer result;
1103
1104 int count = 0;
1105
1106 for (trackingRecHit_iterator hit = track.recHitsBegin(); hit != track.recHitsEnd(); ++hit) {
1107 if ((*hit)->isValid()) {
1108 DetId recoid = (*hit)->geographicalId();
1109 if (recoid.det() == DetId::Tracker)
1110 count++;
1111 }
1112 }
1113 return count;
1114 }