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

卷積型小波變換實現及機械早期故障診斷應用

2013-02-13 06:35:24田福慶李克玉丁慶喜
振動與沖擊 2013年7期
關鍵詞:故障診斷故障信號

羅 榮,田福慶,李克玉,丁慶喜

(海軍工程大學 兵器工程系,武漢 430033)

機械故障診斷領域中通常采用內積型小波變換,并利用Mallat算法實現其快速運算。但Mallat算法用于機械故障診斷存有缺陷:① 該算法引入下采樣環節使分解長度隨分解層數的增加依次減半,從而使概貌信號信息量減少,時間分辨率降低,不利于機械故障特征提取[1]。若對Mallat算法分解后的頻帶進行單支重構,不僅會增加計算量,還會使信號產生畸變[2-3]。②任何實際濾波器組的幅頻特性與理想幅頻特性有一定差距使各濾波器頻帶有一定重疊及Mallat算法特有的上下抽樣環節,共同造成Mallat算法分解后的子帶存在頻率折疊,并非原信號的真實分量,會影響分析結果的可讀性,頻率折疊嚴重時還會造成誤判[4-6]。

為消除Mallat算法固有缺陷對機械故障特征提取的不利影響,更好應用小波變換進行機械故障診斷,文獻[4-5]提出單子帶重構改進算法,用FFT濾波原理強制型去掉多余頻率成分;文獻[6]提出的巴特沃斯小波變換算法基于性能優良的巴特沃斯濾波器組實現無下抽樣環節的小波變換快速算法。本文提出摒棄內積型小波變換,采用卷積型小波變換進行機械故障診斷。基于此,本文首先推導出卷積型小波變換快速分解算法,并給出基于濾波器組的遞歸分解實現方法。與Mallat算法相比,該快速算法不存在隔點取樣過程,只對上一尺度分解結果進行移位,無論信號被分解多少層,每層分解所得逼近信號與細節信號長度始終與原始信號一致。因此該快速算法進行故障特征提取無須單子帶重構過程,計算量少,便于工程實現。更重要的是,由于不存在隔點取樣過程,該快速算法分解所得各子帶信號為信號的真實分量,不存在頻率折疊,有利于提取故障特征。可見,卷積型小波變換完全克服了內積型小波變換不利于機械故障診斷問題。文獻[7]通過詳地比較認為卷積型小波變換對信號的處理效果優于內積型小波變換,不僅頻帶分離效果較內積型小波變換優越,且所提故障特征信息也較內積型小波變換明顯。因此卷積型小波變換更適合機械故障特征提取。

針對滾動軸承早期故障診斷與定量識別難題與共振解調法、沖擊脈沖法的不足,本文提出將卷積型小波變換與共振解調法、沖擊脈沖法相結合的新方法對滾動軸承早期故障進行診斷與定量識別,并給出具體實現過程。仿真實驗與實例分析表明,卷積型小波變換能消除Mallat算法對機械故障特征提取的不利影響,用于機械故障診斷較內積型小波變換優勢明顯。本文所提新方法可有效實現對滾動軸承早期故障診斷與定量識別,具有一定應用價值。

1 卷積型小波變換及快速分解算法

設{Vj}j∈Z為空間L2(R)的正交多分辨率分析,φ(x)∈V0為{Vj}j∈Z的正交尺度函數,φ(x)∈W0為{Vj}j∈Z的正交小波函數,Wj為Vj在Vj-1中的正交補。由于Vj?Vj-1,Wj?Wj-1,所以Vj中元素 2-j/2φ(2-jx)與Wj中元素2-j/2ψ(2-jx)可分別用Vj-1的基表示:

設f(t)∈L2(R),記ψ(t)以s為尺度的擴張函數為:ψs(t)=1/sψ(t/s),則f(t)與ψs(t)的卷積型連續小波變換定義為[8]:

令s=2j,x=k2j得:

可證明Wjf(k)與內積型小波系數僅相差一個常數,故Wjf(k)即f(t)在Wj中的系數。

記φ(t)以s為尺度擴張函數為:φs(t)=1/sφ(t/s),且Ssf(x)=f(t)*φs(t),再令s=2j,x=k2j,則有:

同樣可證明Sjf(k)與內積型小波系數也僅相差一個常數,故Sjf(k)即f(t)在Vj中的系數。

