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

基于交互式多模型粒子濾波的相控陣?yán)走_(dá)自適應(yīng)采樣

2012-09-27 01:41:54郁衛(wèi)華朱曉華
電子設(shè)計(jì)工程 2012年5期
關(guān)鍵詞:模型

郁衛(wèi)華 , 朱 翔 , 朱曉華

(1.南通農(nóng)業(yè)職業(yè)技術(shù)學(xué)院 江蘇 南通 226007;2.南京理工大學(xué) 電子工程與光電技術(shù)學(xué)院,江蘇 南京 210094;3.國網(wǎng)電力科學(xué)研究院 江蘇 南京 210003)

基于交互式多模型粒子濾波的相控陣?yán)走_(dá)自適應(yīng)采樣

郁衛(wèi)華1,2, 朱 翔3, 朱曉華2

(1.南通農(nóng)業(yè)職業(yè)技術(shù)學(xué)院 江蘇 南通 226007;2.南京理工大學(xué) 電子工程與光電技術(shù)學(xué)院,江蘇 南京 210094;3.國網(wǎng)電力科學(xué)研究院 江蘇 南京 210003)

為有效合理利用雷達(dá)資源和解決雷達(dá)測(cè)量值與運(yùn)動(dòng)狀態(tài)間的非線性關(guān)系以及目標(biāo)狀態(tài)本身可能出現(xiàn)的非線性,提出了一種基于交互式多模型粒子濾波(IMMPF)的相控陣?yán)走_(dá)自適應(yīng)采樣目標(biāo)跟蹤方法。將交互式多模型粒子濾波一步預(yù)測(cè)值的后驗(yàn)克拉美羅矩陣代替預(yù)測(cè)協(xié)方差矩陣,通過該矩陣的跡與某一門限值比較來更新采樣周期以適應(yīng)目標(biāo)運(yùn)動(dòng)狀態(tài)的變化。將該方法與基于量測(cè)轉(zhuǎn)換的IMM自適應(yīng)采樣算法進(jìn)行仿真實(shí)驗(yàn),表明了該算法的有效性。

相控陣?yán)走_(dá);交互式多模型;粒子濾波;后驗(yàn)克拉美羅矩陣(PCRB);自適應(yīng)采樣

相控陣?yán)走_(dá)采用電子掃描的方式使得天線波束指向無慣性的捷變,實(shí)現(xiàn)了雷達(dá)跟蹤目標(biāo)時(shí)可根據(jù)目標(biāo)運(yùn)動(dòng)狀態(tài)自適應(yīng)地改變采樣周期[1]。當(dāng)目標(biāo)機(jī)動(dòng)時(shí),采用高采樣率,實(shí)現(xiàn)對(duì)目標(biāo)進(jìn)行更精確跟蹤,以避免跟蹤丟失。反之,當(dāng)目標(biāo)處于平穩(wěn)狀態(tài)時(shí),在一定跟蹤誤差內(nèi),采用一個(gè)相對(duì)低的采樣率。這種依據(jù)目標(biāo)運(yùn)動(dòng)狀態(tài)自適應(yīng)改變采樣周期的手段可在保證基本跟蹤誤差的前提下,更合理有效地利用雷達(dá)資源,充分發(fā)揮相控陣?yán)走_(dá)的多目標(biāo)跟蹤能力[2]。

Bar-Shalom提出的交互式多模型(IMM)算法是目前機(jī)動(dòng)目標(biāo)跟蹤中最有效的算法之一。但由于雷達(dá)測(cè)量值與目標(biāo)運(yùn)動(dòng)狀態(tài)間的高度非線性關(guān)系以及目標(biāo)狀態(tài)本身可能出現(xiàn)的非線性,使得卡爾曼濾波或擴(kuò)展卡爾曼濾波算法性能大大下降。Gordon提出的粒子濾波(PF)作為一種基于貝葉斯估計(jì)的非線性濾波算法,在處理非高斯非線性時(shí)變系統(tǒng)的參數(shù)估計(jì)和狀態(tài)濾波問題方面有獨(dú)到的優(yōu)勢(shì)。融合IMM和PF的新濾波方法,即交互式多模型粒子濾波器(IMMPF),可同時(shí)解決雷達(dá)測(cè)量值與運(yùn)動(dòng)狀態(tài)間的非線性關(guān)系和目標(biāo)機(jī)動(dòng)運(yùn)動(dòng)帶來的跟蹤困難,因此該方法的濾波性能要遠(yuǎn)遠(yuǎn)優(yōu)于傳統(tǒng)的IMM方法。

