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

基于MIKESHE的流域表層土壤含水量時空變化特征分析

2020-04-13 03:07:42雷凱文盧宏瑋
安徽農業科學 2020年6期

雷凱文 盧宏瑋

摘要?基于氣象站點實測資料、MODIS植被產品數據集LAI(leaf area index)及世界土壤數據庫(HWSD)土壤數據,結合氣象插值軟件ANUSPLIN、土壤水分特征曲線擬合軟件RETC、空間信息處理軟件ArcGIS分別對氣象、土壤、植被數據進行計算,使用MIKESHE分布式水文模型對拉薩河流域內2008—2017年3、7月及典型年(豐水年、中水年、枯水年)生長季(3—10月)內的表層1 m土壤含水量進行模擬。結果表明,在空間分布上,拉薩河流域表層土壤含水量自西南向東北呈遞減趨勢,流域中部拉薩河中上游的高山地區土壤含水量相對較低,相較于中下游河谷及西部盆地地區,中部高山地區土壤水分變化更為劇烈。在時間分布上,典型年生長季內不同月份土壤含水量變化量在33%~37%,總體而言,土壤含水量年內變化大于年際變化。拉薩河谷平原、澎波盆地、流域西部念青唐古拉山脈以東的洪積寬谷、羊八井盆地及流域東北部麥地卡濕地植被覆蓋度和土壤水分含量較高、變化幅度小,農業適宜性較強。

關鍵詞?土壤含水量,MODIS,MIKESHE,時空演變,拉薩河流域

中圖分類號?P334+.92文獻標識碼?A

文章編號?0517-6611(2020)06-0050-05

Abstract?Based on meteorological site data,MODIS vegetation product dataset LAI (leaf area index) and soil data of the world soil database (HWSD),soil moisture content in the March and July during 2008-1017 in Lhasa River Basin was simulated and analyzed by the MIKESHE model,a distributed hydrology model.The results showed that in the spatial distribution,the surface soil water content across the Lhasa River Basin gradually decreased from the southwest to northeast,and the water content in the alpine region on the east and west sides of the middle stream was lower.Compared with the valley locate in the middlelower stream,the change of soil moisture in the central alpine area was more obvious.In terms of time distribution,the soil water content of the typical annual growth season changes in different months,varied between 33% and 37%,comparatively,the actual changes in soil moisture content were more obvious during the year.The basins and valley areas in the study area had higher vegetation coverage and higher soil moisture content,i.g.the Libo Basin and the Lhasa Valley Plain locate at middlelower stream of Lhasa River,Hongji Kuangu and Yangbajing Basins in the east of the Nyainqentanglha Mountains in the western part of the basin,and the Medica wetland in the northeastern part of the basin.The resulting agricultural suitability is stronger in these areas.

Key words?Soil moisture,MODIS,MIKESHE,Temporal and spatial variation,Lhasa River Basin

土壤含水量是氣候變化、水文循環、生態平衡過程中的一項重要參數,了解土壤含水量時空分布特征對區域旱澇預報[1-2]、農牧業生產決策[3-4]與水資源調控[5-6]具有重要意義。目前常用的研究土壤含水量分布變化的方法主要包括實地測量法、遙感反演法、分布式水文模型模擬3種。劉雅莉等[7]通過在4種不同地形下設置采樣點并進行長期監測與采樣,分析了微尺度下土壤水分的時空分布格局變化特征及其驅動因子,楊凱悅等[8]利用EM50數據采集器連續記錄監測土壤體積含水量,對高寒沙區人工林在不同生長階段的土壤水分分布特征進行了研究。傳統的實地測量方法雖然可獲得精度更高的數據,但同時也存在數據獲取周期長、人工成本高的缺點,且站點網絡密度的限制使得傳統的站點監測方法在區域尺度上的應用被限制,難以用于大范圍區域土壤含水量的監測[9]。隨著觀測技術與手段的不斷革新發展,基于衛星遙感的土壤含水量反演技術獲得了很大發展,產生了以微波遙感法、熱慣量法、溫度-植被干旱指數法等土壤含水量的反演方法,Kolassa等[10]以ISMN站點監測數據為基礎,對常用的土壤含水量反演產品SMOS、SMAP、MetopA及ESA CCI在不同區域及氣候條件下數據集的精度進行了評價。

使用遙感技術反演估算土壤含水量具有觀測快速、覆蓋區域廣的優點,但同時觀測結果受云雨與復雜地表覆蓋類型的影響,通常難以獲得高精度、連續性強的土壤含水量[11]。由于遙感反演技術在高空間分辨率土壤含水率模擬方面的局限性,在利用衛星遙感獲取的氣象、土地利用類型、LAI及DEM有效數據基礎上,使用分布式水文模型模擬目標區域的土壤含水量在缺少觀測資料的地區具有廣泛的應用前景。如Khan等[12]利用衛星獲取非洲維多利亞湖某流域的氣象、土地利用類型等遙感數據,以此為基礎使用水文模型CREST(coupled routing and excess storage)對該流域的水文過程進行了模擬,Mason等[13]基于TerraSAR-X合成孔徑雷達(SAR)遙感數據,估測英國Tewkesbury城區的洪水淹沒區,研究認為遙感獲取的洪水范圍能夠為水文模型的驗證提供較為可靠的信息。遙感技術為水文模型提供了大量基礎的地表水文信息,為缺乏實測資料區域的水文模擬提供了數據基礎,提高了定量描述流域地貌、水文過程的可能性。基于遙感數據的水文模擬在未來具有廣泛的應用前景[14-15]。

該研究利用氣象站點監測數據、HWSD土壤數據庫數據及MODIS遙感數據,結合分布式水文模型MIKESHE對拉薩河流域表層土壤含水量進行模擬,獲得2008—2017年流域表層土壤含水量,并對模擬結果與流域31個采樣點實測數據進行皮爾遜相關分析,選取各年枯水期與豐水期的空間分布結果,分析拉薩河流域表層土壤水分在不同情景下各年的空間分布狀況及其變化。

1?資料與方法

1.1?研究區概況

拉薩河流域地處雅魯藏布江流域中游左岸,面積33 823 km2,約占雅魯藏布江流域面積的13.5%。干流拉薩河全長約551 km,主要支流包括堆龍曲、墨竹馬曲、麥曲、拉曲、桑曲、和玉年曲等。拉薩河流域氣候以半干旱、半濕潤氣候為主,太陽輻射強烈,日照強烈,降水主要集中于夏季,氣候變化分明,從西南向東北大致分為溫和半干旱區、溫涼半干旱區、濕涼半濕潤區、冷涼半濕潤區4個氣候類型區[16],流域內旱澇并存,土壤水分空間分布存在較明顯差異。流域內主要分布7種土類,分別為高山寒漠土、高山草甸土、亞高山草甸土、灌叢草原土、草甸土、潮土、新積土。中下游河谷區內土壤類型主要為草甸土、草原土,土壤石礫含量較高。流域內主要植被類型為高山沼澤草甸、高山灌叢草甸及少量高山灌叢[17-18]。

該研究區域位于29°20′~31°15′N、90°05′~93°20′E,主要包括拉薩市、墨竹工卡縣、林周縣、當雄縣、嘉黎縣等主要縣市。其中拉薩河河谷區位于流域拉薩河中下游地區,是青藏高原最重要的農業產糧地,主要產出糧食作物包括青稞、小麥、油菜等,表層土壤水分是當地農業生產的重要制約因素。因此,對拉薩河流域內土壤水分分布的研究對該地區農業發展具有重要意義。

1.2?數據來源

氣象降水數據來自青藏高原地區317個氣象站監測數據,使用Anusplin 42插值得到,葉面積指數數據(LAI)來自于美國航天局(NASA)MODIS15產品,時間分辨率為8 d,空間分辨率為500 m,流域土壤性質數據來自聯合國糧農組織(FAO)的HWSD土壤數據庫,使用基于Van-Genuchten模型的土壤水分特征曲線擬合軟件RETC處理后得到,流域內土地利用類型數據來自中國科學院資源環境科學數據中心的基于Landsat-TM/ETM影像的中國多時期土地覆被遙感數據監測數據集(CNLUCC),DEM數據來自中國科學院資源環境數據中心的全國1 km DEM數據產品。

根據多年降雨量的值,采用國內常用的分析枯水年、平水年、豐水年標準:①豐水年,Pi>+0.33W,其中Pi為第i年降水量,W為均方差,②枯水年,Pi<+0.33W,③Pi與+0.33W結果相近時為平水年。以此為依據計算得到2008年為豐水年、2015年為枯水年、2012年為平水年。

1.3?ANUSPLIN氣象數據插值

ANUSPLIN允許引入多元協變量線性子模型,擁有平穩處理二維以上樣條的功能,可引入多個氣象空間要素作為協變量,為考慮多要素下氣象數據的空間插值提供了可能。薄盤光滑樣條原型中擴展了局部薄盤光滑樣條(partial thin plate smoothing splines)部分[19],允許引入海拔、降水、海岸線關系等線性協變量因子,擴展了模型功能。這種模型的表述如下:

式(1)中,xi為獨立變量矢量,f為要估算的關于xi的未知平滑函數,yi為p維獨立協變量矢量,Zi為空間上第i點的因變量,b為yi的p維系數,N為觀測值個數。當p=0時(即協變量不存在時)即為普通薄盤樣條模型,當不存在獨立自變量f(xi)時,模型變為多元線性回歸模型。

由最小二乘估計來確定函數f和系數b:

式(2)中,Jm(f)定義為函數的m階偏導數,ρ為光滑參數(正值),在數據保真度與曲面的粗糙度之間起平衡作用。

1.4?RETC擬合土壤水分特征曲線

土壤水分特征為影響土壤含水量的重要因素,同時,土壤水分特征曲線也是分析土壤水分遷移的重要基本資料[20]。對缺乏資料地區,利用數值模擬方法對土壤水分運動進行預測,對于干旱程度、農牧業適宜性評價等領域具有重要意義[21-22]。在推求土壤水分運動過程及擬合土壤水分特征曲線的研究中,Van Genuchten、Gardner、GardnerRusso、BroodsCorey等土壤水分特征曲線擬合模型被廣泛用于描述水分在土壤中的入滲、遷移、蒸發等過程。其中Van Genuchten模型適用范圍廣、擬合效果好、采用最為廣泛,其水分特征曲線方程為:

式(3)中,θ(h)為土壤含水量,h為壓力水頭(cm),θr為土壤剩余含水量,θs為土壤飽和含水率,α和n為經驗擬合系數(或稱為土壤水分曲線形狀參數),m=1-1/n。

RETC軟件基于Broods-Corey模型和Van Genuchten模型對土壤水分特征曲線進行擬合,常用于非飽和帶土壤水力傳導特性和土壤水分的分析計算[23]。

1.5?MIKESHE分布式水文模型

MIKESHE(MIKE system hydrological european)模型是丹麥水利研究所在SHE模型基礎上研發而成。MIKESHE在濕地修復、洪水預報、水資源管理、水源地保護等流域具有廣泛的應用。MIKESHE主要包括蒸散發(ET)、融雪(SM)、非飽和帶(UZ)、坡面流(OL)、河流湖泊(OC)、飽和帶(SZ)6個模塊,每個模塊模擬不同的水文過程,不同模塊間可以獨立作用也可相互作用[24]。MIKESHE模型可根據需要選擇單一模塊或耦合多個模塊進行建模[25]。該研究選用ET與UZ模塊對拉薩河流域地表土壤含水量進行模擬。近年來,對于缺乏資料地區的水文模擬逐漸成為國際水文水資源研究領域的難點與熱點[25]。這使得水文模型在缺乏資料地區的模擬應用變為可能。

