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

基于OPAQ的城市空氣質量預報系統研究

2016-06-09 08:53:57王淑瑩尹翠芳BinoMaiheuStijnJanssenLisaBlyth
中國環境監測 2016年3期
關鍵詞:模型

王淑瑩, 許 榮, 尹翠芳, Bino Maiheu, Stijn Janssen, Lisa Blyth

1.北京立博威拓環境技術有限公司, 北京 100083 2.中國環境監測總站, 國家環境保護環境監測質量控制重點實驗室, 北京 100012 3.比利時法蘭德斯技術研究院, 摩爾 BE-2400

基于OPAQ的城市空氣質量預報系統研究

王淑瑩1, 許 榮2, 尹翠芳1, Bino Maiheu3, Stijn Janssen3, Lisa Blyth3

1.北京立博威拓環境技術有限公司, 北京 100083 2.中國環境監測總站, 國家環境保護環境監測質量控制重點實驗室, 北京 100012 3.比利時法蘭德斯技術研究院, 摩爾 BE-2400

空氣質量預測在國內的關注度日益提高,傳統的空氣質量預測系統通常運用數值化學傳輸模型,利用物理方程來計算污染物的擴散、沉降和化學反應。而化學傳輸模型的預測準確性很大程度上需要依賴詳細的污染源排放信息和氣象模型的輸出結果。基于統計模型的OPAQ空氣質量預報業務系統,采用人工神經網絡算法,可預測各污染物的日均值或日最大值。并對北京空氣質量預報的結果進行了評價,OPAQ空氣質量預報業務系統對空氣質量預測的準確性較高,能夠利用較低的計算資源得到較為準確的預測結果。與數值預報相比,OPAQ空氣質量預報業務系統不需要大量的基礎數據作為輸入,可彌補數值預報的不足,并成為數值預報的有力補充。

空氣質量預報;統計模型;OPAQ

空氣質量與人類健康之間的關系已經引起了廣泛關注,有不少學者對此做了大量研究[1]。因此,空氣質量預報在污染事件發生前對公眾的預警和制定減災計劃方面尤為重要。例如,在比利時,如果PM10的日均濃度值連續2 d超過70 μg/m3,就會實行高速路限速的措施,且這種措施行之有效[2];或者,在高污染事件發生的城市實行車輛單雙號限行措施。

當前常用的空氣質量預報方法主要有:數值預報、統計預報2種方法。數值預報方法可以模擬從全球到街區等多種尺度。另外,數值預報需要大量的基礎數據作為輸入,除氣象預報數據外,還需要較為準確的污染源排放數據、詳細的地理環境數據、邊界條件數據等,并且需要做大量的計算機數值計算,需要投入較高的人力、物力和技術裝備。同時,由于污染源的污染物排放動態變化較大,難以獲得精確的污染源數據,因此數值預報目前的預報效果往往難以達到理想的效果。而統計預報方法則避免了上述數值預報方法的弊端,該方法只需氣象數據、污染物監測濃度數據即可進行預報。但傳統的統計預報方法只是針對單個站點進行模擬,不能實現整個區域的預報。本研究中,在統計模型的基礎上集成了插值模型,可以進行整個模擬區域的預報。并在應用統計模型過程中,對污染物濃度監測值與不同氣象參數之間的相關性進行了對比。另外,還探討了本研究中所用的空氣質量預報系統與數值預報系統在應用中的互補性。

另一方面,模型評價對模型預測性能也有重要意義。例如,目前,在歐洲空氣質量模型論壇上,有一批專業技術人員在協調制定預測模型的評價標準[3]。尤其是在不同的歐洲成員國分別使用多種不同的空氣質量預測模型的背景下,制定統一的標準非常重要。本研究中也制定了相應的評價方法對模型的預測結果進行了評價。

1 OPAQ空氣質量預報業務系統

法蘭德斯技術研究院(VITO)與比利時環保署合作,并由其資助開發了兩款統計模型用以替代化學傳輸模型,分別是OVL模型,利用人工神經網絡算法預測站點污染物濃度[4];RIO模型,利用地理統計插值技術來繪制污染物濃度分布圖[5]。兩款模型將長期作為比利時官方空氣質量預報業務系統中的一部分,并已在插值繪圖和高污染事件預測方面表現出比化學傳輸模型更為優越的性能[6]。在歐盟環境治理項目中的AirINFORM項目資助下,將OVL模型與RIO模型結合使用,稱之為OPAQ空氣質量預報業務系統(以下簡稱“OPAQ”)。

