96958

Построение и исследование характеристик цифрового эквивалента аналогового фильтра

Курсовая

Коммуникация, связь, радиоэлектроника и цифровые приборы

Исследование характеристик аналоговой цепи. Расчёт цифровой цепи методом Эйлера. Расчёт цифровой цепи методом билинейного преобразования. Расчёт цифровой цепи методом инвариантной импульсной характеристики. Разработка и тестирование алгоритма цифровой фильтрации.

Русский

2015-10-12

1.98 MB

3 чел.

МИНОБРНАУКИ РОССИИ

Федеральное государственное бюджетное образовательное учреждение

высшего профессионального образования

"Московский государственный технический университет радиотехники,

электроники и автоматики"

МГТУМИРЭА

факультет кибернетики

Кафедра автоматических систем

КУРСОВАЯ РАБОТА

по дисциплине

«Цифровая обработка сигналов»

Тема курсовой работы

«Построение и исследование характеристик цифрового эквивалента аналогового фильтра»

Студент группы  КУБ-1-12  

Барашков А.А.   (шифр 121290)

Руководитель курсовой работы

профессор, д.т.н., профессор

Асанов А.З.

Работа представлена к защите

«__»_______2015 г.

(подпись студента)

«Допущен к защите»

«__»_______2015 г.

(подпись руководителя)

Москва 2015

МИНОБРНАУКИ РОССИИ

Федеральное государственное бюджетное образовательное учреждение

высшего профессионального образования

"Московский государственный технический университет радиотехники,

электроники и автоматики"

МГТУМИРЭА

факультет Кибернетика

Кафедра Автоматических систем

Утверждаю

Заведующий

кафедрой_____________

«____» __________201___ г.

ЗАДАНИЕ

на выполнение  курсовой работы

по дисциплине «Цифровая обработка сигналов                                      »

Студент          Барашков Алексей Андреевич                                      Группа   КУБ-1-12       

  1.  Тема  « Построение и исследование характеристик цифрового эквивалента аналогового фильтра».
  2.  Исходные данные:

R1=1000 Ом, R2=4000 Ом, С1=1,5 мкФ, C2=0,5 мкФ.

  1.  Перечень вопросов, подлежащих разработке, и обязательного графического материала:

  1.  Срок представления к защите курсовой работы:до «___» _______2015 г.

Задание на курсовую

работу выдал

«___»______201__г.

Подпись руководителя

Ф.И.О. руководителя

работы

Задание на курсовую

работу получил

«___»______201__г.

Подпись студента –

исполнителя

Ф.И.О. студента -

исполнителя

работы

5. Мониторинг процесса выполнения курсовой работы

№ этапа

Этап курсовой работы выполнил и представил результаты руководителю проекта (работы),

дата и подпись исполнителя

Работу по этапу курсовой  работы принял на рассмотрение,

дата и подпись руководителя

Рекомендации и замечания по
этапу курсовой

работы
выдал
исполнителю,

дата и подпись руководителя

Комментарии
руководителя курсовой проекта работы

1

2

3

4

5

МИНОБРНАУКИ РОССИИ

Федеральное государственное бюджетное образовательное учреждение

высшего профессионального образования

"Московский государственный технический университет радиотехники,

электроники и автоматики"

МГТУМИРЭА

факультет Кибернетика

Кафедра Автоматических систем

Протокол заседания комиссии по защите курсовой работы

от ________________201__г.  №__________

Состав комиссии:______________________________________________________________________

______________________________________________________________________

(должность, ученая степень, ученое звание)

Утверждена распоряжением заведующего кафедрой ________________________________________

(наименование кафедры)

от «___» ________ 201_ г. №________.

Слушали защиту курсовой работы ________________________________________

____________________________________________________________________________________________

по дисциплине _______________________________________________________________________________

студента группы______________________________________________________________________________

    (группа)         (Ф.И.О.)

Во время защиты курсовой работы  были заданы следующие вопросы:

1 _____________________________________________________________________________

2. _____________________________________________________________________________

3. _____________________________________________________________________________

Итоговая (комплексная) оценка выполнения и защиты курсовой работы_____________________

Члены комиссии    ____________________________________________

____________________________________________

(подпись)                                       (Ф.И.О.)


Оглавление

[1] Федеральное государственное бюджетное образовательное учреждение

[2] высшего профессионального образования

[3] "Московский государственный технический университет радиотехники,

[4] электроники и автоматики"

[5] МГТУМИРЭА

[6] Федеральное государственное бюджетное образовательное учреждение

[7] высшего профессионального образования

[8] "Московский государственный технический университет радиотехники,

[9] электроники и автоматики"

[10] МГТУМИРЭА

