陸域重力勘探進(jìn)展

楊亞斌, 荊 磊, 徐夢龍, 韓革命, 邱隆君, 吳新剛,郜曉亮, 郝國江, 孫誠業(yè), 張光之, 蘇振寧

(1.中國地質(zhì)科學(xué)院地球物理地球化學(xué)勘查研究所,廊坊 065000 2.國家現(xiàn)代地質(zhì)勘查工程技術(shù)研究中心,廊坊 065000 3.海洋油氣勘探國家工程研究中心,北京 100028)

0 引言

地球重力場是地球系統(tǒng)物質(zhì)屬性產(chǎn)生的一個(gè)最基本的物理場，其靜態(tài)和動(dòng)態(tài)(時(shí)變)特征反映了地球內(nèi)部和表層物質(zhì)的分布變化和運(yùn)動(dòng)狀態(tài)，重力場觀測數(shù)據(jù)為地球內(nèi)部結(jié)構(gòu)和地球動(dòng)力學(xué)等研究提供基礎(chǔ)觀測信息[1-4]，它是重力勘探的理論基礎(chǔ)[5]。重力勘探是測量與圍巖有密度差異的地質(zhì)體在其周圍引起的重力異常，以確定這些地質(zhì)體存在的空間位置大小和形狀從而對工作地區(qū)的地質(zhì)構(gòu)造做出判斷的一種地球物理勘探方法。

自從十九世紀(jì)三十年代重力被用于地質(zhì)勘探以來，重力勘探經(jīng)過幾十年的發(fā)展，已經(jīng)取得巨大成就。我國1936年首次利用重力扭秤探測鉛鋅礦，之后的很長一段時(shí)間，重力勘探的主要應(yīng)用領(lǐng)域是油氣和固體礦產(chǎn)勘查[6-10]。1967年物化探所完成我國第一個(gè)物探應(yīng)用程序“礦區(qū)規(guī)則網(wǎng)重力測量地形改正”的手編程序[11]。1979年區(qū)域重力調(diào)查列入中國地質(zhì)工作計(jì)劃，開始系統(tǒng)開展區(qū)域重力調(diào)查工作，不斷完善調(diào)查方法與提高我國陸域重力調(diào)查工作程度[12-13]。近些年，隨著重力勘探儀器、觀測方法、測量方式、數(shù)據(jù)處理技術(shù)及解釋方法的發(fā)展，重力勘探更是邁向了一個(gè)新的階段。筆者以重力勘探儀器設(shè)備、觀測方法、數(shù)據(jù)處理與解釋技術(shù)為出發(fā)點(diǎn)，對重力勘探技術(shù)在國內(nèi)地學(xué)的應(yīng)用研究現(xiàn)狀進(jìn)行總結(jié)和展望。

1 重力儀現(xiàn)狀

陸域重力儀包括有相對重力儀和絕對重力儀。

1.1 相對重力儀

最早的重力測量儀器以擺儀和扭秤為主,1890年出現(xiàn)了最早的測量重力變化率的扭秤，哈特利(Hartley)[13]、格拉夫(A.Graf)[13]以垂直彈簧秤為基礎(chǔ)研制出了杠桿彈簧扭秤系統(tǒng)，由此誕生了德國阿斯卡利亞(Askania)GS 型重力儀。在1939年制造出第一臺可以工作的拉科斯特(LaCoste)金屬彈簧重力儀，1947年美國的沃登(Sam.P.Worden)[13]研制成功了自動(dòng)石英彈簧重力儀，此后其他制造商以同樣的原理生產(chǎn)了類似的石英彈簧重力儀，如加拿大 Scintrex 公司的CG型和蘇聯(lián)的GAK型重力儀。目前，Scintrex公司的CG-6型全自動(dòng)石英彈簧重力儀測量精度優(yōu)于0.005×10-5m/s2，成為世界主流的相對重力儀[5,13-14]。

在國內(nèi)，地面相對重力儀主要是由中國地質(zhì)裝備集團(tuán)公司北京奧地探測儀器公司研發(fā)和生產(chǎn)的ZSM系列重力儀，采用機(jī)械調(diào)零式石英彈簧傳感器。1975年，北京地質(zhì)儀器廠生產(chǎn)的ZSM-3型重力儀投入批量生產(chǎn)。目前最新的重力儀已經(jīng)實(shí)現(xiàn)0.001×10-5m/s2讀數(shù)分辨率，重復(fù)性達(dá)到0.005×10-5m/s2，殘余長期漂移<0.1×10-5m/s2，部分指標(biāo)已經(jīng)接近CG-5型重力儀水平。1994年物化探所研制了淺海重力測量系統(tǒng)，填補(bǔ)了我國在此方面的空白[11]。

1.2 絕對重力儀

1946年伏萊(C.Colet)[23]提出應(yīng)用自由落體原理測定絕對重力值，并在1952年實(shí)現(xiàn)了±10.0×10-5m/s2精度絕對重力儀。1980年代后，法拉(J.E.Faller)[24]在之前基礎(chǔ)上制成了JILA型輕便絕對重力儀，此后不斷改進(jìn)，于1993年研制出新一代商業(yè)化可移動(dòng)式FG5型絕對重力儀，精度可達(dá)±0.001×10-5m/s2～0.002×10-5m/s2，目前該產(chǎn)品已升級到FG5-X。

