Оптимальная линейная фильтрация непрерывных сигналов

Рассмотрим теорию оптимальной линейной фильтрации стационарных процессов – фильтр Колмогорова-Винера или оптимальный линейный фильтр (ОЛФ).

Пусть на входе линейного фильтра с передаточной функцией H (j w) действует сумма полезного сигнала s (t ) и помехи n (t ). Отклик фильтра на это действие – восстановленный полезный сигнал (оценка сигнала s (t )). Будем считать, что s (t ) и n (t ) – стационарные взаимнонекоррелированные процессы с известными спектральными плотностями мощности (СПМ) G s (f ) и G n (f ). Нужно найти такую функцию H (j w), которая обеспечивает минимум среднего квадрата ошибки восстановления сигнала

. (22.1)

Иначе говоря, критерием оптимальности фильтра является минимум среднего квадрата ошибки восстановления сигнала. В такой постановке задача была решена А.Н. Колмогоровым (1939 г.) для дискретных случайных последовательностей и Н. Винером (1941 г.) для непрерывных процессов. Поэтому оптимальный (в указанном смысле) линейный фильтр называется фильтром Колмогорова-Винера.

Действующие на входе фильтра сигнал s (t ) и помеха n (t ) проходят через фильтр независимо и создают на выходе фильтра соответственно фильтрованные сигнал s вых (t ) и помеху n вых (t ). С учетом этого ошибка восстановления запишется

Слагаемое Ds (t ) отражает составляющую ошибки, обусловленную линейными искажениями полезного сигнала фильтром. Средний квадрат ошибки запишется

. (22.3)

Величина линейных искажений полезного сигнала фильтром зависит от степени отличия АЧХ фильтра от постоянной величины и степени отличия ФЧХ от линейной зависимости. Средний квадрат шума на выходе фильтра зависит только от АЧХ фильтра. Для того чтобы линейные искажения полезного сигнала были минимальными, ФЧХ фильтра должна быть линейной

j(w) = –wt 0 , (22.4)

где t 0 – задержка сигнала в фильтре. Ясно, что с учетом задержки соотношение (22.1) имеет вид . Это уточнение не влияет на критерий оптимальности, поскольку в системах связи и вещания ожидаемая задержка сигнала в фильтре несущественная.

Перейдем к определению АЧХ фильтра. Для этого определим спектральные плотности мощностей левой и правой частей соотношения (22.2)

Выразим СПМ помехи на выходе фильтра через СПМ помехи n (t ) и искомую АЧХ фильтра:

. (22.6)

По определению СГП ергодичного процесса

, (22.7)

где S D s (w) – амплитудный спектр ошибки за счет линейных искажений Ds (t );

Т – длительность реализации Ds (t ).

Поскольку амплитудный спектр сигнала s в ых (t ) определяется как S (w)×H (w), где S (w) – амплитудный спектр сигнала s (t ), то

Переходя от амплитудных спектров к СПМ, получим

После подстановки соотношений (22.6) и (22.9) в (22.5) получим

Средний квадрат ошибки восстановления (средняя мощность) вычисляется

. (22.11)

Поскольку функция G e (w) ³ 0 на всех частотах, то, обеспечив min G e (w) на всех частотах, достигнем минимума величины . Искомую АЧХ H (w) определим из условия экстремума функции G e (w):

После решения уравнения (22.13) получим выражение для АЧХ фильтра

. (22.14)

На рис. 22.1 иллюстрируется АЧХ ОЛФ, определенная соотношением (22.14).

Из рис. 22.1 видно особенности АЧХ ОЛФ:

На частотах, где G n (f ) = 0, значение АЧХ H (f ) = 1 – в этих областях частот фильтр не должен вносить искажений;

На частотах, где G s (f ) = 0, значение АЧХ H (f ) = 0 – в этих областях частот фильтр должен полностью ослаблять составляющие помехи;

На частотах, на которых G s (f ) = G n (f ), АЧХ H (f ) = 0,5;

На других частотах значения АЧХ определяются вычислениями по формуле (22.14).

Подставим выражение (22.14) в соотношение (22.10) для определения СПМ ошибки:

(22.15)

При подстановке соотношения (22.15) в выражение (22.11) можно вычислить средний квадрат ошибки восстановления сигнала .

Из (22.15) видно, что ошибка равна нулю только в том случае, когда G s (f )G n (f ) = 0, т.е. когда спектры сигнала и помехи не перекрываются (хотя бы один с сомножителей равен нулю).

Перепишем соотношение (22.14) в виде

. (22.16)

Из последнего соотношения видно, что коэффициент передачи оптимального фильтра на каждой из частот тем меньше, чем больше отношение G n (f )/G s (f ) на этой частоте.

Следует отметить, что оптимальные линейные фильтры, обеспечивающие минимум ошибки , существенным образом отличаются от согласованных фильтров, рассмотренных ранее. Если основное назначение рассмотренных здесь фильтров состоит в наилучшем воспроизведении формы сигнала, то задача согласованных фильтров состоит в формировании максимального отношения сигнал/шум в момент отсчета.

При использовании ОЛФ в аналоговых системах связи и вещания выявляется такая особенность. Имеет место высокое отношение спектральных плотностей сигнала и шума: G s (f )/G n (f ) >> 1. Выражение для АЧХ ОЛФ (22.14) в случае полосовых сигналов переходит в следующее

(22.17)

где f min и f max – граничные частоты спектра сигнала. В случае НЧ сигналов

(22.18)

где F max – максимальная частота спектра сигнала.

Таким образом, оптимальные линейные фильтры в системах связи и вещания имеют П-образную АЧХ.

Контрольные вопросы

1. Сформулируйте критерий оптимальности ОЛФ.

2. При каком условии ошибку восстановления сигнала можно свести к нулю?

3. Объясните отличие ОЛФ и согласованного фильтра с точки зрения ослабления помех.

23. Сравнение помехоустойчивости оптимальных
демодуляторов сигналов аналоговых видов модуляции

Мы выяснили, что в условиях слабых помех демодулятор должен содержать: фильтр додетекторной обработки, детектор, фильтр последетекторной обработки. Для того чтобы демодулятор был оптимальным, фильтры должны быть оптимальными. В условиях слабой помехи АЧХ фильтров должны быть П-образными:

– фильтр додетекторной обработки – полосовой фильтр, граничные частоты полосы пропускания которого совпадают с граничными частотами спектра модулированного сигнала;

– фильтр последетекторной обработки – ФНЧ, частота среза которого совпадает с максимальной частотой спектра первичного сигнала F max .

Помехоустойчивость определим в условиях действия АБГШ. Анализ помехоустойчивости состоит в определении выигрыша демодулятора в отношении сигнал/шум

Для определения выигрыша нужно определить 4 величины: P b , P e , P s , P n .

Оптимальныйдемодулятор сигнала балансовой модуляции . Математически сигнал БМ записывается

s БМ (t ) = A 0 b (t )cos2pf 0 t . (23.2)

Полосовой фильтр имеет полосу пропускания DF БМ = 2F max . Для детектирования БМ сигнала необходимо использовать синхронный детектор (рис. 23.1). ФНЧ, являющийся обязательным элементом схемы синхронного детектора, используется как фильтр последетекторной обработки, т.е. частота среза фильтра равна F max , а АЧХ в полосе пропускания постоянная и равна 1.

P s = = = 0,5 P b . (23.3)

Средняя мощность помехи на входе демодулятора P n .

На выходе перемножителя за счет сигнала имеем

A 0 b (t )cos2pf 0 t ×2cos2pf 0 t = A 0 b (t ) + A 0 b (t ) cos2p2f 0 t . (23.4)

