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

改進(jìn)SOA算法在焊縫圖像分割中的應(yīng)用

2016-01-25 06:52:38余勝威丁建明曹中清

余勝威,丁建明,曹中清

(1.西南交通大學(xué) 機(jī)械工程學(xué)院, 四川 成都 610031;

2.西南交通大學(xué) 牽引動(dòng)力國家重點(diǎn)實(shí)驗(yàn)室,四川 成都 610031)

?

改進(jìn)SOA算法在焊縫圖像分割中的應(yīng)用

余勝威1,丁建明2,曹中清1

(1.西南交通大學(xué) 機(jī)械工程學(xué)院, 四川 成都 610031;

2.西南交通大學(xué) 牽引動(dòng)力國家重點(diǎn)實(shí)驗(yàn)室,四川 成都 610031)

摘要:提出一種基于改進(jìn)的自適應(yīng)人群搜索算法(ISOA)的焊縫圖像分割方法,該算法通過對(duì)基本人群搜索算法的高斯隸屬度函數(shù)以及隸屬度參數(shù)計(jì)算的改進(jìn),自適應(yīng)確定搜索步長,實(shí)現(xiàn)對(duì)n尺度圖像進(jìn)行n-1個(gè)閾值尋優(yōu)計(jì)算。實(shí)驗(yàn)結(jié)果表明:對(duì)比PSO,GA和SA算法,ISOA算法具有收斂速度快、穩(wěn)定性強(qiáng)、精度高、全局尋優(yōu)等特點(diǎn),有效地克服了傳統(tǒng)SOA算法易陷入局部最優(yōu)和收斂速度慢等缺陷,可滿足實(shí)際工程需求。

關(guān)鍵詞:多尺度分割;ISOA;類方差;PSO;GA;SA

一種閾值尋優(yōu)的最感性的方法是采用窮舉法,基于Otsu的窮舉法較簡潔,但是計(jì)算量大[8],采用窮舉法對(duì)本文“car測(cè)試圖像”進(jìn)行閾值尋優(yōu),當(dāng)閾值個(gè)數(shù)為2時(shí),CPU耗時(shí)15.41 s,當(dāng)閾值個(gè)數(shù)為3時(shí),CPU耗時(shí)為4 095.01 s,計(jì)算量成千倍增加,采用窮舉法對(duì)n-1個(gè)閾值進(jìn)行計(jì)算,有n(L-n+1)n-1種[7]可能,因此從算法執(zhí)行效率上考慮,窮舉法不適合用于閾值優(yōu)化分析。然而n-1閾值尋優(yōu)可轉(zhuǎn)換為一個(gè)多尺度的函數(shù)優(yōu)化問題。

本文采用生物智能算法進(jìn)行圖像分割。生物智能算法能夠解決傳統(tǒng)算法不能解決的問題,如當(dāng)尋優(yōu)目標(biāo)函數(shù)是離散的、不可微或者含有一些非線性相關(guān)參量時(shí),傳統(tǒng)算法根本無法解決。粒子群算法在搜索域內(nèi),粒子以一定的速度更新當(dāng)前最優(yōu)粒子和最優(yōu)種群。PSO簡單易實(shí)現(xiàn),然而易出現(xiàn)早熟等現(xiàn)象,以致不能全局尋優(yōu)[10]。遺傳算法是一種采用適者生存,優(yōu)勝劣汰的隨機(jī)搜索方法,GA直接對(duì)結(jié)構(gòu)對(duì)象進(jìn)行操作,不存在求導(dǎo)和函數(shù)連續(xù)性的限定,但GA存在收斂性差、耗時(shí)長等缺點(diǎn)[11],給圖像分割閾值全局尋優(yōu)帶來困難。模擬退火算法(SA)是一種啟發(fā)式隨機(jī)搜索算法,通過模擬固體物質(zhì)的退火過程,克服了窮舉法計(jì)算量過大的缺點(diǎn),然而計(jì)算量大,Hammouche等研究表明[12],PSO圖像分割在CPU執(zhí)行時(shí)間,結(jié)果穩(wěn)定性及精確度等方面都優(yōu)于GA,DE,ACO,SA和TS。人群搜索算法(SOA)[13],通過隨機(jī)挑選一定數(shù)量的志愿者,基于一定的啟發(fā)式信息、搜索規(guī)則和策略等,模擬人群搜索/覓食行為。目前該算法已用于函數(shù)優(yōu)化、IIR數(shù)字濾波器設(shè)計(jì)、電力系統(tǒng)無功優(yōu)化等領(lǐng)域[14],與DE,APSO,CFPSO和CLPSO算法比較,SOA算法具有較好的全局尋優(yōu)能力、較快的收斂速度和算法魯棒性[15]。本文提出一種改進(jìn)的自適應(yīng)SOA算法(ISOA),并將該算法應(yīng)用于圖像分割中,分割效果較滿意。本文主要研究性能相對(duì)較好的SOA,PSO,GA和SA算法,同時(shí)也為了驗(yàn)證ISOA能用于多尺度圖像分割。

