pruning image segments' leftovers

Question

As you see in images below, in some of my segmentation result(segmentation done by watershed transformation method), there are some leftovers left. I want to somehow crop the images so that only the rectangles remain. This operation is only based on rectangle shape and it doesn't relate to intensity level.

Answer 1

Solution explanation

I suggest the following approach:

1. generate an initial guess for the 4 corners of the shape according to geometric properties (see code below for further details).

2. create a quadrangular given these 4 corners, by drawing a line between each pair of corresponding corners.

3. find the corners which optimize the Jaccard coefficient of the boundary image and the generated quadrangular map.

In order to spare time, The optimization stage will be done locally. We will try to replace each corner with the best achievable corner in a certain neighborhood. We stop the optimization stage if there is no improvement for each of the 4 corners.

code

%reads image
gray = rgb2gray(imread('Bqx51.png'));
I = gray>0;

%extract boundries
B = bwboundaries(I,8);
B = B{1};
boundriesImage = zeros(size(I));
boundriesImage(sub2ind(size(I),B(:,1),B(:,2))) = 1;

%finds best 4 corners
[ corners ] = optimizeCorners(B);

%generate line mask
linesMask =  drawLines(size(I),corners,corners([2:4,1],:));

%fill holes
rectMask = imfill(linesMask,'holes');

%noise reduction
rectMask = I & rectMask;
rectMask = imopen(rectMask,strel('disk',2));

%calculate result 
result = gray;
result(~rectMask) = 0;

%display results
figure,imshow([gray, 255*ones(size(I,1),1),result]);

corner optimization function

function [ corners] = optimizeCorners(pnts)
%OPTIMIZE4PTS Summary of this function goes here
%   Detailed explanation goes here

Y = pnts(:,1);
X = pnts(:,2);

corners = getInitialGuess(X,Y); 
boundriesIm = zeros(max(Y),max(X));
boundriesIm(sub2ind(size(boundriesIm),pnts(:,1),pnts(:,2))) = 1;

%R represents the search radius
R = 3;

%continue optimizing as long as there is no change in the final result
unchangedIterations = 0;
while unchangedIterations<4

    for ii=1:4
        %optimize corner ii
        currentCorner = corners(ii,:);
        bestCorner = currentCorner;
        bestRes = calcEnergy(boundriesIm,corners);
        cornersToEvaluate = corners;
        candidateInds = sum(((repmat(currentCorner,size(X,1),1)-[Y,X]).^2),2)<(R^2);
        candidateCorners = [Y(candidateInds),X(candidateInds)];
        for jj=length(candidateCorners)
            xx = candidateCorners(jj,2);
            yy = candidateCorners(jj,1);
            cornersToEvaluate(ii,:) = [yy,xx];
            res = calcEnergy(boundriesIm,cornersToEvaluate);
            if res > bestRes
                bestRes = res;
                bestCorner = [yy,xx];
            end
        end
        if isequal(bestCorner,currentCorner)
            unchangedIterations = unchangedIterations + 1;
        else
            unchangedIterations = 0;
            corners(ii,:) = bestCorner;

        end
    end
end

end

Calculate energy function

function res = calcEnergy(boundriesIm,corners)
%calculates the score of the corners list, given the boundries image.
%the result is acutally the jaccard index of the boundries map and the
%lines map
linesMask =  drawLines(size(boundriesIm),corners,corners([2:4,1],:));
res = sum(sum(linesMask&boundriesIm)) / sum(sum(linesMask|boundriesIm));

end

finding initial guess for corners function

function corners = getInitialGuess(X,Y)
%calculates an initial guess for the 4 corners

corners = zeros(4,2);

%preprocessing stage
minYCoords = find(Y==min(Y));
maxYCoords = find(Y==max(Y));
minXCoords = find(X==min(X));
maxXCoords = find(X==max(X));
%top corners
topRightInd = find(X(minYCoords)==max(X(minYCoords)),1,'last');
topLeftInd = find(Y(minXCoords)==min(Y(minXCoords)),1,'last');
corners(1,:) = [Y(minYCoords(topRightInd)) X((minYCoords(topRightInd)))];
corners(2,:) = [Y(minXCoords(topLeftInd)) X((minXCoords(topLeftInd)))];
%bottom corners
bottomRightInd = find(Y(maxXCoords)==max(Y(maxXCoords)),1,'last');
bottomLeftInd = find(X(minYCoords)==min(X(minYCoords)),1,'last');
corners(4,:) = [Y(maxXCoords(bottomRightInd)) X((maxXCoords(bottomRightInd)))];
corners(3,:) = [Y(maxYCoords(bottomLeftInd)) X((maxYCoords(bottomLeftInd)))];


end

drawLine function (taken from the following answer, by @Suever)

function mask = drawLines(imgSize, P1, P2)
%generates a mask with lines, determine by P1 and P2 points

mask = zeros(imgSize);

P1 = double(P1);
P2 = double(P2);

for ii=1:size(P1,1)
    x1 = P1(ii,2); y1 = P1(ii,1);
    x2 = P2(ii,2); y2 = P2(ii,1);

    % Distance (in pixels) between the two endpoints
    nPoints = ceil(sqrt((x2 - x1).^2 + (y2 - y1).^2));

    % Determine x and y locations along the line
    xvalues = round(linspace(x1, x2, nPoints));
    yvalues = round(linspace(y1, y2, nPoints));

    % Replace the relevant values within the mask
    mask(sub2ind(size(mask), yvalues, xvalues)) = 1;
end

Results

first image (before and after):

second image (before and after):

Runtime

Elapsed time is 0.033998 seconds.

Possible improvements / suggestions

1. The energy function may also include constraints which encourage parallel lines to have the similar slopes (in your example they don't have the same slope).

2. Energy function may include constraints which encourage each corner angles to be close to 90 degrees.

3. noise reduction stages (such as imclose) can be done before performing this approach to eliminate small artifacts.

4. It is possible to run the algorithm with several initial guesses and choose the best one.

5. Notice that this solution doesn't estimates the best possible rectangle - it estimates the best quadrangular. The reason is that the input images are not rectangles (the lines are not parallel).