Аналитический сигнал и преобразование Гильберта

hilbert функция находит точный аналитический сигнал для конечного блока данных. Можно также сгенерировать аналитический сигнал с помощью фильтра гильбертового трансформатора с конечной импульсной характеристикой (FIR) для вычисления аппроксимации к мнимой части.

Генерируют последовательность, состоящую из трех синусоид с частотами 203, 721 и 1001 Гц. Последовательность отбирают при частоте 10 кГц в течение примерно 1 секунды. Используйте hilbert функция для вычисления аналитического сигнала. Постройте график между 0,01 и 0,03 секундами.

fs = 1e4;
t = 0:1/fs:1; 

x = 2.5 + cos(2*pi*203*t) + sin(2*pi*721*t) + cos(2*pi*1001*t);

y = hilbert(x);

plot(t,real(y),t,imag(y))
xlim([0.01 0.03])
legend('real','imaginary')
title('hilbert Function')

Вычислите оценки Уэлча спектральных плотностей мощности исходной последовательности и аналитического сигнала. Разделите последовательности на неперекрывающиеся участки длиной 256 с помощью окна Хэмминга. Убедитесь, что аналитический сигнал не имеет мощности на отрицательных частотах.

pwelch([x;y].',256,0,[],fs,'centered')
legend('Original','hilbert')

Используйте designfilt функция для проектирования фильтра FIR трансформатора Гильберта 60-го порядка. Укажите ширину перехода 400 Гц. Визуализация частотной характеристики фильтра.

fo = 60;

d = designfilt('hilbertfir','FilterOrder',fo, ...
       'TransitionWidth',400,'SampleRate',fs); 

freqz(d,1024,fs)

Фильтрация синусоидальной последовательности для аппроксимации мнимой части аналитического сигнала.

hb = filter(d,x);

Групповая задержка фильтра, grd, равно половине порядка фильтра. Компенсировать эту задержку. Удалить первый grd образцы мнимой части и последней grd выборки действительной части и вектора времени. Постройте график результата от 0,01 секунды до 0,03 секунды.

grd = fo/2;
   
y2 = x(1:end-grd) + 1j*hb(grd+1:end);
t2 = t(1:end-grd);

plot(t2,real(y2),t2,imag(y2))
xlim([0.01 0.03])
legend('real','imaginary')
title('FIR Filter')

Оцените спектральную плотность мощности (PSD) приближенного аналитического сигнала и сравните ее с hilbert результат.

pwelch([y;[y2 zeros(1,grd)]].',256,0,[],fs,'centered')
legend('hilbert','FIR Filter')

См. также

designfilt | hilbert