1.1 卷積型小波變換快速分解算法推導

對Sjf(x)=f(t)*φj(t)等號兩邊對x求傅氏變換并將式(3)代入得:

將上式轉換到時域得:

令x=k2j得:

同理對細節系數有:

式(10)、(11)即為卷積型小波變換的快速分解算法,可采用基于濾波器組的遞歸分解方法實現。實現過程為:

對式(10)等號兩邊同時進行z變換得:

同理:

式中:Wjf(k)為Sj-1f(k)通過濾波器1/G(z2j-1)后的輸出,G(z2j-1)由濾波器G(z)每兩點間插入2j-1-1個零所得。卷積型小波變換快速分解算法基于濾波器組的遞歸分解實現方法見圖1。

2 仿真實驗

通過兩個仿真實驗比較兩類小波變換對信號的處理效果。

圖1 卷積型小波變換快速分解算法基于濾波器組的遞歸分解實現方法Fig.1 The recursion realization of the fast decomposition algorithm of convolution type of wavelet transformation based on filter banks

2.1 兩類小波變換對多頻信號分解效果比較

設時間信號y(t)的表達式為:

即信號y(t)包含15 Hz,30 Hz,60 Hz,80 Hz,100 Hz,120 Hz,150 Hz等七個頻率成分。以512 Hz采樣率對此信號采樣512點,分別利用卷積型小波變換快速分解算法與Mallat算法對采樣信號進行3層分解,所用小波均為db40小波,分解頻譜比較如圖2所示。信號y(t)的小波分解得到的各子帶理論上的頻率成分如表1所示。為便于比較,圖2內積型小波分解頻譜為單支重構頻譜。由圖2看出,卷積型小波變換所得各子帶無任何虛假頻率,為信號的真實分量,只因小波濾波器的非理想頻域截止特性造成子帶含相鄰子帶部分頻率成分。此因任何實際濾波器均不具有理想頻域截止特性所致。而Mallat算法分解所得子帶信號,不僅含相鄰子帶部分頻率成分,且存在頻率折疊,出現明顯的虛假頻率成分,如W1子帶 136 Hz,W2子帶 68 Hz,S3子帶34 Hz等。虛假頻率混在子帶信號中可致子帶信號發生畸變。如圖3中內積型小波分解所得近似子帶S3提取的二倍頻信號出現明顯畸變,遠不及卷積型小波變換提取的清晰,此因內積型小波分解所得近似子帶S3中混有頻率34 Hz,幅值0.36的虛假成分所致。此仿真實驗說明Mallat算法用于信號分解會產生嚴重頻率折疊,所得子帶信號并非信號的真實分量,不利于機械故障特征提取。而卷積型小波變換用于信號分解則不會產生任何頻率折疊,完全克服了Mallat算法產生頻率折疊缺陷。

表1 各子帶理論頻率成分Tab.1 Frequencys of bands under ideal conditions

圖2 兩類小波變換所得子帶頻譜比較Fig.2 The comparison between bands spectrums of two types of wavelet transformation

圖3 兩類小波變換所得近似子帶S3比較Fig.3 The comparison between approximate bands of two types of wavelet transformation

2.2 兩類小波變換對多載波調幅信號各調幅波及調制信息提取效果比較

由于軸承故障振動信號大多由復雜的具有多載波與多調制的調幅波,從中提取各調幅波并進行包絡解調即可提取故障特征信息進行故障診斷。設多載波調幅信號為:

對式(15)信號以512 Hz采樣率采樣512點,如圖4所示。① 分別利用卷積型小波變換快速分解算法與Mallat算法對該信號進行2層分解(用db40小波),分解結果及頻譜如圖5、圖6所示。② 分別對兩類小波變換的分解結果進行希爾伯特解調,所得包絡及包絡譜如圖7、圖8所示。比較圖5、圖6可知,卷積型小波變換分解所得各子帶未產生虛假頻率,較完美地提取出各調幅波;而Mallat算法分解所得子帶則有明顯虛假頻率,以致其對各調幅波的提取效果不及卷積型小波變換。比較圖7、圖8可知,卷積型小波變換分解結果包絡基本與各調幅波真實包絡一致,而內積型小波變換分解結果包絡較粗糙并與真實包絡相差較大,甚至包絡譜中出現虛假譜峰。因此,卷積型小波變換對多載波調幅信號各調幅波及調制信息的提取效果好于內積型小波變換。此仿真實驗說明Mallat算法產生的頻率折疊會影響內積型小波變換對各調幅波及調制信息的準確提取,不利于機械特征提取。

