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

Информатика и вычислительная техника

УДК 519.6 DOI: 10.14529/ctcr180401

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

Е.Н. Симонов, Д.В. Аврамов

Южно-Уральский государственный университет, г. Челябинск, Россия

Рассматривается вычислительный алгоритм решения обратной задачи рентгеновской томографии по реконструкции внутренней структуры микрообъектов в ближнем зоне Френеля с применением голографических методов визуализации объемных изображений. Голографиче-ские методы не дают прямого решения задачи по реконструкции внутренней структуры объекта. Они могут только решить задачу объемного отображения некоторой поверхности объекта. Однако, используя данные по поглощению рентгеновского излучения объекта и фазокон-трастных голографических сигналов в ближней зоне Френеля, показана возможность получения объемного голографического изображения внутренних слоев объекта. Решение этой сложной задачи потребовало использование трехмерного (3D) преобразования Радона внутренней функции объекта и двухмерного (2D) преобразования Радона фазоконтрастной голо-графической проекции. Был получен алгоритм реконструкции фазоконтастных томографических изображений внутренней структуры объекта и на основании его разработан вычислительный алгоритм для практической реконструкции объемных томографических изображений внутренней структуры микрообъектов. Результаты исследований были подтверждены математическим моделированием алгоритма реконструкции объемных изображений. Для чего, была разработана математическая модель тестового фантома и для него смоделированы фазо-контрастные проекции с последующей реконструкцией по ним фазоконтрастных томографических изображений томографическими методами с использованием разработанного авторами программного реконструктора. Разработан компактный алгоритм реконструкции фазокон-трастных томографических изображений для экспериментальных установок по исследованию малоразмерных биологических объектов, когда для получения фазоконтрастных проекций осуществляется два массива измерений: измерения яркости от объекта на детекторе в ближней зоне Френеля и измерения от объекта на детекторе в зоне поглощения. Практическое применение разработанного подхода в получении фазоконтрастных томографических изображений томографическими методами может найти для исследования микрообъектов, где требуется высокое разрешение на изображении, чем при применении классических методов КТ.

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

Введение

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

В работах [9-11] была сформулирован задача реконструкции внутренней структуры сечений объекта в фазоконтрастной рентгеновской томографии (показатель преломления) на основе трехмерного преобразования Радона (для микрообъектов).

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

принципу Габора [12, 13]. Для иллюстрации этого принципа рассмотрим облучение объекта монохроматическим когерентным рентгеновским лучом (рис. 1).

Область Френеля Рис. 1. Принцип когерентного получения изображений

Абсорбционно-контрастное изображение наблюдается на плоскости, контактной с объектом. Фазоконтрастные изображения наблюдаются в ближнем поле области преломления Френеля.

Чистый контраст поглощения, который используется в традиционной КТ [14], наблюдается только в плоскости, контактной со снимком (то есть когда расстояние между объектом и плоскостью изображения является нулевым. Фазовый контраст из-за дифракции рентгеновских лучей происходит на всей области дифракции Френеля.

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

Основные соотношения по решению обратной голографической задачи внутренней структуры объекта

Теория [9-11], основанная на фундаментальном отношении между трёхмерным (3D) преобразованием Радона функции объекта и двухмерным (2D) преобразованием Радона фазоконтраст-ной проекции, дает возможность решения обратной голографической задачи по визуализации внутренней структуры объекта с применением алгоритмов реконструкции изображений классической рентгеновской компьютерной томографии. Вводятся базовые основы когерентного построения изображений, выводится формула, описывающая линейное отношение между фазовой функцией и яркостью в ближней зоне дифракции Френеля.

Введем декартову систему координат относительно объекта исследования (хх, х2, х3) (рис. 2).

Рис. 2. Координатная система для объекта и плоскости наблюдения

Детектор снабжён собственной системой координат (х, у), основание которой лежит в плос-

х.

кости (х1, х2), а направление оси у совпадает с направлением оси

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

Мы полагаем, что монохроматическое когерентное волновое поле с длиной волны X облучает объект в направлении, определённом углом обзора 0 .

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

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

х, у ) = Т0( х, у )и, (1)

где и - плоская монохроматическая волна после взаимодействия; Т0(х, у) - комплексная функция передачи объекта на угле 0

Т0 (х у) = м0 (х у) ехР [ф0 (х у)]

- 2 ^0 (х у)

, где

с модулем М 0 (х, у) = ехр

^0(х,у)= | х1,х2,у)5(х-x1cos0-х^т0)dx1dx2,

д2

и фазовая функция 2л г

Ф0 (х,у) = — I /(х1,х2,у)5(х - x1cos0- x2sin0)dx1dx2. X 2

(2)

(3)

Следует учесть, что выражения (2) и (3) описывают двухмерное преобразование Радона в плоскости, перпендикулярной оси х3 . Распределение яркости на расстоянии г от объекта может

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

|2

/0г (х,у) = |йг * *и0,

(4)

где

К (х,у) = ех|Мехр

' хг (х2+у2)

- функция распространения Френеля; двойной звёздочкой обозначена двухмерная свёртка, при

7 2л

этом к = — . X

Задачей восстановления является нахождение функции f (х1, х2, х3) по измеренным значениям /0г (х, у), 0 < 0 < л . Чтобы найти аналитическое решение этой задачи, сделаем некоторые приближения, которые будут соблюдаться в ближнем поле Френеля. Сначала произведём преобразование Фурье для функции распределения Френеля:

Нг (£,= ехр {¡кг)ехр -¡лХг (^2 + )

где и ^ являются пространственными частотами. Для малых Хг частотное представление функции может быть упрощено:

Нг (£, л) = ехр (¡кг)

1 - ¡лХг (^2 +^2 )

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

¡Хг 2 ^ '

/0г ( х, у )■■

1 +—V2 и 4л

0

(5)

где V - лапласиан.

д

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

Ie (х, 7) = V

1 "|^2фе( х, 7)

(6)

По своему происхождению равенство (6) описывает распределение яркости в ближней зоне Френеля, где Xz << ——- .

Равенство (6) устанавливает линейную связь между фазовой функцией и измеренной яркостью (яркостью на детекторе).

Использование равенства (6) для представления яркости помогает установить прямую связь между данными яркости и искомой функцией f (х-, х2, х3) объекта, которую необходимо восстановить. Для этого требуется преобразование Радона, которое будем применять дважды: в плоскости х-у и в трёхмерном пространстве (х1, х2, х3) . В первом случае мы имеем дело с двухмерным

преобразованием Радона, во втором случае - с трёхмерной версией преобразования Радона. Двухмерное преобразование Радона можно представить, как

R2 (s, ю) = J g ( х, у )5( s - х sin ю - у cos ю) dxdy .

R2

Оно описывает интегрирование функции g (х, у) по линиям, определённым равенством х sin ю + у cos ю = s в плоскости х-у, тогда как трёхмерное преобразование Радона

R3 (s,е,ю) = J f (х1,х2,х3 )5(s - (х1 cosю + х2 sinю)sinю - х3 cosю)dх1dх2dх3

R3

описывает интегрирование функции f (х-, х2, х3) по плоскостям. Такая плоскость перпендикулярна единичному вектору n = (cos е sin ю,sin е sin ю,cos ю) и имеет расстояние s до исходной

плоскости облучения объекта (см. рис. 2).

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

ge (х,у) = V (х,у)/(х,у)-1, (7)

где Ied (х, у) - яркость распределения на малом расстоянии d.

Яркость в плоскости, контактной с объектом, равна Ie° = Me (х,у)U| , тогда

82 -1

—R3 (s, e, ю) = — R2e( s, ю). (8)

8s d

Общий алгоритм реконструкции фазоконтрастных томографических изображений внутренней структуры объекта

Используя выше представленные соотношения, можем вывести алгоритм реконструкции. Функцию объекта можно найти по формуле (8), если обратить трёхмерное преобразование Радона.

Точная обратная формула трёхмерного преобразования Радона, открытая Радоном [14], выглядит следующим образом:

-1 г г 82

f ( х1, х2, х3) = —- J sin юd ю1 d e—- R3 ( s, e, ю), (9)

4г ° ° 8s

где, выполняя дифференцирование, s заменяется следующей формулой s' = (х1 cose + х2 sine)sinю + х3 cosю.

Примечательным фактом является то, что эта формула применяет вторую производную к

преобразованию Радона. Это значит, что, подставив выражение (8) в формулу (9), мы сразу перейдём к формуле

1 К К

f ( x1, x2, Х3

)= tj7~2" í SÍn®í d 0r29( ^ ^ (10)

d4к о о

которая представляет решение задачи восстановления фазоконтрастной КТ.

Формула (10) предполагает, что двухмерное преобразование Радона в плоскости x-y должно быть просчитано для каждой двухмерной проекции ge (x, y) . Затем этот результат восстанавливается в трёхмерном пространстве, используя при этом оператор восстановления трёхмерного обратного преобразования Радона. Следует учесть, что эта формула является устойчивой, так как все используемые операторы являются непрерывными и не требуют дифференцирования.

Формула (10) не имеет каких-либо особенностей и может непосредственно использоваться в численной реализации алгоритма восстановления. Тем не менее применение двухмерного преобразования Радона к каждой угловой проекции потребует больших вычислений. Чтобы ускорить эти вычисления, можно использовать быстрое преобразование Фурье. В то же время можно получить более практический алгоритм восстановления, если упростить формулу (10), вычислив интеграл по углу ю . Затем нам следует объединить вычисление двухмерного преобразования Радона и восстановления по углу ю на каждом единичном шаге алгоритма. Чтобы проделать это, для начала перепишем формулу (10), подставляя выражение преобразования Радона для данных ge . Это даст нам

1 К К

f (x1, x2, x3 ) = —TTídeííge(xy)dxdyí5((x1cose + x2sine-x)sinю+(x3 -y)coSю)sinюю (11) d4" 0 0

Интеграл по углу ю может быть выражен умножением и делением интеграла на выражение

yj(x1cos e + x2sin e- x)2 +(x3 - y )2 . Пусть

_i x, cos e + x2 sin e - x

y = tan —-2-,

x3 - y

тогда

x3 - y . x1cos e + x2 sin e- x

COS у = , 2 2' Sln Y= I-i-T

J( XjCOS 0 + x2sln 0-x) +(x3 - y) J( XjCOS 0 + x2sin 0-x) +(x3 - y)

Затем, применяя формулу 5(a(b)) = 5(b -b') /|(da / db)(b')|, которая представляет дельта-

функцию 5(a(b)), и, используя корень b' уравнения a(b) = 0 , получим: п

|s(( x1 cos 0 + x2 sln 0 - x) sin ra + ( x3 - y ) cos ra) sln rad ra =

0

■\j(x1 cos e + x2 sin e - x)2 + (x3 - y)' 0 1 К

"Í5( sin y sin ю + cos у cos ю) sin rad ю =

---cos (ю - y)) sin od ю =

д/(x1cos e + x2sin e- x) +(x3 - y ) 0

1 К -, í§(cos(ю -ю'))sinюdю :

I sin (ю'-у^( x1cos e + x2sin e-x) +( x3 - y ) 0

x3 - y

(x1cos e + x2sin e- x) +(x3 - y ) t Л

где мы использовали ю = у +—.

2

1

2

Подставляя значение интеграла в (11) и используя

q (X, y )=-г^т,

X + y можем записать

1 f

f (xi, x2, x3 ) = -J- J d 0 **ge,

(12)

в котором, следуя свёртке, аргументы х и у функции q (двухмерная свертка) заменяются соответственно на х = x1cos 0 + X2sin 0 и у = Х3 . Равенство (12) представляет алгоритм восстановления в форме фильтрующего восстановления.

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

свёртки д ( х, у ) , что даст

У 1 ^

-ехр[-+ уп)]^ = 21ехР("У 1^1)ехР(-72лУП)Ф = ва л) :

Q & ^) = JJ

X2 + y2

+л2

Функция фильтра в л) имеет бесконечно много нулей вдоль линии ^ = 0 и единственный

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

Дискретная функция фильтра показана на рис. 3.

Рис. 3. Нормализованная дискретная функция Q(Д,i) фильтра восстановления

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

Вычислительный алгоритм реконструкции фазоконтрастных томографических изображений внутренней структуры объекта

Определение искомой функции показателя преломления f (х1, x2, х3) определяется формулой (11)

1 f f

f (x1, x2, x3 ) =-2 JdeJ J ge (x,y)dxdyJs((x1 cos 0 + x2 sin 0 - x)sin ю + (x3 - y)cosю)sin radю,

d4f 0 0

f

где J5((XjCos0 + x2sin0-x)sinю + (x3 -y)cosю^тюdю будет равен

0

y f

' ' ,учитывая, что x = x1 cos0 + x2sin0, y = x3 и ю' = ю+—.

X2 + y2 2

Формула (11) примет вид

1 f С С I V I

f (x1, ^ x3 ) = J d0J J g0 (X y)~2-2 dxdy . (13)

dо X + y

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

1. В формуле (13) двойной интеграл представляет двухмерную свертку функции g0(x, y) и яд-

|y|

ра свертки ' 2 .

X2 + y2

Двухмерную свертку можно определить через двухмерное преобразование Фурье, используя двухмерное быстрое прямое преобразование, а затем обратное преобразование (БПФ):

y I ^2 |

а) двухмерное преобразование Фурье функции 2 2 будет Q = —-- ;

X2 + y2 +П2

б) определение функции g0(X, y).

Вариант 1. Имеется слабоконтрастный объект, у которого коэффициент поглощения ^,(x, y) изменяется незначительно.

Тогда яркость на детекторе при сканировании объекта 70^ (X, y) определяется по формуле (6), как

I0d ( X, y )= I

1 X, y )

где У2ф0 ( х, у) - лапласиан фазовой функции ф0 (х, у) . Подставляя последнее выражение в (7), получим

£0(х, у ) = -М У2ф0 (х, у), (14)

2л

где фазовая функция ф0 (х, у) определяется, как 2л

т

(X,y) = — J f (Xj,X2,y)5(X - Xjcos0- X2sin0)dr1dX2 .

□2

R2

Моделирование g0(x, у) проводится по формуле (14).

Вариант 2. Общий случай: имеется любой объект, у которого коэффициент поглощения ц(х, у) изменяется каким-либо образом.

Тогда g0(x,у) определяется по формуле (7), как g0(х, у) = 10й (х, у) /10° (х, у)-1, где 10(х, у) -

яркость на детекторе после объекта при г = 0, т. е. детектор соприкасается с объектом, а 10^ (х, у) -яркость на детекторе после объекта, находящегося на расстоянии г = й от детектора.

/0(х, у) и 10^ (х, у) определяются экспериментальным путем при измерении яркости на детекторе в каждой его точке (пикселе) с координатами (х, у) при каждом значении угла поворота объекта 0.

Моделирование £е(х,у) по формуле ge(x,у) = 1°с(х,у)/!е0(х,у)-1 требует моделирования яркости/е°(х,у) и IеТ(х,у). Для моделирования значений яркостей 1°(х,у) и (х,у) можно

воспользоваться формулой (х,у) = \к2 * *и° |2.

Функцию распределения Френеля Нг можно представить через фурье-образ

Н2 (£,^) = ехр (¡кг) 1 - ¡пХг + ) Поток волнового поля и° (х, у) полностью определяется комплексной функцией передачи

k=

X

объекта TQ (x, y), которая, в свою очередь, определяется, как

