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

2001?2016 年祁連山地區植被覆蓋度對干旱的響應

2021-04-20 07:17:28劉佳茹王建邦
草業科學 2021年3期
關鍵詞:趨勢區域

劉佳茹,趙 軍,王建邦

(1. 西北師范大學地理與環境科學學院,甘肅 蘭州 730070;2. 蘭州大學資源環境學院,甘肅 蘭州 730030)

植被覆蓋度(fractional vegetation cover,FVC)是表示地表植被數量的參數,能夠在一定程度上反映地表植被覆蓋的生長狀況和趨勢[1]。作為評價氣候變化影響與區域生態系統環境狀況的重要指標[2],FVC在生態系統評價、土地荒漠化、水文、大氣污染等領域廣泛應用[3-6],研究FVC 的方法包括逐點調查法[7]、經驗模型[8]、像元二分模型[9]、人工神經網絡[10]等。

干旱是長時間水分匱缺導致的自然災害,具有發生頻率高、持續時間長、影響范圍廣的特點[11]。在全球氣候變暖的條件下,干旱對農業生產的影響最大[12],呈現出常態化趨勢,極端干旱發生頻率和強度不斷增加,由此導致的損失也隨之增加[13]。標準化降水蒸散發指數(standardized precipitation evaportranspiration index,SPEI )常用來表征氣候的干旱程度。計算多尺度SPEI 有利于判斷不同時間段的干旱狀況,能更有效地監測和評估區域干旱程度[14]。

植被在生長過程中需要良好的氣候條件,因此FVC與區域氣候干旱條件和干旱程度密切相關[15],它們之間的關系成為研究熱點。很多學者利用不同類型的遙感數據分析FVC 在不同時間尺度上的變化趨勢[16-18],結合氣象數據研究FVC 對干旱的響應[19]。大部分研究成果中顯示,FVC 與干旱程度呈正相關關系[20]。干旱半干旱地區FVC 對土地荒漠化極為敏感,因此研究FVC 對干旱程度的響應具有重要意義[21]。

祁連山地區是西北干旱半干旱地區的生態保護功能區及河西多條內陸河的發源地,近年來社會經濟的快速發展使人地關系緊張,人類活動對祁連山地區生態環境的影響加劇,加之祁連山大陸性氣候特征,降水量少,氣候干燥,加劇了干旱的影響。21世紀以來,伴隨著全球氣溫升高等氣候變化,祁連山地區氣溫升高降水量增加,極端氣候事件頻發,冷暖夜日數增加、強降水日數增加,導致冰川消融加快、河川徑流量增加,FVC 增加,干旱程度變小,氣候向暖濕化發展[22]。因此,從不同尺度深入研究SPEI 對FVC 的影響對了解祁連山地區FVC 與干旱之間的關系十分必要。本研究利用MODIS NDVI 產品通過像元二分模型計算得到祁連山地區FVC,結合氣象數據計算得到SPEI,分析植被覆蓋時空變化趨勢以及對多尺度干旱的響應關系,以期了解祁連山地區生態環境狀況,為生態環境保護提供依據。

1 研究區概況與研究方法

1.1 研究區概況

圖1 研究區概況圖Figure 1 Overview of the study area

祁連山地區位于青藏高原、蒙新高原和黃土高原交匯地帶,跨越甘肅省西部和青海省東北部邊境(35.5° ? 40° N,95° ? 104° E),由多條西北至東南走向平行山脈和其形成的寬闊谷地組成(圖1)。祁連山地區地理位置特殊,東起烏鞘嶺,南靠茶卡?共和盆地及柴達木盆地,西至當金山口,北臨河西走廊,東西長約800 km,南北寬約300 km[23-24],地勢險峻,由東北向西南逐漸升高。祁連山地區是典型高原大陸性荒漠氣候,年均溫0.6 ℃,年降水量為400~700 mm,由東向西逐漸減少,夏季降水多且濕潤,冬季降水少且干燥,為西北地區多條內陸河提供水源補給。植被類型多種多樣呈垂直地帶性規律[25]。

