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

基于BJ-1小衛(wèi)星遙感數(shù)據(jù)的礦區(qū)土地覆蓋變化檢測

2011-01-05 07:56:56杜培軍唐偉成柳思聰
自然資源遙感 2011年3期
關(guān)鍵詞:檢測方法

陳 宇,杜培軍,唐偉成,柳思聰

基于BJ-1小衛(wèi)星遙感數(shù)據(jù)的礦區(qū)土地覆蓋變化檢測

陳 宇,杜培軍,唐偉成,柳思聰

(中國礦業(yè)大學(xué)國土環(huán)境與災(zāi)害監(jiān)測國家測繪局重點(diǎn)實(shí)驗(yàn)室,徐州 221116)

為了評價利用北京一號小衛(wèi)星(BJ-1)遙感數(shù)據(jù)監(jiān)測煤礦區(qū)土地利用/地表覆蓋變化的效果,針對其數(shù)據(jù)特點(diǎn),選擇基于圖像信息運(yùn)算和圖像信息變換的直接變化檢測法以及分類后比較法,對徐州東礦區(qū)2007~2008年土地利用/地表覆蓋變化情況進(jìn)行檢測,以對比、評價各種方法在土地利用/地表覆蓋變化檢測中的應(yīng)用效果和BJ-1數(shù)據(jù)的適用性。結(jié)果顯示,變化矢量分析法的檢測精度最高,其后依次為圖像比值法、圖像差值法和多波段主成分分析法。通過變化檢測,確定徐州市東礦區(qū)土地利用/地表覆蓋變化較為集中的幾個區(qū)域,包括東礦區(qū)北部的青山泉礦、韓橋礦、董莊礦和南部的大黃山礦等地區(qū)。

北京一號小衛(wèi)星(BJ-1);土地覆蓋;變化檢測;礦區(qū);變化矢量分析

0 引言

北京一號小衛(wèi)星(BJ-1)是根據(jù)國家“十五”科技攻關(guān)計劃和高技術(shù)研究發(fā)展計劃(863計劃)的安排,在科技部、北京市政府、國家測繪局和21世紀(jì)空間技術(shù)應(yīng)用股份有限公司等部門的共同支持下,與英國薩瑞衛(wèi)星技術(shù)有限公司(SSTL)合作,于2005年10月27日在俄羅斯成功發(fā)射的[1]。目前,BJ-1遙感數(shù)據(jù)在國土資源調(diào)查監(jiān)測、環(huán)境與災(zāi)害監(jiān)測、城市管理建設(shè)等方面已開展了一系列的示范應(yīng)用研究,并取得了可喜的成果[2]。土地利用/地表覆蓋變化檢測是遙感應(yīng)用的熱點(diǎn)之一,尤其是對環(huán)境多變的煤礦區(qū)而言,實(shí)時監(jiān)測土地利用/覆蓋動態(tài)是進(jìn)行全面監(jiān)管和有效治理礦區(qū)生態(tài)環(huán)境的關(guān)鍵。利用遙感圖像進(jìn)行地表動態(tài)變化檢測的方法可分為以下兩類:①非監(jiān)督式的直接變化檢測法;②監(jiān)督式的基于分類后比較的變化檢測法。由于土地利用/地表覆蓋變化受到空間、時間、地物光譜等多種復(fù)雜因素的影響,任何變化檢測方法都有其優(yōu)勢與不足,還沒有一種方法最優(yōu),并適合于所有情況,即使在相同環(huán)境下,不同方法得到的檢測結(jié)果也各不相同,甚至出現(xiàn)爭議性的結(jié)論[3]。

為了探究國產(chǎn)小衛(wèi)星遙感數(shù)據(jù)在煤礦區(qū)土地利用/地表覆蓋變化監(jiān)測中的應(yīng)用效果,本文以BJ-1小衛(wèi)星多光譜遙感影像為數(shù)據(jù)源,采用多種變化檢測法對徐州東礦區(qū)土地利用/地表覆蓋變化進(jìn)行了檢測,并分析、對比了各種檢測方法的優(yōu)缺點(diǎn)和在本研究區(qū)的應(yīng)用效果及BJ-1數(shù)據(jù)的適用性。

1 研究區(qū)及數(shù)據(jù)源概況

