1柒爵、基礎(chǔ)
卷積定理
折疊誤差補(bǔ)零
當(dāng)處理DFT時(shí)，圖像及其變換是周期的赚爵。在周期接近函數(shù)非零部分的持續(xù)周期時(shí)棉胀，對(duì)周期函數(shù)進(jìn)行卷積會(huì)引起相鄰周期的串?dāng)_法瑟。
這種稱(chēng)為 折疊誤差 的串?dāng)_可通過(guò)補(bǔ)零方法來(lái)避免。
假設(shè)函數(shù) f(x,y) 和 h(x,y) 的大小分別為 A x B 和 C x D唁奢。通過(guò)對(duì) f 和 g 補(bǔ)零霎挟，構(gòu)造兩個(gè)大小均為 P x Q 的擴(kuò)展函數(shù)。
按如下方式進(jìn)行選擇可以避免折疊誤差：
下面 paddedsize 的函數(shù)可以滿足上面等式的 P 和 Q 的最小偶數(shù)值麻掸。
function PQ = paddedsize(AB, CD, PARAM)
if nargin == 1
PQ = 2*AB;
elseif nargin == 2 && -ischar(CD)
PQ = AB + CD - 1;
PQ = 2 * ceil(PQ / 2);
elseif nargin == 2
m = max(AB);
P = 2^nextpow2(2*m);
PQ = [P, P];
elseif (nargin == 3) && strcmp(PARAM, 'pwr2')
m = max([AB CD]);
P = 2^nextpow2(2*m);
PQ = [P, P];
else
error('Wrong number of inputs.')
end
示例
[f, revertcalss] = tofloat(f);
sig = 10;
PQ = paddedsize(size(f));
Fp = fft2(f, PQ(1), PQ(2));
Hp = lpfilter('gaussian', PQ(1), PQ(2), 2*sig); %濾波器大小是未填充時(shí)的兩倍
Gp = Hp .* Fp;
gp = ifft2(Gp);
gpc = gp(1:size(f, 1), 1:size(f, 2));
gpc = revertcalss(gpc);
imshow(gp)
2酥夭、DFT 濾波的基本步驟
1、使用函數(shù) tofloat 把輸入圖像轉(zhuǎn)換為浮點(diǎn)圖像：
[f, revertclass] = tofloat(f);
2脊奋、使用函數(shù) paddedsize 獲得填充參數(shù)：
PQ = paddedsize(size(f));
3熬北、得到有填充圖像的傅里葉變換：
F = fft2(f, PQ(1), PQ(2));
4、生成一個(gè)大小為 PQ(1) * PQ(2) 的濾波器函數(shù)诚隙。如果它是居中的讶隐，使用前要令 H = ifftshift(H)。
5久又、使用濾波器乘以該變換：
G = H .* F;
6巫延、獲得G的IFFT：
g = ifft2(G);
7、將左上部的矩形修剪為原始大械叵：
g = g(1:size(f, 1), 1:size(f, 2));
8炉峰、需要時(shí)，將濾波后的圖像轉(zhuǎn)換為輸入圖像的類(lèi)：
g = revertclass(g);
- 用于頻率域?yàn)V波的M函數(shù)
function g = dftfilt(f, H, classout)
[f, revertcalss] = tofloat(f);
F = fft2(f, size(H, 1), size(H, 2))
g = ifft2(H .* F);
g = g(1:size(f, 1), 1:size(f, 2));
if nargin == 2 || strcmp(classout, 'original')
g = revertcalss(g);
elseif strcmp(classout, 'fltpoint')
return
else
error('Undefined class fot the output image.')
end
3犯建、從空間濾波器獲得頻域?yàn)V波器
通常當(dāng)濾波器較小時(shí)讲冠，在計(jì)算上空間濾波要比頻率域?yàn)V波更有效率。
當(dāng)濾波器大約有 32 個(gè)或更多元素時(shí)适瓦，使用 FFT 算法的濾波處理要比空間實(shí)現(xiàn)快竿开。
函數(shù) freqz2 可以將空間濾波器轉(zhuǎn)換為頻率域?yàn)V波器
H = freqz2(h, R, C)
示例
h = fspecial('sobel')'
freqz2(h)
% 頻率域?yàn)V波器
PQ = paddedsize(size(f));
H = freqz2(h, PQ(1), PQ(2));
H1 = ifftshift(H)
函數(shù) dftfilt 用于處理頻率域?yàn)V波
gf = dftfilt(f, H1);
4、在頻率域中直接生成濾波器
- 創(chuàng)建用于實(shí)現(xiàn)頻率域?yàn)V波器的網(wǎng)格數(shù)組
function [U, V] = dftuv(M, N)
u = single(0:(M - 1));
v = single(0:(N - 1));
idx = find(u > M/2);
u(idx) = u(idx) - M;
idy = find(v > N/2);
v(idy) = v(idy) - N;
[V, U] = meshgrid(v, u);
% ----------------------- %
[U, V] = dftuv(8, 5);
DSQ = hypot(U, V)
-
低通（平滑）頻率域?yàn)V波器
理想低通濾波器（ILPF）
巴特沃斯低通濾波器（BLPF）
高斯低通濾波器
示例玻熙，對(duì)圖像 f 應(yīng)用一個(gè)高斯低通濾波器否彩，其中 D0 為所填充圖像寬度的 5%。
[f, revertclass] = tofloat(f);
PQ = paddedsize(size(f));
[U, V] = dftuv(PQ(1), PQ(2));
D = hypot(U, V);
D0 = 0.05 * PQ(2);
F = fft2(f, PQ(1), PQ(2));
H = exp(-(D .^ 2) / (2 * (D0^2)));
g = dftfilt(f, H);
g = revertclass(g)
- 低通濾波器整合函數(shù)
function H = lpfilter(type, M, N, D0, n)
[U, V] = dftuv(M, N);
D = hypot(U, V);
switch type
case 'ideal'
H = single(D <= D0);
case 'btw'
if nargin == 4
n = 1;
end
H = 1 ./ (1 + (D./D0) .^ (2*n));
case 'gaussian'
H = exp(-(D.^2) ./ (2*(D0^2)));
otherwise
error('Unknown filter type.')
end
5嗦随、繪制線框圖和表面圖
繪制二維函數(shù)的最簡(jiǎn)單方法是使用函數(shù) mesh
mesh(H(1:k:end, 1:k:end))
默認(rèn)繪出彩色網(wǎng)格圖命令為
colormap([0 0 0])
關(guān)閉網(wǎng)格和坐標(biāo)軸命令
grid off
axis off
觀察點(diǎn)由函數(shù) view 控制
view(az, el) %az,el 代表方位角和仰角
view(3) %默認(rèn)觀察點(diǎn)
示例列荔，繪制一個(gè)高斯低通濾波器
H = fftshift(test2('gaussian', 500, 500, 50));
mesh(double(H(1:10:500, 1:10:500)))
axis tight
colormap([0 0 0])
- 表面圖
surf(H)
6、高通（銳化）頻率域?yàn)V波器
理想高通濾波器（IHPF）
巴特沃斯高通濾波器（BHPF）
高斯高通濾波器
基于 lpfilter 構(gòu)建一個(gè)高通濾波器函數(shù)
function H = hpfilter(type, M, N, D0, n)
if nargin == 4
n = 1;
end
H1p = lpfilter(type, M, N, D0, n);
H = 1 - H1p;
% ---------------------------%
H = fftshift(hpfilter('ideal', 500, 500, 50));
高頻強(qiáng)調(diào)濾波
高通濾波器偏離了直流項(xiàng)枚尼，將圖像的平均值降低為0.補(bǔ)償方法之一是給高通濾波器加上一個(gè)偏移量贴浙。把偏移量與將濾波器乘以一個(gè)大于1的常數(shù)結(jié)合起來(lái)的方法就稱(chēng)為 高頻強(qiáng)調(diào)濾波。
示例署恍，聯(lián)合使用高頻強(qiáng)調(diào)濾波和直方圖均衡
PQ = paddedsize(size(f));
D0 = 0.05 * PQ(1);
HBW = hpfilter('btw', PQ(1), PQ(2), D0, 2);
H = 0.5 + 2*HBW;
gbw = dftfilt(f, HBW, 'fltpoint');
gbw = gscale(gbw);
ghf = dftfilt(f, H, 'fltpoint');
ghf = gscale(ghf);
ghe = histeq(ghf, 256);
