Steve on Image Processing

Dilation, erosion, and the morphological gradient 53

Posted by Steve Eddins,

The morphological operator dilation acts like a local maximum operator. Erosion acts like a local minimum operator. You can use them together to compute something called the morphological gradient.

Contents

Dilation

The basic form of grayscale image dilation computes, for each image pixel, the maximum value of its neighboring pixels. The neighborhood is defined by the structuring element. For example, this structuring element:

se1 = strel([1 1 1])
 
se1 =
 
Flat STREL object containing 3 neighbors.

Neighborhood:
     1     1     1

 

defines a neighborhood consisting of the pixel itself, together with its left and right neighbors. Whereas this structuring element:

se2 = strel([1; 1; 1])
 
se2 =
 
Flat STREL object containing 3 neighbors.

Neighborhood:
     1
     1
     1

 

defines a neighborhood consisting of the pixel itself, together with the pixels above and below it.

Let's try dilation on everyone's favorite sample MATLAB matrix, the magic square:

m5 = magic(5)
m5 =

    17    24     1     8    15
    23     5     7    14    16
     4     6    13    20    22
    10    12    19    21     3
    11    18    25     2     9

d1 = imdilate(m5, se1)
d1 =

    24    24    24    15    15
    23    23    14    16    16
     6    13    20    22    22
    12    19    21    21    21
    18    25    25    25     9

d1(3,3), which is 20, is the maximum of [m5(3,2), m5(3,3), m5(3,4)], or [6 13 20].

Dilating with se2 computes the 3-pixel local maximum vertically:

d2 = imdilate(m5, se2)
d2 =

    23    24     7    14    16
    23    24    13    20    22
    23    12    19    21    22
    11    18    25    21    22
    11    18    25    21     9

d2(3,3) is the maximum of [m5(2,3), m5(3,3), m5(4,3)].

Erosion

Grayscale image erosion computes the minimum of each pixel's neighborhood.

e1 = imerode(m5, se1)
e1 =

    17     1     1     1     8
     5     5     5     7    14
     4     4     6    13    20
    10    10    12     3     3
    11    11     2     2     2

e1(3,3) is the minimum of [m5(3,2), m5(3,3), m5(3,4)].

Morphological gradient

Dilation and erosion are often used in combination to produce a desired image processing effect. One simple combination is the morphological gradient. P. Soille, in section 3.8 of the second edition of Morphological Image Analysis: Principles and Applications, talks about three kinds of basic morphological gradients:

  • dilated_image - eroded_image
  • original_image - eroded_image
  • dilated_image - original_image

Soille calls the first one the basic morphological gradient, and you compute it this way using MATLAB and the Image Processing Toolbox:

% Read in circuit board image, crop it so it's not too big for the blog
% page, and convert it to grayscale:

rgb = imread('board.tif');
I = rgb2gray(rgb(1:256,1:256,:));

se = strel(ones(3,3));
basic_gradient = imdilate(I, se) - imerode(I, se);

subplot(1,2,1), imshow(I), imcredit('Image courtesy of Alexander V. Panasyuk')

subplot(1,2,2), imshow(basic_gradient, [])

The second form is called the half-gradient by erosion or internal gradient.

internal_gradient = I - imerode(I, se);

subplot(1,2,2), imshow(internal_gradient, [])

"The internal gradient enhances internal boundaries of objects brighter than their background and external boundaries of objects darker than their background. For binary images, the internal gradient generates a mask of the internal boundaries of the foreground image objects." [Soille, page 86]

The third form is called the half-gradient by dilation or external gradient:

external_gradient = imdilate(I, se) - I;

Direction gradients

By using line segments as structuring elements, you can compute directional gradients.

seh = strel([1 1 1]);
sev = strel([1;1;1]);

horizontal_gradient = imdilate(I,seh) - imerode(I,seh);
vertical_gradient = imdilate(I,sev) - imerode(I,sev);

subplot(1,2,1)
imshow(horizontal_gradient, []), title('Horizontal gradient')

subplot(1,2,2)
imshow(vertical_gradient, []), title('Vertical gradient')

Try this on your own grayscale images. You can experiment with larger structuring elements, too, such as strel('disk',7).

I'll be covering ways to combine erosion, dilation, and other fundamental morphological operators in the near future.


Get the MATLAB code

Published with MATLAB® 7.3

