Synthesis of the reverse filter on the criterion of minimum quadratic error of the weight function mismatch of the path «direct - reverse filter»

Abstract


Based on the analysis of known approaches for the synthesis of optimal signal reconstruction algorithms, involving the use of regularizing procedures in solving incorrect inverse problems, and related computational problems, a method for constructing digital filters for solving inverse problems of signal recovery, time series and images using the approximation approach. A method for developing a model of the weight function of the inverse filter, based on the criterion of the minimum quadratic error of the weight function mismatch of the path "direct - inverse filter". The problem of restoration of blurred and defocused images in space systems of remote sensing of the Earth is formulat ed, the possibility of reducing the amount of calculations in the processing of two-dimensional data sets by moving to the solution of the one-dimensional problem is shown. Along with a General approach for the synthesis of the inverse filter, a variant of the algorithm with a number of restrictions imposed on the weight function of the inverse filter is presented. On the basis of the considered method, an approach to the problem of restoring blurred images with a known point scattering function is proposed. Approbation of algorithms on model examples and at processing of the real images received at remote sensing of the Earth is carried out. To quantify the quality of recovery, a relative mean square measure of the difference between the reference and recovered signals (images) was used. The results of testing show that the use of this approach allows to reduce the error of recovery, which gives an advantage in solving problems of approximation and data recovery.

