автореферат диссертации по информатике, вычислительной технике и управлению, 05.13.18, диссертация на тему:Томография по неполным и искаженным данным

доктора технических наук
Лихачев, Алексей Валерьевич
город
Новосибирск
год
2010
специальность ВАК РФ
05.13.18
цена
450 рублей
Диссертация по информатике, вычислительной технике и управлению на тему «Томография по неполным и искаженным данным»

Автореферат диссертации по теме "Томография по неполным и искаженным данным"

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

Лихачев Алексей Валерьевич

Томография по неполным и искажённым данным

05.13.18 "Математическое моделирование, численные методы и комплексы программ"

АВТОРЕФЕРАТ

диссертации на соискание учёной степени доктора технических наук

Новосибирск-2011 9 ИЮН 2011

4849419

Работа выполнена в Учреждении Российской академии наук Институте автоматики и электрометрии Сибирского отделения РАН

Официальные оппоненты: доктор физико-математических наук, профессор

Воскобойников Юрий Евгеньевич

доктор технических наук, профессор Грузман Игорь Семёнович

доктор технических наук Резник Александр Львович

Ведущая организация Всероссийский научно-исследовательский

институт оптико-физических измерений

Защита состоится "21 " 1АУО уУ1\ 2011 г. в 10 часов на заседании диссертационного совета Д 003.005.01

в Учреждении Российской академии наук Институте автоматики и электрометрии СО РАН по адресу 630090, г. Новосибирск, проспект Академика Коптюга, 1.

С диссертацией можно ознакомиться в читальном зале библиотеки ИАЭ СО РАН. Автореферат разослан "20" ,ХК<Х\ 2011 г.

Учёный секретарь диссертационного совета д.ф.-м.н.

Насыров К. А.

Общая характеристика работы.

Актуальность темы.

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

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

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

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

Тем не менее, в ряде случаев это уравнение является слишком грубым приближением. Оно предполагает, что измеряемая интегральная величина полностью определяется стационарным пространственным распределением некоторого скалярного параметра. Обычно им является линейный коэффициент ослабления зондирующего излучения, либо коэффициент эмиссии. Иногда этого оказывается недостаточно, чтобы описать реальную ситуацию. Большой интерес вызывает эмиссионная томография среды, где также присутствует поглощение. В первую очередь это связано с развитием методов медицинской диагностики SPECT (Single Photon Emission Computerized Tomography) и PET (Positron Emission Tomography). Проблемы, связанные с наличием рассеянного излучения, также не могут быть решены в рамках этой модели сбора проекционных данных.

Использование рассматриваемого уравнения томографии подразумевает определённые свойства регистрирующего оборудования. Вообще говоря, они являются идеализированными и на практике могут быть реализованы лишь приближённо. Любой детектор имеет конечные геометрические размеры и угловую апертуру, поэтому принимает излучение из некоторого объёма пространства, а не с бесконечно тонкого луча. Это обстоятельство следует учитывать при разработке методов повышения пространственного разрешения томограмм.

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

Алгоритмы малоракурсной томографии, разработанные при участии автора, были использованы в научных исследованиях: для получения распределения эмиссии микропинчевого разряда (Физический институт имени П.Н. Лебедева, г. Москва); при изучении методом электронно-пучкового зондирования ударной волны, образующейся при гиперзвуковом обтекании твёрдых тел (Институт теоретической и прикладной механики СО РАН, г. Новосибирск); для реконструкции внутренней структуры клеток крови по результатам интерферометрических измерений (Всероссийский научно-исследовательский институт оптико-физических измерений, г. Москва).

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

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

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

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

Созданные в процессе проведённых исследований алгоритмы объединены в многоцелевой томографический комплекс программ с соответствующим интерфейсом, включающим графический пакет.

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

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

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

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

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

Апробация работы.

Основные результаты диссертационной работы докладывались на 9-й Международной конференции "Методы аэрофизических исследований" (Новосибирск, 1998), на Международной конференции "Обратные и некорректно поставленные задачи" (Москва, 1998), на 12-й Международной конференции "Photonics West'98" (Сан Жозе, Калифорния, США, 1998), на 1-м, 2-м и 5-м Международных конгрессах по индустриальной томографии (Букстон, Англия, 1999; Ганновер, Германия, 2001; Берген, Норвегия 2007), на 4-м Сибирском конгрессе по прикладной и индустриальной математике (Новосибирск, 2000), на Ежегодной конференции германского научного общества по проблемам неразрушающего контроля (Инсбрук, Германия, 2000), па 15-й и 16-й Международных конференциях по неразрушающему контролю (Рим, Италия, 2000; Монреаль, Канада 2004), на 6-й Международной научно-технической конференции "Оптические методы исследования потоков" (Москва, 2001), на Международной конференции по вычислительной математике (Новосибирск, 2002), на Международной конференции "Некорректно поставленные и обратные задачи", посвященной 70-ти легию со дня рождения академика М. М. Лаврентьева (Новосибирск, 2002), на 1-й и 2-й Международных конференциях "Автоматизация, управление и информационные технологии" (Новосибирск, 2002, 2005), на 2-м Международном симпозиуме по неразрушающему контролю "NDN in Progress" (Прага, Чехия, 2003), на Международном семинаре "Неразрушающий контроль в гражданском строительстве" (Берлин, Германия, 2003), на Научно-технической конференции "Томография" (Москва, 2005), на Международной конференции "Обратные и некорректные задачи математической физики", посвященной 75-ти летию со дня рождения академика М. М. Лаврентьева (Новосибирск, 2007).

Публикации. По теме диссертации в период с 1998 по 2009 гг. было опубликовано 46 работ, 25 из них в рецензируемых отечественных и зарубежных журналах, в том числе 16 работ в журналах, рекомендованных ВАК. Личный вклад автора.

Автору принадлежат: идея, разработка и исследование алгоритмов двойной фильтрации и повышения контрастности, модифицированной процедуры нелинейной очистки, метода разложения проекций; обоснование алгоритма пополнения проекционных данных; развитие метода разложения в ряд Неймана; получение всех оценок быстродействия и точности, приведённых в работе; разработка моделей сбора проекционных данных в рассеивающих средах и для трансмиссионной томографии при конечной ширине пучка; научная постановка и решение рассмотренных в диссертации

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

Объём и структура работы.

Диссертация состоит из введения, 5 глав, заключения и списка литературы. Материал изложен на 265 страницах, включая 76 рисунков, 3 таблицы и список литературы из 199 наименований.

Основное содержание.

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

Глава 1 носит обзорный характер.

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

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

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

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

Глава 2 посвящена методам малоракурсной томографии.

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

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

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

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

Для анализа ошибок реконструкции был разработан следующий подход. Ошибка в каждой точке восстановленной томограммы, вызванная заменой гатр-фильтра его аппроксимацией, рассматривается как результат обратного проецирования случайной функции от угла сканирования. Средняя ошибка реконструкции оценивается как корень из дисперсии этой случайной величины. Исходя из такого подхода, было объяснено наличие минимума на зависимости ошибки реконструкции от размера носителя аппроксимации гатр-фильтра, которая получается путём регуляризации обобщённой функции 1/г2 [4].

На практике реконструкция производится по дискретным проекциям, фурье-образы которых являются периодическими функциями. Если спектр проекции широкий, то происходит перекрытие соседних периодов. Чтобы уменьшить связанную с этим ошибку при реализации гатр-фильтра часто используется подавление спектра вблизи частоты Найквиста [3]. В диссертации были разработаны два таких фильтра

(ш)={НетрНс°г). и**» i^JH/ИИ"). Нй(а"

о, И>ю* °>

0)

В формулах (1) п - натуральное число; <uv = n/h, где h - шаг дискретизации проекции.

Было показано, что любой из фильтров (1) более точно аппроксимирует гатр-фильтр, по сравнению с фильтром Рамачандрана-Лакшминараянана [5], если выбирать величину параметра а согласно условию

п+Ъ

0<а<2--Я (2)

1 11ГВХ

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

число) является оценкой снизу среднего значения фурье-образа проекции на интервале [-Зюк,-й)к[о]юл,Зй)„]. При выводе оценки (2) функций (1) раскладывались по степеням ее, при этом отбрасывались члены степени выше первой. Предполагалось, что на значения фурье-образа в пределах ценгрального периода влияют только два соседних периода.

В вычислительном эксперименте было получено, что когда проекций меньше 50 степень модуля частоты п = 2 является оптимальной для обоих фильтров (1). Если же их число превосходит 50, то лучшее качество реконструкции обеспечивается при п > 2. Проведено сравнение фильтрации (1) с другими реализациями гатр-фильтра: Рамачандрана-Лакшминараянана [5], Шеппа-Логана [б], Ерохина-Шнайдерова [7], регуляризация обобщённой функции 1/г2 [4]. Получено, что когда проекций менее 30, для гладких фантомов, имеющих непрерывные производные вплоть до второго порядка,

фильтры , „ (со) дают ошибку реконструкции примерно в 1.1 - 1.3 раза меньше,

чем остальные рассмотренные реализации гатр-фильтра.

Согласно [3], т -мерное преобразование Радона имеет общую формулу обращения:

г = О)

Здесь / - совокупность проекционных данных; - оператор обратного проецирования; т - размерность пространства, на котором задана функция I - потенциал Рисса. По определению, действие оператора Г' на некоторую функцию эквивалентно умножению её фурье-образа на модуль частоты в степени у. Вещественный параметр /3 в (3) может принимать любое значение из интервала ]-ш,т[.

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

В диссертации был предложен и исследован алгоритм двумерной томографии, реализующий формулу (3) при р е ]- 2; 2[. Реконструкция этим алгоритмом производится следующим образом. Сначала фильтруются одномерные проекции посредством фильтра |(о|' затем производится обратное проецирование. Полученное в результате двумерное

ах + а»! . Из-за

дискретизации данных, использование разных значений Р приводит, вообще говоря, к разным результатам. В частности, в диссертации было показано, что ошибки реконструкции для алгоритмов с/?=0и/?=1 имеют различные оценки. При получении этих оценок использовался описанный выше подход, представляющий ошибку реконструкции как корень из дисперсии случайной функции угловой переменной.

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

реконструкции от параметра р. На рис.1 приведены зависимости нормированной среднеквадратичной ошибки Л от параметра Р при различном числе проекций М. Кривые 1 - 5 на рис.1 соответствуют А/=30, Л/=40, М=5О, М=70 и М=100. Как видно из рис.1, при М< 50, зависимости Д(Д) имеют два минимума: один при р < 0, другой при р > 0 (р = 0 соответствует алгоритму обратного проецирования с фильтрацией).

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

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

алгебраических уравнений ^_|a¡jgj = fi. Элемент проецирующей матрицы аи обычно

у

определяется как длина пересечения луча, связанного с / -м измерением, с } -м элементом

1

2

4

5

Рис.1. Зависимости ошибки Д от Р при различном числе проекций

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

в) г)

Рис.2. Математический фантом (а) и томограммы, полученные по 25 проекциям: о) фильтрация и обратное проецирование; е) ART; г) предлагаемый метод

Итерационные алгебраические алгоритмы, применяемые в томографии для решения системы уравнений, полученной после дискретизации задачи, как правило, строятся на I основе поиска экстремума какого-либо целевого функционала [3]. В частности, широко распространённый алгоритм ART (Algebraic Reconstruction Technique) минимизирует функционал невязки [3]. Из множества решений, на которых целевой функционал достигает экстремума, часто выбирается то, которое имеет минимальную норму. В связи с этим изображения, полученные такими алгоритмами, могут оказаться переглаженными. 1 Вычислительный эксперимент показал, что при малом числе ракурсов наблюдения границы | объектов на томограммах, восстановленных алгебраическими методами, оказываются размытыми. Алгоритмы обратного проецирования с фильтрацией в тех же условиях дают ! достаточно контрастное изображение, однако на нём сильно проявляются артефакты, вызванные недостатком проекционных данных, см. пример на рис.2.

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

следующем. Выбирается начальное приближение и производится итерация алгебраического алгоритма, при этом получается решение g"l. Функции и ¿в комбинируются по определённым правилам и результат поступает на следующую итерацию, после которой получается решение g'2'. Оно таюке комбинируется с функцией и далее описанные действия повторяются.

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

где выполняется условие g^>(x,y)-gí"> - среднее значение функции gi■"\

- её средняя величина в окрестности точки (л:,)'), £ - некоторая пороговая константа), значение (х,у) остаётся неизменным. В точках, для которых имеет место обратное неравенство, значение (х,у) заменяется на gFB(x.y).

Путём вычислительного эксперимента было показано, что когда число ракурсов составляет 15 - 50, применение предлагаемого метода приводит к томограмме более высокого качества, чем то, которое обеспечивают алгебраический алгоритм и алгоритм обратного проецирования с фильтрацией, использованные для его построения. Оказалось,

что лучшие результаты получаются, когда отношение составляет от 0.05 до 0.15. С

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

В разделе 2.3 разработана и исследована модифицированная процедура нелинейной пороговой очистки томограмм. Её идея состоит в следующем [1]. Пусть в некоторой точке значение проекции равно нулю. Тогда неотрицательная функция должна быть равна нулю всюду на луче, вдоль которого производится интегрирование. Зануляя восстановленную томограмму вдоль всех таких лучей, можно существенно повысить качество реконструкции благодаря тому, что будут уничтожены внешние артефакты, которые, в частности, возникают из-за малого числа ракурсов наблюдения.

Поскольку в реальной ситуации из-за наличия случайного шума значение проекции нигде не равно нулю, вводится пороговая величина е > 0. Томограмма в /-м узле зануляется, если хотя бы одна из проекций в окрестности точки Вр, меньше £ (Вр, -проекция у'-го узла на т-ю плоскость регистрации). После такой процедуры томограмма может оказаться занулённой в области, которую занимает объект, в частности вблизи его границ. Предлагаемая в диссертации модификация заключается в том, чтобы занулять томограмму в у'-м узле только тогда, когда условие /„(в^)<£ выполняется дня

нескольких проекций, число которых обозначим через Мс.

Пусть Р(/т(Вр,) > а) вероятность того, что значение т-й проекции в точке В^ больше величины а. Предположим, что для интересующих нас узлов вероятности Р(/т(Вр,) > е), Р(/„,(В]т) >£)(<!; -- корень из дисперсии шума) одинаковы для всех проекций и обозначим их через Р(/(В]) > £), Р(/[В]) > £) соответственно. Тогда вероятность того, что томограмма не занулится ву'-м узле после процедуры нелинейной очистки с параметрами £, Мс, будет

