精华内容
下载资源
问答
  • 高斯滤波的开始——高斯的计算

    千次阅读 2019-10-01 14:42:35
    高斯滤波的开始——高斯的计算 首先先说说高斯滤波的含义:高斯滤波是一种线性平滑滤波,适用于消除高斯噪声,广泛应用于图像处理的减噪过程。通俗的讲,高斯滤波就是对整幅图像进行加权平均的过程,每一个像素点...

    高斯滤波的开始——高斯核的计算

    首先先说说高斯滤波的含义:高斯滤波是一种线性平滑滤波,适用于消除高斯噪声,广泛应用于图像处理的减噪过程。通俗的讲,高斯滤波就是对整幅图像进行加权平均的过程,每一个像素点的值,都由其本身和邻域内的其他像素值经过加权平均后得到。高斯滤波的具体操作是:用一个模板(或称卷积、掩模)扫描图像中的每一个像素,用模板确定的邻域内像素的加权平均灰度值去替代模板中心像素点的值。(来源于百度)
    对于高斯核的计算,网上有许多,但都讲的马马虎虎,因此,个人写了一个基于matlab,简单,易懂的方法。话不多说,直接上代码;

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    %%                计算高斯核                 %%
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    pi=3.1415926;%高斯公式中Π的值
    sum=0;%计算累计权重的值
    d=1.5;%高斯核标准差
    size=3;%高斯核大小
    center=(size/2)+0.5;%模板的中心,这里取整
    A=1/(2*pi*d*d);%高斯公式指数前面的系数
    for i=1:3
        x2=(i-center)*(i-center);
        for j=1:3
            y2=(j-center)*(j-center);
            B=exp(-(x2+y2)/(2*d*d));
            C(i,j)=A*B;
            sum=sum+C(i,j);
        end
    end
    %整数形式的高斯核,需要进行归一化,即把左上角的值化为1
    %下面进行归一化
    k=1/C(1,1);
    for i=1:3
        for j=1:3
            D(i,j)=C(i,j)*k;%进行归一化后整数形式的高斯核
        end
    end
    %若是小数形式的高斯核,则在不需要归一化,而需要进权重分配,高斯核的每个系数要除以所有系数的和。
    %下面进行权重分配
    for i=1:3
        for j=1:3
            E(i,j)=C(i,j)/sum%进行归一化后整数形式的高斯核
        end
    end
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%截至这里,高斯核以计算完毕%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    下面是输出的结果:标准差=1.5,模板大小为3
    上图为直接输出的高斯核,没有进行归一化和权重分配
    上图是经过归一下后输出的高斯核
    上图是经过权重分配后输出的高斯核
    所谓权重分配,就是把高斯核的每个系数乘以所有系数的和,重新输出。

    展开全文
  • 高斯滤波就是对整幅图像进行加权平均过程,每一个像素点值,都由其本身和邻域内其他像素值经过加权平均后得到。消除图像在数字化过程中产生或者混入噪声。  1.2高斯滤波核计算 二维高斯分布: 假定中心点...
  • 几年前的总结 ,在此记录一下: 一 相关 相关是用来衡量两个信号 "形状" 和 "相位"(起始位置)的相似程度 . 计算相关性的两个信号...相关的计算是满足交换律的, "x 与 h的相关性" 等价于 "h与x的相关性" , 在这...

    几年前的总结 ,在此记录一下:

    一 相关

    相关是用来衡量两个信号  "形状" 和 "相位"(起始位置)的相似程度 . 

    计算相关性的两个信号必须要有相同起始位置和长度,然后对两个信号逐点乘积加和,当两个信号完全一致时相关性是最高的.

    离散信号相关计算定义:

    信号中的直流分量通常对信号表达是没有意义的,反而会产生干扰。所以某些相关性计算要先对信号做减直流处理,如数理统计里用协方差来度量相关性: 

     

    相关的计算是满足交换律的, "x 与 h的相关性" 等价于 "h与x的相关性" , 在这里我们把 h看成特征模板,x看成未知信号(下文卷积也都按这种约定).从这个角度来看

    相关性计算也表达了"模板信号对未知信号加权和"的意思.比如要从一个很长的序列里匹配一段已知信号,相关计算算法可能从任意起始位置截取一段 x(0 ~ n) , 与 h 做加权和,得分最高的位置即是最匹配的信号.

     

    二卷积(通常是指线卷积)

    卷积与相关类似,按字面理解,先卷(与相关的主要区别),再积,离散卷积公式:

    公式如下:

                    先把模板h翻卷(按原点对褶) h(i) -> h(-i)

                    h的零点对齐(平移) 到x的目标点n  h(-i) -> h(n - i)

                    再对两段序列求逐点积,加和

    当卷积模板中心对称时,卷积的结果与相关是一致的!!

    同样的卷积计算满足交换率,即 (x * h)(n)  <==> (h * x)(n) 

    卷积的物理解释:

    a  : 假如一幅只有不同黑色程度的水墨画下雨泡了水。  x(p) 表示水墨画在 p 点的含墨量,h(delta) 表示字画泡水后,单位墨向位置 (p + delta)扩散了多少 ( 实际扩散了 x(p) * h(p + delta)  , 线性时不变系统) , 问题求是这幅画最终每一点的墨迹量。  注意, x 的定义域为整个画的二维坐标, h的定义域为与x,y平面原点中心对称的矩形区域。

     

        显然用卷积的定义求解就很容易, 目标点 P1受到模板h定义的所有邻域点影响,对于邻域点 Q1 , 位移了delta 到 P1 , 则 Q1 + delta = P1 ,即受Q1的影响为 x(Q1) * h (P1 - Q1) ,  p1 点收到总的影响为:

    b  : 假如小明在t时刻脸被挨打强度为x(t),  受单位强度的力挨打后, 脸鼓起的包随时间变化函数为 h(t) , 那问题来了, 求{0 , ... N}任意时刻小明脸鼓起的包强度。  显然 x(t) 和h(t) 的定义域只有时间正轴,包括零点。

    信号里 x 被称为冲激函数, h 为冲激响应,这里冲激函数只会对后面的时刻产生响应。挨打之前脸上就有包那是不负责任的! 同样计算公式如下:

    卷积计算的几种方法:             

    a 按公式 反褶(与零点)、平移(h 0点到目标位置n) , 逐点乘积 再求和

     

    b 多项式方法(竖式计算法则)一种快捷计算方式

    X1 = 

    X2 =  

     

    如图所示,竖式计算时,两个序列要做到齐次项对齐,即定义域零点对齐。

     

    三圆周卷积

    圆周卷积也称作循环卷积,百度百科里的定义如下图:

     

    两个序列x,h (序列下标从0开始) , 作L点的圆周卷积:

            1    首先将x ,h补成L长度的序列(少则补0,多则截断) 

            2    模板序列h逆时针排成一圈,x 顺时针排成一圈 

            3    求 0 点的圆周卷积对应位置相乘加和

            4    求 点n 的圆周圈积,将模板h 顺时针转n下,再逐点乘积求和。

     

    也来个定义:

     

    假设两个序列长度分别为 M,N , L >= M + N - 1 时 圆周卷积和线性卷积值是一致的!! 

     

    圆周卷积求法2:  竖式计算

    X = {1 , 4  , 3 } ,  H = {1 , 2 , 3}

    方法是用H的每一个元素去乘  X 每一组顺时针循环移位的序列,再求和:

                            1 X    {1 , 4 ,  3}  +  2 X   {3 , 1 , 4}  + 3 X {4 , 3 , 1}  =          1          4            3

                                                                                                                              6          2           8         

                                                                                                                              12        9           3

                                                                                                                        ---------------------------------

                                                                                                                             19        15         14

    圆周卷积求法3:  矩阵乘法

    X = {1 , 4  , 3 } ,  H = {1 , 2 , 3}

     

                                             0

    先将 H周期延拓(1,2,3 , 1,2,3 , 1,2,3),沿0点反转 , 再依次求每组顺时针循环右称的序列

    HM = 

        |    1    3     2    |

        |    2    1     3    |

        |    3    2     1    |

    F列向量 =    = {19 , 15 , 14}

    圆周卷积求法4:  DFT 

    经证明圆周卷积与 两个序列先DFT再乘积再逆DFT结果是一致。个人理解圆周卷积只是计算上有DFT 这种快速算法,没有线卷积那么强的物理意义,L = N + M - 1 时二者实现了统一。 L < N + M - 1 时会造成混叠效应。

     

    四 利用圆周卷积来计算线性卷积

    4.1 模板定义域为正轴数据 

    a 将 x & h 扩展成 L = N + M - 1 序列 

    b 计算圆周卷积。

    c 再截掉 最后  2 *floor(M / 2) = M - 1个值 

     

    4.2 模板h定义域包含负值 

    这里举一个 h 中心对称的例子,同样 x 序列长度为N , h 长度为M M为奇数, 且h 定义域为 [ -M / 2 , M/ 2]

    a 将包含负值的定义域平移到0点变成 [0 , M],这样最终卷积的结果会延时(晚)  floor(M/2) 个点!

    b 将 x & h 扩展成 L = N + M - 1 序列 

    c 计算圆周圈积。 

    4 截掉最前面 floor(M / 2) 和最后 (M / 2) 个点

     

    4.3 线卷积matlab 代码 

    function [s] = linear_cov(x , h)
    
    
    N = size(x , 2);
    M = size(h , 2);
    L = N + M - 1;
    
    
    % fprintf(2 , 'M = %d , N = %d , L = %d \n' , M ,  N , L);
    
    
    s = zeros(1 , L);
    
    
    for i=1:N
        for j=1:M
            s(j + i - 1) = s(j + i - 1) + x(i) * h(j);
        end
    end   

    测试:

    X=[4,5,6,7,8,9,10];
    H=[1,2,3];
    s=linear_cov(X,H);
    %s=conv(H,X);

    4.4 圆周卷积matlab 代码

    function [s] = circle_cov(x , h , L)
    
    
    N = size(x , 2);
    M = size(h , 2);
    
    
    if N < L 
        X1 = [x , zeros(1 , L - N)];
    else
        X1 = x(1:L);
    end
    
    
    if M < L
        H1 = [h , zeros(1 , L - M)];
    else
        H1 = h(1:L);
    end
    
    
    s = zeros(1 , L);
    
    
    for i=1:L
        for j=1:L
            s(i) = s(i) + X1(j) * H1(1 + mod(i - j , L));
        end
    end

    测试:

    X=[4,5,6,7,8,9,10];
    H=[1,2,3];
    s=circle_cov(X,H , 9);

     

    9点的圆周卷积结果与线卷积一致。

    五 FFT 快速卷积

    求x, h的卷积 matlab 代码: ( h (x) : x >= 0)

    function [s] = fft_cov(x , h)
    
    
    N = size(x , 2);
    M = size(h , 2);
    L = N + M - 1;
    
    
    X = [x , zeros(1 , L - N)];
    H = [h , zeros(1 , L - M)];
    
    
    f_x = fft(X);
    f_h = fft(H);
    
    
    f_s = f_x .* f_h;
    s   = ifft(f_s);

    如果 h 定义域  x [-a , b] , a < 0 , 卷积输出会延时 a 长度。 对于二维图像卷积适用。

    六 IIR 型高斯滤波器

    参考:

    https://www.cnblogs.com/luo-peng/p/5223910.html

    IIR 型高斯滤波器特点是计算复杂度与 核 大小无关有一定参考意义

    实现代码

    //iir_guass_blur.h
    
    #ifndef __IIR_GUASS_BLUR_H
    #define __IIR_GUASS_BLUR_H
    
    #ifdef __cplusplus
    extern "C"
    {
    #endif
    
    long IIR_guass_filter_create(int w , int h , float sigma);
    
    void IIR_guass_filter_do(long  h_flt , const uint8_t * src_plane , int src_pitch , uint8_t *dst_plane , int dst_pitch);
    
    void IIR_guass_filter_destroy(long  h_flt);
    
    long IIR_guass_filterI_create(int w , int h , float sigma);
    
    void IIR_guass_filterI_do(long  h_flt , const uint8_t * src_plane , int src_pitch , uint8_t *dst_plane , int dst_pitch);
    
    void IIR_guass_filterI_destroy(long  h_flt);
    
    
    
    #ifdef __cplusplus
    }
    #endif
    
    #endif
    
    //iir_guass_blur.c
    
    #include <malloc.h>
    #include <stdlib.h>
    #include <string.h>
    #include <stdint.h>
    #include <assert.h>
    #include <math.h>
    #include "iir_guass_blur.h"
    
    /* *** IIR guass blur formula : ***
     *  
     *  forward pass  : 
            w[n]  =  B * x[n] + (b1 * y[n - 1] + b2 * y[n - 2] + b3 * y[n - 3]) / b0;
      
     *  backword pass :
           out[n] =  B * w[n] + (b1 * out[n + 1] + b2 * out[n + 2] + b3 * out[n + 3]) / b0;
    
     *      B = 1 - (b1 + b2 + b3) / b0; 
     *
    */
    
    typedef struct IIR_guass_coefs
    {
        float   B;
        float   b[4];    
    }IIR_guass_coefs;
    
    typedef struct IIR_guass_flt_context
    {
        IIR_guass_coefs     coeff;
        int                 w;
        int                 h;
        float               sigma;
        float               *w1;
        float               *w2;
        uint8_t             *dst_1pass;
    }IIR_guass_flt_context;
    
    
    static void compute_coefs3(IIR_guass_coefs *c , float sigma)  
    {  
        float q, q2, q3;  
        if (sigma >= 2.5){  
            q = 0.98711 * sigma - 0.96330;  
        }  
        else if ((sigma >= 0.5) && (sigma < 2.5)){  
            q = 3.97156 - 4.14554 * (float) sqrtf ((float) 1 - 0.26891 * sigma);  
        }  
        else {  
            q = 0.1147705018520355224609375;  
        }  
    
        q2 = q * q;  
        q3 = q * q2;  
        c->b[0] = 1/(1.57825+(2.44413*q)+(1.4281 *q2)+(0.422205*q3));  
        c->b[1] = (        (2.44413*q)+(2.85619*q2)+(1.26661 *q3)) *c->b[0];  
        c->b[2] = (                   -((1.4281*q2)+(1.26661 *q3)))*c->b[0];  
        c->b[3] = (                                 (0.422205*q3)) *c->b[0];  
        c->B = 1.0-(c->b[1]+c->b[2]+c->b[3]);  
    }  
    
    long IIR_guass_filter_create(int w , int h , float sigma)
    {
        IIR_guass_flt_context *ctx;
        int max_len;
    
        ctx = (IIR_guass_flt_context *)calloc(1 , sizeof(IIR_guass_flt_context));
        assert(ctx != NULL);
    
        ctx->w      = w;
        ctx->h      = h;
        ctx->sigma  = sigma;
    
        compute_coefs3(&ctx->coeff , sigma);
    
        max_len     = w > h ? w + 3 : h + 3;
    
        ctx->w1     = (float *)malloc(sizeof(float) * max_len);
        ctx->w2     = (float *)malloc(sizeof(float) * max_len);
    
        assert(ctx->w1 != 0);
        assert(ctx->w2 != 0);
    
        ctx->dst_1pass = (uint8_t *)malloc(sizeof(uint8_t) * h * w);
        assert(ctx->dst_1pass != 0);
        
    
        return (long)ctx;
        
    }   
    
    void IIR_guass_filter_onepass( IIR_guass_flt_context *ctx , int width , int height , const uint8_t * src_plane , int src_pitch , uint8_t *dst_plane , int dst_pitch)
    {
        const uint8_t       *p_src;
        uint8_t             *p_dst;
        int n , h;
        float *w1 , *w2 , *b , B;
    
    	p_src       = src_plane;
        p_dst		= dst_plane;
        w1          = ctx->w1;
        w2          = ctx->w2;
        b           = ctx->coeff.b;
        B           = ctx->coeff.B;
    
        for(h = 0 ; h < height ; ++h){
            /*forward process*/
            w1[0] = w1[1] = w1[2] = p_src[0];
            for(n = 3 ; n < width + 3 ; ++n){
                w1[n] = B * p_src[n - 3] + b[1] * w1[n - 1] + b[2] * w1[n - 2] + b[3] * w1[n - 3];
            }
            /*backward process*/
            w2[width] = w2[width + 1] = w2[width + 2] = w1[width + 2];
            for(n = width - 1; n >= 0 ; --n){
                w2[n] = B * w1[n + 3]    + b[1] * w2[n + 1] + b[2] * w2[n + 2] + b[3] * w2[n + 3];
            }
            /*store data*/
            for(n = 0 ; n < width ; ++n){
                p_dst[n * dst_pitch] = w2[n];
            }
            p_src += src_pitch;
            p_dst++;
        }
    }
    	
    void IIR_guass_filter_do(long  h_flt , const uint8_t * src_plane , int src_pitch , uint8_t *dst_plane , int dst_pitch)
    {
        IIR_guass_flt_context *ctx;
    
        ctx = (IIR_guass_flt_context *)h_flt;
    
        /*process horizontal data*/    
    	IIR_guass_filter_onepass(ctx , ctx->w , ctx->h , src_plane , src_pitch , ctx->dst_1pass , ctx->h);
    	
        /*process vertical data*/ 
    	IIR_guass_filter_onepass(ctx , ctx->h , ctx->w , ctx->dst_1pass , ctx->h , dst_plane , dst_pitch);
    	
    }
    
    void IIR_guass_filter_destroy(long  h_flt)
    {
        IIR_guass_flt_context *ctx; 
    
        ctx = (IIR_guass_flt_context *)h_flt;
    
        if (ctx->dst_1pass){
            free(ctx->dst_1pass);
            ctx->dst_1pass = 0;
        }
    
        if (ctx->w1){
            free(ctx->w1);
            ctx->w1 = 0;
        }
    
        if (ctx->w2){
            free(ctx->w2);
            ctx->w2 = 0;
        }
    
        free(ctx);
    }
    
    
    typedef struct IIR_guass_coefsI
    {
        int32_t   B;
        int32_t   b[4];    
    }IIR_guass_coefsI;
    
    typedef struct IIR_guass_fltI_context
    {
        IIR_guass_coefsI    coeff;
        int                 w;
        int                 h;
        float               sigma;
        int32_t             *w1;
        int32_t             *w2;
        uint8_t             *dst_1pass;
    }IIR_guass_fltI_context;
    
    
    static void compute_coefs3_I(IIR_guass_coefsI *c , float sigma)  
    {  
        double q, q2, q3;  
        if (sigma >= 2.5){  
            q = 0.98711 * sigma - 0.96330;  
        }  
        else if ((sigma >= 0.5) && (sigma < 2.5)){  
            q = 3.97156 - 4.14554 * (float) sqrtf ((float) 1 - 0.26891 * sigma);  
        }  
        else {  
            q = 0.1147705018520355224609375;  
        }  
    
        q2 = q * q;  
        q3 = q * q2;  
        c->b[0] = 16384 / (1.57825+(2.44413*q)+(1.4281 *q2)+(0.422205*q3));  
        c->b[1] = (        (2.44413*q)+(2.85619*q2)+(1.26661 *q3))  * c->b[0];  
        c->b[2] = (                   -((1.4281*q2)+(1.26661 *q3))) * c->b[0];  
        c->b[3] = (                                 (0.422205*q3))  * c->b[0];  
        c->B = 16384 - (c->b[1]+c->b[2]+c->b[3]);  
    }  
    
    long IIR_guass_filterI_create(int w , int h , float sigma)
    {
        IIR_guass_fltI_context *ctx;
        int max_len;
    
        ctx = (IIR_guass_fltI_context *)calloc(1 , sizeof(IIR_guass_fltI_context));
        assert(ctx != NULL);
    
        ctx->w      = w;
        ctx->h      = h;
        ctx->sigma  = sigma;
    
        compute_coefs3_I(&ctx->coeff , sigma);
    
        max_len     = w > h ? w + 3 : h + 3;
    
        ctx->w1     = (int32_t *)malloc(sizeof(int32_t) * max_len);
        ctx->w2     = (int32_t *)malloc(sizeof(int32_t) * max_len);
    
        assert(ctx->w1 != 0);
        assert(ctx->w2 != 0);
    
        ctx->dst_1pass = (uint8_t *)malloc(sizeof(uint8_t) * h * w);
        assert(ctx->dst_1pass != 0);
    
        return (long)ctx;
        
    }   
    
    void IIR_guass_filterI_onepass( IIR_guass_fltI_context *ctx , int width , int height , const uint8_t * src_plane , int src_pitch , uint8_t *dst_plane , int dst_pitch)
    {
        const uint8_t       *p_src;
        uint8_t             *p_dst;
        int n , h;
        int32_t *w1 , *w2 , *b , B;
    
    	p_src       = src_plane;
        p_dst		= dst_plane;
        w1          = ctx->w1;
        w2          = ctx->w2;
        b           = ctx->coeff.b;
        B           = ctx->coeff.B;
    
        for(h = 0 ; h < height ; ++h){
            /*forward process*/
            w1[0] = w1[1] = w1[2] = p_src[0];
            for(n = 3 ; n < width + 3 ; ++n){
                w1[n] = B * p_src[n - 3] + b[1] * w1[n - 1] + b[2] * w1[n - 2] + b[3] * w1[n - 3];
            }
            /*backward process*/
            w2[width] = w2[width + 1] = w2[width + 2] = (w1[width + 2] >> 14);
            for(n = width - 1; n >= 0 ; --n){
                w2[n] = (B * (w1[n + 3] >> 14)) + b[1] * w2[n + 1] + b[2] * w2[n + 2] + b[3] * w2[n + 3];
            }
            /*store data*/
            for(n = 0 ; n < width ; ++n){
                p_dst[n * dst_pitch] = (w2[n] >> 14);
            }
            p_src += src_pitch;
            p_dst++;
        }
    }
    	
    void IIR_guass_filterI_do(long  h_flt , const uint8_t * src_plane , int src_pitch , uint8_t *dst_plane , int dst_pitch)
    {
        IIR_guass_fltI_context *ctx;
    
        ctx = (IIR_guass_fltI_context *)h_flt;
    
        /*process horizontal data*/    
    	IIR_guass_filterI_onepass(ctx , ctx->w , ctx->h , src_plane , src_pitch , ctx->dst_1pass , ctx->h);
    	
        /*process vertical data*/ 
    	IIR_guass_filterI_onepass(ctx , ctx->h , ctx->w , ctx->dst_1pass , ctx->h , dst_plane , dst_pitch);
    	
    }
    
    void IIR_guass_filterI_destroy(long  h_flt)
    {
        IIR_guass_fltI_context *ctx; 
    
        ctx = (IIR_guass_fltI_context *)h_flt;
    
        if (ctx->dst_1pass){
            free(ctx->dst_1pass);
            ctx->dst_1pass = 0;
        }
    
        if (ctx->w1){
            free(ctx->w1);
            ctx->w1 = 0;
        }
    
        if (ctx->w2){
            free(ctx->w2);
            ctx->w2 = 0;
        }
    
        free(ctx);
    }
    
    
    

     

    展开全文
  • 针对传统无参考图像质量评价方法计算复杂、难以应用问题,提出一种简单、直接小波高频结构相似性无参考高斯图像质量评价方法。该方法根据自然图像同尺度高频子带间结构相似度(SSIM)随着失真程度增加而降低...
  • 基于CUDA并行计算技术+opencv完成图像高斯滤波和双边滤波,开发版本为VS2019+openCV3.4
  • 图像处理之高斯一阶及二阶导数计算

    万次阅读 多人点赞 2013-11-17 13:09:50
    演示图像中高斯一阶与二阶导数的计算方法,这个在图像的特征提取与处理中十分有用! 一阶导数可以反应出图像灰度梯度的变化情况 二阶导数可以提取出图像的细节同时双响应图像梯度变化情况

    图像处理之高斯一阶及二阶导数计算

     

    图像的一阶与二阶导数计算在图像特征提取与边缘提取中十分重要。一阶与二阶导数的

    作用,通常情况下:

    一阶导数可以反应出图像灰度梯度的变化情况

    二阶导数可以提取出图像的细节同时双响应图像梯度变化情况

    常见的算子有Robot, Sobel算子,二阶常见多数为拉普拉斯算子,如图所示:


    对于一个1D的有限集合数据f(x) = {1…N}, 假设dx的间隔为1则一阶导数计算公式如下:

    Df(x) = f(x+1) – f(x-1) 二阶导数的计算公式为:df(x)= f(x+1) + f(x-1) – 2f(x);

    稍微难一点的则是基于高斯的一阶导数与二阶导数求取,首先看一下高斯的1D与2D的

    公式。一维高斯对应的X阶导数公式:


    二维高斯对应的导数公式:


    二:算法实现

    1.      高斯采样,基于间隔1计算,计算mask窗口计算,这样就跟普通的卷积计算差不多

    2.      设置sigma的值,本例默认为10,首先计算高斯窗口函数,默认为3 * 3

    3.      根据2的结果,计算高斯导数窗口值

    4.      卷积计算像素中心点值。

    注意点计算高斯函数一定要以零为中心点, 如果窗口函数大小为3,则表达为-1, 0, 1

    三:程序实现关键点

    1.      归一化处理,由于高斯计算出来的窗口值非常的小,必须实现归一化处理。

    2.      亮度提升,对X,Y的梯度计算结果进行了亮度提升,目的是让大家看得更清楚。

    3.      支持一阶与二阶单一方向X,Y偏导数计算

    四:运行效果:

    高斯一阶导数X方向效果


    高斯一阶导数Y方向效果


    五:算法全部源代码:

    /*
     * @author: gloomyfish
     * @date: 2013-11-17
     * 
     * Title - Gaussian fist order derivative and second derivative filter
     */
    package com.gloomyfish.image.harris.corner;
    import java.awt.image.BufferedImage;
    
    import com.gloomyfish.filter.study.AbstractBufferedImageOp;
    
    public class GaussianDerivativeFilter extends AbstractBufferedImageOp {
    
    	public final static int X_DIRECTION = 0;
    	public final static int Y_DIRECTION = 16;
    	public final static int XY_DIRECTION = 2;
    	public final static int XX_DIRECTION = 4;
    	public final static int YY_DIRECTION = 8;
    	
    	// private attribute and settings
    	private int DIRECTION_TYPE = 0;
    	private int GAUSSIAN_WIN_SIZE = 1; // N*2 + 1
    	private double sigma = 10; // default
    
    	public GaussianDerivativeFilter()
    	{
    		System.out.println("高斯一阶及多阶导数滤镜");
    	}	
    	
    	public int getGaussianWinSize() {
    		return GAUSSIAN_WIN_SIZE;
    	}
    
    	public void setGaussianWinSize(int gAUSSIAN_WIN_SIZE) {
    		GAUSSIAN_WIN_SIZE = gAUSSIAN_WIN_SIZE;
    	}
    	public int getDirectionType() {
    		return DIRECTION_TYPE;
    	}
    
    	public void setDirectionType(int dIRECTION_TYPE) {
    		DIRECTION_TYPE = dIRECTION_TYPE;
    	}
    
    	@Override
    	public BufferedImage filter(BufferedImage src, BufferedImage dest) {
    		int width = src.getWidth();
            int height = src.getHeight();
    
            if ( dest == null )
                dest = createCompatibleDestImage( src, null );
    
            int[] inPixels = new int[width*height];
            int[] outPixels = new int[width*height];
            getRGB( src, 0, 0, width, height, inPixels );
            int index = 0, index2 = 0;
            double xred = 0, xgreen = 0, xblue = 0;
            // double yred = 0, ygreen = 0, yblue = 0;
            int newRow, newCol;
            double[][] winDeviationData = getDirectionData();
    
            for(int row=0; row<height; row++) {
            	int ta = 255, tr = 0, tg = 0, tb = 0;
            	for(int col=0; col<width; col++) {
            		index = row * width + col;
            		for(int subrow = -GAUSSIAN_WIN_SIZE; subrow <= GAUSSIAN_WIN_SIZE; subrow++) {
            			for(int subcol = -GAUSSIAN_WIN_SIZE; subcol <= GAUSSIAN_WIN_SIZE; subcol++) {
            				newRow = row + subrow;
            				newCol = col + subcol;
            				if(newRow < 0 || newRow >= height) {
            					newRow = row;
            				}
            				if(newCol < 0 || newCol >= width) {
            					newCol = col;
            				}
            				index2 = newRow * width + newCol;
                            tr = (inPixels[index2] >> 16) & 0xff;
                            tg = (inPixels[index2] >> 8) & 0xff;
                            tb = inPixels[index2] & 0xff;
                        	xred += (winDeviationData[subrow + GAUSSIAN_WIN_SIZE][subcol + GAUSSIAN_WIN_SIZE] * tr);
                        	xgreen +=(winDeviationData[subrow + GAUSSIAN_WIN_SIZE][subcol + GAUSSIAN_WIN_SIZE] * tg);
                        	xblue +=(winDeviationData[subrow + GAUSSIAN_WIN_SIZE][subcol + GAUSSIAN_WIN_SIZE] * tb);
            			}
            		}
            		
            		outPixels[index] = (ta << 24) | (clamp((int)xred) << 16) | (clamp((int)xgreen) << 8) | clamp((int)xblue);
            		
            		// clean up values for next pixel
                    newRow = newCol = 0;
                    xred = xgreen = xblue = 0;
                    // yred = ygreen = yblue = 0;
            	}
            }
    
            setRGB( dest, 0, 0, width, height, outPixels );
            return dest;
    	}
    	
    	private double[][] getDirectionData()
    	{
    		double[][] winDeviationData = null;
            if(DIRECTION_TYPE == X_DIRECTION)
            {
            	winDeviationData = this.getXDirectionDeviation();
            }
            else if(DIRECTION_TYPE == Y_DIRECTION)
            {
            	winDeviationData = this.getYDirectionDeviation();
            }
            else if(DIRECTION_TYPE == XY_DIRECTION)
            {
            	winDeviationData = this.getXYDirectionDeviation();
            }
            else if(DIRECTION_TYPE == XX_DIRECTION)
            {
            	winDeviationData = this.getXXDirectionDeviation();
            }
            else if(DIRECTION_TYPE == YY_DIRECTION)
            {
            	winDeviationData = this.getYYDirectionDeviation();
            }
            return winDeviationData;
    	}
    	
    	public int clamp(int value) {
    		// trick, just improve the lightness otherwise image is too darker...
    		if(DIRECTION_TYPE == X_DIRECTION || DIRECTION_TYPE == Y_DIRECTION)
    		{
    			value = value * 10 + 50;
    		}
    		return value < 0 ? 0 : (value > 255 ? 255 : value);
    	}
    	
    	// centered on zero and with Gaussian standard deviation
    	// parameter : sigma
    	public double[][] get2DGaussianData()
    	{
    		int size = GAUSSIAN_WIN_SIZE * 2 + 1;
    		double[][] winData = new double[size][size];
    		double sigma2 = this.sigma * sigma;
    		for(int i=-GAUSSIAN_WIN_SIZE; i<=GAUSSIAN_WIN_SIZE; i++)
    		{
    			for(int j=-GAUSSIAN_WIN_SIZE; j<=GAUSSIAN_WIN_SIZE; j++)
    			{
    				double r = i*1 + j*j;
    				double sum = -(r/(2*sigma2));
    				winData[i + GAUSSIAN_WIN_SIZE][j + GAUSSIAN_WIN_SIZE] = Math.exp(sum);
    			}
    		}
    		return winData;
    	}
    	
    	public double[][] getXDirectionDeviation()
    	{
    		int size = GAUSSIAN_WIN_SIZE * 2 + 1;
    		double[][] data = get2DGaussianData();
    		double[][] xDeviation = new double[size][size];
    		double sigma2 = this.sigma * sigma;
    		for(int x=-GAUSSIAN_WIN_SIZE; x<=GAUSSIAN_WIN_SIZE; x++)
    		{
    			double c = -(x/sigma2);
    			for(int i=0; i<size; i++)
    			{
    				xDeviation[i][x + GAUSSIAN_WIN_SIZE] = c * data[i][x + GAUSSIAN_WIN_SIZE];				
    			}
    		}
    		return xDeviation;
    	}
    	
    	public double[][] getYDirectionDeviation()
    	{
    		int size = GAUSSIAN_WIN_SIZE * 2 + 1;
    		double[][] data = get2DGaussianData();
    		double[][] yDeviation = new double[size][size];
    		double sigma2 = this.sigma * sigma;
    		for(int y=-GAUSSIAN_WIN_SIZE; y<=GAUSSIAN_WIN_SIZE; y++)
    		{
    			double c = -(y/sigma2);
    			for(int i=0; i<size; i++)
    			{
    				yDeviation[y + GAUSSIAN_WIN_SIZE][i] = c * data[y + GAUSSIAN_WIN_SIZE][i];				
    			}
    		}
    		return yDeviation;
    	}
    	
    	/***
    	 * 
    	 * @return
    	 */
    	public double[][] getXYDirectionDeviation()
    	{
    		int size = GAUSSIAN_WIN_SIZE * 2 + 1;
    		double[][] data = get2DGaussianData();
    		double[][] xyDeviation = new double[size][size];
    		double sigma2 = sigma * sigma;
    		double sigma4 = sigma2 * sigma2;
    		// TODO:zhigang
    		for(int x=-GAUSSIAN_WIN_SIZE; x<=GAUSSIAN_WIN_SIZE; x++)
    		{
    			for(int y=-GAUSSIAN_WIN_SIZE; y<=GAUSSIAN_WIN_SIZE; y++)
    			{
    				double c = -((x*y)/sigma4);
    				xyDeviation[x + GAUSSIAN_WIN_SIZE][y + GAUSSIAN_WIN_SIZE] = c * data[x + GAUSSIAN_WIN_SIZE][y + GAUSSIAN_WIN_SIZE];
    			}
    		}
    		return normalizeData(xyDeviation);
    	}
    	
    	private double[][] normalizeData(double[][] data)
    	{
    		// normalization the data
    		double min = data[0][0];
    		for(int x=-GAUSSIAN_WIN_SIZE; x<=GAUSSIAN_WIN_SIZE; x++)
    		{
    			for(int y=-GAUSSIAN_WIN_SIZE; y<=GAUSSIAN_WIN_SIZE; y++)
    			{
    				if(min > data[x + GAUSSIAN_WIN_SIZE][y + GAUSSIAN_WIN_SIZE])
    				{
    					min = data[x + GAUSSIAN_WIN_SIZE][y + GAUSSIAN_WIN_SIZE];
    				}
    			}
    		}
    		
    		for(int x=-GAUSSIAN_WIN_SIZE; x<=GAUSSIAN_WIN_SIZE; x++)
    		{
    			for(int y=-GAUSSIAN_WIN_SIZE; y<=GAUSSIAN_WIN_SIZE; y++)
    			{
    				data[x + GAUSSIAN_WIN_SIZE][y + GAUSSIAN_WIN_SIZE] = data[x + GAUSSIAN_WIN_SIZE][y + GAUSSIAN_WIN_SIZE] /min;
    			}
    		}
    		
    		return data;
    	}
    	
    	public double[][] getXXDirectionDeviation()
    	{
    		int size = GAUSSIAN_WIN_SIZE * 2 + 1;
    		double[][] data = get2DGaussianData();
    		double[][] xxDeviation = new double[size][size];
    		double sigma2 = this.sigma * sigma;
    		double sigma4 = sigma2 * sigma2;
    		for(int x=-GAUSSIAN_WIN_SIZE; x<=GAUSSIAN_WIN_SIZE; x++)
    		{
    			double c = -((x - sigma2)/sigma4);
    			for(int i=0; i<size; i++)
    			{
    				xxDeviation[i][x + GAUSSIAN_WIN_SIZE] = c * data[i][x + GAUSSIAN_WIN_SIZE];				
    			}
    		}
    		return xxDeviation;
    	}
    	
    	public double[][] getYYDirectionDeviation()
    	{
    		int size = GAUSSIAN_WIN_SIZE * 2 + 1;
    		double[][] data = get2DGaussianData();
    		double[][] yyDeviation = new double[size][size];
    		double sigma2 = this.sigma * sigma;
    		double sigma4 = sigma2 * sigma2;
    		for(int y=-GAUSSIAN_WIN_SIZE; y<=GAUSSIAN_WIN_SIZE; y++)
    		{
    			double c = -((y - sigma2)/sigma4);
    			for(int i=0; i<size; i++)
    			{
    				yyDeviation[y + GAUSSIAN_WIN_SIZE][i] = c * data[y + GAUSSIAN_WIN_SIZE][i];				
    			}
    		}
    		return yyDeviation;
    	}
    
    }
    

    国足都战胜亚洲强队印尼了,我还有什么理由不坚持写下去!

    转载请务必注明!!!
    展开全文
  • 图像的高斯模糊

    2016-11-12 20:55:00
    通过手机获取的图像由于多方面的原因或多或少存在一些噪声,即图像的去噪处理。简单的来说就是用 的矩阵在灰度图像一个一个像素移动,以某种逻辑来消除灰度图像中的孤立点,即噪声。去噪的方法主要有均值滤波、中值...

      通过手机获取的图像由于多方面的原因或多或少存在一些噪声,即图像的去噪处理。简单的来说就是用 的矩阵在灰度图像一个一个像素移动,以某种逻辑来消除灰度图像中的孤立点,即噪声。去噪的方法主要有均值滤波、中值滤波、自适应维纳滤波器、形态学噪声滤除器、高斯滤波器。

      图像高斯滤波是在二维空间利用正态分布(高斯函数)计算平滑模版,并利用该模板与灰度图像做卷积运算来达到滤波去噪的目的。若 是正态分布的标准差,n为模板的大小,则模板上 处元素的计算公式为:

          

    根据正态分布的特性, 越大得到的图像越模糊。从理论上讲图像上灰度不为0的点都应该构成卷积矩阵与原图像做变换,也就是说每个像素的计算都得包括整幅图像。然而在实际中计算高斯函数的离散近似时候, 距离外的像素值可以看作不起作用,也就是模板大小为(6*sigma+1)*(6*sigma+1)

     

    。本文实验选择的是 的模板矩阵, sigma取值越大越模糊。由所得的模板矩阵与原灰度图像做卷积运算即可去除图像中的部分噪声。

       因此高斯模糊一般分为以下四步:

    1、获取高斯模板矩阵:

        public float[][] get2DKernalData(int n, float sigma) {
                int size = 2 * n + 1;
                float sigma22 = 2 * sigma * sigma;
                float sigma22PI = (float) Math.PI * sigma22;
                float[][] kernalData = new float[size][size];
                int row = 0;
                for (int i = -n; i <= n; i++) {
                    int column = 0;
                    for (int j = -n; j <= n; j++) {
                        float xDistance = i * i;
                        float yDistance = j * j;
                        kernalData[row][column] = (float) Math
                                .exp(-(xDistance + yDistance) / sigma22) / sigma22PI;
                        column++;
                    }
                    row++;
                }
    
                return kernalData;
            }
    View Code

    2、获取灰度矩阵

    public int[][]getGrayMatrix(Bitmap bt)
        {
            int hdjz[][]=new int[h][w];
            for(int i=0;i<h;i++)
                for(int j=0;j<w;j++)
                {
            int argb=bt.getPixel(j , i);
            int r = (argb >> 16) & 0xFF;
            int g = (argb >> 8) & 0xFF;
            int b = (argb >> 0) & 0xFF;
            int grayPixel = (int) (r + g + b) / 3;
            hdjz[i][j]=grayPixel;
                }
            return hdjz;
        }
    View Code

    3、卷积运算

    public int[][] GS(int[][] hd, int size, float sigma) {
            float[][] gs = get2DKernalData(size, sigma);
            int outmax = 0;
            int inmax = 0;
            for (int x = size; x < w - size; x++)
                for (int y = size; y < h - size; y++) {
                    float hc1 = 0;
                    if (hd[y][x] > inmax)
                        inmax = hd[y][x];
                    for (int k = -size; k < size + 1; k++)
                        for (int j = -size; j < size + 1; j++) {
                            hc1 = gs[size + k][j + size] * hd[y + j][x + k] + hc1;
    
                        }
                    hd[y][x] = (int) (hc1);
                    if (outmax < hc1)
                        outmax = (int) (hc1);
                }
            float rate = inmax / outmax;
    
            for (int x = size; x < w - size; x++)
                for (int y = size; y < h - size; y++) {
                    hd[y][x] = (int) (hd[y][x] * rate);
                }
            
            return hd;
        }    
    View Code

    4、将得到的灰度矩阵创建灰度图

    //由灰度矩阵创建灰度图
        public Bitmap createGrayImage(int[][]hdjz)
        {
            int h=hdjz.length;
            int w = hdjz[0].length;
            Bitmap bt=Bitmap.createBitmap(w, h, image.getConfig());
            for(int i=0;i<h;i++)
                for(int j=0;j<w;j++)
                {
                    int grayValue=hdjz[i][j];
                    int color = ((0xFF << 24)+(grayValue << 16)+(grayValue << 8)+grayValue);
                    bt.setPixel(j, i, color);
                }
            return bt;
        }
    View Code

     

    转载于:https://www.cnblogs.com/bokeofzp/p/6057546.html

    展开全文
  • 图像的高斯模糊是非常经典的图像卷积例子。本质上,图像模糊就是将(灰度)图像 I 和一个高斯核进行卷积操作: 其中 * 表示卷积操作;是标准差为 σ 的二维高斯核,定义为 : 这里提到了卷积操作,那么什么是卷积呢...
  • 在计算机视觉中,有时也简称为高斯函数。高斯函数具有五个重要的性质,这些性质使得它在早期图像处理中特别有用.这些性质表明,高斯平滑滤波器无论在空间域还是...一般来说,一幅图像的边缘方向是事先不知道的,因此,
  • 然后,我们来看一下导数的计算,在离散空间中,我们一般会用f(k) - f(k-1)/(k- (k-1) )来计算导数。 现在我们介绍另外一种基于卷积特性的导数计算方法,由于卷积具有以下性质,[2] [2] 因而,在图像处理中,...
  • 高斯图像模糊算法及其 C 实现

    千次阅读 2014-01-03 12:22:09
    高斯模糊的基本思路是根据二维 正太分布 正态分布 (感谢 xhr 大牛指正错别字) 公式生成一个高斯矩阵, 求新图像中的每一点时, 将高斯矩阵的中心对准旧图像的这一点, 并将所有点根据高斯矩阵上对应的点加权平均....
  • 高斯模糊高斯模糊(英语:Gaussian Blur),也叫高斯平滑,是在Adobe Photoshop、GIMP以及Paint.NET等图像处理软件中广泛使用处理效果,通常用它来减少图像杂讯以及降低细节层次。这种模糊技术生成的图像,其视觉...
  • 高斯模糊高斯模糊(英语:Gaussian Blur),也叫高斯平滑,是在Adobe Photoshop、GIMP以及Paint.NET等图像处理软件中广泛使用处理效果,通常用它来减少图像杂讯以及降低细节层次。这种模糊技术生成的图像,其视觉...
  • 图像的高斯模糊

    2015-04-17 11:03:10
    高斯模糊是一种图像滤波器,它使用正态分布(高斯函数)计算模糊模板,并使用该模板与原图像做卷积运算,达到模糊图像的目的。 N维空间正态分布方程为: 其中,σ是正态分布的标准差,σ值越大,图像越模糊...
  • 很好matlab程序关于高斯噪声详细算法。
  • 将彩色图像转换为NTSC灰度图像进行高斯滤波,采用σ=1,2,3,5。任 选一种像素填补方案。对于σ=1 下结果,与直接调用相关函数结果进行比较(可以简单计算差值图像)。然后,任选两幅图像,比较其他参数条件...
  • 高斯模糊高斯模糊(英语:Gaussian Blur),也叫高斯平滑,是在Adobe Photoshop、GIMP以及Paint.NET等图像处理软件中广泛使用处理效果,通常用它来减少图像杂讯以及降低细节层次。这种模糊技术生成的图像,其视觉...
  • 同时,进一步说明由于这两种矩阵构成元素稀疏性可以简化图像重建过程中投影计算,从而提高重建速度.实验结果表明新测量矩阵均有较好测量效果,在满足一定测量数目要求条件下可以精确重建.最后给出了这两种...
  • 图像金字塔技术在很多层面上都有着广泛应用,很多开源工具也都有对他们建立写了专门函数,比如IPP,比如OpenCV等等,这方面理论文章特别多,我不需要赘述,但是我发现大部多分开源代码实现都不是严格...
  • 基于高斯模型彩色图像反向投影

    千次阅读 2017-06-14 21:08:39
    一:介绍图像反向投影最终目的是获取ROI然后实现对ROI区域标注、识别、测量等图像处理与分析,是计算机视觉与人工智能常见方法之一。图像反向投影通常是彩色图像投影效果会比灰度图像效果要好,原因在于彩色...
  • 高斯模糊高斯模糊(英语:Gaussian Blur),也叫高斯平滑,是在Adobe Photoshop、GIMP以及Paint.NET等图像处理软件中广泛使用处理效果,通常用它来减少图像杂讯以及降低细节层次。这种模糊技术生成的图像,其视觉...

空空如也

空空如也

1 2 3 4 5 ... 20
收藏数 1,912
精华内容 764
关键字:

高斯图像的计算