Сергиенко А.Б. Цифровая обработка сигналов (2002) (1095939), страница 52
Текст из файла (страница 52)
(5.32) »=и+~ )'» «-О Ыа В этих формулах Х» — собственное число, соответствующее я-му собственному вектору ч». В 1101 приводятся сведения о том, что при заданном значении порядка М метод ЕЧ порождает меньше ложных спектральных пиков, чем МП51С, и, как правило, лучше передает форму спектра шума. Еще раз подчеркнем, что псевдоспектры (5,29) и (5.31) не являются оценками истинного спектра плотности мощности, а представляют собой лишь спектральные псевдооценки, позволяющие оценивать частоты синусоидальных или узкополосных составляющих сигнала с разрешением, несколько превосходящим разрешение авторегрессионных методов. Функции спектрального анализа в МАТЫВ Для выполнения спектрального анализа в МАТ1.АВ имеется большое число функ- ций, реализующих разнообразные методы — от обычного ДПФ до М()51С.
288 Функции спектрального анализа в МАТ(АВ Прямое и обратное ДПФ Для выполнения прямого и обратного ДПФ в МАТ1АВ служат функции ГТС и )ттс; 0 у - ттс(х) — вычисляет прямое ДПФ для вектора х; если х — матрица, преобразование производится для каждого ее столбца по отдельности; 0 у = ттс(х, й) — предварительно приводит исхолные данные к размеру М, урезая их или дополняя нулями; 0 х = зттс(у) и х = 1тгс(у, й) — аналогичные варианты вызова лля функции обратного ДПФ. Функции ГТ1 и зтт1 входят в базовую библиотеку МАТБАВ. Вычисления организованы так, что реализуется максимально возможное для каждой длины исходного вектора ускорение вычислений: длина вектора (число строк в матрице) х раскладывается на простые множители, число этих множителей соответствует количеству ступеней БПФ, а сами множители определяют коэффициенты прореживания на разных ступенях БПФ.
Продемонстрируем зависимость числа вычислительных операций от размерности БПФ. Создадим вектор из 128 случайных чисел (на всякий случай — комплексных, чтобы лишить МАТ|АВ возможности что-либо сэкономить на отсутствии мнимой части у исходных данных; впрочем, при использовании алгоритмов БПФ такая экономия очень мала). ВНИМАНИЕ Приводимый ниже пример работает только в младших версиях МАТ|А — до,б.х включительно. Дело в том, что в версии 6 было обновлено ядро МАТ|АВ и функция Норв, подсчитывавшая число выполненных арифметических операций, более не поддерживается. й - 128: х - гапбп(й.
1) + з * гапон(й. 1); Гогй1й т)орз(0): Х В МАТ(АВ б зто не работает! у Т(1(х, й); т)орй) - т1орз: Ж В МАТ(АВ б зто не работает> епб р)ог(Г)ор) По16 оп р1ог((1:й)."2*8. ': ') р1о1((1;й).*1од2(1:й)*б, ': ') п016 Отт Полученный график (рис. 5.14) настолько показателен, что не привести этот пример было бы совершенно непростительно, несмотря на то, что он не работает в последних версиях МАТ|АВ.
Кроме собственно зависимости числа операций от размерности БПФ на графике показаны также верхняя и нижняя границы. Верхняя линия — это число операций, соответствующее прямой формуле ДПФ и равное 8Хз (множитель 8 возникает из-за того, что функция ~1орз подсчитывает вещеплвенные операции). Нижняя а8Е Глава В.
Спектральный анализ линия — число операций, соответствующее максимальному ускорению вычислений и равное 62тГ 1одз(У) (множитель 6 также возникает из-за подсчета функцией Г)срз вещественных операций). Хорошо видно, что график числа операций при БПФ касается верхней границы, если )т'- простое число, и касается нижней границы, если )т' — степень двойки. х 10ь 14 10 О 0 20 40 60 60 100 120 140 Рис. В.! 4.
Зависимость числе арифметических операций от размврнооти БПФ Функция ЛФВЫЙ Элементы вектора, возвращаемого функцией ГГт, соответствуют частотам, равномерно распределенным в диапазоне от нуля и почти до частоты дискретизации. Первый элемент, таким образом, соответствует нулевой частоте, а последний— частоте, меньшей частоты дискретизации на/х/У, где )тт- размер входного и выходного векторов. При выводе спектральных графиков иногда желательно, чтобы нулевая частота находилась в центре, а диапазон отображаемых частот простирался от -~,/2 до/„/2.
Сделать зто позволяет функция тттзй1гт, которая меняет местами половины переданного ей вектора: у - ГгтзП1ГГ(х) Продемонстрируем действие функции тгтзЫ гс на примере двух коротких векторов четной и нечетной длины: » Гтгзп)тгП) 2 3 41) апз = 3 4 1 2 » 1ГГз)т11ГП1 2 3 4 5]) дп5 4 5 ) 2 3 Функции спектрального анализа а ГиАТЫВ Как видите, в случае четной длины действительно происходит перестановка половин входного вектора (при этом первый элемент результирующего вектора соответствует частоте Найквиста). В случае нечетной длины перестановка выполняется так, чтобы первый элемент, соответствующий нулевой частоте, стал средним элементоаг результирующего вектора. Разумеется, функцию Тй~зл1ТС можно использовать не только для вывода спектральных графиков, но и в других случаях, когда требуется поменять местами половины вектора.
Матрица ДПФ В МАТЮКАВ расчет матрицы прямого ДПФ реализуется с помощью функции бугаях. Синтаксис ее вызова следующий: А - от1шГх(л) Здесь и — размерность ДПФ. Матрица обратного ДПФ отличается от матрицы прямого ДПФ комплексным сопряжением и делением на Аг, поэтому она может быть получена следующим образом: А 1лу сол3(бугпгсх(п))!л Блочная фильтрация а частотной области Для реализации блочной фильтрации с помощью БПФ в МАТБАВ используется функция Л~ст111, имеющая следующий синтаксис: у - ТТ1т1'11(л, х. л) Здесь и — импульсная характеристика фильтра, х — входной сигнал, и — размерность БПФ (точнее, размерность БПФ при задании этого параметра определяется так: йттс - 2"пехтрои2(п), то есть л округляется вверх до степени двойки). Выходной параметр у — результат фильтрации.
Третий входной параметр в при вызове можно опускать, тогда размерность БПФ будет выбираться автоматически, исходя из максимальной эффективности вычислений (минимального числа операций). Окна МАТЮКАВ содержит (в пакете 51япа1 Ргосезз1пя) целый ряд стандартных весовых функций. Они возвращают векторы отсчетов, которые могут использоваться в качестве одного из параметров разнообразных функций непараметрического спектрального анализа. Все рассматриваемые ниже функции принимают в качестве параметра требуемую длину вектора (в), которая должна быть целым положительным числом, и возвращают вектор-столбец и. При л=1 все функции возвращают значение 1. Амплитудный спектр весовой функции соответствует частотной характеристике нулевого канала ДПФ при использовании данной весовой функции.
При рас- 288 Глава В. Спектральный анализ смотрении конкретных функций графики их амплитудных спектров будут строиться в логарифмическом масштабе для и = 16. Чтобы обеспечить на нулевой частоте значение спектральной функции, равное единице (О дБ), перед вычислением спектра весовые функции нормируются — делятся на сумму своих отсчетов. Графики спектров будут строиться функцией ггецг (см. раздел «Расчет частотной ь характеристики» главы 4). Поскольку фазовый спектр для всех весовых функций линейно зависит от частоты, его графики не представляют интереса и потому не приводятся. Для повышения наглядности частотная ось градуируется в номерах каналов ДПФ, для этого при вызове функции ггейг указана частота дискретизации, равная длине окна.
Прямоугольное окно Функция Ьохсаг, реализующая «прямоугольное окно», введена в МАТ1.АВ лишь для полноты набора весовых функций, поскольку она соответствует отсутствию взвешивания: и - 'оохсзг(п) Возвращаемый вектор заполнен единицами: и = опез(П,1). Приводить график данной весовой функции не имеет смысла, а ее амплитудный спектр (смещенный по горизонтали) для и - 6 был показан ранее на рис.
5.6. Треугольное окно Функция тг(ап9 реализует треугольное окно: и = тг1дп9(п) Отсчеты треугольного окна рассчитываются по следующей формуле: 1<() ~ —, и+1 2 для нечетного и, и+1 — <и <и 2 2() и+1 2(и — (г+ 1) гл(к) = и+1 21 — 1 1ьк <— 2 для четного и. 2(и — Й+ 1) и и ' 2 ге((г) = При нечетном и треугольное окно является симметричным, его крайние значе- ния (при (з - 1 и (г - и') равны 2/(и + 1), а в середине окна (при (з - (и + 1)/2) достигается единичное значение. При четном п треугольное окно является несимметричным, его можно предста- вить как отсчеты симметричного треугольного импульса, который начинается при (з = 0,5 и заканчивается при 1 = и + 1.
Вершина импульса расположена прн )г = и/2» 0,75, а его амплитуда равна 1 + 1/(2и). При )1 = и/2 + 1 отсчет окна име- ет единичное значение. На рис. 5.15 приведены графики треугольного окна и его амплитудного спектра при и - 16: 289 Функции спектрального анализа в МАТ(А5 » и = Сг)ап9(16): » и - и)зов(и): » р10((и) » Йооге » (ь, 71 Тгеоа(и. 1, [). 16); » р10С(7, 20*10910(аоз(Ь))) » у11а)(Г-60 01) > 9г1О оп 0,12 0,1 о.ов 0.0б о.ог о 2 4 б в то тз 14 тб -10 -20 -30 -ао -50 -70 0 1 2 3 4 б б 7 В Рис. 5.15. Треугольное окно (сверку) и его амплитудный спектр (снизу) Уровень первого бокового лепестка составляет -26,5 дБ.