#pragma once #include "StarInterpolation.hpp" namespace Star { // Describes a simple table from index to value, which operates on bins // corresponding to ranges of indexes. IndexType can be any ordered type, // ValueType can be anything. template class ParametricTable { public: typedef IndexType Index; typedef ValueType Value; ParametricTable(); template explicit ParametricTable(ParametricTable const& parametricFunction); // Construct a ParametricTable with a list of point pairs, which does not // have to be sorted (it will be sorted internally). Throws an exception on // duplicate index values. template explicit ParametricTable(PairContainer indexValuePairs); // addPoint does not need to be called in order, it will insert the point in // the correct ordered position for the given index, and return the position. size_t addPoint(IndexType index, ValueType value); void clearPoints(); size_t size() const; bool empty() const; IndexType const& index(size_t i) const; ValueType const& value(size_t i) const; // Returns true if the values of the table are also valid indexes (true when // the data points are monotonic increasing) bool isInvertible() const; // Invert the table, switching indexes and values. Throws an exception if // the function is not invertible. Will not generally compile unless the // Index and Value types are the same type. void invert() const; // Find the value to the left of the given index. If the index is lower than // the lowest index point, returns the first value. ValueType const& get(IndexType index) const; protected: typedef std::vector IndexList; typedef std::vector ValueList; IndexList const& indexes() const; ValueList const& values() const; private: IndexList m_indexes; ValueList m_values; }; // Extension of ParametricTable that simplifies of all the complex // interpolation code for interpolating an ordered list of points. Useful for // describing a simple 2d or n-dimensional (using VectorN for value type) curve // of one variable. IndexType should generally be float or double, and // ValueType can be any type that can be interpolated. template class ParametricFunction : public ParametricTable { public: typedef ParametricTable Base; ParametricFunction( InterpolationMode interpolationMode = InterpolationMode::Linear, BoundMode boundMode = BoundMode::Clamp); template explicit ParametricFunction(ParametricFunction const& parametricFunction); template explicit ParametricFunction(PairContainer indexValuePairs, InterpolationMode interpolationMode = InterpolationMode::Linear, BoundMode boundMode = BoundMode::Clamp); InterpolationMode interpolationMode() const; void setInterpolationMode(InterpolationMode interpolationType); BoundMode boundMode() const; void setBoundMode(BoundMode boundMode); // Interpolates a value at the given index according to the interpolation and // bound mode. ValueType interpolate(IndexType index) const; // Synonym for interpolate ValueType operator()(IndexType index) const; private: InterpolationMode m_interpolationMode; BoundMode m_boundMode; }; template ParametricTable::ParametricTable() {} template template ParametricTable::ParametricTable( ParametricTable const& parametricTable) { for (size_t i = 0; i < parametricTable.size(); ++i) { m_indexes.push_back(parametricTable.index(i)); m_values.push_back(parametricTable.value(i)); } } template template ParametricTable::ParametricTable(PairContainer indexValuePairs) { if (indexValuePairs.empty()) return; sort(indexValuePairs, [](typename PairContainer::value_type const& a, typename PairContainer::value_type const& b) { return std::get<0>(a) < std::get<0>(b); }); for (auto const& pair : indexValuePairs) { m_indexes.push_back(std::move(std::get<0>(pair))); m_values.push_back(std::move(std::get<1>(pair))); } for (size_t i = 0; i < size() - 1; ++i) { if (m_indexes[i] == m_indexes[i + 1]) throw MathException("Degenerate index values given in ParametricTable constructor"); } } template size_t ParametricTable::addPoint(IndexType index, ValueType value) { size_t insertLocation = std::distance(m_indexes.begin(), std::upper_bound(m_indexes.begin(), m_indexes.end(), index)); m_indexes.insert(m_indexes.begin() + insertLocation, std::move(index)); m_values.insert(m_values.begin() + insertLocation, std::move(value)); return insertLocation; } template void ParametricTable::clearPoints() { m_indexes.clear(); m_values.clear(); } template size_t ParametricTable::size() const { return m_indexes.size(); } template bool ParametricTable::empty() const { return size() == 0; } template IndexType const& ParametricTable::index(size_t i) const { return m_indexes.at(i); } template ValueType const& ParametricTable::value(size_t i) const { return m_values.at(i); } template bool ParametricTable::isInvertible() const { if (empty()) return true; for (size_t i = 0; i < size() - 1; ++i) { if (m_values[i] > m_values[i + 1]) return false; } return true; } template void ParametricTable::invert() const { if (isInvertible()) throw MathException("invert() called on non-invertible ParametricTable"); for (size_t i = 0; i < size(); ++i) std::swap(m_indexes[i], m_values[i]); } template ValueType const& ParametricTable::get(IndexType index) const { if (empty()) throw MathException("get called on empty ParametricTable"); auto i = std::lower_bound(m_indexes.begin(), m_indexes.end(), index); if (i != m_indexes.begin()) --i; return m_values[std::distance(m_indexes.begin(), i)]; } template auto ParametricTable::indexes() const -> IndexList const & { return m_indexes; } template auto ParametricTable::values() const -> ValueList const & { return m_values; } template ParametricFunction::ParametricFunction(InterpolationMode interpolationMode, BoundMode boundMode) : m_interpolationMode(interpolationMode), m_boundMode(boundMode) {} template template ParametricFunction::ParametricFunction( ParametricFunction const& parametricFunction) : Base(parametricFunction) { m_interpolationMode = parametricFunction.interpolationMode(); m_boundMode = parametricFunction.boundMode(); } template template ParametricFunction::ParametricFunction( PairContainer indexValuePairs, InterpolationMode interpolationMode, BoundMode boundMode) : Base(indexValuePairs) { m_interpolationMode = interpolationMode; m_boundMode = boundMode; } template InterpolationMode ParametricFunction::interpolationMode() const { return m_interpolationMode; } template void ParametricFunction::setInterpolationMode(InterpolationMode interpolationType) { m_interpolationMode = interpolationType; } template BoundMode ParametricFunction::boundMode() const { return m_boundMode; } template void ParametricFunction::setBoundMode(BoundMode boundMode) { m_boundMode = boundMode; } template ValueType ParametricFunction::interpolate(IndexType index) const { if (Base::empty()) return ValueType(); if (m_interpolationMode == InterpolationMode::HalfStep) { return parametricInterpolate2(Base::indexes(), Base::values(), index, StepWeightOperator(), m_boundMode); } else if (m_interpolationMode == InterpolationMode::Linear) { return parametricInterpolate2( Base::indexes(), Base::values(), index, LinearWeightOperator(), m_boundMode); } else if (m_interpolationMode == InterpolationMode::Cubic) { // ParametricFunction uses CubicWeights with linear extrapolation (not // configurable atm) return parametricInterpolate4( Base::indexes(), Base::values(), index, Cubic4WeightOperator(true), m_boundMode); } else { throw MathException("Unsupported interpolation type in ParametricFunction::interpolate"); } } template ValueType ParametricFunction::operator()(IndexType index) const { return interpolate(index); } }