数字图像处理实验报告实验四：图像分割

阅读量：

数字图像处理实验报告

实验四：图像分割

一、实验目的

掌握图像分割的核心知识。
熟悉阈值法、基于K均值的聚类算法及其应用、边缘检测技术和区域划分与合并的方法。

二、实验环境

Matlab R2022a

三、实验内容

3.1实验题目

基于以下几种典型的阈值选择方法进行图像分割：迭代法与OTSU算法等；依据其理论基础自行编写相关代码；并允许选用任意合适的图像作为分割对象。
基于K-means聚类算法展开图像分割工作；依据其核心机制独立完成相关代码编写；并允许选用任意合适的图像作为分割对象。
对比分析不同边缘检测器的性能特点：包括Roberts算子、Sobel算子以及拉普拉斯高斯算子；并借助系统提供的相应工具函数辅助完成检测过程。
通过合理的选择与优化设计来实现肺部区域的精确分离；其中将肺部区域标记为白色区域并填充相应颜色信息；同时将背景区域标记为黑色区域以明确区分各部分。

3.2相关知识

基于阈值的图像分割方法是一种将图片分解为多个区域并筛选个人关注部分的技术。该方法通过分析目标与背景在亮度上的差异性来实现分割效果。具体而言：①直方图阈值法通过分析亮度分布选择合适的分割点；②半阈值选择法则仅对背景进行二元化处理；③迭代阈值法则通过逐步优化找到最佳分割点；④Otsu算法利用最大类间方差原理自动确定最佳分割点；⑤最大熵法则选取能最大化信息熵的分割点作为标准。

（4）不断改变t，进行以上步骤，不断判断，得到最大类间方差的最优阈值。
OTSU算法被认为是图像分割中阈值选取的最佳算法，计算简单，不受图像亮度和对比度的影响。因此，使类间方差最大的分割意味着错分概率最小。
4、k-means算法
采用K-means算法实现图像分割时，流程如下：
输入：样本集D，簇的数目k，最大迭代次数N；
输出：簇划分（k个簇，使平方误差最小）；
算法步骤：
（1）为每个聚类选择一个初始聚类中心；
（2）将样本集按照最小距离原则分配到最邻近聚类；
（3）使用每个聚类的样本均值更新聚类中心；
（4）重复步骤（2）、（3），直到聚类中心不再发生变化；
（5）输出最终的聚类中心和k个簇划分；
5、算子原理及特点
（1）Roberts算子：利用局部差分算子寻找边缘，边缘定位精度较高，但容易丢失一部分边缘，同时由于图像没经过平滑处理，因此不具备能抑制噪声能力。对陡峭边缘且含噪声少的图像效果较好。
首先使用以下模板分别对图像进行卷积操作：

接着根据公式对图像的所有点进行处理：

（2）Sobel 算子：先执行加权滤波处理以减少噪声影响，在随后的差分计算中提取边缘信息。由于各区域权重设置不同，在抑制噪声方面存在一定效果但无法完全消除虚假边缘的影响。然而，在实际应用中发现这些算子虽然定位精度较高但在边缘检测结果上存在宽度超过单一像素的现象

（3）LOG 算子：首先通过高斯函数进行滤波处理以实现图像平滑效果；然后采用拉普拉斯算子定位边缘特征；该方法弥补了拉普拉斯算子在抗噪声能力方面的不足；但同时也导致过于平滑的边缘信息丢失；因此尖锐突兀的边缘特征无法被检测出来

利用多种算法对图像进行运算可以获得边缘信息。
（1）主程序采用基于二元阈值的方法完成图像分割任务；接着运用罗伯茨、索贝尔以及拉普拉斯高斯等算子来进行边缘检测；并完成了肺部区域划分的过程。

复制代码

    clear all;close all;clc;
    %% 1、采用阈值法实现图像分割，分别采用两种阈值选取方法实现。
    I=imread('bb.jpg');
    I=rgb2gray(I);
    threshold1=myiteration(I);
    threshold2=myotsu(I);
    J1=im2bw(I,threshold1);
    J2=im2bw(I,threshold2);
    figure,subtitle('基于阈值的图像分割');                %Figure 1
    subplot(131),imshow(I),axis on,title('原图') ; 
    subplot(132),imshow(J1),axis on,title('迭代阈值法');
    subplot(133),imshow(J2),axis on,title('otsu阈值法') ; 
    %% 2、分别用Roberts,Sobel和拉普拉斯高斯算子对图像进行边缘检测，比较三种算子处理的不同之处。
    E1=mysharpen(I,'Roberts');
    E2=mysharpen(I,'Sobel');
    E3=mysharpen(I,'Log');
    figure,subtitle('边缘检测');                          %Figure 2
    subplot(221), imshow(I), title('原图');
    subplot(222), imshow(E1), title('Roberts算子');
    subplot(223), imshow(E2), title('Sobel算子');
    subplot(224), imshow(E3), title('Log算子');
    %% 3、实现肺的分割，结果包括两部分：肺（白色显示）和背景（黑色显示）。 
    Lung=imread('lung.png');
    threshold=myotsu(Lung);
    S1=im2bw(Lung,threshold);%二值化，将肺和胸腔等其他部分分隔开
    S2=~S1;%反色，使肺部变为白色
    X8=mycut(S2);%去掉胸腔外的背景和气管
    figure,subtitle('肺的分割');                          %Figure 3
    subplot(221),imshow(Lung),axis on,title('原图') ; 
    subplot(222),imshow(S1),axis on,title('二值化');
    subplot(223),imshow(S2),axis on,title('反色') ;
    subplot(224),imshow(X8),axis on,title('分割后');

（2）函数myiteration实现了迭代阈值法。

复制代码

    % 函数myiteration:迭代阈值法
    % 输入参数:I:原图像
    % 输出参数:阈值threshold
    % 使用函数:im2double(I):变为双精度
    %         min():最小元素
    %         max():最大元素
    %         find(A):符合条件A的元素
    %         mean():求数组平均值
    %         abs():求绝对值
    function threshold=myiteration(I)
    I=im2double(I);     %变为双精度，即0-1
    T0=0.01;            %参数T0
    T1=(min(I(:))+max(I(:)))/2;%选取图像中最大像素和最小像素的中阈值作为T1
    r1=find(I>T1);%灰度值大于T1的
    r2=find(I<=T1);%其他
    T2=(mean(I(r1))+mean(I(r2)))/2;%新的阈值T2    %mean函数是一个求数组平均值的函数
    while abs(T2-T1)<T0
    T1=T2;
    r1=find(I>T1);%灰度值大于T1的
    r2=find(I<=T1);%其他
    T2=(mean(I(r1))+mean(I(r2)))/2;%新的阈值T2
    end
    threshold=T2;

（3）函数myotsu实现了最大类间差法。

复制代码

    % 函数myotsu:最大类间差法
    % 输入参数:I:原图像
    % 输出参数:阈值threshold
    % 使用函数:size(I):求矩阵大小
    %         rgb2gray():rgb转化为灰度图像
    %         im2double():变为双精度，即0-1
    function threshold=myotsu(I)
    [M,N,i]=size(I);     %得到图像行列值
    if i>1
    I=rgb2gray(I);%若不是灰度图像则先处理为灰度图像
    end
    I=im2double(I);  %变为双精度，即0-1
    number_all=M*N;  %总像素值
    hui_all=0;       %预设图像总灰度值为0
    ICV_t=0;        %预设最大方差为0
    for i=1:M
    for j=1:N
        hui_all=hui_all+I(i,j);%得到图像总灰度值
    end
    end
    all_ave=hui_all*255/number_all;   %图像灰度值的总平均值
     
    %t为某个阈值，把原图像分为A部分（每个像素值>=t）与B部分（每个像素值<t）
    for t=0:255                       %不断试探最优t值
    hui_A=0;                      %不断重置A部分总灰度值
    hui_B=0;                      %不断重置B部分总灰度值
    number_A=0;                   %不断重置A部分总像素
    number_B=0;                   %不断重置B部分总像素
    for i=1:M                     %遍历原图像每个像素的灰度值
        for j=1:N
            if (I(i,j)*255>=t)    %分割出灰度值>=t的像素
                number_A=number_A+1;  %得到A部分总像素
                hui_A=hui_A+I(i,j);   %得到A部分总灰度值
            elseif (I(i,j)*255<t) %分割出灰度值<t的像素
                number_B=number_B+1;  %得到B部分总像素
                hui_B=hui_B+I(i,j);   %得到B部分总灰度值
            end
        end
    end
    PA=number_A/number_all;            %得到A部分像素总数与图像总像素的比列
    PB=number_B/number_all;            %得到B部分像素总数与图像总像素的比列
    A_ave=hui_A*255/number_A;          %得到A部分总灰度值与A部分总像素的比例
    B_ave=hui_B*255/number_B;          %得到B部分总灰度值与B部分总像素的比例
    ICV=PA*((A_ave-all_ave)^2)+PB*((B_ave-all_ave)^2);  %Otsu算法
    if (ICV>ICV_t)                     %不断判断，得到最大方差
        ICV_t=ICV;
        k=t;                           %得到最大方差的最优阈值
    end
    end
    threshold=k/255;                           %输出阈值

（4）函数mycut实现了根据连通性去掉图片胸腔外的背景和气管。

复制代码

    % 函数mycut:根据连通性去掉图片胸腔外的背景和气管
    % 输入参数:I:原图像
    % 输出参数:OUT：处理后的图片数据
    % 使用函数:size(I):求矩阵大小
    %         zeros():建全零矩阵
    function OUT=mycut(I)
    [x,y]=size(I);
    [X8,N] = bwlabel(I,8);%8连通
    % s=regionprops(X8,'Area');
    % X8=ismember(X8,find([s.Area]>=7000&[s.Area]<=40000));
    %% 计算各个连通区域的面积
    Area=zeros(N,1);
    for i=1:x
    for j=1:y
        for k=1:N
            if X8(i,j)==k
                Area(k,1)=Area(k,1)+1;
            end
        end
    end
    end
    %% 将白色部分连通区域的背景部分、气管部分置黑
    for k=1:N
    	if Area(k,1)>=7000 && Area(k,1)<=40000
        for i=1:x
            for j=1:y
                if X8(i,j)==k-1
                X8(i,j)=0;
                end
            end
        end
    	end
    end
    OUT=logical(X8);
    
    （5）函数mysharpen:根据选择的算子对图像进行锐化处理
    % 函数mysharpen:根据选择的算子对图像进行锐化处理
    % 输入参数:I:原图像
    %          modelx:水平边缘算子
    %          modely:垂直边缘算子
    % 输出参数:锐化处理后的图像OUT
    % 使用函数:size(x):求矩阵大小
    %         sqrt(x):开根号
    %         strcmp(choice,''):比较字符串是否相等
    %         myconvolve(I，model):卷积操作
    
    function OUT=mysharpen(I,choice)
    [x,y]=size(I);
    %% 采用Roberts算子进行边缘检测
    if strcmp(choice,'Roberts')
    modelx=[0,0,0;0,-1,0;0,0,1];%水平边缘算子
    modely=[0,0,0;0,0,-1;0,1,0];%垂直边缘算子
    b=myconvolve(I,modelx);
    c=myconvolve(I,modely);
    for i=1:x
        for j=1:y
             b(i,j)=sqrt(b(i,j)^2+c(i,j)^2); 
        end
    end
    OUT=uint8(b);
    end
    %% 采用Sobel算子进行边缘检测
    if strcmp(choice,'Sobel')
    modelx=[-1,-2,-1;0,0,0;1,2,1];%水平边缘算子
    modely=[1,0,-1;2,0,-2;1,0,-1];%垂直边缘算子
    b=myconvolve(I,modelx);
    c=myconvolve(I,modely);
    for i=1:x
        for j=1:y
             b(i,j)=sqrt(b(i,j)^2+c(i,j)^2); 
        end
    end
    OUT=uint8(b);
    end
    %% 采用Prewitt算子进行边缘检测
    if strcmp(choice,'Prewitt')
    modelx=[-1,-1,-1;0,0,0;1,1,1];%水平边缘算子
    modely=[-1,0,1;-1,0,1;-1,0,1];%垂直边缘算子
    b=myconvolve(I,modelx);
    c=myconvolve(I,modely);
    for i=1:x
        for j=1:y
             b(i,j)=sqrt(b(i,j)^2+c(i,j)^2); 
        end
    end
    OUT=uint8(b);
    end
    %% 采用Log算子进行边缘检测
    if strcmp(choice,'Log')
    model=[-2,-4,-4,-4,-2;
           -4,0,8,0,-4;
           -4,8,24,8,-4;
           -4,0,8,0,-4;
           -2,-4,-4,-4,-2];%模板
    b=myconvolve(I,model);
    OUT=uint8(b);
    end

