# pmusic
使用MUSIC算法的伪谱
函数库: TySignalProcessing
# 语法
S, wo, = pmusic(x, p)
S, wo, = pmusic(___; nfft = Integer)
S, wo, = pmusic(___; corr = true)
S, fo, = pmusic(___; nfft = Integer, fs = Number)
S, fo, = pmusic(___; fi = Vector, fs = Number)
S, fo, = pmusic(___; nfft = Integer, fs = Number, nw = Integer, noverlap = Integer)
S, wo, = pmusic(___; freqrange = value)
S, fo, v, e = pmusic(___)
___ = pmusic(___; plotfig = Bool)
# 说明
S, wo, = pmusic(x, p) 实现了多信号分类 (MUSIC) 算法,并返回 S、输入信号 x 的伪频谱估计值以及评估伪频谱的归一化频率 (以rad/sample 为单位) 的矢量wo。可以使用输入参数 p 指定信号子空间维数。 示例
S, wo, = pmusic(___; nfft = Integer) 指定用于估计伪谱的 FFT 的整数长度 nfft。此语法可以包括其他语法输入参数的任何组合。
S, wo, = pmusic(___; corr = true) 使输入信号 x 被解释为相关矩阵而不是信号数据矩阵。对于这种语法,x 必须是一个方阵,并且它的所有特征值必须是非负的。 示例
S, fo, = pmusic(___; nfft = Integer, fs = Number) 以 Hz 为单位提供采样率 fs, 返回指定的频率下计算的伪谱。
S, fo, = pmusic(___; fi) = Vector, fs = Number) 返回指定的频率下计算的伪谱。向量 fi 必须有两个或多个元素,否则函数将其解释为 nfft。
S, fo, = pmusic(___; nfft = Integer, fs = Number, nw = Integer, noverlap = Integer) 通过使用窗口 nw 和重叠长度 noverlap 分割输入数据 x ,并返回伪谱 S。 示例
S, fo, v, e = pmusic(___) 返回噪声特征向量的矩阵 v 以及向量 e 中的相关特征值。
___ = pmusic(___; freqrange= value) 指定要包含在 fo 或 wo 中的频率值范围。 示例
___ = pmusic(___; plotfig = Bool) 在当前图形窗口中绘制伪谱。
# 示例
不指定采样率的 MUSIC 算法计算伪谱
此示例分析信号向量 x,假设信号子空间中存在两个实正弦分量。在这种情况下,信号子空间的维数是 4,因为每个实正弦曲线是两个复指数的和。
using TySignalProcessing
using TyMath
n = 0:199
rng = MT19937ar(1234)
x = cos.(0.257 .* pi .* n) .+ sin.(0.2 .* pi .* n) .+ 0.01 * randn(rng, length(n))
Pxx, w, = pmusic(x, 4; plotfig=true)

指定采样频率和子空间维度
本例分析了相同的信号向量 x,其特征值截止值比最小值高 10%。设置 p(1)=Inf 使信号/噪声子空间决策基于阈值参数 p(2)。使用 nw 参数指定长度为 7 的特征向量,并将采样频率 fs 设置为 8 kHz:
using TyPlot
using TySignalProcessing
using TyMath
rng = MT19937ar(1234)
n = [0:199;]
x = cos.(0.257 .* pi .* n) .+ sin.(0.2 .* pi * n) .+ 0.01 .* randn(rng, length(n))
P, f, = pmusic(x, [Inf, 1.1]; fs=8000, nw=7)
plot(f, 20 * log10.(abs.(P)))
xlabel("Frequency (Hz)")
ylabel("Power (dB)")
title("Pseudospectrum Estimate via MUSIC")
grid("on")

输入相关矩阵
提供一个正定相关矩阵 R,用于估计谱密度。使用默认的 256 个样本。
using TyBase
using TySignalProcessing
using TyMath
R = toeplitz(cos.(0.1 * pi * (0:6))) + 0.1 * eye(7)
pmusic(R, 4; corr=true, plotfig=true)

输入 corrmtx 生成的信号数据矩阵
输入使用 corrmtx 从数据生成的信号数据矩阵 Xm。
using TySignalProcessing
using TyMath
rng = MT19937ar(1234)
n = 0:699
x = cos.(0.257 * pi * (n)) + 0.1 * randn(rng, length(n))
Xm, = corrmtx(x, 7, "modified")
pmusic(Xm, 2; plotfig=true)

使用窗创建信号数据矩阵
使用相同的信号,但让 pmusic 使用其窗口输入参数形成 100 × 7 数据矩阵。此外,指定长度为 512 的 FFT。
using TySignalProcessing
using TyMath
rng = MT19937ar(1234)
n = 0:699
x = cos.(0.257 * pi * (n)) + 0.1 * randn(rng, length(n))
PP, ff, = pmusic(x, 2; plotfig=true, nfft=512, nw=7, noverlap=0, fs=1)

