How to do correctly this for loop for segmentation?

Question

Veronika 2017-3-4

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/328046-how-to-do-correctly-this-for-loop-for-segmentation

评论： Veronika 2017-3-7

在 MATLAB Online 中打开

Dear all,

I have code for segmentation specific parts. Output from this code is this image:

I would like to gradually in this manner do segmentation the other parts. I think, that another for loop could help. I tried this one:

*for k = 1:16*
cnt=0;
for i=1:numel(tmp)
    if numel(tmp{i})>4
    cnt=cnt+1;
    tmp2{cnt}=tmp{i};
   end
end
tmp = tmp2;
*end*

But it didn´t solve my problem. I would like to for loop, which works like this: for k = 1 this output:

for k = 2 This output:

and so on. I hope that you understand me. I attach my code and original image. Thank you for your answers.

2 个评论
显示无隐藏无

Image Analyst 2017-3-4

I don't really understand. What regions do you want? So many times in the past I've done things like finding lungs, skulls, tumors, etc. Search my old posts.

Veronika 2017-3-4

I need segmentation like this

but bloob after bloob. Because I need to filter found objects by size.

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Image Analyst 2017-3-4

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/328046-how-to-do-correctly-this-for-loop-for-segmentation#answer_257297

Here's the plan (which I can't do because you didn't post the original, non-annotated image):

Threshold to find body.
Use bwareafilt() to extract largest blob.
Use that to mask out non-body stuff in the gray level image.
Threshold again at a higher gray level
Use bwareaopen() ror bwareafilt() to get only the blobs of the size you want.
Call bwboundaries() to get the perimeters of the blobs.
Use a for loop and plot() to plot red outlines around the blobs.

Good luck. If you need more help, insert your original image and post your code.

1 个评论
显示 -1更早的评论隐藏 -1更早的评论

Veronika 2017-3-5

在 MATLAB Online 中打开

Okay, I attach original image and code. I think, that I could use this loop

cnt=0;
for i=1:numel(tmp)
    if numel(tmp{i})>4
    cnt=cnt+1;
    tmp2{cnt}=tmp{i};
   end
end
tmp = tmp2;

for my output. Because this for loop I used for one yellow segment i my output, so I thought, that this could be the right way. Thank you for your help.

请先登录，再进行评论。

Answer 2

Image Analyst 2017-3-5

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/328046-how-to-do-correctly-this-for-loop-for-segmentation#answer_257347

在 MATLAB Online 中打开

I didn't understand your comments, so I just wrote the whole thing myself.

clc;    % Clear the command window.
close all;  % Close all figures (except those of imtool.)
clear;  % Erase all existing variables. Or clearvars if you want.
workspace;  % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 15;
%===============================================================================
% Get the name of the image the user wants to use.
baseFileName = 'thorax-mdl.jpg'; % Assign the one on the button that they clicked on.
% Get the full filename, with path prepended.
folder = []; % Determine where demo folder is (works with all versions).
fullFileName = fullfile(folder, baseFileName);
%===============================================================================
% Read in a demo image.
grayImage = imread(fullFileName);
% Get the dimensions of the image.  
% numberOfColorChannels should be = 1 for a gray scale image, and 3 for an RGB color image.
[rows, columns, numberOfColorChannels] = size(grayImage);
if numberOfColorChannels > 1
  % It's not really gray scale like we expected - it's color.
  % Use weighted sum of ALL channels to create a gray scale image.
  grayImage = rgb2gray(grayImage); 
  % ALTERNATE METHOD: Convert it to gray scale by taking only the green channel,
  % which in a typical snapshot will be the least noisy channel.
  % grayImage = grayImage(:, :, 2); % Take green channel.
