Итерационные алгебраические методы реконструкции изображения

Итерационные алгебраические методы реконструкции изображения

Дипломная работа

Итерационные алгебраические методы реконструкции изображения

Студент гр. САУ-05-1

Климов А.М.

Руководитель работы

проф. Литвин О.Г.

2009 р.

РЕФЕРАТ

Пояснительная записка: __ страницы, 31 рисунок, 7 таблиц, 8 источников, 2 приложения.

Бакалаврская аттестационная работа посвящена исследованию итерационных алгебраических методов реконструкции изображения.

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

Предметом исследования является итерационный алгебраический метод ART .

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

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

ПРЕОБРАЗОВАНИЕ РАДОНА, КОЭФФИЦИЕНТ РЕЛАКСАЦИИ, МЕТОД РЕЛАКСАЦИИ, МЕТОД ИТЕРАЦИЙ, АЛГЕБРАИЧЕСКИЙ МЕТОД, ДИСКРЕТИЗАЦИЯ ОБЛАСТИ, ЭЛИЗ, МЕТОД КАЧМАЖА

СОДЕРЖАНИЕ

ВВЕДЕНИЕ

1 СИСТЕМНАЯ МОДЕЛЬ СЛОЖНОЙ ОРГАНИЗАЦИОННОЙ СИСТЕМЫ «НЕВРОЛОГИЧЕСКАЯ ЛЕЧЕБНО-ДИАГНОСТИЧЕСКАЯ КЛИНИКА»

1.1 Вербальное описание системы

1.1.1 Морфологическое описание объекта

1.1.2 Функциональное описание объекта

1.1.3 Информационное описание объекта

1.1.4 Историческое описание объекта

1.2 Формализованное описание системы

1.3 Системный анализ проблематики, связанной с системой «Государственная неврологическая лечебно-диагностическая клиника»

1.3.1 Модель анализа проблемы

2 СОДЕРЖАТЕЛЬНАЯ ПОСТАНОВКА ЗАДАЧИ

3 ФОРМАЛЬНАЯ ПОСТАНОВКА ЗАДАЧИ

4 ВЫБОР И ОБОСНОВАНИЕ ИСПОЛЬЗУЕМОГО МЕТОДА

4.1 Алгебраический метод

4.2 Итерационный метод ART

4.3 Критерии, оценивающие точность восстановления

5 ОБЩАЯ ПОСЛЕДОВАТЕЛЬНОСТЬ РЕШЕНИЯ ЗАДАЧИ

6 ОСОБЕННОСТИ ПРОГРАММНОЙ РЕАЛИЗАЦИИ

6.1 Восстановление функции с носителем в круге

6.2 Восстановление функции с носителем в эллипсе

7 РЕЗУЛЬТАТЫ ИССЛЕДОВАНИЙ

7.1 Результаты решения задачи восстановления функции с носителем в круге

7.2 Результаты решения задачи восстановления функции с носителем в эллипсе

7.3 Анализ результатов

ЗАКЛЮЧЕНИЕ

ВВЕДЕНИЕ

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

Хронология развития вычислительной томографии:

г. - открытие рентгеновских лучей;

г. - преобразование Радона;

г. - рентгенограмма в медицине;

г. - линейная томография, вращательная томография;

г. - РВТ в радиоастрономии;

г. - свёрточный алгоритм;

г. - алгоритм РВТ А. Кормака;

г. - серийный томограф Г. Хаунсфилда;

1977 г. - учебный курс по вычислительной томографии в университете штата Нью-Йорк (г. Буффало);

1979 г. - Нобелевская премия А. Кормаку и Г. Хаунсфилду.

Методы реконструктивной томографии начали активно развиваться в 60-х годах прошлого века, хотя математическая база была, по существу, создана И. Радоном в 1917 г. Однако работа И. Радона, опубликованная в трудах Саксонской академии наук, не попала в поле зрения исследователей и была незаслуженно забыта. В настоящее время предложенные Радоном методы широко используются, а его общая формула обращения не потеряла своего теоретического и прикладного значения.

Реконструктивная томография в медицине стала использоваться после публикации работ А. Кормака в 1963-1964 гг.

Первый серийный рентгеновский томограф был разработан Г. Хаунсфилдом в 1972 г. и выпущен фирмой ЕМI (Великобритания). Хотя в современных томографах алгоритмы Кормака и Хаунсфилда не применяются, в 1979 г. им обоим была присуждена Нобелевская премия по физиологии и медицине. В современных коммерческих томографах, как правило, реализованы различные разновидности более эффективного свёрточного метода. Этот метод в 1961 г. разработала группа математиков, программистов и инженеров для обнаружения неисправностей в атомных реакторах. Позднее было обнаружено, что эта группа повторно открыла метод, разработанный И. Радоном.

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

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

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

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

