Моделирование распространения фотонов методами Монте-Карло - это гибкий, но строгий подход к моделированию переноса фотонов. В этом методе локальные правила переноса фотонов выражаются в виде распределений вероятностей, которые описывают размер шага движения фотона между участками взаимодействия фотона с тканью и углы отклонения траектории фотона, когда происходит событие рассеяния. Это эквивалентно аналитическому моделированию переноса фотонов с помощью уравнения переноса излучения (RTE), которое описывает движение фотонов с помощью дифференциального уравнения. Однако решения RTE в закрытой форме часто невозможны; для некоторых геометрий диффузионное приближениеможет использоваться для упрощения RTE, хотя это, в свою очередь, вносит много неточностей, особенно вблизи источников и границ. Напротив, моделирование методом Монте-Карло можно сделать произвольно точным, увеличив количество отслеживаемых фотонов. Например, посмотрите фильм, в котором моделирование методом Монте-Карло « карандашного луча», падающего на полубесконечную среду, моделирует как начальный поток баллистических фотонов, так и более позднее диффузное распространение.
Метод Монте-Карло обязательно является статистическим и поэтому требует значительного времени вычислений для достижения точности. Кроме того, моделирование методом Монте-Карло позволяет отслеживать несколько физических величин одновременно с любым желаемым пространственным и временным разрешением. Эта гибкость делает моделирование методом Монте-Карло мощным инструментом. Таким образом, будучи неэффективными с вычислительной точки зрения, методы Монте-Карло часто считаются стандартом для моделирования измерений переноса фотонов для многих биомедицинских приложений.
Биомедицинские применения методов Монте-Карло
Биомедицинская визуализация
Оптические свойства биологической ткани предлагают подход к биомедицинской визуализации. Есть много эндогенных контрастов, включая поглощение из крови и меланина и рассеяние от нервных клеток и ядер раковых клеток. Кроме того, флуоресцентные зонды могут быть нацелены на множество различных тканей. Методы микроскопии (включая конфокальную , двухфотонную и оптическую когерентную томографию ) имеют возможность отображать эти свойства с высоким пространственным разрешением, но, поскольку они полагаются на баллистические фотоны, их глубина проникновения ограничена несколькими миллиметрами. Для получения изображений в более глубокие слои тканей, где фотоны были многократно рассеяны, требуется более глубокое понимание статистического поведения большого количества фотонов в такой среде. Методы Монте-Карло обеспечивают гибкую основу, которая используется различными методами для восстановления оптических свойств глубоко внутри ткани. Здесь представлено краткое введение в некоторые из этих методов.
- Фотоакустическая томография В PAT рассеянный лазерный свет поглощается, что приводит к локальному повышению температуры. Это локальное изменение температуры, в свою очередь, генерирует ультразвуковые волны за счет термоупругого расширения, которые регистрируются ультразвуковым преобразователем. На практике варьируются различные параметры настройки (например, длина волны света, числовая апертура преобразователя), и в результате моделирование методом Монте-Карло является ценным инструментом для прогнозирования реакции ткани до применения экспериментальных методов.
- Диффузная оптическая томография DOT - это метод визуализации, который использует массив источников света и детекторов ближнего инфракрасного диапазона для измерения оптических свойств биологических тканей. Можно измерить различные контрасты, включая поглощение оксигемоглобина и дезоксигемоглобина (для функциональной нейровизуализации или обнаружения рака) и концентрацию флуоресцентных зондов. Чтобы восстановить изображение, необходимо знать, каким образом свет проходит от данного источника к данному детектору и как измерение зависит от распределения и изменений оптических свойств (известная как прямая модель). Из-за сильно рассеивающей природы биологической ткани такие пути сложны, а функции чувствительности размыты. Форвардная модель часто создается с использованием методов Монте-Карло.
Лучевая терапия
Цель лучевой терапии - доставить энергию, обычно в форме ионизирующего излучения, к раковой ткани, сохраняя при этом окружающие нормальные ткани. Моделирование методом Монте-Карло обычно используется в лучевой терапии для определения периферической дозы, которую пациент будет испытывать из-за рассеяния как от ткани пациента, так и рассеяния от коллимации выше по потоку в линейном ускорителе.
Фотодинамическая терапия
В фотодинамической терапии (PDT) свет используется для активации химиотерапевтических агентов. Из-за природы ФДТ полезно использовать методы Монте-Карло для моделирования рассеяния и поглощения в ткани, чтобы гарантировать, что соответствующие уровни света доставляются для активации химиотерапевтических агентов.
Реализация фотонного транспорта в рассеивающей среде
Представлена модель метода фотонного Монте-Карло в однородной бесконечной среде. Однако модель легко расширяется для многослойных носителей. Для неоднородной среды необходимо учитывать границы. Кроме того, для полубесконечной среды (в которой фотоны считаются потерянными, если они выходят за верхнюю границу), необходимо уделить особое внимание. Для получения дополнительной информации перейдите по ссылкам внизу страницы. Мы решим проблему, используя бесконечно малый точечный источник (аналитически представленный как дельта-функция Дирака в пространстве и времени). Отклики на произвольную геометрию источника могут быть построены с использованием метода функций Грина (или свертки , если существует достаточная пространственная симметрия). Обязательными параметрами являются коэффициент поглощения, коэффициент рассеяния и фазовая функция рассеяния. (Если рассматриваются границы, также должен быть указан показатель преломления для каждой среды.) Отклики с временным разрешением находятся путем отслеживания общего времени, прошедшего с момента пролета фотона, с использованием длины оптического пути . Ответы на источники с произвольными временными профилями можно затем смоделировать посредством свертки по времени.
В нашей упрощенной модели мы используем следующую технику уменьшения дисперсии, чтобы сократить время вычислений. Вместо того, чтобы распространять фотоны по отдельности, мы создаем пакет фотонов с определенным весом (обычно инициализируемым как единица). Когда фотон взаимодействует в мутной среде, он будет набирать вес из-за поглощения, а оставшийся вес будет рассеиваться на другие части среды. По пути может регистрироваться любое количество переменных, в зависимости от интересов конкретного приложения. Каждый фотонный пакет будет многократно проходить следующие пронумерованные шаги, пока он не будет завершен, отражен или передан. Процесс показан на схеме справа. Любое количество пакетов фотонов может быть запущено и смоделировано до тех пор, пока в результате смоделированных измерений не будет достигнуто желаемое отношение сигнал / шум. Обратите внимание, что поскольку моделирование методом Монте-Карло - это статистический процесс, включающий случайные числа, мы будем использовать переменную ξ повсюду как псевдослучайное число для многих вычислений.
Шаг 1: запуск фотонного пакета
В нашей модели мы игнорируем начальную зеркальную отражательную способность, связанную с попаданием в среду, показатель преломления которой не согласован. Имея это в виду, нам просто нужно установить начальное положение фотонного пакета, а также начальное направление. Удобно использовать глобальную систему координат. Мы будем использовать три декартовых координаты для определения положения, а также три направляющих косинуса для определения направления распространения. Начальные условия запуска будут варьироваться в зависимости от приложения, однако для стержневого луча, инициализированного в начале координат, мы можем установить начальное положение и направляющие косинусы следующим образом (изотропные источники можно легко смоделировать путем рандомизации начального направления каждого пакета):
Шаг 2: выбор размера шага и движение фотонного пакета
Размер шага s - это расстояние, которое пакет фотонов проходит между местами взаимодействия. Существует множество методов выбора размера шага. Ниже приведена основная форма выбора размера шага фотона (полученная с использованием метода обратного распределения и закона Бера – Ламберта ), которую мы используем для нашей однородной модели:
где случайное число и - полный коэффициент взаимодействия (т. е. сумма коэффициентов поглощения и рассеяния).
После выбора размера шага пакет фотонов распространяется на расстояние s в направлении, определяемом направляющими косинусами. Это легко сделать, просто обновив координаты следующим образом:
Шаг 3: Поглощение и рассеяние
Часть веса фотона поглощается в каждом месте взаимодействия. Эта доля веса определяется следующим образом:
где - коэффициент поглощения.
Затем массовая доля может быть записана в виде массива, если распределение поглощения представляет интерес для конкретного исследования. Затем вес фотонного пакета должен быть обновлен следующим образом:
После поглощения фотонный пакет рассеивается. Средневзвешенное значение косинуса угла рассеяния фотона известно как анизотропия рассеяния ( g ), которая имеет значение от -1 до 1. Если оптическая анизотропия равна 0, это обычно указывает на то, что рассеяние является изотропным. Если g приближается к значению 1, это означает, что рассеяние в основном идет в прямом направлении. Чтобы определить новое направление пакета фотонов (и, следовательно, косинусы направления фотонов), нам необходимо знать фазовую функцию рассеяния. Часто используется фазовая функция Хеньи-Гринштейна. Затем угол рассеяния θ определяется по следующей формуле.
И, как правило, предполагается , что полярный угол φ равномерно распределен между 0 и. Исходя из этого предположения, мы можем установить:
На основе этих углов и исходных направляющих косинусов мы можем найти новый набор направляющих косинусов. Новое направление распространения можно представить в глобальной системе координат следующим образом:
Для особого случая
использовать
или же
использовать
C-код:
/ *********************** Индикаторная матрица ********************** Новые направляющие косинусы после рассеяния на угол theta, fi.* mux new = (sin (theta) * (mux * muz * cos (fi) -muy * sin (fi))) / sqrt (1-muz ^ 2) + mux * cos (theta)* muy new = (sin (theta) * (muy * muz * cos (fi) + mux * sin (fi))) / sqrt (1-muz ^ 2) + muy * cos (theta)* muz новый = - sqrt (1-muz ^ 2) * sin (theta) * cos (fi) + muz * cos (theta)* ------------------------------------------------- --------*Вход:* muxs, muys, muzs - косинус направления перед столкновением* mutheta, fi - косинус полярного угла и азимутального угла* ------------------------------------------------- --------*Выход:* muxd, muyd, muzd - косинус направления после столкновения* ------------------------------------------------- --------* /void Indicatrix (двойные мультиплексоры, двойные муйсы, двойные музы, двойные мутхеты, двойные фи, двойные * мультиплексоры, двойные * муйд, двойные * муз){ двойная костета = мутета; двойная синтета = sqrt (1.0-costheta * costheta); // грех (тета) двойной sinfi = sin (fi); двойной cosfi = cos (fi); if (muzs == 1.0) { * muxd = sintheta * cosfi; * muyd = sintheta * sinfi; * музд = костета; } elseif (muzs == -1.0) { * muxd = sintheta * cosfi; * muyd = -sintheta * sinfi; * muzd = -costheta; } еще { двойной деном = sqrt (1.0-muzs * muzs); двойной muzcosfi = muzs * cosfi; * muxd = sintheta * (muxs * muzcosfi-muys * sinfi) / denom + muxs * costheta; * muyd = sintheta * (muys * muzcosfi + muxs * sinfi) / деном + muys * costheta; * муз = -деном * синтета * косфи + муз * костхета; }}
Шаг 4: прекращение действия фотона
Если пакет фотонов испытал много взаимодействий, для большинства приложений вес, оставленный в пакете, не имеет большого значения. В результате необходимо определить средства для прерывания пакетов фотонов достаточно малого веса. В простом методе использовался бы порог, и если вес фотонного пакета ниже порога, пакет считается мертвым. Вышеупомянутый метод ограничен, поскольку он не экономит энергию. Чтобы сохранить постоянную общую энергию, часто используют технику русской рулетки для фотонов ниже определенного порога веса. В этом методе используется постоянная рулетки m, чтобы определить, выживет ли фотон. У фотонного пакета есть один шанс в m на выживание, и в этом случае ему будет присвоен новый вес в мВт, где W - начальный вес (этот новый вес в среднем сохраняет энергию). Во всех остальных случаях вес фотона устанавливается на 0, и фотон прекращается. Математически это выражается ниже:
Графические процессоры (GPU) и быстрое моделирование переноса фотонов методом Монте-Карло
Моделирование миграции фотонов в мутной среде методом Монте-Карло представляет собой задачу с высокой степенью распараллеливания, когда большое количество фотонов распространяется независимо, но в соответствии с одинаковыми правилами и разными последовательностями случайных чисел. Параллельный характер этого особого типа моделирования Монте-Карло делает его очень подходящим для выполнения на графическом процессоре (GPU). Выпуск программируемых графических процессоров положил начало такому развитию, и с 2008 года появилось несколько отчетов об использовании графических процессоров для высокоскоростного моделирования миграции фотонов методом Монте-Карло. [1] [2] [3] [4]
Этот базовый подход сам по себе можно распараллелить, используя несколько связанных между собой графических процессоров. Одним из примеров является «MCML GPU Cluster», который можно загрузить с веб-сайта авторов (Монте-Карло Моделирование переноса света в многослойных мутных средах на основе кластеров GPU): http://bmp.hust.edu.cn/ GPU_Cluster / GPU_Cluster_MCML.HTM
Смотрите также
Ссылки на другие ресурсы Монте-Карло
- Лаборатория оптической визуализации Вашингтонского университета в Сент-Луисе (MCML)
- Медицинский лазерный центр Орегона
- Миграция фотонов Исследование методом Монте-Карло в Лундском университете, Швеция Ускорение моделирования методом Монте-Карло на графическом процессоре и масштабируемый Монте-Карло. Открытый исходный код для скачивания.
- Монте-Карло на основе облаков для переноса света в мутной рассеивающей среде Инструмент можно бесплатно использовать в исследовательской и некоммерческой деятельности.
- Легкий транспорт в тканях как пример моделирования Монте-Карло (с исходным кодом C ++).
Рекомендации
- Wang, LH; У Синь-И (2007). Биомедицинская оптика: принципы и изображения . Вайли.
- Л.-Х. Ванга; С.Л. Жак; L.-Q. Чжэн (1995). «MCML - Монте-Карло моделирование переноса света в многослойных тканях». Компьютерные методы и программы в биомедицине . 47 (2): 131–146. DOI : 10.1016 / 0169-2607 (95) 01640-F .
- Л.-Х. Ванга; С.Л. Жак; L.-Q. Чжэн (1997). «Conv - свертка для ответов на пучок фотонов конечного диаметра, падающий на многослойные ткани» (PDF) . Компьютерные методы и программы в биомедицине . 54 (3): 141–150. DOI : 10.1016 / S0169-2607 (97) 00021-7 .
- С.Л. Жак; Л.-Х. Ван (1995). «Моделирование переноса света в тканях методом Монте-Карло» (PDF) . В AJ Welch; MJC van Gemert (ред.). Оптический термический отклик ткани, облученной лазером . Нью-Йорк: Пленум Пресс. С. 73–100.
- Л.-Х. Ванга; С.Л. Жак (1994). «Оптимизированные радиальные и угловые положения в моделировании Монте-Карло» (PDF) . Медицинская физика . 21 (7): 1081–1083. Bibcode : 1994MedPh..21.1081W . DOI : 10.1118 / 1.597351 . PMID 7968840 .
Встроенные ссылки
- ^ Э. Алерстам; Т. Свенссон; С. Андерссон-Энгельс (2008). «Параллельные вычисления с графическими процессорами для высокоскоростного моделирования миграции фотонов методом Монте-Карло» (PDF) . J. Biomed. Опт . 13 (6): 060504. Bibcode : 2008JBO .... 13f0504A . DOI : 10.1117 / 1.3041496 . PMID 19123645 .
- ^ Q. Fang; Д.А. Боас (2009). «Моделирование миграции фотонов в трехмерных мутных средах методом Монте-Карло, ускоренное графическими процессорами» . Опт. Экспресс . 17 (22): 20178–20190. Bibcode : 2009OExpr..1720178F . DOI : 10.1364 / oe.17.020178 . PMC 2863034 . PMID 19997242 .
- ^ Н. Рен; Дж. Лян; X. Qu; Дж. Ли; Б. Лу; Дж. Тиан (2010). «Моделирование методом Монте-Карло на основе графического процессора для распространения света в сложных гетерогенных тканях» . Опт. Экспресс . 18 (7): 6811–6823. Bibcode : 2010OExpr..18.6811R . DOI : 10.1364 / oe.18.006811 . PMID 20389700 .
- ^ А. Доронин; И. Меглинский (2011). «Онлайновый объектно-ориентированный вычислительный инструмент Монте-Карло для нужд биомедицинской оптики» . Биомед. Опт. Экспресс . 2 (9): 2461–2469. DOI : 10,1364 / boe.2.002461 . PMC 3184856 . PMID 21991540 .