demo_filter_matrices
Illustration of low-pass and high-pass linear filters uing banded filter matrices.
Ivan Selesnick, selesi@nyu.edu, Tandon School of Engineering New York University 2017
Contents
Start
clear printme = @(filename) print('-dpdf', sprintf('figures/demo_filter_matrices_%s', filename));
Define high-pass filter H
N = 100; % N : signal length fc = 0.05; % fc : cut-off frequency (cycles/sample) (0 < fc < 0.5) d = 2; % d : degree of filter is 2d (d = 1, 2, 3) [A, B, B1, D, a, b] = ABfilt(d, fc, N); H = @(x) A\(B*x); % H : high-pass filter L = @(x) x - H(x); % L : low-pass filter
Impulse responses of filters
no = 40; x = zeros(N,1); x(no) = 1; % input signal (impulse) figure(1) clf subplot(2, 1, 1) stem(1:N, L(x), '.') title(sprintf('LPF impulse response (impulse at n = %d)', no)) subplot(2, 1, 2) stem(1:N, H(x), '.') title(sprintf('HPF impulse response (impulse at n = %d)', no)) printme('impulse_responses')
Frequency response of filters
% difference equation coefficients of low-pass filter b_lpf = a - b; a_lpf = a; % Compute frequency response [Hf_hpf, om] = freqz(b, a); [Hf_lpf, om] = freqz(b_lpf, a_lpf); figure(1) clf plot(om/pi/2, abs(Hf_lpf), om/pi/2, abs(Hf_hpf), '--') line(fc, 0.5, 'marker', 'o') title('Frequency responses') legend('Low-pass filter', 'High-pass filter') xlabel('Normalized frequency') ylim([0 1.2]) printme('frequency_responses')
Pole-zero diagram of low-pass filter L
figure(1) clf [hz, hp, hl] = zplane(b_lpf, a_lpf); title('LPF pole-zero diagram') printme('zplane_LPF')
Pole-zero diagram of high-pass filter H
figure(1) clf zplane(b, a); title('HPF pole-zero diagram') printme('zplane_HPF')
Filter noisy data
% Create noisy data n = (1:N)'; x = cos(0.05*pi*n) + 0.2*randn(N,1); % x : noisy data y_lpf = L(x); % Apply LPF to noisy data y_hpf = H(x); % Apply HPF to noisy data figure(1) clf subplot(2, 1, 1) plot(n, x, '.-', n, y_lpf, '-') legend('Input', 'Output') title('Low-pass filtering (smoothing of data)') subplot(2, 1, 2) plot(n, x, '.-', n, y_hpf, '-') legend('Input', 'Output') title('High-pass filtering (removes smooth part)') printme('filtering')