Xref: utzoo comp.graphics:6098 sci.math:7002 sci.math.num-analysis:30
Path: utzoo!attcan!uunet!mcvax!cernvax!cui!ugun2b!ugobs!pfenniger
From: pfenniger@obs.unige.ch
Newsgroups: comp.graphics,sci.math,sci.math.num-analysis
Subject: Re: Need to fit a circle to some points
Message-ID: <221@obs.unige.ch>
Date: 8 Jun 89 12:26:19 GMT
References: <573@lehi3b15.csee.Lehigh.EDU>
Organization: University of Geneva, Switzerland
Lines: 46

In article <573@lehi3b15.csee.Lehigh.EDU>, flash@lehi3b15.csee.Lehigh.EDU (Stephen Corbesero) writes:
> I need an algorithm to fit a set of data points, obtained
> experimentally, to the best-fit circle.  I have searched through many
> book son least squares anaylsis, but can only find simple cases.  I
> started to do the math, but it became rather involved very quickly.
> 

The problem can be stated as follow :
You have a set of points in cartesian coordinates {Xi, Yi} i = 1, ..., n (n>2). 
You want to fit to these points in a least square sense the quadratic form of a 
circle of unknown radius R and center {X0, Y0} : 

    f(X,Y) = (X - X0)^2 + (Y - Y0)^2 - R^2 

           =  a X + b Y + c + X^2 + Y^2  ,

where a = -2 X0, b = -2 Y0, c = X0^2 + Y0^2 - R^2. 
In matrix form the problem is to find Z which minimizes the euclidian norm :

    || A Z - B ||

where 

      [ a ]
  Z = [ b ]
      [ c ]

and 

      [ X1 Y1 1 ]          [ X1^2 + Y1^2 ] 
  A = [ X2 Y2 1 ] ,   B = -[ X2^2 + Y2^2 ]  .
        .......              ...........
      [ Xn Yn 1 ]          [ Xn^2 + Yn^2 ]

The problem is therefore suited to *linear* least squares, since the unknowns
a, b, c in Z enter in a linear way.
Then use any modern linear least squares algorithm (e.g. HFTI in Lawson &
Hanson, 1974, "Solving Least Squares Problems", Prentice-Hall, NJ ; or Linpack
etc.) which just finds this minimum euclidian norm by orthogonal decomposition.
This approach is superior to the traditional inversion of a square normal
matrix. 

Note that fitting ellipses or more complicated forms would use essentially the
same method, since the additional unknown parameters are also linear.

          Daniel Pfenniger