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

靶道空間三點法標定與誤差分析

2016-04-18 07:40:23王唯唐志華南京理工大學瞬態物理國家重點實驗室江蘇南京210094
兵工學報 2016年2期

王唯,唐志華(南京理工大學瞬態物理國家重點實驗室,江蘇南京210094)

?

靶道空間三點法標定與誤差分析

王唯,唐志華
(南京理工大學瞬態物理國家重點實驗室,江蘇南京210094)

摘要:為了避免現有激光準直標定系統在激光光斑對準過程中的困難,簡化標定系統組成,提高標定效率和標定精度,增強標定系統的場地適應性。利用空間內任意不共線的三點分別在標定坐標系和靶道坐標系下的坐標構造過渡坐標系,得到標定坐標系到靶道坐標系的轉換關系。研究標定坐標系下標定點的位置偏差以及靶道坐標系下對應點的測量誤差對最終待測點測量誤差自協方差的影響,通過雅克比矩陣得到相應的誤差模型,導出坐標測量誤差均值和方差的表達式。依據向量范數和矩陣范數之間的關系以及Cauchy-Schwarz不等式研究測量誤差均值與方差的增益系數的理論上限。通過實彈擊發實驗,考察三點標定方法的實踐效果。研究結果表明:最終測量誤差均值向量的模長等于被測點在標定坐標系下均值向量的模長;確定標定坐標系到靶道坐標系轉換關系的任意不共線三點之間距離及其包圍的面積越大,則測量誤差的自協方差增益系數越小;將過渡坐標系原點靠近預計彈道線有利于減小測量誤差的自協方差增益系數;取決于測量方法及其誤差形式,隨著標定坐標系原點遠離靶道坐標系原點,與靶道坐標系中三點坐標測量誤差相關的被測點測量誤差自協方差增益有所不同;基于三點標定方法,使用常見的通用儀器、設備就可以完成標定過程,且經過三點法標定的靶道測量數據質量較高。

關鍵詞:兵器科學與技術;靶道空間標定;靶道實驗;誤差控制;誤差分析

0 引言

在靶道測量過程中,利用彈道靶道配備的多臺套陰影閃光成像系統[1]捕捉彈丸飛行過程中的姿態及空間位置信息。為了能夠正確綜合各個照相站所捕捉到的彈丸在飛行過程中的不同位置信息,準確確定照相站在靶道中的空間位置,提高彈道參數測量結果的準確性,需要對照相站及靶道的空間進行標定[2-5]。現有的標定方法主要依據陰影成像系統和以雙重田字網格為載體的靶道空間基準系統[6-9]建立相機圖片到靶道空間坐標系的轉換關系。本文通過確定空間內任意不共線的三點分別在標定坐標系[5]和靶道坐標系下的坐標,得到標定坐標系到靶道坐標系的轉換關系,據此建立照相站的測量誤差模型,考察標定模型參數[10]對測量誤差的影響和作用,提高了標定方法的場地適應性。通過標定坐標系到靶道坐標系的轉換矩陣研究測量誤差對最終靶道空間測量數據的影響,分析各個誤差來源對測量結果所造成的影響,找出相應的處理方法。最后通過實彈擊發實驗驗證了三點標定法的實際效果。

1 標定模型的建立

1.1標定坐標系的轉換

定義A、B、C為空間中任意不共線三點,則A、B、C三點在標定坐標系中可以表示為: A1=, B1= [χB1yB1zB1]T, C1=.則由A、B、C三點可以確定向量:

由V1和U1可以確定一個以V1為X軸,(V1× U1)×V1為Y軸,V1×U1為Z軸的過渡坐標系OXYZ.若過渡坐標系OXYZ中任意一點的坐標為[χ y z]T,則在標定坐標系中其坐標表示為

式中:M1=表示由過渡坐標系OXYZ轉換到標定坐標系的方向余弦矩陣。

1.2靶道坐標系的轉換

同理,令A、B、C三點在靶道坐標系中的坐標為A2=, B2=, C2=,則由A、B、C三點在靶道坐標系下的坐標可以確定的向量:

相應的過渡坐標系OXYZ中的坐標[χ y z]T在靶道坐標系中的坐標表示為

式中:M2=表示由過渡坐標系OXYZ轉換到靶道坐標系的方向余弦矩陣。聯立(3)式和(6)式可得到標定坐標系到靶道坐標系的轉換關系為

式中: M = M2M1-1表示對標定坐標系中點的坐標旋轉;- M·表示對標定坐標系中點的坐標平移。

2 標定模型的誤差分析

由(7)式可知,影響[χ2y2z2]T的主要因素有:A、B、C三點在標定坐標系下選取時所產生的偏差和;A、B、C三點在靶道坐標系下測量所得的偏差、δB2=和以及坐標[χ1y1z1]T的偏差.于是A、B、C三點在標定坐標和靶道坐標的坐標真值與測量值的關系可以分別表示為=,i = 1,2.其中下標i為1時表示標定坐標系,下標i為2時表示靶道坐標系。

2.1A、B、C的誤差特性

2.1.1A、B和C三點在標定坐標系下的誤差特性

由于A、B、C為空間中任意取定的不共線三點,可以將其偏差Δ1= [δA1δB1δC1]T視為獨立的隨機變量,且均服從期望μ1,方差ε的正態分布。即

式中:I9×9是9階單位方陣。

2.1.2A、B和C三點在靶道坐標系下的誤差特性

靶道坐標系中,為了方便A、B、C三點的測量,取靶道空間觀測點為靶道坐標系原點O2,過原點的水平面為O2X2Z2坐標面,包含基準射向的鉛直面為O2X2Y2坐標面,上述兩坐標面的交線即為X2軸,順射向為正。過原點垂直于X2軸的平面為O2Y2Z2坐標面,該面與O2X2Y2坐標面的交線為Y2軸,向上為正,與O2X2Z2坐標面的交線為Z2軸,Z2軸的正方向使X2、Y2、Z2三坐標軸構成右手螺旋坐標系。點A在靶道坐標系下的位置關系如圖1所示,可以測得點A到靶道坐標系原點O2的距離為LA,O2A在O2X2Z2面內的投影與靶道坐標系X2軸的夾角為αA,以及O2A與O2X2Z2的夾角為βA.同理可得到, 點B到靶道坐標系原點O2的距離為LB,O2B在O2X2Z2面內的投影與靶道坐標系X2軸的夾角為αB,以及O2B與O2X2Z2的夾角為βB.點C到靶道坐標系原點O2的距離為LC,O2C在O2X2Z2面內的投影與靶道坐標系X2軸的夾角為αC,以及O2C與O2X2Z2的夾角βC.

圖1 A點在靶道空間的位置Fig.1 Coordinate of Point A in ballistic range system

令LA=LA+ΔLA、LB=LB+ΔLB、LC=LC+ΔLC、αA=A+ΔαA、αB=B+ΔαB、αC=C+ΔαC、βA=A+ΔβA、βB=B+ΔβB、βC=C+ΔβC,其中ΔLA、ΔLB、ΔLC、ΔαA、ΔαB、ΔαC、ΔβA、ΔβB、ΔβC分別表示真值測量時所產生的偏差,于是A、B、C在靶道坐標系中的坐標表示為

令ΔL= [ΔLAΔLBΔLC]T,Δα= [ΔαAΔαBΔαC]T和Δβ= [ΔβAΔβBΔβC]T分別表示影響A、B和C三點在靶道坐標系下的測量誤差來源。由于測量過程中,測量方法及手段相同,可以將ΔL、Δα和Δβ中的元素分別視為服從期望均為μ2,方差分別為εL、εα和εβ的正態分布的獨立隨機變量。即

式中:I3×3是3階單位方陣。

2.2靶道測量值的誤差特性分析

令Φ1=表示影響[χ2y2z2]T的各個誤差來源,將(7)式按Φ1做級數展開,并忽略各個偏差小量的高次項和耦合項。得到各個偏差小量Φ1對[χ2y2z2]T所造成的影響ΔT為

式中: P3×21=是 [χ2y2z2]T關于Φ1的雅克比矩陣,其中P1、P2、P3和P4分別表示為[χ2y2z2]T關于Δ1、ΔL、Δα和Δβ的雅克比矩陣,除了P1的維數為3×9外,P2、P3和P4的維數均3×3.

對(19)式兩邊取均值,即

將(8)式、(13)式、(14)式和(15)式代入(20)式可得

考察Pi,i =1,2,…,4可知,其行和分別為0,即P1ij=0,i = 1,2,3,Pkij= 0,i = 1,2,3,k = 2, 3,4.于是(21)式右端前4項均為0,因此將(21)式化簡為

(22)式表明靶道坐標下坐標[χ2y2z2]T的測量偏差ΔT的均值向量只取決于對應點[χ1y1z1]T在標定坐標系下的測量偏差ξ的均值向量。由于M是正交矩陣,不改變變換后向量的模長,即得到|E(ΔT) | = |E(ξ) |.

在豆莢成熟、即豆粒將長成圓形,但莢內種子還不很明顯看得見時采收。若采收太遲則豆莢過熟,纖維過多影響品質,一般在盛果期間隔2—3天采一次。

對(19)式兩邊取方差可得

而Var(Φ1)可以表示為

diag([Var(Δ1) Var(ΔL) Var(Δα) Var(Δβ) Var(ξ)]),于是化簡(23)式可得

式中:Ki = PiPT

i,i =1,2,…,4,分別表示ε、εL、εα和 εβ的增益系數矩陣。考察(24)式左端協方差矩陣的主元,即測量值的自協方差。可知,當各個偏差的方差特性一定時,使系數矩陣Ki以及MVar(ξ)MT的主元分別較小時可令測量值的自協方差較小。

用符號tr(·)表示矩陣求跡運算,對于MVar(ξ)MT分項,因為M為正交矩陣,即有tr(MVar(ξ)MT)等于tr(Var(ξ)),而Var(ξ)表示被測點在標定坐標下的坐標偏差,與被測點在標定坐標下的空間位置相關,其具體表達可以參考文獻[3]。以下具體論述減小系數矩陣Ki的跡的方法和途徑。

對于K1,由于P1可以表示為

考察[χ2y2z2]T的(7)式,注意到旋轉矩陣M = M2M1-1,而M2只取決于A、B和C三點在靶道坐標系下的測量值,即?M2/?Δ1= 0.又因為?([χ2Ay2Az2A]T) /?Δ1=0,于是將(25)式化簡為

式中:W1=(M1-1[χ1-χ1A

z1- z1A]T)是 3×9維矩陣。于是K1可以表示為

由于M2為正交矩陣,因此K1是矩陣W1WT1的一個相似變換。而相似變換不改變矩陣的跡,進一步將W1按行分割矩陣并表示為[W11W21W31]T,按照矩陣跡的定義得

式中:符號‖·‖2在下文中除非有特別說明均表示向量的2范數。考察行向量Wi1的具體表達,整理后得

式中:wi1, i = 1,2,3分別是維數為9×3的系數矩陣。由向量范數和矩陣范數的關系可知:式中:符號‖·‖F表示矩陣的Frobenius范數。直接計算矩陣wi1的Frobenius范數,化簡得到‖w11‖2F= 4/ r21、‖w21‖2F= 4/ r21t21、‖w31‖2F= 4/ t21,其中r1= |V1|,t1= |V1×U1|,于是K1矩陣的跡滿足不相等關系:

對于Ki, i =2,3,4,由于P2、P3和P4可以分別表示為

式中:Δ2= [δ2Aδ2Bδ2C]T表示A、B和C三點在靶道坐標系中的偏差。將(32)式、(33)式和(34)式表示為

與前述處理W1的方法類似,將Wj, j =2,3,4表示為向量形式: Wj= [ W1jW2jW3j]T+[ I1I2I3]T,其中I1=[ 1 0 0 0 0 0 0 0 0 ], I2= [ 0 1 0 0 0 0 0 0 0 ]和I3= [0 0 1 0 0 0 0 0 0],且Wij= wij(M· [χ1-χA1y1- yA1z1- zA1]T), i =1,2,3.則

由三角不等式以及Cauchy-Schwarz不等式, (38)式、(39)式和(40)式具有不相等關系:

式中:‖?Δ2/?ΔL‖2、‖?Δ2/?Δα‖2和‖?Δ2/?Δβ‖2分別表示?Δ2/?ΔL、?Δ2/?Δα和?Δ2/ ?Δβ的矩陣2范數。按照矩陣2范數定義有

由向量范數和矩陣范數的關系可知:

式中:i =1,2,3; j =2,3,4.由于M為正交矩陣,只改變原有向量的方向不改變向量模長,則有

注意到V1、U1和V2、U2是A、B和C三點在不同坐標系下的向量表示。并有關系|V2| = |V1| = r1和|V2×U2| = |V1×U1| = t1.進一步計算w1j、w2j和w3j的Frobenius范數,則有

應用(44)式~(49)式,化簡(41)式~(43)式得

由(31)式、(50)式、(51)式和(52)式可知, tr(Ki),i =1,2,…,4的上限與r1和t1的大小以及被測點相對A點的距離相關。增大r1和t1可以使tr(Ki), i =1,2,…,4減小,即削弱了ε、εL、εα和εβ對測量偏差自協方差的影響。由于r1為向量AB的模長,t1是向量AB與向量AC確定的四邊形面積。所以增大A、B兩點距離并增加A、B和C三點所確定四邊形的面積可以使tr(Ki), i =1,2,…,4減小。同時在標定過程中,將標定坐標系內的A點靠近預計彈道線將有利于減小測量誤差的方差。由于分析中3點在靶道坐標中的測量偏差用球坐標形式表達,當標定坐標系遠離靶道坐標系原點時,tr(K3)和tr(K4)增大,即說明當被測點遠離靶道坐標原點時測量誤差增益增大。若可以直接獲得A、B和C三點在靶道坐標系下的測量誤差,那么(44)式、(45)式和(46)式左端均退化為對自身的雅可比矩陣的2范數,即偏差增益系數與被測點知靶道坐標原點的距離無關。

3 標定與擊發實驗

3.1標定實驗

基于現有的靶道陰影成像系統,為了應用三點法完成彈道靶道標定,需要使用全站儀( Leica-Ts11),標定板(定制玻璃刻蝕棋盤格)和標定架(按實際標定環境制作)。其中全站儀用于建立靶道坐標系,標定板安裝在標定架上實現相機鏡頭畸變修正并建立標定坐標系。標定的步驟大致分為鏡頭畸變修正和標定坐標在靶道坐標中位置、姿態的確定兩部分,按照現場條件可以同時或分步進行。只要保持照相系統位置、狀態不變,標定工作亦可在實驗后進行。完成標定后將標定板作為靜態標志物,提取圖像上標定板的特征點,應用標定參數可以計算出標定板上對應特征點在標定坐標中的坐標,經實踐驗證計算得到的特征點坐標與實際坐標的偏差小于0.1 mm,有關標定過程的詳細情況參考文獻[11]。針對文獻[11]使用的標定系統,當被測彈丸在各站標定坐標原點附近時,與圖像識別誤差相關的測量誤差增益系數在0.9左右。

3.2擊發實驗

在經過標定后的靶道進行擊發實驗。發射平臺為7.62 mm線膛彈道槍,彈藥為7.62 mm制式彈,光源控制與時序測量使用自研的控制板卡,其時間分辨率為20 ns,參試陰影成像系統為現有的3號、4號、8號、10號~13號、15號和19號站。完成標定后,與A、B、C三點在靶道坐標下的測量誤差相關的增益系數上限在6左右,與A、B、C三點在標定坐標下的測量誤差相關的增益系數上限在3左右。實驗后提取各個圖像上彈丸特征點(以彈頭為例)的像素坐標,通過標定參數轉換坐標并求解相關方程,即獲得靶道坐標系下的彈頭點坐標。

由于射角較小、接近平射且彈丸飛行的時間和距離均較短,故可以忽略偏流效應,靶道內也不存在橫風擾動,所以彈丸飛行軌跡在水平面內的投影近似為直線,可以使用直線方程擬合。彈丸受到空氣阻力影響,其在豎直面內的軌跡不完全符合拋物線規律,但是在短時間內可以用二次多項式擬合數據來考核測量數據質量。對于本組實驗數據,由于彈道傾角的絕對量始終較小,彈道高變化不大,認為彈丸在水平面內速度分量的大小滿足近似方程υ·+ aυ2=0,該方程具有通解υ= (at + b)-1,其中a和b是待定系數。

圖2給出的是靶道坐標中彈頭坐標在水平面內的投影及其與擬合直線之間的偏差。從圖2中可以看出,射向較靶道坐標X2軸向略偏左,實測數據與直線匹配良好,偏差絕對值主體小于2 mm,最大不超過4 mm.

圖3給出的是彈丸頭部在射向豎直面內的軌跡及其同擬合曲線之間的偏差。從圖3中可以看出,實測數據與二次曲線匹配良好,偏差絕對值主體小于1 mm,最大不超過2 mm.

圖2 不同時刻彈丸頭部位置測量值在水平面內的投影坐標Fig.2 Measured projected coordinates of projectile nose on horizontal plane

圖3 彈丸頭部在射向豎直面內的空間軌跡Fig.3 Measured projected coordinates of projectile noseon vertical plane in the direction of fire

圖4給出的是彈丸頭部速度水平分量的時間歷程及其同擬合曲線的偏差。圖4中時間軸表示的是各站閃光時刻相對于第1次閃光(3號站閃光時間)的相對時間,速度分量是采用站間差分方法獲得的。從圖4中可以看出,數據與擬合曲線偏差較小,速度偏差絕對值小于0.1 m/ s.

圖4 彈丸頭部速度在水平面內的分量Fig.4 Velocity component of projectile nose on horizontal plane

4 結論

本文通過空間任意不共線的三點分別在標定坐標系和靶道坐標系下的坐標,得到標定坐標系到靶道坐標系的轉換關系。通過研究測量值關于各個誤差來源的雅克比矩陣得到了相應的誤差模型,導出坐標測量誤差均值和方差的表達式。研究測量誤差自協方差隨標定點位置的變化規律,基于范數不等式得到減小測量誤差的標定點分布規律。最終通過實驗來驗證標定方法的效果。得出以下結論:

1)靶道坐標下被測點測量誤差的均值向量只取決于被測點在標定坐標系下的均值向量,并且其模長大小不變。

2)在標定過程中,增大A、B兩點距離并增加A、B和C三點所確定的四邊形面積有利于減小被測點測量誤差的自協方差。

3)將標定坐標系內的A點(過渡坐標系原點)靠近預計彈道線,有利于減小被測點測量誤差的自協方差。

4)空間點測量誤差自協方差的增益系數不與被測點到靶道坐標原點的距離直接相關,考慮到測距誤差通常與距離相關,改善三點的測距誤差有利于在整個靶道長度上獲得較為一致的方差增益系數。

5)使用近似簡化的方法對實驗數據進行了粗略的擬合,盡管被擬合函數并不代表真實的彈道(6自由度彈道方程),但反映了測量數據較小的離散情況,并在一定程度上體現了標定方法較好的性能。

參考文獻(References)

[1]顧金良,陳平,夏言,等.數字式靶道陰影照相系統[J].彈道學報, 2009, 21(4):38-41.GU Jin-liang, CHEN ping, XIA Yan, et al.Digital ballistic range shadowgraph system[J].Journal of Ballistics, 2009, 21(4):38-41.(in Chinese)

[2]Zhang Z Y.A flexible new technique for camera calibration[J].IEEE Transactions on Pattern Analysis & Machine Intelligence, 2000, 22(11):1330-1334.

[3]王唯,唐志華,羅紅娥,等.照相站空間基準標定方法及測量誤差分析[J].兵工學報, 2014, 35(9):1414-1418.WANG Wei, TANG Zhi-hua, LUO Hong-e, et al.Calibration of shadowgraphy system and measuring error analysis[J].Acta Armamentarii, 2014, 35(9):1414-1418.(in Chinese)

[4]Tsai R Y.An efficient and accurate camera calibration technique for 3D machine vision[C]∥Proceedings of Conference on Computer Vision & Pattern Recognition Robotics and Automation.Miami, FL, US:IEEE, 1986:364-374.

[5]李竹良,趙宇明.基于單幅圖片的相機完全標定[J].計算機工程, 2013, 39(11):5-8.LI Zhu-liang, ZHAO Yu-ming.Full camera calibration based on single image[J].Computer Engineering, 2013, 39(11):5-8.(in Chinese)

[6]李觀濤.彈道靶道空間基準系統的設計原理[J].彈道學報, 1991,3(2):60-69.LI Guan-tao.Design principle for ballistic range reference frame [J].Journal of Ballistics, 1991,3(2):60-69.(in Chinese)

[7]任國民,李觀濤,張薇.彈道靶道空間坐標測量誤差的初步估計[J].彈道學報,1995, 7(3):57-64.REN Guo-min, LI Guan-tao, ZHANG Wei.Evaluation in measuring error of space coordinate for ballistic range[J].Journal of Ballistics, 1995, 7(3):57-64.(in Chinese)

[8]劉世平,易文俊,顧金良,等.彈道靶道數據判讀與處理方法研究[J].兵工學報,2000,21(3):201-204.LIU Shi-ping, YI Wen-jun, GU Jin-liang, et al.A new method for mage analysis and data extraction of projectiles in flight[J].Acta Armamentarii, 2000, 21(3):201-204.(in Chinese)

[9]馬鎖冬,朱日宏,李建欣,等.瞬態飛行目標三維面形的多視角測量系統[J].中國激光, 2010,37(12):3091-3097.MA Suo-dong, ZHU Ri-hong, LI Jian-xin, et al.A multi-view measurement system for three-dimensional surface distribution of transient moving target [ J].Chinese Journal of Lasers, 2010, 37(12):3091-3097.(in Chinese)

[10]尹洪濤,劉成,李一兵,等.相機標定誤差因素分析[J].信息通信, 2012(1):28-30.YIN Hong-tao, LIU Cheng, LI Yi-bing, et al.Analysis of factors on the error of the camera calibration[J].Information & Communications, 2012(1):28-30.(in Chinese)

[11]唐志華,王唯.一種靶道空間標定方法[C]∥2015年彈道學術會議論文集.南京:中國兵工學會彈道專業委員會, 2015: 184-189.TANG Zhi-hua, WANG Wei.A method of ballistic range calibration[ C]∥Proceedings of the 2015 Conference on Ballistics.Nanjing: Ballistic Professional Committee, China Ordnance Society, 2015: 184-189.(in Chinese)

Three-point Method and Error Analysis of Ballistic Range Calibration

WANG Wei, TANG Zhi-hua
(Science and Technology on Transient Physics Laboratory, Nanjing University of Science and Technology, Nanjing 210094, Jiangsu, China)

Abstract:To avoid the difficulty from laser spot alignment of existing laser collimation system in calibration process, a system composition is simplified, the efficiency and accuracy of calibration are improved, and the adaptability of calibration system is enhanced.A transformation matrix between calibration system and ballistic range system is presented by using the coordinates of three non-collinear points in calibration system and ballistic range system.The mean value and the variance-covariance matrix of the final measuring error are studied by an error model founded by Jacobian matrix between the measuring and various error sources.The self-covariance of measuring error is studied based on the Cauchy-Schwarz inequality and the relationship between vector and matrix norms.An firing experiment of live ammunition is per-book=318,ebook=129formed to verify the three-point calibration method of ballistic range.The research results show that the mean vector length of final measuring error about a point in ballistic range system is equal to the mean vector length of measuring error associated with same point in the calibration system; the larger the distance and the area surrounded by three non-collinear points which links the calibration system and ballistic range system together are, the smaller the self-covariance of the final measuring error is; the self-covariance decreases when the transition coordinate origin which is defined by the three non-collinear points is close to a projected trajectory line; based on the measurement method and the error of the three-point coordinates in ballistic system, the self-covariance of measuring error are changed when the calibration system origin moves away from the ballistic range system origin; and the calibration system based on the three-point method can be build-up by general instruments, and the quality of measured data is better.

Key words:ordnance science and technology; calibration of ballistic range; ballistic range experiment; error control; error analysis

作者簡介:王唯(1979—),男,講師。E-mail:office2007@ njust.edu.cn;唐志華(1991—),男,碩士研究生。E-mail:974724295@ qq.com

基金項目:瞬態物理國家重點實驗室基金項目(9140C300306120C30112)

收稿日期:2015-03-22

DOI:10.3969/ j.issn.1000-1093.2016.02.018

中圖分類號:TB89

文獻標志碼:A

文章編號:1000-1093(2016)02-0317-08

主站蜘蛛池模板: 美美女高清毛片视频免费观看| 久久综合色天堂av| 亚洲福利一区二区三区| 中文字幕免费播放| 国产另类视频| 色婷婷成人网| 午夜爽爽视频| 色亚洲成人| a毛片基地免费大全| 久久夜色精品| 久久精品人人做人人爽电影蜜月| 大香网伊人久久综合网2020| 欧美特黄一级大黄录像| 欧美啪啪精品| 欧美激情视频一区| 伊人久久久久久久久久| 日韩a在线观看免费观看| 久久semm亚洲国产| 黄色不卡视频| 亚洲色图综合在线| 欧美精品高清| 国产91在线|中文| 青草娱乐极品免费视频| 欧美国产三级| 国产精品不卡永久免费| 熟女成人国产精品视频| 丁香五月婷婷激情基地| 久久久噜噜噜久久中文字幕色伊伊| 欧美专区在线观看| 亚洲天堂视频网| 国产麻豆91网在线看| 国产亚洲欧美日韩在线一区| 欧美视频免费一区二区三区| 国产chinese男男gay视频网| 国产一区二区三区免费观看| 国产精品理论片| 在线观看91精品国产剧情免费| 久久精品国产亚洲麻豆| 亚洲综合二区| 国产精品亚洲综合久久小说| 网友自拍视频精品区| 亚洲人成网站色7777| 欧美第九页| 亚洲欧美日韩高清综合678| 人妻丰满熟妇av五码区| 丁香婷婷激情网| 成人福利在线免费观看| 中文字幕日韩丝袜一区| 超碰aⅴ人人做人人爽欧美| 免费看美女毛片| 高清免费毛片| 亚洲国产精品日韩av专区| 在线高清亚洲精品二区| 色综合狠狠操| 最新国产成人剧情在线播放| 青青青国产免费线在| 男女猛烈无遮挡午夜视频| 国产后式a一视频| 成人欧美日韩| 国产精品夜夜嗨视频免费视频| 波多野结衣一区二区三区88| 亚洲天堂视频在线免费观看| 国产va视频| 美女被躁出白浆视频播放| 亚洲人成网站日本片| 久夜色精品国产噜噜| 国产成人三级在线观看视频| 日韩欧美国产三级| 91在线精品麻豆欧美在线| 无码人妻免费| 色偷偷av男人的天堂不卡| 国产欧美性爱网| 国产97色在线| 在线观看免费AV网| 亚洲国产中文精品va在线播放| 久青草免费在线视频| 久久久久亚洲AV成人人电影软件| 国产特一级毛片| 91在线激情在线观看| 亚洲精品久综合蜜| 欧美日韩免费在线视频| 欧美日韩国产精品综合|