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

近10年汾河流域植被覆蓋時空變化研究

2017-12-27 08:22:24韓麗君
關鍵詞:區域模型研究

韓麗君,李 晶

太原師范學院 地理科學學院,山西 晉中 030619

1 引言

1.1 選題背景與依據

土地利用/覆蓋(LUCC)是全球變化研究的核心內容,也是區域生態評價的重要指標,在逐漸升溫的環境變化和可持續發展研究中占據了重要地位。作為土地利用/覆蓋最主要的表征,植被是最具有顯示度的覆蓋類型[1-11]。植被覆蓋可作為生物地球化學過程、水文過程、環境變化的輸入參數,在地學領域的研究中具有基礎地位。

1.2 國內外研究現狀

國內外植被覆蓋提取的方法可分為兩大類:一是地表實測,諸如目估法、樣方法、儀器法等[12]。目估法是在實地或面對相片根據經驗判斷植被覆蓋,此方法簡單易行,但精度受人為因素影響較大。樣方法是將研究區域劃分為若干網格,收割網格內的全部植物,根據生物量和植被覆蓋的比例系數推算植被覆蓋,這種方法具有破壞性,且耗時耗力,多用于小區域點狀研究對象的植被覆蓋,難以大面積推廣。儀器法通常使用相機拍攝植被冠層,利用圖像處理軟件提取植被覆蓋,此方法測量結果具有較高精度,估算效率也較高,為目前最常用的實測方法。二是遙感監測法,諸如回歸模型法、植被指數法、像元分解模型法等。回歸模型法是利用回歸分析得到空間外推公式,對影像波段計算出的植被指數或植被覆蓋度進行空間或時間上的外延擴展,可得到更大范圍或更多時期的植被覆蓋,使用最廣的為線性回歸模型。如Larsson H.利用TM影像反演出阿拉伯地區的NDVI指數,并與植被覆蓋度進行了線性回歸,得到了精度較高的統計模型。Anatoly利用反演的NDVI、ARVI、GreenNDVI三種植被指數與小麥的植被覆蓋進行了線性回歸分析,得出了植被指數對于植被覆蓋變化十分敏感的結論,并在一定程度上消除了大氣的影響。從前人所做的研究來看,回歸模型法存在一定的局限性,盡管對于局部區域的植被覆蓋估算在滿足一定效率的基礎上具有較高的精度,但建立的回歸模型只適用于特定的區域或特定的植被,不宜推廣,同時所用數據多是Landsat、SPOT這種中高分辨率的影像,低分辨率的遙感影像由于混合信息較多不宜進行回歸分析,同時,面對研究范圍較大,時間序列較長的情況時會遇到數據獲取的困難。像元分解模型模型法原理是將像元信息分解成多個組分(如土壤、水、植被等),構建分解模型計算覆蓋。植被指數法是建立植被指數與植被覆蓋的關系,計算植被覆蓋。

最重要的是,作為一種緩變自然現象,植被覆蓋能夠被遙感影像所記錄。同時,由于環境條件是植被覆蓋時空分布和變化的直接驅動力,植被對環境指標具有降低變化幅度等反饋的同時也能表征環境因子的變遷[13]。本研究通過植被覆蓋的時空分布模式來定性與定量評價汾河流域生態環境,為區域的可持續發展提供科學依據。

2 模型與方法

2.1 估算模型

目前已經發展了許多利用遙感影像反演植被覆蓋的模型,較為實用且簡單的方法是利用植被指數估算植被覆蓋。植被指數是遙感影像反演地表植被生長狀況的一個指標[14],常用的植被指數是歸一化植被指數NDVI(Normalized Difference Vegetation Index)[15],它利用了植被在紅光波段的高吸收和近紅外波段的高反射這個特征,計算公式為:

式中ρnir為近紅外波段的反射率ρred為紅光波段的反射率。NDVI值在[-1,1]之間[15],正值表示有植被區,且植被狀況越好值越高,負值或零表明像元為云、水、裸地等無植被區。

