等）。离散化的小波族构成紧标架， 计算小波变换的系数和从小波系数重构信号变得很方便，同样也保证了后面求解优化问题 时算法的收敛性。
[0134] 小波标架的算子表示
[0135]令
[0136]Vm,n= 巧8)
[0137] 则公式（14)和（15)可W写成
[014引本发明定义K为从映射到王2 (化；)的算子，将标架系数C= (C,,w),,wcH映射为 王2 (化）上的信号，即 [0143]
(62)
[0144] 算子K的伴随算子r为从￡2 (吸)映射到f2 (妒)的算子，将王2 (吸)上的信号投影到 标架系数上：
[0145]K*s= <s,Vm，n〉，化3)
[0146] 则公式（19)可W写成
[0147] S=KK*s,化4)
[0148] 算子K为合成算子，r为分析算子，和小波标架{>' 相对应。
[0149] 有效信号能量分布空间的确定
[0150] 当选用合适的小波函数对含噪信号进行小波变换，将其投影到时间尺度域的时 候，有效信号的能量将会分布在较小的子空间V，噪声的能量会扩散到比较大的子空间 V'，甚至整个时间-尺度域。换句话说，有效信号和少数的系数相对应，而噪声几乎分布在 全部的系数上。如果确定了有效信号对应的子空间V，得到有效信号对应的系数，将其它系 数置零，则噪声就会在变换域被压制，信噪比得到提高。
[0151] 为了确定有效信号的能量分布，希望有效信号投影到尽可能少的小波标架原子 上，也就是说希望信号在变换域得到尽可能稀疏的表示，该样就能压制更多的噪声，提高信 噪比。具体来说，一个含噪信号可W表示为
[0152]y=s+n=Kx+n, (65)
[015引其中y表示含噪信号，s为不含噪的有效信号，n为高斯白噪声，K为公式（2。中 的合成算子，X表示变换域的系数。希望得到一个最优的系数X，使得（1)S--KS，S尽可 能稀疏；(2) ||s-4的值小。稀疏性和1。范数最小化问题相对应，是非凸问题，求解难度很 大，常常在一定条件下转化为li范数最小化问题。因此，可W通过求解下列优化问题得到 X:
[0154]
[0155] 其中e和待分析信号的噪声水平有关。根据Elad等人提出的方法，可W用 Lagrange乘子A将该问题转化为W下的无约束问题：
[0156]
[0157] 其中A被称为正则化参数，该样一来，可W直接采用迭代萎缩阔值算法（1ST)来 求解该问题。Daubechies等人指出，采用下列公式所示的迭代方法对系数X进行更新，当迭 代次数足够大的时候，迭代的结果收敛于公式（24)所示优化问题的解。
[0巧引 x(w)=T,[x(k)+K*(y-Kx(k))]，k= 1，2,…，N，化8)
[0159] 其中为阔值函数，定义为
[0160]
[0161] 需要指出的是，由于迭代萎缩阔值算法（1ST)中采用的是软阔值（如公式（27)所 示），此时阔值和正则化参数相等，用A表示。
[0162] 在求解优化问题过程中，需要采取阔值下降策略，即；随着迭代次数的增加，入的 值不断减小。常用的阔值下降策略有线性下降、指数下降等，Gao等人研究发现，采用指数 阔值下降方法能够大大提高算法收敛的速度。
[0163] 在每一次迭代过程中，采用指数阔值下降策略，如公式（28)所示：
[0164]
[0165] 其中
[0166]
[0167] 公式（29)中Amax和Amin分别为正则化参数的最小值和最大值，Pma济Pmin为最 大和最小百分比。
[016引上面的阔值策略只有在迭代次数达到最大值N的时候才会停止，为了在能够及时 中止迭代过程，本发明定义了一个动态停止准则
[0169]
[0170] 其中tolerance为迭代前给定的容许值。在该准则下，当继续迭代不会取得更好 的结果时，迭代过程就会自动中止，W减少不必要的迭代次数。
[0171] 迭代萎缩阔值算法（1ST)的收敛速度很慢，仅为0(l/k)，即具有一阶收敛速度。 Beck和Teboulle提出了快速迭代萎缩阔值算法（FIST)，见公式（31)，该方法使得算法收敛 速度提高为0(l/k2)，即具有二阶收敛速度。
[0172]
[017引其中x?= 0,tW= 1，tW满足W下递推公式：
[0174]
[01巧]表2给出了求取有效信号对应系数的流程（伪代码），通过迭代，本发明得到了优 化问题的解，即与噪声被压制后的信号相对应的系数X。图l(a)-(b)分别画出了不含噪和 信噪比为5地的50化的化cker子波的变换域谱图，可W看到噪声对有效信号的分布造成 的影响。图1(c)画出了经过迭代W后得到的系数构成的谱图，图中噪声被有效压制，准确 刻画了有效信号的能量分布。
[0176] 表2通过迭代求取有效信号对应系数的流程（伪代码）
[0177]
[0178] 图150化的化cker子波的时间尺度谱图，可W看出在含噪信号的小波变换谱图 中，有效信号能量分布受到噪声的影响，经过迭代W后，噪声得到压制，谱图清晰地反映出 有效信号的能量分布，有了系数S之后，本发明可W得到和有效信号对应的解析信号，进一 步计算瞬时属性。
[0179] 本发明的物质基础是地震数据体，采用的逐道处理办法。具体步骤为：
[0180] 步骤1 ;根据公式（17)计算各道地震数据对应的广义Morse标架系数；
[0181] 步骤2;运用表2的流程进行迭代，得到有效信号对应的系数；
[0182] 步骤3 ;由公式似计算解析信号；
[018引步骤4;根据公式（3)-巧）计算瞬时属性。
[0184] 效果分析
[0185] 合成信号算例
[0186] 本发明用=个测试信号来检验不同方法计算瞬时频率的效果。图2(a)-(c)是分 别是50化的化cker子波，用50化的化cker子波与反射系数序列卷积合成的地震记录 和实际地震道数据。首先，本发明用传统的化化6的方法计算=个信号的瞬时频率，如图 2(d)-(f)所示，然后，采用本发明提出的方法，选用0=1，丫 =3的广义Morse小波标架 来计算瞬时频率，如图2(g)-(i)所示，可见两种方法得到的结果几乎相同。
[0187] 图2为测试信号和不同方法计算的瞬时频率，可见在不含噪情况下，两种方法均 能得到准确的瞬时频率计算结果。
[018引然而，当信号被噪声干扰时，HUbed变换法计算的瞬时频率准确性大大降低。本 发明给S个测试信号加上噪声，其信噪比均为5地（图3 (a)-(c))，然后重新用两种方法计 算瞬时属性，如图3(d)-(f)和图3(g)-(i)所不。
[0189] 图3含噪测试信号和不同方法计算的瞬时频率，可见在HUbed变换法的计算结 果中，有效信号的瞬时频率完全被噪声掩盖，而采用本发明方法的计算结果，即使在强噪声 条件下，仍能够清晰地显示出有效信号的瞬时频率
[0190] 在本发明提出的方法中，通过求解优化问题，得到和有效信号对应的系数，然后计 算解析信号和瞬时频率。可见，即使在噪声比较严重的情况下，该方法仍能够准确得到测试 信号的瞬时频率。而传统的HiAed变换方法得到的结果，有效信号的瞬时频率完全被噪 声影响，无法辨识。
[0191] 下面，将本发明提出的方法和基于HUbed变换方法计算的含噪测试信号的瞬时 频率和不含噪信号的瞬时频率画在一起，并且在本发明提出的方法中选用不同的小波标架 进行计算，结果如图4所示。
[0192]图4不同方法计算出来的含噪测试信号的瞬时频率，其中藍线为含噪测试信号的 瞬时频率，红线为不含噪信号的瞬时频率。分别采用HUbed方法和本发明提出的方法计 算，在本发明提出的方法中，分别选用广义Morse标架（GMW)和Morlet小波标架。可W看 出本发明提出的方法具有更好的抗噪性能和精度，同时，采用广义Morse标架的结果