Реконструкция распределений случайных функций в задачах стохастической дифракционной томографии* Текст научной статьи по специальности «Математика»

УДК 519.2

0.В. Шестаков1

РЕКОНСТРУКЦИЯ РАСПРЕДЕЛЕНИЙ СЛУЧАЙНЫХ ФУНКЦИЙ В ЗАДАЧАХ СТОХАСТИЧЕСКОЙ ДИФРАКЦИОННОЙ ТОМОГРАФИИ*

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

Ключевые слова: стохастическая томография, преобразование Радона, проекции, случайные функции, дифракционная теорема.

1. Введение. Задача определения вероятностных характеристик двумерных случайных функций по характеристикам одномерных проекций возникает при исследовании объектов, которые меняются случайным образом во время процесса регистрации проекций. Если для получения проекционных данных используются рентгеновские или 7-лучи, то проекции можно достаточно точно описать преобразованием Радона (см. [1]). Однако в некоторых приложениях (например, биологических) применяется более безвредный тип излучения, например ультразвук или электромагнитные волны. При этом уже нельзя игнорировать волновую природу излучения. В результате возникает проблема восстановления коэффициента преломления, который представляет собой случайную функцию.

В работах [2-5] рассматривается задача определения вероятностных характеристик двумерных случайных функций в предположении, что проекции описываются преобразованием Радона. Показывается, что в общем эта задача характеризуется сильной неоднозначностью, и содержательные результаты удается получить лишь в том случае, когда случайная функция имеет не более чем счетное число состояний. В работе [5] для класса таких функций разрабатывается метод восстановления вероятностных распределений.

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

