Радиотехника и электроника, 2023, T. 68, № 7, стр. 632-649
Реконструкция сигналов линейных стационарных систем с помощью тестовых испытаний
А. В. Новиков-Бородин *
Институт ядерных исследований РАН
117312 Москва, пр. 60-летия Октября, 7а, Российская Федерация
* E-mail: novikov.borodin@gmail.com
Поступила в редакцию 21.03.2022
После доработки 28.06.2022
Принята к публикации 27.09.2022
- EDN: WPIOFD
- DOI: 10.31857/S0033849423070082
Аннотация
Предложены и проанализированы методы тестовых испытаний, предназначенные для математической реконструкции входных сигналов линейных стационарных систем обработки информации из искаженных ими и зашумленных выходных сигналов. Показано, что использование данных тестовых испытаний позволяет осуществить не слепую реконструкцию без определения функций искажений сигналов или аппаратных функций систем обработки, которые, в общем случае, могут принадлежать к классу обобщенных функций. Рассмотрены особенности применения техники регуляризации в методах тестовых испытаний при решении некорректно поставленных и плохо обусловленных задач реконструкции реальных недетерминированных сигналов. Проанализированы критерии выбора тестовых сигналов. Представлены результаты численных экспериментов по восстановлению одномерных сигналов и двумерных изображений при разных уровнях зашумления.
ВВЕДЕНИЕ
Реальные системы обработки информации, такие как системы измерения, усиления, приема-передачи и др. преобразований сигналов и изображений, несовершенны, что приводит к систематическим искажениям обрабатываемых ими данных и к появлению в них шумов и неопределенностей. Уменьшение искажений и неопределенностей в выходных сигналах путем усовершенствования оборудования связано с решением сложных научно-технических задач и с созданием дорогостоящей аппаратуры. Доступной альтернативой аппаратным методам могут стать математические методы восстановления входных сигналов из искаженных и зашумленных выходных на основе математической модели процесса обработки информации. Большинство подобных инверсных задач реконструкции являются некорректно поставленными и плохо обусловленными, универсальных подходов к их решению не существует и разработка эффективных в том или ином случае методов восстановления обрабатываемой информации, является актуальным и востребованным направлением исследований (см., например, [1–11]). Из существующих подходов к решению некорректных задач можно выделить вариационные [12–17] и итерационные [18–25] методы регуляризации, а также методы с использованием вейвлет-преобразований [26–28].
В данной работе рассматриваются линейные стационарные (ЛС) системы обработки информации, имеющие важное значение и широко используемые при измерениях, усилении, приеме-передачи и других преобразованиях сигналов (см., например, [29–37]). При математической реконструкции внешних воздействий $v\left( x \right)$ из откликов $w\left( x \right)$ ЛС-систем, поведение которых характеризуется аппаратной функцией $u\left( x \right)$, используется приближенное вследствие недетерминированности реальных сигналов и неопределенности аппаратной функции уравнение типа свертки $w\left( x \right) \approx \left( {u * v} \right)\left( x \right)$.
Определение аппаратной функции является первоочередной задачей еще на стадии постановки задачи реконструкции, но часто трудности ее определения превосходят проблемы самой реконструкции сигналов. Во многом это связано с тем, что в отличие от сигналов, описываемых основными функциями, аппаратные функции зачастую относятся к классу обобщенных (см., например, [38–40] и примеры в разд. 1), что делает неприемлемыми многие методы их определения, такие, например, как численные или экспериментальные, в которых, по сути, предпринимается попытка аппроксимации обобщенных функций основными. Аналитические методы определения аппаратных функций с помощью математического моделирования работы ЛС-систем возможны далеко не всегда, так как моделирование является крайне непростой задачей в случае сложных ЛС-систем и может быть неточным в случае быстродействующих и прецизионных из-за проблем учета паразитных параметров и взаимосвязей элементов приборов. Подобные проблемы определения аппаратных функций еще на стадии постановки задачи приводят к тому, что математические методы реконструкции сигналов для устранения систематических искажений в сложных системах обработки информации практически не применяются, за исключением отдельных случаев, когда аппаратная функция относится к классу основных или, реже, некоторых разрывных функций.
В данной работе для решения задачи математической реконструкции сигналов, искаженных ЛС-системами, предлагаются методы тестовых испытаний (ТИ), в которых моделирование ЛС-систем и определение их аппаратных функций не требуются, что устраняет многие из указанных выше проблем. Реконструкция сигналов в ТИ-методах осуществляется не из ЛС-уравнения, а из уравнения тестовых испытаний, которое тоже является уравнением типа свертки, но содержит только сигналы из класса основных функций, что значительно упрощает задачу реконструкции и позволяет использовать многие существующие методы решения некорректно поставленных задач. В данной работе для реконструкции сигналов и изображений с помощью тестовых испытаний была использована техника регуляризации по Тихонову [12–17], которая удобна для контроля за процессом регуляризации и для иллюстрации применения предлагаемых ТИ-методов. На основе техники регуляризации рассматриваются варианты восстановления непрерывных сигналов и дискретных данных прямыми и спектральными (с использованием преобразований Фурье) методами регуляризации. Особое внимание в работе уделено практическим вопросам выбора тестовых сигналов, стабилизаторов и параметров регуляризации, напрямую влияющих на точность реконструкции. Представлены примеры восстановления одномерных сигналов и двумерных изображений.
1. ОСНОВНОЕ УРАВНЕНИЕ ТЕСТОВЫХ ИСПЫТАНИЙ
Линейные стационарные системы имеют важное значение и широко используются при измерениях, усилении, приеме-передачи и другой обработке сигналов (см., например, [29–37]). Поведение ЛС-системы в исследуемой области $D$ характеризуется аппаратной функцией $u\left( x \right)$, и ее отклик $w\left( x \right)$ (выходной сигнал) на внешнее воздействие $v\left( x \right)$ (входной сигнал) в случае недетерминированных сигналов приближенно описывается уравнением типа свертки
(1)
$w\left( x \right) \approx \left( {u * v} \right)\left( x \right) = \mathop \smallint \limits_D^{\text{\;}} u\left( {x - \xi } \right)v\left( \xi \right){\text{d}}\xi ,$которое будем называть уравнением линейных стационарных систем или, сокращенно, ЛС-уравнением. Неопределенность ЛС-уравнения для реальных недетерминированных сигналов $w\left( x \right)$, $v\left( x \right)$ и функции $u\left( x \right)$ с шумами и неопределенностями ${{\Delta }}w\left( x \right)$, ${{\Delta }}v\left( x \right)$ и ${{\Delta }}u\left( x \right)$ составляет
При $\eta \left( x \right) = 0$ ЛС уравнение (1) переходит в строгое уравнение свертки $w\left( x \right) = \left( {u * v} \right)\left( x \right)$ для детерминированных сигналов.
Помимо интерпретации функции $u\left( x \right)$ как аппаратной или приборной (в радиотехнике и электронике), в других областях науки и техники в зависимости от основных процессов, описываемых уравнением (1), она интерпретируется как весовая (в теории систем, экспериментальной физике и математической статистике), как функция рассеяния (информационные системы обработки сигналов и изображений) и др. В математике функция $u\left( x \right)$ называется ядром свертки. Далее функцию $u\left( x \right)$ будем называть системной или, подчеркивая специфику данной работы, аппаратной, но понимая ее в общем смысле, т.е. включая все физические проявления линейного стационарного процесса (1) при работе ЛС-систем.
Задача реконструкции воздействий $v\left( x \right)$ из откликов $w\left( x \right)$ математическими методами из приближенного уравнения типа свертки (1) является инверсной, некорректно поставленной и, зачастую, плохо обусловленной, т.е. если пренебречь неопределенностями в (1), то найденное точное решение может сколь угодно много отличаться от реконструируемого воздействия, даже если исходные неопределенности были сколь угодно малы. В общем случае, решение может вообще не существовать или быть неоднозначным.
В существующих методах решения обратных некорректно поставленных задач ищется не точное, а приближенное решение, степень приближения которого к истинному решению критически зависит от неопределенности исходных данных (функций $w\left( x \right)$ и $u\left( x \right)$). Особые затруднения вызывает точность определения аппаратных функций $u\left( x \right)$, которые в отличие от сигналов $w\left( x \right)$ и $v\left( x \right)$, описываемых основными функциями, зачастую относятся к классу обобщенных. Так, если, например, поведение ЛС-системы описывается обыкновенным дифференциальным уравнением n-го порядка с постоянными коэффициентами:
то аппаратная функция такой системы будет равна сумме производных от дельта-функции $~\delta \left( x \right)$:
что следует из представления $v\left( x \right) = \left( {\delta * v} \right)\left( x \right)$ и линейных свойств свертки (подробнее см., например, в [38–40]), т.е. будет обобщенной функцией. В частности, последний член $u\left( x \right) = {{a}_{0}}\delta \left( x \right)$ соответствует аппаратной функции “простейшей” “идеальной” измерительной или усилительной ЛС-системы, не искажающей форму входного сигнала, так как
Принадлежность аппаратных функций к классу обобщенных делает неприемлемыми для ее определения многие существующие подходы. Аналитическое описание аппаратных функций на основе моделирования работы ЛС-систем возможно далеко не всегда, так как является крайне непростой задачей в случае сложных систем и неточной в случае быстродействующих и прецизионных из-за проблем учета паразитных параметров и взаимосвязей элементов приборов. Принципиальные проблемы возникают при попытках определения аппаратных функций численными расчетами или экспериментальным путем, где, по сути, предпринимается попытка аппроксимации обобщенных функций основными, что или невозможно, или приведет к погрешностям, неприемлемым для использования аппаратных функций при дальнейшей реконструкции сигналов из ЛС-уравнений (1).
Например, определить аппаратную функцию $u\left( x \right)$ прямым экспериментом как отклик ЛС-системы на входное воздействие в виде дельта-функции (при $v\left( x \right) = \delta \left( x \right)$: $w\left( x \right) \approx \left( {u * \delta } \right)\left( x \right) = u\left( x \right)$) невозможно из-за невозможности сгенерировать входной сигнал с нулевой длительностью и бесконечной амплитудой, а также из-за ограничений на допустимый диапазон входных сигналов реальных систем. Если же вместо дельта-функции использовать некоторое тестовое воздействие $v\left( x \right) = {{v}_{c}}\left( x \right)$ из допустимого диапазона входных сигналов ЛС-системы, найти отклик на него ${{w}_{c}}\left( x \right)$ и попробовать реконструировать аппаратную функцию $u\left( x \right)$ математическими методами из ЛС-уравнения (1) для этих сигналов: ${{w}_{c}}\left( x \right) \approx \left( {u * {{v}_{c}}} \right)\left( x \right)$, то, поскольку эта задача является некорректно поставленной, методическая погрешность приближенной реконструкции аппаратной функции существующими методами с учетом ее принадлежности к классу обобщенных может быть слишком велика и полученное решение будет непригодно для дальнейшего использования при реконструкции сигналов. Кроме того, требования к тестовым сигналам могут оказаться слишком жесткими и нереализуемыми на практике. Например, чтобы реконструировать аппаратную функцию $u\left( x \right) = a\delta \left( x \right)$ “идеального” усилителя, не искажающего форму входного сигнала ($w\left( x \right) \approx a\left( {\delta * v} \right)\left( x \right) = av\left( x \right)$), частотный диапазон тестовых сигналов должен перекрывать бесконечный диапазон частот, что на практике для реальных сигналов недостижимо.
В методе тестовых испытаний, предлагаемом в данной работе для реконструкции сигналов, также используются тестовое воздействие ${{{v}}_{c}}\left( x \right)$ и отклик ${{w}_{c}}\left( x \right)$, находящиеся в допустимых диапазонах ЛС-систем и, следовательно, связанные уравнением (1): ${{w}_{c}}\left( x \right) \approx \left( {u * {{v}_{c}}} \right)\left( x \right)$, но определение аппаратной функции в ТИ-методах не требуется, что позволяет частично или полностью решить указанные выше проблемы. Действительно, свертывая обе части уравнения (1) с функцией ${{v}_{c}}\left( x \right)$ и учитывая коммутативность и ассоциативность свертки и уравнение связи тестовых сигналов, можно исключить аппаратную функцию и получить основное ТИ-уравнение
(2)
$\left( {w * {{v}_{c}}} \right)\left( x \right) \approx \left( {{{w}_{c}} * v} \right)\left( x \right).{\text{\;\;}}$ТИ-уравнение (2), как и (1), является приближенным уравнением типа свертки, но не содержит аппаратных функций, а только сигналы из класса основных функций, что дает возможность восстановления воздействий $v\left( x \right)$ из него многими существующими методами решения некорректных задач. Неопределенность ТИ-уравнения зависит от неопределенностей ${{\Delta }}w\left( x \right)$, ${{\Delta }}v\left( x \right)$ и ${{\Delta }}{{w}_{c}}\left( x \right)$, ${{\Delta }}{{v}_{c}}\left( x \right)$ в сигналах $w\left( x \right)$, $v\left( x \right)$ и ${{w}_{c}}\left( x \right)$, ${{v}_{c}}\left( x \right)$ и равна:
Тестовые сигналы, содержащие информацию о работе ЛС-систем, можно определить многими способами, например, с помощью численных расчетов, что даже для сложных ЛС-систем не вызывает больших затруднений. При экспериментальном определении тестовых сигналов ЛС-системы рассматриваются как “черный ящик” и их анализ вообще не требуется (достаточно информации о их линейности и стационарности). При этом в тестовых сигналах автоматически учитываются паразитные взаимосвязями элементов устройств, что дает возможность использовать ТИ-методы для реконструкции сигналов быстродействующих и прецизионных ЛС-систем любой сложности.
К сложным ЛС-системам можно отнести системы, состоящие из произвольного числа последовательно и/или параллельно соединенных ЛС-подсистем. Действительно, последовательно соединенные ЛС-системы являются ЛС-системой в силу ассоциативности свертки:
а параллельно соединенные – в силу ее линейности:
Подобные составные ЛС-системы очень распространены. Такими системами являются, например, электронные схемы из любого количества пассивных элементов (резисторов, конденсаторов, индуктивностей, линий связи и т.д.) (см. пример в разд. 6) и активных элементов (операционных усилителей, транзисторов и др.), работающих в линейном режиме.
ТИ-методами можно реконструировать не только сигналы на входе составной ЛС-системы, но и, в силу линейности свертки, сигналы на входе или выходе любой из ее ЛС-подсистем. При реконструкции сигналов на входе ЛС-подсистем тестовые сигналы должны быть определены в соответствующих узлах. При определении тестовых сигналов экспериментальным путем не требуется анализ ни составной системы, ни ее подсистем.
Классификация и обзор существующих методов решения некорректно поставленных задач реконструкции сигналов из класса основных функций из ТИ-уравнений (2) не входит в рамки данной работы. Далее для реконструкции сигналов и изображений будет использоваться техника регуляризации по Тихонову [12–17], которая удобна для контроля за процессом регуляризации и для иллюстрации использования ТИ-методов. Будут рассмотрены варианты восстановления непрерывных сигналов и дискретных данных прямыми и спектральными (с использованием преобразований Фурье) методами регуляризации. Особое внимание будет уделено практическим вопросам выбора тестовых сигналов, стабилизаторов и параметров регуляризации, напрямую влияющих на точность реконструкции.
2. РЕШЕНИЕ УРАВНЕНИЙ ТИПА СВЕРТКИ
Задача реконструкции воздействий из приближенных уравнений типа свертки, таких как уравнение ЛС-систем (1) и уравнение тестовых испытаний (2), является инверсной, некорректно поставленной и, зачастую, плохо обусловленной. В ТИ-методах для реконструкции сигналов используется ТИ-уравнение (2), которое, как и ЛС-уравнение (1) является уравнением типа свертки, но не содержит функций из класса обобщенных, а содержит лишь сигналы, относящиеся к классу основных функций. Если в (1) ограничиться рассмотрением функций из класса основных, то методы решения ЛС-уравнения (1) принципиально не будут отличаться от методов решения ТИ-уравнения (2), и, заменив в известных решениях уравнения (1) функции $w\left( x \right)$ на $\left( {w * {{v}_{c}}} \right)\left( x \right)$ и $u\left( x \right)$ на ${{w}_{c}}\left( x \right)$, можно перейти от решений ЛС-уравнения (1) к решениям ТИ-уравнения (2). Далее эти замены будем называть переходными.
В технике регуляризации по Тихонову регуляризированным решением ${{v}_{r}}\left( x \right)$ уравнения типа свертки (1) считается функция, при которой достигается минимум функционала:
(3)
${{v}_{r}}\left( x \right) = \mathop {\min }\limits_v \left( {{{{\left\| {w - u * v} \right\|}}^{2}} + P\left[ {v,p} \right]} \right).$Здесь $\left\| {f\left( x \right)} \right\|$ – норма функции $f\left( x \right)$, обычно понимаемая как евклидова норма
Стабилизирующий член $P\left[ {v,p} \right]$, где $p = \left\{ {{{p}_{i}}} \right\}$ – параметры регуляризации, обеспечивает устойчивость решения при варьировании исходных данных. При оптимальном выборе стабилизатора $P\left[ {v,p} \right]$ достигается минимум невязки $\left\| {{{v}_{r}}\left( x \right) - v\left( x \right)} \right\|$ регуляризированного решения и реконструируемого воздействия. Вопросы выбора стабилизатора $P\left[ {v,p} \right]$ и параметров регуляризации $p$ будут рассмотрены ниже.
Спектр ${{V}_{r}}\left( \omega \right)$ или преобразование Фурье регуляризированного решения ${{v}_{r}}\left( x \right)$ (3) можно выразить в явном виде через спектры $~W\left( \omega \right)$, $U\left( \omega \right)$, $~V\left( \omega \right)$ функций $w\left( x \right)$, $u\left( x \right)$, $v\left( x \right)$:
(4)
$\begin{gathered} {{V}_{r}}\left( \omega \right) = \frac{{U{\text{*}}\left( \omega \right)W\left( \omega \right)}}{{U{\text{*}}\left( \omega \right)U\left( \omega \right) + Q\left( {\omega ,q} \right)}} = \\ = \frac{{U{\text{*}}\left( \omega \right)W\left( \omega \right)}}{{{{{\left| {U\left( \omega \right)} \right|}}^{2}} + Q\left( {\omega ,q} \right)}}. \\ \end{gathered} $Здесь $Q\left( {\omega ,q} \right)$ – стабилизирующий член с параметрами регуляризации $q = \left\{ {{{q}_{i}}} \right\}$, ${{\left| {U\left( \omega \right)} \right|}^{2}} = U{\text{*}}\left( \omega \right)U\left( \omega \right)$, $U{\text{*}}\left( \omega \right)$ – комплексное сопряжение к $U\left( \omega \right)$ (в многомерном случае сопряженно-транспонированный оператор). Cвязь между стабилизаторами в (3) и (4) можно найти, используя равенство Планшереля
где ${{\Omega }}$ – пространство частот.
Методы (3) и (4) поиска регуляризированного решения эквивалентны в силу линейности Фурье-преобразования, но обладают некоторой спецификой при применении на практике, поэтому далее будем рассматривать оба способа, называя метод (4) спектральным, а (3) – прямым.
Область реконструкции спектра $V\left( \omega \right)$ воздействия v(x) ограничена диапазоном реконструкции
где спектр отклика $W\left( \omega \right)$ может быть выделен на фоне спектра $N\left( \omega \right)$ неопределенностей $\eta \left( x \right)$ уравнения (1). Здесь коэффициент $\gamma > 0$ характеризует методы выделения спектра. Используя равенство Планшереля, можно оценить ограничение $\delta v_{r}^{\eta }$ на точность $\delta {{v}_{r}}$ реконструкции воздействия $v\left( x \right)$ в диапазоне ${{{{\Omega }}}_{\eta }}$:
При стремлении уровня неопределенности в (1) к нулю: $\left\| {\eta \left( x \right)} \right\| \to 0$ (при этом в силу равенства Планшереля норма его спектра также стремится к нулю: $\left\| {N\left( \omega \right)} \right\| \to 0$), для любого $\gamma > 0$ точность $\delta {{v}_{r}}$ будет ограничена полосой пропускания ЛС-системы ${{{{\Omega }}}_{U}} = \left\{ {\omega :U\left( \omega \right) \ne 0} \right\}$, так как
Поскольку для функций из класса основных методы решения ТИ-уравнения (2) принципиально не отличаются от методов решения ЛС-уравнения (1), то, заменяя в (3) функцию w(x) на $\left( {w * {{v}_{c}}} \right)\left( x \right)$ и функцию $u\left( x \right)$ на ${{w}_{c}}\left( x \right)$, получим функционал для регуляризированного решения ТИ-уравнения (2):
(6)
${{v}_{r}}\left( x \right) = \mathop {\min }\limits_{v} \left( {{{{\left\| {w * {{v}_{c}} - {{w}_{c}} * v} \right\|}}^{2}} + P\left[ {v,p} \right]} \right).$Аналогично, заменяя в (4) спектр $W\left( \omega \right)$ на $W\left( \omega \right){{V}_{c}}~\left( \omega \right)$ и спектр $U\left( \omega \right)$ на ${{W}_{c}}\left( \omega \right)$, получим спектр регуляризированного решения:
(7)
${{V}_{r}}\left( \omega \right) = \frac{{W_{c}^{*}\left( \omega \right)W\left( \omega \right){{V}_{c}}\left( \omega \right)}}{{W_{c}^{*}\left( \omega \right){{W}_{c}}\left( \omega \right) + Q\left( {\omega ,q} \right)}}.~~$По аналогии с (5), используя равенство Планшереля, можно оценить ограничение $\delta v_{r}^{c}$ на уровень относительной погрешности $\delta {{v}_{r}}$ реконструкции воздействия $v\left( x \right)$ из ТИ-уравнения (2) в диапазоне реконструкции $\Omega _{\eta }^{c} = \left\{ {\omega :\gamma \left| {{{N}_{c}}\left( \omega \right)} \right| \leqslant \left| {W\left( \omega \right){{V}_{c}}\left( \omega \right)} \right|} \right\}$:
При стремлении уровня неопределенности ТИ-уравнения (2) к нулю $\left\| {{{\eta }_{c}}\left( x \right)} \right\| \to 0$ (в силу равенства Планшереля: $\left\| {{{N}_{c}}\left( \omega \right)} \right\| \to 0$), точность реконструкции будет ограничена не только полосой пропускания ЛС-системы, как в (5), но и частотным диапазоном ${{V}_{c}}\left( \omega \right)$ тестового воздействия, или, поскольку $U\left( \omega \right){{V}_{c}}\left( \omega \right) \to {{W}_{c}}\left( \omega \right)$, – частотным диапазоном ${{W}_{c}}\left( \omega \right)$ тестового отклика. При этом
Отметим, что неопределенность ТИ-уравнения (2) можно уменьшить, уменьшив случайные шумы и неопределенности в тестовых сигналах ${{w}_{c}}\left( x \right)$ и ${{v}_{c}}\left( x \right)$ статистическими методами, например, повторяя тестовые испытания необходимое количество раз при проведении тестовых испытаний экспериментальным путем. Таким же образом можно уменьшить случайные шумы и неопределенности в реконструируемом отклике $w\left( x \right)$, если имеется возможность повторить измерения и вид воздействий при этом не меняется.
Уравнение тестовых испытаний (2) можно использовать не только для реконструкции воздействий, но и для моделирования отклика $w\left( x \right)$ ЛС-систем на заданное воздействие $v\left( x \right)$. При этом в ТИ-уравнении (2) известными будут воздействие $v\left( x \right)$ и тестовые сигналы ${{v}_{c}}\left( x \right)$ и ${{w}_{c}}\left( x \right)$, а неизвестным будет отклик $w\left( x \right)$. Заменяя в (3) $w\left( x \right)$ на $\left( {v * {{w}_{c}}} \right)\left( x \right)$ и $u\left( x \right)$ на ${{v}_{c}}\left( x \right)$, получим функционал для моделирования отклика $w\left( x \right)$:
Аналогично, заменив в (4) $W\left( \omega \right)$ на $V\left( \omega \right){{W}_{c}}\left( \omega \right)$ и $U\left( \omega \right)$ на ${{V}_{c}}\left( \omega \right)$, получим выражение для регуляризированного спектра ${{W}_{r}}\left( \omega \right)$ моделируемого отклика:
Диапазон спектра ${{W}_{r}}\left( \omega \right)$ будет, как и в (8), ограничен диапазоном спектра тестового воздействия и полосой пропускания ЛС-системы ${{\Omega }}_{\eta }^{c}$, и точность моделирования в этом диапазоне будет ограничена погрешностью
3. ВЫБОР ПАРАМЕТРОВ РЕГУЛЯРИЗАЦИИ
Критерием оптимальности выбора стабилизаторов $P\left[ {v,p} \right]$ и Q(ω, q) и параметров регуляризации p и q в выражениях (3), (6) и (4), (7) является минимум невязки $\left\| {{{v}_{r}}\left( x \right) - v\left( x \right)} \right\|$ полученного решения и искомого воздействия, но так как искомое воздействие неизвестно, то напрямую применить этот критерий невозможно. Выбор осуществляют в зависимости от интерпретации процессов (1) и (2) и наличия априорной информации о них. Исследованию этого вопроса в частных случаях различной априорной информации посвящен целый ряд работ (см., например, [1–5, 8, 9, 11, 13, 15, 16, 18, 22]). В данном разделе будут описаны общие закономерности и подходы, позволяющие упростить задачу выбора стабилизаторов и параметров регуляризации на практике применительно к рассматриваемому случаю реконструкции сигналов и изображений ЛС-систем.
Реальные сигналы $w\left( x \right)$, $v\left( x \right)$, ${{w}_{c}}\left( x \right)$ и ${{v}_{c}}\left( x \right)$ в ЛС-системах обычно описываются вещественными функциями, а шумы и неопределенности в каждом из сигналов не коррелируют с ними. В этом случае оптимальный стабилизатор в выражениях (4) и (7), обеспечивающий минимум невязки
не зависит от параметров регуляризации $q = \left\{ {{{q}_{i}}} \right\}$ и соответствует так называемому фильтру Винера:
где ${{\left| {V\left( \omega \right)} \right|}^{2}}$ – энергетический спектр восстанавливаемого сигнала $v\left( x \right)$, а ${{\left| {N\left( \omega \right)} \right|}^{2}}$ – энергетический спектр неопределенностей $\eta \left( x \right)$ или ${{\eta }_{c}}\left( x \right)$ рассматриваемых уравнений (1) или (2). Эти спектры для случайных стационарных процессов понимаются в смысле математического ожидания и, как правило, неизвестны, но из того, что они неотрицательны и для вещественных функций симметричны, следует, что стабилизатор надо искать в виде неотрицательной симметричной функции:
где ${{q}_{i}} \geqslant 0$, ${{q}_{n}} > 0$, а $2n$ – порядок стабилизатора. Используя связь между стабилизаторами $P\left[ {{v},p} \right]$ и $Q\left( {\omega ,q} \right)$ (см. комментарии к (4)) и учитывая преобразование Фурье для производных, получим вид соответствующих стабилизаторов в (3) и (6):
которые, следовательно, подразумевают поиск решения, имеющего непрерывные производные порядка $2n$. Значения параметров регуляризации ${{q}_{i}}$ и ${{p}_{i}}$ определяют “степень сглаженности” $2i$-й производной регуляризированного решения, т.е. при ${{q}_{i}},{{p}_{i}} = 0$ сглаженность $2i$-й производной будет равна нулю, а при увеличении этих параметров ее гладкость будет возрастать. В частности, стабилизаторы нулевого порядка ($2n = 0$): $Q\left( {\omega ,q} \right) = {{q}_{0}}$ и $P\left[ {v,p} \right] = {{p}_{0}}{{\left\| v \right\|}^{2}}$ соответствуют поиску решений, “степень непрерывности” которых определяется значениями ${{q}_{0}}$ и ${{p}_{0}}$.
Если все параметры регуляризации и, следовательно, стабилизирующие члены равны нулю: $~p,q = 0$, то будет найдено “несглаженное”, так называемое обобщенное решение Мура–Пенроуза, которое, в общем случае, неустойчиво. В другом крайнем случае, когда значения параметров регуляризации стремятся к бесконечности $p,q \to \infty $, соответствующее регуляризированное решение будет “сглажено” до нуля: $v_{r}^{{p,q}}\left( x \right) \to 0$, что напрямую следует из выражений (4) и (7).
Поскольку высшие i-тые компоненты стабилизаторов более эффективно, чем низшие, сглаживают функции, подавляя высшие частотные гармоники сигналов, в простейшем случае можно ограничиться рассмотрением стабилизатора вида $Q\left( {\omega ,q} \right) = q{{\left| \omega \right|}^{{2z}}}$, где $q,z > 0$, причем $z$, в общем случае, может быть действительным числом, что в стабилизаторе $P\left[ {v,p} \right]$ будет пониматься как производная “нецелого порядка”:
где $F\left( \omega \right)$ – спектр функции $f\left( x \right)$, $I = \sqrt { - 1} $.
Из приведенных рассуждений о возрастании гладкости решения при увеличении значений параметров регуляризации следует простой способ “визуальной” оценки оптимальных значений параметров регуляризации для реальных сигналов, физический спектр которых, как правило, ограничен, а высокочастотные гармоники в основном соответствуют шумам и неопределенностям. В этом случае при некоторых значениях параметров регуляризации $p$ (или $q$), близких к нулю, решение $v_{r}^{p}\left( x \right)$ будет содержать высокочастотные выбросы, напоминающие случайные, на фоне некоторого низкочастотного тренда, соответствующего физическому сигналу. Высокочастотные выбросы можно сгладить, увеличивая значения $p$ в высших компонентах стабилизаторов (или значение $z$ в простейшем представлении $Q\left( {\omega ,q} \right) = q{{\left| \omega \right|}^{{2z}}}$). Оптимальными можно считать значения, при которых амплитуда выбросов значительно уменьшилась, а форма “основного тренда” не изменилась или изменилась незначительно. Практика работы показывает (см. примеры в разд. 6), что даже значительные отклонения параметров регуляризации от оптимального значения, как правило, не приводят к существенному возрастанию погрешности решения.
Более точные оценки параметров регуляризации следуют из монотонности возрастания невязок $\left\| {w - u * v_{r}^{p}} \right\|$ или $\left\| {{{v}_{c}} * w - {{w}_{c}} * v_{r}^{p}} \right\|$ с ростом значений параметров регуляризации $p$ (или $q$). Действительно, невязки $\left\| {w - u * v_{r}^{p}} \right\|$ и $\left\| {{{v}_{c}} * w - {{w}_{c}} * v_{r}^{p}} \right\|$ монотонно возрастают от $0$ при $p,q = 0$ (когда обобщенное решение не расходится) до $\left\| w \right\|$ и $\left\| {{{v}_{c}} * w} \right\|$ при $p,q \to \infty $ соответственно (так как при $p,q \to \infty $: $v_{r}^{p}\left( x \right) \to 0$) и при оптимальных значениях параметров регуляризации будут соответствовать уровню неопределенности соответствующих уравнений: $\left\| {w - u * v_{r}^{p}} \right\| \approx \left\| \eta \right\|$ для (1) и $\left\| {{{v}_{c}} * w - {{w}_{c}} * v_{r}^{p}} \right\| \approx \left\| {{{\eta }_{c}}} \right\|$ для (2). Отсюда также следуют условия на уровень неопределенности уравнений: $\left\| \eta \right\| \leqslant \left\| w \right\|$ в (1) и $\left\| {{{\eta }_{c}}} \right\| \leqslant \left\| {{{v}_{c}} * w} \right\|$ в (2), при которых реконструкция возможна.
Для поиска оптимального решения можно также использовать и другую априорную информацию о процессах (1) или (2) (как отмечалось выше, см., например, [1–5, 8, 9, 11, 13, 15, 16, 18, 22]). В рассматриваемых случаях обработки сигналов ЛС-систем априорно могут быть известны (или предполагаются известными) некоторые характеристики воздействий по амплитуде, длительности, фронтам, энергии, форме и т.п., которые при оптимальных значениях параметров регуляризации должны быть примерно равны соответствующим характеристикам найденного решения. В общем случае, если приводятся в соответствие несколько априорных данных, для поиска оптимального решения задается целевая функция, включающая сравниваемые характеристики с весовыми коэффициентами, зависящим от значимости характеристик.
4. ВЫБОР ТЕСТОВЫХ СИГНАЛОВ
Характеристики тестовых сигналов напрямую влияют на точность реконструкции сигналов ТИ‑методами, что предъявляет особые требования к их выбору.
Во-первых, как отмечалось при выводе ТИ-уравнения (2), характеристики тестовых сигналов ${{v}_{c}}\left( x \right)$ и ${{w}_{c}}\left( x \right)$ должны удовлетворять допустимым диапазонам для аппаратуры рассматриваемой ЛС-системы по амплитуде, частотному, динамическому диапазону и другим параметрам. Эти требования будем называть аппаратным критерием выбора тестовых сигналов.
Во-вторых, тестовые сигналы должны обеспечивать минимально возможный уровень неопределенности ${{\eta }_{c}}\left( x \right)$ ТИ-уравнения (2), что будем называть критерием неопределенности при их выборе. Если представить ${{\eta }_{c}}\left( x \right)$ в виде
то видно, что она будет минимальна, если тестовые сигналы ${{w}_{c}}\left( x \right)$ и ${{v}_{c}}\left( x \right)$ совпадут или будут пропорциональны реконструируемым сигналам $w\left( x \right)$ и $v\left( x \right)$ с точностью до шумов и неопределенностей. Выполнение этого критерия, по существу, означает реконструкцию воздействий подбором тестовых сигналов, что довольно сложно осуществить на практике. Более реальным является уменьшение шумов и неопределенностей в тестовых сигналах, что при экспериментальном определении можно осуществить статистическими методами, повторяя тестовые измерения необходимое количество раз.
В-третьих, диапазон спектра тестового воздействия должен быть не меньше диапазона спектра реконструируемого воздействия или, если последнее неизвестно, не меньше полосы пропускания ЛС-системы. В противном случае, как было показано ранее (см. (5) и (8)), может возрасти погрешность реконструкции из-за ограничений, связанных с уменьшением частотного диапазона реконструкции
С точки зрения частотного критерия наиболее универсальным было бы использование тестового воздействия в виде дельта-функции, которая имеет равномерный спектр во всем пространстве частот, но использование таких сигналов противоречит аппаратному критерию. Метод тестовых испытаний позволяет в качестве тестового воздействия взять любое приближение дельта-функции, удовлетворяющее аппаратному критерию. Но если, например, в одномерном случае взять воздействие в виде прямоугольного импульса с конечной амплитудой $A$ и длительностью $2T$: ${{v}_{c}}\left( t \right) = d\left( t \right) = A$ при $t \in \left( { - T,T} \right)$, то спектр такого воздействия
будет иметь множество нулей на оси частот при $\omega T = \pi k$, где $k$ – целое, не равное нулю, причем первый нуль при $\omega = {\pi \mathord{\left/ {\vphantom {\pi T}} \right. \kern-0em} T}$, ограничивающий диапазон реконструкции ${{\Omega }}_{\eta }^{c}$, будет обратно пропорционален длительности импульса. Расширить диапазон реконструкции можно, сместив нули спектра в область высших частот за полосу пропускания ЛС-системы или за верхнюю границу частотного диапазона реконструируемого сигнала, но при этом потребуется уменьшить длительность тестового импульса и пропорционально увеличить его амплитуду, чтобы не снижать амплитуду $2AT$ спектра сигнала в области низких частот. Во многих случаях требуемое уменьшение длительности и увеличение амплитуды может снова привести к противоречию с аппаратным критерием.
В качестве приближения к дельта-функции можно попытаться использовать сигналы другой формы, например, гауссовы импульсы
спектр которых
хоть и неравномерен, но не обращается в ноль при любых частотах. Проблему, однако, это не решает, поскольку амплитуда спектра гауссова импульса быстро спадает в области высоких частот, и для расширения диапазона реконструкции, как и в случае с прямоугольным импульсом, потребуется уменьшать его длительность (значение $\sigma $) и увеличивать амплитуду, что также может привести к противоречию с аппаратным критерием.
Решение можно найти, если идти не по пути приближения воздействий к дельта-функции, а по пути использования воздействий специальной формы. Рассмотрим, например, тестовые импульсы в виде разнобедренной трапеции c длительностью $2T$, передним фронтом $b$ и задним $B$:
Спектр этого импульса
можно представить в виде
(9)
$\begin{gathered} R\left( \omega \right) = \frac{{I\omega A}}{{\sqrt {2\pi } }}\frac{{\sin \left( {\frac{{\omega B}}{2}} \right)}}{{\frac{{\omega B}}{2}}}{\text{exp}}\left( { - I\omega \left( {T + \frac{B}{2}} \right)} \right) - \\ ~ - \,\,\frac{{I\omega A}}{{\sqrt {2\pi } }}\frac{{\sin \left( {\frac{{\omega b}}{2}} \right)}}{{\frac{{\omega b}}{2}}}{\text{exp}}\left( {I\omega \left( {T + \frac{b}{2}} \right)} \right).~~ \\ \end{gathered} $В случае равных фронтов ($B = b$) этот спектр
как и прямоугольный импульс, будет иметь множество нулей на оси частот при $\omega \left( {T + {b \mathord{\left/ {\vphantom {b 2}} \right. \kern-0em} 2}} \right) = \pi k$ и ${{\omega b} \mathord{\left/ {\vphantom {{\omega b} 2}} \right. \kern-0em} 2} = \pi m$, где $k,m$ – целые числа, не равные нулю. Причем первый нуль при $\omega = {\pi \mathord{\left/ {\vphantom {\pi {\left( {T + {b \mathord{\left/ {\vphantom {b 2}} \right. \kern-0em} 2}} \right)}}} \right. \kern-0em} {\left( {T + {b \mathord{\left/ {\vphantom {b 2}} \right. \kern-0em} 2}} \right)}}$, ограничивающий диапазон спектра, как и в случае с прямоугольным импульсом, будет обратно пропорционален длительности импульса, что потребует для расширения диапазона реконструкции ${{\Omega }}_{\eta }^{c}$ уменьшать длительность импульса и увеличивать его амплитуду и снова может привести к противоречию с аппаратным критерием. Однако при разных фронтах ($B \ne b$) слагаемые в (9) обратятся в ноль при $\omega b = 2\pi k$ и $\omega B = 2\pi m$, где $k,m$ – целые числа, не равные нулю, и если отношение фронтов ${B \mathord{\left/ {\vphantom {B b}} \right. \kern-0em} b}$ – иррациональное число, то спектр $R\left( \omega \right)$ нигде на оси частот не обратится в ноль.
На рис. 1а показаны графики трапецеидальных импульсов $r\left( t \right)$ с передним фронтом $b = $ 1 нс, длительностью полочки $2T = $ 20 нс, длительностях заднего фронта $B = $ 1, 3π, 10π нс. На рис. 1б – графики модуля спектров $\left| {R\left( \nu \right)} \right|$, $\nu = {\omega \mathord{\left/ {\vphantom {\omega {\left( {2\pi } \right)}}} \right. \kern-0em} {\left( {2\pi } \right)}}$ ГГц, этих импульсов в диапазоне от 0 до 100 МГц и в окрестности $\nu = $ 1 ГГц. При $B = b = 1$ нс спектр равнобедренной трапеции имеет множество нулей на оси частот. При увеличении заднего фронта до $B = ~$ 3π нс (иррациональное отношение ${B \mathord{\left/ {\vphantom {B b}} \right. \kern-0em} b}$) низкочастотные проседания модуля спектра (9) уменьшаются и при $B = $ 10π нс сглаживаются и превращаются в затухающие с частотой осцилляции на фоне основного тренда, определяемого первым слагаемым ${{{\text{sin}}\left( {{{\omega b} \mathord{\left/ {\vphantom {{\omega b} 2}} \right. \kern-0em} 2}} \right)} \mathord{\left/ {\vphantom {{{\text{sin}}\left( {{{\omega b} \mathord{\left/ {\vphantom {{\omega b} 2}} \right. \kern-0em} 2}} \right)} {\left( {{{\omega b} \mathord{\left/ {\vphantom {{\omega b} 2}} \right. \kern-0em} 2}} \right)}}} \right. \kern-0em} {\left( {{{\omega b} \mathord{\left/ {\vphantom {{\omega b} 2}} \right. \kern-0em} 2}} \right)}}$ в (9), нули спектра которого зависят от величины переднего фронта, а не от длительности импульса. В увеличенном масштабе показаны спектры в окрестности первого нуля спектра первого слагаемого в (9) ($\nu = $ 1 ГГц). При иррациональном соотношении фронтов ${B \mathord{\left/ {\vphantom {B b}} \right. \kern-0em} b} = $ 3π, 10π второе слагаемое отлично от нуля и модуль спектра не равен нулю $\left| {R\left( \nu \right)} \right| \ne 0$.
Рис. 1.
Трапецеидальные импульсы $r\left( t \right)$ с разной длительностью заднего фронта 1 (1), $3\pi $ (2) и $10\pi $ нс (3) (а) и графики модуля их спектров $\left| {R\left( \nu \right)} \right|$ в диапазоне частот от 0 до 100 МГц и в окрестности 1 ГГц (на вставке) (б).

Несмотря на то, что при иррациональном соотношении фронтов ${B \mathord{\left/ {\vphantom {B b}} \right. \kern-0em} b}$ модуль спектра нигде не обращается в ноль, при небольших соотношениях фронтов его просадка в области низких частот (~30–50 МГц на рис. 1) может быть значительной, что приведет к сокращению диапазона реконструкции. На практике для уменьшения просадок удобно полагать $B,T \gg b$ и $B \geqslant 3T$ (соотношение длительности заднего фронта и полочки больше, чем 3 : 2: $B \geqslant 3T$).
При воздействии в форме разнобедренной трапеции диапазон реконструкции будет определяться не длительностью импульса, влияющего на амплитуду низкочастотного спектра, а величиной его переднего фронта. Уменьшая длительность переднего фронта разнобедренного трапецеидального импульса, можно расширить диапазон реконструкции вплоть до полосы пропускания ЛС-системы, а увеличивая длительность импульса, – достичь необходимой амплитуды спектра в области низких частот и при этом не вступить в противоречие с аппаратным критерием.
Максимально увеличить амплитуду спектра в области низких частот можно, устремляя задний фронт импульса к бесконечности $B \to \infty $. В пределе получим тестовый сигнал в виде нефинитной справа ступенчатой функции с передним фронтом $b$, спектр которой, согласно (9), будет равен
Уменьшая фронт ступеньки, можно сдвинуть первый нуль спектра ступенчатой функции, определяющий диапазон реконструкции, за полосу пропускания ЛС-системы, при этом амплитуда спектра в области низких частот будет максимально возможной. Отметим, забегая вперед, что использовать в качестве тестового воздействия импульс в виде ступенчатой функции, как и другие нефинитные справа сигналы, можно только при реконструкции с помощью прямых методов регуляризации.
К выводам, аналогичным анализу разнобедренной трапеции, можно прийти, если рассмотреть воздействие в виде импульса с разными экспоненциальными фронтами, такими как сигналы вида
где ${{a}_{0}} \gg a > 0$, но эти импульсы, как и ступенька, не являются финитными справа, что приводит к определенным трудностям при использовании спектральных методов реконструкции. Кроме того, могут возникнуть сложности с генерацией экспоненциальных импульсов при проведении тестовых испытаний экспериментальным путем.
5. ЧИСЛЕННЫЕ РАСЧЕТЫ
В большинстве практических случаев реальные сигналы крайне трудно выразить через аналитические функции, что определяет важность численных расчетов при математической реконструкции сигналов. При проведении численных расчетов сигналы представляют или сразу получают в виде массивов их значений в некоторых дискретных точках интересуемой области. В одномерном случае (многомерный случай будет рассмотрен позже) функции $w\left( x \right)$, $v\left( x \right)$, $u\left( x \right)$ из (1), определенные на вещественной оси $x$, представим в виде множеств ${\mathbf{w}} = \left\{ {{{w}_{i}}} \right\}$, ${\mathbf{v}} = \left\{ {{{v}_{i}}} \right\}$, ${\mathbf{u}} = \left\{ {{{u}_{i}}} \right\}$ их значений ${{w}_{i}} = w\left( {{{x}_{i}}} \right)$, ${{v}_{i}} = v\left( {{{x}_{i}}} \right)$, ${{u}_{i}} = u\left( {{{x}_{i}}} \right)$ в точках ${{x}_{i}} = i\Delta $, где $i$ – целое число, а Δ – шаг разбиения оси $x$.
Согласно (1) значение ${{w}_{k}} = w\left( {{{x}_{k}}} \right)$ в точке ${{x}_{k}} = k\Delta $ можно представить в виде интеграла:
Частичная сумма этого интеграла на введенном разбиении равна
что совпадает с выражением дискретной свертки, где дискретная характеристика системы переопределена как ${{u}_{i}} = u\left( {i\Delta } \right)\Delta $:
(10)
${{w}_{k}} \approx \mathop \sum \limits_i {{u}_{{k - i}}}{{v}_{i}} = \mathop \sum \limits_i {{u}_{i}}{{v}_{{k - i}}}.$Дискретное уравнение свертки (10) является частичной суммой интеграла в (1) и обеспечивает лишь приближенное соответствие между значениями функций ${{w}_{i}} = w\left( {{{x}_{i}}} \right)$, ${{v}_{i}} = v~\left( {{{x}_{i}}} \right)$, ${{u}_{i}} = u\left( {{{x}_{i}}} \right)$ в точках разбиения. Таким образом, уравнение (10) изначально содержит ошибки при дискретизации сигналов, и, в общем случае, задача реконструкции сигналов из него некорректно поставленной. Ошибки при дискретизации можно уменьшить, уменьшая шаг разбиения, т.к. при $\Delta \to 0$ частичная сумма сходится к значению интеграла
Чтобы представить аппаратную функцию из класса обобщенных в дискретном виде, необходимо уменьшать шаг разбиения Δ до нуля, при этом количество дискретных значений сигналов будет стремиться к бесконечности. В ТИ-методах для реконструкции сигналов используется ТИ-уравнение (2), которое не содержит функций из класса обобщенных, а содержит лишь сигналы, относящиеся к классу основных функций, которые вполне можно с необходимой точностью (например, используя теорему Котельникова) представить в дискретном виде.
Далее, как и раньше, ТИ-уравнения (2) в дискретном виде и расчетные выражения для реконструкции дискретных сигналов из них получим с помощью переходных замен из известных выражений для ЛС-уравнений (1), функции в которых относятся к классу основных.
Интерпретируя массивы ${\mathbf{u}}$, ${\mathbf{v}}$ и ${\mathbf{w}}$ как вектора и вводя матрицы Тёплица
с порождающими векторами ${\mathbf{u}}$ и ${\mathbf{v}}$, дискретное уравнение типа свертки (10) можно представить в матричной форме:
Если функции $v\left( x \right)$ и $u\left( x \right)$ из (1) определены на конечных интервалах $\left( {a,c} \right)$ и $\left( {b,d} \right)$, то функция $w\left( x \right)$ будет определена на интервале $\left( {a + b,\,\,c + d} \right)$, и при однородном разбиении оси $x$ размеры $n,m,k$ векторов ${\mathbf{w}}$, $~{\mathbf{v}}$ и ${\mathbf{u}}$ будут конечны и связаны соотношением $n = m + k - 1$. Матрицы ${\mathbf{U}}$ и ${\mathbf{V}}$ в этом случае являются прямоугольными матрицами размером $n \times m$ и $n \times k$ соответственно, а уравнения (11) представляют собой переопределенную прямоугольную систему линейных алгебраических уравнений, содержащую все актуальные значения дискретных функций. Регуляризированное дискретное решение ${{{\mathbf{v}}}_{r}}$ такой системы, соответствующее (3), со стабилизатором первого порядка в дискретном виде дается выражением
(12)
${{{\mathbf{v}}}_{r}} = {{\left[ {{\mathbf{U}}{\text{*}}{\mathbf{U}} + {{p}_{0}}{\mathbf{I}} + {{p}_{1}}{\mathbf{D}}_{1}^{2}} \right]}^{{ - 1}}}{\mathbf{U}}{\text{*}}{\mathbf{w}}.~$Здесь ${\mathbf{U}}{\text{*}}$ – сопряженно-транспонированная матрица, ${\mathbf{I}}$ и ${{\left( \cdot \right)}^{{ - 1}}}$ – единичная и обратная $m \times m$ матрицы, ${{p}_{0}}$ и ${{p}_{1}}$ – неотрицательные параметры регуляризации. Матрица ${\mathbf{D}}_{1}^{2} = {\mathbf{D}}_{1}^{*}{{{\mathbf{D}}}_{1}}$, где ${{{\mathbf{D}}}_{1}} = \left\{ {{{d}_{{ij}}}} \right\}$ – $m \times m$ матрица, соответствует оператору первых конечных разностей с ненулевыми элементами ${{d}_{{ii}}} = 1$ и ${{d}_{{i + 1,i}}} = - 1$. Согласно идеологии методов регуляризации при определенных значениях параметров регуляризации дискретное решение системы (12) существует и устойчиво к вариации исходных данных, а обратная матрица в (12) является симметричной, положительно определенной и достаточно хорошо обусловленной.
Прямоугольные переопределенные системы (11) можно преобразовать в квадратные, дополнив векторы ${\mathbf{v}}$ и ${\mathbf{u}}$ нулями до размера $n$ вектора ${\mathbf{w}}$. Применяя к таким векторам ${\mathbf{w}}$, ${\mathbf{v}}$ и ${\mathbf{u}}$ дискретное $n$-точечное прямое и обратное преобразование Фурье (ДПФ)
($i,j = 0, \ldots ,n - 1$, $I = \sqrt { - 1} $), из (4) для стабилизатора первого порядка получим дискретный регуляризированный спектр ${{V}_{r}} = \left\{ {V_{i}^{r}} \right\}$, $i = 0, \ldots ,n - 1$, восстанавливаемого воздействия ${\mathbf{v}}$:
(13)
$V_{i}^{r} = \frac{{U_{i}^{*}{{W}_{i}}}}{{U_{i}^{*}{{U}_{i}} + {{q}_{0}} + {{q}_{1}}{{i}^{2}}}}~.$Здесь ${{q}_{0}}$ и ${{q}_{1}}$ – неотрицательные параметры регуляризации.
Нефинитные справа функции ${v}\left( x \right)$, $u\left( x \right)$ и $w\left( x \right)$ можно свести к финитным, если для них существуют конечные базовые интервалы: $\left( {a,c} \right)$, $\left( {b,d} \right)$ и $\left( {a + b,c + d} \right)$, $a,b,c,d \ne \infty $, вне которых значениями функций можно пренебречь, т.е. если
Погрешности при сведении функций к финитным можно интерпретировать как дополнительные неточности в исходных данных, что позволяет использовать прямые (12) и спектральные (13) методы регуляризации для решения некорректно поставленных задач реконструкции.
Прямые методы регуляризации (12) можно также применять в случае финитных слева и нефинитных справа функций $u\left( x \right)$, $v\left( x \right)$, $w\left( x \right)$, определенных на интервалах $\left( {a,\infty } \right)$, $\left( {b,\infty } \right)$ и $\left( {a + b,\infty } \right)$, для которых базовых интервалов не существует. Векторы, соответствующие таким функциям, имеют бесконечные размеры, но для конечных векторов ${\mathbf{u}} = \left\{ {{{u}_{i}}} \right\}$, ${\mathbf{v}} = \left\{ {{{{v}}_{i}}} \right\}$ и ${\mathbf{w}} = \left\{ {{{w}_{i}}} \right\}$, $i = 0, \ldots ,n - 1$, содержащих любое количество первых $n$ элементов, можно получить аналогичную (11) конечную систему $n$ линейных алгебраических уравнений и, используя выражения (12), восстановить элементы вектора ${\mathbf{v}} = \left\{ {{{v}_{i}}} \right\}$. Отметим, что использование спектральных методов (13) для решения таких систем приведет к неверным результатам, так как сигналы будут рассматриваться как периодические на данных интервалах, что может не соответствовать условию задачи.
В методе тестовых испытаний вместо уравнения типа свертки (1) рассматривается основное уравнение тестовых испытаний (2), регуляризированное решение которого можно найти с помощью замен функций в выражениях (3) и (4). Аналогично, подставляя в уравнение (11) вместо матрицы Теплица ${\mathbf{U}}$ матрицу ${{{\mathbf{W}}}_{c}}$ с производящим вектором ${{{\mathbf{w}}}_{c}}$ и вместо вектора ${\mathbf{w}}$ – вектор ${{{\mathbf{V}}}_{c}}{\mathbf{w}}$, получим аналогичное (2) основное дискретное уравнение тестовых испытаний:
Соответствующие подстановки в (12) дают решение (14) прямыми методами регуляризации:
(15)
${{{\mathbf{v}}}_{r}} = {{\left[ {{\mathbf{W}}_{c}^{*}{{{\mathbf{W}}}_{c}} + {{p}_{0}}{\mathbf{I}} + {{p}_{1}}{\mathbf{D}}_{1}^{2}} \right]}^{{ - 1}}}{\mathbf{W}}_{c}^{*}{{{\mathbf{V}}}_{c}}{\mathbf{w}},~$а подстановки дискретных спектров соответствующих функций в (13) (при возможности использования интегральных преобразований) дают выражение для дискретного спектра ${{V}_{r}} = \left\{ {V_{i}^{r}} \right\}$ регуляризированного решения ${{{\mathbf{v}}}_{r}}$:
(16)
$V_{i}^{r} = \frac{{W_{i}^{{c*}}V_{i}^{c}{{W}_{i}}}}{{W_{i}^{{c*}}W_{i}^{c} + {{q}_{0}} + {{q}_{1}}{{i}^{2}}}}.$Спектры в (16), как и в (13), получены из векторов, дополненных нулями до размера, одинакового с вектором максимального размера. Оптимальные значения параметров регуляризации ищутся способами, аналогичными описанным ранее для недискретных функций.
Дискретное ТИ-уравнение (14) можно также использовать для моделирования отклика ${\mathbf{w}}$ ЛС-систем на заданное воздействие ${\mathbf{v}}$. В этом случае, заменяя в (12) ${\mathbf{w}}$ на $~{{{\mathbf{W}}}_{c}}{\mathbf{v}}$ и ${\mathbf{U}}$ на ${{{\mathbf{V}}}_{c}}$, получим регуляризированное решение ${{{\mathbf{w}}}_{r}}$ для моделируемого отклика ${\mathbf{w}}$:
(17)
${{{\mathbf{w}}}_{r}} = {{\left[ {{\mathbf{V}}_{c}^{*}{{{\mathbf{V}}}_{c}} + {{p}_{0}}I + {{p}_{1}}{\mathbf{D}}_{1}^{2}} \right]}^{{ - 1}}}{\mathbf{V}}_{c}^{*}{{{\mathbf{W}}}_{c}}{\mathbf{v}}.~~$Аналогичные замены спектров в (12) дают регуляризированный спектр $~{{W}_{r}} = \left\{ {W_{i}^{r}} \right\}$ моделируемого отклика $~{\mathbf{w}}$:
(18)
$W_{i}^{r} = \frac{{V_{i}^{{c*}}W_{i}^{c}{{V}_{i}}}}{{V_{i}^{{c*}}V_{i}^{c} + {{q}_{0}} + {{q}_{1}}{{i}^{2}}}}.$При обработке данных большого объема $n$ прямые методы регуляризации (12), (15) и (17) потребуют выполнения порядка $O(5{{n}^{3}} + 2n)$ операций с матрицами и векторами, а спектральные (13), (16) и (18) с использованием быстрых прямого и обратного преобразований Фурье – порядка $O\left( {7n~{\text{lg}}\left( n \right) + n} \right)$ операций.
Однако спектральные методы можно использовать для реконструкции сигналов далеко не всегда и не на всех интервалах. Дело в том, что при численных расчетах функции рассматриваются на конечных интервалах и по умолчанию предполагается, что при использовании прямых методов регуляризации функции равны нулю вне этих интервалов, а при использовании спектральных – что они периодически продолжены вовне. Если функции финитны и рассматриваемые интервалы включают всю область определения функций, то решения, полученные прямыми и спектральными методами, совпадут. Если функции нефинитны справа, но рассматриваемые интервалы включают их базовые интервалы, то решения, полученные прямыми и спектральными методами, будут примерно равны. Если же функции нефинитны справа и не имеют базовых интервалов, то решения, полученные спектральными методами, будут неверны. Поскольку прямые методы позволяют проводить реконструкцию для любых финитных слева функций (т.е. сигналов, имеющих начало), то с этой точки зрения прямые методы регуляризации более универсальны, чем спектральные.
С физической точки зрения спектральные методы дают установившееся решение для периодических на рассматриваемых интервалах сигналов, а прямые – неустановившееся решение для сигналов, равных нулю вне рассматриваемых интервалов. Если, например, переходные процессы в ЛС-системе долго не затухают, то короткое воздействие приведет к отклику, длительность которого на много порядков может превышать длительность воздействия. При реконструкции прямыми методами регуляризации можно ограничиться рассмотрением интервала, включающего лишь интервал короткого воздействия или даже интересующую его часть, в то время как при использовании спектральных методов необходимо рассматривать весь интервал существования переходных процессов или, по крайней мере, базовый интервал, за пределами которого этими процессами можно пренебречь.
В случае нефинитных справа функций спектральные методы неприменимы. Если функции определены на конечных или имеют базовые интервалы, но интересующий интервал реконструкции меньше базового в $k$ раз, то для реконструкции прямыми методами потребуется $O\left( {{{4{{n}^{3}}} \mathord{\left/ {\vphantom {{4{{n}^{3}}} {{{k}^{3}}~}}} \right. \kern-0em} {{{k}^{3}}~}} + {n \mathord{\left/ {\vphantom {n k}} \right. \kern-0em} k}} \right)$ операций, а для спектральных – $O\left( {7n{\kern 1pt} ~{\text{lg}}\left( n \right) + n} \right)$ операций, и при некоторых $k$ до определенных значений $n$ прямые методы могут быть быстрее спектральных. Например, при $k = {{10}^{3}}$ прямые методы будут быстрее спектральных до $n \approx 1.5 \times {{10}^{5}}$ (${n \mathord{\left/ {\vphantom {n k}} \right. \kern-0em} k} \approx 150$), при $k = {{10}^{4}}$ – до $n \approx 5.2 \times {{10}^{6}}$ (${n \mathord{\left/ {\vphantom {n k}} \right. \kern-0em} k} \approx 520$), а при $k = {{10}^{5}}$ – до $n \approx 18.2 \times {{10}^{7}}$ (${n \mathord{\left/ {\vphantom {n k}} \right. \kern-0em} k} \approx 1820$).
6. РЕКОНСТРУКЦИЯ ОДНОМЕРНЫХ СИГНАЛОВ
Рассмотрим численный эксперимент реконструкции одномерных сигналов ЛС-систем методом тестовых испытаний на примере емкостного делителя напряжения. Принципиальная электрическая схема делителя и сигналы на его входе ${{v}_{1}}\left( t \right)$, ${{v}_{2}}\left( t \right)$, ${{v}_{c}}\left( t \right)$ и выходе ${{w}_{1}}\left( t \right)$, ${{w}_{2}}\left( t \right)$, ${{w}_{c}}\left( t \right)$ представлены на рис. 2. Делитель состоит из емкостей ${{C}_{1}}$ = 100 мкФ, ${{C}_{2}}$ = 100 нФ, сопротивлений ${{R}_{1}}$ = = 50 МОм, ${{R}_{2}}$ = 50 кОм и предназначен для уменьшения амплитуды высоковольтных входных сигналов ($v\left( t \right)$, ${{v}_{c}}\left( t \right)$), поступающих от источника $P$ с внутренним сопротивлением ${{R}_{p}}$ = 10 Ом, до низковольтных сигналов ($w\left( t \right)$, ${{w}_{c}}\left( t \right)$) на выходе устройства. Сигналы передаются по помехозащищенным линиям передачи ${{T}_{1}}$ и ${{T}_{2}}$ (витая пара) с волновым сопротивлением 150 Ом и задержками распространения сигнала 5 и 10 нс соответственно. Делитель представляет собой ЛС-систему, так как является составной системой (см. разд. 1) из пассивных элементов, ${\text{и}}$ сигналы в нем могут быть реконструированы ТИ-методами.
Рис. 2.
Принципиальная электрическая схема емкостного делителя напряжения (а) и графики сигналов на его входе ${{v}_{1}}$, ${{v}_{2}}$, ${{v}_{c}}$ и выходе ${{w}_{1}}$, ${{w}_{2}}$, ${{w}_{c}}$ (б).

