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

        ?

        基于二維特征相似系數(shù)的雷達(dá)輻射源識(shí)別算法

        2020-11-27 07:54:00鑫,王旭,洪
        艦船電子對(duì)抗 2020年5期
        關(guān)鍵詞:特征信號(hào)

        王 鑫,王 旭,洪 偉

        (中國(guó)船舶重工集團(tuán)公司第七二四研究所,江蘇 南京 211153)

        0 引 言

        雷達(dá)輻射源信號(hào)識(shí)別是一項(xiàng)對(duì)采集到的不同型號(hào)雷達(dá)信號(hào)進(jìn)行識(shí)別分類的技術(shù),在現(xiàn)代電子偵察中起到關(guān)鍵作用。隨著電子戰(zhàn)偵察環(huán)境不斷變換,新體制雷達(dá)的不斷投入使用,傳統(tǒng)的目標(biāo)識(shí)別特征參數(shù)如:脈沖重復(fù)周期(PRI)、波達(dá)方向定位(DOA)、二進(jìn)制相移鍵控脈寬(PW)已經(jīng)遠(yuǎn)遠(yuǎn)不能滿足對(duì)雷達(dá)信號(hào)的識(shí)別。

        為了解決上述問題,人們研究發(fā)現(xiàn)對(duì)雷達(dá)輻射源信號(hào)識(shí)別貢獻(xiàn)率最大的特征,即使在其他傳統(tǒng)識(shí)別特征不斷改變的條件下,脈內(nèi)特征的性能十分穩(wěn)定[1-3]。國(guó)內(nèi)外研究人員利用相似系數(shù)、小波變換和脈內(nèi)特征相結(jié)合進(jìn)行識(shí)別,開始對(duì)雷達(dá)輻射源信號(hào)的信號(hào)譜進(jìn)行特征提取與分析[4-7],取得了較好的識(shí)別效果。在文獻(xiàn)[4]中,L. Huadong和H. Jianghong提取了所得信號(hào)的Wigner-Ville分布,將其簡(jiǎn)化為二維空間特征,再利用小波變換對(duì)不同型號(hào)的雷達(dá)輻射源信號(hào)進(jìn)行識(shí)別。文獻(xiàn)[5]中,L. Jian-Dong等人先提取所獲信號(hào)的模糊函數(shù),將其簡(jiǎn)化為二維空間特征,利用相似系數(shù)和KFCM聚類算法進(jìn)行識(shí)別。文獻(xiàn)[6]提取了雷達(dá)信號(hào)的雙譜特征并將其簡(jiǎn)化至二維空間特征,利用相似系數(shù)進(jìn)行信號(hào)識(shí)別。文獻(xiàn)[7]中,J. Han等從二維特征模糊函數(shù)中提取小波,對(duì)雷達(dá)信號(hào)進(jìn)行識(shí)別。

        現(xiàn)有方法可以較好地對(duì)雷達(dá)輻射源信號(hào)進(jìn)行識(shí)別,并且取得了較高的識(shí)別率,但是需要對(duì)信號(hào)多次迭代且處理維數(shù)較高,大大占用了系統(tǒng)資源,降低了信號(hào)識(shí)別的便利性。為了更好地提高不同信噪比條件下的識(shí)別率,對(duì)不同型號(hào)的雷達(dá)輻射源信號(hào)的識(shí)別更為方便可靠,本文提出了一種新的基于二維特征相似系數(shù)的雷達(dá)輻射源信號(hào)的識(shí)別算法。首先,提取雷達(dá)的Wigner 三譜特征,將其簡(jiǎn)化至二維空間特征,之后利用一個(gè)矩陣脈沖序列和一個(gè)三角脈沖序列進(jìn)行雷達(dá)信號(hào)的相似系數(shù)的計(jì)算,從而提取出雷達(dá)信號(hào)的相似系數(shù);最后采用核模糊C均值聚類算法對(duì)不同型號(hào)的雷達(dá)信號(hào)進(jìn)行識(shí)別;最后,通過仿真實(shí)驗(yàn),證明本方法的有效性和可靠性。

        1 雷達(dá)輻射源識(shí)別算法

        1.1 雷達(dá)輻射源的Wigner 三譜及特征相似系數(shù)

        1.1.1 雷達(dá)輻射源的Wigner三譜

        在研究非平穩(wěn)的非高斯信號(hào)的過程中,主要使用的方法是研究原信號(hào)u(t)的K階Winger高階譜。

        文獻(xiàn)[8]介紹了待分析信號(hào)u(t)的高階譜Wku(t,f1,f2,…,fk)的求解:

        (1)

        式中:f為信號(hào)u(t)的頻率;δ為時(shí)差變量,且i=1,2,…,k。

        信號(hào)u(t)的Wigner三譜(WT)即為k=3的Wigner高階譜,求解如下:

        exp(-j2πf2δ2)·exp(-j2πf3δ3)dδ1dδ2dδ3

        (2)

        信號(hào)u(n)的Wigner三譜離散表達(dá)式為:

        exp(-j2π(f1δ1+f2δ+f3δ3))

        (3)

        其中,r4(n,δ1,δ2,δ3)定義如下:

        r4(n,δ1,δ2,δ3)=u*(n-δ)·u(n-δ+δ1)·

        u(n-δ+δ2)·u(n-δ+δ3)

        (4)

        式中:δ=(δ1+δ2+δ3+δ4)/4。

        Wigner三譜的頻率f1=f2=-f3=f,由此可得:

        exp(-j2πf(δ1+δ2-δ3))

        (5)

        在公式(5)的計(jì)算中,根據(jù)奈奎斯特頻率理論,信號(hào)的采樣頻率應(yīng)設(shè)置為大于等于起始信號(hào)頻率的2倍,且頻率軸坐標(biāo)等比例地縮放為原來的0.5倍,以便于計(jì)算與分析。所以,將FFT的長(zhǎng)度設(shè)置為原信號(hào)長(zhǎng)度的4倍,以避免頻譜混疊的狀況出現(xiàn)。

        對(duì)于數(shù)據(jù)長(zhǎng)短的不同,雷達(dá)輻射源的Wigner三譜的計(jì)算方式也不同。對(duì)于短數(shù)據(jù)類型的信號(hào),直接計(jì)算Wigner三譜即可。對(duì)于長(zhǎng)數(shù)據(jù)類型的信號(hào),不妨設(shè)其信號(hào)長(zhǎng)度為N,步驟如下:

        (1) 將信號(hào)長(zhǎng)度N等分為P個(gè)部分,每個(gè)部分采樣點(diǎn)數(shù)為S,即為P=N/S;

        (2) 對(duì)每個(gè)部分的Wigner 三譜[WT(n,f)]進(jìn)行求解;

        (3) 對(duì)每部分的Wigner 三譜計(jì)算進(jìn)行合成,求出雷達(dá)輻射源信號(hào)的Wigner 三譜的平均值:

        (6)

        式中:f=1,…,L,L為FFT的長(zhǎng)度;n=1,…,S;WTp(n,f)為S×L大小的矩陣。

        由上式可知,Wigner三譜是關(guān)于采樣頻率f和信號(hào)長(zhǎng)度n的表達(dá)式,因此采用三維特征圖來描繪Wigner三譜。

        三維特征圖能較完整地反映Wigner三譜的特征,但是完整的三維特征在計(jì)算時(shí)需要大量的計(jì)算,這大大增加了系統(tǒng)的計(jì)算任務(wù),提高了數(shù)據(jù)量,降低了運(yùn)算速度;并且由于數(shù)據(jù)量的增加給目標(biāo)識(shí)別帶來了困難,還是不能很好地利用它來對(duì)雷達(dá)輻射源信號(hào)進(jìn)行甄別。因此,將Wigner三譜簡(jiǎn)化至二維空間特征進(jìn)行分析,可以降低系統(tǒng)的復(fù)雜程度,提高信號(hào)識(shí)別的準(zhǔn)確率。

        1.1.2 簡(jiǎn)化Wigner三譜特征

        三維空間特征圖不能快速有效地反映不同雷達(dá)信號(hào)之間的區(qū)別,因此我們將Wigner三譜簡(jiǎn)化為二維空間特征進(jìn)行觀察。簡(jiǎn)化后的Wigner三譜要滿足2個(gè)條件:(1)簡(jiǎn)化后的二維空間特征仍然可以完整地反映Wigner三譜的特性;(2)計(jì)算量大大降低,可以較快地甄別不同型號(hào)的雷達(dá)信號(hào)。

        對(duì)于三維特征圖定義:沿著X軸、平行于YoZ的平面截取獲得的特征即為Wigner三譜的簡(jiǎn)化,且將其中最大的Wigner 三譜值作為二維空間的特征向量。這種截取手段結(jié)合獲取最大值的方法,可以保證截取到的Wigner三譜較好地反映譜圖特征,并且對(duì)該特征的處理能對(duì)雷達(dá)信號(hào)識(shí)別起到直觀、積極的作用。在對(duì)每一個(gè)雷達(dá)輻射源信號(hào)的二維空間特征進(jìn)行分解時(shí),采樣點(diǎn)數(shù)設(shè)置為S=512,有利于提高算法的處理速度。將連續(xù)波(CW)信號(hào)、非線性調(diào)頻(NLFM)信號(hào)、正交相移鍵控(QPSK)信號(hào)進(jìn)行二維簡(jiǎn)化,分別如圖1(a)、圖1(b)、圖1(c)所示,信號(hào)之間的差異性可以通過二維空間特征圖直觀地反映出來。

        圖1 3種雷達(dá)輻射源信號(hào)的二維空間特征

        1.1.3 提取二維特征的相似系數(shù)

        二維空間特征的采樣點(diǎn)數(shù)為512個(gè),但是對(duì)于實(shí)際計(jì)算過程而言,運(yùn)算量仍然較大,為了保證較低的信號(hào)識(shí)別率,需要對(duì)采樣點(diǎn)數(shù)進(jìn)行進(jìn)一步改進(jìn)。采用相似系數(shù)的方法來降低采樣點(diǎn)數(shù),使其達(dá)到最小化效果。

        2個(gè)一維離散序列{x1(j),j=1,…,N},{x2(j),j=1,…,N}組成相似系數(shù)Cr,且兩序列均為正信號(hào)序列,即(x1(j)≥0,x2(j)≥0),則相似系數(shù)Cr表達(dá)式為:

        (7)

        式中:序列x1(j),x2(j)均不為零,且Cr∈[0,1]。

        由于簡(jiǎn)化后的二維特征得到的Wigner 三譜信息中不可避免地包含了冗余信息的干擾,雷達(dá)信號(hào)的識(shí)別過程難度并未明顯降低。為此,我們?nèi)孕鑼?duì)二維特征進(jìn)行降維運(yùn)算。

        由圖1可知,通過降維可以得到不同雷達(dá)輻射源信號(hào)Wigner三譜的二維特征,它們的特征有著比較大的差別,為了提取這些特征的差異,采用了相似系數(shù)的方法,步驟如下:

        (1) 采用矩形脈沖信號(hào)和三角脈沖信號(hào)為參考信號(hào);

        (2) 利用參考信號(hào),分別計(jì)算出Wigner三譜的二維特征與矩形脈沖序列、三角脈沖序列的相似系數(shù)分別記為Cr1、Cr2;

        (3) 分別用W(j)、T(j)來表示矩形脈沖采樣序列和三角脈沖采樣序列,則它們可以表示為:

        (8)

        (9)

        通過式(8)和式(9),可以用計(jì)算[Cr1,Cr2]的相似系數(shù)來取代計(jì)算采樣點(diǎn)數(shù)多達(dá)512個(gè)的雷達(dá)輻射源信號(hào)。采用相似系數(shù)的方法大大降低了采樣信號(hào)的維數(shù),使雷達(dá)輻射源信號(hào)的識(shí)別過程更加高效、便利。

        1.2 核模糊C均值聚類算法

        1.2.1 模糊C均值聚類算法

        模糊C均值聚類算法(FCM),是在C均值聚類算法的基礎(chǔ)上演變發(fā)展而來。它利用隸屬度的原則,將信號(hào)的數(shù)據(jù)點(diǎn)劃分為不同類別。模糊C均值聚類算法步驟如下:

        (1) 用c個(gè)模糊的組別來劃分n個(gè)向量xi(i=1,2,…,n);

        (2) 找出每個(gè)組別的聚類中心,同時(shí)最小化其他不是類似性指標(biāo)的目標(biāo)函數(shù)。在FCM中,給定任意已知數(shù)據(jù)的情況下,利用模糊函數(shù)進(jìn)行分類,數(shù)據(jù)劃分進(jìn)每組的可能性大小用隸屬度來表示,隸屬度取值為(0,1)之間的任意數(shù),這也是FCM與HCM最顯著的差別所在。隸屬矩陣用U來表示,在進(jìn)行歸一化后,c個(gè)組別的隸屬度之和為恒定值1,表達(dá)為:

        (10)

        則模糊C均值聚類算法的目標(biāo)函數(shù)可用下式表達(dá):

        (11)

        (12)

        其中,在第i個(gè)模糊組中,ci表示聚類中心,隸屬度uij∈[0,1],m∈[1,∞)表示加權(quán)系數(shù),dij=‖ci-xj‖表示從聚類中心ci到第j個(gè)點(diǎn)的歐式距離。由此可得代價(jià)函數(shù):

        (13)

        上式可得式(12)最小化的條件。其中,式(10)的n個(gè)約束條件中,拉格朗日因子用λj,j=1,…,n表示。使式(12)取最小值的條件則可進(jìn)一步用下式表示:

        (14)

        (15)

        由式(14)和式(15)可得,模糊C均值聚類算法就是一個(gè)不斷的迭代更新的過程。在對(duì)數(shù)據(jù)的處理運(yùn)算過程中,聚類中心ci和隸屬度uij的取值非常重要,過程如下:

        (1) 初始化隸屬矩陣U,使每個(gè)樣本集整體之和恒為1,即滿足式(10);

        (2) 按式(12)和式(14)求得目標(biāo)函數(shù)和c個(gè)聚類中心ci,i=1,…,c;

        (3) 判斷目標(biāo)函數(shù)值,若小于閾值,或者函數(shù)值與上一次目標(biāo)函數(shù)值的差值的絕對(duì)值大小小于閾值大小,則迭代結(jié)束;否則,利用式(15)求解一個(gè)新的隸屬矩陣U。重復(fù)步驟(2)。

        上述過程中,可以調(diào)換聚類中心的初始化過程和迭代順序。然而由于初始聚類中心的選擇不同,導(dǎo)致算法得到的解是否是最優(yōu)解得不到保證,因此對(duì)于初始的聚類中心的取值,可以每次取不同的初始值多次運(yùn)算,也可以用其他算法來獲得。

        1.2.2 核模糊C均值聚類算法

        核模糊C均值聚類算法(KFCM)是由模糊C均值聚類算法進(jìn)化而來的基于高斯核函數(shù)的分類算法。利用高維特征空間H來反映原樣本集數(shù)據(jù)X的投影結(jié)果,該映射過程可以表示為φ:x→φ(x),φ(x)∈H,且在該特征空間H中,卷積和歐式距離可以表達(dá)為:

        k(x,y)=<φ(x),φ(y)>=φT(x)φ(y)

        (16)

        (17)

        其中,高斯核函數(shù)表達(dá)式為:

        (18)

        式中:δ表示高斯核函數(shù)的寬度。

        利用式(18)對(duì)(17)進(jìn)行化簡(jiǎn),可得:

        ‖φ(x)-φ(y)‖2=2-k(x,y)

        (19)

        其中,k(x,x)=k(y,y)=1。則核模糊C均值聚類算法的目標(biāo)函數(shù)可寫為:

        (20)

        核模糊C均值聚類算法的流程如下:

        (1) 設(shè)置初始值條件:分類組別數(shù)目最大值為cmax,迭代次數(shù)最大值為kmax,迭代結(jié)束條件ε≥0。

        (2) 將樣本數(shù)據(jù)集X進(jìn)行歸一化處理,求解到高斯核矩陣K=[k(xj,yi)]。

        (3) 將隸屬度矩陣初始化,即為U=[uij]c×n。

        (4) 不斷更新迭代隸屬度矩陣U和聚類中心vi,表達(dá)式如下:

        (21)

        (22)

        (5) 求得第k次迭代更新的代價(jià)函數(shù)值Jφ(k)。

        (6) 進(jìn)行判斷。若滿足兩點(diǎn)之間歐式距離‖v(k+1)-v(k)‖<ε或迭代次數(shù)k>kmax,執(zhí)行(7);否則,令k=k+1,跳至(4)繼續(xù)執(zhí)行。

        (7) 若迭代次數(shù)c

        1.3 算法總體設(shè)計(jì)

        根據(jù)前文所述,基于二維特征相似系數(shù)的信號(hào)識(shí)別算法流程如圖2所示,主要分為雷達(dá)信號(hào)輸入、Wigher二維特征提取、信號(hào)相似系數(shù)提取、核模糊C均值聚類、雷達(dá)輻射源識(shí)別。雷達(dá)信號(hào)輸入部分用于接收經(jīng)過信號(hào)調(diào)理的雷達(dá)信號(hào),Wigher二維特征提取用于對(duì)接收信號(hào)進(jìn)行Wigher三譜的捕獲及二維簡(jiǎn)化,信號(hào)相似性提取用于求解信號(hào)的相似系數(shù)Cr1和Cr2,核模糊C均值聚類用于對(duì)處理好的二維特征和相似系數(shù)進(jìn)行分類,雷達(dá)輻射源識(shí)別根據(jù)核模糊C均值聚類的結(jié)果對(duì)輸入信號(hào)進(jìn)行識(shí)別分類。

        圖2 算法整體流程圖

        輸入的雷達(dá)信號(hào)經(jīng)過Wigher二維特征提取和相似系數(shù)的提取可以更好地凸顯不同信號(hào)的特征,采用核模糊C均值聚類法對(duì)不同信號(hào)特征進(jìn)行有效分類,可以有效解決復(fù)雜電磁環(huán)境中不同型號(hào)雷達(dá)的分類問題,同時(shí)將Wigher三譜轉(zhuǎn)換到二維空間,大大降低了信號(hào)處理的數(shù)據(jù)量,提升了處理速度,有利于硬件系統(tǒng)的實(shí)現(xiàn),從而可以高效快速地提高信號(hào)識(shí)別率和電子偵察效率。

        2 實(shí)驗(yàn)與仿真

        采用6種不同的雷達(dá)輻射源信號(hào),它們分別為:連續(xù)波(CW)雷達(dá)信號(hào)、線性調(diào)頻(LFM)雷達(dá)信號(hào)、移頻鍵控(FSK)雷達(dá)信號(hào)、二進(jìn)制相移鍵控(BPSK)信號(hào)、正交相移鍵控(QPSK)信號(hào)、非線性調(diào)頻(NLFM)信號(hào)。分別以5 dB、10 dB、15 dB、20 dB的信噪比為條件,每個(gè)雷達(dá)輻射源信號(hào)產(chǎn)生100個(gè)脈沖信號(hào)。每個(gè)信號(hào)的采樣點(diǎn)數(shù)設(shè)置為512個(gè),采樣頻率為200 MHz,線性調(diào)頻雷達(dá)信號(hào)帶寬為10 MHz,移頻鍵控雷達(dá)信號(hào)的頻率分別為10 MHz和15 MHz,非線性調(diào)頻信號(hào)采用正弦頻率調(diào)制雷達(dá)信號(hào),其他信號(hào)載頻為15 MHz。對(duì)這些信號(hào)進(jìn)行均值估計(jì),之后分別以信噪比5 dB~20 dB間每隔1.5 dB為條件,計(jì)算信號(hào)的相似系數(shù),計(jì)算結(jié)果如圖3、圖4所示。

        圖3 信噪比與相似系數(shù)Cr1之間的關(guān)系

        圖4 信噪比與相似系數(shù)Cr2之間的關(guān)系

        將連續(xù)波雷達(dá)信號(hào)、線性調(diào)頻雷達(dá)信號(hào)、移頻鍵控雷達(dá)信號(hào)、二進(jìn)制相移鍵控信號(hào)、正交相移鍵控信號(hào)、非線性調(diào)頻信號(hào)分別用R1~R6表示,以相似系數(shù)Cr1和Cr2作為參數(shù),利用核模糊C均值聚類算法對(duì)它們進(jìn)行識(shí)別分類,結(jié)束條件設(shè)定為ε<10-5,得出不同信噪比條件下雷達(dá)輻射源信號(hào)的識(shí)別率,如表1所示。

        表1 雷達(dá)輻射源信號(hào)在不同信噪比條件下的識(shí)別率

        從表1可以看出,利用相似系數(shù)法可以較好地將不同類型的雷達(dá)輻射源信號(hào)進(jìn)行分離識(shí)別。信噪比較低的情況下,樣本會(huì)發(fā)生相互混疊的情況,例如:在信噪比為5 dB時(shí),R4~R6雷達(dá)輻射源信號(hào)的識(shí)別率有所下降;隨著信噪比不斷增大,雷達(dá)信號(hào)的識(shí)別率也逐漸提高;在信噪比為10 dB時(shí),每種雷達(dá)信號(hào)可以有效分離出來,并且總體可以達(dá)到98.5%的識(shí)別率;當(dāng)信噪比達(dá)到15 dB和20 dB時(shí),識(shí)別率大大提高,可以達(dá)到將近100%的識(shí)別率。在實(shí)際運(yùn)用中,即使是在信噪比為5 dB的情況下,對(duì)R6雷達(dá)信號(hào)的識(shí)別率低至76%,這個(gè)識(shí)別率也是相當(dāng)可觀的。

        將本文提出的方法與文獻(xiàn)[10]和文獻(xiàn)[11]中所提出的模糊函數(shù)法和Wigner雙譜法進(jìn)行識(shí)別率的比較,結(jié)果如圖5所示。從圖5中可以看出,本文提出的方法在信噪比為7 dB時(shí),雷達(dá)輻射源信號(hào)的識(shí)別率有較大幅度的提升;當(dāng)信噪比較高時(shí),雷達(dá)輻射源信號(hào)識(shí)別率提高了0.5%以上;與現(xiàn)有方法相比,本文提出的方法可以在復(fù)雜電磁環(huán)境中更好地實(shí)現(xiàn)不同雷達(dá)輻射源信號(hào)的識(shí)別,提高了雷達(dá)系統(tǒng)的識(shí)別性能。

        圖5 本文方法與模糊函數(shù)法和Wigner雙譜法的效果對(duì)比

        本文提出了一種基于二維特征相似系數(shù)的雷達(dá)輻射源識(shí)別算法,用于在復(fù)雜電磁環(huán)境中實(shí)現(xiàn)雷達(dá)輻射源信號(hào)識(shí)別率的提升,將雷達(dá)信號(hào)的Wigner-Ville特征映射到二維空間,使用矩形采樣序列和三角采樣序列提取信號(hào)相關(guān)系數(shù),采用核模糊C均值聚類方法識(shí)別出不同雷達(dá)輻射源信號(hào)。通過實(shí)驗(yàn)證明了對(duì)于常用雷達(dá)信號(hào)(如續(xù)波雷達(dá)信號(hào)、線性調(diào)頻雷達(dá)信號(hào)、移頻鍵控雷達(dá)信號(hào)、二進(jìn)制相移鍵控信號(hào)、正交相移鍵控信號(hào)、非線性調(diào)頻信號(hào)),本方法可以較好地對(duì)它們進(jìn)行分離識(shí)別,隨著信噪比不斷增大,信號(hào)的識(shí)別率也逐漸提高。

        猜你喜歡
        特征信號(hào)
        抓住特征巧觀察
        信號(hào)
        鴨綠江(2021年35期)2021-04-19 12:24:18
        完形填空二則
        新型冠狀病毒及其流行病學(xué)特征認(rèn)識(shí)
        如何表達(dá)“特征”
        不忠誠(chéng)的四個(gè)特征
        孩子停止長(zhǎng)個(gè)的信號(hào)
        抓住特征巧觀察
        基于LabVIEW的力加載信號(hào)采集與PID控制
        一種基于極大似然估計(jì)的信號(hào)盲抽取算法
        亚洲不卡免费观看av一区二区| 久久久诱惑一区二区三区| 无码骚夜夜精品| 亚洲欧美日韩中文在线制服| 欧美 日韩 国产 成人 在线观看 | 国产午夜视频在永久在线观看| 亚洲AV无码国产永久播放蜜芽| 第十色丰满无码| 国产三级在线看完整版| av高清视频在线麻豆免费观看| 蜜桃网站入口可看18禁| 极品一区二区在线视频| 亚洲国产精品一区二区成人片国内 | 韩国三级黄色一区二区| 日本伊人精品一区二区三区| 天天做天天爱夜夜爽毛片毛片| 日本真人做人试看60分钟| 亚洲av无码一区二区三区人妖| 国产在线无码一区二区三区| 不卡高清av手机在线观看| 精品日韩欧美| 丰满少妇又爽又紧又丰满动态视频| 国产一区二区三区四区在线视频| 国产裸体美女永久免费无遮挡| 48久久国产精品性色aⅴ人妻 | 久久精品国产亚洲av热一区| 国产免费观看久久黄av麻豆 | 国产在线视频一区二区三| 水蜜桃在线精品视频网| 国产午夜精品无码| 亚洲av无码一区二区二三区 | 国产精品免费_区二区三区观看| a级大胆欧美人体大胆666| 国产免费久久精品国产传媒| 亚洲VA中文字幕无码毛片春药| 色综合久久久久综合一本到桃花网| 精品国产麻豆免费人成网站| 国产精品狼人久久影院软件介绍| 无码无套少妇毛多18pxxxx| 女人被男人躁得好爽免费视频| 麻豆精产国品|