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

基于高光譜遙感的小麥黃花葉病害等級監測研究

2022-08-09 08:08:34劉夢冉張海艷齊雙麗孫炳劍段劍釗郭天財
麥類作物學報 2022年7期
關鍵詞:模型

劉夢冉,張海艷,齊雙麗,孫炳劍,段劍釗,賀 利,郭天財,馮 偉

(1. 河南農業大學農學院,河南鄭州 450046;2.河南農業大學機電工程學院,河南鄭州 450002;3. 漯河市農科院,河南漯河 462000;4. 河南農業大學植物保護學院,河南鄭州 450002)

小麥黃花葉病是一種低溫型土傳病毒病害,主要通過禾谷多黏菌()介體進行傳播。近年來受全球氣候變化、連年秸稈還田、小麥聯合收割機跨區作業等因素的綜合影響,該病害在我國黃淮及長江中下游麥區的發生呈逐年加重趨勢。據調查,在該病害一般發病田塊小麥產量損失30%左右,嚴重田塊減產可超過50%,甚至絕收。田間實地調查取樣是獲取該病害信息的主要方法,但所得信息病害預防而言具有明顯的滯后性,嚴重影響預報的準確率和有效防治的最佳時機。目前,遙感監測技術發展迅速,并在農業生產中得到廣泛應用。利用遙感技術有望提高小麥黃花葉病的早期預報和動態監測水平,對于有效防控其發生危害、減少農藥投入、保障該區小麥安全生產及可持續發展均具有重要意義。

植物病蟲害發生后葉片形態及內部結構會發生變化,最終改變植物的光譜曲線形狀,從而為遙感技術監測作物病蟲害提供了可能。目前,已有較多研究者嘗試通過分析光譜響應曲線,篩選出作物發生病蟲害后的敏感波段信息。如小麥白粉病敏感波段為580~710 nm;小麥條銹病敏感波段為630~687、740~890 及 976~1 350 nm;小麥蚜蟲敏感光譜波段為 625~635 nm 和680~695 nm;稻葉瘟病敏感波段為550 和680 nm。根據這些敏感波段光譜信息可以建立特征光譜參數與病蟲害嚴重度之間關系,或通過不同波段組合構建相對敏感的植被指數,用于提高作物病蟲害的遙感監測精度。如用DGSR(雙綠簡單比)和DGND(雙綠歸一化差)進行小麥白粉病遙感監測;利用近紅外光譜比值指數NSRI實現對冬小麥條銹病的有效監測;通過構建褐飛虱光譜指數(BPHI)來監測水稻作物的褐飛虱危害;光化學反射指數(PRI)可作為量化黃銹病感染的一個強有力的光譜指標。與其他病蟲害監測研究比較,病毒類黃花葉病害遙感監測相對偏少。研究表明,NDVI可以用來監測大豆的黃花葉病害;利用/比值可對豇豆黃花葉病進行有效監測;和/可用于大豆黃花葉病的監測;/的成像比值可以用來監測木薯黃花葉病。從上述研究結果看,用于植物病蟲害監測的光譜特征位置、敏感光譜指數及分析方法因植物種類、病蟲害類型存在很大差異;雖然利用光譜特征波段及敏感光譜參數可以監測豆類、薯類等作物黃花葉病害,而對低溫型作物而言,小麥黃花葉病的遙感監測預報研究則比較薄弱。

從生長形態、溫度及生育時期的響應看,小麥黃花葉病的發病特征不同于其他病害。當小麥受黃花葉病菌感染后,病株先從心葉葉尖或中部開始褪綠,初為淡綠色至橙黃色條或梭形點,變為黃色或淡綠色不連續的短線條,逐漸擴大成黃綠相間的斑駁或不規則條斑,條斑的中心可引起壞死。拔節后氣溫回升快,當高于20 ℃時,病癥逐漸消退,植株生長纖弱,株形松散,重病株在抽穗前多數分蘗甚至整株枯死。這種病害的特殊屬性決定了其冠層反射光譜的特有變化,以前用于專一指示其他病害的波段和參數將很難適用小麥黃花葉病害的監測。因此,需要加強低溫敏感性小麥黃花葉病的高光譜反射特征研究,明確生育時期對該病害監測的影響,確立病害監測適宜模型,為小麥黃花葉病的及時防控提供技術支持。本研究選取近些年小麥黃花葉病發病日趨嚴重的典型地塊,通過獲取不同等級黃花葉病害的小麥冠層高光譜數據,提取黃花葉病植株的特征光譜信息,明確適宜監測病害的光譜波段及植被指數,確立小麥黃花葉病害等級定量估算模型,以期為遙感大面積監測小麥黃花葉病害提供定量評價方法,增強農業災害的精準化管理水平。

