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

基于Hilbert-Huang變換的制動意圖聚類識別

2016-04-12 01:29:08唐先智楊樹軍
汽車工程 2016年5期
關鍵詞:踏板信號

唐先智,王 波,楊樹軍,馬 雷

(燕山大學車輛與能源學院,秦皇島 066004)

2016094

基于Hilbert-Huang變換的制動意圖聚類識別

唐先智,王 波,楊樹軍,馬 雷

(燕山大學車輛與能源學院,秦皇島 066004)

為了進一步提高電動汽車再生制動系統對駕駛員制動意圖的識別準確率,從而使電動汽車能在制動的過程中回收更多的能量,提出了基于Hilbert-Huang變換(HHT)的電動汽車制動意圖聚類識別方法。建立了HHT的數學模型,基于HHT在時頻域中進一步探尋中等制動和平緩制動兩種制動意圖下制動踏板行程信號的特征。建立了制動踏板行程信號特征提取模型,運用Hilbert局部邊際能量譜得到局部特征能量,從而對信號特征進行提取,獲取信號的特征向量。建立基于模糊C均值聚類算法的制動意圖識別模型,并分別進行了離線實驗和實時實驗。結果表明所提出的基于HHT的制動意圖模糊C均值聚類識別方法能更好地分辨中等制動和平緩制動意圖,提高了識別準確率,并具有較好的實時性。

電動汽車;制動意圖;Hilbert-Huang變換;聚類識別

前言

電動汽車制動控制策略對再生制動和機械制動在制動過程中使用比例的分配依據是駕駛員的制動意圖[1-2]。駕駛員制動意圖識別準確與否直接影響到電動汽車再生制動的能量回收率。因此,如何能夠準確地識別駕駛員制動意圖是電動汽車再生制動技術需要解決的重要問題[3-4]。國外制動意圖識別技術主要應用在電動汽車再生制動系統和制動輔助系統上,研究開展較早,很多汽車企業已經把制動意圖識別技術應用到量產車型中。但是,國外各大汽車企業都沒有公開關于電動汽車制動意圖識別的核心技術細節[5-9]。

國內學者針對電動汽車駕駛員在制動過程中的制動意圖識別方法進行了相關研究。文獻[10]中以制動踏板角速度、車速方差和車輛減速度作為制動意圖識別的參數,利用模糊邏輯對駕駛員制動意圖進行辨識。文獻[11]中確定了以制動踏板角速度識別制動意圖,以制動管路壓力表征目標制動強度。采用模糊辨識算法進行制動意圖識別。文獻[12]和文獻[13]中以制動踏板位移為制動意圖識別的參數。文獻[14]和文獻[15]中以制動踏板位移及變化率并結合車速識別駕駛員的制動意圖。文獻[16]中以制動踏板行程及其變化率為識別參數,運用數理統計法制定識別參數的隸屬函數和模糊推理規則并用神經網絡對識別參數的隸屬函數進行了優化,通過模糊算法可較準確地識別制動意圖。由上述可見,國內制動意圖識別主要是以制動踏板行程及其變化率等時域參量為識別參數,通過邏輯推斷或模糊推理識別駕駛員的制動意圖。這種識別方法平均識別準確率在90%左右,對于緊急制動和中等制動辯識率較高,在95%左右。但在辨識平緩制動和中等制動的過程中,駕駛員對制動踏板的操作特征不明顯,若直接用制動踏板行程及其變化率等時域參數表征制動意圖,則臨界狀態不易區分,易受到測量誤差的干擾,識別準確率低于90%有時甚至低于80%。因此應對識別參數信號進行進一步挖掘,探尋識別參數更深層次的特征,提高識別分辨率。

于是,本文中提出了基于Hilbert Huang變換(Hilbert-Huang transform,HHT)的電動汽車制動意圖識別方法。首先運用HHT對制動踏板行程信號在時頻域內進行數據挖掘,得到制動踏板行程信號在時頻域的特征,然后采用聚類識別算法對制動意圖進行識別,從而進一步提高制動意圖識別的分辨率和準確率。

1 Hilbert-Huang變換

對于制動踏板行程這樣的非平穩數據信號,直接進行Hilbert變換得到的結果很大程度上失去了原有的物理意義。因此數據信號應首先進行經驗模態分解(empirical mode decomposition, EMD),將數據信號分解為平穩的固有模態函數(intrinsic mode function, IMF),即將信號經過經驗模態分解(EMD),使真實存在的不同尺度波動或趨勢逐級分解開來,產生一系列具有不同特征尺度的數據序列,每個序列稱為一個固有模態函數(IMF)。之后再進行Hilbert變換便能夠得到能量在頻率以及時間上的分布規律,從而表征信號的局部特征。