由以上兩仿真實驗可知,卷積型小波變換分解所得子帶是信號的真實分量,無頻率折疊;內積型小波變換所得子帶則存在嚴重頻率折疊。頻率折疊嚴重時必會影響故障特征提取,導致機械故障診斷困難,甚至會造成誤診。故卷積型小波變換較內積型小波變換更適合機械故障診斷。

圖4 多載波調幅信號及其頻譜Fig.4 The amplitude modulation signal and its spectrum

圖5 卷積型小波變換后各子帶信號及頻譜Fig.5 The results of convolution wavelet transformation and their spectrums

圖6 內積型小波變換后各子帶及頻譜Fig.6 The results of traditional wavelet transformation and their spectrums

圖7 卷積型小波變換后各子帶包絡及包絡譜Fig.7 The envelopes of results of convolution wavelet transformation and their spectrums

圖8 內積型小波變換后各子帶包絡及包絡譜Fig.8 The envelopes of results of traditional wavelet transformation and their spectrums

3 滾動軸承早期故障診斷與定量識別新方法提出

滾動軸承運轉中會因各種原因導致其過早損壞,如裝配不當、潤滑不良、水分及異物侵入、腐蝕等。即使安裝、潤滑及使用維護正常,經過一段時間運轉也會出現疲勞剝落及磨損而不能正常工作。因此,對滾動軸承早期故障進行診斷、準確判定其損傷的嚴重程度,預防重大事故發生,成為故障診斷領域中重要研究方向[9]。

共振解調法可將軸承缺陷信息從復雜的調幅振動中分離出來,是目前滾動軸承故障診斷中常用方法之一。共振解調技術關鍵在于帶通濾波的實現。但傳統的共振解調技術中,通常采用單一帶通濾波器濾除混雜在信息中的無用成分以提高信噪比。主要問題有:① 傳統帶通濾波器參數選擇取決于操作者經驗及歷史數據,獲得滿意結果前需多次嘗試,費時費力[10];②滾動軸承發生故障時,通常不止一個載波頻率,用單一帶通濾波器進行包絡分析,難全面提取隱含在振動信號各調制頻帶內的故障信息特征,甚至包絡譜中出現非故障特征頻率的譜峰。因此,傳統的共振解調技術對滾動軸承中晚期故障診斷雖有一定效果,但不適合滾動軸承早期故障診斷[9]。

沖擊脈沖法是另一種簡單適用的滾動軸承故障檢測方法由瑞典SPM公司提出,被公認為對診斷滾動軸承局部損傷故障工程實用性最強的方法[11]。滾動軸承各元件存在缺陷時,在運轉過程中損傷點會反復快速撞擊與之接觸的其它元件表面而發生沖擊,引起脈沖性振動。沖擊脈沖的強弱代表故障的嚴重程度,因此沖擊脈沖的大小可用于對故障程度進行定量識別。沖擊脈沖法即通過檢測沖擊點處的沖擊波幅值判定軸承工作狀態及定量識別故障程度[9,11]。然而,沖擊脈沖值不僅與軸承的油膜厚度、操作程度有關,且與軸承幾何尺寸及轉速有關。為得到衡量各種滾動軸承狀態標準,SPM公司規定了只與軸承工作狀況有關的標準分貝值,即表示沖擊值的增加率,并給出標準分貝值的經驗計算公式[9,12-13]:

式中:dBn為標準分貝值;N為軸轉速(r/min);D為軸承內徑(m);SV為沖擊值(m/s2)。軸承運轉過程中最大標準分貝值dBn,max即可定量識別故障程度。實踐證明,對穩態運行設備,評定各類軸承的工作狀態的dBn,max為:

(1)0 dB≤dBn,max≤21 dB 正常狀態,軸承工作狀態良好;

(2)21 dB≤dBn,max≤35 dB 輕微故障,軸承早期損傷;

(3)35 dB≤dBn,max≤60 dB 嚴重故障,軸承已有明顯損傷。

