计算几何基础与一般模板:POJ 3675&&ZOJ 1450&&POJ 1375
Telescope POJ 3675
一道模板题:
这个模板太难打了,记住就行了吧。。。
就是给你一个圆,求一个多边形在圆里的面积:
分五种情况!!!
因为任何多边形可以分成若干个三角形。。。
那么就判断两点与圆心形成的三角形来求面积。。。
参见:http://hi.baidu.com/billdu/item/703ad4e15d819db52f140b0b(五种情况,带图的哦!亲)
接下来就是模板。。。转载于:http://gzhu-101majia.iteye.com/blog/1128345
1 //多边形与圆点相交面积 poj3675 2 #include <iostream> 3 #include <stdio.h> 4 #include <cmath> 5 using namespace std; 6 7 struct point 8 { 9 double x,y; 10 }a[55]; 11 double r; //半径 12 13 double dist_1point(double x0,double y0) //点到原点距离 14 { 15 return sqrt(x0*x0+y0*y0); 16 } 17 double dist_2point(double x1,double y1,double x2,double y2) //两点距离 18 { 19 return sqrt((x1-x2)*(x1-x2)+(y1-y2)*(y1-y2)); 20 } 21 double dist_line(double x1,double y1,double x2,double y2) //直线到原点距离 22 { 23 double A,B,C,dist; 24 A=y1-y2; 25 B=x1-x2; 26 C=x1*y2-x2*y1; 27 dist=fabs(C)/sqrt(A*A+B*B); //直线到原点距离公式 28 return dist; 29 } 30 double get_cos(double a,double b,double c) //余弦定理求角度 31 { 32 double angel=(b*b+c*c-a*a)/(2*b*c); 33 return angel; 34 } 35 point get_point(double x0,double y0) //计算点与原点的直线与圆相交的交点 36 { 37 double k; 38 point temp; 39 if(x0!=0) //若斜率存在 40 { 41 k=y0/x0; 42 temp.x=fabs(r)/sqrt(1+k*k); //判断是两个点的哪一个 43 if(x0<0) temp.x=-temp.x; //temp.x应该与x0同符号 44 temp.y=k*temp.x; 45 } 46 else //斜率不存在 47 { 48 temp.x=0; 49 if(y0>0) temp.y=r; //判断是两个点的哪一个 50 else temp.y=-r; //temp.y应该与y0同符号 51 } 52 return temp; 53 } 54 int fi(double x1,double y1,double x2,double y2) //关键啊啊啊!!!! 55 { 56 if (x1*y2-x2*y1>0) return 1; //判断是相加还是相减 57 else return -1; 58 } 59 60 double get_area(double x1,double y1,double x2,double y2) //三角剖分 61 { 62 int sign=fi(x1,y1,x2,y2); //判断三角形面积加还是减 63 double s; //总面积 64 double l=dist_line(x1,y1,x2,y2); //l = 直线ab与原点距离 65 double a=dist_1point(x1,y1); //a = 线段a长度 66 double b=dist_1point(x2,y2); //b = 线段b长度 67 double c=dist_2point(x1,y1,x2,y2); //c = 线段c长度 68 if(a==0 || b==0) 69 return 0; //若其中一条边为0,返回面积0 70 71 //第一种情况:三角形的两条边a、b全部短于半径。 72 if(a<=r && b<=r) 73 { 74 s=fabs(x1*y2-x2*y1)/2.0; 75 return s*sign; 76 } 77 78 //第二种情况:a、b两条边长于半径,l也长于半径。 79 else if(a>=r && b>=r && l>=r) 80 { 81 point t1=get_point(x1,y1); 82 point t2=get_point(x2,y2); 83 double d=dist_2point(t1.x,t1.y,t2.x,t2.y); 84 double sita1=acos(get_cos(d,r,r)); 85 double s=fabs(sita1*r*r/2.0); //扇形面积:s=θ*r*r/2 86 return s*sign; 87 } 88 89 //第三种情况:a、b两条边长于半径,但l短于半径,并且垂足落在这条边上。 90 else if(a>=r && b>=r && l<=r && (get_cos(a,b,c)<=0 || get_cos(b,a,c)<=0)) 91 { 92 point t1=get_point(x1,y1); 93 point t2=get_point(x2,y2); 94 double d=dist_2point(t1.x,t1.y,t2.x,t2.y); 95 double sita=acos(get_cos(d,r,r)); 96 s=fabs(sita*r*r/2.0); 97 return s*sign; 98 } 99 100 //第四种情况:a、b两条边长于半径,但l短于半径,且垂足没有落在这条边上。 101 else if(a>=r && b>=r && l<=r && (get_cos(a,b,c)>0 && get_cos(b,a,c)>0)) 102 { 103 double xx1,xx2,yy1,yy2; //点(x1,y1)与(x2,y2)组成的直线与圆的交点 104 if(x1!=x2) //若斜率存在 105 { 106 double k12=(y1-y2)/(x1-x2); 107 double b12=y1-k12*x1; 108 double a0=(1+k12*k12); 109 double b0=(2*k12*b12); 110 double c0=(b12*b12-r*r); 111 //化成一元二次方程,用公式求出两个交点 112 xx1=(-b0+sqrt(b0*b0-4*a0*c0))/(2*a0); 113 yy1=k12*xx1+b12; 114 xx2=(-b0-sqrt(b0*b0-4*a0*c0))/(2*a0); 115 yy2=k12*xx2+b12; 116 } 117 else //若斜率不存在 x1==x2 118 { 119 xx1=x1; 120 xx2=x1; 121 yy1=sqrt(r*r-x1*x1); 122 yy2=-sqrt(r*r-x1*x1); 123 } 124 point t1=get_point(x1,y1); //(x1,y1),(0,0)组成直线与圆的交点 125 point t2=get_point(x2,y2); //(x2,y2),(0,0)组成直线与圆的交点 126 double d1=dist_2point(xx1,yy1,xx2,yy2); //直线1与原点距离 127 double d2=dist_2point(t1.x,t1.y,t2.x,t2.y); //直线2与原点距离 128 double sita1=acos(get_cos(d1,r,r)); //小的扇形弧度 129 double sita2=acos(get_cos(d2,r,r)); //大的扇形弧度 130 double s1=fabs(sita1*r*r/2.0); //小的扇形面积 131 double s2=fabs(sita2*r*r/2.0); //大的扇形面积 132 double s3=fabs(xx1*yy2-xx2*yy1)/2.0; //三角形面积 133 s=s2+s3-s1; //相交面积 134 return s*sign; 135 } 136 137 //第五种情况1:三角形的两条边一条长于半径,另外一条短于半径 138 else if(a>=r && b<=r) //a长于半径,b短于半径 139 { 140 double xxx,yyy; 141 if(x1!=x2) //斜率存在 142 { 143 double k12=(y1-y2)/(x1-x2); 144 double b12=y1-k12*x1; 145 double a0=(1+k12*k12); 146 double b0=(2*k12*b12); 147 double c0=(b12*b12-r*r); 148 //化成一元二次方程,用公式求出两个交点 149 double xx1=(-b0+sqrt(b0*b0-4*a0*c0))/(2*a0); 150 double yy1=k12*xx1+b12; 151 double xx2=(-b0-sqrt(b0*b0-4*a0*c0))/(2*a0); 152 double yy2=k12*xx2+b12; 153 //判断两个交点中的哪一个,应在(x1,x2)两点之间 154 if(x1<=xx1 && xx1<=x2 || x2<=xx1 && xx1<=x1) {xxx=xx1; yyy=yy1;} 155 else {xxx=xx2; yyy=yy2;} 156 } 157 else //斜率不存在 x1==x2 158 { 159 double xx1=x1; 160 double yy1=-sqrt(r*r-x1*x1); 161 double yy2=sqrt(r*r-x1*x1); 162 //判断两个交点中的哪一个,应在(y1,y2)两点之间 163 if(y1<=yy1 && yy1<=y2 || y2<=yy1 && yy1<=y1) {yyy=yy1; xxx=xx1;} 164 else {yyy=yy2; xxx=xx1;} 165 } 166 //判断交点(该点已判断方向) 167 point t1=get_point(x1,y1); 168 double ddd=dist_2point(t1.x,t1.y,xxx,yyy); 169 double sita1=acos(get_cos(ddd,r,r)); 170 double s1=fabs(sita1*r*r/2.0); 171 double s3=fabs(xxx*y2-yyy*x2)/2.0; 172 s=s1+s3; //相交面积 173 return s*sign; 174 } 175 176 //第五种情况2:三角形的两条边一条长于半径,另外一条短于半径,与上述同理!!! 177 else if(a<=r && b>=r) //a短于半径,b长于半径 178 { 179 double xxx,yyy; //与上述同理!!! 180 if(x1-x2!=0) 181 { 182 double k12=(y1-y2)/(x1-x2); 183 double b12=y1-k12*x1; 184 double a0=(1+k12*k12); 185 double b0=(2*k12*b12); 186 double c0=(b12*b12-r*r); 187 double xx1=(-b0+sqrt(b0*b0-4*a0*c0))/(2*a0); 188 double yy1=k12*xx1+b12; 189 double xx2=(-b0-sqrt(b0*b0-4*a0*c0))/(2*a0); 190 double yy2=k12*xx2+b12; 191 if(x1<=xx1 && xx1<=x2 || x2<=xx1 && xx1<=x1) {xxx=xx1; yyy=yy1;} 192 else {xxx=xx2; yyy=yy2;} 193 } 194 else 195 { 196 double yy1=-sqrt(r*r-x1*x1); 197 double yy2=sqrt(r*r-x1*x1); 198 double xx1=x1; 199 if(y1<=yy1 && yy1<=y2 || y2<=yy1 && yy1<=y1) {yyy=yy1; xxx=xx1;} 200 else {yyy=yy2; xxx=xx1;} 201 } 202 point t1=get_point(x2,y2); 203 double ddd=dist_2point(t1.x,t1.y,xxx,yyy); 204 double sita1=acos(get_cos(ddd,r,r)); 205 double s1=fabs(sita1*r*r/2.0); 206 double s3=fabs(xxx*y1-yyy*x1)/2.0; 207 s=s1+s3; 208 return s*sign; 209 } 210 else return 0; 211 } 212 213 int main() 214 { 215 int i,n; 216 double area; 217 while(scanf("%lf",&r)!=EOF) 218 { 219 scanf("%d",&n); 220 for(i=0;i<n;i++) 221 { 222 scanf("%lf%lf",&a[i].x,&a[i].y); 223 } 224 a[n]=a[0]; area=0; 225 for(i=0;i<n;i++) //原点与其中两个点组成三角形来判断 226 { 227 area+=get_area(a[i].x,a[i].y,a[i+1].x,a[i+1].y); 228 } 229 printf("%.2lf\n",fabs(area)); 230 } 231 232 return 0; 233 }
Minimal Circle ZOJ 1450
一道模板题:
给出一系列的点,求一个最小圆将其覆盖的半径与圆心。。。
算法详解:http://wenku.baidu.com/view/584b6d3e5727a5e9856a610d.html(没时间了。。。)
接下来是模板:参见:http://www.cnblogs.com/DreamUp/archive/2010/10/12/1849072.html
1 代码 2 3 #include <stdio.h> 4 #include <math.h> 5 const int maxn = 1005; 6 const double eps = 1e-6; 7 struct TPoint { 8 double x, y; 9 TPoint operator-(TPoint & a) { 10 TPoint p1; 11 p1.x = x - a.x; 12 p1.y = y - a.y; 13 return p1; 14 } 15 }; 16 struct TCircle { 17 double r; 18 TPoint centre; 19 }; 20 struct TTriangle { 21 TPoint t[3]; 22 }; 23 TCircle c; 24 TPoint a[maxn]; 25 double distance(TPoint p1, TPoint p2) { 26 TPoint p3; 27 p3.x = p2.x - p1.x; 28 p3.y = p2.y - p1.y; 29 return sqrt(p3.x * p3.x + p3.y * p3.y); 30 } 31 double triangleArea(TTriangle t) { 32 TPoint p1, p2; 33 p1 = t.t[1] - t.t[0]; 34 p2 = t.t[2] - t.t[0]; 35 return fabs(p1.x * p2.y - p1.y * p2.x) / 2; 36 } 37 TCircle circumcircleOfTriangle(TTriangle t) { 38 //三角形的外接圆 39 TCircle tmp; 40 double a, b, c, c1, c2; 41 double xA, yA, xB, yB, xC, yC; 42 a = distance(t.t[0], t.t[1]); 43 b = distance(t.t[1], t.t[2]); 44 c = distance(t.t[2], t.t[0]); 45 //根据S = a * b * c / R / 4;求半径R 46 tmp.r = a * b * c / triangleArea(t) / 4; 47 48 xA = t.t[0].x; 49 yA = t.t[0].y; 50 xB = t.t[1].x; 51 yB = t.t[1].y; 52 xC = t.t[2].x; 53 yC = t.t[2].y; 54 c1 = (xA * xA + yA * yA - xB * xB - yB * yB) / 2; 55 c2 = (xA * xA + yA * yA - xC * xC - yC * yC) / 2; 56 57 tmp.centre.x = (c1 * (yA - yC) - c2 * (yA - yB)) / 58 ((xA - xB) * (yA - yC) - (xA - xC) * (yA - yB)); 59 tmp.centre.y = (c1 * (xA - xC) - c2 * (xA - xB)) / 60 ((yA - yB) * (xA - xC) - (yA - yC) * (xA - xB)); 61 62 return tmp; 63 } 64 65 TCircle MinCircle2(int tce, TTriangle ce) { 66 TCircle tmp; 67 if (tce == 0) tmp.r = -2; 68 else if (tce == 1) { 69 tmp.centre = ce.t[0]; 70 tmp.r = 0; 71 } else if (tce == 2) { 72 tmp.r = distance(ce.t[0], ce.t[1]) / 2; 73 tmp.centre.x = (ce.t[0].x + ce.t[1].x) / 2; 74 tmp.centre.y = (ce.t[0].y + ce.t[1].y) / 2; 75 } else if (tce == 3) tmp = circumcircleOfTriangle(ce); 76 return tmp; 77 } 78 79 void MinCircle(int t, int tce, TTriangle ce) { 80 int i, j; 81 TPoint tmp; 82 c = MinCircle2(tce, ce); 83 if (tce == 3) return; 84 for (i = 1; i <= t; i++) { 85 if (distance(a[i], c.centre) > c.r) { 86 ce.t[tce] = a[i]; 87 MinCircle(i - 1, tce + 1, ce); 88 tmp = a[i]; 89 for (j = i; j >= 2; j--) { 90 a[j] = a[j - 1]; 91 } 92 a[1] = tmp; 93 } 94 } 95 } 96 97 int main() { 98 //freopen("circle.in", "r", stdin); 99 //freopen("out.txt", "w", stdout); 100 int n,i; 101 while (scanf("%d", &n) != EOF && n) { 102 for (i = 1; i <= n; i++) 103 scanf("%lf%lf", &a[i].x, &a[i].y); 104 TTriangle ce; 105 MinCircle(n, 0, ce); 106 printf("%.2lf %.2lf %.2lf\n", c.centre.x, c.centre.y, c.r); 107 } 108 return 0; 109 }
接下来这道题总算好理解了!!!
思路只要画出图就清楚了。。。
具体思路参见:http://blog.csdn.net/acm_cxlove/article/details/7896110
详细解释看代码:
1 #include<iostream> 2 #include<algorithm> 3 #include<math.h> 4 #include<stdio.h> 5 using namespace std; 6 struct line 7 { 8 double l; 9 double r; 10 }a[505]; 11 double add(double x,double y,double x1,double y1) 12 { 13 return sqrt((x1-x)*(x1-x)+(y1-y)*(y1-y));//两点之间的距离 14 } 15 bool comp1(line i,line j) 16 { 17 return i.l<j.l; 18 } 19 int main() 20 { 21 double L,R,a1,b1,x,y,r; 22 int n,i; 23 while(scanf("%d",&n)!=EOF&&n) 24 { 25 scanf("%lf%lf",&a1,&b1); 26 for(i=0;i<n;i++) 27 { 28 scanf("%lf%lf%lf",&x,&y,&r); 29 double aa1=add(x,y,a1,b1); 30 double ab1=asin(r/aa1); 31 double ab2=asin((a1-x)/aa1); 32 a[i].l=a1-(tan(ab1+ab2)*b1);//图参见链接 33 a[i].r=a1-(tan(ab2-ab1)*b1);//先把每一个圆投影所对应的坐标求出来 34 } 35 sort(a,a+n,comp1);//进行排序,贪心处理 36 L=a[0].l; 37 R=a[0].r; 38 for(i=1;i<n;i++) 39 { 40 if(a[i].l>R)//没有重合部分输出值 41 { 42 printf("%.2lf %.2lf\n",L,R); 43 L=a[i].l; 44 R=a[i].r; 45 } 46 else 47 R=max(R,a[i].r);//取较大值 48 } 49 printf("%.2lf %.2lf\n\n",L,R); 50 } 51 return 0; 52 }
【推荐】编程新体验,更懂你的AI,立即体验豆包MarsCode编程助手
【推荐】凌霞软件回馈社区,博客园 & 1Panel & Halo 联合会员上线
【推荐】抖音旗下AI助手豆包,你的智能百科全书,全免费不限次数
【推荐】博客园社区专享云产品让利特惠,阿里云新客6.5折上折
【推荐】轻量又高性能的 SSH 工具 IShell:AI 加持,快人一步