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

基于自適應(yīng)時頻濾波的變轉(zhuǎn)速齒輪故障特征提取

2018-05-28 02:55:02陳向民晉風華李錄平
振動與沖擊 2018年10期
關(guān)鍵詞:振動故障信號

陳向民,張 亢,晉風華,李錄平

(長沙理工大學 能源與動力工程學院,長沙 410076)

齒輪是旋轉(zhuǎn)機械中的重要動力傳輸部件。長時間、高負荷地運行容易引起齒輪出現(xiàn)局部故障。齒輪出現(xiàn)局部故障時,其振動信號中會產(chǎn)生以轉(zhuǎn)頻及其倍頻為調(diào)制頻率的調(diào)幅調(diào)頻信號[1]。如何從原始振動信號中將包含齒輪局部故障信息的調(diào)幅調(diào)頻信號提取出來是診斷齒輪故障的關(guān)鍵。

對于恒定轉(zhuǎn)速下的齒輪故障信號分離,常用的信號分析方法主要有共振解調(diào)[2]、小波分析[3-4]、EMD/EEMD[5-6]、LMD[7-8]等。這類方法的共同特點是采用類似于帶通濾波的方式對信號進行分析,將包含故障信息的信號成分分離出來,對于恒定轉(zhuǎn)速下的故障信號分析取得了較好的效果。但對于頻率大范圍波動的非平穩(wěn)信號(如機器的啟停、工況的改變等),若直接采用該類方法進行故障信號提取,則分析結(jié)果不太理想[9]。

時頻濾波通過先采用時頻分析方法將信號展開為時間-頻率的聯(lián)合函數(shù)(即時頻分布),再在時頻域?qū)π盘栠M行濾波分析,可有效提取時變非平穩(wěn)信號成分,已成為非平穩(wěn)信號分析的一種重要方法。常用于時頻濾波的時頻分析方法有短時傅里葉變換,Wigner-Ville變換、Gabor變換[10],S變換[11]等,其中,S變換因不含交叉項,且保留了信號的絕對相位信息,因而,近年來被用于非平穩(wěn)信號的時頻濾波分析[12-14]。

時頻濾波的基本原理是根據(jù)分析信號的頻率變化特點,在時頻域選取合適的時頻區(qū)域(即時頻濾波器)進行濾波分析。對于變轉(zhuǎn)速下的齒輪故障信號,其時頻濾波器的設(shè)計主要包括中心頻率和帶寬兩個部分,其中,中心頻率的選取是關(guān)鍵[15]。而對于齒輪故障信號的時頻濾波器的中心頻率一般選取為齒輪的嚙合頻率。因此,如何精確的將嚙合頻率曲線從齒輪故障信號中提取出來對其時頻濾波器的設(shè)計至關(guān)重要。CPP方法是Candes等[16]提出的一種非平穩(wěn)信號頻率曲線估計方法,具有較高的頻率估計精度,且具有很強的抗噪能力,近年來被用于機械故障診斷,取得了較好的效果[17-18]。

基于上述分析,本文將CPP方法與S變換相結(jié)合,提出了基于CPP與S變換的自適應(yīng)時頻濾波方法,并將其應(yīng)用于齒輪故障特征提取中。算法仿真和應(yīng)用實例表明,自適應(yīng)時頻濾波方法可根據(jù)齒輪故障信號的頻率變化特點,自適應(yīng)地改變中心頻率和帶寬,能有效分離出齒輪故障信號成分,并提取齒輪故障特征,且濾取的信號無相位畸變。因此,該方法非常適合于變轉(zhuǎn)速下的非平穩(wěn)信號分析。

1 CPP算法

CPP算法采用的多尺度線調(diào)頻基元函數(shù)庫為

D(haμbμI)={haμbμI(t)}={KaμbμIe[-i(aμt+bμt2)]LI(t)}

(1)

式中:D表示基元函數(shù)庫;haμbμI(t)為多尺度線調(diào)頻基元函數(shù);N為分析信號長度;I表示動態(tài)分析時間段,I=[kN2-j~(k+1)N2-j];k為動態(tài)時間段的序號,k=0,1,…,2j-1;j表示分析尺度系數(shù),j=0,1,…,log2N-1;KaμbμI為歸一化系數(shù),滿足‖haμbμI‖=1;aμ和bμ分別表示頻率偏置系數(shù)和調(diào)頻率,且滿足aμ+2bμt

