Дискретизация. Частота Найквиста и теорема Котельникова

Исходный физический сигнал x (t) является непрерывной функцией времени и называется аналоговым. Числа x_k = x (t_k) называются отсчетами сигнала. Как правило, отсчеты берутся через равные промежутки времени Δ t – период или шаг дискретизации. Величина f_d, обратная периоду дискретизации, называется частотой дискретизации:

Если Δ t выражен в секундах, то частота дискретизации выражается в герцах (Гц).

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

     Если в сигнале x (t) присутствует периодическая составляющая с периодом T <2Δ t, то есть с частотой f >1/2Δ t, то в дискретизованном сигнале она наблюдаться не будет, но в нем останутся ложные периодические составляющие с периодами 2 T, 3 T и т.д. Таким образом, частоту (аргумент f спектральной плотности) имеет смысл анализировать только на интервале

[ -1/(2Δ t), 1/(2Δ t)].

Граничную частоту F = 1/(2Δ t) называют частотой Найквиста. Спектральную плотность в силу симметрии обычно рассматривают только на промежутке [0, 1/(2Δ t)], удваивая ее значения при положительных f для сохранения основных энергетических соотношений. Именно такая спектральная плотность измеряется на практике с помощью прямой фильтрации.

    Если для процесса x (t) заранее известно, что в нем отсутствуют периодические составляющие с частотами большими, чем некоторая величина F, то соответствующий дискретизованный процесс с шагом Δ t ≤1/(2 F) несет в себе всю информацию о частотном составе исходного x (t). Этот результат носит название теоремы Котельникова.

   Рассмотрим теперь ряд измерений { x_k }, из которого тренд, сезонная состав­ляющая и периодические составляющие либо исключены, либо их в нем не было. При этом, если речь идет о сравнении нескольких рядов, требуется, чтобы во всех них выде­ление этих составляющих было произведено одним и тем же способом. В большинстве случаев модель такого ряда { x_k } можно построить, рассматривая его как результат дискретизации некоей реализации стационарного эргодического случайного процесса   x (t):

x_k = x (k D t).

Еще одна важная возможность состоит в том, чтобы перейти в ряде { x_k } к конечным разностям:

y_k = x_{k + 1} - x_k; z_k = y_{k + 1} - y_k и т.д.

Разности достаточно высокого порядка не содержат медленно меняющихся составляющих и их обычно можно изучать методами, применяемыми для ста­ционарных рядов. Этот подход может использоваться в задачах распозна­вания, когда нужно создать некий обобщенный «портрет» ряда и сравни­вать его с «портретами» других рядов, прошедших такую же обработку. Кроме того, многие процедуры прогнозирования исходного ряда можно по­лучать из процедур прогнозирования рядов их конечных разностей.

    Ряд измерений { x_k }, моделью которого является стационарный эргодический случайный процесс с дискретным временем, называется стационарным временным рядом.

7.2. Быстрое преобразование Фурье (БПФ)

Для процесса с дискретным временем { x_k }, x_k = x (k Δ t), k =1,…, n сохраняются все введенные выше определения и результаты с той разницей, что интегралы заменяются суммами, а частоту (аргумент f спектральной плотности) имеет смысл рассматривать только на интервале [0 F_N ], где F_N = 1/(2Δ t) - частота Найквиста. В частности, прямое преобразование Фурье последовательности { x_k } – непараметрическую оценку спектральной плотности F (f) - можно вычислить по формуле

обратное – по формуле

Вычисления по каждой из этих формул требуют n ² операций. В 1966 г. Мак-Кован предложил геометрическую интерпретацию преобразования Фурье как некоторого поворота вектора x в n – мерном пространстве и показал, что его можно разложить на последовательность поворотов вокруг каждой из осей, что требует всего 4 n операций. Этот подход получил название быстрого преобразования Фурье, БПФ. Сейчас такие быстрые версии разработаны для всех основных интегральных преобразований. Благодаря технике БПФ при корреляционном анализе оказалось выгодным находить корреляционную функцию как обратное преобразование Фурье спектральной плотности, а не наоборот, как это делалось раньше. Пример вычисления БПФ в ИМС MatLab приведен в Документе 7.1 и проиллюстрирован на рис.7.1.

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

