段楚峰 張昊楠 匡翠林* 余文坤 阮福明 戴吾蛟

(①中南大學(xué)地球科學(xué)與信息物理學(xué)院，湖南長(zhǎng)沙 410083；②中海油田服務(wù)股份有限公司物探事業(yè)部，天津 300451；③海洋油氣勘探國家工程研究中心，北京100028)

0 引言

海上地震拖纜勘探是由勘探船拖拽槍陣與電纜沿測(cè)線航行，通過槍陣震源激發(fā)模擬地震波，并由拖纜上安裝的檢波器陣列接收海底反射的地震波信號(hào)以獲得海底的構(gòu)造地質(zhì)信息，從而達(dá)到油氣勘探的目的[1-4]。現(xiàn)階段海上地震拖纜勘探一般采用多纜作業(yè)模式以提高施工效率。拖纜定位是海上地震勘探作業(yè)的關(guān)鍵步驟之一，其精度將直接影響地震數(shù)據(jù)的成像精度與可靠性[5-6]。

拖纜具有非剛體特性，一般長(zhǎng)達(dá)數(shù)千米。水下環(huán)境復(fù)雜多變，如何對(duì)拖纜建立合理的數(shù)學(xué)模型并利用特定的算法對(duì)模型進(jìn)行解算，是拖纜獲得高精度纜形和位置的關(guān)鍵[7-8]。

國外較早開展拖纜建模和定位算法研究，Gilbert于1981年利用羅盤觀測(cè)值對(duì)拖纜平面坐標(biāo)之間的多項(xiàng)式系數(shù)進(jìn)行估計(jì)，確定單纜的纜形和位置[9]；又于1982年以拖纜里程為因變量，分別對(duì)單纜的三維坐標(biāo)分量進(jìn)行多項(xiàng)式擬合建模，并利用極大似然估計(jì)確定多項(xiàng)式系數(shù)[10]；Zinn等[11]將拖纜視為若干段相連的圓弧，以羅盤節(jié)點(diǎn)確定每段圓弧的圓心、直徑等參數(shù)，并按照?qǐng)A弧的幾何性質(zhì)進(jìn)行拖纜坐標(biāo)推算；Combier[12]根據(jù)羅盤提供的方位對(duì)單纜的坐標(biāo)進(jìn)行推算，再通過尾標(biāo)的精確位置對(duì)拖纜進(jìn)行定位修正；Goutorbe等[13]利用羅盤觀測(cè)值擬合整條拖纜的方位角，并沿曲線對(duì)坐標(biāo)積分獲得單纜上各點(diǎn)位置。中國于20世紀(jì)90年代開始對(duì)拖纜建模和定位方法進(jìn)行了探討[14]，姚宜斌等[15]先通過近似曲線積分的方法對(duì)拖纜進(jìn)行坐標(biāo)推算，獲得初始纜形，再通過聲學(xué)測(cè)距網(wǎng)絡(luò)的約束平差和擬穩(wěn)平差確定聲學(xué)節(jié)點(diǎn)的精確位置，然后通過分段旋轉(zhuǎn)、縮放的方式使初始纜形與平差后的聲學(xué)節(jié)點(diǎn)相吻合，該方法在多纜定位中取得了較好的效果；易昌華等[16]介紹了弧段模型與曲線積分模型，并對(duì)比分析了二者的適用場(chǎng)景和計(jì)算效率；易昌華[17]提出了利用近似曲線積分派生出前、中、后聲學(xué)網(wǎng)絡(luò)之間的虛擬距離觀測(cè)值，并進(jìn)行聲學(xué)網(wǎng)絡(luò)統(tǒng)一平差的方法。

上述拖纜定位算法大致分為兩類：第一類是對(duì)拖纜建立特定的數(shù)學(xué)模型，并利用拖纜定位觀測(cè)值對(duì)數(shù)學(xué)模型中的參數(shù)進(jìn)行估計(jì)，拖纜的纜形和位置由模型直接導(dǎo)出。此類方法缺點(diǎn)是大部分拖纜模型只針對(duì)單纜定位，在多纜勘探場(chǎng)景下不能較好地利用距離觀測(cè)值。第二類是先利用部分觀測(cè)值進(jìn)行近似曲線積分推算，得到一個(gè)近似的拖纜纜形，再通過聲學(xué)測(cè)距網(wǎng)絡(luò)平差的結(jié)果對(duì)拖纜進(jìn)行調(diào)整。此類方法應(yīng)用較為廣泛，缺點(diǎn)是建立的拖纜模型不夠嚴(yán)密，不能較好地利用羅盤觀測(cè)值，且測(cè)距網(wǎng)絡(luò)平差時(shí)通常需估計(jì)上百個(gè)聲學(xué)節(jié)點(diǎn)的坐標(biāo)，計(jì)算量較大。

