diff --git a/math/minuit2/CMakeLists.txt b/math/minuit2/CMakeLists.txt index c204e61b78ae7..3807757135abb 100644 --- a/math/minuit2/CMakeLists.txt +++ b/math/minuit2/CMakeLists.txt @@ -74,16 +74,15 @@ if(CMAKE_PROJECT_NAME STREQUAL ROOT) Minuit2/MnMinos.h Minuit2/MnParabola.h Minuit2/MnParabolaFactory.h - Minuit2/MnParabolaPoint.h Minuit2/MnParameterScan.h Minuit2/MnPlot.h + Minuit2/MnPoint.h Minuit2/MnPosDef.h Minuit2/MnPrint.h Minuit2/MnScan.h Minuit2/MnSeedGenerator.h Minuit2/MnSimplex.h Minuit2/MnStrategy.h - Minuit2/MnTiny.h Minuit2/MnTraceObject.h Minuit2/MnUserCovariance.h Minuit2/MnUserFcn.h @@ -148,7 +147,6 @@ if(CMAKE_PROJECT_NAME STREQUAL ROOT) src/MnScan.cxx src/MnSeedGenerator.cxx src/MnStrategy.cxx - src/MnTiny.cxx src/MnTraceObject.cxx src/MnUserFcn.cxx src/MnUserParameterState.cxx diff --git a/math/minuit2/inc/Minuit2/InitialGradientCalculator.h b/math/minuit2/inc/Minuit2/InitialGradientCalculator.h index 93dfcbefdac4d..bf6a080f425df 100644 --- a/math/minuit2/inc/Minuit2/InitialGradientCalculator.h +++ b/math/minuit2/inc/Minuit2/InitialGradientCalculator.h @@ -16,30 +16,9 @@ namespace ROOT { namespace Minuit2 { -class MnFcn; class MnUserTransformation; -class MnMachinePrecision; -/** - Class to calculate an initial estimate of the gradient - */ -class InitialGradientCalculator : public GradientCalculator { - -public: - InitialGradientCalculator(const MnFcn &fcn, const MnUserTransformation &par) : fFcn(fcn), fTransformation(par) {} - - FunctionGradient operator()(const MinimumParameters &) const override; - - FunctionGradient operator()(const MinimumParameters &, const FunctionGradient &) const override; - - const MnFcn &Fcn() const { return fFcn; } - const MnUserTransformation &Trafo() const { return fTransformation; } - const MnMachinePrecision &Precision() const; - -private: - const MnFcn &fFcn; - const MnUserTransformation &fTransformation; -}; +FunctionGradient calculateInitialGradient(const MinimumParameters &, const MnUserTransformation &, double errorDef); } // namespace Minuit2 diff --git a/math/minuit2/inc/Minuit2/MnLineSearch.h b/math/minuit2/inc/Minuit2/MnLineSearch.h index c2a36d5ecffc4..630874b4804b7 100644 --- a/math/minuit2/inc/Minuit2/MnLineSearch.h +++ b/math/minuit2/inc/Minuit2/MnLineSearch.h @@ -11,6 +11,7 @@ #define ROOT_Minuit2_MnLineSearch #include "Minuit2/MnMatrix.h" +#include "Minuit2/MnPoint.h" namespace ROOT { @@ -19,7 +20,6 @@ namespace Minuit2 { class MnFcn; class MinimumParameters; class MnMachinePrecision; -class MnParabolaPoint; /** @@ -40,15 +40,15 @@ and Lorenzo Moneta class MnLineSearch { public: - MnParabolaPoint operator()(const MnFcn &, const MinimumParameters &, const MnAlgebraicVector &, double, - const MnMachinePrecision &) const; + MnPoint operator()(const MnFcn &, const MinimumParameters &, const MnAlgebraicVector &, double, + const MnMachinePrecision &) const; #ifdef USE_OTHER_LS - MnParabolaPoint CubicSearch(const MnFcn &, const MinimumParameters &, const MnAlgebraicVector &, double, double, - const MnMachinePrecision &) const; + MnPoint CubicSearch(const MnFcn &, const MinimumParameters &, const MnAlgebraicVector &, double, double, + const MnMachinePrecision &) const; - MnParabolaPoint BrentSearch(const MnFcn &, const MinimumParameters &, const MnAlgebraicVector &, double, double, - const MnMachinePrecision &) const; + MnPoint BrentSearch(const MnFcn &, const MinimumParameters &, const MnAlgebraicVector &, double, double, + const MnMachinePrecision &) const; #endif }; diff --git a/math/minuit2/inc/Minuit2/MnMachinePrecision.h b/math/minuit2/inc/Minuit2/MnMachinePrecision.h index 59d9972f7e254..2ed1c0dece08f 100644 --- a/math/minuit2/inc/Minuit2/MnMachinePrecision.h +++ b/math/minuit2/inc/Minuit2/MnMachinePrecision.h @@ -52,8 +52,12 @@ class MnMachinePrecision { void ComputePrecision(); private: + double One() const; + double Tiny(double epsp1) const; + double fEpsMac; double fEpsMa2; + double fOne = 1.; }; } // namespace Minuit2 diff --git a/math/minuit2/inc/Minuit2/MnParabola.h b/math/minuit2/inc/Minuit2/MnParabola.h index 9c1c11665455f..5ae0ba9791bd2 100644 --- a/math/minuit2/inc/Minuit2/MnParabola.h +++ b/math/minuit2/inc/Minuit2/MnParabola.h @@ -30,118 +30,29 @@ and Lorenzo Moneta class MnParabola { public: - /** - - Constructor that initializes the parabola with its three parameters. - - @param a the coefficient of the quadratic term - @param b the coefficient of the linear term - @param c the constant - - */ - + /// Constructor that initializes the parabola with its three parameters. + /// + /// @param a the coefficient of the quadratic term. + /// @param b the coefficient of the linear term. + /// @param c the constant. MnParabola(double a, double b, double c) : fA(a), fB(b), fC(c) {} - /** - - Evaluates the parabola a the point x. - - @param x the coordinate where the parabola needs to be evaluated. - - @return the y coordinate of the parabola corresponding to x. - - */ - + /// Evaluates the parabola a the point x. double Y(double x) const { return (fA * x * x + fB * x + fC); } - /** - - Calculates the bigger of the two x values corresponding to the - given y Value. - -

- - ???????!!!!!!!!! And when there is none?? it looks like it will - crash?? what is sqrt (-1.0) ? - - @param y the y Value for which the x Value is to be calculated. - - @return the bigger one of the two corresponding values. - - */ - - // ok, at first glance it does not look like the formula for the quadratic - // equation, but it is! ;-) - double X_pos(double y) const { return (std::sqrt(y / fA + Min() * Min() - fC / fA) + Min()); } - // maybe it is worth to check the performance improvement with the below formula?? - // double X_pos(double y) const {return (std::sqrt(y/fA + fB*fB/(4.*fA*fA) - fC/fA) - fB/(2.*fA));} - - /** - - Calculates the smaller of the two x values corresponding to the - given y Value. - -

- - ???????!!!!!!!!! And when there is none?? it looks like it will - crash?? what is sqrt (-1.0) ? - - @param y the y Value for which the x Value is to be calculated. - - @return the smaller one of the two corresponding values. - - */ - - double X_neg(double y) const { return (-std::sqrt(y / fA + Min() * Min() - fC / fA) + Min()); } - - /** - - Calculates the x coordinate of the Minimum of the parabola. - - @return x coordinate of the Minimum. - - */ - + /// Calculate the x coordinate of the Minimum of the parabola. double Min() const { return -fB / (2. * fA); } - /** - - Calculates the y coordinate of the Minimum of the parabola. - - @return y coordinate of the Minimum. - - */ - + /// Calculate the y coordinate of the Minimum of the parabola. double YMin() const { return (-fB * fB / (4. * fA) + fC); } - /** - - Accessor to the coefficient of the quadratic term. - - @return the coefficient of the quadratic term. - - */ - + /// Get the coefficient of the quadratic term. double A() const { return fA; } - /** - - Accessor to the coefficient of the linear term. - - @return the coefficient of the linear term. - - */ - + /// Get the coefficient of the linear term. double B() const { return fB; } - /** - - Accessor to the coefficient of the constant term. - - @return the coefficient of the constant term. - - */ - + /// Get the coefficient of the constant term. double C() const { return fC; } private: diff --git a/math/minuit2/inc/Minuit2/MnParabolaFactory.h b/math/minuit2/inc/Minuit2/MnParabolaFactory.h index 57f7da674f112..89547e42c7160 100644 --- a/math/minuit2/inc/Minuit2/MnParabolaFactory.h +++ b/math/minuit2/inc/Minuit2/MnParabolaFactory.h @@ -10,18 +10,18 @@ #ifndef ROOT_Minuit2_MnParabolaFactory #define ROOT_Minuit2_MnParabolaFactory +#include "Minuit2/MnParabola.h" +#include "Minuit2/MnPoint.h" + namespace ROOT { namespace Minuit2 { -class MnParabola; -class MnParabolaPoint; - class MnParabolaFactory { public: - MnParabola operator()(const MnParabolaPoint &, const MnParabolaPoint &, const MnParabolaPoint &) const; + MnParabola operator()(const MnPoint &, const MnPoint &, const MnPoint &) const; - MnParabola operator()(const MnParabolaPoint &, double, const MnParabolaPoint &) const; + MnParabola operator()(const MnPoint &, double, const MnPoint &) const; }; } // namespace Minuit2 diff --git a/math/minuit2/inc/Minuit2/MnParabolaPoint.h b/math/minuit2/inc/Minuit2/MnParabolaPoint.h deleted file mode 100644 index 21b4ab7e58997..0000000000000 --- a/math/minuit2/inc/Minuit2/MnParabolaPoint.h +++ /dev/null @@ -1,79 +0,0 @@ -// @(#)root/minuit2:$Id$ -// Authors: M. Winkler, F. James, L. Moneta, A. Zsenei 2003-2005 - -/********************************************************************** - * * - * Copyright (c) 2005 LCG ROOT Math team, CERN/PH-SFT * - * * - **********************************************************************/ - -#ifndef ROOT_Minuit2_MnParabolaPoint -#define ROOT_Minuit2_MnParabolaPoint - -namespace ROOT { - -namespace Minuit2 { - -/** - -A point of a parabola. - -

- -????!!!! in reality it is just a general point in two dimensional space, -there is nothing that would indicate, that it belongs to a parabola. -This class defines simply an (x,y) pair!!!! - -@author Fred James and Matthias Winkler; comments added by Andras Zsenei -and Lorenzo Moneta - -@ingroup Minuit - -\todo Should it be called MnParabolaPoint or just Point? - - */ - -class MnParabolaPoint { - -public: - /** - - Initializes the point with its coordinates. - - @param x the x (first) coordinate of the point. - @param y the y (second) coordinate of the point. - - */ - - MnParabolaPoint(double x, double y) : fX(x), fY(y) {} - - /** - - Accessor to the x (first) coordinate. - - @return the x (first) coordinate of the point. - - */ - - double X() const { return fX; } - - /** - - Accessor to the y (second) coordinate. - - @return the y (second) coordinate of the point. - - */ - - double Y() const { return fY; } - -private: - double fX; - double fY; -}; - -} // namespace Minuit2 - -} // namespace ROOT - -#endif // ROOT_Minuit2_MnParabolaPoint diff --git a/math/minuit2/inc/Minuit2/MnPoint.h b/math/minuit2/inc/Minuit2/MnPoint.h new file mode 100644 index 0000000000000..b6272f325266c --- /dev/null +++ b/math/minuit2/inc/Minuit2/MnPoint.h @@ -0,0 +1,52 @@ +// @(#)root/minuit2:$Id$ +// Authors: M. Winkler, F. James, L. Moneta, A. Zsenei 2003-2005 + +/********************************************************************** + * * + * Copyright (c) 2005 LCG ROOT Math team, CERN/PH-SFT * + * * + **********************************************************************/ + +#ifndef ROOT_Minuit2_MnPoint +#define ROOT_Minuit2_MnPoint + +namespace ROOT { + +namespace Minuit2 { + +/** + +A point in x-y. + +@author Fred James and Matthias Winkler; comments added by Andras Zsenei +and Lorenzo Moneta + +@ingroup Minuit + + */ + +class MnPoint { + +public: + /// Initializes the point with its coordinates. + /// + /// @param x the x (first) coordinate of the point. + /// @param y the y (second) coordinate of the point. + MnPoint(double x, double y) : fX(x), fY(y) {} + + /// Get the x (first) coordinate. + double X() const { return fX; } + + /// Get the y (second) coordinate. + double Y() const { return fY; } + +private: + double fX; + double fY; +}; + +} // namespace Minuit2 + +} // namespace ROOT + +#endif // ROOT_Minuit2_MnPoint diff --git a/math/minuit2/inc/Minuit2/MnTiny.h b/math/minuit2/inc/Minuit2/MnTiny.h deleted file mode 100644 index 0410b837903e3..0000000000000 --- a/math/minuit2/inc/Minuit2/MnTiny.h +++ /dev/null @@ -1,32 +0,0 @@ -// @(#)root/minuit2:$Id$ -// Authors: M. Winkler, F. James, L. Moneta, A. Zsenei 2003-2005 - -/********************************************************************** - * * - * Copyright (c) 2005 LCG ROOT Math team, CERN/PH-SFT * - * * - **********************************************************************/ - -#ifndef ROOT_Minuit2_MnTiny -#define ROOT_Minuit2_MnTiny - -namespace ROOT { - -namespace Minuit2 { - -class MnTiny { - -public: - double One() const; - - double operator()(double epsp1) const; - -private: - double fOne = 1.; -}; - -} // namespace Minuit2 - -} // namespace ROOT - -#endif // ROOT_Minuit2_MnTiny diff --git a/math/minuit2/inc/Minuit2/Numerical2PGradientCalculator.h b/math/minuit2/inc/Minuit2/Numerical2PGradientCalculator.h index 2a9daf062f8aa..64ecc9a763825 100644 --- a/math/minuit2/inc/Minuit2/Numerical2PGradientCalculator.h +++ b/math/minuit2/inc/Minuit2/Numerical2PGradientCalculator.h @@ -43,7 +43,6 @@ class Numerical2PGradientCalculator : public GradientCalculator { FunctionGradient operator()(const MinimumParameters &, const FunctionGradient &) const override; - const MnFcn &Fcn() const { return fFcn; } const MnUserTransformation &Trafo() const { return fTransformation; } const MnMachinePrecision &Precision() const; const MnStrategy &Strategy() const { return fStrategy; } diff --git a/math/minuit2/src/CMakeLists.txt b/math/minuit2/src/CMakeLists.txt index 630750baf233b..2344c6e07b808 100644 --- a/math/minuit2/src/CMakeLists.txt +++ b/math/minuit2/src/CMakeLists.txt @@ -64,16 +64,15 @@ set(MINUIT2_HEADERS MnMinos.h MnParabola.h MnParabolaFactory.h - MnParabolaPoint.h MnParameterScan.h MnPlot.h + MnPoint.h MnPosDef.h MnPrint.h MnScan.h MnSeedGenerator.h MnSimplex.h MnStrategy.h - MnTiny.h MnTraceObject.h MnUserCovariance.h MnUserFcn.h @@ -137,7 +136,6 @@ set(MINUIT2_SOURCES MnScan.cxx MnSeedGenerator.cxx MnStrategy.cxx - MnTiny.cxx MnTraceObject.cxx MnUserFcn.cxx MnUserParameterState.cxx diff --git a/math/minuit2/src/FumiliBuilder.cxx b/math/minuit2/src/FumiliBuilder.cxx index f9bf7f3b93caf..9f2d867a29904 100644 --- a/math/minuit2/src/FumiliBuilder.cxx +++ b/math/minuit2/src/FumiliBuilder.cxx @@ -10,7 +10,6 @@ #include "Minuit2/FumiliBuilder.h" #include "Minuit2/FumiliStandardMaximumLikelihoodFCN.h" #include "Minuit2/GradientCalculator.h" -//#include "Minuit2/Numerical2PGradientCalculator.h" #include "Minuit2/MinimumState.h" #include "Minuit2/MinimumError.h" #include "Minuit2/FunctionGradient.h" @@ -20,7 +19,6 @@ #include "Minuit2/MnFcn.h" #include "Minuit2/MnMachinePrecision.h" #include "Minuit2/MnPosDef.h" -#include "Minuit2/MnParabolaPoint.h" #include "Minuit2/MnStrategy.h" #include "Minuit2/MnHesse.h" #include "Minuit2/MnPrint.h" @@ -265,7 +263,6 @@ FunctionMinimum FumiliBuilder::Minimum(const MnFcn &fcn, const GradientCalculato } } - // take a full step //evaluate function only if doing a line search @@ -277,7 +274,7 @@ FunctionMinimum FumiliBuilder::Minimum(const MnFcn &fcn, const GradientCalculato if (doLineSearch && p.Fval() >= s0.Fval()) { print.Debug("Do a line search", fcn.NumOfCalls()); MnLineSearch lsearch; - MnParabolaPoint pp = lsearch(fcn, s0.Parameters(), step, gdel, prec); + auto pp = lsearch(fcn, s0.Parameters(), step, gdel, prec); if (std::fabs(pp.Y() - s0.Fval()) < prec.Eps()) { // std::cout<<"FumiliBuilder: no improvement"< Trafo().Parameter(exOfIn).UpperLimit()) - save2 = Trafo().Parameter(exOfIn).UpperLimit(); + if (trafo.Parameter(exOfIn).HasLimits()) { + if (trafo.Parameter(exOfIn).HasUpperLimit() && save2 > trafo.Parameter(exOfIn).UpperLimit()) + save2 = trafo.Parameter(exOfIn).UpperLimit(); } - double var2 = Trafo().Ext2int(exOfIn, save2); + double var2 = trafo.Ext2int(exOfIn, save2); double vplu = var2 - var; save2 = save1 - werr; - if (Trafo().Parameter(exOfIn).HasLimits()) { - if (Trafo().Parameter(exOfIn).HasLowerLimit() && save2 < Trafo().Parameter(exOfIn).LowerLimit()) - save2 = Trafo().Parameter(exOfIn).LowerLimit(); + if (trafo.Parameter(exOfIn).HasLimits()) { + if (trafo.Parameter(exOfIn).HasLowerLimit() && save2 < trafo.Parameter(exOfIn).LowerLimit()) + save2 = trafo.Parameter(exOfIn).LowerLimit(); } - var2 = Trafo().Ext2int(exOfIn, save2); + var2 = trafo.Ext2int(exOfIn, save2); double vmin = var2 - var; - double gsmin = 8. * Precision().Eps2() * (std::fabs(var) + Precision().Eps2()); + double gsmin = 8. * trafo.Precision().Eps2() * (std::fabs(var) + trafo.Precision().Eps2()); // protect against very small step sizes which can cause dirin to zero and then nan values in grd double dirin = std::max(0.5 * (std::fabs(vplu) + std::fabs(vmin)), gsmin); - double g2 = 2.0 * fFcn.ErrorDef() / (dirin * dirin); + double g2 = 2.0 * errorDef / (dirin * dirin); double gstep = std::max(gsmin, 0.1 * dirin); double grd = g2 * dirin; - if (Trafo().Parameter(exOfIn).HasLimits()) { + if (trafo.Parameter(exOfIn).HasLimits()) { if (gstep > 0.5) gstep = 0.5; } @@ -71,25 +71,13 @@ FunctionGradient InitialGradientCalculator::operator()(const MinimumParameters & gr2(i) = g2; gst(i) = gstep; - print.Trace("Computed initial gradient for parameter", Trafo().Name(exOfIn), "value", var, "[", vmin, ",", vplu, + print.Trace("Computed initial gradient for parameter", trafo.Name(exOfIn), "value", var, "[", vmin, ",", vplu, "]", "dirin", dirin, "grd", grd, "g2", g2); } return FunctionGradient(gr, gr2, gst); } -FunctionGradient InitialGradientCalculator::operator()(const MinimumParameters &par, const FunctionGradient &) const -{ - // Base class interface - return (*this)(par); -} - -const MnMachinePrecision &InitialGradientCalculator::Precision() const -{ - // return precision (is set in transformation class) - return fTransformation.Precision(); -} - } // namespace Minuit2 } // namespace ROOT diff --git a/math/minuit2/src/MnFunctionCross.cxx b/math/minuit2/src/MnFunctionCross.cxx index d3ca08bd7b637..3a66fa46649f2 100644 --- a/math/minuit2/src/MnFunctionCross.cxx +++ b/math/minuit2/src/MnFunctionCross.cxx @@ -11,8 +11,6 @@ #include "Minuit2/FunctionMinimum.h" #include "Minuit2/MnMigrad.h" #include "Minuit2/FCNBase.h" -#include "Minuit2/MnParabola.h" -#include "Minuit2/MnParabolaPoint.h" #include "Minuit2/MnParabolaFactory.h" #include "Minuit2/MnCross.h" #include "Minuit2/MnMachinePrecision.h" @@ -388,8 +386,7 @@ MnCross MnFunctionCross::operator()(std::span par, std::span do { // do parabola fit - MnParabola parbol = MnParabolaFactory()(MnParabolaPoint(alsb[0], flsb[0]), MnParabolaPoint(alsb[1], flsb[1]), - MnParabolaPoint(alsb[2], flsb[2])); + MnParabola parbol = MnParabolaFactory()({alsb[0], flsb[0]}, {alsb[1], flsb[1]}, {alsb[2], flsb[2]}); // aopt = parbol.X_pos(aim); // std::cout<<"alsb1,2,3= "<= maxiter) { // exhausted max number of iterations - return MnParabolaPoint(xvmin, fvmin); + return {xvmin, fvmin}; } print.Trace("after initial 2-point iter:", '\n', " x0, x1, x2:", p0.X(), p1.X(), slam, '\n', " f0, f1, f2:", p0.Y(), p1.Y(), f2); - MnParabolaPoint p2(slam, f2); + MnPoint p2(slam, f2); // do now the quadratic interpolation with 3 points do { @@ -257,7 +256,7 @@ MnParabolaPoint MnLineSearch::operator()(const MnFcn &fcn, const MinimumParamete // std::cout<<"f1, f2, f3= "<= maxiter) { // exhausted max number of iterations print.Debug("Exhausted max number of iterations ",maxiter," return"); - return MnParabolaPoint(xvmin, fvmin); + return {xvmin, fvmin}; } // find worst previous point out of three and replace - MnParabolaPoint p3(slam, f3); + MnPoint p3(slam, f3); if (p0.Y() > p1.Y() && p0.Y() > p2.Y()) p0 = p3; else if (p1.Y() > p0.Y() && p1.Y() > p2.Y()) @@ -307,7 +306,7 @@ MnParabolaPoint MnLineSearch::operator()(const MnFcn &fcn, const MinimumParamete } while (niter < maxiter); print.Debug("f1, f2 =", p0.Y(), p1.Y(), '\n', "x1, x2 =", p0.X(), p1.X(), '\n', "x, f =", xvmin, fvmin); - return MnParabolaPoint(xvmin, fvmin); + return {xvmin, fvmin}; } #ifdef USE_OTHER_LS @@ -316,8 +315,8 @@ MnParabolaPoint MnLineSearch::operator()(const MnFcn &fcn, const MinimumParamete This is used at the beginning when the second derivative is known to be negative */ -MnParabolaPoint MnLineSearch::CubicSearch(const MnFcn &fcn, const MinimumParameters &st, const MnAlgebraicVector &step, - double gdel, double g2del, const MnMachinePrecision &prec) const +MnPoint MnLineSearch::CubicSearch(const MnFcn &fcn, const MinimumParameters &st, const MnAlgebraicVector &step, + double gdel, double g2del, const MnMachinePrecision &prec) const { MnPrint print("MnLineSearch::CubicSearch"); @@ -359,8 +358,8 @@ MnParabolaPoint MnLineSearch::CubicSearch(const MnFcn &fcn, const MinimumParamet double flast = f1; double slam = 1.; - // MnParabolaPoint p0(0., f0); - // MnParabolaPoint p1(slam, flast); + // MnPoint p0(0., f0); + // MnPoint p1(slam, flast); ROOT::Math::SMatrix cubicMatrix; ROOT::Math::SVector cubicCoeff; // cubic coefficients to be found @@ -391,7 +390,7 @@ MnParabolaPoint MnLineSearch::CubicSearch(const MnFcn &fcn, const MinimumParamet // find the minimum need to invert the matrix if (!cubicMatrix.Invert()) { print.Warn("Inversion failed - return"); - return MnParabolaPoint(xvmin, fvmin); + return {xvmin, fvmin}; } cubicCoeff = cubicMatrix * bVec; @@ -496,23 +495,23 @@ MnParabolaPoint MnLineSearch::CubicSearch(const MnFcn &fcn, const MinimumParamet } else if (std::fabs(df2) <= 0) { // gradient is significative different than zero then redo a cubic interpolation // from new point - return MnParabolaPoint(slam, fp); // should redo a cubic interpol. ?? - // niter ++; - // if (niter > maxiter) break; + return {slam, fp}; // should redo a cubic interpol. ?? + // niter ++; + // if (niter > maxiter) break; // MinimumParameters pa = MinimumParameters(st.Vec() + stepNew, fp); // gdel = stepNew(i) - // MnParabolaPoint pp = CubicSearch(fcn, st, stepNew, df, df2 + // MnPoint pp = CubicSearch(fcn, st, stepNew, df, df2 } else - return MnParabolaPoint(slam, fp); + return {slam, fp}; niter++; } while (niter < maxiter); - return MnParabolaPoint(xvmin, fvmin); + return {xvmin, fvmin}; } // class describing Fcn function in one dimension @@ -534,8 +533,8 @@ class ProjectedFcn : public ROOT::Math::IGenFunction { const MnAlgebraicVector &fStep; }; -MnParabolaPoint MnLineSearch::BrentSearch(const MnFcn &fcn, const MinimumParameters &st, const MnAlgebraicVector &step, - double gdel, double g2del, const MnMachinePrecision &prec) const +MnPoint MnLineSearch::BrentSearch(const MnFcn &fcn, const MinimumParameters &st, const MnAlgebraicVector &step, + double gdel, double g2del, const MnMachinePrecision &prec) const { MnPrint print("MnLineSearch::BrentSearch"); @@ -577,8 +576,8 @@ MnParabolaPoint MnLineSearch::BrentSearch(const MnFcn &fcn, const MinimumParamet double x0 = 0; -// MnParabolaPoint p0(0., f0); -// MnParabolaPoint p1(slam, flast); +// MnPoint p0(0., f0); +// MnPoint p1(slam, flast); #ifdef USE_CUBIC ROOT::Math::SMatrix cubicMatrix; @@ -609,7 +608,7 @@ MnParabolaPoint MnLineSearch::BrentSearch(const MnFcn &fcn, const MinimumParamet // find the minimum need to invert the matrix if (!cubicMatrix.Invert()) { print.Warn("Inversion failed - return"); - return MnParabolaPoint(xvmin, fvmin); + return {xvmin, fvmin}; } cubicCoeff = cubicMatrix * bVec; @@ -798,7 +797,7 @@ MnParabolaPoint MnLineSearch::BrentSearch(const MnFcn &fcn, const MinimumParamet if (f0 > fa || f0 > fb) // skip minim 1d try Minuit LS // return (*this) (fcn, st, step, gdel, prec, debug); - return MnParabolaPoint(xvmin, fvmin); + return {xvmin, fvmin}; print.Info("1D minimization using", minType); @@ -810,7 +809,7 @@ MnParabolaPoint MnLineSearch::BrentSearch(const MnFcn &fcn, const MinimumParamet MnPrint::info("result of GSL 1D minimization:", ret, "niter", min.Iterations(), "xmin", min.XMinimum(), "fmin", min.FValMinimum()); - return MnParabolaPoint(min.XMinimum(), min.FValMinimum()); + return {min.XMinimum(), min.FValMinimum()}; } #endif diff --git a/math/minuit2/src/MnMachinePrecision.cxx b/math/minuit2/src/MnMachinePrecision.cxx index 01ca35f8c38e2..ef2b16ecf5838 100644 --- a/math/minuit2/src/MnMachinePrecision.cxx +++ b/math/minuit2/src/MnMachinePrecision.cxx @@ -8,7 +8,6 @@ **********************************************************************/ #include "Minuit2/MnMachinePrecision.h" -#include "Minuit2/MnTiny.h" #include namespace ROOT { @@ -39,8 +38,6 @@ void MnMachinePrecision::ComputePrecision() fEpsMa2 = 2.*std::sqrt(fEpsMac); */ - MnTiny mytiny; - // calculate machine precision double epstry = 0.5; double epsbak = 0.; @@ -49,7 +46,7 @@ void MnMachinePrecision::ComputePrecision() for (int i = 0; i < 100; i++) { epstry *= 0.5; epsp1 = one + epstry; - epsbak = mytiny(epsp1); + epsbak = Tiny(epsp1); if (epsbak < epstry) { fEpsMac = 8. * epstry; fEpsMa2 = 2. * std::sqrt(fEpsMac); @@ -58,6 +55,18 @@ void MnMachinePrecision::ComputePrecision() } } +double MnMachinePrecision::One() const +{ + return fOne; +} + +double MnMachinePrecision::Tiny(double epsp1) const +{ + // evaluate minimal difference between two floating points + double result = epsp1 - One(); + return result; +} + } // namespace Minuit2 } // namespace ROOT diff --git a/math/minuit2/src/MnParabolaFactory.cxx b/math/minuit2/src/MnParabolaFactory.cxx index 78a01e618a437..e23dc8e7ae853 100644 --- a/math/minuit2/src/MnParabolaFactory.cxx +++ b/math/minuit2/src/MnParabolaFactory.cxx @@ -9,14 +9,13 @@ #include "Minuit2/MnParabolaFactory.h" #include "Minuit2/MnParabola.h" -#include "Minuit2/MnParabolaPoint.h" +#include "Minuit2/MnPoint.h" namespace ROOT { namespace Minuit2 { -MnParabola MnParabolaFactory:: -operator()(const MnParabolaPoint &p1, const MnParabolaPoint &p2, const MnParabolaPoint &p3) const +MnParabola MnParabolaFactory::operator()(const MnPoint &p1, const MnPoint &p2, const MnPoint &p3) const { // construct the parabola from 3 points p1,p2,p3 double x1 = p1.X(); @@ -49,7 +48,7 @@ operator()(const MnParabolaPoint &p1, const MnParabolaPoint &p2, const MnParabol return MnParabola(a, b, c); } -MnParabola MnParabolaFactory::operator()(const MnParabolaPoint &p1, double dxdy1, const MnParabolaPoint &p2) const +MnParabola MnParabolaFactory::operator()(const MnPoint &p1, double dxdy1, const MnPoint &p2) const { // construct the parabola from 2 points + derivative at first point dxdy1 double x1 = p1.X(); diff --git a/math/minuit2/src/MnSeedGenerator.cxx b/math/minuit2/src/MnSeedGenerator.cxx index e668c6c3d09a3..f41fe7da50131 100644 --- a/math/minuit2/src/MnSeedGenerator.cxx +++ b/math/minuit2/src/MnSeedGenerator.cxx @@ -20,7 +20,6 @@ #include "Minuit2/MnMachinePrecision.h" #include "Minuit2/MinuitParameter.h" #include "Minuit2/MnLineSearch.h" -#include "Minuit2/MnParabolaPoint.h" #include "Minuit2/MinimumState.h" #include "Minuit2/MnUserParameterState.h" #include "Minuit2/MnStrategy.h" diff --git a/math/minuit2/src/MnTiny.cxx b/math/minuit2/src/MnTiny.cxx deleted file mode 100644 index b7ab518c28543..0000000000000 --- a/math/minuit2/src/MnTiny.cxx +++ /dev/null @@ -1,30 +0,0 @@ -// @(#)root/minuit2:$Id$ -// Authors: M. Winkler, F. James, L. Moneta, A. Zsenei 2003-2005 - -/********************************************************************** - * * - * Copyright (c) 2005 LCG ROOT Math team, CERN/PH-SFT * - * * - **********************************************************************/ - -#include "Minuit2/MnTiny.h" - -namespace ROOT { - -namespace Minuit2 { - -double MnTiny::One() const -{ - return fOne; -} - -double MnTiny::operator()(double epsp1) const -{ - // evaluate minimal difference between two floating points - double result = epsp1 - One(); - return result; -} - -} // namespace Minuit2 - -} // namespace ROOT diff --git a/math/minuit2/src/ModularFunctionMinimizer.cxx b/math/minuit2/src/ModularFunctionMinimizer.cxx index 4fa546c360760..bcdd0b705651e 100644 --- a/math/minuit2/src/ModularFunctionMinimizer.cxx +++ b/math/minuit2/src/ModularFunctionMinimizer.cxx @@ -23,7 +23,6 @@ #include "Minuit2/FCNBase.h" #include "Minuit2/MnStrategy.h" #include "Minuit2/MnLineSearch.h" -#include "Minuit2/MnParabolaPoint.h" #include "Minuit2/MnPrint.h" namespace ROOT { diff --git a/math/minuit2/src/NegativeG2LineSearch.cxx b/math/minuit2/src/NegativeG2LineSearch.cxx index 4fb3f45e93bc7..b06d706479428 100644 --- a/math/minuit2/src/NegativeG2LineSearch.cxx +++ b/math/minuit2/src/NegativeG2LineSearch.cxx @@ -13,7 +13,6 @@ #include "Minuit2/GradientCalculator.h" #include "Minuit2/MnMachinePrecision.h" #include "Minuit2/MnLineSearch.h" -#include "Minuit2/MnParabolaPoint.h" #include "Minuit2/VariableMetricEDMEstimator.h" #include "Minuit2/MnPrint.h" @@ -90,7 +89,7 @@ MinimumState NegativeG2LineSearch::operator()(const MnFcn &fcn, const MinimumSta print.Debug("Iter", iter, "param", i, pa.Vec()(i), "grad2", dgrad.G2()(i), "grad", dgrad.Vec()(i), "grad step", step(i), " gdel ", gdel); - MnParabolaPoint pp = lsearch(fcn, pa, step, gdel, prec); + auto pp = lsearch(fcn, pa, step, gdel, prec); print.Debug("Line search result", pp.X(), "f(0)", pa.Fval(), "f(1)", pp.Y()); diff --git a/math/minuit2/src/Numerical2PGradientCalculator.cxx b/math/minuit2/src/Numerical2PGradientCalculator.cxx index a511ccc09c7b6..e6da20aa7cf42 100644 --- a/math/minuit2/src/Numerical2PGradientCalculator.cxx +++ b/math/minuit2/src/Numerical2PGradientCalculator.cxx @@ -36,8 +36,7 @@ FunctionGradient Numerical2PGradientCalculator::operator()(const MinimumParamete { // calculate gradient using Initial gradient calculator and from MinimumParameters object - InitialGradientCalculator gc(fFcn, fTransformation); - FunctionGradient gra = gc(par); + FunctionGradient gra = calculateInitialGradient(par, fTransformation, fFcn.ErrorDef()); return (*this)(par, gra); } diff --git a/math/minuit2/src/SimplexSeedGenerator.cxx b/math/minuit2/src/SimplexSeedGenerator.cxx index c477b3fc209f4..71f76fec6e177 100644 --- a/math/minuit2/src/SimplexSeedGenerator.cxx +++ b/math/minuit2/src/SimplexSeedGenerator.cxx @@ -33,8 +33,7 @@ operator()(const MnFcn &fcn, const GradientCalculator &, const MnUserParameterSt x(i) = st.IntParameters()[i]; double fcnmin = fcn(x); MinimumParameters pa(x, fcnmin); - InitialGradientCalculator igc(fcn, st.Trafo()); - FunctionGradient dgrad = igc(pa); + FunctionGradient dgrad = calculateInitialGradient(pa, st.Trafo(), fcn.ErrorDef()); MnAlgebraicSymMatrix mat(n); double dcovar = 1.; for (unsigned int i = 0; i < n; i++) diff --git a/math/minuit2/src/VariableMetricBuilder.cxx b/math/minuit2/src/VariableMetricBuilder.cxx index 226378e1205cc..d2b8fcc506c63 100644 --- a/math/minuit2/src/VariableMetricBuilder.cxx +++ b/math/minuit2/src/VariableMetricBuilder.cxx @@ -18,7 +18,6 @@ #include "Minuit2/MnFcn.h" #include "Minuit2/MnMachinePrecision.h" #include "Minuit2/MnPosDef.h" -#include "Minuit2/MnParabolaPoint.h" #include "Minuit2/MnStrategy.h" #include "Minuit2/MnHesse.h" #include "Minuit2/MnPrint.h" @@ -273,7 +272,7 @@ FunctionMinimum VariableMetricBuilder::Minimum(const MnFcn &fcn, const GradientC } } - MnParabolaPoint pp = lsearch(fcn, s0.Parameters(), step, gdel, prec); + auto pp = lsearch(fcn, s0.Parameters(), step, gdel, prec); // <= needed for case 0 <= 0 if (std::fabs(pp.Y() - s0.Fval()) <= std::fabs(s0.Fval()) * prec.Eps()) {