2. Основные математические соотношения дифракционной томографии. В задачах дифракционной томографии рассматривается объект с коэффициентом преломления п(х) = (1 + /(ж))1/2 (х € К2). Объект облучается волной е~ш«/(ж) с частотой к (мы рассматриваем только плоские волны и/(ж) = егк(в'х\ где в € Б1 — направление распространения волны). Зарегистрировав рассеянную волну е~гЫи3(х) за пределами объекта, требуется восстановить функцию /(ж). В данной работе мы будем предполагать, что /(ж) ^ 0 и имеет носителем круг II единичного радиуса с центром в начале координат. Функции, совпадающие всюду за исключением множеств нулевой лебеговой меры, будем считать эквивалентными.

Известно, что и = и/ + и3 удовлетворяет волновому уравнению

Аи + к2(1 + /)« = 0, (1)

где А — оператор Лапласа. Воспользовавшись аппроксимацией Рытова (см. [1]), можно получить приближенное решение п,ц уравнения (1) в форме п,ц = где шц удовлетворяет уравнению

А (и/шд) + к2(и1Шл) = —к/и1. (2)

Решив уравнение (2), получаем (см. [1])

и1(х)юя(х) =-к ! Н0(к \х - у\)}{у)щ(у) (¿у, (3) _ я2

1Факультет ВМиК МГУ, к.ф.-м.н., ассистент, e-mail:oshestakov@cs.msu.su.

* Работа выполнена при финансовой поддержке РФФИ (грант 08-01-00567-а).

где

R1

Здесь 4iHQ — функция Хенкеля первого рода нулевого порядка, a a(t) =

^fk2 — t2. При использовании

аппроксимации Рытова для нахождения коэффициента преломления требуется решить уравнение (3) относительно /(ж) при заданной функции д(в, s) = WR^d+sd-1), где г > 1 — фиксированное число (т. е. шц измеряется за пределами диска радиуса г). При каждом фиксированном в G S1 будем называть функцию g(0,s) проекцией функции /(ж). Преобразовав уравнение (3) (см. [6]), имеем

ik Г pira(t) г

д(е, s) = —e~ikr / eist-~— / e~i{r *)/(r'0 + dr'ds'dt. (4)

4тг J a(t) J

R R2

Далее, переходя в (4) к преобразованиям Фурье, получаем следующее соотношение между двумерным преобразованием Фурье функции /(ж) и одномерным преобразованием Фурье функции g (в, s) по второй переменной:

f((a(t)-k)e + te±) = -(^j ' i^-eir^~a^g(e,t). (5)

Это так называемая дифракционная теорема. При изменении t

в интервале [—k^kj значение

(a(t) - к)0 + te1- = (Vk2 - t2 - к)0 + te1- - кв

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

Соответственно, если в принимает значения из S1, то выражение (5) определяет / в круге радиуса не меньше к. Естественно, для восстановления /(ж) требуется знать ее преобразование Фурье всюду. В общем случае для этого требуется использовать аналитическое продолжение, что представляет собой очень неустойчивую задачу. Однако при некоторых дополнительных предположениях удается разработать устойчивые методы восстановления (см., например, [6]).

3. Стохастические объекты. Перейдем теперь к рассмотрению объектов, меняющихся случайным образом в процессе регистрации проекционных данных. Формально постановка задачи, рассматриваемой в настоящей работе, следующая. Имеется двумерная случайная функция £(ж) — стохастический объект вида £(ж) = /„(ж), где fi(x), /2(ж),... — последовательность интегрируемых функций, определенных в единичном круге U, a v — случайная величина, принимающая целые положительные значения. Вероятностная структура £(ж) полностью определяется набором

(/1(ж),/2(ж),... ;pi,p2, • • ■),

00

где pi = -Р(£(ж) = fi(ж)), i = 1,2,..., Pi = 1- Распределение £(ж) будем обозначать через а

г= 1

распределения проекций (при фиксированном в G S1) через Р^в).

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

Теорема. Пусть случайные функции £(ж) и г](х) имеют описанный выше вид, и

PçW = Рп(в)

для всех в G Л Ç S1, где Л - бесконечное множество, тогда

= Л/т. е. распределение двумерной случайной функции однозначно определяется распределениями любого бесконечного множества проекций.

Доказательство данной теоремы полностью повторяет доказательство аналогичной теоремы в работе [3]. В дальнейшем для простоты будем полагать что Л совпадает с S1.

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

Для простоты изложения в данной работе мы будем рассматривать класс функций, имеющих всего два состояния. Обобщение на любое конечное число состояний очевидно, а для случая счетного числа состояний можно произвести процедуру "усечения" распределений проекций так же, как это делается в [5].

Итак, пусть случайная функция £(ж) принимает значения fi(x) и /г(ж) с вероятностями pi и р2. Предполагается, что известны распределения проекций для всех в € S1, т.е. для каждого в € S1 известны функции дг(в, s), i = 1, 2, являющиеся проекциями функций /¿(ж), / 1. 2. и реализующиеся с вероятностями pi и р2 соответственно. Причем, вообще говоря, заранее неизвестно, какое состояние проекции какому состоянию функции соответствует, т.е. может быть так, что gi(6,s) является проекцией /г(ж), а д2(0,s) — проекцией fi(x). Необходимо разделить функции дг(в, s), i = 1,2, для всех в € S1 на группы так, чтобы каждая группа состояний проекций относилась к одному состоянию случайной функции.

Если pi ф р2, тогда такое разделение можно произвести по вероятностям состояний проекций, т.е. для всех в € S1 значение gi(6,s), которое реализуется с вероятностью рi, мы относим к первой группе, а значение дг(в, s), которое реализуется с вероятностью , — ко второй.

В случае, когда р\ = = ■§, метод группировки проекций основан на некоторых следствиях из дифракционной теоремы (5). Принимая во внимание, что в = в(а) = (cosa, sin а), а й1 = 0±(а) = = (— sin a, cos а), запишем соотношение (5) в следующем виде:

/f\l/2 e¿r(a(t)-fc) л

д(а, t) = I — ) гк-г-f((a(í) — к) cosa — t sin a, (a(t) — к) sin a + t cos a). (6)

\ ¿t f Qj\tJ

Поскольку функция / аналитична, g(a,t) дифференцируема по í в окрестности нуля бесконечное число раз. Обозначим д^т\а, í)|t=o через J(™)(a). Из (6) следует, что функция J^m\a) представляет собой многочлен от cosa и sin a степени не выше п. В частности,

(|)1/2Ш0,0),

/7г \ 1/2 л

