Научная статья на тему 'Оптимальный алгоритм поиска пересечений в задаче Монте-Карло моделирования распространения зондирующего излучения в головном мозге человека'

Оптимальный алгоритм поиска пересечений в задаче Монте-Карло моделирования распространения зондирующего излучения в головном мозге человека Текст научной статьи по специальности «Математика»

CC BY
581
96
i Надоели баннеры? Вы всегда можете отключить рекламу.
Ключевые слова
МОДЕЛИРОВАНИЕ РАСПРОСТРАНЕНИЯ СВЕТА / МЕТОД МОНТЕ-КАРЛО / ПОИСК ПЕРЕСЕЧЕНИЙ / KD-ДЕРЕВЬЯ / BVH-ДЕРЕВЬЯ

Аннотация научной статьи по математике, автор научной работы — Горшков Антон Валерьевич, Коршунова Анна Леонидовна

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

i Надоели баннеры? Вы всегда можете отключить рекламу.

Похожие темы научных работ по математике , автор научной работы — Горшков Антон Валерьевич, Коршунова Анна Леонидовна

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.
i Надоели баннеры? Вы всегда можете отключить рекламу.

OPTIMAL INTERSECTION SEARCH ALGORITHM IN MONTE CARLO SIMULATION OF LIGHT TRANSPORT IN HUMAN BRAIN

We consider the problem of choosing an optimal search algorithm for photon trajectory intersections with the boundaries in Monte Carlo light transport simulation in multilayer biological tissues with complex geometry based on the example of the human head. A review is given of main intersection search algorithms, their strengths and weaknesses, as well as the possibility of their adaptation to the problem. We present the results of computing experiment tests to compare the most appropriate intersection search algorithms, as well as their analysis in terms of selection of the optimal algorithm.

Текст научной работы на тему «Оптимальный алгоритм поиска пересечений в задаче Монте-Карло моделирования распространения зондирующего излучения в головном мозге человека»

Математическое моделирование. Оптимальное управление Вестник Нижегородского университета им. Н.И. Лобачевского, 2012, № 5 (2), с. 73-80

УДК 004.94

ОПТИМАЛЬНЫЙ АЛГОРИТМ ПОИСКА ПЕРЕСЕЧЕНИЙ В ЗАДАЧЕ МОНТЕ-КАРЛО МОДЕЛИРОВАНИЯ РАСПРОСТРАНЕНИЯ ЗОНДИРУЮЩЕГО ИЗЛУЧЕНИЯ В ГОЛОВНОМ МОЗГЕ ЧЕЛОВЕКА

© 2012 г. А.В. Горшков, А.Л. Коршунова

Нижегородский госуниверситет им. Н.И. Лобачевского

anton.v. gorshkov@gmail .com

Поступила в редакцию 10.09.2012

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

Ключевые слова: моделирование распространения света, метод Монте-Карло, поиск пересечений, Ы-деревья, BVH-деревья.

Введение

В настоящее время в медицинских исследованиях, в том числе предклинических, существует потребность в развитии новых неинвазивных и доступных методов диагностики, поскольку используемые традиционные методы (МРТ, КТ, ПЭТ) имеют ряд ограничений, связанных с их небезопасностью, высокими требованиями к инфраструктуре и стоимостью оборудования. Классом наиболее перспективных методов диагностики, которые могут применяться как в сочетании с существующими методами, так и в некоторых случаях вместо них, являются оптические методы. Их основными преимуществами являются неинвазивность, сравнительно невысокая стоимость приборов и широкие функциональные возможности, обусловленные возможной вариативностью параметров зондирующего излучения (длина волны, модуляция, длина импульса и т.д.).

Ряд методов оптической биомедицинской диагностики уже прошел стадию предклиниче-ских испытаний и успешно применяется в клинической практике. Для диагностики биотканей на сравнительно небольших глубинах (до единиц мм) применяются методы, использующие когерентность зондирующего излучения. Для диагностики биотканей на больших глубинах необходимо применять методы, для которых информативным является многократно рассеянное излучение. Таким методом является оптическая диффузионная томография (ОДТ).

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

Актуальность разработки методов оптической диффузионной томографии и спектроскопии подтверждается растущим в последнее время количеством публикаций на эти темы в мировой научной литературе. Основными областями применения ОДТ в настоящее время являются нейробиология, неврология и онкология.

Разработка методов применения ОДТ и ОДС в биомедицинской диагностике требует понимания процессов распространения зондирующего излучения в биологических объектах. Из-за сложной и неоднородной структуры таких объектов использование аналитических подходов к описанию этих процессов затруднительно. В этой связи эффективным выглядит использование компьютерного моделирования распространения зондирующего излучения в объектах со сложной геометрией на основе метода Монте-Карло.

Метод Монте-Карло основан на многократном повторении случайного испытания и последующем обобщении полученных данных. Применительно к распространению излучения оптического диапазона в биотканях он заключа-

ется в многократном моделировании случайных траекторий фотонов в исследуемой среде. На рис. 1 Я - диффузно отраженный фотон, А -поглощенный фотон, Т - прошедший фотон.

Рис. 1. Моделирование случайных траекторий фотонов методом Монте-Карло

Впервые применять метод Монте-Карло для моделирования распространения света в биологических тканях предложил Вильсон (Wilson) [1]. Впоследствии этот подход развивался многими учеными для получения более точных результатов. Так, Прахл (Prahl) ввел в модель учет анизотропии и внутреннего отражения света [2]. Вэнг (Wang) предложил алгоритм моделирования распространения света в многослойных гетерогенных средах с плоскопараллельной геометрией, а также разработал программную систему MCML, реализующую предложенный алгоритм [3]. Позже метод был успешно усовершенствован для моделирования оптического зондирования мозга [4]. При этом использовался воксельный подход для моделирования распространения света в голове человека. Позже этот код был адаптирован для расчета с помощью графических процессоров. Он также применялся для анализа чувствительности флуоресцентных методов имиджинга мозга человека

[5]. Недостатком этого кода является использование кубических вокселей, что не обеспечивает точного описания преломления света на границах слоев. Для преодоления этого недостатка был предложен альтернативный подход, основанный на учете сложной геометрии границ слоев и точном расчете преломления на границах [6, 7].

Моделирование распространения зондирующего излучения в биологических тканях методом Монте-Карло требует существенных вычислительных ресурсов, а в случае учета сложной геометрии границ слоев эта проблема становится еще более актуальной в силу необходимости многократного поиска пересечений фотонов с границами слоев среды. В рамках Нижегородского Центра суперкомпьютерных технологий [8, 9] была разработана программ-мная система для проведения такого моделиро-

вания для систем с распределенной памятью, где использовался простейший алгоритм поиска пересечений, основанный на полном переборе треугольников. Использование вычислительных возможностей кластера позволило существенно сократить время расчетов, но время, необходимое для проведения полезных численных экспериментов, оставалось значительным. Исследование существующего кода с использованием Intel Amplifier выявило, что именно на поиск пересечений тратится порядка 90% от общего времени моделирования.

Общее описание задачи

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

Рис. 2. МРТ головы человека (слева), триангулированная модель коры головного мозга (справа)

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

этапа: тест на пересечение луча с треугольником и выбор набора треугольников для сокращения числа таких тестов.

И если в качестве оптимального алгоритма поиска пересечения луча с треугольником большинство исследователей используют алгоритм, описанный в работе [10], то, в плане сокращения числа тестируемых треугольников, применяемые подходы различны. Так, в работе

[6] предлагается переход из трехмерного пространства в двухмерное путем проецирования триангулированной поверхности на плоскость. Плоскость разбивается на квадраты, выбираются те из них, через которые проходит траектория движения фотона, и тестируются на пересечение только те треугольники, проекции которых пересекаются с выбранными квадратами. Такой подход выглядит обоснованным в случае, когда границы слоев среды представляют собой «примерно» плоские поверхности, например, в задаче моделирования распространения излучения в верхних слоях кожи человека. Е. Маргал-ло-Балбас (Е. Маг§а11о-Ба1Ьа8) предлагает более универсальный способ поиска пересечений на основе применения окто-деревьев [7]. Однако в работе не обоснованы преимущества такого подхода по сравнению с использованием других структур данных.

Постановка задачи

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

- рассматривались биологические ткани, состоящие из нескольких слоев;

- каждый слой обладал одной или несколь-

кими границами, граница представляла собой произвольную триангулированную поверх-

ность;

- осуществлялся поиск пересечения траектории движения фотона с границами текущего слоя с учетом длины шага, на который смещается фотон в заданном направлении;

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

