几何不变矩--Hu矩

た入场券 2023-08-17 16:15 102阅读 0赞

【图像算法】图像特征：

**-------------------------------------------------------------------------------------------------------------------------------**

一 原理

几何矩是由Hu(Visual pattern recognition by moment invariants)在1962年提出的，具有平移、旋转和尺度不变性。 定义如下：

① (**p+q)阶不变矩定义**：

![2011071909102486.png][]

**② 对于数字图像，离散化，定义为**：

![2011071909105531.png][]

**③ 归一化中心矩定义**：

![2011071909121175.png][]

**④Hu矩定义**

![2011071909130059.png][]

![2011071909131361.png][]

**-------------------------------------------------------------------------------------------------------------------------------**

**二 实现(源码)**

**①[自编函数模块][Link 1]C++**

//#################################################################################//
    double M[7] = {0}; //HU不变矩
    bool HuMoment(IplImage* img)
    {
    int bmpWidth = img->width;
    int bmpHeight = img->height;
    int bmpStep = img->widthStep; 
    int bmpChannels = img->nChannels;
    uchar*pBmpBuf = (uchar*)img->imageData;
    
    double m00=0,m11=0,m20=0,m02=0,m30=0,m03=0,m12=0,m21=0; //中心矩 
    double x0=0,y0=0; //计算中心距时所使用的临时变量（x-x'） 
    double u20=0,u02=0,u11=0,u30=0,u03=0,u12=0,u21=0;//规范化后的中心矩
    //double M[7]; //HU不变矩 
    double t1=0,t2=0,t3=0,t4=0,t5=0;//临时变量， 
    //double Center_x=0,Center_y=0;//重心 
    int Center_x=0,Center_y=0;//重心 
    int i,j; //循环变量
    
    // 获得图像的区域重心(普通矩)
    double s10=0,s01=0,s00=0; //0阶矩和1阶矩 
    for(j=0;j<bmpHeight;j++)//y
    {
    for(i=0;i<bmpWidth;i++)//x
    {
    s10+=i*pBmpBuf[j*bmpStep+i];
    s01+=j*pBmpBuf[j*bmpStep+i];
    s00+=pBmpBuf[j*bmpStep+i];
    }
    }
    Center_x=(int)(s10/s00+0.5);
    Center_y=(int)(s01/s00+0.5);
    
    // 计算二阶、三阶矩(中心矩)
    m00=s00; 
    for(j=0;j<bmpHeight;j++) 
    {
    for(i=0;i<bmpWidth;i++)//x 
    { 
    x0=(i-Center_x); 
    y0=(j-Center_y); 
    m11+=x0*y0*pBmpBuf[j*bmpStep+i]; 
    m20+=x0*x0*pBmpBuf[j*bmpStep+i]; 
    m02+=y0*y0*pBmpBuf[j*bmpStep+i]; 
    m03+=y0*y0*y0*pBmpBuf[j*bmpStep+i];
    m30+=x0*x0*x0*pBmpBuf[j*bmpStep+i]; 
    m12+=x0*y0*y0*pBmpBuf[j*bmpStep+i]; 
    m21+=x0*x0*y0*pBmpBuf[j*bmpStep+i]; 
    } 
    }
    
    // 计算规范化后的中心矩: mij/pow(m00,((i+j+2)/2)
    u20=m20/pow(m00,2); 
    u02=m02/pow(m00,2); 
    u11=m11/pow(m00,2);
    u30=m30/pow(m00,2.5); 
    u03=m03/pow(m00,2.5);
    u12=m12/pow(m00,2.5); 
    u21=m21/pow(m00,2.5);
    
    // 计算中间变量
    t1=(u20-u02); 
    t2=(u30-3*u12); 
    t3=(3*u21-u03); 
    t4=(u30+u12);
    t5=(u21+u03);
    
    // 计算不变矩 
    M[0]=u20+u02; 
    M[1]=t1*t1+4*u11*u11; 
    M[2]=t2*t2+t3*t3; 
    M[3]=t4*t4+t5*t5;
    M[4]=t2*t4*(t4*t4-3*t5*t5)+t3*t5*(3*t4*t4-t5*t5); 
    M[5]=t1*(t4*t4-t5*t5)+4*u11*t4*t5;
    M[6]=t3*t4*(t4*t4-3*t5*t5)-t2*t5*(3*t4*t4-t5*t5);
    
    returntrue;
    }

**②调用OpenCV方法**

![复制代码][copycode.gif]

1 //  利用OpenCV函数求7个Hu矩
    2     CvMoments moments;
    3     CvHuMoments hu;
    4     cvMoments(bkImgEdge,&moments,0); 
    5     cvGetHuMoments(&moments, &hu);
    6     cout<<hu.hu1<<"/"<<hu.hu2<<"/"<<hu.hu3<<"/"<<hu.hu4<<"/"<<hu.hu5<<"/"<<hu.hu6<<"/"<<hu.hu7<<"/"<<"/"<<endl;
    7     cvMoments(testImgEdge,&moments,0); 
    8     cvGetHuMoments(&moments, &hu);
    9     cout<<hu.hu1<<"/"<<hu.hu2<<"/"<<hu.hu3<<"/"<<hu.hu4<<"/"<<hu.hu5<<"/"<<hu.hu6<<"/"<<hu.hu7<<"/"<<"/"<<endl;

![复制代码][copycode.gif]

Python调用OpenCV：

#-*-coding:utf-8-*-
    import cv2
    from datetime import datetime
    import numpy as np
    
    def test(img):
        moments = cv2.moments(img)
        humoments = cv2.HuMoments(moments)
        # humoments = no.log(np.abs(humoments)) # 同样建议取对数
        print(humoments)
    
    if __name__ == '__main__':
        t1 = datetime.now()    
        fp = '/home/mamq/images/3.jpg'
        img = cv2.imread(fp)
        img_gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
        test(img_gray)
        print datetime.now() - t1

Python方法：

#-*-coding:utf-8-*-
    import cv2
    from datetime import datetime
    import numpy as np
    np.set_printoptions(suppress=True)
    
    def humoments(img_gray):
        '''
        由于7个不变矩的变化范围很大,为了便于比较,可利用取对数的方法进行数据压缩;同时考虑到不变矩有可能出现负值的情况,因此,在取对数之前先取绝对值
        经修正后的不变矩特征具有平移 、旋转和比例不变性
        '''
        # 标准矩定义为m_pq = sumsum(x^p * y^q * f(x, y))
        row, col = img_gray.shape
        #计算图像的0阶几何矩
        m00 = img_gray.sum()
        m10 = m01 = 0
        #　计算图像的二阶、三阶几何矩
        m11 = m20 = m02 = m12 = m21 = m30 = m03 = 0
        for i in range(row):
            m10 += (i * img_gray[i]).sum()
            m20 += (i ** 2 * img_gray[i]).sum()
            m30 += (i ** 3 * img_gray[i]).sum()
            for j in range(col):
                m11 += i * j * img_gray[i][j]
                m12 += i * j ** 2 * img_gray[i][j]
                m21 += i ** 2 * j * img_gray[i][j]
        for j in range(col):
            m01 += (j * img_gray[:, j]).sum()
            m02 += (j ** 2 * img_gray[:, j]).sum()
            m30 += (j ** 3 * img_gray[:, j]).sum()
        # 由标准矩我们可以得到图像的"重心"
        u10 = m10 / m00
        u01 = m01 / m00
        # 计算图像的二阶中心矩、三阶中心矩
        y00 = m00
        y10 = y01 = 0
        y11 = m11 - u01 * m10
        y20 = m20 - u10 * m10
        y02 = m02 - u01 * m01
        y30 = m30 - 3 * u10 * m20 + 2 * u10 ** 2 * m10
        y12 = m12 - 2 * u01 * m11 - u10 * m02 + 2 * u01 ** 2 * m10
        y21 = m21 - 2 * u10 * m11 - u01 * m20 + 2 * u10 ** 2 * m01
        y03 = m03 - 3 * u01 * m02 + 2 * u01 ** 2 * m01
        # 计算图像的归格化中心矩
        n20 = y20 / m00 ** 2
        n02 = y02 / m00 ** 2
        n11 = y11 / m00 ** 2
        n30 = y30 / m00 ** 2.5
        n03 = y03 / m00 ** 2.5
        n12 = y12 / m00 ** 2.5
        n21 = y21 / m00 ** 2.5
        # 计算图像的七个不变矩
        h1 = n20 + n02
        h2 = (n20 - n02) ** 2 + 4 * n11 ** 2
        h3 = (n30 - 3 * n12) ** 2 + (3 * n21 - n03) ** 2
        h4 = (n30 + n12) ** 2 + (n21 + n03) ** 2
        h5 = (n30 - 3 * n12) * (n30 + n12) * ((n30 + n12) ** 2 - 3 * (n21 + n03) ** 2) + (3 * n21 - n03) * (n21 + n03) \
            * (3 * (n30 + n12) ** 2 - (n21 + n03) ** 2)
        h6 = (n20 - n02) * ((n30 + n12) ** 2 - (n21 + n03) ** 2) + 4 * n11 * (n30 + n12) * (n21 + n03)
        h7 = (3 * n21 - n03) * (n30 + n12) * ((n30 + n12) ** 2 - 3 * (n21 + n03) ** 2) + (3 * n12 - n30) * (n21 + n03) \
            * (3 * (n30 + n12) ** 2 - (n21 + n03) ** 2)
        inv_m7 = [h1, h2, h3, h4, h5, h6, h7]
        inv_m7 = np.log(np.abs(inv_m7))
        return inv_m7
    
    if __name__ == '__main__':
        t1 = datetime.now()
        fp = '/home/mamq/images/3.jpg'
        img = cv2.imread(fp)
        img_gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
        print humoments(img_gray)
        print datetime.now() - t1