[11] Федеральное государственное бюджетное образовательное учреждение

[12] высшего профессионального образования

[13] "Московский государственный технический университет радиотехники,

[14] электроники и автоматики"

[15] МГТУМИРЭА

[16] Оглавление

[17] Исследование характеристик аналоговой цепи

[18] Расчёт цифровой цепи методом Эйлера.

[19] Расчёт цифровой цепи методом билинейного преобразования.

[20] Расчёт цифровой цепи методом инвариантной импульсной характеристики.

[21] Разработка и тестирование алгоритма цифровой фильтрации.


Исследование характеристик аналоговой цепи

Дана аналоговая пассивная цепь:

Рис. 1. Схема аналоговой цепи

Для цепи заданы параметры:

R1=1000 Ом

R2=4000 Ом

С1=1,5 мкФ

C2=0,5 мкФ

Воспользуемся символическим методом расчёта электрических цепей.

Рис. 2. Схема для расчёта аналоговой цепи

Применим метод «двух узлов». Узел a – заземлим, его потенциал будет равен нулю. Найдём потенциал узла c.

Ток, протекающий через резистор R2:

Далее найдём потенциал в точке b. Так как , а ,

Получим частотную передаточную функцию цепи:

Код в Matlab:

clc;clear;

 

R1=1000;

R2=4000;

C1=1.5*10^(-6);

C2=0.5*10^(-6);

 

Wa0=(R2*C2)^2;

Wa1=R2*C2;

Wb0=C1*R1*(R2*C2)^2;

Wb1=C2*R2*(C2*R2+2*C1*R1+C2*R1);

Wb2=2*C2*R2+C1*R1+C2*R1;

 

Wanalog=tf([Wa0 Wa1 0], [Wb0 Wb1 Wb2 1])

Подставив значения параметров R и C, получим:

Так как преобразования Фурье и Лапласа имеют схожий вид, можно выполнить замену p=jω. Тогда передаточная функция аналоговой цепи:

Найдём для этой передаточной функции нули и полюса.

Код в Matlab:

disp('Полюса:');

disp(pole(Wanalog));

disp('Нули:');

disp(zero(Wanalog));

Нули:

p=0; p=-500

Полюса:

p=-1000; p=-500; p=-1000/3

Исходя из найденных значений, можем записать передаточную функцию в виде:

Этой передаточной функции соответствует дифференциальное уравнение:

Получим логарифмические амплитудно-частотные и фазо-частотные характеристики для этого фильтра:

Амплитудно-частотная характеристика:

Фазо-частотная характеристика:

Построим графики этих характеристик:

Код в Matlab:

figure(1);

bode(Wanalog,'g');

grid on;

Рис. 3. Логарифмические частотные характеристики аналогового фильтра

Для проверки предыдущих вычислений можно сравнить полученную выше ЛАЧХ с ЛАЧХ, полученной при непосредственном моделировании электрической цепи в пакете Electronics Workbench:

Рис. 4. Моделирование электрической цепи в пакете Electronics Workbench

Получим импульсную переходную функцию этого фильтра.

Для этого разложим передаточную функцию на простые дроби:

Найдём коэффициенты:

C0=1,5

C1=-0,5

С помощью обратного преобразования Лапласа, найдём ИПФ:

Построим в системе Matlab график ИПФ:

Код в Matlab:

figure(2);

impulse(Wanalog,'g');

grid on;

Рис. 5. График импульсной переходной функции  аналогового фильтра.

Расчёт цифровой цепи методом Эйлера.

 Для получения передаточной функции цифрового фильтра необходимо произвести замену , где Td -  интервал дискретизации.

По ЛАЧХ аналогового фильтра найдём частоту, при которой происходит ослабление амплитуды в 10 раз (-20дБ). Примем её за ширину спектра.

ωmax=6580,668 рад/с

fmax= =1047,346 Гц

Возьмём частоту дискретизации в 5 раз больше ширины спектра:

fd=5fmax=5236,730 Гц

Тогда период дискретизации:

Td=  =0,00019096 с

Выбрав период дискретизации, можно найти коэффициенты передаточной функции и разностного уравнения:

Код Matlab:

clc;clear;

 

Td=0.00019096; %Период дискретизации

 

A=2/3*10^3;

%Полюса аналоговой передаточной функции

p0=1000;

p1=1000/3;

 

u=(1+Td*(p0+p1+Td*p0*p1));

%Коэффициенты передаточной функции цифровой цепи

b0=A*Td/u

b1=-A*Td/u

a1=(-2-Td*(p0+p1))/u

a2=1/u

 

b=[b0 b1];

a=[1 a1 a2];

Полученные коэффициенты:

b0 =  0,1005

b1 = -0,1005

a1 =  -1,7798

a2 = 0,7894

Передаточная функция:

Разностное уравнение:

Найдём нули и полюса этой передаточной функции:

Код Matlab:

[q,p]=tf2zpk(b,a);

disp('Нули');

disp(q);

disp('Полюса');

disp(p);

figure(4);

zplane(b,a);

Нули:

z=0; z=1

Полюса:

z=0,9402; z= 0,8397

Рис. 6. Карта нулей и полюсов цифрового фильтра, полученного методом Эйлера при Td=0,00019096 с

Получим амплитудно-частотные и фазо-частотные характеристики этого фильтра. Для этого в передаточной функции произведём замену .

Амплитудно-частотная характеристика:

Фазо-частотная характеристика:

Код Matlab:

w=logspace(1,5,10000);

figure(1);

Wd=(b0+b1*exp(-1i*w*Td))./(1+a1*exp(-1i*w*Td)+a2*exp(-2*1i*w*Td));

% Построение графиков ЛАЧХ и ЛФЧХ цифровой цепи

subplot(2,1,1), loglog(w,abs(Wd),'b'), grid on, xlabel('w (Rad/s)'), title('MAGNITUDE - |H(w)|');

hold on;

subplot(2,1,2), semilogx(w,180/pi*angle(Wd),'b'), grid on, xlabel('w (Rad/s)'), title('PHASE - arg [H(w)] (deg)');

hold on;

figure(2);

% Построение графиков ЛАЧХ аналоговой и цифровой цепи на одном полотне Wa=A*1i*w./(1i*w+p0)./(1i*w+p1);

loglog(w,abs(Wd),'b');

hold on;

loglog(w,abs(Wa),'g');

hold on;

grid on; xlabel('w (Rad/s)'); title('MAGNITUDE - |H(w)|'); axis([0 10^5 10^(-2) 1]);

Рис. 7.  ЛАЧХ аналогового фильтра (обозначена зелёным) и соответствующего ему цифрового фильтра, полученного методом Эйлера при Td=0,00019096 с (обозначена синим)

Рис. 8.  ЛАЧХ и ЛФЧХ цифрового фильтра, полученного методом Эйлера при Td=0,00019096 с

Рис. 9.  Амплитудно-частотная характеристика (не в логарифмическом масштабе) аналогового фильтра (обозначена зелёным) и соответствующего ему цифрового фильтра, полученного методом Эйлера при Td=0,00019096 с (обозначена синим)

Как видно из полученных графиков (рисунки 7 и 9), частотные характеристики цифрового фильтра по форме почти точно повторяют характеристики аналогового прототипа на участке [0; ], то есть [0; 1,645·104] рад/c. Далее идёт периодическое повторение характеристики цифрового фильтра.

Импульсная переходная функция цифрового фильтра может быть получена за счёт обратного Z-преобразования его передаточной функции. Мы получим её с помощью пакета Matlab.

Код Matlab:

N=0.012/Td;

n=0:(N-1);

h=impz(b,a,N);

figure(3);

title('Impulse Response h(n*Td) - impz');

hold on;

xlabel('n');

ylabel('h(n*Td)');

plot(n,h,'b');

grid on;

Рис. 10.  Импульсная переходная функция цифрового фильтра, полученного методом Эйлера при Td=0,00019096 с

По рисунку 10 видно, что ИПФ цифрового фильтра очень близка по форме к ИПФ аналогового и отличается, в первую очередь, на величину вещественного коэффициента.

Теперь возьмём частоту дискретизации в 10 раз больше ширины спектра:

fd=10fmax=10473,460 Гц

Тогда период дискретизации:

Td=  =0,00009548с

Найдём коэффициенты передаточной функции и разностного уравнения:

b0 =  0,0563

b1 = -0,0563

a1 =  -1,8820

a2 = 0,8847

Передаточная функция:

Разностное уравнение:

Нули:

z=0; z=1

Полюса:

z=0,9692; z= 0,9128

Рис. 11. Карта нулей и полюсов цифрового фильтра, полученного методом Эйлера при Td=0,00009548 с

Получим амплитудно-частотные и фазо-частотные характеристики этого фильтра. Для этого в передаточной функции произведём замену .

Амплитудно-частотная характеристика:

Фазо-частотная характеристика:

Рис. 12.  ЛАЧХ аналогового фильтра (обозначена зелёным) и соответствующего ему цифрового фильтра, полученного методом Эйлера при Td=0,00009548 с  (обозначена синим)

Рис. 13.  ЛАЧХ и ЛФЧХ цифрового фильтра, полученного методом Эйлера при Td=0,00009548 с  