徐州市位于江蘇省西北部,地理坐標(biāo)為116°22'~118°40'E、33°43'~34°58'N,處于蘇、魯、豫、皖四省交界處,交通便利,鐵路、公路四通八達(dá),是全國重要的水路交通樞紐和東西南北經(jīng)濟(jì)聯(lián)系的“十”字路口[4]。徐州市是華東地區(qū)典型的礦業(yè)城市,區(qū)內(nèi)煤炭資源十分豐富,為江蘇省重要的煤炭生產(chǎn)和供應(yīng)基地,也是中國重點(diǎn)煤炭基地之一。

本文的研究區(qū)域?yàn)樾熘菔袞|礦區(qū),包括青山泉礦、大黃山礦、韓橋礦、董莊礦、權(quán)臺礦、旗山礦和鄭莊礦等煤礦。圖1為研究區(qū)2008年BJ-1小衛(wèi)星遙感數(shù)據(jù)的近紅外波段、紅光波段和綠光波段的假彩色RGB合成圖像。

圖1 徐州東礦區(qū)2008年BJ-1衛(wèi)星遙感假彩色合成圖像Fig.1 BJ -1 false-color composite image of east coal mining area in Xuzhou in 2008

試驗(yàn)采用2007年2月2日和2008年2月18日獲取的空間分辨率為32 m的BJ-1多光譜圖像(包括綠光波段、紅光波段和近紅外波段)。由于所選用數(shù)據(jù)的獲取月份相近,因此地物覆蓋類別信息較一致,便于進(jìn)行土地利用/地表覆蓋變化檢測。

數(shù)據(jù)預(yù)處理包括幾何精糾正、圖像裁剪和直方圖匹配。預(yù)處理時,先以2007年圖像為基準(zhǔn),采用Geographic Lat/lon和WGS-84坐標(biāo)系,選擇了33個地面控制點(diǎn),對2008年圖像進(jìn)行幾何配準(zhǔn),配準(zhǔn)后圖像的均方根誤差(RMS)為0.493像元;然后采用二次多項(xiàng)式及最近鄰取值法進(jìn)行重采樣,得到幾何精糾正圖像;接著再依據(jù)研究區(qū)范圍將2幅圖像裁剪為同一區(qū)域;最后以2007年圖像為基準(zhǔn)進(jìn)行直方圖匹配,即得到預(yù)處理后效果最佳的圖像。

2 變化檢測方法

選擇直接變化檢測法和分類后比較法來檢測徐州市東礦區(qū)2007~2008年土地利用/地表覆蓋的變化情況。

2.1 直接變化檢測法

目前,常用的直接變化檢測方法包括基于圖像信息運(yùn)算的方法、基于圖像信息變換的方法和基于模型分析的方法。本文采用前2種方法進(jìn)行實(shí)驗(yàn)[5]。

2.1.1 基于圖像信息運(yùn)算的方法

(1)圖像差值法。該法是應(yīng)用最為廣泛的變化檢測方法。它先通過將配準(zhǔn)后不同時相圖像的對應(yīng)像元的灰度值進(jìn)行相減來生成一幅新的代表不同時相光譜變化的差值圖像;然后用差值圖像的灰度值來表示2幅圖像間的變化程度[6];最后通過選取合適的閾值對差值圖像灰度進(jìn)行分割來確定變化區(qū)域與不變區(qū)域。此方法既可應(yīng)用于單波段圖像,也可用于多波段圖像。

(2)圖像比值法。該法是通過將經(jīng)配準(zhǔn)的不同時相遙感圖像的對應(yīng)波段圖像進(jìn)行逐像元相除,所得到的比值圖像就增強(qiáng)了變化信息[7]。如果像元的亮度值較為接近,則比值為1或近似為1,表明無變化;若像元比值明顯高于或低于1,則認(rèn)為發(fā)生變化,比值的大小取決于發(fā)生變化的性質(zhì)。通常可通過經(jīng)驗(yàn)法判斷比值圖像直方圖上下尾端的閾值來確定變化信息。

(3)植被指數(shù)法。該法是先分別對2個時相圖像求取植被指數(shù),然后對植被指數(shù)進(jìn)行差值處理,進(jìn)而通過選擇合適的閾值來提取出變化信息[6]。

在眾多植被指數(shù)中,歸一化差值植被指數(shù)(NDVI)常被用于土地利用/地表覆蓋變化檢測,利用植被對近紅外波段與紅光波段光譜明顯的響應(yīng)差,可很好地反映地面植被的覆蓋情況。