MIKESHE蒸散發模塊及非飽和帶下滲模塊計算基于以下算法:

(1)Kristensen-Jensen模型。

采用Rutter模型模擬水循環中截留過程,該過程可計算蒸散發量、冠層實際蓄水量和到達地面的凈雨量。使用氣候與植被資料直接計算潛在蒸散發能力。Rutter模型的實質功能是計算冠層蓄水量,其公式如下:

式(4)、(5)中,P是雨強(mm/s),S是冠層蓄水容量(mm),C是冠層實際含水量(mm),K、b是冠層蓄水參數,P1是地表植被覆蓋率,P2是總的葉片面積與植被覆蓋的地表面積之比, t是時間(s),Ep為土壤潛在蒸散發能力。

(2)植被實際蒸騰量。具體計算公式如下:

式(6)中,RDF為根系分布函數,f1(LAI)為葉片面積指數的函數,f2(θ)為土壤含水量函數。

對于發生在表層的土壤蒸發量,包括表層土壤蓄水量的蒸發量和基本蒸發量Epf3(θ)兩部分,其計算公式如下:

公式(8)、(9)中,θp為土壤田間持水量,θW為植物凋萎含水量,θM為土壤剩余含水量,C1、C2、C3為經驗常數。

2?結果與分析

2.1?表層土壤含水量空間變化特征

拉薩河流域按地貌可分為和園區、高山盆地及河谷區3類,流域內氣候類型屬于高原溫帶-寒溫帶半濕潤季風氣候,年內降水分布極不均勻,雨季一般為6月中旬之后開始,8月中旬左右結束,雨季降水量約占全年降水量的90%。利用拉薩河流域2008—2017年各年3月(枯水期)、7月(豐水期)降水、蒸散發、LAI數據,計算得到各年典型月(3、7月)表層土壤體積含水量及2008、2012、2015年的生長季3—9月表層土壤含水量。

由2008—2017年模擬結果可知,拉薩河流域表層土壤含水量自西南向東北呈逐漸減少趨勢,流域內盆地地區平均土壤表層含水量相對較高,包括位于拉薩河中游林周縣境內的澎波盆地、下游拉薩市、墨竹工卡縣中部境內拉薩河谷平原地區,其含水量在雨季為36%~40%,旱季為34%~36%,流域西部當雄縣境內,念青唐古拉山脈以東的中部洪積寬谷盆地、南部羊八井盆地,其含水量在雨季為38%~40%,旱季為34%~36%,流域東北地區嘉黎縣北部麥地卡濕地,其含水量在雨季。流域中部自拉薩河上游兩岸向東西兩側的高山地區平均土壤含水量相對較低,其含水量在雨季為34%~36%,旱季為30%~34%。各典型年旱季和雨季表層土壤平均含水量見圖1。

從拉薩河流域不同年份雨季與旱季的表層土壤含水量變化量空間分布(圖2)可以看出,在各不同年份雨季全流域土壤表層含水量呈增加趨勢,旱季除東北部少量河源濕地地區外全流域土壤表層含水量呈下降趨勢。在雨季,流域內除西部念青唐古拉山脈地區外,全流域土壤呈變濕潤態勢,其中南部河谷地區及西部高山盆地區域內表層土壤含水量上升了2%~4%,流域東北部河源濕地地區表層土壤含水量上升了1%~3%,流域中部山地地區土壤含水量上升了3%~5%,相較而言流域中部山地地區濕潤化程度更為明顯,在旱季,除流域東北部少數區域外,全流域呈變干旱趨勢,其中南部河谷地區土壤含水量下降了1%~2%,東北部河源濕地地區下降基本小于1%,西部高山盆地地區下降了1%~3%,流域中部山地地區下降了1%~4%。總體而言,流域中部山地地區干旱化程度更為顯著。

2.2?表層土壤含水量時間變化特征

