Algoritmo Marching Squares

En el post anterior, comenté que no hay solamente una forma de representar voxels en pantalla. Mostré las formas mas simples y mencioné el algoritmo Marching Cube como una posibilidad para suavizar el renderizado de voxels.

En este post, en lugar de mostrar el algoritmo Marching Cubes voy a mostrar su versión 2D: Marching Squares. El algoritmo Marching Squares permite hacer una analogía directa con Marching Cubes y puede ser comprendido mas fácilmente.

La idea entonces, es que tenemos un campo de voxels 2D. Una posible representación de ese campo de voxels es un simple array bidimensional de números enteros:

  v = [[0,     0,   0, 0],
       [30,  100,  30, 0],
       [100, 100, 100, 0],
       [ 30, 100,  30, 0]];

Una forma de representar esto sería utilizando una escala de grises con 0 = blanco y 100 = negro. Pero para  evitar un resultado pixelado sería necesario un array de mayor dimensión.

Si lo que buscamos es usar un array de tamaño reducido para representar los voxels, pero a su vez un resultado no tan pixelado, podemos generar una "curva de nivel". Esto es, un contorno continuo de este campo usando segmentos de lineas. Es ahi donde entra Marching Squares.

Algoritmo Marching Squares

Este algoritmo permite obtener una "curva de nivel" de una función discreta (en este caso, los voxels). Los pasos a seguir son:

1. Convertir los voxels a un campo de valores 0 o 1 según un isovalor predeterminado.
2. Seleccionar un sub cuadrado 2x2 del campo de valores 0/1 y determinar un caso utilizando una tabla de casos (ver table mas adelante).
3. Utilizar los valores del campo de voxels original junto con el caso determinado en el paso 2 para obtener una linea interpolada.
4. Repetir desde el paso 2 para cada subcuadrado 2x2 del array original.

Utilizando el ejemplo anterior sería:

1. Convertir a campo de valores:

En este paso elijo un isovalor de 50 (punto medio entre 0 y 100) y creo un nuevo array que tiene 0 cuando el valor v[i][j] es menor o igual a 50 y 1 cuando v[i][j] es mayor a 50.

  v_iso = [[0, 0, 0, 0],
           [0, 1, 0, 0],
           [1, 1, 1, 0],
           [0, 1, 0, 0]];

2. Determinar casos:

Viendo cada subarray de 2x2 y comparando con la siguiente tabla de casos, determino el caso.

Tabla de casos. Un circulo negro corresponde a un voxel con valor superior al isovalor, uno blanco corresponde a un valor inferior. Los extremos de la linea resultante tienen que ser interpolados utilizando los valores de los voxels

En nuestro ejemplo, volviendo a v_iso, el primer subarray de 2x2 es

  [[0, 0], 
   [0, 1]

Comparando con la imagen, esto corresponde al caso 2. Ese caso es una linea diagonal de abajo hacia la derecha.

3. Interpolar usando voxels:

Si miramos en la tabla, cada extremo de las lineas celestes tiene una flecha que indica que ese extremo deberá ser ubicado teniendo en cuenta los valores de los extremos del lado. Tomemos el ejemplo encontrado en el paso anterior. Los valores de cada voxel eran 0, 0, 30 y 100, que resultaron en un Caso 2:

En este punto, hay que tener en cuenta que lo que estamos buscando es una linea correspondiente a un valor de 50. Si la distancia vertical y horizontal entre voxels es L y los valores en la posición de cada voxel son los valores dados por el voxel, entonces podemos interpolar las posiciones de cada extremo de la linea usando una proporción.

(100-30)/L = (50-30)/x

que implica que x es L*(50-30)/(100-30) = 0.28L

De la misma manera, para obtener la coordenada Y de la posición del otro extremo, podemos usar una proporción similar. En este caso, dado que los extremos corresponden a 0 y a 100, el valor de y es en el centro exacto de los voxels, o sea y=L*(50-0)/(100-0) = 0.5L

4. Repetir cada paso:

Si repetimos los pasos anteriores para cada subarray de 2x2 obtenemos el siguiente gráfico:

En el dibujo ya están las posiciones interpoladas de cada linea y un sombreado marcando el interior del polígono resultante. El interior del polígono se puede determinar fácilmente modificando la tabla de casos para mostrar cual es el interior en cada caso.

Ejemplo interactivo

La mejor forma de entender el algoritmo es utilizando un ejemplo interactivo y viendo el código. Hice un pequeño widget para poder experimentar. El ejemplo está realizado completamente en HTML5 usando webgl y javascript, por lo que es necesario un navegador moderno. El código fuente está disponible en google code.

Para usar el ejemplo solamente hay que asignar valores en la tabla. Se puede desactivar y activar el renderizado de voxels y de la isolinea.

Mostrar Voxels Mostrar Isolinea

Isovalue

En el  próximo artículo voy a mostrar cómo se generaliza esta idea a un caso 3D con el algoritmo Marching Cubes.