Дозиметрическое планирование дистанционной лучевой терапии на основе метода Монте-Карло тема автореферата и диссертации по физике, 01.04.16 ВАК РФ

На правах рукописи

Чупикин Дмитрий Анатольевич

ДОЗИМЕТРИЧЕСКОЕ ПЛАНИРОВАНИЕ ДИСТАНЦИОННОЙ ЛУЧЕВОЙ ТЕРАПИИ НА ОСНОВЕ МЕТОДА МОНТЕ-КАРЛО

Специальность 01 04 16 - Физика атомного ядра и элементарных частиц

АВТОРЕФЕРАТ диссертации на соискание ученой степени кандидата физико-математических наук

1111111111«»

ООЗ 15еэ г ».я

Москва 2007

Работа выполнена в Государственном образовательном учреждении высшего профессионального образования «МИФИ» - Московском государственном инженерно-физическом институте (техническом университете)

Научный руководитель

Московский государственный инженерно-физический институт, г Москва Доктор физико-математических наук,

Профессор Климанов Владимир Александрович

Научный консультант

ВНИИ экспериментальной физики РАН, г Саров Кандидат физико-математических наук,

доцент Донской Евгений Николаевич

Официальные оппоненты

Обнинский Государственный Технический Университет Атомной Энергетики доктор физико-математических наук,

профессор Андросенко Петр Александрович

НИИ экспериментальной ядерной физики МГУ, г Москва кандидат физико-математических наук,

доцент Кэбин Эдуард Йоханнесович

Ведущая организация

ГНЦ РФ Институт теоретической и экспериментальной физики РАН, г Москва

Защита состоится « 01 » ноября 2007 г в 14_ часов на заседании диссертационного совета Д 501 001 65 при Биологическом факультете МГУ им M В Ломоносова - Московском Государственном Университете им M В Ломоносова по адресу 119992, Россия, Москва, Ленинские горы, д 1, стр 12, Биологический факультет МГУ, ауд

С диссертацией можно ознакомиться в библиотеке Биологического факультета МГУ им М В Ломоносова

Отзыв на автореферат в 2-х экз просим направлять по указанному адресу

Автореферат разослан « » ^^_2007 г

Ученый секретарь диссертационного _ ^

совета, кандидат биологических наук Веселова Т В

ОБЩАЯ ХАРАКТЕРИСТИКА РАБОТЫ

Актуальность работы.

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

В настоящее время в таких передовых системах применяются в основном два метода 1) метод свертки/суперпозиции, использующий дозовые ядра дифференциального тонкого луча, 2) метод тонкого луча в различных модификациях Эти методы пришли на смену эмпирическим и полуэмпирическим методам Они обладают достаточной быстротой и обеспечивают в случае однородности среды высокую точность расчета Однако при наличии гетерогенности в области расчета, а также вблизи поверхности тела пациента погрешность расчета с помощью этих методов сильно возрастает Принято считать, что единственным методом, который в условиях сложной геометрии расчета может обеспечить рекомендуемую международными комиссиями точность, является вероятностный метод Монте-Карло. Важно, что метод не требует записи соответствующего уравнения переноса, фактически здесь нет проблем с формулировкой дискретной модели В обоих случаях реализация вычислительного алгоритма не требует упрощений Метод максимально при-

способлен для использования на ЭВМ Единственным ограничением являются затраты машинного времени для получения результата с необходимой точностью В настоящее время эти проблемы не являются решенными. Но в последнее время и в этом направлении ведутся целенаправленные работы и имеются определенные достижения

Таким образом, метод Монте-Карло является наиболее точным методом расчета доз в лучевой терапии, и поэтому, ускорение расчетов путем модификации метода Монте-Карло является актуальной научно-технической задачей

Пелыо работы является усовершенствование одной из модификаций метода Монте-Карло, а именно, метода РЬ-оценок потока в точке, а также разработка новых моделей источников для использования их в методе Монте-Карло

При этом решаются следующие основные задачи

1 Поиск возможных способов ускорения расчетов поглощенной дозы с использованием метода РЬ-оценок На основе исследований предлагаются методы ускорения расчетов

2 Метод РЬ-оценок потока в точке обобщается на случай расчета доз в неоднородной среде

3 Разработка корректных математических моделей источников излучения радиотерапевтических ускорителей и облучателей с б0Со с целью использования их в методе РЬ-оценок

4 Разработка моделей фотонных источников излучения радиотерапевтических установок для использования их в расчетах различными алгоритмами метода Монте-Карло в лучевой терапии

Теоретическая база исследования Для решения поставленных в работе задач использовалась теория переноса излучения, теория вероятностей и математическая статистика

Научная новизна настоящей работы заключается в следующем

1. Впервые предложены методы ускорения метода РЬ-оценок потока в точке

2 Впервые показана возможность использования метода РЬ-оденок потока в точке при расчетах доз в химически неоднородных средах

3.У совершенствованы модели источников используемых в методе РЬ-оценок потока в точке

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

5 Предложен метод восстановления энергетического спектра дистанционного радиотерапевтического облучателя на основе пространственного распределения энергопоглощений