文中基于交互式多模型粒子濾波(IMMPF)的目標(biāo)跟蹤算法,提出將交互式多模型粒子濾波一步預(yù)測(cè)值的后驗(yàn)克拉美羅矩陣(PCRB)代替預(yù)測(cè)協(xié)方差,選取該矩陣的跡不超過某一門限值的最大采樣間隔作為下一時(shí)刻的采樣間隔,使得雷達(dá)在滿足要求的跟蹤誤差情況下,實(shí)現(xiàn)雷達(dá)資源的有效利用。文章首先分析了IMMPF算法,接著給出了IMMPF的PCRB迭代更新方法,然后提出基于IMMPF的自適應(yīng)采樣算法,最后進(jìn)行了計(jì)算機(jī)仿真實(shí)驗(yàn)。

1 交互式多模型粒子濾波

1.1 粒子濾波

粒子濾波算法(PF)是一種基于蒙特卡洛仿真的近似貝葉斯濾波算法[3],其核心思想是使用一組具有相應(yīng)加權(quán)值的隨機(jī)樣本(粒子)來表示后驗(yàn)概率密度(PDF)[4-5]。它首先在狀態(tài)空間中根據(jù)先驗(yàn)分布產(chǎn)生一組粒子,然后在觀測(cè)的基礎(chǔ)上通過調(diào)節(jié)每個(gè)粒子所對(duì)應(yīng)權(quán)值的大小和位置,來獲得服從實(shí)際分布的粒子,并以該組粒子的均值作為系統(tǒng)狀態(tài)估計(jì)值,最后對(duì)這組粒子權(quán)重采樣以保證權(quán)重的均勻分布。

目標(biāo)跟蹤問題可以描述為在給定一組觀測(cè)的條件下對(duì)系統(tǒng)的狀態(tài)進(jìn)行估計(jì)。在隨機(jī)系統(tǒng)離散的形式下,目標(biāo)運(yùn)動(dòng)的狀態(tài)方程可表述為:

其中xk為狀態(tài)向量,vk為與狀態(tài)無關(guān)的過程噪聲,fk表示系統(tǒng)狀態(tài)轉(zhuǎn)移模型。

相控陣?yán)走_(dá)系統(tǒng)的測(cè)量方程可表述為:

其中nk表示獨(dú)立于系統(tǒng)測(cè)量的熱噪聲,hk表示系統(tǒng)測(cè)量模型。

設(shè)已知系統(tǒng)狀態(tài)的初始概率密度函數(shù)為 p(x0|z0)=p(x0),則狀態(tài)預(yù)測(cè)方程為:

狀態(tài)更新方程為:

其中,

遞推求解式(3)、(4)可得到貝葉斯最優(yōu)估計(jì),但在實(shí)際系統(tǒng)中,很難求出PDF。基于隨機(jī)采樣運(yùn)算的蒙特卡羅方法可將積分運(yùn)算轉(zhuǎn)化為有限樣本點(diǎn)的求和運(yùn)算,可以解決上述問題。如果能從PDF上隨機(jī)抽取一組粒子來逼近,則PDF可用這組粒子來逼近。通常很難直接從PDF上抽樣粒子,但可從一已知的分布上隨機(jī)抽取一組粒子來逼近PDF。

若將重要密度函數(shù)q(x0:k|z1:k)可以寫成如下遞推形式:

把公式(4)、(5)、(8)代入到公式(7)中,于是權(quán)值的遞推形式可表示為:

預(yù)測(cè)、更新和權(quán)值估計(jì)就構(gòu)成了基本的序貫重要抽樣(SIS)法,這是粒子濾波算法的基本框架[6]。粒子濾波依靠重要抽樣,因而要求重要密度分布盡可能地逼近PDF。

1.2 交互式多模型粒子濾波算法

1)隨機(jī)抽取粒子

2)輸入交互

3)模型匹配粒子濾波

4)模型概率更新

5)殘差重抽樣

6)輸出交互

對(duì)m個(gè)模型的各相應(yīng)的粒子群進(jìn)行輸出交互,然后對(duì)所有的粒子求加權(quán)平均和:

2 基于IMMPF的自適應(yīng)采樣