從植被指數到植被覆蓋度的計算采用像元二分模型,它假設一個像元僅由兩種地物(如植被區和無植被區)組成,影像記錄的光譜信息是兩種地物信息的線性加權和,各組分權重是其面積所占比例。

根據上述像元二分模型原理,像元信息P可分解為由植被貢獻的信息Pveg和土壤、裸巖等無植被區貢獻的信息Pbare[16]。

對于由植被區和無植被區構成的混合像元,植被區所占比例即是此像元植被覆蓋度Cveg,與之對應的是非植被區覆蓋:

式中,Pveg為完全有植被覆蓋的純像元所反映的光譜信息,Cveg*Pveg表示混合像元情況下植被所貢獻的光譜信息;同理,Pbare為完全無植被覆蓋的純像元所反映的光譜信息,(1-Cveg)*Pbare表示混合像元情況下無植被區所貢獻的光譜信息。

將上述公式進行方程組求解,可以得到:

像元NDVI值可視為有植被覆蓋所貢獻的NDVIveg和無植被覆蓋所貢獻的NDVIbare的組合,以此作為像元信息P代入公式,可以得到利用歸一化植被指數計算植被覆蓋度的公式:

式中,NDVIbare為完全無植被覆蓋像元的NDVI值,NDVIveg為完全由植被覆蓋像元的NDVI值。

利用這個模型計算植被覆蓋度的關鍵是計算NDVIsoil和NDVIveg。這里有兩種假設:

(1)當區域內可以近似取VFCmax=100%,VFCmin=0%。

公式(1)可變為:VFC=(NDVI-NDVImin)/(NDVImax-NDVImin)

NDVImax和NDVImin分別為區域內最大和最小的NDVI值。由于不可避免存在噪聲,NDVImax和NDVImin一般取一定置信度范圍內的最大值與最小值,置信度的取值主要根據影像實際情況來定。

(2)當區域內不能近似取VFCmax=100%,VFCmin=0%

當有實測數據的情況下,取實測數據中的植被覆蓋度的最大值和最小值作為VFCmax和VFCmin,這兩個實測數據對應影像的NDVI作為NDVImax和NDVImin。

當沒有實測數據的情況下,取一定置信度范圍內的NDVImax和NDVImin。VFCmax和VFCmin根據經驗估算[17]。

一般情況下,可根據影像區域內NDVI分布,以0.5%的置信度的上下限閾值分別代表NDVIveg和NDVIbare。

基于像元二分模型的植被覆蓋計算模型理論簡單,所需參數容易獲取,在實際應用中使用比較廣泛。

2.2 分析方法

本研究在時間和空間兩個維度對研究區植被覆蓋發展變化及其與環境因子的關系進行分析:一是利用空間平均法分析植被覆蓋時間變化及其與環境因子的關系,以發掘植被覆蓋的時間分布模式;二是利用逐像元法分析植被覆蓋空間變化及其與環境因子的關系,以發掘植被覆蓋的空間分布模式。植被與影響因子的關系通過相關系數來研究。考慮到數據獲取等因素的限制,本文計算了植被覆蓋與氣溫、降水量的相關性。對于一個像元來說,10年際的16 d步長溫度(或降水)會形成一個時間序列,同一時期的植被覆蓋也會形成一個同等長度的時間序列,相關系數即可以得到每個像元植被覆蓋與氣象因子的相關性,公式為:。

式中i為16 d植被覆蓋均值,j為16 d氣溫均值或降水量,相關系數為兩個變量集的協方差除以它們的標準差乘積。

為定量化覆蓋程度,等間距將VFC劃分為5級反映覆蓋程度,具體劃分如表1。

表1 植被覆蓋分級標準Table 1 Classification criteria for VFC

3 植被覆蓋信息提取

3.1 研究區概況

汾河流域位于黃河中游,處于 34°N~38 °N、110°E~113°E 之間,行政區域包括忻州、太原、晉中、呂梁、臨汾、運城6個地級市,40個縣。命名河流汾河發源于山西忻州寧武縣管涔山,南向跨越太原后在新絳縣折向西,在河津市匯入黃河,全長709.9 km,流域面積39826 km2,是黃河的第二大支流。河源至太原為上游地區,多為山地、黃土丘陵,此段水土流失嚴重,是洪水、泥沙的來源地;太原至洪洞縣為中游地區,進入晉中平原后坡降變緩,泥沙淤積嚴重。作為山西的母親河,汾河在山西的自然、歷史、文化地位舉足輕重。汾河流域是山西人口密集、工農業集中的區域,人口占全省比例為39%,工業產值占全省46%,耕地面積占全省30%,農業產值占全省64%。歷史時期的統計資料表明,中游地區的洪澇干旱頻發。研究區如圖1所示,右下角為研究區在山西的位置。

圖1 研究區示意圖Fig.1 Location of the studied areas

3.2 數據與處理

本研究利用VFC指數距平來估算植被覆蓋,必須的數據來源于MODIS的NDVI產品。NDVI數據為美國航空航天局TERRA衛星搭載的中分辨率成像光譜儀(MODIS)的多時相合成產品。衛星每天都能獲取同一區域的原始影像,為減少或避免氣象原因導致的無數值區域,其產品數據采用一定的算法進行多天的合成。NDVI16 d合成產品采用了平均值算法,為展開研究工作,所用遙感影像來源于NASA LAADS網站(http://ladsweb.nascom.nasa.gov/data),為2004~2013年MODIS產品數據MOD13A2,為1 km分辨率植被指數16 d合成產品。下載了圖幅號為h26v5和h27v5的MODIS產品,數據為HDF格式,作為產品,已經進行了幾何糾正和輻射糾正,通過鑲嵌,投影,利用研究區矢量邊界裁剪出研究所需的影像。均值合成后的NDVI數據在內部存在無數值區域,利用ArcGIS的Nibble功能修復壞點,形成內部連續的影像,統一為img格式。數據統一為六度分帶的Transverse Mercator平面坐標投影,在這種區域尺度下,能最大程度減少區域內面積或角度屬性的變形。所用參數如表2。

表2 Transverse Mercator投影參數Table 2 Projection parameters of Transverse Mercator

氣象資料來源于國家氣象中心(http://cdc.cma.gov.cn/home.do),利用汾河流域及周邊地區22個氣象站2004~2013年降水、氣溫日值序列,以16 d累積量或平均值,使用反距離加權方法進行插值,分別得到230個時間序列的降水和氣溫面域化數據。數字高程模型DEM來源于ASTER GDEM(http://gdem.ersdac.jspacesystems.or.jp)。

4 植被覆蓋分析

4.1 植被覆蓋時空演化分析

利用2004~2013年16 d合成的MODISNDVI數據,根據上述估算模型,在 ArcGIS、Erdas、Envi、Matlab的協同工作下,結合流域邊界矢量數據,實現汾河流域植被覆蓋的遙感反演。圖2顯示了2004~2013年16 d步長的植被覆蓋時空變化。

圖2 2004~2013年16 d步長的植被覆蓋時空變化Fig.2 Temporal and spatial variations of VFC for 16 d step length from 2004 to 2013

每年從年積日(Day of Year,DOY)第1 d到最后1 d依據16 d步長共有23個時相。綜合10年際的時空分布,從DOY001到DOY081(對應1~3月)這一時期的植被覆蓋變化不大,大部分區域的VFC值在0.4以下,只有少部分區域的植被覆蓋比較可觀,太原西部呂梁山區相對有一個聚集區VFC值在0.4以上。這一時期整體是兩邊的山區(呂梁山脈,太行山脈)VFC值較高,中間的河谷盆地值較低。到后期,VFC值在南部有一定的升高,臨汾盆地相對明顯。DOY097到DOY113這兩個時相相對穩定,從DOY129開始,VFC值在流域邊緣區域有明顯的升高,太原盆地VFC值升高的幅度相對較低,在DOY129開始會形成一個明顯的低值區域,并一直持續到DOY161。這一時期,流域兩邊的山區VFC值增加明顯。在DOY177這一時相,流域北部,特別是太原盆地VFC值增加明顯,這一時相北部的VFC值明顯超過流域南部。從DOY193開始,流域VFC值開始大規模攀升,流域內值基本在0.6以上,兩邊山區的VFC值多在0.8以上,這一狀態一直持續到DOY241和DOY257。這一時期過后,VFC值開始有較為明顯的降低,在太原盆地和臨汾盆地降速尤為顯著。從DOY289開始,流域內VFC值降低到0.4以下,太原盆地的值最低。

4.2 植被覆蓋與氣象因子的關系

利用相關性來分析植被覆蓋與氣象因子的關系,兩種氣象因子對植被覆蓋的影響力和控制力因地理位置的不同而異。圖3顯示了兩種氣象因子與每個像元16 d植被覆蓋時間序列的相關性在空間上的分布。左邊為降水與植被覆蓋的空間相關性,降水與植被覆蓋的相關性呈現出北部高、南部低的特點,整體上相關性偏低;右邊為氣溫與植被覆蓋的空間相關性,氣溫與植被覆蓋的相關性呈現出四周高,中間低的分布,在介休縣附近有一個相關性的低值區域,在太原北部偏西的區域相關性較高,普遍在0.8以上。

圖3 植被覆蓋與氣象因子的關系Fig.3 The relationship between vegetation coverage and meteorological factors

圖4 VFC、溫度、降水的時間序列分布Fig.4 Time series distribution of VFC,temperature and precipitation

為發掘流域VFC的時間分布模式,將流域內VFC、氣溫、降水的均值進行聯合對比分析。如圖4所示,從上到下依次為VFC、氣溫、降水,紅色曲線為最大值,綠色曲線為標準差,藍色曲線為均值,藍綠色曲線為最小值。坐標橫軸為時間序列,每年23個時相,10年一共230個時間序列,橫坐標為各指標值(為清晰顯示做了偏置,觀察變化趨勢即可)。通過統計圖與統計數字,發現VFC值與氣溫、降水在時間分布上有較好的一致性,但在時間上存在一個時相的步長差。VFC峰值出現在DOY209附近,氣溫和降水的峰值出現在DOY225附近。

4.3 植被覆蓋與地形因子的關系

地形因子能影響植被的光熱條件和水分狀況,進而對植被覆蓋產生作用。對區域內數字高程模型進行坡度分析和坡向分析,通過閾值分類,區分出坡度分類圖和坡向分類圖。

海拔分級標準按0~500 m、500~1000 m、1000~1500 m、1500~2000 m、>-2000 m共5個級別。統計得到不同海拔的植被覆蓋均值如圖5所示。

圖5 植被覆蓋的海拔分布Fig.5 Statistics for vegetation cover altitudinal distribution

從圖5可以看出,在海拔1500~2000 m這個區間的VFC值最高,在這個海拔之下,植被的VFC值隨海拔的升高而升高,但幅度較小,植被分布隨海拔高度的變化不大。在海拔大于等于2000 m這個區間的VFC值有明顯的降低,說明這個海拔對對植被的影響較為顯著,限制了植被的生長。通常來說,低海拔地區人類活動較為頻繁,高海拔地區人類活動較少,低海拔區域的VFC值不高可能與此有關。

坡度分級標準為 0°~3°、3°~6°、6°~9°、9°~12°、>=12°共5個類別,將坡度分類圖與植被覆蓋進行疊置分析,計算每個坡度范圍內的植被覆蓋均值,得到植被覆蓋的坡度分布,如圖6。

平均下來,0~3度范圍內植被覆蓋均值為0.3022,3~6度范圍內植被覆蓋均值為0.2199,6~9度范圍內植被均值為0.1787,9~12度范圍內植被覆蓋均值為0.1533,大于等于12度地區植被覆蓋均值為0.1431。坡度越大,VFC值顯著降低,說明坡度是影響植被生長的一個主要因子。其原因在于坡度的大小與水分密切相關,坡度增加會導致水分流失嚴重,致使植被生長狀況不佳。

對數字高程模型進行坡向分級,分級標準為135~225(陽坡)、90~135和225~270(半陽坡)、315~45(陰坡)、45~90和270~315(半陰坡)、無坡向共5類,如圖7所示。將坡向圖與植被覆蓋進行疊置分析。計算每個坡向范圍內的植被覆蓋均值,得到植被覆蓋的坡向分布,如圖7。

圖6 植被覆蓋的坡度分布Fig.6 Statistics for vegetation cover slope distribution

圖7 植被覆蓋的坡向分布Fig.7 Statistics for vegetation cover slope distribution

從圖7中可以看出,陰坡植被覆蓋均值最高,半陰坡植被覆蓋次之,半陽坡植被覆蓋最低。結合上述坡度對于VFC值的影響,說明此區域降水是限制植被生長的一個重要因子,陰坡光照少,有利于水分的涵養,植被茂盛,半陽坡和陽坡因陽光照射導致水分蒸發較多,水分流失嚴重,植被反而不是很理想。

5 總結與展望

5.1 總結

本研究利用MODIS NDVI數據,基于像元二分模型對汾河流域2004~2013年10年際16 d步長的植被覆蓋度進行了估算,對植被覆蓋的時空演變過程進行了全程描述,并結合氣象因子和地形因子給出了分析。研究得出了以下幾點結論:

(1)植被覆蓋在2004~2013年較為穩定,從DOY001到DOY081變化不大,基本為0.4以下,從DOY129開始,VFC值有明顯的升高,在流域邊緣的山區尤為顯著,從DOY193開始,流域VFC值快速升高,流域內VFC值多在0.6以上,DOY257過后,VFC值回落到0.4以下。總體來說,中間的太原盆地和臨汾盆地VFC值較低,兩邊的呂梁山區和太行山區VFC值較高;

(2)對植被覆蓋與氣象因子的空間分布關系與時間分布關系研究表明:空間上,植被覆蓋與降水的相關性北部高、南部低,植被覆蓋與氣溫的相關性邊緣高、中間低;時間上,植被覆蓋與降水、氣溫具有一致的波形,都是夏季高、冬季低,但在相位上存在差異,植被覆蓋的峰值相對于降水與氣溫的峰值提前一個步長(16 d);

(3)通過統計不同海拔、坡度、坡向的VFC值,可以得出流域植被集中分布在500 m以下,坡度小于3°的陰坡地區,符合實際情況。從結果來看,低海拔區域植被覆蓋較為穩定,海拔大于等于2000 m后植被覆蓋有一個明顯的降低;

(4)植被覆蓋是一種復雜地表過程,其成因受多方面因素影響,發展變化過程難以完全由自然因子來控制。

本文作為階段性的成果,特別是在16 d這樣一個高時間分辨率的步長下得到了流域植被覆蓋,相對以往不存在時間分辨率或以季或年為時間分辨率的植被覆蓋計算來說不失為一個有益的嘗試,同時,本研究探索了植被覆蓋變化的自然因子和驅動力,可為植被保持和恢復提供科學依據。

5.2 不足與展望

(1)本文統計了植被覆蓋的海拔、坡度、坡向分布,但影響植被生長的自然因子遠不僅此,地形對于水分的重分配而導致的植被生長不均影響顯著。人為的因子和驅動力因數據獲取困難而未進行相關研究,從本研究的初步結果來看,很有可能是比較重要的一環。可以將有關植被、人口密度、交通網絡、水系分布等圖層進行疊加分析得出人為影響,這一切期待更為開放的數據政策;

(2)本研究使用的VFC基于植被指數,MODIS影像獲取的植被指數在高覆蓋植被區的生長季節具備較高可靠性,在稀疏植被區(覆蓋度小于15%)的非生長季節或無植被區可信度較低,這些區域或時期的植被覆蓋估算需進一步的研究;

(3)植被覆蓋對氣象因子的反饋作用有待研究,如降水、溫度適宜時,植被生物量在生長季會保持上升趨勢,遙感反演的植被指數同步上升,小區域內的氣溫也應保持在一定范圍,如遇干旱,植被受水分脅迫作用,一方面生長受限,植被指數降低,另一方面植被無足夠水分供給蒸騰作用時,部分氣孔關閉,致使小區域內溫度升高;

(4)本文從宏觀尺度研究了流域內的植被覆蓋,但流域內人口空間分布極度不均,工業農業均有可觀規模的分布,特別是農業及其附屬的水利設施對植被覆蓋的影響不可忽視,農作物的收割不可避免帶來VFC值的突降,水利設施的灌溉條件又能極大影響植被覆蓋,這兩者都屬于人類活動,涉及到人與自然的耦合關系。本文僅僅從自然屬性出發,研究了植被覆蓋與降水、溫度的影響,忽視人類活動的影響會致使研究結果存在一定偏差,尤其是近10年城市化提速的情況下,建成區擴展對植被覆蓋的影響應該是拉低作用,這在未來的研究中應予以考慮;

本研究是利用VFC對過去10年植被覆蓋情況進行了研究,給出的只是相對的覆蓋等級,但協同降水、氣溫、地形進行的分析,對既定事實的描述有助于對未來植被覆蓋的恢復,在后續的研究中期待更長時間序列、更多因子的加入。

[1]鐘莉娜,趙文武.基于NDVI的黃土高原植被覆蓋變化特征分析中國水土保持科學,2013,11(5):57-62

[2]王兮之,陸 龍,梁釗雄.青海湟水流域植被覆蓋度時空變化分析[J].水土保持研究,2010,17(6):172-176

[3]樸世龍,方精云.最近18年來中國植被覆蓋的動態變化[J].第四紀研究,2001,21(4):294-302

[4]汪明霞,王衛東.植被覆蓋度的提取方法研究綜述[J].黃河水利職業技術學院學報,2013,25(2):23-26

[5]李 琳,譚炳香,馮秀蘭.北京郊區植被覆蓋變化動態遙感監測—以懷柔區為例[J].農業網絡信息,2008(6):38-41

[6]姜 燁.孫建國,李 慶.基于像元二分模型的植被覆蓋度遙感信息提取[J].科技信息,2013(10):168-170

[7]梁四海,陳 江,金曉媚,等.近21年青藏高原植被覆蓋變化規律[J].地球科學進展,2007,22(1):33-40

[8]陳國榮,劉 健,施聰智,等.基于GIS和遙感技術的植被覆蓋度反演模型[J].北華大學學報:自然科學版,2010,11(3):251-254

[9]秦 偉,朱清科,張學霞,等.植被覆蓋度及其測算方法研究進展[J].西北農林科技大學學報:自然科學版,2006,34(9):163-170

[10]周兆葉,儲少林,王志偉,等.基于NDVI的植被覆蓋度的變化分析—以甘肅省張掖市甘州區為例[J].草葉科學,2008,25(12):23-29

[11]張學珍,朱金峰.1982~2006年中國東部植被覆蓋度的變化[J].氣候與環境研究,2013,18(3):365-374

[12]章文波,符素華,劉寶元.目估法測量植被覆蓋度的精度分析[J].北京師范大學學報:自然科學版,2001,37(3):402-408

[13]龍曉閩,周忠發,張 會,等.基于NDVI像元二分模型植被覆蓋度反演喀斯特石漠化研究——以貴州畢節鴨池示范區為例[J].安徽農業科學,2010(8):4184-4186

[14]李曉香,陳媛麗.基于NDVI的鶴壁市植被覆蓋動態分析研究[J].河南科技2014(7):188

[15]詹起林.基于森林植被指數的快速生態評價[D].上海:上海師范大學,2012

[16]聶 勇,范建容,賀秀斌,等.水土流失遙感調查中植被信息提取與評價指標討論[J].水土保持通報,2007,37(4):10-14

[17]蔣友嚴,韓 濤,王有恒,等.石羊河調水10a來民勤綠洲生態脆弱性變化[J].干旱區研究,2014,31(1):157-162

猜你喜歡
區域模型研究
一半模型
FMS與YBT相關性的實證研究
遼代千人邑研究述論
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
視錯覺在平面設計中的應用與研究
科技傳播(2019年22期)2020-01-14 03:06:54
EMA伺服控制系統研究
3D打印中的模型分割與打包
關于四色猜想
分區域
主站蜘蛛池模板: 亚洲永久色| 成人亚洲视频| 中文字幕 91| 四虎影视国产精品| 国产精品原创不卡在线| 国产手机在线小视频免费观看| 99久久国产综合精品2020| 91年精品国产福利线观看久久| 婷婷综合缴情亚洲五月伊| 亚洲一区色| 国产丝袜啪啪| 久久国语对白| 国产精品亚洲欧美日韩久久| 色欲综合久久中文字幕网| 亚洲人妖在线| 老司机午夜精品网站在线观看| 日韩欧美国产精品| 欧美在线黄| 国产无码网站在线观看| 无码高清专区| 午夜少妇精品视频小电影| 精品少妇人妻一区二区| 激情乱人伦| 在线看片中文字幕| 成人亚洲国产| 精品免费在线视频| 熟女视频91| 超清无码熟妇人妻AV在线绿巨人| 久久天天躁夜夜躁狠狠| 中国美女**毛片录像在线| 99精品视频在线观看免费播放| 婷婷六月综合网| 成人年鲁鲁在线观看视频| 亚洲欧美天堂网| 国产精品亚洲日韩AⅤ在线观看| 成人精品视频一区二区在线| 日韩无码黄色网站| 日韩无码视频网站| 在线免费亚洲无码视频| 亚洲无码视频喷水| 超级碰免费视频91| 欧美三級片黃色三級片黃色1| 亚洲第一页在线观看| 国产玖玖玖精品视频| 午夜福利无码一区二区| 国产肉感大码AV无码| 亚洲一区网站| 欧美日韩中文字幕二区三区| 精品国产91爱| 国产日韩欧美成人| 久久久91人妻无码精品蜜桃HD | 国产人成乱码视频免费观看 | 中日无码在线观看| 午夜福利视频一区| 青青草原国产精品啪啪视频| 亚洲国产成人综合精品2020| 欧美在线观看不卡| 亚洲a级在线观看| 四虎永久免费地址在线网站| 久久久久国色AV免费观看性色| 免费福利视频网站| 久久精品电影| 欧美成人A视频| 超清无码熟妇人妻AV在线绿巨人| 91视频精品| 97视频免费在线观看| 天堂久久久久久中文字幕| a免费毛片在线播放| 色成人亚洲| 18禁高潮出水呻吟娇喘蜜芽| 精品视频一区二区观看| 中文字幕亚洲无线码一区女同| 欧类av怡春院| 97免费在线观看视频| A级毛片无码久久精品免费| 国产va免费精品观看| 国内精品小视频福利网址| 欧美国产在线看| 四虎AV麻豆| 大香网伊人久久综合网2020| 视频二区亚洲精品| 91精品啪在线观看国产|