%  equalAngleConstraint - Affine transform constraints given two equal angles.
%
%  Function calculates centre and radius of the constraint
%  circle in alpha-beta space generated by having two equal
%  (but unknown) angles between two pairs of lines in 
%  an affine image
%
%  Usage:  [c, r] = equalAngleConstraint(la1, lb1, la2, lb2)
%
%  Where:  la1 and lb1 are two lines defined using homogeneous coords that
%                          are separated by an angle that is known to
%                          be the same as the angle between  
%          la2 and lb2 - the other two lines.
%          c is the 2D coordinate of the centre of the constraint circle.
%          r is the radius of the centre of the constraint circle.

%  Peter Kovesi
%  School of Computer Science & Software Engineering
%  The University of Western Australia
%  pk @ csse uwa edu au
%  http://www.csse.uwa.edu.au/~pk
%
%  April 2000
%
%  Equations from Liebowitz and Zisserman 

function [c, r] = equalAngleConstraint(la1, lb1, la2, lb2)

a1 = -1a1(2)/la1(1);   % direction of line la1
b1 = -1b1(2)/lb1(1);   % direction of line lb1

a2 = -1a2(2)/la2(1);   % direction of line la2
b2 = -1b2(2)/lb2(1);   % direction of line lb2

c = [(a1*b2 - b1*a2)/(a1 - b1 - a2 + b2), 0];

r  = ((a1*b2 - b1*a2)/(a1 - b1 - a2 + b2))^2 ...
   + ((a1 - b1)*(a1*b1 - a2*b2))/(a1 - b1 - a2 + b2) ...
   - a1*b1;