How can I find the peak location from the FFT output?

Hi,

Please, how can I find the peak location from the FFT output? I have attached my m file

Answer 1

[Star Strider]

Missing data file.

However this is straightforwazrd —

Fs = 1000;
L = 10;
t = linspace(0,  L*Fs, L*Fs+1)/Fs;

% [1 2 3]'*t

% return

s = sum(sin([1 100:100:400]'*t*2*pi),1).';
t = t(:);

figure
plot(t, s)
grid

Fn = Fs/2;
L = numel(t);
NFFT = 2^nextpow2(L);
FTs = fft(s.*hann(L), NFFT)/L;
Fv = linspace(0, 1, NFFT/2+1)*Fn;
Iv = 1:numel(Fv);

[pks,locs] =  findpeaks(abs(FTs(Iv))*2, 'MinPeakProminence',0.25)

pks = 5×1

0.4825 0.4810 0.4952 0.4952 0.4810

locs = 5×1

17 1639 3278 4916 6555

Results = table(Fv(locs).', pks, 'VariableNames',{'Frequency','Peak Magnitude'})

Results = 5×2 table
    Frequency    Peak Magnitude
    _________    ______________

     0.97656        0.48247    
      99.976        0.48103    
      200.01        0.49516    
      299.99        0.49516    
      400.02        0.48103    

figure
plot(Fv, abs(FTs(Iv))*2)
grid
ylim([min(ylim) 1.5*max(ylim)])
text(Results{:,1}, Results{:,2}, compose('\\leftarrow Freq = %.1f\n    Mag = %.3f',Results{:,[1  2]}), 'Rotation',45)

That approach should also work for your data.

.

[University]

Thank you so much. I have tried this but I received an error:

Fn = Fs1/2; % the Fs1 is my frequency

L = numel(t1);

NFFT = 2^nextpow2(L);

FTs = fft(theta1.*hann(L), NFFT)/L; % theta1 is the signal

Fv = linspace(0, 1, NFFT/2+1)*Fn;

Iv = 1:numel(Fv);

[pks,locs] = findpeaks(abs(FTs(Iv))*2, 'MinPeakProminence',0.25);

%Results = table(Fv(locs).', pks, 'VariableNames',{'Frequency','Peak Magnitude'});

I have also attached the excel file.

[University]

Hi Star,

Thank you so much. I have been able to do it.

[Star Strider]

As always, my pleasure!

[University]

Hi Star,

Sorry to come back. Please how can I perform a spatial frequency. I have the same data but instead of time I have x. Please how can I perform the FFT ? I have checked Matlab documentation but couldn't get how it works. I have attached the excel file.

[Star Strider]

As always, my pleasure!

The Fourier transform will work regardless of the independnet variable units. If theat is time, the frequency units are in cycles/time unit. If spatial, cycles/spatial unit (for example cycles/millimetre).

Try this —

T1 = readtable('Fourier_data101.xlsx')

T1 = 101×2 table
    time     theta 
    ____    _______

       0    0.77481
    0.02    0.35591
    0.04    0.18156
    0.06    0.21909
    0.08    0.21391
     0.1    0.24668
    0.12    0.29114
    0.14    0.31926
    0.16     0.3504
    0.18    0.37767
     0.2    0.39259
    0.22    0.39536
    0.24    0.38582
    0.26    0.36519
    0.28    0.33557
     0.3    0.29983

t1 = T1{:,1};
s1 = T1{:,2};

T2 = readtable('FFT_spatial.xlsx')

T2 = 857×2 table
         x            theta   
    ___________    ___________

         -1e-05    -1.6396e-19
         -1e-05     5.2143e-07
    -9.9543e-06      0.0017766
    -9.9099e-06       0.003522
    -9.9087e-06      0.0035718
    -9.9074e-06      0.0036231
     -9.863e-06      0.0053857
    -9.8196e-06      0.0070574
     -9.817e-06      0.0071556
    -9.8144e-06      0.0072566
    -9.7711e-06      0.0088551
    -9.7292e-06       0.010345
    -9.7251e-06       0.010482
    -9.7212e-06       0.010617
    -9.6782e-06       0.012015
    -9.6357e-06       0.013347

t2 = T2{:,1};
s2 = T2{:,2};

Fs2 = 1/mean(diff(t2))

Fs2 = 4.2800e+07

[s2r,t2r] = resample(s2, t2, Fs2);                          % Resample To Constant Sampling Frequency


[FTs1,Fv1] = FFT1(s1,t1);                                   % Fourier_data101


[pks1,locs1] =  findpeaks(abs(FTs1)*2, 'MinPeakProminence',0.05)

pks1 = 0.1120

locs1 = 7

Results1 = table(Fv1(locs1).', pks1, 'VariableNames',{'Frequency','Peak Magnitude'})

Results1 = 1×2 table
    Frequency    Peak Magnitude
    _________    ______________

     2.3438         0.11201    

figure
plot(Fv1, abs(FTs1)*2)
grid
ylim([min(ylim) 1.5*max(ylim)])
xlabel('Frequency')
ylabel('Magnitude')
title('Fourier\_data101')
text(Results1{:,1}, Results1{:,2}, compose('\\leftarrow Freq = %.1f\n    Mag = %.3f',Results1{:,[1  2]}), 'Rotation',45)

[FTs2,Fv2] = FFT1(s2r,t2r);                                 % FFT_spatial

[pks2,locs2] =  findpeaks(abs(FTs2)*2, 'MinPeakProminence',0.1)

pks2 = 0.2211

locs2 = 15

Results2 = table(Fv2(locs2).', pks2, 'VariableNames',{'Frequency','Peak Magnitude'})

Results2 = 1×2 table
    Frequency     Peak Magnitude
    __________    ______________

    5.8516e+05       0.22111    

figure
plot(Fv2, abs(FTs2)*2)
grid
ylim([min(ylim) 1.5*max(ylim)])
xlabel('Frequency')
ylabel('Magnitude')
title('FFT\_spatial')
text(Results2{:,1}, Results2{:,2}, compose('\\leftarrow Freq = %.1f\n    Mag = %.3f',Results2{:,[1  2]}), 'Rotation',45)

function [FTs1,Fv] = FFT1(s,t)
    s = s(:);
    t = t(:);
    L = numel(t);
    Fs = 1/mean(diff(t));
    Fn = Fs/2;
    NFFT = 2^nextpow2(L);
    FTs = fft((s - mean(s)).*hann(L), NFFT)/sum(hann(L));
    Fv = linspace(0, 1, NFFT/2+1)*Fn;
    Iv = 1:numel(Fv);
    FTs1 = FTs(Iv);
end

That seems appropriate to me, however if you want to identify more peaks, decrese the 'MinPeakProminence' value, or eliminate it altogether if desired. You can also use xlim to decrease the width of the plots, to show more detail (especially at the lower frequency values).

My apologies for not seeing ‘Fourier_data101.xlsx’ earlier. I was off doing other things for a while.

.

[University]

Thank you so mcuh for your assistance. Is possible to get the other half the frequencies (negative parts). See the attached figure to understand what I mean. Also, why do you have to resample the data.

[Star Strider]

As always, my pleasure!

Sure!

It looks like you allready solved that.

Here, I changed ‘FFT1’ to ‘FFT2’ to return two-sided Fourier transforms —

F = openfig('plot2.fig');

T1 = readtable('Fourier_data101.xlsx')

T1 = 101×2 table
    time     theta 
    ____    _______

       0    0.77481
    0.02    0.35591
    0.04    0.18156
    0.06    0.21909
    0.08    0.21391
     0.1    0.24668
    0.12    0.29114
    0.14    0.31926
    0.16     0.3504
    0.18    0.37767
     0.2    0.39259
    0.22    0.39536
    0.24    0.38582
    0.26    0.36519
    0.28    0.33557
     0.3    0.29983

t1 = T1{:,1};
s1 = T1{:,2};

T2 = readtable('FFT_spatial.xlsx')

T2 = 857×2 table
         x            theta   
    ___________    ___________

         -1e-05    -1.6396e-19
         -1e-05     5.2143e-07
    -9.9543e-06      0.0017766
    -9.9099e-06       0.003522
    -9.9087e-06      0.0035718
    -9.9074e-06      0.0036231
     -9.863e-06      0.0053857
    -9.8196e-06      0.0070574
     -9.817e-06      0.0071556
    -9.8144e-06      0.0072566
    -9.7711e-06      0.0088551
    -9.7292e-06       0.010345
    -9.7251e-06       0.010482
    -9.7212e-06       0.010617
    -9.6782e-06       0.012015
    -9.6357e-06       0.013347

t2 = T2{:,1};
s2 = T2{:,2};

Fs2 = 1/mean(diff(t2))

Fs2 = 4.2800e+07

[s2r,t2r] = resample(s2, t2, Fs2);                          % Resample To Constant Sampling Frequency

figure
plot(t2, s2, 'DisplayName','Original Signal')
hold on
plot(t2r, s2r, 'DisplayNAme','Resampled Signal')
hold off
grid
xlabel('Time')
ylabel('Amplitude')
title('FFT\_spatial — Time Domain Plot')
legend('Location','best')

[FTs1,Fv1] = FFT2(s1,t1);                                   % Fourier_data101


[pks1,locs1] =  findpeaks(abs(FTs1)*2, 'MinPeakProminence',0.005)

pks1 = 3×1

0.1120 0.0186 0.1120

locs1 = 3×1

59 65 71

Results1 = table(Fv1(locs1).', pks1, 'VariableNames',{'Frequency','Peak Magnitude'})

Results1 = 3×2 table
    Frequency    Peak Magnitude
    _________    ______________

     -2.3438         0.11201   
           0        0.018623   
      2.3438         0.11201   

figure
plot(Fv1, abs(FTs1)*2)
grid
ylim([min(ylim) 1.5*max(ylim)])
xlabel('Frequency')
ylabel('Magnitude')
title('Fourier\_data101')
xlim(xlim/3)
text(Results1{:,1}, Results1{:,2}, compose('\\leftarrow Freq = %.2f\n    Mag = %.3f',Results1{:,[1  2]}), 'Rotation',45)

[FTs2,Fv2] = FFT2(s2r,t2r);                                 % FFT_spatial

[pks2,locs2] =  findpeaks(abs(FTs2)*2, 'MinPeakProminence',0.01)

pks2 = 8×1

0.2116 0.2211 0.1216 0.1947 0.1947 0.1216 0.2211 0.2116

locs2 = 8×1

497 499 502 512 514 524 527 529

Results2 = table(Fv2(locs2).', pks2, 'VariableNames',{'Frequency','Peak Magnitude'})

Results2 = 8×2 table
     Frequency     Peak Magnitude
    ___________    ______________

    -6.6875e+05       0.21162    
    -5.8516e+05       0.22111    
    -4.5977e+05       0.12159    
         -41797       0.19473    
          41797       0.19473    
     4.5977e+05       0.12159    
     5.8516e+05       0.22111    
     6.6875e+05       0.21162    

figure
plot(Fv2, abs(FTs2)*2)
grid
ylim([min(ylim) 1.5*max(ylim)])
xlabel('Frequency')
ylabel('Magnitude')
title('FFT\_spatial')
xlim(xlim/12)
text(Results2{:,1}, Results2{:,2}, compose('\\leftarrow Freq = %.2f\n    Mag = %.3f',Results2{:,[1  2]}), 'Rotation',45)

function [FTs2,Fv] = FFT2(s,t)
    s = s(:);
    t = t(:);
    L = numel(t);
    Fs = 1/mean(diff(t));
    Fn = Fs/2;
    NFFT = 2^nextpow2(L);
    FTs = fftshift(fft((s - mean(s)).*hann(L), NFFT)/sum(hann(L)));
    Fv = linspace(-Fs/2, Fs/2-Fs/length(FTs), length(FTs));
    Iv = 1:numel(Fv);
    FTs2 = FTs(Iv);
end

Make approopriate changes to the rest of the code to get the result you want. (I do not recoognise the spectrum in ‘plot2.fig’ so I guess it was not a signal already uploaded here.)

EDIT — (26 Oct 2023 at 11:12)

Resmapling is necessary because the fft function requires evenly-sampled (or reasonably close to evenly-sampled) data. (The nufft function does not, however I prefer to use fft with resampled data, since resampling also permits other signal processing techniques such as filtering, that are otherwise impossible.) Looking at ‘t2’, the mean of the differences (‘diff(t)’) is 2.3364e-8 and the standard deviation of the differences is 1.6398e-8. It should be several orders-of-magnitude less than the mean. The sampling time differences range from 2.1069e-12 to 7.0182e-8. That range is simply too much for fft (or other signal processing functions, all of which require regular sampling) to deal with. Without resampling, the results would be inaccurate. Resampling to a common sampling frequency (chosen here as the inverse of the mean of the differences (using the median would also work), makes the subsequent analyses reliable. The time-domain plots would look substantially the same when plotted (added above).

.

[University]

Thank you so much for your assistance. Please what is the function of "hann function" in Matlab?