File indexing completed on 2024-09-07 04:37:23
0001 #include <map>
0002 #include <string>
0003 #include <iomanip>
0004 #include <sstream>
0005 #include <iostream>
0006
0007 #include "TH1.h"
0008
0009 #include "FWCore/Framework/interface/Event.h"
0010 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0011 #include "FWCore/Utilities/interface/InputTag.h"
0012 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0013 #include "FWCore/ServiceRegistry/interface/Service.h"
0014 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0015
0016 #include "DataFormats/PatCandidates/interface/Jet.h"
0017 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0018
0019
0020
0021 static const unsigned int MAXBIN = 8;
0022
0023
0024
0025 static const float BINS[] = {30., 40., 50., 60., 70., 80., 100., 125., 150.};
0026
0027
0028
0029
0030
0031
0032
0033
0034
0035
0036
0037
0038
0039
0040
0041
0042 class PatJetAnalyzer : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0043 public:
0044
0045 explicit PatJetAnalyzer(const edm::ParameterSet& cfg);
0046
0047 ~PatJetAnalyzer() override {}
0048
0049 private:
0050
0051 void analyze(const edm::Event& event, const edm::EventSetup& setup) override;
0052
0053
0054 bool booked(const std::string histName) const { return hists_.find(histName) != hists_.end(); };
0055
0056 void fill(const std::string histName, double value) const {
0057 if (booked(histName))
0058 hists_.find(histName)->second->Fill(value);
0059 };
0060
0061 void print(edm::View<pat::Jet>::const_iterator& jet, unsigned int idx);
0062
0063 private:
0064
0065 std::string corrLevel_;
0066
0067 edm::EDGetTokenT<edm::View<pat::Jet> > jetsToken_;
0068
0069 std::map<std::string, TH1F*> hists_;
0070 };
0071
0072 PatJetAnalyzer::PatJetAnalyzer(const edm::ParameterSet& cfg)
0073 : corrLevel_(cfg.getParameter<std::string>("corrLevel")),
0074 jetsToken_(consumes<edm::View<pat::Jet> >(cfg.getParameter<edm::InputTag>("src"))) {
0075 usesResource(TFileService::kSharedResource);
0076
0077 edm::Service<TFileService> fs;
0078
0079
0080 hists_["mult"] = fs->make<TH1F>("mult", "N_{Jet}", 15, 0., 15.);
0081
0082 hists_["pt"] = fs->make<TH1F>("pt", "p_{T}(Jet) [GeV]", 60, 0., 300.);
0083
0084 hists_["eta"] = fs->make<TH1F>("eta", "#eta (Jet)", 60, -3., 3.);
0085
0086 hists_["phi"] = fs->make<TH1F>("phi", "#phi (Jet)", 60, 3.2, 3.2);
0087
0088 hists_["mass"] = fs->make<TH1F>("mass", "M_{jj} [GeV]", 50, 0., 500.);
0089
0090 for (unsigned int idx = 0; idx < MAXBIN; ++idx) {
0091 char buffer[10];
0092 sprintf(buffer, "jes_%i", idx);
0093 char title[50];
0094 sprintf(title, "p_{T}^{rec}/p_{T}^{gen} [%i GeV - %i GeV]", (int)BINS[idx], (int)BINS[idx + 1]);
0095 hists_[buffer] = fs->make<TH1F>(buffer, title, 100, 0., 2.);
0096 }
0097 }
0098
0099 void PatJetAnalyzer::analyze(const edm::Event& event, const edm::EventSetup& setup) {
0100
0101 edm::Handle<edm::View<pat::Jet> > jets;
0102 event.getByToken(jetsToken_, jets);
0103
0104
0105 for (edm::View<pat::Jet>::const_iterator jet = jets->begin(); jet != jets->end(); ++jet) {
0106
0107
0108
0109
0110 fill("pt", jet->correctedJet(corrLevel_).pt());
0111 fill("eta", jet->eta());
0112 fill("phi", jet->phi());
0113
0114 if (jet->genJet()) {
0115 double resp = jet->correctedJet(corrLevel_).pt() / jet->genJet()->pt();
0116 for (unsigned int idx = 0; idx < MAXBIN; ++idx) {
0117 if (BINS[idx] <= jet->genJet()->pt() && jet->genJet()->pt() < BINS[idx + 1]) {
0118 char buffer[10];
0119 sprintf(buffer, "jes_%i", idx);
0120 fill(buffer, resp);
0121 }
0122 }
0123 }
0124 }
0125
0126 fill("mult", jets->size());
0127
0128 if (jets->size() > 1) {
0129 fill("mass", ((*jets)[0].p4() + (*jets)[1].p4()).mass());
0130 }
0131 }
0132
0133 void PatJetAnalyzer::print(edm::View<pat::Jet>::const_iterator& jet, unsigned int idx) {
0134
0135 std::cout << "[" << idx << "] :: eta=" << std::setw(10) << jet->eta() << " phi=" << std::setw(10) << jet->phi()
0136 << " size: " << jet->availableJECLevels().size() << std::endl;
0137 for (unsigned int idx = 0; idx < jet->availableJECLevels().size(); ++idx) {
0138 pat::Jet correctedJet;
0139 if (jet->availableJECLevels()[idx].find("L5Flavor") != std::string::npos ||
0140 jet->availableJECLevels()[idx].find("L7Parton") != std::string::npos) {
0141 correctedJet = jet->correctedJet(idx, pat::JetCorrFactors::UDS);
0142 } else {
0143 correctedJet = jet->correctedJet(idx, pat::JetCorrFactors::NONE);
0144 }
0145 std::cout << std::setw(10) << correctedJet.currentJECLevel() << " pt=" << std::setw(10) << correctedJet.pt()
0146 << std::endl;
0147 }
0148 }
0149
0150 #include "FWCore/Framework/interface/MakerMacros.h"
0151 DEFINE_FWK_MODULE(PatJetAnalyzer);