Обнаружьте искаженный сигнал в шуме

Присутствие шума часто мешает определять спектральное содержимое сигнала. Анализ частоты может помочь в таких случаях.

Рассмотрите, например, симулированный выход нелинейного усилителя, который вводит искажение третьего порядка.

Входной сигнал является амплитудной модулем синусоидой на 180 Гц, произведенной на уровне 3,6 кГц. Сгенерируйте 10 000 выборок.

N = 1e4;
n = 0:N-1;
fs = 3600;
f0 = 180;
t = n/fs;
y = sin(2*pi*f0*t);

Добавьте белый шум модульного отклонения во вход. Смоделируйте усилитель с помощью полинома третьего порядка. Передайте входной сигнал через усилитель с помощью polyval. Постройте раздел выхода. Поскольку сравнение строит выход чистой синусоиды.

rng default
noise = randn(size(y));

dispol = [0.5 0.75 1 0];
out = polyval(dispol,y+noise);

ns = 300:500; 

plot(t(ns),[out(ns);polyval(dispol,y(ns))])
xlabel('Time (s)')
ylabel('Signals')
axis tight
legend('With white noise','No white noise')

Figure contains an axes object. The axes object contains 2 objects of type line. These objects represent With white noise, No white noise.

Используйте pwelch вычислить и построить спектральную плотность мощности выхода.

[pxx,f] = pwelch(out,[],[],[],fs);

pwelch(out,[],[],[],fs)

Figure contains an axes object. The axes object with title Welch Power Spectral Density Estimate contains an object of type line.

Поскольку усилитель вводит искажение третьего порядка, выходной сигнал, как ожидают, будет иметь:

  • DC (нулевая частота) компонент;

  • Основной компонент с той же частотой как вход, 180 Гц;

  • Две гармоники - частотные составляющие в дважды и три раза частота входа, 360 и 540 Гц.

Проверьте, что выход как ожидалось для кубической нелинейности.

[pks,lox] = findpeaks(pxx,'NPeaks',4,'SortStr','descend');

hold on
plot(f(lox)/1000,10*log10(pks),'or')
hold off

legend('PSD','Frequency Components')

Figure contains an axes object. The axes object with title Welch Power Spectral Density Estimate contains 2 objects of type line. These objects represent PSD, Frequency Components.

components = sort([f(lox) f0*(0:3)'])'
components = 2×4

    0.8789  180.1758  360.3516  540.5273
         0  180.0000  360.0000  540.0000

pwelch работает путем деления сигнала на перекрывающиеся сегменты, вычисления периодограммы каждого сегмента и усреднения. По умолчанию функция использует восемь сегментов с 50%-м перекрытием. Для 10 000 выборок это соответствует 2 222 выборкам на сегмент.

Деление сигнала в более короткие сегменты приводит к большему количеству усреднения. Периодограмма более сглаженна, но имеет более низкое разрешение. Более высокую гармонику нельзя отличить.

pwelch(out,222,[],[],fs)

Figure contains an axes object. The axes object with title Welch Power Spectral Density Estimate contains an object of type line.

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

pwelch(out,4444,[],[],fs)

Figure contains an axes object. The axes object with title Welch Power Spectral Density Estimate contains an object of type line.

Смотрите также

|

Похожие темы