(-J i(fl2 (0, 0) cos a -f'ti (0,0) sin a), /тгу/2 .fl^irk í(n /4 (0,0) cos a+ 4 (0,0) sin a

Ы --к-+

+fí[ti (0, 0) sin2 a - 2 fl[t2 (0, 0) cos a sin a + (0, 0) cos2 a)

и т. д. Эти соотношения похожи на хорошо известные условия Людвига-Хелгасона для преобразования Радона (см. [7]). Воспользуемся ими для группировки проекций.

Если функции gi(a,s), i = 1,2, различны, то найдется номер т, для которого j!>m\a) различаются (это следует из аналитичности преобразований Фурье от функций /¿(ж)). Поскольку функции

J¡m\a) представляют собой тригонометрические многочлены, они либо тождественно равны, либо пересекаются в конечном числе точек. Пусть S задано так, что найдется точка а, в окрестности которой разница между функциями J^m\a), i = 1,2, по модулю больше 8. Посчитав значения J^m\aj), i = 1,2, в 2т + 1 близких точках a:¡ из этой окрестности, мы можем, пользуясь непрерывностью J¡m\a), разделить их на группы, соответствующие каждому состоянию случайной функции (если разности J^7l\a,j) — J^io/j) и J^n\ajjri) — «/^"^(aj+i) имеют разный знак, значит, значения j[m\ctj) и j[m\aj+1) соответствуют различным состояниям случайной функции, и следует поменять местами j[m\aj+i) и J^n\a:jjri)), а затем найти явный вид функций J¡m\a), i = 1,2, от переменной а. Для

JÍ°)(a) = jW(a) =

J(2)(a) =

этого нужно решить следующие две линеиные системы уравнении:

то то

п=1 п=1

гт,а + ^2am,ncos(naj) + ^2Kn,nsin(naj) = Jim\a:j)", i = 0, ...,2m, г = 1,2, (7)

и найти коэффициенты агт п и Ьгт п тригонометрических многочленов, которыми представляются функции J^m\a).

Вывод явного вида функций J^m\a) аналогичен выводу интерполяционного многочлена Лагранжа (подробности можно найти в [8]). Решая (7), имеем

2 то 2то • (а—ак\

л(то)(а) = Е^то)ы Д i = 1>2. (8)

j=0 к=0,кф] 2 >

Далее, вычисляя для каждого а значение J^m\a) как производную m-го порядка по t от gi(a,t) при t = 0, мы относим проекции к той или иной группе в зависимости от того, со значением какой из найденных функций (8) в точке а совпадает это вычисленное значение.

После того как проекции распределены по группам, можно восстановить каждое состояние случайной функции, а значит, и ее распределение обычными методами дифракционной томографии (см., например, [1] или [6]).

СПИСОК ЛИТЕРАТУРЫ

1. Как А. С., Slaney М. Principles of Computerized Tomographic Imaging. IEEE Press, 1988.

2. Аристов В.В., Ушаков Н.Г., Ушакова А.П. Стохастическая томография // Докл. РАН. 1998. 359. № 4. С. 464-466.

3. Ушаков В. Г., Ушаков Н.Г. Восстановление вероятностных характеристик многомерных случайных функций по проекциям // Вестн. Моск. ун-та. Сер. 15. Вычисл. матем. и киберн. 2001. № 4. С. 32-39.

4. Шестаков О.В. О единственности восстановления вероятностных характеристик многомерных случайных функций по вероятностным характеристикам их проекций // Вестн. Моск. ун-та. Сер. 15. Вычисл. матем. и киберн. 2003. № 3. С. 37-41.

5. Shestakov O.V. An algorithm to reconstruct probabilistic distributions of multivariate random functions from the distributions of their projections // J. Math. Sci. 2002. 112. N 2. P. 4198-4204.

6. Klebanov L. В., Rachev S. T. On a special case of the basic problem in diffraction tomography // Stochastic Models. 1996. 12. N 2. P. 181-197.

7. Луис А. К., Наттерер Ф. Математические проблемы реконструктивной вычислительной томографии / / ТИИЭР. 1983. 71. С. 111-125.

8. Гончаров В. J1. Теория интерполирования и приближения функций. М., 1934.

Поступила в редакцию 15.09.08