Лабораторная работа №61

Санкт-Петербургский Государственный

Электротехнический Университет

Кафедра МОЭВМ

Лабораторная работа № 6

«ИССЛЕДОВАНИЕ ПРОЦЕДУР

СПЕКТРАЛЬНОГО АНАЛИЗА СИГНАЛОВ»

Студент: гр. 3341 Митягин С.А.

Преподаватель: Жукова Н.А.

Санкт-Петербург

2007

Цель работы

Целью работы является изучение методов использования алгоритма быстрого преобразования Фурье (БПФ) для решения задач оценки спектральных характеристик различных сигналов с помощью ЭВМ, а также исследование свойств дискретного преобразования Фурье, непосредственно влияющих на качество получаемых спектральных оценок.

2. Постановка задачи

2.1. Содержательная постановка

2.1.1. Спектральное представление детерминированных и случайных процессов обладает рядом полезных свойств, обусловивших его широкое распространение в практических задачах исследования сигналов и динамических систем. В число таких свойств входят:

- высокая информативность (хорошо известно, что временное и частотное представление некоторого процесса связаны взаимнооднозначным соотношением);

- наглядность и возможность физической интерпретации;

- имеющая место в ряде практических случаев независимость соседних частотных отсчетов;

- возможность описания основных характеристик случайных процессов их энергетическим спектром;

- непосредственная связь корелляционных свойств исследуемого процесса с его спектральным представлением;

- наличие эффективных методов цифрового спектрального анализа. Однако поскольку при проведении спектрального анализа на ЭВМ

используется дискретное преобразование Фурье (ДПФ), свойства которого отличаются от свойств непрерывного преобразования Фурье и коэффициентов Фурье, используемых в теоретическом анализе случайных и детерминированных процессов, методика использования ДПФ обладает рядом характерных особенностей.

Типовая последовательность обработки исходного процесса при цифровом спектральном анализе представляется следующим образом. Исследуемая реализация процесса предварительно подвергается дискретизации или квантованию. В результате квантования реализация представляется последовательностью дискретных отсчетов, отстоящих друг от друга на заданный временной интервал. Квантованная реализация разбивается на участки, имеющие одинаковое количество отсчетов. Количество участков и их длина определяется требуемым частотным разрешением и требуемой точностью оценки спектра. С целью уменьшения влияния конечной длины участков на качество получаемых спектральных оценок каждый участок умножается на специальное взвешивающее временное окно. Для каждого взвешенного участка с помощью алгоритма БПФ вычисляется оценка мгновенного энергетического спектра. Полученная последовательность оценок мгновенного энергетического спектра усредняется по частоте и/или амплитуде, давая требуемую оценку энергетического спектра исследуемого процесса.

2.2. Математическая постановка задачи обработки

2.2.1. Функционирование цифровой системы спектрального анализа (ЦССА) можно формально описать следующим образом:

1. На вход ЦССА поступает последовательность отсчетов Y_k, k=0,1,2,... временной функции y(t), взятых в моменты времени 0,T,2T,3T,... . Величина T носит название периода квантования, а сам процесс формирования последовательности по по временной непрерывной функции называется квантованием.

2. Последовательность отсчетов Y_k разбивается на отрезки, состоящие из N отсчетов, где N обычно представляется в виде 2^M (2,4,8,16,32,64,...).

3. Каждый отрезок подвергается дискретному преобразованию Фурье (ДПФ). В результате получается последовательность X_j, j=0,1,..., N-1, состоящая из N комплексных спектральных отсчетов. Для дальнейшей обработки, как правило, используется действительная последовательность, состоящая из N/2 отсчетов, каждый из которых является квадратом модуля соответствующего комплексного коэффициента X_j:

Z_j = (ReX_j)² + (ImXj)².

4. Для улучшения точности полученных спектральных оценок периодограмма подвергается операциям частотного и временного сглаживания, в результате чего получают оценку энергетического спектра процесса Y(t) в виде последовательности действительных коэффициентов.

Рассмотренному явлению наложения соответствует явление переноса частот сущность которого в случае гармонического сигнала может быть сформулирована следующим образом:

- при квантовании гармонического сигнала с частотой f, лежащей в пределах k f_кв(2k+1) f_кв/2, где k - некоторое целое число, спектр квантования сигнала будет в точности соответствовать спектру гармонического сигнала с частотой f-k f_кв;

- при квантовании гармонического сигнала с частотой f_l, лежащего в пределах (2k+1) f_кв/2(k+1) f_кв, где k - некоторое целое число, спектр квантованного сигнала будет в точности соответствовать спектру гармонического сигнала с частотой (k+1) f_кв-f_l.

Кроме явления квантования, влияние которого на свойства ДПФ мы только что рассмотрели, ряд параметров ЦССА существенно зависят от длины анализируемого участка временной последовательности X_k. Характер этой зависимости определяется конечностью длины участка временной реализации. Действительно, рассмотрим случай, когда исходной временной последовательностью является часть функции cos(t), взятая для t=iT, i=0,N-1. полученную последовательность дискретных отсчетов можно рассматривать как результат квантования непрерывной функции, являющейся произведением функции cost на "прямоугольное временное окно" вида:

Oi = (4)

В рассматриваемом случае можно показать, что для косинусоиды с частотой _s квадрат модуля комплексного отсчета ДПФ, соответствующего частоте wi=2p if_кв/N описываются выражением:

P(_s-_i) = |sin[(_s-_i) N/f_кв]/(_s-_i) N/f_кв| (5)

Отметим, что функция (5) имеет центральный "лепесток", соответствующий глобальному максимуму функции, и боковые "лепестки", соответствующие локальным максимумам. С учетом выражения (5) можно дать физическую интерпретацию ДПФ: квадрат модуля комплексного отсчета ДПФ, соответствующего частоте _i, можно интерпретировать как мощность сигнала на выходе полосового фильтра, передаточная характеристика которого описывается выражением (5). Таким образом, первые N/2 полосовых фильтров, передаточная характеристика которых описывается выражением (5), а центральные частоты равны f_i=f_кв/N. Таким образом, функция P(_s,_i) описывает избирательные свойства ДПФ, которые обычно оценивают эффективной полосой пропускания - B_l и степенью подавления боковых лепестков - k_зат.

Величина B_l оценивается на уровне 0.5 от максимального значения коэффициента передачи мощности P_i(). Для ДПФ:

B_l ~ f_кв/N (6)

где f_кв - частота квантования, а N - количество отсчетов временной реализации.

Наличие боковых лепестков может существенно снизить точность спектрального оценивания за счет "проникновения" сигналов из полос частот, соответствующих соседним спектральным отсчетам. Поэтому для увеличения k_зат принимаются специальные меры, приводящие к уменьшению величины боковых лепестков. Наиболее часто для этого используется взвешивание обрабатываемого отрезка реализации функцией временного окна в соответствии с выражением:

X_взвi = Y_i k_окнi

На практике используются самые различные виды функций k_окнi. Общим для всех этих функций является соотношение k_окн = 0, _i<0 и _i.N. Одними из наиболее распространенных функций временного окна являются:

I. Прямоугольное:

Этому окну соответствует P(), определяемое соотношением (5).

2. Треугольное:

3. Хеннинга:

4. Хемминга:

5. 4-х членное окно Блэкмана-Хэрриса:

= 0,35875-0,48829cos(2(i-1)/N))+0,14128cos(4(i-1)/N)) - 0,01168cos(6(i-1)/N)).

На точность оценок спектральных характеристик существенно большее влияние, чем уровень боковых лепестков, оказывает наличие случайных флуктуаций, величину которых можно оценить коэффициентом разброса:

K_разб = _x /m_x,

где m_x и _x - математическое ожидание и среднеквадратическое отклонение оценки величины спектрального отсчета. Для случая, когда оцениваемый спектр равномерен в полосе B_l, коэффициент разброса легко выражается через полосу разрешения B_l и длины реализации или время наблюдения T_Н .

K_разб = 1/().

Из соотношения видно, что понижения K_разб можно добиться двумя способами:

- снижением разрешающей способности путем увеличения величины B_l;

- увеличением времени наблюдения T_Н.

Перечисленные два способа уменьшения коэффицинта разброса соответствуют двум основным практическим методам:

- методу осреднения оценки спектра по частоте и

- методу осреднения по времени.

Изучение частотной характеристики БПФ

Задать параметры:

n=1024 - длина реализации БПФ (N=10, f_дискретизации будет равно 1024 Гц).

f_sin=100 Гц f_sin=100,5 Гц

f_sin=700 Гц f_sin=1100 Гц

f_sin=1700 Гц f_sin=2100 Гц

Выводы: если частота сигнала больше величины 2^(N-1), наблюдаем эффект переноса частот (частоты f_sin=700 Гц, f_sin=1100 Гц). Если частота сигнала не попадает в отсчеты ДПФ, наблюдаем эффект просачивания сигнала через боковые лепестки и уменьшения максимального значения спектра сигнала (f_sin=100,5 Гц). Оно обусловлено некратностью частоты синуса шагу частоты дискретизации.

