Документация

pxx = pwelch(x) возвращает оценку спектральной плотности степени (PSD), pxx, входного сигнала, x, найденный с использованием перекрытой оценки среднего сегмента Welch. Когда x является вектором, он рассматривается как один канал. Когда x является матрицей, PSD вычисляется независимо для каждого столбца и хранится в соответствующем столбце pxx. Если x является реальным, pxx является односторонней оценкой PSD. Если x является комплексным, pxx является двусторонней оценкой PSD. По умолчанию x разделяется на самые длинные возможные сегменты, чтобы получить как близкие, но не превышающие 8 сегментов с 50% перекрытием. Каждый сегмент окончается окном Хэмминга. Модифицированные периодограммы усредняются для получения оценки PSD. Если вы не можете разделить длину x точно в целое число сегментов с 50% перекрытием, x соответственно усечено.

pxx = pwelch(x,window) использует вектор входа или целое число, window, чтобы разделить сигнал на сегменты. Если window является вектором, pwelch разделяет сигнал на сегменты, равные по длине длине window. Измененные периодограммы вычисляются с помощью сегментов сигнала, умноженных на вектор, window. Если window является целым числом, сигнал делится на сегменты длины window. Измененные периодограммы вычисляются с помощью окна Хэмминга длины window.

pxx = pwelch(x,window,noverlap) использует noverlap выборки перекрытия от сегмента к сегменту. noverlap должно быть положительным целым числом, меньшим window если window - целое число. noverlap должно быть положительным целым числом, меньше длины window если window является вектором. Если вы не задаете noverlap, или задайте noverlap пустой, количество перекрывающихся выборок по умолчанию составляет 50% от длины окна.

pxx = pwelch(x,window,noverlap,nfft) задает количество дискретных точек преобразования Фурье (DFT) для использования в оценке PSD. Значение по умолчанию nfft больше 256 или следующая степень 2 больше, чем длина сегментов.

[pxx,w] = pwelch(___) возвращает нормированный вектор частоты, w. Если pxx является односторонней оценкой PSD, w охватывает интервал [0, π] если nfft даже и [0, π), если nfft нечетно. Если pxx является двусторонней оценкой PSD, w охватывает интервал [0,2в).

[pxx,f] = pwelch(___,fs) возвращает вектор частоты, f, в циклах в единицу времени. Частота дискретизации, fs, количество выборок в единицу времени. Если модулем времени является секунды, то f в циклах/с (Гц). Для реальных сигналов, f охватывает интервал [0, fs/ 2] когда nfft является четным и [0, fs/ 2) при nfft нечетно. Для комплексных сигналов, f охватывает интервал [0, fs). fs должен быть пятым входом, pwelch. Чтобы ввести частоту дискретизации и все еще использовать значения по умолчанию предыдущих необязательных аргументов, задайте эти аргументы как пустые [].

[pxx,w] = pwelch(x,window,noverlap,w) возвращает двусторонние оценки PSD Welch на нормализованных частотах, заданных в векторе, w. Векторная w должен содержать как минимум два элемента, поскольку в противном случае функция интерпретирует его как nfft.

[pxx,f] = pwelch(x,window,noverlap,f,fs) возвращает двусторонние оценки PSD Welch на частотах, заданных в векторе, f. Векторная f должен содержать как минимум два элемента, поскольку в противном случае функция интерпретирует его как nfft. Частоты в f в циклах в единицах времени. Частота дискретизации, fs, количество выборок в единицу времени. Если модулем времени является секунды, то f в циклах/с (Гц).

[___] = pwelch(x,window,___,freqrange) возвращает оценку PSD Welch в частотной области значений, заданном freqrange. Допустимые опции для freqrange являются: 'onesided', 'twosided', или 'centered'.

[___] = pwelch(x,window,___,trace) возвращает оценку спектра с максимальным удержанием, если trace задается как 'maxhold' и возвращает оценку спектра с минимальным удержанием, если trace задается как 'minhold'.

[___,pxxc] = pwelch(___,'ConfidenceLevel',probability) возвращает probability × 100% доверительные интервалы для оценки PSD в pxxc.

[___] = pwelch(___,spectrumtype) возвращает оценку PSD, если spectrumtype задается как 'psd' и возвращает спектр степени, если spectrumtype задается как 'power'.

pwelch(___) без выходных аргументов строит график оценки PSD Welch в текущую фигуру окне.

Примеры

Оценка Уэлча с использованием входных параметров по умолчанию

Получите оценку PSD Уэлча входного сигнала, состоящего из синусоиды дискретного времени с угловой частотой $π / 4$ рад/образец с добавкой $N (0, 1)$ белый шум.

Создайте синусоиду с угловой частотой $π / 4$ рад/образец с добавкой $N (0, 1)$ белый шум. Сбросьте генератор случайных чисел для воспроизводимых результатов. Сигнал имеет длину $N_{x} = 320$ выборки.

rng default

n = 0:319;
x = cos(pi/4*n)+randn(size(n));

Получите оценку PSD Welch с помощью окна Хэмминга по умолчанию и длины ДПФ. Длина сегмента по умолчанию является 71 выборкой, и длина ДПФ является 256 точками, получая разрешение частоты $2 π / 256$ рад/образец. Поскольку сигнал является реальным, периодограмма является односторонней, и существует 256/2 + 1 точек. Построение графика подсчета PSD по сварке.

pxx = pwelch(x);

pwelch(x)

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

Повторите расчет.

Разделите сигнал на участки длины $n s c = ⌊ N_{x} / 4.5 ⌋$ . Это действие эквивалентно делению сигнала на самые длинные возможные сегменты, чтобы получить как близкие, но не превышающие 8 сегментов с 50% перекрытием.
Откройте разделы с помощью окна Хэмминга.
Задайте 50% перекрытия между смежными сечениями
Для вычисления БПФ используйте $\max (256, 2^{p})$ точки, где $p = ⌈ \log_{2} n s c ⌉$ .

Проверьте, что эти два подхода дают одинаковые результаты.

Nx = length(x);
nsc = floor(Nx/4.5);
nov = floor(nsc/2);
nff = max(256,2^nextpow2(nsc));

t = pwelch(x,hamming(nsc),nov,nff);

maxerr = max(abs(abs(t(:))-abs(pxx(:))))

maxerr = 0

Разделите сигнал на 8 секций равной длины с 50% перекрытием между секциями. Задайте ту же длину БПФ, что и на предыдущем шаге. Вычислите оценку PSD Welch и проверьте, что она дает тот же результат, что и две предыдущие процедуры.

ns = 8;
ov = 0.5;
lsc = floor(Nx/(ns-(ns-1)*ov));

t = pwelch(x,lsc,floor(ov*lsc),nff);

maxerr = max(abs(abs(t(:))-abs(pxx(:))))

maxerr = 0

Оценка Welch с использованием заданной длины сегмента

Получите оценку PSD Уэлча входного сигнала, состоящего из синусоиды дискретного времени с угловой частотой $π / 3$ рад/образец с добавкой $N (0, 1)$ белый шум.

Создайте синусоиду с угловой частотой $π / 3$ рад/образец с добавкой $N (0, 1)$ белый шум. Сбросьте генератор случайных чисел для воспроизводимых результатов. Сигнал имеет 512 выборки.

rng default

n = 0:511;
x = cos(pi/3*n)+randn(size(n));

Получите оценку PSD Welch, разделяющую сигнал на сегменты 132 выборки по длине. Сегменты сигнала умножаются на Окно Хэмминга 132 выборок длины. Количество перекрывающихся выборок не задано, поэтому установлено на 132/2 = 66. Длина ДПФ составляет 256 точек, что приводит к частотному разрешению $2 π / 256$ рад/образец. Поскольку сигнал является реальным, оценка PSD является односторонней, и существует 256/2 + 1 = 129 точек. Постройте график PSD как функции нормированной частоты.

segmentLength = 132;
[pxx,w] = pwelch(x,segmentLength);

plot(w/pi,10*log10(pxx))
xlabel('\omega / \pi')

Figure contains an axes. The axes contains an object of type line.

Оценка Welch, определяющая перекрытие сегмента

rng default

n = 0:319;
x = cos(pi/4*n)+randn(size(n));

Получите оценку PSD Welch, разделяющую сигнал на сегменты 100 выборки в длину. Сегменты сигнала умножаются на Окно Хэмминга 100 выборок в длину. Количество перекрываемых выборок составляет 25. Длина ДПФ составляет 256 точек, что приводит к частотному разрешению $2 π / 256$ рад/образец. Поскольку сигнал является реальным, оценка PSD является односторонней, и существует 256/2 + 1 точек.

segmentLength = 100;
noverlap = 25;
pxx = pwelch(x,segmentLength,noverlap);

plot(10*log10(pxx))

Figure contains an axes. The axes contains an object of type line.

Оценка Welch с использованием указанной длины ДПФ

rng default

n = 0:319;
x = cos(pi/4*n) + randn(size(n));

Получите оценку PSD Welch, разделяющую сигнал на сегменты 100 выборки в длину. Используйте перекрытие по умолчанию 50%. Задайте длину ДПФ, чтобы она составляла 640 точек, так чтобы частота $π / 4$ рад/образец соответствует интервалу ДПФ (интервал 81). Поскольку сигнал является реальным, оценка PSD является односторонней, и существует 640/2 + 1 точек.

segmentLength = 100;
nfft = 640;
pxx = pwelch(x,segmentLength,[],nfft);

plot(10*log10(pxx))
xlabel('rad/sample')
ylabel('dB / (rad/sample)')

Figure contains an axes. The axes contains an object of type line.

Оценка PSD сигнала Welch с частотой в Герце

Создайте сигнал, состоящий из синусоиды 100 Гц в аддитиве N (0,1) белого шума. Сбросьте генератор случайных чисел для воспроизводимых результатов. Частота дискретизации составляет 1 кГц, и сигнал составляет 5 секунд в длительность.

rng default

fs = 1000;
t = 0:1/fs:5-1/fs;
x = cos(2*pi*100*t) + randn(size(t));

Получите оценку PSD перекрытого сегмента Welch предыдущего сигнала. Используйте длину сегмента 500 выборок с 300 перекрываемыми выборками. Используйте 500 точек ДПФ так, чтобы 100 Гц падали непосредственно на интервал ДПФ. Введите частоту дискретизации, чтобы вывести вектор частот в Гц. Постройте график результата.

[pxx,f] = pwelch(x,500,300,500,fs);

plot(f,10*log10(pxx))

xlabel('Frequency (Hz)')
ylabel('PSD (dB/Hz)')

Figure contains an axes. The axes contains an object of type line.

Спектры с максимальным и минимальным задержками

Создайте сигнал, состоящий из трех шумных синусоидов и щебета, дискретизированный при 200 кГц в течение 0,1 секунды. Частоты синусоидов составляют 1 кГц, 10 кГц и 20 кГц. Синусоиды имеют различные амплитуды и уровни шума. Бесшумный щебет имеет частоту, которая начинается с 20 кГц и увеличивается линейно до 30 кГц во время дискретизации.

Fs = 200e3; 
Fc = [1 10 20]'*1e3; 
Ns = 0.1*Fs;

t = (0:Ns-1)/Fs;
x = [1 1/10 10]*sin(2*pi*Fc*t)+[1/200 1/2000 1/20]*randn(3,Ns);
x = x+chirp(t,20e3,t(end),30e3);

Вычислите оценку PSD Welch и спектры максимального удержания и минимального удержания сигнала. Постройте график результатов.

[pxx,f] = pwelch(x,[],[],[],Fs);
pmax = pwelch(x,[],[],[],Fs,'maxhold');
pmin = pwelch(x,[],[],[],Fs,'minhold');

plot(f,pow2db(pxx))
hold on
plot(f,pow2db([pmax pmin]),':')
hold off
xlabel('Frequency (Hz)')
ylabel('PSD (dB/Hz)')
legend('pwelch','maxhold','minhold')

Figure contains an axes. The axes contains 3 objects of type line. These objects represent pwelch, maxhold, minhold.

Повторите процедуру, на этот раз вычисляя оценки центрированного спектра степени.

[pxx,f] = pwelch(x,[],[],[],Fs,'centered','power');
pmax = pwelch(x,[],[],[],Fs,'maxhold','centered','power');
pmin = pwelch(x,[],[],[],Fs,'minhold','centered','power');

plot(f,pow2db(pxx))
hold on
plot(f,pow2db([pmax pmin]),':')
hold off
xlabel('Frequency (Hz)')
ylabel('Power (dB)')
legend('pwelch','maxhold','minhold')

Figure contains an axes. The axes contains 3 objects of type line. These objects represent pwelch, maxhold, minhold.

Верхняя и нижняя границы 95% -доверие

Этот пример иллюстрирует использование доверия границ с оценкой PSD перекрытого сегмента (WOSA) Welch. Хотя не является необходимым условием для статистической значимости, частоты в оценке Уэлча, где нижняя доверительная граница превышает верхнюю доверительную границу для окружающих оценок PSD, четко указывают на значительные колебания во временных рядах.

Создайте сигнал, состоящий из суперпозиции 100 Гц и 150 Гц синусоид аддитивного белого N (0,1) шума. Амплитуда двух синусоид равна 1. Частота дискретизации составляет 1 кГц. Сбросьте генератор случайных чисел для воспроизводимых результатов .

rng default
fs = 1000;
t = 0:1/fs:1-1/fs;
x = cos(2*pi*100*t)+sin(2*pi*150*t)+randn(size(t));

Получите оценку WOSA с 95% -доверием. Установите длину сегмента равную 200 и перекройте сегменты на 50% (100 выборки). Постройте график оценки PSD WOSA вместе с интервалом доверия и увеличьте частоту необходимой области около 100 и 150 Гц.

L = 200;
noverlap = 100;
[pxx,f,pxxc] = pwelch(x,hamming(L),noverlap,200,fs,...
    'ConfidenceLevel',0.95);

plot(f,10*log10(pxx))
hold on
plot(f,10*log10(pxxc),'-.')
hold off

xlim([25 250])
xlabel('Frequency (Hz)')
ylabel('PSD (dB/Hz)')
title('Welch Estimate with 95%-Confidence Bounds')

Figure contains an axes. The axes with title Welch Estimate with 95%-Confidence Bounds contains 3 objects of type line.

Нижняя доверительная граница в непосредственной близости от 100 и 150 Гц значительно выше верхней доверительной границы за пределами 100 и 150 Гц.

Степень с постоянным током

Создайте сигнал, состоящий из синусоиды 100 Гц в аддитиве $N (0, 1 / 4)$ белый шум. Сбросьте генератор случайных чисел для воспроизводимых результатов. Частота дискретизации составляет 1 кГц, и сигнал составляет 5 секунд в длительность.

rng default

fs = 1000;
t = 0:1/fs:5-1/fs;

noisevar = 1/4;
x = cos(2*pi*100*t)+sqrt(noisevar)*randn(size(t));

Получите спектр степени с постоянным током с помощью метода Уэлча. Используйте длину сегмента 500 выборок с 300 перекрываемыми выборками и длину ДПФ 500 точек. Постройте график результата.

[pxx,f] = pwelch(x,500,300,500,fs,'centered','power');

plot(f,10*log10(pxx))
xlabel('Frequency (Hz)')
ylabel('Magnitude (dB)')
grid

Figure contains an axes. The axes contains an object of type line.

Вы видите, что степень на -100 и 100 Гц близка к ожидаемой степени 1/4 для реальной синусоиды с амплитудой 1. Отклонение от 1/4 обусловлено эффектом аддитивного шума.

Оценка PSD Welch многоканального сигнала

Сгенерируйте 1024 выборки многоканального сигнала, состоящего из трех синусоидов в аддитиве $N (0, 1)$ белый Гауссов шум. Частоты синусоидов $π / 2$ , $π / 3$ , и $π / 4$ рад/образец. Оцените PSD сигнала с помощью метода Уэлча и постройте график.

N = 1024;
n = 0:N-1;

w = pi./[2;3;4];
x = cos(w*n)' + randn(length(n),3);

pwelch(x)

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

Входные параметры

`x` - Входной сигнал
вектор | матрица

Входной сигнал, заданный как строка или вектор-столбец, или как матрица. Если x является матрицей, тогда ее столбцы рассматриваются как независимые каналы.

Пример: cos(pi/4*(0:159))+randn(1,160) является одноканальным вектором-строкой.

Пример: cos(pi./[4;2]*(0:159))'+randn(160,2) является двухканальным сигналом.

Типы данных: single | double
Поддержка комплексного числа: Да

`window` - Окно
целое число | вектор | `[]`

Окно, заданное как строка или столбец или целое число. Если window является вектором, pwelch делит x в перекрывающиеся сегменты длины, равной длине window, а затем умножает каждый сегмент сигнала на вектор, заданный в window. Если window - целое число, pwelch разделяется на сегменты длины, равные целому значению, и используется окно Хэмминга равной длины. Если длина x нельзя разделить точно на целое число сегментов с noverlap количество перекрываемых выборок, x соответственно усечено. Если вы задаете window пустой, окно Хемминга по умолчанию используется для получения восьми сегментов x с noverlap перекрываемые выборки.

Типы данных: single | double

`noverlap` - Количество перекрываемых выборок
положительное целое число | `[]`

Количество перекрывающихся выборок, заданное в виде положительного целого числа, меньшего, чем длина window. Если вы опускаете noverlap или задайте noverlap пустое значение используется для получения 50% перекрытия между сегментами.

`nfft` - Количество точек ДПФ
`max(256,2^nextpow2(length(window)))` (по умолчанию) | целое число | `[]`

Количество точек ДПФ, заданное как положительное целое число. Для действительного входного сигнала, x, оценка PSD, pxx имеет длину (nfft/ 2 + 1), если nfft является четным, и (nfft + 1 )/2, если nfft нечетно. Для комплексного входного сигнала, xоценка PSD всегда имеет длину nfft. Если nfft задан как пустой, значение по умолчанию nfft используется.

Если nfft больше длины сегмента, данные заполнены нулями. Если nfft меньше длины сегмента, сегмент обернут с помощью datawrap чтобы сделать длину равной nfft.

Типы данных: single | double

`fs` - Частота дискретизации
положительная скалярная величина

Частота дискретизации, заданная как положительная скалярная величина. Частота дискретизации является количеством выборок в единицу времени. Если модулем времени является секунды, то частота дискретизации имеет модули измерения Гц.

`w` - Нормированные частоты
вектор

Нормированные частоты, заданные как строка или вектор-столбец с по крайней мере двумя элементами. Нормированные частоты указаны в рад/отсчет.

Пример: w = [pi/4 pi/2]

Типы данных: double

`f` - Частоты
вектор

Частоты, заданные как строка или вектор-столбец с по крайней мере двумя элементами. Частоты указаны в циклах в единицах времени. Единичное время задается частотой дискретизации, fs. Если fs имеет модули измерения проб/секунду, затем f содержит модули Гц.

Пример: fs = 1000; f = [100 200]

Типы данных: double

`freqrange` - Частотная область значений для оценки PSD
`'onesided'` | `'twosided'` | `'centered'`

Частотная область значений для оценки PSD, заданный как один из 'onesided', 'twosided', или 'centered'. Значение по умолчанию является 'onesided' для реальных сигналов и 'twosided' для комплексных сигналов. Частотные области значений, соответствующие каждой опции,

'onesided' - возвращает одностороннюю оценку PSD действительного входного сигнала, x. Если nfft является четным, pxx имеет длину nfft/ 2 + 1 и вычисляется через интервал [0, π] рад/отсчет. Если nfft нечетно, длина pxx is (nfft + 1 )/2 и интервал равен [0, π) рад/отсчету. Когда fs опционально задан, соответствующие интервалы [0, fs/ 2] циклов/единичное время и [0, fs/ 2) циклы/единичное время для четной и нечетной длины nfft соответственно.
'twosided' - возвращает двустороннюю оценку PSD для реального или комплексного входов, x. В этом случае pxx имеет длину nfft и вычисляется через интервал [0,2 π) рад/отсчета. Когда fs опционально задан, интервал равен [0, fs) циклы/единичное время.
'centered' - возвращает центрированную двустороннюю оценку PSD для реального или комплексного входа, x. В этом случае pxx имеет длину nfft и вычисляется по интервалу (- π, π] рад/отсчет для четной длины nfft и (- π, π) рад/отсчет для нечетной длины nfft. Когда fs опционально задано, соответствующие интервалы являются (- fs/2, fs/ 2] циклов/единичное время и (- fs/2, fs/ 2) циклы/единичное время для четной и нечетной длины nfft соответственно.

