DSP Лекция 9 Digital Signal ProcessingDSP Лекция 9 Digital Signal Processing. - презентация

1 DSP Лекция 9 Digital Signal Processing

2 DSP Цифровой спектральный анализ Введение в алгоритмы БПФ с основанием 2 Введение в алгоритмы БПФ с основанием 2 Алгоритм БПФ с прореживанием по частоте Алгоритм БПФ с прореживанием по частоте Алгоритмы БПФ для составного значения N Алгоритмы БПФ для составного значения N Реализация ДПФ на основе цифровой фильтрации Реализация ДПФ на основе цифровой фильтрации Спектральные представления детерминированных сигналов Спектральные представления детерминированных сигналов Разрешение и произведение длительности на ширину спектра Разрешение и произведение длительности на ширину спектра Спектральные представления случайных процессов Спектральные представления случайных процессов Коррелограммный метод оценки СПМ Коррелограммный метод оценки СПМ Периодограммный метод оценки СПМ Периодограммный метод оценки СПМ Использование окон при спектральных измерениях Использование окон при спектральных измерениях Разрешение и произведение «устойчивость длительность ширина полосы» Разрешение и произведение «устойчивость длительность ширина полосы» Классические спектральные оценки тест-последовательности Классические спектральные оценки тест-последовательности

3 DSP Цифровой спектральный анализ Цифровой спектральный анализ (ЦСА) Спектральный анализ-это один из методов обработки сигналов, который позволяет охарактеризовать частотный состав измеряемого сигнала. Преобразование Фурье является математической основой, которая связывает временной или пространственный сигнал (или же некоторую модель этого сигнала) с его представлением в частотной области. Существенный вклад в развитие цифровых методов спектрального анализа внесли эффективные алгоритмы, предназначенные для вычисления дискретного преобразования Фурье, предложенные Д.Кули и Д.Тьюки в 1965 году. Набор алгоритмов, называемых алгоритмами быстрого преобразования Фурье (БПФ), включает разнообразные методы уменьшения времени вычисления дискретного преобразования Фурье (ДПФ). Поскольку вычисление ДПФ является основной операцией в большинстве задач спектрального анализа, то использование БПФ в некоторых встречающихся на практике случаях, позволяю­щее ускорить вычисление ДПФ в 100 и более раз по сравнению с методом прямого вычисления ДПФ, имеет чрезвычайно важное значение и должно рассматриваться как неотъемлемая часть применения методов цифровой обработки сигналов для спектрального анализа. Возможно, именно алгоритмы БПФ более чем какие - либо другие методы существенно расширили область применения методов спектрального анализа как средства обработки сигналов. Поэтому начнем рассмотрение вопросов ЦСА с алгоритмов БПФ, включающих алгоритмы с основанием 2 и прореживанием по времени и по частоте.

4 DSP Цифровой спектральный анализ Введение в алгоритмы БПФ с основанием 2 Напомним, что прямое ДПФ конечной последовательности x(n), 0 n N-1 определяется выражением или в более удобном виде как где а обратное ДПФ имеет вид (5.1) (5.2) Выражения (5.1) и (5.2) различаются только знаком экспоненты от W N и коэффициентом 1/N, поэтому рассуждения, касающиеся вычислительных процедур для (5.1), применимы с очевидными изменениями к (5.2).

5 DSP Из соотношения (5.1) следует, что в случае, когда последовательность х(п) является комплексной, при прямом вычислении N -точечного ДПФ нужно выполнить (N - I) 2 комплексных умножений и N (N 1) комплексных сложений. Таким образом, для достаточно больших N (порядка 1000) прямое вычисление ДПФ требует выполнения чрезмерного количества вычислительных операций. Основная идея БПФ состоит в том, чтобы разбить исходную N -точечную последовательность на две более короткие последовательности, ДПФ которых могут быть скомбинированы таким образом, чтобы получилось ДПФ исходной N -точечной последовательности. Так, например, если N четное, а исходная N -точечная последовательность разбита на две (N/2)- точечные последовательности, то для вычисления искомого N -точечного ДПФ потребуется порядка (N /2) 2 2 = N 2 /2 комплексных умножений, т. е. вдвое меньше по сравнению с прямым вычислением. Здесь множитель (N /2) 2 дает число умножений, необходимое для прямого вычисления (N /2)- точечного ДПФ, а множитель 2 соответствует двум ДПФ, которые должны быть вычислены. Эту операцию можно повторить, вычисляя вместо (N /2)- точечного ДПФ два (N /4)- точечных ДПФ (предполагая, что N/2 четное) и сокращая тем самым объем вычислений еще в два раза. Выигрыш в два раза является приближенным, поскольку не учитывается, каким образом из ДПФ меньшего размера образуется искомое N - точечное ДПФ. Цифровой спектральный анализ

6 DSP Проиллюстрируем описанную методику для N -точечной последовательности {х(п)}, считая, что N равно степени 2. Введем две (N /2)- точечные последовательности {х 1 (п)} и {х 2 (п)} из четных и нечетных членов х(п) соответственно, т. е. (5.3) Цифровой спектральный анализ N -точечное ДПФ последовательности {x(n)} можно записать как (5.4) где учтено, что, и введены обозначения для (N/2)- точечных ДПФ последовательностей x 1 (n) и x 2 (n):

7 DSP Цифровой спектральный анализ Из формулы (5.4) следует, что N -точечное ДПФ X (k) может быть разложено на два (N /2)- точечных ДПФ, результаты которых объединяются согласно (5.4). Если бы ( N /2 )-точечные ДПФ вычислялись обычным способом, то для вычисления N -точечного ДПФ потребовалось бы, очевидно, (N 2 /2 + N) комплексных умножений. При больших N (когда N 2 /2 N ) это позволяет сократить время вычисления на 50%. Поскольку X (k) определено при 0 k N 1, а Х 1 (k) и Х 2 (k) определены при 0 k N/2-1, необходимо доопределить формулу (5.4) для k N/2. Это определение достаточно очевидно и может быть записано следующим образом: (5.5) где учтено, что. Таким образом, вычисление X (k) по X 1 (k) и X 2 (k) можно представить в виде (5.6)

8 DSP На рис.1 с помощью направленного графа представлена последовательность операций при вычислении восьмиточечного ДПФ с использованием двух четырехточечных ДПФ. Незачерненный кружок графа означает операцию сложения/вычитания, причем верхний выход соответствует сумме, а нижний – разности. Стрелка обозначает операцию умножения на значение множителя а, указанного над стрелкой. Входная последовательность х (п) сначала разбивается на две последовательности х 1 (п) и х 2 (п) из четных и нечетных членов х(п), после чего рассчитываются их преобразования Х 1 (k) и Х 2 (k). Затем в соответствии с формулой (5.6) получают X (k). Цифровой спектральный анализ Рисунок 1. Вычисление восьмиточечного ДПФ через два четырехточечных ДПФ.

9 DSP Выражение (5.6) соответствует разбиению исходного N –точечного вычисления ДПФ на два N /2 -точечных вычислений. Если N /2 –четное число, что имеет место всегда, когда N равно степени 2, то можно вычислять каждое N /2 -точечное ДПФ в (5.6) путем разбиения сумм на два N /4 –точечных ДПФ, которые затем объединяются, давая N /2 –точечное ДПФ. Каждая из последовательностей х 1 (п) и х 2 (п) разбивается на две последовательности, состоящие из четных и нечетных членов. Аналогично N/2 -точечные ДПФ могут быть записаны как комбинации двух N /4- точечных ДПФ, т. е. Цифровой спектральный анализ или (5.7) (5.8) где 0 < k < N/4-1, A(k) и В(k) - N/4 -точечные ДПФ соответственно четных и нечетных членов х 1 (п), C(k) и D(k) - N/4 -точечные ДПФ соответственно четных и нечетных членов х 2 (п). На рис.2 показан результирующий направленный граф, в котором четырехточечные ДПФ из рис.1 рассчитываются согласно (5.8).

