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

基于不同階微分高光譜植被指數的牧區草場地上生物量估算

2022-09-29 11:23:40楊震雷張亦然吳宇辰段利民
草地學報 2022年9期
關鍵詞:模型

童 新, 楊震雷, 張亦然, 吳宇辰, 段利民,2,3

(1. 內蒙古農業大學水利與土木建筑工程學院, 內蒙古 呼和浩特 010018; 2. 內蒙古自治區水資源保護與利用重點實驗室, 內蒙古 呼和浩特 010018; 3. 黃河流域內蒙段水資源與水環境綜合治理協同創新中心, 內蒙古 呼和浩特 010018; 4. 西湖大學工學院, 浙江 杭州 310024)

牧區草地生長旺盛期地上生物量(Aboveground biomass,AGB)是草場產量的一個重要指標,它的準確估算與監測是牧草高效管理、草畜供求關系平衡以及放牧制度優化等方面的前提[1-2]。傳統的收獲法(又稱樣方法)能準確估算草地地上生物量,但耗時長、破壞性大、估測范圍小、勞動強度高且無法實時獲取[3-4]。高光譜遙感技術的迅速發展,為快速、高效、定量估算草地地上生物量提供了一個重要的新途徑[5]。高光譜遙感的連續波段數量足、光譜分辨率高,其能在特定光譜區間獲取被測物體連續的反射光譜信息[6]。相較于多光譜遙感的常規寬波段光譜信息,高光譜遙感可以對植被特征進行微弱光譜差異的定量分析,還能進一步豐富植被地上生物量的估算方法,提高估算的精度[7]。

通過遙感估算草地地上生物量通常是利用綠色植被特有的光譜反射特征,選擇不同光譜波段(原始光譜或微分光譜)構建適宜的植被指數,再利用植被指數進行地上生物量的建模與估算。植被指數是不同敏感光譜波段的線性或非線性組合,可以比單一波段反映更多的額外隱含信息[8]。所有植被指數中最常用的有簡單比值植被指數(Simple ratio vegetation index,SRVI)和歸一化植被指數(Normalized difference vegetation index,NDVI),由于它們能較敏感的反映植被地上生物量而被廣泛應用[9-10]。高光譜遙感技術可以將地物光譜分離為更精細的窄波段,因此構建植被指數時對波段的選擇十分重要,合適的波段組合能有效提高研究的精度,同時波段優選也是高光譜遙感的重要研究內容[11-13]。傳統植被指數多是僅利用可見光、近紅外2個固定波段范圍的光譜信息,比如通過對比葉綠素引起的紅光波段強吸收與多重散射引起的近紅外波段高反射,卻沒有考慮其他波段或是在全部可利用波段內進行全面優選,信息量不足難以綜合分析而使得估算模型效果欠佳[14]。另外,牧草冠層在生長旺盛期的葉面積指數和郁閉度大,傳統紅光與近紅外波段組合的植被指數容易出現飽和現象,導致對地上生物量的敏感程度降低從而影響估算精度[15]。高光譜遙感則提供了基于所有可用光譜波段范圍內窄波段反射率構建植被指數的可能性,而不是僅僅局限于紅光波段和近紅外波段[16]。

關于草地地上生物量高光譜遙感估算的研究,大部分還是單獨基于原始或者一階微分光譜,借助常用的SRVI和NDVI 2種植被指數形式進行建模估算。微分光譜分析是一種能夠減少背景信號影響的技術,將反射率光譜視為d次的多項式函數,一階微分反射率則表示該函數隨波長變化的斜率,其消除了所有常數項(d=0)對該函數的影響。土壤光譜信號或陰影通常被認為是波長的線性函數,一階微分又不能完全去除背景信號,若利用二階微分反射率,則所有線性影響(d=1)都將被消除[17]。高宏元等[18]提取和計算了高寒草地原始、一階微分光譜以及植被指數等不同類型特征變量,構建了草地地上生物量估算模型。Bayaraa等[19]在蒙古中北部地區開展野外原位高光譜測量試驗,評價了不同植被指數估算牧草地上生物量的適用性。安海波等[11]利用傳統和所有兩波段組合優選的植被指數,開展了內蒙古天然和人工草地牧草地上生物量的對比估算研究,得出經波段優選的植被指數能顯著提高牧草生物量預測能力的結論。那么利用諸如二階、三階或四階等高階微分光譜,或是其他形式的植被指數,是否可以提高高光譜遙感估算草地地上生物量的精度是值得探究的,尤其是針對內蒙古地區天然草場生長旺盛期地上生物量估算的研究更是少見且值得深入開展。因此,本文針對位于內蒙古科爾沁沙地內部的天然草場,在牧民即將收獲牧草前開展冠層高光譜反射率采集與同期地上生物量收集試驗。對光譜反射率數據進行不同階微分處理,在全波段范圍內挑選最佳波段構建4種高光譜植被指數,建立相應地上生物量估算模型并對比評價各模型精度,為天然草場高效管理提供科學依據。

