|
Текст подпрограммы и версий rsc2c_p.zip |
Тексты тестовых примеров trsc2c_p.zip |
Вычисление двумерных периодических свеpток двух комплексных матриц с применением алгоритма быстрого преобразования Фурье.
Заданы две комплексные матрицы A, B одинаковой размерности N1 * N2 с элементами
An1, n2 , n1 = 0, 1, ..., N1 - 1 , n2 = 0, 1, ..., N2 - 1 , Br1, r2 , r1 = 0, 1, ..., N1 - 1 , r2 = 0, 1, ..., N2 - 1 .
Требуется вычислить периодическую (круговую) свеpтку C этих матриц, т.е.
N1-1 N2-1
(1) ∑ ∑ As1 - r1, s2 - r2 Br1, r2 = Cs1, s2
r1=0 r2=0
s1 = 0, 1, ..., N1 - 1 , s2 = 0, 1, ..., N2 - 1 ,
где элементы матрицы A продолжены периодически для отрицательных значений индексов с периодом N1 в каждом столбце и с периодом N2 в каждой стpоке:
A- n1, n2 = AN1 - n1, n2 , n1 = 1, 2, ..., N1 - 1 , n2 = 0, 1, ..., N2 - 1 , An1, - n2 = An1, N2 - n2 , n1 = 0, 1, ..., N1 - 1 , n2 = 1, 2, ..., N2 - 1 .
Пусть F - оператор двумерного Дискретного Преобразования Фурье, задаваемый формулой (например, для матрицы A):
Am1, m2 = F(An1, n2) =
N1-1 N2-1
= ∑ ∑ An1, n2 * exp [ - 2π i ( n1*m1/N1 + n2*m2/N2)] ,
n1=0 n2=0
m1 = 0, 1, ..., N1 - 1 , m2 = 0, 1, ..., N2 - 1 , i = √-1 .
Аналогично Bm1, m2 = F (Br1, r2), Cm1, m2 = F (Cs1, s2).
По теоpеме о дискретной периодической свеpтке [1]
(2) Am1, m2 * Bm1, m2 = Cm1, m2 , m1 = 0, 1, ..., N1 - 1 , m2 = 0, 1, ..., N2 - 1 .
Искомая матрица С находится путем обратного Дискретного Преобразования Фурье F - 1 от произведения ДПФ матриц A и B:
Cs1, s2 = F - 1 ( Cm1, m2 ) = N1-1 N2-1 = 1/(N1*N2) ∑ ∑ Am1,m2* Bm1,m2* exp[2π i (s1*m1/N1 + s2*m2/N2)], m1=0 m2=0 s1 = 0, 1, ..., N1 - 1 , s2 = 0, 1, ..., N2 - 1 .
Подпрограмма предусматривает работу в двух режимах (в зависимости от значения параметра L).
| 1. | L=1. | Производится вычисление свертки (1) с вычислением ДПФ каждой из двух матриц A, B. |
| 2. | L=2. | Производится вычисление свертки (1) в предположении, что ДПФ матрицы A уже вычислено. Этот режим позволяет избежать повторного получения комплексной матрицы A при многократном вычислении свертки (1) с той же матрицей A. |
Для эффективного вычисления ДПФ и ОДПФ применяется двумерное Быстрое Преобразование Фурье [2]. При этом время вычислений пропорционально N1 N2 (log2 N1 + log2 N2), а объем памяти пропорционален N1 * N2.
Аналогичный алгоритм для вычисления одномерных периодических сверток двух последовательностей описан в подпрограмме AMC1R.
| 1. | Л.Рабинер, Б.Гоулд. Теория и применение цифровой обработки сигналов. M., "Мир", 1978. |
| 2. | В.А.Морозов, Н.Н.Кирсанова, А.Ф.Сысоев. Комплекс алгоритмов быстрого преобразования Фурье дискретных рядов, сб. "Численный анализ на ФОРТРАНЕ", вып.15, Изд - во МГУ, M., 1976. |
procedure RSC2C(var AR :Array of Real; var AI :Array of Real;
var BR :Array of Real; var BI :Array of Real;
var N1 :Integer; var N2 :Integer; L :Integer);
Параметры
| AR,AI - | двумерные вещественные массивы размера N1 * N2, содержащие соответственно вещественные и мнимые части элементов заданной матрицы A; |
| BR, BI - | двумерные вещественные массивы размера N1 * N2, содержащие соответственно вещественные и мнимые части элементов заданной матрицы B; на выходе подпрограммы содержат соответственно вещественные и мнимые части элементов вычисленной матрицы C; |
| N1 - | заданное число стpок матриц A и B, N1 ≥ 4 - целая степень числа два (тип: целый); |
| N2 - | заданное число столбцов матриц A и B, N2 ≥ 4 - целая степень числа два (тип: целый); |
| L - | параметр, определяющий режим использования подпрограммы (тип: целый); |
| L=1 - | вычисление свертки C матриц A, B в общем случае, |
| L=2 - | вычисление свертки C матриц A, B в предположении, что вещественные и мнимые части элементов ДПФ матрицы A уже содержатся соответственно в массивах AR, AI. |
Версии: нет
Вызываемые подпрограммы
| FTFTC - | вычисление двумерного дискретного или обратного дискретного преобразования Фурье комплексной матрицы размера N1 * N2, N1 = 2 J1, N2 = 2 J2 методом быстрого преобразования Фурье. |
Замечания по использованию
| 1. |
B результате работы подпрограммы RSC2C при L = 1 в массивах AR, AI содержатся соответственно вещественные и мнимые части элементов ДПФ матрицы A. | |
| 2. |
При первом обращении к подпрограмме (L = 1) необходимо задать значения параметров AR, AI, BR, BI, N1, N2. При повторном вычислении свертки (1) с той же матрицей A (L = 2) можно менять только параметры BR, BI. | |
| 3. |
Kак следует из формулы (2), сомножители A, B в свеpтке (1) перестановочны, поэтому при необходимости повторного вычисления свертки с той же матрицей B достаточно при L = 1 матрицу B задать в массивах AR, AI, а матрицу A - в массивах BR, BI и далее воспользоваться режимом при L = 2. |
Unit TRSC2C_p;
interface
uses
SysUtils, Math, { Delphi }
LStruct, Lfunc, UtRes_p, RSC2C_p;
function TRSC2C: String;
implementation
function TRSC2C: String;
var
N1,N2,I,J,L :Integer;
AR :Array [0..31] of Real;
АI :Array [0..31] of Real;
BR :Array [0..31] of Real;
BI :Array [0..31] of Real;
label
_1;
begin
Result := '';
N1 := 4;
N2 := 8;
for I:=1 to N1 do
begin
for J:=1 to N2 do
begin
AR[(I-1)+(J-1)*4] := (2.0-I)*(J-4.0);
AI[(I-1)+(J-1)*4] := 0.0;
BR[(I-1)+(J-1)*4] := 6.0-I-J;
_1:
BI[(I-1)+(J-1)*4] := I-J;
end;
end;
L := 1;
RSC2C(AR,AI,BR,BI,N1,N2,L);
Result := Result + #$0D#$0A;
for J:=1 to 8 do
begin
for I:=1 to 4 do
begin
Result := Result + Format('%20.16f ',[BR[(I-1)+(J-1)*4]]) + #$0D#$0A;
end;
end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
for J:=1 to 8 do
begin
for I:=1 to 4 do
begin
Result := Result + Format('%20.16f ',[BI[(I-1)+(J-1)*4]]) + #$0D#$0A;
end;
end;
Result := Result + #$0D#$0A;
UtRes('TRSC2C',Result); { вывод результатов в файл TRSC2C.res }
exit;
end;
end.
Результаты:
| - 16. 24. 48. 56. 48. 24. - 16. - 72. |
BR = | - 8. 32. 56. 64. 56. 32. - 8. - 64. |
| - 16. 24. 48. 56. 48. 24. - 16. - 72. |
| - 40. 0. 24. 32. 24. 0. - 40. - 96. |
| - 16. 24. 48. 56. 48. 24. - 16. - 72. |
BI = | - 24. 16. 40. 48. 40. 16. - 24. - 80. |
| - 16. 24. 48. 56. 48. 24. - 16. - 72. |
| 8. 48. 72. 80. 72. 48. 8. - 48. |