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


 

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

23051. Ознайомлення з основними можливостями пакета програм автоматизованого проектування електронних схем MicroSim PSPICE 8.0 1.35 MB
  Система автоматизованого проектування MicroSim PSPICE використовує один з найбільш вдалих кодів схемотехнічного моделювання SPICE Simulation Program with Integrated Circuit Emphasis який був розроблений на початку 70х років фахівцями Каліфорнійського університету США. Фактично зазначений код став стандартним для моделювання електронних схем і застосовується також у інших відомих системах моделювання схем зокрема MicroCap а вхідний формат мови завдань SPICE підтримується практично усіма пакетами автоматизованого проектування електронних...
23052. Електронний ключ на біполярному транзисторі 482 KB
  Каскад виконує логічну операцію заперечення оскільки високий рівень напруги на вході забезпечує введення транзистора у режим насичення коли напруга на навантаженні буде низькою. При введенні наведеної вище схеми дослідження ключового каскаду застосовуються джерела сталої напруги живлення VCC та імпульсної вхідної напруги VIN. Перелічимо основні параметри даних джерел: Як джерело сталої напруги живлення застосовується стандартна модель VSRC що міститься у бібліотеці source. Основними є такі її параметри: DC стала напруга що її виробляє...
23053. Електронні ключі на МДН-транзисторах 1.07 MB
  Вважайте що напруга живлення дорівнює 10 В амплітуда вхідного імпульсу 10 В тривалість цього імпульсу 500 нс його період 1000 нс. Тривалості фронту і спаду імпульсу задайте дуже малими наприклад по 0. Поясніть зміни у тривалості спаду вихідного імпульсу та рівні напруги логічного нуля на виході. Параметри джерел вважайте такими: напруга живлення 20 В амплітуда вхідного імпульсу 20 В тривалість цього імпульсу 500 нс його період 1000 нс.
23054. Базовий елемент транзисторно-транзисторної логіки (ТТЛ) 1016 KB
  Насправді опором навантаження для виходу ТТЛсхеми звичайно є вхідний опір наступної ТТЛсхеми. Оскільки у реальних ситуаціях на один вихід треба підєднувати досить багато входів важливим є такий параметр схеми як навантажувальна здатність тобто максимальна кількість входів яку можна навантажити на вихід без втрати працездатності схеми. Оскільки транзистори в даній схемі працюють у режимах насичення та відсікання має місце досит значна інерційність схеми потрібен певний час для переведення транзисторів з одного граничного стану в...
23055. Моделювання цифрових логічних схем 178.5 KB
  Перелічимо деякі логічні ІМС 74ї серії: 74x00 базовий елемент 2ІНЕ 74x10 логічний елемент 3ІНЕ 74x20 логічний елемент 4ІНЕ 74x30 логічний елемент 8ІНЕ 74x02 логічний елемент 2АБОНЕ 74x27 логічний елемент 3АБОНЕ 74x08 логічний елемент 2І 74x32 логічний елемент 2АБО 74x04 інвертор логічний елемент НЕ 74x51 логічний елемент 2І2АБОНЕ 74x86 логічний елемент Виключне АБО на 2 входи Пакет OrCAD дозволяє провести суто цифрове моделювання для даного вузла схеми якщо до цього вузла підєднані лише входи та виходи...
23056. Роль та повноваження Ради національної безпеки України в системі забезпечення національної безпеки 40.5 KB
  Роль та повноваження Ради національної безпеки України в системі забезпечення національної безпеки. Це випливає із Закону України Про Раду національної безпеки.подає пропозиції Президентові України щодо .визначення стратегічних національних інтересів України концептуальних підходів та напрямів забезпечення національної безпеки і оборони у політичній економічній соціальній військовій науковотехнологічній екологічній інформаційній та інших сферах проектів державних програм доктрин .
23057. Стратегія національної безпеки України: основні положення та проблеми формування 38 KB
  Стратегія національної безпеки України: основні положення та проблеми формування. ТЕМА: Стратегія національної безпеки Стратегія чітко вивірений шлях і напрям досягнення мети. Стратегія національної безпеки це система державнополітичних рішень головних напрямків діяльності у сфері безпеки послідовна реалізація яких забезпечує досягнення мети головна лінія що дозволяє забезпечити безпеку на певний період спрямована на досягнення середньострокових та довгострокових інтересів. В НАТО змінено пріоритет безпеки на поширення зони...
23058. Загрози національній безпеці України у політичній сфері та шляхи їх нейтралізації 40 KB
  Загрози національній безпеці України у політичній сфері та шляхи їх нейтралізації. Загрози у внутрішньополітичній сфері: порушення з боку органів державної влади та органів місцевого самоврядування Конституції і законів України прав і свобод людини і громадянина в тому числі при проведенні виборчих кампаній недостатня ефективність контролю за дотриманням вимог Конституції і виконанням законів України; можливість виникнення конфліктів у сфері міжетнічних і міжконфесійних відносин радикалізації та проявів екстремізму в діяльності деяких...
23059. Загрози національній безпеці України в економічній сфері та шляхи їх подолання 38 KB
  Загрози в економічній сфері: істотне скорочення внутрішнього валового продукту зниження інвестиційної та інноваційної активності і науковотехнічного та технологічного потенціалу скорочення досліджень на стратегічно важливих напрямах інноваційного розвитку; ослаблення системи державного регулювання і контролю у сфері економіки; нестабільність у правовому регулюванні відносин у сфері економіки в тому числі фінансової фіскальної політики держави; відсутність ефективної програми запобігання фінансовим кризам; зростання...