采用多尺度線性調(diào)頻基函數(shù)對分析信號進行逐段投影,并計算每個時間分析段I內(nèi)的投影系數(shù)和對應(yīng)的線調(diào)頻基元函數(shù)。當分析信號與多尺度線性調(diào)頻基函數(shù)的相似性越高,其投影系數(shù)也就越大,此時,基元函數(shù)對應(yīng)的能量也就越大。因此,需尋求一種動態(tài)連接算法∏,使得所連接基元函數(shù)對應(yīng)的信號在整個分析時間內(nèi)的總能量最大,且連接算法∏應(yīng)覆蓋整個分析時間段,不重疊,即

(2)

此時,對應(yīng)的投影系數(shù)集合和基元函數(shù)集合分別為

(3)

CPP方法中連接算法∏的連接步驟如下:

(1) 初始化。以i表示分析時間段序號,Edi表示第i個分析時間段之前分解信號的總能量,li表示連接到第i個分析時間段的前一分析時間段序號,Eei表示第i個分析時間段投影系數(shù)對應(yīng)的分解信號的能量。初始化時,Edi=0,li=0;

(2) 對于動態(tài)分析時間段集合I={I1,I2,…}中的每一個元素Ii,查找出與其相鄰的所有下一個動態(tài)分析時間段集合{Ij}, 如果

Edi+Eei>Edj

(4)

則有

(5)

連接算法∏可保證在整個分析時間段內(nèi)基元函數(shù)組合所對應(yīng)的信號與分析信號最為相似。基元函數(shù)在動態(tài)分析時間支持區(qū)Ii內(nèi)的瞬時頻率fIi(t)=aμ+2bμt,ti∈Ii,將所有動態(tài)時間段集合I={I1,I2,…}中所對應(yīng)的頻率曲線集合fI={fI1,fI2,…}按時間先后順序連接成線則為信號在整個分析時間段內(nèi)的瞬時頻率估計。

2 S變換

S變換是Stockwell等在研究地球物理數(shù)據(jù)時提出的一種非平穩(wěn)信號分析方法。S變換可以認為是一種加高斯窗、且窗寬與信號頻率成反比的特殊STFT。S變換在信號低頻段具有較高的頻率分辨率,而在高頻段具有較高的時間分辨率,因此,S變換具有多分辨特性,克服了STFT分辨率固定的缺點。并且,S變換也可認為是一種經(jīng)過相位修正的連續(xù)小波變換(Continuous wavelet transform,CWT),即S變換保持了原始信號的絕對相位信息,彌補了CWT缺乏相位信息的不足。

對于給定信號,其S變換與S逆變換分別為

(6)

(7)

S變換是一種線性變化,滿足線性疊加原理,相對于二次型時頻分布(如Wigner-Ville分布,Choi-Williams分布等),不存在交叉項的干擾,即對于信號x(t)=x1(t)+x2(t),其S變換滿足

S[x(t)]=S[x1(t)]+S[x2(t)]

(8)

為了便于S變換的離散化計算,將其寫成關(guān)于頻域信號X(f)的表達式,即

(9)

式中:X(f)為x(t)的傅里葉變換。由此可得到S變換的離散化計算公式為

(10)

式中:T為采樣間隔;N為采樣點數(shù),j,n=0,1,…,N-1。

離散S逆變換為

(11)

式中:k=0,1,…,N-1。

3 變轉(zhuǎn)速齒輪故障特征提取步驟

齒輪出現(xiàn)局部故障時,其振動信號中會出現(xiàn)調(diào)幅調(diào)頻信號,對該調(diào)幅調(diào)頻信號的提取與分析對齒輪的故障診斷具有重要意義。而當齒輪處于變轉(zhuǎn)速下運行時,該調(diào)幅調(diào)頻信號具有與轉(zhuǎn)速相關(guān)的時變特點,此時,傳統(tǒng)的信號濾波方法不再適用。

CPP算法采用分段擬合的思想對信號中頻率呈曲線變化的信號成分進行頻率估計,可準確估計出變轉(zhuǎn)速下非平穩(wěn)信號的頻率曲線;而S變換是一種線性時頻分析方法,不存在交叉項,且保留了信號的絕對相位信息。因此,針對頻率大范圍變化的調(diào)頻調(diào)幅信號的提取,本文將CPP方法與S變換相結(jié)合,提出了基于CPP與S變換的自適應(yīng)時頻濾波方法,并將其用于變轉(zhuǎn)速下的齒輪故障特征提取。本文方法具體計算步驟如下:

H(t,f)={H(ti,f),i=1,2,…,N}