（6）函数myconvolve实现了根据模板对图像进行卷积操作。

复制代码

    % 函数myconvolve:根据模板对图像进行卷积操作
    % 输入参数:I:原图像
    %          modelx:模板
    % 输出参数:卷积操作后的图像数据OUT
    % 使用函数:size(x):求矩阵大小
    %         double(x):增加精度
    %         zeros():为矩阵分配空间
    function OUT=myconvolve(image,model)
    HW=size(image);             %获取原图像的大小
    image=double(image);
    [x,y]=size(model);          %模板大小
    OUT=zeros(HW);              %为输出分配空间
    %% 加边沿
    wall=(x-1)/2;
    nHW=[HW(1)+2*wall,HW(2)+2*wall];%加上墙（补充的边沿）之后的数据区大小
    data=zeros(nHW);%为数据区分配空间
    data(wall+1:wall+HW(1),wall+1:wall+HW(2))=image(:,:);%为数据区赋值
    for i = 1:wall
    data(wall+1:wall+HW(1),i)=data(wall+1:wall+HW(1),wall+1);%左
    data(wall+1:wall+HW(1),i+(nHW(2)-wall))=data(wall+1:wall+HW(1),wall+HW(2));%右
    end
    for i = 1:wall
    data(i,:)=data(wall+1,:);%上
    data(i+(nHW(1)-wall),:)=data(nHW(1)-wall,:);%下
    end
    %% 计算
    for i = wall+1:wall+HW(1)
    	for j = wall+1:wall+HW(2)
        for m = 1:x
            for n = 1:y
                OUT(i-wall,j-wall) = OUT(i-wall,j-wall) + data(i-wall+m-1,j-wall+n-1)*model(m,n);%对整幅图像进行卷积操作
            end
        end
    	end
    end
    
    （7）实现了K-means算法，将其展开实现。先随机选取k个聚类中心，度量每个像素点到每个聚类中心的“距离”并按“距离”分类。对所有像素点的分类完成后，根据每个像素点的类型计算更新其“聚类中心”的灰度级值。
    clear all;clc;
    I=imread("bb.jpg");
    I=im2gray(I);
    k=3;
    Mu=[1,1,1;127,127,127;255,255,255];
    eps=1e-6;
    normMuDiff=[];
    m=size(I,1);
    n=size(I,2);
    I=double(I);
    M=reshape(I,m*n,size(I,3));
    flag=1;
    Cov=repmat(eye(size(I,3)),1,1,k); %k个聚类的协方差
    Label=zeros(m*n,1);
    iter=1;
    figure;
    for iter=1:5
    old_Mu=Mu;
    Num=zeros(1,k);
    ClassedPixs=zeros(m*n,size(I,3),k); %分到k个类中的像素值
    w=zeros(m*n,k); %各像素点对于聚类的权重矩阵
    C=zeros(size(I,3),size(I,3),k);
    for j=1:k
        invCov=inv(Cov(:,:,j));
        C(:,:,j)=invCov/norm(invCov);%使用Mahalanobis距离(即马氏距离)度量
    end
    for i=1:m*n
        dis=zeros(k,1);
        for j=1:k
            dis(j)=(M(i,:)-Mu(j,:))*C(:,:,j)*(M(i,:)-Mu(j,:))'; 
        end
        maxIdx=find(dis==min(dis));
        label=maxIdx(1);
        Label(i)=label;
        Num(label)=Num(label)+1;
        ClassedPixs(Num(label),:,label)=M(i,:);
        w(Num(label),label)=1/(sqrt(min(dis))+1); %避免被0除的情况
    end
    w=w./repmat(sum(w),m*n,1); %权重归一化
    % 更新协方差Cov
    for j=1:k
        PixsThisClass=ClassedPixs(1:Num(j),:,j);
        if Num(j)~=0      
            Cov(:,:,j)=cov(PixsThisClass)+0.01*rand(size(I,3)); %加随机量以防止矩阵奇异
            Mu(j,:)=sum(repmat(w(1:Num(j),j),1,size(I,3)).*PixsThisClass);
        end
    end
    normDiff=norm(old_Mu-Mu);
    if  normDiff<eps
        flag=0;
    end
    normMuDiff=[normMuDiff;normDiff];
    sg=reshape(Label,m,n);
    end
    imshow(mat2gray(sg))

