数字图像处理-图像分割
1.线检测模板检测
使用线检测模板检测右图中的直线和斜线。使用edge函数，分别用Sobel、LoG和Canny边缘检测器检测右图中的直线和斜线，看哪种方法效果更好。
代码：
%使用线检测模板检测右图中的直线和斜线。使用edge函数，分别用Sobel、LoG和Canny边缘检测器检测右图中的直线和斜线，看哪种方法效果更好。
clear;
clc;
%获取图像
f = imread("Fig1034(a)(marion_airport).tif");
%使用edge函数扫描直线和斜线
f_S = edge(f,"Sobel");
f_L = edge(f,"LoG");
f_C = edge(f,"Canny");
%展示原图像和三种算子的效果
subplot(2,2,1);imshow(f);title("原图");
subplot(2,2,2);imshow(f_S);title("使用Sobel算子");
subplot(2,2,3);imshow(f_L);title("使用LoG算子");
subplot(2,2,4);imshow(f_C);title("使用Cannay算子");
效果图：
线检测通过不同的模板检测不同方向的线，通过使用edge函数可以方便的使用不同的算子来对线进行检测，Sobel、LoG和Canny都有各自不同的优点，Sobel算子通过导数来寻找边缘，LoG是先使用高斯拉普拉斯滤波进行平滑后再通过检验交叉零点的方式寻找边缘，Canny通过双阈值检测强弱边缘能够有效对抗噪声。
### 2.Hough变换
使用edge函数，分别用Sobel、LoG和Canny边缘检测器检测右图中的边缘，选出检测效果最好的图像，然后用Hough变换对其做线检测和连接。
代码：
%使用edge函数，分别用Sobel、LoG和Canny边缘检测器检测右图中的边缘，选出检测效果最好的图像，然后用Hough变换对其做线检测和连接。
f = imread("Fig1027(a)(van_original).tif");
f = tofloat(f);
%使用edge函数扫描直线和斜线
f_S = edge(f,"Sobel",0.05);
f_L = edge(f,"LoG",0.003,2.25);
f_C = edge(f,"Canny",[0.04 0.1],1.5);
%展示原图像和三种算子的效果
subplot(2,2,1);imshow(f);title("原图");
subplot(2,2,2);imshow(f_S);title("使用Sobel算子");
subplot(2,2,3);imshow(f_L);title("使用LoG算子");
subplot(2,2,4);imshow(f_C);title("使用Cannay算子");
%由于Log算子的表现最好，使用霍夫变化对其做线检测和连接
H = hough(f_L);
%显示默认的霍夫变换图形
figure;imshow(H,[]);
[H, theta, rho] = hough(f_L);
%显示带有刻度轴的霍夫变换
figure;imshow(H, [], 'XData', theta, 'YData', rho)
axis on, axis normal
xlabel('\theta'), ylabel('\rho');
[H, theta, rho] = hough(f_L, 'Theta',-90:0.2:90-0.2);
imshow(H, [], 'XData', theta, 'YData', rho);
axis on, axis normal;
xlabel('\theta'), ylabel('\rho');
%寻找高计数的累加器单元
peaks = houghpeaks(H, 15);
hold on
plot(theta(peaks(:,2)), rho(peaks(:,1)), 'linestyle', 'none', 'marker', 's', 'color', 'w');
%确定与峰值相关联的有意义线段以及他们的端点
lines = houghlines(f_L, theta, rho, peaks);
figure, imshow(f_L), hold on
for k = 1: length(lines)
xy = [lines(k).point1; lines(k).point2];
plot(xy(:,1), xy(:,2), 'LineWidth', 4, 'Color', [0.8 0.8 0.8]);
end
效果图：
3.阈值处理
使用阈值处理方法分割下图：
代码：
%使用阈值处理方法分割右图
f = imread("Fig1001(d)(noisy_region).tif");
%由于图像存在椒盐噪声，使用中值滤波实现对图像的降噪之后再分割图像
gm = medfilt2(f,"symmetric");
%使用全局阈值和最佳全局阈值去分割图像
T = mean2(gm); %初始值，全局平均
done = false;
%开始迭代，直到上一次的阈值相差不超过0，5
while ~done
count += 1;
g = gm>T;
Tnext = 0.5*(mean(gm(g))+mean(gm(~g)));
done = abs(T-Tnext) < 0.5;
T = Tnext;
end
g = imbinarize(gm,T/255);
%使用最佳全局阈值分割
h = imhist(gm);
[T1, SM] = graythresh(h); %T1在[0, 1]之间
g2 = imbinarize(gm,T1);
subplot(2,2,1);imshow(f);title("原图")
subplot(2,2,2);imshow(gm);title("降噪原图");
subplot(2,2,3);imshow(medfilt2(g,"symmetric"));title("全局赋值分割图");
subplot(2,2,4);imshow(medfilt2(g2,"symmetric"));title("最佳全局阈值分割");
效果图：
基本全局阈值是通过均值的不断迭代使得相差不超过预期值停止，而最佳阈值通过最大类间方差来对图像实现分割，在这幅图像中存在椒盐噪声，使用中值滤波之后再进行分割能够得到得效果更好。
4. 阈值&区域生长
分别使用(1)基本全局阈值(2)Otsu阈值(3)基于边缘改进的全局阈值和(4)区域生长方法分割出下图山的区域，看哪种方法效果更好
基本全局阈值是通过均值的不断迭代使得相差不超过预期值停止，而最佳阈值通过最大类间方差来对图像实现分割，边缘改进得全局阈值是通过对检测出得边缘进行分割而不必考虑非边缘元素，而区域生长是通过将相似像素归在一起直到邻边像素不再一样。而局部可变阈值分割是通过局部得均值和方差来确定局部得阈值。为了能过确定区域生长得种子点，我通过bwlabel函数确定每一个连通区域，将除黑色得连通域外其余点作为种子点。
代码：
clear;
clc;
%分别使用(1)基本全局阈值(2)Otsu阈值(3)基于边缘改进的全局阈值和(4)区域生长方法分割出右图山的区域，看哪种方法效果更好。
%读文件
f = imread("Fig1045(a)(iceberg).tif");
%使用基本全局阈值
T = mean2(f); %全局平均，得到均值
done = false;
%开始迭代，直到两次平均不超过0，5
count = 0;
while ~done
count = count +1;
g = f>T;
Tnext = 0.5 * mean(mean(f(g))+mean(f(~g)));
done = abs(T-Tnext) < 0.5;
T = Tnext;
end
g1 = imbinarize(f,T/255);
fprintf("一共迭代了 %d 次\n",count);
%使用Otsu阈值
h = imhist(f);
[T1,SM] = graythresh(h);
g2 = imbinarize(f,T1);
%基于边缘改进的全局阈值
f = tofloat(f);
%step1:利用sobel算子计算边缘
sx = fspecial("sobel");
sy = sx';
gx = imfilter(f,sx,'replicate');
gy = imfilter(f,sy,'replicate');
grad = sqrt(gx.*gx + gy.*gy); %梯度
%将边缘图像单位化
grad = grad / max(grad(:));
%step2 指定阈值T
h = imhist(grad);
Q = percentile2i(h,0.999); %从图像中找到富含所有像素0.99的阈值
%step3:阈值处理
markerImg = grad > Q;
%step4:仅使用f中对应markerImage中1值像素，计算其直方图
fp = f.* markerImg;
hp = imhist(fp,256);
%step5:使用步骤4的直方图，采用Otsu方法全局分割
hp(1) = 0; %删除背景黑色像素
T = graythresh(hp);
%T*(numel(hp)-1);
g3 = imbinarize(f,T);
%局部分割
g4 = localthresh(f,ones(3),30,1.6,'global'); %使用局部阈值分割
%区域生长法来实现对山的分割
f = tofloat(f);
%找到种子点这幅图像的
[L,n] = bwlabel(f);
S = zeros(size(f));
for i=1:n-1
[r,c] = find(L==i);
S(r(1),c(1)) = 1;
end
[g5,NR,ST,TI] = regiongrow(f,S,0.4);
figure();
imshow(g5);
subplot(2,3,1);imshow(f);title("原图");
subplot(2,3,2);imshow(g1);title("基本全局阈值")
subplot(2,3,3);imshow(g2);title("最佳全局阈值")
subplot(2,3,4);imshow(g3);title("基于边缘改进的")
subplot(2,3,5);imshow(g4);title("基于局部分割的")
subplot(2,3,6);imshow(g5);title("基于区域生长的")
效果图：
5.Otsu阈值和基于边缘改进
分别使用基本全局阈值、Otsu阈值和基于边缘改进的全局阈值方法分割下图，看哪种方法效果相对更好
基本全局阈值是通过均值的不断迭代使得相差不超过预期值停止，而最佳阈值通过最大类间方差来对图像实现分割，边缘改进得全局阈值是通过对检测出得边缘进行分割而不必考虑非边缘元素。基于边缘改进由于使用的是sobel算子，对于噪声的敏感性太强了，在对边缘分割后的图像进行5*5的开操作之后得到的图像显得更加的平滑。
代码：
clear;
clc;
%分别使用基本全局阈值、Otsu阈值和基于边缘改进的全局阈值方法分割右图，看哪种方法效果相对更好
f = imread("FigP1036(blobs).tif");
%f = tofloat(f);
f = medfilt2(f,"symmetric");
%使用基本全局阈值
done = false;
count = 0;
T = mean2(f);
while ~done
count = count+1;
g = f > T;
Tnext = 0.5*mean(mean(f(g))+mean(f(~g)));
done = abs(Tnext-T) < 0.5;
T = Tnext;
end
g1 = imbinarize(f,T/255);
fprintf("迭代了%i次",count);
figure();
%使用otsu阈值进行
h = imhist(f);
[T1,SM] = graythresh(h);
g2 = imbinarize(f,T1);
%使用基于边缘改进的全局阈值方法
f = tofloat(f);
%step1:使用soble算子计算边缘
sx = fspecial("sobel");
sy = sx';
gx = imfilter(f,sx,"replicate");
gy = imfilter(f,sy,"replicate");
grad = sqrt(gx.*gx + gy.*gy)
%将边缘图像单位化
grad = grad / max(grad(:));
%指定阈值T
h = imhist(grad);
Q = percentile2i(h,0.999);
%step3:阈值处理
makeImage = grad > Q;
%step4：仅使用f中非黑元素计算直方图
fp = f.* makeImage;
hp = imhist(fp,256);
%step5:使用步骤4的直方图，采用Otsu方法全局分割
hp(1) = 0; %删除背景黑色像素
T = graythresh(hp);
%T*(numel(hp)-1);
g3 = imbinarize(f,T);
subplot(2,2,1);imshow(f);title("均值降噪后图形");
subplot(2,2,2);imshow(g1);title("基本全局阈值图形");
subplot(2,2,3);imshow(g2);title("Otsu阈值图形");
b = ones(5,5);
subplot(2,2,4);imshow(imopen(g3,b));title("基于边缘改进的阈值");
效果图：
6.统计数量
考虑右侧的图像，它显示了较大圆圈区域包围的小圆圈区域。请将这些圆圈分隔出来，并统计大圆圈和小圆圈的数量
本题主要是考察形态学和分割还有连通的综合应用，对于白背景黑粒子图像，使用结构元闭操作，则大于结构元的粒子会显示出来，小的则删除了。通过这样的方式我们可以得到的是大圆圈的图形，再通过将所有黑像素点置为1，其余像素置为0可通过bwlabel获得多个连通分量从而判定黑色大圆圈的数量。针对得到大圆圈图像使用开操作可删除大圆圈保留小圆圈的内部区域，将这部分内部区域再原图中保留并删除其余部分，最后再使用闭运算得到所有的小圆圈的图形，之后按照统计大圆圈的数量的方式得到小圆圈的数量。
代码：
clear;
clc;
%考虑右侧的图像，它显示了较大圆圈区域包围的小圆圈区域。请将这些圆圈分隔出来，并统计大圆圈和小圆圈的数量
f = imread("FigP0934(blobs_in_circular_arrangement).tif");
figure('name','灰度级形态学'),subplot(1,3,1),imshow(f),title('原图像')
b=strel('disk',30);
g1=imclose(f,b);%对于白背景黑粒子图像，使用结构元闭操作，则大于结构元的粒子会显示出来，小的则删除了。
subplot(1,3,2),imshow(g1,[]),title('半径30结构元闭操作a')
r1 = zeros(size(f)); %得到一个二值图像的初始值
I = g1 == min(g1(:)); %找到黑色元素所在的点
r1(I) = 1; %把所有的黑色素的点置为1
[L1 ,n1] = bwlabel(r1);
fprintf("大圆圈的数量为 %d\n",n1);
b=strel('disk',50);%对于a来说白色为前景，
g2=imopen(g1,b);%使用开操作，删除了小于50的白例子，得到中心的一大块白色
f1 = f;
I = g2 == min(g2(:));
f1(I) = 255;
b=strel('disk',15);
g3 = imclose(f1,b);
r2 = zeros(size(f)); %得到一个二值图像的初始值
I = find(g3 == min(g1(:))); %找到黑色元素所在的点
r2(I) = 1; %把所有的黑色素的点置为1
[L2 ,n2] = bwlabel(r2);
fprintf("小圆圈的数量为 %d",n2);
%g2 = imbinarize(g2,254);
subplot(1,3,3),imshow(g3,[]),title('对a使用半径50结构元开操作')
效果图：
7.请分割出右图图形部分。
代码：
%请分割出右图图形部分。
clear;
clc;
%读入均值为0、标准差为50个灰度的高斯噪声图像
fn = imread('Fig1046(a)(septagon_noisy_shaded).tif');
%直接使用Otsu方法分割
Tn = graythresh(fn);
gn = imbinarize(fn, Tn);
%平滑噪声图像
w = fspecial('average', 5);
fa = imgaussfilt(fn,2);
%fa = imfilter(fn, w, 'replicate');
%使用Otsu方法分割平滑后的图像
Ta = graythresh(fa);
ga = imbinarize(fa, Ta);
%使用5*5局部阈值分割
g = localthresh(fa, ones(5), 10, 1.1, 'global'); %使用局部阈值分割
SIG = stdfilt(fa, ones(5));
%显示效果
subplot(3,3,1); imshow(fn); title('噪声图像');
subplot(3,3,2); imhist(fn); title('噪声图像直方图');
subplot(3,3,3); imshow(gn); title('Otsu直接分割');
subplot(3,3,4); imshow(fa); title('高斯平滑');
subplot(3,3,5); imhist(fa); title('平滑图像直方图');
subplot(3,3,6); imshow(imopen(ga,ones(3,3))); title('Otsu分割平滑图像');
subplot(3,3,7); imshow(SIG, []); title('局部标准差图像');
subplot(3,3,9); imshow(imopen(g,ones(3,3))); title('局部阈值分割');
效果图：
5万+
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