Для поиска пересечения луча с треугольником был использован алгоритм, описанный в работе [10].

Таким образом, целью данного исследования являлся поиск оптимального алгоритма для вы-

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

Существующие алгоритмы поиска пересечений

Задача поиска пересечения траектории фотона с границами аналогична задаче поиска пересечений луча с набором геометрических примитивов, используемой в рейтрейсинге. Отличие заключается в том, что в данном случае нужно искать только пересечения, расположенные не дальше шага свободного пробега фотона от его настоящего положения. Следовательно, возможно применение ускоряющих структур данных, использующихся для повышения производительности трассировки лучей. Как правило, такие структуры организованы в виде иерархий или ЭБ-сеток, организующих геометрические примитивы в ЭБ-пространстве. Ускоряющая структура строится для того, чтобы уменьшить число тестов пересечения лучей и примитивов модели.

Регулярная сетка. Создается однородная ЭБ-сетка, занимающая все ЭБ-пространство модели. Каждый треугольник помещается в ячейки сетки, пересекающие его. При этом во время поиска пересечений можно пробегать только по тем ячейкам, через которые прошел луч.

Достоинствами данной структуры являются:

1) быстрое построение: О(п), где п - количество треугольников;

2) простая процедура прохода луча в поисках пересечений;

3) возможность использования информации о длине шага фотона;

4) хорошая скорость на сложных участках.

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

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

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

Иерархическая сетка лучше, чем регулярная, приближает геометрию, но:

1) переход между узлами вычислительно сложен;

2) сложные регионы сходятся медленно.

Kd-деревья. Kd-дерево [11-13] представляет собой бинарное дерево ограничивающих параллелепипедов, вложенных друг в друга. Каждый нелистовой параллелепипед в kd-дереве разбивается плоскостью, перпендикулярной одной из осей координат, на два дочерних параллелепипеда. Вся геометрия целиком содержится внутри корневого параллелепипеда. Проводя рекурсивное разбиение параллелепипедов, можно добиться, чтобы в каждом листовом узле содержалось небольшое число примитивов.

Достоинства kd-деревьев:

1) временные затраты прослеживания луча: O(log N), где N - число треугольников;

2) быстрое отсечение пустого пространства;

3) простой и эффективный алгоритм поиска пересечений;

4) возможность использования информации о длине шага фотона;

5) нет перекрытий в пространстве между узлами дерева;

6) строго упорядоченное в пространстве перечисление листьев дерева вдоль луча во время поиска пересечения;

7) требуется меньше регистров при реализации для GPU.

Следует отметить, что пункты 6, 7 являются достоинствами относительно ускоряющей структуры BVH, которая будет рассматриваться далее.

Минусы kd-деревьев:

1) построение дерева с использованием SAH медленное, временная сложность алгоритма построения определена как O(Nlog N), где N -число треугольников (однако в данной задаче построение дерева выполняется всего один раз и почти не влияет на итоговую производительность алгоритма моделирования);

2) плохо предсказуемое потребление памяти, так как один примитив может находиться в нескольких листьях дерева одновременно;

3) большая глубина дерева, следовательно, больше шагов при построении и поиске пересечений.

Недостатки 2 и 3 являются недостатками относительно ускоряющей структуры BVH.

BSP-деревья. BSP-деревья (Binary Space Partition) отличаются от kd-деревьев тем, что плоскость разбиения может быть выбрана не перпендикулярно какой-либо оси координат, а произвольно. Фактически, kd-дерево представляет собой частный случай BSP-дерева. Очевидно, что увеличение количества вариантов расположения секущей плоскости ведет к усложнению алгоритма построения дерева. Дополнительно к этому усложняется алгоритм поиска пересечений. На практике, BSP-деревья не лучше для трассировки, чем kd-деревья.

BVH-деревья. Bounding Volume Hierarchy (BVH) - бинарное дерево ограничивающих объемов [14]. В зависимости от используемого ограничивающего объема выделяют несколько типов BVH:

1) Sphere tree. Со сферами просто искать пересечения, но такие деревья медленно приближают геометрию;

2) AABB-деревья используют в качестве ограничивающего объема параллелепипед, выровненный по осям координат (Axis-Aligned Bounding Box). AABB лучше приближают геометрию, чем сферы, с ними просто считать пересечения;

3) OBB-деревья состоят из не выровненных по осям координат параллелепипедов (Oriented Bounding Box), имеют лучшую адаптивность, чем AABB-деревья, но построение дерева и прослеживание луча в них становится сложнее;

