Data approximation

In order to build \(r\), one can either use a pre-defined sequence of control points \(\{ c[k] \}_{k=0,...,M-1}\) or of knots \(\{ n[k] \}_{k=0,...,M-1}\). Alternatively, one can also attempt to retreive the control points that best approximate a set of data points, as in the classical spline approximation setting.

The problem is framed as follows. We consider a set \(\{ p[i] \}_{i=0,...,N-1}\) of \(N\) points to be approximated with the spline model (1) of \(M\) control points. Hereafter, we will assume a periodic spline model, but a similar derivation can easily be done for the non-periodic case.

We obtain an approximation by ensuring that the samples of the spline \(r\) match the data points \(p\), which translates to

(1)\[p[i] = \sum_{k=0}^{M-1}c[k]\phi\left(\frac{Mi}{N}-k\right).\]

Since \(\phi\) is of finite support, we can re-write (1) as

(2)\[\mathbf{\Phi}\mathbf{C} = \mathbf{P},\]

with the basis matrix \(\mathbf{\Phi}\) (size \(N \times M\)), the control point matrix \(\mathbf{C}\) (size \(M \times 1\)), and the data points matrix \(\mathbf{P}\) (size \(N \times 1\)) given by

(3)\[\begin{split}\mathbf{\Phi} = \begin{bmatrix} \phi(0) & \phi(-1) & \dots & \ \phi(-(M-1)) \\ \phi\left(\frac{M}{N}\right) & \phi\left(\frac{M}{N}-1\right) & \dots & \ \phi\left(\frac{M}{N}-(M-1)\right) \\ \vdots & \vdots & \ddots & \vdots \\ \phi\left(\frac{(N-1)M}{N}\right) & \phi\left(\frac{(N-1)M}{N}-1\right) & \dots & \ \phi\left(\frac{(N-1)M}{N}-(M-1)\right) \end{bmatrix}\end{split}\]

(4)\[\begin{split}\mathbf{C} = \begin{bmatrix} c[0] \\ \vdots \\ c[M-1] \end{bmatrix}\end{split}\]

(5)\[\begin{split}\mathbf{P} = \begin{bmatrix} p[0] \\ \vdots \\ p[N-1] \end{bmatrix}.\end{split}\]

The control points \(\mathbf{C}\) can then be retrieved by finding the least-square best solution that minimizes

(6)\[\| \mathbf{P} - \mathbf{\Phi} \mathbf{C} \|^2_2.\]

Boundary conditions

Open splines are padded with additional control points at the ends (Padding). They are located at parameters values \(t = -1, -2, \ldots\) and \(t = M, M+1, \ldots\). Since the data points are only positions on the parameters interval \([0, M-1]\), the additional control points can cause erratic behavior of the spline when fitting noisy data. To control this behaviour, one of the following boundary conditions can be enforced:

• clamped: \(r'(0)=0\) and \(r'(M-1)=0\)

• natural: \(r''(0)=0\) and \(r''(M-1)=0\)

With out loss of generality, we will describe how to fit a spline with the clamped boundary condition.

The spline can be written as:

(7)\[r(t) = \sum_{k=-p}^{M-1+p}c[k]\Phi(t-k),\]

where \(p\) is the amount of padding.

The boundary condition allows us to express the first control point as a combination of the other control points:

\[\begin{split}r'(0) &= 0 \\ 0 &= \sum_{k=-p}^{M-1+p}c[k]\Phi'(-k) \\ c[-p] &= \sum_{k=-p+1}^{M-1+p} -c[k] \frac{\Phi'(-k)}{\Phi'(p)}\end{split}\]

Note: We assume that \(\Phi'(p) \neq 0\).

Plugging this into equation (7) and grouping the summand by control point yields:

(8)\[r(t) = \sum_{k=-p+1}^{M-1+p} c[k] (\Phi(t-k) - \frac{\Phi'(-k)}{\Phi'(p)}\Phi(t+p))\]

Note: \(\Phi(t+p)\) will be zero for most \(t\) since the support of the basis function \(\Phi\) end in the interval \([p, p+1)\).

We can do the same for the last control point, starting from equation (8):

\[\begin{split}r'(M-1) &= 0 \\ 0 &= \sum_{k=-p+1}^{M-1+p} c[k] (\Phi'(M-1-k) - \frac{\Phi'(-k)}{\Phi'(p)}\Phi'(M-1+p)) \\\end{split}\]

Because of the way we choose \(p\) and the fact that \(M > 1\), we know that is outside the support of \(\Phi\) and \(\Phi'(M-1+p)=0\). Therefor we get:

\[\begin{split} 0 &= \sum_{k=-p+1}^{M-1+p} c[k] \Phi'(M-1-k) \\ c[M-1+p] &= \sum_{k=-p+1}^{M-2+p} -c[k] \frac{\Phi'(M-1-k)}{\Phi'(-p)}\end{split}\]

Plugging this into equation (8) yields:

\[r(t) = \sum_{k=-p+1}^{M-2+p} c[k] (\Phi(t-k) - \frac{\Phi'(-k)}{\Phi'(p)}\Phi(t+p) - \frac{\Phi'(M-1-k)}{\Phi'(-p)}\Phi(t-M+1-p))\]

Like equation (1), this equation can be written as a matrix-vector multiplication and can be solved using least-squares.