If I get your question right you have some curve-length distance l
from start point of Hilbert 3D curve and want to get coordinates corresponding to such point.
if you pre-generate the whole 3D Hilbert curve (covering unit cube) as a Polyline then all the sequenced points are at same distance between previous and next point. So you can compute your point using piecewise linear interpolation.
This is how I generate and render 2D/3D Hilbert curve in C++:
//---------------------------------------------------------------------------
#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
//---------------------------------------------------------------------------
I used mine dynamic list template so:
List<double> xxx;
is the same as double xxx[];
xxx.add(5);
adds 5
to end of the list
xxx[7]
access array element (safe)
xxx.dat[7]
access array element (unsafe but fast direct access)
xxx.num
is the actual used size of the array
xxx.reset()
clears the array and set xxx.num=0
xxx.allocate(100)
preallocate space for 100
items
But you can use dynamic or even static 1D array instead as the number of points of Hilbert curve is easily computable (m
at the start of each Hilbert function).
Usage is simple just do something like this:
List<double> pnt;
Hilbert3D(pnt,-0.8,-0.8,+0.8,1.6,n);
Where n
is number of iterations and pnt
is linear list of (x,y,z)
coordinates for each point (3 numbers per point). the start position and initial size is set to cover cube centered around (0,0,0)
with half size 0.8
<-0.8,+0.8>
.
Now just compute unit length between points, index of closest Hilbert curve point to the left and parameter (distance to it) from that just linearly interpolate. Here C++ example:
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 preview:
3D preview: