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

    基于VP/VS波速比模型約束的張渤地震帶深部電性結(jié)構(gòu)研究

    2021-08-03 10:56:56吳萍萍丁志峰譚捍東楊歧焱王興臣
    地球物理學(xué)報(bào) 2021年8期
    關(guān)鍵詞:波速斷裂帶電阻率

    吳萍萍, 丁志峰* , 譚捍東, 楊歧焱, 王興臣

    1 中國地震局地球物理研究所, 北京 100081 2 中國地質(zhì)大學(xué)(北京)地球物理與信息技術(shù)學(xué)院, 北京 100083 3 河北地質(zhì)大學(xué)地球科學(xué)學(xué)院, 石家莊 050031

    0 引言

    華北克拉通中新生代以來經(jīng)歷了獨(dú)特的巖石圈改造和破壞過程,并伴隨著廣泛的巖漿活動和大規(guī)模的構(gòu)造伸展盆地(翟明國等, 2005; 段永紅等, 2016),是中國大陸內(nèi)部地震活動強(qiáng)烈地區(qū)之一(王椿鏞等, 2017).地震活動分布圖(圖1)反映華北地區(qū)的地震活動具有明顯的分帶性,主要有四條重要的地震帶:張家口—渤海地震帶(簡稱張渤地震帶)、河北平原地震帶、汾渭地震帶和鄂爾多斯西緣地震帶.其中張渤地震帶位于華北平原的北部,整體呈現(xiàn)北西向的地震活動帶特征,該區(qū)域分布著一系列雁行排列的北西至北西西向斷裂(圖2),該斷裂帶與多條NNE-NE向斷裂帶交匯,M6.0以上的強(qiáng)震主要群集于NE向斷裂交匯地段.王椿鏞等(2017)統(tǒng)計(jì)該區(qū)域的發(fā)震深度主要集中在8~25 km范圍,在深度大于25km的地震數(shù)量急劇減少.深部介質(zhì)的電阻率對地球內(nèi)部溫度、流體、熔融/半熔融等異常反映靈敏,已有研究表明地震孕育和發(fā)生與地下高導(dǎo)體(例如火山高導(dǎo)巖漿囊、殼內(nèi)滑脫層、含水高導(dǎo)層等區(qū)域)和高阻體的接觸邊界息息相關(guān)(趙國澤等, 2001; 魏文博等, 2010; Zhang et al., 2016; Ye et al., 2018; 葉濤等, 2018).因此通過獲取張渤地震帶殼內(nèi)高導(dǎo)和高阻異常體的分布特征,可為認(rèn)識研究區(qū)的物質(zhì)狀態(tài)、地殼運(yùn)動過程、深部孕震環(huán)境等科學(xué)問題提供深部電性結(jié)構(gòu)依據(jù).

    大地電磁法以天然電磁場為場源,具有探測深度大、對低阻體反映靈敏、不受高阻屏蔽、設(shè)備輕便、成本低等優(yōu)勢,成為研究地球深部電性結(jié)構(gòu)的主要地球物理方法之一.華北克拉通也開展了大量的大地電磁探測研究工作(魏文博等, 2006, 2008, 2010; Ouro-Djobo等, 2018; 詹艷等, 2011; 徐光晶等, 2015),但因在張渤地震帶大地電磁測點(diǎn)數(shù)量少,深部電阻率結(jié)構(gòu)成像分辨率有限.同時(shí)大地電磁法因受限于方法本身的局限性,采集易受電磁干擾、體積效應(yīng)、縱向分辨隨深度增加而降低等因素影響,反演獲取的電阻率模型在深部電阻率經(jīng)常會出現(xiàn)明顯的拖尾現(xiàn)象,對于深部異常體的底界面形態(tài)的分辨率不足(金勝等, 2010).

    已有研究表明基于交叉梯度結(jié)構(gòu)耦合約束(Gallardo and Meju, 2003, 2004)的大地電磁和地震資料的聯(lián)合反演,可以實(shí)現(xiàn)兩種方法的優(yōu)勢互補(bǔ),克服單方法的局限性,提高反演精度(彭淼等, 2013; Bennington et al., 2015; Ogunbo et al., 2018; Peng et al., 2019; Wu et al., 2020).但地震數(shù)據(jù)和大地電磁數(shù)據(jù)的聯(lián)合反演需要考慮兩種不同數(shù)據(jù)的耦合性,在聯(lián)合反演過程中需要權(quán)衡多種權(quán)重參數(shù)(例如,交叉梯度項(xiàng)權(quán)重因子、正則化因子、數(shù)據(jù)權(quán)重因子等),在實(shí)際操作中存在困難.

    目前華北地區(qū)已有很多地震資料的研究成果,包括人工地震深反射地殼結(jié)構(gòu)(段永紅等,2016)、天然地震體波層析成像(黃金莉和趙大鵬, 2005; 楊婷等, 2012; 楊歧焱等, 2018)、地震面波和噪聲面波成像(魯來玉等, 2009; 何正勤等, 2009; Fang et al., 2010; 唐有彩等, 2011; 潘佳鐵等, 2011; 房立華等, 2013)、地殼厚度和泊松比(嵇少丞等, 2009; Lei, 2012)等.這些研究成果為認(rèn)識該區(qū)域深部速度結(jié)構(gòu)提供很好的素材.地下介質(zhì)速度和電阻率模型在一些特殊構(gòu)造環(huán)境(例如局部熔融半熔融異常體、含水的裂隙或斷裂、殼內(nèi)滑脫層等)在三維結(jié)構(gòu)分布上具有相似性(Newman et al., 2008; Berdichevsky et al., 2008),根據(jù)這一特征可以相互驗(yàn)證各自反演結(jié)果的可靠性,這也為本文基于地震學(xué)模型結(jié)構(gòu)耦合約束的大地電磁二維反演算法開發(fā)提供思路.

    彈性參數(shù)VP、VS、VP/VS、泊松比等屬性的空間分布特征對揭示地球內(nèi)部介質(zhì)的非均勻性有重要的意義,影響彈性參數(shù)的因素主要有礦物組成、溫度、壓力和空隙流體等.根據(jù)已有的實(shí)驗(yàn)資料(陳颙等, 2009)可知,不同的地震波屬性對介質(zhì)溫度和成分變化的反應(yīng)靈敏程度不同,例如當(dāng)溫度發(fā)生變化時(shí),S波速度對溫度異常的響應(yīng)能力比P波速度強(qiáng).VP/VS波速比模型(泊松比可以通過波速比換算得到)同時(shí)兼顧了介質(zhì)VP和VS屬性特點(diǎn),可以提供比單純P波和S波速度參數(shù)更為豐富的信息(李志偉等, 2009).在熔融/半熔融體中,波速比值隨著熔融程度的增加而增大(劉瓊林等, 2011),在這個(gè)地質(zhì)環(huán)境下電導(dǎo)率通常也表現(xiàn)出類似的特征.雖然很難定量分析波速比與電導(dǎo)率的物性關(guān)系式,但是可以根據(jù)這兩種物性在這些特殊環(huán)境中表現(xiàn)出來的結(jié)構(gòu)相似性進(jìn)行耦合約束反演,減小反演的多解性.

    因此,為了避免地震資料和大地電磁資料聯(lián)合反演過程中多參數(shù)的權(quán)衡,同時(shí)充分利用已有的地震學(xué)模型,本文提出將VP/VS波速比模型作為固定模型約束,通過交叉梯度結(jié)構(gòu)耦合函數(shù)帶入大地電磁二維反演算法中,實(shí)現(xiàn)基于地震學(xué)模型約束的大地電磁二維反演法,該方法即克服了聯(lián)合反演的局限,同時(shí)又可以有效克服大地電磁法深部電阻率分辨率不足的缺點(diǎn).文中首先介紹大地電磁二維OCCAM(Constable et al.,1987)反演算法中加入VP/VS波速比結(jié)構(gòu)耦合約束策略,通過理論模型合成數(shù)據(jù)檢驗(yàn)算法的可靠性.然后將該方法應(yīng)用于張渤地震帶大地電磁測深剖面,對比分析大地電磁法無約束和有約束反演結(jié)果,檢驗(yàn)約束反演算法在張渤地震帶應(yīng)用的優(yōu)越性和實(shí)用性.根據(jù)張渤地震帶反演結(jié)果設(shè)計(jì)一套理論模型,通過理論模型合成數(shù)據(jù)試算,檢驗(yàn)算法在該研究區(qū)應(yīng)用的有效性.最后結(jié)合該區(qū)域已有的地質(zhì)-地球物理資料,探討張渤地震帶深部電性結(jié)構(gòu)、孕震環(huán)境及其動力學(xué)機(jī)制等問題.

    1 地質(zhì)構(gòu)造背景

    圖1為華北地區(qū)大地構(gòu)造背景及地震分布圖,華北克拉通主要由東、中和西三個(gè)主要構(gòu)造區(qū)組成.西部為穩(wěn)定的鄂爾多斯高原,地震活動性主要集中在地塊的邊緣地區(qū),地塊內(nèi)部表現(xiàn)為穩(wěn)定性特征,在塊體西緣(鄂爾多斯地震帶,圖1中標(biāo)注為D的地震帶)曾在1739年發(fā)生寧夏平羅—銀川M8.0地震.中部地區(qū)為太行山造山帶,位于鄂爾多斯高原和華北平原之間,該區(qū)域構(gòu)造復(fù)雜,地震頻發(fā),曾發(fā)生M8.0以上地震3次(1303年山西洪洞地震、1556年陜西華縣地震、1695年山西臨汾地震).東部地區(qū)華北平原是研究華北克拉通巖石圈活化的熱點(diǎn)地區(qū)(Zhu et al., 2012; 翟明國, 2019),該區(qū)域內(nèi)地震分布呈現(xiàn)兩條共軛地震條帶(圖1),為北東向活動斷裂的河北平原地震帶(圖1中標(biāo)注為B的地震帶)和北西向活動斷裂的張渤地震帶(圖1中標(biāo)注為A的地震帶).本文研究區(qū)為張渤地震帶.

    圖1 華北地區(qū)大地構(gòu)造背景及地震分布圖A 張渤地震帶;B 河北平原地震帶;C 汾渭地震帶;D 鄂爾多斯地震帶. 紅色虛線表示主要斷裂分布.Fig.1 Geological background of the North China and the earthquakes distributionA Zhangbo seismic belt; B Hebei plain seismic belt; C Fenwei seismic belt; D Ordos seismic belt. Red lines donate the distribution of the main faults.

    張渤地震帶西起太行山,橫貫燕山和華北盆地,東入渤海,是中國東部地區(qū)一條重要的北西向地震活動帶,主要發(fā)育NNE和NE向斷裂,部分?jǐn)鄬幼呦驗(yàn)镹WW向.通過GPS資料分析,張渤地震帶以左旋走滑為主,兼有擠壓運(yùn)動.華北平原塊體和燕山塊體相對運(yùn)動時(shí)張渤帶左旋走滑的直接動力源(陳長云, 2016).

    本文大地電磁測線沿著該地震帶布設(shè)(圖2),自東向西經(jīng)過唐山、三河、北京、延慶到張家口地區(qū),主要穿過9條斷裂帶,依次為寧河—昌黎斷裂(F1)、唐山斷裂(F2)、豐田—野雞坨斷裂(F3)、夏墊斷裂(F4)、黃莊—高麗營斷裂(F5)、南口斷裂(F6)、延慶—礬山盆地北緣斷裂(F7)、懷(來)—涿(鹿)盆地北緣斷裂(F8)和張家口斷裂(F9).其中在黃莊—高麗營斷裂(F5)是一條京西山前大斷裂,是華北平原和京西隆起的分界線,走向NE20°~50°,總長約110 km,斷裂面傾向南東,傾角55°~75°.

    圖2 張渤地震帶構(gòu)造背景和MT測點(diǎn)分布圖F1 寧河—昌黎斷裂; F2 唐山斷裂; F3 豐臺—野雞坨斷裂; F4 夏墊斷裂; F5 黃莊—高麗營斷裂; F6 南口斷裂; F7 延慶—礬山盆地北緣斷裂; F8 懷(來)—涿(鹿)盆地北緣斷裂; F9 張家口斷裂. 紅色五角星: 三河平谷M 8.0地震; 藍(lán)色圓M 7.0~M 8.0之間地震; 黃色圓M 6.0~M 7.0之間地震; 綠色圓表示M 5.0~M6.0之間地震. 01~31數(shù)字為大地電磁測點(diǎn)編號.Fig.2 Geological settings of Zhangbo seismic belt and the locations of the MT sitesF1 Ninghe-Changli fault; F2 Tangshan fault; F3 Fengtai-Yejituo fault; F4 Xiadian fault; F5 Huangzhuang-Gaoliying fault; F6 Nankou fault; F7 North margin fault of Yanqing-Fanshan basin; F8 North margin fault of Huailai-Zhoulu basin; F9 Zhangjiakou fault. Red pentagram: Shanhe-Pinggu M8.0 earthquake; blue circles: M 7.0~M 8.0 earthquakes; yellow circles: M6.0~M7.0 earthquakes; green circles: M5.0 ~M6.0 eaqthquakes. The numbers from 01 to 31 are the number of the MT sites.

    測線最東邊為由寧河—昌黎斷裂(F1)、唐山斷裂(F2)、豐臺—野雞坨斷裂(F3)組成的唐山強(qiáng)震區(qū)斷裂構(gòu)造.唐山斷裂帶是1976年唐山7.8級地震的發(fā)震構(gòu)造,人工地震資料展示唐山斷裂帶走向NE30°,在地殼伴生復(fù)雜斷裂構(gòu)造帶,表現(xiàn)為花狀構(gòu)造斷裂特征,主干斷裂沿伸到Moho面(劉保金等, 2011a),屬走滑斷層性質(zhì).夏墊斷裂(F4)是一條NNE向巖石圈尺度的區(qū)域性深斷裂帶,全長28 km,總體走向45°,傾向SE,傾角70~80°,是北京平原歷史上唯一發(fā)生過8級地震的活動斷裂.

    測線最西側(cè)為張家口斷裂帶(F9),屬于首都圈地震重點(diǎn)監(jiān)視防御區(qū),為張家口-渤海地震構(gòu)造帶與山西地震構(gòu)造帶的交匯區(qū)(趙博等, 2011),由北西西向和近東西向兩組、 多段斷裂組成, 北西西向斷層為斷裂主體,傾向南東,傾角60°左右,全長約30 km,屬正斷層性質(zhì).懷涿盆地北緣斷裂(F8)控制著懷來—涿鹿盆地的發(fā)育,總體走向北東,傾向南東,傾角50~70°,是重要的控盆活動正斷裂之一.延慶—礬山盆地北緣斷裂(F7)是一條控制延慶—礬山盆地北緣的主邊界正斷裂帶,總體走向NE,傾向SE,傾角50~80°,全長約102 km.

    2 基于VP/VS波速比模型約束的二維大地電磁反演方法

    2.1 交叉梯度項(xiàng)

    Gallardo和Meju(2003, 2004)通過定義兩種不同物性參數(shù)的梯度的叉乘來定量描述兩種物性的結(jié)構(gòu)相似性,這是本文實(shí)現(xiàn)基于VP/VS波速比模型約束的大地電磁二維反演算法的關(guān)鍵.根據(jù)定義可將電阻率和波速比模型的交叉梯度函數(shù)表達(dá)式寫成:

    (1)

    mr表示電阻率模型,ms表示VP/VS模型.根據(jù)求解旋度公式可將交叉梯度項(xiàng)分解為

    (2)

    (3)

    (4)

    ty=Wcg,s·mr.

    (5)

    Wcg,s為與網(wǎng)格和ms相關(guān)的常數(shù)矩陣. 當(dāng)電阻率反演初始模型mr0為均勻半空間模型時(shí),則Wcg,s·mr0=0,那么公式(5)可以改寫成:

    ty=Wcg,s·(mr-mr0).

    (6)

    取公式(6)的平方和構(gòu)建VP/VS波速比模型的結(jié)構(gòu)耦合約束項(xiàng),則公式可寫成:

    (7)

    2.2 基于VP/VS波速比模型約束的二維大地電磁反演目標(biāo)函數(shù)

    在大地電磁無約束反演目標(biāo)函數(shù)中,通過權(quán)重因子η將VP/VS波速比模型的結(jié)構(gòu)耦合約束項(xiàng)加入到無約束反演目標(biāo)函數(shù)中,形成基于VP/VS波速比結(jié)構(gòu)耦合約束的二維大地電磁反演目標(biāo)函數(shù),具體表示如下:

    (8)

    f(·)是二維大地電磁正演響應(yīng)函數(shù),λ為拉格朗日因子,Cd為數(shù)據(jù)協(xié)方差矩陣,Wm為模型光滑矩陣,dobs為觀測數(shù)據(jù)矩陣.對目標(biāo)函數(shù)φ求梯度,則其梯度表達(dá)式可以寫成:

    (9)

    構(gòu)建

    (10)

    根據(jù)OCCAM反演算法(Constable et al., 1987)的迭代公式,即可獲得基于VP/VS波速比模型約束反演第k次模型更新公式:

    (11)

    公式中λ因子的選取采用自動搜索算法實(shí)現(xiàn),每一次迭代都要求根據(jù)RMS值選取最優(yōu)λ值,具體算法可參考Constable等(1987).圖3為本文開發(fā)的基于VP/VS波速比模型約束的大地電磁二維反演算法流程圖.

    圖3 基于VP/VS波速比模型約束的大地電磁二維反演流程圖Fig.3 The flow chart of 2-D MT inversion algorithm with the VP/VS wave velocity ratio model constraints

    3 理論模型合成數(shù)據(jù)反演試算

    3.1 理論模型設(shè)計(jì)

    為了檢驗(yàn)基于VP/VS波速比模型約束的大地電磁反演算法的可靠性,本文設(shè)計(jì)一套復(fù)雜電阻率模型(圖4a),其中包含有A、B、C1、C2、F電阻率塊體、近垂直構(gòu)造體D和上涌構(gòu)造體E.圖4a中黑色三角形為大地電磁測點(diǎn),基于有限單元法(Constable et al.,1987)正演計(jì)算10-3Hz~1103Hz 區(qū)間42個(gè)頻點(diǎn)的TE和TM模式的視電阻率和相位,分別加入5%的高斯隨機(jī)誤差,作為觀測數(shù)據(jù).

    圖4b為采用的VP/VS波速比固定模型,波速比值范圍為1.5~1.9.對比電阻率模型,波速比模型中設(shè)計(jì)了與電阻率模型結(jié)構(gòu)相似(A、B、C1、F、D、E)和結(jié)構(gòu)不相似(C2)的異常塊體或構(gòu)造體.將VP/VS波速比模型(圖4b)作為固定模型約束,分別采用無約束和有約束大地電磁反演算法獲得電阻率結(jié)構(gòu),對比分析這兩種方法的反演結(jié)果,檢驗(yàn)約束反演算法的可靠性和優(yōu)越性.

    圖4 理論電阻率模型(a)和VP/VS波速比模型(b)(a) 理論電阻率模型; (b) 固定約束的VP/VS波速比地震學(xué)模型. A、B、C1、C2、F表示不同的物性塊體;D為近垂直構(gòu)造體;E為上涌構(gòu)造體.Fig.4 The resistivity model (a) and the VP/VS wave velocity ratio model (b)(a) Synthetic resistivity model; (b) The fixed VP/VS wave ratio model. A, B, C1, C2 and F represent the different physical blocks; D represents the near vertical structure; E represents the upwelling structure.

    3.2 理論模型反演結(jié)果

    圖5a、b分別為無約束和有約束反演結(jié)果.從圖5a中可以發(fā)現(xiàn)大地電磁法在無約束反演情況下,淺部電阻率結(jié)構(gòu)具有較高的分辨率,例如塊體A的異常值和異常邊界都恢復(fù)得與真實(shí)模型相似、C1和C2塊體的上邊界接近于真實(shí)模型邊界、塊體B的電阻率值在淺部恢復(fù)得與真實(shí)值相近.但無約束反演結(jié)果在深部分辨率明顯不足,具體表現(xiàn)在:異常體C1和C2底邊界缺失,形成了與F貫通的兩個(gè)上涌構(gòu)造體;近垂直構(gòu)造體D在淺部有明顯的邊界信息,隨著深度增加,構(gòu)造體的邊界信息減弱;上涌構(gòu)造體E幾乎沒有體現(xiàn)出來.

    圖5 理論模型無約束和有約束的大地電磁二維反演結(jié)果(a) 無約束大地電磁反演結(jié)果; (b) 有約束大地電磁反演結(jié)果; (c) 無約束反演電阻率模型與VP/VS波速比模型的交叉梯度值空間分布; (d) 有約束反演電阻率模型與VP/VS波速比模型的交叉梯度值空間分布. RMS為數(shù)據(jù)擬合差.Fig.5 The results derived from the unconstrained and constrained inversion method(a) The resistivity model derived from the unconstrained inversion method; (b) The resistivity model derived from the constrained inversion method; (c) The spatial distribution of the cross-gradient value between the VP/VS ratio model and the unconstrained inverted resistivity model; (d) The spatial distribution of the cross-gradient value between the VP/VS ratio model and the constrained inverted resistivity model. RMS represents the data residual.

    圖5b為大地電磁法約束反演結(jié)果,從圖中可以看出電阻率反演結(jié)果在深度分辨率有明顯提高.在電阻率模型和VP/VS波速比模型結(jié)構(gòu)相似處電阻率邊界信息有明顯改善,具體表現(xiàn)在:C1和C2塊體由原來的連通狀態(tài)恢復(fù)成兩個(gè)獨(dú)立圈閉的塊體,其中C1的邊界信息與VP/VS波速比模型相似,在C1的異常值和邊界都恢復(fù)得與真實(shí)模型相似,而在C2位置處,因在VP/VS波速比模型中沒有存在C2異常,雖然C2異常值有改善,但邊界信息與真實(shí)情況仍有差距;上涌構(gòu)造E和近垂直構(gòu)造D的邊界信息恢復(fù)得與真實(shí)模型相近.圖5c和d分別為無約束反演結(jié)果、有約束反演結(jié)果與VP/VS波速比模型的交叉梯度值空間分布圖,從圖中可以發(fā)現(xiàn)在塊體C1周邊、近垂直構(gòu)造D、上涌構(gòu)造E的交叉梯度值有明顯降低,這說明經(jīng)過約束反演后的模型,在這幾個(gè)區(qū)域的結(jié)構(gòu)更相似.從反演結(jié)果中不難發(fā)現(xiàn)對于淺部電阻率值和形態(tài)(例如塊體A)在無約束和有約束反演都可以很好地恢復(fù),這說明基于地震學(xué)模型約束的大地電磁反演算法可以在滿足數(shù)據(jù)擬合差的前提下,獲取電阻率模型和地震學(xué)模型結(jié)構(gòu)上更為相似的模型,有效提高大地電磁法深部電阻率分辨率.

    4 張渤地震帶大地電磁測深剖面

    4.1 大地電磁數(shù)據(jù)采集和處理

    圖2為張渤帶大地電磁測點(diǎn)布設(shè)圖,測線跨越張家口、北京、唐山等城市,全長400 km,測線方向?yàn)楸逼?4°.測點(diǎn)布設(shè)采用“+”字形正南北東西向布極方式,采集加拿大鳳凰公司的MTU-5A寬頻帶大地電磁觀測設(shè)備.從該測線中挑選出數(shù)據(jù)質(zhì)量較高的31個(gè)MT測點(diǎn),分別對每個(gè)測點(diǎn)數(shù)據(jù)進(jìn)行濾波、重采樣、頻譜分析等預(yù)處理,基于WinGLink平臺獲取0.001 Hz到320.0 Hz的TE模式和TM模式的視電阻率和相位曲線.同時(shí)采用Berdichevskiy和Rhoplus方法檢驗(yàn)視電阻率和相位曲線的一致性,將處理后的數(shù)據(jù)作為反演數(shù)據(jù).

    圖6中黑色曲線為TE模式的視電阻率和相位,紅色曲線為TM模式視電阻率和相位.根據(jù)地質(zhì)構(gòu)造特征,將測線上MT測點(diǎn)數(shù)據(jù)分成四段,分別為唐山斷裂帶區(qū)域(圖6a,紅色測點(diǎn))、三河—平谷斷裂帶區(qū)域(圖6a,黃色測點(diǎn))、華北平原—京西丘陵過渡區(qū)(圖6a,藍(lán)色測點(diǎn))和京西丘陵區(qū)域(圖6a,黑色測點(diǎn)).從圖中可以看出除了華北平原—京西丘陵過渡帶的視電阻率和相位曲線形態(tài)差異比較大之外,其他三個(gè)區(qū)域的數(shù)據(jù)形態(tài)較為一致.其中唐山斷裂帶區(qū)域視電阻率曲線形態(tài)和京西丘陵區(qū)域曲線形態(tài)相似,但京西丘陵段的最高值出現(xiàn)在0.15 Hz左右,而唐山斷裂帶區(qū)域最高值出現(xiàn)在0.06 Hz左右,這兩段的相位曲線差異比較大,說明地下電性結(jié)構(gòu)山區(qū)和平原存在差別.三河—平谷斷裂帶視電阻率值比較低,最高值約300 Ωm,該區(qū)域的相位曲線形態(tài)與唐山斷裂帶相似.而華北平原—京西丘陵過渡帶的視電阻率和相位曲線形態(tài)復(fù)雜,沒有明顯的規(guī)律,這也間接說明該區(qū)域地下電性結(jié)構(gòu)復(fù)雜的特征.

    圖6 實(shí)測大地電磁資料視電阻率和相位曲線(a) 紅色三角形為Segment 1測點(diǎn)位置; 黃色三角形為Segment 2測點(diǎn)位置; 藍(lán)色三角形為Segment 3測點(diǎn)位置; 黑色三角形為Segment 4測點(diǎn)位置. 黑色曲線為TE模式的視電阻率和相位; 紅色曲線為TM模式的視電阻率和相位.Fig.6 Measured MT apparent resistivity and phase data(a) Red triangles represent the MT sites corresponding to the data of segment 1; Yellow triangles represent the MT sites corresponding to the data of segment 2; Blue triangles represent the MT sites corresponding to the data of segment 3; Black triangles represent the MT sites corresponding to the data of segment 4. Black curves represent the apparent resistivity and phase data of TE mode; red curves represent the apparent resistivity and phase data of TM mode.

    4.2 VP/VS波速比模型

    楊歧焱等(2018)基于華北地區(qū)數(shù)字地震臺網(wǎng)176個(gè)固定臺站,采用近震體波層析成像法獲取了張渤地震帶VP/VS波速比地震學(xué)模型.本文提取了張渤地震帶對應(yīng)大地電磁測線剖面的VP/VS波速比模型(圖7),對該剖面進(jìn)行平滑,去除模型的突變點(diǎn),然后將該模型作為固定約束代入公式(5—11)中的Wcg,s,實(shí)現(xiàn)基于VP/VS波速比模型約束的大地電磁法二維反演.

    圖7 張渤地震帶VP/VS波速比剖面圖TS 唐山; SH 三河; BJ 北京; ZJK 張家口.Fig.7 The VP/VS wave velocity ratio profile along Zhangbo seismic beltTS Tangshan; SH Shanhe; BJ Beijing; ZJK Zhangjiakou.

    4.3 反演網(wǎng)格

    二維大地電磁反演網(wǎng)格剖分中,用均勻和非均勻剖分法分別對目標(biāo)研究區(qū)域和邊界區(qū)域進(jìn)行網(wǎng)格化處理(圖8).在目標(biāo)研究區(qū)采用均勻網(wǎng)格剖分,X軸網(wǎng)格間距為3 km,Z軸網(wǎng)格間距為2 km.在邊界網(wǎng)格采用的是指數(shù)遞增網(wǎng)格間距.圖中8紅色框?yàn)閂P/VS波速比模型約束區(qū),在反演過程中只在紅色框區(qū)域內(nèi)進(jìn)行結(jié)構(gòu)耦合約束.

    圖8 大地電磁正、反演網(wǎng)格剖分黑色線為網(wǎng)格;黑色方框?yàn)榧s束反演區(qū)域.Fig.8 Mesh discretization for MT forward and inversionBlack lines represent the mesh; black rectangle shows the constrained inversion region.

    4.4 張渤地震帶反演結(jié)果

    圖9a為無約束和有約束的大地電磁法反演迭代RMS曲線圖,從圖中可知這兩種方法都具有較好的收斂性.從RMS迭代變化趨勢可以看出,當(dāng)加入VP/VS波速比模型約束項(xiàng)時(shí),RMS值迭代下降的速度變慢,但經(jīng)過多次迭代后,有約束和無約束反演結(jié)果的RMS都收斂到相同水平上.

    圖9 RMS迭代曲線和約束項(xiàng)權(quán)重因子L-Curve曲線(a) 無約束和有約束反演RMS迭代曲線; 空心圈曲線為無約束反演(SP)迭代RMS曲線; 實(shí)心點(diǎn)曲線為有約束反演(CG)迭代RMS曲線. (b) 約束項(xiàng)權(quán)重因子L-Curve曲線; 實(shí)心圈曲線為不同權(quán)重因子的RMS值; 黑色箭頭標(biāo)注為本文選取的權(quán)重因子.Fig.9 Plots of RMS iteration curves and L-curve analysis for selecting constrained weight factors(a) RMS iteration curves for no-constrained and constrained inversion. The hollow circle curve represents the RMS iteration curve for unconstrained inversion. The solid dot curve represents the RMS iteration curve for constrained inversion. (b) L-curve anasys for selcting constrained weight factors; The solid dot curve represents the RMS value for different weighting factors; the black arrow represnets the optimal weighting value.

    圖9b為公式中VP/VS波速比模型約束項(xiàng)權(quán)重因子(公式8中η)的選擇.本文采用“L-Curve”方法,選取了多個(gè)權(quán)重因子進(jìn)行二維大地電磁法約束反演(圖9b),圖中黑色箭頭為本文最終確定的權(quán)重因子.圖10展示了6個(gè)測點(diǎn)無約束和有約束反演模型響應(yīng)曲線,從圖中可以發(fā)現(xiàn)無約束和有約束反演模型響應(yīng)曲線基本上重合,與觀測數(shù)據(jù)有較好的一致性,說明這兩個(gè)模型都可以很好地?cái)M合實(shí)測數(shù)據(jù).

    圖10 無約束和有約束反演模型在測點(diǎn)01、09、16、19、26和30處的正演響應(yīng)紅色點(diǎn): 實(shí)測數(shù)據(jù)(real); 黑色圈: 無約束反演結(jié)果正演響應(yīng)(SP); 藍(lán)色叉: 有約束反演結(jié)果正演響應(yīng)(CG).Fig.10 The synthetic data sets of the unconstrained and constrained inverted model at the sites of 01, 09, 16, 26, and 30.Red dots represent the measured data (real); black circles represent the synthetic data sets of the unconstrained inverted model (SP); blue crosses represent the synthetic data sets of the constrained inverted model (CG).

    圖11展示了大地電磁法二維無約束和有約束反演結(jié)果,圖11a為測線對應(yīng)的地形曲線,可以發(fā)現(xiàn)測線的西部測點(diǎn)主要位于京西丘陵區(qū),測線的東部主要位于華北平原地區(qū),測線上由東向西穿越了唐山斷裂帶、三河—平谷斷裂帶、華北平原—京西丘陵過渡帶、延懷盆地和張家口斷裂帶.

    圖11 大地電磁法無約束和有約束反演結(jié)果(a) 地形剖面圖; (b) 無約束反演電阻率模型; (c) 有約束反演電阻率模型. H1—H5為高阻異常體編號; L1—L3為低阻異常體編號.Fig.11 The resistivity models derived from the unconstrained and constrained inversion method(a) Topographic profile; (b) The resistivity model derived from the unconstrained inversion method; (c) The resistivity model derived from the constrained inversion method. H1—H5 represent the numbered high resistivity anomalies; L1—L3 represent the numbered low resistivity anomalies.

    為了方便對比這兩種方法反演結(jié)果,將圖中5個(gè)高阻體標(biāo)注為H1、H2、H3、H4和H5,3個(gè)低阻區(qū)域標(biāo)注為L1、L2和L3.從圖中可以看出,基于VP/VS波速比模型約束的反演結(jié)果(圖11c)的橫向分辨率和縱向分辨率都有明顯的提高,主要表現(xiàn)有:在唐山強(qiáng)震區(qū)的高阻異常體H5沿伸到下地殼,經(jīng)過約束反演后,高阻異常的拖尾現(xiàn)象有明顯的改善,底部邊界約在20 km左右;從夏墊斷裂帶到京西山前斷裂帶的高阻體異常體H4近垂直地沿伸到地幔,很難識別出該異常體在地殼內(nèi)的異常形態(tài),經(jīng)約束反演后,深部的分辨率有明顯的提升,高阻異常體H4在地殼中表現(xiàn)為圈閉異常形態(tài);在京西丘陵區(qū)域的地殼高阻異常體H1、H2和H3在無約束反演中表現(xiàn)為整體高阻異常,深度延伸到下地殼,經(jīng)約束反演后,可以明顯看出H1、H2和H3這三個(gè)高阻異常體是獨(dú)立分開的三個(gè)高阻異常,同時(shí)這三個(gè)異常體的底界面大約在20 km;經(jīng)過約束反演后,發(fā)現(xiàn)深部高導(dǎo)區(qū)域L1和L2從30~40 km上涌,這個(gè)現(xiàn)象在無約束反演結(jié)果中沒有發(fā)現(xiàn);在夏墊斷裂帶區(qū)域無約束和有約束反演都可以發(fā)現(xiàn)明顯的低阻異常體L3,但低阻異常的形態(tài)、規(guī)模有很大的差異.對比VP/VS波速比模型(圖7)和大地電磁法約束反演結(jié)果(圖11c),在20 km范圍內(nèi),電阻率模型與VP/VS波速比模型的異常體表現(xiàn)出明顯的相似性.而在大于20 km深度,這兩個(gè)模型存在比較大的差異,表現(xiàn)在約30 km處波速比模型出現(xiàn)低波速比異常層,而電阻率模型上未發(fā)現(xiàn)該特征.

    圖12為大地電磁法無約束和有約束反演的電阻率模型與VP/VS波速比模型的交叉梯度值空間分布圖,mean|ty|為約束反演區(qū)域(圖8黑色框)中每個(gè)網(wǎng)格的交叉梯度值的絕對值的平均值.從圖中可以發(fā)現(xiàn)有約束反演結(jié)果與VP/VS波速比模型的交叉梯度值比無約束反演結(jié)果的交叉梯度值有明顯下降(mean|ty|由2.94降到1.5),這說明基于VP/VS波速比模型約束的大地電磁反演可以在滿足數(shù)據(jù)擬合差的前提下,可以獲取與地震學(xué)模型結(jié)構(gòu)上更為耦合的電阻率模型.

    圖12 無約束、有約束反演電阻率模型與VP/VS波速比模型的交叉梯度值分布圖(a) 無約束反演電阻率模型與VP/VS波速比模型的交叉梯度值空間分布; (b) 有約束反演電阻率模型與VP/VS波速比模型的交叉梯度值空間分布. mean|ty|為約束反演區(qū)域(圖8紅色框)中每個(gè)網(wǎng)格的交叉梯度值的絕對值的平均值.Fig.12 The spatial distribution of the cross-gradient value between the VP/VS ratio model and the non-constrained or the constrained inverted resistivity model.(a) The spatial distribution of the cross-gradient value between the VP/VS ratio model and the unconstrained inverted resistivity model; (b) The spatial distribution of the cross-gradient value between the VP/VS ratio model and the constrained inverted resistivity model. mean |ty| represents the average of the absolute value of the cross-gradient values in the common region (red box in Fig.8).

    根據(jù)上述現(xiàn)象,不難看出要使得基于地震學(xué)模型約束反演算法能有效改善大地電磁反演結(jié)果的前提條件是地震學(xué)模型和電阻率模型存在結(jié)構(gòu)相似性.當(dāng)兩種物性結(jié)構(gòu)相似時(shí),基于交叉梯度的約束反演可以有效地加強(qiáng)兩種物性結(jié)構(gòu)耦合關(guān)系,提高結(jié)構(gòu)相似處的分辨率.在結(jié)構(gòu)不相似的區(qū)域,交叉梯度項(xiàng)幾乎不起作用,保留單方法反演的優(yōu)勢.Wu 等(2020)、Peng 等(2019)通過大量的理論模型算例充分地論證這一特征.

    4.5 張渤地震帶約束反演結(jié)果可靠性和有效性檢驗(yàn)

    通過對比無約束、有約束的大地電磁反演結(jié)果,可知在地震學(xué)模型的約束下,張渤地震帶深部電性結(jié)構(gòu)在橫向和縱向的分辨率都有明顯的提升.為了檢驗(yàn)基于VP/VS波速比模型約束的大地電磁反演算法在張渤地震帶應(yīng)用的可靠性和有效性,本文考慮設(shè)計(jì)形態(tài)與反演結(jié)果相似的5個(gè)獨(dú)立的高阻體(圖11中H1—H5)和3個(gè)上涌的低阻構(gòu)造體(圖11中L1—L3),通過理論模型合成數(shù)據(jù)試算,采用與實(shí)測數(shù)據(jù)反演一樣的VP/VS波速比模型(圖7)進(jìn)行約束反演,對比無約束和有約束反演結(jié)果,從理論模型上檢驗(yàn)該算法是否能提高張渤地震帶深部電阻率橫向和縱向分辨率.

    為了能滿足上述5個(gè)獨(dú)立的高阻體和3個(gè)上涌的高導(dǎo)構(gòu)造體的理論模型需求,本文將約束反演結(jié)果(圖11c)中電阻率值大于103.5Ωm的電阻率值設(shè)置為恒定值103Ωm,將小于101.8Ωm的電阻率值設(shè)置為恒值10 Ωm,其他區(qū)域的值設(shè)置為102Ωm.圖13a為本文設(shè)計(jì)的理論電阻率模型,從圖中可以看出理論模型在20km以內(nèi)設(shè)計(jì)有5個(gè)圈閉的高阻異常體(圖13a 中H1—H5),中下地殼3個(gè)低阻隆起(圖13a 中L1—L3).

    測點(diǎn)位置(圖13a黑色三角形)和每個(gè)測點(diǎn)的頻段選擇與實(shí)際觀測一樣,正反演網(wǎng)格也跟上述一樣.正演每個(gè)測點(diǎn)的理論視電阻率和相位曲線,同時(shí)加入與實(shí)測數(shù)據(jù)相同的誤差,形成觀測數(shù)據(jù).

    圖13 理論電阻率模型(a)、無約束(b)和有約束(c)反演結(jié)果(a) 理論電阻率模型; (b) 無約束反演電阻率模型; (c) 有約束反演電阻率模型. H1—H5為高阻異常體編號; L1—L3為低阻異常體編號.Fig.13 Synthetic model (a), unconstrained (b) and constrained (c) inverted models(a) Synthetic model; (b) Unconstrained inverted model; (c) Constrained inverted model. H1—H5 represent the numbered high resistivity anomalies; L1—L3 represent the numbered low resistivity anomalies.

    圖13b和c分別為大地電磁法二維無約束和有約束的反演結(jié)果,對比這兩種方法的反演結(jié)果,可以發(fā)現(xiàn)無約束反演的結(jié)果與上述實(shí)際數(shù)據(jù)無約束反演結(jié)果出現(xiàn)的問題相似,主要表現(xiàn)在:獨(dú)立高阻異常體H1—H4在無約束反演(圖13b)中表現(xiàn)為連接在一起、規(guī)模大的高阻異常;低阻體L1和L2幾乎沒有體現(xiàn)出來,L3的低阻異常體形態(tài)與真實(shí)異常體形態(tài)差別很大;同時(shí)異常體的底界面相比真實(shí)模型有明顯的往下拖尾的現(xiàn)象.從無約束反演結(jié)果可以看出因大地電磁法受限于方法本身體積效應(yīng)、深度分辨率隨深度增加而降低等局限,無約束反演獲得的電阻率模型在橫向和垂向分辨率有限.

    圖13c為基于VP/VS波速比模型約束反演的結(jié)果,從圖中可以明顯的看出淺部高阻異常體(H1—H5)的橫向分辨率和縱向分辨率都有明顯的提高,5個(gè)高阻異常體可以清晰的展示出獨(dú)立異常體的特征,異常體的拖尾現(xiàn)象得到明顯的改善;對于低阻體L1、L2和L3異常體的反演精度有明顯的提高.因此,通過這個(gè)理論模型合成數(shù)據(jù)試算,驗(yàn)證了基于VP/VS波速比模型約束的二維大地電磁反演法在張渤地震帶應(yīng)用可以提高深部電阻率的橫向和縱向分辨率,從而驗(yàn)證了該算法在張渤地震帶應(yīng)用的可靠性和有效性.

    5 討論

    根據(jù)反演得到的張渤地震帶地殼電性結(jié)構(gòu),結(jié)合地質(zhì)資料繪制電性構(gòu)造解釋圖(斷裂信息參考徐錫偉等(2002)).圖14中斷裂F1—F9的深部形態(tài)分別用黑色虛線展示,斷裂地表位置與圖2中斷裂(F1—F9)位置對應(yīng).從整體上來看,電阻率結(jié)構(gòu)表現(xiàn)出明顯的橫向不均勻性,與地表的各個(gè)構(gòu)造相呼應(yīng).對比圖11a地形曲線和構(gòu)造分區(qū),可以發(fā)現(xiàn)測線上斷裂帶在電性結(jié)構(gòu)上表現(xiàn)為明顯的電性分界面或電性梯度帶,在張家口斷裂帶、延慶—懷柔盆地、夏墊斷裂帶和唐山斷裂帶在中上地殼普遍存在高阻異常體,而在中下地殼表現(xiàn)為低阻異常特征.從M≥3.0地震分布特征來看,發(fā)現(xiàn)該區(qū)域的地震活動與斷裂帶的分布息息相關(guān),具有明顯的分段性,以下將分段討論剖面電性結(jié)構(gòu)特征.

    唐山震群主要有三條斷裂帶控制(寧河—昌黎斷裂F1、唐山斷裂F2和豐臺—野雞坨斷裂F3).唐山7.8 級地震發(fā)生在主斷裂帶(F2)上,該區(qū)域?yàn)楦咦鑹K體與低阻塊體的相交地帶、偏向高阻體一側(cè).根據(jù)人工地震反射資料可知唐山斷裂為主斷裂,延伸至Moho面,地殼伴生的斷裂構(gòu)造表現(xiàn)為走滑斷裂(花狀構(gòu)造)特征.圖中可以發(fā)現(xiàn)高阻體的底界面大約在20 km左右,地震的發(fā)生主要集中在20 km以內(nèi),位于高阻區(qū)域或高低阻交界區(qū).人工深反射剖面(劉保金等, 2011a)和天然地震觀測(劉啟元等, 2007; 楊歧焱等, 2018)發(fā)現(xiàn)唐山震區(qū)在深度<20 km內(nèi)呈現(xiàn)高速異常,在約深20 km處存在明顯的低速層(體),該低速層(體)一直延伸到40 km以下的深度.劉國棟等(1983)發(fā)現(xiàn)唐山地震源區(qū)下約20 km深處存在高導(dǎo)異常.嘉世旭和張先康(2005)通過深地震測深資料也發(fā)現(xiàn)唐山強(qiáng)震區(qū)中下地殼存在低速異常.巖漿底侵和置換作用是下地殼和地幔物質(zhì)交換、相互作用的重要形式,推測唐山震區(qū)中下地殼的低速低阻異常體與流體和熱物質(zhì)作用有關(guān),大規(guī)模的幔源物質(zhì)侵入,形成了中上地殼的孕震環(huán)境.

    在三河—平谷斷裂帶區(qū)域,存在明顯的殼淺部到深部低阻特征,中下地殼存在上涌低阻體(L3),淺部和深部的高導(dǎo)體不連續(xù),地震分布主要位于該不連續(xù)層中,說明該區(qū)域深淺高導(dǎo)異常的成因可能不同.鄧前輝等(2001)通過大地電磁法也揭示夏墊斷裂帶中下地殼高導(dǎo)體和地幔上涌現(xiàn)象,淺部和深部的構(gòu)造存在差異,推測震源體的形成和發(fā)震與深部高角度隱伏斷裂有關(guān).對應(yīng)地質(zhì)資料可知,三河—平谷斷裂帶為第四紀(jì)活動隱伏高傾角斷裂帶,地表覆蓋有上新統(tǒng)-第四系沉積層,淺部低阻異常與沉積層對應(yīng).從人工地震深反射剖面(劉保金等, 2011b; 趙金仁等, 2004)可知該斷裂帶是一條深淺共存的斷裂構(gòu)造帶,在Moho存在復(fù)雜的楔形反射帶,說明該區(qū)域地殼發(fā)生過強(qiáng)烈的擠壓、變形,同時(shí)也反映出巖漿活動對地下地殼結(jié)構(gòu)經(jīng)過強(qiáng)烈改造,形成了復(fù)雜的地殼深部結(jié)構(gòu).徐錫偉等(2002)也從人工地震反射剖面發(fā)現(xiàn)Moho面存在反射過渡帶,推測是早期地殼伸展過程中形成斷裂帶,填充了來自上地幔的基性巖墻.從圖14中深部上涌低阻異常(L3)可能反映的是上地幔熱物質(zhì)侵入,Zhu等(2012, 2015)指出西太平洋板塊的俯沖時(shí),俯沖板塊的脫水交代上覆地幔楔使其發(fā)生熔融,俯沖帶后撤后導(dǎo)致巖石圈強(qiáng)烈伸展,導(dǎo)致地幔熱物質(zhì)上涌,形成了中上地殼的孕震環(huán)境.

    黃莊—高麗營斷裂帶(F5)和南口斷裂(F6)為京西丘陵的山前斷裂帶,是大型拆離斷裂帶,在電阻率結(jié)構(gòu)上表現(xiàn)為明顯的電阻率變化梯度帶.在山前斷裂帶到夏墊斷裂帶之間淺部有幾公里的反映沉積層特征的低阻層,在低阻層下方、延伸到下地殼為規(guī)模較大的高阻異常,未見低速體上涌現(xiàn)象.詹艷等(2011)從石家莊地區(qū)大地電磁剖面也發(fā)現(xiàn)太行山山前斷裂從淺地表到下地殼表現(xiàn)為高電阻體,無殼內(nèi)高導(dǎo)層發(fā)育.京西丘陵的山前斷裂帶是華北克拉通東部與中部的邊界,在重力梯度帶、地形、巖石圈厚度等資料顯示在該區(qū)域附近變化顯著(Zhu et al., 2012).規(guī)模大、延伸深的高阻異常(H4)可能是華北克拉通東部和中部地區(qū)構(gòu)造運(yùn)動形成的高密度的輝長巖、鐵鎂質(zhì)巖石或者下地殼麻粒巖(翟明國等, 2005; 翟明國, 2019).華北平原地區(qū)強(qiáng)震主要發(fā)生在NE和NW兩組構(gòu)造帶,黃莊—高麗營斷裂帶(F5)和南口斷裂(F6)山前斷裂帶自新生代以來對北京凹陷起著重要的控制作用(焦青等, 2005;焦青和邱澤華, 2006),是北京地區(qū)值得關(guān)注的強(qiáng)震發(fā)生重點(diǎn)監(jiān)視防御區(qū).

    在南口斷裂(F6)以西區(qū)域?yàn)榫┪髑鹆甑貐^(qū),主要有延慶—懷來盆地和張家口斷裂帶.延懷盆地是由山西裂谷走滑運(yùn)動形成的NE向次級張性構(gòu)造區(qū),地表延慶—礬山次級盆地北緣斷裂(F7)和涿鹿—懷來盆地北緣斷裂(F8)控制盆地構(gòu)造特征.在電阻率結(jié)構(gòu)成像中(圖14),在淺部表現(xiàn)明顯的低阻異常,在中下地殼有明顯的低阻體沿著盆地下方通道侵入的現(xiàn)象.從地震反射剖面(張先康等, 1996)可知在延慶—礬山盆地北緣斷裂正下方存在低速體,其成因可能與該區(qū)域地殼伸展運(yùn)動、幔源物質(zhì)上涌有關(guān),在地表形成盆嶺構(gòu)造.

    圖14 張渤地震帶電性結(jié)構(gòu)解釋圖白色圈為2007年至今M≥3.0地震(由中國地震臺網(wǎng)中心提供). F1 寧河—昌黎斷裂; F2 唐山斷裂; F3 豐臺—野雞坨斷裂; F4 夏墊斷裂; F5 黃莊—高麗營斷裂; F6 南口斷裂; F7 延慶—礬山盆地北緣斷裂; F8 懷(來)—涿(鹿)盆地北緣斷裂; F9 張家口斷裂. F1—F8走向?yàn)镹E或NEE;F9走向?yàn)镹W.Fig.14 The geology interpretation maps of the electrical structure beneath Zhangbo seismic beltWhite circles: M≥3.0 earthquakes since 2007 (provided by China earthquake Networks Center). F1 Ninghe-Changli fault; F2 Tangshan fault; F3 Fengtai-Yejituo fault; F4 Xiadian fault; F5 Huangzhuang-Gaoliying fault; F6 Nankou fault; F7 North margin fault of Yanqing-Fanshan basin; F8 North margin fault of Huailai-Zhoulu basin; F9 Zhangjiakou fault. The strike of F1—F8 is NE or NEE; the strike of F9 is NW.

    張家口斷裂帶位于張渤地震帶和山西地震帶的交匯區(qū)(趙博等, 2011),區(qū)域內(nèi)斷裂相交和切割,形成復(fù)雜的孕震構(gòu)造格局,主要由北西西向和近東西向兩組、 多段斷裂組成, 北西西向斷裂為斷裂帶主體(圖14中F9),與研究區(qū)數(shù)條北北東至北東向的斷裂顯著不同(高戰(zhàn)武等, 2001),是首都圈地震重點(diǎn)監(jiān)視防御區(qū).張學(xué)民等(2012)通過地震體波層析成像發(fā)現(xiàn)太行山北段多構(gòu)造交界的張家口地區(qū)存在局部低速異常,推斷該區(qū)域受到西太平洋板塊俯沖的影響,表現(xiàn)為拉張應(yīng)力場,深部地幔熱物質(zhì)活動及淺部應(yīng)力累積引發(fā)活躍的構(gòu)造活動.地質(zhì)資料顯示張家口地區(qū)出露中生代形成的源于古老下地殼或地幔的堿性花崗巖(楊進(jìn)輝等, 2006; 韋忠良等, 2008; 汪洋和程素華, 2010),從巖石學(xué)、地球化學(xué)角度論證了該區(qū)域經(jīng)歷了大規(guī)模的殼幔相互作用.從本文獲取的電性結(jié)構(gòu)圖上可以發(fā)現(xiàn)在該斷裂帶上深部有低阻體上涌,在張家口中下地殼高導(dǎo)層的規(guī)模比懷來盆地大,這可能暗示張家口斷裂帶的深部物質(zhì)作用更為強(qiáng)烈.

    6 結(jié)論

    本文實(shí)現(xiàn)了基于VP/VS波速比模型約束的大地電磁二維反演算法,通過理論模型合成數(shù)據(jù)試算,檢驗(yàn)了算法的可靠性.將約束反演算法應(yīng)用于張渤地震帶31個(gè)大地電磁實(shí)測資料,研究發(fā)現(xiàn)唐山斷裂帶中上地殼表現(xiàn)為高阻異常,地震主要分布于高阻異常體周邊,在下地殼底部有上涌的低阻異常體,推測該區(qū)域深部低阻異常體與幔源物質(zhì)侵入有關(guān).三河—平谷斷裂帶淺部表現(xiàn)為沉積層的低阻異常特征,深部高導(dǎo)異常與夏墊斷裂帶展布、殼幔物質(zhì)相互作用有關(guān),地震事件主要發(fā)生在淺部和深部高導(dǎo)異常不連續(xù)區(qū)域.太行山山前斷裂電阻率結(jié)構(gòu)上表現(xiàn)為明顯的電阻率變化梯度帶,在與華北平原過渡帶中下地殼存在規(guī)模大、延伸深的高阻異常,推測與該區(qū)域因構(gòu)造運(yùn)動形成高密度的輝長巖、鐵鎂質(zhì)巖石或者下地殼麻粒巖有關(guān).延懷盆地在淺部表現(xiàn)明顯的低阻異常,在中下地殼有明顯的低阻體沿著盆地下方通道侵入的現(xiàn)象,其成因可能是地殼伸展運(yùn)動,幔源物質(zhì)上涌,形成地表的盆嶺構(gòu)造.張家口斷裂帶中下地殼高導(dǎo)異常層比懷來盆地深部高導(dǎo)異常規(guī)模大,這可能暗示張家口斷裂帶的深部物質(zhì)作用更為強(qiáng)烈.

    致謝感謝中國地質(zhì)大學(xué)(北京)王瑜教授和中國地震局地球物理研究所尤惠川研究員對本文的地質(zhì)構(gòu)造解釋提供的幫助.感謝審稿專家和編輯部對本文提出的重要修改意見.

    猜你喜歡
    波速斷裂帶電阻率
    冷凍斷裂帶儲層預(yù)測研究
    基于實(shí)測波速探討地震反射波法超前預(yù)報(bào)解譯標(biāo)志
    依蘭—伊通斷裂帶黑龍江段構(gòu)造運(yùn)動特征
    地震研究(2021年1期)2021-04-13 01:05:12
    三維電阻率成像與高聚物注漿在水閘加固中的應(yīng)用
    吉林地區(qū)波速比分布特征及構(gòu)造意義
    準(zhǔn)噶爾盆地西北緣克-夏斷裂帶構(gòu)造特征新認(rèn)識
    郯廬斷裂帶及兩側(cè)地區(qū)強(qiáng)震異常特征分析
    中國地震(2015年1期)2015-11-08 11:11:26
    隨鉆電阻率測井的固定探測深度合成方法
    海洋可控源電磁場視電阻率計(jì)算方法
    粉煤灰摻量對水泥漿體電阻率與自收縮的影響
    天天躁夜夜躁狠狠躁躁| 久久99一区二区三区| 水蜜桃什么品种好| 在线观看国产h片| 国产精品久久久久久久电影| 十八禁网站网址无遮挡| 久久狼人影院| 久久精品国产鲁丝片午夜精品| 黄色视频在线播放观看不卡| 中文字幕av电影在线播放| 免费av中文字幕在线| 亚洲熟女精品中文字幕| 美女脱内裤让男人舔精品视频| 国产免费又黄又爽又色| 免费高清在线观看日韩| 久久久久人妻精品一区果冻| 99视频精品全部免费 在线| 精品久久蜜臀av无| 90打野战视频偷拍视频| 最近最新中文字幕大全免费视频 | 欧美精品人与动牲交sv欧美| 香蕉国产在线看| 免费不卡的大黄色大毛片视频在线观看| 啦啦啦在线观看免费高清www| 欧美少妇被猛烈插入视频| 欧美精品人与动牲交sv欧美| 久久精品人人爽人人爽视色| 久久人人爽av亚洲精品天堂| 日韩电影二区| 少妇的逼水好多| 亚洲经典国产精华液单| 最近中文字幕2019免费版| 国产成人精品久久久久久| 日本午夜av视频| 精品人妻偷拍中文字幕| 高清黄色对白视频在线免费看| 亚洲精品一二三| 日韩一区二区三区影片| 日本爱情动作片www.在线观看| 男女下面插进去视频免费观看 | www日本在线高清视频| 男女边吃奶边做爰视频| 国产男女内射视频| 中文精品一卡2卡3卡4更新| 美女大奶头黄色视频| 亚洲欧美精品自产自拍| 欧美亚洲 丝袜 人妻 在线| 一区二区av电影网| 哪个播放器可以免费观看大片| 亚洲av男天堂| 少妇的逼好多水| 卡戴珊不雅视频在线播放| 亚洲精品第二区| 91成人精品电影| 久久婷婷青草| 久久久久久久国产电影| 亚洲精品成人av观看孕妇| 青青草视频在线视频观看| 亚洲欧美成人综合另类久久久| 亚洲av国产av综合av卡| 日本免费在线观看一区| 国产永久视频网站| 伊人亚洲综合成人网| 精品久久久精品久久久| 夫妻性生交免费视频一级片| 亚洲av在线观看美女高潮| 国产精品蜜桃在线观看| 99久国产av精品国产电影| 两个人看的免费小视频| 日本猛色少妇xxxxx猛交久久| 国产一区二区在线观看av| 成人国产麻豆网| 日韩在线高清观看一区二区三区| 男女边摸边吃奶| 母亲3免费完整高清在线观看 | av女优亚洲男人天堂| 夜夜骑夜夜射夜夜干| 久久精品久久精品一区二区三区| 人妻人人澡人人爽人人| 最近2019中文字幕mv第一页| 精品午夜福利在线看| 免费日韩欧美在线观看| 90打野战视频偷拍视频| 不卡视频在线观看欧美| 国产淫语在线视频| 欧美3d第一页| 一二三四在线观看免费中文在 | 久久精品久久精品一区二区三区| 久久ye,这里只有精品| 菩萨蛮人人尽说江南好唐韦庄| av.在线天堂| 九九爱精品视频在线观看| 中文欧美无线码| 全区人妻精品视频| 欧美人与善性xxx| 国产精品不卡视频一区二区| 大香蕉久久网| 一区二区三区精品91| 国产亚洲欧美精品永久| 亚洲婷婷狠狠爱综合网| 亚洲成人手机| 欧美日韩亚洲高清精品| 天天影视国产精品| 亚洲国产精品999| 国产午夜精品一二区理论片| 亚洲国产看品久久| 成人国语在线视频| 成年女人在线观看亚洲视频| 91精品三级在线观看| 久久精品熟女亚洲av麻豆精品| 欧美老熟妇乱子伦牲交| 亚洲成av片中文字幕在线观看 | av卡一久久| 最近最新中文字幕免费大全7| 一边摸一边做爽爽视频免费| 亚洲av成人精品一二三区| 一级片免费观看大全| 国产欧美另类精品又又久久亚洲欧美| 狂野欧美激情性xxxx在线观看| 人妻少妇偷人精品九色| 国产淫语在线视频| 国产福利在线免费观看视频| 亚洲国产精品专区欧美| 欧美日韩视频高清一区二区三区二| 欧美 亚洲 国产 日韩一| 男女下面插进去视频免费观看 | 又大又黄又爽视频免费| 精品国产乱码久久久久久小说| 欧美精品人与动牲交sv欧美| 久久久久久久精品精品| 欧美国产精品一级二级三级| 久久 成人 亚洲| 91精品国产国语对白视频| 婷婷色综合www| 国产免费视频播放在线视频| 狂野欧美激情性xxxx在线观看| 多毛熟女@视频| 欧美性感艳星| 最近中文字幕高清免费大全6| 国产精品99久久99久久久不卡 | 97在线人人人人妻| 桃花免费在线播放| 久久久久久久久久久久大奶| 街头女战士在线观看网站| 亚洲精品久久久久久婷婷小说| 久久99精品国语久久久| av在线观看视频网站免费| 亚洲美女黄色视频免费看| 高清av免费在线| 青春草亚洲视频在线观看| 国产成人免费无遮挡视频| av女优亚洲男人天堂| 日本av手机在线免费观看| 毛片一级片免费看久久久久| 日本av手机在线免费观看| 亚洲欧美清纯卡通| 欧美国产精品va在线观看不卡| 久久人人爽人人片av| 亚洲精品中文字幕在线视频| h视频一区二区三区| 插逼视频在线观看| 亚洲av欧美aⅴ国产| 国产免费现黄频在线看| av天堂久久9| 久久人人爽av亚洲精品天堂| 久久97久久精品| 国产1区2区3区精品| 国产精品久久久久成人av| 国产无遮挡羞羞视频在线观看| 一个人免费看片子| videosex国产| 亚洲av日韩在线播放| 亚洲欧洲国产日韩| 天堂8中文在线网| 嫩草影院入口| 亚洲,一卡二卡三卡| 一二三四在线观看免费中文在 | 亚洲av免费高清在线观看| 嫩草影院入口| 欧美丝袜亚洲另类| 人妻 亚洲 视频| 欧美xxⅹ黑人| 2018国产大陆天天弄谢| 国产 精品1| 亚洲国产看品久久| 日本免费在线观看一区| 国产精品无大码| 国产片内射在线| 国产精品久久久久成人av| 高清在线视频一区二区三区| 国产免费现黄频在线看| 视频区图区小说| 激情五月婷婷亚洲| 久久久久网色| 自线自在国产av| 肉色欧美久久久久久久蜜桃| 一本—道久久a久久精品蜜桃钙片| 国语对白做爰xxxⅹ性视频网站| 咕卡用的链子| 一级毛片我不卡| 91精品伊人久久大香线蕉| 如何舔出高潮| 午夜免费观看性视频| 欧美变态另类bdsm刘玥| av一本久久久久| 男女边吃奶边做爰视频| 国产成人a∨麻豆精品| 欧美精品亚洲一区二区| 嫩草影院入口| 女人被躁到高潮嗷嗷叫费观| h视频一区二区三区| 国产成人精品福利久久| 五月开心婷婷网| 九九爱精品视频在线观看| 九色成人免费人妻av| 啦啦啦在线观看免费高清www| 亚洲国产欧美在线一区| 日本av手机在线免费观看| 国精品久久久久久国模美| 51国产日韩欧美| videos熟女内射| 国产色婷婷99| 伦理电影大哥的女人| 日本-黄色视频高清免费观看| 最新的欧美精品一区二区| 美女主播在线视频| 爱豆传媒免费全集在线观看| 欧美亚洲日本最大视频资源| 2021少妇久久久久久久久久久| 女人被躁到高潮嗷嗷叫费观| 国产一区二区激情短视频 | www日本在线高清视频| 少妇猛男粗大的猛烈进出视频| 久久99热这里只频精品6学生| 国产免费一区二区三区四区乱码| 最近的中文字幕免费完整| 内地一区二区视频在线| 一本一本久久a久久精品综合妖精 国产伦在线观看视频一区 | 久久青草综合色| 亚洲成人一二三区av| 亚洲美女搞黄在线观看| 少妇人妻精品综合一区二区| av又黄又爽大尺度在线免费看| 黄片播放在线免费| 亚洲欧美日韩另类电影网站| 97在线人人人人妻| 熟女av电影| 免费高清在线观看视频在线观看| 国产高清国产精品国产三级| 欧美国产精品va在线观看不卡| 国产精品久久久久久精品古装| 青春草视频在线免费观看| 色94色欧美一区二区| 街头女战士在线观看网站| 有码 亚洲区| 中文乱码字字幕精品一区二区三区| 波多野结衣一区麻豆| 午夜av观看不卡| a级片在线免费高清观看视频| 人体艺术视频欧美日本| 超色免费av| 亚洲欧美日韩卡通动漫| 成年av动漫网址| 人体艺术视频欧美日本| 黑人猛操日本美女一级片| 亚洲国产色片| 乱码一卡2卡4卡精品| 久久久久久久久久久免费av| 亚洲伊人久久精品综合| 最近的中文字幕免费完整| 日韩,欧美,国产一区二区三区| 免费在线观看黄色视频的| 亚洲成人手机| 成年av动漫网址| 大香蕉久久网| 日韩精品有码人妻一区| 日本爱情动作片www.在线观看| av女优亚洲男人天堂| 美女中出高潮动态图| 纵有疾风起免费观看全集完整版| 国国产精品蜜臀av免费| 51国产日韩欧美| 97超碰精品成人国产| 欧美 亚洲 国产 日韩一| av在线播放精品| 韩国高清视频一区二区三区| 黄色一级大片看看| 日韩人妻精品一区2区三区| 免费不卡的大黄色大毛片视频在线观看| 久久久久精品性色| 国产在线免费精品| 国产精品免费大片| 黄色 视频免费看| 精品酒店卫生间| 亚洲av成人精品一二三区| 精品国产国语对白av| 国产高清三级在线| 久久免费观看电影| 大片免费播放器 马上看| 80岁老熟妇乱子伦牲交| 国产欧美日韩综合在线一区二区| 国内精品宾馆在线| h视频一区二区三区| 九九在线视频观看精品| 久久久欧美国产精品| 国产欧美日韩一区二区三区在线| 另类亚洲欧美激情| 日本黄色日本黄色录像| 一本—道久久a久久精品蜜桃钙片| 久久99热6这里只有精品| 大陆偷拍与自拍| 亚洲国产精品一区三区| 18禁观看日本| 午夜福利网站1000一区二区三区| 久久久精品94久久精品| 国产精品嫩草影院av在线观看| 国产伦理片在线播放av一区| 又大又黄又爽视频免费| 亚洲国产精品成人久久小说| 国产精品不卡视频一区二区| 少妇人妻 视频| 水蜜桃什么品种好| 国产探花极品一区二区| 国产黄色免费在线视频| 亚洲国产成人一精品久久久| 99热这里只有是精品在线观看| 亚洲欧美成人精品一区二区| 男的添女的下面高潮视频| 国产男人的电影天堂91| 在线免费观看不下载黄p国产| 熟女av电影| 日本wwww免费看| 亚洲av免费高清在线观看| 午夜免费男女啪啪视频观看| 亚洲欧美精品自产自拍| 一本色道久久久久久精品综合| 国产老妇伦熟女老妇高清| 精品一品国产午夜福利视频| 久久国产精品大桥未久av| 十分钟在线观看高清视频www| 亚洲国产最新在线播放| 久久精品国产亚洲av天美| 国产在线免费精品| 99热网站在线观看| 亚洲色图综合在线观看| 99热全是精品| 女人被躁到高潮嗷嗷叫费观| 久久国产亚洲av麻豆专区| 少妇人妻精品综合一区二区| 精品少妇内射三级| 内地一区二区视频在线| 欧美日韩综合久久久久久| 多毛熟女@视频| 九九爱精品视频在线观看| 欧美精品一区二区大全| 在线观看国产h片| 18禁动态无遮挡网站| 精品国产一区二区三区久久久樱花| 免费观看a级毛片全部| 亚洲精品国产av成人精品| 秋霞伦理黄片| 成人手机av| 在线观看免费高清a一片| 精品国产一区二区三区四区第35| 一级片免费观看大全| 一级a做视频免费观看| 免费看av在线观看网站| 大香蕉久久网| 国产成人欧美| 在线观看人妻少妇| 超色免费av| 精品久久久精品久久久| 夜夜爽夜夜爽视频| 国产成人精品久久久久久| 熟女av电影| 国产成人av激情在线播放| 丁香六月天网| 久久99精品国语久久久| 国产高清国产精品国产三级| 超色免费av| 成人国产av品久久久| 久久97久久精品| 亚洲欧美色中文字幕在线| 成人黄色视频免费在线看| 国产精品久久久久久av不卡| 考比视频在线观看| 日韩精品有码人妻一区| 黄色视频在线播放观看不卡| 97超碰精品成人国产| av又黄又爽大尺度在线免费看| 欧美成人午夜免费资源| 91午夜精品亚洲一区二区三区| 日本爱情动作片www.在线观看| 国产深夜福利视频在线观看| 国产av码专区亚洲av| 有码 亚洲区| 久久久久久久大尺度免费视频| videosex国产| 免费不卡的大黄色大毛片视频在线观看| 韩国高清视频一区二区三区| 国产av码专区亚洲av| 亚洲综合色惰| 人人妻人人爽人人添夜夜欢视频| 大香蕉久久成人网| 婷婷色综合www| 午夜免费观看性视频| 精品国产一区二区久久| 丰满饥渴人妻一区二区三| 在线看a的网站| 王馨瑶露胸无遮挡在线观看| 黄色配什么色好看| 伦理电影免费视频| 中文字幕亚洲精品专区| 久久99一区二区三区| 国产综合精华液| 另类亚洲欧美激情| 国产国拍精品亚洲av在线观看| 久久精品久久久久久久性| 国产精品 国内视频| 久久久久网色| 黄色一级大片看看| 国产男女内射视频| 日韩视频在线欧美| freevideosex欧美| 在线观看三级黄色| 大香蕉97超碰在线| 麻豆精品久久久久久蜜桃| 欧美日本中文国产一区发布| 少妇的逼好多水| av女优亚洲男人天堂| 国产精品国产av在线观看| 亚洲精品日韩在线中文字幕| 99久久精品国产国产毛片| 日韩成人伦理影院| 亚洲国产av新网站| 狂野欧美激情性bbbbbb| 中文字幕人妻熟女乱码| 美女xxoo啪啪120秒动态图| 精品一区二区三卡| 亚洲色图 男人天堂 中文字幕 | 日韩一区二区三区影片| 交换朋友夫妻互换小说| 亚洲四区av| 精品酒店卫生间| 夫妻午夜视频| 日本91视频免费播放| 日本免费在线观看一区| 综合色丁香网| 青春草视频在线免费观看| 看非洲黑人一级黄片| 国产精品久久久久久精品电影小说| 国产成人精品福利久久| 久久精品国产亚洲av涩爱| 日韩不卡一区二区三区视频在线| 亚洲五月色婷婷综合| 免费女性裸体啪啪无遮挡网站| 午夜福利视频精品| 国产淫语在线视频| 七月丁香在线播放| 在线亚洲精品国产二区图片欧美| 男女边摸边吃奶| 国产精品嫩草影院av在线观看| 2018国产大陆天天弄谢| av在线老鸭窝| 久久这里只有精品19| 欧美成人精品欧美一级黄| av播播在线观看一区| 免费黄色在线免费观看| 国产有黄有色有爽视频| 日日爽夜夜爽网站| 秋霞在线观看毛片| 欧美日韩亚洲高清精品| 亚洲精品视频女| 亚洲国产色片| 久久人妻熟女aⅴ| 精品国产国语对白av| 精品少妇黑人巨大在线播放| 国国产精品蜜臀av免费| 成人综合一区亚洲| 中国三级夫妇交换| 美女视频免费永久观看网站| 男女啪啪激烈高潮av片| 国产精品熟女久久久久浪| 久久女婷五月综合色啪小说| 最近中文字幕高清免费大全6| 久久久欧美国产精品| 日韩中文字幕视频在线看片| 乱人伦中国视频| 男人操女人黄网站| 久久久精品区二区三区| 婷婷色麻豆天堂久久| 2021少妇久久久久久久久久久| 美女脱内裤让男人舔精品视频| 亚洲精品久久午夜乱码| 免费黄频网站在线观看国产| 草草在线视频免费看| 天天躁夜夜躁狠狠躁躁| 一区在线观看完整版| 另类精品久久| 在线看a的网站| 国语对白做爰xxxⅹ性视频网站| 国产伦理片在线播放av一区| av天堂久久9| 欧美变态另类bdsm刘玥| 日韩三级伦理在线观看| av在线app专区| 精品国产乱码久久久久久小说| 妹子高潮喷水视频| 汤姆久久久久久久影院中文字幕| 9热在线视频观看99| 午夜福利网站1000一区二区三区| www日本在线高清视频| 国产精品.久久久| 久久青草综合色| 国产在线视频一区二区| 黄色毛片三级朝国网站| 成人毛片a级毛片在线播放| 久久久精品区二区三区| 一本—道久久a久久精品蜜桃钙片| 97超碰精品成人国产| 精品卡一卡二卡四卡免费| 22中文网久久字幕| 天堂俺去俺来也www色官网| av片东京热男人的天堂| a级片在线免费高清观看视频| 精品国产国语对白av| 亚洲精品自拍成人| 亚洲美女搞黄在线观看| 丰满乱子伦码专区| 亚洲经典国产精华液单| 少妇 在线观看| 国产麻豆69| 久久精品aⅴ一区二区三区四区 | 国产片内射在线| 国产在线免费精品| 熟女人妻精品中文字幕| 精品国产一区二区久久| 国产成人午夜福利电影在线观看| 成人无遮挡网站| 亚洲精品久久久久久婷婷小说| 日日啪夜夜爽| 久久热在线av| 九九在线视频观看精品| 秋霞在线观看毛片| 十分钟在线观看高清视频www| 亚洲国产精品999| 又黄又粗又硬又大视频| 一级片免费观看大全| 成年av动漫网址| 欧美另类一区| 久久久久精品人妻al黑| 亚洲欧美中文字幕日韩二区| 中文字幕人妻熟女乱码| 一区二区日韩欧美中文字幕 | 热99久久久久精品小说推荐| 国产成人a∨麻豆精品| 精品福利永久在线观看| 日韩制服丝袜自拍偷拍| 国产高清国产精品国产三级| 在线天堂最新版资源| 国产无遮挡羞羞视频在线观看| 最近最新中文字幕免费大全7| 亚洲,欧美,日韩| 欧美日韩国产mv在线观看视频| 只有这里有精品99| 丁香六月天网| 亚洲精品久久午夜乱码| 天天操日日干夜夜撸| 精品99又大又爽又粗少妇毛片| 久久久精品94久久精品| 九九在线视频观看精品| 成人二区视频| 一级a做视频免费观看| 国产精品一区二区在线观看99| 欧美成人午夜精品| 亚洲精品色激情综合| 欧美日本中文国产一区发布| 熟女人妻精品中文字幕| 亚洲国产精品专区欧美| a级片在线免费高清观看视频| 99久久精品国产国产毛片| 美女主播在线视频| 日韩电影二区| 美女内射精品一级片tv| 天天躁夜夜躁狠狠久久av| av在线app专区| 日韩成人av中文字幕在线观看| 人成视频在线观看免费观看| 丝瓜视频免费看黄片| 久久久久人妻精品一区果冻| 亚洲欧美精品自产自拍| 精品99又大又爽又粗少妇毛片| 欧美精品av麻豆av| 亚洲成人av在线免费| 成人二区视频| 国产精品欧美亚洲77777| 黑丝袜美女国产一区| 超碰97精品在线观看| 国产男女内射视频| 精品国产一区二区三区久久久樱花| 久久99热6这里只有精品| av又黄又爽大尺度在线免费看| 看免费成人av毛片| 亚洲精品成人av观看孕妇| 亚洲国产精品专区欧美| 亚洲国产精品999| 王馨瑶露胸无遮挡在线观看| 精品福利永久在线观看| 亚洲五月色婷婷综合| 精品久久蜜臀av无| 亚洲国产精品999| 午夜91福利影院| av网站免费在线观看视频| 国产精品久久久久久久久免| 男女国产视频网站|