1.1 EMD原理與算法

IMF可以直接進行Hilbert變換得到能量在頻率以及時間上的分布規律。而一般的信號往往是復雜的信號,并不滿足IMF的定義。因此文獻[17]和文獻[18]中假設:任何復雜信號都是由許多簡單的、不同的、非正弦的IMF分量(可以是線性的,也可以是非線性的)組成。根據此假設,提出了EMD算法(即Huang變換),該算法有兩個重要作用:一個是去除疊加模態,另一個是使波形更對稱。這是HHT變換的關鍵環節。

首先找出初始信號X(t)上所有的極值點,然后分別對所有極大值點和極小值點進行三次樣插差值,擬合出初始信號的上包絡線Xmax(t)和下包絡線Xmin(t)。取上下包絡線的均值m1(t):

m1(t)=[Xmax(t)+Xmin(t)]/2

(1)

用初始信號減去包絡線均值得到h1(t):

h1(t)=X(t)-m1(t)

(2)

如果h1(t)不滿足IMF所需的條件,則將h1(t)作為初始信號,重復以上步驟得:

h11(t)=h1(t)-m11(t)

(3)

式中:m11(t)為h1(t)的上下包絡線均值。若h11(t)仍不滿足IMF的條件,則重復上述算法k次得:

h1k(t)=h1(k-1)(t)-m1k(t)

(4)

h1k(t)篩選過程終止判據為

(5)

式中SD為相鄰兩個處理結果的標準差。SD的取值不能過小也不能過大。SD過小會造成得到的IMF頻率調制信號幅值恒定,失去信號特征;SD過大會使終止條件過于寬松,篩選結果無法滿足IMF的條件。經驗表明,要保證IMF既是線性的,又具有穩定性并且具有物理特征,SD應取0.2~0.3之間的數值比較合適。

如果式(4)中的h1k(t)符合終止判據的要求,則h1k(t)為1階IMF,用c1(t)表示。用初始信號X(t)減去1階IMF得到信號殘差r1(t):

r1(t)=X(t)-c1(t)

(6)

r1(t)可作為新信號再重復以上的EMD過程,經多次迭代可得所有信號殘差ri(t):

ri(t)=ri-1(t)-ci(t),i=2,3,…,n

(7)

當殘差rn(t)成為單調函數時,即不可能再從中提取IMF分量時,便可終止EMD過程。可見,初始信號X(t)可由n階IMF和殘差rn(t)組成:

(8)

1.2 Hilbert變換與Hilbert譜

對式(8)中的每個IMFci(t)做Hilbert變換:

(9)

構造解析信號:

zi(t)=ci(t)+jH[ci(t)]=ai(t)ejΦi(T)

(10)

(11)

(12)

式中:ai(t)為幅值函數;Φi(t)為相位函數。

根據式(12)可求得瞬時頻率:

(13)

初始信號可表示為

(14)

式中Re表示取實部。

根據式(14)可定義Hilbert譜的表達式:

(15)

式中Hilbert譜H(ω,t)描述了信號幅值隨時間和頻率的變化關系。

1.3 制動踏板行程信號的Hilbert-Huang變換

將制動踏板行程的數據信號輸入計算機,運用MATLAB軟件在計算機中編寫Hilbert-Huang變換的程序。EMD算法(即Huang變換)程序流程如圖1所示。Hilbert變換的流程如圖2所示。

圖3和圖4分別為中等制動和平緩制動時制動踏板行程信號的實驗數據。平緩制動時制動減速度小于1.6m/s2,制動踏板行程一般小于全行程的30%,這種制動的制動強度小,制動緊急程度也很小。對于電動汽車,在這種制動工況下要盡量使用再生制動,回收盡可能多的制動能量。中等制動時制動減速度介于1.6~3m/s2,制動踏板一般介于30%~60%,制動踏板行程變化率一般介于每秒65%~140%,這種制動的制動強度中等,制動緊急程度也為中等。電動汽車在這種制動工況下要適當加大機械制動的比例,確保制動的安全性。各截取如圖3和圖4狹長橢圓圈內所示的兩組實驗數據的一段進行Hilbert-Huang變換。所截取的實驗數據如圖5和圖6所示。分別對截取的兩段實驗數據運用EMD算法進行Huang變換,分解得到的初始數據信號的固有模態函數(即IMF分量)如圖7和圖8所示。分別對中等制動和平緩制動數據信號的IMF分量進行Hilbert變換,得到兩段數據信號的Hilbert譜,如圖9和圖10所示。從圖9和圖10可以看出,中等制動時與平緩制動時的制動踏板行程信號的幅值在頻率和時間上的分布特征有明顯不同。通過適當算法進行特征提取,便可對兩種制動模式進行辨識。

