A Triangle and a Circle
| Time Limit: 5000MS | Memory Limit: 131072K | |||
| Total Submissions: 3001 | Accepted: 679 | Special Judge | ||
Description
Given one triangle and one circle in the plane. Your task is to calculate the common area of these two figures.
Input
The input will contain several test cases. Each line of input describes a test case. Each test case consists of nine floating point numbers, x1, y1, x2, y2, x3, y3, x4, y4 and r, where (x1, y1), (x2, y2) and (x3, y3) are the three vertices of the triangle and (x4, y4) is the center of the circle and r is the radius. We guarantee the triangle and the circle are not degenerate.
Output
For each test case you should output one real number, which is the common area of the triangle and the circle, on a separate line. The result should be rounded to two decimal places.
Sample Input
0 20 10 0 -10 0 0 0 10
Sample Output
144.35
开始一直wa ,后来发现EPS弄得太小了（10 ^ -16 改成 10 ^ -8）
存一下模板
1 #include<iostream> 2 #include<cstdio> 3 #include<cmath> 4 #include<cstring> 5 #include<vector> 6 #include<algorithm> 7 using namespace std; 8 typedef long long LL; 9 10 const double PI = acos(-1.0); 11 const double EPS = 1e-8; 12 13 inline int sgn(double x) { 14 return (x > EPS) - (x < -EPS); 15 } 16 17 struct Point { 18 double x, y; 19 Point() {} 20 Point(double x, double y): x(x), y(y) {} 21 void read() { 22 scanf("%lf%lf", &x, &y); 23 } 24 double angle() { 25 return atan2(y, x); 26 } 27 Point operator + (const Point &rhs) const { 28 return Point(x + rhs.x, y + rhs.y); 29 } 30 Point operator - (const Point &rhs) const { 31 return Point(x - rhs.x, y - rhs.y); 32 } 33 Point operator * (double t) const { 34 return Point(x * t, y * t); 35 } 36 Point operator / (double t) const { 37 return Point(x / t, y / t); 38 } 39 double operator *(const Point &b)const 40 { 41 return x*b.x + y*b.y; 42 } 43 double length() const { 44 return sqrt(x * x + y * y); 45 } 46 Point unit() const { //单位向量 47 double l = length(); 48 return Point(x / l, y / l); 49 } 50 }; 51 double cross(const Point &a, const Point &b) { 52 return a.x * b.y - a.y * b.x; 53 } 54 double sqr(double x) { 55 return x * x; 56 } 57 double dist(const Point &p1, const Point &p2) { 58 return (p1 - p2).length(); 59 } 60 double sdist(Point a,Point b){ 61 return (a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y); 62 } 63 //向量 op 逆时针旋转 angle 64 Point rotate(const Point &p, double angle, const Point &o = Point(0, 0)) { 65 Point t = p - o; 66 double x = t.x * cos(angle) - t.y * sin(angle); 67 double y = t.y * cos(angle) + t.x * sin(angle); 68 return Point(x, y) + o; 69 } 70 Point line_inter(Point A,Point B,Point C,Point D){ //直线相交交点 71 Point ans; 72 double a1=A.y-B.y; 73 double b1=B.x-A.x; 74 double c1=A.x*B.y-B.x*A.y; 75 76 double a2=C.y-D.y; 77 double b2=D.x-C.x; 78 double c2=C.x*D.y-D.x*C.y; 79 80 ans.x=(b1*c2-b2*c1)/(a1*b2-a2*b1); 81 ans.y=(a2*c1-a1*c2)/(a1*b2-a2*b1); 82 return ans; 83 } 84 Point p_to_seg(Point p,Point a,Point b){ //点到线段的最近点 85 Point tmp=p; 86 tmp.x+=a.y-b.y; 87 tmp.y+=b.x-a.x; 88 if(cross(a-p,tmp-p)*cross(b-p,tmp-p)>0) return dist(p,a)<dist(p,b)?a:b; 89 return line_inter(p,tmp,a,b); 90 } 91 void line_circle(Point c,double r,Point a,Point b,Point &p1,Point &p2){ 92 Point tmp=c; 93 double t; 94 tmp.x+=(a.y-b.y);//求垂直于ab的直线 95 tmp.y+=(b.x-a.x); 96 tmp=line_inter(tmp,c,a,b); 97 t=sqrt(sqr(r)-sqr( dist(c,tmp)))/dist(a,b); //比例 98 p1.x=tmp.x+(b.x-a.x)*t; 99 p1.y=tmp.y+(b.y-a.y)*t; 100 p2.x=tmp.x-(b.x-a.x)*t; 101 p2.y=tmp.y-(b.y-a.y)*t; 102 } 103 struct Region { 104 double st, ed; 105 Region() {} 106 Region(double st, double ed): st(st), ed(ed) {} 107 bool operator < (const Region &rhs) const { 108 if(sgn(st - rhs.st)) return st < rhs.st; 109 return ed < rhs.ed; 110 } 111 }; 112 struct Circle { 113 Point c; 114 double r; 115 vector<Region> reg; 116 Circle() {} 117 Circle(Point c, double r): c(c), r(r) {} 118 void read() { 119 c.read(); 120 scanf("%lf", &r); 121 } 122 void add(const Region &r) { 123 reg.push_back(r); 124 } 125 bool contain(const Circle &cir) const { 126 return sgn(dist(cir.c, c) + cir.r - r) <= 0; 127 } 128 bool intersect(const Circle &cir) const { 129 return sgn(dist(cir.c, c) - cir.r - r) < 0; 130 } 131 }; 132 void intersection(const Circle &cir1, const Circle &cir2, Point &p1, Point &p2) { //两圆相交 交点 133 double l = dist(cir1.c, cir2.c); //两圆心的距离 134 double d = (sqr(l) - sqr(cir2.r) + sqr(cir1.r)) / (2 * l); //cir1圆心到交点直线的距离 135 double d2 = sqrt(sqr(cir1.r) - sqr(d)); //交点到 两圆心所在直线的距离 136 Point mid = cir1.c + (cir2.c - cir1.c).unit() * d; 137 Point v = rotate(cir2.c - cir1.c, PI / 2).unit() * d2; 138 p1 = mid + v, p2 = mid - v; 139 } 140 Point calc(const Circle &cir, double angle) { 141 Point p = Point(cir.c.x + cir.r, cir.c.y); 142 return rotate(p, angle, cir.c); 143 } 144 const int MAXN = 1010; 145 Circle cir[MAXN],cir2[MAXN]; 146 bool del[MAXN]; 147 int n; 148 double get_area(Circle* cir,int n) { //多个圆的相交面积 149 double ans = 0; 150 memset(del,0,sizeof(del)); 151 for(int i = 0; i < n; ++i) { 152 for(int j = 0; j < n; ++j) if(!del[j]) { //删除被包含的圆 153 if(i == j) continue; 154 if(cir[j].contain(cir[i])) { 155 del[i] = true; 156 break; 157 } 158 } 159 } 160 for(int i = 0; i < n; ++i) if(!del[i]) { 161 Circle &mc = cir[i]; 162 Point p1, p2; 163 bool flag = false; 164 for(int j = 0; j < n; ++j) if(!del[j]) { 165 if(i == j) continue; 166 if(!mc.intersect(cir[j])) continue; 167 flag = true; 168 intersection(mc, cir[j], p1, p2); //求出两圆的交点 169 double rs = (p2 - mc.c).angle(), rt = (p1 - mc.c).angle(); 170 if(sgn(rs) < 0) rs += 2 * PI; 171 if(sgn(rt) < 0) rt += 2 * PI; 172 if(sgn(rs - rt) > 0) mc.add(Region(rs, PI * 2)), mc.add(Region(0, rt)); //添加相交区域 173 else mc.add(Region(rs, rt)); 174 } 175 if(!flag) { 176 ans += PI * sqr(mc.r); 177 continue; 178 } 179 sort(mc.reg.begin(), mc.reg.end()); //对相交区域进行排序 180 int cnt = 1; 181 for(int j = 1; j < int(mc.reg.size()); ++j) { 182 if(sgn(mc.reg[cnt - 1].ed - mc.reg[j].st) >= 0) { //如果有区域可以合并，则合并 183 mc.reg[cnt - 1].ed = max(mc.reg[cnt - 1].ed, mc.reg[j].ed); 184 } else mc.reg[cnt++] = mc.reg[j]; 185 } 186 mc.add(Region()); 187 mc.reg[cnt] = mc.reg[0]; 188 for(int j = 0; j < cnt; ++j) { 189 p1 = calc(mc, mc.reg[j].ed); 190 p2 = calc(mc, mc.reg[j + 1].st); 191 ans += cross(p1, p2) / 2; // 192 double angle = mc.reg[j + 1].st - mc.reg[j].ed; 193 if(sgn(angle) < 0) angle += 2 * PI; 194 ans += 0.5 * sqr(mc.r) * (angle - sin(angle)); //弧所对应的的面积 195 } 196 } 197 return ans; 198 } 199 double two_cir(Circle t1,Circle t2){ //两个圆的相交面积 200 if(t1.contain(t2)||t2.contain(t1)) return PI * sqr(min(t2.r,t1.r)); 201 if(!t1.intersect(t2)) return 0; 202 double ans=0,len=dist(t1.c,t2.c); 203 double x=(sqr(t1.r)+sqr(len)-sqr(t2.r))/(2*len); 204 double angle1=acos(x/t1.r),angle2=acos((len-x)/t2.r); 205 ans=sqr(t1.r)*angle1+sqr(t2.r)*angle2-len*t1.r*sin(angle1); // 两个扇形 减去一个四边形面积 206 return ans; 207 } 208 double triangle_circle(Point a,Point b,Point o,double r){ 209 double sign=1.0; 210 double ans=0; 211 Point p1,p2; 212 a=a-o;b=b-o; 213 o=Point(0,0); 214 if(sgn(cross(a,b))==0) return 0; 215 if(sdist(a,o)>sdist(b,o)){ 216 swap(a,b); 217 sign=-1.0; 218 } 219 if(sdist(a,o)<r*r+EPS){ //a 在内， b 在外 220 if(sdist(b,o)<r*r+EPS) return cross(a,b)/2.0*sign; 221 line_circle(o,r,a,b,p1,p2); 222 if(dist(p1,b)>dist(p2,b)) swap(p1,p2); 223 double ans1=fabs(cross(a,p1)); 224 double ans2=acos(p1*b/p1.length()/b.length())*r*r; 225 ans=(ans1+ans2)/2.0; 226 if(cross(a,b)<EPS&&sign>0.0||cross(a,b)>EPS&&sign<0.0) return -ans; 227 return ans; 228 } 229 Point tmp=p_to_seg(o,a,b); 230 if(sdist(o,tmp)>r*r-EPS){ //a,b所在直线与圆没有交点 231 ans=acos(a*b/a.length()/b.length())*r*r/2.0; 232 if(cross(a,b)<EPS&&sign>0.0||cross(a,b)>EPS&&sign<0.0) return -ans; 233 return ans; 234 } 235 line_circle(o,r,a,b,p1,p2); 236 double cm=r/(dist(a,o)-r); 237 Point m=Point((o.x+cm*a.x)/(1+cm),(o.y+cm*a.y)/(1+cm) ); 238 double cn=r/(dist(o,b)-r); 239 Point n=Point( (o.x+cn*b.x)/(1+cn) , (o.y+cn*b.y)/(1+cn) ); 240 double ans1 = acos( m*n/m.length()/n.length() )*r*r; 241 double ans2 = acos( p1*p2/p1.length()/p2.length() )*r*r-fabs(cross(p1,p2)); 242 ans=(ans1-ans2)/2.0; 243 if(cross(a,b)<EPS&&sign>0.0||cross(a,b)>EPS&&sign<0.0) return -ans; 244 return ans; 245 } 246 int main() { 247 #ifndef ONLINE_JUDGE 248 freopen("input.txt","r",stdin); 249 #endif // ONLINE_JUDGE 250 Point a,b,c; 251 int cas=0; 252 Circle t1; 253 double d1,d2,d3,d4,d5,d6,d7,d8,d9; 254 while(scanf("%lf%lf%lf%lf%lf%lf%lf%lf%lf",&d1,&d2,&d3,&d4,&d5,&d6,&d7,&d8,&d9)!=EOF){ 255 a.x=d1;a.y=d2; 256 b.x=d3;b.y=d4; 257 c.x=d5;c.y=d6; 258 t1.c.x=d7;t1.c.y=d8;t1.r=d9; 259 double ans=0; 260 ans+=triangle_circle(a,b,t1.c,t1.r); 261 ans+=triangle_circle(b,c,t1.c,t1.r); 262 ans+=triangle_circle(c,a,t1.c,t1.r); 263 printf("%.2f ",fabs(ans)+EPS); 264 } 265 }