Back in July of 2005, one of the early Tips & Techniques articles (since updated) in the Maple Reporter was a comparison of two different approaches to fitting a circle to 3D data points. The impetus for the comparison was Carl Cowen's article on the subject. His approach was algebraic - he used the singular value decomposition to obtain a basis for the best-fit plane, and its normal. Now I knew Carl because he was on the faculty at my alma mater, Purdue University. Purdue is about 100 miles from Rose-Hulman in Indiana, so we ran across each other at some of the local math society meetings during the years I was at RHIT.

Because I knew Carl, and because using the singular value decomposition to perform what I thought of as an analytic task was so intriguing, I just had to explore the work. It took me a few years to get to it, and finally the results appeared in the 2005 reporter article. In the algebraic approach, the data points were translated so their center of mass moved to the origin in the expectation that the plane containing the best-fit circle would pass through the origin.

 When the best-fit circle in this plane was found, its center was given in terms of the basis of this plane. No graph was drawn, the coordinates of the center were not translated back into the original frame, and the work was all done in MATLAB®.  So, or course, the first thing I did was to implement this algebraic solution in Maple.  But initially, I didn't express the coordinates of the center in the original frame, I didn't draw a graph, but I did add some random noise just to see if the algebraic approach was robust.

 Then, my analytic approach was least-squares, with the sum of squares of deviations being the sum of two terms. One term fitted the data to a sphere; and the other, to a plane through the center of that sphere. The result can be seen in Figure 1, taken from the 2005 Reporter article.

It puzzled me that the "best-fit plane" should lie below the cloud of points. I would have expected to find a more equal distribution of points on either side of the plane.  So, finally I decided to investigate, and began by searching the internet for related material. I quickly found [3] in the revised article that's now in the Application Center, and in it I discovered that least squares should be used to fit the points to a plane, then used to minimize the sum of perpendicular distances of each point to a line parallel to the plane's normal. This provides the radius, and the location of the center of the best-fit circle.

The revised fitting functions are

where

  i +  j +  k 

is a unit normal on the plane.  The function  gives the distance of the point  to the plane, and the function  gives the distance of that point to a line parallel to the normal and through , the center of the circle.  Table 1, containing Figure 2, sketches the derivation of the function .

Figure 3 contains the result of the revised analytic approach to finding the best-fit circle.

Rotating the graph in Figure 3 shows that the circle lies in the plane, and the plane is more evenly centered between the noisy data points.

MATLAB is a registered trademark of The MathWorks, Inc.

Download FittingCircles.mw