Find the center of a circle when three points are indicated

Question

I studied this link and encoded accordingly, but received the wrong answer for the example described in the link. While solving the equation, I subtracted equation 2 from equation 1 and equation 3 from equation 2, and then continued on. Please check the link for more details.

My code is:

include<stdio.h> int is_formCircle(float a1,float b1,float a2,float b2,float a3,float b3) { float p1=a2-a1; float p2=a3-a2; float p3=b2-b1; float p4=b3-b2; float alpha=(a1+a2)*(a1-a2) + (b1+b2)*(b1-b2); float beta =(a2+a3)*(a2-a3) + (b2+b3)*(b2-b3); float y1=p1*beta - p2*alpha; float y2=p2*p3 - p1*p4; if(y2==0 || y1==0) return 1; float y=y1/y2; float x1 = 2*p4*y + beta; float x2 = 2*p2; float x = x1/x2; printf("x=%fy=%f\n",x,y); return 0; } int main() { float a1,a2,a3,a4,b1,b2,b3,b4; a1=4.0; b1=1.0; a2=-3.0; b2=7.0; a3=5.0; b3=-2.0; is_formCircle(a1,b1,a2,b2,a3,b3); return 0; } 

MY other code:

 #include<stdio.h> int is_formCircle(float a1,float b1,float a2,float b2,float a3,float b3) { float mid1,mid2,mid3,mid4,m1,m2,D,Dx,Dy,x,y; mid1 = a1+(a2-a1)/2; mid2 = b1+(b2-b1)/2; mid3 = a2+(a3-a2)/2; mid4 = b2+(b3-b2)/2; m1=(b2-b1)/(a2-a1); m2=(b3-b2)/(a3-a2); m1=-1*m1; m2=-1*m2; D=m2-m1; Dx=mid2-(m1*mid1) + (mid3*m2) - mid4; Dy=(m1*(mid3*m2-mid4))-(m2*(mid1*m1-mid2)); x=Dx/D; y=Dy/D; printf("%f %f",x,y); return 0; } int main() { float a1,a2,a3,b1,b2,b3; a1=4.0; b1=1.0; a2=-3.0; b2=7.0; a3=5.0; b3=-2.0; is_formCircle(a1,b1,a2,b2,a3,b3); return 0; } 

Why is my code giving the wrong answer?

Answer 1

I have to say that if you follow the link you provided, this would help keep the variable names the same. We could better understand the algorithm by seeing x1, y1, x2, y2, x3, y3 instead of p1, p2, p3, p4, alpha and beta. In fact, I don’t see much in your algorithm that matches the link. I'm not trying to be as harsh as comments (and if you are worried about switching float to double, that would be a good example for typedef), but debugging algorithms are easiest when you don't need to convert variable names.

I would recommend just using what they give you for h and k in the link, which is done by calculating the determinants of 3x3 matrices. You can find many links for this.

I would perform two functions:

 float calculateH(float x1, float y1, float x2, float y2, float x3, float y3) { float numerator = (x2*x2+y2*y2)*y3 - (x3*x3+y3*y3)*y2 - ((x1*x1+y1*y1)*y3 - (x3*x3+y3*y3)*y1) + (x1*x1+y1*y1)*y2 - (x2*x2+y2*y2)*y1; float denominator = (x2*y3-x3*y2) - (x1*y3-x3*y1) + (x1*y2-x2*y1); denominator *= 2; return numerator / denominator; } float calculateK(float x1, float y1, float x2, float y2, float x3, float y3) { float numerator = x2*(x3*x3+y3*y3) - x3*(x2*x2+y2*y2) - (x1*(x3*x3+y3*y3) - x3*(x1*x1+y1*y1)) + x1*(x2*x2+y2*y2) - x2*(x1*x1+y1*y1); float denominator = (x2*y3-x3*y2) - (x1*y3-x3*y1) + (x1*y2-x2*y1); denominator *= 2; return numerator / denominator; } 

Then your is_formCircle will be simple:

 float is_formCircle(float x1, float y1, float x2, float y2, float x3, float y3) { float h = calculateH(x1, y1, x2, y2, x3, y3); float k = calculateK(x1, y1, x2, y2, x3, y3); printf("x=%fy=%f\n",h,k); } 

There are many ways to optimize this, and there is a chance that I have sealed any of the determinant calculations, but it should make you go.

Answer 2

The solution that was given in the link is a "blind" solution, i.e. you know the equation, the arrow solves it.

However , if you more deeply understand what is behind the scenes , you can:

• Write more readable, reliable, flexible code.
• Debugging is easy.

What happens if you subtract equation 1 from equation 2? You are really trying to find the equation of a line that describes those points that are equidistant from point 1 and point 2. Then you do the same with points 2 and 3. Finally, you will find the intersection between the lines that the center of the circle gives you.

How do you describe a straight line of points equidistant to points 1 and 2? You take a point that is in the middle of two, and move in a direction perpendicular to the direction between points 1 and 2.

If this is not entirely clear, take the paper and draw an example: place points 1,2 and 3, find two lines and find the intersection.

Now that you understand everything, change your code to two functions, on which find the line equidistant between two points, the other that calculates the intersection between two lines ...

---

After your editing, the code looked better, although it was not easy to understand. I think the error is that you decide to intersect the two lines, do not forget that you are in parametric form:

 Dx = (mid4-mid2) - m2*(mid3-mid1); lambda=Dx/D; x = mid1 + lambda*m1; y = mid2 + lambda*1.0; 

Verified graphically using Matlab.