In this post we will experiment with the contents of the section on Filling Holes and Finding Peaks Using Erosion and Dilation of the PluralSight course Building Image Processing Applications Using scikit-image by Janani Ravi.
It seems that this as most of the other sections stand on its own. The idea is to find out that such a library exists, what it is used for, experiment with it, and if the need arises in the future, you will have a basic understanding of how to use the feature.
Like I have mentioned in a previous post, I am also reading a book on NLP. In the book the chapters are related and to some extent connected. It seems that you rely on what you learned in one chapter to enhance or continue in a future one. In this course we are given a broad range of topics which do not seem to be connected to one another. That said, I believe the course is worthwhile and am looking for a next course in the digital image processing subject. That will occur as soon as I complete all the exercises in the NLP book.
Enough of chatting. Let’s get into the main subject of this post.
# **** folder of interest **** cd C:\Documents\_Image Processing\scikit-image-building-image-processing-applications\02\demos # **** open file of interest using VSCode **** (base) C:\Documents\_Image Processing\scikit-image-building-image-processing-applications\02\demos>code ErosionDilationToFillHolesAndFindPeaks.py # **** execute python script of interest **** (base) C:\Documents\_Image Processing\scikit-image-building-image-processing-applications\02\demos>python ErosionDilationToFillHolesAndFindPeaks.py
As in the previous post, we start by getting to our folder of interest. In your case the path will be different. Once there, we open the VSCode IDE with the name of a script that at the time does not exist. Once VSCode opens, please save the blank script. To verify all is well so far, you can run the blank Python script. As we progress in the post, we will be adding code, saving it to the file of interest, and running it.
# **** imports **** import matplotlib.pyplot as plt import numpy as np from skimage import io from skimage import data from skimage.exposure import rescale_intensity from skimage.morphology import reconstruction
As usual we start the script adding the names of the libraries to import. In my opinion this is easier said than done. With experience, one might be able to do so. In general one depends on the IDE to suggest the imports based on the functions being used. The first few lines seem routine. The last two do not. They are associated with morphological operations of erosion and dilation which is the main subject of this post.
# **** read image ****
moon = io.imread('./images/moon.jpg')
# **** display moon image ****
plt.figure(figsize=(8, 8))
plt.imshow(moon, cmap='gray')
plt.title('Moon')
plt.show()
For this post we will be using an image of the moon. We read the image and display it.
# **** let's work with a higher constrast image ****
moon_rescaled = rescale_intensity( moon,
in_range=(30, 200))
# **** display moon_rescaled image ****
plt.figure(figsize=(8, 8))
plt.imshow(moon_rescaled, cmap='gray')
plt.title('Moon Rescaled')
plt.show()
We generate a higher contrast image of the moon. The image is displayed to see the difference between images.
# **** create a seed
# the seed image used for erosion is initialized
# to the maximum value of the original image ****
erosion_seed = np.copy(moon_rescaled)
# **** except at the borders which is the
# starting point of the erosion process ****
erosion_seed[1:-1, 1:-1] = moon_rescaled.max()
# **** create mask ****
mask = moon_rescaled
# **** apply erosion
# eroding inwards from the border fills holes - as the holes are
# by definition surrounded by pixels of higher intensities ****
filled = reconstruction(erosion_seed,
mask,
method='erosion')
# **** display side by side the images ****
fig, ax = plt.subplots( nrows=1,
ncols=2,
figsize=(14, 12),
sharex=True,
sharey=True)
# **** flatten the axes ****
ax = ax.ravel()
# **** display the images ****
ax[0].imshow(moon_rescaled, cmap='gray')
ax[0].set_title('Original Image')
ax[0].axis('off')
ax[1].imshow(filled, cmap='gray')
ax[1].set_title('Erosion Reconstruction')
ax[1].axis('off')
plt.tight_layout()
plt.show()
We create a seed image which we will use to experiment with the erosion process. The result is a mask that we will use in conjunction with the erosion_seed. The images are then displayed side by side for comparison.
# **** create a seed for dilation ****
dilation_seed = np.copy(moon_rescaled)
# **** the seed image used for dilation is initialized
# to the minimum value of the original image
# except at the borders which is the
# starting point of the dilation process
dilation_seed[1:-1, 1:-1] = moon_rescaled.min()
# **** mask ****
mask = moon_rescaled
# **** apply dilation
# highlights bright spots in an image by
# expanding maximal values in local areas ****
highlighted = reconstruction( dilation_seed,
mask,
method='dilation')
# **** display side by side the images ****
fig, ax = plt.subplots( nrows=1,
ncols=2,
figsize=(14, 12),
sharex=True,
sharey=True)
# **** flatten the axes ****
ax = ax.ravel()
# **** display the images ****
ax[0].imshow(moon_rescaled, cmap='gray')
ax[0].set_title('Original Image')
ax[0].axis('off')
ax[1].imshow(highlighted, cmap='gray')
ax[1].set_title('Dilation Reconstruction')
ax[1].axis('off')
plt.tight_layout()
plt.show()
Now we will repeat the previous steps but this time we will make the necessary changes to apply dilation. As before the resulting images are displayed side by side for comparison.
# **** isolate dark regions ~ erode - original ****
holes = filled - moon_rescaled
peaks = highlighted - moon_rescaled
# **** OR issolate bright regions ~ dilate - original ****
# bright = moon_rescaled - highlighted
# peaks = highlighted - moon_rescaled
# **** display side by side the images ****
fig, ax = plt.subplots( nrows=1,
ncols=2,
figsize=(14, 12),
sharex=True,
sharey=True)
# **** flatten the axes ****
ax = ax.ravel()
# **** display the images ****
ax[0].imshow(holes, cmap='gray')
ax[0].set_title('Holes')
ax[0].axis('off')
ax[1].imshow(peaks, cmap='gray')
ax[1].set_title('Peaks')
ax[1].axis('off')
plt.tight_layout()
plt.show()
As you can see we can isolate dark or bright regions in the image. The current code illustrates isolation of dark regions. The commented out lines that follow illustrate the isolation of the bright regions.
Finally we display side by side the two images.
I believe the topic in this post is quite interesting and adds to your toolkit of algorithms a different way to approach different types of images. I wish I could have more time to select a couple of different images and experiment with them. If you do have the time I believe it would be well spent.
If interested in the code for this post you can find it in my GitHub repository ErosionDilationToFillHolesAndFindPeaks.
Keep in mind that one of the best ways to learn is to read, experiment, and then repeat as needed.
Enjoy,
John