一、概述
在图像处理相关的问题中,图像对准是一类典型的问题,也就是要将两幅图严丝合缝地对应起来。通常来讲,两幅图大小不一,一个是模板,一个是母图,也就是要在母图中搜寻定位到与模板图最为接近的区域。
实现的方式有很多,惯常使用的是基于准则匹配的方法和基于特征匹配的方法。基于准则匹配,就是直接地对图的灰度值矩阵进行计算操作,以特定的准则遍历整个母图,找到与目标图(模板图)最相近的子区域;基于特征匹配,就是先提取出图像特征,再基于特征进行操作。这里对基于准则匹配的图像对准基本方法做简单介绍。
二、匹配准则
常见的匹配准则有SAD、MAD、SSD、MSD、NCC。前四种是基于两个矩阵的向量差做运算,NCC是计算两个矩阵的相关系数。事实上,矩阵是一个高阶向量(二阶张量),对两个矩阵向量作差,就得到差向量,对差向量做分析运算,便可在一定程度上获得两个矩阵间的差异性信息。
设
A
=
(
a
i
j
)
A=\left( a_{ij} \right)
A=(aij),
B
=
(
b
i
j
)
B=\left( b_{ij} \right)
B=(bij),
i
=
1
,
2
,
.
.
.
,
M
i=1,2,...,M
i=1,2,...,M,
j
=
1
,
2
,
.
.
.
,
N
j=1,2,...,N
j=1,2,...,N.则差向量
D
=
A
−
B
=
(
a
i
j
−
b
i
j
)
\begin{equation} \begin{aligned} D&=A-B\\&=\left( a_{ij}-b_{ij} \right) \end{aligned} \end{equation}
D=A−B=(aij−bij)
(1) SAD
SAD,绝对误差算法(Sum of Absolute Differences),它是差向量D中各元素的绝对值之和,也就是L1范数,是两个向量间的曼哈顿距离。表达式为 S A D = ∑ i = 1 M ∑ j = 1 N ∣ a i j − b i j ∣ \begin{equation} \begin{aligned}SAD=\sum_{i=1}^{M}{\sum_{j=1}^{N}{\left| a_{ij}-b_{ij} \right|}}\end{aligned} \end{equation} SAD=i=1∑Mj=1∑N∣aij−bij∣
(2) MAD
MAD,平均绝对误差算法(Mean Absolute Differences),它是在SAD基础上进一步求平均值。表达式为 M A D = 1 M × N ∑ i = 1 M ∑ j = 1 N ∣ a i j − b i j ∣ \begin{equation} \begin{aligned}MAD=\frac{1}{M\times N}\sum_{i=1}^{M}{\sum_{j=1}^{N}{\left| a_{ij}-b_{ij} \right|}}\end{aligned} \end{equation} MAD=M×N1i=1∑Mj=1∑N∣aij−bij∣
(3) SSD
SSD,误差平方和算法(Sum of Squared Differences),它是差向量D中各元素的平方和。表达式为 S S D = ∑ i = 1 M ∑ j = 1 N ( a i j − b i j ) 2 \begin{equation} \begin{aligned}SSD=\sum_{i=1}^{M}{\sum_{j=1}^{N}{\left( a_{ij}-b_{ij} \right)^{2}}}\end{aligned} \end{equation} SSD=i=1∑Mj=1∑N(aij−bij)2
(4) MSD
MSD,平均误差平方和算法(Mean Square Differences),它是在SSD的基础上进一步求平均值。表达式为 M S D = 1 M × N ∑ i = 1 M ∑ j = 1 N ( a i j − b i j ) 2 \begin{equation} \begin{aligned}MSD=\frac{1}{M\times N}\sum_{i=1}^{M}{\sum_{j=1}^{N}{\left( a_{ij}-b_{ij} \right)^{2}}}\end{aligned} \end{equation} MSD=M×N1i=1∑Mj=1∑N(aij−bij)2
(5) NCC
NCC,归一化互相关算法(Normalized Cross Correlation)。若将两个矩阵看做两个随机变量,那么NCC就是两个变量之间的皮尔逊相关系数。同时,它也是两个矩阵向量在各自中心化之后彼此间空间夹角的余弦值。它的表达式为
N
C
C
=
∑
i
=
1
M
∑
j
=
1
N
(
a
i
j
−
E
(
A
)
)
(
b
i
j
−
E
(
B
)
)
∑
i
=
1
M
∑
j
=
1
N
(
a
i
j
−
E
(
A
)
)
2
⋅
∑
i
=
1
M
∑
j
=
1
N
(
b
i
j
−
E
(
B
)
)
2
\begin{equation} \begin{aligned}NCC=\frac{\sum_{i=1}^{M}{\sum_{j=1}^{N}{\left( a_{ij}-E(A) \right)\left( b_{ij}-E(B) \right)}}}{\sqrt{\sum_{i=1}^{M}{\sum_{j=1}^{N}{\left( a_{ij}-E(A) \right)^{2}}}}\cdot\sqrt{\sum_{i=1}^{M}{\sum_{j=1}^{N}{\left( b_{ij}-E(B) \right)^{2}}}}}\end{aligned} \end{equation}
NCC=∑i=1M∑j=1N(aij−E(A))2⋅∑i=1M∑j=1N(bij−E(B))2∑i=1M∑j=1N(aij−E(A))(bij−E(B))
易知,ncc值的范围为
[
−
1
,
1
]
[−1,1]
[−1,1],越接近1,两个矩阵越相关;越接近-1,两个矩阵越不相关。
等同于皮尔逊相关系数
皮尔逊相关系数,用以衡量两个变量间的线性相关性。它的表达式为
P
e
a
r
s
o
n
=
C
o
v
(
X
,
Y
)
D
(
X
)
⋅
D
(
Y
)
=
E
(
X
−
E
X
)
(
Y
−
E
Y
)
D
(
X
)
⋅
D
(
Y
)
\begin{equation} \begin{aligned} Pearson&=\frac{Cov\left( X,Y \right)}{\sqrt{D\left( X \right)}\cdot\sqrt{D\left( Y \right)}}\\&=\frac{E\left( X-EX \right)\left( Y-EY \right)}{\sqrt{D\left( X \right)}\cdot\sqrt{D\left( Y \right)}} \end{aligned} \end{equation}
Pearson=D(X)⋅D(Y)Cov(X,Y)=D(X)⋅D(Y)E(X−EX)(Y−EY) 将两个矩阵看做两个随机变量代入,有
P
e
a
r
s
o
n
=
1
M
×
N
∑
i
=
1
M
∑
j
=
1
N
(
a
i
j
−
E
(
A
)
)
(
b
i
j
−
E
(
B
)
)
∑
i
=
1
M
∑
j
=
1
N
(
a
i
j
−
E
(
A
)
)
2
M
×
N
⋅
∑
i
=
1
M
∑
j
=
1
N
(
b
i
j
−
E
(
B
)
)
2
M
×
N
=
∑
i
=
1
M
∑
j
=
1
N
(
a
i
j
−
E
(
A
)
)
(
b
i
j
−
E
(
B
)
)
∑
i
=
1
M
∑
j
=
1
N
(
a
i
j
−
E
(
A
)
)
2
∑
i
=
1
M
∑
j
=
1
N
(
b
i
j
−
E
(
B
)
)
2
=
N
C
C
\begin{equation} \begin{aligned} Pearson&=\frac{\frac{1}{M\times N}\sum_{i=1}^{M}{\sum_{j=1}^{N}{\left( a_{ij}-E(A) \right)\left( b_{ij}-E(B) \right)}}}{\sqrt{\frac{\sum_{i=1}^{M}{\sum_{j=1}^{N}{\left( a_{ij}-E(A) \right)^{2}}}}{M\times N}}\cdot\sqrt{\frac{\sum_{i=1}^{M}{\sum_{j=1}^{N}{\left( b_{ij}-E(B) \right)^{2}}}}{M\times N}}}\\&=\frac{\sum_{i=1}^{M}{\sum_{j=1}^{N}{\left( a_{ij}-E(A) \right)\left( b_{ij}-E(B) \right)}}}{\sqrt{\sum_{i=1}^{M}{\sum_{j=1}^{N}{\left( a_{ij}-E(A) \right)^{2}}}} \sqrt{\sum_{i=1}^{M}{\sum_{j=1}^{N}{\left( b_{ij}-E(B) \right)^{2}}}}}\\&=NCC \end{aligned} \end{equation}
Pearson=M×N∑i=1M∑j=1N(aij−E(A))2⋅M×N∑i=1M∑j=1N(bij−E(B))2M×N1∑i=1M∑j=1N(aij−E(A))(bij−E(B))=∑i=1M∑j=1N(aij−E(A))2∑i=1M∑j=1N(bij−E(B))2∑i=1M∑j=1N(aij−E(A))(bij−E(B))=NCC
等同于余弦距离
余弦距离即空间向量夹角的余弦值,通常用以衡量两个向量间的差异度。它的表达式为
c
o
s
θ
=
<
X
,
Y
>
∣
X
∣
⋅
∣
Y
∣
\begin{equation} \begin{aligned}cos\theta=\frac{<X,Y>}{\left| X \right|\cdot\left| Y \right|}\end{aligned} \end{equation}
cosθ=∣X∣⋅∣Y∣<X,Y> 将两个矩阵向量去中心化后代入,有
c
o
s
θ
=
<
A
−
E
(
A
)
,
B
−
E
(
B
)
>
∣
A
−
E
(
A
)
∣
⋅
∣
B
−
E
(
B
)
∣
=
∑
i
=
1
M
∑
j
=
1
N
(
a
i
j
−
E
(
A
)
)
(
b
i
j
−
E
(
B
)
)
∑
i
=
1
M
∑
j
=
1
N
(
a
i
j
−
E
(
A
)
)
2
∑
i
=
1
M
∑
j
=
1
N
(
b
i
j
−
E
(
B
)
)
2
=
N
C
C
\begin{equation} \begin{aligned} cos\theta&=\frac{<A-E(A),B-E(B)>}{\left| A-E(A) \right|\cdot\left| B-E(B) \right|}\\&=\frac{\sum_{i=1}^{M}{\sum_{j=1}^{N}{\left( a_{ij}-E(A) \right)\left( b_{ij}-E(B) \right)}}}{\sqrt{\sum_{i=1}^{M}{\sum_{j=1}^{N}{\left( a_{ij}-E(A) \right)^{2}}}}\sqrt{\sum_{i=1}^{M}{\sum_{j=1}^{N}{\left( b_{ij}-E(B) \right)^{2}}}}}\\&=NCC \end{aligned} \end{equation}
cosθ=∣A−E(A)∣⋅∣B−E(B)∣<A−E(A),B−E(B)>=∑i=1M∑j=1N(aij−E(A))2∑i=1M∑j=1N(bij−E(B))2∑i=1M∑j=1N(aij−E(A))(bij−E(B))=NCC
三、matlab实现
(1) SAD
clear all;
close all; clc;
%1.读取图片
img_A_dir = '.\data\lena.bmp'; %待寻母图
img_A_raw = imread(img_A_dir);
[r1,c1,d1] = size(img_A_raw);
if d1==3 %灰度化
img_A = rgb2gray(img_A_raw);
else
img_A = img_A_raw;
end
img_B_dir = '.\data\refer.bmp'; %模板图
img_B_raw = imread(img_B_dir);
[r2,c2,d2] = size(img_B_raw);
if d2==3
img_B = rgb2gray(img_B_raw);
else
img_B = img_B_raw;
end
%2.计算SAD矩阵
msad = zeros(r1-r2,c1-c2);
for i = 1:r1-r2
for j = 1:c1-c2
temp = img_A(i:i+r2-1,j:j+c2-1);
msad(i,j) = msad(i,j) + sum(sum(abs(temp - img_B)));
end
end
%3.定位匹配位置
min_sad = min(min(msad));
[x,y] = find(msad == min_sad);
x = x(1); %定位到的第一个位置
y = y(1);
%4.保存结果图
getImg = img_A_raw(x:x+r2-1,y:y+c2-1,1:3);
imwrite(getImg,'.\output\SAD_match.bmp');
fprintf('\n Done. \n');
在这里插入代码片
(2) MAD
clear all;
close all; clc;
%1.读取图片
img_A_dir = '.\data\lena.bmp'; %待寻母图
img_A_raw = imread(img_A_dir);
[r1,c1,d1] = size(img_A_raw);
if d1==3 %灰度化
img_A = rgb2gray(img_A_raw);
else
img_A = img_A_raw;
end
img_B_dir = '.\data\refer.bmp'; %模板图
img_B_raw = imread(img_B_dir);
[r2,c2,d2] = size(img_B_raw);
if d2==3
img_B = rgb2gray(img_B_raw);
else
img_B = img_B_raw;
end
%2.计算MAD矩阵
mmad = zeros(r1-r2,c1-c2);
for i = 1:r1-r2
for j = 1:c1-c2
temp = img_A(i:i+r2-1,j:j+c2-1);
mmad(i,j) = mmad(i,j) + sum(sum(abs(temp - img_B)))/(r2*c2);
end
end
%3.定位匹配位置
min_mad = min(min(mmad));
[x,y] = find(mmad == min_mad);
x = x(1); %定位到的第一个位置
y = y(1);
%4.保存结果图
getImg = img_A_raw(x:x+r2-1,y:y+c2-1,1:3);
imwrite(getImg,'.\output\MAD_match.bmp');
fprintf('\n Done. \n');
(3) SSD
clear all;
close all; clc;
%1.读取图片
img_A_dir = '.\data\lena.bmp'; %待寻母图
img_A_raw = imread(img_A_dir);
[r1,c1,d1] = size(img_A_raw);
if d1==3 %灰度化
img_A = rgb2gray(img_A_raw);
else
img_A = img_A_raw;
end
img_B_dir = '.\data\refer.bmp'; %模板图
img_B_raw = imread(img_B_dir);
[r2,c2,d2] = size(img_B_raw);
if d2==3
img_B = rgb2gray(img_B_raw);
else
img_B = img_B_raw;
end
%2.计算SSD矩阵
mssd = zeros(r1-r2,c1-c2);
for i = 1:r1-r2
for j = 1:c1-c2
temp = img_A(i:i+r2-1,j:j+c2-1);
mssd(i,j) = mssd(i,j) + sum(sum((temp - img_B).^2));
end
end
%3.定位匹配位置
min_ssd = min(min(mssd));
[x,y] = find(mssd == min_ssd);
x = x(1); %定位到的第一个位置
y = y(1);
%4.保存结果图
getImg = img_A_raw(x:x+r2-1,y:y+c2-1,1:3);
imwrite(getImg,'.\output\SSD_match.bmp');
fprintf('\n Done. \n');
(4) MSD
clear all;
close all; clc;
%1.读取图片
img_A_dir = '.\data\lena.bmp'; %待寻母图
img_A_raw = imread(img_A_dir);
[r1,c1,d1] = size(img_A_raw);
if d1==3 %灰度化
img_A = rgb2gray(img_A_raw);
else
img_A = img_A_raw;
end
img_B_dir = '.\data\refer.bmp'; %模板图
img_B_raw = imread(img_B_dir);
[r2,c2,d2] = size(img_B_raw);
if d2==3
img_B = rgb2gray(img_B_raw);
else
img_B = img_B_raw;
end
%2.计算MSD矩阵
mmsd = zeros(r1-r2,c1-c2);
for i = 1:r1-r2
for j = 1:c1-c2
temp = img_A(i:i+r2-1,j:j+c2-1);
mmsd(i,j) = mmsd(i,j) + sum(sum((temp - img_B).^2))/(r2*c2);
end
end
%3.定位匹配位置
min_msd = min(min(mmsd));
[x,y] = find(mmsd == min_msd);
x = x(1); %定位到的第一个位置
y = y(1);
%4.保存结果图
getImg = img_A_raw(x:x+r2-1,y:y+c2-1,1:3);
imwrite(getImg,'.\output\MSD_match.bmp');
fprintf('\n Done. \n');
(5) NCC
clear all;
close all; clc;
%1.读取图片
img_A_dir = '.\data\lena.bmp'; %待寻母图
img_A_raw = imread(img_A_dir);
[r1,c1,d1] = size(img_A_raw);
if d1==3 %灰度化
img_A = rgb2gray(img_A_raw);
else
img_A = img_A_raw;
end
img_B_dir = '.\data\refer.bmp'; %模板图
img_B_raw = imread(img_B_dir);
[r2,c2,d2] = size(img_B_raw);
if d2==3
img_B = rgb2gray(img_B_raw);
else
img_B = img_B_raw;
end
%2.计算NCC矩阵
mNCC = zeros(r1-r2,c1-c2);
for i = 1:r1-r2
for j = 1:c1-c2
temp = img_A(i:i+r2-1,j:j+c2-1);
mean_temp = mean(temp(:)); %temp均值
mean_B = mean(img_B(:)); %img_B均值
inp = sum(sum((temp - mean_temp).*(img_B - mean_B))); %两向量内积
mod1 = sqrt(sum(sum((temp - mean_temp).^2))); %模长1
mod2 = sqrt(sum(sum((img_B - mean_B).^2))); %模长2
ncc = inp / (mod1*mod2);
mNCC(i,j) = mNCC(i,j) + ncc;
end
end
%3.定位匹配位置
max_ncc = max(max(mNCC)); %最大ncc值
[x,y] = find(mNCC == max_ncc);
x = x(1); %定位到的第一个位置
y = y(1);
%4.保存结果图
getImg = img_A_raw(x:x+r2-1,y:y+c2-1,1:3);
imwrite(getImg,'.\output\NCC_match.bmp');
fprintf('\n Done. \n');
End.