将希尔伯特值映射到 3D 点

Question

我有一组希尔伯特值（从Hilbert curve 开始到给定点的长度）。

将这些值转换为 3D 点的最佳方法是什么？原始希尔伯特曲线不是 3D 的，所以我想我必须自己选择我需要的希尔伯特曲线等级。我确实有总曲线长度（即集合中的最大值）。

也许有一个现有的实现？一些允许我使用希尔伯特曲线/值的库？语言并不重要。

Answer 1

不是关于 3D 转换的答案，但这里有一个很好的算法和关于希尔伯特值的讨论Two-dimensional spatial hashing with space-filling curves

来自MIT

4 algorithms for the n-dimensional Hilbert Space-Filling Curve

* A. R. Butz, "Alternative Algorithm for Hilbert's Space-Filling Curve",
  IEEE Trans. Comp., April, 1971, pp 424-426. [Butz 1971]

* S. W. Thomas, "hilbert.c" in the Utah Raster Toolkit circa 1993,
  http://web.mit.edu/afs/athena/contrib/urt/src/urt3.1/urt-3.1b.tar.gz

* D. Moore, Fast Hilbert Curves in C, without Recursion

* J.K.Lawder, Calculation of Mappings Between One and n-dimensional Values Using the Hilbert Space-filling Curve, [JL1_00]

Answer 2

如果我的问题是正确的，那么您与 Hilbert 3D 曲线的起点有一些曲线长度距离l，并且希望获得与该点相对应的坐标。

如果您将整个 3D 希尔伯特曲线（覆盖单位立方体）预先生成为 折线，则所有排序的点在前一个点和下一个点之间的距离相同。因此，您可以使用piecewise linear interpolation 计算您的观点。

这就是我在 C++ 中生成和渲染 2D/3D 希尔伯特曲线的方式：

//---------------------------------------------------------------------------
#ifndef _Hilbert_vector_h
#define _Hilbert_vector_h
//---------------------------------------------------------------------------
#include "list.h"
//---------------------------------------------------------------------------
void Hilbert2D(List<double> &pnt,double x,double y,double z,double a,int n)
    {
    int i,j,m;
    double x0,y0,x1,y1,q;
    for (m=4*3,i=1,j=2;j<=n;j++,i+=i+1) m*=4; a/=i; // m = needed size of pnt[]
    pnt.num=0;
    // init generator
          pnt.add(x); pnt.add(y); pnt.add(z);
    y+=a; pnt.add(x); pnt.add(y); pnt.add(z);
    x+=a; pnt.add(x); pnt.add(y); pnt.add(z);
    y-=a; pnt.add(x); pnt.add(y); pnt.add(z);
    x0=x-0.5*a; // center of generator
    y0=y+0.5*a;
    // iterative subdivision
    for (j=2;j<=n;j++)
        {
        // mirror/rotate 2 qudrants
        x1=x0; y1=y0; m=pnt.num;
        for (i=m;i>=3;)
            {
            i--; z=pnt.dat[i]   ;
            i--; y=pnt.dat[i]-y0;
            i--; x=pnt.dat[i]-x0;
            q=x; x=+y; y=-q;    // z+
            pnt.dat[i+0]=(x1+x);
            pnt.dat[i+1]=(y1-y);
            pnt.dat[i+2]=(   z);
            }
        for (y1+=2.0*a,i=m;i>=3;)
            {
            i--; z=pnt.dat[i]   ;
            i--; y=pnt.dat[i]-y0;
            i--; x=pnt.dat[i]-x0;
            q=x; x=-y; y=+q;    // z-
            pnt.add(x1+x);
            pnt.add(y1+y);
            pnt.add(   z);
            }
        // mirror the rest
        x0+=a; y0+=a; m=pnt.num;
        for (i=m;i>=3;)
            {
            i--; z=pnt.dat[i]   ;
            i--; y=pnt.dat[i]-y0;
            i--; x=pnt.dat[i]-x0;
            pnt.add(x0-x);
            pnt.add(y0+y);
            pnt.add(   z);
            }
        a*=2.0;
        }
/*
        // rotations
        q=x; x=+y; y=-q;    // z+
        q=x; x=-y; y=+q;    // z-
*/
    }