1 СИСТЕМНАЯ МОДЕЛЬ СЛОЖНОЙ ОРГАНИЗАЦИОННОЙ СИСТЕМЫ «НЕВРОЛОГИЧЕСКАЯ ЛЕЧЕБНО-ДИАГНОСТИЧЕСКАЯ КЛИНИКА»

.1       Вербальное описание системы

Цель: получение максимального качества диагноза и лечения заболеваний головного мозга людей.

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

.1.1    Морфологическое описание объекта

Наша система состоит из подсистем:

управляющая подсистема:

− Министерство здравоохранения;

− главный врач;

управляемая подсистема:

− бухгалтерия;

− отдел кадров;

− персонал.

.1.2    Функциональное описание объекта

Метафункция - организация процесса предоставления услуг диагностики с компьютерной томографии, анализа результатов томографии и постановки диагноза.

Функции основных подсистем:

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

Глав. врач (главный врач) - ведущий медицинский работник в медицинском учреждении.

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

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

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

.1.3    Информационное описание объекта

Функции агентства реализуются с помощью таких потоков: финансовые, кадровые, информационные.

Информационные потоки. Из внешней среды в систему:

− Результаты компьютерной томографии;

− Заявки пациентов;

− лицензия.

Кадровые потоки со средой:

− набор необходимых сотрудников.

Финансовые потоки со средой:

− привлечение инвесторов;

− выплата налогов.

1.1.4  Историческое описание объекта

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

В 1963 году американский физик А. Кормак повторно (но отличным от Радона способом) решил задачу томографического восстановления, а в 1969 году английский инженер-физик

Г.Хаунсфилд сконструировал «ЭМИ-сканер» - первый компьютерный рентгеновский томограф, чьи клинические испытания прошли в 1972 году. В 1979 году Кормак и Хаунсфилд «за разработку компьютерной томографии» были удостоены Нобелевской премии по физиологии и медицине.

А в 2003 за изобретение метода магнитно-резонансной томографии, на основе открытия Реймонда Дамадьяна, Нобелевскую премию по физиологии и медицине получили Питер Мэнсфилд и Пол Лотербур.

1.2 Формализованное описание системы

Рисунок 1.1. Модель внешней среды

Рисунок 1.2. Модель «Черный ящик»

Рисунок 1.3. Морфологическая модель

1.3 Системный анализ проблематики, связанной с системой «Государственная неврологическая лечебно-диагностическая клиника»

.3.1 Модель анализа проблемы

В процессе анализа ПС-системы нужно формализовать проблему, выделив ряд неудовлетворённостей в системе.

Эти неудовлетворённости разбиваются на три группы:

а) Нежелательные:

)        Напряженность отношений между коллективом и руководством;

)        Устаревшее оборудование;

)        Текучесть кадров;

)        Неэффективное распределение финансов;

)        Неэффективность распределения рабочего времени;

)        Неэффективное распределение полномочий руководства;

)        Высокая ресурсоемкость;

б) Критические:

)        Низкое качество услуг;

)        Низкий уровень финансирования;

)        Высокий уровень коррупции;

11)     Низкий уровень квалификации сотрудников-управленцев;

в) Желаемые:

)        Низкий уровень квалификации персонала;

)        Низкая заработная плата персонала;

)        Отсутствие системы сбора и накопления данных о результатах деятельности клиники ;

)        Низкий уровень мотивации персонала.

1.3.1.1 Оценка приоритетов неудовлетворенностей

Построим иерархическую структуру этих неудовлетворенностей и применим метод анализа иерархий, используя шкалу Саати.

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

Рисунок 1.4. Матрица парных сравнений и вектор приоритетов для групп неудовлетворенностей

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

Рисунок 1.5. Матрица парных сравнений и вектор приоритетов для желаемых неудовлетворенностей

Рисунок 1.6. Матрица парных сравнений для нежелательных неудовлетворенностей

Рисунок 1.7. Вектор приоритетов для нежелательных неудовлетворенностей

Рисунок 1.8. Матрица парных сравнений и вектор приоритетов для критических неудовлетворенностей

Рисунок 1.9. Иерархическая структура неудовлетворенностей системы

Рисунок 1.10. Глобальный вектор приоритетов неудовлетворенностей

Таблица 1.1 - Глобальные приоритеты неудовлетворенностей

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

Рисунок 1.11. Графическая интерпретация принципа Парето

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

Таблица 1.2 - Нормированные приоритеты наиболее значимых неудовлетворенностей

