I'm having trouble diagnosing and fixing this error. I'm trying to write the OaA algorithm, described in this paper.
import numpy as np
from numpy.fft import fft2 as FFT, ifft2 as IFFT
class convolve(object):
""" contains methods to convolve two images with different convolution
algorithms """
def __init__(self, image_array, kernel):
self.array = image_array
self.kernel = kernel
self.__rangeX_ = self.array.shape[0]
self.__rangeY_ = self.array.shape[1]
self.__rangeKX_ = self.kernel.shape[0]
self.__rangeKY_ = self.kernel.shape[1]
def OaAconvolve(self):
""" faster convolution algorithm, O(N^2*log(n)). """
# solve for the total padding along each axis
diffX = (self.__rangeKX_ -
(self.__rangeX_ - self.__rangeKX_*(
self.__rangeX_ // self.__rangeKX_)))
% self.__rangeKX_
diffY = (self.__rangeKY_ -
(self.__rangeY_ - self.__rangeKY_*(
self.__rangeY_ // self.__rangeKY_)))
% self.__rangeKY_
# each side, i.e. left, right, top and bottom
right = diffX // 2
left = diffX - right
bottom = diffY // 2
top = diffY - bottom
# pad the array [(top, bottom), (left, right)]
self.array = np.pad(self.array, [(left, right), (top, bottom)],
mode='constant', constant_values=0)
# return a list of tuples to partition the array
subsets = [(i*self.__rangeKX_, (i + 1)*self.__rangeKX_,
j*self.__rangeKY_, (j + 1)*self.__rangeKY_)
for i in xrange(self.array.shape[0] // self.__rangeKX_)
for j in xrange(self.array.shape[1] // self.__rangeKY_)]
# set up OaA method
padX = self.__rangeKX_ // 2
padY = self.__rangeKY_ // 2
self.kernel = np.pad(self.kernel, [(padY, padY), (padX, padX)],
mode='constant', constant_values=0)
transformed_kernel = FFT(self.kernel)
# create a blank array of the same size
convolved_image = np.zeros([self.array.shape[0] + 2*padX,
self.array.shape[1] + 2*padY])
print self.array.shape
# transform each partition and OaA on the convolved_image
for tup in tqdm(subsets):
# slice and pad the array subset
subset = np.pad(self.array[tup[0]:tup[1], tup[2]:tup[3]],
[(padY, padY), (padX, padX)],
mode='constant', constant_values=0)
# transform
transformed_subset = FFT(subset)
# multiply the two arrays entrywise and take the IFFT. np.real()
# is used because some residual/negligible imaginary terms are
# left over after the IFFT.
space = np.real(IFFT(transformed_kernel*transformed_subset))
# overlap with indices and add them together to build the image
convolved_image[tup[0]:tup[1] + 2*padX,tup[2]:tup[3] + 2*padY] += space
# crop image and get it back, convolved
return convolved_image[padX + left:padX + left + self.__rangeX_,
padY + bottom:padY + bottom + self.__rangeY_]
But now, when I call it, I get black bars as follows:
Here's an example Gaussian kernel I am using.
[[ 0. 0.02390753 0.03476507 0.02390753 0. ]
[ 0.02390753 0.06241541 0.07990366 0.06241541 0.02390753]
[ 0.03476507 0.07990366 0.10040324 0.07990366 0.03476507]
[ 0.02390753 0.06241541 0.07990366 0.06241541 0.02390753]
[ 0. 0.02390753 0.03476507 0.02390753 0. ]]
I believe I've narrowed the problem down to the multiplication and IFFT
space = np.real(IFFT(transformed_kernel*transformed_subset))
I think it has something to do with the IFFT of discrete Gaussian kernel (for some reason). It's odd because if I just plot
space = np.real(IFFT(transformed_subset))
I get the following (no black bars, and it pieces it back together fine):
And if I plot the opposite, i.e.
space = np.real(IFFT(transformed_kernel))
I again get no black bars, and it appears to place them in the right places.
What am I missing? I've been staring at this for days, editing indices and whatnot, but I can't get rid of this tessellation !!
Aucun commentaire:
Enregistrer un commentaire