(12)

式中:fωi為ti時刻的半帶寬。

(3) 然后采用自適應(yīng)時頻濾波器H(t,f)對信號的時頻分布S(τ,f)進行時頻濾波,得到濾波后的時頻分布S′(τ,f)=S(τ,f)×H(t,f);

(4) 將時頻濾波結(jié)果S′(τ,f)進行S逆變換,即可得到提純后的齒輪故障振動信號S′;

4 算法仿真

為驗證本文方法的有效性,設(shè)置一調(diào)幅調(diào)頻信號sig1(如式(13))來模擬變轉(zhuǎn)速下的齒輪局部故障,模擬齒數(shù)為26,齒輪嚙合頻率被1倍轉(zhuǎn)頻調(diào)制,齒輪局部故障信號如圖1。為模擬強噪聲的干擾,在信號中加入強度為-2 dB的高斯白噪聲,得到的仿真合成信號如圖2。圖2中,齒輪局部故障信號已被完全淹沒。

sig1=(1+cos(2π×(25t+sin(3πt))))×cos(2π×26×(25t+sin(3πt)))

(13)

S=sig1+SN

(14)

圖1 仿真變轉(zhuǎn)速齒輪故障信號Fig.1 Simulation signal of a faulted gear under variable rotational speed

圖2 仿真合成信號Fig.2 Simulation composite signal

圖3 估計的嚙合頻率與實際嚙合頻率對比Fig.3 Comparison of estimated and actual gear mesh frequency

圖4 仿真信號的自適應(yīng)時頻濾波器Fig.4 Adaptive time-frequency filter of simulated signal

采用圖4中的自適應(yīng)時頻濾波器對信號進行濾波,濾取的齒輪故障振動信號如圖5(a)所示。圖5(b)為圖5(a)在時間段0.19~0.22 s的局部放大圖,圖中,濾取的信號與原始齒輪故障信號基本重合,僅在幅值上存在部分差異。

圖5 自適應(yīng)時頻濾波提取的齒輪故障信號及局部放大信號Fig.5 Extracted signal of a faulted gear and its partly enlarged signal using adaptive time-frequency filter

同時采用EEMD方法對圖2所示的合成信號進行分析,得到的第2個IMF如圖6(a)。圖6(b)為圖6(a)在時間段0.19~0.22 s的局部放大圖,圖中,IMF與原始齒輪故障信號不僅存在幅值上的差異,而且存在相位畸變。

圖6 EEMD方法提取的齒輪故障信號及局部放大信號Fig.6 Extracted signal of a faulted gear and its partly enlarged signal using EEMD

對圖5(a)和圖6(a)所示信號分別進行階次譜分析,得到的階次譜如圖7所示,其中,圖7(a)為自適應(yīng)時頻濾波信號的階次譜,圖7(b)為IMF2的階次譜。對比圖7(a)和7(b),雖然兩圖中在階次25.04、26.09和27.04處均出現(xiàn)了峰值,表明齒輪嚙合頻率旁邊出現(xiàn)了1倍轉(zhuǎn)頻調(diào)制現(xiàn)象,但圖7(a)中的噪聲明顯要低,且幅值特性也明顯要好。

圖7 階次譜對比Fig.7 Comparison of order spectrum

5 應(yīng)用實例

為檢驗本文方法在實測齒輪箱振動信號中提取變轉(zhuǎn)速下齒輪故障特征的有效性,在齒輪箱上進行試驗,并采用斷齒齒輪和正常齒輪進行對比分析。試驗齒輪為正齒輪,主動軸與從動軸的齒數(shù)比為55∶75,斷齒故障設(shè)置為主動軸上,整體切割一個齒以模擬齒輪斷齒故障,如圖8所示。

圖8 斷齒齒輪Fig.8 A broken gear

為減少傳遞路徑的影響,將PCB振動加速度傳感器直接置于主動軸的軸承座上。試驗時,齒輪箱處于空載狀態(tài)。試驗采用LMS數(shù)據(jù)采集儀同步采集振動加速度信號和主動軸轉(zhuǎn)速信號,采樣頻率為8 192 Hz。圖9為變轉(zhuǎn)速下拾取的齒輪斷齒故障的振動信號,信號中出現(xiàn)了較多的沖擊成分。

圖9 變轉(zhuǎn)速下齒輪斷齒故障的齒輪箱振動信號Fig.9 Vibration signal of a gearbox with broken gear under variable rotational speed