10 DSP Процесс уменьшения размера ДПФ от L до L/2, где L равно степени 2, может быть продолжен до тех пор, пока не останутся только двухточечные ДПФ. Двухточечное ДПФ, например, A(k), k = 0,1,может быть рассчитано без использования умножений по формулам (5.9)(5.9) Цифровой спектральный анализ Таким образом, восьмиточечное ДПФ (рис.1 и 2) в итоге сводится к алгоритму, описываемому направленным графом, представленным на рис.3. Анализ графа на рис.3 и процедуры последовательного сокращения вдвое размеров преобразований показывает, что на каждом этапе БПФ (т. е. при каждом сокращении размеров ДПФ) необходимо выполнить N/2 комплексных умножений. Поскольку общее количество этапов равно log 2 N, то число комплексных умножений, необходимое для нахождения N -точечного ДПФ, приблизительно равно N/2log 2 N. Слово приблизительно использовано по той причине, что умножения на,, и в действительности сводятся просто к сложениям и вычитаниям комплексных чисел. Так, например, на рис.3 первый этап БПФ содержит только сложения и вычитания комплексных чисел. Даже на втором этапе используются только сложения и вычитания комплексных чисел. Фактически, как следует из направленного графа на рис.3, вместо ожидаемых 12 (т. е. 4log 2 N ) достаточно выполнить всего два нетривиальных умножения. Однако для больших значений N фактическое число нетривиальных умножений хорошо аппроксимируется выражением N /2log 2 N.

11 DSP Цифровой спектральный анализ Рисунок 2. Вычисление восьмиточечного ДПФ через два четырехточечных ДПФ, которые в свою очередь вычисляются через четыре двухточечных ДПФ. Рисунок 3. Восьмиточечное ДПФ, полученное последовательным прореживанием в 2 раза.

12 DSP Описанный выше алгоритм был назван алгоритмом с прореживанием по времени, поскольку на каждом этапе входная (т.е. временная) последовательность разделяется на две обрабатываемые последовательности меньшей длины, т. е. входная последовательность прореживается на каждом этапе. Другая форма алгоритма БПФ (с прореживанием по частоте) будет описана ниже, а сейчас целесообразно обсудить некоторые общие свойства алгоритмов БПФ. Базовая операция алгоритма с прореживанием по времени (так называемая «бабочка») состоит в том, что два входных числа А и В объединяются для получения двух выходных чисел X и Y следующим образом: (5.10) Цифровой спектральный анализ На рис.4 изображен направленный граф базовой операции (5.10). Рисунок 4. Базовая операция алгоритма БПФ.

13 DSP Внимательное рассмотрение направленного графа на рис.3 показывает, что каждый из этапов содержит N/2 базовых операций. В случае когда множитель нетривиальный, для каждой базовой операции необходимо выполнить только одно умножение, поскольку величину можно вычислить и запомнить. Таким образом, структура базовых операций такова, что для выполнения БПФ N -точечной последовательности, размещенной в памяти, достаточно иметь лишь одну дополнительную ячейку памяти. Результаты всех промежуточных этапов БПФ можно размещать в те же ячейки памяти, где находились исходные данные. Поэтому для хранения и входной, и выходной последовательностей можно использовать один и тот же массив ячеек памяти. Алгоритм, в котором для размещения входной и выходной последовательностей используются одни и те же ячейки памяти, называется алгоритмом БПФ с замещением. Еще одной особенностью алгоритма с прореживанием по времени (как, впрочем, и большинства других алгоритмов БПФ) является необходимость такой перестановки элементов входной последовательности, чтобы выходная последовательность X (k) имела естественный (прямой) порядок расположения, т. е. k =0,1,..., N-1. В примере на рис.3 для этого требовался следующий порядок размещения входной последовательности: х(0), х(4), х(2), х(6), х(1), х(5), х(3) и х(7). Характер перестановки элементов входной последовательности может быть описан сравнительно просто. Цифровой спектральный анализ

14 DSP Ниже будет показано, что в случае, когда N является степенью 2, входная последовательность должна быть расположена в памяти в двоично-инверсном порядке для того, чтобы выходная последовательность получалась в прямом порядке. Двоично-инверсный порядок определяется следующим образом. Если записать порядковые номера элементов входной последовательности в двоичном коде, используя L двоичных разрядов, причем N = 2 L, а затем инвертировать порядок следования разрядов, то получаемые при этом числа и будут номерами элементов вход­ной последовательности после их перестановки. Так, для случая N = 8 = 2 3 прямой порядок номеров приведен в табл. 5.1 слева, а двоично-инверсный порядок - справа. Таким образом, для двоичной инверсии входной последовательности необходим соответствующий алгоритм. Простой и эффективный двоично-инверсный счетчик предложен Рейдером. Цифровой спектральный анализ Номер Двоичное представление Двоичная инверсия Двоично- инверсный номер Таблица 5.1.

15 DSP Из сказанного выше ясно, что перестановку входной последовательности можно произвести с замещением, меняя в парах местами числа с прямым и двоично-инверсным номерами и используя для этого лишь одну вспомогательную ячейку памяти. На рис.5 показана схема перестановки данных, представленных в табл Цифровой спектральный анализ Рисунок 5. Перестановка данных с замещением. Простая, составленная на ФОРТРАНе, программа расчета БПФ с прореживанием по времени (авторы ее Кули, Льюис, Уэлч) приведена на рис.6.

16 DSP SUBROUTINE FFT(A,M,N) COMPLEX A(N),U,W,T N=2**M NV2=N/2 NM1=N-1 J=1 DO 7 I=1,NM1 IF(I.GE. J) GO TO 5 T=A(J) A(J)=A(I) A(I)=T 5K=NV2 6IF(K.GE. J) GO TO 7 J=J-K K=K/2 GO TO 6 7J=J+K PI= DO 20 L=1,M LE=2**L LE1=LE/2 U=(1.0,0.) W=CMPLX(COS(PI/LE1),SIN(PI/LE1)) DO 20 J=1,LE1 DO 10 I=J,N,LE IP=I+LE1 T=A(IP)*U A(IP)=A(I)-T 10A(I)=A(I)+T 20U=U*W RETURN END Цифровой спектральный анализ Рисунок 6. Программа БПФ с основанием 2 и прореживанием по времени.

17 DSP Входная комплексная последовательность представляет собой массив А, размером N= 2 м, причем М указывается в операторе обращения к подпрограмме. Все операторы цикла DO 7 предназначены для выполнения двоично-инверсной перестановки элементов входного массива. Остальные операторы используются для вычисления непосредственно БПФ и образуют три вложенных цикла. С помощью первого (внешнего) цикла выполняется М этапов, другой цикл предназначен для выполнения базовых операций в пределах каждого этапа, а третий цикл (внутренний) необходим для вычисления степеней W, используемых при выполнении базовых операций в пределах одного этапа. Для вычисления степеней W используется рекуррентное соотношение, причем синусы и косинусы непосредственно вычисляются лишь для нахождения W по приращениям. Цифровой спектральный анализ