53 CommentsOldest to Newest

Vihang – Are you talking about what regionprops does? For the ConvexHull option of regionprops, it first computes the perimeter pixels. Each perimeter has four corner points. regionprops calls convhull to compute the convex hull of the set of corner points of all the perimeter pixels. Does that help? See regionprops documentation for more information.

Would you please explain dilation and erosion with non-flat structuring elements on grayscale images? Thanks!

Sir i am working on plotting contours of back side of fingers and then find the length and width of finger. how should i find length and width of finger using matlab? i have tried doing it but problem is that image also contains some part of hand that is beneath fingers, so we have to make a reference point and thn only we can find the length. Can u please help me with it

Anubhav—Sounds like an interesting algorithm development problem, although maybe a little vague. But I don’t have any particular suggestions for you.

Pramod—Oversimplifying a bit, dilation is a local max operator on some defined neighborhood of every pixel, and erosion is a local min operator. Sometimes there is also an additive offset applied to each pixel in the neighborhood before the max or min is computed. Details are plentiful in many image processing textbooks, as well as online.

We are processing holograms using your toolbox in our lab. Since holograms contain 3D information, it makes sense to do the image processing in 3D as well. The problem is that the 3D images matrices are huge, that the operations like imdilate and imerode take a long time and that it puts enormous demands on computer memory. The matrices contain however mainly zeros, so we decided to write some codes to do dilation and erosion based on the (x,y,z)-coordinates of the pixels. Our dilation code (fully vectorized) is faster than imdilate for small 3D structure elements, but is quickly outperformed by imdilate for larger ones. I’m now working on the (coordinate-based) erosion code, but I’m a bit clueless how to tackle it in preferably vectorized manner.
Although the concept of dilation and erosion is simple, I’m puzzled as to how MATLAB does this so efficiently. Could you please explain how MATLAB performs dilation and erosion? If you have any ideas about performing these tasks on pixel coordinates only; I’m all ears. Many thanks in advance for your answer.

Michel—I’ve had dilation and erosion algorithms on my blog topics list for a long time, but I’ve been procrastinating because it’s a complex story. There are several different algorithms used in a complementary way by imdilate and imerode: Structuring element decomposition, binary word packing, and a special method just for horizontal and vertical line structuring elements.

In your comparisons, are you using a rectangular structuring element? For that case, our combination of algorithms achieves a running time that is independent of the structuring element size.

For other shapes, you could be seeing the effects of our other structuring element decompositions.

Note also that there’s a lot C code underlying our morphology implementations. See the source code in toolbox/images/images/private/.

I am trying to use morphological operators for one dimensional signal (e.g. sinusoidal signal with random noise).

For the initial start, I am trying to visualize the process manually using a vector representing a signal (eg input signal, x=[1,1.3,1.6,1.4,1.6] and structuring element, g=[0,1,0]). I used the following expression for dilation.

dil(x)=max{x(n-m)+g(m)}, n-length of IP signal,m-length of SE

I get different result when I calculate it manually and when I use matlab to do so.

Is it possible to know how matlab handles this expression? Or what algorithm matlab uses?

Thanks in advance.

Suresh

Suresh—You didn’t say how you were computing the dilation in MATLAB, and your formula is incomplete, but I have a guess about what your problem might be. Perhaps you are using this syntax for imdilate:

  y = imdilate(x, [0 1 0]);

With this syntax, the 2nd argument specifies the domain of a flat, zero-valued structuring element. Since only the center element of your domain is nonzero, your structuring element has only a single element in it, and the output is the same as the input.

To specify and use nonflat structuring element, you need to do this:

  se = strel('arbitrary', [1 1 1], [0 1 0]);
  y = imdilate(x, se);

Steve,
Thank you for the reply. I appreciate your help.
Your suggestion worked but I still get error in one of the element. Here is what I got using matlab,

>> x=[1,2,4,6,-20,3,1];
>> se = strel(‘arbitrary’, [1 1 1], [0 1 0])
>> y = imdilate(x, se)

y =

2 4 6 7 6 4 3

But this is what I get manually,
y = [1 2 4 6 7 6 4]

