Múltiples vistas

1Rayos que deberían cortarse

Una sola imagen fija el rayo en que está un punto del espacio, pero no en qué lugar de ese rayo. Los apuntes de geometría epipolar usaron una segunda imagen para encontrar dónde puede aparecer el punto en ella. Estos apuntes van dos pasos más allá. Dados los píxeles de un punto en dos o más imágenes, encuentran el punto mismo. Y dados sus píxeles en dos imágenes, predicen dónde aparece el punto en una tercera.

Un punto M visto en la imagen i como el píxel m_i está sobre el rayo que sale del centro óptico C_i y pasa por m_i. Está sobre el rayo de cada imagen que lo ve, así que es donde esos rayos se cortan. Encontrar un punto a partir de los rayos de sus imágenes se llama triangulación, y encontrarlo así es la reconstrucción 3D. Aquí esa palabra significa localizar un punto. No es lo que reconstruye la tomografía computarizada, que es de qué está hecho un objeto en cada uno de sus puntos.

Los píxeles reales se miden y las cámaras se estiman, así que los rayos casi nunca se cortan exactamente, ni siquiera dos de ellos. La pregunta pasa a ser qué punto del espacio explica mejor los píxeles.

El ejemplo

El ejemplo de esta página está hecho para ella, con números redondos donde se puede. Una fuente de rayos X, el centro óptico, puede ubicarse en cualquier punto de un círculo horizontal de 1000 mm de radio en torno a un eje vertical Z, y cada imagen que toma es de 2000 × 1600 píxeles. La imagen 1 se toma desde un punto de ese círculo, la imagen 2 desde un cuarto de vuelta más allá, y la imagen 3 desde 175°, casi al otro lado del objeto respecto de la imagen 1. El objeto es una barra plana de 45 mm de largo, en punta en ambos extremos, que está a 12° de la vertical con su punta inferior a la altura de la fuente. Sus dos puntas se marcan en cada una de las tres imágenes y, como pasa con los clics hechos a mano, cada clic queda a unos pocos píxeles de su lugar. Como la barra es conocida, lo que dan las imágenes se puede comparar con dónde están de verdad sus puntas. Las longitudes están en milímetros y las coordenadas de imagen en píxeles.

2Dos ecuaciones por vista

Para una imagen, escribe su cámara P con entradas p_jk, el píxel como m = [x y 1]^T y el punto como M = [X Y Z 1]^T. La proyección λm = PM son tres ecuaciones:

Las incógnitas son X, Y y Z, y también λ, una nueva por cada imagen. La tercera ecuación da λ directamente: es la profundidad de M delante de la cámara, en milímetros en este ejemplo, ya que las tres primeras entradas de la tercera fila de P forman un vector unitario. Al reemplazarla en las dos primeras, λ desaparece:

Por eso cada imagen da dos ecuaciones: aporta tres ecuaciones y una incógnita nueva, λ, y una ecuación se gasta en ella. Las dos que quedan son lineales en X, Y y Z, porque el píxel (x, y) es un número conocido.

Cada una tiene una imagen geométrica. La primera se cumple para todo punto cuya imagen tiene esta x, sea cual sea su y: los puntos que se proyectan sobre una columna de la imagen. Forman un plano que pasa por el centro óptico. La segunda es el plano que pasa por el centro y se proyecta sobre la fila y. El rayo de m es donde se cortan los dos planos. Así, una imagen da una recta, no un punto, y hace falta una segunda imagen para fijar un punto sobre ella.

3Mínimos cuadrados con n vistas

Cada una de las n imágenes da sus dos filas, con su propia cámara y su propio píxel. Apiladas, forman un solo sistema Q[X Y Z]^T = r, donde Q es de 2n × 3 y r tiene 2n componentes, con las dos filas de la imagen i en los lugares 2i − 1 y 2i. Como en los apuntes de calibración de cámaras, los mínimos cuadrados dan el punto que hace más pequeña la suma de los cuadrados de las componentes de Q[X Y Z]^T − r:

El sombrero marca una estimación. Llama M a ese punto. Existe cuando Q tiene rango 3, es decir, cuando los 2n planos no comparten todos una misma recta.

Dos vistas y la restricción epipolar

Dos imágenes dan cuatro ecuaciones para tres incógnitas, una más de las necesarias, así que las cuatro tienen una solución común solo para algunos pares de píxeles. Escríbelas como [Q | −r][X Y Z 1]^T = 0: un sistema de 4 × 4 con lado derecho 0, que tiene una solución distinta de 0 solo cuando det [Q | −r] = 0. Cada fila de [Q | −r] es xp₃ − p₁ o yp₃ − p₂ para uno de los dos píxeles, con p_j las filas de su cámara. Por eso el determinante es lineal en m₁ y lineal en m₂, que es la forma de m₂^TFm₁, y se anula para todo par de puntos correspondientes. Es la restricción epipolar, multiplicada por un factor que depende solo de las cámaras. Para las imágenes 1 y 2, con F = [BC₁]_×BA⁺ como en los apuntes de geometría epipolar, el factor es −6.40·10⁷, sean cuales sean los píxeles. Es el determinante det G de esos apuntes, con λ₁ y λ₂ ya eliminadas.

Así, dos rayos se cortan exactamente cuando sus píxeles cumplen la restricción epipolar. Los clics del ejemplo no la cumplen. En la imagen 2 están a 10.83 y 6.86 píxeles de las líneas epipolares de sus pares en la imagen 1, así que sus rayos no se cortan. Los mínimos cuadrados reparten la diferencia: solo con las imágenes 1 y 2, el M de la punta superior se proyecta a 5.55 píxeles de su clic en la imagen 1 y a 5.41 del de la imagen 2.

Vistas enfrentadas

Los rayos de dos imágenes deben cruzarse con un ángulo claro. Las fuentes de las imágenes 1 y 3 están a 175° una de otra, a ambos lados de la barra, así que los dos rayos que pasan por uno de sus puntos son casi una misma recta: se cortan en unos 172°. A lo largo de esa recta, las dos imágenes apenas ven diferencia, y unos pocos píxeles de error mueven M lejos sobre ella. Solo con estas dos imágenes, la punta inferior queda a 10.99 mm de donde está y la superior a 3.76 mm, y la barra resulta de 48.21 mm de largo. Las imágenes 1 y 2 dan 45.34 mm, y las imágenes 2 y 3 dan 44.43 mm, con ambas puntas a cerca de un milímetro de donde están. En las ecuaciones, las columnas de Q son casi dependientes: sus valores singulares son 11319, 11292 y 786.1, frente a 11325, 8014 y 7987 para las imágenes 1 y 2.

Más vistas

Con tres imágenes, Q es de 6 × 3. Los errores de los clics se compensan en parte, y ningún par de imágenes decide solo el resultado. Las imágenes 1, 2 y 3 juntas dan 45.17 mm, con las puntas a 0.51 y 0.56 mm de donde están.

4Qué minimiza el ajuste

Los mínimos cuadrados hacen pequeño QM − r, escribiendo M por [X Y Z]^T, pero ¿qué es ese residuo? Toma la primera fila de la imagen i. Su lado izquierdo menos su lado derecho es (p₃ · M)x_i − p₁ · M, escribiendo p_j · M por el producto de la fila j de P_i con [X Y Z 1]^T. El primer producto es la profundidad λ_i de M, y el segundo es λ_ix_i, donde m_i = (x_i, y_i) es donde M se proyecta en la imagen i. Así, cada fila del residuo es λ_i(x_i − x_i) o λ_i(y_i − y_i), una profundidad en milímetros por un error en píxeles, y los mínimos cuadrados minimizan

