Steve on Image Processing

Concepts, algorithms & MATLAB

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.



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.

     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.



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)].


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);

imshow(horizontal_gradient, []), title('Horizontal gradient')

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


Comments are closed.

53 CommentsOldest to Newest

Steve replied on : 2 of 53

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.

Alex replied on : 3 of 53

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

Anubhav replied on : 4 of 53

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

Steve replied on : 5 of 53

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

Steve replied on : 7 of 53

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.

Michel replied on : 8 of 53

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.

Steve replied on : 9 of 53

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/.

Suresh replied on : 10 of 53

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.


Steve replied on : 11 of 53

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);
Suresh replied on : 12 of 53

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

Steve replied on : 13 of 53

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)
alice replied on : 14 of 53

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.


Could you please advise how to resolve this problem?

Thanks a lot!

alice replied on : 15 of 53

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.


dewi replied on : 16 of 53

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..

Steve replied on : 17 of 53

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

Suresh replied on : 19 of 53


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?


Steve replied on : 20 of 53

Suresh—See the doc for strel. DEG “specifies the angle (in degrees) of the line, as measured in a counterclockwise direction from the horizontal axis.”

Suresh replied on : 21 of 53


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

Steve replied on : 22 of 53

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.

suhas patel replied on : 25 of 53

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 ???

Steve replied on : 26 of 53

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

N replied on : 27 of 53

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.


Steve replied on : 28 of 53

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);
Student replied on : 29 of 53

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?


praveen replied on : 31 of 53


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.
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.
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.
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.

praveen replied on : 32 of 53

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

Steve replied on : 33 of 53

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.

praveen replied on : 34 of 53

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

Steve replied on : 35 of 53

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.

hira replied on : 36 of 53

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

YCbCr= rgb2ycbcr(im);

% dy= dilate(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.

Steve replied on : 37 of 53

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

hashem replied on : 38 of 53

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

fer replied on : 39 of 53

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?

Steve replied on : 40 of 53

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.

Steve replied on : 43 of 53

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

Ricardo replied on : 46 of 53

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.

Steve replied on : 47 of 53

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

Steve replied on : 49 of 53

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

tom replied on : 51 of 53

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

Steve replied on : 52 of 53

Tom—Image Processing Toolbox functions do not support spatially varying structuring elements.

prajakta replied on : 53 of 53

sir, i need to segment circle, square, thin line, thick line from pcb, i did gradient imdilate and erode but not the the desired output,
i did rgb2gray, thresholding, binarization, then bwmorph with remove, then gradient of dialtion and erosion.
is there any other flow to obtain the result.?
please help