Содержание к диссертации
Введение
Глава 1. Плоские волны в обработке и восстановлении изображений 20
1.1. Восстановление на основе плоских волн 28
1.1.1. Равномерно распределенные направления проекций 32
1.1.2. Проекции как лучевые и полосовые суммы 39
1.1.3. Восстановление по модельным и реальным данным 42
1.2. Алгебраические методы восстановления 47
1.2.1. Вычислительные эксперименты 54
1.3. Информативность набора проекций 60
1.3.1. Декомпозиция текстур в сумму плоских волн 61
1.3.2. Оптимизация схем облучения в радиотерапии 67
1.4. Характеристические конечные элементы 73
1.5. Выводы 75
Глава2.Задача Радонавполосе 76
2.1. Матрицы Грама и проекционные матрицы 81
2.1.1. Структура матриц Грама 92
2.1.2. Вычислительные эксперименты 97
2.2. Веерные плоские волны 101
2.2.1. Пара линейных детекторов 103
2.2.2. Вычислительные эксперименты 110
2.2.3. Пара криволинейных детекторов 111
2.2.4. Вычислительные эксперименты 118
2.3. Интегральные данные по траекториям рассеяния фотонов 121
2.4. Выводы 124
Глава 3. Лучевое преобразование с послойной сверткой 125
3.1. Метод обратного проецирования с послойной фильтрацией 131
3.2. Комптоновское рассеяние в эмиссионной томографии 140
3.2.1. Идеализированная геометрическая модель 144
3.2.2. Статистическое моделирование 148
3.3. Послойная сверточная модель комптоновского рассеяния 154
3.4. Послойные искажения в электронной микроскопии 159
3.4.1. Вычислительные эксперименты 163
3.5. Выводы 170
Глава4. Лучевое сферическое преобразование 171
4.1. Конусы зрения и их трассировка 181
4.2. Гномонические проекции и граничные кубы 192
4.3. Конструктивный алгоритм вычисления многоугольника 197
4.3.1. Вычислительные эксперименты 205
4.4. Численные реконструкционные эксперименты 207
4.4.1. Малоразмерный эксперимент 209
4.4.2. Большеразмерный эксперимент 213
4.4.3. Применение многосеточных методов 214
4.5. Математическое моделирование преобразования Радона 216
4.6. Выводы 223
Заключение 224
Приложениеi 239
Приложение ii
- Проекции как лучевые и полосовые суммы
- Вычислительные эксперименты
- Комптоновское рассеяние в эмиссионной томографии
- Конструктивный алгоритм вычисления многоугольника
Введение к работе
Актуальность работы. Во многих областях науки и техники, таких как медицина, геофизика, электронная микроскопия, кристаллофизика, астрофизика, промышленная дефектоскопия, диагностика плазмы и других, возникает проблема определения внутренней структуры объекта наиболее удобным для восприятия образом, т.е. в виде изображения. Для ее решения во многих случаях неприемлемы прямые методы исследования, связанные с разрушением объекта, и поэтому создаются специальные системы дистанционного получения данных. Физический принцип таких систем состоит в воздействии некоторого физического процесса известной природы (излучение, волновое поле и т.д.) и последующей регистрации отклика этого процесса на объект. Математическая модель используемого физического воздействия является основой построения численных методов обработки регистрируемых данных и разработки программных комплексов, дающих в результате изображение внутренних структур.
В работе исследуется такой вид дистанционно получаемых данных, как интегральные проекции; особенно актуально их применение в томографии. Хотя математические основы были заложены в начале ХХ века в работах по интегральной геометрии (Г. Минковский 1905 г., П. Функ 1913 г., Й. Радон 1917 г.), для полного решения проблемы реконструкции потребовались современные средства соединения цифровой техники, систем программирования, обработки изображений и численных методов обработки экспериментальных данных в единую систему вычислительной томографии. Математический аппарат многомерных интегральных преобразований лежит также в основе цифровой обработки сигналов и изображений. Известная взаимосвязь преобразований Фурье и Радона приводит к применению (лучевого) преобразования Радона в задачах анализа структур на изображениях и созданию новых быстрых многомерных алгоритмов обработки сигналов. Широкое распространение получил термин “обработка изображений в пространстве Радона”. Современный уровень научных исследований требует эффективного математического моделирования и программного обеспечения для исследования влияния множества факторов на качество восстанавливаемых изображений и решения конретных прикладных задач.
Целью диссертационной работы являются: разработка математических моделей, численных методов, алгоритмов и программ для решения задач восстановления и обработки изображений, формулируемых в терминах лучевого преобразования в двух-, трех- и четырехмерных постановках, а также проведение вычислительных экспериментов для проверки созданных моделей и алгоритмов.
В соответствии с поставленными целями в диссертационной работе решены следующие (и сопутствующие им) задачи:
-
Исследование численных методов и алгоритмов обращения двумерного лучевого преобразования с помощью разложения восстанавливаемой функции в суперпозицию плоских волн, и применение созданных новых алгоритмов в задачах томографии и анализа регулярных текстур в обработке изображений.
-
Развитие аналитической и дискретной моделей преобразования Радона в полосе на основе оригинальных методов представления изображений в базисах веерных плоских волн и характеристических конечных элементов для двумерных и трехмерных задач томографии.
-
Моделирование процесса формирования проекций в позитронной эмиссионной томографии с учетом однократного комптоновского рассеяния, в частности в виде нового интегрального преобразования с послойной сверткой; верификация созданного алгоритма методом статистического моделирования с использованием метода Монте-Карло. Разработка численных методов, алгоритмов и программ обращения трехмерного лучевого преобразования с послойными искажениями, имеющими характер свертки.
-
Дискретизация четырехмерного лучевого сферического преобразования в виде квазирегулярной сетки на единичной сфере в четырехмерном пространстве и создание нового конструктивного некомбинаторного быстрого алгоритма вычисления центральных сечений многомерного куба.
-
Апробирование и адаптация созданных алгоритмов, тестирование разработанных программных комплексов и модельных фантомов на сетках больших размеров, необходимых в настоящее время для достижения высокого разрешения томографической реконструции при использовании алгебраических методов.
Методы исследования. Основные результаты работы получены с использованием интегральной геометрии, численных методов и геометрической томографии, а также методов цифровой обработки изображений, математического моделирования, асимптотических методов и статистического моделирования. Для численного моделирования и программной реализации разработанных алгоритмов использовались методы прикладного программирования.
Научная новизна.
Разработан новый алгоритм восстановления изображений по конечному набору проекций в постановке лучевого двумерного преобразования и представления изображений в виде суперпозиции плоских волн произвольных ориентаций. В отличие от существующих методов реконструкции, впервые в алгоритм восстановления включена управляющая матрица, содержащая информацию о геометрии направлений просвечивания.
Впервые теоретически обосновано и в численных экспериментах апробировано использование информативных проекций для анализа изображений регулярных текстур и оптимальных доз просвечивания.
Впервые вычислены проекционные матрицы дискретной прямой задачи Радона в полосе для двумерного и трехмерного случаев, используя геометрию характеристических конечных элементов с носителями в форме призм. Разработан быстрый алгоритм реконструкции по веерным данным на основе нового интегро-дифферен-циального уравнения.
Получена новая асимптотическая формула обращения трехмерного лучевого преобразования с послойными искажениями в виде двумерной свертки. Разработан и апробирован соответствующий численный алгоритм обращения с достаточно точным поведением асимптотики как на высоких, так и на низких частотах.
Впервые разработана математическая модель формирования проекционных данных в позитронной эмиссионной томографии с учетом однократного комптоновского рассеяния. Показано, что при определенных условиях, проекционные данные описываются лучевым преобразованием с послойными искажениями. Проведена верификация модели методом статистического моделирования.
В задаче дискретизации четырехмерного лучевого сферического преобразования впервые разработаны численный метод трассировки окружностей и новый быстрый (без комбинаторного перебора общих точек граней и секущей плоскости) алгоритм вычисления центральных сечений многомерного куба.
Полученные новые теоретические результаты и разработанныепро-граммные комплексы, использующие созданные новые быстрые алгоритмы, впервые позволили провести большие вычислительные томографические эксперименты алгебраических реконструкций с изображениями размером порядка 103 103 103.
Достоверность полученных результатов и выводов подтверждается доказательствами теоретических утверждений, анализом разработанных численных алгоритмов и проведением численных экспериментов. В работе широко используется сравнение реконструкций с эталонным известным тестовым фантомом, а также сравнение с реконструкциями, полученными известными классическими методами. Алгоритм вычисления комптоновских проекций верифицирован математическим моделированием с использованием метода Монте-Карло.
Теоретическая и практическая ценность работы. Диссертация носит теоретический характер с направленностью на задачи, имеющие актуальные физические постановки, сформулированные как математические модели с определенной идеализацией физических явлений. Это обстоятельство требует вычислительного эксперимента для верификации моделей, что в свою очередь мотивирует разработку программных комплексов. Полученные в диссертации результаты и программы могут быть применены на практике в рентгеновской компьютерной томографии, межскважинном сейсмическом просвечивании, электронной микроскопии макромолекул, позитронной эмиссионной томографии, радиационной терапии и анализе текстурных изображений в машинном зрении.
Основные положения, выносимые на защиту
-
Разработанные численные методы и алгоритмы обращения двумерного лучевого преобразования на основе метода суперпозиции плоских волн, учитывающего геометрию просвечивания и информативность проекционных данных, а также основанные на этой методике алгоритмы аппроксимации текстур в обработке изображений.
-
Созданные двух- и трехмерные модели веерного преобразования в задаче Радона в полосе на основе разложений в базисах веерных плоских волн и характеристических конечных элементов.
-
Новое трехмерное интегральное преобразование, моделирующее формирование проекций в позитронной эмиссионной томографии с учетом однократного комптоновского рассеяния, его сведение к лучевому преобразованияю с послойными искажениями в виде свертки, и численные методы, алгоритмы и программы восстановления изображений по проекционным данным нового преобразвания.
-
Разработанные новые алгоритмы численной реализации дискретного четырехмерного лучевого сферического преобразования в задаче итерационного восстановления функции распределения ориентаций кристаллических текстур по данным дифракционных изображений, основанные на новом быстром методе построения сечений многомерного
куба.
5. Программный комплекс томографической реконструкции и визуализации, апробированный в больших вычислительных экспериментах с изображениями функции распределения ориентаций размером порядка 103 103 103 элементов.
Апробация работы. Основные результаты диссертации докладывались и обсуждались на симпозиумах по вычислительной томографии (Самара 1985, Ташкент 1989, Москва 1991, Новосибирск 1993), международной конференции ”Обработка изображений и дистанционные исследования” (Новосибирск 1990), международных конференциях ”Обратные и некорректные задачи математической физики” (Новосибирск 2002, 2007, 2012), Российско-немецком семинаре ”Распознавание образов и понимание изображений” (Нижний Новгород 2011), XXIII Российской конференции по электронной микроскопии (Черноголовка, 2010), IX Международном научном конгрессе и выставке “ИнтерЭкспо Гео-Сибирь-2013”, международном симпозиуме “Вычислительная томография для промышленных приложений” (Берлин 2003), международных симпозиумах “Биомедицинская визуализация” (IEEE Int. Symp. on Biomedical Imaging), (2008, 2010), международных конференциях “Медицинская визуализация” (IEEE Nuclear Sci. Symp. and Medical Imaging Conf.), (1999, 2000, 2002, 2003, 2004, 2005, 2006, 2008), международных конференциях “Обработка изображений” (Int. Conf. on Image Processing, ICIP), (1996, 1997).
На различных стадиях выполнения работа обсуждалась на семинарах, руководимых ведущими специалистами, в российских и зарубежных институтах и университетах: Институт математики им. С.Л. Соболева СО РАН (Романов В.Г., Аниконов Д.С.), Институт автоматики и электрометрии СО РАН (Потатуркин О.И., Киричук В.С.), Институт физики полупроводников им. А.В. Ржанова СО РАН (Латышев А.В.), Институт прикладной математики, Мюнстерский университет, Германия (F. Natterer), Группа цифровой визуализации, Городской университет Нью Йорка, США (G.T. Herman), Группа обработки медицинских изображений, Пенсильванский университет, США (R. Lewitt, S. Matej), Факультет информатики и математического моделирования, Датский технический университет, Дания (P.C. Hansen).
Личный вклад автора. Основные научные и практические результаты диссертации получены автором лично. Вывод формул в работе [3] принадлежит автору, соавторам И.П. Яровенко и И.В. Прохорову принадлежит численная верификация методом статистического моделирования. Доказательства в работе [13] принадлежит автору, а чис-
ленная проверка принадлежит соавтору В.В. Пикалову. Из печатных работ, опубликованных диссертантом в соавторстве, в диссертацию вошли только те результаты, в получении которых он принял непосредственное творческое участие; ему принадлежат ключевые идеи доказательств и большинство программных кодов. Конфликт интересов с соавторами отсутствует.
Публикации. Результаты исследований по теме диссертационной работы опубликованы в 21 печатных работах в изданиях списка ВАК.
Структура и объем работы.
Проекции как лучевые и полосовые суммы
1. Разработанные численные методы и алгоритмы обращения двумерного лучевого преобразования на основе метода суперпозиции плоских волн, учитывающего геометрию просвечивания и информативность проекционных данных, а также основанные на этой методике алгоритмы аппроксимации текстур в обработке изображений.
2. Созданные двух- и трехмерные модели веерного преобразования в задаче Радона в полосе на основе разложений в базисах веерных плоских волн и характеристических конечных элементов.
3. Новое трехмерное интегральное преобразование, моделирующее формирование проекций в позитронной эмиссионной томографии с учетом однократного комптонов-ского рассеяния, его сведение к лучевому преобразованияю с послойными искажениями в виде свертки, и численные методы, алгоритмы и программы восстановления изображений по проекционным данным нового преобразвания.
4. Разработанные новые алгоритмы численной реализации дискретного четырехмерного лучевого сферического преобразования в задаче итерационного восстановления функции распределения ориентаций кристаллических текстур по данным дифракционных изображений, основанные на новом быстром методе построения сечений многомерного куба.
5. Программный комплекс томографической реконструкции и визуализации, апробированный в больших вычислительных экспериментах с изображениями функции распределения ориентаций размером порядка 103 103 103 элементов.
Апробация работы. Основные результаты диссертации докладывались и обсуждались на научных конференциях в России и за рубежом: - симпозиумах по вычислительной томографии (Самара 1985, Ташкент 1989, Москва 1991, Новосибирск 1993), -международной конференции ”Обработка изображений и дистанционные исследования” (Новосибирск 1990), - международных конференциях ”Обратные и некорректные задачи математической физики” (Новосибирск 2002, 2007, 2012), - 8-м Открытом российско-немецком семинаре ”Распознавание образов и понимание изображений” (Нижний Новгород 2011), - XXIII Российской конференции по электронной микроскопии (Черноголовка, Московская область 2010), - международном симпозиуме ”Вычислительная томография для промышленных приложений” (International Symposium on Computed Tomography and Image Processing for Industrial Radiology) (Берлин 2003), - международных симпозиумах ”Биомедицинская визуализация” (IEEE International Symposium on Biomedical Imaging), (2008, 2010), - международных конференциях ”Медицинская визуализация” (IEEE Nuclear Science Symposium and Medical Imaging Conference), (1999, 2000, 2002, 2003, 2004, 2005, 2006, 2008).
На различных стадиях выполнения работа обсуждалась на семинарах, руководимых ведущими специалистами, в российских и зарубежных институтах и университетах: - Институт математики им. С.Л. Соболева СО РАН (Романов В.Г., Аниконов Д.С.) - Институт автоматики и электрометрии СО РАН (Потатуркин О.И., Лихачев А.В.) - Институт физики полупроводников им. А.В. Ржанова СО РАН (Латышев А.В.) - Лаборатория обработи изображений, Институт вычислительной математики и математической геофизики СО РАН (Пяткин В.П.), - Научно-исследовательский институт интроскопии Томского политехнического института (Чахлов Л.В.), - Институт измерений Словацкой академии наук, Братислава, Словакия (I. Bajla), - Институт прикладной математики, Мюнстерский университет, Германия (F. Natterer), - Группа обработки медицинских изображений и сигналов, Гентский университет, Бельгия (I. Lemahieu), -Группа обработки медицинских изображений, Пенсильванский университет, США (R. Lewitt, S. Matej), - Группа цифровой визуализации, Городской университет Нью Йорка, США (G.T.Herman), - Факультет информатики и математического моделирования, Датский технический университет, Копенгаген, Дания (P.C.Hansen).
Личный вклад автора. Основные научные и практические результаты диссертации получены автором лично. Вычисления в разделе 1.2.1, использующие пакет ТОПАЗ, проведены В.В. Пикаловым в совместной статье [128]. Статистическое моделирование в разделе 3.2.2 проведено И.П. Яровенко и И.В. Прохоровым в совместной статье [42]. Вычисления в разделе 4.4.3, использующие графические ускорители, проведены с участием сотрудников Датского технического университета в совместной статье [176]. Из печатных работ, опубликованных диссертантом в соавторстве, в диссертацию вошли только те результаты, в получении которых он принял непосредственное творческое участие.
Публикации. Результаты исследований по теме диссертационной работы опубликованы в 21 печатных работах из списка ВАК.
Структура и объем работы. Диссертация состоит из введения, четырех глав, заключения, списка литературы из 190 наименований и двух приложений. Содержание основного текста диссертации изложено на 238 страницах, содержит 60 иллюстраций и 8 таблиц. Основное содержание работы
Во Введении обоснована актуальность проблем, рассматриваемых в диссертации, определены основные цели и задачи исследования, показана его научная новизна и практическая ценность, сформулированы выносимые на защиту положения и представлен ретроспективный обзор содержания работы.
Главы имеют однотипную структуру: физическая и математическая постановка задачи, обзор существующих в литературе подходов, обоснование предлагаемого решения, алгоритмическая реализация в виде программного комплекса, вычислительные эксперименты на тестовых и реальных данных, выводы. Главы 1 - 4 рассматривают отдельные виды лучевого преобразования (Радона): глава 1 – двумерное лучевое преобразование с параллельной организацией лучей; глава 2 – двумерное и трехмерное лучевое преобразование Радона с веерной организацией лучей в полосе; глава 3 – трехмерное лучевое преобразование с послойными искажениями; глава 4 – четырехмерное лучевое сферическое преобразование.
Вычислительные эксперименты
Функции hi(x cos uji+у sin си,;) являются плоскими волнами и представляют собой ориентированные, или направленные структуры. В англоамериканской литературе широко используется термин “ridge functions” - “хребтовые функции”. Заметим, что в отечественной литературе имеются работы (Вострецов Б.А. и Крейнес М.А., 1961 г.) [17], где доказаны теоремы существования приближения непрерывных функций суперпозициями плоских волн в более общем многомерном случае. Теорема Логана-Шеппа есть утверждение существования плоских волн в выбранных направлениях. В разделе 1.1 будет выведен конструктивный алгоритм вычисления суперпозиции плоских волн Нш(х,у) по произвольному набору полных проекций Rwf. В случае равномерно распределенных углов вычисление становится очень эффективным в силу наличия аналитической формулы для обратных матриц, участвующих в реконструкции. Для минимального нормального решения Нш справедливо равенство из которого следует, что знание нормы Д 2 полезно в оценивании точности восстановления (1.15) по данному набору проекций.
Преобразование Радона является одним из главных математических инструментов вычислительной томографии. История его применения в интегральной геометрии, математической физике, геофизике и других областях науки достаточно полно представлена в [102]. Помимо томографии, преобразование Радона нашло применение в обработке изображений, где оно используется для нахождения линий (границ, разрывов) на цифровых снимках [92], [36], [37], [38]. Линии на изображении могут быть выделены, когда проекционные данные, вычисленные на этом изображении, включают интегралы вдоль прямых, совпадающих с искомыми линейчатыми структурами, поскольку эти интегралы аккумулируют одинаковые (приблизительно) значения вдоль структуры. Присутствие линии проявляется в локальных пиках или провалах на синограммах, или визуализированном изображении преобразования Радона. Та же идея направленной структуры содержится в операторе обратного проецирования и тесно связанном с ним понятии плоской волны. Распознавание линии предполагает ее доминирующее положение на изображении по амплитуде и по длине. Слабые и короткие структуры могут теряться в силу наличия шумов и маскирующего эффекта других структур.
В этой главе плоские волны и их ориентация рассматриваются как индикаторы для определения протяженных структур в обработке изображений. Некоторые из них, такие как трещиноватые дефекты в металлических деталях, могут моделироваться как одна плоская волна. Другие, такие как текстура ткани при контроле текстиля, обладают несколькими выделенными направлениями и объект моделируется как сумма плоских волн. Дефект выражается в нарушении направленности текстильной структуры.
В радиографической дефектоскопии многие характеристики объекта известны – например, физические параметры деталей изделия и его форма. Описания дефектов также известны, среди них протяженные трещины представляют наибольшую опасность и выявляются при дефектоскопии (например, сварочных швов стальных труб), как правило, по малому числу (1-2) проекций. Чувствительность метода дефектоскопии зависит от ориентации возможного трещиноватого дефекта. Если трещина параллельна направлению просвечивания, она обнаруживается легче по сравнению с трещиной, перпендикулярной к направлению излучения. Возникает вопрос о надежности определения анизотропии структуры по малому числу проекций и какое направление просвечивания оптимально, или наиболее информативно.
В приложениях компьютерного зрения [15], таких как автоматический контроль материалов, обнаружение дефектов в текстурах выделилось в отдельную область анализа изображений. Объектами анализа становятся текстиль и другие ткани, денежные банкноты, поверхность дерева, металлов, бумаги. Текстурные дефекты являются местами локального нарушения общей однородности, или структуры. Многие структуры выглядят как сложная суперпозиция протяженных структур [163]. Представляет интерес проблема: будет ли декомпозиция текстуры в сумму определенных (наиболее информативных) плоских волн иметь смысл для структурных (не стохастических) текстур и каковы численные аспекты их выделения. Для рассмотрения подобных проблем оказался полезным полученный в разделе 1.1 алгоритм разложения двумерных изображений в сумму плоских волн. В разделе 1.2 исследуется применение алгебраических методов восстановления изображений для более общих геометрий сбора данных, в частности, для моделирования полосовых сумм (Рис. 1.1 б). Введением в разделе 1.3 количественного критерия качества интегральных данных, называемого информативностью, обеспечивается выделение текстур на изображениях и решаются некоторые известные теоретические и практические задачи томографии и обработки изображений. Метод характеристических элементов, как дискретный аналог метода плоских волн (Рис. 1.1 в) обсуждается в разделе 1.4, создавая базис для рассмотрения в главе 2 задач восстановления по неполным данным. Результаты главы 1 диссертации опубликованы автором в работах
Комптоновское рассеяние в эмиссионной томографии
Для проверки соответствия формулы (3.73) реальному физическому процессу проведено имитационное моделирование на основе метода Монте-Карло [30]. Идея метода основана на представлении переноса излучения в среде в виде случайного процесса - движения фотонов через среду, используя известные законы взаимодействия излучения с веществом. В результате такого моделировния накапливается необходимое количество статистической информации для того, чтобы построить оценку величины в. Такой подход позволяет имитировать сигнал, который регистрирует реальный томограф.
Алгоритм, используемый в работе, состоит из следующих пунктов:
1. Розыгрыш точки аннигиляции фотона. Пусть источник с эмиссионной активностью / распределен в некоторой области Go = D(f) С D(/i). На первом этапе разыгрываем точку г рождения позитрона равномерно внутри Go. Будем считать, что позитрон тут же аннигилирует, так что точки рождения и аннигиляции позитрона совпадают
2. Розыгрыш направления движения фотонов, родившихся в резуль тате аннигиляции.
Здесь мы ограничимся рассмотрением наиболее вероятного случая, когда в результате аннигиляции рождаются два фотона. Для применяемых в медицине фармпрепаратов отклонение родившихся в результате аннигиляции фотонов от антиколлинеарности, как правило, не превышает 10-3 рад. [165]. В связи с этим, будем считать, что аннигиляционные фотоны распространяются строго в диаметрально противоположных направлениях. Вследствие изотропии первоначальных направлений движения фотона, равномерно на сфере разыгрывам направление ш = ( х і,бо 2, з), в котором полетит первый из образовавшихся фотонов. Далее отслеживаем историю фотона, движущегося в направлении со.
3. Розыгрыш свободного пробега фотона в веществе.
Свободный пробег разыгрывается с учетом экспоненциального закона уменьшения с расстоянием интенсивности потока фотонов. Для нахождения длины свободного пробега в кусочно-неоднородной среде применялся алгоритм, приведенный в [30]. Если в результате свободного пробега фотон покинул геометрию томографа, то такая траектория отбрасывается и выполняется переход к п.1. В противном случае разыгрывается тип взаимодействия фотона с веществом.
4. Розыгрыш взаимодействия фотона с веществом.
Для характерного в позитронной эмиссионной томографии диапазона энергий среди всех видов взаимодействия излучения с веществом преобладают фотоэлектрическое поглощение, комптоновское рассеяние и рассеяние по Рэлею [62]. Для определения типа взаимодействия излучения с веществом вычисляются вероятности поглощения фотона, рассеяния по Комптону и где ца - коэффициент поглощения, /ІД, цс - коэффициенты рассеяния по Рэлею и Комптону соответственно.
Здесь и далее по тексту символами а, [5 будем обозначать независимые реализации случайной величины, равномерно распределенной в диапазоне [0,1].
Если а Ра, то фотон поглотился, и такая траектория дальше не отслеживается. Если а Ра, то фотон рассеется и необходимо определить тип рассеяния. При /3 Рг фотон рассеется по закону Рэлея, в противном случае фотон рассеется по закону Комптона.
5. Рассеяние фотона по Рэлею. Если фотон рассеялся по Рэлею, то разыгрываем новое направление движения u/ = (OJ[,UJ2,UJ 3), которое связано с первоначальным направлением ш следующими соотношениями: азимутальный угол рассеяния ф = 2тг(3. Формулы (3.77)-(3.79) для нахождения нового направления справедливы при х з ±1. В противном случае После нахождения нового направления движения возвращаемся назад к п. 3.
6. Рассеяние фотона по Комптону. При рассеянии фотона по закону Комптона происходит изменение не только направления, но и энергии фотона. Вероятность рассеяному фотону иметь энергию от Е до Е в зависимости от значения случайного числа а определяется из уравнения где Emm = Е/{1 + 2E/EQ) - минимальная энергия, которую может приобрести фотон в результате рассеяния, dac/dE - дифференциальное сечение комптоновского рассеяния в интервал энергий от Е до Е + dE, определяемое комбинацией сечения Кляйна - Нишины - Тамма и функции некогерентного рассеяния [62], позволяющей учитывать влияние связи электронов в атоме вещества на процесс рассеяния.
Решение уравнения (3.82) - достаточно трудоемкий процесс. Для ускорения в расчетах использовалась заранее насчитанная двумерная таблица, содержащая решения уравнения (3.82) с заданной точностью на равномерных сетках по Е и а. В промежуточных точках применялась линейная интерполяция. Найденное значение энергии рассеянного фотона позволяет определить косинус угла рассеяния: + Ео/Е - Е0/Е : где EQ - энергия покоя электрона. Далее, разыгрывая равновероятно в диапазоне от 0 до 27Г азимутальный угол рассеяния и применяя формулы (3.77)-(3.80), находим новое направление движения фотона. После этого возвращаемся к п. 3.
7. Регистрация фотона детектором. Если в результате трассировки фотон попал в детектор, запоминаем номер детектора и переходим к отслеживанию фотона в направлении —со. В случае, когда оба фотона, родившиеся в результате аннигиляции, попали соответственно в детекторы А и , и их энергия попадает в заданный диапазон [а,а + Да], зависящий от в, оценка функции в увеличивается на величину интенсивности активного вещества в точке аннигиляции позитрона. Итоговая оценка величины fB имеет вид:
Конструктивный алгоритм вычисления многоугольника
Будем говорить, что абсолютная синусоида \sx(t)(t)\ мажорирует остальные абсолютные синусоиды \sk(t)\ в точке t. Если отобразить все абсолютные синусоиды совместно на графике с абсциссой t, то результатом мажорирования будет кривая /І, состоящая из мажорирующих участков отдельных абсолютных синусоид:
Назовем эту кривую мажорирующей, или огибающей. Тогда можно видеть, что если в точке t абсолютная синусоида с номером K(t) мажорирует остальные абсолютные синусоиды, то точка S(t) = Pcost + Qsint окружности C(P,Q) принадлежит конусу зрения с номером K(t), если sK(t)(t) 0, и конусу зрения с номером {-K{t)) если sm{t) 0. Покажем, что случай n(t) = 0 невозможен при любых t.
Утверждение 4.2. Функция ц является строго положительной.
Доказательство. Доказательство проводится от противного. Предположим, в точке to произошло мажорирование нулем всех абсолютных синусоид. Тогда \pkcosto + sin to = 0 в точке to для всех синусоид к = 1,..., iV. Возведя в квадрат и сложив эти N равенств, получаем
Конструктивный алгоритм вычисления вершин состоит из следующих построений. Введем систему координат (р, q) и отобразим на ней точки с координатами (pk,qk),k = 1,N, являщимися компонентами векторов Р = (php2,...,pN) и Q = (qhq2,...,qN). Обозначим эти точки Тк = [РкАк]) а их центрально симметричные антиподы Т_к = {—pk)—qk)-Когда угловой параметр пробегает t Є [0, 27г), координаты точек окружности С(Р, Q) описываются синусоидами где ak определяется из равенств cos(a ) = Pkl\Jp\ + Qk и sin(o ) = Як/у/Рк + Як. Мы видим, что \OTk\\cos{t-ak)\ есть модуль синусоиды sk{t), где \ОТк\ = y/Pk + qjj. обозначает длину вектора ОТк. Перейдем от абсолютных синусоид sfc() к их эквивалентному представлению на плоскости {p,q) в виде пучков окружностей. Введем для этого окружности А±к, построенные на отрезках ОТ±к как на диаметрах. Следующее утверждение является очевидным следствием определения окружностей.
Утверждение 4.3. Геометрическое место точек с полярными координатами (r,t) где г = \OTk\\cos(t - ак)\ и t Є [0,2тг), составляют две центрально симметричные окружности Ак и А_к) т.е.
Следовательно, объединение окружностей Ак U А_к является эквивалентной геометрической моделью абсолютной синусоиды \sk{t)\, а полярная координата г точки (r,t) Є Ак U А_к) лежащей на этих окружностях, равна значению абсолютной синусоиды \sk(t)\. Обозначим объедине 200 ние окружностей А±к как
Таким образом, мы получили двумерную интерпретацию абсолютных синусоид \sk\ в виде пучка окружностей А. Очевидно, что мажорирующие части абсолютных синусоид соответствуют дугам окружностей, составляющим замкнутую кривую, огибающую множество А. Из утверждения 4.2 следует, что начало координат не принадлежит огибающей множества А, и значит, любой луч с началом в начале координат пересекает какую-то из окружностей множества А в точке, отличной от начала координат. Вычислим выпуклую оболочку множества точек {Т±к}
Очевидно, точки Тк и Т_к принадлежат (или не принадлежат) множеству Н совместно. Сформулируем основной результат.
Теорема 4.4. Конус зрения с номером ±к имеет ненулевое пересечение с окружностью C(P,Q) = {Рcos(t) + Qsin(J), t Є [0,2тг)} тогда и только тогда, когда точка Т±к принадлежит выпуклой оболочке Н.
Доказательство. Доказываемое утверждение эквивалентно тому, что окружность Ак имеет хотя бы одну мажорирующую точку тогда и только тогда, когда Тк Є Н. Пусть Тк Є Н для некоторого к. Покажем, что окружность Ак имеет мажорирующую точку. Проведем линию С ортогонально ОТк (Рис. 4.8). Обозначим полуплоскость, определямую линией , и содержащую начало координат О = (0, 0) Є С0 и саму линию С С С0. Выберем из центрально симметричного выпуклого многоугольника Н точ-ки Тг Є Н и Tj Є Н - ближайшие соседние с точкой Тк слева и справа (смежные точки). В силу выпуклости Н, возможны только два различных
Рассмотрим первый случай (4.56). Построим окружности Аг и А3 с диаметрами ОТг и OTj соответственно (Рис. 4.8 (а)). Заметим, что углы ТгОТк и TkOTj не могут одновременно быть тупыми, иначе точки антиподы T_j и T-j окажутся ближе к точке Тк. Пусть угол T,OTfe - острый. Через точки Тг и Тк проведем линию ТгТк. Обозначим X О точку пересечения окружностей Аг и Ак. Имеем X є ГгГ& и ОХ _1_ ТгТк. Из Afe с следует, что X Є С0. Пусть точка S - пересечение окружности Aj и диаметра ОТк. Покажем, что \OS\ \ОТк\. Для этого продолжим линию ОТг до пересечения с линией С в точке Q. Из подобия треугольников QTkO и TtSO и неравенства \ОТг\ \OQ\ следует, что \OS\ \OTk\. Они подобны так как ZQTkO = ZTnSO = 7г/2. По построению имеем \ОТп\ \OQ\ и поэтому OS! ОТ . Из этого заключаем, что точки дуги ХТк С Ак мажорируют точки дуги XS С Aj. Аналогично для точки Y ф О, Y = Ak П Aj, дуга TkY мажорирует часть дуги OY окружности Aj. При этом точки X и У разделены линией ОТ& и поэтому дуга XTkY состоит из точек окружности Ak, мажорирующих обе окружности Aj и Aj. Следовательно, конус зрения, соответствующий кубу Fjfe, посещается окружностью C{P,Q). В предельном случае, когда TnTj Є С, имеем X = Tk = Y и тогда единственная точка Тк является мажорирующей.
Во втором случае (4.57), смежные точки ТІ и 7} из Н разделены линией С (Рис. 4.8 (б)). Точка 7} должна располагаться между линиями С и TjTk и не может находиться выше линии T{Tk, так как это нарушило бы условие Tk Є Н. Заметим, что точка Т& попадает внутрь окружности Aj, построенной на диаметре ОТг Точка Е пересечения окружностей Ак и А3 находится на дуге ХТк) следовательно, дуга ХЕ состоит из мажорирующих точек. Отсюда следует, что на окружности Ак существует мажорирующая дуга, хотя и не содержащая точку Tk. Тогда конус зрения с основанием F& является посещаемым. Мы показали, что принадлежность точки Tk выпуклой оболочке Н всегда указывает на то, что некоторая дуга окружности Ak является мажорантной и, таким образом, конус зрения с номером к является посещаемым.
