精通图像特征提取

Github个人博客：https://joeyos.github.io

在现实生活中，我们要交往一个人，首先要了解这个人的特征，相貌和性格；同样，在基于数字图像的模式识别中，我们也要提取图像的特征。对图像特征的提取与研究提供了一种具有统计意义的图像内容表达，正所谓：“透过现象看本质”。

• 面特征：图像金字塔，图像矩特征
• 线特征：边缘检测
• 局部区域特征：斑点特征检测
• 点特征：角点检测
• 不变点特征：尺度不变特征提取
• 尺度不变特征变换SIFT算法、解压偶数页鲁棒性尺度不变特征提取SURF算法

图像金字塔

图像金字塔：图像经过一个低通滤波器进行平滑，然后对这个平滑图像进行抽样，抽样比例一般在水平和垂直方向上都为1/2,从而得到一系列尺寸缩小、分辨率降低的图像。

% 读入图像并灰度化
I = imread('test.jpg');
I = rgb2gray(I);
% 生成高斯滤波器的核
w = fspecial('gaussian',3,0.5);
size_a = size(I);
% 进行高斯滤波
g = imfilter(I,w,'conv','symmetric','same');
% 降采样
t = g(1:2:size_a(1),1:2:size_a(2));
% 显示结果
figure(1);imshow(I);
figure(2);imshow(t);

• 高斯金字塔分解：低通滤波和降采样
• 拉普拉斯金字塔分解：相当于高斯金字塔与其高一级图像内插放大后图像的差，相当于带通滤波。

function [ pyr ] = genPyr( img, type, level )
%   功能： 产生图像的高斯金字塔或拉普拉斯金字塔
%   输入：img-灰度图像；
%         type-变换的类型：'gauss' or 'laplace'；
%         level-分解层数
%   输出：PYR-1*LEVEL 元胞数组

pyr = cell(1,level);
pyr{1} = im2double(img);
%%%%%%图像的金字塔化%%%%%%

%计算图像的高斯金字塔
for p = 2:level
    pyr{p} = pyr_reduce(pyr{p-1});
end
if strcmp(type,'gauss'), return; end
 
% 调整图像尺寸
for p = level-1:-1:1
    osz = size(pyr{p+1})*2-1;
    pyr{p} = pyr{p}(1:osz(1),1:osz(2),:);
end

% 计算图像的拉普拉斯金字塔
for p = 1:level-1
    pyr{p} = pyr{p}-pyr_expand(pyr{p+1});
end
 
end

function [ imgout ] = pyr_reduce( img )
% 功能：计算图像的高斯金字塔
% 输入：img-灰度图像；
% 输出：imgout-分解后的金子塔图像数组

% 生成高斯核
kernelWidth = 5; 
cw = .375; 
ker1d = [.25-cw/2 .25 cw .25 .25-cw/2];
kernel = kron(ker1d,ker1d');
 
img = im2double(img);
sz = size(img);
imgout = [];
 
% 产生图像的高斯金字塔
for p = 1:size(img,3)
    img1 = img(:,:,p);
    imgFiltered = imfilter(img1,kernel,'replicate','same');
    imgout(:,:,p) = imgFiltered(1:2:sz(1),1:2:sz(2));
end
end

function [ imgout ] = pyr_expand( img )
% 功能：图像金字塔扩张
% 输入：img-待扩张的图像
% 输出：imgout-扩张后的图像

kw = 5; 
cw = .375; 
ker1d = [.25-cw/2 .25 cw .25 .25-cw/2];
kernel = kron(ker1d,ker1d')*4;
 
ker00 = kernel(1:2:kw,1:2:kw); 
ker01 = kernel(1:2:kw,2:2:kw); 
ker10 = kernel(2:2:kw,1:2:kw); 
ker11 = kernel(2:2:kw,2:2:kw); 
 
img = im2double(img);
sz = size(img(:,:,1));
osz = sz*2-1;
imgout = zeros(osz(1),osz(2),size(img,3));
 
for p = 1:size(img,3)
    img1 = img(:,:,p);
img1ph = padarray(img1,[0 1],'replicate','both');
img1pv = padarray(img1,[1 0],'replicate','both');    
    img00 = imfilter(img1,ker00,'replicate','same');
    img01 = conv2(img1pv,ker01,'valid'); 
    img10 = conv2(img1ph,ker10,'valid');
    img11 = conv2(img1,ker11,'valid');
    
    imgout(1:2:osz(1),1:2:osz(2),p) = img00;
    imgout(2:2:osz(1),1:2:osz(2),p) = img10;
    imgout(1:2:osz(1),2:2:osz(2),p) = img01;
    imgout(2:2:osz(1),2:2:osz(2),p) = img11;
end
end

%% 金字塔分解函数实现
close all;clear all;clc;
img=imread('qingdao.jpg');
n=3;
%[ pyr ] = genPyr( img, 'gauss', n );%高斯金字塔分解
[ pyr ] = genPyr( img, 'laplace', n );%拉普拉斯金字塔分解
for i=1:n
    figure(i);
    imshow(pyr{i});
end

图像的矩特征

由于被识别的图像与原图像相比一般有很大程度的失真，如平移、旋转和其他变化。不变矩是一种高度浓缩的图像特征，具有平移、灰度、尺度、旋转不变性，因此矩和矩函数被广泛用于图像的模式识别、图像分类、目标识别和场景分析中。M.K.Hu在1961年首次提出不变矩的概念，并将几何矩用于图像描述，但其低阶几何矩与图像的整体特征有关，不包含太多的图像细节信息，而高阶几何矩易受噪声影响，因此，很难用几何矩恢复图像。Zernike矩能够很容易的构造图像的任意高阶矩，并能使用较少的矩来重建图像。

• 计算Hu矩

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];   

• 计算Zernike矩
• 首先确定图像函数的大小N*N，从而确定N值
• 进一步确定r和delta的范围
• 利用Zernike多项式的快速递推性质计算各阶Rnm,Cnm和Snm
• 对Cnm和Snm求模，求得|Znm|

function        [A_nm,zmlist,cidx,V_nm]        = zernike(img,n