Midiendo el rango dinámico de un sensor de imagen

En este artículo vamos a desgranar paso a paso la metodología clásica para calcular el Rango dinámico de un sensor digital de imagen, basándonos en mediciones precisas de ruido hechas sobre capturas RAW de parches de color uniforme.

Se define el rango dinámico (RD) de un sensor como el intervalo de niveles de luminosidad (típicamente medido en pasos EV) en los que dicho sensor es capaz de capturar información útil. El RD está limitado en su extremo superior por la propia saturación del sensor y en el extremo inferior por el nivel de ruido máximo aceptable.

Vamos a echar mano de una serie de capturas RAW amablemente compartidas por Hugo Rodríguez procedentes de su Olympus OM-1, una para cada ISO de la cámara. La carta de parches de color usada por Hugo es la que Bill Claff (Photons to Photos) propone aquí para calcular el RD, y consiste en proyectar sobre el monitor una cuadrícula de 11x7 parches de color magenta:

El color magenta no es caprichoso, el astuto de mi tocayo lo ha escogido sabedor de que el canal G siempre obtiene mayores niveles en el RAW. Con parches magenta (pobres en verdes), los niveles RAW del canal G se hacen más parejos a los R y B proporcionando un set de muestras más uniforme. La toma a ISO12800 tenía este aspecto:

Que la toma esté desenfocada no solo no es un error sino que es imprescindible para no capturar ningún tipo de textura (en este caso de los píxeles del monitor) que sería interpretada como ruido adicional, falseando todos los cálculos de relación S/N.

También es importante que no exista ningún gradiente de iluminación causado por viñeteo de la óptica o por reflejos en la pantalla, el cual igualmente incrementaría la medición de ruido. Siendo el área útil que aprovecharemos en cada parche pequeña en relación a la imagen total, estos gradientes sería raro que supusieran un problema.

~~~

El siguiente paso es disponer de los niveles de negro y de saturación RAW del sensor bajo estudio. En cuanto al primero, seguramente por adecuación a la electrónica del conversor A/D y a la aritmética digital con enteros, los fabricantes sitúan deliberadamente el nivel de negro (valor medio que tendría un RAW obtenido haciendo una foto con el objetivo tapado) en las potencias de 2:

• Sensores de 10 bits (1.024 niveles): nivel de negro 64 o 128
• Sensores de 12 bits (4.096 niveles): nivel de negro 256
• Sensores de 14 bits (16.384 niveles): nivel de negro 512 o 1.024

Con el comando dcraw -v -D -t 0 -4 -T *.DNG extraemos en formato TIFF los datos de todos los archivos RAW. Siendo la Olympus OM-1 una cámara de 12 bits lo esperable sería encontrar el nivel de negro en el valor 256, lo que comprobamos con el histograma RAW de la toma ISO800:

Comentamos que el rango dinámico del sensor está limitado en el extremo derecho por la saturación, por lo que necesitamos también conocer en qué valor satura para referir todos los demás niveles RAW a dicha saturación. Podría pensarse que los sensores saturan en el fondo de su escala de bits pero no siempre es así, aunque éste sí lo hace:

Que la toma ISO25600 sature en 4.095 no garantiza que las de ISOs inferiores lo hagan, pero ningún RAW por debajo de éste tiene píxeles saturados así que asumiremos este valor de saturación para todos. El error en cualquier caso será limitado; desviaciones en este valor son menos determinantes que en el de negro.

Los niveles RAW los leemos como valores enteros y los normalizamos al rango 0..1 en coma flotante, permitiendo valores normalizados negativos:

 BLACK=256  # Olympus OM-1 black level

 SAT=4095  # Olympus OM-1 sat level

 img=readTIFF("iso25600.tiff", as.is=TRUE)

 img=img-BLACK

 img=img/(SAT-BLACK)

~~~

Ahora vamos a realizar la corrección de perspectiva de la captura original. Esta corrección no es imprescindible pero tener los parches alineados facilita mucho la tarea de leer los datos en ellos.

