MATLAB etude demonstrates spectrum for multitone signal with different amplitudes

unit_spectrum_et1.m

% Spectral analysis etudes
% This etude demonstrates spectrum for multitone signal with different
% amplitudes.
% code by [email protected]

clear

% PARAMETERS
% -------------------------------------------------------------------------

% Sample rate (Hz)
Fs = 48e3;
% Time resolution (secs)
dt = 1/Fs;
% Number of signal samples
NSIG = 2^10;
% Number of FFT input samples (zero padded signal samples)
NFFT = 2^14;
% Actual frequency resolution (Hz)
df = Fs/NSIG;
% Interpolated frequency resolution (Hz)
dfi = Fs/NFFT;

% Central frequency (Hz)
f0 = 100*df;
% Frequency difference between neighbor tones (Hz)
f1 = 20*df;
% Number of tone components in signal
nf = 5;
% Tones frequencies (Hz)
fsig = (f0-(nf-1)/2*f1):f1:(f0+nf/2*f1);
% Tones amplitudes 0-Pk (Volts)
%asig = ones(size(fsig)); % uniform
asig = (1:numel(fsig)) / numel(fsig); % stairs

% GENERATOR
% -------------------------------------------------------------------------

% Generate tone signals
s = zeros(nf, NSIG);
t = (0:NSIG-1)*dt;
for i = 1:nf
    s(i,1:NSIG) = asig(i)*sin(2*pi*fsig(i)*t);
end

% PROCESSOR
% -------------------------------------------------------------------------

% Calculate FFT magnitude
x = zeros(1,NFFT);
x(1:NSIG) = sum(s,1);
X = fft(x, NFFT);
X = X(2:NFFT/2); % X(1) is DC
X = abs(X);

% Multiply by 2 for 0-Pk amplitude
X = X/NSIG*2;

% RESULTS
% -------------------------------------------------------------------------

% Interpolated spectrum frequnecies (Hz)
f = (1:NFFT/2-1)*dfi;
% convert Hz to kHz
f_khz = f*1e-3;
fsig_khz = fsig*1e-3;

% plot spectrum
figure; stem(f_khz, X, '.');
fb = max(fsig_khz) - min(fsig_khz);
fc = mean(fsig_khz);
xlim([fc-fb fc+fb]);
ylim([0 1.1*max(asig)]);
xlabel('kHz'); ylabel('V 0-Pk');
grid on; grid minor;