м\

'с )(М-МС)!МС!

(4)

Р(е,Мс,М,Л = -(1 -Р(/(В1)>е))1

Здесь М - полное число проекций. Будем считать, что вероятность того, что искомая функция в ] -м узле не равна нулю есть > . После процедуры нелинейной

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

Р„ {е,Мс,М,]) = Р(е,МсМ,]) + Р(/{Б^> <?)-

-2 р{е,мс,м,])р(/(в,)>£). (5)

Если предположить, что вероятность Р^/^В^ > ^ стремится к единице при »0 и уменьшается с ростом Е,, то, исходя из уравнения (5), можно заключить, что при низком уровне шума лучшего результата очистки можно достичь, выбирая Мс > 1.

Вычислительный эксперимент показал, что процедура нелинейной пороговой очистки эффективнее простого зануления отрицательных значений на томограмме. Было получено, что, когда число проекций больше, чем 30, при одинаковой интенсивности шума в данных нелинейная очистка с Мс > 1 обеспечивает меньшую среднеквадратичную ошибку Л, чем с Мс = 1. При любом е €[0.25|; 0.75|] зависимости А {Мс) достигают минимума при Мс> 1.

Причём большим значениям е соответствуют большие Мс.

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

Пусть скорость исследуемого стационарного потока имеет единственную компоненту 1>г=о(х, у, г). В момент времени 1 = 0« плоскости г = 0 молекулы переводятся в возбуждённое состояние. Будем полагать, что уменьшение числа возбуждённых молекул происходит главным образом за счёт их перехода в основное состояние путём спонтанного излучения. При этих условиях было получено следующее соотношение:

(б)

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

за время измерения через окрестность точки (х,у,г), ДК - объём этой окрестности.

Функция п[х,у,г) может быть найдена путём томографической реконструкции по интенсивности флюоресценции, зарегистрированной двумерными детекторами, расположенными параллельно направлению потока под различными азимутальными углами. Если п(х, у, ¿) известна, скорость в каждой точке определяется из уравнения (6).

Разработанный метод диагностики потоков исследовался посредством вычислительного эксперимента. Было получено, что при увеличении параметра /3 (т.е. при уменьшении среднего времени пребывания молекул в возбуждённом состоянии) точность восстановление поля скорости падает.

В главе 3 рассматриваются задачи томографии с ограниченным доступом к исследуемому объекту.

Раздел 3.1 посвящен томографии с ограниченным углом обзора объекта. Для двумерной томографии с ограниченным диапазоном углов в диссертации был разработан и исследован алгорипъч генерации проекций. Пусть проекции измерены в интервале углов [о,<р0]. По ним реконструируется начальное приближение, в него вносится априорная информация, известная об искомом распределении. От полученной функции вычисляются проекции в тех направлениях, где они отсутствуют. Затем производится реконструкция по набору, состоящему из известных и вычисленных по начальному приближению данных. В решение вновь вносится априорная информация, и далее описанные действия повторяются. В результате приходим к итерационному процессу:

(7)

Здесь через "Г1 обозначен оператор, соответствующий алгоритму решения задачи двумерной томографии; через Ф<п) - операторы априорной информации для каждой итерации; через Р^ 2Л[ - оператор, вычисляющий проекции от функции в интервале

углов ]<р0,2я[; через/- имеющиеся данные. Выражение / + '' означает полный

набор данных, включающий существующие и вычисленные проекции.

В качестве априорной информации рассматривались неотрицательность искомой функции и её принадлежность определённым классам гладкости. Поэтому в операторы ф'"' входили процедуры сглаживания и зануления отрицательных значений.

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

9т={р-\)\, = + (8)

м V. .'^0 )

Здесь предполагается, что имеется М проекций, равномерно распределённых в угловом диапазоне [0,р0] с шагом . Каждая из них определена на интервале [— г0, г0 ] в Л" точках с

шагом дискретизации А.

Для случая, когда оператор Т"1 соответствует алгоритму обратного проецирования с фильтрацией, была выведена оценка средней суммарной ошибки восстановления каждой проекции по начальному приближению :

--Х, 1 -у- (/,. со* <р„) + /, МП ( I, соч,<рт)

мх

8(р 41 гп' ' гт 81

//

9т = %+тИ9, I, = -г0 + Ш. (9)

В (9) чертой сверху обозначено усреднение по угловой переменной в диапазоне

/и' - проекционные данные после гатр-фильтрации; М\ - количество вычисленных

проекций, которые равномерно распределены в угловом интервале ]<р0,2тг[ с шагом И^.

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

Можно считать, что в предлагаемом алгоритме приближения ¿г''"1 и получаются при реконструкции по набору данных, в котором проекции в угловом интервале 2/г[ известны с ошибкой Д^0' и Д^' соответственно. Естественно предположить, что если будет выполняться условие то ошибка решения полученного после первой

итерации, будет меньше ошибки начального приближения .

В рассмотренном обосновании алгоритма генерации проекций используются оценки Др' и Д«, которые могут быть получены из имеющихся проекционных данных. Обычно оценка точности реконструкции по данным, заключённым в ограниченном диапазоне углов, производится, исходя из других принципов, основанных на анализе фурье-образа искомой функции. При этом для её получения требуется нетривиальная априорная информация.

Было предложено прекращать итерационный процесс (7), когда ошибка восстановления недостающих проекций будет больше величины самих проекций. Если предположить, что проекции во всех направлениях восстанавливаются с приблизительно одинаковой интегральной ошибкой, то для проверки этого условия можно брать известные данные из углового интервала [0,<ро]. В результате приходим к следующему критерию останова:

—У мЬ

^ М ! N \

Ч>м=(т-Щ> + (10)

Посредством вычислительного эксперимента было показано, что для данных, известных в ограниченном диапазоне углов, алгоритм генерации проекций позволяет повысить качество реконструкции по сравнению с алгоритмами обратного проецирования с фильтрацией в случае, когда выполняется условие А^' < . В частности при <р„е[90°; 180°] он обеспечивает среднеквадратичную ошибку в 2 - 3 раза меньшую, чем алгоритм Шеппа-Логана. Пример представлен на рис.3.

Во всех проведенных расчётах зависимости среднеквадратичной ошибки реконструкции от числа итераций п имели минимум, положение которого зависит от количества проекций М, величины угла <Ро и уровня шума. При увеличении М, минимум сдвигается в сторону больших п. Так для М= 100 минимальная ошибка достигается на 8 -12-й итерации, а для М=500 на 50 - 60-й. Минимум перемещается к началу кривой при уменьшении интервала углов, в котором известны проекционные данные, и при повышении уровня шума. Номера итераций, полученные по критерию (10), всегда оказывались несколько больше, чем те, на которых достигался минимум среднеквадратичной ошибки.

Однако отклонение между этими номерами не превышало 15%.

'''чаимЗСТОЗКА'''

а)

б)

Рис.3. Результаты реконструкции математического фантома, 500 проекций в угловом интервале 90°: а) алгоритм генерации проекций, Л=0.278; б) алгоритм Шеппа-Логана, А=0.742

«)

Рис.4. Математический фантом («) и результаты реконструкции, 500 проекций: б) алгоритм обратного проецироваиия с фильтрацией; в) метод генерации проекций

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

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

Алгоритмы тестировались путём вычислительного эксперимента и реконструкции по реальным данным. На рис.4,а изображена часть математического фантома, отражающего основные черты продольного сечения полой трубы с метками на её внешней и внутренней поверхности. На рис.4,б и рис.4,в приведены результаты его реконструкции алгоритмом обратного проецирования с фильтрацией и методом генерации проекций соответственно. Угол обзора центральной части фантома составлял примерно 70°.

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

Путём численного моделирования было проведено методическое исследование возникающей при этом задачи восстановления функции трёх переменных по набору из нескольких десятков двумерных проекций, содержащихся в ограниченном диапазоне по азимутальному и полярному углам. Реконструкция проводилась алгебраическими алгоритмами ART [3], MART (Multiplicative Algebraic Reconstruction Technique) [9], MENT (Maximum entropy) [10]. Кроме того расчёты велись методом Гершберга-Папулиса [1], а также комбинированным алгоритмом cART, полученным в результате объединения ART с методом Гершберга-Папулиса [11]. Для экваториальной схемы сбора данных (т.е. когда все плоскости регистрации параллельны оси Z ) получено, что, если проекции заключены в узком интервале углов, составляющем 20° - 50°, cART и MENT обеспечивают ошибку реконструкции в 1.5 - 2 раза меньшую, чем остальные алгоритмы.

В условиях ограниченного угла сканирования проведено исследование задачи обнаружения включений на постоянном фоне в зависимости от их размера и контрастности. При этом использовался следующий критерий качества реконструкции: c/ = |gi'-gjj(/||gJ/-gJi|j, который чувствителен к малым изменениям размера или амплитуды включения. Здесь gM и g - математический фантом и результат реконструкции соответственно; gM - среднее значение фантома. Усреднение ведётся по той области, где

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

а) б)

Рис.5. Томограммы клеток крови человека; а) эритроцит; б) лимфоцит

Результаты описанного выше исследования послужили основой для реконструкции клеток крови - эритроцитов и лимфоцитов - по данным, зарегистрированным на интерференционном микротомографе во Всероссийском научно-исследовательском институте оптико-физических измерений (г. Москва). Из снятых ингерферограмм были получены экваториальные проекции. Для эритроцита их было 18, а для лимфоцита - 17, при этом диапазон по углу составлял [-18°; 16°] и [-24°; 40°]. На рис.5,а даны две изоповерхности распределения показателя преломления внутри эритроцита, восстановленного алгоритмом сАКТ. На рис.5,б приведён лимфоцит в разрезе. Расчет произведён методом Гершберга-Папулиса.

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

В работе [12] были сформулированы условия, которые должны накладываться на траекторию источника, чтобы в пределе бесконечного числа измерений была возможна точная томографическая реконструкция. Такие траектории будуг в дальнейшем называться полными. Основным из условий полноты является следующее: любая плоскость, пересекающая носитель реконструируемой функции, должна также пересекать траекторию источника хотя бы в одной точке. При разработке томографического оборудования по различным причинам не всегда удаётся реализовать полную траекторию источника. В связи с этим в диссертации исследуются различные массы траекторий.

В силу своей простоты часто используется следующая схема регистрации проекционных данных. Источник движется в плоскости г = 0 по окружности, охватывающей объект. Двумерный детектор, регистрирующая поверхность которого перпендикулярна плоскости г = 0, движется синхронно с ним с противоположной стороны объекта, см. рис.6,а. Эта траектория не является полной, поэтому для реконструкции по такому набору данных используются приближённые методы. В работе сравниваются два из них: алгоритм Фельдкампа [13] и послойное восстановление, которое осуществлялось алгоритмом фурье-синтеза [3].

Обычно в послойных алгоритмах для реконструкции слоя с координатой г = г0 берутся проекционные данные, имеющие на детекторе ту же координату. Проведённое исследование показало, что выбор данных с координатой г = (г, + ъ)г0/»; обеспечивает меньшую ошибку. Здесь гг и гл -расстояния от центра траектории до источника и детектора соответственно.

Рис.6. Схемы сбора проекционных данных для трехмерной томографической реконструкции: и) круговая траектория источника; б) спиральная траектория источника;

1 — источник, 2 — исследуемый объект, 3 - двумерный детектор

Были сделаны оценки времени реконструкции рассматриваемыми методами. Пусть имеется Мпроекций, каждая из которых задана в ЫхИ узлах, искомая функция вычисляется на сетке NxNxN узлов, N равно целой степени числа 2, Если в методе фурье-синтеза используется алгоритм быстрого преобразования Фурье, то при больших N и М он оказывается быстрее алгоритма Фельдкампа в !"*0(А?71о£2 И) раз. Эту оценку

подтверждает таблица 1, где приведены данные для N = М.

Таблица 1. Время реконструкции методом послойного Фурье-синтеза и алгоритмом Фельдкампа

N Послойный Фурье-синтез (с) Алгоритм Фельдкампа (с) N1 Хо^Ы

128 3 22 7.33 18.29

256 32 840 26.25 32

512 311 16272 52.32 56.89

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

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

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

Для полных траекторий источника имеет место формула обращения, которая была выведена в [12]. С вычислительной точки зрения она является неэффективной, поскольку для реконструкции функции в одной точки требуется пятикратное интегрирование. В работе [4] формула обращения была приведена к виду, более удобному для построения на её основе численных алгоритмов:

1 2 л л 1