采用CPP方法對圖9所示信號進行分析,估計出的嚙合頻率曲線如圖10中的虛線所示,圖中,頻率逐漸增加,表明齒輪箱處于升速階段。圖10中的實線為實測嚙合頻率曲線(即光電式轉(zhuǎn)速傳感器測取的主動軸的轉(zhuǎn)頻與齒數(shù)的乘積)。圖10中,兩條曲線基本重合。

圖10 估計的嚙合頻率與實測嚙合頻率對比Fig.10 Comparison of estimated and measured gear mesh frequency

根據(jù)圖10中估計的嚙合頻率曲線(虛線)設(shè)計自適應(yīng)時頻濾波器,得到的自適應(yīng)時頻濾波器的時頻分布如圖11,黑色表示阻帶,白色表示通帶。圖11中可看出,濾波器的中心頻率隨著信號頻率的變化而改變。

圖11 斷齒齒輪的自適應(yīng)時頻濾波器Fig.11 Adaptive time-frequency filter of a broken gear

采用圖11所示的自適應(yīng)時頻濾波器對信號進行濾波分析,濾取的齒輪故障振動信號如圖12。圖12中存在一定的調(diào)制現(xiàn)象,但無法確定故障類型。進一步對濾取的信號(如圖12)進行階次分析,得到的階次譜如圖13。圖13中,在階次53、55.01、57.01處出現(xiàn)明顯的峰值,表示齒輪的嚙合頻率被2倍所轉(zhuǎn)頻調(diào)制,與齒輪斷齒的轉(zhuǎn)頻調(diào)制現(xiàn)象相符,驗證了本文方法提取變轉(zhuǎn)速下齒輪故障特征有效性。

Fig.12 自適應(yīng)時頻濾波提取的齒輪斷齒故障信號Fig.12 Extracted signal of a broken gear using adaptive time-frequency filter

圖13 斷齒齒輪自適應(yīng)時頻濾波振動信號的階次譜Fig.13 Order spectrum of filtered signal of a broken gear using adaptive time-frequency filter

為增加對比,采用EEMD方法對圖9所示信號進行分析,得到的第1個IMF如圖14。對IMF1進行階次分析,得到的階次譜如圖15。圖15中,在階次53、55.01、57.01處也出現(xiàn)了峰值,表明出現(xiàn)了齒輪局部故障,但與圖13對比可知,圖15中在高頻和低頻段均存在了較多的噪聲,其效果要遜色于圖13。

圖14 第1個IMFFig.14 IMF1

圖15 第1個IMF的階次譜Fig.15 Order spectrum of IMF1

同時,采用正常齒輪來進行試驗,以便進行對比分析。圖16為正常齒輪在變轉(zhuǎn)速下采集到的齒輪箱振動信號,圖中,信號的幅值較圖9要小,且沖擊成分較少。

圖16 變轉(zhuǎn)速下正常齒輪箱振動信號Fig.16 Vibration signal of a normal gearbox under variable rotational speed

采用本文方法對圖16所示正常齒輪箱振動信號進行分析,自適應(yīng)時頻濾波后的信號如圖17。對圖17所示信號進行階次譜分析,得到的階次譜如圖18。圖18中,在嚙合頻率階次55處有一個明顯的峰值,但該峰值相對于齒輪斷齒故障的幅值要小很多,且嚙合頻率附近的調(diào)制邊頻帶并不明顯,即不存在轉(zhuǎn)頻調(diào)制現(xiàn)象,與正常齒輪的振動特性相符。

圖17 正常齒輪箱自適應(yīng)時頻濾波得到的濾波信號Fig.17 Filtered signal of a normal gearbox using adaptive time-frequency filter

圖18 正常齒輪箱濾波信號的階次譜Fig.18 Order spectrum of filtered signal of a normal gearbox

6 結(jié) 論

針對變轉(zhuǎn)速下的齒輪故障信號分離與故障特征提取,提出了基于CPP與S變換的自適應(yīng)時頻濾波方法,算法仿真和應(yīng)用實例驗證了本文方法的有效性和優(yōu)越性。本文主要結(jié)論如下:

(1) 自適應(yīng)時頻濾波方法可根據(jù)信號的頻率變化特點自適應(yīng)地改變中心頻率和帶寬,能有效濾取頻率呈曲線變化的調(diào)幅調(diào)頻信號,具有較好的信號自適應(yīng)性,非常適合于變轉(zhuǎn)速下非平穩(wěn)信號的分析。

(2) 與EEMD方法進行了對比分析,結(jié)果表明,自適應(yīng)時頻濾波方法濾取的信號無相位畸變,實現(xiàn)了零相位濾波。