В результате сравнения графиков на рисунках 7 и 12, 8 и 13 можно сделать вывод, что при уменьшении периода дискретизации различия между частотными характеристиками цифрового фильтра и его аналогового прообраза становятся несколько меньше. Также у частотных характеристик цифрового фильтра увеличивается период.

Импульсная переходная функция:

Рис. 14.  Импульсная переходная функция цифрового фильтра, полученного методом Эйлера при Td=0,00009548 с  

Из сравнения ИПФ при разных периодах дискретизации (рис. 10 и рис. 14) видно, что при уменьшении Td в несколько раз,  во столько же раз уменьшается амплитуда ИПФ.

Расчёт цифровой цепи методом билинейного преобразования.

 Для получения передаточной функции цифрового фильтра необходимо произвести замену , где Td -  интервал дискретизации.

Примем период дискретизации:

Td=0,00019096 с

Выбрав период дискретизации, можно найти коэффициенты передаточной функции и разностного уравнения:

Код Matlab:

clc;clear;

 

Td=0.00019096; %Период дискретизации

 

A=2/3*10^3;

%Полюса аналоговой передаточной функции

p0=1000;

p1=1000/3;

 

u=(2/Td+p0)*(2/Td+p1);

%Коэффициенты передаточной функции цифровой цепи

b0=2*A/Td/u

b1=0

b2=-2*A/Td/u

a1=(-8/Td^2+2*p0*p1)/u

a2=(-2/Td+p0)*(-2/Td+p1)/u

 

b=[b0 b1 b2];

a=[1 a1 a2];

Полученные коэффициенты:

b0 =  0,0563

b1 = 0

b2 = -0,0563

a1 =  -1,7640

a2 = 0,7747

Передаточная функция:

Разностное уравнение:

Найдём нули и полюса этой передаточной функции:

Код Matlab:

[q,p]=tf2zpk(b,a);

disp('Нули');

disp(q);

disp('Полюса');

disp(p);

figure(4);

zplane(b,a);

Нули:

z=-1; z=1

Полюса:

z=0,9383; z= 0,8257

Рис. 15. Карта нулей и полюсов цифрового фильтра, полученного методом билинейного преобразования при Td=0,00019096 с

Получим амплитудно-частотные и фазо-частотные характеристики этого фильтра. Для этого в передаточной функции произведём замену .

Амплитудно-частотная характеристика:

Фазо-частотная характеристика:

Код Matlab:

w=logspace(1,5,10000);

Wd=(b0+b2*exp(-2*1i*w*Td))./(1+a1*exp(-1i*w*Td)+a2*exp(-2*1i*w*Td));

% Построение графиков ЛАЧХ и ЛФЧХ цифровой цепи

figure(1);

subplot(2,1,1), loglog(w,abs(Wd),'b'), grid on, xlabel('w (Rad/s)'), title('MAGNITUDE - |H(w)|');

hold on;

subplot(2,1,2), semilogx(w,180/pi*angle(Wd),'b'), grid on, xlabel('w (Rad/s)'), title('PHASE - arg [H(w)] (deg)');

hold on;

% Построение графиков ЛАЧХ аналоговой и цифровой цепи на одном полотне

figure(2);

Wa=A*1i*w./(1i*w+p0)./(1i*w+p1);

loglog(w,abs(Wd),'b');

hold on;

loglog(w,abs(Wa),'g');

hold on;

grid on; xlabel('w (Rad/s)'); title('MAGNITUDE - |H(w)|'); axis([0 10^5 10^(-2) 1]);

Рис. 16.  ЛАЧХ аналогового фильтра (обозначена зелёным) и соответствующего ему цифрового фильтра, полученного методом билинейного преобразования при Td=0,00019096 с (обозначена синим)

Рис. 17.  ЛАЧХ и ЛФЧХ цифрового фильтра, полученного методом билинейного преобразования при Td=0,00019096 с

Полученная АЧХ цифрового фильтра, синтезированного методом билинейного преобразования, больше соответствует АЧХ аналогового прототипа в области частот от 10 до 10^4 рад/c, чем методом Эйлера при том же периоде дискретизации.

Импульсная переходная функция цифрового фильтра может быть получена за счёт обратного Z-преобразования его передаточной функции. Мы получим её с помощью пакета Matlab.

Код Matlab:

N=0.012/Td;

n=0:(N-1);

h=impz(b,a,N);

figure(3);

title('Impulse Response h(n*Td) - impz');

hold on;

xlabel('n');

ylabel('h(n*Td)');

plot(n,h,'b');

grid on;

Рис. 18.  Импульсная переходная функция цифрового фильтра, полученного методом билинейного преобразования при Td=0,00019096 с

