Изучение алгоритма цифровой сверки

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

x=[2 1 3];

h=[4 1 3 4];

y=conv(x,h)

y =

8 6 19 14 13 12.

Интеграл свертки осуществляет сглаживание входных сигналов и поэтому может быть использован для проектирования низкочастотных фильтров. В этом случае задают весовые коэффициенты функции w(nT_П). В пакете Matlab дискретная фильтрация выполняется командой filter:

y=filter(b,a,x);

где b – вектор нерекурсивной части фильтра (числителя функции передачи); а – вектор коэффициентов рекурсивной части фильтра (знаменателя функции передачи), х – входной сигнал. Возвращаемой величиной является вектор отсчетов входного сигнала фильтра.

a=[2 2 2 1];

b=[4 1 3]

y =

4.0000 -1.0000 6.5000

Между функциями filter и conv существуют два различия. Во-первых, функция convтребует задания импульсной характеристики фильтра в виде вектора. Поэтому импульсная характеристика должна иметь конечную длину. Фактически это означает, что с помощью функции conv можно реализовать только нерекурсивные фильтры. С помощью функции filter можно реализовать и рекурсивные фильтры, т.е. фильтры с бесконечной импульсной характеристикой. Во-вторых, функция filter возвращает результат, длина которого равна длине входного сигнала.

Повторим пример по расчету фильтра, используя функцию filter. Поскольку используется импульсная характеристика фильтра, то в передаточной функции фильтра (команда filter) в качестве коэффициента в нерекурсивной части фильтра должна быть задана его импульсная характеристика (вектор h), а в рекурсивной части фильтра вектор адолжен быть представлен скаляром, равным единице:

z= filter(h,1,x)

    z =

8 6 19.

Видно, что результат фильтрации совпадает в трех первых пунктах со значениями, полученными с помощью функции conv. Если бы входной сигнал х имел четыре отсчета, то совпадение наблюдалось бы в четырех точках, то есть векторz не содержит затухающих «хвостов». Различия в вычислениях объясняется тем, что функция conv предполагает завершенность входного сигнала, а функция filter завершенность входного сигнала не требует и позволяет, при необходимости, продолжать обработку поступающих на вход фильтра новых отсчетов.

Чтобы получить результат, совпадающий с результатом функции conv, входной вектор х нужно дополнить необходимым количеством нулей:

z1=filter(h,1,[x,0,0,0])

z1 =

8 6 19 14 13 12.

При фильтрации сигналов в ряде случаев возникает необходимость сохранения фазовых соотношений в исходном сигнале. Для этого требуется фильтр с линейной ФЧХ или иначе с постоянной фазовой задержкой. Это достигается в нерекурсивных фильтрах с симметричным набором коэффициентов b_i = b_m-i, которые определяются конечной импульсной характеристикой. Нерекурсивные фильтры обладают бесконечной импульсной характеристикой и поэтому принципиально не могут иметь линейную ФЧХ.

Для реализации фильтрации с нулевым фазовым сдвигом в Matlab разработана функция filtfilt. Ее название символизирует смысл выполняющих ею действий: сначала входной сигнал фильтруется как обычно, а затем производится второй проход фильтрации, при котором сигнал обрабатывается от конца к началу. Если заданный фильтр характеризуется передаточной функцией W( ), импульсной характеристикой h(nT_П) и АЧХ – W(ω), то фильтр, эквивалентный двум проходам обрабатываемого сигнала (прямой и обратный проход), будет характеризоваться передаточной функцией W(z)W(z^-1), АЧХ |W(ω)|² и нулевой ФЧХ. Импульсная характеристика фильтра с нулевым фазовым сдвигом представляет собой корреляционную функцию исходной импульсной характеристики:

,

Синтаксис вызова функцииfiltfilt:

y=filtfilt(b,a,x).

Приведем пример, демонстрирующий компенсацию фазового сдвига с помощью функции filtfilt. Создадим сигнал в виде прямоугольного импульса, дополненного с обеих сторон нулями:

S=[zeros(50,1);ones(100,1);zeros(50,1)].

В качестве фильтра используем ФНЧ Чебышева первого рода с достаточно сильными пульсациями АЧХ в полосе пропускания (3 дБ) и низкой частотой среза (0,05 частоты Найквиста), чтобы выходной импульс был заметно искажен [b,a]=cheby1(4,2,0,05).

Обработайте сигнал с помощью filter и filtfilt и постройте графики входного и выходного сигналов

s=[zeros(50,1);ones(100,1);zeros(50,1)]; %Входной сигнал

[b,a]=cheby1(4,2,0.05);      %Коэффициенты фильтра Чебышева

b =

1.0e-004 *

0.0589 0.2356 0.3535 0.2356 0.0589

a =

1.0000 -3.8634 5.6223 -3.6525 0.8936

s1=filter(b,a,s);           %Результат фильтрации

s2=filtfilt(b,a,s);              %Результат фильтрации

figure(1);

plot(s);

hold on;

plot(s1,'--'),grid on;      %Графики выходных сигналов

plot(s2,'-.'),grid on;      %Графики выходных сигналов

Hold on.

Интеграл свертки используется при проектировании цифровых фильтров. Для этого требуется задать весовые коэффициенты или, что одно и то же, импульсную характеристику, которую можно получить путем подачи на вход фильтра единичного сигнала, дополненного некоторым количеством нулей filter(b,a,[1,zeros(1,n)]).

n=[2 1 3 4 1 3 4]

ans =

1.0e-003 *

0.0059 0.0463 0.1812

По известным коэффициентам передаточных функций коэффициенты импульсной характеристики определяются как n=impz(b,a,n), гдеb, a – полиномы числителя и знаменателя; n – вектор отсчетов импульсной характеристики. Если он не указан, то Matlab автоматически выбирает число отсчетов в зависимости от поведения импульсной характеристики. На рис. 8.2 представлена импульсная характеристика для фильтра Баттерворта пятого порядка.

n=impz(b,a,n)

n =

0.0002

0.0000

0.0005

0.0010

0.0000

0.0005

0.0010

[b,a]=butter(4,0.05); %Фильтр Баттерворта пятого порядка

%с частотой среза 0.2

b =

1.0e-003 *

0.0312 0.1250 0.1874 0.1250 0.0312

a =

1.0000 -3.5897 4.8513 -2.9241 0.6630

figure(2); impz(b,a),grid on %Импульсная характеристика фильтра.

В m-файле задан дискретный интервал времени и по пяти точкам определен интеграла свертки для входного сигнала, состоящего из двух гармонических сигналов: x₁ и x₂.

t=[0.1:0.1:0.8];     %Временной интервал входного сигнала

f=0.1;             %Частота задающего сигнала

x1=5*sin(2*pi*f*t); %Задающий сигнал

f1=4;              %Частота сигнала помехи

x2=0.3*sin(10*pi*f1*t); %Сигнал помехи

x3=x1+x2           %Сумма сигналов

x3 =

0.3140 0.6267 0.9369 1.2434 1.5451 1.8406 2.1289 2.4088

x22=[2 2 2 1 1 2 1 3]; %Импульсная характеристика фильтра

y1=filter(x3,1,x22)/5 %Выходной сигнал фильтра

y1 =

0.1256 0.3762 0.7510 1.1856 1.6783 2.2900 2.9554 3.7975