1 材料和方法

1.1 試驗區概況

試驗于2016-2017年、2018-2019年和2020-2021年在河南省漯河市源匯區(33°57'N,114°02'E)進行,土壤類型為砂姜黑土。由于該區秋季作物騰茬早、陰雨天較多和光照不足,小麥聯合收割機跨區作業且高茬收割,輪作倒茬困難,大部分品種對黃花葉病的抗性差,導致在自然條件下當地小麥黃花葉病發病和蔓延迅猛。

1.2 病害調查

本試驗田間數據采集通過農戶田塊隨機調查方式進行。試驗區小麥種植品種主要包括百農207、周麥27、西農979、周麥22、鄭麥379和矮抗58。依據往年該區域地塊的發病情況,按照重度(4級)、中度(3級)、中輕度(2級)、輕度(1級)和健康(0級) 5個等級調查該病害田間實際發病程度。被調查地塊的黃花葉病在越冬期已出現零星發生,至返青期發病連片、癥狀特征明顯,因此田間取樣選擇在小麥返青期、拔節前期和拔節后期進行。 2016-2017年每個時期每個發病等級取5個樣點,2018-2019年每個時期每個發病等級取4個樣點,2020-2021年每個時期每個發病等級取6個樣點。小麥黃花葉病嚴重度分級標準參照孫炳劍等方法并作適當改進(表1)。

表1 小麥黃花葉病發病等級及發病癥狀

圖1 小麥黃花葉病不同發病等級癥狀圖

1.3 冠層光譜測量

小麥冠層光譜數據采集與病害調查同步,采用美國Analytical Spectral Device(ASD)公司生產的Field Spec Hand Held手持便攜式光譜分析儀進行。光譜測定波段值為325~1 075 nm,光譜采樣間隔為1.6 nm,光譜分辨率為3.5 nm。冠層光譜測定選擇在天氣晴朗無風的中午 (10:00-14:00)進行,傳感器探頭垂直向下,光譜儀視場角為25°,距離冠層頂垂直高度約1.0 m,冠層視場范圍面積為0.44 m。測量過程中及時進行標準白板校正,每個采樣點記錄10個光譜,以其平均值作為該點的光譜反射率。

1.4 生長特征參數調查

與光譜采集同步同位,每個小區選有代表性的植株 10 株,分別測定株高、葉片水分含量、色素含量及氮含量。

1.4.1 葉片含水量測定

將植株帶回實驗室稱取葉片鮮重,再放入烘箱 105 ℃殺青 30 min,然后在 80 ℃下烘干至恒重并稱取干重。葉片含水量=(WF-WD)/ WF×100%,WF和WD分別為葉片鮮重和干重。

1.4.2 葉片氮含量測定

將烘干至恒重的葉片取出研磨成粉末狀,稱取0.200 g放入K1100 全自動凱氏定氮儀測定葉片全氮含量。

1.4.3 色素含量測定

葉片去葉脈剪碎稱 0.200 g()置于50 mL()棕色瓶中用95%濃度的無水乙醇定容,避光24~48 h至葉綠素完全被提取后使用日立 U-2800型紫外分光光度計測出葉綠素 a和葉綠素 b 在 665 nm 和 649 nm 波段的吸光值(),再根據Lichtenthaler法算出葉綠素和類胡蘿卜素的濃度。公式如下:

Ca=13.95-6.88

Cb=24.96-7.32

Chl=(Ca+Cb) ×

C=(1 000-2.05Ca-114.8Cb)/245×

式中Ca、Cb、Chl和C分別為葉綠素 a、葉綠素 b、葉綠素和類胡蘿卜素的濃度。

1.5 數據計算處理

本研究采用Savitaky-Golay平滑法對原始光譜進行平滑,在Matlab語言環境下編程實現了光譜指數與黃花葉病病情等級間相關系數運算,選擇10個與病情等級最顯著相關的光譜特征參數列于表2,并將其輸入到偏最小二乘法(PLSR)模型中,進一步提高監測精度和穩定性。PLSR結合了主成分分析和多重線性分析兩種回歸方法,通過最小化交叉驗證的均方根誤差(RMSECV)來優化因子數量,可以規避模型潛在的過擬合問題。以2016-2017年和2018-2019年試驗資料為基礎,建立小麥黃花葉病病情等級監測模型。利用2020-2021年試驗獨立資料對模型精度進行驗證,采用決定系數()和均方根誤差(RMSE)進行評價。

