Search code examples
matlabfftspectrum

how to find the peak(r,theta) from DFT spectrum in matlab?

Suppose i have a matrix M I apply DFT on it and shift it to the centre

img_fft = fft2(double(M));
img_spec = abs(img_fft);
img_shift = fftshift(img_spec);

Then according to my question, how could I find the peak of the spectrum ,i.e. (u,v) of the peak? Thanks a lot!

Solution

  • Use find to determine where in the image the maximum value occurs. However, the question says to ignore the DC value, which will very likely be the largest component in your spectrum so far. As such, what you should do is find the maximum value, which will should be the DC value currently, and set it to NaN. Once you do this, the next maximum value is what you're searching for to answer the question.

    As such, do this first:

    idx = find(img_shift == max(img_shift(:)));
img_shift(idx) = NaN;

    If you set the maximum value here to NaN, the next time you call find, it will disregard this location and skip this location when searching for the maximum value. Now, call find again:

    [u,v] = find(img_shift == max(img_shift(:)));

    Here, u would be the horizontal location of your image while v would be the vertical position. This should now properly return the maximum peak that is not the DC value. However, this isn't with respect to the centre as you have shifted the image, so once you're done with the above, you need to do:

    [rows, cols] = size(img_fft);

u = u - floor(rows/2);
v = v - floor(cols/2);

    This will now shift the origin to the centre of the image so this should now find where the maximum peak is with respect to the centre of the image.

    You can now use u and v to determine r and theta according to your question. In other words:

    r = sqrt(u^2 + v^2);
theta = atan2(v, u);

    A small note according to your comment is that you are encountering two values that share the maximum. This is because the spectrum is symmetric. After you do find, you will get two answers. As such, either answer will give you the right result. You can either choose any one of the r,theta pairs, or you can modify your find statement so that it only returns one answer. In other words, do this:

    [u,v] = find(img_shift == max(img_shift(:)), 1);

    Note the second parameter of 1 at the end. This means only return one result.