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

基于單極化TerraSAR-X影像提取建筑區研究

2016-05-25 00:37:04丹,盧剛,陳
地理與地理信息科學 2016年1期
關鍵詞:特征建筑

蔣 丹 丹,盧 剛,陳 成

(1.中國礦業大學環境與測繪學院,江蘇 徐州 221116;2.江蘇省測繪工程院,江蘇 南京 210013)

基于單極化TerraSAR-X影像提取建筑區研究

蔣 丹 丹1,盧 剛2,陳 成2

(1.中國礦業大學環境與測繪學院,江蘇 徐州 221116;2.江蘇省測繪工程院,江蘇 南京 210013)

基于多時相單極化TerraSAR-X影像,提取后向散射系數特征影像、變差函數紋理特征影像、干涉特征影像,并對其波段組合,在基于eCognition軟件多尺度分割的基礎上,選擇訓練樣本對象,充分挖掘樣本對象的光譜、紋理特征信息,通過建立特征分離性指標,降低特征維數,最后應用面向對象的分類方法實現了建筑區的自動提取,提取精度達93.50%。這一研究表明特征組合影像具有較大的分割優勢,特征的優化選擇在減少特征冗余的同時維持了較高的信息提取精度,且彌補了eCognition軟件紋理特征計算慢的不足。

TerraSAR-X影像;建筑區提取;影像分割;降低特征維數;面向對象分類

0 引言

近些年,隨著城市擴張和城市建設的迅猛發展,運用遙感技術手段進行城市地理國情監測具有重要的應用價值。建筑區作為城市地理國情監測的一個重要的地理空間要素,實現城市建筑區的提取在城市規劃、土地利用分析、災害監測等方面具有重要意義[1]。在多云、多雨、多霧的天氣狀況下,難以獲取清晰的光學影像,而合成孔徑雷達(Synthetic Aperture Radar,SAR)具有獨特的全天候、全天時數據獲取能力,數據質量幾乎不受天氣狀況與云層覆蓋的影響,這就從一定程度上彌補了光學影像的缺陷[2]。隨著高分辨率SAR衛星(TerraSAR-X、COSMO-SkyMed、ALOS-2、Sentinel-A)等的成功發射,獲取高分辨率SAR影像變得越來越方便,SAR影像也開始廣泛應用于城市土地利用調查中,建筑區作為城市結構中一個典型的地物類型,利用高分辨率SAR影像進行建筑區提取的研究得到關注。

目前,許多學者基于SAR影像進行了建筑區提取研究。Ban等[3]結合多景時間序列的Radarsat-1幅度圖像與紋理信息,應用最大似然分類器和人工神經網絡分類器對城市目標進行了分類。吳樊等[4]利用灰度共生矩陣計算高分辨率SAR圖像的紋理特征,采用非監督聚類分析的方法提取居民區。趙凌君等[5]提出一種基于變差函數紋理特征的高分辨率SAR圖像建筑區提取方法。徐佳等[6]提出一種綜合利用灰度和紋理特征的高分辨率星載SAR圖像建筑區提取方法。張舞燕等[7]基于SAR的相干系數圖像對硬目標相干系數較高的特點,進行了城市建成區邊界的提取。但以上方法較少使用SAR影像非灰度特征參與影像分割并進行面向對象的建筑區提取;另外,在挖掘SAR影像對象特征時,以上研究很少進行特征的優選,容易引起一定的數據特征冗余。

針對以上問題,本文以常州市武進區為研究區,提出了一種基于單極化TerraSAR-X影像提取建筑區的方法。在常規的多尺度分割方法的基礎上,嘗試特征影像輔助分割,通過優選樣本對象特征,采用面向對象的分類方法,實現建筑區輪廓的提取,為城市建筑區的變化監測提供一定的技術參考。

1 研究區概況與數據源

1.1 研究區概況

研究區位于江蘇省常州市武進區東南部,地處長江三角洲太湖平原西北部,屬于我國工業化和城市化發展較快的區域,地理位置為31°36' ~ 31°39'N 、119°29' ~ 120°1'E,該區域緊鄰太湖,自然條件優越;地勢較緩,以平原為主;主要土地利用類型有建筑用地、道路、植被、水體等。該區域的建筑以高密度低矮房屋為主,零星分布有多層及高層房屋、獨立房屋等。

1.2 數據源及預處理

本研究選用的是研究區的TerraSAR-X影像、DEM數據、資源三號多光譜數據。其中,TerraSAR-X影像選用的是TSX SAR L1B多時相單視斜距復影像,具體數據說明如表1所示;DEM數據是“十一五”期間生產的江蘇省5 m分辨率的數字高程模型;資源三號多光譜影像數據作為目視解譯檢驗樣本的參考數據。本文數據的預處理主要是利用SARScape 5.1軟件完成對多時相TerraSAR-X影像的后向散射系數圖的生成,其基本過程包括:多視處理、影像配準、多時相濾波、地理編碼與輻射定標、影像裁剪。其中地理編碼與輻射定標的過程要輔以預處理后的參考DEM數據參與,參考DEM數據的預處理主要是對原始的DEM數據先進行上采樣至3 m分辨率,然后進行坐標系的轉換,生成WGS-84坐標系下的DEM數據。

表1 使用的TerraSAR-X數據

Table 1 TerraSAR-X data applied

獲取時間成像模式極化方式軌道方向入射角(°)距離向分辨率(m)方位向分辨率(m)2013.01.06StripMapHHDescending21.483.213.302013.01.17StripMapHHDescending21.463.213.302013.03.02StripMapHHDescending21.453.213.302013.07.01StripMapHHDescending21.463.213.30

2 研究方法

本研究主要包括影像特性分析與提取、影像分割、特征提取與優選、精度評價等,首先在預處理TerraSAR-X影像的基礎上,從影像的特性出發,提取影像的特征影像并對其波段組合;然后基于eCognition 8.7軟件進行多尺度分割,將影像分割為一個個互不重疊的同質對象,選擇訓練樣本,提取影像對象特征并進行特征優化組合,選擇合適的監督分類方法實現研究區建筑區的自動提取;最后以研究區樣本點的目視解譯結果為參考數據,對建筑區的提取結果進行精度評價,主要技術流程如圖1所示。

圖1 技術流程

Fig.1 Methodological flowchart adopted in this study

3 TerraSAR-X影像特性分析與提取

3.1 后向散射特性分析

隨著TerraSAR-X影像獲取時相的改變,觀測地物目標參數隨之改變,從而影響雷達后向散射的性能,利用這些變化或者不變的特性可以將不同的地物類型區分開來。研究區的建筑區在2013年1-7月時段保持穩定的形態,而對于其他地物很難保持時間上的相干性,所以本研究選取3個時相(2013-01-06、2013-07-01、2013-03-02)的水平極化的TerraSAR-X影像進行紅、綠、藍彩色組合顯示。本研究分別選取不同地物類型的感興趣區SAR影像樣本,求取后向散射系數的平均值(圖2),建筑區的后向散射系數最高,植被次之,水體最低,并且建筑區因其穩定的結構與介電常數使后向散射特性隨時相變化不大,具有一定的時相穩定性,因此后向散射特性可以作為建筑區、非建筑區區分的一個可選特征。

圖2 不同地物的平均后向散射系數

Fig.2 The mean backscatter coefficient of different land cover types

3.2 紋理特性分析

3.2.1 變差函數理論 變差函數,又稱為半變差函數(Semivariogram Function),定義二維平面內的實值區域化變量{f(x),x∈D},即定義于二維空間R2內子集D的實值隨機過程,區域化變量f(x)與f(x+d)(同時包含兩點的距離與方向信息)兩點之差的方差的一半,即:

(1)

式中:兩個區域化變量的值僅僅與兩個點間的空間距離相關,γ(d)為變差(或半變差)函數,即指區域化變量f(x)在空間上距離為d的位置x和x+d兩處的差值方差的一半。

對于離散的柵格數據,實驗變差函數定義為:

(2)

式中:P(d)指觀測窗口內距離為d的點對數目,估計值γ*(d)是實驗變差函數,變差函數刻畫了空間樣本間的相關性,可反映圖像數據的隨機性與結構性[8]。

3.2.2 變差函數計算方法 變差函數用于紋理分析的計算,確定步長d、窗口寬度M、計算方向,利用式(2)計算窗口M內所有間距為d的點對的變差函數值,取其平均值作為計算窗口中心點的變差函數值,最后遍歷整幅影像得到影像的變差函數紋理特征圖。對于影像而言,一般計算4個點對方向(0°、45°、90°、135°)的變差函數均值[5]。以窗口M=5,像素間距d=1為例,圖3給出了0°、45°、90°、135°方向變差函數紋理特征計算的示意圖。圖3中虛線框為計算窗口的大小,窗口寬度為M,黑色實心方塊代表當前的像素(x0,y0),即計算窗口的中心,黑色箭頭代表的是參與計算的像素點對數目。

圖3 變差函數紋理特征計算示意

Fig.3 Calculation of semivariogram function texture feature

3.2.3 變差函數分析建筑區的周期性紋理 建筑區在高分辨率城區SAR影像上呈現有一定周期規律性的紋理,而道路、植被、水體等區域的紋理相對均勻。在預處理后的2013年1月6日的單極化TerraSAR-X影像上選取建筑區、道路、植被、水體,分別計算其變差函數值。圖4顯示了建筑區、道路、植被、水體樣本在單極化高分辨率TerraSAR-X影像上的變差函數曲線。從圖4中可以發現,建筑區屬于強散射的不均勻區域,具有很強的非相似性,其變差函數曲線要遠遠高于植被區、道路、水體,并且隨著間距d的增加周期性地出現峰值與谷點,類似于正弦函數曲線的周期性;道路的變差函數值要低于建筑區,并且會上下波動;植被區的變差函數曲線隨著間距d的增加,先緩慢上升后趨于平穩;水體因為鏡面反射回波很弱,像素的灰度值較低,隨著間距d的增加,變差函數曲線逐漸趨于直線。因此,SAR影像的變差函數紋理特征可以作為區分建筑區、非建筑區的一個可選特征。從地物類別的可分離程度方面考慮,變差函數曲線在d=5的時候達到第一個峰值,此時建筑區與非建筑區(植被、水體、道路)的變差函數數值差異最大,對建筑區與非建筑區的區分性較強,則選擇變差函數的變程d=5。窗口內點對數目太小容易影響變差函數值的準確性,窗口寬度至少為3d~5d,本文選擇窗口寬度M=5d=25,此時的窗口圖像能夠突顯建筑區明暗相間的紋理特征。

圖4 不同地物的變差函數曲線

Fig.4 Semivariogram function of different land cover types

3.3 相干性分析

雷達數據包含有強度與相位信息,雷達的相位信息主要應用在利用雷達干涉技術提取數字高程模型、地表形變監測等。干涉相干作為一個關鍵的量值,被越來越多地應用在土地利用分類、目標識別與檢測中。復數SAR影像之間通常采用相干系數度量影像間的相似程度以及干涉條紋圖的質量,它反映了地面各種散射體的重要信息[9,10]。城市環境中,相干系數圖像對硬目標(建筑物、道路)具有相干系數高的特點,對非硬目標(植被、水體)具有相干系數低的特點,因此可以利用不同地物類型相干系數圖像的特點進行類別的區分。

選取實驗區時間基線為12 d的兩景水平極化的TerraSAR-X影像單視復數(SLC)數據(2013-01-06為主影像, 2013-01-17為從影像)在SARScape5.1軟件中進行干涉處理,生成相干系數圖。建筑區的相干性比較穩定,相干系數較大。為了更直觀地分析不同的地物類型的相干性,本文選取了典型地物樣本(建筑區、道路、植被、水體)進行相干系數的統計分析(圖5)。從圖5可以發現,建筑區的相干系數明顯高于非建筑區(道路、植被、水體),平均相干系數值接近0.8,因此相干系數特征也可以作為區分建筑區與非建筑區的一個可選特征。

圖5 相干系數統計

Fig.5 Interferometry coherence statistics

4 實驗與分析

4.1 影像分割

將提取得到的后向散射特征圖、變差函數紋理特征圖、相干性圖進行特征波段組合,得到研究區含有5個波段的影像圖,各波段圖層有:3個時相的后向散射系數圖層(2013-01-06、2013-07-01、2013-03-02)、變差函數紋理值圖層、干涉系數圖層。選用eCognition8.7軟件對該影像數據進行多尺度分割試驗,利用集成在該軟件中的工具(Estimation of Scale Parameter,ESP)確定分割參數,該工具是通過分割后影像對象的平均局部方差(Local Variance,LV)與相鄰尺度平均局部方差變化率(Rate of Change,ROC)曲線提取最優分割參數的參考值[11]。本文多尺度實驗參數設置為:影像各波段權重設為1,初始的分割尺度(Scale)為5,分割尺度的遞增步長為5,循環次數為40,光譜因子(Colour)與形狀因子(Shape)的權重比為0.9∶0.1,光滑度因子(Smoothness)與緊致度(Compactness)因子的權重比為0.5∶0.5。多尺度分割完成得到的ROC和LV曲線如圖6所示,一般情況下,最優分割尺度出現在ROC-LV變化劇烈的地方,其參考值出現在局部方差的峰值且局部方差變化率開始趨于下降的地方[12]。特征組合影像表現出了較大的分割優勢,從圖6中可以看出最佳分割尺度的可能值為70、105、180,分割的局部效果如圖7(見封3)所示。從建筑區提取的角度,當分割尺度Scale為70時,建筑區被分割得過于零碎,存在過分割現象;當分割尺度Scale為105時,建筑區分割較為完整,為建筑區的最優分割尺度;當分割尺度Scale為180時,建筑區與周圍的植被聚合在一起,存在少量欠分割現象。

圖6 ROC和LV變化曲線

Fig.6 Graph of ROC and LV with a increasing scale parameter

4.2 特征提取與優選

在多尺度分割的基礎上,研究區的影像被分割成1 267個對象,選擇對象總量的3%作為樣本,包括建筑區、道路、植被、水體4類樣本。eCognition軟件提供了很多對象特征,包括光譜特征、紋理特征、形狀特征、鄰域特征、語義對象關系等,本文提取的特征有光譜特征和紋理特征。光譜特征包括最大差分、亮度以及影像各波段的均值和標準差;紋理特征包括基于灰度共生矩陣(Gray Level Co-occurrence Matric,GLCM)的均值、方差、對比度、非相似度、同質性、角二階矩、熵、相關性,以及基于灰度差分矢量(Grey Level Difference Vector,GLDV)的角二階矩、熵、對比度、均值,這些紋理特征的計算方向包括0°、45°、90°、135°、全方向。利用選擇的樣本對象,基于提取的72個特征進行面向對象分類的一個重要步驟就是對這些特征的篩選,但是直接在eCognition軟件對紋理特征的計算耗費時間太長,本文采用導出這些樣本對象的特征方法,利用改進的SEaTH(Separability and Thresholds)算法進行特征的優選,其核心思想包括兩個方面:其一是去除特征之間的相關性,其二是綜合考慮類內距離和類間距離的特征選擇[13]。

(1)去除特征之間的相關性。通過設置相關系數閾值,調節入選特征的數目。閾值設置越小,則入選的特征個數越少;反之,則入選的特征個數越大。當兩個特征間的相關系數大于閾值時,則舍棄與其他特征具有較強的相關性且方差較小的特征。從圖像信息量的角度,方差越大,則說明特征中含有的分類信息越多;方差越小,則說明特征中含有的分類信息越少[14]。通過反復試驗,確定本次實驗的相關系數閾值為0.95。