//---------------------------------------------------------------------------
void Hilbert3D(List<double> &pnt,double x,double y,double z,double a,int n)
    {
    int i,j,m;
    double x0,y0,z0,x1,y1,z1,q;
    for (m=8*3,i=1,j=2;j<=n;j++,i+=i+1) m*=8; a/=i; // m = needed size of pnt[]
    pnt.num=0;
    // init generator
          pnt.add(x); pnt.add(y); pnt.add(z);
    z-=a; pnt.add(x); pnt.add(y); pnt.add(z);
    x+=a; pnt.add(x); pnt.add(y); pnt.add(z);
    z+=a; pnt.add(x); pnt.add(y); pnt.add(z);
    y+=a; pnt.add(x); pnt.add(y); pnt.add(z);
    z-=a; pnt.add(x); pnt.add(y); pnt.add(z);
    x-=a; pnt.add(x); pnt.add(y); pnt.add(z);
    z+=a; pnt.add(x); pnt.add(y); pnt.add(z);
    x0=x+0.5*a; // center of generator
    y0=y-0.5*a;
    z0=z-0.5*a;
    // iterative subdivision
    for (j=2;j<=n;j++)
        {
        // mirror/rotate qudrants
        x1=x0; y1=y0; z1=z0; m=pnt.num;
        for (i=m;i>=3;)
            {
            i--; z=pnt.dat[i]-z0;
            i--; y=pnt.dat[i]-y0;
            i--; x=pnt.dat[i]-x0;
            q=y; y=-z; z=+q;    // x-
            pnt.dat[i+0]=(x1+x);
            pnt.dat[i+1]=(y1+y);
            pnt.dat[i+2]=(z1-z);
            }
        for (z1-=2.0*a,i=m;i>=3;)
            {
            i--; z=pnt.dat[i]-z0;
            i--; y=pnt.dat[i]-y0;
            i--; x=pnt.dat[i]-x0;
            q=z; z=+x; x=-q;    // y+
            q=y; y=+z; z=-q;    // x+
            pnt.add(x1-x);
            pnt.add(y1+y);
            pnt.add(z1+z);
            }
        for (x1+=2.0*a,i=m;i>=3;)
            {
            i--; z=pnt.dat[i]-z0;
            i--; y=pnt.dat[i]-y0;
            i--; x=pnt.dat[i]-x0;
            q=y; y=+z; z=-q;    // x+
            pnt.add(x1+x);
            pnt.add(y1+y);
            pnt.add(z1+z);
            }
        for (z1+=2.0*a,i=m;i>=3;)
            {
            i--; z=pnt.dat[i]-z0;
            i--; y=pnt.dat[i]-y0;
            i--; x=pnt.dat[i]-x0;
            q=z; z=+x; x=-q;    // y+
            pnt.add(x1-x);
            pnt.add(y1-y);
            pnt.add(z1+z);
            }
        // mirror octants
        x0+=a; y0+=a; z0-=a; m=pnt.num;
        for (i=m;i>=3;)
            {
            i--; z=pnt.dat[i]-z0;
            i--; y=pnt.dat[i]-y0;
            i--; x=pnt.dat[i]-x0;
            pnt.add(x0+x);
            pnt.add(y0-y);
            pnt.add(z0+z);
            }
        a*=2.0;
        }
/*
        // rotations
        q=z; z=+x; x=-q;    // y+
        q=z; z=-x; x=+q;    // y-
        q=y; y=+z; z=-q;    // x+
        q=y; y=-z; z=+q;    // x-
        q=x; x=+y; y=-q;    // z+
        q=x; x=-y; y=+q;    // z-
*/
    }
//---------------------------------------------------------------------------
void pnt_draw2(List<double> &pnt)   // piecewise linear
    {
    int i;
    glBegin(GL_LINE_STRIP);
    for (i=0;i<pnt.num;i+=3) glVertex3dv(pnt.dat+i);
    glEnd();
    }
