Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-12-13 13:14:14

0001 /* Example analizer to use L7 flavor JetCorrector services
0002    Applies different parton corrections randomly
0003    original from F.Ratnikov (UMd)  Nov 16, 2007 - Adapted deom A.Santocchia Mar 01, 2008
0004 */
0005 
0006 #include <string>
0007 
0008 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0009 #include "FWCore/Utilities/interface/InputTag.h"
0010 
0011 class PartonJetCorrectionExample : public edm::one::EDAnalyzer<> {
0012 public:
0013   explicit PartonJetCorrectionExample(const edm::ParameterSet& fParameters);
0014   ~PartonJetCorrectionExample() override {}
0015   void analyze(const edm::Event&, const edm::EventSetup&) override;
0016 
0017 private:
0018   edm::InputTag mInput;
0019   std::string m_gJ_CorrectorName;
0020   std::string m_qJ_CorrectorName;
0021   std::string m_bJ_CorrectorName;
0022   std::string m_bT_CorrectorName;
0023 };
0024 
0025 #include "FWCore/Framework/interface/Event.h"
0026 #include "DataFormats/Common/interface/Handle.h"
0027 #include "FWCore/Framework/interface/EventSetup.h"
0028 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0029 #include "DataFormats/JetReco/interface/CaloJetCollection.h"
0030 #include "JetMETCorrections/Objects/interface/JetCorrector.h"
0031 
0032 using namespace std;
0033 using namespace reco;
0034 
0035 PartonJetCorrectionExample::PartonJetCorrectionExample(const edm::ParameterSet& fConfig)
0036     : mInput(fConfig.getParameter<edm::InputTag>("src")),
0037       m_gJ_CorrectorName(fConfig.getParameter<std::string>("gJetCorrector")),
0038       m_qJ_CorrectorName(fConfig.getParameter<std::string>("qJetCorrector")),
0039       m_bJ_CorrectorName(fConfig.getParameter<std::string>("bJetCorrector")),
0040       m_bT_CorrectorName(fConfig.getParameter<std::string>("bTopCorrector")) {}
0041 
0042 void PartonJetCorrectionExample::analyze(const edm::Event& fEvent, const edm::EventSetup& fSetup) {
0043   // get all correctors
0044   const JetCorrector* gJetCorrector = JetCorrector::getJetCorrector(m_gJ_CorrectorName, fSetup);
0045   const JetCorrector* qJetCorrector = JetCorrector::getJetCorrector(m_qJ_CorrectorName, fSetup);
0046   const JetCorrector* bJetCorrector = JetCorrector::getJetCorrector(m_bJ_CorrectorName, fSetup);
0047   const JetCorrector* bTopCorrector = JetCorrector::getJetCorrector(m_bT_CorrectorName, fSetup);
0048   const JetCorrector* corrector = nullptr;
0049 
0050   // get input jets (supposed to be MC corrected already)
0051   edm::Handle<CaloJetCollection> jets;
0052   fEvent.getByLabel(mInput, jets);
0053   // loop over jets
0054   for (unsigned ijet = 0; ijet < jets->size(); ++ijet) {
0055     const CaloJet& jet = (*jets)[ijet];
0056     std::cout << "PartonJetCorrectionExample::analize-> jet #" << ijet;
0057     if (ijet % 4 == 0) {  // assume it is gluon from diJet
0058       std::cout << ": use gJ corrections" << std::endl;
0059       corrector = gJetCorrector;
0060     } else if (ijet % 4 == 1) {  // assume it is light quark from diJet
0061       std::cout << ": use qJ corrections" << std::endl;
0062       corrector = qJetCorrector;
0063     } else if (ijet % 4 == 2) {  // assume it is b quark from diJet
0064       std::cout << ": use bJ corrections" << std::endl;
0065       corrector = bJetCorrector;
0066     } else {  // assume it is b quark from ttbar
0067       std::cout << ": use bT corrections" << std::endl;
0068       corrector = bTopCorrector;
0069     }
0070     // get selected correction for the jet
0071     double correction = corrector->correction(jet);
0072     // dump it
0073     std::cout << "  jet pt/eta/phi: " << jet.pt() << '/' << jet.eta() << '/' << jet.phi()
0074               << " -> correction factor: " << correction << ", corrected pt: " << jet.pt() * correction << std::endl;
0075   }
0076 }
0077 
0078 #include "FWCore/Framework/interface/MakerMacros.h"
0079 DEFINE_FWK_MODULE(PartonJetCorrectionExample);