Why the brightness of the image got reduced after applying gaussian filter?

Question

I just learned how to apply a gaussian filter from scratch to a Grayscale image in python from this blog:
http://www.adeveloperdiary.com/data-science/computer-vision/applying-gaussian-smoothing-to-an-image-using-python-from-scratch/

Now I want to apply the Gaussian filter to a 3 channel(RGB) image.
For that, I implemented the code but the output that I'm getting is a blurred dull image with very low brightness. Also, the edges of the image are not blurred properly.

Here is my code:

# import libraries
import cv2
import numpy as np
import os
import matplotlib.pyplot as plt
%matplotlib inline
import cv2

# loading image
img_orig = cv2.imread('home.jpg')

# convert GBR image to RGB image
img_orig = cv2.cvtColor(img_orig, cv2.COLOR_BGR2RGB)

# Gaussian function
def dnorm(x, mu, sd):
    return 1 / (np.sqrt(2 * np.pi) * sd) * np.exp(-((x-mu)/sd)** 2 / 2)

# function for making gaussian kernel
def gaussian_kernel(kernel_size, mu = 0):
    # initializing mu and SD
    sd = np.sqrt(kernel_size)

    # creating 1D kernel
    kernel_1D = np.linspace(-(kernel_size // 2), kernel_size // 2, kernel_size)

    # normalizing 1D kernel
    for i in range(kernel_size):
        kernel_1D[i] = dnorm(kernel_1D[i], mu, sd)

    # creating 2D kernel
    kernel_2D = np.outer(kernel_1D, kernel_1D)
    kernel_2D /= kernel_2D.max()

    return kernel_2D

Here is what the 11 X 11 kernel looks like:

# Covolution function with zero padding
def convolution(image, kernel):
    # find row and column of 3 channel (RGB) image
    img_row, img_col, img_channel = image.shape

    kernel_size = kernel.shape[0]
    padding_width = (kernel_size - 1) // 2

    #initialize output image
    output = np.zeros(image.shape, dtype = np.uint8)

    # initialize padded image with zeros
    padded_img = np.zeros((img_row + 2*padding_width, img_col + 2*padding_width, img_channel), dtype = np.uint8)

    # copy orignal image inside padded image
    padded_img[padding_width : padding_width + img_row, padding_width : padding_width + img_col] = image

    # average pixel values using gaussian kernel
    for i in range(img_row):
        for j in range(img_col):
            # average each pixel's R channel value
            output[i, j, 0] = np.sum(padded_img[i : i+kernel_size, j : j+kernel_size, 0] * kernel) // (kernel_size * kernel_size)

            # average each pixel's G channel value
            output[i, j, 1] = np.sum(padded_img[i : i+kernel_size, j : j+kernel_size, 1] * kernel) // (kernel_size * kernel_size)

            # average each pixel's B channel value
            output[i, j, 2] = np.sum(padded_img[i : i+kernel_size, j : j+kernel_size, 2] * kernel) // (kernel_size * kernel_size)

    return output

def gaussian_filter(image, kernel_size = 3):
    # initialize mu
    mu = 0

    # create gaussian kernel
    kernel = gaussian_kernel(kernel_size, mu)

    # apply convolution to image
    conv_img = convolution(image, kernel)

    # return blurred image
    return conv_img

Testing code for Gaussian filter:

plt.figure(figsize = (7, 5))
print('orignal image')
plt.imshow(img_orig)
plt.show()

plt.figure(figsize = (7, 5))
print('blurred image')
plt.imshow(gaussian_filter(img_orig, 11))
plt.show()

Output:

Comparying with openCV GaussianBlur:

print('openCV blurred image')
plt.imshow(cv2.GaussianBlur(img_orig, (11,11), 0))
plt.show()

Output:

My Questions are:
1) why I am getting a dull image as output.
2) Is the above implementation of the Gaussian filter for RGB image WRONG? If it is wrong how can I make it correct?
3) Why the edges are not properly blurred (see the blackish shadow at the edges)?
4) The above implementation of the Gaussian filter is taking a very long time to execute as compared to the OpenCV GaussianBlur. How can I make it Time efficient?

Answer 1

Two things are wrong that cause the image intensity to not be preserved: you first normalize the kernel by dividing by its maximum value, then in the convolution you divide by the number of samples in the kernel. Instead of these two normalizations, normalize just once, when you create the kernel, by dividing by the sum of the kernel values. This makes the sum of the kernel weights equal 1, and causes the convolution to preserve the average image intensity. Note that the convolution computes a local weighted average; in a weighted average we need the weights to add to 1 to avoid a bias.

The dark edges are caused by the padding: you pad with zeros (black), which mixes in with the values at the edges of the image in the convolution. OpenCV likely uses a different bounday condition, or padding of the image. Typical options involve mirroring values, or simply extending the edge values out.

Finally, the main reason your code is slow is that you are using a loop in Python. Python is an interpreted language and hence slow. You could use Numba to speed up the loops (it’s a just-in-time compiler for Python), or simply use the convolution in NumPy, which is implemented in a compiled language.

The other reason your code is slow (which won’t matter much until you fix the first) is that you are not making use of the separability of the Gaussian. You build a 2D Gaussian by multiplying two 1D Gaussian, but instead you could apply two 1D convolutions in sequence. For your example of a 11x11 kernel, the computational cost is reduced from 11*11=121 multiplications and additions to 11+11=22 multiplications and additions. The larger the kernel, the better the speed gain.

Answer 2

Hi I'm a new learner as well, so please correct me if I'm wrong. Apart from Cris Luengo's answer, which I think it's great, I also made two observations:

1. in the calculation of output[i, j, 0], output[i, j, 1] and output[i, j, 2], the output from the right is a float, to be assigned to the value on the left, which is np.uint8, what it does is ignore the decimal part, i.e. 1.9 -> 1. so I think it might be better to use round() before assigning the calculated value, so that 1.9 -> 2

2. according to the definition of convolution, which I also just learnt, you need to flip the kernel upside down and left to right, before the calculation, this is to keep convolution associative. As you know, Gaussian kernel is separable, which means you need to convolve twice, vertically and horizontally, associative make sure the order doesn't matter. But of course, a Gaussian kernel is symmetrical, so it doesn't matter if you flip it nor not.