Hilbert curve

Hilbert curve 는 좌표계를 mapping하는 방법이다.

참고한 코드는 https://people.sc.fsu.edu/~jburkardt/py_src/hilbert_curve/hilbert_curve.html 이곳이며, [WIKIPEDIA]에도 같은 C코드가 있다. 위 사이트의 python코드에서 xy2d 함수의 xcopy와 ycopy 앞에 abs 는 typo이며 삭제해야 한다. Hilbert curve를 generate하여 visualization해 주는 사이트는 http://bit-player.org/extras/hilbert/hilbert-construction.html 이다.

 

Hilbert curve coordinate mapping은 총 3개의 함수로 구성되어 있다. 1) 1차원에서 2차원으로 mapping (d2xy), 2) 2차원에서 1차원으로 mapping (xy2d), 3) 1)과 2)함수의 종속 함수인 (rot)가 있다. 어떠한 사이즈의 Hilbert curve가 주어지건 사분면으로 나누어 진다. 가장 기본이 되는 틀은 m=2일 때인 2x2 이며, 코드에서 정의한 s=1이다. (0,0)에 해당하는 값은 x와 y의 Transpose가 일어나며, (1,0)에 해당하는 값은 180도 Rotation과 Transpose가 일어난다.

 

[ 2x2 Hilbert curve unit matrix (rx,ry) ]

1) d2xy ( m, d ) 함수는 1차원에서 2차원으로의 mapping이며, 1차원의 수, d,가 들어오면 1차적으로 d 값이 2x2 사이즈인 unit matrix에서 어느 곳에 해당하는지 먼저 찾게 된다. 방법은 다음과 같다.
1-1) 1차원에서 값인 d를 2로 나눈 몫은 짝수와 홀수가 나온다. 짝수는 왼쪽 그림의 0, 1에 해당하는 rx=0인 위치에 해당되고, 홀수는 rx=1인 2, 3에 위치하게 된다.

1-2) rx=0일때, d가 짝수면 (0,0), 홀수면 ry=1인 (0,1)가 된다. rx=1인 경우, d가 짝수면 ry=1을 갖도록하는 (1,1), 홀수이면 ry=0을 갖는 (1,0)에 대응된다. 이것으로 d값이 unit matrix에 어느 위치에 해당하는지 알았다. 위치가 정해졌으니 (x, y)가 어디에 대응하는지, 즉 (x,y) 변환이 없이 갈지, transpose를 해야할지, rotation 후에 transpose를 해야할지 판단이 가능하다.

1-3) 초기 (x,y)=(0,0)이며 s=1이다. 좌표의 transpose나 rotation 이후 (x',y')=(x+s*rx,y+s*ry)가 된다. s=1일때 값을 알았으니 원래 사이즈인 s를 찾아가 보자. s=2일때는 d를 4로 나눈 몫을 d라고 생각하고 1-2)과정을 진행하여 unit matrix에서 어느 위치에 해당하는지 확인하고 (x,y)를 변환해 준다.

1-4) 2*s가 실제 size n이 될 때까지 반복하면 된다.

 

2) xy2d ( m, x, y ) 함수는 1)과 반대로 d=0 초기값으로 설정하고, 실제 사이즈, nxn에서의 주어진 (x,y)에 해당하는 s^2=(n/2)^2값을 d에 더해주고 s=1일때까지 위치에서의 s제곱값을 계속 더해주면 (x,y)에 해당하는 d값을 구할 수 있다.

 

* rot ( s, x, y, rx, ry ) 실제 코드 코멘트에서 Reflect 는 Rotation이라고 표현하는게 맞고, 이것은 input (x, y)가 sxs사이즈 matrix에서 180도 돌린 위치로 R(x, y)를 바꾼다. 예를들어 s=4 이면 (x, y)가 속한 4x4 내에서 180도 돌린 위치가 된다. Flip은 Transpose이며 T(x, y)=(y,x)이다.

 

[ Example of finding the coordinate (x,y) of 51 ]

 

Ex) 51이 어떤 (x,y)를 가지는지 찾아보자.

1) s=1이고 (x,y)=(0,0)에서 시작한다.

  - t=51을 2로 나눈 몫을 2로 나눈 나머지는 1로, rx=1이다.

  - rx가 1 일때는 51은 홀수여서 ry=0을 가지게 된다. 만약 짝수이면 ry=1을 가지게 된다.

  - (rx,ry)=(1,0) 이기 때문에 R, T 변환이 일어나고 (x',y')=(x+s*rx,y+s*ry)가 된다. 따라서 (0,0) -> (1,0)이 된다.

2) s=2이고 (x,y)=(1,0)에서 시작하고, t=51//4=12이다.

  - t=12을 2로 나눈 몫을 2로 나눈 나머지는 0으로, rx=0이다.

  - rx가 0 일때는 12는 짝수여서 ry=0을 가지게 된다.

  - (rx,ry)=(0,0) 이기 때문에 T 변환이 일어나고 (x',y')=(0+2*0,1+2*0)=(0,1)이 된다.

3) s=4이고 (x,y)=(0,1)에서 시작하고, t=12//4=3이다.

  - t=3을 2로 나눈 몫을 2로 나눈 나머지는 1로, rx=1이다.

  - rx가 1 일때는 3은 홀수여서 ry=0을 가지게 된다.

  - (rx,ry)=(1,0) 이기 때문에, 4x4 사이즈(0:3,0:3)의 중심에서 180도로 (x,y) Rotation과 Transpose 변환이 일어나고 (x',y')=(2+4*1,3+4*0)=(6,3)이 된다.

 

CODE : Generate Hilbert curve

def d2xy ( m, d ): 
  n = 2 ** m 
  x = 0; y = 0; t = d; s = 1;

  while ( s < n ): 
    rx = ( ( t // 2 ) % 2 ) 
    if ( rx == 0 ): 
      ry = ( t % 2 ) 
    else: 
      ry = ( ( t ^ rx ) % 2 ) 
    x, y = rot ( s, x, y, rx, ry ) 
    x = x + s * rx 
    y = y + s * ry 
    t = ( t // 4 ) 

    s = s * 2 

  return x, y


def xy2d ( m, x, y ): 
  xcopy = x; ycopy = y;
  d = 0 
  n = 2 ** m 
  s = ( n // 2 )

  while ( 0 < s ): 
    if ( 0 <  ( xcopy & s ) ): 
      rx = 1 
    else: 
      rx = 0 

    if ( 0 < ( ycopy & s ) ): 
      ry = 1 
    else: 
      ry = 0

    d = d + s * s * ( ( 3 * rx ) ^ ry ) 
    xcopy, ycopy = rot ( s, xcopy, ycopy, rx, ry) 
    s = ( s // 2 ) 

  return d 


def rot ( s, x, y, rx, ry ): 
  if ( ry == 0 ): 

    #  180 Rotation in s*s matrix
    if ( rx == 1 ): 
      x = s - 1 - x 
      y = s - 1 - y 

    #  Transpose 
    t = x 
    x = y 
    y = t 

  return x, y

// : 몫
% : 나머지
& : AND 연산자
^ : XOR 연산자

 

CODE : 2D array generator

def xy2d_test ( ):
  n = 2 ** m
  a=len(str(n*n))+1

  print ( '==' )
  print ( '     '+' '*a, end = '' )
  for x in range ( 0, n ):
    print(( "{0:>"+str(a)+"}").format(x), end = '' )
  print ( '' )

  print ( '' )
  for y in range ( n - 1, -1, -1 ):
    print(( "  {0:>"+str(a)+"}:  ").format(y), end = '' )
    for x in range ( 0, n ):
      d = xy2d ( m, x, y )
      print(( "{0:>"+str(a)+"}").format(d), end = '' )
    print ( '' )
  print()

 

 

Reference

- https://people.sc.fsu.edu/~jburkardt/py_src/hilbert_curve/hilbert_curve.html

- http://bit-player.org/extras/hilbert/hilbert-construction.html