AnyMath 文档

普威尔奇

Welch功率的谱密度的估计。

库::`工程师`

语法

函数调用

  • pwelch（___) -绘制SPM或Welch功率谱的估计。

争论

输入参数

# x-输入信号

+ 向量资料 | 矩阵

Details

指定为矢量或矩阵的输入信号。 如果 x -矩阵，函数 pwelch 将其列视为独立通道。

数据类型

漂浮物32, 漂浮64, Int64</无翻译> 支持复数::是

# 赢-窗口

+ 空数组 （默认情况下）| 一个整数 | 向量资料

Details

设置为空数组的窗口 地板32[], 浮标64[] 或 Int64[]，正整数或向量。

功能 pwelch 分割输入信号 [参数:x] 划分为段并应用窗口以将每个段的计数递增地乘以窗口值。 窗口值取决于此参数的设置方式。:

  • 地板32[], 浮标64[] 或 Int64[] -该函数应用这样长度的汉明窗口进行分割 [参数:x] 上 8 与重叠的段 [参数:noverlap] 计数;

  • 一个正整数 -功能划分 [参数:x] 成段长度 赢 并应用相同长度的汉明窗;

  • 向量资料 -功能划分 [参数:x] 成长度等于向量长度的段 赢，并应用参数中指定的窗口 赢.

如果长度为 [参数:x] 不允许您获得整数个重叠段 [参数:noverlap] 计数，功能 pwelch 相应地截断信号 [参数:x].

有关可用窗口的列表，请参阅库 光谱分析.

默认情况下，汉明窗口的旁瓣在 42.5 dB，因此频谱分量低于该水平（相对于最大分量）它们可以被隐藏。 选择不同的窗口允许您更改分辨率（例如，使用矩形窗口）和旁瓣衰减（例如，使用矩形窗口）之间的比率。
数据类型

漂浮物32, 漂浮64, Int64</无翻译>

# nOverlap-重叠计数的数量

+ 空数组 (默认)| 一个整数

Details

设置为空数组的重叠样本数 地板32[], 浮标64[] 或 Int64[] 或非负整数。

  • 地板32[], 浮标64[] 或 Int64[] -该函数选择一个值，以确保段重叠 50%.

    如果未指定段长度，则函数设置此参数的值。 ，在哪里 -输入信号的长度，以及 表示向下取整函数。 这个动作相当于将信号分成尽可能长的段，以便得到尽可能多的，但不超过八个重叠的段。 50%.

  • 一个非负整数 -相邻段按此参数中指定的样本数重叠。

**如果 [参数：赢] -标量，值 诺弗拉普 它应该少一些 [参数：赢].

**如果 [参数：赢] -向量，值 诺弗拉普 必须短于长度 [参数：赢].

# freqSpec-频率响应

+ 空数组 （默认情况下）| 一个整数 | 向量资料

Details

以下列方式之一指定的频率响应:

  • 地板32[], 浮标64[] 或 Int64[] -该函数使用等于最大值的DFT点的数量计算频谱估计值 256 和二的最近幂，等于或大于窗口长度;

  • 一个正整数 -函数使用此自变量中指定的DFT点数计算光谱估计值:

如果 freqSpec 更长的段长，功能 pwelch 向段数据添加零; 如果 freqSpec 小于段长，函数对段进行模分割 freqSpec 并总结所得片段;

  • 向量资料，包含至少两个元素，—该函数计算在该参数中指定的频率下的频谱估计值:

如果设置了采样率 [参数:fs]，函数考虑值 freqSpec 作为循环频率在相同的单位 [参数:fs]; 如果未设置采样频率，则函数假定值 freqSpec -这些是rad/count中的归一化频率。

数据类型

漂浮物32, 漂浮64, Int64</无翻译>

# Fs-采样 费率+ 标量,标量

Details

的采样率，设置为正标量。 它等于每单位时间的计数次数。 如果时间以秒为单位，则采样率以Hz表示。

如果参数 财政司司长 设置为 地板32[], 浮标64[] 或 Int64[]，功能 *pwelch 接受输入信号的采样频率 [参数:x] 平等 1 赫兹。 如果参数 财政司司长 未指定，功能 *pwelch 接受采样频率等于 2π rad/倒计时。

# freqRange-用于估计SPM或功率谱的频率范围

+ "onesided" | "twosided" | "居中"

Details

SPM估计值或功率谱的频率范围，定义为 "onesided", "twosided" 或 "居中".

功能 pwelch 申报表 [参数:pxx] 随着行数和频率范围取决于值 freqRange,freqRange，dft点数的奇偶性在 [参数：freqspec] 这取决于参数是否设置。 [参数:fs].

freqRange

freqSpec

行数 pxx

频率范围 pxx

财政司司长 未指定

财政司司长 设置

"onesided" （默认情况下，如果 [参数:x] -一个实数）

偶数

freqSpec/2 + 1

[0, π] rad/倒计时

[0，Fs/2] 每单位时间的周期

奇数

（freqSpec + 1) / 2

[0, π) rad/倒计时

[0，Fs/2) 每单位时间的周期

"twosided" （默认情况下，如果 x -复数）

偶数或奇数

freqSpec

[0, 2π) rad/倒计时

[0，Fs) 每单位时间的周期

"居中"

偶数

freqSpec

(−π, π] rad/倒计时

(-Fs/2,Fs/2] 每单位时间的周期

奇数

(−π, π) rad/倒计时

(-Fs/2,Fs/2) 每单位时间的周期

以下情况不支持此参数 [参数：freqspec] 它被定义为循环或归一化频率的向量。 *对于值 "onesided" 功能 *pwelch 将功率乘以 2 在除零频率和奈奎斯特频率之外的所有频率上保持全功率。 *价值 "onesided" 不支持复杂值 [参数:x].
数据类型

字符串</无翻译>

# spectrumType-功率谱的缩放

+ "psd" (默认)| "权力"

Details

功率谱的缩放，设置为以下值之一:

  • "psd" -功能 pwelch 返回频谱功率密度;

  • "权力" -该函数将每个SPM估计值按窗口的等效噪声频带缩放，并返回每个频率处功率的估计值。

该表显示了SPM估计值与返回到的功率谱的估计值之间的缩放比例 [参数:pxx]，与指定的输入信号 [参数:x]，窗口向量 [参数：赢]，重叠的长度 [参数:noverlap]、频率响应 [参数：freqspec] 和采样率 [参数:fs].

采样率 缩放比例

财政司司长 套装

 
psd = pwelch(x, win, nOverlap, freqSpec, Fs, "psd")
pow = pwelch(x, win, nOverlap, freqSpec, Fs, "power")

然后 战俘 相当于 psd*enbw(win,Fs)，在哪里 赢 -这是窗口向量。

财政司司长 未指定

 
psd = pwelch(x, win, nOverlap, freqSpec, "psd")
pow = pwelch(x, win, nOverlap, freqSpec, "power")

然后 战俘 相当于 psd*enbw(win,2*π).

# 跟踪-跟踪模式

+ "卑鄙" (默认)| "马克斯霍尔德" | "minhold"

Details

跟踪模式，由以下值之一设置:

  • "卑鄙" -功能 pwelch 通过平均每个频率范围内所有段的功率谱的估计值，返回每个输入通道的功率谱的估计值;

  • "马克斯霍尔德" -该函数返回每个输入通道的最大保留频谱，在每个频率范围内按段选择功率谱的最大值;

  • "minhold" -该函数返回每个输入通道的最小保留谱，在每个频率范围内按段选择功率谱的最小值。

输入参数 «名称-值»

将可选参数对指定为 名称，值，在哪里 姓名 -参数的名称，以及 价值 -适当的值。 类型参数 «名称-值» 它们应该放在其他参数之后，但对的顺序无关紧要。 您可以指定多个对 «名称-值».

使用逗号分隔名称和值，以及 姓名 将其用引号括起来，或者使用等号分隔名称和值，并且 姓名 指定不带引号。

# ConfidenceLevel-SPM评估的置信区间

+ 0.95 （默认情况下）| 标量,标量

Details

SPM估计的复盖概率，以范围内的标量给出 (0, 1).

如果指定 ConfidenceLevel=p （或 "ConfidenceLevel"，p), 功能 pwelch 在输出参数中返回 [参数:pxxc] 置信区间与水平的下限和上限 p×100%为真SPM。

# Parent-目标父容器

+ 空数组 (默认)/对象 轴心 /对象 事务委员会

Details

目标父容器，设置为对象 轴心 或 事务委员会.

如果指定 父=h （或 "家长"，h），功能 pwelch 在指定的目标父容器上绘制SPM或Welch功率谱的估计，而不考虑是否存在输出参数。

# out-输出数据的类型

