标记控制的分水岭算法

1、和分水岭算法相关的资料

[1] 标记符控制的分水岭算法原理及matlab实现_CSDN
[2] 图像处理——分水岭算法
[3]

2、《数字图像处理：原理与实践》

基于标记控制的分水岭分割方法有以下基本步骤：
1.计算分割函数。图像中较暗的区域是要分割的对象。
2.计算前景标志。这些是每个对象内部连接的斑点像素。
3.计算背景标志。这些是不属于任何对象的像素。
4.修改分割函数，使其仅在前景和后景标记位置有极小值。
5.对修改后的分割函数做分水岭变换计算。

（1）读入彩色图，转为为灰度图

rgb = imread('tudou500.jpg');
I = rgb2gray(rgb);

（2）使用Sobel算子边缘滤波

hy = fspecial('sobel');  % 水平 y 轴
hx = hy';  % 竖直 x
Iy = imfilter(double(I), hy, 'replicate');
Ix = imfilter(double(I), hx, 'replicate');
gradmag = sqrt(Ix.^2 + Iy.^2);
figure;
subplot(221);imshow(rgb);title('original');
subplot(222);imshow(gradmag, []);title('梯度幅值图像');
subplot(223);imshow(Iy, []);title('水平边缘检测');
subplot(224);imshow(Ix, []);title('垂直边缘检测');

（2-1）直接对梯度幅值图像gradmag进行watershed()操作

L = watershed(gradmag);
Lrgb = label2rgb(L);
figure;
subplot(211);imshow(I);title('original');
subplot(212);imshow(Lrgb);title('对梯度幅值直接做分水岭分割');

结果，过度分割达不到我们想要的效果。

（3）对图像进行增强

se = strel('disk',13);
Ie = imerode((I),se);
% Ie = imerode(imcomplement(I),se);  % 腐蚀
% Ie = imcomplement(Ie);
Iobr = imreconstruct(Ie,I);  % 重建
figure
subplot(231);imshow(I);title('original I');
subplot(232);imshow(Ie);title('腐蚀 Ie');
subplot(233);imshow(Iobr);title('基于开的重建 Iobr');

Iobrd = imdilate(Iobr, se);
Iobrdr = imreconstruct(imcomplement(Iobrd), imcomplement(Iobr));
Iobrdr = imcomplement( Iobrdr );
subplot(234);imshow(I);title('original I');
subplot(235);imshow(Iobrd);title('Iobr + dilate = Iobrd');
subplot(236);imshow(Iobrdr);title('基于重建的开闭操作 Iobrd+r=Iobrdr');

（4）获取前景标记

fgm = imregionalmax(Iobrdr);
It1 = rgb(:,:,1);
It2 = rgb(:,:,2);
It3 = rgb(:,:,3);
% 将fgm区域变成红色，其他区域不变
It1(fgm) = 255; It2(fgm) = 0;  It3(fgm) = 0;  
I2 = cat(3,It1, It2, It3);
figure
subplot(121);imshow(fgm,[]);title('局部极大值图像');
subplot(122);imshow(I2,[]);title('局部极大值叠加图像');

（5）极大值结果修正

se2 = strel(ones(12,12));
fgm2 = imcomplement(imclose(imcomplement(fgm), se2));
% fgm2 = imclose(fgm, se2);
fgm3 = imerode(fgm2, se2);
% fgm3 = imerode(fgm3, se2);
fgm3 = imdilate(fgm3, se2);
fgm3 = imdilate(fgm3, se2);
fgm3 = imerode(fgm3, se2);
% fgm3 = imopen(fgm2, se2);
fgm4 = bwareaopen(fgm3, 300);
figure
subplot(221);imshow(fgm);title('fgm');
subplot(222);imshow(fgm2);title('fgm2');
subplot(223);imshow(fgm3);title('fgm3');
subplot(224);imshow(fgm4);title('fgm4');

（6）使用分水岭分割watershed()

BW = im2bw(Iobrdr, graythresh(Iobrdr));
D = bwdist(BW);
DL = watershed(D);
bgm = DL == 0;  % DL == 0 确定矩阵DL等于0的值所在位置，=0的位置置为1
gradmag2 = imimposemin(gradmag, bgm | fgm4);


%% 分水岭分割
L = watershed(gradmag2);
% first way show
Lrgb = label2rgb(L, 'jet', 'w', 'shuffle');
figure
subplot(121);imshow(Lrgb);title('分水岭结果显示');

% second way show
subplot(122);imshow(rgb, []);title('分水岭结果显示（叠加）');
hold on;
himage = imshow(Lrgb);
set(himage, 'AlphaData', 0.3)

---

（上周就做完的实验，但是换方向，拖到了周六/(ㄒoㄒ)/~~）
end