The example quoted above (plus manual solution) is taken from a paper [X. Lin, H. Liu, H. Weng, P. Liu, B. Wang, and Z. Q. Bo, Senior Member, IEEE “A Dual-Window Transient Energy Ratio-Based Adaptive Single-Phase Reclosure Criterion for EHV Transmission Line”, IEEE Transactions on Power Delivery, Vol. 22, NO. 4, pp.2080-2086, October 2007].

Also I have some papers, (not quoted here to make shorter) which say that Structuring Elements can be semi-circular, triangular or even sinusoidal with real numbers as its element. But Matlab documentation says SE consists of only 0 and 1. Is there any way to deal with this case in matlab?

Thank you
Suresh

Suresh—It looks to me like your result is simply shifted by one position. Note that the Image Processing Toolbox is for the center of the neighborhood matrix is the origin. Your paper and your manual computations are probably using a different convention.

MATLAB documentation does not say that the structuring element consists of only 0s and 1s. The structuring element neighborhood (or domain) input is 0s and 1s. As I mentioned in my previous response, to create a nonflat structuring element you should use this syntax:

se = strel('arbitrary', nhood, height)

Hello, Steve,

I tried to do binary expansion with bwmorph. The test.m file contains
mask = zeros(12, 12);
mask(5:7, 5:7) = 1;
mask = bwmorph(logical(mask), ‘dilate’);

It works well in matlab script. There was also no error when I compiled test.m with mcc -m -B sgl -T link:exe test.m. However, executing !test.exe yields
“Warning: Duplicate directory name: c:\matlab6p5\toolbox\matlab\iofun.
ERROR: Reference to unknown function ‘checknargin’ from FEVAL in stand-alone mode.

EXITING ”

Could you please advise how to resolve this problem?

Thanks a lot!
Alice

Thanks, Steve

Yes, it was the compiler’s problem. I have tried the solution on Matlab v7.1 R14 and it works!

Thanks for your help.

Alice

Sir, coulf you tell me about how erotion and dilation work in segmentation of an image??

I’m making a face detection and now I’ve segmented the RGB image using the intensity. I want to know what method you would like to recommended that gives optimum automated threshold value? is it using dilation/erotion, using Otsu,or something else?

Thank you so much..

Dewi—The function graythresh, which uses Otsu’s method, does a pretty good job of automated threshold selection.

Steve,

Could you explain what the “DEG” in flat line SE mean in the following function?

SE = STREL(‘line’,LEN,DEG)

How does the angle affect the dilation or erosion?

Thanks,
Suresh

Steve,

Thank you for the reply.

I am using dilation/erosion for filtering noise in sinusoidal signal. A flat line SE with DEG=0 works for some signals but for others DEG=90 has to be used.

I could not figure out what effect does this angle has on dilation/erosion. Does the fundamental definition of ‘local minima’ and ‘local maxima’ for erosion and dilation apply in case where the angle is not 0 degree?

I could find algorithms for case where DEG=0 but not for other angles. Could you suggest me any book/paper or website, where I could find some help?

Thanks in advance
Suresh

Suresh—The shape of the structuring element controls which neighborhood pixels participate in the local min / local max operation. For example, dilating by a 1-by-5 structuring element (line at 0 degrees) is finding the local max of a pixel and its 4 nearest horizontal neighbors. Dilating by a 5-by-1 structuring element (line at 90 degrees) is finding the local max of a pixel and its 5 nearest vertical neighbors. Most image processing textbooks cover morphological operations. You could try Digital Image Processing by Gonzalez and Woods, or Digital Image Processing Using MATLAB by Gonzalez, Woods, and Eddins.

i want to know use of morphological operator to find convexhull…
in matlab convex hull is find out using qhull technique.
so is it possible to find out convex hull using morphological operator ???

Suhas—Try the P. Soille book, Morphological Image Analysis. I believe it contains a description of computing the convex hull using morphology.

Hi Steve,

I would like to use imerode on a 3D stack of segmentations (binary masks) of medical images of the brain. I see in some of your other posts that imerode should work with 3D data. Can I use strel to make a 3D structuring element? If so, can you give an example? I don’t see mention of anything but 2d images in the help documentation.

Also, can you think of a way to handle the fact that my in-plane resolution is different than the slice thickness? For example in my case the dimensionality of each voxel is ~1x1x3 mm.

My problem is to identify voxels that are less than a certain distance from the brain surface.

Thanks

N—Sure, you can use imerode to do 3-D erosion. This call makes a 3-D structuring element:

se = strel(ones(3, 3, 3));

You can handle the spatial dimensions yourself by adjusting the size of the structuring element in the third dimension. For example, you could make a 9-by-9-by-3 structuring element:

se = strel(9, 9, 3);

Dear Steve,

My teacher is currently giving your document as if it is his for a lecture in my class. He doesn’t really understand it all that well. Do you think you can explain in a bit more detail for me?

Thanks,
Student

steve,

i have problem in choosing flat structuring element to perform morphological opening on a sinusoidal signal(1-D signal).
i will expalain my problem clearly.
my main aim is to decompose the sinusoidal signal in to three small segments.
in my problem, i have 3 levels to decompose signal in to
3 small segments.
first i sampled the sinusoidal signal and collected the sampled values in to an array.

at level-1: if A is original signal, i will apply morphological opening on A by using flat structuring element of length one sample.
R1=opening(A);
S1=A-R1; where S1 is the first segment.

at level 2: i will apply morphological opening on R1 by using flat structuring element of length two samples.
R2=opening(R1);
S2=R1-R2; where S2 is the second segment.

at level 3: i will apply morphological opening on R2 by using flat structuring element of length three samples.
R3=opening(R2);
S3=R3-R2; where S3 is the third segment.

S1,S2,S3 are the my needed segments.
if i am using structuring element of length one sample at level-1, i am getting opening output is as same as input.
as a result S1 is becoming zero (it is not correct for my application).

for above my application, what type of structuring element is suitable by considering above conditions.
please kindly reply me.
thank you.

steve
in the above comment at level -3: S3=R2-R3 is the correct,
previously mentioned S3= R3-R2 is wrong.

Praveen—Sounds to me like the difficulty is in the problem definition itself. Flat dilation with a single-element structuring element (se = 1) doesn’t change the image at all. It follows directly from the definition of dilation that the output image would equal the input image.

ok steve,
can you suggest me which type of structuring element is suitable for my application to decompose signal in to segments with out using structuring element of length one sample.
thank you

Praveen—I find the problem statement too vague. Without more information about the purpose of the decomposition, or what properties the decomposed segments are supposed to have, I have no idea what to suggest.

hi steve

i want to dilate an image to propagate the high values and erode for the low ones. could you please suggest me how to do it. the image is the Y component of YCbCr. after performing dilation and erosion i have to divide the images some thing like this

im=imread('image058.jpg');
YCbCr= rgb2ycbcr(im);
y=YCbCr(:,:,1);

% dy= dilate(y);
%ey=erode(y);
% map= dy/ey;

the division gives error
please guide me how i could perform the dialtion and erosion and how to divide both images.

Hira—Use imdilate and imerode to dilate and erode your image. Use the ./ operator instead of / to do element-wise division.

Hi
At first I convert RGB image in uint8 class to double class with im2double function.
then I apply erosion to an image with ball structure, the resulting pixels have negative values. how can I convert them to positive equal?
Thanks alot

Sorry, I’m not sure this is the right place to put this but…………
I’m trying to use “imreconstruct” command but I keep on getting the same error message:

“???MARKER pixels must be <= MASK pixels”

could naybody help me out?
Thanks!!

Fer—With the sparse information you have provided, all I can say is that you must have marker pixels that are greater than your mask pixels. :-) Perhaps the imreconstruct doc will help.

Danish—The Image Processing Toolbox supports morphology on arrays (rectlinear grids) only. You’ll need to develop your own implementations for other geometries.

Hello Steve. I am currently working on computing the gradient of a synthetic image to which noise has been added. By computing the normal gradient, when I calculate the gaussian curvature of the noisy image, the results are really bad.

So I have been trying to compute the morphological gradient instead, that has been giving me better results. However, I’ve noticed that when computing the directional gradients, by using larger structuring elements, the result is a lot more resistant to noise.

My question is, is there a disadvantage to using large structuring elements, and why is it more robust to noise? Thank you.

Ricardo—Large structuring elements take longer to compute and can cause detailed structures of interest (as opposed to variations caused by noise) to vanish.

Tom—Sure. imdilate and imerode support arbitrarily shaped, multidimensional structuring elements. See the doc for strel.

i have another doubt,is it possible to create a spatially variant structuring element by knowing the orientation of the pixels in an image

These postings are the author's and don't necessarily represent the opinions of MathWorks.