根據(jù)IMMPF算法可以看出,粒子濾波過程中并不能直接求出預(yù)測(cè)協(xié)方差矩陣,因此無法建立預(yù)測(cè)協(xié)方差矩陣與采樣間隔的關(guān)系,導(dǎo)致無法直接利用IMMPF進(jìn)行自適應(yīng)采樣。在參數(shù)估計(jì)理論中,對(duì)任何無偏估計(jì)量的方差存在一個(gè)下限。盡管存在許多這樣的限,但克拉美羅下限是最容易確定且最有效的一種。文獻(xiàn)[8]中給出了非線性濾波下后驗(yàn)克拉美羅界(PCRB)及其迭代更新公式,因此將預(yù)測(cè)值代入PCRB中就可以得到該預(yù)測(cè)值的估計(jì)質(zhì)量。將預(yù)測(cè)協(xié)方差矩陣用預(yù)測(cè)值的PCRB矩陣來替代,通過建立該矩陣與采樣間隔的關(guān)系,就可以實(shí)現(xiàn)基于IMMPF的自適應(yīng)采樣。

2.1 交互式多模型粒子濾波下的后驗(yàn)克拉美羅界(PCRB)

假設(shè)在k時(shí)刻考慮目標(biāo)跟蹤問題,根據(jù)上述IMMPF濾波算法通過雷達(dá)測(cè)量值z(mì)1:k來估計(jì)目標(biāo)軌跡x0:k,那么目標(biāo)軌跡的貝葉斯信息矩陣(Bayesian Information Matrix BIM,其逆矩陣為PCRB)可以表示為:

根據(jù)文獻(xiàn)[8]中PCRB更新迭代公式,單模型下更新k+1時(shí)刻的BIM為的遞推公式為:

其中Q為過程噪聲協(xié)方差矩陣,F(xiàn)為系統(tǒng)狀態(tài)轉(zhuǎn)移矩陣,Γk+1為測(cè)量誤差所引起的協(xié)方差矩陣。因此交互式多模型下更新k+1時(shí)刻的Jk+1的遞推公式表示為

其中uj為模型j的模型更新概率,Qj為模型j的過程噪聲協(xié)方差矩陣,F(xiàn)j為模型j的線性系統(tǒng)狀態(tài)轉(zhuǎn)移矩陣或非線性狀態(tài)轉(zhuǎn)移矩陣的一階導(dǎo)數(shù)。當(dāng)采樣周期變化時(shí),則Fj也隨之變化,從而Jk+1變化。Γk+1為雷達(dá)測(cè)量誤差所引起的協(xié)方差矩陣,其定義為

由式(20)、(21)可知,目標(biāo)軌跡的 PCRB矩陣與狀態(tài)模型和觀測(cè)模型都相關(guān)的,通過后者可以算出Γk+1。但由于模型采用粒子濾波時(shí),則需要計(jì)算每一粒子的Γk+1,計(jì)算量很大。根據(jù)文獻(xiàn)[9],Γk+1可表示為

其中Ξk+1是由測(cè)量值z(mì)k+1來估計(jì)狀態(tài)值xk+1的標(biāo)準(zhǔn)的費(fèi)舍信息矩陣(Fisher Information Matrix,F(xiàn)IM)。在自適應(yīng)采樣過程中,由于zk+1在k時(shí)刻是未知的,因此采用預(yù)測(cè)值代替量測(cè)值來計(jì)算。

2.2 自適應(yīng)采樣原則

目標(biāo)狀態(tài)預(yù)測(cè)誤差能很大程度上反映目標(biāo)跟蹤的情況,因此常用預(yù)測(cè)協(xié)方差門限法來實(shí)現(xiàn)雷達(dá)系統(tǒng)的自適應(yīng)采樣控制。本文依據(jù)該思想,用下一時(shí)刻預(yù)測(cè)PCRB矩陣的跡與給定門限值比較,選擇最佳下一時(shí)刻雷達(dá)采樣間隔。

具體算法是,當(dāng)計(jì)算得到的下一時(shí)刻預(yù)測(cè)PCRB矩陣的跡小于給定的門限值時(shí),則雷達(dá)的采樣間隔增加一個(gè)步進(jìn),直到PCRB矩陣的跡超過門限值時(shí),則此時(shí)雷達(dá)的采樣間隔的前一步進(jìn)即為最終輸出的下一個(gè)時(shí)刻的雷達(dá)采樣間隔。反之,當(dāng)計(jì)算得到的下一時(shí)刻預(yù)測(cè)PCRB的跡超過給定的門限值時(shí),則雷達(dá)采樣間隔減少一個(gè)步進(jìn),直到預(yù)測(cè)PCRB的跡小于門限值時(shí),同樣此時(shí)雷達(dá)的采樣間隔的前一步進(jìn)為最終輸出的下一個(gè)時(shí)刻的雷達(dá)采樣間隔。

假設(shè)當(dāng)前時(shí)刻為k,采樣間隔為Tk,門限值為Pth算法可表示如下:

這樣,對(duì)于穩(wěn)定跟蹤的目標(biāo),目標(biāo)跟蹤精度較高,則采樣間隔可以增大,有效地節(jié)約了雷達(dá)資源。當(dāng)目標(biāo)處于機(jī)動(dòng)狀態(tài)時(shí),雷達(dá)的跟蹤精度下降,迫使目標(biāo)狀態(tài)的預(yù)測(cè)PCRB矩陣增大,超過門限值,則采樣間隔相應(yīng)減小以提高雷達(dá)跟蹤精度。

3 仿真實(shí)驗(yàn)與結(jié)果分析

為了驗(yàn)證基于IMMPF的相控陣?yán)走_(dá)自適應(yīng)采樣方法的有效,本文選取典型的目標(biāo)機(jī)動(dòng)環(huán)境進(jìn)行仿真,并將仿真結(jié)果與基于量測(cè)轉(zhuǎn)換的IMM自適應(yīng)采樣方法[10]進(jìn)行比較。假定目標(biāo)運(yùn)動(dòng)在X-Y平面內(nèi)運(yùn)動(dòng)120 s,高度保持不變,起始位置(2 500 m,8 000 m)。在t=0 s開始以100 m/s沿X軸勻速運(yùn)動(dòng),在t=21 s時(shí)以1.5 g加速度機(jī)動(dòng)轉(zhuǎn)彎,在t=36 s時(shí)勻加速運(yùn)動(dòng),在t=66 s時(shí)又恢復(fù)勻速運(yùn)動(dòng)。目標(biāo)整個(gè)過程的運(yùn)動(dòng)軌跡如圖1所示。仿真過程中假設(shè)相控陣?yán)走_(dá)的測(cè)距誤差Pr=100 m,測(cè)角誤差 ρθ=0.03 rad。 自適應(yīng)采樣間隔為 0.8 s、1.0 s、1.2 s、1.4 s、1.6 s、1.8 s和2.0 s。采用勻速、勻加速和常速率轉(zhuǎn)彎3種模型,模型先驗(yàn)概率 Mo=[1/3,1/3,1/3],Markov模型轉(zhuǎn)移概率

圖1 目標(biāo)的運(yùn)動(dòng)軌跡Fig.1 The trajectory of target

圖2 平均更新時(shí)間1.2 s下IMMPF和IMM算法X軸和Y軸跟蹤誤差Fig.2 Tracking error of X-axis and Y-axis using IMMPF and IMM algorithm on average update time 1.2 s

圖3 相同門限下IMMPF和IMM的更新時(shí)間Fig.3 Time interval using IMMPF and IMM algorithm on the same threshold

平均更新時(shí)間1.2 s下IMMPF和IMM算法X軸和Y軸跟蹤誤差如圖2所示。由圖可見,IMMPF算法的跟蹤誤差要明顯小于基于量測(cè)轉(zhuǎn)換的IMM算法。在初始非機(jī)動(dòng)時(shí)間段,由于PF比量測(cè)轉(zhuǎn)換在非線性量測(cè)方程下具有更好的性能,因此IMMPF在初始階段的跟蹤誤差就小于IMM。當(dāng)目標(biāo)突然開始機(jī)動(dòng)時(shí),由于IMMPF能從粒子群中選擇能較好反映目標(biāo)運(yùn)動(dòng)狀態(tài)的粒子并進(jìn)行復(fù)制,而且各模型間的粒子群還能存在交互,因而IMMPF可較快的適應(yīng)機(jī)動(dòng),導(dǎo)致IMM跟蹤誤差與IMMPF相比進(jìn)一步擴(kuò)大。

相同門限下IMMPF和IMM算法的更新時(shí)間如圖3所示。由于IMMPF本身的跟蹤性能很好,在目標(biāo)機(jī)動(dòng)段通過減小更新時(shí)間來保持跟蹤質(zhì)量,在非機(jī)動(dòng)段增大更新時(shí)間。而IMM由于其跟蹤性能有限,在非機(jī)動(dòng)段無法增大更新時(shí)間,限制了其自適應(yīng)采樣的性能。在整個(gè)時(shí)間段中,IMMPF的平均更新時(shí)間為1.38 s,IMM的平均更新時(shí)間1.2 s。在相同時(shí)間內(nèi),IMMPF所需的采樣次數(shù)比IMM要少15%,因而IMMPF具有比IMM更好的自適應(yīng)采樣性能。

4 結(jié) 論

雷達(dá)采樣間隔時(shí)間,對(duì)相控陣?yán)走_(dá)多目標(biāo)跟蹤性能影響很大,因此自適應(yīng)的采樣間隔不僅合理有效地利用了雷達(dá)資源,同時(shí)對(duì)確保跟蹤的連續(xù)性,可靠性和跟蹤精度有重要意義。本文提出的基于IMMPF的自適應(yīng)采樣算法能有效地根據(jù)目標(biāo)運(yùn)動(dòng)狀態(tài)自適應(yīng)更新采樣間隔,使得系統(tǒng)在目標(biāo)處于非機(jī)動(dòng)的穩(wěn)定狀態(tài)時(shí)用較大的采樣間隔,而在目標(biāo)處于機(jī)動(dòng)狀態(tài)時(shí)采用較小的采樣間隔,以確保目標(biāo)跟蹤不丟失。

[1]張光義,趙玉潔.相控陣?yán)走_(dá)技術(shù)[M].北京,電子工業(yè)出版社,2006.

[2]王峰,張洪才,潘泉.相控陣?yán)走_(dá)采樣周期自適應(yīng)策略研究[J].系統(tǒng)仿真學(xué)報(bào),2003,15(9):1230-1233.

WANG Feng,ZHANG Hong-cai,PAN Quan.A study on adaptive sampling period of phased array radar[J].Acta Simulata Systematica Sinica,2003,15(9):1230-1233.

[3]王鑫,胡昌華,暴飛虎.基于貝葉斯原理的粒子濾波算法[J].彈箭與制導(dǎo)學(xué)報(bào),2006,26(2):269-271.

WANG Xin,HU Chang-hua,BAO Fei-hu.Particle filtering method based on Bayesian theorem[J].Journal of Projectiles,Rockets,Missiles and Guidance,2006,26(2):269-271.

[4]Arulampalam S,Maskell S,Gordon N,et al.A tutorial on particle filters for online nonlinear/non-Gaussian Bayesian tracking[J].IEEE Trans.on Signal Processing,2002,50(2):174-188.

[5]楊德貴,張長城,鄭濤.粒子濾波算法的性能分析與改進(jìn)[J].電光與控制,2008,15(5):72-75.

YANG De-gui,ZHANG Chang-cheng,ZHENG Tao.Performance analysis and improvement of particle filter algorithm[J].Electronics Optics&Control,2008,15(5):72-75.

[6]張琪,胡昌華,喬玉坤.基于權(quán)值選擇的粒子濾波算法研究[J].控制與決策,2008,23(1):117-120.

ZHANG Qi,HU Chang-hua,QIAO Yu-kun.Particle filter algorithm based on weight selected[J].Control and Decision,2008,23(1):117-120.

[7]鄧小龍,謝劍英,楊煜普.基于交互式多模型的粒子濾波算法[J].系統(tǒng)仿真學(xué),2005,17(10):2360-2380.

DENG Xiao-long,XIE Jian-ying,YANG Yu-pu.Particle filter based on interacting multiple model[J].Acta Simulata Systematica Sinica,2005,17(10):2360-2380.

[8]Tichavsky P,Muravchik C H,Nehorai A.Posterior cramerrao bounds for discrete-time nonlinear filtering[J].IEEE Trans.on Signal Processing,1998,46(5):1386-1395.

[9]Hurtado M,Nehorai A.Adaptive path design of a moving radar[M].Sensor, Signal and Information Processing Workshop,2008.

[10]周文輝.相控陣?yán)走_(dá)及組網(wǎng)跟蹤系統(tǒng)資源管理技術(shù)研究[D].長沙:國防科技大學(xué),2004.

Adaptive sampling method for phased array radar based on interacting multiple model particle filter

