当前位置：文档之家› 把点集凸包化Granham-Scan算法(使用水平序)概要

把点集凸包化Granham-Scan算法(使用水平序)概要

把点集凸包化Granham-Scan算法（使用水平序）

#include

using nameaace std;

const int M=100000+5;

struct Point

{

double x,y;

} p[M];

double dis(Point A,Point B)

{

return sqrt((B.x-A.x)*(B.x-A.x)+(B.y-A.y)*(B.y-A.y));

}

bool cmp(Point a,Point b)

{

if (a.x

if (a.x>b.x) return false;

if (a.y

return false;

}

double Xdet(Point A,Point B,Point C)

{

double x1,x2,y1,y2;

x1=B.x-A.x; y1=B.y-A.y; x2=C.x-A.x; y2=C.y-A.y;

return x1*y2-x2*y1; //大于0在左手边，逆时针

}

bool bo[M];

Point pp[M];

int stack[M]; //from 1 to t

void Gram_Scan(Point *p,int &n) //p从1-n,把点集土包化 n*log(n); {

int i,t;

sort(p+1,p+1+n,cmp);

for(t=0,i=1;i<=n;i++){ //合并排序

if (i>1 && p[i].x==p[i-1].x && p[i].y==p[i-1].y) continue;

p[++t]=p[i];

}

n=t; t=0;

memset(bo+1,true,n*sizeof(bo[0]));

if (n>0){

stack[++t]=1;

bo[stack[t]]=false;

}

if (n>1){

stack[++t]=2;

bo[stack[t]]=false;

}

if (n>2){

for(i=3;i

if(bo[i] && Xdet(p[stack[t-1]],p[stack[t]],p[i])>=0) stack[++t]=i,bo[i]=false;

else{

while(t>=2&& Xdet(p[stack[t-1]],p[stack[t]],p[i])<0) bo[stack[t]]=true,t--;

stack[++t]=i;

bo[stack[t]]=false;

}

for(i=n;i>=1;i--)

if (bo[i] && Xdet(p[stack[t-1]],p[stack[t]],p[i])>=0) stack[++t]=i,bo[i]=false;

else{

while(t>=2&&Xdet(p[stack[t-1]],p[stack[t]],p[i])<0) bo[stack[t]]=true,t--;

stack[++t]=i;

bo[stack[t]]=false;

}

t--;

}

for(i=1;i<=t;i++) pp[i]=p[stack[i]];

memcpy(p+1,pp+1,t*sizeof(Point));

n=t;

}

求凸点集的面积

double area(Point *p,int n)

{

double sum=0;

int i;

p[n+1]=p[1];

for(i=1;i<=n;i++)

sum+=p[i].x*p[i+1].y-p[i].y*p[i+1].x;

return sum/2.0;

}

const double EPS=1e-9;

const double maxn=1e6-1;

const double PI=acos(-1.);

const int M=100+5;

const double offset=11213;

#define zero(x) (((x)>0?(x):-(x))

struct Point{

double x,y;

} p[M];

struct LINESEG{

Point st,ed;

};

struct LINE{

double a,b,c;

};

double dist(Point a,Point b){

return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y));

}

double multiply(Point a,Point b,Point c){

return (a.x-c.x)*(b.y-c.y)-(b.x-c.x)*(a.y-c.y);

}

double dotmultiply(Point a,Point b,Point c){

return (a.x-c.x)*(b.x-c.x)+(a.y-c.y)*(b.y-c.y);

}

