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

基于自適應聯邦濾波的列控車載設備測速測距算法

2022-02-16 06:56:44楊志杰程劍鋒李一楠
中國鐵道科學 2022年1期
關鍵詞:檢測模型

張 淼,楊志杰,程劍鋒,李一楠

(1.中國鐵道科學研究院 研究生部,北京 100081;2.中國鐵道科學研究院集團有限公司 通信信號研究所,北京 100081)

測速測距是列控車載設備的關鍵技術之一。列控系統行車許可的生成必須依賴列車的速度和位置,因此列控車載設備實時、準確地采集和計算列車運行速度和走行距離尤為重要[1]。考慮到轉速傳感器實現較簡單、受環境與地形限制小、測量精度與可靠性高等優點,傳統的測速測距系統普遍采用轉速傳感器技術。為了保證測距測速系統的互補和冗余關系,一般將2 個轉速傳感器分別安裝在不同輪軸處,車輪旋轉1 周轉速傳感器會產生固定的脈沖個數,結合單位時間內的脈沖計數和車輪直徑,即可得到列車運行速度;再對速度積分,即可得到列車走行距離。

測速測距系統的主要功能需求為計算當前列車的加速度、運行速度和走行距離,利用列車運動模型來描述列車的加速度和運行速度隨時間變化的規律,但其建模過程中會不可避免地產生系統狀態噪聲,轉速傳感器也會受外界環境影響而產生測量噪聲。這2 種噪聲帶來的偏差常被近似為零均值高斯白噪聲,故采用對白噪聲有很好抑制效果的卡爾曼濾波,能提供統計意義上的最優估計。

已有文獻對測速測距系統核心算法進行研究,但多存在不足。文獻[2]基于卡爾曼濾波,針對高速列車測速測距系統建立了恒速和恒加速的列車運動模型,并進行仿真驗證,但其提出的恒速模型對加速度變化不敏感,對列車速度變化的感知存在延時效應,且恒加速模型在處理加速度變化大的臨界點時會出現“超調現象”。文獻[3]在建立列車運動模型時綜合考慮了牽引力、制動力、運行阻力、坡度和彎道阻力等因素,但實踐中測速測距系統無法直接獲得上述變量的取值,而是通過加速度計或模型方法直接獲取列車運行加速度進而得到列車運行速度。文獻[4]介紹了當輪對受軌面條件、環境氣候和列車運行狀況等影響時,列車會發生空轉/滑行,影響計算結果準確度,隨著空轉/滑行時間的增加,僅依靠單一的加速度計會發生漂移,導致累計誤差增大。文獻[5]基于滑模控制的卡爾曼濾波融合轉速傳感器和加速度計信息,雖然實現了空轉/滑行誤差降低、定位精度提升,但僅對仿真數據進行分析說明。文獻[6]利用轉速傳感器和雷達,建立了空轉/滑行檢測及誤差校正計算模型,然而雷達的使用條件苛刻,既要保證安裝精度,又會受到天氣、軌面狀況等因素影響,故障率偏高。

本文在僅采用轉速傳感器技術的基礎上,提出1 種基于自適應聯邦濾波的列控車載設備測速測距算法,先建立自適應參數列車運動模型,對每個轉速傳感器的采集數據進行卡爾曼濾波處理;再利用檢測方法判定列車空轉/滑行情況,通過空轉/滑行速度校正模型得到虛擬傳感器的輸出,即空轉/滑行校正速度;然后根據實際情況動態調整信息分配系數,利用聯邦濾波算法將自適應參數列車運動模型和空轉/滑行時速度校正模型的計算結果進行數據融合,得到列車運行速度和走行距離。通過列車正常運行場景和制動時輪對滑行場景下的仿真分析,驗證算法應用于列控車載設備測速測距系統的準確性和安全性。

1 自適應參數列車運動模型

在列車運行過程中,司機會施加牽引或制動,列車也會因為彎道和坡度變化而改變加速度,故列車在實際運行中很難保持恒加速狀態。為能更好地跟蹤列車運動趨勢,在速度及加速度變化較大時仍能實時追蹤列車運行速度,故選取自適應參數機動目標模型估計方法[7],建立自適應列車運動模型,由狀態方程和測量方程組成。