те (x,y) = MQ (x, y)exp [7фе (x, y)]:

где

1

- ^ ^е (x, у )

M е (x, y ) = exp ^е (x,y) = J x1,x2,y)5(x-x1 cosе-x2sinе)dx1dx2 ,

R2

Фе (x,y) = — J f (x1,x2,y)5(x - x1 cosе - x2 sinе)dx1dx2;

в) полученное в варианте 1 или варианте 2 двухмерное преобразование Фурье ge(x, у) обозначим Ge .

Умножим его на Q , возьмем от произведения Ge(^,Q двухмерное обратное

преобразование Фурье, результат будет равен (* (* ge (х, у) 2| у 1 2

х2 + у2

dxdy, формула (13).

Обозначим этот двойной интеграл как В° (х, у).

п

2. В формуле (13) интеграл | В° (х, у)С e можно определить любым из известных методов.

Вычислительный эксперимент

Реализация вычислительного алгоритма осуществлялась в среде МА^АВ, которая эффективно осуществляет матричные вычисления интегралов, дифференциалов, прямого и обратного одномерного и двухмерного преобразования Фурье, одномерной и двухмерной свертки, одномерного и двухмерного преобразования Радона с задаваемыми одномерными, двухмерными и трехмерными матрицами переменных.

Для того чтобы оценить подход, использовалась модель фантома (рис. 4). Фазовый объект был представлен в виде цилиндра диаметром 300 мкм и высотой 300 мкм, диаметр вставок 30 мкм.

Рис. 4. Компьютерная модель фантома. Показаны пять его сечений

Значения яркости серо-белого для объекта были -5 10- для цилиндра и -10-10, -4-10 и -2-10-7 для вставок.

Значения яркостей 7°(х, у) и (х, у) моделировались по варианту 1 вычислительного алгоритма.

Длина волны была принята равной X = 0,1 нм, размер пикселя детектора - А = 0,1 мкм.

Параметры переменных для реконструкции изображений:

- угол сканирования объекта 0 изменялся 0 < 0 < л с дискретностью Л0 = 0,5°, т. е. было сделано 360 проекций;

- матрица восстанавливаемого изображения 512x512 пикселей;

- центр изображения находится в центре координат хОу, размер изображения по горизонтали (по оси х) равен 512 пикселей (±256 пикселей), т. е. каждый пиксель изображения равен 300/512 мкм размера объекта, количество пикселей детектора в одном 300/512/0,1 = 3000/512. Параметр х изменяется с дискретностью Лх = 300/512 мкм, т. е. с дискретностью 3000/512 пикселей детектора;

- аналогично по параметру у.

В вычислительном эксперименте были смоделированы фа-зоконтрастные данные на расстоянии й = 3 см от объекта. Для одного из сечений объекта результаты представлены на рис. 5, где показаны значения яркости между 0,7 и 1,1.

Структуры на рис. 5 появляются в результате фазового контраста.

Восстановленные изображения с параметрами реконструкции, указанными выше, показаны на рис. 6. Для восстановления сечений объекта по формуле (13) использовался разработанный авторами реконструктор томографического изображения [14, 15].

Рис. 6. Восстановленные сечения объекта

пикселе изображения равно

Рис. 5. Яркость данных (х, у)

для одного из сечений объекта при d = 3 см

Сравнивая изображения на рис. 4 и 6, можно видеть, что предлагаемый метод восстанавливает изображение достаточно хорошо.

Выводы

1. Проведено усовершенствование теории и алгоритма фазоконтрастной компьютерной томографии [9-11]. Рассмотрены варианты восстановления изображения для слабоконтрастных объектов и общий случай - для любых объектов. Более детально рассмотрены алгоритмы реконструкции фазоконтрастных томографических изображений внутренней структуры объекта для практической разработки программных реконструкторов.

2. Установлена связь между трёхмерным преобразованием Радона функции преломления объекта и двухмерным преобразованием Радона фазоконтрастных проекций, отражающих соотношение яркости на детекторе в ближней зоне Френеля (размер объекта гораздо меньше Аг) и яркости на детекторе в зоне поглощения (г = 0). Эта связь дает возможность осуществить реконструкцию фазконтрастных изображений томографическими методами, принятыми в рентгеновской КТ.

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

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

ний: измерения яркости от объекта на детекторе в ближней зоне Френеля и измерения от объекта на детекторе в зоне поглощения.

5. Практическое применение разработанного подхода в получении фазоконтрастных томографических изображений томографическими методами может найти там, где требуется более высокое разрешение на изображении, чем при применении классических методов КТ.

Литература/References

1. 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. Rev. Sci. Instrum., 1995, 66, pp. 5486-5492. DOI: 10.1063/1.1146073

2. Davis T., Gao D., Gureyev T., Stevenson A., Wilkins S. Phase Contrast Imaging of Weakly Absorbing Materials Using Hard X-Rays. Nature, 1995, 373, pp. 595-598. DOI: 10.1038/373595a0

3. Wilkins S.W., Gureyev T.E., Gao D., Pogany A., Stevenson A.W. Phase-Contrast Imaging Using Polychromatic Hard X-Rays. Nature, 1996, 384, pp. 335-338. DOI: 10.1038/384335a0

4. 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. Appl. Phys. Let., 1999, 75, pp. 2912-2914. DOI: 10.1063/1.125225

5. Anastasio M.A., Shi D., De Carlo F., Pan X. Analytic Image Reconstruction in Local Phase-Contrast Tomography. Phys. Med. Biol., 2004, 49, pp. 121-144. DOI: 10.1088/0031-9155/49/1/009

6. Groso A., Stampanoni M., Abela R., Schneider P., Linga S., Muller R. Phase Contrast Tomography: an Alternative Approach. Appl. Phys. Let., 2006, 88, p. 214104. DOI: 10.1063/1.2207221

7. Groso A., Abela R., Stampanoni M. Implementation of a Fast Method for High Resolution Phase Contrast Tomography. Optics Express, 2006, vol. 14, iss. 18, pp. 8103-8110. DOI: 10.1364ЮЕ.14.008103

8. Gureyev T.E., Paganin D.M., Myers G.R., Nesterets Ya.I., Wilkins S. Phase-and-Amplitude Computer Tomography. Appl. Phys. Let., 2006, 89, p. 034102. DOI: 10.1063/1.2226794

9. Bronnikov A.V. Theory of Quantitative Phase-Contrast Computed Tomography. J. Opt. Soc. Am. A, 2002, vol. 19, no. 3, March, pp. 472-480. DOI: 10.1364/JOSAA.19.000472

10. Bronnikov A.V. Reconstruction Formulas in Phase-Contrast Tomography. Optics Communications, 1999, 171, pp. 239-244. DOI: 10.1016/S0030-4018(99)00575-1

11. Bronnikov A.V. Theory of Quantitative Phase-Contrast Computed Tomography. J. Opt. Soc. Am. A., 2002, 19, pp. 472-480. DOI: 10.1364/JOSAA.19.000472

12. Ярославский Л.П., Мерзляков Н.С. Цифровая голография. М.: Наука, 1982. 219 с. [Yaroslavskiy L.P., Merzlyakov N.S. Tsifrovaya golografiya [Digital Holography]. Moscow, Nauka Publ., 1982. 219 p.]

13. Gabor D. A New Microscopic Principle. Nature, 1948, 161, pp. 777-778. DOI: 10.1038/161777a0

14. Симонов Е.Н. Рентгеновская компьютерная томография: моногр. Снежинск: Изд-во РФЯЦ-ВНИИТФ, 2002. 364 с. [Simonov E.N. Rentgenovskaya komp'yuternaya tomografiya: mono-grafiya [X-Ray Computer Tomography: Monograph]. Snezhinsk, RFNC VNIITF Publ., 2002. 364 p.]

15. Симонов Е.Н. Физика визуализации изображений в рентгеновской компьютерной томографии. Челябинск: Издат. центр ЮУрГУ, 2014. 505 с. Simonov E.N. Fizika vizualizatsii izobra-zheniy v rentgenovskoy komp'yuternoy tomografii [Physics of Visualization of Images in a X-Ray Computer Tomography]. Chelyabinsk, South Ural St. Univ. Publ., 2014. 505 p.

Симонов Евгений Николаевич, д-р техн. наук, профессор кафедры техники и технологии Политехнического института, Южно-Уральский государственный университет, г. Челябинск; e.n.simonov@yandex.ru.

Аврамов Денис Витальевич, аспирант, Южно-Уральский государственный университет, г. Челябинск; avramov.91@mail.ru.

Поступила в редакцию 15 марта 2018 г

DOI: 10.14529/ctcr180401

DEVELOPMENT OF COMPUTATIONAL ALGORITHM

FOR RECONSTRUCTION OF VOLUMETRIC IMAGES IN X-RAY

COMPUTED TOMOGRAPHY USING HOLOGRAPHIC METHODS

E.N. Simonov, e.n.simonov@yandex.ru,

D.V. Avramov, avramov.91@mail.ru

South Ural State University, Chelyabinsk, Russian Federation

A computational algorithm for solving the inverse problem of x-ray tomography for reconstruction of the internal structure of micro objects in the near Fresnel zone using holographic methods of visualization of three-dimensional images is considered. Holographic methods do not provide a direct solution to the problem of reconstruction of the internal structure of the object. They can only solve the problem of volumetric display of some object surface. However, using data on the absorption of x-ray radiation of the object and phase-contrast holographic signals in the near Fresnel zone, the possibility of obtaining a three-dimensional holographic image of the inner layers of the object is shown. The solution of this complex problem required the use of three-dimensional (3D) Radon transformation of the internal function of the object and two-dimensional (2D) Radon transformation of phase-contrast holographic projection. An algorithm for reconstruction of phase-contrast tomographic images of the internal structure of the object was obtained and based on it a computational algorithm for the practical reconstruction of three-dimensional tomographic images of the internal structure of micro objects was developed. The results of the research were confirmed by mathematical modeling of the algorithm for reconstruction of three-dimensional images. For this purpose, a mathematical model of the test phantom was developed and phase-contrast projections were simulated for it with subsequent reconstruction of phase-contrast tomographic images by tomographic methods using the program reenactor developed by the authors. A compact algorithm for reconstruction of phase-contrast tomographic images for experimental facilities for the study of small biological objects, when two arrays of measurements are carried out to obtain phase-contrast projections: measurements of brightness from the object on the detector in the near Fresnel zone and measurements from the object on the detector in the absorption zone. The practical application of the developed approach in obtaining phase-contrast tomographic images by tomographic methods can be found for the study of micro objects where high resolution is required in the image than in the application of classical CT methods.

Keywords: x-ray computed tomography, holography, two-dimensional and three-dimensional radon transformation.

Received 15 March 2018

ОБРАЗЕЦ ЦИТИРОВАНИЯ

Симонов, Е.Н. Разработка вычислительного алгоритма реконструкции объемных изображений в рентгеновской компьютерной томографии с применением голографических методов / Е.Н. Симонов, Д.В. Аврамов // Вестник ЮУрГУ. Серия «Компьютерные технологии, управление, радиоэлектроника». -2018. - Т. 18, № 4. - С. 5-15. DOI: 10.14529/ctcr180401

FOR CITATION

Simonov E.N., Avramov D.V. Development of Computational Algorithm for Reconstruction of Volumetric Images in X-Ray Computed Tomography Using Holographic Methods. Bulletin of the South Ural State University. Ser. Computer Technologies, Automatic Control, Radio Electronics, 2018, vol. 18, no. 4, pp. 5-15. (in Russ.) DOI: 10.14529/ctcr180401