?

隴東黃土高原旱區短時強降水的時空分布特征及地形影響研究

2023-01-14 05:20楊麗杰曹彥超劉維成徐麗麗張洪芬孫子茱
干旱氣象 2022年6期
關鍵詞:旱區隴東黃土高原

楊麗杰,曹彥超,劉維成,徐麗麗,張洪芬,孫子茱

(1.甘肅省慶陽市氣象局,甘肅 慶陽 745000;2.蘭州中心氣象臺,甘肅 蘭州 730020)

引 言

我國短時強降水通常是指1 h 降水量大于等于20.0 mm 或3 h降水量大于等于50.0 mm 的強降水事件[1],它由中小尺度天氣系統直接驅動[2],其間常伴有雷暴大風、冰雹等天氣,具有突發性強、落區分散、發展迅猛等特點,極易誘發山洪、滑坡、泥石流等災害,造成嚴重人員傷亡和財產損失。伴隨著全球氣候變暖,我國多地短時強降水事件呈增多趨勢[3-5],給防災減災工作帶來巨大挑戰,也成為預報業務中的重點和難點問題之一。

我國幅員遼闊、地形地貌復雜,短時強降水時空分布極不均勻??臻g上,以遼東半島南部—燕山—陰山—河套—兩廣地區一線為界,分界線的東部、南部為短時強降水高發區[6],且區域性特征明顯,尤其在山區、丘陵[7-8]、城區[9-10]等較小空間范圍內更為凸顯[11-13]。時間上,短時強降水的月、候變化總體呈單峰型,日變化則表現出單峰型、雙峰型、多峰型和持續活躍型等多種特征[14]。這些研究為短時強降水預報預警打下了堅實基礎,但結論的普適性有待驗證,因此還需開展更多特定區域短時強降水的時空分布特征研究。

作為影響降水分布的重要因素[15-17],地形對降水影響具有明顯的“空間尺度效應”,大尺度地形可以直接影響大氣環流,中尺度地形主要通過渦旋運動和背風波影響降水,而小尺度地形主要通過迎風坡影響降水。此外,地形對層云和對流云降水的影響亦有不同,前者地形的動力作用更重要,后者地形的熱力作用更明顯[18]??梢?,地形對降水的影響極其復雜,準確描述地形對降水的影響機制,是改進模式預報性能、提升預報預警能力的關鍵。然而,受當前實況觀測數據分辨率低、數值模式對地形刻畫不夠精準[19]、地形對天氣系統影響的非線性[20]等客觀因素限制,地形對降水影響機制尚存在許多亟待探討的問題。

隴東黃土高原地處西北半干旱半濕潤氣候過渡帶,地勢西北高東南低,東西北三面環山,形成開口向南的“簸箕”狀地形[21]。該區域地表植被稀疏、土質疏松、溝壑縱橫,在短時強降水沖擊下下墊面承災能力較差。如2021年8月18日23:00(北京時,下同)至19日01:00,慶陽市鎮原縣閆溝站連續3 h出現短時強降水,累積雨量達170.0 mm,引發暴雨洪澇災害,直接經濟損失達968.35 萬元人民幣。因此,加強該地區短時強降水研究有著迫切的現實需求。隨著遙感技術的快速發展,高分辨率的數字高程模型(digital elevation model,DEM)數據得到廣泛應用,利用GIS技術可從中提取坡度、坡向等地形因子,為短時強降水影響機制分析提供了重要支撐[22-23]。本文利用2013—2020年隴東黃土高原半干旱半濕潤區分布較為稠密的302個區域自動氣象觀測站逐時降水數據和DEM 數據,在短時強降水精細化分布特征分析基礎上,初步探討地理、地形因子與短時強降水之間的關系,并結合常規氣象觀測資料和再分析數據,剖析地形對2021年8月18—19日極端短時強降水的影響機制,以期提升當地短時強降水的預報預警能力。

1 資料與方法

隴東黃土高原旱區位于甘肅省最東端,行政區域包括慶陽市和平涼市的崆峒區、靈臺縣、崇信縣、涇川縣和華亭縣共13 個縣(區)[21],西依六盤山,東接子午嶺,北起白于山分水嶺,南至涇河,海拔高度882~2737 m,地勢起伏大[圖1(a)]。其中,涇河及其支流在境內構成“樹枝”狀河網,切割出塬、殘塬、峁梁丘陵、掌地、丘陵和侵蝕堆積區6 種地貌單元[24],地形特征極其復雜。該區域處在季風邊緣地帶,年平均降水量300~600 mm,夏季降水頻發。除東部的子午嶺林區和西南部的六盤山林區外,地表植被稀疏[圖1(b)],土壤蓄水能力差,夏季洪澇災害頻發。

圖1 隴東黃土高原旱區海拔高度、河網及氣象站點(a)和2013—2020年生長季平均歸一化植被指數(b)空間分布Fig.1 The spatial distribution of elevation,river network,meteorological stations (a) and average normalized difference vegetation index (NDVI) in growing season from 2013 to 2020 (b) in the Loess Plateau arid region of eastern Gansu

所用資料包括氣象大數據云平臺提供的2013—2020年隴東黃土高原旱區302個地面區域自動氣象觀測站4—10月逐小時降水數據、地理空間數據云(https://www.gscloud.cn)提供的航天飛機雷達地形測繪任務(shuttle radar topography mission,SRTM)DEM數據(分辨率為90 m×90 m)、甘肅省慶陽市氣象局歸檔的常規氣象觀測數據和歐洲中期天氣預報中心(https://cds.climate.copernicus.eu)的第五代ERA5逐小時再分析數據(分辨率為0.25°×0.25°)。隴東黃土高原旱區包括的市、縣(區)行政邊界是基于國家測繪地理信息局標準地圖服務網站下載的審圖號為GS(2019)3082號的標準地圖制作,底圖無修改。

本文將雨強大于等于20.0 mm·h-1作為短時強降水事件的篩選條件,并根據甘肅省短時強降水等級標準(20.0~29.9、30.0~49.9、50.0 mm·h-1及以上)[25]進行分級,分別統計不同時段、不同等級短時強降水日數、次數和極值的時空分布。短時強降水日數定義為統計時段內研究區域出現短時強降水的日數,將24 h 內研究區域至少1 站出現1 次及以上的短時強降水事件記為1個短時強降水日數。單站短時強降水次數定義為統計時段內某站點出現短時強降水事件的時次;統計時段內研究區域各單站短時強降水次數之和定義為研究區短時強降水次數。短時強降水極值定義為統計時段內所有短時強降水事件小時雨量的最大值。經統計,2013—2020年隴東黃土高原旱區共出現2115次短時強降水,通過氣候極值、內部一致性、時間一致性和空間一致性[26]等對這些樣本進行質量控制,得到信度較高的1653個樣本,占比78.2%。

通過交叉驗證法和t檢驗[27],選取精度最高(表略)的反距離權重插值法獲取各要素的空間分布。

采用傅抱璞[28]提出的基于爬坡機制的地形降水模型,地形抬升速度公式如下:

式中:vg(m·s-1)為地形抬升速度;v(m·s-1)為風速;c1、c2分別為坡度、坡向影響系數;α(°)為地形坡度;θ(°)為盛行風向;β(°)為地形坡向。當c1>0(即0°<α≤45°)時,坡度對降水起正反饋,且α=45°時影響最大,反之起負反饋;當c2>0(即∣θ-β∣<90°)時為迎風坡,對降水起正反饋,反之為背風坡,對降水起負反饋。

利用ArcGIS 軟件,對DEM 數據進行空間重采樣,并提取各氣象站點的坡向、坡度、海拔高度3 個地形因子。統計發現,夏季環縣、鎮原及崆峒、華亭西部部分鄉鎮的主導風向為西西北風,其余地區為東南風或偏南風。根據上述公式,計算各站點的坡度和坡向影響系數,并同經度、緯度、海拔高度因子作為自變量,通過主成分分析和線性回歸分析,探討地形的動力抬升作用與短時強降水之間的關系。

2 短時強降水時間分布特征

2.1 年際變化

2013—2020年,隴東黃土高原旱區短時強降水共出現220 d、1653 次,短時強降水發生日數和次數占比存在明顯的年際變化,且趨勢基本同步[圖2(a)],2018年占比最大,分別為20.5%和23.2%,這可能是副熱帶高空急流和西太平洋副熱帶高壓(簡稱“西太副高”)明顯偏北引起的水汽輸送偏強造成的[29];2014年受強厄爾尼諾事件影響,西太副高偏西、偏南,水汽向北輸送困難[30],短時強降水日數和次數占比最小,分別為6.4%和5.1%。另外,2013、2016、2018年短時強降水占比高于平均值,為多發年,其余年份為少發年。結合逐年雨強箱線圖[圖2(b)]發現,近8 a 短時強降水雨強的25%分位數接近22.0 mm·h-1,除2015、2016、2020年外75%分位數均未超過31.0 mm·h-1,即雨強主要集中在22.0~31.0 mm·h-1。對比發現,雨強與日數和次數占比之間無明顯對應關系,短時強降水最多的2018年雨強并非最強。整體而言,2015、2016、2020年雨強較大,但2013、2016、2019年降水的極端性更強。

圖2 2013—2020年隴東黃土高原旱區短時強降水日數、次數占比(a)和雨強(b)的年際變化Fig.2 The inter-annual changes of proportion of days,times (a) and rainfall intensity (b) of short-term heavy rainfall in the Loess Plateau arid region of eastern Gansu from 2013 to 2020

2.2 月變化

統計發現,隴東黃土高原旱區短時強降水在4—10月均有發生,但主要集中在6—8月??紤]到樣本數量及季節背景,在分析雨強月變化時將4—5月、9—10月分別作為整體與6、7、8月進行對比。

