MATLAB Answers

# Why the amplitude of theFFT trasform is not equal as the amplitude of the input signal?

5 views (last 30 days)
Arianna Senesi on 26 Apr 2021
Answered: Mathieu NOE on 27 Apr 2021
close all
clear all
clc
sine = [];
rate = 100000;
fbase=1300;
freq = fbase+100;
sinesine1 = 0.25*sin(linspace(0, 2*pi,rate/freq));
freq = fbase+200;
sinesine2 = 0.25*sin(linspace(0,2*pi,rate/freq));
freq = fbase+300;
sinesine3 =0.25*sin(linspace(0, 2*pi,rate/freq));
freq = fbase+400;
sinesine4 = 0.25*sin(linspace(0, 2*pi,rate/freq));
freq =fbase+500;
sinesine5 = 0.25*sin(linspace(0, 2*pi,rate/freq));
n=200;
for i = 1: n
sine = [sine sinesine1(1:end-1)];
end
for i = 1: n
sine = [sine sinesine2(1:end-1)];
end
for i = 1: n
sine = [sine sinesine3(1:end-1)];
end
for i = 1: n
sine = [sine sinesine4(1:end-1)];
end
for i = 1: n
sine = [sine sinesine5(1:end-1)];
end
Fs = rate;
L = length(sine); % Length of signal
T = 1/Fs; % Sampling period
t = (0:L-1)*T; % Time vector
Z = fft(sine);
P2= abs(Z)./L;
P1 = P2(1:L/2+1);
P1(2:end-1) = 2*P1(2:end-1);
f = Fs*(0:(L/2))/L;
figure('Name','Spettro 1400-1800','NumberTitle','off');
plot(f,P1)
axis([fbase fbase+1100 0 max(P5(2:end))+0.01*max(P5(2:end))])
##### 4 CommentsShowHide 3 older comments
Arianna Senesi on 26 Apr 2021
Where does it give you errors? because it works for me. thanks anyway

Sign in to comment.

### Accepted Answer

Mathieu NOE on 27 Apr 2021
hi again
I was not happy to see that your code needed a very high sampling rate to display correctly the 5 frequencies
this is by no mean required by Shannon's theorem
I modified the way the data are generated and now we can have the right display even with Fs = 5000;
also one single for loop needed vs five
close all
clearvars
clc
sine = [];
% rate = 100000;
% fbase= 1300;
% freq = fbase+100;
% sinesine1 = 0.25*sin(linspace(0, 2*pi,rate/freq));
% freq = fbase+200;
% sinesine2 = 0.25*sin(linspace(0,2*pi,rate/freq));
% freq = fbase+300;
% sinesine3 =0.25*sin(linspace(0, 2*pi,rate/freq));
% freq = fbase+400;
% sinesine4 = 0.25*sin(linspace(0, 2*pi,rate/freq));
% freq =fbase+500;
% sinesine5 = 0.25*sin(linspace(0, 2*pi,rate/freq));
% n=200;
% for i = 1: n
% sine = [sine sinesine1(1:end-1)];
% end
% for i = 1: n
% sine = [sine sinesine2(1:end-1)];
% end
% for i = 1: n
% sine = [sine sinesine3(1:end-1)];
% end
% for i = 1: n
% sine = [sine sinesine4(1:end-1)];
% end
% for i = 1: n
% sine = [sine sinesine5(1:end-1)];
% end
% Fs = rate;
% my suggestion
periods = 200;
Fs = 5e3;
dt = 1/Fs;
fbase= 1300;
freq = fbase + [100 200 300 400 500];
for ci = 1:length(freq)
omega = 2*pi*freq(ci);
t = (0:dt:periods/freq(ci)-dt);
sine = [sine 0.25*sin(omega*t)];
end
% T = 1/Fs; % Sampling period
% t = (0:L-1)*T; % Time vector
% L = length(sine); % Length of signal
% Z = fft(sine);
% P2= abs(Z)./L;
% P1 = P2(1:L/2+1);
% P1(2:end-1) = 2*P1(2:end-1);
% f = Fs*(0:(L/2))/L;
% figure('Name','Spettro 1400-1800','NumberTitle','off');
% plot(f,P1)
% axis([fbase fbase+1100 0 max(P5(2:end))+0.01*max(P5(2:end))])
%%%%%%%%%% PEAK HOLD AVERAGING FFT SPECTRUM DEMO %%%%%%%%%%%%%%%%%
NFFT = 250;
Overlap = 0.75;
[freq,fft_spectrum] = myfft_peak(sine(:), Fs, NFFT, Overlap); % data must be column oriented (samples x channels)
figure(1),plot(freq,fft_spectrum,'b');grid
title(['Peak Hold FFT Spectrum / Fs = ' num2str(Fs) ' Hz / Delta f = ' num2str(freq(2)-freq(1)) ' Hz ']);
xlabel('Frequency (Hz)');ylabel(' Amplitude')
%%%%%%%%%% SPECTROGRAM DEMO %%%%%%%%%%%%%%%%%
[S,F,T] = myspecgram(sine, Fs, NFFT, Overlap);
figure(2);
imagesc(T,F,(abs(S)))
colorbar('vert');
set(gca,'YDir','Normal')
xlabel('Time (secs)')
ylabel('Freq (Hz)')
title('Short-time Fourier Transform spectrum')
colormap('jet');
function [freq_vector,fft_spectrum] = myfft_peak(signal, Fs, nfft, Overlap)
% FFT peak spectrum of signal (example sinus amplitude 1 = 0 dB after fft).
% Linear averaging
% signal - input signal,
% Fs - Sampling frequency (Hz).
% nfft - FFT window size
% Overlap - buffer percentage of overlap % (between 0 and 0.95)
[samples,channels] = size(signal);
if channels > samples
signal = signal'; % flip dimensions
end
% fill signal with zeros if its length is lower than nfft
if samples<nfft
s_tmp = zeros(nfft,channels);
s_tmp((1:samples),:) = signal;
signal = s_tmp;
samples = nfft;
end
% window : hanning
window = hanning(nfft);
window = window(:);
% compute fft with overlap
offset = fix((1-Overlap)*nfft);
spectnum = 1+ fix((samples-nfft)/offset); % Number of windows
% % for info is equivalent to :
% noverlap = Overlap*nfft;
% spectnum = fix((samples-noverlap)/(nfft-noverlap)); % Number of windows
% main loop
fft_spectrum = 0;
for i=1:spectnum
start = (i-1)*offset;
sw = signal((1+start):(start+nfft),:).*(window*ones(1,channels));
fft_spectrum = max(fft_spectrum,(abs(fft(sw))*4/nfft)); % X=fft(x.*hanning(N))*4/N; % hanning only
end
% fft_spectrum = fft_spectrum/spectnum; % to do linear averaging scaling
% one sidded fft spectrum % Select first half
if rem(nfft,2) % nfft odd
select = (1:(nfft+1)/2)';
else
select = (1:nfft/2+1)';
end
fft_spectrum = fft_spectrum(select,:);
freq_vector = (select - 1)*Fs/nfft;
end
function [fft_specgram,freq_vector,time] = myspecgram(signal, Fs, nfft, Overlap)
% FFT peak spectrogram of signal (example sinus amplitude 1 = 0 dB after fft).
% signal - input signal,
% Fs - Sampling frequency (Hz).
% nfft - FFT window size
% Overlap - buffer overlap % (between 0 and 0.95)
signal = signal(:);
samples = length(signal);
% fill signal with zeros if its length is lower than nfft
if samples<nfft
s_tmp = zeros(nfft,1);
s_tmp((1:samples),:) = signal;
signal = s_tmp;
samples = nfft;
end
% window : hanning
window = hanning(nfft);
window = window(:);
% compute fft with overlap
offset = fix((1-Overlap)*nfft);
spectnum = 1+ fix((samples-nfft)/offset); % Number of windows
% % for info is equivalent to :
% noverlap = Overlap*nfft;
% spectnum = fix((samples-noverlap)/(nfft-noverlap)); % Number of windows
% main loop
fft_specgram = [];
for ci=1:spectnum
start = (ci-1)*offset;
sw = signal((1+start):(start+nfft)).*window;
fft_specgram = [fft_specgram abs(fft(sw))*4/nfft]; % X=fft(x.*hanning(N))*4/N; % hanning only
end
% one sidded fft spectrum % Select first half
if rem(nfft,2) % nfft odd
select = (1:(nfft+1)/2)';
else
select = (1:nfft/2+1)';
end
fft_specgram = fft_specgram(select,:);
freq_vector = (select - 1)*Fs/nfft;
% time vector
% time stamps are defined in the middle of the buffer
time = ((0:spectnum-1)*offset + round(offset/2))/Fs;
end
##### 0 CommentsShowHide -1 older comments

