How can I smooth elements of a two-dimensional array with differing gaussian functions in python?

astro_ash picture astro_ash · Nov 5, 2015 · Viewed 14.2k times · Source

How could I smooth the x[1,3] and x[3,2] elements of the array,

x = np.array([[0,0,0,0,0],[0,0,0,1,0],[0,0,0,0,0],[0,0,1,0,0],[0,0,0,0,0]])   

with two two-dimensional gaussian functions of width 1 and 2, respectively? In essence I need a function that allows me to smooth single "point like" array elements with gaussians of differing widths, such that I get an array with smoothly varying values.

Answer

Flabetvibes picture Flabetvibes · Nov 9, 2015

I am a little confused with the question you asked and the comments you have posted. It seems to me that you want to use scipy.ndimage.filters.gaussian_filter but I don't understand what you mean by:

[...] gaussian functions with different sigma values to each pixel. [...]

In fact, since you use a 2-dimensional array x the gaussian filter will have 2 parameters. The rule is: one sigma value per dimension rather than one sigma value per pixel.

Here is a short example:

import matplotlib.pyplot as pl
import numpy as np
import scipy as sp
import scipy.ndimage

n = 200 # widht/height of the array
m = 1000 # number of points
sigma_y = 3.0
sigma_x = 2.0

# Create input array
x = np.zeros((n, n)) 
i = np.random.choice(range(0, n * n), size=m)
x[i / n, i % n] = 1.0

# Plot input array
pl.imshow(x, cmap='Blues', interpolation='nearest')
pl.xlabel("$x$")
pl.ylabel("$y$")
pl.savefig("array.png")

# Apply gaussian filter
sigma = [sigma_y, sigma_x]
y = sp.ndimage.filters.gaussian_filter(x, sigma, mode='constant')

# Display filtered array
pl.imshow(y, cmap='Blues', interpolation='nearest')
pl.xlabel("$x$")
pl.ylabel("$y$")
pl.title("$\sigma_x = " + str(sigma_x) + "\quad \sigma_y = " + str(sigma_y) + "$")
pl.savefig("smooth_array_" + str(sigma_x) + "_" + str(sigma_y) + ".png")

Here is the initial array:

enter image description here

Here are some results for different values of sigma_x and sigma_y:

enter image description here

enter image description here

enter image description here

enter image description here

This allows to properly account for the influence of the second parameter of scipy.ndimage.filters.gaussian_filter.

However, according to the previous quote, you might be more interested in the assigement of different weights to each pixel. In this case, scipy.ndimage.filters.convolve is the function you are looking for. Here is the corresponding example:

import matplotlib.pyplot as pl
import numpy as np
import scipy as sp
import scipy.ndimage

# Arbitrary weights
weights = np.array([[0, 0, 1, 0, 0],
                    [0, 2, 4, 2, 0],
                    [1, 4, 8, 4, 1],
                    [0, 2, 4, 2, 0],
                    [0, 0, 1, 0, 0]],
                   dtype=np.float)
weights = weights / np.sum(weights[:])
y = sp.ndimage.filters.convolve(x, weights, mode='constant')

# Display filtered array
pl.imshow(y, cmap='Blues', interpolation='nearest')
pl.xlabel("$x$")
pl.ylabel("$y$")
pl.savefig("smooth_array.png")

And the corresponding result:

enter image description here

I hope this will help you.