基于小波变换的重力异常计算方法及装置与流程

1.本发明涉及大地测量技术领域，尤其涉及一种基于小波变换的重力异常计算方法及装置。


背景技术：

2.重力异常是地表、地壳、地幔及深部构造的密度不均匀分布的叠加产生的共同影响。随着现代大地测量技术的发展，对于油、气、矿藏等自然资源勘探，国家高程基准的构建所采用的基础资料——重力异常信号的获取手段，已经由传统的地面测量获取逐渐发展成为空间手段，比如航空重力测量等方式。空间手段获得的重力资料具有覆盖面广、测量效率高和节省人力资源开销的优点。但在实际应用过程中，需要采用一定的数据处理方法将重力资料向下延拓至地面，才方便使用。
3.重力异常向下延拓是经典的不适定问题，向下延拓的方法可分为空间域和频率域进行分类。在空间域，由于观测矩阵条件数大，求逆容易受到误差的影响，使得延拓计算结果不稳定。为此，为了改善法方程求解不稳定的问题，涌现了大量的针对法方程正则化等的多种改进算法。在频率域进行向下延拓时，同样存在对误差的放大效应，通常采用滤波的方法抑制误差。实际上在抑制误差的同时，也抑制了有用的信号。同时，现有方法常通过将观测值的多尺度边缘视为地面延拓值的多尺度边缘，并据此来判断虚假边缘，实现向下延拓。但是多尺度边缘约束只适用于延拓高度较低的情况，对卫星高度的延拓则不适用。
4.为了得到依据地面高精度的重力异常数据的产品，需要针对以上不足之处进行改进完善。


技术实现要素：

5.本发明提供一种基于小波变换的重力异常计算方法及装置，用以解决现有技术中已有向下延拓方法法方程条件数过大、频率域对误差向下延拓方向的不足的缺陷，利用小波变换的方法实现空中重力异常向地面延拓。
6.本发明提供一种基于小波变换的重力异常计算方法，包括：获取空中重力异常信号；对所述空中重力异常信号进行粗差探测，获得无粗差的空中重力异常信号；对所述无粗差的空中重力异常信号进行小波变换，获得小波细节系数和小波近似值；对所述小波细节系数进行误差消除获得小波细节系数估值，根据所述小波细节系数估值和所述小波近似值进行小波反变换，获得重构的空中重力异常信号；采用点对点的方法将所述重构的空中重力异常信号延拓至地面，获得地面重力异常信号。
7.根据本发明提供的一种基于小波变换的重力异常计算方法，所述对所述无粗差的空中重力异常信号进行小波变换，获得小波细节系数和小波近似值，包括：根据正交性、对称性和消失矩的特征选取小波函数，根据所述无粗差的空中重力
异常信号选取阶数，采用所述选取的小波函数对所述无粗差的空中重力异常信号在所述阶数进行小波变换，获得小波细节系数和小波近似值。
8.根据本发明提供的一种基于小波变换的重力异常计算方法，所述选取的小波函数为bior3.5小波函数，所述阶数为4阶。
9.根据本发明提供的一种基于小波变换的重力异常计算方法，所述对所述空中重力异常信号进行粗差探测，获得无粗差的空中重力异常信号，包括：对所述空中重力异常信号进行粗差探测，删除粗差，获得无粗差的空中重力异常信号。
10.根据本发明提供的一种基于小波变换的重力异常计算方法，所述对所述小波细节系数进行误差消除获得小波细节系数估值，包括：设置阈值，对大于阈值的所述小波细节系数置零，获得误差消除后的所述小波细节系数估值。
11.根据本发明提供的一种基于小波变换的重力异常计算方法，采用点对点的方法将所述重构的空中重力异常信号延拓至地面，获得地面重力异常信号，包括：获取所述空中重力异常信号在地面上点对点投影位置处的原始地面重力异常信号；根据所述空中重力异常信号和所述原始地面重力异常信号获取改正项；根据所述改正项和所述重构的空中重力异常信号获得地面重力异常信号。
12.本发明还提供一种基于小波变换的重力异常计算装置，包括：空中采集模块，获取空中重力异常信号；粗差消除模块，对所述空中重力异常信号进行粗差探测，获得无粗差的空中重力异常信号；小波变换模块，用于对所述无粗差的空中重力异常信号进行小波变换，获得小波细节系数和小波近似值；还用于对小波细节系数进行误差消除获得小波细节系数估值，根据所述小波细节系数估值和所述小波近似值进行小波反变换，获得重构的空中重力异常信号；对地延拓模块，用于采用点对点的方法将所述重构的空中重力异常信号延拓至地面，获得地面重力异常信号。
13.根据本发明提供的一种基于小波变换的重力异常计算装置，所述对地延拓模块包括：投影信号获取模块，用于获取所述空中重力异常信号在地面上点对点投影位置处的原始地面重力异常信号；改正项计算模块，用于根据所述空中重力异常信号和所述原始地面重力异常信号获取改正项；输出模块，用于根据所述改正项和所述重构的空中重力异常信号获得地面重力异常信号。
14.本发明还提供一种电子设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述程序时实现如上述任一种所述基于小波变换的重力异常计算方法的步骤。
15.本发明还提供一种非暂态计算机可读存储介质，其上存储有计算机程序，该计算机程序被处理器执行时实现如上述任一种所述基于小波变换的重力异常计算方法的步骤。
16.本发明提供的一种基于小波变换的重力异常计算方法及装置，通过小波变换将重力异常数据的粗差与误差采用一步法融合处理，提高了数据处理精度，采用了点对点的向下延拓方案，避免矩阵求逆对误差的放大效应。本发明可以将航空重力测量或卫星重力异常数据延拓到地面，形成覆盖均匀、范围广、精度高的重力异常模型，促进自然资源的物理勘探和国家高程基准的建立。
附图说明
17.为了更清楚地说明本发明或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作一简单地介绍，显而易见地，下面描述中的附图是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
18.图1是本发明实施例提供的基于小波变换的重力异常计算方法的流程示意图；图2是本发明实施例提供的点对点延拓的流程示意图；图3（a）是本发明实施例提供的不同方法获得的1891测线重力异常产品对比图；图3（b）是本发明实施例提供的不同方法获得的2790测线重力异常产品对比图；图3（c）是本发明实施例提供的不同方法获得的3660测线重力异常产品对比图；图4是本发明实施例提供的基于小波变换的重力异常计算装置的结构示意图；图5是本发明实施例提供的电子设备的结构示意图。
具体实施方式
19.为使本发明的目的、技术方案和优点更加清楚，下面将结合本发明中的附图，对本发明中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
20.重力异常是探测地下油、气、矿产资源指示器，是探测地球地质构造的信息载体，是国家高程基准建设的重要基础。传统地面重力测量是获取重力异常的重要手段。航空重力测量和卫星重力测量逐渐成为获取重力异常的发展趋势，但高于地面获取的重力数据不便于使用。因此需要通过重力异常向下延拓的方法来获取地面重力异常数据。
21.重力异常向下延拓是经典的不适定问题，向下延拓的方法可分为空间域和频率域进行分类。公式（1）即为频率域的向下延拓表达式，该方法表现出对数据高频部分的振荡，向下延拓后的数据通常无法直接使用。
22.（1）是实测航空重力异常，是待求的地面重力异常，是测量高度，、是傅里叶变换后、方向对应的波数。 表示对的傅里叶变换，表示对变量的傅里叶逆变换。公式（1）是向下延拓的频率域表达式。
[0023] 在空间域一般表示为积分形式，离散后可以表示为：
（2）a是系数阵。由于航空重力数据向下延拓固有的特点，矩阵a的条件数较大，在观测数据中有极小的误差时，延拓到地面的重力异常会受到严重的影响。
[0024]
基于上述问题，本发明实施例提出了一种基于小波变换的重力异常计算方法，如图1所述，方法包括如下步骤：步骤101、获取空中重力异常信号；步骤102、对所述空中重力异常信号进行粗差探测，获得无粗差的空中重力异常信号；步骤103、对所述无粗差的空中重力异常信号进行小波变换，获得小波细节系数和小波近似值；对所述小波细节系数进行误差消除获得小波细节系数估值，根据所述小波细节系数估值和所述小波近似值进行小波反变换，获得重构的空中重力异常信号；步骤104、采用点对点的方法将所述重构的空中重力异常信号延拓至地面，获得地面重力异常信号。
[0025]
需要说明的是，所述空中重力异常信号中需包括表示重力异常数据位置的经度、纬度和大地高，除此之外还应该有重力异常数据观测时刻和重力异常值的信息。
[0026]
本发明实施例提供的基于小波变换的重力异常计算方法首次完整提出了使用小波变换的方法进行了重力异常向下延拓，形成地面可用的重力异常产品的解决方案。通过小波变换将重力异常数据的粗差与误差采用一步法融合处理，提高了数据处理精度，采用了点对点的向下延拓方案，避免矩阵求逆对误差的放大效应。本发明还可以实现从卫星高度处的重力异常向地面延拓，扩展了小波变换进行向下延拓的空间范围。
[0027]
在本发明的至少一个实施例中，所述对所述无粗差的空中重力异常信号进行小波变换，获得小波细节系数和小波近似值，包括：根据正交性、对称性和消失矩的特征选取小波函数，根据所述无粗差的空中重力异常信号选取阶数，采用所述选取的小波函数对所述无粗差的空中重力异常信号在所述阶数进行小波变换，获得小波细节系数和小波近似值。
[0028]
需要说明的是，根据无粗差的空中重力异常信号需进行多尺度分析的特征，本发明实施例中选择正交性、对称性和消失矩作为特征来选取小波函数，具体要求如下所示：（1）小波函数应具备完备正交性的要求，该要求也是多尺度分析的基本要求。
[0029]
（2）小波函数应该是正交的。正交小波可保证多尺度分析的低阶细节不变性。
[0030]
（3）小波函数应具有对称性。该特征与信号的重构效果有直接的关系。若对称性好，重构信号的失真度就小，能较好的逼近输入信号。
[0031]
（4）小波函数的消失矩应尽量高，这个要求与小波函数逼近光滑函数时的收敛率有直接的关系。函数的高阶变化可以用高阶消失矩小波函数进行分析。
[0032]
针对重力异常数据，在分解的层次的选择过多或过少都是不合适的，因此本发明分析了各阶小波变换所具有的地球物理含义，从而根据所述无粗差的空中重力异常信号选取阶数：1阶小波变换重力细节d1主要反映上地壳浅层密度不均匀体的分布；2阶小波细节d2主要反映了中上地壳浅层密度不均匀体分布；3阶小波变换细节d3中，在d1和d2中弱小的
异常逐渐消失，合并显示为连成较大范围的异常；4阶小波细节d4中，连成片的特点更为明显，它反映了区域地质构造特征。即浅部异常变化特征少，深部异常信号突出；5阶小波细节的特征主要表现了盆地和造山带的轮廓与构造；6阶小波细节表现的特征只是将5阶细节由小波伸缩得到的，仅仅增大了异常特征的波长，在更大尺度上表现了异常特征。
[0033]
在本发明的至少一个实施例中，所述选取的小波函数为bior3.5小波函数，所述阶数为4阶。
[0034]
由于在实际的数据处理过程中，完全正交对称小波是不存在的。在处理中尽量选用具有较好对称性，消失矩较高的正交小波函数。优选的，在后续步骤中选择bior3.5的小波函数进行重力异常数据的向下延拓。进行更高阶的分解，能反映出地球内部更深的结构，但会偏离了利用小波变换提升重力异常产品延拓到地面形成可用产品的初衷。对于航空重力测量而言，将重力异常数据分解到4解细节是合适的。优选的，将小波分解确定为4阶。
[0035]
在本发明的至少一个实施例中，所述对所述空中重力异常信号进行粗差探测，获得无粗差的空中重力异常信号，包括：对所述空中重力异常信号进行粗差探测，删除粗差，获得无粗差的空中重力异常信号。在本实施例中取数据的3倍标准差作为约束进行粗差剔除。
[0036]
在本发明的至少一个实施例中，所述对所述小波细节系数进行误差消除获得小波细节系数估值，包括：设置阈值，对大于阈值的所述小波细节系数置零，获得误差消除后的所述小波细节系数估值。
[0037]
需要说明的是，误差消除可以选择合适的阈值对细节信号进行软处理，减少数据中存在的误差。所述阈值为小波细节系数3倍的中误差。
[0038]
在本发明的至少一个实施例中，采用点对点的方法将所述重构的空中重力异常信号延拓至地面，获得地面重力异常信号，包括：步骤201、获取所述空中重力异常信号在地面上点对点投影位置处的原始地面重力异常信号；需要说明的是，所述空中重力异常信号和原始地面重力异常信号根据全球重力场模型获得，且均需要进行粗差处理。
[0039]
步骤202、根据所述空中重力异常信号和所述原始地面重力异常信号获取改正项，；步骤203、根据所述改正项和所述重构的空中重力异常信号获得地面重力异常信号；需要说明的是，对航空重力测线的重力异常添加改正项，实现观测数据向下延拓，即。即为地面重力异常成果。
[0040]
在本发明的至少一个实施例中，利用点对点的向下延拓方法，将延拓到地面。即是需要的重力异常数据，可以进行后的应用研究。其实现步骤如下：（1）利用全球重力场模型，计算航空重力测线位置处的重力异常异常模型值。对其
实行步骤
②
、步骤
③
中去粗差处理，得到测线位置处的重力异常模型。
[0041]
（2）利用全球重力场模型，计算航空重力测线在地面上投影位置处的重力异常异常模型值。对其实行步骤
②
、步骤
③
中去粗差处理，得到测线位置处的重力异常模型。
[0042]
（3）基于和计算，测线高度处与对应地面点位置处两点之间的重力异常改正数。即。
[0043]
（4）对航空重力测线的重力异常添加改正项，实现观测数据向下延拓，即。即为地面重力异常成果。
[0044] 本发明实施例首次基于小波变换的方法，借助于全球重力场模型，求取重力异常的改正数，实现空间数据向下延拓。
[0045]
对比例本发明公开了一种重力异常计算方法，旨在根据建立的观测值与预测值的方差与协方差函数模型，利用配置法由大地水准面转计算重力异常。在计算过程中，考虑了地形对大地水准面及重力异常的影响。方法包括如下步骤：步骤、利用已有延拓方法对观测信号，比如航空重力数据，做向下延拓，其中，已有延拓方法如快速傅里叶变换（fft）、最小二乘配置等；步骤、计算观测信号和延拓信号的多尺度边缘；步骤
③
、修改延拓信号的多尺度边缘，即只保留与观测信号的多尺度边缘位置相同的边缘；步骤
④
、利用修改后的多尺度边缘重构延拓值。
[0046]
本发明对比例的缺陷表现在：（1）在步骤中，在不考虑粗差的情况下，使用fft技术或配置法等方法进行向下延拓时，会造成误差迅速扩大，甚至会超过有用的信号。若数据中存在粗差，此步骤的不足将会表现的更突出。
[0047]
（2）方法的适应范围有限，尤其不能应用于卫星重力的向下延拓，而且在多尺度边缘分析时存在阈值和尺度等几个参数的选取问题上，需要开展进一步的研究工作。
[0048]
实施例1在上述实施例的基础上，本发明实施例公开了一种基于小波变换的重力异常计算方法，如图4所示，在我国某海域进行了航空重力向下延拓，并与船载重力测量数据进行了对比，图中网格线表示航空重力测量的飞行轨迹，粗实线表示测试的三条测线1891、2790与3660。方法包括如下步骤：步骤、收集空中重力异常数据；步骤、根据小波函数的性质和重力异常信号之间的关系，选择合适的小波函数；步骤
③
、对所述重力异常数据进行粗差探测，并将粗差剔除，得到无粗差的重力数据；步骤
④
、对无粗差的重力数据，利用合适的小波函数进行小波变换
，得到重力异常信号的细节和近似值；步骤
⑤
、选择合适的阈值对细节信号进行软处理，减少数据中存在的误差，得到细节；需要说明的是，所述软处理采用一个连续光滑的函数，使得在阈值附近的数据变化平缓，不存在突跳。
[0049]
步骤
⑥
、根据和恢复重力异常信号；步骤
⑦
、利用点对点的向下延拓方法，将延拓到地面。其中，即是需要的重力异常数据，与船载重力数据进行对比，可评价其精度。
[0050]
本发明实施例中的基于小波变换的重力异常计算方法，无论输入数据（航空重力测量数据或卫星重力数据）是否存在粗差，都可以通过引入小波变换的方法和对其进行多尺度分析，避免引入不可度量的误差。本发明实施例采用点对点的向下延拓方案，避免向下延拓问题的法方程求逆，增强求解的稳定性和计算效率。同时，本发明还能发挥小波变换的优势，将该优势方法引入到重力数据处理领域。
[0051]
如图3（a）、图3（b）和图3（c）所示，是本发明实施例提供的不同方法获得的重力异常产品对比图，该图中的各测线处采用航空重力数据依据不同方法获得的重力异常数值模型的具体数值如下表：由上述数据分析可得，本发明实施例的基于小波变换的重力异常计算方法可使航空重力测量数据向地面延拓时计算结果稳定。计算效率提高，对计算内存要求也不高。将有关模型数据准备后，一次性导入某测区全部测线数据，依据本发明方案进行点对点延拓，避免频率算法的边界效应、放大高频误差，且不需要额外开辟计算机内存。
[0052]
下面对本发明提供的基于小波变换的重力异常计算装置进行描述，下文描述的基于小波变换的重力异常计算装置与上文描述的基于小波变换的重力异常计算方法可相互对应参照。如图4所示，本发明实施例的基于小波变换的重力异常计算装置，包括：空中采集模块401，获取空中重力异常信号；粗差消除模块402，对所述空中重力异常信号进行粗差探测，获得无粗差的空中重力异常信号；小波变换模块403，用于对所述无粗差的空中重力异常信号进行小波变换，获得小波细节系数和小波近似值；还用于对小波细节系数进行误差消除获得小波细节系数估值，根据所述小波细节系数估值和所述小波近似值进行小波反变换，获得重构的空中重力异常信号；对地延拓模块404，用于采用点对点的方法将所述重构的空中重力异常信号延拓至地面，获得地面重力异常信号。
[0053]
本发明实施例提供的基于小波变换的重力异常计算装置首次完整提出了使用小波变换的方法进行了重力异常向下延拓，形成地面可用的重力异常产品的解决方案。通过小波变换将重力异常数据的粗差与误差采用一步法融合处理，提高了数据处理精度，采用了点对点的向下延拓方案，避免矩阵求逆对误差的放大效应。本发明还可以实现从卫星高度处的重力异常向地面延拓，扩展了小波变换进行向下延拓的空间范围。
[0054]
在本发明的至少一个实施例中，所述对地延拓模块404包括：投影信号获取模块，用于获取所述空中重力异常信号在地面上点对点投影位置处的原始地面重力异常信号；改正项计算模块，用于根据所述空中重力异常信号和所述原始地面重力异常信号获取改正项；输出模块，用于根据所述改正项和所述重构的空中重力异常信号获得地面重力异常信号。
[0055]
图5示例了一种电子设备的实体结构示意图，如图5所示，该电子设备可以包括：处理器(processor)510、通信接口(communications interface)520、存储器(memory)530和通信总线540，其中，处理器510，通信接口520，存储器530通过通信总线540完成相互间的通信。处理器510可以调用存储器530中的逻辑指令，以执行基于小波变换的重力异常计算方法，该方法包括：获取空中重力异常信号；对所述空中重力异常信号进行粗差探测，获得无粗差的空中重力异常信号；对所述无粗差的空中重力异常信号进行小波变换，获得小波细节系数和小波近似值；对所述小波细节系数进行误差消除获得小波细节系数估值，根据所述小波细节系数估值和所述小波近似值进行小波反变换，获得重构的空中重力异常信号；采用点对点的方法将所述重构的空中重力异常信号延拓至地面，获得地面重力异常信号。
[0056]
此外，上述的存储器530中的逻辑指令可以通过软件功能单元的形式实现并作为独立的产品销售或使用时，可以存储在一个计算机可读取存储介质中。基于这样的理解，本发明的技术方案本质上或者说对现有技术做出贡献的部分或者该技术方案的部分可以以软件产品的形式体现出来，该计算机软件产品存储在一个存储介质中，包括若干指令用以使得一台计算机设备（可以是个人计算机，服务器，或者网络设备等）执行本发明各个实施例所述方法的全部或部分步骤。而前述的存储介质包括：u盘、移动硬盘、只读存储器（rom，read-only memory）、随机存取存储器（ram，random access memory）、磁碟或者光盘等各种可以存储程序代码的介质。
[0057]
另一方面，本发明还提供一种计算机程序产品，所述计算机程序产品包括计算机程序，计算机程序可存储在非暂态计算机可读存储介质上，所述计算机程序被处理器执行时，计算机能够执行上述各方法所提供的基于小波变换的重力异常计算方法，该方法包括：获取空中重力异常信号；对所述空中重力异常信号进行粗差探测，获得无粗差的空中重力异常信号；对所述无粗差的空中重力异常信号进行小波变换，获得小波细节系数和小波近似值；对所述小波细节系数进行误差消除获得小波细节系数估值，根据所述小波细节系数估
值和所述小波近似值进行小波反变换，获得重构的空中重力异常信号；采用点对点的方法将所述重构的空中重力异常信号延拓至地面，获得地面重力异常信号。
[0058]
又一方面，本发明还提供一种非暂态计算机可读存储介质，其上存储有计算机程序，该计算机程序被处理器执行时实现以执行上述各方法提供的基于小波变换的重力异常计算方法，该方法包括：获取空中重力异常信号；对所述空中重力异常信号进行粗差探测，获得无粗差的空中重力异常信号；对所述无粗差的空中重力异常信号进行小波变换，获得小波细节系数和小波近似值；对所述小波细节系数进行误差消除获得小波细节系数估值，根据所述小波细节系数估值和所述小波近似值进行小波反变换，获得重构的空中重力异常信号；采用点对点的方法将所述重构的空中重力异常信号延拓至地面，获得地面重力异常信号。
[0059]
以上所描述的装置实施例仅仅是示意性的，其中所述作为分离部件说明的单元可以是或者也可以不是物理上分开的，作为单元显示的部件可以是或者也可以不是物理单元，即可以位于一个地方，或者也可以分布到多个网络单元上。可以根据实际的需要选择其中的部分或者全部模块来实现本实施例方案的目的。本领域普通技术人员在不付出创造性的劳动的情况下，即可以理解并实施。
[0060]
通过以上的实施方式的描述，本领域的技术人员可以清楚地了解到各实施方式可借助软件加必需的通用硬件平台的方式来实现，当然也可以通过硬件。基于这样的理解，上述技术方案本质上或者说对现有技术做出贡献的部分可以以软件产品的形式体现出来，该计算机软件产品可以存储在计算机可读存储介质中，如rom/ram、磁碟、光盘等，包括若干指令用以使得一台计算机设备（可以是个人计算机，服务器，或者网络设备等）执行各个实施例或者实施例的某些部分所述的方法。
[0061]
最后应说明的是：以上实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的精神和范围。