ФНЧ пропускает первый компонент, а второй ослабляет. С учетом умножения на 1/A 0 на выходе демодулятора получим b (t P b .

Помеху на входе синхронного детектора (как полосовой процесс) представим квадратурными составляющими

n (t ) = N C (t )cos2pf 0 t + N S (t )sin2pf 0 t . (23.5)

Мощность помехи делится поровну между квадратурными составляющими, мощность каждой из них P n /2. Квадратурная составляющая помехи не проходит через синхронный детектор, и на выходе демодулятора получим

e(t ) = N C (t )/A 0 . (23.6)

Поскольку = P n /2, а левая часть равенства равняется 0,5 , то = P n , а

P e = P n / . (23.7)

g БМ = = 2. (23.8)

Оптимальныйдемодулятор сигнала однополосной модуляции . Математически сигнал ОМ записывается

s ОМ (t ) = A 0 b (t )cos2pf 0 t ± A 0 sin2pf 0 t . (23.9)

Оптимальный демодулятор сигнала ОМ выполняется по схеме, приведенной на рис. 23.1. Полосовой фильтр имеет полосу пропускания DF ОМ = F max .

Средняя мощность модулированного сигнала

P s = = P b . (23.10)

Средняя мощность помехи на входе синхронного детектора P n .

Синхронный детектор не пропускает квадратурную составляющую сигнала ОМ, поэтому на основе анализа демодуляции сигнала БМ на выходе демодулятора сигнала ОМ получим сигнал b (t ). Его средняя мощность равняется P b . Прохождение помехи через синхронный детектор проанализирован выше и получено значение мощности помехи на выходе демодулятора (23.7).

Определим выигрыш демодулятора

g ОМ = = 1. (23.11)

Оптимальныйдемодулятор сигнала амплитудной модуляции на базесинхронного детектора . Математически сигнал АМ записывается

s АМ (t ) = A 0 (1 + m АМ ×b (t ))cos2pf 0 t . (23.12)

Схема демодулятора сигнала АМ приведена на рис. 11.5. Полосовой фильтр имеет полосу пропускания DF АМ = 2F max . Исходя из приведенного выше анализа, очевидно, что на выходе ФНЧ за счет сигнала получим A 0 + A 0 m АМ ×b (t ). Фильтр верхних частот ослабляет постоянную составляющую A 0 и пропускает вторую составляющую A 0 m АМ ×b (t ).

Средняя мощность модулированного сигнала

P s = = 0,5 (1 + m АМ P b ). (23.13)

Для дальнейшего анализа удобно учесть (см. модуль 1), что P b = 1/ , где K A – коэффициент амплитуды сигнала b (t ).

Прохождение шума через синхронный детектор проанализирован выше. Мощность шума на выходе демодулятора сигнала АМ

P e = . (23.14)

Определим выигрыш демодулятора

= . (23.15)

Демодулятор сигнала амплитудной модуляции на базедетектора огибающей . Схема такого демодулятора приведена на рис. 23.3. Она отличается от схемы демодулятора на рис. 23.2 типом амплитудного детектора – синхронный детектор заменен детектором огибающей с целью упрощения схемы демодулятора. Выходной сигнал детектора огибающей пропорциональный огибающей входного сигнала .

Поэтому помеху на выходе демодулятора создают как косинусная, так и синусная составляющие. Мощность помехи e(t ) будет вдвое большей, чем в случае синхронного детектора: Выигрыш демодулятора будет вдвое меньшим:

(23.16)

Оптимальныйдемодулятор сигнала фазовой модуляции . Математически сигнал ФМ записывается

s ФМ (t ) = A 0 cos(2pf 0 t + Δj д ∙b (t )), (23.17)

где Δj д – девиация фазы сигнала, которую часто называют индексом фазовой модуляции m ФМ.

Схема оптимального демодулятора сигнала ФМ приведена на рис. 23.4: ФНЧ1 и ФНЧ2 имеют частоты среза F max (m ФМ + 1); ФНЧ3 – фильтр последетекторной обработки с частотой среза F max ; АЧХ фильтров постоянная в полосах пропускания и равна 1; ФД – фазовый детектор.

Средняя мощность модулированного сигнала

. (23.18)

Средняя мощность исходного сигнала P b = 1/ .

Средняя мощность помехи на входе ФД P n . Прохождение помехи через ФД анализируют при отсутствии модулирующего сигнала, т.е. b (t ) º 0, и s (t ) = A 0 cos2pf 0 t . Помеху представляют квадратурными составляющими в виде (23.5). Тогда

В условиях слабой помехи |N c (t )| << A 0 , и помеха на выходе ФД N s (t )/(A 0 m ФМ), а ее мощность равна P n /(A 0 m ФМ) 2 .

Полоса пропускания ФНЧ3 в раз меньше, чем полоса пропускания ФНЧ2 и ширина спектра помехи N s (t ).

Поскольку спектр помехи равномерный, то мощность помехи e(t ) в m ФМ + 1 раз меньше, чем мощность помехи N s (t )/(A 0 m ФМ) и определяется

. (23.20)

Определим выигрыш демодулятора

Оптимальныйдемодулятор сигнала частотной модуляции. Математически сигнал ЧМ записывается

где Δf д – девиация частоты. Для последующего изложения удобно использовать индекс ЧМ m ЧМ = Δf д /F max .

Схема оптимального демодулятора сигнала ЧМ приведена на рис. 23.5. Схема отличается от схемы демодулятора ФМ сигнала (рис. 23.4) наличием дифференциатора; ЧД – частотный детектор.

Средняя мощность модулированного сигнала

. (23.23)

Средняя мощность сигнала на выходе демодулятора P b =1/ .

Средняя мощность помехи на входе ЧД P n .

Прохождение помехи через ЧД аналогично прохождению помехи через ФД, следует рассмотреть прохождение помехи N s (t )/(A 0 2pΔf д) через дифференциатор. Поскольку помеха N s (t ) – квазибелый шум в полосе частот F max (m ЧМ + 1), то спектральная плотность мощности этой помехи на входе дифференциатора

Поскольку передаточная функция дифференциатора j ω, то спектральная плотность мощности помехи на выходе дифференциатора

Определим мощность помехи e(t )

Определим выигрыш демодулятора

При анализе мы выявили, что спектральная плотность мощности помехи e(t ) имеет параболическую зависимость – формула (23.24):

G e (f ) = kf 2 , (23.26)

где k – коэффициент пропорциональности. Эта особенность спектра часто учитывается при разработке систем передачи с ЧМ.

Сравнениеаналоговых систем передачи. Основными параметрами, по которым сравниваются системы передачи, является выигрыш демодулятора в отношении сигнал/шум g и коэффициент расширения полосы частот при модуляции a = ΔF s /F max. Для рассмотренных методов модуляции эти параметры сведены в табл. 23.1.

Проведем сравнение числовых значений параметров при некоторых исходных данных: K А = 5; m ЧМ = m ФМ = 6; m АМ = 1.

Вычисления дают: g АМ = 0,038; g БМ = 2; g ОМ = 1; g ЧМ = 60,5; g ФМ = 20,2.

Сравнение числовых значений выигрыша показывает, что самой низкая помехоустойчивость присуща системе передачи с АМ: выигрыш g АМ << 1, что логически назвать проигрышем. Однако АМ используется в системах радиовещания, где этот недостаток компенсируется простотой демодулятора на основе детектора огибающей (огромное количество более простых радиоприемников и один радиопередатчик с большей мощностью, чем при использовании БМ и ОМ).

Наибольшая помехоустойчивость присуща системе передаче с ЧМ. «Платой» за высокую помехоустойчивость является широкая полоса частот сигнала. Так, при F max = 3,4 кГц ΔF ЧМ = 47,6 кГц, в то время как полоса частот сигнала ОМ ΔF ОМ = 3,4 кГц.

Таблица 23.1 – Основные параметры аналоговых систем передачи

Метод модуляции	g	a	Примечания
АМ			Синхронный детектор
		Детектор огибающей
БМ
ОМ
ЧМ	×a	2(m ЧМ + 1)	r вх > r пр
ФМ	×a	2(m ФМ + 1)	r вх > r пр

Порогпомехоустойчивости демодулятора сигнала ЧМ. Из соотношения для выигрыша демодулятора ЧМ (23.25) вытекает, что, чем больше индекс m ЧМ, тем больше выигрыш (правда, ценой увеличения полосы частот сигнала). Может показаться, что это дает возможность работать демодулятору со слабым сигналом (низким отношением сигнал/шум). Но, когда отношение сигнал/шум на входе демодулятора r вх меньше порогового отношения сигнал/шум r пр, то выигрыш демодулятора резко уменьшается (рис. 23.6). Такое явление резкого уменьшения величины выигрыша называют порогом помехоустойчивости приема сигнала ЧМ.

Пороговое отношение сигнал/шум r пр несколько зависит от значения m ЧМ (рис. 23.6). Считают, что демодулятор по схеме стандартного частотного детектора характеризуется ориентировочным значением r пр = 10. Область значений r вх, когда r вх < r пр, – это нерабочая область.

Были предложены так называемые порогопонижающие схемы демодуляторов сигналов ЧМ, которые получили названия:

Демодулятор со следящим фильтром;

Демодулятор с обратной связью по частоте;

Демодулятор на основе синхронно-фазового детектора.

Схемы этих демодуляторов описаны в специальной литературе. Демодуляторы, которые выполнены по таким схемам, характеризуются пороговым отношением сигнал/шум r пр = 5...7 дБ (в зависимости от исходных данных на систему передачи). Снижение r пр разрешает:

1) работать демодулятору с более низким отношением сигнал/шум;

