poj3528 三维凸包面积
2012-10-06 19:57
411 查看
#include<iostream> #include<cstring> #include<cmath> #include<algorithm> using namespace std; #define N 600 #define EPS 1e-8 struct tpoint { double x,y,z; tpoint(){} tpoint(double _x,double _y,double _z):x(_x),y(_y),z(_z){} tpoint operator-(const tpoint p) { return tpoint(x-p.x,y-p.y,z-p.z); } tpoint operator *(const tpoint p) { return tpoint(y*p.z-z*p.y,z*p.x-x*p.z,x*p.y-y*p.x);//叉积 } double operator ^(const tpoint p) { return x*p.x+y*p.y+z*p.z;//点积 } }; struct fac { int a,b,c;//凸包一个面上的三个点的编号 bool ok;//该面是否是最终凸包的面 }; struct t3dhull { int n;//初始点数 tpoint ply ;//初始点 int trianglecnt;//土包上三角形数 fac tri[6*N];//凸包三角形创建的面 int vis ;//点i到j是属于哪个面的 //两点长度 double dist(tpoint a) { return sqrt(a.x*a.x+a.y*a.y+a.z*a.z); } //三角形面积 *2 double area(tpoint a,tpoint b,tpoint c) { return dist((b-a)*(c-a)); } //四面体有向体积*6 double volume(tpoint a,tpoint b,tpoint c,tpoint d) { return (b-a)*(c-a)^(d-a); } //正:点在面的铜像 double ptoplane(tpoint &p,fac &f) { tpoint m=ply[f.b]-ply[f.a],n=ply[f.c]-ply[f.a],t=p-ply[f.a]; return (m*n)^t; } void deal(int p,int a,int b) { int f=vis[a][b]; fac add; if(tri[f].ok) { if((ptoplane(ply[p],tri[f]))>EPS) dfs(p,f); else { add.a=b; add.b=a; add.c=p; add.ok=1; vis[p][b]=vis[a][p]=vis[b][a]=trianglecnt; tri[trianglecnt++]=add; } } } void dfs(int p,int cnt) { tri[cnt].ok=0; deal(p,tri[cnt].b,tri[cnt].a); deal(p,tri[cnt].c,tri[cnt].b); deal(p,tri[cnt].a,tri[cnt].c); } void construct() { int i,j; trianglecnt=0; if(n<4)return; bool tmp=true; for(i=1;i<n;i++) { if((dist(ply[0]-ply[i]))>EPS) { swap(ply[1],ply[i]); tmp=false; break; } } if(tmp)return; tmp=true; for(i=2;i<n;i++) { if((dist((ply[0]-ply[1])*(ply[1]-ply[i])))>EPS) { swap(ply[2],ply[i]); tmp=false; break; } } if(tmp)return; tmp=true; for(i=3;i<n;i++) { if(fabs((ply[0]-ply[1])*(ply[1]-ply[2])^(ply[0]-ply[i]))>EPS) { swap(ply[3],ply[i]); tmp=false; break; } } if(tmp)return; fac add; for(i=0;i<4;i++) { add.a=(i+1)%4; add.b=(i+2)%4; add.c=(i+3)%4; add.ok=1; if((ptoplane(ply[i],add))>0) swap(add.b,add.c); vis[add.a][add.b]=vis[add.b][add.c]=vis[add.c][add.a]=trianglecnt; tri[trianglecnt++]=add; } for(i=4;i<n;i++) { for(j=0;j<trianglecnt;j++) { if(tri[j].ok && (ptoplane(ply[i],tri[j]))>EPS) { dfs(i,j); break; } } } int cnt=trianglecnt; trianglecnt=0; for(i=0;i<cnt;i++) { if(tri[i].ok) tri[trianglecnt++]=tri[i]; } } double area() { double ret=0; for(int i=0;i<trianglecnt;i++) ret+=area(ply[tri[i].a],ply[tri[i].b],ply[tri[i].c]); return ret/2.0; } }hull; int main() { while(scanf("%d",&hull.n)!=EOF) { int i; for(i=0;i<hull.n;i++) scanf("%lf%lf%lf",&hull.ply[i].x,&hull.ply[i].y,&hull.ply[i].z); hull.construct(); printf("%.3lf\n",hull.area()); } return 0; }
相关文章推荐
- 三维凸包求其表面积(POJ3528)
- 三维凸包求其表面积(POJ3528)
- POJ3528 (三维凸包模板)
- POJ3528最大化最小值+三维凸包
- UVA 11769 All Souls Night 的三维凸包要求的表面面积
- poj3528,hdu3662三维凸包模板题
- poj3528--Ultimate Weapon(三维凸包)
- UVA 11769 All Souls Night 三维凸包求表面积
- poj3528 Ultimate Weapon 求三维凸包的表面积
- HDU1616计算几何-凸包+点在多边形内+多边形的面积
- hdu4449Building Design(三维凸包+平面旋转)
- 三维凸包求内部一点到表面的最近距离(HDU4266)
- 二维凸包模板(凸包重心,周长,面积,直径,最大三角形,最小环绕矩形)
- pku 3348 Cows(求凸包面积)
- BZOJ 1069 SCOI 2007 最大土地面积 凸包+旋转卡壳
- poj 2079 Triangle 凸包+旋转卡壳 求最大三角形面积
- Cows(凸包面积)
- uva109求凸包面积,判断点是不是在凸包内
- SWUST0249 (凸包面积)
- 三维模型凸包提取