考慮到測速測距系統需要在線處理大量數據,此處選用二階矩陣自適應參數列車運動模型,以列車運行速度v和加速度a為狀態變量,以便在較短的處理時間內得到較優的濾波結果。計算時,以機動頻率α和加速度方差δ2作為模型的自適應參數,先利用狀態方程預測列車運行狀態和協方差,在獲得預測結果的基礎上,結合實際測量值更新列車運行狀態,估計當前狀態下的列車運行速度v和加速度a;再根據列車加速度估計值,實時修正自適應參數α和δ2,并不斷更新列車運動模型和系統狀態噪聲,最終實現模型計算結果與列車實際運動狀態趨于一致。

1.1 狀態方程

由于轉速傳感器無法直接觀測量到列車運行加速度,故采用列車運行過程中前2 個周期估計加速度的均值作為狀態方程的加速度,取列車運行時的速度v和加速度a作為狀態變量,則狀態向量x=(v,a)T,構造狀態方程可表示為

式中:k為當前時刻;T為k?1 至k時刻的測速測距采樣周期;xk?1和xk分別為k?1 時刻和k時刻下的狀態向量;A為狀態轉移矩陣,Ak,k?1描述了系統狀態從k?1 至k時刻的變化規律;U為輸入控制矩陣,Uk,k?1引入加速度參量來描述k?1 至k時刻外力的作用;為k?2時刻至k?1時刻的列車估計加速度均值;Wk?1為k?1時刻下的系統狀態噪聲,假設其是均值為零、方差為Qk?1的白噪聲矢量。

利用系統的自適應參數可對對稱矩陣Qk?1實現實時調整,矩陣元素分別用q11,q12,q21和q22表示,即

其中,

式中:E為期望;rk?1(0)為k?1 時刻列車加速度估計值的自相關系數;rk?1(1)為k?1時刻與k?2時刻列車加速度估計值的相關系數。

1.2 測量方程

測速測距系統中,利用已知的轉速傳感器齒數M和車輪直徑D,可得到轉速傳感器每轉過1 個齒對應的列車走行距離為若1 個測速測距采樣周期T內采集到的脈沖為n個,那么轉速傳感器的采樣頻率f=測量速度v=則列控車載設備測速測距系統的測量方程可表示為

式中:Zk為k時刻下的測量值,在取值上與fk相等;fk為轉速傳感器在k時刻的脈沖頻率;Hk為測量矩陣;Vk為k時刻下的轉速傳感器測量噪聲,假設其是均值為零、方差為Rk的白噪聲矢量,且與系統狀態噪聲Wk相互獨立。

2 空轉/滑行時列車速度校正模型

受列車運行狀況、軌面條件和環境氣候等因素影響,列車運行過程中輪對容易發生空轉/滑行,此時轉速傳感器測得的列車運行速度與實際情況間會出現較大偏差,并且這一偏差無法僅憑轉速傳感器自身進行有效修正[8]。

為解決該問題,設計列車速度校正模型,并將其視為1 個虛擬傳感器,與安裝在不同輪軸上的2個轉速傳感器共同構成列控車載設備測速測距系統。列車速度校正模型可快速檢測出轉速傳感器所在輪對的運行狀態,當2 個輪對均出現空轉/滑行動作時,便由速度校正模型獲取并計算空轉/滑行狀態下的列車運行速度和走行距離,以便獲取的結果更為準確。

2.1 空轉/滑行檢測方法

按照列車輪軸的動力形式(是否帶有牽引力或制動力),車輪可劃分為帶制動的動力輪、帶制動的從動輪和不帶制動的從動輪這3類,其中帶制動的動力輪較易發生空轉和滑行,帶制動的從動輪較易發生滑行,不帶制動的從動輪不易發生空轉和滑行[9],因此空轉和滑行檢測需要分車輪情況具體討論。

常規的檢測空轉/滑行的方法主要有滑移率法、速度差法和加/減速度檢測法[10],這3 種方法既可單獨使用也可組合使用,使用時根據不同車型的牽引和制動性能選取相應判據參數。為能檢測到較微弱的列車空轉/滑行情況,在這3 種方法的基礎上,本文再引入加/減速度變化率的檢測方法[11],結合列車輪軸類型檢測輪對是否處于空轉/滑行狀態。

1)變量定義

當列車以速度v、加速度a運行時,第i個轉速傳感器所在輪對的運行速度為vi(i=1,2),加速度為ai,相應的加速度變化率(加加速度)為ji=則當前測速測距采樣周期下2 個輪對的速度差為Δv=|v1?v2|。

為描述車輪的滑行程度,定義車軸對應車輪的滑移率λi為列車在制動時車輪的抱死程度,則λi=當車輪完全抱死純滑動時取λi=1。同理,為描述車輪的空轉程度,定義車軸對應車輪的滑轉率λri為列車加速起步時車輪的滑轉程度,即λri=×100% ,當車輪滑轉完全不前進時取λri=1。

由于不同車型的牽引和制動性能并不相同,考慮通過調整門限值的形式來適應不同車型的空轉/滑行情況,設速度差門限為Δvlim,最大加速度門限為amax,最大減速度門限為amin,滑移率門限為λlim,滑轉率門限為λrlim。

2)空轉檢測

對于較易發生空轉的帶制動動力輪,當a>0時首先對每個轉速傳感器所在輪對分別進行檢測;經過測速測距采樣周期T后,若ai+jiT>amax或λi>λlim這2 個條件任一成立,則檢測到此轉速傳感器所在輪對發生空轉;若此時還滿足Δv>vlim的條件,由于2 個轉速傳感器所在輪對互為參考,無法判定到底哪個輪對發生空轉,考慮安全性,檢測結果為這2個轉速傳感器所在輪對同時發生空轉。

3)滑行檢測

對于較易發生滑行的帶制動動力輪和帶制動從動輪,當a<0 時首先對每個轉速傳感器所在輪對分別進行檢測;經過測速測距采樣周期T后,若ai+jiTλrlim這2 個條件任一成立,則檢測到此轉速傳感器所在輪對發生滑行;若此時還滿足Δv>vlim的條件,由于2 個轉速傳感器所在輪對互為參考,無法判定到底哪個軸發生滑行,考慮安全性,檢測結果為這2 個轉速傳感器所在輪對同時發生滑行。

2.2 列車速度校正模型

檢測到空轉/滑行后,需要通過列車運行速度校正模型獲得準確且安全的列車運行速度和走行距離。校正時,根據不同車型車輛的空轉/滑行程度和黏著恢復情況,設定相應的空轉/滑行持續校正時間。為簡單起見,以amax和amin分別作為速度校正模型的加速度和減速度。

當檢測到2 個轉速傳感器所在輪對同時發生滑行時,記其在上個測速測距采樣周期的速度為vbefore,加速度為abefore。在vbefore的基礎上,需計算k時刻下的列車校正速度vcorrect(k),即

在滑行持續校正時間內,實時監測車輪輪軌黏著恢復情況,將2 個轉速傳感器所在輪對的測量值與校正后的速度進行比較,以vi>vcorrect時記作時刻n,此時車輪輪軌黏著已經恢復,

校正速度不應小于車輪輪軌黏著恢復后的速度,故增大列車校正速度,即

待校正速度大于2 個轉速傳感器的測量速度,即vcorrect>vi后繼續以式(4)推算列車運行速度,校正距離即為校正速度的積分。

列車減速滑行運行過程中,當車輪抱死時測量速度v1和v2同時為零,此時測速測距系統缺乏有效的測量值,因列車運行速度不可能突變,綜合考慮可用性和安全性,取較小的默認減速度作為校正減速度并繼續推算列車運行速度,使列車運行速度平穩下降直至列車完全停止(即校正速度為零)。

空轉時情況與滑行時相反,先以速度校正模型的加速度求解速度,當車輪輪軌黏著恢復后,校正速度不應大于車輪輪軌黏著恢復后的速度,故減小列車校正速度。校正結果同樣成立,不再展開。

列車速度校正模型基于加速度二重積分法[12?13]獲得計算速度和距離,同時考慮了轉速傳感器實時的測量值,這樣通過計算速度與測量速度的不斷比較,能夠不斷修正校正速度和校正距離,使獲得的列車運行速度和走行距離符合車輪輪軌實際黏著情況。

3 基于聯邦濾波的數據融合算法

3.1 測速測距數據融合結構

基于聯邦濾波的融合算法具有容錯性能好、融合結果精確、計算量小等優點[14],且聯邦濾波結構能夠保證數據先分散處理、再全局融合,通過時間更新和測量更新進行算法遞推,以融合后的全局最優估計值作為系統的最終輸出結果,因此在列控車載設備測速測距系統中采用基于聯邦濾波的融合算法。

設計列控車載設備測速測距系統由2 個轉速傳感器和1 個虛擬傳感器(空轉/滑行速度校正模型)共3 個部分組成,其算法框圖如圖1所示。圖中:箭頭表示數據流的走向;為全局狀態估計;Pf為協方差矩陣;(i′=1,2,3,i′=3 時代表虛擬傳感器)為轉速傳感器對應子濾波器的狀態估計,由速度和加速度估計值構成;Pi′為子濾波器狀態估計的誤差協方差矩陣;βi′為3 個子濾波器的信息分配系數,

圖1 基于聯邦濾波的測速測距數據融合結構

列車實際運行過程中,測距測速系統難免會受到外界的電磁干擾,測速齒輪箱體和轉速傳感器也有可能在鋼軌銜接處受到振動的影響,這些情況都將導致轉速干擾信號的產生。為避免輪速發生突變造成空轉/滑行的誤判,在數據融合處理前,先要對測量結果進行數據預處理[15]。

由圖1可知:2 個轉速傳感器輸出的數據分別經過數據預處理去除干擾信號,輸出轉速傳感器有效測量值;校正模型作為虛擬傳感器,輸出轉速傳感器所在輪對發生空轉/滑行時的列車校正速度;將測量值和校正值分別輸入3 個相互獨立、并行運行的子濾波器,其中2 個轉速傳感器對應的子濾波器1 和子濾波器2 均采用自適應參數列車運動模型,得到轉速傳感器所在輪對速度和加速度的估計值,與子濾波器3 得到的校正速度估計值共同描述列車運行的動態過程;最終由主濾波器根據列車空轉/滑行情況和轉速傳感器的濾波估計值動態調整信息分配系數,主濾波器數據融合得到的狀態估計值即為列車運行速度和加速度。

在聯邦濾波算法中,將系統噪聲Qf和狀態協方差Pf通過信息分配系數βi′分配到各子濾波器,信息分配系數βi′滿足信息守恒原理并能隨著系統狀態的動態變化進行相應調整,以增加系統的安全性和魯棒性,即

主濾波器和各子濾波得到的狀態估計(k)與協方差矩陣Pi′(k)對應融合,得到全局狀態估計及協方差矩陣,即

其中,

3.2 動態調整信息分配系數方案

確立信息分配系數方案是建立基于聯邦濾波的列車測速測距系統的關鍵[16]。對信息分配系數β值的動態調整方案需要滿足以下速度處理安全原則。

(1)當檢測到某1 個轉速傳感器所在輪對發生空轉/滑行時,為減少空轉/滑行輪對的采集數據影響,將此輪對的信息分配系數取值調小,另1個輪對的則相應調大。

(2)當檢測到2 個轉速傳感器所在輪對同時發生空轉/滑行時,采用空轉/滑行校正策略得到列車運行速度,2 個輪對的信息分配系數取值均調小,列車校正速度的信息分配系數相應調大。

(3)當檢測到2 個轉速傳感器所在輪對均正常但存在速度差時,從控車安全性角度出發,速度值較大輪對的車軸信息分配系數取值較大,速度值較小的信息分配系數取值也相應較小。

在滿足上述安全原則的基礎上,設進入子濾波器的數據為1個包含速度和狀態信息的向量,對于2個轉速傳感器,若檢測到傳感器所在輪對發生空轉/滑行則狀態值Si取0,否則Si取1;對于虛擬傳感器,狀態值此時信息分配系數βi′(k)的動態表達式為

由于存在分配系數的倒數運算,βi′(k)不能取0,故計算時取極小值。

根據不同的運行情況動態,利用式(8)調整得到相應的信息分配系數后,通過基于聯邦濾波的數據融合算法,可將自適應參數列車運動模型和列車速度校正模型的計算結果進行濾波融合,輸出列車運行速度;再對速度積分,即得到列車走行距離。

4 測速測距算法的仿真及結果分析

為驗證基于自適應聯邦濾波的列控車載設備測速測距算法的應用性能,利用接口型式實驗采集轉速傳感器測量數據和MATLAB 仿真軟件進行仿真分析。分別取型式實驗中列車正常運行場景和制動時輪對滑行場景下的測速測距系統監測數據,分析列控車載設備測速測距算法是否滿足準確性和安全性要求。