bool online(LINESEG L,Point q){

return fabs(multiply(L.ed,q,L.st))

LINE makeline(Point p1,Point p2){

LINE tl;

tl.a=p2.y-p1.y; tl.b=p1.x-p2.x;

tl.c=p1.y*p2.x-p1.x*p2.y;

if (tl.a

return tl;

}

inline double MAX(double x,double y){return x>y?x:y;}

inline double MIN(double x,double y){return x

bool intersect(LINESEG u,LINESEG v){

return MAX(u.st.x,u.ed.x)>=MIN(v.st.x,v.ed.x)

&& MAX(v.st.x,v.ed.x)>=MIN(u.st.x,u.ed.x)

&& MAX(u.st.y,u.ed.y)>=MIN(v.st.y,v.ed.y)

&& MAX(v.st.y,v.ed.y)>=MIN(u.st.y,u.ed.y)

&& multiply(v.st,u.ed,u.st)*multiply(u.ed,v.ed,u.st)>=0

&& multiply(u.st,v.ed,v.st)*multiply(v.ed,u.ed,v.st)>=0;

}

bool intersect_A(LINESEG u,LINESEG v){

return intersect(u,v) && !online(u,v.st) && !online(u,v.ed)

&& !online(v,u.st) && !online(u,v.ed);

}

bool intersect_L(LINESEG u,LINESEG v){

return multiply(u.st,v.ed,v.st)*multiply(v.ed,u.ed,v.st)>=0;

}

bool intersection_line(LINE L1,LINE L2,Point &inter){ double det=L1.a*L2.b-L2.a*L1.b;

if (fabs(det)

inter.x=(L2.c*L1.b-L1.c*L2.b)/det;

inter.y=(L2.a*L1.c-L1.a*L2.c)/det;

return true;

}

bool intersection_lineseg(LINESEG u,LINESEG v,Point &inter){ LINE L1,L2;

L1=makeline(u.st,u.ed);

L2=makeline(v.st,v.ed);

if (intersection_line(L1,L2,inter)) return online(u,inter) && online(v,inter);

else return false;

}

bool intersection_LLS(LINE L1,LINESEG u,Point &inter){ LINE L2;

L2=makeline(u.st,u.ed);

if (intersection_line(L1,L2,inter))

return online(u,inter);

else return false;

}

int Inside(Point q,int n)

{

Point q2;

const int on_edge=0;

int i=0,count;

p[n]=p[0];

while (i

for

(count=0,i=0,q2.x=rand()+offset,q2.y=rand()+offset;i

if (zero(multiply(q,p[i],p[i+1])) && (p[i].x-q.x)*(p[i+1].x-q.x)

else if (zero(multiply(q,q2,p[i]))) break;

else if (multiply(q,p[i],q2)*multiply(q,p[i+1],q2)<-EPS &&

multiply(p[i],q,p[i+1])*multiply(p[i],q2,p[i+1])<-EPS) count++;

return count&1;

}

bool cmp(Point p,Point q){

if (p.x

if (p.x>q.x) return false;

if (p.y

return false;

}

//线段与多边形相交定义为至少有1个点在多边形内，返回true; bool LINESEG_intersect_Polygon(LINESEG LS1,int n) {

LINESEG LS2;

Point mid;

int i;

int j,top;

Point stack[M];

p[n]=p[000];

for(i=0;i

LS2.st=p[i]; LS2.ed=p[i+1];

if (intersect_A(LS1,LS2)) return true;

}

top=0;

stack[top++]=LS1.st;

stack[top++]=LS1.ed;

for(i=0;i

if (online(LS1,p[i])) stack[top++]=p[i];

sort(stack,stack+top,cmp);

stack[top]=stack[0];

for(j=0;j

mid.x=(stack[j].x+stack[j+1].x)/2;

mid.y=(stack[j].y+stack[j+1].y)/2;

if (Inside(mid,n)) return true;

}

return false;

}

const int M=15000;

struct Line

{

double a,b,c;

}L[M];

struct Point

{

double x,y;

} p[M],pp[M];

//由两点化一般式

Line LineChange(Point P1,Point P2)

{

Line K;

K.a=P2.y-P1.y;

//(y2-y1)*x+(x1-x2)*y-(x1-x2)*y1-x1*(y2-y1)=0; counterclock

K.b=P1.x-P2.x;

K.c=-P1.x*K.a-K.b*P1.y;

return k;

}

//求两直线交点

bool GetSegcross(Line K1,Line K2,Point &pp)

{

double det=K1.a*K2.b-K2.a*K1.b;

if (fabs(det)<1e-8) return 0;

pp.x=(K1.b*K2.c-K2.b*K1.c)/det;

pp.y=-(K1.a*K2.c-K2.a*K1.c)/det;

return 1;

}

//p在直线上返回0.00

double p_on_line(Line K,Point p)

{

return K.a*p.x+K.b*p.y+K.c;

}

//求两点的垂直平分线

void MidLine(Point P1,Point P2,Line &K)

{

K.a=2*(P1.x-P2.x);

K.b=2*(P1.y-P2.y);

K.c=(-K.a*(P1.x+P2.x)-K.b*(P1.y+P2.y))/2;

}

//求任意点对的最小距离，分治nlogn.

double nearest(Point *p, int n) //p should already be sorted by x

{

if(n==2) return dist(p[0], p[1]);

if(n==1) return INF;

int m=p[n/2].x,i,j;

double left,right,tmp,ret,t;

left=nearest(p,n/2); right=nearest(&p[n/2],n-n/2);

ret=tmp =left

for(i=n/2-1; i>=0 && p[i].x>m-tmp_nearest;i--)

for(j=n/2; j

t=dist(p[i],p[j]);

if(t

}

return ret;

}

Euler的任意四面体体积公式（已知边长求体积）

已知4点坐标求体积（其中四个点的坐标分别为（x1,y1,z1）,（x2,y2,z2）,（x3,y3,z3）,

（x4,y4,z4））

1111

1234

(1/6)*.

1234

x x x x

y y y y

z z z z

222222

36.

p q n p r m

p q n q r l

V q

p r m q r l

+-+-