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

2015-2020 年山西省空氣質量時空變化特征及其驅動因素分析

2023-10-27 06:22:54何向陽張鳳英何立環馬廣文林蘭鈺方德昆
環境科學研究 2023年10期
關鍵詞:風速污染模型

何向陽,張鳳英,何立環,馬廣文,林蘭鈺,方德昆*

1.中國環境科學研究院環境健康風險評估與研究中心,北京 100012

2.中國環境監測總站,北京 100012

近年來,很多學者對我國空氣質量不同空間尺度下的時空特征開展了研究[1-3],多以社會經濟因素(人口聚集[4]、城市化水平[5]、城市綠化覆蓋率[6]等)和自然因素(風速[7]、氣溫[8]、降水[9]、地形[10]等)作為潛在影響因素,采用線性回歸模型[11]、地理加權回歸(GWR)模型[12]、地理探測器[13]和小波分析[14]等方法開展研究.回歸模型可以建立多影響因子與因變量之間的定量關系,Fang 等[11]發現,城市化對我國的空氣質量起到了負面作用,人口、城市化率、汽車密度、第二產業占比等因素均對空氣質量產生重要影響.小波變換可以對非平穩時間序列進行有效分析,能夠識別時間序列中的不同頻率成分并揭示其周期性變化,Tan 等[15]結合成渝地區2015-2021 年空氣質量監測數據和小波分析發現,短周期內空氣質量指數與氣溫和日照指數均呈明顯正相關,與溫度和降雨量均呈明顯負相關.區域空氣質量變化特征和驅動因素研究可為環境管理提供科學依據和支撐.

山西省是我國重要的煤炭生產基地和能源消耗大省,以煤炭、冶金、機械等為主要產業,是近年來我國空氣質量較差的地區之一,太原、臨汾、運城等城市在168 個城市空氣質量排名中多次位列倒數20 名內[16].山西省空氣污染受季節、氣象條件和區域交通影響,且城市和農村地區的污染特征存在明顯差異[17].研究[18]發現,山西省空氣污染主要來自工業和交通排放,且區域間傳輸和擴散對污染的影響較大.但這些研究多聚焦于汾渭平原(山西省南部)和京津冀及周邊地區覆蓋的山西省北部城市[19-21],或僅對單一污染物的時空特征開展分析[22],對山西省重污染地區(如臨汾市)或對重要污染物(O3)以及綜合分析全省空氣質量時空特征的研究較少.

因此,本文以2015-2020 年為研究期,基于山西省空氣質量監測結果,結合社會經濟、氣象等相關數據,使用地理加權回歸模型和小波分析等研究方法,開展空氣質量時空變化特征及驅動因素分析,以期為改善山西省空氣質量提供科學支撐.

1 材料與方法

1.1 研究區概況

山西省(34°34′N~40°44′N、110°14′E~114°33′E)位于我國華北地區(見圖1),共轄11 個地級市,轄區面積約15.67×104km2,地勢起伏不平,山區面積占80.1%.氣候類型屬溫帶大陸性季風氣候,四季分明、雨熱同期,從西北部到東南部由半干旱氣候過渡到半濕潤氣候[22-23].

圖1 山西省海拔及空氣質量監測點位Fig.1 Air quality monitoring points and altitude in Shanxi Province

1.2 數據來源和處理

空氣質量數據來源于國家生態環境監測網,包括山西省11 個地級市的六項污染物監測結果(SO2、NO2、CO、PM2.5和PM10濃度為24 h 平均值,O3濃度為日最大1 h 平均值)以及每日空氣質量指數(AQI)數據;社會經濟數據來源于《山西省統計年鑒-2021 年》[24],包括人口數量、人口密度、城鎮化率、GDP、人均GDP、第一產業生產值、第一產業占比、第二產業生產值、第二產業占比、民用汽車數量、工業用電量、城市綠化面積和城市綠化覆蓋率等;氣象數據中風速數據和平均溫度數據來源于美國國家海洋和大氣管理局(NOAA)下設的國家環境信息中心(NCEI)開源獲取的中國地面氣候資料日值數據集V3.0,降水量和平均相對濕度數據來源于中國地面氣候資料日值數據集V3.0.

1.3 分析方法