1.1 OPAQ所采用的預測模型介紹

OPAQ中的OVL模型采用人工神經網絡算法中的多層感知方法,其中一個隱含層代表一個輸出層。預測的結果是污染物的日均值或日最大值。假如模型在當天早上運行,則當天00:00 至模型開始運行最近時刻的監測數據以及未來幾天的氣象數據作為OVL的輸入數據。對每個監測站點配置幾個不同的模型進行訓練,不同模型的區別在于所采用的輸入參數不同。所預測當天早上的監測數據及行星邊界層高度(BLH)作為預測未來第N天的重要輸入參數。行星邊界層高度可認為是理查德系數(Ri)超過0.5時的高度[7]。其中Ri是浮力(由垂直溫廓線確定)與慣性力(由大氣湍流確定)的比率。如果Ri非常小(甚至是負的),湍流非常強烈,并且足以打破穩定的溫廓線將顆粒物帶至垂直方向,則從底部到行星邊界層高度之間,顆粒物可以隨湍流在此之間擴散。研究表明行星邊界層高度在基于人工神經網絡算法的空氣質量預報系統中有非常重要的作用[8]。另外,還使用了風速與風向、溫度、濕度及云量等氣象參數提高預測的準確性[9-11]。對每個監測站點都配置幾種模型進行訓練,最后選擇一個最優的模型作為日常預報模型。

空氣質量監測值通常符合的正態分布,高濃度值對應的數據量明顯小于其他濃度值所對應的數據量。這表明,高濃度事件發生概率較小。盡管如此,依然期望對高濃度事件進行準確預測。對此,研究中將人工神經網絡算法的訓練樣品進行重采樣,使其呈均勻分布,這樣可以給高濃度事件賦予和其他濃度相同的權重值,從而提高高濃度事件預測的準確性。

相對于原始的神經網絡輸出而言,OPAQ具有動態偏差修正功能。將日常預測值及實際監測值存儲在臨時數據庫中,通過該數據庫,可以將過去的預測誤差統計出來,再通過不同的方法,如計算過去預測誤差的均值或加權平均,計算得到當天的預測誤差。然后將預測誤差應用到原始的模型輸出結果上,這就是所謂的實時修正。同時,對于每個站點、每個預測時間采用何種修正方法都需要進行優化。

1.2 OPAQ預測模型與數值模型的互補性

OPAQ的機器自學習和統計特性可認為是復雜的數值模型預報系統的補充。這2種方法都有各自的優點,但也也都存在一些不足。有學者曾經做過統計模型和數值模型的比較[12],數值化學傳輸模型通過擴散和化學反應方程很好地描述了污染物的物理化學過程,并且充分利用了三維氣象模型的輸出結果。因此,也可以模擬污染物的長距離傳輸,而這對于污染事件起因的研究是非常重要的。由于化學傳輸模型需要使用污染物排放數據,且輸出結果會隨著排放數據的變化而發生相應的變化。因此,化學傳輸模型可以用于情景分析。然而,如果沒有準確的排放數據將會影響空氣質量預報的準確性。很多化學傳輸模型進行評價、預測時都需要進行校正或數據同化,對比空氣質量模擬結果和實際監測值來提高模型的性能。此外,化學傳輸模型對計算機的計算能力也有較高要求,尤其是在設置較高空間分辨率的情況下。因此,通常運行化學傳輸模型需要較昂貴的硬件設備。另外,要運行基于數值模型的空氣質量預報系統以及解釋其輸出結果也需要在大氣科學方面有較深背景的專業人員。

與之相反,OPAQ在安裝和部署方面相對簡單易學。首先,與數值模型預報系統相比對硬件和人員的要求相對較低。其次,它只需要空氣質量監測數據,無需污染物的排放數據。OPAQ或其他統計預報方法同樣也有一些不足之處,例如,統計模型不能很好地模擬過去沒有發生過的情景,統計方法的特性也決定了其不能很好地解釋污染的過程。盡管如此,當訓練神經網絡模型時,為避免過度擬合化,須將統計模型的性能進行一般化處理。相對于復雜的化學傳輸模型預報方法,統計預報方法的預測準確性與其相當,甚至要優于前者,對于高污染事件的預測性能也是如此[12-13]。

