999精品在线视频,手机成人午夜在线视频,久久不卡国产精品无码,中日无码在线观看,成人av手机在线观看,日韩精品亚洲一区中文字幕,亚洲av无码人妻,四虎国产在线观看 ?

基于波場外推和閾值截斷的鬼波壓制方法

2018-04-09 05:53:23馬繼濤王艷冬陳小宏王建花
石油地球物理勘探 2018年2期
關鍵詞:方法

馬繼濤 王艷冬 陳小宏 王建花

(①中國石油大學(北京)地球物理與信息工程學院物探系,北京 102249;②中海油研究總院有限責任公司,北京 100028)

1 引言

海上地震數據采集中,拖纜必須固定在一定的深度,這會使采集到的數據中存在虛反射,一般將這種虛反射稱為鬼波。鬼波導致數據頻譜發生陷波現象,極大地影響了數據的頻帶寬度。因此鬼波壓制是海上地震數據拓頻處理的必要環節。

在采集設備方面,技術裝備在過去的幾十年來已經取得了很大進步,如固體電纜降低了海洋環境噪聲,新的水鳥控制技術能夠保證不同深度的電纜始終位于同一個垂直平面內[1];在采集技術方面,開發了各種在信號采集過程中壓制鬼波的方法,如上下源采集技術[2],該技術采用兩個垂直的、沉放深度不同的震源,可采集兩個震源的數據,將其進行合并,可降低震源鬼波的影響。雙檢采集技術在同一個深度放置兩個不同的檢波器,速度檢波器測量速度分量,壓力檢波器測量壓力分量,速度分量可彌補壓力分量的陷波點[3,4]。拖纜深度可以適當增加,但存在一個極限值[5],若拖纜放置過深,第一個填充的水檢陷波點對應陸檢所要填充的陷波頻率,而在這個頻率陸檢會有過多的噪聲導致雙檢合并效果變差。上下纜采集技術是在不同深度位置放置兩條拖纜,較深拖纜測量低頻信息,較淺拖纜填補深部拖纜的陷波點[6]。結合兩條拖纜采集的數據可壓制鬼波,得到更加豐富的低頻信息。然而這種采集方式要求高低纜位于同一個垂直平面內,對電纜的定位精度要求高,增加了采集難度。變深度纜采集技術采用一套固體電纜以變深度拖曳方式進行數據采集,因檢波器深度不同,鬼波的陷頻點分散化,可以利用不同沉放深度鬼波特征的差異,用專門的反鬼波和成像處理技術對鬼波進行壓制,獲取低頻和高頻信號,最大限度拓展原始地震數據的頻帶寬度[7]。前人針對海洋地震勘探中的鬼波問題也展開了相應的采集方法調研與研究。吳志強[8]綜合闡述了海洋寬頻帶地震勘探技術的進展,介紹了海洋地震勘探中存在的虛反射問題及人們針對這些問題開發出的一系列采集技術;張振波等[9]對雙檢和上下纜采集數據進行了聯合處理和分析,改善了地震資料品質;張振波等[10]介紹了斜纜采集技術在珠江口盆地的應用,對斜纜采集系統、激發參數以及處理技術流程和應用效果進行了闡述和分析;劉春成等[11]提出了一種“犁式”電纜地震采集方法,對采集參數進行了分析,并應用于實際采集,拓寬了地震頻帶。

在海洋地震數據鬼波壓制處理方面也有針對性的研究。Fred等[12]早在1989年提出了雙檢合并鬼波壓制技術,該方法認為,對于雙檢采集得到的數據,由于兩種檢波器的自身結構和工作原理不同,二者對上行波響應極性相同,而對于鳴震等下行波的響應極性相反,因此二者對于鬼波具有極性相反的記錄,將兩種數據進行合并會有效壓制鬼波、增強有效信號。然而該方法由于陸檢資料信噪比低,頻率存在差異,以及海底反射系數變化等因素,簡單疊加兩種數據無法得到滿意的鬼波壓制效果,必須求取準確的校準算子對陸檢數據進行校正。對上下纜采集的海洋地震數據,Nick等[6]于2007年提出針對性的聯合鬼波壓制方法,該方法將一次波認定為上行波、鬼波認定為下行波,利用波場外推對接收到的地震波場進行波場分離,得到上行波場,并用其替代原來的總波場,以達到消除檢波器鬼波的目的;針對斜纜采集得到的地震數據,Robert[13]提出利用偏移和鏡像偏移數據進行聯合反褶積的鬼波壓制方法,該方法首先生成相對海平面對稱的鏡像數據,再對數據和鏡像數據分別進行疊前時間偏移,最后利用偏移后的兩個數據道集進行聯合反褶積,達到鬼波壓制的目的;許自強等[14,15]、丁洪波[16]、王沖等[17]也針對聯合反褶積算法展開研究,對該算法進行了優化完善,改進了該方法的應用效果; Wang等[18,19]還提出分別在F-XY域和τ-p域自動壓制斜纜采集地震數據的疊前鬼波的方法,該方法首先生成鏡像數據,隨后在疊前道集中利用數據及其鏡像數據,通過迭代反演壓制鬼波,并不需要精確的纜深信息,可用于多種形式的采集數據,大大提高了鬼波壓制精度; 宋建國等[20]基于頻率域高分辨率Radon變換擬合鬼波,進而將其自適應減去,取得了良好的鬼波壓制效果; Lu等[21,22]認為地震信號呈超高斯分布,在時空域可以利用地震信號的超高斯性壓制鬼波,用海水表面的反射系數和鬼波時差模擬鬼波,在時空域用二維掃描估計反射系數和時差,并設計濾波器濾除鬼波; 王芳芳等[23]研究了基于逆散射理論的鬼波壓制技術,該方法不依賴于速度模型,適用于多維復雜介質,能在頻率—波數域較好地壓制鬼波,恢復地震有效信息,有效保護一次波,但該方法在處理低信噪比地震數據時對震源和檢波點深度信息有一定依賴性,且計算量較大; 赫建偉等[24]研究了基于波場延拓最小二乘法去除鬼波的技術,該方法基于上下纜采集數據波場與上行波和鬼波算子之間的關系,利用最小二乘法去除鬼波并將剩余噪聲最小化,提取上下纜數據中的上行波; 張威等[25]基于LSMR算法壓制斜纜數據鬼波,該方法基于平面波傳播理論,建立了τ-p域拖纜觀測的總波場與海水表面上行波波場之間的關系,利用最小平方殘差算法精確求解線性τ-p方程,得到鬼波壓制后的波場,在提高計算效率的同時改善了鬼波壓制效果; Berkhout等[26]研究了基于混合采集地震數據解混算法的鬼波壓制技術,該方法將含有鬼波(如震源鬼波)的地震數據視為由混合地震震源形成的地震數據,一個震源位于真實震源位置,另一個震源位于鏡像震源位置,因此可以將鬼波的壓制問題視為多源數據的分離問題,進而應用解混算法對鬼波進行壓制。

本文針對海洋地震數據的鬼波壓制問題,基于Berkhout等[26]的理論研究,給出了一種基于波場外推、閾值截斷的鬼波壓制算法。該算法首先將地震數據進行兩次波場外推,然后對外推結果求和并做閾值截斷處理,生成一次波數據;利用生成的一次波數據進行二次波場外推,得到鬼波記錄并將其從原數據中自適應減去。如此迭代處理,在迭代過程中逐步改進鬼波壓制效果。該算法是完全數據驅動,可以壓制震源一側的鬼波,也可以壓制檢波點一側的鬼波。本文利用模擬數據和實際斜纜數據驗證了算法的有效性。

2 方法原理

帶有鬼波的地震數據,可以視為真實震源產生的波場與虛擬震源產生的波場之和,假定真實震源位置為+zs,虛擬震源位置為-zs,檢波器所在位置為z0,則利用Berkhout[27]數據矩陣的概念,可將帶鬼波的地震數據表示為

P-(z0;±zs)=P-(z0;+zs)+P-(z0;-zs)

(1)

而+zs處真實震源和-zs處虛擬震源產生的波場,可以視為由位于水表面z0處震源產生的波場外推得到

P-(z0;+zs)=P-(z0;z0)F-(z0;+zs)

(2)

P-(z0;-zs)=P-(z0;z0)R∩(z0;z0)×

W+(z0;-zs)

(3)

式中:F-和W+分別為在水層中的反向波場外推算子和正向波場外推算子;R∩為海水表面的反射系數。若震源位置、水速等參數已知,則可以計算出F-、W+和R∩。

依據本文方法進行鬼波壓制的流程如圖1所示。

利用本算法進行鬼波壓制處理需遵循以下步驟。

(1)對原含有鬼波的數據進行兩次波場外推,一次外推是將位于+zs處震源產生的波場(一次波波場)外推至水表面,第二次外推是將位于-zs處的虛擬震源產生的波場(鬼波波場)外推至水表面。依據式(2)和式(3),這兩次外推(第一次迭代)可以表示為

(4)

=-[P-(z0;-zs)](0)F-(-zs;z0)R-1

(5)

式(4)將原本位于+zs處的震源外推至水表面z0處,此時數據中位于-zs處的虛擬震源也受到算法的影響,被外推至-2zs處;而式(5)將虛擬震源由-zs處外推至海水表面z0處,真實震源受到算法的影響,被外推至+2zs處。

(2)將外推結果求和并做閾值截斷處理。將式(4)與式(5)相加并除以2,這時震源位于海水表面z0處產生的波場,其振幅強度應為位于-2zs和+2zs處震源產生波場振幅強度的2倍。若對此結果進行閾值截斷處理,可以得到位于海水表面處的一次波波場。相加結果可以表示為

[P-(z0;-zs)](k)F-(-zs;z0)R-1}

(6)

式中上標k代表循環迭代次數(k=0,1,2,…)。將式(2)和式(3)的波場外推表達式代入上式,可以得到第一次迭代產生的結果

(7)

(3)將閾值截斷得到的位于海水表面處震源產生的波場外推,即將位于z0處的震源分別外推至+zs和-zs處,得到新的一次波波場和鬼波波場,即

(8)

R∩(z0;z0)W+(z0;-zs)

(9)

(4)將波場延拓后得到的一次波和鬼波分別從輸入數據中自適應減去,分別得到鬼波記錄和一次波記錄

(10)

(11)

式中A(k)(-zs)、A(k)(+zs)為自適應相減的算子,對波場外推中由于實際數據水速、水表面反射系數等因素引起的外推波場與實際數據之間的差異進行調整,達到鬼波壓制的目的。

(5)最后將得到的一次波波場、鬼波波場從原輸入數據中減去,得到波場的殘差

ΔP-(z0;±zs)=P-(z0;±zs)-

(12)

圖1 基于閾值截斷和波場外推的迭代鬼波壓制方法流程

對于檢波器一側鬼波的壓制,可以遵循類似的迭代步驟。如果同一個炮記錄中,檢波器深度是變化的,也可以利用本文的算法壓制檢波器鬼波,變深度纜鬼波壓制與常規拖纜鬼波壓制唯一的區別在于波場外推算法不同。

3 數據實例

下面利用一個模擬數據和實際斜纜數據驗證本文方法的正確性與有效性。

3.1 模擬數據

模擬數據來自兩個水平界面組成的模型,檢波點深度為5m,炮點深度為50m,觀測系統采用單邊接收。利用有限差分正演算法模擬單炮記錄,經過直達波切除后,如圖2所示。

圖2 含有鬼波的炮記錄

首先基于式(4)和式(5),利用水速和波場外推算法對帶有鬼波記錄的炮記錄做波場外推,分別將炮記錄從真實震源位置(50m)外推至水表面處(0m),從虛擬震源位置(-50m)外推至水表面處(0m);然后利用式(6)將得到的兩個波場外推記錄相加并除以2。需要注意的是,為保證方法效果,在進行波場外推時,要選取一個相對較小的孔徑進行外推,孔徑隨著迭代次數的增加而增大,在最后一次迭代時,外推孔徑達到最大,即所有的炮記錄都參與到波場外推中。波場外推相加后的結果如圖3左所示;對該結果進行閾值截斷處理,得到的結果如圖3右所示。

閾值截斷的結果就是水表面(0m)處的記錄。利用式(8)和式(9)對該記錄做進一步的波場外推,分別得到真實震源位置產生的一次波記錄和虛擬震源位置產生的鬼波記錄(圖4);利用式(10)和式(11)將此次波場外推的結果分別從輸入數據中自適應減掉,得到第一次迭代的鬼波和一次波(圖5)。

第一次迭代的結果明顯不滿足算法的終止條件,因此利用圖5的一次波和鬼波進行第二次迭代處理,同樣經過波場外推、求和、閾值截斷、波場外推,得到第二次迭代的鬼波和一次波。如此反復,第30次迭代得到的求和后的波場及閾值截斷后的結果如圖6所示,對閾值截斷結果進行波場外推后的結果如圖7所示,將波場外推結果從輸入數據中自適應減掉的結果如圖8所示。可以看到,經30次迭代已經達到了較好的鬼波壓制效果,多次波對應的虛反射也得到了較好的壓制。

圖3 原始數據波場外推求和結果(左)及閾值截斷后的結果(右)

圖4 閾值截斷結果波場外推產生的一次波(左)和鬼波(右)記錄

圖5 第1次迭代自適應相減得到的一次波(左)和鬼波(右)記錄

3.2 實際數據

為進一步驗證方法的應用效果,利用某區實際采集得到的斜纜地震數據對方法的有效性進行了驗證。斜纜數據采集中,每個檢波點的深度都是不同的,在波場延拓時,需要將每個檢波點深度分別做延拓,將其延拓至水表面處,同時也需要對檢波點的鏡像位置分別做延拓,將其延拓至水表面處。由于該數據的波場延拓算法需運行多次(次數與數據的道數有關),因此斜纜數據鬼波壓制算法計算速度較慢。圖9左為該斜纜數據經過相關去噪處理后的單炮記錄,圖9右為本文方法30次迭代后的處理結果,圖10為鬼波壓制前、后的頻譜。從炮記錄中,尤其是箭頭所指的位置可以明顯看出,本文方法很好地壓制了檢波點鬼波。從頻譜圖(圖10)中也可以看出,本文方法較好地補充了原數據中的陷波帶。

圖11和圖12展示了鬼波壓制前、后該斜纜數據的兩個共炮檢距剖面,圖13為鬼波壓制前、后共炮檢距剖面的自相關譜。從共炮檢距剖面(圖11和圖12)箭頭所指位置可以看出,一次波下方緊跟的鬼波被較好地壓制;鬼波壓制后自相關譜的旁瓣消失,也說明了本文方法對鬼波壓制的效果較好。

圖6 第30次迭代波場外推求和結果(左)和進行閾值截斷后的結果(右)

圖8 第30次迭代自適應相減得到的一次波(左)和鬼波(右)記錄

圖9 斜纜數據炮記錄(左)及本文算法進行鬼波壓制后的結果(右)

圖10 斜纜數據炮記錄頻譜(左)及本文算法進行鬼波壓制結果的頻譜(右)

圖11 鬼波壓制前(左)、后(右)的共炮檢距剖面(炮檢距為1350m)

圖12 鬼波壓制前(左)、后(右)的共炮檢距剖面(炮檢距為3105m)

圖13 鬼波壓制前(左)、后(右)的共炮檢距剖面的自相關譜

4 結論

本文給出了一種基于波場外推和閾值截斷的鬼波壓制迭代方法,方法基于波場外推理論,通過對外推求和后的波場數據進行閾值截斷處理,得到一次波模型數據,然后再次進行波場外推得到鬼波;之后將鬼波從原始數據中自適應減去,得到鬼波壓制結果。應用模擬數據和實際變深度纜數據驗證了方法的有效性,并得到以下認識:

(1)為得到最優的鬼波壓制效果,本文方法需進行迭代處理,且隨著迭代次數增加,所選用截斷的閾值逐漸減小為零或接近于零,波場外推孔徑逐漸增大;

(2)本文方法可以進行炮點鬼波的處理,也可以進行檢波點鬼波的處理;

(3)本文方法可以進行檢波器深度變化的數據檢波點鬼波壓制,該數據與常規數據鬼波壓制的唯一不同之處在于波場外推算法的不同;

(4)算法的計算量較大,實際數據處理時應考慮并行運算。

[1]吳志強,閆桂京,童思友等.海洋地震采集技術新進展及對我國海洋油氣地震勘探的啟示.地球物理學進展,2013,28(6):3056-3065.

Wu Zhiqiang,Yan Guijing,Tong Siyou et al.New advances in marine seismic data acquisition technology and enlightenment to domestic offshore petroleum seismic exploration.Progress in Geophysics,2013,28(6):3056-3065.

[2]Moldoveanu N.Vertical source array in marine seismic exploration.SEG Technical Program Expanded Abstracts,2000,19:53-56.

[3]Posthumus B J.Deghosting using a twin streamer configuration.Geophysical Prospecting,1993,41(3):267-286.

[4]Carlson D H,Long A,S?llner W et al.Increased resolution and penetration from a towed dual-sensor streamer.First Break,2007,25(12):71-77.

[5]Tenghamn R,Vaage S and Borresen C.A dual-sensor towed marine streamer:Its viable implementation and initial results.SEG Technical Program Expanded Abstracts,2007,26:989-993.

[6]Moldoveanu N,Combee L,Egan M et al.Over/under towed-streamer acquisition:A method to extend seismic bandwidth to both higher and lower frequencies.The Leading Edge,2007,26(1):41-58.

[7]Soubaras R and Lafet Y.Variable-depth streamer acquisition:broadband data for imaging and inversion.Geophysics,2013,78(2):WA27-WA39.

[8]吳志強.海洋寬頻帶地震勘探技術新進展.石油地球物理勘探,2014,49(3):421-430.

Wu Zhiqiang.New advances in marine broadband seismic exploration.OGP,2014,49(3):421-430.

[9]張振波,軒義華,劉賓等.雙檢與上下纜地震數據聯合成像.石油地球物理勘探,2014,49(5):884-891.

Zhang Zhenbo,Xuan Yihua,Liu Bin et al.Joint imaging of dual-sensor and over-under streamer data.OGP,2014,49(5):884-891.

[10]張振波,李東方.斜纜寬頻地震勘探技術在珠江口盆地的應用.石油地球物理勘探,2014,49(3):451-456.

Zhang Zhenbo,Li Dongfang.Variable depth streamer seismic acquisition and processing in Pearl River Mouth Basin.OGP,2014,49(3):451-456.

[11]劉春成,劉志斌,陶杰等.海洋“犁式”電纜采集技術研究與實踐.石油地球物理勘探,2016,51(6):1069-1074.

Liu Chuncheng,Liu Zhibin,Tao Jie et al.Marine plow streamer acquisition.OGP,2016,51(6):1069-1074.

[12]Fred J B and Joe I S.Attenuation of water-column reverberations using pressure and velocity detectors in a water-bottom cable.SEG Technical Program Expanded Abstracts,1989,8:653-656.

[13]Robert S.Deghosting by joint deconvolution of a migration and a mirror migration.SEG Technical Program Expanded Abstracts,2010,29:3406-3410.

[14]許自強,方中于,顧漢明等.海上變深度纜數據最優化壓制鬼波方法及其應用.石油物探,2015,54(4):404-413.

Xu Ziqiang,Fang Zhongyu,Gu Hanming et al.The application of optimal deghosting algorithm on marine variable-depth streamer data.GPP,2015,54(4):404-413.

[15]許自強,李添才,王用軍等.傾斜電纜地震資料處理關鍵技術及其效果分析.中國海上油氣,2015,27(6):10-18.

Xu Ziqiang,Li Tiancai,Wang Yongjun et al.Key technology for slant streamer seismic data processing and its application analysis.China Offshore Oil and Gas,2015,27(6):10-18.

[16]丁洪波.海上地震資料鬼波機理與消除方法研究[學位論文].山東青島:中國石油大學(華東),2013.

Ding Hongbo.Research on the Mechanism of Ghost and Deghosting Method for Marine Seismic Data[D].China University of Petroleum(East China),Qing-dao,Shandong,2013.

[17]王沖,顧漢明,許自強等.頻率慢度域自適應迭代反演算法壓制海上傾斜纜鬼波方法及其應用.地球物理學報,2016,59(12):4677-4689.

Wang Chong,Gu Hanming,Xu Ziqiang et al.The application of adaptive iteration inversion algorithm to deghost for marine variable-depth streamer data in frequency-slowness domain.Chinese Journal of Geophysics,2016,59(12):4677-4689.

[18]Wang P and Peng C.Premigration deghosting for marine towed streamer data using a bootstrap approach.SEG Technical Program Expanded Abstracts,2012,31:1-5.

[19]Wang P,Ray S,Peng C et al.Premigration deghosting for marine streamer data using a bootstrap approach in Tau-p domain.SEG Technical Program Expanded Abstracts,2013,32:4221-4225.

[20]宋建國,宮云良,李珊.斜纜數據頻率域高分辨率Radon變換與虛反射消除方法.應用地球物理,2015,12(4):564-572.

Song Jianguo,Gong Yunliang,Li Shan.High resolution frequency domain Radon transform and variable-depth streamer data deghosting.Applied Geophysics,2015,12(4):564-572.

[21]Lu W K,Xu Z Q,Fang Z Y et al.Deghosting for variable depth streamer based on super-Gaussianity.SEG Technical Program Expanded Abstracts,2015,34:4545-4549.

[22]Lu W K,Xu Z Q,Fang Z Q et al.Receiver deghosting in the t-x domain based on super-Gaussianity.Journal of Applied Geophysics,2017,36:134-141.

[23]王芳芳,李景葉,陳小宏.基于逆散射級數法的鬼波壓制方法.地球物理學報,2013,56(5):1628-1636.

