Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:32:11

0001 #include "Validation/EventGenerator/interface/TTbar_GenJetAnalyzer.h"
0002 #include "Validation/EventGenerator/interface/DQMHelper.h"
0003 
0004 TTbar_GenJetAnalyzer::TTbar_GenJetAnalyzer(const edm::ParameterSet& iConfig)
0005     : jets_(iConfig.getParameter<edm::InputTag>("jets")),
0006       genEventInfoProductTag_(iConfig.getParameter<edm::InputTag>("genEventInfoProductTag")) {
0007   genEventInfoProductTagToken_ = consumes<GenEventInfoProduct>(genEventInfoProductTag_);
0008   jetsToken_ = consumes<std::vector<reco::GenJet> >(jets_);
0009 }
0010 
0011 TTbar_GenJetAnalyzer::~TTbar_GenJetAnalyzer() {}
0012 
0013 void TTbar_GenJetAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0014   using namespace edm;
0015 
0016   // --- the MC weights ---
0017   Handle<GenEventInfoProduct> evt_info;
0018   iEvent.getByToken(genEventInfoProductTagToken_, evt_info);
0019   if (!evt_info.isValid())
0020     return;
0021   weight = evt_info->weight();
0022 
0023   // Gather information in the GenJet collection
0024   edm::Handle<std::vector<reco::GenJet> > jets;
0025   iEvent.getByToken(jetsToken_, jets);
0026 
0027   if (!jets.isValid())
0028     return;
0029   // loop Jet collection and fill histograms
0030   int njets = 0;
0031   for (std::vector<reco::GenJet>::const_iterator jet_it = jets->begin(); jet_it != jets->end(); ++jet_it) {
0032     ++njets;
0033 
0034     hists_["jetPtAll"]->Fill(jet_it->pt(), weight);
0035     hists_["jetEtaAll"]->Fill(jet_it->eta(), weight);
0036 
0037     if (njets == 1) {
0038       hists_["jetPt1"]->Fill(jet_it->pt(), weight);
0039       hists_["jetEta1"]->Fill(jet_it->eta(), weight);
0040     }
0041     if (njets == 2) {
0042       hists_["jetPt2"]->Fill(jet_it->pt(), weight);
0043       hists_["jetEta2"]->Fill(jet_it->eta(), weight);
0044     }
0045     if (njets == 3) {
0046       hists_["jetPt3"]->Fill(jet_it->pt(), weight);
0047       hists_["jetEta3"]->Fill(jet_it->eta(), weight);
0048     }
0049     if (njets == 4) {
0050       hists_["jetPt4"]->Fill(jet_it->pt(), weight);
0051       hists_["jetEta4"]->Fill(jet_it->eta(), weight);
0052     }
0053   }
0054 }
0055 
0056 void TTbar_GenJetAnalyzer::bookHistograms(DQMStore::IBooker& i, edm::Run const&, edm::EventSetup const&) {
0057   DQMHelper dqm(&i);
0058   i.setCurrentFolder("Generator/TTbar");
0059   hists_["jetPtAll"] =
0060       dqm.book1dHisto("TTbar_jetPtAll", "pt", 1000, 0., 1000., "P_{t}^{All-Jets} (GeV)", "Number of Events");
0061   hists_["jetPt1"] =
0062       dqm.book1dHisto("TTbar_jetPt1", "pt", 1000, 0., 1000., "P_{t}^{1st-Jet} (GeV)", "Number of Events");
0063   hists_["jetPt2"] =
0064       dqm.book1dHisto("TTbar_jetPt2", "pt", 1000, 0., 1000., "P_{t}^{2nd-Jet} (GeV)", "Number of Events");
0065   hists_["jetPt3"] =
0066       dqm.book1dHisto("TTbar_jetPt3", "pt", 1000, 0., 1000., "P_{t}^{3rd-Jet} (GeV)", "Number of Events");
0067   hists_["jetPt4"] =
0068       dqm.book1dHisto("TTbar_jetPt4", "pt", 1000, 0., 1000., "P_{t}^{4th-Jet} (GeV)", "Number of Events");
0069 
0070   hists_["jetEtaAll"] = dqm.book1dHisto("TTbar_jetEtaAll", "eta", 100, -5., 5., "#eta^{All-Jets}", "Number of Events");
0071   hists_["jetEta1"] = dqm.book1dHisto("TTbar_jetEta1", "eta", 100, -5., 5., "#eta^{1st-Jet}", "Number of Events");
0072   hists_["jetEta2"] = dqm.book1dHisto("TTbar_jetEta2", "eta", 100, -5., 5., "#eta^{2nd-Jet}", "Number of Events");
0073   hists_["jetEta3"] = dqm.book1dHisto("TTbar_jetEta3", "eta", 100, -5., 5., "#eta^{3rd-Jet}", "Number of Events");
0074   hists_["jetEta4"] = dqm.book1dHisto("TTbar_jetEta4", "eta", 100, -5., 5., "#eta^{4th-Jet}", "Number of Events");
0075 }