1.2 研究數據

本研究選取2001?2016 年MODIS NDVI 數據,空間分辨率為250 m,時間分辨率為16 d,來源于美國國家航天局(https://modis.gsfc.nasa.gov/)。采用最大值合成法,排除云、大氣和太陽高度角的干擾,得到月NDVI。選取植被生長季5 月 ? 9 月的NDVI 求平均值得到年NDVI。氣象數據選取祁連山地區及周邊35 個氣象站點2001?2016 年月平均氣溫和降水量數據,數據來源于中國氣象科學數據網(https://cdc.cma.gov.cn)。

1.3 研究方法

1.3.1 計算植被覆蓋度

采用像元二分模型[15,26]計算植被覆蓋度(FVC)。假設一個像元的地表由有植被覆蓋部分與無植被覆蓋部分組成,傳感器觀測到的光譜信息也由這兩個組分因子線性加權合成,各因子的權重是各自的面積在像元中所占的比率[27]。由于不同的影像在計算過程中的權重也不相同,本研究采用在5%置信度范圍內累計頻率為95%對應的NDVI 作為NDVIv,累計頻率為5%對應的NDVI 作為NDVIs值[9]。公式如下:

式中:NDVI 為像元植被覆蓋部分和非植被覆蓋部分的加權平均值,NDVIs為無植被覆蓋像元或純裸土像元的NDVI,NDVIv為全植被覆蓋區域像元NDVI。

1.3.2 計算標準化降水蒸散發指數

評估和監測干旱程度的方法有標準化降水指數(standardized precipitation index, SPI)[28]、標準化降水蒸散發指數(SPEI)等,SPEI 結合了SPI 多尺度和溫度變化影響[29],是研究氣候變暖條件下干旱的有效工具。采用Thornthwaite 方法,公式如下[30-31]:

1)計算潛在蒸散量(potential evapotranspiration,PET):

式中:PET 為潛在蒸散量;Ti為月平均溫度;i 為月;H 為年熱量指數;A 為常數,由H 決定,A = 0.492 +0.179H ? 0.000 077 1H2+ 0.000 000 675H3。

2)計算水平衡(D):

式中:Pi為月降水量;PETi為月潛在蒸散量。

3)采用三參數的log-logistic 概率分布對Di數據序列進行正態化,計算相應的SPEI 指數:

參數α、β、γ 計算公式如下:

式中:Γ 為階乘函數,ω0、ω1、ω2為數據序列Di的概率加權矩。

式中:N 為參與計算的月份數量,ωs為概率加權矩。

對累計概率密度進行標準化:

式中:常數c0= 2.515 517,c2= 0.802 853,c3= 0.010 328,d1= 1.432 788,d2= 0.189 269,d3= 0.001 208。

根據國家氣象局制定的SPEI 干旱等級劃分標準對祁連山區域干旱等級進行劃分,SPEI > ?0.5 為無旱,?1 < SPEI ≤ ?0.5 為輕旱,?1.5 < SPEI ≤ ?1 為中旱,?2 < SPEI ≤ ?1.5 為重旱,SPEI ≤ ?2 為特旱。

1.3.3 Theil-Sen 趨勢分析和Mann-Kendall 檢驗

利用非參數化趨勢度(Sen)[32]分析方法計算植被覆蓋度變化趨勢,并通過Mann-Kendall 統計檢驗方法[33]檢驗變化趨勢顯著性[34-35]。

1) Sen 趨勢分析方法不需要樣本服從一定的分布,且不受異常值干擾,具有一定的規避性。公式如下:

式中:β 為植被覆蓋度變化趨勢;i、j 為時間序列;FVCH(j)、FVCH(i)分別為第j 年和第i 年植被覆蓋度。

2) Mann-Kendall 檢驗是判斷Sen 趨勢顯著性。

設定時間序列:FVCH(i),i = 1,2,···,n,標準化檢驗統計量Z 定義為:

S 統計量定義為:

式中:FVCH(i)和FVCH(j)分別為第i 年和第j 年植被覆蓋度。

在給定顯著性水平α 下,當|Z| > Z(1?α)/2時,表示研究序列在α 水平上存在顯著變化。若通過置信度95%的顯著性檢驗視為顯著,若通過置信度90%的顯著性檢驗則視為弱顯著。

1.3.4 Hurst 指數

通過Hurst 指數[36-37]來分析祁連山地區FVC 未來演變趨勢。R/S(重標極差)分析方法是一種定量描述時間序列信息的長程相關思想的分析方法,其原理是根據時間序列{FVCH(t)}(t = 1,2,3,···,n)定義均值序列:

1)累計離差:

2)極差:

3)標準差:

考慮比值R(τ)/S(τ)≌R/S,若存在R/S∝τH,說明時 間 序 列{FVCH(t)},t = 1,2,3,···,n 存 在Hurst 現象。H 為Hurs 指數,當0.5 < H < 1,表明該時間序列為持續性序列,且未來變化趨勢與過去一致;當H =0.5,表明時間序列為隨機序列,不存在長期相關性;當0 <H < 0.5,表明該時間序列為反持續性,且未來變化趨勢與過去趨勢相反。

1.3.5 FVC 與SPEI 相關性分析

基于像元尺度,進行FVC 與SPEI 的相關性分析,公式如下:

式中:R 為FVC、SPEI 兩變量的相關系數;FVCi為第i 年植被覆蓋度;SPEIi為第i 年年平均標準化降水蒸散發指數。

2 結果分析

2.1 祁連山地區植被覆蓋度時空變化趨勢

2.1.1 祁連山地區植被覆蓋度時間變化分析

2001?2016 年祁連山地區FVC 整體呈增加趨勢,增長速率為每10 年0.019,速度緩慢。2001 年FVC 最低,為0.420 4;2011 年FVC 最大,為0.469 5。2001?2016 年FVC 動 態 變 化 波 動 明 顯,2002?2008年波動幅度較小,2008?2016 年波動較大,16 年間FVC 平均值為0.452 2 (圖2)。2001?2009年與16 年間植被覆蓋FVC 平均值差值為?0.031 7~0.001 2,表明在此期間植被生長狀況雖有所好轉,但較為緩慢,可能與受人類活動影響以及氣溫和降水有關,與2001?2009 年SPEI 指數變化趨勢相近;除2014 年外,2010?2016 年FVC 均高于16 年間平均值,表明FVC 長勢開始好轉,與國家建立自然保護區以及退耕還林還草等相關政策的實施有關。從FVC 季節變化來看(圖2),季節變化明顯,夏季FVC 最高,秋季和春季次之。

圖2 2001?2016 年植被覆蓋度年際和季節變化趨勢Figure 2 The change of inter-annual and seasonal fractional vegetation cover from 2001 to 2016

2.1.2 祁連山地區植被覆蓋度空間變化分析

圖3 2001?2016 年植被覆蓋度空間動態Figure 3 Spatial distribution and variation of fractional vegetation cover from 2001 to 2016

2001?2016 年祁連山地區FVC 整體格局呈現出東南高西北低的分布特征,FVC 區域差異性明顯,達坂山、大通山、青海南山以及青海湖周邊FVC 較高,黨河南山、野馬南山、柴達木山以及宗務隆山FVC 較低(圖3)。由2001?2016 年祁連山地區FVC的變化趨勢可以看出,77.13%的區域植被覆蓋度呈改善趨勢,最大改善速率為0.04·a?1,22.87%的區域植被覆蓋度呈退化趨勢,最大退化速率為0.06·a?1。為檢驗FVC 的變化趨勢,本研究計算了2001?2016年祁連山地區FVC 的Sen 趨勢度,并進行了MK 檢驗。將MK 檢驗結果劃分為顯著變化(Z ≥ 1.96 或Z ≤ ?1.96)、弱顯著變化(?1.96 < Z ≤ ?1.65 或1.65 <Z < 1.96)和無顯著變化(?1.65 < Z ≤ 1.65),Sen 趨勢度結果劃分改善區域(Sen ≥ 0)和退化區域(Sen <0)。將Sen 趨勢分析和MK 檢驗分級結果進行疊加,得到像元尺度上FVC 變化趨勢數據,并將結果分為6 種變化類型(表1),可以看出,得到改善的區域占總面積的77.13%,發生退化的區域占22.87%。

