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

利用不同植被指數估算植被覆蓋度的比較研究

2012-12-27 06:41:02沈潤平楊曉月
自然資源遙感 2012年4期
關鍵詞:模型

徐 爽,沈潤平,楊曉月

(1.南京信息工程大學氣象災害省部共建教育部重點實驗室,南京 210044;2.南京信息工程大學遙感學院,南京 210044)

利用不同植被指數估算植被覆蓋度的比較研究

徐 爽1,2,沈潤平1,楊曉月1,2

(1.南京信息工程大學氣象災害省部共建教育部重點實驗室,南京 210044;2.南京信息工程大學遙感學院,南京 210044)

選用蔬菜地和草地2種植被類型,利用ASD光譜儀實測二者在不同覆蓋度下的光譜響應,分析了歸一化植被指數(NDVI)、差值植被指數(DVI)、比值植被指數(RVI)、修正植被指數(MVI)、修改型土壤調節植被指數(MSAVI)以及全球環境監測植被指數(GEMI)等6種植被指數所用的最佳波段及其組合,進而研究了利用像元二分模型估算植被覆蓋度時的不同植被指數的表現。結果表明,與蔬菜地植被指數相關系數較高的波段組合為620~740 nm譜段和780~900 nm譜段內波段的組合,與草地植被指數相關系數較高的波段組合為620~750 nm譜段和760~900 nm譜段內波段的組合,相關系數均達0.8以上;在高光譜數據構建的植被指數和模擬衛星數據構建的植被指數中,用DVI和MSAVI估算植被覆蓋度,平均總體精度分別達到83.7%和79.5%,與其他4種植被指數相比,這2種指數更適合于利用像元二分模型進行植被覆蓋度的估算。

植被覆蓋度;植被指數;像元二分模型

0 引言

植被作為陸地生態系統的主要成分之一,對氣候變化具有重要的調節作用[1]。植被覆蓋度(vegetation coverage)是指植被在地面的垂直投影面積占統計區總面積的百分比[2],是衡量地表植被狀況的一個重要指標,也是影響土壤侵蝕與水土流失的主要因子,對于區域環境變化和監測研究具有重要意義[3]。利用遙感數據估算植被覆蓋度是測定區域植被覆蓋度的主要手段之一,目前應用的遙感測量方法較多,其中應用最廣泛的是回歸模型法和像元分解模型法。Wittich等[4]和刁兆巖等[5]都曾建立了歸一化差值植被指數(normalized difference vegetation index,NDVI)和植被覆蓋度的回歸模型,并對所選研究區植被覆蓋度進行了估算;但回歸模型法一般只適用于特定區域和特定植被類型,并且需要大量的實測數據建立模型,有一定的局限性,不利于普遍推廣應用。像元分解法則具有理論基礎較好,不依賴于實測數據,不受地域限制,可削弱大氣、土壤背景和植被類型等影響,較易于推廣等優點。其中,像元二分模型(dimidiate pixel model)法得到廣泛的認可,具體應用時,需要首先構建植被指數。

目前,應用NDVI來估算植被覆蓋度的研究較多[6-8],有關其估算精度問題,Jiang 等[9]對 NDVI代入模型所存在的誤差問題進行了公式推導;Small[10]研究發現,在中等植被覆蓋度時利用TM數據計算的NDVI隨著葉面積指數(leaf area index,LAI)的增加逐漸趨于飽和,并且會高估植被密度較高區域周圍零散分布的植被區的覆蓋度,這種現象主要是由于陰影對NDVI的影響造成的。由于不同植被指數的構建有著不同的理論依據和應用價值,利用像元二分模型法開展不同植被指數應用效果的比較尚缺乏系統、深入的研究。

本文基于地表高光譜實測數據,構建了歸一化植被指數(NDVI)、差值植被指數(difference vegetation index,DVI)、比值植被指數(ratio vegetation index,RVI)、修正植被指數(modified vegetation index,MVI)、修改型土壤調節植被指數(modified soiladjusted vegetation index,MSAVI)及全球環境監測植被指數(global environment monitoring index,GEMI)等6種植被指數,采用像元二分模型對植被覆蓋度進行估算,對比分析不同植被指數估算植被覆蓋度的效果,為植被覆蓋度的遙感估算提供理論依據。

1 數據源及其預處理

1.1 實驗數據獲取

選用位于 N 32.203°,E 118.705°處,覆蓋度達100%的青菜地為代表的蔬菜地和草地。光譜測定日期為2011-10-14,天氣狀況良好,晴朗無云。采用ASD Field Spec Pro FRTM野外光譜輻射儀進行采集,光譜采集范圍350~2500 nm,采樣間隔1.4 nm(350~1000 nm)和2 nm(1000~2500 nm),重采樣間隔1 nm。測定光譜時,傳感器探頭垂直向下,視場角25°,距離冠層頂1.6 m;每次采集10條光譜,最后取光譜平均值。所測光譜數據經處理后,選擇適當譜帶構建植被指數,以便利用像元二分模型進行植被覆蓋度估算。

植被覆蓋度通過數碼相機拍攝獲得。為了保證植被覆蓋度測量的準確性,并保證與光譜測量范圍相對應,在拍攝前放置一鐵圈于草地周圍,鐵圈大小盡量與光譜儀的視場角范圍一致;拍攝時相機垂直向下,保證整個鐵圈(樣圓)內的草地在拍攝范圍內(測量光譜時將鐵圈移走)。為研究不同植被覆蓋度下光譜曲線的變化,在利用光譜儀完成樣區內植被覆蓋度100%下的光譜測量后,分7次在該樣圓內進行人工均勻去除一定量植被,調節植被覆蓋度,直至樣圓內所有植被全部去除;每去除一次植被后,進行一次光譜測量和植被覆蓋度拍照。對蔬菜地也做同樣處理,對不同覆蓋度下的蔬菜地進行了光譜采集和覆蓋度測量。

1.2 數據預處理

在Photoshop和ENVI軟件下完成植被覆蓋度的求算。首先在Photoshop軟件下進行樣圓裁切,然后在ENVI軟件下進行監督分類,求算植被所占百分比,得到樣地的植被覆蓋度。草地植被覆蓋度分別為 100%,78%,62%,44%,37%,25%,18% 和0%;蔬菜地的植被覆蓋度分別為100%,85%,74% ,62%,58%,52%,42%,39%,33%,25%,10%和0%。利用ViewSpecPro軟件剔除2個樣地光譜數據中的異常數據,并取光譜數據的平均值得到不同覆蓋度下的植被光譜曲線。

2 植被覆蓋度估算

2.1 不同植被指數構建

首先要確定有效譜段范圍和波段組合,波段選擇方法主要有逐步回歸選擇法、簡單波段自相關選擇法和主成分系數權重選擇法[11-12]等。比較所有波段組合下求算的植被指數與植被覆蓋度的相關系數大小,來確定構建植被指數的較優波段。由于運算量較大,所以通過VBA編程實現運算。另外,根據光譜響應函數模擬衛星傳感器波段的反射率,利用高光譜數據模擬美國陸地衛星Landsat7 TM數據和我國環境減災衛星HJ-1A CCD數據,分別計算NDVI,DVI,RVI,MVI,MSAVI和 GEMI等 6 種植被指數,用以分析代入像元二分模型的適宜性。6種植被指數的計算公式見表1。

表1 植被指數一覽表Tab.1 Vegetation indices

2.2 基于像元二分模型的植被覆蓋度估算

像元二分模型[19]是一種簡單實用的遙感估算模型,它假設一個像元的地表由有植被覆蓋地表與無植被覆蓋地表共同組成,而傳感器觀測到的光譜信息也由這2個組分的因子線性加權合成,各因子的權重是各自的面積在像元中所占比例,如植被的權重可以看作是植被覆蓋度。