(4)變化矢量分析法。變化矢量是描述從時相1到時相2光譜向量變化的方向和大小的矢量,每個像元可生成一個具有變化方向和變化強(qiáng)度(大小)2個特征的變化向量[3]。當(dāng)某一像元的變化強(qiáng)度值超出一定范圍時,則認(rèn)為該像元發(fā)生變化;否則,認(rèn)為該像元未發(fā)生變化。

2.1.2 基于圖像信息變換的方法

(1)主成分分析法。主成分分析(PCA)變換檢測法[7]可分為主成分差異法、差異主成分分析法和多波段主成分分析法。

主成分差異法是將2個時相的多光譜圖像先分別進(jìn)行PCA變換,然后將變換后得到的2個第一主成分進(jìn)行差值處理,再對得到的差值圖像通過設(shè)定閾值來提取變化信息。

差異主成分分析法是先將2個時相多光譜圖像的各個波段對應(yīng)圖像分別進(jìn)行相減,然后對差值圖像進(jìn)行PCA變換,變化信息主要集中于第一主成分。

多波段主成分分析法是將2個時相圖像的波段先合成形成一個2倍于原圖像波段數(shù)的新圖像,然后對該圖像進(jìn)行PCA變換。變換結(jié)果的前幾個主成分集中了2個圖像的不變化信息,而后面的幾個主成分則反映2個圖像的變化信息,因此可以通過提取后幾個主成分進(jìn)行合成來得到變化信息。

(2)對應(yīng)成分分析法。對應(yīng)分析(CV)是采用V2統(tǒng)計量來表示波段內(nèi)各個像素與波段間的相關(guān)關(guān)系[8]。設(shè)遙感圖像數(shù)據(jù)有n個波段,每個波段有r個像素,則可構(gòu)造一個矩陣Xrxn。矩陣X中的每一個元素都除以矩陣元素數(shù)值的總和,即構(gòu)造出新矩陣Z,新矩陣元素為

再構(gòu)造一個r×1大小的向量p,使它的每個元素代表矩陣Z的每行元素之和,同時構(gòu)造一個n×1大小的向量q,使它的每個元素代表矩陣Z的每列元素之和。

V2統(tǒng)計量可以描述矩陣Z中行列變量之間的獨(dú)立性,V2統(tǒng)計量為

矩陣U為

若計算出U的特征值和特征向量,則對應(yīng)的分析變換系數(shù)矩陣即為由U的特征向量組成的矩陣。

對應(yīng)分析法的處理過程與主成分分析法類似,即首先分別對2個時相的多光譜圖像進(jìn)行CA變換,然后將變換后得到的2個第一分量進(jìn)行差值處理,再提取變化信息。研究中采用IDL語言編程實(shí)現(xiàn)對應(yīng)成分變換。

2.2 分類后比較法

分類后比較法也是目前應(yīng)用較廣泛的一種遙感變化檢測方法,其原理是對2個不同時相的圖像先采用相同的算法分別進(jìn)行單獨(dú)分類,然后通過對分類結(jié)果進(jìn)行對比來檢測變化信息[9]。

這種方法可以避免直接變化檢測法對2個不同時相圖像的成像條件要基本一致的條件,以及2個圖像間的輻射校正、匹配等問題,對用多源遙感數(shù)據(jù)進(jìn)行變化檢測也具有較好的效果[7];同時,通過選擇合適的分類方法還可以克服地形因素引起的偽變化[3];另外,此方法還可以直接提供變化前后的地物類別信息和變化信息。但是此方法也有一定的局限,它要求2次地物影像分類必須遵循統(tǒng)一的分類標(biāo)準(zhǔn),而變化分析的精度又依賴于地物影像分類的精度,這就會受到利用2個圖像分別進(jìn)行地物影像分類所帶來的誤差影響,從而不可避免地夸大了地表變化的程度[7]。

3 實(shí)驗(yàn)和分析

以BJ-1遙感數(shù)據(jù)為基礎(chǔ),采用上述各種變化檢測方法分別進(jìn)行土地利用/地表覆蓋變化檢測實(shí)驗(yàn),得到的徐州市東礦區(qū)土地利用/地表覆蓋變化結(jié)果如圖2所示。

圖2 9種變化檢測方法的檢測結(jié)果Fig.2 Results of nine change detection methods

圖2 (i)為分類后比較法的檢測結(jié)果,圖中僅選擇變化較明顯且相對重要的區(qū)域顯示。