計算得到2008—2017年各年平均含水量及各年生長季內土壤含水量變化量如圖3所示,結果顯示,各年雨季降水量充沛,雨季、旱季表層土壤含水量變化趨勢與年際間降水變化趨勢大致相同,在植物生長期內各年7月末8月初平均表層土壤含水量最大,各年3、4月表層土壤含水量最小,在枯水年7月由于植被蒸騰作用茂盛而降雨量稀少,導致水分支出失衡,故表層土壤呈干旱化趨勢。

從圖4可以看出,典型年內3—10月流域內表層土壤含水量變化趨勢隨枯水年、平水年、豐水年基本呈上升趨勢。流域內3—5月為旱季,降水較少,此時土壤水分變化較小,為33%~34%,7—8月為雨季,降水充沛,土壤含水量主要依靠降水補給,變化較大,為35%~37%。總體而言,土壤含水量年內變化大于年際變化。

2.3?土壤含水量與氣象要素及LAI的相關性分析

土壤含水量受多種影響因子如氣象、植被覆蓋類型、高程、土壤成分、土層厚度等共同影響,作用機理復雜,各影響因子間亦存在相互促進或削弱的作用。例如降水、蒸發、氣溫等氣象因子直接影響流域內土壤的補給水量,地表植被覆蓋類型的不同影響土壤蒸發量、植物蒸騰量、冠層截留水量及根系截留水量,高程及坡度影響流域內壤中流的走向[22],壤中流在重力作用下向海拔較低的地區匯流,故往往海拔較低的地區如流域內的河谷、盆地地區的土壤水分含量相對較高,這也與該研究所得結果相吻合,土壤類型、土壤成分(黏土、礫石組分)、土層厚度等因素則會影響土壤的蓄水能力與導水能力。土壤水分的分布結果是各類影響因子共同作用的結果,該研究主要探索氣象、植被覆蓋對土壤水分含量分布的影響,利用SPSS 18統計分析軟件對拉薩河流域2008—2017年土壤水分含量與降水、蒸發、葉面積指數(LAI)進行皮爾遜相關分析(Pearson),結果表明(表1),拉薩河流域內土壤水分含量與降水、LAI呈正向關關系,置信度較高,故可認為降水及植被覆蓋度是流域內土壤水分含量的重要作用因子。

3?結論與討論

(1)拉薩河流域表層土壤含水量自西南向東北呈逐漸減少趨勢,流域內盆地地區平均土壤表層含水量相對較高,流域南部河谷、西部盆地及東北部濕地地區土壤水分含量高于流域中部山地地區,各年旱季中,中部山地地區干旱化趨勢強于流域中河谷、盆地地區,各年雨季中,流域中部山地地區濕潤化趨勢強于流域中河谷、盆地地區。拉薩河流域河谷及盆地地區表層土壤砂石含量相對中部山地較高,植被覆蓋度更高。因此造成的土壤蓄水、導水能力的不同,形成了土壤水分空間變化的差異。總體而言,高海拔的流域中部地區土壤濕度相對較低,更適宜安排耐旱作物或品種。

(2)拉薩河流域各年雨季土壤平均含水量在33.5%~35.5%。干旱年份雨季土壤水分含量增加較少。各典型年內3—4月土壤含水量相對較低,變化幅度較小,5—6月降水充沛、植被相對茂盛的年份其土壤補給水量充沛、植被截留水量增加,相應的土壤水分含量升高,在7月末至8月初為植被最為茂盛的階段,由于蒸騰作用劇烈導致土壤水分略有下降,8月后降水減少,氣溫下降、植被覆蓋度降低,蒸騰減弱,故流域內土壤水分含量略有上升后持續降低。總體而言,土壤水分含量的年內變化相較于年際變化更為劇烈,可將流域內表層土壤水分的旱澇差異作為安排旱澇作物的依據。