Wang Fangfang,Li Jingye,Chen Xiaohong.Deghos-ting method based on inverse scattering series.Chinese Journal of Geophysics,2013,56(5):1628-1636.

[24]赫建偉,顧漢明,吳耀樂等.波場延拓最小二乘法去虛反射技術在上下纜采集數據處理中的應用.石油地球物理勘探,2015,50(3):424-430.

He Jianwei,Gu Hanming,Wu Yaole et al.Ghost elimi-nation with extrapolation least square for seismic data acquired by over/under towed-streamer.OGP,2015,50(3):424-430.

[25]張威,韓立國,李洪建等.基于LSMR算法的斜纜數據鬼波壓制方法.石油地球物理勘探,2017,52(3):434-441.

Zhang Wei,Han Liguo,Li Hongjian et al.Deghosting of variable depth streamer data based on LSMR.OGP,2017,52(3):434-441.

[26]Berkhout A J and Blacquière G.Deghosting by echo-deblending.Geophysical Prospecting,2016,64(2):406-420.

[27]Berkhout A J.Seismic Migration:Imaging of Acoustic Energy by Wave Field Extrapolation.Elsevier,1982.

猜你喜歡
方法
中醫特有的急救方法
中老年保健(2021年9期)2021-08-24 03:52:04
高中數學教學改革的方法
河北畫報(2021年2期)2021-05-25 02:07:46
化學反應多變幻 “虛擬”方法幫大忙
變快的方法
兒童繪本(2020年5期)2020-04-07 17:46:30
學習方法
用對方法才能瘦
Coco薇(2016年2期)2016-03-22 02:42:52
最有效的簡單方法
山東青年(2016年1期)2016-02-28 14:25:23
四大方法 教你不再“坐以待病”!
Coco薇(2015年1期)2015-08-13 02:47:34
賺錢方法
捕魚
主站蜘蛛池模板: 一边摸一边做爽的视频17国产| 日韩中文字幕亚洲无线码| 久久精品免费看一| 免费在线播放毛片| 日韩精品免费在线视频| 毛片久久网站小视频| 国产九九精品视频| 亚洲AV无码一区二区三区牲色| 成人免费一级片| 亚洲国语自产一区第二页| 免费Aⅴ片在线观看蜜芽Tⅴ | 热久久这里是精品6免费观看| 欧美性久久久久| 国产免费怡红院视频| 9999在线视频| 国产视频大全| 久久久噜噜噜| 国产96在线 | 中文字幕波多野不卡一区| 亚洲欧美成人网| 精品国产一区二区三区在线观看 | 欧美一区精品| 国产特级毛片| 另类综合视频| 免费大黄网站在线观看| 国产91在线|中文| 成人毛片免费在线观看| 国产经典在线观看一区| 久久青草精品一区二区三区| 亚洲—日韩aV在线| 国产高清在线观看91精品| 亚洲人视频在线观看| 国产成人凹凸视频在线| 婷婷中文在线| 亚洲欧洲美色一区二区三区| 国产视频自拍一区| 91激情视频| 小13箩利洗澡无码视频免费网站| 国禁国产you女视频网站| 国产欧美日韩综合一区在线播放| 一本大道香蕉高清久久| 亚洲欧美极品| 尤物在线观看乱码| 精品国产香蕉在线播出| 在线毛片网站| 青草91视频免费观看| 日韩在线欧美在线| aaa国产一级毛片| 国产又粗又猛又爽视频| 人妻精品全国免费视频| 亚洲视频三级| 国产欧美日韩另类精彩视频| 久久婷婷人人澡人人爱91| 一级福利视频| 91无码国产视频| www中文字幕在线观看| 亚洲精品国产自在现线最新| 欧美一区二区福利视频| 男女男免费视频网站国产| 狠狠色综合久久狠狠色综合| 亚洲av中文无码乱人伦在线r| 免费观看无遮挡www的小视频| 大陆精大陆国产国语精品1024| 九色视频在线免费观看| 中文字幕欧美日韩| 亚洲区一区| 9啪在线视频| 亚洲区第一页| 亚洲天堂网在线观看视频| 欧美成人手机在线观看网址| 亚洲欧美日韩动漫| 亚洲成a人在线播放www| 青草视频久久| 久久精品人妻中文系列| 在线国产91| 国产 日韩 欧美 第二页| 欧美色综合久久| 亚洲无码精品在线播放| 波多野结衣无码中文字幕在线观看一区二区 | 国模沟沟一区二区三区| 亚洲精品综合一二三区在线| 欧美一级黄片一区2区|