(3) 由于通帶內(nèi)噪聲和頻域截斷的影響,自適應(yīng)時頻濾波方法濾取的信號與原始信號存在一定的幅值誤差,因此,如何減少上述兩因素對幅值的影響有待進一步改進。

參 考 文 獻

[1] 丁康,李巍華,朱小勇.齒輪及齒輪箱故障診斷實用技術(shù)[M].北京:機械工業(yè)出版社,2005.

[2] 王平,廖明夫.滾動軸承故障診斷的自適應(yīng)共振解調(diào)技術(shù)[J].航空動力學報,2005,20(4): 606-612.

WANG Ping, LIAO Mingfu.Adaptive demodulated resonance technique for the rolling bearing fault diagnosis[J].Journal of Aerospace Power, 2005, 20(4): 606-612.

[3] JENA D P, SAHOO, S, PANIGRAHI S N.Gear fault diagnosis using active noise cancellation and adaptive wavelet transform[J].Measurement, 2014, 47(1): 356-372.

[4] YAN R Q, GAO R X, CHEN X F.Wavelet for fault diagnosis of rotary machines: A review with applications[J].Signal Processing, 2014, 96(5): 1-15.

[5] CHENG Junsheng, YU Dejie, YANG Yu.The application of energy operator demodulation approach based on EMD in machinery fault diagnosis[J].Mechanical Systems and Signal Processing, 2007, 21(2): 668-677.

[6] LEI Yaguo, LIN Jing, HE Zhengjia, et al.A review on empirical mode decomposition in fault diagnosis of rotating machinery[J].Mechanical System and Signal Processing, 2013, 35(1/2): 108-126.

[7] CHENG Junsheng, YANG Yu.A rotating machinery fault diagnosis method based on local mean decomposition[J].Digital Signal Processing, 2012, 22(2): 356-366.

[8] 王衍學,何正嘉,訾艷陽,等.基于LMD的時頻分析方法及其機械故障診斷應(yīng)用研究[J].振動與沖擊,2012,31(9): 9-12.

WANG Yanxue, HE Zhengjia, ZI Yanyang, et al.Several key issues of local mean decomposition method used in mechanical fault diagnosis[J].Journal of Vibration and Shock, 2012, 31(9): 9-12.

[9] 彭富強,于德介,劉堅.一種基于多尺度線調(diào)頻基的稀疏信號分解方法[J].振動工程學報,2010,23(3): 333-338.

PENG Fuqiang, YU Dejie, LIU Jian.Sparse signal decomposition method based on multi-scale chirplet[J].Journal of Vibration and Engineering, 2010, 23(3): 333-338.

[10] 騰偉,安宏文,馬志勇,等.基于時頻濾波的汽輪機半速渦動故障成分提取[J].振動與沖擊,2015,34(3): 178-182.

TENG Wei, AN Hongwen, MA Zhiyong, et al.Semi-speed oil whirl fault component extraction in a stream turbine based on time-frequency filtering[J].Journal of Vibration and Shock, 2015, 34(3): 178-182.

[11] STOCKWELL R G, MANSINHA L, LOWE R P.Localization of the complex spectrum: the S transform[J].IEEE Transactions on Signal Processing, 1996, 44(4): 998-1001.

[12] 郭遠晶,魏燕定,周曉軍,等.S變換時頻譜SVD降噪的沖擊特征提取方法[J].振動工程學報,2014,27(4): 621-628.

GUO Yuanjing, WEI Yanding, ZHOU Xiaojun, et al.Impact feature extracting method based on S transform time-frequency spectrum denoised by SVD[J].Journal of Vibration and Engineering, 2014, 27(4): 621-628.

[13] 張云強,張培林,李兵.分數(shù)階時頻局部二值模式譜在齒輪故障診斷中的應(yīng)用[J].振動與沖擊,2015,34(11): 122-127.

ZHANG Yunqiang, ZHANG Peilin, LI Bing.Application of fractional order time-frequency local binary pattern spectra in gear fault diagnosis[J].Journal of Vibration and Shock, 2015,34(11): 122-127.

[14] 劉建敏,劉遠宏,江鵬程,等.基于包絡(luò)S變換時頻圖像提取齒輪故障特征[J].振動與沖擊,2014,33(1): 165-169.

LIU Jianmin, LIU Yuanhong, JIANG Pengcheng, et al.Extraction of gear fault features based on the envelope and time-frequency image of S transformation[J].Journal of Vibration and Shock, 2014, 33(1): 165-169.