ИПФ фильтра, полученного методом БП, несколько отличается по форме от ИПФ прототипа для начальных отсчётов.

Теперь возьмём другой период дискретизации

Td=0,00009548с

Найдём коэффициенты передаточной функции и разностного уравнения:

b0 =  0,0299

b1 = 0

b2 = -0,0299

a1 =  -1,8775

a2 = 0,8804

Передаточная функция:

Разностное уравнение:

Нули:

z=-1; z=1

Полюса:

z=0,9687; z= 0,9089

Рис. 19. Карта нулей и полюсов цифрового фильтра, полученного методом билинейного преобразования при Td=0,00009548 с

Получим амплитудно-частотные и фазо-частотные характеристики этого фильтра. Для этого в передаточной функции произведём замену .

Амплитудно-частотная характеристика:

Фазо-частотная характеристика:

Рис. 20.  ЛАЧХ аналогового фильтра (обозначена зелёным) и соответствующего ему цифрового фильтра, полученного методом билинейного преобразования при Td=0,00009548 с  (обозначена синим)

Рис. 21.  ЛАЧХ и ЛФЧХ цифрового фильтра, полученного методом билинейного преобразования при Td=0,00009548 с  

Импульсная переходная функция:

Рис. 22.  Импульсная переходная функция цифрового фильтра, полученного методом билинейного преобразования при Td=0,00009548 с  

При уменьшении периода дискретизации у фильтра, полученного методом билинейного преобразования, характеристики изменились так же, как и у фильтра, полученного методом Эйлера.

Расчёт цифровой цепи методом инвариантной импульсной характеристики.

Ранее была получена импульсная переходная функция аналогового фильтра:

В этом методе синтеза ИПФ ЦФ определяется через дискретные значения ИПФ аналогового прототипа:

 Восстановим передаточную функцию по ИПФ:

Примем период дискретизации:

Td=0,00019096 с

Выбрав период дискретизации, можно найти коэффициенты передаточной функции и разностного уравнения:

Код Matlab:

clc;clear;

 

Td=0.00019096; %Период дискретизации

 

A=2/3*10^3;

C0=1.5;

C1=-0.5;

%Полюса аналоговой передаточной функции

p0=1000;

p1=1000/3;

 

%Коэффициенты передаточной функции цифровой цепи

b0=A*Td*(C0+C1)

b1=A*Td*(-C0*exp(-p1*Td)-C1*exp(-p0*Td))

a1=-exp(-p1*Td)-exp(-p0*Td)

a2=exp(-(p0+p1)*Td)

 

b=[b0 b1];

a=[1 a1 a2];

Полученные коэффициенты:

b0 =  0,1273

b1 = -0,1266

a1 =  -1,7645

a2 = 0,7752

Передаточная функция:

Разностное уравнение:

Найдём нули и полюса этой передаточной функции:

Код Matlab:

[q,p]=tf2zpk(b,a);

disp('Нули');

disp(q);

disp('Полюса');

disp(p);

figure(4);

zplane(b,a);

Нули:

z=0; z=0,9944

Полюса:

z=0,9383; z= 0,8262

Рис. 23. Карта нулей и полюсов цифрового фильтра, полученного методом инвариантной импульсной характеристики при Td=0,00019096 с

Получим амплитудно-частотные и фазо-частотные характеристики этого фильтра. Для этого в передаточной функции произведём замену .

Амплитудно-частотная характеристика:

Фазо-частотная характеристика:

Код Matlab:

w=logspace(1,5,10000);;

Wd=(b0+b1*exp(-1i*w*Td))./(1+a1*exp(-1i*w*Td)+a2*exp(-2*1i*w*Td));

% Построение графиков ЛАЧХ и ЛФЧХ цифровой цепи figure(1)

subplot(2,1,1), loglog(w,abs(Wd),'b'), grid on, xlabel('w (Rad/s)'), title('MAGNITUDE - |H(w)|');

hold on;

subplot(2,1,2), semilogx(w,180/pi*angle(Wd),'b'), grid on, xlabel('w (Rad/s)'), title('PHASE - arg [H(w)] (deg)');

hold on;

% Построение графиков ЛАЧХ аналоговой и цифровой цепи на одном полотне

figure(2);

Wa=A*1i*w./(1i*w+p0)./(1i*w+p1);

loglog(w,abs(Wd),'b');

hold on;

loglog(w,abs(Wa),'g');

hold on;

grid on; xlabel('w (Rad/s)'); title('MAGNITUDE - |H(w)|'); axis([0 10^5 10^(-2) 1]);

Рис. 24.  ЛАЧХ аналогового фильтра (обозначена зелёным) и соответствующего ему цифрового фильтра, полученного методом инвариантной импульсной характеристики при Td=0,00019096 с (обозначена синим)