針對(duì)以上存在的問題，本文建立了一種基于曲線積分的拖纜數(shù)學(xué)模型，通過模型參數(shù)描述拖纜的形狀和位置，并基于模型充分利用拖纜定位網(wǎng)絡(luò)中的觀測(cè)值進(jìn)行參數(shù)估計(jì)，從而實(shí)現(xiàn)多纜定位和精度評(píng)估，最后通過仿真與實(shí)測(cè)數(shù)據(jù)對(duì)算法進(jìn)行驗(yàn)證。

1 基于曲線積分的拖纜數(shù)學(xué)模型

在進(jìn)行海上地震勘探時(shí)，拖纜在水下呈現(xiàn)“光滑的曲線”形態(tài)，其上安裝有深度控制傳感器，用于控制拖纜的沉放深度，實(shí)際作業(yè)時(shí)拖纜各處檢波器的深度差一般在1m之內(nèi)。由于拖纜長(zhǎng)度通常為3～10km，檢波器深度差可忽略不計(jì)。以單條拖纜為例，其數(shù)學(xué)描述如圖1所示，其中α(s)為拖纜任意位置s處切線方位角的函數(shù)。

圖1 拖纜數(shù)學(xué)模型示意圖

設(shè)拖纜位置s0處坐標(biāo)為(x0,y0)，以該點(diǎn)為基準(zhǔn)點(diǎn)，則拖纜上任意點(diǎn)s的平面坐標(biāo)(x,y)可表示為

(1)

如圖1所示，拖纜上任意點(diǎn)的坐標(biāo)(x,y)與其對(duì)應(yīng)位置s之間存在如下微分關(guān)系

(2)

將式(2)代入式(1)可得拖纜上任意位置s處的坐標(biāo)

(3)

通過該模型，若能確定(x0,y0)和α(s)，即可通過曲線積分導(dǎo)出拖纜上任意里程處點(diǎn)的坐標(biāo)，確定整條拖纜的纜形。式(3)中x0、y0可看作拖纜的位置參數(shù)，函數(shù)α(s)中的參數(shù)可看作拖纜的形態(tài)參數(shù)。

2 基于曲線積分模型的多纜定位算法

2.1 拖纜導(dǎo)航定位設(shè)備及定位模型

海上地震拖纜勘探中使用的導(dǎo)航定位設(shè)備(圖2)有：①勘探船上安裝的差分全球定位系統(tǒng)(Differential Global Positioning System，DGPS)和電羅經(jīng)。其中，DGPS用于獲取船的絕對(duì)位置，電羅經(jīng)用于獲取船的姿態(tài)。②槍陣與拖纜尾部掛載了GPS相對(duì)定位系統(tǒng)(Relative Global Positioning System，RGPS)浮標(biāo)，用于測(cè)量其與拖船上RGPS參考站之間的基線向量，確定它們相對(duì)于拖纜的位置(有時(shí)拖纜頭部也掛載RGPS)。③水下拖纜安裝有聲學(xué)設(shè)備和羅盤。其中，聲學(xué)設(shè)備用于測(cè)定聲學(xué)節(jié)點(diǎn)之間的距離，拖纜上的聲學(xué)設(shè)備與槍陣、尾標(biāo)上的聲學(xué)設(shè)備一起組成聲學(xué)測(cè)距網(wǎng)絡(luò)，羅盤用于確定拖纜在各羅盤節(jié)點(diǎn)處的切線方位[18-22]。

圖2 拖纜導(dǎo)航定位系統(tǒng)示意圖

海上地震勘探作業(yè)時(shí)，DGPS和RGPS的觀測(cè)值精度較高，RGPS可將拖船的位置基準(zhǔn)傳遞到槍陣與拖纜尾標(biāo)，位于槍陣與拖纜尾標(biāo)上的聲學(xué)節(jié)點(diǎn)因此可獲得較高的位置精度，可作為拖纜定位網(wǎng)絡(luò)中的已知點(diǎn)。拖纜的α(s)是未知的，受勘探船修舵、海浪以及洋流的影響，拖纜在水下的纜形隨之變化，因此，可以使用一個(gè)多項(xiàng)式對(duì)α(s)進(jìn)行描述

α(s)=a0+a1s+a2s2+…+ansn

(4)

式中：n為多項(xiàng)式階數(shù)；α0，α1，…，αn為多項(xiàng)式系數(shù)，即待估計(jì)的拖纜形態(tài)參數(shù)。

多項(xiàng)式曲線具有良好的幾何特性，可以較準(zhǔn)確地?cái)M合出不同的拖纜形狀。拖纜基準(zhǔn)點(diǎn)的選取對(duì)于最終的拖纜定位結(jié)果有直接影響，且基準(zhǔn)點(diǎn)最好位于拖纜上。盡管拖纜尾標(biāo)具有較高的位置精度，但由于其漂浮于水面之上，和拖纜之間有一定的高差，將其作為拖纜基準(zhǔn)點(diǎn)進(jìn)行積分計(jì)算會(huì)帶來較大的系統(tǒng)誤差。因此，可選取拖纜尾端的聲學(xué)節(jié)點(diǎn)作為拖纜基準(zhǔn)點(diǎn)，它與尾標(biāo)之間存在聲學(xué)距離觀測(cè)值，將其坐標(biāo)(x0,y0)作為待估計(jì)的拖纜位置參數(shù)，與拖纜形態(tài)參數(shù)一起參與平差，進(jìn)行整體估計(jì)。這樣既可使拖纜基準(zhǔn)點(diǎn)獲得較高的位置精度，又可避免積分過程中產(chǎn)生系統(tǒng)誤差。

2.2 誤差方程的建立

拖纜上的定位傳感器羅盤和聲學(xué)測(cè)距安裝位置已知。羅盤觀測(cè)值和聲學(xué)距離觀測(cè)值的誤差方程推導(dǎo)如下。

2.2.1 羅盤觀測(cè)值誤差方程

羅盤的觀測(cè)值為拖纜該節(jié)點(diǎn)處切線的磁方位角，拖纜上位于s處的羅盤方位觀測(cè)值觀測(cè)方程為

As=α(s)+m+εa

(5)

式中：As為s處的羅盤觀測(cè)值；m為磁偏角改正量；εa為方位觀測(cè)隨機(jī)誤差。該觀測(cè)方程為線性，結(jié)合式(4)、式(5)可寫為誤差方程形式

(6)

2.2.2 聲學(xué)距離觀測(cè)值誤差方程

聲學(xué)設(shè)備的距離觀測(cè)值分為同纜觀測(cè)值和異纜觀測(cè)值。其中，同纜距離觀測(cè)值兩端對(duì)應(yīng)的聲學(xué)節(jié)點(diǎn)位于同一條拖纜，異纜觀測(cè)值的聲學(xué)節(jié)點(diǎn)則位于不同拖纜上。為不失一般性，先對(duì)異纜觀測(cè)情況下的誤差方程進(jìn)行推導(dǎo)，同纜觀測(cè)可看作異纜觀測(cè)的一種特殊情況。設(shè)G纜與H纜為相鄰的兩條拖纜，i點(diǎn)位于G纜里程SG,i處，j點(diǎn)位于H纜里程SH,j處，i、j兩點(diǎn)間的距離觀測(cè)值為

(7)

式中：xG,i與yG,i、xH,j與yH,j分別為i、j兩點(diǎn)的平面坐標(biāo)；εd為距離觀測(cè)值隨機(jī)誤差。對(duì)式(7)進(jìn)行線性化，得到

(8)

(9)

將式(9)代入式(8)，得

(10)

式(10)即為距離觀測(cè)誤差方程的一般形式。對(duì)于同纜觀測(cè)的情況，式(10)可簡(jiǎn)化為

(11)

2.3 拖纜參數(shù)估計(jì)

據(jù)式(6)、式(10)可確定設(shè)計(jì)矩陣A和閉合差向量L。隨機(jī)模型P由觀測(cè)值先驗(yàn)信息確定，目前一般采用經(jīng)驗(yàn)定權(quán)的方法。根據(jù)最小二乘估計(jì)確定待估參數(shù)

