• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    利用協(xié)同反演方法反演地震序列滑動(dòng)分布

    2021-08-03 11:10:36王樂洋谷旺旺孫龍翔
    地球物理學(xué)報(bào) 2021年8期
    關(guān)鍵詞:主震殘差滑動(dòng)

    王樂洋 , 谷旺旺, 孫龍翔

    東華理工大學(xué)測繪工程學(xué)院, 南昌 330013

    0 引言

    大地測量反演是研究地球科學(xué)的關(guān)鍵手段,通常采用大地測量數(shù)據(jù)反演斷層面滑動(dòng)的大小與分布、應(yīng)力分布以及地震機(jī)制等(Xu et al., 2017; Wang et al., 2018b; Wang and Gu, 2020).隨著大地測量觀測技術(shù)的不斷發(fā)展以及觀測手段的增加,越來越多的觀測資料用于滑動(dòng)分布反演中,其中GPS和InSAR以其高時(shí)間和空間分辨率成為滑動(dòng)分布反演的重要數(shù)據(jù).隨著數(shù)據(jù)種類的增加,聯(lián)合不同數(shù)據(jù)進(jìn)行反演成為研究的方向.由于大地測量反演問題通常是病態(tài)的,導(dǎo)致不同數(shù)據(jù)反演結(jié)果存在差異;顯而易見,能夠同時(shí)解釋多類數(shù)據(jù)的模型更有說服力.由于不同類數(shù)據(jù)的優(yōu)勢與劣勢各不相同,如何揚(yáng)長避短彌補(bǔ)不同數(shù)據(jù)之間的劣勢,發(fā)揮多類數(shù)據(jù)的優(yōu)勢成為研究的熱點(diǎn).其中,聯(lián)合反演方法有效將多類數(shù)據(jù)進(jìn)行融合,統(tǒng)一反演最優(yōu)模型.而聯(lián)合反演應(yīng)用的關(guān)鍵問題在于相對權(quán)比的確定,它體現(xiàn)了各類數(shù)據(jù)在聯(lián)合反演中的貢獻(xiàn)程度.目前已有大量學(xué)者對此進(jìn)行相關(guān)研究,其中Xu等(2009)總結(jié)目前較為常見的定權(quán)方法:將相對權(quán)比作為未知參數(shù)統(tǒng)一求解;根據(jù)先驗(yàn)信息確定;將不同類數(shù)據(jù)視為等權(quán)處理;采用赫爾默特方差分量估計(jì)法(Helmert Variance Component Estimation, HVCE)確定.Xu等(2009)表明HVCE法能夠得到合理的相對權(quán)比.Xu等(2006)研究了病態(tài)問題中方差分量估計(jì)法的適用性,提出對參數(shù)的殘差進(jìn)行偏差改正的方差分量估計(jì)法.Xu(2009)研究了利用正則化方法求解的參數(shù)是有偏的,通過偏差改正以獲得更精確的驗(yàn)后方差.王樂洋等(2012)提出了確定聯(lián)合反演相對權(quán)比的兩步法,并針對不同目標(biāo)函數(shù)對相對權(quán)比的影響進(jìn)行了討論.Huang等(2013)采用方差縮減法確定聯(lián)合反演的相對權(quán)比和正則化參數(shù).Wang等(2018a)采用U曲線法確定正則化參數(shù),判別函數(shù)最小化法確定相對權(quán)比,實(shí)驗(yàn)結(jié)果表明該方法能夠在均方誤差意義下優(yōu)于方差分量估計(jì)法.許才軍等(2016)提出采用方差分量估計(jì)法同時(shí)確定相對權(quán)比和正則化參數(shù),拓展了方差分量估計(jì)法的應(yīng)用范圍.Goldberg等(2020)根據(jù)不同相對權(quán)比下GPS與InSAR的均方根誤差確定兩類數(shù)據(jù)的相對權(quán)比.以上方法均將不同觀測數(shù)據(jù)聯(lián)合到一個(gè)模型下進(jìn)行求解,通過權(quán)比控制不同數(shù)據(jù)反演所占權(quán)重.考慮到相對權(quán)比的確定對反演的結(jié)果會(huì)產(chǎn)生不同的影響,而相對權(quán)比的確定并沒有一個(gè)統(tǒng)一標(biāo)準(zhǔn).在地震序列反演中,由于部分?jǐn)?shù)據(jù)為單次地震獲取的同震形變場,而其他數(shù)據(jù)包含了整個(gè)地震序列的同震形變,這將導(dǎo)致這部分?jǐn)?shù)據(jù)無法用于單次地震的反演.Li等(2020)表明采用整體形變數(shù)據(jù)進(jìn)行反演會(huì)造成反演結(jié)果存在折衷,并不能反演實(shí)際滑動(dòng)分布.針對此問題,Goldberg等(2020)根據(jù)地質(zhì)信息勘測到單次地震的破裂情況,利用多種數(shù)據(jù)進(jìn)行聯(lián)合反演,并根據(jù)數(shù)據(jù)的均方根誤差確定相對權(quán)比.Wang等(2020)通過單次地震數(shù)據(jù)作為先驗(yàn)信息,通過迭代的方式將包含整個(gè)地震序列的同震形變進(jìn)行分離,從而反演出單次地震的滑動(dòng)分布.以上方法均將不同數(shù)據(jù)聯(lián)合到同一個(gè)方程中進(jìn)行求解,由于InSAR數(shù)據(jù)為兩次地震共同產(chǎn)生的形變,為了得到分離的形變數(shù)據(jù),考慮采用協(xié)同反演方法進(jìn)行反演.

    協(xié)同反演方法在地球物理聯(lián)合反演中具有廣泛的應(yīng)用,Lines等(1988)首次提出協(xié)同反演方法,并利用地表地震數(shù)據(jù)、聲波測井、地表重力數(shù)據(jù)和井眼重力儀數(shù)據(jù)進(jìn)行了反演.隨后大量學(xué)者對此進(jìn)行了深入研究.Anderson等(1998)采用順序聯(lián)合反演法利用地震和重力資料反演速度和密度參數(shù),為深度偏移成像提供準(zhǔn)確的速度模型.Paasche等(2012)擴(kuò)展了基于模糊C-means聚類分析和單一數(shù)據(jù)集反演算法的區(qū)域協(xié)同反演方法,用于部分位于同一模型區(qū)域的數(shù)據(jù)的協(xié)同反演.McMillan和Oldenburg(2014)采用協(xié)同反演方法反演了具有邊界約束的電磁數(shù)據(jù),得到了具有改進(jìn)分辨率的一致的3D電阻率模型.Takam-Takougang等(2015)提出了一種協(xié)同反演方法,利用地震和大地電磁數(shù)據(jù)聯(lián)合反演聲阻抗.Le等(2016)測試了不同聯(lián)合反演策略,驗(yàn)證了所提出的協(xié)同反演方法與其他反演手段一致的電導(dǎo)率分布.認(rèn)為在并行計(jì)算的幫助下,地震的大地電磁數(shù)據(jù)的協(xié)同反演可以實(shí)現(xiàn)自動(dòng)化.Moorkamp(2017)討論了聯(lián)合反演和約束反演方法的基本原理,并討論了聯(lián)合反演方法中不同耦合方法的特點(diǎn).Gon?alves和Leite(2019)提出了一種利用復(fù)雜地質(zhì)結(jié)構(gòu)疊后地震反射和重力數(shù)據(jù)偏移2D協(xié)同反演方法,通過減少變量數(shù)量提高了算法的計(jì)算效率.Singh等(2019)提出了一種利用模糊C-means聚類協(xié)同反演方法反演出更可靠的電阻率和密度模型.

    協(xié)同反演方法在地質(zhì)勘探中應(yīng)用較多,但并沒有應(yīng)用于滑動(dòng)分布反演中.在大地測量反演中,聯(lián)合反演均采用同步聯(lián)合反演,即將多類數(shù)據(jù)聯(lián)合到一個(gè)模型中進(jìn)行解算,采用不同的相對權(quán)比區(qū)分每類數(shù)據(jù)的權(quán)重影響;而并未考慮到協(xié)同反演方法,即順序聯(lián)合反演,利用單一類型數(shù)據(jù)進(jìn)行反演,將反演結(jié)果作為另一類數(shù)據(jù)的先驗(yàn)信息進(jìn)行反演,如此反復(fù)迭代.由于順序聯(lián)合反演方法具有無需考慮相對權(quán)比的優(yōu)勢,本文將順序聯(lián)合反演方法引入到大地測量滑動(dòng)分布反演中,一方面可以有效地獲得地震序列的滑動(dòng)分布情況,另一方面避免了相對權(quán)比確定問題.

    1 協(xié)同反演

    1.1 地震同震滑動(dòng)分布反演基本原理

    在同震滑動(dòng)分布反演中,地表同震形變位移與斷層滑動(dòng)之間可以通過線性函數(shù)連接,其表達(dá)式為

    d=Gm,

    (1)

    其中,d表示觀測向量(如InSAR視線向位移、GPS位移場等);G表示格林函數(shù)矩陣,利用均勻彈性半空間矩形位錯(cuò)模型(Okada, 1985)確定;m為滑動(dòng)參數(shù).

    利用最小二乘法對式(1)進(jìn)行求解時(shí),由于格林函數(shù)矩陣的復(fù)共線性,造成了求解過程中法矩陣病態(tài),通常對滑動(dòng)參數(shù)添加一定的平滑約束,使系數(shù)矩陣滿秩,且能夠避免滑動(dòng)解出現(xiàn)非物理因素的震蕩.一般采用拉普拉斯二階平滑算子對斷層單元進(jìn)行平滑約束,可表示為

    Hm=0,

    (2)

    聯(lián)合式(1)和式(2),得到滑動(dòng)分布反演公式為

    (3)

    利用最小二乘法求解式(3),即可獲得滑動(dòng)參數(shù)解.滑動(dòng)參數(shù)解可表示為

    m=(GTPG+α2HTH)-1GTPd,

    (4)

    其中,P為觀測數(shù)據(jù)的權(quán)陣,α為平滑因子,文中采用方差分量估計(jì)法進(jìn)行確定.

    1.2 方差分量估計(jì)法基本原理

    根據(jù)式(1)和(2),得到GPS和InSAR數(shù)據(jù)聯(lián)合反演滑動(dòng)分布公式(Xu et al., 2017; Wang et al., 2018a)為

    (5)

    利用方差分量估計(jì)法對式(5)進(jìn)行迭代求解,可以得到GPS和InSAR數(shù)據(jù)的相對權(quán)比和平滑因子,當(dāng)僅為GPS或者InSAR數(shù)據(jù)時(shí),采用方差分量估計(jì)法可以用于確定平滑因子.

    1.3 協(xié)同反演方法基本原理

    利用多類數(shù)據(jù)進(jìn)行協(xié)同反演時(shí),首先利用一類數(shù)據(jù)進(jìn)行滑動(dòng)分布反演,將反演的滑動(dòng)參數(shù)作為另一類數(shù)據(jù)反演的先驗(yàn)信息進(jìn)行約束.假設(shè)本文僅考慮利用GPS和InSAR數(shù)據(jù)進(jìn)行協(xié)同反演.由于在協(xié)同反演中進(jìn)行了迭代運(yùn)算,故利用何種數(shù)據(jù)首先進(jìn)行滑動(dòng)分布反演對結(jié)果并沒有影響.由于本文考慮將協(xié)同反演方法應(yīng)用于地震序列中,假設(shè)一個(gè)地震序列包含了兩次地震,第一次地震稱為前震,第二次地震為主震.由于GPS數(shù)據(jù)的時(shí)間分辨率非常高,能夠得到單次地震的同震形變,故本文首先分別利用兩次GPS數(shù)據(jù)進(jìn)行滑動(dòng)分布反演,得到滑動(dòng)分布反演公式為

    (6)

    利用最小二乘法求解式(6),得到利用GPS數(shù)據(jù)反演兩次地震的滑動(dòng)分布解m1、m2,其中平滑因子采用方差分量估計(jì)法確定.

    根據(jù)式(6)可以得到前震和主震的滑動(dòng)分布,將GPS數(shù)據(jù)反演結(jié)果作為InSAR數(shù)據(jù)反演的先驗(yàn)信息,可以得到InSAR數(shù)據(jù)的滑動(dòng)分布反演公式為:

    (7)

    利用最小二乘法求解式(7),得到利用GPS數(shù)據(jù)反演兩次地震的滑動(dòng)分布解m′1、m′2.

    根據(jù)式(7)得到InSAR數(shù)據(jù)約束下的滑動(dòng)分布,將反演結(jié)果作為下一次GPS數(shù)據(jù)反演的先驗(yàn)信息,可以得到附有先驗(yàn)信息的滑動(dòng)分布反演公式為

    (8)

    利用式(7)和(8)進(jìn)行迭代,最終能夠得到擬合兩類數(shù)據(jù)的滑動(dòng)模型.

    本文給出利用GPS和InSAR數(shù)據(jù)進(jìn)行協(xié)同反演的迭代步驟:

    (1)分別利用前震和主震GPS數(shù)據(jù)進(jìn)行滑動(dòng)分布反演,反演公式如式(6),利用方差分量估計(jì)法確定正則化參數(shù),采用最小二乘法求解滑動(dòng)分布,得到前震與主震的滑動(dòng)分布分別為m1和m2;

    (2)將斷層滑動(dòng)分布m1和m2作為先驗(yàn)信息,采用最小二乘法求解式(7),得到InSAR數(shù)據(jù)約束下反演的滑動(dòng)分布分別為m′1和m′2;

    (3)將步驟(2)中求得的滑動(dòng)分布m′1和m′2作為先驗(yàn)信息,根據(jù)式(8)利用GPS數(shù)據(jù)反演前震和主震滑動(dòng)分布,得到兩次地震的滑動(dòng)分布為m″1和m″2;

    (4)計(jì)算步驟(2)和步驟(3)分別求得的滑動(dòng)分布解m′1和m′2與m″1和m″2的擬合殘差,若殘差小于閾值,迭代終止;否則重復(fù)步驟(2)—(4).

    協(xié)同反演迭代流程圖如圖1所示.由于設(shè)置不同的閾值,導(dǎo)致迭代的結(jié)果并不相同,當(dāng)閾值設(shè)置較小時(shí),反演得到的形變值殘差較小,但這伴隨著反演的滑動(dòng)分布可能出現(xiàn)非物理震蕩,基于此我們根據(jù)殘差與平滑度之間的折衷選取閾值的大小.

    圖1 協(xié)同反演迭代流程圖Fig.1 Iterative flow chart of cooperation inversion method

    2 模擬實(shí)驗(yàn)

    為了驗(yàn)證協(xié)同反演方法在地震序列中的可行性,本文做了以下模擬實(shí)驗(yàn).在模擬實(shí)驗(yàn)中,模擬了一對正交共軛斷層,其斷層面幾何參數(shù)如表1所示,斷層面分布如圖2所示.模擬兩次地震累計(jì)產(chǎn)生的形變量,并施加觀測誤差N~(0,12cm2),圖3為InSAR數(shù)據(jù)降采樣后的數(shù)據(jù)分布圖.同時(shí)模擬了兩次地震GPS3方向形變數(shù)據(jù)如圖4和圖5所示,并給形變點(diǎn)施加觀測誤差,其中水平方向施加N~(0,32mm2)的觀測誤差,垂直方向施加N~(0,52mm2)的觀測誤差.在模擬實(shí)驗(yàn)中,由于加入的誤差是標(biāo)準(zhǔn)正態(tài)分布的,且誤差較小,故我們設(shè)置了較小的閾值進(jìn)行迭代,閾值為10-5m.

    圖2 模擬實(shí)驗(yàn)的滑移分布結(jié)果與殘差分布(a) 聯(lián)合反演方法反演的滑移分布結(jié)果; (b) 協(xié)同反演方法反演的滑移分布結(jié)果; (c) 聯(lián)合反演方法的殘差分布; (d) 協(xié)同反演方法的殘差分布.Fig.2 The results of simulated slip distribution inversion and their residuals(a) The results of slip distribution inverted by the joint inversion method; (b) The results of slip distribution inverted by the cooperation inversion method; (c) The results of residuals by the joint inversion method; (d) The results of residuals by the cooperation inversion method.

    表1 模擬地震震源參數(shù)Table 1 The real source parameters of the simulated earthquake

    在模擬實(shí)驗(yàn)中,分別采用聯(lián)合反演和協(xié)同反演方法進(jìn)行滑動(dòng)分布反演,其中聯(lián)合反演中利用方差分量估計(jì)法確定的相對權(quán)比為1∶ 0.9697∶ 1.6039∶0.0004,在協(xié)同反演中正則化參數(shù)采用方差分量估計(jì)法進(jìn)行確定,正則化參數(shù)分別為0.0583和0.0136.圖2為兩種方法反演的滑動(dòng)分布及其殘差分布,其中(a)和(b)分別為聯(lián)合反演和協(xié)同反演方法反演的滑動(dòng)分布結(jié)果,(c)和(d)分別為兩種方法對應(yīng)的殘差分布結(jié)果.在模擬實(shí)驗(yàn)中,利用聯(lián)合反演和協(xié)同反演方法得到的地震滑動(dòng)分布各參數(shù)結(jié)果見表2.

    表2 模擬實(shí)驗(yàn)反演地震滑動(dòng)分布結(jié)果Table 2 The results of slip distribution inversion in simulation experiment

    根據(jù)表2反演結(jié)果可以看出,兩種方法均能很好反演出模擬的滑動(dòng)分布,且數(shù)據(jù)的擬合殘差僅有微小差別,其中采用方差分量估計(jì)法反演結(jié)果正演的GPS數(shù)據(jù)的擬合殘差分別為3.5 mm和3.4 mm,InSAR數(shù)據(jù)的擬合殘差為9.1 mm;而采用協(xié)同反演方法反演結(jié)果正演的GPS數(shù)據(jù)的擬合殘差分別為3.8 mm和3.6 mm,InSAR數(shù)據(jù)的擬合殘差為9.0 mm.但兩種方法反演結(jié)果略有不同,其中采用方差分量估計(jì)法確定相對權(quán)比的聯(lián)合反演方法反演出的斷層S2及F2的最大滑動(dòng)量超過了模擬設(shè)置的滑動(dòng)量,而其他斷層的最大滑動(dòng)量與模擬值相比較小,這很有可能是反演不同的斷層之間的滑動(dòng)存在一定的折中,這與Li等(2020)反演結(jié)果相一致.由于模擬實(shí)驗(yàn)中,模擬的觀測數(shù)據(jù)所加入的隨機(jī)誤差較小,反演的結(jié)果與模擬值相差不大.對于協(xié)同反演而言,各個(gè)斷層的滑動(dòng)量與模擬值較為接近,僅前震F3斷層反演的滑動(dòng)量較小,主要是模擬中F3斷層的滑動(dòng)位置處于斷層較深處,且其滑動(dòng)量較小,故難以反演出其真實(shí)滑動(dòng)大小.總體而言,兩種方法的反演結(jié)果較為一致,且其數(shù)據(jù)的擬合殘差較小,僅有微小差異,故能夠驗(yàn)證協(xié)同反演方法應(yīng)用于滑動(dòng)分布反演的有效性.

    圖2展示了利用方差分量估計(jì)法和協(xié)同反演方法反演的主震與前震的滑動(dòng)分布以及兩種反演結(jié)果與模擬滑動(dòng)之間的殘差分布.從圖中可以看出整個(gè)地震產(chǎn)生了多個(gè)斷層的滑動(dòng),故設(shè)置了不用的斷層來反演每個(gè)斷層的滑動(dòng)情況;對比兩種方法的反演結(jié)果,主要滑動(dòng)區(qū)的滑動(dòng)保持一致,僅在非滑動(dòng)區(qū)存在較小的差異,這是平滑約束所引入的,無法通過本文方法消除;且設(shè)置的最大滑動(dòng)量和反演的滑動(dòng)量之間較一致,這也與模擬數(shù)據(jù)所加入的誤差有關(guān),模擬數(shù)據(jù)的誤差為正態(tài)分布且誤差量級較小,所以通過本文兩種方法能夠有效地反演出每個(gè)斷層上的滑動(dòng),且整體誤差較小,故通過模擬實(shí)驗(yàn)?zāi)軌蝌?yàn)證協(xié)同反演的有效性.從圖中兩種方法的殘差圖可以看出,在模擬滑動(dòng)區(qū)域,兩種方法都能夠較好地反演出實(shí)際滑動(dòng),僅在S1斷層的地表處顯示出部分滑動(dòng)低于實(shí)際滑動(dòng),這主要是平滑約束造成的.而在斷層的深處,兩種方法都出現(xiàn)了部分偽滑移,主要是因?yàn)楸疚脑O(shè)置的斷層相比于實(shí)際滑動(dòng)區(qū)域較大.

    圖3展示了兩種方法反演的滑動(dòng)分布模型正演的InSAR形變場以及模擬形變.對比前兩列,兩種方法反演的形變場與模擬形變場保持較好的一致性,最后一列為兩種方法的殘差圖,從殘差的形變尺度上可以看出,最大形變誤差在3 cm以內(nèi),表明了兩種反演方法的精度是非常高的.

    圖3 InSAR觀測形變值與模擬形變值對比及殘差圖(a)和(d) 表示InSAR觀測數(shù)據(jù); (b)和(d) 分別為采用聯(lián)合反演和協(xié)同反演方法得到的模擬值; (c)和(f) 分別表示兩種方法反演的殘差分布.Fig.3 The observed InSAR displacements, the modeled displacements and their residuals(a,d) represent InSAR displacements, respectively; (b,d) The modeled displacements obtained by the joint inversion and the cooperation inversion method, respectively; (c,f) represent the residuals obtained by the joint inversion and the cooperation inversion method, respectively.

    圖4展示了前震GPS形變圖的模擬形變量和采用方差分量估計(jì)法和協(xié)同反演方法正演的GPS形變場;圖5展示了主震GPS形變圖的模擬形變量和兩種方法正演的GPS形變場;從圖4和圖5的結(jié)果可以看出,本文兩種反演的模型的正演形變場與模擬形變場非常吻合,由于模擬實(shí)驗(yàn)中所加入的誤差是正態(tài)分布且誤差較小,所以反演模型與模擬值較為一致.

    圖4 前震GPS觀測形變值與模擬形變值及殘差圖(a) 藍(lán)色箭頭和紅色箭頭分別表示觀測的和聯(lián)合反演模擬的GPS水平位移; (b) 藍(lán)色箭頭表示聯(lián)合反演模擬的GPS水平位移的殘差; (c) 藍(lán)色箭頭和紅色箭頭分別表示觀測的和協(xié)同反演模擬的GPS水平位移; (d) 藍(lán)色箭頭表示協(xié)同反演模擬的GPS水平位移的殘差; (a)和(c) 中每個(gè)測站的彩色圓圈表示測得的(外圈)和預(yù)測的(內(nèi)圈)GPS垂直位移; (b)和(d)中每個(gè)測站的彩色圓圈表示垂直形變的殘差.Fig.4 The foreshock GPS observation displacements, the modeled displacements and their residuals(a) The results with blue and red arrows represent the observed and simulated GPS horizontal displacements by joint inversion method, respectively; (b) The results with the blue arrow represent the residuals of the GPS horizontal displacements; (c) The results with the blue and red arrow represent the observed and simulated GPS horizontal displacements of the cooperation inversion method, respectively; (d) The results with the blue arrow represent the residuals of the GPS horizontal displacements; The colored circles at each station in (a) and (c) denotes measured (outer circle) and predicted (inner circle) GPS vertical displacements. The colored circles at each station in (b) and (d) denote the residuals of vertical displacements.

    圖5 主震GPS觀測形變值與模擬形變值及殘差圖(a) 藍(lán)色箭頭和紅色箭頭分別表示觀測的和聯(lián)合反演模擬的GPS水平位移; (b) 藍(lán)色箭頭表示聯(lián)合反演模擬的GPS水平位移的殘差; (c) 藍(lán)色箭頭和紅色箭頭分別表示觀測的和協(xié)同反演模擬的GPS水平位移; (d) 藍(lán)色箭頭表示協(xié)同反演模擬的GPS水平位移的殘差; (a)和(c)中每個(gè)測站的彩色圓圈表示測得的(外圈)和預(yù)測的(內(nèi)圈)GPS垂直位移; (b)和(d)中每個(gè)測站的彩色圓圈表示垂直形變的殘差.Fig.5 The main shock GPS observation displacements, the modeled displacements and their residuals(a) The results with the blue and red arrow represent the observed and simulated horizontal GPS displacements by joint inversion method, respectively; (b) The results with the blue arrow represent the residuals of the GPS horizontal displacements; (c) The results with the blue and red arrow are the observed and simulated GPS horizontal displacement by cooperation inversion method, respectively; (d) The blue arrow represents the residuals of the horizontal GPS displacements; The colored circles at each station in (a) and (c) denote the measured (outer circle) and predicted (inner circle) GPS vertical displacements; The colored circles at each station in (b) and (d) denote the residual of vertical displacements.

    綜上所述,采用協(xié)同反演方法能夠有效地反演出多斷層序列地震的滑動(dòng)分布,且反演精度較高,能夠應(yīng)用于地震序列滑動(dòng)分布反演研究中.

    3 Ridgecrest地震序列

    2019年7月4日,加利福尼亞州當(dāng)?shù)貢r(shí)間17時(shí)33分,在加州東部西爾斯谷西南發(fā)生MW6.4地震(35.705°N,117.504°W),在此之后發(fā)生了一系列的余震,34 h后發(fā)生了MW7.1地震(35.770°N,117.599°W).在Ridgecrest地震序列中,MW6.4地震視為前震,MW7.1地震為主震.該地震序列被大量大地測量數(shù)據(jù)記錄到,為研究該地震序列提供了數(shù)據(jù)支持.本文從Li等(2020)中獲取Ridgecrest地震序列76個(gè)GPS點(diǎn)位形變數(shù)據(jù),其中前震與主震的形變被單獨(dú)記錄到;從Xu等(2020)中獲取Ridgecrest地震序列T65升軌數(shù)據(jù)1908個(gè),T66升軌數(shù)據(jù)1976個(gè).根據(jù)已有文獻(xiàn)公布的斷層參數(shù)信息(Li et al., 2020; Wang et al., 2020),我們設(shè)置了部分?jǐn)鄬訁?shù)如表3所示,根據(jù)斷層參數(shù)構(gòu)建斷層面,將斷層面均勻剖分成2 km×2 km的矩形單元,并利用協(xié)同反演方法反演Ridgecrest地震序列滑動(dòng)分布,其中正則化參數(shù)采用方差分量估計(jì)確定,分別為0.0761和0.0339.根據(jù)圖6擬合殘差和平滑度之間的折中曲線,我們得到閾值大小為0.251 m.

    表3 Ridgecrest地震序列的滑動(dòng)分布參數(shù)Table 3 Fault geometry and source parameters of the Ridgecrest earthquake sequence

    表4 Ridgecrest地震序列主震與前震的滑動(dòng)分布參數(shù)Table 4 The slip parameters of the main shock and foreshock of the Ridgecrest earthquake sequence

    反演的滑動(dòng)分布模型與GPS觀測數(shù)據(jù)擬合較好,其中MW6.4前震GPS數(shù)據(jù)的水平和垂直方向的擬合殘差分別為1.7 mm和6.2 mm;MW7.1主震GPS數(shù)據(jù)的水平和垂直方向的擬合殘差分別為5.7 mm和10.2 mm.圖8、9分別展示了前震和主震GPS觀測數(shù)據(jù)與反演模型的正演形變.對于InSAR數(shù)據(jù)來說,擬合殘差較大,尤其在斷層線附近,誤差高達(dá)數(shù)十厘米,Li等(2020)表明這些誤差的來源可能是模型中表面破裂軌跡的簡化、衛(wèi)星軌道誤差、大氣延遲和非彈性形變等.圖10展示了T065和T066軌道的InSAR形變、反演模型的正演形變以及殘差大小.根據(jù)反演的滑動(dòng)模型,主震破裂了斷層F1和F3,且F1斷層的破裂幾乎是純右旋走滑的,F(xiàn)3斷層以右旋走滑為主伴有少量的傾滑.F1斷層出現(xiàn)兩個(gè)主滑動(dòng)區(qū),其中最大滑動(dòng)量達(dá)到4.26 m.主震破裂高滑動(dòng)區(qū)在地下10 km以上區(qū)域,隨著深度的增加,滑動(dòng)逐漸減少.根據(jù)本文反演的滑動(dòng)分布模型,計(jì)算出MW7.1主震釋放的地震矩為4.49×1019N·m,假設(shè)剪切模量為30 GPa,對應(yīng)于MW7.07,與USGS公布結(jié)果較為一致.MW6.4前震破裂了F2和F3斷層,且這兩個(gè)斷層具有共軛斷層結(jié)構(gòu)的特征,兩條斷層幾乎相互垂直.MW6.4前震的滑動(dòng)分布如圖7所示,主要破裂發(fā)生在F2斷層,表現(xiàn)為幾乎純左旋走滑特征,其最大滑動(dòng)量達(dá)到1 m,而F3斷層對應(yīng)滑動(dòng)較小.針對F3斷層,MW6.4主震破裂區(qū)域離F2斷層更近,且破裂深度更深.前震反演的滑動(dòng)分布模型對應(yīng)總的地震矩為5.28×1018N·m,對應(yīng)于MW6.45,略高于USGS給出的結(jié)果.Li等(2020)在利用GPS數(shù)據(jù)進(jìn)行單獨(dú)反演時(shí),反演的F2斷層滑動(dòng)量較小,在進(jìn)行聯(lián)合反演時(shí)F2斷層上的滑動(dòng)量有所提高.而聯(lián)合反演后F3斷層的滑動(dòng)為兩次地震的疊加,其最大滑動(dòng)量達(dá)到~3.95 m,相較于本文結(jié)果來說偏大.總體而言,Li等(2020)在僅利用GPS數(shù)據(jù)進(jìn)行反演時(shí),由于GPS數(shù)據(jù)大多分布在遠(yuǎn)場,約束力不足,其反演結(jié)果較為平滑;而加入InSAR數(shù)據(jù)后,模型得到很好的約束,模型反演的最大滑動(dòng)量達(dá)到~5.5m,累積釋放地震矩為4.93×1019,約為MW7.1.Wang等(2020)利用迭代的方式分離了主震與前震的InSAR形變場,并聯(lián)合GPS數(shù)據(jù)反演了前震與主震的滑動(dòng)分布.滑動(dòng)分布結(jié)果顯示,F(xiàn)2斷層最大滑動(dòng)量達(dá)到2 m,F(xiàn)1斷層最大滑動(dòng)量達(dá)到1 m,釋放地震矩為5.12×1018,約為MW6.44.對于MW7.1主震,本文滑動(dòng)分布結(jié)果與Wang等(2020)較為相似,存在兩個(gè)較大滑動(dòng)區(qū)域,最大滑動(dòng)量達(dá)到~5m,略大于本文結(jié)果.Goldberg等(2020)提出了一種新的運(yùn)動(dòng)學(xué)滑動(dòng)分布反演方法,有效地反演出兩次地震的滑動(dòng)分布情況,并分析了兩次地震之間的靜、動(dòng)態(tài)庫倫應(yīng)力觸發(fā)關(guān)系.Goldberg等(2020)反演結(jié)果表明,MW6.4前震最大滑動(dòng)量達(dá)到~1 m,釋放地震矩4.37×1018,約為MW6.36;MW7.1主震最大滑動(dòng)量達(dá)到~4 m,釋放地震矩4.41×1019,約為MW7.03.Liu等(2019)聯(lián)合大地測量和地震波數(shù)據(jù)反演了2019 Ridgecrest地震序列,其反演結(jié)果顯示MW6.4前震最大滑動(dòng)量達(dá)到1.1 m,釋放地震矩5.05×1018,約為MW6.4;MW7.1主震最大滑動(dòng)量達(dá)到3.7 m,釋放地震矩5.0×1019,約為MW7.03.以上文獻(xiàn)采用不同的觀測數(shù)據(jù)以及不同的方法對2019 Ridgecrest地震序列進(jìn)行了同震滑動(dòng)分布反演,雖然滑動(dòng)分布存在差異,但其最大破裂區(qū)域較為一致;盡管不同學(xué)者設(shè)置了不同的斷層數(shù)據(jù)進(jìn)行反演,最終反演矩震級與USGS公布結(jié)果較為一致.

    圖6 確定迭代閾值的折中曲線圖Fig.6 Graph of compromise curve for determining iteration threshold

    圖7 MW6.4前震和MW7.1主震滑動(dòng)分布結(jié)果(a)和(b)分別表示MW6.4前震F2和F3斷層的滑移分布; (c)和(d)分別表示MW7.1主震F1和F3斷層的滑移分布.Fig.7 The slip distribution of the MW6.4 foreshock and the MW7.1 main shock(a) and (b) represent the slip distribution of faults F2 and F3 for the MW6.4 foreshock, respectively; (c) and (d) represent the slip distribution of the F1 and F3 faults of the MW7.1 main shock, respectively.

    圖8 MW6.4前震GPS形變場及其殘差(a) 中藍(lán)色和黃色箭頭分別表示GPS水平方向觀測值和模擬值,每個(gè)測站的彩色圓圈表示垂直方向觀測值(外圈)和模擬值(內(nèi)圈); (b) 中藍(lán)色箭頭表示水平方向殘差,彩色圓圈表示垂直方向殘差.Fig.8 GPS deformation measurements of MW6.4 foreshock and its residuals(a) The results with the blue and yellow arrows represent the observed and simulated GPS horizontal displacements, respectively, and colored circles at each station represent the observed (outer circle) and predicted (inner circle) GPS vertical displacements; The blue arrow in (b) denotes the residual in the horizontal direction, and the colored circle indicates the residuals in the vertical direction. The black line represents the fault traces.

    圖9 MW7.1主震GPS形變場及其殘差(a) 藍(lán)色和青色箭頭分別表示不同比例GPS水平方向觀測值,黃色和紅色箭頭分別表示不同比例模擬值,每個(gè)測站的彩色圓圈表示垂直方向觀測值(外圈)和模擬值(內(nèi)圈); (b) 藍(lán)色箭頭表示水平方向殘差,彩色圓圈表示垂直方向殘差.Fig.9 GPS deformation measurements of MW7.1 main shock and its residuals(a) The results with the blue and cyan-blue arrows represent different scales of GPS horizontal displacements, the results with the yellow and red arrows represent different scales of predicted displacements, and the colored circles at each station denote measured (outer circle) and predicted (inner circle) vertical displacements; (b) The results with the blue arrow represent the residuals in the horizontal direction, and the results with the colored circle represent the residuals in the vertical direction. The black line represents the fault traces.

    圖10 InSAR形變值、模型值以及殘差(a)和(d)分別為ALOS-2上升軌道T065和上升軌道T066的LOS形變場; (b)和(e)分別為利用協(xié)同反演方法得到的模擬形變場; (c)和(f)分別為對應(yīng)軌道的殘差;黑線表示斷層跡線.Fig.10 The observed InSAR displacements, the modeled displacements and their residuals(a) and (d) are the LOS deformation fields of the ALOS-2 ascending track T065 and T066, respectively; (b) and (e) are the predicted displacements obtained by the cooperation inversion method, respectively; (c) and (f) are the residuals of the corresponding tracks, respectively; The black line represents the fault traces.

    綜上所述,采用協(xié)同反演方法可以有效地分離出該地震序列InSAR同震形變,并聯(lián)合GPS數(shù)據(jù)反演了該地震序列的滑動(dòng)分布,與已有研究進(jìn)行了對比分析,證明了本文方法在同震滑動(dòng)分布反演中的有效性.

    4 結(jié)論

    本文采用的協(xié)同反演方法為反演序列地震滑動(dòng)分布提供了一種新的思路.本文通過模擬實(shí)驗(yàn),利用聯(lián)合反演方法和協(xié)同反演方法進(jìn)行序列地震滑動(dòng)分布反演,驗(yàn)證了本文方法的有效性;此外本文方法無需確定多類數(shù)據(jù)的相對權(quán)比,針對序列地震InSAR數(shù)據(jù)包含多次地震形變,利用本文方法能夠有效地反演出單次地震事件的滑動(dòng)分布.并將本文方法應(yīng)用于Ridgecrest地震序列滑動(dòng)分布反演,結(jié)果表明:MW7.1主震滑動(dòng)主要集中在0~10 km范圍,最大滑動(dòng)量達(dá)到4.26 m,位于地表下0~5 km.累積釋放地震矩4.49×1019N·m,對應(yīng)矩震級約為MW7.07.MW6.4前震破裂了一對共軛斷層,其滑動(dòng)主要集中在0~10 km范圍,且東北向斷層破裂更大,最大滑動(dòng)量達(dá)到1 m,位于地表下2~4 km處.累積釋放地震矩5.28×1018N·m,對應(yīng)矩震級約為MW6.45.反演結(jié)果均與USGS公布結(jié)果較為一致.根據(jù)模擬實(shí)驗(yàn)與Ridgecrest地震序列反演結(jié)果表明,采用協(xié)同反演方法能夠針對InSAR無法分離出單次地震同震形變場的情況下,利用GPS數(shù)據(jù)反演的結(jié)果作為先驗(yàn)信息,從而分離出單次地震的InSAR形變數(shù)據(jù),使InSAR數(shù)據(jù)能夠有效用于單次地震的滑動(dòng)分布反演中,一定程度上突破了InSAR數(shù)據(jù)在時(shí)間上的局限性.

    致謝感謝審稿專家對本文提出的寶貴意見.文中部分圖件是采用開源軟件GMT(Generic Mapping Tools)繪制,在此表示感謝.感謝許光煜博士提出的寶貴建議和幫助.

    猜你喜歡
    主震殘差滑動(dòng)
    基于雙向GRU與殘差擬合的車輛跟馳建模
    基于殘差學(xué)習(xí)的自適應(yīng)無人機(jī)目標(biāo)跟蹤算法
    基于遞歸殘差網(wǎng)絡(luò)的圖像超分辨率重建
    一種新型滑動(dòng)叉拉花鍵夾具
    Big Little lies: No One Is Perfect
    多塔斜拉橋在主震-余震序列波下地震位移研究
    平穩(wěn)自相關(guān)過程的殘差累積和控制圖
    河南科技(2015年8期)2015-03-11 16:23:52
    滑動(dòng)供電系統(tǒng)在城市軌道交通中的應(yīng)用
    龍卷流旋轉(zhuǎn)與地震成因
    一種基于變換域的滑動(dòng)聚束SAR調(diào)頻率估計(jì)方法
    大片免费播放器 马上看| 日本黄色片子视频| 校园人妻丝袜中文字幕| 777米奇影视久久| 少妇人妻 视频| 欧美bdsm另类| 最新中文字幕久久久久| 哪个播放器可以免费观看大片| 国产精品一区二区性色av| 亚洲欧美精品自产自拍| 成人亚洲精品一区在线观看 | 国产精品国产三级国产av玫瑰| 欧美精品一区二区大全| 中文字幕久久专区| 伊人久久国产一区二区| 日韩 亚洲 欧美在线| 日本黄色片子视频| 乱码一卡2卡4卡精品| 午夜免费观看性视频| 高清视频免费观看一区二区| 看免费成人av毛片| 亚洲欧美精品自产自拍| 亚洲欧美精品专区久久| 久久久成人免费电影| 久久99热这里只有精品18| 国产欧美另类精品又又久久亚洲欧美| 国产黄a三级三级三级人| 少妇裸体淫交视频免费看高清| 国产免费一区二区三区四区乱码| 国产片特级美女逼逼视频| 亚洲国产最新在线播放| 一个人观看的视频www高清免费观看| 欧美日韩视频高清一区二区三区二| 一区二区av电影网| 亚洲av成人精品一二三区| 草草在线视频免费看| 亚洲精品成人久久久久久| 婷婷色麻豆天堂久久| 免费看不卡的av| 18+在线观看网站| 午夜视频国产福利| 联通29元200g的流量卡| 别揉我奶头 嗯啊视频| 免费黄频网站在线观看国产| av在线亚洲专区| 成人免费观看视频高清| 精品国产一区二区三区久久久樱花 | 最新中文字幕久久久久| av又黄又爽大尺度在线免费看| 女的被弄到高潮叫床怎么办| 国产在线一区二区三区精| 欧美潮喷喷水| 天天一区二区日本电影三级| 欧美xxⅹ黑人| 五月伊人婷婷丁香| 一级a做视频免费观看| 成人无遮挡网站| 免费在线观看成人毛片| 在线观看人妻少妇| 在线精品无人区一区二区三 | 寂寞人妻少妇视频99o| 永久网站在线| 国产爱豆传媒在线观看| 成人国产麻豆网| 狂野欧美激情性xxxx在线观看| 一级毛片我不卡| 成人二区视频| 午夜精品一区二区三区免费看| 人人妻人人爽人人添夜夜欢视频 | 欧美zozozo另类| 久久99热6这里只有精品| 欧美潮喷喷水| 夜夜看夜夜爽夜夜摸| 18禁动态无遮挡网站| 午夜视频国产福利| 免费观看av网站的网址| 亚洲成人中文字幕在线播放| 亚洲国产欧美人成| 美女cb高潮喷水在线观看| 国产一级毛片在线| 日韩电影二区| .国产精品久久| 日韩欧美精品v在线| 午夜激情久久久久久久| 卡戴珊不雅视频在线播放| 国产精品一区www在线观看| 中文字幕亚洲精品专区| 国产免费又黄又爽又色| 中文资源天堂在线| 国产有黄有色有爽视频| 一区二区三区精品91| 日韩三级伦理在线观看| av在线播放精品| 日韩 亚洲 欧美在线| 又爽又黄a免费视频| 91狼人影院| 22中文网久久字幕| 99久久精品热视频| 校园人妻丝袜中文字幕| 亚洲欧美一区二区三区国产| 在线亚洲精品国产二区图片欧美 | 国产女主播在线喷水免费视频网站| 少妇被粗大猛烈的视频| 日本猛色少妇xxxxx猛交久久| 久久6这里有精品| 人妻少妇偷人精品九色| 中文字幕人妻熟人妻熟丝袜美| 亚洲av成人精品一区久久| 性色avwww在线观看| 亚洲成人久久爱视频| 狂野欧美白嫩少妇大欣赏| 亚洲av免费高清在线观看| 老师上课跳d突然被开到最大视频| 国产淫语在线视频| 亚洲av日韩在线播放| 18禁裸乳无遮挡动漫免费视频 | 国产精品久久久久久精品古装| 欧美人与善性xxx| 深夜a级毛片| 日韩av免费高清视频| 王馨瑶露胸无遮挡在线观看| 欧美3d第一页| 国产成人午夜福利电影在线观看| 欧美日韩一区二区视频在线观看视频在线 | h日本视频在线播放| 国产一区二区在线观看日韩| 国产av国产精品国产| 狠狠精品人妻久久久久久综合| 欧美潮喷喷水| 国产精品一区二区性色av| 亚洲自偷自拍三级| 最近手机中文字幕大全| 青春草国产在线视频| 99久久精品一区二区三区| 成人毛片60女人毛片免费| 啦啦啦在线观看免费高清www| 久久午夜福利片| 亚洲最大成人手机在线| 国产 一区 欧美 日韩| 欧美人与善性xxx| 神马国产精品三级电影在线观看| 两个人的视频大全免费| 国产高潮美女av| 亚州av有码| 免费观看的影片在线观看| 亚洲精品国产色婷婷电影| 国产免费一区二区三区四区乱码| 热re99久久精品国产66热6| 成人午夜精彩视频在线观看| 三级经典国产精品| 少妇被粗大猛烈的视频| 国产男人的电影天堂91| 99久国产av精品国产电影| 精品人妻熟女av久视频| 亚洲国产精品国产精品| 国产白丝娇喘喷水9色精品| av天堂中文字幕网| 2021少妇久久久久久久久久久| 美女高潮的动态| 亚洲av成人精品一区久久| 精品国产一区二区三区久久久樱花 | 欧美日韩亚洲高清精品| 波多野结衣巨乳人妻| 久久热精品热| 亚洲内射少妇av| 男插女下体视频免费在线播放| 男人狂女人下面高潮的视频| 国产在线男女| 你懂的网址亚洲精品在线观看| 精品国产乱码久久久久久小说| 丰满乱子伦码专区| 性插视频无遮挡在线免费观看| 亚洲色图综合在线观看| 久久精品国产亚洲网站| 联通29元200g的流量卡| 午夜福利在线观看免费完整高清在| 日本与韩国留学比较| 日韩中字成人| 有码 亚洲区| 好男人视频免费观看在线| 国产69精品久久久久777片| 精品久久久噜噜| 女人久久www免费人成看片| 亚洲av中文av极速乱| 亚洲真实伦在线观看| 一级片'在线观看视频| 国产免费一区二区三区四区乱码| 国产精品久久久久久精品古装| 亚洲av在线观看美女高潮| 青春草视频在线免费观看| 黄色配什么色好看| 亚洲精品乱码久久久久久按摩| 神马国产精品三级电影在线观看| 九九久久精品国产亚洲av麻豆| 夫妻午夜视频| 丰满少妇做爰视频| 欧美xxxx性猛交bbbb| 观看免费一级毛片| 国产成人午夜福利电影在线观看| 麻豆成人av视频| 伦精品一区二区三区| 亚洲成人av在线免费| 国产成人精品久久久久久| 国产永久视频网站| 日产精品乱码卡一卡2卡三| 亚洲成人中文字幕在线播放| 亚洲伊人久久精品综合| 中国国产av一级| 校园人妻丝袜中文字幕| 最后的刺客免费高清国语| 99热全是精品| 久久亚洲国产成人精品v| 伦理电影大哥的女人| 欧美区成人在线视频| 青春草亚洲视频在线观看| 午夜福利视频1000在线观看| 真实男女啪啪啪动态图| 亚洲在久久综合| 成年版毛片免费区| 国产精品久久久久久久电影| 国内精品美女久久久久久| 男女边摸边吃奶| 精品国产乱码久久久久久小说| 久久精品综合一区二区三区| 一级毛片我不卡| 亚洲久久久久久中文字幕| 日韩成人伦理影院| 国产永久视频网站| 日韩制服骚丝袜av| 校园人妻丝袜中文字幕| 99热全是精品| 久久这里有精品视频免费| 午夜福利高清视频| 日韩伦理黄色片| 激情五月婷婷亚洲| 人妻少妇偷人精品九色| 国产亚洲av片在线观看秒播厂| 尾随美女入室| 久久久久久久久久久免费av| 久久久久精品久久久久真实原创| 啦啦啦中文免费视频观看日本| 热re99久久精品国产66热6| 80岁老熟妇乱子伦牲交| 午夜福利网站1000一区二区三区| 日韩大片免费观看网站| 色视频在线一区二区三区| 狂野欧美激情性bbbbbb| 中文字幕人妻熟人妻熟丝袜美| 亚洲成色77777| 22中文网久久字幕| 美女国产视频在线观看| 日韩成人av中文字幕在线观看| 国产精品国产三级国产专区5o| 亚洲丝袜综合中文字幕| 人人妻人人爽人人添夜夜欢视频 | av免费在线看不卡| 久久99热这里只频精品6学生| 国产精品嫩草影院av在线观看| 亚洲人成网站在线播| 欧美人与善性xxx| av国产免费在线观看| 五月玫瑰六月丁香| 免费av毛片视频| 丰满少妇做爰视频| 人妻一区二区av| 男女下面进入的视频免费午夜| 欧美xxxx黑人xx丫x性爽| 涩涩av久久男人的天堂| 国内精品美女久久久久久| 亚洲最大成人av| 中文字幕亚洲精品专区| 天堂中文最新版在线下载 | 国语对白做爰xxxⅹ性视频网站| 国产色爽女视频免费观看| av播播在线观看一区| 国语对白做爰xxxⅹ性视频网站| 赤兔流量卡办理| 在线观看一区二区三区| 久久久久九九精品影院| 99久久精品一区二区三区| 日韩不卡一区二区三区视频在线| 日本一二三区视频观看| 免费看不卡的av| 亚洲天堂av无毛| 美女内射精品一级片tv| 又爽又黄a免费视频| 国产老妇伦熟女老妇高清| 国产精品蜜桃在线观看| 大香蕉久久网| 神马国产精品三级电影在线观看| 日本猛色少妇xxxxx猛交久久| 国内精品宾馆在线| 小蜜桃在线观看免费完整版高清| 亚洲精品,欧美精品| 偷拍熟女少妇极品色| 毛片一级片免费看久久久久| 欧美成人一区二区免费高清观看| 久热久热在线精品观看| 岛国毛片在线播放| 成人黄色视频免费在线看| 精品国产乱码久久久久久小说| 精品午夜福利在线看| 又粗又硬又长又爽又黄的视频| 亚洲最大成人手机在线| 男女无遮挡免费网站观看| 边亲边吃奶的免费视频| 欧美国产精品一级二级三级 | 日本色播在线视频| 久久久久精品久久久久真实原创| 亚洲国产欧美人成| 日韩免费高清中文字幕av| 国产精品久久久久久久久免| 午夜视频国产福利| 你懂的网址亚洲精品在线观看| av免费观看日本| 免费大片18禁| 秋霞在线观看毛片| 高清视频免费观看一区二区| 精华霜和精华液先用哪个| 2021少妇久久久久久久久久久| 欧美一级a爱片免费观看看| 欧美激情久久久久久爽电影| 亚洲欧美精品专区久久| 黄色怎么调成土黄色| 欧美成人午夜免费资源| 日本黄大片高清| 少妇的逼好多水| 精品人妻熟女av久视频| 看黄色毛片网站| 午夜福利网站1000一区二区三区| 免费看a级黄色片| 日韩制服骚丝袜av| 亚洲国产最新在线播放| 国产一区有黄有色的免费视频| 国内精品美女久久久久久| 成人特级av手机在线观看| 亚洲av不卡在线观看| 国产久久久一区二区三区| 97超碰精品成人国产| 亚洲欧美精品专区久久| 六月丁香七月| 久久久久久久久久久免费av| 国产精品一区二区在线观看99| 色吧在线观看| 亚洲欧美中文字幕日韩二区| www.av在线官网国产| 国产欧美另类精品又又久久亚洲欧美| 毛片女人毛片| 99久久人妻综合| 成人亚洲精品一区在线观看 | 久久精品国产鲁丝片午夜精品| www.av在线官网国产| 中文欧美无线码| a级毛片免费高清观看在线播放| 男人爽女人下面视频在线观看| 中文字幕人妻熟人妻熟丝袜美| 国产白丝娇喘喷水9色精品| 精华霜和精华液先用哪个| 偷拍熟女少妇极品色| 三级国产精品片| 国产精品不卡视频一区二区| 久久热精品热| 男女边吃奶边做爰视频| 国内精品宾馆在线| tube8黄色片| 天堂中文最新版在线下载 | 嫩草影院入口| 老师上课跳d突然被开到最大视频| 欧美最新免费一区二区三区| 成人国产麻豆网| 午夜老司机福利剧场| 男女无遮挡免费网站观看| 国产精品无大码| 亚洲精品国产av成人精品| 看黄色毛片网站| 99久久中文字幕三级久久日本| 精品一区二区三区视频在线| 国产午夜精品一二区理论片| 在线看a的网站| 大片电影免费在线观看免费| 丝袜脚勾引网站| 色5月婷婷丁香| 搡女人真爽免费视频火全软件| 嫩草影院精品99| 狂野欧美白嫩少妇大欣赏| a级一级毛片免费在线观看| 午夜福利高清视频| 一级毛片电影观看| 久久久久性生活片| 91精品一卡2卡3卡4卡| 国产亚洲5aaaaa淫片| 狂野欧美激情性xxxx在线观看| 99热全是精品| 国产精品99久久久久久久久| 免费观看a级毛片全部| 午夜福利视频精品| 久久久久久久久久久丰满| 国内精品美女久久久久久| 欧美日韩视频高清一区二区三区二| 日日啪夜夜爽| 男女边摸边吃奶| 一区二区三区免费毛片| 在线播放无遮挡| 高清日韩中文字幕在线| 在线免费观看不下载黄p国产| 国产精品人妻久久久久久| 日韩av不卡免费在线播放| 18禁在线播放成人免费| 国产一级毛片在线| 女人久久www免费人成看片| 国产一区亚洲一区在线观看| 少妇熟女欧美另类| 成年女人看的毛片在线观看| 日韩强制内射视频| 又爽又黄a免费视频| 亚洲精品日韩av片在线观看| 在线免费十八禁| 一个人看的www免费观看视频| 国产高清国产精品国产三级 | 蜜桃亚洲精品一区二区三区| 我的老师免费观看完整版| 少妇的逼好多水| 男人狂女人下面高潮的视频| 插逼视频在线观看| 免费观看无遮挡的男女| 国产欧美日韩精品一区二区| 久久国产乱子免费精品| 国产高清有码在线观看视频| 久久精品国产亚洲av天美| 亚洲欧美精品自产自拍| 视频中文字幕在线观看| 大又大粗又爽又黄少妇毛片口| 国产女主播在线喷水免费视频网站| 超碰av人人做人人爽久久| 国产日韩欧美在线精品| 久久精品久久久久久噜噜老黄| 男人爽女人下面视频在线观看| 黄色配什么色好看| 18禁裸乳无遮挡动漫免费视频 | 黄色一级大片看看| 久久精品综合一区二区三区| 亚洲人成网站在线播| 成年女人看的毛片在线观看| 国产亚洲最大av| 直男gayav资源| 人妻夜夜爽99麻豆av| 又粗又硬又长又爽又黄的视频| 国产高清有码在线观看视频| 国产综合懂色| 99热全是精品| 国产伦在线观看视频一区| 精华霜和精华液先用哪个| 一级毛片我不卡| 80岁老熟妇乱子伦牲交| 婷婷色综合大香蕉| 国产亚洲av嫩草精品影院| 联通29元200g的流量卡| 深夜a级毛片| 亚洲精品国产av成人精品| 看黄色毛片网站| 日韩精品有码人妻一区| 国产精品熟女久久久久浪| 国产伦理片在线播放av一区| 99久久精品热视频| 最近2019中文字幕mv第一页| 岛国毛片在线播放| 一级片'在线观看视频| 又爽又黄无遮挡网站| 一级毛片黄色毛片免费观看视频| 亚洲精品国产成人久久av| 熟妇人妻不卡中文字幕| 国产爱豆传媒在线观看| 日日摸夜夜添夜夜爱| 精品一区二区三卡| 一级毛片黄色毛片免费观看视频| 亚洲国产精品专区欧美| 2021少妇久久久久久久久久久| 中国美白少妇内射xxxbb| 日本色播在线视频| 啦啦啦中文免费视频观看日本| 久久久久久久国产电影| 91久久精品电影网| av又黄又爽大尺度在线免费看| 99热6这里只有精品| 插逼视频在线观看| 卡戴珊不雅视频在线播放| 国产精品不卡视频一区二区| 五月玫瑰六月丁香| 久久久久久久久久久丰满| 少妇人妻一区二区三区视频| 久久久久网色| 欧美激情久久久久久爽电影| av在线天堂中文字幕| 大又大粗又爽又黄少妇毛片口| 菩萨蛮人人尽说江南好唐韦庄| 免费看日本二区| 亚洲av国产av综合av卡| av女优亚洲男人天堂| 全区人妻精品视频| 国产免费视频播放在线视频| 麻豆成人av视频| 成人漫画全彩无遮挡| 国产色婷婷99| 国产午夜精品久久久久久一区二区三区| 人妻少妇偷人精品九色| 亚洲欧美一区二区三区黑人 | 亚洲婷婷狠狠爱综合网| 日韩,欧美,国产一区二区三区| 听说在线观看完整版免费高清| 午夜视频国产福利| 久久久久久国产a免费观看| 国产欧美日韩精品一区二区| 九色成人免费人妻av| 黄色欧美视频在线观看| 97超碰精品成人国产| 特级一级黄色大片| 大片免费播放器 马上看| 国产精品久久久久久久电影| 日韩亚洲欧美综合| 大陆偷拍与自拍| 日韩欧美精品免费久久| 国产高清不卡午夜福利| 最近手机中文字幕大全| 天堂中文最新版在线下载 | 国产亚洲av片在线观看秒播厂| 日韩大片免费观看网站| 国产色爽女视频免费观看| 亚洲av福利一区| 成人美女网站在线观看视频| 午夜激情久久久久久久| 一级爰片在线观看| 五月玫瑰六月丁香| 内地一区二区视频在线| 99久久精品热视频| 久久99热6这里只有精品| 街头女战士在线观看网站| 免费黄网站久久成人精品| 久久久久久久久大av| 国产亚洲av片在线观看秒播厂| 日韩制服骚丝袜av| 亚洲人成网站高清观看| 精品午夜福利在线看| 午夜视频国产福利| 黄色怎么调成土黄色| 一级毛片aaaaaa免费看小| 2018国产大陆天天弄谢| a级一级毛片免费在线观看| 综合色av麻豆| 欧美日本视频| av免费观看日本| a级毛片免费高清观看在线播放| 国产日韩欧美在线精品| 免费不卡的大黄色大毛片视频在线观看| 亚洲不卡免费看| 日韩伦理黄色片| 亚洲精品国产av蜜桃| 免费观看a级毛片全部| 国产精品伦人一区二区| av在线观看视频网站免费| 精品人妻熟女av久视频| 日韩一本色道免费dvd| 亚洲成人久久爱视频| 国产永久视频网站| 国产美女午夜福利| 国产永久视频网站| 男人和女人高潮做爰伦理| 一级毛片黄色毛片免费观看视频| 天天躁日日操中文字幕| 亚洲在久久综合| 欧美性猛交╳xxx乱大交人| 亚洲图色成人| 午夜福利网站1000一区二区三区| 精品一区在线观看国产| 国产欧美日韩一区二区三区在线 | a级毛片免费高清观看在线播放| 久久精品国产a三级三级三级| 久久久午夜欧美精品| 男人和女人高潮做爰伦理| 成人高潮视频无遮挡免费网站| 高清日韩中文字幕在线| 国产黄色视频一区二区在线观看| 国产精品99久久久久久久久| 美女国产视频在线观看| 国产精品一区www在线观看| 一个人看视频在线观看www免费| 少妇人妻久久综合中文| 免费黄色在线免费观看| 大陆偷拍与自拍| 久久综合国产亚洲精品| 中文在线观看免费www的网站| 亚洲天堂av无毛| 2021少妇久久久久久久久久久| 国产精品国产av在线观看| 亚洲欧美成人综合另类久久久| 久久久色成人| 免费看光身美女| 国产男人的电影天堂91| 国产精品久久久久久精品古装| 精品国产露脸久久av麻豆| 嫩草影院精品99| 在线观看人妻少妇| 国产黄频视频在线观看| 一本久久精品| 中文资源天堂在线| 黑人高潮一二区| 国产爽快片一区二区三区| 只有这里有精品99| 久久精品久久久久久久性| 国产黄频视频在线观看| 女人被狂操c到高潮| 深夜a级毛片| 日韩电影二区| 欧美少妇被猛烈插入视频| 91精品一卡2卡3卡4卡| 观看美女的网站|