當空氣質量預測用于預警或是否采取消減措施(比如機動車限速或其他措施)時,就需要評價預測模型是否能夠預測超標值。為了量化模型的預測性能,使用了預測成功指數(SI)這個指標。SI基于散布圖來計算,見圖1。

圖1 預測值與監測值的散布圖

圖1中橫坐標為監測值,縱坐標為預測值,單位均為μg/m3,100 μg/m3為超標臨界值。利用預測到的超標天數總量與預測到的未超標天數總量分別與總的超標天數及總的未超標天數的比值相加再減去1得到,其中減1是為使指數在-1~1。詳見式(1):

(1)

式中:SI為預測成功指數,N3為預測到的超標天數,N4為未預測到的超標天數,N1為預測到的達標天數,N2為未預測到的達標天數。

圖2展示了統計與數值兩種模型的預測性能對比,以PM10為例,比較了兩種模型的預測值與監測值的相關性系數(r)及SI兩種預測性能。圖中淺灰色方塊和灰色三角形分別代表OPAQ中OVL的不同模型的預測結果,深灰色菱形代表化學傳輸模型CHIMERE的預測結果。參與評價的數據為比利時2008年11月—2009年4月3種模型運行的預測數據[6]。結果表明,統計模型在預測高污染事件的發生方面具有較高的潛力。

圖2 統計模型與數值化學傳輸模型就PM10日均值的預測性能對比

2 OPAQ應用于北京空氣質量預報的評價

2.1 模型配置

為北京建立的OPAQ空氣質量模型,采用中國環境監測總站發布的北京萬壽西宮、定陵、東四、天壇、農展館、官園、海淀區萬柳、順義新城、懷柔鎮、昌平鎮、奧體中心及古城12個國控點的小時監測數據,美國國家環境預報中心(NCEP)下載的再分析氣象數據(FNL)進行人工神經網絡算法訓練。FNL數據基于全球預報系統(GFS),與GFS數據相比,FNL數據中使用了更多的實際觀測數據。因此,對于統計模型來說,所使用的GFS數據和FNL數據具有一致性,使用FNL數據作為歷史數據進行訓練,GFS數據作為實時預測數據進行預報,所采用的氣象參數有2 m高處空氣溫度,2 m高處相對濕度,10 m高處緯向和經向風速,BLH(行星邊界層高度),1 000 hPa處分別與 975、 950、925、900、850 hPa處的溫度差。然而,FNL和GFS數據的空間精度較低,可以使用當地氣象模型數據進行改善(見4結語部分)。

首先,對PM2.5、PM10、NO2的日均濃度值和O3的8 h滑動平均最大濃度值與氣象數據進行了相關性分析。相關性分析使用污染物濃度的對數值,有些氣象數據也使用其對數值,如BLH,見圖3。

圖3 北京東四站點NO2監測值與各氣象參數相關性分析

圖3中,從上到下,從左至右,分別為NO2日均濃度值的對數值log(1 + NO2)分別與各個氣象參數的相關性系數,這些氣象參數分別為2 m高處空氣溫度(2 m Temp),2 m高處相對濕度(2 m RH),10 m高處緯向(10 m V)和經向風速(10 m U),行星邊界層高度(BLH),1 000 hPa處分別與 975、950、925、900、850 hPa處的溫度差(分別用IS975、IS950、IS925、IS900、IS850表示)。另外,圖中還標出了相關性系數值。研究在之前HOOYBERGHS J等[4]所采用參數的基礎上,OPAQ所用的OVL模型還增加了反映大氣穩定性的垂直溫度信息,分別為975~850 hPa高處分別與1 000 hPa高處的溫度差,提高了模型性能。

HOOYBERGHS J等研究發現,NO2、PM10、PM2.5與BLH呈顯著負相關,O3與BLH呈正相關,這與NO2、NO和O3之間的光化學平衡有關。北京東四站點O3與氣象參數的相關性分析結果(圖4)與HOOYBERGHS J等的研究結論一致,尤其是O3與2 m處溫度日均值相關性顯著。

O3濃度在多數站點與1 000 hPa和850、900 hPa處溫差相關性顯著,而與近地面處的溫差相關性較小。目前,對于O3濃度與較高處的溫差相關性更顯著的原因還不清楚。

另外,分析表明(圖3),當風向為北風時,污染物的濃度明顯低于風向為南風時(經向風速分別為負和正時)。同樣,當風向為西風時污染物濃度較低,而風向為東風時污染物濃度較高。這可能與城市北部、西部多山區和郊區分布有關。

