Physics/src/SMPDQM.cc

0001 #include "DQM/Physics/src/SMPDQM.h"
0002
0003 using namespace std;
0004 using namespace reco;
0005
0006 struct SortByPt
0007
0008 {
0009   bool operator()(const TLorentzVector& a, const TLorentzVector& b) const { return a.Pt() > b.Pt(); }
0010 };
0011 SMPDQM::SMPDQM(const edm::ParameterSet& iConfig) {
0012   //now do what ever initialization is needed
0013   muons_ = consumes<reco::MuonCollection>(iConfig.getParameter<edm::InputTag>("muonCollection"));
0014   pvs_ = consumes<edm::View<reco::Vertex>>(iConfig.getParameter<edm::InputTag>("pvs"));
0015
0016   elecs_ = consumes<reco::GsfElectronCollection>(iConfig.getParameter<edm::InputTag>("elecCollection"));
0017   jets_ = consumes<edm::View<reco::PFJet>>(iConfig.getParameter<edm::InputTag>("jets"));
0018
0019   for (edm::InputTag const& tag : iConfig.getParameter<std::vector<edm::InputTag>>("mets"))
0020     mets_.push_back(consumes<edm::View<reco::MET>>(tag));
0021 }
0022
0023 SMPDQM::~SMPDQM() {
0024   // do anything here that needs to be done at desctruction time
0025   // (e.g. close files, deallocate resources etc.)
0026 }
0027 // ------------ method called for each event  ------------
0028 void SMPDQM::bookHistograms(DQMStore::IBooker& bei, edm::Run const&, edm::EventSetup const&) {
0029   bei.setCurrentFolder("Physics/SMP");
0030
0031   NPV = bei.book1D("NPV", "Number of primary vertices", 40, 0., 80.);
0032   MET = bei.book1D("MET", "MET", 100, 0.0, 200);
0033   METphi = bei.book1D("METphi", "#phi(MET)", 50, -3.14, 3.14);
0034
0035   pt_muons = bei.book1D("pt_muons", "p_{T}(muons)", 40, 0., 200.);
0036   eta_muons = bei.book1D("eta_muons", "#eta(muons)", 50, -5., 5.);
0037   phi_muons = bei.book1D("phi_muons", "#phi(muons)", 32, -3.2, 3.2);
0038   muIso_CombRelIso03 = bei.book1D("muIso_CombRelIso03", "Iso_{rel}^{#mu}", 20, 0., 1.);
0039   Nmuons = bei.book1D("Nmuons", "Number of muons", 20, 0., 10.);
0040   isGlobalmuon = bei.book1D("isGlobalmuon", "isGlobalmuon", 2, 0, 1);
0041   isTrackermuon = bei.book1D("isTrackermuon", "isTrackermuon", 2, 0, 1);
0042   isStandalonemuon = bei.book1D("isStandalonemuon", "isStandalonemuon", 2, 0, 1);
0043   isPFmuon = bei.book1D("isPFmuon", "isPFmuon", 2, 0, 1);
0044   muIso_TrackerBased03 = bei.book1D("muIso_TrackerBased03", "Iso_{trk03}^{#mu}", 20, 0, 10);
0045
0046   Nelecs = bei.book1D("Nelecs", "Number of electrons", 20, 0., 10.);
0047   HoverE_elecs = bei.book1D("HoverE_elecs", "HoverE", 50, 0., 1.);
0048   pt_elecs = bei.book1D("pt_elecs", "p_{T}(elecs)", 40, 0., 200.);
0049   eta_elecs = bei.book1D("eta_elecs", "#eta(elecs)", 50, -5., 5.);
0050   phi_elecs = bei.book1D("phi_elecs", "#phielecs)", 32, -3.2, 3.2);
0051   elIso_cal = bei.book1D("elIso_cal", "Iso_{cal}^{el}", 21, -1., 20.);
0052   elIso_trk = bei.book1D("elIso_trk", "Iso_{trk}^{el}", 21, -2., 40.);
0053   elIso_CombRelIso = bei.book1D("elIso_CombRelIso", "Iso_{rel}^{el}", 20, 0., 1.);
0054
0055   PFJetpt = bei.book1D("PFJetpt", "p_{T}(jets)", 100, 0.0, 100);
0056   PFJeteta = bei.book1D("PFJeteta", "#eta(jets)", 50, -2.5, 2.5);
0057   PFJetphi = bei.book1D("PFJetphi", "#phi(jets)", 50, -3.14, 3.14);
0058   PFJetMulti = bei.book1D("PFJetMulti", "jet multiplicity", 5, -0.5, 4.5);
0059   PFJetRapidity = bei.book1D("PFJetRapidity", "y(jets)", 50, -6.0, 6.0);
0060   mjj = bei.book1D("mjj", "m_{jj}", 100, 0, 1000);
0061   detajj = bei.book1D("detajj", "#Delta#etajj", 20, 0, 5);
0062
0063   dphi_lepMET = bei.book1D("dphi_lepMET", "#Delta#phi(lep,MET)", 60, -3.2, 3.2);
0064   mass_lepMET = bei.book1D("mass_lepMET", "m(lep,MET)", 200, 0, 200);
0065   pt_lepMET = bei.book1D("pt_lepMET", "p_{T}(lep,MET)", 200, 0, 200);
0066   detall = bei.book1D("detall", "#Delta#etall", 20, -5, 5);
0067   dphill = bei.book1D("dphill", "#Delta#phill", 20, -6.4, 6.4);
0068   mll = bei.book1D("mll", "mll", 200, 0, 200);
0069   etall = bei.book1D("etall", "#Delta#etall", 60, -6, 6);
0070   ptll = bei.book1D("ptll", "p_{T}ll", 200, 0, 200);
0071   mjj = bei.book1D("mjj", "mjj", 100, 0, 1000);
0072   detajj = bei.book1D("detajj", "#Delta#etajj", 20, 0, 5);
0073
0074   dphi_lepjet1 = bei.book1D("dphi_lepjet1", "#Delta#phi(lep,jet1)", 60, -3.2, 3.2);
0075
0076   dphi_lep1jet1 = bei.book1D("dphi_lep1jet1", "#Delta#phi(lep1,jet1)", 60, -3.2, 3.2);
0077   dphi_lep2jet1 = bei.book1D("dphi_lep2jet1", "#Delta#phi(lep2,jet1)", 60, -3.2, 3.2);
0078 }
0079 void SMPDQM::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0080   std::vector<TLorentzVector> recoPFJets;
0081   recoPFJets.clear();
0082   TLorentzVector imet;
0083   imet.Clear();
0084   std::vector<TLorentzVector> selected_recoPFJets;
0085   selected_recoPFJets.clear();
0086   std::vector<TLorentzVector> selected_lep;
0087   selected_lep.clear();
0088
0089   for (std::vector<edm::EDGetTokenT<edm::View<reco::MET>>>::const_iterator met_ = mets_.begin(); met_ != mets_.end();
0090        ++met_) {
0091     edm::Handle<edm::View<reco::MET>> met;
0092     if (!iEvent.getByToken(*met_, met))
0093       continue;
0094     if (met->begin() != met->end()) {
0095       MET->Fill(met->begin()->et());
0096       METphi->Fill(met->begin()->phi());
0097       imet.SetPtEtaPhiM(met->begin()->et(), 0., met->begin()->phi(), 0.0);
0098     }
0099   }
0100
0101   // Muons
0102
0103   edm::Handle<edm::View<reco::Vertex>> pvs;
0104   if (!iEvent.getByToken(pvs_, pvs)) {
0105     return;
0106   }
0107
0108   unsigned int pvMult = 0;
0109
0110   for (edm::View<reco::Vertex>::const_iterator pv = pvs->begin(); pv != pvs->end(); ++pv) {
0111     if (pv->position().Rho() < 2 && abs(pv->position().z()) <= 24. && pv->ndof() > 4 && !pv->isFake()) {
0112       pvMult++;
0113     }
0114   }
0115   NPV->Fill(pvMult);
0116
0117   edm::Handle<reco::MuonCollection> muons;
0118   iEvent.getByToken(muons_, muons);
0119   reco::MuonCollection::const_iterator mu;
0120   if (!muons.failedToGet()) {
0121     Nmuons->Fill(muons->size());
0122
0123     for (mu = muons->begin(); mu != muons->end(); ++mu) {
0124       if (mu->pt() < 3.0)
0125         continue;
0126       TLorentzVector Mu;
0127       Mu.SetPtEtaPhiM(mu->pt(), mu->eta(), mu->phi(), 0.0);
0128       selected_lep.push_back(Mu);
0129       pt_muons->Fill(mu->pt());
0130       eta_muons->Fill(mu->eta());
0131       phi_muons->Fill(mu->phi());
0132       isGlobalmuon->Fill(mu->isGlobalMuon());
0133       isTrackermuon->Fill(mu->isTrackerMuon());
0134       isStandalonemuon->Fill(mu->isStandAloneMuon());
0135       isPFmuon->Fill(mu->isPFMuon());
0136
0137       reco::MuonIsolation muIso03 = mu->isolationR03();
0138       double muonCombRelIso = 1.;
0139
0140       muonCombRelIso = (muIso03.emEt + muIso03.hadEt + muIso03.hoEt + muIso03.sumPt) / mu->pt();
0141
0142       muIso_TrackerBased03->Fill(muIso03.sumPt / mu->pt());
0143       muIso_CombRelIso03->Fill(muonCombRelIso);
0144
0145     }  //size of muons
0146
0147   }  // muons
0148
0149   // electrons
0150
0151   edm::Handle<reco::GsfElectronCollection> elecs;
0152   iEvent.getByToken(elecs_, elecs);
0153   reco::GsfElectronCollection::const_iterator elec;
0154
0155   if (!elecs.failedToGet()) {
0156     Nelecs->Fill(elecs->size());
0157
0158     for (elec = elecs->begin(); elec != elecs->end(); ++elec) {
0159       if (elec->pt() < 5.0)
0160         continue;
0161       TLorentzVector El;
0162       El.SetPtEtaPhiM(elec->pt(), elec->eta(), elec->phi(), 0.0);
0163       selected_lep.push_back(El);
0164
0165       HoverE_elecs->Fill(elec->hcalOverEcal());
0166       pt_elecs->Fill(elec->pt());
0167       eta_elecs->Fill(elec->eta());
0168       phi_elecs->Fill(elec->phi());
0169
0170       reco::GsfTrackRef track = elec->gsfTrack();
0171       reco::GsfElectron::IsolationVariables elecIso = elec->dr03IsolationVariables();
0172
0173       double elecCombRelIso = 1.;
0174
0175       elecCombRelIso = (elecIso.ecalRecHitSumEt + elecIso.hcalRecHitSumEt[0] + elecIso.tkSumPt) / elec->pt();
0176       elIso_CombRelIso->Fill(elecCombRelIso);
0177       elIso_cal->Fill(elecIso.ecalRecHitSumEt);
0178       elIso_trk->Fill(elecIso.tkSumPt);
0179     }
0180
0181   }  // electrons
0182   // jets
0183
0184   edm::Handle<edm::View<reco::PFJet>> jets;
0185   if (!iEvent.getByToken(jets_, jets)) {
0186     return;
0187   }
0188
0189   for (edm::View<reco::PFJet>::const_iterator jet = jets->begin(); jet != jets->end(); ++jet) {
0190     if (jet->pt() < 15.0)
0191       continue;
0192     TLorentzVector ijet;
0193     ijet.SetPtEtaPhiM(jet->pt(), jet->eta(), jet->phi(), jet->mass());
0194     recoPFJets.push_back(ijet);
0195   }
0196
0197   std::sort(recoPFJets.begin(), recoPFJets.end(), SortByPt());
0198   std::sort(selected_lep.begin(), selected_lep.end(), SortByPt());
0199
0200   for (unsigned int i = 0; i < recoPFJets.size(); i++) {
0201     bool goodjet = false;
0202     for (unsigned int j = 0; j < selected_lep.size(); j++) {
0203       if (recoPFJets[i].DeltaR(selected_lep[j]) > 0.4) {
0204         goodjet = true;
0205         continue;
0206       } else {
0207         goodjet = false;
0208         break;
0209       }
0210     }
0211     if (goodjet) {
0212       TLorentzVector temp;
0213       temp.Clear();
0214       temp.SetPtEtaPhiM(recoPFJets[i].Pt(), recoPFJets[i].Eta(), recoPFJets[i].Phi(), recoPFJets[i].M());
0215       selected_recoPFJets.push_back(temp);
0216     }
0217   }
0218
0219   std::sort(selected_recoPFJets.begin(), selected_recoPFJets.end(), SortByPt());  // for safety
0220   int njet = 0;
0221   for (unsigned int k = 0; k < selected_recoPFJets.size(); k++) {
0222     if (k > 4)
0223       break;
0224     else {
0225       njet++;
0226       PFJetpt->Fill(selected_recoPFJets.at(k).Pt());
0227       PFJeteta->Fill(selected_recoPFJets.at(k).Eta());
0228       PFJetphi->Fill(selected_recoPFJets.at(k).Phi());
0229       PFJetRapidity->Fill(selected_recoPFJets.at(k).Rapidity());
0230     }
0231   }
0232   PFJetMulti->Fill(njet);
0233
0234   // now we have selected jet and lepton collections
0235
0236   if (selected_lep.size() > 1) {
0237     detall->Fill(selected_lep[0].Eta() - selected_lep[1].Eta());
0238     dphill->Fill(selected_lep[0].DeltaPhi(selected_lep[1]));
0239     mll->Fill((selected_lep[0] + selected_lep[1]).M());
0240     ptll->Fill((selected_lep[0] + selected_lep[1]).Pt());
0241     etall->Fill((selected_lep[0] + selected_lep[1]).Eta());
0242     if (!selected_recoPFJets.empty()) {
0243       dphi_lep1jet1->Fill(selected_recoPFJets[0].DeltaPhi(selected_lep[0]));
0244       dphi_lep2jet1->Fill(selected_recoPFJets[0].DeltaPhi(selected_lep[1]));
0245     }
0246   }
0247
0248   else if (selected_lep.size() == 1) {
0249     dphi_lepMET->Fill(selected_lep[0].DeltaPhi(imet));
0250     mass_lepMET->Fill((selected_lep[0] + imet).M());
0251     pt_lepMET->Fill((selected_lep[0] + imet).Pt());
0252     if (!selected_recoPFJets.empty()) {
0253       dphi_lepjet1->Fill(selected_recoPFJets[0].DeltaPhi(selected_lep[0]));
0254     }
0255   }  // W case
0256
0257   else {
0258     // std::cout << "zero lepton case" << endl;
0259   }
0260   if (selected_recoPFJets.size() > 1) {
0261     detajj->Fill(abs(selected_recoPFJets[0].Eta() - selected_recoPFJets[1].Eta()));
0262     mjj->Fill((selected_recoPFJets[0] + selected_recoPFJets[1]).M());
0263   }
0264
0265 }  // analyze
0266
0267 //define this as a plug-in
0268 //DEFINE_FWK_MODULE(SMPDQM);
0269
0270 //dilepton eta,pt  and phi using lorentz vectors
0271 //first five jets
0272 // dphi(lepton jet)
0273 //vbf detajj, mjj