Практическая значимость результатов работы состоит прежде всего в том, что предложенные в работе быстрые модификации метода РЬ-оценок с корректным учетом неоднородностей, а также новые математические модели источников дистанционной лучевой терапии, позволили создать на основе метода Монте-Карло методику расчета доз в лучевой терапии за приемлемое для клиник расчетное время

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

Основные положения, выносимые на защиту:

1 Способы ускорения расчетов при использовании метода РЬ-оценок

2 Методика учета гетерогенностей при расчете доз в лучевой терапии при использовании метода РЬ-оценок

3 Новые модели фотонных источников излучения для использования их в методе РЬ-оценок

4 Методика реконструкции эффективного спектра тормозного излуче-

ния на основании обработки пространственных дозных распределений, измеренных при вводе в эксплуатацию радиотерапевтических ЛУЭ и кобальтовых облучателей

Апробация работы Основные материалы и положения работы докладывались и обсуждались на следующих научных конференциях и семинарах Научная конференция в Московском инженерно-физическом институте «Научная сессия МИФИ 2004, 2005» - Москва, 2004г, 2005 г

II Евразийский конгресс по медицинской физике и инженерии «Медицинская физика 2005», Москва, 2005 г Публикации.

Материалы диссертации опубликованы в 5 печатных работах 1 статье в научном журнале из списка ВАК , 4 научных трудах и тезисах конференций 1 статья принята к изданию в журнале «Атомная энергия»

Структура и объем диссертации Диссертация состоит из введения, трех глав, заключения, списка литературы и приложения Общий объем составляет 131 страницу печатного текста, включая 50 рисунков Список литературы включает 42 наименования, из них 30 иностранных

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

производителей систем дозиметрического планирования лучевой терапии, у которых запланирована (или уже введена) реализация метода Монте-Карло.

Показаны существующие методы ускорения расчетов методом Монте-Карло В главе рассмотрены наиболее часто встречающиеся способы ускорения расчетов, показаны их условия применения и уменьшение времени расчета вследствие их применения Использование данных способов напрямую влияет на ускорение расчетов методом Монте-Карло, однако их применение также косвенно влияет на ускорение расчетов с помощью метода РЬ-оценок

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

Вторая глава настоящей работы посвящена исследованию метода РЬ-оценок Суть РЬ-оценок состоит в том, что внутри традиционной схемы метода Монте-Карло используются новые оценки для вычисления вклада в рассматриваемые функционалы При этом преодолеваются два существенных недостатка, присущие методу Монте-Карло медленная сходимость метода и значительная флуктуация результатов расчетов При использовании РЬ-оценок скорость счета повышается и исчезают флуктуации в результатах расчетов, так как методу РЬ-оценок присуще такое замечательное свойство как внутренне сглаживание расчетных данных Расчеты с использованием РЬ-оценок также можно комбинировать с другими быстрыми приближенными методами расчета доза, основанными либо на использовании экспериментальных данных, либо на методе суперпозиции/свертки и применять их только для вычисления только в тех точках, где приближенные методы дают недопустимую погрешность

Расчет поглощенной дозы с использование РЬ-оценок делится на два этапа

• построение дерева траекторий в бесконечной однородной водной среде

• РЬ обработка дерева траекторий При построении дерева траектории используются следующие правила частица источника - фотон, фазовые координаты для всех частиц источника (?„,.£„) Непосредственно в точке источника происходит взаимодействие частицы с веществом, те построение траектории происходит для источника первых столкновений

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

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

1 Импульсные фазовые координаты основной и подобной траектории совпадают Е"т1°г = Е^пм"у = Е,

2 Оптические длины каждого звена основной и подобной траектории совпадают

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

В случае воксельной геометрии уравнение поиска расстояния до узловой точки ДГ( выглядит следующим образом

N

Ы1

N - число пересеченных вокселей, рк- плотность водной среды в к-м вокселе, 1к -расстояние пройденное в k-м вокселе, вдоль направления (-и,)

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

^ h{?rlar\^t / f8^, (2)

Н source

С^-переводной множитель из МэВ/г в пГр, w, -вес звена на траектории (Ё) 1

Функционал = -Е£а{Е) если текущее звено - трек гамма-кванта, для Р

которого перенос заряженных частиц не учитывается, 0 если текущее звено -трек гамма-кванта, для которого перенос заряженных частиц учитывается, у/(Е) 1

г—=---(Д) если текущее звено - трек заряженной частицы М, -длина

v, р ах

текущего звена траектории, Esmrce = ^Ер {E)dE -средняя энергия источника в МэВ. Ssourcc - площадь виртуального источника Функция источника h{p("""'ar, Еа) представляется в виде

к(гГ'аг ,£„)=£(£„)ехр

f Z&)

\

р

Рн-/> J

(3)

х(£0)- полное сечение взаимодействия первичного гамма-кванта со средой, ¿„расстояние в -^-от точки первого взаимодействия г/"""" до границы

см

воксельной геометрии в направлении обратном движению первичного гам-

_ пГр см1 ма-кванта Размерность дозы есть —--

г МэВ

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

hetero -Kjhetero

nequal _ г_^"e

^ ът water 5

• pke""° — плотность вещества в неоднородности

• iV^'"'" - число электронов в единице объема вещества в неоднородности

• N"°"r - число электронов в единице объема воды

NovaXYZ -MCPL

Рисунок 1 Глубинное распределение поглощенной дозы от источника б МУ соответственно в воде Расположение химической неоднородности {2 см < г< 5 см}, состав неоднородности - алюминий

Во-вторых, для учета поглощения излучения веществами с Т отличными от воды, необходимо в формулу (2) ввести следующий коэффициент

к = ^

мГШ {)

• Ма ~ массовый коэффициент поглощения вещества в детекторной точке

• мГ"г ~ массовый коэффициент поглощения в воде

• Е0- энергия источника Для энергетически распределенных источников используется средняя энергия спектра

По этой методике рассчитаны глубинные дозовые распределения, результаты сравнивались с результатами, посчитанными по программе N0-\zXYZ (рис 1)

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

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

Ниже приведена зависимость относительного времени расчета дозы от выбора граничного уровня взаимодействия

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

УРОВЕНЬ ВЗАИМОДЕЙСТВИЯ

Рисунок 2 Зависимость относительного времени расчета при разных граничных уровнях взаимодействия для источника излучения с энергией 5 МэВ

Пусть имеется три подряд идущих электронных вершины с координатами {Е,,?} 1=1, 2, 3, которые являются вершинами двух подряд идущих электронных треков Расстояния между вершинами определяется по формуле Ц Для того чтобы исключить из выборки среднюю вершину, те

рассматривать два последовательных трека как один, необходимо выполнение двух условий

а) в точке 2 не происходит образования тормозных фотонов и 5-электронов,

б) длины звеньев удовлетворяют следующему отношению

/12+/23 <к /13, (6) где к - коэффициент сглаживания

Поглощенная энергия на сглаженном треке /,3 вычисляется исходя из формулы-

(¿Е\

/„+ -

V ¿х

к.

(7)

Были выполнены расчеты дозовых распределений для различных источников излучения с применением описанной методики для к=[1 1,10] с шагом 0,1 Зависимость времени расчета от к показана на рис 3

1,0

,0,9

(г1

8

а.

I

Л

§ 0,7

н §

§ н

О., 0,6

0,5 -г

\

\ \ V I

9 •-<

1,02 1,04 1,06

Критерий сглаживания

1,08

1,1

Рисунок 3 Зависимость времени расчета от критерия сглаживания для энергии Е= 6МУ

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

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

der г

(8)

в'

dcran

dCl

.. dCl

- комптоновское дифференциальное сечение рассеяния на угол

-комптоновское сечение рассеяния на угол в, т е. угол рассея-

ния

47

45

Ч

ф

х я'

I

43

41

I

О 39

о с

37

35

-10

—— Бе' — — см ! смещения гщение по углу

.............

I i !

-5 0 5

Расстояние от оси источника

10

Рисунок 4 Дозовые профили, рассчитанные с использованием технологии смещения угла и без нее Источник 1,25 МэВ, &НЭ 80 см, размер поля 10x10 см

С применением описанной технологии учета расходимости первичного излучения были рассчитаны профили дозовых распределений от различных источников излучения Для полей 10x10 см и ББВ 100 см различие между профилями, рассчитанными с использованием сдвига и без него, незначительны Разница растет с уменьшением вББ, т к поток энергии все более и более расходится В качестве примера на рис 4 приведены дозовые профили от источника 1,25 МэВ, 88Б 80 см, размера поля 10x10 Как и следовало ожидать, применение описанной технологии «поднимает» дозу на краях плато дозового профиля

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

Рассмотрим сначала метод восстановления спектров тормозного излучения, который осуществляется на основе комбинированного использования нескольких глубинных дозных профилей в плоскостях перпендикулярных оси падения пучка Если измерения Рк проведены в плоскостях на глубинах к = 1,2 и известны соответствующие функции отклика Рк (Е ), где Е - значения энергий искомого спектра, то для каждого к - го уравнения Р* = 1<р(Е) Рк(Е)йЕ или в дискретном представлении

Р* =^<рк{Е„) Р^Е„)АЕ„ (9)

Дозовые профили для моноэнергетических фотонов - ядра Рк{Еп) уравнения (9) - рассчитывались методом Монте-Карло помощью кода БРМ в диапазоне энергий 0 1-20 МэВ для геометрических параметров, указанных ниже. Статистика расчета ядер составляла 50-Ю6 историй, статистическая погрешность результатов не превосходила 1-2% Исходные профили Рк так же рассчитывались по той же программе

Расчеты проводились для трех плоскостей к = 1,2,3 (водный фантом 30x30x30 см3, глубины 5см, 15см, 25см, БББ = 100 см и площадь облучения Б = 10x10 см2), полученные для первичного спектра начального приближения и для каждой итерации у три спектра усреднялись срку ={д>\ +<р*+<р*)/з и итерационный спектр ¡р) являлся исходным для у+1 итерации Естественным обоснованием этой процедуры является независимость формы восстановленного спектра от используемого глубинного дозного профиля с индексом к После проведения расчета комплекса итераций для всех трех глубин дозных профилей определялось относительное отклонение Щк) 10 центральных точек исходного профиля Р(к,х) от соответствующих значений 8(к,1), вычисляемых в ходе итерационного процесса

7 Р(к,1)

Рисунок 5 Результаты восста/оеления спектра тормозного излучения с граничной энергией 6 МУ (нормированы по пощади)

Выбиралась величина максимального отклонения Я(к)тах и с учетом знака отклонения изменялся параметр начального приближения в сторону уменьшения или увеличения (на 10%), формировалось новое начальное приближение и вся расчетная процедура повторялась до достижения минимума значения К(к)тах Результаты восстановления спектра показаны на рис 5

Суть другого метода моделирования источника состоит в следующем Мы имеем набор координат частиц (х,у,Е,и,у,м), которые должна получить программа метода Монте-Карло от программы моделирования вторичного источника, однако для определения направляющих косинусов доста-

точно двух углов в и ср (полярный и азимутальный углы, соответственно) Итак, необходимо сформировать функции распределения вероятности для пяти переменных (х,у,Е,в,<р) В настоящем случае обработки файлов фазового пространства в функции распределения вероятности порядок в иерархии переменных подразумевает, что для данного интервала переменной х суммируются все попавшие в интервал фотоны (независимо от у, Е, в, ср), для пары интервалов переменных (х,у) суммируются все попавшие в интервал фотоны (независимо от Е, в, /р), для тройки интервалов (х,у,Е) суммируются все попавшие в интервал фотоны (независимо от в, <р), и так далее

Для каждого считанного фотона определяются номера ячеек, номер энергетической группы Выполняется преобразование направляющих косинусов (и,у/w) —* (и ', V)') для уменьшения анизотропии по углу

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

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

1Е*0

1Е-1-

1Е-2-

1Р-1-1Е-4-

0.0 01 02 03 8 4 05 06 0.7 08 09 10 11 12

Е.МеУ

Рисунок 6 "Дифференциальные" энергетические распределения фотонов установки Пост, идущих через ячейки (iX.il9 поверхности запоминания фазового пространства черный - (гХ.гУ) = (1,1), красный-(¡ХД) = (4,4), синий-(гХ,гУ) = (5,5)

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

ЗАКЛЮЧЕНИЕ

В настоящей работе исследован метод РЬ-потока в точке и предложены его улучшения, направленные на увеличение скорости и точности расчетов Отметим, что скорость счета с использованием РЬ-оценок имеет почти линейную зависимость от количества детекторных точек. Это позволяет практически сразу рассчитать дозу в любой критической точке за время от 0,04 до 1,5 сек, в зависимости от условий расчета Использование обычных схем в подобных ситуациях практически невозможно, т к зависимость скорости от числа детекторных точек незначительна Эта особенность делает метод уникальным для проведения расчетов связанных с небольшим количеством точек. Обычно такие расчеты проводятся с исследовательскими или научными целями, где нужно получить глубинные дозовые распределения или профили дозовых распределений

Были предложены и исследованы методы ускорения расчетов поглощенной дозы методом РЬ-оценок Эти методы могут быть использованы как отдельно друг от друга, так и совместно Выигрыш по времени при применении этих методик может составлять до 2,5 раз Отметим, что расчет с использованием вышеописанных методик не вносит существенных изменений в сам метод РЬ-оценок А это значит, что использование этих двух способов ускорения не изменяет ограничения на область применения метода РЬ-оценок Метод возможно комбинировать с любыми другими методами расчета дозы

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

Полученная в результате данного исследования информация по восстановлению эффективного действующего спектра тормозного излучения на ос-

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

Метод РЬ-оценок с вышеописанными модификациями реализован в программном модуле МСРЬ, который является динамически подключаемой библиотекой Подобная реализация дает возможность включать модуль МСРЬ в другие программы, в том числе и в системы планирования

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

ВЫВОДЫ

1 Предложены методы ускорения метода РЬ-оценок потока в точке, основанные на обработке дерева траекторий

2 Показана возможность использования метода РЬ-оценок потока в точке, при расчетах доз в средах с неоднородностью по атомному номеру (для г<13)

3 Метод РЬ-оценок с вышеописанными модификациями реализован в программном модуле МСРЬ, который является динамически подключаемой библиотекой Подобная реализация дает возможность включать модуль МСРЬ в другие программы, в том числе и в системы планирования

4 Предложен метод восстановления энергетического спектра дистанционного радиотерапевтического облучателя на основе пространственного распределения энергопоглощений

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

ПУБЛИКАЦИИ

1 Чупикин Д А "Моделирование дозового распределения в лабиринте медицинского ускорителя методом Монте-Карло", //Сборник научных трудов научной сессии «МИФИ-2004», Москва, 2004

2 Климанов В А, Смирнов В В, Чупикин Д А, "Реализация быстрой схемы расчета дозового распределения, основанной на методе Монте-Карло" // Сборник научных трудов научной сессии «МИФИ-2005», Москва, 2005

3 Журов Ю В, Климанов В А, Петров Д Э, Смирнов В В, Чупикин Д А, "Современное состояние применения численных методов теории переноса при планировании гамма/электронной терапии" // Сборник материалов П-го Евразийского конгресса по медицинской физике и инженерии «Медицинская физика - 2005», Москва, 2005

4 Климанов В А, Чупикин Д А, Донской Е Н "Применение метода РЬ-оценок потока в точке для расчета доз методом Монте-Карло" // Сборник материалов П-го Евразийского конгресса по медицинской физике и инженерии «Медицинская физика - 2005», Москва, 2005

5 Климанов В А, Чупикин Д А, Донской Е Н "Применение метода РЬ-оценок потока в точке для расчета доз методом Монте-Карло" // Медицинская физика №25, Москва, 2005

Автор выражает благодарность Смирнову В В за ряд ценных советов и всестороннюю помощь в работе и подготовке диссертации, Кононову И Н за разъяснения некоторых сложных моментов в работе и помощь в подготовке диссертации, Журову Ю В за критические замечания и за полезные предложения по их устранению, Донскому Е Н за удаленные консультации Особую благодарность хочется выразить научному руководителю Климанову Владимиру Александровичу за бесчисленное множество ценных и полезных советов, а также за неоценимую помощь в подготовке рукописи диссертации

Подписано в печать 14 09 2007 г Исполнено 14 09 2007 г г Печать трафаретная

Заказ № 748 Тираж 80 экз

Типография «11-й ФОРМАТ» ИНН 7726330900 115230, Москва, Варшавское ш, 36 (495) 975-78-56 \vw\v аШоге£ега1 ги

ВВЕДЕНИЕ.

ГЛАВА I. РАЗЛИЧНЫЕ СХЕМЫ И МОДИФИКАЦИИ МЕТОДА МОНТЕ-КАРЛО, ИСПОЛЬЗУЕМЫЕ В ЛУЧЕВОЙ ТЕРАПИИ.

1.1 Схемы метода Монте-Карло, используемые для переноса заряженных частиц.

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

1.2.1 VMC.

1.2.2 VMC++.

1.2.3 ММС.

1.2.4 PEREGRINE.

1.2.5 SMC.

1.2.6 MCDOSE.

1.2.7 DPM.

1.3 Существующие методы повышения эффективности расчетов при расчете поля фотонов в лучевой терапии.

1.3.1 Изменение длины свободного пробега.

1.3.2 Предварительное вычисление источника первых взаимодействий

1.3.3 Русская рулетка и расщепление.

1.3.4 Экспоненциальное преобразование.

1.3.5 Увеличение числа вторичных частиц.

1.3.6 Разделение задач, использование заранее рассчитанных результатов.

1.4 Методы увеличения эффективности расчетов, используемые в лучевой терапии для электронов.

1.4.1 Поглощение электронов внутри области интереса.

1.4.2 Исключение электронов по пробегу.

1.5 Макроскопические теории учета гетерогенностей при вычислении дозы.

1.5.1 Теория Брэгга-Грея.

1.5.2 Теория Спенсера-Аттикса.

1.5.3 Учет больших полостей в поле фотонного излучения.

1.5.4 Теория Берлина для фотонного излучения.

1.6 Выводы.

ГЛАВА II. МЕТОД PL-ОЦЕНОК ПОТОКА В ТОЧКЕ.

2.1 Построение дерева траекторий в однородной бесконечной водной среде.

2.2 PL обработка дерева траектории.

2.2.1 Построение обратной подобной траектории.

2.2.2 Вычисление вклада в дозу в детекторной точке от текущего звена траектории.

2.3 Методы ускорения расчетов при использовании PL-оценок.

2.3.1 Керма-приближение.

2.3.2 Сглаживание электронных треков.

2.4 Учет дивергенции первичного излучения от немононаправленных источников.

2.4.1 Комптоновское (некогерентное) рассеяние.

2.4.2 Зависимость углов сдвига от местоположения точки первого взаимодействия.

2.4.3 Зависимость весового коэффициента от углов рассеяния и сдвига.

2.4.4 Зависимость весового коэффициента от энергии начальных фотонов.

2.5 Учет гетерогенностей в методе PL оценок.

2.5.1 Гетерогенности 1-го рода.

2.5.2 Гетерогенности П-го рода.

2.6 Выводы.

ГЛАВА III. МОДЕЛИРОВАНИЕ ИСТОЧНИКОВ ИЗЛУЧЕНИЯ РАДИОТЕРАПЕВТИЧЕСКИХ УСТАНОВОК.

3.1 Реконструкция эффективного спектра тормозного излучения на основании обработки пространственных дозных распределений

3.2 Обработка файлов фазового пространства радиотерапевтических установок.

3.3 Моделирование источников в методе PL-оценок.

3.3.1 Моделирование энергетического спектра начальных частиц.

3.3.2 Моделирование интенсивности.

3.3.3 Моделирование расходящегося источника.

3.4 Выводы.

Актуальность работы.

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

В настоящее время в таких передовых системах применяются в основном два метода: 1) метод свертки/суперпозиции, использующий дозовые ядра дифференциального тонкого луча; 2) метод тонкого луча в различных модификациях. Эти методы пришли на смену эмпирическим и полуэмпирическим методам. Они обладают достаточной быстротой и обеспечивают в случае однородности среды высокую точность расчета. Однако при наличии гетерогенности в области расчета, а также вблизи поверхности тела пациента погрешность расчета с помощью этих методов сильно возрастает. Принято считать, что единственным методом, который в условиях сложной геометрии расчета может обеспечить рекомендуемую международными комиссиями точность, является вероятностный метод Монте-Карло. Важно, что метод не требует записи соответствующего уравнения переноса, фактически здесь нет проблем с формулировкой дискретной модели. В обоих случаях реализация вычислительного алгоритма не требует упрощений. Метод максимально приспособлен для использования на ЭВМ. Единственным ограничением являются затраты машинного времени для получения результата с необходимой точностью. В настоящее время эти проблемы не являются решенными. Но в последнее время и в этом направлении ведутся целенаправленные работы и имеются определенные достижения.

Таким образом, метод Монте-Карло является наиболее точным методом расчета доз в лучевой терапии, и поэтому, ускорение расчетов путем модификации метода Монте-Карло является актуальной научно-технической задачей.

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

При этом решаются следующие основные задачи:

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

2. Метод PL-оценок потока в точке обобщается на случай расчета доз в неоднородной среде.

3. Разработка корректных математических моделей источников излучения радиотерапевтических ускорителей и облучателей с б0Со с целью использования их в методе PL-оценок.

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

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

Научная новизна настоящей работы заключается в следующем:

1. Впервые предложены методы ускорения метода PL-оценок потока в точке.

2.Впервые показана возможность использования метода PL-оценок потока в точке при расчетах доз в химически неоднородных средах.

3.У совершенствованы модели источников используемых в методе PL-оценок потока в точке.

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

5. Предложен метод восстановления энергетического спектра дистанционного радиотерапевтического облучателя на основе пространственного распределения энергопоглощений.

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

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

Основные положения, выносимые на защиту:

1. Способы ускорения расчетов при использовании метода PL-оценок.

2. Методика учета гетерогенностей при расчете доз в лучевой терапии при использовании метода PL-оценок.

3. Новые модели фотонных источников излучения для использования их в методе PL-оценок.

4. Методика реконструкции эффективного спектра тормозного излучения на основании обработки пространственных дозных распределений, измеренных при вводе в эксплуатацию радиотерапевтических ЛУЭ и кобальтовых облучателей.

