6 Summary of the algorithm

Given the high diffusion of stereo in research and applications, we have endeavoured to make our algorithm as easily reproducible and usable as possible. To this purpose, we give below our working MATLAB code, which is simple and compact (22 lines). The comments provided should make it understandable without knowledge of MATLAB. Moreover, a ``rectification kit'' including code, data sets, and instructions on how to use the algorithm can be downloaded from

  ftp://taras.dimi.uniud.it/pub/sources/rectif_m.tar.gz

   function [T1,T2,Pn1,Pn2] = rectify(Po1,Po2)
   % RECTIFY computes the rectified projection matrices
   % Pn1, Pn2, and the rectifying transformation of
   % the retinal plane T1 and T2 (in homog. coords.)
   % Po1 and Po2 are the original projection  matrices.

   % focal length (extp(a,b) is external product of vectors a,b)
   au = norm(extp(Po1(1,1:3)', Po1(3,1:3)'));
   av = norm(extp(Po1(2,1:3)', Po1(3,1:3)'));

   % optical centers
   c1 = - inv(Po1(:,1:3))*Po1(:,4);
   c2 = - inv(Po2(:,1:3))*Po2(:,4);

   % unit vectors of retinal planes 
   fl = Po1(3,1:3)';   fr = Po2(3,1:3)';
   nn = extp(fl,fr);

   % select the first 3 components
   L = [1 1 1 0 ]';

   % lagra(A,L,k) computes the minimum of the norm of Ax,
   % with the constraint "squared norm of x(L) equal to k"

   A = [ [c1' 1]' [c2' 1]' [nn'  0]' ]'
   a3 = lagra(A,L,1) ;

   A = [ [c1' 1]' [c2' 1]' [a3(1:3)' 0]' ]';
   a2 = lagra(A,L,av^2-1) ;

   A = [ [c1' 1]' [a2(1:3)' 0]' [a3(1:3)' 0]' ]';
   a1 = lagra(A,L,au^2-1) ;

   A = [ [c2' 1]' [a2(1:3)' 0]' [a3(1:3)' 0]' ]';
   b1 = lagra(A,L,au^2-1);

   % adjustment transformation
   T = [
     1 0 0
     0 1 0
     0 0 1 ];

   % rectifying projection matrices
   Pn1 = T * [ a1 a2 a3 ]';     Pn2 = T * [ b1 a2 a3 ]';

   % rectifying image transformations
   T1 = Pn1(1:3,1:3)* inv(Po1(1:3,1:3));
   T2 = Pn2(1:3,1:3)* inv(Po2(1:3,1:3));

This MATLAB code receives in input two projection matrices, and returns two image transformation matrices (to be applied to image points) and two rectifying projection matrices (to be used to perform reconstruction from rectified images).

Reconstruction.

Reconstruction is performed through Eq. (1). Given two conjugate points, and , and the two projection matrices and , we can write the overconstrained linear system

where

and gives the position of the 3-D point projected to the conjugate points. Notice that different sizes or centers can be obtained, if required, by pre-multiplying both rectifying projection matrices by a suitable matrix , set to identity in our MATLAB code, of the same structure as (3).