J LJ(ß,Wfi

Щ

sinddödtp (jj)

Ц(с'(л),я)дЯ

Здесь С(А) = (Я),Су(Я),С1(Л)) - уравнение траектории источника; С'(Л) - производная

от С(А) по А. Функция получается из проекционных данных, путём деления на

длину вектора ß, ßeR\ Через S(n) обозначена окружность, по которой пересекаются единичная сфера и плоскость, проходящая через начало координат ортогонально вектору й , имеющему координаты (sin б cos (р, sin б sin (р, COS0); dß - мера на S{n); -

оператор дифференцирования в направлении вектора п ; - скалярное произведение. Значение параметра А при фиксированных г и п определяется из уравнения (r,n)-(c(X),nj. Существование его решения, такого что (с'{к),п^ О, гарантируется

выполнением условия полноты траектории.

Формула (11) для реконструкции функции в одной точке требует всего три интегрирования, что является существенным преимуществом по сравнению с формулой из работы [12]. Тем не менее, она до сих пор не нашла практического применения. Для спиральной траектории источника разработаны специализированные алгоритмы, которые обладают высокой вычислительной эффективностью, однако не могут применяться для траекторий другой формы.

В диссертации на основании формулы (11) был разработан алгоритм реконструкции, который может применяться как для полных, гак и для неполных траекторий. Из трёх рассматриваемых траекторий — спирали, двух взаимно перпендикулярных окружностей и трёх параллельных окружностей — первые две являются полными, последняя - неполной. Кроме реализации формулы (11), в вычислительном эксперименте использовался алгебраический алгоритм ART. Для реконструкции по данным, полученным при спиральном движении источника, также использовались два алгоритма, разработанные для такой формы траектории. Оба строятся по следующей схеме: производится фильтрация двумерных проекций, затем интегрирование отфильтрованных проекционных данных по траектории источника. Отличие состоит в том, что в одном используется одномерная фильтрация [14], в другом - двумерная [15], которая обеспечивает более высокую точность.

Если функцию g(?) вычислять по формуле (11) отдельно в каждом узле, то количество операций для её реконструкции в NxNxN узлах по М проекциям, заданным на сетке JVxJV, оценивается следующим образом: Nl MN + NN ¿К . Здесь через N й

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

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

оценивается как ЛГ2/В'МЫ2 (ЛЙТ + Л^). При больших./Vимеет место И^/Ы^р-К.

Таблица 2. Время реконструкции для разных алгоритмов и разных траекторий

"~-\Алгоритм Траектория Спиральный ID фильтр (с) Спиральный 2D фильтр (с) Реализация формулы (11) (с) ART (с)

Спираль 832 1250 32425 3492

Перпендикулярные окружности - - 32944 3518

Параллельные окружности - - 32837 3540

Сравнение быстродействия алгоритмов реконструкции, ориентированных на спиральную томографию, с алгоритмом, реализующим (11), который может использоваться для любой траектории источника, даёт отношение числа операций порядка М/Ыя .

В вычислительном эксперименте на области реконструкции была задана сетка 257x257x257 узлов. Двумерные проекции определялись на сетке 257*257 >"ЗЛов, их число варьировалось от 50 до 500. При реконструкции алгоритмом, основанным на формуле обращения (11), бралось N в = 120x60 направлений вектора п (через 3° по азимутальному и полярному углу), для вычисления внутреннего интеграла в (11) использовались 257 значений подынтегральной функции.

Некоторые результаты вычислительного эксперимента представлены в таблттах 2, 3. Количество двумерных проекций равнялось 240. В таблице 2 приведены времена реконструкции, а в таблице 3 - нормированные среднеквадратичные ошибки. Количество итераций при реконструкции алгоритмом ART равнялось пяти (в таблице 2 указано суммарное время для всех итераций).

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

Алгоритм Траектория —. Спиральный ID фильтр (с) Спиральный 2D фильтр (с) Реализация формулы (И) (с) ART (с)

Спираль 0.154 0.082 0.075 0.119

Перпендикулярные окружности - - 0.077 0.126

Параллельные окружности - - 0.228 0.223

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

был получен как для реализации формулы (11), так и для ART. Таким образом, свойство полноты траектории играет определяющую роль для качества реконструкции. Её форма оказывает значительно меньшее влияние. В частности, это видно из таблицы 3.

Когда число проекций составляет несколько сотен, алгоритм, реализующий общую формулу (11), при реконструкции по спиральным данным на 5% - 10% превосходит по точности специализированный алгоритм, использующий двумерную фильтрацию проекций. Причём точность алгоритма, основанного на (11), может быть повышена путём увеличения числа направлений вектора п .

Время счёта для алгоритма, реализующего формулу (11), оказалось примерно в 26 и в 39 раз больше времени счёта для спиральных алгоритмов с двумерной и одномерной фильтрацией соответственно. При этом величина N^/M равнялась 30. Отношение времён реконструкции спиральными алгоритмами с двумерной и одномерной фильтрацией составило 1.5. Время, затрачиваемое на одну итерацию ART, приблизительно в 1.2 раза меньше того, которое нужно на реконструкцию спиральным алгоритмом с одномерной фильтрацией. Однако для получения хорошего восстановления посредством ART надо сделать от нескольких единиц до нескольких десятков итераций в зависимости от числа проекций.

В разделе 3.3 промоделирована задача малоракурсной томографии, ориентированная на экспериментальное изучение обтекания твёрдых тел потоками жидкости или газа. Рассматривается следующая математическая постановка. Функция gixy^) реконструируется внутри ограниченной области ilcR3 (supp(g)cQ) по набору своих двумерных проекций. При этом для любого луча, пересекающего область П, с П, значение проекции в соответствующей точке приравнивается к нулю. Множество всех таких точек называется тенью, а область ii, - непрозрачным телом.

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

Для решения рассматриваемой задачи были использованы алгебраические методы. При построении системы уравнений не учитывались лучи, пересекающие непрозрачное тело, так как значения проекций в точках, попадающих в тень, обращаются в нуль, в то время как интегралы or g(x,y,z), взятые вдоль соответствующих прямых, вообще говоря, нулю не равны. В результате посредством алгоритма ART по шести проекциям, искажённым тенью непрозрачного тела, качественно верно был восстановлен сложный объект, моделирующий картину обтекания вращающегося конуса струями горячего газа.

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

В разделе 4.1 для решения задачи реконструкции по данным, искажённым аддитивным фоном, был разработан метод разложения проекций.

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

Рассматривается двумерная задача. Пусть неискажённые проекционные данные /о связаны с искомой функцией g посредством двумерного преобразования Радона. Будем называть /а истинными проекционными данными. Предположим, что в результате измерений известны данные / = /о+Х, где % ~ аддитивный фон. Пусть Т"1 - оператор, соответствующий какому-либо алгоритму обращения преобразования Радона. В результате его применения к истинным данным /о будет получена функция ge, которая, вообще говоря, отличается от функции g из-за того, что набор данных, по которому производится реконструкция, конечен. Обозначим через Р проецирующий оператор, переводящий функцию в набор неискажённых проекционных данных. Имеет место

Рге = РТ-1/о = /о+<5. (и)

где 8- некоторая функция, определённая на том же множестве, что и /а. Выражение (12) показывает, насколько хорошо решение £ согласуется с истинными проекционными данными. Как правило, 8ф 0, но её можно считать достаточно малой. Если положить в (12) 5=0, то /о инвариантная функция оператора РТ"1. Если РТ"1 линейный оператор в гильбертовом пространстве Я, то при дополнительных условиях, которые обычно выполняются на практике, имеет место статистическая эргодическая теорема [16]:

л—ню и 4 ' '

___, (13)

где /" - функция, инвариантная относительно оператора РТ"1, / - любая функция, принадлежащая гильбертову пространствуя.

Возьмём в (13) в качестве / имеющиеся искажённые проекционные данные/=/о+Х-Инвариантную функцию оператора РТ"1, к которой при этом будет сходиться (13), будем считать оценкой истинных проекционных данных /о. На основе формулы (13) в диссертации предложен следующий алгоритм поиска /о:

йп) +Я(")(РГ1Г/))), Г'-0, „ = 1,2,.... (14,

В (14) А(п)

— числовой параметр, Фр 1 и ^ — операторы априорной информации для проекционных данных и фона соответственно, которые могут изменяться в зависимости от номера итерации. После того, как оценка истинных проекционных данных будет определена в итерационном процессе (14), по ней производится томографическая реконструкция искомой функции g.

При / -/ц + Х выражение (13) сходится к /" =/0 + Х" • Таким образом, истинные проекционные данные восстанавливаются с точностью до инвариантной функции оператора РТ"1, зависящей от фона. В качестве критерия останова итерационного процесса (14) было принято неравенство Цу,""1' —Л|| -котоР°е приводится к виду

аз)

IIИ /¡=0 II

Для вычисления правой части (15) в диссертации была получена рекуррентная формула.

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

переменный фон Хт (О = (Р , где /Зт и ут - случайные числа

с равномерным распределением в интервале [—1;1].

В проведённом вычислительном эксперименте применение метода разложения проекций позволило уменьшить среднеквадратичную ошибку реконструкции фантома в 2.5 - 5 раз. Зависимости ошибки реконструкции проекций, фона и фантома имели минимумы, приходящиеся на близкие итерации. Помер итерации, полученный по правилу останова (15), оказался также близким к ним. При заданных операторах априорной информации и фиксированных значениях параметров положение минимумов зависит от величины фона, числа проекций и уровня шума. При увеличении амплитуды фона, а также уровня шума, минимумы сдвигаются ближе к началу кривых.

Рис.7. Реконструкция по данным, содержащим фон; а) метод разложения проекций, Д=0.299; б) алгоритм Шенпа-Логана, Д=0.879

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

На рис.7 представлены томограммы. Для реконструкции использовались 500 проекций, искажённых аддитивным фоном с/? = 0.5. Шум в данных отсутствовал.

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

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

/= jg(/?)exp -|rj(p')dp'

dp

(16)

-со у р

где ^(р), Г] (г) - распределения коэффициентов эмиссии и абсорбции соответственно.

Существующие в настоящее время методы реконструкции функции g(f) по данным

(16) требуют знания распределения коэффициента абсорбции Т]{г). Для двумерного случая известны формулы, выражающие через данные (16), зарегистрированные в

параллельной и веерной схемах сканирования [17, 18]. В трёхмерном случае используются алгебраические методы. При этом элементы проецирующей матрицы строятся, исходя из соотношения (16).

В диссертации для решения рассматриваемой задачи был развит метод, предложенный ранее в [19], который основывается на разложении в ряд Неймана. Неизвестное распределение эмиссии аппроксимируется следующей суммой:

^'^-'[Т^Я^Е^Т'ГЛ]. (17)

Здесь Т"1 - оператор решения задачи эмиссионной томографии без поглощения, т.е. "Г"1 реконструирует функцию по набору её интегралов вдоль прямых. Через Р,, обозначен оператор, вычисляющий проекционные данные по уравнению (16); Ф1"' - оператор априорной информации для распределения эмиссии; Е - единичный оператор; А -положительный вещественный параметр; / - совокупность проекций с поглощением. Поскольку каждый из членов суммы (17), взятый в отдельности, не отвечает какому-либо реальному распределению параметров, то оператор Ф*-"' применяется ко всей сумме каждый раз после добавления очередного члена.

При л = оо в правой части (17) стоит сумма ряда Неймана оператора РцТ"1, которая

равна в области сходимости обратному оператору. Следовательно, = (Р„Т~1) ' /ч, и в

случае существования обратного оператора Р получаем Я1*'=РЧ_1/,-

Ряд Неймана сходится внутри единичного шара в гильбертовом пространстве при условии | А | |Е — РПГ1)| < 1. Исходя из этого, было предложен следующий критерий останова алгоритма (17):

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

Разработанные методы исследовались путём вычислительного эксперимента. Было проведено их сравнение с алгебраическими алгоритмами реконструкции ART и MART, при помощи которых решалась система уравнений, построенных с учётом поглощения. Моделировалась трёхмерная задача.

Было установлено, что зависимость среднеквадратичной ошибки от числа членов суммы (17) имеет минимум. Число, при котором этот минимум достигается, и число, полученное по критерию (18), отличаются не более чем на 10%. В проведённых

вычислениях эти числа принимали значения 15 - 30 в зависимости от условий эксперимента.

При небольшом поглощении алгебраические алгоритмы дают лучший результат. Однако с увеличением поглощения преимущество переходит к разложению в ряд Неймана. При этом для всех методов точность реконструкции падает при увеличении оптической толщины объекта. Зависимость качества реконструкции от количества проекций для разложения в ряд Неймана более ярко выражена, чем для алгебраических алгоритмов. Кроме того, оказалось, что при большом поглощении ошибка реконструкции алгоритмом ART больше, чем алгоритмом MART.

Расчёты по проекционным данным, содержащим случайный шум, показали, что при включении в итерационный процесс сглаживающих процедур все алгоритмы достаточно

устойчивы. Наилучшие результаты были получены для регуляризующего фильтра 1/а)®)2 с вычислением параметра регуляризации а на каждой итерации по критерию невязки.

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

На рис.8 сравниваются различные методы реконструкции. Здесь приведены зависимости среднеквадратичной ошибки Д от параметра щ, (усредненная оптическая толщина вдоль координатных осей). Число ракурсов 25. Кривая 1 -реконструкции в предположении, что поглощение отсутствует (алгоритм ART); 2 -алгоритм ART с учётом поглощения; 3 - разложение в ряд Неймана; 4 - метод разложения проекций. Видно, что при q0< 5 метод разложения проекций уступает по точности другим алгоритмам, учитывающим поглощение. Но для его применения не требуется знать распределение абсорбции.

В разделе 4.3 исследуется томография рассеивающих сред.

Решение некоторых томографических задач осложняется из-за засветки детектора рассеянным излучением. Для улучшения качества реконструкции по таким проекционным данным часто используется метод компенсации. Он состоит в следующем. Производится математическое моделирование процесса засветки при условиях, максимально приближенных к тем, при которых производилась регистрация проекционных данных. Затем вычисленное распределение интенсивности рассеянного излучения по детектору вычитается из зарегистрированного сигнала [20].

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

А

Рис.8. Сравнение методов эмиссионной томографии среды частично поглощающей излучение

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

Чтобы построить оператор Ра, были выведены уравнения эмиссионной и трансмиссионной (когда объект исследуется путём зондирования излучением от внешнего источника) томографии в рассеивающей среде. Эти уравнения имеют различную форму.

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

Г Г „ Л Л

