Geometría epipolar

1Dos imágenes de un punto

Una sola imagen no puede decir a qué distancia está un punto. Un píxel m₁ fija el rayo en que está su punto del espacio, el rayo que sale del centro óptico C₁ y pasa por m₁, pero no en qué lugar de ese rayo. Toma ahora una segunda imagen desde otro lugar. Un punto M del espacio aparece en ambas, como m₁ en la primera y m₂ en la segunda. Dos píxeles así se llaman puntos correspondientes: imágenes del mismo punto del espacio.

Conociendo solo m₁, ¿dónde puede estar m₂? M está en algún lugar del rayo de m₁, y cada punto de ese rayo es un candidato. La segunda cámara ve cada candidato desde su propio centro óptico C₂. Una cámara lleva una recta del espacio a una recta de la imagen, así que las imágenes de todos los candidatos están sobre una recta de la segunda imagen. Esa recta es la línea epipolar de m₁, ℓ₂.

Esta es la restricción epipolar: si m₁ y m₂ son puntos correspondientes, entonces m₂ está sobre la línea epipolar de m₁. La restricción es necesaria, pero no suficiente. Todos los píxeles de la línea la cumplen, y solo uno de ellos es la imagen de M. Una tercera imagen puede señalar cuál, como muestran los apuntes de múltiples vistas. Aun así, convierte la búsqueda de m₂ en toda la segunda imagen en una búsqueda a lo largo de una línea. En una imagen de N × N píxeles, eso son N candidatos en vez de N².

2El plano epipolar y los epipolos

Los dos centros ópticos y M están en un mismo plano, el plano epipolar π. Ambos rayos están en él: el rayo de m₁ une C₁ y M, y el rayo de m₂ une C₂ y M. Cada plano de imagen corta a π en una recta, y esas dos rectas son las líneas epipolares. Con R₁ y R₂ para los dos planos de imagen,

Así, la restricción funciona en ambos sentidos. m₁ está sobre ℓ₁, la línea epipolar de m₂ en la primera imagen.

La recta que pasa por los dos centros ópticos es la línea base. Otro punto M′ del espacio tiene su propio plano epipolar π′, con otras líneas epipolares. Pero todo plano epipolar contiene a C₁ y a C₂, así que todos contienen la línea base. Por eso todas las líneas epipolares de una imagen pasan por un mismo punto, donde la línea base cruza su plano de imagen. Ese punto es el epipolo de la imagen:

La línea base pasa por ambos centros, así que e₁ es donde la primera cámara ve a C₂, y e₂ es donde la segunda cámara ve a C₁.

El ejemplo

El ejemplo de esta página es una serie de imágenes de rayos X, de 2688 × 2208 píxeles cada una, de un objeto que gira 2° en torno a un eje vertical entre una imagen y la siguiente. La fuente de rayos X hace el papel del centro óptico, como en los apuntes de transformaciones 3D. Vista desde el objeto, la fuente recorre un círculo horizontal. Las imágenes 1 y 82 están a 162° una de otra en ese círculo, así que sus fuentes casi se miran de frente. De ahí salen dos cosas.

• El plano del círculo contiene ambos centros ópticos, así que es un plano epipolar. Los dos epipolos están sobre la recta en que ese plano corta las imágenes, que es casi una fila de píxeles, en v ≈ 1285.
• Cada fuente está casi justo enfrente de la otra, así que cada imagen ve la otra fuente, cerca de su borde: los dos epipolos caen dentro de sus imágenes.

3La matriz fundamental a partir de las dos cámaras

Ahora en ecuaciones. Llama A y B a las matrices de proyección generales de las dos cámaras, como en los apuntes de calibración de cámaras. Un punto M = [X Y Z 1]^T, en las mismas coordenadas para ambas cámaras, aparece en

El objetivo es la línea epipolar ℓ₂ a partir de m₁, A y B solamente.

Dos puntos del rayo

ℓ₂ es la imagen del rayo de m₁, y una recta queda fijada por dos de sus puntos, así que bastan dos puntos del rayo. Uno es C₁. Todos los rayos de la primera cámara parten de ahí, así que A lo lleva al vector cero, AC₁ = 0, y se puede encontrar a partir de A solamente, como en los apuntes de calibración de cámaras.

El otro punto sale de m₁. Cualquier M con AM = m₁ se proyecta en m₁, así que está sobre el rayo. Tres ecuaciones no fijan cuatro incógnitas, pero la pseudo-inversa de A elige una solución:

AA^T es una matriz de 3 × 3, invertible porque A tiene rango 3. Entonces AA⁺ = AA^T(AA^T)⁻¹ = I, así que AM⁺ = m₁. El producto en el otro orden, A^TA, no serviría: es de 4 × 4 pero sigue teniendo rango 3, así que no tiene inversa.

M⁺ es un punto distinto de C₁, porque A lo lleva a m₁, no a 0. De todos los puntos del rayo, es aquel cuyo vector [X Y Z W]^T es más corto, lo que lo hace perpendicular a C₁ como vector de cuatro números. Esa elección es cosa del álgebra: cualquier otro punto del rayo daría la misma recta.