Таким образом, после применения принципа Парето получим иерархическую структуру неудовлетворенностей, представленную на рисунке 1.12

Рисунок 1.12. Иерархическая структура неудовлетворенностей после применения принципа Парето

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

Рисунок 1.13. Исходный направленный граф

Таблица 1.3 - Матрица достижимости

Рисунок 1.14. Иерархическая структура графа неудовлетворенностей

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

1.3.1.2 Оценка требуемого состояния

Выберем показатели, которые позволят оценить состояние ПС-системы:

) качество услуг;

) финансовые затраты.

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

Основываясь на этих показателях можно сформулировать такие логически возможные сценарии:

1)  Сц.1 - низкое качество услуг (вес 0,185);

2)      Сц.2 - повышение финансирования (вес 0,173);

)        Сц.3 - приемлемое качество и стоимость (вес 0,642).

Таблица 1.4 - Калибровочная таблица

Рисунок 1.15. Иерархическая структура сценариев решения проблемы

Рисунок 1.16. Вектор приоритетов сценариев решения проблемы

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

Рисунок 1.17. Иерархическая модель прямого анализа проблемы

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

Рисунок 1.18. Иерархическая модель оценки вектора приоритетов по критерию «выгода»

Рисунок 1.19. Вектор приоритетов по критерию «выгода»

Рисунок 1.20. Иерархическая модель оценки вектора приоритетов по критерию «ущерб»

Рисунок 1.21. Вектор приоритетов по критерию «ущерб»

Таблица 1.5 - Оценка вектора приоритетов по критерию «выгода/ущерб»

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

В данном случае это сценарий внедрение новых технологий в производство.

2 СОДЕРЖАТЕЛЬНАЯ ПОСТАНОВКА ЗАДАЧИ

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

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

.1 Физическая модель

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

Рисунок 2.1

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

Исходное соотношение:

где:

–          I=I(x) - интенсивность излучения;

–          f (x)∆x - количество лучей, поглощенных на промежутке ∆x;

здесь:

–          ∆І - приращение интенсивности;

–          І0 - начальная интенсивность луча;

–          І0 = I(d) - интенсивность после прохождения луча;

–          - относительное уменьшение интенсивности излучения.

Перейдем от (2.1) к дифференциальной форме:

Уравнение (2.2) - это стационарное уравнение переноса излучения в чисто поглощающей неоднородной среде, описывающее процесс ослабления излучения в веществе, и представляет собой баланс частиц или энергии.

Проинтегрируем левую и правую части (2.2):

Обозначив , получаем:

Согласно постановке задачи I0 и I1 заданные величины, значит правая часть соотношений (2.3), (2.4) - также задана.

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

 (2.2)

Далее будем рассматривать данный интеграл для функции от двух переменных:

 (2.3)

Задача состоит в том, чтобы по совокупности проекций gk восстановить функцию f (x,y) с использованием интерлинации функциц.

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

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

Процесс решения: используется алгебраический метод восстановления с использованием итерационного метода ART.

Исходными данными нашей задачи являются:

проекционные данные gk, полученные с использованием тестовых задач;

значение числа разбивки по осям для дискретизации области;

число направлений и число сечений на перпендикулярных направлениях для формирования информации о различных направлениях сканирования;

параметр релаксации щ.

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

Область применения:

в медицине;

в программном обеспечении компьютерных томографов.

Требуемые ресурсы: наличие ЭВМ, обладающих большой памятью и высокой скоростью.

3 ФОРМАЛЬНАЯ ПОСТАНОВКА ЗАДАЧИ

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

Рисунок 3.1. Круговая геометрия измерений с параллельными проекциями

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

Такие измерения повторяются, пока система "Источник-Детектор" не повернётся на угол . Строго говоря, достаточно поворота на угол , т. к. затем результаты измерений станут повторяться. По полученному набору одномерных проекций необходимо восстановить двумерное распределение f(х,у). Поскольку система "Источник-Детектор" вращается вокруг объекта, такую схему измерений называют круговой геометрией измерений, а т. к. для получения следующего отсчёта в проекции пучок смещается параллельно предыдущему положению, проекции называют параллельными проекциями.

Рисунок 3.2. Неподвижная (х, у) и вращающаяся (х’, у’) системы координат

Для математического описания связи проекций с искомым распределением f(х,у) наряду с неподвижной системой j координат (x, у) введём вращающуюся систему координат (х’, у’) (рис. 1.3):

 (3.1)

Обозначим через (x’,y’) распределение линейного коэффициента ослабления в системе координат (x’,y’), повёрнутой относительно неподвижной системы координат (х,у) на угол :

  (3.2)

В частности (x’,y’) = f(x’,y’)

