CC BY
50
УДК 535.31
ИССЛЕДОВАНИЕ МАТРИЧНОГО МЕТОДА ВЫЧИСЛЕНИЯ ВЫТЯНУТЫХ СФЕРОИДАЛЬНЫХ ФУНКЦИЙ НУЛЕВОГО ПОРЯДКА
© 2001 С.Н. Хонина
Институт систем обработки изображений РАН, г. Самара
Проводится исследование матричного метода вычисления вытянутых сфероидальных функций нулевого порядка, реализованного с использованием программного обеспечения Matlab. Точность вычисления определяется по выполнению основных свойств вытянутых сфероидальных функций, таких как инвариантность к интегральному преобразованию с sinc-ядром и преобразованию Фурье. Предложен алгоритм правильного упорядочивания функций по собственным значениям. Определены границы работоспособности данной реализации метода.
Введение
Вытянутые сфероидальные функции (ВСФ) нулевого порядка уп(с,х) являются собственными функциями интегрального уравнения с 8шс-ядром [1-4] и преобразования Фурье [5], что позволяет эффективно использовать их во многих областях науки и техники: теория синтеза антенн, восстановление объектов по изображению, сверхразрешение, теория резонаторов [1-10].
Подходы к вычислению ВСФ произвольного порядка отличаются сложностью и трудоемкостью: представление через ряды по другим, более простым базисам [11-13], решение дифференциального уравнения методом ортогональной дифференциальной прогонки с переносом граничных условий из особых точек [14]. Для расчета ВСФ нулевого порядка можно дополнительно воспользоваться присущими только им свойствами и решать интегральное уравнения на собственные значения, например, итерационно [15].
Учитывая трудоемкость упомянутых подходов, в одном случае [16,17] и чувствительность к начальному приближению, в другом [15], желательно иметь простой и устойчивый метод расчета ВСФ нулевого порядка. В работах [8,10] предлагается вместо интегрального уравнения на собственные значения решать матричное уравнение, полученное путем дискретизации.
В данной работе проводится исследование матричного метода вычисления ВСФ нулевого порядка, реализованного с исполь-
зованием программного обеспечения Ма1;1аЪ. Работоспособность оценивается по выполнению основных свойств рассчитанных функций: ортогональность, инвариантность к интегральному преобразованию с 8Іпс-ядром и преобразованию Фурье. Правильность упорядочивания функций определяется по убыванию собственных значений. Целью работы является определение границ работоспособности метода.
Основные свойства ВСФ нулевого порядка
ВСФ произвольного порядка являются решением скалярного волнового уравнения в вытянутых сфероидальных координатах. После разделения переменных трехмерного волнового уравнения в этих координатах возникает дифференциальное уравнение [18]:
—(1 - x') МЙ +
dx dx
2 2 X-c x -
m
1 - x 2
y(x) = 0, (1)
которое имеет непрерывное вещественное решение у(х), ограниченное при любых х. С точностью до константы решение единственно. Здесь и далее, с - произвольное положительное число, т - неотрицательное целое.
Обозначим специальные значения с, при которых существуют непрерывные решения через %тп, п = т, т+1... . Эти собственные значения могут быть так обозначены,
поскольку 0<с <с <... . Соответствую' /*т,т /*т,т+1 *
щие им решения обозначим ут и(х). Собственные функции Ути(х) четные, если п-т - четное число и нечетные, если п-т - нечетное
число. Также они имеют ровно п-т нулей на открытом интервале |х|<1.
Функции ут и(х) известны [11-13] как, вытянутые сфероидальные волновые функции порядка т. Однако особый интерес представляют ВСФ нулевого порядка. Благодаря ряду замечательных свойств [1-5] они играют важную роль во многих прикладных задачах [1-10].
Функции у0 и(х), которые далее будем обозначать просто Уи(х) удовлетворяют интегральному уравнению [18]:
. . . г б1{с У - х)]
ЯпУп (у) = ]----г-Х)йх' (2)
-Г Р (у - х )
где Я - собственные числа:
п
1 > Я0 >1 > ... > 0 . (3)
При этом собственные числа интегрального преобразования определяют количество энергии соответствующей собственной функции, концентрирующейся на данном ограниченном интервале. Собственные числа близкие к единице показывают, что данная ВСВФ имеет за пределами данного интервала малую долю энергии. Если рассматривать 1 как функцию от п, то собственные числа имеют характерную особенность поведения: до некоторого п=2с/р они имеют значения близкие к единице, а после - резко спадают до нуля.
ВСФ нулевого порядка являются также собственными функциями преобразования Фурье на ограниченном интервале [5,7,18]:
бесконечном интервале [5, 7, 18].
Исследование матричного метода
Рассмотрим модифицированное интегральное уравнение (2):
Т /2
КУп у)= | Уп Х)
2с
Біп----- (у - X)
Т
dx, (6)
-Т/2 Р(У - X)
где параметр с = к ЖТ, Т - ширина интервала в объектной области, Ж - ширина интервала в частотной области.
Пусть N - число точек дискретизации, тогда получим дискретный вариант уравнения (6)
ЛУп ук )= Еуп Хх)
. і 2с Бі:і— (ук -хх) Т
а
пУп у)= |еХР1СХУп Х)ЙХ
к = 0,Ы , (7)
который можно записать в матричном виде:
ЯпУп = 2 -Уп . (8)
Решать уравнение (8) можно с помощью прикладных пакетов программ, например, Ма1;1аЪ. Заметим, что матрица 5 квадратная, и число точек дискретизации N определяет также общее количество вычисляемых функций.
Проверка свойства ортогональности.
Для проверка свойства ортогональности функций, рассчитанных матричным методом (7)- (8) вычислялись следующие величины:
с I 2
—К .
2Р 1 '
(4)
В [1] предложено следующее правило нормировки:
Т /2
|Уп ХУт Х)ІХ = ^тп
-Т/2
|Уп ХУт Х^Х =^шп •
(5)
Также ВСФ нулевого порядка являются собственными функциями для 8тс-преобра-зования (2) и преобразования Фурье (4) на
ХУп хк )ут хк) ____
А пт = . п,т = 0,м .(9)
ХУп2 хк ) к = 0
Критерием выполнения свойства ортогональности (5) может служить величина:
м м
С =ХХ|А пт |+Х|1-Д пт |> (10)
п =0т ^п п = 0
которая в идеале д.б. равна нулю. В таблице 1 приведены значения величины Б для десяти первых функций (в ур. (10) М=9) с различным параметром с.
1=0
1
Таблица 1. Выполнение свойства ортогональности для М = 9, N = 101
с=2 с=5 с=10 с=20
Б 1,0463-10-14 1,4003 10-14 9,6979-10-15 1,9848-10-14
Из таблицы 1 можно сделать вывод о выполнении свойства ортогональности для полученных собственных функций.
Проверка выполнения свойства инвариантности к преобразованию с $тс-ядром.
Выполнение этого свойства оценивалось по среднеквадратичному отклонению
"выходной" уп (у) функции после преобразования (7) от "входной" у (х):
N —1
Хк)-Уп Хк)]2
(11)
соответственно. Значения Яп в таблицах 2 и 3 вычислялись исходя из ур. (2) как отношение максимумов функций уп (у) и уп(х) на одинаковых интервалах, то есть
, тах уп у)}
[-Т /2Т /2]_____________
тах Уп (х)}-
[-Т /2,Т /2] п
(12)
где функции уп (у), уп(х) предварительно выровнены по максимуму.
Результаты для некоторых значений параметра с=10, 20 приведены в таблице 2 и 3,
Из таблиц 2, 3 видно, что не выполняется свойство о не возрастании собственных чисел 1 . Ранжировка полученных функций была проведена в соответствии с (3), что скорректировало порядок, установленный встроенной функцией Ма1!аЪ.
Из таблиц 2 и 3 также видно, что свойство (2) нарушается, начиная с некоторого номера п, зависящего от параметра с, при чем пс > п0. Так, например,
- для с=2: п0=1, пс=8,
- для с=5: п0=3, п =11,
Таблица 2. Значения 1п и ошибки 5 (11) для с=10, N=101, М=18
Номер функции (Ма^аЬ) Скорректированный номер функции Значение Ошибка, б
6 0 0,99999 1,1781-10-15
5 1 0,99998 7,7543 -10-16
4 2 0,99993 7,2443-10-16
3 3 0,9982 6,5013-10-16
2 4 0,9778 7,3094-10-16
0 5 0,8419 7,0206-10-16
1 6 0,4678 1,0162-10-15
7 7 0,1254 2,7448-10-15
8 8 0,0173 1,0685-10-14
9 9 0,0016 7,1481 -10-14
10 10 0,0001 1,2566-10-12
11 11 5,7821-10-6 1,8201 -10-11
12 12 2,6109-10-7 4,2788-10-10
13 13 9,9312-10-9 1,5854-10-8
14 14 3,2293-10-10 3,7692-10-7
15 15 9,0837-10-12 1,6731-10-5
16 16 2,2325-10-13 4,3169-10-4
17 17 4,8268-10-15 0,0204
18 18 1,5834-10-16 2,6302
Таблица 3. Значения 1п и ошибки 5 (11) для с = 20, N = 101, М = 27
Номер функции (ММаЪ) Скорректированный номер функции Значение 1 Ошибка, 8
16 0 0,999999 9,7149-10-15
15 1 0,999998 8,6180-10-16
14 2 0,999996 1,4720-10-15
13 3 0,999992 1,0531 -10-15
12 4 0,99999 1,2603-10-15
11 5 0,99997 8,5340-10-16
9 6 0,99993 1,2798-10-15
8 7 0,9999 8,0873-10-16
6 8 0,9998 5,4342-10-16
5 9 0,9981 1,0142-10-15
3 10 0,9840 2,2465-10-15
1 11 0,9019 1,3170-10-15
0 12 0,6365 2,2044-10-15
4 13 0,0621 1,1603-10-14
2 14 0,2667 2,6452-10-15
7 15 0,0095 2,7672-10-14
10 16 0,0011 2,3090-10-13
17 17 1,0868 10-4 1,7534-10-12
18 18 9,157110-6 2,9481-10-11
19 19 6,7550-10-7 4,1404-10-10
20 20 4,4129-10-8 5,2743-10-9
21 21 2,5752-10-9 1,0669-10-7
22 22 1,3519-10-10 1,7211 -10-6
23 23 6,4209-10-12 3,0271 -10-5
24 24 2,7739-10-13 8,1059-10-4
25 25 1,0941 -10-14 0,0195
26 26 4,9774-10-16 0,4351
27 27 3,9142-10-16 0,9099
- для с=10: п0=6, пс=16,
- для с=20: п0=12, пс=24.
Индексы п , п0 определяют границы применимости метода:
- для п < пп значения 1 » 1, и величина
0 п 5
ошибки 5 <<1% и имеет значение отличное от нуля в четырнадцатом порядке;
- для п0 < п < пс значение ошибки 5 приемлемо и порядок ее не более 1%, причем числа 1 имеют значения отличные от нуля в порядке меньше или равном пятнадцати;
- для п > пс величина ошибки 5 становится не приемлемой, а числа 1 имеют значения отличные от нуля в порядке больше или равном шестнадцати;
Следует отметить, что в прикладных задачах наибольший интерес представляют случаи, когда 1 » 1, когда рассматриваемый метод демонстрирует хорошие результаты.
Проверка выполнения свойства инвариантности к преобразованию Фурье.
Выполнение этого свойства оценивалось по среднеквадратичному отклонению "выходной" уп (у) функции после преобразования (4) от "входной" Уп(х).
На рис. 1, 2 приведены графики абсолютных значений собственных функций Уп(х) и их Фурье-образов уп (у) для различных параметров. Из рисунков 1 и 2 видно, что при
аЬв(1тО
а)
аЬэ(1т1)
б)
аЬэ(1т|)
в)
Рис. 1. Графики абсолютный значений собственный функций уп(х) (аЪв(/п)) и их Фурье-образов уп (у) (аЪя^п)) для N=101, с=2, п=0 (а), п=1 (б), п=3 (в)
значениях 1 » 1 энергия ВСФ сконцентрирована на интервале [-Т /2, Т /2], и, наоборот, при значениях 1 » 0 вне этого интервала. Однако, вид функции внутри интервала после преобразования Фурье не меняется с точностью до множителя.
Из таблицы 4 видно, что ошибка 5 становится не приемлемой только, если собственные значения 1 намного меньше еди-
п
ницы, при значениях 1 близких к единице значение ошибки стремится к нулю.
Недостатки метода Выше было показано, что полученные собственные функции удовлетворяют основ-
аЬэ(1т|)
а)
аЬэ^п)
аЬэ(Рп)
в)
Рис. 2. Графики абсолютныгх значений собственныгх функций у (х) (аЪ$([п)) и их Фурье-образов уп (у) (аЪв(¥п)) для N=101, с=10, п=6 (а), п=1 (б), п=8 (в)
ным свойствам ВСФ при 1п » 1.
Однако, из таблиц 2-4 видно несоответствие номеров собственных функций порядку собственных чисел. Правильно отранжи-ровать ВСФ можно, используя (12) и (3).
При уменьшении 1п наблюдается не только нарушение свойств ВСФ, но и искажение вида функций. Например, собственные функции на рис.3. Нужно при этом отметить фильтрующее свойство преобразования Фурье - после его применения функции приняли стандартный вид. Преобразование с Бшс-ядром такого эффекта не дало. Однако "выглаживание" рассчитанных функций с помощью (4) снижает ошибку в выполнении
Таблица 4. Значения отклонения абсолютного значения ВСФ и их Фурье-образа для М = 9, N = 101
Номер функции (МаЙаЪ) с=2 с=10
Собственные числа, 1 Ошибка, 8 Собственные числа, 1 Ошибка, 8
0 0,8929 0,0046 0,8394 0,0419
1 0,3678 0,0047 0,4719 0,0537
2 0,0395 0,0169 0,9796 0,0374
3 0,0014 0,0439 0,9984 0,0307
4 0,0000 0,0804 0,9988 0,0228
5 0,0000 0,1244 0,9996 0,0154
6 0,0000 0,1743 0,9999 0,0099
7 0,0000 0,2286 0,1422 0,0728
8 0,0000 0,2840 0,0225 0,1253
9 0,0000 0,2805 0,0024 0,1866
(2) только при 1 » 1, в остальных случаях ошибка даже увеличивается.
Заключение
В данной работе реализован и исследован матричный метод вычисления ВСФ нулевого порядка с использованием программного обеспечения МайаЬ. Исследование метода заключается в проверке получаемых функций на выполнение основных свойств. Результатом данной проверки являются следующие выводы:
1. Данный метод позволяет получать удовлетворительную точность только в некотором диапазоне изменения индекса п, обеспечивающие 1 » 1. Верхней границей диапазона можно считать п„ = 2с/я.
2. Выполнение свойства инвариантности к интегральному преобразованию с 8тс-ядром обеспечивается с более высокой точностью, чем свойство инвариантности к преобразованию Фурье. В первом случае отклонение составляет менее 1% даже для функций с 1п » 10-13. Во втором случае отклонение было меньше 5% только для п<п0.
3. Недостатком метода также является неправильное упорядочивание вычисляемых функций, которую нужно корректировать в соответствии с убыванием 1п.
4. При 1 » 0 наблюдается не только нарушение свойств ВСФ, но и искажение вида функций. "Выгладить" такие функции можно с помощью преобразования Фурье, однако ошибка при этом снижается только при 1 »1.
-С 5
а) б)
Рис. 3. Графики абсолютных значений собственных функций Уп(х) (фп) и их Фурье-образов уп (у) (Рп) для N = 101, с = 2, п = 9 (а), п = 10 (б)
Тем не менее, учитывая простоту реализации, рассмотренный метод может успешно применяться в практических задачах, использующих ВСФ с собственными чисел близкими к единице.
Работа поддержана Российским фондом фундаментальных исследований (гранты 0015-96114, 00-01-00031).
СПИСОК ЛИТЕРАТУРЫ
1. Slepian D., Pollak H.O. Prolate spheroidal wave functions. Fourier Analysis and Uncertainty
- I // Bell Syst. Tech. J. V.40. 1961.
2. Landau H.J., Pollak H. O. Prolate spheroidal wave functions. Fourier Analysis and Uncertainty
- II // Bell Syst. Tech. J. V.40. 1961.
3. Slepian D., Pollak H. O. Prolate spheroidal wave functions, Fourier analysis and uncertainty - III: The dimension of essentially time-and band-limited signals // Bell Syst. Tech. J. V.41. 1962.
4. ХургинЯ.И., Яковлев В.П. Методы теории целых функций в радиофизике, теории связи и оптике. М.: Физматгиз, 1962.
5. De Santis P., Palma C. Degrees of freedom of aberrated images // Opt. Acta. V.23. 1976. №9.
6. Bertero M., Viano G.A., De Mol C. Resolution beyond the diffraction limit for regularized object restoration // Opt. Acta. V.27. 1980. №3.
7. Bertero M., Pike E.R. Resolution in diffraction-limited imaging, a singular value analysis. I. The case of coherent illumination // Opt. Acta. V.29. 1982. №6.
8. Barbosa L.C. A maximum-energy concentration spectral window // IBM J. Res. Develop. V.30. 1986. №3.
9. Rogers M.E., Erkkila J.H. Resonator mode analysis using linear prolate functions // Appl.
Opt. V.22. 1983. №13.
10. Latham W.P., Tilton M.L. Calculation of prolate functions for optical analysis // Appl. Opt. V.26. 1987. №13.
11. Фламмер К. Таблицы волновых сфероидальных функций. М.: ВЦ АН СССР, 1962.
12. Комаров И. В., Пономарев Л. П., Славя-нов С. Ю. Сфероидальные и кулоновские сферои-дальные функции. М.: Наука, 1976.
13.Лоуэн А. Сфероидальные волновые функции / Справочник по специальным функциям под ред. М. Абрамовица и И. Сти-ган. М.: Наука, 1979.
14. Абрамов А.А., Дышко А.Л., Конюхова Н.Б., Пак ТВ., Парийский Б.С. Вычисление вытянутых сфероидальных функций решением соответствующих дифференциальных уравнений // Ж. вычисл. матем. и матем. физ. Т.24. 1984. № 1.
15.Хонина С.Н., КотлярВ.В. Вытянутые сфероидальные функции в дифракционной оптике // Материалы международной молодежной школы по оптике, лазерной физике и биофизике. Саратов: Сарат. ун-т, 2000.
16. Хонина С.Н. Приближение сфероидальных волновых функций конечными рядами // Компьютерная оптика. Самара, 1999. Вып.19.
17. Волотовский С.Г., Казанский Н.Л., Хонина С.Н. Анализ и разработка методов вычисления собственных значений вытянутых сфероидальных функций нулевого порядка // Труды 5-ой международной конференции РОАИ-5-2000. Т.4. Самара, 2000.
18. D. Slepian, Some asymptotic expansions for prolate spheroidal wave functions // J. Math. & Phys. V.44. 1965.
STUDIES INTO A MATRIX TECHNIQUE FOR CALCULATING PROLATE SPHEROIDAL FUNCTIONS OF ZERO ORDER
© 2001 S.N. Khonina
Image Processing Systems Institute of Russian Academy of Science, Samara
A matrix technique for calculating prolate spheroidal functions of zero order that is implemented using the Matlab tools is studied. The calculational accuracy is assessed on the basis of fulfillment of the major properties of the prolate spheroidal functions: invariance to the integral transform with sinc-kernel and to the Fourier transform. An algorithm for properly ordering the functions using their eigenvalues is proposed. Applicability limits of the proposed technique implementation are determined.
