Левин Г.Г., Вишняков Г.Н. - Оптическая томография (1989) (1032160), страница 7
Текст из файла (страница 7)
для определения значения функции 7 в точке х необходимо знать интегралы по плоскостям, проходящим через точку х, и по бес- конечно близким плоскостям. Это следует непосредственно из формулы (1.27), так как для вычисления по ней функции 7 не- обходимо двукратно продифференцировать все проекции и сумми- ровать значения производных в точке х. В двумерном случае, как следует из (1.23), для определения значения 7 в точке х необходимо знать интегралы по всем пря- мым, так как для вычисления внутреннего интеграла по р требу- ется задание проекции для всех р. Это существенное отличие и определяет алгоритмы вычисления функций )(х) по их проекциям.
12. МЕТОДЫ ВОССТАНОВЛЕНИЯ ДВУМЕРНЫХ ТОМОГРАММ ПО ОДНОМЕРНЫМ ПРОЕКЦИЯМ Рассмотрим методы восстановления томограмм двумерных объектов, описываемых функцией 7(х,у) по их одномерным параллельным проекциям (1.6). 28 1 Г 1'(г, ф) = — — ~ — 2) Г'(гсоэ(ф — ф), р)Ы<р+ 4«' о + Г Г У(Р 6+3(-Р, ~) „„ [р — г соэ (<р — ф))' (!.28) Для ее вычисления необходима следующая последовательность операций: 1) в одном канале каждая проекция «сглаживается», т. е. выполняется операция одномерной свертки Г(р,ф) с функцией р-'. В другом канале проекции остаются без изменения. Далее все операции над проекциями одинаковы в обоих каналах; 2) каж. дая проекция поворачивается в плоскости ху на соответствующий угол ф; 3) все повернутые проекции «растягиваются» в направлении, перпендикулярном оси р, т.
е. реализуется переход от одномерной функции к двумерной Г(гсоэ (ф — ф), Ч~). Эта операция носит название обратного проецирования; 4) в обоих каналах преобразованные проекции суммируются (интегрируются по ф); 5) вычисляется разность изображений, полученных в каждом из каналов. Достоинством данного алгоритма можно назвать то, что на всех этапах вычисления томограммы преобразуемые функции положительные, и поэтому первые четыре операции могут быть выполнены как в когерентной, так и в некогерентной оптической системе. Недостатком данного алгоритма является двухканальность обработки проекций. 1.2.2.
Метод фурье-синтеза Выполним одномерное преобразование Фурье над обеими частями уравнения (1.6): В этом случае направление проецирования определяется одним параметром — углом между осями р и х. Обзор алгоритмов решения уравнения (1.6), которые широко используются в медицинской компьютерной томографии, приведен в [11, 23, 39]. Специфика применения этих алгоритмов в физических исследованиях, в том числе в оптических и голографических, отражена в [22, 37).
В настоящем параграфе рассмотрим алгоритмы, основанные на интегральных преобразованиях, и определим набор тех операций над проекциями, которые необходимы для восстановления томограмм. Реализация указанного класса алгоритмов в оптических системах наиболее перспективна, так как одно из основных достоинств оптических процессоров — простота выполнения инте- . гральных преобразований. 1.2.1. Инверсное преобразование Радона Формула обращения Радона (1.23) может быть преобразована к виду [3) Р(К т)= ] г (Р у)е ~~гзР =-Ц~(х, у)ехр] — !р(хсоз в+ + у и!и ~)] г(л1у, (!.29) В (1.29) фурье-сопряженная к р частотная координата р может изменяться в пределах ( — со, +со), а переменная ф является параметром, который может принимать любые значения из диапазона углов (О, и) или ( — я/2, +я/2).
Этот факт отражен в (1 29) путем записи переменных р и ~р через точку с запятой В полярной системе координат р, у правую часть (1 29) можно записать в виде двумерного преобразования Фурье Р(р, ф). Спектры проекций также можно записать в полярной системе координат в виде следующей функции: Р( ) Р(Р! Ч~) О < 'Р(я Р~~ О (! 30) Р( — р; юр — я), я ~ ср(2я, р)~0. Тогда выражение (1.29) будет соответствовать теореме о центральном слое: Р(р, ч) = Р(р, у), р~О„О~ у < 2я, (1.31) где Р— фурье-образ 1. Отсюда следует, что искомая томограмма ).(х,у) получается простым обратным преобразованием Фурье выражения (1.31): 3»»» У (х, у) =- — ! Р(р» в) ехр ]Бр(хсоз в + у з1п <р)] рЫрс6р = (2„,)г й~~ !»» = — ) ) Р (р; р) ] р ! екр ]гр(х соз у + у з!и ~р)] дрйр.
(1.32) 2яз Таким образом, в методе фурье-синтеза вначале формируется двумерный спектр томограммы из одномерных фурье-образов проекций по полярной сетке отсчетов, а затем выполняется обратное двумерное преобразование Фурье также в полярной системе координат. При численной реализации метода фурье-синтеза используют алгоритмы быстрого преобразования Фурье (БПФ). Для этоп це' ли необходима промежуточная операция интерполяции значений спектров проекций с полярной сетки отсчетов на прямоугольную.
Возможна также следующая модификация метода фурье-синтеза. Из выражений (1.29), (1.31) следует, что значения одномерного фурье-образа любой ф-й проекции совпадают со значениями двумерного фурье-спектра искомой томограммы вдоль прямой рз!п (р — ф) =О, проходящей через начало координат частотной плоскости рф под углом ф к базовой линии Если все спектры проекций в частотной плоскости р~р повернуть на соответствующий угол и накопить, т е. просуммировать, то в результате сформируется двумерный спектр ЗО 5(р, Ф)= ~ Р(р; Ф)6(р з(п(Ф вЂ” Ф))НФ. к)1 Для вычисления интеграла применим следующее свойство б-функций (1): Ь(а(1)) = ~' 3(4 — 4 )/ ( а'(1„) ), а(1„)=0, а'(1) =ба!И, а'(1„) чЬО.
Тогда ~п ' 1 5(р. Ф)= ) Р(р' Ф)! р 1 '3(Ф вЂ” ФИФ = — Р(р; Ф)= !р1 (1.33) 1.2.3. Метод суммирования грильтрованных обратных проекций Подставим в (1.32) выражение (1.29) для фурье-образа проекции и, поменяв порядок интегрирования по р и р, получим 1 С ~(х, у) = ~ ~ у (р, <р) ~ ехр [(р(хсоз р+уз(по — р)1Х Х ~ р1Ф( 'г. (1.37) 31 = — '~(р,.) (1.34) р Из теоремы о центральном слое (1.31) следует, что 5Ь. 'г) = (1/р) РЬ* 1). (1.35) Поэтому искомую томограмму можно восстановить из произведения р5(р, ~р) двумерным обратным преобразованием Фурье: й:с ю г (х у) = — Г ) 5 (р, Ф) р ехр (1 р (х соз <р+ у з(п р)1ьдрсАр.
(1.36) (2")з ' Для восстановления томограммы по формуле (1.36) необходимо выполнить следующие операции: 1) вычислить фурье-образы проекций г" (р; Ф); 2) повернуть полученный одномерный спектр в частотной плоскости р~г на соответствующий угол ф. В двумерной плоскости такой одномерный фурье-спектр можно описать функпней Г(р; Ф)й(рз1п (~р — Ф)); 3) накопить все спектры в плоскости ра. Результат накопления описывается функцией 5(р, Ч~); 4) умножить полученный спектр 5(р, Ч~) на двумерный фильтр р; 6) выполнить двумерное обратное преобразование Фурье После выполнения первых трех операций в фурье-плоскости одномерные спектры будут расположены в виде некоторой радиальной структуры-звезды (см, рис.
В.б,в) Поэтому плотность отсчетов в двумерном спектре по координате р обратно пропорциональна р. Это видно также из выражения (1.34) Отсюда становится ясным физический смысл четвертой операции. Она выравнивает интегральную плотность отсчетов в спектре по всей частотной плоскости. Введем функцию ] в /з(р) = — ~ ! р ! ехр(1рр)е(р.
2н (1.38) Тогда (1.37) примет вид е ~(х, р) = — ~T(хсоз р+рз,.пф, ф),(ф, о (1.39) где ар ру Рнс 13 Фильтру«адан функция (а) н одномерный р-фильтр (о) у(р, ф) = ~ у(д, р)й(д — р) оп=у(р, ф)(рй(р) (140) — фильтрованная проекция. Фильтрующая функция Ь(р) может быть записана в виде обобщенной функции Р(1/р»), т. е.
/т (р) = — (1/2 з) Р [1/Ф!. (1,41) На рис 1.3 схематично изображена фильтрующая функция Ь(р) и ее частотный фильтр !р!. Видно, что, подставив данную функцию в (1.40), из (139) получим формулу инверсного преобразования Радона (1 23) Для восстановления томограмм по формуле (1.39) необходимы следующие операции 1) одномерная фильтрация каждой проекции (1 40) фильтрами (1 38), (1 41), называемая р-фильтрацией проекций; 2) поворот проекций в плоскости хр на соответствующий угол гр, 3) замена координаты р на хсозф+уз!и ф.
Эта операция называется обратным проецированием, так как она в некотором смысле обратна операции получения проекций, т. е. операции прямого проецирования. Каждая одномерная проекция при этом как бы «растягивается» в направлении, перпендикулярном оси р, т е становится двумерной функцией /(хсозф+рз1пф); 4) суммирование всех фильтрованных обратных проекций. Отметим особенность данного алгоритма, которая может оказаться существенной при его реализации в оптической системе Суть ее заключается в том, что при фильтрации проекции фильтром (1.38) значения спектра на нулевых частотах становятся равными нулю, т. е.
г(0, ~р) =О Но, с другой стороны, Г(О, р)= ( У(р, р)ар=О. Следовательно, фильтрованные проекции 7 должны быть биполярными Поэтому при их суммировании отрицательные участки одних фильтрованных обратных проекций компенсируются положительными значениями других обратных проекций, что приводит в конечном счете к восстановлению самой томограммы 1.2.4. Метод фильтрации суммарноео изображения Возможен также другой подход к решению задачи восстановления томограммы, Который впервые был предложен в [8). Применив к обеим частям выражения (1.34) обратное преобразование Фурье, получим Ь« з(х, у) = — ) ) 5(р, ~р) ехр[юр(хсозр+уяприрФройр= (2) оо й ~ — — ) ) Р(р; ср) ехр11р(х соз р+ у 31п р))о(робр = о 1 = — ) у' (х соз р + у зш р, р) аьр.
о Видно, что функция з(х, у) образована путем суммирования по углу обратных проекций 11(х сов р+уз)ну), поэтому ее называют суммарным или интегральным изображением. Учитывая выраже. ние (1 35) и ок 1 (2я)' о о ~ ~ — ехр1)р(х соа р+ у яп р))рара'р = „$~'х2+у2 можно получить следующую связь суммарного изображения с исходной функцией: 1 о""1 з(х, у)= — ~ ( — 7 (р, р)ехр (р(х спз р+у япр)рф1йр = (2п)хо 3' Р 1 =~(х, у)®(3 (1,42) „)/'хз + у2 где 88 — значок двумерной свертки. Таким образом, суммарное изображение з(х,у) представляет собой низкочастотный вариант томограммы, так как фильтр р ' 3 — 1137 33 ослабляет ее выа>кие пространственные частоты, Это приводит к «размытию» мелких деталей на суммарном изображении. Как следует из (1.31), для устранения этого размытия спектр суммарного изображения 5(р,ф) необходимо умножить на двумерный р-фильтр. Восстанавливающее действие этого фильтра основано на усилении высоких пространственных частот суммарного изображения, что приводит к подчеркиванию его мелких деталей.
В пространственной области алгоритму (1.31) соответствует выражение Г (х, у) = з (х, у) ® ® Р ( — 114«з (х«+ уз)зя), так как обратное двумерное преобразование Фурье от р имеет вид обобщенной функции Р[( — 1/4п') (х»+у») — з~»). В заключение отметим, что процесс восстановления томограмм является линейным. Этот факт можно использовать для анализа алгоритмов восстановления томограммы следующим образом. Предположим, что исследуемый объект описывается двумерной б-функцией, поэтому фильтрованная проекция (1.40) совпадает с самим фильтром Ь(р).