1N尺度閾值分割模型

多尺度圖像分割是一種有效的圖像分割方法,而n尺度閾值尋優(yōu)的魯棒性對(duì)于圖像分割則是一個(gè)難點(diǎn)。本節(jié)提供一種更精確的解決方法,下面介紹一些基本的表達(dá)式[9]。

給定一個(gè)最大灰度值L(本文直接設(shè)為256),對(duì)于一副灰度圖像而言,閾值取值為{0,1,2,…,L-1},定義

其中i表示具體的灰度值,0≤i≤L-1,N表示圖像的像素個(gè)數(shù),hi表示i灰度值的個(gè)數(shù),pi表示概率值,則灰度的平均值為:

二維閾值尺度問題能夠拓展到n維閾值尺度,考慮n-1個(gè)不同閾值尺度tj,j=1,2,3,…,n-1,具體的表達(dá)式如下:

(1)

其中:x表示H×W圖像的寬度坐標(biāo)值;y表示H×W圖像的高度坐標(biāo)值,其像素二維函數(shù)為f(x,y),待分割圖像將被分成D1,D2,…,Dn類,組成F(x,y)。

最簡單并且最快速的閾值尋優(yōu)方法是,直接找出不同類之間最大方差對(duì)應(yīng)的灰度值,一般定義如下:

(2)

其中:j代表一個(gè)具體的類;wj表示第j類的概率;μj為第j類的均值;具體的D1,D2,…,Dn類的wj和μj表示如下:

(3)

(4)

從而一個(gè)n尺度閾值問題降維為一個(gè)函數(shù)尋優(yōu)問題,尋找閾值tj,使不同類間方差目標(biāo)函數(shù)最大,目標(biāo)函數(shù)如下:

(5)

相應(yīng)的算法核心代碼如下:

fitR(j)=sum(probR(1:xR(j,1)))*(sum((1:xR(j,1)).*probR(1:xR(j,1))/sum(probR(1:xR(j,1)))) - sum((1:Lmax).*probR(1:Lmax)) )^2;

forjlevel=2:level-1

fitR(j)=fitR(j)+sum(probR(xR(j,jlevel-1)+1:xR(j,jlevel)))*(sum((xR(j,jlevel-1)+1:xR(j,jlevel)).*probR(xR(j,jlevel-1)+1:xR(j,jlevel))/sum(probR(xR(j,jlevel-1)+1:xR(j,jlevel))))- sum((1:Lmax).*probR(1:Lmax)))^2;

end

fitR(j)=fitR(j)+sum(probR(xR(j,level-1)+1:Lmax))*(sum((xR(j,level-1)+1:Lmax).*probR(xR(j,level-1)+1:Lmax)/sum(probR(xR(j,level-1)+1:Lmax)))-sum((1:Lmax).*probR(1:Lmax)))^2;

fitBestR(j)=fitR(j);

%xR為種群個(gè)體,level為分割尺度,Lmax=256;