Es el error en píxeles ponderado por el cuadrado de la profundidad, en mm·px: pesa más una imagen cuya fuente está más lejos del punto. La cantidad con un significado directo es el error de reproyección, la distancia en píxeles de cada m_i marcado a la proyección m_i del punto encontrado, con λ_im_i = P_iM. El mejor punto en ese sentido minimiza

Ninguna fórmula da ese punto, porque x_i y y_i dividen por la profundidad, que depende de M. Se encuentra iterando, con un método de gradiente por ejemplo, a partir de M. Cuando todas las profundidades se parecen, ambas respuestas casi coinciden. Aquí van de 985 a 1020 mm: minimizar J mueve las puntas 0.026 y 0.024 mm, y cambia el largo en 0.028 mm.

5Una tercera vista a partir de dos

Ahora tres imágenes, con cámaras A, B y C. Dados los puntos correspondientes m₁ y m₂ en las dos primeras, ¿dónde aparece su punto del espacio en la tercera? La geometría epipolar ya lo responde. m₃ está sobre la línea epipolar de m₁ en la imagen 3, y sobre la de m₂:

donde F₁₃ es la matriz fundamental de las imágenes 1 y 3, a partir de A y C, y F₂₃ la de las imágenes 2 y 3, a partir de B y C. Así, m₃ es donde se cruzan las dos líneas, y el punto donde se cruzan dos rectas es su producto cruz, como en los apuntes de coordenadas homogéneas:

También es una prueba. m₁, m₂ y m₃ son correspondientes exactamente cuando m₃ es ese punto: la condición es necesaria y suficiente.

Cuando las dos líneas son una

El cruce no está definido cuando ℓ₁₃ y ℓ₂₃ son la misma línea. ℓ₁₃ es donde la imagen 3 corta el plano epipolar de M, C₁ y C₃, y ℓ₂₃ es donde corta el plano de M, C₂ y C₃. Ambos planos contienen M y C₃, y son el mismo plano cuando además contiene a C₁ y C₂. Eso ocurre de dos maneras:

• Los tres centros ópticos están sobre una recta, como cuando una cámara toma las tres imágenes moviéndose en línea recta. Es frecuente, y entonces ocurre para todo punto.
• Los centros no están sobre una recta, pero M está sobre el plano que pasa por los tres.

Cerca de cualquiera de los dos casos, las líneas se cruzan con un ángulo pequeño, y un error pequeño en m₁ o m₂ mueve su cruce lejos. Y aun cuando las líneas se cruzan bien, el método da dos pasos, dos líneas y luego su cruce, en vez de una sola fórmula.

En el ejemplo

Aquí todas las fuentes están en un mismo círculo, así que el plano que pasa por tres de ellas es el plano del círculo, Z = 0. Pasa por todos los centros ópticos, así que cada cámara lo ve de canto, como una recta, y es la misma recta en todas las imágenes: la fila v = 800, a la altura del eje de cada cámara.

La punta inferior de la barra está sobre ese plano. Al predecirla en la imagen 3 a partir de sus clics en las imágenes 1 y 2, las dos líneas epipolares se cruzan en apenas 0.46°, y su cruce, (298.48, 797.28), cae a 842.56 píxeles de donde se marcó la punta en la imagen 3, (1141, 805). Justo sobre esa fila, un píxel de error en m₂ mueve el cruce más de cien píxeles, y desde los píxeles verdaderos de la punta las dos líneas son la misma fila v = 800, sin ningún cruce único. El plano explica también por qué las imágenes 1 y 3 solas dejan esta punta a 10.99 mm de su lugar: para un punto sobre él, las filas de ambas imágenes dan el mismo plano, Z = 0, y solo sus columnas, casi paralelas, lo ubican a lo largo de los rayos. La punta superior está 44 mm por encima del plano. Ahí las líneas se cruzan en 37.93°, y su cruce cae a 25.30 píxeles de su clic.

6El tensor trifocal