18 DSP Цифровой спектральный анализ Алгоритм БПФ с прореживанием по частоте Другая распространенная форма алгоритма БПФ (при условии, что N равно степени 2) - так называемый алгоритм БПФ с прореживанием по частоте. В этом варианте алгоритма БПФ входная последовательность {х(п)} разбивается на две последовательности, содержащие по N/2 отсчетов каждая следующим образом: первая последовательность {х 1 (п)} состоит из первых (N/2) отсчетов {х(п)}, а вторая {х 2 (п)} из остальных (N/2) отсчетов {х(п)}, т, е (5.11) При таком разбиении N -точечное ДПФ последовательности х(п) можно записать в виде где учтено, что Запишем выражения раздельно для четных и нечетных отсчетов ДПФ: (5.12) (5.13)

19 DSP Из выражений (5.12) и (5.13) видно, что четные и нечетные отсчеты ДПФ можно получить из (N /2)-точечных ДПФ последовательностей f(n) и g(n), равных Цифровой спектральный анализ (5.14) Таким образом, снова вычисление N -точечного ДПФ удалось свести к вычислению двух (N /2)- точечных ДПФ. На рис. 7 эта методика иллюстрируется для случая N = 8. Рисунок 7. Переход от восьмиточечного ДПФ к двум четырехточечным ДПФ при прореживании по частоте.

20 DSP Описанную методику можно применить повторно и представить каждое из (N/2)- точечных ДПФ в виде комбинации двух (N/4)- точечных ДПФ. На рис. 8 и 9 показан переход от четырех­точечных ДПФ (рис.9) к двухточечным ДПФ с последующим прямым вычислением двухточечных ДПФ. Цифровой спектральный анализ Рисунок 8. Переход от четырехточечных ДПФ на рис.7 к двухточечным ДПФ. Рисунок 9. Полный направленный граф восьмиточечного ДПФ с замещением и прореживанием по частоте.

21 DSP Сравнение алгоритмов, иллюстрированных на рис. 5.3 и 5.9, позволяет выявить два очевидных различия между ними. Во-первых, при прореживании по времени порядок следования входных отсчетов двоично-инверсный, а выходных прямой и наоборот при прореживании по частоте (рис.9). Второе отличие заключается в несколько ином выполнении базовой операции (см. рис.10 и 4): при прореживании по частоте комплексное умножение выполняется после сложения - вычитания. Цифровой спектральный анализ Рисунок 10. Базовая операция алгоритма БПФ с прореживанием по частоте. Легко заметить и сходство между алгоритмами с прореживанием по времени и по частоте. В обоих случаях при вычислении ДПФ требуется около N Iog 2 N операций, вычисления могут быть проведены с замещением и должно быть предусмотрено выполнение двоичной инверсии. Отметим еще одно сходство между вычислениями, соответствующим базовым операциям («бабочкам»), и между рис.3 и 9 т.е. то, что рис.9 можно получить из рис.3 изменив направление стрелок (сигналов) и поменяв вход и выход. Как уже отмечалось, для вычисления обратного ДПФ (5.2) можно использовать БПФ-алгоритм, если разделить результат на N и использовать вместо степеней степени. В приложении 5.1 приведена программа на ФОРТРАНе, реализующая вычисления прямого и обратного ДПФ по основанию 2 с прореживанием по частоте.

22 DSP Приложение 5.1 SUBROUTINE FFT2(X,N,INV) C программа быстрого преобразования Фурье CX- комплексный массив отсчетов функции или спектра CN- размер массива X, равный целой степени двойки CINV- указатель преобразования (=0-прямое, =1-обратное) C COMPLEX X(1),W,T,CMPLX ITER=0 IREM=N 10IREM=IREM/2 IF(IREM.EQ.0) GO TO 20 ITER=ITER+1 GO TO 10 20CONTINUE SIGN=-1. IF(INV.EQ.1) SIGN=1. NXP2=N DO 50 IT=1,ITER NXP=NXP2 NXP2=NXP/2 DEL=FLOAT(NXP2) WPWR= /DEL DO 40 M=1,NXP2 ARG=FLOAT(M-1)*WPWR W=CMPLX(COS(ARG),SIGN*SIN(ARG)) Цифровой спектральный анализ

23 DSP Приложение 5.1 DO 40 MXP=NXP,N,NXP J1=MXP-NXP+M J2=J1+NXP2 T=X(J1)-X(J2) X(J1)=X(J1)+X(J2) 40X(J2)=T*W 50CONTINUE N2=N/2 N1=N-1 J=1 DO 65 I=1,N1 IF(I.GE.J) GO TO 55 T=X(J) X(J)=X(I) X(I)=T 55K=N2 60IF(K.GE.J) GO TO 65 J=J-K K=K/2 GO TO 60 65J=J+K IF(INV.NE.1) GO TO 75 DO 70 I=1,N 70X(I)=X(I)/FLOAT(N) 75CONTINUE RETURN END Цифровой спектральный анализ

24 DSP Цифровой спектральный анализ Алгоритмы БПФ для составного значения N Рассмотренные выше алгоритмы предполагали, что N является целой степенью 2, т.е. N=2 L. В более общем случае эффективное вычисление ДПФ связано с представлением N в виде сомножителей: (5.15) Если N представлено в виде произведения одинаковых сомножителей r, то соответствующий алгоритм называют алгоритмом БПФ по основанию r. Понятие смешанное основание означает, что не все сомножители N одинаковы. N -точечные алгоритмы по основанию 2 особенно просты для реализации, поэтому на практике выгодно всегда иметь дело с последовательностями длины N=2 L. Это можно сделать во многих случаях, просто дополняя последовательность конечной длины нулями, если это необходимо. В случаях, когда невозможно выбрать N=2 L, приходится рассматривать более общую ситуацию (5.15). Рассмотрим применение принципа прореживания по времени, когда N является произведением сомножителей не все из которых равны 2. Пусть q 1 =p 2 p 3 …p L, так что N=p 1 q 1. Если N является степенью 2, можно выбрать p 1 =2, а q 1 = N/2. Используя прореживание по времени, можно разложить x(n) на две последовательности длины N/2, состоящие из четных и нечетных выборок соответственно, как это мы делали выше. Если N=p 1 q 1, то можно разделить входную последовательность на p 1 последовательностей длины q 1, так, что каждая p 1 -я выборка попадает в одну последовательность.

25 DSP Например, если p 1 =3, а q 1 = 4, так что N=12, то можно разложить x(n) на три последовательности длины 4, причем первая последовательность состоит из выборок x(0), x(3), x(6), x(9), вторая – из x(1), x(4), x(7), x(10), а третья – из x(2), x(5), x(8), x(11). В общем случае можно записать X(k) в виде Цифровой спектральный анализ (5.16) Внутренние суммы можно представить как q 1 -точечные ДПФ: потому что, как легко проверить,, при N=p 1 q 1. Таким образом, (5.16) представляет X(k) в виде p 1 ДПФ последовательностей длины q 1. Чтобы определить число комплексных умножений и сложений для вычисления ДПФ по (5.16), будем считать, что q 1 -точечные ДПФ получаются путем прямого вычисления. Из (5.16) видно, что нужно рассчитать p 1 q 1 -точечных ДПФ. Поэтому общее число требуемых комплексных сложений и умножений равно p 1 q 1 2.

