考慮到時誤差的地震定位算法及其在四川地區(qū)2001—2008年地震定位的應(yīng)用

萬永革 盛書中 程萬正 張致偉 武 曄趙曉燕 卜玉菲 薛志芳 劉金麗

1)防災(zāi)科技學(xué)院，河北三河燕郊開發(fā)區(qū) 065201

2)四川省地震局，成都 610041

3)甘肅省地震局，蘭州 730000

4)河北省地震局秦皇島中心臺，秦皇島 066100

5)河北省曲周縣實驗中學(xué)，邯鄲 057250

考慮到時誤差的地震定位算法及其在四川地區(qū)2001—2008年地震定位的應(yīng)用

萬永革1)盛書中1)程萬正2)張致偉2)武 曄1)趙曉燕1)卜玉菲3)薛志芳4)劉金麗5)

1)防災(zāi)科技學(xué)院，河北三河燕郊開發(fā)區(qū) 065201

2)四川省地震局，成都 610041

3)甘肅省地震局，蘭州 730000

4)河北省地震局秦皇島中心臺，秦皇島 066100

5)河北省曲周縣實驗中學(xué)，邯鄲 057250

精確的地震位置對于地震活動性、地震層析成像和地殼應(yīng)力場反演具有相當(dāng)重要的意義，對于地震速報也具有重要的應(yīng)用價值。將觀測到時的不確定性、臺站高程、地震震源深度進(jìn)行約束的同時，根據(jù)反演理論給出了地震震源位置精確估計和誤差估計的方法。該算法聯(lián)合考慮Pg波、Sg波、Pn波和Sn波的到時進(jìn)行反演，數(shù)據(jù)量的增加可以增強地震位置的準(zhǔn)確性，并可同時應(yīng)用于地方震和區(qū)域地震。采用模擬數(shù)據(jù)對該地震定位算法進(jìn)行檢驗發(fā)現(xiàn)，該算法在觀測數(shù)據(jù)的不確定性不等時明顯優(yōu)于其他方法。將該算法應(yīng)用于四川地區(qū)2001—2008年間的地震定位，得到的地震位置更加符合地震的叢集性并集中于斷裂帶附近。這些結(jié)果為四川地區(qū)的地震活動性、斷層構(gòu)造以及地震層析成像研究打下了基礎(chǔ)，并且為汶川地震之前的地震活動前兆研究也提供了有益幫助。

地震定位 不確定性 臺站高程 四川地區(qū)

0 引言

地震定位是地震學(xué)中最基本的問題之一，對于研究地震活動構(gòu)造、地球內(nèi)部結(jié)構(gòu)、震源的幾何構(gòu)造等地震學(xué)基本問題都有重要意義。另外，基于快速準(zhǔn)確的地震定位的地震速報，對于震后的減災(zāi)、救災(zāi)工作也具有重要的實際意義。因此，地震學(xué)家一直在不斷改進(jìn)或提出新的定位方法。

地震定位算法有很多種(Lee et al.，1975;趙仲和，1983)。一些研究者采用S波與P波的到時差和虛波速度確定地震的位置，如Inglada法。該方法的優(yōu)點在于不需要給定初始震中位置，但采用了與實際相差較大的均勻地殼模型(徐果明等，1982)。最常用的地震定位算法為Geiger法(Geiger，1912)。這種方法可采用較為復(fù)雜的地殼模型，甚至三維不均勻速度結(jié)構(gòu)(Wu et al.，2008)，其好處是能夠根據(jù)與實際相近的地殼模型反演震源位置，但必須給出與真實震源位置相差不大的初始震源位置和發(fā)震時刻。還有的研究者是采用遺傳算法(萬永革等，1995，1997)和單純形法(趙珠等，1994)，這類方法的好處是不需要求理論走時對震源參數(shù)的偏導(dǎo)數(shù)，缺點是搜索較慢，且難以估計參數(shù)的不確定性。還有一些研究者采用相對定位法，如主事件定位法(Spence et al.，1980;周仕勇等，1999)、雙差定位法(Waldhauser et al.，2000)。主事件定位法必須給出較為準(zhǔn)確的主事件位置，并且其他事件與主事件的空間距離要比震源到臺站的距離小很多。雙差定位法可以得到事件之間的相對位置，其絕對位置的準(zhǔn)確性也能得到保證，但該方法不能對零散地震進(jìn)行定位。