(3)拉薩河流域尺度上土壤水分含量與降水、LAI呈正相關關系(P<0.05),與蒸發呈負相關關系(P>0.1),可能是因為降水及植物覆被直接影響土壤水分補給量,關系較為簡單,而蒸發受到地表植被蒸騰作用及土壤貯水、導水能力的影響,作用更為復雜,故而無明顯相關關系。之后的研究中可進一步加強對土壤水分及其主要影響因子的長時間序列監測,以揭示更穩定的土壤水分變化規律,這對調整布局農作物生產具有積極意義。

參考文獻

[1] 郭燾,于紅博,馬梓策,等.基于MODIS的土壤含水量時空變化及干旱化程度分析[J].水土保持研究,2019,26(4):185-189.

[2] BRONSTERT A,CREUTZFELDT B,GRAEFF T,et al.Potentials and constraints of different types of soil moisture observations for flood simulations in headwater catchments[J].Natural hazards,2012,60(3):879-914.

[3] 張淑芳,柴守璽,常磊,等.冬小麥不同播期對土壤水分及產量的影響[J].中國農學通報,2018,34(29):13-19.

[4] 連玉珍,曹麗花,劉合滿,等.西藏農牧交錯帶農田與相鄰草地土壤含水量空間分布特征[J].西南農業學報,2019,32(5):1092-1097.

[5] 高海峰,白軍紅,王慶改,等.霍林河下游典型洪泛區濕地土壤pH值和土壤含水量分布特征[J].水土保持研究,2011,18(1):268-271.

[6] BROCCA L,TULLO T,MELONE F,et al.Catchment scale soil moisture spatialtemporal variability[J].Journal of hydrology,2012,422/423:63-75.

[7] 劉雅莉,杜劍卿,李鋒,等.微尺度下城市公園人造綠地土壤水分的時空分異格局及其驅動機制[J].生態學報,2019,39(18):6794-6802.

[8] 楊凱悅,賈志清,張立恒,等.高寒沙區典型人工林土壤水分空間分布特征研究[J].干旱區資源與環境,2019,33(5):88-94.

[9] 汪倩倩,汪權方,王新生,等.地面資料稀缺區域的農田土壤水分微波與光學遙感協同反演方法研究[J].中國農學通報,2018,34(36):117-123.

[10] KOLASSA J,REICHLE R H,LIU Q,et al.Estimating surface soil moisture from SMAP observations using a Neural Network technique[J].Remote sensing of environment,2017,204:43-59.

[11] WANG L L,QU J J.Satellite remote sensing applications for surface soil moisture monitoring:A review[J].Frontiers of earth science in China,2009,3(2):237-247.

[12] KHAN S I,HONG ?Y,WANG J H,et al.Satellite remote sensing and hydrologic modeling for flood inundation mapping in Lake Victoria Basin:Implications for hydrologic prediction in ungauged basins[J].IEEE Transactions on Geoscience & Remote Sensing,2010,49(1):85-95.

[13] MASON D C,SPECK R,DEVEREUX B,et al.Flood detection in urban areas using terraSAR-X[J].IEEE Transactions on Geoscience & Remote Sensing,2010,48(2):882-894.

[14] GASHAW T,TULU T,ARGAW M,et al.Modeling the hydrological impacts of land use/land cover changes in the Andassa watershed,Blue Nile Basin,Ethiopia[J].Science of the total environment,2018,619/620:1394-1408.

[15] 徐宗學,程濤,洪思揚,等.遙感技術在城市洪澇模擬中的應用進展[J].科學通報,2018,63(21):2156-2166.

[16] GUO D L,FAN J,MI M X.Scaledependency of spatial variability of surface soil moisture under different land use types in Heihe Oasis,China[J].Chinese journal of applied ecology,2013,24(5):1199-1208.

[17] 李筱金,徐琳瑜.拉薩河流域水電開發帶來的河岸帶土壤特征變化研究[J].環境科學與技術,2015,38(5):148-156.

[18] 王慶海,李翠,龐卓,等.中國草地主要有毒植物及其防控技術[J].草地學報,2013,21(5):831-841.

