Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:02:18

0001 //
0002 // Original Author:  Fedor Ratnikov Nov 9, 2007
0003 // $Id: JetCorrectorParameters.cc,v 1.20 2012/03/01 18:24:53 srappocc Exp $
0004 //
0005 // Generic parameters for Jet corrections
0006 //
0007 #include "CondFormats/JetMETObjects/interface/JetCorrectorParameters.h"
0008 #include "CondFormats/JetMETObjects/interface/JetCorrectorParametersHelper.h"
0009 #include <iostream>
0010 #include <iomanip>
0011 #include <fstream>
0012 #include <sstream>
0013 #include <cstdlib>
0014 #include <algorithm>
0015 #include <cmath>
0016 #include <iterator>
0017 
0018 //------------------------------------------------------------------------
0019 //--- JetCorrectorParameters::Definitions constructor --------------------
0020 //--- takes specific arguments for the member variables ------------------
0021 //------------------------------------------------------------------------
0022 JetCorrectorParameters::Definitions::Definitions(const std::vector<std::string>& fBinVar,
0023                                                  const std::vector<std::string>& fParVar,
0024                                                  const std::string& fFormula,
0025                                                  bool fIsResponse) {
0026   for (unsigned i = 0; i < fBinVar.size(); i++)
0027     mBinVar.push_back(fBinVar[i]);
0028   for (unsigned i = 0; i < fParVar.size(); i++)
0029     mParVar.push_back(fParVar[i]);
0030   mFormula = fFormula;
0031   mIsResponse = fIsResponse;
0032   mLevel = "";
0033 }
0034 //------------------------------------------------------------------------
0035 //--- JetCorrectorParameters::Definitions constructor --------------------
0036 //--- reads the member variables from a string ---------------------------
0037 //------------------------------------------------------------------------
0038 JetCorrectorParameters::Definitions::Definitions(const std::string& fLine) {
0039   std::vector<std::string> tokens = getTokens(fLine);
0040   if (!tokens.empty()) {
0041     if (tokens.size() < 6) {
0042       std::stringstream sserr;
0043       sserr << "(line " << fLine << "): less than 6 expected tokens:" << tokens.size();
0044       handleError("JetCorrectorParameters::Definitions", sserr.str());
0045     }
0046     unsigned nvar = getUnsigned(tokens[0]);
0047     unsigned npar = getUnsigned(tokens[nvar + 1]);
0048     for (unsigned i = 0; i < nvar; i++)
0049       mBinVar.push_back(tokens[i + 1]);
0050     for (unsigned i = 0; i < npar; i++)
0051       mParVar.push_back(tokens[nvar + 2 + i]);
0052     mFormula = tokens[npar + nvar + 2];
0053     std::string ss = tokens[npar + nvar + 3];
0054     if (ss == "Response")
0055       mIsResponse = true;
0056     else if (ss == "Correction")
0057       mIsResponse = false;
0058     else if (ss == "Resolution")
0059       mIsResponse = false;
0060     else if (ss.find("PAR") == 0)
0061       mIsResponse = false;
0062     else {
0063       std::stringstream sserr;
0064       sserr << "unknown option (" << ss << ")";
0065       handleError("JetCorrectorParameters::Definitions", sserr.str());
0066     }
0067     mLevel = tokens[npar + nvar + 4];
0068   }
0069 }
0070 //------------------------------------------------------------------------
0071 //--- JetCorrectorParameters::Record constructor -------------------------
0072 //--- reads the member variables from a string ---------------------------
0073 //------------------------------------------------------------------------
0074 JetCorrectorParameters::Record::Record(const std::string& fLine, unsigned fNvar) : mMin(0), mMax(0) {
0075   mNvar = fNvar;
0076   // quckly parse the line
0077   std::vector<std::string> tokens = getTokens(fLine);
0078   if (!tokens.empty()) {
0079     if (tokens.size() < 3) {
0080       std::stringstream sserr;
0081       sserr << "(line " << fLine << "): "
0082             << "three tokens expected, " << tokens.size() << " provided.";
0083       handleError("JetCorrectorParameters::Record", sserr.str());
0084     }
0085     for (unsigned i = 0; i < mNvar; i++) {
0086       mMin.push_back(getFloat(tokens[i * 2]));
0087       mMax.push_back(getFloat(tokens[i * 2 + 1]));
0088     }
0089     unsigned nParam = getUnsigned(tokens[2 * mNvar]);
0090     if (nParam != tokens.size() - (2 * mNvar + 1)) {
0091       std::stringstream sserr;
0092       sserr << "(line " << fLine << "): " << tokens.size() - (2 * mNvar + 1) << " parameters, but nParam=" << nParam
0093             << ".";
0094       handleError("JetCorrectorParameters::Record", sserr.str());
0095     }
0096     for (unsigned i = (2 * mNvar + 1); i < tokens.size(); ++i)
0097       mParameters.push_back(getFloat(tokens[i]));
0098   }
0099 }
0100 std::ostream& operator<<(std::ostream& out, const JetCorrectorParameters::Record& fBin) {
0101   for (unsigned j = 0; j < fBin.nVar(); j++)
0102     out << fBin.xMin(j) << " " << fBin.xMax(j) << " ";
0103   out << fBin.nParameters() << " ";
0104   for (unsigned j = 0; j < fBin.nParameters(); j++)
0105     out << fBin.parameter(j) << " ";
0106   out << std::endl;
0107   return out;
0108 }
0109 //------------------------------------------------------------------------
0110 //--- JetCorrectorParameters constructor ---------------------------------
0111 //--- reads the member variables from a string ---------------------------
0112 //------------------------------------------------------------------------
0113 JetCorrectorParameters::JetCorrectorParameters(const std::string& fFile, const std::string& fSection) {
0114   std::ifstream input(fFile.c_str());
0115   std::string currentSection = "";
0116   std::string line;
0117   std::string currentDefinitions = "";
0118   while (std::getline(input, line)) {
0119     if (line.empty())
0120       continue;
0121     std::string section = getSection(line);
0122     std::string tmp = getDefinitions(line);
0123     if (!section.empty() && tmp.empty()) {
0124       currentSection = section;
0125       continue;
0126     }
0127     if (currentSection == fSection) {
0128       if (!tmp.empty()) {
0129         currentDefinitions = tmp;
0130         continue;
0131       }
0132       Definitions definitions(currentDefinitions);
0133       if (!(definitions.nBinVar() == 0 && definitions.formula().empty()))
0134         mDefinitions = definitions;
0135       Record record(line, mDefinitions.nBinVar());
0136       bool check(true);
0137       for (unsigned i = 0; i < mDefinitions.nBinVar(); ++i)
0138         if (record.xMin(i) == 0 && record.xMax(i) == 0)
0139           check = false;
0140       if (record.nParameters() == 0)
0141         check = false;
0142       if (check)
0143         mRecords.push_back(record);
0144     }
0145   }
0146   if (currentDefinitions.empty())
0147     handleError("JetCorrectorParameters", "No definitions found!!!");
0148   if (mRecords.empty() && currentSection.empty())
0149     mRecords.push_back(Record());
0150   if (mRecords.empty() && !currentSection.empty()) {
0151     std::stringstream sserr;
0152     sserr << "the requested section " << fSection << " doesn't exist!";
0153     handleError("JetCorrectorParameters", sserr.str());
0154   }
0155   std::sort(mRecords.begin(), mRecords.end());
0156   valid_ = true;
0157 
0158   if (mDefinitions.nBinVar() <= MAX_SIZE_DIMENSIONALITY) {
0159     init();
0160   } else {
0161     std::stringstream sserr;
0162     sserr << "since the binned dimensionality is greater than " << MAX_SIZE_DIMENSIONALITY
0163           << " the SimpleJetCorrector will default to using the legacy binIndex function!";
0164     handleError("JetCorrectorParameters", sserr.str());
0165   }
0166 }
0167 //------------------------------------------------------------------------
0168 //--- initializes the correct JetCorrectorParametersHelper ---------------
0169 //------------------------------------------------------------------------
0170 void JetCorrectorParameters::init() {
0171   std::sort(mRecords.begin(), mRecords.end());
0172   helper = std::make_shared<JetCorrectorParametersHelper>();
0173   helper->init(mDefinitions, mRecords);
0174 }
0175 //------------------------------------------------------------------------
0176 //--- returns the index of the record defined by fX ----------------------
0177 //------------------------------------------------------------------------
0178 int JetCorrectorParameters::binIndex(const std::vector<float>& fX) const {
0179   int result = -1;
0180   unsigned N = mDefinitions.nBinVar();
0181   if (N != fX.size()) {
0182     std::stringstream sserr;
0183     sserr << "# bin variables " << N << " doesn't correspont to requested #: " << fX.size();
0184     handleError("JetCorrectorParameters", sserr.str());
0185   }
0186   unsigned tmp;
0187   for (unsigned i = 0; i < size(); ++i) {
0188     tmp = 0;
0189     for (unsigned j = 0; j < N; j++)
0190       if (fX[j] >= record(i).xMin(j) && fX[j] < record(i).xMax(j))
0191         tmp += 1;
0192     if (tmp == N) {
0193       result = i;
0194       break;
0195     }
0196   }
0197   return result;
0198 }
0199 //------------------------------------------------------------------------
0200 //--- returns the index of the record defined by fX (non-linear search) --
0201 //------------------------------------------------------------------------
0202 int JetCorrectorParameters::binIndexN(const std::vector<float>& fX) const {
0203   if (helper->size() > 0) {
0204     return helper->binIndexN(fX, mRecords);
0205   } else {
0206     return -1;
0207   }
0208 }
0209 //------------------------------------------------------------------------
0210 //--- returns the neighbouring bins of fIndex in the direction of fVar ---
0211 //------------------------------------------------------------------------
0212 int JetCorrectorParameters::neighbourBin(unsigned fIndex, unsigned fVar, bool fNext) const {
0213   int result = -1;
0214   unsigned N = mDefinitions.nBinVar();
0215   if (fVar >= N) {
0216     std::stringstream sserr;
0217     sserr << "# of bin variables " << N << " doesn't correspond to requested #: " << fVar;
0218     handleError("JetCorrectorParameters", sserr.str());
0219   }
0220   unsigned tmp;
0221   for (unsigned i = 0; i < size(); ++i) {
0222     tmp = 0;
0223     for (unsigned j = 0; j < fVar; j++)
0224       if (fabs(record(i).xMin(j) - record(fIndex).xMin(j)) < 0.0001)
0225         tmp += 1;
0226     for (unsigned j = fVar + 1; j < N; j++)
0227       if (fabs(record(i).xMin(j) - record(fIndex).xMin(j)) < 0.0001)
0228         tmp += 1;
0229     if (tmp < N - 1)
0230       continue;
0231     if (tmp == N - 1) {
0232       if (fNext)
0233         if (fabs(record(i).xMin(fVar) - record(fIndex).xMax(fVar)) < 0.0001)
0234           tmp += 1;
0235       if (!fNext)
0236         if (fabs(record(i).xMax(fVar) - record(fIndex).xMin(fVar)) < 0.0001)
0237           tmp += 1;
0238     }
0239     if (tmp == N) {
0240       result = i;
0241       break;
0242     }
0243   }
0244   return result;
0245 }
0246 //------------------------------------------------------------------------
0247 //--- returns the number of bins in the direction of fVar ----------------
0248 //------------------------------------------------------------------------
0249 unsigned JetCorrectorParameters::size(unsigned fVar) const {
0250   if (fVar >= mDefinitions.nBinVar()) {
0251     std::stringstream sserr;
0252     sserr << "requested bin variable index " << fVar << " is greater than number of variables "
0253           << mDefinitions.nBinVar();
0254     handleError("JetCorrectorParameters", sserr.str());
0255   }
0256   unsigned result = 0;
0257   float tmpMin(-9999), tmpMax(-9999);
0258   for (unsigned i = 0; i < size(); ++i)
0259     if (record(i).xMin(fVar) > tmpMin && record(i).xMax(fVar) > tmpMax) {
0260       result++;
0261       tmpMin = record(i).xMin(fVar);
0262       tmpMax = record(i).xMax(fVar);
0263     }
0264   return result;
0265 }
0266 //------------------------------------------------------------------------
0267 //--- returns the vector of bin centers of fVar --------------------------
0268 //------------------------------------------------------------------------
0269 std::vector<float> JetCorrectorParameters::binCenters(unsigned fVar) const {
0270   std::vector<float> result;
0271   for (unsigned i = 0; i < size(); ++i)
0272     result.push_back(record(i).xMiddle(fVar));
0273   return result;
0274 }
0275 //------------------------------------------------------------------------
0276 //--- prints parameters on screen ----------------------------------------
0277 //------------------------------------------------------------------------
0278 void JetCorrectorParameters::printScreen() const {
0279   std::cout << "--------------------------------------------" << std::endl;
0280   std::cout << "////////  PARAMETERS: //////////////////////" << std::endl;
0281   std::cout << "--------------------------------------------" << std::endl;
0282   std::cout << "Number of binning variables:   " << definitions().nBinVar() << std::endl;
0283   std::cout << "Names of binning variables:    ";
0284   for (unsigned i = 0; i < definitions().nBinVar(); i++)
0285     std::cout << definitions().binVar(i) << " ";
0286   std::cout << std::endl;
0287   std::cout << "--------------------------------------------" << std::endl;
0288   std::cout << "Number of parameter variables: " << definitions().nParVar() << std::endl;
0289   std::cout << "Names of parameter variables:  ";
0290   for (unsigned i = 0; i < definitions().nParVar(); i++)
0291     std::cout << definitions().parVar(i) << " ";
0292   std::cout << std::endl;
0293   std::cout << "--------------------------------------------" << std::endl;
0294   std::cout << "Parametrization Formula:       " << definitions().formula() << std::endl;
0295   if (definitions().isResponse())
0296     std::cout << "Type (Response or Correction): "
0297               << "Response" << std::endl;
0298   else
0299     std::cout << "Type (Response or Correction): "
0300               << "Correction" << std::endl;
0301   std::cout << "Correction Level:              " << definitions().level() << std::endl;
0302   std::cout << "--------------------------------------------" << std::endl;
0303   std::cout << "------- Bin contents -----------------------" << std::endl;
0304   for (unsigned i = 0; i < size(); i++) {
0305     for (unsigned j = 0; j < definitions().nBinVar(); j++)
0306       std::cout << record(i).xMin(j) << " " << record(i).xMax(j) << " ";
0307     std::cout << record(i).nParameters() << " ";
0308     for (unsigned j = 0; j < record(i).nParameters(); j++)
0309       std::cout << record(i).parameter(j) << " ";
0310     std::cout << std::endl;
0311   }
0312 }
0313 //------------------------------------------------------------------------
0314 //--- prints parameters on file ----------------------------------------
0315 //------------------------------------------------------------------------
0316 void JetCorrectorParameters::printFile(const std::string& fFileName) const {
0317   std::ofstream txtFile;
0318   txtFile.open(fFileName.c_str());
0319   txtFile.setf(std::ios::right);
0320   txtFile << "{" << definitions().nBinVar() << std::setw(15);
0321   for (unsigned i = 0; i < definitions().nBinVar(); i++)
0322     txtFile << definitions().binVar(i) << std::setw(15);
0323   txtFile << definitions().nParVar() << std::setw(15);
0324   for (unsigned i = 0; i < definitions().nParVar(); i++)
0325     txtFile << definitions().parVar(i) << std::setw(15);
0326   txtFile << std::setw(definitions().formula().size() + 15) << definitions().formula() << std::setw(15);
0327   if (definitions().isResponse())
0328     txtFile << "Response" << std::setw(15);
0329   else
0330     txtFile << "Correction" << std::setw(15);
0331   txtFile << definitions().level() << "}"
0332           << "\n";
0333   for (unsigned i = 0; i < size(); i++) {
0334     for (unsigned j = 0; j < definitions().nBinVar(); j++)
0335       txtFile << record(i).xMin(j) << std::setw(15) << record(i).xMax(j) << std::setw(15);
0336     txtFile << record(i).nParameters() << std::setw(15);
0337     for (unsigned j = 0; j < record(i).nParameters(); j++)
0338       txtFile << record(i).parameter(j) << std::setw(15);
0339     txtFile << "\n";
0340   }
0341   txtFile.close();
0342 }
0343 
0344 namespace {
0345   const std::vector<std::string> labels_ = {
0346       "L1Offset",
0347       "L2Relative",
0348       "L3Absolute",
0349       "L4EMF",
0350       "L5Flavor",
0351       "L6UE",
0352       "L7Parton",
0353       "L1JPTOffset",
0354       "L2L3Residual",
0355       "Uncertainty",
0356       "L1FastJet",
0357       "UncertaintyAbsolute",
0358       "UncertaintyHighPtExtra",
0359       "UncertaintySinglePionECAL",
0360       "UncertaintyFlavor",
0361       "UncertaintyTime",
0362       "UncertaintyRelativeJEREC1",
0363       "UncertaintyRelativeJEREC2",
0364       "UncertaintyRelativeJERHF",
0365       "UncertaintyRelativeStatEC2",
0366       "UncertaintyRelativeStatHF",
0367       "UncertaintyRelativeFSR",
0368       "UncertaintyPileUpDataMC",
0369       "UncertaintyPileUpOOT",
0370       "UncertaintyPileUpPtBB",
0371       "UncertaintyPileUpBias",
0372       "UncertaintyPileUpJetRate",
0373       "UncertaintySinglePionHCAL",
0374       "UncertaintyRelativePtEC1",
0375       "UncertaintyRelativePtEC2",
0376       "UncertaintyRelativePtHF",
0377       "UncertaintyRelativeSample",
0378       "UncertaintyPileUpPtEC",
0379       "UncertaintyPileUpPtHF",
0380       "L1RC",
0381       "L1Residual",
0382       "UncertaintyAux3",
0383       "UncertaintyAux4",
0384   };
0385 
0386   const std::vector<std::string> l5Flavors_ = {"L5Flavor_bJ",
0387                                                "L5Flavor_cJ",
0388                                                "L5Flavor_qJ",
0389                                                "L5Flavor_gJ",
0390                                                "L5Flavor_bT",
0391                                                "L5Flavor_cT",
0392                                                "L5Flavor_qT",
0393                                                "L5Flavor_gT"};
0394 
0395   const std::vector<std::string> l7Partons_ = {"L7Parton_gJ",
0396                                                "L7Parton_qJ",
0397                                                "L7Parton_cJ",
0398                                                "L7Parton_bJ",
0399                                                "L7Parton_jJ",
0400                                                "L7Parton_qT",
0401                                                "L7Parton_cT",
0402                                                "L7Parton_bT",
0403                                                "L7Parton_tT"};
0404 }  // namespace
0405 std::string JetCorrectorParametersCollection::findLabel(key_type k) {
0406   if (isL5(k))
0407     return findL5Flavor(k);
0408   else if (isL7(k))
0409     return findL7Parton(k);
0410   else
0411     return labels_[k];
0412 }
0413 
0414 std::string JetCorrectorParametersCollection::findL5Flavor(key_type k) {
0415   if (k == L5Flavor)
0416     return labels_[L5Flavor];
0417   else
0418     return l5Flavors_[k / 100 - 1];
0419 }
0420 
0421 std::string JetCorrectorParametersCollection::findL7Parton(key_type k) {
0422   if (k == L7Parton)
0423     return labels_[L7Parton];
0424   else
0425     return l7Partons_[k / 1000 - 1];
0426 }
0427 
0428 void JetCorrectorParametersCollection::getSections(std::string inputFile, std::vector<std::string>& outputs) {
0429   outputs.clear();
0430   std::ifstream input(inputFile.c_str());
0431   while (!input.eof()) {
0432     char buff[10000];
0433     input.getline(buff, 10000);
0434     std::string in(buff);
0435     if (in[0] == '[') {
0436       std::string tok = getSection(in);
0437       if (!tok.empty()) {
0438         outputs.push_back(tok);
0439       }
0440     }
0441   }
0442   std::cout << "Found these sections for file: " << std::endl;
0443   copy(outputs.begin(), outputs.end(), std::ostream_iterator<std::string>(std::cout, "\n"));
0444 }
0445 
0446 // Add a JetCorrectorParameter object, possibly with flavor.
0447 void JetCorrectorParametersCollection::push_back(key_type i, value_type const& j, label_type const& flav) {
0448   std::cout << "i    = " << i << std::endl;
0449   std::cout << "flav = " << flav << std::endl;
0450   if (isL5(i)) {
0451     std::cout << "This is L5, getL5Bin = " << getL5Bin(flav) << std::endl;
0452     correctionsL5_.push_back(pair_type(getL5Bin(flav), j));
0453   } else if (isL7(i)) {
0454     std::cout << "This is L7, getL7Bin = " << getL7Bin(flav) << std::endl;
0455     correctionsL7_.push_back(pair_type(getL7Bin(flav), j));
0456   } else if (flav.empty()) {
0457     corrections_.push_back(pair_type(i, j));
0458   } else {
0459     std::cout << "***** NOT ADDING " << flav << ", corresponding position in JetCorrectorParameters is not found."
0460               << std::endl;
0461   }
0462 }
0463 
0464 // Access the JetCorrectorParameter via the key k.
0465 // key_type is hashed to deal with the three collections
0466 JetCorrectorParameters const& JetCorrectorParametersCollection::operator[](key_type k) const {
0467   collection_type::const_iterator ibegin, iend, i;
0468   if (isL5(k)) {
0469     ibegin = correctionsL5_.begin();
0470     iend = correctionsL5_.end();
0471     i = ibegin;
0472   } else if (isL7(k)) {
0473     ibegin = correctionsL7_.begin();
0474     iend = correctionsL7_.end();
0475     i = ibegin;
0476   } else {
0477     ibegin = corrections_.begin();
0478     iend = corrections_.end();
0479     i = ibegin;
0480   }
0481   for (; i != iend; ++i) {
0482     if (k == i->first)
0483       return i->second;
0484   }
0485   throw cms::Exception("InvalidInput") << " cannot find key " << static_cast<int>(k)
0486                                        << " in the JEC payload, this usually means you have to change the global tag"
0487                                        << std::endl;
0488 }
0489 
0490 // Get a list of valid keys. These will contain hashed keys
0491 // that are aware of all three collections.
0492 void JetCorrectorParametersCollection::validKeys(std::vector<key_type>& keys) const {
0493   keys.clear();
0494   for (collection_type::const_iterator ibegin = corrections_.begin(), iend = corrections_.end(), i = ibegin; i != iend;
0495        ++i) {
0496     keys.push_back(i->first);
0497   }
0498   for (collection_type::const_iterator ibegin = correctionsL5_.begin(), iend = correctionsL5_.end(), i = ibegin;
0499        i != iend;
0500        ++i) {
0501     keys.push_back(i->first);
0502   }
0503   for (collection_type::const_iterator ibegin = correctionsL7_.begin(), iend = correctionsL7_.end(), i = ibegin;
0504        i != iend;
0505        ++i) {
0506     keys.push_back(i->first);
0507   }
0508 }
0509 
0510 // Find the L5 bin for hashing
0511 JetCorrectorParametersCollection::key_type JetCorrectorParametersCollection::getL5Bin(std::string const& flav) {
0512   std::vector<std::string>::const_iterator found = find(l5Flavors_.begin(), l5Flavors_.end(), flav);
0513   if (found != l5Flavors_.end()) {
0514     return (found - l5Flavors_.begin() + 1) * 100;
0515   } else
0516     return L5Flavor;
0517 }
0518 // Find the L7 bin for hashing
0519 JetCorrectorParametersCollection::key_type JetCorrectorParametersCollection::getL7Bin(std::string const& flav) {
0520   std::vector<std::string>::const_iterator found = find(l7Partons_.begin(), l7Partons_.end(), flav);
0521   if (found != l7Partons_.end()) {
0522     return (found - l7Partons_.begin() + 1) * 1000;
0523   } else
0524     return L7Parton;
0525 }
0526 
0527 // Check if this is an L5 hashed value
0528 bool JetCorrectorParametersCollection::isL5(key_type k) { return k == L5Flavor || (k / 100 > 0 && k / 1000 == 0); }
0529 // Check if this is an L7 hashed value
0530 bool JetCorrectorParametersCollection::isL7(key_type k) { return k == L7Parton || (k / 1000 > 0); }
0531 
0532 // Find the key corresponding to each label
0533 JetCorrectorParametersCollection::key_type JetCorrectorParametersCollection::findKey(std::string const& label) const {
0534   // First check L5 corrections
0535   std::vector<std::string>::const_iterator found1 = find(l5Flavors_.begin(), l5Flavors_.end(), label);
0536   if (found1 != l5Flavors_.end()) {
0537     return getL5Bin(label);
0538   }
0539 
0540   // Next check L7 corrections
0541   std::vector<std::string>::const_iterator found2 = find(l7Partons_.begin(), l7Partons_.end(), label);
0542   if (found2 != l7Partons_.end()) {
0543     return getL7Bin(label);
0544   }
0545 
0546   // Finally check the default corrections
0547   std::vector<std::string>::const_iterator found3 = find(labels_.begin(), labels_.end(), label);
0548   if (found3 != labels_.end()) {
0549     return static_cast<key_type>(found3 - labels_.begin());
0550   }
0551 
0552   // Didn't find default corrections, throw exception
0553   throw cms::Exception("InvalidInput") << " Cannot find label " << label << std::endl;
0554 }
0555 
0556 //#include "FWCore/Framework/interface/EventSetup.h"
0557 //#include "FWCore/Framework/interface/ESHandle.h"
0558 //#include "FWCore/Framework/interface/ModuleFactory.h"
0559 #include "FWCore/Utilities/interface/typelookup.h"
0560 
0561 TYPELOOKUP_DATA_REG(JetCorrectorParameters);
0562 TYPELOOKUP_DATA_REG(JetCorrectorParametersCollection);