從圖3看出,隴東黃土高原旱區短時強降水月變化呈單峰型特征,并表現出急增驟減、盛夏活躍強盛的特點。4—5月我國季風雨帶位于長江以南,隴東黃土高原區短時強降水少且增加緩慢,雨強最小,主要為22.0~26.0 mm·h-1。主汛期6—8月短時強降水共出現195 d、1597 次,占比分別為88.6%和96.6%,為全年多發時段,且雨強大、極端性強。其中,7月發生日數最多為79 d、占比35.9%,極端性最強,雨強主要在22.0~32.0 mm·h-1之間;而8月次數最多為775 次、占比46.9%,雨強主要在22.0~33.0 mm·h-1之間,較7月略強。這種變化趨勢可能與東亞季風的強度變化、西太副高的西伸北抬等有關。8月以后,隨著西太副高東退,短時強降水迅速減少,強度也隨之減弱。

圖3 2013—2020年隴東黃土高原旱區短時強降水日數、次數占比(a)和雨強(b)月變化Fig.3 The monthly changes of proportion of days,times (a) and rainfall intensity (b) of short-term heavy rainfall in the Loess Plateau arid region of eastern Gansu from 2013 to 2020

2.3 日變化

從圖4看出,隴東黃土高原旱區的短時強降水在任何時次都可能出現,其次數占比的日變化呈多峰型特征,傍晚至前半夜(17:00至次日00:00)最為活躍,而04:00—12:00 出現較少。其中,活躍時段內短時強降水共出現938 次、占比56.8%,并分別于00:00、17:00和22:00出現峰值,對應次數(占比)分別為143(8.7%)、125(7.6%)、130(7.9%)次,且00:00和22:00雨強的箱體范圍明顯大于其他時次,說明離散性更強。值得注意的是,在短時強降水相對不活躍的時段中,06:00 的次數占比和強度均較大,且極值(2016年8月15日06:00 環縣八珠雨強76.7 mm·h-1)僅次于22:00 和23:00,也是日常業務中需關注的時段。

圖4 2013—2020年隴東黃土高原旱區短時強降水次數占比(a)和雨強(b)逐時變化Fig.4 The hourly variation of proportion of times (a) and rainfall intensity (b) of short-term heavy rainfall in the Loess Plateau arid region of eastern Gansu from 2013 to 2020

3 短時強降水空間分布特征

2013—2020年隴東黃土高原旱區短時強降水空間分布極不均勻,總體表現為東南多、西北少,全區普遍累積發生5~7 次,有3 個相對集中的高發區,分別在崆峒東部—靈臺、鎮原東北部—慶城北部—華池中部—合水西部和寧縣中南部—正寧北部[圖5(a)]。其中,崆峒東部—靈臺短時強降水最為活躍,四十里鋪、東王溝村等站點最多超過16次,鎮原東北部—慶城北部—華池中部—合水西部局地多達13~15 次,而范圍相對最小的寧縣中南部—正寧北部高發區,亦有個別站點超過15次。從短時強降水的小時雨量極值分布[圖5(b)]來看,隴東黃土高原旱區短時強降水的小時雨量極值呈東北和西南大、中間小的空間特征,極值普遍在31.0~40.0 mm,50.0 mm 以上極值集中在慶城東部—合水西部,而環縣東北部、平涼中東部亦有部分站點極值突破60.0 mm,但極值大值區的范圍整體較次數明顯偏小。

進一步對比不同等級的短時強降水次數空間分布[圖5(c)和圖5(d)],發現隨著雨強增大,短時強降水次數驟減,20.0~29.9 mm·h-1的降水次數與總次數的空間分布特征最為相似,普遍發生4~5 次;30.0~49.9 mm·h-1的短時強降水普遍發生1~2 次,高發區集中在慶陽市中部慶城—鎮原交界及華池東北部等地,達6次,而環縣西北部、慶城東北部、鎮原中部、崆峒中部及靈臺東部也有部分站點超過5次;50.0 mm·h-1及以上短時強降水的局地性進一步凸顯,8 a間31站出現的34次50 mm·h-1及以上短時強降水中,僅有平涼東部和慶陽市慶城、合水交界處和環縣西北部的3站出現2次,其余站均出現1次。

圖5 2013—2020年隴東黃土高原旱區短時強降水次數(a)、小時雨量極值(b)以及20.0~29.9(c)、30.0~49.9和50.0 mm·h-1及以上(d)降水次數的空間分布(紅色方框為短時強降水多發區域,綠色方框為掌地地貌)Fig.5 The spatial distributions of times (a) and hourly precipitation extremum (b) of short-term heavy rainfall and times of short-term heavy rainfall with rainfall intensity of 20.0-29.9 (c),30.0-49.9 and greater than or equal to 50.0 mm·h-1 (d) in the Loess Plateau arid region of eastern Gansu from 2013 to 2020(The red boxes are the high incidence areas of short-term heavy rainfall,and the green box is the palm landform)

結合地形和河網分布,發現短時強降水的高發區主要分布在河谷地勢升高區或河流交匯處,如六盤山東側西北—東南向的涇河河谷(紅框A)、西北轉南方向的環江—馬蓮河河谷(紅框B)、子午嶺西側上坡帶的東北—西南向四郎河河谷(紅框C),它們共同構成了開口向南的喇叭口地形。此外,30.0 mm·h-1以上的短時強降水在環縣西北部的黃土掌地區(綠框)也有分布。

與甘肅短時強降水最多的隴南市[4]對比發現,兩地的短時強降水時間變化特征較為相似,夜雨特征明顯,但隴東黃土高原旱區短時強降水的雨強偏弱、發生次數相對較少、局地性更強。這除了與隴南地處溫帶、亞熱帶半濕潤氣候區[31],水汽輸送條件更好外,還與其高山峻嶺與峽谷、盆地相間的復雜地形有一定關系。

4 地形、地理因子與短時強降水之間的關系

采用線性回歸分析方法探究地形、地理因子對隴東黃土高原旱區短時強降水的影響。由于地形、地理因子之間彼此存在共線性問題,首先利用主成分分析法進行降維,發現排名靠前的3個主成分是:由經度和緯度貢獻最大的反映地理位置的因子P1;由坡度影響系數、緯度和海拔高度貢獻最大的反映地形強迫抬升的因子P2;主要由坡向影響系數貢獻的迎風坡因子P3(表略)。這3個主成分的方差膨脹因子(variance inflation factor,VIF)等于1,共線性較低,作為自變量是合理的。然后,分別以短時強降水的次數和小時雨量極值作為因變量進行線性回歸分析,并進行F檢驗,其結果僅短時強降水次數的回歸方程通過顯著性檢驗,表明地形、地理因子對短時強降水的小時雨量極值分布無顯著影響,而對于短時強降水次數分布僅地理位置的影響顯著(通過α=0.01的顯著性檢驗)(表1)。這與張寧等[10]開展的蘭州市短時強降水研究結論不同,可見地形對短時強降水的影響復雜,不同地區的影響機制差異較大。

表1 隴東黃土高原旱區短時強降水次數主成分的線性回歸分析Tab.1 The linear regression analysis of principal components of short-term heavy rainfall times in the Loess Plateau arid region of eastern Gansu

通常,海拔高度和坡向、坡度主要通過迎風坡的強迫抬升影響降水分布[28]。然而,隴東黃土高原旱區中南部地形坡度小,動力抬升有限,而在坡度較大的東北部,迎風坡與背風坡交替分布且范圍相當,對降水的反饋互相抵消,致使地形對降水無明顯動力抬升作用。對比坡向和坡度影響系數空間分布(圖6)發現,僅子午嶺沿線及環縣東部等地坡向和坡度影響系數對降水有正反饋。因此,地形強迫抬升并非是隴東黃土高原旱區短時強降水的主要影響機制。

圖6 隴東黃土高原旱區坡向(a)和坡度(b)影響系數的空間分布Fig.6 The spatial distributions of aspect (a) and slope (b) influence coefficients in the Loess Plateau arid region of eastern Gansu

5 地形對短時強降水的影響機制

地形通過強迫抬升作用影響短時強降水外,還可以通過邊界層輻合引起的上升運動、陸面空間差異導致的局地熱動力作用不均勻等[32]對降水產生影響。下面結合一次極端短時強降水事件剖析隴東黃土高原旱區地形的影響機制。

2021年8月18—19日,隴東黃土高原旱區發生極端短時強降水事件,共有8 站過程累積雨量突破100.0 mm,并伴有35 次短時強降水;降水最多的鎮原閆溝站24 h 累積雨量達195.6 mm,其中18日23:00 至19日01:00 的3 h 累積雨量為170.0 mm(占比86.9%),突破建站以來的歷史極值,19日00:00小時雨量為86.3 mm,為建站以來歷史次強的短時強降水[圖7(a)]。結合地形和河流分布,發現此次短時強降水集中在空間尺度較小的河谷喇叭口地形區[圖7(b)],由中小尺度天氣系統產生。

圖7 2021年8月18—19日極端強降水期間短時強降水次數和雨量極值逐時變化(a)、次數空間分布(b)和18日21:00(c)、19日03:00(d)地面風場(風矢,單位:m·s-1)及輻合線(藍色線)(彩色填充區為海拔高度,單位:m)以及18日16:00(e)、19日00:00(f)沿35.6°N的假相當位溫(彩色填充區,單位:K)緯向垂直剖面(黑色陰影為地形)Fig.7 The hourly evolutions of times and precipitation extremum (a) and the spatial distribution of times (b) of short-term heavy rainfall during the extremely strong precipitation process from 18 to 19 August,the surface wind field (wind vectors,Unit:m·s-1)and convergence lines (blue lines)(the color shaded for elevation,Unit:m) at 21:00 BST 18 August (c) and 03:00 BST 19 August (d),and the zonal vertical section of pseudo-equivalent potential temperature (the color shaded,Unit:K)(the black shaded for terrain) along 35.6°N at 16:00 BST 18 August (e) and 00:00 BST 19 August (f) 2021

此次過程發生在有利的天氣背景下,18日20:00(圖略),隴東黃土高原處于500 hPa 西太副高邊緣強盛的西南氣流控制下,其上空200 hPa高空急流分流區的輻散抽吸作用和700 hPa 切變線為降水提供動力抬升條件。距離最近的平涼探空站觀測的對流有效位能為1362.9 J·kg-1、700 hPa 比濕為12 g·kg-1、溫度露點差為0.8 ℃、0 ℃層高度和抬升凝結高度分別為5675.0 m 和814.9 m,且中高層有干冷空氣侵入,環境場處于高能、高濕的對流不穩定狀態,暖云層深厚,有利于短時強降水發生。之后,西太副高東退,西風槽前分裂的冷空氣下滲,與低層切變線東側的暖濕氣流在河谷地區相遇,導致此次極端短時強降水事件的發生,地面輻合線是其主要的中尺度影響系統。

地面逐小時加密觀測數據顯示,18日21:00[圖7(c)]開始,在隴東黃土高原西部河網密布的喇叭口地形區和北部河流上游山坡分別有一條中尺度輻合線發展,并持續到19日03:00[圖7(d)],與短時強降水活躍時段一致。其間,輻合線A 略有東移,觸發的短時強降水主要集中在中部的茹河、蒲河河谷,而輻合線B則由緯向逐漸轉為經向分布,觸發的短時強降水主要分布在其前側的環江、元城川河谷[圖7(b)]。