f = J s{p)exp -J^o,{Р')4Р' +К{">Р)

dp.

(19)

V V у у

Здесь 1 - интенсивность излучения, зарегистрированная в точке детектора, через которую проходит прямая, вдоль которой производится интегрирование; g - исходное

распределение эмиссии; оь/ - плотность рассеяния, отнесённая к единице длины. Через обозначено распределение коэффициента эмиссии, индуцированное благодаря рассеянию. Поверхность детектора перпендикулярна вектору п .

Разобьём исследуемую область на кубические вокселы. В пределах каждого воксела будем считать функции ^ и Оо/ постоянными. Значение функции в точке В, принадлежащей рассматриваемой прямой, определяется путём суммирования вкладов от всех вокселов:

( (

А, г.*тя

arceos

V

1 —

п-

Г.-Г.

\r.-r.

(20)

J)

где Vo — объём воксела; гл - радиус-вектор точки А, являющейся центром воксела; а0 (гв) -плотность рассеяния в точке В\ <70 - средняя плотность рассеяния; Oi(/3) - угловое распределение рассеяния (угол /3 отсчитывается от направления п ). В численном моделировании предполагалось, что в основном происходит фотоэлектрическое взаимодействие фотонов с веществом, для которого имеет место <7, (J3) = 1 + cos2 ¡3.

В случае многократного рассеяния для получения полного сигнала на детекторе суммировались вклады, обусловленные рассеянием каждого порядка. В результате первичного рассеяния генерируется эффективная эмиссия £ffi, величина которой в каждой точке вычисляется по формуле (20). Эмиссия eai будет порождать эмиссию и так далее. В разработанной компьютерной программе ограничение на количество учитываемых порядков рассеяния устанавливалось путём введения предельной величины г^. Все значения индуцированной эмиссии меньшие Сщш приравнивались нулю.

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

По результатам вычислительного эксперимента были сделаны следующие выводы. При не очень больших значениях ст0 (приблизительно до 0.5) преимущество алгоритма (17) по сравнению с методами реконструкции, не учитывающими рассеяние, тем существеннее, чем больше количество проекций. С увеличением <Т0 для алгоритма (17) имеет место ухудшение качества восстановления, как и в задаче эмиссионной томографии в поглощающей среде при увеличении оптической толщины. Разложение в ряд Неймана более чувствительно к занижению величины плотности рассеяния, чем к её завышению. Метод разложения проекций позволяет получить результаты лучшие, чем те, которые дают стандартные томографические алгоритмы, не учитывающие рассеяния, но уступает по точности разложению в ряд Неймана.

Рис.9. Математический фантом (а) и результаты реконструкции по 5 проекциям с рассеянием: о) ART без учёта рассеяния, Д = 0.960; к) метод разложения проекций, Д = 0.891

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

В этом разделе также исследовались проблемы, связанные с рассеянным излучением, возникающие в случае зондирования твёрдых тел. Было получено выражение для засветки плоского детектора однократно рассеянными фотонами. Интенсивность излучения, рассеянного с луча, принадлежащего оси X, перпендикулярной плоскости регистрации, в точке детектора («о, Vo) выражается следующим образом:

) = /oJc70/(x,0,0)exp -J g, (>',0,0]ск'- jg(x,у,z)dp

а V а I

(l + cos2(£))cos3(j8)

•d

+ xf

-dx

(21)

где а, Ь - точки пересечения луча с поверхностью объекта; гд - расстояние от детектора до начала координат; g+Oa¡ - распределение полного линейного коэффициента ослабления, включающего как поглощение g, так и рассеяние ст0/, 10 - интенсивности источника зондирующего излучения. Предполагается, что в основном происходит фотоэлектрическое взаимодействие излучения и среды. Интегрирование во втором члене показателя

экспоненты производится вдоль отрезка, соединяющего текущую точку на луче (х,0,0) с точкой, имеющей координату на детекторе (ио,М>); Р угол между этим отрезком и осьюХ

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

Для решения рассматриваемой задачи использовались ряд Неймана и метод разложения проекций. Проецирующий оператор строился, исходя из уравнения (21), т.е. с учётом только однократного рассеяния. Предполагалось, что число фотонов, рассеянных в малом объёме, пропорционально с коэффициентом /числу фотонов, поглощённых в этом же объёме. Было получено, что применение алгоритмов, учитывающих рассеяние, позволяет понизить ошибку реконструкции в 1.1-1.8 раза (разложение проекций) и в 1.9 -2.5 раза (ряд Неймана). Исследование влияния многократно рассеянных фотонов показало, что они начинают играть заметную роль, когда величина /становится больше чем 10 3.

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

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

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

Поскольку в рассматриваемом случае одной из проблем является размытие мелких деталей томограммы, представляется обоснованным применять методы оценки искомой функции, усиливающие градиенты. В диссертации для этой цели проведено обобщение алгоритма Вайнберга [21] на трёхмерный случай. Разработанная модификация метода Вайнберга состоит в следующем. Каждая двумерная проекция, полученная в конических пучках, дифференцируется двумерным оператором Лапласа. После этого производится процедура обратного проецирования для лучевого приближения задачи. Фильтрация посредством оператора Лапласа соответствует умножению фурье-образа проекции на квадрат модуля частоты, в то время как гатр-фнльтрация - это умножение на модуль частоты. Таким образом, при использовании метода Вайнберга получается не томограмма, а контрастное изображение, на котором лучше проявляется мелкомасштабная структура.

В вычислительном эксперименте было проведено сравнение алгебраических алгоритмов и метода разложения в ряд Неймана. Реконструировалась непрерывная функция трёх переменных по 6-ти проекциям в конических пучках. Было получено, что алгебраические алгоритмы ART и MART обеспечивают среднеквадратичную ошибку в 1.6 раза меньшую, чем разложение в ряд Неймана. Однако они оказались менее устойчивыми по отношению к случайному шуму. Это объясняется тем, что в случае широких пучков

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

Вычисление проецирующей матрицы в приближении пучков конечной ширины приводит к резкому увеличению времени счёта. В проведённом вычислительном эксперименте оно оказалось в 300 - 500 раз больше (в зависимости от параметров пучка), чем для лучевого приближения. В связи с этим было проведено исследование возможности использования лучевого приближения для детекторов, имеющих конечные параметры при реконструкции алгебраическими методами. В случае цилиндрических пучков оно обеспечивает неплохую точность реконструкции до величины диаметра пучка 4-7 длин стороны воксела. Для конических пучков уже при угле полураскрыва в 3° оно заметно уступает приближению пучков конечной ширины.

в) г)

Рис. 10. Примеры реконструкций по проекциям в конических пучках

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

В разделе 5.2 рассматривается рентгеновская томография в широких пучках. Модель детектора берётся такой же, как и в разделе 5.1. Считается, что выходное пятно источника имеет круглую форму, а интенсивность излучения распределена равномерно по его площади. В предположении, что для распределения линейного коэффициента ослабления рентгеновского излучения T](x,y,z) имеет место ri(x,y,z) = ri0+g(x,y,z), причём

J g(x,y,z)dl1& 1 (L(y,z) - любой луч, лежащий внутри пучка), показано, что из Чу,')

проекционных данных могут быть получены интегралы от функции g(x,y,z) по объёмам ^ прямых усечённых конусов. Таким образом, при этих условиях задача трансмиссионной ! томографии с широким зондирующим пучком описывается таким же уравнением, что и эмиссионная задача, исследованная в разделе 5.1. Следовательно, для её решения можно применять разработанные там методы. |

