1.概念

高斯分布

圖像濾波之高斯濾波介紹

圖像處理算法|高斯濾波

高斯濾波（Gaussian filter）包含很多種，包括低通、高通、帶通等，在圖像上說的高斯濾波通常是指的高斯模糊（Gaussian Blur），是一種高斯低通濾波。通常這個算法也可以用來模糊圖像，提供模糊濾鏡。也可以用來過濾自然界的高斯白噪聲。

高斯分布（正態分布）是一個常見的連續概率分布，正態分布的數學期望值或期望值μ muμ等于位置參數，決定了分布的位置，其方差σ 2 sigma^2σ2的開平方或者標準差σ sigmaσ等于尺度參數，決定了分布的幅度。正態分布的概率密度函數曲線呈鐘形，所以又被稱為鐘形曲線。我們常說的標準正態分布是位置參數μ = 0 , 方差 σ 2 = 1 mu=0,方差sigma^2=1μ=0,方差σ2=1的正態分布。
??若隨機變量X XX服從一個位置參數為μ 、方差為 σ 2 mu、方差為sigma^2μ、方差為σ2的正態分布，可以記為X N ( μ , σ 2 ) X~N(mu,sigma^2)XN(μ,σ2),其概率密度函數為：
g ( x ) = 1 2 π σ e ( ? ( x ? μ ) 2 2 σ 2 ) (1) g(x)=frac{1} {sqrt{2pi} sigma }e^{(-frac{{(x-mu)}^{2} }{2sigma^{2}})} ag{1}g(x)=2πσ1e(?2σ2(x?μ)2)(1)
??高斯濾波器是一種根據高斯函數的形狀來選擇權值的線性平滑濾波器，對于抑制服從正態分布的噪聲非常有效，一維零均值高斯函數為：
g ( x ) = 1 2 π σ e ? x 2 2 σ 2 (2) g(x)=frac{1} {sqrt{2pi} sigma }e^{-frac{x^{2} }{2sigma^{2}}} ag{2}g(x)=2πσ1e?2σ2x2(2).
??其中，高斯分布參數（σ sigmaσ）決定了高斯函數的寬度，一維高斯函數的圖形如下圖：

??二維高斯高斯分布的函數為：
g ( x , y ) = 1 2 π σ 2 e ? x 2 + y 2 2 σ 2 (3) g(x,y)=frac{1}{2pisigma^2}e^{frac{-x^2+y^2}{2sigma^2}} ag{3}g(x,y)=2πσ21e2σ2?x2+y2(3)
??二維高斯分布的圖像為：

2.高斯濾波性質

高斯函數具有五個重要的性質，這些性質使得它在早期圖像處理中特別有用．這些性質表明，高斯平滑濾波器無論在空間域還是在頻率域都是十分有效的低通濾波器，且在實際圖像處理中得到了工程人員的有效使用．高斯函數具有五個十分重要的性質，它們是：
??1.二維高斯函數具有旋轉對稱性，即濾波器在各個方向上的平滑程度是相同的．一般來說，一幅圖像的邊緣方向是事先不知道的，因此，在濾波前是無法確定一個方向上比另一方向上需要更多的平滑．旋轉對稱性意味著高斯平滑濾波器在后續邊緣檢測中不會偏向任一方向．
??2.高斯函數是單值函數．這表明，高斯濾波器用像素鄰域的加權均值來代替該點的像素值，而每一鄰域像素點權值是隨該點與中心點的距離單調增減的．這一性質是很重要的，因為邊緣是一種圖像局部特征，如果平滑運算對離算子中心很遠的像素點仍然有很大作用，則平滑運算會使圖像失真．
??3.高斯函數的傅立葉變換頻譜是單瓣的．正如下面所示，這一性質是高斯函數付立葉變換等于高斯函數本身這一事實的直接推論．圖像常被不希望的高頻信號所污染(噪聲和細紋理)．而所希望的圖像特征（如邊緣），既含有低頻分量，又含有高頻分量．高斯函數付立葉變換的單瓣意味著平滑圖像不會被不需要的高頻信號所污染，同時保留了大部分所需信號．
??4.高斯濾波器寬度(決定著平滑程度)是由參數σ表征的，而且σ和平滑程度的關系是非常簡單的．σ越大，高斯濾波器的頻帶就越寬，平滑程度就越好．通過調節平滑程度參數σ，可在圖像特征過分模糊(過平滑)與平滑圖像中由于噪聲和細紋理所引起的過多的不希望突變量(欠平滑)之間取得折衷．
??5.由于高斯函數的可分離性，較大尺寸的高斯濾波器可以得以有效地實現．二維高斯函數卷積可以分兩步來進行，首先將圖像與一維高斯函數進行卷積，然后將卷積結果與方向垂直的相同一維高斯函數卷積．因此，二維高斯濾波的計算量隨濾波模板寬度成線性增長而不是成平方增長．

