/*********************************************/
* *
* 计算几何学库函数 *
* *
* *
/*********************************************/
#include <stdlib.h>
#define infinity 1e20
#define EP 1e-10
struct TPoint{
float x,y;
};
struct TLineSeg{
TPoint a,b;
};
//求平面上两点之间的距离
float distance(TPoint p1,TPoint p2)
{
return(sqrt((p1.x-p2.x)*(p1.x-p2.x)+(p1.y-p2.y)*(p1.y-p2.y)));
}
/********************************************************
返回(P1-P0)*(P2-P0)的叉积。
若结果为正，则<P0,P1>在<P0,P2>的顺时针方向；
若为0则<P0,P1><P0,P2>共线；
若为负则<P0,P1>在<P0,P2>的在逆时针方向;
可以根据这个函数确定两条线段在交点处的转向,
比如确定p0p1和p1p2在p1处是左转还是右转，只要求
(p2-p0)*(p1-p0)，若<0则左转，>0则右转，=0则共线
*********************************************************/
float multiply(TPoint p1,TPoint p2,TPoint p0)
{
return((p1.x-p0.x)*(p2.y-p0.y)-(p2.x-p0.x)*(p1.y-p0.y));
}
//确定两条线段是否相交
int intersect(TLineSeg u,TLineSeg 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))&&
(multiply(v.a,u.b,u.a)*multiply(u.b,v.b,u.a)>=0)&&
(multiply(u.a,v.b,v.a)*multiply(v.b,u.b,v.a)>=0));
}
//判断点p是否在线段l上
int online(TLineSeg l,TPoint p)
{
return( (multiply(l.b,p,l.a)==0)&&( ((p.x-l.a.x)*(p.x-l.b.x)<0 )||( (p.y-l.a.y)*(p.y-l.b.y)<0 )) );
}
//判断两个点是否相等
int Euqal_Point(TPoint p1,TPoint p2)
{
return((abs(p1.x-p2.x)<EP)&&(abs(p1.y-p2.y)<EP));
}
//一种线段相交判断函数，当且仅当u,v相交并且交点不是u,v的端点时函数为true;
int intersect_A(TLineSeg u,TLineSeg v)
{
return((intersect(u,v))&&
(!Euqal_Point(u.a,v.a))&&
(!Euqal_Point(u.a,v.b))&&
(!Euqal_Point(u.b,v.a))&&
(!Euqal_Point(u.b,v.b)));
}
/*===============================================
判断点q是否在多边形Polygon内，
其中多边形是任意的凸或凹多边形，
Polygon中存放多边形的逆时针顶点序列
================================================*/
int InsidePolygon(int vcount,TPoint Polygon[],TPoint q)
{
int c=0,i,n;
TLineSeg l1,l2;
l1.a=q;
l1.b=q;
l1.b.x=infinity;
n=vcount;
for (i=0;i<vcount;i++)
{
l2.a=Polygon[i];
l2.b=Polygon[(i+1)%n];
if ( (intersect_A(l1,l2))||
(
(online(l1,Polygon[(i+1)%n]))&&
(
(!online(l1,Polygon[(i+2)%n]))&&
(multiply(Polygon[i],Polygon[(i+1)%n],l1.a)*multiply(Polygon[(i+1)%n],Polygon[(i+2)%n],l1.a)>0)
||
(online(l1,Polygon[(i+2)%n]))&&
(multiply(Polygon[i],Polygon[(i+2)%n],l1.a)*multiply(Polygon[(i+2)%n],Polygon[(i+3)%n],l1.a)>0)
)
)
) c++;
}
return(c%2!=0);
}
/********************************************/
* *
* 计算多边形的面积 *
* *
* 要求按照逆时针方向输入多边形顶点 *
* 可以是凸多边形或凹多边形 *
* *
/********************************************/
float area_of_polygon(int vcount,float x[],float y[])
{
int i;
float s;
if (vcount<3) return 0;
s=y[0]*(x[vcount-1]-x[1]);
for (i=1;i<vcount;i++)
s+=y[i]*(x[(i-1)]-x[(i+1)%vcount]);
return s/2;
}
/*====================================================================================
寻找凸包 graham 扫描法
PointSet为输入的点集；
ch为输出的凸包上的点集，按照逆时针方向排列;
n为PointSet中的点的数目
len为输出的凸包上的点的个数；
设下一个扫描的点PointSet[i]=P2,当前栈顶的两个点ch[top]=P1,ch[top-1]=P0,
如果P1P2相对于P0P1在点P1向左旋转(共线也不行)，则P0,P1一定是凸包上的点；
否则P1一定不是凸包上的点，应该将其出栈。
比如下面左图中点1就一定是凸包上的点，右图中点1就一定不是凸包上的点，因为
可以连接点0,2构成凸包的边
2
|
| _____2
1 1
/ /
____0/ ____0/
====================================================================================*/
void Graham_scan(TPoint PointSet[],TPoint ch[],int n,int &len)
{
int i,j,k=0,top=2;
TPoint tmp;
//选取PointSet中y坐标最小的点PointSet[k]，如果这样的点右多个，则取最左边的一个
for(i=1;i<n;i++)
if ((PointSet[i].y<PointSet[k].y)||((PointSet[i].y==PointSet[k].y)&&(PointSet[i].x<PointSet[k].x)))
k=i;
tmp=PointSet[0];
PointSet[0]=PointSet[k];
PointSet[k]=tmp; //现在PointSet中y坐标最小的点在PointSet[0]
for (i=1;i<n-1;i++) //对顶点按照相对PointSet[0]的极角从小到大进行排序，极角相同的按照距离PointSet[0]从近到远进行排序
{ k=i;
for (j=i+1;j<n;j++)
if ( (multiply(PointSet[j],PointSet[k],PointSet[0])>0)
||
( (multiply(PointSet[j],PointSet[k],PointSet[0])==0)
&&(distance(PointSet[0],PointSet[j])<distance(PointSet[0],PointSet[k])) )
) k=j;
tmp=PointSet[i];
PointSet[i]=PointSet[k];
PointSet[k]=tmp;
}
ch[0]=PointSet[0];
ch[1]=PointSet[1];
ch[2]=PointSet[2];
for (i=3;i<n;i++)
{
while (multiply(PointSet[i],ch[top],ch[top-1])>=0) top--;
ch[++top]=PointSet[i];
}
len=top+1;
}
扫一扫
评论
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
查看更多评论
添加红包
