本技术属于海洋地震勘探领域，具体涉及一种多震源炮点预测方法、系统、终端及可读存储介质。
背景技术：
1、海洋地震勘探作业过程为：勘探船拖拽着震源阵列沿着设计测线(前绘测线)航行并在设计的炮点位置激发地震波，检波器接收经过海底岩层分界面反射或折射的地震波，通过对地震波资料进行处理和解释可以确定海底地震构造，从而达到探测海洋油气资源的目的。其中，实际响炮位置与设计炮点位置越接近，地震数据采集精度越高，因此，需要实现一种高精度、高稳定的炮点预测方法。目前海洋地震勘探炮点预测方法是：基于勘探船沿着测线直线匀速行走的作业场景，使用多项式模型、常速度模型、常加速度模型等定位模型用于估计勘探船和震源阵列的运动状态。然而，上述模型过于理想化无法适用于复杂海况环境下震源阵列导航定位。同时，炮点预测过程中将震源阵列当成刚体使用其几何中心点(mgp，mean gun position)进行炮点预测，但实际作业中震源之间并没有硬件装置连接且受到风浪影响震源阵列会发生倾斜，因此将mgp当成响炮参考点进行炮点预测方法会影响地震数据采集精度。
2、因此，亟需一种可以克服上述技术缺陷的炮点预测方法。
技术实现思路
1、本技术的目的是，提供一种多震源炮点预测方法、系统、终端及可读存储介质，能够提高海洋地震勘探炮点预测精度。
2、本技术提供的技术方案为：
3、第一方面，本技术提供一种多震源炮点预测方法，包括以下步骤：
4、步骤1：获取勘探船和震源阵列中各个震源上的数据采集装置采集的观测数据；
5、步骤2：对观测数据进行预处理，去除粗差数据并进行时间同步；
6、步骤3：根据预处理后的观测数据计算当前时刻震源阵列中各个震源的坐标；分别将当前时刻震源的坐标作为观测向量，输入基于卡尔曼滤波的震源导航定位模型，得到各个震源的空间状态信息；其中，所述基于卡尔曼滤波的震源导航定位模型具有根据输入的观测向量，计算相应震源的空间状态信息的功能；所述震源的空间状态信息包括震源的坐标；
7、步骤4：基于各个震源的空间状态信息以及设计炮点的位置，确定下次响炮对应的设计炮点及下次响炮震源。
8、在一种可能的实现方式中，所述步骤1中，所述勘探船和震源阵列中各个震源上的数据采集装置，包括：勘探船上的区域增强全球导航卫星系统(regionalglobalnavigation satellite system，简称rgnss)基准站、勘探船上的差分全球卫星导航系统(differential global navigation satellite system，简称dgnss)和电罗经，以及震源阵列中各个震源上的rgnss流动站；
9、所述步骤2，包括：基于几何结构约束和数理统计特性对观测数据进行粗差探测，去除粗差数据；
10、所述步骤2，还包括：对于各一震源上的rgnss流动站，若当前时刻不存在有效观测数据，则基于历史观测数据进行观测数据时间同步，获取当前时刻的观测数据；其中，所述基于历史观测数据进行观测数据时间同步，获取当前时刻的观测数据，包括：
11、若当前时刻两侧均存在有效观测数据，则使用线性插值方法获取当前时刻的观测数据；
12、若当前时刻左侧或右侧不存在有效观测数据时，则通过对历史观测数据进行多项式拟合外推得到当前时刻的观测数据，多项式拟合公式为：
13、y＝a0+a1t+a2t2+lantn
14、式中：y为t时刻对应的观测数据，t为时刻，an为多项式系数，n为多项式阶数；
15、通过对历史观测数据进行多项式拟合得到时刻与观测数据之间的拟合关系，确定多项式拟合公式中的多项式系数；再将当前时刻代入，得到当前时刻的观测数据。
16、在一种可能的实现方式中，所述步骤3中，根据预处理后的观测数据计算当前时刻震源阵列中各个震源的坐标，包括：
17、
18、式中：(xref,yref)为勘探船上的rgnss基准站坐标，(xdgnss,ydgnss)为勘探船上的dgnss坐标；α为电罗经采集的勘探船方位角，(δx,δy)为勘探船上的rgnss基准站与勘探船上的dgnss之间的相对偏移值；(x,y)为k时刻的震源的坐标，l和β分别为rgnss基准站和震源上rgnss流动站之间的方位和斜距。
19、在一种可能的实现方式中，所述步骤3中，所述基于卡尔曼滤波的震源导航定位模型，包括：
20、xk/k-1＝φk/k-1xk-1
21、
22、xk＝xk/k-1+kk(zk-hkxk/k-1)
23、pk＝(i-kkhk)pk/k-1
24、式中，xk为k时刻的状态向量，即k时刻的震源的空间状态信息；xk-1为k-1时刻的状态向量；pk和pk-1分别为xk和xk-1的协方差矩阵；zk为k时刻观测向量，即利用步骤2预处理后的观测数据计算得到的k时刻的震源的坐标，xk/k-1为k-1时刻对k时刻状态向量预测值；φk/k-1为k-1时刻对k时刻的状态转移矩阵；pk/k-1为xk/k-1的协方差矩阵，qk-1为k-1时刻的过程噪声矩阵的协方差矩阵；kk为k时刻的增益矩阵，hk为系数矩阵；i为单位矩阵；rk为k时刻的观测噪声矩阵vk的协方差矩阵。
25、在一种可能的实现方式中，所述方法还包括：
26、对k时刻的观测噪声矩阵vk的协方差矩阵rk进行调整，得到然后将作为rk代入到kk的计算公式中，进而求解各个震源的空间状态信息；
27、其中，对k时刻的观测噪声矩阵vk的协方差矩阵rk进行调整，得到计算公式为：
28、vk/k-1＝zk-hkxk/k-1
29、
30、式中，为调整后的k时刻的观测噪声矩阵vk的协方差矩阵；vk/k-1和分别为卡尔曼滤波的新息及其协方差矩阵，vi为第i个观测数据的标准化残差；ri为第i个观测数据对应的降权因子，k0和k1为相应的临界值，esp为极小值。
31、在一种可能的实现方式中，针对上述实施例中的多震源船，先根据震源阵列中各个震源与设计炮点的位置关系分别找到与各个震源最接近的下一设计炮点的炮号(即沿着设计测线方向前进最接近的设计炮点的编号)，再根据与各个震源最接近的下一设计炮点的炮号的大小关系及前绘给定的设计炮点的炮号增量(沿着设计测线方向设计炮点的编号的变化量)，判断出下次响炮对应的设计炮点(即预测炮点)及下次响炮震源，判断方法为：
32、设与震源ei最接近的下一设计炮点的炮号为mi，其中i＝1,2,...,n，n为震源阵列中震源的数量；若炮号增量大于0，则找出{mi,i＝1,2,...,n}中的最小值，该最小值对应的炮点即为下次响炮对应的设计炮点，该最小值对应的震源即为下次响炮震源；若炮号增量小于0，则找出{mi,i＝1,2,...,n}中的最大值，该最大值对应的炮点即为下次响炮对应的设计炮点，该最大值对应的震源即为下次响炮震源。
33、在一种可能的实现方式中，所述步骤3中，所述震源的空间状态信息，还包括震源阵列中各个震源的速度；
34、所述方法还包括：步骤5：根据下次响炮对应的设计炮点与下次响炮震源的相对位置关系，以及下次响炮震源的速度，计算剩余响炮距离ds及剩余响炮时间dt，用于枪控系统根据剩余响炮时间判断是否激发地震波。
35、在一种可能的实现方式中，所述步骤5中，根据下次响炮对应的设计炮点与下次响炮震源的相对位置关系，计算剩余响炮距离ds及剩余响炮时间dt，包括：
36、计算下次响炮对应的设计炮点相对于测线起点的沿测线方向距离(distancealong-line，简称da)dapre；
37、dapre＝sprecos(θ-α)
38、式中，spre和θ分别为下次响炮对应的设计炮点与测线起点连线的距离和方位角，α为测线方位角；
39、计算下次响炮震源相对于测线起点的沿测线方向距离dasource：
40、dasource＝ssourcecos(θ-α)
41、式中，ssource和λ分别为下次响炮震源与测线起点连线的距离和方位角，根据步骤3得到的各个震源的空间状态信息与测线起点的位置信息计算得到；α为测线方位角；
42、基于dapre和dasource，计算剩余响炮距离ds及剩余响炮时间dt并发送给枪控系统，包括：
43、ds＝da-da
44、ds＝dapre-dasource
45、dt＝|ds/vda|
46、式中：vda为下次响炮震源沿测线方向速度，根据步骤3得到的各个震源的空间状态信息获取。
47、第二方面，本技术提供一种多震源炮点预测系统，包括：
48、观测数据获取模块，用于获取勘探船和震源阵列中各个震源上的数据采集装置采集的观测数据；
49、数据预处理模块，用于对观测数据进行预处理，去除粗差数据并进行时间同步；
50、震源空间状态信息计算模块，用于根据预处理后的观测数据计算当前时刻震源阵列中各个震源的坐标；分别将当前时刻震源的坐标作为观测向量，输入基于卡尔曼滤波的震源导航定位模型，得到各个震源的空间状态信息；其中，所述基于卡尔曼滤波的震源导航定位模型具有根据输入的观测向量，计算相应震源的空间状态信息的功能；所述震源的空间状态信息包括震源的坐标；
51、炮点预测模块，用于基于各个震源的空间状态信息以及设计炮点的位置，确定下次响炮对应的设计炮点及下次响炮震源。
52、在一种可能的实现方式中，所述炮点预测模块，还用于根据下次响炮对应的设计炮点与下次响炮震源的相对位置关系，计算剩余响炮距离ds及剩余响炮时间dt。
53、第三方面，本技术提供一种终端，包括：存储器和处理器；
54、所述存储器，用于存储计算机程序；
55、所述处理器，用于调用所述计算机程序，以执行如上所述的多震源炮点预测方法。
56、第四方面，本技术提供一种计算机可读存储介质，所述计算机可读存储介质中存储有计算机程序，所述计算机程序在电子设备上运行时，使得所述电子设备实现如上所述的多震源炮点预测方法。
57、第五方面，本技术提供一种计算机程序产品，包括计算机程序，所述计算机程序在电子设备上运行时，使得所述电子设备实现如上所述的多震源炮点预测方法。
58、本技术上述第二到第五方面的具体实现方式可以参考上述第一方面的实现方式，在此不进行赘述。
59、有益效果：
60、1.本技术提供的一种多震源炮点预测方法，其相较于现有技术，基于卡尔曼滤波的震源导航定位模型更能够真实描述震源实际运动轨迹，同时能够抵抗粗差对炮点预测精度的影响。
61、2.本技术采用左右侧震源独立预测的方法，相较于传统采用震源阵列几何中心预测炮点，更能够满足海洋地震勘探需要，能够保证海洋地震勘探炮点预测精度需求。
62、3.本技术使用几何结构约束和数理统计相结合的粗差探测方法，更能够为炮点预测提供“干净”的数据。
1.一种多震源炮点预测方法，其特征在于，包括以下步骤：
2.根据权利要求1所述的方法，其特征在于，所述步骤1中，所述勘探船和震源阵列中各个震源上的数据采集装置，包括：勘探船上的rgnss基准站、dgnss和电罗经，以及震源阵列中各个震源上的rgnss流动站；
3.根据权利要求1所述的方法，其特征在于，所述步骤3中，所述基于卡尔曼滤波的震源导航定位模型，包括：
4.根据权利要求1所述的方法，其特征在于，所述方法还包括：
5.根据权利要求1所述的方法，其特征在于，针对上述实施例中的多震源船，先根据震源阵列中各个震源与设计炮点的位置关系分别找到与各个震源最接近的下一设计炮点的炮号，再根据与各个震源最接近的下一设计炮点的炮号的大小关系及设计给定的设计炮点的炮号增量，判断出下次响炮对应的设计炮点及下次响炮震源，判断方法为：
6.根据权利要求1所述的方法，其特征在于，所述步骤3中，所述震源的空间状态信息，还包括震源阵列中各个震源的速度；
7.一种多震源炮点预测系统，其特征在于，包括：
8.一种终端，其特征在于，包括：存储器和处理器；
9.一种计算机可读存储介质，其特征在于，所述计算机可读存储介质中存储有计算机程序，所述计算机程序在电子设备上运行时，使得所述电子设备实现如权利要求1至6中任一项所述的方法。
10.一种计算机程序产品，包括计算机程序，其特征在于，所述计算机程序在电子设备上运行时，使得所述电子设备实现如权利要求1至6中任一项所述的方法。