4.1 列車正常運行場景

選取列車進、出站及在區間中正常運行的場景,將每個輪對現場采集速度疊加噪聲后作為測量值、參考值(即真實值)分別與采用自適應參數列車運動模型得到的濾波結果進行對比,如圖2所示。圖中:速度的采樣周期取0.05 s;采樣點個數為0 時表示采樣開始時刻,采樣點個數與采樣周期的乘積表示當前采樣時刻。由圖2可知:采樣點4 000個左右,即采樣開始200 s時(4 000×0.05 s=200 s),列車加速度開始發生較大變化,從加速度為零變為加速度為正;采樣點11 500個左右,即采樣開始575 s 時(11 500×0.05 s=575 s),列車加速度開始從零變負;濾波結果和速度參考值變化趨勢一致,這表示自適應參數列車運動模型能有效濾除噪聲;濾波結果對加速度變化敏感,在速度變化的過程中也不存在滯后現象,滿足測速測距對準確性的要求。

圖2 列車正常運行仿真結果

4.2 列車制動時輪對滑行場景

以涂灑減磨液的形式制造滑行場景。考慮列車制動過程中可能發生的2 種情況:一是單個轉速傳感器所在輪對發生滑行;二是2 副輪對同時發生滑行。分情況測量2 個轉速傳感器所在輪對的運行速度數據,采用基于自適應聯邦濾波的測速測距算法模型進行仿真。單個轉速傳感器所在輪對發生滑行時的算法仿真結果如圖3所示。圖中:傳感器1 和傳感器2 分別為2 個安裝在不同輪軸上的轉速傳感器;紅線和綠線分別表示不同轉速傳感器所在輪對的測量值通過列車運動模型得到的估計速度;黑線表示列車運行速度。由圖3可知:采樣點750 個左右,即采樣開始37.5 s 時(750×0.05 s=37.5 s)可視之為分界點;37.5 s之前列車正常運行,信息分配系數權重相當,這表明聯邦濾波融合結果同時考慮了2個轉速傳感器所在輪對的測量速度;37.5 s之后,利用空轉/滑行檢測方法檢測出傳感器2 所在輪對發生滑行,于是減少對其信息分配系數的權重,聯邦濾波融合結果主要基于傳感器1所在輪對的速度,因此列車運行速度基本沒有受到滑行的影響,滿足列控車載設備測速測距對準確性和安全性的要求。

圖3 單個轉速傳感器所在輪對滑行場景下的仿真效果

2 個轉速傳感器所在輪對同時發生滑行時的算法仿真結果如圖4所示。由圖4可知:采樣點12 個左右,即采樣開始1.6 s 時(12×0.05 s=1.6 s),利用滑行檢測方法檢測出2 個轉速傳感器所在輪對同時發生滑行,此時主要由列車速度校正模型計算列車運行速度;采樣點40個左右,即采樣開始2.0 s時(40×0.05 s=2.0 s),傳感器1 所在輪對的速度超過了校正后的列車運行速度,此時轉速傳感器所在輪對黏著恢復,校正速度由4.2 m·s?1調整為4.7 m·s?1,在滑行持續校正時間內繼續利用滑行速度校正模型計算校正速度;采樣點150 個左右,即采樣開始7.5 s時(150×0.05 s=7.5 s),2 個轉速傳感器所在輪對的速度均為零(滑行期內若2 個轉速傳感器的測量速度為零,則視之為車輪抱死),此時列車還在減速運行,為保證測速測距的安全,在原有校正速度的基礎上,取較小的減速度計算列車運行速度直至列車完全停止,整個速度計算過程滿足列控測速測距對準確性和安全性的要求。

圖4 2個轉速傳感器所在輪對滑行場景下仿真效果

可見,無論是列車正常運行場景還是制動時輪對滑行場景,無論滑行的是單個轉速傳感器所在軸的輪對還是2 副輪對,本文算法均能夠自適應地跟蹤列車運行速度,通過實時檢測轉速傳感器所在軸的空轉/滑行情況,在安全導向的前提下得到符合車輛運行實際情況的列車運行速度,兼顧了列控車載測速測距系統對準確性和安全性的要求。

5 結 論

