remove DC component of real time signal

103 次查看(过去 30 天)
I am using a interface for acquiring a signla in matlab.
The signal is sampled at 500Hz. I have to remove the DC component in real time.
How can I do that?
I know that simulink has the DC blocker block, but I need to accomplish this inside a matlab routine!

回答(1 个)

Star Strider
Star Strider 2023-11-14
The D-C component is the mean of the signal, so for a recorded signal, simply subtract the mean.
Otherwise, since the D-C offsset is defined as having a frequency of zero (rad/s, Hz, or anything else), a highpass filter (or bandpass filter) with its passband (or lower passband) slightly above zero should remove the D-C offset. The resulting signal should then have a zero D-C offset.
  5 个评论
L
L 2023-11-14
编辑:L 2023-11-14
sampling is 500Hz (2ms). Matlab receives a vector of 10 points at each 20ms. I can save up to 30 samples, if I want.
My signal of interest is actually between 0.25 and 4 hz.
So I could band pass it.. but if I could just remove the DC will also help.
What is important is to keep the signal without any introduced delays.
Star Strider
Star Strider 2023-11-14
This is a difficult filter to realise because the lower transition region is too steep. I can make it work with a lower attenuation (35 dB rather than 50 dB) for ‘Rs’ and then it distorts the filtered signal significantly less. I encourage you to experiment with various passband, stopband, and atttenuation values to get the result you want.
Fs = 500; % Sampling Frequency (Hz)
Fn = Fs/2; % Nyquist Frequency (Hz)
Wp = [0.25 4.0]/Fn; % Passband Frequency (Normalised)
Ws = [0.01 5.0]/Fn; % Stopband Frequency (Normalised)
Rp = 1; % Passband Ripple
Rs = 35; % Passband Ripple (Attenuation)
[n,Wp] = ellipord(Wp,Ws,Rp,Rs); % Elliptic Order Calculation
[z,p,k] = ellip(n,Rp,Rs,Wp); % Elliptic Filter Design: Zero-Pole-Gain
[sos,g] = zp2sos(z,p,k); % Second-Order Section For Stability
figure
freqz(sos, 2^16, Fs) % Filter Bode Plot
set(subplot(2,1,1), 'XLim',[0 10]) % Zoom To Region-Of-Interest
set(subplot(2,1,2), 'XLim',[0 10])
L = 10*30;
t = linspace(0, Fs*L, Fs*L+1).'/Fs;
signal = sum(sin(2*pi*t.*[0.1:0.1:5 10:10:240]),2);
filtered_signal = filtfilt(sos,g,signal);
figure
plot(t, signal, 'DisplayName','Original')
hold on
plot(t, filtered_signal, 'DisplayName','Filtered')
hold off
xlabel('Time (s)')
ylabel('Amplitude')
legend('Location','best')
[FTsignal,Fv] = FFT1(signal, t);
[FTfiltered_signal,Fv] = FFT1(filtered_signal,t);
figure
plot(Fv, abs(FTsignal)*2, 'DisplayName','Original')
hold on
plot(Fv, abs(FTfiltered_signal)*2, 'DisplayName','Filtered')
hold off
grid
xlim([0 5])
xlabel('Frequency (Hz)')
ylabel('Magnitude')
legend('Location','best')
function [FTs1,Fv] = FFT1(s,t)
s = s(:);
t = t(:);
L = numel(t);
Fs = 1/mean(diff(t));
Fn = Fs/2;
NFFT = 2^nextpow2(L);
FTs = fft((s - mean(s)).*hann(L), NFFT)/sum(hann(L));
Fv = linspace(0, 1, NFFT/2+1)*Fn;
Iv = 1:numel(Fv);
FTs1 = FTs(Iv);
end
.

请先登录,再进行评论。

类别

Help CenterFile Exchange 中查找有关 Audio Processing Algorithm Design 的更多信息

Community Treasure Hunt

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

Start Hunting!

Translated by