+ :剧情 （默认情况下）| :数据 | :全部

Details

输出数据的类型:

  • :剧情 -函数返回一个图形;

  • :数据 -函数返回数据;

  • :全部 -函数返回数据和图形。

对于此参数，名称和值用等号(=)分隔。

输出参数

# pxx-SPM分数或功率谱

+ 向量资料 | 矩阵

Details

作为实非负列向量或矩阵返回的SPM估计值或功率谱。

*每列 pxx的 它表示SPM或相应列的功率谱的估计。 [参数:x] 取决于参数的设置方式 [参数：光谱类型]. *SPM估计值的测量单位是每单位频率的输入信号值的平方。 例如，如果输入数据是 [参数:x] 以伏特为单位，采样率 [参数:fs] -在赫兹里，阻力是 1 欧姆，那么SPM估计值将以瓦特/赫兹表示。

*功率谱的测量单位是输入信号的平方。 例如，如果输入数据是 [参数:x] 它们以伏特为单位设置，电阻为 1 欧姆，那么SPM估计将以瓦特表示。

# f-与SPM评估相关的频率

+ 向量资料

Details

与SPM评估相关联的频率作为实列向量返回。

*如果参数 [参数:fs] 套装, f 它包含以Hz为单位的循环频率。 *如果参数 [参数:fs] 未指定, f 包含以rad/count为单位的归一化频率。

# pxxc-置信区间

+ 矩阵

Details

作为实矩阵返回的置信区间的边界。

*矩阵 pxxc 有尽可能多的行 [参数:pxx]. *矩阵 pxxc 它有两倍多的列 [参数:pxx]. 奇数列包含置信区间的下限。 偶数列包含置信区间的上限。

+ 因此, pxxc[m,2*n − 1] —这是置信区间的下限，并且 pxxc[m,2*n] -估计值对应的置信区间的上限 <参数:pxx>[m,n].

*复盖置信区间的概率由自变量的值决定 [参数：ConfidenceLevel].

数据类型

漂浮物32, 漂浮64, Int64</无翻译>

例子：

使用默认输入数据进行Welch估计

Details

我们得到输入信号的Welch SPM估计值，该信号是具有角频率的离散正弦曲线 rad/计数和加性白噪声 .

让我们创建一个角频率的正弦曲线 rad/计数和加性白噪声 . 重置随机数发生器的状态以获得可再现的结果。 信号的长度为 数数。

using Random
Random.seed!(0)

n = 0:319
x = cos.(π/4 * n) + randn(length(n))

我们使用Hamming窗口和默认DFT长度获得Welch SPM的估计值。 默认段长度为 71 的计数，并且dft的长度包含 256 分，这给出了频率分辨率 rad/倒计时。 由于信号是真实的，SPM估计是片面的，并且有 256 / 2 + 1 分。

import EngeeDSP.Functions: pwelch

pxx, f, pxxc = pwelch(x, out = :data)
([0.1863630724000894; 0.33918260379618265; … ; 0.24304012258269603; 0.11708917663794487;;], [0.0; 0.02454369260617026; … ; 3.117048960983623; 3.141592653589793;;], [0.08146885100079476 0.7719778085293696; 0.1818050422612173 0.8436289786035511; … ; 0.130271774739574 0.6044994291550672; 0.05118568052389489 0.48502229985464423])

让我们建立一个Welch的SPM评估图。

pwelch(x)

pwelch 1

让我们重复计算。

*将信号分成长度的部分 . 这个动作相当于将信号分成尽可能长的段，以便得到尽可能多的，但不超过8个重叠的段。 50%. *我们使用汉明窗口段。 *指定重叠 50 相邻段之间的％。 *要计算FFT，我们使用 点在哪里 .

让我们确保这两种方法产生相同的结果。

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

import EngeeDSP.Functions: hamming

pxxt = pwelch(x, hamming(nsc), nov, nff, out = :data)

maxerr = maximum(abs.(abs.((pxxt[1])[:]) - abs.(pxx[:])))
0.0

将信号分成8个相等的重叠段 50 之间的％。 我们将指定与上一步相同的FFT长度。 让我们计算Welch的SPM分数，并确保结果与之前的估计相匹配。

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

pxxt8 = pwelch(x, lsc, floor(ov * lsc), nff, out = :data)

maxerr8 = maximum(abs.(abs.((pxxt8[1])[:]) - abs.(pxx[:])))
0.0

给定段长度的韦尔奇分数

Details

