JetMETObjects/src/METCorrectorParameters.cc

0001 // Generic parameters for MET corrections
0002 //
0003 #include "CondFormats/JetMETObjects/interface/METCorrectorParameters.h"
0004 #include "CondFormats/JetMETObjects/interface/Utilities.h"
0005 #include <iostream>
0006 #include <iomanip>
0007 #include <fstream>
0008 #include <sstream>
0009 #include <cstdlib>
0010 #include <algorithm>
0011 #include <cmath>
0012 #include <iterator>
0013 //#include <string>
0014
0015 //------------------------------------------------------------------------
0016 //--- METCorrectorParameters::Definitions constructor --------------------
0017 //--- takes specific arguments for the member variables ------------------
0018 //------------------------------------------------------------------------
0019 METCorrectorParameters::Definitions::Definitions(const std::vector<std::string>& fBinVar,
0020                                                  const std::vector<std::string>& fParVar,
0021                                                  const std::string& fFormula) {
0022   for (unsigned i = 0; i < fBinVar.size(); i++)
0023     mBinVar.push_back(fBinVar[i]);
0024   for (unsigned i = 0; i < fParVar.size(); i++)
0025     mParVar.push_back(fParVar[i]);
0026   mFormula = fFormula;
0027 }
0028 //------------------------------------------------------------------------
0029 //--- METCorrectorParameters::Definitions constructor --------------------
0030 //--- reads the member variables from a string ---------------------------
0031 //------------------------------------------------------------------------
0032 METCorrectorParameters::Definitions::Definitions(const std::string& fLine) {
0033   std::vector<std::string> tokens = getTokens(fLine);
0034   // corrType N_bin binVa.. N_var var... formula
0035   if (!tokens.empty()) {
0036     if (tokens.size() < 6) {
0037       std::stringstream sserr;
0038       sserr << "(line " << fLine << "): less than 6 expected tokens:" << tokens.size();
0039       handleError("METCorrectorParameters::Definitions", sserr.str());
0040     }
0041     // No. of Bin Variable
0042     std::cout << "Definitions===========" << std::endl;
0043     ptclType = getSigned(tokens[0]);
0044     std::cout << tokens[0] << "\t";
0045     unsigned nBinVar = getUnsigned(tokens[1]);
0046     std::cout << tokens[1] << "\t";
0047     // No of Parameterization Variable
0048     unsigned nParVar = getUnsigned(tokens[nBinVar + 2]);
0049     std::cout << tokens[2] << "\t";
0050     for (unsigned i = 0; i < nBinVar; i++) {
0051       mBinVar.push_back(tokens[i + 2]);
0052       std::cout << tokens[i + 2] << "\t";
0053     }
0054     for (unsigned i = 0; i < nParVar; i++) {
0055       mParVar.push_back(tokens[nBinVar + 3 + i]);
0056       std::cout << tokens[nBinVar + 3 + i] << "\t";
0057     }
0058     mFormula = tokens[nParVar + nBinVar + 3];
0059     std::cout << tokens[nParVar + nBinVar + 3] << std::endl;
0060   }
0061 }
0062 //------------------------------------------------------------------------
0063 //--- METCorrectorParameters::Record constructor -------------------------
0064 //--- reads the member variables from a string ---------------------------
0065 //------------------------------------------------------------------------
0066 METCorrectorParameters::Record::Record(const std::string& fLine, unsigned fNvar) : mMin(0), mMax(0), mParameters(0) {
0067   mNvar = fNvar;
0068   // quckly parse the line
0069   std::vector<std::string> tokens = getTokens(fLine);
0070   if (!tokens.empty()) {
0071     if (tokens.size() < 6) {
0072       std::stringstream sserr;
0073       sserr << "(line " << fLine << "): "
0074             << "three tokens expected, " << tokens.size() << " provided.";
0075       handleError("METCorrectorParameters::Record", sserr.str());
0076     }
0077     std::cout << "Record ===============" << std::endl;
0078     for (unsigned i = 0; i < mNvar; i++) {
0079       mMin.push_back(getFloat(tokens[i * 2]));
0080       mMax.push_back(getFloat(tokens[i * 2 + 1]));
0081       std::cout << tokens[i * 2] << "\t";
0082       std::cout << tokens[i * 2 + 1] << "\t";
0083     }
0084     unsigned nParam = getUnsigned(tokens[2 * mNvar]);
0085     std::cout << tokens[2 * mNvar] << "\t";
0086     if (nParam != tokens.size() - (2 * mNvar + 1)) {
0087       std::stringstream sserr;
0088       sserr << "(line " << fLine << "): " << tokens.size() - (2 * mNvar + 1) << " parameters, but nParam=" << nParam
0089             << ".";
0090       handleError("METCorrectorParameters::Record", sserr.str());
0091     }
0092     for (unsigned i = (2 * mNvar + 1); i < tokens.size(); ++i) {
0093       mParameters.push_back(getFloat(tokens[i]));
0094       std::cout << tokens[i] << "\t";
0095     }
0096     std::cout << std::endl;
0097   }
0098 }
0099 //------------------------------------------------------------------------
0100 //--- JetCorrectorParameters constructor ---------------------------------
0101 //--- reads the member variables from a string ---------------------------
0102 //------------------------------------------------------------------------
0103 METCorrectorParameters::METCorrectorParameters(const std::string& fFile, const std::string& fSection) {
0104   std::ifstream input(fFile.c_str());
0105   std::string currentSection = "";
0106   std::string line;
0107   std::string currentDefinitions = "";
0108   while (std::getline(input, line)) {
0109     //std::cout << " Line of parameters " << line << std::endl;
0110     std::string section = getSection(line);
0111     std::string tmp = getDefinitions(line);
0112     if (!section.empty() && tmp.empty()) {
0113       currentSection = section;
0114       continue;
0115     }
0116     if (currentSection == fSection) {
0117       if (!tmp.empty()) {
0118         currentDefinitions = tmp;
0119         continue;
0120       }
0121       Definitions definitions(currentDefinitions);
0122       if (!(definitions.nBinVar() == 0 && definitions.formula().empty()))
0123         mDefinitions = definitions;
0124       Record record(line, mDefinitions.nBinVar());
0125       bool check(true);
0126       for (unsigned i = 0; i < mDefinitions.nBinVar(); ++i)
0127         if (record.xMin(i) == 0 && record.xMax(i) == 0)
0128           check = false;
0129       if (record.nParameters() == 0)
0130         check = false;
0131       if (check)
0132         mRecords.push_back(record);
0133     }
0134   }
0135   if (currentDefinitions.empty())
0136     handleError("METCorrectorParameters", "No definitions found!!!");
0137   if (mRecords.empty() && currentSection.empty())
0138     mRecords.push_back(Record());
0139   if (mRecords.empty() && !currentSection.empty()) {
0140     std::stringstream sserr;
0141     sserr << "the requested section " << fSection << " doesn't exist!";
0142     handleError("METCorrectorParameters", sserr.str());
0143   }
0144   std::sort(mRecords.begin(), mRecords.end());
0145   valid_ = true;
0146 }
0147 //------------------------------------------------------------------------
0148 //--- returns the index of the record defined by fX ----------------------
0149 //------------------------------------------------------------------------
0150 /*
0151 int METCorrectorParameters::binIndex(const std::vector<float>& fX) const
0152 {
0153   int result = -1;
0154   unsigned N = mDefinitions.nVar();
0155   if (N != fX.size())
0156     {
0157       std::stringstream sserr;
0158       sserr<<"# bin variables "<<N<<" doesn't correspont to requested #: "<<fX.size();
0159       handleError("METCorrectorParameters",sserr.str());
0160     }
0161   unsigned tmp;
0162   for (unsigned i = 0; i < size(); ++i)
0163     {
0164       tmp = 0;
0165       for (unsigned j=0;j<N;j++)
0166         if (fX[j] >= record(i).xMin(j) && fX[j] < record(i).xMax(j))
0167           tmp+=1;
0168       if (tmp==N)
0169         {
0170           result = i;
0171           break;
0172         }
0173     }
0174   return result;
0175 }
0176 */
0177 //------------------------------------------------------------------------
0178 //--- returns the neighbouring bins of fIndex in the direction of fVar ---
0179 //------------------------------------------------------------------------
0180 /*
0181 int METCorrectorParameters::neighbourBin(unsigned fIndex, unsigned fVar, bool fNext) const
0182 {
0183   int result = -1;
0184   unsigned N = mDefinitions.nVar();
0185   if (fVar >= N)
0186     {
0187       std::stringstream sserr;
0188       sserr<<"# of bin variables "<<N<<" doesn't correspond to requested #: "<<fVar;
0189       handleError("METCorrectorParameters",sserr.str());
0190     }
0191   unsigned tmp;
0192   for (unsigned i = 0; i < size(); ++i)
0193     {
0194       tmp = 0;
0195       for (unsigned j=0;j<fVar;j++)
0196         if (fabs(record(i).xMin(j)-record(fIndex).xMin(j))<0.0001)
0197           tmp+=1;
0198       for (unsigned j=fVar+1;j<N;j++)
0199         if (fabs(record(i).xMin(j)-record(fIndex).xMin(j))<0.0001)
0200           tmp+=1;
0201       if (tmp<N-1)
0202         continue;
0203       if (tmp==N-1)
0204         {
0205           if (fNext)
0206             if (fabs(record(i).xMin(fVar)-record(fIndex).xMax(fVar))<0.0001)
0207               tmp+=1;
0208           if (!fNext)
0209             if (fabs(record(i).xMax(fVar)-record(fIndex).xMin(fVar))<0.0001)
0210               tmp+=1;
0211         }
0212       if (tmp==N)
0213         {
0214           result = i;
0215           break;
0216         }
0217     }
0218   return result;
0219 }
0220 */
0221 //------------------------------------------------------------------------
0222 //--- returns the number of bins in the direction of fVar ----------------
0223 //------------------------------------------------------------------------
0224 /*
0225 unsigned METCorrectorParameters::size(unsigned fVar) const
0226 {
0227   if (fVar >= mDefinitions.nVar())
0228     {
0229       std::stringstream sserr;
0230       sserr<<"requested bin variable index "<<fVar<<" is greater than number of variables "<<mDefinitions.nVar();
0231       handleError("METCorrectorParameters",sserr.str());
0232     }
0233   unsigned result = 0;
0234   float tmpMin(-9999),tmpMax(-9999);
0235   for (unsigned i = 0; i < size(); ++i)
0236     if (record(i).xMin(fVar) > tmpMin && record(i).xMax(fVar) > tmpMax)
0237       {
0238         result++;
0239         tmpMin = record(i).xMin(fVar);
0240         tmpMax = record(i).xMax(fVar);
0241       }
0242   return result;
0243 }
0244 */
0245 //------------------------------------------------------------------------
0246 //--- returns the vector of bin centers of fVar --------------------------
0247 //------------------------------------------------------------------------
0248 /*
0249 std::vector<float> METCorrectorParameters::binCenters(unsigned fVar) const
0250 {
0251   std::vector<float> result;
0252   for (unsigned i = 0; i < size(); ++i)
0253     result.push_back(record(i).xMiddle(fVar));
0254   return result;
0255 }
0256 */
0257 //------------------------------------------------------------------------
0258 //--- prints parameters on screen ----------------------------------------
0259 //------------------------------------------------------------------------
0260 void METCorrectorParameters::printScreen() const {
0261   std::cout << "--------------------------------------------" << std::endl;
0262   std::cout << "////////  PARAMETERS: //////////////////////" << std::endl;
0263   std::cout << "--------------------------------------------" << std::endl;
0264   std::cout << "Number of binning variables:   " << definitions().nBinVar() << std::endl;
0265   std::cout << "Names of binning variables:    ";
0266   for (unsigned i = 0; i < definitions().nBinVar(); i++)
0267     std::cout << definitions().binVar(i) << " ";
0268   std::cout << std::endl;
0269   std::cout << "--------------------------------------------" << std::endl;
0270   std::cout << "Number of parameter variables: " << definitions().nParVar() << std::endl;
0271   std::cout << "Names of parameter variables:  ";
0272   for (unsigned i = 0; i < definitions().nParVar(); i++)
0273     std::cout << definitions().parVar(i) << " ";
0274   std::cout << std::endl;
0275   std::cout << "--------------------------------------------" << std::endl;
0276   std::cout << "Parametrization Formula:       " << definitions().formula() << std::endl;
0277   std::cout << "--------------------------------------------" << std::endl;
0278   std::cout << "------- Bin contents -----------------------" << std::endl;
0279   for (unsigned i = 0; i < size(); i++) {
0280     for (unsigned j = 0; j < definitions().nBinVar(); j++)
0281       std::cout << record(i).xMin(j) << " " << record(i).xMax(j) << " ";
0282     std::cout << record(i).nParameters() << " ";
0283     for (unsigned j = 0; j < record(i).nParameters(); j++)
0284       std::cout << record(i).parameter(j) << " ";
0285     std::cout << std::endl;
0286   }
0287 }
0288 //------------------------------------------------------------------------
0289 //--- prints parameters on file ----------------------------------------
0290 //------------------------------------------------------------------------
0291 void METCorrectorParameters::printFile(const std::string& fFileName) const {
0292   std::ofstream txtFile;
0293   txtFile.open(fFileName.c_str());
0294   txtFile.setf(std::ios::right);
0295   txtFile << "{" << definitions().nBinVar() << std::setw(15);
0296   for (unsigned i = 0; i < definitions().nBinVar(); i++)
0297     txtFile << definitions().binVar(i) << std::setw(15);
0298   txtFile << definitions().nParVar() << std::setw(15);
0299   for (unsigned i = 0; i < definitions().nParVar(); i++)
0300     txtFile << definitions().parVar(i) << std::setw(15);
0301   txtFile << std::setw(definitions().formula().size() + 15) << definitions().formula() << std::setw(15);
0302   txtFile << "}"
0303           << "\n";
0304   for (unsigned i = 0; i < size(); i++) {
0305     for (unsigned j = 0; j < definitions().nBinVar(); j++)
0306       txtFile << record(i).xMin(j) << std::setw(15) << record(i).xMax(j) << std::setw(15);
0307     txtFile << record(i).nParameters() << std::setw(15);
0308     for (unsigned j = 0; j < record(i).nParameters(); j++)
0309       txtFile << record(i).parameter(j) << std::setw(15);
0310     txtFile << "\n";
0311   }
0312   txtFile.close();
0313 }
0314
0315 namespace {
0316   const std::vector<std::string> labels_ = {"MiniAod"};
0317   const std::vector<std::string> MiniAodSource_ = {"MiniAod_ShiftX", "MiniAod_ShiftY"};
0318
0319 }  //namespace
0320
0321 std::string METCorrectorParametersCollection::findLabel(key_type k) {
0322   std::cout << "findLabel with key_type: " << k << std::endl;
0323   if (isMiniAod(k)) {
0324     std::cout << "is MiniAod" << std::endl;
0325     return findMiniAodSource(k);
0326   }
0327   return labels_[k];
0328 }
0329 std::string METCorrectorParametersCollection::findMiniAodSource(key_type k) {
0330   if (k == MiniAod)
0331     return labels_[MiniAod];
0332   else
0333     return MiniAodSource_[k - MiniAod * 100 - 1];
0334 }
0335 void METCorrectorParametersCollection::getSections(std::string inputFile, std::vector<std::string>& outputs) {
0336   outputs.clear();
0337   std::ifstream input(inputFile.c_str());
0338   while (!input.eof()) {
0339     char buff[10000];
0340     input.getline(buff, 10000);
0341     std::string in(buff);
0342     if (in[0] == '[') {
0343       std::string tok = getSection(in);
0344       if (!tok.empty()) {
0345         outputs.push_back(tok);
0346       }
0347     }
0348   }
0349   std::cout << "Found these sections for file: " << std::endl;
0350   copy(outputs.begin(), outputs.end(), std::ostream_iterator<std::string>(std::cout, "\n"));
0351 }
0352
0353 // Add a METCorrectorParameter object.
0354 void METCorrectorParametersCollection::push_back(key_type i, value_type const& j, label_type const& source) {
0355   std::cout << "i    = " << i << std::endl;
0356   std::cout << "source = " << source << std::endl;
0357   if (isMiniAod(i)) {
0358     std::cout << "This is MiniAod, getMiniAodBin = " << getMiniAodBin(source) << std::endl;
0359     correctionsMiniAod_.push_back(pair_type(getMiniAodBin(source), j));
0360   } else {
0361     std::cout << "***** NOT ADDING " << source << ", corresponding position in METCorrectorParameters is not found."
0362               << std::endl;
0363   }
0364 }
0365
0366 // Access the METCorrectorParameter via the key k.
0367 // key_type is hashed to deal with the three collections
0368 METCorrectorParameters const& METCorrectorParametersCollection::operator[](key_type k) const {
0369   collection_type::const_iterator ibegin, iend, i;
0370   if (isMiniAod(k)) {
0371     ibegin = correctionsMiniAod_.begin();
0372     iend = correctionsMiniAod_.end();
0373     i = ibegin;
0374   }
0375   for (; i != iend; ++i) {
0376     if (k == i->first)
0377       return i->second;
0378   }
0379   throw cms::Exception("InvalidInput") << " cannot find key " << static_cast<int>(k)
0380                                        << " in the METC payload, this usually means you have to change the global tag"
0381                                        << std::endl;
0382 }
0383
0384 // Get a list of valid keys. These will contain hashed keys
0385 // that are aware of all three collections.
0386 void METCorrectorParametersCollection::validKeys(std::vector<key_type>& keys) const {
0387   keys.clear();
0388   for (collection_type::const_iterator ibegin = correctionsMiniAod_.begin(),
0389                                        iend = correctionsMiniAod_.end(),
0390                                        i = ibegin;
0391        i != iend;
0392        ++i) {
0393     keys.push_back(i->first);
0394   }
0395 }
0396
0397 METCorrectorParametersCollection::key_type METCorrectorParametersCollection::getMiniAodBin(std::string const& source) {
0398   std::vector<std::string>::const_iterator found = find(MiniAodSource_.begin(), MiniAodSource_.end(), source);
0399   if (found != MiniAodSource_.end()) {
0400     return (found - MiniAodSource_.begin() + 1) + MiniAod * 100;
0401   } else
0402     return MiniAod;
0403 }
0404
0405 bool METCorrectorParametersCollection::isMiniAod(key_type k) {
0406   return k == MiniAod || (k > MiniAod * 100 && k < MiniAod * 100 + 100);
0407 }
0408
0409 #include "FWCore/Utilities/interface/typelookup.h"
0410
0411 TYPELOOKUP_DATA_REG(METCorrectorParameters);
0412 TYPELOOKUP_DATA_REG(METCorrectorParametersCollection);