Вычислительный эксперимент проводился с фантомом, имеющим особенности, характерные для задач дефектоскопии. Его сечение плоскостью х = 0 представлено на рис. 10,а. На рис. 10,6 - рис. 10,г показаны сечения той же плоскостью реконструированных [ функций. Использовались 25 проекций в конических пучках. Средний диаметр пучка составлял 8 размеров стороны воксела. Разбиение области реконструкции 129x129x129 вокселов. Рисунок 10,(5 соответствует алгоритму ART в лучевом приближении (Д=0.227), рис.10,в - алгоритму ART в приближении конечной ширины пучка (Д=0.129), рис.10,г - 1 разложению в ряд Неймана (Д=0.171).

Среднеквадратичная ошибка реконструкции путём разложения в ряд Неймана оказалась приблизительно в 1.3 раза больше, чем для алгебраических методов. Тем не 1 менее, мелкомасштабные высококонтрастные структуры, размер которых в несколько раз меньше диаметра пучка, восстанавливаются лучше алгоритмом, основанным на разложении в ряд Неймана. В частности, это видно из сравнения рис. 10,в и рис. 10,г. ^

По исследованиям, проведённым в разделах 5.1, 5.2, были сделаны следующие выводы. Несмотря на то, что в вычислительном эксперименте метод разложения в ряд Неймана показал большую среднеквадратичную ошибку, чем алгебраические алгоритмы, при решении рассматриваемой задачи он имеет ряд преимуществ перед последними. Есть 1 основания предполагать, что при большом числе проекций он будет давать томограммы более высокого качества. Уже для 25 проекций мелкомасштабная структура,

восстановленная с его помощью, выглядит предпочтительнее. Кроме того, разложение в ряд Неймана | оказалось более устойчивым по ; отношению к случайным шумам, чем j алгебраические методы. |

В случае, когда представляют интерес только форма и размер I мелкомасштабной структуры, но не её амплитуда, может быть использован метод Вайнберга, основное преимущество которого заключается в его быстродействии. [ Однако метод Вайнберга очень чувствителен к высокочастотному шуму.

В разделе 5.3 исследуются проблемы, связанные с размытием

Рис. 11. Функция импульсного отклика и её аппроксимация гауссовой функцией

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

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

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

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

В процессе вычислительного эксперимента было установлено, что применение процедуры деконволюции приводит к уменьшению ошибки в оценке длины трещины на 25 - 30%. С ростом шума оценка длины трещины увеличивается. Причём это явление наблюдается как для размытых, так и для лучевых проекционных данных. При уровне шума выше чем 2% использование деконволюции посредством деления фурье-образов оказалось неэффективным из-за сложности определения параметра регуляризации.

а) б)

Рис.12. Реконструкции сечения стальной плитки с бороздами

Был разработан и опробован метод для экспериментального определения аппаратной функции томографического оборудования. Предположим, что измеряется одномерная проекция /(/) от объекта в форме ступеньки. Если источник излучения расположен строго над её срезом и находится достаточно далеко, то можно считать, что в лучевом

приближении проекция вблизи среза приблизительно равна в -функции. Тогда размытая проекция будет являться свёрткой в -функции с аппаратной функцией оборудования К{1). Взяв от неё производную, с учётом свойства в -функции, получим:

* rrj j 00

-/(/)= \-e{l-V)K(V)dr = \ 8{l-r)K{l')dl'=K{l). (22)

Предложенная методика была реализована на практике. Эксперимент проводился в Государственном институте материаловедения и прочности (г. Берлин, Германия). От стальной ступеньки высотой 5 мм было измерено 100 проекций. Профиль, полученный после усреднения, сглаживался медианным фильтром. Первая производная вычислялась разностным методом. Оценка ядра представлена на рис.11 (штриховая линия). По оси абсцисс указано количество пикселов, размер пиксела составляет 200 микрон. Амплитуда ядра нормирована на единицу. Экспериментальное ядро аппроксимировалось элементарными функциями, параметры которых находились методом наименьших квадратов. Лучшая аппроксимация была достигнута при помощи гауссовой функции. Она представлена на рис. 11 в виде сплошной линии.

Полученная оценка ядра была использована для деконволюции измеренных данных. Объектом реконструкции являлась стальная плитка с нанесёнными на её поверхность бороздами. Па рис.12,а показана томограмма, полученная без деконволюции проекционных данных. Для рис. 12,6 была произведена деконволюция с гауссовой функцией, аппроксимирующей экспериментально определённое ядро. Глубина борозд в рассматриваемом сечении составляла 0.65 мм, 0.63 мм и 0.56 мм (справа налево). Глубины борозд, измеренные по томограммам, составили 0.87 мм, 0.85 мм, 0.78 мм и 0.74 мм, 0.73 мм, 0.67 мм для рис.12,а и рис.12,б соответственно.

В заключении сформулированы основные результаты диссертации.

1. Предложен алгоритм томографической реконструкции, основанный на общей формуле обращения преобразования Радона, выраженной через потенциалы Рисса. Показано, что когда число проекций меньше 50-ти, он обеспечивает среднеквадратичную ошибку в 1.1-1.15 раз меньше, чем алгоритмы обратного проецирования с фильтрацией.

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

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

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

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

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

Список основных публикаций по теме диссертации.

[1] Лихачев А.В., Пикалов В.В. Сверхразрешение в трехмерной эмиссионной томографии при конечных параметрах детектора. // Оптика и спектроскопия. 1998. Т.85, № 3. С.490-497.

[2] Лихачёв А.В., Пикалов В.В. Синтезированный алгоритм трехмерной томографии. // Математическое моделирование. 1998. Т.10, № 1. С.73-85.

[3] Лихачёв А.В., Пикалов В.В. Трехмерная томография в диагностике газовых потоков при наличии непрозрачного тела. // Прикладная механика и техническая физика. 1998. Т.39, № 1. С. 174-180.

[4] Лихачев А.В., Маслов А.А. Миронов С.Г. Пикалов В.В. Электронно-пучковая томография плотности газа при гиперзвуковом обтекании тел. // ЖТФ. 1998. Т.68, № 4. С.125-133.

[5] Likhachov A.V., Pickalov V.V. Three-dimensional tomography with finite aperture beams. // Nuclear Instruments & Methods in Physics Research (A). 1998. V.405, N 2-3. P.506-510.

[6] Levin G.G., Vishnyakov G.N., Zakarin C.S., Likhachov A.V., Pickalov V.V., Rozinets G.I., Novoderzhkina J.K., Streletskaya E.A. Three-dimensional limited-angle microtomography of blood cells: experimental results. // Proc. SPIE. 1998. V.3261. P.159-164.

[7] Pickalov V., Likhachov A., Fuchs G., Ohl A., Roettig G. and Bandlow I. Region-of-interest tomography of fine structures in plasmas // Proceedings of the 1998ICPP & 25th EPS Conf. on Contr. Fusion and Plasma Physics. Prague. ECA. 1998. V. 22C. P.1570-1573.

[8] Likhachov A.V. and Pickalov V.V. New iterative algorithm for the SPECT diagnostics. H 1998 IEEE Nuclear Science Symposium and Medical Imaging Conference, Toronto, Canada. 1998. Abstract Book. Toronto, 1998. P.106.

[9] Лихачев A.B., Пикалов В.В. Метод подвижных сеток в алгоритмах трехмерной томографии. // Обрати, it некоррект. поставленные задачи. Москва, 1998. Тез. докл. М: МГУ, 1998. С.52.

[10] Лихачев А.В., Пикалов В.В. Определение неизвестного фона в проекционных данных для 3D томографии. // Обрати, и некоррект. поставленные задачи. Москва, 1998. Тез. докл. М: МГУ, 1998. С.51.

[11] Вишняков Г.Н., Левин Г.Г., Лихачев А.В., Пикалов В.В. Фазовая томография трёхмерных биологических мпкрообъектов: численное моделирование и экспериментальные результаты. // Оптика и спектроскопия. 1999. Т.87, № 3. С.413-419.

[12] Лихачёв А.В., Пикалов В.В. Подсеточное сглаживание в алгебраических алгоритмах трехмерной томографии. // Математическое моделирование. 1999. Т.11, № 8. С.79-90.

[13] Likhachov A.V., Pickalov V.V., Chugunova N.V., Baranov V.A. Development of iterative algorithms for industrial tomography. // Proceedings of 1st World Congress on Industrial Process Tomography, Buxton, Greater Manchester, 1999. P.463-469.

[14] Лихачев A.B., Пикалов В.В. Трёхмерная эмиссионная томография оптически плотных объектов при известном поглощении. // Оптика и спектроскопия. 2000. Т.88, № 3. С.740-749.

[15] Лихачев А.В., Пикалов В.В. Трехмерная трансмиссионная томография при рассеянии зондирующего пучка. // Четвертый сибирский конгресс по прикладной и

индустриальной математике (ИНПРИМ-2000). Тезисы докладов. Часть II. Новосибирск: Изд-во ИМ СОРАН, 2000, С.90-91.

[16] Лихачев А.В., Пикапов В.В. Новый метод в ROI томографии // Обратные и некорректно поставленные задачи. Москва, МГУ, 2000. Тез. докл. М: МАКС Пресс, 2000, С.47.

[17] Redmer В., Ewert U., Likhachov A.V., Pickalov V.V., Onel Y., Li Zheng, Baranov V.A. Sensitive Detection of Planar Defects by a Mechanised Radiometric Weld Inspection System. // Proceedings of the 15th World Conference on Non-Destructive Testing. October 2000, Rome. Paper No.IDN370.

[18] Pickalov V.V., Likhachov A.V. Iteration algorithm to correct absorption in PET. // IEEE Trans. Nucl. Sci. 2001. V.48, N 1, Pt.l. P.82-88.

[19] Пикапов B.B., Лихачев A.B. Визуализация численных расчетов трехмерных течений газа на нерегулярной сетке. // Оптические методы исследования потоков: Тез. докл. VI научно-техн. конф. (2001, МЭИ). М.: Издательство МЭИ, 2001. С.244-247.

[20] Likhachov A.V., Pickalov V.V., Ewert U., Redmer В. Influence of unsharpness on crack length evaluation with computerized tomography methods. // Proceedings of the 2nd World Congress on Industrial Process Tomography Hannover, Germany, 2001. Hannover, 2001. P.711-718.

[21] Likhachov A.V., Pickalov V.V. Scattered radiation problem in three-dimensional transmission tomography. // Proceedings of 2nd World Congress on Industrial Process Tomography. Hannover, Germany, 2001. Hannover, 2001. P.742-750.

[22] Лихачёв A.B., Пикалов В.В. Новый метод определения неизвестного аддитивного фона в проекционных данных в задаче трехмерной томографии. // ЖВММФ. 2002. Т.42, № 3. С.85-97.

[23] Лихачев А.В., Пикалов В.В. Трехмерная эмиссионная томография рассеивающей плазмы. II Оптика и спектроскопия. 2002. Т.92, № 6. С.988-999.

[24] Likhachov A.V., Pickalov V.V. Modification of Feldkamp algorithm for bifocal tomography. // Proc. 1ASTED Int. Conf, 2002, Novosibirsk Russia. Anaheim: ACTA Press, 2002. P.474-479.

[25] Likhachov A.V., Pickalov V.V. Influence of approximation of 3D tomography integral equation on inversion accuracy. // Proc. Intern. Conf. Computational Mathematics (ICCM 2002). Pt.ll. Novosibirsk: ICMMG Publ., 2002. P.607-613.

[26] Likhachov A.V., Pickalov V.V. Tomography algorithms for 3D bifocal geometry. II Intern. Conf. Ill-posed and Inverse Problems Novosibirsk, 2002. Novosibirsk: Sobolev Inst. Press, 2002. P. 106.

[27] Лихачев A.B. Измерение скорости в одномерном стационарном потоке методами эмиссионной томографии. // Прикладная механика и техническая физика. 2003. Т.44, № 2. С.83-91.

[28] Likhachov A.V., Pickalov V.V., Ewert U., Redmer В. Resolution enhancement in computerized tomography by deconvolution methods. II Proceedings of the 2nd Workshop "NDT in Progress''. International Meeting of NDT Experts, 2003, Prague. Brno: Bmo Univ. Technol., 2003. P.137-150.

[29] Redmer В., Ewert U., Radel Ch., Likhachov A.V., Vengrinovich V., Solotarev S. Planar tomography under limited view conditions for crack detection in austenitic welds. // Proceedings of the 2nd Workshop "NDT in Progress". International Meeting of NDT Experts, October 6-8, 2003, Prague. Brno: Brno Univ. Technol., 2003. P.321.

[30] Redmer В., Weise F„ Ewert U., Likhatchev A. Location of reinforcement in structures by different methods of gamma-radiography. II Proceedings of the International Symposium NonDestructive Testing in Civil Engineering, Berlin, 2003. No.v020. Berlin: BAM, 2003. P.l-6.

[31] Пикалов В.В., Лихачев A.B. Применение метода Гершберга-Папулиса в трехмерной доплеровской томографии. // Вычислительные методы и программирование. 2004. Т.5, № 2. С.27-34.

[32] Пикалов В.В., Лихачев A.B. Сравнение алгоритмов спиралыюй томографии. // Вычислительные методы и программирование. 2004. Т.5, № 2. С.51-64.

[33] Vishnyakov G.N., Levin G.G., Minaev V.L., Pickalov V.V., Likhachov A.V. Tomographie interference microscopy of living cells. // Microscopy and Analysis (UK). 2004. V.18, N 1. P. 15-17.

[34] Redmer В., Ewert U., Rädel С., Zscherpcl U., Likhatchov A.V., Pickalov V.V., Neundorf B. Mobile 3D-X-ray tomography for analysis of planar defects in welds by "TOMOCAR". // Proceedings of the 16th World Conference on Nondestructive Testing, Montreal, Canada,

2004. Montreal, 2004. P.l-8.

[35] Пикалов B.B., Лихачев A.B. Сравнение проекционных схем с прямолинейной и криволинейной траекториями источника в трехмерной томографии. // Автометрия.

2005. Т.41, № 5. С.74-80.

[36] Likhachov A.V., Shaposhnikova E.V., Trofimov O.E. On the differentiation and integration of ray transform on the unite sphere. // Proc. 2-nd IASTED - ACIT. Novosibirsk, Russia, 2005. ACTA Press. P.290-295.

[37] Шапошникова E.B., Лихачев A.B., Трофимов O.E. Программная реализация методов трехмерной томографической реконструкции в конусе лучей. // Науч.-техн. конф. "Томография". Москва, 2005. М: Из-во "Спектр", С.6.

[38] Лихачев A.B. Сравнение алгоритма Фельдкампа с алгоритмом синтеза Фурье для трёхмерной томографии. //Автометрия. 2006. Т.42, № 1. С.88-101.

[39] Лихачёв A.B. Нелинейная пороговая очистка томограмм. // Сибирский журнал индустриальной математики. 2006. Т.9, № 5. С.102-110.

[40] Лихачёв A.B. Исследование 1/z2 фильтрации в алгоритмах томографии. // Автометрия. 2007. Т.43, № 3. С.57-64.

[41] Likhachov A.V. Trofimov О. Е., Vengrinovich V.L., Zolotarev S.A. Virtual projections technique for limited angle tomography. // Proc. 5-th World Congress on Industrial Process Tomography. 2007, Bergen, Norway. Bergen: Bergen University Press, 2007. P.1038-1045.

[42] Трофимов O.E., Лихачев A.B. Сравнение некоторых алгоритмов томографической реконструкции в конусе лучей. // Сибирский журнал индустриальной математики. 2008. Т.11, № 3. С. 126-134.

[43] Лихачёв A.B. Регуляризующая фильтрация проекций в алгоритмах двумерной томографии. // Сибирский журнал вычислительной математики. 2008. Т. 11, № 2. С. 187200.

[44] Лихачев A.B. Алгоритм двойной фильтрации для двумерной томографии. // Математическое моделирование. 2009. Т.21, № 8. С.21-29.

[45] Лихачев A.B. Повышение контрастности малоракурсных томограмм, полученных алгебраическими алгоритмами реконструкции. // Вычислительные технологии. 2009. Т.14, № 3. С.38-48.

[46] Лихачев A.B. Алгоритм генерации проекций в задачах томографии с ограниченным диапазоном углов. // Автометрия. 2009. Т.45, No.l. С.83-91.

Цитируемая литература.

[1] Пикалов В.В., Мельникова Т.С. Томография плазмы. Новосибирск: Наука. Сиб. издат. фирма, 1995.

[2] Филонин О.В. Малоракурсная томография. Самара: СНЦ РАН, 2006.

3] Наттерср Ф. Математические аспекты компьютерной томографии. М.: Мир, 1990. (перевод с английского).

4] Лаврентьев М.М., Зеркаль С.М., Трофимов О.Е. Численное моделирование в томографии и условно-корректные задачи. Новосибирск: Изд-во ИДМИ НГУ, 1999.

5] Ramachandran G.N., Lakshminarayanan A.V. Three-dimensional reconstruction from radiographs and electron micrographs: application of convolutions instead of Fourier transforms. // Proc. Nat. Acad. Sci. U.S. 1971. V.68. P.2236-2240.

6] Shepp L.A. and Logan B.F. The Fourier reconstruction of a head section. II IEEE Trans. Nucl. Sci. 1974. V.21, No.3. P.21-43.

7] Ерохин B.A., Шнейдеров B.C. Трехмерная реконструкция (машинная томография). Моделирование на ЭВМ. // Препринт No.23, ЛНИВЦ; Ленинград, 1981.

8] Пикалов В.В., Преображенский Н.Г. Реконструктивная томография в газодинамике и физике плазмы. Новосибирск: Наука, 1987.

9] Ценсор Я. Методы реконструкции изображений, основанные на разложении в конечные ряды. // ТИИЭР. 1983. Т.71, No.3. С.148-160.

10] Minerbo G.N., Sanderson J.G., van Hulsteyn D.B., Lee P. Three-Dimensional Reconstruction of the X-ray emission in laser imploded targets. // Applied Optics. 1980. V.19, No. 10. P. 17231728.

11] Likhachov A.V., Pickalov V.V. Frequency Filtration in the Algebraic Algorithms of 3D-Tomography. // Proc. 1995 Intern. Meeting on Fully 3D Image Reconstruction in Radiology and Nuclear Medicine. Grenoble: Leti, 1995. P.103-107.

12] Tuy H.K. An inversion formula for cone-beam reconstruction. // SIAM J. Applied Mathematics. 1983. V.43, No.3. P.546-552.

13] Feldkamp L.A., Davis L.C., Kress J.W. Practical cone-beam algorithm. // J. Opt. Soc. Amer. A. 1984. V.l.No.6. P.612-619.

14] Yan X., and Leahy R.M. Cone beam tomography with circular, elliptical and spiral orbits. // Phys. Med. Biol. 1992. V.37, No.3. P.493-506.

15] Kudo H., Noo F., Defrise M. Cone-beam filtred-backprojeetion algorithm for truncated helical data. // Phys. Med Biol. 1998. V.43, No.10. P.2885-2909.

16] Рисс Ф., Секефальви-Надь Б. Лекции по функциональному анализу. М: Мир, 1979 (перевод с французского).

17] Novikov R. G. An inversion formula for the attenuated x-ray transformation. II Ark Mat. 2002. V.40. P.145-167.

18] Bukhgeim A.A., Kazantsev S.G. The attenuated Radon transform inversion formula for divergent beam geometry. // Proc. Int. Conf. Automation, Control, and Information Technology (IASTED 2002, Novosibirsk Russia). Anaheim: ACTA Press, 2002. P. 460-465

19] Лихачев A.B., Пикалов B.B. Трёхмерная эмиссионная томография оптически плотных объектов при известном поглощении. // Оптика и спектроскопия. 2000. Т.88, № 3. С.740-749.

20] Zaidi Н., Koral К. F. Scatter modelling and compensation in emission tomography. И Eur. J. Nucl. Med. Mol. Imaging. 2004. V.31, No.5. P.761-782.

21] Вайнберг Э.И., Казак И.А., Курозаев В.П. Реконструкция внутренней пространственной структуры объектов по интегральным проекциям в реальном масштабе времени. II ДАН СССР. 1981. T.257,No.l. С.89-94.

Автореферат докторской диссертации: Подписано в печать 18.03.2011 Формат 60x84/16. Объём 2,2 усл. печ. л. Тираж 100 экз. Заказ №3аказ № 864.

Отпечатано ЗАО РИЦ «Прайс-курьер» ул. Кутателадзе, 4г, т. 330-7202

Оглавление автор диссертации — доктора технических наук Лихачев, Алексей Валерьевич

Глава 1.

Введение 2010 год, диссертация по информатике, вычислительной технике и управлению, Лихачев, Алексей Валерьевич

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

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

При помощи томографии исследуются объекты различной природы. В соответствии с этим физические принципы, лежащие в основе систем сбора проекционных данных, могут быть весьма разнообразны. Наиболее известным приложением является медицинская, рентгеновская томография. Здесь измеряемой величиной является интенсивность проникающего излучения, прошедшего через тело пациента. Интенсивность излучения- также регистрируется в задачах дефектоскопии, диагностики плазмы и во многих других случаях. Однако существуют методы, оперирующие с данными совсем иной природы. Это, например, ЯМР-томография, электротомография, геотомография и другие.

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

