File Exchange

image thumbnail

FFT-based convolution

version 1.7.0.2 (5.78 KB) by Bruno Luong
Discrete convolution using FFT method

28 Downloads

Updated 20 Sep 2018

View License

As opposed to Matlab CONV, CONV2, and CONVN implemented as straight forward sliding sums, CONVNFFT uses Fourier transform (FT) convolution theorem, i.e. FT of the convolution is equal to the product of the FTs of the input functions.
In 1-D, the complexity is O((na+nb)*log(na+nb)), where na/nb are respectively the lengths of A and B.

Optional arguments to control the dimension(s) along which convolution is carried out.

Slightly less accurate than sliding sum convolution.

Good usage recommendation:
In 1D, this function is faster than CONV for nA, nB > 1000.
In 2D, this function is faster than CONV2 for nA, nB > 20.
In 3D, this function is faster than CONVN for nA, nB > 5.

Cite As

Bruno Luong (2019). FFT-based convolution (https://www.mathworks.com/matlabcentral/fileexchange/24504-fft-based-convolution), MATLAB Central File Exchange. Retrieved .

Comments and Ratings (37)

Refaat Gabr

Per-Ivar

Ryan Muir

I'm running Matlab 2017a. It installs successfully compiling with MinGW, but crashes matlab when I try to run it.

Thanks Bruno, it works and still very helpful.

Bruno Luong

Vijendhar, I don't have access to MacOS.
You can replace the inplaceprod(A,B) with
A = A.*B;
it works with more memory required.

Hello Bruno. I used convnfft function in my office computer, Windows 2010, MATLAB 2015a. It worked very well and I saved a lot of time. Thanks for great work.

Unfortunately, I get error with my PC, MacOS, MATLAB 2018a. It keeps saying an error,

Undefined function or variable 'inplaceprod'.

Error in convnfft (line 166)
inplaceprod(A,B);

Any help? I did some research and got to know that the MEX file should be compiled for the mexmaci64 platform to run on Apple Mac. Is it possible for you to provide MEX file of inplaceprod in mexmaci64 platform? I'm actually new to Matlab, please correct my question if I'm wrong and provide a solution for the error. Thanks

Jered Wells

Thank you for updating the mex call. This is a wonderful tool! The folks at Mathworks would be wise to implement FFT-based math options in their core convolution functions.

Jered Wells

Bruno, up until MATLAB R2017a, this function works fine, but it breaks when running in MATLAB 2018a. Any ideas? No errors get thrown, but image A just seems to get "post" zero padded with a border with same dimensions as B.

Sorry, I'm obviously doing something wrong.

First, this is way slower than convn (but given the size of the matrices it shouldn't be). And second I'm getting different answers.

Help?

i=200; pwdc=30;
stimImg=randn(i, i, 240);

hemo=randn(30, 1);
tic
for r=1:1
tmp1 = convnfft(stimImg, hemo(:), 'valid');
size(tmp1)
end
toc
tic
for r=1:1
tmp2 = convn(stimImg, hemo(:), 'valid');
size(tmp2)
end
toc

round(sum(abs(tmp1(:))-abs(tmp2(:))))

Excellent function Bruno. I have a suggestion. Would it be possible to do the convolution with NaNs? Also, sometimes when using the 'same' argument for the shape it would be interesting to add NaNs to the edges instead of zeros. This would eliminate artifacts one gets near the edges of the convolution.

wtmlma

AP

This is what I experimented with the code.

>> tic,A = convnfft(rand(300,300,300), ones(5,5,5), 'same');toc
Elapsed time is 8.061082 seconds.

>> tic,A = convn(rand(300,300,300), ones(5,5,5), 'same');toc
Elapsed time is 2.085360 seconds.

>> 8.061082/2.085360
ans =
3.8656

Ian

I am running 2014a on a machine with 192Gb of RAM and 20 cores. I am trying to convolute two vectors, one with 3,060,663 elements, the other with 693. The built-in conv took 0.06 seconds. convnfft filled the memory and then crashed the machine.

I downloaded your code and tried to install launching the installation function from Matlab command line. I use a Mac running Maverick (10.9.4)
I got the following error that I am copying in the following.
I verified that clang is installed:
f4230:~ mauede$ which clang
/usr/bin/clang
I would greatly appreciate your help.Thank you so much.
Regards,
Maura Monville

>> convnfft_install
-> mexopts.sh sourced from directory (DIR = $MATLAB/bin)
FILE = /Applications/MATLAB_R2013b.app/bin/mexopts.sh
----------------------------------------------------------------
-> MATLAB = /Applications/MATLAB_R2013b.app
-> CC = xcrun -sdk macosx10.7 clang
-> CC flags:
CFLAGS = -fno-common -arch x86_64 -isysroot -mmacosx-version-min=10.7 -fexceptions
CDEBUGFLAGS = -g
COPTIMFLAGS = -O2 -DNDEBUG
CLIBS = -L/Applications/MATLAB_R2013b.app/bin/maci64 -lmx -lmex -lmat -lstdc++
arguments =
-> CXX = xcrun -sdk macosx10.7 clang++
-> CXX flags:
CXXFLAGS = -fno-common -fexceptions -arch x86_64 -isysroot -mmacosx-version-min=10.7
CXXDEBUGFLAGS = -g
CXXOPTIMFLAGS = -O2 -DNDEBUG
CXXLIBS = -L/Applications/MATLAB_R2013b.app/bin/maci64 -lmx -lmex -lmat -lstdc++
arguments =
-> FC = gfortran
-> FC flags:
FFLAGS = -fexceptions -m64 -fbackslash
FDEBUGFLAGS = -g
FOPTIMFLAGS = -O
FLIBS = -L/Applications/MATLAB_R2013b.app/bin/maci64 -lmx -lmex -lmat -L -lgfortran -L -lgfortranbegin
arguments =
-> LD = xcrun -sdk macosx10.7 clang
-> Link flags:
LDFLAGS = -arch x86_64 -Wl,-syslibroot, -mmacosx-version-min=10.7 -bundle -Wl,-exported_symbols_list,/Applications/MATLAB_R2013b.app/extern/lib/maci64/mexFunction.map
LDDEBUGFLAGS = -g
LDOPTIMFLAGS = -O
LDEXTENSION = .mexmaci64
arguments =
-> LDCXX =
-> Link flags:
LDCXXFLAGS =
LDCXXDEBUGFLAGS =
LDCXXOPTIMFLAGS =
LDCXXEXTENSION =
arguments =
----------------------------------------------------------------

-> xcrun -sdk macosx10.7 clang -c -I/Applications/MATLAB_R2013b.app/extern/include -I/Applications/MATLAB_R2013b.app/simulink/include -DMATLAB_MEX_FILE -fno-common -arch x86_64 -isysroot -mmacosx-version-min=10.7 -fexceptions -O2 -DNDEBUG "inplaceprod.c"

xcodebuild: error: SDK "macosx10.7" cannot be located.
xcrun: error: unable to find utility "clang", not a developer tool or in PATH

mex: compile of ' "inplaceprod.c"' failed.

Unable to complete successfully.

Error in convnfft_install (line 17)
mex(mexopts{:},'inplaceprod.c');

Matt Taylor

This function is indeed faster than CONV, but as soon as I attempted to use it on larger data sets, Matlab produced an 'out of memory' error, whereas CONV can cope just fine with the same datasets (albeit taking longer).

FYI if I run the 'memory' command my output is as follows:
Maximum possible array: 11862 MB (1.244e+10 bytes) *
Memory available for all arrays: 11862 MB (1.244e+10 bytes) *
Memory used by MATLAB: 820 MB (8.597e+08 bytes)
Physical Memory (RAM): 8011 MB (8.400e+09 bytes)

So the problem definitely isn't my hardware

Excellent job! Nicely documented and elegant code and to the point!

Works much faster than conv2 for full case, and also faster than conv2 with option 'valid', which misteriously makes conv2 35x faster with a 500x500 matrix with a 400x400 one (makes me suspect that conv2 + 'valid' does not just extract the mid part but saves computations).

Bruno Luong

Young, You must change the directory where convfft is installed (including the c code) to install it.

Young Gyu

Hi Bruno,

First of all, thanks for the nice job. Unfortunately, I'm having problem with running it. When I try to run convnfft_install, it keeps saying
--------------------------------------
C:\PROGRA~1\MATLAB\R2013B\BIN\MEX.PL: Error: 'inplaceprod.c' not found.
Error in convnfft_install (line 17)
mex(mexopts{:},'inplaceprod.c');
--------------------------------------

Do you have any suggestion?

Thanks,
YG

Petr

Bruno, sorry, I forgot to mention that I used 2012b. Moreover, speed difference might depend on inputs and I'm lucky with the inputs I use?

Bruno Luong

Petr, I just test with 2012a, the recommendation stands.

Petr

Hi Bruno, as for usage recommendations, are they updated for new MATLAB releases? In case of my application for 2D, both nA and nB are greater than 20 (about 200*300 each). However, MATLAB conv2 takes almost the same time as convnfft (even slightly faster).

Luke Phai

Thanks Bruno.

Bruno Luong

The GPU acceleration depends on user's hardware. It is impossible to give reliable number without what is your computer setup. A test I run long ago shows an acceleration between 3-5 times.

Bruno Luong

Sorry, Luc -> Luke

Luke Phai

Excellent work Bruno. Many thanks.
Quick question. How much faster it would be with the GPU jacket enabled?
I try GPUmat using fft2 to programme similar code but it turns out to be slower.
I am thinking to get MATLAB Parallel Computing Toolbox to run the GPU if it is a lot faster.
I hope it is.
Could you please give me some idea? Say
A = rand(1000,1000);
B = rand(1000,1000);
tic;C=convnfft(A, B, 'same', [1, 2],'false');toc
given me> 0.213153 seconds without GPU

tic;C=convnfft(A, B, 'same', [1, 2],'true');toc
time???
Thanks

Nikolay S.

Hello again. Apparently it's Matlab filter2 function fault: "Given a matrix X and a two-dimensional FIR filter h, filter2 rotates your filter matrix 180 degrees to create a convolution kernel."
Why the rotation is needed? hack knows...

Nikolay S.

Hi Bruno. Excellent contribution and elegant code.
A few small comments:
1) The speed up is smaller then the one you state- as conv is optimized both by both Мatlab and by CPU vendor (Intel in my case).
2) The convolution shift in your (and btw mine) is different form the one resulting from filter2 function. Try running the following code:

img=rgb2gray(imread('peppers.png'));
filt=zeros(50);
filt(1,1)=1;

convFFT=convnfft(img, filt, 'same', [1, 2]);
regConv=filter2(filt, img);

figure;
subplot(1, 2, 1);
imshow(uint8(convFFT));
title('FFT based convolution');
subplot(1, 2, 2);
imshow(uint8(regConv));
title('Matlab regular convolution');

Jonathan

Edwin

Moreover, when I checked the result form this code, there are some different between this and convn for two 3d matrix, the 2 input matrix are both positive.
the result from this code has negative values while convn does not.

Bruno Luong

Hi Michael,
May be the inplace times is no longer necessary for recent Matlab.

I remember implement that from a user request.

Thanks for the useful feedback.

Hi Bruno, are you sure your inplaceprod() is (still) useful?

I'm pretty sure, MATLAB does A=A.*B in-place itself. I just compared my memory usage for very large A/B with both methods and there was no difference.
This post is from 2007: http://blogs.mathworks.com/loren/2007/03/22/in-place-operations-on-data/

In a heterogeneous environment, it is useful to avoid mex code for such small tasks. If you are a non-privileged user on a compute server it is really a mess when compiling fails due to compiler version, libraries or whatever.

A minor notice is that (i)fftn is faster than for-loops around 1D (i)fft calls. At least as long as the input and output are of the same size. So I got at least a little speed gain by replacing
for dim=dims
A = ifft(A,[],dim);
end
with
A = ifftn( A );
for the MATLAB ifft case.

Thank you for the code!

Romesh

I think this says it all...

>> tic;C = convn(Vs,Vs);toc;
Elapsed time is 473.103412 seconds.
>> tic;C2 = convnfft(Vs,Vs);toc;
Elapsed time is 1.351315 seconds.
>> max(max(max(abs(C-C2))))
ans =
5.2208e-15

Thanks so much for this!

Arun

Well written (IMHO).

Alexander

Awesome function! My code runs 60x faster now (thanks to your GPU support).

Eric

Updates

1.7.0.2

Installation script use now right mex compilation options for R2018a or later

1.7.0.1

Make inplaceprod compatible with interleaved complex

1.7.0.0

Add the syntax conv2fft(H1, H2, A, ...)

1.6.0.0

Option allows to disable padding to next power-two. Mex implement inplace product that saves about 1/3 memory. These two enhancement might be useful when perform convolution with very large arrays.

1.5.0.0

GPU unable by default + changes in help section

1.4.0.0

GPU/Jacket capable

1.1.0.0

correct bug when ndims(A)<ndims(B)

MATLAB Release Compatibility
Created with R2009a
Compatible with any release
Platform Compatibility
Windows macOS Linux
Acknowledgements

Inspired: conv2fft_reuse, Matching pursuit for 1D signals