**MATLAB方法：**

invariable_moment(imread('lena.jpg'));
    
    function inv_m7 = invariable_moment(in_image)
    % 功能：计算图像的Hu的七个不变矩
    % 输入：in_image-RGB图像
    % 输出：inv_m7-七个不变矩
     
    % 将输入的RGB图像转换为灰度图像   
    image=rgb2gray(in_image);     
    %将图像矩阵的数据类型转换成双精度型
    image=double(image);      
    %%%=================计算 、 、 =========================
    %计算灰度图像的零阶几何矩 
    m00=sum(sum(image));     
    m10=0;
    m01=0;
    [row,col]=size(image);
    for i=1:row
        for j=1:col
            m10=m10+i*image(i,j);
            m01=m01+j*image(i,j);
        end
    end
    %%%=================计算 、 ================================
    u10=m10/m00;
    u01=m01/m00;
    %%%=================计算图像的二阶几何矩、三阶几何矩============
    m20 = 0;m02 = 0;m11 = 0;m30 = 0;m12 = 0;m21 = 0;m03 = 0;
    for i=1:row
        for j=1:col
            m20=m20+i^2*image(i,j);
            m02=m02+j^2*image(i,j);
            m11=m11+i*j*image(i,j);
            m30=m30+i^3*image(i,j);
            m03=m03+j^3*image(i,j);
            m12=m12+i*j^2*image(i,j);
            m21=m21+i^2*j*image(i,j);
        end
    end
    %%%=================计算图像的二阶中心矩、三阶中心矩============
    y00=m00;
    y10=0;
    y01=0;
    y11=m11-u01*m10;
    y20=m20-u10*m10;
    y02=m02-u01*m01;
    y30=m30-3*u10*m20+2*u10^2*m10;
    y12=m12-2*u01*m11-u10*m02+2*u01^2*m10;
    y21=m21-2*u10*m11-u01*m20+2*u10^2*m01;
    y03=m03-3*u01*m02+2*u01^2*m01;
    %%%=================计算图像的归格化中心矩====================
            n20=y20/m00^2;
            n02=y02/m00^2;
            n11=y11/m00^2;
            n30=y30/m00^2.5;
            n03=y03/m00^2.5;
            n12=y12/m00^2.5;
            n21=y21/m00^2.5;
    %%%=================计算图像的七个不变矩======================
    h1 = n20 + n02;                      
    h2 = (n20-n02)^2 + 4*(n11)^2;
    h3 = (n30-3*n12)^2 + (3*n21-n03)^2;  
    h4 = (n30+n12)^2 + (n21+n03)^2;
    h5 = (n30-3*n12)*(n30+n12)*((n30+n12)^2-3*(n21+n03)^2)+(3*n21-n03)*(n21+n03)*(3*(n30+n12)^2-(n21+n03)^2);
    h6 = (n20-n02)*((n30+n12)^2-(n21+n03)^2)+4*n11*(n30+n12)*(n21+n03);
      h7 = (3*n21-n03)*(n30+n12)*((n30+n12)^2-3*(n21+n03)^2)+(3*n12-n30)*(n21+n03)*(3*(n30+n12)^2-(n21+n03)^2);
     
    inv_m7= [h1 h2 h3 h4 h5 h6 h7];

