Хабр, жаркий августовский привет тебе от отдела компьютерной томографии компании Smart Engines! Раньше мы тебе рассказывали о задаче поиска положения смещенной и наклоненной оси вращения объекта в компьютерной томографии. Мы обещали рассказать о нами разработанном методе решения этой задачи, и вот, мы здесь! Мы вернулись к тебе с опубликованной статьей о нашем методе и с полученным патентом РФ!

Введение

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

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

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

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

Побороть часть недостатков существующих методов удалось нашему отделу компьютерной томографии – мы разработали быстрый и устойчивый к шумам автоматический метод определения положения оси вращения объекта исключительно по набору его проекционных данных, полученных в параллельной или конусной моделях распространения рентгеновских лучей. Метод был нами запатентован в РФ, описание метода и исследование его свойств опубликовано в международном научном журнале. Ниже мы и хотели бы немного рассказать нашему читателю о разработанном методе и показать, как он нам не раз помогал улучшить качество реконструкции в Smart Tomo Engine!

Описание метода

Предложенный метод решает задачу определения координат положения оси вращения объекта: ее сдвига x вдоль детектора и угла наклона вдоль детектора (в предположении отсутствия наклона оси ортогонально плоскости детектора). Входные данные метода – набор проекций, полученных при обращении объекта не обязательно на полный круг. Проекционные данные могут быть получены как в параллельной, так и в конусной моделях распространения рентгеновских лучей. Метод устойчив к шумам в проекционных данных и к небольшим выходам объекта из поля видимости детектора. Он не требует дополнительной настройки пользователем и является полностью автоматическим.

Ключевой инструмент нашего метода автоматического определения параметров положения оси вращения – (арифметически) усредненное проекционное изображение (сокращенно, УПИ) :

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

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

Важное свойство УПИ – осесимметричность (рис. 3). Мы предполагаем, что УПИ – осесимметричное изображение, причем его ось симметрии совпадает с проекцией оси вращения объекта на плоскость окна детектора при просвечивании рентгеновским излучением.  Именно свойство симметричности УПИ позволяет извлечь из него информацию об истинном расположении оси вращения. Оно дает возможность ассоциировать ось симметрии УПИ с осью вращения объекта.

Можно уточнить, что, имея в виду круговую схему сканирования с конусной моделью распространения рентгеновских лучей, УПИ на самом деле является осесимметричным изображением в случае отсутствия сдвигов оси вращения. В случае наличия сдвигов УПИ можно считать приближенно осесимметричным изображением: несимметричность УПИ становится заметной при сдвигах, превышающих несколько десятков-сотен пикселей детектора. Это мы экспериментально проверили. В параллельной геометрической схеме при рассмотрении множества допустимых движений оси вращения в плоскости, параллельной плоскости окна детектора, УПИ всегда строго осесимметрично.

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

В предположении осесимметричности ось симметрии УПИ совпадает с проекцией оси вращения объекта на плоскость окна детектора (рис. 4): величина угла наклона оси в плоскости, параллельной плоскости окна детектора, равна величине угла между осью симметрии УПИ и вертикалью, а сдвиг оси вращения отличается от сдвига оси симметрии УПИ от центра изображения в количество раз, задаваемое коэффициентом увеличения конусной геометрической схемы, или совпадает со сдвигом оси симметрии УПИ в параллельной геометрической схеме. Таким образом, задача определения скалярных параметров положения оси вращения объекта – величин сдвига и наклона в плоскости, параллельной плоскости окна детектора, – сведена к задаче нахождения уравнения оси симметрии УПИ, в частности, ее угла наклона и сдвига в горизонтальном направлении (от центра изображения).

Как можно найти уравнение оси симметрии УПИ, а точнее величины ее угла наклона и сдвига ? Заметим, что искомые величины и таковы, что при повороте УПИ на градусов и последующем сдвиге на пикселей полученное трансформированное УПИ имеет строго вертикальную ось симметрии, проходящую через центр изображения. Следовательно, наша цель – действиями группы сдвигов и поворотов преобразовать УПИ в изображение с вертикальной осью симметрии, проходящей через центр изображения. Или, другими словами, наша цель состоит в минимизации величины – функции оценки несимметрии УПИ объекта относительно прямой, образующей угол наклона с вертикалью и сдвинутой на пикселей от центра изображения:

где – оператор поворота изображения на градусов вокруг его центра, – оператор зеркального отражения изображения в горизонтальном направлении, – оператор сдвига изображения на в горизонтальном направлении, – операция композиции отображений, а – стандартная -норма Минковского. Вместо нормы Минковского может быть вычислена иная матричная норма.

Ось, наклоненная на градусов и сдвинутая горизонтально относительно центра изображения на пикселей, является осью симметрии изображения , если при повороте на градусов, горизонтальном отзеркаливании и последующем сдвиге на пикселей также получим изображение , повернутое на градусов. Функция оценки несимметрии изображения относительно оси, повернутой на градусов и сдвинутой от центра на пикселей, есть – норма разности между повернутым на градусов изображением и его зеркально отраженной и сдвинутой на пикселей копией . Чем меньше значение , тем более уверенными мы можем быть в том, что – ось симметрии изображения (отсчитывая от центра изображения). Тем самым, действительно, искомые величины угла наклона и сдвига оси симметрии характеризуется тем, что пара параметров доставляет минимум функции (рис. 5). Значения параметров положения оси вращения, - угол между осью и вертикалью и сдвиг в плоскости, параллельной плоскости окна детектора, - напрямую выражаются через и : значение угла наклона оси равно , а значение сдвига оси отличается от в количество раз, равное коэффициенту увеличения геометрической схемы. Итак, задача определения параметров положения оси вращения сведена в свою очередь к минимизации функционала .

Минимизация функционала в нашем методе производится каскадным образом: перебираются значения параметра из нескольких измельченных равномерных сеток значений (вокруг значения 0), для каждого из которых перебираются значения по некоторой равномерной сетке (вокруг значения 0) значений с последующим вычислением ; после каждой такой итерации расчетные сетки измельчаются и центрируются, то есть помещаются так, что их центром являются значения и , при которых был достигнут минимум целевого функционала на предыдущей итерации.  Искомые значения параметров положения оси вращения и – доставляют минимум функционала . 

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

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

Результаты

Наш метод был неоднократно протестирован как на синтетических, так и на реальных данных. Результаты работы метода на отдельных пакетах синтетических и реальных данных представлены на рис. 7. На рисунках изображены реконструкции центральных слоев синтетических данных (фантом с геометрическими фигурами) и реальных данных (тубус с шахматами) (сравните с реконструкциями без коррекции, помещенными слева на рис. 7). Для реконструкции использовался алгоритм Feldkamp-Davis-Kress.

Синтетические проекционные снимки (пакет 1024x1024x720) симулированы с различными значениями угла наклона и сдвига оси вращения (табл. 1). На рисунке 7 – синтетические данные симулированы с углом наклона оси вращения величины 9.20 градуса в плоскости, параллельной плоскости окна детектора, и сдвигом величины 2.83 пикселя. Результаты работы нашего метода на наборах синтетических проекционных снимков оформлены в виде таблицы (табл. 1). Для большинства пакетов синтетических данных абсолютная погрешность каждого из параметров положения оси вращения составляет не более нескольких сотых долей градуса или пикселя. Так, например, для набора данных на рис. 7 найденное предложенным методом значение угла наклона оси составляет 9.20 градуса, найденное значение сдвига оси равно 2.81 пикселя (абсолютная погрешность 0.02 пикселя). С результатами других экспериментов читатель может ознакомиться в работе.

Разработанный нами метод спас качество не одной томографической реконструкции (А вы видели наш youtube-канал? Для почти каждой из реконструкций мы использовали изложенный метод определения оси). И наш метод за борьбу с артефактами в реальных экспериментах уже не раз был представлен к награде. Среди его подвигов, к примеру, коррекция артефактов смещенной оси на пакете реальных данных тубуса с шахматами (пакет 2064x1548x1024): на рис. 7 визуально очевидно улучшение качества реконструкции после применения метода. Найденная величина угла наклона оси вращения равна 0.33 градуса, найденная величина сдвига оси составляет -7.09. Ручной подбор параметров положения оси вращения подтверждает точность результатов метода: оптимальное значение сдвига между -7.00 и -8.00 пикселями, в то время как оптимальный угол наклона - между 0.00 и 0.50 градусов.

А вот случай, в котором наш метод помог побороть беду, а именно починить фонарик! И все благодаря предложенному нами методу! 

До встречи с методом коррекции оси было непонятно, почему фонарик не работает (рис. 8). Казалось бы, все соединения исправны, и фонарик должен гореть (рис. 9)!

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

Заключение

Преимущество метода в том, что он избегает попарного сравнения множества отдельно взятых проекционных снимков (вместе с этим попарное сравнение – ядро ряда классических методов). На УПИ случайные шумы не так явно выражены, как на отдельных проекциях, – шумы усредняются в результате суммирования всех проекций и становятся менее выраженными на усредненном снимке. Последнее выступает аргументом в пользу использования именно УПИ вместо единичных, иногда очень шумных, проекционных снимков. Вместе с отмеченным, не менее важен тот факт, что предложенный метод применим не к единственной параллельной геометрической схеме: он может использован как в параллельной, так и конусной геометрической схемах.

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

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