Рис. 25.  ЛАЧХ и ЛФЧХ цифрового фильтра, полученного методом инвариантной импульсной характеристики при Td=0,00019096 с

АЧХ цифрового фильтра, синтезированного методом инвариантной импульсной характеристики, имеет существенно большие различия с АЧХ аналогового прототипа, чем у фильтров, полученных методами Эйлера и билинейного преобразования.

Импульсная переходная функция цифрового фильтра может быть получена за счёт обратного Z-преобразования его передаточной функции. Мы получим её с помощью пакета Matlab.

Код Matlab:

N=0.012/Td;

n=0:(N-1);

h=impz(b,a,N);

figure(3);

title('Impulse Response h(n*Td) - impz');

hold on;

xlabel('n');

ylabel('h(n*Td)');

plot(n,h,'b');

grid on;

Рис. 26.  Импульсная переходная функция цифрового фильтра, полученного методом инвариантной импульсной характеристики при Td=0,00019096 с

ИПФ цифрового фильтра,  полученного методом ИИХ отличается от ИПФ аналогового прототипа на вещественный множитель, равный периоду дискретизации.

Теперь возьмём другой период дискретизации

Td=0,00009548с

Найдём коэффициенты передаточной функции и разностного уравнения:

b0 =  0,0637

b1 = -0,0636

a1 =  -1,8776

a2 = 0,8805

Передаточная функция:

Разностное уравнение:

Нули:

z=0; z=0,9985

Полюса:

z=0,9687; z= 0,9089

Рис. 27. Карта нулей и полюсов цифрового фильтра, полученного методом инвариантной импульсной характеристики при Td=0,00009548 с

Получим амплитудно-частотные и фазо-частотные характеристики этого фильтра. Для этого в передаточной функции произведём замену .

Амплитудно-частотная характеристика:

Фазо-частотная характеристика:

Рис. 28.  ЛАЧХ аналогового фильтра (обозначена зелёным) и соответствующего ему цифрового фильтра, полученного методом инвариантной импульсной характеристики при Td=0,00009548 с  (обозначена синим)

Рис. 29.  ЛАЧХ и ЛФЧХ цифрового фильтра, полученного методом инвариантной импульсной характеристики при Td=0,00009548 с  

Импульсная переходная функция:

Рис. 30.  Импульсная переходная функция цифрового фильтра, полученного методом инвариантной импульсной характеристики при Td=0,00009548 с  

При уменьшении периода дискретизации у фильтра, полученного методом билинейного преобразования, характеристики изменились так же, как и у фильтра, полученного методом Эйлера.

Разработка и тестирование алгоритма цифровой фильтрации.

Создадим алгоритм цифровой фильтрации. Для этого воспользуемся данными, полученным с помощью метода Эйлера для периода дискретизации  Td=0,00019096 с.

b0 =  0,1005

b1 = -0,1005

a1 =  -1,7798

a2 = 0,7894

Передаточная функция:

Разностное уравнение:

 С помощью пакета Matlab Simulink составим структурную схему:

Рис. 31.  Структурная схема алгоритма фильтрации

Этой схеме соответствует следующая функция:

Код Matlab:

function V=EuFilter(U)

b0=0.1005;

b1=-0.1005;

a1=-1.7798;

a2=0.7894;

length=size(U,2);

V=zeros(1,length);

V(1)=b0*U(1);

V(2)=b0*U(2)+b1*U(1)-a1*V(1);

for n=3:1:length

   V(n)=b0*U(n)+b1*U(n-1)-a1*V(n-1)-a2*V(n-2);

end

end

В качестве тестовых сигнала используем синусоиды с разными частотами. Пропустим их через фильтр, построим графики входных и выходных сигналов, их спектров.

Чтобы легче было проанализировать полученные результаты, приведём АЧХ и ФЧХ цифровой цепи:

Код Matlab:

figure(5);

w=50:1:6580;

Wd=(b0+b1*exp(-1i*w*Td))./(1+a1*exp(-1i*w*Td)+a2*exp(-2*1i*w*Td));

subplot(2,1,1), plot(w,abs(Wd),'b');

hold on;

subplot(2,1,2), plot(w,180/pi*angle(Wd),'b');

hold on;

subplot(2,1,1), grid on, xlabel('w (Rad/s)'), title('MAGNITUDE - |H(w)|');

subplot(2,1,2), grid on, xlabel('w (Rad/s)'), title('PHASE - arg [H(w)] (deg)');

Рис. 32.  АЧХ и ФЧХ фильтра

