Matlab’s built-in function for Butterworth filtering has stability issues.
The following code snippet yields an unstable filter. The cutoff frequency is 4 Hz. The sample rate is 20000 Hz. The filter is a high pass, sixth order.
fh=4;
sr=20000;
order=6;
Wn=fh/(sr/2);
[bb,aa] = butter(order,Wn,’high’);
flag = isstable(bb,aa);
if(flag==1)
fprintf(‘\n filter is stable \n’);
else
fprintf(‘\n filter is unstable \n’);
end
(Note that the single quotation marks need to be fixed if you copy and paste this code into Matlab).
The same filtering operation can be performed in a stable manner using the following functions with a staged cascade implementation. Butterworth_filter.zip (right mouse click and save link as)
Digital Filtering, Butterworth, Main Function
Butterworth_filter_function.m
Supporting Functions
th_weighted_filter_coefficients.m
th_weighted_apply_filter.m
References:
An Introduction to the Filtering of Digital Signals: filter.pdf
– Tom Irvine