Remove 50 Hz noise
38 次查看(过去 30 天)
显示 更早的评论
I am recording 6 channels from a force transducer, but there is a 50 Hz noise I want to remove. I was considering to use a simple solution as iirnotch, but I don't have this function on my Matlab R2015a. How can I create a filter (filtfilt and fir1?) to eliminate this 50 Hz without a signal processing toolbox?
Thank you
16 个评论
Rik
2018-6-25
I added the 'no toolbox' tag, because this is almost trivial otherwise. My Google search (for 'bandstop filter') didn't turn up something either. There are several contributors more skilled in signal processing than me, so they might be able to help here.
Paola
2018-6-26
编辑:Walter Roberson
2018-6-26
20000 Hz. Precisely my fs is:
fs = Fdata.streams.Wav1.fs that is defined my acquisition system.
I have already tried usign a butterworth filter, but of course, it doesn't help so much.
ForceData = double(ForceData);
[B,A]=butter(2,50/fs/2);
for ch=1:6
ForceData(ch,:) = filtfilt(B, A, ForceData(ch,:));
end
Paola
2018-6-26
编辑:Paola
2018-6-26
What am I doing wrong?
According to the "cut off frequency" paragraph: "If Wn is scalar, then butter designs a lowpass or highpass filter with cutoff frequency Wn." So apparently for a bandstop I should use a two element vector. In this case what value should I put, to remove only the 50 Hz? And how do I choose the order?
ForceData = double(ForceData);
[B,A]=butter(2,[50]/(fs/2),'stop');
for ch=1:6
ForceData(ch,:) = filtfilt(B, A,ForceData(ch,:));
end
Thank you
Walter Roberson
2018-6-26
You cannot remove only 50 Hz exactly with a digital filter of that type. Instead you effectively construct a lowpass filter that starts to drop down near 50 Hz, cutting off high frequencies, "plus" (in some sense) a highpass filter that starts coming up near 50 Hz, cutting off low frequencies. The downwards dip heading into 50 Hz and the upwards dip heading after 50 Hz together diminish 50 Hz, along with diminishing frequencies approaching 50 Hz on other side.
The higher the order that you use, the sharper the cutoff can be, but also the longer the filter representation comes out when expressed as a convolution, so the processing rate is higher for lower orders.
Paola
2018-6-26
编辑:Paola
2018-6-26
Thank you for your clear explanation.
ex:
ForceData = double(ForceData);
[B,A]=butter(2,[49 51]/(fs/2),'stop');
for ch=1:6
ForceData(ch,:) = filtfilt(B, A,ForceData(ch,:));
end
Is it correct in this way or should I express the cutoff frequencies in radians? I will try different orders to check what could be the best for me.
Thanks again
Walter Roberson
2018-6-26
I think that looks okay, but I hardly do any filter work so I am not sure.
Paola
2018-6-27
Hi, I tried this code but it doesn't change much. The ideal would be a notch but I don't have DSP toolbox and I am not able to calculate it by myself. Any suggestion to talk to a filter expert?
Walter Roberson
2018-6-27
Tests I did with butter suggest that it is not very good at narrow band suppression :( With low orders it still left a fair bit of signal in the target range, and with orders from about 5 up, it created quite large transients.
Paola
2018-6-28
Definitely! I am trying to write down a notch but I find it quite difficult since I am not an expert...
Star Strider
2018-6-28
I thought you said that you do not have the Signal Processing Toolbox! I agree with Rik Wisselink that with it this is essentially trivial. However, you appear to be using its functions, so I will add this:
Fs = 20000;
Fn = Fs/2;
Wp = [47 53]/Fn;
Ws = [49 51]/Fn;
Rp = 1;
Rs = 50;
[n,Wn] = buttord(Wp,Ws,Rp,Rs);
[z,p,k] = butter(n,Wn,'stop');
[sos,g] = zp2sos(z,p,k);
figure
freqz(sos, 2^14, Fs)
set(subplot(2,1,1), 'XLim',[0 100])
set(subplot(2,1,2), 'XLim',[0 100])
The freqz call lets you see what the filter is doing in the frequency domain. There is also no need for the loop, because the filter functions operate column-wise. You simply need to transpose your matrix.
This should work:
ForceDataFilt = filtfilt(sos, g, ForceData'); % Transpose & Filter
ForceDataFilt = ForceDataFilt';
I transposed ‘ForceDataFilt’ so it will match the rest of your code. Experiment with the filter passband and stopband frequencies to get the result you want. Remember that for a notch or bandstop filter, the stopband frequencies (here ‘Ws’) are within the passband limits, ‘Wp’. The passband and stopband attenuation values are ‘Rp’ and ‘Rs’ respectively.
Paola
2018-6-28
Actually, I don't have the DSP toolbox and I have an old version of the Signal Processing toolbox.
回答(1 个)
Kouichi C. Nakamura
2020-2-20
编辑:Walter Roberson
2020-2-21
This looks like an answer.
Remove the 60 Hz Hum from a Signal
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Digital Filter Analysis 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!发生错误
由于页面发生更改,无法完成操作。请重新加载页面以查看其更新后的状态。
您也可以从以下列表中选择网站:
如何获得最佳网站性能
选择中国网站(中文或英文)以获得最佳网站性能。其他 MathWorks 国家/地区网站并未针对您所在位置的访问进行优化。
美洲
- América Latina (Español)
- Canada (English)
- United States (English)
欧洲
- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)
- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- United Kingdom(English)
亚太
- Australia (English)
- India (English)
- New Zealand (English)
- 中国
- 日本Japanese (日本語)
- 한국Korean (한국어)