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

        ?

        基于提升小波的心電信號(hào)P、T波檢測(cè)快速算法

        2013-07-07 03:04:32司玉娟郎六琪程延偉李賀佳臧國(guó)華
        關(guān)鍵詞:信號(hào)檢測(cè)

        姚 成,司玉娟,郎六琪,程延偉,李賀佳,臧國(guó)華

        (1.吉林大學(xué) 通信工程學(xué)院,長(zhǎng)春130022;2.長(zhǎng)春工程技術(shù)學(xué)院,長(zhǎng)春130117;3.吉林大學(xué) 珠海學(xué)院,廣東珠海519041)

        目前,隨著人們對(duì)心臟健康水平關(guān)注的不斷 加強(qiáng),心電信號(hào)(Electrocardiograph,ECG)的檢測(cè)分析技術(shù)也得到了迅猛的發(fā)展。但要完成對(duì)心電信號(hào)的分析診斷,首先要將心電信號(hào)的主要成分P波、P-R間期、QRS波群、ST段、T段、QT間期和U波等心臟病變?cè)\斷依據(jù)提取和識(shí)別出來(lái)。在這些波形及參數(shù)中,由于P波和T波頻率低,幅度小,通常都混有干擾和噪聲,因此檢測(cè)的難度較大。但檢測(cè)、定位P波和T波,直接影響對(duì)心電信號(hào)診斷的精度。

        當(dāng)前對(duì)于P波和T波檢測(cè)的研究成果層出不窮。較早發(fā)展起來(lái)的P、T波檢測(cè)方法有局部極值和局部變換法[1]、斜率域值法[2]、函數(shù)擬合法[3-4]等。這些方案的算法相對(duì)簡(jiǎn)單,計(jì)算量不大。但當(dāng)心電信號(hào)中混有復(fù)雜噪聲時(shí),算法的識(shí)別精度明顯下降,算法的抗噪效果不好。不適于在一些簡(jiǎn)單的心電分析儀中使用。對(duì)此,近年來(lái)出現(xiàn)基于小波變換的自適應(yīng)QRS-T對(duì)消P波檢測(cè)算法[5]、形態(tài)數(shù)學(xué)等的方法。算法主要是在P、T波的頻域內(nèi)對(duì)其進(jìn)行識(shí)別。在頻域中分析提高了算法整體抗噪聲的性能,但算法計(jì)算量大,需要高性能的計(jì)算機(jī)處理系統(tǒng)來(lái)實(shí)現(xiàn),難于在簡(jiǎn)單的處理器及硬件上實(shí)現(xiàn)。

        本文在研究現(xiàn)有檢測(cè)方法之后,針對(duì)目前計(jì)算機(jī)的性能和算法復(fù)雜度,提出了一種應(yīng)用提升小波去噪并結(jié)合差分法的復(fù)合檢測(cè)P、T波的方法。該方法在保證檢測(cè)精度的同時(shí),大大降低了算法的計(jì)算復(fù)雜度,使算法更容易在硬件上實(shí)現(xiàn)。

        1 提升小波變換及小波去噪原理

        1.1 提升小波原理

        小波分析帶有明顯的傅里葉變換對(duì)信號(hào)進(jìn)行頻域分析的能力,同時(shí)也拓展了其對(duì)信號(hào)進(jìn)行時(shí)域分析的能力。在小波分析算法的基礎(chǔ)上,1996年,Swedens提出了不依賴于傅里葉變化的小波提升算法(liftingscheme),其較經(jīng)典小波方法提供了一種快速實(shí)現(xiàn)方法,與經(jīng)典的Mallat算法相比,運(yùn)算量減少一半。同時(shí)能夠?qū)崿F(xiàn)小波變換的原位(in-place)計(jì)算,節(jié)省存儲(chǔ)單位,能夠?qū)崿F(xiàn)整數(shù)小波變換。其特性更適合于快速算法的構(gòu)造和硬件平臺(tái)的應(yīng)用。提升方案實(shí)現(xiàn)過(guò)程如下。

        正向小波變換的提升實(shí)現(xiàn)算法(預(yù)測(cè)步驟由奇序列預(yù)測(cè)偶序列開(kāi)始),其提升過(guò)程為

        步驟1 懶小波變換

        式中:x為信號(hào);l=0,1,…,N/2 -1。

        步驟2 提升與對(duì)偶提升

        式中:l=0,1,…,N/2 -1。

        步驟3 比例變換

        式中:s={s0,s1,…,sN/2-1};d={d0,d1,…,dN/2-1},最后得到的s和d分別為小波分解的低頻分量和高頻分量。

        算法 逆向小波變換的提升實(shí)現(xiàn)算法

        步驟1 比例變換

        步驟2 提升與對(duì)偶提升

        步驟3 逆懶小波變換

        算法分析表明,當(dāng)輸入數(shù)據(jù)量很大時(shí),提升算法比Mallat算法的計(jì)算量減小一半。

        1.2 小波閾值去噪原理

        小波閾值去噪的基本算法流程如下:

        ①分解過(guò)程:選定一種小波,對(duì)信號(hào)進(jìn)行N層小波分解;

        ②作用閾值過(guò)程:對(duì)分解得到的各層系數(shù)選擇一個(gè)閾值,并對(duì)細(xì)節(jié)系數(shù)作閾值處理;

        ③重建過(guò)程:降噪處理后的系數(shù)通過(guò)小波重建恢復(fù)原始信號(hào)。

        小波閾值收縮去噪法的理論是依據(jù)小波變換。特別是正交小波變換具有很強(qiáng)去數(shù)據(jù)相關(guān)性的特性,它能使信號(hào)的能量在小波域集中在一些大的小波系數(shù)中,而噪聲的能量卻分布于整個(gè)小波域內(nèi)。因此,經(jīng)小波分解后,信號(hào)的小波系數(shù)幅值要大于噪聲的系數(shù)幅值,可認(rèn)為幅值比較大的小波系數(shù)一般以信號(hào)為主,而幅值比較小的系數(shù)在很大程度上是噪聲。于是,采用設(shè)定閾值的辦法可把信號(hào)系數(shù)保留,而使大部分噪聲系數(shù)減小至零。

        2 結(jié)合去噪和P、T波識(shí)別的算法

        本文設(shè)計(jì)的復(fù)合檢測(cè)算法依據(jù)心電信號(hào)波形的特征和頻域特點(diǎn),首先利用選定的DB4小波進(jìn)行提升小波變換,對(duì)信號(hào)進(jìn)行多層分解;再對(duì)提升小波系數(shù)進(jìn)行加權(quán)閾值去噪。最后將去噪后的提升小波分解系數(shù)重構(gòu)至第四層,進(jìn)行差分運(yùn)算。在差分后得到的數(shù)據(jù)上在R波位置之前、后查詢模極值對(duì),確定P、T波的位置區(qū)間,在查詢模極值對(duì)過(guò)程中算法結(jié)合了相應(yīng)的可變跟隨閾值的方法,其目的是使算法更好地適應(yīng)心電信號(hào)個(gè)體差異和病變心電信號(hào)突然的變化。再將相應(yīng)的位置區(qū)間映射到去噪重構(gòu)信號(hào)上,進(jìn)行最后波峰點(diǎn)確定。

        2.1 P、T波小波分解層數(shù)的確定

        在心電信號(hào)中,T波在0~8Hz上聚集了近92%的能量,其波峰能量在0~8Hz之間;P波的頻譜帶寬在0~20Hz范圍內(nèi)聚集了90%的能量,其波峰的能量主要集中在3~12Hz。在本文實(shí)驗(yàn)中信號(hào)的采樣率為360Hz,經(jīng)過(guò)提升小波分解后,第四層的低頻小波系數(shù)頻帶范圍為0~22.5 Hz,第五層低頻小波系數(shù)頻帶范圍為0~11.25 Hz。因此,根據(jù)心電信號(hào)頻譜分析,確定在第四層來(lái)實(shí)現(xiàn)P、T波的檢測(cè)。

        2.2 提升小波加權(quán)閾值去噪

        在提升小波去噪過(guò)程中,根據(jù)心電信號(hào)的頻譜特征,選定DB4小波進(jìn)行8層分解。對(duì)各層高頻系數(shù)進(jìn)行閾值處理是去噪算法的核心過(guò)程。各層閾值的求取也是直接影響去噪效果的關(guān)鍵。本算法中將噪聲視為白噪聲,采用根據(jù)噪聲強(qiáng)度σ[6-7]和各層閾值加權(quán)系數(shù)k確定各層去噪閾值thrk。即:式中:n為需要進(jìn)行閾值處理信號(hào)長(zhǎng)度;σ為噪聲強(qiáng)度。且σ,式中:d(k)為提升小波分解得到的各尺度細(xì)節(jié)系數(shù)[8];k為每層小波去噪閾值的權(quán)重系數(shù)。算法根據(jù)實(shí)驗(yàn)結(jié)果,采取的權(quán)重系數(shù)為將第1層細(xì)節(jié)系數(shù)和第8層概貌系數(shù)置零。直接去掉心電信號(hào)中90~180Hz頻率成分和基線漂移。從第2層至第4層,其權(quán)重系數(shù)依次為0.8,0.6,0.4,第 5 層至第 8 層全部為 0.1。最大限度保留小波分解系數(shù)中心電信號(hào)有用成分,減小重構(gòu)信號(hào)的失真。

        2.3 算法流程

        Step1對(duì)采樣信號(hào)利用DB4小波進(jìn)行8層小波提升分解;

        Step2求取對(duì)各層細(xì)節(jié)系數(shù)處理的閾值(閾值采用加權(quán)處理,由2至8層,閾值的權(quán)重遞減);

        Step3對(duì)各層細(xì)節(jié)系數(shù)應(yīng)用相應(yīng)加權(quán)閾值進(jìn)行軟閾值去噪處理(其中對(duì)第1層高頻和第8層低頻系數(shù)直接置零);

        Step4由第8層開(kāi)始逐層重構(gòu)出去噪后信號(hào),其中將重構(gòu)得到第四層低頻系數(shù)CA4保存;

        Step6在差分結(jié)果R波頂點(diǎn)前后相應(yīng)時(shí)間窗口內(nèi),求識(shí)別模極大值對(duì)的正負(fù)閾值;

        Step7在差分結(jié)果R波頂點(diǎn)前后相應(yīng)時(shí)間窗口內(nèi),求識(shí)別模極值對(duì),并記錄極值位置;

        Step8在模極值序列中采用漏檢回溯策略進(jìn)行漏檢查詢;

        Step9在模極值序列中按照誤檢查詢策略進(jìn)行誤檢查詢;

        Step10根據(jù)模極值對(duì)的位置,在去噪重構(gòu)信號(hào)(Xdiff)中的相對(duì)模極值位置,進(jìn)行P、T波頂點(diǎn)定位;

        Step11在去噪信號(hào)中P、T波頂點(diǎn)前后一定時(shí)間窗口內(nèi)查找差分值零點(diǎn),即為起點(diǎn)和終點(diǎn)。

        2.4 可變跟隨閾值

        由于心電信號(hào)具有較大的個(gè)體差異,其信號(hào)的質(zhì)量也由于測(cè)試環(huán)境和手段變化而不同。同一個(gè)體的心電信號(hào)有時(shí)候也會(huì)有較大的變化,因此如果使用固定閾值,必然不能較好地適應(yīng)變化的測(cè)試群體,也不能適應(yīng)個(gè)體信號(hào)幅值的突變。因此本文提出了一種可變跟隨閾值的方法用于提高小波模極大值對(duì)的檢出率,使之適應(yīng)心電信號(hào)突變的檢測(cè)。

        具體方法為:在進(jìn)行檢測(cè)的過(guò)程中,當(dāng)檢測(cè)到滿足閾值的小波系數(shù)極值點(diǎn)后,為保證對(duì)心電信號(hào)突變具有良好的檢測(cè)效果,設(shè)置閾值跟隨當(dāng)前檢測(cè)點(diǎn)幅值,即新閾值更新為檢測(cè)極值點(diǎn)幅值與原閾值和的平均值0.6倍,如式(5)、式(6)

        式中:nmax為新檢測(cè)到的小波系數(shù)極大值點(diǎn)幅值。

        2.5 基于差分P、T波識(shí)別

        2.5.1 P 波檢測(cè)算法

        圖1 標(biāo)準(zhǔn)心電圖Fig.1 Standardelectrocardiogram

        ①搜索P波頂點(diǎn)。

        P波位于Q波時(shí)間軸左側(cè)。根據(jù)P波寬度與RR間期(心電周期)的寬度比例,確定其搜索范圍。步驟如下。

        步驟1確定P波頂點(diǎn)搜索范圍:由于PR間期范圍為0.12~0.2s,即相比 RR間期寬度的比例為0.12~0.2,則設(shè)定P波頂點(diǎn)的搜索范圍為:R波頂點(diǎn)坐標(biāo)減去PR段最小寬度至R波減去PR段的長(zhǎng)度。式中:Rmean為平均心電周期。

        步驟3根據(jù)模極值對(duì)相對(duì)位置判斷P波是正向還是倒置。

        步驟4在Xdiff中對(duì)應(yīng)的模極值對(duì)區(qū)間內(nèi),根據(jù)P波方向搜索幅度極大值或幅度極小值作為P波頂點(diǎn)。

        ②搜索P波起始點(diǎn)。

        搜索策略是從P波頂點(diǎn)位置開(kāi)始,向時(shí)間零點(diǎn)方向搜索一個(gè)P波寬度的區(qū)域,搜索區(qū)間內(nèi)查找第一個(gè)拐點(diǎn)作為P波的起始點(diǎn)。其方法還是在去噪信號(hào)差分值中尋找轉(zhuǎn)折點(diǎn),即連續(xù)兩個(gè)點(diǎn)基本平滑,第三個(gè)點(diǎn)的差分值變化較大。具體步驟如下。

        步驟2在Xdiff信號(hào)上指定搜索區(qū)間內(nèi)查找信號(hào)突變點(diǎn),即零點(diǎn)。

        ③搜索P波終點(diǎn)。

        搜索策略是從P波頂點(diǎn)位置開(kāi)始,向QRS波方向搜索一個(gè)P波寬度的區(qū)域,搜索區(qū)間內(nèi)查找第一個(gè)拐點(diǎn)作為P波的終點(diǎn)。具體過(guò)程與搜索P波起始點(diǎn)基本相同。

        2.5.2 T 波檢測(cè)算法

        T波檢測(cè)算法首先要找到T波的頂點(diǎn),之后在T波頂點(diǎn)的基礎(chǔ)上尋找T波的起始點(diǎn)和終止點(diǎn)。其搜索的策略與P波的檢測(cè)基本相同,不同之處是搜索區(qū)間不同。

        ①T波頂點(diǎn)搜索。

        T波搜索的方法與P波搜索方法基本相同。根據(jù)RR間期和T波信號(hào)時(shí)寬性質(zhì)設(shè)定查找T波頂點(diǎn)的搜索區(qū)間。在小波系數(shù)差分結(jié)果中相應(yīng)區(qū)間中查找模極值對(duì)。在模極值對(duì)指定的區(qū)間內(nèi),在中按照T波正向和倒置兩種情況確定T波頂點(diǎn)。其中搜索T波頂點(diǎn)范圍為0.32 ~0.44s。QRS 波的間期為 0.06~0.10s。因此設(shè)定 T波查詢范圍,起點(diǎn)為,終點(diǎn)為

        ②T波終點(diǎn)搜索。

        ③T波起始點(diǎn)搜索。

        查找到T波的頂點(diǎn)以后,在T波終點(diǎn)開(kāi)始向S波終點(diǎn)方向順序考察信號(hào)查找到第一個(gè)拐點(diǎn)即作為T(mén)波起始點(diǎn)。如果波形變化劇烈,查不到拐點(diǎn),則我們將假設(shè)T波為對(duì)稱波,利用T波終點(diǎn)對(duì)T波起始點(diǎn)進(jìn)行估計(jì)。

        3 算法硬件實(shí)現(xiàn)可行性分析

        本文構(gòu)造算法中采用的小波提升和差分法,由于提升小波變換結(jié)構(gòu)簡(jiǎn)單,且算法執(zhí)行速度快,在心電監(jiān)護(hù)儀等硬件平臺(tái)上只需利用加法器和乘法器即可實(shí)現(xiàn),且逆變化可直接反轉(zhuǎn)實(shí)現(xiàn),很容易在FPGA上實(shí)現(xiàn);差分運(yùn)算只需進(jìn)行簡(jiǎn)單的減法運(yùn)算。另外,去噪和識(shí)別算法涉及的閾值處理,只需利用硬件資源中的比較器實(shí)現(xiàn)。因此,本文提出的算法在硬件算法實(shí)現(xiàn)中只涉及加法、移位和比較操作,便于硬件開(kāi)發(fā)語(yǔ)言描述,易于在基于FPGA等大規(guī)模集成電路器件的硬件平臺(tái)上實(shí)現(xiàn)。由此,本文提出算法的主要用途是將小波提升算法進(jìn)行改進(jìn)并移植到心電監(jiān)護(hù)等硬件平臺(tái)中,也可利用大規(guī)模集成電路芯片對(duì)算法進(jìn)行封裝,增強(qiáng)算法在不同硬件平臺(tái)上的可移植性。

        4 實(shí)驗(yàn)結(jié)果分析

        為驗(yàn)證本文算法的有效性,實(shí)驗(yàn)采用MITBIH心電數(shù)據(jù)庫(kù)中11個(gè)典型的數(shù)據(jù)記錄與文獻(xiàn)[5]中的方法進(jìn)行了對(duì)比試驗(yàn),結(jié)果如表1所示。在試驗(yàn)過(guò)程中,對(duì)于R波位置的定位直接采用心電數(shù)據(jù)庫(kù)注釋文件中的數(shù)據(jù),避免了由于R波位置不準(zhǔn)確帶來(lái)算法檢測(cè)上的偏差。圖2是本文算法對(duì)100.dat心電數(shù)據(jù)進(jìn)行P、T波檢測(cè)的結(jié)果,圖中標(biāo)出了各波的起點(diǎn)、頂點(diǎn)和終點(diǎn)。

        圖2 100.dat心電數(shù)據(jù)P、T波檢測(cè)結(jié)果Fig.2 Theresultsof100.datECGdataP,T-wavedetection

        表1 P波檢測(cè)結(jié)果統(tǒng)計(jì)Table1 TheresultsofPwavedetectionstatistics

        以上兩種方法都是首先根據(jù)注釋文件確定QRS波,再進(jìn)行P波或P、T波的位置查找。從表1可知,兩種方法對(duì)于心電記錄 100,101,103,119等信號(hào)質(zhì)量好,P波明顯時(shí)兩種方法的正確檢測(cè)率都很高。對(duì)于105,205中噪聲干擾較為嚴(yán)重,其部分P波波形不明顯。由于本文算法在利用提升小波進(jìn)行去噪后的相應(yīng)頻域中進(jìn)行差分識(shí)別,故相比較文獻(xiàn)[5],算法準(zhǔn)確率有所提升??傮w上,文獻(xiàn)[5]中方法相對(duì)于本文提出的方法P波檢測(cè)準(zhǔn)確率平均高出0.08%,而本文對(duì)于T波的檢測(cè)方法與P波檢測(cè)策略基本相同,與P波檢測(cè)的準(zhǔn)確率相近。故文中只對(duì)P波識(shí)別準(zhǔn)確率進(jìn)行了對(duì)比試驗(yàn)。

        5 結(jié)論

        本文分析了心電信號(hào)中噪聲頻域特點(diǎn),利用小波的時(shí)-頻域分析能力,將小波去噪算法與差分法相結(jié)合,對(duì)P、T波進(jìn)行識(shí)別。通過(guò)與文獻(xiàn)[5]中方法進(jìn)行對(duì)比試驗(yàn),本文算法在信號(hào)比較干凈時(shí),與文獻(xiàn)[5]算法的檢測(cè)準(zhǔn)確率基本相當(dāng),達(dá)到99.6%以上。在信號(hào)被噪聲污染嚴(yán)重時(shí),本文算法相比文獻(xiàn)[5]的準(zhǔn)確率有近0.6%的提高。同時(shí),由于采用提升小波對(duì)信號(hào)進(jìn)行分解、去噪、重構(gòu),使本文算法比利用傳統(tǒng)小波分析的算法在運(yùn)算速度上提高了一倍。采用對(duì)CA4進(jìn)行差分識(shí)別,將參與運(yùn)算的數(shù)據(jù)變?yōu)樵瓉?lái)信號(hào)數(shù)據(jù)量的1/16。因此,本文提出的算法不但繼承了小波分析對(duì)心電信號(hào)處理的特性,且對(duì)P/T波識(shí)別的準(zhǔn)確性還略有提高,并提高了執(zhí)行速度。由算法特性分析可知,提出的算法不但適用于高性能計(jì)算機(jī)平臺(tái)應(yīng)用,也同樣適用于硬件電路及大規(guī)模集成電路的封裝,使算法面向硬件平臺(tái)的應(yīng)用具有良好的通用性和可移植性。

        [1]范曉東,朱澤煌.心電特征點(diǎn)定位算法[J].北京生物醫(yī)學(xué)工程,1996,15(1):15-18.FanXiao-dong,ZhuZe-huang.An algorithm for locating ECG key points[J].Beijing Biomedical Engineering,1996,15(1):15-18.

        [2]JoengHee-Kyo,KimKwang-Keun,HwangSun-Chul,et al.A new algorithm for P-wave detection inthe ECG signal[C]∥Engineering in Medicineand Biology Society,Proceedings of the Annual International Conference ofthe IEEE Engineering,1989:42.

        [3]楊振野,李玲華,林家瑞.一種基于函數(shù)逼近原理的心電圖P波識(shí)別的研究[J].生物醫(yī)學(xué)工程學(xué)雜志,1998,15(2):120-122.YangZhen-ye,LiLing-hua,LinJia-rui.Approach to recognition of ECG p waves based on approximating functions[J].Journal of Biomedical Engineering,1998,(2):120-122.

        [4] ZhouLi-gao,LiHui-jun,HuGuang-shu.Real-time base-line drift correction and P-wave detection of ECG signal[C]∥Engineering in Medicineand Biology Society,Proceedings of the Annual International Conference of the IEEE Engineering,1989:49.

        [5]季虎,孫即祥,王春光.基于小波變換的自適應(yīng)QRST對(duì)消P波檢測(cè)算法[J].電子與信息學(xué)報(bào),2007,29(8):1868-1871.JiHu,SunJi-xiang,WangChun-guang.An adaptive QRS-T cancellation based on wavelet transform for pwave detection[J].Journal of Electronics & Information Technology,2007,29(8):1868-1871.

        [6]PrasadVV,SiddaiahP,PrabhakaraRaoB.An ew wavelet based method for denoising of biological signals[J].IJCSNS International Journal of Computer Science and Network Security,2008,18(1):238-244.

        [7]Umamaheswara ReddyG,Muralidhar M,Varadarajan S.ECG de-noising using improved thresholding based on wavelet transforms[J].JCSNS International Journal of Computer Scienceand Network Security,2009,9(9):221-225.

        [8]AssalehK.Extraction of fetal electrocardiogram using adaptive neuro-fuzzy inference system[J].IEEE Trans Biomed Eng,2007,54(1):59-68.

        猜你喜歡
        信號(hào)檢測(cè)
        “不等式”檢測(cè)題
        “一元一次不等式”檢測(cè)題
        “一元一次不等式組”檢測(cè)題
        信號(hào)
        鴨綠江(2021年35期)2021-04-19 12:24:18
        “幾何圖形”檢測(cè)題
        “角”檢測(cè)題
        完形填空二則
        孩子停止長(zhǎng)個(gè)的信號(hào)
        小波變換在PCB缺陷檢測(cè)中的應(yīng)用
        基于LabVIEW的力加載信號(hào)采集與PID控制
        国产视频激情视频在线观看| 99精品热6080yy久久| 国产高清天干天天视频| 国产女人av一级一区二区三区| 精品国产天堂综合一区在线| 乱人伦中文无码视频| 91精品全国免费观看青青| 日本高清在线一区二区三区| 久久777国产线看观看精品| 亚洲码国产精品高潮在线| 无国产精品白浆免费视| 美女一区二区三区在线观看视频| 国产乱精品女同自线免费| 久久香蕉国产线看观看精品yw| 国产精品美女| 亚洲日本一区二区在线观看| 一区二区三区视频在线观看免费| 无码国产福利av私拍| 高清国产日韩欧美| 日本道免费一区日韩精品| 国产自拍精品视频免费| 人妻丰满熟妇av无码区不卡| 久久久久久久98亚洲精品| 女同中文字幕在线观看| 久久精品国产成人午夜福利| 日本熟妇色xxxxx欧美老妇| 国产午夜精品久久久久99| 中文字幕亚洲中文第一| 日本一卡二卡3卡四卡免费观影2022| 处破痛哭a√18成年片免费| 亚洲AV无码日韩综合欧亚| 日韩av天堂一区二区三区在线| 特级精品毛片免费观看| 国产精品视频牛仔裤一区| 亚洲av一区二区三区网站| 国产激情久久久久影院小草| 大陆极品少妇内射aaaaa| 男人天堂AV在线麻豆| 国产精品大片一区二区三区四区| 亚洲中文字幕无码中文字在线| 一级呦女专区毛片|