3d curve and curvature

[William Higgs]

Hello,

I have a series of x,y,z data points and would like to create a 3d curve plot
then calculate the curvature of this line at discrete locations.

Can anyone offer some assistance with this.

Thanks

bill

[Penny Anderson]

William Higgs wrote in message ...

Hi Bill.

This response is from Carl de Boor, author of The Spline Toolbox:

Assume you have the spline toolbox for matlab and that you have listed your
3d-points in order as the columns of the matrix points (which then has
three rows, and as many columns as you have points). Then the command
c = cscvn(points);
will produce (the description of) a 3d-curve that passes through those
points
in the specified order.
Its parameter interval is given by its `basic interval', i.e.,
interv = fnbrk(c,'interval').
Next, the curvature of the resulting curve t |--> c(t) into 3-space can be
computed as
kappa(t) := norm{Dc(t) x D^2 c(t)}/norm{Dc(t)}^3
with Dc and D^2 c the first and second derivative of that curve wrto its
parameter, with a x b the cross product of the two 3-vectors a and b ,
and norm{a} the Euclidean length of the vector a..
So, to compute the curvature at 21 equi-parametric points on this curve, do
the following (which uses the commands fnder and fnval from the spline
toolbox):
t = linspace(interv(1), interv(2),21);
dc = fnder(c); dct = fnval(dc,t); ddct = fnval(fnder(dc),t);
kappa = sqrt(sum(cross(dct,ddct).^2)./(sum(dct.^2)).^3);
... and there are the curvature values.

=========================================
Penny Anderson mailto:pe...@mathworks.com
The MathWorks, Inc. http://www.mathworks.com
3 Apple Hill Drive Phone: (508) 647-7506
Natick MA 01760-2098 Fax: (508) 647-7002
=========================================

[Jonathan]

Using the above approach, I get

Error using cross (line 46)
A and B must have at least one dimension of length 3.

Error in AoCentreline (line 38)
kappa = sqrt(sum(cross(dct,ddct).^2)./(sum(dct.^2)).^3);

Have I applied this incorrectly? I just copy/pasted directly from the above response.

[Jonathan]

The problem is fixed if one uses

c = cscvn(points');

In other words, for N points, the argument to cscvn must have the coordinates arranged in N columns, rather than N rows.