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

基于高光譜數(shù)據(jù)的水稻葉片凈光合速率反演

2020-03-17 08:50:16方圣輝
關(guān)鍵詞:水稻區(qū)域模型

劉 聰 彭 漪 方圣輝

(武漢大學(xué) 遙感信息工程學(xué)院,武漢 430079)

光合作用是綠色植物生長過程中的重要生理活動,也是生態(tài)系統(tǒng)碳循環(huán)的基礎(chǔ)[1]。凈光合速率(Pn, Net photosynthefic rate, μmol/(m2·s))是衡量植被光合作用能力的直接指標(biāo),對于植被生產(chǎn)力的研究而言,其可與總初級生產(chǎn)力近似等同,凈光合速率值越高,表明葉片的結(jié)構(gòu)和功能狀況越好[2]。傳統(tǒng)研究方法以實(shí)地采樣測量為主,耗時(shí)費(fèi)力,而遙感技術(shù)具有快速、動態(tài)、大尺度的特點(diǎn),可以為植被光合作用研究提供更為便捷的手段。水稻作為我國第一大糧食作物,研究其不同氮肥水平下的凈光合速率快速反演,對于水稻的田間管理和生產(chǎn)力評估有著重要意義。

植被生產(chǎn)力是太陽發(fā)出的輻射被植被冠層截獲的結(jié)果[3]。目前,采用遙感手段利用光能利用率模型對植被生產(chǎn)力進(jìn)行大尺度的定量觀測已成為熱點(diǎn),其基礎(chǔ)為Monteith所提出的光能利用率方程[4-5]:

GPP=fAPAR×PAR×LUE

(1)

式中:fAPAR(Fraction of absorbed photosynthetically active radiation)為植被吸收的光合有效輻射占入射光合有效輻射的比例;PAR(Photosynthetically active radiation)為光合有效輻射;LUE(Light use efficiency)為光能利用率,即植被將所吸收的光合有效輻射轉(zhuǎn)化為有機(jī)碳的效率[6]。該模型的關(guān)鍵在于fAPAR和LUE的估算,相關(guān)研究表明,其均與冠層葉綠素含量密切相關(guān)[7-8],基于此Gitelson等[9-11]提出了基于冠層葉綠素含量的植被生產(chǎn)力模型:

GPP∝Chlcanopy×PAR

(2)

葉是綠色植物光合作用的主要器官,利用遙感手段快速無損估測植被葉片葉綠素含量已相當(dāng)成熟[7]。在葉片尺度上,可使用葉綠素相關(guān)植被指數(shù)作為葉綠素含量的替代來對葉片凈光合速率進(jìn)行反演。相關(guān)研究有:陳俊英等[12]和雷利琴[13]使用冠層光譜反射率計(jì)算得到的植被指數(shù)分別對棉花和油菜葉片凈光合速率進(jìn)行反演,驗(yàn)證了該方法的可行性,但精度有限;孫少波等[14]發(fā)現(xiàn)使用經(jīng)小波變換后的光譜反射率所構(gòu)建的植被指數(shù)可以提高毛竹葉片凈光合速率反演結(jié)果精度。葉片是植被進(jìn)行碳匯釋氧的主要場所,考慮到目前葉片尺度上的植被生產(chǎn)力反演研究較少,且多未考慮光合有效輻射對植被光合作用的影響,因此,本研究將引入光合入射有效輻射建模,對比引入前后模型精度差異,旨在建立不同氮水平下,葉片尺度水稻凈光合速率快速反演模型。

1 材料與方法

1.1 試驗(yàn)材料

于武漢大學(xué)雜交水稻國家重點(diǎn)實(shí)驗(yàn)室(18°32′28″N, 110°03′10″E)的試驗(yàn)田進(jìn)行。該試驗(yàn)田位于海南省陵水黎族自治縣光坡鎮(zhèn)官坡村,屬于熱帶季風(fēng)性氣候,光照充足,年日照時(shí)數(shù)2 400 h左右,多年平均氣溫25 ℃,年降水量1 700 mm左右,降雨主要集中于每年5—10月,終年無霜。

研究區(qū)如圖1所示,試驗(yàn)田總面積約為 518 m2,包含12個裂區(qū),每個裂區(qū)長10 m、寬3 m,種植珞優(yōu)9348和豐兩優(yōu)4號2個品種的水稻。2個品種的水稻均設(shè)置4個氮肥梯度:N0、N120、N180和N240,分別代表純氮0、120、180和240 kg/hm2,裂區(qū)間田埂上鋪有0.4 m寬黑色塑料厚膜用于隔水隔肥,同時(shí)對各施氮水平進(jìn)行3次重復(fù),共計(jì)24塊試驗(yàn)小區(qū)。珞優(yōu)9348和豐兩優(yōu)4號水稻均于2017年12月10日育秧,2018年1月5日進(jìn)行移栽,移栽密度為300 000 株/hm2,田間水稻生長周期約為6個月。

在本研究中,考慮到大田中儀器架設(shè)和測量的不便,選擇了施氮水平分別為N0、N180、N240的小區(qū)各2個,共計(jì)6個。自2018年3月15日起,到4月5日,對2個品種,3個氮水平下的水稻葉片進(jìn)行了光譜和光合數(shù)據(jù)采集。

1.2 光合數(shù)據(jù)測量

光合數(shù)據(jù)采集使用的儀器為美國LI-COR公司生產(chǎn)的LI-6400XT便攜式光合作用儀,該儀器可以測量植被凈光合速率Pn、入射光強(qiáng)PARin以及氣孔導(dǎo)度Gs、蒸騰速率Tr、胞間CO2濃度Ci等光合相關(guān)參數(shù)。

數(shù)據(jù)采集自3月16日開始,持續(xù)21 d,均選擇晴朗無風(fēng)無云天氣條件,有效數(shù)據(jù)共計(jì)采集14次。在每個小區(qū)隨機(jī)選擇3株長勢均一且良好、無病蟲害影響的水稻,依次編號標(biāo)記,控制葉室內(nèi)溫度恒定且適宜,測定完全展開葉片的凈光合速率Pn和光合相關(guān)參數(shù),每個點(diǎn)待曲線穩(wěn)定后記錄5次,取均值作為該點(diǎn)的光合數(shù)據(jù),共計(jì)測得有效光合數(shù)據(jù)247條。

圖1 研究區(qū)域示意圖
Fig.1 The location of study area

1.3 光譜數(shù)據(jù)測量

光譜數(shù)據(jù)測定使用的儀器為美國ASD公司研制的FieldSpec Pro FR2500背掛式野外高光譜儀,該儀器波段范圍為350~2 500 nm,其中在350~1 000 nm 區(qū)間采樣間隔為1.4 nm,經(jīng)重采樣處理后所得光譜數(shù)據(jù)波長間隔1 nm。

在每次光合數(shù)據(jù)測量完畢后,需盡快按編號對相應(yīng)葉片的光譜數(shù)據(jù)進(jìn)行測定,每個點(diǎn)同樣記錄5條光譜值,取其平均值作為該點(diǎn)的光譜反射率數(shù)據(jù)。另外,每次光譜測定前需使用白板進(jìn)行校正,以保證光譜數(shù)據(jù)的準(zhǔn)確性。

1.4 葉綠素相關(guān)植被指數(shù)

基于遙感手段的植被葉綠素反演方法已相當(dāng)成熟,在中低植被覆蓋度的條件下,歸一化植被指數(shù)(Normalized difference vegetation index, NDVI)作為最常用的植被指數(shù)[15],可以很好地對葉綠素含量進(jìn)行反演,然而在高葉綠素含量條件下,由于近紅外波段反射率遠(yuǎn)高于紅波段,NDVI會趨近于飽和[16-18];具體地,當(dāng)葉綠素含量>2 g/m2時(shí),NDVI失效[7]。Huete等[19]提出的增強(qiáng)型植被指數(shù)(Enhanced vegetation index, EVI)對于中高植被覆蓋區(qū)敏感度更優(yōu)[19],同時(shí)在對于作物生產(chǎn)力預(yù)測中,表現(xiàn)出了更優(yōu)的預(yù)測精度[20],但易出現(xiàn)異常值,受觀測角度影響大。在此基礎(chǔ)上,Jiang等[21]提出了完全基于近紅外和紅光波段的EVI2。除此之外,減少飽和性影響的植被指數(shù)還有寬動態(tài)植被指數(shù)(Wide dynamic range vegetation index, WDRVI)、葉綠素植被指數(shù)(Chlorophyll index, CI)、歸一化紅邊差值指數(shù)(Normalized difference red edge index,NDRE)等[17,22-23]。

針對植被葉片色素含量和植被生產(chǎn)力評估,Gitelson等[8]提出了三波段模型,并基于此得到綠波段葉綠素指數(shù)(Green chlorophyll index, CIgreen)和紅邊葉綠素指數(shù)(Rededge chlorophyll index, CIrededge),其形式簡單,且算法系數(shù)不受植被品種所影響。另外,陸地葉綠素指數(shù)(Terrestrial chlorophyll index,MTCI)為MERIS傳感器葉綠素產(chǎn)品的核心算法[24],在凈初級生產(chǎn)力的估算方面也有著一定應(yīng)用[25]。

在本研究中,將對比上述8個植被指數(shù)在不同氮水平下水稻葉片光合速率反演中的表現(xiàn),相應(yīng)植被指數(shù)計(jì)算公式見表1。其中,近紅外波段(ρNIR)取780~800 nm反射率均值,紅波段(ρred)取660~680 nm 反射率均值,綠波段(ρgreen)取540~560 nm反射率均值,紅邊波段(ρrededge)取710~730 nm反射率均值,ρ750與ρ705分別表示750與705 nm處光譜反射率,WDRVI中常數(shù)值α取0.2。

1.5 建模與精度驗(yàn)證

結(jié)合式(2),在葉片尺度上,分別以葉綠素相關(guān)植被指數(shù)VIchl-related和葉綠素相關(guān)植被指數(shù)與入射光合有效輻射的乘積VIchl-related×PARin為自變量對凈光合速率Pn進(jìn)行回歸建模,并比較引入PARin前后模型預(yù)測精度。回歸建模過程使用Matlab R2014a軟件進(jìn)行。

為了避免品種與施氮水平的影響,將數(shù)據(jù)按品種與氮水平進(jìn)行分組討論。對于各組數(shù)據(jù),隨機(jī)選取2/3作為建模集,余下1/3作為驗(yàn)證集。采用決定系數(shù)R2、均方根誤差(Root mean squared error,RMSE)與噪聲當(dāng)量(Noise equivalent, NE)NEΔPn分別作為模型精度與敏感度的評判標(biāo)準(zhǔn),NEΔPn計(jì)算公式如下[26]:

表1 植被指數(shù)計(jì)算公式Table 1 Formula of vegetation index

NEΔPn=RMSE{VI×PARin VS.Pn}/ d(VI×PARin)/d(Pn)

(3)

式中:RMSE{VI×PARin VS.Pn}為回歸模型VI×PARin VS.Pn為均方根誤差,d(VI×PARin)/d(Pn)為VI×PARin VS.Pn的一階導(dǎo)數(shù)。NEΔPn值越小,模型靈敏度越高,效果越佳;反之模型靈敏度低,效果越差。

2 結(jié)果與分析

對2個品種的水稻,以N0、N180、N240分別代表純氮肥0、180和240 kg/hm2,對表1中的8種葉綠素相關(guān)植被指數(shù),以VIchl-related作為自變量,凈光合速率Pn作為因變量,分別建立線性和非線性模型,選取各組結(jié)果中決定系數(shù)R2最高的模型并記錄(表2)。

表2 植被指數(shù)反演水稻葉片凈光合速率結(jié)果Table 2 The inversion result of net photosynthetic rate based only on VIs of rice leaves

由表2可知,在N0區(qū)域中,植被指數(shù)NDRE、CI、CIgreen、CIrededge和MTCI模型決定系數(shù)R2均高于0.34,且表現(xiàn)為在0.01水平雙側(cè)顯著相關(guān);但在N180和N240區(qū)域中,8種植被指數(shù)擬合模型結(jié)果的R2均低于0.27。同時(shí),不同氮水平區(qū)域間模型精度差異明顯,并不存在單一植被指數(shù)可用于不同組別數(shù)據(jù)反演。因此,僅使用葉綠素相關(guān)植被指數(shù)來對葉片凈光合速率進(jìn)行反演是不可行的。究其原因,除了葉片葉綠素含量外,光合作用的強(qiáng)弱與環(huán)境入射光強(qiáng)密切相關(guān),對于同一片葉來說,入射光合有效輻射不同,光合速率隨之不同[27]。

考慮光合入射有效輻射,以VIchl-related×PARin作為自變量,凈光合速率Pn作為因變量進(jìn)行建模,同樣選取各組結(jié)果中決定系數(shù)R2最高的模型,并記錄RMSE(表3和4)。

由表3和4可知,對于品種珞優(yōu)9348,葉綠素植被指數(shù)CI都給出了最佳的擬合結(jié)果,R2均達(dá)到了0.87以上,RMSE不超過1.13 μmol/(m2·s);驗(yàn)證集中,R2高于0.75,RMSE不超過1.63 μmol/(m2·s)。除了CI外, NDRE、CIgreen、CIrededge同樣給出了不錯的擬合結(jié)果。對于品種豐兩優(yōu)4號,N0區(qū)域中,NDVI表現(xiàn)最優(yōu),訓(xùn)練集和驗(yàn)證集中R2分別為0.77和0.80,RMSE分別為1.01和2.08 μmol/(m2·s);其次為CI。N180區(qū)域中,NDRE、WDRVI、CI和CIrededge擬合效果最好,R2均為0.83,RMSE不超過1.40 μmol/(m2·s),驗(yàn)證集中R2高于0.78,RMSE不超過1.81 μmol/(m2·s)。N240區(qū)域中,NDVI和CI擬合效果最好,R2均為0.88,RMSE不超過1.85 μmol/(m2·s),驗(yàn)證集中,R2均為0.80,RMSE不超過1.75 μmol/(m2·s)。

表3 不同氮含量區(qū)域珞優(yōu)9348葉片凈光合速率反演模型Table 3 Net photosynthetic rate retrieval models of Luoyou 9348 leaves in different nitrogen area

表4 不同氮含量區(qū)域豐兩優(yōu)4號葉片凈光合速率反演模型Table 4 Net photosynthetic rate retrieval models of Fengliangyou 4 leaves in different nitrogen area

比較不同氮水平區(qū)域結(jié)果,可以發(fā)現(xiàn)N180和N240區(qū)域的擬合結(jié)果都優(yōu)于N0區(qū)域,且該現(xiàn)象在豐兩優(yōu)4號中尤為明顯,其N0區(qū)域中,R2均不超過0.80,而驗(yàn)證集RMSE多為2~4 μmol/(m2·s),明顯高于N180與N240區(qū)域,說明氮肥脅迫會對反演模型精度造成影響。

從整體上看,除了WDRVI外,各VIchl-related×PARin用以預(yù)測水稻葉片凈光合速率都得到了穩(wěn)定的結(jié)果,訓(xùn)練集R2多為0.70~0.90,驗(yàn)證集R2為0.60~0.89,均表現(xiàn)為在0.01水平雙側(cè)顯著相關(guān),RMSE多低于2 μmol/(m2·s),滿足反演精度要求。其中,以植被指數(shù)CI給出的結(jié)果最優(yōu)(圖2),其次為NDVI和NDRE。

a)珞優(yōu)9348-N0區(qū)域 Luoyou 9348-N0 Area; b)珞優(yōu)9348-N180區(qū)域 Luoyou 9348-N180 Area; c)珞優(yōu)9348-N240區(qū)域 Luoyou 9348-N240 Area; d)豐兩優(yōu)4號-N0區(qū)域 Fengliangyou-N0 Area; e)豐兩優(yōu)4號-N180區(qū)域 Fengliangyou 4-N180 Area; f)豐兩優(yōu)4號-N240區(qū)域 Fengliangyou 4-N240 Area 下同。The same below。
圖2 基于葉綠素植被指數(shù)CI的水稻葉片凈光合速率Pn反演結(jié)果
Fig.2 The relationship between CI×PARin and net photosynthetic rate of rice leaves

為比較基于植被指數(shù)CI與NDRE、NDVI模型的靈敏度,作出其NEΔPn曲線如圖3所示。對于線性模型來說,其敏感度NEΔPn為一常量;在大部分情況下,NEΔPn值會隨著凈光合速率Pn的升高而升高,表現(xiàn)為預(yù)測靈敏度降低。對比上述3種植被指數(shù)模型的NEΔPn,多數(shù)情況下,當(dāng)凈光合速率Pn達(dá)到一定數(shù)值后,NDVI模型預(yù)測敏感度將低于CI和NDRE,即NDVI會在一定程度上達(dá)到飽和,但該飽和程度尚不致其失效。對比CI與NDRE可知,大多數(shù)同等條件下,植被指數(shù)CI都表現(xiàn)的更加靈敏,這也是表3和4中基于CI的模型決定系數(shù)R2多優(yōu)于NDRE的原因之一。綜合上述分析結(jié)果,基于植被指數(shù)CI的葉片凈光合速率模型結(jié)果最優(yōu),接下來的討論分析將圍繞CI展開。

