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

    大興安嶺諾敏河火山群遠(yuǎn)震P波衰減研究

    2021-12-30 08:28:38劉翰林吳慶舉
    地球物理學(xué)報(bào) 2021年1期
    關(guān)鍵詞:臺(tái)站火山反演

    劉翰林, 吳慶舉*

    1 中國(guó)地震局地球物理研究所, 北京 100081 2 中國(guó)地震局地球物理研究所, 地震觀測(cè)與地球物理成像重點(diǎn)實(shí)驗(yàn)室, 北京 100081

    0 引言

    中國(guó)東北地區(qū)位于中亞造山帶(Central Asian Orogenic Belt)東端,北部是西伯利亞克拉通,南部有華北克拉通,東部為西太平洋俯沖帶.該地區(qū)是一個(gè)盆地、造山帶和火山豐富發(fā)育的地質(zhì)構(gòu)造環(huán)境復(fù)雜區(qū)域,經(jīng)歷了晚二疊紀(jì)至早三疊紀(jì)時(shí)期與華北克拉通的碰撞(Eizenh?fer et al., 2014)、晚侏羅紀(jì)至早白堊紀(jì)時(shí)期與西伯利亞克拉通的碰撞(van Der Voo et al., 2015),以及中生代時(shí)期的斷裂運(yùn)動(dòng)(Hu et al., 1998; Feng et al., 2010; 葛榮蜂等,2010).中國(guó)東北在中生代斷裂運(yùn)動(dòng)期,火山活動(dòng)豐富發(fā)育,且在新生代時(shí)期,火山活動(dòng)延續(xù)并以松遼盆地為中心向盆地四周邊緣移動(dòng)(Liu et al., 2001).如長(zhǎng)白山火山、鏡泊湖火山、五大連池火山、哈拉哈河—淖爾河火山、錫林浩特—阿巴嘎火山等.這些火山群落屬于遠(yuǎn)離板塊邊界的板內(nèi)火山,從而使中國(guó)東北成為研究新生代板內(nèi)火山的理想場(chǎng)所.世界上多數(shù)火山活動(dòng)發(fā)生在板塊邊界上并能用板塊理論進(jìn)行解釋(Morgan, 1968),而板內(nèi)火山活動(dòng)卻有著更為復(fù)雜多樣的成因機(jī)制.近年來(lái),各學(xué)者對(duì)中國(guó)東北地區(qū)開(kāi)展了大量的地球物理研究并提出了系列板內(nèi)火山成因機(jī)制模型.遠(yuǎn)震體波成像研究(Lei and Zhao, 2005; Zhao, 2004; Huang and Zhao, 2006; Zhao, 2009; 雷建設(shè)等, 2018; Wei et al., 2018)在地幔過(guò)渡帶發(fā)現(xiàn)了大范圍的高速異常,推測(cè)為滯留在地幔過(guò)渡帶的俯沖太平洋板塊,并在長(zhǎng)白山等火山群下方發(fā)現(xiàn)了延伸到地幔過(guò)渡帶的低速異常,推測(cè)滯留俯沖板塊上方發(fā)生脫水熔融上涌導(dǎo)致了地表火山活動(dòng).環(huán)境噪聲、面波成像研究(Guo et al., 2016; Liu et al., 2017)則試圖用一個(gè)小尺度局部地幔對(duì)流模型來(lái)解釋各火山的成因,認(rèn)為長(zhǎng)白山下方來(lái)自轉(zhuǎn)換帶或更深的熱物質(zhì)上涌可能會(huì)引發(fā)松遼盆地西側(cè)小尺度的局部熱地幔對(duì)流上涌,從而導(dǎo)致阿巴嘎、阿爾山地區(qū)的火山活動(dòng).

    本文研究區(qū)位于中國(guó)東北北部,研究圍繞著隱匿在大興安嶺中的諾敏河火山展開(kāi)(圖1).有地球化學(xué)研究結(jié)果表明,諾敏河火山巖與五大連池火山巖的相似鉀質(zhì)屬性,體現(xiàn)了兩座火山群的親緣性,并與阿爾山及中國(guó)東部大多數(shù)的新生代鈉質(zhì)火山巖區(qū)別開(kāi)來(lái)(樊祺城等, 2012; 李霓等, 2012; 趙勇偉等, 2013).由于大興安嶺被大面積森林覆蓋,惡劣的野外環(huán)境導(dǎo)致地震觀測(cè)不便,諾敏河火山周邊的數(shù)據(jù)匱乏,東北地區(qū)眾多的地球物理學(xué)研究對(duì)該火山群并沒(méi)有很好的覆蓋,該地區(qū)深部地球物理結(jié)構(gòu)有待約束.2015年7月至2017年5月,中國(guó)地震局地球物理研究所在國(guó)家自然科學(xué)基金委的資助下圍繞諾敏河火山布設(shè)了43套寬頻地震儀進(jìn)行為期兩年的觀測(cè),臺(tái)間距約為30~50 km(圖1).在該項(xiàng)目的支持下,張風(fēng)雪和吳慶舉(2019)在該地區(qū)開(kāi)展遠(yuǎn)震體波成像研究,結(jié)果顯示諾敏河與五大連池下方200~300 km深度范圍有水平展布的低速異常,并認(rèn)為此低速異??赡苁莾蓚€(gè)火山群的共有地幔巖漿房.另外,SKS分裂研究在諾敏河火山下方得到多個(gè)無(wú)效分裂結(jié)果(強(qiáng)正陽(yáng)和吳慶舉,2019),接收函數(shù)研究結(jié)果表明諾敏河火山下方有較薄的地殼和較高的泊松比(謝振新等,2018),均顯示該火山區(qū)下方可能存在熱地幔上涌.

    新的地球物理觀測(cè)數(shù)據(jù)補(bǔ)充,為諾敏河火山地區(qū)的深部結(jié)構(gòu)提供了一定的有效約束,但仍需要新的研究手段,為以上提到的各向異性、地震波速度層析成像等傳統(tǒng)手段的研究成果提供更多的觀測(cè)輔證,以進(jìn)一步探索該地區(qū)板內(nèi)火山活動(dòng)乃至于整個(gè)東北地區(qū)新生代板內(nèi)火山活動(dòng)深部成因.地震波衰減在推斷地下結(jié)構(gòu)中溫度、孔隙流體、部分熔融、巖石組分以及主要構(gòu)造特征時(shí)扮演著十分重要的角色(Dziewonski et al., 1982),可豐富我們對(duì)深部結(jié)構(gòu)、火山構(gòu)造活動(dòng)及巖石圈演化的綜合理解(Adams and Humphreys, 2010; Zhao et al., 2010; He et al., 2017; Bezada, 2017; Byrnes et al., 2019).通常高衰減意味著高溫、高含水量的巖石成分,并與之對(duì)應(yīng)的是地震波低速異常、高地形的年輕造山帶或火山地區(qū),而低衰減則意味著低溫、低含水量的巖石成分,與之對(duì)應(yīng)的則是地震波高速異常、低地形的盆地或古老克拉通塊體.本研究中,我們利用上文提到的諾敏河流動(dòng)臺(tái)陣43套寬頻帶地震儀兩年時(shí)間的觀測(cè)數(shù)據(jù),挑選遠(yuǎn)震事件中的直達(dá)P波震相,采用Bezada(2017)提出的時(shí)間域波形匹配方法計(jì)算遠(yuǎn)震P波衰減因子,并反演了諾敏河火山及周邊地區(qū)的2D相對(duì)衰減結(jié)構(gòu)模型,為該地區(qū)的深部地球物理結(jié)構(gòu)提供了進(jìn)一步約束.

    圖1 研究區(qū)地質(zhì)概況及臺(tái)站分布 SP:西伯利亞板塊;BR:貝加爾裂谷;CAOB:中亞造山帶;NCC:華北克拉通;HB:海拉爾盆地;SLB:松遼盆地;GXR:大興安嶺; LXR:小興安嶺;NSGL:南北重力梯度帶;紅色三角形為火山位置;黃色正方形為臺(tái)站位置;灰色粗實(shí)線為南北重力梯度線.Fig.1 Map of tectonic setting and the distribution of seismic stations SP: Siberian Plate; BR: Baikal Rift; CAOB: Central Asian Orogenic Belt; NCC: North China Craton; HB: Hailar Basin; SLB: Songliao Basin; GXR: Da Hinggan; LXR: Xiao Hinggan; NSGL: North-South Gravity Line; red triangles are volcanoes; yellow squares represent seismic stations; thick grey solid line is North-South Gravity Line.

    1 數(shù)據(jù)與方法

    地震波衰減主要由巖石黏滯性導(dǎo)致的本征吸收以及巖石地震波速度橫向不均勻性導(dǎo)致的散射衰減組成(Akinci et al., 1995),可以用品質(zhì)因子參數(shù)Q(Knopoff, 1964)來(lái)量化,定義如下:

    (1)

    其中,ΔE是地震波每個(gè)周期的能量變化,Emax為每個(gè)周期中的最大彈性能量.該參數(shù)與地震波信號(hào)的關(guān)系(Aki and Richards, 2002)定義如下:

    (2)

    式中A(ω)是信號(hào)頻譜,A0(ω)為未被衰減的信號(hào)頻譜,ω為角頻率,x和c(ω)分別為傳播距離和相速度.

    以上為地震波衰減和品質(zhì)因子Q的基本定義,通常我們用于地震波衰減分析中的是衰減因子參數(shù)t*,該參數(shù)是一個(gè)單位為秒(s)的量化參數(shù),代表全射線路徑上所有巖石Q值對(duì)地震波的作用累積的衰減總和,定義如下(Stein and Wysession, 2003):

    (3)

    式中V與Q分別代表射線路徑所經(jīng)過(guò)的速度模型和品質(zhì)因子模型,dt為走時(shí)微分.

    最早由Teng(1968)提出譜比法,將兩個(gè)臺(tái)站的振幅譜在去儀器響應(yīng)后作比,以消除震源項(xiàng)和近震源衰減項(xiàng),從而得到兩個(gè)臺(tái)站下方的Δt*.后人基于譜比法提出了共譜法(Halderman and Davis,1991)、虛擬源法(Boyd and Sheehan,2005)以及參考譜法(Adams and Humphreys,2010)等方法,不再選取單一的參考臺(tái)站,而通過(guò)尋找平均參考譜,同時(shí)計(jì)算各個(gè)臺(tái)站之間的Δt*.頻率域方法在衰減研究中應(yīng)用廣泛,但仍存在一定局限性.頻譜計(jì)算和譜比計(jì)算對(duì)窗口選擇較為敏感(Solomon,1973),且散射和聚焦效應(yīng)對(duì)地震波形振幅的影響具有頻率依賴性(Cafferky and Schmandt,2015;Dalton et al., 2008).本研究中,采用由Bezada(2017)提出的時(shí)間域波形匹配搜索方法來(lái)探索研究區(qū)下方的本征衰減特征,該方法與Adams 和 Humphreys (2010)在黃石火山地區(qū)應(yīng)用的時(shí)間域衰減計(jì)算方法類似,并在摩洛哥、中阿巴拉契亞山脈等地區(qū)的衰減研究應(yīng)用中,均取得可靠的研究結(jié)果(Bezada,2017; Byrnes et al.,2019).

    1.1 數(shù)據(jù)挑選

    選取30°~90°的遠(yuǎn)震事件數(shù)據(jù).基于巖石圈與軟流圈Q值顯著差異(Dziewonski and Anderson, 1981),為了盡量減小接收信號(hào)里所包含的近震源衰減項(xiàng),使其所包含的衰減信息主要來(lái)自于臺(tái)站下方,我們選取震源深度在250 km以下的事件.為了得到高信噪比的波形數(shù)據(jù),事件震級(jí)下限被設(shè)為Mb5.5,同時(shí)為了減小震源方向性效應(yīng)對(duì)衰減計(jì)算結(jié)果的影響,得到脈沖性的直達(dá)P波震相,震級(jí)上限被設(shè)為Mb7.0(Hwang et al., 2011).綜上,我們共計(jì)挑選了17個(gè)遠(yuǎn)震事件(圖2).通過(guò)去均值、去趨勢(shì)、去儀器響應(yīng)等數(shù)據(jù)處理后,對(duì)波形進(jìn)行0.02~3 Hz的帶通濾波,挑選地震記錄垂向分量中SNR(信噪比)大于3、波形相似并具有脈沖性的直達(dá)P波震相進(jìn)行Δt*的計(jì)算.

    圖2 事件分布圖 藍(lán)色五角星為地震臺(tái)陣位置,彩色圓圈為地震事件位置. 其中圓圈顏色尺度為事件震源深度,半徑大小為事件震級(jí).Fig.2 Distribution of teleseismic events used in present paper Blue star denotes array, and colorful circles are events. Color and radius of the circles scale the focus depth and magnitude respectively.

    1.2 波形匹配與Δt*搜索

    時(shí)間域波形匹配法第一步是求取衰減為零或衰減極弱的震源時(shí)間函數(shù)(Bezada, 2017).如圖3a所示,各臺(tái)站接收到的相似的波形中,衰減顯著的波形因更多的高頻成分被衰減導(dǎo)致脈沖時(shí)間變長(zhǎng)而呈現(xiàn)為較寬的震相,衰減微弱的波形則反之,呈現(xiàn)為較窄的震相.我們選取被衰減最少的幾個(gè)臺(tái)站波形作為源波形(圖3),將其疊加以減小尾波能量,得到參照源(Bezada, 2017).

    第二步為生成理論波形并搜索Δt*值(Bezada, 2017).在此,我們引入Azimi等(1968)提出的衰減項(xiàng):

    (4)

    式中i為復(fù)數(shù)的虛部單位,ω0為參照頻率常數(shù).ω0的選取只改變理論波形的到時(shí)并不改變其形狀,故ω0的選取對(duì)衰減因子計(jì)算結(jié)果并無(wú)影響(Bezada, 2017).通過(guò)對(duì)參照源添加衰減項(xiàng)來(lái)合成理論波形,給定Δt*越大脈沖越寬,Δt*越小脈沖越窄(圖4).將理論與觀測(cè)波形作振幅標(biāo)準(zhǔn)化處理后再進(jìn)行波形擬合,計(jì)算波形振幅差的L2范數(shù)作為擬合偏差,同時(shí)對(duì)Δt*進(jìn)行步長(zhǎng)為0.01 s的網(wǎng)格搜索(圖4),最終選取使得擬合偏差最小的Δt*作為搜索結(jié)果(Bezada, 2017).圖4展示的擬合結(jié)果為歸一化處理后的波形.

    最后一步即對(duì)擬合震相進(jìn)行質(zhì)量評(píng)估,保留擬合效果較好的結(jié)果.為盡可能保留更多的高質(zhì)量觀測(cè)數(shù)據(jù),我們并未采用統(tǒng)一的評(píng)估閾值對(duì)結(jié)果進(jìn)行自動(dòng)篩選,而對(duì)搜索結(jié)果進(jìn)行人工視覺(jué)篩選(圖4).最終,我們保留了300個(gè)Δt*搜索結(jié)果,并將其分別對(duì)每個(gè)事件進(jìn)行去均值處理,再用作成像反演.

    1.3 反演

    我們利用觀測(cè)Δt*數(shù)據(jù)反演研究區(qū)衰減結(jié)構(gòu)的橫向變化(Hwang et al., 2011; Lawrence et al.,2006; Bezada, 2017).考慮到觀測(cè)數(shù)據(jù)針對(duì)于每個(gè)事件去均值處理,但每個(gè)事件均值并不相等,另外對(duì)于每個(gè)臺(tái)站其近地表效應(yīng)也存在差異,則觀測(cè)到的Δt*為衰減結(jié)構(gòu)差異、事件均值差異及臺(tái)站近地表效應(yīng)差異的作用總和,于是有(Bezada, 2017):

    表1 圖3、4及7中所展示四個(gè)事件的詳細(xì)信息Table 1 Details of 4 events presented in Fig.3, 4 and 7

    圖3 源波形選取及參照源疊加示意圖(展示事件的信息詳見(jiàn)表1中Evt.1) (a) 選取源波形示意; (b) 疊加生成參照源及選取擬合窗口示意.Fig.3 Source Traces and Reference Source (Details of presented event are shown in Table 1-Evt.1) (a) Shows the selection of Source Traces; (b) Shows the stacking and selection of fitting window on the Reference Source.

    圖4 波形擬合及Δt*質(zhì)量評(píng)估示意圖 圖中展示的計(jì)算結(jié)果來(lái)自于圖3中同一事件.圖中黑色實(shí)線為觀測(cè)波形; (a) 中紅色實(shí)線; (b) 中紅色虛線為理論合成波形; (a) 中紅色實(shí)線部分為擬合時(shí)窗;XM31臺(tái)站波形為被篩除的擬合效果較差波形.Fig.4 Waveform matching and Δt* quality control The presented event is the same with it in Fig.3. Black solid lines represent observed waveforms; red solid lines in panel (a) and red dashed lines in panel (b) are synthetic waveforms; red solid lines in panel (b) represent fitting window; result of station XM31 is removed due to the bad matching of waveforms.

    圖5 理論測(cè)試輸入模型及測(cè)試結(jié)果 圖中圓圈為臺(tái)站.(a1—3)分別為三個(gè)理論測(cè)試輸入模型; (b1—3) 分別展示了對(duì)三個(gè)輸入模型加入高斯隨機(jī)噪聲, 并按事件去均值處理后再按臺(tái)站疊加求均值的Δt*結(jié)果; (c1—3) 分別為三個(gè)理論測(cè)試輸出結(jié)果.Fig.5 Input and output models of three synthetic tests Circles in the figure denote stations. Panels of (a1—3) show the three input models respectively; panels of (b1—3) show the average Δt* values per station taken from the input models, with the mean per event removed and Gaussian noise added; panels of (c1—3) represent three output models respectively.

    圖6 二維Δt*反演結(jié)果模型及模型標(biāo)準(zhǔn)差 (a) 二維Δt*模型,圖中黑色圓點(diǎn)為臺(tái)站;紫色三角形為火山;HA1與HA2分別為兩個(gè)高衰減觀測(cè)區(qū)域; LA1為低衰減觀測(cè)區(qū)域;兩條藍(lán)色線段為圖8中兩條剖面所在位置.(b)模型標(biāo)準(zhǔn)差,圖中黑色圓點(diǎn)為臺(tái)站.Fig.6 Maps of 2D Δt* model and model standard deviation (a) Black dots are stations and purple triangles are volcanoes; HA1 and HA2 represent two regions of higher Δt* observations; LA1 represent the region of lower Δt* observations; blue lines show the locations of vertical profiles in Fig.8; (b) The model standard deviation and the black dots are stations.

    圖7 三個(gè)事件Δt*搜索結(jié)果及部分臺(tái)站波形擬合展示(三個(gè)事件的信息詳見(jiàn)表1中Evt.2、Evt.3及Evt.4) 左側(cè)為每個(gè)事件各臺(tái)站Δt*的計(jì)算結(jié)果,黑色空心圓代表未接收到該事件信號(hào)或波形被篩除的臺(tái)站; 右側(cè)展示了左側(cè)圖中標(biāo)注臺(tái)站的波形擬合結(jié)果.Fig.7 Map of observations of Δt* values from three events and waveforms matching from some of the stations(Details of the three events are shown in Table 1-Evt.2, Evt.3 and Evt.4.) Figures on the left side show the Δt* of stations estimated from each event, and black circles are stations which did not receive the signal from the event or whose searching results are removed due to the bad matching; figures on the right side are the waveform-matching results from stations pointed out in the left figures.

    圖8 Δt*與速度模型對(duì)比垂直剖面圖 (a)和(b)兩條剖面位置如圖6所示.(a)、(b)從上至下依次為地形海拔、Δt*值、Tao 等(2018)速度模型的P波速度擾動(dòng)、S波速度擾動(dòng)、 張風(fēng)雪和吳慶舉(2019)速度模型的P波速度擾動(dòng).圖中粉色三角為火山.NM為諾敏河火山,其余縮寫(xiě)簡(jiǎn)稱均與圖1、6中一致.Fig.8 Vertical profiles of the comparison between Δt* and velocity Locations of sections (a) and (b) are showed in Fig.6. From top to bottom, panels (a) and (b) show the topography elevation, Δt* values, Vp and Vs perturbations from model of Tao et al. (2018), Vp perturbations from model of Zhang and Wu (2019). Pink triangle is volcano in the figure. NM represents Nuomin River volcano, and other abbreviations are the same with those in theFig.1,6.

    GΔm+e+s=Δd,

    (5)

    式中Δd為觀測(cè)數(shù)據(jù),e為獨(dú)立待求解參數(shù)事件項(xiàng),代表事件均值差異,s為獨(dú)立待求解參數(shù)臺(tái)站項(xiàng),代表臺(tái)站近地表效應(yīng)差異(Bezada, 2017).GΔm代表固有衰減差異,Δm為二維相對(duì)衰減結(jié)構(gòu),G為距離加權(quán)因子,定義如下(Bezada, 2017):

    (6)

    其中,i為觀測(cè)數(shù)據(jù)序號(hào),j為模型網(wǎng)格序號(hào),k為臺(tái)站序號(hào),rjk為第j個(gè)網(wǎng)格到第k個(gè)臺(tái)站的距離,rmax為預(yù)設(shè)最大距離(Bezada, 2017).

    通常衰減研究中,波形擬合偏差或頻譜擬合偏差被作為觀測(cè)Δt*不確定性,用作成像反演中的權(quán)重因子(Boyd and Sheehan, 2005; Adams and Humphreys, 2010).本研究中,我們采用非線性貝葉斯-蒙特卡洛反演方法解上述問(wèn)題,將觀測(cè)數(shù)據(jù)不確定性作為獨(dú)立待求解參數(shù)進(jìn)行反演,并且相比于Bezada(2017)采用的線性最小二乘解法,該方法不需要選取平滑和阻尼來(lái)控制模型變化,針對(duì)于觀測(cè)數(shù)據(jù)量稀少的反演問(wèn)題亦不會(huì)導(dǎo)致過(guò)載的計(jì)算量(Byrnes et al., 2019).貝葉斯-蒙特卡洛反演方法由Bodin等(2012a,b)引入面波成像研究,Byrnes等(2019)將其應(yīng)用于衰減成像研究,獲得了中阿巴拉契亞山脈下方的上地幔衰減結(jié)構(gòu),并有效地估計(jì)了觀測(cè)數(shù)據(jù)的噪聲水平.

    貝葉斯-蒙特卡洛反演通過(guò)Voronoi圖隨機(jī)生成初始模型m(模型的網(wǎng)格數(shù)量、網(wǎng)格空間分布、網(wǎng)格賦值)以及數(shù)據(jù)不確定性σd,基于貝葉斯理論,該模型的概率函數(shù)如下:

    P∝p(m|σd)p(d|m,σd),

    (7)

    式中p(m|σd)為模型的先驗(yàn)概率,p(d|m,σd)為評(píng)估模型預(yù)測(cè)數(shù)據(jù)與觀測(cè)數(shù)據(jù)擬合質(zhì)量的概率函數(shù)(Bodin et al., 2012a,b; Byrnes et al., 2019).通過(guò)對(duì)初始模型進(jìn)行如下隨機(jī)擾動(dòng)得到新模型m′:1.移動(dòng):隨機(jī)選取一個(gè)模型網(wǎng)格改變其空間位置;2.改變:隨機(jī)選取一個(gè)模型網(wǎng)格改變其模型值;3.σd:隨機(jī)改變數(shù)據(jù)不確定性σd的大小;4.誕生:隨機(jī)生成一個(gè)新的網(wǎng)格并賦值;5.消亡:隨機(jī)刪除一個(gè)模型網(wǎng)格(Bodin et al., 2012a,b; Byrnes et al., 2019).然后對(duì)新模型m′與初始模型m之間構(gòu)造如下判定函數(shù):

    (8)

    式中q(m|m′)代表從新模型跳躍到舊模型的發(fā)生概率,q(m′|m)則反之.反演通過(guò)馬爾可夫鏈(Markov Chain)進(jìn)行迭代,可多條鏈平行迭代.每條鏈中當(dāng)新模型被判定函數(shù)接受,則將新模型作為新的起點(diǎn)進(jìn)行下一步迭代;若否,則將舊模型進(jìn)行新的隨機(jī)擾動(dòng),直到被接受(Bodin et al., 2012a, b; Byrnes et al., 2019).若數(shù)據(jù)擬合偏差不顯著增大,迭代將使得反演結(jié)果向更簡(jiǎn)單的模型逼近,而當(dāng)更復(fù)雜的模型不能顯著改善數(shù)據(jù)擬合情況時(shí),數(shù)據(jù)不確定性將會(huì)向某一確定值逼近(Bodin et al., 2012a,b; Byrnes et al., 2019).

    1.4 理論測(cè)試

    為驗(yàn)證反演方法可行性,以及觀測(cè)數(shù)據(jù)對(duì)研究區(qū)橫向衰減差異的恢復(fù)能力,我們給定三個(gè)不同的初始模型分別進(jìn)行理論測(cè)試.初始模型一:Δt*值在121°E—126°E范圍內(nèi)隨經(jīng)度線性增加,由0至0.4 s;Δt*值在48°N—52°N范圍內(nèi)隨緯度線性增加,由0至0.4 s;初始模型三:以經(jīng)度124°E為界,西區(qū)Δt*為0 s,東區(qū)Δt*為0.4 s(圖5a1—3).

    理論Δt*數(shù)據(jù)均依從觀測(cè)數(shù)據(jù)集中的臺(tái)站-事件對(duì)關(guān)系進(jìn)行輸入,并加入20%噪聲(標(biāo)準(zhǔn)差為0.08 s高斯隨機(jī)噪聲),對(duì)理論數(shù)據(jù)按每個(gè)事件進(jìn)行去均值處理(圖5b1—3).如圖5b1—3所示,對(duì)于按事件去均值處理后的Δt*,若直接按臺(tái)站疊加求均值結(jié)果與初始輸入模型存在一定偏差,不能很好反映真實(shí)Δt*橫向擾動(dòng).

    測(cè)試結(jié)果顯示(圖5c1—3),我們使用的觀測(cè)數(shù)據(jù)和反演方法對(duì)隨經(jīng)、緯度增加的線性梯度與分區(qū)突變的陡峭梯度均有較好的形狀恢復(fù),但橫向衰減擾動(dòng)幅度較輸入值有細(xì)微差距,三個(gè)測(cè)試反演結(jié)果中分別得到的三個(gè)數(shù)據(jù)不確定性為0.0815 s、0.0777 s、0.0733 s,與輸入的隨機(jī)噪聲值大致相符.綜上,本研究中采用的反演方法及觀測(cè)數(shù)據(jù)集可獲得較為可靠的研究區(qū)橫向衰減差異.

    2 結(jié)果與討論

    對(duì)300條Δt*觀測(cè)數(shù)據(jù)進(jìn)行貝葉斯-蒙特卡洛反演,采用10條平行鏈進(jìn)行迭代,每條鏈迭代1×106次,最后將保留模型疊加取平均得到反演結(jié)果(圖6a),Δt*橫向變化大致在-0.15~0.15 s范圍內(nèi).并將反演模型概率分布標(biāo)準(zhǔn)差作為模型不確定性(圖6b),模型不確定性在臺(tái)站覆蓋區(qū)域較好處不超過(guò)0.05 s.反演得到數(shù)據(jù)不確定性均值為0.08 s(與模型不確定性不同,數(shù)據(jù)不確定性為反演得到的觀測(cè)數(shù)據(jù)噪聲水平).

    2.1 相對(duì)衰減結(jié)構(gòu)基本特征

    如圖6所示,研究區(qū)的相對(duì)衰減表現(xiàn)出較大幅度橫向擾動(dòng).通常來(lái)說(shuō),火山活動(dòng)與巖石圈減薄、熱地幔上涌、部分熔融等現(xiàn)象聯(lián)系緊密,并伴隨著低速異常、高衰減異常的觀測(cè)(Shapiro et al., 2000).本研究區(qū)的火山活動(dòng)與我們的高衰減觀測(cè)結(jié)果呈現(xiàn)較好的空間一致性(圖6).諾敏河火山大致位于高衰減區(qū)域HA1的中心,Δt*最高達(dá)約0.15 s,且Δt*圍繞著火山群向四周逐漸減小(圖6、8).另一個(gè)高衰減區(qū)域HA2位于研究區(qū)最北端,Δt*最高值與HA1接近(圖6、8),主要表現(xiàn)在XM01與XM02這兩個(gè)臺(tái)站的觀測(cè)(圖7),可能與研究區(qū)北部空間位置接近的小古里河火山活動(dòng)有關(guān).位于諾敏河火山東部、科洛河火山西部,緊鄰兩個(gè)火山群的LA1區(qū)域觀測(cè)到顯著低衰減異常,Δt*最低至約-0.15 s(圖6、8),主要表現(xiàn)為XM17、XM10、XM06、XM04這四個(gè)臺(tái)站的觀測(cè)(圖7).

    選取三個(gè)事件的Δt*觀測(cè)結(jié)果成圖,并展示分別位于高衰減區(qū)域和低衰減區(qū)域的部分臺(tái)站波形(圖7).如圖7所示,三個(gè)事件分別在三個(gè)區(qū)域中的Δt*觀測(cè)結(jié)果與我們的反演結(jié)果基本相符.位于HA1的XM12、XM15、XM31、XM40臺(tái)站與位于HA2的XM01、XM02臺(tái)站具有較大的Δt*值及較寬的脈沖波形,而位于LA1的XM04、XM06、XM10、XM17臺(tái)站觀測(cè)到較小的Δt*值及較窄的脈沖波形.

    2.2 與前人結(jié)果對(duì)比及討論

    通常遠(yuǎn)震射線中,下地幔的衰減遠(yuǎn)小于上地幔(Booth et al., 1974),并且對(duì)比Zhao等(2010)對(duì)中國(guó)東北地區(qū)的Lg波衰減研究,其結(jié)果在本文覆蓋區(qū)域并未顯示如此明顯大幅的橫向變化,由于Lg波衰減主要受地殼影響,我們認(rèn)為本研究所得到的橫向衰減變化主要來(lái)自于軟流圈和巖石圈地幔的固有衰減.

    對(duì)二維Δt*模型作兩條垂直剖面(a)、(b),橫穿HA1、HA2及LA1,與前人速度觀測(cè)結(jié)果進(jìn)行對(duì)比(圖8a—b).圖8a—b中分別采用Tao等(2018)通過(guò)全波形反演得到的東亞地區(qū)P、S波速度模型FWEA18(Full Waveform inversion of East Asia in 2018)以及張風(fēng)雪和吳慶舉(2019)的遠(yuǎn)震P波成像結(jié)果在兩條剖面處的速度異常成圖.Tao等(2018)及張風(fēng)雪和吳慶舉(2019)的速度模型在兩條剖面下方均顯示分層低速異常:100~200 km深度低速異常及200~400 km低速異常,橫向貫穿HA1、HA2及LA1衰減異常區(qū)域,速度異常幅度略有差異,最高處達(dá)約4%.另在Tao等(2018)模型中LA1、HA2區(qū)域下方,淺層低速異常上覆厚度為約100 km高速異常,而該高速異常在張風(fēng)雪和吳慶舉(2019)的模型中只表現(xiàn)在LA1區(qū)域下方(圖8a).

    通常來(lái)說(shuō),溫度是影響上地幔巖石礦物速度的主要因素,低速異??赡馨凳局邷禺惓?,而高速異常可能代表較冷的構(gòu)造體(Karato,1993; Cammarano et al., 2003; Goes et al., 2000).假設(shè)只受熱力學(xué)條件控制,溫度對(duì)Q與對(duì)速度的影響類似,因此衰減研究中通常將速度與衰減的負(fù)相關(guān)現(xiàn)象主要?dú)w因于溫度因素(Lawrence et al., 2006; Hwang et al., 2011).例如Hwang等(2011)的全球體波衰減研究中觀測(cè)到在古老的、冷的大陸克拉通地區(qū)的低衰減、高速異常,以及熱的、年輕的洋中脊地區(qū)的高衰減、低速異常.

    本研究中HA1、HA2下方觀測(cè)到高衰減異常,對(duì)應(yīng)著前人速度模型中同區(qū)域下方上地幔兩層低速異常,暗示著該區(qū)域下方存在高溫異常.另外張風(fēng)雪和吳慶舉(2019)提出諾敏河火山下方兩層低速異??赡芘c地表火山活動(dòng)有關(guān),且認(rèn)為深層低速異常為諾敏河與五大連池火山群的共有地幔巖漿房,可能與晚中生代巖石圈拆沉事件導(dǎo)致的軟流圈熱物質(zhì)上涌有關(guān).強(qiáng)正陽(yáng)和吳慶舉(2019)在該研究區(qū)的SKS分裂研究顯示,HA1區(qū)域和HA2區(qū)域均觀測(cè)到純無(wú)效分裂結(jié)果,并歸因于熱地幔物質(zhì)的上涌侵蝕了巖石圈殘存化石各向異性.謝振新等(2018)在該研究區(qū)的接收函數(shù)研究顯示,HA1區(qū)域和HA2區(qū)域均具有較薄的地殼厚度與較高的泊松比,并認(rèn)為這與富鉀巖石圈拆沉觸發(fā)的熱物質(zhì)上涌有關(guān).綜上,熱地幔物質(zhì)上涌可能是HA1、HA2區(qū)域觀測(cè)到高衰減異常的主要原因.

    LA1下方觀測(cè)到顯著低衰減異常,前人速度模型在該區(qū)域下方有厚度不一(約100 km)的高速異常體(Tao et al., 2018; 張風(fēng)雪和吳慶舉,2019),且強(qiáng)正陽(yáng)和吳慶舉(2019)在LA1區(qū)域中XM17臺(tái)站觀測(cè)到延遲時(shí)間大于1.0 s的SKS分裂結(jié)果,暗示其下方可能有未直接受到熱地幔物質(zhì)上涌侵蝕而殘存巖石圈.Zhang等(2014)的S波接收函數(shù)結(jié)果顯示大興安嶺下方巖石圈厚度約為140~160 km,據(jù)PREM模型(Dziewonski and Anderson, 1981)中巖石圈、軟流圈的Qp值進(jìn)行衰減估算,150 km巖石圈厚度差異導(dǎo)致的Δt*約為0.08 s,盡管該估算具有一定不確定性,但其結(jié)果僅約為該研究中衰減觀測(cè)差值的三分之一(約0.3±0.05 s),因而巖石圈厚度差異并不足以解釋LA1區(qū)域低衰減觀測(cè).另外LA1區(qū)域下方同樣對(duì)應(yīng)著兩層低速異常(Tao et al., 2018; 張風(fēng)雪和吳慶舉,2019),LA1區(qū)域較之HA1區(qū)域具有更薄的地殼和更高的泊松比(謝振新等,2018),且SKS分裂研究顯示LA1區(qū)域亦觀測(cè)到純無(wú)效分裂結(jié)果(強(qiáng)正陽(yáng)和吳慶舉,2019).這暗示著LA1區(qū)域下方同樣可能存在熱地幔物質(zhì)上涌帶來(lái)的高溫異常,那是否與本研究中在該區(qū)域下方觀測(cè)到的顯著低衰減異常矛盾呢?

    除卻溫度因素外,含水量是影響巖石礦物Q值的另一重要因素.部分巖石物理研究結(jié)果表明,上地幔條件下含水量(主要以氫離子形態(tài))對(duì)Q值的影響能達(dá)到75%以上,對(duì)速度影響卻很小(Karato,1995,2003; Karato and Jung,1998).盡管含水量對(duì)Q值影響的程度至今仍尚存爭(zhēng)議,但近年來(lái)觀測(cè)到低速異常對(duì)應(yīng)著低衰減異常的現(xiàn)象,已不能基于熱力學(xué)條件進(jìn)行很好地解釋了(Hammond and Humphreys, 2000; Humphreys et al., 2003; Lawrence et al., 2006).例如Yang等(2007)以及Yang和Forsyth(2008)分別在東太平洋隆起和南加州地區(qū)的衰減研究中發(fā)現(xiàn),若只基于熱力學(xué)條件控制,通過(guò)速度異常預(yù)測(cè)的Q值比觀測(cè)Q值偏小,而加入含水量因素后擬合效果得到顯著改善.Lawrence等(2006)與Adams和Humphreys(2010)相繼在黃石火山下方觀測(cè)到顯著低衰減異常,認(rèn)為部分熔融相的介入會(huì)導(dǎo)致地幔脫水,從而使得殘留地幔橄欖巖Q值增加,衰減減小.Bezada(2017)在摩洛哥火山活動(dòng)地區(qū)下方觀測(cè)到低衰減異常,亦歸因于部分熔融導(dǎo)致的殘留地幔礦物脫水.類似地,我們推測(cè)LA1區(qū)域下方可能存在著部分熔融導(dǎo)致的殘留地幔礦物脫水,進(jìn)而導(dǎo)致衰減減小.

    綜上,LA1區(qū)域觀測(cè)到的低衰減異??赡苡蓛蓚€(gè)原因共同導(dǎo)致,一方面是該區(qū)域下方殘留未被侵蝕的巖石圈,另一方面是局部區(qū)域可能存在部分熔融導(dǎo)致的殘留地幔礦物脫水.

    3 結(jié)論

    利用圍繞著大興安嶺諾敏河火山布設(shè)的43個(gè)流動(dòng)臺(tái)站所接收到的17個(gè)遠(yuǎn)震事件,我們采用時(shí)間域波形匹配方法計(jì)算了300條直達(dá)P波震相的Δt*,并采用貝葉斯蒙特卡洛方法進(jìn)行反演,獲得到了中國(guó)東北諾敏河火山及周邊地區(qū)的二維Δt*模型.結(jié)果表明研究區(qū)高衰減區(qū)域與火山活動(dòng)位置分布呈現(xiàn)較好的空間一致性.熱地幔物質(zhì)上涌可能是導(dǎo)致諾敏河火山區(qū)域以及北部臨近小古里河火山區(qū)域高衰減異常的主要原因.諾敏河火山東部、科洛河火山西部區(qū)域呈觀測(cè)到顯著低衰減異常,與高衰減區(qū)域最大Δt*差值達(dá)約0.3±0.05 s,可能是兩個(gè)原因共同作用的結(jié)果:該區(qū)域下方殘留未被侵蝕的較厚巖石圈,以及局部區(qū)域可能存在的部分熔融導(dǎo)致的殘留地幔礦物脫水.

    致謝感謝國(guó)家自然科學(xué)基金項(xiàng)目(41674094,41474074)的資助.感謝參與該項(xiàng)目野外流動(dòng)臺(tái)站架設(shè)、維護(hù)及數(shù)據(jù)收集的全體人員.感謝審稿專家的寶貴修改意見(jiàn),使得本文更加嚴(yán)謹(jǐn)、流暢.

    猜你喜歡
    臺(tái)站火山反演
    中國(guó)科學(xué)院野外臺(tái)站檔案工作回顧
    反演對(duì)稱變換在解決平面幾何問(wèn)題中的應(yīng)用
    氣象基層臺(tái)站建設(shè)
    西藏科技(2021年12期)2022-01-17 08:46:38
    海底火山群
    有趣的火山圖
    火山
    基于低頻軟約束的疊前AVA稀疏層反演
    基于自適應(yīng)遺傳算法的CSAMT一維反演
    我是火山
    基層臺(tái)站綜合觀測(cè)業(yè)務(wù)管理之我見(jiàn)
    西藏科技(2015年6期)2015-09-26 12:12:13
    听说在线观看完整版免费高清| 人妻丰满熟妇av一区二区三区| 午夜福利在线观看吧| 亚洲自偷自拍三级| 国产乱人视频| 国产伦在线观看视频一区| 国产精品伦人一区二区| 嫁个100分男人电影在线观看| 色综合站精品国产| 精品久久久久久久久亚洲 | 男女那种视频在线观看| 久久精品综合一区二区三区| 真人做人爱边吃奶动态| 深爱激情五月婷婷| 少妇高潮的动态图| 亚洲成av人片免费观看| 久久久久久九九精品二区国产| 桃色一区二区三区在线观看| 一区福利在线观看| 又黄又爽又刺激的免费视频.| 亚洲av日韩精品久久久久久密| 91狼人影院| 久久久久九九精品影院| 日韩精品青青久久久久久| 国产成人啪精品午夜网站| 一本久久中文字幕| 欧美一区二区亚洲| 88av欧美| 搞女人的毛片| 亚洲人成网站在线播| 18禁裸乳无遮挡免费网站照片| 99热只有精品国产| 97碰自拍视频| 韩国av一区二区三区四区| 可以在线观看的亚洲视频| 国产亚洲欧美在线一区二区| 日本五十路高清| 国产伦在线观看视频一区| 中文在线观看免费www的网站| av视频在线观看入口| 日韩欧美一区二区三区在线观看| 成年版毛片免费区| 婷婷丁香在线五月| 两个人视频免费观看高清| 最新中文字幕久久久久| 亚洲经典国产精华液单 | 国产精品久久久久久久电影| 久久久久久久亚洲中文字幕 | 欧美区成人在线视频| 人人妻人人看人人澡| 国内久久婷婷六月综合欲色啪| 欧美色欧美亚洲另类二区| 国产精品嫩草影院av在线观看 | 中文字幕人成人乱码亚洲影| 一区二区三区高清视频在线| 网址你懂的国产日韩在线| 色5月婷婷丁香| 欧美潮喷喷水| 一本精品99久久精品77| 男人狂女人下面高潮的视频| 日韩亚洲欧美综合| 性插视频无遮挡在线免费观看| 国产精品自产拍在线观看55亚洲| 丰满人妻一区二区三区视频av| 一级a爱片免费观看的视频| 精品人妻偷拍中文字幕| 国产精品亚洲一级av第二区| 国产精品乱码一区二三区的特点| 真人做人爱边吃奶动态| 91av网一区二区| 亚洲欧美日韩东京热| 国产精品电影一区二区三区| 直男gayav资源| 校园春色视频在线观看| 欧美中文日本在线观看视频| 非洲黑人性xxxx精品又粗又长| 国产亚洲精品久久久com| 丝袜美腿在线中文| 日韩欧美免费精品| 少妇的逼好多水| 国产精品一区二区三区四区久久| 九九久久精品国产亚洲av麻豆| 日本与韩国留学比较| 日本免费一区二区三区高清不卡| 亚洲中文字幕一区二区三区有码在线看| 亚洲国产精品合色在线| 亚洲国产日韩欧美精品在线观看| 脱女人内裤的视频| 亚洲在线自拍视频| 精品久久久久久久久久久久久| 成年女人毛片免费观看观看9| 亚洲成a人片在线一区二区| 国产成人aa在线观看| 亚洲av免费高清在线观看| 午夜福利视频1000在线观看| 色精品久久人妻99蜜桃| 毛片一级片免费看久久久久 | 色噜噜av男人的天堂激情| 丁香欧美五月| 色在线成人网| 亚洲真实伦在线观看| 人人妻,人人澡人人爽秒播| 男人舔奶头视频| 中文字幕av成人在线电影| 国产高潮美女av| 小说图片视频综合网站| 在线观看66精品国产| 白带黄色成豆腐渣| 又黄又爽又免费观看的视频| 少妇熟女aⅴ在线视频| 中文字幕av在线有码专区| 国产在线男女| 国产探花在线观看一区二区| 精品人妻一区二区三区麻豆 | 亚洲 欧美 日韩 在线 免费| av福利片在线观看| 成人国产综合亚洲| 少妇熟女aⅴ在线视频| 亚洲欧美日韩东京热| 757午夜福利合集在线观看| 人妻制服诱惑在线中文字幕| 日本熟妇午夜| 又黄又爽又免费观看的视频| 少妇高潮的动态图| 婷婷精品国产亚洲av| 在线观看66精品国产| 欧美日韩综合久久久久久 | 亚洲天堂国产精品一区在线| 精品人妻一区二区三区麻豆 | 91午夜精品亚洲一区二区三区 | 一进一出抽搐gif免费好疼| 精品人妻偷拍中文字幕| 日韩精品青青久久久久久| 99国产极品粉嫩在线观看| 日本与韩国留学比较| 美女免费视频网站| 一级a爱片免费观看的视频| 午夜福利欧美成人| 午夜a级毛片| 久久婷婷人人爽人人干人人爱| 他把我摸到了高潮在线观看| 久久中文看片网| 精品乱码久久久久久99久播| 97超视频在线观看视频| 国产欧美日韩精品亚洲av| 少妇的逼好多水| 90打野战视频偷拍视频| 听说在线观看完整版免费高清| 天美传媒精品一区二区| 免费搜索国产男女视频| 亚洲精品影视一区二区三区av| 成人无遮挡网站| 国产精华一区二区三区| 一级作爱视频免费观看| 热99re8久久精品国产| 国产亚洲欧美98| 校园春色视频在线观看| 日韩欧美国产在线观看| 国产精品一区二区免费欧美| 日韩大尺度精品在线看网址| 中文字幕久久专区| 99精品在免费线老司机午夜| 亚洲狠狠婷婷综合久久图片| 久久精品国产清高在天天线| 精品久久久久久,| 久久99热这里只有精品18| 亚洲自偷自拍三级| 国产精品爽爽va在线观看网站| 国产黄色小视频在线观看| 欧美成人a在线观看| 真人一进一出gif抽搐免费| 国内精品久久久久久久电影| 欧美成人免费av一区二区三区| 中文字幕人妻熟人妻熟丝袜美| 亚洲在线自拍视频| 国产精品不卡视频一区二区 | 日本免费a在线| 毛片一级片免费看久久久久 | 欧美日韩中文字幕国产精品一区二区三区| 亚洲无线在线观看| 狠狠狠狠99中文字幕| 69av精品久久久久久| 最近视频中文字幕2019在线8| 亚洲国产色片| 在线播放国产精品三级| 欧美日韩福利视频一区二区| 国产精品野战在线观看| 欧美+亚洲+日韩+国产| 黄色视频,在线免费观看| 老司机午夜福利在线观看视频| 99久久九九国产精品国产免费| 亚洲欧美日韩无卡精品| 天堂动漫精品| 久久6这里有精品| 精品日产1卡2卡| 亚洲经典国产精华液单 | 精品久久久久久久末码| 热99re8久久精品国产| 最近最新中文字幕大全电影3| 91狼人影院| 久久精品人妻少妇| 免费看美女性在线毛片视频| 国产欧美日韩一区二区精品| 十八禁国产超污无遮挡网站| 国产爱豆传媒在线观看| 女生性感内裤真人,穿戴方法视频| 欧美黄色淫秽网站| 成人国产综合亚洲| 久久精品国产清高在天天线| 啪啪无遮挡十八禁网站| 精品人妻一区二区三区麻豆 | 国产精华一区二区三区| 又黄又爽又刺激的免费视频.| 人妻夜夜爽99麻豆av| 国产黄a三级三级三级人| 热99re8久久精品国产| 免费人成在线观看视频色| 动漫黄色视频在线观看| 精品久久国产蜜桃| 桃色一区二区三区在线观看| 国产精品乱码一区二三区的特点| 一级a爱片免费观看的视频| 一个人免费在线观看电影| 国内精品久久久久久久电影| 久久久国产成人免费| 精品人妻偷拍中文字幕| av在线蜜桃| 人妻制服诱惑在线中文字幕| 欧美日韩国产亚洲二区| 91麻豆av在线| 亚洲欧美激情综合另类| 国产亚洲欧美98| 老司机午夜福利在线观看视频| 一级av片app| 欧美另类亚洲清纯唯美| 国产在线男女| 国产黄片美女视频| 中文字幕av成人在线电影| 少妇丰满av| 一二三四社区在线视频社区8| 日韩欧美免费精品| 99在线视频只有这里精品首页| 久久草成人影院| 中文字幕久久专区| 国产白丝娇喘喷水9色精品| 大型黄色视频在线免费观看| 久久亚洲精品不卡| 少妇的逼好多水| 性插视频无遮挡在线免费观看| 午夜激情欧美在线| 日韩欧美在线二视频| 亚洲午夜理论影院| 99视频精品全部免费 在线| 久久欧美精品欧美久久欧美| 亚洲激情在线av| 免费大片18禁| av女优亚洲男人天堂| 亚洲男人的天堂狠狠| 少妇的逼水好多| 在线a可以看的网站| 亚洲av成人不卡在线观看播放网| 色吧在线观看| 日日干狠狠操夜夜爽| 日韩 亚洲 欧美在线| 中文在线观看免费www的网站| 欧美+日韩+精品| 亚洲中文字幕一区二区三区有码在线看| 亚洲成人中文字幕在线播放| 亚洲国产欧洲综合997久久,| 在线观看美女被高潮喷水网站 | 日本撒尿小便嘘嘘汇集6| 精品不卡国产一区二区三区| 欧美午夜高清在线| 韩国av一区二区三区四区| 一区二区三区高清视频在线| www.999成人在线观看| 欧美黑人巨大hd| 亚洲人成网站在线播| 看片在线看免费视频| 亚洲精品一卡2卡三卡4卡5卡| 日本黄色视频三级网站网址| 成人国产一区最新在线观看| 亚洲成人中文字幕在线播放| 色吧在线观看| 亚洲av日韩精品久久久久久密| 91在线精品国自产拍蜜月| 男人的好看免费观看在线视频| 日韩 亚洲 欧美在线| 成熟少妇高潮喷水视频| 真实男女啪啪啪动态图| 天堂av国产一区二区熟女人妻| 特级一级黄色大片| 精品国产亚洲在线| 色噜噜av男人的天堂激情| 国产精品自产拍在线观看55亚洲| 亚洲最大成人手机在线| 亚洲男人的天堂狠狠| 欧美乱色亚洲激情| 人妻夜夜爽99麻豆av| 国产精品,欧美在线| 色播亚洲综合网| 久久精品国产99精品国产亚洲性色| 丁香六月欧美| 亚洲av二区三区四区| 美女免费视频网站| 欧美性猛交╳xxx乱大交人| 久久午夜亚洲精品久久| 亚洲不卡免费看| 又黄又爽又免费观看的视频| 全区人妻精品视频| 999久久久精品免费观看国产| 身体一侧抽搐| 两个人的视频大全免费| av在线蜜桃| 神马国产精品三级电影在线观看| 中文字幕av成人在线电影| 亚洲精品在线观看二区| 亚洲真实伦在线观看| 91av网一区二区| 国产伦精品一区二区三区视频9| 久久久久久久久大av| 国产成年人精品一区二区| 国产精品自产拍在线观看55亚洲| 国产中年淑女户外野战色| 色综合婷婷激情| 亚洲 国产 在线| 麻豆国产97在线/欧美| 欧美色视频一区免费| 精品久久久久久久久久久久久| 亚洲专区国产一区二区| eeuss影院久久| 久久久久久久精品吃奶| 亚洲久久久久久中文字幕| 久久久久国产精品人妻aⅴ院| 亚洲成人中文字幕在线播放| 69人妻影院| 一卡2卡三卡四卡精品乱码亚洲| 欧美绝顶高潮抽搐喷水| 久久国产乱子免费精品| 美女xxoo啪啪120秒动态图 | 悠悠久久av| 免费高清视频大片| 男人和女人高潮做爰伦理| 日本一本二区三区精品| 身体一侧抽搐| 国产成人福利小说| 99热6这里只有精品| 91九色精品人成在线观看| 午夜影院日韩av| 亚洲最大成人av| 美女 人体艺术 gogo| 日韩欧美在线乱码| 国产野战对白在线观看| 欧美性猛交╳xxx乱大交人| 色视频www国产| 内射极品少妇av片p| 五月伊人婷婷丁香| 啪啪无遮挡十八禁网站| 亚洲第一电影网av| 久久精品91蜜桃| 国产高清激情床上av| 一边摸一边抽搐一进一小说| 黄色丝袜av网址大全| 日本黄色视频三级网站网址| 十八禁网站免费在线| 亚洲人成电影免费在线| 看黄色毛片网站| 精品久久久久久,| 亚洲精品影视一区二区三区av| 精品人妻一区二区三区麻豆 | 久久久久久大精品| 色视频www国产| 欧美日韩瑟瑟在线播放| 国产伦精品一区二区三区视频9| 性色avwww在线观看| 啦啦啦观看免费观看视频高清| 日韩欧美 国产精品| 亚洲最大成人av| 中出人妻视频一区二区| 亚洲国产精品合色在线| 精品一区二区免费观看| 久久九九热精品免费| 国产伦精品一区二区三区四那| 97碰自拍视频| 亚洲av第一区精品v没综合| 日本三级黄在线观看| 精品国产亚洲在线| 久久久久久久午夜电影| 99久久精品国产亚洲精品| 久久精品国产亚洲av涩爱 | 国产国拍精品亚洲av在线观看| 三级毛片av免费| 一个人看视频在线观看www免费| 九九久久精品国产亚洲av麻豆| 日韩成人在线观看一区二区三区| 欧美性猛交╳xxx乱大交人| 男人狂女人下面高潮的视频| h日本视频在线播放| 亚洲自偷自拍三级| 18禁黄网站禁片午夜丰满| 日日摸夜夜添夜夜添av毛片 | 色在线成人网| 91在线精品国自产拍蜜月| 午夜久久久久精精品| 我的老师免费观看完整版| 极品教师在线免费播放| 亚洲欧美精品综合久久99| 午夜两性在线视频| 国产精品日韩av在线免费观看| 桃色一区二区三区在线观看| 亚洲av二区三区四区| 国产av不卡久久| 麻豆av噜噜一区二区三区| 免费无遮挡裸体视频| 免费高清视频大片| 亚洲一区高清亚洲精品| av天堂中文字幕网| 国产av不卡久久| 久久草成人影院| 亚洲中文字幕日韩| or卡值多少钱| 日本在线视频免费播放| 国产免费一级a男人的天堂| 久久久精品欧美日韩精品| 日韩欧美在线二视频| 成人av一区二区三区在线看| 男女之事视频高清在线观看| 国产精品永久免费网站| 一区二区三区高清视频在线| 观看美女的网站| 日本一本二区三区精品| 自拍偷自拍亚洲精品老妇| 69av精品久久久久久| 精品国内亚洲2022精品成人| 亚洲成人中文字幕在线播放| 亚洲最大成人手机在线| 欧美日韩国产亚洲二区| 国产91精品成人一区二区三区| 亚洲人成伊人成综合网2020| 国产毛片a区久久久久| 久久精品国产清高在天天线| 亚洲人成伊人成综合网2020| 国产色婷婷99| 欧美日韩黄片免| 欧美日韩中文字幕国产精品一区二区三区| 午夜福利成人在线免费观看| 1000部很黄的大片| 岛国在线免费视频观看| 十八禁国产超污无遮挡网站| 少妇裸体淫交视频免费看高清| 欧美极品一区二区三区四区| 久久久久久久久久成人| 国产一区二区在线av高清观看| 老司机午夜十八禁免费视频| 两性午夜刺激爽爽歪歪视频在线观看| 91av网一区二区| 啪啪无遮挡十八禁网站| 国产高清三级在线| 男人舔奶头视频| 亚洲av熟女| 婷婷丁香在线五月| 午夜福利高清视频| 中国美女看黄片| 全区人妻精品视频| 性色avwww在线观看| 日韩国内少妇激情av| 国产欧美日韩精品亚洲av| 国产精品伦人一区二区| 夜夜看夜夜爽夜夜摸| 制服丝袜大香蕉在线| 99国产精品一区二区三区| 黄片小视频在线播放| 高潮久久久久久久久久久不卡| 国产精品综合久久久久久久免费| 亚洲成人久久爱视频| 日本成人三级电影网站| 18禁在线播放成人免费| 亚洲av成人精品一区久久| www日本黄色视频网| 国产精品美女特级片免费视频播放器| 精品一区二区三区人妻视频| 好男人在线观看高清免费视频| 亚洲色图av天堂| 琪琪午夜伦伦电影理论片6080| 特大巨黑吊av在线直播| 人人妻人人澡欧美一区二区| 最新在线观看一区二区三区| 日本与韩国留学比较| 一级黄片播放器| 精品人妻1区二区| 老司机午夜十八禁免费视频| 亚洲av成人不卡在线观看播放网| 可以在线观看的亚洲视频| 国产毛片a区久久久久| 精品久久国产蜜桃| 日本三级黄在线观看| 欧美又色又爽又黄视频| 99riav亚洲国产免费| 亚洲一区二区三区不卡视频| 88av欧美| 人妻久久中文字幕网| 亚洲专区中文字幕在线| 国产精品久久视频播放| 在线播放无遮挡| 国产成年人精品一区二区| 中出人妻视频一区二区| 老熟妇仑乱视频hdxx| 中亚洲国语对白在线视频| 中文字幕av在线有码专区| 国产91精品成人一区二区三区| 最新在线观看一区二区三区| 亚洲中文字幕一区二区三区有码在线看| 日日干狠狠操夜夜爽| 一二三四社区在线视频社区8| 性欧美人与动物交配| 中文亚洲av片在线观看爽| 性色av乱码一区二区三区2| 夜夜躁狠狠躁天天躁| 国产av不卡久久| a级一级毛片免费在线观看| 国产精品自产拍在线观看55亚洲| 国产大屁股一区二区在线视频| 亚洲国产精品999在线| 国产淫片久久久久久久久 | 伦理电影大哥的女人| 亚洲无线在线观看| 久久人妻av系列| 在线观看一区二区三区| 夜夜夜夜夜久久久久| 少妇丰满av| 十八禁国产超污无遮挡网站| 亚洲精品色激情综合| 日韩欧美在线二视频| 久久国产乱子伦精品免费另类| 成人av在线播放网站| 亚洲经典国产精华液单 | 窝窝影院91人妻| 国产精品一及| 精品午夜福利视频在线观看一区| 男人舔奶头视频| 1024手机看黄色片| 国产老妇女一区| 久久久久久久精品吃奶| 欧美日本视频| 精品久久久久久久久av| 中文字幕人妻熟人妻熟丝袜美| 国产真实伦视频高清在线观看 | 最近在线观看免费完整版| 麻豆成人av在线观看| 丝袜美腿在线中文| 国产成+人综合+亚洲专区| 亚洲精品粉嫩美女一区| 乱码一卡2卡4卡精品| 久久这里只有精品中国| 精品日产1卡2卡| 亚洲精品日韩av片在线观看| xxxwww97欧美| 97碰自拍视频| 日韩有码中文字幕| 两个人的视频大全免费| 久久久久国内视频| 又黄又爽又免费观看的视频| 制服丝袜大香蕉在线| 中文资源天堂在线| 90打野战视频偷拍视频| 可以在线观看毛片的网站| 美女黄网站色视频| av黄色大香蕉| 久久久成人免费电影| 亚洲第一欧美日韩一区二区三区| 2021天堂中文幕一二区在线观| 美女 人体艺术 gogo| 国产精品久久久久久久电影| 69人妻影院| 91在线观看av| 久久精品影院6| 别揉我奶头 嗯啊视频| 国产精品影院久久| 久久久久九九精品影院| 中文字幕久久专区| 国产乱人视频| 久久久久九九精品影院| 亚洲,欧美,日韩| 国产精品三级大全| 欧美中文日本在线观看视频| 特级一级黄色大片| 国产精品日韩av在线免费观看| 无人区码免费观看不卡| 一本综合久久免费| 久久久久久国产a免费观看| 99久久成人亚洲精品观看| 日韩欧美精品v在线| 免费搜索国产男女视频| 国产成+人综合+亚洲专区| 床上黄色一级片| 丁香六月欧美| 久久久久久久久大av| 91麻豆av在线| 欧美又色又爽又黄视频| 最近在线观看免费完整版| 一级黄色大片毛片| 免费人成视频x8x8入口观看| 特级一级黄色大片| 琪琪午夜伦伦电影理论片6080| 国产欧美日韩精品一区二区| 亚洲av.av天堂| 日韩欧美在线二视频| 丁香欧美五月| 日本黄色片子视频| 日韩人妻高清精品专区| 亚洲最大成人手机在线| 三级男女做爰猛烈吃奶摸视频| 69人妻影院| 国产亚洲精品综合一区在线观看| 欧美一区二区国产精品久久精品|