Image Filtering and Frequency Analysis

This blog explores image filtering techniques, frequency analysis, and multiresolution blending. We'll cover convolution operations, edge detection, hybrid image creation, and advanced blending techniques using Gaussian and Laplacian pyramids.

1.1 Convolution

Here are the two implementations of 2D convolution:

def conv2d_4for(image, kernel): # use odd sized kernel padded_shape = (image.shape[0] + 2 * (kernel.shape[0] // 2), image.shape[1] + 2 * (kernel.shape[1] // 2)) padded_image = np.zeros(padded_shape) padded_image[kernel.shape[0] // 2: -kernel.shape[0] // 2 + 1, kernel.shape[1] // 2: -kernel.shape[1] // 2 + 1] = image output_shape = (padded_image.shape[0] - kernel.shape[0] + 1, padded_image.shape[1] - kernel.shape[1] + 1) output = np.zeros(output_shape) for i in range(output_shape[0]): for j in range(output_shape[1]): for ki in range(kernel.shape[0]): for kj in range(kernel.shape[1]): output[i, j] += padded_image[i+ki, j+kj] * kernel[ki, kj] return output def conv2d_2for(image, kernel): # use odd sized kernel padded_shape = (image.shape[0] + 2 * (kernel.shape[0] // 2), image.shape[1] + 2 * (kernel.shape[1] // 2)) padded_image = np.zeros(padded_shape) padded_image[kernel.shape[0] // 2: -kernel.shape[0] // 2 + 1, kernel.shape[1] // 2: -kernel.shape[1] // 2 + 1] = image output_shape = (padded_image.shape[0] - kernel.shape[0] + 1, padded_image.shape[1] - kernel.shape[1] + 1) output = np.zeros(output_shape) for i in range(output_shape[0]): for j in range(output_shape[1]): output[i, j] = np.sum(padded_image[i:i+kernel.shape[0], j:j+kernel.shape[1]] * kernel) return output

We see that the two implementations as well as the scipy.signal.convolve2d function (when set to 'same' mode) all give the same output, but the speeds are significantly different with scipy taking a few seconds, 2d for loop in the tens of seconds, and 4d for loop in the multiple minutes when ran on test images.

Here are some images that we used this function on:

Box filter result
Box Filter
D_x convolution
D_x Convolution
D_y convolution
D_y Convolution
Grayscale image
Grayscale Image

1.2 Finite Difference Operator

Here is the original cameraman image:

Cameraman original
Original Cameraman Image

Here are the results when we run a convolution with d_x and d_y which are the finite difference operators on the cameraman image:

Cameraman d_x convolution
Cameraman D_x Convolution
Cameraman d_y convolution
Cameraman D_y Convolution

This doesn't show much so instead, we check different possible thresholds for the edges to make a binary plot. The threshold that we found was most useful was about 15. This was chosen as a higher threshold did not show the tower in the background whereas a lower threshold showed little specs in the sky which is shown in the following images:

Cameraman d_x binary
Cameraman D_x Binary (threshold=15)
Cameraman d_y binary
Cameraman D_y Binary (threshold=15)

1.3 Gaussian Filter

Here are the results when we run a convolution with the gaussian filter:

Gaussian convolution
Gaussian Convolution Result

We further check the partial derivatives of the gaussian filter:

Gaussian partial x
Gaussian Partial X
Gaussian partial y
Gaussian Partial Y

We see that the partial derivatives of the gaussian filter are very similar to the finite difference operators, which is expected.

We also verify that this can be done with a single convolution, which is possible since convolutions are associative, using something like the following code:

single_conv_filter_d_x = scipy.signal.convolve2d(two_d_gaussian, d_x, mode='same') single_conv_filter_d_y = scipy.signal.convolve2d(two_d_gaussian, d_y, mode='same') camera_man_gaussian_partial_x_once = scipy.signal.convolve2d(camera_man, single_conv_filter_d_x, mode='same') camera_man_gaussian_partial_y_once = scipy.signal.convolve2d(camera_man, single_conv_filter_d_y, mode='same') # camera_man_gaussian_partial_x is the result after doing two convolution assert np.allclose(camera_man_gaussian_partial_x_once, camera_man_gaussian_partial_x) assert np.allclose(camera_man_gaussian_partial_y_once, camera_man_gaussian_partial_y)

We also show the output of the binary plots here:

Gaussian partial x binary
Gaussian Partial X Binary
Gaussian partial y binary
Gaussian Partial Y Binary

We see that the output is much smoother than the original because the gaussian filter helps reduce the noise in the image.

2. Frequency Analysis

2.1 Getting Low and High Frequencies to Sharpen Images

To get the low and high frequencies, we can use a gaussian filter to get the lower frequencies as the higher frequencies will be lost to the filter and then subtract it from the original to get the remaining frequencies, which will be the high frequencies.

The following shows the above process on the Taj Mahal:

Taj frequency analysis
Taj Mahal Frequency Analysis - Low and High Frequency Separation

We can see the output of various values of alpha:

Taj sharpening
Taj Mahal Sharpening with Different Alpha Values

We can see that the output is much sharper when alpha is larger.

2.2 Hybrid Image Creation

To get a hybrid image, the steps go as follows:

  1. Align the two images
  2. Get the low and high frequencies of the two images
  3. Combine the low frequency of one image with the high frequency of the other image

To get the low and high frequencies, we can use a gaussian filter to get the lower frequencies as the higher frequencies will be lost to the filter and then subtract it from the original to get the remaining frequencies, which will be the high frequencies.

Here is the output of the hybrid image creation with various steps included:

Hybrid FFT analysis
Hybrid Image Creation Process with FFT Analysis

You can see that the FFT of the hybrid is equal to the sum/average of the FFT of the two images.

Here are some other hybrid images created from different source pairs:

Tchaikovsky & Kreisler Hybrid

Original source images:

Tchaikovsky original
Tchaikovsky (Low Frequencies)
Kreisler original
Kreisler (High Frequencies)
Tchaikovsky Kreisler hybrid
Tchaikovsky & Kreisler Hybrid Result

Tiger & Lion Hybrid

Original source images:

Tiger original
Tiger (Low Frequencies)
Lion original
Lion (High Frequencies)
Tiger Lion hybrid
Tiger & Lion Hybrid Result

2.3 & 2.4 Multiresolution Blending

Here is the process of the multiresolution blending using Gaussian and Laplacian stacks:

Laplacian stack analysis
Optimized Laplacian Stack Analysis - Orange, Apple, and Oraple Frequency Decomposition

Mango & Pear Multiresolution Blend

Original source images:

Mango original
Mango (Left Side)
Pear original
Pear (Right Side)
Mango pear blend
Mango and Pear Multiresolution Blending Result

The multiresolution blending technique allows for seamless combination of images by working in different frequency bands. High frequencies are blended sharply while low frequencies transition smoothly, creating natural-looking results without visible seams.