`spectrumtype` - Масштабирование спектра степени
`'psd'` (по умолчанию) | `'power'`

Масштабирование спектра степени, заданное как одно из 'psd' или 'power'. Опускание spectrumtype, или определение 'psd', возвращает спектральную плотность степени. Определение 'power' масштабирует каждую оценку PSD по эквивалентной шумовой полосе окна. Используйте 'power' опция для получения оценки степени на каждой частоте.

`trace` - Режим трассировки
`'mean'` (по умолчанию) | `'maxhold'` | `'minhold'`

Режим трассировки, заданный как один из 'mean', 'maxhold', или 'minhold'. Значение по умолчанию является 'mean'.

'mean' - возвращает оценку спектра Welch для каждого входного канала. pwelch вычисляет оценку спектра Уэлча в каждом частотном интервале путем усреднения оценок спектральной степени всех сегментов.
'maxhold' - возвращает максимально удерживаемый спектр каждого входного канала. pwelch вычисляет максимально удерживаемый спектр в каждом частотном интервале путем сохранения максимального значения среди оценок спектральной степени всех сегментов.
'minhold' - возвращает минимально удерживаемый спектр каждого входного канала. pwelch вычисляет спектр с минимальным удержанием в каждом частотном интервале путем сохранения минимального значения среди оценок спектра степени всех сегментов.

`probability` - Доверительный интервал для оценки PSD
`0.95` (по умолчанию) | скаляром в области значений (0,1)

Вероятность покрытия для истинного PSD, заданная в виде скаляра в области значений (0,1). Выход, pxxc, содержит нижнюю и верхнюю границы probability × оценку 100% интервала для истинного PSD.

Выходные аргументы

`pxx` - Оценка PSD
вектор | матрица

Оценка PSD, возвращенная как действительный, неотрицательный вектор-столбец или матрица. Каждый столбец pxx - оценка PSD соответствующего столбца x. Модули оценки PSD указаны в квадрате величины модулей из данных временных рядов на модуль частоту. Для примера, если входные данные в вольтах, оценка PSD в модули от квадратов В на модуль частоту. Для временных рядов в вольтах, если вы принимаете сопротивление 1 Ом и задаете частоту дискретизации в герцах, оценка PSD в ваттах на герц.

Типы данных: single | double

`w` - Нормированные частоты
вектор

Нормированные частоты, возвращенные как реальный вектор-столбец. Если pxx является односторонней оценкой PSD, w охватывает интервал [0, π], если nfft чётный и [0, π), если nfft нечетно. Если pxx является двусторонней оценкой PSD, w охватывает интервал [0,2 π). Для оценки PSD с центром DC, w охватывает интервал (- π, π] для четных nfft и (- π, π) для нечетных nfft.

Типы данных: double

`f` - Циклические частоты
вектор

Циклические частоты, возвращенные как реальный вектор-столбец. Для односторонней оценки PSD, f охватывает интервал [0, fs/ 2] когда nfft является четным и [0, fs/ 2) при nfft нечетно. Для двусторонней оценки PSD, f охватывает интервал [0, fs). Для оценки PSD с центром DC, f охватывает интервал (- fs/2, fs/ 2] циклов/единичное время для четной длины nfft и (- fs/2, fs/ 2) циклы/единичное время для нечетной длины nfft.

Типы данных: double | single

`pxxc` - Доверительные границы
матрица

Доверительные границы, возвращенные как матрица с вещественными элементами. Размер строки матрицы равен длине оценки PSD, pxx. pxxc имеет в два раза больше столбцов, чем pxx. Столбцы с нечетным номером содержат нижние границы доверительных интервалов, а столбцы с четным номером - верхние границы. Таким образом, pxxc(m,2*n-1) - нижняя доверительная граница и pxxc(m,2*n) - верхняя доверительная граница, соответствующая оценке pxx(m,n). Вероятность покрытия доверительных интервалов определяется значением probability вход.

Типы данных: single | double

Подробнее о