Загрузка и визуализация нестационарного непрерывного сигнала, состоящего из синусоидальных волн с отчетливым изменением частоты. Дискретизация сигнала осуществляется на частоте 250 Гц.

fs = 250;
load nonstatdistinct
t = (0:length(nonstatdistinct)-1)/fs;
plot(t,nonstatdistinct)
xlabel('Time (s)')
ylabel('Signal')
axis tight

Использование ewt для получения мультирезолюционного анализа (MRA) сигнала.

mra = ewt(nonstatdistinct);

Используйте компоненты MRA с hht функция и график спектра Гильберта.

hht(mra,fs)

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

Создайте нестационарный непрерывный сигнал, состоящий из синусоидальных волн с отчетливым изменением частоты. Дискретизация сигнала осуществляется на частоте 1000 Гц.

Fs = 1000;
t = 0:1/Fs:4;
x1 = sin(2*pi*50*t) + sin(2*pi*200*t);
x2 = sin(2*pi*25*t) + sin(2*pi*100*t) + sin(2*pi*250*t);
x = [x1 x2] + 0.1*randn(1,length(t)*2);
t1 = (0:length(x)-1)/Fs;
plot(t1,x)
xlabel('Time (s)')
ylabel('Amplitude')
title('Signal')

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

[mra,~,wfb,info] = ewt(x);
Fs*info.PeakFrequencies

ans = 5×1

  249.9375
  200.0750
  100.1000
   50.1125
   25.1187

Fs*info.FilterBank.Passbands 

ans = 5×2

  223.6941  500.0000
  141.5896  223.6941
   70.8573  141.5896
   35.4911   70.8573
         0   35.4911

Постройте график спектра величины сигнала и банка фильтров. Местоположения peaks определяют полосы пропускания фильтра.

f = 0:Fs/length(x):Fs-1/length(x);
plot(f,wfb)
ylabel('Magnitude')
grid on
yyaxis right
plot(f,abs(fft(x)),'k--','linewidth',1.5)
ylabel('Magnitude')
xlabel('Hz')

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

Загрузка сигнала ЭКГ. Дискретизация сигнала производится на частоте 180 Гц.

load wecg

Использование ewt для получения мультирезолюционного анализа (MRA) сигнала и соответствующих коэффициентов анализа. Используйте четыре самых больших peaks, чтобы определить полосы пропускания фильтра.

mp = 4;
[mra,cfs] = ewt(wecg,'MaxNumPeaks',mp);

Постройте график сигнала и компонентов MRA.

fs = 180;
subplot(mp+1,1,1)
t = (0:length(wecg)-1)/fs;
plot(t,wecg)
title('MRA of Signal')
ylabel('Signal')
axis tight
for k=1:mp
    subplot(mp+1,1,k+1)
    plot(t,mra(:,k))
    ylabel(['MRA ',num2str(k)])
    axis tight
end
xlabel('Time (s)')

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

max(abs(wecg-sum(mra,2)))

ans = 8.8818e-16

Проверьте сохранение энергии коэффициентов EWT-анализа.

cfsenergy = sum(sum(abs(cfs).^2));
[cfsenergy norm(wecg,2)^2]

ans = 1×2

  298.2759  298.2759