Fit countours to shape

Yo, image processing experts, any way to automatically fit the 4 boundary lines of this shape:

This what I want:

I believe this should be easy.

Any ideas? This should be done in less than 1 second.

[John D'Errico]

It is amazing how what the human brain thinks is easy, even trivial, is not so trivial at all when done by computer. I would note that you only THINK it should be easy. Perhaps you have been watching too many TV shows, where computers can do anything? Where any problem can be solved using now more than a couple of keystrokes?

[André]

So this is not easy?

[André]

Wow, so many diverse answers. Thank you a lot. I will evaluate all them and choose the one I think is best.

Personally, I tought of using the Hough transform to detect the almost vertical lines, but I am reluctant of using it as it is a somewhat computationally demanding algorithm.

[John D'Errico]

NO. It is not at all trivial to do what your brain thinks is easy. Obviously, if it was as easy as you seem to think, you would already have written the code.

What are the issues?

1. Some of the sides of that "polygon" are curved. Straight lines are one thing. But a general curved arc is less easy.
2. Those curved lines were constructed so they seem to ignore SOME of the bumps in the black regions. So sort of a least squares fit, but that would involve a fit with deviations in both x and y. And that alone is a little more difficult, but doable, but when that involves a curved line, things go all to hell.
3. There is stuff inside the region that you have blithely chosen to ignore. Again, you know what you want to see. But writing code for what you want to see is not as easy as you think.
4. The curves for each side intersect in a point that is sometimes unconstrained by your data. Essentially, you are extrapolating those curved arcs to where they intersect. But extrapolation is a difficult thing to accomplish (at least, doing extrapolation well by any definition of the word "well".)
5. And you want code that does all of this automatically in less than a second. Sigh.
6. I'm sure there is a number 6, but we can stop here.

Again, if it is easy to do, then start writing! Surely you know how to solve all of the issues I described? As I said, I think we all see too many TV shows and movies where the computers are godlike. That may be coming one day, sooner than I like, but...

My recommendation is an easy one. Draw the lines yourself. Take advantage of the processor in your brain that knows exactly what you want to see. You can use a tool like ginput to pick points off the figure, then connect them in a polyfon. Easy, peasy.

The result will be just as good as any code you will write. Yes, it will not be as fast as you want. Life sucks.

Answer 1

[Mathieu NOE]

arrgh - I don't have time to finish my work - it's getting very late here ...

and I am not an expert in image processing also so I guess someone will post an answer with 3 lines of code and I will look a bit dumb.

anyway this is the result so far and my code if that interest you

basically I could fit a parabola on the top and bottom boundary

what I still need to do is find a way to fit a straight line on the peaks of the magenta and cyan data segments, connect all four boundaries

but now it's time for me to leave and I'm coming back only in january so happy christmas !!!!!!!!!!!!!

A = imread('image.png');
A = double(A(:,:,1));
A = flipud(A);
[m,n] = size(A);

% transform the surrounding blank frame into black - ortherwise function top_bottom_boundary will not work fine...
ll = 10; % left / right / bottom / top lines thickness in pixels to be converted.
B = zeros(size(A));
B(ll:end-ll,ll:end-ll) = A(ll:end-ll,ll:end-ll);

%% main code
[y,x] = find(B>0.5*256); % find black dots coordinates : threshold is set at 50% of 256
[x3,y3,x4,y4] = top_bottom_boundary(x,y);  % top and bottom boundary
[y3b,x3b,y4b,x4b] = top_bottom_boundary(y,x);  %left and right  boundary


% top line fit
yfit3 = parabola_fit(x3,y3);
ylow_top = min(yfit3);

% bottom line fit
yfit4 = parabola_fit(x4,y4);
yup_bot = max(yfit4);

% remove point of right boundary that are beyong the parabolas
ind = (y3b>ylow_top) | (y3b<yup_bot);
x3b(ind) = [];
y3b(ind) = [];
% remove points that are not on the "top of the hill" 
ind = x3b>(min(x3b) + 25);
x3b(ind) = [];
y3b(ind) = [];

% remove point of left boundary that are beyong the parabolas
ind = (y4b>ylow_top) | (y4b<yup_bot);
x4b(ind) = [];
y4b(ind) = [];
% remove points that are not on the "top of the hill" 
ind = x4b<(max(x4b) - 25);
x4b(ind) = [];
y4b(ind) = [];

figure(1)
imagesc(B);
set(gca,'YDir','normal');
colormap('gray');
% colorbar('vert');
hold on
plot(x3b,y3b, '-c',x4b,y4b, '-m');% left and right  boundary transformed into linear segments ? 
plot(x3, yfit3, '-r', x4, yfit4, '-g','linewidth',5); % top and bottom boundary transformed with parabola fit
hold off

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function yfit = parabola_fit(x,y)
    % PARABOLA FIT
    m = length(x);
    %Set up the appropriate matrix A to find the best-fit parabola of the form y=C+Dx+Ex^2.  The 
    %first column of A will contain all 1's, using the ones() command.  The second column of A 
    %contains x values that are stored in X.  The third column of A contains the squared x values
    %that are stored in X.  Elementwise multiplication of X by itself, using .* operator, will 
    %produce the desired values for the third column.  
    A = [ones(m,1) x x.*x];
    A_transposeA = A.' * A;
    A_transposeY = A.' * y;
    %backslash operation to solve the overdetermined system.
    Soln2 = A_transposeA\A_transposeY;
    %x values to use for plotting the best-fit parabola.  
    yfit = Soln2(1) + Soln2(2)*x + Soln2(3)*x.*x; %
    % 
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [x3,y3,x4,y4] = top_bottom_boundary(x,y)
% based on FEX : https://fr.mathworks.com/matlabcentral/answers/299796-tight-boundary-around-a-set-of-points
    %split data into classes and find max/min for each class
    class_label = unique(x);
    upper_boundary = zeros(size(class_label));
    lower_boundary = zeros(size(class_label));
    for idx = 1:numel(class_label)
        class = y(x == class_label(idx));
        upper_boundary(idx) = max(class);
        lower_boundary(idx) = min(class);
    end
    left_boundary = y(x == class_label(1));
    right_boundary = y(x == class_label(end));

%     % left_boundary
%     x1 = class_label(1)*ones(size(left_boundary));
%     y1 = left_boundary;
% 
%     % right_boundary
%     x2 = class_label(end)*ones(size(right_boundary));
%     y2 = right_boundary;
    % top boundary
    x3 = class_label;
    y3 = upper_boundary;
    % bottom boundary
    x4 = class_label;
    y4 = lower_boundary;
end

[André]

Excellent method. I had tought of using a quadratic fit, but it seemed challenging. You helped me a lot. I also used a linear fit for the (almost) vertical borders. Thanks a lot.

[Mathieu NOE]

My pleasure and happy new year !!!

Answer 2

[Cris LaPierre]

I made an attempt using the Image Segmenter app. I first thresholded the image, then eroded the mask to get a roughly square shape, filled in the holes, then dilated the mask using a square shape that touched all 4 sides.

Once I had that, I borrowed some code from the Boundary Tracing in Images example to generate an outline.

The solution is defniitely customized for this particular image, but may work for similar images.

img = imread('https://www.mathworks.com/matlabcentral/answers/uploaded_files/1821389/image.png');
imshow(img)

% Threshold image with global threshold
BW = imbinarize(im2gray(img));

% Erode mask with square
width = 100;
se = strel('square', width);
BW = imerode(BW, se);

% Fill holes
BW = imfill(BW, 'holes');

% Dilate mask with square
width = 100;
se = strel('square', width);
BW = imdilate(BW, se);

boundaries = bwboundaries(BW);
b = boundaries{1};
hold on
plot(b(:,2),b(:,1),"g",LineWidth=3);

●

[André]

Excellent approach. I will try it with other images. However, I have to chose the next answer because using a fitting algorithm may be more robust when used with different images.

Answer 3

[André]

Thank you everybody fot the excellent answers. Unfortunately I could choose only one. If anyone still have different ideas, feel free to post, from the high quality of the answers so far, it would be nice to view all them.

The shape I was segmenting is the human crystalline by the way.