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

ВЕСТНИК САНКТ-ПЕТЕРБУРГСКОГО УНИВЕРСИТЕТА

Сер. 4. 2009. Вып. 1

УДК 535.42

Д. А. Кузьмин, Ю. А. Толмачёв

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

Введение. Успешное развитие нейтронной голографии [1, 2] ставит на повестку дня анализ общих свойств метода и, прежде всего, его разрешающей способности. Основной задачей нейтронной голографии является исследование пространственной структуры вещества на атомном и субатомном уровнях расстояний. Если на первых этапах в качестве главной была сформулирована проблема анализа структуры отдельной кристаллической ячейки и измерения расстояний между ближайшими атомами, то сейчас явно звучит задача измерения расстояний до вторых-третьих «соседей», анализа смещений атомов, вызванных малыми возмущениями решётки и определения положения устойчивых комплексов частиц [3]. Вместе с тем, форма аппаратной функции и величина погрешности восстановления положения атомных ядер, рассевающих нейтроны, в явном виде не рассматривалась. Из самых общих волновых свойств частиц следует, что наиболее вероятная ошибка составляет примерно половину длины волны де-Бройля для данной энергии нейтрона, однако, реальная оценка этого важнейшего параметра осталась за скобками исследований.

В представленной работе сделана попытка устранить данную неопределённость в применении к варианту постановки опыта, когда источник опорной нейтронной волны является точечным и находится внутри кристаллической ячейки. Исследована сравнительно простая и ясная по своей физической сущности компьютерная модель эксперимента, конечной целью которого является получение информации о «форме» точечного объекта, восстановленной в результате анализа голограммы. Такая постановка задачи имеет прямое отношение и к проблемам распознавания образов, которые в теории оптической обработки информации, после 60-70-х гг. ХХ в. вновь выходят на первый план в связи с развитием роботизированных систем управления в технике и проблемами анализа сложных белковых молекул в биофизике [3].

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

Постановка задачи. Математическая модель голографического опыта предполагает, что имеется первичный узконаправленный поток монокинетических нейтронов. Первичный поток облучает точечный центр, некогерентно рассеивающий падающий пучок, который создаёт сферическую опорную волну вида Фо = ехр(ИгН)/\Н| (рис. 1). Интересующий нас объект - рассеивающая точка - находится на расстоянии й

© Д. А. Кузьмин, Ю. А. Толмачёв, 2009

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

Пренебрегая затуханием опорной волны на малом расстоянии до рассевающей точки, представим объектную волну в виде:

ei(\k\\d\+\k\\(d-R )|)

Ф = а--------^^,

\d — R\

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

Эксперимент по рассеянию нейтронов ставится всегда так, что R ^ d. Неравенство это выполнено с большим запасом, так как величина d имеет порядок 10-10 м, а R - порядка 1 м. В наших дальнейших расчётах и рассуждениях будем предполагать, что детектор (детекторы) располагаются на сфере, так что R = const.

Вероятность регистрации нейтрона Р в данной точке сферы описывается обычной интерференционной формулой:

P = Р0 [1 + 2a(kd(1 — cos 0))] + O(a), (1)

где 8 - угол между направлением на детектор и направлением на объект, Ро - интенсивность, соответствующая «засветке» регистрирующей системы от источника в отсутствие интерференции двух рассеянных сферических волн, О(а) - поправки более высоких порядков малости, которыми в дальнейшем мы будем пренебрегать. При выводе соотношения (1) предполагалось, что а ^ 1.

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

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

Решение обратной задачи - восстановление положения рассевающего центра по измеренной голограмме - позволяет обнаружить искомый объект в данной точке пространства только с некоторой вероятностью, которую мы будем полагать пропорциональной квадрату модуля соответствующей волновой функции. Физическая модель соответствующего расчёта выглядит следующим образом: сходящаяся в начало координат (в опорный источник) сферическая волна «освещает» голограмму и рассеивается всеми её точками, каждая из которых создаёт свою сферическую расходящуюся волну с амплитудой, пропорциональной величине Р/Ро в данной точке. Складываясь в пространстве, вторичные волны формирую три новые сходящиеся волны: а) обратную опорной, б) создающую правильное действительное изображение интересующего нас объекта и в) - ложное его изображение.

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

з да,! =

3(р)К, !] = 4п2

п

/е—гк(К—р соб 8)