7、实验结果
（1）下图1为迭代和otsu两种阈值分割法的结果。

图1、阈值分割法
（2）下图2为K-means聚类算法实现图像分割的结果。

图表2采用k-means聚类算法进行分类。（4）本研究中所用的下标为(4)的图表展示了利用罗伯茨、索贝尔以及拉普拉斯高斯滤波器对图像进行边缘检测后的结果。

图像3展示了基于R-O算法以及拉普拉斯高斯滤波器的边缘检测效果。（在图4中展示了肺部的处理结果，在主色调区域（对应白色显示）与背景区域（对应黑色显示）之间进行了区分。）

图4、肺的分割

四、实验心得

熟悉了Matlab系统函数的功能细节，并对图像阈值法、K-means聚类方法、边缘提取及区域生长和分裂合并等方法的工作原理进行了深入的学习和理解；
通过实践迭代阈值法和OTSU算法来确定阈值参数，进一步掌握了如何利用连通区域分析来处理医学图像中的细节问题。

全部评论 (0)

还没有任何评论哟~

数字图像处理实验报告实验四：图像分割

数字图像处理实验报告实验四：图像分割一、实验目的 1、理解图像分割的基本概念。 2、掌握阈值法、Kmeans聚类方法、边缘提取及区域生长和分裂合并等方法进行图像分割。

matlab数字图像处理实验报告,数字图像处理实验报告(完整版)..doc

数字图像处理实验报告完整版. 数字图像处理实验一MATLAB数字图像处理初步1．利用imread函数读取一幅图像，假设其名为.tif，存入一个数组中； 2．利用whos命令提取该读入图像flower...

数字图像处理实验报告实验五：形态学图像处理

实验五形态学图像处理一、实验目的 1、了解形态学的基本理论和方法。 2、掌握对图像进行开、闭、膨胀、腐蚀的方法。二、实验环境 MatlabR2022a 三、实验内容 3.1实验题目 1、任选一种结...

医学图像处理matlab实验报告,matlab-数字图像处理实验报告(五份

《matlab数字图像处理实验报告五份》由会员分享，可在线阅读，更多相关《matlab数字图像处理实验报告五份13页珍藏版》请在人人文库网上搜索。 1、数字图像处理实验报告实验一图像的增强一.实验目的...

matlab数字图像处理实验报告

MATLAB数字图像处理仿真系统实验一常用MATLAB图像处理命令一、实验目的 1、熟悉并掌握MATLAB工具的使用； 2、实现图像的读取、显示、代数运算和简单变换。

Matlab数字图像处理实验报告

尘挂菠塘码惹忧良袍啥恢贱迸待梅夸甲陇挡集解拂违希菜瓷趟堆坠冕珍叶眶臣峙堵嘉幼酒俘颧栅贷榜六瘪塑港腰朋治策昏档浦顿承柴享瓷佣尧浇蔬汗推揖杯赘吐鸯获鲜纬坡眩要揍圈件刺战宴完蛮已茬权闸滚磐傣楞砰季涅甄好引保...

数字图像处理实验报告实验一：图像处理基本操作

数字图像处理实验报告实验一：图像处理基本操作一、实验目的 1、熟悉并掌握MATLAB工具的使用； 2、实现图像的读取、显示、存储、平移、镜像、放大、缩小及旋转操作； 3、掌握常用的插值方法，并了解...

数字图像处理 matlab 报告总结,matlab 数字图像处理实验报告(五份)

等原因，成像后的图像与原景物图像相比，会产生比例失调，甚至扭曲。以上图像退化现象称之为几何失真有几何畸变的图像，不但视觉效果不好，而且在对图像进行定量分析时提取的形状、距离、面积等数据也不准确。

数字图像处理实验四

%分别在灰度图像中加入一定量的高斯噪声和椒盐噪声，然后采用3×3的均值滤波器和3×3中值滤波器分别对噪声图像进行处理， %给出两种处理方法的峰值信噪比PSNR；仿效“中值滤波”的方法，对原图像分别进行...

图像处理-图像锐化实验报告

实验报告课程名：图像处理创新实验内容：图像锐化学生姓名专业班级学号指导教师常化文实验日期2022年10月25日一、实验目的掌握matlab环境下，函数（function）文件的编写方...

是否确定退出登录?

数字图像处理实验报告 实验四：图像分割