分析時段為2015 年1 月1 日-2020 年12 月31 日,空氣質量監測數據和氣象數據有效天數均為2 192 d.為滿足時間序列分析模型要求,在數據預處理階段通過Python 軟件中“sklearn”包的“Random Forest Regressor”函數補全空氣質量數據中的缺失值.空氣質量數據統計和可視化展示在R 4.2.1 和ArcGIS 10.8 軟件中完成;地理加權回歸(GWR)和多尺度地理加權回歸(MGWR)分析利用GWmodelS 軟件[25]完成;臨汾市和大同市PM2.5、O3濃度與氣象因子間的Pearson 相關性分析在R 4.2.1 軟件中完成;最后在MATLAB R2021b 軟件中利用小波分析探究臨汾市和大同市PM2.5、O3濃度與風速的周期性及相關關系.

1.3.1 地理加權回歸模型(GWR)

GWR 模型是普通線性回歸模型的拓展,引入數據的空間位置信息,進而對觀測點進行回歸分析求出對應的空間回歸系數[26-27].GWR 模型的一般表達式如下:

式中,yi為觀測點i的因變量值,(ui,vi) 表示觀測點i的位置坐標(如經緯度),β0(ui,vi) 為觀測點i的截距,xik(k=1,2,···,m)表示觀測點i第k個自變量值,βk(ui,vi)表示觀測點i第k個自變量的回歸參數,εi為誤差項[28].

1.3.2 多尺度地理加權回歸模型(MGWR)

GWR 與MGWR 均被認為是廣義相加模型,但與GWR 在相同的空間尺度進行操作不同,MGWR讓恒定帶寬(bw)假設變得更為靈活,允許因變量和不同自變量之間的條件關系在不同空間尺度上變化[29-31].MGWR 的計算公式:

式中,βbwk(ui,vi)表示觀測點i第k個自變量確定最優帶寬條件下的回歸參數.

1.3.3 小波分析

小波方差分析可判別時間序列的周期性變化規律,對于污染因子的時間序列x,小波方差檢驗可判斷序列的主周期[14],計算公式:

在空氣質量的影響因素分析中經常需要量化兩個非平穩信號之間的統計關系,為了量化兩個非平穩信號之間的關系,可以計算小波交叉譜和小波相干[32].引入另一個時間序列y(t) 后,小波交叉譜可由Wx,y(f,τ)=Wx(f,τ)(f,τ)得到,其中*表示復共軛,而小波相干則可以定義為各信號頻譜歸一化后的交叉譜[33],計算公式:

式中,Rx,y(f,τ)表示序列x和序列y的小波相干系數,<>表示時間和尺度上的平滑算子.通過式(4)計算可以確定Rx,y(f,τ)的范圍為[0,1].

2 結果與討論

2.1 山西省空氣質量時空特征

山西省2015-2020 年六項空氣污染物(SO2、NO2、O3、CO、PM2.5、PM10)濃度逐日變化趨勢如圖2所示,其中,2016-2017 年六項污染物濃度均呈上升趨勢;2017 年之后SO2、CO 和PM2.5濃度均呈下降趨勢,NO2和O3濃度在2017 年后均呈較強的季節性特征,未呈現明顯年際變化趨勢,PM10濃度僅在2020年后呈下降趨勢.從年周期上看,山西省六項污染物濃度日均值變化特征主要呈單峰型,除O3濃度的單峰出現在夏季外,其余污染物濃度峰值均出現在冬季,峰值過后濃度均開始下降.這是由于O3污染一般發生在夏季,而冬季采暖期SO2等污染物排放量增加,加之不利氣象條件,污染物擴散條件較差,因此在時間上O3污染與其余污染物造成的污染不同步.

圖2 2015-2020 年山西省六項污染物濃度逐日變化情況Fig.2 Daily concentration of six pollutants in Shanxi Province from 2015 to 2020

