您的位置:首页 > 其它

二维计算几何模板

2014-03-06 12:28 766 查看
Ps:还不完整,可能有错

#include<cstdio>
#include<cstring>
#include<cmath>
#include<algorithm>
using namespace std;
//-----------------基础定义-----------------
struct point
{
double x,y;
point(double x=0,double y=0):x(x),y(y){}//构造
};

typedef point vector;

//各种符号重载
vector operator + (vector a,vector b) {return vector(a.x+b.x,a.y+b.y);}

vector operator - (point a,point b) {return vector(a.x-b.x,a.y-b.y);}

vector operator * (vector a,double p) {return vector(a.x*p,a.y*p);}

vector operator / (vector a,double p) {return vector(a.x/p,a.y/p);}

bool operator < (const point& a,const point& b) {return a.x<b.x||(a.x==b.x&&a.y<b.y);}

//高精度比较
const double eps=1e-10;
int dcmp(double x) {if(fabs(x)<eps) return 0; else return x<0?-1:1;}

bool operator == (const point &a,const point &b)
{ return dcmp(a.x-b.x)==0&&dcmp(a.y-b.y)==0;}
//----------------------------------------

//-----------------向量-------------------
//点积
double dot(vector a,vector b) { return a.x*b.x+a.y*b.y;}

//向量模
double length(vector a) { return sqrt(dot(a,a));}

//两向量夹角
double angle(vector a,vector b) {return acos(dot(a,b)/length(a)/length(b));}

//叉积
double cross(vector a,vector b) {return a.x*b.y-a.y*b.x;}

//三点平行四边形面积,除以2为三角形面积
double area2(point a,point b,point c) {return cross(b-a,c-a);}

//向量旋转,rad>0逆时针,rad<0顺时针
vector rotate(vector a,double rad)
{return vector(a.x*cos(rad)-a.y*sin(rad),a.x*sin(rad)+a.y*cos(rad));}

//向量逆时针旋转90度
vector normal(vector a) {double l=length(a);return vector(-a.y/l,a.x/l);}
//---------------------------------------

//-----------------直线------------------
//直线交点,调用前确保有交点,叉积判断
point GetLineIntersection(point p,vector v,point q,vector w)
{
vector u=p-q;
double t=cross(w, u)/cross(v, w);
return p+v*t;
}

//点到直线距离
double DstanceToLine(point p,point a,point b)
{
vector v1=b-1,v2=p-a;
return fabs(cross(v1, v2))/length(v1);
}

//求p点在直线ab上的投影
point GetLineProjection(point p,point a,point b)
{
vector v=b-a;
return a+v*(dot(v,p-a))/dot(v,v);
}
//---------------------------------------

//---------------线段---------------
//点到线段的距离
double distancetosegma(point p,point a,point b)
{
if(a==b) return length(p-a);
vector v1=b-a,v2=p-a,v3=p-b;
if(dcmp(cross(v1, v2))<0) return length(v2);
else if(dcmp(dot(v1,v3))>0) return length(v1);
else return fabs(cross(v1, v2))/length(v1);
}

//线段规范相交,交点不含端点
bool segmentproperintersection(point a1,point a2,point b1,point b2)
{
double c1=cross(a2-a1, b1-a1),c2=cross(a2-a1, b2-a1),
c3=cross(b2-b1, a1-b1),c4=cross(b2-b1, a2-b1);
return dcmp(c1)*dcmp(c2)<0&&dcmp(c3)*dcmp(c4)<0;
}

//p点是否在线段ab上,<=0表示p点可能在端点上
bool onsegment(point p,point a1,point a2)
{return dcmp(cross(a1-p, a2-p))==0&&dcmp(dot(a1-p, a2-p))<=0;}
//二者结合判断线段相交