Sign in to comment.

### More Answers (1)

Mathieu NOE on 26 Apr 2021
hello
FYI see code below
it iplement a peak hold fft averaging technique; also the spectrogram principle is coded (so you can see what matlab is doing)
both methods show the right amplitude
I reduced your sampling freq that was too high for the demo purpose
close all
clearvars
clc
sine = [];
rate = 10000;
fbase= 1300;
freq = fbase+100;
sinesine1 = 0.25*sin(linspace(0, 2*pi,rate/freq));
freq = fbase+200;
sinesine2 = 0.25*sin(linspace(0,2*pi,rate/freq));
freq = fbase+300;
sinesine3 =0.25*sin(linspace(0, 2*pi,rate/freq));
freq = fbase+400;
sinesine4 = 0.25*sin(linspace(0, 2*pi,rate/freq));
freq =fbase+500;
sinesine5 = 0.25*sin(linspace(0, 2*pi,rate/freq));
n=200;
for i = 1: n
sine = [sine sinesine1(1:end-1)];
end
for i = 1: n
sine = [sine sinesine2(1:end-1)];
end
for i = 1: n
sine = [sine sinesine3(1:end-1)];
end
for i = 1: n
sine = [sine sinesine4(1:end-1)];
end
for i = 1: n
sine = [sine sinesine5(1:end-1)];
end
Fs = rate;
% T = 1/Fs; % Sampling period
% t = (0:L-1)*T; % Time vector
% L = length(sine); % Length of signal
% Z = fft(sine);
% P2= abs(Z)./L;
% P1 = P2(1:L/2+1);
% P1(2:end-1) = 2*P1(2:end-1);
% f = Fs*(0:(L/2))/L;
% figure('Name','Spettro 1400-1800','NumberTitle','off');
% plot(f,P1)
% axis([fbase fbase+1100 0 max(P5(2:end))+0.01*max(P5(2:end))])
%%%%%%%%%% PEAK HOLD AVERAGING FFT SPECTRUM DEMO %%%%%%%%%%%%%%%%%
NFFT = 500;
Overlap = 0.75;
[freq,fft_spectrum] = myfft_peak(sine(:), Fs, NFFT, Overlap); % data must be column oriented (samples x channels)
figure(1),plot(freq,fft_spectrum,'b');grid
title(['Peak Hold FFT Spectrum / Fs = ' num2str(Fs) ' Hz / Delta f = ' num2str(freq(2)-freq(1)) ' Hz ']);
xlabel('Frequency (Hz)');ylabel(' Amplitude')
%%%%%%%%%% SPECTROGRAM DEMO %%%%%%%%%%%%%%%%%
[S,F,T] = myspecgram(sine, Fs, NFFT, Overlap);
figure(2);
imagesc(T,F,(abs(S)))
colorbar('vert');
set(gca,'YDir','Normal')
xlabel('Time (secs)')
ylabel('Freq (Hz)')
title('Short-time Fourier Transform spectrum')
colormap('jet');
function [freq_vector,fft_spectrum] = myfft_peak(signal, Fs, nfft, Overlap)
% FFT peak spectrum of signal (example sinus amplitude 1 = 0 dB after fft).
% Linear averaging
% signal - input signal,
% Fs - Sampling frequency (Hz).
% nfft - FFT window size
% Overlap - buffer percentage of overlap % (between 0 and 0.95)
[samples,channels] = size(signal);
if channels > samples
signal = signal'; % flip dimensions
end
% fill signal with zeros if its length is lower than nfft
if samples<nfft
s_tmp = zeros(nfft,channels);
s_tmp((1:samples),:) = signal;
signal = s_tmp;
samples = nfft;
end
% window : hanning
window = hanning(nfft);
window = window(:);
% compute fft with overlap
offset = fix((1-Overlap)*nfft);
spectnum = 1+ fix((samples-nfft)/offset); % Number of windows
% % for info is equivalent to :
% noverlap = Overlap*nfft;
% spectnum = fix((samples-noverlap)/(nfft-noverlap)); % Number of windows
% main loop
fft_spectrum = 0;
for i=1:spectnum
start = (i-1)*offset;
sw = signal((1+start):(start+nfft),:).*(window*ones(1,channels));
fft_spectrum = max(fft_spectrum,(abs(fft(sw))*4/nfft)); % X=fft(x.*hanning(N))*4/N; % hanning only
end
% fft_spectrum = fft_spectrum/spectnum; % to do linear averaging scaling
% one sidded fft spectrum % Select first half
if rem(nfft,2) % nfft odd
select = (1:(nfft+1)/2)';
else
select = (1:nfft/2+1)';
end
fft_spectrum = fft_spectrum(select,:);
freq_vector = (select - 1)*Fs/nfft;
end
function [fft_specgram,freq_vector,time] = myspecgram(signal, Fs, nfft, Overlap)
% FFT peak spectrogram of signal (example sinus amplitude 1 = 0 dB after fft).
% signal - input signal,
% Fs - Sampling frequency (Hz).
% nfft - FFT window size
% Overlap - buffer overlap % (between 0 and 0.95)
signal = signal(:);
samples = length(signal);
% fill signal with zeros if its length is lower than nfft
if samples<nfft
s_tmp = zeros(nfft,1);
s_tmp((1:samples),:) = signal;
signal = s_tmp;
samples = nfft;
end
% window : hanning
window = hanning(nfft);
window = window(:);
% compute fft with overlap
offset = fix((1-Overlap)*nfft);
spectnum = 1+ fix((samples-nfft)/offset); % Number of windows
% % for info is equivalent to :
% noverlap = Overlap*nfft;
% spectnum = fix((samples-noverlap)/(nfft-noverlap)); % Number of windows
% main loop
fft_specgram = [];
for ci=1:spectnum
start = (ci-1)*offset;
sw = signal((1+start):(start+nfft)).*window;
fft_specgram = [fft_specgram abs(fft(sw))*4/nfft]; % X=fft(x.*hanning(N))*4/N; % hanning only
end
% one sidded fft spectrum % Select first half
if rem(nfft,2) % nfft odd
select = (1:(nfft+1)/2)';
else
select = (1:nfft/2+1)';
end
fft_specgram = fft_specgram(select,:);
freq_vector = (select - 1)*Fs/nfft;
% time vector
% time stamps are defined in the middle of the buffer
time = ((0:spectnum-1)*offset + round(offset/2))/Fs;
end
##### 2 CommentsShowHide 1 older comment
Mathieu NOE on 27 Apr 2021
Ok
so let's put back the rate at 100 kHz, then the 5 peaks will show up but you'll have to increase the NFFT to 10 times higher if you'd like to distinguish the relatively close frequencies

Sign in to comment.

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!