由圖3 可見,2015-2020 年,山西省11 個地級市六項污染物年均值的變化趨勢大致相同,其中SO2、CO 和PM2.5年均濃度均明顯下降,但各市NO2、O3和PM10年均濃度差值有明顯差異.相較于2015 年,各市SO2年均濃度變化范圍為-48~-11 μg/m3,其中太原市下降最多,從2015 年的65 μg/m3降至2020 年的17 μg/m3,低于GB 3095-2012《環境空氣質量標準》一級濃度限值(20 μg/m3)[34],這受益于近年來太原市多次開展秋冬季大氣污染綜合治理攻堅等專項行動以及全省SO2污染物減排工作的大力推進[35].CO 年均濃度變化范圍為-1.8~-0.2 mg/m3,呂梁市從2015 年CO 年均濃度最高(2.5 mg/m3)的城市,變為2020 年CO 年均濃度最低(0.7 mg/m3)的城市.PM2.5年均濃度變化范圍為-11~-1 μg/m3,其中太原市PM2.5年均濃度僅從2015 年的55 μg/m3降至2020 年的54 μg/m3.除朔州市、呂梁市、臨汾市外,其余城市近年來PM10年均濃度均下降,這可能是呂梁等市長期存在結構性污染(焦煤產業作為工業主體)問題所致[36].NO2年均濃度變化范圍為-3~15 μg/m3,其中長治市、晉城市、運城市和朔州市均降低,其余城市均上升.O3年均濃度變化范圍為-2~67 μg/m3,僅朔州市降低,其余城市均上升,其中晉城市O3年均濃度從2015 年的60 μg/m3升至127 μg/m3,增加了111.1%,這與當地電力源、焦化源、其他工業源和柴油交通源逐漸增加有關[18].

山西省11 個地級市AQI 值逐日特征(見圖4)表明,冬季各市空氣質量存在明顯空間差異,大同市和呂梁市的空氣質量明顯優于其他城市,而運城市和臨汾市的空氣質量較差.除大同市外,其余各市2016 年冬季空氣質量均較差,這可能是在采暖期的高污染排放情況下,伴隨極端氣候造成的不利氣象條件所致[37].隨著山西省秋冬季空氣污染治理攻堅行動的不斷推進以及近年來相對有利的氣象條件,2017 年后各市空氣質量在冬季有一定改善.在新型冠狀病毒感染疫情造成社會經濟活動減少的背景下,2020 年冬季各市空氣質量較2016 年同期有明顯改善[38].

結合山西省11 個地級市空氣質量等級分布(見圖5)來看,除臨汾市在2017-2019 年優良天數比例低于50%外,其余各市在2015-2020 年的空氣質量優良天數比例均大于50%,其中大同市每年的空氣質量優良天數都超過300 d,2016 年優良天數最多,達334 d,且無重度污染及嚴重污染天氣出現,這也是近年來大同市政府開展包括清潔取暖改造、燃煤鍋爐綜合治理等多項精準治氣措施帶來的成果[39].

圖5 2015-2020 年山西省11 個地級市空氣質量等級分布情況Fig.5 Air quality grade distribution of 11 cities in Shanxi Province from 2015 to 2020

2015-2020 年,臨汾市重度污染和嚴重污染天氣出現次數最多,分別達124 和41 d,其中2016 年臨汾市全年的重度污染天數達30 d,嚴重污染天數為14 d,與臨汾市的煤炭消費量巨大、散煤燃燒排放多、焦化企業排放超標以及集中供熱脫硫裝置數量不夠等有關[40].

綜合來看,近年來大同市空氣質量在山西省內相對較好,而臨汾市的空氣污染防治情況依舊嚴峻.

2.2 社會經濟影響(基于GWR 和MGWR 模型)

山西省近年來SO2濃度持續下降,但煤炭開發和利用使得冬季SO2污染形勢嚴峻,且各市NO2污染程度有所加重.據統計,2020 年全國SO2年均濃度達GB 3095-2012《環境空氣質量標準》二級標準的10個城市以及NO2年均濃度超過GB 3095-2012 二級標準的6 個城市中,山西省所轄城市占比分別為1/5和1/3[16].為探究社會經濟因子與SO2和NO2濃度的相關關系,運用GWR 和MGWR 模型對其相關關系的空間異質性進行研究.選擇2020 年各市的人口數量、人口密度、城鎮化率、GDP、人均GDP、第一產業生產值、第一產業占比、第二產業生產值、第二產業占比、民用汽車數量、工業用電量、城市綠化面積和城市綠化覆蓋率共13 個因子作為候選社會經濟因子.根據最小二乘(OLS)模型回歸結果中VIF 值(VIF 值<7.5)剔除多余因子,排除多重共線性問題,最終選擇人口密度(104人/km2)、城鎮化率(%)、人均GDP(元/人)、工業用電量(108kW?h)、城市綠化覆蓋率(%)五個因子.

