Primitivas y operaciones booleanas

Las primitivas de una SDF no se memorizan, se deducen. Cada una sale de responder a la misma pregunta —cuál es el punto más cercano de esta forma— con álgebra vectorial elemental, y entender la deducción es lo que te permite escribir la tuya cuando ninguna del catálogo encaja. Las operaciones que las combinan son todavía más simples, y son la razón de que valga la pena todo esto.

Tres primitivas y su deducción

Esfera. Ya la tienes: distancia al centro menos el radio.

float sdEsfera( vec3 p, float r ) {
  return length( p ) - r;
}

Caja. Esta es la interesante, porque su forma canónica confunde hasta que se ve de dónde sale. Una caja centrada en el origen con semilados b:

float sdCaja( vec3 p, vec3 b ) {
  vec3 q = abs( p ) - b;
  return length( max( q, 0.0 ) ) + min( max( q.x, max( q.y, q.z ) ), 0.0 );
}

El abs(p) reduce los ocho octantes a uno solo por simetría. q es entonces, componente a componente, cuánto sobresale el punto de la caja en cada eje: positivo si se sale por ese eje, negativo si está dentro del rango.

El primer término, length( max( q, 0.0 ) ), es la distancia por fuera. Al recortar las componentes negativas a cero te quedas solo con los ejes por los que el punto se sale, y su módulo es exactamente la distancia al punto más cercano de la caja: si se sale por un eje, la distancia a la cara; por dos, la distancia a la arista; por tres, la distancia al vértice. Los tres casos salen de la misma expresión, que es lo elegante del asunto.

El segundo término, min( max( q.x, max( q.y, q.z ) ), 0.0 ), es la distancia por dentro, y solo es distinto de cero cuando todas las componentes de q son negativas. En ese caso max de las tres es la que menos negativa es, es decir, la distancia a la cara más cercana, con el signo correcto. Y cuando el punto está fuera, este término se anula porque el min con cero lo mata.

Los dos términos son mutuamente excluyentes: uno u otro es cero. Sumarlos en lugar de ramificar evita un condicional en la GPU.

Toro. Un toro de radio mayor t.x y radio menor t.y alrededor del eje Y:

float sdToro( vec3 p, vec2 t ) {
  vec2 q = vec2( length( p.xz ) - t.x, p.y );
  return length( q ) - t.y;
}

La deducción tiene dos pasos. Primero se colapsa el problema tridimensional a uno bidimensional: length( p.xz ) - t.x es la distancia con signo desde el punto al círculo central, medida en el plano horizontal. Junto con la altura p.y, forma un punto en un plano de sección. Y en ese plano, el toro es simplemente una circunferencia de radio t.y centrada en el origen, así que se resuelve con la fórmula de la esfera. Este patrón —reducir a un problema de menor dimensión y resolver ahí— es el que genera casi todas las primitivas de revolución.

Las booleanas exactas

float opUnion( float a, float b )        { return min( a, b ); }
float opInterseccion( float a, float b ) { return max( a, b ); }
float opResta( float a, float b )        { return max( a, -b ); }

La unión es el mínimo porque un punto está fuera de la unión solo si está fuera de las dos formas, y la distancia a la más cercana es la que manda. La intersección es el máximo por el argumento simétrico. Y la resta es el máximo de a con el campo negado de b, porque negar una SDF invierte dentro y fuera: convierte el sólido en el complemento del sólido.

Estas tres operaciones no son aproximaciones. Para la unión, el resultado es una SDF exacta siempre que las dos entradas lo sean. Para la intersección y la resta, el resultado es una cota inferior: cerca de las aristas donde las dos superficies se cortan, el máximo subestima la distancia real. Eso es seguro para el algoritmo de marcha, solo que ligeramente más lento junto a esas aristas.

Este es el punto donde conviene detenerse a apreciar lo que acaba de pasar. Restar un cilindro de un cubo en una malla de triángulos es una operación de geometría computacional con casos degenerados, problemas de precisión numérica y bibliotecas dedicadas. Aquí es max( a, -b ).

La unión suave

Casi todo lo que hace bonito al raymarching viene de aquí. La unión con min produce una arista viva donde las dos formas se encuentran. La unión suave interpola entre ellas en una banda de anchura k, generando una transición orgánica que en mallas es un problema abierto.

La versión polinómica, que es la que se usa por defecto por su bajo coste:

// Union suave polinomica. k controla el radio de la fusion.
float smin( float a, float b, float k ) {
  float h = clamp( 0.5 + 0.5 * ( b - a ) / k, 0.0, 1.0 );
  return mix( b, a, h ) - k * h * ( 1.0 - h );
}

Léela por partes. h es un factor de mezcla que vale uno cuando a domina claramente, cero cuando domina b, y transita linealmente entre ambos en la banda de anchura k. mix( b, a, h ) es la interpolación lineal entre las dos distancias. Y el término final, - k * h * ( 1.0 - h ), es una parábola que vale cero en los extremos y k/4 en el centro: es lo que hunde el campo en la zona de encuentro y crea el filete cóncavo.

Las variantes de intersección y resta suaves salen por el mismo camino:

float smax( float a, float b, float k ) {
  return -smin( -a, -b, k );
}

float sResta( float a, float b, float k ) {
  return smax( a, -b, k );
}

Modificadores sobre cualquier forma

Estos operan sobre el campo o sobre el dominio y valen para cualquier primitiva.

Redondeo. Restar una constante al campo desplaza la superficie hacia fuera de forma uniforme, lo que redondea todas las aristas con ese radio:

float d = sdCaja( p, vec3( 0.5 ) ) - 0.08;   // cubo con aristas de radio 0.08

Ojo: la caja resultante mide 0.5 + 0.08 de semilado, no 0.5. Compensa restando el radio a las dimensiones si quieres conservar el tamaño exterior.

Vaciado. El valor absoluto del campo menos un grosor convierte un sólido en una cáscara:

float d = abs( sdEsfera( p, 1.0 ) ) - 0.02;   // cascara esferica de 0.04 de grosor

Repetición infinita. Aplicando módulo al dominio antes de evaluar, una primitiva se convierte en una retícula infinita a coste cero:

vec3 q = p - c * round( p / c );   // c es el periodo por eje
float d = sdEsfera( q, 0.3 );

Con round en lugar de mod la celda queda centrada, que es lo que casi siempre quieres. Y la repetición limitada, que suele ser más útil:

vec3 q = p - c * clamp( round( p / c ), -lim, lim );
float d = sdEsfera( q, 0.3 );

Elongación. Estirar una primitiva sin deformarla, insertando un tramo recto:

vec3 q = p - clamp( p, -h, h );   // h es el vector de elongacion
float d = sdEsfera( q, 0.3 );     // ahora es una capsula/cajon redondeado

Y así se monta la escena entera, que en raymarching es una sola función:

float mapa( vec3 p ) {

  float suelo  = p.y + 1.0;
  float bola   = sdEsfera( p - vec3( 0.0, 0.0, 0.0 ), 0.7 );
  float bloque = sdCaja( p - vec3( 0.9, -0.4, 0.0 ), vec3( 0.45 ) ) - 0.05;
  float anillo = sdToro( p - vec3( -1.0, 0.1, 0.0 ), vec2( 0.5, 0.15 ) );

  float cuerpo = smin( bola, bloque, 0.35 );
  cuerpo = min( cuerpo, anillo );

  return min( cuerpo, suelo );

}