3.高斯濾波的原理與實現

3.1 高斯模板的生成

高斯濾波既能用來過濾高斯噪音，也可用來做高斯模糊。
??要模糊一張圖像，可以直接用均值濾波來做簡單的模糊，但是這樣做顯然不大合理，因為圖像是連續的，離卷積核中心的點關系更加密切，越遠的點關系越疏遠，這個時候就需要加權平均。明顯的離中心點越近的像素點權重越大。而正態分布顯然是一種可取的權重分配方式，又由于圖像是二維的，所以需要使用二維的高斯函數。所以高斯濾波的本質是利用高斯函數來生成高斯核（高斯卷積模板）來對圖像進行卷積操作。
??理論上高斯分布在所有定義域上都有非負值，這就需要一個無限大的卷積核，但是實際上，僅需要取均值的三倍標準差（即3 σ 3sigma3σ）內的值，以外的部分去掉即可。
??高斯濾波最重要的就是找到高斯模板然后進行卷積，以3X3高斯模板為例，假設中心點的坐標為（0，0），根據二維高斯函數g ( x , y ) g(x,y)g(x,y)，還需要設定σ sigmaσ的值，假定σ = 0.8 sigma=0.8σ=0.8(這個值不宜過大，否則就會變成均值濾波)，可以根據坐標值來算出對應的高斯模板。
??假定中心點的坐標為（0，0），那么距離它最近的八個點的坐標如下：

??將坐標以及設定的σ sigmaσ帶入二維高斯函數中，可以得到對應點的坐標的權重如下：

為了防止高斯濾波后的圖像偏亮或者偏暗，我們還需要對圖像進行歸一化，這九個點的權重和等于0.9125991，因此需要分別對這九個數除以0.9125991，最終得到的高斯模板為：

??有了高斯模板就可以對圖像進行卷積了，但是在FPGA中對于小數的運算不友好，于是我們將這個3X3模板擴大了16倍，得到了近似的整數模板，卷積完成后再除以16來做定浮點數的近似計算。16倍后的模板高斯卷積模板為：

??3X3的高斯模板生成與定浮點數的MATLAB算法實現為：

clear;
clc;
close all;
sigma = 0.8 ;
A = exp(-(1+1)/(2*sigma*sigma))/(2*pi*sigma*sigma);
B = exp(-(1+0)/(2*sigma*sigma))/(2*pi*sigma*sigma);
C = exp(-(0+0)/(2*sigma*sigma))/(2*pi*sigma*sigma);
D = A*4 + B*4 + C;

gauss_double = [A,B,A;B,C,B;A,B,A];
gauss_normal = gauss_double / sum(sum(gauss_double));
gauss_integer = floor(gauss_normal/gauss_normal(1,1));

3.2 高斯算法的FPGA實現

生成卷積模板后的FPGA實現與仿真參考前面幾章，這里只貼高斯算法的FPGA實現部分：

module gassin_filter#(
    parameter DW = 8
)(
    input   wire                clk         ,
    input   wire                rst_n       ,

    input   wire                matrix_de   ,
    input   wire    [DW-1:0]    matrix11    ,
    input   wire    [DW-1:0]    matrix12    ,
    input   wire    [DW-1:0]    matrix13    ,
    input   wire    [DW-1:0]    matrix21    ,
    input   wire    [DW-1:0]    matrix22    ,
    input   wire    [DW-1:0]    matrix23    ,
    input   wire    [DW-1:0]    matrix31    ,
    input   wire    [DW-1:0]    matrix32    ,
    input   wire    [DW-1:0]    matrix33    ,

    output  wire                gassin_data_de     ,
    output  wire    [DW-1:0]    gassin_data 
);

//gassin_filter
//  [1      2       1]
//  [2      4       2]
//  [1      2       1]
reg     [1:0]           matrix_de_r ;

