Why so serious? --[NKU]schindlerlee

2010-06-16 12:04:40 pku3528 三维凸包


三维凸包的求解如果有四点共面,那么就比较繁琐了。
增量法,每次添加一个点,然后更新凸包。
关于增量法的图示看这里:
http://www.cse.unsw.edu.au/~lambert/java/3d/hull.html

现在的主要问题就是,一直所有面的法向量,如何更新凸包。
法向量点乘一个向量,可以得到这个向量在法向量上的偏移,如果这个偏移是正的,说明这个向量相关的点在面的上方,负的是下方,所以只要找到一条边关联的两个面的点积是一正一负即可。

注意这里的所有面都是由三个点组成,且是有序的,每个边只被记录两次,一次顺时针,一次逆时针。
这里不能使用int,数值太大会越界。

 1 
 2 const int N = 1024;
 3 struct F {
 4     int v[4];
 5     F(){}
 6     F(int a,int b,int c){ v[0= a, v[1= b, v[2= c, v[3= v[0];}
 7 };
 8 
 9 struct point_t {
10     double x, y, z;
11     point_t (){}
12     point_t (double a, double b, double c){x = a, y = b, z = c;}
13 };
14 point_t operator +(point_t a, point_t b) { return point_t (a.x + b.x, a.y + b.y, a.z + b.z);}
15 point_t operator -(point_t a, point_t b) { return point_t (a.x - b.x, a.y - b.y, a.z - b.z);}
16 
17 double dot_mul(point_t a, point_t b) { return a.x * b.x + a.y * b.y + a.z * b.z;}
18 point_t cross_mul(point_t a,point_t b)
19 {
20     return point_t (a.y*b.z - b.y*a.z ,
21                     a.z*b.x - b.z*a.x ,
22                     a.x*b.y - b.x*a.y);
23 }
24 
25 double area(point_t a, point_t b)
26 {
27   point_t v = cross_mul(a, b);
28   return sqrt(0.0 + v.x * v.x + v.y * v.y + v.z * v.z) / 2;
29 }
30 
31 point_t p[N];
32 int vis[N][N], n;
33 
34 int main()
35 {
36   int i, j, k;
37   scanf("%d"&n);
38   for (i = 0;i < n;i++) {
39       scanf("%lf %lf %lf"&p[i].x, &p[i].y, &p[i].z);
40   }
41   vector<F> cur;
42   cur.pb(F(012)), cur.pb(F(210));
43   for (i = 3;i < n;i++) {
44       vector<F> next;
45       for (j = 0;j < cur.size();j++) {
46           F f = cur[j];
47           point_t v = cross_mul(p[f.v[1]] - p[f.v[0]],
48                                 p[f.v[2]] - p[f.v[1]]);
49           int val = dot_mul(v, p[i] - p[f.v[1]]) < 0 ? -1 : 1;
50           if (val < 0) {
51               next.pb(f);
52           }
53           for (k = 0;k < 3;k++) {
54               if (vis[f.v[k+1]][f.v[k]] == 0) {
55                   vis[f.v[k]][f.v[k+1]] = val;
56               }else { //http://www.cppblog.com/schindlerlee
57                   if (vis[f.v[k+1]][f.v[k]] != val) {
58                       if (val > 0) {
59                           next.pb(F(f.v[k], f.v[k+1], i));
60                       }else {
61                           next.pb(F(f.v[k+1], f.v[k], i));
62                       }
63                   }
64                   vis[f.v[k+1]][f.v[k]] = 0;
65               }
66           }
67       }
68       cur = next;
69   }
70   double ans = 0;
71   for (i = 0;i < cur.size();i++) {
72       F f = cur[i];
73       ans += area(p[f.v[0]] - p[f.v[1]], p[f.v[2]] - p[f.v[1]]);
74   }
75   printf("%.3f\n", ans);
76   return 0;
77 }


posted on 2010-06-16 12:02 schindlerlee 阅读(2006) 评论(2)  编辑 收藏 引用 所属分类: 解题报告

Feedback

# re: 2010-06-16 12:04:40 pku3528 三维凸包 2010-06-18 13:37 AngelClover

你这是啥复杂度  回复  更多评论   

# re: 2010-06-16 12:04:40 pku3528 三维凸包 2010-06-18 16:03 schindlerlee

@AngelClover
O(n^2)  回复  更多评论   


只有注册用户登录后才能发表评论。
网站导航: 博客园   IT新闻   BlogJava   知识库   博问   管理