2 特征提取與聚類識別

從前面的分析可以看出,通過對不同制動意圖下制動踏板行程信號的Hilbert-Huang變換,進一步挖掘了該信號的時頻域特征。下面從信號能量的角度,對中等制動和平緩制動時制動踏板行程信號的特征進行提取,并運用模糊C均值聚類算法對這兩種制動意圖進行識別。

2.1 特征提取

由于制動踏板行程信號為電壓信號,因此可以把制動踏板行程信號的平方值當作信號的能量。對它進行Hilbert-Huang變換,便得到Hilbert能量譜,記為H2(ω,t)。根據能量守恒定律,Hilbert-Huang變換前后制動踏板行程信號能量守恒,如式(16)所示。使H2(ω,t)對時間進行積分,可得制動踏板信號的Hilbert邊際能量譜,它表達了每個頻率在整個時間長度內能量的累積,描述了累積的能量隨頻率分布的情況,如式(17)所示。同理,求出制動踏板行程信號的每個IMF分量的Hilbert能量譜Hi2(ω,t),使Hi2(ω,t)對時間進行積分,可得到每個IMF分量的Hilbert邊際能量譜,即制動踏板信號的Hilbert局部邊際能量譜Ei(ω),如式(18)所示。使Ei(ω)對頻率進行積分,得到局部特征能量,如式(19)所示。由此可得制動踏板行程信號的特征向量,如式(20)所示。特征向量的每個元素為每個IMF分量的局部特征能量與總特征能量的比值。

(16)

(17)

(18)

(19)

T=[E1/E,E2/E,…,En/E]

(20)

2.2 制動意圖識別

模糊C均值聚類算法的核心思想是找到樣本的聚類中心,使樣本到聚類中心的加權距離平方和最小。該算法的目標函數為

(21)

式中:μij為第j個樣本隸屬于第i類制動意圖的程度;dij為樣本點xj到第i類制動意圖聚類中心的歐式距離,即dij(xj,zi)=‖xj-zi‖;Z為聚類中心,即Z=(z1,z2,…,zc);m為加權指數;C為制動意圖的分類數,表示將樣本集{x1,x2,…,xn}分為C類制動意圖;U為初始隸屬度矩陣。

模糊C均值聚類算法的迭代方程如式(22)和式(23)所示,迭代終止條件如式(24)所示。

(22)

(23)

‖Z(p)-Z(p+1)‖<ε

(24)

式中:C為制動意圖分類數,2≤C≤n,n為樣本個數;ε為迭代停止閾值;Z(0)為初始聚類中心;Z(p)為第p次迭代的聚類中心,p為迭代次數;m為加權指數。

3 實驗驗證

3.1 離線驗證

選取中等制動和平緩制動的制動踏板行程信號樣本,均取各自的前3階IMF分量構造特征向量。以50個標準樣本特征向量中每個元素的均值為初始聚類中心,如表1所示。加權指數m取2,聚類算法迭代終止閾值為10-6。

表1 樣本的初始聚類中心

中等制動和平緩制動各選擇10組制動踏板行程信號作為檢測樣本,兩種制動意圖共計20個檢測樣本。檢測樣本1~10為中等制動,11~20為平緩制動。在MATLAB軟件環境下,對基于HHT的制動意圖模糊C均值聚類識別算法進行離線驗證,識別結果如表2所示。表3為基于模糊推理的制動意圖識別結果,識別方法見文獻[16]。

表2 基于HHT的制動意圖聚類識別算法識別結果

表3 制動意圖模糊推理識別算法識別結果

從表2中可以看出,識別算法將制動踏板行程信號的檢測樣本分別歸入了兩類制動意圖。表中樣本編號下面標記橫線的樣本為識別錯誤的樣本,樣本13和樣本5識別錯誤,識別準確率為90%。從表3中可以看出,采用相同的檢測樣本,樣本1,5,6,13和14識別錯誤,識別準確率為75%。

模糊推理識別法識別參數模糊化后的隸屬函數如圖11和圖12所示。模糊推理規則如表4所示。

