I've made a simple first order IIR highpass filter with a zero at z = 1 and a pole at z = 0.9. Its frequency response looks like this:
Now, I filter a DC signal using this filter. Here's the MATLAB code I use to do it:
b = [1 -1]; % Zero at z = 1
a = [1 -0.9]; %Pole at z = 0.9
figure(1)
freqz(b, a)
t = 1:100;
x(1:length(t)) = 1; % Constant function
y = filter(b, a, x);
figure(2)
plot(t, x)
xlabel('Time');
ylabel('Input Signal');
figure(3)
plot(t, y)
xlabel('Time');
ylabel('Output Signal');
As my filter is highpass, I expect the DC to become zero, or atleast become severely attenuated. However, the output signal I get looks like this:
From my understanding, this exponential output is a transient produced because I haven't set the initial conditions correctly. Sure enough, setting x[-1] = 1 solves the problem. However, this works only for this particular input DC signal. For any general input signal, how do I set the initial conditions so that transients aren't produced?
Edit : I'm aware that the filtfilt() function does forward-backward filtering with transient minimization, but I really want to port the filter to an embedded platform, so I need to understand how transient removal works. Thanks in advance for the help!
Edit 2 : As suggested by Kuba Ober, I tried setting x[-1] as the value that it actually should have been. It works fine for a DC input, but here's what happened for a sinusoidal input:
clc; clear all;
p = 0.9
a = [1 -p]
b = [1 -1]
n = 1:100; % Samples
f = 0.2; % Frequency in Hz
Fs = 10; % Sampling rate in samples per second
t = n/Fs; % Time axis
x = sin(2*pi*f*t);
% Filter with the appropriate initial conditions
y = filter(b, a, x, filtic(b, a, [], [sin(2*pi*f*0)]));
figure(1)
plot(t, x)
xlabel('Time');
ylabel('Input Signal');
figure(2)
plot(t, y)
xlabel('Time');
ylabel('Output Signal');
Here's the input signal :
And here's the output :
The first peak is visibly smaller than the second, which indicates some transients being present. I'm not entirely sure about this, but I think the reason it doesn't work is because just setting x[-1] is not enough, I also need to set y[-1]. The problem here is that there's no way to find out what y[-1] actually should be.
Edit 3 : Let me provide a little more info on the problem I'm working on. I'm trying to use filters to remove noise from ECG (Electrocardiogram) signals in an embedded platform. Here's a typical ECG signal, after filtering:
Here's what an ECG signal looks before filtering:
Note the DC offset in the signal before filtering. For filtering, I need a notch filter to remove high frequency power line noise and a highpass filter to remove the DC and the low frequency "drifting" of the signal.
The filters I use need to be linear phase, since the time domain morphology of an ECG signal is very important for diagnosis. However, my filter doesn't need to be real-time, as I'm doing the processing offline after acquiring the ECG signal from the patient. So, for implementing nonlinear phase IIR filters, I'm currently using forward-backward zero phase filtering.
One opinion that's shared by @Matt L. and @Royi is that transients are unavoidable in real-time filtering and that I should use a longer input signal and crop off the first few seconds of the filtered output instead. This is something I'd like to avoid, as acquiring long ECGs from a living patient is somewhat difficult. Also, I do not have to filter in real-time, so any technique of transient removal that hinges on knowing the entire signal in advance is perfectly all right. Any help is appreciated!
x[-1]
to whatever value it actually had? As in - take one more sample before the samples that you'll filter, and use its value? $\endgroup$filter(b, a, x, -1)
. That is absolutely correct, but it's not the same as setting x[-1] = -1. The fourth argument is actually the initial conditions for the delays in the direct-form II implementation. To actually set x[-1] = -1, you need to usefilter(b, a, x, filtic(b, a, [], [-1])
. This doesn't quite work, though - only setting x[-1] = 1 works for removing a DC with amplitude 1. link $\endgroup$