Для интенсивности I(x’,y’) излучения, прошедшего через объект вдоль

оси y’, в соответствии с (1.2) получим

  (3.3)

При этом было учтено, что за пределами объекта f(x,y)≡0 проекцией R(x’,y’) называют такую величину

 (3.4)

Таким образом, получим следующее выражение для проекции:

  (3.5)

Соотношение (3.5) называется преобразованием Радона двумерной функции f(x,y). Отметим, что если функция  находится по , то это прямое преобразование Радона, а если  по , то обратное.

Если  отлична от нуля внутри круга радиуса r, то (3.5) можно записать в виде

   (3.6)

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

 (3.7)

Происхождение представления (3.8) становится понятным, если вспомнить, что д-функция будет равна нулю везде, кроме прямой y’= x cos() + y sin().

Получено выражение для проекции с помощью определенного интеграла:

 , ,  (3.8)

где

,

gk - заданное число.

Требуется: восстановить функцию f(x,y).

Отметим основные этапы решения задачи восстановления функции по её проекциям.

) преобразование криволинейного интеграла к определенному, т.е. представим (3.8) в виде

  (3.9)

) приближенное решение, полученное на основе интерлинации, имеет вид:

 (3.10)

где

Получили вид приближенного решения, в котором общее количество неизвестных .

4)      подставив в (3.9) получим систему линейных алгебраических уравнений имеем:

 (3.11)

) решение системы отыскиваем с помощью итерационного алгебраического метода ART:

 (3.12)

4 ВЫБОР И ОБОСНОВАНИЕ ИСПОЛЬЗУЕМОГО МЕТОДА

4.1     Алгебраический метод

Пусть функция f(x) = f(x, y) описывает некоторое распределение плотностей в каком-либо выделенном сечении объекта. Основная задача вычислительной томографии состоит в восстановлении функции f(x) по набору экспериментально полученных проекций:

 (4.1)

которые представляют собой линейные интегралы от искомого распределения вдоль прямых L:. Здесь  - угол сканирования,  - дельта-функция.

На практике, как правило, проекции  заданы не для всех значений  и , а только для конечного их числа. Существует целый ряд практических задач, для которых число дискретизаций  по 0 весьма ограничено (от 3 до 5). Задачи такого типа относятся к задачам малоракурсной томографии и являются одними из наиболее трудно решаемых. Задача может быть поставлена следующим образом: по заданному конечному набору проекций  функции двух переменных  получить наилучшую оценку этой функции.

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

Пусть требуется восстановить двумерную функцию f(x)=f(x,y), заданную в области D Ì R2. Предположим, что область восстановления D заключена в квадрат К, который разбит на п равных маленьких квадратиков, называемых элизами. Пронумеруем все элизы от 1 до п. При этом примем основное ограничение, которое заключается в том, что восстанавливаемая функция f(x) принимает постоянное значение fj внутри j-го элиза, т. е. функцию f (x) заменяем дискретизированным выражением

  (4.2)

где

если (х) Î j-му элизу;

в противном случае. (4.3)

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

 (4.4)

Тогда  - проекция функции f(х) вдоль луча Li.

Применяя операторы  к равенству (4.2) и учитывая их непрерывность и линейность, получаем систему линейных алгебраических уравнений

где , i = 1, ..., m; j = 1, ..., n.

Если семейство базисных функций {bj} задается формулой (4.3), то

длина пересечения i-го луча с j-м элизом.

Матрицу коэффициентов обозначим А=(), вектор изображений - f=(f1, f2, ..., fn), вектор проекций - R=(R1, R1, ..., Rт). Тогда решение задачи сводится к решению системы линейных алгебраических уравнений вида

Af = R.   (4.5)

При этом вектор R задан заведомо с некоторой погрешностью.

Стоит отметить, что вид системы (4.5) зависит от конкретного выбора системы базисных функций bi и набора функционалов Ri. Существуют другие способы выбора сетки разбиения области D (а значит, и базисных функций bi). Функционалы  выбираются не только в виде (4.4), но и с учетом реальной длины лучей и с использованием кусочно-постоянных функций. Кроме того, постановка задачи не зависит от геометрии лучей и легко формулируется для трехмерного случая.

4.2     Использование операторов интерлинации

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

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

 ;

разбиения Е2 на четырехугольники. Введем следующие обозначения.

Оператор О1 является оператором аппроксимации f(x,y) кусочно-постоянными функциями по x. Если y=const, то  находится из условия наилучшей аппроксимации f(x,y) в полосе , yÎE. Аналогично, оператор О2 является оператором аппроксимации f(x,y) кусочно-постоянными функциями по y.

