File indexing completed on 2024-09-07 04:34:37
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021 #include <algorithm> // std::sort
0022 #include <memory>
0023 #include <string>
0024 #include <vector> // std::vector
0025
0026
0027 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0028 #include "DataFormats/TrackReco/interface/Track.h"
0029 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0030 #include "DataFormats/VertexReco/interface/Vertex.h"
0031 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0032 #include "FWCore/Framework/interface/Event.h"
0033 #include "FWCore/Framework/interface/EventSetup.h"
0034 #include "FWCore/Framework/interface/Frameworkfwd.h"
0035 #include "FWCore/Framework/interface/MakerMacros.h"
0036 #include "FWCore/Framework/interface/MakerMacros.h"
0037 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0038 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0039 #include "FWCore/ServiceRegistry/interface/Service.h"
0040 #include "FWCore/Utilities/interface/InputTag.h"
0041 #include "Geometry/Records/interface/GlobalTrackingGeometryRecord.h"
0042 #include "RecoVertex/AdaptiveVertexFit/interface/AdaptiveVertexFitter.h"
0043 #include "RecoVertex/VertexTools/interface/VertexDistance3D.h"
0044 #include "TrackingTools/Records/interface/TransientTrackRecord.h"
0045 #include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
0046 #include "Alignment/OfflineValidation/interface/SmartSelectionMonitor.h"
0047 #include "DataFormats/Common/interface/TriggerResults.h" // Classes needed to print trigger results
0048 #include "FWCore/Common/interface/TriggerNames.h"
0049
0050
0051
0052
0053 using reco::TrackCollection;
0054
0055 class JetHTAnalyzer : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0056 public:
0057 explicit JetHTAnalyzer(const edm::ParameterSet&);
0058 ~JetHTAnalyzer() override = default;
0059
0060 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0061 static bool mysorter(reco::Track i, reco::Track j) { return (i.pt() > j.pt()); }
0062
0063 private:
0064 void beginJob() override;
0065 void analyze(const edm::Event&, const edm::EventSetup&) override;
0066 void endJob() override;
0067
0068
0069 const edm::InputTag pvsTag_;
0070 const edm::EDGetTokenT<reco::VertexCollection> pvsToken_;
0071
0072 const edm::InputTag tracksTag_;
0073 const edm::EDGetTokenT<reco::TrackCollection> tracksToken_;
0074
0075 const edm::InputTag triggerTag_;
0076 const edm::EDGetTokenT<edm::TriggerResults> triggerToken_;
0077
0078 const int printTriggerTable_;
0079 const double minVtxNdf_;
0080 const double minVtxWgt_;
0081
0082 const std::vector<double> profilePtBorders_;
0083 const std::vector<int> iovList_;
0084
0085
0086 edm::Service<TFileService> outfile_;
0087 TH1F* h_ntrks;
0088 TH1F* h_probePt;
0089 TH1F* h_probeEta;
0090 TH1F* h_probePhi;
0091 TH1F* h_probeDxy;
0092 TH1F* h_probeDz;
0093 TH1F* h_probeDxyErr;
0094 TH1F* h_probeDzErr;
0095
0096 SmartSelectionMonitor mon;
0097
0098
0099 static constexpr double cmToum = 10000;
0100 };
0101
0102
0103
0104
0105 JetHTAnalyzer::JetHTAnalyzer(const edm::ParameterSet& iConfig)
0106 : pvsTag_(iConfig.getParameter<edm::InputTag>("vtxCollection")),
0107 pvsToken_(consumes<reco::VertexCollection>(pvsTag_)),
0108 tracksTag_(iConfig.getParameter<edm::InputTag>("trackCollection")),
0109 tracksToken_(consumes<reco::TrackCollection>(tracksTag_)),
0110 triggerTag_(iConfig.getParameter<edm::InputTag>("triggerResults")),
0111 triggerToken_(consumes<edm::TriggerResults>(triggerTag_)),
0112 printTriggerTable_(iConfig.getUntrackedParameter<int>("printTriggerTable")),
0113 minVtxNdf_(iConfig.getUntrackedParameter<double>("minVertexNdf")),
0114 minVtxWgt_(iConfig.getUntrackedParameter<double>("minVertexMeanWeight")),
0115 profilePtBorders_(iConfig.getUntrackedParameter<std::vector<double>>("profilePtBorders")),
0116 iovList_(iConfig.getUntrackedParameter<std::vector<int>>("iovList")) {
0117
0118 usesResource(TFileService::kSharedResource);
0119 }
0120
0121
0122
0123
0124
0125
0126 void JetHTAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0127 using namespace edm;
0128
0129 const auto& vertices = iEvent.get(pvsToken_);
0130 const reco::VertexCollection& pvtx = vertices;
0131
0132 edm::Handle<reco::TrackCollection> tracks = iEvent.getHandle(tracksToken_);
0133
0134
0135 const int runNumber = iEvent.id().run();
0136 std::string iovString = "iovNotFound";
0137
0138
0139 if (runNumber >= iovList_.at(0)) {
0140 for (std::vector<int>::size_type i = 1; i < iovList_.size(); i++) {
0141 if (iovList_.at(i) > runNumber) {
0142 iovString = Form("iov%d-%d", iovList_.at(i - 1), iovList_.at(i) - 1);
0143 break;
0144 }
0145 }
0146 }
0147
0148
0149 if (printTriggerTable_) {
0150 const auto& triggerResults = iEvent.get(triggerToken_);
0151
0152 const edm::TriggerNames& triggerNames = iEvent.triggerNames(triggerResults);
0153 for (unsigned i = 0; i < triggerNames.size(); i++) {
0154 const std::string& hltName = triggerNames.triggerName(i);
0155 bool decision = triggerResults.accept(triggerNames.triggerIndex(hltName));
0156 std::cout << hltName << " " << decision << std::endl;
0157 }
0158 }
0159
0160 int counter = 0;
0161 for (reco::VertexCollection::const_iterator pvIt = pvtx.begin(); pvIt != pvtx.end(); pvIt++) {
0162 const reco::Vertex& iPV = *pvIt;
0163 counter++;
0164
0165 if (iPV.isFake())
0166 continue;
0167 reco::Vertex::trackRef_iterator trki;
0168
0169 const math::XYZPoint pos_(iPV.x(), iPV.y(), iPV.z());
0170
0171
0172 if (iPV.ndof() < minVtxNdf_ || (iPV.ndof() + 3.) / iPV.tracksSize() < 2 * minVtxWgt_)
0173 continue;
0174
0175 reco::TrackCollection allTracks;
0176 for (trki = iPV.tracks_begin(); trki != iPV.tracks_end(); ++trki) {
0177 if (trki->isNonnull()) {
0178 reco::TrackRef trk_now(tracks, (*trki).key());
0179 allTracks.push_back(*trk_now);
0180 }
0181 }
0182
0183
0184 std::sort(allTracks.begin(), allTracks.end(), mysorter);
0185 uint ntrks = allTracks.size();
0186 h_ntrks->Fill(ntrks);
0187
0188 for (uint tracksIt = 0; tracksIt < ntrks; tracksIt++) {
0189 auto tk = allTracks.at(tracksIt);
0190
0191 double dxyRes = tk.dxy(pos_) * cmToum;
0192 double dzRes = tk.dz(pos_) * cmToum;
0193
0194 double dxy_err = tk.dxyError() * cmToum;
0195 double dz_err = tk.dzError() * cmToum;
0196
0197 float trackphi = tk.phi();
0198 float tracketa = tk.eta();
0199 float trackpt = tk.pt();
0200
0201 h_probePt->Fill(trackpt);
0202 h_probeEta->Fill(tracketa);
0203 h_probePhi->Fill(trackphi);
0204
0205 h_probeDxy->Fill(dxyRes);
0206 h_probeDz->Fill(dzRes);
0207 h_probeDxyErr->Fill(dxy_err);
0208 h_probeDzErr->Fill(dz_err);
0209
0210 mon.fillHisto("dxy", "all", dxyRes, 1.);
0211 mon.fillHisto("dz", "all", dzRes, 1.);
0212 mon.fillHisto("dxyerr", "all", dxy_err, 1.);
0213 mon.fillHisto("dzerr", "all", dz_err, 1.);
0214
0215 mon.fillProfile("dxyErrVsPt", "all", trackpt, dxy_err, 1.);
0216 mon.fillProfile("dzErrVsPt", "all", trackpt, dz_err, 1.);
0217
0218 mon.fillProfile("dxyErrVsPhi", "all", trackphi, dxy_err, 1.);
0219 mon.fillProfile("dzErrVsPhi", "all", trackphi, dz_err, 1.);
0220
0221 mon.fillProfile("dxyErrVsEta", "all", tracketa, dxy_err, 1.);
0222 mon.fillProfile("dzErrVsEta", "all", tracketa, dz_err, 1.);
0223
0224
0225 for (std::vector<double>::size_type i = 0; i < profilePtBorders_.size(); i++) {
0226 if (trackpt < profilePtBorders_.at(i))
0227 break;
0228 mon.fillProfile("dxyErrVsPtWide", "all", i, dxy_err, 1.);
0229 mon.fillProfile("dzErrVsPtWide", "all", i, dz_err, 1.);
0230 }
0231
0232
0233 mon.fillHisto("dxy", iovString, dxyRes, 1.);
0234 mon.fillHisto("dz", iovString, dzRes, 1.);
0235 mon.fillHisto("dxyerr", iovString, dxy_err, 1.);
0236 mon.fillHisto("dzerr", iovString, dz_err, 1.);
0237
0238 mon.fillProfile("dxyErrVsPt", iovString, trackpt, dxy_err, 1.);
0239 mon.fillProfile("dzErrVsPt", iovString, trackpt, dz_err, 1.);
0240
0241 mon.fillProfile("dxyErrVsPhi", iovString, trackphi, dxy_err, 1.);
0242 mon.fillProfile("dzErrVsPhi", iovString, trackphi, dz_err, 1.);
0243
0244 mon.fillProfile("dxyErrVsEta", iovString, tracketa, dxy_err, 1.);
0245 mon.fillProfile("dzErrVsEta", iovString, tracketa, dz_err, 1.);
0246
0247
0248 for (std::vector<double>::size_type i = 0; i < profilePtBorders_.size(); i++) {
0249 if (trackpt < profilePtBorders_.at(i))
0250 break;
0251 mon.fillProfile("dxyErrVsPtWide", iovString, i, dxy_err, 1.);
0252 mon.fillProfile("dzErrVsPtWide", iovString, i, dz_err, 1.);
0253 }
0254
0255 if (std::abs(tracketa) < 1.) {
0256 mon.fillHisto("dxy", "central", dxyRes, 1.);
0257 mon.fillHisto("dz", "central", dzRes, 1.);
0258 mon.fillHisto("dxyerr", "central", dxy_err, 1.);
0259 mon.fillHisto("dzerr", "central", dz_err, 1.);
0260
0261 mon.fillProfile("dxyErrVsPt", "central", trackpt, dxy_err, 1.);
0262 mon.fillProfile("dzErrVsPt", "central", trackpt, dz_err, 1.);
0263
0264 mon.fillProfile("dxyErrVsPhi", "central", trackphi, dxy_err, 1.);
0265 mon.fillProfile("dzErrVsPhi", "central", trackphi, dz_err, 1.);
0266
0267
0268 for (std::vector<double>::size_type i = 0; i < profilePtBorders_.size(); i++) {
0269 if (trackpt < profilePtBorders_.at(i))
0270 break;
0271 mon.fillProfile("dxyErrVsPtWide", "central", i, dxy_err, 1.);
0272 mon.fillProfile("dzErrVsPtWide", "central", i, dz_err, 1.);
0273 }
0274 }
0275
0276 }
0277 }
0278
0279 mon.fillHisto("nvtx", "all", counter, 1.);
0280 }
0281
0282
0283 void JetHTAnalyzer::beginJob() {
0284 h_ntrks = outfile_->make<TH1F>("h_ntrks", "n. trks;n. of tracks/vertex;n. vertices", 100, 0, 100);
0285 h_probePt = outfile_->make<TH1F>("h_probePt", "p_{T} of probe track;track p_{T} (GeV); tracks", 100, 0., 500.);
0286 h_probeEta = outfile_->make<TH1F>("h_probeEta", "#eta of the probe track;track #eta;tracks", 54, -2.8, 2.8);
0287 h_probePhi = outfile_->make<TH1F>("h_probePhi", "#phi of probe track;track #phi (rad);tracks", 100, -3.15, 3.15);
0288
0289 h_probeDxy =
0290 outfile_->make<TH1F>("h_probeDxy", "d_{xy}(PV) of the probe track;track d_{xy}(PV);tracks", 200, -100, 100);
0291 h_probeDz = outfile_->make<TH1F>("h_probeDz", "d_{z}(PV) of the probe track;track d_{z}(PV);tracks", 200, -100, 100);
0292 h_probeDxyErr = outfile_->make<TH1F>(
0293 "h_probeDxyErr", "error on d_{xy}(PV) of the probe track;track error on d_{xy}(PV);tracks", 100, 0., 100);
0294 h_probeDzErr = outfile_->make<TH1F>(
0295 "h_probeDzErr", "error on d_{z}(PV) of the probe track;track error on d_{z}(PV);tracks", 100, 0., 100);
0296
0297 mon.addHistogram(new TH1F("nvtx", ";Vertices;Events", 50, 0, 50));
0298 mon.addHistogram(new TH1F("dxy", ";d_{xy};tracks", 100, -100, 100));
0299 mon.addHistogram(new TH1F("dz", ";d_{z};tracks", 100, -100, 100));
0300 mon.addHistogram(new TH1F("dxyerr", ";d_{xy} error;tracks", 100, 0., 200));
0301 mon.addHistogram(new TH1F("dzerr", ";d_{z} error;tracks", 100, 0., 200));
0302 mon.addHistogram(new TProfile("dxyErrVsPt", ";track p_{T};d_{xy} error", 100, 0., 200, 0., 100.));
0303 mon.addHistogram(new TProfile("dzErrVsPt", ";track p_{T};d_{z} error", 100, 0., 200, 0., 100.));
0304 mon.addHistogram(new TProfile("dxyErrVsPhi", ";track #varphi;d_{xy} error", 100, -M_PI, M_PI, 0., 100.));
0305 mon.addHistogram(new TProfile("dzErrVsPhi", ";track #varphi;d_{z} error", 100, -M_PI, M_PI, 0., 100.));
0306 mon.addHistogram(new TProfile("dxyErrVsEta", ";track #eta;d_{xy} error", 100, -2.5, 2.5, 0., 100.));
0307 mon.addHistogram(new TProfile("dzErrVsEta", ";track #eta;d_{z} error", 100, -2.5, 2.5, 0., 100.));
0308
0309
0310 int nBins = profilePtBorders_.size();
0311 mon.addHistogram(
0312 new TProfile("dxyErrVsPtWide", ";track p_{T} wide bin;d_{xy} error", nBins, -0.5, nBins - 0.5, 0.0, 100.0));
0313 mon.addHistogram(
0314 new TProfile("dzErrVsPtWide", ";track p_{T} wide bin;d_{z} error", nBins, -0.5, nBins - 0.5, 0.0, 100.0));
0315 }
0316
0317
0318 void JetHTAnalyzer::endJob() { mon.Write(); }
0319
0320
0321 void JetHTAnalyzer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0322 edm::ParameterSetDescription desc;
0323 desc.setComment("JetHT validation analyzer plugin.");
0324 desc.add<edm::InputTag>("vtxCollection", edm::InputTag("offlinePrimaryVerticesFromRefittedTrks"));
0325 desc.add<edm::InputTag>("triggerResults", edm::InputTag("TriggerResults", "", "HLT"));
0326 desc.add<edm::InputTag>("trackCollection", edm::InputTag("TrackRefitter"));
0327 desc.addUntracked<int>("printTriggerTable", false);
0328 desc.addUntracked<double>("minVertexNdf", 10.);
0329 desc.addUntracked<double>("minVertexMeanWeight", 0.5);
0330 desc.addUntracked<std::vector<double>>("profilePtBorders", {3, 5, 10, 20, 50, 100});
0331 desc.addUntracked<std::vector<int>>("iovList", {0, 500000});
0332 descriptions.addWithDefaultLabel(desc);
0333 }
0334
0335
0336 DEFINE_FWK_MODULE(JetHTAnalyzer);
