Tarango

Untitled

Aug 3rd, 2016
235
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
C++ 3.65 KB | None | 0 0
  1. // 简单多边形与圆面积交
  2. #include <cstdio>
  3. #include <cstring>
  4. #include <cmath>
  5. #include <algorithm>
  6. using namespace std;
  7. const double eps = 1e-8;
  8. const double pi = acos(-1.0);
  9. const int maxn = 103;
  10. struct Point {
  11.     double x, y;
  12.     Point(){}
  13.     Point (double tx,double ty)
  14.     {
  15.         this->x = tx;
  16.         this->y = ty;
  17.     }
  18.     Point operator - (const Point &t) const {
  19.         Point res;
  20.         res.x = x - t.x;
  21.         res.y = y - t.y;
  22.         return res;
  23.     }
  24. }p[maxn];
  25. double multi(Point &o, Point &a, Point &b) {                            // 点积
  26.     return (a.x-o.x)*(b.x-o.x) + (a.y-o.y)*(b.y-o.y);
  27. }
  28. double cross(Point &o, Point &a, Point &b) {                            // 叉积
  29.     return (a.x-o.x)*(b.y-o.y) - (b.x-o.x)*(a.y-o.y);
  30. }
  31. double cp(Point &a, Point &b) {                                
  32.     return a.x*b.y - b.x*a.y;
  33. }
  34. double angle(Point &a, Point &b) {                                              // 两向量夹角
  35.     double ans = fabs((atan2(a.y, a.x) - atan2(b.y, b.x)));
  36.     return ans > pi+eps ? 2*pi-ans : ans;
  37. }
  38. void cir_line(Point ct, double r, Point l1, Point l2, Point& p1, Point& p2) {   // 直线与圆
  39.     double a1, a2, b1, b2, A, B, C, t1, t2;
  40.     a1 = l2.x - l1.x; a2 = l2.y - l1.y;
  41.     b1 = l1.x - ct.x; b2 = l1.y - ct.y;
  42.     A = a1*a1 + a2*a2;
  43.     B = (a1*b1 + a2*b2)*2;
  44.     C = b1*b1 + b2*b2 - r*r;
  45.     t1 = (-B - sqrt(B*B - 4.0*A*C))/2.0/A;
  46.     t2 = (-B + sqrt(B*B - 4.0*A*C))/2.0/A;
  47.     p1.x = l1.x + a1*t1; p1.y = l1.y + a2*t1;
  48.     p2.x = l1.x + a1*t2; p2.y = l1.y + a2*t2;
  49. }
  50. double cir_polygon(Point ct, double R, Point *p, int n) {                   // 圆与简单多边形
  51.     Point o, a, b, t1, t2;
  52.     double sum=0, res, d1, d2, d3, sign;
  53.     o.x = o.y = 0;
  54.     p[n] = p[0];
  55.     for (int i=0; i < n; i++) {
  56.         a = p[i]-ct;
  57.         b = p[i+1]-ct;
  58.         sign = cp(a,b) > 0 ? 1 : -1;
  59.         d1 = a.x*a.x + a.y*a.y;
  60.         d2 = b.x*b.x + b.y*b.y;
  61.         d3 = sqrt((a.x-b.x)*(a.x-b.x) + (a.y-b.y)*(a.y-b.y));
  62.         if (d1 < R*R+eps && d2 < R*R+eps) { //两个点都在圆内
  63.             res = fabs(cp(a, b));
  64.         }
  65.         else if (d1 < R*R-eps || d2 < R*R-eps) { //一个点在圆内
  66.             cir_line(o, R, a, b, t1, t2);
  67.             if ((a.x-t2.x)*(b.x-t2.x) < eps && (a.y-t2.y)*(b.y-t2.y) < eps) {
  68.                 t1 = t2;
  69.             }
  70.             if (d1 < d2)
  71.                 res = fabs(cp(a, t1)) + R*R*angle(b, t1);
  72.             else
  73.                 res = fabs(cp(b, t1)) + R*R*angle(a, t1);
  74.         }
  75.         else if (fabs(cp(a, b))/d3 > R-eps) { // 两个点都在园外,且线段与圆之多只有一个交点
  76.             res = R*R*angle(a, b);
  77.         }
  78.         else {  // 线段与圆有两个交点
  79.             cir_line(o, R, a, b, t1, t2);
  80.             if (multi(t1, a, b) > eps || multi(t2, a, b) > eps) {
  81.                 res = R*R*angle(a, b);
  82.             }
  83.             else {
  84.                 res = fabs(cp(t1, t2));
  85.                 if (cross(t1, t2, a) < eps)
  86.                     res += R*R*(angle(a, t1) + angle(b, t2));
  87.                 else
  88.                     res += R*R*(angle(a, t2) + angle(b, t1));
  89.             }
  90.         }          
  91.         sum += res * sign;
  92.     }
  93.     return fabs(sum)/2.0;
  94. }
  95. int main()
  96. {
  97.     int t,ca=1;
  98.     scanf("%d",&t);
  99.     while(t--)
  100.     {
  101.         Point cen;
  102.         double R;
  103.         scanf("%lf%lf%lf",&cen.x,&cen.y,&R);
  104.         Point p[4]; double a,b,c,d;
  105.         scanf("%lf%lf%lf%lf",&a,&b,&c,&d);
  106.         p[0] = Point(a,b);
  107.         p[1] = Point(c,b);
  108.         p[2] = Point(c,d);
  109.         p[3] = Point(a,d);
  110.         printf("Case %d: %lf\n",ca++,cir_polygon(cen,R,p,4));
  111.     }
  112.     return 0;
  113. }
Advertisement
Add Comment
Please, Sign In to add comment