Sobre algunos problemas reales en visión por ordenador que conducen a sistemas de ecuaciones algebraicos

SOBRE ALGUNOS PROBLEMAS REALES EN VISION POR ORDENADOR QUE CONDUCEN A SISTEMAS DE ECUACIONES ALGEBRAICOS Luis Alvarez y Javier Sanchez Departamento de Informatica y Sistemas Universidad de Las Palmas de G.C. Campus Universitario de Tara 35017, Las Palmas Email: f lalvarez/jsanchez g @dis.ulpgc.es WebPage: http://serdis.dis.ulpgc.es/~lalvarez/ Abstract In this pap er, we present some real problems which appear in computer vision which yields to nonlinear system of algebraic equations. We study the problem of camera calibration. Roughly sp eaking, camera calibration consists in lo oking at the camera p osition in the 3 ; D world using as information the pro jection of a 3 ; D scene in a 2 ; D plane (the photogram). The problem is quite dierentwhenweuseasingleview or several views (stereo vision) of the 3 ; D scene. We will show in this pap er how these problems yields to nonlinear algebraic system of equations. Intro duccion En esta comunicacion, vamos a estudiar el problema de la calibracion de imagenes, haciendo esp ecial enfasis en los asp ectos geometricos y algebraicos. El ob jetivo de este traba jo es presentar algunos problemas que aparecen en la teora, y que dan lugar al estudio de sistemas de ecuaciones algebraicos no-lineales. Consideramos EACA-99, un congreso de algebra computacional y aplicaciones, un marco idoneo para presentar estos problemas, y estimular una discusion sobre ellos. La mayor parte de los resultados que se van a presentar, se encuentran en el excelente libro de Olivier Faugeras 3 ]. Presentaremos el problema de la calibracion cuando se utiliza una sola camara y cuando se utilizan dos camaras, como veremos, la forma de ab ordar el problema es bastante distinta,yenamb os casos, en ciertas situaciones, llegamos a sistemas de ecuaciones algebraicos. El mo delo que vamos adoptar para estudiar la calibracion de camaras es el denominado "pinhole mo del" que asume un mo delo proyectivo simple para el funcionamiento de la camara. Es decir, la camara viene representada por un plano R en 3 ; D y un fo co 1 Figure 1: Proyeccion de un punto 3-D. C =( C x C y C z ). La proyeccion de un punto M =( x y  z ) de la escena 3 ; D en la camara viene dada p or la interceccion de la recta MC con el plano R  (vease gura 1). 1 Calibracion de Camaras utilizando una sola vista. Calibrar una camara signica encontrar los parametros que determinan la proyeccion de la escena 3 ; D en el plano que determina la camara. Hay dos tip os fundamentales de parametros. Los parametros intrnsecos y los parametros extrnsecos. Una imagen captada por una camara digital viene dada p or una matriz de valores I ( i j )(niveles de gris). I ( i j ) representa la intensidad de luz recibida p or el pixel ( i j )delacamara. Los pixels representan las peque~nas celdas (o captores) en los que esta dividido el negativo de la camara. Los parametros intrnsecos de la camara son la distancia fo cal f el tama ~no del pixel ( p x p y )y la p osicion en la camara ( c x c y )de la interceccion de la recta p erp endicular al plano de la camara que pasa p or el fo co con la propia camara. Si p or ejemplo, el negativodelacamara estacentrado resp ecto a la p osicion del fo co, entonces c x = N  umer o de pixel s hor iz ontal es 2 c y = N  umer o de pixel s v er tical es 2 anivel proyectivo los parametros f y( p x p y ) estan relacionados entre s, es decir, si tomamos una constante k , la proyeccion generada por una camara de parametros f y ( p x p y ) es la misma que la generada por los parametros kf y ( kp x kp y ) : Por tanto, en lugar de los parametros f y( p x p y ) se utilizan los parametros  x = p x =f y  y = p y =f : Concluyendo, una camara tiene 4 parametros intr nsecos:  x  y c x c y : Como su nombre indica, los parametros 2 intrnsecos no dep enden de la p osicion 3 ; D que o cupa la camara en el universo 3 ; D: En muchas situaciones estos parametros se pueden cono cer "a priori" teniendo en cuenta las caractersticas tecnicas de la camara. Los parametros extrnsecos de la camara hacen referencia a la p osicion 3 ; D de la camara con resp ecto a la escena que estamos contemplando. Para ello se ja un sistema de referencia 3 ; D en la escena, y los parametros extrnsecos de la camara vienen dados p or un vector de traslacion t =( t x t y t z ) y una matriz de rotacion R =( r ij ) con resp ecto a dicho sistema de referencia. Dado que una rotacion dep ende de 3 parametros (el vector unitario sobre el que se rota y el angulo de rotacion), el n umero total de parametros extrnsecos es 6 : Para establecer la manera en que los parametros de la camara determinan la proyeccion de un punto M (3 ; D ) que se proyecta en un punto m (2 ; D ) resulta mas conveniente por diversos motivos traba jar con co ordernadas proyectivas. El resultado fundamental se establece en el siguiente teorema: Teorema 1.1. (Faugeras 3]) Sea un punto f M = ( M 1 M 2 M 3 M 4 ) 2 P 3 , entonces la proyeccion e m =( m 1 m 2 m 3 ) 2P 2 de f M en la camara que tiene por parametros intrnsecos  x  y c x c y ypor parametros extrnsecos t y R viene dada por la siguiente expresion: 0 @ m 1 m 2 m 3 1 A = 0 @  x r 1 + c x r 3  x t x + c x t z  y r 2 + c y r 3  y t y + c y t z r 3 t z 1 A 0 B B @ M 1 M 2 M 3 M 4 1 C C A (1) donde r i presenta el vector columna i ;  esimo de la matriz de rotacion R: A la matriz 4 x 3 que determina la proyeccion se le denomina e P es decir tenemos e m = e P f M: Notese que como aplicacion lineal proyectiva, la matriz e P esta denida modulo la multiplicacion por una constante. Para calibrar la camara sup ondremos que existe un conjunto de puntos f M i de los cuales cono cemos sus proyecciones en la camara e m i  el problema fundamental de la calibracion consiste en encontrar una matriz P de 4 x 3elementos que verique e m i = e P f M i 8 i =1  :: N (2) donde N es el numero de puntos. En la practica, para obtener el conjunto de puntos ( f M i  e m i )se utiliza, lo que se denomina, un calibrador. Un calibrador es un ob jeto como el representado en la imagen de la gura 2, del cual cono cemos con exactitud las co ordenadas 3 ; D de algunos puntos. Por ejemplo, en este caso, cono cemos (p orque lo hemos medido fsicamente), las co ordenadas de los puntos que corresp onden a las esquinas de los rectangulos negros. Por tanto los puntos f M i se elegiran entre las esquinas de los rectangulos, y los puntos e m i seran las co ordenadas de dicho punto en la imagen. Notese que la matriz e P posee 12 elementos, y que para cada par ( f M i  e m i ) la ecuacion (2) determina, una vez eliminado el parametro proyectivo, dos ecuaciones que involucran a 3 Figure 2: Imagen de un calibrador. los elementos de e P: que vienen dadas p or 3 X i =1 p 1 i M i + p 14 ; m 1  3 X i =1 p 3 i M i + p 34 ! =0 (3) 3 X i =1 p 2 i M i + p 14 ; m 2  3 X i =1 p 3 i M i + p 34 ! =0 Cuando N (el n umero de pares ( f M i  e m i ) cono cidos) es mayor o igual que 6, y dichos puntos estan en p osicion general (es decir no son ciertas conguraciones geometricas esp eciales de puntos) entonces es p osible recup erar la matriz e P a partir de las relaciones (2). La tecnica para calcular e P es muy sencilla. Si denotamos p =( p 11 p 12 :::: p 34 )elvector de 12 elementos dados p or las las de matriz e P entonces las relaciones (3) pueden escribirse como Ap =0  donde A es una matriz de 2 Nx 12 : Como p esta denido modulo una constante proyectiva, p o demos sup oner sin p erdida de generalidad que k p k =1 : El vector p puede interpretarse como el vector donde alcanza el mnimo la energa E ( p ) dada p or: E ( p )= k Ap k 2 = p T A T Ap con la restriccion k p k =1 : Un sencillo calculo nos lleva a que el mnimo de E ( p ) corresp onde al autovector del autov alor mas p eque ~no de la matriz A T A: Por otro lado, tambien por tecnicas standard una vez cono cida la matriz e P es posible calcular a partir de ella las conguraciones p osibles de parametros intrnsecos y extrnsecos. Cuando N es inferior a 6  p or ejemplo N =5  estan tecnicas standard no son op erativas, yno dan ning un resultado debido a que tendramos 10 ecuaciones lineales y 12 incognitas (los elementos de e P ). En 4 este caso, lo que sucede es que el espacio de autovectores de la matriz A T A aso ciado al autovalor mas p eque ~no tiene dimension mayor que 1  yportanto no p o demos determinar p de forma unica. Sin embargo, formalmente N =5 debera ser suciente para obtener alg un tip o de resultado p orque la matriz e P dep ende de 4 parametros intrnsecos y 6 parametros extrnsecos, es decir un total de 10 parametros, p or tanto, 10 ecuaciones seran en principio sucientes para encontrar p osibles valores de los parametros intrnsecos y extrnsecos Veamos como cuando N =5  p o demos llegar a un sistema algebraico de 11 ecuaciones y11 incognitas para determinar los parametros intrnsecos y extrnsecos. Para ello, vamos a expresar la matriz de rotation R algebraicamen te. Es bien cono cido que una matriz de rotacion puede expresarse como R = 0 @ s 2 + l 2 ; m 2 ; n 2 2( lm ; sn )2( nl + sm ) 2( lm + sn ) s 2 ; l 2 + m 2 ; n 2 2( mn ; sl ) 2( nl ; sm )2( mn + sl ) s 2 ; l 2 ; m 2 + n 2 1 A (4) con la condicion adicional de que s 2 + l 2 + m 2 + n 2 = 1 : Si sustituimos esta ecuacion en la expresion de e P dada en (1) obtenemos e P expresada algebraicamente en terminos de s l  m n t x t y t z c x c y  x y  y : (11 incognitas). Al sustituir los valores de P en las ecuaciones (3) obtenemos 2 ecuaciones algebraicas. Si N = 5 tenemos 10 ecuaciones algebraicas. La undecima ecuacion viene dada p or s 2 + l 2 + m 2 + n 2 =1 : Lo que determina un sistema algebraico de ecuaciones de 11 ecuaciones y 11 incognitas. Otra manera de formular el sistema de ecuaciones es el siguiente: Para que una matriz e P determine una proyeccion correctamente (es decir e P pueda expresarse como en (1) a partir de los parametros intrnsecos y extrnsecos) tiene que cumplir 2 condiciones (vease 3 ]) p 2 31 + p 2 32 + p 2 33 =1 0 @ 0 @ p 11 p 12 p 13 1 A ^ 0 @ p 31 p 32 p 33 1 A 1 A  0 @ 0 @ p 21 p 22 p 23 1 A ^ 0 @ p 31 p 32 p 33 1 A 1 A =0 donde ^ indica el pro ducto vectorial y  el pro ducto escalar. Por tanto, en el caso N =5  si a~nadimos estas 2 ecuaciones a las 10 ecuaciones lineales que determinan (3) obtenemos un sistema algebraico de 12 ecuaciones y 12 incognitas. 2 Calibracion de Camaras utilizando dos vistas. En el caso de tener 2 vistas sobre una escena, como se muestra en la gura 3, no es necesario utilizar un calibrador para obtener informacion 3 ; D sobre la escena, el par de imagenes stereo es suciente para recup erar la geometra 3 ; D . Sup ondremos en este caso que cono cemos los parametros intrnsecos de la camara, de hecho dichos parametros pueden normalizarse de tal forma que c x = c y =0 y  x =  y =1 : Vamos a centrar nuestro sistema de referencia 3 ; D en una de las imagenes. En este caso, entendemos p or calibrar las camaras, encontrar la rotacion R y traslacion t que existe 5 Figure 3: Ejemplo de un par de imagenes stereo. Figure 4: Proyeccion de un punto 3-D en dos planos. entre las p osiciones de las dos camaras. La informacion de la que partimos es un conjunto de N pares de puntos m i =( u i v i ), m 0 i =( u 0 i v 0 i ) que representan la proyeccion del mismo punto ( x i y i z i ) en las dos imagenes. Notese que para cualquier par de puntos m i y m 0 i las rectas Cm i y C 0 m 0 i estan en el mismo plano, esto establece una condicion, denominada condicion de epip olaridad, que relaciona m i m 0 i t y R: Un calculo directo nos llevaa queesta condicion puede expresarse como: T e m 0 i Q e m i =0 : (5) donde T e m 0 i =( u 0 i v 0 i  1)  T e m i =( u i v i  1), y Q viene dada p or la matriz Q = 0 @ ( r 13 t y ; r 12 t z ) ( r 11 t z ; r 13 t x ) ( r 12 t x ; r 11 t y ) ( r 23 t y ; r 22 t z ) ( r 21 t z ; r 23 t x ) ( r 22 t x ; r 21 t y ) ( r 33 t y ; r 32 t z ) ( r 31 t z ; r 33 t x ) ( r 32 t x ; r 31 t y ) 1 A (6) 6 El meto do standard de calibracion en este caso consiste en tomar N pares de puntos ( m i m 0 i )(con N  8), montar un sistema lineal de ecuaciones para encontrar los elementos de una matriz Q que verique (5) para to dos los pares de puntos ( m i m 0 i ) : Para ello se utiliza el denominado algoritmo de 8 puntos intro ducido en 6]. Este algoritmo se basa, como en el problema anterior de una sola camara, en el calculo del mnimo de una energa que se resuelve calculando el autovector aso ciado al autovalor mas p eque ~no de una cierta matriz. Una vez calculado Q se calculan las p osibles rotaciones R y traslaciones t que cumplen la relacion (6). Notese que la relacion (5) es homogenea en t ypor tanto Q solo puede determinarse modulo un factor de escala. Normalmente se normaliza k t k = 1 : Por tanto para recup erar de manera exacta la p osicion 3 ; D de los los puntos cuyas proyecciones son ( m i m 0 i ) es necesario alguna informacion adicional (p or ejemplo cono cer la distancia exacta entre dos puntos concretos 3 ; D de la escena ), si no p oseemos alguna informacion adicional de este tip o, entonces la informacion 3 ; D que recup eramos es modulo un factor de escala. Es decir la forma de los ob jetos es recono cible, p ero no su tama ~no exacto. Notese que para cualquier matriz Q de 3 x 3, no existe una rotacion R y traslacion t tal que Q pueda expresarse como (6). El resultado fundamentaleneste sentido es el siguiente: Teorema 2.1. (Faugeras 3]) Dada una matriz Q ,existe una rotacion R y una traslacion t tal que Q pueda expresarse como (6) si solo si j Q j =0 1 2 ; tr az a ; T QQ  2 ; tr az a ; T QQ  2 =0 Notese que con resp ecto a los elementos de Q ,la primera ecuacion es un p olinomio de grado 3  la segunda ecuacion es un p olinomio de grado 4 : El n umero de parametros libres necesarios para determinar la rotacion y traslacion es 5:3para la rotacion, y 2 para la traslacion (teniendo en cuenta la condicion k t k =1) : Por tanto resulta razonable p ensar que con 5 pares de puntos ( m i m 0 i )sera suciente para determinar p osibles valores para la rotacion R y traslacion t: Si N = 5 (n umero de pares ( m i m 0 i )) p o demos plantear un sistema algebraico de 8 ecuaciones y 9 incognitas dado p or las ecuaciones lineales (5) para los 5 pares de puntos ( m i m 0 i ) las dos ecuaciones dadas p or el teorema anterior, y la ecuacion adicional P ij q 2 ij =1 : A p esar de que hay menos ecuaciones que incognitas, el n umero de soluciones esp erado del sistema es nito. En una consulta p ersonal realizada al Profesor Faugeras sobre esta cuestion, este nos indico que el sistema tiene un n umero nito de soluciones debido a que la segunda ecuacion del teorema anterior determina implcitamente 2 condiciones sobre Q: Por ejemplo, una ecuacion del tip o ( x 2 + y 2 + z 2 ) 2 +( xy ; 1) 2 = 0 es una unica ecuacion que determina 2 condiciones sobre las variables x y  z : Este problema tambien fue ab ordado en 4] donde se presenta un algoritmo de 5 puntos para calcular las p osibles conguraciones de rotaciones y traslaciones. El algoritmo se basa en primer lugar en calcular los denominados epip olos (puntosdeinterceccion entre la recta CC 0 y los planos de las camaras), y a partir de ahi recup erar las rotaciones y traslaciones p osibles. En dicho articulo se demuestra que el n umero de conguraciones p osibles es 10 : 7 Otra formulacion alternativa del problema es utilizar la descrip cion de una rotacion, utilizando cuaterniones, dada por la expresion (4), sustituyendo dicha expresion en las ecuaciones (6) y (5) obtenemos la ecuacion: 0 B B B B B B B B B B B B @ u 0 i u i u 0 i v i u 0 i v 0 i u i v 0 i v i v 0 i u i v i 1 1 C C C C C C C C C C C C A T  0 B B B B B B B B B B B B @ t y (2 nl +2 sm ) ; t z (2 lm ; 2 sn ) t z ; s 2 + l 2 ; m 2 ; n 2  ; t x (2 nl +2 sm ) t x (2 lm ; 2 sn ) ; t y ; s 2 + l 2 ; m 2 ; n 2  t y (2 mn ; 2 sl ) ; t z ; s 2 ; l 2 + m 2 ; n 2  t z (2 lm +2 sn ) ; t x (2 mn ; 2 sl ) t x ; s 2 ; l 2 + m 2 ; n 2  ; t y (2 lm +2 sn ) t y ; s 2 ; l 2 ; m 2 + n 2  ; t z (2 mn +2 sl ) t z (2 nl ; 2 sm ) ; t x ; s 2 ; l 2 ; m 2 + n 2  t x (2 mn +2 sl ) ; t y (2 nl ; 2 sm ) 1 C C C C C C C C C C C C A =0 Si N =5  la expresion anterior determina 5 ecuaciones cuyas incognitas son 7 : s l  m n t x t y y t z  estas 5 ecuaciones se completan con las dos ecuaciones s 2 + l 2 + m 2 + n 2 =1 t 2 x + t 2 y + t 2 z =1 construyendo as un sistema algebraico de 7 ecuaciones y 7 incognitas. Este problema ha sido resuelto utilizando tecnicas algebraicas por Emeris (vease 1] y 2]) simplicando el sistema de la siguiente forma: Se realiza un cambio de variable utilizando cuaterniones. Los cuaterniones son elementos del espacio IR 4  que se expresan como _ s =( s 0  s ) donde s 0 es un numero real, y s 2 IR 3  y donde se ha denido una ley de comp osicion interna de la siguiente forma: ( s 0  s )  ( d 0  d )=( s 0 d 0 ; sd T s 0 d + d 0 s + s ^ d ) donde s ^ d dene el pro ducto vectorial usual en IR 3 : El espacio IR 3 se sumerge en en conjunto de los cuaterniones, como los elementos ( s 0  s ) tales que s 0 =0 : De manera similar alos numeros complejos, se dene el conjugado de un cuaterni on como ( s 0  s )=( s 0  ; s )  cumpliendose la propiedad ( s 0  s )  ( s 0  s )=( s 2 0 + k s k 2  0 ) : Volviendo a nuestro problema, vamos a denotar p or _ q =( s l  m n ) el cuaternion que determina la matriz de rotacion R y por _ t =(0 t x t y t z ) el cuaternion que determina el vector de traslacion t: Vamos a sustituir ahora _ t ,por el cuaternion _ d denido como _ d = _ t  _ q: De tal manera que ahora las nuevas incognitas son _ d y _ q: Realizando ahora algunas manipulaciones algebraicas, y sup oniendo que q 0 = d 0 = 1 se obtiene que la relacion de epip olaridad T e m 0 i Q e m i = 0  en las nuevas incognitas se escrib e como ; e m T i q ; e m 0 T i d  + e m 0 T i e m i +( e m i ^ q ) T e m 0 i +( e m i ^ q ) ; d ^ e m 0 i  + e m T i ; d ^ e m 0 i  =0 para i =1  :: 5 : La gran venta ja de estas ecuaciones es que son bilineales resp ecto a e m i y e m 0 i : Ademas a estas 5 ecuaciones tenemos que a ~nadir una ecuacion mas que determina que 8 cuando recup eremos el cuaternion _ t ,este verique que t 0 =0  ello determina la condicion adicional 1 ; dt T =0 a partir de este sistema de ecuaciones, Emeris desarrolla un meto do puramente algebraico que determina de manera exacta el conjunto de soluciones p osibles utilizando tecnicas en el contexto de la teora de variedades toricas y la eliminacion disp ersa. En el caso de tener 3 o mas imagenes de una escena, tambien aparecen de forma natural condiciones de tip o algebraico determinadas por la geometra epip olar de las diferentes camaras. Un estudio sobre estas condiciones ha sido desarrollado p or Faugeras y Mourrain en 5 ]. Existen otros problemas interesantes en Robotica, como la manipulacion de un brazo articulado de rob ot con 6 grados de lib ertad (vease 7 ]), que conducen de forma natural ala resolucion de sistemas de ecuaciones algebraicos. En este caso, el problema consiste en dada una p osicion determinada en el espacio, encontrar las p osibles conguraciones de los elementos del brazo articulado del rob ot, para que el rob ot alcance dicha p osicion en el espacio. Es un problema que surge de forma natural en las cadenas de monta je de las industrias que utilizan rob ots articulados. Agradecimientos Este traba jo ha sido parcialmente nanciado p or el proyecto PIB95-1225 del M.E.C. y p or la accion integrada Hispano-Francesa HF98-0098. Agradecemos al Comite Cientco del EACA'99, p or p oner en nuestro cono cimiento la existencia de los traba jos 1], 2 ], 7] y 5 ], que nosotros descono camos, y cuya inclusion ha supuesto una sensible mejora en la presentaci on de algunos de los resultados expuestos en este traba jo. References 1] I. Emiris "Sparse Elimination and Applications in Kinematics", Ph.D. Thesis U.C. Berkeley , diciembre 1994. 2] I. Emiris "A general Solver Based on Sparse Resultants: Numerical Issues & Kinematic Applications", Informe Tecnico N o 3110, INRIA , 1997. 3] O.Faugeras, \3-D computer vision. A geometric viewp oint," MIT Press , 1993. 4] O.Faugeras and S.Maybank \Motion from p oint matches: multiplicity of solutions," International Journal of Computer Vision ,Vol. 4(3) pp 225-246, 1990. 5] O.Faugeras and B.Mourrain \On the geometry and algebra of the point and line corresp ondences b etween N images" Informe Tecnico N o 2665, INRIA , 1995. 9