将希尔伯特值映射到 3D 点
Posted
技术标签:
【中文标题】将希尔伯特值映射到 3D 点【英文标题】：Mapping Hilbert values to 3D points 【发布时间】：2010-10-04 16:39:29 【问题描述】：我有一组希尔伯特值（从Hilbert curve 开始到给定点的长度）。
将这些值转换为 3D 点的最佳方法是什么？原始希尔伯特曲线不是 3D 的，所以我想我必须自己选择我需要的希尔伯特曲线等级。我确实有总曲线长度（即集合中的最大值）。
也许有一个现有的实现？一些允许我使用希尔伯特曲线/值的库？语言并不重要。
【问题讨论】：
我最近遇到了Smooth Hilbert curves，成功后我决定也尝试 3D。你的问题很容易直接解决，见我的回答（全部用 C++ 编码） 【参考方案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]
【讨论】：
【参考方案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 预览：
【讨论】：
以上是关于将希尔伯特值映射到 3D 点的主要内容，如果未能解决你的问题，请参考以下文章
如何使用 std::map 将 bool 映射到 3d 点结构？