一种监测承压水储量变化的卫星重力协同正演模拟方法

本发明涉及重力卫星与地下水科学交叉领域，具体而言，涉及一种监测承压水储量变化的卫星重力协同正演模拟方法。

背景技术：

1、地下水作为地球上重要的淡水储备之一，由浅层和深层地下水两部分组成，在生态平衡的维护、人类社会可持续发展的支撑以及生命健康的保障方面扮演着重要角色。特别是在干旱、半干旱地区和水资源稀缺地带，地下水的存在为人类提供了一种相对稳定的水源，填补了地表水供应不足的空白，满足了农业、工业、城市等多个领域的水需求。在水资源供需紧张的背景下，地下水的利用被视为缓解水资源压力的关键策略。

2、深层地下水(承压水)，相较于浅层地下水，具有更高的稳定性。然而，在高耗水量的地区，随着人类活动的不断扩展和水资源利用的进一步加剧，对深层地下水的长期过度开采使地下水正面临着不容忽视的紧迫压力，部分地区呈现出不同程度的地面沉降现象。因此，定量分析承压含水层系统中可开采的地下水总储量至关重要。深入探究深层地下水的空间分布、补给机制，以及与浅层地下水的相互关系，将有助于更加有效地实施保护和合理利用的策略，为未来的水资源管理和可持续发展提供可靠的科学基础。

3、目前，对于深层地下水的研究可采取的方法主要有利用传统监测井水位数据监测、地球物理探测方法、通过模拟地下水流动与补给情况构建数值模型以及利用热红外、微波遥感技术实现间接监测等几种。前两者对于小尺度空间的监测效果较佳，但监测点的数量及空间分布会直接影响研究结论，后两者虽然可以实现大尺度的监测，但构建模型需要输入大量的先验知识与参数，并且监测结果的不确定性较大。结合上述各方法的优势与局限，本文提出一种将监测井水位数据与遥感技术相结合的深层地下水储监测方法。

4、近年来，随着美德联合研制的grace(gravity recovery and climateexperiment)和grace-fo重力卫星的成功发射和持续发展，国内外学者在监测陆地水储量和地下水储量变化等方面进行了长期而大量的研究工作。作为目前能够有效监测地表以下水储量变化的遥感卫星，它们克服了传统地面观测空间局限性的问题，为地下水科学提供了宝贵的数据资源，该卫星目前已成为全球和大区域水储量变化监测的有效工具。

5、在对grace/grace-fo数据进行处理的过程中，截断球谐系数、滤波等操作会导致信号产生不同程度的向内或向外的泄漏，从而产生泄露误差。因此，在得出结论前需对所得数据进行校正，以保证恢复正常信号。在目前常用的校正方法中，加法校正、尺度因子校正方法的有效程度均严重依赖于水文模型数据的数量与质量，对于迭代正演模拟方法来说，尽管其克服了前者对于先验知识的过度依赖问题，但其恢复信号空间分辨率的能力高度又依赖于真实信号的空间分布的均匀程度，因此对地下水来说，尤其是承压水，传统的迭代正演模拟方法也很难实现高质量、高精度的信号恢复。

技术实现思路

1、本发明提供一种监测承压水储量变化的卫星重力协同正演模拟方法，用以解决上述现有技术存在的问题。

2、为达到上述目的，本发明提供了一种监测承压水储量变化的卫星重力协同正演模拟方法，其包括：

3、s1：获取预设研究区边界及经纬度信息，并根据误差可能影响的范围确定预设背景区的边界及经纬度信息；

4、s2：基于预设背景区的边界及经纬度信息分别获取月尺度的监测井水位数据、背景区给水度经验值、grace/grace-fo月尺度的重力卫星时变重力场数据以及非地下水组分数据，非地下水组分数据包括土壤含水量数据、雪水当量数据、冠层水量数据、地表水储量数据；

5、s3：将月尺度的grace/grace-fo重力卫星时变重力场数据截断至一定阶数的球谐系数并进行处理，扣除一定时间段的平均重力场得到含有泄漏误差的陆地水储量变化月时间序列，同时，对其余非地下水各组分数据展开为球谐系数并采取同样的数据处理方式，得到非地下水组分变化月时间序列，从含有泄漏误差的陆地水储量变化月时间序列中扣除非地下水组分变化月时间序列，得到含有泄露误差的地下水储量变化月时间序列；