Апробация работы. Основные материалы и положения работы докладывались и обсуждались на следующих научных конференциях и семинарах:

Научная конференция в Московском инженерно-физическом институте «Научная сессия МИФИ 2004,2005». - Москва, 2004г, 2005 г.

II Евразийский конгресс по медицинской физике и инженерии «Медицинская физика 2005», Москва, 2005 г.

Публикации.

Материалы диссертации опубликованы в 5 печатных работах: 1 статье в научных журналах, 4 научных трудах и тезисах конференций. 1 статья принята к изданию в журнале «Атомная энергия».

Структура и объем диссертации. Диссертация состоит из введения, трех глав, заключения, списка литературы и приложения. Общий объем составляет 131 страницу печатного текста, включая 50 рисунков. Список литературы включает 42 наименования, из них 30 иностранных.

3.4 Выводы

По результатам проведенных исследований можно сделать следующие выводы: Создан новый метод реконструкции спектра тормозного излучения ЛУЭ на основе обработки типовых дозовых распределений в водном фантоме, используя метод направленного расхождения и

MCPL • Этими

I I I I 1JLI 4- -J — —ы-£+- i i к i IIII 1 1 1 IIII, JLLI. —i—I—1—1— —1—!—1- —1—-IIII IIII 4lLi -4—1—44 — -4—1—4-4 — IIII IIII ±JL± ^44-44- -44-44- llll i i /I i -ti/rt- "Tl/TT" ~тлггт~ i / i i — 1—1-—1—Г8*» ~ГТ~1"Г" ~гт~гг~ IIII IIII ^f-i-M- "~ГТ 1 i 1 i ' IIII -1-htl- ^vi i i IIII llll