26 DSP Внешняя сумма в (5.16) получается путем умножения q 1 -точечных ДПФ на коэффициенты и сложения результатов. Так как двойное суммирование в (5.16) выполняется для N значений k, то для объединения p 1 q 1 -точечных ДПФ требуется N(p 1 -1) комплексных сложений и умножений. Следовательно, общее число комплексных умножений и сложений для вычисления ДПФ по (5.16) равно N(p 1 -1)+ p 1 q 1 2. Теперь q 1 -точечное ДПФ может быть разложено аналогичным образом. В частности, если представить q 1 в виде q 1 =p 2 q 2, то q 1 -точечные последовательности во внутренней сумме (5.16) могут быть разбиты на p 2 последовательностей, каждая из которых состоит из q 2 точек, так, что внутренняя сумма в (5.16) может быть заменена на двойную сумму тем же способом, с которого мы начали. Тогда число операций требуемых для расчета q 1 -точечных ДПФ в (5.16) вместо q 1 2 станет равным q 1 (p 2 -1)+ p 2 q 2 2. Следовательно, общее число комплексных умножений и сложений для вычисления ДПФ станет равным Цифровой спектральный анализ N(p 1 -1)+ N(p 2 -1)+ p 1 p 2 q 2 2. Если продолжить эту процедуру, разлагая далее q 2 -точечные ДПФ, то, в конце концов, общее число комплексных умножений и сложений для вычисления ДПФ станет равным N(p 1 + p 2 + …+p L –L). (5.17) Из (5.17) видно, что лучше производить разложение на максимально возможное число сомножителей и лучше выбирать простые сомножители.

27 DSP Цифровой спектральный анализ Реализация ДПФ на основе цифровой фильтрации Покажем, что ДПФ можно рассматривать как отклик цифрового фильтра и найдем структуру фильтра, определяющего один спектральный отсчет X(k) в точке k =2 k/N : (5.18) В такой форме записи (5.18) можно рассматривать X(k) как выходной сигнал y k (N-1) фильтра с импульсной характеристикой h k (n) в момент (N-1). Из (5.18) импульсная характеристика h k (n ) определится соотношением: (5.19) Таким образом, цифровой фильтр – это КИХ-фильтр (комплексный), передаточная функция которого имеет вид: (5.20) Эта структура (5.20) уже встречалась – один канал фильтра с частотной выборкой. Для вычисления всех значений X(k) (k=0, 1, …, N-1) необходимо N таких фильтров – частотных каналов анализатора спектра. Каждый частотный канал анализатора спектра представляет собой комплексный резонатор – когерентный накопитель отсчетов комплексной гармоники

28 DSP На рис. 11 представлена структурная схема анализатора спектра в виде гребенки фильтров. Достоинством такого выполнения является возможность вычисления «скользящего» X n (k) спектра по N предшествующим текущему моменту n отсчетам: x(n), x(n-1), x(n-2), …, x(n-N+1). Заметим, что при скользящем анализе на основе БПФ для каждого нового входного отсчета приходится вычислять полное БПФ, т.е. N/2 log 2 N базовых операций на отсчет. При анализе с помощью гребенки фильтров (рис.11) число комплексных умножений на отсчет равно N. Таким образом, в большинстве представляющих интерес случаев гребенка фильтров оказывается эффективнее скользящего БПФ. Определим амплитудно-частотную характеристику одного канала анализатора спектра, вычислив H k (z) на единичной окружности: Цифровой спектральный анализ

29 DSP Цифровой спектральный анализ Рисунок 11. Анализатор спектра в виде гребенки фильтров. X n (0) X n (1) X n (k) X n (N-1)

30 DSP На рис.12 изображены графики функций D N ( - k ) для 16 –точечного ДПФ, реализованного гребенкой фильтров в соответствии с рис.11. Частотные характеристики фильтров с четными номерами k показаны на рис.12 вверху, а с нечетными номерами k – внизу, при этом для всех фильтров, за исключением 8-го, изображены лишь главные лепестки. Видно, что скользящее ДПФ эквивалентно довольно грубому набору фильтров с АЧХ, имеющими относительно большие боковые лепестки и существенным перекрытием между соседними фильтрами. При поступлении на вход анализатора спектра комплексной гармоники с частотой = k =2 k/N, все отсчеты X(i) 0, (i k), а X(k) =N. Если же 2 k/N, то отличны от нуля все X(i), т.е. происходит «растекание» (просачивание) энергии входной гармоники по всем N спектральным отсчетам. Цифровой спектральный анализ Рисунок 12. Анализатор спектра в виде гребенки фильтров.

31 DSP Чтобы улучшить частотную характеристику канальных фильтров анализатора можно ввести весовую обработку входной последовательности x(n) умножением ее на весовую функцию (окно) w(n). Аналогично проделанной выше процедуре, можно показать, что АЧХ канальных фильтров с весовой обработкой совпадает с амплитудно- частотным спектром окна, сдвинутого на частоту k. Выбор функций окна обсудим ниже в специальном разделе. На основе цифровых фильтров предложен более экономный по числу операций, чем прямой алгоритм вычисления X(k), получивший название алгоритма Герцеля. Заметим, что в этих алгоритмах нет необходимости вычислять все N значений X(k), т.е. в общем случае можно вычислять X(k) для любых M значений k. Поэтому, при M

32 DSP Цифровой спектральный анализ Спектральные представления детерминированных сигналов Напомним определения понятий, связанных со спектральными представлениями детерминированных сигналов непрерывного и дискретного времен. Для непрерывных во времени детерминированных сигналов понятие спектральной плотности амплитуды связано с парой преобразований Фурье (НВПФ): (5.21) Энергия E сигнала определяется соотношением где спектральная плотность энергии (СПЭ), описывающая распределение энергии сигнала по частоте. Для сигналов с ограниченной мощностью вводят спектральную плотность мощности, описывающую распределение мощности (энергии, рассеиваемой в единицу времени) по частоте.

33 DSP В случае, когда мы располагаем лишь отсчетами (выборками) непрерывной функции времени взятыми с периодом T, понятие спектральной плотности амплитуды последовательности отсчетов определяется парой дискретно-временного преобразования Фурье (ДВПФ) Цифровой спектральный анализ (5.22) Заметим, что определяет периодическую ( с периодом 1/T ) функцию частоты, которая является периодическим продолжением спектральной плотности сигнала. Если имеет ограниченную локализацию по оси частот, например, в интервале, а период отсчетов удовлетворяет соотношению то спектральные полосы функции не перекрываются, поэтому Следовательно, спектр дискретной последовательности повторяет форму спектра функции, которая может быть восстановлена по ее отсчетам. В противном случае имеет место эффект наложения спектров и восстановление функции по ее отсчетам оказывается невозможным.

34 DSP Цифровой спектральный анализ Заметим также, что связано с z-преобразованием последовательности Поэтому пару соотношений ДВПФ можно записать с введением шкалы частот, нормированной к частоте дискретизации в виде (5.23) Для сигнала, ограниченного по длительности интервалом 0-NT секунд, спектральные представления иногда связывают с непрерывно-временным рядом Фурье (НВРФ): (5.24) который по существу на всей временной оси представляет периодически (с периодом NT ) продолженную функцию, т.е. Для а т.е. коэффициенты ряда Фурье суть отсчеты непрерывного по частоте спектра в дискретной последовательности точек на оси частот.

35 DSP Наконец для последовательности конечной длины в N отсчетов спектральные представления связываются с дискретно - временным рядом Фурье (ДВРФ): Цифровой спектральный анализ (5.25) Обе последовательности и периодичны с периодом N и на интервале в один период представляют соответственно совокупность отсчетов спектра в точках и отсчетов сигнала в точках. Соотношения (5.25) называют еще парой дискретного преобразования Фурье (ДПФ), их отличия от традиционно используемых, вызванные наличием множителей T и 1/T, связаны с желанием обеспечить корректность масштабов при вычислении энергии и мощности. Теорема об энергии для ДПФ имеет вид Поэтому спектральная плотность энергии примет вид а соотношение для спектральной плотности мощности (СПМ) примет форму (5.26) (5.27)

