Текст подпрограммы и версий
rsc1r_p.zip 		Тексты тестовых примеров
trsc1r_p.zip

Подпрограмма:  RSC1R (модуль RSC1R_p)

Назначение

Построение автокорреляционных и взаимных корреляционных функций двух случайных процессов по оценкам их спектров.

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

Пусть SXm, SYm, SXYm - заданные значения автоспектров и взаимного спектра двух случайных процессов  X, Y начастотах  λm = m * Δλ, где  m = 0, 1, ..., NC - 1,  Δλ = 2π/(N * DT) - шаг сетки по частотной оси, N = 2 (NC - 1), DT - шаг сетки по оси времени. Доопределим (внутри подпрограммы) функции SXm, SYm, SXY1m = Re SXYm симметрично относительно точки  m = NC - 1, а функцию SXY2m = Im SXYm антисимметрично относительно  m = NC - 1. Применяя обратное дискретное преобразование Фурье к  N значениям спектров, получаем соответствующие оценки для автокорреляционных и взаимной корреляционной функций [1]: 

                                   N-1
   CXs  =  1/(N*DT)    ∑    SXm exp( i λm ts)  =
                                  m=0
                                                                         N-1   
                                                   =  1/(N*DT)    ∑   SXm exp(2π i m s / N) ,
                                                                         m=0
                                   N-1
   CYs  =  1/(N*DT)    ∑    SYm exp( i λm ts)  =
                                  m=0
                                                                         N-1   
                                                   =  1/(N*DT)    ∑   SYm exp(2π i m s / N) ,
                                                                         m=0 
                                      N-1
   CXYs  =  1/(N*DT)    ∑    SXYm exp( i λm ts)  =
                                     m=0
                                                                         N-1   
                                                   =  1/(N*DT)    ∑   SXYm exp(2π i m s / N) ,
                                                                         m=0

где   ts = S * DT,   s = 0, 1, ..., N - 1,   i = √ -1.

Построенные оценки CXs, CYs, являются вещественными функциями симметричными относительно центра S = N/2 = NC - 1. Эти оценки "циклические" [1, стp.357], поэтому рекомендуется применять данный метод для быстро затухающих корреляций.

Полное описание реализованного алгоритма содержится в статье [2] (подпрограмма CORREL).

1.  Дж.Бендат, А.Пирсол, Измерение и анализ случайных процессов, Изд - во "Мир", M., 1974.
2.  М.В.Арефьева, Корреляционный и спектральный анализ стационарных случайных процессов (часть 1), сб. "Численный анализ на ФОРТРАНе", вып.15, Изд - во МГУ, M., 1976.

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

procedure RSC1R(var SX :Array of Real; var SY :Array of Real;
                var SXY1 :Array of Real; var SXY2 :Array of Real;
                NC :Integer; DT :Real; var CX :Array of Real;
                var CY :Array of Real; var CXY :Array of Real;
                var R :Array of Real);


   Параметры
SX - одномерный массив длины NC, содержащий заданные значения автоспектра процесса  X (тип: вещественный);
SY - одномерный массив длины NC, содержащий заданные значения автоспектра процесса  Y (тип: вещественный);
SXY1 - одномерный массив длины NC, содержащий заданные значения вещественной части взаимного спектра процессов  X, Y (тип: вещественный);
SXY2 - одномерный массив длины NC, содержащий заданные значения мнимой части взаимного спектра процессов  X, Y (тип: вещественный);
NC - количество заданных значений автоспектров и взаимного спектра процессов  X, Y,  NC = 2n + 1, где n ≥ 1 - целое число (тип: целый);
DT - заданный шаг сетки по оси времени, DT > 0 (тип: вещественный);
CX - одномерный массив длины NC, содержащий вычисленные значения автокорреляционной функции процесса  X (тип: вещественный);
CY - одномерный массив длины NC, содержащий вычисленные значения автокорреляционной функции процесса  Y (тип: вещественный);
CXY - одномерный массив длины 2 NC - 1, содержащий в первых (2 NC - 2) элементах вычисленные значения взаимной корреляционной функции процессов  X, Y (тип: вещественный);
R - одномерный массив длины 2 NC - 1, используемый в подпрограмме как рабочий (тип: вещественный). 
   Версии:  нет 
   Вызываемые подпрограммы
FTF1C - подпрограмма вычисления дискретного или обратного дискретного преобразования Фурье комлексного ряда длины, равной степени двух, методом быстрого преобразования Фурье. 
   Замечания по использованию
  1. 

C целью экономии числа используемых массивов результаты можно получать на месте исходной информации, а именно, допустимы совпадения параметров: CX = SX,  CY = SY.

  2.  B частных случаях применения подпрограммы возможны совмещения массивов: SX = SY = SXY1 = SXY2,  CX = CY - при построении автокорреляционной функции одного случайного процесса  X и SY = SXY1 = SXY2 - при построении только автокорреляцонных функций двух случайных процессов  X, Y без вычисления их взаимной корреляции; в этих случаях массив CXY используется в подпрограмме как рабочий.

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

Unit TRSC1R_p;
interface
uses
SysUtils, Math, { Delphi }
LStruct, Lfunc, UtRes_p, RSC1R_p;

function TRSC1R: String; 

implementation

function TRSC1R: String;
var
NC,_i :Integer;
DT :Real;
СХ :Array [0..2] of Real;
CY :Array [0..2] of Real;
CXY :Array [0..4] of Real;
R :Array [0..4] of Real;
const
SX :Array [0..2] of Real = ( 25.0,2.0,1.0 );
SY :Array [0..2] of Real = ( 25.0,2.0,1.0 );
SXY1 :Array [0..2] of Real = ( 25.0,-2.0,-1.0 );
SXY2 :Array [0..2] of Real = ( 0.0,0.0,0.0 );
begin
Result := '';
NC := 3;
DT := 1.0;
RSC1R(SX,SY,SXY1,SXY2,NC,DT,CX,CY,CXY,R);
Result := Result + #$0D#$0A;
for _i:=0 to 2 do
 begin
  Result := Result + Format('%20.16f ',[CX[_i]]);
  if ( ((_i+1) mod 3)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
for _i:=0 to 2 do
 begin
  Result := Result + Format('%20.16f ',[CY[_i]]);
  if ( ((_i+1) mod 3)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + #$0D#$0A;
for _i:=0 to 4 do
 begin
  Result := Result + Format('%20.16f ',[CXY[_i]]);
  if ( ((_i+1) mod 4)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
UtRes('TRSC1R',Result);  { вывод результатов в файл TRSC1R.res }
exit;
end;

end.

Результаты:

      CX    =  (7.5,  6.,  5.5)
      CY    =  (7.5,  6.,  5.5)
      CXY  =  (5.,  6.5,  7.,  6.5)