% probR為Pi概率值。

隨著灰度級(jí)的增加,目標(biāo)函數(shù)尋優(yōu)將變得更復(fù)雜,計(jì)算量大大增加;因此,選取哪種算法進(jìn)行優(yōu)化模型分析,成為我們需要考慮的問題。

最近生物智能算法被廣泛應(yīng)用于目標(biāo)優(yōu)化分析,它是一種高效的智能算法。

2算法分析

隨著智能算法的日新月異,算法的改進(jìn)顯得越來越重要,其計(jì)算速度、穩(wěn)定性、全局尋優(yōu)能力是當(dāng)前研究的熱點(diǎn)問題,改進(jìn)的自適應(yīng)人群搜索算法(ISOA)立足簡化SOA算法[15],通過對(duì)隸屬度函數(shù)以及慣性權(quán)值的改進(jìn),簡化了算法結(jié)構(gòu),縮短了算法執(zhí)行時(shí)間,提高了算法穩(wěn)定性以及全局尋優(yōu)能力等。

2.1 搜索步長的確定

基本SOA算法,采用個(gè)體在當(dāng)前種群中的相對(duì)位置排序,進(jìn)行隸屬度計(jì)算,增大了計(jì)算復(fù)雜度;而ISOA不確定推理行為的隸屬度隨迭代次數(shù)線性遞增:

(6)

uij=ui+rand·(1-ui),j=1,…,D

(7)

其中,ui為第i代個(gè)體對(duì)應(yīng)目標(biāo)函數(shù)值φ的隸屬度,本文采用線性隸屬函數(shù),在最佳位置有最大隸屬度值umax=1.0,最差位置有最小隸屬度umin=0.011 1;iter和itermax分別為當(dāng)前迭代次數(shù)和最大迭代次數(shù);uij為j維搜索空間第i代個(gè)體對(duì)應(yīng)目標(biāo)函數(shù)值φ的隸屬度;D為搜索空間維數(shù)。

由不確定推理可得步長:

(8)

其中,αij為j維搜索空間的搜索步長;δij為高斯隸屬函數(shù)參數(shù),計(jì)算如下:

(9)

(10)

式中:pg(t)為種群全局最優(yōu)值;rands(·)為[-1, 1]隨機(jī)變量;ω是慣性權(quán)值,Wmin和Wmax分別為相應(yīng)的最小最大權(quán)值。

2.2 搜索方向的確定

采用第i代個(gè)體的利己方向、利他方向以及預(yù)動(dòng)方向隨機(jī)加權(quán)平均,確定搜索方向如式(11):

(11)

式中:pi(t)為種群個(gè)體最優(yōu)值;xi(t-1)為上一代個(gè)體;sign(·)為符號(hào)函數(shù);?1和?2為[0,1]隨機(jī)數(shù)。

2.3 個(gè)體位置的更新

確定搜索方向和步長后,進(jìn)行位置更新:

Δxij(t+1)=αij(t)dij(t)

(12)

xij(t+1)=xij(t)+Δxij(t+1)

(13)

相應(yīng)的算法核心代碼如下:

% maxgen=itermax,t=iter迭代次數(shù)

% 慣性權(quán)值,Wmax=0.9,Wmin=0.1

W=Wmax-t*(Wmax-Wmin)/maxgen;

% 隸屬度值Umax=1.0,Umin=0.011 1

u=Umin+t*(Umax-Umin)/maxgen;

U=u+(1-u)*rand; % 隸屬度

T= sqrt(-log(U)); % 搜索步長

% 確定利己方向,gbest為pi(t)

Diego(i,:) = (gbest(i,:) - pop(i,:));

% 確定利他方向,zbest為pg(t)

Dialt(i,:) = (zbest - pop(i,:));

% 確定預(yù)動(dòng)方向,pop_1為上一代個(gè)體

Dipro(i,:) = (pop(i,:) - pop_1(i,:));

% 確定經(jīng)驗(yàn)梯度方向