Минимальное подавление амплитуды наблюдается при частоте ω=544,5 рад/с. В этом случае выходная амплитуда будет составлять 0,478 от входной. ФЧХ пересекает линию 0º при частоте ω=577,64 рад/с. При этой частоте выходной сигнал не будет сдвинут относительно входного.


 Тестирование фильтра:

  1.  Примем частоту ω=40 рад/с

Код Matlab:

clc;clear;

Td=0.00019096; % Период дискретизации

N=2^13; % Количество отсчётов

t=0:Td*(N-1)*Td;

w=40;

U=sin(w*t); % Тестовый сигнал

V=EuFilter(U); % Фильтрация тестового сигнала

figure(1);

plot(t(1:1024),U(1:1024),'r');

hold on;

plot(t(1:1024),V(1:1024),'g');

grid on, xlabel('t (s)'), title('Signal');

 

Fd=1/Td; % Частота дискретизации 

SpectrU=fft(U,N);%Быстрое преобразование Фурье

SpectrU=2*SpectrU./N; % Нормировка спектра по амплитуде

SpectrU(1)=SpectrU(1)/2; % Нормировка постоянной составляющей в спектре

W=0:2*pi*Fd/N:1000;

figure(2);

subplot(2,1,1), plot(W,abs(SpectrU(1:length(W))),'r'), grid on, xlabel('w (Rad/s)'), title('MAGNITUDE');

hold on;

subplot(2,1,2), plot(W,angle(SpectrU(1:length(W)))*180/pi,'r'), grid on, xlabel('w (Rad/s)'), title('PHASE (deg)');

hold on;

SpectrV=fft(V,N); %Быстрое преобразование Фурье

SpectrV=2*SpectrV./N; % Нормировка спектра по амплитуде

SpectrV(1)=SpectrV(1)/2; % Нормировка постоянной составляющей в спектре

subplot(2,1,1), plot(W,abs(SpectrV(1:length(W))),'g');

subplot(2,1,2), plot(W,angle(SpectrV(1:length(W)))*180/pi,'g');


Рис. 33.  Тестовый сигнал U=sin(40·t) (обозначен красным), и сигнал V, полученный в результате фильтрации (обозначен зелёным)

Рис. 34.  Амплитудный и фазовый спектры тестового сигнала U=sin(40·t) (обозначены красным), и сигнала V, полученный в результате фильтрации (обозначены зелёным). Получены с помощью быстрого преобразования Фурье.

По спектрам определим амплитуды и фазы входного и выходного сигналов:

Амплитуда входного: A(U)≈1;

Амплитуда выходного: A(V)≈0,08;

Соотношение амплитуд:

Фаза входного: φ(U)≈-90º;

Фаза выходного: φ(V)≈-10º;

Разность фаз: φ(V)-φ(U)=-10º+90º=80º

Полученные результаты согласуются с АЧХ и ФЧХ фильтра для этой частоты. Она находится левее основной полосы пропускания.

  1.  Примем частоту ω=545 рад/с

Рис. 35.  Тестовый сигнал U=sin(545·t) (обозначен красным), и сигнал V, полученный в результате фильтрации (обозначен зелёным)

Рис. 36.  Амплитудный и фазовый спектры тестового сигнала U=sin(545·t) (обозначены красным), и сигнала V, полученный в результате фильтрации (обозначены зелёным).

По спектрам определим амплитуды и фазы входного и выходного сигналов:

Амплитуда входного: A(U)≈1;

Амплитуда выходного: A(V)≈0,473;

Соотношение амплитуд:

Фаза входного: φ(U)≈-90º;

Фаза выходного: φ(V)≈-87º;

Разность фаз: φ(V)-φ(U)=-87º+90º=3º

Полученные результаты согласуются с АЧХ и ФЧХ фильтра для этой частоты. Она находится в середине основной полосы пропускания.

  1.  Примем частоту ω=1000 рад/с

Рис. 37.  Тестовый сигнал U=sin(1000·t) (обозначен красным), и сигнал V, полученный в результате фильтрации (обозначен зелёным)

Рис. 38.  Амплитудный и фазовый спектры тестового сигнала U=sin(1000·t) (обозначены красным), и сигнала V, полученный в результате фильтрации (обозначены зелёным).

По спектрам определим амплитуды и фазы входного и выходного сигналов:

Амплитуда входного: A(U)≈1;

Амплитуда выходного: A(V)≈0,4;

Соотношение амплитуд:

Фаза входного: φ(U)≈-90º;

Фаза выходного: φ(V)≈-114º;

Разность фаз: φ(V)-φ(U)=-114º+90º=-24º

Полученные результаты согласуются с АЧХ и ФЧХ фильтра для этой частоты. Она находится в правой части основной полосы пропускания.

Подадим на вход фильтра равномерный белый шум:

