Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-10-25 09:36:19

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     auto value = f.evaluate(emptyV, emptyV);
0818     CPPUNIT_ASSERT(std::abs(value - (-(-2.36997 + 0.413917 * std::log(208.)) / 208.)) / value < 5.0E-16);
0819   }
0820 
0821   {
0822     //For Jet energy corrections
0823     reco::FormulaEvaluator f("2*TMath::Erf(4*(x-1))");
0824 
0825     std::vector<double> x = {1.};
0826 
0827     std::vector<double> xValues = {1., 2., 3.};
0828     std::vector<double> emptyV;
0829 
0830     auto func = [](double x) { return 2 * TMath::Erf(4 * (x - 1)); };
0831 
0832     for (auto const xv : xValues) {
0833       x[0] = xv;
0834       CPPUNIT_ASSERT(compare(f.evaluate(x, emptyV), func(x[0])));
0835     }
0836   }
0837 
0838   {
0839     //For Jet energy corrections
0840     reco::FormulaEvaluator f("2*TMath::Landau(2*(x-1))");
0841 
0842     std::vector<double> x = {1.};
0843 
0844     std::vector<double> xValues = {1., 2., 3.};
0845     std::vector<double> emptyV;
0846 
0847     auto func = [](double x) { return 2 * TMath::Landau(2 * (x - 1)); };
0848 
0849     for (auto const xv : xValues) {
0850       x[0] = xv;
0851       CPPUNIT_ASSERT(compare(f.evaluate(x, emptyV), func(x[0])));
0852     }
0853   }
0854 
0855   {
0856     //From SimpleJetCorrector
0857     reco::FormulaEvaluator f("([0]+([1]/((log10(x)^2)+[2])))+([3]*exp(-([4]*((log10(x)-[5])*(log10(x)-[5])))))");
0858 
0859     std::vector<double> x = {1.};
0860 
0861     std::vector<double> v = {1., 4., 2., 0.5, 2., 1.};
0862 
0863     std::vector<double> xValues = {1., 10., 100.};
0864 
0865     auto func = [&v](double x) {
0866       return (v[0] + (v[1] / (((std::log(x) / std::log(10)) * (std::log(x) / std::log(10))) + v[2]))) +
0867              (v[3] *
0868               std::exp(-1. * (v[4] * ((std::log(x) / std::log(10.) - v[5]) * (std::log(x) / std::log(10.) - v[5])))));
0869     };
0870 
0871     for (auto const xv : xValues) {
0872       x[0] = xv;
0873       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
0874     }
0875   }
0876 
0877   {
0878     //From SimpleJetCorrector
0879     reco::FormulaEvaluator f("[0]*([1]+[2]*TMath::Log(x))");
0880 
0881     std::vector<double> x = {1.};
0882 
0883     std::vector<double> v = {1.3, 4., 2.};
0884 
0885     std::vector<double> xValues = {1., 10., 100.};
0886 
0887     auto func = [&v](double x) { return v[0] * (v[1] + v[2] * std::log(x)); };
0888 
0889     for (auto const xv : xValues) {
0890       x[0] = xv;
0891       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
0892     }
0893   }
0894 
0895   {
0896     //From SimpleJetCorrector
0897     reco::FormulaEvaluator f("[0]+([1]/((log10(x)^[2])+[3]))");
0898 
0899     std::vector<double> x = {1.};
0900 
0901     std::vector<double> v = {1.3, 4., 1.7, 1.};
0902 
0903     std::vector<double> xValues = {1., 10., 100.};
0904 
0905     auto func = [&v](double x) { return v[0] + (v[1] / ((std::pow(log(x) / log(10.), v[2])) + v[3])); };
0906 
0907     for (auto const xv : xValues) {
0908       x[0] = xv;
0909       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
0910     }
0911   }
0912 
0913   {
0914     //From SimpleJetCorrector
0915     reco::FormulaEvaluator f("max(0.0001,1-y*([0]+([1]*z)*(1+[2]*log(x)))/x)");
0916 
0917     std::vector<double> v = {.1, 1., .5};
0918 
0919     std::vector<double> p = {1.3, 5., 10.};
0920 
0921     std::vector<double> xValues = {1., 10., 100.};
0922 
0923     auto func = [&p](double x, double y, double z) {
0924       return std::max(0.0001, 1 - y * (p[0] + (p[1] * z) * (1 + p[2] * std::log(x))) / x);
0925     };
0926 
0927     for (auto const xv : xValues) {
0928       v[0] = xv;
0929       CPPUNIT_ASSERT(compare(f.evaluate(v, p), func(v[0], v[1], v[2])));
0930     }
0931   }
0932   {
0933     reco::FormulaEvaluator f("(-2.36997+0.413917*TMath::Log(x))/x-(-2.36997+0.413917*TMath::Log(208))/208");
0934 
0935     std::vector<double> x = {1.};
0936 
0937     std::vector<double> v;
0938 
0939     auto func = [](double x) {
0940       return (-2.36997 + 0.413917 * std::log(x)) / x - (-2.36997 + 0.413917 * std::log(208)) / 208;
0941     };
0942 
0943     std::vector<double> xValues = {.1, 1., 10., 100.};
0944     for (auto const xv : xValues) {
0945       x[0] = xv;
0946       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
0947     }
0948   }
0949 
0950   {
0951     reco::FormulaEvaluator f(
0952         "TMath::Max(0.,1.03091-0.051154*pow(x,-0.154227))-TMath::Max(0.,1.03091-0.051154*TMath::Power(208.,-0.154227)"
0953         ")");
0954     std::vector<double> x = {1.};
0955 
0956     std::vector<double> v;
0957 
0958     std::vector<double> xValues = {.1, 1., 10., 100.};
0959 
0960     auto func = [](double x) {
0961       return std::max(0., 1.03091 - 0.051154 * std::pow(x, -0.154227)) -
0962              std::max(0., 1.03091 - 0.051154 * std::pow(208., -0.154227));
0963     };
0964 
0965     for (auto const xv : xValues) {
0966       x[0] = xv;
0967       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
0968     }
0969   }
0970 
0971   {
0972     reco::FormulaEvaluator f("[2]*([3]+[4]*TMath::Log(max([0],min([1],x))))");
0973 
0974     std::vector<double> x = {1.};
0975 
0976     std::vector<double> v = {1., 4., 2., 0.5, 2., 1., 1., -1.};
0977     std::vector<double> xValues = {.1, 1., 10., 100.};
0978 
0979     auto func = [&v](double x) { return v[2] * (v[3] + v[4] * std::log(std::max(v[0], std::min(v[1], x)))); };
0980     for (auto const xv : xValues) {
0981       x[0] = xv;
0982       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
0983     }
0984   }
0985   {
0986     //From SimpleJetCorrector
0987     reco::FormulaEvaluator f(
0988         "((x>=[6])*(([0]+([1]/((log10(x)^2)+[2])))+([3]*exp(-([4]*((log10(x)-[5])*(log10(x)-[5])))))))+((x<[6])*[7])");
0989 
0990     std::vector<double> x = {1.};
0991 
0992     std::vector<double> v = {1., 4., 2., 0.5, 2., 1., 1., -1.};
0993 
0994     std::vector<double> xValues = {.1, 1., 10., 100.};
0995 
0996     auto func = [&v](double x) {
0997       return ((x >= v[6]) * ((v[0] + (v[1] / (((std::log(x) / std::log(10)) * (std::log(x) / std::log(10))) + v[2]))) +
0998                              (v[3] * std::exp(-1. * (v[4] * ((std::log(x) / std::log(10.) - v[5]) *
0999                                                              (std::log(x) / std::log(10.) - v[5]))))))) +
1000              ((x < v[6]) * v[7]);
1001     };
1002 
1003     for (auto const xv : xValues) {
1004       x[0] = xv;
1005       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1006     }
1007   }
1008 
1009   {
1010     reco::FormulaEvaluator f(
1011         "(TMath::Max(0.,1.03091-0.051154*pow(x,-0.154227))-TMath::Max(0.,1.03091-0.051154*TMath::Power(208.,-0.154227))"
1012         ")+[7]*((-2.36997+0.413917*TMath::Log(x))/x-(-2.36997+0.413917*TMath::Log(208))/208)");
1013 
1014     std::vector<double> x = {1.};
1015 
1016     std::vector<double> v = {1., 4., 2., 0.5, 2., 1., 1., -1.};
1017     std::vector<double> xValues = {.1, 1., 10., 100.};
1018 
1019     auto func = [&v](double x) {
1020       return (std::max(0., 1.03091 - 0.051154 * std::pow(x, -0.154227)) -
1021               std::max(0., 1.03091 - 0.051154 * std::pow(208., -0.154227))) +
1022              v[7] * ((-2.36997 + 0.413917 * std::log(x)) / x - (-2.36997 + 0.413917 * std::log(208)) / 208);
1023     };
1024 
1025     for (auto const xv : xValues) {
1026       x[0] = xv;
1027       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1028     }
1029   }
1030 
1031   {
1032     //From SimpleJetCorrector
1033     reco::FormulaEvaluator f("100./3.*0.154227+2.36997");
1034 
1035     std::vector<double> x = {1.};
1036 
1037     std::vector<double> v = {1., 4., 2., 0.5, 2., 1., 1., -1.};
1038     std::vector<double> xValues = {.1, 1., 10., 100.};
1039 
1040     auto func = [](double x) { return 100. / 3. * 0.154227 + 2.36997; };
1041 
1042     for (auto const xv : xValues) {
1043       x[0] = xv;
1044       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1045     }
1046   }
1047 
1048   {
1049     //From SimpleJetCorrector
1050     reco::FormulaEvaluator f(
1051         "[2]*([3]+[4]*TMath::Log(max([0],min([1],x))))*1./([5]+[6]*100./"
1052         "3.*(TMath::Max(0.,1.03091-0.051154*pow(x,-0.154227))-TMath::Max(0.,1.03091-0.051154*TMath::Power(208.,-0."
1053         "154227)))+[7]*((-2.36997+0.413917*TMath::Log(x))/x-(-2.36997+0.413917*TMath::Log(208))/208))");
1054 
1055     std::vector<double> x = {1.};
1056 
1057     std::vector<double> v = {1., 4., 2., 0.5, 2., 1., 1., -1.};
1058     std::vector<double> xValues = {.1, 1., 10., 100.};
1059 
1060     auto func = [&v](double x) {
1061       return v[2] * (v[3] + v[4] * std::log(std::max(v[0], std::min(v[1], x)))) * 1. /
1062              (v[5] +
1063               v[6] * 100. / 3. *
1064                   (std::max(0., 1.03091 - 0.051154 * std::pow(x, -0.154227)) -
1065                    std::max(0., 1.03091 - 0.051154 * std::pow(208., -0.154227))) +
1066               v[7] * ((-2.36997 + 0.413917 * std::log(x)) / x - (-2.36997 + 0.413917 * std::log(208)) / 208));
1067     };
1068 
1069     for (auto const xv : xValues) {
1070       x[0] = xv;
1071       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1072     }
1073   }
1074 
1075   {
1076     //Tests that pick proper evaluator for argument of function
1077     reco::FormulaEvaluator f("exp([4]*(log10(x)-[5])*(log10(x)-[5]))");
1078     std::vector<double> x = {10.};
1079 
1080     std::vector<double> v = {0.88524, 28.4947, 4.89135, -19.0245, 0.0227809, -6.97308};
1081     std::vector<double> xValues = {10.};
1082 
1083     auto func = [&v](double x) {
1084       return std::exp(v[4] * (std::log(x) / std::log(10) - v[5]) * (std::log(x) / std::log(10) - v[5]));
1085     };
1086 
1087     for (auto const xv : xValues) {
1088       x[0] = xv;
1089       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1090     }
1091   }
1092 
1093   {
1094     reco::FormulaEvaluator f(
1095         "max(0.0001,[0]+[1]/(pow(log10(x),2)+[2])+[3]*exp(-1*([4]*((log10(x)-[5])*(log10(x)-[5])))))");
1096 
1097     std::vector<double> x = {10.};
1098 
1099     std::vector<double> v = {0.88524, 28.4947, 4.89135, -19.0245, 0.0227809, -6.97308};
1100     std::vector<double> xValues = {10.};
1101 
1102     auto func = [&v](double x) {
1103       return std::max(
1104           0.0001,
1105           v[0] + v[1] / (std::pow(std::log(x) / std::log(10), 2) + v[2]) +
1106               v[3] * std::exp(-1 * v[4] * (std::log(x) / std::log(10) - v[5]) * (std::log(x) / std::log(10) - v[5])));
1107     };
1108 
1109     for (auto const xv : xValues) {
1110       x[0] = xv;
1111       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1112     }
1113   }
1114 
1115   {
1116     std::vector<double> x = {425.92155818};
1117     std::vector<double> v = {0.945459, 2.78658, 1.65054, -48.1061, 0.0287239, -10.8759};
1118     std::vector<double> xValues = {425.92155818};
1119 
1120     reco::FormulaEvaluator f("-[4]*(log10(x)-[5])*(log10(x)-[5])");
1121     auto func = [&v](double x) { return -v[4] * (std::log10(x) - v[5]) * (std::log10(x) - v[5]); };
1122 
1123     for (auto const xv : xValues) {
1124       x[0] = xv;
1125       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1126     }
1127   }
1128 
1129   {
1130     reco::FormulaEvaluator f("max(0.0001,[0]+[1]/(pow(log10(x),2)+[2])+[3]*exp(-[4]*(log10(x)-[5])*(log10(x)-[5])))");
1131 
1132     std::vector<double> x = {10.};
1133 
1134     std::vector<double> v = {0.88524, 28.4947, 4.89135, -19.0245, 0.0227809, -6.97308};
1135     std::vector<double> xValues = {10.};
1136 
1137     auto func = [&v](double x) {
1138       return std::max(
1139           0.0001,
1140           v[0] + v[1] / (std::pow(std::log(x) / std::log(10), 2) + v[2]) +
1141               v[3] * std::exp(-v[4] * (std::log(x) / std::log(10) - v[5]) * (std::log(x) / std::log(10) - v[5])));
1142     };
1143 
1144     for (auto const xv : xValues) {
1145       x[0] = xv;
1146       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1147     }
1148   }
1149 
1150   {
1151     reco::FormulaEvaluator f(
1152         "[2]*([3]*([4]+[5]*TMath::Log(max([0],min([1],x))))*1./([6]+[7]*100./"
1153         "3.*(TMath::Max(0.,1.03091-0.051154*pow(x,-0.154227))-TMath::Max(0.,1.03091-0.051154*TMath::Power(208.,-0."
1154         "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))))"
1155         "/max(30.,min(6500.,x)))-(1+0.04432-1.304*pow(208.,-0.4624)+(0+1.724*TMath::Log(208.))/208.))))");
1156 
1157     std::vector<double> v = {55, 2510, 0.997756, 1.000155, 0.979016, 0.001834, 0.982, -0.048, 1.250};
1158 
1159     std::vector<double> x = {100};
1160 
1161     auto func = [&v](double x) {
1162       return v[2] *
1163              (v[3] * (v[4] + v[5] * TMath::Log(std::max(v[0], std::min(v[1], x)))) * 1. /
1164               (v[6] +
1165                v[7] * 100. / 3. *
1166                    (TMath::Max(0., 1.03091 - 0.051154 * std::pow(x, -0.154227)) -
1167                     TMath::Max(0., 1.03091 - 0.051154 * TMath::Power(208., -0.154227))) +
1168                v[8] *
1169                    ((1 + 0.04432 - 1.304 * std::pow(std::max(30., std::min(6500., x)), -0.4624) +
1170                      (0 + 1.724 * TMath::Log(std::max(30., std::min(6500., x)))) / std::max(30., std::min(6500., x))) -
1171                     (1 + 0.04432 - 1.304 * std::pow(208., -0.4624) + (0 + 1.724 * TMath::Log(208.)) / 208.))));
1172     };
1173 
1174     CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1175   }
1176 
1177   {
1178     std::vector<std::pair<std::string, double>> formulas = {
1179         {"(1+0.04432+(1.724+100.))-1", 101.76832},
1180         {"(1+(1.724+100.)+0.04432)-1", 101.76832},
1181         {"((1.724+100.)+1+0.04432)-1", 101.76832},
1182         {"(1+0.04432+1.724/100.)-1", .06156},
1183         {"(1+1.724/100.+0.04432)-1", .06156},
1184         {"(1.724/100.+1+0.04432)-1", .06156},
1185         {"(1+0.04432+(1.724/100.))-1", (1 + 0.04432 + (1.724 / 100.)) - 1},
1186         {"(1+(1.724/100.)+0.04432)-1", (1 + 0.04432 + (1.724 / 100.)) - 1},
1187         {"((1.724/100.)+1+0.04432)-1", (1 + 0.04432 + (1.724 / 100.)) - 1},
1188         {"0.997756*(1.000155*(0.979016+0.001834*TMath::Log(max(55.,min(2510.,100.))))*1./(0.982+-0.048*100./"
1189          "3.*(TMath::Max(0.,1.03091-0.051154*pow(100.,-0.154227))-TMath::Max(0.,1.03091-0.051154*TMath::Power(208.,-0."
1190          "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."
1191          ",100.))))/max(30.,min(6500.,100.)))-(1+0.04432-1.304*pow(208.,-0.4624)+(0+1.724*TMath::Log(208.))/208.))))",
1192          0.997756 *
1193              (1.000155 * (0.979016 + 0.001834 * TMath::Log(std::max(55., std::min(2510., 100.)))) * 1. /
1194               (0.982 +
1195                -0.048 * 100. / 3. *
1196                    (TMath::Max(0., 1.03091 - 0.051154 * std::pow(100., -0.154227)) -
1197                     TMath::Max(0., 1.03091 - 0.051154 * TMath::Power(208., -0.154227))) +
1198                1.250 * ((1 + 0.04432 - 1.304 * std::pow(std::max(30., std::min(6500., 100.)), -0.4624) +
1199                          (0 + 1.724 * TMath::Log(std::max(30., std::min(6500., 100.)))) /
1200                              std::max(30., std::min(6500., 100.))) -
1201                         (1 + 0.04432 - 1.304 * std::pow(208., -0.4624) + (0 + 1.724 * TMath::Log(208.)) / 208.))))}};
1202 
1203     std::vector<double> x = {};
1204     std::vector<double> v = {};
1205     for (auto const& form_val : formulas) {
1206       reco::FormulaEvaluator f(form_val.first);
1207 
1208       CPPUNIT_ASSERT(compare(f.evaluate(x, v), form_val.second));
1209     }
1210   }
1211 
1212   //tests for JER
1213   {
1214     reco::FormulaEvaluator f("[0]+[1]*exp(-x/[2])");
1215 
1216     std::vector<double> v = {0.006467, 0.02519, 77.08};
1217     std::vector<double> x = {100.};
1218 
1219     auto func = [&v](double x) { return v[0] + v[1] * std::exp(-x / v[2]); };
1220 
1221     CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1222   }
1223 
1224   {
1225     reco::FormulaEvaluator f("max(0.0001,1-y/x*([1]*(z-[0])*(1+[2]*log(x/30.))))");
1226 
1227     std::vector<double> v = {1.4, 0.453645, -0.015665};
1228     std::vector<double> x = {157.2, 0.5, 23.2};
1229 
1230     auto func = [&v](double x, double y, double z) {
1231       return std::max(0.0001, 1 - y / x * (v[1] * (z - v[0]) * (1 + v[2] * std::log(x / 30.))));
1232     };
1233 
1234     CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0], x[1], x[2])));
1235   }
1236   {
1237     reco::FormulaEvaluator f("max(0.0001,1-y*[1]*(z-[0])*(1+[2]*log(x/30.))/x)");
1238 
1239     std::vector<double> v = {1.4, 0.453645, -0.015665};
1240     std::vector<double> x = {157.2, 0.5, 23.2};
1241 
1242     auto func = [&v](double x, double y, double z) {
1243       return std::max(0.0001, 1 - y / x * (v[1] * (z - v[0]) * (1 + v[2] * std::log(x / 30.))));
1244     };
1245 
1246     CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0], x[1], x[2])));
1247   }
1248   {
1249     reco::FormulaEvaluator f("max(0.0001,1-y*([1]*(z-[0])*(1+[2]*log(x/30.)))/x)");
1250 
1251     std::vector<double> v = {1.4, 0.453645, -0.015665};
1252     std::vector<double> x = {157.2, 0.5, 23.2};
1253 
1254     auto func = [&v](double x, double y, double z) {
1255       return std::max(0.0001, 1 - y * (v[1] * (z - v[0]) * (1 + v[2] * std::log(x / 30.))) / x);
1256     };
1257 
1258     CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0], x[1], x[2])));
1259   }
1260 
1261   {
1262     reco::FormulaEvaluator f("sqrt([0]*abs([0])/(x*x)+[1]*[1]*pow(x,[3])+[2]*[2])");
1263 
1264     std::vector<double> v = {1.326, 0.4209, 0.02223, -0.6704};
1265     std::vector<double> x = {100.};
1266 
1267     auto func = [&v](double x) {
1268       return std::sqrt(v[0] * std::abs(v[0]) / (x * x) + v[1] * v[1] * std::pow(x, v[3]) + v[2] * v[2]);
1269     };
1270 
1271     CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1272   }
1273 
1274   {
1275     reco::FormulaEvaluator f("sqrt([0]*[0]/(x*x)+[1]*[1]/x+[2]*[2])");
1276 
1277     std::vector<double> v = {2.3, 0.20, 0.009};
1278     std::vector<double> x = {100.};
1279 
1280     auto func = [&v](double x) { return std::sqrt(v[0] * v[0] / (x * x) + v[1] * v[1] / x + v[2] * v[2]); };
1281 
1282     CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1283   }
1284 
1285   {
1286     //This was previously evaluated wrong
1287     std::array<double, 3> xs = {{224., 225., 226.}};
1288     std::vector<double> x = {0.};
1289     std::vector<double> v = {-3., -3., -3., -3., -3., -3.};
1290 
1291     reco::FormulaEvaluator f(
1292         "([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)");
1293 
1294     auto func = [&v](double x) {
1295       return (v[0] + v[1] * x + v[2] * (x * x)) * (x < 225) +
1296              (v[0] + v[1] * 225 + v[2] * (225 * 225) + v[3] * (x - 225) + v[4] * ((x - 225) * (x - 225)) +
1297               v[5] * ((x - 225) * (x - 225) * (x - 225))) *
1298                  (x > 225);
1299     };
1300 
1301     for (auto x_i : xs) {
1302       x[0] = x_i;
1303       CPPUNIT_ASSERT(compare(f.evaluate(x, v), func(x[0])));
1304     }
1305   }
1306 
1307   {
1308     auto t = []() { reco::FormulaEvaluator f("doesNotExist(2)"); };
1309 
1310     CPPUNIT_ASSERT_THROW(t(), cms::Exception);
1311   }
1312 
1313   {
1314     auto t = []() { reco::FormulaEvaluator f("doesNotExist(2) + abs(-1)"); };
1315 
1316     CPPUNIT_ASSERT_THROW(t(), cms::Exception);
1317   }
1318 
1319   {
1320     auto t = []() { reco::FormulaEvaluator f("abs(-1) + doesNotExist(2)"); };
1321 
1322     CPPUNIT_ASSERT_THROW(t(), cms::Exception);
1323   }
1324 
1325   {
1326     auto t = []() { reco::FormulaEvaluator f("abs(-1) + ( 5 * doesNotExist(2))"); };
1327 
1328     CPPUNIT_ASSERT_THROW(t(), cms::Exception);
1329   }
1330 
1331   {
1332     auto t = []() { reco::FormulaEvaluator f("( 5 * doesNotExist(2)) + abs(-1)"); };
1333 
1334     CPPUNIT_ASSERT_THROW(t(), cms::Exception);
1335   }
1336 
1337   {
1338     auto t = []() { reco::FormulaEvaluator f("TMath::Exp(2)"); };
1339 
1340     CPPUNIT_ASSERT_THROW(t(), cms::Exception);
1341   }
1342 
1343   {
1344     //this was previously causing a seg fault
1345     auto t = []() { reco::FormulaEvaluator f("1 + 2 * 3 + 5 * doesNotExist(2) "); };
1346 
1347     CPPUNIT_ASSERT_THROW(t(), cms::Exception);
1348   }
1349 
1350   {
1351     //Make sure spaces are shown in exception message
1352     try {
1353       reco::FormulaEvaluator f("1 + 2#");
1354     } catch (cms::Exception const& e) {
1355       auto r =
1356           "An exception of category 'FormulaEvaluatorParseError' occurred.\n"
1357           "Exception Message:\n"
1358           "While parsing '1 + 2#' could not parse beyond '1 + 2'\n";
1359       CPPUNIT_ASSERT(std::string(r) == e.what());
1360     }
1361   }
1362 }