36 DSP Цифровой спектральный анализ Таким образом, ДПФ (5.25) позволяет преобразовать N отсчетов непрерывного во времени сигнала со спектральной плотностью в N отсчетов его спектральной плотности. Для получения большего количества отсчетов спектра L>N, т.е. для частот, следует доопределить исходную последовательность, заданную на интервале, необходимым числом нулевых отсчетов, полагая,. Алгоритм ДПФ, реализуемый одним из вариантов БПФ, является основным при реализации спектрального анализа сигналов на основе цифровых вычислительных устройств. Рис.13 показывает как влияет дополнение нулями исходной 16-точечной последовательности данных, содержащих три синусоиды, с соответствующим увеличением объема ДПФ на воспроизведение ее спектра.

37 DSP Цифровой спектральный анализ Рисунок 13. Интерполяция за счет дополнения нулями: а - модуль ДВПФ 16-точечной записи данных, содержащих три синусоиды; б - модуль ДПФ той же последовательности данных без дополнения нулями (неопределенности не разрешены); в - модуль ДПФ той же последовательности после двукратного увеличения ее отсчетов за счет дополнения нулями (неопределенности разрешены, так как различимы все три синусоиды); г -модуль ДПФ той же последовательности данных после восьмикратного увеличения числа ее отсчетов за счет дополнения нулями.

38 DSP Цифровой спектральный анализ Разрешение и произведение длительности на ширину спектра Одним из важных параметров устройства спектрального оценивания является его разрешающая способность, под которой понимается способность разрешать (раздельно измерять) спектральные отклики двух синусоидальных сигналов, близких по частоте и амплитуде. Считается, что спектральные пики разрешены, если провал между ними составляет величину не менее 3 дб. Обычно полагают, что разнесение частот двух синусоид не может быть меньше эквивалентной ширины полосы частот спектра окна, через которое наблюдаются отрезки этих синусоид. Поэтому говорят, что разрешение в герцах приближенно равно величине, обратной времени наблюдения в секундах, поскольку. Обосновать это утверждение можно на основе следующих соображений. Сигнал x(n) не может быть одновременно быть ограниченным по длительности и по ширине полосы спектра. Однако его можно охарактеризовать некоторым интервалом T e секунд, в котором сосредоточена большая часть его энергии при представлении во временной области, и некоторым интервалом B e герц, в котором сосредоточена большая часть его энергии при представлении в частотной области. Для количественного описания временной концентрации энергии дискретно-временной последовательности отсчетов сигнала и соответствующей частотной концентрации ее ДВПФ предложены несколько различных мер. Определим две из них.

39 DSP Эквивалентная длительность T e сигнала x(n) определяется как Цифровой спектральный анализ (5.28) т. е. иными словами, как «площадь» этого дискретного сигнала, поделенная на его центральное значение. Эквивалентная длительность сигнала равна длительности сигнала с прямоугольной огибающей, высота которого равна значению x(n) в начале координат, т. е. x(0), а площадь равна площади исходного сигнала; см. рис.14а и 14б. Заметим, что в общем случае Т е не будет целочисленно кратным интервалу отсчетов Т. Эквивалентная ширина полосы В е дискретно-временного преобразования Фурье X(f) сигнала x(n) определяется аналогично как (5.29) см. рис.14 в и 14 г. Эти две меры временной и частотной концентраций применимы лишь к действительнозначным симметричным сигналам с максимальным значением в начале координат, интеграл от которых имеет конечное ненулевое значение. Эти условия точно выполняются для весовых функций (окон), которые рассматриваются в связи с обсуждением классических спектральных оценок.

40 DSP Цифровой спектральный анализ Рисунок 14. а – последовательность отсчетов исходного сигнала; б – эквивалентный сигнал с прямоугольной огибающей; в – ДВПФ исходного сигнала; г – эквивалентное ДВПФ прямоугольной формы.

41 DSP Из соотношений ДВПФ с очевидностью следует, что Цифровой спектральный анализ откуда получаем, что произведение длительности и ширины полосы равно (5.30) а это означает, что эквивалентная длительность сигнала и эквивалентная ширина его преобразования являются взаимно обратными величинами. Именно это равное единице произведение, в котором Т е полагается равным интервалу наблюдения, и положено в основу эмпирического правила, определяющего частотное разрешение. Возможно другое определение произведения длительности и ширины полосы, основанное на понятии среднеквадратичной длительности (ширины). Среднеквадратичная ширина является мерой дисперсии (среднеквадратичного отклонения) некоторой функции от ее среднего значения.

42 DSP Для того чтобы получить возможность оперировать с более широким классом функций, чем это возможно при использовании эквивалентных длительности (5.28) и ширины полосы (5.29), можно использовать определения среднеквадратичной ширины для функций вида |х(n)| 2 и X(f) 2. Такой подход позволяет получать значимые среднеквадратичные величины ширины для комплексных функций, осциллирующих функций и функций, которые имеют нулевую интегральную площадь. На рис. 15 – 18 показано, как влияет на величину частотного разрешения длина N окна данных на примере вычисления ДПФ объемом 512 точек суммы двух комплексных синусоид с равными амплитудами и близкими частотами 0,2 и 0,24. Длина N окна данных выбиралась равной 10, 30, 100 и 300 отсчетов соответственно и обрабатываемая реализация дополнялась необходимым количеством нулевых отсчетов до длины 512 отсчетов, весовая обработка не применялась. На рис.15 синусоиды не разрешены, поскольку разнос их частот составляет величину f=0,04, а частотное разрешение 1/N=1/10. На рис.16 синусоиды уже разрешены, поскольку частотное разрешение 1/N=1/ Рис.17 и 18 иллюстрируют спектры разрешенных синусоид при N=100 и N=300 соответственно. Цифровой спектральный анализ

43 DSP Цифровой спектральный анализ Рисунок 15. Спектр суммы двух комплексных синусоид с равными амплитудами и близкими частотами 0,2 и 0,24, N=10.

44 DSP Цифровой спектральный анализ Рисунок 16. Спектр суммы двух комплексных синусоид с равными амплитудами и близкими частотами 0,2 и 0,24, N=30.

45 DSP Цифровой спектральный анализ Рисунок 17. Спектр суммы двух комплексных синусоид с равными амплитудами и близкими частотами 0,2 и 0,24, N=100.

46 DSP Цифровой спектральный анализ Рисунок 18. Спектр суммы двух комплексных синусоид с равными амплитудами и близкими частотами 0,2 и 0,24, N=300.

47 DSP Цифровой спектральный анализ Спектральные представления случайных процессов Дискретный случайный процесс можно рассматривать как некоторую совокупность или ансамбль действительных или комплексных дискретных временных последовательностей, каждую из которых можно было бы наблюдать как результат проведения некоторого эксперимента. При фиксации временного индекса в сечениях случайного процесса образуются случайные величины, при задании распределений плотностей вероятностей которых можно определить соответствующие моменты распределений. Среднее значение (момент 1-го порядка) определяется выражением Автокорреляция случайного процесса в два различных момента времени и определяется выражением где * - знак комплексного сопряжения. Автокорреляция центрированного случайного процесса, т.е. с удаленным средним значением, называется автоковариацией:

48 DSP Случайный процесс называется стационарным в широком смысле, если его среднее значение не зависит от времени, а автокорреляционная последовательность (АКП) зависит только от разности индексов времени m=n 1 -n 2. Спектральная плотность мощности (СПМ) стационарного случайного процесса определяется как дискретно - временное преобразование Фурье автокорреляционной последовательности Цифровой спектральный анализ (5.31) СПМ дискретного случайного процесса периодически (с периодом 1/T Гц) размножает по оси частот СПМ исходного непрерывного во времени процесса в предположении, что она локализована в интервале частот Обратное ДВПФ имеет вид (5.32) поэтому при m=0 имеем

49 DSP Цифровой спектральный анализ Последнее соотношение оправдывает трактовку как распределение средней мощности процесса по частоте. Пару преобразований Фурье (5.31) и (5.32) называют соотношениями Винера - Хинчина для случая дискретного времени. Поскольку то СПМ должна быть действительной и неотрицательной функцией. Если АКП – действительная функция, то и СПМ можно представить в форме (5.33) а это означает четность СПМ как функции частоты, т.е. P xx (f)= P xx (-f). Если стационарный случайный процесс обладает свойством эргодичности, то процедуру усреднения по ансамблю при вычислении среднего значения и автокорреляции можно заменить усреднением по времени с использованием одной единственной реализации процесса. Допущение об эргодичности позволяет ввести еще одно определение СПМ: (5.34)

50 DSP Цифровой спектральный анализ Эта эквивалентная (5.32) форма СПМ получается путем статистического усреднения квадрата модуля ДВПФ взвешенной совокупности данных, поделенного на длину записи данных, при увеличении числа отсчетов данных до бесконечности. Метод определения СПМ на основе оценки автокорреляционной последовательности называют косвенным (коррелограммным), поскольку случайный процесс непосредственно не используется для оценки СПМ. Определение СПМ по формуле (5.34) называют прямым методом (периодограммным), т.к. процесс непосредственно используется в расчете СПМ.

51 DSP Цифровой спектральный анализ Коррелограммный метод оценки СПМ. Рассмотрим подробнее практическую реализацию коррелограммного метода оценивания СПМ. Автокорреляционная последовательность эргодического процесса как предел среднего по времени определяется соотношением (5.35) На практике эта последовательность, как правило, не известна и поэтому должна оцениваться по имеющейся конечной последовательности данных. Предположим, что мы располагаем N отсчетами данных, тогда одним из вариантов соотношения для оценки АКП может быть, например, следующее: (5.36) применимое при При отрицательных значениях индекса значения АКП получаются из свойства сопряженной симметрии. Таким образом, (2N-1) значений АКП могут быть оценены по N отсчетам данных. Дискретная последовательность является несмещенной оценкой истинной АКП, поскольку

52 DSP При увеличении N значение дисперсии оценки стремится к нулю, поэтому это статистически состоятельная оценка АКП. Цифровой спектральный анализ Другой вариант соотношения для оценки АКП имеет вид (5.37) Эта оценка отличается от только нормирующим множителем поэтому при конечном N она является смещенной, поскольку Для типичных приложений средний квадрат ошибки будет, как правило, больше для оценки, чем для оценки. По этой и ряду других причин предпочтение часто отдается смещенной оценке АКП. Суммирование, предусматриваемое в выражениях (5.36) и (5.37) можно записать как линейную свертку Следовательно, для вычисления оценок дискретной автокорреляции можно применить эффективные в вычислительном отношении процедуры, такие как «быстрая свертка», в которых используются алгоритмы БПФ.

53 DSP Цифровой спектральный анализ Коррелограммный метод оценивания СПМ есть подстановка в соотношение (5.31) конечной последовательности значений оценки автокорреляции (коррелограммы) на место бесконечной последовательности неизвестных истинных значений АКП. Так, при подстановке несмещенной оценки АКП с максимальными значениями временного сдвига получается оценка СПМ в виде (5.38) Максимальный индекс временного сдвига L обычно ограничивают величиной. Причина выбора такого максимального значения – стремление устранить большие значения дисперсии, связанные с оценками автокорреляции при больших временных сдвигах. Среднее значение оценки определяется соотношением (5.39) где прямоугольная функция (окно) единичной амплитуды

54 DSP Цифровой спектральный анализ Соотношение (5.39) есть ДВПФ произведения последовательностей АКП и, которое выражается через свертку ДВПФ сомножителей: (5.40) где ДВПФ прямоугольного окна Таким образом среднее значение оценки отличается от истинного спектра, поскольку это свертка истинного спектра и спектра прямоугольного окна. Для уменьшения эффекта просачивания (растекания) спектра из-за влияния прямоугольного окна и уменьшения смещения оценки в соотношение (5.38) следует ввести весовую функцию (корреляционное окно),, симметричное относительно начала отсчета. Тогда оценка СПМ, получаемая на основе коррелограммного метода, принимает форму (5.41) Среднее значение этой оценки будет сверткой истинного спектра и спектра окна :

55 DSP Цифровой спектральный анализ Правильный выбор окна позволит уменьшить растекание спектра и его смещение. Для вычисления оценки СПМ, определяемой на сетке из K частот f i =i/KT, где 0 i K-1, можно использовать алгоритм БПФ. Обычно значение K=2 M >>L, а это требует доопределения, при L

56 DSP Цифровой спектральный анализ Периодограммный метод оценки СПМ. Как уже отмечалось выше, для эргодического случайного процесса эквивалентным определением СПМ является следующая дискретно- временная форма: При наличии конечного множества данных, и единственной реализации это соотношение преобразуется в СПМ выборки или периодограмму: (5.42) Несложно показать, что среднее значение оценки такое же, как и у коррелограммной оценки, дисперсия оценки не стремится к нулю при, а среднеквадратичное значение оценки при любом N сравнимо со средним значением. Следовательно, выборочный спектр не является состоятельной (статистически устойчивой) оценкой СПМ. На рис.19 – 22 показаны периодограммы, вычисленные по (5.42), для суммы комплексной синусоиды с амплитудой, равной 1,0 и частотой 0,2 и белого нормального (гауссова) шума с дисперсией, равной 1,0 при N=50, 100, 200, 500 соответственно. Из сравнения рис. 19 – 22 видно что с ростом N дисперсия спектральной оценки не уменьшается, а значения периодограммы начинают быстрее флюктуировать.

57 DSP Цифровой спектральный анализ Рисунок 19. Периодограмма суммы комплексной синусоиды с амплитудой, равной 1,0 и частотой 0,2 и белого нормального (гауссова) шума с дисперсией, равной 1,0 при N=50.

58 DSP Цифровой спектральный анализ Рисунок 20. Периодограмма суммы комплексной синусоиды с амплитудой, равной 1,0 и частотой 0,2 и белого нормального (гауссова) шума с дисперсией, равной 1,0 при N=100.

59 DSP Цифровой спектральный анализ Рисунок 21. Периодограмма суммы комплексной синусоиды с амплитудой, равной 1,0 и частотой 0,2 и белого нормального (гауссова) шума с дисперсией, равной 1,0 при N=200.

60 DSP Цифровой спектральный анализ Рисунок 22. Периодограмма суммы комплексной синусоиды с амплитудой, равной 1,0 и частотой 0,2 и белого нормального (гауссова) шума с дисперсией, равной 1,0 при N=500.

