Is there a equivalent of scipy.signal.deconvolve for 2D arrays?

These functions using fftn, ifftn, fftshift and ifftshift from the SciPy's fftpack package seem to work:

from scipy import fftpackdef convolve(star, psf):    star_fft = fftpack.fftshift(fftpack.fftn(star))    psf_fft = fftpack.fftshift(fftpack.fftn(psf))    return fftpack.fftshift(fftpack.ifftn(fftpack.ifftshift(star_fft*psf_fft)))def deconvolve(star, psf):    star_fft = fftpack.fftshift(fftpack.fftn(star))    psf_fft = fftpack.fftshift(fftpack.fftn(psf))    return fftpack.fftshift(fftpack.ifftn(fftpack.ifftshift(star_fft/psf_fft)))star_conv = convolve(star, psf)star_deconv = deconvolve(star_conv, psf)f, axes = plt.subplots(2,2)axes[0,0].imshow(star)axes[0,1].imshow(psf)axes[1,0].imshow(np.real(star_conv))axes[1,1].imshow(np.real(star_deconv))plt.show()

The image in the left bottom shows the convolution of the two Gaussian functions in the upper row, and the reverse of the effects of convolution is shown in the bottom right.

Note that deconvolving by division in the fourier domain isn't really useful for anything but demonstration purposes; any kind of noise, even numerical, may render your outcome completely unusable. One may regularize the noise in various ways; but in my experience, an RL iteration is easier to implement, and in many ways more physically justifiable.