為了厘清短時強降水前后河谷地形對溫濕條件和局地環流的影響,沿降水最強的鎮原閆溝站(107.5°E,35.6°N)做假相當位溫的緯向垂直剖面。18日白天,隴東黃土高原區處于西太副高控制下,天氣晴好。由于太陽輻射加熱不均,白天山坡較山谷升溫迅速,16:00[圖7(e)]六盤山(106.5°E)近地面800 hPa 附近假相當位溫超過360 K,較其東側茹河河谷(107.8°E)近地層850 hPa高出6 K,使得河谷中以西南或東南向的谷風為主,與沿地面低壓倒槽東北側進入河谷的風方向一致,故白天無輻合線活動;夜間,六盤山東側山坡較茹河河谷降溫迅速,假相當位溫梯度逐漸轉為與白天相反方向。19日00:00[圖7(f)],六盤山東側山坡近地層假相當位溫降至350 K 左右,比茹河河谷近地層356 K 低了6 K,使得山上的偏西風不斷下滑,形成與白天相反的近地層環流(山風),下滑的偏西氣流在河谷附近與倒槽東側的西南暖濕氣流交匯,形成地面輻合線,導致短時強降水增強。由此可見,下墊面熱力差異產生的山谷風環流誘發了地面中尺度輻合線并穩定維持,是地形影響此次極端短時強降水事件的主要機制。

6 結 論

(1)近8 a隴東黃土高原旱區短時強降水的雨強主要集中在22.0~31.0 mm·h-1,其日數、次數和強度具有明顯的年際變化,2013、2016、2018年相對多發;月變化呈單峰型分布,7月日數占比(35.9%)最多、極端性最強,8月次數占比(46.9%)最多、強度最強;任何時次都可能出現短時強降水,日變化呈多峰型特征,17:00 至次日00:00 最為活躍,次數占比56.8%,且強度和極端性最強,而06:00 的次數占比和強度驟升,在業務中需警惕。

(2)隴東黃土高原旱區短時強降水次數和小時雨量極值空間分布極不均勻,次數東南多、西北少,8 a 全區普遍為5~7 次,崆峒東部—靈臺最為頻發,其次為鎮原東北部—慶城北部—華池中部—合水西部,寧縣中南部—正寧北部也有小范圍的頻發區域;小時雨量極值則表現為東北和西南大、中間小,普遍為31.0~40.0 mm,50.0 mm以上的高值區集中在慶城東部—合水西部。

(3)隨著降水強度增大,短時強降水次數驟減。20.0~29.9 mm·h-1的短時強降水次數與總次數空間格局一致,整體也是西北少、東南多,8 a 普遍為4~5次;30.0 mm·h-1以上的普遍為1~2 次,其中30.0~49.9 mm·h-1的高發區在慶陽市中部,而50.0 mm·h-1及以上的高發區分散在崆峒的個別站點。結合地形及河網分布,發現隴東黃土高原旱區短時強降水主要集中在河谷喇叭口地形的地勢升高區或河流交匯處,而環縣西北部的掌地也是30.0 mm·h-1以上強降水高發區。

(4)地理和地形因子對隴東黃土高原旱區短時強降水的小時雨量極值分布無明顯影響,而對降水次數分布的顯著影響則主要由地理位置造成??傮w而言,地形的動力抬升對短時強降水影響不大,僅子午嶺沿線及環縣東南部等坡度較大的迎風坡對降水有一定的正反饋。

(5)下墊面熱力差異導致的山谷風環流及其誘發的地面中尺度輻合線是隴東黃土高原旱區喇叭口河谷地形區產生短時強降水的重要原因。

猜你喜歡
旱區隴東黃土高原
加氣對西北旱區膜下滴灌棉花生長與水分利用效率的影響
基于根系層水分狀態的旱區凈灌溉需水量模型構建和應用
隴東能源大數據中心
不同種植區隴東苜蓿營養價值的比較研究
寧南旱區溝壟集雨結合覆蓋對土壤水熱肥與馬鈴薯產量的影響
基于Brutsaert模型的黃土高原蒸發規律研究
行于大道——憶從黃土高原走向高峰的劉文西
隴東黃土高原軟棗獼猴桃栽培技術
淺談有機培肥措施對土壤肥力及作物生長的影響
91香蕉高清国产线观看免费-97夜夜澡人人爽人人喊a-99久久久无码国产精品9-国产亚洲日韩欧美综合