中國計(jì)量科學(xué)研究院從1960年代開始，先后研制了NIM型系列絕對重力測量儀器[15-17]；清華大學(xué)研制了T-2等號型絕對重力儀[18]等。此外，中國地震局研制了新型激光干涉絕對重力儀，上述儀器精度可達(dá)0.005×10-8m/s2～0.010×10-8m/s2左右[19-20]。2002年開始，浙江工業(yè)大學(xué)、浙江大學(xué)、華中科技大學(xué)等單位研發(fā)基于激光干涉法的自由落體原理，原子絕對重力儀。近年來，絕對重力儀在小型化、儀器的測量精度和測量靈敏度等性能得以較快提升[21-22]。

2 方法技術(shù)現(xiàn)狀

2.1 測地工作

重力勘探的前身是研究地球形狀的重力測量學(xué)，它與測量學(xué)緊密相連，測繪技術(shù)的發(fā)展也促進(jìn)著重力勘探技術(shù)進(jìn)步。

2.1.1 坐標(biāo)系統(tǒng)與高程系統(tǒng)

重力勘探工作必須用到坐標(biāo)和高程資料，而大地坐標(biāo)系統(tǒng)和高程系統(tǒng)是一個(gè)國家最基本、最重要的兩個(gè)大地測量基準(zhǔn)。我國重力勘探初期，采用了1954年北京坐標(biāo)系[25]和1956年黃海高程系統(tǒng)[29]，1990年又啟用了1980西安坐標(biāo)系[25]和1985國家高程基準(zhǔn)[29]。從2008年7月起，中國啟用了新一代國家大地坐標(biāo)系-2000國家大地坐標(biāo)系[25](China Geodetic Coordinate System 2000， CGCS 2000)。

重力勘探工作所采用的坐標(biāo)系統(tǒng)與高程系統(tǒng)也在不斷更新，最新的行業(yè)標(biāo)準(zhǔn)《區(qū)域重力調(diào)查規(guī)范》(DZ/T 0082-2021)[27]也提出了新的區(qū)域重力調(diào)查技術(shù)要求：統(tǒng)一采用2000國家重力基本網(wǎng)系統(tǒng)；統(tǒng)一采用2000國家大地坐標(biāo)系和大地高程基準(zhǔn)；統(tǒng)一采用CGCS2000橢球面正常重力公式計(jì)算正常重力值。

2.1.2 點(diǎn)位測量技術(shù)

重力勘探工作中測地工作的主要任務(wù)是布設(shè)測網(wǎng),確定測點(diǎn)和各級重力基點(diǎn)的平面位置和高程。其工作方法，已由地形圖定點(diǎn)，水準(zhǔn)、氣壓測高和航片定點(diǎn)測高等傳統(tǒng)方法技術(shù)發(fā)展為全球?qū)Ш叫l(wèi)星系統(tǒng)(GNSS)進(jìn)行三維坐標(biāo)測量。

通常區(qū)域重力調(diào)查測點(diǎn)平面位置和高程使用GNSS接收機(jī)測定，采用包括快速靜態(tài)測量、實(shí)時(shí)動(dòng)態(tài)測量(RTK)以及其他滿足精度要求的測量方法[27-28]?？刂泣c(diǎn)的聯(lián)測方法為載波相對定位測量(靜態(tài)定位)，測量控制點(diǎn)應(yīng)由已知GNSS控制點(diǎn)引出，布網(wǎng)方式根據(jù)需要可采用三角網(wǎng)、單三角形、附合導(dǎo)線網(wǎng)等，測網(wǎng)至少應(yīng)該包含2個(gè)已知GNSS控制點(diǎn)。

2.2 重力測量工作

2.2.1 重力基準(zhǔn)網(wǎng)現(xiàn)狀

全球重力基準(zhǔn)是將全球相對重力測量納入統(tǒng)一系統(tǒng)的基礎(chǔ)。其發(fā)展經(jīng)歷了三個(gè)主要的系統(tǒng)：維也納重力系統(tǒng)、波茨坦重力系統(tǒng)和國際重力標(biāo)準(zhǔn)網(wǎng)1971[3]。

1956年在前蘇聯(lián)國家重力測量隊(duì)的協(xié)助下，我國建立了第一個(gè)國家重力基準(zhǔn)網(wǎng)，由21個(gè)基本重力點(diǎn)和82個(gè)一等重力點(diǎn)組成，基本點(diǎn)精度優(yōu)于±0.32×10-5m/s2，一等點(diǎn)精度優(yōu)于±0.40×10-5m/s2，稱為1957國家重力基準(zhǔn)[26]。1985年我國完全自主地重新建立了國家重力基準(zhǔn)，由6個(gè)基準(zhǔn)重力點(diǎn)、46個(gè)基本重力點(diǎn)和5個(gè)引點(diǎn)組成，精度為±0.015×10-5m/s2，稱為1985國家重力基本網(wǎng)[26]。2000年左右我國又完成了全國分布的295個(gè)重力基本點(diǎn)的高精度測量和數(shù)據(jù)處理，建立精度為±0.0074×10-5m/s2的2000國家重力基本網(wǎng)[26]。2017年～2019年建立了新一代重力基準(zhǔn)網(wǎng)(尚未正式發(fā)布)。

2.2.2 重力測量

野外重力數(shù)據(jù)測量為基點(diǎn)-測點(diǎn)-基點(diǎn)的單程閉合方法?；c(diǎn)網(wǎng)采用雙程往返(精度較高儀器)或三程小循環(huán)法(精度較低儀器)[27]。

2.2.3 地形改正

在重力勘探數(shù)據(jù)精度評價(jià)中，地形起伏的影響是個(gè)不可忽略的問題[30]。

1)近區(qū)地形改正(0 m～50 m)。近區(qū)地形改正值的獲取，已從目估法發(fā)展到儀器測量法。用于近區(qū)地形改正的儀器主要有森林羅盤儀[31]、激光測距儀[32-35]和高精度GNSS[31,34]。目前，已有學(xué)者進(jìn)行無人機(jī)激光掃描獲取地形數(shù)據(jù)，利用棱柱公式進(jìn)行近區(qū)地形改正試驗(yàn)，另外，還有一些學(xué)者利用大比例尺高精度DEM數(shù)據(jù)進(jìn)行近區(qū)地形改正試驗(yàn)。

2)中區(qū)地形改正(50 m～2 000 m)。中區(qū)地形改正技術(shù)發(fā)展很快，以DEM高程數(shù)據(jù)基礎(chǔ)，利用共用點(diǎn)法、表面積分法和三棱錐體元剖分法計(jì)算地形改正值[36]。此外，對湖泊水體的改正也有研究試驗(yàn)[33]。

3)遠(yuǎn)區(qū)地形改正(2 000 m～166 700 m)。遠(yuǎn)一區(qū)(2 000 m～20 000 m)地形改正采用平面公式計(jì)算地形改正值，使用1 km×1 km高程節(jié)點(diǎn)網(wǎng)按狹義地形改正(相對測點(diǎn))和有限中間層(指20 km以內(nèi))模式計(jì)算；遠(yuǎn)二區(qū)(20 km～166.7 km)地形改正因包含中間層改正又稱為廣義地形改正,采用5′×5′平均高程用計(jì)算機(jī)完成[38-39]。

遠(yuǎn)區(qū)地形改正計(jì)算涉及大量的高程數(shù)據(jù)，方法技術(shù)進(jìn)展緩慢，影響重力勘探數(shù)據(jù)總精度進(jìn)一步提高，急需高程數(shù)據(jù)的更新和方法研究。

2.3 數(shù)據(jù)預(yù)處理

2.3.1 正常場改正公式

在重力勘探的工作中以往用過的正常場改正公式有1901～1909赫爾墨特公式、1930卡西尼公式、(IAG)推薦的1980公式。隨著2000 國家大地坐標(biāo)系(CGCS2000)的啟用，重力正常場改正公式也須更新[40-41]。最新的行業(yè)標(biāo)準(zhǔn)《區(qū)域重力調(diào)查規(guī)范》(DZ/T 0082-2021)也提出了新的采用2000橢球面正常重力公式計(jì)算正常重力值[27]。

2.3.2 布格改正

高度改正、中間層改正計(jì)算公式所用的高程由1985高程基準(zhǔn)改進(jìn)為CGCS2000橢球大地高，減去了正高和大地高之間物質(zhì)產(chǎn)生的重力值，使重力勘探工作使用的布格異常更加清晰、準(zhǔn)確地反映了地質(zhì)目標(biāo)意義[42-44]。

2.4 數(shù)據(jù)處理

野外獲取的重力數(shù)據(jù)經(jīng)過零點(diǎn)漂移校正、固體潮校正、緯度校正(正常場改正)、布格校正(高度改正和中間層改正)和地形校正等預(yù)處理，得到反映地下密度異常體分布的布格重力異常[5,27,45]。

由于布格重力異常是地下所有的密度不均勻體引起的疊加異常，重力數(shù)據(jù)處理的目的是從疊加異常中分離或突出目標(biāo)體引起的重力異常，使其產(chǎn)生的信息更易于識別和辨認(rèn)。重力數(shù)據(jù)處理和反演是后期地質(zhì)解釋的關(guān)鍵[46-48]。

2.4.1 異常分離

為了獲得勘探目標(biāo)地質(zhì)體引起的異常，通常需要對重力異常進(jìn)行異常分離。能否正確提取出目標(biāo)異常，決定著反演結(jié)果是否可靠。正是由于異常分離的重要性，多年來很多地球物理學(xué)者都致力于尋找更優(yōu)的異常分離方法[5,48,49]。

國內(nèi)、外的地球物理工作者[50-83]提出了多種區(qū)域和剩余異常分離方法，主要包括圖形法、 趨勢分析法、濾波法(向上延拓法、匹配濾波、小波變換)、統(tǒng)計(jì)法、剝皮法等方法。

圖形法是最早使用的方法，缺點(diǎn)是主觀性比較大，常用于剖面數(shù)據(jù)處理，對于大量的平面數(shù)據(jù)難以實(shí)現(xiàn)。

Agocs[50]提出了基于最小二乘擬合的趨勢分析方法,之后大量關(guān)于多項(xiàng)式擬合的方法及其改進(jìn)方法被提出，這些方法包括采用正交多項(xiàng)式擬合，有變階滑動(dòng)趨勢分析法、篩選-趨勢分析法[53-55]。趨勢分析法的結(jié)果好壞直接取決于所選取數(shù)學(xué)模型與實(shí)際區(qū)域異常的逼近程度。

濾波是一種很寬泛的概念，Dean[51]提出了重磁處理中的線性濾波理論，認(rèn)為二階導(dǎo)數(shù)、解析延拓、圓滑等可以被認(rèn)為是濾波。隨著計(jì)算機(jī)的廣泛應(yīng)用，快速傅立葉變換被引入到地球物理數(shù)據(jù)處理中，學(xué)者們設(shè)計(jì)了各種理想濾波器，推導(dǎo)了基于傅里葉變換和采樣理論的新的卷積濾波系數(shù)，波數(shù)域?yàn)V波方法逐漸成為主流[56-60]。常用的濾波主要有：

1)垂向二階導(dǎo)數(shù)法常被用來提取埋藏淺的局部異常。Evjen[52]利用垂向一階導(dǎo)數(shù)與垂向二階導(dǎo)數(shù)處理并解釋了重力資料。之后，眾多學(xué)者提出了關(guān)于垂向二階導(dǎo)數(shù)在空間域的近似計(jì)算公式[61-63]。目前，導(dǎo)數(shù)的計(jì)算在波數(shù)域中被廣泛的研究與應(yīng)用。

