Sunday 29 March 2015

Design Matlab filter with cutoff frequency

Design a simple set of filter coefficients, apply it to a test signal:
(http://developer.mbed.org/handbook/Matlab-FIR-Filter)

320 samples at 48kHz for a sum of two sinusoids (1000Hz + 15000 Hz)

sample_rate = 48000;
nsamples = 256;

F = [1 15] * 1000;
A = [1 0.5];

% Time vector - use colon operator to generate integer vector of sample
% numbers
t = (0:nsamples-1) / sample_rate;

% Test signal - use matrix notation to compose it with single expression
signal = A * sin(2*pi*F'*t);



FIR filter

Assume a lowpass filter with cutoff frequency of 6 kHz.
The expectation is this should filter out the 15 kHz component from the test signal

% Choose filter cutoff frequency (6 kHz)
cutoff_hz = 6000;

% Normalize cutoff frequency (wrt Nyquist frequency)
nyq_freq = sample_rate / 2;
cutoff_norm = cutoff_hz / nyq_freq;

% FIR filter order (i.e. number of coefficients - 1)
order = 28;

% Create lowpass FIR filter through a direct approach: provide
% (normalized) cutoff frequency and filter order (assumed as known).
% fir1 takes care of designing the filter by imposing the constraints in
% the frequency domain and transforming back to time using a given window
% (the dafault used here is the Hamming window).
% For more advanced requirements see e.g. firpmord and firpm
% NOTE: fir1, firpmord and firpm all require Signal Processing Toolbox
fir_coeff = fir1(order, cutoff_norm);

% Analyse the filter using the Filter Visualization Tool
fvtool(fir_coeff, 'Fs', sample_rate)

% Filter the signal with the FIR filter
filtered_signal = filter(fir_coeff, 1, signal);


Plot the Signals
Also align filtered signal with original and discard transient samples (in this case the first order samples)


% Group delay as a scalar, in number of samples
group_delay = median(grpdelay(fir_coeff));

% Group delay in seconds
group_delay_s = group_delay / sample_rate;


% Plot the original signal...
figure(1)
plot(t, signal)
% ...and allow adding more plots
hold on

% Align and plot the filtered signal
% (On top of existing one)
plot(t-group_delay_s, filtered_signal, 'r-')

% Align and plot only the d6esired part of the filtered signal (discarding
% the transient)
plot(t(order:end)-group_delay_s, filtered_signal(order:end), ...
    'g', 'LineWidth', 4);

grid on
hold off

No comments:

Post a Comment