CC BY
161
противном случае необходим пересчет. Разобьем отрезок [х р ,х р+1 ] фазовой траектории пополам. Можно записать систему уравнений
u(tр+і) = m-1БТК(1 р )Rx(tр),
u(t р+і) = m-1Б ТКа р )Rx(t рр), (10)
где в левой части величина управлящего воздействия соответствует области насыщения.
Подставив в систему уравнений (10) желаемые значения фазовых координат и решив ее, получим значения
параметров —, і є [1,2] критерия качества (3). Проделав аналогичную процедуру для каждого отрезка разбиения фазовой траектории, получим закон изме-
J*i
нения параметров —, і є [1,2] критерия качества (3)
в процессе отработки синтезированной СУ требуемого задания.
7. Принцип реализации синтезированного алгоритма управления
Пусть хж(to) = x(to), хж(t1) = x(t1). Тогда получим хж (t) = eAiKtxж (t0) или xж (t0) = (eA t)-1xж (t) .
Поэтому при x(t) = x ж (t) для открытой области получим L(t) = -m^Б^Д^е^)-1х ,
■t- = (A - m-1BBTK(t)R(eAiKt)-1)x .
8. Заключение
В классе задач АКОР исследована процедура структурного и параметрического синтеза гармонического осциллятора с демпфированием как ОУ. Проведенное исследование для рассматриваемого класса задач позволило получить следующие новые результаты, имеющие научное и прикладное значение:
- в рамках единого подхода получено решение задачи структурного синтеза для двух возможных состояний (устойчивое и неустойчивое) рассматриваемого ОУ;
- реализовано аналитическое решение задачи структурного синтеза, что позволило определить размерности управляющих параметров и установить их связь с параметрами рассматриваемого ОУ;
- процедура параметрического синтеза реализуется на основе пакета прикладных программ “линейная алгебра”.
Практическая значимость результатов исследования определяется возможностью их использования в качестве математического обеспечения при проектировании СУ для рассматриваемого класса задач динамического синтеза.
Литература: 1. Мандельштам Л.И. Лекции по теории колебаний. М.: Наука, 1978. 470с. 2. Лойцянский Л.Г, Лурье А.И. Курс теоретической механики, т.2. М.: Госте-хиздат, 1954. 595с. 3. Моисеев Н.Н. Асимптотические методы в нелинейной механике. М.: Наука, 1969. 380с. 4. Атанс М., Фалб П Оптимальное управление. М.: Машиностроение, 1968.764с. 5. Андронов А.А., Витт А.А., Хай-кин В.Э. Теория колебаний. М.: Наука, 1981.568с. 6. Магнус К. Колебания. Введение в исследование колебательных систем. М.: Мир, 1982.303с. 7. Божко А.Е. Синтез оптимального управления колебательными системами. Киев: Наук. думка, 1990. 164с. 8. Радиевский А.Е. Динамический синтез в общей проблематике современной теории управления // Радиоэлектроника и информатика. 2004. №3. С.70-75. 9. Радиевский А.Е. Функциональноаналитический метод синтеза детерминированного регулятора // Оптимизация электромеханических систем с упругими элементами. Харьков: ИМиС, 1995. С.137-148. 10. Попов Е.П. Прикладная теория процессов управления в нелинейных системах. М.: Наука, 1973.584с. 11. Летов А.М. Динамика полета и управление. М.:Наука, 1969.360с. 12. Бесекерский В.А.,Попов Е.П. Теория систем автоматического регулирования. М.: Наука, 1972. 768с. 13. Моисеев Н.Н. Математические задачи системного анализа. М.: Наука, 1981.488с.
Поступила в редколлегию 14.06.2006
Рецензент: д-р техн.наук, проф. Кузнецов Б.И.
Радиевский Анатолий Евгеньевич, канд.техн.наук, с.н.с., заведующий лабораторией Харьковского НИИ комплексной автоматизации Минтопэнерго Украины. Научные интересы: математическая теория экстремальных задач, динамические задачи многокритериальной оптимизации. Адрес: Украина, 61003, Харьков, пер. Кузнечный, 2, тел. 731-35-67, 731-41-80.
УДК517.373:517.443:517.444:519.6
ОБ ОДНОМ МЕТОДЕ ВЫЧИСЛЕНИЯ 3D ДИСКРЕТНОГО ОБРАТНОГО ПРЕОБРАЗОВАНИЯ РАДОНА
ТЕВЯШЕВ А.Д., СМИРНОВА В.С.________________
Описывается численный метод нахождения обратного 3D дискретного преобразования Радона. Проводится оценка погрешности вычисления по данному методу. Разработанные и описанные методы и результаты могут быть использованы в реконструктивной томографии для восстановления объекта в пространстве R3 из некоторого набора его проекций.
1. Введение
В данной статье рассматривается один из методов нахождения 3D дискретного обратного преобразования Радона, который основан на использовании обратного псевдополярного преобразования Фурье.
Целью исследования является разработка метода нахождения 3D дискретного обратного преобразования Радона, который обеспечивает высокую точность при больших размерностях изображений.
Для достижения цели исследов ания были поставлены и реализованы следующие задачи:
1) изучение одного из численных методов нахождения 3D дискретного обратного преобразования Радона;
РИ, 2006, № 3
37
2) реализация метода 3D дискретного обратного преобразования Радона на ЭВМ;
3) локализация мест наибольшего накопления погрешности рассмотренного метода;
4) устранение погрешности путем модификации метода;
5) проведение сравнительного анализа погрешностей исходного и разработанного метода.
2. Определение 3D непрерывного преобразования Радона
Рассмотрим кусочно-непрерывную и в общем случае нелинейную функцию f (x) = f (x, y, z), определенную
3
в вещественном пространстве R , и плоскость, которая задана вектором нормали а и расстоянием s до начала координат. Тогда 3D непрерывное преобразование Радона функции f(X) для этой плоскости определяется выражением [1]:
ТО ТО ТО
fa ,s) = J J J f(x )S(xT а - s)dx =
—TO —TO —TO
ТОТОТО
= J J Jf(x,y,z) -x , (1)
—ТО —ТО —ТО
x S(x sin 9 cos ф + ysin 9 sin ф + zcos 9 — s)dxdydz
где x = [x,y,z]T, a = [sin9cosф^ш9sinф,cos9]T, 5 -
дельта-функция Дирака, которая удовлетворяет условиям:
5(x) = 0, x Ф 0,
ТО
J 5(x)dx = 1.
—ТО
Каждой паре (a,s) соответствует плоскость в про-
3
странстве R , а непрерывное преобразование Радона
r 3
функции f(x) из R ставит в соответствие множе-
3
ству всевозможных плоскостей в R множество ее интегралов по этим плоскостям.
В общем случае, для того, чтобы интеграл (1) сходился для почти всех значений а и s, достаточно потребовать, чтобы функция f(x) была абсолютно интегри-
3
руема по всему пространству R [2].
3. Связь 3D непрерывного преобразования Радона с преобразованием Фурье
1D непрерывное преобразование Фурье g(ra) скалярной функции g(x) определяется как интеграл [3]:
ТО
g(ra) = J g(x)e—2ni“xdx.
—ТО
Соответственно, обратное преобразование Фурье задается следующим образом:
ТО
g(x) = 2- Jg(ra)e2ni“xdra 2п _
—ТО
Обозначим 1D преобразование Фурье проекции Ж (a, s) по параметру s через 5Rf (а, £):
ТО
5Ща,О = Jf а,s)e—2ni^ds .
—ТО
Теорема о центральном сечении устанавливает соотношение между 3D преобразованием Фурье функции
f (x) и 1D преобразованием Фурье ее преобразования Радона. Согласно этой теореме [2]:
5Rf(cx, |) = f(^5t) = f(| sin 9 cos ф, \ sin 9 sin ф, \ cos 9) ,(2)
где
то то то T
f(^) = J J Jf(x)e—2ni^(x ■a)dX =
—TO —TO —TO
TO TO TO
= J J Jf(x,y,z)e
—2ni^(xsin 9 cos ф+ysin 9 sin ф+zcos 9)
dxdydz.
—TO —TO —TO
Соотношение (2) позволяет вычислить непрерывное преобразование Радона, используя 3D преобразование Фурье.
4. 3D псевдополярная сетка
3D псевдополярной сеткой называется дискретное множество “ pp , которое в декартовой прямоугольной системе координат (x, y, z) задается следующим образом [2]:
12 3
“ pp = “ pp ^ “ pp ^ “ pp :
где
~1 3vu 3wu . .
“pp = !<3u---ЦІ--^T>: u-v-w 6 {-V2-
2 ,, 3vu „ 3wu,
“pp ={(—u>v>w 6 {-n/2> “pp = {(— 37T,3u): u,v,w є {-n/2,
,n/2>>;
,n/2>>,
,n/2>>.
pp " n/2’ n/2
Множества “pp , “pp , “pp называются секторами псевдополярной сетки “ pp или псевдополярными секторами.
Геометрически множество “pp подобно «песочным
часам» с центральной осью Ox и высотой 3n+1 (рис. 2
1). Аналогично, “ pp лежит вдоль оси Oy (рис. 2),
3
а “pp - вдоль оси Oz (рис. 3). Объединившись в “ pp, они образуют куб с центром в точке начала координат и ребром, равным 3n +1.
Тройка величин (u, v, w) называется псевдополярной системой координат. Координаты v и w можно интерпретировать как «псевдоуглы», а координату и - как «псевдорадиус».
38
РИ, 2006, № 3
5. 3D псевдополярное преобразование Фурье
В декартовой системе координат (x,y,z) введем множество Q ch , которое представляет собой сетку размерности n х n х n:
Qch = {(x,y,z): x,y,z є {-n/2,...,n/2}}, где n - четное.
Множество значений функции f (X) = f (x, y, z), заданной на дискретном множестве точек Q ch , будем называть 3D дискретным изображением I(x, y, z).
Определим понятие псевдополярного преобразования Фурье для 3D дискретного изображения I(x, y, z), заданного на множестве точек Qch . С этой целью введем тригонометрический полином I(u,v,w) вида [1]:
n/2-1 n/2-1 n/2-1 „
I(u,v,w) = E E E I(x,y,z)e-2ni(ux+vy+wz)/m
x=-n/2 y=-n/2 z=-n/2
(3)
где m = 3n +1.
Псевдополярным преобразованием Фурье PPI(u,v, w) дискретного 3D изображения I(x, y, z) будем называть преобразование, которое определяется как значение тригонометрического полинома (3) в узлах псевдополярной сетки:
PPI(u,v,w)
PP1I(u,v,w) = I(3u,-
3vu 3wu
n/2 n/2
ч ї, 3vu „ 3wu,
PP2I(u,v,w) = I(-— ,3u,---—),
n/2 n/2
3vu 3wu
PP3I(u,v,w) = ї(--— ,-
n/2 ’ n/2
3u)
(4)
где u,v,w є {-n/2,...,n/2}.
Применяя псевдополярноге преобразование Фурье получаем преобразование Фурье на псевдополярной
сетке Q pp .
6. 3D дискретное преобразование Радона
Пусть X = {x j, j = 0,...,n -1} - последовательность чисел длины n. 1D дискретным преобразованием Фурье конечной длины n для последовательности X называется сумма [3]:
n-1
yk = E xje j=0
-2nikj/n
k = 0,
,n -1,
где Y = {yk,k = 0,...,n-1}.
(5)
В операторной форме выражение (5) принимает вид:
Y = F(X), где F - оператор 1D дискретного преобразования Фурье. Аналогично определяется обратное преобразование Фурье:
xj
1 n^1 2nikj/n
- E yke '
nk=0
j = 0,..., n -1.
В операторной форме это выражение (13) принимает вид:
X = F-1(Y), где F-1 - оператор 1D дискретного обратного преобразования Фурье.
РИ, 2006, № 3
39
Используя псевдополярное преобразование Фурье, можно дать определение дискретного преобразования Радона. Далее будем использовать операторную форму записи.
Пусть PP - оператор 3D псевдополярного преобразования Фурье: PP(I(x,y,z)) = PPI(u,v,w), где
x,y,z є {-n/2,...,n/2}, u,v, w є {-n/2,..., n/2}.
Известно [1], что 3D дискретное преобразование Радона изображения I(x,y,z) может быть получено путем применения 1D преобразования Фурье к его 3D псевдополярному преобразованию Фурье:
R(I(x, y, z)) = F-1 o PP(I(x, y, z)). (6)
Здесь R - оператор 3D дискретного преобразования Радона.
7. 3D дискретное обратное преобразование Радона
Для восстановления изображения I(x, y, z) из его 3D дискретного преобразования Радона (6) необходимо применить к нему обратное преобразование Радона. Для этого нужно найти оператор R -1, обратный оператору R : R-1 o R(I(x,y,z)) = I(x,y,z).
Применив к обеим частям равенства (6) последовательно слева оператор 3D псевдополярного обратного преобразования Фурье PP-1 и оператор 1D дискретного преобразования Фурье F , получим:
PP-1o Fo R(I(x, y, z)) =
= PP -1 o Fo F-1 o PP(I(x, y, z)) = I(x, y, z) .
Сравнивая полученное выражение с (6), имеем
R-1 (I(x, y, z)) = PP-1 o F(I(x, y, z)) .
Следовательно, для восстановления изображения I(x, y, z) по известному его 3D дискретному преобразованию Радона (6) необходимо к правой части равенства (6) последовательно применить 1D преобразование Фурье и 3D псевдополярное обратное преобразование Фурье. Выражение для 1D дискретного преобразования Фурье известно (5), а задача восстановления изображения I(x, y, z) в этом случае сводится к нахождению 3D псевдополярного обратного преобразования Фурье.
В данной статье рассматривается один из возможных методов нахождения обратного 3D псевдо-полярного преобразования Фурье. Этот метод состоит из двух этапов: первый - это отображение значений тригонометрического полинома (3) из одного множества точек на другое (метод «чистки лука»), второй этап -это восстановление исходного изображения по известным значениям тригонометрического полинома (инвертирование десятичных частот).
7.1. Метод «чистки лука» для 3D-случая
Исходными данными для применения метода «чистки лука» является 3D псевдополярное преобразование Фурье (4) дискретного изображения I(x,y,z), т.е. значения тригонометрического полинома (3) в точках псевдополярного множества Qpp . На первом этапе восстановления изображения I(x, y, z) из (4) необходимо найти значения тригонометрического полинома
(3) в точках множества Qch :
Qch = {(3u,3v,3w): u,v,w = - n/2,...,n/2}.
Точки множества Qc^ , аналогично точкам множества Q pp , заполняют куб с центром в начале координат и ребром, равным 3n +1. Однако в отличие от множества Qpp , точки множества Qch распределены внутри куба равномерно. Точки множеств Qch и Q pp совпадают только на поверхности куба и в точке начала координат.
Метод «чистки лука» состоит из n/2 +1 шагов и основан на последовательном нахождении значений тригонометрического полинома (3) в точках сечений данного куба.
На каждом шаге через точки множеств Q ch и Q pp проводится шесть секущих плоскостей, каждая из которых параллельна одной из граней куба. На первом шаге эти плоскости проходят через грани куба. На каждом последующем шаге секущие плоскости параллельно сдвигаются в сторону центра куба на 3 единицы. Формально это описывается следующим образом.
На j-м шаге (j = 1,...,n/2 +1) необходимо найти значения тригонометрического полинома I(u,v,w) (3) в точках множества ©k :
©k =©x и ©-x и ©k и ©-y и ©k и ©-z,
где
©k = {(3k,3y,3z): y,z = -k,...,k},
©-x = {(-3k,3y,3z): y,z = -k,...,k} ,
©y = {(3x,3k,3z): x,z = -k,...,k} ,
©-y = {(3x,-3k,3z): x,z = -k,...,k} ,
©k = {(3x,3y,3k): x,y = -k,...,k} ,
©-z = {(3x,3y,-3k): x,y = -k,...,k} . k = n/2,...,0, k = n/2 - j +1.
Процедура нахождения значений тригонометрического полинома I(u,v,w) в точках множества ©k , k = n/2,...,0, сводится к последовательному нахождению значений I(u,v,w) в точках его шести подмножеств © x, © - x, © y, © - y, © у, © - z.
РИ, 2006, № 3
40
На первом шаге рассматриваем множество @k при
k = n/2, которое состоит из точек, лежащих на поверхности куба. Поскольку на поверхности куба точки множеств Q ch и Q pp совпадают, значения тригонометрического полинома I(u, v, w) в этих точках известны. Поэтому точки множества @k при k = n/2 можно считать восстановленными.
На j -м шаге необходимо найти значения тригонометрического полинома I(u, v, w) в точках множества @k при k = n/2 - j +1. Рассмотрим процесс нахождения
значений I(u,v,w) в точках подмножества @ X множества @k . Множество точек в сечении куба можно рассматривать как 2D изображение (рис. 4).
г
О о оо о оо о оо о оо о о
. g °о О <0 8 °о О 8 •
■Ф ■о О О О О О О О О о
% 8 о° о °о 8 о° о °о 8 г
0 О оо оо о оо о оо о о
. 8 °о о о° 8 °о о о° 8 •
Ф О о о о о о о о о о
% 8 о° о °о 8 о° о °о 8 $
о О оо оо О оо о оо О о
У
Рис. 4. Сечение куба при n = 8, j = 2
В сечении куба условно можно выделить три типа точек. Первый - это точки множества Qch , в которых
значения полинома I(u,v,w) были найдены на предыдущем шаге. Второй тип - это точки множества Q, в которых необходимо найти значения полинома I(u, v,w). Третий тип - это точки псевдо-полярного множества Qpp , в которых значения полинома I(u, v, w) изначально известны.
Рассмотрим точки, выделенные на рис. 4 (рис. 5).
__________6__________
О__О оооооооооооо о
1 а
Рис. 5. Точки из рис.4
Необходимо найти значения полинома I(u,v,w) в серых точках по известным значениям II(u, v, w) в черных и белых точках. Формально эта задача описывается следующим образом.
Обозначим через фррХ множество точек:
фkpx = (№-^ -~^): y,z = -п/2,...,п/2}
n/2 n/2 •
k
Множество фррХ состоит из точек третьего типа, а
значит, значения полинома II(u, v, w) в этих точках известны.
Введем множество у ррХ :
і 3kz
уррх = {(3k,3y,- —): y,z = - п/2,...,п/2} .
Значения полинома I(u,v, w) в точках множества
k
фррх :
I(3k,-?%.-Щ) =
n/2 n/2
n/2-1 n/2-1 n/2-1
= S S S I(u,v,w)e
u=-n/2 v=-n/2 w=-n/2
, 3ky 3kz
-2ni(3ku---y- v--— w)/m
n/2 n/2
или
, . .3ky
3ky 3kz n/2-1 2m—— v/m
I(3k,-%-^ = S ave П2 , (7)
V2 V2 v = - n/2
где
n 2-1 n 2-1
av = S S I(u,v,w)e
u = -n/2 w = -n/2
- 2nl(3ku - ^-^w)/ m n2
Аналогично, значения полинома I(u,v, w) в точках множества у ррХ:
I(3k,3y,-^) = n/2
n 2-1 n 2-1 n 2-1
= S S S I(u,v,w)e
u=-n/2 v=-n/2 w=-n/2
3kz
-2nl(3ku+3yv--—w)/m
n2
n/2-1 2 -3 /
= S ave-2ni3yv/n
(8)
v=-n/2
Выражения (7) и (8) представляют собой две системы линейных алгебраических уравнений:
n/2 -1
fj = S aveivyj
v = -n/2
n 2-1
gl = S a ve
lvXl
v=-n/2
(9)
(10)
3kj
где yj = 2n—j-/m, xi =-2n3l/m, n/2
fj=*(*•-j^)' gi=I(3M--^2)-
j,z = -n/2,...,n/2, l = -k,...,k .
3kz,
РИ, 2006, № 3
41
Из системы (9) находим значения a v и подставляем их в систему (10). В результате находим gj, соответствующие значениям I(u,v,w) на уjppx . Значения gj также могут быть вычислены по приближенной формуле, основанной на интерполяции полиномом Лагранжа [4]:
k 1
gi = ci Z fjdj(—----— - 0
n/2-1
gl = Z aweivxl
v=-n/2
(15)
где
j=-k tg((xl - yjV2)
n/2
cl = Пsin((xl - yj)/2) j=-n/2
k 1
dj p?ksin((yj - ypV2),
p* j
l = -k,...,k, j = - n/2,...,n/2.
(11)
(12)
(13)
В результате получается множество точек, показанное на рис. 6. Аналогично рассматриваются точки, расположенные вдоль оси Oz.
• о О О О О О О •
• 8 8 8 8 8 8 8 •
о о ❖ о о о ❖
• 8 8 8 8 8 8 8 •
• О О О О О О О •
• 8 8 8 8 8 8 8 •
о о ❖ о о о ❖
• 8 8 8 8 8 8 8 •
• О О О О О О О •
Еш-
Рис. 6. Множество точек
3kz4
Далее представим I(3k,3y,-----—) и I(3k,3y,3z) соот-
n/2
ветственно в виде:
, „ . 3kz ,
3kz у 2-1 2ni——w / m
I(3k,3y,- —) = Z a we n/2 ,
n/2 w = -n/2
n/2-1
I(3k,3y,3z) = Z awe2ni3zw/m
w=-n/2 ’
n/2-1 n/2-1
где aw = Z Z I(u,v,w)e-2ni(3ku+3yv)/m .
u=-n/2 v=-n/2
Отсюда имеем две системы линейных алгебраических уравнений:
fj = П/Г a welvyj
v = -n/2
(14)
здесь y: = 2n-^f^/m, xl =-2n3l/m, j n/2 j
fj = I(3k,3y,-І, gl = I(3k,3y,3l),
J n/2
j,y = -n|2,...,n|2, \ = -k,...,k.
Для нахождения gl из системы (14) находим значения aw и подставляем их в систему (15), либо используем приближенные формулы (11)-(13). В результате находим gl , соответствующие значениям I(3k,3y,3l), которые представляют собой I(u,v,w) на ©X. Аналогично находятся значения полинома Fu,v, w) в остальных точках множества ©k .
7.2. Инвертирование десятичных частот
Значения полинома (3) на декартовом множестве Q.ch можно также найти как
n/2-1 n/2-1 У 2-1
I(u,v,w) = Z Z Z I(x,y,z)e-2™3(ux+uy+wz)/m x =-n 2 y=-n 2 z=-n 2
на множестве Q ch = {(x,y,z): x,y,z = - n/2,..., n/2}.
Введем одномерный оператор Fq : Cn ^ Cn+1 такой, что
n 2-1
(FDx)(k) = Z x(u)e-2niu3k/m
u=-n 2
k = -n/2,...,n/2, m = 3n +1.
Для всех
x,y,z = -n/2,...,n/2 - 1,u,v,w = -n/2,...,n/2
значения I могут быть получены следующим образом:
I(x,y,w) = FD(I(x,y,z))(w),
I(x, v,w) = Fq (I(x, y, w))(v) ,
I(u,v,w) = Fq(I(x,v,w))(u) .
Это значит, что 3D преобразование Фурье, заданное на декартовой сетке, может быть получено после применения 1D преобразования Фурье вдоль каждой координаты. Таким образом, зная оператор, обратный
Fq , можно найти обратное 3D преобразование Фурье. Нахождение оператора Fq-1 описано в [1]. В нашем случае для m=3n+1:
(!&■№) = Zy<u)e2""3t/m , k = _nl2,...,nl2 -1,
u=-n 2
42
РИ, 2006, № 3
* ~ т,
где Fd - оператор, сопряженный Fq ,
(FDFD)k,l
n/2
Еy(u)e
2niu3
(k-l)
m
u=-n/2
,k,l = - n/2,
,n/2 -1.
* —I *
Поскольку y = Fqx ^ x = (FqFq) FQy , можно найти оператор, обратный Fq , и восстановить I(x,y,z) из I(u,v, w).
8. Оценка погрешности
Тестирование данного метода вычисления обратного преобразования Радона проводилось на трехмерных изображениях размерности n х n х n, где n принимало значения степеней двойки. Элементами этих изображений являлись случайные числа, имеющие равномерное распределение на интервале [0, 1 ]. Для каждого изображения сначала вычислялись значения прямого преобразования Радона на псевдополярной сетке, а затем, согласно изложенному выше алгоритму, вычислялось обратное преобразование Радона. В программе было реализовано два метода вычисления значений тригонометрического полинома (3). Метод вычисления значений полинома (3), основанный на решении системы линейных алгебраических уравнений, будем называть точным методом, а метод, основанный на формулах (11)-(13) - приближенным методом. Полученные в результате двух преобразований изображения сравнивались с исходным. Параметры ПЭВМ, на котором производились вычисления: процессор Atlon 1600, 512 Mb ОЗУ. Оценка относительной погрешности реконструкции изображения проводилась по формуле:
Вычисление обратного преобразования Радона было реализовано последовательным применением прямого одномерного преобразования Фурье, а затем обратного псевдополярного трехмерного преобразования Фурье. В свою очередь, для обратного псевдополярного трехмерного преобразования Фурье был применен алгоритм, состоящий из двух шагов. Первый шаг состоял в отображении преобразования Фурье, вычисленного на псевдополярной сетке, на декартово множество (метод «чистки лука»). На втором шаге из полученного преобразования Фурье методом инвертирования десятичных частот восстанавливалось исходное изображение I(x, y, z).
Было проведено исследование зависимости относительной погрешности прямого и обратного преобразования Радона от размерности и порядка исходных данных. В качестве исходных данных были использованы изображения размерности n х n х n, где n = 2i,i = 2...7 .
Научная новизна проведенных исследований состоит в разработке и реализации численного метода нахождения 3D дискретного обратного преобразования Радона, который дает значительно меньшую погрешность преобразования по сравнению с рассмотренными существующими методами.
Практическая значимость рассматриваемой задачи состоит в том, что преобразование Радона применяется в самых различных задачах математической физики и прикладной математики. Наиболее широкое изучение аспектов проблемы реконструкции изображений обусловлено практическим применением в медицине.
max | I(x, y, z) - I(x, y, z) |
x,y,z
max | I(x, y, z) |
x,y,z
где I - исходное изображение, I - восстановленное изображение. Результаты приведены в таблице.
n Т очный метод Приближенный метод
E t, с E t, с
4 1.2454e-14 0.11 1.0153e-15 0.98
8 2.9173e-13 0.380 1.7643e-15 1.64
16 5.9752e-11 2.250 2.5664e-14 7.31
32 2.5513e-10 17.740 1.2676e-12 39.01
64 1.0717e-09 254 1.8581e-09 292.74
128 2.1649e-08 4640 1.0156 5734
9. Выводы
Для вычисления прямого трехмерного дискретного преобразования Радона был применен алгоритм, состоящий из O(NHogN) шагов, где N=n3 - число пикселей в изображении.
В качестве перспектив развития данной темы можно рассматривать изучение 3D непрерывного прямого и обратного преобразования Радона и их приложений.
Литература: 1. Averbuch A., Shkolnisky Yoel. 3D Fourier based discrete Radon transform. Applied and Computational Harmonic Analysis. 2003. Vol. 15. P.33-69. 2. Гельфанд И.М., Граев М.И., Виленкин Н.Я. Интегральная геометрия и связанные с ней вопросы теории представлений. М.: Физматгиз, 1962. 3. Смоленцев Н.К. Основы теории вейвлетов. Вейвлеты в MATLAB. М: ДМК, 2005. 4. DuttA. and Rokhlin V. Fast Fourier transforms for nonequispaced data, II. Applied and Computatational Harmonic Analysis, 1995. Vol. 2. P. 85-100.
Поступила в редколлегию 17.06.2006 Рецензент: д-р физ.-мат. наук, проф. Литвин О.Н.
Тевяшев Андрей Дмитриевич, д-р техн. наук, профессор, зав кафедрой ПМ ХНУРЭ. Научные интересы: системный анализ и теория оптимального стохастического управления. Увлечения и хобби: теннис, горные лыжи, рыбалка. Адрес: Украина, 61166, Харьков, пр. Ленина, 14, тел. 70-21-436.
Смирнова Виктория Сергеевна, аспирантка кафедры ПМ ХНУРЭ. Научные интересы: математическое моделирование нестационарных режимов транспорта и распределения природного газа. Адрес: Украина, 61166, Харьков, пр. Ленина, 14, тел. 70-21-436.
РИ, 2006, № 3
43
