[Minuit2] Migrate InitialGradientCalculator to free function

Author: guitargeek

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
 

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; }

math/minuit2/src/HessianGradientCalculator.cxx

@@ -29,8 +29,7 @@ namespace Minuit2 {
 FunctionGradient HessianGradientCalculator::operator()(const MinimumParameters &par) const
 {
    // use initial gradient as starting point
-   InitialGradientCalculator gc(fFcn, fTransformation);
-   FunctionGradient gra = gc(par);
+   FunctionGradient gra = calculateInitialGradient(par, fTransformation, fFcn.ErrorDef());
 
    return (*this)(par, gra);
 }

math/minuit2/src/InitialGradientCalculator.cxx

@@ -22,13 +22,13 @@ namespace ROOT {
 
 namespace Minuit2 {
 
-FunctionGradient InitialGradientCalculator::operator()(const MinimumParameters &par) const
+/// Initial rough estimate of the gradient using the parameter step size.
+FunctionGradient
+calculateInitialGradient(const MinimumParameters &par, const MnUserTransformation &trafo, double errorDef)
 {
-   // initial rough  estimate of the gradient using the parameter step size
-
    assert(par.IsValid());
 
-   unsigned int n = Trafo().VariableParameters();
+   unsigned int n = trafo.VariableParameters();
    assert(n == par.Vec().size());
 
    MnPrint print("InitialGradientCalculator");
@@ -38,58 +38,46 @@ FunctionGradient InitialGradientCalculator::operator()(const MinimumParameters &
    MnAlgebraicVector gr(n), gr2(n), gst(n);
 
    for (unsigned int i = 0; i < n; i++) {
-      unsigned int exOfIn = Trafo().ExtOfInt(i);
+      unsigned int exOfIn = trafo.ExtOfInt(i);
 
       double var = par.Vec()(i);
-      double werr = Trafo().Parameter(exOfIn).Error();
-      double save1 = Trafo().Int2ext(i, var);
+      double werr = trafo.Parameter(exOfIn).Error();
+      double save1 = trafo.Int2ext(i, var);
       double save2 = save1 + werr;
-      if (Trafo().Parameter(exOfIn).HasLimits()) {
-         if (Trafo().Parameter(exOfIn).HasUpperLimit() && save2 > 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;
       }
       gr(i) = grd;
       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

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);
 }

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++)