表4 模糊推理規則

運用模糊推理法識別正確的檢測樣本,如樣本11,12,15,16,17,18,19和20為平緩制動樣本。制動踏板行程變化率相對于制動踏板行程信號的測量難度更大,容易對識別結果產生干擾。其中樣本11,12,17,18和20的制動踏板行程為“小”,這樣從模糊推理規則可以看出無論制動踏板變化率為“小”、“中”或“大”,識別結果都為平緩制動。制動踏板變化率信號不會對識別結果產生負面影響,因此這5個樣本相對容易識別。樣本15,16和19也為平緩制動,制動踏板行程為“中”,制動踏板行程變化率為“小”并且小于0.4,此時制動踏板行程變化率很小,信號穩定,易于測量并且此時制動踏板行程變化率隸屬于“小”的隸屬度為100%,因此對于這3個樣本模糊識別方法也能夠準確識別。

對于識別錯誤的樣本1,5和6,其制動踏板行程為“中”,制動踏板行程變化率為“小”,識別結果為平緩制動,這與實際不符。從模糊推理規則可以看出,只有制動踏板行程為“中”時,推理的結果才有可能是中等制動。由于隸屬函數與模糊推理規則是通過大量的數理統計并經過反復優化得到,相對穩定,因此導致識別結果錯誤的原因應是制動踏板行程變化率信號測量有誤,經模糊化后為“小”,實際應為“中”。可見,為更好地區分平緩制動和中等制動而引入的識別參數制動踏板行程變化率,由于其測量難度相對較大,測量值會有誤差,有時反而會影響模糊推理結果。樣本2,3,4,7,8,9,10,13和14的識別結果可以同理分析。文獻[16]中所做的后續研究工作在文獻[19]中加入了另一個識別參數車速,但仍沒有解決上述問題。

基于HHT的制動意圖模糊C均值聚類識別算法的識別參數只有制動踏板行程一個參數,并且易于測量。排除了制動踏板行程變化率這一不易測量準確參數的干擾,并省去了制動踏板角速度傳感器,節省了成本。把對制動強度和制動緊急程度的辨識引入頻域范圍,在時頻域內進一步挖掘制動踏板行程信號的時頻特征,運用Hilbert能量譜提取信號的特征向量并進行聚類識別。識別準確率提高了15%。

基于HHT的制動意圖模糊C均值聚類識別方法中的EMD算法,容易受到信號中干擾成分的影響,會導致信號極值點的偏差,從而對極值點的包絡擬合也產生偏差。尤其是信號有間斷或脈沖時,容易產生模式混疊現象。因此在實際對信號進行EMD分解的過程中,同一頻率的成分會被分到多個IMF中,也可能在一個IMF中包含多個分量的成分。從而使每個IMF局部特征能量改變,導致特征向量和初始聚類中心的偏移,使某些識別結果產生偏差。如表2所示,基于HHT的制動意圖聚類識別算法識別結果中仍有兩個樣本5和13識別錯誤。今后的工作是進一步提高EMD算法的抗干擾能力,進一步提高制動意圖的識別準確率。

3.2 實時驗證

實時驗證采用的實驗用車為解放牌混合動力公交車,其型號為CA6120URH1。制動踏板行程信號的在線采集、數據處理和制動意圖的在線識別則選用與MATLAB軟件兼容的通用控制器dSPACE/MicroAutoBox1401/1504。由于CA6120URH1車的CAN總線上,制動踏板行程為開關量,即其數值為0或1,不能直接用來進行制動意圖識別,因此須自制制動踏板行程傳感器對制動踏板行程信號進行采集,并輸入到MicroAutoBox中進行在線數據處理和制動意圖在線識別。在線識別之前須通過程序下載接口將基于HHT的制動意圖模糊C均值聚類識別算法的程序代碼下載到MicroAutoBox中。制動踏板行程實時信號如圖13所示,它是通過駕駛員按要求(即一次中等制動,一次平緩制動)對制動踏板進行操作并在線采集得到的。圖13中標記了制動踏板行程實時信號與其相應的制動意圖。制動意圖的在線識別結果如圖14所示,可見該算法可以準確地在線識別駕駛員的制動意圖。圖15示出制動意圖在線識別響應時間,可以看出,從進入制動模式到識別出第一個中等制動意圖,在線識別時間為0.4s,符合中等制動意圖與平緩制動意圖識別的響應速度要求。

4 結論