4) k-DOP-деревья (Discrete Oriented Polytope) состоят из выпуклых многогранников с k гранями (AABB является 6-DOP).

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

Алгоритм построения BVH аналогичен алгоритму построения kd-дерева, только после выбора плоскости разбиения нужно ограничить каждое подмножество примитивов (слева и справа от секущей плоскости) AABB.

Достоинства BVH:

1) временные затраты прослеживания луча: O(log N), где N - число треугольников;

2) самое быстрое среди SAH-деревьев отсечение пустого пространства, поскольку узлы ограждаются от пустого пространства не отдельными секущими плоскостями, а целыми шестигранниками (в случае AABB), охваты-

вающими каждый узел;

3) возможность использования информации о длине шага фотона;

4) меньше глубина дерева (по сравнению с kd-деревьями);

5) предсказуемое потребление памяти, так как каждый примитив находится лишь в одном узле дерева.

Недостатки:

1) построение SAH-дерева медленное

(O(Nlog N), где N - количество треугольников);

2) возможны перекрытия в пространстве между узлами, не строго упорядоченное в пространстве перечисление листьев дерева вдоль луча во время поиска пересечения;

3) алгоритм поиска пересечений сложнее, чем для kd-дерева.

При построении таких ускоряющих структур, как kd-деревья и иерархии ограничивающих объемов (BVH), на каждом шаге необходимо выбирать плоскость разбиения. От способа выбора разделяющей плоскости напрямую зависит качество дерева, наилучшим критерием которого является SAH (Surface Area Heuristic)

- эвристика площадей поверхности [15]. Для повышения качества дерева, нужно на каждом шаге построения выбирать плоскость разбиения так, чтобы минимизировать SAH.

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

Так же часто используется так называемый биновый (binned SAH) алгоритм выбора плоскости разбиения [16]. Количество перебираемых плоскостей в данном алгоритме небольшое, за счет этого построение дерева происходит быстрее. Метод прост в реализации и позволяет создавать дерево высокого качества.

Выбор оптимального алгоритма

После анализа достоинств и недостатков рассмотренных выше ускоряющих структур применительно к задаче моделирования распространения зондирующего излучения, для дальнейшего экспериментального исследования были выбраны kd-деревья и BVH. Алгоритм прослеживания луча в kd-дереве проще и позволяет завершить поиск, как только найдено пересечение с примитивом (поскольку нет перекрытий между узлами). В то же время BVH, как прави-

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

Для проведения экспериментов использовались три группы тестов с геометриями, аппроксимирующими границы тканей головного мозга человека:

1) плоскопараллельная геометрия (plane) -6 плоскостей, каждая состоит из двух треугольников (всего 12 треугольников);

2) квазисферическая геометрия (sphere) -5 концентрических квазиполусфер, ограниченных плоскостью, сферы аппроксимировались разным числом треугольников. Квазиполусфера представляет собой полусферу, на поверхности которой встречаются выпуклые места и впадины (рис. 3). Такая поверхность не описывается аналитически, не является гладкой или выпуклой, и может считаться поверхностью достаточно общего вида;

3) геометрия коры головного мозга человека (cortex), полученная посредством данных МРТ (см. рис. 2), состоит из 686 348 треугольников.

Рис. З. Квазиполусфера

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

В рамках каждого теста 10 000 раз производился поиск пересечения фотона с границами слоя, текущее положение и направление движения фотона выбирались случайным образом. Тестирование проводилось с использованием алгоритмов, способных учитывать величину смещения фотона в заданном направлении (шаг фотона). Для kd-деревьев дополнительно вводились ограничения по высоте дерева (такие ограничения необходимы для ограничения объема памяти, потребляемой в процессе построения дерева).

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

- процессор Intel Core i5 2410M 2.З0 ГГц;

- 4 ГБ оперативной памяти;

- операционная система Microsoft Windows 7 x64;

- компилятор MS Visual Studio 2008.

Таблица 1

Время (в секундах) работы тестового приложения при использовании BVH

Тестовые наборы, количество треугольников Выбор разделяющей плоскости методом бинов (16 бинов) Выбор разделяющей плоскости методом бинов (32 бина)

Построение дерева Поиск без учета шага Поиск с учетом шага Построение дерева Поиск без учета шага Поиск с учетом шага

plane, 12 0 0.015 0 0 0.016 0

sphere, 4 962 0 0.078 0.046 0 0.093 0.062

sphere, 49 502 0.093 0.109 0.046 0.093 0.109 0.062

sphere, 4 995 002 13.728 0.281 0.062 14.399 0.327 0.078

cortex, 686 348 1.607 0.062 0.015 1.7 0.063 0.015

Таблица 2

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

Тестовые наборы, количество треугольников Ограничение на глубину: 15 Ограничение на глубину: 20

Построение дерева Поиск без учета шага Поиск с учетом шага Построение дерева Поиск без учета шага Поиск с учетом шага

plane, 12 0 0.16 0 0 0 0

sphere, 4 962 0.078 0.125 0.078 0.125 0.141 0.078

sphere, 49 502 0.375 0.203 0.125 0.624 0.234 0.14

sphere, 4 995 002 35.49 3.463 2.574 49.796 1.31 0.842

cortex, 686 348 3.182 0.577 0.421 4.805 0.343 0.219

Таблица 3

Время (в секундах) работы тестового приложения для к^деревьев, построенных биновым методом с количеством бинов, равным 16

Тестовые наборы, количество треугольников Ограничение на глубину: 15 Ограничение на глубину: 20

Построение дерева Поиск без учета шага Поиск с учетом шага Построение дерева Поиск без учета шага Поиск с учетом шага

plane, 12 0 0 0.015 0 0.016 0

sphere, 4 962 0.047 0.172 0.078 0.172 0.234 0.109

sphere, 49 502 0.14 0.218 0.14 0.514 0.327 0.187

sphere, 4 995 002 8.299 3.182 2.309 14.087 1.185 0.749

cortex, 686 348 1.248 0.686 0.468 2.698 0.437 0.25

Таблица 4

Время (в секундах) работы тестового приложения для к^деревьев, построенных биновым методом с количеством бинов, равным 32

Тестовые наборы, количество треугольников Ограничение на глубину: 15 Ограничение на глубину: 20

Построение дерева Поиск без учета шага Поиск с учетом шага Построение дерева Поиск без учета шага Поиск с учетом шага

plane, 12 0 0.015 0 0 0 0

sphere, 4 962 0.047 0.171 0.093 0.218 0.234 0.14

sphere, 49 502 0.156 0.234 0.14 0.639 0.343 0.187

sphere, 4 995 002 8.331 3.327 2.464 14.508 1.326 0.842

cortex, 686 348 1.263 0.624 0.436 2.98 0.452 0.249

■ BVH (16 бинов)

□ BVH (32 бина)

0KD-дерево (по границам треугольников, глубина 15) а KD-дерево (по границам треугольников, глубина 20)

□ KD-дерево (16 бинов, глубина 15)

□ KD-дерево (16 бинов, глубина 20)

□ KD-дерево (32 бинов, глубина 15)

□ KD-дерево (32 бинов, глубина 20)

Рис. 4. Сравнение реализаций алгоритма поиска пересечений с учетом шага фотона

Результаты проведенных экспериментов (рис. 4) показывают, что в задаче моделирования распространения зондирующего излучения в головном мозге человека более предпочтительным является использование БУИ-дерева, причем для его построения достаточно 16 бинов. Алгоритм на основе данной структуры на большинстве тестовых данных обеспечивает более быстрый поиск пересечений, чем алгоритм на основе кё-деревьев, а также позволяет эффективно учитывать длину свободного пробега фотона. Единственным недостатком БИУ, выявленным в рамках данных тестов, может являться большое время построения дерева. Однако построение кё-деревьев выполняется быстрее только при небольшой их высоте, что не позволяет быстро осуществлять поиск пересечений, к тому же в рассматриваемой задаче построение дерева выполняется один раз, в то время как поиск может осуществляться миллиарды раз.

Следует отметить, что, несмотря на превосходство БУИ-деревьев в данных тестах, кё-деревья могут являться им хорошей альтернативой, например, при переносе алгоритма на графические процессоры.

Заключение

В работе были исследованы алгоритмы поиска пересечений, наиболее подходящие для задачи Монте-Карло моделирования распространения зондирующего излучения в головном мозге человека. Описаны основные подходы к созданию таких алгоритмов, указаны их достоинства и недостатки с учетом поставленной за-

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

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

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