(1 + 2а ссЦЫ — кс1 соэ 8])-—-эш 6с1д

о

п 2

р е—ък(К—р соб 8)

эш 6 <¿8 + А

К

/е—гк(Н—р соб 8)

со в[Ы — к(1 сое 0]--—-вт 0сЮ;

о

3(р)[К, !] = 4п2|3 + 32\2.

2

п

0 20 40 60 80 100 120 140 160 180

0 20 40 60 80 100 120 140 160 180

Рис. 2. Голограммы одной точки в угловых координатах при d = 4 А (а), d = 2, 5 А (б)

Здесь р - текущая координата на луче, совпадающем по направлению с направлением на объект; R, d - параметры прямой задачи.

Интеграл (2) может быть приближённо вычислен аналитически при упрощении его вида путём поворота сферы так, чтобы направление оси симметрии голограммы совпало с направлением на объект. Задачу решаем с точностью до малых порядка R-2, оставив возведение в квадрат при переходе от амплитуды волновой функции к вероятности на заключительный этап.

Рассмотрим оба слагаемых в соотношении (2):

п

«X-WO fe,l.p,r„«anm= 2 ЛЯ^Р=сопв1(р) (3)

R J R кр

о

п

exp(—¿kR) Г e<kd-kdcos е) — e-i(kd-kd cos е) exp(—¿kR)

j2 = 2----------aj ------------------------------С = (A - B)-------------,

о

-ikd

A — _1_____________ (e-ik{p+d) _ e;ifc(p+d)\ .

(p + d)ik V J '

(4)

eikd , , д =_____________ I p_ifc(p-rf) — P*fc(p-d) )

(p - d)ik v J

C = eikp cos 0 sin 0d0.

Окончательно для суммы Ji + J2 имеем:

