Result of a hybrid image.
The goal of this assignment is to understand a local feature matching algorithm using techniques.
The assignment contains several parts:
In this assignment, we are taught to use harris corner detector to choose interested points. The picture could be divided into 3 regions: flat, edge and corner, and it is corners that are most interested. Here we use second moment matrix to find corners.
In this assignment, the feature width is set to 16X16 , and divded to 16 4X4 cell which detects 8 orientations. So that each discriptor contains 4X4X8 dimentsion histogram that discribe the gradient of pixels around the interested point. Since the connected components will generate severl interested point, we need to find the local maxima in the connected components. In this assignment, function bwconnecomp is used to find connected components and generate one local maxima interested point.
In this assignment, distance between two features of different interested point pair is calculated to eliminate wrong pairing caused by having similar pattern. Therefore pair with least distance will be chosen.
[m,n,L] = size(image);
image = double(image);
g1 = fspecial('gaussian',[7,7], 1);
g2 = fspecial('gaussian',[8,8],0.1);
[dx,dy]=gradient(g1);
%%Build second moment matrix
for i = 1:L
Ix(:,:,i) = imfilter(image(:,:,i),dx);
Iy(:,:,i) = imfilter(image(:,:,i),dy);
Ix2(:,:,i) = imfilter(Ix(:,:,i).^2, g2);
Iy2(:,:,i) = imfilter(Iy(:,:,i).^2, g2);
Ixy(:,:,i) = imfilter(Ix(:,:,i).*Iy(:,:,i), g2);
end
%%Measure corner response
alpha = 0.05;
cim = (Ix2.*Iy2 - Ixy.^2) - alpha*(Ix2 + Iy2).^2;
thresh = 0.000001;
%%Transform to binary image
cimMark = cim>thresh;
list=[];
%%Use bwconnecomp function to generate local maxima
for k = 1:L
connect(k)=bwconncomp(cimMark(:,:,k),8);
cim2 = cim(:,:,k);
[height,length]=size(connect(k).PixelIdxList);
for i = 1:length
cell = cell2mat(connect(k).PixelIdxList(i));
a = -Inf;
b = 1;
for j = 1:size(cell)
if cim2(cell(j))>a
a = cim2(cell(j));
b = cell(j);
end
end
if ~any(list==b)
list = [list b];
end
end
end
y = ceil(list./m);
x = mod(list,m);
delMask = (ceil(x - (feature_width / 2)) < 1 | ceil(x + (feature_width / 2 - 1)) > m | ceil(y - feature_width / 2) < 1 | ceil(y + feature_width / 2 - 1) > n);
x(delMask) = [];
y(delMask) = [];
x = round(x);
y = round(y);
[M,N,L]=size(image);
g1 = fspecial('gaussian',[feature_width,feature_width],1);
[dx, dy] = gradient(g1);
grad_imX = imfilter(image,dx);
grad_imY = imfilter(image,dy);
n = length(x);
features = zeros(n,128);
for k = 1:n
his_int{k} = cell(4);
for i = 1:16
his_int{k}{i}=zeros(1,8);
end
%%Generate gradients
for l = 1:L
grad_interestX = grad_imX(round(x(k)-feature_width/2):round(x(k)+feature_width/2-1),(round(y(k)-feature_width/2):round(y(k)+feature_width/2-1)));
grad_interestY = grad_imY(round(x(k)-feature_width/2):round(x(k)+feature_width/2-1),(round(y(k)-feature_width/2):round(y(k)+feature_width/2-1)));
Mag = sqrt(grad_interestX.^2+grad_interestY.^2);
[m_grad,n_grad,L_grad]=size(grad_interestX);
%% Assign magnitude to different orientation according to the angle. For convinience, magnitude will be add to both orientation if the angle is located between them.
for i = 1:m_grad
for j = 1:n_grad
if grad_interestY(i,j)==0&&grad_interestX==0
bins = [1 5];
end
if grad_interestX(i,j)==0&&grad_interestY(i,j)>0
bins = [3];
end
if grad_interestX(i,j)==0&&grad_interestY(i,j)<0
bins = [7];
end
if grad_interestX(i,j)>0
if grad_interestY(i,j)/grad_interestX(i,j) ==0
bins = [1];
elseif grad_interestY(i,j)/grad_interestX(i,j)>0 && grad_interestY(i,j)/grad_interestX(i,j) <1
bins = [1 2];
elseif grad_interestY(i,j)/grad_interestX(i,j)==1
bins = [2];
elseif grad_interestY(i,j)/grad_interestX(i,j)>1
bins = [2 3];
elseif grad_interestY(i,j)/grad_interestX(i,j)<-1
bins = [7 8];
elseif grad_interestY(i,j)/grad_interestX(i,j)==-1
bins = [8];
elseif grad_interestY(i,j)/grad_interestX(i,j)>-1 && grad_interestY(i,j)/grad_interestX(i,j)<0
bins = [1 8];
end
end
if grad_interestX(i,j)<0
if grad_interestY(i,j)/grad_interestX(i,j) ==0
bins = [5];
elseif grad_interestY(i,j)/grad_interestX(i,j)>0 && grad_interestY(i,j)/grad_interestX(i,j) <1
bins = [5 6];
elseif grad_interestY(i,j)/grad_interestX(i,j)==1
bins = [6];
elseif grad_interestY(i,j)/grad_interestX(i,j)>1
bins = [6 7];
elseif grad_interestY(i,j)/grad_interestX(i,j)<-1
bins = [3 4];
elseif grad_interestY(i,j)/grad_interestX(i,j)==-1
bins = [4];
elseif grad_interestY(i,j)/grad_interestX(i,j)>-1 && grad_interestY(i,j)/grad_interestX(i,j)<0
bins = [4 5];
end
end
%%Assign magnitude
his_int{k}{ceil(i*4/feature_width),ceil(j*4/feature_width)}(bins) = his_int{k}{ceil(i*4/feature_width),ceil(j*4/feature_width)}(bins)+Mag(i,j);
end
end
end
for und = 1 : 16
features(k, (und - 1) * 8 + 1 : und * 8) = his_int{k}{und};
end
features(k, :) = features(k, :) .* (1 / norm(features(k, :)));
end
[num,area]=size(features1);
[num2,area2]=size(features2);
matches = zeros(num,2);
confidences = zeros(num,1);
for i = 1:num
dis = Inf;
dis_2 = Inf;
point1 = 1;
for j = 1:num2
L = sqrt(sum((features1(i,:)-features2(j,:)).^2));
if L < dis
point1 = j;
dis_2 = dis;
dis = L;
end
end
matches(i,1) = i;
matches(i,2) = point1;
confidences(i)=(dis_2/dis);
end
|
|
|