(1)基于HHT的制動意圖模糊C均值聚類識別算法可在時頻域內對制動踏板行程信號的特征進行進一步的挖掘,提取其信號的特征向量,能對中等制動和平緩制動進行準確地辨識。離線實驗證明,本算法相比于單純在時域內的模糊推理算法,減少了識別參數,減小了參數測量誤差對識別結果的影響,識別準確率提高了15%,進一步提高了對中等制動和平緩制動意圖的分辨能力。

(2)基于HHT的制動意圖模糊C均值聚類識別算法具有較好的實時性。實時實驗證明,本算法可實時對制動意圖進行準確的在線識別,為制動意圖識別技術在電動汽車再生制動領域中更好的實時應用奠定理論基礎。

[1] 初亮, 王彥波, 姚亮,等. 制動能量回收系統的制動力矩協調控制仿真[J]. 華南理工大學學報(自然科學版), 2014, 42(4):137-142.

[2] HAN J, PARK Y, PARK Y. Cooperative Regenerative Braking Control for Front-wheel-drive Hybrid Electric Vehicle Based on Adaptive Regenerative Brake Torque Optimization Using Under-steer Index[J]. International Journal of Automotive Technology, 2014, 15(6):989-1000.

[3] 趙軒, 馬建, 汪貴平. 基于制動駕駛意圖辨識的純電動客車復合制動控制策略[J]. 交通運輸工程學報, 2014(4):64-75.

[4] 王玉海, 宋健, 李興坤. 駕駛員意圖與行駛環境的統一識別及實時算法[J]. 機械工程學報, 2006, 42(4):206-212.

[5] WASEKURA M, WANG C, LORENZ R D. A Transient Core Loss Analysis of Multiple-gap Inductor Designed for the 2010 Prius[C].Energy Conversion Congress and Exposition (ECCE), 2014 IEEE:2177 - 2181.

[6] MILLER F P, VANDOME A F, Mcbrewster J. Honda EV Plus[M].Alphascript Publishing, 2010.

[7] HIROSE T, TANIGUCHI T, HATANO T, et al. A Study on the Effect of Brake Assist Systems (BAS)[J]. SAE International Journal of Passenger Cars: Mechanical Systems, 2009(1):729-735.

[8] SCHICK B, BüTTNER R, BALTRUSCHAT K, et al. Evaluation Methods for the Function and Quality of Driver Assistance Systems with Active Brake Control[J]. ATZ Worldwide, 2007, 109(5):14-18.

[9] Bernhard Eichhorn M, Steffen K?nig, Thorsten Ullrich. Electronic Brake Control for Greater Active Safety[J]. ATZ worldwide, 2014, 116(9):50-53.

[10] 孫逸神. 基于模糊邏輯的制動意圖離線識別方法研究[J]. 北京汽車, 2009(6): 21-23.

[11] 張元才, 余卓平, 徐樂,等. 基于制動意圖的電動汽車復合制動系統制動力分配策略研究[J]. 汽車工程, 2009, 31(3):244-249.

[12] 王英范,寧國寶,余卓平.乘用車駕駛員制動意圖識別參數的選擇[J].汽車工程, 2011,33(3): 213-216.

[13] 李玉芳,吳炎花.電-液復合制動電動汽車制動感覺一致性及實現方法[J].中國機械工程, 2012,23(4):488-492.

[14] 林逸,沈沉,王軍,等.汽車線控制動技術及發展[J]. 汽車技術, 2005(12):1-3.

[15] 馬朝永,化北,王震,等.電子制動踏板感覺模擬器研究[J].電子測量技術, 2011,34(7):28-32.

[16] 王慶年,孫磊,唐先智,等. HEV 制動意圖識別的研究[J].汽車工程, 2013,35(9): 769-774.

[17] HUANG N E, SHEN Z, LONG S R, et al. The Empirical Mode Decomposition and the Hilbert Spectrum for Nonlinear and Non-stationary Time Series Analysis[J]. Proceedings of the Royal Society A Mathematical Physical & Engineering Sciences, 1998, 454:903-995.

[18] BATTISTA B, KNAPP C, MCGEE T, et al. Application of the Empirical Mode Decomposition and Hilbert-Huang Transform to Seismic Reflection Data[J]. Geophysics, 2007, 72(2):H29-H37.

[19] 王慶年, 王俊, 陳慧勇,等. 混合動力車輛中的加速與制動意圖識別[J]. 吉林大學學報(工學版), 2014, 44(2):281-286.

Cluster Identification of Braking Intention Based on Hilbert-Huang Transform

Tang Xianzhi, Wang Bo, Yang Shujun & Ma Lei