1 材料與方法

1.1 研究區域

研究區位于內蒙古自治區通遼市科左后旗的內蒙古農業大學阿古拉生態水文試驗區,地處科爾沁沙地東南緣,地理坐標為122°33′00″~122°41′00″ E,43°18′48″~43°21′24″ N,面積約55 km2。該區屬半干旱大陸性季風氣候,春季多風、夏季炎熱、秋季干燥、冬季寒冷。多年平均氣溫為6.6℃,多年平均降水量為389 mm,多年平均蒸發量(Ф20 cm口徑蒸發皿)1 412 mm[20-21]。研究區地貌特點為沙丘(流動、半流動、固定)、草甸相間分布,中部有一狹長小型淺水湖泊,名為王巴哈嘎湖。地勢南北高,中間低,且整體西高東低,海拔為184~235 m。

環繞王巴哈嘎湖的草甸天然草場為當地牧民收獲牧草的主要區域,草場的優勢植被物種為蘆葦(Phragmitesaustralis)、羊草(Leymuschinensis)和節節草(EquisetumramosissimumDesf.)。在王巴哈嘎湖北側草甸草場上選擇兩處采樣區域進行高光譜和地上生物量數據采集,采樣點根據前期野外實地踏勘與采樣區域附近實際情況確定(試驗路線需征得牧民同意)。具體來說,針對西側采樣區域1,牧民允許從最北端的起點開始,向南行進100 m距離,隨后可向東和向西各行進100 m;東側采樣區域2則是從中心為起點,可向東和向西各行進160 m、向南和向北各行進120 m。為了使獲取的地上生物量數據盡量滿足兩個采樣區域內的生物量空間差異性,實際采樣過程中,沿著牧民允許的行進路線,有針對性地選擇不同植被茂密程度的采樣點。采樣區域1的采樣點整體呈現倒“T”型,共布設30個采樣點;采樣區域2的采樣點整體呈現“十”字型,也布設30個采樣點(圖1)。兩處采樣區域累計布設60個采樣點,為保證高光譜數據與地上生物量數據空間匹配,采樣點均用高精度GPS定位并在地面插小紅旗標記。

圖1 研究區及試驗點所在位置Fig.1 Location of the study area and the sampling sites

1.2 草地冠層高光譜數據采集

(1)光譜測量時間

分別于2019年7月31日10∶40和12∶10對兩個采樣區域布設的各30個采樣點進行高光譜反射率數據采集,采樣過程分別持續1.5 h左右。當天天氣晴朗,光照穩定,微風且無云。

(2)光譜測量儀器

草地冠層高光譜反射率數據采用美國ASD(Analytical Spectral Devices,Inc.,Boulder,CO,USA)公司生產的Field Spec4可見光/近紅外便攜式地物高光譜儀進行采集。儀器視場角為25°,波段范圍為350~2 500 nm,其中350~1 000 nm范圍的光譜采樣間隔為1.4 nm,光譜分辨率為3 nm,1 001~2 500 nm范圍的光譜采樣間隔為2 nm,光譜分辨率為10 nm。隨儀器附帶的RS3軟件可對350~2 500 nm光譜范圍的原始反射率、輻照度和輻亮度光譜數據進記錄、存儲和顯示。

(3)光譜測量

野外光譜測量前,先對光譜儀進行預熱、優化處理,隨后進行參考板檢驗和比對(1 min內完成),且在測量過程中根據環境情況改變及時進行參考板校正。測量方式為試驗人員手持,光譜傳感器探頭垂直向下,距采樣點處冠層頂高度75 cm左右。在儀器輸出光譜設置項中,每條光譜的平均采樣數為10,每個采樣點又進行2次重復測量,因此每個采樣點共采集記錄20次高光譜反射率數據。

