curve fitting and image segmentation

Hi all, I have a grayscale image-'FLIR0888.jpg' (attatched within) I have segmented this image into right and left breasts by fitting a polynomial curve and finding the intersection point of the two fitted curves. Now, below is the image of right breast and the curve fitted is shown in blue. How do i segment the grayscale image along this blue curve? I want the part of the image below this blue curve to be gone. How do I approach this? I have attatched the code and the original image below, thanks in advance.

clc;
clear;
close all;
a=rgb2gray(imread('FLIR0888.jpg'));
ed=edge(a,'canny',0.2,0.5);
bw1=bwareaopen(ed,10);
bw=imdilate(bw1,true(3));
se = strel('disk',4);
I=imclose(bw,se);
BW=imclearborder(I);
[B,L]= bwboundaries(BW,'noholes');
[~,I] = sort(cellfun(@length,B),'descend');
BB= B(I);
x1 = BB{1}(:, 2);y1 = BB{1}(:, 1);
x2 = BB{2}(:, 2);y2 = BB{2}(:, 1);
P1=polyfit(y1,x1,3);
P2=polyfit(y2,x2,3);
yy = linspace( 1, size(L,1));
p1=polyval( P1, yy );
p2=polyval( P2, yy);
x_intersect = fzero(@(x) polyval(P1-P2,x),3);
y_intersect = polyval(P1,x_intersect);
segimgR=a(:,1:y_intersect);
[m,n]=size(segimgR);
segimgl=a(:,y_intersect:end);
segimgL=imresize(segimgl,[m,n]);
figure()
subplot(2,2,1);imshow(a);title('original image');
subplot(2,2,2);imshow(a, 'border', 'tight' );
hold on
plot(p1,yy);
plot(p2,yy);
hold off
subplot(2,2,3);imshow(segimgR);
subplot(2,2,4);imshow(segimgL);
figure()
subplot(1,2,1);imshow(segimgR);
hold on
plot(p1,yy);
hold off

Answer 1

[Thiago Henrique Gomes Lobato]

In this case, for each column you have there's a different line point, so you have to check for each one individually and do the threshold accordingly. Here is a implementation that does this, if you ever want another direction it is equivalent:

% Loop for all pixels in one direction
MaskUP = zeros(size(L));
AllDirectionIndex = 1:size(L,2);
DivisionPoints=polyval( P1, 1:size(L,1) );
for idx=1:size(L,1)
AboveIndex   = AllDirectionIndex<DivisionPoints(idx);
MaskUP(idx,AboveIndex) = 1;
end

figure,imshow( uint8(MaskUP.*double(a) ))

[Apoorva Maiya]

Thank you so much sir, i was able to segment both right and left images seperately.

[Cansu Kelebek]

Hi, Could you please help me to separate the left and right half of breast? There is code below which I tried to seperate both left and right of breast. @Apoorva Maiya@Thiago Henrique Gomes Lobato

clc; clear; close all;
a=rgb2gray(imread('DINAMIC-FRONTAL1.jpg'));
en= imsharpen(a,'Radius',2,'Amount',1);
B = imgaussfilt(en,1.4);
ed=edge(B,'canny',0.3,0.5);
figure();imshow(ed);title('canny edge');
bw1=bwareaopen(ed,10);
se = strel('disk',4);
bw=imdilate(bw1,se);
figure();imshow(bw);
[y x] = find( bw );
right = find(min(a));
as=findpeaks(y);
xr = x(right);
yr = y(right);
xl = x(as);
yl = y(as);
pr = polyfit( yr, xr, 2);  %// fit 2rd deg poly
pl = polyfit( yl, xl, 2 ); 
yy = linspace( 1, size(bw,1), 50 );
figure; imshow(a, 'border', 'tight' );
hold all
plot( polyval( pr, yy ), yy, '.-', 'LineWidth', 1 );
plot( polyval( pl, yy ), yy, '.-', 'LineWidth', 1 );

Answer 2

[Image Analyst]

@Cansu Kelebek, try this to find the dividing line between the 2 breasts:

% Demo by Image Analyst
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 short g;
format compact;
fontSize = 20;

grayImage = rgb2gray(imread('DINAMIC-FRONTAL.jpeg'));
en= imsharpen(grayImage,'Radius',2,'Amount',1);
B = imgaussfilt(en,1.4);
ed=edge(B,'canny',0.3,0.5);
subplot(2, 2, 1);
imshow(ed);
title('Initial Canny Edge', 'FontSize', fontSize);
bw1=bwareaopen(ed,10);
se = strel('disk', 4);
bw=imdilate(bw1,se);
subplot(2, 2, 2);
imshow(bw);
title('Canny Edge Dilated', 'FontSize', fontSize);

% Non-working code from poster commented out.
% [y x] = find( bw );
% right = find(min(a));
% as=findpeaks(y);
% xr = x(right);
% yr = y(right);
% xl = x(as);
% yl = y(as);
% pr = polyfit( yr, xr, 2);  %// fit 2rd deg poly
% pl = polyfit( yl, xl, 2 ); 
% yy = linspace( 1, size(bw,1), 50 );
% subplot(2, 2, 3);
% imshow(a, 'border', 'tight' );
% hold all
% plot( polyval( pr, yy ), yy, '.-', 'LineWidth', 1 );
% plot( polyval( pl, yy ), yy, '.-', 'LineWidth', 1 );

