计算3D网格边界框的对角线长度的简单方法是什么?

问题描述 投票:0回答:1

我想计算3D网格的边界框的对角线长度。使用C ++迭代顶点,并搜索X坐标的(min,max),Y坐标的(min,max)和Z坐标的(min,max)。但是,我不知道如何利用这些获得的最小值/最大值来计算边界框的对角线长度。有什么帮助吗?

3d geometry mesh
1个回答
0
投票

为简单起见,我们考虑将n 3D点(点云)的列表作为输入(而不是网格),足以用于多边形网格。

网格的“对角线”就是网格中两个最远点之间的直线。这可以通过简单的O(n^2)蛮力搜索(2个嵌套循环来记住最远的点)轻松地计算出来。还有一些更快的方法可以利用点的排序。这里是蛮力的例子:

line pointcloud::diagonal()
    {
    int i,j;
    line l,ll;
    l=line(vec3(0.0,0.0,0.0),vec3(0.0,0.0,0.0)); // empty line
    for (i=0;i<pnt.num-1;i++)                    // O(n^2) search through all point pairs
     for (j=i+1;j<pnt.num-1;j++)
        {
        ll=line(pnt.dat[i],pnt.dat[j]);          // prepare line
        if (l.l<ll.l) l=ll;                      // compare sizes and remember the longer one
        }
    return l;
    }

有关linepointcloud类实现的更多信息,请阅读下面的链接和OBB的源代码。

但是从评论中我感觉到您需要3D OBB(定向边界框),而不只是对角线。您现在拥有的只是AABB(与轴对齐的边界框),它不会给您网格对角线(除非它的幸运方向与AABB对角线匹配)。

当心AABB和OBB对角线与网格对角线不同!!

有很多方法可以使用特征向量,凸包等从蛮力(〜O(n^6))更快地计算OBB ...

我设法将2D OBB approximation移植到3D中。

想法是一样的。将最大距离存储在“所有”(m)可能的方向/角度(在2D中覆盖整个球面而不是圆形),以将数据从n减少到m。然后只搜索计算出的数据以获取最小边界体积(而不是2D中的面积)。

我将Cone to box collision用于测试并作为起点。

算法:

  1. 计算枢轴点p0

    它必须在OBB的内部。通常,AABB的中心或平均点就足够了。

  2. 计算每个可能方向的距离

    存在无限可能的方向,因此我们需要将其限制为mm越大,计算速度越慢,但精度更高。为了快速存储和获取这些值,我使用了cube_map

    它是一个二维纹理,覆盖单位立方体的表面(6个正方形),并通过方向矢量而不是纹理坐标进行寻址。

    我实现了2个函数,用于在纹理数据(存储为1D数组)中的indexdirection向量之间进行转换。有关更多信息,请参见示例中的cube_map ...

    dp在某个方向p0上的距离dir的计算如下:

    d = dot( p-p0 , dir )
    

    因此生成m个可能的方向,并针对点源列表中所有点的每个计算距离,记住最大的一个,然后将其存储到cube_map。这是O( m * n )

    这里存储一帧距离的示例(cube_map的内容):

    distances in directions

  3. 寻找最小边界体积

    仅生成某个坐标系(覆盖半球)的所有m旋转。您不需要覆盖整个领域,因为另一半只是否定...

    现在,通过获取两个方向上沿其3个轴的距离并计算形成的框的体积,并记住最小的框(轴,距离和体积)来计算每个体积。由于舍入和非线性问题,在cube_map中可能存在统一的数据,从而导致volume = 0(如果cube_map在开始时清零),因此请忽略此类体积。

    在此之后,您应该使用OBB代替。这是OBB的几个旋转位置的预览:

    OBB preview

    这有点跳跃,因为对于这种对称形状,有无限数量的有效OBB,并且在不同的旋转中,可以首先在搜索中找到不同的一个。]]

  4. 提高准确性

  5. 简单地搜索附近发现的几个旋转<< OBB

    并记住最小的旋转。这可以递归完成。但是,我实在太懒惰了,因为OBB结果的当前状态对我来说足够了。

    此处为C ++ / GL源(其余内容可在上面的链接中找到:

]//--------------------------------------------------------------------------- class pointcloud { public: // cfg List<vec3> pnt; pointcloud() {} pointcloud(pointcloud& a) { *this=a; } ~pointcloud() {} pointcloud* operator = (const pointcloud *a) { *this=*a; return this; } //pointcloud* operator = (const pointcloud &a) { ...copy... return this; } void reset(){ pnt.num=0; } void add(vec3 p){ pnt.add(p); } void add(point p){ pnt.add(p.p0); } void compute(){}; void draw() { glBegin(GL_POINTS); for (int i=0;i<pnt.num;i++) glVertex3fv(pnt.dat[i].dat); glEnd(); } }; //--------------------------------------------------------------------------- template<class T,int N> class cube_map { public: int n,nn,sz; float fn2; T map[6*N*N]; cube_map() { n=N; nn=N*N; sz=6*nn; fn2=0.5*float(n); } cube_map(cube_map& a) { *this=a; } ~cube_map() {} cube_map* operator = (const cube_map *a) { *this=*a; return this; } //cube_map* operator = (const cube_map &a) { ...copy... return this; } vec3 ix2dir(int ix) { float x,y; vec3 dir=vec3(0.0,0.0,0.0); if ((ix<0)||(ix>=sz)) return dir; x=ix%n; ix/=n; x/=fn2; x--; y=ix%n; ix/=n; y/=fn2; y--; if (ix==0){ dir.y=x; dir.z=y; dir.x=-1.0; } if (ix==1){ dir.y=x; dir.z=y; dir.x=+1.0; } if (ix==2){ dir.x=x; dir.z=y; dir.y=-1.0; } if (ix==3){ dir.x=x; dir.z=y; dir.y=+1.0; } if (ix==4){ dir.x=x; dir.y=y; dir.z=-1.0; } if (ix==5){ dir.x=x; dir.y=y; dir.z=+1.0; } return normalize(dir); } int dir2ix(vec3 dir) { int ix=0,x=0,y=0; float a=0.0,b; b=fabs(dir.x); if (a<b){ a=b; if (dir.x<0) ix=0; else ix=1; } b=fabs(dir.y); if (a<b){ a=b; if (dir.y<0) ix=2; else ix=3; } b=fabs(dir.z); if (a<b){ a=b; if (dir.z<0) ix=4; else ix=5; } dir/=a; dir+=vec3(1.0,1.0,1.0); dir*=fn2; if (ix==0){ x=dir.y; y=dir.z; } if (ix==1){ x=dir.y; y=dir.z; } if (ix==2){ x=dir.x; y=dir.z; } if (ix==3){ x=dir.x; y=dir.z; } if (ix==4){ x=dir.x; y=dir.y; } if (ix==5){ x=dir.x; y=dir.y; } ix=(ix*nn)+(y*n)+(x); if ((ix<0)||(ix>=sz)) ix=0; return ix; } void set(vec3 dir,T &a){ map[dir2ix(dir)]=a; } T get(vec3 dir ){ return map[dir2ix(dir)]; } void clear(T &a){ for (int i=0;i<sz;i++) map[i]=a; } }; //--------------------------------------------------------------------------- class OBB // Oriented Bounding Box { public: // computed vec3 p0; // center vec3 u,v,w; // basis half vectors (p0 origin) OBB() {} OBB(OBB& a) { *this=a; } ~OBB() {} OBB* operator = (const OBB *a) { *this=*a; return this; } //OBB* operator = (const OBB &a) { ...copy... return this; } void compute(pointcloud &pcl) { const int N=24; int i,j,k,na=6*N,nb=2*N; cube_map<float,N> map; mat4 m,ma; vec3 o,p,q,pp0; int a,b; float da,db,d,dd,e,ee,V,VV; p0=vec3(0.0,0.0,0.0); u=vec3(0.0,0.0,0.0); v=vec3(0.0,0.0,0.0); w=vec3(0.0,0.0,0.0); if (pcl.pnt.num<=0) return; // init constants and stuff da=2.0*M_PI/float(na ); db= M_PI/float(nb-1); // compute avg point for (j=0;j<pcl.pnt.num;j++) p0+=pcl.pnt.dat[j]; p0/=pcl.pnt.num; // [compute perpendicular distances] // fill whole surface of cubemap for (map.clear(0.0),i=0;i<map.sz;i++) { // cube map index to 3D direction p=map.ix2dir(i); // compute max distance from p0 in direction p d=dot(pcl.pnt.dat[0]-p0,p); for (j=1;j<pcl.pnt.num;j++) { dd=dot(pcl.pnt.dat[j]-p0,p); if (d<dd) d=dd; } // store it in cube map for latter map.map[i]=d; } // [pick the smallest volume OBB combination] V=1e300; pp0=p0; // try half of "all" rotations (the other one is just negation) ma=mat4 // unit matrix -> unrotated coordinate system ( 1.0,0.0,0.0,0.0, 0.0,1.0,0.0,0.0, 0.0,0.0,1.0,0.0, 0.0,0.0,0.0,1.0 ); for ( a=0;a<na;a+=2,ma=lrotz(ma,da)) for (m=lroty(ma,float(-0.5*M_PI)),b=0;b<nb;b++,m=lroty(m,db)) { // get OBB per orientation of m p.x=map.get(-m[0].xyz); q.x=map.get(+m[0].xyz); p.y=map.get(-m[1].xyz); q.y=map.get(+m[1].xyz); p.z=map.get(-m[2].xyz); q.z=map.get(+m[2].xyz); o=p+q; VV=fabs(o.x*o.y*o.z); if ((V>VV)&&(VV>1e-6)) { V=VV; u=m[0].xyz; v=m[1].xyz; w=m[2].xyz; o*=0.5; pp0=p0+(u*(o.x-p.x))+(v*(o.y-p.y))+(w*(o.z-p.z)); u*=o.x; v*=o.y; w*=o.z; } } p0=pp0; } void draw() { const vec3 p[8]= { p0-u-v-w, p0+u-v-w, p0+u+v-w, p0-u+v-w, p0-u-v+w, p0+u-v+w, p0+u+v+w, p0-u+v+w, }; const int ix[24]= { 0,1,1,2,2,3,3,0, 4,5,5,6,6,7,7,4, 0,4,1,5,2,6,3,7, }; glBegin(GL_LINES); for (int i=0;i<24;i++) glVertex3fv(p[ix[i]].dat); glEnd(); } }; //---------------------------------------------------------------------------
希望我没有忘记复制某些东西……我想使代码尽可能地简单,因此它不是很优化,还有很多改进的余地。使用的数学基于

