Signal processing

rng(1)
fs = 500;                          % sample rate (Hz)
tt = 0:1/fs:1;
s = sin(2*pi*12*tt) + 0.6*sin(2*pi*47*tt) + 0.4*randn(size(tt));

subplot(2, 1, 1)
plot(tt, s)
title('noisy signal'), xlabel('time (s)')

subplot(2, 1, 2)
[f, A] = spectrum(s, fs);
plot(f, A, 'LineWidth', 1.5)
xlim([0 100])
title('amplitude spectrum'), xlabel('frequency (Hz)')

% A 9 Hz sine sampled at 10 Hz looks like a 1 Hz sine
tt = linspace(0, 1, 1000);
f0 = 9; fs = 10;
ts = 0:1/fs:1;
plot(tt, sin(2*pi*f0*tt), 'Color', [0.75 0.75 0.75]), hold on
stem(ts, sin(2*pi*f0*ts))
plot(tt, sin(2*pi*(f0 - fs)*tt), 'r--', 'LineWidth', 2)
hold off
legend('9 Hz signal', 'samples at 10 Hz', 'apparent signal: 1 Hz')
xlabel('time (s)')
title('Aliasing: sampling below the Nyquist rate (2 × 9 = 18 Hz)')

% A strong tone leaks across the spectrum and hides a weak one — a window helps
fs = 1000; n = 256; tt = (0:n-1)/fs;
x = sin(2*pi*101.3*tt) + 0.01*sin(2*pi*160*tt);
w = 0.5 - 0.5*cos(2*pi*(0:n-1)/(n-1));         % Hann window
[f, A1] = spectrum(x, fs);
[f, A2] = spectrum(2 * x .* w, fs);            % ×2 makes up for the window's average of ½
semilogy(f, A1 + 1e-7, 'LineWidth', 1.5), hold on
semilogy(f, A2 + 1e-7, 'LineWidth', 1.5), hold off
xlim([0 300]), ylim([1e-5 2])
xline(160, 'k:', 'weak 160 Hz tone')
legend('rectangular (no window)', 'Hann window')
xlabel('frequency (Hz)'), ylabel('amplitude')
title('Spectral leakage')