function sift_matching()
% Read the images
Ia = single(rgb2gray(imread('cathedral1.bmp')));
Ib = single(rgb2gray(imread('cathedral2.bmp')));
% Use SIFT to compute feature matches
[fa, da] = vl_sift(Ia);
[fb, db] = vl_sift(Ib);
[matches, scores] = vl_ubcmatch(da, db);
% Sort by the scores to get the best 100 matches, if there are that many
[sortedScores, indices] = sort(scores);
sortedMatches = matches(:, indices);

ind1 = 0;   ind2 = 0;
imshow(appendimages(Ia,Ib), []); hold on;
for i = 1 : min(100, length(scores))
    ind1 = sortedMatches(1, i);
    ind2 = sortedMatches(2, i);
    plot(fa(1, ind1), fa(2, ind1), 'ro'); hold on;
    plot(size(Ia,2)+fb(1, ind2), fb(2, ind2), 'ro'); hold on;
    plot([fa(1, ind1) size(Ia,2)+fb(1, ind2)], ...
         [fa(2, ind1) fb(2, ind2)], 'g-');
end

end