Известия РАН. Серия физическая, 2019, T. 83, № 2, стр. 198-202

ВВЕДЕНИЕ

Метод компьютерной томографии позволяет восстанавливать 3D-структуру объекта без его физического разрушения [1]. Это достигается путем математической обработки набора измеренных томографических 2D-проекций. Основным этапом обработки является процедура реконструкции, которая переводит данные из пространства измерения в пространство объекта. В индустриальных томографических комплексах для реконструкции в основном используются интегральные методы. Основным их достоинством является высокое быстродействие. Однако методы не позволяют использовать априорные знания о восстанавливаемом изображении, что затрудняет их применение в нестандартных условиях измерения, т.е. при малом времени экспозиции и при небольшом числе проекционных углов, которые, к тому же могут быть распределены неравномерно. В таких случаях используется алгебраический подход [1, 2]. Интегральные выражения, которые связывают результат измерения с пространственным распределением коэффициента ослабления рентгеновского излучения, в алгебраическом подходе заменяются на лучевые суммы, и далее решается система алгебраических уравнений. Следует обратить внимание, что количество неизвестных в таких системах уравнений может достигать десятков миллиардов. Один из подходов к решению такой системы – свести решение системы к решению оптимизационной задачи. В такой задаче оптимизируемое выражение может состоять из двух слагаемых – суммарной невязки по проекциям и регуляризирующего члена, описывающего априорные знания о восстанавливаемом объекте. На примере использования измеренных в полихроматической моде проекций в работе, продемонстрировано влияние регуляризации на результат реконструкции. Измерения проводились на томографическом комплексе, созданном и функционирующем во ФНИЦ “Кристаллография и фотоника” РАН.

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

Пусть объект зондируется полихроматическим излучением, спектр которого описывается функцией ${{I}_{0}}\left( E \right).$ В модели используется приближение параллельной схемы измерения (рис. 1). Такое приближение правомерно, если расстояние источник – образец достаточно велико. В эксперименте, данные которого проанализированы в статье, расстояние источник – образец 1.2 м, включая распространение излучения по вакуумному пути. Рассмотрим 2D-задачу, т.е. для одного из горизонтальных сечений объекта. Система “источник–детектор” неподвижна. Детектор позиционно-чувствительный. Объект укреплен на гониометре. Величина $\varphi $ описывает текущий проекционный угол, $r$ определяет ячейку детектора, для которой будет записано выражение, связывающее величину регистрируемого сигнала и распределение коэффициента ослабления (рис. 1).

Согласно закону Бугера–Ламберта–Бера наблюдается экспоненциальное ослабление зондирующего излучения при его прохождении через объект:

Здесь $I\left( {\varphi ,r} \right)$ – это регистрируемое ячейкой детектора излучение, прошедшее через объект, повернутый на угол $\varphi .$ Положение ячейки детектора задает $r.$ описывает траекторию регистрируемого ячейкой луча. $\mu \left( {l,E} \right)$ – коэффициент линейного ослабления объектом излучения на энергии $E$ в текущей точке $l$ на луче $L\left( {r,\varphi } \right)$. Обратим внимание, что записанное выражение содержит набор неизвестных распределений коэффициентов линейного ослабления, т.е. каждой энергии энергетического спектра $E$ соответствует свое распределение $\mu .$ Чтобы применить алгоритмы реконструкции, созданные математиками для восстановления изображений в монохроматической моде (т.е. для случаев, когда для зондирования используется только одна энергия), необходимо скорректировать зарегистрированные проекции. Мы применили метод гамма-коррекции проекций, в котором используется концепция радоновского инварианта для автоматического выбора параметра коррекции [3]. После выполнения процедуры коррекции над проекциями $\quad$обратная задача формулируется следующим образом: требуется восстановить пространственное распределение $\mu {\text{*,}}$ используя выражение:

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

ИСПОЛЬЗОВАНИЕ РЕГУЛЯРИЗИРУЮЩИХ ПРОЦЕДУР

При использовании алгебраических методов дискретизация проводится на первом шаге. Восстанавливаемый объем опишем воксельной моделью. Перепишем в матричной форме выражение (2):

где A – так называемая проекционная матрица, описывающая лучи, в направлении которых проводится суммирование, x – вектор восстанавливаемых значений $\mu {\text{*}}$ объема, y – вектор лучевых сумм с компонентами $\ln \left( {\frac{{\int {dE{{I}_{0}}\left( E \right)} }}{{I{{{\left( {\varphi ,r} \right)}}^{\gamma }}}}} \right).$ Для решения системы уравнений (3) используется метод численной оптимизации. Минимизируемый функционал имеет вид:

Здесь второе слагаемое представляет линейную комбинацию регуляризирующих членов, учитывающих свойства изображения. Обратим внимание, что для изменения вклада невязки производится взвешивание слагаемых с коэффициентами $\alpha $ и $1 - \alpha ,$ где $\alpha $ лежит в диапазоне от 0 до 1. Мы решаем задачу в такой постановке, поскольку допускаем, что так называемая идеальная геометрия, заложенная в модели формирования лучевой суммы, может отличаться от геометрии проведенного измерения, например, из-за погрешностей юстировки. Для решения оптимизационной задачи был использован метод градиентного спуска.

Такой подход позволяет использовать априорные знания об используемом объекте. В нашем случае это была информация о том, что томографируемый объект содержал области сильного поглощения, т.е. в наборе сигналов, регистрируемых пикселями детектора, содержались такие сигналы, величина которых была на уровне шума. Оценить поглощение на траекториях лучей, приходящих в эти пиксели, невозможно, однако известно, что оно больше предела чувствительности регистрирующего прибора. Для таких пикселей в системе уравнений (2) равенства были заменены на неравенства. И задача решалась в постановке с мягкими ограничениями [4].

ОБСУЖДЕНИЕ РЕЗУЛЬТАТОВ

Объектом томографического исследования, для которого проводилось сравнение процедур реконструкции, был молочный зуб с нарушенным покровом эмали, физиологическое рассасывание корней которого не завершено. Не рассосавшиеся корни зуба образовали внутреннюю полость, в которую была помещено металлическое свинцовое включение, моделирующее сильно поглощающую пломбу или имплант. Фотография зуба представлена на рис. 2.

Параметры эксперимента следующие. Вакуумный путь 1.2 м от рентгеновской трубки с Mo-анодом (40 кэВ, 20 мА) до образца, расстояние образец–детектор 0.05 м. Время измерения 5 с на кадр. Собрано 400 проекций с шагом поворота объекта угла 0.5 градуса. Детектор XIMEA-xiRay 11 Mpix использован для регистрации излучения.

Результат реконструкции полного объема методом SIRT [5] приведен на рис. 3.

В области сильнопоглощающего включения наблюдаются артефакты, известные как “metal-like”. Исследования, целью которых является уменьшение такого рода артефактов, не прекращаются и сегодня [6]. Ниже представлены несколько горизонтальных сечений восстановленного разными методами объема. На рис. 4 приведены результаты реконструкции разными методами. Реконструкция проводилась методом FBP (Filtered Back Projection) [2] (рис. 4а), методом SIRT (Simultaneous Iterative Reconstruction Technique) [5] (рис. 4б), методом SIRT с включенной процедурой регуляризации по первой норме градиента (рис. 4в), методом SIRT в постановке задачи с мягкими ограничениями (рис. 4г) и методом SIRT в постановке задачи с мягкими ограничениями и включенной процедурой регуляризации по первой норме градиента (рис. 4д).

Хотелось бы обратить внимание, что методы FBP (рис. 4a) и SIRT без учета мягких ограничений (рис. 4б и 4в) порождают яркую границу вокруг сильно поглощающего включения. Решение задачи в условиях мягких ограничений (рис. 4г и 4д) выравнивает плотность включения. Добавление к мягким ограничениям регуляризирующего члена (рис. 4д) проявляет и границы полости, в которую было добавлено включение.

Все численные реализации использованных методов выполнены авторами данной статьи.

ЗАКЛЮЧЕНИЕ

В работе представлены результаты томографического эксперимента, проведенного на томографическом комплексе, созданном и функционирующем во ФНИЦ “Кристаллография и фотоника” РАН. Зондировался молочный зуб с нарушенным покровом эмали, физиологическое рассасывание корней которого не завершено. Не рассосавшиеся корни зуба образовали внутреннюю полость, в которую была помещена частица Pb, моделирующая сильно поглощающее включение. Обсуждаются результаты реконструкции с использованием интегрального и алгебраического подходов. Показано, что использование регуляризирующих процедур в алгебраическом методе реконструкции позволяет повысить качество восстанавливаемого изображения даже в случае наличия в зондируемом объекте одновременно областей высокой и низкой оптической плотности.

Исследование влияния рентгенооптических свойств объектов при полихроматическом зондировании на динамику изменения параметров процедур предобработки томографических проекций и анализ восстановленных изображений проводились при поддержке Министерства науки и высшего образования в рамках выполнения работ по Государственному заданию ФНИЦ “Кристаллография и фотоника” РАН. Исследования, связанные с выбором параметров процедуры регуляризации в зависимости от локальных рентгенооптических свойств зондируемых объектов, выполнялись при частичной финансовой поддержке РФФИ (проект № 17-29-03492).