基本原理
f ( x , y ) ∗ h ( x , y ) = ∑ m = 0 M − 1 ∑ n = 0 N − 1 f ( m , n ) h ( x − m , y − n ) f(x,y)\ast h(x,y) = \sum_{m=0}^{M-1}\sum_{n=0}^{N-1}f(m,n)h(x-m,y-n) f(x,y)∗h(x,y)=m=0∑M−1n=0∑N−1f(m,n)h(x−m,y−n)
上述公式为二维循环卷积,这是极为重要的一个公式,是因为:
f
(
x
,
y
)
∗
h
(
x
,
y
)
⇔
F
(
u
,
v
)
H
(
u
,
v
)
f(x,y)\ast h(x,y) \Leftrightarrow F(u,v)H(u,v)
f(x,y)∗h(x,y)⇔F(u,v)H(u,v)
在频率域的乘积通过IDFT变换为空间域的卷积,而这就是线性滤波的基础。从而两个离散函数的卷积会产生周期性问题,需要填0来进行解决
对于二维图像而言:假设f(x,y) h(x,y)
为A×B C×D
像素图像阵列,我们对其进行填充,按如下公式:
P
≥
A
+
C
−
1
P \geq\ A+C-1
P≥ A+C−1
Q ≥ B + D − 1 Q\geq\ B+D-1 Q≥ B+D−1
当Q和P为2次幂时,执行速度较快,故而我们需要代码来得到P和Q的最小偶数值
function PQ = paddedsize(AB,CD,PARAM)
% 根据关系式:P>=A+C-1,Q>=B+D-1
% PQ = paddedsize(AB) 意味着只有一个参数AB,为一个含两个元素的矩阵,那么
% PQ = 2*AB
%
% PQ = paddedsize(AB,CD) 意味着两个参数,AB CD为含两个参数的矩阵,那么根据
% 关系式就可以得到最小的PQ = AB+CD-1,此时我们还需要判断第二个参数是CD还是
% PARAM,通过ischar进行判断
%
% PQ = paddedsize(AB,'PWR2') 'PWR2'为一个可算参数,当输入时,本函数的尺寸
% 等于最接近2的整数幂的方形图像。通过nextpow2可以实现
%
% PQ = paddedsize(AB,CD,'PWR2')同理
if nargin == 1
PQ = 2*AB;
elseif nargin == 2 && ~ischar(CD)
% ceil(X)将X的每个元素四舍五入到大于或等于该元素的最接近整数(返回的都是偶数)。
% P Q为偶数时非常的快
PQ = AB + CD -1;
PQ = 2 * ceil(PQ/2)
elseif nargin == 2
m = max(AB);
P = 2^nextpow2(2*m);
% nextpow2返回2的更高次幂的指数,注意是指数,故而前面是2^。P返回的是
PQ = [P,P];
elseif nargin == 3
m = max([AB CD]); % 最大维数
P = 2^nextpow2(2*m);
PQ = [P,P];
else
error('Wrong number of inputs')
end
有无填充效果展示
% lpfilter函数并没有在路径当中,需要通过addpath向其中进行添加
% addpath 'F:\Mycode\MY_MATLAB_CODE\冈萨雷斯数字图像处理第三版matlab源代码'
I = imread('F:\Mycode\MY_MATLAB_CODE\DIP3E_Original_Images_CH04\Fig0432(a)(square_original).tif');
[M,N] = size(I);
[f,revertclass] = tofloat(I);
F = fft2(f);
sig = 10;
H = lpfilter('gaussian',M,N,sig);
G = H.*F;
g = ifft2(G);
g = revertclass(g);
imshow(g);
我们可以看到,上面是有模糊的黑色边,这正是低通滤波后的结果,但是两边并没有:正是由于没有填充0