2)解析延拓法是將獲得的測區(qū)重力異常，利用相關(guān)的數(shù)學(xué)算法變換到更高或更低的觀測面，以實(shí)現(xiàn)區(qū)域異常與局部異常的分離。最初由Peters[59-64]提出，他認(rèn)為延拓可以對不同波長的信號進(jìn)行增強(qiáng)或衰減。之后，眾多學(xué)者對方法進(jìn)行完善與發(fā)展[65-70]。

3)匹配濾波，位場信號的各頻段成分不相關(guān)時(shí)，擬合對數(shù)功率譜并計(jì)算直線斜率和截距，確定匹配濾波因子，實(shí)現(xiàn)深源場和淺源場的分離。1970年Spector[71]等首先根據(jù)磁異常的功率譜設(shè)計(jì)了匹配濾波器來分離不同深度場源產(chǎn)生的異常。之后，學(xué)者們不斷完善、推進(jìn)該方法[72-78]。

4)利用小波多尺度分析方法，將重力異常分解到不同尺度空間中，尺度大小表示了異常所反映的地質(zhì)體規(guī)模和埋深的大小。湯井田[79]將小波變換引入重磁數(shù)據(jù)處理。Fedi[80]利用小波變換的多尺度分析，用不同分量的組合來代表區(qū)域場和局部場。隨著技術(shù)推進(jìn)，利用小波變換把重力異常進(jìn)行多重分解，在一定程度上將頻譜互相重疊的異常提取出來，獲得地殼三維密度結(jié)構(gòu)的更多細(xì)節(jié)[82-83]。

剝皮法是由Hammer[83]提出的，即在淺層的構(gòu)造已經(jīng)很清楚的前提下，通過正演獲得其異常，再從觀測異常中減去這一部分異常，從而獲得深層的信息。有大量的地質(zhì)資料或鉆孔、地震資料的控制，效果很好[84-86]。

2.4.2 邊界識別

邊界識別技術(shù)主要是對異常體的水平位置進(jìn)行定位。在研究地質(zhì)目標(biāo)體的橫向不均勻性特別是地質(zhì)目標(biāo)體邊緣位置時(shí)，重、磁位場有其獨(dú)特的優(yōu)勢。位場源邊界識別方法包括：梯度運(yùn)算類、數(shù)理統(tǒng)計(jì)類、其它類等，以前兩者發(fā)展及應(yīng)用最為廣泛[79]。

梯度運(yùn)算類是通過對位場進(jìn)行梯度計(jì)算(或組合)，放大高頻成份增強(qiáng)邊界信息的識別，該類方法發(fā)展快速、應(yīng)用廣泛[88]。Evjen[89]提出利用垂向?qū)?shù)識別地質(zhì)體邊界信息；Cordell[90]認(rèn)為梯度模的極大值對目標(biāo)體的邊界有較好的指示效應(yīng)；王萬銀等[91]利用n階垂向?qū)?shù)峰值及其比值對邊界信息的有效提??；之后，其他學(xué)者進(jìn)一步提高了該類方法的邊界識別效果[92-94]。Nabighina[95-96]認(rèn)為針對二維位場數(shù)據(jù)的解析信號與磁化方向無關(guān)；Nabighina與Roest[96-97]提出三維解析信號能提供有效的場源邊界信息，Miller和Singh等[98]將解析信號相位的概念引入了邊緣識別，提出傾斜角(Tilt Angle)法；隨后，學(xué)者們對傾斜角法展開研究[99-104]；Wijns等[105]首次利用解析信號振幅對總水平導(dǎo)數(shù)進(jìn)行歸一化，提出了θ圖法；周文納等[106]于幾何意義上對傾斜角法與Theta map法實(shí)現(xiàn)了統(tǒng)一；陳安國等[107]基于Theta圖法進(jìn)行了改進(jìn)，實(shí)現(xiàn)了深部和淺部邊界的均衡。

數(shù)理統(tǒng)計(jì)運(yùn)算方法，壓制噪聲干擾并且結(jié)果相對穩(wěn)定。楊高印等[108]提出小子域?yàn)V波方法進(jìn)行異常的分離。之后，其他學(xué)者進(jìn)行了該方法的改進(jìn)[109-111]；Cooper等[112]首次提出了歸一化標(biāo)準(zhǔn)差；張恒磊等[113]提出各向異性標(biāo)準(zhǔn)化方差法；王彥國等[114]提出歸一化差分法，在邊界定位及識別能力方面體現(xiàn)出優(yōu)勢性；徐夢龍等[115]提出各方向均方差相關(guān)系數(shù)法；Du等[116]提出基于滑動(dòng)窗口的垂向?qū)?shù)平均值與其標(biāo)準(zhǔn)差的相關(guān)系數(shù)法。

其他類方法包括多尺度信號分析、結(jié)果圖像的檢測、加強(qiáng)等，也在不斷發(fā)展。朱保健[117]利用小波模極大值進(jìn)行位場的邊界信息加強(qiáng)；楊文采、孫艷云、付麗華等[118-120]利用譜矩的位場幾何特征對弧形特征的有效刻畫；夏玲燕等[121]將非線性濾波技術(shù)與斜導(dǎo)數(shù)及其水平導(dǎo)數(shù)進(jìn)行有效結(jié)合，實(shí)現(xiàn)弱信息加強(qiáng)及線性構(gòu)造的提?。悔w希剛等[122]利用圖像細(xì)化技術(shù)對斷裂帶信息的窄化處理；肖鋒等[123]利用Canny邊界檢測算子對斜導(dǎo)數(shù)法結(jié)果圖中的邊界信息提??；陳濤[124]提出Harris 濾波器技術(shù)，在提取邊界的同時(shí)抑制噪聲；張雙喜等[125]結(jié)合二維經(jīng)驗(yàn)?zāi)B(tài)分解與小波模極大值，對位場數(shù)據(jù)多尺度分解及各尺度異常邊界檢測；張志厚等[126]提出了基于深度學(xué)習(xí)的位場邊界識別方法。

2.4.3 重力三維反演技術(shù)

重力反演計(jì)算就是將異常數(shù)據(jù)轉(zhuǎn)化為描述地球內(nèi)部場源密度變化及分布的數(shù)據(jù)。針對不同的地質(zhì)信息，重磁反演主要分為形態(tài)反演和物性反演。形態(tài)反演是在已知場源物性參數(shù)的基礎(chǔ)上，利用多面體來擬合地下地質(zhì)體的幾何形態(tài)分布，從而刻畫地下地質(zhì)體的界面。吳文鸝[127]在研究三角形多面體模型重、磁異常三維正演和反演技術(shù)的基礎(chǔ)上,吸取人機(jī)交互正演模擬的優(yōu)點(diǎn),實(shí)現(xiàn)了三角形多面體模型重、磁異常數(shù)據(jù)三維人機(jī)聯(lián)作反演。物性反演解釋技術(shù)是對地下空間地質(zhì)體的物性參數(shù)進(jìn)行定量計(jì)算。物性反演由于具有反演條件寬松、容易操作且能夠模擬復(fù)雜地質(zhì)體的優(yōu)點(diǎn)，成為近年來位場反演的主要發(fā)展方向[48，128-130]。

三維反演方法可以分為非線性化反演和線性化反演兩類，其中線性反演的迭代優(yōu)化方法有最小二乘法、牛頓法、正則化方法、最速下降法和共軛梯度法等。牛頓法收斂速度快但有Hessian矩陣及其逆矩陣的大規(guī)模計(jì)算問題。最速下降法計(jì)算量少但收斂速度慢，而共軛梯度法既解決了收斂速度慢的問題又避免了大型矩陣向量的計(jì)算和儲存問題，是一種求解大型最小化問題的反演算法[131-136]。這些線性迭代優(yōu)化算法收斂速度快，但容易陷入局部極小值、計(jì)算耗時(shí)等困境。隨著數(shù)學(xué)優(yōu)化方法的發(fā)展，多種全局隨機(jī)優(yōu)化的非線性反演算法被引入到地球物理反演中，并得到成功應(yīng)用[48,128,137-138]。非線性化反演方法包括遺傳算法、神經(jīng)網(wǎng)絡(luò)算法、模擬退火算法、最小二乘方法等。1983年Silva等[139]將遺傳算法應(yīng)用于磁法反演中；Montesinos等[140]采用遺傳算法進(jìn)行重力三維反演；姚長利等[128-129,159]改進(jìn)方法，提高了非線性反演計(jì)算速度；Eslam 等[141]采用神經(jīng)網(wǎng)絡(luò)算法進(jìn)行重力實(shí)測數(shù)據(jù)反演，且準(zhǔn)確獲得了異常的分布位置，該算法計(jì)算具有并行性，但收斂效率低，易存在局部極小化問題；馬國慶[142]利用模擬退火法完成重磁異常物性反演，該方法有效地避免了大型矩陣的計(jì)算，可以克服局部極小化問題，但搜索效率變低；Tarantola 等[143]利用阻尼最小二乘方法求解非線性反演問題，從而得到穩(wěn)定的解，但其計(jì)算速度慢。最小二乘方法既可以應(yīng)用于非線性反演問題的計(jì)算，還可以應(yīng)用于線性反演問題計(jì)算中。

同其他地球物理反演方法一樣，由于復(fù)雜的地質(zhì)情況，三維重力場數(shù)據(jù)信息不足和場源的等效性問題會(huì)致使重力反演研究存在兩個(gè)重要難題，即解的不穩(wěn)定性和多解性，解決此問題成為重磁反演研究的重點(diǎn)和難點(diǎn)。為解決此問題，很多學(xué)者開展了約束反演研究[83，144-146]。當(dāng)約束條件采用不同的加權(quán)因子以及不同的計(jì)算策略時(shí)會(huì)得到不同意義上的解，主要的約束方式有：①深度約束，用于抵消核函數(shù)隨深度自然衰減，消除其近地表權(quán)重過大而使反演密度分布不符合真實(shí)異常源的情況，進(jìn)而提高深度分辨率[131-132,147-148]；②聚焦約束,可以清晰地刻畫異常體的邊界特征[149-152]，方便后期的處理解釋；③物性邊界約束，為了獲得更合理的反演物性分布，在物性反演過程中需要加入地質(zhì)體物性上下限約束條件，強(qiáng)制將反演密度值限制在一定的范圍內(nèi)[145,150]；④結(jié)構(gòu)約束，可用于不同地球物理屬性之間的聯(lián)合約束，有交叉梯度約束[153-155]和總梯度約束(summative gradients)[156-158]等；⑤地質(zhì)統(tǒng)計(jì)約束，可以將區(qū)域地質(zhì)特征規(guī)律信息和地質(zhì)學(xué)家對地質(zhì)情況的認(rèn)識加入到模型中。