牟磊育等(2006)研究了原始的Geiger法，給出了聯(lián)合采用Inglada法和Geiger法對區(qū)域臺網(wǎng)地震到時數(shù)據(jù)的精確定位，編寫了MATLAB程序。雖然該方法是一種較為理想的區(qū)域臺網(wǎng)定位法，但該程序沒有利用Pn波到時，沒有考慮各種震相到時拾取準(zhǔn)確性的差別，并且不能給出定位參數(shù)的不確定范圍。本文在他們工作的基礎(chǔ)上增加了Pn波到時的利用，并且考慮了不同震相到時拾取的不同誤差，利用反演理論給出了估計震源參數(shù)的標(biāo)準(zhǔn)差。將我們研制的方法應(yīng)用于四川地區(qū)2001—2008年的地震定位，得到了較為精確的地震位置。

1 定位算法及其改進(jìn)

1.1 考慮觀測誤差的Geiger法反演地震震源參數(shù)

4個參數(shù)描述了發(fā)震時間和震源位置，我們把這些參數(shù)稱為模型，定義模型矢量

采用Geiger法是通過考慮目標(biāo)位置的擾動把問題線性化

這里m0是猜想的最佳位置，可采用Inglada法求出粗略的震源發(fā)震時間和位置(徐果明等，1982)。m是離m0很小距離的新的位置?？梢杂锰├占墧?shù)展開式的第1項來近似表示在m處所預(yù)測的時間

為了使這些殘差最小，我們試圖用下列方式求解模型的修正量Δm

這里G 是偏導(dǎo)數(shù)Gij=?ti/?mj的矩陣，i=1，2，…，n，j=1，…，4。 為了得到位置的調(diào)整值 Δm，用標(biāo)準(zhǔn)的最小二乘法得到方程(4)的最佳擬合。下一步，我們使m0處于m0+Δm的位置，重復(fù)這一過程，直到定位收斂。假如初始猜想的位置離實際位置不是太遠(yuǎn)，迭代過程一般能相當(dāng)快地收斂。

考慮觀測誤差時，式(6)可以寫為

其中，δd為觀測誤差，對于Pg波，Sg波和Pn波可以有不同的誤差，Cd為觀測誤差的協(xié)方差矩陣。利用最大似然解估計公式可得(Jackson et al.，1985;Shen et al.，1996，萬永革等，2004，2008a，b)

解的協(xié)方差矩陣可表示為

該矩陣對角線元素的開方即為解的標(biāo)準(zhǔn)差。

平均殘差定義為

其中m為數(shù)據(jù)個數(shù)，即方程的個數(shù)。

1.2 Pn，Sn波走時及其偏導(dǎo)數(shù)

設(shè)地球平層模型的地殼層的序列編號為l，最下一層為地幔層m，地殼層厚度為dl，速度為Vl，入射角為θl，地幔層的速度為Vm。震源位于j層中，dz為震源層中震源距該層頂部的距離。根據(jù)該地球模型，只有震源層的入射角滿足時才會出現(xiàn)首波，在其他層中的入射角為

臨界觀測到首波的距離為

震中距為Δ的首波的走時為

Pn波走時的偏導(dǎo)數(shù)

1.3 臺站高程的校正

在平原或高原地區(qū)，臺站高程相差不大，對定位沒有太大影響。但如果在地形起伏比較大的地區(qū)，臺站高程往往相差較大，如果不考慮臺站高程會導(dǎo)致定位結(jié)果誤差較大。牟磊育等(2006)的程序因用于平原地區(qū)而沒有采用臺站高程校正。本研究在計算某一臺站走時時，采用地殼模型第1層厚度增加臺站高程值的方法來校正臺站高程的影響，使得定位結(jié)果更為準(zhǔn)確。

1.4 地震震源深度的約束

其他定位方法一般沒有對深度進(jìn)行約束。只是限定如果深度不滿足一定的范圍(如空中地震或深度特別大地震)，則強制給出某一深度值。這樣處理的結(jié)果會導(dǎo)致震源參數(shù)不滿足最小殘差的要求。我們在每一次迭代過程中判斷修正后的震源深度是否在合理的范圍之內(nèi)，如果不合理，只將參數(shù)增加到合理深度的邊界，進(jìn)行下一步迭代。下一步迭代還有可能使得震源深度向較為合理的方向進(jìn)展。最終得到殘差最小且深度在合理范圍內(nèi)的震源參數(shù)。

1.5 地震深度的精確求解

由于地震震相數(shù)據(jù)觀測數(shù)據(jù)有限，并且地震臺站大多位于同一水平面上，對深度的分辨率往往不高。為了使深度有更好的約束，我們分別在初始水平位置的基礎(chǔ)上設(shè)置每一層的中點作為深度的初始深度進(jìn)行迭代求解，選擇殘差最小的解為最優(yōu)解。

2 數(shù)值實驗

為了測試以上的定位方法，假設(shè)一個由9個臺站組成的虛擬臺網(wǎng)，虛擬臺網(wǎng)之形狀為邊長300km的正方形(圖1)。虛擬臺網(wǎng)下的速度模型按照地球物理較為全面的觀測結(jié)果得到的Crust5.1(Mooney et al.，1998)給定。具體參數(shù)見表 1。

圖1 虛擬臺站及虛擬地震位置Fig.1 Locations of the virtual stations and earthquake.臺站的編號同表1

假設(shè)在虛擬臺網(wǎng)內(nèi)發(fā)生一個地震，發(fā)震時間 0時刻，位于橫軸 75.071km，縱軸50km處，深度為19km(圖1)。為簡單起見，虛擬臺網(wǎng)的高程全部設(shè)為0，根據(jù)虛擬臺網(wǎng)內(nèi)臺站的位置坐標(biāo)、速度模型、震中位置和發(fā)震時間，計算出虛擬臺網(wǎng)內(nèi)各臺站的Pg、Sg和Pn的到時，分別隨機產(chǎn)生以0.5，2.5，1.5幅度的隨機干擾，得到模擬 Pg、Sg和 Pn的到時(表2)，然后采用我們修訂的程序進(jìn)行定位計算。初始位置X，Y，Z，t分別設(shè)定為50，30，18，0。當(dāng)采用相同的到時拾取誤差時，得到的地震震源 X，Y，Z，t分別為 76.0866，50.1371，18.1681，0.1544;采用與虛擬隨機擾動相同的誤差進(jìn)行反演得到的 X，Y，Z，t分別為 75.6551km、50.1251km、17.6933km、0.1366s，得到的該參數(shù)的標(biāo)準(zhǔn)差為 1.4504km、1.5413km、4.2750km 和0.1958s。采用不同的到時拾取誤差得到的結(jié)果更加接近虛擬設(shè)置的地震位置，可見引用相對到時標(biāo)準(zhǔn)差的重要性。

3 2001—2008年的四川地區(qū)地震定位實例

采用我們的程序?qū)λ拇ǖ貐^(qū)2001—2008年的地震進(jìn)行了研究，對震相進(jìn)行核對，得到1,831個地震的位置(圖2)。原來所有地震震相的平均殘差減少為原來的0.35%。定位后的位置更加集中在斷裂帶上，并具有更強的叢集性，反映了小震分布的基本特征。重新定位后的小震叢集區(qū)主要為龍門山斷裂帶、鮮水河斷裂帶及其西南區(qū)域、馬邊-鹽津斷裂帶附近及華鎣山斷裂帶北部和東南區(qū)域。很多研究者采用雙差定位法來確定地震的叢集特性(Waldhauser et al.，2000)，其重定位是將地震組成地震簇而精確確定地震之前的相對位置，由于去除了不能組成地震簇的離散地震，因此能反映地震的叢集特性從而確定斷層的形狀(萬永革等，2008a;王福昌等，2008)。而本研究沒有去除離散地震，也明顯表現(xiàn)出地震的叢集特性，說明了研究結(jié)果的客觀性。趙珠等(1997)采用該地區(qū)早期較長時間的地震資料確定的地震位置也表現(xiàn)為地震叢集在龍門山斷裂帶附近的特性，與本文的研究結(jié)果有類似之處，表明了本文處理方法的正確性。

在深度上，地震目錄列出的某些震源深度達(dá)60km左右(已達(dá)不易發(fā)生地震的上地幔)，重新定位后的震源深度主要集中在8km以上的上地殼(圖3)，下地殼的地震較少。這與Takagi等(1997)的結(jié)果類似。另外，由于這期間大部分為小震，均發(fā)生在淺部的4～8km的范圍內(nèi)，深部缺少地震，似乎可以稱之為深部“地震空區(qū)”，這也許是地震預(yù)測研究人員要注意的現(xiàn)象。

4 結(jié)論與討論

與其他地震定位方法對比，我們的方法有如下優(yōu)勢:1)可以利用Pn波到時數(shù)據(jù);2)可以對數(shù)據(jù)給予不同的觀測誤差，從而確定地震定位的不確定性;3)考慮了臺站高程;4)可以約束震源深度的搜索范圍;5)通過對各個層中分別設(shè)置震源來去除由于地殼結(jié)構(gòu)各層速度差別較大而造成非最優(yōu)解的情況。實驗結(jié)果表明，在數(shù)據(jù)具有較大誤差的情況下，我們的程序仍然可以得到較為精確的震源位置。將該算法用于2001年至汶川地震之前的四川地區(qū)地震記錄，得到了殘差更小的震源位置。

圖2 定位的地震震源位置Fig.2 The hypocentral locations of earthquakes located.

將該算法用于2001年至汶川地震之前的四川地區(qū)1831個小震震源精確定位，結(jié)果更為精確地顯示了地震的叢集性和在斷裂帶附近發(fā)生的特性，并且震源深度集中在淺部4～8km的范圍內(nèi)。這些結(jié)果為進(jìn)一步研究汶川地震之前的地震活動特征及分析地震前兆現(xiàn)象打下了基礎(chǔ)。

本文提出了考慮P波、S波和Pn波到時拾取誤差來進(jìn)行地震定位的方法。實際上地震波到時的讀取和地震定位由不同的人員進(jìn)行，僅根據(jù)地震觀測報告無法得到確切的到時讀取的誤差。因此，本文采用直達(dá)P波到時的精度為0.5、S波到時拾取誤差為1.5s、Pn波到時的誤差為1.0s的原則來進(jìn)行定位。更為嚴(yán)格的地震定位需要給定不同臺站不同拾取震相不同到時的誤差，這樣才能更好地利用到時數(shù)據(jù)的拾取誤差。因此，為了更好確定地震位置，建議在地震觀測報告中增加震相拾取的誤差參數(shù)。

圖3 1831個地震的深度分布Fig.3 The depth distribution of the 1831 earthquakes.

我們編寫的MATLAB程序可向感興趣的讀者提供。

致謝 審稿人對本文提出了建設(shè)性修改意見，特此致謝。

牟磊育，趙仲和，張偉，等.2006.用INGLADA與GEIGER方法實現(xiàn)近震精定位［J］.中國地震，22(3):294—302.

MOU Lei-yu，ZHAO Zhong-he，ZHANGWei，et al.2006.Using the INGLADA and GEIGER method to realize regional earthquake location accurately［J］.Earthquake Research in China，22(3):294—302(in Chinese).

萬永革，李鴻吉.1995.遺傳算法在確定震源位置中的應(yīng)用［J］.地震地磁觀測與研究，16(6):1—7.

WAN Yong-ge，LI Hong-ji.1995.The preliminary study on the seismic hypocenter location using genetic algorithms［J］.Seismological and Geomagnetic Observation and Research，16(6):1—7(in Chinese).

萬永革，李清河，李鴻吉，等.1997.用遺傳算法確定三維橫向不均勻介質(zhì)中的近震震源位置［J］.西北地震學(xué)報，19(2):7—14.

WAN Yong-ge，LIQing-he，LIHong-ji，et al.1997.Hypocentral locations of the near earthquakes in 3-D lateral heterogeneousmedium determined by the genetic algorithms［J］.Northwestern Seismological Journal，19(2):7—14(in Chinese).

萬永革，王敏，沈正康，等.2004.利用GPS和水準(zhǔn)測量資料反演2001年昆侖山口西8.1級地震的同震滑動分布［J］.地震地質(zhì)，26(3):393—404.

WAN Yong-ge，WANG Min，SHEN Zheng-kang，et al.2004.Co-seismic slip distribution of the 2001 west of Kunlun Mountain Pass earthquake inverted by GPSand leveling data［J］.Seismology and Geology，26(3):393—404(in Chinese).

萬永革，沈正康，刁桂苓，等.2008a.利用小震分布和區(qū)域應(yīng)力場確定大震斷層面參數(shù)方法及其在唐山地震序列中的應(yīng)用［J］.地球物理學(xué)報，51(3):793—804.

WAN Yong-ge，SHEN Zheng-kang，DIAO Gui-ling，et al.2008a.An algorithm of fault parameter determination using distribution of small earthquakes and parameters of regional stress field and its application to Tangshan earthquake sequence［J］.Chinese JGeophys，51(3):569—583(in Chinese).

萬永革，沈正康，王敏，等.2008b.根據(jù)GPS和InSAR數(shù)據(jù)反演2001年昆侖山口西地震同震破裂分布［J］.地球物理學(xué)報，51(4):1074—1084.

WAN Yong-ge，SHEN Zheng-kang，WANG Min，et al.，2008b.Coseismic slip distribution of the 2001 West Kunlun Mountain Pass earthquake constrained using GPSand InSAR data［J］.Chinese JGeophysics，51(4):753—764.

萬永革，沈正康，盛書中，等.2009.2008年汶川大地震對周圍斷層的影響［J］.地震學(xué)報，31(2):128—139.

WAN Yong-ge，SHEN Zheng-kang，SHENG Shu-zhong，etal.2009.The influence of2008Wenchuan earthquake on surrounding faults［J］.Acta Seismologica Sinica，31(2):128—139(in Chinese).

王福昌，萬永革，胡順田.2008.粒子群算法在主震斷層面參數(shù)估計中的應(yīng)用［J］.地震研究，31(2):149—154.

WANG Fu-chang，WAN Yong-ge，HU Shun-tian.2008.Application of particle swarm optimization to the estimation of mainshock fault plane parameters［J］.Journal of Seismological Research，31(2):149—154(in Chinese).

徐果明，周惠蘭.1982.地震學(xué)原理［M］.北京:科學(xué)出版社.

XU Guo-ming，ZHOU Hui-lan.1982.Principle of Seismology［M］.Science Press，Beijing(in Chinese).

趙珠，丁志峰，易桂喜，等.1994.西藏地震定位:一種使用單純形優(yōu)化的非線性方法［J］.地震學(xué)報，16(2):212—219.

ZHAO Zhu，DING Zhi-feng，YIGui-xi，et al.1994.Location of Tibetan earthquakes:A nonlinear approach by a simplex optimized technique［J］.Acta Seismologica Sinica，16(2):212—219(in Chinese).

趙珠，范軍，鄭斯華，等.1997.龍門山斷裂帶地殼速度結(jié)構(gòu)和震源位置的精確修定［J］.地震學(xué)報，19(6):615—622.

ZHAO Zhu，F(xiàn)AN Jun，ZHENG Si-hua，et al.1997.Crustal structure and accurate hypocenter determination along the Longmenshan Fault zone［J］.Acta Seismologica Sinica，19(6):615—622(in Chinese).

趙仲和.1983.多重模型地震定位程序及其在北京臺網(wǎng)的應(yīng)用［J］.地震學(xué)報，5(2):242—254.

ZHONG Zhong-he.1983.An earthquake location program withmultiple velocitymodel and its application in the Beijing seismic network［J］.Acta Seismologica Sinica，5(2):242—254(in Chinese).

周仕勇，許忠淮，韓京，等.1999.主地震定位方法分析以及1997年新疆伽師震群高精度定位［J］.地震學(xué)報，21(3):258—265.

ZHOU Shi-yong，XU Zhong-huai，HAN Jing，et al.1999.Analysis on the master eventmethod and precise location of 1997 Jiashi strong earthquake swarm in western China［J］.Acta Seismologica Sinica，21(3):258—265(in Chinese).

Geiger L.1912.Probability method for the determination of earthquake epicenters from arrival time only ［J］.Bull St Louis Univ，8:60—71.

Jackson D D，Matsu'ura M.1985.A Bayesian approach to nonlinear inversion［J］.JGeophys Res，90(B1):581—591.

Lee W H K，Lahr JC.1975.A computer program for determining hypocenter，magnitude，and firstmotion pattern of local earthquakes［J］.U SGeol Surv Open-file Rep，75—311.

Mooney W D，Laske G，Masters G.1998.A global crustalmodel at5°x 5°［J］.JGeophys Re，103:727—747.

Shen Z K，Jackson D，Ge B.1996.Crustal deformation across and beyond the Los Angeles Basin from geodeticmeasurements［J］.JGeophys Res，101(B12):27957—27980.

SpenceW.1980.Relative epicenter determination using P-wave arrival-time differences［J］.Bull Seism Soc Amer，70(1):171—183.

Takagi A，Hasegawa A，Umino N.1977.Seismic activity in the Northeastern Janpan Arc［J］.JPhys Earth，25(Suppl):95—104.

Waldhauser F，Ellsworth W L.2000.A double-difference earthquake location algorithm:Method and application to the Northern Hayward Fault，California［J］.Bull Seism Soc Amer，90(6):1353—1368.

Wu Y M，Chang C H，Zhao L，et al.2008.A comprehensive relocation of earthquakes in Taiwan from 1991 to 2005［J］.Bull Seism Soc Amer，98(3):1471—1481.

EARTHQUAKE LOCATION METHODW ITH ARRIVAL TIME UNCERTAINTY CONSIDERED AND ITS APPLICATION TO LOCATION OF EARTHQUAKES FROM 2001 TO 2008 IN SICHUAN AREA

WAN Yong-ge1)SHENG Shu-zhong1)CHENGWan-zheng2)ZHANG Zhi-wei2)WU Ye1)ZHAO Xiao-yan1)BU Yu-fei3)XUE Zhi-fang4)LIU Jin-li5)

1)Institute of Disaster Prevention，Yanjiao，Hebei 065201，China

2)Sichuan Earthquake Administration，Chengdu 610041，China

3)Gansu Earthquake Administration，Lanzhou 730000，China

4)Qinhuangdao Center Station，Hebei Earthquake Administration，Qinhuangdao 066100，China

5)Quzhou Experimental High School，Hebei Province，Quzhou 057250，China

The precisely located earthquake catalogue is important to seismicity，seismic tomography and crustal stress inversion studies.It also has greatapplication value in rapid report of an earthquake that just occurred.By considering the arrival time uncertainty，and the constraints on station elevation and seismic depth，we propose a relatively accuratemethod to estimate hypocentral location and its uncertainty based on inversion theory.Ourmethod can combine the arrival times of Pg wave，Sg wave，Pn wave and Sn wave in hypocenter location，so it increases the location accuracy by involvingmore data;and it can be also used in local and regional earthquake location simultaneously.In order to test our location method，we located earthquakes by using the simulated data with different uncertainty of Pg，Sg，Pn，Sn arrivals.The result shows that the location determined by using ourmethod ismore accurate than that by using othermethod.We apply it to earthquakes occurring in the period from 2001 to 2008 in Sichuan area，and obtained amore clustered hypocentral distribution convergent to the fault zones.The result provides a solid foundation for studies of seismicity，geometry of the active faults and seismic tomography in Sichuan region.It is also helpful to study the seismicity precursors before theWenchuan earthquake.

earthquake location，uncertainty，station elevation，Sichuan area

P315.5

A

0253-4967(2012)01-0001-10

10.3969/j.issn.0253-4967.2012.01.001

2011-05-03收稿，2011-07-30改回。

國家自然科學(xué)基金(41074072、40874022)、中國地震局教師科研基金(20100101)、中央高?；究蒲袠I(yè)務(wù)費專項(ZY20110101)和地震行業(yè)專項(200808053，201208009)共同資助。

萬永革，男，1967年生，2001年在中國地震局地球物理研究所獲得博士學(xué)位，2002—2006年在中國地震局地質(zhì)研究所進(jìn)行博士后合作研究，研究員，主要從事地震學(xué)和地球動力學(xué)方面的研究，電話:010-61597607;E-mail:wanyg217217@vip.sina.com;wanyongge8@gmail.com。