Часто случается так, что по некоторым объективным обстоятельствам оказывается возможным провести только небольшое количество измерений, которых, с точки зрения существующей теории, недостаточно для качественного восстановления внутренней структуры объекта. Среди причин, в силу которых приходится использовать неполные данные (это понятие будет уточнено по ходу изложения в главах 2, 3), можно назвать следующие: ограниченный доступ к исследуемой области, сложность и высокая стоимость оборудования, вредное влияние проникающего излучения. Отсюда видно, что томография при недостатке проекционных данных представляет значительный практический интерес. Этой теме посвящено большое количество работ. В них, как. правило, предлагаются алгоритмы, позволяющие повысить качество реконструкции в некоторых частных случаях. При этом обычно, используются свойства симметрии искомых функций. Однако универсальных методов до сих пор не создано. Поэтому исследование подобных задач является весьма актуальным.

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

Тем не менее, в ряде случаев описанное выше уравнение является слишком грубым приближением. Оно предполагает, что измеряемая интегральная величина полностью определяется стационарным пространственным распределением некоторого скалярного параметра. При этом учитывается- всего лишь один процесс взаимодействия излучения с веществом: поглощение, или испускание. Иногда этого оказывается недостаточно, чтобы описать реальную ситуацию. Большой интерес вызывает эмиссионная томография среды, где также присутствует поглощение. В первую очередь это связано с развитием методов медицинской диагностики SPECT и PET, см. раздел 1.1.2. Кроме того, подобная задача возникает при изучении оптически плотной плазмы. Проблемы, связанные с наличием рассеянного излучения, также не могут быть решены в рамках этой модели сбора проекционных данных.

Использование рассматриваемого уравнения томографии подразумевает определённые свойства регистрирующего оборудования. Вообще говоря, они являются идеализированными и на практике могут быть реализованы лишь приближённо. Любой детектор имеет конечные геометрические размеры и угловую апертуру, поэтому принимает излучение из некоторого объёма пространства, а не с бесконечно тонкого луча. Это обстоятельство следует учитывать при разработке методов повышения пространственного разрешения томограмм.

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

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

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

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

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

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

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

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

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

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

3. Предложена двухпараметрическая модификация процедуры нелинейной пороговой очистки томограмм: Для её анализа, разработан метод вычисления вероятности' ошибки при определении наличия артефакта в рассматриваемой точке томограммы по статистическим характеристикам шума в проекционных данных.

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

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

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

Практическая ценность и реализация научных результатов работы.

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

Алгоритмы малоракурсной томографии, разработанные при участии автора, были использованы в научных исследованиях: для получения распределения эмиссии микропинчевого разряда (Физический институт имени П.Н. Лебедева, г. Москва); при изучении методом электронно-пучкового зондирования ударной волны, образующейся при гиперзвуковом обтекании твёрдых тел (Институт теоретической и прикладной механики СО РАН, г. Новосибирск); для реконструкции внутренней структуры клеток крови по результатам интерферометрических измерений (Всероссийский научно-исследовательский институт оптико-физических измерений, г. Москва).

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

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

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

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

Созданные в процессе проведённых исследований алгоритмы объединены в многоцелевой томографический комплекс программ с соответствующим интерфейсом, включающим графический пакет.

Значительная часть вошедших в диссертацию исследований была поддержана грантами российских и международных научных фондов, в том числе гранты Российского фонда фундаментальных исследований: N0.93-05-14137, N0.95-02-03615, 06-01-81000-Бела, N0.03-07-90060-6, N0.07-07-00251-а; грант Научного общества Нидерландов N0.713-097; грант N0.1101000 в области фундаментальных наук от Международного научного фонда.

На защиту выносятся следующие положения.

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

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

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

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

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

Апробация работы.

Основные результаты диссертационной работы докладывались на 9-й Международной конференции "Методы аэрофизических исследований" (Новосибирск, 1998), на Международной конференции "Обратные и некорректно поставленные задачи" (Москва, 1998), на 12-й Международной конференции "Photonics West'98" (Сан Жозе, Калифорния, США, 1998), на 1-м, 2-м и 5-м Международных конгрессах по индустриальной томографии (Букстон, Англия, 1999; Ганновер, Германия, 2001; Берген, Норвегия 2007), на 4-м Сибирском конгрессе по прикладной и индустриальной математике (Новосибирск, 2000), на Ежегодной конференции германского научного общества по проблемам неразрушающего контроля (Инсбрук, Германия, 2000), на 15-й и 16-й Международных конференциях по неразрушающему контролю (Рим, Италия, 2000; Монреаль, Канада 2004), на 6-й Международной научно-технической конференции "Оптические методы исследования потоков" (Москва, 2001), на Международной конференции по вычислительной математике (Новосибирск, 2002), на Международной конференции "Некорректно поставленные и обратные задачи", посвящённои 70-ти летию со дня рождения* академика М: М!. Лаврентьева (Новосибирск,

2002), на 1-й и 2-й Международных конференциях "Автоматизация, управление и информационные технологии" (Новосибирск, 2002, 2005), на 2-м Международном симпозиуме по неразрушающему контролю "NDN. in Progress" (Прага, Чехия, 2003), на Международном семинаре "Неразрушающий контроль в гражданском строительстве" (Берлин, Германия,

2003), на Научно-технической конференции "Томография" (Москва, 2005), на Международной конференции "Обратные и некорректные задачи математической физики", посвященной 75-ти летию со дня рождения академика М. М. Лаврентьева (Новосибирск, 2007). Лично автору принадлежат: идея, разработка и исследование алгоритмов двойной фильтрации и повышения контрастности, модифицированной процедуры нелинейной очистки, метода разложения проекций; обоснование алгоритма пополнения проекционных данных; развитие метода разложения в ряд Неймана; получение всех оценок быстродействия и точности, приведённых в работе; разработка моделей сбора проекционных данных в рассеивающих средах и для трансмиссионной томографии при конечной ширине пучка; научная постановка и решение рассмотренных в диссертации задач диагностики одномерного потока и экспериментального определения аппаратной функции оборудования; программная реализация исследованных алгоритмов; планирование и проведение вычислительных экспериментов.

Публикации.

По теме диссертации опубликовано 46 работ, 25 из них в рецензируемых журналах, в том числе 16 работ в журналах, рекомендованных ВАК.

Заключение диссертация на тему "Томография по неполным и искаженным данным"

Основные результаты, полученные в диссертации, состоят в следующем.

1. Предложен алгоритм томографической реконструкции, основанный на общей формуле обращения преобразования Радона, выраженной через потенциалы Рисса. Показано, что когда число проекций меньше 50-ти, он обеспечивает среднеквадратичную ошибку в 1.1 - 1.15 раз меньше, чем алгоритмы обратного проецирования с фильтрацией.

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

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

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

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

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

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

Заключение.

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

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

Библиография Лихачев, Алексей Валерьевич, диссертация по теме Математическое моделирование, численные методы и комплексы программ

1. Oldendorf W.H. Isolated flying-spot detection of radiodensity discontinuities; displaying the internal structural pattern of a complex object. // IRE Trans. Bio-Med. Electron. 1961. V.BME-8. P.68-72.

2. Cormack A.M. Representation of a function by its line integrals, with some radiological applications. II J. Appl. Phys. 1963. V.34. P. 1722-1727.

3. Hounsfield G.N. Computerized transverse axial scanning tomography: Part 1, description of the system. I/ Br. J. Radiol. 1973. V.46. P. 1016-1022.

4. Хермен Г.Т. Восстановление изображений по проекциям. Основы реконструктивной томографии. М.: Мир, 1983 (перевод с английского).

5. Kniipfer W., Hell Е., Mattern D. Novel X-ray detectors for medical imaging. // Nuclear Physics B. 1999. V.78. P.610-615.

6. Kotter E., Langer M. Digital radiography with large-area flat-panel detectors. // Eur. Radiol. 2002. № 12. P.2562-2570.

7. Verdonck В., Bloch I., Maitre H. Spiral CT: A survey on accuracy, resolution and artifacts. II Appl. Sign. Proc. 1998. V.5, № 1. P.3-23.

8. Garvey C. J., Hanlon R. Computed tomography in clinical practice. // Br. Med. J. 2002. V.324. P.1077-1080.

9. Del Sole A., Falini A., Ravasi L. et al. Anatomical and biochemical investigation of primary brain tumours. II Eur. J. Nucl Med. 2001. V.28. P.1851-1872.

10. John C., Eisner E., Mtiller A. et al. Computed tomography in acute cerebral ischemia. // Radiologe. 1997. V.37. P.853-858.

11. Hansell D.M. Computed tomography of diffuse lung disease: functional correlates. U Eur. Radiol. 2001. № 9. P. 1666-1680.

12. Hughes S.W., D'Arcy T.J., Maxwell D.J. et al. In vitro estimation of foetal liver volume using ultrasound, x-ray computed tomography and magnetic resonance imaging. /1 Physiol. Measur. 1997. V.18, № 4. P.401-410.

13. Ollinger J.M., Fessler J.A. Positron emission tomography. // IEEE Signal Processing Magazine. 1997. V.14, № 1. P.43-55.

14. Ussov W.Y., Riannel J.E., Barysheva E.V. et al. Role of single photon emission tomography with 99mTc-MIBI in diagnosis of metastatic widespread in breast cancer. // European Journal of Cancer. 1998. V.34. P.S19-S20.

15. Majumder D. D., Bhattacharya M. Frank H. George Research Award Winning Paper. Cybernetic approach to medical technology: application to cancer screening and other diagnostics. // Kybernetes. 2000. V.9, № 7-8. P.871-895.

16. Moses W.W. Trends in PET imaging. // Nuclear Instruments and Methods in Physics Research. Section A. 2001. V.471, № 1-2. P.209-214.

17. Del Guerra A., Di Domenico G., Fantini A. et al. A dedicated system for breast cancer study with combined SPECT-CT modalities. // Nuclear Instruments and Methods in Physics Research. Section A. 2003. V.497, № 1. P. 129-134.

18. Kallen K., Burtscher I.M., Holtas S. et al. 201Thallium SPECT and 1H-MRS compared with MRI in chemotherapy monitoring of high-grade malignant astrocytomas. II Journal of N euro-Oncology. 2000. V.46, № 2. P.173-185.

19. Grosu A.-L., Feldmann H.J., Dick S. et al. Implications of IMT-SPECT for postoperative radiotherapy planning in patients with gliomas. // International Journal of Radiation Oncology *Biology*Physics. 2002. V.54, № 3. P.842-854.

20. Yao Z., Liu X. J., Shi R. et al. A comparison of 99mTc-MIBI myocardial SPET with electron beam computed tomography in the assessment of coronary artery disease. II Eur. J. Nucl. Med. 1997. V.24, № 9. P.l 115-1120.

21. Tillfors M., Furmark Т., Marteinsdottir I. Cerebral blood flow in subjects with social phobia during stressful speaking tasks: a PET study. // American Journal of Psychiatry. 2001. V.158, № 8. P.1220-1226.

22. Videbech P., Ravnkilde В., Pedersen Т.Н. et al. The Danish PET/depression project: clinical symptoms and cerebral blood flow. A regions-of-interest analysis. // Acta Psychiatrica Scandinavica. 2002. V.106, № 1. P.35-44.

23. Пикалов B.B., Мельникова T.C. Томография плазмы. Новосибирск: Наука. Сиб. издат. фирма, 1995.

24. Minerbo G.N., Sanderson J.G., van Hulsteyn D.B. et al. Three-Dimensional Reconstruction of the X-ray emission in laser imploded targets. // Applied Optics. 1980. V.19, № 10. P.1723-1728.

25. Chen Y.-W., Miyanaga N., Yamanaka M. et al. Three-Dimensional imaging of laser imploded targets. II J. Appl. Phys. 1990. V.68, № 4. P.1483-1488.

26. Chen Y.-W., Yamanaka M., Miyanaga N. et al. Three-dimensional reconstruction of laser-irradiated targets using URA coded aperture cameras. // Optics Communications. 1989. V.71, № 5. P.249-254.

27. Escande D.F., Cappello S., D'Angelo F. et al. Single helicity: a new paradigm for the reversed field pinch. // Plasma Phys. & Contr. Fusion. 2000. V.42,№ 12B. P.243-253.

28. Stehle C., Busquet M. et al. On Stark broadening as a tool for diagnostics of high density plasmas. // Laser and Particle Beams. 2005. V.23, № 3. P.357-363.

29. Antony M., Weiseny H., Dutchy M. J. et al. X-ray tomography on the TCV tokamak. IIPlasma Phys. & Contr. Fusion. 1996. V.38, № 11. P.1849-1878.

30. Peysson Y., Coda S., Imbeaux F. Hard X-ray CdTe tomography of tokamak fusion plasma. // Nuclear Instruments and Methods in Physics Research Section• A. 2001. V.458, № i2. P.1849-1878.

31. Преображенский Н.Г. (Отв. ред.) Инверсия Абеля и её обобщения. Новосибирск: ИТПМ, 1978.

32. Yamaguchi S., Igami Н., Tanaka Н. et al. Observation of sawtooth crashes by a multi-toroidally positioned soft x-ray computer tomography system in the WT-3 tokamak. // Plasma Phys. & Contr. Fusion. 2004. V.46, № 8. P. 1163-1180

33. Bates D.R. A suggestion regarding the use rockets to vary the amount of atmospheric sodium. II J. Geophys. Res. 1951. V.55, № 3. P.347-349.

34. Алпатов В.В., Левин Г.Г., Пикалов В.В. и др. Оптическая томография искусственных образований в околоземной среде. // Космические исследования. 1993. Т.31, Вып.1. С.121-134.

35. Williamson J.H., Clarke М.Е. Construction of electron distribution function from laser scattering spectra. // J. Plasma Physics. 1971. V.6, № 1. P.211-221.

36. Kinsey J.L. Fourier transform Doppler spectroscopy: A new means of obtaining velocity-angle distributions in scattering experiments. // J. of Chemical Physics. 1977. V.66, № 6. P.2560-2565.

37. Koslover R., McWilliams R. Measurement of multidimensional ion velocity distributions by optical tomography. // Rev. Sci. Instrum. 1986. V.57, № 10. P.2441-2448.r

38. Kontrym-Sznajd G., Samsel-Czekala M. New reconstruction method of electron momentum density from Compton profiles. // Appl. Phys. (A). 2000. V.70, № 1. P.89-92.

