MATLAB Answers

0

Using FFT in for-loop is extremely slow - how to accelerate?

Asked by Apoorva Ayyagari on 5 Nov 2019
Latest activity Commented on by Apoorva Ayyagari on 19 Nov 2019 at 23:05
Below is the code that I have set up to run FFT on img_norm, a 4D data set. This code takes nearly 1 hour to run partly because I need "n" in the fft function to be 10,000 for my analysis. Any suggestions or thoughts on how I can speed up this calculation would be greatly appreciated! Thank you.
xdim = 128;
ydim = 44;
zdim = 25;
newspacing = 10000;
ft = zeros(xdim, ydim, zdim, newspacing);
ftshift = zeros(xdim, ydim, zdim, newspacing);
ftshiftmag = zeros(xdim, ydim, zdim, newspacing);
phase = zeros(xdim, ydim, zdim, newspacing);
for x = 1:xdim
for y = 1:ydim
for z = 1:zdim
if mask(x,y,z) == 1
ft(x,y,z,:) = fft(squeeze(img_norm(x,y,z,:)),newspacing); %calculate fourier transform of time series in each voxel
ftshift(x,y,z,:) = fftshift(ft(x,y,z,:)); %shift fourier transorm of each time series to be centered around 0
ftshiftmag(x,y,z,:) = abs(ftshift(x,y,z,:)); %get the magnitude of the shifted fourier transform
phase(x,y,z,:) = angle(ftshift(x,y,z,:)); %calculate the phase angle
end
end
end
end

  1 Comment

Your results should be consuming 60 GB in double floating point. Do you really have that much RAM?

Sign in to comment.

Products

1 Answer

Answer by Matt J
on 5 Nov 2019
Edited by Matt J
on 5 Nov 2019
 Accepted Answer

To conserve memory, I only store the results for x,y,z inside the mask. It would be straightforward to re-embed them in 4D arrays if you really have sufficient RAM.
I=img_norm(mask(:),:);
ft=fft(I,newspacing,2);
ftshift=fftshift(ft,2);
ftshiftmag=abs(ftshift);
phase=angle(ftshift);

  1 Comment

Sign in to comment.