Back to home page

Project CMSSW displayed by LXR

 
 

    


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

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