% HOMOGRAPHY2D - computes 2D homography % % Usage: H = homography2d(x1, x2) % H = homography2d(x) % % Arguments: % x1 - 3xN set of homogeneous points % x2 - 3xN set of homogeneous points such that x1<->x2 % % x - If a single argument is supplied it is assumed that it % is in the form x = [x1; x2] % Returns: % H - the 3x3 homography such that x2 = H*x1 % % This code follows the normalised direct linear transformation % algorithm given by Hartley and Zisserman "Multiple View Geometry in % Computer Vision" p92. % Copyright (c) 2003-2005 Peter Kovesi % School of Computer Science & Software Engineering % The University of Western Australia % pk at csse uwa edu au % http://www.csse.uwa.edu.au/~pk % % Permission is hereby granted, free of charge, to any person obtaining a copy % of this software and associated documentation files (the "Software"), to deal % in the Software without restriction, subject to the following conditions: % % The above copyright notice and this permission notice shall be included in % all copies or substantial portions of the Software. % % The Software is provided "as is", without warranty of any kind. % May 2003 - Original version. % Feb 2004 - Single argument allowed for to enable use with RANSAC. % Feb 2005 - SVD changed to 'Economy' decomposition (thanks to Paul O'Leary) function H = homography2d(varargin) [x1, x2] = checkargs(varargin(:)); % Attempt to normalise each set of points so that the origin % is at centroid and mean distance from origin is sqrt(2). [x1, T1] = normalise2dpts(x1); [x2, T2] = normalise2dpts(x2); % Note that it may have not been possible to normalise % the points if one was at infinity so the following does not % assume that scale parameter w = 1. Npts = length(x1); A = zeros(3*Npts,9); O = [0 0 0]; for n = 1:Npts X = x1(:,n)'; x = x2(1,n); y = x2(2,n); w = x2(3,n); A(3*n-2,:) = [ O -w*X y*X]; A(3*n-1,:) = [ w*X O -x*X]; A(3*n ,:) = [-y*X x*X O ]; end [U,D,V] = svd(A,0); % 'Economy' decomposition for speed % Extract homography H = reshape(V(:,9),3,3)'; % Denormalise H = T2\H*T1; %-------------------------------------------------------------------------- % Function to check argument values and set defaults function [x1, x2] = checkargs(arg); if length(arg) == 2 x1 = arg{1}; x2 = arg{2}; if ~all(size(x1)==size(x2)) error('x1 and x2 must have the same size'); elseif size(x1,1) ~= 3 error('x1 and x2 must be 3xN'); end elseif length(arg) == 1 if size(arg{1},1) ~= 6 error('Single argument x must be 6xN'); else x1 = arg{1}(1:3,:); x2 = arg{1}(4:6,:); end else error('Wrong number of arguments supplied'); end