2) увеличить выигрыш, так как

и если допустить уменьшение отношения сигнал/шум r вх, то можно увеличить индекс m ЧМ, а увеличение m ЧМ приводит к увеличению выигрыша.

Литература

Основная

1. Стеклов В.К., Беркман Л.Н. Теорія електричного зв’язку: Підручник для ВНЗ під редакцією В.К. Стеклова. – К.: Техніка, 2006.

2. Теория электрической связи: Учебник для вузов / А.Г. Зюко, Д.Д. Кловский, В.И. Коржик, М.В. Назаров; Под ред. Д.Д. Кловского. – М.: Радио и связь, 1998.

Условие оптимальности фильтра. Фильтр Колмогорова-Винера является оптимальным фильтром формирования из входного сигнала x(t) выходного сигнала z(t) при известной форме полезного сигнала s(t), который содержится во входном сигнале в сумме с шумами. В качестве критерия его оптимизации используется среднее квадратическое отклонение сигнала y(t) на выходе фильтра от заданной формы сигнала z(t). Подставим уравнение свертки (12.2.1) в раскрытой форме весового суммирования в выражение (12.2.2") и получим отклонение e 2 выходного сигнала y(k) = h(n)③x(k-n) от заданной формы выходного сигнала z(k) по всем точкам массива данных:

. (12.3.1)

В частном случае воспроизведения формы полезного сигнала в качестве функции z(k) принимается функция s(k). Минимум выражения (12.3.1) определяет значения коэффициентов h(n) оптимального фильтра. Стоит сказать, что для нахождения их значений продифференцируем выражение (12.3.1) по коэффициентам фильтра и приравняем полученные уравнения нулю. В итоге получаем:

где - корреляционная функция входного сигнала, - взаимная корреляционная функция сигналов z(k) и x(k). Отсюда:

h n R(m-n) = B(m), n = m = 0,1,2, ... , M. (12.3.2)

В краткой форме символической записи:

h(n) ③ R(m-n) = B(m). (12.3.3)

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

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

m=0: h o R(0)+ h 1 R(1)+ h 2 R(2)+ h 3 R(3)+ ...+ h M R(M) = B(0), (12.3.3")

m=1: h o R(1)+ h 1 R(0)+ h 2 R(1)+ h 3 R(2)+ ...+ h M R(M-1) = B(1),

m=2: h o R(2)+ h 1 R(1)+ h 2 R(0)+ h 3 R(1)+ ...+ h M R(M-2) = B(2),

..............................................................................................................

m=M: h o R(M)+ h 1 R(M-1)+ h 2 R(M-2)+ .... + h M R(0) = B(M).

Решение данной системы уравнений относительно значений h i дает искомый оператор фильтра.

При расчете каузальных (односторонних) фильтров выходной сигнал z(t) следует задавать со сдвигом вправо по оси координат таким образом, чтобы значимая часть функции взаимной корреляции B(m) полностью располагалась в правой части системы уравнений (12.3.3").

Отметим, что R(m) = R s (m)+R q (m), где R s - функция автокорреляции сигнала, R q - функция автокорреляции шума, а B(m) = B zs (m)+B zq (m), где B zs - функция взаимной корреляции сигналов z(k) и s(k), B zq - функция взаимной корреляции сигнала z(k) и помех q(k). Подставляя данные выражения в (12.3.3), получаем:

h(n) ③ = B zs (m)+B zq (m). (12.3.4)

Частотная характеристика фильтра находится преобразованием Фурье левой и правой части уравнения (12.3.4):

H(w) = W zs (w)+W zq (w),

H(w) = / , (12.3.5)

где W s (w) ó R s (m) и W q (w) ó R q (m) - энергетические спектры (плотности мощности) сигнала и помех, W zs (w) ó B zs (m) - взаимный энергетический спектр входного и выходного сигналов, W zq (w) ó B zq (m) - взаимный энергетический спектр выходного сигнала и помех.

В геофизической практике обычно имеет место статистическая независимость полезного сигнала, а, следовательно, и сигнала z(k), от шумов, при этом B zq = 0 и фильтр называют оптимальным по сглаживанию шумов при заданной форме выходного сигнала:

H(w) = W zs (w) / , (12.3.6)

Фильтр (12.3.6) оптимален в том смысле, что максимизирует отношение мощности сигнала к мощности шума по всему интервалу сигнала, но не в каждой индивидуальной точке.

Выражения (12.3.5-12.3.6), как правило, всегда имеют решение. При этом это не означает возможность формирования фильтром любой заданной формы выходного сигнала. Из чисто практических соображений можно сразу предполагать, что если спектр заданного сигнала z(t) больше значимой части спектра полезного сигнала s(t), то оператор фильтра попытается сформировать требуемые высокие частоты заданного сигнала из незначимых частот спектра полезного сигнала, что может потребовать огромных коэффициентов усиления на этих частотах, под действие которых попадут и частотные составляющие шумов. Результат такой операции непредсказуем. Эти практические соображения можно распространить и на все частотные соотношения входного и выходного сигналов линейных фильтров: значимые гармоники спектров выходных сигналов должны формироваться из значимых гармоник спектров входных сигналов.

В случае если заданная форма сигнала z(k) совпадает с формой полезного сигнала s(k), то B(m) = B ss = R s и фильтр называют фильтром воспроизведения полезного сигнала :

H(w) = W s (w) / , (12.3.7)

Выражения (12.3.6-12.3.7) достаточно наглядно демонстрируют физический смысл формирования передаточной функции фильтра. При воспроизведении сигнала частотная функция взаимной корреляции входного сигнала с выходным W zs (плотность взаимной мощности) повторяет частотную функцию автокорреляции W s (плотность мощности сигнала). Плотность мощности статистических шумов W q распределена по частотному диапазону равномерно, в отличие от плотности мощности сигнала W s , которая, в зависимости от формы сигнала, может занимать любые частотные интервалы спектрального диапазона. Частотная передаточная функция фильтра воспроизведения сигнала формируется отношением W s (w)/. На частотах, где сосредоточена основная энергия сигнала, имеет место W s (w)>>W q (w) и H(w) Þ 1 (как минимум, больше 0.5). Там, где значение W s (w) становится меньше W q , коэффициент передачи фильтра становится меньше 0.5, и в пределе H(w)=0 на всех частотах, где полностью отсутствуют частотные составляющие сигнала.

Аналогичный процесс имеет место и при формировании произвольного сигнала z(t) на выходе фильтра, только в этом случае из частот входного сигнала устанавливаются на выделение и усиление частоты взаимной мощности входного и выходного сигнала, необходимые для формирования сигнала z(t), причем коэффициент на этих частотах может быть много больше 1, а подавляться могут не только шумы, но и частоты основного сигнала, если их нет в сигнале z(t).

Τᴀᴋᴎᴍ ᴏϬᴩᴀᴈᴏᴍ, оптимальные фильтры учитывают особенности спектрального состава сигналов и способны формировать передаточные функции любой сложности на выделение полезных частот сигналов из любых диапазонов спектра с максимальных подавлением шумов на всех частотах спектрального диапазона, не содержащих полезных сигналов, при этом границы усиления-подавления устанавливаются автоматически по заданному уровню шумов.

Задание мощности шумов. Следует внимательно относиться к заданию функции шумов Wq(w). При полной неопределенности шума крайне важно, как минимум, выполнять оценку его дисперсии s 2 и распространять на весь частотный диапазон с соответствующей нормировкой на его величину (2Wq(w) dw = s 2), ᴛ.ᴇ. считать его белым шумом. При известной функции полезного сигнала в зарегистрированной реализации значение дисперсии шумов в первом приближении может быть оценено по разности дисперсий реализации и функции полезного сигнала. Можно выполнить и выделение шумов из входного сигнала в отдельный шумовой массив, к примеру, вейвлетным преобразованием. При этом использовать выделенный шум непосредственно для вычисления функции Wq(w) допустимо только по достаточно представительной реализации при условии стационарности и эргодичности шума. В противном случае функция Wq(w) будет отображать только распределение шумов в зарегистрированной реализации сигнала, а соответственно фильтр будет оптимален только к этой реализации, что не гарантирует его оптимальности к любой другой реализации. Но для обработки единичной зарегистрированной реализации сигнала такой метод не только вполне допустим, но и может существенно повысить точность формирования выходного сигнала.

Эффективность фильтра. Из выражений (12.3.5-12.3.7) следует, что с позиции минимального искажения полезного сигнала при максимальном подавлении шумов фильтр Колмогорова-Винера эффективен в тем большей степени, чем больше отношение сигнал/шум на входе фильтра. В пределе, при W q (w)<>W s (w) имеем H(w) Þ 0 и сигнал будет сильно искажен, но никакой другой фильтр лучшего результата обеспечить не сможет.

Пример. Расчет оптимального фильтра воспроизведения сигнала. Выполняется в среде Mathcad.

Форма входного сигнала считается известной и близка к гауссовой. На входной сигнал наложен статистический шум с равномерным распределением мощности по всему частотному диапазону (белый шум), некоррелированный с сигналом, и функцию Wzq принимаем равной нулю. Стоит сказать, что для наглядного просмотра связи параметров фильтра с параметрами сигнала задаем модели сигналов в двух вариантах:

K:= 1000 k:= 0 .. K A:= 50

s1 k:= A·exp[-0.0005·(k-500) 2 ] s2 k:= A·exp[-0.00003·(k-500) 2 ] Ü информационные сигналы

Q:= 30 q k:= rnd(Q) – Q/2 x1 k:= s1 k + q k x2 k:= s2 k + q k Ü входные сигналы

Рис. 12.3.1. Модельные сигналы.

В качестве выходных сигналов задаем те же самые функции s1 и s2. Быстрым преобразованием Фурье вычисляем спектры сигналов и формируем спектры плотности мощности:

S1:= CFFT(s1) S2:= CFFT(s2) Q:= CFFT(q) Ü спектры сигналов

Ü спектры мощности

Ds1:= var(s1) Ds2:= var(s2) Dx1:= var(x1) Dx2:= var(x2) Dq:= var(q) Ü дисперсии

Ds1 = 124.308 Ds2 = 310.264 Dx1 = 202.865 Dx2 = 386.78 Dq = 79.038 Ü информация

mean(Wq) = 0.079 Wq1:= (Dx1 – Ds1)/(K+1) Wq1 = 0.078 Ü информация

Wq2:= (Dx2 – Ds2)/(K+1) Wq2 = 0.076 Ü информация

Wq k:= Wq1 Ü замена на постоянное распределение

Для воспроизведения сигналов вычисления функций Wzs не требуется, т.к. Wzs = Ws. Вычисление Wq также имеет только информативный характер.

Передаточные функции оптимальных фильтров (приведены на рис. 12.3.2):

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

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

Как следует из рисунка 12.3.2, для плавных монотонных функций, основная энергия которых сосредоточена в низкочастотной части спектра, передаточные функции оптимальных фильтров, по существу, представляют собой низкочастотные сглаживающие фильтры с автоматической подстройкой граничной частоты пропускания под основные частоты входного сигнала. Операторы фильтров можно получить обратным преобразованием Фурье:

h1:= ICFFT(H1)/(K+1) h2:= ICFFT(H2)/(K+!) Ü обратное преобразование Фурье

Рис. 12.3.3. Импульсные отклики фильтров.

Оператор фильтра, в принципе, бесконечен. В данном случае, при использовании БПФ максимальное число отсчетов равно К/2 = 500. Усечение размеров оператора может выполняться по типовым методам с применением весовых функций (в расчете операторы нормируются к 1, весовые функции не применяются).

N1:= 160 n1:= 0 .. N1 N2 ;= 500 n2:= 0 .. N2 Ü размеры и нумерация операторов

hm1:= h1 0 + 2·h1 n 1 hm1=0.988 h1:= h1/hm1 Ü нормировка

hm2:= h2 0 + 2·h2 n 2 hm2=1.001 h2:= h2/hm2 Ü нормировка

Ü свертка

Рис. 12.3.4. Проверка действия оптимальных фильтров.

Коэффициент усиления дисперсии помех Þ Kd:= (h 0) 2 + 2·h n Kd1=0.021 Kd2= 0.0066

Среднеквадратическое отклонение воспроизведения сигнала:

e1= 1.238 e2 = 0.701

Проверка действия оператора фильтра приведена на рис. 12.3.4.

Особую эффективность оптимальный фильтр имеет при очистке от шумов сигналов, имеющих достаточно сложный спектральный состав. Оптимальный фильтр учитывает конфигурацию спектра сигнала и обеспечивает максимальное подавление шумов, в том числе внутри интервала основного частотного диапазона сигнала. Это наглядно видно на рис. 12.3.5 для сигнала, близкого к прямоугольному, спектр которого имеет кроме основной низкочастотной части затухающие боковые осцилляции. Расчет фильтра выполнялся по методике, приведенной в примере 1.

Рис. 12.3.5. Оптимальная фильтрация сигнала со сложным спектральным составом.

Рис. 12.3.6. Оптимальная фильтрация радиоимпульса.

На рис. 12.3.6 приведен пример фильтрации оптимальным фильтром радиоимпульса. Основной пик спектра радиоимпульса находится в области несущей частоты, а боковые полосы определяются формой модулирующего сигнала, в данном случае – прямоугольного импульса. На графике модулей сигнала и передаточной функции фильтра можно видеть, что оптимальный фильтр превратился в полосовой фильтр, при этом учитывается форма боковых полос спектра сигнала.

Фильтры прогнозирования и запаздывания. В случае если в правой части уравнения (12.3.3) желаемым сигналом задать входной сигнал со сдвигом на величину kDt, то при этом B(m) = R(m+k) и уравнение принимает вид:

h(n) ③ R(m-n) = R(m+k). (12.3.8)

При k > 0 фильтр принято называть фильтром прогнозирования и вычисляет будущие значения сигнала по его предшествующим значениям. При k < 0 фильтр является фильтром запаздывания. Реализация фильтра заключается в решении соответствующих систем линейных уравнений для каждого заданного значения k. Фильтр может использоваться для интерполяции геофизических полей, в том числе в наперед заданные точки, а также для восстановления утраченных данных.

Вольсков Дмитрий Геннадьевич, кандидат технических наук, доцент кафедры «Самолётостроение» ИАТУ УлГТУ. Имеет монографию, научные статьи в журналах ВАК, методические пособия. Поступила 11.03.2016 г. УДК 621.391 К. К. ВАСИЛЬЕВ ДИСКРЕТНЫЙ ФИЛЬТР ВИНЕРА Рассмотрена классическая задача синтеза и анализа фильтра Колмогорова – Винера. Предложен способ синтеза оптимального реализуемого фильтра в дискретном времени. Приведены примеры решения задач синтеза и анализа оптимальных фильтров для случайных последовательностей с экспоненциальной корреляционной функцией. Ключевые слова: случайная последовательность, дисперсия ошибки, корреляционная функция, линейная оценка, фильтр Введение Во многих случаях полезная информация заключена в последовательности значений x1 , x2 ,..., xk изменяющегося в дискретном времени параметра. Для извлечения этой информации используются наблюдения z1 , z2 ,..., zk ,... , являющиеся функциями полезного параметра и помех. При этом на основе наблюдений необходимо дать наилучшую в определённом смысле оценку xˆk = xˆk (z1 , z2 ,..., zk ,...) значений изменяющегося параметра . К этому классу относятся, например, задачи оценивания изменяющихся параметров сигналов и помех в радиолокации, радионавигации и радиосвязи . Решение рассматриваемой задачи для дискретного времени впервые было дано А. Н. Колмогоровым . Построение оптимальных фильтров для непрерывного времени после фундаментального исследования Н. Винера представлено в большом числе работ. При синтезе непрерывного фильтра Винера проблема реализуемости решается достаточно просто на основе выбеливания наблюдений и факторизации передаточной функции. Вместе с тем для дискретного времени методика синтеза реализуемого фильтра, по-видимому, забыта в связи с появлением методов калмановской и квазилинейной фильтрации . В настоящей работе приводятся результаты синтеза и анализа оптимальных нереализуемого и реализуемого линейных фильтров в дискретном времени. Оценки = xˆk Оценивание в дискретном времени xˆk (z1 , z2 ,...,= zk ,...) xˆk (zi , i ∈ Dk) изменяющегося параметра функциями наблюдений { } { zi , i ∈ Dk } , сделанных на интервалах дискретного времени ных задачах для оценивания xk xk являются Dk . В различ- могут использоваться наблюдения, выполненные как до момента i = k , так и после этого момента. Если область Dk содержит моменты дискретного времени i ≤ k , то нахождение оценок вида xˆk = xˆk (z1 , z2 ,..., zk) называется фильтрацией. Оценивание xˆk на основе наблюдений z1 , z2 , ..., zk , zk +1 , ... , zk + m , m ≥ 1, полученных до и после момента i = k , называется интерполяцией, или сглаживанием. Если же i < k , то оценивание будущего значения xˆk xˆ= = xˆk (z1 , z2 , ... , zk −m), m ≥ 1, является экстраполяцией, или проk (z1 , z2 , ... , zi) гнозированием значения xk на основе предшествующих наблюдений. времени Васильев К. К., 2016 Вестник УлГТУ 1/2016 47 Способ получения данных об информационном параметре, как правило, известен и поэтому может быть дано и математическое описание связи наблюдений и значений параметра. Наиболее общим способом такого описания является совместная плотность распределения вероятностей w(z1 , z2 ,..., zk ,..., x1 , x2 ,..., xk ,...) = w(x1 , x2 ,..., xk ,...) w(z1 , z2 ,..., zk ,... / x1 , x2 ,..., xk ,...). При этом w(x1 , x2 ,..., xk ,...) определяет динамику изменения информационной случайной последовательности (СП) x1 , x2 ,..., xk ,... Для гауссовских СП полная информация о динамике СП x1 , x2 ,..., xk ,... может быть задана с помощью корреляционной функции (КФ) Bx (i, j) = M {(xi − mi)(x j − m j)} , i, j = 1, 2,... В рассматриваемых задачах описание процесса изменения параметра xi в дискретном времени xi , i = 1,2,..., k , сочетающий про- имеет принципиальное значение. Правильный выбор модели СП стоту математического представления и адекватность реальным физическим явлениям, зачастую представляет довольно сложную проблему и всегда требует тщательного анализа . Во многих реальных системах осуществляются наблюдения = zi h= (xi , ni), i 1,2,... , зависящие только от текущего значения параметра ностях ni , i = 1,2,... xi и случайной погрешности ni . При независимых погреш- с известными распределениями совместное распределение факторизуется k w(z1 , z2 ,..., zk / x1 , x2 ,..., xk) = ∏ w(zi / xi), и полная информация о способе получения данных i =1 w(zi / xi), i = 1,2,... h(xi , ni), i = 1,2,... . Заметим, что при изменяющихся случайных параметрах наиболее информативным о значении xk , как правило, оказывается наблюдение zk = h (xk , nk) в этот же момент времени. Вклад других наблюдений в результирующую оценку xˆk = xˆk (z1 , z2 ,..., zk ,...) должен зависеть от их числа и временных свойств СП x1 , x2 ,..., xk ,... модели наблюдения и характеристик помех ni , i = 1,2,... . содержится в или функциях Оптимальный линейный фильтр Пусть реализация СП x1 , x2 ,..., xk ,... представляет собой значения изменяющегося в дискретном времени информационного параметра. Будем предполагать, что до проведеи КФ ния наблюдений известны среднее значение = mi M= X (ti) , i 1,2,... { } Bx (i, j) = M {(X (ti) − mi)(X (t j) − m j)} , i, j = 1, 2,... такой СП. Для описания взаимодействия zi = h (xi , ni) скалярных параметров xi и гауссовских некоррелированных помех ni известными предполагаются КФ Bz (i, j) = M { zi z j } и взаимная КФ Bzx (i, k) = M { zi xk } . Ограничиваясь линейными оценками, можно учесть всё многообразие представленных априорных данных с помощью выбора весовых коэффициентов xˆk = {g } : ∑g i∈Dk i ,k i ,k zi , где Dk – область моментов времени выполненных наблюдений. Рассмотрим решение задачи построения оптимального в смысле минимума дисперсии ошибки = σ ε2k M {(xˆk − xk) 2 } линейного алгоритма оценивания xˆk = ∑ gi ,k zi изменяющегося параметра i∈Dk xk . Такой алгоритм нахождения оптимальных оценок c помощью весового суммирования наблюде- ний называют фильтром Колмогорова–Винера, или дискретным фильтром Винера . 48 Вестник УлГТУ 1/2016 После элементарных преобразований получим следующее выражение для дисперсии ошибки алгоритма с произвольными весовыми коэффициентами      2  2 M  ∑ g i ,k zi − xk   = σ ε2k = σ xk + M ∑   i∈Dk где    σ xk2 = M { xk2 } . i∈Dk Для поиска оптимальных весовых коэффициентов ошибки оценивания σ ε2k , продифференцируем σ ε2k {g } : i ,k ∑g j∈Dk i ,k  g j ,k zi z j − 2 xk ∑ g i ,k zi  , i∈Dk  { g , i ∈ D } , минимизирующих дисперсию по { g , i ∈ D } и приравняем производные i ,k k i ,k k нулю. Получим следующую систему уравнений:   M (∑ g j ,k z j − xk) zi  = 0, i ∈ Dk , j D ∈ k   или M {(xˆk − xk) zi } = 0, i ∈ Dk , показывающих, что ошибка оптимального оценивания ε= xˆk − xk и каждое из используемых k наблюдений { zi , i ∈ Dk } должны быть некоррелированны. Это условие называют принципом ортогональности ошибки и наблюдений, или леммой об ортогональном проектировании . После вычисления математического ожидания полученная система линейных уравнений B (i, j) ∑g = j∈Dk j ,k Bzx (i, k), i ∈ Dk z в качестве коэффициентов содержит значения КФ имной КФ Bz (i, j) = M { zi z j } (1) СП z1 , z2 ,..., zk ,... и вза- Bzx (i, k) = M { zi xk } СП z1 , z2 ,..., zk ,... и СП x1 , x2 ,..., xk ,... Систему уравнений (1) называют уравнениями Винера-Хопфа для дискретного времени. С учётом (1) легко находится и минимально достижимая в условиях рассматриваемой задачи дисперсия ошибки оценивания: 2 σ= σ xk2 − ∑ gi ,k Bzx (i, k) . εk i∈Dk При небольшом числе элементов в области индексов найти оптимальные весовые коэффициенты ходимо будет находить {g i ,k , i ∈ Dk } {g i ,k Dk можно решить систему уравнений (1) и, i ∈ Dk } . Однако для нестационарных СП для каждого необ- k -го шага оценивания, что может вызвать значи- тельные вычислительные проблемы даже для скалярной СП x1 , x2 ,..., xk ,... Решение значительно упрощается для стационарных СП, наблюдаемых на фоне аддитивных помех. В этом случае z= xk + nk , k Bzx (i, k)= M { zi xk }= M { xi xk }= Bx (i − k), Bz (i, j)= M { zi z j }= Bx (i − j) + 1, если i = j , + M {ni n j }= Bx (i − j) + σ 2δ K (i − j), δ K (i − j)=  0, если i ≠ j. Таким образом, система уравнений (1) преобразуется к виду: σ 2 gi + Вестник УлГТУ 1/2016 ∑g j∈D j Bx (i −= j) Bx (i), i ∈ D . (2) 49 Для бесконечного интервала наблюдений D: (−∞ < i < ∞) можно воспользоваться теоремой о свертке и тогда, после z -преобразования, получим H (z)(σ 2 + F (z)) = F (z) или где = H (z) F (z) / (σ 2 + F (z)) , ∞ F (z) = ∑ j =−∞ Bx (j) z − j , H (z) = димо воспользоваться обратным ∞ ∑gz j =−∞ −j j . Чтобы найти весовые коэффициенты, необхо- z -преобразованием: 1 gj = H (z) z j −1dz ,  ∫ 2π i C (3) где C − единичная окружность на плоскости комплексного переменного, а минимально достижимая дисперсия ошибки ∞ σ ε2k = σ x2 − ∑ g j Bx (j) = σ 2 g0 . j =−∞ Рассмотрим точное решение задачи нахождения коэффициентов дискретного фильтра Винера для экспоненциальной КФ F (z) = ∞ Bx (j) = σ x2 ρ j информационной СП. Тогда ∞ ∞ j=0 j =1 2 σ x2 ρ z − j σ x2 ∑ (ρ / z) j + σ x= = ∑ ∑ (ρ z) j + σ x2 j j = −∞ Весовая характеристика оптимального фильтра gj = где 1 H (z) z j −1dz =  ∫ 2π i C q 1 + q 2 + 2q (1 + ρ 2) / (1 − ρ 2) α= (q (1 − ρ) + (1 + ρ)) / 2 ρ , 2 2 1− ρ2 . (1 − ρ z)(1 − ρ z −1) (α − α 2 −1) j , легко находится, например, с помощью вычетов подынте- гральной функции. Анализ весовых коэффициентов g j , j = 0, ± 1, ± 2, ... , показывает, что синте- зированная структура оценивания является физически нереализуемой, так как для нахождения оптимальных оценок xˆk = ∞ ∑gz j =−∞ j j −k необходимо использовать не только предшествующие xˆk наблю- дения z1 , z2 , ..., zk , но и все последующие наблюдения zk +1 , zk + 2 , ... Можно реализовать приближение к оптимальной оценке с помощью запоминания части наблюдений и взвешивания в скользящем окне xˆk = N ∑gz j =− N j j −k , содержащем 2N + 1 элемент. Такая реализация приведёт к задержке получения оценок на N тактовых интервалах, а также к потерям в эффективности фильтрации за счёт потери информации, содержащихся в отброшенных наблюдениях. Для оценки величины таких потерь можно сравнить дисперсию ошибки фильтрации в скользящем окне с минимально достижимой величиной 2 g0 σ 2 = σ ε2 σ= q 1 + q 2 + 2q (1 + ρ 2) / (1 − ρ 2) . (4) Оптимальный реализуемый дискретный фильтр Для получения физически реализуемого алгоритма оценивания можно воспользоваться известным методом декорреляции наблюдений в непрерывном времени , который можно сформулировать следующим образом. Если бы на вход фильтра поступала последовательность некоррелированных наблюдений, то весовую характеристику можно было бы представить в виде двух компонент 50 Вестник УлГТУ 1/2016 g j , j = 0, 1, 2, ... и g j , j =−1, − 2, ... , каждая их которых позволяла бы получить две незави∞ симые оценки параметра. Оставляя только одну из них, можно найти структуру xˆk = ∑ g j z j −k оп- j =0 тимального реализуемого фильтра. Воплощение этой идеи возможно с помощью преобразования СП z= xk + nk , имеющих спектр σ 2 + F (z) = Ψ (z)Ψ∗* (z) , декоррелирующим фильтk ром с передаточной функцией H= 1 / Ψ (z) . Действительно, после прохождения СП наблюдений 0 z= xk + nk через декоррелирующий фильтр получим СП ηk с равномерным энергетическим k наблюдений спектром σ 2 + F (z) = (σ 2 + F (z)) H 0 (z) = Ψ (z)Ψ * (z) H 0 (z) H 0* (z) = 1 Bη (j) = δ K (j). Заметим также, что при этом не происходит потери информации о значении 2 и КФ параметра ηk xk , поскольку исходная СП z= xk + nk может быть в принципе восстановлена из СП k с помощью обратного преобразования. Таким образом, появляется возможность решения задачи оптимальной реализуемой фильтрации по схеме, представленной на рис. 1. Рис. 1. Оптимальное оценивание на основе декорреляции ∞ После декоррелирующего преобразования η j = ∑ g вi z j −i необходимо найти передаточную ха- i =0 рактеристику Hη (z) реализуемого линейного фильтра, преобразующего некоррелированную СП η j , j = 0, 1, 2, ... ∞ xˆk = ∑ gη j ηk − j . в оценку При этом коэффициенты j =0 gη j , j = 0,1,... опти- мального линейного преобразования находятся из уравнений Винера-Хопфа (2), записанных применительно к заданным условиям в виде ∞ ∑ gη i =0 i Bη (i= − j) Bη x (= j), j 0,1,... Поскольку для некоррелированной СП η j , j = 0, 1, 2, ... КФ Bη (i − j)= δ K (i − j), то систе- ма распадается на отдельные уравнения, являющиеся её решениями:  ∞  ∞ gη j = Bη x (j) = M { xkη k − j } = M  xk ∑ g вi zk − j −i  = ∑ g вi Bx (i + j) , j = 0,1,... =  i 0=  i 0 После этого с помощью фильтра (рис. 1) z -преобразования можно найти передаточную функцию линейного = Hη (z) а затем H (z) = и характеристику H 0 (z) Hη (z) . ∞ = gη j z − j ∑ =j 0 ∞ ∞ ∑ ∑g =j 0=i 0 оптимального вi Bx (i + j)z − j , реализуемого дискретного фильтра Винера: Покажем возможности синтеза реализуемого фильтра на примере рассмотренной задачи фильтрации стационарной СП z= xk + nk k j xk с КФ Bx (j) = σ x2 ρ . Для этого представим спектр наблюдений в виде двух комплексно-сопряжённых сомножителей: Вестник УлГТУ 1/2016 51 (1 − ρ 2) α (1 − β z −1) α (1 − β z) σ + F (z) = σ +σ = , (1 − ρ z −1)(1 − ρ z) (1 − ρ z −1) (1 − ρ z) 2 2 где коэффициенты β = ρσ / α и (1 + ρ 2) σ x2 2 2 2  2 2 2  , α 0,5σ (1 − ρ)  q + = + q + 1 + 2q (1 + ρ) / (1 − ρ)= , q 2 2 (1) ρ σ −   2 2 2 x находятся из условия тождественности полиномов в числителях левой и правой частей представленного разложения. Таким образом, передаточная характеристика выбеливающего фильтра имеет следующий вид: 1 (1 − ρ z −1) H0 = . = Ψ (z) α (1 − β z −1) Весовая функция такого фильтра может быть найдена с помощью интегрирования: 1 / α , если j = 0,  1 j −1 = H (z) z dz  0 j −1 ∫ 2π i  ((β − ρ) / α) β , если j ≥ 1. C = g вj Для рассматриваемого примера σ x2 (ρ − β) j gη j ∑ g вi Bx (= i + j) ∑ g вiσ= ρ ρ , j 0,1,..., = = αβ q =i 0=i 0 Hη (z) = σ x2 (ρ − β) / αβ q(1 − ρ z −1) ∞ ∞ 2 x i+ j , и поэтому (1 − β / ρ) H (z) = H 0 (z) Hη (z) = , gj = (1 − β / ρ) β j , j = 0, 1, 2, ... . −1 (1 − β z) Полученное выражение для передаточной функции позволяет представить процесс фильтрации в виде следующего рекуррентного соотношения: (1 − β z −1) xˆk = (1 − β / ρ) zk или β xˆk −1 + (1 − β / ρ) zk . xˆ= k Таким образом, появляется возможность значительно сократить объём вычислений за счёт замены ∞ весового суммирования xˆk = ∑ g j z j −k на эквивалентные с точки зрения достижения минимальной j =0 дисперсии ошибки рекуррентные вычисления с минимальным числом операций на каждом k-м шаге оценивания. Величина дисперсии ошибки σ ε2 = σ 2 g 0 = σ 2 (1 − β 2q)= σ2 ρ (1 + q) 1 + 1 + 4q ρ 2 / (1 − ρ 2)(1 + q) 2 () (5) для физически реализуемых алгоритмов фильтрации оказывается больше, чем (4). Например, при малых q и (1 − ρ) для нереализуемого фильтра σ ε2 / σ x2  1 / 1 + 2q / (1 − ρ) , а применение (5) σ ε2 / σ x2  2 / (1 + 1 + 2q / (1 − ρ)) . Сравнение этих формул показывает, что при медленном изменении информационной СП (1 − ρ) << q дисперсия ошиб- приводит к следующему выражению: ки (5) в 2 раза больше, чем для нереализуемого фильтра. Это объясняется использованием в реализуемом алгоритме вдвое меньшего числа наблюдений. 52 Вестник УлГТУ 1/2016 Рис. 2. Дисперсия ошибки реализуемого фильтра Заключение Основным результатом работы является завершенное изложение теории синтеза и анализа реализуемого и нереализуемого оптимального линейного фильтра в дискретном времени. Представленная методика и рассмотренные примеры будут полезны преподавателям, аспирантам и студентам, изучающим статистическую теорию оптимального приёма сигналов и теорию стохастического управления. СПИСОК ЛИТЕРАТУРЫ 1. Колмогоров А. Н. Интерполирование и экстраполирование стационарных случайных последовательностей // Изв. АН СССР. Сер. Математика. – 1941. − Т. 5, №1. − С. 3−14. 2. Wiener N. Extrapolation, Interpolation and Smoothing of Stationary Time Series. – N. Y. : MIT Press/John Wiley, 1964. − 171 p. 3. Тихонов В. И., Харисов В. Н. Статистический анализ и синтез радиотехнических устройств и систем: учебное пособие для вузов. – М. : Радио и связь, 2004. – 608 с. 4. Перов А. И. Статистическая теория радиотехнических систем: учебное пособие для вузов. – М. : Радиотехника, 2003. – 400 с. 5. Сейдж Э. П., Мелс Дж. Теория оценивания и её применение в связи и управлении / Пер. с англ.; под ред. Б. Р. Левина. – М. : Связь, 1976. – 495 с. 6. Васильев К. К., Служивый М. Н. Математическое моделирование систем связи: учебное пособие. – Ульяновск: УлГТУ, 2010. – 170 с. Васильев Константин Константинович, доктор технических наук, профессор, заведующий кафедрой «Телекоммуникации» УлГТУ. Поступила 14.03.2016 г. Вестник УлГТУ 1/2016 53

Линейную фильтрацию широко используют в системах передачи информации для обработки сигналов несмотря на то, что во многих случаях необходима нелинейная обработка. Объясняется это прежде всего простотой реализации линейных фильтров, которые сравнительно легко синтезируются, и существованием развитой теории их построения, чего нельзя сказать о нелинейных фильтрах.

Линейные фильтры являются неотъемлемой частью любого приёмного устройства. С их помощью осуществляется как додетекторная, так и последетекторная обработка сигналов. С помощью линейных фильтров сигналы часто разделяются в многоканальных системах передачи. Требования к этим

фильтрам могут быть весьма различными в зависимости от их назначения. Здесь рассмотрим теорию оптимальной линейной фильтрации.

Пусть сигнал на входе линейного фильтра с импульсной характеристикой представляет сумму переданного сигнала и помехи

Требуется найти такую функцию которая минимизирует среднеквадратическую ошибку

где оценка сигнала на выходе фильтра. Здесь считаем, что время запаздывания сигнала в фильтре а среднее значение берется по ансамблям сигналов и помех Будем полагать, что стационарные взаимно некоррелированные процессы с известными . В такой постановке задача была решена независимо друг от друга Колмогоровым 1939 г.) и Винером 1942 г.), и поэтому оптимальный (в указанном смысле) линейный фильтр называют фильтром Колмогорова-Винера (ФКВ). Требование физической реализуемости фильтра, как известно, сводится к тому, что импульсная характеристика фильтра должна удовлетворять условию для всех Это ограничение учитывается в записи

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

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

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

Докажем справедливость условия (8.60). Пусть импульсная характеристика оптимального фильтра, удовлетворяющего условию (8.60), импульсная характеристика любого другого линейного фильтра. Отклики фильтров соответственно обозначим через Тогда

И если справедливо (8.60), то

Следовательно,

Очевидно, что последнее выражение будет минимальным, когда что и доказывает справедливость условия (8.60) для оптимальной фильтрации. Геометрический смысл этого условия состоит в том, что случайный вектор должен быть строго ортогональной проекцией на линейное подпространство, порождаемое случайным вектором

Представим условие (8.60) в виде для всех из у. Отсюда с учётом или

В том случае, когда сигнал и помеха некоррелированы, (8.61) принимает вид

Это основное интегральное уравнение теории линейной фильтрации называется уравнением Винера-Хопфа. Его решением является искомая функция минимизирующая средний квадрат ошибки

Не следует пугать оптимальные линейные фильтры, определяемые (8.61) или (8.62), с согласованными фильтрами, рассмотренными в § 5.4. Если основное назначение рассматриваемых здесь фильтров состоит в наилучшем воспроизведении неизвестной формы сигнала, то задача согласованных фильтров заключается в формировании максимально возможного пика сигнала известной формы в момент отсчёта на фоне шума.

Уравнение (8.62) легко решается для нереализуемых фильтров, т.е. когда Для этого случая, применив преобразование Фурье к обеим частям (8.62), получим в частотной области

Отсюда коэффициент передачи оптимального линейного ФКВ

Докажем, что дисперсия ошибки при оптимальной нереализуемой линейной фильтрации в общей постановке

СПМ для случайного процесса

Записав можно видеть, что обеспечивается лишь при ФЧХ оптимального фильтра стало быть, при передаточной функции фильтра

Нетрудно видеть, что минимальная дисперсия ошибки

Легко заметить, что ошибка только в том случае, когда т.е. когда спектры сигнала и помехи не перекрываются. Во всех других случаях оптимальный фильтр пропускает различные частоты с тем меньшим весом, чем больше отношение при данной частоте.

Мы не будем здесь обсуждать вопросы реализации фильтра, характеристики которого приближаются к характеристике нереализуемого ФКВ, поскольку известен всегда реализуемый вариант фильтра, оптимальный по среднеквадратическому критерию, определяемый методом переменных состояния (см. ниже § 8.7).

Результаты оптимальной фильтрации можно существенно улучшить, если применить так называемое предыскажение сигнала с последующей его коррекцией на приеме. Сущность метода предыскажения состоит в том, что на передающей стороне сигнал пропускается через фильтр с коэффициентом передачи Полученный таким образом видоизменённый сигнал передаётся по каналу. На приёмной стороне включён другой фильтр Характеристики фильтров и выбираются так, чтобы обеспечить минимум среднеквадратической ошибки. Расчёты показывают, что предыскажение даёт тем больший выигрыш, чем меньше относительная ширина полосы перекрытия спектров сигнала и помехи. Предыскажение позволяет перераспределить мощность полезного сигнала в полосе частот канала так, чтобы обеспечить лучшие условия согласования источника сигнала с каналом (в общем случае полезно стремиться к тому, чтобы сумма спектральных плотностей мощности сигнала и мощности помехи была постоянной в пределах полосы частот канала). Это означает, что предыскажение можно рассматривать как некоторое "линейное кодирование" непрерывного сигнала, позволяющее уменьшить ошибку и улучшить использование пропускной способности канала.

Линейное предыскажение широко используется в современных системах связи. Характерными в этом отношении являются системы, в которых используется частотная модуляция. Согласно (8.41) плотность мощности шума на выходе ЧМ демодулятора увеличивается пропорционально квадрату частоты, так что верхние частотные составляющие сообщения подвержены шумам сильнее, чем нижние. Метод предыскажения и последующая коррекция позволяют снизить шум на верхних частотах и тем самым создать примерно одинаковые условия для передачи как нижних, так и верхних частот сообщения.

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

Условие оптимальности фильтра. Фильтр Колмогорова-Винера является оптимальным фильтром формирования из входного сигнала x(t) выходного сигнала z(t) при известной форме полезного сигнала s(t), который содержится во входном сигнале в сумме с шумами. В качестве критерия его оптимизации используется среднее квадратическое отклонение сигнала y(t) на выходе фильтра от заданной формы сигнала z(t). Подставим уравнение свертки (12.2.1) в раскрытой форме весового суммирования в выражение (12.2.2") и получим отклонение  2 выходного сигнала y(k) = h(n)③x(k-n) от заданной формы выходного сигнала z(k) по всем точкам массива данных:

В частном случае воспроизведения формы полезного сигнала в качестве функции z(k) принимается функция s(k). Минимум выражения (12.3.1) определяет значения коэффициентов h(n) оптимального фильтра. Для нахождения их значений продифференцируем выражение (12.3.1) по коэффициентам фильтра и приравняем полученные уравнения нулю. В итоге получаем:

где
- корреляционная функция входного сигнала,
- взаимная корреляционная функция сигналов z(k) и x(k). Отсюда:

h n R(m-n) = B(m), n = m = 0,1,2, ... , M. (12.3.2)

В краткой форме символической записи:

h(n) ③ R(m-n) = B(m). (12.3.3)

m=0: h o R(0)+ h 1 R(1)+ h 2 R(2)+ h 3 R(3)+ ...+ h M R(M) = B(0), (12.3.3")

m=1: h o R(1)+ h 1 R(0)+ h 2 R(1)+ h 3 R(2)+ ...+ h M R(M-1) = B(1),

m=2: h o R(2)+ h 1 R(1)+ h 2 R(0)+ h 3 R(1)+ ...+ h M R(M-2) = B(2),

..............................................................................................................

m=M: h o R(M)+ h 1 R(M-1)+ h 2 R(M-2)+ .... + h M R(0) = B(M).

Решение данной системы уравнений относительно значений h i дает искомый оператор фильтра.

h(n) ③ = B zs (m)+B zq (m). (12.3.4)

Частотная характеристика фильтра находится преобразованием Фурье левой и правой части уравнения (12.3.4):

H() = W zs ()+W zq (),

H() = / , (12.3.5)

где W s ()  R s (m) и W q ()  R q (m) - энергетические спектры (плотности мощности) сигнала и помех, W zs ()  B zs (m) - взаимный энергетический спектр входного и выходного сигналов, W zq ()  B zq (m) - взаимный энергетический спектр выходного сигнала и помех.

H() = W zs () / , (12.3.6)

Выражения (12.3.5-12.3.6), как правило, всегда имеют решение. Однако это не означает возможность формирования фильтром любой заданной формы выходного сигнала. Из чисто практических соображений можно сразу предполагать, что если спектр заданного сигнала z(t) больше значимой части спектра полезного сигнала s(t), то оператор фильтра попытается сформировать требуемые высокие частоты заданного сигнала из незначимых частот спектра полезного сигнала, что может потребовать огромных коэффициентов усиления на этих частотах, под действие которых попадут и частотные составляющие шумов. Результат такой операции непредсказуем. Эти практические соображения можно распространить и на все частотные соотношения входного и выходного сигналов линейных фильтров: значимые гармоники спектров выходных сигналов должны формироваться из значимых гармоник спектров входных сигналов.

Если заданная форма сигнала z(k) совпадает с формой полезного сигнала s(k), то B(m) = B ss = R s и фильтр называют фильтром воспроизведения полезного сигнала :

H() = W s () / , (12.3.7)

Выражения (12.3.6-12.3.7) достаточно наглядно демонстрируют физический смысл формирования передаточной функции фильтра. При воспроизведении сигнала частотная функция взаимной корреляции входного сигнала с выходным W zs (плотность взаимной мощности) повторяет частотную функцию автокорреляции W s (плотность мощности сигнала). Плотность мощности статистических шумов W q распределена по частотному диапазону равномерно, в отличие от плотности мощности сигнала W s , которая, в зависимости от формы сигнала, может занимать любые частотные интервалы спектрального диапазона. Частотная передаточная функция фильтра воспроизведения сигнала формируется отношением W s ()/. На частотах, где сосредоточена основная энергия сигнала, имеет место W s ()>>W q () и H()  1 (как минимум, больше 0.5). Там, где значение W s () становится меньше W q , коэффициент передачи фильтра становится меньше 0.5, и в пределе H()=0 на всех частотах, где полностью отсутствуют частотные составляющие сигнала.

Таким образом, оптимальные фильтры учитывают особенности спектрального состава сигналов и способны формировать передаточные функции любой сложности на выделение полезных частот сигналов из любых диапазонов спектра с максимальных подавлением шумов на всех частотах спектрального диапазона, не содержащих полезных сигналов, при этом границы усиления-подавления устанавливаются автоматически по заданному уровню шумов.

Задание мощности шумов. Следует внимательно относиться к заданию функции шумов Wq(). При полной неопределенности шума необходимо, как минимум, выполнять оценку его дисперсии  2 и распространять на весь частотный диапазон с соответствующей нормировкой на его величину (2Wq() d =  2), т.е. считать его белым шумом. При известной функции полезного сигнала в зарегистрированной реализации значение дисперсии шумов в первом приближении может быть оценено по разности дисперсий реализации и функции полезного сигнала. Можно выполнить и выделение шумов из входного сигнала в отдельный шумовой массив, например, вейвлетным преобразованием. Однако использовать выделенный шум непосредственно для вычисления функции Wq() допустимо только по достаточно представительной реализации при условии стационарности и эргодичности шума. В противном случае функция Wq() будет отображать только распределение шумов в зарегистрированной реализации сигнала, а соответственно фильтр будет оптимален только к этой реализации, что не гарантирует его оптимальности к любой другой реализации. Но для обработки единичной зарегистрированной реализации сигнала такой метод не только вполне допустим, но и может существенно повысить точность формирования выходного сигнала.

Эффективность фильтра. Из выражений (12.3.5-12.3.7) следует, что с позиции минимального искажения полезного сигнала при максимальном подавлении шумов фильтр Колмогорова-Винера эффективен в тем большей степени, чем больше отношение сигнал/шум на входе фильтра. В пределе, при W q ()<>W s () имеем H() 0 и сигнал будет сильно искажен, но никакой другой фильтр лучшего результата обеспечить не сможет.

Пример. Расчет оптимального фильтра воспроизведения сигнала. Выполняется в среде Mathcad .

Форма входного сигнала считается известной и близка к гауссовой. На входной сигнал наложен статистический шум с равномерным распределением мощности по всему частотному диапазону (белый шум), некоррелированный с сигналом, и функцию Wzq принимаем равной нулю. Для наглядного просмотра связи параметров фильтра с параметрами сигнала задаем модели сигналов в двух вариантах:

K:= 1000 k:= 0 .. K A:= 50

s1 k:= A·exp[-0.0005·(k-500) 2 ] s2 k:= A·exp[-0.00003·(k-500) 2 ]  информационные сигналы

Q:= 30 q k:= rnd(Q) – Q/2 x1 k:= s1 k + q k x2 k:= s2 k + q k  входные сигналы

Рис. 12.3.1. Модельные сигналы.

S1:= CFFT(s1) S2:= CFFT(s2) Q:= CFFT(q)  спектры сигналов

 спектры мощности

Ds1:= var(s1) Ds2:= var(s2) Dx1:= var(x1) Dx2:= var(x2) Dq:= var(q)  дисперсии

Ds1 = 124.308 Ds2 = 310.264 Dx1 = 202.865 Dx2 = 386.78 Dq = 79.038  информация

mean(Wq) = 0.079 Wq1:= (Dx1 – Ds1)/(K+1) Wq1 = 0.078  информация

Wq2:= (Dx2 – Ds2)/(K+1) Wq2 = 0.076  информация

Wq k:= Wq1  замена на постоянное распределение

Передаточные функции оптимальных фильтров (приведены на рис. 12.3.2):

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

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

h1:= ICFFT(H1)/(K+1) h2:= ICFFT(H2)/(K+!)  обратное преобразование Фурье

Рис. 12.3.3. Импульсные отклики фильтров.

N1:= 160 n1:= 0 .. N1 N2 ;= 500 n2:= 0 .. N2  размеры и нумерация операторов

hm1:= h1 0 + 2·h1 n 1 hm1=0.988 h1:= h1/hm1  нормировка

hm2:= h2 0 + 2·h2 n 2 hm2=1.001 h2:= h2/hm2  нормировка

 свертка

Рис. 12.3.4. Проверка действия оптимальных фильтров.

Коэффициент усиления дисперсии помех  Kd:= (h 0) 2 + 2·h n Kd1=0.021 Kd2= 0.0066

Среднеквадратическое отклонение воспроизведения сигнала:

e1= 1.238 e2 = 0.701

Проверка действия оператора фильтра приведена на рис. 12.3.4.

Рис. 12.3.5. Оптимальная фильтрация сигнала со сложным спектральным составом.

Рис. 12.3.6. Оптимальная фильтрация радиоимпульса.

Фильтры прогнозирования и запаздывания. Если в правой части уравнения (12.3.3) желаемым сигналом задать входной сигнал со сдвигом на величину kt, то при этом B(m) = R(m+k) и уравнение принимает вид:

h(n) ③ R(m-n) = R(m+k). (12.3.8)

При k > 0 фильтр называется фильтром прогнозирования и вычисляет будущие значения сигнала по его предшествующим значениям. При k < 0 фильтр является фильтром запаздывания. Реализация фильтра заключается в решении соответствующих систем линейных уравнений для каждого заданного значения k. Фильтр может использоваться для интерполяции геофизических полей, в том числе в наперед заданные точки, а также для восстановления утраченных данных.

Курсовая работа 18-07. Разработка программы расчета оптимального фильтра приема данных в коде Манчестер-II и фильтрации цифровых данных.

Курсовая работа 19-07. Разработка программы расчета фильтра синхронизации приема данных в коде Манчестер-II и фильтрации цифровых данных.