File indexing completed on 2025-03-23 15:57:38
0001 #include "DQM/Physics/src/SingleTopTChannelLeptonDQM.h"
0002 #include "DataFormats/BTauReco/interface/JetTag.h"
0003 #include "DataFormats/Common/interface/Handle.h"
0004 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
0005 #include "DataFormats/JetReco/interface/CaloJet.h"
0006 #include "DataFormats/JetReco/interface/PFJet.h"
0007
0008 #include "DataFormats/Math/interface/deltaR.h"
0009 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
0010 #include "DataFormats/VertexReco/interface/Vertex.h"
0011 #include "JetMETCorrections/Objects/interface/JetCorrectionsRecord.h"
0012 #include <iostream>
0013 #include <memory>
0014
0015 #include "FWCore/Framework/interface/ConsumesCollector.h"
0016 #include "FWCore/Framework/interface/EDConsumerBase.h"
0017 #include "FWCore/Utilities/interface/EDGetToken.h"
0018 #include "FWCore/Utilities/interface/InputTag.h"
0019
0020 using namespace std;
0021 namespace SingleTopTChannelLepton {
0022
0023
0024
0025 static const unsigned int MAXJETS = 4;
0026
0027
0028 static const double WMASS = 80.4;
0029
0030 MonitorEnsemble::MonitorEnsemble(const char* label,
0031 const edm::ParameterSet& cfg,
0032 const edm::VParameterSet& vcfg,
0033 edm::ConsumesCollector&& iC)
0034 : label_(label),
0035 elecSelect_(nullptr),
0036 pvSelect_(nullptr),
0037 muonIso_(nullptr),
0038 muonSelect_(nullptr),
0039 jetIDSelect_(nullptr),
0040 jetlooseSelection_(nullptr),
0041 jetSelection_(nullptr),
0042 includeBTag_(false),
0043 lowerEdge_(-1.),
0044 upperEdge_(-1.),
0045 logged_(0) {
0046
0047 edm::ParameterSet sources = cfg.getParameter<edm::ParameterSet>("sources");
0048 muons_ = iC.consumes<edm::View<reco::Muon>>(sources.getParameter<edm::InputTag>("muons"));
0049 elecs_ = iC.consumes<edm::View<reco::GsfElectron>>(sources.getParameter<edm::InputTag>("elecs"));
0050 jets_ = iC.consumes<edm::View<reco::Jet>>(sources.getParameter<edm::InputTag>("jets"));
0051 for (edm::InputTag const& tag : sources.getParameter<std::vector<edm::InputTag>>("mets"))
0052 mets_.push_back(iC.consumes<edm::View<reco::MET>>(tag));
0053 pvs_ = iC.consumes<edm::View<reco::Vertex>>(sources.getParameter<edm::InputTag>("pvs"));
0054
0055
0056 if (cfg.existsAs<edm::ParameterSet>("elecExtras")) {
0057
0058
0059
0060
0061 edm::ParameterSet elecExtras = cfg.getParameter<edm::ParameterSet>("elecExtras");
0062
0063
0064 if (elecExtras.existsAs<std::string>("select")) {
0065 elecSelect_ = std::make_unique<StringCutObjectSelector<reco::GsfElectron>>(
0066 elecExtras.getParameter<std::string>("select"));
0067 }
0068
0069 if (elecExtras.existsAs<std::string>("rho")) {
0070 rhoTag = elecExtras.getParameter<edm::InputTag>("rho");
0071 }
0072
0073
0074 if (elecExtras.existsAs<edm::ParameterSet>("electronId")) {
0075 edm::ParameterSet elecId = elecExtras.getParameter<edm::ParameterSet>("electronId");
0076 electronId_ = iC.consumes<edm::ValueMap<float>>(elecId.getParameter<edm::InputTag>("src"));
0077 eidCutValue_ = elecId.getParameter<double>("cutValue");
0078 }
0079 }
0080
0081 if (cfg.existsAs<edm::ParameterSet>("pvExtras")) {
0082 edm::ParameterSet pvExtras = cfg.getParameter<edm::ParameterSet>("pvExtras");
0083
0084
0085 if (pvExtras.existsAs<std::string>("select")) {
0086 pvSelect_ =
0087 std::make_unique<StringCutObjectSelector<reco::Vertex>>(pvExtras.getParameter<std::string>("select"));
0088 }
0089 }
0090
0091 if (cfg.existsAs<edm::ParameterSet>("muonExtras")) {
0092 edm::ParameterSet muonExtras = cfg.getParameter<edm::ParameterSet>("muonExtras");
0093
0094
0095 if (muonExtras.existsAs<std::string>("select")) {
0096 muonSelect_ =
0097 std::make_unique<StringCutObjectSelector<reco::Muon>>(muonExtras.getParameter<std::string>("select"));
0098 }
0099
0100
0101 if (muonExtras.existsAs<std::string>("isolation")) {
0102 muonIso_ =
0103 std::make_unique<StringCutObjectSelector<reco::Muon>>(muonExtras.getParameter<std::string>("isolation"));
0104 }
0105 }
0106
0107
0108
0109 if (cfg.existsAs<edm::ParameterSet>("jetExtras")) {
0110 edm::ParameterSet jetExtras = cfg.getParameter<edm::ParameterSet>("jetExtras");
0111
0112
0113 if (jetExtras.existsAs<std::string>("jetCorrector")) {
0114 jetCorrector_ =
0115 iC.consumes<reco::JetCorrector>(edm::InputTag(jetExtras.getParameter<std::string>("jetCorrector")));
0116 }
0117
0118 if (jetExtras.existsAs<edm::ParameterSet>("jetID")) {
0119 edm::ParameterSet jetID = jetExtras.getParameter<edm::ParameterSet>("jetID");
0120 jetIDLabel_ = iC.consumes<reco::JetIDValueMap>(jetID.getParameter<edm::InputTag>("label"));
0121 jetIDSelect_ =
0122 std::make_unique<StringCutObjectSelector<reco::JetID>>(jetID.getParameter<std::string>("select"));
0123 }
0124
0125
0126
0127 if (jetExtras.existsAs<std::string>("select")) {
0128 jetSelect_ = jetExtras.getParameter<std::string>("select");
0129 jetSelection_ = std::make_unique<StringCutObjectSelector<reco::PFJet>>(jetSelect_);
0130 jetlooseSelection_ = std::make_unique<StringCutObjectSelector<reco::PFJet>>(
0131 "chargedHadronEnergyFraction()>0 && chargedMultiplicity()>0 && chargedEmEnergyFraction()<0.99 && "
0132 "neutralHadronEnergyFraction()<0.99 && neutralEmEnergyFraction()<0.99 && "
0133 "(chargedMultiplicity()+neutralMultiplicity())>1");
0134 }
0135
0136
0137
0138
0139 includeBTag_ = jetExtras.existsAs<edm::ParameterSet>("jetBTaggers");
0140 if (includeBTag_) {
0141 edm::ParameterSet btagCSV =
0142 jetExtras.getParameter<edm::ParameterSet>("jetBTaggers").getParameter<edm::ParameterSet>("cvsVertex");
0143 btagCSV_ = iC.consumes<reco::JetTagCollection>(btagCSV.getParameter<edm::InputTag>("label"));
0144 btagCSVWP_ = btagCSV.getParameter<double>("workingPoint");
0145 }
0146 }
0147
0148
0149 if (cfg.existsAs<edm::ParameterSet>("triggerExtras")) {
0150 edm::ParameterSet triggerExtras = cfg.getParameter<edm::ParameterSet>("triggerExtras");
0151 triggerTable_ = iC.consumes<edm::TriggerResults>(triggerExtras.getParameter<edm::InputTag>("src"));
0152 triggerPaths_ = triggerExtras.getParameter<std::vector<std::string>>("paths");
0153 }
0154
0155
0156
0157
0158 if (cfg.existsAs<edm::ParameterSet>("massExtras")) {
0159 edm::ParameterSet massExtras = cfg.getParameter<edm::ParameterSet>("massExtras");
0160 lowerEdge_ = massExtras.getParameter<double>("lowerEdge");
0161 upperEdge_ = massExtras.getParameter<double>("upperEdge");
0162 }
0163
0164
0165
0166
0167
0168 verbosity_ = STANDARD;
0169 if (cfg.existsAs<edm::ParameterSet>("monitoring")) {
0170 edm::ParameterSet monitoring = cfg.getParameter<edm::ParameterSet>("monitoring");
0171 if (monitoring.getParameter<std::string>("verbosity") == "DEBUG")
0172 verbosity_ = DEBUG;
0173 if (monitoring.getParameter<std::string>("verbosity") == "VERBOSE")
0174 verbosity_ = VERBOSE;
0175 if (monitoring.getParameter<std::string>("verbosity") == "STANDARD")
0176 verbosity_ = STANDARD;
0177 }
0178
0179 directory_ = cfg.getParameter<std::string>("directory");
0180 }
0181
0182 void MonitorEnsemble::book(DQMStore::IBooker& ibooker) {
0183
0184 std::string current(directory_);
0185 current += label_;
0186 ibooker.setCurrentFolder(current);
0187
0188
0189
0190 hists_["pvMult_"] = ibooker.book1D("PvMult", "N_{good pvs}", 50, 0., 100.);
0191
0192 hists_["muonPt_"] = ibooker.book1D("MuonPt", "pt(#mu TightId, TightIso)", 40, 0., 200.);
0193
0194 hists_["muonMult_"] = ibooker.book1D("MuonMult", "N_{loose}(#mu)", 10, 0., 10.);
0195
0196 hists_["muonMultTight_"] = ibooker.book1D("MuonMultTight", "N_{TightIso,TightId}(#mu)", 10, 0., 10.);
0197
0198 hists_["elecPt_"] = ibooker.book1D("ElecPt", "pt(e TightId, TightIso)", 40, 0., 200.);
0199
0200 hists_["jetMult_"] = ibooker.book1D("JetMult", "N_{30}(jet)", 10, 0., 10.);
0201
0202 hists_["jetMultLoose_"] = ibooker.book1D("JetMultLoose", "N_{30,loose}(jet)", 10, 0., 10.);
0203
0204
0205 hists_["metPflow_"] = ibooker.book1D("METPflow", "MET_{Pflow}", 50, 0., 200.);
0206
0207 hists_["massW_"] = ibooker.book1D("MassW", "M(W)", 60, 0., 300.);
0208
0209 hists_["massTop_"] = ibooker.book1D("MassTop", "M(Top)", 50, 0., 500.);
0210
0211 hists_["MTWm_"] = ibooker.book1D("MTWm", "M_{T}^{W}(#mu)", 60, 0., 300.);
0212
0213 hists_["mMTT_"] = ibooker.book1D("mMTT", "M_{T}^{t}(#mu)", 50, 0., 500.);
0214
0215
0216 hists_["MTWe_"] = ibooker.book1D("MTWe", "M_{T}^{W}(e)", 60, 0., 300.);
0217
0218 hists_["eMTT_"] = ibooker.book1D("eMTT", "M_{T}^{t}(e)", 50, 0., 500.);
0219
0220
0221 triggerBinLabels(std::string("trigger"), triggerPaths_);
0222
0223 if (verbosity_ == STANDARD)
0224 return;
0225
0226
0227
0228
0229 hists_["muonEta_"] = ibooker.book1D("MuonEta", "#eta(#mu TightId, TightIso)", 30, -3., 3.);
0230
0231 hists_["muonRelIso_"] = ibooker.book1D("MuonRelIso", "Iso_{Rel}(#mu TightId) (#Delta#beta Corrected)", 50, 0., 1.);
0232
0233 hists_["muonPhi_"] = ibooker.book1D("MuonPhi", "#phi(#mu TightId, TightIso)", 40, -4., 4.);
0234
0235 hists_["elecEta_"] = ibooker.book1D("ElecEta", "#eta(e tightId, TightIso)", 30, -3., 3.);
0236
0237 hists_["elecRelIso_"] = ibooker.book1D("ElecRelIso", "Iso_{Rel}(e TightId)", 50, 0., 1.);
0238
0239 hists_["elecPhi_"] = ibooker.book1D("ElecPhi", "#phi(e tightId, TightIso)", 40, -4., 4.);
0240
0241 hists_["elecMultTight_"] = ibooker.book1D("ElecMultTight", "N_{TightIso,TightId}(e)", 10, 0., 10.);
0242
0243
0244 hists_["jet1Eta_"] = ibooker.book1D("Jet1Eta", "#eta_{30,loose}(jet1)", 60, -3., 3.);
0245
0246 hists_["jet1Pt_"] = ibooker.book1D("Jet1Pt", "pt_{30,loose}(jet1)", 60, 0., 300.);
0247
0248 hists_["jet2Eta_"] = ibooker.book1D("Jet2Eta", "#eta_{30,loose}(jet2)", 60, -3., 3.);
0249
0250 hists_["jet2Pt_"] = ibooker.book1D("Jet2Pt", "pt_{30,loose}(jet2)", 60, 0., 300.);
0251
0252 hists_["jet3Eta_"] = ibooker.book1D("Jet3Eta", "#eta_{30,loose}(jet3)", 60, -3., 3.);
0253
0254 hists_["jet3Pt_"] = ibooker.book1D("Jet3Pt", "pt_{30,loose}(jet3)", 60, 0., 300.);
0255
0256 hists_["jet4Eta_"] = ibooker.book1D("Jet4Eta", "#eta_{30,loose}(jet4)", 60, -3., 3.);
0257
0258 hists_["jet4Pt_"] = ibooker.book1D("Jet4Pt", "pt_{30,loose}(jet4)", 60, 0., 300.);
0259
0260 hists_["muonDelZ_"] = ibooker.book1D("MuonDelZ", "d_{z}(#mu)", 50, -25., 25.);
0261
0262 hists_["muonDelXY_"] = ibooker.book2D("MuonDelXY", "d_{xy}(#mu)", 50, -0.1, 0.1, 50, -0.1, 0.1);
0263
0264 hists_["muonDxy_"] = ibooker.book1D("MuonDxy", "d_{xy}(#mu)", 100, -0.05, 0.05);
0265
0266
0267 hists_["muonDelXY_"]->setAxisTitle("x [cm]", 1);
0268 hists_["muonDelXY_"]->setAxisTitle("y [cm]", 2);
0269 hists_["muonDxy_"]->setAxisTitle("d_{xy} [cm]", 1);
0270
0271 if (verbosity_ == VERBOSE)
0272 return;
0273
0274
0275
0276
0277 hists_["muonChHadIso_"] = ibooker.book1D("MuonChHadIsoComp", "ChHad_{IsoComponent}(#mu TightId)", 50, 0., 5.);
0278
0279
0280 hists_["muonNeHadIso_"] = ibooker.book1D("MuonNeHadIsoComp", "NeHad_{IsoComponent}(#mu TightId)", 50, 0., 5.);
0281
0282
0283 hists_["muonPhIso_"] = ibooker.book1D("MuonPhIsoComp", "Photon_{IsoComponent}(#mu TightId)", 50, 0., 5.);
0284
0285
0286 hists_["elecChHadIso_"] = ibooker.book1D("ElectronChHadIsoComp", "ChHad_{IsoComponent}(e tightId)", 50, 0., 5.);
0287
0288
0289 hists_["elecNeHadIso_"] = ibooker.book1D("ElectronNeHadIsoComp", "NeHad_{IsoComponent}(e tightId)", 50, 0., 5.);
0290
0291
0292 hists_["elecPhIso_"] = ibooker.book1D("ElectronPhIsoComp", "Photon_{IsoComponent}(e tightId)", 50, 0., 5.);
0293
0294
0295 hists_["jetMultBCSVM_"] = ibooker.book1D("JetMultBCSVM", "N_{30}(CSVM)", 10, 0., 10.);
0296
0297 hists_["jetBCSV_"] = ibooker.book1D("JetDiscCSV", "BJet Disc_{CSV}(JET)", 100, -1., 2.);
0298
0299
0300
0301
0302
0303
0304
0305
0306
0307 hists_["eventLogger_"] = ibooker.book2D("EventLogger", "Logged Events", 9, 0., 9., 10, 0., 10.);
0308
0309
0310 hists_["eventLogger_"]->getTH1()->SetOption("TEXT");
0311 hists_["eventLogger_"]->setBinLabel(1, "Run", 1);
0312 hists_["eventLogger_"]->setBinLabel(2, "Block", 1);
0313 hists_["eventLogger_"]->setBinLabel(3, "Event", 1);
0314 hists_["eventLogger_"]->setBinLabel(4, "pt_{30,loose}(jet1)", 1);
0315 hists_["eventLogger_"]->setBinLabel(5, "pt_{30,loose}(jet2)", 1);
0316 hists_["eventLogger_"]->setBinLabel(6, "pt_{30,loose}(jet3)", 1);
0317 hists_["eventLogger_"]->setBinLabel(7, "pt_{30,loose}(jet4)", 1);
0318 hists_["eventLogger_"]->setBinLabel(8, "M_{W}", 1);
0319 hists_["eventLogger_"]->setBinLabel(9, "M_{Top}", 1);
0320 hists_["eventLogger_"]->setAxisTitle("logged evts", 2);
0321 return;
0322 }
0323
0324 void MonitorEnsemble::fill(const edm::Event& event, const edm::EventSetup& setup) {
0325
0326 edm::Handle<edm::TriggerResults> triggerTable;
0327
0328
0329
0330
0331
0332
0333
0334
0335
0336 edm::Handle<edm::View<reco::Vertex>> pvs;
0337
0338 if (!event.getByToken(pvs_, pvs))
0339 return;
0340 const reco::Vertex& Pvertex = pvs->front();
0341 unsigned int pvMult = 0;
0342 for (edm::View<reco::Vertex>::const_iterator pv = pvs->begin(); pv != pvs->end(); ++pv) {
0343 if (!pvSelect_ || (*pvSelect_)(*pv))
0344 pvMult++;
0345 }
0346 fill("pvMult_", pvMult);
0347
0348
0349
0350
0351
0352
0353
0354
0355
0356 edm::Handle<edm::View<reco::GsfElectron>> elecs;
0357 reco::GsfElectron e;
0358 edm::Handle<double> _rhoHandle;
0359 event.getByLabel(rhoTag, _rhoHandle);
0360 if (!event.getByToken(elecs_, elecs))
0361 return;
0362
0363
0364 edm::Handle<edm::ValueMap<float>> electronId;
0365 if (!electronId_.isUninitialized()) {
0366 if (!event.getByToken(electronId_, electronId)) {
0367 return;
0368 }
0369 }
0370
0371 unsigned int eMult = 0, eMultIso = 0;
0372 std::vector<const reco::GsfElectron*> isoElecs;
0373 for (edm::View<reco::GsfElectron>::const_iterator elec = elecs->begin(); elec != elecs->end(); ++elec) {
0374
0375 if (electronId_.isUninitialized()) {
0376 if (!elecSelect_ || (*elecSelect_)(*elec)) {
0377 double el_ChHadIso = elec->pfIsolationVariables().sumChargedHadronPt;
0378 double el_NeHadIso = elec->pfIsolationVariables().sumNeutralHadronEt;
0379 double el_PhIso = elec->pfIsolationVariables().sumPhotonEt;
0380 double absEta = std::fabs(elec->superCluster()->eta());
0381
0382
0383 double eA = 0;
0384 if (absEta < 1.000)
0385 eA = 0.1703;
0386 else if (absEta < 1.479)
0387 eA = 0.1715;
0388 else if (absEta < 2.000)
0389 eA = 0.1213;
0390 else if (absEta < 2.200)
0391 eA = 0.1230;
0392 else if (absEta < 2.300)
0393 eA = 0.1635;
0394 else if (absEta < 2.400)
0395 eA = 0.1937;
0396 else if (absEta < 5.000)
0397 eA = 0.2393;
0398
0399 double rho = _rhoHandle.isValid() ? (float)(*_rhoHandle) : 0;
0400 double el_pfRelIso = (el_ChHadIso + max(0., el_NeHadIso + el_PhIso - rho * eA)) / elec->pt();
0401
0402
0403 if (eMult == 0) {
0404 fill("elecRelIso_", el_pfRelIso);
0405 fill("elecChHadIso_", el_ChHadIso);
0406 fill("elecNeHadIso_", el_NeHadIso);
0407 fill("elecPhIso_", el_PhIso);
0408 }
0409 ++eMult;
0410
0411 if (!((el_pfRelIso < 0.0588 && absEta < 1.479) || (el_pfRelIso < 0.0571 && absEta > 1.479)))
0412 continue;
0413
0414
0415 if (eMultIso == 0) {
0416 e = *elec;
0417 fill("elecPt_", elec->pt());
0418 fill("elecEta_", elec->eta());
0419 fill("elecPhi_", elec->phi());
0420 }
0421 ++eMultIso;
0422 }
0423 }
0424 }
0425
0426 fill("elecMultTight_", eMultIso);
0427
0428
0429
0430
0431
0432
0433
0434
0435
0436
0437 unsigned int mMult = 0, mTight = 0, mTightId = 0;
0438
0439 edm::Handle<edm::View<reco::Muon>> muons;
0440 edm::View<reco::Muon>::const_iterator muon;
0441
0442 reco::Muon mu;
0443
0444 if (!event.getByToken(muons_, muons))
0445 return;
0446
0447 for (edm::View<reco::Muon>::const_iterator muon = muons->begin(); muon != muons->end(); ++muon) {
0448
0449 if (muon->isGlobalMuon()) {
0450 fill("muonDelZ_", muon->innerTrack()->vz());
0451 fill("muonDelXY_", muon->innerTrack()->vx(), muon->innerTrack()->vy());
0452
0453
0454 if (muon->muonBestTrack().isNonnull()) {
0455 double dxy = muon->muonBestTrack()->dxy(Pvertex.position());
0456 fill("muonDxy_", dxy);
0457 }
0458
0459
0460 if ((!muonSelect_ || (*muonSelect_)(*muon))) {
0461 mMult++;
0462 double chHadPt = muon->pfIsolationR04().sumChargedHadronPt;
0463 double neHadEt = muon->pfIsolationR04().sumNeutralHadronEt;
0464 double phoEt = muon->pfIsolationR04().sumPhotonEt;
0465 double pfRelIso = (chHadPt + max(0., neHadEt + phoEt - 0.5 * muon->pfIsolationR04().sumPUPt)) /
0466 muon->pt();
0467
0468 if (!(muon->isGlobalMuon() && muon->isPFMuon() && muon->globalTrack()->normalizedChi2() < 10. &&
0469 muon->globalTrack()->hitPattern().numberOfValidMuonHits() > 0 && muon->numberOfMatchedStations() > 1 &&
0470 muon->innerTrack()->hitPattern().numberOfValidPixelHits() > 0 &&
0471 muon->innerTrack()->hitPattern().trackerLayersWithMeasurement() > 5 &&
0472 fabs(muon->muonBestTrack()->dxy(Pvertex.position())) < 0.2 &&
0473 fabs(muon->muonBestTrack()->dz(Pvertex.position())) < 0.5))
0474 continue;
0475
0476 if (mTightId == 0) {
0477 fill("muonRelIso_", pfRelIso);
0478 fill("muonChHadIso_", chHadPt);
0479 fill("muonNeHadIso_", neHadEt);
0480 fill("muonPhIso_", phoEt);
0481 }
0482 mTightId++;
0483
0484 if (!(pfRelIso < 0.15))
0485 continue;
0486
0487 if (mTight == 0) {
0488 mu = *(muon);
0489 fill("muonPt_", muon->pt());
0490 fill("muonEta_", muon->eta());
0491 fill("muonPhi_", muon->phi());
0492 }
0493 mTight++;
0494 }
0495 }
0496 }
0497 fill("muonMult_", mMult);
0498 fill("muonMultTight_", mTight);
0499
0500
0501
0502
0503
0504
0505
0506
0507
0508
0509 edm::Handle<reco::JetTagCollection> btagEff, btagPur, btagVtx, btagCSV;
0510 if (includeBTag_) {
0511 if (!event.getByToken(btagCSV_, btagCSV))
0512 return;
0513 }
0514
0515
0516
0517 const reco::JetCorrector* corrector = nullptr;
0518 if (!jetCorrector_.isUninitialized()) {
0519
0520 edm::Handle<reco::JetCorrector> correctorHandle = event.getHandle(jetCorrector_);
0521 if (correctorHandle.isValid()) {
0522 corrector = correctorHandle.product();
0523 } else {
0524 edm::LogVerbatim("SingleTopTChannelLeptonDQM")
0525 << "\n"
0526 << "-----------------------------------------------------------------"
0527 "-------------------- \n"
0528 << " No JetCorrector available from Event:\n"
0529 << " - Jets will not be corrected.\n"
0530 << "-----------------------------------------------------------------"
0531 "-------------------- \n";
0532 }
0533 }
0534
0535 std::vector<reco::Jet> correctedJets;
0536 std::vector<double> JetTagValues;
0537 reco::Jet TaggedJetCand;
0538 unsigned int mult = 0, multLoose = 0, multCSV = 0;
0539 vector<double> bJetDiscVal;
0540
0541 edm::Handle<edm::View<reco::Jet>> jets;
0542 if (!event.getByToken(jets_, jets)) {
0543 return;
0544 }
0545
0546 for (edm::View<reco::Jet>::const_iterator jet = jets->begin(); jet != jets->end(); ++jet) {
0547 bool isLoose = false;
0548
0549
0550
0551
0552
0553
0554
0555
0556
0557
0558 if (dynamic_cast<const reco::PFJet*>(&*jet)) {
0559 reco::PFJet sel = dynamic_cast<const reco::PFJet&>(*jet);
0560 if ((*jetlooseSelection_)(sel))
0561 isLoose = true;
0562 sel.scaleEnergy(corrector ? corrector->correction(*jet) : 1.);
0563 if (!(*jetSelection_)(sel))
0564 continue;
0565 }
0566
0567 reco::Jet monitorJet = *jet;
0568
0569 ++mult;
0570 monitorJet.scaleEnergy(corrector ? corrector->correction(*jet) : 1.);
0571
0572 if (isLoose) {
0573 unsigned int idx = jet - jets->begin();
0574 correctedJets.push_back(monitorJet);
0575 if (includeBTag_) {
0576
0577 edm::RefToBase<reco::Jet> jetRef = jets->refAt(idx);
0578 fill("jetBCSV_", (*btagCSV)[jetRef]);
0579 if ((*btagCSV)[jetRef] > btagCSVWP_) {
0580 if (multCSV == 0) {
0581 TaggedJetCand = monitorJet;
0582 bJetDiscVal.push_back((*btagCSV)[jetRef]);
0583 } else if (multCSV == 1) {
0584 bJetDiscVal.push_back((*btagCSV)[jetRef]);
0585 if (bJetDiscVal[1] > bJetDiscVal[0])
0586 TaggedJetCand = monitorJet;
0587 }
0588 ++multCSV;
0589 }
0590
0591
0592 JetTagValues.push_back((*btagCSV)[jetRef]);
0593 }
0594
0595
0596 if (multLoose == 0) {
0597 fill("jet1Pt_", monitorJet.pt());
0598 fill("jet1Eta_", monitorJet.eta());
0599 };
0600 if (multLoose == 1) {
0601 fill("jet2Pt_", monitorJet.pt());
0602 fill("jet2Eta_", monitorJet.eta());
0603 };
0604 if (multLoose == 2) {
0605 fill("jet3Pt_", monitorJet.pt());
0606 fill("jet3Eta_", monitorJet.eta());
0607 };
0608 if (multLoose == 3) {
0609 fill("jet4Pt_", monitorJet.pt());
0610 fill("jet4Eta_", monitorJet.eta());
0611 };
0612 multLoose++;
0613 }
0614 }
0615 fill("jetMult_", mult);
0616 fill("jetMultLoose_", multLoose);
0617 fill("jetMultBCSVM_", multCSV);
0618
0619
0620
0621
0622
0623
0624
0625
0626
0627
0628 reco::MET mET;
0629 for (std::vector<edm::EDGetTokenT<edm::View<reco::MET>>>::const_iterator met_ = mets_.begin(); met_ != mets_.end();
0630 ++met_) {
0631 edm::Handle<edm::View<reco::MET>> met;
0632 if (!event.getByToken(*met_, met))
0633 continue;
0634 if (met->begin() != met->end()) {
0635 unsigned int idx = met_ - mets_.begin();
0636 if (idx == 0) {
0637 fill("metPflow_", met->begin()->et());
0638 mET = *(met->begin());
0639 }
0640 }
0641 }
0642
0643
0644
0645
0646
0647
0648
0649
0650
0651
0652 Calculate eventKinematics(MAXJETS, WMASS);
0653 double wMass = eventKinematics.massWBoson(correctedJets);
0654 double topMass = eventKinematics.massTopQuark(correctedJets);
0655 if (wMass >= 0 && topMass >= 0) {
0656 fill("massW_", wMass);
0657 fill("massTop_", topMass);
0658 }
0659
0660 if ((lowerEdge_ == -1. && upperEdge_ == -1.) || (lowerEdge_ < wMass && wMass < upperEdge_)) {
0661 if (!triggerTable_.isUninitialized())
0662 fill(event, *triggerTable, "trigger", triggerPaths_);
0663 if (logged_ <= hists_.find("eventLogger_")->second->getNbinsY()) {
0664
0665
0666 fill("eventLogger_", 0.5, logged_ + 0.5, event.eventAuxiliary().run());
0667 fill("eventLogger_", 1.5, logged_ + 0.5, event.eventAuxiliary().luminosityBlock());
0668 fill("eventLogger_", 2.5, logged_ + 0.5, event.eventAuxiliary().event());
0669 if (!correctedJets.empty())
0670 fill("eventLogger_", 3.5, logged_ + 0.5, correctedJets[0].pt());
0671 if (correctedJets.size() > 1)
0672 fill("eventLogger_", 4.5, logged_ + 0.5, correctedJets[1].pt());
0673 if (correctedJets.size() > 2)
0674 fill("eventLogger_", 5.5, logged_ + 0.5, correctedJets[2].pt());
0675 if (correctedJets.size() > 3)
0676 fill("eventLogger_", 6.5, logged_ + 0.5, correctedJets[3].pt());
0677 fill("eventLogger_", 7.5, logged_ + 0.5, wMass);
0678 fill("eventLogger_", 8.5, logged_ + 0.5, topMass);
0679 ++logged_;
0680 }
0681 }
0682 if (multCSV != 0 && mTight == 1) {
0683 double mtW = eventKinematics.tmassWBoson(&mu, mET, TaggedJetCand);
0684 fill("MTWm_", mtW);
0685 double MTT = eventKinematics.tmassTopQuark(&mu, mET, TaggedJetCand);
0686 fill("mMTT_", MTT);
0687 }
0688
0689 if (multCSV != 0 && eMultIso == 1) {
0690 double mtW = eventKinematics.tmassWBoson(&e, mET, TaggedJetCand);
0691 fill("MTWe_", mtW);
0692 double MTT = eventKinematics.tmassTopQuark(&e, mET, TaggedJetCand);
0693 fill("eMTT_", MTT);
0694 }
0695 }
0696 }
0697
0698 SingleTopTChannelLeptonDQM::SingleTopTChannelLeptonDQM(const edm::ParameterSet& cfg)
0699 : vertexSelect_(nullptr),
0700 beamspot_(""),
0701 beamspotSelect_(nullptr),
0702 MuonStep(nullptr),
0703 PFMuonStep(nullptr),
0704 ElectronStep(nullptr),
0705 PFElectronStep(nullptr),
0706 PvStep(nullptr),
0707 METStep(nullptr) {
0708
0709 CaloJetSteps.clear();
0710 PFJetSteps.clear();
0711
0712 edm::ParameterSet presel = cfg.getParameter<edm::ParameterSet>("preselection");
0713 if (presel.existsAs<edm::ParameterSet>("trigger")) {
0714 edm::ParameterSet trigger = presel.getParameter<edm::ParameterSet>("trigger");
0715 triggerTable__ = consumes<edm::TriggerResults>(trigger.getParameter<edm::InputTag>("src"));
0716 triggerPaths_ = trigger.getParameter<std::vector<std::string>>("select");
0717 }
0718 if (presel.existsAs<edm::ParameterSet>("beamspot")) {
0719 edm::ParameterSet beamspot = presel.getParameter<edm::ParameterSet>("beamspot");
0720 beamspot_ = beamspot.getParameter<edm::InputTag>("src");
0721 beamspot__ = consumes<reco::BeamSpot>(beamspot.getParameter<edm::InputTag>("src"));
0722 beamspotSelect_ =
0723 std::make_unique<StringCutObjectSelector<reco::BeamSpot>>(beamspot.getParameter<std::string>("select"));
0724 }
0725
0726 std::vector<edm::ParameterSet> sel = cfg.getParameter<std::vector<edm::ParameterSet>>("selection");
0727
0728 for (unsigned int i = 0; i < sel.size(); ++i) {
0729 selectionOrder_.push_back(sel.at(i).getParameter<std::string>("label"));
0730 selection_[selectionStep(selectionOrder_.back())] =
0731 std::make_pair(sel.at(i),
0732 std::make_unique<SingleTopTChannelLepton::MonitorEnsemble>(
0733 selectionStep(selectionOrder_.back()).c_str(),
0734 cfg.getParameter<edm::ParameterSet>("setup"),
0735 cfg.getParameter<std::vector<edm::ParameterSet>>("selection"),
0736 consumesCollector()));
0737 }
0738 for (std::vector<std::string>::const_iterator selIt = selectionOrder_.begin(); selIt != selectionOrder_.end();
0739 ++selIt) {
0740 std::string key = selectionStep(*selIt), type = objectType(*selIt);
0741 if (selection_.find(key) != selection_.end()) {
0742 using std::unique_ptr;
0743
0744 if (type == "muons") {
0745 MuonStep = std::make_unique<SelectionStep<reco::Muon>>(selection_[key].first, consumesCollector());
0746 }
0747 if (type == "muons/pf") {
0748 PFMuonStep = std::make_unique<SelectionStep<reco::Muon>>(selection_[key].first, consumesCollector());
0749 }
0750 if (type == "elecs") {
0751 ElectronStep = std::make_unique<SelectionStep<reco::GsfElectron>>(selection_[key].first, consumesCollector());
0752 }
0753 if (type == "elecs/pf") {
0754 PFElectronStep = std::make_unique<SelectionStep<reco::GsfElectron>>(selection_[key].first, consumesCollector());
0755 }
0756 if (type == "pvs") {
0757 PvStep = std::make_unique<SelectionStep<reco::Vertex>>(selection_[key].first, consumesCollector());
0758 }
0759 if (type == "jets") {
0760 JetSteps.push_back(std::make_unique<SelectionStep<reco::Jet>>(selection_[key].first, consumesCollector()));
0761 }
0762 if (type == "jets/pf") {
0763 PFJetSteps.push_back(std::make_unique<SelectionStep<reco::PFJet>>(selection_[key].first, consumesCollector()));
0764 }
0765 if (type == "jets/calo") {
0766 CaloJetSteps.push_back(
0767 std::make_unique<SelectionStep<reco::CaloJet>>(selection_[key].first, consumesCollector()));
0768 }
0769 if (type == "met") {
0770 METStep = std::make_unique<SelectionStep<reco::MET>>(selection_[key].first, consumesCollector());
0771 }
0772 }
0773 }
0774 }
0775 void SingleTopTChannelLeptonDQM::bookHistograms(DQMStore::IBooker& ibooker, edm::Run const&, edm::EventSetup const&) {
0776 for (auto selIt = selection_.begin(); selIt != selection_.end(); ++selIt) {
0777 selIt->second.second->book(ibooker);
0778 }
0779 }
0780 void SingleTopTChannelLeptonDQM::analyze(const edm::Event& event, const edm::EventSetup& setup) {
0781 if (!triggerTable__.isUninitialized()) {
0782 edm::Handle<edm::TriggerResults> triggerTable;
0783 if (!event.getByToken(triggerTable__, triggerTable)) {
0784 return;
0785 }
0786 if (!accept(event, *triggerTable, triggerPaths_)) {
0787 return;
0788 }
0789 }
0790 if (!beamspot__.isUninitialized()) {
0791 edm::Handle<reco::BeamSpot> beamspot;
0792 if (!event.getByToken(beamspot__, beamspot)) {
0793 return;
0794 }
0795 if (!(*beamspotSelect_)(*beamspot)) {
0796 return;
0797 }
0798 }
0799
0800
0801 unsigned int nPFJetSteps = -1;
0802 unsigned int nCaloJetSteps = -1;
0803 for (std::vector<std::string>::const_iterator selIt = selectionOrder_.begin(); selIt != selectionOrder_.end();
0804 ++selIt) {
0805 std::string key = selectionStep(*selIt), type = objectType(*selIt);
0806 if (selection_.find(key) != selection_.end()) {
0807 if (type == "empty") {
0808 selection_[key].second->fill(event, setup);
0809 }
0810 if (type == "presel") {
0811 selection_[key].second->fill(event, setup);
0812 }
0813 if (type == "elecs" && ElectronStep != nullptr) {
0814 if (ElectronStep->select(event)) {
0815 selection_[key].second->fill(event, setup);
0816 } else {
0817 break;
0818 }
0819 }
0820 if (type == "elecs/pf" && PFElectronStep != nullptr) {
0821 if (PFElectronStep->select(event, "electron")) {
0822 selection_[key].second->fill(event, setup);
0823
0824 } else {
0825 break;
0826 }
0827 }
0828 if (type == "muons" && MuonStep != nullptr) {
0829 if (MuonStep->select(event)) {
0830 selection_[key].second->fill(event, setup);
0831 } else {
0832 break;
0833 }
0834 }
0835 if (type == "muons/pf" && PFMuonStep != nullptr) {
0836 if (PFMuonStep->select(event, "muon")) {
0837 selection_[key].second->fill(event, setup);
0838 } else {
0839 break;
0840 }
0841 }
0842
0843
0844
0845
0846
0847
0848
0849
0850
0851
0852 if (type == "jets/pf") {
0853 nPFJetSteps++;
0854 if (PFJetSteps[nPFJetSteps]) {
0855 if (PFJetSteps[nPFJetSteps]->select(event, setup)) {
0856 selection_[key].second->fill(event, setup);
0857 } else {
0858 break;
0859 }
0860 }
0861 }
0862 if (type == "jets/calo") {
0863 nCaloJetSteps++;
0864 if (CaloJetSteps[nCaloJetSteps]) {
0865 if (CaloJetSteps[nCaloJetSteps]->select(event, setup)) {
0866 selection_[key].second->fill(event, setup);
0867 } else {
0868 break;
0869 }
0870 }
0871 }
0872 if (type == "met" && METStep != nullptr) {
0873 if (METStep->select(event)) {
0874 selection_[key].second->fill(event, setup);
0875 } else {
0876 break;
0877 }
0878 }
0879 }
0880 }
0881 }