Skip to content

Commit cce19a8

Browse files
committed
[Math] Avoid using std::map in Math::FitResult and fix logic
Avoid using `std::map` in Math::FitResult and fix the logic to set whether a parameter is fixed or bound. Before, these flags were only changed from `false` to `true`, but not the other way around, which can also happen if the user releases a parameter in the configuration object. A unit test that covers this case was also implemented. Closes #20703. (cherry picked from commit d5085b1)
1 parent 309e8ec commit cce19a8

File tree

4 files changed

+111
-33
lines changed

4 files changed

+111
-33
lines changed

math/mathcore/inc/Fit/FitResult.h

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -354,8 +354,8 @@ class FitResult {
354354
std::shared_ptr<ROOT::Math::IMultiGenFunction> fObjFunc; ///<! objective function used for fitting
355355
std::shared_ptr<IModelFunction> fFitFunc; ///<! model function resulting from the fit.
356356
std::shared_ptr<FitData> fFitData; ///<! data set used in the fit
357-
std::map<unsigned int, bool> fFixedParams; ///< list of fixed parameters
358-
std::map<unsigned int, unsigned int> fBoundParams; ///< list of limited parameters
357+
std::vector<bool> fFixedParams; ///< if parameters are fixed
358+
std::vector<unsigned int> fBoundParams; ///< if parameters are limited
359359
std::vector<std::pair<double,double> > fParamBounds; ///< parameter bounds
360360
std::vector<double> fParams; ///< parameter values. Size is total number of parameters
361361
std::vector<double> fErrors; ///< errors

math/mathcore/src/FitResult.cxx

Lines changed: 30 additions & 31 deletions
Original file line numberDiff line numberDiff line change
@@ -79,19 +79,19 @@ FitResult::FitResult(const FitConfig & fconfig) :
7979

8080
// get parameter values and errors (step sizes)
8181
unsigned int npar = fconfig.NPar();
82+
fParamBounds.resize(npar);
83+
fFixedParams.resize(npar);
84+
fBoundParams.resize(npar);
8285
for (unsigned int i = 0; i < npar; ++i ) {
8386
const ParameterSettings & par = fconfig.ParSettings(i);
8487
fParams[i] = par.Value();
8588
fErrors[i] = par.StepSize();
8689
fParNames[i] = par.Name();
87-
if (par.IsFixed() ) fFixedParams[i] = true;
88-
else fNFree++;
89-
if (par.IsBound() ) {
90-
double lower = (par.HasLowerLimit()) ? par.LowerLimit() : - std::numeric_limits<double>::infinity() ;
91-
double upper = (par.HasUpperLimit()) ? par.UpperLimit() : std::numeric_limits<double>::infinity() ;
92-
fBoundParams[i] = fParamBounds.size();
93-
fParamBounds.push_back(std::make_pair(lower,upper));
94-
}
90+
fFixedParams[i] = par.IsFixed();
91+
if (!par.IsFixed() ) fNFree++;
92+
double lower = par.HasLowerLimit() ? par.LowerLimit() : - std::numeric_limits<double>::infinity() ;
93+
double upper = par.HasUpperLimit() ? par.UpperLimit() : std::numeric_limits<double>::infinity() ;
94+
fParamBounds.emplace_back(lower,upper);
9595
}
9696
std::cout << "create fit result from config - nfree " << fNFree << std::endl;
9797
}
@@ -155,17 +155,17 @@ void FitResult::FillResult(const std::shared_ptr<ROOT::Math::Minimizer> & min, c
155155

156156
// check for fixed or limited parameters
157157
unsigned int nfree = 0;
158-
if (!fParamBounds.empty()) fParamBounds.clear();
158+
fParamBounds.resize(npar);
159+
fFixedParams.resize(npar);
160+
fBoundParams.resize(npar);
159161
for (unsigned int ipar = 0; ipar < npar; ++ipar) {
160162
const ParameterSettings & par = fconfig.ParSettings(ipar);
161-
if (par.IsFixed() ) fFixedParams[ipar] = true;
162-
else nfree++;
163-
if (par.IsBound() ) {
164-
double lower = (par.HasLowerLimit()) ? par.LowerLimit() : - std::numeric_limits<double>::infinity() ;
165-
double upper = (par.HasUpperLimit()) ? par.UpperLimit() : std::numeric_limits<double>::infinity() ;
166-
fBoundParams[ipar] = fParamBounds.size();
167-
fParamBounds.push_back(std::make_pair(lower,upper));
168-
}
163+
fFixedParams[ipar] = par.IsFixed();
164+
fBoundParams[ipar] = par.IsBound();
165+
if (!par.IsFixed() ) nfree++;
166+
double lower = par.HasLowerLimit() ? par.LowerLimit() : - std::numeric_limits<double>::infinity() ;
167+
double upper = par.HasUpperLimit() ? par.UpperLimit() : std::numeric_limits<double>::infinity() ;
168+
fParamBounds.emplace_back(lower,upper);
169169
}
170170
// check if nfree (from FitConfig) and fNFree (from minimizer) are consistent
171171
if (nfree != fNFree ) {
@@ -343,25 +343,24 @@ int FitResult::Index(const std::string & name) const {
343343
return -1; // case name is not found
344344
}
345345

346-
bool FitResult::IsParameterBound(unsigned int ipar) const {
347-
return fBoundParams.find(ipar) != fBoundParams.end();
346+
bool FitResult::IsParameterBound(unsigned int ipar) const
347+
{
348+
return ipar < fBoundParams.size() ? fBoundParams[ipar] : false;
348349
}
349350

350-
bool FitResult::IsParameterFixed(unsigned int ipar) const {
351-
return fFixedParams.find(ipar) != fFixedParams.end();
351+
bool FitResult::IsParameterFixed(unsigned int ipar) const
352+
{
353+
return ipar < fFixedParams.size() ? fFixedParams[ipar] : false;
352354
}
353355

354-
bool FitResult::ParameterBounds(unsigned int ipar, double & lower, double & upper) const {
355-
std::map<unsigned int, unsigned int>::const_iterator itr = fBoundParams.find(ipar);
356-
if (itr == fBoundParams.end() ) {
357-
lower = -std::numeric_limits<Double_t>::infinity();
358-
upper = std::numeric_limits<Double_t>::infinity();
359-
return false;
356+
bool FitResult::ParameterBounds(unsigned int ipar, double &lower, double &upper) const
357+
{
358+
constexpr double inf = std::numeric_limits<double>::infinity();
359+
if (ipar < fParamBounds.size()) {
360+
lower = fParamBounds[ipar].first;
361+
upper = fParamBounds[ipar].second;
360362
}
361-
assert(itr->second < fParamBounds.size() );
362-
lower = fParamBounds[itr->second].first;
363-
upper = fParamBounds[itr->second].second;
364-
return true;
363+
return lower != -inf || upper != inf;
365364
}
366365

367366
std::string FitResult::ParName(unsigned int ipar) const {

math/mathcore/test/CMakeLists.txt

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -100,3 +100,5 @@ ROOT_ADD_GTEST(testKNNDensity testKNNDensity.cxx LIBRARIES Core MathCore Hist)
100100
if(clad)
101101
ROOT_ADD_GTEST(CladDerivatorTests CladDerivatorTests.cxx LIBRARIES Core MathCore)
102102
endif()
103+
104+
ROOT_ADD_GTEST(testFitter testFitter.cxx LIBRARIES Core MathCore)

math/mathcore/test/testFitter.cxx

Lines changed: 77 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,77 @@
1+
// Tests for the ROOT::Fit::Fitter
2+
3+
#include <Fit/Fitter.h>
4+
#include <Math/Functor.h>
5+
6+
#include <gtest/gtest.h>
7+
8+
namespace {
9+
10+
double gauss(double *x, double *p)
11+
{
12+
double A = p[0];
13+
double mu = p[1];
14+
double sigma = p[2];
15+
double xx = x[0];
16+
return A * exp(-(xx - mu) * (xx - mu) / (2 * sigma * sigma));
17+
}
18+
19+
} // namespace
20+
21+
/// Check if we can fix and release parameters, and that this will be correctly
22+
/// reflected in the FitResult.
23+
/// Covers https://github.com/root-project/root/issues/20703
24+
TEST(Fitter, FitAndReleaseParams)
25+
{
26+
// --- Create data ---
27+
const int n = 50;
28+
double x[n];
29+
double y[n];
30+
std::vector<double> initVals{1.0, 2.0, 1.0};
31+
for (int i = 0; i < n; ++i) {
32+
x[i] = i * 0.1;
33+
y[i] = gauss(&x[i], initVals.data());
34+
}
35+
36+
// --- Model functor ---
37+
ROOT::Math::Functor f(
38+
[&](const double *p) {
39+
double chi2 = 0.0;
40+
for (int i = 0; i < n; ++i) {
41+
double yi = gauss(&x[i], (double *)p);
42+
double diff = y[i] - yi;
43+
chi2 += diff * diff;
44+
}
45+
return chi2;
46+
},
47+
3);
48+
49+
ROOT::Fit::Fitter fitter;
50+
51+
std::vector<bool> fixState{false, false, false};
52+
fitter.Config().SetParamsSettings(initVals.size(), initVals.data());
53+
54+
fitter.Config().ParSettings(0).SetName("A");
55+
fitter.Config().ParSettings(1).SetName("mu");
56+
fitter.Config().ParSettings(2).SetName("sigma");
57+
58+
// Repeatedly run fits using the SAME fitter
59+
for (int iter = 0; iter < 6; ++iter) {
60+
bool fix = (iter % 2 == 0);
61+
int iparam = (iter / 2) % initVals.size();
62+
63+
fixState[iparam] = fix;
64+
if (fix) {
65+
fitter.Config().ParSettings(iparam).Fix();
66+
} else {
67+
fitter.Config().ParSettings(iparam).Release();
68+
}
69+
70+
fitter.FitFCN(f, nullptr, 3);
71+
const ROOT::Fit::FitResult &res = fitter.Result();
72+
73+
for (unsigned int i = 0; i < initVals.size(); ++i) {
74+
EXPECT_EQ(res.IsParameterFixed(i), fixState[i]);
75+
}
76+
}
77+
}

0 commit comments

Comments
 (0)