Back to home page

Project CMSSW displayed by LXR

 
 

    


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

0001 // Implementation of class LXXXCorrector.
0002 // Generic LX jet corrector class.
0003 
0004 #include "JetMETCorrections/Algorithms/interface/LXXXCorrector.h"
0005 #include "CondFormats/JetMETObjects/interface/FactorizedJetCorrectorCalculator.h"
0006 
0007 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0008 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0009 #include "FWCore/ParameterSet/interface/FileInPath.h"
0010 #include "DataFormats/JetReco/interface/Jet.h"
0011 #include "DataFormats/JetReco/interface/CaloJet.h"
0012 #include "DataFormats/JetReco/interface/JPTJet.h"
0013 
0014 using namespace std;
0015 
0016 //------------------------------------------------------------------------
0017 //--- LXXXCorrector constructor ------------------------------------------
0018 //------------------------------------------------------------------------
0019 LXXXCorrector::LXXXCorrector(const JetCorrectorParameters& fParam, const edm::ParameterSet& fConfig) {
0020   string level = fParam.definitions().level();
0021   if (level == "L2Relative")
0022     mLevel = 2;
0023   else if (level == "L3Absolute")
0024     mLevel = 3;
0025   else if (level == "L4EMF")
0026     mLevel = 4;
0027   else if (level == "L5Flavor")
0028     mLevel = 5;
0029   else if (level == "L7Parton")
0030     mLevel = 7;
0031   else if (level == "L2L3Residual")
0032     mLevel = 8;
0033   else
0034     throw cms::Exception("LXXXCorrector") << " unknown correction level " << level;
0035   vector<JetCorrectorParameters> vParam;
0036   vParam.push_back(fParam);
0037   mCorrector = new FactorizedJetCorrectorCalculator(vParam);
0038 }
0039 //------------------------------------------------------------------------
0040 //--- LXXXCorrector destructor -------------------------------------------
0041 //------------------------------------------------------------------------
0042 LXXXCorrector::~LXXXCorrector() { delete mCorrector; }
0043 //------------------------------------------------------------------------
0044 //--- Returns correction for a given 4-vector ----------------------------
0045 //------------------------------------------------------------------------
0046 double LXXXCorrector::correction(const LorentzVector& fJet) const {
0047   // L4 correction requires more information than a simple 4-vector
0048   if (mLevel == 4) {
0049     throw cms::Exception("Invalid jet type") << "L4EMFCorrection is applicable to CaloJets only";
0050     return 1;
0051   }
0052 
0053   FactorizedJetCorrectorCalculator::VariableValues values;
0054   values.setJetEta(fJet.eta());
0055   values.setJetE(fJet.energy());
0056   values.setJetPt(fJet.pt());
0057   values.setJetPhi(fJet.phi());
0058 
0059   return mCorrector->getCorrection(values);
0060 }
0061 //------------------------------------------------------------------------
0062 //--- Returns correction for a given jet ---------------------------------
0063 //------------------------------------------------------------------------
0064 double LXXXCorrector::correction(const reco::Jet& fJet) const {
0065   double result = 1.;
0066   // L4 correction applies to Calojets only
0067   if (mLevel == 4) {
0068     const reco::CaloJet& caloJet = dynamic_cast<const reco::CaloJet&>(fJet);
0069     FactorizedJetCorrectorCalculator::VariableValues values;
0070     values.setJetEta(fJet.eta());
0071     values.setJetPt(fJet.pt());
0072     values.setJetEMF(caloJet.emEnergyFraction());
0073     result = mCorrector->getCorrection(values);
0074   } else
0075     result = correction(fJet.p4());
0076   return result;
0077 }