POINTS2CIRCLE

determine the circle through 3 points (v2.0)

Description

Given three points, this function determines the center and the radius of the circle passing through these three points.

It does so analytically using matrix minors.

Although written in R14, it is likely to work in most (older) versions as well.

Current version 2.0, march 2008

This submission was inspired by a cssm post by Peter Bone on 3-6-2008, and his solution (File ID 19083)

• MATLAB release: MATLAB 7 (R14)

[Peter Bone]

Nothing wrong with a bit of friendly competition!

[John D'Errico]

This submission and its cousin, calc_circle (FEX 19083) were submitted together, in the spirit of friendly competition. Lets take an extra careful look at these two, in that same spirit. How do they compare?

The most obvious question is, does this function work? I'll create a set of three points on a known circle, here with center [1 2], and radius 1.3.

center = [1,2];
theta = [.1 .2 .98]; % radians
radius = 1.3;
circ = @(center,rad,t) center + rad*[cos(t),sin(t)];
p1 = circ(center,radius,theta(1));
p2 = circ(center,radius,theta(2));
p3 = circ(center,radius,theta(3));

xyr = points2circle(p1,p2,p3)
xyr =
            1 2 1.3

So, yes, it does indeed seem to work properly, returning the correct center and radius for my test points.

How fast is this code? I've used Steve's code, timeit (find it on the fex!) for the comparison. It gives a good estimate for the time required.

timeit(@() points2circle(p1,p2,p3))
ans =
   0.00047171

timeit(@() calc_circle(p1,p2,p3))
ans =
   0.00017641

Note that points2circle is quite fast. There really is not a huge amount of work to do, so the time should not be a serious factor, UNLESS you have this function call deep inside nested loops and you need maximum speed. In that event, calc_circle (Peter's code) will give you almost 3x the throughput.

Pure speed is not all that matters however. While it is important, we also need to look at other factors. Does this code detect collinear data?

p1 = [0 0]; p2 = [1 1]; p3 = [3 3];

xyr = points2circle(p1,p2,p3)
Warning: No solution! Points may be on a straight line (colinear).
> In points2circle at 70
xyr =
   NaN NaN NaN

As I might hope, it finds the problem and returns a reasonable result, with a warning message. Kudos to points2circle here, since the warning message is done properly. See that we can turn off this warning message as below:

warning('off','points2circle:ColinearPoints')
xyr = points2circle(p1,p2,p3)
xyr =
   NaN NaN NaN

Next, I looked at the interface. This code, points2circle, has a more flexible interface, allowing the user to provide their points in either of two logical forms.

Points2circle has error checks on its arguments. I found an H1 line. The help is well done, in a style that Matlab users will find completely standard. They will also find an excellent example of use for the code.

An interesting point to note is whether this code survives a common problem. What happens if your data is not composed of doubles?

p1 = rand(1,2);p2 = rand(1,2);p3 = rand(1,2);

This one works:
xyr = points2circle(p1,p2,p3)
xyr =
      0.31796 0.3303 0.3232

Single data is fine:
xyr = points2circle(single(p1),single(p2),single(p3))
xyr =
      0.31796 0.3303 0.3232

But no test is done for non-floating point data.

xyr = points2circle(uint8(p1),uint8(p2),uint8(p3))
??? Undefined function or method 'det' for input arguments of type 'uint8'.

Error in ==> points2circle>local_minordet at 82
md = det(M(~([1:4]==i),~([1:4]==j))) ;

Error in ==> points2circle at 66
M11 = local_minordet(M,1,1) ;

An often overlooked issue is, can points2circle survive complex data?
xyr = points2circle(p1+i,p2+i,p3+i)
xyr =
      0.31796 + 1i 0.3303 + 1i 0.3232 + 6.8702e-16i

Note that it got the center right, but the radius was a complex number. (Calc_circle had its own problems on these test, returning an incorrect result for the unit8 data with no error generated. But it did get the radius correct, since it used norm to compute the radius of the circle.)

Overall, this code gets a rating of 4 for me. The uint8 problem and complex radius pushed me below a 5 rating. If your data is required to be single or double precision, then you should test for that and return an error message! (The author of points2circle knows that of course.) On the other hand, the excellent help and flexibility of the code, plus niceties like a proper use of warning were all points in the favor of points2circle.

With a little cleanup, I'd rate this a 5.

[Duane Hanselman]

Returns complex valued radius for some real valued data, e.g.,
theta=[0 pi/10 pi/3]';
r=.21;
xc=1e3;
yc=-2e6;
x=r*sin(theta)+xc;
y=r*cos(theta)+yc;
xyr = points2circle([x y])

[Urs (us) Schwarz]

jos, unlike other (very verbose) referees, i won't rate your submission for now - because the snippet may - indeed(!) - yield erroneous (rating=1)/unexpected results (that were not caught by others, btw); and i do know that YOU can do MUCH better...
us