Sus imágenes en la segunda vista

B proyecta ambos puntos en la segunda imagen:

El primero es el epipolo, la imagen de C₁. El segundo es algún punto de ℓ₂.

La recta que pasa por ellos

La recta que pasa por dos puntos es su producto cruz, como en los apuntes de coordenadas homogéneas:

El producto cruz con un vector fijo u es lineal en el otro factor, así que es un producto de matriz por vector, u × v = [u]_×v, con la matriz antisimétrica

Multiplicarla por v = [v₁ v₂ v₃]^T da las tres componentes de u × v; la primera es u₂v₃ − u₃v₂. Con u = BC₁, la recta queda como una matriz por m₁:

La matriz F de 3 × 3 es la matriz fundamental. Depende solo de A y B, no de m₁, m₂ ni M: una sola matriz da la línea epipolar de cada píxel de la primera imagen.

La restricción epipolar

Un punto está sobre una recta cuando su producto punto es cero, ℓ^Tm = 0, como en los apuntes de coordenadas homogéneas. Así, m₂ está sobre ℓ₂ = Fm₁ exactamente cuando

Esa es la restricción epipolar en una sola ecuación, que cumple todo par de puntos correspondientes.

4Propiedades de F

Líneas epipolares en la primera imagen

Trasponer la restricción da m₁^TF^Tm₂ = 0, que dice que m₁ está sobre la recta F^Tm₂. Así,

F lleva puntos de la primera imagen a rectas de la segunda, y F^T va en sentido contrario. Intercambiar las cámaras, con B como la primera y A como la segunda, da F^T, salvo un factor de escala.

Escala

Para cualquier número k ≠ 0, kF da las mismas rectas y la misma restricción. Así, F es homogénea, como A y B.

Su determinante es 0

[u]_×u = u × u = 0, así que [u]_× lleva un vector distinto de 0 a 0, y es singular. F la tiene como factor, así que det F = det [BC₁]_× · det BA⁺ = 0. Una matriz homogénea de 3 × 3 tiene ocho números libres, y det F = 0 quita uno más: F tiene siete grados de libertad.

Los epipolos

Todas las líneas epipolares de la segunda imagen pasan por e₂, así que e₂^TFm₁ = 0 para todo m₁. Eso solo puede ocurrir si e₂^TF es la fila cero, y lo mismo vale para el otro epipolo:

Fe₁ = 0 también tiene una imagen. El rayo de e₁ es la propia línea base, que pasa por C₂, así que la segunda cámara ve todo el rayo como un solo punto, e₂, no como una recta. Y el vector cero no representa ninguna recta.

Estas ecuaciones dan los epipolos a partir de F solamente. F^Te₂ = 0 dice que e₂ es perpendicular a todas las columnas de F, así que es el producto cruz de dos columnas de F, dividido por su tercera componente para leerlo. Del mismo modo, e₁ es el producto cruz de dos filas de F.

Dos imágenes desde un mismo centro óptico

Si C₂ = C₁, entonces BC₁ = 0 y F = 0: no hay línea base, ni epipolo, ni línea epipolar. Ambas imágenes muestran entonces los mismos rayos, y se relacionan, en cambio, por una homografía, como en los apuntes de transformaciones 2D.

5La misma F a partir de un determinante

La misma F sale de las dos ecuaciones de proyección sin dibujar ninguna recta. Pasa todos sus términos a la izquierda, AM − λ₁m₁ = 0 y BM − λ₂m₂ = 0. Son seis ecuaciones con seis incógnitas, las cuatro componentes de M, λ₁ y λ₂, y forman un solo sistema Gv = 0:

Aquí G es de 6 × 6: A y B llenan sus cuatro primeras columnas, y m₁ y m₂ las dos últimas. Si m₁ y m₂ se corresponden, el punto M del que vienen da una solución v ≠ 0. Un sistema cuadrado con lado derecho 0 tiene una solución distinta de 0 solo cuando su determinante es 0, así que det G = 0.

Desarrolla det G por sus dos últimas columnas, las únicas que contienen los píxeles. Elegir la componente j de m₁ elimina la fila j de A, y elegir la componente i de m₂ elimina la fila i de B. Cada término es el producto de esas dos componentes por el determinante de 4 × 4 de las filas que quedan: dos filas de A y dos de B. Así, det G es lineal en m₁ y lineal en m₂, que es justo la forma m₂^TFm₁, con

Para las dos imágenes del ejemplo, esta F es igual a [BC₁]_×BA⁺ salvo un factor de escala, hasta once cifras significativas. Además muestra por qué la restricción tiene esta forma y no otra: det G es lineal en cada píxel, sea cual sea la posición de las cámaras.

6Los puntos medidos no caen en la línea

Los píxeles reales se miden y A y B se estiman, así que un par de puntos correspondientes casi nunca cumple m₂^TFm₁ = 0 exactamente: m₂ queda cerca de su línea epipolar, no sobre ella. Por eso la restricción se usa con una tolerancia. m₂ se acepta como posible correspondiente cuando su distancia a ℓ₂ = Fm₁ = [a b c]^T es menor que cierta d₀:

Por qué esta es la distancia: (a, b) es perpendicular a la recta ax + by + c = 0. Mover un punto una distancia t en la dirección del vector unitario (a, b)/√(a² + b²) cambia ax + by + c en t√(a² + b²), y el valor es 0 sobre la recta. Así, el valor en m₂, que es m₂^TFm₁, dividido por √(a² + b²), es la distancia de m₂ a la recta. Una recta escalada para que a² + b² = 1 da la distancia directamente, con un signo que indica el lado.

La tolerancia depende de qué tan bien se conocen los píxeles y las cámaras: un píxel o dos para esquinas encontradas por un programa, más para puntos marcados a mano. Con una tercera imagen, la prueba se vuelve mucho más estricta, como muestran los apuntes de múltiples vistas.

7Estimar F a partir de puntos correspondientes

Hasta aquí, F salió de A y B. Cuando las cámaras no están calibradas, F se puede estimar solo a partir de las imágenes, con pares de puntos que se sabe que se corresponden, como la misma esquina encontrada en ambas. Cada par m₁ = [x₁ y₁ 1]^T, m₂ = [x₂ y₂ 1]^T cumple m₂^TFm₁ = 0. Al multiplicar,

F_ij multiplica la componente i de m₂ y la componente j de m₁, y los píxeles son números conocidos, así que la ecuación es lineal en las nueve entradas de F. Apilándolas fila por fila en f = [F₁₁ F₁₂ … F₃₃]^T, cada par da una fila:

Ocho pares

F es homogénea, así que se puede fijar una entrada, como se hizo con a₃₄ en la calibración: F₃₃ = 1. El último término de cada fila pasa al lado derecho, y las ocho incógnitas que quedan, f′ = [F₁₁ … F₃₂]^T, cumplen

donde cada fila de Q es [x₁x₂ y₁x₂ x₂ x₁y₂ y₁y₂ y₂ x₁ y₁] y cada entrada de r es −1. Ocho pares dan un sistema cuadrado, que se resuelve como f′ = Q⁻¹r. Con más pares se ajusta por mínimos cuadrados, como en la calibración: f′ = (Q^TQ)⁻¹Q^Tr.

Fijar F₃₃ = 1 falla solo si el verdadero F₃₃ es 0. F₃₃ es m₂^TFm₁ para m₁ = m₂ = [0 0 1]^T, así que eso ocurre cuando el píxel (0, 0) de la segunda imagen está sobre la línea epipolar del píxel (0, 0) de la primera. La elección también se puede evitar, como con las homografías y la calibración: mantén las nueve incógnitas y resuelve Qf = 0 con ‖f‖ = 1, cuya respuesta es el vector singular derecho de Q con el menor valor singular.

Por qué ocho, si F tiene siete grados de libertad

El método lineal trata las ocho razones de F como independientes. No puede imponer además det F = 0, que no es una condición lineal en las entradas de F, así que necesita ocho ecuaciones. Y una estimación a partir de píxeles medidos no es exactamente singular.

Qué pares

Ocho pares no bastan si sus puntos del espacio están todos en un mismo plano. Los píxeles de un plano en dos imágenes se relacionan por una homografía, como en los apuntes de transformaciones 2D, y entonces Q tiene rango 6 a lo más: toda una familia de matrices cumple cada ecuación. Los puntos deben repartirse en profundidad, fuera de un único plano.

Qué tan sensible es

En píxeles, las columnas de Q tienen tamaños muy distintos: x₁x₂ llega a 7·10⁶ en imágenes de este tamaño, mientras x₁ no pasa de 3000. Las entradas de F difieren en el sentido contrario, y los errores pequeños en los píxeles se amplifican. Con exactamente ocho pares, la estimación hace pasar cada línea epipolar exactamente por su m₂, incluso uno erróneo, y las líneas de todos los demás puntos pagan el costo. Más pares, ajustados por mínimos cuadrados, permiten que el ajuste note el erróneo: su línea ya no lo sigue, y las demás se mueven menos.

En la figura de abajo, los pares son puntos del objeto, proyectados en las imágenes 1 y 82 y redondeados a píxeles enteros. El redondeo por sí solo deja las líneas de otros puntos del objeto a 1.10 píxeles de sus verdaderos correspondientes en promedio con ocho pares, y a 0.47 con dieciséis. Mover 30 píxeles hacia abajo el punto del par 3 en la imagen 82 deja clara la diferencia. Con ocho pares, su línea lo sigue exactamente, y las líneas de los demás puntos se desvían 22.38 píxeles en promedio. Con dieciséis, su línea queda a 21.12 píxeles de él, y las demás se desvían solo 2.90.

8Resumen

Las ideas clave junto a sus equivalentes en NumPy. En código, A, B y F son mat_a, mat_b y mat_f, y los píxeles son vectores homogéneos como m1 = np.array([x1, y1, 1]).