61 DSP Предложен ряд модификаций соотношения (5.42), суть которых сводится к сглаживанию оценки путем введения псевдоусреднения по некоторому ансамблю. Наибольшую популярность в настоящее время получил метод Уэлча, идея которого состоит в образовании псевдоансамбля периодограмм за счет разбиения наблюдаемых данных на перекрывающиеся сегменты и обработки каждого сегмента с помощью окна данных. Если запись комплексных данных, разбита на P сегментов по D отсчетов в каждом со сдвигом S отсчетов между соседними сегментами, то максимальное число сегментов P есть целая часть величины (N-D)/S+1. После взвешивания окном сегмент с номером будет содержать отсчеты Цифровой спектральный анализ Выборочный спектр P - го сегмента определяется выражением (5.43) где -ДВПФ сегмента, - энергия окна. Среднее значение сегментных периодограмм определяет оценку СПМ: (5.44)

62 DSP Уэлч, в частности, предложил использовать окно Ханна и 50%-ное перекрытие сегментов, которое обеспечивало очень эффективные реализации его метода на основе алгоритма БПФ. Кроме того, при 50%- ном перекрытии сегментов все данные используются дважды, за исключением D/2 отсчетов на каждом конце исходной N -точечной последовательности данных, а это выравнивает обработку большинства отсчетов данных, поскольку те отсчеты, которые имели малые веса на одном сегменте, получают большие веса на следующем сегменте. Анализ поведения периодограммы Уэлча для гауссовских процессов показал, что минимальная дисперсия для окна Ханна достигается при 65%-ном перекрытии, при этом величина дисперсии увеличивается приблизительно на 8% при использовании 50%-ного перекрытия сегментов. Дисперсия периодограммы Уэлча примерно обратно пропорциональна числу сегментов, т.е. Цифровой спектральный анализ (5.45) в предположении независимости сегментов (хотя перекрытие сегментов приводит, конечно, к некоторой их взаимозависимости). Благодаря перекрытию по заданной записи данных можно сформировать большее число сегментов и уменьшить величину дисперсии. На рис.23 представлена периодограмма Уэлча суммы комплексной синусоиды с амплитудой, равной 1,0 и частотой 0,2 и белого нормального (гауссова) шума с дисперсией, равной 1,0 при N=500, D=50, S=25 ( 19 сегментов).

63 DSP Цифровой спектральный анализ Рисунок 23. Периодограмма Уэлча суммы комплексной синусоиды с амплитудой, равной 1,0 и частотой 0,2 и белого нормального (гауссова) шума с дисперсией, равной 1,0 при N=500, D=50, S=25 (19 сегментов).

64 DSP Цифровой спектральный анализ Использование окон при спектральных измерениях. Один из важных вопросов, который является общим для всех классических методов спектрального оценивания, связан с применением функций окна. Заметим, что конечную запись данных из N отсчетов можно представить как некоторую часть исходной бесконечной последовательности, видимую через прямоугольное окно, в форме произведения: где При этом мы полагаем, что все ненаблюдаемые отсчеты равны нулю независимо от того, так ли это на самом деле или нет. ДВПФ взвешенной окном последовательности, выраженное через преобразования последовательности и прямоугольного окна, равно свертке этих преобразований: Где- ядро Дирихле, представляющее ДВПФ прямоугольной функции. ДВПФ наблюдаемой конечной последовательности является искаженной версией ДВПФ бесконечной последовательности.

65 DSP Влияние прямоугольного окна на дискретно-временную синусоиду с частотой f 0 иллюстрирует рис. 24, из которого видно, что острые спектральные пики ДВПФ исходной синусоидальной последовательности расширились из-за воздействия ДВПФ окна. При этом минимальная ширина спектральных пиков конечной последовательности ограничена шириной, определяемой главным лепестком ДВПФ окна, и не зависит от исходных данных. Боковые лепестки ДВПФ окна, называемые растеканием [просачиванием (spectral leakage)] будут изменять амплитуды соседних спектральных пиков, приводя к смещению спектральных оценок. Аналогичные искажения будут наблюдаться и в случае несинусоидальных сигналов. Просачивание приводит не только к появлению амплитудных ошибок, но может также маскировать присутствие слабых сигналов и, следовательно, затруднять их обнаружение. Цифровой спектральный анализ

66 DSP Цифровой спектральный анализ Рисунок 24. Иллюстрация эффекта растекания спектра из-за взвешивания данных: а – исходная синусоидальная последовательность; б – модуль периодического ДВПФ синусоидальной последовательности; в – взвешенная синусоидальная последовательность; г - модуль ДВПФ взвешенной синусоидальной последовательности.

67 DSP Предложен ряд других функций окна, имеющих меньший уровень боковых лепестков их спектра, чем в случае прямоугольного окна. Однако, как правило, это достигается ценой расширения главного лепестка спектра окна. Следовательно, окно должно выбираться с учетом компромисса между шириной главного лепестка и уровнем подавления боковых лепестков. В таблице 5.2 представлены сравнительные характеристики некоторых оконных функций. Из всех приведенных окон самый узкий главный лепесток имеет частотная характеристика прямоугольного окна, однако у него самый высокий уровень боковых лепестков. На рис.25 представлены типичные дискретно-временные функции окна (слева) и логарифм модуля их ДВПФ (справа) для прямоугольного окна, треугольного окна, окна Ханна и окна Хэмминга. В таблице 5.2 приведены определения некоторых используемых N - точечных функций окна симметричных относительно n=0. Именно такие окна используются в качестве корреляционных окон в коррелограммном методе оценки СПМ. При этом используются только окна нечетной длины с w(0)=1. Для окон данных (в периодограммном методе оценки СПМ) применяются окна четной и нечетной длины, начинающиеся с n=0 и имеющие симметрию относительно точки (N-1)/2. При этом у окон четной длины точка симметрии находится посередине между отсчетами. Цифровой спектральный анализ

68 DSP Цифровой спектральный анализ ОкноДискретно- временная функция w(n), -N/2 n N/2 Максимальный уровень боковых лепестков (дб) Ширина полосы по уровню 0,5 мощности относительно 1/NT Прямоугольное1-13,30,89 Треугольное (Бартлетта)) -26,51,28 Ханна-31,51,44 Хэмминга-431,3 Таблица 5.2. Характеристики оконных функций

69 DSP Цифровой спектральный анализ Рисунок 25. Типичные дискретно-временные функции окна (слева) и логарифм модуля их ДВПФ (справа): а-б -прямоугольное окно; в-г -треугольное окно; д-е -окно Ханна; ж-з- окно Хэмминга.

70 DSP Цифровой спектральный анализ Разрешение и произведение «устойчивость длительность ширина полосы». Для детерминированного сигнала с длительностью секунд произведение длительности на ширину полосы удовлетворяет условию. Поэтому разрешение спектральной плотности энергии, полученной по конечной последовательности отсчетов детерминированного сигнала, равно B e герц, т.е. приблизительно равно величине, обратной интервалу записи данных. Для того чтобы определить разрешение спектральной плотности мощности, получаемой по конечной последовательности отсчетов случайного сигнала, произведение длительности на ширину полосы необходимо несколько видоизменить, с тем чтобы учесть случайный характер сигнала, влияющий на статическое качество спектральной оценки. Одной из подходящих модификаций этого произведения является произведение устойчивость длительность ширина полосы, где - статистический показатель качества, определяемый как отношение дисперсии оценки СПМ к квадрату математического ожидания этой оценки (5.46)

71 DSP Цифровой спектральный анализ Показатель качества -это инвертированное отношение сигнал/шум (SNR), которое непосредственно связано со статистической устойчивостью спектральной оценки. Значения, много меньшие единицы, соответствуют гладким спектральным оценкам с малой дисперсией. Значения, много большие единицы, соответствуют весьма зашумленным спектральным оценкам с большой дисперсией. Эквивалентная ширина полосы B e замеряется теперь величиной B s - эффективной статистической шириной полосы, которая в случае спектрального окна определяется выражением Величина характеризует ширину полосы эквивалентного прямоугольного окна с тем же значением отношения дисперсии к квадрату среднего значения на его выходе, как и у окна в том случае когда на его вход подан белый шумовой процесс. Однако обе эти величины можно связать соотношением, где выбирается в зависимости от окна, так для прямоугольного окна =0.8, для окна Ханна-1.3 и для окна Хемминга Ширина полосы позволяет грубо судить о разрешении спектральной оценки

72 DSP Показано, что для гауссовского случайного процесса условие Цифровой спектральный анализ (5.47) применимо ко всем классическим процедурам спектрального оценивания. Если, то, как правило. С другой стороны, если значение произведения выбрано значительно большим единицы, то для всех классических спектральных оценок будет выполнятся приближенное равенство (5.48) Произведение «устойчивость длительность ширина полосы» устанавливает взаимосвязь между тремя фундаментальными параметрами, от которых зависят характеристики спектральных оценок. Для заданной записи данных длительностью T e невозможно получить оценки, которые одновременно обладают высоким разрешением (малыми значениями B S ) и высокой устойчивостью (малыми значениями ). Например, если необходимо, чтобы дисперсия спектральной оценки составляла десятую часть от ее среднего значения, то Q=0,1 и достижимое разрешение будет равно B e 10/T e герц. Это значение в 10 раз больше эмпирического значения 1/T e для разрешения детерминированного сигнала..

73 DSP Этот пример подчеркивает важность принятия такого компромиссного решения, которое не нарушало бы баланс между статистической устойчивостью и разрешением, которые требуются от всех методов спектрального оценивания. Заметим, что увеличение числа отсчетов за счет увеличения частоты дискретизации при неизменном интервале T e не будет влиять на максимально достижимое разрешение, так как на него влияет только длина записи данных, но не число отсчетов. На рис.26 и 27 показаны периодограммы Уэлча суммы двух комплексных синусоид с амплитудами, равными 1,0 и частотами 0,2 и 0,24 и белого нормального (гауссова) шума с дисперсией, равной 1,0 при N=500, D=50, S=25 ( 19 сегментов) и при N=500, D=250, S=125 ( 3 сегмента) соответственно, которые иллюстрируют соотношение между устойчивостью и разрешением. Рис. 26 соответствует спектральной оценке с большей устойчивостью, но меньшим разрешением, а рис. 27 соответствует спектральной оценке с меньшей устойчивостью, но большим разрешением. Цифровой спектральный анализ

74 DSP Цифровой спектральный анализ Рисунок 26. Периодограмма Уэлча суммы двух комплексных синусоид с амплитудами, равными 1,0 и частотами 0,2 и 0,24 и белого нормального (гауссова) шума с дисперсией, равной 1,0 при N=500, D=50, S=25 (19 сегментов).

75 DSP Цифровой спектральный анализ Рисунок 27. Периодограмма Уэлча суммы двух комплексных синусоид с амплитудами, равными 1,0 и частотами 0,2 и 0,24 и белого нормального (гауссова) шума с дисперсией, равной 1,0 при N=500, D=250, S=125 (3 сегмента).

76 DSP Цифровой спектральный анализ Классические спектральные оценки тест-последовательности. Для того, чтобы охарактеризовать поведение каждой спектральной оценки при анализе короткой последовательности данных Марпл предложил использовать искусственно синтезированную последовательность данных с заранее известными свойствами. Эта тест-последовательность данных содержит 64 отсчета некоторого комплекснозначного процесса, состоящего из 4 -х комплексных синусоид и комплексного аддитивного окрашенного шума. Истинный спектр тест-последовательности показан на рис.28. Рисунок 28. Истинный спектр процесса, соответствующего тест-последовательности данных.

77 DSP Частоты двух синусоид выбраны очень близкими друг к другу( 0,2 и 0,21 ) с тем, чтобы проверить разрешающую способность той или иной спектральной оценки. Мощность двух более слабых синусоид с частотами 0,1 и –0,15 на 20 дб меньше мощности двух более сильных. Эти слабые синусоиды используются для испытания способности спектральных оценок обеспечить обнаружение слабых компонент сигналов на фоне сильных сигналов. Окрашенный шумовой процесс был сформирован пропусканием двух независимых белых шумов через идентичные фильтры с частотными характеристиками типа приподнятой косинусоиды для раздельного получения действительной и мнимой составляющих шумового процесса в тест-последовательности. Высота линий, характеризующих синусоидальные составляющие выбрана так, чтобы отображать мощность каждой синусоиды относительно мощности окрашенного шума. Следовательно, мощность каждой из сильных синусоид превышает мощность шума, а мощность каждой из 2-х слабых синусоид меньше мощности шума. Приведем примеры классических периодограммного и коррелограммного методов спектрального оценивания СПМ тест-последовательности. Цифровой спектральный анализ

78 DSP На рис. 29 приведены две периодограммы тест-последовательности Марпла. Периодограмма этой последовательности, соответствующая максимальному разрешению и минимальной устойчивости, показана на рис.29а. Она получена для случая только одного сегмента (полной записи данных) без применения окна. Заметим, что близкие по частоте синусоиды с разносом 0,01 не разрешены, поскольку предельное разрешение в данном случае составляет величину 1/64=0,0156. Слабая синусоида с частотой 0,1 маскируется боковыми лепестками сильных синусоид. На рис. 29 б показана гладкая усредненная периодограмма с низким разрешением, которая получена разбиением записи данных на три сегмента по 32 отсчета в каждом, взвешенным окном Хэмминга, с перекрытием сегментов, равным 16 отсчетов. На этой периодограмме отчетливо видны слабые синусоиды. Коррелограммные оценки СПМ тест-последовательности показаны на рис. 30 а, б для 16 и 32 значений автокорреляции соответственно и взвешивания окном Хэмминга. Графики на рис. 30 дают полное представление о взаимосвязи между разрешением и гладкостью оценки. Цифровой спектральный анализ

79 DSP Цифровой спектральный анализ Рисунок 29. Периодограммные оценки СПМ для 64- точечной тест- последовательности: а- один сегмент, окно не применялось; б – три сегмента, окно Хэмминга. а) б)

80 DSP Цифровой спектральный анализ Рисунок 30. Коррелограммные оценки СПМ для 64-точечной тест- последовательности: а – 16 значений автокорреляции, окно Хэмминга; б – 32 значения автокорреляции, окно Хэмминга. а) б)

81 DSP Заметим, что если не известно как выбирать компромиссное решение относительно устойчивости и разрешения оценок, то рекомендуется начать с выбора максимальной устойчивости (наиболее гладкой оценки), а это значит – с минимального разрешения, а затем постепенно увеличивать длину сегментов (в случае периодограммы) или число корреляционных членов (в случае коррелограммы). Эта процедура будет увеличивать дисперсию оценки и уменьшать ширину спектральных пиков (повышать разрешение). Следует продолжать указанную процедуру до тех пор, пока достигнутое разрешение не позволит передать наиболее тонкие и важные детали спектра и дальнейшее увеличение разрешения к существенным изменениям в спектральных оценках уже не приводит. Цифровой спектральный анализ