根據像元二分模型的原理,通過遙感傳感器所獲信息S可以表達為由綠色植被組分所貢獻的信息Sv和由土壤組分所貢獻的信息Ss組成,即

設像元中有植被覆蓋的面積比例為fc(即該像元的植被覆蓋度),則無植被覆蓋地表(土壤覆蓋)的面積比例為1-fc;設全部由植被所覆蓋的純像元所得的遙感信息為Sveg,則混合像元中的植被組分所貢獻的信息Sv可以表示為Sveg與fc的乘積,即

同理,設全部由土壤覆蓋的純像元所得的遙感信息為Ssoil,則混合像元中的土壤成分所貢獻的信息Ss可以表示為Ssoil與1-fc的乘積,即

將式(2)(3)代入式(1),可得

對式(4)進行變換,可得

式中Ssoil與Sveg是模型的2個參數。本文中Ssoil即為植被覆蓋度為0%時的植被指數信息;Sveg即為植被覆蓋度為100%時的植被指數信息;S則為所求覆蓋度下的植被指數信息。將 NDVI代入式(5)[20],可得

將表示其他5種植被指數的變量(DVI,RVI,MVI,MSAVI和GEMI)也分別代入式(5)進行植被覆蓋度估算。

3 結果與分析

3.1 不同覆蓋度下植被光譜特征

圖1為不同植被覆蓋度下蔬菜地和草地的光譜反射率曲線圖。

圖1 不同植被覆蓋度下植被反射光譜曲線Fig.1 Reflectance spectra of vegetation in different levels of vegetation coverage

從圖1可以看出,隨著植被覆蓋度的增大,可見光譜段內藍光波段反射率有減小趨勢,綠光波段反射率有增大趨勢(但不是很明顯),紅光波段反射率減小較為明顯;近紅外譜段反射率則逐漸增大,并且增加幅度較大。裸土光譜反射率在紅光波段內明顯高于有植被覆蓋時的反射率,在近紅外譜段則明顯低于有植被覆蓋時的反射率。

3.2 高光譜數據構建植被指數估算覆蓋度比較

3.2.1 植被指數最佳波段選擇

因植被光譜反射率隨植被覆蓋度變化的規律在1000 nm后表現不明顯,因此本文選取1000 nm前的數據進行了植被指數構建。由于在可見光和近紅外譜段植被光譜隨植被覆蓋度變化的規律性明顯,因此首先對這2個譜段內所有波段兩兩組合,求算其植被指數,然后計算植被指數與植被覆蓋度的相關系數,依照相關系數的大小來確定植被覆蓋度估算較優的植被指數波段組合。確定的結果是:蔬菜地植被指數相關系數較高的波段組合為620~740 nm譜段和780~900 nm譜段內波段的組合,草地植被指數相關系數較高的波段組合為620~750 nm譜段和760~900 nm譜段內波段的組合,植被指數的相關系數均達0.8以上。本文選取相關系數最高的最佳波段組合構建 NDVI,DVI,RVI,MVI,MSAVI和GEMI(表2),再分別代入像元二分模型估算對應的植被覆蓋度。

表2 不同植被指數最佳波段組合Tab.2 Optimized band combination of different vegetation indices

3.2.2 不同植被指數估算精度比較

將相關系數最高的波段組合構建計算得到的NDVI,DVI,RVI,MVI,MSAVI和 GEMI分別代入像元二分模型,估算植被覆蓋度,并與實測值比較(圖2)。結果表明,對于蔬菜地,用DVI和MSAVI代入像元二分模型估算的覆蓋度與實測值較接近,在圖2(a)中兩者散點較接近于1∶1直線,用MSAVI的估算結果略好于DVI;對于草地(圖2(b)),同樣用DVI和MSAVI的估算結果較接近1∶1直線。而RVI,MVI和NDVI這3種植被指數存在高估植被覆蓋度的現象,對于蔬菜地植被覆蓋度的高估現象則更加嚴重(尤其以用RVI估算植被覆蓋度時高估現象最為嚴重)。具體估算精度見表3。