面對大范圍、大數(shù)據(jù)量的重力數(shù)據(jù)，反演效率低的問題非常突出。為此國內(nèi)外學(xué)者從不同角度進(jìn)行了研究，主要包括：①降維度處理方法，通過減少維度從而減少存儲空間和計(jì)算量，有隨機(jī)子域反演、小波壓縮、基于多項(xiàng)式的反演[159-161]；②利用對稱性的處理方法，幾何格架法利用重力場正演核矩陣的平移不變的性質(zhì)減少計(jì)算量[128]。在此基礎(chǔ)上荊磊利用快速傅里葉變換(FFT)進(jìn)一步實(shí)現(xiàn)了具有空間域計(jì)算精度和頻率域計(jì)算速度的快速算法[162]。

3 重力數(shù)據(jù)庫與圖件編制

3.1 重力數(shù)據(jù)庫

1981年物化探所研發(fā)了“區(qū)域重力資料整理處理成圖自動(dòng)化”軟件系統(tǒng)[11]。原地礦部物化探局組織研發(fā)我國第一代制基于M-160型計(jì)算機(jī)的重力數(shù)據(jù)庫，數(shù)據(jù)庫由高程子庫、重力數(shù)據(jù)子庫和程序子庫組成。高程子庫存放全國及全球的高程數(shù)據(jù)。可以對高程數(shù)據(jù)檢索、拼接、更新和進(jìn)行從一公里到全球的地改及均衡改正的計(jì)算。重力子庫存放重力測點(diǎn)的平面座標(biāo)、高程、測點(diǎn)原始重力值和近中區(qū)地改值等五個(gè)基本數(shù)據(jù)，存放基點(diǎn)的座標(biāo)和重力值、觀測時(shí)間及施測單位的編號。可以通過程序控制對重力數(shù)據(jù)進(jìn)行檢索、更新及增補(bǔ)。程序子庫包括庫管理程序、資料整理程序和數(shù)據(jù)處理程序三類[163]。1997年由原地礦部區(qū)域重力調(diào)查方法技術(shù)中心升級遷移到PC計(jì)算機(jī)上，形成windows系統(tǒng)下可視化管理數(shù)據(jù)庫[38]。2006年中國地質(zhì)調(diào)查局發(fā)展研究中心再次升級為RGIS2006版。具備對全國范圍從野外到室內(nèi)的重力數(shù)據(jù)整理的數(shù)據(jù)預(yù)處理、時(shí)間域及頻率域數(shù)據(jù)處理、數(shù)據(jù)庫可視化管理、GIS制圖、通用GIS數(shù)據(jù)轉(zhuǎn)換等的功能,適合于重力調(diào)查及重磁綜合解釋工作[39，164]。

3.2 重力圖件的編制

地球物理圖件對于區(qū)域地質(zhì)構(gòu)造研究和資源調(diào)查戰(zhàn)略部署等都發(fā)揮了重要作用，具有重要價(jià)值[165-167]。隨著重力工作的持續(xù)推進(jìn)，重力工作者在不同階段編制出不同比例尺的全國和區(qū)域性重力圖件。

1978年物化探所對所編1°×1°布格重力異常圖解釋認(rèn)為，大興安嶺-太行山-武陵山的北北東向梯級帶及青藏高原周緣梯級帶反映了地殼厚度變化，有可能是深部斷裂存在的一種標(biāo)志。一些學(xué)者對1:4 000 000布格重力異常圖進(jìn)行系統(tǒng)解釋，將全國重力場劃分為22個(gè)異常區(qū)和18個(gè)異常帶，為現(xiàn)代中國大地構(gòu)造劃分提供了地球物理依據(jù)[168-170]。左愚等編制1:12000 000重力系列圖[171]。1999年孫文珂等利用最新實(shí)測重力數(shù)據(jù)編制1:5 000 000中國及毗鄰海區(qū)重力系列圖，并進(jìn)行了說明。將布格重力異常圖分為5個(gè)區(qū)，分布對應(yīng)濱西太平洋構(gòu)造域的內(nèi)帶及其外帶、特提斯-喜馬拉雅構(gòu)造域、古亞洲構(gòu)造域和三大構(gòu)造復(fù)合地帶，每個(gè)區(qū)的異常幅值反映了地殼厚度及密度的變化[167]。李志等利用最新完成的重力資料編制了青藏高原及鄰區(qū)重力系列圖，并進(jìn)行了大地構(gòu)造分區(qū)、推斷斷裂、圈定巖體、盆地形態(tài)等研究[172]。張明華等利用不同來源且基準(zhǔn)未知的重力數(shù)據(jù)資料，采用背景一致性原則并進(jìn)行局部調(diào)平、拼接，編制了1∶5 000 000比例尺的中國海陸及鄰域布格重力異常圖，并依據(jù)重力場特征、斷裂帶及造山帶的相互關(guān)系，將重力場劃分為38個(gè)一級異常區(qū)和部分二級異常區(qū)[173]。

4 重力勘探標(biāo)準(zhǔn)與規(guī)范

重力規(guī)范一般提出重力調(diào)查任務(wù)、基本技術(shù)要求、技術(shù)設(shè)計(jì)、重力儀調(diào)節(jié)及性能測試、野外工作、資料整理、精度評價(jià)、基礎(chǔ)圖件編繪、數(shù)據(jù)處理與解釋、成果報(bào)告編寫與資料提交等方面的技術(shù)要求。所有規(guī)范都是經(jīng)過一段時(shí)間的工作實(shí)踐后，制訂出技術(shù)規(guī)定，最后升級為規(guī)范，并隨著相關(guān)技術(shù)進(jìn)步而不斷修訂完善。如1978年代有計(jì)劃地開展了區(qū)域重力調(diào)查工作后，首先制訂了《區(qū)域重力調(diào)查技術(shù)規(guī)定》[174]。1993年升級為行標(biāo)《區(qū)域重力調(diào)查規(guī)范》(DZ/T 0082-93)[175],2021升級為新版《區(qū)域重力調(diào)查規(guī)范》(DZ/T 0082-2021)[23]。