-4/4 — 44-+/-t h + -"tl-ГТ" 1 [II IIII ТТГГ и ii~i~ --1-4-4-1— ■ 1 H i ■ "TtTT' -4-1--t 1 -t—t —r-rri- —44-44-~ T ~t— Г T

ТГГГ A M M 1 1 l :ttFb IIII L 1 1 ! IIII 1 Li 1 —f—1—t—t— - 1 4-4 1 j IIII 1 1 1 1 J 1 JJ - 4 —1— 44 — 4 1 4-4 IIII llll 14 L4 —44—4444 44 IIII iiii it-1-1—t-- /Mil IIII IIII —1—1—I—1— —1—1—I—1— 1 1 i ! IIII IIII 1 1 ! 1 I ! i 1 IIII IIII IIII —M~-1— II 1 1 llll llll -t-t-ht- llll llll предварительно рассчитанные методом Монте-Карло дозовые профили для пучков моноэнергетических фотонов. Разработанный метод позволяет также приближенно учитывать зависимость спектра от расстояния до геометрической оси пучка.

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

Дозовые распределения в водном фантоме, рассчитанные методом Монте-Карло с использованием энергетического спектра тормозного излучения и использованием математической модели источника хорошо совпали с экспериментальными данными и расчетами других авторов.

Развита методика моделирования головки радиотерапевтических установок и источников фотонного излучения при проведении расчетов дозовых распределений методом Монте-Карло с использованием PL-оценок потока в точке.

ЗАКЛЮЧЕНИЕ

В настоящей работе исследован метод PL-потока в точке и предложены его улучшения, направленные на увеличение скорости и точности расчетов. Отметим, что скорость счета с использованием PL-оценок имеет почти линейную зависимость от количества детекторных точек. Это позволяет практически сразу рассчитать дозу в любой критической точке за время от 0,04 до 1,5 сек., в зависимости от условий расчета. Использование обычных схем в подобных ситуациях практически невозможно, т.к. зависимость скорости от числа детекторных точек незначительна. Эта особенность делает метод уникальным для проведения расчетов связанных с небольшим количеством точек. Обычно такие расчеты проводятся с исследовательскими или научными целями, где нужно получить глубинные дозовые распределения или профили дозовых распределений.

Были предложены и исследованы методы ускорения расчетов поглощенной дозы методом PL-оценок. Эти методы могут быть использованы как отдельно друг от друга, так и совместно. Выигрыш по времени при применении этих методик может составлять до 2,5 раз. Отметим, что расчет с использованием вышеописанных методик не вносит существенных изменений в сам метод PL-оценок. А это значит, что использование этих двух способов ускорения не изменяет ограничения на область применения метода PL-оценок. Метод возможно комбинировать с любыми другими методами расчета дозы.

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

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

Метод PL-оценок с вышеописанными модификациями реализован в программном модуле MCPL, который является динамически подключаемой библиотекой. Подобная реализация дает возможность включать модуль MCPL в другие программы, в том числе и в системы планирования.

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

1. Климанов В.А., Крянев А.В. Постановка задач оптимизации планирования лучевой терапии. // Медицинская физика, 2000, №7, стр 62-68

2. A. Kryanev, V. Klimanov, S. Klimanov, D. Rubinsky, A. Zrajun. Profiles Beams Optimization Problems for Remote Radiation Therapy as Multicrite-rion Problem. // Chicago-2000 World Congress Proceeding, 2000.,p.83

3. Klimanov V.A., Kryanev A.V., Rubinsky D.A. Numeric Solution for Radiation Therapy Dose Planning Optimization Problem Based on the Pencil Beam Algorithm and Large-Scaled Elements Methods.// Physica Medica, 1999,15, p. 166.

4. Rubinsky D.A., Klimanov V.V., Klimanov S.G., Kryanev A.V. Profiles beams optimization problem for remote radiation therapy as multi-criterion problem. // Proceeding of Biological Engineering and Computing, Croatia, 2001, p. 32-33

5. Klimanov V.V., Klimanov S.G., Kryanev A.V. Formulation and numerical solution of the radiation intensity profile optimization problem as multi-criteria problem using physical and biological objective function// Med. Phys., 2001,11.

6. Климанов B.A., Климанов С.Г., Крянев A.B., Беляков А.И., Головкин Ю.В. Влияние погрешностей во входных данных на величину вероятности контроля над опухолью при оптимальном плане облучения. // Сборник научных трудов МИФИ-2002., стр. 45

7. Климанов В.А., Климанов С.Г., Крянев А.В. Выбор оптимального расположения заданного количества портов облучения с использованием физических целевых функций. // Сборник научных трудов МИФИ-2002., стр. 46.

8. J. Sempau, S.J. Wilderman and A.F. Bielajew DPM, a fast, accurate Monte Carlo code optimized for photon and electron radiotherapy treatment planning dose calculations. // Phys. Med. Biol., 2000,45, p. 2263-2291.

9. Донской E. H. PL-оценка потока через поверхность и их применение в расчетах транспорта фотонов и заряженных частиц методом Монте-Карло.// ВАНТ, Сер: математическое моделирование физических процессов, № 2, стр. 25-29, 1993

10. Kawrakow, М. Fippel Investigation of variance reduction techniques for Monte Carlo photon dose calculation using XVMC. // Phys. Med. Biol, 2000, 45, p. 2163-2184

11. F. Bielajew, D. W. O. Rogers, and A. E. Nahum. Monte Carlo simulation of ion chamber response to 60Co Resolution of anomalies associated with interfaces.// Phys. Med. Biol., 1985,30, p. 419-428.

12. D. W. O. Rogers, A. F. Bielajew, Calculated buildup curves for photon with energies up to 60Co. // Med. Phys., 1985,12, p. 738-744.

13. D. W. O. Rogers and A. F. Bielajew. The use of EGS for Monte Carlo calculations in medical physics. // National Research Council of Canada, (Ottawa, Canada K1A 0R6), Report PXNR-2692, 1984.

14. Ahnesjo A, Andreo P and Brahme A Calculation and application of point spread functions for treatment planning with high energy photon beams // Acta Oncologica, 1987,26, p. 49-55

15. M. Rrmar, D. Nikolic, P. Krstonosic, A simple method for bremsstrahlung spectra reconstruction from transmission measurements// Med. Phys., 2002, 29, p. 932-938.

16. A. Ahnesjo, P. Andreo, Determination of effective bremsstrahlung spectra and electron contamination for photon dose calculations// Phys. Med. Biol., 1989, 34, p. 1451-1464.

17. Тараско M.3., Крамер-Агеев E.A., Тихонов Е.Г. Применение метода направленного расхождения для восстановления спектра быстрых нейтронов. // В сб. "Вопросы дозиметрии и защиты от излучений". Вып.П., М., Атомиздат, 1970, с. 125.

18. В.Ф.Баранов. Дозиметрия электронного излучения.// М., Атомиздат,1974.

19. J. Sempau, S.J. Wilderman and A.F. Bielajew, DPM, a fast, accurate Monte Carlo code optimized for photon and electron radiotherapy treatment planning dose calculations, Phys. Med. Biol., 2000, 45, 2263-2291.

20. M.J. Berger, S.M. Seltzer Bremsstrahlung and Photoneutrons from Thick Tungsten and Tantalum Targets//Phys. Rev., 1970, 2, p. 621-626.

21. A. Bielajew Fundamentals of the Monte-Carlo method for neutral and charged particles transport // University of Michigan, 2000.

22. А. Ф. Аккерман Моделирование траекторий заряженных частиц в веществе//М. Энергоатомиздат, 1991.

23. W.R. Nelson, Н. Hirayama, D. W. О. Rogers The EGS4 Code System // Stanford Linear Accelerator Center, Report SLAC-265,1985.

24. J.A. Halbleib et al ITS Version 3.0: The Integrated TIGER Series of Coupled Electron/Photon Monte Carlo Transport Codes // Sandia Report SAND91-1634 1992

25. J.F. Briesmeister (Editor) MCNP A general Monte Carlo N-particle transport code // LANL (Los Alamos, NM), Report LA-12625-M, 1993

26. F. Salvat et al PENELOPE, an algorithm and computer code for Monte Carlo simulation of electron-photon showers // University of Barcelona preprint, 1996

27. H. Neuenschwander, E.J. Born A macro Monte Carlo method for electron beam dose calculations.// Phys. Med. Biol., 1992,37, p. 107-125

28. H. Neuenschwander et al. MMC-a high-performance Monte Carlo code for electron beam treatment planning.// Phys. Med. Biol., 1995,40, p. 543-574.

29. Kawrakow, M. Fippel, K. Friedrich 3D electron dose calculation using a Voxel based Monte Carlo algorithm (VMC) // Med. Phys., 1996, 23, p. 445 -457

30. Kawrakow Improved modeling of multiple scattering in the Voxel Monte Carlo model. //Med. Phys. 1997,24, p.505-517.

31. M. Fippel Fast Monte Carlo dose calculation for photon beams based onthe VMC electron algorithm. //Med. Phys., 199,26, p. 1466-1475.

32. Kawrakow, M. Fippel Investigation of variance reduction techniques for Monte Carlo photon dose calculation using XVMC. // Phys. Med. Biol.,2000,45, p. 2163-2183.

33. P.J. Keall, P.W. Hoban Super-Monte Carlo: A 3-D electron beam dose calculation algorithm // Med. Phys., 1996, 23, p. 2023

34. C.L. Hartmann-Siantar et al LLNL's PEREGRINE project. // Proceedings of the XH-th Conference on the Use of Computers in Radiotherapy, Madison, Wisconsin, May 27-30,1997, p. 19-22

35. C.-M. Ma, J.S. Li, T. Pawlicki, S. B. Jiang, J. Deng MCDOSE A Monte Carlo dose calculation tool for radiation therapy treatment planning. // Proceedings of the ХШ-th Conference on the Use of Computers in Radiotherapy, Heidelberg, 2000, p. 123 - 125

36. J. Sempau, S.J. Wilderman, A. F. Bielajew: Phys. Med. Biol.,2000, 45, p. 2263

37. Kawrakow, M. Fippel VMC++, a fast MC algorithm for radiation treatment planning. // Proceedings of the XHI-th Conference on the Use of Computers in Radiotherapy, Heidelberg, 2000,p. 126-128

38. Донской E. H., Климанов B.A., Чупикин Д.А. Применение метода PL-оценок потока в точке для расчета доз методом Монте-Карло. // Медицинская физика, 2005, №25., стр. 15-22

39. Кейз К., Цвайфель П., Линейная теория переноса.// М.- Мир, 1972.

40. Кольчужкин A.M., Учайкин В.В. Введение в теорию прохождения частиц через вещество. // М,- Атомиздат, 1978.

41. Спанье Дж., Гелбард Э. Метод Монте-Карло и задачи переноса нейтронов.//М.-Атомиздат, 1972.

42. Донской Е. Н. Методика и программа PATIENTPL для расчета дозовых распределений в теле пациента методом Монте-Карло с использованием PL-оценок потока в точке.// Отчет по проекту МНТЦ №107999.