亚洲免费av电影一区二区三区,日韩爱爱视频,51精品视频一区二区三区,91视频爱爱,日韩欧美在线播放视频,中文字幕少妇AV,亚洲电影中文字幕,久久久久亚洲av成人网址,久久综合视频网站,国产在线不卡免费播放

        ?

        一種自適應(yīng)Morlet小波包絡(luò)解調(diào)的弱故障檢測(cè)方法

        2015-10-24 03:49:54侯新國(guó)楊忠林
        船電技術(shù) 2015年10期
        關(guān)鍵詞:特征故障信號(hào)

        牛 超,侯新國(guó),楊忠林

        (海軍工程大學(xué)電氣工程學(xué)院,武漢 430033)

        一種自適應(yīng)Morlet小波包絡(luò)解調(diào)的弱故障檢測(cè)方法

        牛 超,侯新國(guó),楊忠林

        (海軍工程大學(xué)電氣工程學(xué)院,武漢 430033)

        為自適應(yīng)實(shí)現(xiàn)Morlet小波與故障沖擊特征成分的最優(yōu)匹配,采用基于Shannon小波熵的方法優(yōu)化帶寬參數(shù)設(shè)計(jì)了最優(yōu)Morlet小波。針對(duì)最佳尺度的求取難題,利用譜峭度與小波熵均能敏感反映沖擊性的特性,提出了基于峭熵比的最佳尺度的求取方法。基于此,提出基于自適應(yīng)Morlet小波包絡(luò)解調(diào)的弱故障特征提取方法。仿真實(shí)驗(yàn)與實(shí)例分析表明:該方法克服了傳統(tǒng)包絡(luò)解調(diào)需要人為設(shè)定帶通濾波器參數(shù)的不足,減少人工干預(yù),能自動(dòng)有效地從強(qiáng)噪背景中提取微弱故障特征,具有一定的工程應(yīng)用價(jià)值。

        Morlet小波 自適應(yīng)包絡(luò)解調(diào) shannon小波熵 譜峭度 峭熵比

        0 引言

        包絡(luò)解調(diào)法是電機(jī)軸承故障診斷中最常使用的方法之一[1]。在對(duì)故障振動(dòng)信號(hào)進(jìn)行包絡(luò)解調(diào)前,通常要對(duì)故障振動(dòng)信號(hào)進(jìn)行帶通濾波。然而,帶通濾波的實(shí)現(xiàn)需要人為選擇濾波中心頻率與帶寬,應(yīng)用極為不便[2]。將小波應(yīng)用于振動(dòng)信號(hào)的包絡(luò)解調(diào)分析是一種常用方法[3]。Morlet小波是在故障診斷領(lǐng)域內(nèi)常用的一種小波函數(shù),這不僅由于Morlet小波是左右兩邊都按指數(shù)衰減的余弦信號(hào),可以與故障振動(dòng)信號(hào)實(shí)現(xiàn)較好的匹配,而且Morlet小波具有非常好的帶通濾波特性[4-5]。文獻(xiàn)[6]提出基于Morlet小波的包絡(luò)檢波算法提取淹沒(méi)在噪聲中的故障沖擊成分和沖擊成分所在頻帶,取得一定效果,但該方法沒(méi)有考慮Morlet小波的優(yōu)化,不能自適應(yīng)實(shí)現(xiàn)對(duì)沖擊成分的最優(yōu)匹配,影響了包絡(luò)解調(diào)效果。利用最小Shannon小波熵自適應(yīng)優(yōu)化帶寬參數(shù)設(shè)計(jì)與沖擊特征成分相匹配的Morlet小波,將含有被調(diào)制信號(hào)的信息最豐富的所在頻段從原信號(hào)頻帶中分離出來(lái),是優(yōu)化小波變換算法的一種新思路,但最佳變換尺度的求取是基于Morlet小波分析的難點(diǎn)。文獻(xiàn)[4][7]采用小波變換系數(shù)矩陣的奇異值比譜(SingularValue Ratio,sVR)求取最佳尺度,該方法在利用SVR探測(cè)信號(hào)的周期分量時(shí),存在截?cái)嗾`差導(dǎo)致其經(jīng)常失效。由于譜峭度和Shannon小波熵均能反映強(qiáng)噪聲背景下沖擊信號(hào)的特征,本文提出基于譜峭度與Shannon小波熵比值即峭熵比(Kurtosis Entropy Ratio, KER)的方法求取最佳變換尺度,以彌補(bǔ)最佳尺度求取方法的缺陷。在此基礎(chǔ)上,提出基于自適應(yīng)Morlet小波包絡(luò)解調(diào)的弱故障特征提取方法,以有效提取強(qiáng)噪聲背景下的微弱故障特征。該方法首先依據(jù)信號(hào)的特點(diǎn)自適應(yīng)地設(shè)計(jì)最優(yōu)Morlet小波,并利用最優(yōu)Morlet連續(xù)小波變換對(duì)信號(hào)進(jìn)行全頻段帶通濾波;然后依據(jù)峭熵比自適應(yīng)地獲取最優(yōu)頻帶;最后對(duì)最優(yōu)頻帶小波系數(shù)取實(shí)部,并對(duì)其求模,即可獲得最優(yōu)頻帶的包絡(luò)實(shí)現(xiàn)包絡(luò)解調(diào)。仿真實(shí)驗(yàn)分析表明:該方法克服了傳統(tǒng)包絡(luò)解調(diào)需要人為設(shè)定帶通濾波器參數(shù)的不足,能有效地從強(qiáng)噪背景中提取弱故障特征。

        1 Morlet小波及其包絡(luò)解調(diào)原理分析

        Morlet小波是高斯包絡(luò)下的復(fù)指數(shù)函數(shù),其通常的定義式為:

        其頻譜為:

        式中fc是小波的中心頻率;fb是帶寬參數(shù),決定了波形振蕩衰減的快慢程度和母小波的帶寬,fb越小,波形衰減得越快,帶寬越寬,直至逼近一個(gè)脈沖信號(hào);fb越大,則波形衰減得越慢,帶寬越窄,最后退化成一個(gè)余弦信號(hào)。因此,在應(yīng)用Morlet小波進(jìn)行故障診斷時(shí),需要優(yōu)化fb,以實(shí)現(xiàn)小波與沖擊特征成分的最優(yōu)匹配。

        由式(2),經(jīng)傅里葉逆變換知小波ψ(t)的實(shí)部ψr(t)和虛部ψi(t)分別是t的實(shí)偶函數(shù)與實(shí)奇函數(shù),因此它們的頻譜ψr(f)和ψi(f)分別是實(shí)偶函數(shù)和純虛奇函數(shù);由于Morlet小波頻譜函數(shù)僅含有正頻率部分,故Morlet小波是解析小波,于是可推知ψr(f)和ψi(f)滿(mǎn)足下列關(guān)系:

        式中sgn(f)為符號(hào)函數(shù)。上式是一個(gè)時(shí)間函數(shù)Hilbert變換式的等價(jià)頻域表示,即式(3)表明ψi(t)是ψr(t)的Hilbert變換。

        由式(5)知,Wxr(a, b)本質(zhì)是信號(hào)x( t)通過(guò)一個(gè)頻率響應(yīng)為Ψr(aω )濾波器輸出。

        其中, H[·]表示Hilbert變換運(yùn)算,因此有:

        故信號(hào)x( t)的Morlet小波變換系數(shù)的實(shí)虛部也是一對(duì)Hilbert變換,所以Morlet小波變換系數(shù)Wx( a, b)的求模即可提取到Wxr(a, b)的包絡(luò),從而實(shí)現(xiàn)對(duì)信號(hào)的包絡(luò)解調(diào)。

        從以上分析可知,Morlet小波通過(guò)實(shí)虛部?jī)纱螏V波構(gòu)建了信號(hào)的解析信號(hào),實(shí)現(xiàn)將帶通濾波與包絡(luò)解調(diào)兩個(gè)過(guò)程合二為一。

        2 基于Shannon小波熵自適應(yīng)設(shè)計(jì)最優(yōu)Morlet小波的方法

        在突變信號(hào)的檢測(cè)中,為了突出特征成分,可根據(jù)信號(hào)特征,自適應(yīng)地選擇最優(yōu)帶寬參數(shù)fb以使Morlet基小波與特征成分具有最大的相似性。Shannon小波熵[5,8]是一種評(píng)價(jià)相似的標(biāo)準(zhǔn),其定義為:

        Shannon小波熵值的大小反映了概率分布的均勻性,最不確定的概率分布(等概率分布)具有最大的熵值。當(dāng)Shannon小波熵最小時(shí),對(duì)應(yīng)的Morlet基小波就是與特征成分最相似的小波。根據(jù)這一原理,可自適應(yīng)設(shè)計(jì)出最優(yōu)Morlet小波。

        以仿真信號(hào)為例驗(yàn)證基于Shannon小波熵自適應(yīng)設(shè)計(jì)最優(yōu)Morlet小波有效性,設(shè)信號(hào)x( t)為:

        以1 Hz的頻率對(duì)x( t)采樣1024點(diǎn),其時(shí)域波形如圖1所示。fb的取值范圍為[1, 100],步長(zhǎng)為0.5,經(jīng)計(jì)算可得到fb與Shannon小波熵之間的關(guān)系如圖2所示。由圖2知,fb=15時(shí),具有最小的小波熵,對(duì)應(yīng)的基小波就是與特征成分最匹配的Morlet小波,其小波變換時(shí)頻圖如圖3(b)所示。作為對(duì)比,分別令fb取5和90,再作Morlet連續(xù)小波變換,時(shí)頻圖分別如圖3(a)和3(c)所示。

        圖1 仿真信號(hào)時(shí)域波形

        圖2 fb與Shannon小波熵之間的關(guān)系曲線(xiàn)

        圖3 fb取不同值時(shí)的Morlet小波變換時(shí)頻圖

        從圖3可知,當(dāng)fb取5時(shí),信號(hào)的時(shí)間分辨率稍好于fb取15時(shí)的分辨率,但其頻率分辨率比f(wàn)b取15時(shí)的低,信號(hào)的高頻分量與低頻分量沒(méi)有分開(kāi)。當(dāng)fb取95時(shí),信號(hào)的頻率分辨率稍好于fb取15時(shí)的分辨率,但其時(shí)間分辨率比f(wàn)b取15時(shí)的低,以致信號(hào)的兩個(gè)高頻分量沒(méi)有被分開(kāi)。fb取15時(shí),信號(hào)的時(shí)頻聚集性最好,3個(gè)分量被清晰的區(qū)分開(kāi)來(lái)。綜合考慮時(shí)頻分辨率,fb取15時(shí)為最優(yōu),因此,基于Shannon小波熵自適應(yīng)設(shè)計(jì)最優(yōu)Morlet小波的方法是有效的。

        3 基于峭熵比的最佳尺度的求取方法

        尺度為a的Morlet小波ψa(t)的中心頻率為fc/a。在對(duì)信號(hào)作小波變換時(shí),將各尺度的小波的中心頻率fc/a預(yù)先歸一化成一組固定值,信號(hào)的采樣率與固定值的乘積即為相應(yīng)尺度的小波系數(shù)的中心頻率。故尺度a不僅控制濾波的頻帶范圍,且濾波效果主要由帶寬參數(shù)fb和尺度a決定。因此,設(shè)計(jì)出最優(yōu)Morlet小波后,還須求取最佳尺度a,以將含有被調(diào)制信號(hào)信息最豐富的那個(gè)頻段從強(qiáng)噪聲背景信號(hào)中分離出來(lái)。

        譜峭度(Spectral Kurtosis:SK)的本質(zhì)為每根譜線(xiàn)處的峭度值,其計(jì)算公式為[9]:

        由前述,Shannon小波熵值的大小反映了小波系數(shù)分布的均勻性。故障沖擊越強(qiáng)的頻帶,即是含有被調(diào)制信號(hào)的信息最豐富的頻帶,其信號(hào)幅值分布越不均勻,其譜峭度也越大,Shannon小波熵值越小,故Shannon小波熵值同譜峭度一樣,也能反映信號(hào)中的沖擊及其所在的頻帶。但實(shí)際信號(hào)分析表明,譜峭度的最大值與Shannon小波熵的最小值不一定在同一頻帶取到,也就是說(shuō)基于譜峭度求取的最佳尺度與基于小波熵求取的最佳尺度不一定是同一尺度。鑒于此,綜合考慮譜峭度與Shannon小波熵,本文提出一個(gè)新的沖擊性評(píng)價(jià)指標(biāo)——峭熵比,即譜峭度與Shannon小波熵的比值來(lái)衡量各個(gè)頻帶的沖擊性。顯然,沖擊性越強(qiáng)的頻帶,其譜峭度越大,其Shannon小波熵值越小,因此其峭熵比也就越大?;诖?,本文提出基于峭熵比的最佳尺度的求取方法,其基本步驟如下:1)先對(duì)信號(hào)進(jìn)行最優(yōu)Morlet小波變換;2)再用小波系數(shù)代替式(12)中的H(t, f )計(jì)算各個(gè)尺度的譜峭度;3)再將各個(gè)尺度的小波系數(shù)的實(shí)部處理成一個(gè)概率分布序列并計(jì)算其Shannon小波熵,第j個(gè)尺度的Shannon小波熵計(jì)算公式為:

        再依次計(jì)算各個(gè)尺度峭熵比,當(dāng)峭熵比取得最大值時(shí)的尺度即為最佳尺度。

        4 自適應(yīng)Morlet小波包絡(luò)解調(diào)的弱故障特征提取方法

        為了有效地提取強(qiáng)噪聲背景下的微弱故障特征和克服傳統(tǒng)包絡(luò)解調(diào)需要人為確定帶通濾波器參數(shù)的不足,在前文研究基礎(chǔ)上提出基于自適應(yīng)Morlet小波包絡(luò)解調(diào)的弱故障特征提取方法。該方法首先根據(jù)信號(hào)的特征成分采用Shannon小波熵方法優(yōu)化帶寬參數(shù),自適應(yīng)地設(shè)計(jì)出最優(yōu)Morlet小波,并對(duì)信號(hào)進(jìn)行最優(yōu)Morlet連續(xù)小波變換,實(shí)現(xiàn)在全頻帶范圍內(nèi)對(duì)信號(hào)進(jìn)行帶通濾波。然后,依據(jù)峭熵比自適應(yīng)地求取最佳尺度并提取最佳尺度的小波系數(shù),從而實(shí)現(xiàn)對(duì)信號(hào)的最優(yōu)帶通濾波。最后,對(duì)最佳尺度的小波系數(shù)取模即可實(shí)現(xiàn)對(duì)最優(yōu)頻帶的包絡(luò)解調(diào),提取到最優(yōu)頻帶的包絡(luò),再對(duì)包絡(luò)求頻譜,即可得到包絡(luò)譜,從而進(jìn)一步提取到明顯的微弱故障特征。

        5 仿真實(shí)例

        以仿真信號(hào)驗(yàn)證基于自適應(yīng)Morlet小波包絡(luò)解調(diào)的故障診斷方法的有效性。用類(lèi)似于沖擊信號(hào)的10個(gè)高斯原函數(shù)組成的信號(hào)疊加強(qiáng)噪聲模擬實(shí)際故障振動(dòng)信號(hào)(信噪比為-12.4dB),由下式給出:

        式(14)中,ti=200, 400 …,2000,n(t)為強(qiáng)噪聲。用fs=2048 Hz的采樣率對(duì)x(t)進(jìn)行采樣,采樣點(diǎn)數(shù)為2048點(diǎn),采樣信號(hào)及其頻譜如圖4所示。由圖4時(shí)域圖可知強(qiáng)背景噪聲完全掩蓋了沖擊信號(hào),因此對(duì)其故障沖擊特征的提取是具有相當(dāng)難度的。首先依據(jù)最小Shannon小波熵原理自適應(yīng)設(shè)計(jì)最優(yōu)Morlet小波。設(shè)母小波中心頻率fc=0.5,各小波歸一化的中心頻率范圍[0.005, 0.5],步長(zhǎng)為0.005,fb的范圍[1, 100],步長(zhǎng)為0.5。依次計(jì)算各個(gè)fb對(duì)應(yīng)的Shannon小波熵得到fb與Shannon小波熵之間的關(guān)系如圖5所示。由圖5可知,fb=17時(shí),具有最小的小波熵,其對(duì)應(yīng)的基小波就是與特征成分最匹配的最優(yōu)Morlet小波。設(shè)計(jì)出了最優(yōu)Morlet小波,然后基于峭熵比的方法求取最佳尺度。各尺度的峭熵比與尺度的關(guān)系如圖6所示。由圖6可知,峭熵比在中心頻率為256 Hz的頻帶取得最大值,因此基于峭熵比求取方法求得的最優(yōu)變換尺度為fs?fc/256=4。此尺度的中心頻率即是各尺度的中心頻率離振動(dòng)固有頻率250 Hz最近的頻率,所以其中必然包含了最為豐富的故障沖擊調(diào)制信息。

        圖4 含強(qiáng)噪聲的模擬故障信號(hào)

        圖5 fb與Shannon小波熵的關(guān)系

        設(shè)計(jì)出最優(yōu)Morlet小波和求取最優(yōu)變換尺度后,下面對(duì)原始信號(hào)進(jìn)行最優(yōu)Morlet連續(xù)小波變換并提取最佳尺度的小波系數(shù)。尺度為4的小波系數(shù)的實(shí)部及其頻譜如圖7所示。從圖7可以看出,最優(yōu)Morlet小波有非常好的帶通濾波性能,將包含故障調(diào)制信息最多的頻帶從包含強(qiáng)噪聲的原始信號(hào)頻譜中完美分離出來(lái),充分的降低了最優(yōu)頻帶中的噪聲以致其小波系數(shù)中清晰地出現(xiàn)了間隔約為0.1s的沖擊,這是非常有利于提取微弱故障特征的。對(duì)最佳尺度的小波系數(shù)取模即可得小波系數(shù)的包絡(luò),再對(duì)包絡(luò)求傅里葉變換即可得相應(yīng)的包絡(luò)譜。尺度為4的小波系數(shù)包絡(luò)和包絡(luò)譜如圖8所示。從圖8中可知,尺度為4的小波系數(shù)的包絡(luò)譜在故障沖擊頻率10 Hz及其倍頻處出現(xiàn)很明顯的譜線(xiàn),這說(shuō)明自適應(yīng)Morlet小波包絡(luò)解調(diào)已相當(dāng)成功的提取到了隱藏在強(qiáng)噪聲背景下的微弱故障沖擊特征。

        圖6 各尺度的峭熵比與尺度的關(guān)系

        圖7 最佳尺度的小波實(shí)部系數(shù)及其頻譜

        圖8 最佳尺度的小波系數(shù)實(shí)部包絡(luò)及包絡(luò)譜

        從以上分析過(guò)程可知:1)最優(yōu)Morlet小波有非常好的帶通濾波性能和降噪能力,這非常有利于弱故障特征的提取;2)自適應(yīng)Morlet小波包絡(luò)解調(diào)能根據(jù)信號(hào)特點(diǎn)自適應(yīng)設(shè)計(jì)最優(yōu)Morlet小波并在全頻帶內(nèi)自適應(yīng)選擇最優(yōu)頻帶進(jìn)行包絡(luò)解調(diào),無(wú)須人工設(shè)定濾波中心頻率與帶寬,便于用計(jì)算機(jī)進(jìn)行故障自動(dòng)識(shí)別。

        以上仿真表明,自適應(yīng)Morlet小波包絡(luò)解調(diào)的弱故障提取方法能有效的隱藏在強(qiáng)噪聲背景下的微弱故障特征。

        6 結(jié)論

        綜合本文可得如下結(jié)論:

        1)根據(jù)信號(hào)自身的沖擊特征,采用基于Shannon小波熵的方法優(yōu)化了帶寬參數(shù)設(shè)計(jì)了最優(yōu)Morlet小波,自適應(yīng)地實(shí)現(xiàn)其與沖擊特征成分的最優(yōu)匹配。

        2)利用譜峭度與小波熵均能敏感反映沖擊性的特性,提出基于峭熵比的新方法來(lái)求取最優(yōu)變換尺度。該方法綜合考慮了譜峭度與小波熵兩種指標(biāo),能方便有效的求出最優(yōu)變換尺度,具有一定的工程應(yīng)用價(jià)值。

        3)提出基于自適應(yīng)Morlet小波包絡(luò)解調(diào)的弱故障特征提取方法。仿真實(shí)驗(yàn)分析表明:該方法無(wú)須人工設(shè)定濾波中心頻率與帶寬,能自動(dòng)有效地從強(qiáng)噪背景中提取微弱故障特征,具有一定的工程應(yīng)用價(jià)值。

        [1] Y Yang, D J Yu, Js Cheng. A fault diagnosis approach for roller bearing based on IMF envelopespectrum andsVM[J]. Measurement, 2007, 40(9): 943-950.

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

        [3] Yuh-Taysheen. On thestudy of applying Morlet wavelet to the Hilbert transform for the envelope detection of bearingVibrations[J]. Mechanicalsystems andsignal Processing, 2009, 23(5): 1518- 1527.

        [4] 程發(fā)斌, 湯寶平, 鐘佑明. 基于最優(yōu)Morlet小波和SVD的濾波消噪方法及故障診斷的應(yīng)用[J]. 振動(dòng)與沖擊, 2008, 27(2) : 91-94.

        [5] JIN L,LIANGSHENG Q U.Feature extraction based on morlet wavelet and its application for mechanical fault diagnosis[J].sound andVibration, 2000, 234(1): 135-148.

        [6] 何嶺松, 李巍華, 用Morlet小波進(jìn)行包絡(luò)檢波分析[J],振動(dòng)工程學(xué)報(bào),2002,15(1): 119-122.

        [7] 蔣永華, 湯寶平, 董紹江. 自適應(yīng)Morlet小波降噪方法及在軸承故障特征提取中的應(yīng)用[J]. 儀器儀表學(xué)報(bào), 2010, 31(12) : 2712-2717.

        [8] NIKOLAOU N G,ATONIADIS I A.Demodulation ofVibrationsignals generated by defects in rolling element bearings using complexshifted morlet wavelets[J]. Mechanicalsystems andsignal Processing, 2002,16(4): 677- 694.

        [9]sawalhi N, Randall R B, Endo H. The enhancement of fault detection and diagnosis in rolling element bearings using minimum entropy deconvolution combined withspectral kurtosis [ J]. Mechanicalsystems andsignal Processing, 2007, 21(6): 2616- 2633.

        The Method to Extract Weak Fault Feature Based on Adaptive Morlet Wavelet Envelope Demodulation

        Niu Chao, Hou Xinguo, Yang Zhonglin
        (College of Electrical Engineering, Naval University of Engineering, Wuhan 430033, China)

        In order to achieve adaptive optimal match with the impact component,shannon wavelet entropy is used to optimize bandwidth parameter of the Morlet wavelet to design the optimal Morlet wavelet. Aimed at the problem of obtaining the optimalscale, the new approach which is based on kurtosis entropy ratio to acquire the optimalscale is presented. In addition, the new method to extract weak fault feature based on adaptive morlet wavelet envelope demodulation is proposed in the paper. Thesimulation results and analysis on the incipient fault data of ball bearingshow that the method overcomes the defect that the parameters of band-pass filter areselected by experience of the user in conventional envelope demodulation, decreases the dependence on user, and can automatic and effectively extract the weak fault feature.so, the method hassome applicationValue.

        Morlet wavelet; adaptive envelope demodulation;shannon wavelet entropy; kurtosis entropy ratio;spectral kurtosis

        TM712

        A

        1003-4862(2015)10-0026-05

        2015-03-14

        牛超(1982-),碩士研究生。研究方向:信號(hào)檢測(cè)與處理,設(shè)備的故障診斷等。

        猜你喜歡
        特征故障信號(hào)
        信號(hào)
        鴨綠江(2021年35期)2021-04-19 12:24:18
        完形填空二則
        故障一點(diǎn)通
        如何表達(dá)“特征”
        不忠誠(chéng)的四個(gè)特征
        基于FPGA的多功能信號(hào)發(fā)生器的設(shè)計(jì)
        電子制作(2018年11期)2018-08-04 03:25:42
        抓住特征巧觀(guān)察
        奔馳R320車(chē)ABS、ESP故障燈異常點(diǎn)亮
        基于LabVIEW的力加載信號(hào)采集與PID控制
        故障一點(diǎn)通
        亚洲免费av电影一区二区三区| 波多野42部无码喷潮在线| 毛片大全真人在线| 亚洲天堂资源网| av在线男人的免费天堂| 日韩精品免费一区二区三区观看| 亚洲国产精品一区二区www| 人妻在卧室被老板疯狂进入国产 | 欧美大黑帍在线播放| 二区久久国产乱子伦免费精品 | 漂亮丰满人妻被中出中文字幕| 少妇粉嫩小泬喷水视频www| 青草福利在线| 亚洲素人日韩av中文字幕| 丝袜美腿福利一区二区| 国产探花在线精品一区二区| 国产午夜精品理论片| 一本之道加勒比在线观看| 99re66在线观看精品免费| 成av免费大片黄在线观看| 亚洲精品成人av观看| 久久国产精品国语对白| 欧美日韩精品乱国产| 在线播放亚洲第一字幕| 3亚洲日韩在线精品区| 加勒比东京热一区二区| 又大又紧又粉嫩18p少妇| 国产亚洲欧美在线观看的| 中国av一区二区三区四区| 黑人巨大精品欧美| 精品国产一区二区三区av 性色| 日韩中文字幕精品免费一区| 亚洲成人激情深爱影院在线 | 日本中文字幕婷婷在线| 99久久久无码国产精品6| 日韩一区二区超清视频| 亚洲av天堂一区二区| 亚洲综合成人婷婷五月网址| 玩弄放荡人妻一区二区三区| 久久亚洲精品国产精品婷婷| 天堂视频在线观看一二区|