1.3 草地地上生物量數據獲取

高光譜反射率數據采集后立即對采樣點處的地上生物量進行獲取,地上生物量獲取點對應光譜采樣點,收割牧草的樣方范圍則根據光譜儀視場角在地面的投影面積來確定(0.3~0.5 m2)。齊地收割樣方內所有牧草,現場稱其鮮重并裝入紙袋,隨后帶回實驗室置于烘箱內,105℃殺青30 min,再以65℃烘干至恒重,稱其干質量后計算地上生物量。

1.4 數據處理與分析

光譜儀自帶的軟件程序可將野外原位試驗采集的光譜數據重采樣至逐波段高光譜反射率曲線,共輸出2 151個波段。每個采樣點共獲取了20條高光譜曲線,異常光譜剔除后對余下各光譜波段取均值作為各采樣點最終光譜反射率曲線。草地地上生物量在可見光至近紅外波段內與光譜反射率關系密切[16],且所用光譜儀在400~900 nm波段范圍的信噪比佳,故本文選取該波段范圍高光譜反射率(共501個波段)進行研究。對初始光譜數據進行平滑去噪、一階至四階光譜微分變化處理,微分光譜計算式如下:

(1)

式中:DR(λi)為λi波段處高一階(一階至四階)的微分反射率,R(λi)和R(λj)是在λi波段和λj波段處低一階(原始至三階)的反射率,Δλ為微分窗口大小,Δλ=λi-λj(λi>λj)。由于過大的微分窗口會導致重要光譜特征的丟失,同時窗口小于光譜分辨率又會導致信號失真。本研究選取的400~900 nm波段范圍光譜分辨率為3 nm,綜合考慮下,微分窗口Δλ大小選為5 nm。

將原始、一階至四階高光譜反射率與草場地上生物量進行逐波段相關性分析,尋找與地上生物量相關性強的波段。除了SRVI和NDVI,還借助土壤調節植被指數(Soil adjusted vegetation index,SAVI)、增強型植被指數(Enhanced vegetation index,EVI)這2種常用的植被指數形式,將400~900 nm范圍內原始與不同階微分逐波段光譜反射率代入其公式中(表1)。根據相應指數形式,選擇2個或3個單一波段的反射率構建植被指數,計算其與草地地上生物量的相關性并繪制相關性矩陣圖,同時挑選出由最佳波段組合構成的不同階(原始相當于零階)最佳高光譜植被指數。

表1 本研究選用的植被指數Table 1 Equations of vegetation indexes in this study

1.5 模型建立與驗證

利用一元線性回歸法分別建立不同階高光譜最佳植被指數與草場地上生物量的回歸估算模型。數據集劃分時,考慮采樣點數量與樣本代表性,對60個地上生物量數據進行降序排列,從前兩個數據開始,采用隔一選二的方式,選擇訓練集樣本,剩下的數據為驗證集樣本。共選取40個訓練樣本建立回歸模型,剩余20個樣本對模型進行精度驗證。模型精度采用決定系數(Coefficient of determination,R2)和均方根誤差(Root mean square error,RMSE)來評價,其中R2越接近于1,RMSE越小,所構建模型的精度越高。本文所有數據處理、分析計算、模型建立與繪制圖形等都在MATLAB R2016a軟件環境下完成。

2 結果與分析

2.1 草地冠層不同階高光譜反射率和地上生物量的關系

對草地地上生物量與逐波段原始、一階至四階高光譜冠層光譜反射率進行線性相關性分析,得到不同階光譜反射率與地上生物量相關系數隨波長變化的規律(圖2和圖3)。由圖2可知,對于生長旺盛期,草地地上生物量與單一波段原始高光譜冠層反射率相關性整體偏低,只在660~680 nm(紅光波段)、750~900 nm(近紅外波段)兩個波段范圍形成相關性稍好高值平臺。由圖3可知,光譜微分階數越高,相關性曲線的波動越劇烈。一階至四階微分高光譜反射率與草地地上生物量的相關性曲線沒有出現明顯高相關性高值平臺,但從400~900 nm波段范圍的單一波段相關性極值來看,微分處理后高光譜反射率優于原始反射率。

圖2 草地原始高光譜反射率與地上生物量的相關性曲線Fig.2 Correlograms of the correlation coefficient (R) between the original hyperspectral reflectance and AGB

圖3 草地一至四階高光譜微分反射率與地上生物量的相關性曲線Fig.3 Correlograms of the correlation coefficient (R) between the first-,second-,third- and fourth order derivative hyperspectral reflectance and AGB

2.2 不同階高光譜植被指數與地上生物量的關系

以簡單比值植被指數(SRVI)、歸一化植被指數(NDVI)、土壤調節植被指數(SAVI)和增強型植被指數(EVI)4種植被指數形式對400~900 nm范圍逐個光譜反射率進行組合,構建基于原始光譜反射率、一階至四階微分光譜反射率的5類不同植被指數,并繪制相關性矩陣圖分析其與草地地上生物量的相關性,圖4至圖8分別給出了5類植被指數與草地地上生物量的相關性矩陣。由于NDVI和SAVI的相關性矩陣圖沿著對角線對稱且正負相關性相反,僅繪制對角線右下側的圖形。

對于基于原始高光譜反射率構建的SRVI,NDVI和SAVI三種植被指數,如果不考慮SRVI的非對稱性,SRVI,NDVI和SAVI的矩陣圖形式相對一致(圖4)。而基于一階至四階微分高光譜反射率構建的3種植被指數相關性矩陣圖,相似性隨著微分階數的增加而不斷減弱,僅基于一階微分高光譜反射率的SRVI和NDVI矩陣圖相似性尚佳(圖5),說明隨著光譜微分階數的增加,全波段范圍內構建的不同植被指數間的差異性也在不斷增加。基于原始高光譜反射率的植被指數與草地地上生物量相關性矩陣圖中,存在較高相關性的區域:波段組合在730~900 nm與700~740 nm的區域為SRVI,NDVI和SAVI的高相關性區域,波段組合為735~755 nm,755~785 nm與650~670 nm的區域為EVI的高相關性區域(圖4)。基于微分高光譜反射率構建的植被指數與地上生物量相關性矩陣圖呈現明顯的斑塊化(圖5,圖6),且斑塊化程度隨著階數的增加而增大,并逐漸破碎化(圖7,圖8),因此越是高階微分光譜,矩陣圖內成片的高相關性區域越不明顯。

圖4 基于原始高光譜反射率的4種植被指數與牧草地上生物量的相關性矩陣圖Fig.4 Correlation matrix plots of the four original hyperspectral reflectance vegetation indexes and AGB

以原始高光譜反射率、一階至二階微分高光譜反射率構建的4種植被指數與草地地上生物量的相關系數絕對值的最大值都在0.63以上,三階微分高光譜反射率植被指數與地上生物量的相關性較前者有所下降,四階情形下相關性更低一些。其中,基于一階微分高光譜反射率并以SRVI形式構建的植被指數(兩波段分別為712 nm和749 nm)與地上生物量的相關性最好(相關系數為—0.68)。原始與一至四階微分高光譜反射率植被指數的最佳波段組合見表2。

圖7 基于三階高光譜微分反射率的4種植被指數與牧草地上生物量的相關性矩陣圖Fig.7 Correlation matrix plots of the four third-order hyperspectral derivative reflectance vegetation indexes and AGB

圖8 基于四階高光譜微分反射率的4種植被指數與牧草地上生物量的相關性矩陣圖Fig.8 Correlation matrix plots of the four fourth-order hyperspectral derivative reflectance vegetation indexes and AGB

由表2可知,基于原始光譜反射率與一階至三階微分光譜反射率構建的不同最佳植被指數的主要組成波段都以紅邊、近紅外波段為主,且SRVI與NDVI的兩波段組合具有高度的相似性,SAVI與EVI的前兩個組成波段有著高度相似性。基于四階微分光譜反射率的最佳植被指數的波段組成則較分散,藍、綠、紅光以及紅邊、近紅外波段都有涉及。整體上看,就與草地地上生物量的相關性而言,基于一階微分光譜反射率構建的最佳植被指數最高,二階微分光譜反射率略微高于原始光譜反射率,三階、四階微分反射率依次遞減。

表2 最優原始與一至四階微分反射率植被指數及其波段組合Table 2 Optimal vegetation indexes against AGB based on the original and different order derivative reflectance

2.3 草地地上生物量的最優植被指數估算模型對比

