连续小波变换，Java 实现
4.31/5 (5投票s)
本文介绍了一个执行连续小波变换的 Java 示例应用程序。
sha256sum CwtSrc.zip: 1e3d72421b8d345a0b2f6f90cc209f8e6a8aa13f17c5bb6bc1a6fdbf5721729f
sha256sum CwtRunnable.zip: c3d872322f4b7803c22c7a1264dcd4493ae75eecc197c44e470d1d39e14223b6
引言
这是我之前一篇关于离散小波变换的文章的配套文章。在本文中，我提供了一个使用连续小波变换来探索一维信号的应用程序。该方法可以进一步开发，以转换具有更高维度（如图像）的信号。还有其他方法可以将小波功能添加到 Java 中，例如使用 Matlab/Scilab-Wavelab 包装器或开源库。为了说明，我直接在 Java 方法中展示了计算的基础知识。
背景
在傅里叶信号分析中，任何信号都可以理解为叠加的复数波形分量的线性组合，每个分量都是静止的。信号的离散傅里叶变换揭示了分量的频率，这些分量在分辨率允许的情况下可以相加重构原始信号。但是，由于信号的分量是频率的平稳函数，因此确定的所有频率都被解释为在信号的整个持续时间内都存在。如果信号具有需要仔细检查的局部结构，例如气候历史中的温度趋势、地震记录中的事件、经济中的崩盘等，该怎么办？
一种解决方案是迭代地对信号中的子采样时间窗口执行一系列短期傅里叶变换。另一种方法是检查信号与一个可以迭代缩放和平移的短时函数——小波——的卷积。尺度为s的 CWT 是尺度为s的平移小波与信号的卷积。多个尺度的 CWT 一起构成多分辨率分析 (MRA)。
小波的一个例子：高斯分布的第 6 导数
我的第一篇关于离散小波变换 (DWT) 的文章演示了离散小波变换的过程，该过程通过连续减半的小波尺度来迭代执行变换。连续小波变换 (CWT) 允许相邻尺度之间进行更精细的分离。术语连续具有误导性，因为分析是在离散采样的信号上，在离散分离的尺度上进行的。对于某些应用，CWT 比 DWT 更能揭示结构。例如，DWT 文章中的一个兴趣点是揭示的坎托阶梯函数的自相似性。这可以通过 CWT 更好地演示。
由于其相对的初创性和适应机会，小波分析实现受到非标准化、有时冲突的特定应用处理方式的影响。小波本身有多个同义词。变换的参数以不同的方式表达。特征量的术语各不相同。一些处理不归一化尺度，而另一些使用不同的（有时是任意的）归一化方法。数据呈现方式差异巨大。要理解本文中的实现如何适应这种分散性，以下参考文献非常重要且易于获取
我的 CWT 实现的开发最初基于对他们出版物的研究以及对他们wavelet.f Fortran 代码中提供的 Fortran 函数的算法分析。上述作者的贡献超出了我文章的范围，特别是对于他们演示的地球物理应用中对结果进行有意义解释所需的统计分析而言。
本文演示的目的是为有兴趣使用 CWT 强大功能的程序员提供一个切入点。它提供了一个环境，用于探索几个函数的多分辨率分析，其中一些函数是解析的且具有可验证的启示。
尝试可运行的 jar 文件
要运行演示，请获取上面的 CwtRunnable.zip 下载，解压缩，找到 CwtRunnable.jar 文件，然后运行它。在 Windows 上，我只需双击 CwtRunnable.jar 文件即可运行它。在 Linux 上，我在终端中键入 java -jar <thePath>/CwtRunnable.jar。显然，您需要安装 Java，最好是最新版本并包含所有更新。
这是直观的初始窗口。首先，尝试默认设置，然后进行实验。
主入口点：绘图生成器
点击生成表格将在一个窗口中列出函数。从这里可以绘制、变换、保存等。
函数的表格窗口及其图
点击图将显示点击点的坐标在文本框中。
表格中的 CWT 菜单将打开一个窗口，用于输入变换小波和相关参数。当您点击“确定”时，将创建一个多分辨率表格。
CWT 参数选择器和 MRA 表格
现在我们可以以图形方式检查数据。首先，选择同步绘图菜单项。顾名思义，它在同一图表上同时绘制多分辨率分析的分量。这是 Chirp 示例的同步绘图
Chirp 示例的同步 MRA 绘图
绘图菜单呈现一个选择器表，您可以使用它来选择性地隐藏/显示分量。使用 [Ctrl] 上/下键缩放 x 轴，使用 [Alt] 上/下键缩放 y 轴。同样，点击图表会显示坐标。
MRA 表格菜单中的尺度图菜单选项提供了文献中最常见的表示方法。
解压的演示中包含两个逗号分隔的数据文件。从绘图生成器窗口，点击从 csv 打开 xy，然后导航到下载文件夹中的TestSignal.csv文件。从数据表格窗口绘制函数。
该信号使用以下分量以线性组合的形式合成
- 从 t = 0 到 0.255875，dt = 0.000125
- 分量 1：sin(2pi500t)
- 分量 2：sin(2pi1000t)
- 分量 3：0.1sin(2pi8t)
- 分量 4：一个幅度为 5 的脉冲，位于 t = 0.15 和 0.154
获取 CWT 并从 MRA 表格菜单中选择尺度图 > 模值。信号的分量规格告诉我们，在 f = 8、f = 500 和 f = 1000 附近应该有主要的暗带。在 0.150 和 0.154 处也应该有两个脉冲。
当您尝试这些示例时，某些尺度图需要进行对比度调整，该功能由鼠标滚轮控制。
从 MRA 表格窗口的菜单中选择重构可以现实地检查您选择的小波和参数选择是否适合有效分析您的信号。为了最大程度地重构，请在参数输入窗口中使用较小的值作为最小尺度和较大的值作为总尺度数。较低的相对残差表示良好的重构。尽管减小最小尺度并增加变换中的尺度数量可以提高重构质量，但它不一定能提高 MRA 数据本身的可用性。在这里，我们必须考虑影响锥。
CWT 实现涉及使用卷积定理迭代地将缩放的小波与信号进行卷积。卷积是两个分量的傅里叶变换的点积的傅里叶逆变换。信号的有限、离散采样在序列的开始和结束处存在突然的边界，并且与边界重叠的卷积会受到误差的污染。这种影响可以通过一个在信号开始和结束处具有脉冲而在其他地方为零的函数来可视化。这是该函数 MRA 尺度图。
灰色区域按比例显示了 MRA 中受误差影响的区域（在影响锥内）。尺度越大的变换（越往上）越会落入灰色区域。有各种方法可以图形化地补偿误差（例如，余弦阻尼），但它们需要对信号进行人工失真/外推。本文开头处的“魔鬼阶梯函数”的尺度图故意进行了余弦阻尼处理以生成图形。如果可能，可以收集更长时间的信号，并选择性地忽略 MRA 中早期和晚期的特征。还有其他方法可以管理误差，例如对称扩展，但适用性因应用而异。
与合成函数不同，自然系统产生的信号是嘈杂的。必须以预定义的置信度统计确定 MRA 观测值的重要性。提出的问题是：我们是否可以合理地确信一个观测到的特征实际上不是噪声引起的偏移？回答这个问题需要理解噪声在所研究系统中如何分布以及相关统计知识。这是特定于应用程序的，并且由于此演示侧重于精确的合成信号，因此不涉及显著性检验。
使用代码
每个类的头文件中都包含一个许可证：Copyright Mark Bishop, GNU v3。其他贡献已注明。作者不对所提供信息的准确性、完整性、安全性或有用性作任何保证，也不声明其使用不会侵犯私人所有权。
CwtSrc.zip 下载是一个 zip 格式的 Eclipse 项目文件夹。如果您使用 Eclipse，请将其解压缩到您选择的工作区，然后在 IDE 中创建一个名为 CwtDemo 的项目。这将填充项目并启用测试运行。使用其他 IDE 或选择其他方式创建 Eclipse 项目时，您需要在构建路径中链接到 jel.jar 库才能使用示例函数表达式。它位于 jel-2.0.1 文件夹中。我提交了一篇关于使用 JEL 表达式库的文章，网址为www.codeproject.com/Articles/737497/Java-Function-Compiler。JEL 表达式是演示的一部分，但不是使用 CWT 特定类所必需的。
计算类是获得 CWT 多分辨率分析所必需的唯一类。其余类和 JEL 表达式库是演示的辅助功能，提供了图形界面、表格、对话框、绘图等。
项目中包含以下类
- 计算类
- ComplexCalc.java (复数计算的静态方法)
- ComplexNumber.java (包含单个复数的对象)
- CWT.java (本文的主要主题)
- Fft.java (正向和反向傅里叶变换方法)
- Gamma.java (计算伽马函数)
- Hermite.java (计算埃尔米特多项式的方法)
- MatrixOps.java (此处用于数组管理)
- 交互式 GUI 类
- MainForm.java (主绘图生成器窗口)
- MraForm.java (MRA 表格窗口)
- XYForm.java (函数表格窗口，t 和 f(t))
- PlotFrame.java (函数和 MRA 图)
- PlotVisiblePicker.java (同步绘图的可见性选择器)
- ParameterInputCwt.java (用于输入变换参数的对话框窗口)
- Scalogram.java (尺度图表示)
- 其他实用程序
- Examples.java (几个有趣的示例函数)
- FileOps.java
- ImagePanel.java
- Signal.java (用于管理 xy 数据的助手)
- StringUtils.java
- ConsoleTest.java (每个小波的有限重构测试例程)
- VariableProvider.java (jel 所需)
最重要的一类是 CWT.java。在此类中，一个方法负责返回多分辨率数据：cWT()。
/**
* Continuous wavelet transform for a real, discretely sampled function y =
* f(t) of length n where n is an integer power of 2
*
* @param y
* The time-sampled signal to transform.
* @param dt
* The sampling increment
* @param mother
* The wavelet to use in the transform (Morlet, DOG, Paul)
* @param param
* wavelet parameter (Morlet param = wave number, Paul param =
* order, DOG param = m::m-th derivative of Gaussian)
* @param s0
* The smallest desired scale in the transform (try s0=dt for
* Morlet; s0=dt/4 for Paul)
* @param dj
* The desired spacing between successive scales (try 0.25)
* @param jtot
* The desired number of scales to be computed
* @return An ArrayList[Object] wspcmf[6] where: xspcmf[0] is a
* Complex[n][jtot] multi-resolution matrix, wspcmf[1] is a
* double[jtot] vector of scales, wspcmf[2] is a double[jtot] vector
* of periods, wspcmf[3] is a double[n] cone-of-influence vector,
* wspcmf[4] is the signal mean, and wspcmf[5] is the Fourier
* factor.
**/
public static ArrayList<Object> cWT(double[] y, double dt, Wavelet mother, double param, double s0, double dj, int jtot)
- 输入参数
- y[]: 用于变换的信号
- mother: 小波（Morlet、高斯导数 (DOG) 或 Paul）
- param: 小波参数（Morlet：波数，Paul：阶数，DOG：m，其中我们使用高斯分布的第 m 阶导数）
- s0: MRA 中要包含的最小尺度
- dj: 尺度之间的间隔
- jtot: 请求的总尺度数
- 输出：一个名为 wspcmf 的 ArrayList<Object>，包含以下内容
- xspcmf[0] 是一个 ComplexNumber[n][jtot] 多分辨率矩阵
- wspcmf[1] 是一个 double[jtot] 尺度的向量
- wspcmf[2] 是一个 double[jtot] 周期的向量
- wspcmf[3] 是一个 double[n] 影响锥向量（演示中未使用）
- wspcmf[4] 是信号的平均值
- wspcmf[5] 是傅里叶因子
该方法首先创建一个名为 yNoDC 的 y 的副本，并从中删除了信号的平均值。接下来，其傅里叶变换作为 yfft。创建了波数向量后，主循环调用 getDaughter() 方法，该方法在迭代尺度上返回相应小波函数的傅里叶变换（参见 Torrence & Compo，表 1）。然后，通过获得 daughter 和 yfft 的乘积的傅里叶逆变换，循环在尺度上执行卷积。
每个尺度的 CWT 可以与尺度集、周期或周期的倒数（频率）一起制成表格或共绘。在演示中使用了倒数周期。傅里叶因子和周期向量包含在返回对象中以便于使用，也可以根据需要在一个消耗方法中计算。可以处理影响锥，以设想信号边界效应污染 MRA 的区域。
CWT.java 中的 cwtReconstruct() 方法从 MRA 矩阵的实部和尺度向量重构原始信号。它需要信号的平均值（在变换期间从信号中减去）以及采样增量 dt。其他输入参数必须与原始变换中使用的值匹配。cwtReconstruct() 方法调用 getEmpiricalFactors() 方法，该方法通过 delta 函数校准（参见 Torrence & Compo）确定小波的重构因子。请查看我的文章，题为选定小波的连续小波变换重构因子，网址为：mark-bishop.net/signals/CWTReconstructionFactors.pdf。显然，重构因子可以硬编码，而不是在每次执行重构时计算。
cwtReconstruct.java 类还包含一个方便的方法 getSelectedParamChoices()，它提供了演示中三种小波的一些可能的小波参数。DOG 5 的选择作为脑筋急片被包含在内。最后，该类有一个枚举 Wavelet，包含三种小波选择。
一些参考资料
Marie Farge, Wavelet Transforms and Their Applications to Turbulence, Annu. Rev. Fluid Mechanics, 24:395-457, 1992
S. D. Meyers, Kelly, Obrien; An Introduction to Wavelet Analysis in Oceanography and Meteorology: With Application to the Dispersion of Yanai Waves, Mon. Wea. Rev., 121, 2858-2866, 1993。
E. Oran Brigham, The Fast Fourier Transform, Prentice-Hall, 1974。
Paul Addison, The Illustrated Wavelet Transform Handbook, Taylor & Francis, 2002。
Jaideva Goswami, Andrew Chan, Fundamentals of Wavelets, Wiley, 1999。
http://users.rowan.edu/~polikar/WAVELETS/WTtutorial.html
历史
在此处保持您所做的任何更改或改进的实时更新。