File indexing completed on 2024-04-06 12:10:12
0001 #include "DQMServices/Examples/interface/DQMExample_Step1.h"
0002 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0003
0004
0005 #include "DataFormats/Math/interface/LorentzVector.h"
0006 #include "DataFormats/Math/interface/deltaPhi.h"
0007 #include "DataFormats/Math/interface/deltaR.h"
0008 #include "TLorentzVector.h"
0009
0010 #include <cmath>
0011 #include <cstdio>
0012 #include <iomanip>
0013 #include <iostream>
0014 #include <sstream>
0015 #include <string>
0016
0017
0018
0019
0020
0021 DQMExample_Step1::DQMExample_Step1(const edm::ParameterSet &ps) {
0022 edm::LogInfo("DQMExample_Step1") << "Constructor DQMExample_Step1::DQMExample_Step1 " << std::endl;
0023
0024
0025 theElectronCollection_ = consumes<reco::GsfElectronCollection>(ps.getParameter<edm::InputTag>("electronCollection"));
0026 theCaloJetCollection_ = consumes<reco::CaloJetCollection>(ps.getParameter<edm::InputTag>("caloJetCollection"));
0027 thePfMETCollection_ = consumes<reco::PFMETCollection>(ps.getParameter<edm::InputTag>("pfMETCollection"));
0028 theConversionCollection_ =
0029 consumes<reco::ConversionCollection>(ps.getParameter<edm::InputTag>("conversionsCollection"));
0030 thePVCollection_ = consumes<reco::VertexCollection>(ps.getParameter<edm::InputTag>("PVCollection"));
0031 theBSCollection_ = consumes<reco::BeamSpot>(ps.getParameter<edm::InputTag>("beamSpotCollection"));
0032 triggerEvent_ = consumes<trigger::TriggerEvent>(ps.getParameter<edm::InputTag>("TriggerEvent"));
0033 triggerResults_ = consumes<edm::TriggerResults>(ps.getParameter<edm::InputTag>("TriggerResults"));
0034 triggerFilter_ = ps.getParameter<edm::InputTag>("TriggerFilter");
0035 triggerPath_ = ps.getParameter<std::string>("TriggerPath");
0036
0037
0038 ptThrL1_ = ps.getUntrackedParameter<double>("PtThrL1");
0039 ptThrL2_ = ps.getUntrackedParameter<double>("PtThrL2");
0040 ptThrJet_ = ps.getUntrackedParameter<double>("PtThrJet");
0041 ptThrMet_ = ps.getUntrackedParameter<double>("PtThrMet");
0042 }
0043
0044
0045
0046
0047 DQMExample_Step1::~DQMExample_Step1() {
0048 edm::LogInfo("DQMExample_Step1") << "Destructor DQMExample_Step1::~DQMExample_Step1 " << std::endl;
0049 }
0050
0051
0052
0053
0054
0055 void DQMExample_Step1::dqmBeginRun(edm::Run const &, edm::EventSetup const &) {
0056 edm::LogInfo("DQMExample_Step1") << "DQMExample_Step1::beginRun" << std::endl;
0057 }
0058
0059
0060
0061
0062 void DQMExample_Step1::bookHistograms(DQMStore::IBooker &ibooker_, edm::Run const &, edm::EventSetup const &) {
0063 edm::LogInfo("DQMExample_Step1") << "DQMExample_Step1::bookHistograms" << std::endl;
0064
0065
0066 bookHistos(ibooker_);
0067 }
0068
0069
0070
0071
0072 void DQMExample_Step1::analyze(edm::Event const &e, edm::EventSetup const &eSetup) {
0073 edm::LogInfo("DQMExample_Step1") << "DQMExample_Step1::analyze" << std::endl;
0074
0075
0076
0077
0078 edm::Handle<reco::VertexCollection> vertexHandle;
0079 e.getByToken(thePVCollection_, vertexHandle);
0080 if (!vertexHandle.isValid()) {
0081 edm::LogError("DQMClientExample") << "invalid collection: vertex"
0082 << "\n";
0083 return;
0084 }
0085
0086 int vertex_number = vertexHandle->size();
0087 reco::VertexCollection::const_iterator v = vertexHandle->begin();
0088
0089 math::XYZPoint PVPoint(-999, -999, -999);
0090 if (vertex_number != 0)
0091 PVPoint = math::XYZPoint(v->position().x(), v->position().y(), v->position().z());
0092
0093 PVPoint_ = PVPoint;
0094
0095
0096
0097
0098 edm::Handle<reco::PFMETCollection> pfMETCollection;
0099 e.getByToken(thePfMETCollection_, pfMETCollection);
0100 if (!pfMETCollection.isValid()) {
0101 edm::LogError("DQMClientExample") << "invalid collection: MET"
0102 << "\n";
0103 return;
0104 }
0105
0106
0107
0108 edm::Handle<reco::GsfElectronCollection> electronCollection;
0109 e.getByToken(theElectronCollection_, electronCollection);
0110 if (!electronCollection.isValid()) {
0111 edm::LogError("DQMClientExample") << "invalid collection: electrons"
0112 << "\n";
0113 return;
0114 }
0115
0116 float nEle = 0;
0117 int posEle = 0, negEle = 0;
0118 const reco::GsfElectron *ele1 = nullptr;
0119 const reco::GsfElectron *ele2 = nullptr;
0120 for (reco::GsfElectronCollection::const_iterator recoElectron = electronCollection->begin();
0121 recoElectron != electronCollection->end();
0122 ++recoElectron) {
0123
0124 if (MediumEle(e, eSetup, *recoElectron)) {
0125 if (!ele1 && recoElectron->pt() > ptThrL1_)
0126 ele1 = &(*recoElectron);
0127
0128 else if (!ele2 && recoElectron->pt() > ptThrL2_)
0129 ele2 = &(*recoElectron);
0130 }
0131
0132 if (recoElectron->charge() == 1)
0133 posEle++;
0134 else if (recoElectron->charge() == -1)
0135 negEle++;
0136
0137 }
0138
0139 nEle = posEle + negEle;
0140
0141
0142
0143
0144 edm::Handle<reco::CaloJetCollection> caloJetCollection;
0145 e.getByToken(theCaloJetCollection_, caloJetCollection);
0146 if (!caloJetCollection.isValid()) {
0147 edm::LogError("DQMClientExample") << "invalid collection: jets"
0148 << "\n";
0149 return;
0150 }
0151
0152 int nJet = 0;
0153 const reco::CaloJet *jet1 = nullptr;
0154 const reco::CaloJet *jet2 = nullptr;
0155
0156 for (reco::CaloJetCollection::const_iterator i_calojet = caloJetCollection->begin();
0157 i_calojet != caloJetCollection->end();
0158 ++i_calojet) {
0159
0160 if (ele1)
0161 if (Distance(*i_calojet, *ele1) < 0.3)
0162 continue;
0163
0164 if (ele2)
0165 if (Distance(*i_calojet, *ele2) < 0.3)
0166 continue;
0167
0168 if (i_calojet->pt() < ptThrJet_)
0169 continue;
0170
0171 nJet++;
0172
0173 if (!jet1)
0174 jet1 = &(*i_calojet);
0175
0176 else if (!jet2)
0177 jet2 = &(*i_calojet);
0178 }
0179
0180
0181
0182
0183
0184
0185 edm::Handle<edm::TriggerResults> hltresults;
0186 e.getByToken(triggerResults_, hltresults);
0187
0188 if (!hltresults.isValid()) {
0189 edm::LogError("DQMClientExample") << "invalid collection: TriggerResults"
0190 << "\n";
0191 return;
0192 }
0193
0194 bool hasFired = false;
0195 const edm::TriggerNames &trigNames = e.triggerNames(*hltresults);
0196 unsigned int numTriggers = trigNames.size();
0197
0198 for (unsigned int hltIndex = 0; hltIndex < numTriggers; ++hltIndex) {
0199 if (trigNames.triggerName(hltIndex) == triggerPath_ && hltresults->wasrun(hltIndex) && hltresults->accept(hltIndex))
0200 hasFired = true;
0201 }
0202
0203
0204 edm::Handle<trigger::TriggerEvent> triggerEvent;
0205 e.getByToken(triggerEvent_, triggerEvent);
0206 if (triggerEvent.failedToGet()) {
0207 edm::LogError("DQMClientExample") << "invalid collection: TriggerEvent"
0208 << "\n";
0209 return;
0210 }
0211
0212 reco::Particle *ele1_HLT = nullptr;
0213 int nEle_HLT = 0;
0214
0215 size_t filterIndex = triggerEvent->filterIndex(triggerFilter_);
0216 trigger::TriggerObjectCollection triggerObjects = triggerEvent->getObjects();
0217 if (!(filterIndex >= triggerEvent->sizeFilters())) {
0218 const trigger::Keys &keys = triggerEvent->filterKeys(filterIndex);
0219 std::vector<reco::Particle> triggeredEle;
0220
0221 for (size_t j = 0; j < keys.size(); ++j) {
0222 trigger::TriggerObject foundObject = triggerObjects[keys[j]];
0223 if (abs(foundObject.particle().pdgId()) != 11)
0224 continue;
0225
0226 triggeredEle.push_back(foundObject.particle());
0227 ++nEle_HLT;
0228 }
0229
0230 if (!triggeredEle.empty())
0231 ele1_HLT = &(triggeredEle.at(0));
0232 }
0233
0234
0235
0236
0237
0238
0239 h_vertex_number->Fill(vertex_number);
0240
0241
0242 h_pfMet->Fill(pfMETCollection->begin()->et());
0243
0244
0245 h_eMultiplicity->Fill(nEle);
0246 h_jMultiplicity->Fill(nJet);
0247 h_eMultiplicity_HLT->Fill(nEle_HLT);
0248
0249
0250 if (ele1) {
0251 h_ePt_leading->Fill(ele1->pt());
0252 h_eEta_leading->Fill(ele1->eta());
0253 h_ePhi_leading->Fill(ele1->phi());
0254 }
0255 if (ele1_HLT) {
0256 h_ePt_leading_HLT->Fill(ele1_HLT->pt());
0257 h_eEta_leading_HLT->Fill(ele1_HLT->eta());
0258 h_ePhi_leading_HLT->Fill(ele1_HLT->phi());
0259 }
0260
0261 if (jet1) {
0262 h_jPt_leading->Fill(jet1->pt());
0263 h_jEta_leading->Fill(jet1->eta());
0264 h_jPhi_leading->Fill(jet1->phi());
0265 }
0266
0267
0268 if (ele1 && ele1_HLT && deltaR(*ele1_HLT, *ele1) < 0.3 && hasFired == true) {
0269 h_ePt_leading_matched->Fill(ele1->pt());
0270 h_eEta_leading_matched->Fill(ele1->eta());
0271 h_ePhi_leading_matched->Fill(ele1->phi());
0272
0273 h_ePt_leading_HLT_matched->Fill(ele1_HLT->pt());
0274 h_eEta_leading_HLT_matched->Fill(ele1_HLT->eta());
0275 h_ePhi_leading_HLT_matched->Fill(ele1_HLT->phi());
0276
0277 h_ePt_diff->Fill(ele1->pt() - ele1_HLT->pt());
0278 }
0279 }
0280
0281
0282
0283
0284
0285
0286
0287
0288
0289 void DQMExample_Step1::bookHistos(DQMStore::IBooker &ibooker_) {
0290 ibooker_.cd();
0291 ibooker_.setCurrentFolder("Physics/TopTest");
0292
0293 h_vertex_number = ibooker_.book1D("Vertex_number", "Number of event vertices in collection", 40, -0.5, 39.5);
0294
0295 h_pfMet = ibooker_.book1D("pfMet", "Pf Missing E_{T}; GeV", 20, 0.0, 100);
0296
0297 h_eMultiplicity = ibooker_.book1D("NElectrons", "# of electrons per event", 10, 0., 10.);
0298 h_ePt_leading_matched = ibooker_.book1D("ElePt_leading_matched", "Pt of leading electron", 50, 0., 100.);
0299 h_eEta_leading_matched = ibooker_.book1D("EleEta_leading_matched", "Eta of leading electron", 50, -5., 5.);
0300 h_ePhi_leading_matched = ibooker_.book1D("ElePhi_leading_matched", "Phi of leading electron", 50, -3.5, 3.5);
0301
0302 h_ePt_leading = ibooker_.book1D("ElePt_leading", "Pt of leading electron", 50, 0., 100.);
0303 h_eEta_leading = ibooker_.book1D("EleEta_leading", "Eta of leading electron", 50, -5., 5.);
0304 h_ePhi_leading = ibooker_.book1D("ElePhi_leading", "Phi of leading electron", 50, -3.5, 3.5);
0305
0306 h_jMultiplicity = ibooker_.book1D("NJets", "# of electrons per event", 10, 0., 10.);
0307 h_jPt_leading = ibooker_.book1D("JetPt_leading", "Pt of leading Jet", 150, 0., 300.);
0308 h_jEta_leading = ibooker_.book1D("JetEta_leading", "Eta of leading Jet", 50, -5., 5.);
0309 h_jPhi_leading = ibooker_.book1D("JetPhi_leading", "Phi of leading Jet", 50, -3.5, 3.5);
0310
0311 h_eMultiplicity_HLT = ibooker_.book1D("NElectrons_HLT", "# of electrons per event @HLT", 10, 0., 10.);
0312 h_ePt_leading_HLT = ibooker_.book1D("ElePt_leading_HLT", "Pt of leading electron @HLT", 50, 0., 100.);
0313 h_eEta_leading_HLT = ibooker_.book1D("EleEta_leading_HLT", "Eta of leading electron @HLT", 50, -5., 5.);
0314 h_ePhi_leading_HLT = ibooker_.book1D("ElePhi_leading_HLT", "Phi of leading electron @HLT", 50, -3.5, 3.5);
0315
0316 h_ePt_leading_HLT_matched = ibooker_.book1D("ElePt_leading_HLT_matched", "Pt of leading electron @HLT", 50, 0., 100.);
0317 h_eEta_leading_HLT_matched =
0318 ibooker_.book1D("EleEta_leading_HLT_matched", "Eta of leading electron @HLT", 50, -5., 5.);
0319 h_ePhi_leading_HLT_matched =
0320 ibooker_.book1D("ElePhi_leading_HLT_matched", "Phi of leading electron @HLT", 50, -3.5, 3.5);
0321
0322 h_ePt_diff = ibooker_.book1D("ElePt_diff_matched", "pT(RECO) - pT(HLT) for mathed candidates", 100, -10, 10.);
0323
0324 ibooker_.cd();
0325 }
0326
0327
0328
0329
0330
0331 double DQMExample_Step1::Distance(const reco::Candidate &c1, const reco::Candidate &c2) { return deltaR(c1, c2); }
0332
0333 double DQMExample_Step1::DistancePhi(const reco::Candidate &c1, const reco::Candidate &c2) {
0334 return deltaPhi(c1.p4().phi(), c2.p4().phi());
0335 }
0336
0337
0338 double DQMExample_Step1::calcDeltaPhi(double phi1, double phi2) {
0339 double deltaPhi = phi1 - phi2;
0340 if (deltaPhi < 0)
0341 deltaPhi = -deltaPhi;
0342 if (deltaPhi > 3.1415926) {
0343 deltaPhi = 2 * 3.1415926 - deltaPhi;
0344 }
0345 return deltaPhi;
0346 }
0347
0348
0349
0350
0351
0352 bool DQMExample_Step1::MediumEle(const edm::Event &iEvent,
0353 const edm::EventSetup &iESetup,
0354 const reco::GsfElectron &electron) {
0355
0356 edm::Handle<reco::ConversionCollection> conversions_h;
0357 iEvent.getByToken(theConversionCollection_, conversions_h);
0358
0359 bool isMediumEle = false;
0360
0361 float pt = electron.pt();
0362 float eta = electron.eta();
0363
0364 int isEB = electron.isEB();
0365 float sigmaIetaIeta = electron.sigmaIetaIeta();
0366 float DetaIn = electron.deltaEtaSuperClusterTrackAtVtx();
0367 float DphiIn = electron.deltaPhiSuperClusterTrackAtVtx();
0368 float HOverE = electron.hadronicOverEm();
0369 float ooemoop = (1.0 / electron.ecalEnergy() - electron.eSuperClusterOverP() / electron.ecalEnergy());
0370
0371 int mishits = electron.gsfTrack()->hitPattern().numberOfLostHits(reco::HitPattern::MISSING_INNER_HITS);
0372 int nAmbiguousGsfTracks = electron.ambiguousGsfTracksSize();
0373
0374 reco::GsfTrackRef eleTrack = electron.gsfTrack();
0375 float dxy = eleTrack->dxy(PVPoint_);
0376 float dz = eleTrack->dz(PVPoint_);
0377
0378 edm::Handle<reco::BeamSpot> BSHandle;
0379 iEvent.getByToken(theBSCollection_, BSHandle);
0380 const reco::BeamSpot BS = *BSHandle;
0381
0382 bool isConverted = ConversionTools::hasMatchedConversion(electron, *conversions_h, BS.position());
0383
0384
0385 if ((pt > 12.) && (fabs(eta) < 2.5) &&
0386 (((isEB == 1) && (fabs(DetaIn) < 0.004)) || ((isEB == 0) && (fabs(DetaIn) < 0.007))) &&
0387 (((isEB == 1) && (fabs(DphiIn) < 0.060)) || ((isEB == 0) && (fabs(DphiIn) < 0.030))) &&
0388 (((isEB == 1) && (sigmaIetaIeta < 0.010)) || ((isEB == 0) && (sigmaIetaIeta < 0.030))) &&
0389 (((isEB == 1) && (HOverE < 0.120)) || ((isEB == 0) && (HOverE < 0.100))) &&
0390 (((isEB == 1) && (fabs(ooemoop) < 0.050)) || ((isEB == 0) && (fabs(ooemoop) < 0.050))) &&
0391 (((isEB == 1) && (fabs(dxy) < 0.020)) || ((isEB == 0) && (fabs(dxy) < 0.020))) &&
0392 (((isEB == 1) && (fabs(dz) < 0.100)) || ((isEB == 0) && (fabs(dz) < 0.100))) &&
0393 (((isEB == 1) && (!isConverted)) || ((isEB == 0) && (!isConverted))) && (mishits == 0) &&
0394 (nAmbiguousGsfTracks == 0))
0395 isMediumEle = true;
0396
0397 return isMediumEle;
0398 }