表1列舉了各站點污染物監測值與氣象參數的相關性系數的均值。氣象參數從左至右依次為2 m高處空氣溫度,2 m高處相對濕度,10 m高處緯向和經向風速,BLH,1 000 hPa處溫度分別與 975、950、925、900、850 hPa處溫度差。

為每項污染物配置5組模型,計算模擬預測值和監測值的均方根誤差(RMSE),然后為每個站點每個預測天(預測當天至未來第5 d)選取最優模型和實時修正方法。相較HOOYBERGHS J等的研究,配置的模型加入了風速、星期屬性(是否工作日)、溫度、濕度參數后,OPAQ預測結果得到了明顯改善。新增的850 hPa高度逆溫強度提高了預報結果的準確性,但加入其他高度處逆溫強度對模型改善效果不顯著。

圖4 北京東四站點O3監測值與各氣象參數相關性分析

表1 北京各站點PM2.5、PM10、NO2濃度日均值和O38 h滑動平均濃度日最大值分別與氣象參數日均值相關性系數的平均值

污染物2mTemp2mRH10mV10mUBLHIS975IS950IS925IS900IS850PM100.150.190.42-0.33-0.24-0.21-0.16-0.010.150.32PM2.50.100.380.47-0.44-0.41-0.22-0.17-0.040.160.40NO2-0.060.110.37-0.28-0.480.030.100.230.370.46O30.670.060.18-0.150.46-0.39-0.31-0.26-0.33-0.45

2.2 模型評價

使用2015年以前所有可獲得的歷史數據對OPAQ預測模型進行訓練,并用2015年1月1日—4月30日的數據進行檢驗。以下重點討論O3和 PM2.5的預報情況。圖5、圖6分別為北京東四站點PM2.5和O3的預報與監測值對比時間序列圖。圖5中監測值表示PM2.5監測日均值,預測值表示OPAQ原始預報值,修正值表示采用實時修正方法的值。圖中用黑色虛直線標出超標值(75 μg/m3)。

圖5 北京東四站點PM2.5的預報與監測值對比時間序列圖

圖6 北京東四站點O3的預報與監測值對比時間序列圖

圖6中監測值表示O3的8 h滑動平均監測濃度最大值,預測值表示OPAQ原始預報值,修正值表示實時修正值。圖中用黑色虛直線標出超標值(160 μg/m3)。從上至下分別為當天至未來第5 d的時間序列圖。

由圖5、圖6可見,OPAQ預測模型很好地預測了O3和 PM2.5的濃度。且當監測值和原始預報值差距較大時,采用實時修正方法得到的結果更好。

圖7為模擬預測值與實際監測值的相關性系數統計箱型圖,橫坐標為預報未來天數(用當天至第5 d表示)。預報的每個天數均作了原始預報結果(左側箱圖)和實時修正后結果(右側箱圖)的對比。

圖7 模擬預測值與實際監測值的相關系數統計箱型圖

圖8 PM2.5其他檢驗指標的統計箱型圖

由圖7可見,O3每個預報天的模擬預報結果與實際監測值的相關性系數都約為0.8,說明氣象預報結果較好。而對于PM2.5,未來第1 d的相關系數突然下降,約0.6~0.7,未來第2 d以后相關性系數相當。實時修正結果相關性(右側箱圖)和原始人工神經網絡算法結果相關性(左側箱圖)略有不同,可能是因為實時修正進行了偏移修正。RMSE的結果也類似,如圖8中a圖所示,實時修正結果明顯低于原始預報結果。

圖8還展示了超標天預報的結果統計,設定超標濃度為75 μg/m3。PM2.5的未來第1 d預測成功指數約55~60,實時修正結果(右側箱圖)有明顯提高。正確預警指數(FCF)為預報出正確的超標天數除以總超標天數,錯誤預警指數(FFA)為預報錯誤的超標天數(實際不超標)除以總超標天數。

由圖8可見,FFA約為10%~20%,實時修正結果稍差一些。而FCF統計中,實時修正結果有明顯提高,約70%~80%。

3 OPAQ空氣質量預報業務系統在國內的應用案例

針對OPAQ在其他幾個城市的應用情況,對其預測準確率做了相應的統計,如表2所示。表中的預測準確率是根據中國環境監測總站編寫的《環境空氣質量預報預警方法技術指南》第二章第六節中關于級別準確分的描述所做的統計,統計了完全準確天數與級別準確天數占所有預報天數的比例,其中完全準確是指預報值±10個點(或±10%)與實況相符;表中的級別準確是指預報不跨級且級別與實況相符。

