2020-2021学年——图像图形编程实践实验5_彩色图像的切割

彩色图像的切割

实验目的

1. 理解图像分割的基本概念；
2. 理解彩色图像的基本概念；
3. 掌握进行欧氏距离和绝对值距离进行图像分割的基本方法；

实验设备

PC机、matlab2018b

实验原理

  指定的颜色区域的物体。给定一个感兴趣的有代表性彩色的彩色样点集，可得到我们希望分割的颜色的“平均”估计。令这个平均彩色用RGB向量a来表示。分割的目的是将给定图像中的每个RGB像素分类，即在指定的区域内是否有一种颜色。为了执行这一比较，有一个相似性度量是必要的。最简单的度量之一是欧氏距离。令z表示RGB空间中的任意一点。如果它们之间的距离小于特定的阙值D_0。.则称z与a是相似的。z和a间的欧氏距离由下式给出:

$$D(z,a)=\Vert z-a\Vert =[(z-a{)}^T(z-a){]}^{\frac{1}{2}}=[(z_R-a_R{)}^2+(z_G-a_G{)}^2+(z_B-a_B{)}^2{]}^{\frac{1}{2}}\text{\hspace{1em}(5.1)}$$

  Kmeans算法是最常用的聚类算法，主要思想是:在给定K值和K个初始类簇中心点的情况下，把每个点(亦即数据记录)分到离其最近的类簇中心点所代表的类簇中，所有点分配完毕之后，根据一个类簇内的所有点重新计算该类簇的中心点(取平均值)，然后再迭代的进行分配点和更新类簇中心点的步骤，直至类簇中心点的变化很小，或者达到指定的迭代次数。

问题描述

  编程实现彩色图像的切割。设计一个通用的方法，基于RGB空间，从下列的图像中，分割出指定的目标（如草莓）。
1.特征：RGB向量（3维）；
2.计算特征的距离—定义草莓的空间（球：欧式距离，立方体：绝对值距离），分析这两种不同的距离的差异。
  欧式距离：中心点，半径阈值R
  绝对值距离：中心点，距离阈值R
  中心点定义：[255,0,0] ； 随机截取一个草莓区域，这个区域的RGB向量的均值作为中心点，阈值，可以用方差表示。

计算方法

  选用图片如下图所示：

图5-1

1.欧氏距离
1.1 不同区域取样
  (1)读取实验用图像记为img：使用imread函数将图像读入Matlab；
  (2)将图像转为double类型；
  (3)分别在图像上截取样本草莓红色部分、叶子部分、种子部分如下图所示：

图5-2

图5-3

图5-4

  (4)将截取的图像拆分为R、G、B三层，分别计算出草莓红色部分、叶子部分、种子部分每一层的均值average，按照公式5.1计算计算截取图像中每一点与均值的欧式距离：
再计算样本标准差作为阈值$D_0$如下：

$$D_0=\sqrt{\frac{{\sum }_{i=1}^N(x_i-\mu {)}^2}{N}}\text{\hspace{1em}(5.2)}$$

  其中，N为元素个数，为均值，为每个元素的值。
  (5)创建一个和原图像长宽一样的矩阵记为img_gray，将原图像拆分为R、G、B三层，遍历图像计算每一个点与样本均值计算欧式距离，若距离小于等于$D_0$，将图像img_gray相应位置置1，若距离大于等于$D_0$，将图像img_gray相应位置置为0，切割出来的草莓红色部分、叶子部分、种子部分如下图所示：

图5-5

图5-6

图5-7

(6)将切割出来的草莓红色部分、叶子部分、种子部分相加得到如下图像记为img_grayAll：

图5-8

(7)在彩色图像上做如下处理，将原图拆分为R、G、B三层，遍历R、G、B三层图像，在(i,j)点，若img_grayAll(i,j)为0，则R、G、B三层图像(i,j)点置0，否则不做变化，最后将R、G、B三层合并为一张图像记为img_colour如下图所示：

图5-9

1.1.2 取样整个草莓
  (1)读取实验用图像记为img：使用imread函数将图像读  入Matlab；
  (2)将图像转为double类型；
  (3)在图像上截取包含草莓红色部分、叶子部分、种子部分如下图所示：

图5-10

  余下步骤与1.1中的(4),(5),(6),(7)相同，最后得到图像如下所示：

图5-11

2.曼哈顿距离
  (1)读取实验用图像记为img：使用imread函数将图像读入Matlab；
  (2)将图像转为double类型；
  (3)在图像上截取包含草莓红色部分、叶子部分、种子部分如下图所示：

图5-12

  (4)将截取的图像拆分为R、G、B三层，计算截取部分每一层的均值average，按照如下公式计算计算截取图像中每一点与均值的曼哈顿距离：

$$D=\sum _{k=1}^n\left|x_{1k}-x_{2k}\right|\text{\hspace{1em}(5.3)}$$

  其中，n为数据维度，$x_{1k}$代表第一个元素的第k维的值，$x_{2k}$代表第二个元素的第k维的值。
再计算阈值$D_0$如下所示：

$$D_0=\frac{{\sum }_{i=1}^N|x_i-\mu |}{N}\text{\hspace{1em}(5.4)}$$

  其中，N为元素个数，为均值，为每个元素的值。
  (5)创建一个和原图像长宽一样的矩阵记为img_gray，将原图像拆分为R、G、B三层，遍历图像计算每一个点与样本均值计算曼哈顿距离，若距离小于等于$D_0$，将图像img_gray相应位置置1，若若距离大于等于$D_0$，将图像img_gray相应位置置0，切割出来的草莓如下图所示：