沖擊脈沖法可有效判定軸承工作狀態及準確定量識別故障,卻不能確定故障具體位置。卷積型小波變換與內積型小波變換相同,為較精細的小波分解方法,具有良好時頻局域化功能及優良帶通濾波功能,無須重構過程,便于工程實現,分解的各子帶信號不存在頻率折疊,具有較內積型小波變換更好的頻帶分離特性。因此,為彌補共振解調法與沖擊脈沖法不足,實現故障的準確定位與定量識別,本文將卷積型小波變換引入共振解調法與沖擊脈沖法中,提出將卷積型小波變換與共振解調法、沖擊脈沖法相結合的新方法對滾動軸承早期故障進行診斷與定量識別,具體過程如下:

(1)將原始振動信號按卷積型小波變換快速分解算法進行分解;

(2)將分解所得各子帶進行Hilbert包絡解調得各子帶包絡,對包絡進行快速FFT變換得各個子帶包絡譜,即可從中提取滾動軸承早期故障特征,從而定位故障到具體元件;

(3)將各子帶包絡譜幅值作為沖擊值代入式(16)計算標準分貝值,提取各子帶故障特征頻率處對應的標準分貝值;

(4)選取各子帶中同一故障特征頻率對應標準分貝值的最大值dBn,max,對滾動軸承故障進行定量識別。

4 實例分析

據實驗臺實驗數據,實驗時用電火花制造支撐軸承外圈點蝕損傷模擬早期微弱損傷故障,為支撐軸承座加墊片模擬平行不對中故障。支撐軸承為6205深溝球軸承,幾何參數為:軸承節徑39.04 mm,滾動體直徑7.94 mm,滾動體9個,接觸角0°。實驗時軸轉頻約29 Hz。據上述參數計算出外圈故障特征頻率為103 Hz。實驗數據采樣頻率8 kHz,數據長度為8 192點,數據單位為10 m/s2。

4.1 兩類小波變換對外圈故障特征提取效果比較及故障定量識別

圖9為所測故障振動信號及頻譜。對此故障振動信號分別進行4層卷積型小波分解與內積型小波分解,所用小波均為db40。卷積型小波分解所得細節子帶W1,W2,W3的包絡譜與內積型小波分解所得細節子帶W1,W2,W3的包絡譜比較見圖10。由圖10知,卷積型小波分解所得細節子帶W1,W2,W3的包絡譜在外圈故障特征頻率103 Hz及2,3,4倍頻處均出現明顯譜線,而內積型小波分解所得細節子帶W1,W2,W3的包絡譜僅在外圈故障特征頻率103 Hz及2倍頻處出現明顯譜線。說明卷積型小波提取的外圈故障特征好于內積型小波。卷積型小波分解所得細節子帶W1,W2,W3的包絡譜分貝值見圖11。由圖11知,卷積型小波分解各子帶外圈故障特征頻率及倍頻對應的分貝值最大為25.61 dB。說明滾動軸承外圈存在早期損傷故障。因此本文提出的卷積小波包絡解調方法可有效診斷并定量識別滾動軸承早期故障。

圖9 振動信號及頻譜Fig.9 The vibration signals and its spectrum

4.2 兩類小波變換對不對中故障特征提取效果比較

對此故障信號分別進行6層卷積型小波分解與內積型小波分解,所用小波均為db40。卷積型小波與內積型小波分解所得第6層近似子帶S6及頻譜見圖12、圖13。由圖12知,卷積型小波分解所得子帶S6為基頻分量29.3 Hz與二次諧波分量57.62 Hz的疊加。說明經卷積型小波“去粗存精,層層剝離”分解,所得近似信號即“骨架”較清楚地提取到不對中故障的兩倍頻振動現象。由圖13可知,內積型小波提取的兩倍頻振動現象明顯不及卷積型小波提取到的清晰。因為內積型小波分解所得子帶S6中存在虛假頻率成分67.38 Hz。此虛假頻率經分析發現由Mallat算法中隔點插零環節由頻率折疊所致。67.38 Hz即為二次諧波分量57.62Hz關于S6子帶理論上最大頻率62 Hz的映像,該映像經重構低通濾波器后未完全濾掉而留在S6子帶中。由此可見,小波濾波非理想特性及Mallat算法上下抽樣環節共同造成了Mallat算法中的頻率折疊;而卷積型小波無上下抽樣環節,不會產生頻率折疊,且分解后的頻帶無須重構。故卷積型小波變換無上下抽樣環節使其在機械故障診斷中較內積型小波變換更具優勢。

圖10 兩類小波變換后各子帶包絡譜比較Fig.10 The comparison between bands envelope spectrums of two types of wavelet transformation

圖11 卷積型小波變換后各子帶包絡譜分貝值圖Fig.11 The dB values of bands envelope spectrums of convolution wavelet transformation

圖12 卷積型小波變換后子帶S6及其頻譜Fig.12 The band S6 of convolution wavelet transformation and its spectrum

圖13 內積型小波變換后子帶S6及其頻譜Fig.13 The band S6 of traditional wavelet transformation and its spectrum

5 結論

(1)與內積型小波變換相比,卷積型小波變換無須重構過程,便于工程實現,分解的各子帶不存在頻率折疊,具有較內積型小波變換更好的頻帶分離特性,可消除Mallat算法存在頻率折疊等固有缺陷對機械故障診斷的不利影響。

(2)采用卷積型小波變換進行機械故障診斷,推導出卷積型小波變換快速分解算法,給出基于濾波器組的遞歸分解實現方法。

(3)提出的將卷積型小波變換與共振解調法、沖擊脈沖法相結合的新方法及具體實現過程,彌補了滾動軸承早期故障診斷與定量識別難題及共振解調法、沖擊脈沖法的不足。

(4)卷積型小波變換能消除Mallat算法存在頻率折疊等固有缺陷對機械故障診斷的不利影響,用于機械故障診斷較內積型小波變換更有優勢。本文所提新方法可有效實現對滾動軸承早期故障診斷與定量識別,具有一定工程應用價值。

[1]萬書亭,呂路勇,何玉靈.基于提升模式非抽樣小波變換的滾動軸承故障診斷方法研究[J].振動與沖擊,2009,28(1):170-173.

WAN Shu-ting,Lü Lu-yong,HE Yu-ling.Fault diagnosis method of rolling bearing based on undecimated wavelet transformation of lifting scheme[J].Journal of Vibration and Shock,2009,28(1):170-173.

[2] Zhao X Z,Ye B Y.Convolution wavelet packet transform and its applications to signal processing[J].Digital Signal Processing,2010,20(5):1352-1364.

[3]張兆偉,陳維義,林春生.一種對稱小波重構算法及其在突變信號檢測中的應用[J].武漢理工大學學報(交通科學與工程版),2005,29(6):685-688.

ZHANG Zhao-wei,CHEN Wei-yi,LIN Chun-sheng.Symmetric wavelets based non-distortion decomposition and reconstruction algorithm[J].Journal of Wuhan University of Technology(Transportation Science &Engineering),2005,29(6):685-688.

[4] Yang J G,Park S T.An anti-aliasing algorithm for discrete wavelet transform [J].Mechanical Systems and Signal Processing,2003,17(5):945-954.

[5]楊建國.小波分析及其工程應用[M].北京:機械工業出版社,2005:47-55.

[6]何嶺松,王峻峰.巴特沃斯小波變換算法在故障診斷中的應用[J].華中理工大學學報,2000,28(10):31-33.

HE Ling-song WANG Jun-feng.Butterworth based on discrete wavelettransform and itsapplication in machine fault diagnosis[J].J.Huazhong Univ.of Sci.& Tech,2000,28(10):31-33.

[7]趙學智,陳統堅,葉邦彥,等.內積型和卷積型小波變換對信號處理效果的研究[J].機械工程學報,2004,40(3):55-60.

ZHAO Xue-zhi,CHEN Tong-jian,YE Bang-yan,et al.Research on signal processing effect of both inner product type of wavelet transformation and convolution type of wavelet transformation[J].Chinese Journal of Mechanical Engineering,2004,40(3):55-60.

[8]秦前清,楊宗凱.實用小波分析[M].西安:西安電子科技大學出版社,1995.

[9]何正嘉,訾艷陽,張西寧.現代信號處理及工程應用[M].西安:西安交通大學出版社,2007.

[10]蘇文勝,王奉濤,張志新,等.EMD降噪和譜峭度法在滾動軸承早期故障診斷中的應用[J].振動與沖擊,2010,3(29):18-21.

SU Wen-sheng,WANG Feng-tao,ZHANG Zhi-xin,et al.Application of EMD denoising and spectral kurtosis in early fault diagnosis of rolling element bearings[J].Journal of Vibration and Shock,2010,3(29):18-21.

[11]黃 海,王曉萍,張衛東.沖擊脈沖法評價滾動軸承故障的系數自修正方法[J].機械工程學報,1999,35(2):85-87.

HUANG hai,WANG Xiao-ping,ZHANG Wei-dong.Parameter self-modification method for SPM using in rolling bearing fault diagnosis[J].Chinese Journal of Mechanical Engineering,1999,35(2):85-87.

[12]姜洪開,王仲生.第二代小波包構造及發動機微弱損傷識別[J].北京航空航天大學學報,2007,33(7):777-780.

JIANG Hong-kai,WANG Zhong-sheng.Second generation wavelet packet construction and aircraft engine weak damage identification[J].Journal of Beijing University of Aeronautics and Astronautics,2007,33(7):777-780.

[13]徐玉秀,原陪新,王占國.沖擊脈沖法在電機軸承故障診斷中的應用[J].沈陽工業大學學報,2000,22(6):453-455.

XU Yu-xiu,YUAN Pei-xin,WANG Zhan-guo.Application of shock impulse method on motor bearing fault diagnosis[J].Journal of Shenyang University of Technology,2000,22(6):453-455.

猜你喜歡
故障診斷故障信號
信號
鴨綠江(2021年35期)2021-04-19 12:24:18
完形填空二則
故障一點通
基于FPGA的多功能信號發生器的設計
電子制作(2018年11期)2018-08-04 03:25:42
奔馳R320車ABS、ESP故障燈異常點亮
基于LabVIEW的力加載信號采集與PID控制
因果圖定性分析法及其在故障診斷中的應用
故障一點通
江淮車故障3例
基于LCD和排列熵的滾動軸承故障診斷
主站蜘蛛池模板: 久久网综合| 黄色网在线| 天天干天天色综合网| 天堂亚洲网| 国产精品密蕾丝视频| 5388国产亚洲欧美在线观看| 日韩a在线观看免费观看| 国产又爽又黄无遮挡免费观看 | 一区二区三区国产精品视频| 亚洲性一区| 亚洲手机在线| 欧美色伊人| 国产成人久久综合777777麻豆 | 亚洲第七页| 国产99在线观看| 亚洲成人网在线观看| 一级毛片在线免费视频| 毛片免费高清免费| 久草视频一区| 99热在线只有精品| 亚洲一区色| 日韩亚洲高清一区二区| 国产成人精品优优av| 永久成人无码激情视频免费| 九九热免费在线视频| 久久久91人妻无码精品蜜桃HD| 亚洲国产AV无码综合原创| 国产精品亚洲五月天高清| 青青国产成人免费精品视频| 在线一级毛片| 日本免费新一区视频| 亚洲天堂久久久| 国产色偷丝袜婷婷无码麻豆制服| 色综合久久88| 日韩国产亚洲一区二区在线观看| 国产91久久久久久| 欧美 国产 人人视频| 日本精品视频一区二区| 久久精品国产在热久久2019| 亚洲男人天堂久久| 制服丝袜亚洲| 日韩麻豆小视频| 国产剧情伊人| 国产黄色片在线看| 久久这里只精品国产99热8| 97精品国产高清久久久久蜜芽| 国产精品手机视频| 亚洲天堂日韩在线| 六月婷婷精品视频在线观看| 国国产a国产片免费麻豆| 国产极品美女在线播放| 亚洲区欧美区| 一本色道久久88| 欧美一级特黄aaaaaa在线看片| 亚洲香蕉伊综合在人在线| 国产精品无码影视久久久久久久| 国产全黄a一级毛片| 又猛又黄又爽无遮挡的视频网站| 成人字幕网视频在线观看| 久久这里只有精品免费| 成人福利一区二区视频在线| 午夜影院a级片| 伦精品一区二区三区视频| 青青操视频免费观看| 日韩第一页在线| 国产视频自拍一区| 久久久亚洲国产美女国产盗摄| 男人天堂伊人网| 福利视频一区| 欧美日韩综合网| 91麻豆久久久| 女同国产精品一区二区| 干中文字幕| 国产福利观看| 狂欢视频在线观看不卡| 国产美女精品人人做人人爽| 欧美成人h精品网站| v天堂中文在线| 国产综合精品一区二区| 一区二区三区四区精品视频 | 色九九视频| 欧美中文字幕第一页线路一|