2d convolution using python and numpy

Maybe it is not the most optimized solution, but this is an implementation I used before with numpy library for Python:

def convolution2d(image, kernel, bias):    m, n = kernel.shape    if (m == n):        y, x = image.shape        y = y - m + 1        x = x - m + 1        new_image = np.zeros((y,x))        for i in range(y):            for j in range(x):                new_image[i][j] = np.sum(image[i:i+m, j:j+m]*kernel) + biasreturn new_image

I hope this code helps other guys with the same doubt.

Regards.

Edit [Jan 2019]

@Tashus comment bellow is correct, and @dudemeister's answer is thus probably more on the mark. The function he suggested is also more efficient, by avoiding a direct 2D convolution and the number of operations that would entail.

Possible Problem

I believe you are doing two 1d convolutions, the first per columns and the second per rows, and replacing the results from the first with the results of the second.

Notice that numpy.convolve with the 'same' argument returns an array of equal shape to the largest one provided, so when you make the first convolution you already populated the entire data array.

One good way to visualize your arrays during these steps is to use Hinton diagrams, so you can check which elements already have a value.

Possible Solution

You can try to add the results of the two convolutions (use data[:,c] += .. instead of data[:,c] = on the second for loop), if your convolution matrix is the result of using the one dimensional H_r and H_c matrices like so:

Another way to do that would be to use scipy.signal.convolve2d with a 2d convolution array, which is probably what you wanted to do in the first place.

Since you already have your kernel separated you should simply use the sepfir2d function from scipy:

from scipy.signal import sepfir2dconvolved = sepfir2d(data, H_r, H_c)

On the other hand, the code you have there looks all right ...