Перекрестный спектр и квадратичная когерентность

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

Создайте двухмерный временной ряд. Отдельные ряды состоят из двух синусоидальных волн с частотами 100 и 200 Гц. Серии встроены в дополнительный белый гауссов шум и дискретизированы на частоте 1 кГц. Синусоидальные волны в серии x имеют амплитуды, равные 1. Синусоидальная волна 100 Гц в y-ряду имеет амплитуду 0,5, а синусоидальная волна 200 Гц в y-ряду имеет амплитуду 0,35. Синусоидальные волны с частотой 100 Гц и 200 Гц в Y-рядах имеют фазовое отставание, соответственно, на $$\mathrm{\delta /4}$$ радиан и $$\mathrm{\lambda /2}$$ радиан. Y-ряд можно рассматривать как искаженный шумом выходной сигнал линейной системы с входом x. Установите генератор случайных чисел в значения по умолчанию для воспроизводимых результатов.

rng default

Fs = 1000;
t = 0:1/Fs:1-1/Fs;

x = cos(2*pi*100*t) + sin(2*pi*200*t) + 0.5*randn(size(t));
y = 0.5*cos(2*pi*100*t - pi/4) + 0.35*sin(2*pi*200*t - pi/2) + 0.5*randn(size(t));

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

Чтобы предотвратить получение оценки когерентности в квадрате величины, которая одинакова 1 для всех частот, необходимо использовать усредненную оценку когерентности. Подходящими являются методы усреднения перекрывающихся сегментов Welch (WOSA) и multitaper. mscohere реализует оценщик WOSA.

Установите длину окна равной 100 образцам. Эта длина окна содержит 10 периодов синусоидальной волны 100 Гц и 20 периодов синусоидальной волны 200 Гц. Используйте перекрытие 80 образцов в окне Хэмминга по умолчанию. Введите частоту дискретизации в явном виде, чтобы получить выходные частоты в Гц. Постройте график квадратичной когерентности. Квадратичная когерентность больше 0,8 при 100 и 200 Гц.

[Cxy,F] = mscohere(x,y,hamming(100),80,100,Fs);

plot(F,Cxy)
title('Magnitude-Squared Coherence')
xlabel('Frequency (Hz)')
grid

Получение перекрестного спектра x и y с помощью cpsd. Используйте те же параметры для получения перекрестного спектра, который использовался в оценке когерентности. Пренебрегайте перекрестным спектром, когда когерентность мала. Постройте график фазы перекрестного спектра и укажите частоты со значительной когерентностью между двумя раз. Отметьте известные фазовые задержки между синусоидальными компонентами. При частоте 100 Гц и 200 Гц фазовые запаздывания, оцененные из перекрестного спектра, близки к истинным значениям.

[Pxy,F] = cpsd(x,y,hamming(100),80,100,Fs);

Pxy(Cxy < 0.2) = 0;

plot(F,angle(Pxy)/pi)
title('Cross Spectrum Phase')
xlabel('Frequency (Hz)')
ylabel('Lag (\times\pi rad)')
grid

См. также

cpsd | mscohere | pwelch