Можно рассчитать оптимальные фильтры деконволюции, метрика приближения которых меньше, чем у усеченных фильтров деконволюции. Для получения общего уравнения оптимальной деконволюции будем считать, что число коэффициентов оператора hn равно M+1, a число коэффициентов инверсного оператора hn-1 равно N+1.
Принцип оптимизации. Выходная функция приближения при использовании уравнения свертки (13.1.2) с ограничением числа членов оператора фильтра:
F = Е2 = [do(k)-sk]2. sk = hn-1 hk-n. (13.3.1)
Чтобы определить минимум функции, приравняем нулю частные производные от Е по неизвестным коэффициентам фильтра:
dF/dhj-1 = -2hk-j [do(k) - hn-1 hk-n] = 0. (13.3.2)
hk-j hn-1 hk-n = hk-j do(k) = h-j. (13.3.3)
hn-1 hk-n hk-j = hn-1 aj-n = h-j, j = 0,1,2,..., N, (13.3.4)
где aj-n - функция автоковариации импульсной реакции h(n). Учитывая также, что hn = 0 при n<0 и aj = a-j (функция автоковариации является четной функцией), окончательное решение определяется следующей системой линейных уравнений:
a0 h0-1 + a1 h1-1 + a2 h2-1 + a3 h3-1 +... + aN hN-1 = h0 (13.3.5)
a1 h0-1 + a0 h1-1 + a1 h2-1 + a2 h3-1 +... + aN-1 hN-1 = 0
a2 h0-1 + a1 h1-1 + a0 h2-1 + a1 h3-1 +... + aN-2 hN-1 = 0
|
|
...........................................
aN h0-1 + aN-1 h1-1 + aN-2 h2-1 + aN-3 h3-1 +... +a0 hN-1 = 0
Пример расчета оптимального фильтра деконволюции.
Повторим инверсию оператора, приведенного на рис. 13.2.1, N=6.
Рис. 13.3.1. |
Значения оптимального инверсного оператора в сопоставлении с усеченным:
h-1(n) = {4.56, -6.033, 2.632, 0.417, -0.698, -0.062, 0.267} – прямой расчет по (13.2.2).
h-1(n) = {4.557, -6.026, 2.633, 0.397, -0.693, -0.009, 0.145} – расчет по (13.3.5).
Значения свертки инверсных операторов с прямыми и метрики приближения:
Оператор по (13.2.2) – рис. 13.3.1(А): sn= {1, 0, 0, 0, 0, 0, 0, 0.005, 0.031, 0.027, 0.013, 0.004, 0.001, 0, 0,…}. E=0.044.
Оператор по (13.3.5) – рис. 13.3.1(В): sn= {0.999, <0.001, 0.002, -0.003, -0.003, 0.013, -0.008, -0.012, 0.011, 0.013, 0.007, 0.002, <0.001, 0, 0, …}. E=0.027.
Метрика приближения оптимального оператора в 1.6 раза меньше усеченного.
Как видно на рис. 13.3.1, оптимизация инверсного оператора заключается в центрировании ошибок приближения и с распределением по интервалу суммарной длины прямого и инверсного оператора.
Уравнение оптимальной инверсии. Оптимальный инверсный фильтр может быть получен непосредственно с использованием z-образов импульсной реакции и автоковариационной функции прямого фильтра. Если для прямого фильтра мы имеем передаточную функцию H(z), то z-образ автоковариационной функции фильтра (как z-отображение спектральной плотности мощности) представляет собой произведение:
A(z) = H(z)H*(z), (13.3.6)
где H*(z)- функция, комплексно сопряженная с H(z). Переносим H(z) в левую часть формулы и для функций диракоидного типа заменяем выражение 1/H(z) = H-1(z):
А(z)H-1(z) = H*(z). (13.3.7)
Запишем последнее равенство в развернутом виде:
(a-Nz-N+... +a-1z-1+a0+a1z1+... +aNzN)(h0-1+h1-1z1+h2-1z2+... +hN-1zN) =
= h0*+h1*z-1+h2*z-2+... +hN*zN. (13.3.8)
|
|
В выражении (13.3.8) сумма коэффициентов при одинаковых степенях z в левой части равенства должна быть равна коэффициенту при соответствующей степени z в правой части равенства, что позволяет составить следующую систему из N уравнений для коэффициентов при степенях z0, z1, z2,..., zN:
a0 h0-1 + a-1 h1-1 + a-2 h2-1 + a-3 h3-1 +... + a-N hN-1 = h0* (13.3.9)
a1 h0-1 + a0 h1-1 + a-1 h2-1 + a-2 h3-1 +... + a-N-1 hN-1 = 0
a2 h0-1 + a1 h1-1 + a0 h2-1 + a5 h3-1 +... + a-N-2 hN-1 = 0
...........................................
aN h0-1 + aN-1 h1-1 + aN-2 h2-1 + aN-3 h3-1 +... +a0 h-N -1 = 0
В случае вещественных фильтров, когда ai = a-i и h0* = h0, уравнение (13.3.9) идентично уравнению (13.3.5).
Уравнение Левинсона. Практический способ расчета оптимальных инверсных фильтров по уравнению (13.3.9) предложен в 1947 году Н.Левинсоном.
Перепишем уравнение (13.3.9) в матричной форме:
(13.3.10)
Так как коэффициенты инверсного фильтра достаточно определить с точностью до произвольного масштабного множителя, приведем ho-1 к 1, a функцию автоковариации переведем в функцию коэффициентов корреляции делением обеих частей уравнения на ао. Обозначая Ai = ai/ao, Wi = hi-1/ho-1 и V = ho*/(ho-1ao) = hoho*/ao, получаем:
(13.3.11)
где для значений W и V введен индекс j номеров предстоящих итераций по циклу вычисления коэффициентов фильтра.
При нулевой итерации (N=0, j=0) имеем только одно уравнение:
. (13.3.12)
Благодаря проведенной нормировке решения уравнения (13.3.12) не требуется:
А0= 1, V0= 1, W00= 1.
Составим уравнение для двучленного фильтра (N=1, j=1):
(13.3.13)
Перепишем уравнение (13.3.12) в прямой форме:
А0 W00 = V0. (13.3.14)
Запишем вспомогательную систему, для чего к уравнению (13.3.14) добавим вторую строку с новой постоянной Ej:
A0 W00 + A1·0 = V0,
A1 W00 + A0 ·0 = E1.
В матричной форме:
(13.3.15)
Реверсируем уравнение (13.3.15):
(13.3.16)
Вычтем (13.3.16) из (13.3.15) с неопределенным множителем Rj:
(13.3.17)
Из верхней строки уравнения (13.3.16) можно получить значение Е1:
Е1= A1W00. (13.3.18)
Уравнение (13.3.13) можно сделать равнозначным уравнению (13.3.17), если правую часть нижней строки уравнения (13.3.17) приравнять к правой части нижней строки уравнения (13.3.13):
E1 - R1V0 = 0, R1 = E1/V0. (13.3.19)
При этом из правых частей верхних строк уравнений (13.3.13, 13.3.17) будем иметь:
V1 = V0 - R1E1. (13.3.20)
Приравнивая друг другу левые части уравнений (13.3.13,13.3.17), получаем:
W01 = W00 - R1·0 = W00 = 1.
W11 = 0 - R1W00 = -R1W00. (13.3.21)
Этим заканчивается первая итерация. Аналогично, для второй итерации:
(13.3.22)
(13.3.23)
(13.3.24)
(13.3.25)
Из верхней строки уравнения (13.3.24):
Е2 = A1W11+A2W01.
Из правых частей нижней и верхней строк уравнений (13.3.22, 13.3.25):
R2 = E2/V1,
V2 = V1 - R2E2.
Новые коэффициенты из левых частей уравнений (13.3.22, 13.3.25):
W02 = W01 - R2 0= 1,
W12 = W11 - R2W11,
W22 = 0 - R2W01.
Анализ расчетов позволяет вывести следующие рекуррентные формулы:
Ej = AiWj-i,j , j = 1,2,...,M. (13.3.26)
Rj = Ej/Vj-1,
Vj = Vj-1 - RjEj,
Wi,j = Wi,j-1 - RjWj-1,j-1, i = 0,1,.., j.
Подпрограммы решения уравнений для ЭВМ приведены в литературе /12,22/.
Курсовая работа 21-07. Разработка программы расчета оптимального фильтра деконволюции сигналов и фильтрации цифровых данных.