圖2 植被指數估算覆蓋度與實測覆蓋度對比Fig.2 Comparison between observed and estimated coverage using the six different vegetation indices

表3 高光譜數據構建不同植被指數估算覆蓋度誤差對照表Tab.3 Evaluation of different methods to derive vegetation coverage

由表3可以看出,在蔬菜地用DVI估算植被覆蓋度的總體精度為88.0%,用MSAVI估算的總體精度達到88.3%,估算總體精度均高于另外4種植被指數,且MAE和MRE也均低于其他4種植被指數。可見,用DVI和MSAVI代入像元分解模型效果較好。同樣,在草地情況下也是DVI和MSAVI這2種植被指數估算總體精度較高,平均誤差較低,代入像元二分模型的結果優于另外4種植被指數。

3.3 模擬衛星數據構建植被指數估算植被覆蓋度比較

利用光譜響應函數模擬美國陸地衛星Landsat7 TM數據和我國環境減災衛星HJ-1A CCD波段數據,分別計算以上6種植被指數,代入像元二分模型估算植被覆蓋度,并與實測值比較(圖3—圖4),結果表明,用DVI估算植被覆蓋度的效果最好。

圖3 模擬TM數據構建的植被指數估算覆蓋度與實測覆蓋度對照Fig.3 Comparison between observed and estimated coverage using the six different vegetation indices by simulated TM data

圖4 模擬HJ-1A CCD數據構建的植被指數估算覆蓋度與實測覆蓋度對照Fig.4 Comparison between observed and estimated coverage using the six different vegetation indices by simulated CCD data

從圖3—圖4可知,DVI估算結果最接近1∶1直線,無論是對于草地還是蔬菜地估算的植被覆蓋度均較接近實測值;其次為MSAVI;另外4種植被指數對植被覆蓋度的高估較為嚴重,RVI偏離最大,其他依次是MVI,NDVI和GEMI。具體估算精度見表4。

表4 模擬衛星數據構建不同植被指數估算覆蓋度誤差對照表Tab.4 Evaluation of different methods to derive vegetation coverage

從表4可以看出,無論是對于草地還是蔬菜地,用DVI估算植被覆蓋度的精度都最高,平均總體精度達到80.6%,平均絕對誤差和平均相對誤差均低于其他植被指數,代入像元二分模型的估算結果最好;用MSAVI估算的精度也較高,僅次于DVI,因此模擬衛星數據中也是DVI和MSAVI兩種植被指數更適宜代入像元二分模型。

綜上所述,用6種不同植被指數估算植被覆蓋度,以DVI表現最優。對于蔬菜地與草地2種植被類型,無論是高光譜數據構建的植被指數,還是模擬衛星數據構建的植被指數,用DVI估算的植被覆蓋度總體精度都較高,高光譜數據與模擬衛星數據2種方法估算的植被覆蓋度平均總體精度達83.7%,平均絕對誤差只有4.05%,平均相對誤差為16.25%;其次為MSAVI,平均總體精度達到79.5%,在利用模擬衛星數據情況下總體精度低于DVI,但在利用高光譜數據構建情況下,總體精度略高于DVI;用GEMI估算的植被覆蓋度精度居中,平均總體精度為70.8%;用RVI估算的植被覆蓋度平均總體精度僅為40.09%,平均絕對誤差達24.06%,平均相對誤差則高達50.91%;用MVI和NDVI估算的植被覆蓋度平均總體精度分別為53.27%和53.19%。因此,用RVI,MVI與NDVI代入模型估算精度都較低,代入像元二分模型估算植被覆蓋度的結果不理想。

4 結論

1)本文選用蔬菜地和草地2種植被類型,通過實測不同覆蓋度下植被光譜,構建6種植被指數,估算植被覆蓋度。研究表明,蔬菜地各植被指數相關系數較高的波段組合為620~740 nm譜段和780~900 nm譜段內波段的組合,草地各植被指數相關系數較高的波段組合為620~750 nm譜段和760~900 nm譜段內波段的組合。

2)取相關系數最高波段組合構建 NDVI,DVI,RVI,MVI,MSAVI和 GEMI,代入像元二分模型估算植被覆蓋度,發現DVI和MSAVI表現較好。將模擬TM數據和HJ-1A CCD數據構建的6種植被指數代入像元二分模型,估算植被覆蓋度,DVI和MSAVI同樣表現優于其他植被指數。目前很多研究用NDVI代入像元二分模型進行植被覆蓋度估算,本研究表明NDVI估算精度較低,代入像元二分模型表現不理想。另外RVI,MVI植被指數估算精度同樣較低,代入模型也并不理想。

[1]Qi J,Cabot F,Moran M S,et al.Biophysical Parameter Estimations Using Multidirectional Spectral Measurements[J].Remote Sensing Of Environment,1995,54(1):71 -83.

[2]Gitelson A A,Kaufman Y J,Stark R,et al.Novel Algorithms for Remote Estimation of Vegetation Fraction[J].Remote Sensing of Environment,2002,80(1):76 -87.

[3]黨 青,楊武年.近20年成都市植被覆蓋度動態變化檢測及原因分析[J].國土資源遙感,2011(4):121 -125.Dang Q,Yang W N.Dynamic Supervision and Reason Analysis of Vegetation Coverage Changes of Chengdu in the Past 20 Years[J].Remote Sensing for Land and Resource,2011(4):121 - 125(in Chinese with English Abstract).

[4]Wittich K P,Hansing O.Area-averaged Vegetative Cover Fraction Estimated from Satellite Data[J].International Journal of Biometerology,1995,38(4):209 -215.

[5]刁兆巖,徐立榮,馮朝陽,等.呼倫貝爾沙化草原植被覆蓋度估算光譜模型[J].干旱區資源與環境,2012,26(2):139 -144.Diao Z Y,Xu L R,Feng C Y,et al.The Ground Spectral Model for Estimating Vegetation Coverage on Desertified Grassland,Hulunbeier,Inner Mongolia,China[J].Journal of Arid Land Resources and Environment,2012,26(2):139 -144(in Chinese with English Abstract).

[6]李登科,范建忠,王 娟.陜西省植被覆蓋度變化特征及其成因[J].應用生態學報,2010,21(11):2896 -2903.Li D K,Fan J Z,Wang J.Change Characteristics and Their Causes of Fractional Vegetation Coverage(FVC)in Shaanxi Province[J].Chinese Journal of Applied Ecology,2010,21(11):2896 -2903(in Chinese with English Abstract).

[7]劉 瑞,王世新,周 藝,等.基于遙感技術的縣級區域環境質量評價模型研究[J].中國環境科學,2012,32(1):181 -186.Liu R,Wang S X,Zhou Y,et al.Ecological Environment Condition Evaluation Mode of County Region Based on Remote Sensing Techniques[J].China Environmental Science,2012,32(1):181 - 186(in Chinese with English Abstract).

[8]廖春華,張顯峰,孫 權,等.基于HJ-1高光譜數據的植被覆蓋度估測方法研究[J].遙感信息,2011(5):65-70.Liao C H,Zhang X F,Sun Q,et al.Fractional Vegetation Cover Estimation Using HJ -1 Spaceborne Hyperspectral Data[J].Remote Sensing Information,2011(5):65 -70(in Chinese with English Abstract).

[9]Jiang Z Y,Huete A R,Chen J,et al.Analysis of NDVI and Scaled Difference Vegetation Index Retrievals of Vegetation Fraction[J].Remote Sensing of Environment,2006,101(3):366 -378.

[10]Small C.Estimation of Urban Vegetation Abundance by Spectral Mixture Analysis[J].International Journal of Remote Sensing,2001,22(7):1305 -1334.

[11]Huang Z,Turner B J,Dury S J,et al.Estimating Foliage Nitrogen Concentration from HyMAP Data Using Continuum Removal Analysis[J].Remote Sensing of Environment,2004,93(1/2):18 -29.

[12]Thenkabail P S,Enclona E A,Ashton M S,et al.Accuracy Assessments of Hyperspectral Waveband Performance for Vegetation Analysis Applications[J].Remote Sensing of Environment,2004,91(3/4):354-376.

[13]Rouse J W,Haas R H,Schell J A,et al.Monitoring Vegetation Systems in the Ggreat Plains with ERTS[C]//Third Earth Resources Technology Satellite - 1 Symposium.Washington.DC.NASA,1973,1:309 -317.

[14]Jordan C F.Derivation of Leaf-area Index from Quality of Light on the Forest Floor[J].Ecology,1969,50(4):663 -666.

[15]Pearson R L,Miller D L.Remote Mapping of Standing Crop Biomass for Estimation of the Productivity of the Shortgrass Prairie[C]//Proceedings of the English International Symuposium on Remote Sensing of Environment,1972,2:1375 -1381.

[16]McDaniel K C,Haas R H.Assessing Mesquite- grass Vegetation Condition from Landsat[J].Photogrammetric Engineering and Remote Sensing,1982,48(3):441 -450.

[17]Qi J,Chenbouni A,Huete A R,et al.A Modified Soil Adjusted Vegetation Index[J].Remote Sensing of Environment,1994,48(2):119-126.

[18]Pinty B,Verstraete M M.GEMI:A Non-linear Index to Monitor Global Vegetation from Satellites[J].Vegetation,1992,101(1):15 -20.

[19]Zribi M,Le Hegarat- Mascle S,Taconet O,et al.Derivation of Wild Vegetation Cover Density in Semi-arid Regions:ERS2/SAR Evaluation[J].International Journal of Remote Sensing,2003,24(6):1335 -1352.

[20]李 麗,童立強,李小慧.基于植被覆蓋度的石漠化遙感信息提取方法研究[J].國土資源遙感,2010(2):59-62.Li L,Tong L Q,Li X H.The Remote Sensing Information Extraction Method Based on Vegetation Coverage[J].Remote Sensing for Land and Resources,2010(2):59 -62(in Chinese with English Abstract).

A Comparative Study of Different Vegetation Indices for Estimating Vegetation Coverage Based on the Dimidiate Pixel Model

XU Shuang1,2,SHEN Run - ping1,YANG Xiao - yue1,2
(1.Key Laboratory of Meteorological Disasters,Ministry of Education,Nanjing University of Information Science and Technology,Nanjing 210044,China;2.College of Remote Sensing,Nanjing University of Information Science and Technology,Nanjing 210044,China)

ASD Field Spec Pro FRTM spectroradiometer was used to measure the spectral response of the vegetable and grass at different vegetation coverage levels.The data were applied to calculate six vegetation indices,i.e.,NDVI(normalized difference vegetation index),DVI(difference vegetation index),RVI(ratio vegetation index),MVI(modified vegetation index),MSAVI(modified soil adjusted vegetation index)and GEMI(global environment monitoring index).Then the best combination of spectral bands was analyzed.Furthermore,the performance of different vegetation indices was investigated when they were used to estimate the vegetation coverage by using the dimidiate pixel model.The results show that,for the green vegetable,the best combinations of bands in the spectral region from 620 to 740 nm and from 780 to 900 nm have the best correlation with the vegetation index,whereas for the grass,the best combinations of bands are from 620 to 750 nm and from 760 to 900 nm,with the correlation coefficients of the two cases being all larger than 0.8.The bands of Landsat7 and HJ-1A CCD1 simulated according to the spectral response function were employed to calculate the six vegetation indices.The average overall accuracy for estimating the vegetation fraction by DVI and MSAVI is 83.7%and 79.5%respectively,indicating that they are superior to the other four vegetation indices as the input of vegetation index for the dimidiate pixel model.

vegetation coverage;vegetation index;dimidiate pixel model

10.6046/gtzyyg.2012.04.16

TP 79

A

1001-070X(2012)04-0095-06

2012-02-15;

2012-04-08

國家重點基礎研究發展計劃(973計劃)項目(編號:2010CB950701-1,2005CB121108-6)和江蘇省高校“青藍工程”項目共同資助。

徐 爽(1985-),女,碩士研究生,主要從事遙感建模與應用研究。E-mail:kening_xsh@sina.com。

沈潤平(1963 -),男,教授,博士生導師。E -mail:rpshen@nuist.edu.cn。

(責任編輯:李 瑜)

猜你喜歡
模型
一半模型
一種去中心化的域名服務本地化模型
適用于BDS-3 PPP的隨機模型
提煉模型 突破難點
函數模型及應用
p150Glued在帕金森病模型中的表達及分布
函數模型及應用
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
3D打印中的模型分割與打包
主站蜘蛛池模板: 波多野结衣久久高清免费| 看看一级毛片| 免费人成黄页在线观看国产| 国产啪在线| 久久午夜夜伦鲁鲁片无码免费| 亚洲成年网站在线观看| 亚洲国产亚洲综合在线尤物| 青草国产在线视频| 草草线在成年免费视频2| 亚洲swag精品自拍一区| 亚洲国产精品无码AV| 国产精品午夜福利麻豆| 亚洲国产AV无码综合原创| 四虎AV麻豆| 午夜福利在线观看入口| 欧美一道本| 亚洲第一天堂无码专区| 欧美精品成人一区二区在线观看| 亚洲天堂免费| 国产高清在线丝袜精品一区 | 宅男噜噜噜66国产在线观看| 国产丝袜无码一区二区视频| 亚洲精品成人片在线观看 | 色综合色国产热无码一| 无码一区中文字幕| 久久动漫精品| 日韩A级毛片一区二区三区| 精品少妇三级亚洲| 成人在线观看一区| 日韩高清中文字幕| 久久国产免费观看| 夜夜操国产| 国产欧美精品一区二区 | 国产在线97| 欧美www在线观看| 日本午夜三级| 99久久精品国产综合婷婷| 久久99久久无码毛片一区二区 | 国产视频一二三区| 大陆国产精品视频| 一级成人欧美一区在线观看| www精品久久| 亚洲精品高清视频| 久久综合成人| 国产天天色| 欧美精品在线视频观看| 欧美精品导航| 在线精品欧美日韩| 毛片视频网址| 免费无码AV片在线观看国产| 久久久久免费看成人影片 | 四虎永久免费地址在线网站 | 国产91久久久久久| 性色一区| av色爱 天堂网| 国产精品久久自在自线观看| 亚洲欧美日韩中文字幕在线| h视频在线播放| 日韩精品免费一线在线观看| 久草视频精品| 伦伦影院精品一区| 国产精品大尺度尺度视频| 人妻中文字幕无码久久一区| 亚洲色无码专线精品观看| 国产亚洲欧美另类一区二区| 色欲综合久久中文字幕网| 亚洲中文精品人人永久免费| 午夜久久影院| 亚洲欧美自拍中文| 99re经典视频在线| 一区二区无码在线视频| 亚洲欧美在线综合图区| 亚洲香蕉伊综合在人在线| 第九色区aⅴ天堂久久香| 香蕉视频在线观看www| 国产第一色| 国产精品lululu在线观看| 91精品人妻互换| 日韩AV无码一区| 亚洲AV成人一区二区三区AV| 亚洲AV人人澡人人双人| 国产成人乱码一区二区三区在线|