Диссертация (1145332), страница 10
Текст из файла (страница 10)
Так, при реконструкции 256 изображений размером 1024×1024 пикселей каждое, при использовании просвечиваний по 40401прямым, время вычислений на упомянутом выше процессоре увеличится более чем в 1000 раз. Таким образом, использование технологии параллельныхвычислений при реконструкции структуры трехмерного объекта представляется весьма уместным.Эффективная параллелизация алгоритмов, решающих задачу компьютерной томографии, может осуществляться для алгебраических методов обращения преобразования Радона [20], а также Фурье-алгоритма [29].Остановимся на обсуждении перспектив применения параллельных вычислений при реализации алгоритма свертки и обратной проекции. Вначалеобсудим перспективы реализации параллельных вычислений с использованием графических процессоров на основе технологии CUDA. Будем осуществлять параллелизацию при реконструкции изображения в каждом сечении.Пусть радиационное просвечивание проводится на основе параллельной схемы сканирования, которая характеризуется тем, что для каждого фиксиро-63ванного направления облучения источники излучения и детекторы располагаются на параллельных прямых (см.
рисунок 1.3). Пусть для фиксирован-Рисунок 1.3 — Параллельная схема сканирования.ного направления число линий, по которым осуществляется радиографическое облучение, равно 2q + 1. Число направлений, по которым осуществляется облучение, равно p. В результате облучения получаем двумерный массивзначений интегралов gj,k , где j есть индекс, характеризующий направлениеоблучения (j = 1, ..., p), а k определяет номер линии для фиксированного направления (k = −q, ..., q). Реконструкцию среды в отдельном сечении будемосуществлять на основе алгоритма свертки и обратной проекции, которыйсостоит из двух шагов.Шаг 1.
Вычисление свертки.Для всех индексов j = 1, ..., p и k = −q, ..., q вычисляем значения свертки:vj,k = hqXwb (h(k − i))gj,i .i=−qЗдесь h есть расстояние между двумя ближайшими параллельными линиями,wb есть фильтр с заданной спектральной шириной b. Например, b можновзять равным π/h, если эта величина меньше p.
В этом случае wb (h(k − i)) =1/π 2 h2 (1 − 4(k − i)2 ).64Шаг 2. Вычисление обратной проекции.Для каждой точки (r1 , r2 , r3 ) при r3 = const вычисляется искомая функция µ(r):p2π X((1 − u)vj,k + uvj,k+1 ),µ(r) =p j=1где u = s/h − k, s = r1 cos γ + r2 sin γ, γ = jπ/p, k ≤ s/h ≤ k + 1.Обсудим некоторые результаты тестирования параллельной версии алгоритма, реализованной на основе технологии CUDA. Расчеты были проведеныс помощью графической карты nVidia GeForce GTX 660M (CUDA Capability3.0).Исходя из того, что каждая свертка рассчитывается независимо от других, вычисление сверток можно выполнять на графических процессорах впараллельных потоках (нитях). В таблице 1.1 приведены результаты тестирования первого шага алгоритма для разных размеров блока (blockDim.x) иразного количества блоков (gridDim.x).
Для оптимизации кода были введенывременные переменные для хранения некоторых промежуточных результатов вычислений, что уменьшило время выполнения в 4,226 раза. Кроме того,для уменьшения вычислений внутри функции все арифметические операции,не зависящие от индексов, были рассчитаны на основном процессоре и помещены в константную память устройства.Если количество просвечиваний, равное 2q + 1, на одном повороте меньше, чем максимальное допустимое количество нитей в блоке, то вычислениесверток для каждого направления можно поместить в отдельный блок.
Приэтом для уменьшения времени работы программы требуется выполнить следующие оптимизационные действия.1) Сделать блоки с данными просвечивания одной проекции кратными128 битам (см. рисунок 1.4).2) Исходя из того, что размер варпа равен 32, определить размер блокакак 32s, где 32(s − 1) < 2q + 1 ≤ 32s.3) При работе блока с очередной проекцией копировать сначала данныепросвечивания в этой проекции из глобальной памяти в разделяемую память.65Таблица 1.1 — Время вычисления операции свертки в зависимости от числанитей/блоков№12345678910111213141516Кол-во нитей в блоке8163264128256512102410241024512256128643216Кол-во блоков111111112448163264128Время работы (сек.)2,331191,304240,7499460,4410250,2571050,1803250,1080490,0549560,02750240,02752490,02752380,02758350,02757930,02756840,05408810,100913Рисунок 1.4 — Схема хранения данных.66Параллельная реализация обратного проецирования производится аналогично первому шагу, за исключением того, что использование разделяемойпамяти здесь уже невозможно, так как расчеты значения в определенномпикселе требуют значений для всех углов просвечивания.Относительно задачи восстановления трехмерной структуры объекта понабору сечений следует отметить, что при последовательном восстановлениипоследовательности 2D изображений, время, затрачиваемое на реконструкцию, может уменьшаться за счет встроенных алгоритмов оптимизации.
Такпри восстановлении 32 сечений фантома Кормака [112] получены результаты,представленные в таблице 1.2. Восстановление изображения в одном сеченииосуществлялось в 160000 пикселях по данным просвечиваний по 10201 прямым.Таблица 1.2 — Время реконструкции в последовательности сечений№ сечения12345678910111213141516Время (сек.)0,06859530,06611150,05457590,03703710,03596080,03174360,02917940,02921950,02905450,02873040,02869610,0287040,02789760,02777890,02771460,0276469№ сечения17181920212223242526272829303132Время (сек.)0,02748420,02684620,0272620,02664650,02615850,02623160,02620620,02568610,02539580,02544550,02543910,02508390,0248420,02539340,02482620,0242855Представленные результаты демонстрируют эффективность параллелизации алгоритма решения задачи томографии на основе технологии CUDA.Учитывая также относительно недорогую стоимость GPU устройств, под-67держивающих технологию CUDA, дальнейшее развитие этого направленияпредставляется перспективным.Рассмотрим теперь иной подход в параллелизации алгоритма реконструкции структуры трехмерной среды.
Будем осуществлять параллелизацию понабору сечений на основе технологии MPI. Предложенный алгоритм заключается в следующем: выделим расчетные процессы и управляющий процесс,не использующий большие вычислительные мощности.Опишем алгоритм работы управляющего процесса.Шаг 1.
Обнуляются счетчики количества срезов и расчетных процессов.Шаг 2. Если значение счетчика количества срезов меньше либо равнообщего количества срезов, то оно инкрементируется.Шаг 3. Ожидается запрос от одного из расчетных процессов. В ответ емуотправляется значение счетчика количества срезов.Шаг 4. Если значение счетчика количества срезов больше общего количества срезов, то инкрементируется значение счетчика расчетных процессов.Шаг 5.
Если значение счетчика расчетных процессов равно общему количеству расчетных процессов, то выводится общее время работы расчетныхпроцессов, иначе осуществляется переход к шагу 2.Алгоритм работы расчетного процесса.Шаг 1. Отправляется запрос управляющему процессу, в ответ приходитзначение счетчика количества срезов.Шаг 2. Если значение счетчика количества срезов больше общего количества срезов, то работа процесса завершается, иначе осуществляется переходк шагу 3.Шаг 3.
С помощью алгоритма свертки и обратной проекции восстанавливаются данные среза с порядковым номером равным значению счетчикаколичества срезов. Осуществляется переход к шагу 1.Данный алгоритм позволяет динамически распределить нагрузку между процессами в зависимости от их быстродействия. В качестве исходныхданных использовались зашумленные результаты сканирования 256 сеченийфантома Кормака, состоящих из 101 × 101 измерений каждое, находящихсяв 256 файлах соответственно.
Вычисления были осуществлены на класте-68ре Дальневосточного федерального университета на основе технологии MPIпри использовании 32-х двухядерных процессоров 2 × Intel(R) Xeon(R) CPUE5345 @ 2.33GHz. В таблице 1.3 приведены результаты работы программы.Здесь помимо экспериментальных данных приведен расчет доли параллельной части кода программы, основанный на законе Амдала:Sp =1,1 − α + α/pгде Sp – ускорение, получаемое при использовании p процессоров, α – доляпараллельных вычислений в программе.Таблица 1.3 — Эффективность параллелизации алгоритмаКол-во ядер1248163264Время работы (сек.)482,135241,738121,00562,109431,412617,22159,94997Ускорение11,994453,984427,7626715,348627,996148,4559α10,99720,99870,99560,99720,99540,9949На основе полученных экспериментальных данных можно сделать вывод,что данная параллельная реализация алгоритма реконструкции структурытрехмерной среды достаточно хорошо масштабируется.
Это позволяет существенно уменьшить время обработки результатов сканирования, что являетсявесьма актуальным, например, в области медицинской томографии.Основные результаты и выводыНесмотря большое количество работ по исследованию краевых задач дляуравнения переноса с энергетической зависимостью, ряд вопросов представляет интерес для дальнейших исследований. К задачам, представляющимпрактический интерес, можно отнести задачи компьютерной томографии,формулируемые как обратные задачи для уравнения переноса излучения.69Одной из основных характеристик среды, дающих достаточно полное представление о ее строении, является коэффициент ослабления (коэффициентполного взаимодействия). Поэтому, наибольший теоретический и практический интерес связан с определением именно этой характеристики.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
