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

        ?

        粒子濾波算法改進及其應(yīng)用研究

        2015-05-30 10:48:04王永翔陳國初
        計算機時代 2015年8期
        關(guān)鍵詞:粒子濾波

        王永翔 陳國初

        摘 要: 粒子濾波算法是一種基于貝葉斯估計和非參數(shù)化蒙特卡羅模擬的新型算法。詳細介紹了粒子濾波算法的基本原理,闡述了遞推貝葉斯估計方法、序貫重要性采樣以及序貫重要性重采樣算法;分析了粒子濾波主要存在的問題及關(guān)鍵改進之處,然后介紹粒子濾波算法的應(yīng)用,并對算法未來發(fā)展進行了展望。

        關(guān)鍵詞: 貝葉斯估計; 序貫重要性采樣; 重采樣; 粒子濾波

        中圖分類號:TP301.6 文獻標志碼:A 文章編號:1006-8228(2015)08-01-04

        Research on development of particle filtering algorithm

        Wang Yongxiang, Chen Guochu

        (College of Electrical Engineering, Shanghai Dianji University, Shanghai 200240, China)

        Abstract: Particle filtering is a new algorithm based on Bayesian estimation and non-parametric Monte Carlo simulation. This paper introduces in detail the basic principle of particle filtering algorithm, describes the recursive Bayesian estimation method, sequential importance sampling and sequential importance resample algorithm, analyzes the main problems and the key improvements on the particle filter. Then, the application and the future development of particle filtering algorithm is discussed.

        Key words: Bayesian estimation; sequential importance sampling; resample; particle filter

        0 引言

        許多復(fù)雜的系統(tǒng)都存在著濾波處理、參數(shù)估計等問題,尤其在通信、自動控制、金融、彈道武器等領(lǐng)域濾波處理顯得更加重要。1960年,卡爾曼濾波器(KF)[1]被提出,其作為一種線性、無偏的最小方差狀態(tài)估計器得到很廣泛應(yīng)用。然而,由于卡爾曼濾波對系統(tǒng)要求比較苛刻,其只能處理一些線性并且噪聲服從高斯分布的線性系統(tǒng),對于非線性、非高斯噪聲的復(fù)雜系統(tǒng),其往往達不到濾波要求。盡管后來有些學(xué)者提出了擴展Kalman濾波器(EKF)[2]、無跡Kalman濾波器(UKF)[3]等,處理的系統(tǒng)也往往是高斯噪聲分布的弱非線性的,對于強非線性的系統(tǒng),其濾波精度及穩(wěn)定性等仍然達不到要求。

        20世紀50年代,Hammersley等人提出了一種序貫重要性采樣(SIS)蒙特卡羅方法[4],通過離散的隨機樣本來逼近概率分布。該算法主要是從一個重要性分布中抽取一個帶有權(quán)值的樣本集合來代表所研究的系統(tǒng)狀態(tài)的后驗概率密度函數(shù),并利用這些樣本及其權(quán)值對系統(tǒng)狀態(tài)進行估計。然而,由于算法本身存在樣本退化和計算復(fù)雜度高等問題,該算法很長一段時間沒有得到較大發(fā)展。直到1993年,Gordon等人提出了重采樣概念[5],將其引入到SIS中產(chǎn)生了基本粒子濾波算法—序貫重要性重采樣算法(SIR)。新算法克服了原算法中的樣本退化問題,同時計算機技術(shù)的快速發(fā)展也為SIR算法的計算提供了有力保證。

        本文首先通過狀態(tài)空間模型詳細介紹遞推貝葉斯估計,粒子濾波算法的基本原理,闡述序貫重要性采樣以及重要性重采樣算法,然后介紹粒子濾波常見問題及關(guān)鍵改進之處,最后介紹粒子濾波算法的應(yīng)用并對算法未來發(fā)展進行了展望。

        1 基本粒子濾波算法

        1.1 遞推貝葉斯估計方法

        粒子濾波(Particle Filter,PF)的思想是基于貝葉斯估計原理利用粒子集來表示概率的蒙特卡洛(Monte Carlo)模擬方法,它可以用在任何形式的狀態(tài)空間模型上。假設(shè)狀態(tài)空間模型如下:

        狀態(tài)模型:

        xk=f(xk-1,uk-1) ⑴

        觀測模型:

        yk=h(xk,vk) ⑵

        其中,f和h為已知的狀態(tài)函數(shù)和觀測函數(shù),xk是系統(tǒng)在k時刻的狀態(tài)量,uk和vk分別為零均值的系統(tǒng)噪聲和觀測噪聲,yk是k時刻狀態(tài)的觀測值,uk和vk相互獨立并獨立于過去和當前狀態(tài)。

        用統(tǒng)計描述法描述上式方程,系統(tǒng)模型對應(yīng)為狀態(tài)轉(zhuǎn)移概率密度p(xk|xk-1),觀測模型對應(yīng)為觀測似然概率密度p(yk|xk),并假定狀態(tài)的初始概率密度函數(shù)p(x0|y0)=p(x0)已知。設(shè)已知k-1時刻的概率密度p(xk-1|y1:k-1),且狀態(tài)xk服從一階馬爾可夫過程,則由第一步預(yù)測可得k時刻的先驗概率密度函數(shù):

        其中,系統(tǒng)的轉(zhuǎn)移概率密度:

        假定p(uk-1|xk-1)=p(uk-1),則:

        其中,δ(·)是狄拉克函數(shù),若δ(·)函數(shù)存在,xk-1和uk-1已知,那么狀態(tài)值xk可以通過⑴式得到。然后,若yk測量可知,則可由貝葉斯準則實現(xiàn)第二步的更新,得后驗濾波概率密度:

        其中,似然概率密度:

        通常是一個歸一化常數(shù)。

        由式⑶~⑹構(gòu)成了貝葉斯遞推估計過程,整個過程實現(xiàn)了由k-1時刻的p(xk-1|y1:k-1)到k時刻的p(xk|y1:k)預(yù)測更新,理論上解決了求解后驗概率密度的問題。但是由于式⑶的積分在實際中很難求解,并且整個過程中運用了許多假設(shè),在實際應(yīng)用中受到了很大的限制,所以對于非線性非高斯分布的系統(tǒng)貝葉斯估計法還不能給出很好的解決方法。

        1.2 序貫重要性采樣(SIS)

        針對貝葉斯估計無法解決的積分問題,1954年Hammersley等人提出了SIS算法。根據(jù)Monte Carlo方法從后驗概率密度中抽取N個獨立同分布樣本,用經(jīng)驗分布逼近狀態(tài)的概率密度函數(shù),使

        然而,在實際中p(x0:k|y1:k)往往是非標準的多變量形式,并且它很難寫成組合形式的解析分布函數(shù),從樣本p(x0:k|y1:k)中抽樣通常是非常困難的甚至沒法實現(xiàn)。所以,常常從被引入的一個容易實現(xiàn)采樣的重要分布函數(shù)q(x0:k|y1:k)中抽樣。為了便于采樣,首先對q(x0:k|y1:k)進行分解:

        當觀測值到來時,將新的數(shù)據(jù)值加入到舊樣本中得到新的樣本。這樣,從q(x0:k|y1:k)抽取N個獨立同分布樣本并加權(quán)求和來近似表示p(x0:k|y1:k),得:

        其中,為歸一化后的權(quán)值。

        假設(shè):

        得:

        結(jié)合式⑽可得權(quán)值的遞推公式為:

        為算法簡便,設(shè)q(xk|x0:k-1,y1:k)=q(xk|xk-1,y1:k),這樣重要性密度函數(shù)就只由xk-1和yk決定,每一步只需估計出p(xk|y1:k),權(quán)值遞推式變?yōu)椋?/p>

        相應(yīng)的近似后驗概率密度函數(shù):

        1.3 序貫重要性重采樣算法

        SIS算法解決了Bayesian估計無法解決的積分問題,以樣本均值代替了積分運算實現(xiàn)了遞推貝葉斯估計。但是,SIS算法卻存在著嚴重的樣本退化問題。當空間維數(shù)或時間增加時,經(jīng)過很多次迭代之后,較少量的粒子權(quán)重會變得很大,而絕大多數(shù)的粒子權(quán)重會變得很小,權(quán)值的更新大部分時間都浪費在這些權(quán)重小甚至不起作用的粒子上。

        1993年,Gordon等人提出的重采樣算法有效的解決了樣本退化問題。重采樣的根本目的是去除那些權(quán)值小的粒子,保留并復(fù)制那些權(quán)值較大的粒子[6],從式⒃中抽取相互獨立的權(quán)值相等的粒子,可以設(shè)。如圖1,形象地說明了重要性采樣與重采樣過程。

        圖1 重要性采樣與重采樣過程

        粒子退化程度一般采用有效粒子數(shù)Neff來衡量

        有效的粒子個數(shù)Neff越小,意味著粒子退化越厲害,當Neff小于一定的閾值NT時,然后才進行重采樣。

        基本粒子濾波算法就是在序貫重要性采樣(SIS)基礎(chǔ)上加入了重采樣步驟,具體步驟如下。

        第一步:初始化(k=0)

        從初始分布p(x0)中采樣N個樣本:{,i=1,2,…,N},權(quán)值:=1/N,令k=1。

        第二步:重要性采樣

        ⑴ 采樣得N個k時刻樣本:{,i=1,2,…,N}

        ⑵ 觀測值到來時,計算每個粒子權(quán)值,并由下式更新權(quán)值

        ⑶ 權(quán)值的歸一化

        第三步:判斷是否重采樣

        計算有效粒子數(shù)Neff=1/,若Neff?NT則執(zhí)行第四步;否則執(zhí)行第五步。

        第四步:重采樣

        依據(jù)權(quán)值大小采樣新粒子,使其滿足pr(=)=,設(shè)定=,=1/N。

        第五步:估計輸出

        后驗概率密度估計輸出:

        狀態(tài)估計輸出:

        第六步:令k=k+1,返回第二步。

        2 關(guān)鍵問題及改進之處

        2.1 粒子濾波存在的關(guān)鍵問題

        ⑴ 重要性密度函數(shù)選擇問題

        重要性密度函數(shù)實現(xiàn)了粒子濾波算法中的抽樣問題,其選擇往往是以最小化重要性權(quán)值為標準。實驗證明,為了有效地限制粒子退化,一般情況下是選擇重要性密度函數(shù)作為先驗概率密度函數(shù),即q(xk|,y1:k)=p(xk|,y1:k)。然而由于這種方法沒有考慮當前時刻的觀測值,使現(xiàn)有的模型處理當前狀態(tài)時會產(chǎn)生很大的誤差,不能將概率密度函數(shù)真實分布準確的表示出來。

        ⑵ 樣本退化問題

        樣本退化(或粒子退化)現(xiàn)象是粒子濾波技術(shù)中最嚴重的一個問題,選擇好的重要性密度函數(shù)和重采樣是解決這個問題的關(guān)鍵。重采樣算法的提出對樣本退化問題有了一定的抑制作用,但同時也產(chǎn)生了另外一個問題—粒子耗盡問題。由于在重采樣過程中,隨著采樣次數(shù)的增加,較少的權(quán)值大的粒子被復(fù)制保留,大多數(shù)權(quán)值小的粒子被遺棄,這樣導(dǎo)致粒子多樣性被逐漸削減。

        ⑶ 樣本貧乏問題

        樣本貧乏問題往往是由重采樣過程導(dǎo)致的,重采樣后的少數(shù)粒子不能真正的反映出后驗概率密度的特征。尤其是在模型參數(shù)不變或過程噪聲較小時,樣本貧乏問題顯得尤為突出,因為若重采樣后的粒子數(shù)僅剩幾個或變成一個點時,就會導(dǎo)致樣本枯竭而嚴重影響濾波性能。

        ⑷ 濾波實時性問題

        實時性也是粒子濾波中的一個重要問題。粒子濾波雖然解決了卡爾曼算法不能解決的非線性、非高斯等問題,但由于粒子濾波是源于Monte Carlo的思想隨著粒子數(shù)的不斷增加,算法的計算復(fù)雜度也逐漸變高,這就會嚴重影響算法的實時性。

        2.2 粒子濾波的改進之處

        ⑴ 基于重要性密度函數(shù)選擇的改進

        重要性密度函數(shù)的選擇,其根本目的是限制粒子的退化,提高算法的準確度。SIR粒子濾波器選擇先驗概率密度函數(shù)p(xk|)作為重要性密度函數(shù)使得算法簡潔易求,但由于其選擇沒有考慮最新來的觀測值,使得算法精度不是很高。文獻[7]將EKF方法結(jié)合到PF算法中構(gòu)造重要性密度函數(shù),提出了擴展Kalman粒子濾波器(EKF-PF),不僅考慮了先驗概率分布,也結(jié)合了當前的觀測值,使得重要性密度函數(shù)越來越接近真實的后驗密度函數(shù)分布。Merwe等人[8]將估計精度優(yōu)于EKF的UKF應(yīng)用于粒子濾波算法中,提出了無味粒子濾波器(UPF),獲得更好的采樣,使濾波算法性能得到進一步的改善。此外,對重要性密度函數(shù)的改進,還提出了高斯-厄米特粒子濾波、高斯和粒子濾波等等。

        ⑵ 基于重采樣技術(shù)的改進

        基于重采樣改進PF算法中最典型的是1999年P(guān)itt等人提出的輔助粒子濾波算法,通過引入一個輔助變量來獲得一個考慮當前狀態(tài)觀測值的重要性概率密度函數(shù)。正則化粒子濾波器[9]也是在重采樣階段通過核密度函數(shù)進行抽樣的改進算法。此外,除了對重采樣本身改進之外,也有在重采樣之前、重采樣之后和免去重采樣等地方改進。文獻[10]提出的進化粒子濾波器,將進化規(guī)則融合到PF算法中,即在重采樣前先對未更新粒子進行變異處理,通過競爭選擇出似然密度函數(shù)大的粒子。文獻[11]通過在重采樣之后對粒子進行轉(zhuǎn)移,然后提出了改進的馬爾可夫蒙特卡羅粒子濾波算法,其通過Markov轉(zhuǎn)換核函數(shù)使p(x0:k|y1:k)產(chǎn)生新的粒子序列來逼近目標的真實后驗概率分布。文獻[12]直接免去重采樣步驟,將高斯分布近似概率密度p(xk|y1:k)和p(xk|y1:k-1),提出了高斯粒子濾波器消除了重采樣帶來的樣本匱乏問題。

        ⑶ 其他改進方法

        基于重要性密度函數(shù)和重采樣的改進是粒子濾波算法中最常見的兩種改進方法,這兩種方法也是針對粒子濾波的有效性進行改進。除此之外,也有對粒子的實時性進行改進的算法,如改變粒子濾波的樣本數(shù)量或改變系統(tǒng)狀態(tài)的維數(shù)等。文獻[13]提出的邊緣化粒子濾波和文獻[14]通過KLD采樣樣本分布來自適應(yīng)改變粒子的數(shù)量方法,都使算法的實時性得到了很大提高。

        3 粒子濾波的應(yīng)用

        隨著科學(xué)技術(shù)的快速發(fā)展以及各領(lǐng)域各系統(tǒng)的復(fù)雜度的提高,濾波處理技術(shù)也不斷得到發(fā)展。粒子濾波以其可以處理非線性非高斯等復(fù)雜系統(tǒng)的優(yōu)勢廣泛的被應(yīng)用于各個領(lǐng)域,如目標跟蹤、系統(tǒng)識別與參數(shù)估計、金融數(shù)據(jù)分析、故障診斷等等。

        3.1 目標跟蹤

        衛(wèi)星導(dǎo)航系統(tǒng)、車輛、雷達跟蹤等領(lǐng)域都涉及到目標定位與跟蹤問題,粒子濾波技術(shù)在這些復(fù)雜的非線性系統(tǒng)的處理中得到廣泛的應(yīng)用。傳統(tǒng)的粒子濾波在目標跟蹤中,其模型的樣本數(shù)控制往往是比較困難的,文獻[15]將狀態(tài)估計和系統(tǒng)估計分解開來計算濾波的粒子,使結(jié)果精度得到很大提高。文獻[16]將粒子濾波技術(shù)結(jié)合徑向神經(jīng)網(wǎng)絡(luò)最佳逼近特性應(yīng)用于直升機運動跟蹤中,在跟蹤精度與實時性中都顯示出了較好地效果,很好解決了卡爾曼等傳統(tǒng)濾波器無法解決的非線性難題。

        3.2 參數(shù)估計與系統(tǒng)識別

        現(xiàn)代大多數(shù)復(fù)雜的非線性系統(tǒng)都存在著參數(shù)估計問題,由于系統(tǒng)的數(shù)據(jù)是隨機的或未知的,通過粒子濾波器對系統(tǒng)離散變量進行模型建立來估計系統(tǒng)的參數(shù)是一個有效地處理方法。由于粒子濾波技術(shù)可以處理線性的和非線性系統(tǒng),所以可將系統(tǒng)分為線性和非線性兩部分來處理,不但擴展了粒子濾波在單一的非線性系統(tǒng)中的應(yīng)用,而且也省去采用其他濾波器處理線性部分的麻煩。

        3.3 金融數(shù)據(jù)分析

        模糊隨機系統(tǒng)的估值問題是金融領(lǐng)域數(shù)據(jù)分析中常見的問題,金融數(shù)據(jù)是隨機波動不確定的,粒子濾波在處理這些數(shù)據(jù)估計時起到了一定的作用。文獻[17]就將粒子濾波技術(shù)應(yīng)用在股票交易數(shù)據(jù)預(yù)測中,根據(jù)交互式的粒子濾波原理分別從跟蹤和濾波角度對數(shù)據(jù)進行預(yù)測,達到很好的預(yù)測效果。

        3.4 故障診斷

        故障診斷是各類復(fù)雜系統(tǒng)中的一個重要研究分支,對于線性系統(tǒng)的故障診斷研究目前已經(jīng)取得了較好的成果,而對于非線性和強非線性系統(tǒng)的故障診斷仍然很不成熟。粒子濾波技術(shù)的較快發(fā)展對這方面的研究提供了解決方法。如一些模型尚不完備的復(fù)雜動態(tài)系統(tǒng)的故障診斷處理等,文獻[18]提出了一種有效的PF算法,通過把粒子集的規(guī)格化因子和最大后驗概率估計狀態(tài)的信度兩個值提取出來,并在此基礎(chǔ)上設(shè)出檢測未知故障模式閾值邏輯,達到了很好地預(yù)報故障效果。

        除了以上介紹的領(lǐng)域外,粒子濾波在圖像處理、計算機視覺、生物學(xué)、海洋學(xué)等領(lǐng)域都有相關(guān)性的應(yīng)用。

        4 總結(jié)與展望

        粒子濾波技術(shù)在處理非線性系統(tǒng)時具有獨特的優(yōu)勢,近些年其被研究的趨勢也在逐漸上升。在詳細闡述粒子濾波基本原理后,對算法中常見的問題以及一些改進的方法作了分析介紹,并舉例說明了其在各種領(lǐng)域中的廣泛應(yīng)用。但作為一種新的濾波處理方法,粒子濾波技術(shù)仍處于不斷的發(fā)展中,理論方法和實際應(yīng)用還有待進一步深入研究。

        ⑴ 首先,粒子濾波算法本身還有待改進,如算法的收斂性、重要性密度函數(shù)選擇、重采樣部分,以及基于系統(tǒng)觀測模型和狀態(tài)模型準確建模等。此外,怎樣減少計算量、削減粒子退化,以及進一步提高算法的準確性和實時性也仍然是未來的研究方向。

        ⑵ 其次,粒子濾波算法是基于Monte Carlo模擬的實現(xiàn)遞推貝葉斯估計的數(shù)學(xué)算法,可以將其與其他智能算法相結(jié)合來對研究目標進行處理。雖然已經(jīng)有不少學(xué)者根據(jù)其他智能算法對粒子濾波某個部分進行優(yōu)化,但在組合算法的性能上還有待進一步提高。怎樣增強粒子濾波等組合算法的穩(wěn)定性和實用性也是值得更深一步研究。

        ⑶ 最后,盡管粒子濾波技術(shù)已應(yīng)用于某些領(lǐng)域,但就其目前的發(fā)展還不是太成熟,僅僅部分的應(yīng)用于某些領(lǐng)域。并且如今大多數(shù)學(xué)者對粒子濾波的實際研究往往是基于仿真模擬性的研究,怎樣將粒子濾波技術(shù)擴展到更多的應(yīng)用新領(lǐng)域還有待于進一步探索。

        參考文獻:

        [1] Kalman R E. A new approach to linear filtering and prediction

        problem. Transactions of the ASME-Journal of Basic Engineering,1960.82:35-45

        [2] JAZW INSKIAH. Stochastic Processes and Filter ing Theory[M].

        New York: Academic Press,1970.

        [3] S. J. Julier, J. K. Uhimann, H. F. Durrant-Whyte. A new method

        for nonlinear transformation of means a and estimators, IEEE T ans Automatic Control,2000.45(3):477-482

        [4] Hammersley J M, Morton K W. Poor man's Monte Carlo[J].Journal

        of the Royal Statistical Society B,1954.16(1):23-38

        [5] Gordon N J, Salmond D J, Smith A F M. Novel approach to

        nonlinear/non-Gaussian Bayesian state estimation[J]. IEE proceedings F: Radar and Signal Processing,1993.140(2):107-113

        [6] 趙梅,張三同等.輔助粒子濾波算法及仿真舉例[J].北京交通大學(xué)學(xué)

        報,2006.30(2):24-28

        [7] Ristic B, Arulampalam S, Gordon N. Beyond the Kalman filter:

        Particle filters for tracking applications[M].Boston: Artech House, 2004.Boston: Artech House,2004.

        [8] Merwe R V D, Doucet A, Freitas, et al. The unscented particle filter[J].

        Advances in Neural Information Processing Systems,2001.96(390):584-590

        [9] Musso C, Oudjane N, Legland F. Improving regularized particle

        filter[M].Doucet A, de Freitas J F G,Gordon N J.Sequential Monte Carlo Methods in Practical, New York:Springer-Verlag,2001:247-272

        [10] 莫以為, 蕭德云. 進化粒子濾波算法及其應(yīng)用[J] .控制理論與應(yīng)用,

        2005.22(2):269-272

        [11] SPALL J C. Estimation via Markov Chain Monte Carlo[J].IEEE

        Control Systems Magazine,2003.23(2):34-45

        [12] BOLIC M, DJURIC P M, HONG SANGJIN. Resampling

        Algorithms and Architectures for Distributed Particle Filters[J]. IEEE Trans Signal Processing,2005.53(7):2442-2450

        [13] Schon T, Gustafsson F, Nordlund P J. Marginalized particle

        filters for mixed linear/nonlinear state-space models[J]. IEEE Transactions on: Signal Processing,2005.53(7):2279-2289

        [14] Fox D. Adapting the sample size in particle filters through

        KLD-sample[J]. The International Journal of Robotics Research,2003.22(12):985-1003

        [15] 鑒福升,徐躍民等.改進的多模型粒子濾波機動目標跟蹤算法[J].控

        制理論與應(yīng)用,2010.27(8):1012-1016

        [16] 李春鑫,王孝通等.改進粒子濾波算法在目標跟蹤中的應(yīng)用[J].控制

        工程,2009.16(5):575-577

        [17] Papavasiliou A. Adaptive particle filters with applications[R].

        Princeton: University of Princeton,2002.

        [18] 段琢華,蔡自興,于金霞.不完備多模型混合系統(tǒng)故障診斷的粒子濾

        波算法[J].自動化學(xué)報,2008.34(5):581-587

        猜你喜歡
        粒子濾波
        基于粒子濾波算法的視頻目標跟蹤研究
        數(shù)字濾波技術(shù)的現(xiàn)狀與發(fā)展
        可見光與紅外CPCT彩色圖像融合研究與實現(xiàn)
        OCV處于平臺期的汽車鋰電池SOC估算的研究
        粒子濾波法在移動物體定位中的應(yīng)用
        基于智能粒子濾波的多傳感器信息融合算法
        基于重采樣粒子濾波的目標跟蹤算法研究
        航空兵器(2016年5期)2016-12-10 17:12:24
        復(fù)雜場景中的視覺跟蹤研究
        一種針對海面弱小目標的檢測方案
        交互多模型的Rao—Blackwellized粒子濾波算法在多目標跟蹤中的應(yīng)用
        航空兵器(2014年4期)2015-01-15 02:31:43
        国产成人AV无码精品无毒| 日韩产的人妻av在线网| 91国内偷拍一区二区三区| 精品国产av一区二区三区| 国产精品无码制服丝袜| 精品国产一区二区三区av天堂| 欧美猛少妇色xxxxx猛交| 国产精品白浆一区二小说| 欧美人与动牲交片免费| 激情五月开心五月啪啪| 国产精品妇女一区二区三区| 亚洲av午夜国产精品无码中文字| 亚洲av电影天堂男人的天堂| 精品三级久久久久久久电影| 国产精品天天看大片特色视频| av网站韩日在线观看免费| 亚洲精品国产电影| 青青草97国产精品免费观看| 亚洲网站地址一地址二| 亚洲欧洲日产国码久在线观看| 99视频偷拍视频一区二区三区| 亚洲人成网站18禁止| 国产莉萝无码av在线播放| 91老司机精品视频| 亚洲av午夜成人片精品| 男女后入式在线观看视频| 日本妇人成熟免费2020| 18禁无遮拦无码国产在线播放| 日本VA欧美VA精品发布| 欧美日本视频一区| 综合图区亚洲偷自拍熟女| 亚洲一区二区国产一区| 熟女无套高潮内谢吼叫免费| 99国产精品无码| 国产欧美日韩在线观看 | 久久精品视频在线看99| 国产欧美日韩网站| 亚洲天堂一二三四区在线| 日本男人精品一区二区| 无码不卡av东京热毛片| 久久精品国产亚洲av瑜伽|