表2 經典植被指數的計算方法

式中,分別為預測值和實測值,為樣本數。

2 結果與分析

2.1 小麥黃花葉病發病等級與株高及生理指標的相關性

小麥株高及葉片水分、氮素、胡蘿卜素和葉綠素含量在不同時期存在差異,其與黃花葉病病情等級的相關性也因時期而不同(表3)。除拔節后期的葉片水分含量外,各時期株高及四個生理指標與病害等級均呈顯著或極顯著負相關;從相關程度看,返青期和拔節前期各指標與病害等級的相關性均較好,整體優于拔節后期,其中拔節前期表現最佳,相關系數均超過0.80。

表3 不同時期小麥黃花葉病發病等級與形態生理指標的相關性

2.2 小麥黃花葉病不同發病等級下反射率變化特征

小麥冠層光譜反射率隨黃花葉病發病等級的增加呈現規律性變化(圖2)。在返青期和拔節前期,在560~720 nm波段范圍內,葉綠素對光的吸收減少,正常植株的反射率最高;隨著發病等級的加重,小麥葉片變黃,葉綠素含量減少,可見光波段反射率逐漸增加;而在近紅外反射區域800~900 nm波段間,細胞結構對光的反射開始增強,黃花葉病會破壞小麥葉片組織結構,隨病害等級的增加,反射率降低。在拔節后期,隨著氣溫升高,冬小麥開始快速生長,黃花葉病逐漸被遏制,不同病害程度的反射率在可見光區域差異變小,而在近紅外波段差異逐漸加大,病害等級越高,反射率越低。

“Health”為發病等級0級;“1”為發病等級1級;“2” 為發病等級2級;“3”為發病等級3級;“4”為發病等級4級。

2.3 病害等級及相關生理指標與反射率相關性

在返青期和拔節前期,小麥黃花葉病的病害等級與560~720 nm反射率呈正相關,相關系數均達0.80以上,而株高、葉片含水量、氮含量和色素含量與反射率呈負相關,除返青期葉片氮含量外相關性均較好(> -0.75);而在800~900 nm波段,病害等級、株高、葉片含水量、氮含量和色素含量與冠層光譜反射率的相關性相比可見光波段表現出相反趨勢,且相關性較差(-0.5<<0.3)(圖3)。這說明小麥黃花葉病病害等級及相關指標的敏感波段在返青期和拔節前期主要集中在560~720 nm。在拔節后期由于小麥快速生長及病菌被抑制,病害等級及相關指標的敏感波段主要集中在800~900 nm,病情等級與該波段光譜反射率的相關系數低于-0.60,而株高、葉片含水量、氮含量和色素含量與反射率均呈正相關 (=0.2~0.6),而可見光波段的相關系數隨波段位置變化較大,且相關性偏差。

圖3 小麥黃花葉病害等級及相關生理指標與冠層反射率間關系

2.4 植被指數與黃花葉病等級之間的線性回歸模型及驗證結果

2.4.1 植被指數與黃花葉病等級之間的線性回歸關系

以表征不同形態、生理指標的10個表現較好的常規植被指數為基礎,分析其與小麥黃花葉病害等級之間的定量關系,篩選指示病害等級的植被指數。表4列出了表現較好的回歸方程及其決定系數。在返青期,與生物量、色素和氮含量相關的植被指數的較高,EVI-1、WDRVI-3、mND705、VOG-1、PRI和NDDA的為0.70~0.74,均達到極顯著水平,其中以表征色素含量的mND705和VOG-1的最高(=0.74);在拔節前期,植被指數NDVI、EVI-1、WDRVI-3、mND705、PSNDb和NDDA的為0.80~ 0.82,均達到極顯著水平,其中以表征葉綠素b含量的PSNDb的最高(=0.82);在拔節后期,mND705、VOG-1和NDDA的為0.69~ 0.72,均達到極顯著水平,其中mND705的最高(=0.72)。綜合來看,植被指數與病害等級的相關性總體上好于原始光譜及株高、葉片含水量、氮含量和色素含量,其中表征色素變化的光譜指數mND705與病害等級的關系在3個發病時期表現最穩定且密切,其次為VOG-1和NDDA,因此可選擇植被指數對小麥黃葉病進行預測。