end
% Display the image.
subplot(2, 3, 1);
imshow(grayImage, []);
axis on;
caption = sprintf('Original Gray Scale Image');
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
drawnow;
hp = impixelinfo();
% Let's compute and display the histogram.
[pixelCount, grayLevels] = imhist(grayImage);
subplot(2, 3, 2); 
bar(grayLevels, pixelCount); % Plot it as a bar chart.
grid on;
title('Histogram of original image', 'FontSize', fontSize, 'Interpreter', 'None');
xlabel('Gray Level', 'FontSize', fontSize);
ylabel('Pixel Count', 'FontSize', fontSize);
xlim([0 grayLevels(end)]); % Scale x axis manually.
% Set up figure properties:
% Enlarge figure to full screen.
set(gcf, 'Units', 'Normalized', 'OuterPosition', [0 0 1 1]);
% Get rid of tool bar and pulldown menus that are along top of figure.
set(gcf, 'Toolbar', 'none', 'Menu', 'none');
% Give a name to the title bar.
set(gcf, 'Name', 'Demo by ImageAnalyst', 'NumberTitle', 'Off') 
drawnow;
% Binarize the image by thresholding.
bodyMask = grayImage > 128;
% Get rid of white surround.
bodyMask = imclearborder(bodyMask);
% Extract the largest blob only.
bodyMask = bwareafilt(bodyMask, 1);
subplot(2, 3, 3);
imshow(bodyMask);
axis on;
title('Binary Image of Whole Body', 'fontSize', fontSize);
drawnow;
% Find the ribs mask by thresholding and ANDing with the body mask
% to make sure we don't get anything from inside the lungs.
ribsMask = bodyMask & (grayImage > 220);
% Fill holes.
ribsMask = imfill(ribsMask, 'holes');
% Find areas.
labeledImage = bwlabel(ribsMask);
props = regionprops(labeledImage, 'Area');
allAreas = sort([props.Area])
% Get rid of areas below 100 pixels.
ribsMask = bwareaopen(ribsMask, 100);
% Display the image.
subplot(2, 3, 4);
imshow(ribsMask, []);
axis on;
caption = sprintf('Ribs Mask');
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
drawnow;
% bwboundaries() returns a cell array, where each cell contains the row/column coordinates for an object in the image.
% Plot the borders of all the coins on the original grayscale image using the coordinates returned by bwboundaries.
subplot(2, 3, 5);
imshow(grayImage);
title('Outlines, from bwboundaries()', 'FontSize', fontSize); 
axis image; % Make sure image is not artificially stretched because of screen's aspect ratio.
hold on;
boundaries = bwboundaries(ribsMask);
numberOfBoundaries = size(boundaries, 1);
for k = 1 : numberOfBoundaries
  thisBoundary = boundaries{k};
  plot(thisBoundary(:,2), thisBoundary(:,1), 'r', 'LineWidth', 2);
end
hold off;

5 个评论
显示 3更早的评论隐藏 3更早的评论

Veronika 2017-3-6

编辑：Veronika 2017-3-6

在 MATLAB Online 中打开

This your segmentation I already have, but it dosn´t solve my problem. I would like to create FEM model from this segmentation. But some blobs are too small (1 pixel) for reconstruction, so I think, that I have to choose only bloobs bigger than 1 pixel. This is my code:

clc;    % Clear the command window.
close all;  % Close all figures (except those of imtool.)
clear;  % Erase all existing variables. Or clearvars if you want.
workspace;  % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 15;
%===============================================================================
% Get the name of the image the user wants to use.
baseFileName = 'thorax-mdl.jpg'; % Assign the one on the button that they clicked on.
% Get the full filename, with path prepended.
folder = []; % Determine where demo folder is (works with all versions).
fullFileName = fullfile(folder, baseFileName);
%===============================================================================
% Read in a demo image.
grayImage = imread(fullFileName);
% Get the dimensions of the image.  
% numberOfColorChannels  = 1; 
[rows, columns, numberOfColorChannels] = size(grayImage);
if numberOfColorChannels > 1
  % It's not really gray scale like we expected - it's color.
  % Use weighted sum of ALL channels to create a gray scale image.
grayImage = rgb2gray(grayImage); 
  % ALTERNATE METHOD: Convert it to gray scale by taking only the green channel,
  % which in a typical snapshot will be the least noisy channel.