(1)提出1 種由自適應參數列車運動模型、空轉/滑行時列車速度校正模型和基于聯邦濾波的融合算法3 部分組成的列控車載設備測速測距算法,僅使用2 個轉速傳感器,即可實現對轉速傳感器所在輪對當前測量速度的實時、有效修正,克服了以往測速測距系統中轉速傳感器因所在輪對空轉/滑行引起的速度測量精度下降問題,更符合列車實際運行情況。

(2)建立自適應參數列車運動模型,2 個轉速傳感器采集得到的速度脈沖信號經過預處理后,可分別得到其所在輪對的估計速度和加速度,在列車運行加速度變化較大時,仍能夠實現對轉速傳感器所在輪對速度的實時追蹤。

(3)設計速度校正模型并將其視為1 個虛擬傳感器,能夠在檢測確定2 個轉速傳感器所在輪對空轉/滑行情況的基礎上,通過校正速度與測量速度的不斷比對與修正,與自適應參數列車運動模型共同描述列車運行的動態過程。

(4)利用基于聯邦濾波的融合算法,設計動態調整信息分配系數方案,實現自適應參數列車運動模型和列車速度校正模型的計算結果的有效融合,使得到的運行速度和走行距離能夠反映出列車實際運行狀態。

猜你喜歡
檢測模型
一半模型
“不等式”檢測題
“一元一次不等式”檢測題
“一元一次不等式組”檢測題
“幾何圖形”檢測題
“角”檢測題
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
3D打印中的模型分割與打包
小波變換在PCB缺陷檢測中的應用
主站蜘蛛池模板: 精品天海翼一区二区| 国产成在线观看免费视频| 国产va免费精品观看| 精品视频一区二区观看| 亚洲区欧美区| 国产在线观看第二页| 九九热这里只有国产精品| 99久久精品免费看国产免费软件 | 亚洲av中文无码乱人伦在线r| 日本不卡在线播放| 波多野结衣在线一区二区| 成人va亚洲va欧美天堂| 四虎精品黑人视频| 国产免费好大好硬视频| 欧美特级AAAAAA视频免费观看| 成年午夜精品久久精品| 看国产毛片| 青青操国产| 青草视频久久| 国产69精品久久| 国产高潮视频在线观看| 久久国产成人精品国产成人亚洲| 中文字幕 日韩 欧美| 国产jizzjizz视频| 国产真实乱了在线播放| 香蕉在线视频网站| 扒开粉嫩的小缝隙喷白浆视频| 亚洲中文字幕国产av| 成人午夜免费观看| 午夜视频在线观看免费网站| 精品国产免费人成在线观看| 激情综合激情| 日本黄色不卡视频| 国产欧美亚洲精品第3页在线| 777国产精品永久免费观看| 亚洲av成人无码网站在线观看| 国产精品露脸视频| 亚洲国产AV无码综合原创| 欧美精品1区2区| 国产香蕉一区二区在线网站| 成人av手机在线观看| 美女国产在线| 毛片在线播放a| 国产婬乱a一级毛片多女| 夜夜拍夜夜爽| 国产在线观看人成激情视频| 日本成人福利视频| 亚洲精品波多野结衣| 国产欧美日韩va另类在线播放| 91国内视频在线观看| 国产精品视频白浆免费视频| 中国黄色一级视频| 日韩东京热无码人妻| 欧美成a人片在线观看| 欧美一区二区丝袜高跟鞋| 毛片基地美国正在播放亚洲| 99在线视频精品| 免费aa毛片| 精品欧美一区二区三区久久久| 中文字幕亚洲精品2页| 久久鸭综合久久国产| 波多野结衣一区二区三区88| 亚洲区视频在线观看| 无码av免费不卡在线观看| 亚洲福利视频一区二区| 1024你懂的国产精品| 国产凹凸视频在线观看| 国产丝袜精品| 日韩欧美国产精品| 欧美在线黄| 久久久成年黄色视频| 久久不卡国产精品无码| 嫩草影院在线观看精品视频| 国产精品xxx| 麻豆精品久久久久久久99蜜桃| 午夜不卡福利| 日韩精品成人在线| 九色在线视频导航91| 91成人在线观看| 老熟妇喷水一区二区三区| 国产精品漂亮美女在线观看| 午夜精品久久久久久久无码软件|