表4 不同生育時期植被指數與黃花葉病等級之間的定量關系

2.4.2 小麥黃花葉病等級監測模型精度的檢驗

利用2020-2021年度數據對建模效果較好的三個植被指數mND705,VOG-1和NDDA模型進行精度檢驗。結果(圖4)顯示,在返青期, 三個模型的驗證分別為0.68、0.66和0.65,RMSE分別為0.97、0.86和0.75;在拔節前期,驗證分別為0.78、0.67和0.75,RMSE分別為0.79、0.83和0.88;在拔節后期,驗證分別為0.59、0.52和0.51,RMSE分別為0.98、1.11和1.05。總體上以mND705模型的預測效果最好,精度最高,誤差較小。

a:返青期;b:拔節前期;c:拔節后期。

2.5 基于植被指數的小麥黃花葉病PLSR模型預測效果

為了進一步提高小麥黃花葉病害分級模型的監測精度以及穩定性,利用2016-2017和2018-2019年度數據,基于所選取的10個植被指數建立PLSR模型。圖5展示了返青期、拔節前期和拔節后期的PLSR模型所需的因子數及交叉驗證的均方根誤差(RMSE)值,以上三個時期的最佳因子數分別為前三因子、前二因子以及前九 因子。

圖5 基于偏最小二乘回歸的小麥黃花葉病害分級模型交叉驗證的均方根誤差

采用2020-2021年度獨立試驗數據對模型精度進行了檢驗(圖6)。總體而言,三個時期的PLSR模型的預測精度均高于植被指數模型。在返青期,PLSR訓練集=0.94,RMSE=0.32;驗證集=0.93, RMSE=0.29,PLSR模型預測精度比植被指數模型提高幅度較大(按照計算增幅27.02%~36.76%)。在拔節前期,PLSR模型預測精度較植被指數模型進一步提升,訓練集和驗證集分別為0.97和0.96,RMSE分別為0.24和0.26。在拔節后期,PLSR模型精度較植被指數也有一定程度的改善,訓練集和驗證集分別為0.84和0.81,RMSE分別為0.56和 0.61。由此可見,前兩個時期的PLSR模型精度較高(>0.90),而拔節后期較低。

Train:訓練集數據點;Test:驗證集數據點。

3 討 論

3.1 小麥黃花葉病的病征及反射率

小麥受黃花葉病感染后,病株葉片先從心葉葉尖或中部開始褪綠,初為淡綠色至橙黃色條或梭形點,變為黃色或淡綠色不連續的短線條,逐漸擴大成黃綠相間的斑駁或不規則條斑,條斑的中心可引起壞死。小麥黃花葉病的病害特征會隨著氣溫的上升而消失,當氣溫高于20 ℃時,病癥消退恢復生長,但此時已造成產量損失。不同病害侵染特征及對溫度的反應存在差異,其光譜曲線存在不同之處。如隨小麥白粉病和全蝕病病情嚴重度加重,可見光波段(550~750 nm)反射率升高,近紅外波段的反射率下降。隨著銹病對小麥脅迫程度的加重,紅邊發生“藍移”,黃邊發生“紅移”;水稻受到稻葉瘟危害后,隨著侵染加重光譜反射率,在綠光波段和近紅外區域呈現下降趨勢,紅光波段呈現上升趨勢,近紅光波段的紅邊區域“藍移”;棉花受到枯萎病侵染后,冠層光譜反射率相對于健康對照在綠光和近紅外波段下降。而本研究表明,小麥遭受黃花葉病后可見光波段光譜反射率逐漸增加,出現紅谷抬升現象,而在近紅外反射區域800~900 nm波段間,細胞結構對光的反射開始增強,隨病害等級增加,反射率降低。拔節后期,隨著氣溫升高,黃花葉病逐漸遏制,冬小麥恢復生長,反射率在可見光區域差異變小,在近紅外波段差異加大,病害等級越高,反射率越低。可見,小麥黃花葉病侵染后冠層光譜反射率變化不同于其他病害,利用反射光譜進行該病害監測需要考慮這一特點。

3.2 小麥黃花葉病的敏感波段

