Pioneer21/Lecture03

This page serves as a supplement to the third group meeting of the Summer 2021 Pioneer program for image processing.

Python Requirements

The programs from this lecture will require several modules:

• NumPy - https://numpy.org/
• Matplotlib - https://matplotlib.org/ (and specifically matplotlib.pyplot)
• scikit-image - https://scikit-image.org/
• SciPy - https://www.scipy.org/ (and specifically scipy.signal)

The start of any script we use in this lecture will be:

import numpy as np
import matplotlib.pyplot as plt
import skimage as ski
import scipy.signal as sig

scikit-image commands

The main command for this session is:

• DATA = ski.data.NAME() where DATA is a variable name you give to store the information from a built-in image and NAME is the name of a built-in image from scikit-image; the built-in images can be found in the Data section of the General Examples page.

Starting on July 26, 2021, we will be using AXIS.imshow(DATA) instead of ski.io.imshow(DATA) to display an image. It looks like ski.io.imshow is going away.

Depending on the nature of an image, its DATA variable will be different shapes and contain different data types. Here is how scikit-image and imshow() interpret things:

• A 2D-array will either be black and white, grayscale, or it will be mapped to a colormap
• A 3D-array will be interpreted as an RGB-based color image
• For a 2D-array, the following data types and sizes might be useful:
  • An array of boolean (True or False) values or integer values with only integers 0 or 1 will produce an image with values at the extreme edges of the colormap. Assume that the code:

    bwimage = ski.data.binary_blobs()
    

    has been run already for each of the cases below. By default in Python with matplotlib this is purple and gold.
    • Default case with booleans

      fig, ax = plt.subplots(num=1, clear=True)
      ax.imshow(bwimage)
      

      You can add some keyword arguments to the command to get black and white as follows:
    • Black and white with booleans

      fig, ax = plt.subplots(num=1, clear=True)
      ax.imshow(bwimage, cmap=plt.cm.gray)
      

    • Black and white and no axis ticks with booleans

      fig, ax = plt.subplots(num=1, clear=True)
      ax.imshow(bwimage, cmap=plt.cm.gray)
      ax.set_axis_off()
      

      Note that for each of the above, you can replace the array of booleans with an array of integers 0 or 1 by replacing the bwimage line with:

      bwimage = np.random.choice([0,1], (50,100))
      

      or

      bwimage = np.random.randint(0, 2, (50, 100))
      

  • An array of integer values will produce an image with values mapped to the default Viridis (purple - green - blue - yellow) colormap. Assume that the code:

    grayimage = ski.data.coins()
    

    has been run already for each of the cases below.
    • Default case with integers

      fig, ax = plt.subplots(num=1, clear=True)
      ax.imshow(grayimage)
      

      You can add some keyword arguments to the command to get shades of gray as above:
    • Gray with integers

      fig, ax = plt.subplots(num=1, clear=True)
      ax.imshow(grayimage, cmap=plt.cm.gray)
      

    • Gray and no axis ticks with integers

      fig, ax = plt.subplots(num=1, clear=True)
      ax.imshow(grayimage, cmap=plt.cm.gray)
      ax.set_axis_off()
      

• Note 1 - the addition of the colormap is required for using the axis version of imshow versus the scikit-image version of imshow.
• Note 2 - by default, imshow will stretch the range of integers to cover the whole map. If you want to specifically map to a certain range - for instance, have 0 as pure black and 255 as pure white, you need to add vmin and vmax keyword arguments to explicitly give the values that map to the low and high end of the volormap. For instance:

  fig, ax = plt.subplots(num=1, clear=True)
  grayimage = ski.data.coins()
  ax.imshow(grayimage, cmap=plt.cm.gray, vmin=0, vmax=255))
  

  will produce the exact same image that ski.io.imshow(grayimage) would have produced.

• • An array of floating-point values will also produce an image with values mapped to the default Viridis (purple - green - blue - yellow) colormap. You can add some keyword arguments to the command to get shades of gray as above. Assume that the code:

    x, y = np.meshgrid(np.linspace(-2, 2, 201), np.linspace(-1, 1, 101))
    z = np.exp(-np.sqrt(x**2+y**2))*np.cos(2*np.pi*x)*np.sin(4*np.pi*y)
    

    has been run already for each of the cases below.
    • Default case with floats

      fig, ax = plt.subplots(num=1, clear=True)
      ax.imshow(z)
      

    • Gray with floats

      fig, ax = plt.subplots(num=1, clear=True)
      grayimage = ski.data.coins()
      ax.imshow(grayimage, cmap=plt.cm.gray)
      

    • Gray and no axis ticks with floats

      fig, ax = plt.subplots(num=1, clear=True)
      grayimage = ski.data.coins()
      ax.imshow(grayimage, cmap=plt.cm.gray)
      ax.set_axis_off()
      

• Note - once again, by default, imshow will stretch the range of floats to cover the whole map. If you want to specifically map to a certain range - for instance, have 0 as pure black and 0.2 as pure white, you need to add vmin and vmax keyword arguments to explicitly give the values that map to the low and high end of the colormap. Anything below vmin will be mapped to the low color and anything above vmax will be mapped to the high. For instance:

  fig, ax = plt.subplots(num=1, clear=True)
  ax.imshow(z, vmin=0, vmax=255))
  

  will produce an image with yellow rounded rectangles on a purple background with greenish-blue edges the further you get away from the center. The large all-yellow rectangles contain many values that are above 0.2 and the purple parts are all the negative parts of the array.

• For a 3D array, the third dimension should have three layers representing the red, green, and blue levels for each pixel. The entries either need to be integers between 0 and 255 or floating point numbers between 0 and 1. If the array contains values outside of those limits, strange things happen.

Links

• 2D Convolution information: HTTP://www.songho.ca/dsp/convolution/convolution2d_example.HTML
• Kernels: https://en.wikipedia.org/wiki/Kernel_(image_processing)
• Test Card for later use: https://en.wikipedia.org/wiki/Test_card

Examples

The following sections will contain examples of different images. You may want to run them in your own version of Python to make sure you are getting the same answers. These are provided so you can compare what you get with what we think you should get.

Example 1: Black & White Images

a = np.array([
    [1, 0, 1, 0, 0],
    [1, 0, 1, 0, 1],
    [1, 1, 1, 0, 0],
    [1, 0, 1, 0, 1],
    [1, 0, 1, 0, 1]])
fig, ax = plt.subplots(num=1, clear=True)
ax.imshow(a, cmap=plt.cm.gray)
ax.set_axis_off()
fig.tight_layout()
fig, ax = plt.subplots(num=2, clear=True)
aplot2 = ax.imshow(a, cmap=plt.cm.gray)
fig.colorbar(aplot2)
fig.tight_layout()

Example 2: Simple Grayscale Images

bx, by = np.meshgrid(range(255), range(255))
fig, ax = plt.subplots(num=1, clear=True)
ax.imshow(bx, cmap=plt.cm.gray)
fig.tight_layout()
fig, ax = plt.subplots(num=2, clear=True)
ax.imshow(by, cmap=plt.cm.gray)
fig.tight_layout()

Example 3: Less Simple Grayscale Images

x, y = np.meshgrid(np.linspace(0, 2*np.pi, 201),
                   np.linspace(0, 2*np.pi, 201));
z = np.cos(x)*np.cos(2*y);
fig, ax = plt.subplots(num=1, clear=True)
zplot = ax.imshow(z, cmap=plt.cm.gray)
ax.axis('equal')
fig.colorbar(zplot)
fig.tight_layout()

Notice how the use of ax.axis('equal') will make the image look like a square since it is 201x201 but also caused the display to be filled with some whitespace as a result of the figure size versus the image size. Also note that the values from -1 to 1 were automatically mapped to the full colormap.

Example 4: Building an Image

rad = 100
delta = 10
x, y =np.meshgrid(range((-3*rad-delta),(3*rad+delta)+1),
                  range((-3*rad-delta),(3*rad+delta)+1))
rows, cols = x.shape
dist = lambda x, y, xc, yc: np.sqrt((x-xc)**2+(y-yc)**2)
venn_img = np.zeros((rows, cols, 3));
venn_img[:,:,0] = (dist(x, y, rad*np.cos(0), rad*np.sin(0)) < 2*rad);
venn_img[:,:,1] = (dist(x, y, rad*np.cos(2*np.pi/3), rad*np.sin(2*np.pi/3)) < 2*rad);
venn_img[:,:,2] = (dist(x, y, rad*np.cos(4*np.pi/3), rad*np.sin(4*np.pi/3)) < 2*rad);
fig, ax = plt.subplots(num=1, clear=True)
ax.imshow(venn_img)
ax.axis('equal')
fig.tight_layout()

Example 5: Exploring Colors

nr = 100
nc = 200
x, y = np.meshgrid(np.linspace(0, 1, nc),
                   np.linspace(0, 1, nr))
other = 0.5;
palette = np.zeros((nr, nc, 3))
palette[:,:,0] = x
palette[:,:,1] = y
palette[:,:,2] = other
fig, ax = plt.subplots(num=1, clear=True)
ax.imshow(palette)
ax.axis('equal')
fig.tight_layout()

Example 6: 11x11 Blurring

x = ski.data.coins()
h = np.ones((11, 11))/11**2;
y = sig.convolve2d(x, h, 'same');
fig, ax = plt.subplots(num=1, clear=True)
ax.imshow(x, cmap=plt.cm.gray)
ax.axis('equal')
ax.set(title='Original')
fig.tight_layout()
fig, ax = plt.subplots(num=2, clear=True)
ax.imshow(y, cmap=plt.cm.gray)
ax.axis('equal')
ax.set(title='Blurred')
fig.tight_layout()

Example 7: Basic Edge Detection and Display

x, y = np.meshgrid(np.linspace(-1, 1, 200),
                   np.linspace(-1, 1, 200));
z1 = (.7<np.sqrt(x**2+y**2)) * (np.sqrt(x**2+y**2)<.9)
z2 = (.3<np.sqrt(x**2+y**2)) * (np.sqrt(x**2+y**2)<.5)
zimg = 100*z1+200*z2;
fig, ax = plt.subplots(num=1, clear=True)
ax.imshow(zimg, cmap=plt.cm.gray)
ax.axis('equal')
ax.set(title='Original')
fig.tight_layout()

hx = np.array([[1, -1]])
edgex = sig.convolve2d(zimg, hx, 'same')
fig, ax = plt.subplots(num=2, clear=True)
ax.imshow(edgex, cmap=plt.cm.gray)
ax.axis('equal')
ax.set(title='Vertical Edges')
fig.tight_layout()

hy = np.array([[1], [-1]])
edgey = sig.convolve2d(zimg, hy, 'same')
fig, ax = plt.subplots(num=3, clear=True)
ax.imshow(edgey, cmap=plt.cm.gray)
ax.axis('equal')
ax.set(title='Horizontal Edges')
fig.tight_layout()

edges = np.sqrt(edgex**2 + edgey**2)
fig, ax = plt.subplots(num=4, clear=True)
ax.imshow(edges, cmap=plt.cm.gray)
ax.axis('equal')
ax.set(title='Edges')
fig.tight_layout()