6、s4：根据监测井类型分离潜水井水位变化月时间序列与承压水井水位变化月时间序列，之后将潜水井水位数据与s2中获取的给水度经验值根据格网的经纬度对应相乘，再将其展开为球谐系数，采取与s3中相同的数据处理方式，得到潜水储量变化月时间序列；

7、s5：在地下水储量变化月时间序列中扣除潜水储量变化月时间序列，得到含有泄露误差的0.5°格网上的承压水储量变化月时间序列；

8、s6：根据s4得到的承压水井水位变化月时间尺度数据以及s5得到的含有泄露误差的承压水储量变化月时间序列，通过给定符合水文地质条件范围内的任意值作为弹性释水系数的模拟初始值，运行预设的迭代正演模拟算法；

9、s7：基于预设的迭代正演模拟算法，获得预设研究区地下水储量变化月时间序列反演值所对应的弹性释水系数值，将弹性释水系数值与s4得到的承压水井水位变化月时间尺度数据根据空间位置相乘，得到预设研究区校正后的承压水储量变化月时间序列。

10、在本发明的一实施例中，步骤s2中，背景区给水度经验值是以预设研究区给水度经验值的空间分布图为基础，通过面积加权平均算法计算背景区内0.5°格网尺度上的给水度值。

11、在本发明的一实施例中，步骤s3中，对grace/grace-fo重力卫星时变重力场数据进行以下处理：

12、使用截断至60阶的球谐数据；

13、采用卫星激光测距测得的数值对c20项进行替换；

14、对球谐数据作半径为300km的高斯滤波并结合p3m10去相关处理，即对10次及以上的球谐系数按照奇、偶数阶分别用3阶多项式拟合一次并扣除。

15、在本发明的一实施例中，步骤s4中，采取反距离加权插值法对潜水井水位变化月时间序列在背景区范围内插值，根据近邻点的加权均值来估计未知点，以插值点与已知的样本点的距离为权重来加权平均，离插值点越近的样本点赋予的权重越大，以获取在背景区范围上空间连续的0.5°格网尺度上的潜水井水位变化月时间序列；

16、将潜水井水位变化月时间序列与步骤s2中得到的0.5°格网尺度上的给水度值，根据格网的经纬度对应相乘，得0.5°格网尺度上的潜水储量变化月时间序列，计算公式如下：

17、sgwsa＝sgwla×sy

18、式中：sgwsa为潜水储量变化月时间序列，sgwla为潜水井水位变化月时间序列，sy为给水度经验值；

19、将潜水储量变化月时间序列其展开为球谐系数，采取与步骤s3中相同的数据处理方式进行正演模拟，得到正演模拟后的潜水储量变化月时间序列。

20、在本发明的一实施例中，步骤s5的计算公式如下：

21、dgwsash＝gwsash-sgwsash

22、式中：dgwsash为待校正的承压水储量变化月时间序列，gwsash为含有泄露误差的地下水储量变化月时间序列，sgwsash为正演模拟后的潜水储量变化月时间序列。

23、在本发明的一实施例中，步骤s6运行预设的迭代正演模拟算法时，通过比较迭代模拟后的承压水储量变化月时间序列与含有泄露误差的承压水储量变化月时间序列的差值与给定阈值的大小判断是否停止迭代，计算公式如下：

24、dgwsai＝dgwla×μi

25、式中：dgwsai为第i次迭代所得出的承压水储量变化月时间序列，dgwla为研究区承压水井水位变化月时间序列，μi为第i次迭代模拟的弹性释水系数，

26、μi+1＝μi±λ

27、式中：μi+1为第i+1次迭代所需的弹性释水系数，λ为预设的迭代正演模拟算法所设定的调节步长，此处加减步长的判断依据为第i次迭代模拟后的承压水储量变化月时间序列与含有泄露误差的承压水储量变化月时间序列差值的正负情况，直至差值小于给定阈值，整个迭代正演模拟算法结束。

28、在本发明的一实施例中，步骤s7，基于最小二乘拟合的谐波分析算法对其进行信号分解，以求取周年、半周年振幅、相位和变化趋势，计算公式如下：

29、

30、式中：δh(t)为月时间序列；t为时间；a和b分别为常数项和趋势项；ai和分别为振幅和相位(i＝1，为周年；i＝2，为半周年)；ε(t)为拟合残差。

31、本发明提供的监测承压水储量变化的卫星重力协同正演模拟方法通过联合重力卫星和监测井水位数据进行迭代校正，将宏观监测技术与微观观测技术相结合，集遥感技术以及传统水文监测方法的优势于一体，以获取更加高质量的承压水层的信号，为有关承压含水层的研究提供了一种新思路。