植物葉肉細胞中葉綠素含量等生化成分的不同,導致不同波段呈現出不同特征的光譜曲線。400~700 nm是植物的強吸收波段,葉綠素強烈吸收藍光和紅光,其反射率在480和680 nm附近非常低,在綠光范圍內有弱反射。在680~750 nm范圍內,反射率急劇上升,該區域是葉綠素在紅波段的強吸收到近紅外波段多次散射形成的高反射平臺的過渡波段,與植物單位面積葉綠素含量有關。700~900 nm出現“紅外高平臺”,此波段范圍的反射率主要受葉肉細胞結構和冠層結構的影響,反射率較高。植株感病后由病原體導致的植物生理代謝發生變化,出現病害癥狀,進而改變植物的光譜反應模式。前人對小麥病蟲害的敏感波段開展了較多研究,白粉病的敏感波段主要集中在580~710 nm范圍,而條銹病的敏感波段主要為630~687、740~890 及 976~1 350 nm,蚜蟲的敏感光譜波段范圍為 625~635和 680~695 nm,以上這些小麥病蟲害的敏感反應波段在生育時期間是一致的,較少因時期而變化。本研究結果顯示,小麥黃花葉病在不同生育時期的敏感波段有所差異,在返青期和拔節前期的敏感波段為560~720 nm,相關性較高;而在拔節后期,敏感波段主要集中在800~900 nm,且相關性較差,這可能與黃花葉病發病特征受氣溫影響較大有關。本研究提取的小麥黃花葉病敏感波段與其他作物黃花葉病的敏感波段也不盡相同。如豇豆黃花葉病的敏感波段為571、593、705和721 nm;大豆黃花葉病的敏感波段為445、688和750 nm。這可能是因為豆類的生育時期主要為氣溫較高的夏季,高溫高濕環境有利于豆類該病害發生;而小麥拔節前期的氣溫較低,低溫有利于小麥黃花葉病的發展,當拔節后期氣溫回升,對該病害的發生產生抑制,小麥生長開始恢復。因此,病害不同,作物類型不同,監測的敏感波段也不同。

3.3 小麥黃花葉病的適宜監測模型

植被指數由多波段組合構成,可以突出反映植被信息,減少背景等噪音影響,與監測指標的相關性優于單波段。在小麥返青期,黃花葉病導致新生葉片褪綠,葉片上呈現橙黃色梭形點或黃綠相間的不規則狀條斑,發病等級與色素含量的相關性最高。在本研究中,用來指示植株色素變化的mND705和VOG-1可以較好地指示黃花葉在返青期的發病狀態。在拔節前期,小麥植株發病癥狀明顯,平均株高、葉片含水量、葉片氮含量、葉綠素和類胡蘿卜素等生理指標變化較大,發病等級與各植被指數的相關性均較好。PSNDb可以很好指示色素濃度,EVI-1在監測和評估植被數量和時空變化方面具有很好的敏感性,WDRVI-3用來監測作物初級生產力效果較好,NDDA可以準確地估算小麥氮含量,以上植被指數均可以用來很好監測小麥拔節前期的黃花葉病情況。在拔節后期,隨著氣溫回升,病株癥狀快速消退并恢復生長,但此時小麥已產生了不可逆的損失,此時各植被指數的監測精度都較低。從不同時期的植被指數表現看,黃花葉病侵染后小麥返青期主要表現為葉色變黃,因此可以選擇使用色素類植被指數;拔節前期各種癥狀明顯,采用遙感手段監測病害的植被指數較多,包括氮素、色素及株高類植被指數等;而隨著發病時期推進,溫度升高,拔節后期的癥狀處于消失狀態,植株處于恢復狀態,此時監測精度不高。綜合來說,在三個時期表現均較好的mND705適宜監測小麥黃花葉病。這與Prabhakar等利用/對豇豆、Gazala等采用/對大豆黃花葉病進行監測的研究結果不一致,但在紅邊波段信息利用方面有一定相似之處。可見,作物類型不同,發病特點不一樣,所適宜的植被指數也不相同。另外,本研究中植株水分狀況與病害等級間相關性極顯著,但本試驗所用儀器的波段較窄、近紅外噪音較大,導致水分類植被指數(如WI)監測效果不好,在以后的研究中可考慮采用分辨率高、波動范圍寬的遙感監測儀器,加強水分類植被指數(如NDWI)在該病害監測中的應用。

