function [c,r] = circle(x1,y1,x2,y2,x3,y3)

N1 = det ([x2*x2 + y2*y2 - (x1*x1 + y1*y1), y2-y1;
           x3*x3 + y3*y3 - (x1*x1 + y1*y1), y3-y1]);

N2 = det ([x2-x1, x2*x2 + y2 * y2 - (x1*x1 + y1*y1);
           x3-x1, x3*x3 + y3 * y3 - (x1*x1 + y1*y1)]);

D = 2 * det ([x2-x1, y2-y1;
              x3-x1, y3-y1]);

if D == 0,
  c = [0;0];
  r = 0
else
  c = [N1/D;N2/D];
  r = norm([x1;y1] - c);
end