Получить многоаппаратную PSD-оценку входного сигнала, состоящего из дискретно-временной синусоиды с угловой частотой $$\mathrm{\delta /4}$$ рад/образец с аддитивным N (0,1) белым шумом.

Создайте синусоидальную волну с угловой частотой $$\mathrm{\delta /4}$$ рад/образец с аддитивным N (0,1) белым шумом. Сигнал имеет длину 320 выборок. Получите многоаппаратную оценку PSD, используя произведение времени на полбандаж по умолчанию длиной 4 и DFT. По умолчанию количество точек DFT равно 512. Поскольку сигнал является действительным, оценка PSD является односторонней, и в оценке PSD имеются точки 512/2 + 1.

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

Постройте график многоаппаратной оценки PSD.

pmtm(x)

Сформировать 2048 выборок двухканального сигнала, внедренного в аддитивный N (0,1) белый гауссов шум.

Используйте метод multitaper, чтобы оценить PSD сигнала за 1024-интервал выборки от 0.1δ рад/образец до 0.4δ рад/образец. Используйте 13 синусоидальных конусов, взвешенных одинаково.

n = (0:2047)';

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

w = linspace(0.1,0.4,1024);

ntp = 13;
pmtm(x,ntp,'Tapers','sine',w*pi)

Повторите вычисления, но теперь весите 13 конусов в линейном порядке убывания. Вы можете разместить 'Tapers','sine' пара имя-значение в любом месте после x в вызове функции.

pmtm(x,(ntp:-1:1)/sum(1:ntp),w*pi,'Tapers','sine')

Повторите вычисления, но теперь используйте 13 Slepian сужений и укажите произведение полуполосы времени 7,5.

nw = 7.5;

pmtm(x,{nw,ntp},w*pi)

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

fs = 2e3;

pmtm(x,{nw,ntp},w*(fs/2),fs)

Получите многоаппаратную оценку PSD с указанным продуктом с полубандовой шириной.

Создайте синусоидальную волну с угловой частотой $$\mathrm{\delta /4}$$ рад/образец с аддитивным N (0,1) белым шумом. Сигнал имеет длину 320 выборок. Получите многоаппаратную оценку PSD с результатом полуполосы времени 2,5. Пропускная способность разрешения равна $$[\mathrm{-}\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}\mathrm{}2.5\lambda /320$$] рад/выборка. По умолчанию количество точек DFT равно 512. Поскольку сигнал является действительным, оценка PSD является односторонней, и в оценке PSD имеются точки 512/2 + 1.

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

Получить многоаппаратную PSD-оценку входного сигнала, состоящего из дискретно-временной синусоиды с угловой частотой $$\mathrm{\delta /4}$$ рад/образец с аддитивным N (0,1) белым шумом. Используйте длину DFT, равную длине сигнала.

Создайте синусоидальную волну с угловой частотой $$\mathrm{\delta /4}$$ рад/образец с аддитивным N (0,1) белым шумом. Сигнал имеет длину 320 выборок. Получите многоаппаратную оценку PSD с произведением полуполосы времени 3 и длиной DFT, равной длине сигнала. Поскольку сигнал имеет действительное значение, односторонняя оценка PSD возвращается по умолчанию с длиной, равной 320/2 + 1.

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

Получите многоаппаратную оценку PSD сигнала, дискретизированного на частоте 1 кГц. Сигнал представляет собой синусоидальную волну 100 Гц в аддитивном N (0,1) белом шуме. Длительность сигнала составляет 2 с. Используйте произведение полуполосы времени 3 и длины ДПФ, равное длине сигнала.

fs = 1000;
t = 0:1/fs:2-1/fs;
x = cos(2*pi*100*t)+randn(size(t));
[pxx,f] = pmtm(x,3,length(x),fs);

Постройте график многоаппаратной оценки PSD.

pmtm(x,3,length(x),fs)

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

