Документ взят из кэша поисковой машины. Адрес
оригинального документа
: http://num-anal.srcc.msu.su/lib_na/cat/r_htm_c/rsa1r_c.htm
Дата изменения: Fri Dec 4 15:26:15 2015 Дата индексирования: Sun Apr 10 00:59:10 2016 Кодировка: Windows-1251 |
Текст подпрограммы и версий rsa1r_c.zip , rsa1d_c.zip |
Тексты тестовых примеров trsa1r_c.zip , trsa1d_c.zip |
Вычисление прямой и обратной дискретной периодической (круговой) свертки двух одномерных массивов вещественных чисел.
Пусть известны N значений sj = s ( tj ) функции сигнала s ( t ) на равномерной сетке tj = j Δt, j = 0, 1, 2, ..., N - 1 (где Δt - шаг сетки, а N равняется целой степени двух). Предполагается, что сигнал sj является периодическим с периодом N (отсюда, в частности, следует, что значения sj известны и для отрицательных значений j, т.е. s- j = sN - j, j = 1, 2, ..., N - 1). Значения sj задаются в вещественном одномерном массиве ARRAY длины N: ARRAY ( J ) = sj ( J = 1, 2, ..., N; j = 0, 1, 2, ..., N - 1).
Пусть известны также M значений rk = r ( tk ) импульсной функции отклика r ( t ) на той же равномерной сетке tk = k Δt, k = - (M - 1)/2, ..., 0, ..., (M - 1)/2 (где Δt - тот же шаг сетки, что и для функции сигнала, а M - нечетное натуральное число такое, что M < N). Значения rk задаются в M первых элементах одномерного вещественного массива RESP длины N следующим образом:
RESP(1) = r0 RESP( I ) = rk , I = 2, 3, ..., (M+1)/2 ; k = 1, 2, ..., (M-1)/2 RESP( I ) = rk , I = (M+3)/2, ..., M ; k = (M-1)/2, ..., - 1
Подпрограмма rsa1r_c имеет два режима работы, задаваемых при обращении к ней значением параметра IREG. В первом режиме ( IREG = 1) выполняется прямая дискретная периодическая свертка массивов ARRAY и RESP (L = 1, 2, ..., N):
N (r*s)L = ∑ ARRAY(L - J) × RESP( J ) j =1
где * - знак свертки, а
× - знак умножения.
Предварительно в массив RESP добавляются нули следующим образом:
RESP( I ) = RESP( I ) , I = 1, 2, ..., (M+1)/2 RESP( I ) = 0.0 , I = (M+3)/2, ..., N - (M-1)/2 RESP(N+1-I ) = RESP(M+1-I ) , I = 1, 2, ..., (M-1)/2
Для вычисления указанной свертки используется быстрое дискретное преобразование Фурье в соответствии с теоремой о дискретной периодической (круговой) свертке. Полученные значения (r * s)L помещаются в первых N компонентах массива ANS.
Во втором режиме ( IREG = - 1) выполняется восстановление значений sj функции сигнала по известным значениям rk импульсной функции ответа, задаваемым в массиве RESP, и значениям соответствующей свертки, задаваемым в массиве ARRAY. Полученные значения sj помещаются в первых N элементах массива ANS.
Р.Отнес, Л.Эноксон. Прикладной анализ временных рядов. Изд - во "Мир", 1982.
int rsa1r_c (real *array, integer *n, real *resp, integer *m, complex *ans, complex *fft, integer *ireg, integer *iflag)
Параметры
array - | вещественный одномерный массив длины n, содержащий либо значения sj ( ireg = 1), либо значения свертки ( ireg = - 1); |
n - | длина массива array, равная целой степени двух (тип: целый); |
resp - | вещественный одномерный массив длины n, в первых m компонентах которого задаются значения rk; |
m - | количество заданных значений rk (тип: целый); |
ans - | вещественный одномерный массив длины 2n, в первых n элементах которого в результате работы подпрограммы размещаются либо компоненты вектора свертки (ireg = 1), либо восстановленные значения sj (ireg = - 1); |
fft - | комплексный одномерный массив длины n, используемый в подпрограмме в качестве рабочего; |
ireg - | задает режим работы подпрограммы (тип: целый); при этом: |
ireg= 1 - | когда вычисляется свертка; |
ireg= -1 - | когда восстанавливается функция сигнала sj |
iflag - | целая переменная, служащая для сообщения о характере окончания работы подпрограммы во втором режиме ( ireg = - 1); при этом: |
iflag=0 - | когда восстановлены все значения sj; |
iflag=i - | когда i - ая компонента преобразования Фурье импульсной функции отклика rk равна нулю; это означает, что значение si не может быть восстановлено из - за потери информации в векторе свертки; значения sj в этом случае подпрограммой не восстанавливаются |
Версии
rsa1d_c - | вычисление прямой и обратной дискретной периодической (круговой) свертки двух одномерных массивов вещественных чисел в режиме удвоенной точности; при этом параметры array, resp и ans должны иметь тип double, параметр fft - тип doublecomplex (это означает, что данная версия существует только для тех компиляторов, входной язык которых допускает комплексные данные удвоенной точности) |
Вызываемые подпрограммы
fta2r_c - fta2d_c | одновременное выполнение прямого быстрого дискретного преобразования Фурье двух одномерных массивов вещественных чисел, каждый из которых имеет длину n, равной целой степени двух, в режимах одинарной и удвоенной точности; используются в подпрограммах rsa1r_c и rsa1d_c соответственно |
fta3r_c - fta3d_c | выполнение прямого и обратного быстрых преобразований Фурье одномерного массива вещественных чисел длины 2n, где n равняется целой степени двух, в режимах одинарной и удвоенной точности; используются в подпрограммах rsa1r_c и rsa1d_c соответственно |
Замечания по использованию
В подпрограммах rsa1r_c и rsa1d_c проверка того, что значение n должно быть целой степенью двух, а m должно быть нечетным числом, меньшим n, не производится |
int main(void) { /* Initialized data */ static float array[16] = { 0.f,0.f,.5f,0.f,0.f,-.5f,0.f,0.f,.5f,.5f,.5f, .5f,.5f,0.f,0.f,0.f }; /* Local variables */ static float resp[16]; extern int rsa1r_c(float *, int *, float *, int *, float *, complex *, int *, int *); static int m, n, iflag, i; static complex fft[16]; static float ans[32]; n = 16; m = 3; resp[0] = .7f; resp[1] = .4f; resp[2] = .4f; rsa1r_c(array, &n, resp, &m, ans, fft, &c__1, &iflag); printf("\n %5i \n", iflag); for (i = 0; i <= 28; i += 4) { printf("\n %16.7e %16.7e %16.7e %16.7e \n", ans[i], ans[i+1], ans[i+2], ans[i+3]); } return 0; } /* main */ Результаты: ans = (- 0.596046e - 7, 0.2, 0.35, 0.2, - 0.2, - 0.35, - 0.2, 0.2, 0.55, 0.75, 0.75, 0.75, 0.55, 0.2, - 0.745058e - 8, 0.521541e - 07)