Численный расчет выходных ${{w}_{1}}\left( t \right)$, ${{w}_{2}}\left( t \right)$, ${{w}_{c}}\left( t \right)$ и входных ${{v}_{1}}\left( t \right)$, ${{v}_{2}}\left( t \right)$, ${{v}_{c}}\left( t \right)$ сигналов делителя проводился с помощью симулятора аналоговых/цифровых схем Micro-Cap 12.2.0.4 от Spectrum Software. Принципиальная электрическая схема делителя была взята из стандартной библиотеки этой программы. Для удобства сравнения сигналов отклики ${{w}_{1}}\left( t \right)$, ${{w}_{2}}\left( t \right)$, ${{w}_{c}}\left( t \right)$ показаны с постоянным масштабным коэффициентом, компенсирующим падение напряжения на делителе. Частота среза делителя составляет ≈1 МГц, что значительно меньше диапазонов спектров входных сигналов, поэтому в процессе обработки высокочастотные спектры сигналов искажаются и форма выходных сигналов существенно отличается от формы входных.
Восстановлению ТИ-методами подлежат входные трапецеидальный ${{v}_{1}}\left( t \right)$ и гауссов ${{v}_{2}}\left( t \right)$ сигналы, показанные на рис. 2. Исходными данными для их реконструкции служат выходные сигналы ${{w}_{1}}\left( t \right)$ и ${{w}_{2}}\left( t \right)$ и тестовые сигналы ${{w}_{c}}\left( t \right)$ и ${{v}_{c}}\left( t \right)$. Согласно описанной ранее процедуре выбора тестовых сигналов в качестве тестового воздействия был использован финитный импульс в виде разнобедренной трапеции с передним фронтом 1 нс, длительностью полочки 20 нс и задним фронтом 10π нс. Реконструкцию проводили прямыми ТИ‑методами, так как выбранный интервал реконструкции от 0 до 200 нс не является базовым для рассматриваемых сигналов и применять спектральные методы в этом случае нельзя (см. комментарии в разд. 5). При реконструкции использован стабилизатор нулевого порядка (с ${{p}_{1}}$ = 0 в (15)). Результаты реконструкции ${\mathbf{v}}_{1}^{r}$ и ${\mathbf{v}}_{2}^{r}$ входных сигналов ${{{\mathbf{v}}}_{1}}$ и ${{{\mathbf{v}}}_{2}}$ из выходных ${{w}_{1}}$ и ${{w}_{2}}$ в точках дискретизации с шагом $\Delta t$ = 1.0 нс представлены на рис. 3а и 3б.
Рис. 3.
Реконструкция ${\mathbf{v}}_{1}^{r}$ входного трапецеидального сигнала ${{{\mathbf{v}}}_{1}}$ из выходного сигнала ${{{\mathbf{w}}}_{1}}$ (a) и реконструкция ${\mathbf{v}}_{2}^{r}$ входного гауссова сигнала ${{{\mathbf{v}}}_{2}}$ из выходного сигнала ${{{\mathbf{w}}}_{2}}$ (б) прямыми численными ТИ-методами.

