2016-04-17 83 views
3

过滤不需要的频率我正在做一个带通滤波器,我已经创建了基于正弦波一些不需要的频率的信号:带通滤波器未能在MATLAB

Fs = 8e3; % Sampling Rate 
fn = Fs/2; % Nyquist frequency 
L = 1e3; % Length of signal 
T = 1/Fs; % Sampling period 
t = T*linspace(0,L,Fs); % Time domain 

% Frequencies in Hz 
f1 = 1500; 
f2 = 700; 
f3 = 2500; 
f4 = 3500; 

% Signal 
x = 6*sin(2*pi*f1*t); 

% Noise 
noise = 3*sin(2*pi*f2*t)... 
     + 2*sin(2*pi*f3*t)... 
     + 1*sin(2*pi*f4*t); 

x_noise = x + noise; 

然后创建一个巴特沃思带通滤波器:

[b,a] = butter(10,[1000 2000]/fn,'bandpass'); 

在时间和频率空间中的信号,与所述带通响应(与freqz)看起来像这样:

Fig. 1 Signal with corruption | Fig. 2 Frequency with bandpass response

我会一直在这里上算了一下,简单地做

xf = filter(b,a,x_noise); 

就已经产生了非常类似于原始信号的东西,但很可惜,我得到的是really far from the filtered signal with a high response far from the bandpass

我在这里做错了什么?

下面是完整的代码:

clear all 
Fs = 8e3; % Sampling Rate 
fn = Fs/2; % Nyquist frequency 
L = 1e3; % Length of signal 
T = 1/Fs; % Sampling period 
t = T*linspace(0,L,Fs); % Time domain 

% Frequencies in Hz 
f1 = 1500; 
f2 = 700; 
f3 = 2500; 
f4 = 3500; 

% Signal 
x = 6*sin(2*pi*f1*t); 

% Noise 
noise = 3*sin(2*pi*f2*t)... 
     + 2*sin(2*pi*f3*t)... 
     + 1*sin(2*pi*f4*t); 

x_noise = x + noise; 

subplot(221); 
idx = 1:round(length(t)/30); 
plot(t(idx),x(idx),t(idx),x_noise(idx)); 
xlabel('Time (s)'); ylabel('Signal Amplitudes'); 
legend('Original signal','Noisy signal'); 

% Frequency space 
f = fn*linspace(0,1,L/2+1); 
X = fft(x_noise)/L; 

[b,a] = butter(10,[1000 2000]/fn,'bandpass'); 
h = abs(freqz(b,a,floor(L/2+1))); 

subplot(222); 
plot(f,abs(X(1:L/2+1)),f,h*max(abs(X))); 
xlabel('Freq (Hz)'); ylabel('Frequency amplitudes'); 
legend('Fourier Transform of signal','Filter amplitude response'); 

% Filtered signal 
xf = filter(b,a,x_noise); 
subplot(223) 
plot(t(idx),xf(idx)); 
xlabel('Time (s)'); ylabel('Signal Amplitudes'); 
legend('Filtered signal'); 

% Filtered in frequency space 
Xf = abs(fft(xf)/L); 
subplot(224); 
plot(f,Xf(1:L/2+1),f,h*5e-6); 
xlabel('Freq (Hz)'); ylabel('Frequency amplitudes'); 
legend('Fourier Transform of filtered signal','Bandpass'); 

回答

1

你的时间变量t是错误的,e.g 1/(t(2)-t(1))应该给Fs但事实并非如此。

尝试,而不是:

t = T*(0:L-1); 
+0

现在只是传递的整个信号[原创与“过滤”() –

+0

那不是我所看到的。看看频率域(你自己的情节) –

+0

[this](http://i.imgur.com/Fo21pxS.png)是我得到 –