(12)

(13)

3 實(shí)驗(yàn)驗(yàn)證與分析

為驗(yàn)證本文算法的應(yīng)用效果，設(shè)計(jì)了如下的仿真實(shí)驗(yàn)與實(shí)測(cè)數(shù)據(jù)實(shí)驗(yàn)。

3.1 仿真實(shí)驗(yàn)

分別選取正弦函數(shù)和對(duì)數(shù)函數(shù)作為拖纜方位角函數(shù)的真值，并據(jù)此生成了兩套不同纜形的拖纜仿真數(shù)據(jù)(圖3)。每套仿真數(shù)據(jù)均為6纜配置，纜間距為100m，每條拖纜長(zhǎng)為7km，且槍陣與拖纜尾標(biāo)的位置已知。測(cè)線走向?yàn)檎?。拖纜上包含有羅盤觀測(cè)值和聲學(xué)距離觀測(cè)值，同一拖纜上的羅盤間距為300m，聲學(xué)測(cè)距網(wǎng)絡(luò)為前—后網(wǎng)絡(luò)配置。其中，仿真數(shù)據(jù)1(圖3a)模擬拖纜多段彎曲的情況，仿真數(shù)據(jù)2(圖3b)模擬拖纜整體彎曲的情況。

圖3 仿真數(shù)據(jù)拖纜配置示意圖

仿真數(shù)據(jù)1和數(shù)據(jù)2采用的方位角函數(shù)分別為

(14)

(15)

方案1：采用本文基于曲線積分的多纜定位算法對(duì)拖纜參數(shù)進(jìn)行估計(jì)，并導(dǎo)出拖纜的纜形和位置信息。其中，拖纜方位角函數(shù)選用7階多項(xiàng)式，這里的7階為經(jīng)驗(yàn)值。

方案2：采用傳統(tǒng)的測(cè)邊網(wǎng)平差算法進(jìn)行計(jì)算，即先利用羅盤觀測(cè)值對(duì)拖纜纜形進(jìn)行近似曲線積分推算，再使用測(cè)邊網(wǎng)平差的結(jié)果對(duì)拖纜纜形進(jìn)行調(diào)整。

將仿真數(shù)據(jù)1和仿真數(shù)據(jù)2的解算結(jié)果與實(shí)際纜形進(jìn)行對(duì)比，如圖4、圖5所示。由圖可見，使用本文算法得到的拖纜纜形在整體上與實(shí)際纜形更加吻合，而使用測(cè)邊網(wǎng)平差算法得到的拖纜在空間位置上存在更大的偏差，且其纜形多處不平滑。這是由于測(cè)邊網(wǎng)平差算法是對(duì)拖纜上離散的聲學(xué)節(jié)點(diǎn)進(jìn)行較為精確的位置估計(jì)，位于聲學(xué)節(jié)點(diǎn)之間的拖纜纜形只能通過羅盤觀測(cè)值進(jìn)行確定，而羅盤觀測(cè)值沒有參與平差，加入至其中的隨機(jī)誤差未能被較好地消除或減弱。其次，測(cè)邊網(wǎng)平差算法中的拖纜需要通過分段調(diào)整以吻合平差后的聲學(xué)節(jié)點(diǎn)位置，導(dǎo)致其各段的連接處出現(xiàn)纜形過渡不平滑的現(xiàn)象，不符合拖纜實(shí)際物理特性。本文算法使用一系列的拖纜參數(shù)描述拖纜的纜形和位置，并且將羅盤觀測(cè)值和聲學(xué)距離觀測(cè)值共同納入平差計(jì)算，直接對(duì)拖纜參數(shù)進(jìn)行估計(jì)，再通過曲線積分導(dǎo)出拖纜任意里程處的坐標(biāo)，更好地利用了拖纜上的不同類型的定位觀測(cè)值，估計(jì)得到的纜形更加平滑，精度更高。

圖4 仿真數(shù)據(jù)1兩種方案拖纜纜形對(duì)比

圖5 仿真數(shù)據(jù)2兩種方案拖纜纜形對(duì)比

計(jì)算各拖纜位于任意里程處的位置偏差，并將其位置偏差投影到沿測(cè)線方向(Direction Along Line，DA)和垂直測(cè)線方向(Direction Cross Line，DC)，這是因?yàn)榭碧阶鳂I(yè)中更關(guān)心拖纜上檢波器陣列在這兩個(gè)方向上的位置偏差。由于拖纜數(shù)量較多，這里僅給出各拖纜在DA與DC方向上的平均偏差DA與DC對(duì)比。如圖6、圖7所示，由圖可見，本文算法的拖纜定位精度整體上明顯高于測(cè)邊網(wǎng)平差算法。由于仿真數(shù)據(jù)均采用前—后聲學(xué)網(wǎng)絡(luò)配置(圖3)，拖纜前中部缺少聲學(xué)距離觀測(cè)，測(cè)邊網(wǎng)平差算法的定位誤差在拖纜中部迅速增加，而本文算法在拖纜各處的定位精度均較為穩(wěn)定。

圖6 仿真數(shù)據(jù)1兩種方案拖纜位置偏差對(duì)比

圖7 仿真數(shù)據(jù)2兩種方案拖纜位置偏差對(duì)比

表1給出了兩種方案計(jì)算結(jié)果與實(shí)際纜形之間的平均點(diǎn)位偏差、DA方向平均偏差以及DC方向平均偏差。結(jié)果表明，方案1的定位精度顯著優(yōu)于方案2。其中，方案1的DA方向平均偏差均小于1m，DC方向平均偏差均少于3m。針對(duì)仿真數(shù)據(jù)1，方案1相較于方案2，點(diǎn)位偏差減小了63.8%，DA方向偏差減小了79.5%，DC方向偏差減小了61.1%。針對(duì)仿真數(shù)據(jù)2，方案1相較于方案2點(diǎn)位偏差減小了84.1%，DA方向偏差減小了77.4%，DC方向偏差減小了84.5%。

表1 不同方案偏差結(jié)果統(tǒng)計(jì)

3.2 實(shí)測(cè)數(shù)據(jù)實(shí)驗(yàn)

實(shí)測(cè)數(shù)據(jù)來自于2019年中國南海某測(cè)線拖纜定位數(shù)據(jù)，該測(cè)線為東南走向，方位角為149°，全長(zhǎng)38km。如圖8所示，該次作業(yè)為6纜配置，拖纜長(zhǎng)約7km，纜間距約為100m，每條拖纜掛載檢波器564個(gè)，間距為12.5m；掛載羅盤25個(gè)，間距約300m；掛載聲學(xué)設(shè)備15個(gè)，采用前—后聲學(xué)網(wǎng)絡(luò)配置。該套數(shù)據(jù)包含1600個(gè)響炮歷元，響炮間隔為12～13s。使用本文算法對(duì)該套數(shù)據(jù)進(jìn)行逐個(gè)歷元解算，拖纜方位角函數(shù)根據(jù)經(jīng)驗(yàn)選擇為7階多項(xiàng)式，并將解算結(jié)果與國外拖纜勘探導(dǎo)航定位數(shù)據(jù)處理軟件SPRINT的解算結(jié)果進(jìn)行對(duì)比分析。

圖8 實(shí)測(cè)數(shù)據(jù)拖纜配置示意圖

將SPRINT軟件的解算結(jié)果作為參考值，計(jì)算本文算法得到的拖纜檢波器坐標(biāo)與SPRINT解算結(jié)果之間的差值，將坐標(biāo)偏差投影到DA、DC方向，并統(tǒng)計(jì)每個(gè)歷元6條電纜上共3384個(gè)檢波器的平均DA、DC偏差。由表2和圖9可以看出，本文算法的結(jié)果與SPRINT總體上非常接近，且各歷元的檢波器平均DA、DC偏差一直維持在一個(gè)較小的水平。這表明本文算法在應(yīng)用于實(shí)測(cè)數(shù)據(jù)時(shí)可以獲得較為正確、可靠的結(jié)果，且算法具有良好的穩(wěn)定性。與SPRINT的解算結(jié)果相比，本文算法1600個(gè)響炮的歷元平均DA偏差為0.79m，其中，99.6%的歷元DA偏差小于2.00m，80.6%的歷元DA偏差小于1m，標(biāo)準(zhǔn)差小于0.50m；平均DC偏差為1.73m，其中，99.8%的歷元DC偏差小于3.00m，85.4%的歷元DC偏差小于2.00m，標(biāo)準(zhǔn)差小于0.50m。