Код Matlab:

U=rand(1,N);

Рис. 39.  Белый шум до фильтрации (обозначен красным) и после (обозначен зелёным)

Рис. 40.  Спектр белого шума до фильтрации (обозначен красным) и после (обозначен зелёным)

Подадим на вход сигнал, состоящий из трёх гармоник на частотах 40 рад/с, 545 рад/с 5000 рад/с и нормального белого шума.

Код Matlab:

clc;clear;

Td=0.00019096; % Период дискретизации

N=2^13; % Количество отсчётов

t=0:Td:(N-1)*Td;

w1=40;

w2=545;

w3=5000;

U=sin(w1*t)+sin(w2*t)+sin(w3*t)+randn(1,N);

V=EuFilter(U);

figure(1);

plot(t(1:1024),U(1:1024),'r');

hold on;

plot(t(1:1024),V(1:1024),'g','LineWidth',2);

grid on, xlabel('t (s)'), title('Signal');

 

Fd=1/Td; % Частота дискретизации 

SpectrU=fft(U,N);%Быстрое преобразование Фурье

SpectrU=2*SpectrU./N; % Нормировка спектра по амплитуде

SpectrU(1)=SpectrU(1)/2; % Нормировка постоянной составляющей в спектре

W=0:2*pi*Fd/N:5500;

figure(2);

subplot(2,1,1), plot(W,abs(SpectrU(1:length(W))),'r'), grid on, xlabel('w (Rad/s)'), title('MAGNITUDE');

hold on;

subplot(2,1,2), plot(W,angle(SpectrU(1:length(W)))*180/pi,'r'), grid on, xlabel('w (Rad/s)'), title('PHASE (deg)');

hold on;

SpectrV=fft(V,N); %Быстрое преобразование Фурье

SpectrV=2*SpectrV./N; % Нормировка спектра по амплитуде

SpectrV(1)=SpectrV(1)/2; % Нормировка постоянной составляющей в спектре

subplot(2,1,1), plot(W,abs(SpectrV(1:length(W))),'g','LineWidth',2);

subplot(2,1,2), plot(W,angle(SpectrV(1:length(W)))*180/pi,'g','LineWidth',2);

Рис. 41.  Тестовый сигнал U(t)= sin(40·t) + sin(545·t) +sin(5000·t)+n(t)  (обозначен красным), и сигнал V, полученный в результате фильтрации (обозначен зелёным)

Рис. 42.  Амплитудный спектр тестового сигнала U(t) до фильтрации (обозначен красным) и после (обозначен зелёным)

После прохождения фильтра существенное ослабление претерпела первая гармоника (с частотой 40 рад/с), находящаяся левее основной полосы пропускания фильтра. Также существенно подавлена третья гармоника (с частотой 5000 рад/с), находящаяся правее основной полосы. Кроме того,  на выходе наблюдается меньшее отношение шума к полезному сигналу.

PAGE   \* MERGEFORMAT43


 

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

77409. Понятие брака. Порядок его заключения 43.5 KB
  Понятие признаки правовая природа брака. Брак является одним из наиболее древних институтов. Основные признаки брака в РФ: любой брак это союз.
77413. Реорганизация юридических лиц 66.5 KB
  В учебной литературе вопросы реорганизации носят описательный характер. Проблемы возникающие из реорганизации практически не освещаются в научной литературе. В этом и есть принципиальное отличие реорганизации от ликвидации Реорганизация – процедура совокупность юридических действий определяющая переход в порядке правопреемства прав и обязанностей от одного или нескольких юридических лиц правопредшественников к другому или другим юридическим лицам правопреемникам связанный с прекращением правопредшественников и созданием...
77414. Ликвидация коммерческих организации 23.27 KB
  При добровольной учредители или орган юридического лица это решение может приниматься и по истечении срока на который создавалось юридическое лицо по достижению цели в случае если суд признал недействительной государственную регистрацию юридического лица. К ликвидационной комиссии переходят правомочия управления коммерческой организацией она выступает в суде от имени ликвидируемого юридического лица. Этот баланс утверждается органами юридического лица по соглашению с органами государственной регистрации. Ликвидация юридического лица 1.
77416. Полное товарищество 19.53 KB
  Закон специально подчеркивает право участника полного товарищества знакомиться со всей документацией по ведению дел товарищества даже в случае когда такой участник в соответствии с учредительным договором не уполномочен на ведение дел товарищества а отказ от этого права или его ограничение в том числе по соглашению участников товарищества ничтожны. Кроме того он вправе возмездно или безвозмездно передать свою долю в складочном капитале товарищества или ее часть как другому товарищу так и третьему лицу не участвующему в товариществе....