(2)綜合考慮類內距離和類間距離的特征選擇。根據不同類別的樣本特征值計算類內距離,依據同類樣本特征值計算類內距離[15]。假設去除相關性的特征子集為Fs=(f1,f2,…,fs),對象類別為Ct=(C1,C2,…,Ct),從各類別中選取的樣本對象個數為Nt=(n1,n2,…,nt)。以C1和C2兩類為例,對其特征選擇過程進行詳細描述。

①計算類間距離J。計算C1和C2兩類樣本對象的某個特征(如fj,j=1,…,s)的均值mi與方差σi,i=1,2,代入式(3)和式(4),得到特征fj對應的C1與C2的類間距離J,選取類間距離較大的前10個特征構成特征子集,從而進行下一步的篩選。

(3)

(4)

②特征歸一化處理。由于所選的特征值在數量級上表現不同,需要對上一步所選的特征進行歸一化處理,使其取值范圍為[0,1],計算公式如下:

(5)

③計算加權類內距離D。遍歷C1與C2類中的所有樣本,計算每個樣本和其他同類樣本的某個特征值(如fj,j=1,…,s)的距離并累加,分別獲取到C1與C2的類內距離d1和d2,如式(6)所示,依據類別C1與C2的樣本數目賦予類內距離d1和d2相對應的權重,從而得到加權類內距離D,如式(7)所示。

(6)

D=(n1d1+n2d2)/(n1+n2)

(7)

④構建特征篩選指標Jf。根據類間距離J大、加權類內距離D小的準則,構建篩選指標Jf,如式(8)所示。Jf值越大,說明特征之間的分離性越好;反之,說明特征之間的分離性越差。對Jf值進行降序排列,選擇排在前3個的若干特征參與分類。

(8)

通過試驗,最終得到建筑區樣本對象與其他類樣本對象進行類別區分的5個優選特征:GLCM對比度(45°)、GLCM標準差(135°)、最大差分、變差函數紋理均值、干涉系數均值(表2)。

表2 優選特征指標值

Table 2 Selected best feature indexes

兩兩類別優選特征類間距離(J)加權類內距離(D)篩選指標(Jf)建筑區-道路最大差分1.90400.027569.2364變差函數紋理均值1.65890.051532.2117GLCM對比度(45°)1.47730.047031.4319建筑區-植被最大差分1.99990.0143139.8531 變差函數紋理均值1.99860.029467.9796干涉系數均值1.90680.053135.9096建筑區-水體變差函數紋理均值1.99930.023186.5498最大差分1.94410.027171.7380GLCM標準差(135°)1.86530.032557.3938

4.3 建筑區提取與精度評價

根據優選出的特征構建特征空間,選擇支持向量機(SVM)的面向對象分類方法,將研究區劃分為建筑區與非建筑區,分類結果如圖8a(圖8見封3)所示,其中紅色區域表示建筑區。根據研究區的地理范圍,平均劃分為400個矩形樣方,取樣方的中心點為檢驗樣本點,檢驗樣本點分布如圖8b所示。以研究區的資源三號衛星多光譜影像為參考數據,對這些檢驗樣本點進行人工目視解譯,確認檢驗樣本的類別屬性,建立混淆矩陣,計算相關精度評價指標評價研究區的建筑區結果,如表3所示。

表3 總體分類精度統計

Table 3 Overall classification accuracy statistics

用戶類參考類建筑區非建筑區合計用戶精度建筑區12071270.94非建筑區192542730.93合計139261400生產者精度0.860.97總體精度=93.50%Kappa系數=0.85

5 結語

本文針對高分辨率單極化TerraSAR-X影像,從影像的基本特性出發,組合影像特征波段,在面向對象的多尺度分割技術的基礎上,選擇訓練樣本,充分挖掘影像對象的光譜、紋理、干涉信息,建立分離性指標,確立優選特征,采用支持向量機的面向對象分類方法,實現了研究區建筑區輪廓的自動提取,分類的總體精度為93.50%,Kappa系數為0.85。

本研究的提取方法一方面可以滿足一般性研究分析的精度要求,另一方面有助于城市建筑區的變化監測,同時也為其他高分辨率SAR影像進行城市建筑區提取提供一定的參考依據。但本研究只是針對單極化TerraSAR-X影像對以高密度低矮房屋為主的建筑區進行了提取,對于更加復雜的建筑區的提取方法還有待完善,同時多極化的TerraSAR-X影像進行建筑區提取還有待進一步研究。

致謝:江蘇省測繪工程院提供了實驗數據,此致謝忱!

[1] FENG T,ZHAO J.Review and comparison:Building extraction methods using high-resolution images[A].Information Science and Engineering (ISISE),2009 Second International Symposium onIEEE[C]2009.419-422.

[2] 溫禮,程博,柴淵,等.SAR遙感數據監測土地利用變化的研究[J].測繪科學,2014(6):65-69.

[3] BAN Y,WU Q.Radarsat SAR data for land-use/land-cover classification in the rural-urban fringe of the Greater Toronto area[A].8th Conference on Geographic Information Science[C].Estoril,2005.43-50.

[4] 吳樊,王超,張紅.基于紋理特征的高分辨率SAR影像居民區提取[J].遙感技術與應用,2005,20(1):148-152.

[5] 趙凌君,高貴,匡綱要.基于變差函數紋理特征的高分辨率SAR圖像建筑區提取[J].信號處理,2009,25(9):1433-1442.

[6] 徐佳,陳媛媛,黃其歡,等.綜合灰度與紋理特征的高分辨率星載SAR圖像建筑區提取方法研究[J].遙感技術與應用,2012,27(5):692-698.

[7] 張舞燕,劉清臣,孟秀軍.基于SAR相干系數圖像的城市邊界提取[J].測繪與空間地理信息,2014,37(5):56-59.

[8] 王燕紅,程博,尤淑撐,等.基于改進變差函數的高分辨率SAR圖像建筑區提取[J].遙感信息,2014,29(2):1-6.

[9] 舒寧.微波遙感原理(修訂版)[M].武漢:武漢測繪科技大學出版社,2003.

[10] 王超,張紅,劉智.星載合成孔徑雷達干涉測量[M].北京:科學出版社,2002.

[12] 陳鐵橋,陳杰,梅小明,等.基于多層次場聚類的高分辨率遙感影像分割方法[J].地理與地理信息科學,2013,29(6):10-13.

[13] 余曉敏,湛飛并,廖明生,等.利用改進SEaTH算法的面向對象分類特征選擇方法[J].武漢大學(信息科學版),2012,37(8):921-924.

[14] 張秀英,馮學智,江洪.面向對象分類的特征空間優化[J].遙感學報,2009,13(4):664-669.

[15] 王賀,陳勁松,余曉敏.面向對象分類特征優化選取方法及其應用[J].遙感學報,2013,17(4):822-829.

Built-up Areas Extraction Based on Single Polarization TerraSAR-X Images

JIANG Dan-dan1,LU Gang2,CHEN Cheng2

(1.School of Environment Science and Spatial Informatics,China University of Mining and Technology,Xuzhou 221116;2.Jiangsu Province Surveying & Mapping Engineering Institute,Nanjing 210013,China)

This paper shows a object-oriented method that the extraction of urban built-up areas based on multi-temporal single polarization TerraSAR-X images.Firstly,the backscattering coefficient images,variogram texture image and interferometric coherent image were combined,and the combined image was segmented with a multiresolution segmentation algorithm based on eCognition software.And then choosing the train samples,extracting the spectral and textural features of sample objects,sample objects separability index was established in order to reduce feature dimensions.Finally,the built-up areas were extracted automatically with a object-oriented classification approach,and total accuracy is up to 93.50%.The results indicates that combined image has great advantage of segmentation,optimal selection of features reduces redundant features,maintains the high accuracy of information extraction,and makes up for the inadequacy that texture characteristics calculation is slow based on eCognition software.

