亚洲免费av电影一区二区三区,日韩爱爱视频,51精品视频一区二区三区,91视频爱爱,日韩欧美在线播放视频,中文字幕少妇AV,亚洲电影中文字幕,久久久久亚洲av成人网址,久久综合视频网站,国产在线不卡免费播放

        ?

        GPS網(wǎng)共模誤差的空間響應(yīng)分析

        2016-11-17 02:32:44龔國棟賀小星花向紅丁凌航
        關(guān)鍵詞:區(qū)域

        龔國棟 賀小星,3 花向紅,2 舒 穎 丁凌航

        1 武漢大學(xué)測繪學(xué)院,武漢市珞喻路129號(hào),430079 2 武漢大學(xué)災(zāi)害監(jiān)測與防治研究中心,武漢市珞喻路129號(hào),430079 3 華東交通大學(xué)土木建筑學(xué)院,南昌市雙港東大街808號(hào),330013

        ?

        GPS網(wǎng)共模誤差的空間響應(yīng)分析

        龔國棟1賀小星1,3花向紅1,2舒 穎1丁凌航1

        1 武漢大學(xué)測繪學(xué)院,武漢市珞喻路129號(hào),430079 2 武漢大學(xué)災(zāi)害監(jiān)測與防治研究中心,武漢市珞喻路129號(hào),430079 3 華東交通大學(xué)土木建筑學(xué)院,南昌市雙港東大街808號(hào),330013

        對共模誤差的空間變化進(jìn)行分析。以美國南加州及周邊區(qū)域的GPS臺(tái)站網(wǎng)中3個(gè)站臺(tái)為中心,通過劃分9個(gè)不同的尺度空間,采用基于主成分分析(PCA)的空間濾波方法,對每個(gè)空間尺度GPS坐標(biāo)序列共模誤差進(jìn)行提取,探討不同空間尺度下GPS網(wǎng)共模誤差的空間響應(yīng)特性,為子區(qū)域的劃分提供相應(yīng)的理論基礎(chǔ)。

        GPS坐標(biāo)時(shí)間序列;共模誤差;PCA空間濾波;空間尺度

        共模誤差(common mode error,CME)是GPS時(shí)間序列數(shù)據(jù)誤差的主要來源之一,極大地影響GPS解的精度與可靠性[1-2]。為了準(zhǔn)確地提取CME,國內(nèi)外學(xué)者作了一些研究工作。斯克里普斯軌道與永久陣列中心(Scripps orbit and permanent array center,SOPAC)將整個(gè)GPS觀測網(wǎng)(板塊邊界觀測網(wǎng)plate boundary observation,PBO)劃分為若干個(gè)子區(qū)域,再對各子區(qū)分別進(jìn)行濾波。但該方法存在一定的局限性,對如何進(jìn)行子區(qū)域的劃分,不同的區(qū)域劃分是否存在差異,沒有進(jìn)行深入的分析探討[3]。另外,有學(xué)者指出共模誤差在空間上不是等值的,即共模誤差的空間響應(yīng)并不是完全一致的[4]??紤]到CME的物理起源尚不明了,同時(shí)目前對其時(shí)空響應(yīng)機(jī)制缺乏深入的研究,使得CME的分離存在一定的局限性,尤其是在大尺度范圍的GPS網(wǎng)中,給準(zhǔn)確、可靠地分離出CME造成了一定的困難。本文針對共模誤差分離方法中存在的上述問題,以美國南加州及周邊區(qū)域的一個(gè)GPS網(wǎng)為例,采用主成分分析的空間濾波方法,通過對不同尺度下的GPS網(wǎng)進(jìn)行共模誤差分離處理,分析共模誤差的空間分布特性,探討共模誤差空間響應(yīng)機(jī)制及其最佳分離策略。

        1 GPS時(shí)間序列獲取及其模型

        為了研究不同尺度下共模誤差的時(shí)空變化,本文選取了ITRF2008框架下美國南加州及周邊區(qū)域的一個(gè)GPS網(wǎng)為例,分析不同尺度下共模誤差的空間分布及其變化規(guī)律。以3個(gè)GPS站(BEMT,PIN1,MONP)為中心,通過控制GPS網(wǎng)的尺度,不斷增加GPS網(wǎng)的大小,使其從小區(qū)域網(wǎng)逐漸擴(kuò)展到大尺度的網(wǎng)形,根據(jù)距離劃分多尺度GPS網(wǎng)區(qū)域,最大尺度的GPS網(wǎng)包含39個(gè)站,其站點(diǎn)分布見圖1。通過對臺(tái)站的時(shí)間序列進(jìn)行處理,對不同尺度下的GPS網(wǎng)的共模誤差空間變化進(jìn)行探討。

        圖1 GPS網(wǎng)站點(diǎn)分布Fig.1 Distribution of GPS satellites

        通過對圖1中GPS站2001~2006年的連續(xù)觀測數(shù)據(jù)進(jìn)行計(jì)算處理與分析,GPS坐標(biāo)時(shí)間序列獲取步驟如下[5]:

        1)首先采用GAMIT軟件對站點(diǎn)的原始觀測數(shù)據(jù)進(jìn)行處理,選取的站點(diǎn)數(shù)據(jù)時(shí)間跨度為2001~2006年,數(shù)據(jù)要求觀測時(shí)間大于2.5 a,數(shù)據(jù)缺失程度小(本文數(shù)據(jù)平均缺失的比率為1.46%)。以每天24 h觀測數(shù)據(jù)為基本單位,解算出各站點(diǎn)的三維坐標(biāo)及其相關(guān)參數(shù),形成單天解。

        2) 在統(tǒng)一的基準(zhǔn)下(ITRF2008),結(jié)合SOPAC的單天解序列(GIPSY解算結(jié)果),通過公共站點(diǎn)和衛(wèi)星,采用QOCA[6](quasi-observation combination analysis)軟件對單天解進(jìn)行聯(lián)合解算(GAMIT與GIPSY的權(quán)比為1∶2∶4),獲得臺(tái)站網(wǎng)的坐標(biāo)及其速度解等參數(shù),得到GPS坐標(biāo)時(shí)間序列。在QOCA軟件進(jìn)行聯(lián)合平差過程中,給出三坐標(biāo)分量方向(E、N、U)的殘差限制值分別為30 mm、30 mm、100 mm。一旦觀測值的殘差超過這個(gè)值,就認(rèn)為是過于弱的觀測,予以剔除[6]。解算后的GPS原始坐標(biāo)時(shí)間序列見圖2a(限于篇幅,以PIN1站點(diǎn)為例)。

        已有的研究表明,GPS單站、單分量位置間序列通常滿足模型[7-9]:

        y(ti)=a+bti+csin(2πti)+dcos(2πti)+

        esin(4πti)+fcos(4πti)+

        (1)

        式中,ti(i=1,…,n)為GPS站點(diǎn)單日歷元,以a為單位;a為站點(diǎn)位置(為序列的平均值);b為線速度;c、d和e、f分別為年周期和半年周期項(xiàng)的系數(shù);g為發(fā)生在歷元Tg處的由于各種原因引起的階躍(或稱偏移);H為Heaviside階梯函數(shù),這里假定發(fā)生偏移的時(shí)刻Tg已知;τj為粘滯系數(shù)。根據(jù)上述經(jīng)典模型對GPS時(shí)間序列建立模型,以便對時(shí)間序列進(jìn)行后續(xù)分析。

        2 PCA空間濾波方法

        為了評估共模誤差對GPS坐標(biāo)殘差時(shí)間序列的影響,采用空間濾波改正前后GPS時(shí)間序列的加權(quán)均方根(weighted root mean square,WRMS)進(jìn)行比較[10]。若空間濾波后GPS坐標(biāo)殘差時(shí)間序列WRMS減小,即表明非線性變化幅度減弱。

        2.1 區(qū)域?yàn)V波數(shù)學(xué)模型

        假設(shè)共模誤差在某一區(qū)域均勻分布,將單日解誤差作為權(quán)重因子[9],區(qū)域疊加濾波法的共模誤差分量εi可以表示為:

        (2)

        式中,εi為第i個(gè)臺(tái)站的共模誤差值;v為殘差坐標(biāo)時(shí)間序列;δi,j為單日坐標(biāo)解得的殘差;j為參與計(jì)算共模誤差的臺(tái)站數(shù)。

        2.2 PCA數(shù)學(xué)模型

        本文以GPS臺(tái)站的時(shí)間序列為例,各臺(tái)站的時(shí)間序列排列起來形成一個(gè)n×m(n>m,n為觀測數(shù)或歷元數(shù),m為觀測類型)的數(shù)據(jù)矩陣X,其協(xié)方差陣為CX,則CX=XTX。數(shù)據(jù)矩陣如下:

        (3)

        (4)

        作如下定義:

        (5)

        (6)

        (7)

        (8)

        (9)

        令Λ=ΣTΣ,則有:

        CX=XTX=(UΣVT)TUΣVT=VΛVT

        (10)

        即V構(gòu)成X的正交基底。矩陣X按照KLE展開可得:

        (11)

        ak(ti)可由下式求出:

        X=AV?A=XV-1?A=XVT

        (12)

        (13)

        式中,ak(t)是第k個(gè)主成分,vk(x)是對應(yīng)主成分的響應(yīng)特征向量,分別代表時(shí)間特征和空間響應(yīng)。

        3 GPS時(shí)間序列算例分析

        原始序列中包含著速度項(xiàng)、階躍(offset)等非構(gòu)造信號(hào),往往需要對其進(jìn)行去趨勢、去周期項(xiàng)等處理。通過對單日解觀測時(shí)間序列建立時(shí)間序列模型,對原始時(shí)間序列去除速度項(xiàng)、階躍、指數(shù)和對數(shù)衰減等項(xiàng),最后得到殘差時(shí)間序列如圖2(b)。從圖2可知,臺(tái)站殘差序列波動(dòng)比較大,坐標(biāo)分量中誤差相對較大,且坐標(biāo)序列存在周期性波動(dòng)趨勢,在U方向(垂向分量)尤為明顯。這說明GPS殘差時(shí)間序列中存在與時(shí)間相關(guān)的噪聲;從圖3對應(yīng)的剩余頻譜圖可知,在低頻處的譜能量較大,頻譜呈傾斜(斜率趨近于-1~0之間),說明噪聲中包含著有色噪聲。

        圖2 站臺(tái)PIN1Fig.2 GPS satellite site PIN1

        圖3 PIN1站殘差序列頻譜圖Fig.3 PIN1 station residuals spectrum

        本案例以圖1GPS臺(tái)站網(wǎng)中3個(gè)站臺(tái)(BEMT,PIN1,MONP)為中心,劃分9個(gè)尺度空間(空間尺度分別為100km、200km、420km、550km、660km、810km、1 000km、1 660km、2 000km),GPS時(shí)間序列處理過程中采用相同的處理策略,然后對殘差序列采取主成分分析(PCA)的空間濾波方法對共模誤差進(jìn)行分離。PCA濾波后各GPS臺(tái)站的殘差坐標(biāo)時(shí)間序列(以200km區(qū)域網(wǎng)中PIN1站點(diǎn)為例)如圖2(c)。對GPS殘差時(shí)間序列進(jìn)行PCA主分量分離之前,考慮到環(huán)境負(fù)載會(huì)對GPS坐標(biāo)時(shí)間序列WRMS產(chǎn)生影響,本文通過采用mload程序?qū)PS站作負(fù)載改正,對大氣、積雪和土壤水、海洋非潮汐4項(xiàng)負(fù)荷效應(yīng)進(jìn)行負(fù)載改正,最后分別對不同尺度下的GPS網(wǎng)進(jìn)行處理,濾波處理過程中選取PCA的前4個(gè)分量進(jìn)行空間濾波。

        圖4表明,空間濾波前BEMT、PIN1、MONP3個(gè)站點(diǎn)E方向的WRMS分別為1.665mm、2.468mm、1.966mm;N方向的WRMS分別為1.233mm、1.532mm、1.381mm;U方向的WRMS分別為4.213mm、5.795mm、4.544mm。濾波后的WRMS-距離曲線顯示,WRMS值總體上隨距離的增加而增加,表明隨著距離增加臺(tái)站之間相互關(guān)系減弱,大尺度空間的共模誤差更加難以提取處理。在距離小于420km時(shí)共模誤差的濾波效果較明顯,在這個(gè)區(qū)間中隨著距離的增加,WRMS的增加明顯,曲線斜率大,到達(dá)420km~600km區(qū)間的時(shí)候WRMS的增加變緩慢。距離在大于600km區(qū)間時(shí)E方向和N方向的WRMS幾乎不再變化,1 000km距離以后U方向的WRMS也不再變化??臻g濾波后在100km尺度空間上提取效果最好,BEMT、PIN1、MONP3個(gè)站點(diǎn)E、N、U方向的WRMS減少最明顯,E方向WRMS減少分別為39.7%、93.2%、71.8%;N方向WRMS減少分別為40.2%、97%、61.3%;E方向WRMS減少分別為56.3%、86.4%、37.4%。圖4中各條曲線超過420km后WRMS的曲線斜率基本不再發(fā)生變化,且濾波效果較差。根據(jù)這一曲線的變化趨勢,筆者選取200km、420km、600km、1 000km這4個(gè)區(qū)域的所有測站進(jìn)行研究。

        圖4 BEMT、PIN1、MONP濾波前后WRMS-距離曲線Fig.4 BEMT, PIN1, MONP WRMS- distance curve before and after the filter

        表1說明,在200km空間尺度下區(qū)域?yàn)V波后,各站點(diǎn)在E、N、U 3個(gè)坐標(biāo)分量上WRMS平均減少百分比分別為49.89%、47.64%、42.92%;時(shí)間序列殘差的WRMS平均值減少0.982mm、0.792mm、1.973mm。表2說明了400km空間尺度下區(qū)域?yàn)V波以后各站點(diǎn)在N、E、U 3個(gè)方向上的加權(quán)均方根均有減少。圖5(a)、5(b)、5(c)WRMS減少百分比曲線表明,在E、N、U 3個(gè) 坐標(biāo)分量上WRMS平均減少約44.51%、39.72%、37.05%。

        表1 200 km區(qū)域空間濾波前后殘差坐標(biāo)時(shí)間序列WRMS對比

        表2 400 km區(qū)域空間濾波前后殘差坐標(biāo)時(shí)間序列WRMS對比

        對尺度大于600 km的空間進(jìn)行區(qū)域?yàn)V波,通過表3、表4分析得出各站點(diǎn)在N、E、U3個(gè)方向上的加權(quán)均方根均有減少,但減少幅度明顯小于420 km區(qū)域上的各點(diǎn)。這些站點(diǎn)中MNMC、MASW、HUNT、RNCH在E、N或U方向未濾波的時(shí)候WRMS有超過10,在這種情況下判定這些點(diǎn)本地效應(yīng)嚴(yán)重。

        表3 600 km區(qū)域空間濾波前后殘差坐標(biāo)時(shí)間序列WRMS對比

        表4 1 000 km區(qū)域空間濾波前后殘差坐標(biāo)時(shí)間序列WRMS對比

        根據(jù)地震記錄及SOPAC發(fā)布的階躍日志(ftp://sopac-ftp.ucsd.edu/pub/gamit/setup/siteOffsets.txt)可知,2003-12-22 San Simeon 6.5級(jí)地震和2004-09-28 Park Field 6.0級(jí)地震后,上述站點(diǎn)發(fā)生不同程度的階躍,與上述WRMS異常判斷相符合。在利用PCA提取共模誤差時(shí),由于階躍的影響,使得主分量中包含本地效應(yīng),影響共模誤差提取,因此在進(jìn)行統(tǒng)計(jì)分析時(shí)應(yīng)該剔除這些階躍站點(diǎn)。圖5(d)、5(e)、5(f)為濾波前后WRMS減少百分比曲線,結(jié)果表明,在E、N、U3個(gè)坐標(biāo)分量上,WRMS平均減少約21.44%、6.88%、26.26%。根據(jù)兩種尺度空間上所有站點(diǎn)WRMS的比較發(fā)現(xiàn),各個(gè)方向平均WRMS減少的百分比,420 km尺度下均比1 000 km尺度下好很多,而且600 km、1 000 km空間區(qū)域?yàn)V波中WRMS的減少很多都在20%以下,420 km的時(shí)候大部分點(diǎn)的WRMS減少都在30%~40%,200 km的時(shí)候點(diǎn)的WRMS減少都在40%以上。結(jié)合圖3可以得出,空間尺度大于600 km的時(shí)候,PCA區(qū)域疊加濾波方法效果不明顯。

        圖5 400 km與1 000 km各站點(diǎn)濾波前后WRMS比較Fig.5 The site WRMS of 400km and 1 000 km before and after the filter

        4 結(jié) 語

        本文采用主成分分析的時(shí)空濾波方法,對不同空間尺度的GPS網(wǎng)進(jìn)行共模誤差分離,探討共模誤差的空間變化規(guī)律。結(jié)果表明,當(dāng)GPS網(wǎng)在200 km尺度范圍內(nèi)時(shí),PCA空間濾波效果最為明顯;當(dāng)GPS網(wǎng)的尺度增大到400 km~500 km時(shí),在這個(gè)尺度區(qū)間中PCA濾波仍可以較好地提取出共模誤差;當(dāng)GPS網(wǎng)尺度大于600 km時(shí),PCA空間濾波方法效果不是很理想,應(yīng)該選取其他空間濾波方法(例如相關(guān)系數(shù)加權(quán)濾波)、分塊區(qū)域?yàn)V波方法等。本文研究結(jié)果為大尺度GPS網(wǎng)共模誤差的提取提供了一些參考,尤其是為區(qū)域尺度的劃分提供了依據(jù)。

        [1] 田云鋒.GPS坐標(biāo)時(shí)間序列中的異常高頻周期性噪聲[J].測繪科學(xué), 2011, 36(1):26-28(Tian Yunfeng, Anomalous High Frequency Seasonal Noises in GPS Positions Time Series[J]. Science of Surveying and Mapping, 2011,36(1):26-28)

        [2] Wdowinski S, Bock Y, Zhang J, et al. Southern California Permanent GPS Geodetic Array: Spatial Filtering of Daily Position for Estimating Coseismic and Post Seismic Displacements Induced by the 1992 Landers Earthquake[J]. Journal of Geophysical Research: Solid Earth(1978-2012), 1997, 102(B8):18 057-18 070

        [3] 田云峰, 沈正康.GPS觀測網(wǎng)絡(luò)中共模分量的相關(guān)加權(quán)疊加濾波[J].地震學(xué)報(bào).2011,33(2):198-208(Tian Yunfeng, Shen Zhengkang. Correlation Weighted Stacking Filtering of Common-Mode Component in GPS Observation Network[J].Acta Seismologica Sinica, 2011,33(2):198--208)

        [4] 楊博, 占偉, 劉志廣, 等. 多核函數(shù)法對GNSS 大空間域共模誤差的識(shí)別[J].測繪科學(xué)技術(shù)學(xué)報(bào).2014, 31(2):127-132(Yang Bo, Zhan Wei,Liu Zhiguang,et al.The Identification of GNSS Large Space Domain Common-Mode Error Based on Multi-Kernel Function Method[J]. Journal ofGeomatics Science and Technology. 2014, 31(2):127-132)

        [5] 王敏, 張祖勝, 許明元,等.2000國家GPS大地控制網(wǎng)的數(shù)據(jù)處理和精度評估[J].地球物理學(xué)報(bào), 2005,48(4):817-823(Wang Min, Zhang Zusheng, Xu Mingyuan, et al. Data Processing and Accuracy Analysis of National 2000 GPS Geodetic Control Network[J].Chinese J Geophys, 2005,48(4):817-823)

        [6] Dong Danan. Spatiotemporal Filtering Using Principal Component Analysis and Karhunen-Loeve Expansion Approaches for Regional GPS Network Analysis[J]. Journal of Geophysical Research:Solid Earth,2006,111(B3)

        [7] 黃立人.GPS臺(tái)站坐標(biāo)分量時(shí)間序列的噪聲特性分析[J].大地測量與地球動(dòng)力學(xué), 2006,26(2):31-33(Huang Liren. Noise Properties in Time Series of Coordinate Component at GPS Fiducial Stations[J]. Journal of Geodesy and Geodynamics, 2006,26(2):31-33)

        [8] 胡守超,伍吉倉,孫亞峰.區(qū)域GPS網(wǎng)三種時(shí)空濾波方法的比較[J].大地測量與地球動(dòng)力學(xué),2009,29(3):95-99(Hu shouchao, Wu Jicang, Sun Yafeng. Comparison Among Three Spatiotemporal Filtering Methods for Regional GPS Networks Analysis[J]. Journal of Geodesy and Geodynamics, 2009,29(3):95-99)

        [9] Nikolaidis R. Observation of Geodetic and Seismic Deformation with the Global Positioning System[D], San Diego: University of California, 2002

        [10]Jiang W P, Li Z, Dam T V, et al. Comparative Analysis of Different Environmental Loading Methods and Their Impacts on the GPS Height Time Series [J]. Journal of Geodesy,2013, 87(7): 687-703

        Spatial Response Analysis for GPS Network CME

        GONGGuodong1HEXiaoxing1,3HUAXianghong1,2SHUYing1DINGLinghang1

        1 School of Geodesy and Geomatics, Wuhan University, 129 Luoyu Road,Wuhan 430079, China 2 Hazard Monitoring and Prevention Research Center, Wuhan University,129 Luoyu Road,Wuhan 430079, China 3 School of Civil Engineering and Architecture, East China Jiaotong University, 808 East-Shuanggang Street,Nanchang 330013, China

        In GPS regional network coordinate time series, common mode errors (CME) are widespread. In order to extract CME, a sub-region is divided to get them superimposed and filtered. However, there is no principle to rely on for dividing sub-regions. According to the CME mentioned above, regarding three stations in southern California as centers, we divide the space into nine different spaces within different scales, carrying out the PCA space filtering method, to get the CME in each space scale. After analyzing and comparing the CME in each space scale, we can get the spatial response characteristics in CME, which provides a principle for the dividing sub-regions in the scale field.

        GPS coordinate time series; CME; PCA; spatial filtering; spatial space

        National Natural Science Foundation of China, No.41674005,41374011,41464001; Open Fund of Key Laboratory for Digital Land of Jiangxi Province, No.DLLJ201605.

        HUA Xianghong, professor,PhD supervisor,majors in application and process of GNSS data, analyzing and processing the 3D laser scanning data, E-mail: xhhua@sgg.whu.edu.cn.

        2015-10-28

        項(xiàng)目來源:國家自然科學(xué)基金(41674005,41374011,41464001);江西省數(shù)字國土重點(diǎn)實(shí)驗(yàn)室開放研究基金(DLLJZ01605)。

        龔國棟,碩士生,主要從事GNSS數(shù)據(jù)處理及精密工程測量研究,E-mail:1264392842@qq.com。

        花向紅,教授,博士生導(dǎo)師,主要從事GNSS數(shù)據(jù)處理及應(yīng)用、三維激光掃描數(shù)據(jù)處理與質(zhì)量評價(jià)研究,E-mail:xhhua@sgg.whu.edu.cn。

        10.14075/j.jgg.2016.11.003

        1671-5942(2016)011-0951-07

        P228

        A

        About the first author:GONG Guodong, postgraduate,majors in high-rate GPS data processing and its application on engineering,E-mail:1264392842@qq.com.

        猜你喜歡
        區(qū)域
        分割區(qū)域
        探尋區(qū)域創(chuàng)新的密碼
        科學(xué)(2020年5期)2020-11-26 08:19:22
        基于BM3D的復(fù)雜紋理區(qū)域圖像去噪
        軟件(2020年3期)2020-04-20 01:45:18
        小區(qū)域、大發(fā)展
        商周刊(2018年15期)2018-07-27 01:41:20
        論“戎”的活動(dòng)區(qū)域
        區(qū)域發(fā)展篇
        區(qū)域經(jīng)濟(jì)
        關(guān)于四色猜想
        分區(qū)域
        公司治理與技術(shù)創(chuàng)新:分區(qū)域比較
        狠狠色婷婷久久综合频道日韩| 538在线视频| 日本高清无卡一区二区三区| 日韩av一区二区三区精品久久| 日本av亚洲中文字幕| 亚洲av综合一区二区在线观看| 肉体裸交137日本大胆摄影| 国产va免费精品高清在线| 欧美性受xxxx黑人xyx性爽 | 久久99精品久久久久久9蜜桃| 中文无码一区二区不卡αv| 日韩精品成人一区二区三区| 中文字幕人成人乱码亚洲| 国产粉嫩美女一区二区三| 美利坚合众国亚洲视频| 亚洲悠悠色综合中文字幕| 曰韩无码无遮挡a级毛片| 国产熟女露脸大叫高潮| 国产午夜无码视频免费网站| 亚洲精品成人av观看| 国产精品综合色区av| 久久久亚洲av午夜精品| 国产精品国产三级国产av中文| 亚洲av成人片色在线观看高潮| 亚洲av日韩综合一区在线观看| 亚洲欧美日韩高清专区一区| 国产精品激情综合久久| 国产精品成人黄色大片 | 欧美不卡一区二区三区| 国产熟妇高潮呻吟喷水| 精品2021露脸国产偷人在视频| 国产亚洲精品日韩香蕉网 | 少妇性俱乐部纵欲狂欢少妇| 专干老肥熟女视频网站300部| 亚洲av无码乱码国产精品fc2| 思思99热| 国产在线观看一区二区三区av| 26uuu在线亚洲欧美| 无码一区二区三区亚洲人妻 | 久久老熟女一区二区三区| 伊人中文字幕亚洲精品乱码|