您的位置:首页 > 其它

计算几何基础(模板)

2015-08-16 15:52 393 查看
1.多边形面积计算

double S(Point p[],int n)
{
double ans = 0;
p
= p[0];
for(int i=1;i<n;i++)
ans += cross(p[0],p[i],p[i+1]);
if(ans < 0) ans = -ans;
return ans / 2.0;
}


2.求凸包

bool cmp(Point A,Point B)
{
double k = cross(MinA,A,B);
if(k<0) return 0;
if(k>0) return 1;
return dist(MinA,A)<dist(MinA,B);
}

void Graham(Point p[],int n)
{
for(int i=1;i<n;i++)
if(p[i].y<p[0].y || (p[i].y == p[0].y && p[i].x < p[0].x))
swap(p[i],p[0]);
MinA = p[0];
p
= p[0];
sort(p+1,p+n,cmp);
stack[0] = p[0];
stack[1] = p[1];
top = 2;
for(int i=2;i<n;i++)
{
while(top >= 2 && cross(stack[top-2],stack[top-1],p[i])<=0) top--;
stack[top++] = p[i];
}
}


3.任意多边形求重心

Point Gravity(Point p[],int n)
{
Point O,t;
O.x = O.y = 0;
t.x = t.y = 0;
p
= p[0];
double A = 0;
for(int i=0; i<n; i++)
A += cross(O,p[i],p[i+1]);
A /= 2.0;
for(int i=0; i<n; i++)
{
t.x += (p[i].x + p[i+1].x) * cross(O,p[i],p[i+1]);
t.y += (p[i].y + p[i+1].y) * cross(O,p[i],p[i+1]);
}
t.x /= 6*A;
t.y /= 6*A;
return t;
}


4.求线段交点的坐标

bool Segment_crossing(Segment u,Segment v)   /** 判断线段是否相交 */
{
return((max(u.a.x,u.b.x)>=min(v.a.x,v.b.x))&&
(max(v.a.x,v.b.x)>=min(u.a.x,u.b.x))&&
(max(u.a.y,u.b.y)>=min(v.a.y,v.b.y))&&
(max(v.a.y,v.b.y)>=min(u.a.y,u.b.y))&&
(cross(v.a,u.b,u.a)*cross(u.b,v.b,u.a)>=0)&&
(cross(u.a,v.b,v.a)*cross(v.b,u.b,v.a)>=0));
}

/**求直线交点的坐标,如果没有交点返回NULL,否则返回交点p的地址*/
Point* CrossPoint(Segment u,Segment v)
{
Point p;
if(Segment_crossing(u,v))
{
p.x=(cross(v.b,u.b,u.a)*v.a.x-cross(v.a,u.b,u.a)*v.b.x)/(cross(v.b,u.b,u.a)-cross(v.a,u.b,u.a));
p.y=(cross(v.b,u.b,u.a)*v.a.y-cross(v.a,u.b,u.a)*v.b.y)/(cross(v.b,u.b,u.a)-cross(v.a,u.b,u.a));
return &p;
}
return NULL;
}


5.三角形外接圆的半径与圆心

Point Circle_Point(Point A,Point B,Point C)
{
double a=dist(B,C);
double b=dist(A,C);
double c=dist(A,B);
double p=(a+b+c)/2.0;
double S=sqrt(p*(p-a)*(p-b)*(p-c));
R=(a*b*c)/(4*S);    //三角形外接圆的半径为R

double t1=(A.x*A.x+A.y*A.y-B.x*B.x-B.y*B.y)/2;
double t2=(A.x*A.x+A.y*A.y-C.x*C.x-C.y*C.y)/2;
Point center;
center.x=(t1*(A.y-C.y)-t2*(A.y-B.y))/((A.x-B.x)*(A.y-C.y)-(A.x-C.x)*(A.y-B.y));
center.y=(t1*(A.x-C.x)-t2*(A.x-B.x))/((A.y-B.y)*(A.x-C.x)-(A.y-C.y)*(A.x-B.x));
return center;
}


6.旋转卡壳求凸包的直径,也就是平面最远点对,p[]为凸包的点集

double rotating_calipers(Point p[],int n)
{
int k = 1;
double ans = 0;
p
= p[0];
for(int i=0;i<n;i++)
{
while(fabs(cross(p[i],p[i+1],p[k])) < fabs(cross(p[i],p[i+1],p[k+1])))
k = (k+1) % n;
ans = max(ans, max(dist(p[i],p[k]),dist(p[i+1],p[k])));
}
return ans;
}


7.求凸包的宽度

double rotating_calipers(Point p[],int n)
{
int k = 1;
double ans = 0x7FFFFFFF;
p
= p[0];
for(int i=0;i<n;i++)
{
while(fabs(cross(p[i],p[i+1],p[k])) < fabs(cross(p[i],p[i+1],p[k+1])))
k = (k+1) % n;
double tmp = fabs(cross(p[i],p[i+1],p[k]));
double d   = dist(p[i],p[i+1]);
ans = min(ans,tmp/d);
}
return ans;
}


点的向量,叉乘,点乘

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(){}
Point(double _x,double _y)
{
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;
}
};
struct Line
{
Point s,e;
Line(){}
Line(Point _s,Point _e)
{
s = _s;e = _e;
}
};
//判断线段相交
bool inter(Line l1,Line l2)
{
return
max(l1.s.x,l1.e.x) >= min(l2.s.x,l2.e.x) &&
max(l2.s.x,l2.e.x) >= min(l1.s.x,l1.e.x) &&
max(l1.s.y,l1.e.y) >= min(l2.s.y,l2.e.y) &&
max(l2.s.y,l2.e.y) >= min(l1.s.y,l1.e.y) &&
sgn((l2.s-l1.s)^(l1.e-l1.s))*sgn((l2.e-l1.s)^(l1.e-l1.s)) <= 0 &&
sgn((l1.s-l2.s)^(l2.e-l2.s))*sgn((l1.e-l2.s)^(l2.e-l2.s)) <= 0;
}
double dist(Point a,Point b)
{
return sqrt((b-a)*(b-a));
}
内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息
标签: