COMPGROUPS.NET | Search | Post Question | Groups | Stream | About | Register

### Fitting a curve into a 3D point cloud

• Follow

```I have data that is made up of a 3d point cloud (x,y,z coordinates). It essentially looks like a curved tube. My problem is that I can't find a way to calculate the center-line of this "tube of points".  I can only find information on fitting surfaces to 3d point clouds using least square methods.

Any ideas?

Help would be greatly appreciated!
```
 0

```"Stan " <zufallsacc@yahoo.com> wrote in message <i7vk5t\$pir\$1@fred.mathworks.com>...
> I have data that is made up of a 3d point cloud (x,y,z coordinates). It essentially looks like a curved tube. My problem is that I can't find a way to calculate the center-line of this "tube of points".  I can only find information on fitting surfaces to 3d point clouds using least square methods.
>

%Simulated data on the parametric line [1;1;1]+t*[0 0 1]

N=10;
XYZ=bsxfun(@plus, [ 1 1 1], (-N:N).'*[0 0 1]);

%Fit to a line xyz0 + t*direction

xyz0=mean(XYZ);
A=bsxfun(@minus,XYZ,xyz0); %center the data

[U,S,V]=svd(A);

xyz0,
direction=cross(V(:,end),V(:,end-1) ),
```
 0

```"Matt J " <mattjacREMOVE@THISieee.spam> wrote in message <i7vnqi\$rus\$1@fred.mathworks.com>...

> %Fit to a line xyz0 + t*direction
>
>
> xyz0=mean(XYZ);
> A=bsxfun(@minus,XYZ,xyz0); %center the data
>
>
> [U,S,V]=svd(A);
>
> xyz0,
> direction=cross(V(:,end),V(:,end-1) ),
=================

Forget it. I hadn't noticed in your post that the tube was curved. Do you have a model for the shape of the tube, i.e., an equation it is supposed to obey?
```
 0

```"Matt J " <mattjacREMOVE@THISieee.spam> wrote in message <i7vr8q\$f0m\$1@fred.mathworks.com>...
> "Matt J " <mattjacREMOVE@THISieee.spam> wrote in message <i7vnqi\$rus\$1@fred.mathworks.com>...
>
> > %Fit to a line xyz0 + t*direction
> >
> >
> > xyz0=mean(XYZ);
> > A=bsxfun(@minus,XYZ,xyz0); %center the data
> >
> >
> > [U,S,V]=svd(A);
> >
> > xyz0,
> > direction=cross(V(:,end),V(:,end-1) ),
> =================
>
>
> Forget it. I hadn't noticed in your post that the tube was curved. Do you have a model for the shape of the tube, i.e., an equation it is supposed to obey?

No, I don't have a model right now. I thought about just projecting the points onto the xy, xz and yz planes and then do 3 2D curve fits with polynomials of a sufficiently high order and thus get an approximation for the center line shape, but I guess this is not really elegant..and won't work in more complicated cases involving loops..
```
 0

```"Stan " <zufallsacc@yahoo.com> wrote in message <i819b0\$rg6\$1@fred.mathworks.com>...
> "Matt J " <mattjacREMOVE@THISieee.spam> wrote in message <i7vr8q\$f0m\$1@fred.mathworks.com>...
> > "Matt J " <mattjacREMOVE@THISieee.spam> wrote in message <i7vnqi\$rus\$1@fred.mathworks.com>...
> >
> > > %Fit to a line xyz0 + t*direction
> > >
> > >
> > > xyz0=mean(XYZ);
> > > A=bsxfun(@minus,XYZ,xyz0); %center the data
> > >
> > >
> > > [U,S,V]=svd(A);
> > >
> > > xyz0,
> > > direction=cross(V(:,end),V(:,end-1) ),
> > =================
> >
> >
> > Forget it. I hadn't noticed in your post that the tube was curved. Do you have a model for the shape of the tube, i.e., an equation it is supposed to obey?
>
>
> No, I don't have a model right now. I thought about just projecting the points onto the xy, xz and yz planes and then do 3 2D curve fits with polynomials of a sufficiently high order and thus get an approximation for the center line shape, but I guess this is not really elegant..and won't work in more complicated cases involving loops..

No. That will work terribly. And high order polynomials
are not a good choice for the fit.

Are you willing to assume the tube has constant radius,
or can that vary? Is the tube perfectly circular in cross-
section at any point, if you look orthogonally to the
tube axis?

These are things you need to specify, if you will choose
a method for deriving a "model" for this.

John
```
 0

```"John D'Errico" <woodchips@rochester.rr.com> wrote in message <i81kik\$mkp\$1@fred.mathworks.com>...
> "Stan " <zufallsacc@yahoo.com> wrote in message <i819b0\$rg6\$1@fred.mathworks.com>...
> > "Matt J " <mattjacREMOVE@THISieee.spam> wrote in message <i7vr8q\$f0m\$1@fred.mathworks.com>...
> > > "Matt J " <mattjacREMOVE@THISieee.spam> wrote in message <i7vnqi\$rus\$1@fred.mathworks.com>...
> > >
> > > > %Fit to a line xyz0 + t*direction
> > > >
> > > >
> > > > xyz0=mean(XYZ);
> > > > A=bsxfun(@minus,XYZ,xyz0); %center the data
> > > >
> > > >
> > > > [U,S,V]=svd(A);
> > > >
> > > > xyz0,
> > > > direction=cross(V(:,end),V(:,end-1) ),
> > > =================
> > >
> > >
> > > Forget it. I hadn't noticed in your post that the tube was curved. Do you have a model for the shape of the tube, i.e., an equation it is supposed to obey?
> >
> >
> > No, I don't have a model right now. I thought about just projecting the points onto the xy, xz and yz planes and then do 3 2D curve fits with polynomials of a sufficiently high order and thus get an approximation for the center line shape, but I guess this is not really elegant..and won't work in more complicated cases involving loops..
>
> No. That will work terribly. And high order polynomials
> are not a good choice for the fit.
>
> Are you willing to assume the tube has constant radius,
> or can that vary? Is the tube perfectly circular in cross-
> section at any point, if you look orthogonally to the
> tube axis?
>
> These are things you need to specify, if you will choose
> a method for deriving a "model" for this.
>
> John

Thanks for the input, I was already suspecting that it might be a terrible idea ;).

Yes the tube has a given radius and is circular in cross-section.

Another idea I have is partitioning the tube and using RCA on each section, but I wonder if there is a better solution

Stan
```
 0

```Hi all,

I have basically the same problem, perhaps slightly simpler in that:

1. My pts are of a curve with noise in 3D (rather than surface pts of a cylinder, but it's a similar idea).
2. I have known anchor pts at either end (ie, the ends of my curve are constrained)
3. My curve is of reasonably known shape (human rib - basically a skewed "C")
4. I only want a "reasonable" approximation.

For my reasonable approximation, I'm thinking of the following rules:
1. The final curve should be of, say, 10 equally spaced control points
2. The maximum angle between adjacent curve control points is, say, 30 degrees

I was thinking of approaching the problem by starting with 10 equally spaced control points from my start-to-end anchor pts. Then I would:

1. Assign each of my sample pts to its nearest control pt
2. Shift each control pt to the mean of the sample pts in its domain. Any control pts without sample pts in its domain get linearly distributed between its neighbours which do.
3. Re-sample the control points to be equally spaced
4. Repeat from step 1 until a certain number of iterations are reached

Any thoughts on this? Did the OP manage to get towards a solution for their problem?

Thanks,
Sven.
```
 0

```"Sven" <sven.holcombe@gmail.deleteme.com> wrote in message <i92rq8\$80r\$1@fred.mathworks.com>...
> Hi all,
>
> I have basically the same problem, perhaps slightly simpler in that:
>
> 1. My pts are of a curve with noise in 3D (rather than surface pts of a cylinder, but it's a similar idea).
> 2. I have known anchor pts at either end (ie, the ends of my curve are constrained)
> 3. My curve is of reasonably known shape (human rib - basically a skewed "C")
> 4. I only want a "reasonable" approximation.
>
> For my reasonable approximation, I'm thinking of the following rules:
> 1. The final curve should be of, say, 10 equally spaced control points
> 2. The maximum angle between adjacent curve control points is, say, 30 degrees
>
> I was thinking of approaching the problem by starting with 10 equally spaced control points from my start-to-end anchor pts. Then I would:
>
> 1. Assign each of my sample pts to its nearest control pt
> 2. Shift each control pt to the mean of the sample pts in its domain. Any control pts without sample pts in its domain get linearly distributed between its neighbours which do.
> 3. Re-sample the control points to be equally spaced
> 4. Repeat from step 1 until a certain number of iterations are reached
>
> Any thoughts on this? Did the OP manage to get towards a solution for their problem?
>
> Thanks,
> Sven.

Why not simply fitting the curve with three parametrized splines, each applied respectively on x, y, and z coordinates?

Bruno
```
 0