表2 OPAQ在國內應用案例的預測級別準確率匯總

注:“—”表示無數據。

4 結語

1)通過OPAQ統計預報系統與數值預報系統的對比,探討了與數值預報系統的互補性。結果表明,OPAQ統計預報系統與數值預報各有優劣,在應用中可以互為補充。

2)用相關性系數及預測成功指數兩種指標對統計和數值模型進行預測結果的評價。結果表明,統計模型在預測高污染事件的發生方面具有較高的潛力。

3)與北京監測站預報結果的初步比對表明,OPAQ對空氣質量預測的準確性較高。O3每個預報天的模擬預報結果與實際監測值的相關性系數都約為0.8,而對于PM2.5,未來第1 d至第5 d的相關系數在0.6~0.7,而且實時修正后的結果在監測值和原始預報值有差距時更好。

4)OPAQ對污染等級的預測準確性較高,預測未來24 h的級別準確率在56.1%~70%。但是,還有待與數值模型進行預測性能的對比。

5)基于人工神經網絡算法的統計預報系統能夠利用較低的計算資源來較為準確地預測空氣質量。另外,如果使用能夠很好地反映當地天氣狀況的氣象模型參數,可以對OPAQ做進一步的改進,并提高其預測準確率。

[1] World Health Organization.Review of evidence on health aspects of air pollution-REVIHAAP Project 309[R].2013.

[2] LEFEBVRE W,FIERENS F,TRIMPENEERS E,et al.Modeling the effects of a speed limit reduction on traffic-related elemental carbon (EC) concentrations and population exposure to EC[J].Atmospheric Environment,2011,45:197-207.

[3] THUNIS P,PEDERZOLI A,PERNIGOTTI D.Perfor-mance criteria to evaluate air quality modeling applic-ations[J].Atmospheric Environment,2012,59:476-482.

[4] HOOYBERGHS J,MENSINK C,DUMONT G,et al.A neural network forecast for daily average PM concentrations in Belgium[J].Atmospheric Environment,2005,39:3 279-3 289.

[5] JANSSEN S,DUMONT G,FIERENS F,et al.Spatial interpolation of air pollution measurements using CORINE land cover data[J].Atmospheric Environment,2008,42:4 884-4 903.

[6] JANSSEN S,DUMONT G,FIERENS F,et al.Land use to characterize spatial representativeness of air quality monitoring stations and its relevance for model validation[J].Atmospheric Environment,2012,59:492-500.

[7] VOGELEZANG D H P,HOLTSLAG A A M.Evaluation and model impacts of alternative boundary-layer height formulations[J].Boundary-Layer Meteorology,1996,81:245-269.

[8] DUTOT A L,RYNKIEWICZ J,STEINER F E,et al.A 24-h forecast of ozone peaks and exceedance levels using neural classifiers and weather predictions[J].Environmental Modelling and Software,2007,22:1 261-1 269.

[9] 程念亮,李云婷,邱啟鴻,等.2013年北京市PM2.5重污染日時空分布特征研究[J].中國環境監測,2015,31(3):38-42.

[10] 陳強,梅琨,朱慧敏,等.鄭州市PM2.5濃度時空分布特征及預測模型研究[J].中國環境監測,2015,31(3):105-112.

[11] 沈路路,王聿絢,段雷.神經網絡模型在O3濃度預測中的應用[J].環境科學,2011,32(8):2 231-2 235.

[12] ZHANG Y,BOCQUET M,MALLET V,et al.Real-time air quality forecasting,part I:History,techniques,and current status[J].Atmospheric Environment,2012,60:632-655.

[13] 張靜,李旭祥,許先意,等.大氣環境數據分析預測方法對比研究[J].中國環境監測,2010,26(6):66-69.

Study on Urban Air Quality Forecasting with OPAQ

WANG Shuying1,XU Rong2,YIN Cuifang1,BINO Maiheu3,STIJN Janssen3,LISA Blyth3

1.Beijing LIBOVITO Environmental Technology Company, Beijing 100083, China 2.State Enrironmental Protection Key Laboratory of Quality Control in Environmental Monitoring, China National Environmental Monitoring Centre, Beijing 100012, China 3.Flemish Institute for Technological Research, Mol BE-2400, Belgium

