File indexing completed on 2024-09-07 04:37:18
0001
0002 #include "FWCore/Framework/interface/Frameworkfwd.h"
0003 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0004 #include "FWCore/Framework/interface/ESHandle.h"
0005 #include "FWCore/Framework/interface/Event.h"
0006 #include "FWCore/Framework/interface/EventSetup.h"
0007 #include "FWCore/Framework/interface/MakerMacros.h"
0008 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0009
0010 #include "DataFormats/Candidate/interface/Candidate.h"
0011 #include "DataFormats/Candidate/interface/CandidateFwd.h"
0012 #include "DataFormats/JetReco/interface/Jet.h"
0013 #include "DataFormats/JetReco/interface/JetCollection.h"
0014 #include "DataFormats/JetReco/interface/CaloJetCollection.h"
0015 #include "DataFormats/JetReco/interface/JetFloatAssociation.h"
0016 #include "DataFormats/Common/interface/Ref.h"
0017
0018 #include "DataFormats/JetMatching/interface/JetFlavour.h"
0019 #include "DataFormats/JetMatching/interface/JetFlavourMatching.h"
0020
0021 class printGenJetRatio : public edm::one::EDAnalyzer<> {
0022 public:
0023 typedef reco::JetFloatAssociation::Container JetBCEnergyRatioCollection;
0024
0025 explicit printGenJetRatio(const edm::ParameterSet&);
0026 ~printGenJetRatio() {}
0027 void analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup);
0028
0029 private:
0030 edm::EDGetTokenT<JetBCEnergyRatioCollection> sourceBratioToken_;
0031 edm::EDGetTokenT<JetBCEnergyRatioCollection> sourceCratioToken_;
0032 edm::Handle<JetBCEnergyRatioCollection> theBratioValue;
0033 edm::Handle<JetBCEnergyRatioCollection> theCratioValue;
0034 };
0035
0036
0037 #include <memory>
0038 #include <string>
0039 #include <iostream>
0040 #include <vector>
0041
0042 using namespace std;
0043 using namespace reco;
0044 using namespace edm;
0045
0046 printGenJetRatio::printGenJetRatio(const edm::ParameterSet& iConfig) {
0047 sourceBratioToken_ = consumes<JetBCEnergyRatioCollection>(iConfig.getParameter<InputTag>("srcBratio"));
0048 sourceCratioToken_ = consumes<JetBCEnergyRatioCollection>(iConfig.getParameter<InputTag>("srcCratio"));
0049 }
0050
0051 void printGenJetRatio::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0052 cout << "[printGenJetRatio] analysing event " << iEvent.id() << endl;
0053
0054 try {
0055 iEvent.getByToken(sourceBratioToken_, theBratioValue);
0056 iEvent.getByToken(sourceCratioToken_, theCratioValue);
0057 } catch (std::exception& ce) {
0058 cerr << "[printJetFlavour] caught std::exception " << ce.what() << endl;
0059 return;
0060 }
0061
0062 cout << "-------------------- GenJet Bratio ------------------------" << endl;
0063 for (JetBCEnergyRatioCollection::const_iterator itB = theBratioValue->begin(); itB != theBratioValue->end(); itB++) {
0064 const Jet& jetB = *(itB->first);
0065 float cR = 0;
0066 for (JetBCEnergyRatioCollection::const_iterator itC = theCratioValue->begin(); itC != theCratioValue->end();
0067 itC++) {
0068 if (itB->first == itC->first)
0069 cR = itC->second;
0070 }
0071 printf("printGenJetRatio] (pt,eta,phi) jet = %7.3f %6.3f %6.3f | bcRatio = %7.5f - %7.5f \n",
0072 jetB.et(),
0073 jetB.eta(),
0074 jetB.phi(),
0075 itB->second,
0076 cR);
0077 }
0078 }
0079
0080 DEFINE_FWK_MODULE(printGenJetRatio);