Введение В настоящее время при решении задач обработки и интерпретации экспери- ментальных данных часто возникает необходимость рассмотрения обратной за- дачи, заключающейся в восстановлении неизвестного входного воздействия по результатам регистрации откликов на выходе средств измерения. В большинстве случаев это задача компенсации искажающего действия аппаратной функции, обеспечивающая улучшение разрешающей способности различного рода изме- рительных приборов и систем [1, 2]. В случае, когда для обработки доступна только часть искаженного сигнала, без начальных условий, задача становится недоопределенной и, соответственно, некорректной [3, 4]. Примером таких задач являются задачи восстановления смазанных и расфо- кусированных изображений [5, 6, 7]. Так, в современных космических системах дистанционного зондирования Земли (КСДЗЗ) изображение формируется с по- мощью устройств с зарядовой связью, работающих в режиме временной задерж- ки и накопления оптического сигнала. Для правильной работы таких приборов необходимо, чтобы скорость космического аппарата была точно согласована с периодом опроса светочувствительной матрицы. На практике такое равенство может нарушаться из-за ошибки вычисления скорости спутника [8]. В результате изображение подстилающей поверхности оказывается смазанным вдоль траекто- рии движения летательного аппарата. Конструктивные особенности светочув- ствительных элементов позволяют получить параметры функции рассеяния точ- ки. Способы определения параметров смаза представлены в [9, 10]. Полученные в результате несоответствия скоростей искажения имеют одну пространствен- ную составляющую, что позволяет перейти от двумерной задачи к одномерной и существенно снизить объем вычисляемых данных. Решение двумерной задачи восстановления смазанного изображения представлено, например, в [11]. Постановка задачи Каждая строка полученного с помощью КСДЗЗ смазанного изображения может быть представлена как свертка значений строки исходного изображения х(m) с известной функцией рассеяния точки h0 N0 1 xCM (m)  h0 (i)x(m  i) , (1) i0 здесь N0 - величина весовой функции фильтра - представляет собой величину смаза. Значение h0 для всех строк одинаково. Задача реконструкции полученного с помощью (1) изображения сводится к нахождению функции x’(m), в некотором роде близкой к x(m), по имеющимся значениям xCM(m). Отсутствие начальных значений xCM(m) переводит данную за- дачу в класс обратных некорректно поставленных. Существующие в настоящее время известные подходы к синтезу оптималь- ных алгоритмов реконструкции изображений либо требуют для своей реализации априорной информации, которая не всегда доступна, либо сталкиваются с вы- числительными проблемами, связанными с некорректностью обратных задач и необходимостью использования регуляризующих процедур [12, 13]. В данной статье рассматривается подход, основанный на построении модели весовой функции обратного (восстанавливающего) фильтра (рис. 1). Построение аппроксимационных моделей на базисе стохастических функций при решении задач восстановления сигналов и изображений описано в [14]. В представляемом подходе в качестве критерия адекватности модели ис- пользуется минимум квадратической погрешности рассогласования весовой функции тракта «прямой - обратный фильтр». Синтез обратных КИХ-фильтров на основе критерия моментов представлен в [15]. Решением задачи восстановления изображения будем считать нахождение функции h(), представляющей собой весовую функцию обратного фильтра, поз- воляющего получить с помощью операции свертки оценку восстановленного изображения х’(m): N 1 x'(m)  h(i)xCM (m  i) . (2) i0 h h0 Прямой (искажающий) фильтр Обратный (восстанавлива- ющий) фильтр Исходный сигнал Искаженный сигнал Восстановленный сигнал Весовая функция тракта H Рис. 1. Тракт «прямой - обратный фильтр» Рассмотрим случай, когда прямой и обратный фильтры являются КИХ- фильтрами. Весовыми функциями этих фильтров, соответственно, будут h0(i), i=0,…,N0-1 и h(i), i=0,…,N-1. Тогда весовая функция тракта «прямой - обратный фильтр» будет равна i H (i)  h0 (v)h(i  v) . (3) v0 i  0,...,N  N0  2 В идеальном случае должно быть выполнено условие H (i)  (i) где δ(i) - символ Кронекера. i  0,...,N  N0  2 , (4) Это условие принципиально выполнено быть не может, потому что обрат- ный фильтр должен быть БИХ-фильтром. Поэтому всегда будет иметь место по- грешность, квадратичное значение которой будет равно N  N0 1    i 0 (H (i)  (i))2  (H (0) 1)2  N  N0 1 H (i)2 . (5) i 1 При заданном прямом фильтре оно будет зависеть от весовой функции h(i) обратного фильтра и от N. Связь погрешности восстановления сигнала с квадратичной погрешностью Неискаженный сигнал x(m) и его восстановленное значение x’(m) связаны между собой соотношением x'(m)  N  N0 2 H (i)x(m  i) . (6) i0 Погрешность восстановления равна (m)  x'(m)  x(m) , или с учетом (6) (m)  N  N0 2 N  N0 2 H (i)x(m  i)  x(m)  i0  H (i)x(m  i)  (H (0) 1)x(m)  i1 (7) где N  N0 2  i x(m  i) i0   H (0) 1 при i  0 .  i  H (i) при i  1 В соответствии с неравенством Коши - Буняковского из (7) находим, что или 2 (m)  N  N0 2  i 0 2 (i)  N  N0 2  x2 (m  i) , i 0   N  N0 2 N  N0 2 2 (m)   (H (0) 1)2   H 2 (i)  x2 (m  i) .  Принимая во внимание (5), получим i 1  i 0 2 (m)   N  N0 2  x2 (m  i) . (8) i 0 Отсюда следует, что значение ε должно быть как можно меньшим. Синтез обратного фильтра по минимуму квадратичной погрешности Найдем решение данной задачи в частотной области. Введем в рассмотрение частотные характеристики фильтров   W0 ( j)     W ( j)    N  N0 2 N0 1 h0 (i) exp( ji) i0 N0 1 h(i) exp( ji) i0 (9) H ( j)   H (i) exp( ji)  W0 ( j)W ( j) i0 С учетом этих соотношений квадратичная погрешность (5) примет вид     2  W ( j)W0 ( j) 1W ( j)W0 ( j) 1d    (10) Для обеспечения минимума этой погрешности значения весовой функции обратного фильтра должны исходить из условия d dh(k) Из (8) с учетом (9) находим  0, k  0,...,N 1 . (11)     d  N 1   2     2 h(v)  W0 ( j) cos((k  v))d  W0 ( j)exp( jk)d dh(k)   v0  2    2       (12) С учетом этого соотношения из (10) находим, что min  1 H (0) . (13) Справедливы соотношения:      2 N0 1  2  W0 ( j) cos(k)d  h0 (v)h0 (v  k)          vk  W0 ( j) exp( jk)d  h0 (k)  2      С учетом этих соотношений формула (12) примет вид d  N 1   2 h(v)E0 ( k  v )  h0 (k), где dh(k)    v0  N0 1 E0 (k)  h0 (v)h0 (v  k), k  0,...,N0 1 . (14) vk В итоге для определения значений h(k) необходимо решить систему уравне- ний: N 1 h(v)E0 ( k  v )  h0 (k), k  0,...,N 1 . (15) v0 Алгоритм решения этой системы уравнений таков:  g(0,0)  1  (0)  E0 (0)  1 k p1   g(k,0)   (k 1)  g(k 1, p 1)E0 ( p), k  1,...,N0  2  N 1  1 0 (16) g(k,0)   (k 1)  g(k 1, p 1)E0 ( p), k  N0 1,...,N 1  p1 g(k,v)  g(k 1,v 1)  g(k 1, k 1  v)g(k,0), v  1,...,k 1  g(k, k)  1 (k )  (k 1)(1  g 2 (k,0), k  1,...,N 1 N 1 g(q, k)g(q,0) h(k)  h0 (0)  qk (q) , k  0,...,N 1 (17) Вычисляя значение минимальной погрешности из (3), имеем: H (0)  h(0)h0 (0) . Подставив сюда h(0) из (17), получаем: 2 N 1 g 2 (q,0) H (0)  h0 (0)  q0 (q) . Тогда в соответствии с (13) минимальное значение погрешности будет равно 2 N 1 g 2 (q,0) min  1  h0 (0)  q0 Для справки приведем еще такие формулы: (q) . (18) N  N0 2 2 N 1N 1    k 0 H (k )   h(i)h( j)E0 ( i  j ) i0 j0  N  N0 2  d  H 2 (k )   k 0   dh(v) N 1 2 h( j)E0 ( v  j ) j0 (19)  N 1N 1    1  2H (0)   h(i)h( j)E0 ( i  j )  i0 j0    Синтез обратного фильтра по минимуму квадратичной погрешности при внесении ряда ограничений, наложенных на его весовую функцию Исходя из принципов физической реализуемости весовая функция синтезируемого фильтра должна удовлетворять условиям N 1 (i,v)h(i)  (v),v  0,..., p . (20) i0 Так как функции H(i) и h(i) однозначно взаимосвязаны, то условия, нало- женные на функцию H(i), могут быть сведены к условиям, наложенным на h(i). Например, если функция H(i) должна удовлетворять условиям N  N0 2 (i,v)H (i)  (v),v  0,..., p i0 то эти условия могут быть записаны в виде (20), где N0 1 , (21)  (i, v)   h0 (k) (i  k, v), v  0,..., p . (22) k 0 Задачу синтеза обратного фильтра по минимуму квадратичной погрешности с учетом ограничений физической реализуемости, наложенных на его весовую функцию, будем решать методом Лагранжа. Для этого рассмотрим функцию p N 1 f (h(0),...,h(N 1))    2 Av (i,v)h(i) , где Av - коэффициенты Лагранжа. v0 i0 Значения h(i) будем определять из условий минимума функции f(h(0),…,h(N-1)): df (h(0),...,h(N 1))  0, k  0,...,N 1. dh(k) В результате на основании рассмотренных ранее выводов для обеспечения этих условий получим следующую систему уравнений: N 1 N h(i)E0 ( k  i )  h0 (k)   Av(k,v),k  0,...,N 1. i0 v0 Решив ее, получим N 1 g(q,i)  p  h(i)    g(q,0)h0 (0)   Ak D(q, k) . (23) Здесь qi (q)    k 0  q D(q, k)  g(q,v)(v, k) ; (24) v0 g(q,k), Ψ(q) - рассмотренные ранее величины. Подставив h(i) из (23) в (20), получим следующую систему уравнений для определения значений коэффициентов Лагранжа: p N 1 D(q,v)D(q, k) N 1D(q,v)  Ak  (q)   (q) g(q,0)h0 (0)  (v),v  0,..., p . k 0 q0 q0 Для простоты дальнейших выкладок введем обозначения: N 1 D(q, v)D(q, k) f (v, k)   q0 N 1 D(q, v) (q) (25) F (v)   q0 (q) g(q,0)h0 (0)  (v) В этом случае система уравнений примет вид p  Ak f (v, k)  F (v) . (26) k 0 Алгоритм решения этой системы уравнений будет таким. Сначала по форму- лам (24) и (25) определим значения f(v,k), v=0,…,p; k=0,…,v и F(v), v=0,…,p. За- тем выполним вычисления по формулам:  k 1 E(v, k)  f (v, k)  C(k, q)E(v, q)   C(v, k)   E(v, k) E(k, k) q0 (27) k  0,...,v; v  0,..., p     v1 (v)  F (v)  C(v, k)(k)   v  0,..., p k 0 A  (k)  k E(k, k) p  AvC(v, k),k  p,...,0 vk 1 Подставив найденные значения Ak в (23), получим искомые значения весовой функции обратного фильтра. Далее по формуле (2) получим значение восстанов- ленного изображения х’(m). Апробация результатов Для апробации алгоритма было взято изображение, полученное в процессе дистанционного зондирования Земли. Благодаря конструктивным особенностям регистрирующего аппаратуры изображение имеет 1024 градации серого, что в 4 раза выше значений, принятых в распространенных форматах хранения графиче- ских файлов на персональных компьютерах. Из тестового изображения был взят тестовый фрагмент 480 на 285 пикселей, над которым и проводились экспери- менты (рис. 2, а). Каждая строка тестового изображения была обработана фильтром (1) с весо- вой функцией: 1 h0  N . 0 Таким образом, был выполнен смаз изображения вдоль горизонтальной оси на различное количество пикселей N0 = 3, 4, 5, …, 10. Пример полученного сма- занного изображения приведен на рис. 2, б. Рис. 2. Эталонное изображение Рис. 3. Изображение, смазанное на 10 пикселей Далее были предприняты попытки восстановления по алгоритмам (16) и (27). Для количественной оценки качества восстановления использовалась от- носительная среднеквадратическая погрешность (ОСП), вычисляемая по форму- ле ОСП  M 1  (x '( j)  x( j))2 , j 0  M 1 x( j)2 j 0 (28) где х(j) - значение пикселя строки эталонного изображения; х’(j) - значение пикселя строки восстановленного изображения; M - количество пикселей в строке. Данный способ оценки возможен благодаря наличию и эталонного, и вос- становленного изображения. Полученная величина показывает степень отклоне- ния результата от исходного изображения. Другие способы оценки качества вос- становления изображений приведены, например, в [16]. Значения ОСП смазанных изображений без восстановления приведены в табл. 1. Таблица 1 ОСП невосстановленных изображений N0 3 4 5 6 7 8 9 10 ОСП 0,0505 0,0565 0,0672 0,0743 0,0833 0,0913 0,0974 0,1041 Восстановление изображений по алгоритму (16) и (17) без наложенных ограничений. Анализ формул (14) - (19) показывает, что результат восстановле- ния зависит от значений весовых функций прямого и обратного фильтров: N и N0. Некоторые значения погрешности восстановления приведены в табл. 2. Таблица 2 Результаты восстановления по алгоритму (16), 17) N0 N ОСП N0 N ОСП N0 N ОСП 3 3 0,044 5 5 0,0617 7 7 0,062 3 5 0,043 5 9 0,0494 7 13 0,056 3 7 0,043 5 13 0,0602 7 19 0,060 4 4 0,061 6 6 0,0593 8 8 0,075 4 7 0,053 6 11 0,0543 8 15 0,066 4 10 0,060 6 16 0,0539 8 22 0,077 Графики зависимости ОСП от N0 показаны на рис. 4. При увеличении длины весовой функции обратного фильтра N процесс восстановления становится более неустойчивым (рис. 5), что не позволяет значительно увеличивать N. Рис. 4. Зависимость ОСП от N0 и N: 1 - значение ОСП при N0 = 8; 2 - значение ОСП при N0 = 4; 3 - значение ОСП при N0 = 5; 4 - значение ОСП при N0 = 7; 5 - значение ОСП при N0 = 6; 6 - значение ОСП при N0 = 3 Рис. 5. Фрагмент строки эталонного и восста- новленного изображений при N0 = 3: x(m), x’(m) - значения пикселей; m - номер пикселя; 1 - восстановленное изображение при N = 13; 2 - эталонное изображение Восстановление изображений по алгоритму (27) c наложенными ограниче- ниями. Для выполнения условий (21) и (22) были выбраны следующие соотно- шения: H (0)  h(0)h0 (0) 1 ; N  N0 1 H (i)iv  (v),v  0,..., p 1; i0  (i, v)   h (k)(i  k) . N0 1 v 0 k 0 Анализ алгоритма построения весовой функции обратного фильтра (20) - (27) показывает, что результат восстановления зависит от значений весовых функций прямого и обратного фильтров: N и N0, а также p из (22). При проведе- нии экспериментов с различными значениями N0, N и p наименьшие значения ОСП получались при N= 2*p-1. Зависимость ОСП от значений N и p показана на рис. 6. Рис. 6. Зависимость ОСП от N: 1 - N0 = 6; 2 - N0 = 5; 3 - N0 = 4; 4 - N0 = 3 Рис. 7. Зависимость ОСП от порядка модели с различными N0: 1 - N0 = 8; 2 - N0 = 5; 3 - N0 = 4; 4 - N0 = 7; 5 - N0 = 6; 6 - N0 = 3 Некоторые полученные результаты представлены в табл. 3. Таблица 3 Результаты восстановления по алгоритму (27) N0 P ОСП N0 P ОСП N0 P ОСП 3 2 0,034 5 2 0,044 7 2 0,055 3 3 0,044 5 3 0,046 7 3 0,068 3 4 0,054 5 4 0,065 7 4 0,067 4 2 0,029 6 2 0,046 8 2 0,063 4 3 0,042 6 3 0,062 8 3 0,074 4 4 0,065 6 4 0,069 8 4 0,068 Выводы На основе анализа известных подходов к синтезу оптимальных алгоритмов реконструкции сигналов, основанных на использовании регуляризующих проце- дур при решении некорректных обратных задач, и связанных с этим вычисли- тельных проблем предложен метод построения цифровых фильтров для решения обратных задач восстановления сигналов, временных рядов и изображений с использованием аппроксимационного подхода. Сформулирована постановка задачи восстановления смазанных и расфоку- сированных изображений в космических системах дистанционного зондирования Земли, показана возможность снижения объема вычислений при обработке дву- мерных массивов данных путем перехода к решению одномерной задачи. Разработаны и исследованы алгоритмы синтеза обратного фильтра для вос- становления сигналов с известной весовой функцией прямого фильтра на основе критерия минимума квадратической погрешности рассогласования весовой функции всего тракта «прямой - обратный фильтр» для реконструкции смазан- ных изображений. Проведена апробация алгоритмов на модельных примерах и при обработке реальных изображений, полученных при дистанционном зондировании Земли. Для количественной оценки качества восстановления использовалась относи- тельная среднеквадратическая мера различия эталонного и восстановленного сигналов (изображений). Приведенные результаты апробации (табл. 1, 2, 3) пока- зывают, что использование данного подхода позволяет уменьшить погрешность восстановления и аппроксимации, что дает преимущество при решении задач аппроксимации и восстановления данных.