Air quality forecast has been and continues to be a growing concern in China. Traditional air quality forecasting systems employ deterministic chemical transport models in which pollutant dispersion, deposition and chemical reactions are computed from the governing physical equations. Detailed pollutant emission information and meteorological model output areneeded to run the chemical transport models. An alternative statistical forecasting system OPAQ is discussed that is based on the BP neural network algorithm. The results of the OPAQ forecast system for a nuhPaer of monitoring stations in Beijing show that forecast accuracy of OPAQ is high and it can take advantage of lower computing resources to predict air quality. Compared with deterministic chemical transport model, it doesn’t need multiple input data and can be a complementary approach of traditional method.

air quality forecasting;statistical model;OPAQ

2015-09-08;

2015-12-31

環保公益性行業專項“京津冀地區大氣重污染過程應急方案研究”(201309071)、“京津冀城市大氣邊界層過程對重污染形成的影響研究”(201409001-03);科技部科技支撐計劃環境領域項目“大氣復合污染區域聯合預測預報關鍵技術研究”(2014BAC22B04)、“京津冀空氣監測預報及防控技術研究與示范”(2014BAC06B04)、中科院先導項目“大氣灰霾追因與控制專項數值模式與協同控制方案課題”(XDB05030200)

王淑瑩(1982-),女,山西霍州人,碩士,工程師。

許 榮

X830.3

A

1002-6002(2016)03- 0013- 08

10.19316/j.issn.1002-6002.2016.03.02

猜你喜歡
模型
一半模型
一種去中心化的域名服務本地化模型
適用于BDS-3 PPP的隨機模型
提煉模型 突破難點
函數模型及應用
p150Glued在帕金森病模型中的表達及分布
函數模型及應用
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
3D打印中的模型分割與打包
主站蜘蛛池模板: 理论片一区| 国产高清在线观看| 欧美一级视频免费| 午夜国产精品视频| 亚洲欧美日韩天堂| 久久青草精品一区二区三区| 午夜综合网| 999国内精品久久免费视频| 国产精品区视频中文字幕| 亚洲国产天堂在线观看| 国产打屁股免费区网站| 亚洲第一精品福利| 国产人人射| 激情爆乳一区二区| 成AV人片一区二区三区久久| 国产欧美精品午夜在线播放| 久久a级片| 国产精品成人一区二区| 99色亚洲国产精品11p| 婷婷亚洲最大| 91精品aⅴ无码中文字字幕蜜桃 | 久久黄色一级视频| 国产成人在线小视频| 免费A级毛片无码免费视频| 国产精品免费电影| 永久在线精品免费视频观看| 99久久人妻精品免费二区| 亚洲国产亚综合在线区| 亚洲熟女中文字幕男人总站| 亚洲伊人天堂| 91精品伊人久久大香线蕉| 99在线视频免费观看| 久久精品国产一区二区小说| 国产区网址| 欧美一级高清免费a| 久久青草免费91线频观看不卡| 国产清纯在线一区二区WWW| 天天操天天噜| 欧洲在线免费视频| 国产精品永久不卡免费视频| 亚洲一区二区三区国产精华液| 久久久久青草大香线综合精品 | 丁香婷婷激情综合激情| 国产精品久久久久久久久| 国产福利2021最新在线观看| 久久久久久久久18禁秘| 米奇精品一区二区三区| 天堂在线亚洲| 亚洲一区波多野结衣二区三区| 欧美α片免费观看| 中文字幕亚洲精品2页| 国产毛片久久国产| 亚洲日韩精品综合在线一区二区| 亚洲精品欧美日本中文字幕 | 在线综合亚洲欧美网站| 久久国产高潮流白浆免费观看| av在线手机播放| 免费毛片网站在线观看| 欧美a级在线| 久久久精品国产SM调教网站| 日本三区视频| 91精品视频播放| 久久久久国产精品免费免费不卡| 亚洲精品福利视频| 91精品亚洲| 国产精品福利尤物youwu | 人妻21p大胆| 国产理论精品| 伊人久久精品无码麻豆精品| 国产免费精彩视频| 亚洲精品黄| 好吊日免费视频| 伊人久久大香线蕉影院| 九色在线观看视频| 国产91无毒不卡在线观看| 亚洲最新在线| 51国产偷自视频区视频手机观看| 国产视频欧美| 日韩欧美在线观看| 亚洲国产成人无码AV在线影院L| 在线观看视频一区二区| 成人精品亚洲|