matlab - How can I detect the maximum-sized rectangle that I can draw onto the mask?

Question

I'm making an image processing project and I have stuck in one the project's steps. Here is the situation;

This is my mask:

and I want to detect the maximum-sized rectangle that can fit into this mask like this.

I'm using MATLAB for my project. Do you know any fast way to accomplish this aim. Any code sample, approach or technique would be great.

EDIT 1 : The two algorithms below are works with lot's of the cases. But both of them give wrong results in some difficult cases. I'am using both of them in my project.

Answer 1

This approach starts with the entire image and shrinks each border in turn pixel-by-pixel until it finds an acceptable rectangle.

It takes ~0.02 seconds to run on the example image, so it's reasonably fast.

EDIT: I should clarify that this isn't meant to be a universal solution. This algorithm relies on the rectangle being centered and having roughly the same aspect ratio as the image itself. However, in the cases where it is appropriate, it is fast. @DanielHsH offered a solution which they claim works in all cases.

---

The code:

clear; clc;
tic;
%% // read image
imrgb= imread('box.png');
im = im2bw(rgb2gray(imrgb));    %// binarize image
im = 1-im;                      %// convert "empty" regions to 0 intensity
[rows,cols] = size(im);

%% // set up initial parameters
ULrow = 1;       %// upper-left row        (param #1)
ULcol = 1;       %// upper-left column     (param #2)
BRrow = rows;    %// bottom-right row      (param #3)
BRcol = cols;    %// bottom-right column   (param #4)

parameters = 1:4;   %// parameters left to be updated
pidx = 0;           %// index of parameter currently being updated

%% // shrink region until acceptable
while ~isempty(parameters); %// update until all parameters reach bounds

    %// 1. update parameter number
    pidx = pidx+1;
    pidx = mod( pidx-1, length(parameters) ) + 1;
    p = parameters(pidx);   %// current parameter number

    %// 2. update current parameter
    if p==1; ULrow = ULrow+1; end;
    if p==2; ULcol = ULcol+1; end;
    if p==3; BRrow = BRrow-1; end;
    if p==4; BRcol = BRcol-1; end;

    %// 3. grab newest part of region (row or column)
    if p==1; region = im(ULrow,ULcol:BRcol); end;
    if p==2; region = im(ULrow:BRrow,ULcol); end;
    if p==3; region = im(BRrow,ULcol:BRcol); end;
    if p==4; region = im(ULrow:BRrow,BRcol); end;

    %// 4. if the new region has only zeros, stop shrinking the current parameter
    if isempty(find(region,1))
        parameters(pidx) = [];
    end

end

toc;
params = [ULrow ULcol BRrow BRcol]
area = (BRrow-ULrow)*(BRcol-ULcol) 

---

The results for this image:

Elapsed time is 0.027032 seconds.

params =

    10    25   457   471


area =

      199362

Code to visualize results:

imrgb(params(1):params(3),params(2):params(4),1) = 0;
imrgb(params(1):params(3),params(2):params(4),2) = 255;
imrgb(params(1):params(3),params(2):params(4),3) = 255;
imshow(imrgb);

Another example image: