51128

Фільтрація сигналів

Лабораторная работа

Информатика, кибернетика и программирование

Мета роботи: набути навичок проектування цифрових фільтрів, задавання специфікації фільтрів залежно від властивостей сигналів, які треба фільтрувати; набути навичок реалізації дискретної фільтрації сигналів у середовищі MatLAB.

Украинкский

2014-02-06

889.81 KB

10 чел.

Національний технічний університет України

«Київський політехнічний інститут»

Факультет електроніки

Лабораторна робота № 5

з дисципліни «Теорія сигналів»

«Фільтрація сигналів»

Виконав:  студент 3-го курсу

групи ДП-92

 Лонтковський С.А.

Київ – 2011

Мета роботи: набути навичок проектування цифрових фільтрів, задавання специфікації фільтрів залежно від властивостей сигналів, які треба фільтрувати; набути навичок реалізації дискретної фільтрації сигналів у середовищі MatLAB.

Порядок роботи:

1. Сформувати вектор відліків часу тривалістю 5 с для частоти дискретизації 128 Гц. Сформувати прямокутний імпульс в момент часу 3 с тривалості 0.1 с. Спроектувати ФНЧ Батерворта для позбавлення сигналу від шуму (функції buttord, butter, filter). Оцінити відношення сигнал/шум на вході та виході фільтра.

2. Сформувати вектор відліків часу тривалістю 1 с для частоти дискретизації 128 Гц. Сформувати сигнали ділянки синусоїди частотою 10 Гц амплітуди 1 В. Додати випадковий сигнал з нульовим середнім значенням амплітуди 2 В. Спроектувати ФНЧ Чебишова І-го роду, ФВЧ Чебишова ІІ-роду, СФ Кауера для позбавлення сигналу від шуму. Оцінити відношення сигнал/шум на вході та виході фільтра.

4. Для оцифрованих сигналів електрокардіограми, електроенцефалограми, прочитаної з файлу, а також ЕЕГ здорової та хворої людини, сигналів артеріального та внутрішньочерепного тиску та плетизмограми спроектувати фільтри для позбавлення від мережевої перешкоди 50 Гц.

5. Для звукових сигналів, які отримані з різною частотою дискретизації, виконати розділення на три спектральні діапазони: до 450 Гц; від 450 Гц до 1 кГц; від 1 кГц до 4 кГц з використанням фільтрів Кауера. Прослухати отримані сигнали, зробити висновки.

close all;

clear all;

clc;

 

%завдання 1------------------------------------------------

fs=128;

T=5;

fn=fs/2;

Wp=15/fn;Ws=20/fn;

Rp=3;

Rs=25;

f=[0:1/T:fs-1/T];

t=[0:1/fs:T-1/fs];

x=rectpuls(t-3,0.1);

ssx=abs(fft(x)/length(x));

figure()

subplot(2,2,1);plot(t,x);grid on;

title('Вхідний сигнал X');xlabel('Час, c');ylabel('Значення');

subplot(2,2,2);plot(f,ssx);grid on;

title('Амплітудний спектр сигнала x');

xlabel('Частота, Гц');ylabel('Значення');

N=rand(size(x))-0.5;

subplot(2,2,3);plot(t,N);grid on;title('Шум ');

xlabel('Час, c');ylabel('Значення');

ssN=abs(fft(N))/length(N);

subplot(2,2,4);plot(f,ssN);grid on;

title('Амплітудний спектр шума');

xlabel('Частота, Гц');ylabel('Значення');

Nx=x+N;

figure()

subplot(2,2,1);plot(t,Nx);grid on;title('Зашумленний вхідний сигнал');

xlabel('Час, с');ylabel('Значення');

ssNx=abs(fft(Nx)/length(Nx));

subplot(2,2,2);plot(f,ssNx);grid on;title('Амплітудний спектр Nx');

xlabel('Частота, Гц');ylabel('Значення');

[N,W]=buttord(Wp,Ws,Rp,Rs);

[b,a]=butter(N,W);

h=freqz(b,a,f,fs);

subplot(2,1,2);plot(f,abs(h));grid on;title('Амплітудний спектр Фільтра');

xlabel('Частота, Гц');ylabel('Значення');

y=filter(b,a,Nx);

figure()

subplot(2,1,1);plot(t,y);grid on;title('Відфільтрований сигнал');

xlabel('Час, с');ylabel('Значення');

ssy=abs(fft(y)/length(y));

subplot(2,1,2);plot(f,ssy);xlabel('Частота, Гц');ylabel('Значення');

title('Амплітудний спектр відфільтрованого сигналу');grid on;

 

%завдання 2------------------------------------------------

fs=128;

T=1;

fn=fs/2;

Wp=15/fn;

Ws=20/fn;

Rp=3;

Rs=70;

Rs3=150;

t=[0:1/fs:T-1/fs];

f=[0:1/T:fs-1/T];

x=sin(2*pi*10*t);

figure()

subplot(2,2,1);plot(t,x);grid on;title('Вхідний синусоїдальний сигнал');

xlabel('Час, с');ylabel('Значення');

ssx=abs(fft(x)/length(x));

subplot(2,2,2);plot(t,ssx);grid on;

title('Амплітудний спектр синусоїдального сигнала');

xlabel('Частота, Гц');ylabel('Значення');

r=2*rand(1,length(t));

subplot(2,2,3);plot(t,r);grid on;title('Випадковий сигнал');

xlabel('Час, с');ylabel('Значення');

ssr=abs(fft(r)/length(r));

subplot(2,2,4);plot(f,ssr);grid on;

title('Амплітудний спектр випадкового сигналу');

xlabel('Частота, Гц');ylabel('Значення');

rx=x+r;

figure()

subplot(2,2,1);plot(t,rx);grid on;

title('Вхідний і випадковий сигнал');

xlabel('Час, с');ylabel('Значення');

ssrx=abs(fft(rx)/length(rx));

subplot(2,2,2);plot(ssrx);title('Амплітудний спектр сигнала rs');

xlabel('Частота, Гц');ylabel('Значення');grid on;

[N1,Wn1]=cheb1ord(Wp,Ws,Rp,Rs);

[b,a]=cheby1(N1,Rp,Wn1,'low');

h=freqz(b,a,f,fs);

subplot(2,1,2);plot(f,abs(h));grid on;title('Чебешов 1-го порядку');

xlabel('Частота, Гц');ylabel('Значення');

y1=filter(b,a,rx);

figure()

subplot(2,2,1);plot(t,y1);grid on;

title('Перший відфільтрований сигнал');

xlabel('Час, с');ylabel('Значення');

ssy1=abs(fft(y1)/length(y1));

subplot(2,2,2);plot(f,ssy1);grid on;

title('Амплітудний спектр');xlabel('Частота, Гц');ylabel('Значення');

[N2,Wn2]=cheb2ord(Wp,Ws,Rp,Rs);

[b2,a2]=cheby2(N2,Rp,Wn2,'high');

y2=filter(b2,a2,rx);

subplot(2,2,3);plot(t,y2);grid on;

title('Другий відфільтрований сигнал');

xlabel('Час, с');ylabel('Значення');

ssy2=abs(fft(y2)/length(y2));

subplot(2,2,4);plot(f,ssy2);grid on;

title('Амплітудний спектр');xlabel('Частота, Гц');ylabel('Значення');

h2=freqz(b2,a2,f,fs);

figure()

subplot(2,2,1);plot(f,abs(h2));grid on;title('Чебешов 2-го порядка');

xlabel('Частота, Гц');ylabel('Значення');

[N3,Wn3]=ellipord([15/fn 17/fn],[13/fn 20/fn],Rp,Rs3);

[b3,a3]=ellip(N3,Rp,Rs,Wn3);

h3=freqz(b3,a3,f,fs);

subplot(2,2,2);plot(f,abs(h3));title('СФ Кауера');grid on;

xlabel('Частота, Гц');ylabel('Значення');

y3=filter(b3,a3,rx);

subplot(2,2,3);plot(t,y3);title('Третій відфільтрований сигнал');

xlabel('Час, с');ylabel('Значення');grid on;

