Transcript B-splines - Elderlab - Elder Laboratory: Human & Computer

Splines
Vida Movahedi
January 2007
Outline
•
•
•
•
•
•
Piecewise polynomials
Splines
Two approaches
B-splines
Multivariate functions
Examples
What problem splines solve?
• In the simplest situation, one is given points (ti,
yi) and is looking for a piecewise polynomial
function f that satisfies f(ti)=yi , all i, more or less.
– Exact fit  Interpolation
– Approximate fit  least squares approximation or
smoothing splines
• Or we have a complex function g that we want to
approximate with a piecewise polynomial
function f that satisfies f(ti)=g(ti) , all i, more or
less.
Linear Interpolation
Cubic Interpolation
Cubic Smoothing
Why polynomials?
• Easy evaluation, differentiation, and integration
• If done properly, the work needed for data fitting
grows only linearly with the number of data
points.
Definition of Polynomial of order n
p(x)  a 1  a 2 x    a n x
n
n -1
  a j x j-1
j 1
• Note that a polynomial of order n has degree<n
• The set of all polynomials of order n forms a
linear space
• <n: linear space of polynomials of order n
Why piecewise polynomials?
n
p(x)  a1  a 2 x    a n x n -1   a j x j-1   n
j 1
• Jackson’s Theorem: If g has r continuous derivatives
on [a..b] and n> r+1, then
ba
ba
(r )
g  C ( r ) [a..b]  (n  r  1)  dist( g ,  n )  constr 
)
 .w( g ;
n

1
2
(
n

1

r
)


w( g ; h) : max| g ( x)  g ( y) |: x, y [a..b],| x  y | h
r
• The only way to make the error small is to make (b-a)/
(n-1) small  partitioning [a..b]  using piecewise
polynomial approximation.
• It is usually much more efficient to make (b-a) small
than to increase n. This is the justification for piecewise
polynomial approximation.
Why not higher order polynomials?
• We expect the error between the function g and the
polynomial approximation pn on n sites to decrease
when n increases. If the sites are uniformly spaced, it
can be shown that this is not true and the interpolation
error increases with n for some examples.
•
It is shown that choosing the sites as the zeros of
the Chebyshev polynomial of degree n can lead to much
decreasing interpolation error with increasing n. For
some examples, the error, while decreasing, decreases
far too slowly.
What is a spline?
• An interval [a..b] is subdivided into sufficiently small
intervals [j.. j+1], with a=1<…<l+1=b,
• On each such interval, a polynomial pj of relatively low
degree can provide a good approximation to g.
• This can even be done in such a way that the polynomial
pieces blend smoothly, i.e. so that the resulting
composite function s(x) that equals pj(x) for x[j j+1] , all
j , has several continuous derivatives.
• Any such smooth piecewise polynomial function is called
a spline.
1.4
1
0.9
1.2
0.8
1
0.7
0.6
0.8
0.5
0.6
0.4
0.3
0.4
0.2
0.2
0.1
0
0
0.5
1
1.5
2
2.5
Continuous 1st derivative
3
3.5
0
0
0.5
1
1.5
2
2.5
Continuous 2nd derivative
3
3.5
Two major approaches
• The ppform: Use of piecewise polynomials
• The B-form: Use of Basis functions
ppform
• More efficient when evaluating the spline extensively
• The ppform of a polynomial spline of order k provides a
description in terms of its breaks 1,…,l+1 and the local
polynomial coefficients cji of its l pieces.
k
p j ( x)   ( x   j ) k i .c ji , j  1 : l
i 1
• For linear approximation
p j ( x)  ( x   j ).c j1  c j 2  ax  b
g(x)=sin(x)
Example
>> x=0:0.6:pi; y=sin(x);
>> figure; plot (x,y, 'o')
>> sp=spapi(2,x,y);
>> hold on; fnplt(sp);
>> pp=sp2pp(sp);
>>fnbrk(pp)
The input describes a ppform
breaks(1:l+1)
Columns 1 through 5
0 0.6000 1.2000 1.8000
Column 6
3.0000
coefficients(d*l,k)
0.9411
0
0.6123 0.5646
0.0697 0.9320
-0.4973 0.9738
-0.8906 0.6755
pieces number l
5
order k
2
dimension d of target
1
1
0.9
Breaks
(interpolation)
0.8
0.7
0.6
0.5
2.4000
0.4
0.3
Interpolating
spline of order 2
0.2
0.1
0
0
0.5
1
1.5
2
2.5
3
3.5
The ppform provides the coefficients for
each of the polynomial pieces between the
breaks. In the above example the
coefficients for each of the lines is provided.
The B-form
• Useful during construction of a spline
• The B-form
describes a spline as a weighted sum
n
f (t )   B j ,k (t ).a j
j 1
• of B-splines of the required order k, with their number
n>k-1+l
• Each Bj,k is defined on an interval [tj..tj+k] and is zero
elsewhere.
• tj are called knots and are provided based on the
smoothness required (breaks, but sometimes repeated!)
• B-splines are functions that:
n
B
j 1
j ,k
( x)  1, x  [t k ..t n1 ]
Example
1
0.9
0.8
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0
0
0.5
1
1.5
2
2.5
3
3.5
B-splines of order 2
(uniform spacing between knots)
Example (cont.)
>> fnbrk(sp)
The input describes a B-form
knots(1:n+k)
Columns 1 through 7
0
0 0.6000 1.2000 1.8000 2.4000 3.0000
Column 8
3.0000
coefficients(d,n)
0 0.5646 0.9320 0.9738 0.6755 0.1411
number n of coefficients
6
Note that in ppform, we had 6
order k
breaks, and 10 coefficients.
2
dimension d of target
1
Control Points
n
f (t )   B j ,k (t ).a j
j 1
• aj are called the control points
• If B are B-splines  aj are de Boor control points
• If B are Bezier splines  aj are Bezier control points
Bezier
Control
Points
de Boor
Control Points
for B-splines
Why B-splines?
• Composed of (n-k+2) curves of k-order joined Ck-2
continuously at knot values (t0,…,tn+k)  smoothness
and continuity
• Each point is affected by k control points
• Each control point affects k segments  local
controllability
• Inside convex hull  Built-in boundedness  does not
shoot to infinity (unlike Bezier)
• Affine invariance
• Any degree of detail desired  Generative nature
Uniform B-splines
• If the knots are equi-spaced, for a given order k,
the splines are simply shifted copies of one
another.
Linear B-splines (k=2)
Cubic B-splines (k=4)
knot multiplicity
• Knots are a sequence of breaks, sometimes
repeated!
• The rule is
Knot multiplicity + condition multiplicity = order
• For example, for a B-spline of order k=3
Simple knot  smoothness_cond=2  continuity
of function and first derivative
Double knot  smoothness_cond= 1  just
continuity
Triple knot  No condition  function can even
be discontinuous
Multivariate Functions
• If f is a function of x, and g is a function of y, then their
tensor-product p(x,y):=f(x)g(y) is a function of x and y.
• For example for a bivariate spline of order h in x and k
in y, we have:
• Knot sequences
– s= (s1,…, sm+h)
– t= (t1,…,tn+k)
• Coefficients (aij:i=1..m, j=1..n)
m
n
f ( x, y)   B( x | si ,...,si  h ) B( y | t j ,..,t j  k )aij
i 1 j 1
Chord-Length Method
• The knots must be non-decreasing!
• So, what should we do for closed curves?!
• [1,4,7,10,6,3]  not acceptable
• t: new variable , increasing as we follow the contour
• x and y will be functions of t
ti  ti 1  K
|| Qi  Qi 1 ||
m
 || Q
j 2
j
 Q j 1 ||
Representing Image Contours
50
100
150
200
20
40
60
80
100
120
140
Contour from Paco’s grouping algorithm
160
Interpolation
unit step: 50
csapi:Interpolation
Cubic Spline Interpolation
csapi: Cubic Spline
spapi:
Order k=2 Interpolation
spapi: Order k=2
Interpolation
Order k=3 Interpolation
spapi: Orderspapi:
k=3 Interpolation
spapi:
Order k=4 Interpolation
spapi: Order k=4
Interpolation
spapi:
Order k=3 Interpolation
spapi: Order k=3
Interpolation
spapi:
Order k=4 Interpolation
spapi: Order k=4
Interpolation
Interpolation
unit step: 10
Cubic Spline Interpolation
csapi: Cubic Spline csapi:
Interpolation
spapi: Order k=2 Interpolation
spapi: Order k=2 Interpolation
spapi: Order k=3spapi:
Interpolation
Order k=3 Interpolation
spapi: Order k=4spapi:
Interpolation
Order k=4 Interpolation
spapi: Order k=3 Interpolation
spapi: Order k=3 Interpolation
spapi: Order k=4 Interpolation
spapi: Order k=4 Interpolation
Smoothing
unit step: 50
csaps:
Cubic Smoothing Spline
csaps: Cubic Smoothing
Spline
spaps: Tolerance
spaps: Tolerance Smoothing
(tol=2) Smoothing (tol=2)
spap2:
Least Square Smoothing
spap2: Least Square
Smoothing
spaps: Tolerance
Smoothing
& Original Lines
spaps: Tolerance Smoothing
& Original
Lines
Least Square Smoothing
spap2: Least Squarespap2:
Smoothing
spaps: Tolerance
Smoothing
& Original Lines
spaps: Tolerance Smoothing
& Original
Lines
Smoothing
unit step: 10
csaps: Cubic Smoothing
SplineSmoothing Spline
csaps: Cubic
spap2: Least Square Smoothing
spap2: Least Square Smoothing
spap2: Least Square Smoothing
spap2: Least Square Smoothing
spaps: Tolerance Smoothing
(tol=2) Smoothing (tol=2)
spaps: Tolerance
spaps: Tolerance Smoothing & Original Lines
spaps: Tolerance Smoothing & Original Lines
spaps: Tolerance Smoothing & Original Lines
spaps: Tolerance Smoothing & Original Lines
Note the corners
At least one sample
from each line?
Unit step: 20
Cubic Spline Interpolation
csapi: Cubic Splinecsapi:
Interpolation
spapi: Order k=2 Interpolation
spapi: Order k=2 Interpolation
spapi: Order k=3 Interpolation
spapi: Order k=3 Interpolation
spapi: Order k=4 Interpolation
spapi: Order k=4 Interpolation
spapi: Order k=3 Interpolation
spapi: Order k=3 Interpolation
spapi: Order k=4 Interpolation
spapi: Order k=4 Interpolation
Questions?
• Which method is better? Interpolation or
smoothing?
• What order? Linear, cubic, …?
• How can we calculate the error? Should we
compare with the original lines? Or with image?
References
1.
Matlab Spline Toolbox
2.
Carl de Boor, A Practical Guide to Splines, Springer, 2001.
3.
‘An Interactive Introduction to Splines’, http://www.ibiblio.org/enotes/Splines/Intro.htm
4.
Cohen, F.S.; Wang, J.Y. , “Part I: Modeling Image Curves Using
Invariant 3-D Object Curve Models- A Path to 3-D Recognition and
Shape Estimation from Image Contours”, IEEE Transactions on
Pattern Analysis and Machine Intelligence, vol. 16, no. 1, Jan. 1994.
5.
Miroslaw Kuc, COSC 4080.03 Project report: “Representing Image
Contours with Splines”, Jan. 2005
6.
Branka Otasevic, COSC 4080.03 Project report, Winter 1999.