File indexing completed on 2024-04-06 12:19:23
0001 #include <iostream>
0002 #include <sstream>
0003 #include <istream>
0004 #include <fstream>
0005 #include <iomanip>
0006 #include <string>
0007 #include <cmath>
0008 #include <functional>
0009
0010 #include "JetMETCorrections/MCJet/plugins/CaloMCTruthTreeProducer.h"
0011 #include "JetMETCorrections/MCJet/plugins/JetUtilMC.h"
0012 #include "FWCore/Framework/interface/Event.h"
0013 #include "DataFormats/Common/interface/Handle.h"
0014 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0015 #include "FWCore/Framework/interface/EventSetup.h"
0016 #include "FWCore/Framework/interface/ESHandle.h"
0017 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0018 #include "DataFormats/JetReco/interface/CaloJet.h"
0019 #include "DataFormats/JetReco/interface/GenJet.h"
0020 using namespace edm;
0021 using namespace reco;
0022 using namespace std;
0023
0024
0025
0026 CaloMCTruthTreeProducer::CaloMCTruthTreeProducer(edm::ParameterSet const& cfg) {
0027 jets_ = consumes<CaloJetCollection>(edm::InputTag(cfg.getParameter<std::string>("jets")));
0028 genjets_ = consumes<GenJetCollection>(edm::InputTag(cfg.getParameter<std::string>("genjets")));
0029 gen_ = consumes<GenEventInfoProduct>(edm::InputTag("generator"));
0030 histogramFile_ = cfg.getParameter<std::string>("histogramFile");
0031 }
0032
0033 void CaloMCTruthTreeProducer::beginJob() {
0034 file_ = new TFile(histogramFile_.c_str(), "RECREATE");
0035 mcTruthTree_ = new TTree("mcTruthTree", "mcTruthTree");
0036
0037 mcTruthTree_->Branch("ptJet", &ptJet_, "ptJet_/F");
0038 mcTruthTree_->Branch("ptGen", &ptGen_, "ptGen_/F");
0039 mcTruthTree_->Branch("ptHat", &ptHat_, "ptHat_/F");
0040 mcTruthTree_->Branch("emfJet", &emfJet_, "emfJet_/F");
0041 mcTruthTree_->Branch("etaJet", &etaJet_, "etaJet_/F");
0042 mcTruthTree_->Branch("etaGen", &etaGen_, "etaGen_/F");
0043 mcTruthTree_->Branch("phiJet", &phiJet_, "phiJet_/F");
0044 mcTruthTree_->Branch("phiGen", &phiGen_, "phiGen_/F");
0045 mcTruthTree_->Branch("dR", &dR_, "dR_/F");
0046 mcTruthTree_->Branch("rank", &rank_, "rank_/I");
0047 }
0048
0049 void CaloMCTruthTreeProducer::endJob() {
0050 if (file_ != nullptr) {
0051 file_->cd();
0052 mcTruthTree_->Write();
0053 }
0054 file_ = nullptr;
0055 }
0056
0057 void CaloMCTruthTreeProducer::analyze(edm::Event const& event, edm::EventSetup const& iSetup) {
0058 edm::Handle<GenJetCollection> genjets;
0059 edm::Handle<CaloJetCollection> jets;
0060 edm::Handle<GenEventInfoProduct> hEventInfo;
0061 CaloJetCollection::const_iterator i_jet, i_matched;
0062 GenJetCollection::const_iterator i_genjet;
0063 event.getByToken(genjets_, genjets);
0064 event.getByToken(jets_, jets);
0065 event.getByToken(gen_, hEventInfo);
0066 ptHat_ = hEventInfo->binningValues()[0];
0067 float rr;
0068 int njet(0);
0069 if (!jets->empty() && !genjets->empty()) {
0070 for (i_genjet = genjets->begin(); i_genjet != genjets->end(); i_genjet++) {
0071 float rmin(99);
0072 for (i_jet = jets->begin(); i_jet != jets->end(); i_jet++) {
0073 rr = radius(i_genjet, i_jet);
0074 if (rr < rmin) {
0075 rmin = rr;
0076 i_matched = i_jet;
0077 }
0078 }
0079 ptGen_ = i_genjet->pt();
0080 etaGen_ = i_genjet->eta();
0081 phiGen_ = i_genjet->phi();
0082 ptJet_ = i_matched->pt();
0083 etaJet_ = i_matched->eta();
0084 phiJet_ = i_matched->phi();
0085 emfJet_ = i_matched->emEnergyFraction();
0086 dR_ = rmin;
0087 rank_ = njet;
0088 mcTruthTree_->Fill();
0089 njet++;
0090 }
0091 }
0092 }
0093
0094 CaloMCTruthTreeProducer::~CaloMCTruthTreeProducer() {
0095 delete file_;
0096 delete mcTruthTree_;
0097 }
0098