GLSL

,因此您可以使用GLM。我使用了自己的库,以便在需要时可以在上面的链接中找到vec(但需要将其作为〜220KByte的代码生成),但是它与GLSL和GLM完全匹配,因此您可以使用cn。但是mat4使用的某些功能在GLM中不存在,因此以防万一:
template <class T> class _mat4 { public: _vec4<T> col[4]; // columns!!! _mat4(T a00,T a01,T a02,T a03,T a04,T a05,T a06,T a07,T a08,T a09,T a10,T a11,T a12,T a13,T a14,T a15) { col[0]=vec4(a00,a01,a02,a03); // x axis col[1]=vec4(a04,a05,a06,a07); // y axis col[2]=vec4(a08,a09,a10,a11); // z axis col[3]=vec4(a12,a13,a14,a15); // origin } _mat4() { col[0]=vec4(1,0,0,0); col[1]=vec4(0,1,0,0); col[2]=vec4(0,0,1,0); col[3]=vec4(0,0,0,1); } _mat4(const _mat4& a) { *this=a; } ~_mat4() {} // operators (matrix math) _mat4* operator = (const _mat4 &a) { for (int i=0;i<4;i++) col[i]=a.col[i]; return this; } // =a[][] _vec4<T>& operator [](const int i){ return col[i]; } // a[i] _mat4<T> operator * (_mat4<T>&m) // =a[][]*m[][] { _mat4<T> q; int i,j,k; for (i=0;i<4;i++) for (j=0;j<4;j++) for (q.col[i][j]=0,k=0;k<4;k++) q.col[i][j]+=col[k][j]*m.col[i][k]; return q; } _mat4<T> operator * (_vec4<T>&v) // =a[][]*v[] { _vec4<T> q; int i,j; for (i=0;i<4;i++) for (q.dat[i]=0,j=0;j<4;j++) q.dat[i]+=col[i][j]*v.dar[j]; return q; } _mat4<T> operator * (T &c) // =a[][]*c { _mat4<T> q; int i,j; for (i=0;i<4;i++) for (j=0;j<4;j++) q.dat[i]=col[i][j]*c; return q; } _mat4<T> operator / (T &c) // =a[][]/c { _mat4<T> q; int i,j; for (i=0;i<4;i++) for (j=0;j<4;j++) q.dat[i]=divide(col[i][j],c); return q; } _mat4<T> operator *=(_mat4<T>&m){ this[0]=this[0]*m; return *this; }; _mat4<T> operator *=(_vec4<T>&v){ this[0]=this[0]*v; return *this; }; _mat4<T> operator *=(const T &c){ this[0]=this[0]*c; return *this; }; _mat4<T> operator /=(const T &c){ this[0]=this[0]/c; return *this; }; // members void get(T *a) { int i,j,k; for (k=0,i=0;i<4;i++) for (j=0;j<4;j++) a[k]=col[i].dat[j]; } }; //--------------------------------------------------------------------------- template <class T> _mat4<T> transpose(const _mat4<T> &m) { _mat4<T> q; int i,j; for (i=0;i<4;i++) for (j=0;j<4;j++) q.col[i][j]=m.col[j][i]; return q; } //--------------------------------------------------------------------------- template <class T> _mat4<T> inverse(const _mat4<T> &m) { T p[3] _mat4<T> q; T x,y,z; int i,j; // transpose rotation for (i=0;i<3;i++) for (j=0;j<3;j++) q.col[i][j]=m.col[j][i]; // copy projection for (i=0;i<4;i++) q.col[i][3]=m.col[i][3]; // convert origin: new_pos = - new_rotation_matrix * old_pos for (i=0;i<3;i++) for (p[i]=0,j=0;j<3;j++) p[i]+=q.col[j][i]*col[3][j]; for (i=0;i<3;i++) q.col[3][i]=-p[i]; return q; } //--------------------------------------------------------------------------- template <class T> _mat4<T> lrotx(_mat4<T> &m,T ang) { T c=cos(ang),s=sin(ang); mat4 r=mat4( 1, 0, 0, 0, 0, c, s, 0, 0,-s, c, 0, 0, 0, 0, 1); r=m*r; return r; }; //--------------------------------------------------------------------------- template <class T> _mat4<T> lroty(_mat4<T> &m,T ang) { T c=cos(ang),s=sin(ang); mat4 r=mat4( c, 0,-s, 0, 0, 1, 0, 0, s, 0, c, 0, 0, 0, 0, 1); r=m*r; return r; }; //--------------------------------------------------------------------------- template <class T> _mat4<T> lrotz(_mat4<T> &m,T ang) { T c=cos(ang),s=sin(ang); mat4 r=mat4( c, s, 0, 0, -s, c, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1); r=m*r; return r; }; //--------------------------------------------------------------------------- template <class T> _mat4<T> grotx(_mat4<T> &m,T ang){ return inverse(lrotx(inverse(m),ang)); }; template <class T> _mat4<T> groty(_mat4<T> &m,T ang){ return inverse(lroty(inverse(m),ang)); }; template <class T> _mat4<T> grotz(_mat4<T> &m,T ang){ return inverse(lrotz(inverse(m),ang)); }; //--------------------------------------------------------------------------- typedef _mat4<float > mat4; typedef _mat4<double> dmat4; typedef _mat4<bool > bmat4; typedef _mat4<int > imat4; typedef _mat4<DWORD > umat4; //--------------------------------------------------------------------------- mat4 GLSL_math_test4x4; //---------------------------------------------------------------------------
为了理解它或自己写,我建议看看:

List<double> xxx;double xxx[];相同xxx.add(5);5添加到列表的末尾xxx[7]访问数组元素(安全)xxx.dat[7]访问数组元素(不安全但快速的直接访问)xxx.num是数组的实际使用大小xxx.reset()清除数组并设置xxx.num=0xxx.allocate(100)100个项目预分配空间

现在将结果放入OBB

仅由其中心p0和一半向量u,v,w描述的方框。因此,要获取点云PCL

OBB

,只需计算:
OBB obb; pointcloud PCL; PCL.reset(); PCL.add(...); // here feed points into PCL obb.compute(PCL);
仅此而已。
© www.soinside.com 2019 - 2024. All rights reserved.