ssy3=abs(fft(y3)/length(y3));

subplot(2,2,4);plot(ssy3);grid on;

title('Амплiтудний спектр');

xlabel('Частота, Гц');ylabel('Значення');

 

%завдання 4------------------------------------------------

load('ecg_10')

d=d/max(abs(d));

fs=400;

fn=fs/2;

t=[0:1/fs:(length(d)-1)/fs];

f=[0:1/t(length(d)):fn];

ssd=fft(d)/length(d);

assd=abs(ssd);

Ws=[49.9/fn 50.1/fn];Wp=[49/fn 51/fn];

[N,Wn]=ellipord(Wp,Ws,3,160);

[B,A]=ellip(N,0.5,30,Wn,'stop');

modh=abs(freqz(B,A,length(f),fs));

y=filter(B,A,d);

ssy=fft(y)/length(y);

ay=abs(ssy);

subplot(4,1,1);plot(f,modh(1:length(f)));grid on;

xlabel('Частота,Гц');ylabel('Значення, мВ'); title('ЗФ фільтр');

subplot(3,2,3);plot(t,d);xlabel('Час,сек');ylabel('Значення, мВ'); title('ЕКГ');grid on;

subplot(3,2,4);plot(f,assd(1:length(f)));xlabel('Частота,Гц');ylabel('Значення, мВ');

title('Амплітудний спектр ЕКГ');grid on;

subplot(3,2,5);plot(t,y),xlabel('Час,сек');ylabel('Значення, мВ'); grid on;

subplot(3,2,6);plot(f,ay(1:length(f)));xlabel('Час,Гц');

ylabel('Значення, мВ'); grid on;

 

load('eeg_healthy_10')

sig=sig/max(abs(sig));

fs=256;

fn=fs/2;

t=[0:1/fs:(length(sig)-1)/fs];

f=[0:1/t(length(sig)):fn];

ssd=fft(sig)/length(sig);

assd=abs(ssd);

Ws=[49.9/fn 50.1/fn];Wp=[49/fn 51/fn];

[N,Wn]=ellipord(Wp,Ws,3,160);

[B,A]=ellip(N,0.5,30,Wn,'stop');

modh=abs(freqz(B,A,length(f),fs));

y=filter(B,A,sig);

ssy=fft(y)/length(y);

ay=abs(ssy);

figure

subplot(4,1,1);plot(f,modh(1:length(f)));xlabel('Частота,Гц');

ylabel('Значення, мВ');title('ЗФ фільтр');grid on;

subplot(3,2,3);plot(t,sig);xlabel('Час,сек');ylabel('Значення, мВ');

title('ЕЕГ здорової');grid on;

subplot(3,2,4);plot(f,assd(1:length(f)));xlabel('Час,Гц');

ylabel('Значення, мВ');title('амплітудний спектр ЕЕГ');grid on;

subplot(3,2,5);plot(t,y);xlabel('Час,сек');ylabel('Значення, мВ');grid on;

subplot(3,2,6);plot(f,ay(1:length(f)));xlabel('Час,Гц');

ylabel('Значення, мВ');grid on;

 

load('eeg_sick_10')

sig=sig/max(abs(sig));

fs=256;

fn=fs/2;

t=[0:1/fs:(length(sig)-1)/fs];

f=[0:1/t(length(sig)):fn];

ssd=fft(sig)/length(sig);

assd=abs(ssd);

Ws=[49.9/fn 50.1/fn];Wp=[49/fn 51/fn];

[N,Wn]=ellipord(Wp,Ws,3,160);

[B,A]=ellip(N,0.5,30,Wn,'stop');

modh=abs(freqz(B,A,length(f),fs));

y=filter(B,A,sig);

ssy=fft(y)/length(y);

ay=abs(ssy);

figure

subplot(4,1,1);plot(f,modh(1:length(f)));xlabel('Частота,Гц');

ylabel('Значення, мВ');title('ЗФ фільтр');grid on;

subplot(3,2,3);plot(t,sig);xlabel('Час,сек');ylabel('Значення, мВ');

title('ЕЕГ Хворої');grid on;

subplot(3,2,4);plot(f,assd(1:length(f)));xlabel('Час,Гц');

ylabel('Значення, мВ');title('амплітудний спектр ЕЕГ');grid on;

subplot(3,2,5);plot(t,y);xlabel('Час,сек');ylabel('Значення, мВ');grid on;

subplot(3,2,6);plot(f,ay(1:length(f)));xlabel('Час,Гц');

ylabel('Значення, мВ');grid on;

 

 

fid=fopen('TBI_ICP.txt', 'r');

d=fscanf(fid,'%f');

fs=125;

fn=fs/2;

t=[0:1/fs:(length(d)-1)/fs];

f=[0:1/t(length(d)):fn];

ssd=fft(d)/length(d);

assd=abs(ssd);

Ws=[49.9/fn 50.1/fn];Wp=[49/fn 51/fn];

[N,Wn]=ellipord(Wp,Ws,3,160);

[B,A]=ellip(N,0.5,30,Wn,'stop');

modh=abs(freqz(B,A,length(f),fs));

y=filter(B,A,d);

ssy=fft(y)/length(y);

ay=abs(ssy);

subplot(4,1,1);plot(f,modh(1:length(f)));grid on;

xlabel('Частота,Гц');ylabel('Значення'); title('ЗФ фільтр');

subplot(3,2,3);plot(t,d);xlabel('Час,сек');ylabel('Значення');

title('Внутрішньочерепний тиск');grid on;

subplot(3,2,4);plot(f,assd(1:length(f)));xlabel('Частота,Гц');ylabel('Значення');

title('Амплітудний спектр');grid on;

subplot(3,2,5);plot(t,y),xlabel('Час,сек');ylabel('Значення'); grid on;

subplot(3,2,6);plot(f,ay(1:length(f)));xlabel('Час,Гц');

ylabel('Значення'); grid on;

 

 

fid=fopen('TBI_ABP.txt', 'r');

d=fscanf(fid,'%f');

fs=125;

fn=fs/2;

t=[0:1/fs:(length(d)-1)/fs];

f=[0:1/t(length(d)):fn];

ssd=fft(d)/length(d);

assd=abs(ssd);

Ws=[49.9/fn 50.1/fn];Wp=[49/fn 51/fn];

[N,Wn]=ellipord(Wp,Ws,3,160);

[B,A]=ellip(N,0.5,30,Wn,'stop');

modh=abs(freqz(B,A,length(f),fs));

y=filter(B,A,d);

ssy=fft(y)/length(y);

ay=abs(ssy);

subplot(4,1,1);plot(f,modh(1:length(f)));grid on;

xlabel('Частота,Гц');ylabel('Значення'); title('ЗФ фільтр');

subplot(3,2,3);plot(t,d);xlabel('Час,сек');ylabel('Значення');

title('Артеріальний тиск');grid on;

subplot(3,2,4);plot(f,assd(1:length(f)));xlabel('Частота,Гц');ylabel('Значення');

title('Амплітудний спектр');grid on;

subplot(3,2,5);plot(t,y),xlabel('Час,сек');ylabel('Значення'); grid on;

subplot(3,2,6);plot(f,ay(1:length(f)));xlabel('Час,Гц');

ylabel('Значення'); grid on;

 

 

fid=fopen('vavreschuk','r');

k=fread(fid,[9,inf],'int16');

d=k(8,:);

fs=155.1;

fn=fs/2;

t=[0:1/fs:(length(d)-1)/fs];

f=[0:1/t(length(d)):fn];

ssd=fft(d)/length(d);

assd=abs(ssd);

Ws=[49.9/fn 50.1/fn];Wp=[49/fn 51/fn];

[N,Wn]=ellipord(Wp,Ws,3,160);

[B,A]=ellip(N,0.5,30,Wn,'stop');

modh=abs(freqz(B,A,length(f),fs));

y=filter(B,A,d);

ssy=fft(y)/length(y);

ay=abs(ssy);

subplot(4,1,1);plot(f,modh(1:length(f)));grid on;

xlabel('Частота,Гц');ylabel('Значення, мВ'); title('ЗФ фільтр');

subplot(3,2,3);plot(t,d);xlabel('Час,сек');ylabel('Значення, мВ');

title('Плетизмограма');grid on;

subplot(3,2,4);plot(f,assd(1:length(f)));xlabel('Частота,Гц');ylabel('Значення, мВ');

title('Амплітудний спектр');grid on;

subplot(3,2,5);plot(t,y),xlabel('Час,сек');ylabel('Значення, мВ'); grid on;

subplot(3,2,6);plot(f,ay(1:length(f)));xlabel('Час,Гц');

ylabel('Значення, мВ'); grid on;

 

 

%завдання 5------------------------------------------------

[x,fs,bits]=wavread('8kHz.wav');

fn=fs/2;

t=[0:1/fs:(length(x)-1)/fs];

f=[0:1/t(length(x)):fn];

ssx=fft(x)/length(x);

assx=abs(ssx);

figure()

subplot(2,1,1);plot(t,x);xlabel('Час, сек');ylabel('Значення');

title('Вхідний сигнал');grid on;

subplot(2,1,2);plot(f,assx(1:length(f)));xlabel('Частота, Гц');

ylabel('Значення');title('Амплітудний спектр вхідного сигналу');grid on;

[N,W]=cheb2ord(450/fn,460/fn,5,10);

[B,A]=cheby2(N,30,W);

modh=abs(freqz(B,A,length(f),fs));

y=filter(B,A,x);

ssy=fft(y)/length(y);

assy=abs(ssy);

wavplay(y,fs)

figure()

subplot(3,1,1);plot(f,modh);xlabel('Частота, Гц');ylabel('Значення');

title('АЧХ ФНЧ(до 450Гц)');grid on;

subplot(3,1,2);plot(t,y);xlabel('Час, сек');ylabel('Значення');

title('Відфільтрований сигнал');grid on;

subplot(3,1,3);plot(f,assy(1:length(f)));xlabel('Частота, Гц');grid on;

ylabel('Значення');title('Амплітудний спектр відфільтрованого сигналу');

 

[N,W]=cheb2ord([450/fn 1000/fn],[445/fn 1005/fn],5,6);

[B,A]=cheby2(N,50,W);

modh=abs(freqz(B,A,length(f),fs));

y=filter(B,A,x);

ssy=fft(y)/length(y);

assy=abs(ssy);

wavplay(y,fs)

figure()

subplot(3,1,1);plot(f,modh);xlabel('Частота, Гц');ylabel('Значення');

title('АЧХ ПФ(від 450Гц до 1000Гц)');grid on;

subplot(3,1,2);plot(t,y);xlabel('Час, сек');ylabel('Значення');

title('Відфільтрований сигнал');grid on;

subplot(3,1,3);plot(f,assy(1:length(f)));xlabel('Частота, Гц');grid on;

ylabel('Значення');title('Амплітудний спектр відфільтрованого сигналу');

 

[N,W]=cheb2ord([1000/fn 3980/fn],[990/fn 3990/fn],5,6);

[B,A]=cheby2(N,40,W);

modh=abs(freqz(B,A,length(f),fs));

y=filter(B,A,x);

ssy=fft(y)/length(y);

assy=abs(ssy);

wavplay(y,fs)

figure()

subplot(3,1,1);plot(f,modh);xlabel('Частота, Гц');ylabel('Значення');

title('АЧХ ПФ(від 450Гц до 1000Гц)');grid on;

subplot(3,1,2);plot(t,y);xlabel('Час, сек');ylabel('Значення');

title('Відфільтрований сигнал');grid on;

subplot(3,1,3);plot(f,assy(1:length(f)));xlabel('Частота, Гц');grid on;

ylabel('Значення');title('Амплітудний спектр відфільтрованого сигналу');

[x,fs,bits]=wavread('44.1kHz.wav');

fn=fs/2;

t=[0:1/fs:(length(x)-1)/fs];

f=[0:1/t(length(x)):fn];

ssx=fft(x)/length(x);

assx=abs(ssx);

figure()

subplot(2,1,1);plot(t,x);xlabel('Час, сек');ylabel('Значення');

title('Вхідний сигнал');grid on;

subplot(2,1,2);plot(f,assx(1:length(f)));xlabel('Частота, Гц');

ylabel('Значення');title('Амплітудний спектр вхідного сигналу');grid on;

[N,W]=cheb2ord(450/fn,460/fn,5,10);

[B,A]=cheby2(N,30,W);

modh=abs(freqz(B,A,length(f),fs));

y=filter(B,A,x);

ssy=fft(y)/length(y);

assy=abs(ssy);

wavplay(y,fs)

figure()

subplot(3,1,1);plot(f,modh);xlabel('Частота, Гц');ylabel('Значення');

title('АЧХ ФНЧ(до 450Гц)');grid on;

subplot(3,1,2);plot(t,y);xlabel('Час, сек');ylabel('Значення');

title('Відфільтрований сигнал');grid on;

subplot(3,1,3);plot(f,assy(1:length(f)));xlabel('Частота, Гц');grid on;

ylabel('Значення');title('Амплітудний спектр відфільтрованого сигналу');

 

[N,W]=cheb2ord([450/fn 1000/fn],[445/fn 1005/fn],5,10);

[B,A]=cheby2(N,50,W);

modh=abs(freqz(B,A,length(f),fs));

y=filter(B,A,x);

ssy=fft(y)/length(y);

assy=abs(ssy);

wavplay(y,fs)

figure()

subplot(3,1,1);plot(f,modh);xlabel('Частота, Гц');ylabel('Значення');

title('АЧХ ПФ(від 450Гц до 1000Гц)');grid on;

subplot(3,1,2);plot(t,y);xlabel('Час, сек');ylabel('Значення');

title('Відфільтрований сигнал');grid on;

subplot(3,1,3);plot(f,assy(1:length(f)));xlabel('Частота, Гц');grid on;

ylabel('Значення');title('Амплітудний спектр відфільтрованого сигналу');

 

[N,W]=cheb2ord([1000/fn 3980/fn],[990/fn 3990/fn],5,8);

[B,A]=cheby2(N,40,W);

modh=abs(freqz(B,A,length(f),fs));

y=filter(B,A,x);

ssy=fft(y)/length(y);

assy=abs(ssy);

wavplay(y,fs)

figure()

subplot(3,1,1);plot(f,modh);xlabel('Частота, Гц');ylabel('Значення');

title('АЧХ ПФ(від 450Гц до 1000Гц)');grid on;

subplot(3,1,2);plot(t,y);xlabel('Час, сек');ylabel('Значення');

title('Відфільтрований сигнал');grid on;

subplot(3,1,3);plot(f,assy(1:length(f)));xlabel('Частота, Гц');grid on;

ylabel('Значення');title('Амплітудний спектр відфільтрованого сигналу');

Графіки.

1)

2)

4)

5)

Ч    Частота 8кГц.

Ч  Частота 44.1кГц.


 

А также другие работы, которые могут Вас заинтересовать

60807. Урок по моделированию многоэтажного здания в 3ds Max 6.4 MB
  В данном уроке мы рассмотрим способ моделирования высокополигонального современного многоэтажного здания в 3ds Max. Чтобы ясно себе представлять будущую модель, следует пользоваться...
60809. Модификатор EDIT MESH (редактирование сетки) 173.5 KB
  Создание яблока Рассмотрим пример создания модели яблока при помощи модификатора Edit Mesh Редактирование поверхности. 2 сформировав углубление в месте крепления корешка яблока. Не меняя настройки плавного выделения переместите выделенные...
60810. Назначение и настройка модификаторов в 3ds max 733.5 KB
  Кнопка Закрепить стек позволяет зафиксировать меню стека на экране таким образом что оно не исчезнет если снять выделение с объекта или даже выделить другой объект. Кнопка Показывать конечный результат показывает конечный...
60814. Зимние праздники 45.5 KB
  Цель урока: Формирование у учащихся более полного представления о народной праздничной культуре русского народа. Задачи: 1. Обобщить полученные знания по русским зимним праздникам.