4
$\begingroup$
def fft_convolution(image, kernel):
    imageC = image.copy()
    kSize = (kernel.size // 2)+1
    kernelShape = tuple(ti//2 for ti in kernel.shape)
    centery ,centerx = kernelShape
    imagePad = np.pad(image, kSize, mode = "edge")
    imageC = np.fft.fft2(imagePad)
    kernel= np.fft.fft2(kernel, imagePad.shape)
    output = np.real(np.fft.ifft2(np.multiply(imageC, kernel)))
    return output[kSize+centery:-kSize+centery, kSize+centerx:-kSize+centerx]

So my code currently takes an image M x N and kernel P x K Pad on all sides the kernelSize // 2 + 1 and then make the kernel return from fft2 as the correct shape

I then crop using the shape of the kernel from my pad size, I dont really know why, something to do with circular convolutions maybe.

But this is really confusing and makes no sense to me especially when we get to cropping where I had to deal with artifacts of cylic convolution.

I am wondering what is the correct method to apply cropping and padding perfectly, taking into account i want my image to stay the same size and the edges must be interpolated rather than using constant zero values.

$\endgroup$
3
  • $\begingroup$ What would the corresponding code in spatial domain look like? "the edges must be interpolated rather than using constant zero values."? Extrapolated? From what? $\endgroup$
    – Oscar
    Nov 1, 2021 at 14:54
  • $\begingroup$ Basically there is an image and kernel in spatial domain, i must prepare for the frequency domain multiplicaiton by convolution theorem. There is something I must do to the kernel in terms of alignign the axis which i am not sure about $\endgroup$ Nov 1, 2021 at 16:25
  • 1
    $\begingroup$ Related to dsp.stackexchange.com/questions/79145. $\endgroup$
    – Royi
    Nov 21, 2021 at 17:22

1 Answer 1

0
$\begingroup$

I created a MATLAB function which is basically conv2() applied in Frequency Domain:

function [ mO ] = ImageConvFrequencyDomain( mI, mH, convShape )
% ----------------------------------------------------------------------------------------------- %
% [ mO ] = ImageConvFrequencyDomain( mI, mH, convShape )
% Applies Image Convolution in the Frequency Domain.
% Input:
%   - mI                -   Input Image.
%                           Structure: Matrix.
%                           Type: 'Single' / 'Double' (Single Channel).
%                           Range: (-inf, inf).
%   - mH                -   Filtering Kernel.
%                           Structure: Matrix.
%                           Type: 'Single' / 'Double'.
%                           Range: (-inf, inf).
%   - convShape         -   Convolution Shape.
%                           Sets the convolution shape.
%                           Structure: Scalar.
%                           Type: 'Single' / 'Double'.
%                           Range: {1, 2, 3}.
% Output:
%   - mI                -   Output Image.
%                           Structure: Matrix (Single Channel).
%                           Type: 'Single' / 'Double'.
%                           Range: (-inf, inf).
% References:
%   1.  MATLAB's 'conv2()' - https://www.mathworks.com/help/matlab/ref/conv2.html.
% Remarks:
%   1.  A
% TODO:
%   1.  
%   Release Notes:
%   -   1.0.000     29/04/2021  Royi Avital     [email protected]
%       *   First release version.
% ----------------------------------------------------------------------------------------------- %

CONV_SHAPE_FULL     = 1;
CONV_SHAPE_SAME     = 2;
CONV_SHAPE_VALID    = 3;

numRows     = size(mI, 1);
numCols     = size(mI, 2);

numRowsKernel = size(mH, 1);
numColsKernel = size(mH, 2);

switch(convShape)
    case(CONV_SHAPE_FULL)
        numRowsFft  = numRows + numRowsKernel - 1;
        numColsFft  = numCols + numColsKernel - 1;
        firstRowIdx = 1;
        firstColIdx = 1;
        lastRowIdx  = numRowsFft;
        lastColdIdx = numColsFft;
    case(CONV_SHAPE_SAME)
        numRowsFft  = numRows + numRowsKernel;
        numColsFft  = numCols + numColsKernel;
        firstRowIdx = ceil((numRowsKernel + 1) / 2);
        firstColIdx = ceil((numColsKernel + 1) / 2);
        lastRowIdx  = firstRowIdx + numRows - 1;
        lastColdIdx = firstColIdx + numCols - 1;
    case(CONV_SHAPE_VALID)
        numRowsFft = numRows;
        numColsFft = numCols;
        firstRowIdx = numRowsKernel;
        firstColIdx = numColsKernel;
        % The Kernel when transformed is shifted (Namely its (0, 0) is top
        % left not middle).
        lastRowIdx  = numRowsFft;
        lastColdIdx = numColsFft;
end

mO = ifft2(fft2(mI, numRowsFft, numColsFft) .* fft2(mH, numRowsFft, numColsFft), 'symmetric');
mO = mO(firstRowIdx:lastRowIdx, firstColIdx:lastColdIdx);


end


It is fully compatible and validated.
The full code is available on my StackExchange Signal Processing Q74803 GitHub Repository (Look at the SignalProcessing\Q74803 folder).

If you want to apply image filter with Constant / Mirror or Replicate boundary mode just pad the image before the function and use CONV_SHAPE_VALID as convShape. Without padding the boundary condition applied is Periodic.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.