1 #include <cmath>
 2 #include <cstring>
 3 #include <cstdio>
 4 #include <cfloat>
 5
 6 // 内切圆半径
 7 // 内切圆半径r=2S/(a+b+c)，其中S是三角形面积，a、b、c是三角形三边。
 8 // 另外S=根号下p(p-a)(p-b)(p-c)，其中p=(a+b+c)/2
 9 float TrinagleInCircle(float xA, float yA, float zA,
10                        float xB, float yB, float zB,
11                        float xC, float yC, float zC)
12 {
13     float ab = (xA - xB)*(xA - xB) + (yA - yB)*(yA - yB) + (zA - zB)*(zA - zB);
14     float bc = (xC - xB)*(xC - xB) + (yC - yB)*(yC - yB) + (zC - zB)*(zC - zB);
15     float ca = (xA - xC)*(xA - xC) + (yA - yC)*(yA - yC) + (zA - zC)*(zA - zC);
16
17     ab = sqrtf(ab);
18     bc = sqrtf(bc);
19     ca = sqrtf(ca);
20
21     float p = (ab+bc+ca)/2;
22
23     float s = p*(p-ab)*(p-bc)*(p-ca);
24     if (s < FLT_EPSILON)
25     {
26         return 0.0f;
27     }
28     s = sqrtf(s);
29
30     if (ab+bc+ca < FLT_EPSILON)
31     {
32         return 0.0f;
33     }
34
35     float r = 2*s/(ab+bc+ca);
36
37     return r;
38 }
39
40 // 外接圆半径
41 // 已知三角形三边长a,b,c ，及其外接圆的半径R
42 // s=a*b*c/(4*R)      (多半用于求外接圆半径 R=a*b*c/(4*s)
43 float TrinagleOutCircle(float xA, float yA, float zA,
44                         float xB, float yB, float zB,
45                         float xC, float yC, float zC)
46 {
47     float ab = (xA - xB)*(xA - xB) + (yA - yB)*(yA - yB) + (zA - zB)*(zA - zB);
48     float bc = (xC - xB)*(xC - xB) + (yC - yB)*(yC - yB) + (zC - zB)*(zC - zB);
49     float ca = (xA - xC)*(xA - xC) + (yA - yC)*(yA - yC) + (zA - zC)*(zA - zC);
50
51     ab = sqrtf(ab);
52     bc = sqrtf(bc);
53     ca = sqrtf(ca);
54
55     float p = (ab+bc+ca)/2;
56
57     float s = p*(p-ab)*(p-bc)*(p-ca);
58     if (s < FLT_EPSILON)
59     {
60         return 0.0f;
61     }
62     s = sqrtf(s);
63
64     float r = ab*bc*ca/(4*s);
65
66     return r;
67 }
68
69 int ratio (int   x1,int   y1,int   z1,int   x2,int   y2,int   z2,int   x3,int   y3,int   z3)
70 {
71     float r = TrinagleInCircle(float(x1), float(y1), float(z1), float(x2), float(y2), float(z2), float(x3), float(y3), float(z3));
72     float R = TrinagleOutCircle(float(x1), float(y1), float(z1), float(x2), float(y2), float(z2), float(x3), float(y3), float(z3));
73     if (R < FLT_EPSILON || r < FLT_EPSILON)
74     {
75         return 0;
76     }
77     return (int)(1000*R*R/r/r);
78 }