Текст подпрограммы и версий ( Фортран )
rsa1r.zip , rsa1d.zip 		Тексты тестовых примеров ( Фортран )
trsa1r.zip , trsa1d.zip
Текст подпрограммы и версий ( Си )
rsa1r_c.zip , rsa1d_c.zip 		Тексты тестовых примеров ( Си )
trsa1r_c.zip , trsa1d_c.zip
Текст подпрограммы и версий ( Паскаль )
rsa1r_p.zip , rsa1e_p.zip 		Тексты тестовых примеров ( Паскаль )
trsa1r_p.zip , trsa1e_p.zip

Подпрограмма:  RSA1R

Назначение

Вычисление прямой и обратной дискретной периодической (круговой) свертки двух одномерных массивов вещественных чисел.

Математическое описание

Пусть известны  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 имеет два режима работы, задаваемых при обращении к ней значением параметра 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.

Использование 

    SUBROUTINE  RSA1R (ARRAY, N, RESP, M, ANS, FFT, IREG,
                                             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 - вычисление прямой и обратной дискретной периодической (круговой) свертки двух одномерных массивов вещественных чисел в режиме удвоенной точности; при этом параметры ARRAY, RESP и ANS должны иметь тип DOUBLE PRECISION, параметр FFT - тип COMPLEX * 16 (это означает, что данная версия существует только для тех компиляторов, входной язык которых допускает комплексные данные удвоенной точности). 
   Вызываемые подпрограммы
       FTA2R -
       FTA2D   		одновременное выполнение прямого быстрого дискретного преобразования Фурье двух одномерных массивов вещественных чисел, каждый из которых имеет длину  N, равной целой степени двух, в режимах одинарной и удвоенной точности; используются в подпрограммах RSA1R и RSA1D соответственно.
       FTA3R -
       FTA3D   		выполнение прямого и обратного быстрых преобразований Фурье одномерного массива вещественных чисел длины 2N, где  N равняется целой степени двух, в режимах одинарной и удвоенной точности; используются в подпрограммах RSA1R и RSA1D соответственно. 
   Замечания по использованию
  В подпрограммах RSA1R и RSA1D проверка того, что значение  N должно быть целой степенью двух, а  M должно быть нечетным числом, меньшим  N, не производится.

Пример использования 

       REAL  ARRAY(16), RESP(16), ANS(32) 
       COMPLEX  FFT(16) 
       DATA  ARRAY /0.0, 0.0, 0.5, 0.0, 0.0, - 0.5, 0.0, 0.0,  
      *                           0.5, 0.5, 0.5, 0.5, 0.5, 0.0, 0.0, 0.0/ 
       N = 16 
       M = 3 
       RESP(1) = 0.7 
       RESP(2) = 0.4 
       RESP(3) = 0.4 
       CALL  RSA1R (ARRAY, N, RESP, M, ANS, FFT, 1, IFLAG) 

Результаты:

       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)