Ji + J3 = sip(^ ikR

A =

eikp _ e-ikp

0—ikd

В = -------- a (e-tK(P+a> - e

( -ik(p+d.) _ J-k(p+d.)\ .

p + d

(5)

f 1

e—ikd , ,

С =-------a (e-'iHp-d) - eik(p-dA

pd

р — d

Первая дробь (А) в выражении (5) имеет смысл постоянной засветки, и её фазовая зависимость не представляет для нас интереса. Второе и третье слагаемое (В и С) несут информацию о положении объекта и «двойника». Нетрудно заметить, что соответствующей группировкой второго и третьего слагаемого в выражении (5) можно выделить два множителя: ехр(—2ikd — 2гкр)/(р + d) и ехр(2*Ы — 2гкр)/(р — d). Один из них даёт истинное изображение, а второй - изображение двойника. Процедура фильтрации двойников основывается на различии в поведении фаз действительного изображения и двойника при изменении параметра d.

Используя полученные результаты, перепишем выражения (2):

J - 4n2|A - B — C|2, sin kp

А = 2-B = а C = а

kp

2i(cos kd — i sin kd)sink(p + d) (6)

p + d

2i(cos kd + i sin kd) sin k(p — d)

р — d

Это соотношение удобно применить для оценки размеров восстановленного изображения нулевого порядка, для этого достаточно положить d = 0, и учитывая, что а ^ 1, оставить только первое слагаемое. Получаем в первом приближении сферически симметричное относительно начала координат распределение вероятности. Нормируя его в максимуме на единицу, имеем:

(7>

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

Численное моделирование эксперимента по определению координат рассеивающего центра. Если в отсутствие шумов можно однозначно указать положение интересующей нас точки, то при численном моделировании необходимо полностью построить в пространстве «тело неопределённости» и лишь затем принять решение о значении наиболее вероятных координат объекта. Точка, в которой будет найдено максимальное значение интеграла (2) по сфере, будет соответствовать положению объекта с точностью до шага пространственного сканирования. Величина этого шага не может быть бесконечно малой и имеет ограничение, связанное с характерной шириной функции, описывающей пространственное распределение интенсивности вблизи положения объекта. По порядку величины результирующая погрешность должна быть близка к половине длины волны излучения (в данном случае - волны де-Бройля), которая используется при записи и восстановлении. Однако в конечном итоге она определяется отношением сигнал/шум в голограмме.

В проведённом компьютерном моделировании исследовался процесс восстановления изображения одной точки, расположенной на расстоянии d = 4 А (1-й эксперимент), и d = 2, 5 А (2-й эксперимент) от источника. Источник опорной сферической волны X = 0, 7 А находился в начале координат. Голограмма, полученная при расчётах, имеет вид таблицы, столбцы которой содержат значения углов (ф, 8) с угловым шагом 4е и значение интенсивности Р. На рис. 2 центральная часть голограммы по структуре напоминает систему интерференционных полос равной толщины, но при этом надо помнить искажения, вносимые проекцией Меркатора: происходит сильное растяжение полос в полярных областях сферы регистрации, причем «полюс» задаётся для каждого нового положения точки-объекта в пространстве соответствующим направлением радиуса-вектора из опорного источника на объект. В нашем случае положение точки-объекта мы априорно знаем. Задача состоит в том, чтобы изучить форму «тела неопределённости» в восстановленном изображении. Для этого в процессе восстановления,

прежде всего, фиксируется интенсивность 3> непосредственно в точке, соответствующей положению объекта. На полученное таким образом значение нормируются все остальные результаты восстановления «интенсивности» изображения. Затем мы можем найти точки пространства, в которых J будет близко (например, с относительной погрешностью не более 10 %) к заданной наперёд величине 3/3>.

Расчёты проводились в сферической системе координат, шаг восстановления составлял 0,13е по сетке (8, ф). Совокупность полученных значений 3/3> представлялась в сферических координатах (р, 8, ф). Их набор даёт информацию о форме «тела» объекта для фиксированного значения интенсивности. Процедура повторялась с изменённым значением интенсивности от 0,0053> до 0, 983> с заранее выбранным шагом. Окончательным результатом является семейство трёхмерных картинок для выбранных значений 3/3>. Рассматривая их сечение в любой плоскости, можно судить о полученной форме распределения 3/3> в окрестности нашего точечного объекта.

Первым этапом восстановления изображения точки был расчёт 3/3> на сфере малого радиуса, центрированной на начало координат. Величина радиуса сферы изменялась в диапазоне, включавшем в себя значение d, соответствующее положению объекта. В первой серии расчётов, например, начальное значение радиуса было взято равным 3,7 А, конечное - 4,3 А. Во втором, соответственно - 2 А и 3 А, шаг изменения радиуса сферы составлял 0,02 А. Использование столь малого значения позволяет увеличить набор конечных данных, что повышает информативность численного эксперимента.

В результате расчётов был получен 31 файл с данными для первого эксперимента и 51 - для второго. Каждый файл соответствует текущему значению радиуса сканирования р и представляет собой таблицу. Таблицы содержат значения углов (8, ф), интенсивности (3/3>) и р. Для каждого фиксированного радиуса шаг восстановления по углу определяется средней разрешающей способностью. Определение средней разрешающей способности производится для расстояния К = d с помощью специальной подпрограммы. Полученная величина шага делилась на 20 и, таким образом, каждый радиус р, как параметр, определяет сетку по углам (8, ф) с шагом приблизительно 0,13е.

При дальнейшей обработке из семейства таблиц выделялись строки, которые соответствовали заданному значению интенсивности с погрешностью не хуже 5 %. Обработанные строки записывались в новую таблицу, хранящую значения углов ф, 8 и расстояния, центрированного уже на положение объекта (обозначим его 5). Таким образом получалось новое семейство таблиц, отсортированных уже по интенсивности (набор таблиц, соответствующих значениям 3/30 от 0,005 до 0,98). На концах указанного интервала точность была увеличена. Приведённые ниже рисунки, имеющие смысл вероятности обнаружить объект, построены с использованием таких данных.

Форма показанных на рисунках тел соответствует аппаратной функции процесса восстановления объекта. Изображать фигуру «тела вероятности» будем в координатах (8, ф, 5) не забывая, что за начало отсчёта координат взято уже положение объекта, а не опорного источника. Ожидаемой формой является эллипсоид вращения, близкий к сфере. Трёхмерная картина для аппаратной функции, соответствующей точечному объекту, изображена на рис. 3 для уровня 0,1Как показывают результаты численных экспериментов для двух различных положений объекта, форма тела неопределённости действительно близка к сферической.

На основе совокупности полученных данных была построена вероятность обнаружить точку (в нашем случае - величину относительной интенсивности) как функцию расстояния от точки наблюдения до истинного её положения (рис. 4). Данные для двух

- 8

- 10

Рис. 3. Типичный вид пространственного распределения вероятности для 3/30 = 0,1: хорошо видна слоистая структура тела плотности вероятности

значений расстояния между опорным источником и точечным центром ^ = 4,0 А и d = 4, 0 А) совпали в пределах ошибок моделирования. Таким образом, можно сделать вывод, что форма пространственного распределения вероятности обнаружить точечный объект не изменяется в зависимости от его начального положения. Предположив, согласно традиции, что уверенное разрешение пары соседних точек равной интенсивности происходит, если соответствующие распределения интенсивности пересекаются на уровне 0,4 от максимума (провал 20 % интенсивности в суммарном распределении согласно критерию Рэлея), получаем оценку разрешающей способности 0,17-0,18 А, что в пределах ошибок соответствует половине длины волны использованного излучения и полностью согласуется с выводами теории дифракции. Отметим, что в реальном эксперименте положение рассевающего центра не известно, и статистический разброс точек, составляющих «тело неопределённости», будет несколько больше. Тем не менее, корректная статистическая обработка совокупности полученных значений оценки

1,1 “I 1,0 -

0,9 -

0,8 -

0,7 -

0,6 -

0,5 -

0,4 -

0,3 -

0,2 “ 0,1 -

0,0 -

- 0,1 “——i——i——i——i——i——i——i——i——i——i——i——i——i——i——i—|

- 0,1 0,0 0,1 0,2 0,3 0,4 0,5 0,6 0,7 0,8 0,9 1,0 1,1 1,2 1,3 1,4

Рис. 4- Зависимость J/J0 от расстояния:

точки - численный эксперимент, кривая - расчёт по формуле (7)

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

Результаты исследования с очевидностью продемонстрировали осциллирующий характер зависимости вероятности обнаружения объекта в заданной точке от расстояния между ним и источником опорной волны, как это следует из соотношения (6). В качестве примера на рис. 5 приведена структура тела неопределённости, полученная для значения P(р) = 0, 005 (столь малое значение выбрано нами только для того, чтобы сделать более яркими наблюдаемые эффекты). Из этого рисунка хорошо видно, что пространственная структура плотности вероятности действительно имеет характер концентрических слоёв, образованных побочными максимумами функции (sin ар/ар)2. При наличии пространственной периодичности в исследуемой системе, даже при некогерентном рассеянии, в некоторых точках возможно сложение подобных (медленно убывающих по плотности) слоёв, что может привести к образованию ложных изображений. Для уменьшения влияния побочных максимумов необходимо введение специальной процедуры - аподизации - хорошо известной в оптике.

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

10 8 6 4 2 0 - 2

- 4

- 6 - 8

- 10

Рис. 5. Образование множества слоёв при J/Jo = 0,005

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

Осциллирующий характер спада вероятности может приводить к дополнительным систематическим ошибкам при исследовании периодических структур. Этот вывод имеет важнейшее значение для решения задач исследования строения кристаллической ячейки вещества, а также при определении координат априорно известных комплексов частиц. Отметим, что ширина аппаратной функции - это лишь первоначальная оценка точности определения положения рассеивающего центра. Результирующая погрешность будет определяться и параметрами накопленной за время эксперимента статистики, и наличием в ближайшем окружении (на расстоянии порядка ширины тела неопределённости) иных частиц. Подчеркнём также, что в процесс восстановления исследователь непосредственно имеет дело с амплитудой вероятности, а не квадратом её модуля, соответственно, при когерентном рассеянии в описанной нами процедуре побочные максимумы имеют величину приблизительно 20 % от центрального, что усложняет анализ и может приводить к серьёзным качественным ошибкам.

Литература

1. Шарков И. Г., Крекснер Г., Пастор А. А., Толмачёв Ю. А., Торок Д., Чер Л. Численное моделиорование эксперимента по нейтронной голографии // Вестн. С.-Петерб. ун-та. Сер. 4: Физика, химия. 2001. Вып. 3. С. 124-128.

2. Sur B, Rogge R. B, Hammond R. P., Anghel V. N. P., Katsaras J. Atomic structure holography using thermal neutrons // Nature. 2001. Vol. 414. P. 525-527.

3. Пастор А. А., Толмачёв Ю. А., Шарков И. Г. Опознавание образов как метод анализа голограмм сложных молекул // Вестн. С.-Петерб. ун-та. Сер. 4: Физика, химия. Вып. 1. 2004. С. 82-90.

Принято к публикации 10 июня 2008 г.