Suppose you have an image stored in the variable called im. It is also assumed that im is in double format [0-1]. Let's step through each part one at a time, with accompanying code.
(1) Select two thresholds T1 and T2.
Let's say T1 = 0.1 and T2 = 0.7. We need to make sure that these thresholds are at the lower and higher end of the intensity spectrum.
T1 = 0.1;
T2 = 0.7;
(2) Partition the image into three type regions :
R1, containing all pixels with gray values below T1.
R2, containing all pixels with gray values between T1 and T2, inclusive.
R3, containing all pixels with gray values above T2.
where R1=core region, R2= fringe(intermediate, transition) region, R3=background.
We will create three binary maps that will tell us where in the image each of these conditions are satisfied.
R1 = im < T1;
R2 = im >= T1 & im <= T2;
R3 = im > T2;
(3) Visit each pixel in region R2. If the pixel has a neighbor in region R1, then reassign the
pixel to region R1.
I'm going to assume 8-connectedness, which means that a neighbour would be N, NE, E, SE, S, SW, W, or NW, using cardinal directions. What we can do is take each overlapping neighbourhood of pixels for R2 and place them into columns. Do the same for R1. We can use the function called im2col to help us do that. 8-connectedness checks for neighbours within a 3 x 3 pixel neighbourhood. What's cool about im2col is that the for the output of this function, middle row will contain the centre of each neighbourhood. We will definitely need this when we want to reconstruct our regions back. When we extract the centre of each neighbourhood, it's just a matter of checking to see if:
- If the centre of each neighbourhood in region
R2 evaluates to true
- Any pixels in the corresponding neighbourhood of
R1 evaluates to true, then these are all neighbours of R2.
%//Zero-pad the regions so we can check the borders
R1 = padarray(R1, [1 1], 'replicate');
R2 = padarray(R2, [1 1], 'replicate');
%//Transform into columns
B1 = im2col(R1, [3 3]);
B2 = im2col(R2, [3 3]);
%//Extract size for later
[rows cols] = size(R1);
%// Extract centre of each neighbourhood
middleB2 = B2(5,:);
% // Find the indices of those neighbourhoods that have 1 in the centre
windowsHaving1 = find(middleB2 == 1);
% // Access the corresponding neighbourhoods in R1
% // If ANY of them have a pixel of 1, then we know
% // these locations from R2 need to go back to R1
% // In that case, all you really have to do is take whichever
% // pixels are true, and set them to true in R1
% // See which neighbourhoods in `R1` have any pixels that are 1
% // given the corresponding region in R2
finalColumns = any(B1(:,windowsHaving1), 1);
% // Grab the indices of these columns
finalColumnsIndex = windowsHaving1(finalColumns);
% // Reasign those pixels that are neighbours to R1
B1(5, finalColumnsIndex) = 1;
B2(5, finalColumnsIndex) = 0;
%// Restructure back - Recall this image is padded
%// We took the rows and columns of the padded image
%// and so we need to subtract each dimension by 2
R1 = reshape(B1(5,:), rows-2, cols-2);
R2 = reshape(B2(5,:), rows-2, cols-2);
(4) Repeat step 3 until no pixels are reassigned.
Already performed.
(5) Reassign any pixels left in region R2 to R3.
That's pretty straight forward. Any pixels that are true in R2, simply turn these on in R3.
R3 = R2 | R3;
As such, your final double thresholded images are left in R1 and R3. If you want to create an integer map that tells you which pixels are in what region, you can do something like:
map = zeros(size(im));
map(R1) = 1;
map(R3) = 2;
Any pixels that are 0 did not meet the criteria after you went through the double threshold process. As a test, let's take a look at this on a regular image. I'll be using cameraman.tif that is in the MATLAB path. As such, load in your image like so:
im = im2double(imread('cameraman.tif'));
After, run through my code. These are the results I get. White pixels belong to the corresponding region, while black pixels don't.
Good luck!