其中,差值法和比值法均先對各單波段分別進(jìn)行變化檢測,然后將3個波段的檢測結(jié)果進(jìn)行疊加,故最終所得結(jié)果為彩色圖像,它反映了3個波段的綜合變化信息。分類后比較法是先根據(jù)研究區(qū)土地覆蓋的實(shí)際情況,對2個時相的圖像分別選擇訓(xùn)練樣本;然后采用基于CART的決策樹分類方法,將徐州東礦區(qū)的土地覆蓋類型分為農(nóng)田、建筑用地、水體和林地4類;最后選擇測試樣本進(jìn)行分類精度評價,得到的2007年的總分類精度和Kappa系數(shù)分別為84.2%和0.821;2008年的總分類精度和Kappa系數(shù)分別為85.9%和0.837。

對不同變化檢測方法得到的結(jié)果進(jìn)行精度評價,通常采用類似于分類精度評價的方法,即先將變化檢測結(jié)果作為分類圖像,包括“變化”和“不變化”2個類別;然后在原圖像上通過選擇“變化”類別和“不變化”類別的測試樣本,即可構(gòu)造混淆矩陣,并得到總體檢測精度、Kappa系數(shù) 、虛檢率及漏檢率等精度指標(biāo)(見表1)。

表1 不同變化檢測方法檢測結(jié)果精度對比Tab.1 Precision comparison between different change detection methods

通過對各種變化檢測方法的檢測結(jié)果進(jìn)行統(tǒng)計分析,得到的變化區(qū)域面積(單位為km2)以及占地區(qū)總面積的百分比如圖3所示。

圖3 不同變化檢測方法檢測到的變化區(qū)域面積百分比Fig.3 Area percentage of the changed area detected by different methods

由表1和圖3可以看出:

(1)變化矢量分析法的變化檢測精度最高,虛檢率最低,但漏檢率相對偏高,這使得用該方法檢測到的變化區(qū)域面積百分比最低。

(2)圖像差值法和圖像比值法也具有較好的檢測精度,但二者檢測到的變化區(qū)域面積明顯高于其他方法,其原因有二:首先,此2種方法均先利用單波段數(shù)據(jù)分別進(jìn)行差值和比值運(yùn)算,然后再對3個波段圖像的檢測結(jié)果進(jìn)行合成,由于其包括了3個波段圖像的變化信息,因而能充分反映每個波段圖像的地表信息變化,但因綠光波段圖像的質(zhì)量不佳,會對檢測結(jié)果產(chǎn)生影響;其次,這2種方法的漏檢率相對偏低而虛檢率則相對偏高,從而使得檢測到的變化區(qū)域的面積偏大。

(3)在主成分分析的3種方法中,多波段主成分分析法的精度最高,虛檢率和漏檢率也最低,但在所有方法中漏檢率偏高。另2種主成分分析方法的檢測精度十分接近,檢測到的變化區(qū)域面積比例也都在5.9%附近,說明它們在一定程度上檢測效果相當(dāng),可選擇使用;但是與其他方法相比,這2種方法的變化檢測效果并不理想,檢測精度較低,漏檢率和虛檢率都偏高。

(4)對應(yīng)成分分析法的檢測精度相對較低,虛檢率和漏檢率較高。雖然對應(yīng)分析是將多光譜圖像的波譜空間轉(zhuǎn)換為成分空間,前幾個成分集中了大部分信息,且能排除噪聲干擾,但同樣會使部分變化信息被遺漏。此方法對于徐州市東礦區(qū)土地利用/地表覆蓋變化檢測的適用性相對較差。

(5)分類后比較法雖然可以檢測出土地利用/地表覆蓋變化的性質(zhì)和位置(這是其他檢測方法不具有的特點(diǎn)),但在本實(shí)驗(yàn)中,其變化檢測的精度最低。這主要是受2個時相圖像分類精度的影響,導(dǎo)致其土地利用/地表覆蓋變化檢測精度最低,虛檢率最高,檢測出的變化區(qū)域面積百分比最大。在檢測結(jié)果圖(圖2(i))上可以看出,有很多虛檢部分為道路或建筑物輪廓線。此方法對于徐州市東礦區(qū)的土地利用/地表覆蓋變化檢測的適用性較差。