2001?2016 年祁連山地區植被覆蓋度改善的區域遠遠大于植被覆蓋度退化的區域(圖4),顯著改善的區域零星分布在祁連山邊緣區域及高海拔地區,這是由于降水的明顯增加及祁連山自然保護區的設立,同時大力開展退耕還林、封山育林等保護措施,還有高海拔地區人為影響較小,利于植被生長。弱顯著改善區域主要分布在青海南山、黨河南山、疏勒南山及走廊南山。弱顯著退化區域主要分布在各城市區域周邊及青海湖和大通河附近,表明城市擴張及人類活動對植被覆蓋度有影響。顯著退化區域零星分布在祁連山中東段。

表1 植被覆蓋度變化趨勢分類Table 1 stype of trend of fractional vegetation cover change

2.2 祁連山地區植被覆蓋度未來趨勢特征

祁連山地區植被覆蓋度的Hurst 指數均值為0.471 2,其中Hurst 指數小于0.5 的區域面積比例為63.36%,大于0.5 的區域面積比例為36.64%,表明祁連山地區植被覆蓋度的負向特征較強。通過將植被覆蓋度變化趨勢結果和Hurst 指數結果進行疊加得到耦合信息。將結果劃分為持續退化、由退化到改善、持續改善、由改善到退化4 種類型。

持續改善的面積比例為28.10%,主要分布在大雪山、冷龍嶺附近及青海省烏蘭縣西北部的柴達木山;持續退化的面積比例為8.59%,零散分布在各城市周邊區域、宗務隆山、拉脊山附近和烏鞘嶺以南部分區域;由退化到改善的面積比例為14.08%,零星分布在祁連山中東段的門源縣境內、達坂山以及拉脊山部分地區;由改善到退化的面積比例最大,為49.23%,主要分布在祁連山西段魚卡河、柴達木山西南部和疏勒河上游區域以及青海湖以西(圖5)。由改善到退化的區域需要政府部門持續關注并加大力度進行保護。

圖4 2001?2016 年年均植被覆蓋度變化趨勢Figure 4 Trend of inter-annual fractional vegetation coverage from 2001 to 2016

圖5 植被覆蓋度未來變化空間分布Figure 5 Spatial distribution of future changes in fractional vegetation cover

2.3 祁連山地區SPEI 時空變化趨勢

2.3.1 不同時間尺度的祁連山地區SPEI 年際變化分析

圖6 祁連山地區多時間尺度SPEI 長期動態特征Figure 6 Long-term dynamic characteristics of multi-time scale SPEI in Qilian Mountains

根據2001?2016 年祁連山地區不同尺度的SPEI變化趨勢(圖6)可以看出,1、3、6、12 個月尺度的SPEI 均呈現下降趨勢,且下降速度與時間尺度呈正比,干旱化趨勢增強。根據不同時間尺度SPEI 年際變化來看,時間尺度越短,干濕交替周期越頻繁,時間尺度越大,干濕變化越平緩,周期越大。1 和3 個月尺度的SPEI對短期降水和溫度都較為敏感,SPEI 波動周期較短,體現了祁連山地區旱澇無常的特性;6 個月尺度的SPEI 干旱持續半年之久的頻率較大,數值波動也較大。1 個月尺度上,2001?2002 和2013?2014 年干旱頻率較高且干旱持續時間較長,2013 年全年8 個月時間處于輕度干旱狀態;3 個 月 尺 度 上, 2001?2003 年, 2006?2007、 2008?2009 和2012?2013 年干旱頻率也較高,2006 年全年9 個月處于干旱狀態,2011 年全年10 個月處于干旱狀態;6 個月尺度上,2001?2002、2006?2007 和2011?2012 年干旱頻率較高且持續時間較長,2013 年全年處于輕度和中度干旱狀態;12 個月尺度上,2006?2007 和2008?2010 年,2011?2012 和2013?2014 年干旱頻率較高且持續時間較長,2009 年全年處于輕度干旱狀態。由此得出,根據不同時間尺度的 SPEI 反映的干旱程度及周期,可以判斷不同干旱類型的干濕交替狀況。1 個月和3 個月的SPEI 對短期氣溫和降水較敏感,3 個月SPEI 可以反映區域的氣象干旱,如空氣濕度;6 個月SEPI 可以反映區域農業生態干旱,如土壤水分含量情況;12 個月SPEI 可以反映水文干旱,如河流徑流量狀況。

2.3.2 季節、年尺度的祁連山地區SPEI 年際變化

季節尺度上的SPEI 可以反映季節干旱變化。選取35 個氣象站點的2 月、5 月、8 月的SPIE 值分別代表春、夏、秋3 個季節的干旱情況。2001?2016 年祁連山地區16 年間有9 年春季處于干旱狀態,2001 年為中度干旱,其他8 年為輕度干旱;夏季有一半年份處于干旱狀態,2001 年也為中度干旱,其他年份為輕度干旱;秋季也有一半年份處于干旱狀態,均為輕度干旱;2006 年、2011 年和2013 年春、夏、秋3 個季節均處于輕度干旱(圖7)。秋季SPEI年際變化呈現下降趨勢,干旱加強,春季和夏季SPEI呈現上升趨勢,干旱度降低。年尺度干旱變化如圖8所示,2001?2016 年祁連山地區SPEI 呈線性下降趨勢,下降速率為每10 年0.059,干旱趨勢增強,16 年間接近一半年份處于干旱狀態,2013 年為最干旱年為中度干旱(結果通過α = 0.05置信度檢驗)。

圖7 祁連山地區季節尺度SPEI 長期動態特征Figure 7 Long-term seasonal dynamic characteristics of SPEI

圖8 2001?2016 年祁連山地區年尺度SPEI 年際變化Figure 8 Inter-annual variation of SPEI from 2001 to 2016

2.3.3 祁連山地區SPEI 動態分析

2001?2016 年祁連山絕大多數地區處于濕潤狀態,中段和東段部分地區處于輕度干旱狀態(圖9a)。中段、西段西北部和東段西部SPEI呈增加趨勢,干旱化程度降低,濕潤度增加;剩余區域SPEI 呈下降趨勢,趨于干旱,干旱化程度增加(圖9b)??傮w上來看,2001?2016 年祁連山地區干旱化程度呈降低趨勢,說明氣候條件往好的方向發展。

圖9 2001?2016 年祁連山地區標準化降水蒸散發指數(SPEI)空間動態Figure 9 Spatial distribution and variation of standardized precipitation evapotranspiration index (SPEI)in Qilian Mountains from 2001 to 2016

2.4 植被覆蓋度對干旱指數的響應

2.4.1 年尺度SPEI 與植被覆蓋度的相關分析

2001?2016 年祁連山地區年均FVC 和年均SPEI相關性分析表明(圖10a),整個研究區平均相關系數為0.039,正相關區域占55.2%,其中通過P < 0.05和P < 0.01 檢驗的區域面積占比分別為18.8%和2.3%,負相關區域占44.8%,其中通過P < 0.05 和P < 0.01 檢驗的區域面積占比分別為4%和0.8%。由此可以看出,SPEI 與FVC 之間具有一定的相關性,即年尺度上FVC 與SPEI 之間呈正相關關系,SPE 值越大,干旱程度越低,FVC 越高,通過顯著性水平檢驗的區域越大。

2.4.2 季節尺度SPEI 與FVC 的相關分析

圖10 2001?2016 年祁連山地區不同時間尺度植被覆蓋度(FVC)與標準化蒸散發指數(SPEI)的相關關系Figure 10 Correlation between fractional vegetation cover(FVC) and standardized precipitation evapotranspiration(SPEI)index at different times in Qilian Mountains from 2001 to 2016

對2001?2016 年季節尺度上祁連山地區FVC和SPEI 進行相關性分析得出,祁連山地區春季、夏季、秋季的平均相關系數分別為0.056、0.212 和0.058,呈正相關關系。春季(圖10b) FVC 與SPEI 正相關區域占62.7%,其中通過P < 0.05 和P < 0.01 檢驗的區域面積占比分別為2.9%和1.1%,負相關區域占37.3%,其中通過P < 0.05 和P < 0.01 檢驗的區域面積占比分別為1.3%和0.4%;夏季(圖10c) FVC與SPEI 正相關區域占76.6%,其中通過P < 0.05 和P < 0.01 檢驗的區域面積占比分別為12.1%和4.6%,負相關區域占23.4%,其中通過P < 0.05 和P < 0.01 檢驗的區域面積占比分別為1%和0.3%;秋季(圖10d)FVC 與SPEI 正相關區域占59.2%,其中通過P < 0.05和P < 0.01 檢驗的區域面積占比分別為4%和1%,負相關區域占40.8%,其中通過P < 0.05 和P < 0.01檢驗的區域面積占比分別為2.3%和0.7%。季節尺度上夏季祁連山地區FVC 對干旱的響應程度最為明顯,春季和秋季次之,說明FVC 與SPEI 之間顯著性越高,FVC 對干旱的響應也就越明顯。

3 討論與結論

陳京華[38]對2000?2014 年祁連山植被指數NDVI 與降水和氣溫進行了研究,結果顯示夏季NDVI 與降水和氣溫的相關性較高,春季和秋季次之。本研究通過對不同時間尺度的FVC 和SPEI 進行分析,在年尺度上,FVC 與SPEI 呈正相關關系,表明SPEI 越大,FVC 越高;在季節尺度上,夏季相關性最大,其次是春季、秋季,FVC 對干旱的響應程度與年尺度相同,FVC 與SPEI 相關系數越大,表明FVC 對干旱的響應也就越明顯。本研究在季節尺度上得到的結果與以往研究結果基本一致。了解FVC 對干旱的響應程度,對祁連山地區的生態環境保護具有一定的實用性。但是,本研究只是分析了FVC 對SPEI 的響應,沒有考慮不同土地利用類型的FVC 對SPEI 的響應,而且只考慮了氣象因素對FVC 的影響,沒有綜合人為因素以及社會經濟因素對FVC 的影響;還有關于SPEI 的計算,受數據的限制性采用了比較簡單實用的Thornthwaite 方法,在后續的研究中可以采用更為精確的FAO Penman-Monteith 方法提高SPEI 計算精度,以便更好地反映祁連山地區FVC 對干旱的響應程度。

本研究得出以下結論:

1) 2001?2016 年 祁 連 山 地 區FVC 年 均 值 為0.452 2,整體上屬于中覆蓋度,分布特征為東南高西北低,區域差異明顯,整體上呈增加趨勢。季節尺度上,夏季FVC 最高,秋季和春季次之。16 年間,FVC 改善區域遠大于退化區域,FVC 未來由改善到退化的面積為49.23%,主要分布在祁連山西段魚卡河、柴達木山西南部和疏勒河上游區域以及青海湖以西,在后續的發展中要加大力度進行保護。