Vitaly I Batishchev

Samara State Technical University

244, Molodogvardeyskaya st., Samara, 443100, Russian Federation (Dr. Sci. (Techn.)), Professor

Igor I Volkov

Samara State Technical University

244, Molodogvardeyskaya st., Samara, 443100, Russian Federation (Ph.D. (Techn.)), Associate Professor

A G Zolin

Samara State Technical University

Email: zolin.a.g@gmail.com
244, Molodogvardeyskaya st., Samara, 443100, Russian Federation (Ph.D. (Techn.)), Associate Professor

  1. Клебанов Я.М., Карсаков А.В., Хонина С.Н., Давыдов А.Н., Поляков К.А. Компенсация аберраций волнового фронта в телескопах космических аппаратов с регулировкой температурного поля телескопа // Компьютерная оптика. - 2017. - Т. 41. - № 1. - С. 30-36.
  2. Tokovinin A. DONUT: measuring optical aberrations from a single extrafocal image / A. Tokovinin, S. Heathcote // Publications of the Astronomical Society of the Pacific. - 2006. - Vol. 118(846). - Pp. 1165-1175.
  3. Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. - М.: Наука, 1979.
  4. Engl H.W., Hanke M., Neubauer A. Regularization of Inverse Problems. - Dordrecht, Kluwer Acad. Publ., 1996. 322 p.
  5. Василенко Г.И., Тараторин A.М. Восстановление изображений. - М.: Радио и связь, 1986. - 304 с.
  6. Сизиков B.C., Римских М.В., Мирджамолов Р.К. Реконструкция смазанных и зашумленных изображений без использования граничных условий // Оптический журнал. - 2009. - Т. 76. - № 5. - С. 38-46.
  7. Гонсалес Р., Вудс Р. Цифровая обработка изображений: 3-е изд., испр. и доп. - М.: Техносфера, 2012. - 1104 с.
  8. Кузнецов П.К., Семавин В.И., Солодуха А.А. Алгоритм компенсации скорости смаза изображения подстилающей поверхности, получаемого при наблюдении Земли из космоса // Вестник Самарского государственного технического университета. Сер. Технические науки. - 2005. - Вып. 37. - С. 150-157.
  9. Егошкин Н.А., Еремеев В.В. Коррекция смаза изображений в системах космического наблюдения земли // Цифровая обработка сигналов. - 2010. - № 4. - С. 28-32.
  10. Кузнецов П.К., Мартемьянов Б.В., Мятов Г.Н., Юдаков А.А. Методика вычисления оценок параметров смаза изображений, получаемых целевой аппаратурой КАН типа «Ресурс» // 7-я Международная научно-техническая конференция «К.Э. Циолковский - 160 лет со дня рождения. Космонавтика. Радиоэлектроника. Геоинформатика». Рязань, 2017. - С. 344-350.
  11. Волков И.И., Золин А.Г. Решение двумерной обратной задачи восстановления смазанного изображения // Вестник Самарского государственного технического университета. Сер. Технические науки. - 2013. - Вып. 39. - С. 223-226.
  12. Батищев В.И., Волков И.И., Золин А.Г. Синтез фильтров для восстановления смазанных изображений с использованием методов регуляризации // Проблемы управления и моделирования в сложных системах (ПУМСС-2013): Труды XV Международной конференции, ИПУСС РАН. - Самара, 2013. - С. 615-618.
  13. Батищев В.И., Золин А.Г., Косарев Д.Н., Романеев А.Е. Аппроксимационный подход к решению обратных задач анализа и интерпретации экспериментальных данных // Вестник Самарского государственного технического университета. Сер. Технические науки. - 2006. - Вып. 40. - С. 57-65.
  14. Батищев В.И., Волков И.И., Золин А.Г. Использование стохастического базиса в задачах восстановления сигналов и изображений // Автометрия. - 2017. - № 4. - С. 127-134.
  15. Батищев В.И., Волков И.И., Золин А.Г. Синтез цифровых КИХ-фильтров для решения задач восстановления сигналов с использованием критерия моментов // Вестник Самарского государственного технического университета. Сер. Технические науки. - 2012. - Вып. 36. - С. 98-105.
  16. Avcibas I., Sankur B., Sayood K. Statistical evaluating of image quality measures // Journal of Electronic Imaging. - April 2002. - Vol. 11, № 2. - Рp. 206-223.

Views

Abstract - 9

PDF (Russian) - 2

Cited-By


PlumX

Dimensions

Refbacks

  • There are currently no refbacks.

Copyright (c) 2018 Samara State Technical University

Creative Commons License
This work is licensed under a Creative Commons Attribution 4.0 International License.

This website uses cookies

You consent to our cookies if you continue to use our website.

About Cookies