OLS 模型和GWR 模型的分析結果如表1 所示,根據R2和調整R2的結果可以判斷,兩種污染物應用GWR 模型分析的擬合效果均比OLS 模型好,但修正的Akaike 信息準則(AICc)值升高說明GWR 模型的擬合優度較低.

表1 OLS 和GWR 模型分析結果Table 1 Analysis results of OLS and GWR models

為了進一步探究變量估計的尺度差異性,選擇MGWR 模型進行分析,通過迭代選擇最優帶寬.相較于GWR 模型,MGWR 模型AICc 值有小幅提高,但提升了R2值,可以更好地對因變量進行解釋(見表2).而調整R2的降低說明由于自由度改變造成變量的擬合優度略有下降[41].

表2 MGWR 模型分析結果Table 2 Analysis results of MGWR model

應用MGWR 模型分析SO2濃度與社會經濟因素的關系,R2達0.934,SO2濃度與人口密度呈顯著負相關,與城鎮化率呈顯著正相關.這可能是由于城鎮化率高的地區,其能源及產業結構附帶的SO2排放源較多,從而導致SO2濃度較高[42];而人口密度大的地區,對于SO2污染治理投入較大,可以有效降低SO2濃度[43].

結果顯示,NO2濃度與人口密度、城鎮化率和人均GDP 均呈正相關,與城市綠化覆蓋率呈顯著負相關.這是由于空氣中的NO2多源自城市內機動車尾氣排放,人口密度大、城鎮化率高及人均GDP 較高的城市其機動車數量越多,增加了NO2排放量[44].城市綠化覆蓋率高對空氣中氮氧化物的吸收降解產生較大作用,可以有效降低NO2濃度.NO2濃度與工業用電量呈負相關,這與近年來山西省移動源大氣污染排放中氮氧化物排放量較高有關,說明工業生產可能不是近年來山西省城市NO2污染的主要來源[45-46].

MGWR 模型結果中帶寬越大,代表空間異質性越低、空間關系越平穩[29].本研究中人均GDP 和城市綠化覆蓋率與SO2濃度相關關系的空間異質性較高,人口密度和人均GDP 與NO2濃度相關關系的空間異質性較高.模型統計結果中回歸系數的分布(見圖6)表明:人均GDP 與山西省南部城市(運城市、臨汾市和晉城市)SO2濃度呈負相關;城鎮綠化率對山西省北部城市(大同市、朔州市和忻州市)SO2濃度的影響較大,越靠南其影響越小;人均GDP 對山西省南部城市NO2濃度的影響較大,越靠北其影響越小;人口密度對山西省北部城市NO2濃度的影響較大,越靠南其影響越小.

圖6 MGWR 模型回歸系數分布Fig.6 Distribution of MGWR regression coefficient

2.3 氣象因素對空氣質量的影響

為探究氣象因子對山西省空氣質量變化的影響,選取山西省空氣質量最差(臨汾市)和最好(大同市)的城市,結合風速、平均氣溫、降水量和平均相對濕度作為氣象因子,探究近年來濃度明顯降低的PM2.5以及污染明顯加重的O3與氣象因素的相關關系.

Pearson 分析結果(見表3)顯示,臨汾市和大同市氣象因子與PM2.5濃度的相關性一致,除平均相對濕度外,其余氣象因子與PM2.5濃度均呈負相關,表明風速增強、氣溫升高和降雨增加均有利于緩解PM2.5污染.臨汾市平均氣溫與PM2.5濃度的相關性強于大同市,這與兩市所處緯度有關,臨汾市位于晉南地區,其日平均氣溫基本高于大同市,氣溫升高利于混合層高度的提升,加速了顆粒物擴散,從而降低顆粒物污染程度[47].兩市平均相對濕度與PM2.5濃度的正相關性均較弱.

表3 2015-2020 年臨汾市和大同市PM2.5、O3 濃度與氣象因子的Pearson 相關系數Table 3 Pearson correlation coefficient between PM2.5,O3 and meteorological factors in Linfen City and Datong City from 2015 to 2020