Если x=const, тогда hj(x) находится из условия наилучшей аппроксимации f(x,y) в полосе , хÎE.

Введем следующие операторы:

Значения  найдем из условия наилучшей аппроксимации f числом f(оij,h ij) в

Лемма 3.1 Пусть функция, r=1,2 или  и является функцией с ограниченной вариацией. Тогда операторы Onm обладают свойствами

,

где

Доказательство. Свойства (3.25) и (3,26) вытекают из того, что

Поэтому

Свойство (3,27) вытекает из того, что

Свойства (3,29) выполняются для всех дифференцируемых функций и для непрерывных функций с ограниченной вариацией.

Лемма 1 доказана.

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

Следствие 2. Заменяя функции  кусочно-постоянными функциями одной переменной с той же самой оценкой погрешности

получим оператор

Получим значения для gi (x)

Получим значения для Gi (y)

со следующими свойствами:

Следствие 3. Оператор

имеет следующие свойства:

Если , r=1,2 или  и является функцией с ограниченной вариацией, тогда

Доказательство. Для погрешности  можно написать равенство

Отсюда вытекает неравенство

Применяя оценки 3 и 4 к правой части полученного выражения, придем к оценке (3,42).

Следствие 3 доказано.

Если m=n, тогда оператор  имеет погрешность  (он использует  постоянных ); приближение оператором  имеет погрешность . То есть оператор  (он использует  постоянных ) имеет ту же погрешность, как и оператор :

В следующих пунктах отмечаются преимущества указанного метода.

Количество неизвестных

Использование интерлинации функций при построении приближенного решения, а именно представление приближенного решения  в виде:

привело к появлению 2n3+n2 постоянных , которые являются неизвестными. Следовательно оператор  использует O(n3) постоянных-неизвестных. Оператор  имеет погрешность .

Использование оператора  - классическое представление приближенного решения - приводит к появлению n4 постоянных , которые являются неизвестными. Следовательно оператор  использует O(n4) постоянных-неизвестных. Оператор  имеет погрешность .

Обобщая сказанное, делаем вывод, что использование оператора  требует нахождение O(n3) неизвестных, в то время как использование оператора  требует нахождения O(n4) неизвестных для приближения решения с той же самой погрешностью .

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

Для иллюстрации указанного факта приводим следующую таблицу:

Таблица 1

n                 - неизвестных

Сравнения показывают, что для достижения одной и той же точности, при использовании оператора , можно брать меньшее количество уравнений. Например, для n=9 количество неизвестных в классическом методе в 4 раза больше.

В силу того, что система должна быть переопределенной, а для n=9 неизвестных 1539 (для случая с интерлинацией) и 6561(для классического метода), и следует брать число уравнений больше, чем число неизвестных, то ясно, что в методе с интерлинацией этих уравнений будет меньше.

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

Дискретизация области

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

Для  - нерегулярная сетка: разбивка на квадраты со стороной  и прямоугольники со сторонами ,  и , , вытянутые вдоль оси Ox и Oy соответственно. Узлы сетки располагаются в центрах квадратов и прямоугольников.

Для  - регулярная сетка: разбивка на квадраты со стороной . Узлы сетки располагаются в центрах квадратов.

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

,

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

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

.

Относительно точного совпадения в указанных центрах . Значит,

4.2     Итерационный метод ART

реконструкция алгебраический функция восстановление

Среди специфически томографических методов одним из наиболее популярных является метод ART (algebraic reconstruction technique), для которого разработано много модификаций. В первоначальном варианте метода ART вычисления проводят по формуле

    (4.6)

где aiТ = [(ai)1, (ai)2, ..., (ai)n] - i-я строка матрицы В, записанная в виде вектора (столбца) (ai)j≡aij, а  - релаксационный множитель, находящийся в диапазоне 0 <  < 2.

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

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

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

  (4.7)

В полностью ограниченном варианте дополнительно учитывают и информацию о максимальном значении Q, которое не может превысить ни одна из компонент искомого вектора:

(4.8)

Геометрически схема метода ART соответствует последовательному ортогональному проецированию на гиперплоскости, определяемые каждым из уравнений решаемой СЛАУ.

4.3 Критерии, оценивающие точность восстановления

В качестве меры близости восстановленного изображения  к истинному  используются такие критерии:

) среднеквадратическая погрешность:

 (4.9)

) абсолютная погрешность:

  (4.10)

) максимальная по модулю погрешность:

 (4.11)

5 ОБЩАЯ ПОСЛЕДОВАТЕЛЬНОСТЬ РЕШЕНИЯ ЗАДАЧИ

Будем рассматривать проекционные данные по разным направлениям.

Дана система интегральных уравнений:

 , .  (5.1)

