CC BY
30
УДК 519.6
DOI: 10.14529/ mmp190309
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ РЕКОНСТРУКЦИИ ОБЪЕМНЫХ ИЗОБРАЖЕНИЙ В РЕНТГЕНОВСКОЙ КОМПЬЮТЕРНОЙ ТОМОГРАФИИ С ПРИМЕНЕНИЕМ ГОЛОГРАФИЧЕСКИХ МЕТОДОВ
Е.Н. Симонов1, А.В. Прохоров1, А.В. Акинцева1
1 Южно-Уральский государственный университет, г. Челябинск, Российская Федерация
Рассматривается вычислительный алгоритм решения обратной задачи рентгеновской томографии по реконструкции внутренней структуры микрообъектов в ближней зоне Френеля с применением голографических методов визуализации объемных изображений. Голографические методы не дают прямого решения задачи по реконструкции внутренней структуры объекта. Они могут только решить задачу объемного отображения некоторой поверхности объекта. Однако, используя данные по поглощению рентгеновского излучения объекта и фазоконтрастных голографических сигналов в ближней зоне Френеля, показана возможность получения объемного голографическо-го изображения внутренних слоев объекта. Решение этой сложной задачи потребовало использование трехмерного (3D) преобразования Радона внутренней функции объекта и двумерного (2D) преобразования Радона фазоконтрастной голографической проекции. Получен алгоритм реконструкции фазоконтрастных томографических изображений внутренней структуры объекта, и на основании его разработан вычислительный алгоритм для практической реконструкции объемных томографических изображений внутренней структуры микрообъектов. Результаты исследований подтверждены математическим моделированием алгоритма реконструкции объемных изображений, для чего была разработана математическая модель тестового фантома, и для него смоделированы фазоконтрастные проекции с последующей реконструкцией по ним фазо-контрастных томографических изображений томографическими методами с использованием разработанного авторами программного реконструктора.
Ключевые слова: 'рентгеновская компьютерная томография; голография; двумерное и трехмерное преобразование Радона.
Введение
В работах [1-8] рассматривались различные подходы получения фазоконтрастных изображений некоторых слоев объекта исследования в рентгеновских лучах, используя физические основы отражения и преломления по принципам Лауэ и Брэгга. В этих работах также проводится попытка связать фазоконтрастные изображения с го-лографическими. Однако, разработка метода получения внутренней структуры в любой точке объекта по фазоконтрастным изображениям (проекциям) не была достигнута. В работах [9, 10] сформулирована задача реконструкции внутренней структуры сечений объекта в фазоконтрастной рентгеновской томографии (показатель преломления) на основе трехмерного преобразований Радона (для микрообъектов).
Фазоконтрастные изображения могут быть получены применением осевой голографии по принципу Габора [11, 12]. Для иллюстрации этого принципа, рассмотрим облучение объекта монохроматическим когерентным рентгеновским лучом (рис. 1). Абсорбционно-контрастное изображение наблюдается на плоскости, контактной с
объектом. Фазоконтрастные изображения наблюдаются в ближнем поле области преломления Френеля. Чистый контраст поглощения, который используется в традиционной КТ [13], наблюдается только в плоскости, контактной со снимком (то есть когда расстояние между объектом и плоскостью изображения является нулевым). Фазовый контраст из-за дифракции рентгеновских лучей происходит на всей области дифракции Френеля. Поэтому, для того чтобы получить фазоконтрастное изображение, достаточно отодвинуть детектор от объекта.
РИ
Рис. 1. Принцип когерентного получения изображений
1. Основные соотношения по решению обратной
голографической задачи внутренней структуры объекта
Теория [9, 10], основанная на фундаментальном отношении между трехмерным (3В) преобразованием Радона функции объекта и двумерным (2В) преобразованием Радона фазоконтрастной проекции, дает возможность решения обратной голо-графической задачи по визуализации внутренней структуры объекта с применением алгоритмов реконструкции изображений классической рентгеновской компьютерной томографии. Вводятся базовые основы когерентного построения изображений, выводится формула, описывающая линейное отношение между фазовой функцией и яркостью в ближней зоне дифракции Френеля.
Введем декартову система координат относительно объекта исследования (х\, х2,х3) (рис. 2). Детектор снабжен собственной системой координат (х,у), основание которой лежит в плоскости (х1, х2), а направление оси у совпадает с направлением оси х3. Смешанный фазо-амплитудный объект можно описать коэффициентами линейного поглощения ^ (х1,х2, х3) и функцией / (х1,х2,х3) = п (х1,х2, х3) — 1, где п (х1,х2,х3) представляет действительную часть показателя преломления. Мы полагаем, что монохроматическое когерентное волновое поле с длиной волны Л облучает объект в направлении, определенном углом обзора в. Реализация подобной геометрии может быть достигнута вращением объекта вокруг оси х3.
Поток поля монохроматической волны объекта можно описать следующим образом:
ив (х,у) = Т (х,у) иг, (1)
где иг - плоская монохроматическая волна после взаимодействия, Т (х, у) - комплексная функция передачи объекта на угле в.
Контактное Область Френеля
поле
Рис. 2. Координатная система для объекта и плоскости наблюдения
Тв (х, у) = Мв (х, у) ехр [гфв (х, у)] , с модулем Мв (х, у) = ехр [-\цв (х, у)] , где
цв (x,y)= / ц (xi, x2,y) S (x — xi cos 0 — x2 sin 0) dx\dx2,
R2
и фазовая функция Рв (x,y)
2тг
т
f (xi, x2, y) S (x — xi cos 0 — x2 sin 0) dxidx2.
(2)
(3)
R2
Следует учесть, что выражения (2) и (3) описывают двумерное преобразование Радона в плоскости, перпендикулярной оси х3. Распределение яркости на расстоянии г от объекта может быть представлено следующим образом:
i| (x,y ) = |2
(4)
п
^ + у'
, ±± . . ехр (%к£)
где обозначение - двумерная свертка, п? (х, у) = ---ехр
глг
функция распространения Френеля, при этом к = - волновое число.
Задачей восстановления является нахождение функции f (х\,х2,х3) по измеренным значениям (х, у), 0 < в < п. Чтобы найти аналитическое решение этой задачи, сделаем некоторые приближения, которые будут соблюдаться в ближнем поле Френеля. Сначала произведем преобразование Фурье для функции распределения Френеля:
Н (С, п) = ехр (гкг) ехр \-in\z (С2 + п2)] ,
где С и п являются пространственными частотами. Для малых Лг частотное представление функции может быть упрощено:
Н (С, п) = ехр (¿Ь) [1 - гпЛг (С2 + п2)] , где мы сохранили только линейные части ряда Маклорена для экспоненциальной функции. Таким образом, яркость на малом расстоянии г становится равной:
II (x, y)
2 . тт
1 + -7-1 Ue
(5)
где V2 - лапласиан.
2
При незначительном изменении показателя поглощения объекта ^ (х,у), что соответствует для малоконтрастных объектов, распределение яркости на значительно малом расстоянии й от объекта определится следующим образом:
1 - (х,у)
Id (x,y) = I0
(6)
2п
По своему происхождению равенство (6) описывает распределение яркости в ближней зоне Френеля, где Az <С • Равенство (6) устанавливает линейную связь между фазовой функцией и измеренной яркостью (яркостью на детекторе). Использование равенства (6) для представления яркости помогает установить прямую связь между данными яркости и искомой функцией f (xi, x2,x3) объекта, которую необходимо восстановить. Для этого требуется преобразование Радона, которое будем применять дважды: в плоскости x—y и в трехмерном пространстве (x1,x2, x3). В первом случае мы имеем дело с двумерным преобразованием Радона, во втором случае - с трехмерной версией преобразования Радона.
Двумерное преобразование Радона можно представить как
R2 (s, w) = У g (x, y) 8 (s — x sin w — y cos w) dxdy.
R2
Оно описывает интегрирование функции g (x,y) по линиям, определенным равенством x sin ш + y cos w = s в плоскости x — y, тогда как трехмерное преобразование Радона
R3 (s, 0, w) = У f (x1, x2, x3) 8 (s — (x1 cos w + x2 sin w) sin w — x3 cos w) dx1 dx2dx3
R3
описывает интегрирование функции f (x1, x2, x3) по плоскостям. Такая плоскость перпендикулярна единичному вектору n = (cos 0 sin w,sin0 sin w, cos w) и имеет расстояние s до исходной плоскости облучения объекта (рис. 1). Мы можем сформулировать следующий результат. Пусть
ge (x,y) = Id (x,y) /10 (x,y) — 1, (7)
где Ifi (x,y) - яркость распределения на малом расстоянии d, яркость в плоскости, контактной с объектом, равна I0 = |Me (x, y) U¿|2 , тогда:
d2 —1
— Н3(з,е,ш) = —Н2в(з,ш). (8)
2. Общий алгоритм реконструкции фазоконтастных
томографических изображений внутренней структуры объекта
Используя представленные выше соотношения, можем вывести алгоритм реконструкции. Функцию объекта можно найти по формуле (8), если обратить трехмерное преобразование Радона. Точная обратная формула трехмерного преобразования Радона, выглядит следующим образом [13]:
п п
— 1 Í' ^52
/ (х1,х2,х3) = — / siníüduj / —R3(s,9,uj)d9, (9)
где, выполняя дифференцирование, s заменяется следующей формулой:
s' = (xi cos 9 + x2 sin 9) sin ш + x3 cos ш.
Примечательным фактом является то, что эта формула применяет вторую производную к преобразованию Радона. Это значит, что, подставив выражение (8) в формулу (9), мы сразу перейдем к формуле
f (Xi ,Х2,Хз)
1
d4n2
sin шdш / R2e (s', ш) d9,
(10)
0 0
которая представляет решение задачи восстановления фазоконтрастной КТ.
Формула (10) предполагает, что двумерное преобразование Радона в плоскости х — у должно быть просчитано для каждой двумерной проекции д$ (х, у). Затем этот результат восстанавливается в трехмерном пространстве, используя при этом оператор восстановления трехмерного обратного преобразования Радона. Следует учесть, что эта формула является устойчивой, так как все используемые операторы являются непрерывными и не требуют дифференцирования.
Формула (10) не имеет каких-либо особенностей и может непосредственно использоваться в численной реализации алгоритма восстановления. Тем не менее, применение двумерного преобразования Радона к каждой угловой проекции потребует больших вычислений. Чтобы ускорить эти вычисления, можно использовать быстрое преобразование Фурье. В то же время, можно получить более практический алгоритм восстановления, если упростить формулу (10), вычислив интеграл по углу ш. Затем нам следует объединить вычисление двумерного преобразования Радона и восстановления по углу ш на каждом единичном шаге алгоритма. Чтобы проделать это, для начала перепишем формулу (10), подставляя выражение преобразования Радона для данных д$. Это даст нам
f (Xi ,Х2,Хз):
1
d4n2.
d9 J J ge (x, y) dxdyy ó I (x1 cos 9+x2 sin 9 — x) sin ш+ о ^
+ (x3 — y) cos ш ) sin шdш.
(11)
Интеграл по углу ш может быть выражен умножением и делением интеграла на выражение
\J{х\ cos 9 + х2 sin 9 — х)2 + (хз — у)2.
Пусть
тогда
Y = tan
cos y =
sin y
_1 x1 cos 9 + x2 sin 9 — x
x3 — y
x3 — y
(x1 cos 9 + x2 sin 9 — x)2 + (x3 — y)2
x1 cos 9 + x2 sin 9 — x (x1 cos 9 + x2 sin 9 — x)2 + (x3 — y)2
n
П
П
Затем, применяя формулу 8 (а (Ь)) = 8 (Ь — Ь') / (Ь')|, которая представля-
ет дельта-функцию 8 (а (Ь)), и используя корень Ь' уравнения а (Ь) = 0, получим:
J 8 ((xi cos в + х2 sin в — х) sin ш + (х3 — y) cos ш) sin ш^ш 0
п
1
8 (sin y sin ш + cos y cos ш) sin ш^ш
y''(xi COS в + X2 sin в — xf + (х3 — y)2 о
п
= —. = / 8 (cos (ш — y)) sinwdu; =
22 д/ (x1 cos в + x2 sin в — х) + (х3 — y) 0
п
1 г
8 (cos (ш — ш')) sin ш^ш
(sin (а;' — y) I у''(xi cos в + х2 sin в — х)2 + (х3 — у)
|хз — У1
2
0
(x1 cos в + х2 sin в — х) + (х3 — y)2
п
где мы использовали ш = 7 + 2 •
Подставляя значение интеграла в (11) и используя
q(x,y) = }v} 2,
х2 + У2
можем записать:
п
f (xi,x2,x3) = J q**ged9, (12)
о
в котором, следуя свертке, аргументы х и y функции q**g$ (двумерная свертка) заменяются соответственно на х = х1 cos в + х2 sin в и y = х3. Равенство (12) представляет алгоритм восстановления в форме фильтрующего восстановления.
Этот интеграл свертки можно рассчитать в области Фурье, что также позволит сократить время вычислений применением быстрого преобразования Фурье. Соответствующая функция фильтра Q (С, n) может быть найдена применением двумерного преобразования Фурье к ядру свертки q (х,у), что даст:
ехр [—¿2тг (х£ + уг/)] dxdy = ^ х2 + y2 2
1С I
Q(€,V)= II 2 , 2 ехР Н2?г « + У^)] dxdy = - I ехр (—у |£|) ехр (—12тгуг]) dy
Qa (С,П) =
£2 + п2'
Функция фильтра Q (£, п) имеет бесконечно много нулей вдоль линии £ = 0, и единственный полюс в начале координат. Чтобы точно аппроксимировать полюс на дискретной решетке, мы должны использовать малый интервал дискретизации в частотной области. Этого можно достичь, обнулив дискретные наборы данных, что является стандартным способом решения схожих проблем в восстановлении традиционной КТ. Дискретная функция фильтра показана на рис. 3.
Следует учесть, что Q (£, п) является функцией переноса низкочастотного фильтра, и поэтому не требуется дополнительной регуляризации для подавления высокочастотных компонентов зашумленных данных.
п
Рис. 3. Нормализованная дискретная функция п) фильтра восстановления
3. Вычислительный алгоритм реконструкции
фазоконтрастных томографических изображений внутренней структуры объекта
Определение искомой функции показателя преломления f (xi,x2,x3) определяет-
п
ся формулой (11), где / ó ((xi cos в+x2 sin в — x) sin ш +(x3 — y) cos ш) sin wdw будет равен
JV+ 2 1 учитывая, ЧТО X = X\ COS в + Ж2 SÍI10, у = Жз И Cü' = Cü + Формула (11) примет вид
f (xi,x2,x3) =
1
d4n2
de ge (x,y)
|y|
x2 + y2
dxdy.
(13)
Рассмотрим выражение (13).
I. В формуле (13) двойной интеграл представляет двумерную свертку функции дв(х, у) и ядра свертки 2. Двумерную свертку можно определить через двумерное
X +у
преобразование Фурье, используя двумерное быстрое прямое преобразование, а затем обратное преобразование (БПФ):
а) двумерное преобразование Фурье функции Л1 2 будет Qi^rj) =
X +у
_ И
б) определение функции до (ж, у): Вариант 1. Имеется слабоконтрастный объект, у которого коэффициент поглощения ^(х, у) изменяется незначительно. Тогда яркость на детекторе при сканирова-
нии
объекта (x, y) определяется по формуле (6),
как
rd _ т0
/e = /e
i ^d 2 Í \
~ 27TV WWiy)
где (х,у) - лапласиан фазовой функции (х, у). Подставляя последнее выражение в (7), получим
ж, у) =
(14)
п
где фазовая функция (ж, у) определяется как
2П [
'•Рв(х, у) = — / ¡(х1,х2,у)5{х — ХхСОъв — Х2$шв)с1х1с1х2-
в2
Моделирование дв (ж, у) проводится по формуле (14).
Вариант 2. Общий случай: имеется любой объект, у которого коэффициент поглощения ^(ж, у) изменяется каким-либо образом. Тогда дв (ж, у) определяется по формуле (7), как дв(ж,у) = /?(ж,у)/10(ж,у) — 1, где 10(ж,у) яркость на детекторе после объекта при г = 0, т.е. детектор соприкасается с объектом, а /^(ж,у) яркость на детекторе после объекта, находящегося на расстоянии г = й от детектора. /0(ж,у) и /^(ж,у) определяются экспериментальным путем при измерении яркости на детекторе в каждой его точке (пикселе) с координатами (ж, у) при каждом значении угла поворота объекта 0.
Моделирование дв(ж, у) по формуле дв(ж, у) = /^(ж, у)/10(ж, у) — 1 требует моделирования яркости /0(ж,у) и /^(ж,у). Для моделирования значений яркостей /0(ж,у) и /^(ж,у) можно воспользоваться формулой /|(ж,у) = |Л./*ив|2. Функцию распределения Френеля кг можно представить через фурье-образ Иг(£,п) = ещ>{гкг) [1 — + г]2)] Поток волнового поля 11в(х,у) полностью опреде-
ляется комплексной функцией передачи объекта Тв (ж, у), которая, в свою очередь, определяется, как
Те (ж, у) = Мв (ж, у) ехр [г^в (ж, у)],
где
Ме (ж, у) = ехр
^в (ж, у) = у Мжъ ж2, у)^(ж — ж1 сов 0 — ж2 вШ 0)йж1йж2,
в2
2П Г
<рв(х,у) = — / /(ж1, ж2, у)^(ж — Ж1 сое 0 — ж2 вт 6|)б?ж1б?ж2.
в2
в) полученное в варианте 1 или варианте 2 двумерное преобразование Фурье дв (ж, у) обозначим Ов (С,п). Умножим его на возьмем от произведения Овдвумерное обратное преобразование Фурье, результат будет равен
'ж2 + у
Обозначим этот двойной интеграл, как Вв (ж, у)
дв(х,у)—--йхйу, формула (13). (15)
II. В формуле (13) интеграл / Вв (ж, у)й0 можно определить любым из известных
о
методов.
4. Вычислительный эксперимент
Реализация вычислительного алгоритма осуществлялась в среде Ма11аЬ, которая эффективно осуществляет матричное вычисление интегралов, дифференциалов, прямого и обратного одномерного и двумерного преобразования Фурье, одномерной
и двумерной свертки, одномерного и двумерного преобразования Радона с задаваемыми одномерными, двумерными и трехмерными матрицами переменных. Для того, чтобы оценить подход, использовалась модель фантома (рис. 4). Фазовый объект был представлен в виде цилиндра диаметром 300 мкм и высотой 300 мкм, диаметр вставок 30 мкм.
Рис. 4. Компьютерная модель фантома. Показаны пять его сечений
Значения яркости серо-белого для объекта были —5 * 10-7 для цилиндра и
- 10 * 10-7, —4 * 10-7 и —2 * 10-7 для вставок. Значения яркостей /0(ж,у) и !^(х,у) моделировались по варианту 1 вычислительного алгоритма.
Параметры переменных для реконструкции изображений:
- длина волны была принята равной Л = 0,1нм, размер пикселя детектора -А = 0,1мкм;
- угол сканирования объекта 9 изменялся 0 < 9 < п с дискретностью А9 = 0, 50, т.е. было сделано 360 проекций;
- матрица восстанавливаемого изображения 512 х 512 пикселей;
- центр изображения находится в центре координат хОу, размер изображения по горизонтали (по оси х) равен 512 пикселей (±256 пикселей), т.е. каждый пиксель изображения равен 300/512 мкм размера объекта, количество пикселей детектора в одном пикселе изображения равно 300/512/0,1=3000/512. Параметр х изменяется с дискретностью Ах = 300/512 мкм, т.е. с дискретностью 3000/512 пикселей детектора;
- аналогично по параметру у.
В вычислительном эксперименте были смоделированы фазоконтрастные данные на расстоянии d=3 см от объекта. Для одного из сечений объекта результаты представлены на рис. 5, где показаны значения яркости между 0,7 и 1,1. Структуры на рис. 5 появляются в результате фазового контраста. Восстановленные изображения с параметрами реконструкции, указанными выше, показаны на рис. 6. Для восстановления сечений объекта по формуле (13) использовался разработанный авторами реконструктор томографического изображения [13, 14]. Сравнивая изображения на рис. 4 и 6, можно видеть, что предлагаемый метод восстанавливает изображение достаточно хорошо.
Выводы
1. Проведено усовершенствование теории и алгоритма фазоконтрастной компьютерной томографии [9, 10]. Рассмотрены варианты восстановления изображения для слабоконтрастных объектов и общий случай - для любых объектов. Более детально рассмотрены алгоритмы реконструкции фазоконтрастных томографических изображений внутренней структуры объекта для практической разработки программных реконструкторов.
Рис. 5. Яркость данных (х, у) для одного из сечений объекта при й =3 см
Рис. б. Восстановленные сечения объекта
2. Установлена связь между трехмерным преобразованием Радона функции преломления объекта и двумерным преобразованием Радона фазоконтрастных проекций, отражающих соотношение яркости на детекторе в ближней зоне Френеля (размер объекта гораздо меньше Лz) и яркости на детекторе в зоне поглощения (z = O). Эта связь дает возможность осуществить реконструкцию фазоконтрастных изображений томографическими методами, принятыми в рентгеновской KT.
3. Моделирование объекта и фазоконтрастных проекций, последующая реконструкция по ним фазоконтрастных томографических изображений томографическими методами показало хорошую схожесть результатов. При реконструкции использовался разработанный авторами реконструктор томографического изображения.
4. Разработан компактный алгоритм реконструкции фазоконтрастных томографических изображений для экспериментальных установок по исследованию малоразмерных биологических объектов, когда для получения фазоконтрастных проекций осуществляется два массива измерений: измерения яркости от объекта на детекторе в ближней зоне Френеля и измерения от объекта на детекторе в зоне поглощения.
5. Практическое применение разработанного подхода в получении фазоконтраст-ных томографических изображений томографическими методами может найти там, где требуется более высокое разрешение на изображении, чем при применении классических методов KT.
Литература
1. Anastasio, M.A. Analytic Image Reconstruction in Local Phase-Contrast Tomography У M.A. Anastasio, D. Shi, F. De Carlo, X. Pan ^ Physics in Medicine and Biology. - 2004. -V. 49, № 1. - P. 121-144.
2. Cloetens, P. Holotomography: Quantitative Phase Tomography with Micrometer Resolution Using Hard Synchrotron Radiation X-Rays У P. Cloetens, W. Ludwig, J. Baruchel, D. Van Dyke, J. Van Landuyt, J.P. Guigay, M. Schlenker У У Applied Physics Letters. - 1999. - V. 75, № 19. - P. 2912-2914.
3. Davis, T. Phase-Contrast Imaging of Weakly Absorbing Materials Using Hard X-Rays У T. Davis, D. Gao, T. Gureyev, A. Stevenson, S. Wilkins ^ Nature. - 1995. - V. 373, № 6515. - P. 595-598.
4. Groso, A. Implementation of a Fast Method for High Resolution Phase Contrast Tomography / A. Groso, R. Abela, M. Stampanoni // Optics Express. - 2006. - V. 14, № 18. -P. 8103-8110.
5. Groso, A. Phase Contrast Tomography: an Alternative Approach / A. Groso, M. Stampanoni, R. Abela, P. Schneider, S. Linga, R. Muller // Applied Physics Letters. - 2006. - V. 88, № 21. -P. 214104.
6. Gureyev, T.E. Phase-and-Amplitude Computer Tomography / T.E. Gureyev, D.M Paganin, G.R. Myers, Ya.I. Nesterets, S. Wilkins // Applied Physics Letters. - 2006. - V. 89, № 3. -P. 034102.
7. Snigirev, A. On the Possibilities of X-Ray Phase Contrast Microimaging by Coherent High-Energy Synchrotron Radiation / A. Snigirev, I. Snigireva, V. Kohn, S. Kuznetsov, I. Schelokov // Review of Scientific Instruments. - 1995. - V. 66, № 12. - P. 5486-5492.
8. Wilkins, S.W. Phase-Contrast Imaging Using Polychromatic Hard X-Rays / S.W. Wilkins, T.E. Gureyev, D. Gao, A. Pogany, A.W. Stevenson // Nature. - 1996. - V. 384, № 6607. -P. 335-338.
9. Bronnikov, A.V. Reconstruction Formulas in Phase-Contrast Tomography / A.V. Bronnikov // Optics Communications. - 1999. - V. 171, № 4-6. - P. 239-244.
10. Bronnikov, A.V. Theory of Quantitative Phase-Contrast Computed Tomography / A.V. Bronnikov // Journal of the Optical Society of America A. - 2002. - V. 19, № 3. -P. 472-480.
11. Ярославский, Л.П. Цифровая голография / Л.П. Ярославский, H.P. Мерзляков. - М.: Наука, 1982.
12. Gabor, D. A New Microscopic Principle / D. Gabor // Nature. - 1948. - V. 161, № 4098. -P. 777-778.
13. Симонов, Е.Н. Рентгеновская компьютерная томография / Е.Н. Симонов. - Снежинск: РФЯЦ-ВНИИТФ, 2002.
14. Симонов, Е.Н. Физика визуализации изображений в рентгеновской компьютерной томографии / Е.Н. Симонов. - Челябинск: Издат. центр ЮУрГУ, 2014.
Евгений Николаевич Симонов, доктор технических наук, профессор, кафедра «Техника, технология и строительство:», Институт открытого и дистанционного образования, Южно-Уральский государственный университет (г. Челябинск, Российская Федерация), e.n.simonov@yandex.ru.
Александр Владимирович Прохоров, кандидат технических наук, доцент, заведующий кафедрой «Современные образовательные технологии», Институт открытого и дистанционного образования, Южно-Уральский государственный университет (г. Челябинск, Российская Федерация), prokhorovav@susu.ru.
Александра Викторовна Акинцева, кандидат технических наук, старший преподаватель кафедры «Техника, технология и строительство», Институт открытого и дистанционного образования, Южно-Уральский государственный университет (г. Челябинск, Российская Федерация), akintsevaav@susu.ru.
Поступила в редакцию 8 ноября 2018 г.
MSC 68 DOI: 10.14529/mmp190309
MATHEMATICAL MODELLING OF RECONSTRUCTION OF VOLUMETRIC IMAGES IN X-RAY COMPUTED TOMOGRAPHY USING HOLOGRAPHIC METHODS
E.N. Simonov1, A.V. Prokhorov1, A.V. Akintseva1
1 South Ural State University, Chelyabinsk, Russian Federation
E-mails: e.n.simonov@yandex.ru, prokhorovav@susu.ru, akintsevaav@susu.ru
We consider a computational algorithm for solving the inverse problem of x-ray tomography to reconstruct the internal structure of micro objects in the near Fresnel zone using holographic methods of visualization of three-dimensional images. Holographic methods do not give a direct solution to the problem on reconstruction of the internal structure of the object. These methods can only solve the problem on volumetric mapping of some object surface. However, using the data on the absorption of both x-ray radiation of the object and phase contrast holographic signals in the near Fresnel zone, we show the possibility to obtain a volumetric holographic image of the inner layers of the object. In order to solve this complex problem, it is necessary to use a three-dimensional (3D) radon transformation of the internal function of the object and a two-dimensional (2D) radon transformation of the phase contrast holographic projection. We obtain an algorithm to reconstruct phase-contrast tomographic images of the internal structure of an object. Based on the algorithm, we construct a computational algorithm for the practical reconstruction of volumetric tomographic images of the internal structure of microobjects. The results of the research were confirmed by the mathematical modelling of the algorithm to reconstruct three-dimensional images. To this end, we develop a mathematical model of the test phantom, and simulate the phase contrast projections for the test phantom. Then, we develop a software in order to reconstruct the phase contrast tomographic images by tomographic methods on the basis of the obtained phase contrast projections.
Keywords: x-ray computed tomography; holography; two-dimensional and three-dimensional radon transformation.
References
1. Anastasio M.A., Shi D., De Carlo F., Pan X. Analytic Image Reconstruction in Local Phase-Contrast Tomography. Physics in Medicine and Biology, 2004, vol. 49, no. 1, pp. 121-144. DOI: 10.1088/0031-9155/49/1/009
2. Cloetens P., Ludwig W., Baruchel J., Van Dyke D., Van Landuyt J., Guigay J.P., Schlenker M. Holotomography: Quantitative Phase Tomography with Micrometer Resolution Using Hard Synchrotron Radiation X-Rays. Applied Physics Letters, 1999, vol. 75, no. 19, pp. 2912-2914. DOI: 10.1063/1.125225
3. Davis T., Gao D., Gureyev T., Stevenson A., Wilkins S. Phase-Contrast Imaging of Weakly Absorbing Materials Using Hard X-Rays. Nature, 1995, vol. 373, no. 6515, pp. 595-598. DOI: 10.1038/373595a0
4. Groso A., Abela R., Stampanoni M. Implementation of a Fast Method for High Resolution Phase Contrast Tomography. Optics Express, 2006, vol. 14, no. 18, pp. 8103-8110. DOI: 10.1364/0E.14.008103
5. Groso A., Stampanoni M., Abela R., Schneider P., Linga S., Muller R. Phase Contrast Tomography: an Alternative Approach. Applied Physics Letters, 2006, vol. 88, no. 21, p. 214104. DOI: 10.1063/1.2207221
6. Gureyev T.E., Paganin D.M, Myers G.R., Nesterets Ya.I., Wilkins S. Phase-and-Amplitude Computer Tomography. Applied Physics Letters, 2006, vol. 89, no. 3, p. 034102. DOI: 10.1063/1.2226794
7. Snigirev A., Snigireva I., Kohn V., Kuznetsov S., Schelokov I. On the Possibilities of X-Ray Phase Contrast Microimaging by Coherent High-Energy Synchrotron Radiation. Review of Scientific Instruments, 1995, vol. 66, no. 12, pp. 5486-5492. DOI: 10.1063/1.1146073
8. Wilkins S.W., Gureyev T.E., Gao D., Pogany A., Stevenson A.W. Phase-Contrast Imaging Using Polychromatic Hard X-Rays. Nature, 1996, vol. 384, no. 6607, pp. 335-338. DOI: 10.1038/384335a0
9. Bronnikov A.V. Reconstruction Formulas in Phase-Contrast Tomography. Optics Communications, 1999, vol. 171, no. 4-6, pp. 239-244.DOI: 10.1016/S0030-4018(99)00575-1
10. Bronnikov A.V. Theory of Quantitative Phase-Contrast Computed Tomography. Journal of the Optical Society of America A, 2002, vol. 19, no. 3, pp. 472-480. DOI: 10.1364/JOSAA.19.000472
11. Yaroslavsky L.P., Merzlyakov N.S. Cifrovaya golografiya [Digital Holography]. Moscow, Nauka, 1982. (in Russian)
12. Gabor D. A New Microscopic Principle. Nature, 1948, vol. 161, no. 4098, pp. 777-778. DOI: 10.1038/161777a0
13. Simonov E.N. Rentgenovskaya komp'yuternaya tomografiya [X-Ray Computed Tomography]. Snezhinsk: Publishing House of the Russian Federal Nuclear Center-VNIITF, 2002. (in Russian)
14. Simonov E.N. Fizika vizualizacii izobrazhenij v rentgenovskoj komp'yuternoj tomografii [Physics of Image Visualization in X-Ray Computed Tomography]. Chelyabinsk: Publishing Centre SUSU, 2014. (in Russian)
Received November 8, 2018
