CC BY
33
УДК 621.372.542
АЛГОРИТМЫ СПЛАЙН - АППРОКСИМАЦИИ ДВУМЕРНЫХ АНАЛИТИЧЕСКИХ ДИСКРЕТИЗИРОВАННЫХ СИГНАЛОВ
© 2004 П.К. Ланге
Самарский государственный технический университет
Рассмотрена задача построения цифровых фильтров с параболической сплайн - аппроксимацией двумерных дискретизированных сигналов. Определены погрешности аппроксимации такими фильтрами двумерного Гауссового пика.
В настоящее время в практике аналитических измерений все более широко применяются комбинированные аналитические приборы (газовый хроматограф + масс-спектрометр, газовый хроматограф + ИК-Фурье спектрометр, жидкостный хроматограф + УФ спектрометр и т.д.). В состав таких приборов, кроме хроматографа, входит и спектро-анализатор, определяющий спектр в дискретных точках хроматограммы.
Аналогичным образом организован и аналитический прибор, состоящий из дифференциального термического анализатора и спектроанализатора.
Выходной сигнал такого комбинированного прибора является двумерным, определенным относительно двух независимых переменных - например, времени и оптической длины волны. Спектрограмма, полученная с помощью комбинированного анализатора, является функцией двух переменных, определена в трехмерном пространстве и представляет собой совокупность двумерных аналитических пиков (рис.1).
В случае, когда комбинированный анализатор представляет собой хроматограф + УФ спектрометр, сечение такой спектрограммы рядом плоскостей, параллельных плоскости xOz, т.е плоскости рисунка, представляет собой ряд хроматограмм, снятых при разных длинах оптического спектра.
Аналогично этому сечение рядом плоскостей, перпендикулярных плоскости рисунка, представляет собой ряд оптических спектров, снятых в ряде точек хроматограммы, то есть в разные моменты времени.
При преобразовании и обработке как
Рис. 1. Двумерная спектрограмма: ось х - время, ось у - длина волны спектра, z ■ значение аналитического сигнала
одномерных, так и двумерных аналитических сигналов с целью их сжатия и определения информационных параметров актуальными являются задачи аппроксимации дискретных аналитических данных.
Рассмотрим аппроксимацию двумерных аналитических данных двумерными сплайнами на основе методов цифровой фильтрации.
В общем случае двумерный параболический сплайн на участке дискретизации (х* У* хп+1, Уп+1) определяется двумерной параболой
/ (х, у) = аих2 + 2а14х + а22у2 + 2а24у + 2а12ху + а44
(1)
Задача определения коэффициентов а в выражении (1) существенно упрощается, если в качестве аппроксимирующей функции применить выражение, все члены которого зависят только от какой - либо одной переменной.
В этом смысле рассмотрим возможность том последнего соотношения получим выра-использования "усеченной" двумерной пара- жения болы (1), положив а12 = 0. Такая парабола описывается выражением
[п1 +1] = а2 [п1 ] + а1 [п1 ] + а0 [п1 ]; а [ п1+1] = 2а2 [ п1 ] + а1 [ п1 ];
Л (х= У) = (а2х' + а1Х + а0 ) + (ь2у2 + Ь1У + ь0 ) = Ла (х) + Л (у)
(2)
Первая и вторая производные параболической функции (2) определяются выражениями
Л (х У)
йх
= 2а2 [п1 ] х + а1 [п1 ];
й У (х, у)
йх2
йЛ (х У)
= 2а2 [ п ];
йу
= 2Ь2 [ п2 ] У + Ь1 [ п2 ];
(3)
Ь0 [п2 + 1] = Ь2 [п2 ] + Ь1 [п2 ] + Ь0 [п2 ] ; (5)
Ь1 [П2 + 1] = 2Ь2 [П2 ] + Ь1 [П2 ] ;
Применив к системе (5) z - преобразование [1], получаем систему уравнений
а2 []+а1 []+ (1 -^К [] =0 ;
2а2 [ 21 ] + (1 - ¿>1 [ 71 ] = 0 ;
Ь2 [¿2] + ¿1 [¿2]+ (1 - 12)Ь0 [¿2] = 0 ; (6) 2Ь2 [¿2] + (1 - 12)ЬХ [¿2] = 0 ;
й У (х, У )
йу2
= 2Ь2 [ п2 ];
Рассмотрим возможность определения коэффициентов аппроксимирующих функций с использованием методов цифровой фильтрации.
Пусть передаточная функция цифрового фильтра, формирующего младшие коэффициенты а0[п;] , Ь0[п2] двумерной парабо-
ют вид
Л (х, У)
йх
( х, У )
= а1 [ п1 ];
. й2 Л (х, У)
х=0
йх
х=0
йу
= ¿1 [ п2 ];
. й2 Л (х, У)
У=0
йу2
= 2а2 [п,]; = 2Ь2 [п2 ];
У=0
Здесь п1 определяет дискретизацию по параметру х, а п2 - по параметру у.
При х = 0 , У = 0 (на границах интервалов дискретизации) эти выражения принима- лы (2) на каждом участке дискретизации
определяется выражением
F[ъ, п2] = Fa [п1 ] + Fь [п] , (7) состоящим из суммы функций
т
Ра [п1 ] = а0 [п1 ]= Е 4х [п1 + i] ,
i=-т т
^ [п2 ] = Ь0 [п2 ] = Е 1 [п2 + 1 ] ,
] =-т
Эти функции определены на нечетном числе (2т+1) дискретных значений сигнала.
Такой фильтр является разделимым [2], что в значительной степени упрощает его проектирование и позволяет для анализа фильтра привлечь методы проектирования одномерных цифровых фильтров.
Для оценки частотных свойств фильтра необходимо выполнить z - преобразование выражения (7).
Используя дискретное z - преобразование (7), получаем:
(4)
Пусть дискретные интервалы измеряются в относительных единицах:
хв = х [п1 +1] - х [п1 ] = ^ Ув = У [п2 + 1] - У [п2 ] = 1.
Тогда выполняются равенства
а01 п
[п1 +1] = /(х У)|х=1; а1 [п1 +1] = Ь0 [п2 +1] = 1(x, у) у=1; Ь1 [п2 +1] =
(x, у )
йх йЛ (х У)
йу
У=1
В связи с тем, что на границах дискретных интервалов параболическая сплайн -функция не должна иметь разрывов по нулевой и первой производным, из (3) и (4) с уче-
а.
[ ^ ] = Fa [ ^ ] , Ь [ ] = Fь [ ]
Подставляя это уравнение в систему (6),
определяем систему
a
[ Z ] + а1 [ Z ] = (Z - 1) Fa [ Z ] ;
Ь2 [] + Ь1 [] = (- l)Fb [Z2 ] ;
2a2 [ z ]+(i- z1)a1 [z ] = 0 ;
2b2 [Z2] + (1 - z2)\ [Z2] = 0;
Коэффициенты a , b определяются решением этой системы:
a [ z ] = 2 Fa [ z1 ] ; b, = 2 Fb [ z2 ] z, +1 z +1
(Z -1)
(z2 -1)2
1 Т Т 2
Л* = Т1 ПКХ^пюу)"7(х,У)] У,
Т 0 0
(11)
где Т - период сигнала.
Используя замену переменных интегрирования в (11) с учетом соотношения
2п
а =
T '
окончательно получаем: 1
1
a2 hЬ1^ F [Z] ; b2 [z2 Fb [z2].
(8)
Необходимо выбирать такую функцию F[z], чтобы она без остатка делилась на двучлен (z+1), иначе в противном случае появится фазовая погрешность фильтра [2].
Оценим качество аппроксимации гармонического сигнала (двумерной синусоиды) аппроксимирующим сплайн-фильтром [2]
f ( х, y ) = ( sin сох )•( sin coy ) . (9)
В связи с тем, что аппроксимирующие фильтры позволяют определить значение сигнала для произвольного момента времени внутри интервалов дискретизации, выходной сигнал такого фильтра можно считать квазинепрерывным, и для определения его эффективности использовать методы оценки погрешности аппроксимации непрерывных функций.
При малой погрешности аппроксимации гармонического сигнала, что обычно имеет место на практике, выходной сигнал такого фильтра можно также считать близким к синусоидальному:
z(х, y) « AX (с) sin [ах - cpX (с)] х
X Ar (c)sin [cy - <рг (y)] , (10)
где А, р - значения амплитудно - частотной характеристики (АЧХ) и фазо - частотной характеристики (ФЧХ) фильтра на частоте а .
Тогда значение средне - квадратичной погрешности аппроксимации определяется выражением
Аа(ю) =4-АХИ4(ю)()(со*Фт)+4 ;
(12)
Отсюда видно, что минимальное значение погрешности А а достигается при минимальных значениях ФЧХ (р(ю)« 0), а
также при минимальном отличии АЧХ от единицы. При нулевых фазовых погрешностях выражение для средне-квадратичной погрешности аппроксимации принимает вид
Л (с) Р-AX (с)A (с)]
A a (с) =-4-
(13)
В связи с тем, что выходной сигнал цифрового фильтра представляет собой в данном случае дискретные значения гармонического сигнала (9), то выражение (13) может быть применено и для оценки средне-квадратичной погрешности аппроксимации гармонического сигнала
/[п1,п2] = (Бтет*)•(Бтетп2) . (14)
где та =
= 2я/
N
относительная угловая
частота, N - число дискретных участков на периоде функции /[п1,п2].
Будем искать такую функцию цифрового фильтра F[z] , чтобы его фазовая погрешность была бы равна нулю:
р (ет) = 0 ,
что выполняется при симметрии их весовых функций (рис. 2).
Такой фильтр описывается выражениями
2
2
ЬК
(п-1) п (п+1) (п+2) k
Рис. 2. Симметричная весовая функция цифрового фильтра
а0 [п1] = а0х[п1] + а1х[п +1] +... + а1х[п1 -1] + А2х[п1 -2] +... =
= а0х [«1 ] + ЕА( х [п + i] + х [п1 - i])
Ьс к ] = в„ у [п2 ] + в у [п +1] +... + У [п2 -1] + В2 у [п - 2] +... =
= ВУ [п2 ] + ЕВ, (У [п2 + i] + У [п2 - i])
(15)
Z - преобразование этого выражения имеет вид
а0 [ ¿1 ] = х [ ¿1 ] Ь0 [ ¿2 ] = У [ ¿2 ]
т
А + Е А, (4+¿г)
i=1 _
т
В0 + Е В, (+г-)
i=1
(16)
К [
[ * ] = ах^= А Ф (* + ^ );
[ ] = МФ В0 + £в, ()
У1121 ,=1
= еШ1 ; ¿2 = еш2
Выражение для частотной характеристики такого фильтра имеет вид
,=1
ка (М ) = А0 +Е А,
т
= А0 + 2 Е а,
1 ш, г - 1 ш, г е^ 1 + е 1
(18)
г=1
г=l
К (л 2 ) = В0 +Е в,
т
= В0 + 2 Е в,
е1 Ш + е -1Ш2/
(19)
/=1
Следовательно, при выбранной весовой функции (рис. 2) АЧХ фильтра определяются выражениями (18) и (19).
С целью определения коэффициентов двумерной сплайн - аппроксимации разложим АЧХ (18) и (19) в степенной ряд в окрестности точек ш1 = 0 ; ш2 = 0.
Эти разложения имеют вид
Ка (л) к (0)+кА (0) л + кА1 (0)ш
+. .
Дискретная частотная характеристика такого фильтра определяется выражением
Кв (п)*Кв (0)+кВ (0+КВ11 (0+. . .
(20)
где К ( 0), К1 ( 0), К11 ( 0) - значения функций К (ш) и ее производных в точке
та« 0.
Следовательно, минимальное значение средне-квадратичной погрешности Дя , определяемой выражением (13) для выбранной (17) симметричной функции (15) фильтра, будет достигаться при выполнении соотношений
КА (0) = 1; КВ (0) = 1; КА1 (0) = 0; КВ1 (0) = 0; (21)
К а11 (0 ) = 0; Квп (0 ) = 0;
Выражения (21) позволяют определить систему уравнений относительно искомых к°т°р°е для разделимого фильтра делится на коэффициентов аппроксимирующего цифро-
Качество аппроксимации в дискретных точках можно оценить по частотной и фазовой характеристикам цифрового фильтра, получаемых из функций (17) при подстановке в него соотношения [2]
2 = е1Ш1 . е1Ш2,
соотношения
вого фильтра.
В частности, из (18) и (19) с учетом (21) получаем:
т
Ка ( 0) = А + 2Е А = 1;
/=1
к4 (0)=-2Е 4 81пш1г ш = 02Е 4 ш1г=
,=1 1
т
к А1 (0) = -2Е'2 А С03Ш1гШ
/=1 т1 тт
«-2Е А +Е А,4 Ш12 = 0;
,=1 ,=1
(22)
т
Кв (0) = В0 + 2ЕВг = 1;
,=0
т
КВ (0) = -2 Е/В/ в1пш2
г=1
т
К11 (0) = -2Е/2 В
/=1
тт
«-2Е/2 В +Е В/ 4 Ш22 = 0;
/=1 /=1
«-2Е/В/ Ш2 / = 0;
=0 /'=1
/ С0Б ш2 /
Ш2 = 0
(23)
На основании (22) и (23) составим систему уравнений относительно искомых коэффициентов А,, В :
4 + 2Е А■ = 1; В0 + 2Е В/ = 1;
/=1 г=1
т т
Е г2 А, = 0; Е г2 В, = 0;
,=1
г=1
Е г4 А, = 0; Е г4 В, = 0;
г'=1
г=1
Z - преобразование этих функций имеет вид
а, [ ¿1 ] = (А г12+аа+40+аа-1+а2 гг2) х [ ¿1 ];
ь0 [ 22 ] = (в2222 + в122 + в0 + в122-1 + в222^ ) у [ 22 ] ; (2б)
Из (24) для фильтра с нечетным числом членов получаем:
А0 + 2 А1 + 2 А2 = 1 ;
В0 + 2 В1 + 2 В2 = 1 ;
А1 + 4 А2 = 0 ; (27)
В1 + 4 В2 = 0 ;
Третье соотношение для определения коэффициентов А , В определяется из условия делимости выражений (8) на двучлены ^+1), (^+1).
Это условие позволяет определить два недостающих уравнения:
4 -2А1 + 2А2 = 0; В0 -2В1 + 2В2 = 0.
(28)
(24)
Решение системы (24) определяет коэффициенты двумерной сплайн - аппроксимации дискретных значений сигнала.
Двумерный пятиточечный параболический сплайн - фильтр, например, определяется выражениями
а [п1 ] = а2х[п1- 2]+А1х[п1 -1]+4)х[п1 ]+А1х[п1+1]+а2х[п1 + 2]; (25)
Ь [п ] = В2У [и2 - 2] + В1У [п2 -1] + В0У [и2 ] + В1У [п2 +1] + ВгУ [и2 + 2] Л '
Решение системы (27) и (28) позволяет определить коэффициенты А. , В. :
4 = В0 = 1C16 ; А1 = В1 = /1б ; А2 = В2 = - /1б
Таким образом, функции, определяющие младшие коэффициенты пятиточечного параболического двумерного сплайн - фильтра, имеют вид
а0 [п] = /1б(_х[И1 _2] + 4х[п1 -1] + 10х[п1] + 4х[^ + 1]-х[п1 + 2]) ; Ь0 ["2 ] = Ц6 (-У ["2 - 2] + 4У [п -1] +10У [п ] + 4У [пг +1] - у [п + 2]).
(29)
Выражения для коэффициентов а1[п1] , а2[п;] , ЬДп2] и Ь2[п2] определяются из (8):
а1 [п1 ] = 18 (х [п - 2]- 6х [п1 -1] + 6х [п1 +1]- х [п1 + 2]) ;
ь1 [п2 ] = /8 (у [п2 - 2]-6 у [п2 -1]+6 у [п2 +1]-у [п2 + 2]);
а2 [п1 ] = Х6 (-х [п1 - 2] + 7х [п1 -1] - 6х1 [п ]
-6 х [п1 +1]+ / х [п1 + 2] - х [п1 + 3]);
Ь2 [п2 ] = /{6 (- У [ п2 - 2] + 7 У [п2 - 1]-6 У1 [п2 ]-- 6 У [ п2 +1] + 7 у [п2 + 2]-у [п + 3]).
(30)
На практике свойства цифровых фильтров удобно иллюстрировать их амплитудно -
8 %
Частотная характеристика разработанного параболического пятиточечного сплайн -фильтра в зависимости от относительных частот гармонического сигнала
Рис. 3. Зависимость максимального значения
погрешности аппроксимации двумерного гармонического сигнала от числа интервалов дискретизации на периоде сигнала
частотными характеристиками, переходными функциями, а также аппроксимациями измерительных сигналов конкретной формы.
f = >N1; f2 = }N2 ;
где N1 , N2 - число дискретных отсчетов на периодах аппроксимируемых сигналов х , у позволяет определить зависимость максимального значения погрешности аппроксимации фильтром двумерного гармонического сигнала от числа интервалов дискретизации на периоде сигнала (рис. 3).
Рассмотрим пример аппроксимации реального аналитического сигнала разработанным выше сплайн - фильтром.
Такой двумерный сигнал (рис. 1), представляет собой совокупность двумерных Га-уссианов, каждый из которых описывается
7.5
5
2.5
0
n
А
Рис. 5. Поверхность, изображающая погрешности аппроксимации двумерного Гауссиана параболическим сплайном
выражением
f ( X y ) = exp
(x - X0) +(y - Y0)
2g2x 2GY
(31)
где Х(, У0 - координаты вершины Гауссиана, ахс от - параметры ширины Гауссиана
по осям X , У Аппроксимация функции (31) двумерным параболическим пятиточечным сплайном, описываемым выражениями (29) и (30) на 8 дискретных участках по осям X, У позволяют определить погрешность ее аппроксимации.
Абсолютная погрешность А аппроксимации Гауссиана параболическим сплайном представлена поверхностью, изображенной на рис. 4. Анализ этой погрешности показывает, что ее максимальное значение не пре-
вышает 4%.
Таким образом, разработанный двумерный параболический пятиточечный сплайн-фильтр может быть с успехом использован для аппроксимации реального двумерного аналитического сигнала на ограниченном числе интервалов дискретизации, что позволяет осуществить сжатие измерительной информации, и характеризуется удовлетворительной для практических нужд погрешностью аппроксимации.
СПИСОК ЛИТЕРАТУРЫ
1. Цыпкин Я.З. Теория линейных импульсных систем. М.: Физматгиз, 1963.
2. Рабинер Л., Гоулд Б. Теория и применение цифровой обработки сигналов, М.: Мир, 1978
ALGORITHMS FOR SPLINE-APPROXIMATION OF TWO DIMENSIONAL ANALYTICAL SAMPLING DATA
© 2004 P.K. Lange Samara State Technical University
The task of construction of digital filters with parabolic spline - approximation of two dimensional analytical sampling data is considered. The errors of approximation of two dimensional Gaussian Peak by such filters are determined.
