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 ).