基于原始、一階至四階微分5類光譜反射率數據,利用最佳波段組合構建的最佳SRVI,NDVI,SAVI和EVI植被指數,建立了共20個最佳植被指數草地地上生物量線性估算模型,各模型及評價結果見表3。由表3可知,基于原始光譜反射率數據構建的4個植被指數模型估算效果相近,決定系數R2都在0.54以上,均方根誤差RMSE不超過231.13 g·m-2。基于一階和二階微分光譜反射率的植被指數模型中,SRVI和NDVI模型的精度高于SAVI和EVI模型,除一階微分光譜反射率NDVI模型外,其余3個模型的決定系數都達到了0.60,RMSE均小于219.74 g·m-2。整體來看,基于原始和一階微分光譜反射率的4類最佳植被指數模型精度相當,基于二階微分光譜反射率的最佳SRVI和NDVI模型精度最優(決定系數分別為0.69和0.70,RMSE分別為196.60 g·m-2和196.65 g·m-2)。而隨著光譜微分階數的提高,模型的精度呈下降趨勢。圖9給出了精度最高的基于二階微分光譜反射率最佳SRVI與NDVI模型的估算及評價結果,由圖可知二階微分SRVI和NDVI模型估算效果十分接近,在低地上生物量區域模擬精度較好且對地上生物量高值區域存在一定的低估。這是由于全波段范圍優選構建的二階微分植被指數能夠在一定程度避免傳統紅光與近紅外波段組合的植被指數容易出現的飽和現象,但是通過單一植被指數構建的線性模型精度還有進一步上升的空間。

圖9 基于二階微分光譜反射率的最優SRVI與NDVI模型估算結果Fig.9 Results of the models using the optimal SRVI (a) and NDVI (b) based on the second-order derivative reflectance

表3 最優原始與一至四階微分反射率植被指數模型估算結果Table 3 Results of the models using the optimal vegetation indexes based on the original and different order derivative reflectance

3 討論

草地生長旺盛期單一波段冠層原始高光譜反射率與地上生物量的相關性較弱(圖2),這主要是受草地植被種類與冠層立體結構綜合干擾的結果,且地上生物量又是一個受多因素影響的植被參數,其變化不能僅靠原始冠層反射率的變化來體現[26]。通過微分轉換后,不同階微分高光譜敏感波段反射率與草地地上生物量的關系要好于原始高光譜敏感波段反射率(圖3),說明光譜微分變換可削弱背景噪聲,強化并突出高光譜特征信息[27-29]。

利用所有波段范圍內原始、一階至四階微分高光譜反射率,挑選出的20個由最佳波段構建的最佳植被指數與草地地上生物量的相關性比最佳單波段反射率高,且組成最佳植被指數的波段以紅邊和近紅外波段占比最大,分別占所有已挑選波段的49%與33%(表2)。高宏元等[18]對高寒草地地上生物量進行高光譜估算時所選取的特征波段也主要集中在紅邊和近紅外這兩個波段范圍。由紅光過渡到近紅外的紅邊波段是描述植被生長狀況的重要指示波段,有研究表明紅邊波段所包含的光譜信息能表征80%以上的包括地上生物量在內的多種植物理化參數[30-31]。從表2還可以發現,20個最佳植被指數中僅有1個植被指數(基于四階微分光譜的NDVI最佳波段組合為近紅外與紅光波段)的最佳波段屬于表1所列傳統波段的范疇,說明在全波段范圍內優選構建的植被指數與草地地上生物量的相關性要好于利用傳統波段構建的植被指數。

基于原始高光譜反射率的4個最佳植被指數估算模型精度相近,植被指數形式復雜的SAVI,EVI模型略好于形式簡單的SRVI,NDVI模型。而基于一階至三階微分光譜反射率的最佳植被指數估算模型卻相反,形式更為簡單的SRVI,NDVI模型精度優于SAVI,EVI模型。此外,基于原始、一階至四階的SRVI和NDVI兩類模型,隨著微分階數的提高,模型精度先增加后減少,轉折點出現在基于二階微分光譜的植被指數模型。對比前者,SAVI和EVI模型精度卻在一直下降。雖然前面提到,光譜微分變換可以削弱背景噪聲影響,在一定程度上強化光譜特征,但是光譜轉換過程會導致潛在敏感信息的丟失[32],且隨著階數的提高,高頻噪聲會被逐漸放大,信噪比降低[33],繼而影響了估算模型的精度。

