1c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// This file is part of Eigen, a lightweight C++ template library 2c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// for linear algebra. 3c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// 4c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// Copyright (C) 20010-2011 Hauke Heibel <hauke.heibel@gmail.com> 5c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// 6c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// This Source Code Form is subject to the terms of the Mozilla 7c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// Public License v. 2.0. If a copy of the MPL was not distributed 8c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 9c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 10c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath#ifndef EIGEN_SPLINE_H 11c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath#define EIGEN_SPLINE_H 12c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 13c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath#include "SplineFwd.h" 14c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 15c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathnamespace Eigen 16c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath{ 17c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** 18c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \ingroup Splines_Module 19c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \class Spline class 20c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \brief A class representing multi-dimensional spline curves. 21c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 22c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * The class represents B-splines with non-uniform knot vectors. Each control 23c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * point of the B-spline is associated with a basis function 24c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \f{align*} 25c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * C(u) & = \sum_{i=0}^{n}N_{i,p}(u)P_i 26c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \f} 27c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 28c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \tparam _Scalar The underlying data type (typically float or double) 29c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \tparam _Dim The curve dimension (e.g. 2 or 3) 30c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \tparam _Degree Per default set to Dynamic; could be set to the actual desired 31c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * degree for optimization purposes (would result in stack allocation 32c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * of several temporary variables). 33c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath **/ 34c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <typename _Scalar, int _Dim, int _Degree> 35c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath class Spline 36c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 37c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath public: 38c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef _Scalar Scalar; /*!< The spline curve's scalar type. */ 39c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath enum { Dimension = _Dim /*!< The spline curve's dimension. */ }; 40c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath enum { Degree = _Degree /*!< The spline curve's degree. */ }; 41c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 42c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** \brief The point type the spline is representing. */ 43c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename SplineTraits<Spline>::PointType PointType; 44c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 45c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** \brief The data type used to store knot vectors. */ 46c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename SplineTraits<Spline>::KnotVectorType KnotVectorType; 47c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 48c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** \brief The data type used to store non-zero basis functions. */ 49c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename SplineTraits<Spline>::BasisVectorType BasisVectorType; 50c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 51c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** \brief The data type representing the spline's control points. */ 52c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename SplineTraits<Spline>::ControlPointVectorType ControlPointVectorType; 53c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 54c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** 55c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \brief Creates a spline from a knot vector and control points. 56c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \param knots The spline's knot vector. 57c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \param ctrls The spline's control point vector. 58c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath **/ 59c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <typename OtherVectorType, typename OtherArrayType> 60c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Spline(const OtherVectorType& knots, const OtherArrayType& ctrls) : m_knots(knots), m_ctrls(ctrls) {} 61c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 62c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** 63c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \brief Copy constructor for splines. 64c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \param spline The input spline. 65c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath **/ 66c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <int OtherDegree> 67c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Spline(const Spline<Scalar, Dimension, OtherDegree>& spline) : 68c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_knots(spline.knots()), m_ctrls(spline.ctrls()) {} 69c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 70c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** 71c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \brief Returns the knots of the underlying spline. 72c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath **/ 73c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const KnotVectorType& knots() const { return m_knots; } 74c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 75c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** 76c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \brief Returns the knots of the underlying spline. 77c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath **/ 78c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const ControlPointVectorType& ctrls() const { return m_ctrls; } 79c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 80c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** 81c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \brief Returns the spline value at a given site \f$u\f$. 82c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 83c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * The function returns 84c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \f{align*} 85c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * C(u) & = \sum_{i=0}^{n}N_{i,p}P_i 86c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \f} 87c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 88c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \param u Parameter \f$u \in [0;1]\f$ at which the spline is evaluated. 89c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \return The spline value at the given location \f$u\f$. 90c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath **/ 91c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath PointType operator()(Scalar u) const; 92c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 93c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** 94c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \brief Evaluation of spline derivatives of up-to given order. 95c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 96c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * The function returns 97c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \f{align*} 98c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \frac{d^i}{du^i}C(u) & = \sum_{i=0}^{n} \frac{d^i}{du^i} N_{i,p}(u)P_i 99c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \f} 100c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * for i ranging between 0 and order. 101c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 102c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \param u Parameter \f$u \in [0;1]\f$ at which the spline derivative is evaluated. 103c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \param order The order up to which the derivatives are computed. 104c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath **/ 105c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename SplineTraits<Spline>::DerivativeType 106c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath derivatives(Scalar u, DenseIndex order) const; 107c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 108c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** 109c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \copydoc Spline::derivatives 110c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * Using the template version of this function is more efficieent since 111c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * temporary objects are allocated on the stack whenever this is possible. 112c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath **/ 113c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <int DerivativeOrder> 114c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename SplineTraits<Spline,DerivativeOrder>::DerivativeType 115c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath derivatives(Scalar u, DenseIndex order = DerivativeOrder) const; 116c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 117c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** 118c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \brief Computes the non-zero basis functions at the given site. 119c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 120c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * Splines have local support and a point from their image is defined 121c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * by exactly \f$p+1\f$ control points \f$P_i\f$ where \f$p\f$ is the 122c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * spline degree. 123c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 124c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * This function computes the \f$p+1\f$ non-zero basis function values 125c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * for a given parameter value \f$u\f$. It returns 126c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \f{align*}{ 127c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * N_{i,p}(u), \hdots, N_{i+p+1,p}(u) 128c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \f} 129c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 130c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \param u Parameter \f$u \in [0;1]\f$ at which the non-zero basis functions 131c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * are computed. 132c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath **/ 133c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename SplineTraits<Spline>::BasisVectorType 134c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath basisFunctions(Scalar u) const; 135c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 136c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** 137c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \brief Computes the non-zero spline basis function derivatives up to given order. 138c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 139c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * The function computes 140c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \f{align*}{ 141c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \frac{d^i}{du^i} N_{i,p}(u), \hdots, \frac{d^i}{du^i} N_{i+p+1,p}(u) 142c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \f} 143c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * with i ranging from 0 up to the specified order. 144c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 145c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \param u Parameter \f$u \in [0;1]\f$ at which the non-zero basis function 146c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * derivatives are computed. 147c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \param order The order up to which the basis function derivatives are computes. 148c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath **/ 149c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename SplineTraits<Spline>::BasisDerivativeType 150c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath basisFunctionDerivatives(Scalar u, DenseIndex order) const; 151c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 152c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** 153c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \copydoc Spline::basisFunctionDerivatives 154c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * Using the template version of this function is more efficieent since 155c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * temporary objects are allocated on the stack whenever this is possible. 156c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath **/ 157c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <int DerivativeOrder> 158c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename SplineTraits<Spline,DerivativeOrder>::BasisDerivativeType 159c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath basisFunctionDerivatives(Scalar u, DenseIndex order = DerivativeOrder) const; 160c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 161c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** 162c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \brief Returns the spline degree. 163c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath **/ 164c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath DenseIndex degree() const; 165c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 166c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** 167c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \brief Returns the span within the knot vector in which u is falling. 168c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \param u The site for which the span is determined. 169c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath **/ 170c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath DenseIndex span(Scalar u) const; 171c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 172c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** 173c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \brief Computes the spang within the provided knot vector in which u is falling. 174c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath **/ 175c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath static DenseIndex Span(typename SplineTraits<Spline>::Scalar u, DenseIndex degree, const typename SplineTraits<Spline>::KnotVectorType& knots); 176c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 177c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** 178c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \brief Returns the spline's non-zero basis functions. 179c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 180c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * The function computes and returns 181c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \f{align*}{ 182c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * N_{i,p}(u), \hdots, N_{i+p+1,p}(u) 183c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \f} 184c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 185c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \param u The site at which the basis functions are computed. 186c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \param degree The degree of the underlying spline. 187c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \param knots The underlying spline's knot vector. 188c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath **/ 189c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath static BasisVectorType BasisFunctions(Scalar u, DenseIndex degree, const KnotVectorType& knots); 190c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 191c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 192c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath private: 193c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath KnotVectorType m_knots; /*!< Knot vector. */ 194c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath ControlPointVectorType m_ctrls; /*!< Control points. */ 195c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath }; 196c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 197c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <typename _Scalar, int _Dim, int _Degree> 198c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath DenseIndex Spline<_Scalar, _Dim, _Degree>::Span( 199c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename SplineTraits< Spline<_Scalar, _Dim, _Degree> >::Scalar u, 200c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath DenseIndex degree, 201c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const typename SplineTraits< Spline<_Scalar, _Dim, _Degree> >::KnotVectorType& knots) 202c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 203c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // Piegl & Tiller, "The NURBS Book", A2.1 (p. 68) 204c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (u <= knots(0)) return degree; 205c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const Scalar* pos = std::upper_bound(knots.data()+degree-1, knots.data()+knots.size()-degree-1, u); 206c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return static_cast<DenseIndex>( std::distance(knots.data(), pos) - 1 ); 207c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 208c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 209c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <typename _Scalar, int _Dim, int _Degree> 210c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename Spline<_Scalar, _Dim, _Degree>::BasisVectorType 211c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Spline<_Scalar, _Dim, _Degree>::BasisFunctions( 212c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename Spline<_Scalar, _Dim, _Degree>::Scalar u, 213c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath DenseIndex degree, 214c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const typename Spline<_Scalar, _Dim, _Degree>::KnotVectorType& knots) 215c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 216c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename Spline<_Scalar, _Dim, _Degree>::BasisVectorType BasisVectorType; 217c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 218c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const DenseIndex p = degree; 219c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const DenseIndex i = Spline::Span(u, degree, knots); 220c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 221c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const KnotVectorType& U = knots; 222c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 223c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath BasisVectorType left(p+1); left(0) = Scalar(0); 224c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath BasisVectorType right(p+1); right(0) = Scalar(0); 225c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 226c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath VectorBlock<BasisVectorType,Degree>(left,1,p) = u - VectorBlock<const KnotVectorType,Degree>(U,i+1-p,p).reverse(); 227c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath VectorBlock<BasisVectorType,Degree>(right,1,p) = VectorBlock<const KnotVectorType,Degree>(U,i+1,p) - u; 228c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 229c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath BasisVectorType N(1,p+1); 230c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath N(0) = Scalar(1); 231c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for (DenseIndex j=1; j<=p; ++j) 232c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 233c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Scalar saved = Scalar(0); 234c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for (DenseIndex r=0; r<j; r++) 235c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 236c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const Scalar tmp = N(r)/(right(r+1)+left(j-r)); 237c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath N[r] = saved + right(r+1)*tmp; 238c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath saved = left(j-r)*tmp; 239c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 240c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath N(j) = saved; 241c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 242c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return N; 243c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 244c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 245c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <typename _Scalar, int _Dim, int _Degree> 246c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath DenseIndex Spline<_Scalar, _Dim, _Degree>::degree() const 247c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 248c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (_Degree == Dynamic) 249c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return m_knots.size() - m_ctrls.cols() - 1; 250c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath else 251c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return _Degree; 252c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 253c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 254c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <typename _Scalar, int _Dim, int _Degree> 255c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath DenseIndex Spline<_Scalar, _Dim, _Degree>::span(Scalar u) const 256c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 257c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return Spline::Span(u, degree(), knots()); 258c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 259c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 260c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <typename _Scalar, int _Dim, int _Degree> 261c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename Spline<_Scalar, _Dim, _Degree>::PointType Spline<_Scalar, _Dim, _Degree>::operator()(Scalar u) const 262c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 263c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath enum { Order = SplineTraits<Spline>::OrderAtCompileTime }; 264c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 265c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const DenseIndex span = this->span(u); 266c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const DenseIndex p = degree(); 267c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const BasisVectorType basis_funcs = basisFunctions(u); 268c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 269c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const Replicate<BasisVectorType,Dimension,1> ctrl_weights(basis_funcs); 270c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const Block<const ControlPointVectorType,Dimension,Order> ctrl_pts(ctrls(),0,span-p,Dimension,p+1); 271c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return (ctrl_weights * ctrl_pts).rowwise().sum(); 272c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 273c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 274c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /* --------------------------------------------------------------------------------------------- */ 275c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 276c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <typename SplineType, typename DerivativeType> 277c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void derivativesImpl(const SplineType& spline, typename SplineType::Scalar u, DenseIndex order, DerivativeType& der) 278c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 279c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath enum { Dimension = SplineTraits<SplineType>::Dimension }; 280c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath enum { Order = SplineTraits<SplineType>::OrderAtCompileTime }; 281c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath enum { DerivativeOrder = DerivativeType::ColsAtCompileTime }; 282c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 283c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename SplineTraits<SplineType>::Scalar Scalar; 284c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 285c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename SplineTraits<SplineType>::BasisVectorType BasisVectorType; 286c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename SplineTraits<SplineType>::ControlPointVectorType ControlPointVectorType; 287c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 288c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename SplineTraits<SplineType,DerivativeOrder>::BasisDerivativeType BasisDerivativeType; 289c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename BasisDerivativeType::ConstRowXpr BasisDerivativeRowXpr; 290c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 291c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const DenseIndex p = spline.degree(); 292c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const DenseIndex span = spline.span(u); 293c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 294c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const DenseIndex n = (std::min)(p, order); 295c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 296c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath der.resize(Dimension,n+1); 297c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 298c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // Retrieve the basis function derivatives up to the desired order... 299c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const BasisDerivativeType basis_func_ders = spline.template basisFunctionDerivatives<DerivativeOrder>(u, n+1); 300c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 301c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // ... and perform the linear combinations of the control points. 302c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for (DenseIndex der_order=0; der_order<n+1; ++der_order) 303c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 304c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const Replicate<BasisDerivativeRowXpr,Dimension,1> ctrl_weights( basis_func_ders.row(der_order) ); 305c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const Block<const ControlPointVectorType,Dimension,Order> ctrl_pts(spline.ctrls(),0,span-p,Dimension,p+1); 306c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath der.col(der_order) = (ctrl_weights * ctrl_pts).rowwise().sum(); 307c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 308c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 309c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 310c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <typename _Scalar, int _Dim, int _Degree> 311c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename SplineTraits< Spline<_Scalar, _Dim, _Degree> >::DerivativeType 312c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Spline<_Scalar, _Dim, _Degree>::derivatives(Scalar u, DenseIndex order) const 313c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 314c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename SplineTraits< Spline >::DerivativeType res; 315c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath derivativesImpl(*this, u, order, res); 316c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return res; 317c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 318c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 319c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <typename _Scalar, int _Dim, int _Degree> 320c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <int DerivativeOrder> 321c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename SplineTraits< Spline<_Scalar, _Dim, _Degree>, DerivativeOrder >::DerivativeType 322c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Spline<_Scalar, _Dim, _Degree>::derivatives(Scalar u, DenseIndex order) const 323c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 324c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename SplineTraits< Spline, DerivativeOrder >::DerivativeType res; 325c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath derivativesImpl(*this, u, order, res); 326c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return res; 327c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 328c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 329c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <typename _Scalar, int _Dim, int _Degree> 330c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename SplineTraits< Spline<_Scalar, _Dim, _Degree> >::BasisVectorType 331c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Spline<_Scalar, _Dim, _Degree>::basisFunctions(Scalar u) const 332c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 333c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return Spline::BasisFunctions(u, degree(), knots()); 334c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 335c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 336c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /* --------------------------------------------------------------------------------------------- */ 337c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 338c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <typename SplineType, typename DerivativeType> 339c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void basisFunctionDerivativesImpl(const SplineType& spline, typename SplineType::Scalar u, DenseIndex order, DerivativeType& N_) 340c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 341c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath enum { Order = SplineTraits<SplineType>::OrderAtCompileTime }; 342c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 343c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename SplineTraits<SplineType>::Scalar Scalar; 344c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename SplineTraits<SplineType>::BasisVectorType BasisVectorType; 345c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename SplineTraits<SplineType>::KnotVectorType KnotVectorType; 346c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename SplineTraits<SplineType>::ControlPointVectorType ControlPointVectorType; 347c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 348c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const KnotVectorType& U = spline.knots(); 349c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 350c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const DenseIndex p = spline.degree(); 351c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const DenseIndex span = spline.span(u); 352c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 353c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const DenseIndex n = (std::min)(p, order); 354c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 355c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath N_.resize(n+1, p+1); 356c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 357c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath BasisVectorType left = BasisVectorType::Zero(p+1); 358c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath BasisVectorType right = BasisVectorType::Zero(p+1); 359c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 360c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Matrix<Scalar,Order,Order> ndu(p+1,p+1); 361c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 362c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath double saved, temp; 363c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 364c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath ndu(0,0) = 1.0; 365c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 366c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath DenseIndex j; 367c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for (j=1; j<=p; ++j) 368c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 369c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath left[j] = u-U[span+1-j]; 370c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath right[j] = U[span+j]-u; 371c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath saved = 0.0; 372c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 373c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for (DenseIndex r=0; r<j; ++r) 374c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 375c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /* Lower triangle */ 376c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath ndu(j,r) = right[r+1]+left[j-r]; 377c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath temp = ndu(r,j-1)/ndu(j,r); 378c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /* Upper triangle */ 379c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath ndu(r,j) = static_cast<Scalar>(saved+right[r+1] * temp); 380c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath saved = left[j-r] * temp; 381c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 382c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 383c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath ndu(j,j) = static_cast<Scalar>(saved); 384c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 385c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 386c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for (j = p; j>=0; --j) 387c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath N_(0,j) = ndu(j,p); 388c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 389c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // Compute the derivatives 390c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath DerivativeType a(n+1,p+1); 391c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath DenseIndex r=0; 392c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for (; r<=p; ++r) 393c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 394c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath DenseIndex s1,s2; 395c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath s1 = 0; s2 = 1; // alternate rows in array a 396c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath a(0,0) = 1.0; 397c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 398c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // Compute the k-th derivative 399c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for (DenseIndex k=1; k<=static_cast<DenseIndex>(n); ++k) 400c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 401c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath double d = 0.0; 402c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath DenseIndex rk,pk,j1,j2; 403c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath rk = r-k; pk = p-k; 404c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 405c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (r>=k) 406c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 407c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath a(s2,0) = a(s1,0)/ndu(pk+1,rk); 408c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath d = a(s2,0)*ndu(rk,pk); 409c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 410c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 411c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (rk>=-1) j1 = 1; 412c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath else j1 = -rk; 413c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 414c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (r-1 <= pk) j2 = k-1; 415c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath else j2 = p-r; 416c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 417c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for (j=j1; j<=j2; ++j) 418c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 419c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath a(s2,j) = (a(s1,j)-a(s1,j-1))/ndu(pk+1,rk+j); 420c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath d += a(s2,j)*ndu(rk+j,pk); 421c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 422c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 423c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if (r<=pk) 424c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 425c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath a(s2,k) = -a(s1,k-1)/ndu(pk+1,r); 426c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath d += a(s2,k)*ndu(r,pk); 427c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 428c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 429c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath N_(k,r) = static_cast<Scalar>(d); 430c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath j = s1; s1 = s2; s2 = j; // Switch rows 431c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 432c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 433c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 434c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /* Multiply through by the correct factors */ 435c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /* (Eq. [2.9]) */ 436c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath r = p; 437c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for (DenseIndex k=1; k<=static_cast<DenseIndex>(n); ++k) 438c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 439c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for (DenseIndex j=p; j>=0; --j) N_(k,j) *= r; 440c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath r *= p-k; 441c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 442c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 443c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 444c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <typename _Scalar, int _Dim, int _Degree> 445c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename SplineTraits< Spline<_Scalar, _Dim, _Degree> >::BasisDerivativeType 446c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Spline<_Scalar, _Dim, _Degree>::basisFunctionDerivatives(Scalar u, DenseIndex order) const 447c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 448c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename SplineTraits< Spline >::BasisDerivativeType der; 449c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath basisFunctionDerivativesImpl(*this, u, order, der); 450c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return der; 451c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 452c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 453c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <typename _Scalar, int _Dim, int _Degree> 454c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <int DerivativeOrder> 455c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename SplineTraits< Spline<_Scalar, _Dim, _Degree>, DerivativeOrder >::BasisDerivativeType 456c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Spline<_Scalar, _Dim, _Degree>::basisFunctionDerivatives(Scalar u, DenseIndex order) const 457c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 458c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename SplineTraits< Spline, DerivativeOrder >::BasisDerivativeType der; 459c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath basisFunctionDerivativesImpl(*this, u, order, der); 460c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return der; 461c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 462c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath} 463c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 464c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath#endif // EIGEN_SPLINE_H 465