При этом на спектре вместо единственного истинного значения амплитуды виден набор значений с меньшей амплитудой, максимальное из которых будет зависеть от соотношения частоты сигнала и шага частоты дискретизации.

Изучение разрешающей способности БПФ по частоте.

Задать параметры:

N=1024 - длина реализации БПФ (n=10) ,

f1=300, f2=310 проверить разрешение по частоте.

Подобрать длину реализации сигнала ,при которой сохраняется разрешение. Объяснить результаты. Надо учесть ,что при моделировании для заданного N длина реализации сигнала уменьшается путем замены части отсчетов нулями. Реальное разрешение по частоте равно при этом 1/T, где T- длина сигнала.

f1=200, f2=210

f1=200, f2=205

При уменьшении длины реализации сигнала разрешение по частоте уменьшается.

Так как по Теореме Котельникова сигнал можно восстановить с большой точностью, если максимальная частота, на которой спектр непрерывной функции отличен от нуля, не превышает величину f_кв/2.

Чем ближе по частоте расположены сигналы, тем больше должна быть длина их реализации, чтобы была возможность их разделить.

Применение окон для спектрального анализа

Задать параметры:

n=1024 - длина реализации БПФ (N=10) ,

M=20 - ширина окна.

Необходимо учесть , что при таком выборе М и N 1 бин будет соответство-

вовать 1024/20 =50 точек. ( здесь использована интерполяция для расширения

графического изображения окна).

Нажать кнопку "Без синуса"

Выбрать тип окна, нажать кнопку "Выполнить БПФ"

Включить " шкала логарифмическая"

Включить "пропускать через окно"

(При включенной кнопке "пропускать через окно" программа работает с

заданным окном , при выключении кнопки - можно рассчитать спектр

синусоидального сигнала без окна.)

Просмотреть частотные характеристики различных окон. Определить

ширину главного лепестка и скорость спада боковых лепестков (для первых

пяти лепестков ) для всех окон.

Провести сравнение спектра синусоидального сигнала при использовании

БПФ без окна (при этом выключается "пропускать через окно") и с

использованием всех окон ( установить "с синусом", " частота синуса",

n=10 , ширина окна M=100)

Прямоугольное окно:

Сигнал без окна:

Прямоугольное окно:

Без синуса:

С синусом:

Треугольное окно:

Без синуса:

С синусом:

Окно Хеннинга:

Без синуса:

С синусом:

Окно Хемминга:

Без синуса:

С синусом:

Окно Блекмена-Херриса:

Без синуса:

С синусом:

тип окна	Мах боковой лепесток, дб	Скорость спада боковых лепест дб/ окт	Полоса осн. лепестка на уровне 3дб, бин
Прямоугольное	-13	-6	0,9
Треугольное	-27	-12	1,3
Хеннинга	-32	-12	1,4
Хемминга	-67	-20	1,8
Блэкмана-Хэрриса (4 члена)	-92	-6	1,9

Как правило, чем сильнее применяемое окно подавляет боковые лепестки, тем больше у него величина полосы пропускания.

Выводы: наличие боковых лепестков может существенно снизить точность спектрального оценивания за счет «проникновения» сигналов из полос частот, соответствующих соседним спектральным отчетам. Наиболее эффективно боковые лепестки подавляются при использовании окна Блэкмана-Хэрриса. При этом платой за большую скорость убывания боковых лепестков является расширение центрального лепестка спектра окна, поэтому на практике необходимо выбирать компромисс между скоростью спада боковых лепестков и шириной центрального лепестка.

Исследование влияние накопления на характеристики спектрального анализа.

Количество точек реализации n= 1024 (частота дискретизации 1024 Гц).

Мощность шума здесь-С.К.О (sigma)

Амплитуда полезного сигнала (SIN)-A

Отношение сигнал/помеха q=A*A/(2*sigma*sigma) вычисляется

в программе. Это отношение сигнал/помеха во входном сигнале. После

спектрального анализа отношение в каждом фильтре БПФ становится q1=q*N.

Если сделать M накоплений по времени, то отношение q2=q1*(M)^0.5. Cигнал на

спектре визуально можно обнаружить при отношении сигнал/помеха не менее

8-10.

Sigma=10

Теоретический расчет:

Применение усреднения по времени уменьшает коэффициент разброса в раз.

Тогда чтобы отношение сигнал/помеха был равен 10:

Подтверждение:

Sigma=25

Теоретический расчет:

Sigma=50

Теоретический расчет: