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

介質(zhì)目標的太赫茲波近場散射特性計算

2015-07-12 13:59:39程志華謝擁軍馬曉東毛煜茹
電子與信息學報 2015年4期
關鍵詞:方法

程志華謝擁軍 馬曉東 毛煜茹 畢 博

(北京航空航天大學電子信息工程學院 北京 100191)

介質(zhì)目標的太赫茲波近場散射特性計算

程志華*謝擁軍 馬曉東 毛煜茹 畢 博

(北京航空航天大學電子信息工程學院 北京 100191)

該文基于廣義的Kirchhoff阻抗邊界條件和物理光學法,對太赫茲頻段介質(zhì)體近場散射特性進行了研究,給出了介質(zhì)體近場散射計算公式。針對在太赫茲頻段由于波長較短引起的計算量大幅提升的問題,采用以面元為計算單位、以像素為遮擋判斷單位的太赫茲頻段介質(zhì)體近場散射的快速計算方法,該方法在保證計算精度的基礎上,大大降低了計算復雜度和時間。計算了圓柱體和鴨嘴形介質(zhì)體在不同距離下的雷達散射截面,并且分析了電磁場與物體相互作用后,相位項在不同距離、不同頻率下對介質(zhì)體雷達散射截面的影響。

雷達散射截面;圖形電磁學;太赫茲;近場

1 引言

太赫茲(THz)波通常是指介于0.1~10.0 THz之間的電磁波,在電磁波譜中太赫茲波處于紅外和微波之間。由于太赫茲頻段的空間分辨率和時間分辨率頻率都很高,與微波和X射線成像相比,太赫茲成像不僅可以顯示更精密的圖像,還可以獲得頻域內(nèi)的相位信息,因此太赫茲成像技術是太赫茲應用的關鍵技術。太赫茲成像是一個典型的電磁逆散射問題,是通過接收目標散射數(shù)據(jù)重建描寫介質(zhì)或目標特性的基本函數(shù),同時,由于太赫茲頻段頻率較高,近場距離范圍大幅度提升,成像中散射模型應為近場散射計算模型。因此,太赫茲頻段介質(zhì)目標近場散射特性的研究是十分重要和有意義的。

太赫茲頻段目標散射特性研究可分為實驗測量、數(shù)值計算兩種手段。在實驗測量上,美國麻省LOWELL大學毫米波實驗室利用1.56 THz源在緊縮場中對粗糙面圓柱體的目標散射特性進行了研究[1]。新澤西理工學院物理系對太赫茲頻段危險物品的散射成像進行了研究[2]。天津大學太赫茲研究中心搭建了以0.2 THz 返波管振蕩器源、熱釋電探測器、小型自動旋轉(zhuǎn)光學平臺等組成的太赫茲波目標散射特性實驗測試系統(tǒng),并對粗糙銅面的散射特性等進行了研究[3]。首都師范大學太赫茲實驗室研制了太赫茲數(shù)字全息成像系統(tǒng),該系統(tǒng)可以對太赫茲電磁波的振幅、相位、頻率及偏振等全部光學信息的3維空間分布進行精確測量[4]。

在數(shù)值計算領域,東南大學毫米波實驗室對太赫茲波頻段的金屬目標體的遠場目標特性計算進行了研究,對比分析了多種電磁計算方法,并且獨立自主地開發(fā)了準確高效的目標特性和電磁兼容仿真分析軟件[5]。本課題組對太赫茲頻段金屬目標的近場散射特性計算進行了研究,并取得了一定成果。介質(zhì)體近場散射特性的研究較少見報道,遠場條件下,目標的散射特性已得到廣泛的研究[6?11],介質(zhì)目標體計算方法可以分為3類:一是以等效流為基礎的方法[12];二是近似邊界條件法(ABC)[13];三是數(shù)值方法,如基于積分方程的矩量法、快速多極子法、體積分方程法等,基于微分方程的有限元法、時域有限差分法等,高頻方法有射線追蹤法(SBR)、幾何光學法(GO)、物理光學法(PO)、幾何繞射理論(GTD)等。基于上述方法,近場條件下,目標體RCS的計算思路主要有兩種:一種是將物體表面剖分成小面元,用遠場RCS計算方法分別計算各小面元的RCS值,然后通過相位疊加得到總的散射場,是從部分到整體的思路[14]。另外一種是以物理光學法為基礎,通過散射場積分方程在無遠場近似的情況下推導得出的方法,是從整體到部分的思路[15]。在近場區(qū)域,不同距離下的近場散射分析更具實際意義,目標體照射區(qū)域?qū)㈦S著距離的不同而變化,隨之相互遮擋情況也會產(chǎn)生變化,從而導致大量的時間用在相互遮擋的判斷上,因此,采用上述的高頻方法仍需要耗費極大的內(nèi)存和極長的計算時間。1993年RIUS等人[16]提出了基于PO的圖形電磁學(GRECO)方法,該方法充分利用了計算機硬件的優(yōu)勢,由計算機硬件完成最困難、最費時的遮擋和消影工作,這使得它具有很大的工程應用價值,但是GRECO方法只適用于遠場區(qū)域的電磁散射計算,并且GRECO方法以提取像素顏色強度獲得參數(shù)的本質(zhì)決定其無法適用于太赫茲頻段近場的計算。如何把GRECO方法在遮擋判斷上的優(yōu)勢應用在太赫茲頻段近場散射計算中是十分有意義和價值的。

本文利用基于廣義的Kirchhoff的阻抗邊界條件對介質(zhì)體近場RCS計算公式進行了推導。基于PO和GRECO方法,提出一種以面元為計算單位、以像素為遮擋判斷單位的太赫茲波介質(zhì)體近場RCS的快速計算方法。該方法在保證計算精度的基礎上,大大降低了遮擋判斷的計算復雜度和時間。以圓柱體和鴨嘴形介質(zhì)體在不同距離下的雷達散射截面的計算為例,驗證了本文方法的有效性和準確性。

2 理論分析

2.1 介電性質(zhì)

德拜方程是討論一般介質(zhì)極化弛豫特性的重要關系式,它給出了介電常數(shù)和頻率的關系,可表示為

其中ε(0)為靜態(tài)介電常數(shù),ε(∞)為高頻節(jié)點常數(shù),τD為弛豫時間,當介質(zhì)體的弛豫時間的倒數(shù)遠大于電磁波頻率時,由式(1)可知,介質(zhì)介電常數(shù)可近似認為是ε(0)。太赫茲頻段,頻率在1011~1013Hz 之間,而大部分介質(zhì)體的弛豫時間大于皮秒級,因此在太赫茲頻段,大部分介質(zhì)體的介電常數(shù)將是與介質(zhì)體弛豫時間、高頻介電常數(shù)和靜態(tài)介電常數(shù)相關的變量。材料的弛豫時間和介電常數(shù)無法通過理論推導求出,需要通過實驗獲得,目前針對太赫茲頻段的物質(zhì)特性的研究大都通過對吸收系數(shù)和折射系數(shù)展開,通過折射系數(shù)和吸收系數(shù)可以求得介質(zhì)的介電常數(shù)。

2.2 近場介質(zhì)體散射計算理論

由Stratton-Chu積分公式,當散射體內(nèi)無電磁源時,無遠場近似的散射場為

其中Es為目標體散射電場,積分面S為光照區(qū)表面,ds為表面單位單元,其單位法向矢量為n, s為散射場方向單位矢量,J為電流密度,E和H分別為總電場和總磁場,Z0為波阻抗,φ=ejkr/r, φ1= ejkr/r2, r為源與面元ds的距離,其單位矢量為r。假設物體表面的曲率半徑遠大于波長,考慮接收單位的極化方向,由散射場積分方程可得

其中i為入射場方向單位矢量,iH為入射磁場。

由于是球面波入射,因此,

其中hi為入射磁場方向單位矢量,由式(3),式(4)可得

近場雷達散射截面公式可定義為

其中R為源點到目標體中心的距離。

利用廣義的Kirchhoff的阻抗邊界條件,介質(zhì)體表面電流和磁流可以表示為

由式(5),式(6)結(jié)合式(7),式(8)化簡后得到介質(zhì)體近場RCS的計算表達式為

R⊥,R//分別為介質(zhì)表面在垂直極化和水平極化時的反射系數(shù),對非磁性介質(zhì)進行考慮,μ1=μ2=μ0,其可表示為

其中ε1, ε2分別為空氣和媒質(zhì)介電系數(shù)。

3 太赫茲頻段近場散射快速計算方法

由式(9),式(10)可知,該計算公式是以面元為計算單位。在太赫茲頻段,為保證計算結(jié)果的準確性,剖分面元將會是千萬級的數(shù)量,甚至更高,以面元為遮擋判斷單位計算量十分巨大,因此,物理光學法在太赫茲頻段的RCS計算中工程應用性很低。以像素為計算單位的圖形電磁學方法,在判斷遮擋上具有很大的優(yōu)勢,但是,圖形電磁學中,計算所需參數(shù)是以像素為單位通過對光強值的提取來得到的,一般情況下顏色分配模式為RGBA模式,硬件為R, G, B, A成分保留一定數(shù)量的位平面。目前顯卡支持的最高顏色位數(shù)為32位,R, G, B和A各占8位,8 位的R成分可以存儲0~255之間的整數(shù)[17]。位平面中的0,1,…,255將對應于0/255=0.0, 1/255,…,255/255=1的顏色值,因此,通過光照提取出來的參數(shù)值都將是1/255=0.0039的倍數(shù),而太赫茲頻段頻率范圍在0.1~10.0 THz之間,最低頻率0.1 THz時,波長為0.003 m,對比通過光照提取的參數(shù)值的分辨率(1/255),可知以像素為計算單位,通過提取像素顏色強度來獲得所需參數(shù)的方法無法滿足太赫茲頻段的計算需求。同時,在近場的計算中,由于是球面波照射,照射區(qū)域和相位因素不能通過平面波光源和深度緩存得到,因此傳統(tǒng)的圖形電磁學中提取參數(shù)方法和遮擋判斷方法在太赫茲目標散射的計算中完全失效。本文通過以面元為計算單位保證太赫茲頻段散射特性計算精度,以像素為遮擋判斷單位保證遮擋判斷的快速性,最終通過建立像素與面元的對應關系從而達到快速計算太赫茲目標近場散射值的目的。由式(9),式(10)可知,所需計算的參數(shù)包括:各面元到光源的距離,各面元法向與入射波夾角的余弦值。

按照目標體實際尺寸,用FEMAP對其進行建模并進行剖分,剖分文件中包括各面元的法相矢量和位置信息。每一個面元的信息包括3個頂點的坐標(x1,y1,z1),(x2,y2,z2),(x3,y3,z3)和面元的法相矢量n,點源的坐標(x,y,z)可以通過光源的位置設置得到。根據(jù)坐標和法向可以得到各面元到點源的距離和入射波與面元法向的夾角余弦值。

在近場條件下,照射區(qū)域?qū)㈦S距離的變化而變化,隨著遮擋情況也會產(chǎn)生變化。通過OPENGL將光源設置為位置性光源,光源位置信息通過模型視圖矩陣進行變化,并以視覺坐標的形式存儲。光照區(qū)域?qū)㈦S光源位置的不同自動進行光照區(qū)域判斷,結(jié)合深度緩存自動進行遮擋處理,從而得到準確的光照區(qū)域。利用OPENGL將消隱后的目標圖像顯現(xiàn)在計算機屏幕上,如圖1所示。

圖1 隨距離變化的光照區(qū)域變化

光照區(qū)域的面元顏色是以像素為單位進行存儲的,無法直接判斷面元是否處于光照區(qū)域。通過OPENGL工具函數(shù)庫中的gluUnproject函數(shù)將屏幕上光照區(qū)像素與剖分體面元一一對應,從而以像素為光照區(qū)域判斷單位,當像素有顏色時,由對應的面元進行計算。在計算過程中會出現(xiàn)單個像素對應多個面元的情況,在這種情況下通過利用OPENGL的深度緩存判定離剪裁面最近的面元為計算單位,放棄其他面元的計算。為了防止出現(xiàn)單個像素對應多個面元的情況,采用遠場散射計算中的分屏顯示方法,保證像素對應的實際的尺寸小于或等于面元剖分精度。

以半徑0.254 m球體為例,入射頻率300 GHz,剖分精度為1/3波長(m),所需要的分屏數(shù)為

其中max(L), max(W)分別為目標體最大長度和最大寬度;pl, pw分別別為顯示區(qū)域像素的長度和寬度。該方法的計算量為

射線法是通過判斷各剖分面元中點與源點形成的射線與其他面元是否相交,然后保留與源點最近的面元進行計算。DirectX中采用的就是此方法,其具體計算方法是:通過求解三角形面元中兩條邊的權重u,v,如果u+v≤1則判定射線與該面元存在交點,其計算量為:n(n?1)×[N(u)+M(v)+1],其中n為面元數(shù),N(u), M(v)分別為u,v的計算量。

表1為同樣以半徑0.254 m球體為目標體,入射頻率為300 GHz,剖分精度為0.0003 m時,采用射線法判斷遮擋與本文方法的計算量對比。從表1可以看出,采用本文方法大幅度減少了計算量。

表1 射線法與本文方法遮擋判斷的計算量對比

4 計算結(jié)果及分析

為驗證本文方法的有效性,計算了邊長為0.3 m的介質(zhì)立方體在距離為0.04 m時的近場RCS,入射波長為1 m, εr=4, ur=1,入射方向為與立方體表面法向夾角5°~40°的范圍,計算間隔5°。文獻[18]中采用矩量法(MOM)對其進行了計算,兩者計算結(jié)果如圖2所示,從中可以看出當角度較小時與文獻[18]中計算結(jié)果與本文計算結(jié)果吻合較好,當角度增大時誤差略微增大(由棱邊引起),最大誤差小于1 dB,因此驗證了本文方法的正確性和有效性。

圖2 介質(zhì)立方體近場RCS本文計算結(jié)果與文獻計算結(jié)果對比

以介質(zhì)體為含碳納米結(jié)構的高分子復合材料為例[19],分別計算了圓柱體和同比尺寸縮小后的鴨嘴型散射體,計算頻率為100 GHz, 200 GHz和300 GHz,其介電常數(shù)如表2所示。

表2 含碳納米結(jié)構的高分子復合材料介電常數(shù)表

(1)介質(zhì)圓柱體計算距離為10~200 m,剖分精度為0.0005 m,正Z軸入射,垂直極化,模型如圖3所示,計算結(jié)果如圖4所示。

圖3 圓柱介質(zhì)體模型

從圖4可以看出當目標體與源距離較近時,同一距離下,不同頻率的近場RCS值并沒有像遠場RCS一樣隨頻率的增大而增大,隨著距離的增加,頻率與近場RCS值逐漸成單調(diào)趨勢。圖5為固定介電常數(shù)為12.45,圓柱介質(zhì)體在頻率1~200 GHz,距離分別為10 m, 50 m和200 m時的近場RCS值。可以看出,當距離較近時,近場RCS值隨頻率的增大波動較強;當距離較遠時,近場RCS值隨頻率的增加波動降低。由式(16),式(17)可知,當距離不變、頻率發(fā)生變化時,只有系數(shù)項4π/λ2與相位項e2jkr(其中k=2π/λ)發(fā)生變化,其中系數(shù)項與頻率成單調(diào)正比趨勢,結(jié)合圖5中的計算結(jié)果,可以認為在近場散射計算中,當距離越近時,相位項對近場散射值的影響越大;隨著距離的增加,相位項對近場散射值的影響越來越小,到達遠場條件時,散射值與頻率成正比,本文計算結(jié)果與遠場計算結(jié)果吻合。

(2)同比尺寸縮小后的鴨嘴型散射體,計算距離為1~50 m,正Z軸入射,水平極化,剖分精度0.0005 m,模型如圖6所示,計算結(jié)果如圖7所示。

圖4 圓柱介質(zhì)體10~200 m近場RCS

圖5 1~200 GHz固定距離圓柱介質(zhì)體近場RCS

圖6 鴨嘴型天線罩介質(zhì)體模型

圖7 鴨嘴型散射體不同頻率近場RCS

從圖7中同樣可以看出,當距離越遠時相位項對近場RCS值的影響越小,在近距離一定范圍內(nèi),近場RCS值并不隨頻率增加而增大,當?shù)竭_遠場條件時與遠場計算結(jié)果吻合。

5 結(jié)論

本文研究了太赫茲頻段介質(zhì)目標體的近場散射特性,基于廣義的Kirchhoff阻抗邊界條件和物理光學法對介質(zhì)體近場散射計算公式進行了推導,采用以面元為計算單位,以像素為遮擋判斷單位的太赫茲近場散射的快速計算方法。該方法利用OPENGL引入點光源和深度緩存結(jié)合的光照區(qū)域判斷方法,由計算機硬件自動完成耗時的遮擋和消隱以及建立面元與像素之間的對應關系,從而快速準確地計算了太赫茲頻段介質(zhì)目標的近場雷達散射截面,并且分析了電磁場與物體相互作用后,相位項在不同距離、不同頻率下對介質(zhì)體雷達散射截面的影響。

[1] Danylov A A, Goyette T M, Waldman J, et al.. Terahertz inverse synthetic aperture radar (ISAR) imaging with a quantum cascade laser transmitter[J]. Optics Express, 2010, 18(15): 16264-16272.

[2] Federici F, Schulkin B, Huang F, et al.. THz imaging and sensing for security applications-explosives, weapons and drugs[J]. Semiconductor Science and Technology, 2005, 20(2): S266-S280.

[3] 楊洋, 姚建銓, 張鏡水, 等. 粗糙銅表面對低頻太赫茲波的散射實驗[J]. 紅外與毫米波學報, 2013, 32(1): 36-39.

Yang Yang, Yao Jian-quan, Zhang Jing-shui, et al.. Terahertz scattering on rough copper surface[J]. Journal of Infrared and Millimeter Waves, 2013, 32(1): 36-39.

[4] 鄭顯華, 王新柯, 孫文峰, 等. 太赫茲數(shù)字全息術的研發(fā)與應用[J]. 中國激光, 2014, 41(2): 1-10.

Zheng Xian-hua, Wang Xin-ke, Sun Wen-feng, et al.. Developments and applications of the terahertz digital holography[J]. Chinese Journal of Lasers, 2014, 41(2): 1-10.

[5] Li Zhuo, Cui Tie-jun, Zhong Xing-jian, et al.. Electromagnetic scattering characteristics of PEC targets in the terahertz regime[J]. IEEE Antennas and Propagation Magazine, 2009, 51(1): 39-50.

[6] Solomon R C, Leong H, and Antar Y M. Forward scattering effects in RCS of complex targets in the 3-20 MHz high frequency range[C]. Antennas and Propagation Society International Symposium, Boston, 2008: 1-4.

[7] 李曉峰, 謝擁軍, 楊銳. 半空間電大導體目標散射的高頻分析方法[J]. 電子與信息學報, 2009, 31(5): 1268-1270.

Li Xiao-feng, Xie Yong-jun, and Yang Rui. High-frequency analysis on scattering from conductive targets with electrically large size in half space[J]. Journal of Electronics & Information Technology, 2009, 31(5): 1268-1270.

[8] Gillion E, Rochefort E, and Claverie J. Improvement of RCS estimation of large targets by using near-field approach[C]. IEEE Radar Conference (RADAR), Ottawa, ON, Canada, 2013: 1-5.

[9] Altin N and Yazgan E. RCS prediction using fast ray tracing in Plücker coordinates[C]. 2013 7th European Conference on Antennas and Propagation (EuCAP), Gothenburg, 2013: 284-288.

[10] Lipuma D, Me?ric S, and Gillard R. RCS enhancement of flattened dihedral corner reflector using reflect array approach[J]. Electronics Letters, 2013, 49(2): 152-154.

[11] Li Xiao-feng, Xie Yong-jun, and Yang Rui. High-frequency method for scattering from coated targets with electrically large size in half space[J]. IEEE Transactions on Antennas and Propagation, 2009, 57(2): 181-186.

[12] Huddleston P L, Medgyesi-Mitschang L N, and Putnam J M. Combined field integral equation formulation for scattering by dielectrically coated conducting bodies[J]. IEEE Transactions on Antennas and Propagation, 1986, 34(4): 510-520.

[13] Senior T B A and Volakis J L. Generalized impedance boundary condition in scattering[J]. Proceedings of the IEEE, 1991, 79(10): 1413-1420.

[14] Li Jing and Li Xiang-jun. An efficient algorithm for near field RCS of electrically large dynamic targets[C]. International Conference on Information Science and Technology, Nanjing, 2011: 413-415.

[15] Pouliguen P, Damiens J F, et al.. RCS computation in near field[C]. Proceedings of the International Conference Day on Diffraction, Petersburg, 2006: 252-265.

[16] Rius J M, Ferrando M, and Jofre L. High frequency RCS of complex radar targets in real time[J]. IEEE Transactions on Antennas and Propagation, 1993, 41(9): 1308-1319.

[17] Shreniner D. The Khronos OpenGL ARB Working Group.

李軍, 徐波, 等譯. OPENGL編程指南[M]. 第7版, 北京: 機械工業(yè)出版社, 2010: 118-119.

[18] 賀昌輝. 矩量法分析目標的近場電磁散射[D]. [碩士論文], 華中師范大學, 2006.

He Chang-hui. The analysis of near-field EM scattering from targets by moment methods[D]. [Master dissertation], Central China Normal University, 2006.

[19] Zhuravleva V A, Suslyaeva V I, and Dunaevskiia G E. Complex permittivity of polymer composites containing carbon nanostructures in frequency range 0.17-1.1[C]. International Conference on THz Infrared, Millimeter and Terahertz Waves (IRMMW-THz), Wollongong, NSW, 2012: 1-2.

程志華: 男,1982年生,博士生,研究方向為計算電磁學、太赫茲目標特性、天線工程、電磁兼容等.

謝擁軍: 男,1968年生,教授,研究方向為計算微波與計算電磁學、微波通信、天線工程、電磁兼容等.

毛煜茹: 女,1984年生,博士生,研究方向為計算電磁學、天線工程、電磁兼容等.

馬曉東: 男, 1989年生,碩士生,研究方向為目標特性、天線工程、電磁兼容等.

畢 博: 男,1993年生,碩士生,研究方向為目標特性等.

Near-field Electromagnetic Scattering Characteristics of Dielectric Targets in the Terahertz Regime

Cheng Zhi-hua Xie Yong-jun Ma Xiao-dong Mao Yu-ru Bi Bo
(School of Electronic and Information Engineering, Beihang University, Beijing 100191, China)

The near-field scattering characteristics of dielectric targets are studied in the terahertz band based on the generalized Kirchhoff impedance boundary conditions and the physical optics method. The formula of the near field scattering for dielectric targets is deduced. In the light of the increase in calculation amount caused by the shorter wavelength, a fast computational method using surface element as the calculating unit and pixel as occlusion judgment unit is proposed for the near field scattering computation in terahertz band. The method ensures the calculation accuracy and reduces greatly the computational complexity and the time consumption of occlusion judging. The calculation of the near-field Radar Cross Section (RCS) produced by a dielectric cylinder and a duck mouth scatterer is performed in terahertz band. Meanwhile, the effect of the phase on the near field RCS in different distances and frequency is analyzed.

Radar Cross Section (RCS); GRaphical-Electromagnetic COmputing (GRECO); Terahertz; Near field

TN011

: A

:1009-5896(2015)04-1002-06

10.11999/JEIT140807

2014-06-20收到,2014-08-15改回

*通信作者:程志華 asakaka1982@163.com

猜你喜歡
方法
中醫(yī)特有的急救方法
中老年保健(2021年9期)2021-08-24 03:52:04
高中數(shù)學教學改革的方法
河北畫報(2021年2期)2021-05-25 02:07:46
化學反應多變幻 “虛擬”方法幫大忙
變快的方法
兒童繪本(2020年5期)2020-04-07 17:46:30
學習方法
用對方法才能瘦
Coco薇(2016年2期)2016-03-22 02:42:52
最有效的簡單方法
山東青年(2016年1期)2016-02-28 14:25:23
四大方法 教你不再“坐以待病”!
Coco薇(2015年1期)2015-08-13 02:47:34
賺錢方法
捕魚
主站蜘蛛池模板: 久久人体视频| 日韩黄色精品| 国产精品无码在线看| 国产成年女人特黄特色毛片免 | 国产内射一区亚洲| 在线看片免费人成视久网下载| 99热这里只有成人精品国产| 国产一级无码不卡视频| 91久久偷偷做嫩草影院电| 97成人在线观看| 99视频精品全国免费品| 国产视频久久久久| 日韩免费毛片| 色久综合在线| 国内精自视频品线一二区| 99爱视频精品免视看| 国产精品页| P尤物久久99国产综合精品| 香蕉伊思人视频| 久久婷婷人人澡人人爱91| 日本www色视频| 免费在线a视频| 亚洲色图欧美视频| 午夜视频www| 91福利在线观看视频| 午夜国产精品视频| 91免费国产高清观看| 亚洲精品国产综合99| 福利在线不卡| 日本亚洲欧美在线| 国产无人区一区二区三区| 最新国产你懂的在线网址| 欧美成人A视频| 欧美一区二区啪啪| 国产aⅴ无码专区亚洲av综合网| 国产第三区| 国产黄在线观看| 国产菊爆视频在线观看| 国产jizzjizz视频| 中文字幕在线看| 午夜电影在线观看国产1区| 国产成在线观看免费视频| 小说 亚洲 无码 精品| 亚洲va视频| 国产av色站网站| 亚洲AV电影不卡在线观看| 六月婷婷精品视频在线观看| 久久国产亚洲欧美日韩精品| 精品人妻无码中字系列| 日韩亚洲高清一区二区| 国内精品小视频福利网址| 国产在线小视频| 在线播放国产99re| 欧美不卡视频在线观看| 亚洲欧美日韩中文字幕一区二区三区| 精品国产自在在线在线观看| 国产一区二区色淫影院| 欧美黄网站免费观看| 免费观看欧美性一级| 最新日本中文字幕| JIZZ亚洲国产| 免费看黄片一区二区三区| 国产va欧美va在线观看| 伊人久久影视| 六月婷婷综合| 久久国产毛片| 婷婷综合在线观看丁香| 青青操国产视频| 日本成人不卡视频| 国产欧美性爱网| 亚洲天堂成人在线观看| 日韩亚洲综合在线| 澳门av无码| 国内精品久久人妻无码大片高| 97影院午夜在线观看视频| 免费国产黄线在线观看| 亚洲欧美综合在线观看| 激情在线网| 欧美a在线看| 国产国语一级毛片| 91免费在线看| 欧美一级在线|