CollegeofVehiclesandEnergy,YanshanUniversity,Qinhuangdao066004

For further enhancing identification accuracy rate on the braking intention of driver so as to make electric vehicle recover more energy, a cluster identification method of braking intention for electric vehicles is proposed based on Hilbert-Huang transform (HHT). A mathematical model for HHT is built, and based on HHT the brake pedal signal features under both moderate and gentle braking intentions are further explored in time/frequency domains. An extraction model for brake pedal travel signal feature is set up, Hilbert marginal energy spectrum is used to acquire local feature energy, and hence the signal features are extracted with corresponding eigenvectors obtained. Then the braking intention identification model is established based on fuzzy C-means clustering algorithm, and both offline and real-time experiments are performed respectively. The results show that the identification method proposed can better distinguish moderate and gentle braking intentions, raise identification rate with better real-time performance.

electric vehicles; braking intention; Hilbert-Huang transform; cluster identification

*國家自然科學基金青年基金(51505414)、河北省高等學校科學研究項目(Z2015081)和燕山大學博士基金(B794)資助。

原稿收到日期為2015年12月24日,修改稿收到日期為2016年3月3日。

猜你喜歡
踏板信號
單踏板不可取
車主之友(2022年6期)2023-01-30 07:58:16
信號
鴨綠江(2021年35期)2021-04-19 12:24:18
完形填空二則
淺談延音踏板在鋼琴演奏中的用法
黃河之聲(2019年23期)2019-12-17 19:08:43
淺談汽車制動踏板的型面設計
孩子停止長個的信號
淺談鋼琴踏板的運用
黃河之聲(2017年13期)2017-01-28 13:30:17
論鋼琴踏板的正確使用
基于LabVIEW的力加載信號采集與PID控制
一種基于極大似然估計的信號盲抽取算法
主站蜘蛛池模板: 一级香蕉视频在线观看| 国产一区二区三区日韩精品| 国产理论精品| 久久精品国产91久久综合麻豆自制| 色综合久久无码网| 欧美精品在线免费| 国产精品女同一区三区五区| 在线色综合| 伊人久久大香线蕉影院| 亚洲一级毛片在线观| 色久综合在线| 亚洲天堂日韩在线| 国产精品污污在线观看网站| 免费a级毛片18以上观看精品| 18黑白丝水手服自慰喷水网站| 97在线免费视频| 一级做a爰片久久免费| 精品视频免费在线| 亚洲人成色在线观看| YW尤物AV无码国产在线观看| 被公侵犯人妻少妇一区二区三区 | 一区二区三区国产精品视频| 九色综合视频网| 国产69精品久久久久孕妇大杂乱| 国产亚洲精品自在线| 国产va欧美va在线观看| 亚洲美女一区| 91精品免费高清在线| 亚洲综合精品香蕉久久网| 久久青草免费91线频观看不卡| 久久精品国产在热久久2019| 国产乱码精品一区二区三区中文| 午夜视频免费试看| 人妻夜夜爽天天爽| 中文字幕色站| 日韩无码视频网站| 97国产精品视频人人做人人爱| 国产精品成人免费视频99| 激情在线网| 欧洲亚洲一区| 国产在线观看成人91| 欧美一级在线| 99伊人精品| 亚洲日韩精品无码专区| 国产va在线观看免费| 狠狠色成人综合首页| 国产白浆视频| 久久久久青草线综合超碰| 日韩a级毛片| 欧美色视频日本| 97国产在线播放| 中文字幕va| 欧美笫一页| 四虎免费视频网站| 男女性色大片免费网站| 欧美一级黄色影院| 无码一区中文字幕| 人妻21p大胆| 欧美精品另类| 成人免费午间影院在线观看| 亚洲人成网址| 欧美在线综合视频| 无码精油按摩潮喷在线播放| 91探花在线观看国产最新| 五月婷婷综合在线视频| 久热精品免费| 亚洲免费三区| 国产91无码福利在线| 国产成人av一区二区三区| 欧美成a人片在线观看| 中文字幕不卡免费高清视频| 国产成人亚洲综合A∨在线播放| 国产精品林美惠子在线观看| 九九九精品成人免费视频7| 国精品91人妻无码一区二区三区| 日本午夜精品一本在线观看| 国产欧美日韩资源在线观看| 人妻中文久热无码丝袜| 国产成人在线小视频| 亚洲天堂伊人| 亚洲色成人www在线观看| 欧洲亚洲欧美国产日本高清|