Тема 12. Оптимальные линейные цифровые фильтры (стр. 3 из 5)
Эффективность фильтра. Из выражений (12.3.5-12.3.7) следует, что с позиции минимального искажения полезного сигнала при максимальном подавлении шумов фильтр Колмогорова-Винера эффективен тем в большей степени, чем больше отношение сигнал/шум на входе фильтра. В пределе, при Wq(w)<<Ws(w) имеем H(w)Þ 1 и фильтр воспроизводит входной (или заданный) сигнал без искажений (помех нет, исправлять нечего). Отметим также, что помеха, коррелированная с полезным сигналом, как это следует из (12.3.5), используется фильтром для повышения точности воспроизведения сигнала. С другой стороны, при Wq(w)>>Ws(w) имеем H(w) Þ 0 и сигнал будет сильно искажен, но никакой другой фильтр не обеспечит лучшего результата.
Пример. Расчет оптимального фильтра воспроизведения сигнала. Выполняется в среде Mathcad.
Форма входного сигнала считается известной и близка к гауссовой. На входной сигнал наложен статистический шум с равномерным распределением мощности по всему частотному диапазону (белый шум), некоррелированный с сигналом, и функцию Wzq принимаем равной нулю. Для наглядного просмотра связи параметров фильтра с параметрами сигнала задаем модели сигналов в двух вариантах:
K := 1000 k := 0 .. K A := 50
s1k := A·exp[-0.0005·(k-500)2] s2k := A·exp[-0.00003·(k-500)2] Ü информационные сигналы
Q := 30 qk := rnd(Q) – Q/2 x1k := s1k + qk x2k := s2k + qk Ü входные сигналы
Рис. 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 Ü информация
Wqk := 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 := h10 + 2·
h1n1 hm1=0.988 h1 := h1/hm1 Ü нормировкаhm2 := h20 + 2·
h2n2 hm2=1.001 h2 := h2/hm2 Ü нормировка 
Ü свертка 
Рис. 12.3.4. Проверка действия оптимальных фильтров.
Коэффициент усиления дисперсии помех Þ Kd := (h0)2 + 2·
hn 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. Фильтр может использоваться для интерполяции геофизических полей, в том числе в наперед заданные точки, а также для восстановления утраченных данных.
12.4. Оптимальные фильтры сжатия сигналов.
Условие оптимальности. Фильтр сжатия сигнала x(t), по существу, представляет собой фильтр формирования сигнала z(t) с эффективной шириной длительности, меньшей по сравнению с эффективной шириной длительности полезного сигнала s(t) во входном сигнале x(t). Расчет оптимального фильтра сжатия может выполняться непосредственно по выражениям (12.3.3).
В предельном случае сжатия сигнала до импульса Кронекера положим, что z(k)=d(k) при статистической независимости сигнала и шума. Отсюда:
Bsz(m) = d(m) ③ s(k+m) = s(-m).
h(n) ③ (Rs(m-n)+Rq(m-n)) = s(-m). (12.4.1)
H(w) = S*(w) / (|S(w)|2+Wq(w)). (12.4.2)
При некоррелированной помехе с дисперсией s2 система уравнений для определения значений коэффициентов h(n):
ho(R(0)+s2)+ h1R(1)+ h2R(2)+ h3R(3)+ ...+ hMR(M) = s(0), (12.4.3)
hoR(1) + h1R(0)+ h2R(1)+ h3R(2)+ ...+ hMR(M-1) = 0,
hoR(2) + h1R(1)+ h2R(0)+ h3R(1)+ ...+ hMR(M-2) = 0,
. . . . . . . . . . . . .
hoR(M) + h1R(M-1)+ h2R(M-2)+ ....... + hMR(0) = 0.
При расчете коэффициентов фильтра значение s(0) обычно принимается произвольным, чаще всего равным площади сигнала s(t). Тем самым делается попытка полного сжатия площади сигнала до весовой функции Кронекера, что возможно только для сигналов со спектром в главном диапазоне до частоты Найквиста.
Рис. 12.4.1. Сжатие гладких сигналов с разным уровнем шумов.
Для гладких и монотонных функций со спектром в низкочастотной части главного диапазона сжатие до импульса Кронекера невозможно, и в зависимости от уровня шумов фильтр поднимает насколько возможно высокие частоты сигнала, учитывая значение уровня шумов. При этом нарушаются условия нормировки площади оператора фильтра к 1, о чем можно судить по значению передаточной функции H(w) при w=0, которое становится меньше 1, и при обратном преобразовании H(w) Þ h(m) оператор h(m) требует нормировки к 1. Все эти факторы можно наглядно видеть на рис. 12.4.1.
На рисунке приведены три сигнала с одной и той же базовой функцией, на которую наложены шумы разного уровня. При малом уровне шумов (сигнал х1) фильтр в максимальной степени использует высокие частоты сигнала (|H1| >>1 на этих частотах), сохраняя устойчивость работы фильтра при достаточно удовлетворительном (хотя и больше 1) коэффициенте усиления дисперсии помех при максимально возможном сжатии сигнала. При повышении уровня шумов (сигналы х2 и х3) подъем высоких частот сигнала уменьшается, и сжатие сигнала соответственно также уменьшается, предпочтение отдается максимальному подавлению шумов.