高光譜遙感可提供比多光譜更精細的連續窄波段信息,但是如何從眾多信息中挑選出有效信息也是個不小的挑戰。Fava等[13]研究表明,相較于波寬(Full width half maximum,FWHM)而言,波段選擇的微小偏移會帶來草地地上生物量估算的較大差異。本研究借助地面遙感手段,對比討論了利用不同階微分光譜、不同最佳植被指數對草地地上生物量的估算,優選出的最佳波段組合能夠為航天衛星多光譜與航空無人機高光譜影像的波段選取提供參考依據。未來應該在此基礎上,探尋地面高光譜波段和航天或航空遙感影像波段間的對應關系,進而實現“空-天-地”協同的區域草地地上生物量模型構建和估算。

4 結論

微分高光譜敏感波段反射率與草地地上生物量的關系要好于原始高光譜敏感波段反射率,全波段優選構建的最佳植被指數能進一步提高相關性。最佳植被指數的組合波段中,紅邊和近紅外是重點考慮波段,兩者能占所有已挑選波段的80%以上。過高的微分階數會降低估算模型的精度,利用原始、一階至四階光譜微分構建的最優SRVI,NDVI,SAVI和EVI模型中,基于二階微分光譜的SRVI和NDVI模型精度最好(R2分別為0.69和0.70,RMSE分別為196.60 g·m-2和196.65 g·m-2),表現出良好的準確性與便捷性,不失為估算草地地上生物量的一種方便高效方法。

猜你喜歡
模型
一半模型
一種去中心化的域名服務本地化模型
適用于BDS-3 PPP的隨機模型
提煉模型 突破難點
函數模型及應用
p150Glued在帕金森病模型中的表達及分布
函數模型及應用
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
3D打印中的模型分割與打包
主站蜘蛛池模板: 国产不卡国语在线| 热久久综合这里只有精品电影| 国产成人福利在线| 国产女人18毛片水真多1| 亚洲精品爱草草视频在线| 亚洲人成网站在线播放2019| www精品久久| 91久久国产综合精品女同我| 免费高清毛片| 亚洲制服丝袜第一页| 成年人视频一区二区| 精品久久久久久成人AV| 精品1区2区3区| 五月天婷婷网亚洲综合在线| 1级黄色毛片| 亚洲三级色| 亚洲成在线观看| 日本不卡在线| 亚洲成A人V欧美综合天堂| 亚洲男人天堂网址| 国产精品福利社| 毛片在线播放a| 久久精品无码国产一区二区三区| 国产成人精品一区二区三在线观看| 亚洲天堂视频在线观看| 蜜臀AV在线播放| 制服丝袜一区二区三区在线| 一级片一区| 91成人免费观看| 三上悠亚一区二区| 国产在线一区视频| 国产成人综合欧美精品久久| 国产成人精品日本亚洲77美色| 欧美第一页在线| 麻豆精选在线| 精品国产www| 国产精品免费露脸视频| 色久综合在线| 欧美一级高清视频在线播放| 亚洲va视频| 亚洲第一视频网| 国产区精品高清在线观看| 国产午夜人做人免费视频中文| 无码免费的亚洲视频| 老司国产精品视频91| 国产一区二区三区免费观看| 好紧太爽了视频免费无码| 成人午夜视频免费看欧美| 72种姿势欧美久久久久大黄蕉| 久久亚洲国产一区二区| 女人18毛片一级毛片在线| 国产精品任我爽爆在线播放6080| 九色视频在线免费观看| 欧美日韩福利| 欧美激情福利| 欧美另类图片视频无弹跳第一页| 日韩免费毛片视频| 国产成熟女人性满足视频| 激情综合网址| 在线亚洲精品福利网址导航| 波多野结衣在线se| 国产在线精品人成导航| 热久久综合这里只有精品电影| 91久久偷偷做嫩草影院| 乱人伦99久久| 麻豆国产在线观看一区二区| 3344在线观看无码| 久草视频精品| 一区二区三区四区精品视频 | 国产在线啪| 无码av免费不卡在线观看| 毛片在线播放网址| 91网址在线播放| 国产免费久久精品44| 亚洲乱亚洲乱妇24p| 日韩免费视频播播| 亚洲日本在线免费观看| 又污又黄又无遮挡网站| 91精品视频播放| 特级aaaaaaaaa毛片免费视频| 91视频区| 日韩在线2020专区|