Задание правой части gk осуществляется с использованием тестовых задач, где функция f(x,y) задана и gk подсчитана по следующей формуле:

.  (5.2)

Применение интерлинации функций обусловило дискретизацию области с нерегулярной сеткой.

Проводим дискретизацию области, разделив её на n2 квадратов, n3 прямоугольников первого типа («вертикальных») и n3 прямоугольников второго типа («горизонтальных»). Общее количество элизов в данном случае n2+2n3.

Рисунок 5.1. Квадрат с11

Рисунок 5.1. Дискретизация области

Тут aij - центры «горизонтальных» прямоугольников; bij - центры «вертикальных» прямоугольников; сij - центры квадратов.

Используем вид решения, построенный на основе операторов интерлинации:

  (5.3)

Учтем, что

.   (5.4)

Подставив (5.3) в (5.4) получим

 (5.5)

Имеем систему уравнений для определения неизвестных коэффициентов Aiµ, Bнj, Cij.

Количество неизвестных определяется дискретизацией области и равно n3+2n2.

-значения в центрах элизов, Aiµ - значения в центрах «горизонтальных» прямоугольников, Bнj - значения в центрах «вертикальных» прямоугольников.

Получаем систему:

. (5.7)

Эта система неопределенна (недоопределена, если , переопределена, если ).

Её решение осуществляется таким вариантом итерационно-алгебраического метода:

 (5.8)

где  - i-я строка матрицы А, записанная в виде вектора (столбца), а щ - релаксационный множитель, находящийся в диапазоне .

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

 (5.9)

Оценку погрешности данного метода будем производить по такому критерию:

 (5.10)

6 ОСОБЕННОСТИ ПРОГРАММНОЙ РЕАЛИЗАЦИИ

Реализацию данной задачи производили в пакете Mathcad. Благодаря этому удалось реализовать метод в матричном виде. При реализации область сканирования размещаем в круге с радиусом , что определяет использование формулы (3.5) в алгоритме.

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

.1 Восстановление функции с носителем в круге

Пусть задано уравнение окружности радиуса r=0.5, помещенной в квадрат [-1, 1]  [-1, 1]:

В пакете MathCad данная задача будет иметь вид

Далее данную тестовую функцию присвоим функции :

.

Рисунок 6.1. График и линии уровня искомой функции f(x, y)

= 10 - дискретизация области: разбиение области [-1, 1]  [-1, 1] сеткой на n2 элизов.= 32 - число прямых вдоль каждого направления.= 32 - число направлений сканирования.

∆ =  - шаг сетки.

Функция, задающая направления имеет вид:

, .

Разбиения вдоль направлений:

, .

Центры элизов:

, .

Матрица, представляющая значения функции в центрах элизов:

.

Внутри каждого элиза значение искомой функции постоянно.

Сформируем правую часть системы уравнений - вектор B:

Зададим функцию-индикатор, принимающую значение 1, если точка (x, y) принадлежит элизу , и равна 0 в противном случае:

Сформируем проекционную матрицу А размером MN:

Полученная система линейных уравнений AX=B является неопределенной (в данном случае переопределенной), так как число переменных не равно числу уравнений системы. Решение данной системы - вектор размером , представляющий значения восстановленной функции в центрах элизов.

Для решения этой системы используем итерационно-алгебраический метод.

 - релаксационный множитель, находящийся в диапазоне 0 << 2.

NM = N· M - количество строк матрицы А.

h = cols(A) - количество столбцов матрицы А.

 - i-я строка матрицы А, записанная в виде столбца.h = 0 - нулевой вектор-столбец размера h.= 10 - количество итераций.

Получение решения задачи в виде вектора итерационно-алгебраическим методом ART:

= FF(J,щ).

Процедура преобразования вектора в матрицу:

Погрешность данного метода ищем в виде:

.2 Восстановление функции с носителем в эллипсе

Пусть задано уравнение эллипса полуосями a=0.5, b=0.3, помещенного в квадрат [-1, 1]  [-1, 1]:

Рисунок 6.2. График и линии уровня искомой функции f(x, y)

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

7 РЕЗУЛЬТАТЫ ИССЛЕДОВАНИЙ

.1 Результаты решения задачи восстановления функции с носителем в круге

Рассмотрим результаты решения задачи восстановления функции с носите лем в круге:

при

Таблица 7.1 - Результаты исследований при щ= 0.05, J=10

Рисунок 7.1. Истинный и восстановленный графики функции

Рисунок 7.2. Истинные и восстановленные линии уровня функции

.2 Результаты решения задачи восстановления функции с носителем в эллипсе.

