Спектральный анализ и его приложения Том 2 - Дженкинс Г.
Скачать (прямая ссылка):
6. Напечатать сглаженные выборочные спектральные оценки SPEC(I); 1 = 0, NF, ширину полосы частот В = 4/(3M-DELTA) и число степенен свободы D = (8*N)/(3«•M) для окна Тьюки и соответствующих значений N, M и DELTA.
7. Построить на одном и том же графике логарифмы сглаженных спектров в зависимости от частоты для всех использованных значений М.
Для построения графика последовательности чисел LOGSPEC(I) нужно найти максимальное из этих чисел, которое мы обозначим MLOG. Далее мы пользовались следующей методикой построения графика. Находим ближайшую к MLOG степень десяти D так, чтобы D >- MLOG, и строим логарифм спектра в диапазоне от D — 4 до D. Если, например, максимальное значение спектра было 2, так что MLOG = 0,303 и, следовательно, D = 1, то логарифм спектра строился в диапазоне от —3 до 1, что соответствовало значениям спектра от 0,001 до 10. Диапазон значений 104 можно считать подходящим для большинства целей, так как если требуется еще больший диапазон, то, вероятно, целесообразней расфильтровать данные, чтобы получить лучшие выборочные спектральные оценки на тех частотах, где мощность мала. Значение LOGSPEC(I) = —100 автоматически строится на графике ниже самой нижней линии, ограничивающей выбранный диапазон.
Алгоритм
Для вычисления SPEC(I)
SetC = COS VO = O., Vl=O.
Dol, K = M-I, 1
V2 = 2. * С * Vl - VO + W (К) * COV (К) VO = Vl 1 V1=V2.
SPEC (I) = 2. * DELTA * (COV (0) + 2. * (Vl * С - VO)).
Пример. Рассмотрим приведенный в разд. 7.1.1 пример, для которого М = 3, DELTA = LO, NF = 8 и COV(O)=I., W(I)•COV(I) = 0.430, W(2)•COV(2) = 0.065.
З Зак. ІІ78
66
Г лава 7
Тогда для 1 = 0, C = COS (0) = 1.,
VO = 0., Vl =0., и, проходя через цикл ,,do", получаем К = 2, V2 = 2 (1) (0)-0 + 0.065 = 0.065 VO = O Vl = 0.065
K=I, V2 = 2(1)(0.065)-0 + 0.430 = 0.560 VO = 0.065 Vl =0.560
И затем SPEC (0) = 2- 1 [1 +2(0.560- 1 -0.065)] = 3.980, что совпадает с величиной, приведенной в табл. 7.1. Для I = 1
С = COS (я/8) = 0.924, V0 = 0, Vl=O.
Проходя через цикл ,,do", получаем
К = 2, V2 = 2 (0.924) (0)-0 + 0.065 = 0.065
VO = O
Vl =0.065;
K=I, V2 = 2 (0.924) (0.065) - 0 + 0.430 = 0.550 VO = 0.65 Vl = 0.550.
Тогда SPEC(I) = 2(1) (1 +2(0.550) (0.924) —0.065)) = 3.772. Этот алгоритм, хотя работает и не так быстро, как БПФ, тем не менее имеет относительно высокую скорость, высокую точность и требует вычисления косинуса только один раз на каждую точку по частоте.
Л ИТЕРАТУPА
I. Cooley J. W., T u key J. W., An algorithm for the machine calculation of complex Fourier series, Math, of Computation, 19, 90, 297 (1965).
ПРИЛОЖЕНИЕ П7.2
ДИСПЕРСИЯ ОЦЕНОК НАКЛОНА
Дискретное время. Выборочная оценка наклона Oi в рассмотренной выше (разд. 7.4.2) модели
Yt ^ Q1X1+ Zt (П7.2.1)
Примеры одномерного спектрального анализа
67
имеет вид
n
Є, = S ЩУі, (П7.2.2)
t-i
где Wt — (2/N a) cos (nt/b). Действуя так же, как и при выводе формулы (5.2.9), получим дисперсию соответствующей оценки
n n
Var [O1] = 22 wtwryYY (t - г). (П7.2.3)
Если теперь предположим, что наблюдения производятся через единичные интервалы времени, то
V2
Vw (*)= jryr(f)e^kdf. (П7.2.4)
-V2
Подставляя (П7.2.4) в (П7.2.3), получаем
V2
Var [G1]= [ \W(f)\2VyY(f)df, (П7.2.5)
где
n
W (/) = 2 a^eW. (П7.2.6)
Теперь предположим, что wt = (2/Na) cos (nt/b), как в (7.4.1). Тогда
д7/А= 1 ( 1 -exp H(N+ 1)2я(/- 1/26)] . 1 -exp [/(Af+ 1)2я(/ + 1/26)] ) ^ ^ Afa 1 1-ехр[/2я(/-1/26)] I - ехр [;2я (/+ 1/26)] ( '
Подставляя квадрат модуля этого выражения в (П7.2.5), получаем
'/2
1 Г Г sin2 (N+ 1) я(/- 1/26)
L 1J Na2 J L Л7 sin2 я (/ — 1/26) ' -V2
і sin2 (Af+ 1)я (/ + 1/26)1 г /fuf 4 ЛГ sin2 я (/+1/26) J 1 УУ (Г) аГ
плюс члены с перекрестными произведениями.
Используя тот факт, что выражение в квадратных скобках стремится к
б(/-і) + б(/ + ^)
при N -* оо. а члены с перекрестными произведениями имеют порядок 1/N2, мы получаем
1іт/УУаг[в,] = -^ГикШ. (П7.2.7)
68
Глава 7
Для конечных значений отсюда следует
Var [вДъ^ГууШ. (П7.2.8)
TVa2 уу \2Ь
Поскольку для модели (П7.2.1) имеет место равенство ГУу(1/26) = rzz(l/26), то формула (П7.2.8) эквивалентна (7.4.2). Таким образом, дисперсия минимальна тогда, когда частота возмущения 1/2Ь равна частоте, на которой спектр шума достигает минимума.
Непрерывное время. В случае дискретного времени изменения характеристик процесса между регулировками были пренебрежимо малы. Если на входе имеется непрерывное синусоидальное возмущение х (t) = a cos 2nfot, то модель (П7.2.1) нужно видоизменить следующим образом:
Y (t) = 9,aG (Z0) cos 2пУ + Z (г),
где G (/о)—значение функции усиления на частоте /о. Дальнейший вывод вполне аналогичен случаю дискретного времени. Окончательный результат имеет вид
4 г (/) Var [в,] ~у^~§гщ-
Следовательно, дисперсия минимальна тогда, когда частота возмущения /о соответствует максимальному отношению сигнал/шум G2(fo)/IYy(/o).
ПРИЛОЖЕНИЕ П7.3
АЛГОРИТМ БЫСТРОГО ПРЕОБРАЗОВАНИЯ ФУРЬЕ
Недавним новшеством в спектральном анализе является алгоритм быстрого преобра зования Фурье (БПФ). С помощью этого алгоритма дискретное преобразование Фурье вычисляется гораздо быстрее, чем с помощью прямого метода, приведенного в разд. 2.1.2, и с той же самой точностью. Так, используя прямой метод для вычисления дискретного преобразования Фурье ряда из N членов, потребовалось бы приблизительно N2 операций, в то время как БПФ требует лишь 2N log2 А' операций. Экономия времени вычислений может быть очень велика, если нужно проводить анализ Фурье длинных рядов. Например, для вычисления с помощью БПФ коэффициентов Фурье ряда из N = 8192 членов [1] требовалось около 5 сек на вычислительной машине IBM 8094, в то время как для прямого метода нужно было около 30 мин.