Получите многоаппаратную оценку PSD сигнала, дискретизированного на частоте 1 кГц. Сигнал представляет собой синусоидальную волну 100 Гц в аддитивном N (0,1) белом шуме. Длительность сигнала составляет 2 с. Используйте произведение полуполосы времени, равное 3, и длину ДПФ, равную длине сигнала. Используйте 'unity' возможность придания равного веса в среднем каждой из отдельных конусных прямых спектральных оценок.

fs = 1000;
t = 0:1/fs:2-1/fs;
x = cos(2*pi*100*t)+randn(size(t));
[pxx,f] = pmtm(x,3,length(x),fs,'unity');

Постройте график многоаппаратной оценки PSD.

pmtm(x,3,length(x),fs,'unity')

В этом примере рассматриваются концентрации в частотной области последовательностей DPSS. В этом примере производится многоаппаратная PSD-оценка входного сигнала путем предварительного вычисления последовательностей Слепиана и выбора только тех, у которых более 99% энергии сосредоточено в разрешающей полосе.

Сигнал представляет собой синусоидальную волну 100 Гц в аддитивном N (0,1) белом шуме. Длительность сигнала составляет 2 секунды.

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

Установите для продукта с полуполосой пропускания значение 3.5. Для длины сигнала 2000 отсчетов и интервала дискретизации 0,001 секунд это приводит к ширине полосы разрешения [-1,75 1,75] Гц. Вычислите первые 10 слепиановых последовательностей и проверьте их частотные концентрации в заданной разрешающей способности.

[e,v] = dpss(length(x),3.5,10);
lv = length(v);

stem(1:lv,v,'filled')
ylim([0 1.2])
title('Proportion of Energy in [-w,w] of k-th Slepian Sequence')

Определите количество слеповских последовательностей с концентрациями энергии более 99%. Используя выбранные последовательности DPSS, получите многоаппаратную оценку PSD. Набор 'DropLastTaper' кому false для использования всех выбранных конусов.

hold on
plot([1 lv],0.99*[1 1])
hold off

idx = find(v>0.99,1,'last')

idx = 5

[pxx,f] = pmtm(x,e(:,1:idx),v(1:idx),length(x),fs,'DropLastTaper',false);

Постройте график многоаппаратной оценки PSD.

pmtm(x,e(:,1:idx),v(1:idx),length(x),fs,'DropLastTaper',false)

Получение многоаппаратной PSD-оценки синусоидальной волны 100 Гц в аддитивном шуме N (0,1). Данные дискретизируются при частоте 1 кГц. Используйте 'centered' возможность получения PSD с центром постоянного тока.

fs = 1000;
t = 0:1/fs:2-1/fs;
x = cos(2*pi*100*t)+randn(size(t));
[pxx,f] = pmtm(x,3.5,length(x),fs,'centered');

Постройте график оценки PSD с центром постоянного тока.

pmtm(x,3.5,length(x),fs,'centered')

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

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

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

Получите многоаппаратную оценку PSD с границами достоверности 95%. Постройте график оценки PSD вместе с доверительным интервалом и увеличьте масштаб интересующей области частот около 100 и 150 Гц.

[pxx,f,pxxc] = pmtm(x,3.5,length(x),fs,'ConfidenceLevel',0.95);

plot(f,10*log10(pxx))
hold on
plot(f,10*log10(pxxc),'r-.')
xlim([85 175])
xlabel('Hz')
ylabel('dB')
title('Multitaper PSD Estimate with 95%-Confidence Bounds')

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

Генерируют 1024 выборки многоканального сигнала, состоящего из трех синусоид в аддитивном $$N\mathrm{(}\mathrm{}0,1$$) белом гауссовом шуме. Частоты синусоид составляют $$\mathrm{}$$δ/2, $$\mathrm{}$$λ/3 и $$\mathrm{}$$λ/4 рад/образец. Оцените PSD сигнала с помощью многоаппаратного метода Томсона и постройте его график.

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

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

pmtm(x)