Prev: damage measurement in a micrograph
Next: save timeseries to file
From: Eric on 9 Feb 2010 13:20
This post is a bit stale, but here's my take:

There are times when the fast Fourier transform is not necessarily the fastest way to get the data you're looking for. If you want a high-resolution point-spread function (PSF) from a wavefront, one's first thought is you have to pad the heck out of the wavefront and then FFT. In this case you're using a fast algorithm on a large data set.

However, you know a priori that that the PSF is going to be zero outside of a very small region in the middle of the array. So if you perform the algorithm described above you spend a lot of time calculating data you don't really need: all of the data outside of the PSF.

Often times a better implementation is to use matrix multiplication techniques and only calculate the central portion of the Fourier transform. See Matlab's dftmtx function in the Signal Processing Toolbox and "Fast computation of Lyot-style coronagraph propagation", R. Soummer, L. Pueyo, A. Sivaramakrishnan, and R. J. Vanderbei, Optics Express, Vol 15 No 24, 15935-15951 (2007).

The idea here is to use a slow algorithm on a small data set. If I'm only interested in the central 256x256 pixels of a 5120x5120 Fourier transform plane, the computation time on my machine is 39.6 milliseconds. The time to take FFT2 of a 5120x5120 array is 2902 milliseconds. Thus the speed-up is a factor of 73. Of course you have to make sure that the PSF is fully contained in the 256x256 region, I'm just trying to give an example.

Furthermore, you can calculate the central region of a PSF that you could not calculate with FFT2 due to memory constraints. I've calculated the central 1024x1024 region of a 51200x51200 on a 32-bit machine with no problem. Clearly this would not be possible with FFT2.

-Eric
 | 
Pages: 1
Prev: damage measurement in a micrograph
Next: save timeseries to file