//---------------------------------------------------------------------------
void pnt_draw4(List<double> &pnt)   // piecewise cubic
    {
    int i,j;
    double  d1,d2,t,tt,ttt,*p0,*p1,*p2,*p3,a0[3],a1[3],a2[3],a3[3],p[3];
    glBegin(GL_LINE_STRIP);
    for (i=-3;i<pnt.num;i+=3)
        {
        j=i-3; if (j>pnt.num-3) j=pnt.num-3; if (j<0) j=0; p0=pnt.dat+j;
        j=i  ; if (j>pnt.num-3) j=pnt.num-3; if (j<0) j=0; p1=pnt.dat+j;
        j=i+3; if (j>pnt.num-3) j=pnt.num-3; if (j<0) j=0; p2=pnt.dat+j;
        j=i+6; if (j>pnt.num-3) j=pnt.num-3; if (j<0) j=0; p3=pnt.dat+j;
        for (j=0;j<3;j++)
            {
            d1=0.5*(p2[j]-p0[j]);
            d2=0.5*(p3[j]-p1[j]);
            a0[j]=p1[j];
            a1[j]=d1;
            a2[j]=(3.0*(p2[j]-p1[j]))-(2.0*d1)-d2;
            a3[j]=d1+d2+(2.0*(-p2[j]+p1[j]));
            }
        for (t=0.0;t<=1.0;t+=0.1)   // single curve patch/segment
            {
            tt=t*t;
            ttt=tt*t;
            for (j=0;j<3;j++) p[j]=a0[j]+(a1[j]*t)+(a2[j]*tt)+(a3[j]*ttt);
            glVertex3dv(p);
            }
        }
    glEnd();
    }
//---------------------------------------------------------------------------
#endif
//---------------------------------------------------------------------------

我使用了我的动态列表模板，所以：

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

但是您可以使用动态甚至静态 1D 数组代替，因为 希尔伯特曲线 的点数很容易计算（每个希尔伯特函数开头的m )。

使用很简单，只需这样做：

List<double> pnt;
Hilbert3D(pnt,-0.8,-0.8,+0.8,1.6,n); 

其中n 是迭代次数，pnt 是每个点的 (x,y,z) 坐标的线性列表（每个点 3 个数字）。开始位置和初始大小设置为覆盖以(0,0,0) 为中心的立方体，一半大小为0.8 <-0.8,+0.8>。

现在只需计算点之间的单位长度、最靠近左边的希尔伯特曲线点的索引以及线性插值的参数（到它的距离）。这里以C++为例：

if (pnt.num>=6)
    {
    int i;
    double x,y,z,t,l,dl;
    dl=sqrt(                                                    // base distance between points
            ((pnt[0]-pnt[3])*(pnt[0]-pnt[3]))
           +((pnt[1]-pnt[4])*(pnt[1]-pnt[4]))
           +((pnt[2]-pnt[5])*(pnt[2]-pnt[5]))
            );
    l=double(Form1->sb_t->Position)/double(Form1->sb_t->Max);   // <0,1>
    l*=dl*double((pnt.num/3)-1);                                // <0,Hilbert_curve_lenght>
    i=floor(l/dl); t=(l-(double(i)*dl))/dl; i*=3;               // index in pnt[] and single line segment paramerer
    x=pnt[i+0]+(pnt[i+3]-pnt[i+0])*t;                           // linear interpolation
    y=pnt[i+1]+(pnt[i+4]-pnt[i+1])*t;
    z=pnt[i+2]+(pnt[i+5]-pnt[i+2])*t;
    glColor3f(0.0,1.0,0.0); t=0.05;                             // render of marker
    glBegin(GL_LINES);
    glVertex3d(x-t,y-t,z); glVertex3d(x+t,y+t,z);
    glVertex3d(x+t,y-t,z); glVertex3d(x-t,y+t,z);
    glVertex3d(x,y-t,z-t); glVertex3d(x,y+t,z+t);
    glVertex3d(x,y-t,z+t); glVertex3d(x,y+t,z-t);
    glVertex3d(x-t,y,z-t); glVertex3d(x+t,y,z+t);
    glVertex3d(x+t,y,z-t); glVertex3d(x-t,y,z+t);
    glEnd();
    }

2D 预览：

3D 预览：