Di(i,:) = sign(W* Dipro(i,:) +rand*Diego(i,:) +rand*Dialt(i,:));

% 確定高斯函數(shù)的參數(shù),N_PAR=D

C(i,:) =W*abs(zbest.*rands(1,N_PAR));

Bc(i,:) =C(i,:)*T;% 確定搜索步長

pop(i,:)=pop(i,:)+Di(i,:).*Bc(i,:); % 更新位置

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

電腦配置為AMD A8-4500M(1.9 GHz),2 GB可用內(nèi)存,軟件MATLAB(2014a)。圖像大小均為200×400,圖1為不同圖像及其直方圖。

圖1 不同的圖形及對(duì)應(yīng)的直方圖Fig.1 Histogram of different image

初始化粒子速度均設(shè)定為0,粒子的位置在可行域內(nèi)隨機(jī)產(chǎn)生,搜索空間取決于圖像灰度值的范圍,也就是說如果圖像是8位色圖,那么粒子的尋值空間為0~255。ISOA,PSO,GA和SA都是參數(shù)化的算法,為了增強(qiáng)對(duì)比效果,選取迭代次數(shù)為50,種群數(shù)為30(等效于SA迭代次數(shù)(鏈長)為30)。其中,PSO算法粒子速度滿足[vmin,vmax]=[-5, 5],GA算法交叉概率為pc=0.8,變異概率為pm=0.2,SA初始溫度T=1 000,降溫速率q=0.5,本文中閾值尺度2≤n≤5,選取適應(yīng)度估計(jì)值、CPU計(jì)算時(shí)間以及算法穩(wěn)定性3個(gè)指標(biāo)進(jìn)行算法對(duì)比研究。

3.1 適應(yīng)度值估計(jì)

因?yàn)樗惴ㄊ请S機(jī)產(chǎn)生個(gè)體的,對(duì)于每一個(gè)尺度n,分別采用ISOA,PSO,GA和SA算法運(yùn)行20次,統(tǒng)計(jì)各算法得到的適應(yīng)度均值和標(biāo)準(zhǔn)偏差值,如表1所示。

如表1所示,閾值個(gè)數(shù)分別為2,3,4和5,尺度對(duì)應(yīng)為3,4,5和6。盡管對(duì)于每一尺度,各算法計(jì)算的最優(yōu)適應(yīng)度值是相近的,然而隨著閾值個(gè)數(shù)的增加,相對(duì)標(biāo)準(zhǔn)偏差值也增大。相比較,ISOA算法穩(wěn)定且不易陷入局部最優(yōu)。

表1 不同圖像分割算法得到的適應(yīng)度平均值±STD

3.2 CPU處理時(shí)間

衡量算法有效性的一個(gè)重要指標(biāo),就是CPU計(jì)算時(shí)間[9],CPU計(jì)算時(shí)間越小,算法越能夠被考慮來應(yīng)用于實(shí)際控制器,算法執(zhí)行時(shí)間過長,則算法本身有待改進(jìn)。表2為不同算法的處理時(shí)間對(duì)比。

表2 不同圖像分割算法的CPU處理時(shí)間

如表2所示,ISOA算法表現(xiàn)更好的性能,用時(shí)最短,其次是PSO,再次是GA,最后是SA算法。Hammouche等研究表明[12],PSO圖像分割是在CPU執(zhí)行時(shí)間,結(jié)果穩(wěn)定性、結(jié)果精確度等方面優(yōu)于GA和SA。

為了直觀的觀察ISOA算法的分割結(jié)果,圖2為不同尺度下的ISOA算法分割結(jié)果。

圖2 圖像分割結(jié)果(從左到右表示閾值個(gè)數(shù)為2,3,4,5)Fig.2 Image segmentation results(n=2,3,4,5 for image from left to right)

如圖2所示,隨著閾值個(gè)數(shù)增加,圖像分割更加細(xì)致;當(dāng)閾值個(gè)數(shù)為5時(shí),即將圖像分為6部分,圖像顯得更加平滑。

3.3 算法全局尋優(yōu)分析

考慮到生物智能算法都具有隨機(jī)性特點(diǎn),即每一次運(yùn)行結(jié)果是不完全相同的,因此有必要進(jìn)行算法的全局尋優(yōu)能力分析。

為了定義算法的全局尋優(yōu)能力,本文采用標(biāo)準(zhǔn)偏差來衡量,如式(14)所示。

(14)

式中:STD是標(biāo)準(zhǔn)偏方差;σi是第i次算法運(yùn)行得到的適應(yīng)度值;μ是σ的平均值;N是算法執(zhí)行的次數(shù),在此取N=20。STD值越小,則表示算法的穩(wěn)定性越強(qiáng)。閾值取值分別為2,3,4,5,針對(duì)每一個(gè)閾值,算法均運(yùn)行N次,統(tǒng)計(jì)數(shù)據(jù)得如表1所示結(jié)果。

為了提高表1的可讀性,圖3表示采用不同算法,對(duì)不同尺度下的不同圖像的適應(yīng)度值標(biāo)準(zhǔn)偏差計(jì)算統(tǒng)計(jì)圖。

如圖3所示,GA是最不穩(wěn)定的,ISOA是相對(duì)最穩(wěn)定的,綜合各項(xiàng)指標(biāo),ISOA是一種性能最佳的圖像分割算法。

4結(jié)論

1)提出了一種ISOA算法的焊縫圖像分割算法,ISOA算法通過改進(jìn)高斯隸屬函數(shù)以及自適應(yīng)慣性權(quán)值,達(dá)到性能改進(jìn)的目的,實(shí)驗(yàn)結(jié)果表明,該算法具有快速收斂、全局尋優(yōu)能力等特點(diǎn)。

圖3 不同的圖像不同的STD值(從左到右表示2,3,4,5個(gè)閾值)Fig.3 STD values for image segmentation results(n=2,3,4,5 for image from left to right)

2)將ISOA算法成功應(yīng)用于焊縫圖像分割,實(shí)現(xiàn)了Otsu多級(jí)閾值下的焊縫圖像分割,對(duì)比于PSO,Ga和SA算法,基于適應(yīng)度值,STD和CPU計(jì)算時(shí)間性能指標(biāo),試驗(yàn)結(jié)果表明,采用ISOA算法性能更好,不易陷入局部最優(yōu),且隨著閾值個(gè)數(shù)的增加,分割效果越好。

3)一方面,焊縫圖像分割的自動(dòng)化能夠很好地實(shí)現(xiàn)焊縫的定位,且快速檢測(cè)表面有明顯缺陷的焊縫,另一方面,能夠降低人工目測(cè)檢測(cè)周期長,操作困難等缺陷。由于焊縫圖像受外界干擾影響較大,就需要保證算法的適應(yīng)性,ISOA算法根據(jù)用戶設(shè)定的分類數(shù),自動(dòng)進(jìn)行圖像的閾值最大化分割,并且分割結(jié)果能夠很好的自動(dòng)識(shí)別焊縫,因此基于ISOA算法在T型焊縫圖像分割的成功應(yīng)用,具有廣泛的應(yīng)用前景,例如PID參數(shù)整定、遙感圖像分割、視頻分析等領(lǐng)域。

參考文獻(xiàn):

[1] 李娜,李元香.基于自適應(yīng)粒子群算法和數(shù)據(jù)場的圖像二維閾值分割[J].計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào),2014,24(5):628-635.

LI Na, LI Yuanxiang.Image segmentation with two-dimension threshold based on adaptive particle swarm optimization and data field[J].Journal of Computer-Aided Design & Computer Graphics, 2014, 24(5):628-635.

[2] 余勝威,丁建明.轉(zhuǎn)向架構(gòu)架焊縫表面質(zhì)量檢測(cè)研究[J].鐵道科學(xué)與工程學(xué)報(bào),2015,12(2):402-406.

YU Shengwei, DING Jianming.Quality detection reaearch on bogie frame welding surface[J].Journal of Railway Science and Engineering, 2015,12 (2):402-406.

[3] 梁硼.X射線焊縫圖像缺陷自動(dòng)提取與識(shí)別技術(shù)研究[D].南京:南京航空航天大學(xué),2012.

LIANG Peng.Research on technology of automatic extraction and identification for weld defects in X-Ray image[D].Nanjing: Nanjing University of Aeronautics and Astronautics, 2012.

[4] 余剛.鋼制對(duì)接焊縫缺陷超聲相控陣檢測(cè)圖像特征與識(shí)別[D].南昌:南昌航空大學(xué),2012.

YU Gang.Defect image characterization and recognition of steel butt welds based on ultrasonic phased array[D].Nanchang:Nanchang Hangkong University, 2012.

[5] ZHANG Jing,JIN Yanxia.Genetic algorithm for weld image edge detection[J].IEEE,2010:511- 513.

[6] Xianghua Hou,Honghai Liu.Welding image edge detection and identification research based on canny operator[J].IEEE,2012:250-253.

[7] Brink A D.Minimum spatial entropy threshold selection[J].IEE Proceedings on Vision Image and Signal Processing, 1995, 142(3):128-132.

[8] Kulkarni R V, Venayagamoorthy G K.Bio-inspired algorithms for autonomous deployment and localization of sensor nodes[J].IEEE Transactions, 2010, 40(6):663-675.

[9] Pedram Ghamisi, Micael S Couceiro, Jón Atli Benediktsson, et al.An efficient method for segmentation of images based on fractional calculus and natural selection[J].Expert Systems with Applications, 2012, 39(16): 12407-12417.

[10] 胡衛(wèi)東,曹文貴.基于粒子群BP網(wǎng)絡(luò)混合算法的邊坡穩(wěn)定性評(píng)價(jià)[J].鐵道科學(xué)與工程學(xué)報(bào),2015,2(1):66-71.

HU Weidong, CAO Wengui.Slope stablitity evaluation based on hybrid algorithm of particle swarm optimization and BP neural network[J].Journal of Railway Science and Engineering, 2015,2(1):66-71.

[11] 鄧連波,史峰,莫輝輝.物流配送車輛路徑問題多代競爭遺傳算法[J].鐵道科學(xué)與工程學(xué)報(bào),2015,2(5):75-79.

DENG Lianbo, SHI Feng, MO Huihui.Muilt-generation compete genetic algorithm for logistics distribution vehicle routing problem[J].Journal of Railway Science and Engineering, 2015,2(5):75-79.

[12] Hammouche K, Diaf M, Siarry P.A comparative study of various meta-heuristic techniques applied to the multilevel thresholding problem[J].Engineering Applications of Artificial Intelligence, 2010, 23(5):676-688.

[13] 戴朝華,陳維榮,朱云芳,等.IIR數(shù)字濾波器設(shè)計(jì)的搜尋者優(yōu)化算法[J].西南交通大學(xué)學(xué)報(bào), 2009, 44(6): 871-876.

DAI ChaoHua, CHEN Weirong, ZHU Yunfang, et al.IIR digital filter design via seeker optimization algorithm[J].Journal of Southwest Jiaotong University, 2009, 44(6): 871-876.

[14] Dai Chaohua, Chen Weirong, Zhu Yunfang, et al.Reactive power dispatch considering voltage stability with seeker optimization algorithm[J].Electric Power System Research, 2009, 79(10): 1462-1471.

[15] 余勝威,曹中清.基于人群搜索算法的PID控制器參數(shù)優(yōu)化[J].計(jì)算機(jī)仿真,2014,31(9): 347-350.

YU Shengwei, CAO Zhongqing.Optimization parameters of PID controller paraters based on seeker optimization algorithm[J].Computer Simulation,2014,31(9): 347-350.

(編輯陽麗霞)

Improved seeker optimization algorithm application in weld image segmentation

YU Shengwei1, DING Jianming2, CAO Zhongqing1

(1.College of Mechanical Engineering, Southwest Jiaotong University, Chengdu 610031, China;

2.State Key Laboratory of Traction Power, Southwest Jiaotong University, Chengdu 610031, China)

Abstract:This paper present's a new novel method which is called improved seeker optimization algorithm .Gaussian membership function and its parameters of traditional SOA algorithm are optimized, and the step length is determined by uncertainty reasoning based on adaptive principle.Besides, improved SOA is able to decide the n-1 optimal for n-level threshold on a given image.Compared with PSO, GA and SA algorithm, testing results show that improved SOA enhances the performance in terms of convergence speed, stability, solution accuracy and global optimality.ISOA could also greatly overcome the shortcomings of traditional methods, such as local optima and slow convergence speed.Hence, improved SOA can be widely employed in practical engineering.

Key words:multi-scale segmentation; ISOA; class variance; PSO; GA; SA

通訊作者:丁建明(1981-),男,四川巴中人,博士,從事機(jī)械設(shè)計(jì)、故障診斷研究;E-mail:fdingjianming@126.com

基金項(xiàng)目:國家自然科學(xué)基金重點(diǎn)資助項(xiàng)目(61134002)

收稿日期:2015-06-11

中圖分類號(hào):TP391.4

文獻(xiàn)標(biāo)志碼:A

文章編號(hào):1672-7029(2015)06-1471-07

主站蜘蛛池模板: 99热这里只有精品在线观看| 91精品国产91欠久久久久| 九色视频在线免费观看| 国产成人精品亚洲日本对白优播| 国产亚洲精品无码专| 国产精品福利社| 色一情一乱一伦一区二区三区小说| 一区二区影院| av天堂最新版在线| 一边摸一边做爽的视频17国产| 国语少妇高潮| 夜夜拍夜夜爽| 亚洲午夜片| 六月婷婷激情综合| 91亚瑟视频| 狼友视频国产精品首页| 先锋资源久久| 中国黄色一级视频| 国产jizzjizz视频| 亚洲免费成人网| 久久99国产综合精品1| 四虎精品国产永久在线观看| 亚洲人成网7777777国产| 日韩午夜福利在线观看| AV无码一区二区三区四区| 精品国产香蕉在线播出| 91久久精品国产| 国产免费黄| 欧美国产精品不卡在线观看| 久久国产精品波多野结衣| 69av在线| 精品第一国产综合精品Aⅴ| 亚洲浓毛av| 亚洲制服丝袜第一页| 狼友视频一区二区三区| 国产亚洲精品在天天在线麻豆 | 免费A级毛片无码免费视频| 91福利在线观看视频| 国产福利一区视频| 欧美另类一区| 亚洲一区二区三区国产精品 | 亚洲成在人线av品善网好看| 在线看AV天堂| 国产91视频观看| 国产一二三区在线| 国产精品成| 91热爆在线| 无码精品国产dvd在线观看9久| 亚洲无卡视频| 亚洲日韩图片专区第1页| 中文成人在线视频| 成人免费一区二区三区| 中国美女**毛片录像在线| 永久免费av网站可以直接看的| 精品成人免费自拍视频| 香蕉在线视频网站| 91久久国产综合精品| 国产真实乱子伦视频播放| 香蕉视频在线精品| 国产欧美在线观看一区| www.狠狠| 亚洲免费黄色网| 国产在线视频欧美亚综合| 国产精品99久久久久久董美香| www亚洲天堂| 在线观看免费国产| 天天色天天综合网| 亚洲无码高清一区| 久久黄色影院| 四虎影视8848永久精品| 国产成人禁片在线观看| 草草线在成年免费视频2| 午夜无码一区二区三区| 国产欧美中文字幕| 毛片网站在线播放| 国产成a人片在线播放| 久久国产精品电影| 免费亚洲成人| 四虎精品黑人视频| 国产又粗又爽视频| 天天综合网色| 日韩经典精品无码一区二区|