File indexing completed on 2025-03-31 22:27:14
0001 #ifndef MuonTrackValidatorBase_h
0002 #define MuonTrackValidatorBase_h
0003
0004
0005
0006
0007
0008
0009 #include <memory>
0010
0011 #include "FWCore/Framework/interface/Event.h"
0012 #include "FWCore/Framework/interface/ESHandle.h"
0013 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0014 #include "FWCore/Framework/interface/ConsumesCollector.h"
0015
0016 #include "MagneticField/Engine/interface/MagneticField.h"
0017 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0018
0019 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0020 #include "DataFormats/TrackReco/interface/Track.h"
0021 #include "SimDataFormats/TrackingAnalysis/interface/TrackingParticle.h"
0022 #include "SimDataFormats/TrackingAnalysis/interface/TrackingParticleFwd.h"
0023 #include "SimDataFormats/PileupSummaryInfo/interface/PileupSummaryInfo.h"
0024
0025 #include "DQMServices/Core/interface/DQMStore.h"
0026 #include "FWCore/ServiceRegistry/interface/Service.h"
0027
0028 #include "SimTracker/Common/interface/TrackingParticleSelector.h"
0029 #include "CommonTools/RecoAlgos/interface/CosmicTrackingParticleSelector.h"
0030
0031 #include <iostream>
0032 #include <sstream>
0033 #include <string>
0034 #include <TH1F.h>
0035 #include <TH2F.h>
0036
0037 struct HistoParams {
0038 HistoParams(edm::ParameterSet muonHistoParameters) {
0039 mindR = muonHistoParameters.getParameter<double>("mindR");
0040 maxdR = muonHistoParameters.getParameter<double>("maxdR");
0041 nintdR = muonHistoParameters.getParameter<int>("nintdR");
0042 minEta = muonHistoParameters.getParameter<double>("minEta");
0043 maxEta = muonHistoParameters.getParameter<double>("maxEta");
0044 nintEta = muonHistoParameters.getParameter<int>("nintEta");
0045 minPt = muonHistoParameters.getParameter<double>("minPt");
0046 maxPt = muonHistoParameters.getParameter<double>("maxPt");
0047 nintPt = muonHistoParameters.getParameter<int>("nintPt");
0048 useLogPt = muonHistoParameters.getUntrackedParameter<bool>("useLogPt", false);
0049 minNHit = muonHistoParameters.getParameter<double>("minNHit");
0050 maxNHit = muonHistoParameters.getParameter<double>("maxNHit");
0051 nintNHit = muonHistoParameters.getParameter<int>("nintNHit");
0052
0053 minDTHit = muonHistoParameters.getParameter<double>("minDTHit");
0054 maxDTHit = muonHistoParameters.getParameter<double>("maxDTHit");
0055 nintDTHit = muonHistoParameters.getParameter<int>("nintDTHit");
0056
0057 minCSCHit = muonHistoParameters.getParameter<double>("minCSCHit");
0058 maxCSCHit = muonHistoParameters.getParameter<double>("maxCSCHit");
0059 nintCSCHit = muonHistoParameters.getParameter<int>("nintCSCHit");
0060
0061 minRPCHit = muonHistoParameters.getParameter<double>("minRPCHit");
0062 maxRPCHit = muonHistoParameters.getParameter<double>("maxRPCHit");
0063 nintRPCHit = muonHistoParameters.getParameter<int>("nintRPCHit");
0064
0065
0066 minNTracks = muonHistoParameters.getParameter<int>("minNTracks");
0067 maxNTracks = muonHistoParameters.getParameter<int>("maxNTracks");
0068 nintNTracks = muonHistoParameters.getParameter<int>("nintNTracks");
0069 minFTracks = muonHistoParameters.getParameter<int>("minFTracks");
0070 maxFTracks = muonHistoParameters.getParameter<int>("maxFTracks");
0071 nintFTracks = muonHistoParameters.getParameter<int>("nintFTracks");
0072 minLayers = muonHistoParameters.getParameter<double>("minLayers");
0073 maxLayers = muonHistoParameters.getParameter<double>("maxLayers");
0074 nintLayers = muonHistoParameters.getParameter<int>("nintLayers");
0075 minPixels = muonHistoParameters.getParameter<double>("minPixels");
0076 maxPixels = muonHistoParameters.getParameter<double>("maxPixels");
0077 nintPixels = muonHistoParameters.getParameter<int>("nintPixels");
0078 minPhi = muonHistoParameters.getParameter<double>("minPhi");
0079 maxPhi = muonHistoParameters.getParameter<double>("maxPhi");
0080 nintPhi = muonHistoParameters.getParameter<int>("nintPhi");
0081 minDxy = muonHistoParameters.getParameter<double>("minDxy");
0082 maxDxy = muonHistoParameters.getParameter<double>("maxDxy");
0083 nintDxy = muonHistoParameters.getParameter<int>("nintDxy");
0084 minDz = muonHistoParameters.getParameter<double>("minDz");
0085 maxDz = muonHistoParameters.getParameter<double>("maxDz");
0086 nintDz = muonHistoParameters.getParameter<int>("nintDz");
0087 minRpos = muonHistoParameters.getParameter<double>("minRpos");
0088 maxRpos = muonHistoParameters.getParameter<double>("maxRpos");
0089 nintRpos = muonHistoParameters.getParameter<int>("nintRpos");
0090 minZpos = muonHistoParameters.getParameter<double>("minZpos");
0091 maxZpos = muonHistoParameters.getParameter<double>("maxZpos");
0092 nintZpos = muonHistoParameters.getParameter<int>("nintZpos");
0093 minPU = muonHistoParameters.getParameter<double>("minPU");
0094 maxPU = muonHistoParameters.getParameter<double>("maxPU");
0095 nintPU = muonHistoParameters.getParameter<int>("nintPU");
0096
0097 ptRes_rangeMin = muonHistoParameters.getParameter<double>("ptRes_rangeMin");
0098 ptRes_rangeMax = muonHistoParameters.getParameter<double>("ptRes_rangeMax");
0099 ptRes_nbin = muonHistoParameters.getParameter<int>("ptRes_nbin");
0100 etaRes_rangeMin = muonHistoParameters.getParameter<double>("etaRes_rangeMin");
0101 etaRes_rangeMax = muonHistoParameters.getParameter<double>("etaRes_rangeMax");
0102 etaRes_nbin = muonHistoParameters.getParameter<int>("etaRes_nbin");
0103 phiRes_rangeMin = muonHistoParameters.getParameter<double>("phiRes_rangeMin");
0104 phiRes_rangeMax = muonHistoParameters.getParameter<double>("phiRes_rangeMax");
0105 phiRes_nbin = muonHistoParameters.getParameter<int>("phiRes_nbin");
0106 cotThetaRes_rangeMin = muonHistoParameters.getParameter<double>("cotThetaRes_rangeMin");
0107 cotThetaRes_rangeMax = muonHistoParameters.getParameter<double>("cotThetaRes_rangeMax");
0108 cotThetaRes_nbin = muonHistoParameters.getParameter<int>("cotThetaRes_nbin");
0109 dxyRes_rangeMin = muonHistoParameters.getParameter<double>("dxyRes_rangeMin");
0110 dxyRes_rangeMax = muonHistoParameters.getParameter<double>("dxyRes_rangeMax");
0111 dxyRes_nbin = muonHistoParameters.getParameter<int>("dxyRes_nbin");
0112 dzRes_rangeMin = muonHistoParameters.getParameter<double>("dzRes_rangeMin");
0113 dzRes_rangeMax = muonHistoParameters.getParameter<double>("dzRes_rangeMax");
0114 dzRes_nbin = muonHistoParameters.getParameter<int>("dzRes_nbin");
0115
0116 usetracker = muonHistoParameters.getParameter<bool>("usetracker");
0117 usemuon = muonHistoParameters.getParameter<bool>("usemuon");
0118 do_TRKhitsPlots = muonHistoParameters.getParameter<bool>("do_TRKhitsPlots");
0119 do_MUOhitsPlots = muonHistoParameters.getParameter<bool>("do_MUOhitsPlots");
0120
0121 if (useLogPt) {
0122 minPt = log10(std::max(0.01, minPt));
0123 maxPt = log10(maxPt);
0124 }
0125 }
0126
0127 int minNTracks, maxNTracks, nintNTracks;
0128 int minFTracks, maxFTracks, nintFTracks;
0129 double mindR, maxdR;
0130 int nintdR;
0131 double minEta, maxEta;
0132 int nintEta;
0133 double minPt, maxPt;
0134 int nintPt;
0135 bool useLogPt;
0136 double minNHit, maxNHit;
0137 int nintNHit;
0138 double minDTHit, maxDTHit;
0139 int nintDTHit;
0140 double minCSCHit, maxCSCHit;
0141 int nintCSCHit;
0142 double minRPCHit, maxRPCHit;
0143 int nintRPCHit;
0144 double minLayers, maxLayers;
0145 int nintLayers;
0146 double minPixels, maxPixels;
0147 int nintPixels;
0148 double minPhi, maxPhi;
0149 int nintPhi;
0150 double minDxy, maxDxy;
0151 int nintDxy;
0152 double minDz, maxDz;
0153 int nintDz;
0154 double minRpos, maxRpos;
0155 int nintRpos;
0156 double minZpos, maxZpos;
0157 int nintZpos;
0158 double minPU, maxPU;
0159 int nintPU;
0160
0161 double ptRes_rangeMin, ptRes_rangeMax;
0162 int ptRes_nbin;
0163 double etaRes_rangeMin, etaRes_rangeMax;
0164 int etaRes_nbin;
0165 double phiRes_rangeMin, phiRes_rangeMax;
0166 int phiRes_nbin;
0167 double cotThetaRes_rangeMin, cotThetaRes_rangeMax;
0168 int cotThetaRes_nbin;
0169 double dxyRes_rangeMin, dxyRes_rangeMax;
0170 int dxyRes_nbin;
0171 double dzRes_rangeMin, dzRes_rangeMax;
0172 int dzRes_nbin;
0173 bool usetracker, usemuon;
0174 bool do_TRKhitsPlots, do_MUOhitsPlots;
0175 };
0176
0177 class MuonTrackValidatorBase {
0178 typedef dqm::legacy::DQMStore DQMStore;
0179 typedef dqm::reco::MonitorElement MonitorElement;
0180
0181 public:
0182
0183 MuonTrackValidatorBase(const edm::ParameterSet& pset, edm::ConsumesCollector iC) : MuonTrackValidatorBase(pset) {
0184 bsSrc_Token = iC.consumes<reco::BeamSpot>(bsSrc);
0185
0186 if (label_tp_refvector)
0187 tp_refvector_Token = iC.consumes<TrackingParticleRefVector>(label_tp);
0188 else
0189 tp_Token = iC.consumes<TrackingParticleCollection>(label_tp);
0190
0191 pileupinfo_Token = iC.consumes<std::vector<PileupSummaryInfo>>(label_pileupinfo);
0192
0193 histoParameters.reserve(label.size());
0194 for (unsigned int www = 0; www < label.size(); ++www) {
0195 track_Collection_Token[www] = iC.consumes<edm::View<reco::Track>>(label[www]);
0196
0197 histoParameters.push_back(HistoParams(muonHistoParameters[www]));
0198 }
0199 }
0200
0201 MuonTrackValidatorBase(const edm::ParameterSet& pset)
0202 : label(pset.getParameter<std::vector<edm::InputTag>>("label")),
0203 bsSrc(pset.getParameter<edm::InputTag>("beamSpot")),
0204 label_tp(pset.getParameter<edm::InputTag>("label_tp")),
0205 label_tp_refvector(pset.getParameter<bool>("label_tp_refvector")),
0206 label_pileupinfo(pset.getParameter<edm::InputTag>("label_pileupinfo")),
0207 associators(pset.getParameter<std::vector<std::string>>("associators")),
0208 out(pset.getParameter<std::string>("outputFile")),
0209 parametersDefiner(pset.getParameter<std::string>("parametersDefiner")),
0210 ignoremissingtkcollection_(pset.getUntrackedParameter<bool>("ignoremissingtrackcollection", false)),
0211 muonHistoParameters(pset.getParameter<std::vector<edm::ParameterSet>>("muonHistoParameters")) {
0212 histoParameters.reserve(label.size());
0213 for (size_t www = 0; www < label.size(); ++www) {
0214 histoParameters.emplace_back(muonHistoParameters[www]);
0215 }
0216 }
0217
0218 virtual ~MuonTrackValidatorBase() noexcept(false) {}
0219
0220 template <typename T>
0221 void fillPlotNoFlow(MonitorElement* h, T val) {
0222 h->Fill(std::min(std::max(val, ((T)h->getTH1()->GetXaxis()->GetXmin())), ((T)h->getTH1()->GetXaxis()->GetXmax())));
0223 }
0224
0225 void doProfileX(TH2* th2, MonitorElement* me) {
0226 if (th2->GetNbinsX() == me->getNbinsX()) {
0227 TProfile* p1 = (TProfile*)th2->ProfileX();
0228 p1->Copy(*me->getTProfile());
0229 delete p1;
0230 } else {
0231 throw cms::Exception("MuonTrackValidator") << "Different number of bins!";
0232 }
0233 }
0234
0235 void doProfileX(MonitorElement* th2m, MonitorElement* me) { doProfileX(th2m->getTH2F(), me); }
0236
0237
0238 double getEta(double eta) { return eta; }
0239
0240
0241 double getPt(double pt) { return pt; }
0242
0243 void BinLogX(TH1* h) {
0244 TAxis* axis = h->GetXaxis();
0245 int bins = axis->GetNbins();
0246
0247 float from = axis->GetXmin();
0248 float to = axis->GetXmax();
0249 float width = (to - from) / bins;
0250 float* new_bins = new float[bins + 1];
0251
0252 for (int i = 0; i <= bins; i++) {
0253 new_bins[i] = TMath::Power(10, from + i * width);
0254 }
0255 axis->Set(bins, new_bins);
0256 delete[] new_bins;
0257 }
0258
0259 protected:
0260 std::vector<edm::InputTag> label;
0261 edm::InputTag bsSrc;
0262 edm::InputTag label_tp;
0263 bool label_tp_refvector;
0264 edm::InputTag label_pileupinfo;
0265 std::vector<std::string> associators;
0266 std::string out;
0267 std::string parametersDefiner;
0268 std::vector<edm::EDGetTokenT<edm::View<reco::Track>>> track_Collection_Token;
0269 edm::EDGetTokenT<reco::BeamSpot> bsSrc_Token;
0270 edm::EDGetTokenT<TrackingParticleCollection> tp_Token;
0271 edm::EDGetTokenT<TrackingParticleRefVector> tp_refvector_Token;
0272 edm::EDGetTokenT<std::vector<PileupSummaryInfo>> pileupinfo_Token;
0273 edm::ESHandle<MagneticField> theMF;
0274
0275 bool ignoremissingtkcollection_;
0276
0277 std::vector<edm::ParameterSet> muonHistoParameters;
0278
0279 std::vector<HistoParams> histoParameters;
0280
0281
0282 std::vector<MonitorElement*> h_assoc_coll, h_simul_coll, h_reco_coll, h_assoc2_coll;
0283
0284 std::vector<MonitorElement*> h_tracks, h_fakes, h_nhits, h_charge;
0285 std::vector<MonitorElement*> h_recoeta, h_assoceta, h_assoc2eta, h_simuleta, h_misideta;
0286 std::vector<MonitorElement*> h_recopT, h_assocpT, h_assoc2pT, h_simulpT, h_misidpT;
0287 std::vector<MonitorElement*> h_recohit, h_assochit, h_assoc2hit, h_simulhit, h_misidhit;
0288 std::vector<MonitorElement*> h_recophi, h_assocphi, h_assoc2phi, h_simulphi, h_misidphi;
0289 std::vector<MonitorElement*> h_recodxy, h_assocdxy, h_assoc2dxy, h_simuldxy, h_misiddxy;
0290 std::vector<MonitorElement*> h_recodz, h_assocdz, h_assoc2dz, h_simuldz, h_misiddz;
0291 std::vector<MonitorElement*> h_recopu, h_assocpu, h_assoc2pu, h_simulpu, h_misidpu;
0292 std::vector<MonitorElement*> h_recodR, h_assocdR, h_assoc2dR, h_simuldR, h_misiddR;
0293
0294 std::vector<MonitorElement*> h_assocpTB, h_simulpTB;
0295 std::vector<MonitorElement*> h_assocphiB, h_simulphiB;
0296
0297 std::vector<MonitorElement*> h_assocpTO, h_simulpTO;
0298 std::vector<MonitorElement*> h_assocphiO, h_simulphiO;
0299
0300 std::vector<MonitorElement*> h_assocpTE, h_simulpTE;
0301 std::vector<MonitorElement*> h_assocphiE, h_simulphiE;
0302
0303 std::vector<MonitorElement*> h_assocRpos, h_simulRpos, h_assocZpos, h_simulZpos;
0304 std::vector<MonitorElement*> h_etaRes;
0305
0306 std::vector<MonitorElement*> h_nchi2, h_nchi2_prob, h_losthits;
0307 std::vector<MonitorElement*> h_nmisslayers_inner, h_nmisslayers_outer, h_nlosthits;
0308 std::vector<MonitorElement*> h_assochi2, h_assochi2_prob;
0309 std::vector<MonitorElement*> h_assocFraction, h_assocSharedHit;
0310
0311
0312 std::vector<MonitorElement*> nRecHits_vs_nSimHits;
0313 std::vector<MonitorElement*> h_PurityVsQuality;
0314 std::vector<MonitorElement*> chi2_vs_nhits, etares_vs_eta;
0315 std::vector<MonitorElement*> ptres_vs_phi, chi2_vs_phi, nhits_vs_phi, phires_vs_phi;
0316
0317 std::vector<MonitorElement*> nhits_vs_eta, nDThits_vs_eta, nCSChits_vs_eta, nRPChits_vs_eta, nGEMhits_vs_eta,
0318 nME0hits_vs_eta;
0319 std::vector<MonitorElement*> chi2_vs_eta, nlosthits_vs_eta;
0320 std::vector<MonitorElement*> nTRK_LayersWithMeas_vs_eta, nPixel_LayersWithMeas_vs_eta;
0321
0322 std::vector<MonitorElement*> dxyres_vs_eta, ptres_vs_eta, dzres_vs_eta, phires_vs_eta, thetaCotres_vs_eta;
0323 std::vector<MonitorElement*> dxyres_vs_pt, ptres_vs_pt, dzres_vs_pt, phires_vs_pt, thetaCotres_vs_pt;
0324
0325 std::vector<MonitorElement*> dxypull_vs_eta, ptpull_vs_eta, dzpull_vs_eta, phipull_vs_eta, thetapull_vs_eta;
0326 std::vector<MonitorElement*> ptpull_vs_phi, phipull_vs_phi, thetapull_vs_phi;
0327 std::vector<MonitorElement*> h_dxypulleta, h_ptpulleta, h_dzpulleta, h_phipulleta, h_thetapulleta;
0328 std::vector<MonitorElement*> h_ptpullphi, h_phipullphi, h_thetapullphi;
0329 std::vector<MonitorElement*> h_ptpull, h_qoverppull, h_thetapull, h_phipull, h_dxypull, h_dzpull;
0330 };
0331
0332 #endif