[19] HONG Y,NIX H A,HUTCHINSONB M F,et al.Spatial interpolation of monthly mean climate data for china[J].International journal of climatology,2005,25(10):1369-1379.

[20] 劉貫群,孫蓓蓓,朱良超,等.內蒙孿井灌區包氣帶土壤水分特征研究[J].中國海洋大學學報(自然科學版),2011,41(Z1):107-112.

[21] 王全九,邵明安,鄭紀勇.土壤中水分運動與溶質遷移[M].北京:水利水電出版社,2007.

[22] 陳雪,宋婭麗,王克勤,等.基于Van Genuchten模型的等高反坡階下土壤水分特征[J].水土保持研究,2019,26(5):45-52.

[23] MA L,HE C G,BIAN H F,et al.MIKE SHE modeling of ecohydrological processes:Merits,applications,and challenges[J].Ecological engineering,2016,96:137-149.

[24] QIU Y,FU B J,WANG J,et al.Spatial prediction of soil moisture content using multiple-linear regressions in a gully catchment of the Loess Plateau,China[J].Journal of arid environments,2010,74(2):208-220.

[25] 胡相明,趙艷云,程積民,等.云霧山自然保護區環境因素對土壤水分空間分布的影響[J].生態學報,2008,28(7):2964-2971.

主站蜘蛛池模板: 小说区 亚洲 自拍 另类| 亚洲第一区在线| 国产成人毛片| 精品国产自在在线在线观看| 久久精品最新免费国产成人| 色综合天天综合| 最近最新中文字幕在线第一页| 中文字幕免费视频| 亚洲视频在线网| 欧美一级夜夜爽www| 亚洲综合片| 亚洲国产第一区二区香蕉| 国产亚洲成AⅤ人片在线观看| a色毛片免费视频| 女人爽到高潮免费视频大全| 美臀人妻中出中文字幕在线| 黄色网页在线观看| 久久精品国产精品国产一区| 国产内射一区亚洲| 国产麻豆91网在线看| 国产亚洲一区二区三区在线| 伊人成人在线视频| 最新无码专区超级碰碰碰| 美女一级毛片无遮挡内谢| 2021国产乱人伦在线播放| 91亚洲免费| 国产精品久久久久久久伊一| 亚洲天堂伊人| 99手机在线视频| 欧美α片免费观看| 久久久91人妻无码精品蜜桃HD| 亚洲免费播放| 欧美中文字幕第一页线路一| 国产日韩欧美在线视频免费观看| 999精品视频在线| 99久久免费精品特色大片| 一区二区无码在线视频| 国产AV无码专区亚洲A∨毛片| 亚洲一区二区黄色| 国产亚洲精品资源在线26u| 91小视频在线观看免费版高清| 激情综合网址| 国产欧美又粗又猛又爽老| 亚洲区第一页| 国产成本人片免费a∨短片| 午夜精品影院| 久久青草精品一区二区三区 | 日韩欧美国产区| 一级片一区| 日本影院一区| 欧美高清三区| a亚洲视频| 国产欧美日韩综合一区在线播放| 四虎国产精品永久一区| 国产不卡国语在线| 婷婷色丁香综合激情| 亚洲女人在线| 99久久精彩视频| 亚洲码在线中文在线观看| 日韩A∨精品日韩精品无码| 欧美综合成人| 国产无码高清视频不卡| 国产精品自在在线午夜| 中文字幕不卡免费高清视频| 亚洲精品麻豆| 国产成人精品一区二区三区| 国产亚洲精品97AA片在线播放| 国产精品毛片一区| 国产精品手机在线播放| 18禁色诱爆乳网站| 又粗又硬又大又爽免费视频播放| 美女啪啪无遮挡| 99久久人妻精品免费二区| 亚洲精品片911| 毛片卡一卡二| 国产精品免费电影| 国产v精品成人免费视频71pao | 国产成人精品亚洲日本对白优播| 国产福利影院在线观看| 自拍亚洲欧美精品| 欧美日韩福利| 天天躁夜夜躁狠狠躁躁88|