Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-06-04 06:14:12

0001 #include "Utilities/Testing/interface/CppUnit_testdriver.icpp"
0002 #include "cppunit/extensions/HelperMacros.h"
0003 
0004 #include "FWCore/Utilities/interface/FileInPath.h"
0005 #include "CondFormats/JetMETObjects/interface/JetCorrectorParameters.h"
0006 #include "CondFormats/JetMETObjects/interface/FactorizedJetCorrector.h"
0007 
0008 #include <iostream>
0009 #include <fstream>
0010 #include <vector>
0011 #include <string>
0012 #include <iomanip>
0013 #include <cstdlib>
0014 
0015 #include "TBenchmark.h"
0016 
0017 using std::cout;
0018 using std::endl;
0019 using std::flush;
0020 using std::ofstream;
0021 using std::setw;
0022 using std::string;
0023 using std::vector;
0024 
0025 class testJetCorrectorParameters : public CppUnit::TestFixture {
0026   CPPUNIT_TEST_SUITE(testJetCorrectorParameters);
0027   CPPUNIT_TEST(benchmarkBinIndex1D);
0028   CPPUNIT_TEST(benchmarkBinIndex3D);
0029   CPPUNIT_TEST(compareBinIndex1D);
0030   CPPUNIT_TEST(compareBinIndex3D);
0031   CPPUNIT_TEST_SUITE_END();
0032 
0033 public:
0034   testJetCorrectorParameters() {}
0035   ~testJetCorrectorParameters() {}
0036   void setUp();
0037   void tearDown();
0038   void setupCorrector(bool is3D);
0039   void destroyCorrector();
0040   void generateFiles();
0041 
0042   void benchmarkBinIndex(bool is3D);
0043   void benchmarkBinIndex1D() { benchmarkBinIndex(false); }
0044   void benchmarkBinIndex3D() { benchmarkBinIndex(true); }
0045   void compareBinIndex1D();
0046   void compareBinIndex3D();
0047 
0048   inline void loadbar3(
0049       unsigned int x, unsigned int n, unsigned int w = 50, unsigned int freq = 100, string prefix = "") {
0050     if ((x != n) && (x % (n / freq) != 0))
0051       return;
0052     float ratio = x / (float)n;
0053     int c = ratio * w;
0054 
0055     cout << prefix << std::fixed << setw(8) << std::setprecision(0) << (ratio * 100) << "% [";
0056     for (int x = 0; x < c; x++)
0057       cout << "=";
0058     for (unsigned int x = c; x < w; x++)
0059       cout << " ";
0060     cout << "] (" << x << "/" << n << ")\r" << flush;
0061   }
0062 
0063 private:
0064   string filename;
0065   TBenchmark* m_benchmark;
0066   JetCorrectorParameters* L1JetPar;
0067   vector<JetCorrectorParameters> vPar;
0068   FactorizedJetCorrector* jetCorrector;
0069   vector<float> fX;
0070   vector<float> veta;
0071   vector<float> vrho;
0072   vector<float> vpt;
0073 };
0074 
0075 ///registration of the test so that the runner can find it
0076 CPPUNIT_TEST_SUITE_REGISTRATION(testJetCorrectorParameters);
0077 
0078 void testJetCorrectorParameters::setUp() {
0079   m_benchmark = new TBenchmark();
0080   veta = {-5.191, -4.889, -4.716, -4.538, -4.363, -4.191, -4.013, -3.839, -3.664, -3.489, -3.314, -3.139,
0081           -2.964, -2.853, -2.65,  -2.5,   -2.322, -2.172, -2.043, -1.93,  -1.83,  -1.74,  -1.653, -1.566,
0082           -1.479, -1.392, -1.305, -1.218, -1.131, -1.044, -0.957, -0.879, -0.783, -0.696, -0.609, -0.522,
0083           -0.435, -0.348, -0.261, -0.174, -0.087, 0,      0.087,  0.174,  0.261,  0.348,  0.435,  0.522,
0084           0.609,  0.696,  0.783,  0.879,  0.957,  1.044,  1.131,  1.218,  1.305,  1.392,  1.479,  1.566,
0085           1.653,  1.74,   1.83,   1.93,   2.043,  2.172,  2.322,  2.5,    2.65,   2.853,  2.964,  3.139,
0086           3.314,  3.489,  3.664,  3.839,  4.013,  4.191,  4.363,  4.538,  4.716,  4.889,  5.191};
0087   for (unsigned int irho = 1; irho <= 50; irho++) {
0088     vrho.push_back(irho);
0089   }
0090   vpt = {1,   2,   3,   4,    5,    6,    7,    8,    9,    10,   11,   12,   13,   14,   15,
0091          17,  20,  23,  27,   30,   35,   40,   45,   57,   72,   90,   120,  150,  200,  300,
0092          400, 550, 750, 1000, 1500, 2000, 2500, 3000, 3500, 4000, 4500, 5000, 5500, 6000, 6500};
0093 
0094   generateFiles();
0095 }
0096 
0097 void testJetCorrectorParameters::tearDown() {
0098   if (m_benchmark != nullptr)
0099     delete m_benchmark;
0100 }
0101 
0102 void testJetCorrectorParameters::setupCorrector(bool is3D) {
0103   string path = "CondFormats/JetMETObjects/data/";
0104   try {
0105     string to_try = path + ((is3D) ? "testJetCorrectorParameters_3D_L1FastJet_AK4PFchs.txt"
0106                                    : "testJetCorrectorParameters_1D_L1FastJet_AK4PFchs.txt");
0107     edm::FileInPath strFIP(to_try);
0108     filename = strFIP.fullPath();
0109   } catch (edm::Exception const& ex) {
0110     throw ex;
0111   }
0112 
0113   L1JetPar = new JetCorrectorParameters(filename);
0114   vPar.push_back(*L1JetPar);
0115   jetCorrector = new FactorizedJetCorrector(vPar);
0116 }
0117 
0118 void testJetCorrectorParameters::destroyCorrector() {
0119   delete L1JetPar;
0120   delete jetCorrector;
0121 }
0122 
0123 void testJetCorrectorParameters::generateFiles() {
0124   string path = std::getenv("CMSSW_BASE");
0125   path += "/src/CondFormats/JetMETObjects/data/";
0126   string name1D = "testJetCorrectorParameters_1D_L1FastJet_AK4PFchs.txt";
0127   string name3D = "testJetCorrectorParameters_3D_L1FastJet_AK4PFchs.txt";
0128 
0129   ofstream txtFile1D;
0130   txtFile1D.open((path + name1D).c_str());
0131   CPPUNIT_ASSERT(txtFile1D.is_open());
0132   txtFile1D << "{1 JetEta 3 Rho JetPt JetA 1 Correction L1FastJet}" << endl;
0133   for (unsigned int ieta = 0; ieta < veta.size() - 1; ieta++) {
0134     txtFile1D << std::setw(15) << veta[ieta] << std::setw(15) << veta[ieta + 1] << std::setw(15) << "6" << std::setw(15)
0135               << "0" << std::setw(15) << "200" << std::setw(15) << "1" << std::setw(15) << "6500" << std::setw(15)
0136               << "0" << std::setw(15) << "10" << std::setw(15) << endl;
0137   }
0138   txtFile1D.close();
0139 
0140   vector<float> ptbins = {4,  8,  10,  13,  16,  19,  22,  25,  30,  35,   40,  45,
0141                           60, 75, 100, 160, 230, 330, 440, 600, 820, 1100, 1400};
0142   ofstream txtFile3D;
0143   txtFile3D.open((path + name3D).c_str());
0144   CPPUNIT_ASSERT(txtFile3D.is_open());
0145   txtFile3D << "{3 JetEta Rho JetPt 1 JetPt 1 Correction L1FastJet}" << endl;
0146   for (unsigned int ieta = 0; ieta < veta.size() - 1; ieta++) {
0147     for (unsigned int irho = 0; irho < vrho.size() - 1; irho++) {
0148       for (unsigned int ipt = 0; ipt < ptbins.size() - 1; ipt++) {
0149         txtFile3D << std::setw(15) << veta[ieta] << std::setw(15) << veta[ieta + 1] << std::setw(15) << vrho[irho]
0150                   << std::setw(15);
0151         if (vrho[irho + 1] == 50)
0152           txtFile3D << "200" << std::setw(15);
0153         else
0154           txtFile3D << vrho[irho + 1] << std::setw(15);
0155         txtFile3D << ptbins[ipt] << std::setw(15);
0156         if (ipt + 1 == ptbins.size() - 1)
0157           txtFile3D << "6500" << std::setw(15);
0158         else
0159           txtFile3D << ptbins[ipt + 1] << std::setw(15);
0160         txtFile3D << "2" << std::setw(15) << ptbins[ipt] << std::setw(15) << ptbins[ipt + 1] << std::setw(15) << endl;
0161       }
0162     }
0163   }
0164   txtFile3D.close();
0165 }
0166 
0167 void testJetCorrectorParameters::benchmarkBinIndex(bool is3D) {
0168   float oldCPU = 0, newCPU = 0, oldReal = 0, newReal = 0;
0169   unsigned int ntests = (is3D) ? 100000 : 1000000;
0170   if (is3D)
0171     fX = {5.0, 50.0, 100.0};
0172   else
0173     fX = {5.0};
0174   setupCorrector(is3D);
0175 
0176   cout << endl << "testJetCorrectorParameters::benchmarkBinIndex NTests = " << ntests << endl;
0177   cout << "testJetCorrectorParameters::benchmarkBinIndex benchmarking binIndex with file " << filename << " ... "
0178        << flush;
0179   m_benchmark->Reset();
0180   m_benchmark->Start("event");
0181   for (unsigned int i = 0; i < ntests; i++) {
0182     L1JetPar->binIndex(fX);
0183   }
0184   m_benchmark->Stop("event");
0185   cout << "DONE" << endl
0186        << "testJetCorrectorParameters::benchmarkBinIndex" << endl
0187        << "\tCPU time = " << m_benchmark->GetCpuTime("event") / double(ntests) << " s" << endl
0188        << "\tReal time = " << m_benchmark->GetRealTime("event") / double(ntests) << " s" << endl;
0189   oldCPU = m_benchmark->GetCpuTime("event") / double(ntests);
0190   oldReal = m_benchmark->GetRealTime("event") / double(ntests);
0191 
0192   cout << "testJetCorrectorParameters::benchmarkBinIndex benchmarking binIndexN with file " << filename << " ... "
0193        << flush;
0194   m_benchmark->Reset();
0195   m_benchmark->Start("event");
0196   for (unsigned int i = 0; i < ntests; i++) {
0197     L1JetPar->binIndexN(fX);
0198   }
0199   m_benchmark->Stop("event");
0200   cout << "DONE" << endl
0201        << "testJetCorrectorParameters::benchmarkBinIndex" << endl
0202        << "\tCPU time = " << m_benchmark->GetCpuTime("event") / double(ntests) << " s" << endl
0203        << "\tReal time = " << m_benchmark->GetRealTime("event") / double(ntests) << " s" << endl;
0204   newCPU = m_benchmark->GetCpuTime("event") / double(ntests);
0205   newReal = m_benchmark->GetRealTime("event") / double(ntests);
0206 
0207   cout << "testJetCorrectorParameters::benchmarkBinIndex" << endl
0208        << "\tCPU speedup = " << oldCPU / newCPU << endl
0209        << "\tReal speedup = " << oldReal / newReal << endl;
0210 
0211   CPPUNIT_ASSERT(oldCPU / newCPU >= 1.0);
0212   //CPPUNIT_ASSERT(oldReal/newReal >= 1.0); //this might fail not due to longer L1JetPar->binIndexN(fX) execution
0213   if (oldReal / newReal >= 1.0)
0214     cout << "newReal value increased oldReal, which might be due to system load" << endl;
0215 
0216   destroyCorrector();
0217 }
0218 
0219 void testJetCorrectorParameters::compareBinIndex1D() {
0220   fX = {0.0};
0221   float eta = -9999;
0222   int oldBin = -1, newBin = -1;
0223   setupCorrector(false);
0224 
0225   cout << endl << "testJetCorrectorParameters::compareBinIndex1D" << endl;
0226   for (unsigned int ieta = 0; ieta < veta.size() - 1; ieta++) {
0227     loadbar3(ieta + 1, veta.size() - 1, 50, 10, "\tProgress:");
0228     eta = (veta[ieta] + veta[ieta + 1]) / 2.0;
0229     fX = {eta};
0230     oldBin = L1JetPar->binIndex(fX);
0231     newBin = L1JetPar->binIndexN(fX);
0232     if ((oldBin < 0 && newBin >= 0) || (oldBin >= 0 && newBin < 0)) {
0233       cout << "ERROR::testJetCorrectorParameters::compareBinIndex1D Unable to find the right bin for (eta)=(" << eta
0234            << ")" << endl
0235            << "\t(oldBin,newBin)=(" << oldBin << "," << newBin << ")" << endl;
0236       CPPUNIT_ASSERT(oldBin >= 0 && newBin >= 0);
0237     } else if (oldBin != newBin) {
0238       cout << "ERROR::testJetCorrectorParameters::compareBinIndex1D oldBin!=newBin (" << oldBin << "!=" << newBin
0239            << ") for (eta)=(" << eta << ")" << endl;
0240       CPPUNIT_ASSERT(oldBin == newBin);
0241     }
0242 
0243     jetCorrector->setJetEta(eta);
0244     jetCorrector->setRho(50.0);
0245     jetCorrector->setJetA(0.5);
0246     jetCorrector->setJetPt(100.0);
0247     CPPUNIT_ASSERT(jetCorrector->getCorrection() >= 0);
0248   }
0249   destroyCorrector();
0250   cout << endl
0251        << "testJetCorrectorParameters::compareBinIndex1D All bins match between the linear and non-linear search "
0252           "algorithms for 1D files."
0253        << endl;
0254 }
0255 
0256 void testJetCorrectorParameters::compareBinIndex3D() {
0257   fX = {0.0, 0.0, 0.0};
0258   float eta = -9999, rho = -9999, pt = -9999;
0259   int oldBin = -1, newBin = -1;
0260   setupCorrector(true);
0261 
0262   cout << endl << "testJetCorrectorParameters::compareBinIndex3D" << endl;
0263   for (unsigned int ieta = 0; ieta < veta.size() - 1; ieta++) {
0264     for (unsigned int irho = 0; irho < vrho.size() - 1; irho++) {
0265       for (unsigned int ipt = 0; ipt < vpt.size() - 1; ipt++) {
0266         loadbar3(ieta * ((vrho.size() - 1) * (vpt.size() - 1)) + irho * (vpt.size() - 1) + ipt + 1,
0267                  (veta.size() - 1) * (vrho.size() - 1) * (vpt.size() - 1),
0268                  50,
0269                  100,
0270                  "\tProgress");
0271         eta = (veta[ieta] + veta[ieta + 1]) / 2.0;
0272         rho = (vrho[irho] + vrho[irho + 1]) / 2.0;
0273         pt = (vpt[ipt] + vpt[ipt + 1]) / 2.0;
0274         fX = {eta, rho, pt};
0275         oldBin = L1JetPar->binIndex(fX);
0276         newBin = L1JetPar->binIndexN(fX);
0277         if ((oldBin < 0 && newBin >= 0) || (oldBin >= 0 && newBin < 0)) {
0278           cout << "ERROR::testJetCorrectorParameters::compareBinIndex3D Unable to find the right bin for (eta,rho,pt)=("
0279                << eta << "," << rho << "," << pt << ")" << endl
0280                << "\t(oldBin,newBin)=(" << oldBin << "," << newBin << ")" << endl;
0281           CPPUNIT_ASSERT(oldBin >= 0 && newBin >= 0);
0282         } else if (oldBin != newBin) {
0283           cout << "ERROR::testJetCorrectorParameters::compareBinIndex3D oldBin!=newBin (" << oldBin << "!=" << newBin
0284                << ") for (eta,rho,pt)=(" << eta << "," << rho << "," << pt << ")" << endl;
0285           CPPUNIT_ASSERT(oldBin == newBin);
0286         }
0287 
0288         jetCorrector->setJetEta(eta);
0289         jetCorrector->setRho(rho);
0290         jetCorrector->setJetA(0.5);
0291         jetCorrector->setJetPt(pt);
0292         CPPUNIT_ASSERT(jetCorrector->getCorrection() >= 0);
0293       }
0294     }
0295   }
0296   destroyCorrector();
0297   cout << endl
0298        << "testJetCorrectorParameters::compareBinIndex3D All bins match between the linear and non-linear search "
0299           "algorithms for 3D files."
0300        << endl;
0301 }