重力勘探規(guī)范由兩大行業(yè)規(guī)范構(gòu)成，包括地質(zhì)礦產(chǎn)行業(yè)規(guī)范和石油天然氣行業(yè)規(guī)范。

地質(zhì)礦產(chǎn)行業(yè)重力規(guī)范由自然資源部發(fā)布，按工作比例尺分為三種，現(xiàn)階段執(zhí)行的分別標(biāo)準(zhǔn)是《大比例尺重力勘查規(guī)范》(DZ/T 0171-2017)[176]、《重力調(diào)查技術(shù)規(guī)范(1∶50 000)》(DZ/T 0004-2015)[177]、《區(qū)域重力調(diào)查規(guī)范》(DZ/T 0082-2021)[23]。

石油天然氣行業(yè)規(guī)范由國家能源局發(fā)布，現(xiàn)階段執(zhí)行的是《陸上重力勘查技術(shù)規(guī)程》(SY/T 5819-2010)版[178]。

5 重力勘探技術(shù)在地質(zhì)工作中應(yīng)用

重力勘探方法具有輕便、快捷和經(jīng)濟(jì)的特點(diǎn)，能快速獲取大面積高精度重力數(shù)據(jù)，在地質(zhì)工作領(lǐng)域有著廣泛的應(yīng)用。

5.1 在能源勘探方面應(yīng)用

重力資料在確定沉積盆地邊界、劃分構(gòu)造單元、確定基底起伏和蓋層構(gòu)造、進(jìn)行油氣遠(yuǎn)景評價(jià)等方面發(fā)揮了重要作用。在大慶、勝利、大港、任丘、長慶、二連等油氣田的發(fā)現(xiàn)中，重力勘探都做出了貢獻(xiàn)[9]。隨著聯(lián)合反演、三維建模技術(shù)的不斷發(fā)展，重力勘探在能源資源勘探領(lǐng)域繼續(xù)發(fā)揮重要作用。1999年在華北地區(qū)綜合應(yīng)用重力與地震勘探發(fā)現(xiàn)了信安鎮(zhèn)北潛山，獲得了高產(chǎn)油氣流。2000年應(yīng)用重力、地震聯(lián)合反演解釋，在武清凹陷斜坡碼頭潛山發(fā)現(xiàn)了工業(yè)氣流。之后相繼發(fā)現(xiàn)了信安鎮(zhèn)潛山、碼頭潛山等一系列潛山，勘探取得了重大進(jìn)展[179-183]。

5.2 在礦產(chǎn)勘查方面應(yīng)用

重力勘探在探測礦產(chǎn)資源時(shí)，利用來自礦化作用、蝕變帶和特定主巖或構(gòu)造的位場響應(yīng)，可以查找與成礦有關(guān)的地層、構(gòu)造、巖漿巖、蝕變巖和礦化帶等控礦因素，研究金屬礦床賦存的巖體或構(gòu)造以推斷礦體的位置[184-189]。在直接找礦方面，重力勘探發(fā)現(xiàn)了加拿大魁北克西部以黃鐵礦和閃鋅礦為主的大型礦床[190]。國內(nèi)，在冀東鐵礦利用重力資料對隱伏鐵礦或含鐵建造進(jìn)行了系統(tǒng)的解釋和推斷，推斷9條與隱伏鐵礦床或含鐵建造有關(guān)的重力異常帶，部分已被鉆探工程所證實(shí)[191]。在間接找礦方面，重力勘探方法運(yùn)用更為廣泛。在美國加利福尼亞州礦區(qū)，利用地面重力和航磁數(shù)據(jù)，構(gòu)建了2.5維重力模型，圈定潛在的拆離斷層和礦化巖體分布[192]。

在加拿大拉布拉多，利用重力數(shù)據(jù)約束反演來預(yù)測和描繪的大型鎳、銅、鈷礦床的地下特征和展布，為礦體的開采提供了非常有用的地下信息[193]。重力勘探是我國20世紀(jì)60年代～70年代的鉻鐵礦和富鐵礦的普查中的主要方法。1990年代，高精度重力測量發(fā)現(xiàn)了新疆小熱泉子銅礦Ⅰ號礦床重力異常，并確定了孔位，鉆孔的驗(yàn)證取得了當(dāng)年找礦的重要突破[194]。在德興銅礦勘探中，利用重力資料研究與侵入巖體成礦相關(guān)礦種中提取與導(dǎo)礦相關(guān)的導(dǎo)礦巖體與構(gòu)造[195]。在黑龍江省下嘎來奧伊河礦勘查中，分析與矽卡巖型金屬礦有關(guān)的重、磁異常，總結(jié)找礦模型并圈定靶區(qū)，鉆探驗(yàn)證情況良好[196]。劉天佑等[197]在西藏朗縣秀溝工區(qū)完成了高精度重磁勘探，推斷出鉻鐵礦與礦化體的重磁遠(yuǎn)景異常，并且已經(jīng)得到證。潘彤等[198]根據(jù)柴達(dá)木盆地南緣和北緣不同物性分布及成礦地質(zhì)條件，優(yōu)選重磁、電法、地震等有效的方法進(jìn)行探測，總結(jié)分析其找礦模型，取得了一定的找礦效果。

5.3 在地質(zhì)調(diào)查方面應(yīng)用

5.3.1 工作程度

