Application Center - Maplesoft

# Classroom Tips and Techniques: Fitting Circles in Space to 3-D Data

You can switch back to the summary page by clicking here.

Classroom Tips and Techniques

Fitting Circles in Space to 3-D Data

Robert J. Lopez

Emeritus Professor of Mathematics and Maple Fellow
© Maplesoft, a division of Waterloo Maple Inc., 2010

Introduction

An algebraic solution was given in [1] for the problem of fitting a circle to three-dimensional data points not lying in the -plane.  In [2], this solution is crafted as a Matlab project whose steps we articulate in Table 1.

 1. Translate all  points by their centroid so that Π, the plane containing the fitting circle, will pass through (or nearly through) the origin. 2. Obtain the 3 × 3 factor U from the singular-value decomposition of , the 3 ×  matrix whose columns are vectors from the origin to the translated points. 3. The first two columns in  are , an orthonormal basis for Π. 4. The third column in  is , a unit vector normal to Π. 5. The columns in the 3 ×  matrix  are the components of the vectors  with respect to the basis vectors .  Since n is orthogonal to Π, the entries in the third row of  are small. 6. A least-squares fit to a circle is made in the plane Π by treating the first two entries in each column of  as , a point in an -plane. Table 1   Essential steps in an algebraic solution to the problem of fitting a circle in space to 3-D data

We will implement the algebraic solution in Maple, then provide an alternate analytic solution based on the direct minimization of the sum of squares of deviations. The analytic technique, suggested by [3], simultaneously fists a plane to the data, while finding the radius of the circle by minimizing the sum of squares of deviations of each point from a line parallel to the normal to the fitting plane.

 Initializations Load packages and define initial data by clicking on the icon below.

Algebraic Solution

Each of the ten columns of the following matrix are the Cartesian coordinates of data points taken from [2].

To each data point, we add random noise of magnitude less than 1.5. The perturbations for each coordinate are in the matrix

so that the noisy data are contained in the matrix

The coordinates of the centroid of these points appear as the components of the following vector.

The coordinates of the centroid is subtracted from each point; the columns of the resulting matrix

are vectors ,, from the origin to each of the translated points.

The singular value decomposition of this matrix provides the factor  and the vector S of singular values.  The matrix  contains orthonormal bases for the column space of  and for its orthogonal complement.

The third component of

=

the vector of singular values, is a measure of how close the translated points are to the plane Π.  As the magnitude of the noise in the data increases, this third singular value gets larger.  Hence, the algebraic approach advocated in [1] is less robust than the analytic approach sketched in the next section.

If the translated points lie on a plane through the origin, the column space of  will be two-dimensional, that is, the plane Π.  The orthogonal complement will be spanned by a vector orthogonal to the plane Π.  Consequently, the first two columns of  comprise the orthonormal basisfor the plane Π, and the third column of  namely,

is a unit vector normal to Π.

That  is an orthogonal matrix is verified by the calculation , implemented in Maple via

The columns of the product  that is,

are the components of the  along  and n.  Thus, the first two rows of  are taken as - and -coordinates of points in an -plane spanned by  and   In this plane, a least-squares fit of a circle to the  pairs represented by these rows is effected by writing

in the form

where

(as suggested in [2]).  This strategy leads to the column vectors

and to the matrix A and vector b

forming the system Ac = b.

The least-squares solution for c is then

Consequently, in the translated and tilted plane, the fitting circle has

for the radius and center, respectively.  Again, note that the center is in coordinates relative to the orthonormal basis . The coordinates of this center, in the space of the original data, are then

In the presence of noise, the third row in the matrix  does not represent a zero component along the normal to the plane spanned by  Thus, an essential condition for the validity of the algebraic solution is not met. The theory requires that the data lie in a plane for the algebraic method illustrated here to be applicable. Yet, as we will see, the algebraic solution is close to the analytical solution obtained below.

Analytic Solution

An analytic solution of the problem of fitting a circle to points in three-dimensions begins with the following definitions.

Set equal to zero, the first expression defines the plane through the point  and having

i +  j +  k

for its normal of length one. The second expression is the distance of the point  to the line parallel to the normal, and passing through . Figure 1, using concepts from multivariate calculus, is a basis for obtaining this function.

 When the point  is projected onto the line whose direction is given by the unit vector , a right triangle is formed. Elementary trigonometry gives , so . Because  has length one, we therefore can write       The desired result follows from Figure 1   Distance from the point  to the line with direction

The noisy data are stored in the matrix

=

To obtain an estimate of where the center of the sphere might lie, we form lists of -, -, and -coordinates

then compute the average of the extreme values for each coordinate.

The sum of squares of the deviations of the data points from both the plane and the line is given by

and the minimum of this sum is given by

where we constrain the normal vector to have unit length.  Thus, the center and radius of the sphere are

while the normal to the plane containing the fitting circle is

The equation of the plane containing the fitting circle is then

Figure 2 shows this plane, the data points, and the computed center.

 Click the icon below to generate Figure 2. Figure 2   The data points, the computed center, and the least-squares plane containing the data points

Obtaining a graph of the fitting circle requires considerably more work.

In a spherical coordinate system centered at , and defined by

the equation of the plane containing the analytically-found circle becomes

from which we obtain  as

To get a four-quadrant version of this result, write

which specializes to

The parametric representation of the circle in space is then

Figure 3 shows the fitted points, the fitting circle, the plane in which this circle lies, and the center of the fitting circle.

 Click the icon below to generate Figure 3. Figure 3   The least-squares circle and the noisy data to which it has been fit.

References

 [1] Carl C. Cowen, "A Project on Circles in Space."  In Resources for Teaching Linear Algebra.  Edited by David Carlson, Charles R. Johnson, David Lay, A. Duane Porter, Ann Watkins, and William Watkins.  Washington, D.C.: MAA, 1996. [2] ATLAST: Computer Exercises for Linear Algebra, Steven Leon, Eugene Herman, Richard Faulkenberry, Prentice-Hall, Inc., Upper Saddle River, NJ, 1996. [3] Craig M. Shakarji, "Least-Squares Fitting Algorithms of the NIST Algorithm Testing System," Journal of Research of the National Institute of Standards and T2chnology, Volume 103, Number 6, November-December 1998

Legal Notice: © Maplesoft, a division of Waterloo Maple Inc. 2010. Maplesoft and Maple are trademarks of Waterloo Maple Inc. This application may contain errors and Maplesoft is not liable for any damages resulting from the use of this material. This application is intended for non-commercial, non-profit use only. Contact Maplesoft for permission if you wish to use this application in for-profit activities.