臨汾市和大同市平均氣溫與O3濃度均呈顯著正相關,O3作為二次污染物,其濃度隨溫度的升高而增加[48].其余氣象因子與O3濃度的相關性均較弱,且在兩市中存在明顯差異,臨汾市風速與O3濃度呈顯著正相關,而大同市風速對O3濃度無顯著影響;大同市降水量與O3濃度呈顯著正相關,而臨汾市降水量對O3濃度無顯著影響;臨汾市平均相對濕度與O3濃度呈顯著負相關,而大同市平均相對濕度對O3濃度無顯著影響.

風速影響著空氣污染物的稀釋、擴散和轉化過程,為進一步探究風速對山西省空氣質量的影響,利用小波分析方法研究臨汾市和大同市2015-2020年PM2.5和O3濃度與風速的周期性變化以及二者的關系.

由圖7 可見,臨汾市和大同市PM2.5濃度、O3濃度、風速均存在約1 a 左右的周期變化特征.兩市PM2.5濃度均存在約半年的長周期,但所處的時間范圍跨度不同.兩市PM2.5濃度均呈冬季高、夏季低的特征,且在冬季存在多個短周期.2017 年之后臨汾市PM2.5濃度的短周期頻次比大同市多,表明臨汾市PM2.5污染較嚴重、周期性特征更加復雜.兩市O3濃度在夏季和秋季較高,且在這兩季存在多個短周期,短周期的天數基本小于30 d,這是由于北方地區的氣象因素導致的,特別是太陽輻射和西南季風在季節性變化中起關鍵作用[49].整體上,大同市的風速強度高于臨汾市,且存在更多的短周期,這是因為大同市與內蒙古自治區相鄰,受內陸地區季風影響明顯.

圖7 臨汾市和大同市PM2.5 濃度、O3 濃度、風速的小波圖Fig.7 Wavelet diagram of PM2.5 concentration,O3 concentration and wind speed in Linfen City and Datong City

通過交叉小波和小波相干分析發現,臨汾市和大同市PM2.5和O3兩種污染物濃度與風速在不同周期中存在不同的相關關系(見圖8).從年周期來看,兩市PM2.5濃度變化均滯后于風速,這是因為PM2.5污染多發在冬季,而山西省春季大風天氣較多;O3濃度變化則與之相反,這是由于近地面O3污染多發生在夏、秋兩季.在低于30 d 的共有短周期中,PM2.5濃度與風速多呈負相關性,這是因為風速會影響PM2.5的擴散,較快的風速可以通過水平輸送和稀釋擴散效應降低PM2.5濃度.

圖8 臨汾市和大同市PM2.5 濃度、O3 濃度與風速的交叉小波圖、小波相干圖Fig.8 Cross wavelet and wavelet coherence diagram of PM2.5 concentration,O3 concentration and wind speed in Linfen City and Datong City

相較于PM2.5,兩市O3濃度與風速的共有短周期較少,其中大同市明顯多于臨汾市,在共有短周期中O3濃度與風速多呈負相關.在O3濃度與風速的共有長周期中,2015 年臨汾市O3濃度與風速呈正相關,2019-2020 年大同市O3濃度與風速呈負相關,這可能是由于氣溫等其他氣象因素的差異對O3污染造成了不同效應,使其在不同地區出現差異性[50].

3 結論

a)近年來山西省空氣質量持續向好,SO2、CO和PM2.5濃度均呈下降趨勢,但NO2、O3和PM10濃度仍處于高位.山西省11 個地級市冬季的空氣質量在空間上存在明顯差異,大同市和呂梁市的空氣質量明顯優于其他城市,而運城市和臨汾市的空氣質量較差.2020 年冬季各市空氣質量較2016 年同期有明顯改善.

b)不同社會經濟因子與SO2濃度和NO2濃度的相關關系存在明顯差異,且不同社會經濟因子與這兩種污染物濃度的相關關系空間異質性表現不同.其中,山西省人均GDP 與山西省南部城市SO2濃度呈負相關;城鎮綠化率對山西省北部城市SO2濃度影響較大,且越靠南其影響越小;人均GDP 對山西省南部城市NO2濃度的影響較大,且越靠北其影響越小;人口密度對山西省北部城市NO2濃度的影響較大,且越靠南其影響越小.