**-------------------------------------------------------------------------------------------------------------------------------**

**三 相似性准则**

**①法一**

//  计算相似度1
    double dbR =0; //相似度
    double dSigmaST =0;
        double dSigmaS =0;
        double dSigmaT =0;
        double temp =0;  
        {for(int i=0;i<7;i++)
        {
            temp = fabs(Sa[i]*Ta[i]);
            dSigmaST+=temp;
            dSigmaS+=pow(Sa[i],2);
            dSigmaT+=pow(Ta[i],2);
        }}
        dbR = dSigmaST/(sqrt(dSigmaS)*sqrt(dSigmaT));

**②法二**

1 //  计算相似度2
     2 double dbR2 =0; //相似度
     3 double temp2 =0;
     4 double temp3 =0;  
     5     {for(int i=0;i<7;i++)
     6     {
         
     7         temp2 += fabs(Sa[i]-Ta[i]);
     8         temp3 += fabs(Sa[i]+Ta[i]);
     9     }}
    10     dbR2 =1- (temp2*1.0)/(temp3);

**-------------------------------------------------------------------------------------------------------------------------------**

*Author:         [SKySeraph][]*

*Email/GTalk: [zgzhaobo@gmail.com][zgzhaobo_gmail.com]    QQ：452728574*

*From：         [http://www.cnblogs.com/skyseraph/][SKySeraph]*

** ****-------------------------------------------------------------------------------------------------------------------------------**

作者：[skyseraph][SKySeraph]   
出处：[http://www.cnblogs.com/skyseraph/][SKySeraph]   
更多精彩请直接访问SkySeraph个人站点：[http://skyseraph.com//][http_skyseraph.com]   
Email/GTalk: zgzhaobo@gmail.com   
本文版权归作者和博客园共有，欢迎转载，但未经作者同意必须保留此段声明，且在文章页面明显位置给出原文连接，否则保留追究法律责任的权利。

![wechat_sk.jpg][]

转载于:https://www.cnblogs.com/sddai/p/11430732.html

[2011071909102486.png]: /images/20230810/c98b1895f4a4488c9a45fb66be55b618.png
[2011071909105531.png]: /images/20230810/3cb8aaf113414eac8efcc14ba4e4fcd9.png
[2011071909121175.png]: /images/20230810/c61e5e8a20cc48678162f6a0c3020581.png
[2011071909130059.png]: /images/20230810/c6c74df6c9a64896987eb66c27e2e5d9.png
[2011071909131361.png]: /images/20230810/db7a9f1d0d8145949b4b2b92f7845201.png
[Link 1]: http://www.cnblogs.com/skyseraph/archive/2011/03/29/1998432.html
[copycode.gif]: /images/20230810/ebaaca87a7cf44d29fece23ebb37ecf8.png
[SKySeraph]: http://www.cnblogs.com/skyseraph/
[zgzhaobo_gmail.com]: http://gmail.com/
[http_skyseraph.com]: http://skyseraph.com/
[wechat_sk.jpg]: http://7xo4q8.com1.z0.glb.clouddn.com/skyseraph/2016/wechat_sk.jpg