# 输入参数
x - 输入信号向量 | 矩阵
输入信号,指定为向量或矩阵。如果 x 是一个向量,那么它被视为信号的一个观测值。如果 x 是矩阵,则 x 的每一行表示信号的单独观察。例如,在数组处理中,每一行是传感器数组的一个输出,因此 x'*x 是相关矩阵的估计值。
提示
您可以使用 corrmtx 的输出生成 x。
数据类型: Float
是否支持复数: 是
p - 子空间维度实正整数 | 二元向量
子空间维度,指定为实正整数或双元素向量。如果 p 是实正整数,则将其视为子空间维数。如果 p 是一个双元素向量,则 p 的第二个元素表示一个阈值,该阈值乘以信号相关矩阵的最小估计特征值 λmin。将低于阈值
提示
如果 pmusic 的输入是实的正弦曲线,则将 p 的值设置为正弦曲线数的两倍。如果输入是复正弦曲线,则将 p 设置为等于正弦曲线的数量。
是否支持复数: 是
nfft - DFT 点数256(默认值)|整数
DFT 点数,指定为正整数。
fs - 采样率1(默认值)|正标量
采样率,指定为以 Hz 为单位的正标量。采样率默认为 1 Hz。
fi - 输入频率(Hz)向量
输入频率,指定为向量。在向量中指定的频率下计算伪谱。
nw - 矩形窗的长度2*p(1) (默认值) | 非负整数
矩形窗口的长度,指定为非负整数。
noverlap - 重叠样本数nw-1(默认值) | 非负整数
重叠样本数,指定为小于窗口长度的非负整数。
提示
在语法中包含 "corr" 时,将忽略参数 nw 和 noverlap。
freqrange - 伪谱估计的频率范围"half" | "whole" | "centered"
伪谱估计的频率范围,指定为 "half"、"whole" 或 "centered" 之一。
"half":返回实际输入信号 x 的一半频谱。如果 nfft 为 偶数,则 S 的长度为 nfft/2+1,并在区间 [0, π] 上计算。如果 nfft 为奇数,则 S 的长度为 (nfft+1)/2,频率间隔为 [0,π)。当指定 fs 时,偶数和奇数 nfft 的间隔分别为 [0, fs/2) 和 [0, fs/2];
"whole":返回实数或复数输入 x 的整个频谱。在这种情况下,S 的长度为 nfft,并在间隔 [0,2π) 上计算。如果指定 fs,则频率间隔为 [0,fs);
"centered":返回实数或复数输入 x 的居中全谱。在这种情况下,S 的长度为 nfft,对于偶数 nfft,在区间 (–π,π) 上计算,对于奇数 nfft,则在区间 (-π,π] 上计算。指定 fs 时,偶数和奇数 nfft 的频率间隔分别为 (–fs/2,fs/2] 和 (–fs/2,fs/2)。
提示
您可以将参数 freqrange 或 corr=true 放在输入参数列表中 p 之后的任何位置。
# 输出参数
S - 伪谱估计向量
伪谱估计,作为向量返回。使用与输入数据 x 相关联的相关矩阵的特征向量的估计来计算伪谱。
wo - 输出归一化频率向量
输出归一化频率,指定为向量。S 和 wo 具有相同的长度。通常,FFT 的长度和输入 x 的值决定了计算的 S 的长度和相应的归一化频率的范围。该表指示了 S 和 wo 的长度以及第一语法的对应归一化频率的范围。
FFT 长度为 256 的 S 特性 (默认值)
| 输入数据类型 | S 和 wo 的长度 | 相应归一化频率的范围 |
|---|---|---|
| 实数 | 129 | [0, π] |
| 复数 | 256 | [0, 2π) |
如果指定了 nfft,下表显示了 S 和 wo 的长度以及 wo 的频率范围。
S 和频率向量特性
| 输入数据类型 | nfft 偶数或奇数 | w 的范围 | 相应归一化频率的范围 |
|---|---|---|---|
| 实数 | 偶数 | (nfft/2)+1 | [0, π] |
| 实数 | 奇数 | (nfft+1)/2 | [0, π) |
| 复数 | 偶数或奇数 | nfft | [0, 2π) |
fo - 输出频率向量
输出频率,以向量形式返回。fo 的频率范围取决于 nfft、fs 和输入 x 的值。S 和 fo 的长度与上述 S 和频率向量特性相同。如果指定了 nfft 和 fs,下表指示了 fo 的频率范围。
指定 fs 的 S 和频率向量特性
| 输入数据类型 | nfft 偶数或奇数 | f 的范围 |
|---|---|---|
| 实数 | 偶数 | [0, fs/2] |
| 实数 | 奇数 | [0, fs/2) |
| 复数 | 偶数或奇数 | [0, fs) |
此外,如果还指定了 nw 和 noverlap,则在公式化用于估计相关矩阵特征值的矩阵之前,对输入数据 x 进行分段和加窗。数据的分割取决于 nw、noverlap 和 x 的形式。下表描述了对生成的窗口段的注释。
窗口数据取决于 x 和 nwn
| x 的形式 | nw 的形式 | 窗口化数据 |
|---|---|---|
| 向量 | 标量 | 长度为 nwn |
| 向量 | 系数的向量 | 长度为 length(nwn) |
| 矩阵 | 标量 | 数据未窗口化。 |
| 矩阵 | 系数的向量 | length(nwn) 必须与 x 的 列长度相同,并且不使用 noverlap。 |
v - 噪声特征向量矩阵
噪声特征向量,以矩阵形式返回。v 的列跨越维度大小 (v,2) 的噪声子空间。信号子空间的维数是 size(v,1)-size(v,2)。
e - 估计特征值向量
相关矩阵的估计特征值,作为向量返回。