c)臨汾市和大同市氣象因子與PM2.5濃度的相關性一致,除平均相對濕度外,風速、降水量、平均氣溫與PM2.5濃度均呈負相關;兩市平均氣溫與O3濃度均呈顯著正相關,其余氣象因子與O3濃度相關性較弱,且在兩市中存在明顯差異.結合小波分析發現,兩市PM2.5濃度、O3濃度、風速均存在1 a 左右的周期變化特征,兩市PM2.5濃度和O3濃度與風速在不同的共有周期中存在不同的相關關系,其中,在短周期中兩市PM2.5濃度與風速多呈負相關,長周期中PM2.5濃度變化滯后于風速;短周期中O3濃度與風速的相關性不一致,長周期中風速變化滯后于O3濃度變化.

d)為確保山西省空氣質量持續改善,相關管理部門需根據不同地區空氣質量狀況、地方產業結構等經濟運行情況,進一步開展空氣污染溯源工作,并有針對性地制定相應污染治理措施,以減少人為因素對空氣質量的影響;同時,在極端氣候事件頻發的背景下,相關管理部門還需要預防不利氣象條件引發的空氣污染事件.

猜你喜歡
風速污染模型
一半模型
基于Kmeans-VMD-LSTM的短期風速預測
基于最優TS評分和頻率匹配的江蘇近海風速訂正
海洋通報(2020年5期)2021-01-14 09:26:54
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
堅決打好污染防治攻堅戰
當代陜西(2019年7期)2019-04-25 00:22:18
堅決打好污染防治攻堅戰
3D打印中的模型分割與打包
基于GARCH的短時風速預測方法
對抗塵污染,遠離“霾”伏
都市麗人(2015年5期)2015-03-20 13:33:49
主站蜘蛛池模板: 色婷婷综合激情视频免费看 | 无遮挡国产高潮视频免费观看| 东京热高清无码精品| 精品国产免费观看| 国产精品观看视频免费完整版| 国产区在线看| 久久久久久尹人网香蕉| 毛片久久网站小视频| 国产欧美日韩va另类在线播放| 国产精品女熟高潮视频| 婷婷色中文网| 一本一道波多野结衣av黑人在线| 伊人五月丁香综合AⅤ| 久久久精品无码一区二区三区| 美女被狂躁www在线观看| 亚洲人成网站18禁动漫无码| 亚洲无码四虎黄色网站| 亚洲中文字幕手机在线第一页| 国产精品久久久久无码网站| 91欧美在线| 成人精品视频一区二区在线 | 一级黄色欧美| 欧美色视频在线| 特黄日韩免费一区二区三区| 国产一级在线播放| 亚洲日本www| 亚洲无码视频图片| 国产亚洲高清在线精品99| 国产原创演绎剧情有字幕的| 精品国产电影久久九九| 日韩 欧美 国产 精品 综合| 青草午夜精品视频在线观看| 亚洲精品动漫在线观看| 亚洲国产亚综合在线区| 在线中文字幕网| 国产中文一区二区苍井空| 国产系列在线| 欧美日本在线一区二区三区| 亚洲欧美国产高清va在线播放| 九九免费观看全部免费视频| 欧美日韩国产在线观看一区二区三区| 国产福利在线免费| 欧美A级V片在线观看| 少妇精品网站| www.亚洲国产| 欧美日韩专区| 亚洲aⅴ天堂| 中文字幕欧美成人免费| 伊人色综合久久天天| 极品性荡少妇一区二区色欲| 中文字幕乱码中文乱码51精品| 亚洲国产看片基地久久1024 | 人人91人人澡人人妻人人爽| 欧美日本在线| 亚洲无码高清一区| 成年网址网站在线观看| 国产精品一老牛影视频| 国产精品尤物在线| 怡红院美国分院一区二区| 97成人在线视频| 精品少妇人妻av无码久久| 丁香亚洲综合五月天婷婷| 综合色天天| 亚洲国产精品无码AV| 国产一二视频| 欧美成人亚洲综合精品欧美激情| 国产情侣一区| AV天堂资源福利在线观看| 国产成人久视频免费| 日韩人妻少妇一区二区| 精品少妇人妻一区二区| 91精品伊人久久大香线蕉| 久久动漫精品| 婷五月综合| 国产激情无码一区二区免费| 青青青视频免费一区二区| 伊人久久久大香线蕉综合直播| 久久99国产精品成人欧美| 国产一区在线观看无码| 亚洲啪啪网| 少妇精品网站| 亚洲第一精品福利|