2D FFT Convolution Questions

[teeter]

I am currently using VSL's convolution functions to do 2D image convolution in C. Typically, I have 3 arrays and about 50 kernels. I implement convolution using the following pseudo code:

For each array from 1 to 3
{
For each kernel from 1 to 50
{
Call vslsConvExecX(...);
}
}

Typically for convolution using FFT, a forward FFT is applied to both the array and the kernel. The array and the kernel are multiplied in frequency space and a backward FFT is applied to the result. In the above example, it would seem like forward FFT is applied 3 times for each kernel.

To eliminate this redundancy, I considered implementing the following pseudo code:

1. Apply forward FFT to all 3 arrays.
2. For each kernel from 1 to 50
{
Apply forward FFT to the kernel
For each array from 1 to 3
{
Multiply the array and kernel in frequency space. Store results in a separate array.
Apply backward FFT to obtain the convolution results.
}
}

Here are my questions:

1. Is my approach a reasonable way to speed up my program?

2. I implemented the new algorithm using the DFTI calls in MKL, but the resultant images appear garbled. I may have set some of the DFTI parameters wrongly. Here is the way I initialized the descriptor:

status = DftiCreateDescriptor(&descriptor, DFTI_SINGLE, DFTI_REAL, 2, length);
status = DftiSetValue(descriptor, DFTI_ORDERING, DFTI_BACKWARD_SCRAMBLED);
status = DftiSetValue(descriptor, DFTI_PACKED_FORMAT, DFTI_PACK_FORMAT);
status = DftiSetValue(descriptor, DFTI_BACKWARD_SCALE, 1.0f / float(width * height));
status = DftiCommitDescriptor(descriptor);

For debugging purposes, I applied forward FFT to an array and applied backward FFT to get the same array back. I repeated this with the kernel and got the same kernel back too. But if I multiple the array and kernel in frequency space, I get a garbled output array after the backward FFT. What am I doing wrong?

Thanks.

[Dmitry_B_Intel]

Hi,

You do in-place real FFT which happened to be broken (a known issue), that is likely the reason for you see the garbage in the output. I would suggest, as a workaround, trying out-of-placetransform.

Thanks,
Dima

[apolo74]

Hi guys,

I see that this discussion has more than 3 years but I hope you can help me. I'm trying to implement more or less the same application using MKL: I have to create a convolution between a 3D matrix and 3 different kernels (one for each dimension) The size of the kernels can be different between each other and all of them are smaller than their respective dimension of the main 3D matrix.

By following the examples I have managed to create a Fourier transform of a 1D matrix(vector) just to start learning how MKL works. My 1D example has size 7 and my kernel has size 3, so I'm confused on how to do the Fourier multiplication. Should I pad the kernel with zeros??? and then:

- apply Fourier to both objects

- usevcMul(n,a,b,y); to multiply them element by element (NOT_IN_PLACE)

- Transform back FROM fourier

thanks for your help and suggestions,

Boris

PS: Is it important to control the size of these objects??? i.e. should they always be power of 2???

[Victor_Gladkikh]

Hi Boris

I recommend use convolution routines for multidimensional case with the fixed first operand vector. You can compute convolution between your 3D matrix (3D matrix is the fixed first operand in your case) and 3 kernels. This functionality provided by MKL.

Victor

[apolo74]

Hi Victor and thanks for replying and for the tip. I just noticed that it's possible to select the method, FFT or Direct convolution so I guess that they take care of my problem automatically. However I'm not sure if this tool will help me to convolve a 3D matrix with a 1D vector in the right dimension. Let's say that my 3D matrix represents space with coordinates X, Y and Z... and I have 3 different kernels to convolve that space with, a different kernel for each dimension in space. So I would like to convolve kernelX with my matrix but only in the X dimension, kernelY with the 3D matrix but only in the Y dimension and so on. Do you know if this would be possible? or am I forced to create a 3D matrix multiplying the 3 kernels???

Boris

[Victor_Gladkikh]

Boris,

MKL Convolution/Correlation routineshave?shape input parameters. As I understand you can use next shapes for kernels:
xkernels (N,1,1), where N is length of Xkernel

ykernels (1,M,1), where M is length of Ykernel

zkernels (1,1,K), where K is length of Zkernel

Victor

[apolo74]

Yes, it's perfect.

Thanks for your help,

Boris