Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 14:21:38

0001 //
0002 // Original Author:  Fedor Ratnikov Feb. 16, 2007
0003 //
0004 // Correction which chains other corrections
0005 //
0006 
0007 #include "JetMETCorrections/Objects/interface/ChainedJetCorrector.h"
0008 
0009 /// get correction using Jet information only
0010 double ChainedJetCorrector::correction(const LorentzVector& fJet) const {
0011   LorentzVector jet = fJet;
0012   double result = 1;
0013   for (size_t i = 0; i < mCorrectors.size(); ++i) {
0014     double scale = mCorrectors[i]->correction(jet);
0015     jet *= scale;
0016     result *= scale;
0017   }
0018   return result;
0019 }
0020 
0021 /// apply correction using Jet information only
0022 double ChainedJetCorrector::correction(const reco::Jet& fJet) const {
0023   std::unique_ptr<reco::Jet> jet(dynamic_cast<reco::Jet*>(fJet.clone()));
0024   double result = 1;
0025   for (size_t i = 0; i < mCorrectors.size(); ++i) {
0026     double scale = mCorrectors[i]->correction(*jet);
0027     jet->scaleEnergy(scale);
0028     result *= scale;
0029   }
0030   return result;
0031 }
0032 
0033 /// apply correction using all event information
0034 double ChainedJetCorrector::correction(const reco::Jet& fJet,
0035                                        const edm::Event& fEvent,
0036                                        const edm::EventSetup& fSetup) const {
0037   std::unique_ptr<reco::Jet> jet(dynamic_cast<reco::Jet*>(fJet.clone()));
0038   double result = 1;
0039   for (size_t i = 0; i < mCorrectors.size(); ++i) {
0040     double scale = mCorrectors[i]->correction(*jet, fEvent, fSetup);
0041     jet->scaleEnergy(scale);
0042     result *= scale;
0043   }
0044   return result;
0045 }
0046 /// apply correction using all event information and reference to the raw jet
0047 double ChainedJetCorrector::correction(const reco::Jet& fJet,
0048                                        const edm::RefToBase<reco::Jet>& fJetRef,
0049                                        const edm::Event& fEvent,
0050                                        const edm::EventSetup& fSetup) const {
0051   std::unique_ptr<reco::Jet> jet(dynamic_cast<reco::Jet*>(fJet.clone()));
0052   double result = 1;
0053   for (size_t i = 0; i < mCorrectors.size(); ++i) {
0054     double scale = mCorrectors[i]->correction(*jet, fJetRef, fEvent, fSetup);
0055     jet->scaleEnergy(scale);
0056     result *= scale;
0057   }
0058   return result;
0059 }
0060 
0061 /// if correction needs event information
0062 bool ChainedJetCorrector::eventRequired() const {
0063   for (size_t i = 0; i < mCorrectors.size(); ++i) {
0064     if (mCorrectors[i]->eventRequired())
0065       return true;
0066   }
0067   return false;
0068 }
0069 
0070 /// if correction needs jet reference
0071 bool ChainedJetCorrector::refRequired() const {
0072   for (size_t i = 0; i < mCorrectors.size(); ++i) {
0073     if (mCorrectors[i]->refRequired())
0074       return true;
0075   }
0076   return false;
0077 }