Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 12:50:18

0001 #include <cppunit/extensions/HelperMacros.h>
0002 #include "CommonTools/Utils/src/formulaConstantEvaluator.h"
0003 #include "CommonTools/Utils/src/formulaParameterEvaluator.h"
0004 #include "CommonTools/Utils/src/formulaVariableEvaluator.h"
0005 #include "CommonTools/Utils/src/formulaBinaryOperatorEvaluator.h"
0006 #include "CommonTools/Utils/src/formulaFunctionOneArgEvaluator.h"
0007 #include "CommonTools/Utils/src/formulaFunctionTwoArgsEvaluator.h"
0008 #include "CommonTools/Utils/src/formulaUnaryMinusEvaluator.h"
0009 
0010 #include "CommonTools/Utils/interface/FormulaEvaluator.h"
0011 
0012 #include "FWCore/Utilities/interface/Exception.h"
0013 
0014 #include <algorithm>
0015 #include <cmath>
0016 #include "TMath.h"
0017 
0018 class testFormulaEvaluator : public CppUnit::TestFixture {
0019   CPPUNIT_TEST_SUITE(testFormulaEvaluator);
0020   CPPUNIT_TEST(checkEvaluators);
0021   CPPUNIT_TEST(checkFormulaEvaluator);
0022   CPPUNIT_TEST_SUITE_END();
0023 
0024 public:
0025   void setUp() {}
0026   void tearDown() {}
0027   void checkEvaluators();
0028   void checkFormulaEvaluator();
0029 };
0030 
0031 namespace {
0032   bool compare(double iLHS, double iRHS) {
0033     return std::fabs(iLHS) != 0 ? (std::fabs(iLHS - iRHS) < 1E-6 * std::fabs(iLHS))
0034                                 : (std::fabs(iLHS) == std::fabs(iRHS));
0035   }
0036 }  // namespace
0037 
0038 CPPUNIT_TEST_SUITE_REGISTRATION(testFormulaEvaluator);
0039 
0040 void testFormulaEvaluator::checkEvaluators() {
0041   using namespace reco::formula;
0042   {
0043     ConstantEvaluator c{4};
0044 
0045     CPPUNIT_ASSERT(c.evaluate(nullptr, nullptr) == 4.);
0046   }
0047 
0048   {
0049     ParameterEvaluator pe{0};
0050 
0051     double p = 5.;
0052 
0053     CPPUNIT_ASSERT(pe.evaluate(nullptr, &p) == p);
0054   }
0055 
0056   {
0057     VariableEvaluator ve{0};
0058 
0059     double v = 3.;
0060 
0061     CPPUNIT_ASSERT(ve.evaluate(&v, nullptr) == v);
0062   }
0063 
0064   {
0065     auto cl = std::unique_ptr<ConstantEvaluator>(new ConstantEvaluator(4));
0066     auto cr = std::unique_ptr<ConstantEvaluator>(new ConstantEvaluator(3));
0067 
0068     BinaryOperatorEvaluator<std::minus<double>> be(std::move(cl), std::move(cr), EvaluatorBase::Precedence::kPlusMinus);
0069 
0070     CPPUNIT_ASSERT(be.evaluate(nullptr, nullptr) == 1.);
0071   }
0072 
0073   {
0074     auto cl = std::unique_ptr<ConstantEvaluator>(new ConstantEvaluator(4));
0075 
0076     FunctionOneArgEvaluator f(std::move(cl), [](double v) { return std::exp(v); });
0077 
0078     CPPUNIT_ASSERT(f.evaluate(nullptr, nullptr) == std::exp(4.));
0079   }
0080 
0081   {
0082     auto cl = std::unique_ptr<ConstantEvaluator>(new ConstantEvaluator(4));
0083     auto cr = std::unique_ptr<ConstantEvaluator>(new ConstantEvaluator(3));
0084 
0085     FunctionTwoArgsEvaluator f(std::move(cl), std::move(cr), [](double v1, double v2) { return std::max(v1, v2); });
0086 
0087     CPPUNIT_ASSERT(f.evaluate(nullptr, nullptr) == 4.);
0088   }
0089 
0090   {
0091     auto cl = std::unique_ptr<ConstantEvaluator>(new ConstantEvaluator(4));
0092 
0093     UnaryMinusEvaluator f(std::move(cl));
0094 
0095     CPPUNIT_ASSERT(f.evaluate(nullptr, nullptr) == -4.);
0096   }
0097 }
0098 
0099 void testFormulaEvaluator::checkFormulaEvaluator() {
0100   {
0101     reco::FormulaEvaluator f("5");
0102 
0103     std::vector<double> emptyV;
0104 
0105     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 5.);
0106   }
0107 
0108   {
0109     reco::FormulaEvaluator f("3+2");
0110 
0111     std::vector<double> emptyV;
0112 
0113     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 5.);
0114   }
0115 
0116   {
0117     reco::FormulaEvaluator f(" 3 + 2 ");
0118 
0119     std::vector<double> emptyV;
0120 
0121     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 5.);
0122   }
0123 
0124   {
0125     reco::FormulaEvaluator f("3-2");
0126 
0127     std::vector<double> emptyV;
0128 
0129     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 1.);
0130   }
0131 
0132   {
0133     reco::FormulaEvaluator f("3*2");
0134 
0135     std::vector<double> emptyV;
0136 
0137     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 6.);
0138   }
0139 
0140   {
0141     reco::FormulaEvaluator f("6/2");
0142 
0143     std::vector<double> emptyV;
0144 
0145     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 3.);
0146   }
0147 
0148   {
0149     reco::FormulaEvaluator f("3^2");
0150 
0151     std::vector<double> emptyV;
0152 
0153     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 9.);
0154   }
0155 
0156   {
0157     reco::FormulaEvaluator f("4*3^2");
0158 
0159     std::vector<double> emptyV;
0160 
0161     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 36.);
0162   }
0163   {
0164     reco::FormulaEvaluator f("3^2*4");
0165 
0166     std::vector<double> emptyV;
0167 
0168     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 36.);
0169   }
0170 
0171   {
0172     reco::FormulaEvaluator f("1+2*3^4+5*2+6*2");
0173 
0174     std::vector<double> emptyV;
0175 
0176     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 1 + 2 * (3 * 3 * 3 * 3) + 5 * 2 + 6 * 2);
0177   }
0178 
0179   {
0180     reco::FormulaEvaluator f("1+3^4*2+5*2+6*2");
0181 
0182     std::vector<double> emptyV;
0183 
0184     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 1 + 2 * (3 * 3 * 3 * 3) + 5 * 2 + 6 * 2);
0185   }
0186 
0187   {
0188     reco::FormulaEvaluator f("3<=2");
0189 
0190     std::vector<double> emptyV;
0191 
0192     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 0.);
0193   }
0194   {
0195     reco::FormulaEvaluator f("2<=3");
0196 
0197     std::vector<double> emptyV;
0198 
0199     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 1.);
0200   }
0201   {
0202     reco::FormulaEvaluator f("3<=3");
0203 
0204     std::vector<double> emptyV;
0205 
0206     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 1.);
0207   }
0208 
0209   {
0210     reco::FormulaEvaluator f("3>=2");
0211 
0212     std::vector<double> emptyV;
0213 
0214     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 1.);
0215   }
0216   {
0217     reco::FormulaEvaluator f("2>=3");
0218 
0219     std::vector<double> emptyV;
0220 
0221     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 0.);
0222   }
0223   {
0224     reco::FormulaEvaluator f("3>=3");
0225 
0226     std::vector<double> emptyV;
0227 
0228     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 1.);
0229   }
0230 
0231   {
0232     reco::FormulaEvaluator f("3>2");
0233 
0234     std::vector<double> emptyV;
0235 
0236     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 1.);
0237   }
0238   {
0239     reco::FormulaEvaluator f("2>3");
0240 
0241     std::vector<double> emptyV;
0242 
0243     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 0.);
0244   }
0245   {
0246     reco::FormulaEvaluator f("3>3");
0247 
0248     std::vector<double> emptyV;
0249 
0250     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 0.);
0251   }
0252 
0253   {
0254     reco::FormulaEvaluator f("3<2");
0255 
0256     std::vector<double> emptyV;
0257 
0258     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 0.);
0259   }
0260   {
0261     reco::FormulaEvaluator f("2<3");
0262 
0263     std::vector<double> emptyV;
0264 
0265     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 1.);
0266   }
0267   {
0268     reco::FormulaEvaluator f("3<3");
0269 
0270     std::vector<double> emptyV;
0271 
0272     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 0.);
0273   }
0274 
0275   {
0276     reco::FormulaEvaluator f("2==3");
0277 
0278     std::vector<double> emptyV;
0279 
0280     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 0.);
0281   }
0282   {
0283     reco::FormulaEvaluator f("3==3");
0284 
0285     std::vector<double> emptyV;
0286 
0287     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 1.);
0288   }
0289 
0290   {
0291     reco::FormulaEvaluator f("2!=3");
0292 
0293     std::vector<double> emptyV;
0294 
0295     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 1.);
0296   }
0297   {
0298     reco::FormulaEvaluator f("3!=3");
0299 
0300     std::vector<double> emptyV;
0301 
0302     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 0.);
0303   }
0304 
0305   {
0306     reco::FormulaEvaluator f("1+2*3");
0307 
0308     std::vector<double> emptyV;
0309 
0310     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 7.);
0311   }
0312 
0313   {
0314     reco::FormulaEvaluator f("(1+2)*3");
0315 
0316     std::vector<double> emptyV;
0317 
0318     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 9.);
0319   }
0320 
0321   {
0322     reco::FormulaEvaluator f("2*3+1");
0323 
0324     std::vector<double> emptyV;
0325 
0326     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 7.);
0327   }
0328 
0329   {
0330     reco::FormulaEvaluator f("2*(3+1)");
0331 
0332     std::vector<double> emptyV;
0333 
0334     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 8.);
0335   }
0336 
0337   {
0338     reco::FormulaEvaluator f("4/2*3");
0339 
0340     std::vector<double> emptyV;
0341 
0342     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 6.);
0343   }
0344 
0345   {
0346     reco::FormulaEvaluator f("1-2+3");
0347 
0348     std::vector<double> emptyV;
0349 
0350     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 2.);
0351   }
0352 
0353   {
0354     reco::FormulaEvaluator f("(1+2)-(3+4)");
0355     std::vector<double> emptyV;
0356     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == -4.);
0357   }
0358 
0359   {
0360     reco::FormulaEvaluator f("3/2*4+1");
0361 
0362     std::vector<double> emptyV;
0363     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 3. / 2. * 4. + 1);
0364   }
0365 
0366   {
0367     reco::FormulaEvaluator f("1+3/2*4");
0368     std::vector<double> emptyV;
0369     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 1 + 3. / 2. * 4.);
0370   }
0371 
0372   {
0373     reco::FormulaEvaluator f("1+4*(3/2+5)");
0374     std::vector<double> emptyV;
0375     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 1 + 4 * (3. / 2. + 5.));
0376   }
0377 
0378   {
0379     reco::FormulaEvaluator f("1+2*3/4*5");
0380     std::vector<double> emptyV;
0381     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 1 + 2. * 3. / 4. * 5);
0382   }
0383 
0384   {
0385     reco::FormulaEvaluator f("1+2*3/(4+5)+6");
0386     std::vector<double> emptyV;
0387     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 1 + 2. * 3. / (4 + 5) + 6);
0388   }
0389 
0390   {
0391     reco::FormulaEvaluator f("100./3.*2+1");
0392 
0393     std::vector<double> emptyV;
0394 
0395     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 100. / 3. * 2. + 1);
0396   }
0397 
0398   {
0399     reco::FormulaEvaluator f("100./3.*(4-2)+2*(3+1)");
0400 
0401     std::vector<double> emptyV;
0402 
0403     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 100. / 3. * (4 - 2) + 2 * (3 + 1));
0404   }
0405 
0406   {
0407     reco::FormulaEvaluator f("2*(3*4*5)/6");
0408 
0409     std::vector<double> emptyV;
0410 
0411     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 2 * (3 * 4 * 5) / 6);
0412   }
0413 
0414   {
0415     reco::FormulaEvaluator f("2*(2.5*3*3.5)*(4*4.5*5)/6");
0416 
0417     std::vector<double> emptyV;
0418 
0419     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 2 * (2.5 * 3 * 3.5) * (4 * 4.5 * 5) / 6);
0420   }
0421 
0422   {
0423     reco::FormulaEvaluator f("x");
0424 
0425     std::vector<double> emptyV;
0426     std::array<double, 1> v = {{3.}};
0427 
0428     CPPUNIT_ASSERT(f.evaluate(v, emptyV) == 3.);
0429   }
0430 
0431   {
0432     reco::FormulaEvaluator f("y");
0433 
0434     std::vector<double> emptyV;
0435     std::array<double, 2> v = {{0., 3.}};
0436 
0437     CPPUNIT_ASSERT(f.evaluate(v, emptyV) == 3.);
0438   }
0439 
0440   {
0441     reco::FormulaEvaluator f("z");
0442 
0443     std::vector<double> emptyV;
0444     std::array<double, 3> v = {{0., 0., 3.}};
0445 
0446     CPPUNIT_ASSERT(f.evaluate(v, emptyV) == 3.);
0447   }
0448 
0449   {
0450     reco::FormulaEvaluator f("t");
0451 
0452     std::vector<double> emptyV;
0453     std::array<double, 4> v = {{0., 0., 0., 3.}};
0454 
0455     CPPUNIT_ASSERT(f.evaluate(v, emptyV) == 3.);
0456   }
0457 
0458   {
0459     reco::FormulaEvaluator f("[0]");
0460 
0461     std::vector<double> emptyV;
0462     std::array<double, 1> v = {{3.}};
0463 
0464     CPPUNIT_ASSERT(f.evaluate(emptyV, v) == 3.);
0465   }
0466 
0467   {
0468     reco::FormulaEvaluator f("[1]");
0469 
0470     std::vector<double> emptyV;
0471     std::array<double, 2> v = {{0., 3.}};
0472 
0473     CPPUNIT_ASSERT(f.evaluate(emptyV, v) == 3.);
0474   }
0475 
0476   {
0477     reco::FormulaEvaluator f("[0]+[1]*3");
0478 
0479     std::vector<double> emptyV;
0480     std::array<double, 2> v = {{1., 3.}};
0481 
0482     CPPUNIT_ASSERT(f.evaluate(emptyV, v) == 10.);
0483   }
0484 
0485   {
0486     reco::FormulaEvaluator f("log(2)");
0487 
0488     std::vector<double> emptyV;
0489 
0490     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == std::log(2.));
0491   }
0492   {
0493     reco::FormulaEvaluator f("TMath::Log(2)");
0494 
0495     std::vector<double> emptyV;
0496 
0497     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == std::log(2.));
0498   }
0499 
0500   {
0501     reco::FormulaEvaluator f("log10(2)");
0502 
0503     std::vector<double> emptyV;
0504 
0505     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == std::log(2.) / std::log(10.));
0506   }
0507 
0508   {
0509     reco::FormulaEvaluator f("exp(2)");
0510 
0511     std::vector<double> emptyV;
0512 
0513     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == std::exp(2.));
0514   }
0515 
0516   {
0517     reco::FormulaEvaluator f("pow(2,0.3)");
0518 
0519     std::vector<double> emptyV;
0520 
0521     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == std::pow(2., 0.3));
0522   }
0523 
0524   {
0525     reco::FormulaEvaluator f("TMath::Power(2,0.3)");
0526 
0527     std::vector<double> emptyV;
0528 
0529     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == std::pow(2., 0.3));
0530   }
0531 
0532   {
0533     reco::FormulaEvaluator f("TMath::Erf(2.)");
0534 
0535     std::vector<double> emptyV;
0536 
0537     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == TMath::Erf(2.));
0538   }
0539 
0540   {
0541     reco::FormulaEvaluator f("erf(2.)");
0542 
0543     std::vector<double> emptyV;
0544 
0545     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == std::erf(2.));
0546   }
0547 
0548   {
0549     reco::FormulaEvaluator f("TMath::Landau(3.)");
0550 
0551     std::vector<double> emptyV;
0552 
0553     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == TMath::Landau(3.));
0554   }
0555 
0556   {
0557     reco::FormulaEvaluator f("max(2,1)");
0558 
0559     std::vector<double> emptyV;
0560 
0561     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 2);
0562   }
0563 
0564   {
0565     reco::FormulaEvaluator f("max(1,2)");
0566 
0567     std::vector<double> emptyV;
0568 
0569     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 2);
0570   }
0571 
0572   {
0573     reco::FormulaEvaluator f("TMath::Max(2,1)");
0574 
0575     std::vector<double> emptyV;
0576 
0577     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 2);
0578   }
0579 
0580   {
0581     reco::FormulaEvaluator f("TMath::Max(1,2)");
0582 
0583     std::vector<double> emptyV;
0584 
0585     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 2);
0586   }
0587 
0588   {
0589     reco::FormulaEvaluator f("cos(0.5)");
0590 
0591     std::vector<double> emptyV;
0592 
0593     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == std::cos(0.5));
0594   }
0595 
0596   {
0597     reco::FormulaEvaluator f("TMath::Cos(0.5)");
0598 
0599     std::vector<double> emptyV;
0600 
0601     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == TMath::Cos(0.5));
0602   }
0603 
0604   {
0605     reco::FormulaEvaluator f("sin(0.5)");
0606 
0607     std::vector<double> emptyV;
0608 
0609     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == std::sin(0.5));
0610   }
0611 
0612   {
0613     reco::FormulaEvaluator f("TMath::Sin(0.5)");
0614 
0615     std::vector<double> emptyV;
0616 
0617     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == TMath::Sin(0.5));
0618   }
0619 
0620   {
0621     reco::FormulaEvaluator f("tan(0.5)");
0622 
0623     std::vector<double> emptyV;
0624 
0625     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == std::tan(0.5));
0626   }
0627 
0628   {
0629     reco::FormulaEvaluator f("TMath::Tan(0.5)");
0630 
0631     std::vector<double> emptyV;
0632 
0633     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == TMath::Tan(0.5));
0634   }
0635 
0636   {
0637     reco::FormulaEvaluator f("acos(0.5)");
0638 
0639     std::vector<double> emptyV;
0640 
0641     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == std::acos(0.5));
0642   }
0643 
0644   {
0645     reco::FormulaEvaluator f("TMath::ACos(0.5)");
0646 
0647     std::vector<double> emptyV;
0648 
0649     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == TMath::ACos(0.5));
0650   }
0651 
0652   {
0653     reco::FormulaEvaluator f("asin(0.5)");
0654 
0655     std::vector<double> emptyV;
0656 
0657     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == std::asin(0.5));
0658   }
0659 
0660   {
0661     reco::FormulaEvaluator f("TMath::ASin(0.5)");
0662 
0663     std::vector<double> emptyV;
0664 
0665     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == TMath::ASin(0.5));
0666   }
0667 
0668   {
0669     reco::FormulaEvaluator f("atan(0.5)");
0670 
0671     std::vector<double> emptyV;
0672 
0673     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == std::atan(0.5));
0674   }
0675 
0676   {
0677     reco::FormulaEvaluator f("TMath::ATan(0.5)");
0678 
0679     std::vector<double> emptyV;
0680 
0681     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == TMath::ATan(0.5));
0682   }
0683 
0684   {
0685     reco::FormulaEvaluator f("atan2(-0.5, 0.5)");
0686 
0687     std::vector<double> emptyV;
0688 
0689     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == std::atan2(-0.5, 0.5));
0690   }
0691 
0692   {
0693     reco::FormulaEvaluator f("TMath::ATan2(-0.5, 0.5)");
0694 
0695     std::vector<double> emptyV;
0696 
0697     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == TMath::ATan2(-0.5, 0.5));
0698   }
0699 
0700   {
0701     reco::FormulaEvaluator f("cosh(0.5)");
0702 
0703     std::vector<double> emptyV;
0704 
0705     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == std::cosh(0.5));
0706   }
0707 
0708   {
0709     reco::FormulaEvaluator f("TMath::CosH(0.5)");
0710 
0711     std::vector<double> emptyV;
0712 
0713     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == TMath::CosH(0.5));
0714   }
0715 
0716   {
0717     reco::FormulaEvaluator f("sinh(0.5)");
0718 
0719     std::vector<double> emptyV;
0720 
0721     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == std::sinh(0.5));
0722   }
0723 
0724   {
0725     reco::FormulaEvaluator f("TMath::SinH(0.5)");
0726 
0727     std::vector<double> emptyV;
0728 
0729     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == TMath::SinH(0.5));
0730   }
0731 
0732   {
0733     reco::FormulaEvaluator f("tanh(0.5)");
0734 
0735     std::vector<double> emptyV;
0736 
0737     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == std::tanh(0.5));
0738   }
0739 
0740   {
0741     reco::FormulaEvaluator f("TMath::TanH(0.5)");
0742 
0743     std::vector<double> emptyV;
0744 
0745     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == TMath::TanH(0.5));
0746   }
0747 
0748   // std::acosh and std::atanh are using delta fabs instead of equality because gcc compute the value differently at compile time.
0749   {
0750     reco::FormulaEvaluator f("acosh(2.0)");
0751 
0752     std::vector<double> emptyV;
0753 
0754     CPPUNIT_ASSERT(fabs(f.evaluate(emptyV, emptyV) - std::acosh(2.0)) < 1e-9);
0755   }
0756 
0757   {
0758     reco::FormulaEvaluator f("TMath::ACosH(2.0)");
0759 
0760     std::vector<double> emptyV;
0761 
0762     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == TMath::ACosH(2.0));
0763   }
0764 
0765   {
0766     reco::FormulaEvaluator f("asinh(2.0)");
0767 
0768     std::vector<double> emptyV;
0769 
0770     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == std::asinh(2.0));
0771   }
0772 
0773   {
0774     reco::FormulaEvaluator f("TMath::ASinH(2.0)");
0775 
0776     std::vector<double> emptyV;
0777 
0778     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == TMath::ASinH(2.0));
0779   }
0780 
0781   {
0782     reco::FormulaEvaluator f("atanh(0.5)");
0783 
0784     std::vector<double> emptyV;
0785 
0786     CPPUNIT_ASSERT(fabs(f.evaluate(emptyV, emptyV) - std::atanh(0.5)) < 1e-9);
0787   }
0788 
0789   {
0790     reco::FormulaEvaluator f("TMath::ATanH(0.5)");
0791 
0792     std::vector<double> emptyV;
0793 
0794     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == TMath::ATanH(0.5));
0795   }
0796 
0797   {
0798     reco::FormulaEvaluator f("max(max(5,3),2)");
0799 
0800     std::vector<double> emptyV;
0801 
0802     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 5);
0803   }
0804 
0805   {
0806     reco::FormulaEvaluator f("max(2,max(5,3))");
0807 
0808     std::vector<double> emptyV;
0809 
0810     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == 5);
0811   }
0812 
0813   {
0814     reco::FormulaEvaluator f("-(-2.36997+0.413917*TMath::Log(208))/208");
0815     std::vector<double> emptyV;
0816 
0817     CPPUNIT_ASSERT(f.evaluate(emptyV, emptyV) == -(-2.36997 + 0.413917 * std::log(208.)) / 208.);
0818   }
0819 
0820   {
0821     //For Jet energy corrections
0822     reco::FormulaEvaluator f("2*TMath::Erf(4*(x-1))");
0823 
0824     std::vector<double> x = {1.};
0825 
0826     std::vector<double> xValues = {1., 2., 3.};
0827     std::vector<double> emptyV;
0828 
0829     auto func = [](double x) { return 2 * TMath::Erf(4 * (x - 1)); };
0830 
0831     for (auto const xv : xValues) {
0832       x[0] = xv;
0833       CPPUNIT_ASSERT(compare(f.evaluate(x, emptyV), func(x[0])));
0834     }
0835   }
0836 
0837   {
0838     //For Jet energy corrections
0839     reco::FormulaEvaluator f("2*TMath::Landau(2*(x-1))");
0840 
0841     std::vector<double> x = {1.};
0842 
0843     std::vector<double> xValues = {1., 2., 3.};
0844     std::vector<double> emptyV;
0845 
0846     auto func = [](double x) { return 2 * TMath::Landau(2 * (x - 1)); };
0847 
0848     for (auto const xv : xValues) {
0849       x[0] = xv;
0850       CPPUNIT_ASSERT(compare(f.evaluate(x, emptyV), func(x[0])));
0851     }
0852   }
0853 
0854   {
0855     //From SimpleJetCorrector
0856     reco::FormulaEvaluator f("([0]+([1]/((log10(x)^2)+[2])))+([3]*exp(-([4]*((log10(x)-[5])*(log10(x)-[5])))))");
0857 
0858     std::vector<double> x = {1.};
0859 
0860     std::vector<double> v = {1., 4., 2., 0.5, 2., 1.};
0861 
0862     std::vector<double> xValues = {1., 10., 100.};
0863 
0864     auto func = [&v](double x) {
0865       return (v[0] + (v[1] / (((std::log(x) / std::log(10)) * (std::log(x) / std::log(10))) + v[2]))) +
0866              (v[3] *
0867               std::exp(-1. * (v[4] * ((std::log(x) / std::log(10.) - v[5]) * (std::log(x) / std::log(10.) - v[5])))));
0868     };
0869 
0870     for (auto const xv : xValues) {
0871       x[0] = xv;
0872       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
0873     }
0874   }
0875 
0876   {
0877     //From SimpleJetCorrector
0878     reco::FormulaEvaluator f("[0]*([1]+[2]*TMath::Log(x))");
0879 
0880     std::vector<double> x = {1.};
0881 
0882     std::vector<double> v = {1.3, 4., 2.};
0883 
0884     std::vector<double> xValues = {1., 10., 100.};
0885 
0886     auto func = [&v](double x) { return v[0] * (v[1] + v[2] * std::log(x)); };
0887 
0888     for (auto const xv : xValues) {
0889       x[0] = xv;
0890       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
0891     }
0892   }
0893 
0894   {
0895     //From SimpleJetCorrector
0896     reco::FormulaEvaluator f("[0]+([1]/((log10(x)^[2])+[3]))");
0897 
0898     std::vector<double> x = {1.};
0899 
0900     std::vector<double> v = {1.3, 4., 1.7, 1.};
0901 
0902     std::vector<double> xValues = {1., 10., 100.};
0903 
0904     auto func = [&v](double x) { return v[0] + (v[1] / ((std::pow(log(x) / log(10.), v[2])) + v[3])); };
0905 
0906     for (auto const xv : xValues) {
0907       x[0] = xv;
0908       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
0909     }
0910   }
0911 
0912   {
0913     //From SimpleJetCorrector
0914     reco::FormulaEvaluator f("max(0.0001,1-y*([0]+([1]*z)*(1+[2]*log(x)))/x)");
0915 
0916     std::vector<double> v = {.1, 1., .5};
0917 
0918     std::vector<double> p = {1.3, 5., 10.};
0919 
0920     std::vector<double> xValues = {1., 10., 100.};
0921 
0922     auto func = [&p](double x, double y, double z) {
0923       return std::max(0.0001, 1 - y * (p[0] + (p[1] * z) * (1 + p[2] * std::log(x))) / x);
0924     };
0925 
0926     for (auto const xv : xValues) {
0927       v[0] = xv;
0928       CPPUNIT_ASSERT(compare(f.evaluate(v, p), func(v[0], v[1], v[2])));
0929     }
0930   }
0931   {
0932     reco::FormulaEvaluator f("(-2.36997+0.413917*TMath::Log(x))/x-(-2.36997+0.413917*TMath::Log(208))/208");
0933 
0934     std::vector<double> x = {1.};
0935 
0936     std::vector<double> v;
0937 
0938     auto func = [](double x) {
0939       return (-2.36997 + 0.413917 * std::log(x)) / x - (-2.36997 + 0.413917 * std::log(208)) / 208;
0940     };
0941 
0942     std::vector<double> xValues = {.1, 1., 10., 100.};
0943     for (auto const xv : xValues) {
0944       x[0] = xv;
0945       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
0946     }
0947   }
0948 
0949   {
0950     reco::FormulaEvaluator f(
0951         "TMath::Max(0.,1.03091-0.051154*pow(x,-0.154227))-TMath::Max(0.,1.03091-0.051154*TMath::Power(208.,-0.154227)"
0952         ")");
0953     std::vector<double> x = {1.};
0954 
0955     std::vector<double> v;
0956 
0957     std::vector<double> xValues = {.1, 1., 10., 100.};
0958 
0959     auto func = [](double x) {
0960       return std::max(0., 1.03091 - 0.051154 * std::pow(x, -0.154227)) -
0961              std::max(0., 1.03091 - 0.051154 * std::pow(208., -0.154227));
0962     };
0963 
0964     for (auto const xv : xValues) {
0965       x[0] = xv;
0966       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
0967     }
0968   }
0969 
0970   {
0971     reco::FormulaEvaluator f("[2]*([3]+[4]*TMath::Log(max([0],min([1],x))))");
0972 
0973     std::vector<double> x = {1.};
0974 
0975     std::vector<double> v = {1., 4., 2., 0.5, 2., 1., 1., -1.};
0976     std::vector<double> xValues = {.1, 1., 10., 100.};
0977 
0978     auto func = [&v](double x) { return v[2] * (v[3] + v[4] * std::log(std::max(v[0], std::min(v[1], x)))); };
0979     for (auto const xv : xValues) {
0980       x[0] = xv;
0981       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
0982     }
0983   }
0984   {
0985     //From SimpleJetCorrector
0986     reco::FormulaEvaluator f(
0987         "((x>=[6])*(([0]+([1]/((log10(x)^2)+[2])))+([3]*exp(-([4]*((log10(x)-[5])*(log10(x)-[5])))))))+((x<[6])*[7])");
0988 
0989     std::vector<double> x = {1.};
0990 
0991     std::vector<double> v = {1., 4., 2., 0.5, 2., 1., 1., -1.};
0992 
0993     std::vector<double> xValues = {.1, 1., 10., 100.};
0994 
0995     auto func = [&v](double x) {
0996       return ((x >= v[6]) * ((v[0] + (v[1] / (((std::log(x) / std::log(10)) * (std::log(x) / std::log(10))) + v[2]))) +
0997                              (v[3] * std::exp(-1. * (v[4] * ((std::log(x) / std::log(10.) - v[5]) *
0998                                                              (std::log(x) / std::log(10.) - v[5]))))))) +
0999              ((x < v[6]) * v[7]);
1000     };
1001 
1002     for (auto const xv : xValues) {
1003       x[0] = xv;
1004       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1005     }
1006   }
1007 
1008   {
1009     reco::FormulaEvaluator f(
1010         "(TMath::Max(0.,1.03091-0.051154*pow(x,-0.154227))-TMath::Max(0.,1.03091-0.051154*TMath::Power(208.,-0.154227))"
1011         ")+[7]*((-2.36997+0.413917*TMath::Log(x))/x-(-2.36997+0.413917*TMath::Log(208))/208)");
1012 
1013     std::vector<double> x = {1.};
1014 
1015     std::vector<double> v = {1., 4., 2., 0.5, 2., 1., 1., -1.};
1016     std::vector<double> xValues = {.1, 1., 10., 100.};
1017 
1018     auto func = [&v](double x) {
1019       return (std::max(0., 1.03091 - 0.051154 * std::pow(x, -0.154227)) -
1020               std::max(0., 1.03091 - 0.051154 * std::pow(208., -0.154227))) +
1021              v[7] * ((-2.36997 + 0.413917 * std::log(x)) / x - (-2.36997 + 0.413917 * std::log(208)) / 208);
1022     };
1023 
1024     for (auto const xv : xValues) {
1025       x[0] = xv;
1026       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1027     }
1028   }
1029 
1030   {
1031     //From SimpleJetCorrector
1032     reco::FormulaEvaluator f("100./3.*0.154227+2.36997");
1033 
1034     std::vector<double> x = {1.};
1035 
1036     std::vector<double> v = {1., 4., 2., 0.5, 2., 1., 1., -1.};
1037     std::vector<double> xValues = {.1, 1., 10., 100.};
1038 
1039     auto func = [](double x) { return 100. / 3. * 0.154227 + 2.36997; };
1040 
1041     for (auto const xv : xValues) {
1042       x[0] = xv;
1043       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1044     }
1045   }
1046 
1047   {
1048     //From SimpleJetCorrector
1049     reco::FormulaEvaluator f(
1050         "[2]*([3]+[4]*TMath::Log(max([0],min([1],x))))*1./([5]+[6]*100./"
1051         "3.*(TMath::Max(0.,1.03091-0.051154*pow(x,-0.154227))-TMath::Max(0.,1.03091-0.051154*TMath::Power(208.,-0."
1052         "154227)))+[7]*((-2.36997+0.413917*TMath::Log(x))/x-(-2.36997+0.413917*TMath::Log(208))/208))");
1053 
1054     std::vector<double> x = {1.};
1055 
1056     std::vector<double> v = {1., 4., 2., 0.5, 2., 1., 1., -1.};
1057     std::vector<double> xValues = {.1, 1., 10., 100.};
1058 
1059     auto func = [&v](double x) {
1060       return v[2] * (v[3] + v[4] * std::log(std::max(v[0], std::min(v[1], x)))) * 1. /
1061              (v[5] +
1062               v[6] * 100. / 3. *
1063                   (std::max(0., 1.03091 - 0.051154 * std::pow(x, -0.154227)) -
1064                    std::max(0., 1.03091 - 0.051154 * std::pow(208., -0.154227))) +
1065               v[7] * ((-2.36997 + 0.413917 * std::log(x)) / x - (-2.36997 + 0.413917 * std::log(208)) / 208));
1066     };
1067 
1068     for (auto const xv : xValues) {
1069       x[0] = xv;
1070       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1071     }
1072   }
1073 
1074   {
1075     //Tests that pick proper evaluator for argument of function
1076     reco::FormulaEvaluator f("exp([4]*(log10(x)-[5])*(log10(x)-[5]))");
1077     std::vector<double> x = {10.};
1078 
1079     std::vector<double> v = {0.88524, 28.4947, 4.89135, -19.0245, 0.0227809, -6.97308};
1080     std::vector<double> xValues = {10.};
1081 
1082     auto func = [&v](double x) {
1083       return std::exp(v[4] * (std::log(x) / std::log(10) - v[5]) * (std::log(x) / std::log(10) - v[5]));
1084     };
1085 
1086     for (auto const xv : xValues) {
1087       x[0] = xv;
1088       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1089     }
1090   }
1091 
1092   {
1093     reco::FormulaEvaluator f(
1094         "max(0.0001,[0]+[1]/(pow(log10(x),2)+[2])+[3]*exp(-1*([4]*((log10(x)-[5])*(log10(x)-[5])))))");
1095 
1096     std::vector<double> x = {10.};
1097 
1098     std::vector<double> v = {0.88524, 28.4947, 4.89135, -19.0245, 0.0227809, -6.97308};
1099     std::vector<double> xValues = {10.};
1100 
1101     auto func = [&v](double x) {
1102       return std::max(
1103           0.0001,
1104           v[0] + v[1] / (std::pow(std::log(x) / std::log(10), 2) + v[2]) +
1105               v[3] * std::exp(-1 * v[4] * (std::log(x) / std::log(10) - v[5]) * (std::log(x) / std::log(10) - v[5])));
1106     };
1107 
1108     for (auto const xv : xValues) {
1109       x[0] = xv;
1110       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1111     }
1112   }
1113 
1114   {
1115     std::vector<double> x = {425.92155818};
1116     std::vector<double> v = {0.945459, 2.78658, 1.65054, -48.1061, 0.0287239, -10.8759};
1117     std::vector<double> xValues = {425.92155818};
1118 
1119     reco::FormulaEvaluator f("-[4]*(log10(x)-[5])*(log10(x)-[5])");
1120     auto func = [&v](double x) { return -v[4] * (std::log10(x) - v[5]) * (std::log10(x) - v[5]); };
1121 
1122     for (auto const xv : xValues) {
1123       x[0] = xv;
1124       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1125     }
1126   }
1127 
1128   {
1129     reco::FormulaEvaluator f("max(0.0001,[0]+[1]/(pow(log10(x),2)+[2])+[3]*exp(-[4]*(log10(x)-[5])*(log10(x)-[5])))");
1130 
1131     std::vector<double> x = {10.};
1132 
1133     std::vector<double> v = {0.88524, 28.4947, 4.89135, -19.0245, 0.0227809, -6.97308};
1134     std::vector<double> xValues = {10.};
1135 
1136     auto func = [&v](double x) {
1137       return std::max(
1138           0.0001,
1139           v[0] + v[1] / (std::pow(std::log(x) / std::log(10), 2) + v[2]) +
1140               v[3] * std::exp(-v[4] * (std::log(x) / std::log(10) - v[5]) * (std::log(x) / std::log(10) - v[5])));
1141     };
1142 
1143     for (auto const xv : xValues) {
1144       x[0] = xv;
1145       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1146     }
1147   }
1148 
1149   {
1150     reco::FormulaEvaluator f(
1151         "[2]*([3]*([4]+[5]*TMath::Log(max([0],min([1],x))))*1./([6]+[7]*100./"
1152         "3.*(TMath::Max(0.,1.03091-0.051154*pow(x,-0.154227))-TMath::Max(0.,1.03091-0.051154*TMath::Power(208.,-0."
1153         "154227)))+[8]*((1+0.04432-1.304*pow(max(30.,min(6500.,x)),-0.4624)+(0+1.724*TMath::Log(max(30.,min(6500.,x))))"
1154         "/max(30.,min(6500.,x)))-(1+0.04432-1.304*pow(208.,-0.4624)+(0+1.724*TMath::Log(208.))/208.))))");
1155 
1156     std::vector<double> v = {55, 2510, 0.997756, 1.000155, 0.979016, 0.001834, 0.982, -0.048, 1.250};
1157 
1158     std::vector<double> x = {100};
1159 
1160     auto func = [&v](double x) {
1161       return v[2] *
1162              (v[3] * (v[4] + v[5] * TMath::Log(std::max(v[0], std::min(v[1], x)))) * 1. /
1163               (v[6] +
1164                v[7] * 100. / 3. *
1165                    (TMath::Max(0., 1.03091 - 0.051154 * std::pow(x, -0.154227)) -
1166                     TMath::Max(0., 1.03091 - 0.051154 * TMath::Power(208., -0.154227))) +
1167                v[8] *
1168                    ((1 + 0.04432 - 1.304 * std::pow(std::max(30., std::min(6500., x)), -0.4624) +
1169                      (0 + 1.724 * TMath::Log(std::max(30., std::min(6500., x)))) / std::max(30., std::min(6500., x))) -
1170                     (1 + 0.04432 - 1.304 * std::pow(208., -0.4624) + (0 + 1.724 * TMath::Log(208.)) / 208.))));
1171     };
1172 
1173     CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1174   }
1175 
1176   {
1177     std::vector<std::pair<std::string, double>> formulas = {
1178         {"(1+0.04432+(1.724+100.))-1", 101.76832},
1179         {"(1+(1.724+100.)+0.04432)-1", 101.76832},
1180         {"((1.724+100.)+1+0.04432)-1", 101.76832},
1181         {"(1+0.04432+1.724/100.)-1", .06156},
1182         {"(1+1.724/100.+0.04432)-1", .06156},
1183         {"(1.724/100.+1+0.04432)-1", .06156},
1184         {"(1+0.04432+(1.724/100.))-1", (1 + 0.04432 + (1.724 / 100.)) - 1},
1185         {"(1+(1.724/100.)+0.04432)-1", (1 + 0.04432 + (1.724 / 100.)) - 1},
1186         {"((1.724/100.)+1+0.04432)-1", (1 + 0.04432 + (1.724 / 100.)) - 1},
1187         {"0.997756*(1.000155*(0.979016+0.001834*TMath::Log(max(55.,min(2510.,100.))))*1./(0.982+-0.048*100./"
1188          "3.*(TMath::Max(0.,1.03091-0.051154*pow(100.,-0.154227))-TMath::Max(0.,1.03091-0.051154*TMath::Power(208.,-0."
1189          "154227)))+1.250*((1+0.04432-1.304*pow(max(30.,min(6500.,100.)),-0.4624)+(0+1.724*TMath::Log(max(30.,min(6500."
1190          ",100.))))/max(30.,min(6500.,100.)))-(1+0.04432-1.304*pow(208.,-0.4624)+(0+1.724*TMath::Log(208.))/208.))))",
1191          0.997756 *
1192              (1.000155 * (0.979016 + 0.001834 * TMath::Log(std::max(55., std::min(2510., 100.)))) * 1. /
1193               (0.982 +
1194                -0.048 * 100. / 3. *
1195                    (TMath::Max(0., 1.03091 - 0.051154 * std::pow(100., -0.154227)) -
1196                     TMath::Max(0., 1.03091 - 0.051154 * TMath::Power(208., -0.154227))) +
1197                1.250 * ((1 + 0.04432 - 1.304 * std::pow(std::max(30., std::min(6500., 100.)), -0.4624) +
1198                          (0 + 1.724 * TMath::Log(std::max(30., std::min(6500., 100.)))) /
1199                              std::max(30., std::min(6500., 100.))) -
1200                         (1 + 0.04432 - 1.304 * std::pow(208., -0.4624) + (0 + 1.724 * TMath::Log(208.)) / 208.))))}};
1201 
1202     std::vector<double> x = {};
1203     std::vector<double> v = {};
1204     for (auto const& form_val : formulas) {
1205       reco::FormulaEvaluator f(form_val.first);
1206 
1207       CPPUNIT_ASSERT(compare(f.evaluate(x, v), form_val.second));
1208     }
1209   }
1210 
1211   //tests for JER
1212   {
1213     reco::FormulaEvaluator f("[0]+[1]*exp(-x/[2])");
1214 
1215     std::vector<double> v = {0.006467, 0.02519, 77.08};
1216     std::vector<double> x = {100.};
1217 
1218     auto func = [&v](double x) { return v[0] + v[1] * std::exp(-x / v[2]); };
1219 
1220     CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1221   }
1222 
1223   {
1224     reco::FormulaEvaluator f("max(0.0001,1-y/x*([1]*(z-[0])*(1+[2]*log(x/30.))))");
1225 
1226     std::vector<double> v = {1.4, 0.453645, -0.015665};
1227     std::vector<double> x = {157.2, 0.5, 23.2};
1228 
1229     auto func = [&v](double x, double y, double z) {
1230       return std::max(0.0001, 1 - y / x * (v[1] * (z - v[0]) * (1 + v[2] * std::log(x / 30.))));
1231     };
1232 
1233     CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0], x[1], x[2])));
1234   }
1235   {
1236     reco::FormulaEvaluator f("max(0.0001,1-y*[1]*(z-[0])*(1+[2]*log(x/30.))/x)");
1237 
1238     std::vector<double> v = {1.4, 0.453645, -0.015665};
1239     std::vector<double> x = {157.2, 0.5, 23.2};
1240 
1241     auto func = [&v](double x, double y, double z) {
1242       return std::max(0.0001, 1 - y / x * (v[1] * (z - v[0]) * (1 + v[2] * std::log(x / 30.))));
1243     };
1244 
1245     CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0], x[1], x[2])));
1246   }
1247   {
1248     reco::FormulaEvaluator f("max(0.0001,1-y*([1]*(z-[0])*(1+[2]*log(x/30.)))/x)");
1249 
1250     std::vector<double> v = {1.4, 0.453645, -0.015665};
1251     std::vector<double> x = {157.2, 0.5, 23.2};
1252 
1253     auto func = [&v](double x, double y, double z) {
1254       return std::max(0.0001, 1 - y * (v[1] * (z - v[0]) * (1 + v[2] * std::log(x / 30.))) / x);
1255     };
1256 
1257     CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0], x[1], x[2])));
1258   }
1259 
1260   {
1261     reco::FormulaEvaluator f("sqrt([0]*abs([0])/(x*x)+[1]*[1]*pow(x,[3])+[2]*[2])");
1262 
1263     std::vector<double> v = {1.326, 0.4209, 0.02223, -0.6704};
1264     std::vector<double> x = {100.};
1265 
1266     auto func = [&v](double x) {
1267       return std::sqrt(v[0] * std::abs(v[0]) / (x * x) + v[1] * v[1] * std::pow(x, v[3]) + v[2] * v[2]);
1268     };
1269 
1270     CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1271   }
1272 
1273   {
1274     reco::FormulaEvaluator f("sqrt([0]*[0]/(x*x)+[1]*[1]/x+[2]*[2])");
1275 
1276     std::vector<double> v = {2.3, 0.20, 0.009};
1277     std::vector<double> x = {100.};
1278 
1279     auto func = [&v](double x) { return std::sqrt(v[0] * v[0] / (x * x) + v[1] * v[1] / x + v[2] * v[2]); };
1280 
1281     CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1282   }
1283 
1284   {
1285     //This was previously evaluated wrong
1286     std::array<double, 3> xs = {{224., 225., 226.}};
1287     std::vector<double> x = {0.};
1288     std::vector<double> v = {-3., -3., -3., -3., -3., -3.};
1289 
1290     reco::FormulaEvaluator f(
1291         "([0]+[1]*x+[2]*x^2)*(x<225)+([0]+[1]*225+[2]*225^2+[3]*(x-225)+[4]*(x-225)^2+[5]*(x-225)^3)*(x>225)");
1292 
1293     auto func = [&v](double x) {
1294       return (v[0] + v[1] * x + v[2] * (x * x)) * (x < 225) +
1295              (v[0] + v[1] * 225 + v[2] * (225 * 225) + v[3] * (x - 225) + v[4] * ((x - 225) * (x - 225)) +
1296               v[5] * ((x - 225) * (x - 225) * (x - 225))) *
1297                  (x > 225);
1298     };
1299 
1300     for (auto x_i : xs) {
1301       x[0] = x_i;
1302       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1303     }
1304   }
1305 
1306   {
1307     auto t = []() { reco::FormulaEvaluator f("doesNotExist(2)"); };
1308 
1309     CPPUNIT_ASSERT_THROW(t(), cms::Exception);
1310   }
1311 
1312   {
1313     auto t = []() { reco::FormulaEvaluator f("doesNotExist(2) + abs(-1)"); };
1314 
1315     CPPUNIT_ASSERT_THROW(t(), cms::Exception);
1316   }
1317 
1318   {
1319     auto t = []() { reco::FormulaEvaluator f("abs(-1) + doesNotExist(2)"); };
1320 
1321     CPPUNIT_ASSERT_THROW(t(), cms::Exception);
1322   }
1323 
1324   {
1325     auto t = []() { reco::FormulaEvaluator f("abs(-1) + ( 5 * doesNotExist(2))"); };
1326 
1327     CPPUNIT_ASSERT_THROW(t(), cms::Exception);
1328   }
1329 
1330   {
1331     auto t = []() { reco::FormulaEvaluator f("( 5 * doesNotExist(2)) + abs(-1)"); };
1332 
1333     CPPUNIT_ASSERT_THROW(t(), cms::Exception);
1334   }
1335 
1336   {
1337     auto t = []() { reco::FormulaEvaluator f("TMath::Exp(2)"); };
1338 
1339     CPPUNIT_ASSERT_THROW(t(), cms::Exception);
1340   }
1341 
1342   {
1343     //this was previously causing a seg fault
1344     auto t = []() { reco::FormulaEvaluator f("1 + 2 * 3 + 5 * doesNotExist(2) "); };
1345 
1346     CPPUNIT_ASSERT_THROW(t(), cms::Exception);
1347   }
1348 
1349   {
1350     //Make sure spaces are shown in exception message
1351     try {
1352       reco::FormulaEvaluator f("1 + 2#");
1353     } catch (cms::Exception const& e) {
1354       auto r =
1355           "An exception of category 'FormulaEvaluatorParseError' occurred.\n"
1356           "Exception Message:\n"
1357           "While parsing '1 + 2#' could not parse beyond '1 + 2'\n";
1358       CPPUNIT_ASSERT(std::string(r) == e.what());
1359     }
1360   }
1361 }