我们得到输入信号的Welch SPM估计值，该信号是具有角频率的离散正弦曲线 rad/计数和加性白噪声 .

让我们创建一个角频率的正弦曲线 rad/计数和加性白噪声 . 重置随机数发生器的状态以获得可再现的结果。 信号包含 512 数数。

using Random
Random.seed!(0)

n = 0:511
x = cos.(π/3 * n) + randn(length(n))

我们通过将信号分成长度段来获得Welch的SPM估计值 132 倒计时。 信号段乘以长度为Hamming窗口 132 倒计时。 未指定重叠计数的数量，因此将其设置为 132 / 2 = 66. DFT的长度为 256 分，这给出了频率分辨率 2π / 256 rad/倒计时。 由于信号是真实的，SPM估计是片面的，并且有 256 / 2 + 1 = 129 分。 让我们将SPM绘制为归一化频率的函数。

import EngeeDSP.Functions: pwelch

segmentLength = 132
pxx, w, pxxc = pwelch(x, segmentLength, out = :data)

plot(w/π, 10 * log10.(pxx),
     xlabel = "Normalized Frequency (× π rad/sample)",
     title = "Welch Power Spectral Density Estimate")

pwelch 2

具有给定分段重叠的Welch分数

Details

我们得到输入信号的Welch SPM估计值，该信号是具有角频率的离散正弦曲线 rad/计数和加性白噪声 .

让我们创建一个角频率的正弦曲线 rad/计数和加性白噪声 . 重置随机数发生器的状态以获得可再现的结果。 信号的长度为 320 数数。

using Random
Random.seed!(0)

n = 0:319
x = cos.(π/4 * n) + randn(length(n))

我们通过将信号分成长度段来获得Welch的SPM估计值 100 数数。 信号段乘以长度为Hamming窗口 100 数数。 重叠计数的数目为 25. DFT的长度为 256 分，这给出了频率分辨率 2π / 256 rad/倒计时。 由于信号是真实的，SPM估计是片面的，并且有 256 / 2 + 1 分。

import EngeeDSP.Functions: pwelch

segmentLength = 100
noverlap = 25
pxx, f, pxxc = pwelch(x, segmentLength, noverlap, out = :data)

plot(10 * log10.(pxx),
     xlabel = "Frequency Samples",
     title = "Welch Power Spectral Density Estimate")

pwelch 3

具有给定DFT长度的Welch估计

Details

我们得到输入信号的Welch SPM估计值，该信号是具有角频率的离散正弦曲线 rad/计数和加性白噪声 .

让我们创建一个角频率的正弦曲线 rad/计数和加性白噪声 . 重置随机数发生器的状态以获得可再现的结果。 信号的长度为 320 数数。

using Random
Random.seed!(0)

n=0：319
x=cos。(π/4*n)+randn(长度(n))

我们通过将信号分成长度段来获得Welch的SPM估计值 100 默认情况下具有重叠的计数 50%. 我们指定DFT的长度等于 640 分，使频率 π/4 Rad/计数对应于DFT的第81个bin。 由于信号是真实的，SPM估计是片面的，并且有 640 / 2 + 1 分。

import EngeeDSP.Functions: pwelch

segmentLength = 100
nfft = 640
pxx, f, pxxc = pwelch(x, segmentLength, [], nfft, out = :data)

plot(10 * log10.(pxx),
     xlabel = "Frequency Points",
     title = "Welch Power Spectral Density Estimate")

pwelch 4

使用频率为Hz的Welch方法估计SPM信号

Details

让我们创建一个由频率为 100 加性白噪声中的Hz . 重置随机数发生器的状态以获得可再现的结果。 采样率为 1 kHz，信号的持续时间 — 5 几秒钟。

using Random
Random.seed!(0)

Fs = 1000
t = 0:1/Fs:5-1/Fs
x = cos.(2π * 100 * t) + randn(length(t))

我们获得了前一个信号在重叠段上平均的Welch SPM的估计值。 使用段长度 500 重叠计数 300 数数。 我们使用 500 DFT点，使频率 100 Hz直接落入DFT仓中。 让我们输入采样频率以获得以Hz为单位的频率矢量。 让我们绘制结果图。

import EngeeDSP.Functions: pwelch

pxx, f, pxxc = pwelch(x, 500, 300, 500, Fs, out = :data)

plot(f, 10 * log10.(pxx),
     xlabel = "Frequency (Hz)",
     ylabel = "PSD (dB/Hz)")