При анализе сигнала s (t) на фоне аддитивного шума ξ(t) большое значение имеет безразмерный параметр snr (signal to noise ratio) – отношение сигнал/шум. В случае дискретного времени его определяют как отношение дисперсий:

Эту величину часто выражают в децибелах (дБ):

дБ.

Документ 7.1. Построение спектральной плотности с помощью БПФ

Спектральная плотность гармонического сигнала с частотой f 0 на фоне аддитивного шума. Файл sp _ my _0. m clear; clc; del=0.005; Nyq=1/2/del; n=200; T=n*del; t=0:del:T; f0=30;     %Несущая частота (Гц) fi=pi/3/f0; s=0.5;     %СКО аддитивного шума y 0= sin (2* pi * f 0* t + fi); y=y0+s*randn(size(y0)); snr = sqrt (2)/0.5; %отношение сигнал/шум SNR =20* log (snr); %отношение сигнал/шум в Дб S 0= abs (fft (y)).^2; %спектральная плотность мощности f=linspace(0,Nyq,n/2+1); subplot(2,1,1); plot(t,y,'LineWidth',2); grid; title(' Синусоида с шумом ','FontName','Courier New Cyr'); xlabel(' Время (с)','FontName','Courier New Cyr'); ylabel(' Сигнал (мВ)','FontName','Courier New Cyr'); axis([0 T min(y)*1.1 max(y)*1.1]); subplot(2,1,2); plot(f,S0(1:n/2+1)*22'LineWidth',2); grid; title(' Спектральная плотность мощности ','FontName','Courier New Cyr'); xlabel(' Частота (Гц)','FontName','Courier New Cyr'); ylabel(' СПМ (мВ ^2)','FontName','Courier New Cyr'); axis([0 Nyq 0 max(S0)*1.1]);

 

    Для случайного ряда длиной n БПФ дает оценку СПМ S (f), значения которой представляют собой практически не коррелированные случайные величины с дисперсией порядка S ²(f), которая, очевидно, не убывает с ростом n. Основной выход состоит в том, чтобы вычислить оценки для набора интервалов, а затем усреднить. В ИМС MatLab широко используется функция spectrum, в которой реализован алгоритм Уэлча. В этой процедуре данный отрезок временного ряда разбивается на M (по умолчанию M =8) частично пересекающихся отрезков. В каждом отрезке удаляется линейный тренд, для каждого из остатков вычисляется БПФ с использованием одного из сглаживающих окон (по умолчанию используется окно Хеннинга) и результаты усредняются. Удаление линейного тренда приводит к небольшому смещению, поэтому результаты такого анализа могут несколько отличаться от результатов работы других спектроанализаторов.

Рис.7.1. СПМ синусоидального сигнала с частотой 30 Гц, наблюдаемого на фоне аддитивного шума. Отношение сигнал/шум snr = 2.8284=20.7944 дБ

Контрольные вопросы

1. Частота дискретизации

2. Частота Найквиста

3. Теорема Котельникова

4. Быстрое преобразование Фурье (БПФ)

5. БПФ смеси синусоид

6. Отношение «сигнал-шум». Децибелы

Задания на лабораторную работу № 7

1. Сформировать линейную систему с экспоненциально-убывающей весовой функцией

2. Сформировать выборку из стационарного процесса Юла и подать ее на вход системы

3. Получить процесс на выходе системы.

4. Оценить СПМ процессов на входе и выходе

5. Найти передаточную функцию системы

6. Проверить на примере выполнение основной теоремы теории линейных систем

ПАРАМЕТРИЧЕСКИЙ АНАЛИЗ