Работа выполнена при поддержке ФЦП «Научные и научно-педагогические кадры инновационной России» (проект№ 02.740.11.0839) и гранта Президента РФ МК-1652.2012.2 при организационной поддержке Лаборатории информационных технологий ННГУ.

Список литературы

1. Wilson B.C., Adam G.A Monte Carlo model for the absorption and flux distributions of light in tissue // Med. Phys. 1983. 10. P. 824-830.

2. Prahl S.A., Keijzer M., Jacques S.L., Welch A.J. A Monte Carlo model of light propagation in tissue // Proc. SPIE IS. 1989. 5. P. 102-111.

3. Wang L.V., Jacques S.L., Zheng L.Q. MCML-Monte Carlo modeling of light transport in multi-layered tissues // Comput. Meth. Prog. Biol. 1995. 47. P. 131-146.

4. Boas D.A., Culver J.P., Stott J.J., Dunn A.K. Three dimensional Monte Carlo code for photon migration through complex heterogeneous media including the adult human head // Opt. Express. 2002. 10. P. 159-170.

5. Tian P.F., Devor A., Sakadzic S., et al. Monte Car-

lo simulation of the spatial resolution and depth sensitivity of two-dimensional optical imaging of the brain // J. Biomed. Opt. 2011. 16. P. 13.

6. Ren N., Liang J., Qu X., et al. GPU-based Monte

Carlo simulation for light propagation in complex heterogeneous tissues // Opt. Express. 2010. 18.

P. 6811-6823.

7. Margallo-Balbas E., French P.J. Shape based Monte Carlo code for light transport in complex heterogeneous tissues // Opt. Express. 2007. 15. P. 14086-14098.

8. Гергель В.П., Стронгин Р.Г. Опыт Нижегородского университета по подготовке специалистов в области суперкомпьютерных технологий // Вестник Нижегородского университета им. Н.И. Лобачевского. 2010. № 3 (1). С. 191-199.

9. Баркалов К.А., Гергель В.П., Гергель А.В. и др. Организация и проведение Всероссийской школы по суперкомпьютерным технологиям // Открытое и дистанционное образование. 2010. №2. С. 24-29

10. Moller T., Trumbore B. Fast, minimum storage ray-triangle intersection // J. Graphics Tools. 1997. 2. P. 21-28.

11. Wald I. Realtime ray tracing and interactive glo-

bal illumination // PhD thesis, Saarland University, 2004.

12. Foley T., Sugerman J. KD-tree acceleration structures for a GPU raytracer // In Proceedings of the ACM SigGraph/Eurographics Conference on Graphics hardware. 2005. P. 15-22.

13. Horn D., Sugerman J., Houston M., Hanrahan P. Interactive kd-tree GPU raytracing // Proceedings of the symposium on Interactive 3D graphics and games on fast rendering. 2007. P. 167-174.

14. Ernst M., Greiner G. Early split clipping for bounding volume hierarchies // Proceedings of the IEEE Symposium on Interactive Ray Tracing. 2007. P. 73-78.

15. Wald I. On fast construction of SAH-based bounding volume hierarchies // Proceedings of the Eurographics Symposium on Interactive Ray Tracing. 2007. P. 33-40.

16. Боголепов Д.К., Сопин Д.П., Турлапов В.Е., Ульянов Д.Я. Построение SAH BVH деревьев для трассировки лучей на GPU в реальном времени // Высокопроизводительные параллельные вычисления на кластерных системах: Докл. XI Всерос. конф. Н. Новгород, 2-3 нояб. 2011. С. 301-306.

OPTIMAL INTERSECTION SEARCH ALGORITHM IN MONTE CARLO SIMULATION OF LIGHT TRANSPORT IN HUMAN BRAIN

A.V. Gorshkov, A.L. Korshunova

We consider the problem of choosing an optimal search algorithm for photon trajectory intersections with the boundaries in Monte Carlo light transport simulation in multilayer biological tissues with complex geometry based on the example of the human head. A review is given of main intersection search algorithms, their strengths and weaknesses, as well as the possibility of their adaptation to the problem. We present the results of computing experiment tests to compare the most appropriate intersection search algorithms, as well as their analysis in terms of selection of the optimal algorithm.

Keywords: light transport simulation, Monte-Carlo method, intersection search, kd-trees, BVH trees.

i Надоели баннеры? Вы всегда можете отключить рекламу.