Una forma de evitar las líneas epipolares parte de las tres proyecciones mismas. Elige una matriz H de 4 × 4 cuyas tres primeras filas sean A, y cambia de coordenadas a M̃ = HM. Como A = [I | 0]H,

y λ₁m₁ = ÃM̃, λ₂m₂ = B̃M̃, λ₃m₃ = C̃M̃: las mismas proyecciones, en coordenadas donde la primera cámara es lo más simple posible. H debe ser invertible. Agregar la fila [0 0 0 1] bajo A sirve siempre que C₁ sea un punto finito, porque entonces el bloque izquierdo de 3 × 3 de A es invertible.

El rayo de m₁

λ₁m₁ = [I | 0]M̃ dice que las tres primeras componentes de M̃ son λ₁m₁. Con la cuarta igual a 1,

A medida que λ₁ recorre todas las profundidades, M̃ recorre el rayo de m₁. Encontrar el punto es encontrar λ₁.

La profundidad a partir de una coordenada

Proyecta este M̃ en la imagen 2 con B̃. Escribe b̃_j para las tres primeras entradas de la fila j de B̃, como fila, y b̃_j4 para la última. La primera coordenada del píxel es x₂ = (λ₁b̃₁m₁ + b̃₁₄) / (λ₁b̃₃m₁ + b̃₃₄), y al despejar λ₁:

Es la profundidad a la que el rayo de m₁ se proyecta sobre la columna x₂ de la imagen 2: donde cruza el plano que pasa por C₂ y cuya imagen es esa columna, como en la sección 2. El mismo paso en la imagen 3, con las filas c̃_k de C̃ y la columna x₃, da λ₁ = (c̃₁₄ − c̃₃₄x₃) / (x₃c̃₃m₁ − c̃₁m₁).

Cuatro trilinealidades

Si m₁, m₂ y m₃ son correspondientes, las dos profundidades son el mismo λ₁. Al igualarlas, multiplicar y agrupar los términos en 1, x₂, x₃ y x₂x₃, se obtiene una ecuación lineal en m₁, en x₂ y en x₃:

Cada T^jk es aquí una columna de tres números, T^jk = [T₁^jk T₂^jk T₃^jk]^T, con

para i, j, k = 1, 2, 3. Usar y₂ en vez de x₂, o y₃ en vez de x₃, da tres ecuaciones más del mismo tipo. Las cuatro son las trilinealidades:

Los 27 números T_i^jk forman el tensor trifocal de las tres imágenes. Como F, dependen solo de las cámaras, no de los puntos, así que se calculan una vez para todos los puntos. Y como F, son homogéneos: otra cuarta fila para H multiplica los 27 por un mismo número, lo que no cambia ninguna de las cuatro ecuaciones. Las cuatro ecuaciones son independientes, y se cumplen solo para puntos correspondientes: cuando se cumplen las cuatro, m₁, m₂ y m₃ son correspondientes.

Las mismas ecuaciones a partir de un determinante

Como con F en los apuntes de geometría epipolar, las mismas ecuaciones salen de las proyecciones escritas como un solo sistema, Gv = 0. Aquí G es de 9 × 7, con las filas de Ã, B̃ y C̃ en sus cuatro primeras columnas y m₁, m₂ y m₃ en las tres últimas, y v = [M̃; −λ₁; −λ₂; −λ₃]. Los puntos correspondientes dan una solución v ≠ 0, así que toda matriz de 7 × 7 formada por 7 de las 9 filas de G tiene determinante 0. Dejar fuera dos filas de una cámara deja solo las otras dos cámaras, y da su restricción epipolar, como m₃^TF₂₃m₂ = 0. Dejar fuera una fila de cada una de dos cámaras conserva las tres, y esos determinantes dan las trilinealidades.

7Transferencia con el tensor

Conociendo m₁ y m₂, las dos primeras trilinealidades son lineales en x₃ y en y₃. Al despejarlas, dan directamente el punto en la imagen 3:

Por qué funciona: x₂ fija λ₁, que es donde el rayo de m₁ cruza el plano que pasa por C₂ y la columna x₂, y ese punto del espacio, proyectado por C, es m₃. Es una reconstrucción a partir de una imagen y media, seguida de una proyección en la tercera. No se cruzan líneas epipolares, así que el plano de los tres centros no hace daño. Para la punta inferior de la barra, la transferencia da (1152.43, 803.88), a 11.49 píxeles del clic en la imagen 3, donde el cruce de las líneas epipolares quedaba a 842.56 píxeles. Lo que queda es el error propio de los clics: en las imágenes 1 y 2, están a 6.86 píxeles de la línea epipolar del otro.

Su propio punto débil

Si la línea epipolar de m₁ en la imagen 2 es la propia columna x₂, el plano que pasa por C₂ y esa columna contiene todo el rayo de m₁, y x₂ no puede fijar λ₁. Entonces hacen el trabajo las dos últimas trilinealidades, que usan y₂. En el ejemplo, las líneas epipolares son casi filas, ya que los epipolos están en la fila v = 800, casi siempre lejos a los lados de las imágenes, así que la forma con x₂ es la que conviene. La forma con y₂ falla en esa fila, igual que el cruce de líneas epipolares. Para imágenes cuyas fuentes están enfrentadas, como la 1 y la 3, los epipolos caen dentro de las imágenes, en (1349.29, 800) en la imagen 1. Si son las dos primeras imágenes, la forma con x₂ falla entonces para píxeles m₁ cercanos a la columna que pasa por ese epipolo, cuyas líneas epipolares también son columnas.

¿Es una tercera correspondencia?

m₁, m₂ y m₃ son correspondientes si m₁ y m₂ cumplen la restricción epipolar, y los tres cumplen las dos primeras trilinealidades, o las dos últimas. Las dos primeras usan x₂ pero nunca y₂, y es la restricción epipolar entre m₁ y m₂ la que verifica y₂. Con puntos medidos, ambas pruebas llevan una tolerancia:

con F₁₂m₁ = [a b c]^T, como en los apuntes de geometría epipolar. La prueba epipolar va primero porque la transferencia nunca mira y₂. Emparejados al revés, el clic de la punta superior en la imagen 1 con el de la punta inferior en la imagen 2, los clics quedan a 350.93 píxeles de la línea epipolar del otro. Pero la barra se ve vertical en la imagen 2, así que ahí sus dos clics tienen casi el mismo x₂, y la transferencia cae en (1068.56, 447.49), a 12.86 píxeles del clic de la punta superior en la imagen 3, como si los tres correspondieran.

8Cuatro vistas y más

Con una cuarta imagen, las mismas pruebas se encadenan. m₁, m₂, m₃ y m₄ son correspondientes si m₁, m₂, m₃ son correspondientes y m₂, m₃, m₄ también: m₂ y m₃, compartidos, fijan un punto del espacio, así que los cuatro son sus imágenes. Eso son tres verificaciones: la restricción epipolar para m₁ y m₂, dos trilinealidades para las imágenes 1, 2 y 3, y dos para las imágenes 2, 3 y 4, o en la práctica la distancia de m₄ al punto transferido a partir de m₂ y m₃. Cuatro vistas tienen tensores propios, 81 números con 16 ecuaciones, pero la cadena los hace innecesarios. Cadenas así siguen un objeto detectado a lo largo de una secuencia de imágenes de rayos X, donde cada imagen nueva debe estar de acuerdo con las dos anteriores.

9Resumen

Las ideas clave junto a sus equivalentes en NumPy. En código, las cámaras son mat_a, mat_b y mat_c, o p para cualquiera de ellas, los píxeles son vectores homogéneos como m1 = np.array([x1, y1, 1]), y mat_q y mat_r son Q y r.