39. Samsel-Czekala M., Kontrym-Sznajd G. et al. Electron momentum density in CU0.9AI0.1. II Appl. Phys. (A). 2003. V.76, № 1. P.87-92.

40. Matsumoto К., Mennickent R.E. On the secondary star of the galactic supersoft X-ray source RXJ0925.7-4759 and Doppler tomography of the emission-lines. // Astron. Astrophys. 2000. V.356, № 2. P.579-584.

41. Torres M.A.P., Casares J., Martinez-Pais I.G., Charles P.A. Rotational broadening and Doppler tomography of the quiescent X-ray nova Centaurus X-4. // Mon. Not. R. Astron. Soc. 2002. V.334, № 1. P.233-240.

42. Challenor P.G., Cipollini G., Cromwell D. Use of 3D Radon transform to examine the properties of oceanic Rossby Waves. // J. Atmospheric and Oceanic Technology. 2001. V.18, № 6. P. 1558-1566.

43. Lindmo Т., Smithies D.J., Chen Z., Nelson J.S., Milner Т.Е. Accuracy and noise in optical Doppler tomography studied by Monte Carlo simulations. // Phys. Med. Biol 1998. V. 43, № 10. P.3045-3064.

44. Yang V.X.D., Gordon M.L., Мок A. et al. Improved phase-resolved optical Doppler tomography using Kasai velocity estimator and histogram segmentation. // Optics Communications. 2002. V.208, № 4-6. P.209-214.

45. Левин Г.Г., Вишняков Г.Н. Оптическая томография. М: Радио и связь, 1989.

46. Мишин Г.И. (Отв. ред.) Оптические методы исследований в баллистическом эксперименте. JI: Наука, 1979.

47. Вест Ч. Голографическая интерферометрия. М: Мир, 1982 (перевод с английского).

48. Vukicevic D., Jäger Н., Philipp Н. et al. Tomographie reconstruction of the temperature distribution in a convective heat flow using multidirectional holographic interferometry. // Appl. Optics. 1989. V.28, № 8. P. 1508-1516.

49. Snyder R., Hesselink L. Optical tomography for flow visualization of the density field around a revolving helicopter rotor blade. // Appl. Optics. 1988. V.20, № 23. P.3650-3656.

50. Kittelton J.K., Yu Y.H. Reconstruction of a three-dimensional transonic rotor flow field from holographic interferogram data. // AIAA Paper. 1985, 85-0370.

51. Modarress D., Tan H., Trolinger J. Tomographic reconstruction of three-dimensional flow over airfoils. II AIAA Paper. 1985, 85-0479.

52. Timmerman B.H., Watt D.W. Tomographic high-speed digital holographic interferometry. II Meas. Sei. Technol. 1995, V.6, № 9. P.1270-1277.

53. Watt D.W., Timmerman B.H., Biyanston-Cross P.J. Quantitative visualization of high-speed 3D,turbulent flow structures using holographic interferometric tomography. // Optics & Laser Technology. 1999, V.31, № 1. P.53-65.

54. Апресян JT.А., Кравцов Ю.А. Теория переноса излучения. М.: Наука, 1983.

55. Funk Р. Über eine Geometrische Anwendung der Abelschen Integralgleichung. II Math Ann. 1916. V.77. P.129-135.

56. Йон Ф. Плоские волны и сферические средние в применении к дифференциальным уравнениям в частных производных. М.: Иностранная Литература, 1958 (перевод с английского).

57. Helgason S. Differential operators on homogeneous spaces. 11 Acta Math. 1959. V.102. P.239-299.

58. Helgason S. A duality in integral geometry; some generalizations of the Radon transform. II Bull. Amer. Math. Soc. 1964. V.70. P.435-436.

59. Наттерер Ф. Математические аспекты компьютерной томографии. М.: Мир, 1990. (перевод с английского).64., Shepp L.A., Logan B.F. The Fourier reconstruction of a head section. // IEEE Trans. Nucl. Sei. 1974. V.21, № 3. P.21-43.

60. Ramachandran G.N., Lakshminarayanan A.V. Three-dimensional reconstruction from radiographs and electron micrographs: application of convolutions instead of Fourier transforms. // Proc. Nat. Acad. Sci. U.S. 1971. V.68. P.2236-2240.

61. Herman G.T., Rowland S.W. Three methods for reconstructing objects from x-rays: a comparative study. // Cumput. Graph. Img. Proc. 1973. № 2. P.151-178.

62. Ерохин B.A., Шнейдеров B.C. Трехмерная реконструкция (машинная томография). Моделирование на ЭВМ. // Препринт № 23, ЛНИВЦ, 1981.

63. Луитт P.M. Алгоритмы реконструкции с использованием интегральных преобразований. // ТИИЭР. 1983. Т.71, № 3. С.125-147.

64. Лаврентьев М.М., Зеркаль С.М., Трофимов О.Е. Численное моделирование в томографии и условно-корректные задачи. Новосибирск: Изд-во ИДМИ НГУ, 1999.

65. Лихачёв А.В. Исследование l/z2 фильтрации в алгоритмах томографии. // Автометрия. 2007. Т.43, № 3. С.57-64.

66. Орлов С.С. Теория трехмерной реконструкции. // Кристаллография. 1975. Т.20,№З.С.511-515.

67. Ra J.B., Cho Z.N. Generalized true three-dimensional reconstruction algorithms. II Proc. IEEE. 1981. V.69. P.668-670.

68. Пикалов B.B., Преображенский Н.Г. Реконструктивная томография в газодинамике и физике плазмы. Новосибирск: Наука, 1987.

69. Вишняков Г.Н. Восстановление томограмм трёхмерных объектов по двумерным проекциям. // Оптика и спектроскопия. 1988.Т.65,№ 3.C.677-683.

70. Кириллов А.А. Об одной проблеме И.М. Гельфанда. // ДАН. 1961. Т. 137, № 2. С.276-277.

71. Tuy Н.К. An inversion formula for cone-beam reconstruction. // SIAM J. Applied Mathematics. 1983. V.43, № 3. P.546-552.

72. Likhachov A.V., Shaposhnikova E.V., Trofimov O.E. On the differentiation and integration of ray transform on the unite sphere. // Proc. 2-nd IASTED -ACIT. Novosibirsk, Russia, 2005. ACTA Press. P.290-295.

73. Minerbo G.N. Convolution reconstruction from cone-beam projection data. // IEEE Trans. Nucl. Sci. 1979; V. 26, № 2. P.2682-2684.

74. Hamaker C., Smith K.T., Solmon D.C. et al. The divergent beam X-ray transform. II Rocky Mount. J. of Math. 1980. V.10, № 1. P.253-283.

75. Smith B.D. Image reconstruction from cone-beam projections: necessary and sufficient conditions and reconstruction methods. // IEEE Trans. Med. Image. 1985. №4. P.14-28.

76. Smith B.D. Cone-Beam tomography: recent advances and tutorial review. // Optical Engineering. 1990. V.29, № 5. P.524-534.

77. Трофимов O.E. К задаче восстановления функции трех переменных по ее интегралам вдоль прямых, пересекающих заданную кривую. // Автометрия. 1991. Т.27, № 2. С.55-61.

78. Zeng G.L., Gullberg G.T. A cone-beam tomography algorithm for orthogonal circle-and-line orbit. II Phys. Med. Biol. 1992. V.37, № 3. P.563-577.

79. Taguchi K., Zeng G.L., Gullberg G.T. Cone-beam image reconstruction using spherical harmonics. II Phys. Med. Biol. 2001. V.46, № 6. P. 127-138.

80. Grangeat P. Mathematical framework of cone-beam 3D reconstruction via the first derivative of the Radon transform. // Proc. conf. Mathematical methods in tomography. Oberwolfach, Germany. 1990, P.66-97.

81. Feldkamp L.A., Davis L.C., Kress J.W. Practical cone-beam algorithm. // J. Opt. Soc. Amer. A. 1984. V.l, № 6. P.612-619.

82. Shih A., Wang G., Cheng P.C. Fast algorithm for X-ray cone-beam microtomography. //Microsc Microanal. 2001 V.7. № 1. P.13-23.

83. Likhachov A.V., Pickalov V.V. Modification of Feldkamp algorithm for bifocal tomography. // Proc. IASTED Int. Conf., 2002, Novosibirsk, Russia. Anaheim: ACTA Press, 2002. P.474-479.

84. Gregor J., Gleason S.S. et al. Fast Feldkamp reconstruction based on focus of attention and distributed computing. // IJIST. 2002. V.12, № 6. P.229-234.

85. Xiao S., Bresler Y., Munson D. C. Jr. Fast Feldkamp algorithm for cone-beam computer tomography. // ICIP. 2003. № 2. P.819-822.

86. Лихачёв A.B. Сравнение алгоритма Фельдкампа с алгоритмом синтеза Фурье для трёхмерной томографии. // Автометрия. 2006.Т.42,№ 1.С.88-101.

87. Defrise М., Kinahan P. Е. et al. Exact and approximate rebinning algorithms for 3-D PET data. II IEEE Trans. Med. Imag. 1997. V.16, № 2. P.145-158.

88. Ни H. Multi-slice helical CT: Scan and reconstruction. // Med. Phys. 1999. V.26, № 1. P.5-18.

89. Noo F., Defrise M., Clackdoyle R. Single-slice rebinning method for helical cone-beam СТ. II Phys. Med. Biol. 1999. V.44, № 2. P. 561-570.

90. Defrise M., Noo F., Kudo H. Improved two-dimensional rebinning of helical cone-beam computerized tomography data using John's equation. // Inverse Problems. 2003. V.l9. P.41-54.

91. Tam K.C., Samarasekera S., Sauer F. Exact cone beam CT with spiral scan. // Phys. Med. Biol. 1998. V.43, № 4. P.1015-1024.

92. Kudo H., Noo F., Defrise M. Cone-beam filtred-backprojection algorithm for truncated helical data. И Phys. Med. Biol. 1998. V.43, № 10. P.2885-2909.

93. Defrise M., Noo F., Kudo H. A solution to the long-object problem in helical cone-beam tomography. II Phys. Med. Biol. 2000. V.45, № 3. P.623-643.

94. Tretiak O.J., Metz C. The exponential Radon transform. // SI AM J. Appl Math. 1980. V.39, № 2. P.341-354.

95. Inouye T., Kose K., Hasegawa A. Image reconstruction algorithm for single photon emission computed tomography with uniform attenuation. // Phys. Med. Biol. 1989. V.34, № 3. P.299-304.

96. Metz C.E., Pan X.A unified analysis of exact methods of inverting the 2D exponential Radon transform, with implications for noise control in SPECT. // IEEE Trans. Med: Imag. 1995. V.14. P.643-658.

97. Kuchment P., Shneiberg I. Some inversion formulae in the single photon emission computed tomography. II Appl. Anal. 1994. V.53. P.221-231.

98. Wagner J.-M., Noo F., Clackdoyle R. Exact inversion of the exponential x-ray transform for rotating slant-hole (RSH) SPECT. II Phys. Med. Biol. 2002. V.47, № 15. P.2713-2726.

99. Wagner J.-M., Noo F., Clackdoyle R. et al. Attenuation correction for rotating slant-hole (RSH) SPECT using exact rebinning. II IEEE Trans. Nucl. Sci. 2003. V.50. P.110-116.

100. Kunyansky L.A. Inversion of the 3D exponential parallel-beam transform and the Radon transform with angle-dependent attenuation. // Inverse Problems. 2004. V.20, № 5. P. 1455-1478.

101. Bellini S., Piacentini M., Cafforio C. et al. Compensation of tissue absorption in emission tomography. // IEEE Trans. Acoust. Speech Signal Process. 1979. V.27, № 3. P.213-218.

102. Arbuzov E.V., Bukhgeim A.L. Kazantsev S.G. 1998 Two-dimensional tomography problems and the theory of ^-analytic functions. // Siberian Adv. Math. V.8. P. 1-20.

103. Natterer F. Inversion of the attenuated Radon transform. // Inverse Problems. 2001. V.17,№ l.P. 113-119.

104. Novikov R.G. An inversion formula for the attenuated x-ray transformation. //

105. Ark. Mat. 2002. V.40. P. 145-167.

106. Bukhgeim A.A., Kazantsev S.G. The attenuated Radon transform inversion formula for divergent beam geometry. // Proc. Int. Conf. Automation, Control, and Information Technology (IASTED 2002, Novosibirsk, Russia). Anaheim: ACTA Press, 2002. P. 460-465

107. Likhachov A.V., Pickalov V.V. Influence of approximation of 3D tomography integral equation on inversion accuracy. // Proc. Intern. Conf. Computational Mathematics (ICCM2002). Pt.II. Novosibirsk: ICMMG Publ., 2002. P.607-613.

108. Herman G.T., Lent A. Iterative Reconstruction Algorithms. // Comput. Biol, and Med. 1976. V.6. P.273-294.

109. Ценсор Я. Методы реконструкции изображений, основанные на разложении в конечные ряды. // ТИИЭР. 1983. Т.71, № 3. С. 148-160.

110. Herman, G.T. Image reconstruction from projections. // Real-Time Imaging. 1995. V.l.P.3-18.

111. Kaczmarz S. Angenaehrte Aufloesung von Systemen Linearer Gleichunger. // Bull. Int. Acad. Pol. Sci. Lett. 1937. V.A. P.355-357.

112. Gordon R, Bender R., Herman G.T. Algebraic reconstruction techniques (ART) for tree-dimensional electron microscopy and X-ray photography. //

113. J. Teor. Biol. 1970. V.29. P.471-481.

114. Minerbo G. MENT: A maximum entropy algorithm for reconstructing a source from projection data. // Comput. Graph. Image Proc. 1979. V.10. P.48-68.

115. Minerbo G. Maximum entropy reconstruction from cone-beam projection data. // Comput. Biol, and Med. 1979. V.9. P.29-37.

116. Alenius S., Ruotsalainen U. Bayesian image reconstruction for emission tomography based on median root prior. // Eur. J. Nucl. Med. 1997. V.24, № 3. P.258-265.