總體而言,9種檢測方法對于利用BJ-1小衛(wèi)星遙感數(shù)據(jù)來對徐州市東礦區(qū)土地利用/地表覆蓋變化進(jìn)行檢測都有一定的適用性,其中變化矢量分析法、圖像差值法、比值法及多波段主成分分析法的檢測效果相對較好;從檢測到的變化區(qū)域面積比例來看,圖像差值和比值法檢測到的變化區(qū)域面積較大,但虛檢率較高;變化矢量分析法和多波段主成分分析法的虛檢率雖然較低,但漏檢率較高。所以,對檢測方法的選擇還要依據(jù)具體情況而定。

4 結(jié)論

(1)本文基于徐州市東礦區(qū)BJ-1小衛(wèi)星遙感數(shù)據(jù),采用9種變化檢測方法進(jìn)行了土地利用/地表覆蓋變化檢測實(shí)驗(yàn),并評價了各種方法對該地區(qū)土地利用/地表覆蓋變化檢測的適用性。結(jié)果表明,變化矢量分析法的檢測精度最高,其次為圖像比值法、圖像差值法和多波段主成分分析法。其他方法的適用性相對較差,植被指數(shù)法主要對植被的變化較為敏感,而對其他地物的敏感度較低;分類后比較法不僅受分類精度影響較大,并且存在較嚴(yán)重的虛檢現(xiàn)象。

(2)通過對各種方法檢測結(jié)果的分析比較,并參照其他資料,最終確定了徐州市東礦區(qū)2007~2008年土地利用/地表覆蓋變化主要集中的幾個區(qū)域,主要包括東礦區(qū)北部的青山泉礦、韓橋礦、董莊礦和南部的大黃山礦等區(qū)域。

(3)隨著新的土地利用/地表覆蓋變化檢測方法的不斷涌現(xiàn),集成多種變化檢測方法的優(yōu)勢和特點(diǎn),以提高整體檢測精度,將是下一步研究的主要方向。

[1] 童慶禧,衛(wèi) 征.北京一號小衛(wèi)星及其數(shù)據(jù)應(yīng)用[J].航天器工程,2007,16(2):1 -5.

[2] 魏 軍,張萍萍,紀(jì)中奎.“北京一號”小衛(wèi)星數(shù)據(jù)在北京森林生態(tài)系統(tǒng)資源動態(tài)分析中的應(yīng)用[J].生態(tài)學(xué)報,2008,28(7):3465-3469.

[3] 李月臣,楊 華,劉春霞,等.土地覆蓋變化遙感監(jiān)測方法[J].水土保持研究,2006,13(1):209 -216.

[4] 周興東,于勝文,趙長勝.徐州市土地利用/覆蓋動態(tài)變化分析[J].測繪通報,2008(1):33 -36.

[5] 杜培軍,胡召玲,郭達(dá)志,等.工礦區(qū)陸面演變監(jiān)測分析與調(diào)控治理研究[M].北京:地質(zhì)出版社,2005.

[6] 張振龍,曾志遠(yuǎn),李 碩,等.遙感變化檢測方法研究綜述[J].遙感信息,2005(5):64-66.

[7] 李向軍.遙感土地利用變化檢測方法探討[D].北京:中國科學(xué)院遙感應(yīng)用研究所,2006.

[8] 龔 潔,張景雄,申邵洪.基于對象的對應(yīng)分析在高分辨率遙感影像變化檢測中的應(yīng)用[J].武漢大學(xué)學(xué)報,2009,34(5):544-547.

[9] Coppin P,Jonckheere I,Nackaerts K,et al.Digital Change Detection Methods in Ecosystem Monitoring:a Review[J].International Journal of Remote Sensing,2004,25(9):1565 -1569.

Land Cover Change Detection in Coal Mining Area Using BJ-1 Small Satellite Remote Sensing Data

CHEN Yu,DU Pei-jun,TANG Wei-cheng,LIU Si-cong
(Key Laboratory for Land Environment and Disaster Monitoring of SBSM,China University of Mining and Technology,Xuzhou 221116,China)

In order to evaluate the performance of monitoring land cover change in mining areas by BJ-1 small satellite remote sensing data,the authors made an experimental and comparative study of several change detection methods,with the east coal mining district of Xuzhou City in Jiangsu Province as the study area.Direct change detection methods based on image information operation and image information transformation as well as the postclassification comparison method were adopted in the experiment.The results show that the change vector analysis can attain the highest accuracy,followed by the method of image ratio method,image differencing and multi- band principal component analysis.Through the change detection process,several areas with great land cover change were detected,such as Qingshanquan,Hanqiao and Dongzhuang mining areas in northern east mining district and Dahuangshan mining area in southern east mining district.

