# wsstridge
小波同步压缩的时频脊
函数库: TyWavelet
# 语法
fridge, = wsstridge(sst)
fridge, iridge = wsstridge(sst)
___ = wsstridge(sst, penalty)
___ = wsstridge(___, f)
___ = wsstridge(___, Name = Value)
# 说明
fridge, = wsstridge(sst) 从小波同步压缩变换(sst)中提取每个样本周期的最大能量时频脊。sst 输入是 wsst 的输出。每个小波脊都是一个单独的信号模式。
fridge, iridge = wsstridge(sst) 返回 iridge 中 sst 的行索引。行索引是每个样本处的最大时频脊。使用 iridge 使用 iwsst 沿着时间-频率脊重建信号模式。示例
___ = wsstridge(sst, penalty) 将频率仓之间的平方距离乘以惩罚值。您可以包含以前语法中的任何输出参数。
___ = wsstridge(___, f) 基于 f 输入频率向量,以每单位时间的周期为单位返回最大能量时频脊。f 是 wsst 的频率输出。f 输入和 fridge 输出具有相同的单位。
___ = wsstridge(___, Name = Value) 返回具有由一个或多个 Name = Value 对参数指定的附加选项的时间-频率脊。示例
# 示例
从 Chirp 信号中提取时频脊
获得二次 Chirp 的小波同步压缩变换,并提取 fridge 中的最大时频脊和 iridge 中的相关行索引。
加载 Chirp 信号并获得其同步压缩变换。
using TyWavelet
using TyBase
using TyPlot
pkg_dir = pkgdir(TyWavelet)
source_path = pkg_dir * "/examples/Resources/quadchirp.mat"
load(source_path)
sst, f = wsst(quadchirp)
提取最大时频脊。
fridge, iridge = wsstridge(sst)
绘制同步压缩变换。
h = pcolor(tquad, f, abs.(sst))
h.set_edgecolor("flat")
title("Synchrosqueezed Transform")
叠加最大能量频率脊的图。
hold("on")
plot(tquad, fridge; linewidth=2)
title("Synchrosqueezed Transform with Overlaid Ridge")

从多分量信号中提取时频脊
从多分量信号中提取两种最高能量模式。
获得并绘制小波同步压缩变换。
using TyWavelet
using TyBase
using TyPlot
pkg_dir = pkgdir(TyWavelet)
source_path = pkg_dir * "/examples/Resources/multicompsig.mat"
load(source_path)
sig = sig1 + sig2
sst, F = wsst(sig, sampfreq)
contour(vec(t), F, abs.(sst))
xlabel("Time")
ylabel("Hz")
grid("on")
title("Synchrosqueezed Transform of Two-Component Signal")
使用 10 的惩罚值,提取两个最高能量模式并绘制结果图。
fridge, iridge = wsstridge(sst, 10, F; NumRidges=2)
hold("on")
plot(t, fridge, "k"; linewidth=2)