Несмотря на то, что дополнительные случайные шумы в исходные сигналы на рис. 3 не добавлялись, но решаемая задача реконструкции является некорректно поставленной, так как исходные данные содержат ошибки дискретизации (см. комментарии к (10)), а также ошибки вычислений и округлений при представлении данных. Погрешность восстановления трапецеидального импульса ${\mathbf{v}}_{1}^{r}$ (рис. 3а) в точках выборки составила
а гауссова импульса ${\mathbf{v}}_{2}^{r}$ (рис. 3б) –
При добавлении случайных шумов в исходные выходные и тестовые сигналы увеличивается общая неопределенность ${{{\mathbf{\eta }}}_{c}}$ ТИ-уравнений (2) и (14), что влияет на погрешность реконструкции. На рис. 4а представлен пример восстановления ${\mathbf{v}}_{1}^{r}$ трапецеидального входного сигнала ${{{\mathbf{v}}}_{1}}$ из выходных данных ${{{\mathbf{w}}}_{1}}$ со случайными шумами $\Delta {{{\mathbf{w}}}_{1}}$ ($\delta {{w}_{1}} = {{\left\| {\Delta {{{\mathbf{w}}}_{1}}} \right\|~} \mathord{\left/ {\vphantom {{\left\| {\Delta {{{\mathbf{w}}}_{1}}} \right\|~} {\left\| {{{{\mathbf{w}}}_{1}}} \right\|~}}} \right. \kern-0em} {\left\| {{{{\mathbf{w}}}_{1}}} \right\|~}}$ = 4.5%), а на рис. 4б – пример восстановления ${\mathbf{v}}_{2}^{r}$ гауссова входного сигнала ${{{\mathbf{v}}}_{2}}$ из выходных данных ${{{\mathbf{w}}}_{2}}$ со случайными шумами $\Delta {{{\mathbf{w}}}_{2}}$ ($\delta {{w}_{2}} = {{~\left\| {\Delta {{{\mathbf{w}}}_{2}}} \right\|~} \mathord{\left/ {\vphantom {{~\left\| {\Delta {{{\mathbf{w}}}_{2}}} \right\|~} {\left\| {{{{\mathbf{w}}}_{2}}} \right\|}}} \right. \kern-0em} {\left\| {{{{\mathbf{w}}}_{2}}} \right\|}}$ = 4.9%). В этих примерах шумы в тестовых данных были равны $\delta {{w}_{c}} = {{\left\| {\Delta {{{\mathbf{w}}}_{c}}} \right\|~} \mathord{\left/ {\vphantom {{\left\| {\Delta {{{\mathbf{w}}}_{c}}} \right\|~} {\left\| {{{{\mathbf{w}}}_{c}}~} \right\|}}} \right. \kern-0em} {\left\| {{{{\mathbf{w}}}_{c}}~} \right\|}}$ = 2.9% и $\delta {{v}_{c}} = {{\left\| {\Delta {{{\mathbf{v}}}_{c}}} \right\|~} \mathord{\left/ {\vphantom {{\left\| {\Delta {{{\mathbf{v}}}_{c}}} \right\|~} {\left\| {{{{\mathbf{v}}}_{c}}} \right\|~}}} \right. \kern-0em} {\left\| {{{{\mathbf{v}}}_{c}}} \right\|~}}$ = 1.5%. Погрешность восстановления трапецеидального и гауссова входных сигналов составила соответственно
Рис. 4.
Восстановление ${\mathbf{v}}_{1}^{r}$ трапецеидального входного сигнала ${{{\mathbf{v}}}_{1}}$ из выходных данных ${{{\mathbf{w}}}_{1}}$ с шумами ${{\Delta }}{{{\mathbf{w}}}_{1}}$ (а) и восстановление ${\mathbf{v}}_{2}^{r}$ гауссова входного сигнала ${{{\mathbf{v}}}_{2}}$ из выходных данных ${{{\mathbf{w}}}_{2}}$ с шумами ${{\Delta }}{{{\mathbf{w}}}_{2}}$ (б) прямыми численными ТИ-методами.

