|
Документ взят из кэша поисковой машины. Адрес
оригинального документа
: http://num-anal.srcc.msu.ru/lib_na/cat/r_htm_c/rss1r_c.htm
Дата изменения: Mon Mar 2 17:20:46 2015 Дата индексирования: Sun Apr 10 01:46:11 2016 Кодировка: Windows-1251 |
|
Текст подпрограммы и версий rss1r_c.zip |
Тексты тестовых примеров trss1r_c.zip |
Построение автопериодограмм и взаимных периодограмм двух вещественных случайных процессов по одному отрезку их реализаций.
Пусть заданы реализации двух вещественных, стационарных, центрированных случайных процессов Xq, Yq, q = 0, 1, ..., M - 1 на сетке по оси времени tq = q * DT с шагом DT > 0. Рассмотрим отрезки этих реализаций длины N: 1 ≤ N ≤ M, левые границы которых a, b: 0 ≤ a ≤ M - N, 0 ≤ b ≤ M - N, вообще говоря, различны, т.е. отрезки вида: Xa + s, Yb + s, s = 0, 1, ..., N - 1, и вычислим их дискретные преобразования Фурье Am, Bm:
N-1
Am = ∑ Xa + s exp(- i λm ts ) =
s=0
N-1
= ∑ Xa + s exp(- 2π i m s / N) ,
s=0
N-1
Bm = ∑ Yb + s exp(- i λm ts ) =
s=0
N-1
= ∑ Yb + s exp(- 2π i m s / N) ,
s=0
где λm = m * Δλ, m = 0, 1, ..., N - 1, Δλ = 2π/(N * DT), i = √-1.
Соответствующие автопериодограммы и взаимные периодограммы строятся по формулам [1]:
PXm = (DT / N) | Am | 2 , PYm = (DT / N) | Bm | 2 ,
PXYm = (DT / N) A*m Bm , m = 0, 1, ..., N -1 ,
где * - знак комплексного сопряжения, при этом функции РХm, PYm и PXY1m = Re PXYm симметричны относительно точки m = N/2, а функция PXY2m = Im PXYm антисимметрична относительно m = N/2.
Периодограммы, полученные по одному отрезку реализаций, являются при N → +∞ асимптотически несмещенными оценками спектров случайных процессов, но эти оценки асимптотически несостоятельны.
Полное описание реализованного алгоритма содержится в статье [2] (подпрограмма PERSEG).
| 1. | Дж.Бендат, А.Пирсол, Измерение и анализ случайных процессов, Изд - во "Мир", M., 1974. |
| 2. | М.В.Арефьева, Корреляционный и спектральный анализ стационарных случайных процессов (часть 1), сб. "Численный анализ на ФОРТРАНе", вып.15. Изд - во МГУ, M., 1976. |
int rss1r_c (real *x1, real *y1, integer *n, real *dt,
real *px, real *py, real *pxy1, real *pxy2)
Параметры
| x1, y1 - | одномерные мссивы длины n + 1, первые n элементов которых содержат значения исходных реализаций случайных процессов Xa + s, Yb + s, s = 0, 1, ..., n - 1 (тип: вещественный); |
| n - | количество заданных значений каждого из двух исходных рядов, n ≥ 4 - целая степень числа два (тип: целый); |
| dt - | заданный шаг сетки по оси времени, dt > 0 (тип: вещественный); |
| px - | одномерный массив длины n/2 + 1, содержащий вычисленные значения автопериодограммы процесса X (тип: вещественный); |
| py - | одномерный массив длины n/2 + 1, содержащий вычисленные значения автопериодограимы процесса Y (тип: вещественный); |
| pxy1 - | одномерный массив длины n/2 + 1, содержащий вычисленные значения вещественной части взаимной периодограммы процессов X, Y (тип: вещественный); |
| pxy2 - | одномерный массив длины n/2 + 1, содержащий вычисленные значения мнимой части взаимной периодограммы процессов X, Y (тип: вещественный). |
Версии: нет
Вызываемые подпрограммы
| ftf1c_c - | подпрограмма вычисления дискретного или обратного дискретного преобразования Фурье комплексного ряда длины, равной степени двух, методом быстрого преобразования Фурье. |
Замечания по использованию
| 1. |
После выхода из подпрограммы rss1r_c значения исходных случайных процессов в массивах x1, y1 не сохраняются. | |
| 2. |
B частном случае построения автопериодограммы одного случайного процесса X в массив y1 следует заслать нули, при этом допустимы совпадения параметров: px = py = pxy1 = pxy2. | |
| 3. | B частном случае построения только автопериодограмм двух случайных процессов X, Y без вычисления их взаимной периодограммы допустимы совпадения параметров: p = pxy1 = pxy2. |
int main(void)
{
/* Initialized data */
static float x1[5] = { 1.f,2.f,3.f,4.f,0.f };
static float y1[5] = { 4.f,3.f,2.f,1.f,0.f };
/* Local variables */
static float pxy1[3], pxy2[3];
extern int rss1r_c(float *, float *, int *, float *, float *,
float *, float *, float *);
static int n;
static float dt, px[3], py[3];
n = 4;
dt = 1.f;
rss1r_c(x1, y1, &n, &dt, px, py, pxy1, pxy2);
printf("\n %16.7e %16.7e %16.7e \n", px[0], px[1], px[2]);
printf("\n %16.7e %16.7e %16.7e \n", py[0], py[1], py[2]);
printf("\n %16.7e %16.7e %16.7e \n", pxy1[0], pxy1[1], pxy1[2]);
printf("\n %16.7e %16.7e %16.7e \n", pxy2[0], pxy2[1], pxy2[2]);
return 0;
} /* main */
Результаты:
px = (25., 2., 1.) , py = (25., 2., 1.) ,
pxy1 = (25., - 2., - 1.) , pxy2 = ( 0., 0., 0.)