BJ-1 micro-satellite(BJ-1);Land use/cover;Change detection;Mining area;Change vector analysis

TP 79:X 835

A

1001-070X(2011)03-0146-05

2010-11-02;

2011-03-22

國家自然科學(xué)基金項(xiàng)目(編號:40871195)、中國地質(zhì)調(diào)查局地質(zhì)調(diào)查工作項(xiàng)目(編號:1212011088034)、江蘇省“333工程高層次人才培養(yǎng)計劃”科研項(xiàng)目(編號:2009-32)、江蘇省自然科學(xué)基金項(xiàng)目(編號:BK2010182)及江蘇省大學(xué)生實(shí)踐創(chuàng)新訓(xùn)練計劃項(xiàng)目(編號:S200913)共同資助。

陳 宇(1988-),女,中國礦業(yè)大學(xué)碩士研究生,主要從事資源環(huán)境遙感方面的研究。

杜培軍,Email:dupjrs@gmail.com;dupjrs@126.com。

(責(zé)任編輯:丁 群)

猜你喜歡
檢測方法
“不等式”檢測題
“一元一次不等式”檢測題
“一元一次不等式組”檢測題
“幾何圖形”檢測題
“角”檢測題
學(xué)習(xí)方法
小波變換在PCB缺陷檢測中的應(yīng)用
用對方法才能瘦
Coco薇(2016年2期)2016-03-22 02:42:52
四大方法 教你不再“坐以待病”!
Coco薇(2015年1期)2015-08-13 02:47:34
賺錢方法
主站蜘蛛池模板: 色网站在线免费观看| 啪啪永久免费av| 免费看av在线网站网址| 六月婷婷精品视频在线观看| P尤物久久99国产综合精品| 伊人大杳蕉中文无码| 亚洲成在线观看| 日韩第八页| 99久久精品视香蕉蕉| 色网在线视频| 国产精品网曝门免费视频| 国产精品成| 亚洲国产日韩在线观看| 91国内外精品自在线播放| 亚洲三级片在线看| 成·人免费午夜无码视频在线观看| 国产手机在线观看| 美女无遮挡拍拍拍免费视频| 精品无码一区二区在线观看| 久久精品视频亚洲| 毛片网站观看| 亚洲欧美国产视频| 一本色道久久88| 成人午夜久久| 青青草原国产免费av观看| 久久午夜影院| 亚洲精品你懂的| 国产又色又爽又黄| 日韩精品成人在线| 国产精品分类视频分类一区| 欧美另类图片视频无弹跳第一页| 精品少妇人妻一区二区| 久久亚洲国产最新网站| 亚洲一级毛片在线观播放| 日韩精品毛片| 无码'专区第一页| 亚洲v日韩v欧美在线观看| 激情综合网激情综合| 呦女精品网站| 亚洲精品天堂在线观看| 91精品啪在线观看国产91九色| 114级毛片免费观看| 久久精品国产精品国产一区| 国产精品手机视频一区二区| 97国产精品视频人人做人人爱| 色有码无码视频| 亚洲国产亚洲综合在线尤物| 亚洲网综合| 国产精品夜夜嗨视频免费视频| 88国产经典欧美一区二区三区| 国产一级毛片在线| 国产97视频在线| 成人亚洲天堂| 亚洲精品你懂的| 亚洲欧美不卡视频| 潮喷在线无码白浆| 久久99精品久久久久纯品| 亚洲一区二区成人| 天天视频在线91频| 亚洲最大看欧美片网站地址| 黄色三级毛片网站| 欧美成人日韩| 一区二区日韩国产精久久| 啊嗯不日本网站| 成人第一页| 色综合日本| 天天婬欲婬香婬色婬视频播放| 国产偷国产偷在线高清| 中文字幕 日韩 欧美| 久久精品国产在热久久2019| 国产黄色片在线看| 无码中文字幕精品推荐| 国产麻豆91网在线看| 亚洲三级视频在线观看| 精品一区二区三区水蜜桃| 亚洲第一黄片大全| 天堂中文在线资源| 国产亚洲欧美另类一区二区| 色有码无码视频| 性欧美在线| 高清无码手机在线观看| 99re在线免费视频|