что даже меньше уровня неопределенностей в исходном выходном сигнале ${{{\mathbf{w}}}_{2}}$ ($\delta {{w}_{2}}$ = 4.9%) из-за сглаживания шумов при регуляризации.
На рис. 3 и 4 приведены результаты реконструкции сигналов при оптимальных значениях параметра регуляризации ${{p}_{0}}$ в (15), соответствующие минимуму невязки $\left\| {{{{\mathbf{v}}}_{r}}{\text{\;}} - {\mathbf{v}}} \right\|$ реконструированного ${{{\mathbf{v}}}_{r}}$ и истинного ${\mathbf{v}}$ воздействия на интервале реконструкции $T$ от 0 до 150 нс. Но поскольку истинное воздействие изначально неизвестно, на практике применяются описанные ранее методы выбора параметров регуляризации. На рис. 5а представлены результаты реконструкции $v_{1}^{r}$, $v_{1}^{{r - }}$ и $v_{1}^{{r + }}$ трапецеидального воздействия ${{v}_{1}}$ при оптимальном значении параметра регуляризации ${{p}_{0}} = 0.155$ ($v_{1}^{r}$) и неоптимальных – на два порядка меньше и больше оптимального: $0.01{{p}_{0}}$ ($v_{1}^{{r - }}$) и $100~{{p}_{0}}$ ($v_{1}^{{r + }}$) соответственно. При значениях параметра регуляризации меньше оптимального в реконструкции ($v_{1}^{{r - }}$) на всем интервале даже до начала воздействия присутствуют высокочастотные выбросы, напоминающие случайные. С увеличением параметра регуляризации амплитуда выбросов начинает уменьшаться и при превышении оптимальных значений начинает сглаживаться основной импульс, что приводит к появлению сигнала перед началом воздействия. На практике оптимальным значением параметра регуляризации можно считать значения, при которых амплитуда высокочастотных выбросов в реконструируемом сигнале минимальна и форма “основного тренда” не искажена. Погрешности при таком “визуальном” определении даже при больших отклонениях от оптимальных значений, как правило, незначительны. Так, при оптимальном значении параметра регуляризации ${{p}_{0}}$ (см. рис. 5а) погрешность реконструкции сигнала $v_{1}^{r}$ составляет $\delta v_{1}^{r} = {{\left\| {{\mathbf{v}}_{1}^{r} - {{{\mathbf{v}}}_{1}}~} \right\|} \mathord{\left/ {\vphantom {{\left\| {{\mathbf{v}}_{1}^{r} - {{{\mathbf{v}}}_{1}}~} \right\|} {\left\| {{{{\mathbf{v}}}_{1}}} \right\|}}} \right. \kern-0em} {\left\| {{{{\mathbf{v}}}_{1}}} \right\|}} = $ = 6.7%. При двукратном уменьшении значения параметра регуляризации погрешность реконструкции составила 7.3%, а при трехкратном – 8.2%. При двукратном увеличении параметра регуляризации погрешность реконструкции составила 7.0%, а при трехкратном – 7.3%.
Рис. 5.
Результаты реконструкции ${\mathbf{v}}_{1}^{{r - }}$, ${\mathbf{v}}_{1}^{r}$ и ${\mathbf{v}}_{1}^{{r + }}$ трапецеидального воздействия ${{v}_{1}}$ при оптимальном ${{p}_{0}}$ (${\mathbf{v}}_{1}^{r}$) и неоптимальных $0.01{{p}_{0}}$ (${\mathbf{v}}_{1}^{{r - }}$) и $100{{p}_{0}}$ (${\mathbf{v}}_{1}^{{r + }}$) значениях параметра регуляризации (а) и графики погрешности $\delta {{v}_{r}} = {{\left\| {{{{\mathbf{v}}}_{r}} - {\mathbf{v}}} \right\|} \mathord{\left/ {\vphantom {{\left\| {{{{\mathbf{v}}}_{r}} - {\mathbf{v}}} \right\|} {\left\| {\mathbf{v}} \right\|}}} \right. \kern-0em} {\left\| {\mathbf{v}} \right\|}}$ реконструкции ${{{\mathbf{v}}}_{r}}$ воздействия ${\mathbf{v}}$ из отклика ${\mathbf{w}}$ при различных уровнях шумов $\delta w = {{\left\| {\Delta {\mathbf{w}}} \right\|} \mathord{\left/ {\vphantom {{\left\| {\Delta {\mathbf{w}}} \right\|} {\left\| {\mathbf{w}} \right\|}}} \right. \kern-0em} {\left\| {\mathbf{w}} \right\|}}$ в нем, наиболее ожидаемых значений $\langle \delta {{v}_{r}}\rangle $ и эмпирических пределов $\delta {{v}_{{{\text{min}}}}}$ и $\delta {{v}_{{{\text{max}}}}}$ погрешности реконструкции (б).

При определении оптимального значения параметров регуляризации также можно использовать различные априорные сведения об ЛС-системе, обрабатываемых сигналах и шумах в них. Например, может быть известна форма части импульса воздействия, скажем, длительность переднего фронта или длительность полочки реконструируемого трапецеидального сигнала ${{v}_{1}}\left( t \right)$. В этих случаях при определении параметра регуляризации ${{p}_{0}}$ можно использовать критерий минимума невязки ${{\left\| {{{{\mathbf{v}}}_{r}}{\text{\;}} - {\mathbf{v}}} \right\|}_{{{{T}_{K}}}}}$ на соответствующей части ${{T}_{K}}$ интервала реконструкции $T$. В рассмотренных примерах реконструкции трапецеидального сигнала такие ограничения интервала невязки ${{T}_{K}}$ от 0 до 30 нс при сравнении переднего фронта приведут к увеличению погрешности $\delta v_{1}^{r}$ на всем интервале реконструкции с 2.2 до 2.4% на рис. 3а и с 6.687 до 6.693% на рис. 4а, а при сравнении длительности при ${{T}_{K}}$ от 30 до 80 нс – на рис. 3а с 2.2 до 2.3% и на рис. 4а с 6.7 до 7.5%.
Погрешность реконструкции зависит не только от уровня, но и от распределения шумов и неопределенностей в исходных данных, что вызывает значительный разброс погрешности реконструкции при одинаковом уровне шумов и неопределенностей. На рис. 5б показаны значения погрешности $\delta {{v}_{r}} = {{\left\| {{{{\mathbf{v}}}_{r}}{\text{\;\;}} - {\mathbf{v}}} \right\|} \mathord{\left/ {\vphantom {{\left\| {{{{\mathbf{v}}}_{r}}{\text{\;\;}} - {\mathbf{v}}} \right\|} {\left\| {\mathbf{v}} \right\|}}} \right. \kern-0em} {\left\| {\mathbf{v}} \right\|}}$ реконструкции ${{{\mathbf{v}}}_{r}}$ воздействия ${\mathbf{v}}$ в разных численных экспериментах при различных уровнях шумов $\delta w = {{\left\| {\Delta {\mathbf{w}}} \right\|} \mathord{\left/ {\vphantom {{\left\| {\Delta {\mathbf{w}}} \right\|} {\left\| {\mathbf{w}} \right\|}}} \right. \kern-0em} {\left\| {\mathbf{w}} \right\|}}$ в отклике ${\mathbf{w}}$ и фиксированном уровне неопределенностей в тестовых данных
Для наглядности показаны наиболее ожидаемые значения $\langle \delta {{v}_{r}}\rangle $ и примерные эмпирические пределы $\delta {{v}_{{{\text{min}}}}}$ и $\delta {{v}_{{{\text{max}}}}}$ погрешности реконструкции. Замедление роста погрешности реконструкции при увеличении уровня неопределенности ТИ-уравнения объясняется тем, что при возрастании амплитуды шума увеличиваются высокочастотные составляющие его спектра, которые более эффективно сглаживаются при регуляризации.
Программа реконструкции сигналов на основе выражения (15) со стабилизатором нулевого порядка (при ${{p}_{1}} = 0$) была написана на языке Maple 2017.3 Waterloo Maple Inc. и оформлена в виде подпрограммы, входными данными для которой являются векторы ${\mathbf{w}}$, ${{{\mathbf{w}}}_{c}}$, ${{{\mathbf{v}}}_{c}}$ и параметр регуляризации ${{p}_{0}}$. Время работы подпрограммы при реконструкции сигналов на рис. 3 и 4 (размерность векторов и матриц равна 500) при расчетах на компьютере HP 255 G7 с процессором AMD Ryzen 3, 2.5 ГГц в среде Maple 2017 составило 0.2…0.4 с.
7. ВОССТАНОВЛЕНИЕ МНОГОМЕРНЫХ ДАННЫХ И ИЗОБРАЖЕНИЙ
Многомерные данные и изображения обычно представляют собой массивы дискретных значений обрабатываемых сигналов в некоторых точках дискретизации в рассматриваемых областях. Для больших, но конечных объемов данных предпочтительны методы с использованием интегральных преобразований (см. оценки в разд. 5). В случае реконструкции ${{{\mathbf{V}}}_{r}} = \left\{ {v_{{ij \ldots }}^{r}} \right\}$ многомерного массива данных ${\mathbf{V}} = \left\{ {{{v}_{{ij \ldots }}}} \right\}$ из массива искаженных данных ${\mathbf{W}} = \left\{ {{{w}_{{ij \ldots }}}} \right\}$ с использованием тестовых данных ${{{\mathbf{V}}}_{c}} = \left\{ {v_{{ij \ldots }}^{{\text{c}}}} \right\}$ и ${{{\mathbf{W}}}_{c}} = \left\{ {w_{{ij \ldots }}^{{\text{c}}}} \right\}$ уравнение (15) со стабилизатором первого порядка имеет вид
(19)
$V_{{ij \ldots }}^{r} = \frac{{W_{{ij \ldots }}^{{c*}}V_{{ij \ldots }}^{c}{{W}_{{ij \ldots }}}}}{{W_{{ij \ldots }}^{{c*}}W_{{ij \ldots }}^{c} + {{q}_{0}} + {{q}_{1}}\left( {{{i}^{2}} + {{j}^{2}} + \ldots } \right)}},$где
– многомерные спектры массивов ${{{\mathbf{V}}}_{r}}$, ${\mathbf{V}}$, ${\mathbf{W}}$, ${{{\mathbf{V}}}_{c}}$ и ${{{\mathbf{W}}}_{c}}$. Многомерные спектры исходных данных можно получить с помощью многомерного n-точечного ДПФ. При этом, чтобы достичь минимальной погрешности восстановления, исходные массивы необходимо преобразовать в массивы размера $n \times n \times \cdots $, где $n$ – наибольший размер для всех массивов, заполняя нулями добавленные элементы. Например, если массив ${\mathbf{W}}$ имеет размер $k \times n \times m \times \cdots $ и $n$ является наибольшим числом для всех массивов, то все массивы, включая ${\mathbf{W}}$, должны быть преобразованы в массивы размера $n \times n \times \cdots $ с нулевыми добавленными элементами.
На рис. 6 представлен пример восстановления ${{{\mathbf{V}}}_{r}}$ цветного изображения ${\mathbf{V}}$ из размытого изображения ${\mathbf{W}}$ в формате jpeg по формуле (19). Ядро искажения ${\mathbf{U}}$ было выбрано случайным образом в виде асимметричной прямоугольной матрицы размером 50 × 70 и полагалось одинаковым для всех трех цветовых слоев, поэтому при реконструкции достаточно было использовать тестовые изображения ${{{\mathbf{W}}}_{c}}$ и ${{{\mathbf{V}}}_{c}}$ в градациях серого.
Рис. 6.
Восстановление ${{{\mathbf{V}}}_{r}}$ изображения ${\mathbf{V}}$ из изображения ${\mathbf{W}}$ при его размытии с ядром ${\mathbf{U}}$ спектральными ТИ-методами с использованием тестовых изображений ${{{\mathbf{V}}}_{c}}$ и ${{{\mathbf{W}}}_{c}}$.

В общем случае, если искажения в цветовых слоях изображений различны, тестовые изображения должны быть в цветовом формате реконструируемого изображения. Каждый цветовой слой восстанавливаемого изображения ${\mathbf{V}}$ составлял 500 × 700 пикселей, а слои искаженного и увеличенного за счет размытия изображения ${\mathbf{W}}$ – 549 × 769 пикселей. Тестовое изображение ${{{\mathbf{V}}}_{c}}$ имело размер 201 × 201 пикселей, а размытое тестовое изображение ${{{\mathbf{W}}}_{c}}$ – 250 × 270 пикселей.
Согласно алгоритму восстановления (16) все изображения ${\mathbf{W}}$, $~{{{\mathbf{W}}}_{c}}$ и ${{{\mathbf{V}}}_{c}}$ сначала преобразовывались в изображения максимального размера 769 × 769 пикселей с заполнением добавленных элементов нулями. Затем с помощью двумерного n‑точечного ДПФ определялись их многомерные спектры $W$, $~{{W}_{c}}$ и ${{V}_{c}}$ размером 769 × 769 и был найден спектр регуляризированного решения (16). Изображение ${{{\mathbf{V}}}_{r}}$ было восстановлено из спектра ${{V}_{r}}$ с помощью обратного n-точечного ДПФ. При реконструкции использовался стабилизатор нулевого порядка. Исходный размер (500 × 700 пикселей) восстановленного изображения определяли по разнице размеров тестовых изображений. Погрешность реконструкции (без дополнительных шумов в исходных данных) составила
Точность реконструкции сильно зависит от используемого тестового изображения (см. разд. 4). Так, при использовании тестового изображения в виде шахматной доски размером 201 × 201 пикселя погрешность реконструкции составит $\delta {{{\mathbf{V}}}_{r}}$ = 0.3%, а при размере 80 × 80 – $\delta {{{\mathbf{V}}}_{r}}$ = 0.6%. Погрешность реконструкции изображений оказалась намного меньше, чем при реконструкции одномерных сигналов, что связано с использованием ядра искажений ${\mathbf{U}}$ из класса основных функций, а не обобщенных, как в примере, приведенном в разд. 6.
Программа реконструкции изображений на основе выражения (19) со стабилизатором нулевого порядка (${{p}_{1}} = 0$) была написана на языке Maple 2017.3 Waterloo Maple Inc. в виде подпрограммы, входными данными для которой являются изображения ${\mathbf{W}}$, $~{{{\mathbf{W}}}_{c}}$, ${{{\mathbf{V}}}_{c}}$ и параметр регуляризации ${{p}_{0}}$. Согласно оценкам, представленным в разд. 5, число операций $N\left( n \right)$ при реконструкции изображения размером $n \times n$ составляет порядка $N\left( n \right)\sim O\left( {14{{n}^{2}}~{\text{lg}}\left( n \right) + {{n}^{2}}} \right)$. Время реконструкции цветного (три слоя) jpeg-изображения 500 × 700 пикселей (максимальный размер $n = 769$), представленного на рис. 6, при расчетах на компьютере HP 255 G7 с процессором AMD Ryzen 3, 2.5 ГГц в среде Maple 2017 составило 127.1 с, что соответствует времени расчета каждого из слоев $t = {{127.1} \mathord{\left/ {\vphantom {{127.1} 3}} \right. \kern-0em} 3} = 42.4$ с. Время реконструкции другого jpeg-изображения в градациях серого (имеющего один слой) размером 201 × 201 пикселей (максимальный размер ${{n}_{2}} = 270$) составило ${{t}_{2}} = 4.453$ с. Пересчет времени расчета большего изображения из данных для меньшего дает близкое значение
что свидетельствует о хорошей точности оценок числа операций при реконструкции.
Аналогичным образом можно обрабатывать массивы данных любой размерности. Как и в одномерном случае, метод тестовых испытаний позволяет провести реконструкцию без детального анализа каналов и устройств получения и передачи информации.
ЗАКЛЮЧЕНИЕ
В работе предложены и проанализированы математические методы реконструкции внешних воздействий из откликов на эти воздействия линейных стационарных систем обработки информации с помощью проведения тестовых испытаний. Использование данных тестовых испытаний позволяет осуществить не слепую реконструкцию сигналов без определения функций их искажений или аппаратных функций систем обработки. Реконструкция сигналов осуществляется из уравнения тестовых испытаний, которое является уравнением типа свертки и описывает работу ЛС-систем, но при этом не содержит аппаратных функций, которые, в общем случае, могут принадлежать к классу обобщенных, а только сигналы из класса основных функций. Такой подход значительно упрощает задачу реконструкции и позволяет использовать многие существующие методы решения некорректно поставленных и плохо обусловленных задач для восстановления реальных недетерминированных сигналов.
Рассмотрены особенности применения в методах тестовых испытаний техники регуляризации при решении некорректно поставленных задач реконструкции. Проанализированы критерии выбора тестовых сигналов, погрешности реконструкции. Представлены результаты численных экспериментов по восстановлению одномерных сигналов и двумерных изображений при разных уровнях зашумления.
Предложенные математические методы восстановления сигналов позволяют повысить точность, разрешающую способность и расширить область применения ЛС-систем, что является доступной альтернативой аппаратным методам, связанным с трудноразрешимыми научно-техническими задачами и дорогостоящими технологическими решениями. Рассмотренные математические методы восстановления, в отличие от аппаратных, не ограничены физическими возможностями реализации устройств, что обусловливает широкие возможности их практического применения.
Список литературы
Mueller J.L., Siltanen S. Linear and Nonlinear Inverse Problems with Practical Applications. Philadelphia: SIAM, 2012. V. 10.
Hansen P.C. Discrete Inverse Problems: Insight and Algorithms. Fundamentals of Algorithms. Philadelphia: SIAM, 2010.
Kaipio J., Somersalo E. Statistical and Computational Inverse Problems. N.Y.: Springer, 2010.
Chen K. Matrix Preconditioning Techniques and Applications. Cambridge: Cambridge Univ. Press, 2005.
Tarantola A. Inverse Problem Theory and Methods for Model Parameter Estimation. Philadelphia: SIAM, 2005.
Vito E.D., Rosasco L., Caponnetto A. et al. // J. Machine Learning Research. 2005. P. 883.
Ben-Israel A., Greville T.N.E. Generalized Inverses: Theory and Applications. 2nd ed. N.Y.: Springer, 2003. https://doi.org/10.1007/b97366
Zhang X., Burger M., Bresson X., Osher S. // SIAM J. Imaging Sci. 2010. V. 3. № 3. P. 253. https://doi.org/10.1137/090746379
Afonso M.V., Bioucas-Dias J.M., Figueiredo M.A.T. // IEEE Trans. 2011. V. IP-20. № 3. P. 681.
Gonzalez C.R., Woods R.E. Digital Image Processing. 3rd ed. Upper Saddle River: Pearson Prentice Hall, 2008.
Василенко Г.И., Тараторин А.М. Восстановление изображений. М.: Радио и связь, 1986.
Benning M., Burger M. // 2018. arXiv:1801.09922v1 [math.NA].
Kazufumi I., Bangti J. Inverse Problems: Tikhonov Theory and Algorithms. Singapore: World Scientific., 2014.
Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Численные методы решения некорректных задач. М.: Наука, 1990.
Тихонов А.Н., Гончарский А.В., Степанов В.В., Ягола А.Г. Регуляризирующие алгоритмы и априорная информация. М.: Наука, 1983.
Иванов В.К., Васин В.В., Танана В.П. Теория линейных некорректных задач и её приложения. М.: Наука, 1978.
Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. М.: Наука, 1974.
Матысик О.В. Итерационная регуляризация некорректных задач. Saarbrücken: LAP LAMBERT Acad. Publ., 2015.
Самарский А.А., Вабищевич П.Н. Численные методы решения обратных задач математической физики. М.: УРСС, 2004.
Vogel C.R. Computational Methods for Inverse Problems. Philadelphia: SIAM, 2002.
Gilyazov S.F., Goldman N.L. Regularization of Ill-posed Problems by Iteration Methods. Dordrecht: Kluwer Acad. Publ., 2000.
Шлома А.М. // Журн. вычисл. математики и матем. физики. 1996. Т. 36. № 3. С. 15.
Vaidyanathan P.P., Chen T. // IEEE Trans. 1995. V. SP-43. № 5. P. 1090.
Денисов А.М. Введение в теорию обратных задач. М.: МГУ, 1994.
Бакушинский А.Б., Гончарский А.В. Итеративные методы решения некорректных задач. М.: Наука, 1989.
Chen S.S., Donoho D.L., Saunders M.A. // SIAM Rev. 2001. V. 43. № 1. P. 129.
Percival D.B., Walden A.T. Wavelet Methods for Time Series Analysis. Cambridge: Cambridge Univ. Press, 2000.
Mallat S. A Wavelet Tour of Signal Processing. 2nd ed. San Diego: Academic Press, 1999.
Domínguez A. // IEEE Pulse. 2015. V. 6. № 1. P. 38.
von zur Gathen J., Gerhard J. Modern Computer Algebra. 3-rd ed. Cambridge: Cambridge Univ. Press, 2013.
Crutchfield S. The Joy of Convolution. Web Applet. N.Y.: Johns Hopkins Univ., 2010. http:// www.jhu.edu/signals/convolve/.
Hespanha J.P. Linear System Theory. Princeton: Princeton Univ. Press, 2009.
Phillips C.I., Parr J.M., Riskin E.A. Signals, Systems and Transforms. Upper Saddle River: Prentice Hall, 2007.
Баскаков С.И. Радиотехнические цепи и сигналы. М.: Высш. школа, 2005.
Uludag A.M. // J. Mathematical Analysis and Appl. 1998. V. 227. № 2. P. 335. https://doi.org/10.1006/jmaa.1998.6091
Sobolev V.I. Convolution of Functions, Encyclopedia of Mathematics. Helsinki: EMS Press, 2001.
Напалков В.В. Уравнения свертки в многомерных пространствах. М.: Наука, 1982.
Колмогоров А.Н., Фомин С.В. Элементы теории функций и функционального анализа. М.: Наука, 2004. 7-е изд.
Владимиров В.С. Уравнения математической физики. М.: Наука, 1981.
Прасолов В.В. Задачи и теоремы линейной алгебры. М.: Наука, 1996.
Дополнительные материалы отсутствуют.
Инструменты
Радиотехника и электроника







