diff options
author | André Nusser <andre.nusser@googlemail.com> | 2020-04-25 23:03:57 +0200 |
---|---|---|
committer | André Nusser <andre.nusser@googlemail.com> | 2020-04-25 23:03:57 +0200 |
commit | 4c616cf4a001f3b4adda87889a98e403c2e8ab8a (patch) | |
tree | e57bde3e6bd00be4f1554f99d0f670df53ee1b6a /src | |
parent | 8d22262afd0f91b49cb046d090b186c4305d470b (diff) |
Just leave the Steffen method in and clean up powermap.
Diffstat (limited to 'src')
-rw-r--r-- | src/powermap.cc | 78 | ||||
-rw-r--r-- | src/powermap.h | 7 |
2 files changed, 4 insertions, 81 deletions
diff --git a/src/powermap.cc b/src/powermap.cc index ca17a60..914b5c1 100644 --- a/src/powermap.cc +++ b/src/powermap.cc @@ -179,84 +179,10 @@ void Powermap::updateSpline() spline_needs_update = false; } +// This follows the monotone cubic spline algorithm of Steffen, from: +// "A Simple Method for Monotonic Interpolation in One Dimension" std::vector<float> Powermap::calcSlopes(Powers const& X, Powers const& Y) { - // return calcSlopesFritschCarlson(X, Y); - return calcSlopesSteffen(X, Y); -} - -std::vector<float> Powermap::calcSlopesFritschCarlson(Powers const& X, Powers const& Y) -{ - Powers deltas(X.size()-1); - Powers m(X.size()); - - // 1 - for (std::size_t i = 0; i < deltas.size(); ++i) { - deltas[i] = (Y[i+1]-Y[i])/(X[i+1]-X[i]); - } - - // 2 - m[0] = deltas[0]; - for (std::size_t i = 1; i < deltas.size(); ++i) { - m[i] = (deltas[i-1] + deltas[i])/2.; - } - m.back() = deltas.back(); - - // 3 - std::vector<bool> ignore(deltas.size(), false); - for (std::size_t i = 0; i < deltas.size(); ++i) { - if (deltas[i] == 0) { - m[i] = 0; - m[i+1] = 0; - ignore[i] = true; - } - } - - for (std::size_t i = 0; i < deltas.size(); ++i) { - if (ignore[i]) { - continue; - } - - // 4 - auto alpha = m[i]/deltas[i]; - auto beta = m[i+1]/deltas[i]; - assert(alpha >= 0.); - assert(beta >= 0.); - - // 5 - // TODO: expose this parameter for testing both - bool const option1 = false; - if (option1) { - if (alpha > 3 || beta > 3) { - m[i] = 3*deltas[i]; - } - } - else { - auto const a2b2 = alpha*alpha + beta*beta; - - // hard change to enforce monotonicity - // if (a2b2 > 9) { - // auto const tau = 3./sqrt(a2b2); - // m[i] = tau*alpha*deltas[i]; - // m[i+1] = tau*alpha*deltas[i]; - // } - - // soft change to enforce monotonicity - Power threshold = 4.5; // must be >= 0 and < 9. - if (a2b2 >= threshold) { - auto const l = std::min(1., (a2b2-threshold)/(9.-threshold)); - auto const tau = 3./sqrt(a2b2); - m[i] = (1-l)*m[i] + l*tau*alpha*deltas[i]; - m[i+1] = (1-l)*m[i+1] + l*tau*alpha*deltas[i]; - } - } - } - - return m; -} - -std::vector<float> Powermap::calcSlopesSteffen(Powers const& X, Powers const& Y) -{ Powers m(X.size()); Powers d(X.size()-1); diff --git a/src/powermap.h b/src/powermap.h index b942722..2cb7efb 100644 --- a/src/powermap.h +++ b/src/powermap.h @@ -63,15 +63,12 @@ private: bool shelf; // spline parameters (deterministically computed from the input parameters) - std::array<float, 5> m; bool spline_needs_update; + std::array<float, 5> m; + const Power eps = 1e-3; void updateSpline(); std::vector<float> calcSlopes(Powers const& X, Powers const& P); - std::vector<float> calcSlopesFritschCarlson(Powers const& X, Powers const& P); - std::vector<float> calcSlopesSteffen(Powers const& X, Powers const& P); Power clamp(Power in, Power min, Power max) const; - - const Power eps = 1e-3; }; |