TerraSAR-X images;built-up areas extraction;image segmentation;reduce feature dimensions;object-oriented classification

2015-05-25;

2015-08-24

國家高技術研究發展計劃(863計劃)項目(2013AA122003)

蔣丹丹(1990-),女,碩士研究生,研究方向為遙感數據處理。E-mail:jddcumt@163.com

10.3969/j.issn.1672-0504.2016.01.012

P237

A

1672-0504(2016)01-0060-06

猜你喜歡
特征建筑
抓住特征巧觀察
《北方建筑》征稿簡則
北方建筑(2021年6期)2021-12-31 03:03:54
關于建筑的非專業遐思
文苑(2020年10期)2020-11-07 03:15:36
新型冠狀病毒及其流行病學特征認識
建筑的“芯”
現代裝飾(2020年6期)2020-06-22 08:43:12
山居中的石建筑
現代裝飾(2020年4期)2020-05-20 08:55:08
如何表達“特征”
不忠誠的四個特征
當代陜西(2019年10期)2019-06-03 10:12:04
抓住特征巧觀察
聽,建筑在訴說
主站蜘蛛池模板: 久久久久九九精品影院| 色男人的天堂久久综合| 热99re99首页精品亚洲五月天| 国产精品一区二区在线播放| 免费中文字幕一级毛片| 亚洲一区二区视频在线观看| 91麻豆国产视频| www.日韩三级| 青草国产在线视频| 欧美日韩导航| 精品人妻系列无码专区久久| 美女扒开下面流白浆在线试听| 女人毛片a级大学毛片免费| 欧洲亚洲欧美国产日本高清| 午夜国产小视频| 91精品网站| 精品国产自在在线在线观看| 国产精品亚洲专区一区| 国产精品主播| 亚洲综合网在线观看| 国产成人高清在线精品| 91免费国产高清观看| 999国产精品| 午夜三级在线| 国产第一页第二页| 午夜a级毛片| 久青草国产高清在线视频| 久青草免费在线视频| 伊人五月丁香综合AⅤ| 免费网站成人亚洲| 一级毛片在线免费视频| 黄色a一级视频| 免费一级无码在线网站| 国产亚洲欧美日韩在线一区| 亚洲日韩在线满18点击进入| 一级毛片免费的| 香蕉eeww99国产在线观看| 女人一级毛片| 小说区 亚洲 自拍 另类| 欧美日韩精品在线播放| 国产玖玖视频| 人妻精品久久无码区| 国产亚洲精品自在久久不卡 | 日本一区二区三区精品视频| AV在线麻免费观看网站| 免费99精品国产自在现线| 露脸国产精品自产在线播| 亚洲国产日韩视频观看| 夜精品a一区二区三区| 麻豆国产精品一二三在线观看| 久久综合婷婷| 国产精品女同一区三区五区| 妇女自拍偷自拍亚洲精品| 国产凹凸一区在线观看视频| 国产成人亚洲无码淙合青草| 中文一区二区视频| 婷婷综合色| 在线观看免费国产| 久久人人97超碰人人澡爱香蕉 | 最新精品久久精品| 91破解版在线亚洲| 欧美啪啪精品| 97青青青国产在线播放| 久久人人妻人人爽人人卡片av| 亚洲毛片网站| 美女高潮全身流白浆福利区| 精品无码一区二区三区电影| 亚洲国产综合精品一区| 亚洲丝袜第一页| 亚洲成人77777| 亚洲天堂网2014| 伊人网址在线| 亚洲AV无码乱码在线观看裸奔| 欧美色99| 亚洲精品无码抽插日韩| 国产精品网拍在线| 性色一区| 国产九九精品视频| 91极品美女高潮叫床在线观看| 99在线视频网站| 91极品美女高潮叫床在线观看| 国产二级毛片|