[15] WU Chunyan, LIU Jian, PENG Fuqiang, et al.Gearbox fault diagnosis using adaptive zero phase time-varying filter based on multi-scale chirplet sparse signal decomposition[J].Chinese Journal of Mechanical Engineering, 2013, 26(4): 831-838.

[16] CANDES E J, CHARLTON P R, HELGASON H.Detecting highly oscillatory signals by chirplet path pursuit[J].Applied and Computational Harmonic Analysis, 2008, 24(1): 14-40.

[17] LUO Jiesi, YU Dejie, LIANG Ming.Gear fault detection under time-varying rotating speed via joint application of multiscale chirplet path pursuit and multiscale morphology analysis[J].Structural Health Monitoring, 2012, 11(5): 526-537.

[18] PENG Fuqiang, YU Dejie, LUO Jiesi.Sparse signal decomposition method based on multi-scale chirplet and its application to the fault diagnosis of gearboxes[J].Mechanical System and Signal Processing, 2011, 25(2): 549-557.

猜你喜歡
振動故障信號
振動的思考
科學大眾(2023年17期)2023-10-26 07:39:14
信號
鴨綠江(2021年35期)2021-04-19 12:24:18
完形填空二則
振動與頻率
天天愛科學(2020年6期)2020-09-10 07:22:44
故障一點通
基于FPGA的多功能信號發(fā)生器的設(shè)計
電子制作(2018年11期)2018-08-04 03:25:42
中立型Emden-Fowler微分方程的振動性
奔馳R320車ABS、ESP故障燈異常點亮
基于LabVIEW的力加載信號采集與PID控制
故障一點通
主站蜘蛛池模板: 色欲不卡无码一区二区| 成色7777精品在线| 久久免费看片| 亚洲无码精品在线播放| 欧美成人综合视频| 日本高清免费一本在线观看 | 综合社区亚洲熟妇p| 国产精品夜夜嗨视频免费视频| 中文字幕66页| 玖玖精品视频在线观看| 色欲国产一区二区日韩欧美| 三上悠亚精品二区在线观看| 精品福利视频网| 一级毛片在线播放免费| 视频二区亚洲精品| 精品1区2区3区| 高清不卡一区二区三区香蕉| 欧美一级高清片久久99| 亚洲无码高清免费视频亚洲| 亚洲精品无码AⅤ片青青在线观看| 亚洲91精品视频| 亚洲男人的天堂久久香蕉网| 美女视频黄频a免费高清不卡| 欧美一级色视频| 一级做a爰片久久毛片毛片| 制服丝袜国产精品| 91丝袜美腿高跟国产极品老师| 伊大人香蕉久久网欧美| 亚洲 欧美 偷自乱 图片| 国产福利免费观看| 色婷婷成人网| 伊人色综合久久天天| 国产正在播放| 成人免费黄色小视频| 亚洲色图欧美在线| 亚洲国产看片基地久久1024| 呦系列视频一区二区三区| 久久香蕉国产线看观看式| 国产麻豆福利av在线播放| 另类重口100页在线播放| 欧美国产日产一区二区| 亚洲AV无码精品无码久久蜜桃| 99草精品视频| 四虎成人精品| 无码AV动漫| 国产精品人成在线播放| 亚洲av中文无码乱人伦在线r| 91视频首页| 国产精品13页| 香蕉综合在线视频91| 四虎影视永久在线精品| 国产精品无码影视久久久久久久 | 久久久久人妻一区精品| 久久精品中文字幕少妇| 久久精品一品道久久精品| 粗大猛烈进出高潮视频无码| 国产导航在线| 青青草欧美| 免费人欧美成又黄又爽的视频| 无码aaa视频| 国产chinese男男gay视频网| 国产在线视频欧美亚综合| 久久黄色视频影| 国产污视频在线观看| 无码福利日韩神码福利片| 亚洲一区二区日韩欧美gif| 欧美成人在线免费| 久久人人妻人人爽人人卡片av| 男女精品视频| 免费高清毛片| 国产原创演绎剧情有字幕的| 亚洲美女AV免费一区| 亚洲三级成人| 九九精品在线观看| 中文字幕无码电影| 亚洲日韩欧美在线观看| 亚洲男人的天堂网| 国产AV无码专区亚洲精品网站| 欧美日本在线一区二区三区| 国产欧美日韩一区二区视频在线| 中国精品久久| 一级毛片免费不卡在线|