Свертка для оптических откликов широкого луча в рассеивающих средах

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

Общие формулы свертки

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

Реакциями методов переноса фотонов могут быть физические величины, такие как поглощение , флюенс , отражение или пропускание . Учитывая конкретную физическую величину G(x,y,z) от карандашного луча в декартовом пространстве и коллимированного источника света с профилем луча S(x,y) , отклик широкого луча можно рассчитать, используя следующие 2- Формула свертки D:

C(x,y,z)=\int _{-\infty }^{\infty }\int _{-\infty }^{\infty }\ G(x-x',y-y',z)S(x',y')\,dx'\,dy'.\qquad (1)

Подобно 1D свертке, 2D свертка коммутативна между G и S с заменой переменных. $x''=x-x'\,$ и $y''=y-y'\,$ :

C(x,y,z)=\int _{-\infty }^{\infty }\!\int _{-\infty }^{\infty }\ G(x'',y'',z)S(x-x'',y-y'')\,dx''\,dy''.\qquad (2)

Поскольку ответ широкого луча $C(x,y,z)\,$ имеет цилиндрическую симметрию, его интегралы свертки можно переписать как:

C(r,z)=\int _{0}^{\infty }\ S(r')r'\left[\int _{0}^{2\pi }\ G\left({\sqrt {r^{2}+r'\,^{2}-2rr'cos\phi '}},z\right)\,d\phi '\right]dr'\qquad (3)

C(r,z)=\int _{0}^{\infty }G(r'',z)r''\left[\int _{0}^{2\pi }S\left({\sqrt {r^{2}+r''\,^{2}-2rr''cos\phi ''}}\right)\,d\phi ''\right]dr''\qquad (4)

где $r'={\sqrt {x'^{2}+y'^{2}}}$ . Поскольку внутреннее интегрирование уравнения 4 не зависит от z , его необходимо рассчитать только один раз для всех глубин. Таким образом, эта форма ответа с широким лучом более выгодна в вычислительном отношении.

Общие профили балок

Гауссов пучок

Для гауссова луча профиль интенсивности определяется выражением

S(r')=S_{0}\exp \left[-2\left({\frac {r'}{R}}\right)^{2}\right].\qquad (5)

Здесь R обозначает ${\tfrac {1}{e^{2}}}\,$ радиус луча, а S ₀ обозначает интенсивность в центре луча. S ₀ связана с полной мощностью P ₀ соотношением

S_{0}={\frac {2P_{0}}{\pi R^{2}}}.\qquad (6)

Подставив уравнение 5 в уравнение. 4, получаем

C(r,z)=2\pi S(r)\int _{0}^{\infty }G(r'',z)\exp \left[-2\left({\frac {r''}{R}}\right)^{2}\right]I_{0}\left({\frac {4rr''}{R^{2}}}\right)r''\,dr'',\qquad (7)

где I ₀ нулевого порядка — модифицированная функция Бесселя .

Цилиндрическая балка

Для балки цилиндра радиуса R функция источника принимает вид

S(r')={\begin{cases}S_{0},&{\text{if }}r'\leq R\\\,0,&{\text{if }}r'>R\end{cases}}\qquad (8)

где S ₀ обозначает интенсивность внутри пучка. S ₀ связана с полной мощностью пучка P ₀ соотношением

S_{0}={\frac {P_{0}}{\pi R^{2}}}.\qquad (9)

Подставив уравнение 8 в уравнение. 4, получаем

C(r,z)=2\pi S_{0}\int _{0}^{\infty }G(r'',z)I_{\phi }(r,r'')r''\,dr'',\qquad (10)

где

I_{\phi }(r,r'')={\begin{cases}1,&{\mbox{if }}R\geq r+r'',\\{\tfrac {1}{\pi }}\cos ^{-1}\left({\tfrac {r^{2}+r''^{2}-R^{2}}{2rr''}}\right),&{\mbox{if }}\left|r-r''\right|\leq R<r+r'',\\0,&{\mbox{if }}R<\left|r+r''\right|.\end{cases}}\qquad (11)

Ошибки в числовой оценке

Первые взаимодействия

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

C(r,z)=G_{1}(0,z){\frac {\delta (r)}{2\pi r}}+G_{2}(r,z),\qquad (12)

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

C(r,z)=G_{1}(0,z)S(r)+2\pi S_{0}\int _{0}^{\infty }G_{2}(r'',z)\,exp\left[-2\left({\frac {r''-r}{R}}\right)^{2}\right]I_{0e}\left({\frac {4rr''}{R^{2}}}\right)r''\,dr''.\qquad (13)

Для балки-цилиндра мы имеем

C(r,z)=G_{1}(0,z)S(r)+2\pi S_{0}\int _{0}^{\infty }G_{2}(r'',z)I_{\phi }(r,r'')r''\,dr''.\qquad (14)

Ошибка усечения

Для балки-цилиндра верхние пределы интегрирования могут быть ограничены r _max , так что r ≤ r _max − R . Таким образом, ограниченное покрытие сетки в направлении r не влияет на свертку. Для надежной свертки физических величин при r в ответ на луч в цилиндре мы должны гарантировать, что r _max в методах транспортировки фотонов достаточно велико, чтобы r ≤ r _max − R. выполнялось Для гауссова луча не существует простых верхних пределов интегрирования, поскольку теоретически он простирается до бесконечности. При r >> R гауссовский луч и цилиндр с одинаковыми R и S ₀ имеют сравнимые результаты свертки. Следовательно, r ≤ r _max − R можно приближенно использовать и для гауссовских пучков.

Реализация свертки

Для реализации дискретной свертки используются два распространенных метода: определение свертки и быстрое преобразование Фурье (БПФ и ОБПФ) согласно теореме о свертке . Для расчета оптического отклика широкого луча импульсная характеристика карандашного луча свертывается с функцией луча. Как показано уравнением 4, это двумерная свертка. Чтобы вычислить реакцию светового луча на плоскости, перпендикулярной оси z, функция луча (представленная матрицей b × b ) свертывается с импульсной характеристикой в этой плоскости (представленной матрицей a × a ). Обычно a больше b . Эффективность расчета этих двух методов во многом зависит от b , размера светового луча.

При прямой свертке матрица решения имеет размер ( a + b − 1) × ( a + b − 1). Вычисление каждого из этих элементов (кроме тех, что находятся вблизи границ) включает умножения b × b и сложения b × b − 1, поэтому временная сложность равна O [( a + b ) ²б ²]. При использовании метода БПФ основными этапами являются БПФ и ОБПФ матриц ( a + b − 1) × ( a + b − 1), поэтому временная сложность равна O[( a + b ) ² журнал( а + б )]. Сравнивая O[( a + b ) ²б ²] и O[( a + b ) ² log( a + b )], очевидно, что прямая свертка будет быстрее, если b намного меньше a , но метод БПФ будет быстрее, если b относительно велико.

См. также

Ссылки на другие ресурсы Монте-Карло

Ссылки

Л.-Х. Ван и Х.-И. Ву. Биомедицинская оптика: принципы и визуализация. Уайли 2007.
Л.-Х. Ван, С.Л. Жак и Л.-К. Чжэн, «Моделирование транспорта фотонов в многослойных тканях методом Монте-Карло», Компьютерные методы и программы в биомедицине 47, 131–146 (1995).
Л.-Х. Ван, С.Л. Жак и Л.-К. Чжэн, «Свертка для ответов на фотонный пучок конечного диаметра, падающий на многослойные ткани», Компьютерные методы и программы в биомедицине 54, 141–150 (1997). Скачать статью .