117. Seret A. Median root prior and ordered subsets in Bayesian image reconstruction of single-photon emission tomography. И Eur. J. Nucl. Med. 1998. V.25, № 3. P.215-219

118. Siltanen S., Kolehmainen V., Jarvenpaa S. et* al. Statistical inversion for medical x-ray tomography with few radiographs: I. General theory. // Phys. Med. Biol 2003. V.48, № ю. P.1437-1463

119. Niinimaki K., Kolehmainen V., Siltanen S. Bayesian multiresolution method for local tomography in dental x-ray imaging. // Phys. Med. Biol 2007. V.52, №22. P.6663-6678.

120. Тихонов A.H., Арсенин В.Я. Методы решения некорректных задач. М.: Наука, 1986.

121. Тихонов А.Н., Арсенин В.Я., Тимонов А.А. Математические задачи компьютерной томографии. М: Наука, 1987.

122. Бакушинский А.Б., Гончарский А.В. Итеративные методы решения некорректных задач. М.: Наука, 1989.

123. Iwama N., Yoshida Н., Takimoto Н. et al. Phillips-Tikhonov regularization of plasma image reconstruction with the generalized cross validation. // Appl. Phys. Lett. 1989. Y.54, № 6. P.502-504.

124. Bhatia M., Karl W.C., Willsky A.S. A wavelet-based method for multiscale tomographic reconstruction. IIIEEE Trans. Med. Imag. 1996.V.15,№ l.P.92-101.

125. Gourion D., Noll D. The inverse problem of emission tomography. // Inverse Problems. 2002. V.18, № 5. P.1435-1460.

126. Kalifa J., Laine A., Esser P.D. Regularization in tomographic reconstruction using thresholding estimators. // IEEE Trans. Med. Imag. 2003. V.22, № 3. P.351-359.

127. Zhang-O'Connor Y., Fessler J. A. Fast predictions of variance images for fan-beam transmission tomography with quadratic regularization. // IEEE Trans. Med. Imag. 2007. V.26, № 3. P. 335-346.

128. Филонин O.B. Малоракурсная томография. Самара: СНЦРАН, 2006.

129. Konovalov A.B., Kiselev A.N., Vlasov V.V. Spatial resolution of few-view computed tomography using algebraic reconstruction techniques. // Pattern Recognition and Image Analysis. 2006; V. 16; № 2. P.249-255".

130. Пикапов B.B., Казанцев > Д.И. Свойства регуляризованного алгоритма Гершберга-Папулиса в задаче веерной томографии. // Вычислительные технологии. 2008. Т.13, № 6. С.121-133.

131. Лихачёв A.B. Регуляризующая фильтрация проекций в алгоритмах двумерной томографии. // Сибирский журнал вычислительной математики. 2008. Т.11, № 2. С. 187-200.

132. Морозов В.А. О принципе невязки при решении операторных уравнений методом регуляризации. ИЖВММФ. 1968. Т.8, № 2. С.295-309.

133. Ягола А.Г. О выборе параметра регуляризации по обобщённому принципу невязки. // ДАН СССР. 1979. Т.245, № 1. С.37-39.

134. Воскобойников Ю.Е. Устойчивые методы и алгоритмы параметрической идентификации. Новосибирск: Изд-во НГАСУ, 2006.

135. Василенко В.А. Сплайн-функции: теория, алгоритмы, программы. Новосибирск: Наука, 1983.

136. Лихачёв A.B., Пикалов В.В. Синтезированный алгоритм трехмерной томографии. // Математическое моделирование. 1998. Т. 10, № 1. С.73-85.

137. Лихачёв A.B. Алгоритм двойной фильтрации для двумерной томографии. // Математическое моделирование. 2009. Т.21, № 8. С.21-29.

138. Лихачёв A.B. Повышение контрастности малоракурсных томограмм, полученных алгебраическими алгоритмами реконструкции. // Вычислительные технологии. 2009. Т. 14, № 3. С.38-48.

139. Лихачёв A.B. Нелинейная пороговая очистка томограмм. // Сибирский журнал индустриальной математики. 2006. Т.9, № 5. С.102-110.

140. Лихачёв A.B. Измерение скорости в одномерном стационарном потокеметодами эмиссионной томографии. // Прикладная механика и техническая физика. 2003. Т.44, No.2. С.83-91.

141. Milanfar P., Karl W.C., Willsky A.S. A moment-based variational approach to tomographic reconstruction. II IEEE Trans. Img. Proc. 1996. V.5, № 3.P.459-470.

142. Delaney A.H., Bresler Y. Globally convergent edge-preserving regularized" reconstruction: An application to limited-angle tomography. /I IEEE Trans. Image Process. 1998. V.7, № 2. P.204-221.

143. Quintoy E. T. Exterior and limited-angle tomography in non-destructive evaluation. I/ Inverse Problems. 1998. V.14. P.339-353.

144. Hanson K.M., Wecksung G.W. Bayesian approach to limited-angle reconstruction in computed tomography. // Journal of the Optical Society of America. 1983. V.73. P.1501-1509.

145. Persson M., Bone D., Elmqvist H. Total variation norm for three-dimensional iterative reconstruction in limited view angle tomography. // Phys. Med. Biol. 2001. V.46, № 3, P.853-866

146. Rantala M., Yanska S. et al. Wavelet-based reconstruction for limited-angle X-ray tomography. II IEEE Trans. Med. Image. 2006. V.25, № 2. P.210-217.

147. Defrise M., De Mol C. A regularized iterative algorithmfor limited-angle inverse Radon transform. // Optica Acta. 1983. V.30, № 4. P.403-408.

148. Вишняков Г.Н., Гильман Г.А., Левин Г.Г. Восстановление томограмм при ограниченном числе проекций. Итерационные методы. // Оптика и спектроскопия. 1985. Т.58, № 2. С.406-413.

149. Melnikova T.S., Pickalov V.V. Computer-aided plasma tomography. // High temperature dust laden jets in plasma technology. Proc. Inter. Workshop 1988, Novosibirsk. Utrecht: VSP, 1989. P.257-282.

150. Бронников A.B., Воскобойников Ю.Е. Преображенский Н.Г. Итерационные алгоритмы в задачах томографии полупрозрачных сред. // Препринт № 18, ИТПМ; Новосибирск, 1989.

151. Вишняков Г.Н., Левин Г.Г., Лихачёв A.B. и др. Фазовая томография трёхмерных биологических микрообъектов: численное моделирование и экспериментальные результаты. // Оптика и спектроскопия. 1999. Т.87, №3. С.413-419.

152. Лихачёв A.B. Алгоритм пополнения проекционных данных в задачах томографии с ограниченным диапазоном углов обзора. // Автометрия. 2009. Т.45, № 1. С.83-91.

153. Vishnyakov G., Levin G., Zakarin С. Interferometric computed-microtomography of 3D phase objects. // Proc. SPIE. 1997. V.2984. P.67-71.

154. Levin G.G., Vishnyakov G.N., Zakarin C.S. et al. Three-dimensional limited-angle microtomography of blood cells: experimental results. // Proc. SPIE. 1998. V.3261. P.159-164.

155. Vishnyakov G.N., Levin G.G., Minaev V.L., Pickalov V.V., Likhachov A.V. Tomographic interference microscopy of living cells. // Microscopy and Analysis (UK). 2004. V.18, № 1. P.15-17.

156. Redmer В., Ewert U., Likhachov A.V. et al. Sensitive Detection of Planar Defects by a Mechanised Radiometric Weld Inspection System. // Proc. 15th World Conf. Non-Destructive Testing. 15-21 October 2000, Rome. № IDN370.

157. Ewert U., Röbbel J., Bellon. C. et al. Digital Laminography. // Materialpruefung. 1995. V.37, № 6. P.218-224.

158. Пикалов В.В., Лихачёв А.В. Сравнение проекционных схем с прямолинейной и криволинейной траекториями источника в трехмерной томографии. II Автометрия. 2005. Т.41, № 5. С.74-80.

159. Трофимов О.Е., Лихачёв А.В. Сравнение некоторых алгоритмов томографической реконструкции в конусе лучей. // Сибирский журнал индустриальной математики. 2008. Т. 11, № 3. С.126-134.

160. Пикалов В.В., Лихачёв А.В. Сравнение алгоритмов спиральной томографии. И Вычислительные методы и программирование. 2004. Т.5, №2. С.51-64.

161. Yan X., and Leahy R.M. Cone beam tomography with circular, elliptical and spiral orbits. IIPhys. Med. Biol. 1992. V.37, № 3. P.493-506.

162. Лихачёв A.B., Пикалов B.B. Трехмерная томография в диагностике газовых потоков при наличии непрозрачного тела. // Прикладная механика и техническая физика. 1998. Т.39, № 1. С.174-180.

163. Лихачев А.В., Маслов А.А. Миронов С.Г. и др. Электронно-пучковая томография плотности газа при гиперзвуковом обтекании тел. // ЖТФ. 1998. Т.68, № 4. С.125-133.

164. Лихачёв А.В., Пикалов В.В. Новый метод определения неизвестного аддитивного фона в проекционных данных в задаче трехмерной томографии. IIЖВММФ. 2002. Т.42, № 3. С.85-97.

165. Likhachov А.У., Pickalov V.V. Scattered radiation problem in three-dimensional transmission tomography. // Proc. 2nd World Congress on Industrial Process Tomography. Hannover, Germany, 29-31 Aug 2001. P.742-750.

166. Pickalov V.V., Likhachov A.V. Iteration algorithm to correct absorption in PET. //IEEE Trans. Nucl. Sci. 2001. V.48, № 1, Pt.l. P.82-88.

167. Рисс Ф., Секефальви-Надь Б. Лекции по функциональному анализу. М: Мир, 1979 (перевод с французского).

168. Lewitt R.M., Muehllehner G. Accelerated iterative reconstruction for positron emission tomography based on the EM algorithm for maximum likelihood estimation. ¡/IEEE Trans. Med. Imaging. 1986. № 5. P.16-22.

169. Zeng G.L., Gullberg G.T., Tsui B.M.W. et al. Three-dimensional iterative reconstruction algorithms with attenuation and geometric point response correction. И IEEE Trans. Nucl. Sci. 1991. V.38, № 2. P.693-702.

170. Pelegrini M., Buvat I., Benali H. et al. A spline-regularized minimal residual algorithm for iterative attenuation correction in SPECT. // Phys. Med. Biol. 1999. V.44, № 10. P.2623-2642.

171. Lalush D.S., Frey E.C., Tsui B.M.W. Fast maximum entropy approximation in SPECT using the RBI-MAP algorithm. // IEEE Trans. Med. Imag. 2000. V.19, № 4. P.286-294.

172. Лихачёв A.B., Пикалов В.В. Трёхмерная эмиссионная томография оптически плотных объектов при известном поглощении. // Оптика и спектроскопия. 2000. Т.88, № 3. С.740-749.

173. Barrett Н.Н., Swindell W. Radiological Imaging, Volume 2. New York: Academic Press, 1981.

174. Левин Г.Г., Старостенко О.В. О возможности томографических исследований рассеивающих сред. // Линейные и нелинейные задачи вычислительной томографии. Новосибирск, 1985. С 86-99.

175. Аниконов Д.С., Ковтанюк А.Е., Прохоров И.В. Использование уравнения переноса в томографии. М.: Логос, 2000.

176. Любимов В.В. Основы флуоресцентной лазерной томографии сильнорассеивающих сред. // Оптика и спектроскопия. 2000. Т.88, № 2. С.321-324.

177. Кравценюк О.В., Любимов В.В. Особенности статистических характеристик траекторий фотонов в сильнорассеивающей среде вблизи поверхности объекта. // Оптика и спектроскопия. 2000. Т.88, №4.С.670-676.

178. Levin C.S., Dahlbom М., Hoffman E.J. A Monte Carlo correction for the effect of Compton scattering in 3-D PET brain imaging. // IEEE Trans. Nucl. Sci. 1995. V.42. P.l 181-1188.

179. Beekman F.J., Kamphuis C., Frey E.C. Scatter compensation methods in 3D iterative SPECT reconstruction: a simulation study. // Phys. Med. Biol. 1997. V.42. P.1619—1632.

180. Werling A., Bublitz O., Doll J. et al. Fast implementation of the single scatter simulation algorithm and its use in iterative image reconstruction of PET data. // Phys. Med. Biol. 2002. V.47 P.2947-2960.

181. Zaidi H., Koral K. F. Scatter modelling and compensation in emission tomography. II Eur. J. Nucl. Med. Mol. Imaging. 2004. V.31, № 5. P.761-782.

182. Лихачев A.B., Пикалов В.В. Трехмерная эмиссионная томография рассеивающей плазмы. // Оптика и спектроскопия. 2002.Т.92,№6.С.988-999.

183. Пикалов В.В., Чугунова Н.В. Широкоапертурная томография эмиссионных объектов. // Оптика и сие/ст/юскога/я.2000.Т.88,№2.С.325-329;

184. Вайнберг Э.И., Казак И.А., Курозаев В.П. Реконструкция внутренней пространственной структуры объектов по интегральным проекциям в реальном масштабе времени. И ДАН СССР. 1981. Т.257, № 1. С.89-94.

185. Лихачёв А.В., Пикал ов В.В. Сверхразрешение в трехмерной эмиссионной томографии при конечных параметрах детектора. // Оптика и спектроскопия. 1998. Т.85, № 3. С.490-497.

186. Likhachov A.V., Pickalov V.V. Three-dimensional tomography with finite aperture beams. // Nuclear Instruments and Methods in Physics Research Section A. 1998. V.405, № 2-3. P.506-510.

187. Likhachov A.V., Pickalov V.V., Ewert U. et al. Influence of unsharpness on crack length evaluation with computerized tomography methods. // Proc. 2nd World Congress on Industrial Process Tomography Hannover, Germany, 29-31 Aug 2001. P.711-718.

188. Bertero M., Bindi D., Boccacci P., et al. A novel blind-deconvolution method with an application to seismology. // Inverse Problems. 1998. V.14. P.815-833.

189. Burger M., Scherzer O. Regularization method for blind deconvolution and blind source separation problems. // Math. Control Signal Systems. 2001. V.14. P.358-383.