Previamente a la corrección debemos deshacer la matriz de Bayer para procesar de forma individual cada uno de los 4 canales RAW (R, G₁, G₂, B). En nuestro caso escogemos el canal R del RAW asumiendo que la respuesta de los demás será equivalente. Superado el efecto del filtro Bayer no hay motivo para pensar que los fotocaptores de un color puedan tener un rendimiento diferente a los demás (más sobre esto al final del artículo).

Por simplicidad y también por rigor, la corrección la hacemos sin interpolar valores. Realizar un muestreo nearest neighbour implicará tanto desechar algunos píxeles origen como duplicar otros. Dado que la selección es aleatoria ninguno de los efectos afectará a la precisión del cálculo posterior.

Para la corrección de perspectiva se ha utilizado el método visto en 'Transformación trapezoidal de imágenes con R (I). Algoritmo', pero optimizado en C++ con Rcpp reduciendo el tiempo de cálculo por un factor 240!!!

~~~

El siguiente paso es el meollo de todo el ejercicio donde recorremos individualmente los parches calculando en cada uno:

• El nivel de Señal útil: promedio de los valores RAW
• El nivel de Ruido: desviación estándar de los valores RAW
• La Relación S/N: cociente entre los anteriores

Desechamos los parches en los que el nivel de Señal resulte un valor negativo, la relación S/N sea inferior a -10dB (zonas muy ruidosas que no aportan nada) o el número de píxeles cercanos a saturar (>90% en escala lineal) llegue al 1% del total. Este análisis de los parches también lo hemos acelerado en C++.

En la toma ISO25600 por ejemplo descartamos los primeros 3 parches, usando el dato en los restantes 74 del total de 77 parches:

A partir de las parejas de Señal y relación S/N obtenidas, cada una proveniente de un parche, ya podemos intuir la forma de las Curvas de relación S/N en un scatterplot. En el eje X representamos los niveles de Señal en escala logarítmica respecto a la saturación (EV) y en el eje Y los valores de relación S/N igualmente en escala logarítmica (dB) (hacer clic para ver en alta resolución):

~~~

Ahora solo queda calcular el rango dinámico para cada ISO como el margen desde la saturación (0EV) hasta el nivel de Señal en que la relación S/N caiga a cierto valor mínimo. Para ello es preciso establecer un criterio umbral de relación S/N mínima aceptable, y aquí usaremos los dos baremos habituales:

• RD fotográfico (práctico): pasos desde saturación hasta S/N=12dB
• RD ingenieril (teórico): pasos desde saturación hasta S/N=0dB

Para poder automatizar este cálculo, parametrizamos curvas splines cúbicas que aproximen suavemente, es decir sin necesariamente pasar por ellos, las series de puntos correspondientes a cada valor ISO. Usamos la relación S/N como variable independiente y el nivel de Señal como variable dependiente, ya que los umbrales de 12dB y 0dB son valores de relación S/N de entrada para los que queremos obtener sus correspondientes niveles de Señal (el RD).

Las splines serán más efectivas si las calculamos en el dominio logarítmico para ambos ejes, tanto el del nivel de Señal (EV) como el de relación S/N (dB). La transformación genera curvas más suaves que evitan los saltos bruscos de la segunda derivada que penalizan el cálculo de splines.

Sobre la misma gráfica se marcan los dos criterios de rango dinámico y los valores obtenidos de rango dinámico como intersección entre las curvas y cada uno de los criterios umbral (hacer clic para ver en alta resolución):

Las curvas hasta ISO800 tienen un comportamiento parejo, el cual cambia en ISOs superiores a ISO800 donde mejora el rendimiento en las sombras profundas. Esto puede constatarse en la gráfica de RD para esta cámara de Bill Claff que precisamente presenta un salto en el RD tras ISO800, ver aquí. En sus propias palabras: "(...) the Olympus System OM-1 looks to have High Conversion Gain (HCG) starting at ISO1000 (and Noise Reduction starting at ISO16000)".

La siguiente tabla resume el RD calculado para cada ISO y criterio, indicando además las muestras o puntos de la curva (parches) que participaron del total de 77 parches:

La ISO invariancia plena se alcanza a ISO1600. Esto quiere decir que a efectos prácticos ningún valor por encima de ISO1600 tiene utilidad disparando en RAW con esta cámara si buscamos mejorar el ruido, cuando haciéndolo podríamos perder altas luces:

~~~

Vamos a terminar el ejercicio haciendo un estudio de sensibilidad del método descrito a tres variables que hemos elegido en el proceso:

• El canal RAW analizado (R, G₁, G₂, B).
• El nivel de negro del sensor usado para linealizar y normalizar los niveles RAW.
• El tamaño de los parches sobre los que hemos hecho las mediciones de ruido.

Sobre el canal RAW analizado, todo el ejercicio se ha hecho con el canal R, casualmente el que mayores niveles RAW obtuvo en los parches magenta. El resultado sería diferente usando los otros canales?. Aquí puede verse la superposición del cálculo para los 4 canales RAW (R, G₁, G₂, B), arrojando todos ellos curvas, y por tanto valores de RD, muy similares (hacer clic para ampliar):

Tiene mucha más influencia, y por eso debemos resaltarlo como punto crítico del método, el nivel de negro escogido (BLACK=256 en nuestro caso). El cálculo de RD a los ISOs más bajos puede fluctuar del orden de 0,1EV solo por usar el nivel de negro adyacente (BLACK=255) (hacer clic para ampliar):

Esta sensibilidad disminuye conforme aumentamos el ISO o el número de bits del sensor, al reducirse el error relativo que supone desplazar BLACK en una unidad. En conclusión es más complicado medir con precisión el RD de un sensor cuanto menor sea el ISO o sus archivos RAW se codifiquen con menos bits.

En el siguiente histograma de la captura ISO200 vemos que solo se tienen 58 valores tonales por canal para codificar toda la información de la escena. Mover en una sola unidad el nivel de negro tiene un gran impacto a la hora de situar correctamente los niveles de exposición respecto a la saturación, lo que trastoca todos los puntos de las gráficas:

El último parámetro de sensibilidad a evaluar es el tamaño de los parches considerados. Cuanto más grandes tendremos más muestras y por tanto mayor robustez estadística, pero también más peligro de interpretar como ruido efectos indeseables como la invasión de parches cercanos por mal alineamiento o por el desenfoque de la óptica, o directamente incorporar gradientes de luminosidad en la captura (viñeteos, reflejos,...).

Comparando parches de tamaños bastante dispares (separación entre parches de 100 píxeles vs 160 píxeles, suponiendo los primeros 4 veces más muestras) vemos que la diferencia es mínima, lo que valida las capturas usadas en este aspecto (hacer clic para ampliar):

~~~

Como guinda del ejercicio elaboramos sobre el nivel de negro de la cámara estudiada. Hugo me ha enviado dos archivos RAW adicionales a ISO200: uno completamente quemado donde comprobamos que el nivel de saturación usado de 4.095 era correcto, y otro con el objetivo tapado, la cámara en frío y alta velocidad de obturación para minimizar el ruido, cuyo histograma RAW es así:

El histograma lo dice todo: el mejor candidato entero para situar el nivel de negro del sensor es 255, no 256 como hemos usado en todo el tutorial. Pero podemos ir un poco más lejos y fijarnos en que el histograma RAW tiene más píxeles con valor 254 que 256. Esto hace pensar que el verdadero nivel de negro analógico que genera el sensor no caería exactamente en ninguno de los valores enteros que arroja el conversor A/D, sino que debería estar algo por debajo de 255.

Siendo así un valor de negro aún más preciso podría ser la media de los valores RAW (BLACK=254.85). Estoy preparando un artículo donde se explica por qué ese promedio tan inmediato de obtener resulta ser una buena estimación del punto preciso dónde caería el valor de negro analógico a la entrada del conversor A/D.

Usar valores con decimales tiene sentido si el procesado se va a hacer con software que trabaje en coma flotante, como es nuestro cálculo de RD. Dudo que los reveladores comerciales estándar pudieran aprovechar tal nivel de refinamiento ya que a efectos de la imagen resultante la diferencia respecto a usar enteros tampoco se iba a notar demasiado.

Finalmente las curvas más precisas que podemos calcular (BLACK=254.85):

~~~

Repositorio con el código R: GitHub.