File indexing completed on 2024-04-06 12:02:18
0001 #include "CondFormats/JetMETObjects/interface/JetCorrectionUncertainty.h"
0002 #include "CondFormats/JetMETObjects/interface/SimpleJetCorrectionUncertainty.h"
0003 #include "CondFormats/JetMETObjects/interface/JetCorrectorParameters.h"
0004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0005 #include "Math/PtEtaPhiE4D.h"
0006 #include "Math/Vector3D.h"
0007 #include "Math/LorentzVector.h"
0008 #include <vector>
0009 #include <string>
0010
0011
0012 JetCorrectionUncertainty::JetCorrectionUncertainty() {
0013 mJetEta = -9999;
0014 mJetPt = -9999;
0015 mJetPhi = -9999;
0016 mJetE = -9999;
0017 mJetEMF = -9999;
0018 mLepPx = -9999;
0019 mLepPy = -9999;
0020 mLepPz = -9999;
0021 mIsJetEset = false;
0022 mIsJetPtset = false;
0023 mIsJetPhiset = false;
0024 mIsJetEtaset = false;
0025 mIsJetEMFset = false;
0026 mIsLepPxset = false;
0027 mIsLepPyset = false;
0028 mIsLepPzset = false;
0029 mAddLepToJet = false;
0030 mUncertainty = new SimpleJetCorrectionUncertainty();
0031 }
0032
0033 JetCorrectionUncertainty::JetCorrectionUncertainty(const std::string& fDataFile) {
0034 mJetEta = -9999;
0035 mJetPt = -9999;
0036 mJetPhi = -9999;
0037 mJetE = -9999;
0038 mJetEMF = -9999;
0039 mLepPx = -9999;
0040 mLepPy = -9999;
0041 mLepPz = -9999;
0042 mIsJetEset = false;
0043 mIsJetPtset = false;
0044 mIsJetPhiset = false;
0045 mIsJetEtaset = false;
0046 mIsJetEMFset = false;
0047 mIsLepPxset = false;
0048 mIsLepPyset = false;
0049 mIsLepPzset = false;
0050 mAddLepToJet = false;
0051 mUncertainty = new SimpleJetCorrectionUncertainty(fDataFile);
0052 }
0053
0054 JetCorrectionUncertainty::JetCorrectionUncertainty(const JetCorrectorParameters& fParameters) {
0055 mJetEta = -9999;
0056 mJetPt = -9999;
0057 mJetPhi = -9999;
0058 mJetE = -9999;
0059 mJetEMF = -9999;
0060 mLepPx = -9999;
0061 mLepPy = -9999;
0062 mLepPz = -9999;
0063 mIsJetEset = false;
0064 mIsJetPtset = false;
0065 mIsJetPhiset = false;
0066 mIsJetEtaset = false;
0067 mIsJetEMFset = false;
0068 mIsLepPxset = false;
0069 mIsLepPyset = false;
0070 mIsLepPzset = false;
0071 mAddLepToJet = false;
0072 mUncertainty = new SimpleJetCorrectionUncertainty(fParameters);
0073 }
0074
0075 JetCorrectionUncertainty::~JetCorrectionUncertainty() { delete mUncertainty; }
0076
0077 void JetCorrectionUncertainty::setParameters(const std::string& fDataFile) {
0078
0079 delete mUncertainty;
0080 mUncertainty = new SimpleJetCorrectionUncertainty(fDataFile);
0081 }
0082
0083 float JetCorrectionUncertainty::getUncertainty(bool fDirection) {
0084 float result;
0085 std::vector<float> vx, vy;
0086 vx = fillVector(mUncertainty->parameters().definitions().binVar());
0087 vy = fillVector(mUncertainty->parameters().definitions().parVar());
0088 result = mUncertainty->uncertainty(vx, vy[0], fDirection);
0089 mIsJetEset = false;
0090 mIsJetPtset = false;
0091 mIsJetPhiset = false;
0092 mIsJetEtaset = false;
0093 mIsJetEMFset = false;
0094 mIsLepPxset = false;
0095 mIsLepPyset = false;
0096 mIsLepPzset = false;
0097 return result;
0098 }
0099
0100
0101
0102 std::vector<float> JetCorrectionUncertainty::fillVector(const std::vector<std::string>& fNames) {
0103 std::vector<float> result;
0104 for (unsigned i = 0; i < fNames.size(); i++) {
0105 if (fNames[i] == "JetEta") {
0106 if (!mIsJetEtaset) {
0107 edm::LogError("JetCorrectionUncertainty::") << " jet eta is not set";
0108 result.push_back(-999.0);
0109 } else {
0110 result.push_back(mJetEta);
0111 }
0112 } else if (fNames[i] == "JetPt") {
0113 if (!mIsJetPtset) {
0114 edm::LogError("JetCorrectionUncertainty::") << " jet pt is not set";
0115 result.push_back(-999.0);
0116 } else {
0117 result.push_back(mJetPt);
0118 }
0119 } else if (fNames[i] == "JetPhi") {
0120 if (!mIsJetPhiset) {
0121 edm::LogError("JetCorrectionUncertainty::") << " jet phi is not set";
0122 result.push_back(-999.0);
0123 } else {
0124 result.push_back(mJetPhi);
0125 }
0126 } else if (fNames[i] == "JetE") {
0127 if (!mIsJetEset) {
0128 edm::LogError("JetCorrectionUncertainty::") << " jet energy is not set";
0129 result.push_back(-999.0);
0130 } else {
0131 result.push_back(mJetE);
0132 }
0133 } else if (fNames[i] == "JetEMF") {
0134 if (!mIsJetEMFset) {
0135 edm::LogError("JetCorrectionUncertainty::") << " jet emf is not set";
0136 result.push_back(-999.0);
0137 } else {
0138 result.push_back(mJetEMF);
0139 }
0140 } else if (fNames[i] == "LepPx") {
0141 if (!mIsLepPxset) {
0142 edm::LogError("JetCorrectionUncertainty::") << " lepton px is not set";
0143 result.push_back(-999.0);
0144 } else {
0145 result.push_back(mLepPx);
0146 }
0147 } else if (fNames[i] == "LepPy") {
0148 if (!mIsLepPyset) {
0149 edm::LogError("JetCorrectionUncertainty::") << " lepton py is not set";
0150 result.push_back(-999.0);
0151 } else {
0152 result.push_back(mLepPy);
0153 }
0154 } else if (fNames[i] == "LepPz") {
0155 if (!mIsLepPzset) {
0156 edm::LogError("JetCorrectionUncertainty::") << " lepton pz is not set";
0157 result.push_back(-999.0);
0158 } else {
0159 result.push_back(mLepPz);
0160 }
0161 }
0162
0163 else {
0164 edm::LogError("JetCorrectionUncertainty::") << " unknown parameter " << fNames[i];
0165 result.push_back(-999.0);
0166 }
0167 }
0168 return result;
0169 }
0170
0171
0172
0173 float JetCorrectionUncertainty::getPtRel() {
0174 typedef ROOT::Math::LorentzVector<ROOT::Math::PtEtaPhiE4D<float> > PtEtaPhiELorentzVector;
0175 typedef ROOT::Math::DisplacementVector3D<ROOT::Math::Cartesian3D<float> > XYZVector;
0176 PtEtaPhiELorentzVector jet;
0177 XYZVector lep;
0178 jet.SetPt(mJetPt);
0179 jet.SetEta(mJetEta);
0180 jet.SetPhi(mJetPhi);
0181 jet.SetE(mJetE);
0182 lep.SetXYZ(mLepPx, mLepPy, mLepPz);
0183 float lj_x = (mAddLepToJet) ? lep.X() + jet.Px() : jet.Px();
0184 float lj_y = (mAddLepToJet) ? lep.Y() + jet.Py() : jet.Py();
0185 float lj_z = (mAddLepToJet) ? lep.Z() + jet.Pz() : jet.Pz();
0186
0187 float lj2 = lj_x * lj_x + lj_y * lj_y + lj_z * lj_z;
0188 float pTrel2 = -999.0;
0189 if (lj2 > 0) {
0190 float lep2 = lep.X() * lep.X() + lep.Y() * lep.Y() + lep.Z() * lep.Z();
0191
0192 float lepXlj = lep.X() * lj_x + lep.Y() * lj_y + lep.Z() * lj_z;
0193
0194 float pLrel2 = lepXlj * lepXlj / lj2;
0195
0196 pTrel2 = lep2 - pLrel2;
0197 } else
0198 edm::LogError("JetCorrectionUncertainty") << " not positive lepton-jet momentum: " << lj2;
0199 return (pTrel2 > 0) ? std::sqrt(pTrel2) : 0.0;
0200 }
0201
0202
0203
0204 void JetCorrectionUncertainty::setJetEta(float fEta) {
0205 mJetEta = fEta;
0206 mIsJetEtaset = true;
0207 }
0208
0209 void JetCorrectionUncertainty::setJetPt(float fPt) {
0210 mJetPt = fPt;
0211 mIsJetPtset = true;
0212 }
0213
0214 void JetCorrectionUncertainty::setJetPhi(float fPhi) {
0215 mJetPhi = fPhi;
0216 mIsJetPhiset = true;
0217 }
0218
0219 void JetCorrectionUncertainty::setJetE(float fE) {
0220 mJetE = fE;
0221 mIsJetEset = true;
0222 }
0223
0224 void JetCorrectionUncertainty::setJetEMF(float fEMF) {
0225 mJetEMF = fEMF;
0226 mIsJetEMFset = true;
0227 }
0228
0229 void JetCorrectionUncertainty::setLepPx(float fPx) {
0230 mLepPx = fPx;
0231 mIsLepPxset = true;
0232 }
0233
0234 void JetCorrectionUncertainty::setLepPy(float fPy) {
0235 mLepPy = fPy;
0236 mIsLepPyset = true;
0237 }
0238
0239 void JetCorrectionUncertainty::setLepPz(float fPz) {
0240 mLepPz = fPz;
0241 mIsLepPzset = true;
0242 }
0243