§ 23. Итерационный метод восстановления фазы

Из изложенного в § 22 следует один очень важный вывод: если может быть восстановлен приемлемый вариант фазы то он должен соответствовать правильной форме изображения . В данном параграфе мы изложим наиболее эффективный в настоящее время метод вычисления фазы Так же как и в случае многих других интересных приложений математики в физике, пока неизвестно, как строго «доказать», что данный метод «правилен». Но «наука» основана на эксперименте, а эксперименты могут проводиться и в вычислительной лаборатории, не только в физической. Если установлено, то некоторая вычислительная процедура дает приемлемые результаты при определенных условиях, а также при различных других условиях, то приходится сделать вывод, что она и вообще верна, раз нет данных, опровергающих этот вывод. Конечно, теоретическое обоснование нужно искать, но неудачи в поисках не должны обескураживать.

Первым этапом итерационного восстановления фазы по заданным первичным отсчетам модуля спектра является оценка размеров кадра Напомним (см. § 12), что фурье-образ величины, которая задана только в регулярно расположенных точках растра, является периодической функцией. С учетом этого по теореме об автокорреляции [формула (7.10)] получим

где кадр в плоскости изображения, в пределах которого повторяется оценка автокорреляции (это обозначение согласуется с обозначением, введенным в последнем абзаце § 10). Выберем (на основе вычисленного уровня зашумленности данных) такое пороговое значение с, что спектр можно будет считать пренебрежимо малым при

Рассчитав спектр на основе первичных отсчетов модуля спектра по его форме можно определить, действительно ли первичные отсчеты удовлетворяют неравенствам (20.12). Дело в том, что спектр должен удовлетворять условию (23.2) на всей внешней границе кадра Конечно, нужно надеяться, что величины несколько больше величин следовательно, условие (23.2) будет выполняться всюду в пределах кадра овне некоего прямоугольного кадра, значительно меньшего, чем кадр . Этот меньший кадр является наилучшей имеющейся оценкой кадра Назовем предварительно обработанной автокорреляцией неотрицательную часть спектра существующую в пределах кадра [все отрицательные части спектра нужно удалить как связанные с каким-либо видом шума, поскольку изображение по определению имеет неотрицательные значения].

Размеры кадра обязательно вдвое меньше размеров кадра в силу теоремы об ограничении протяженности (8.15). Таким образом, величины можно вычислить непосредственно по первичным отсчетам

Кратко остановимся на случае, когда на некоторых частях внешней границы кадра т. е. не выполняются условия (20.12). В этом случае нужна дальнейшая обработка спектра Его необходимо умножить на функцию окна (см. § 15), которая равна единице на большей части внутренней области кадра Убывает до нуля (с нулевым наклоном) на внешней границе. Идея в том, чтобы сделать изображение, соответствующее кадру 0 у, как можно более близким к автокорреляции, в то же время в минимальной степени искажая Следовательно, далее можно воспользоваться методами § 15. Неотрицательная часть этого взвешенного варианта принимается за так что кадры и должны считаться одинаковыми. Хотя изображение видимо, существенно отличается от автокорреляции это, вероятно, самая лучшая оценка автокорреляции истинного изображения.

Обозначим через предварительно обработанный спектр интенсивностей

Кроме того, введем обозначения

где произвольные целые числа, а вариант изображения восстановимый по (эти обозначения согласуются с обозначениями § 14). Величина «фаза в формуле (23.5) фактически представляет собой величину, которая может быть восстановлена по Величины будем называть модифицированными первичными отсчетами величин При подходе к задаче восстановления фазы, описываемом в данном параграфе, вычисляются величины называемые здесь модифицированными первичными отсчетами фазы. Величины это теперь фактически исходные данные рассматриваемой фазовой задачи в частотной плоскости.

Для нахождения точных значений приходится обращаться к итерационным методам. При этом можно в более или менее прямой форме использовать такие алгоритмы, как алгоритм Герхберга (§ 11) или любой из алгоритмов максимальной энтропии (§ 15). Причины эффективности этих алгоритмов пока не ясны, а потому при выборе того или иного алгоритма можно руководствоваться лишь своим опытом. Но какой бы алгоритм не использовался, необходимы начальные оценки значений Оценки первичных отсчетов величины определяются следующим образом:

Способы получения величин рассматриваются ниже. Важно, что функция оценки первичных отсчетов будет отлична от нуля всюду на кадре если не считать маловероятного случая, когда значения у составляют множество значений фаз, точно соответствующее форме изображения

В том случае, если форма изображения является «сосредоточенной», т. е. содержит отдельную четко ограниченную деталь, имеется несколько итерационных процедур восстановления фазы, которые дают удовлетворительные результаты. Существенно, что спектры таких изображений характеризуются «сильной интерференцией»,

т. е. большая часть частотной плоскости, где модуль спектра заметно отличен от нуля, покрыта пересекающимися системами резких полос. Это напоминает случаи элементарного и простого голографическою подходов со смешенным опорным участком в которых из-за широкого разделяющего промежутка между опорным участком и неизвестной частью изображения в спектре возникают резко выраженные частые полосы.

Очень трудный случай восстановления фазы когда форма изображения содержит неяркие детали, наложенные на сравнительно яркий и более или менее однородный фон. В этом случае в частотной плоскости имеет место лишь слабая интерференция. Спектры, как правило, состоят из больших центральных лепестков и небольшой боковой «ряби». Чтобы применить алгоритмы восстановления фазы в таких ситуациях, вначале нужно выполнить процедуру устранения вуали. Сначала определяют «край» центрального лепестка [т. е. замкнутую кривую в частотной плоскости, окружающую центральный лепесток, которую видит глаз там, где последний переходит в остальную часть модуля спектра либо визуально (что обычно наиболее эффективно благодари очень высокой способности человеческого мозга распознавать образы), либо с использованием пороговых алгоритмов. Далее путем устранения всякой ряби вокруг центрального лепестка модуля спектра получают его сглаженный вариант Вариант величины лишенный вуали, вычисляется как

где коэффициент к выбирается так, чтобы значение примерно в 1,5 раза превышало наибольший соседний максимум величины Это дает вместо почти безынтерференционной картины спектральное распределение, обнаруживающее сильную интерференцию в том смысле, что высота «гребней» и глубина «впадин» модуля спектра почти не изменяется в той части частотной плоскости, где величина значительна.

Подчеркнем, что процедура «устранения вуали» не имеет детального теоретического обоснования. Но она делает свое дело! К тому же результаты мало зависят от выбранного значения коэффициента к. Сглаженная функция вводится для того, чтобы сохранить в любую «рябь», имеющуюся в центральном лепестке Для упрощения обозначений дадим другое определение величины

Далее восстановимым изображением аппроксимируется

упомянутая выше неяркая деталь. Форма изображения вуали реконструируется из описываемым ниже методом реконструкции из После этого неяркая деталь, т. е. , может быть наложена на иуаль.

Следующий этап состоит в расчете величин

где коэффициенты определяются формулами (22.14). Укажем, что величины, первоначально введенные в формулах (22.19) и (22.20), были переопределены по двум причинам: чтобы не изобретать новых обозначений и потому, что это пригодится в дальнейшем. Величины (23.10) следует рассматривать теперь как фактические и промежуточные отсчеты модуля а не Все это нужно для того, чтобы можно было получить величины с помощью простого рекурсивного алгоритма, описанного в § 22 [текст, относящийся к формулам (22.21)-(22.30)]. Прежде чем выполнять этот алгоритм, нужно выбрать (конечное) значение для неотрицательного целого числа введенного в § 22 [см. начало абзаца, содержащего формулу (22.6)]. Оказывается, что совершенно не требуется большое число Как показывает вычислительная практика, обычно лучше всего выбирать

При фактические отсчеты величины идентичны ее модифицированным первичным отсчетам, что, конечно, удобно в вычислительном отношении.

Величины первоначально введенные в формуле (22.19), определяются теперь как оценки фазы, даваемые после замены на рекурсивным алгоритмом, описанным в § 22. При мы имеем, очевидно,

Такой метод получения начальных оценок фазы назовем грубой оценкой фазы. Вычислительная практика показывает, что от применения или неприменения этого метода может полностью зависеть успешность итерационного подхода, описываемого ниже, в том случае, когда исходные данные обработки выводятся из сильно зашумленных данных измерения.

Наиболее эффективной при восстановлении фаз двумерных спектров итерационной процедурой явилась модификация алгоритма Герхберга (§ 11), носяшая название алгоритма Файнапа. Далее описываются два варианта такого алгоритма.

Основой алгоритма Файнапа является алгоритм коррекции ошибок. Для сто реализации необходимо знать начальную функцию фазы

или начальное изображение. Сам Файнап использовал псевдослучайные начальные значения. Преимущество псевдослучайного начального изображения состоит в том, что его кадр может быть взят одного размера с кадром размеры которого определяются по виду функции визуально, как говорилось ранее в данном параграфе. Наш опыт показывает, что можно исходить и из некой начальной функции фазы, которую мы обозначим здесь через Ниже при описании алгоритма коррекции ошибок предполагается, что начальной величиной является функция фазы.

Псевдослучайная функция очень хорошо подходит для обработки данных, синтезированных и необработанных данных измерения в случае сосредоточенного изображения с низким уровнем зашумленности. В общем же случае более эффективной представляется грубая оценка фазы. Кроме того, при использовании псевдослучайной функции фазы итерационный процесс иногда оказывается расходящимся (или сходящимся столь медленно, что фактически оказывается бесполезным). Какой бы ни был вид взятой функции первая оценка изображения вычисляется по формуле

Мы хотим, чтобы оценка изображения существовала на кадре равном кадру но мы не можем сказать априори, где поместить его на плоскости изображения. Мы оцениваем положение центра калра как положение точки в которой максимальная величина

Затем положительная часть изображения в пределах кадра -принимается за текущую оценку изображения Таким образом, мы ограничили оценку изображения переменным кадром. Кадр называется переменным, поскольку значения а и варьируются до достижения максимального значения величины Далее на основе фурье-образа изображения строится новая оценка фазовой функции

Подставив новую оценку (23.14) снова в формулу (23.12), получаем другой вариант оценки чем устанавливается итерационный цикл. Отметим, что, вообше говоря, положение кадра иное на кажлой итерации.

Чтобы оборвать итерационный процесс, необходимо иметь объективный критерий обрывания. Мы выбираем для этого некое

пороговое значение с, как и пороговое значение в условии (23.2), исходя из оценки уровня зашумленности. После каждой итерации мы берем заданные фактические амплитуды отсчетов и последние оценки их фаз и рассчитываем соответствующие величины для промежуточных отсчетов методом интерполяции функциями (см. § 10 и 11), записывая эти рассчитанные величины как а х для некоторого множества значений индексов промежуточных отсчетов. Итерационный процесс обрывается, когда выполняется условие

Напомним, что фактические амплитуды отсчетов никогда не изменяются, поскольку это исходные данные, а потому их незачем включать в условие (23.15), тем более что их включение привело бы к ненужным затратам машинного времени.

Хотя алгоритм коррекции ошибок идеально подходит для ранних этапов процесса итераций, он обычно оказывается непригодным для улучшения оценки до окончательной формы изображения, которая достигается, само собой разумеется, уже после того, как выполняется критерий обрывания итерационного процесса (23.15).

В гибридном алгоритме ввода-вывода предыдущая оценка восстановимого изображения играет роль исходной оценки («драйвера»): текущая оценка функции получается путем определенного изменения. Чтобы показать, в чем состоит эта модификация, мы обозначим через множество точек, принадлежащих кадру которые лежат за пределами кадра а через множества точек (принадлежащих кадру , в которых функция имеет соответственно второму индексу неотрицательные и отрицательные значения. Текущая оценка изображения вычисляется по формуле

где неотрицательная константа, называемая коэффициентом усиления.

В конце итерационного процесса используется фиксированный и увеличенный кадр занимающий область причем , обычно примерно на 10% больше величин

Эти итерации выгодно проводить с использованием циклов Файнапа. Каждый такой цикл состоит из итераций коррекции ошибок и итераций гибридного алгоритма ввода-вывода; суммарное число итераций как правило, берется равным 50, но, конечно, может быть значительно большим или меньшим. Значения параметров

выбираются эмпирически (см. литературу, указанную в вводной части главы), так же как решается вопрос о том, с каких итераций (коррекции ошибок или гибридного алгоритма ввода-вывода) должен начинаться каждый цикл Файнапа.