圖3 植被指數(shù)CI, NDRE, NDVI噪聲當(dāng)量NEΔPn曲線圖
Fig.3 The NEΔPnof models CI, NDRE, NDVI×PARin VS.Pnof rice leaves

從Pn∝VIchl-related×PARin模型出發(fā),根據(jù)Ciganda等[7]的結(jié)論,綠色植物葉片的葉綠素可由CIrededge計(jì)算得到:Chl=(37.904+1 353.7×CIrededge)/1 000。為了探究水稻葉片光能利用能力與葉綠素含量間的關(guān)系,以播種后天數(shù)DOY為橫坐標(biāo),作出兩品種水稻Pn/PARin與葉片葉綠素含量間的關(guān)系折線圖(圖4),考慮到單位和尺度上的不同,首先需要對數(shù)據(jù)進(jìn)行歸一化處理。

圖4 凈光合速率/入射光合有效輻射與葉片葉綠素含量間變化關(guān)系折線圖
Fig.4 Line chart betweenPn/PARin and rice leaves chlorophyll content

由圖4可知,2個品種的水稻,不論所施氮肥濃度高低,Pn/PARin與葉片葉綠素含量在變化趨勢上都具有明顯一致性,即葉綠素含量的升高或降低,與水稻葉片對光能利用能力的高低一致,這與Gitelson等[3]和Peng等[28]的研究結(jié)論一致。基于上述分析和光能利用率模型公式(1),在葉片尺度上,用植被指數(shù)CI作為fAPAR的替代,可以得光能利用率LUE=Pn/(CI×PARin)。針對訓(xùn)練集數(shù)據(jù),選取CI分別計(jì)算珞優(yōu)9348和豐兩優(yōu)4號水稻葉片的線性模型,同時(shí),為了比較低氮和高氮的區(qū)別,增加高氮組別N180+240來對應(yīng)低氮組N0,分別記錄斜率、截距、R2和RMSE,其中,斜率即代表了葉片對光能的利用能力。

由表5可知,對于2個水稻品種,N0區(qū)域斜率明顯低于其他區(qū)域,而N180、N240與N180+240彼此之間非常近似。即低氮條件下,2個品種的水稻葉片對光能的利用能力有著明顯受限,氮肥濃度升高,光能利用率上升;當(dāng)?shù)屎? 180 kg/hm2后,水稻葉片光能利用能力變化不明顯。因此,在對不同施氮水平的水稻凈光合速率進(jìn)行反演時(shí),可考慮直接按低氮、高氮濃度進(jìn)行分組,以簡化計(jì)算。

表5 不同氮水平下的水稻葉片凈光合速率線性模型Table 5 Linear models of rice in different nitrogen content area

品種間橫向?qū)Ρ龋P托甭氏嘟2粎^(qū)分品種,僅區(qū)分高低氮濃度,用建模集數(shù)據(jù)對水稻葉片凈光合速率進(jìn)行建模,用驗(yàn)證集數(shù)據(jù)進(jìn)行檢驗(yàn),所得結(jié)果如圖5和6所示。由圖5可知,N180+240區(qū)域的擬合結(jié)果R2為0.89,優(yōu)于N0區(qū)域的0.80,這與前述區(qū)分品種反演所得結(jié)果一致。與表5相比,不區(qū)分品種的條件下,模型仍保持了較高的預(yù)測精度, RMSE均不超過1.37 μmol/(m2·s),圖6說明了區(qū)分高低氮肥含量條件下,凈光合速率的實(shí)測值與預(yù)測值之間存在有良好的線性關(guān)系,均勻分布在 1∶1 線附近,R2均高于0.86,RMSE不超過1.54 μmol/(m2·s)。因此,基于葉綠素相關(guān)植被指數(shù)CI的葉片模型可以很好地反演不同品種水稻葉片的凈光合速率,即對水稻品種不敏感。

圖5 不同氮水平區(qū)域水稻葉片凈光合速率最佳模型
Fig.5 BestPnretrieval models of rice leaves in different nitrogen content area

圖6 不同氮水平區(qū)域凈光合速率實(shí)測值VS.預(yù)測值
Fig.6 MeasuredPnVS. predictedPnin different nitrogen content area

3 結(jié) 論

本研究以2個品種的水稻葉片光譜為基礎(chǔ),結(jié)合光合入射有效輻射數(shù)據(jù),分別使用8種葉綠素相關(guān)植被指數(shù)VIchl-related和其與入射光合有效輻射的乘積VIchl-related×PARin來反演水稻葉片凈光合速率,對不同氮肥水平條件下的數(shù)據(jù)進(jìn)行對比分析,得出如下結(jié)論:

1)光合作用的強(qiáng)弱與植被葉片所能截獲的光能密切相關(guān),因此,不考慮光合入射有效輻射而僅使用葉綠素相關(guān)植被指數(shù)來反演凈光合速率是不可行的。

2)對于不同氮肥條件的水稻葉片,使用葉綠素相關(guān)植被指數(shù)與光合入射有效輻射的乘積來反演其凈光合速率是可行的,但需區(qū)分高、低氮濃度;以基于葉綠素植被指數(shù)CI的模型結(jié)果最好,可滿足反演精度需求,且對水稻品種不敏感。

3)遙感技術(shù)在植被凈光合速率與生產(chǎn)力研究中具有獨(dú)特的優(yōu)勢,本研究方法預(yù)測精度較高,可以彌補(bǔ)傳統(tǒng)實(shí)地調(diào)查耗時(shí)費(fèi)力的劣勢,對作物生產(chǎn)力評估和估產(chǎn)具有潛在應(yīng)用價(jià)值。

猜你喜歡
水稻區(qū)域模型
一半模型
什么是海水稻
水稻種植60天就能收獲啦
軍事文摘(2021年22期)2021-11-26 00:43:51
重要模型『一線三等角』
一季水稻
文苑(2020年6期)2020-06-22 08:41:52
重尾非線性自回歸模型自加權(quán)M-估計(jì)的漸近分布
水稻花
文苑(2019年22期)2019-12-07 05:29:00
3D打印中的模型分割與打包
關(guān)于四色猜想
分區(qū)域
主站蜘蛛池模板: 国产网站免费看| 新SSS无码手机在线观看| 福利国产微拍广场一区视频在线| a毛片在线| 亚洲综合精品香蕉久久网| 欧美一区福利| 秋霞一区二区三区| 国产精品制服| 欧美高清三区| 国产99视频精品免费观看9e| 国产精品成人免费视频99| 国产好痛疼轻点好爽的视频| 91毛片网| 91福利免费| 99青青青精品视频在线| 亚洲精品日产AⅤ| 久久久精品无码一二三区| 欧美国产日产一区二区| 一级看片免费视频| 青青草原国产免费av观看| 色婷婷在线影院| 亚洲成人福利网站| 国产国产人成免费视频77777| 一级福利视频| 久久国产精品麻豆系列| 91九色视频网| 无码区日韩专区免费系列 | 日本久久免费| 国产亚洲精品自在线| 国产高清毛片| 国产午夜精品鲁丝片| 91原创视频在线| 国产欧美网站| 波多野结衣二区| 国产成人亚洲精品无码电影| 国产成人av一区二区三区| 久久综合亚洲鲁鲁九月天| 波多野结衣一区二区三区AV| 国产理论最新国产精品视频| 亚洲综合天堂网| 午夜爽爽视频| 亚洲免费成人网| 中文字幕人妻av一区二区| 国产sm重味一区二区三区| 91小视频在线播放| 91免费在线看| 伊人久久大线影院首页| 毛片在线看网站| 国产精品一线天| 欧美天堂久久| 精品无码国产一区二区三区AV| 国产成人精品高清在线| 国产污视频在线观看| 黄色一级视频欧美| 久久综合色播五月男人的天堂| 免费毛片网站在线观看| 麻豆国产在线不卡一区二区| 国产成人调教在线视频| 亚洲狼网站狼狼鲁亚洲下载| 爱做久久久久久| 国产美女主播一级成人毛片| 午夜视频www| 欧美激情视频一区| 国产精品福利尤物youwu| 亚洲欧美日韩高清综合678| 青青草原国产av福利网站| 国产福利一区视频| 亚洲成网777777国产精品| 亚洲AV人人澡人人双人| 午夜国产精品视频| 欧美劲爆第一页| 国产91高跟丝袜| 亚洲永久色| 欧美成人日韩| 日本AⅤ精品一区二区三区日| 色哟哟色院91精品网站| 青青国产视频| 99国产精品国产高清一区二区| 日本国产精品一区久久久| 国产高清无码麻豆精品| 一本大道无码日韩精品影视| 成人国产精品一级毛片天堂 |