```"Bruno Luong" <b.luong@fogale.findmycountry> wrote in message <i92utj\$q9e\$1@fred.mathworks.com>...
> "Sven" <sven.holcombe@gmail.deleteme.com> wrote in message <i92rq8\$80r\$1@fred.mathworks.com>...
> > Hi all,
> >
> > I have basically the same problem, perhaps slightly simpler in that:
> >
> > 1. My pts are of a curve with noise in 3D (rather than surface pts of a cylinder, but it's a similar idea).
> > 2. I have known anchor pts at either end (ie, the ends of my curve are constrained)
> > 3. My curve is of reasonably known shape (human rib - basically a skewed "C")
> > 4. I only want a "reasonable" approximation.
> >
> > For my reasonable approximation, I'm thinking of the following rules:
> > 1. The final curve should be of, say, 10 equally spaced control points
> > 2. The maximum angle between adjacent curve control points is, say, 30 degrees
> >
> > I was thinking of approaching the problem by starting with 10 equally spaced control points from my start-to-end anchor pts. Then I would:
> >
> > 1. Assign each of my sample pts to its nearest control pt
> > 2. Shift each control pt to the mean of the sample pts in its domain. Any control pts without sample pts in its domain get linearly distributed between its neighbours which do.
> > 3. Re-sample the control points to be equally spaced
> > 4. Repeat from step 1 until a certain number of iterations are reached
> >
> > Any thoughts on this? Did the OP manage to get towards a solution for their problem?
> >
> > Thanks,
> > Sven.
>
> Why not simply fitting the curve with three parametrized splines, each applied respectively on x, y, and z coordinates?
>
> Bruno

I don't think that I can do this immediately on my input points as aren't necessarily ordered and since the curve can be "C" shaped, they can double back in any given direction, so ordering by their coordinates can give unpredictable results. Below is some sample input - if I'm wrong and it's actually quite simple to fit 3 splines, can you show me how? Perhaps a good way to order my kind of input points would be by their distance from the start/end locations. I'm working on a function of the form
function cPts = fitCurveTo3DPts(XYZ, stPt, endPt, n)

I'll post here if/when it's done, but I'd appreciate any input you have.

Thanks,
Sven.

stPt = [53   136  -222];
endPt = [150   265  -112];
n = 10;
XYZ = [122.6333  277.2984 -113.9278
125.6526  275.7797 -113.9278
127.1622  275.7797 -113.9278
127.1622  275.7797 -111.9149
128.6719  274.2609 -111.9149
130.1815  274.2609 -111.9149
134.7105  271.2234 -111.9149
137.7298  271.2234 -111.9149
139.2394  271.2234 -111.9149
149.8070  265.1484 -111.9149
143.7684  268.1859 -109.9021
44.1314  149.7234 -214.5722
41.1121  154.2797 -210.5464
39.6025  155.7984 -208.5335
41.1121  154.2797 -208.5335
38.0928  160.3547 -206.5206
36.5832  172.5047 -198.4691
39.6025  184.6547 -186.3918
41.1121  189.2109 -182.3660
41.1121  190.7297 -180.3531
42.6218  193.7672 -178.3402
42.6218  196.8047 -176.3273
42.6218  199.8422 -174.3144
44.1314  202.8797 -172.3015
44.1314  204.3984 -170.2887
45.6411  205.9172 -170.2887
45.6411  208.9547 -168.2758
45.6411  210.4734 -166.2629
92.4403  271.2234 -128.0180
95.4596  272.7422 -126.0052
96.9692  272.7422 -123.9923
98.4789  274.2609 -123.9923
99.9885  275.7797 -123.9923
101.4982  275.7797 -121.9794
104.5175  275.7797 -121.9794
107.5368  277.2984 -121.9794
107.5368  275.7797 -119.9665
109.0464  277.2984 -119.9665
110.5561  277.2984 -117.9536
112.0657  278.8172 -117.9536
113.5754  280.3359 -115.9407
69.7955  254.5172 -138.0825
72.8148  257.5547 -138.0825
71.3051  257.5547 -136.0696
72.8148  259.0734 -136.0696
75.8341  260.5922 -136.0696
75.8341  260.5922 -134.0567
78.8534  262.1109 -134.0567
80.3631  263.6297 -134.0567
78.8534  263.6297 -132.0438
80.3631  265.1484 -132.0438
83.3824  266.6672 -132.0438
83.3824  266.6672 -130.0309
84.8920  268.1859 -130.0309
87.9113  269.7047 -130.0309
86.4017  268.1859 -128.0180
87.9113  269.7047 -128.0180
90.9306  271.2234 -128.0180
92.4403  271.2234 -126.0052
54.6990  231.7359 -152.1727
56.2086  233.2547 -150.1598
57.7183  236.2922 -150.1598
57.7183  237.8109 -148.1469
59.2279  239.3297 -146.1340
60.7376  240.8484 -146.1340
60.7376  242.3672 -144.1211
62.2472  243.8859 -144.1211
62.2472  245.4047 -142.1082
65.2665  248.4422 -142.1082
66.7762  249.9609 -140.0954
68.2858  252.9984 -140.0954
68.2858  252.9984 -138.0825
47.1507  213.5109 -164.2500
48.6604  216.5484 -162.2371
48.6604  218.0672 -160.2242
50.1700  219.5859 -160.2242
50.1700  221.1047 -158.2113
51.6797  222.6234 -158.2113
51.6797  224.1422 -156.1985
53.1893  225.6609 -156.1985
53.1893  228.6984 -154.1856
38.0928  172.5047 -198.4691
39.6025  172.5047 -196.4562
39.6025  174.0234 -194.4433
50.1700  139.0922 -224.6366
53.1893  136.0547 -224.6366
50.1700  139.0922 -222.6237
53.1893  136.0547 -222.6237
50.1700  139.0922 -220.6108
51.6797  139.0922 -218.5979];
```
 0

