[Minuit2] Only compute global correlation matrix when needed

README/ReleaseNotes/v638/index.md

@@ -62,6 +62,11 @@ The following people have contributed to this new version:
 
 ## Math
 
+### Minimizer interface
+
+* The function `double ROOT::Math::Minimizer::GlobalCC(unsigned int ivar)` was changed to `std::vector<double> ROOT::Math::Minimizer::GlobalCC()`, always returning the full list of global correlations.
+  This change was needed because global correlations are not unconditionally computed and cached in the minimizer anymore. Only computing them when calling `GlobalCC()` avoids unneeded overhead when global correlations are not required.
+
 ### Minuit2
 
 * Behavior change: building ROOT using `minuit2_omp=ON` option no longer enables OpenMP parallelization by default. One has to call now additionally GradientCalculator::SetParallelOMP().

math/mathcore/inc/Math/Minimizer.h

@@ -283,7 +283,7 @@ class Minimizer {
       return ( tmp < 0) ? 0 : CovMatrix(i,j) / std::sqrt( tmp );
    }
 
-   virtual double GlobalCC(unsigned int ivar) const;
+   virtual std::vector<double> GlobalCC() const;
 
    virtual bool GetMinosError(unsigned int ivar , double & errLow, double & errUp, int option = 0);
    virtual bool Hesse();

math/mathcore/src/FitResult.cxx

@@ -209,13 +209,7 @@ void FitResult::FillResult(const std::shared_ptr<ROOT::Math::Minimizer> & min, c
       // minos errors are set separately when calling Fitter::CalculateMinosErrors()
 
       // globalCC
-      fGlobalCC.reserve(npar);
-      for (unsigned int i = 0; i < npar; ++i) {
-         double globcc = min->GlobalCC(i);
-         if (globcc < 0) break; // it is not supported by that minimizer
-         fGlobalCC.push_back(globcc);
-      }
-
+      fGlobalCC = min->GlobalCC();
    }
 
 }
@@ -280,16 +274,8 @@ bool FitResult::Update(const std::shared_ptr<ROOT::Math::Minimizer> & min, const
       }
 
       // update global CC
-      if (fGlobalCC.size() != npar) fGlobalCC.resize(npar);
-      for (unsigned int i = 0; i < npar; ++i) {
-         double globcc = min->GlobalCC(i);
-         if (globcc < 0) {
-            fGlobalCC.clear();
-            break; // it is not supported by that minimizer
-         }
-         fGlobalCC[i] = globcc;
-      }
-
+      fGlobalCC = min->GlobalCC();
+      fGlobalCC.resize(npar); // pad with zeros
    }
    return true;
 }

math/mathcore/src/Minimizer.cxx

@@ -167,10 +167,9 @@ bool Minimizer::GetHessianMatrix(double *hMat) const
  * other parameters which is most strongly correlated with i.
  * Minimizer must overload method if implemented
  */
-double Minimizer::GlobalCC(unsigned int ivar) const
+std::vector<double> Minimizer::GlobalCC() const
 {
-   MATH_UNUSED(ivar);
-   return -1;
+   return {};
 }
 
 /**

math/minuit/inc/TMinuitMinimizer.h

@@ -170,7 +170,7 @@ class TMinuitMinimizer  : public ROOT::Math::Minimizer {
    int CovMatrixStatus() const override;
 
    ///global correlation coefficient for variable i
-   double GlobalCC(unsigned int ) const override;
+   std::vector<double> GlobalCC() const override;
 
    /// minos error for variable i, return false if Minos failed
    bool GetMinosError(unsigned int i, double & errLow, double & errUp, int = 0) override;

math/minuit/src/TMinuitMinimizer.cxx

@@ -730,16 +730,26 @@ int TMinuitMinimizer::CovMatrixStatus() const {
    return fMinuit->fISW[1];
 }
 
-double TMinuitMinimizer::GlobalCC(unsigned int i) const {
-   // global correlation coefficient for parameter i
-   if (!fMinuit) return 0;
-   if (!fMinuit->fGlobcc) return 0;
-   if (int(i) >= fMinuit->fNu) return 0;
-   // get internal number in Minuit
-   int iin = fMinuit->fNiofex[i];
-   // index in TMinuit starts from 1
-   if (iin < 1) return 0;
-   return fMinuit->fGlobcc[iin-1];
+std::vector<double> TMinuitMinimizer::GlobalCC() const
+{
+   std::vector<double> out;
+   // global correlation coefficients
+   if (!fMinuit)
+      return out;
+   if (!fMinuit->fGlobcc)
+      return out;
+   out.resize(fMinuit->fNu);
+   for (int i = 0; i < fMinuit->fNu; ++i) {
+      // get internal number in Minuit
+      int iin = fMinuit->fNiofex[i];
+      // index in TMinuit starts from 1
+      if (iin < 1) {
+         out.clear();
+         break;
+      }
+      out[i] = fMinuit->fGlobcc[iin - 1];
+   }
+   return out;
 }
 
 bool TMinuitMinimizer::GetMinosError(unsigned int i, double & errLow, double & errUp, int ) {

math/minuit2/inc/Minuit2/Minuit2Minimizer.h

@@ -216,7 +216,7 @@ class Minuit2Minimizer : public ROOT::Math::Minimizer {
       is most strongly correlated with i.
       If the variable is fixed or const the return value is zero
     */