结果分析

  在欧氏距离分割过程中，分别使用草莓各部分相加得到图像不如使用整个草莓分割的效果好。
  欧氏距离与曼哈顿距离相比，欧氏距离分割效果较好，可能是因为欧氏距离计算的是直线距离，在低维数据上使用欧氏距离效果比较好。
  在Kmeans分割过程中，2分类效果较好，保留了较多的草莓红色部分。
  欧氏距离分割方法与Kmeans分割方法相比，欧氏距离分割方法保留了更多的草莓红色部分，Kmeans分割方法保留了更多的绿色叶子。
  欧氏距离与曼哈顿距离分割方法更依赖于样本点的选择，Kmeans分割方法则要求提前知道类别数目。

实验代码

function img_2 = ex5_MyDivision_1(img,sample,D_0_choose,D_0)
% img：输入图像，sample：样本点，D_0_choose：是否选择默认阈值，非0代表使用自定义阈值D_0，D_0：阈值
% img_2:欧氏距离分割的图像
img = double(img);
sample = double(sample);
[m,n] = size(img(:,:,1));
sample_R = sample(:,:,1);
sample_G = sample(:,:,2);
sample_B = sample(:,:,3);
img_R = img(:,:,1);
img_G = img(:,:,2);
img_B = img(:,:,3);
% 均值
average = [mean(mean(sample_R)),mean(mean(sample_G)),mean(mean(sample_B))];
% 默认阈值为方差
if D_0_choose == 0
    [p,q] = size(sample(:,:,1));
    std_column = zeros(p,q);
    for i = 1:p
        for j = 1:q
            std_column(i,j) = (sample_R(i,j) - average(1))^2 + (sample_G(i,j) - average(2))^2 + (sample_B(i,j) - average(3))^2;
        end
    end
    D_0 = sqrt(sum(sum(std_column)) / p / q);
end
img_2 = zeros(m,n);
% 遍历图像计算
for i = 1:m
    for j = 1:n
        D = sqrt((img_R(i,j) - average(1))^2 + (img_G(i,j) - average(2))^2 + (img_B(i,j) - average(3))^2);
        if D <= D_0
            img_2(i,j) = 1;
        else
            img_2(i,j) = 0;
        end
    end
end
function img_2 = ex5_MyDivision_2(img,sample,D_0_choose,D_0)
% img：输入图像，sample：样本点，D_0_choose：是否选择默认阈值，非0代表使用自定义阈值D_0，D_0：阈值
% img_2:曼哈顿距离分割的图像
img = double(img);
sample = double(sample);
[m,n] = size(img(:,:,1));
sample_R = sample(:,:,1);
sample_G = sample(:,:,2);
sample_B = sample(:,:,3);
img_R = img(:,:,1);
img_G = img(:,:,2);
img_B = img(:,:,3);
% 均值
average = [mean(mean(sample_R)),mean(mean(sample_G)),mean(mean(sample_B))];
% 计算默认阈值
if D_0_choose == 0
    [p,q] = size(sample(:,:,1));
    ManDist = zeros(p,q);
    for i = 1:p
        for j = 1:q
            ManDist(i,j) = abs(sample_R(i,j) - average(1)) + abs(sample_G(i,j) - average(2)) + abs(sample_B(i,j) - average(3));
        end
    end
    D_0 = sum(sum(ManDist)) / p / q;
end
 
img_2 = zeros(m,n);
% 遍历图像计算
for i = 1:m
    for j = 1:n
        D = (abs(img_R(i,j) - average(1)) + abs(img_G(i,j) - average(2)) + abs(img_B(i,j) - average(3)));
        if D <= D_0
            img_2(i,j) = 1;
        else
            img_2(i,j) = 0;
        end
    end
end
function [mu,mask] = ex5_MyKmeans(ima,k)
% ima:输入的灰度图像 K:分类数
% mu:均值类向量  mask:分类后的图像
ima=double(ima);
copy=ima;
ima=ima(:);
mi=min(ima);%找到最小值
ima=ima-mi+1;
s=length(ima);%有多少灰度级
%计算图像灰度直方图
m=max(ima)+1;%图像最大灰度值
h=zeros(1,m);
hc=zeros(1,m);
for i=1:s
    if  (ima(i)>0)
        h(ima(i))=h(ima(i))+1;%灰度值i累加
    end
    
end
ind =find(h);
h1=length(ind);
%初始化质心
mu=(1:k)*m/(k+1);
%start process
while(true)
    oldmu=mu;
    %现有的分类·
    for i=1:h1
        c=abs(ind(i)-mu);
        cc=find(c==min(c));
        hc(ind(i))=cc(1);
    end
    %重新计算均值
    for i=1:k
        a=find(hc==i);
        mu(i)=sum(a.*h(a))/sum(h(a));
    end
    if(mu==oldmu) 
        break;
    end
end
    %计算生成分类后的图像
s=size(copy);
mask=zeros(s);
for i=1:s(1)
    for j=1:s(2)
        c=abs(copy(i,j)-mu);
        a=find(c==min(c));
        mask(i,j)=a(1);
    end
end
mu=mu+mi-1;