分享一下最近学习的图像分类方面知识，整体的思路如下（之前的汇报ppt里截的）
把这个过程拆分几个部分共同学习一下吧
1.Otsu法原理
最大类间方差法OTSU是一种自适应的全局阈值确定的方法，根据灰度阈值T将图像分成背景和目标两部分，像素点数占比为ω_0和ω_1，平均灰度值为μ_0和μ_1
由已知，
图像平均灰度
类间方差
代入式(2-1)、(2-2)可得目标函数为
背景和目标之间的类间方差最大时，目标和背景的错分概率最小
2.原始Otsu算法与改进算法的实现
改进算法思想是在原始Otsu算法基础上求取最大类间、类内方差比值最小的情况。
在分割之前利用了一些灰度变换的方法做图像的预处理，详细内容可以看链接
% 原始的与改进的OTSU算法
clc,clear,close all
f = imread('C:\Users\HS\Desktop\duck.jpg');
f = rgb2gray(f);f1 = f;
figure,imshow(f);title('原图的灰度化');
set(gcf,'color','w');set(gca,'linewidth',1,'fontsize',16);
figure,histogram(f);title('原图的灰度直方图');
set(gcf,'color','w');set(gca,'linewidth',1,'fontsize',16);
g = imadjust(f,[0.1,0.85],[],0.5);%mamma<1时，变换凸显高灰度区域
ng = imnoise(g,'gaussian',0,0.1);
[frest,psnr] = get_winer(g,ng); % 维纳滤波
figure();imshow(frest); title(['维纳滤波去噪,PSNR=',num2str(psnr)]);
set(gcf,'color','w');set(gca,'linewidth',1,'fontsize',16);
figure,histogram(frest);title('0.1高斯噪声滤波后的灰度直方图');
axis([0,255,0,7000]);
set(gca,'xtick',0:50:255);
set(gca,'ytick',0:1000:7000);
set(gcf,'color','w');set(gca,'linewidth',1,'fontsize',16);
%% 正式开始：传统的OTSU算法
[thresh,SM] = graythresh(frest); % OTSU算法，自适应选取二值化阈值
BWimg = im2bw(frest,thresh);%im2bw的结果是白底黑目标
figure;
imshow(BWimg);title('OTSU算法二值化结果');
set(gcf,'color','w');set(gca,'linewidth',1,'fontsize',16);
figure,histogram(BWimg);title('OTSU算法二值化结果的灰度直方图');
set(gcf,'color','w');set(gca,'linewidth',1,'fontsize',16);
%% 最大类间、类内比值法
[m,n] = size(f);
th = zeros(1,254); % 准备存放比值
f = double(f); % 必须转换成double格式！！！！否则不能计算像素！！！
% fi = 0;
% for i=1:m
% for j=1:n
% fprintf('fi=%d,',fi);
% fprintf('f(%d,%d)=%d\n',i,j,f(i,j));
% fi = fi+f(i,j);
%
% end
% end
for T = 1:254
s1 = 0;
s2 = 0;
c1 = 0;
c2 = 0;
q1 = 0;
q2 = 0;
u1 = 0;
u2 = 0;
p1 = 0;
p2 = 0; % 初始化
for i = 1:m
for j = 1:n
if f(i,j)>=T
s1 = f(i,j)+s1; % 这一类的像素点灰度值求和
c1 = c1+1; % 这一类的像素点数量
else
s2 = f(i,j)+s2;
c2 = c2+1;
end
end
end
u1 = s1/c1; % 这一类的的灰度均值u1
u2 = s2/c2;
p1 = c1/(c1+c2); % 这一类的分布概率
p2 = c2/(c1+c2);
u = (s1+s2)/(c1+c2); % 总体均值
x = p1*(u-u1)*(u-u1)-p2*(u-u2)*(u-u2); % 类间方差
for i = 1:m
for j = 1:n
if f(i,j)>=T
q1 = q1+(f(i,j)-u1)*(f(i,j)-u1); % 这一类的灰度方差
else
q2 = q2+(f(i,j)-u2)*(f(i,j)-u2);
end
end
end
y = p1*q1+p2*q2; % 类内方差
c = x/y; % 类间、类内方差的比值
th(T) = c;
end
d = find(th==max(th(:))); % 寻找最大比值对应的位置，即T的值
Th = d; % 最优阈值
b = zeros(m,n);
for i = 1:m
for j = 1:n
if f(i,j)>=Th
b(i,j) = 0;
else
b(i,j) = 255;
end
end
end
figure,imshow(b);title('最大类间、类内方差比值法所得二值化图像');
set(gcf,'color','w');set(gca,'linewidth',1,'fontsize',16);
figure,histogram(b);title('新Otsu法的灰度直方图');
set(gcf,'color','w');set(gca,'linewidth',1,'fontsize',16);
figure,subplot(231); imshow(f1);title('原图灰度化结果');
set(gcf,'color','w');set(gca,'linewidth',1,'fontsize',16);
subplot(232); imshow(BWimg);title('OTSU算法二值化结果');
set(gcf,'color','w');set(gca,'linewidth',1,'fontsize',16);
subplot(233); imshow(b);title('改进的OTSU算法二值化结果');
set(gcf,'color','w');set(gca,'linewidth',1,'fontsize',16);
subplot(234);histogram(f1);title('原图的灰度直方图');
set(gcf,'color','w');set(gca,'linewidth',1,'fontsize',16);
subplot(235);histogram(BWimg);title('OTSU算法二值化结果的灰度直方图');
set(gcf,'color','w');set(gca,'linewidth',1,'fontsize',16);
subplot(236);histogram(b);title('改进的OTSU算法二值化结果的灰度直方图');
set(gcf,'color','w');set(gca,'linewidth',1,'fontsize',16);
结果是这样的
似乎原始Otsu算法已经可以实现把目标及其影子分割出来了呢，改进的算法可以把影子也剔除掉
3.基于OTSU方法的灰度图像二值化，提取ROI及其形状特征
对文件夹下图像统一操作——Otsu分割、形态学操作、提取roi，最后计算形状特征
说明几点：
1.势必会提取到多个ROI，对单亮点模型，我设置了只保留外接矩形面积最大的那一个
2.二值化的结果必须是前景为白色，背景为黑色，否则下一步瘪瘪出错
3.形状特征是用以分类的主要特征，包括但不限于我这里计算的四维–面积、周长、离心率、zernike矩
clc,clear,close all
file_path = 'C:\Users\HS\Desktop\图片素材\滤波后';% 获取指定文件夹下所有子文件的路径
SCALE = [];% 存放形状特征
ALL_EM = [];% 存放分割率
ALL_num = [];% 存放ROI数量
p=genpath(file_path);
%这些路径存在于字符串p中，以';'分割
length_p=size(p,2); % 字符串p的长度
path={}; % 建立一个元胞数组，数组的每个单元中包含一个目录
temp=[];
for i=1:length_p % 寻找分隔符';'，一旦找到，则将路径temp写入path数组中
if p(i)~=';'
temp=[temp p(i)];
else
temp=[temp '\'];%在路径的后面加上\
path=[path;temp];
temp=[];
end
end
clear p length_p temp
%获得了指定文件夹及其所有子文件夹，接下来就是逐一文件夹读取图像，利用上面的path数组，将路径取出来
file_num=size(path,1);%获取子文件夹的个数
for i=1:file_num
file_path = path{i};
img_path_list=dir(strcat(file_path,'*.jpg')); %获取每个文件夹下的图像列表
img_num=length(img_path_list);%获取文件夹下图像的数量
if img_num > 0
for j=1:img_num
image_name=img_path_list(j).name;%获取图像名
%new_image_name = image_name(1:end-4);%去掉原图名中的.jpg
SIFT_imageFile = strcat(file_path,image_name);
image=imread(SIFT_imageFile);
%image=imread(strcat(file_path,image_name));%读取文件夹下的每一帧图像
fprintf('第%d个文件夹中第%d个图片 %s\n',i,j,strcat(file_path,image_name));%显示正在读取的文件夹和其中图像数量
%% 正式操作
% 根据滤波后图像二值化提取尺度特征scale
[EM,num,scale] = getScaleFeatures(image);
SCALE = [SCALE;scale];
ALL_EM = [ALL_EM,EM]; % 窥探图片分割率如何
ALL_num = [ALL_num,num];% 分别划分了几个ROI
end
end
end
fprintf('完成！记得保存好每一类目标的特征矩阵！！\n');
其中调用了几个小函数
%% 基于OTSU方法的灰度图像二值化，提取尺度特征
function [EM,num,scale] = getScaleFeatures(frest)
% 对调整好的灰度图片frest，返回尺度特征矩阵scale、分割率EM、ROI数量num
[thresh,EM] = graythresh(frest); % OTSU算法，自适应选取二值化阈值
% if EM>0.7
% fprintf('EM=%f, 适合分割\n',EM);
% else
% fprintf('EM=%f, 恐怕不适合分割\n',EM);
% end
% fprintf('\n');
% 依次替换灰度值，自定义目标为黑色、背景为白色
[width,height]=size(frest);
BWimg = ones(width,height);
T = 255*thresh;
for i=1:width
for j=1:height
if(frest(i,j)>T)
BWimg(i,j)= 0;
else
BWimg(i,j)= 255;
end
end
end
% BWimg = getNewOtsu(f); % 改用一种新的Otsu算法得到BWimg
%% 分割、提取ROI，统计标注连通域
% 1.创建平面结构化元素，闭运算，开运算
se = strel('disk',5,4); %创建一个指定半径5的平面圆盘形的结构元素,
% N缺省值为4，圆盘形结构元素由一组N(或N+2)个周期线结构元素来近似
% 创建一个指定领域的平面结构化元素
% 用于膨胀腐蚀及开闭运算等操作的结构元素对象
BWimg = imclose(BWimg,se);%消弭狭窄的间断，消除小的孔洞
BWimg = imopen(BWimg,se);%平滑物体轮廓，断开狭窄的间断和消除细小的突出物
%
% figure;
% imshow(BWimg);title('形态学操作后的图像');
% set(gcf,'color','w');set(gca,'linewidth',1,'fontsize',16);
% 2.统计标注连通域
%使用外接矩形框选连通域，并使用形心确定连通域位置
[mark_image,num] = bwlabel(BWimg,4); % 找二值图像中的连通区域
status = regionprops(mark_image,'BoundingBox'); % 最小外接矩形
centroid = regionprops(mark_image,'Centroid');% 每个区域的质心（重心）
figure;
imshow(mark_image);title('Otsu方法标记后的图像');
set(gcf,'color','w');set(gca,'linewidth',1,'fontsize',16);
for i=1:num
rectangle('position',status(i).BoundingBox,'edgecolor','r','LineWidth',2);
%参数说明：position绘制的为二维图像（他是通过对角的两点确定矩形框）
%edgecolor 指边缘图像，r表示变换为红色。
%facecolor 指内部填充颜色。
text(centroid(i,1).Centroid(1,1)-15,centroid(i,1).Centroid(1,2)-15, num2str(i),'Color', 'r','FontSize',15)
%这个是为绘制出来的矩形框图标记数字
%S(i) = status(i).BoundingBox(2,1)*status(i).BoundingBox(2,1);
end
AREA = zeros(1,num);% 选择外接矩形最大的一个
for ii = 1:num
AREA(ii) = status(ii).BoundingBox(3)*status(ii).BoundingBox(4);
end
[maxAREA,bestNum] = max(AREA);
% 选择观察目标
% 注意:regionprops计算周长和面积，计算的是白色区域的周长和面积！！！
% 黑白颜色颠倒后会出错！！！
copy_mark_image = mark_image;
number1 = bestNum; % 必要适合选择保留最大的区域！！！！！！！！！！
image_part1 = (copy_mark_image == number1); %进行区域的选择，例如只保留3
% 求面积，对二值图像而言就是求目标像素点数
round_area = regionprops(image_part1,'Area');
% fprintf('round_area面积/目标像素点 = %f\n', round_area.Area);
% 求周长
girth = regionprops(image_part1,'Perimeter');
% fprintf('s.Perimeter周长 = %f\n', girth.Perimeter);
%求离心率 oval.Eccentricity
oval = regionprops(image_part1,'Eccentricity');%离心率 0 < e < 1之间，e越小，越像圆。
% fprintf('oval.Eccentricity离心率 = %f\n', oval.Eccentricity);
% 求不变矩
image_part1 = imresize(image_part1,[256 256]); % 拉伸为正方形！！！！
image_part1 = logical(not(image_part1));
[Z, A, Phi] = Zernikmoment(image_part1,4,2);
% fprintf('不变矩 = %f\n', A);
scale = []; % 形状特征！！！
scale = [scale,round_area.Area,girth.Perimeter,oval.Eccentricity,A];
end
以及用来计算正方形图像ROI（需要提前用imresize修改ROI的尺寸）的不变矩的
function [Z, A, Phi] = Zernikmoment(p,n,m)
% Function to find the Zernike moments for an N x N binary ROI，
% 计算泽尔尼克矩（必须是正方形的ROI）
% p = input image N x N matrix (N should be an even number)
% 输入图像N x N矩阵(N应为偶数)
% n = The order of Zernike moment (scalar)
% Zernike矩(标量)的阶数
% m = The repetition number of Zernike moment (scalar)
% Zernike矩重复次数(标量)
% and
% Z = Complex Zernike moment Zernike矩
% A = Amplitude of the moment 振幅(我需要的)
% Phi = phase (angle) of the mement (in degrees) 相位（°）
%
% Example:
% 1- calculate the Zernike moment (n,m) for an oval shape,
% 2- rotate the oval shape around its centeroid,
% 3- calculate the Zernike moment (n,m) again,
% 4- the amplitude of the moment (A) should be the same for both images
% 5- the phase (Phi) should be equal to the angle of rotation
N = size(p,1);
x = 1:N; y = x;
[X,Y] = meshgrid(x,y);
R = sqrt((2.*X-N-1).^2+(2.*Y-N-1).^2)/N;
Theta = atan2((N-1-2.*Y+2),(2.*X-N+1-2));
R = (R<=1).*R;
Rad = radialpoly(R,n,m); % get the radial polynomial
Product = p(x,y).*Rad.*exp(-1i*m*Theta);
Z = sum(Product(:)); % calculate the moments
cnt = nnz(R)+1; % count the number of pixels inside the unit circle
Z = (n+1)*Z/cnt; % normalize the amplitude of moments
A = abs(Z); % calculate the amplitude of the moment
Phi = angle(Z)*180/pi; % calculate the phase of the mement (in degrees)
我的结果如下，有效提取出鸭子的影子呵呵，其实我想要的是鸭子本鸭~我的算法对这样的图只能这样了呵呵
我的方法是用来探测其他东西的，只是拿着鸭子的图片来说明一下哦，如果是想要去除背景实现抠图功能，建议试一下这个免费网站
结果很不错的哦
4.参考文章
借鉴了很多前人经验，列出如下
详解-OTUS(大津法-最大类间方差)原理
基于最大类间、类内方差比法图像分割（有两个小错误，我在程序里面已经修改过了）
zernike矩
基于Matlab的形状识别与计算图形周长，面积，圆周率
这类方法对图像分类很有用但已经很老了，大家可以多思考提出一些新的方法
6129
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
