script7

clear
close all
load filelist

img1 = rgb2gray(imread(files{1}));
img2 = rgb2gray(imread(files{2}));
points1 = detectHarrisFeatures(img1);
points2 = detectHarrisFeatures(img2);

whos points1 points2
  Name            Size            Bytes  Class           Attributes

  points1      1421x1             17164  cornerPoints              
  points2      1291x1             15604  cornerPoints              

[features1,valid_points1] = extractFeatures(img1,points1);
[features2,valid_points2] = extractFeatures(img2,points2);
whos features1 features2 valid_points1 valid_points2
  Name                  Size            Bytes  Class             Attributes

  features1             1x1             91064  binaryFeatures              
  features2             1x1             81912  binaryFeatures              
  valid_points1      1421x1             17164  cornerPoints                
  valid_points2      1278x1             15448  cornerPoints                

[indexPairs,matchmetric] = matchFeatures(features1,features2);
whos indexPairs matchmetric
  Name              Size            Bytes  Class     Attributes

  indexPairs       45x2               360  uint32              
  matchmetric      45x1               180  single              

Retrieve the locations of the corresponding points for each image.

matchedPoints1 = valid_points1(indexPairs(:,1),:);
matchedPoints2 = valid_points2(indexPairs(:,2),:);
figure(1);
showMatchedFeatures(img1,img2,matchedPoints1,matchedPoints2,'Montage');
npts = size(matchedPoints1,1);
title([num2str(npts) ' Candidate point matches']);
legend('Matched points 1','Matched points 2');
Warning: Image is too big to fit on screen; displaying at 13% 
figure(2);
pts1 = double(matchedPoints1.Location);
show_points(img1,pts1);
pts2 = double(matchedPoints2.Location);
Warning: Image is too big to fit on screen; displaying at 25% 
[F, inliersIndex] = estimateFundamentalMatrix(pts1,pts2);
format shortg
F
format

dev = zeros(npts,1);
for k = 1:npts
    dev(k) = [pts2(k,1:2) 1]*F*[pts1(k,1:2) 1]';
end

ndx = 1:npts;
[ndx' inliersIndex dev]
F =

   2.2777e-07   1.5967e-06    -0.001249
  -1.7029e-06   8.7924e-08    0.0026363
   0.00031063   -0.0026733      0.99999


ans =

    1.0000         0   -0.0005
    2.0000    1.0000    0.0002
    3.0000         0   -0.0004
    4.0000    1.0000    0.0002
    5.0000    1.0000   -0.0002
    6.0000    1.0000    0.0002
    7.0000         0   -0.0003
    8.0000    1.0000   -0.0001
    9.0000         0   -0.0005
   10.0000         0    0.0008
   11.0000    1.0000   -0.0001
   12.0000    1.0000   -0.0002
   13.0000    1.0000    0.0003
   14.0000    1.0000    0.0003
   15.0000         0    0.0172
   16.0000    1.0000    0.0007
   17.0000    1.0000    0.0004
   18.0000         0    0.0019
   19.0000    1.0000    0.0005
   20.0000         0   -0.6873
   21.0000    1.0000   -0.0001
   22.0000    1.0000    0.0000
   23.0000    1.0000    0.0004
   24.0000    1.0000    0.0006
   25.0000    1.0000    0.0008
   26.0000    1.0000   -0.0002
   27.0000         0   -0.0006
   28.0000    1.0000    0.0001
   29.0000    1.0000    0.0009
   30.0000    1.0000    0.0001
   31.0000         0   -0.0004
   32.0000         0   -0.1759
   33.0000         0   -0.1260
   34.0000         0    1.4491
   35.0000         0   -0.0586
   36.0000         0   -0.0347
   37.0000    1.0000    0.0014
   38.0000         0    0.0083
   39.0000         0   -0.1145
   40.0000    1.0000    0.0029
   41.0000         0   -0.0337
   42.0000         0   -0.0108
   43.0000         0   -0.7671
   44.0000         0   -0.7269
   45.0000         0   -0.6820

kdx = find(inliersIndex>0);
[mx, mdx]  = sort(abs(dev(kdx)),'descend');
n = kdx(mdx(1));
fprintf('number of inliers %d\n',length(kdx));
fprintf('maximum inlier dev(%d) = %g\n',n,abs(dev(n)))
number of inliers 23
maximum inlier dev(40) = 0.0029295