-   double GlobalCC(unsigned int i) const override;
+   std::vector<double> GlobalCC() const override;
 
    /**
       get the minos error for parameter i, return false if Minos failed

math/minuit2/inc/Minuit2/MnUserParameterState.h

@@ -37,8 +37,11 @@ class MnUserParameterState {
 public:
    /// default constructor (invalid state)
    MnUserParameterState()
-      : fValid(false), fCovarianceValid(false), fGCCValid(false), fCovStatus(-1), fFVal(0), fEDM(0), fNFcn(0),
-        fParameters(MnUserParameters()), fCovariance(MnUserCovariance()), fIntParameters(std::vector<double>()),
+      : fValid(false),
+        fCovStatus(-1),
+        fParameters(MnUserParameters()),
+        fCovariance(MnUserCovariance()),
+        fIntParameters(std::vector<double>()),
         fIntCovariance(MnUserCovariance())
    {
    }
@@ -61,7 +64,7 @@ class MnUserParameterState {
    // user external representation
    const MnUserParameters &Parameters() const { return fParameters; }
    const MnUserCovariance &Covariance() const { return fCovariance; }
-   const MnGlobalCorrelationCoeff &GlobalCC() const { return fGlobalCC; }
+   MnGlobalCorrelationCoeff GlobalCC() const;
 
    // hessian (inverse of covariance matrix)
    MnUserCovariance Hessian() const;
@@ -78,7 +81,6 @@ class MnUserParameterState {
 
    bool IsValid() const { return fValid; }
    bool HasCovariance() const { return fCovarianceValid; }
-   bool HasGlobalCC() const { return fGCCValid; }
 
    double Fval() const { return fFVal; }
    double Edm() const { return fEDM; }
@@ -151,16 +153,14 @@ class MnUserParameterState {
 
 private:
    bool fValid;
-   bool fCovarianceValid;
-   bool fGCCValid;
+   bool fCovarianceValid = false;
    int fCovStatus; // covariance matrix status
-   double fFVal;
-   double fEDM;
-   unsigned int fNFcn;
+   double fFVal = 0.;
+   double fEDM = 0.;
+   unsigned int fNFcn = 0;
 
    MnUserParameters fParameters;
    MnUserCovariance fCovariance;
-   MnGlobalCorrelationCoeff fGlobalCC;
 
    std::vector<double> fIntParameters;
    MnUserCovariance fIntCovariance;

math/minuit2/src/Minuit2Minimizer.cxx

@@ -849,21 +849,27 @@ double Minuit2Minimizer::Correlation(unsigned int i, unsigned int j) const
    return 0;
 }
 
-double Minuit2Minimizer::GlobalCC(unsigned int i) const
+std::vector<double> Minuit2Minimizer::GlobalCC() const
 {
    // get global correlation coefficient for the parameter i. This is a number between zero and one which gives
    // the correlation between the i-th parameter  and that linear combination of all other parameters which
    // is most strongly correlated with i.
 
-   if (i >= fDim)
-      return 0;
+   std::vector<double> out;
+   MnGlobalCorrelationCoeff globalCC = fState.GlobalCC();
    // no info available when minimization has failed or has some problems
-   if (!fState.HasGlobalCC())
-      return 0;
-   if (fState.Parameter(i).IsFixed() || fState.Parameter(i).IsConst())
-      return 0;
-   unsigned int k = fState.IntOfExt(i);
-   return fState.GlobalCC().GlobalCC()[k];
+   if (!globalCC.IsValid())
+      return out;
+   out.resize(fDim);
+   for (unsigned int i = 0; i < fDim; ++i) {
+      if (fState.Parameter(i).IsFixed() || fState.Parameter(i).IsConst())
+         out[i] = 0;
+      else {
+         unsigned int k = fState.IntOfExt(i);
+         out[i] = globalCC.GlobalCC()[k];
+      }
+   }
+   return out;
 }
 
 bool Minuit2Minimizer::GetMinosError(unsigned int i, double &errLow, double &errUp, int runopt)

math/minuit2/src/MnPrint.cxx

@@ -293,7 +293,7 @@ std::ostream &operator<<(std::ostream &os, const MnUserParameterState &state)
       os << state.Covariance();
    else
       os << "matrix is not present or not valid";
-   if (state.HasGlobalCC())
+   if (state.GlobalCC().IsValid())
       os << "\n  Global correlation coefficients: " << state.GlobalCC();
 
    os.precision(pr);

math/minuit2/src/MnUserParameterState.cxx

@@ -20,15 +20,11 @@ namespace Minuit2 {
 // construct from user parameters (before minimization)
 //
 MnUserParameterState::MnUserParameterState(std::span<const double> par, std::span<const double> err)
-   : fValid(true), fCovarianceValid(false), fGCCValid(false), fCovStatus(-1), fFVal(0.), fEDM(0.), fNFcn(0),
-     fParameters(MnUserParameters(par, err)),
-     fIntParameters(par.begin(), par.end())
+   : fValid(true), fCovStatus(-1), fParameters(MnUserParameters(par, err)), fIntParameters(par.begin(), par.end())
 {
 }
 
-MnUserParameterState::MnUserParameterState(const MnUserParameters &par)
-   : fValid(true), fCovarianceValid(false), fGCCValid(false), fCovStatus(-1), fFVal(0.), fEDM(0.), fNFcn(0),
-     fParameters(par)
+MnUserParameterState::MnUserParameterState(const MnUserParameters &par) : fValid(true), fCovStatus(-1), fParameters(par)
 {
    // construct from user parameters (before minimization)
 
@@ -46,13 +42,7 @@ MnUserParameterState::MnUserParameterState(const MnUserParameters &par)
 // construct from user parameters + errors (before minimization)
 //
 MnUserParameterState::MnUserParameterState(std::span<const double> par, std::span<const double> cov, unsigned int nrow)
-   : fValid(true),
-     fGCCValid(false),
-     fCovStatus(-1),
-     fFVal(0.),
-     fEDM(0.),
-     fNFcn(0),
-     fIntParameters(par.begin(), par.end())
+   : fValid(true), fCovStatus(-1), fIntParameters(par.begin(), par.end())
 {
    // construct from user parameters + errors (before minimization) using std::vector for parameter error and    // an
    // std::vector of size n*(n+1)/2 for the covariance matrix  and n (rank of cov matrix)
@@ -68,13 +58,7 @@ MnUserParameterState::MnUserParameterState(std::span<const double> par, std::spa
 }
 
 MnUserParameterState::MnUserParameterState(std::span<const double> par, const MnUserCovariance &cov)
-   : fValid(true),
-     fGCCValid(false),
-     fCovStatus(-1),
-     fFVal(0.),
-     fEDM(0.),
-     fNFcn(0),
-     fIntParameters(par.begin(), par.end())
+   : fValid(true), fCovStatus(-1), fIntParameters(par.begin(), par.end())
 {
    // construct from user parameters + errors (before minimization) using std::vector (params) and MnUserCovariance
    // class
@@ -90,7 +74,7 @@ MnUserParameterState::MnUserParameterState(std::span<const double> par, const Mn
 }
 
 MnUserParameterState::MnUserParameterState(const MnUserParameters &par, const MnUserCovariance &cov)
-   : fValid(true), fGCCValid(false), fCovStatus(-1), fFVal(0.), fEDM(0.), fNFcn(0), fParameters(par)
+   : fValid(true), fCovStatus(-1), fParameters(par)
 {
    // construct from user parameters + errors (before minimization) using
    // MnUserParameters and MnUserCovariance objects
@@ -110,8 +94,7 @@ MnUserParameterState::MnUserParameterState(const MnUserParameters &par, const Mn
 //
 //
 MnUserParameterState::MnUserParameterState(const MinimumState &st, double up, const MnUserTransformation &trafo)
-   : fValid(st.IsValid()), fCovarianceValid(false), fGCCValid(false), fCovStatus(-1), fFVal(st.Fval()), fEDM(st.Edm()),
-     fNFcn(st.NFcn())
+   : fValid(st.IsValid()), fCovStatus(-1), fFVal(st.Fval()), fEDM(st.Edm()), fNFcn(st.NFcn())
 {
    //
    // construct from internal parameters (after minimization)
@@ -165,8 +148,6 @@ MnUserParameterState::MnUserParameterState(const MinimumState &st, double up, co
       fIntCovariance =
          MnUserCovariance({st.Error().InvHessian().Data(), st.Error().InvHessian().size()}, st.Error().InvHessian().Nrow());
       fCovariance.Scale(2. * up);
-      fGlobalCC = MnGlobalCorrelationCoeff(st.Error().InvHessian());
-      fGCCValid = fGlobalCC.IsValid();
 
       assert(fCovariance.Nrow() == VariableParameters());
 
@@ -237,7 +218,6 @@ void MnUserParameterState::Add(const std::string &name, double val, double err)
    if (fParameters.Add(name, val, err)) {
       fIntParameters.push_back(val);
       fCovarianceValid = false;
-      fGCCValid = false;
       fValid = true;
    } else {
       // redefine an existing parameter
@@ -262,7 +242,6 @@ void MnUserParameterState::Add(const std::string &name, double val, double err,
    if (fParameters.Add(name, val, err, low, up)) {
       fCovarianceValid = false;
       fIntParameters.push_back(Ext2int(Index(name), val));
-      fGCCValid = false;
       fValid = true;
    } else { // Parameter already exist - just set values
       int i = Index(name);
@@ -335,7 +314,6 @@ void MnUserParameterState::Fix(unsigned int e)
       fIntParameters.erase(fIntParameters.begin() + i, fIntParameters.begin() + i + 1);
    }
    fParameters.Fix(e);
-   fGCCValid = false;
 }
 
 void MnUserParameterState::Release(unsigned int e)
@@ -346,7 +324,6 @@ void MnUserParameterState::Release(unsigned int e)
       return;
    fParameters.Release(e);
    fCovarianceValid = false;
-   fGCCValid = false;
    unsigned int i = IntOfExt(e);
    if (Parameter(e).HasLimits())
       fIntParameters.insert(fIntParameters.begin() + i, Ext2int(e, Parameter(e).Value()));
@@ -378,7 +355,6 @@ void MnUserParameterState::SetLimits(unsigned int e, double low, double up)
    // set limits for parameter e (external index)
    fParameters.SetLimits(e, low, up);
    fCovarianceValid = false;
-   fGCCValid = false;
    if (!Parameter(e).IsFixed() && !Parameter(e).IsConst()) {
       unsigned int i = IntOfExt(e);
       if (low < fIntParameters[i] && fIntParameters[i] < up)
@@ -395,7 +371,6 @@ void MnUserParameterState::SetUpperLimit(unsigned int e, double up)
    // set upper limit for parameter e (external index)
    fParameters.SetUpperLimit(e, up);
    fCovarianceValid = false;
-   fGCCValid = false;
    if (!Parameter(e).IsFixed() && !Parameter(e).IsConst()) {
       unsigned int i = IntOfExt(e);
       if (fIntParameters[i] < up)
@@ -410,7 +385,6 @@ void MnUserParameterState::SetLowerLimit(unsigned int e, double low)
    // set lower limit for parameter e (external index)
    fParameters.SetLowerLimit(e, low);
    fCovarianceValid = false;
-   fGCCValid = false;
    if (!Parameter(e).IsFixed() && !Parameter(e).IsConst()) {
       unsigned int i = IntOfExt(e);
       if (low < fIntParameters[i])
@@ -425,7 +399,6 @@ void MnUserParameterState::RemoveLimits(unsigned int e)
    // remove limit for parameter e (external index)
    fParameters.RemoveLimits(e);
    fCovarianceValid = false;
-   fGCCValid = false;
    if (!Parameter(e).IsFixed() && !Parameter(e).IsConst())
       fIntParameters[IntOfExt(e)] = Value(e);
 }
@@ -549,6 +522,21 @@ void MnUserParameterState::SetPrecision(double eps)
    fParameters.SetPrecision(eps);
 }
 
+MnGlobalCorrelationCoeff MnUserParameterState::GlobalCC() const
+{
+   if (!fCovarianceValid) {
+      return MnGlobalCorrelationCoeff{}; // invalidate
+   }
+   int n = fIntCovariance.Nrow();
+   MnAlgebraicSymMatrix invHessian{static_cast<unsigned int>(n)};
+   for (int i = 0; i < n; ++i) {
+      for (int j = 0; j <= i; ++j) {
+         invHessian(i, j) = fIntCovariance(i, j);
+      }
+   }
+   return MnGlobalCorrelationCoeff{invHessian};
+}
+
 } // namespace Minuit2
 
 } // namespace ROOT