Оптическая когерентная томография (ОКТ) – неинвазивная методика исследования внутренней структуры объекта, основанная на принципах низкокогерентной интерферометрии. ОКТ нашла свое применение во многих прикладных и научно-исследовательских задачах. Использование оптического излучения ближнего инфракрасного диапазона (длина волны источника, как правило, составляет 700–1500 нм) с малой когерентностью позволяет визуализировать внутреннюю структуру исследуемого образца с высоким пространственным разрешением. В связи с этим одной из наиболее перспективных областей использования ОКТ является биомедицинская диагностика, где в качестве основных приложений данной техники можно выделить офтальмологию, дерматологию и кардиологию [1, 7, 8].
При исследовании верхних слоев кожи пространственное аксиальное разрешение составляет 3–15 мкм с глубиной когерентного зондирования 1–2 мм, что обусловлено высоким значением показателя анизотропии g ~ 0,9–0,99 и коэффициента рассеяния среды μs = 70–200 см–1 [6]. С увеличением глубины залегания объекта получение его достоверного структурного изображения с помощью данной методики затрудняется в связи с многократным рассеиванием оптического излучения, появлением спекл-шумов, особенностями применяемого источника излучения и т.д. Для усовершенствования современных ОКТ систем необходимо полное понимание технических особенностей данной методики и закономерностей, лежащих в основе взаимодействия оптического излучения с биологическими тканями.
На данный момент методика моделирования транспорта фотонов на основе метода Монте-Карло является стандартным способ исследования взаимодействия оптического излучения и биологических тканей. Помимо фотонного транспорта, методика Монте-Карло используется также во многих других алгоритмах моделирования физических процессов благодаря ее высокой точности [2, 5]. Часто используемый алгоритм и программное обеспечение (MCML), реализующее данную методику, можно использовать для источников с высокой и с низкой когерентностью [9]. Данный алгоритм позволяет проводить моделирование фотонного транспорта в стационарном режиме в многослойной полубесконечной среде, в то время как реальные биологические объекты, как правило, имеют более сложную неоднородную пространственную структуру. Также MCML является довольно медленным приложением, работающим в однопоточном режиме, что не позволяет использовать потенциал современных многоядерных процессоров.
Алгоритм, лежащий в основе MCML, может быть использован для создания методики моделирования ОКТ [11]. Недостатком данной методики является то, что границы слоев задаются в аналитическом виде, что позволяет описать их структуру только в приближенном виде.
При дерматологических исследованиях особый интерес представляет возможность визуализации подкожных структур, таких как кровеносные сосуды [3, 4]. В некоторых случаях визуализация таких объектов затрудняется в связи с крайне высоким показателем рассеяния крови μs = 600–1000 см–1. Теоретические возможности ОКТ при исследовании подобных структур могут быть проверены путем соответствующего моделирования.
Целью данной работы является создание методики моделирования изображений оптической когерентной томографии, позволяющей учитывать спекл-структуру В-скана, особенности биологических объектов, содержащих кровеносные сосуды, и высокие значения показателя рассеяния.
Материалы и методы исследования
Для того чтобы получить более детальное описание структуры, предлагается перейти к вокселной модели формирования среды для моделирования и к соответствующему характеру распространения фотонов в подобной среде. При таком подходе объект задается в пространстве в качестве набора трехмерных элементов (вокселей), каждый из которых имеет свои соответствующие оптические свойства (коэффициент рассеяния μs, коэффициент поглощения μa, показатель анизотропии g и показатель преломления n). Размер ребра вокселя при моделировании определяется пространственным разрешением используемой ОКТ системы, что отражает по сути минимальный размер объекта, который теоретически можно визуализировать с помощью данной методики, и, как правило, составляет несколько микрометров. Общее число вокселей по какому-либо направлению N определяется исходя из размеров исследуемого образца l по этому направлению:
где d – размер ребра вокселя. В общем случае количество и длина ребер вокселей по осям x и y могут отличаться от соответствующих значений по оси z. При моделировании in vivo исследований размеры образца по оси z определяются теоретической глубиной когерентного зондирования методики. Важным является также задание объектов, лежащих ниже данной глубины, так как интенсивность многократно рассеянных фотонов также вносит определенный вклад в формирование сигнала (паразитная составляющая ОКТ сигнала).
Формирование трехмерного массива оптических свойств, представляющих исследуемый объект может производиться как в аналитическом виде, так и непосредственно с помощью экспериментально полученных ОКТ изображений внутренней структуры данного объекта. Во втором случае элементам массива оптических свойств, которые используются при моделировании, ставятся в соответствие пиксели структурного ОКТ изображения, с подчеркнутыми контурами, на котором выделяются определенные зоны с характерными оптическими свойствами (кровь, эпидермис, дерма и т.д.), что может быть выполнено с помощью специализированного программного алгоритма идентификации слоев на структурном ОКТ изображении. При моделировании на основе единичного B-скана распределение оптических характеристик по всем сечениям на оси y в упрощенном виде считается таким же, как и на этом B-скане (например, для показателя преломления 
), либо корректируется в соответствии с предполагаемой структурой объекта.
Дистанция свободного пробега фотонного пакета определяется согласно следующей формуле:
где μt = μa + μs – полный коэффициент взаимодействия воксела, в котором находится фотон, а ξ – случайное число, равномерно распределенное между 0 и 1 [9].
Если показатели преломления текущего воксела и воксела, в который переходит фотонный пакет, различны, то вероятность пересечения границы определяется с помощью формул Френеля:
где αi – угол падения, а αt – угол преломления. Значение R затем сравнивается со случайной величиной ξ ∈ [0, 1), и если ξ > R, то фотон переходит в следующий сегмент, а если ξ ≤ R, то фотон отражается от границы.
Интенсивность излучения I(z), принимаемая датчиком при моделировании, описывается с помощью следующей зависимости [11]:
где I0 – константа, определяемая ОКТ системой; W – статистический вес фотонного пакета, выходящий из ткани внутри области нахождения датчика; Li – оптический путь, пройденный фотонным пакетом; λ – длина волны источника; z – оптический путь в опорном плече интерферометра; lc – длина когерентности:
где Δλ – полная ширина на уровне половинной амплитуды (FWHM). Множитель 
представляет спекл-структуру ОКТ изображения, который может быть удален, если она не учитывается при моделировании. Для проверки адекватности модели получаемые результаты могут сравниваться с экспериментально полученными А-сканами и структурными изображениями (В-сканы).
Представленная модель реализована в виде многопоточного приложения с помощью языка C#. Количество потоков определяется согласно общему количеству ядер центрального процессора. Для расчетов использовался четырехъядерный процессор Intel Core i5-4670. Каждый поток производит вычисления одинакового числа фотонов и заданного количества А-сканов. Разница процедуры вычислений на разных потоках определяется тем, что генератор случайных чисел, используемый при моделировании (стандартная функция Random.NextDouble()), принимает различные значения ключа начального состояния генерации в зависимости от времени активации потока.
Результаты исследования и их обсуждение
Структура объекта, используемого при моделировании, представлена на рис. 1. Оптические свойства слоев, используемые в моделировании, приведены в таблице. Они соответствуют экспериментальным значениям для источника излучения с длиной волны λ = 1000 нм [6].
Результаты моделирования представлены на рис. 2. Изображение, соответствующее моделированию с учетом спекл-структуры (рис. 2, а), характеризуется тем, что зона, расположения подкожного сосуда имеет несколько большую интенсивность, по сравнению с дермой, и в целом сосуд может быть с трудом локализован в точке его реального расположения, в то время как нижняя стенка сосуда полностью размывается. В случае моделирования без учета спекл-структуры (рис. 2, б) сосуд может быть непосредственно локализован в точке его расположения. При этом нижняя стенка сосуда имеет более выраженную границу, но непосредственно локализовать ее местоположение можно лишь приблизительно. В целом полученные изображения соответствуют ранее полученным экспериментальным данным, как с точки зрения отражения спекл-структуры ОКТ сигнала, так и визуализации подкожного кровеносного сосуда [4].
Рис. 1. Структура верхних слоев кожного покрова человека, используемая при моделировании: 1 – воздух; 2 – роговой слой; 3 – эпидермис; 4 – дермис; 5 – стенка сосуда; 6 – кровь
Оптические свойства слоев кожи и стенки сосуда
|
Среда |
μs (см–1) |
μa (см–1) |
g |
n |
|
Роговой слой |
300 |
0,15 |
0,95 |
1,47 |
|
Эпидермис |
90 |
0,02 |
0,85 |
1,34 |
|
Дермис |
70 |
0,07 |
0,9 |
1,34 |
|
Стенка сосуда [10] |
2,8 |
0,9 |
0,9 |
1,4 |
|
Кровь |
650 |
5 |
0,98 |
1,37 |
а б
Рис. 2. Результат моделирования структурного изображения методом Монте-Карло. Показаны верхние слои кожи и подкожного кровеносного сосуда с учетом спеклов (а) и без их учета (б)
Выводы
Предложенная методика моделирования ОКТ методом Монте-Карло на основе воксельной геометрии позволяет производить моделирование исследования среды со сложной пространственной структурой. Адекватность представленного алгоритма исследована с помощью моделирования ОКТ исследования внутренней структуры верхних слоев кожи человека, содержащих кровеносный сосуд, при этом получено хорошее соответствие результатов моделирования с ранее полученными экспериментальными данными [8]. В тоже время нерешённым остается вопрос чёткой визуализации нижней стенки сосуда. При использовании данной методики моделирования этого удалось добиться лишь приблизительно, в то время как экспериментальные данные свидетельствуют о возможности ее непосредственной визуализации в области их расположения, что будет являться объектом дальнейшей работы.