Рассмотрим результаты решения задачи восстановления функции с носителем в эллипсе:

при , .

Таблица 7.2 - Результаты исследований при щ= 0.05, J=10

Рассмотрим результат восстановления функции с носителем в эллипсе при n=20, N=40, M=40, щ= 0.05 и J=10.

Рисунок 7.3. Истинный и восстановленный графики функции

Рисунок 7.4. Истинные и восстановленные линии уровня функции

.3 Анализ результатов

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

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

Однако при больших значениях n, N, M и J время вычислений сильно возрастает. Для получения результатов, изображенных на рисунках, достаточно мощному компьютеру потребовалось около 10 минут.

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

7. БЕЗОПАСНОСТЬ ЖИЗНИ И ДЕЯТЕЛЬНОСТИ ЧЕЛОВЕКА

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

.1 Анализ условий труда

Рассматривается помещение лаборатории, которое располагается на первом этаже четырехэтажного кирпичного здания. Размеры помещения: ширина 6м; длина 8м; высота 3,2м; площадь 48м2; объем 153,6м3.

В лаборатории постоянно работает 4 человека: системный администратор, дизайнер с ПК и два программиста. На каждого человека приходится по 12м2 площади и 38,4 м3 объема, что соответствует нормам (не менее 6м2 и 20м3) согласно ДНАОП 0.00-1.31-99.

Помещение оборудовано шестью компьютерами типа IBM PC, Pentium 2.1ГГц 512 Мб ОЗУ.

Для анализа условий труда используем систему «Человек-Чашина-Чреда» с точки зрения задач охраны труда.

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

Поделим элемент «Человек» на три функциональные части:

Ч1 - человек, выполняющий целенаправленные действия (программирование, системное администрирование и т.п.);

Ч2 - это человек, оказывающий непосредственное влияние на окружающую среду;

Ч3 - психофизиологическое состояние человека (настроение, удовлетворенность работой, психическая и физическая нагрузка).

Элемент «машина» можно также разделить на три части:

М1 - элемент, который выполняет основную техническую функцию (компиляция программы);

М2 - элемент функции аварийной защиты (зануление, заземление, изоляция, предохранители);

М3 - элемент влияния на окружающую среду (выделение тепла, шума, электромагнитное излучение).

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

Рисунок 1.1 - Схема связей в системе “Ч-М-С”.

Выделим следующие элементы данной системы:

1)      Ч1М1: воздействие человека на управление (мелкие стереотипные движения: ввод данных, нажатие на клавиши);

)        М3Ч3: зрительная информация о состоянии ЭВМ (нагрузка на зрительные анализаторы);

)        М1ПТ: влияние машины на предмет труда: на ней разрабатывается ПО (электромагнитные излучения);

)        ПТМ1: информация о предмете труда (программном продукте), получаемая машиной;

)        СЧ3: влияние среды на состояние организма человека (повышенная температура, шум, микроклимат);

)        Ч2С: влияние человека как биологического объекта на среду (тепловыделение, влаговыделение, потребление кислорода);

)        СМ2: влияние среды на состояние машины (срабатывание аварийных систем)

)        М3С: влияние машины на среду (шум, электромагнитные излучения);

)        ПТЧ3: влияние технологического процесса на человека (эмоциональное и интеллектуальное напряжение).

)        C ↔ СРМ: влияние параметров среды на окружающую среду соседних рабочих мест (повышенная температура, шум, микроклимат и т.п.).

Согласно ГОСТ 12.0.003-74 выделим следующие ОВПФ:

Опасный: повышенное значение напряжения в электрической цепи, замыкания которой может пройти через тело человека; недостаточная освещенность рабочей зоны; опасные факторы пожара.

Вредные факторы представлены в таблице 1.1

Таблица 7.1 - Карта условий труда

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

В лаборатории возможны такие аварийные ситуации, как короткое замыкание, возгорание.

.2 Техника безопасности

Лабораторию, согласно ПУЭ, можно отнести к классу помещений без повышенной опасности поражения человека электрическим током, поскольку это сухое помещение с нормальной температурой воздуха. Питание ПК производится от трехфазной сети переменного тока с глухозаземленной нейтралью напряжением 380/220В и частотой импульсов тока 50Гц.

В соответствии с требованиями ГОСТ 12.1.030-81 для защиты от поражения электрическим током применяется зануление. Для снижения напряжения прикосновения до срабатывания защиты предусмотрено повторное заземление нулевого провода. Визуальный контроль изоляции производится ежедневно, а инструментальное измерение сопротивления изоляции - один раз в год. Для предупреждения воздействие статического электричества ежедневно проводится общее и местное увлажнение воздуха и проветривание помещения. На видном месте установлен аварийный выключатель, полностью выключающий электроснабжение помещения кроме освещения.

