DC gain for FIR high-pass filter using fir1

28 次查看(过去 30 天)
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

采纳的回答

Paul
Paul 2022-11-18
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
  3 个评论
Paul
Paul 2022-11-18
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
Masao 2022-11-20
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!

请先登录,再进行评论。

更多回答(0 个)

产品


版本

R2022b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by