pwelch 5

置信区间的上限和下限

Details

此示例说明了将置信边界与平均重叠Welch分段的SPM估计方法一起使用。 虽然这不是统计显着性的必要条件，但Welch估计中的频率，其中置信下限超过SPM周围估计的置信上限，清楚地表明时间序列中的显着波动。

让我们创建一个由频率正弦的叠加组成的信号 100 Hz和 150 加性白噪声中的Hz . 两个正弦波的幅值为 1. 采样率为 1 千赫。 重置随机数发生器的状态以获得可再现的结果。

using Random
Random.seed!(0)

Fs = 1000
t = 0:1/Fs:1-1/Fs
x = cos.(2π * 100 * t) + sin.(2π * 150 * t) + randn(length(t))

我们将得到一个具有置信水平的估计 95%. 让我们将段长设置为 200 计数和重叠 — 50％（即 100 计数）。

import EngeeDSP.Functions: hamming, pwelch

L = 200
win = hamming(L)
nOverlap = 100

pxx, f, pxxc = pwelch(x, win, nOverlap, 200, Fs, ConfidenceLevel = 0.95, out = :data)

让我们将SPM估计与置信区间一起绘制并放大感兴趣的频域。 100 和 150 赫兹。 置信区间的下限在紧邻 100 和 150 Hz明显高于该区域外置信区间的上限。

plot(f, 10 * log10.(pxx),
     xlim = (25, 250),
     xlabel = "Frequency (Hz)",
     ylabel = "PSD (dB/Hz)",
     title = "Welch Estimate with 95%-Confidence Bounds")

plot!(f, 10 * log10.(pxxc),
      linestyle = :dashdot,
      color = "#DD5400")

pwelch 6

以零频率为中心的功率谱

Details

让我们创建一个由频率为 100 加性白噪声中的Hz . 重置随机数发生器的状态以获得可再现的结果。 采样率为 1 kHz，信号的持续时间 — 5 几秒钟。

using Random
Random.seed!(0)

Fs = 1000
t = 0:1/Fs:5-1/Fs

noisevar = 1/4
x = cos.(2π * 100 * t) + sqrt(noisevar) * randn(length(t))

我们使用Welch方法使用长度段获得中心功率谱 500 计数，重叠 300 计数和DFT长度 500 分。

import EngeeDSP.Functions: pwelch

pxx，f，pxxc=pwelch（x，500，300，500，Fs，"居中"，"功率"，out=：数据）

让我们绘制结果图。 频率下的功率 −100 和 100 Hz接近预期功率 1/4 对于一个真正的正弦与振幅 1. 偏离 1/4 这是由于加性噪声的影响。

plot(f, 10 * log10.(pxx),
     xlabel = "Frequency (Hz)",
     ylabel = "Magnitude (dB)")

pwelch 7

评估多通道信号的SPM

Details

生成 1024 由三个频率正弦曲线组成的多通道信号的采样 , 和 加性白高斯噪声中的rad/计数 . 让我们评估Welch的SPM并构建一个图表。

N = 1024
n = 0:N-1

w = π ./ [2, 3, 4]

using Random
Random.seed!(0)

x = cos.(w * collect(n)')' + randn(length(n), 3)

import EngeeDSP.Functions: pwelch

pwelch(x)

pwelch 8

此外

通过平均重叠焊缝段进行SPM估计

Details

周期图不是对最广泛意义上的静止过程的真实谱功率密度的一致估计。 Welch方法通过将输入信号分割成重叠段来减小周期图方差。

Welch方法为每个段计算修改的周期图，然后对这些估计进行平均以获得频谱功率密度的估计。 由于该过程在最广泛的意义上是静止的，并且Welch方法使用输入信号的各个段的功率谱密度的估计，因此修改的周期图是真实功率谱密度的近似不相关

通常，分段乘以窗口函数，例如汉明窗口，使得韦尔奇方法减少到平均修改的周期图。 因为段通常重叠，段的开始和结束处的数据值，由一个段中的窗口平滑，远离相邻段的末端。 这可以防止因使用窗口模式而导致的信息丢失。

文学作品

  1. Hayes，Monson H. _Statistical Digital Signal Processing and Modeling._纽约：John Wiley＆Sons，1996。

  2. 斯托伊卡，佩特里和兰道夫*摩西。 _信号的光谱分析。_上鞍河,NJ:普伦蒂斯霍尔,2005.