%   grayImage = grayImage(:, :, 2); % Take green channel.
end
% Display the image.
eq_grayImage = histeq(grayImage);
subplot(2, 3, 1);
imshow(grayImage, []);
axis on;
caption = sprintf('Original Gray Scale Image');
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
drawnow;
hp = impixelinfo();
% Let's compute and display the histogram.
[pixelCount, grayLevels] = imhist(grayImage);
subplot(2, 3, 2); 
bar(grayLevels, pixelCount); % Plot it as a bar chart.
grid on;
title('Histogram of original image', 'FontSize', fontSize, 'Interpreter', 'None');
xlabel('Gray Level', 'FontSize', fontSize);
ylabel('Pixel Count', 'FontSize', fontSize);
xlim([0 grayLevels(end)]); % Scale x axis manually.
% Set up figure properties:
% Enlarge figure to full screen.
set(gcf, 'Units', 'Normalized', 'OuterPosition', [0 0 1 1]);
% Get rid of tool bar and pulldown menus that are along top of figure.
set(gcf, 'Toolbar', 'none', 'Menu', 'none');
% Give a name to the title bar.
set(gcf, 'Name', 'Demo by ImageAnalyst', 'NumberTitle', 'Off') 
drawnow;
% Binarize the image by thresholding.
bodyMask = grayImage > 128;
% Get rid of white surround.
bodyMask = imclearborder(bodyMask);
% Extract the largest blob only.
bodyMask = bwareafilt(bodyMask, 1);
subplot(2, 3, 3);
imshow(bodyMask);
axis on;
title('Binary Image of Whole Body', 'fontSize', fontSize);
drawnow;
% Find the ribs mask by thresholding and ANDing with the body mask
% to make sure we don't get anything from inside the lungs.
ribsMask = bodyMask & (grayImage > 220);
% Fill holes.
ribsMask = imfill(ribsMask, 'holes');
% Find areas.
labeledImage = bwlabel(ribsMask);
props = regionprops(labeledImage, 'Area');
allAreas = sort([props.Area])
% Get rid of areas below 100 pixels.
ribsMask = bwareaopen(ribsMask, 100);
% Display the image.
subplot(2, 3, 4);
imshow(ribsMask, []);
axis on;
caption = sprintf('Ribs Mask');
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
drawnow;
% bwboundaries() returns a cell array, where each cell contains the row/column coordinates for an object in the image.
% Plot the borders of all the coins on the original grayscale image using the coordinates returned by bwboundaries.
figure(4)
imshow(grayImage);
title('Outlines, from bwboundaries()', 'FontSize', fontSize); 
axis image; % Make sure image is not artificially stretched because of screen's aspect ratio.
hold on;
set(gcf, 'Units', 'Normalized', 'OuterPosition', [0 0 1 1]);
% Odstranění panelu nástrojů a rozbalovacího menu.
set(gcf, 'Toolbar', 'none', 'Menu', 'none');
set(gcf, 'Name', 'Segmentace páteře', 'NumberTitle', 'Off') 
boundaries = bwboundaries(ribsMask);
numberOfBoundaries = size(boundaries, 1);
for k = 1 : numberOfBoundaries
  thisBoundary = boundaries{k};
  plot(thisBoundary(:,2), thisBoundary(:,1), 'r', 'LineWidth', 2);
end 
hold off;
%Práh pro vytvoření binárního obrazu okolí
thresholdValue = 180;
binaryImage_okoli = eq_grayImage > thresholdValue;
% Odstranění okolí.
binaryImage_okoli = imclearborder(binaryImage_okoli);
% Vyplnění otvorů.
binaryImage_okoli = imfill(binaryImage_okoli, 'holes');
% Vymazání menších otvorů.
binaryImage_okoli = bwareaopen(binaryImage_okoli, 750);
figure(3)
subplot (2, 3, 1)
imshow(binaryImage_okoli, []);
axis on;
caption = sprintf('Binární obraz okolí');
title(caption, 'FontSize', fontSize, 'Interpreter', 'None');
run C:/Users/ced0005/eidors-v3.8-ng/eidors/startup.m
if 2==exist ('Pater_model.mat')
load ('Pater_model')
disp('Model páteře byl načten.')
else
shape = [];
tmp = bwboundaries(binaryImage_okoli);
shape.thorax = tmp{1};
tmp = bwboundaries(ribsMask);
shape.spine = tmp{1};
size(grayImage)   
% cnt=0;
% for i=1:numel(tmp)
%     if numel(tmp{i})>4
%     cnt=cnt+1;
%     tmp2{cnt}=tmp{i};
%    end
% end
% tmp = tmp2;
% interpolace hraničních bodů, fourier
Npts1 = 35;      % počet bodů interpolovaných hranic
Npts2 = 25;      % počet bodů interpolovaných hranic
Npts3 = 25;       % počet bodů interpolovaných hranic
Ncomps = 20;     % počet fourierových komponentů
shape.thorax = fourier_fit(fourier_fit(shape.thorax, Ncomps), linspace(0, 1, Npts1));
shape.spine = fourier_fit(fourier_fit(shape.spine, Ncomps), linspace(0, 1, Npts3));
figure(7);
imshow(grayImage);
hold on
plot(shape.thorax(:,2), shape.thorax(:,1), 'o-b')
plot(shape.spine(:,2), shape.spine(:,1), 'o-y')
legend('Okolí hrudníku', 'Páteř + žebra');

As you can see, if I put as input ribsMask, the output of function fourier_fit is only one rib, but I need to all ribs and spine. Can you advise me?

Image Analyst 2017-3-6

I don't know what that code does. What is its purpose? And I don't understand what the code after my code does because I can't read that language. So either convert to English so I can read and understand what's going on, or try yourself by just adjusting the parameters to bwareaopen() or all those other parameters you introduced like Ncomps, etc.

Veronika 2017-3-7

I tried change the parameters, but without success. I attach english version of my code and original image. I hope, that it will be enough understanding for you.

请先登录，再进行评论。