```"Sven" <sven.holcombe@gmail.deleteme.com> wrote in message <i94i0j\$beo\$1@fred.mathworks.com>...
> "Bruno Luong" <b.luong@fogale.findmycountry> wrote in message <i92utj\$q9e\$1@fred.mathworks.com>...
> > "Sven" <sven.holcombe@gmail.deleteme.com> wrote in message <i92rq8\$80r\$1@fred.mathworks.com>...
> > > Hi all,
> > >
> > > I have basically the same problem, perhaps slightly simpler in that:
> > >
> > > 1. My pts are of a curve with noise in 3D (rather than surface pts of a cylinder, but it's a similar idea).
> > > 2. I have known anchor pts at either end (ie, the ends of my curve are constrained)
> > > 3. My curve is of reasonably known shape (human rib - basically a skewed "C")
> > > 4. I only want a "reasonable" approximation.
> > >
> > > For my reasonable approximation, I'm thinking of the following rules:
> > > 1. The final curve should be of, say, 10 equally spaced control points
> > > 2. The maximum angle between adjacent curve control points is, say, 30 degrees
> > >
> > > I was thinking of approaching the problem by starting with 10 equally spaced control points from my start-to-end anchor pts. Then I would:
> > >
> > > 1. Assign each of my sample pts to its nearest control pt
> > > 2. Shift each control pt to the mean of the sample pts in its domain. Any control pts without sample pts in its domain get linearly distributed between its neighbours which do.
> > > 3. Re-sample the control points to be equally spaced
> > > 4. Repeat from step 1 until a certain number of iterations are reached
> > >
> > > Any thoughts on this? Did the OP manage to get towards a solution for their problem?
> > >
> > > Thanks,
> > > Sven.
> >
> > Why not simply fitting the curve with three parametrized splines, each applied respectively on x, y, and z coordinates?
> >
> > Bruno
>
> I don't think that I can do this immediately on my input points as aren't necessarily ordered and since the curve can be "C" shaped, they can double back in any given direction, so ordering by their coordinates can give unpredictable results. Below is some sample input - if I'm wrong and it's actually quite simple to fit 3 splines, can you show me how? Perhaps a good way to order my kind of input points would be by their distance from the start/end locations. I'm working on a function of the form
> function cPts = fitCurveTo3DPts(XYZ, stPt, endPt, n)
>
> I'll post here if/when it's done, but I'd appreciate any input you have.

Ok, so here's what I've come up with. I think it's reasonably concise and seems to be working well for my type of input. Currently I'm hard-coding 5 iterations based purely on my own testing. I imagine it could be generalised to exit based on some degree of stability in the resulting control points. I can certainly imagine the type of twisting or complicated input for which my function below would handle quite poorly. Any suggestions would be great.

Cheers,
Sven.

function cPts = fitCurveTo3DPts(XYZ, stPt, endPt, n)
MAX_ITERATIONS = 5;
completedIterations = 0;

% Initialise control pts linearly between start/end anchors
cPts = interp1(0:1, [stPt; endPt], linspace(0,1,n));

while completedIterations < MAX_ITERATIONS
% Get the nearest-neighbour cntrl pt for each of the sample points
sqDists = cellfun(@(x)sum(bsxfun(@minus, XYZ, x).^2,2), num2cell(cPts,2),'UniformOutput', false);
[~, nnIdxs] = min(cat(2,sqDists{:}),[],2);
% Keep the anchors, update inner cPts to the mean of their linked input pts
for i = 2:n-1
cPts(i,:) = mean(XYZ(nnIdxs==i,:),1);
end
% Handle any cPts that didn't have linked pts so their mean became NaN
goodIdxs = find(~isnan(cPts(:,1)));
% Re-spread the control points out linearly
cPtCumSumDists = cumsum([0; sqrt(sum(diff(cPts,1).^2,2))]);
cPts = interp1(cPtCumSumDists, cPts, linspace(0, cPtCumSumDists(end), n));
completedIterations = completedIterations + 1;
end
```
 0
Reply sven.holcombe8218 (103) 10/13/2010 3:56:04 PM

9 Replies
538 Views

Similiar Articles:

7/22/2012 4:36:03 PM