Инструктажи по технике безопасности проводятся в соответствии с НПАОП 0.00-4.12-05.

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

.3 Производственная санитария и гигиена труда

В лаборатории возможно присутствие нескольких ОВПФ.

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

Каждое рабочее место в лаборатории соответствует требованиям ДНАОП 0.00-1.31-99. Рабочие места расположены относительно световых проемов так, что естественный свет падает с левой стороны. При этом выдержаны следующие расстояния:

от стен со световыми проемами до рабочего места - 1м;

между боковыми поверхностями видеотерминалов - 1,5м;

между тыльной поверхностью одного видеотерминала и экраном - 2,5м;

высота рабочей зоны (столы с видеотерминалами) - 0,75м.

Для искусственного освещения лаборатории используются люминесцентные лампы, с высокой световой отдачей (до 75 лм/Вт и более), продолжительным сроком службы (до 10000 ч), малой яркостью светящейся поверхности и близким к естественному спектральным составом излучаемого света. Для освещения лаборатории применяются двухламповые потолочные светильники типа ЛСП12, с мощностью - 2х40Вт, габаритами 1,2Ч0,02Ч0,005м и типом кривой силы света - Г1. В каждый светильник установлены по две лампы ЛБ 40, световой поток которых - 3120лм.

Поскольку светильники крепятся к потолку, то их высота над полом практически равна высоте помещения =3,2м. Высота рабочей поверхности hр=0.75м.

Расчетным уравнением метода является:

,

где N - число светильников в ряду;

Ен - нормируемая минимальная освещенность, лк;

Кз - коэффициент запаса, учитывающий запыление светильников и износ источников света в процессе эксплуатации;- площадь пола помещения, м2;- коэффициент неравномерности освещения;- число рядов светильников;

Фсв - световой поток, излучаемый светильником, лм.

Т.к. в помещении ведутся работы с дисплеями, и расчетной документацией, то, согласно нормам ДНАОП-0.00-1.31-99, Ен=400лк. Лаборатория является помещением общественного здания, поэтому при условии чистки светильников не реже двух раз в год Кз=1.5. При оптимальном расположении светильников для люминесцентных ламп

.

Коэффициенты отражения светового потока от потолка, стен и пола, исходя из их состояния, соответственно равны: rп=50%, rc=30%, rпола=10% для данного помещения. Определим индекс помещения:

,

где А - длина помещения, В - ширина помещения, h - высота подвеса светильников над рабочей поверхностью.

Коэффициент использования светового потока з зависит от типа светильников, коэффициентов отражения светового, а также индекса помещения i. Согласно приложению 51 справочника [5] з = 0,64. Поток, излучаемый светильником: Фсв=2·3120=6240 лм.

Находим оптимальное расстояние между светильниками:

м,

где л=0.91 выбирается в зависимости от КСС.

Расстояние от рядов до стен рассчитываем по такой формуле:

Располагаем светильники вдоль длинной стороны помещения. Число рядов  определяем из соотношения:

Определяем количество светильников в ряду:

шт.

При использовании светильников освещенность составляет:

 Лк.

Расстояние между светильниками:

 м.

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

Рисунок 7.1. Схема расположения светильников

.4 Пожарная профилактика

Степень огнестойкости здания - II, согласно ДБН В.1.1-7-2002. Помещение находится в здании категории В по взрывопожарной опасности.

Пожар может быть вызван такими причинами: короткое замыкание, плохие контакты в местах соединения проводов, небрежное отношение с огнем. В помещении предусмотрена система автоматической пожарной сигнализации с датчиками, реагирующими на появление дыма ФНП-1 из расчета 2 на 20м2 (в помещении установлено 6 датчиков ФНП-1) (по ГОСТ 12.1.004-91) .Предусматривается установка первичных средств пожаротушения: два углекислотных огнетушителя ОУ-5 (из расчета 1 огнетушитель на 3 ПК) и ящик с песком по ГОСТ 12.4.009-83.

Организационно-технические мероприятия по пожарной безопасности:

включение вопросов пожарной профилактики во все инструктажи по технике безопасности;

запрет курения в помещении;

назначение ответственного за пожарную безопасность;

контроль изоляции и состояние состояния электропроводки;

размещение на видном месте плана эвакуации при пожаре;

ограничение массы и объема горючих веществ (использовать в лаборатории расходных материалов, необходимых на одну рабочую смену; не допускать захламление рабочих мест).

Вынужденная эвакуация при пожаре производится через рабочий выход шириной 1,5 м согласно ДБН В.1.1-7-2002.

Рисунок 7.2.