# 输入参数
sst - 同步压缩变换矩阵
同步压缩变换,指定为矩阵。sst 是一个时频矩阵,是 wsst 的输出。
penalty - 频段缩放惩罚0(默认值) | 非负标量
频段缩放惩罚,指定为非负标量。此输入通过将惩罚值乘以频带间距离的平方来惩罚频率的变化。在提取多条脊线或在加性噪声中提取单一调制成分时,请使用惩罚项。惩罚项可以防止时频平面上最高能量区域突然变化时出现的频率跳跃。
f - 同步压缩变换频率向量
同步压缩变换的行对应的同步压缩变换频率,它是 wsst 的向量输出。频率向量中的元素数量等于 sst 输入中的行数量。
# Name - Value 参数
以 Name1=Value1,...,NameN=ValueN 的形式指定可选的参数对,其中 Name 是参数名,Value 是相应的值。名称-值参数必须出现在其他参数之后,但参数对的顺序并不重要。
示例: NumRidges = 3。
NumRidges - 最高能量时频脊的数量1(默认) | 正整数
要提取的最高能量时间频率脊的数量,由 "NumRidges "和一个正整数组成。如果该整数大于 1,wsstridge 会通过移除先前计算的脊和每个脊仓两侧的默认或指定 "NumFrequencyBins",迭代确定最大能量时间频率脊。
NumFrequencyBins - 要移除的频带数4(默认值) | 正整数
提取多个脊时,要从同步压缩变换 sst 中删除的频率仓数,指定为由 "NumFrequencyBins" 和一个正整数组成。此整数必须小于或等于四舍五入 round(size(sst,1)/4)。只有提取多个小波脊时,才能指定要删除的频率仓数。在提取最高能量的时频脊之后,wsstridge 在每个时间步长去除与 iridge 指数相对应的 sst 值。能量沿着在 iridge 索引两侧延伸指定数量的频率仓的时频脊被去除。如果扩展的时间-频率脊的索引在任何时间步长超过频率仓的数量,则 wsstridge 截断第一个或最后一个频率仓处的去除区域。若要指定 "NumFrequencyBins",必须指定 "NumRidges"。
# 输出参数
fridge - 时频脊频率矩阵
时间频率脊频率,以矩阵形式返回。fridge 是一个 N × nr 矩阵,其中 N 是 sst 中时间样本(列)的数量,nr 是脊的数量。矩阵的第一列包含 sst 中最大能量时频脊的频率,随后各列按能量递减顺序包含时频脊的频率。默认情况下,fridge 包含以每个采样周期为单位的频率。
iridge - 时间频率脊指数矩阵
sst 的时频脊行列指数,以矩阵形式返回。iridge 是一个 N × nr 矩阵,其中 N 是 sst 中时间样本(列)的数量,nr 是脊的数量。矩阵的第一列包含 sst 中最大能量时频脊的索引,随后各列按能量递减顺序包含时频脊的索引。
# 算法
该函数使用一种受惩罚的前向后向贪婪算法,从时间频率矩阵中提取最大能量脊。该算法通过最小化每个时间点的
下面的示例说明了时频脊算法使用的惩罚是频带间距离的两倍。具体来说,元素 (j,k) 和 (m,n) 之间的距离定义为 (j-m)
假设你有一个矩阵:
1 4 4 2 2 2 5 5 4更新 (1,2) 元素的值如下。
a. 保持第一个时间点的数值不变。从矩阵中的 (1,2) 元素开始运算,该元素表示第二个时间点的第一个频带,频带值为 4。根据第一列数值与 (1,2) 元素的距离对其进行惩罚。对第一列进行惩罚的结果是
original value + penalty × distance 1 + 2 × 0 = 1 2 + 2 × 1 = 4 5 + 2 × 4 = 13 1 4 4 2 13 5第一列的最小值为 1,位于第 1 仓。
b. 将第 1 列中的最小值加到当前分隔值 4 中,(1,2) 的更新值就变成了 5,即来自分隔 1 的值。
更新第 2 栏其余元素的值如下。
使用与步骤 2a 相同的方法,利用惩罚因子重新计算第 1 栏的原始值。使用与步骤 2b 相同的方法获取剩余的第二列值。例如,当更新(2,2)元素时,该元素的二进制值为 2,对该列应用惩罚因子后得到的结果是
original value + penalty × distance 1 + 2 × 1 = 3 2 + 2 × 0 = 2 5 + 2 × 1 = 7将最小值 2 添加到当前二进制数值中。(2,2) 的更新值变为 4。更新 (3,2) 元素后,矩阵为
1 5(1) 4 2 4(2) 2 5 9(2) 4只更新了第二列。下标表示数值来自前一列中的分仓索引。
对第三列重复步骤 2。但现在,惩罚将应用于更新后的第二列。例如,更新 (1,3) 元素时,惩罚为
5 + 2 × 0 = 5 4 + 2 × 1 = 6 9 + 2 × 4 = 17第一分区中的最小值 5 被添加到(1,3)分区值中。更新第三列中的所有值后,最终矩阵为
1 5(1) 9(1) 2 4(2) 6(2) 5 9(2) 10(2)从矩阵的最后一列开始,找出最小值。在矩阵中回溯时间,从当前分隔点回到上一个时间点该分隔点的原点。跟踪构成脊的路径的分区索引。该算法通过使用原点仓而不是最小值仓来平滑过渡。在本例中,脊指数为 2、2、2,与步骤 1 中矩阵第 2 行的正弦波能量路径相匹配。
如果要提取多条脊线,算法会从时频矩阵中移除第一条脊线,然后重复这一过程。
# 参考文献
[1] Daubechies, I., J. Lu, and H.-T. Wu. "Synchrosqueezed wavelet transforms: an empirical mode decomposition-like tool." Applied and Computational Harmonic Analysis. Vol. 30, Number 2, 2011, pp. 243–261.
[2] Thakur, G., E. Brevdo, N. S. Fučkar, and H.-T. Wu. "The Synchrosqueezing algorithm for time-varying spectral analysis: Robustness properties and new paleoclimate applications." Signal Processing. Vol. 93, Number 4, 2013, pp. 1079–1094.