Хермен Г. Восстановление изображений по проекциям: Основы реконструктивной томографии

Подождите немного. Документ загружается.

208

ГЛАВА II

ж) если X*** > 2, то х^

к + Х)

располагается по другую сторону от гипер-

плоскости H

и на большем расстоянии от нее, нежели оценка х

{к

Для дальнейшего анализа релаксационного метода для неравенств вновь

обратимся к рис. 11.1, на котором представлены два неравенства (т.е.

Р = 2), а также векторы п

и п

и скаляры q

{

n q

, определенные выраже-

ниями (11.12) и (11.13) соответственно. Предположим, что для всех к пара-

метр Х

(

*> = 1. Необходимо также выбрать вектор начальных значений, ко-

торый мы примем равным

(0)

•й

(1114)

тогда легко проверить, что

*'"-й - *'*-[!]

Поскольку оценка лг

(2)

принадлежит одновременно двум гиперплоскостям

A/j и N

, все остальные величины х

(к)

при к > 2 будут одинаковыми и рав-

ными х

К Следовательно, метод обеспечивает сходимость для х* = х

лежащих на гиперплоскости N.

Сходимость, имевшая место в предыдущем примере, называется конеч-

ной сходимостью, поскольку оценка х

(Аг)

остается постоянной после первых

нескольких итераций. Однако это не означает, что релаксационный метод

для неравенств обеспечивает конечную сходимость при выборе величины

(

*> в соответствии с неравенством (11.10).

Теперь перейдем к изучению систем уравнений. Предположим, что зада-

ны У-мерные векторы a

и вещественные числа Ь

{

(1 ^ i ^ /). Пусть также

^. = {*!<",, х> = />,), (11.16)

L= f]L

. (11.17)

Можно подметить, что L также можно представить в виде пересечения

множества полуплоскостей. В самом деле, если положить Р = 2/ и для лю-

бого /, лежащего в интервале 1 < / < /, определить величины

^2,

1 = M<-</„.Y> < -Ь\

(11.18)

то

= /V

i г N

(11.20)

АЛГЕБРАИЧЕСКИЕ АЛГОРИТМЫ РЕКОНСТРУКЦИИ ИЗОБРАЖЕНИЙ 209

и аналогично тому, как N определяется формулой (11.7), так и величина N

определяется (11.17) при условии, что N, задано соотношениями (11.18) и

(11.19).

Следовательно, релаксационный метод можно применить к системе

неравенств для нахождения элементов множества L. Здесь мы предполага-

ем,

ЧТО <<!,., Я,> > 0 при 1 ^ / < /.

Заметим, однако, что

= //i„ (11.21)

где И

— граничная гиперплоскость полуплоскости N

как это определено

равенством (11.11). Следовательно, суммарный эффект (2к

—

1)-й и 2/г-й

итераций состоит в перемещении (если оно вообще имеется) в направлении,

перпендикулярном гиперплоскости L

. Эти две итерации можно объеди-

нить в одну операцию и получить следующий алгоритм для релаксационно-

го метода решения систем уравнений:

(0)

произвольное

\-*

'= х-*» + c

lk)

. (11.22)

Из аналчза результатов, полученных при использовании релаксационно-

го метода для неравенств, легко видеть, что если для всех значений к имеем

а \1

к)

удовлетворяет соотношению (11.10), то релаксационный метод для

равенств дает последовательность оценок дг

(0)

, дг

л*

... , которая сходит-

ся к вектору множества L при условии, что L непустое.

Алгоритм, описываемый выражениями (11.1) и (11.2), является частным

случаем релаксационного метода для равенств при Х

(А:)

= 1 и любом к.

Этот случай иллюстрируется рис. 11.1. Полагая, что а

}

= л

= и

^1 = Я\ и Ь

= Я

определяются значениями (11.12) и (11.13), можно пока-

зать,

что если в качестве начального значения взято х^ согласно (11.14), то

мы получим оценки дг

(1)

и xS

\ соответствующие значениям (11.15). С этой

точки зрения последовательность, получаемая с помощью релаксационного

метода для равенств, отличается от последовательности, полученной ранее

обсуждавшимся релаксационным методом для неравенств. Это обусловле-

но тем, что оценка х

(2)

не принадлежит множеству L

, и поэтому в даль-

нейшем необходимо брать значения ближе и ближе к элементу множества

который в данном случае определяется единственной точкой пересечения

* двух линий L

{

и L

. С помощью данной геометрической модели можно

прийти к выводу, что последовательность векторов получается путем опу-

скания перпендикуляров либо на линию L

либо на L

(рис. 11.1).

В общем случае множество L может содержать более чем один эле-

мент. При несколько более тщательном выборе значения \

(())

можно полу-

210

ГЛАВА tl

чить гарантию того, что релаксационный метод для равенств будет давать

сходимость именно к тому элементу множества L, который удовлетворяет

критерию оптимальности, а именно критерию минимума нормы (разд. 6.4).

Чтобы сделать это, введем множество 5 векторов, представляющее собой

множество всевозможных линейных комбинаций из a

, т. е.

S \х\х = У Pjdj для некоторых вещественных значений

/3,.

| (11.24)

Последующее утверждение иногда называют теоремой о минимальной

норме. Если L — непустое множество, то в L Г\ S существует один и толь-

ко один его элемент х*; кроме того, для всех элементов х множества L,

отличных от х*, справедливо неравенство

Их*|| < ||х||. (11.25)

Согласно соотношению (6.37), ИдгИ

= <дг. дг>.

При этом х* — элемент с минимальной нормой на множестве L.

В разд. 6.4 уже обсуждался вопрос о том, почему выбор элемента с мини-

мальной нормой считается полезным при реконструкции изображений.

Из условий (11.22) видно, что если лг

* выбрано элементом множества

5, то и величина **** также будет элементом множества 5 при всех к. Это

приводит к тому, что, согласно основным положениям линейно'й алгебры,

предел х* последовательности оценок дг

\ х

\ дг

(2)

... также принадлежит

5. Поскольку предельное значение х* также принадлежит множеству L

(благодаря сходимости релаксационного метода для равенств), то и само

значение х* лежит в L П 5. Следовательно, как следует из теоремы для

минимума нормы, значение х* является элементом с минимальной нормой

в L.

В примере, представленном на рис. 11.1, S представляет собой множе-

ство всевозможных двумерных векторов. Следовательно, х

(0)

принадлежит

множеству 5, однако это значение уже выбрано. Поскольку в данном прос-

том примере имеется лишь одно решение, необходимо, чтобы это было

именно решение с минимальной нормой. Однако приведенный пример не-

достаточно типичен. В том случае, когда число неизвестных превышает

число уравнений, число решений будет, конечно, больше и задача будет со-

стоять

выборе значения

дг*

для получения решения с минимальной нормой.

Имеется ряд разновидностей релаксационного метода, которые обеспе-

чивают получение решения с минимальной нормой для системы нера-

венств, однако последние более сложны по своей сущности, чем упомяну-

тые выше методы, и поэтому в данной книге не рассматриваются. В следу-

ющем разделе мы приступим к алгоритму, реализация которого аналогич-

на релаксационному методу нахождения решения с минимальной нормой

для системы неравенств.

АЛГЕБРАИЧЕСКИЕ АЛГОРИТМЫ РЕКОНСТРУКЦИИ ИЗОБРАЖЕНИЙ 211

11.3.

АДДИТИВНЫЙ АЛГЕБРАИЧЕСКИЙ АЛГОРИТМ

РЕКОНСТРУКЦИИ

В данном разделе будут рассмотрены вопросы применения релаксацион-

ных методов к реконструкции изображения. Указанные методы называют

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

процессе единственной итерации текущая оценка обновляется путем добав-

ления к ней скалярной величины, кратной транспонированной строке мат-

рицы проекций [соотношения (11.1) и

(11.2)].

Наиболее простая реализация алгоритма предполагает использование

релаксационного метода решения системы уравнений, описываемого выра-

жениями (11.22) и (11.23) при a

= г, и b

= у

Основываясь на результатах,

полученных в предыдущем разделе, мы получаем последовательность век-

торов х®\ х^

\ х^

\ ... , которые сходятся к величине оценки дг*, так что

Rx* = у, при условии что она существует. Проблема состоит в том, что

из-за существующего соотношения между вектором изображения х и векто-

ром измеренных данных у [выражение (6.24)] оценка дг* может не сущест-

вовать, а даже если такая

дг*

и существует, то не пригодна для реконструк-

ции в дискретном варианте. С этой точки зрения вызывает приятное удив-

ление тот факт, что даже такой простой подход дает приемлемое качество

реконструкции, особенно в тех случаях, когда параметр релаксации выбран

достаточно малым (например, 0,05), что иллюстрируется примерами

разд.

11.5.

Один из способов применить теорию, полученную в предыдущем разде-

ле,

к решению проблемы реконструкции изображений состоит в использо-

вании методов решения системы неравенств [выражение

(6.42)].

Большое

число исследований было проведено именно в этом направлении; в частно-

сти,

были разработаны близкие к алгебраическим алгоритмам реконструк-

ции процедуры нахождения решения системы неравенств с минимальной

нормой. Указанные процедуры требуют более сложного математического

описания, чем то, что было дано в рамках соотношения (11.1). В дополне-

ние к последовательности У-мерных векторов дг

(0)

, х

(1)

, дг

(2)

... в этом случае

получают и используют последовательность 7-мерных векторов и®\ i/'\

и<

>... .

В данной книге будет рассмотрен иной подход, имеющий, однако, сход-

ное назначение, а именно аддитивный алгебраический алгоритм рекон-

струкции путем нахождения байесовской оценки (разд. 6.4) при некоторых

ограничивающих предположениях.

На основе разд. 6.4 мы сделаем следующие предположения: пусть вели-

чины X и Е — многомерные гауссовские случайные функции с нулевым

средним $i

и с дисперсиями У

пропорциональными единичным мат-

рицам соответствующей размерности. Другими словами, мы предполага-

ем,

что компоненты выборок X

—

не коррелированы между собой, а

каждая из них является выборкой одной и той же гауссовской случайной

функции. Предполагается также, что компоненты выборки Е некоррелиро-

ваны и каждая из них является выборкой одной и той же гауссовской слу-

чайной функции с нулевым средним значением.

212

ГЛАВА II

Величины t

мы используем для обозначения диагональных элементов

матрицы V

, as

— для обозначения диагональных значений матрицы У

Положим г = t/s. В соответствии с выражением (6.27) байесовская оценка

представляет собой вектор

дг,

который минимизирует выражение вида

\\у - Ях||

+ ||х-М

. (11.26)

Заметим, что малые значения г свидетельствуют о том, что априорная ин-

формация о наблюденных значениях вектора изображения играет более

важную роль, нежели измеренные данные, тогда как при больших г имеет

место обратная ситуация.

Теперь мы приступим к следующим операциям. Рассмотрим уравнение

Rx + e = ycI + J неизвестными, а именно с неизвестными компонентами

векторов х и е. Таким образом, образуется совместная система уравнений,

поскольку для любого х имеется решение уравнения е = у

—

Rx. Хотя к

упомянутой системе уравнений можно применить обычные методы нахож-

дения значений дг, минимизирующих выражение (11.26), однако в данном

случае потребуется немного более сложный подход.

[;]■

Обозначим вектор-столбцы размерности / + / величиной I, где и

имеет / компонент, a z имеет / компонент. Мы также используем обозна-

чение \Е rR] для матрицы размерностью / х (/ + /), первые / столбцов

которой образуют единичную матрицу Е размерностью / ж /, а последние

/ столбцов — матрицу R, каждый элемент которой кратен г.

Система уравнений вида

(:Н-

[£ гЯ] =г(у-ад, (11.27)

является совместной, поскольку, рассматривая произвольный /мерный

вектор £ и полагая

й = r(y - Rn

- R2), (11.28)

делаем вывод, что решение I I удовлетворяет системе (11.27)

>ai_

Причина введения системы уравнений (11.27) состоит в том, что если

м* и z* — векторы, для которых I I является решением системы урав-

нений (11.27) с минимальной нормой, и если также

г* = **

** +

/<А;,

(11.29)

то вектор х* минимизирует величину (11.26).

Для проверки этого утверждения рассмотрим произвольный /-мерный

АЛГЕБРАИЧЕСКИЕ АЛГОРИТМЫ РЕКОНСТРУКЦИИ ИЗОБРАЖЕНИЙ 213

вектор

х.

Пусть

= * - Их

(П.ЗО)

определяется согласно соотношению (11.28), тогда

U =

r(y-RZ).

(11

.3i)

Отсюда следует,

что

\\у - М\\

+ ||* - М

= \\Ч

+ Р||

||и*||

||z*||

(11.32)

[:••]-

поскольку

|—

решение системы (11.27)

минимальной нормой.

Учитывая,

что м*, z*

х*

удовлетворяют уравнениям (11.27)

(11.29),

получаем

и* = tiy - Rx*).

(11.33)

что в сочетании

выражениями (11.29) и (11.32) приводит

неравенству вида

||>;

- /Щ

+ ||* - М

-> г

\\у -

Ях*||

+ ||х* - М

(П/34)

Поскольку

— произвольный /-мерный вектор,

то из

этого следует,

что

х*

минимизирует (11.26).

Из приведенного анализа следует также,

что

любой метод, позволяю-

шли получать решение системы (11.27)

минимальной нормой, автомати-

чески дает

нам

вектор, минимизирующий значение (11.26). Одним

из

спосо-

бов нахождения решения

минимальной нормой совместной системы урав-

нений является релаксационный метод решения уравнений. Заметим,

что,

применяя итерацию (11.22)

системе уравнений (11.27), получим

где

означает транспозицию #-й строки

матрице

(которая оказывается

тождественной

/му

столбцу, поскольку

является единичной матрицей)

1 + г-II г,

Отметим также,

что

если

определено выражением (11.24),

то

нулевой

вектор содержится

в S.

Следовательно, единственным способом, обеспечи-

вающим сходимость

решению

минимальной нормой системы (11.27)

[в

рамках релаксационного метода решения системы уравнений

сочетании

процедурой итерации, аналогичной используемой

(11.35)], является выбор

(0)

и z

(0)

кач

естве нулеяых векторов соответствующей размерности.

214

ГЛАВА 11

— и

—

Их*

{к+1

= и

{к)

+ с

{к

Для всех значений к положим

(к)

= z

(k)

^ (11.37)

Если последовательность z®\ z^

\ z^

\ ... сходится кг*, то последователь-

ность х^°\ х^

\ х&\ ... сходится к величине х*, определяемой условием

(11.29) и минимизирующей значение (11.26).

На практике нет необходимости вводить в алгоритм вычислений вели-

чины £***. Используя выражения (11.35) — (11.37), а также условие, что и^

Z*®

— нулевые векторы, можно получить иной алгоритм. Последователь-

ность оценок х*°\ х

(1)

, х*

\ ... получается при их сходимости к байесовской

оценке х* [при условии что параметры релаксации X*** удовлетворяют со-

отношению (11.10)1, т.е.

(0)

есть /-мерный нулевой вектор ,

Их,

= и

{к)

+ с<%,

= x

(fc)

+ rc

{k)

, (11.38)

где

^>_

«,Ky

ifc

-<^,x*>»-<>

l+r»|rj|* • <»■*>

Заметим, что указанный алгоритм нельзя получить в рамках соотноше-

ния (11.1), однако его применение гораздо сложнее, чем метод описывае-

мый выражением (11.2), поскольку в первом из них требуется дополнитель-

ное определение последовательности /-мерных векторов и*

\ На к-й итера-

ции требуется использовать или заменить лишь один элемент последова-

тельности и^

\ а именно ее #^-Й элемент. Поскольку индекс i

определяется

круговой перестановкой (/

= 1, i

= 2 и т.д.), компоненты вектора и^

можно занести во вспомогательный блок памяти ЭВМ с последовательной

адресацией (доступом) и вводить в основную память ЭВМ один за другим

(или по несколько значений за один раз), когда это потребуется. Аналогич-

ные замечания относятся и к элементам вектора измерений у. Как уже от-

мечалось в разд. 11.1, в наших приложениях значения г

1к

обычно полностью

не запоминаются, а положение и величины ненулевых элементов r

вычис-

ляются по мере необходимости. Следовательно, алгоритм, описываемый

соотношениями (11.38) и (11.39), так же как и простой алгебраический алго-

ритм реконструкции, рассмотренный в разд. 11.1, обладает высокой эффек-

тивностью по использованию объема памяти. Легко показать, что близки

и вычислительные характеристики этих алгоритмов.

Алгоритм, описываемый соотношениями (11.38) и (11.39), является ти-

пичным примером аддитивного алгебраического алгоритма реконструкции.

Оставляя подробное рассмотрение этого алгоритма до следующего раздела,

АЛГЕБРАИЧЕСКИЕ АЛГОРИТМЫ РЕКОНСТРУКЦИИ ИЗОБРАЖЕНИЙ 215

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

реконструкции как процедуру, которая дает последовательность оценок

х®\ х^

\ х&\ ... , сходящуюся к решению с минимальной нормой при

Vi£<r

x>£d

(Н.40)

где 1 < г < / [выражение (6.40)J, так что

(0)

- /-мерный нулевой вектор,

{к+1)

(0)

- 7-мерный нулевой вектор, х

{к +1)

где в свою очередь

с*' = mid{<\ (S

- <r

, x*»»/||rj|

(}',■„

- <r

, x«>y||rj|

}, (11.42)

a mid(w, i;, w) означает среднее значение (медиану) трех величин и, v и и>.

Алгоритм, описываемый выражениями (11.41) и (11.42), в литературе

получил название ART-4 и отличается от других модификаций алгебраиче-

ских алгоритмов, имеющих другие особенности сходимости.

11.4.

НЕКОТОРЫЕ ИСКУССТВЕННЫЕ ПРИЕМЫ

Опыт применения алгебраических алгоритмов реконструкции показыва-

ет,

что эффективность итерационных процедур при реконструкции изобра-

жений часто можно повысить посредством введения некоторых операций

обработки векторов изображений между различными итерациями. Указан-

ные операции в литературе получили название искусственных приемов

(«трюков»).

Рассмотрим подробнее итерационную процедуру в алгебраическом алго-

ритме реконструкций, которая описывается выражением (11.1). Пусть

—

функция, отображающая У-мерные векторы в векторы той же размерно-

сти.

Тогда метод итераций, применяемый к выражению (11.1) в сочетании

с серией искусственных приемов т^, дает последовательность оценок х®\

*Р\ х<

\ ... , определяемых как

{k+l)

= oL

(x<

X (11.43)

№+|

> = т

к+1

(*<*

>). (11.44)

Искусственные приемы особенно полезны, когда они сочетаются с априор-.

ной информацией о пространстве ожидаемых векторов изображений. Иног-

да эти приемы можно использовать для ускорения сходимости к величине

вектора изображения, удовлетворяющего конкретным критериям опти-

мальности. В других случаях они действительно приводят к сходимости к

вектору изображения без оптимизации конкретной функции, что тем не ме-

нее дает лучшие результаты с точки зрения меры расстояний между изо-

бражениями, рассмотренной в разд. 5.1. Последнее имеет место, например,

в том случае, когда ожидаемые дискретные изображения обладают одним

общим свойством, которое нельзя описать простой функцией, но которое

(к)

+ c«V.

(И.41)

216

ГЛАВА П

можно тем не менее выявить при использовании того или иного искусст-

венного приема.

При последующем анализе эвристическое обоснование для всех исполь-

зуемых искусственных приемов будет основано на предположении, что ба-

зисное изображение дискретно и состоит из п х п элизов и определено в со-

ответствии с выражением (6.17).

Прием, который получил название селективное сглаживание, неодно-

кратно использовался нами ранее. Во многих случаях изображения содер-

жат области, в пределах которых значения с большой точностью можно

считать постоянными, а также области со значительными вариациями этих

же величин. Селективное сглаживание позволяет получать изображения та-

кого типа следующим образом.

Пусть i;j, t>

, ... , v

обозначают плотности изображения вр-м элизе и в

его окрестности, как это показано на следующей матрице:

Пусть r, w,,

— вещественные неотрицательные величины, называе-

мые пороговыми значениями и сглаживающими весовыми коэффициента-

ми соответственно. После селективного сглаживания плотность изображе-

ния вр-м элизе будет равна

*Vi + w

Х?=_

fi»i +

> Z?=6 Аъ

где

fl при|г, - fj < Г,

[О для остальных t.

Если

p-Vi

элиз находится на краю изображения и поэтому для некоторых /

величина v

не определена, то мы полагаем f

= 0 для соответствующих

значений /.

На рис. 8.10 и 10.5 можно заметить эффект однократного применения

процедуры селективного сглаживания (в предыдущих разделах именовав-

шейся нелинейным сглаживанием) к данным, полученным с помощью свер-

точного алгоритма для параллельного и веерного пучков соответственно.

Данные табл. 8.3 и 10.2 показывают увеличение точности меры расстояний

между изображениями при использовании рассмотренной процедуры. Во

всех указанных экспериментах принималось / = 0,004, w

{

= 9, w

= 4 и

В тех случаях, когда данный прием применяется в сочетании с алгебраи-

ческими алгоритмами реконструкции, функции

в выражении (11.44) для

селективного сглаживания обычно берутся достаточно редко, например

(11.45)

(11.46)

АЛГЕБРАИЧЕСКИЕ АЛГОРИТМЫ РЕКОНСТРУКЦИИ ИЗОБРАЖЕНИЙ 217

лишь тогда, когда к становится кратным числу измерений /. Для осталь-

ных значений к функцию т

выбирают единичной, не изменяющей величи-

ны вектора изображения.

В противоположность этому искусственный прием ограничения обычно

применяется в алгебраических алгоритмах на каждой итерации и оказыва-

ется оправданным в тех случаях, когда имеется априорная информация о

диапазоне, внутри которого лежат компоненты допустимых векторов изо-

бражений. Например, линейный коэффициент поглощения всегда неотрица-

телен (при любой энергии), а в медицинской диагностике обычно можно

предполагать, что этот коэффициент ограничен сверху величиной коэффи-

циента поглощения плотной костной ткани.

Искусственный прием ограничения можно ввести в алгоритмы, исполь-

зующие разложения в ряд, различными способами: либо в качестве состав-

ной части в систему неравенств (11.5), либо в качестве искусственного прие-

ма при проведении итераций.

Например, если нам известно, что при 1 ^

j'^.

А<х

. </А (11.47)

то полезен следующий прием: положим

(*) = х, (11.48)

где при 1 ^ j ^ У, отсюда

А при Xj < К

Xj при X < Xj < /i, (11.49)

[

fi ПРИ ^ < Xj.

Подобный прием легко может сочетаться с алгебраическими алгоритма-

ми.

Чтобы показать это, рассмотрим релаксационный метод решения си-

стемы неравенств. Последующий анализ показывает справедливость этого

утверждения.

Если множество /V, как это определено условием (11.7), содержит по

меньшей мере один вектор х, компоненты которого удовлетворяют соот-

ношению (11.47), то приводимый ниже алгоритм дает последовательность

оценок х®\ х^\ х^

\ ... , которые сходятся к элементу множества /V, эле-

менты которого удовлетворяют соотношению (11.47), т.е.

(0)

произвольное

ikrl)

lk%

(11.50)

= T

fc+1

(Y<

k+n

где о^ определяется в соответствии с соотношением (11.8) при значениях

параметров \

(/с

\ удовлетворяющих (11.10), а т

определяется условием

(11.48).