1950年代，由于石油勘探的需要我國大規(guī)模開展重力工作。經(jīng)過70多年的工作，截至2021年底，全國完成1∶1 000 000區(qū)域重力調(diào)查833×104km2,基本實(shí)現(xiàn)我國陸域1∶1 000 000區(qū)域重力調(diào)查的覆蓋；1∶250 000區(qū)域重力調(diào)查完成587.5×104km2，1∶250 000區(qū)域重力調(diào)查完成全國陸域國土可調(diào)查面積的69.5%，另有115×104km2屬高山、荒漠屬1∶250 000地面重力調(diào)查不可工作區(qū)，見圖1。

圖1 全國1:250 000區(qū)域重力調(diào)查工作程度圖Fig.1 1:250 000 regional gravity survey progress map of China

5.3.2 地質(zhì)調(diào)查成果

1)區(qū)域地質(zhì)成果。孫文珂[199]在區(qū)域重力調(diào)查列入中國地質(zhì)工作計(jì)劃的40年后，編寫了《中國區(qū)域重力調(diào)查40年回顧》。他提出①西拉木倫河深大斷裂是板塊結(jié)合帶的祖先線，板塊的多次聚散、開合，導(dǎo)致極為復(fù)雜的地質(zhì)構(gòu)造現(xiàn)象出現(xiàn)；②認(rèn)為秦嶺-桐柏-大別山同屬一個(gè)構(gòu)造帶；③長江中下游構(gòu)造巖漿巖帶是由5個(gè)構(gòu)造巖漿帶構(gòu)成的以銅鐵為主的成礦帶；④揚(yáng)子陸塊與華夏褶皺帶之間的欽州-錢塘結(jié)合帶；⑤對西南三江褶皺帶及相鄰深部構(gòu)造與成礦關(guān)系新認(rèn)識?！度珖V產(chǎn)資源潛力評價(jià)》項(xiàng)目利用全國重力資料，結(jié)合地質(zhì)構(gòu)造與航磁異常特征分析，綜合推斷了我國地質(zhì)構(gòu)造的縱橫格局，研究對一些基礎(chǔ)地質(zhì)問題取得新認(rèn)識。例如：大興安嶺-太行-武陵巨型重力梯度帶是太平洋板塊俯沖的前緣響應(yīng)；二連浩特-賀根山斷裂混雜巖帶是晚古生代大陸邊緣增生的構(gòu)造邊界[200]。楊亞斌等[201]對2010年～2014年度全國區(qū)域重力調(diào)查成果進(jìn)行了梳理，編寫了《重點(diǎn)成礦帶1:250 000區(qū)域重力調(diào)查成果文集》。對大興安嶺梯級帶、郯廬斷裂、古特提斯縫合帶位置提出新認(rèn)識；探討均衡重力異常與地震關(guān)系；分析巖漿巖與礦產(chǎn)關(guān)系等。

2)地球動(dòng)力學(xué)成果。隨著地球科學(xué)研究的不斷深入，重力資料研究深部構(gòu)造及地球動(dòng)力學(xué)的水平也在不斷提高。王謙[202]身利用重力資料研究汶川地震。根據(jù)均衡重力異常圖看出，龍門山、邛崍山及其鄰近地區(qū)存在十分明顯的大面積很均衡的區(qū)域，認(rèn)為這可能是該處地殼不穩(wěn)定運(yùn)動(dòng)的重要根源之一。楊文采[203]利用重力資料分析和反演，刻劃青藏高原地殼分層的三維密度結(jié)構(gòu)，取得的6個(gè)等效層密度擾動(dòng)圖件，為研究地殼構(gòu)造和物質(zhì)運(yùn)動(dòng)提供了重要佐證。

5.4 在其他方面的應(yīng)用

重力在城市地質(zhì)調(diào)查、地?zé)嵴{(diào)查、工程勘查等方面也發(fā)揮了很重要的作用，取得了不少成果，如秦始皇陵的探測[204]。

6 展望

由于重力儀性能的提高，測量技術(shù)的發(fā)展、數(shù)據(jù)精度和處理方法、技術(shù)及解釋水平的提高，重力勘探在地質(zhì)工作中發(fā)揮著越來越重要作用，不僅在傳統(tǒng)的成礦預(yù)測、圈定巖體和構(gòu)造等方面提供強(qiáng)有力的支撐，也將不斷拓寬服務(wù)領(lǐng)域。目前和未來的發(fā)展趨勢包括以下幾個(gè)方面：

1)隨著便攜式原子絕對重力儀進(jìn)入實(shí)用階段，實(shí)測重力數(shù)據(jù)精度將進(jìn)一步提高；地形高程進(jìn)一步加密、數(shù)據(jù)精度進(jìn)一步提高，地形改正精度也將提高一個(gè)量級，從而提高布格異?？偩?。

2)測量參數(shù)由單一的加速度標(biāo)量測量，向多方向矢量測量發(fā)展，矢量數(shù)據(jù)獲取、處理、解釋將得到發(fā)展、完善。

3)隨著場分離、弱信號提取、聯(lián)合反演技術(shù)發(fā)展，將不斷提高三維反演的準(zhǔn)確度與精確度，進(jìn)一步提升地質(zhì)解釋水平。

4)重力勘探測量方式將由傳統(tǒng)的靜態(tài)測量，向靜態(tài)測量與動(dòng)態(tài)監(jiān)測并重發(fā)展。

5)服務(wù)領(lǐng)域不斷拓展，將由傳統(tǒng)的地質(zhì)勘查領(lǐng)域擴(kuò)展到地質(zhì)災(zāi)害監(jiān)測、深部地下水監(jiān)測等領(lǐng)域。