CC BY
42
ОПРЕДЕЛЕНИЕ ПАРАМЕТРОВ РУСЛА СОСУДОВ С ИСПОЛЬЗОВАНИЕМ ТРЕХМЕРНОГО ЛОКАЛЬНОГО ВЕЕРНОГО ПРЕОБРАЗОВАНИЯ
А.В. Куприянов, Н.Ю. Ильясова, А.Г. Храмов, А.М. Малафеев, О.А. Титова Институт систем обработки изображений РАН, Самарский государственный аэрокосмический университет
Аннотация
В статье рассматривается трехмерное локальное веерное преобразование в задаче выделения структуры древовидного объекта. Данное преобразование является модификацией лучевого преобразования, которое приводится к виду, удобному для анализа древовидных структур. На основе предложенного преобразования разрабатывается метод определения направления русла ветви сосуда. Также в статье описывается метод выделения структуры сосуда. Проводятся экспериментальные исследования предложенных методов и алгоритмов.
Введение
Одна из возможных задач при анализе изображений - это обработка изображений объектов, характеризующихся наличием древовидных структур.
Древовидная структура (рис. 1) характеризуется наличием вершин и ребер, соединяющих эти вершины.
Чтобы избежать этого, необходимо по рентгеновским снимкам выделить структуру исходного объекта. Для этого нужно найти разветвления и направления русла сосуда. Предварительно по имеющимся изображениям проекций восстанавливается пространственное распределение интенсивности /1/. При восстановлении используется лучевое преобразование. Затем на основе анализа пространственной интенсивности выделяется структура древовидного объекта. Одной из важнейших задач при этом явля-
В качестве вершин выступают узлы, начальная точка сосуда и конечные точки его ветвей. Начальная точка или корень - точка, откуда на изображении начинается сосуд. Будем называть узлом ту часть сосуда, где происходит ветвление. Ветвью называется множество точек сосуда, соединяющих вершины. Конечная точка - точка, где ветвь заканчивается.
Примером таких изображений могут быть изображения различного рода трещин, русла рек, сеть дорог и другие. Важным примером такой структуры является сердечно-сосудистая система. Для диагностирования ее состояния получают рентгеновские проекции сердца, по которым специалист ставит диагноз. Однако оценка врача носит субъективный характер и может быть ошибочной, так как он судит о трехмерном объекте по двухмерным проекциям (рис. 2), то есть часть информации о структуре может быть потеряна.
-у >
i
\
ется определение направлений русла ветвей. Предложенный в статье метод основан на методе, описанном в статье /2/, и является его расширением для трехмерного случая.
Трехмерное локальное веерное преобразование
Основной задачей, возникающей при формировании структуры древовидного объекта, является определение направления в точке ветви.
Рассмотрим метод определения направления в точке, основанный на модификации лучевого преобразования /3/. Лучевое преобразование позволяет проводить эффективный анализ пространственного распределения интенсивности. Модифицируем преобразование так, чтобы алгоритм хорошо работал в местах разветвления.
Запишем лучевое преобразование (рис. 3а), с помощью которого проводится анализ окрестности выбранной точки, в общем случае:
P( х,в) =
| / ^
+ гв)йг.
(1)
у©
у®
\
а)
Рис. 3. Модификации преобразования.
Это преобразование представляет собой интеграл функции / по прямой, проходящей через точку
х в направлении в. Здесь / - это восстановленное
пространственное распределение интенсивности, х -заданная точка, в которой определяются направления ветвления, а направление в задается углами.
Предлагается ограничить область интегрирования одним направлением (рис. 3Ь), то есть интегрировать по лучу. Действительно, если в точке посчитать лучевые преобразования по разным направлениям и выделить одно, которое проходит по ветви, то невозможно точно определить, в какую из двух противоположных сторон от точки по прямой интегрирования идет ветвь. Веерное преобразование представляет собой интеграл функции / по лучу с
Ь) с)
а) лучевое, Ь) веерное, с) локальное веерное
Данное преобразование будем называть локальным веерным преобразованием (ЛВП). Оно представляет собой интеграл функции / по лучу с
началом в точке х и направлением в внутри сферической окрестности радиусом г.
Дискретизация переменных ЛВП Так как обработка данных с помощью вычислительных машин имеет дискретный характер, происходит дискретизация переменных, являющихся аргументами полученного преобразования.
Имеется пространственное распределение интенсивности I (х, у, 2), где (х, у, 2) - координаты точки в пространстве. Координаты являются дис-
I = Г7 т = 1М
началом в точке х и направлением
Б(х,в) = |/(х + .
(2)
Дополнительно ограничим область интегрирования относительно небольшой окрестностью, размеры которой следует выбирать в зависимости от кривизны ветвей так, чтобы внутри выбранной окрестности ветвь можно было достаточно хорошо аппроксимировать отрезком прямой (рис. 3с). Это позволяет точнее определить направление ветвления, так как в случае совпадения направления интегрирования с ветвью вся область интегрирования будет лежать целиком на ветви.
Кроме того, в окрестности разветвлений по направлению ветви в результате преобразования получится более выраженный пик, чем для лучевого преобразования, так как интегрирование не будет производиться по фону изображения. Пусть выбранная окрестность имеет радиус г . Тогда формула (2) перепишется следующим образом:
г
^ (~х,в, г) = | / (х + в)Ж.
(3)
где
кретными: х1, ут, 2
п = 1, N. Здесь Ь, М, N - размеры области, ограничивающей сосуд.
Направление, в котором проводится интегрирование, задается вектором в. Суммирование значений интенсивности происходит по этому направлению с увеличением расстояния от точки х , для которой производится ЛВП, то есть по лучу. Таким образом, естественно будет перейти к сферической системе координат. Тогда положение любой точки в пространстве будет определяться сферическими координатами (в, р, г), где в - угол между проекцией вектора, соединяющего заданную точку с началом координат, на плоскость ох2 и осью 2 , р -
угол между проекцией вектора, соединяющего заданную точку с началом координат, на плоскость оху и осью х , а г - радиус, то есть расстояние, на котором от центра координат взята точка (в, р, г). В
нашем случае пространственное распределение интенсивности было восстановлено для координат в декартовой системе. Поэтому, при интегрировании, когда определяются сферические координаты точки, в которой нужно найти интенсивность, необходимо
перейти от сферических координат к декартовым. Формулы перехода выглядят так:
X = х0 + Г С05(р) * 5Ш(в),
•у = Уо + Г вт(р) * 5т(в), (4)
2 = 2о + Г СОв(#),
где (х0, у0,2о) - координаты точки, в которой проводится ЛВП.
При вышеуказанном переходе получаемые координаты не всегда являются целочисленными, а значения интенсивности есть только для целых значений координат. На данном этапе проводится интерполяция.
В простейшем случае проводится интерполяция по ближайшему соседу. Можно также попытаться повысить точность алгоритма путем применения более сложной схемы интерполяции.
Теперь преобразование (3) можно переписать в следующем виде:
Г
^ (х;в,р; г ) = v / *( хо + г со5(р) 5ш(в),
<=0 (Э)
у0 + < 5Ш(р)5Ш(в), 20 + < С05(в)),
где / (х, у, 2) - значение интенсивности в точке с координатами, найденными в результате интерполяции. Преобразование проводится в точке х = (х0, у0,20) в направлении, определяемом углами (в,р).
Зададим шаги интегрирования и Д 2, соответственно, для углов в и р. Эти шаги зависят от требуемой точности алгоритма и выбираются либо априорно, исходя из соображений об исследуемом объекте, либо определяются в процессе работы.
Обозначим число направлений по углу в К1,
а по р - К2. Тогда преобразование (5) примет вид:
Г
^(х; в1 ,р ; г) = v /*(хо + < соэ(р,) 5ш(в,.),
(6)
уо + < sm(рj )5Ш(в,.), 2 о + < СОБ(в,.)),
где I = 1,Кх, ] = 1,К2; х = (хо,уо,2о).
Преобразование, приведенное в формуле (6), представляет собой трехмерное локальное веерное преобразование.
Алгоритм определения ветвления сосуда
Предлагается применить разработанный метод, чтобы определить направления, вдоль которых происходит ветвление сосуда или идет продолжение сосуда. Для этого в каждой точке считается ЛВП, по полученным значениям определяются направления ветвления, а затем используется алгоритм трехмерной трассировки, предложенный в статье [4]. Суть этого алгоритма в том, что на каждом шаге определяются направления «оптимального движения», в которых предположительно сосуд продолжается, и затем сосуд отслеживается в найденных направле-
ниях. Легко видеть, что недостатком такого подхода является большой объем вычислений, но в результате работы алгоритма получается структура исследуемого объекта.
Для того, чтобы определить направление в точке через, которую проходит сосуд или направление ветвления, (если сосуд ветвится), используется полученное преобразование. Считается ЛВП в зафиксированной точке во всех направлениях (рис. 4а, изображение матрицы значений ЛВП) для заданного значения радиуса ограничивающей окрестности. Изображение матрицы представляет собой развертку поверхности сферы, на которой точками отмечены два направления обнаруженных в соответствующей точке пространства.
• I
а) Ь)
Рис. 4. Примеры набора значений ЛВП: а) изображение матрицы значений ЛВП, Ь) направления русла сосуда, найденные после выделения минимумов в матрице
Известно, что на рентгеновских проекциях точки сосуда имеют меньшую яркость, чем точки фона, так что восстановленная интенсивность также темнее там, где находится объект. Поэтому в тех направлениях, где сосуд продолжается, посчитанное преобразование будет меньше, чем в тех, где присутствует фон. Для того чтобы выделить направления ветвления, набор значений ЛВП представляется как двумерная матрица М. ], индексы которой
представляют собой углы, задающие направление в пространстве. Значение матрицы М представляет
собой значение ЛВП, посчитанного в направлении (5, V). В полученной матрице находим минимумы.
Однако сосуд имеет некоторый диаметр, поэтому минимумы, определяющие близкие направления, на самом деле соответствуют одной ветви сосуда. Поэтому после выделения минимумов нужно произвести сегментацию изображения, полученного после выделения минимумов изображения, найти центры полученных областей (рис. 4Ь), которые определяют направления русла сосуда в трехмерном пространстве. Полученные значения углов определяют, в каких направлениях в заданной точке продолжается или ветвится сосуд.
Предложенные методы и алгоритмы были исследованы на изображениях проекций сердечнососудистой системы. Сначала по имеющимся проекциям было восстановлено пространственное рас-
пределение интенсивности. Затем в узле были посчитаны направления ветвления сосуда (рис. 5).
sl=(161°, 56%, s2=(72°, 79%, s3=(86°, 179%
Рис. 5. Изображение восстановленного пространственного распределения интенсивности с найденными направлениями русла сосуда
Заключение
Был разработан метод, позволяющий определить направления ветвей древовидного объекта, используя ЛВП. В [5] решается аналогичная задача. Авторы прослеживают указанную ветвь с помощью трассировки. Недостатками предлагаемого метода можно назвать необходимость явно указывать начальную и конечную точки ветви и высокую чувствительность к шумам.
Предлагаемый в настоящей статье метод не требует указания ветви и, как показали исследования, является устойчивым к шумам.
В качестве области применения можно назвать диагностирование состояния сердечно-сосудистой системы по рентгеновским снимкам. В настоящее время разработанные алгоритмы внедряются в компьютерную систему диагностики заболеваний сердца. Дальнейшее развитие работы связано с разработкой эффективных алгоритмов локального веерного преобразования. В частности, желательно сократить время работы алгоритма.
Благодарность
Работа выполнена при поддержке Министерства образования РФ, Администрации Самарской области и Американского фонда гражданских исследований и развития (CRDF Project SA-014-02) в рамках российско-американской программы "Фундаментальные исследования и высшее образование" («BRHE»), а также Российского фонда фундаментальных исследований (грант № 03-01-00642) .
Литература
1. Бейтс Р., Мак-Доннел М. Восстановление и реконструкция изображений // Пер.с англ. М.: Мир,
1989. 336 с.
2. Баранов В.Г., Храмов А.Г. Дискретное веерное преобразование Радона в задаче выделения центров ветвей сетчатых структур // Компьютерная оптика, М., МЦНТИ, 2001.В 23.
3. Наттерер Ф. Математические аспекты компьютерной томографии // Пер. с англ. М.: Мир,
1990. 288 с.
4. Ильясова Н. Ю., Корепанов А.О., Куприянов А.В., Баранов В.Г., Храмов А.Г. Анализ структуры сосудистой системы сердца методом трассировки изображений проекций // Компьютерная оптика, М., МЦНТИ, 2001. В 23.
5. Branchevsky S.L., Durasov A.B., Ilias-ova N.Yu., Ustinov A.V. Methods for estimating geometric parameters of retinal vessels using diagnostic images of fundus // Proceedings SPIE, 1998. Vol. 3348. P. 316-325.
