博客
关于我
强烈建议你试试无所不能的chatGPT,快点击我
计算几何
阅读量:4286 次
发布时间:2019-05-27

本文共 6003 字,大约阅读时间需要 20 分钟。

此文章是一些自己未接触过的计算几何知识,不是太明白。。。。。

 &&     本题代码来自 大牛  :

/**    最小圆覆盖 (随机增量法)    求给你n个点,让你用一个最小的圆覆盖所有的点,求这个圆的坐标&&半径**/#include 
#include
#include
#include
#include
using namespace std;const double eps = 1e-12;const int len = 100010;struct Point{ double x,y;} p[len];struct Line{ Point a,b;};int dbcmp(double n){ return n < -eps ? -1 : n > eps;}double dis(Point a, Point b){ return sqrt((a.x-b.x) * ( a.x-b.x) + ( a.y-b.y) * ( a.y-b.y));}//求两直线的交点Point intersection(Line u,Line v){ Point ret=u.a; double t=((u.a.x-v.a.x)*(v.b.y-v.a.y)-(u.a.y-v.a.y)*(v.b.x-v.a.x)) /((u.a.x-u.b.x)*(v.b.y-v.a.y)-(u.a.y-u.b.y)*(v.b.x-v.a.x)); ret.x+=(u.b.x-u.a.x)*t; ret.y+=(u.b.y-u.a.y)*t; return ret;}//三角形外接圆圆心(外心)Point center(Point a,Point b,Point c){ Line u,v; u.a.x=(a.x+b.x)/2; u.a.y=(a.y+b.y)/2; u.b.x=u.a.x+(u.a.y-a.y); u.b.y=u.a.y-(u.a.x-a.x); v.a.x=(a.x+c.x)/2; v.a.y=(a.y+c.y)/2; v.b.x=v.a.x+(v.a.y-a.y); v.b.y=v.a.y-(v.a.x-a.x); return intersection(u,v);}void min_cir(Point * p, int n, Point & cir, double &r){ random_shuffle(p, p + n); cir = p[0]; r = 0; for(int i = 1; i < n; ++i) { if(dbcmp(dis(p[i],cir) -r) <= 0)continue; cir = p[i]; r = 0; for(int j =0; j < i ; ++j) { if(dbcmp(dis(p[j], cir) -r) <= 0 )continue; cir.x = (p[i].x + p[j].x) /2; cir.y = (p[i].y + p[j].y) /2; r = dis( cir, p[j]); for(int k =0; k < j; ++k) { if(dbcmp( dis(p[k], cir) -r) <= 0) continue; cir = center(p[i], p[j], p[k]); r = dis( cir, p[k]); } } }}int main(){ int n; while( ~scanf("%d", &n), n) { for (int i = 0; i < n; ++i) { scanf("%lf%lf", &p[i].x, &p[i].y); } Point cir; double r; min_cir(p, n, cir, r); printf("%.2lf %.2lf %.2lf\n", cir.x, cir.y, r); }}

题意:求两凸包之间的距离;摘自 kuangbin 大牛:

#include 
#include
#include
#include
#include
using namespace std;const double eps = 1e-8;int sgn(double x){ if(fabs(x) < eps)return 0; if(x < 0)return -1; else return 1;}struct Point{ double x,y; Point(double _x = 0.0,double _y = 0.0){ x = _x;y = _y; } Point operator -(const Point &b)const{ return Point(x - b.x, y - b.y); } double operator ^(const Point &b)const{ return x*b.y - y*b.x; } double operator *(const Point &b)const{ return x*b.x + y*b.y; } void input(){ scanf("%lf%lf",&x,&y); }};struct Line{ Point s,e; Line(){} Line(Point _s,Point _e){ s = _s; e = _e; }};double dist(Point a,Point b){ return sqrt((a-b)*(a-b));}Point NearestPointToLineSeg(Point P,Line L){ Point result; double t = ((P-L.s)*(L.e-L.s))/((L.e-L.s)*(L.e-L.s)); if(t >= 0 && t <= 1){ result.x = L.s.x + (L.e.x - L.s.x)*t; result.y = L.s.y + (L.e.y - L.s.y)*t; } else{ if(dist(P,L.s) < dist(P,L.e)) result = L.s; else result = L.e; } return result;}/* * 求凸包,Graham算法 * 点的编号0~n-1 * 返回凸包结果Stack[0~top-1]为凸包的编号 */const int MAXN = 10010;Point list[MAXN];int Stack[MAXN],top;//相对于list[0]的极角排序bool _cmp(Point p1,Point p2){ double tmp = (p1-list[0])^(p2-list[0]); if(sgn(tmp) > 0)return true; else if(sgn(tmp) == 0 && sgn(dist(p1,list[0]) - dist(p2,list[0])) <= 0) return true; else return false;}void Graham(int n){ Point p0; int k = 0; p0 = list[0]; //找最下边的一个点 for(int i = 1;i < n;i++){ if( (p0.y > list[i].y) || (p0.y == list[i].y && p0.x > list[i].x) ){ p0 = list[i]; k = i; } } swap(list[k],list[0]); sort(list+1,list+n,_cmp); if(n == 1){ top = 1; Stack[0] = 0; return; } if(n == 2){ top = 2; Stack[0] = 0; Stack[1] = 1; return ; } Stack[0] = 0; Stack[1] = 1; top = 2; for(int i = 2;i < n;i++){ while(top > 1 && sgn((list[Stack[top-1]]-list[Stack[top-2]])^(list[i]-list[Stack[top-2]])) <= 0) top--; Stack[top++] = i; }}//点p0到线段p1p2的距离double pointtoseg(Point p0,Point p1,Point p2){ return dist(p0,NearestPointToLineSeg(p0,Line(p1,p2)));}//平行线段p0p1和p2p3的距离double dispallseg(Point p0,Point p1,Point p2,Point p3){ double ans1 = min(pointtoseg(p0,p2,p3),pointtoseg(p1,p2,p3)); double ans2 = min(pointtoseg(p2,p0,p1),pointtoseg(p3,p0,p1)); return min(ans1,ans2);}//得到向量a1a2和b1b2的位置关系double Get_angle(Point a1,Point a2,Point b1,Point b2){ Point t = b1 - ( b2 - a1 ); return (a2-a1)^(t-a1);}//旋转卡壳,求两个凸包的最小距离double rotating_calipers(Point p[],int np,Point q[],int nq){ int sp = 0, sq = 0; for(int i = 0;i < np;i++) if(sgn(p[i].y - p[sp].y) < 0) sp = i; for(int i = 0;i < nq;i++) if(sgn(q[i].y - q[sq].y) > 0) sq = i; double tmp; double ans = 1e99; for(int i = 0;i < np;i++){ while(sgn(tmp = Get_angle(p[sp],p[(sp+1)%np],q[sq],q[(sq+1)%nq])) < 0 ) sq = (sq + 1)%nq; if(sgn(tmp) == 0) ans = min(ans,dispallseg(p[sp],p[(sp+1)%np],q[sq],q[(sq+1)%nq])); else ans = min(ans,pointtoseg(q[sq],p[sp],p[(sp+1)%np])); sp = (sp+1)%np; } return ans;}double solve(Point p[],int n,Point q[],int m){ return min(rotating_calipers(p,n,q,m),rotating_calipers(q,m,p,n));}Point p[MAXN],q[MAXN];int main(){ int n,m; while(scanf("%d%d",&n,&m)==2) { if(n == 0 && m == 0)break; for(int i = 0;i < n;i++) list[i].input(); Graham(n); n = top; for(int i = 0;i < n;i++) p[i] = list[Stack[i]]; for(int i = 0;i < m;i++) list[i].input(); Graham(m); m = top; for(int i = 0;i < m;i++) q[i] = list[Stack[i]]; printf("%.5lf\n",solve(p,n,q,m)); } return 0;}

转载地址:http://qcsgi.baihongyu.com/

你可能感兴趣的文章
MySql数据类型
查看>>
MySql简单sql使用
查看>>
"未能加载文件或程序集“MySql.Data, Version=6.9.3.0”或它的某一个依赖项。
查看>>
CodeFirst for MySql
查看>>
Code Frist for Mysql 实例
查看>>
Visual Studio 开源控件扩展 NuGet 常用命令及常用组件
查看>>
mysql局域网访问设置
查看>>
UEditor 编辑器跨域上传解决方法
查看>>
VisualSVN Server搭建SVN服务器
查看>>
AngularJs directive指令详解
查看>>
AngularJs directive-scope
查看>>
AngularJs directive-link实例
查看>>
Js实现Base64编码、解码
查看>>
AngularJs directive-scope双向绑定方法处理-实例2
查看>>
AngularJs Ajax分页控件
查看>>
LocalDB数据库修改排序规则,修复汉字变问号
查看>>
C# Json序列化工具--Newtonsoft.Json简介和使用
查看>>
EntityFramework中Json序列化的循环引用问题解决--Newtonsoft.Json
查看>>
AngularJs----ng-class
查看>>
Bootstrap3 datetimepicker控件的使用
查看>>