2) 2001?2016 年祁連山地區不同尺度的SPEI 均呈現出下降趨勢,大部分區域處于濕潤狀態,總體上祁連山地區干旱化程度呈降低趨勢。季節尺度上,祁連山地區秋季SPEI 呈下降趨勢,干旱加強,春季和夏季SPEI 成上升趨勢,干旱度降低。

3)年尺度上,祁連山地區FVC 與SPEI 整體上呈正相關關系。季節尺度上,夏季FVC 與SPEI 相關性最大,春季和秋季次之。說明FVC 與SPEI 之間顯著性越高,FVC 對干旱的響應也就越明顯。

猜你喜歡
趨勢區域
永久基本農田集中區域“禁廢”
今日農業(2021年9期)2021-11-26 07:41:24
趨勢
第一財經(2021年6期)2021-06-10 13:19:08
分割區域
初秋唇妝趨勢
Coco薇(2017年9期)2017-09-07 21:23:49
SPINEXPO?2017春夏流行趨勢
關于四色猜想
分區域
“去編”大趨勢
中國衛生(2015年7期)2015-11-08 11:09:38
基于嚴重區域的多PCC點暫降頻次估計
電測與儀表(2015年5期)2015-04-09 11:30:52
趨勢
汽車科技(2015年1期)2015-02-28 12:14:44
主站蜘蛛池模板: 国产特一级毛片| 免费a级毛片视频| 国禁国产you女视频网站| 成人一级黄色毛片| 成人久久18免费网站| 久草视频精品| 激情在线网| 国产视频入口| 99久久精品国产综合婷婷| 香蕉网久久| 欧美成人在线免费| 久久天天躁狠狠躁夜夜2020一| 亚洲欧美精品在线| 亚洲第一成网站| 亚洲精品无码高潮喷水A| 欧美成人影院亚洲综合图| 波多野结衣一区二区三视频| 秘书高跟黑色丝袜国产91在线 | 免费无码一区二区| 精品久久高清| 成人福利在线免费观看| 国产偷倩视频| 四虎永久免费地址| 婷婷午夜影院| 为你提供最新久久精品久久综合| 成人久久精品一区二区三区| 任我操在线视频| 午夜人性色福利无码视频在线观看| 激情無極限的亚洲一区免费| 亚洲综合在线最大成人| 亚洲男人天堂久久| 国产伦片中文免费观看| 永久免费av网站可以直接看的 | 国产亚洲精久久久久久久91| 亚洲中文字幕97久久精品少妇| 国产91视频免费| 日本a级免费| 再看日本中文字幕在线观看| 天堂亚洲网| 亚洲精品高清视频| 亚洲成肉网| 精品视频一区二区三区在线播| 欧美自慰一级看片免费| 国产97视频在线| 亚洲精品你懂的| 呦视频在线一区二区三区| 久久人妻系列无码一区| 国产乱人伦偷精品视频AAA| 亚洲永久视频| 国产精品视频3p| 找国产毛片看| 亚洲精品无码抽插日韩| 亚洲永久精品ww47国产| 国产真实自在自线免费精品| 亚洲精品国产首次亮相| 狠狠色丁香婷婷综合| 亚洲精品男人天堂| 亚洲视频在线网| 91视频精品| 欧洲在线免费视频| 性做久久久久久久免费看| 男人的天堂久久精品激情| 日韩精品无码免费专网站| 精品国产一区91在线| 亚洲综合色区在线播放2019| 91精品国产丝袜| 中文字幕亚洲另类天堂| 国产乱子伦视频三区| 精品国产91爱| 国产人成午夜免费看| 波多野结衣一二三| 欧美亚洲一区二区三区导航| a级毛片免费播放| 亚洲愉拍一区二区精品| 狼友视频一区二区三区| 亚洲一级无毛片无码在线免费视频 | 无码免费试看| 国产成a人片在线播放| 91在线播放国产| 奇米影视狠狠精品7777| 国产精品人莉莉成在线播放| 91色在线视频|