YU Wei-hua1,2, ZHU Xiang3, ZHU Xiao-hua2
(1.Nantong Agricultural Vocational Technology College,Nantong226007,China;2.School of Electronic Engineering and Optoelectronic Technology,Nanjing University of Science and Technology,Nanjing210094,China;3.State Grid Electric Power Research Institute,Nanjing210003,China)

In order to effectively utilize the resources of radar and solve the nonlinear relationship between radar measurement and target motion state,an adaptive sample target tracking algorithm for phased array radar based on interacting multiple model particle filter (IMMPF) is proposed.This algorithm first predicts Posterior Cramer-Rao Bound (PCRB) matrix of the target state,then updates the sample interval adapted to changing target dynamics by comparing the trace of the predicted PCRB with a certain threshold.Performances of constant and adaptive data rates are compared.Simulation results demonstrate the effectiveness of this algorithm.

phased array radar; interacting multiple model; particle filter; PCRB; adaptive sampling

TN958

A

1674-6236(2012)05-0029-04

2011-12-01稿件編號(hào):201111148

郁衛(wèi)華(1966—),女,江蘇南通人,副教授。研究方向:雷達(dá)信號(hào)處理及應(yīng)用電子技術(shù)教學(xué)。

猜你喜歡
模型
一半模型
一種去中心化的域名服務(wù)本地化模型
適用于BDS-3 PPP的隨機(jī)模型
提煉模型 突破難點(diǎn)
函數(shù)模型及應(yīng)用
p150Glued在帕金森病模型中的表達(dá)及分布
函數(shù)模型及應(yīng)用
重要模型『一線三等角』
重尾非線性自回歸模型自加權(quán)M-估計(jì)的漸近分布
3D打印中的模型分割與打包
主站蜘蛛池模板: 尤物午夜福利视频| 萌白酱国产一区二区| 亚洲一级无毛片无码在线免费视频| www亚洲天堂| 亚洲福利视频一区二区| 国产第二十一页| 亚洲国产AV无码综合原创| 国产亚洲视频播放9000| 国产成人综合在线观看| 在线免费观看a视频| 国产精品无码翘臀在线看纯欲| 嫩草影院在线观看精品视频| 狠狠色狠狠色综合久久第一次| 91免费国产高清观看| 香蕉久人久人青草青草| 在线观看免费AV网| 欧美第一页在线| 性网站在线观看| 亚洲码在线中文在线观看| 欧美翘臀一区二区三区| 日韩国产一区二区三区无码| 午夜视频在线观看免费网站| 手机在线看片不卡中文字幕| 71pao成人国产永久免费视频| 亚洲天堂视频在线免费观看| 日韩在线影院| 老司机午夜精品网站在线观看| 免费在线a视频| 欧美性猛交xxxx乱大交极品| 尤物特级无码毛片免费| 精品成人一区二区| 精品色综合| 夜夜爽免费视频| 午夜天堂视频| 国产精品成人一区二区不卡| 最新国产你懂的在线网址| 国产综合另类小说色区色噜噜 | 国产XXXX做受性欧美88| 国产真实乱子伦精品视手机观看| h视频在线播放| 国产美女免费| 国产精品久久久久婷婷五月| 2021天堂在线亚洲精品专区| 国产激情无码一区二区三区免费| 亚洲色无码专线精品观看| 青草精品视频| 91在线丝袜| 老司机午夜精品网站在线观看| 久久国产高潮流白浆免费观看| 精品无码日韩国产不卡av| 欧美色视频日本| 一级毛片免费不卡在线视频| 色综合久久综合网| 国产午夜在线观看视频| 国产精品无码影视久久久久久久 | 亚洲美女视频一区| 亚洲欧美日韩综合二区三区| 老色鬼久久亚洲AV综合| 免费看美女毛片| 无码不卡的中文字幕视频| 丰满的熟女一区二区三区l| 国产精品视频猛进猛出| 毛片免费高清免费| 中文字幕精品一区二区三区视频| 99在线视频精品| 夜夜操天天摸| 国产高清无码第一十页在线观看| 色久综合在线| 精品五夜婷香蕉国产线看观看| 亚洲av日韩av制服丝袜| 91在线高清视频| 手机在线免费不卡一区二| 日韩精品一区二区深田咏美| 成人精品在线观看| 午夜小视频在线| 亚洲香蕉在线| 亚洲精品卡2卡3卡4卡5卡区| 久久国产黑丝袜视频| 国产日本欧美亚洲精品视| 日韩乱码免费一区二区三区| 18禁影院亚洲专区| 国产小视频a在线观看|