雖然植被指數與單波段相比具有明顯優勢,但不同植被指數反映的信息不同,單獨建模很容易導致其他信息丟失。因此,多變量模型可以保證對所有信息的充分利用,達到對目標物綜合評價的最佳效果。PLSR不僅將共線變量簡化為幾個不相關的變量,而且使相關變量的協方差最大化,近年來在植物生長監測與評價中得到廣泛應用。本試驗將植被指數與PLSR算法相結合,在返青期和拔節前期均取得較好效果(>0.90),充分顯示了該模型在小麥黃花葉病害等級監測方面的優越性,在小麥生產中推薦使用該方法及時準確監測黃花葉病的發病狀況。

及早監測小麥黃花葉病的發病等級,有利于早期防控。本研究發現,不論是反射率還是植被指數與病情等級的相關性均在返青期和拔節前期表現較好,在拔節后期表現較差。這主要是因為黃花葉病的發生與土壤溫度和濕度密切相關,小麥生長進入返青期和拔節期,溫度維持在5~17 ℃范圍,有利于禾谷多黏菌介體傳播病毒。進入拔節后期,溫度明顯升高,病毒侵染明顯受到抑制,發病較輕的病株癥狀快速消退,感病嚴重的病株小分蘗也開始快速恢復生長,因病害導致的穗粒數顯著下降,造成的產量損失已無法挽回。與植被指數模型相比,PLSR模型明顯提升了返青期監測精度,在黃花葉病盛發期間的監測精度均高于0.9。因此,及時監測小麥返青至拔節前期的黃花葉病害所呈現的表型特征,有利于病害嚴重度等級評估及精確防控,對保證國家糧食安全及綠色可持續發展具有重要意義和價值。

猜你喜歡
模型
一半模型
一種去中心化的域名服務本地化模型
適用于BDS-3 PPP的隨機模型
提煉模型 突破難點
函數模型及應用
p150Glued在帕金森病模型中的表達及分布
函數模型及應用
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
3D打印中的模型分割與打包
主站蜘蛛池模板: 国产青青草视频| 亚洲无码一区在线观看| 国产第二十一页| YW尤物AV无码国产在线观看| 国产精品13页| 国产高清不卡视频| 一本无码在线观看| 亚洲欧美色中文字幕| 亚洲综合极品香蕉久久网| 九九视频免费在线观看| 美女无遮挡免费网站| 亚洲无码A视频在线| 亚洲永久色| 亚洲第一区欧美国产综合| 日韩人妻精品一区| 在线观看国产小视频| 99无码熟妇丰满人妻啪啪| 国产在线观看精品| 美女无遮挡拍拍拍免费视频| 国产激情第一页| 日韩在线永久免费播放| a毛片基地免费大全| 国产永久在线视频| 亚洲欧美不卡中文字幕| 国产欧美日本在线观看| 久久国产拍爱| 精品视频91| 91在线丝袜| 久久综合色播五月男人的天堂| 亚洲精品另类| 欧美综合一区二区三区| 999精品色在线观看| 亚洲日韩国产精品无码专区| 亚洲丝袜第一页| 亚洲中文字幕在线一区播放| 久久国产黑丝袜视频| 99国产在线视频| 青青草一区二区免费精品| 国产高清精品在线91| 亚洲国产欧美自拍| 国产日韩欧美精品区性色| 亚洲h视频在线| 国产精鲁鲁网在线视频| 国产91无毒不卡在线观看| 国产大全韩国亚洲一区二区三区| 国产爽妇精品| 在线观看网站国产| 国产午夜一级毛片| 天堂网国产| 亚洲最大情网站在线观看 | 日日拍夜夜嗷嗷叫国产| 天天综合色网| 国产亚洲精品自在久久不卡| 国产精品吹潮在线观看中文| 亚洲自偷自拍另类小说| 女人毛片a级大学毛片免费| 亚洲精品第五页| 制服丝袜 91视频| 日本国产精品| 2022国产91精品久久久久久| 亚洲高清在线播放| 成人毛片免费在线观看| 国产91精品最新在线播放| 92午夜福利影院一区二区三区| 亚洲综合日韩精品| 亚洲视频免| 香蕉99国内自产自拍视频| 国产资源免费观看| 亚洲一区精品视频在线| 国产交换配偶在线视频| 全部无卡免费的毛片在线看| 国产精品无码在线看| 天天摸天天操免费播放小视频| 国产成人综合久久| 国产内射一区亚洲| av手机版在线播放| 婷婷综合色| 91国语视频| 日韩国产综合精选| 思思热在线视频精品| 成年人久久黄色网站| 亚洲精品麻豆|