reg     [DW+2:0]        one_line    ;
reg     [DW+2:0]        two_line    ;
reg     [DW+2:0]        three_line  ;
reg     [DW+4:0]        sum_matrix  ;

always @(posedge clk)begin
    if(rst_n == 0)begin
        matrix_de_r     <=  0;
    end 
    else begin
        matrix_de_r     <=  {matrix_de_r[0],matrix_de};
    end
end 


always @(posedge clk)begin
    if(rst_n==0)begin
        one_line    <=  0   ;
        two_line    <=  0   ;
        three_line  <=  0   ;
    end 
    else if(matrix_de)begin
        one_line    <=  matrix11    +   2*matrix12  +   matrix13    ;
        two_line    <=  2*matrix21  +   4*matrix22  +   2*matrix23  ;
        three_line  <=  matrix31    +   2*matrix32  +   matrix33    ;
    end
    else begin
        one_line    <=  0   ;
        two_line    <=  0   ;
        three_line  <=  0   ;
    end  
end 

always @(posedge clk)begin
    if(rst_n == 0)begin
        sum_matrix  <=  0;
    end 
    else if(matrix_de_r[0])begin
        sum_matrix  <=  one_line + two_line + three_line ;
    end
    else begin
        sum_matrix  <=  0;
    end
end

assign  gassin_data     =   sum_matrix[DW+4:4]  ;
assign  gassin_data_de  =   matrix_de_r[1]      ;

endmodule

3.3 高斯算法的MATLAB實現以及驗證

clear;
clc;
close all;
sigma = 0.8 ;
A = exp(-(1+1)/(2*sigma*sigma))/(2*pi*sigma*sigma);
B = exp(-(1+0)/(2*sigma*sigma))/(2*pi*sigma*sigma);
C = exp(-(0+0)/(2*sigma*sigma))/(2*pi*sigma*sigma);
D = A*4 + B*4 + C;

gauss_double = [A,B,A;B,C,B;A,B,A];
gauss_normal = gauss_double / sum(sum(gauss_double));
gauss_integer = floor(gauss_normal/gauss_normal(1,1));

GRAY = imread('../img/gray.bmp');
[row,col] = size(GRAY);
gassin_padding  =   zeros(row+2,col+2);
gassin_result   =   zeros(row,col);

for i = 1:row
    for j = 1:col
        gassin_padding(i+1,j+1) = GRAY(i,j);
    end
end

for i = 1:row+2
    gassin_padding(i,1) = gassin_padding(i,2);
    gassin_padding(i,col+2) = gassin_padding(i,col+1);
end

for i = 1:col+2
   gassin_padding(1,i) = gassin_padding(2,i);
   gassin_padding(row+2,i) = gassin_padding(row+1,i);
end

for i = 2:row+1
    for j = 2:col+1
        matrix11 = gassin_padding(i-1,j-1);
        matrix12 = gassin_padding(i-1,j);
        matrix13 = gassin_padding(i-1,j+1);
        
        matrix21 = gassin_padding(i,j-1);
        matrix22 = gassin_padding(i,j);
        matrix23 = gassin_padding(i,j+1);
        
        matrix31 = gassin_padding(i+1,j-1);
        matrix32 = gassin_padding(i+1,j);
        matrix33 = gassin_padding(i+1,j+1);

        matrix = [matrix11,matrix12,matrix13;matrix21,matrix22,matrix23;matrix31,matrix32,matrix33];
        gassin_mult = matrix.* gauss_integer;

        sum_gassin_matrix = sum(sum(gassin_mult()));
        gassin_result(i-1,j-1) = sum_gassin_matrix/16;
    end
end

a = textread('../data/gassin_filter.txt','%s');
IMdec1 = hex2dec(a);

IM1 = reshape(IMdec1,col,row);
fpga_Y = uint8(IM1)';

b = textread('../data/pre.txt','%s');

subplot(1,3,1)
matlab_Y = uint8(floor(gassin_result));
imshow(matlab_Y),title('MATLAB gassin算法圖像');
subplot(1,3,2)
imshow(fpga_Y),title('FPGA gassin算法圖像');
subplot(1,3,3)
imshow(GRAY),title('原圖像');

sub = matlab_Y - fpga_Y;

min_sub = min(min(sub));
max_sub = max(max(sub));