function LindHallCompositeHomography
im1 = imread('lindhall_left.png');
im2 = imread('lindhall_right.png');

f = 1224;
px = size(im1,2)/2;
py = size(im1,1)/2;

K = [f 0 px; 0 f py; 0 0 1];

im_warped = zeros(2000,4000,3);

u2 = [1275 1095;
      1291 812;
      400 666;
      359 1054];

u1 = [3564 1205;
      3525 817;
      2624 896;
      2629 1184];
  
H = ComputeHomography(u1,u2);
H = H/H(3,3);
R = inv(K) * H * K;
detR = det(R)
R = 1/detR^(1/3) * R
det(R)

v = [2200 1080;
     2200 920;
     1800 920;
     1800 1080];

u1 = [u1 ones(4,1)];
u2 = [u2 ones(4,1)];
v = [v ones(4,1)];





figure(5)
clf;
plot3([0 0], [0 0], [0 1], 'b-');
hold on
plot3([0 0], [0 1], [0 0], 'g-');
hold on
plot3([0 R(3,1)], [0 R(3,2)], [0 R(3,3)], 'b-');
hold on
plot3([0 R(2,1)], [0 R(2,2)], [0 R(2,3)], 'g-');
axis equal
grid on
% return;

H1 = ComputeHomography(v, u1);
H2 = ComputeHomography(v, u2);

im_warped1 = ImageWarpingGivenTargetImage(im1, H1, im_warped);
im_warped2 = ImageWarpingGivenTargetImage(im2, H2, im_warped);

im = 0.5*im_warped1 + 0.5*im_warped2;

figure(1)
clf;
imshow(im);

imwrite(im, 'lindhall_composite_homography.png');

% figure(2)
% clf;
% imshow(im1);
% hold on
% plot(u1(:,1), u1(:,2), 'rx');
% 
% figure(3)
% clf;
% imshow(im2);
% hold on
% plot(u2(:,1), u2(:,2), 'rx');
return;
 


