Wednesday, April 1st, 2026

Last class, we started to discuss salt and pepper noise.

import numpy as np
import matplotlib.pyplot as plt

Project 4: Image denoising

We looked at loading in image files as 2D numpy arrays storing grayscale values. For whichever image file we want to work with, we'll read it in and convert it (if necessary) to grayscale values scaled between 0 and 1.

rgb_img = plt.imread('everton.png')[:,:,:3]    # Load in the image file and remove any transparency channel information
img = np.mean(rgb_img, axis=2)                 # Average the RGB values together for each pixel to obtain grayscale values

if img.dtype == int:                           # Check if the image contains integer data
    img = img / 255                            # If so, convert to float data

plt.imshow(img, cmap='gray', vmin=0, vmax=1)

<matplotlib.image.AxesImage at 0x15f428f0ec0>

Adding salt and pepper noise

We started too look at using the np.random.random function to select which pixels will be hit by salt and pepper noise.

nrows, ncols = img.shape

noisy_img = img.copy()

salt_array = np.random.random((nrows, ncols))
pepper_array = np.random.random((nrows, ncols))

for row in range(nrows):
    for col in range(ncols):
        if salt_array[row, col] > .9:
            noisy_img[row, col] = 1
        if pepper_array[row, col] > .9:
            noisy_img[row, col] = 0

plt.imshow(noisy_img, cmap='gray', vmin=0, vmax=1)

<matplotlib.image.AxesImage at 0x15f44862210>

Note: the implementation above gives a preference to pepper noise, since it's possible for salt pixels to be overwritten as pepper (but not the reverse). Instead, let's use a single random array, but select select different subintervals between 0 and 1 to identify salt/pepper noise:

nrows, ncols = img.shape

noisy_img = img.copy()

random_array = np.random.random((nrows, ncols))

for row in range(nrows):
    for col in range(ncols):
        if random_array[row, col] > .9:
            noisy_img[row, col] = 1
        if random_array[row, col] < .1:
            noisy_img[row, col] = 0

plt.imshow(noisy_img, cmap='gray', vmin=0, vmax=1)

<matplotlib.image.AxesImage at 0x15f4400ca50>

Note: In the example code above, we iterated through every row and column and decided for each entry whether or not to add salt and pepper noise. This is perfect case to use Boolean arrays.

Idea: create a salt_mask Boolean array containing True in each entry that we wish to hit with salt noise. That is, salt_mask should contain True everywhere that the random_array is less than, say, 0.1.

random_array = np.random.random((8,8))

bool_mask = random_array > .5

bool_mask

array([[ True, False,  True, False, False, False,  True,  True],
       [ True, False,  True,  True,  True, False,  True, False],
       [False,  True,  True, False, False,  True, False, False],
       [False, False, False,  True, False,  True, False, False],
       [False,  True,  True,  True,  True,  True,  True, False],
       [False, False,  True, False,  True,  True,  True,  True],
       [ True,  True, False, False,  True, False,  True, False],
       [False,  True,  True, False,  True,  True,  True, False]])

noisy_img[salt_mask] = 1
noisy_img[pepp_mask] = 0

Exercise. Use the above sample code to write a function sp_noise(img, noise) that adds salt and pepper noise to an image. Its first argument img should be a 2-dimensional numpy array representing the image and the second argument noise should be the fraction of pixels that are to be replaced by noise (for example, with noise = 0.05 about 5% of pixels should be noise, consisting in roughly equal parts of white and black pixels). The function should return a 2-dimensional numpy array representing the original image with noise added.

The mean and median filters

A description of the mean filter is provided on the project page. Let's try to apply the mean filter.

• For right now, we will use a 3 by 3 grid centered centered at each pixel to compute the mean.
• For right now, we will ignore all pixels at the edges of the image.

For now, I will use the pre-noised image file that was downloaded from the course webpage. If you have a function sp_noise function, you should generate your own noisy image to use.

noisy_img = np.mean(plt.imread('noisy_img.png')[:,:,:3], axis=2)
plt.imshow(noisy_img, cmap='gray', vmin=0, vmax=1)

<matplotlib.image.AxesImage at 0x15f4406ead0>

Again, let's make a copy filtered_img of the noisy array that we will filter so that we do not modify the noisy array. This way, we will be able to compare how our filtered image looks with how the noisy image looked.

filtered_img = noisy_img.copy()

If we are considering the pixel in the ith row and jth column (i.e. noisy_img[i,j], we want a 3 by 3 slice centered on the [i,j]th pixel. That is, we want a slice that includes the i-1, i, and i+1 rows and the j-1, j, and j+1 columns.

i = 50
j = 75

grid = noisy_img[i-1:i+2, j-1:j+2]

filtered_img[i,j] = np.mean(grid)

plt.imshow(grid, cmap='gray', vmin=0, vmax=1)

<matplotlib.image.AxesImage at 0x15f455da5d0>

Now let's go pixel by pixel and take the mean of each 3 by 3 grid to create our filtered image array. As a reminder, for now we will ignore any problematic pixels along the edge of the image.

nrows, ncols = noisy_img.shape

filtered_img = noisy_img.copy()
for i in range(1,nrows-1):
    for j in range(1,ncols-1):
        grid = noisy_img[i-1:i+2, j-1:j+2]
        filtered_img[i,j] = np.mean(grid)

plt.imshow(filtered_img, cmap='gray', vmin=0, vmax=1)

<matplotlib.image.AxesImage at 0x15f4547d310>

What about the median filter? The np.median function will compute the median of an array.

nrows, ncols = noisy_img.shape

filtered_img = noisy_img.copy()
for i in range(1,nrows-1):
    for j in range(1,ncols-1):
        grid = noisy_img[i-1:i+2, j-1:j+2]
        filtered_img[i,j] = np.median(grid)

plt.imshow(filtered_img, cmap='gray', vmin=0, vmax=1)
plt.axis('off')

(np.float64(-0.5), np.float64(699.5), np.float64(494.5), np.float64(-0.5))

Exercise: Use the code above to write a function simple_mean_filter that applies the mean filter to non-edge pixels using a 3 by 3 grid.

Exercise: Modify the simple_mean_filter function to write a simple_median_filter that instead applies the median filter.

Note: We could combine these two functions into a single simple_filter function that takes in a noisy image array (noisy_img) along with some function (filter_func) that will be applied to each grid.

The filter_func input should be a function that takes in an 2D array and returns a float (e.g. np.mean or np.median).

def simple_filter(noisy_img, filter_func=np.median):
    ...
    ...
    for ...
        for ...
            grid = ...
            filtered_img[i,j] = filter_func(grid)

    return filtered_img

filtered_img = simple_filter(noisy_img, np.mean)