Recalcular normales II: la vía numérica

Cuando la deformación no se puede derivar a mano —un mapa de alturas, un FBM de seis octavas, un domain warping— queda el método que funciona siempre: mover el vértice, mover también dos vecinos imaginarios, y reconstruir la normal a partir de los tres puntos resultantes. Cuesta tres evaluaciones en vez de una, funciona con cualquier función, y tiene un único parámetro que hay que elegir con criterio en vez de a ojo.

Diferencias finitas hacia delante

La idea completa cabe en seis líneas. Si desplazar( p ) es la función que mueve un punto, evalúala en el punto y en dos vecinos separados un épsilon en cada eje del plano tangente original. Los dos vectores que van del punto desplazado a los vecinos desplazados son las tangentes de la superficie deformada, y su producto vectorial es la normal.

uniform float uTiempo;
uniform float uAmplitud;
uniform float uEpsilon;

varying vec3 vNormal;
varying vec3 vPosicionVista;

// Aqui va snoise y fbm de las lecciones del nivel anterior

float altura( vec2 xy ) {
  return fbm( xy * 0.35 + vec2( 0.0, uTiempo * 0.05 ) ) * uAmplitud;
}

void main() {
  float e = uEpsilon;

  vec3 p  = position;
  vec3 px = position + vec3( e, 0.0, 0.0 );
  vec3 py = position + vec3( 0.0, e, 0.0 );

  p.z  += altura( p.xy );
  px.z += altura( px.xy );
  py.z += altura( py.xy );

  vec3 tangenteX = px - p;
  vec3 tangenteY = py - p;
  vec3 normalObjeto = normalize( cross( tangenteX, tangenteY ) );

  vNormal = normalMatrix * normalObjeto;
  vec4 posicionVista = modelViewMatrix * vec4( p, 1.0 );
  vPosicionVista = posicionVista.xyz;
  gl_Position = projectionMatrix * posicionVista;
}

Comprueba la orientación mentalmente con la superficie sin deformar: tangenteX vale ( e, 0, 0 ), tangenteY vale ( 0, e, 0 ), y su producto vectorial es ( 0, 0, e² ), que normalizado da ( 0, 0, 1 ). Es la normal correcta del PlaneGeometry sin tocar. Si en tu caso sale invertida, cambia el orden de los factores del cross.

Fíjate también en que las dos evaluaciones adicionales usan position original, no p ya desplazado. Evaluar la función sobre coordenadas ya movidas mezclaría dos espacios y daría una normal que no corresponde a nada.

Elegir el épsilon

Éste es el único parámetro y tiene un valor óptimo que se puede calcular, no adivinar.

Con un épsilon demasiado grande, las tangentes atraviesan detalle de la función y la normal resultante es un promedio suavizado. El terreno se ilumina como si fuera más liso de lo que es. No es un error catastrófico y a veces es lo que quieres.

Con un épsilon demasiado pequeño, la resta px - p opera sobre dos números casi idénticos y la cancelación catastrófica se lleva casi todos los bits significativos. Con un flotante de treinta y dos bits, restar dos valores que coinciden en seis cifras decimales deja dos cifras útiles, y la normal se convierte en ruido. El síntoma es un moteado de píxeles brillantes que aparece y desaparece.

El óptimo de una diferencia finita hacia delante está en la raíz cuadrada del épsilon de máquina multiplicada por la escala del problema. Para flotante de treinta y dos bits el épsilon de máquina vale unos 1.19 por diez elevado a menos siete, y su raíz cuadrada es aproximadamente 3.4 por diez elevado a menos cuatro. Multiplicado por la escala del terreno:

// Para un plano de 20 unidades de lado
const escala = 20;
material.uniforms.uEpsilon.value = 3.4e-4 * escala;   // unos 0.007

En la práctica, cualquier valor entre una milésima y una centésima de la extensión del objeto funciona. Lo que hay que evitar son los extremos: por debajo de 1e-4 en unidades absolutas hay ruido, y por encima de la separación entre vértices hay suavizado.

Y hay un uso deliberado del épsilon grande que conviene conocer: fijarlo exactamente en el espaciado entre vértices hace que la normal describa la malla real en lugar de la función continua, y por tanto que el sombreado y la silueta sean coherentes. Con PlaneGeometry( 20, 20, 256, 256 ) ese espaciado es 20 / 256, unos 0.078.

El coste

Tres evaluaciones en vez de una. Con una función barata da igual; con un FBM de seis octavas de simplex son dieciocho evaluaciones de snoise por vértice.

En una malla de 256 por 256 son 66049 vértices, así que 1.2 millones de evaluaciones de snoise por cuadro. Y si la malla proyecta sombra sobre dos luces, tres veces eso.

Frente a los cientos de millones de invocaciones de fragmento de las que hablábamos en el nivel 30, sigue siendo poco. La cuenta que hay que hacer es ésta: mientras el número de vértices sea mucho menor que el de píxeles cubiertos, el vertex shader puede permitirse ser caro. El punto en el que eso deja de ser cierto es cuando los triángulos bajan del tamaño de unos pocos píxeles, y ahí el problema ya no es el coste del método sino que la teselación está mal dimensionada.

Dos optimizaciones reales cuando el método pesa:

// A) Menos octavas para la normal que para la altura.
// El detalle fino apenas cambia la orientacion y si multiplica el coste.
float alturaDetalle( vec2 xy ) { return fbmSeisOctavas( xy ) * uAmplitud; }
float alturaGruesa( vec2 xy )  { return fbmTresOctavas( xy ) * uAmplitud; }

// B) Reutilizar el valor central. Con diferencias hacia delante ya lo haces:
// tres evaluaciones dan dos derivadas. Con diferencias centradas serian cuatro
// y solo ganarias precision de segundo orden, que aqui no compensa.

La tercera vía y la comparativa

Derivadas de pantalla

Hay una forma de obtener normales sin evaluar la función ni una sola vez de más: calcularlas en el fragment shader a partir de las derivadas de la posición interpolada.

// Fragment shader
varying vec3 vPosicionVista;

void main() {
  vec3 N = normalize( cross( dFdx( vPosicionVista ), dFdy( vPosicionVista ) ) );
  ...
}

Coste: dos instrucciones. Y produce la normal geométrica exacta del triángulo, porque dentro de un triángulo la posición interpolada es lineal y sus derivadas son constantes.

Ahí está a la vez su virtud y su límite. Es exacta para la malla, así que la iluminación y la silueta coinciden siempre. Pero es plana por triángulo: el resultado es sombreado facetado, con cada polígono visiblemente distinto de su vecino. Con una teselación muy densa eso puede pasar por suave; con una normal es inaceptable.

Tres usos donde es la respuesta correcta. Cuando quieres facetado a propósito, que es un recurso estético frecuente. Cuando la malla es tan densa que el facetado no se distingue. Y como respaldo de depuración: comparar el sombreado con derivadas de pantalla contra el sombreado con tu normal calculada es la forma más rápida de saber si tu normal está mal, porque la primera nunca miente sobre la geometría.

Y recuerda las dos condiciones de las derivadas: solo existen en el fragment shader, y exigen control de flujo uniforme. Calcularlas dentro de una rama divergente da resultados indefinidos.

Las tres vías en una tabla