% Take the 2 largest blobs.
bw = bwareafilt(bw, 2);
subplot(2, 2, 3);
imshow(bw);
title('2 Largest Blobs', 'FontSize', fontSize);

% Find the bounding boxes
props = regionprops(bw, 'BoundingBox');
bb = vertcat(props.BoundingBox);
xLeft = bb(:, 1);
widths = bb(:, 3);
% Find the right edge of the left blob.
x1 = xLeft(1) + widths(1);
% Find the midpoint between the left edge of the right blob and the right edge of the left blob.
xMiddle = (x1 + xLeft(2)) / 2;
% Put a red line between them
xline(xMiddle, 'Color', 'r', 'LineWidth', 3);

subplot(2, 2, 4);
imshow(grayImage);
title('2 Largest Blobs', 'FontSize', fontSize);
% Put a red line between them
xline(xMiddle, 'Color', 'r', 'LineWidth', 3);
title('Original with Dividing Line', 'FontSize', fontSize);

[Cansu Kelebek]

Thanks for your help but i can’t appropriate edge detection for different breast in size and shape. Could you please help me for other images? @Image Analyst

[Image Analyst]

Why can't you just say the middle is in the middle?

xMiddle = size(grayImage, 2) / 2;

What's wrong with that? It might not be exactly where you'd place it but it could be good enough for what comes next. So, let's say you had the middle. What comes next?

[Cansu Kelebek]

Actually i didn’t get what you mean. I tried vertical projection profile method to find lower and upper borders of breast but i couldn’t get appropriate result. Also I tried to detech branch of breast but i couldn’t get acceptable result for every image due to inappropriate edge detection of all images. And I need one code for every breast thermogram

[Cansu Kelebek]

Also I have some info about that but i couldn’t apply it.

Conditions of detecting 4 curves of breast, circle fitting and breast region: Since the breast boundaries are parabolic in shape, a 2nd order polynomial fitting algorithm was used to detect the bottom breast boundaries. A parabolic curve was first fitted to each object, then two curves with maximum coefficient of determination candidate for the lower breast boundaries. If the difference between x-coordinates of the lowest peaks of two curves was less than 20 and the difference of their y-coordinates was greater than 250, they were selected as lower breast boundaries (X–X ′ < 20 and Y–Y ′ > 250). Otherwise, the algorithm searched the other curves. The next step is to determine the top of the breast. In Qi and co-workers [19,21] study, the region with the highest curvature in the left and right body edges was used to estimate the position of the armpit. Often, poor edge detection indicated major points of concavity elsewhere in the body edge, and the point of maximum concavity in the image was found to be lower or higher than what was expected. Besides, choosing the fitted parabola as lower breast boundaries may be lead to crop the region of the breast and cause information loss because of the poor lower breast edge in many cases. This work proposed to fit a circle to each lower edge object for delimiting the breasts in the region of the breasts in the image. The circle-fitting algorithm was taken from Kasa method [22]. The same procedure was repeated for both sides. After detection of the breast regions, a middle separation line was created to separate between left and right segments along the point of contact of the two breast curves.

[Image Analyst]

If this group's algorithm isn't working for you then you need to just use drawpoints to manually pick some points on the lower breast boundary. Then fit your circle. Obviously this algorithm won't work for some images, like where the woman is fairly flat chested. I assume you're going to have to have something that works for all images. So if you just want to get moving, I'd suggest doing it manually. By this I mean you want to analyze the images somehow and the main point of your research is not to develop a superior breast finding algorithm to the published ones.

[Cansu Kelebek]

Actually i couldn’t apply that algorithm including curve fitting (x-x’<20,y-y’>250) and circle fitting. I don’t know how to apply this conditions below

Since the breast boundaries are parabolic in shape, a 2nd order polynomial fitting algorithm was used to detect the bottom breast boundaries. A parabolic curve was first fitted to each object, then two curves with maximum coefficient of determination candidate for the lower breast boundaries. If the difference between x-coordinates of the lowest peaks of two curves was less than 20 and the difference of their y-coordinates was greater than 250, they were selected as lower breast boundaries (X–X ′ < 20 and Y–Y ′ > 250).

[Cansu Kelebek]

But i have to do it automatically not manually.

[Image Analyst]

Why? Are you wanting to develop and publish your own algorithm that is better than Qi's?

[Cansu Kelebek]

That’s my thesis project so i have to seperate left half and right half of breast automatically. Maybe I’ll use it for some images that as far as i can.

[Image Analyst]

Oh, ok. Well if it's the main subject of your thesis then I can't do it for you. A thesis is something you have to do on your own, with guidance from your professor. But good luck with it.

[Cansu Kelebek]

Nope it’s just a small part of my thesis. But I couldn’t pass this part.