//求两线段交点,线段共线(cross()==0)不行
point GetSegmentIntesection(point a1,point a2,point b1,point b2)
{
double s1=fabs(cross(b1-a1, a2-a1)),
s2=fabs(cross(a2-a1, b2-a1));
point p;
p.x=(s1*b2.x+s2*b1.x)/(s1+s2);
p.y=(s1*b2.y+s2*b1.y)/(s1+s2);
return p;
}
//-----------------------------------

//--------------n边形-----------------
//极角排序
bool cmpPolarange(point a,point b)
{
point c(1,1);//按c点排序
return dcmp(cross(a-c, b-c))>0;
}

//n边形有向面积,逆时针?
double polygenarea(point *p,int n)
{
double area=0;
for(int i=1;i<n-1;i++)
area+=cross(p[i]-p[0], p[i+1]-p[0]);
return area/2;
}

//判断点是否在多边形内,凹凸都适合,逆时针
int isPointInPolygon(point p,point *poly,int n)
{
int wn=0,i;
for(i=0;i<n;i++)
{
if(onsegment(p, poly[i], poly[(i+1)%n])) return -1;//在边界上
int k=dcmp(cross(poly[(i+1)%n]-poly[i], p-poly[i]));
int d1=dcmp(poly[i].y-p.y);
int d2=dcmp(poly[(i+1)%n].y-p.y);
if(k>0&&d1<=0&&d2>0) wn++;
if(k<0&&d2<=0&&d1>0) wn--;
}
if(wn!=0) return 1;//在内部
return 0;//不在内部
}
/////////////////////////////
//Andrew算法求凸包,确保没有重复点
bool cmpXY(point a,point b)
{
if(a.x==b.x) return a.y<b.y;
else return a.x<b.x;
}
//若让ch不保存凸包边上的点则将两个<=改为<
int ConvexHull(point *p,int n,point *ch)
{
sort(p,p+n,cmpXY);
int m=0,i;
for(i=0;i<n;i++)
{
while(m>1&&cross(ch[m-1]-ch[m-2], p[i]-ch[m-2])<=0) m--;
ch[m++]=p[i];
}
int k=m;
for(i=n-2;i>=0;i--)
{
while(m>k&&cross(ch[m-1]-ch[m-2], p[i]-ch[m-2])<=0) m--;
ch[m++]=p[i];
}
if(n>1) m--;
return m;
}
/////////////////////////////
//------------------------------------

//-------------半平面------------------
//有向直线
struct line
{
point p;
vector v;
double ang;//向量v的极角
line(){}
line(point p,vector v):p(p),v(v){ang=atan2(v.y,v.x);}
bool operator < (const line & l) const//极角排序
{return ang<l.ang;}
};
//点p在直线l的左边
bool OnLeft(line l,point p) {return cross(l.v, p-l.p)>0;}
//直线a,b的交点
point GetIntersection(line a,line b)
{
vector u=a.p-b.p;
double t=cross(b.v, u)/cross(a.v, b.v);
return a.p+a.v*t;
}
//求半平面交
int HalfplaneIntersection(line *l,int n,point *poly)
{
sort(l,l+n);
int first,last;//双端队列指针
point *p=new point
;
line *q=new line
;
q[first=last=0]=l[0];
int i;
for(i=1;i<n;i++)
{
while(first<last&&!OnLeft(l[i], p[last-1])) last--;
while(first<last&&!OnLeft(l[i], p[first])) first++;
q[++last]=l[i];
if(fabs(cross(q[last].v, q[last-1].v))<eps)
{
last--;
if(OnLeft(q[last], l[i].p)) q[last]=l[i];
}
if(first<last) p[last-1]=GetIntersection(q[last-1], q[last]);
}
while(first<last&&!OnLeft(q[first], p[last-1])) last--;
if(last-first<=1) return 0;
p[last]=GetIntersection(q[last], q[first]);

int m=0;
for(i=first;i<=last;i++) poly[m++]=p[i];
return m;
}
//------------------------------------

int main()
{
return 0;
}
内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息
标签:  计算几何 二维 模板