DC gain for FIR high-pass filter using fir1

Dear signal processing experts,

I'm analysing data with huge low frequency artifacts like DC offset, decay, trend, and so on, so I'm trying to remove these slow component with FIR filter.

However, if I use fir1 in matlab, the low frequency component, especially the DC component is not removed completely. Actually, when I check the cofficients:

>> b = fir1(500,5/750,'high');
>> sum(b)
ans =
0.0039

So there's about 0.39% of the original DC component still remain in the filtered data, which is still very large.

This is quite in contrast with low-pass filter:

>> b = fir1(500,5/750,'low');
>> 1-sum(b)
ans =
0

In order to remove DC component by FIR filter, I thought upon two options:

1. (1) Apply low-pass filter and subtract the filtered signal from the original.
2. (2) Force filter coefficient to zero-mean:

>> b = fir1(500,5/750,'high');
>> b = b-mean(b);
>> sum(b)
ans =
   5.5511e-17     % this is small enough

My questions are:

• Are the methods (1) and (2) above acceptable?
• What is the best way to handle this (remove slow components, esp DC offset)?

Thank you in advance for your kind attention and advice.

Best,

Masao

Answer 1

[Paul]

Hi Masao,

Let's try the proposed methods.

High-pass fir1 filter

bhigh = fir1(500,5/750,'high');
[h0,w] = freqz(bhigh,1,0:.0001:.1);

Low pass fir1 filter

blow = fir1(500,5/750,'low');
[hlow,w] = freqz(blow,1,w);

Response for Method 2, subtract the mean from the high-pass coefficients

h2 = freqz(bhigh - mean(bhigh),1,w);

Method 1, response of filter that is equivalent to subtracting the low-pass filtered signal from the original signal, i..e, H1(z) = 1 - Hlow(z)

h1 = freqz([1-blow(1) -blow(2:end)],1,w);

plot(w,abs([h0 ; h1 ; h2]))
xlim([0 .1])
legend('h0','h1','h2')

Method 1 is not desirable.

Replotting h0 and h2 in dB and zoming in a bit more

plot(w,mag2db(abs([h0 ; inf*h1 ; h2])))
xlim([0 .02])
legend('h0','h1','h2')

Should probably explore more across the entire frequency band from 0 - pi and compare the phase resposnes as well, if phase matters in your application.

Having said that, if you know that your signal has mean, linear, etc. components, using detrend on the data before filtering might be useful

[Masao]

Thank you, now I see how to use freqz.

bhigh = fir1(500,5/750,'high');
[h0,w] = freqz(bhigh,1,0:.0001:.1);
h2 = freqz(bhigh - mean(bhigh),1,w);
plot(w,angle(h0),'+',w,angle(h2),'-')
xlim([0,0.02])
legend('h0','h2')

It seems phase distortion is similar with and without subtraction of the mean, so is it reasonably acceptable to use the method (2)?

The slow artifact has some fluctuation, so subtraction of segmental liniar fit may be one way but it introduced high frequency artifacts at the segment ends and was not suitable.

[Paul]

Whether or not method 2 is acceptable is up to you, based on your problem and requirements. I'd suggest carefully comparing the response of h0 and h2 across the entire frequency range of interest.

[Masao]

Thanks, I checked the higher freqency range:

bhigh = fir1(500,5/750,'high');
[h0,w] = freqz(bhigh,1,0:.0001:pi);
h2 = freqz(bhigh - mean(bhigh),1,w);
plot(w,abs([h0 ;h2]))
xlim([0.05 1])
legend('h0','h2')

Ripple amplitude is slightly smaller with h2 than h0, this is a good behavior!

plot(w,angle(h0)-angle(h2))

And the difference is phase angle is within acceptable range. Thank you!