圖9 檢波器平均位置偏差隨時(shí)間變化

表2 檢波器位置偏差結(jié)果統(tǒng)計(jì) m

計(jì)算各拖纜上處于相同里程處檢波器的平均DA、DC偏差，并統(tǒng)計(jì)1600個(gè)響炮歷元的平均值，如圖10所示。可以看出，檢波器平均DA偏差沿拖纜里程變化平緩，拖纜上99.1%的檢波器平均DA偏差優(yōu)于1.00m；檢波器平均DC偏差在拖纜前中部與拖纜尾部都出現(xiàn)了峰值，拖纜上97.3%的檢波器平均DC偏差優(yōu)于3.00m，70.1%的檢波器平均DC偏差優(yōu)于2.00m。由于該條測(cè)線作業(yè)時(shí)采用了前—中后聲學(xué)網(wǎng)絡(luò)配置(圖8)，拖纜前中部存在約3.7km長(zhǎng)的部分沒有聲學(xué)距離觀測(cè)，導(dǎo)致檢波器DC偏差在拖纜前中部會(huì)出現(xiàn)一個(gè)峰值，仿真實(shí)驗(yàn)的結(jié)果也驗(yàn)證了這一點(diǎn)。而在拖纜尾部出現(xiàn)的DC偏差的快速增長(zhǎng)，可能是由于SPRINT軟件在定位解算時(shí)對(duì)拖纜尾部采取了某些不同的處理。

圖10 檢波器平均位置偏差隨拖纜里程變化

計(jì)算羅盤方位觀測(cè)值與聲學(xué)距離觀測(cè)值各個(gè)歷元的殘差RMS值并進(jìn)行統(tǒng)計(jì)，結(jié)果如表3所示。其中，羅盤方位觀測(cè)值殘差RMS小于0.6°，聲學(xué)距離觀測(cè)值殘差RMS小于1.3m，說明本文算法的建模效果較好，所建立的拖纜模型內(nèi)符精度較高。

表3 定位觀測(cè)值殘差RMS統(tǒng)計(jì)

4 結(jié)論

高精度拖纜定位算法是海上地震勘探作業(yè)的重要保障，是后續(xù)地震波數(shù)據(jù)高精度成像的關(guān)鍵，且該技術(shù)長(zhǎng)期以來受到國外的封鎖。為解決現(xiàn)有定位算法數(shù)學(xué)模型不嚴(yán)密、定位觀測(cè)值利用不充分的問題，本文針對(duì)水下拖纜光滑的物理特性，建立了一種基于嚴(yán)密曲線積分的拖纜數(shù)學(xué)模型，使用一系列位置參數(shù)與形態(tài)參數(shù)描述拖纜，并基于該模型提出了一種多纜定位嚴(yán)密平差算法，該算法能充分利用不同類型的定位觀測(cè)值估計(jì)拖纜參數(shù)。

仿真實(shí)驗(yàn)表明了本文算法顯著優(yōu)于傳統(tǒng)的測(cè)邊網(wǎng)定位算法，驗(yàn)證了新算法在理論上的可行性、優(yōu)越性。實(shí)測(cè)數(shù)據(jù)實(shí)驗(yàn)表明了本文算法的定位結(jié)果與國外成熟商業(yè)軟件SPRINT的定位結(jié)果較為接近，驗(yàn)證了新算法在實(shí)際應(yīng)用中的有效性、可靠性。綜合各實(shí)驗(yàn)數(shù)值結(jié)果，本文提出的基于曲線積分的多纜定位算法在測(cè)線和垂直測(cè)線方向上的定位精度分別小于1m和3m，滿足現(xiàn)代海上地震勘探對(duì)于拖纜定位精度的需求。

需要指出的是，文中拖纜方位角多項(xiàng)式階數(shù)是一個(gè)經(jīng)驗(yàn)值，并未考慮海洋流態(tài)環(huán)境的影響。階數(shù)是否合理以及如何顧及海洋流態(tài)環(huán)境優(yōu)化階數(shù)，還需結(jié)合實(shí)際數(shù)據(jù)開展進(jìn)一步的研究與驗(yàn)證。