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

一種基于外推信息的慣導/陸基單站導航修正方案

2018-11-21 03:34:56左啟耀
導航定位與授時 2018年6期
關鍵詞:信息

王 勛,左啟耀,陳 亮,馬 超,任 鵬

(1.北京自動化控制設備研究所,北京 100074; 2.中國空間技術研究院,北京 100094)

0 引言

作為一種傳統的導航方式,陸基導航早在20世紀50、60年代已運用于軍民用航空領域, 如羅蘭系統、塔康系統等,具有不受時間、天候限制,測量定位精度不隨時間漂移,定位數據更新率高,接收設備簡單、價格低廉,可靠性高等諸多優點[1]。而隨著全球定位系統(Global Positioning System,GPS)技術的成功應用,陸基導航技術發展曾幾乎停滯。

近年來,各國對GPS拒止環境下導航技術的需求日益增加,陸基導航以其極高的可靠性、抗干擾性能,而被一些國內外學者重啟研究工作。Mauricio A等[2]提出了基于陸基網絡群的導航方法,該方法定位精度高、連續性好,但算法復雜,對導航終端硬件資源要求高,且至少需要4個陸基站同時運行,成本較高。Li等[3]用數學布站的方法將陸基導航測量值引入到彈道末段,但該方法僅在主動段使用陸基導航信息,在自由段進行彈道遞推,其導航誤差受制于彈道模型的準確性。Zhu等[4]通過觀測不同時刻的陸基單站測量結果, 經數據平滑等預處理后, 求出慣導定位和測速的修正參數,該方法僅利用一個陸基站,修護成本低,但陸基信號覆蓋區域有限,導彈飛出覆蓋區域后,存在導航誤差低或陸基信號缺失的問題,且布站固定,戰時陸基站容易被摧毀。

為此,本文利用陸基單站測距信息修正慣導,在載體發射前根據目標坐標等信息預先裝訂方案軌跡信息,或者在發射后對載體進行實時外推,將此外推軌跡與預先裝訂的方案軌跡進行比較,計算出軌跡偏差,根據偏差控制彈上的修正機構進行距離或方向修正。此外,某些載體(以彈道導彈為例)的落點精度主要取決于主動段關機點的速度誤差[5],載體縱向誤差系數遠大于橫向誤差系數,即主動段關機后同樣大小的速度誤差引起的縱向誤差較橫向誤差大,且由于主動段過載主要在射面內,導致關機時射面內的速度誤差較大,所以載體的縱向落點偏差一般大于橫向偏差。這時可結合陸基單站測距信息,進一步減小縱向落點偏差,提高載體速度和方向修正精度。

1 陸基測距原理

陸基單站系統測距原理是根據電磁波恒速直線傳播的特點,采用詢問-應答方式,通過測定電磁波發射點到接收點的傳播時間進行測距。在陸基導航系統中,測距機安裝在載體上,應答站設置在地面固定位置。如圖1所示,其測距工作原理為:載體上的測距機發出成對的問詢脈沖,應答站接收到脈沖后,經過一個固定的時延ts后,發出成對的應答脈沖。應答信號測距機接收到后,將發出問詢信號和收到應答信號之間所經過的時間tu減去固定時延ts,便可算出載體和應答站的距離r。

根據定位原理,距離計算公式如下:

(1)

式中,tu為載體上的測量設備測得的發射和接收信號之間的時間間隔;ts為地面應答器的固定延時;c為光速。

2 基于主動段樣條濾波算法的軌跡信息獲取方法

主動段運動模型是確定載體(如彈道導彈)軌跡外推參數的關鍵[6]。由于主動段受力情況十分復雜,主動段的推力和控制力通常不能精確已知,變質量過程也很難準確描述,導致載體的運動模型非常復雜,難以建立精確的運動模型。

對載體主動段建模需考慮影響軌跡的主要因素:重力、推力、關機時間、質量變化率、攻角、氣動阻力、自旋以及地球自轉偏向力等。這些因素均體現在載體的運動方程中:

∑F(t)=m(t)a(t)

(2)

式中,∑F(t)為作用于載體的合力;m(t)、a(t)分別為載體的實時質量和加速度。

給定初始條件,結合已知載體的作用力,求解方程式(2)能夠得到載體在飛行過程中任意時刻的位置、速度、加速度等信息。而實際上很難獲得上述幾個參數的完全解,所以在實際建模中必須進行一定的折衷,采用更加簡單有效的方法。考慮到載體軌跡在時序上是相關的,其軌道是滿足一定光滑特性的空間曲線,可表示為時間的函數f(t)。此外,參數化的目標是用一組形式簡單且容易計算的函數來高精度的逼近f(t),將軌跡估計問題轉化為這些函數的參數估計問題,能夠達到簡化計算和提高精度的目的。

載體主動段內軌跡可由3次樣條函數表示[9],為了避免狀態耦合導致的濾波計算量大、估計精度下降以及實時性差等問題,考慮利用解耦濾波器在一定范圍內進行狀態解耦,并可以根據解耦模型進行濾波估計。

3 單站測距信息修正慣導方法

修正方法僅利用一個地面站在不同時刻的測量結果,經數據平滑等預處理后,求出慣導定位的修正系數;建立精確的定位模型,解算出每段時間間隔的慣導修正系數,對慣導進行修正[9]。該方案實現簡單,但需要連續觀測,實時性較差。且對于高動態載體,連續測量時間過短會增加觀測量的相關性,引起觀測系數矩陣發生奇異,從而導致估計錯誤或發散。而考慮到慣導隨時間的漂移主要集中在距離通道,所以從修正距離通道出發,利用地面站跟蹤測量的測距數據,精確求出慣導在任一時刻tj的位置(Xj,Yj,Zj)的修正值(δXj,δYj,δZj)。

3.1 建立測距誤差方程

設地面站P的精確坐標為(XP,YP,ZP),慣導在tj時刻的測定結果為(Xj,Yj,Zj);tj時刻陸基導航系統對載體的測距結果為Rj。陸基導航設備測量的地面站與目標之間的距離與實際距離的誤差Vj,可表示為:

(3)

(4)

設tj時刻,慣導測定的搭載Mj坐標為(Xj,Yj,Zj),則(Xj,Yj,Zj)可表示為

(5)

(6)

由此可知,觀測方程的所有修正參數均為(δX0,δY0,δZ0),所以式(3)寫為

Vj=

(7)

將式(7)進行泰勒級數展開,線性化后可得

(8)

式中:

j=0,1,2,…,n

將式(8)寫成:

Vj=ljδX0+mjδY0+njδZ0-Lj

(9)

式(9)即為測距觀測的誤差方程。

3.2 誤差修正參數解算

根據誤差方程式(9),可采用最小二乘法解算修正參數(δX0,δY0,δZ0)。設觀測時刻為t0,t1,…,tn,共n+1個采樣時刻,可組成n+1個觀測誤差方程,表示為矩陣形式為:

V=AδX-L

(10)

式中:

V=[V0,V1,…,Vn]T,δX=[δX0,δY0,δZ0]T,

根據最小二乘方法,可組成法方程,即

ATAδX-ATL=0

(11)

求解法方程式(11),可得

δX=(ATA)-1ATL

(12)

式(12)的計算結果即為慣導修正參數的結果。

3.3 實時偏差求解與補償

載體在主動段關機后轉入自由飛行段,自由飛行段所處空間的大氣非常稀薄,載體所受的空氣動力近似為零,此時載體自由飛行段軌跡可看作橢圓軌道[8]。設載體所處的位置M與目標T的地心矢徑分別為rM、rT,矢徑間地心角為β,則經過載體位置M與目標T的橢圓軌道簇如圖2所示。根據等射程線原理,若位置M點的速度矢量的端點在雙曲線AB上,載體即可命中目標T,橢圓軌跡方程可描述為:

(13)

式中,P為橢圓的半通徑,e為橢圓的偏心率,ξM為軌跡上某點的遠地點角,P、e、ξM均為待定參數;β為矢徑rM、rT之間的地心角。

利用等射程線的特性,沿著等射程線的方向改變速度不會引起載體射程的改變,但等射程線法線方向速度的變化將會改變載體射程。結合單站測距系統、慣導信息及方案軌跡信息,可考慮采用距離測量信息精確估計載體等射程線法線方向的速度,以達到減小載體縱向射程的目的。

4 基于外推信息的慣導/單站測距修正方案

4.1 修正方案原理及實現步驟

4.1.1 方案原理

為了準確地估算飛行軌跡,以便依據軌跡偏差量對其實施軌跡修正,結合軌跡運動模型,應用主動段樣條濾波算法,對地面測控站測量得到的一段飛行軌跡參數進行濾波,進而外推軌跡實時參數。陸基測距系統融合SINS導航信息進行誤差修正,修正結果送入狀態估計器求解法向速度,進而得到實時軌跡的待修正量,從而執行載體速度和方向修正。圖3所示為所提出的基于外推信息的慣導/單站測距修正方案框圖。

圖3中,方案由軌跡探測解算模塊、陸基單站測距系統和補償模塊3個部分構成。軌跡探測解算模塊載入目標坐標信息,構建主動段運動樣條模型,采用解耦濾波器在一定的精度范圍內進行狀態解耦,利用地面測控信息,并通過遞推濾波解算軌跡參數;陸基單站測距系統在工作區域內,一方面結合修正SINS模塊提供的位置精度閾值選擇切換可用陸基站,另一方面為SINS修正模塊提供站點位置信息和一維測距數據,同時,SINS輸出姿態信息進一步反饋修正陸基窄波束天線指向,以提高陸基信號收發信干比。補償模塊利用陸基單站測量地面站至載體的距離ρ,并融合修正后的慣導信息ωe、pg、Vg進行多觀測器的狀態估計,得到相對精確的徑向速度VL,進而得出等射程線法向速度Vd,并結合探測解算遞推的軌跡偏差完成縱向速度誤差和方向的修正。

4.1.2 實現步驟

根據上述分析,基于外推信息的陸基單站/慣導測距修正方案的執行步驟可歸結如下:

Step1:軌跡信息遞推解算

利用主動段樣條濾波算法模型,在某一確定的時間間隔內,載體主動段加速度的變化率為常數的假定,采用分段多項式表示載體主動段運動。設定狀態變量為位置、速度和加速度,使用解耦模型解耦后,對其進行單通道獨立濾波,并且在狀態向量中引入調節項,以保證其具有優越的機動跟蹤性能。以x方向為例,樣條濾波的狀態方程為:

(14)

同樣,樣條濾波需要的地面測控站測量得到的觀測方程為:

(15)

Step2:陸基單站測距信息慣導修正

利用無線電精確單站測距信息,對一定區域內的載體進行慣導誤差修正,推導出單站測距誤差方程,并用最小二乘法解算出慣導修正數據?,F設一段時間內的相對距離誤差觀測值為向量V=[V0,V1,…,Vn]T,其中,Vj,j=0,1,…,n為測量的地面站與目標之間的距離與實際距離的誤差序列,則n+1個觀測誤差方程構成如下向量形式:

V=AδX-L

(16)

式中:

V=[V0,V1,…,Vn]T

δX=[δX0,δY0,δZ0]T

L=[L0,L1,…,Ln]T

令三維矩陣B=ATA,根據最小二乘理論可知,B-1矩陣非奇異是保證數據有效收斂的必要條件[10]。此外,工程中的測距誤差、修正時間間隔及載體與地面站的相對方位均會影響算法求解的收斂性。由此,采用合理的修正時間間隔、合理的布局方案及高精度的測距設備是算法收斂的前提,利用陸基單站測距信息輔助修正慣導,以達到補償慣導誤差的目的。

Step3:慣導修正參數輔助徑向速度求解

利用陸基單站測距系統測量得到的地面站至載體的距離信息,融合搭載慣導系統輸出的導航信息得到的徑向速度。而為了計算方便,將測距信息與導航解算轉至同一發射系下,發射系下慣導修正后輸出位置Pg和速度Vg,而加速度轉換為:

(17)

(18)

式中:

將式(18)兩邊求一階導數,并轉至發射坐標系下可得加速度為

(19)

利用單站測距信息得到的地面站距載體的距離,融合慣導系統輸出的導航信息進行多觀測器的狀態估計,可以得到高精度的徑向速度VL。

Step4:等射程線法向速度估計

根據圖3所示的橢圓軌跡簇,則軌跡方程可以描述為[14-15]

(20)

將式(20)代入式(13)可得

(21)

式中,K為動量矩矢量的模,且有K=rMVM·cosθH。

整理可得:

(22)

又橢圓的半通徑P為

P=K2/fM

(23)

式中,f為萬有引力常數,M為地球質量。

由式(22)和式(23)不難得到等射程線方程為:

(24)

將式(24)極坐標轉到軌道直角坐標系,令

(25)

利用式(25)的轉換關系可以得到當地坐標系下等射程線方程為:

(26)

整理可得:

(27)

式(27)兩邊同時對u求導數,可得:

(28)

由式(28)可得等射程線上任一點的斜率為:

(29)

那么等射程線法線方向為:

(30)

則等射程線法線方向單位矢量轉到發射坐標系中可表示為:

(31)

Step5:實時修正量求解與偏差補償

(32)

(33)

等射程線法線方法的速度Vd為:

(34)

利用速度增量的計算公式計算發射慣性系下的速度增量Vr,將此增量轉換到發射坐標系為

(35)

將速度增量投影至等射程線法線方向為

(36)

利用單站測距信息融合慣導修正信息,估計等射程線法線方向的速度和該方向的速度增量,便可利用此增量進行一步修正,以修正縱向落點偏差,達到提高命中精度的目的。

5 仿真驗證

為了驗證基于外推信息的慣導/陸基單站測距修正方案,設計了修正方案的仿真平臺,并對其進行仿真驗證。利用該平臺,通過與純慣導系統進行導航性能對比,首先在定位精度、快速性等方面對修正方案的導航性能進行研究,依次從陸基/慣導位置修正精度、速度偏差及速度修正精度等3個方面進行分析與驗證。

5.1 仿真條件設定

(1)慣性元件仿真參數

? 初始對準誤差:方位失準角誤差為10′(3σ),水平調平誤差為2′(3σ);

? 陀螺儀:常值漂移誤差為0.03(°)/h(3σ),各個方向安裝偏差均為15″(3σ),驅動白噪聲均方差為0.02(°)/h;

? 加速度計:零偏穩定性5×10-5g(3σ),標度因數穩定性:5×10-5(3σ),驅動白噪聲均方差為10-5g。

? 數據輸出頻率:慣性器件的數據輸出頻率為100Hz。

(2)陸基系統仿真參數

? 測距系統:測距接收機標準偏差為20.0m(1σ),測距應答站標準差為20.0m(1σ);總體均方差為28.3m(1σ);

? 數據輸出頻率:陸基導航器件的數據更新率為15Hz;

(3)地面測控站仿真參數

地面測控站坐標為(117.5°E,39.3°N),雷達測量距離標準差為10.0m(1σ),測量仰角和方位角噪聲標準差為1.0°(1σ),采樣時間為1.0s,雷達探測時間段為載體發射后2s~160s。

5.2 仿真驗證平臺設計

首先設計載體仿真軌跡,如圖4和表1所示。載體發射點坐標為(116.0°E,40.0°N),向正東發射,發射仰角為90°(垂直發射),射程為2217.0km,軌跡頂點高度為371.4km,飛行時間為501.0s。前20s為垂直上升段,60s主動段轉彎結束,160s發動機關機,501.0s載體落地,落地坐標為(141.2°E,40.0°N)。陸基單站測距系統工作段設定在載體上升段,工作時間為發射后62s~160s,即SINS/陸基組合修正工作段,一維測距關機點速度修正時刻為關機點,即160s。SINS系統單獨工作時的位置誤差如圖5所示。

表1 模擬的軌跡參數

圖5所示為組合導航工作時段單獨SINS的位置誤差,在初始設定的加速度計、陀螺儀誤差模型和初始對準誤差的條件下,由于慣性元件(加速度計、陀螺儀)的誤差以及慣導系統自身的誤差發散特性,SINS單獨工作時,位置誤差有隨時間的增加而逐漸積累的趨勢。根據仿真數據,在上升段x方向位置誤差約為28m,y方向位置誤差約為55m,z方向位置誤差約為30m。

圖6和圖7所示分別為軌跡、地面站相對布局圖和軌跡、地球相對位置關系圖,載體發射點坐標為(116.0°E,40.0°N),地面站坐標為(117.5°E,39.3°N),載體真實軌跡、組合濾波生成軌跡和純SINS解算軌跡基本重合,紅色標記軌跡為62s~150s的組合導航軌跡。該投影圖能夠直觀地反映出組合導航前后軌跡對比效果和載體相對于地面站A的實時相對位置。

5.3 仿真結果及分析

5.3.1 樣條濾波軌跡信息結果與分析

在主動段,其3個坐標方向的運動軌跡是滿足一定光滑特性的空間曲線,3個方向的加加速度大小變化很小,幾乎是一個常數,x方向大小為-0.069(m/s2)/s,y方向大小為0.15(m/s2)/s,z方向大小為0.21(m/s2)/s。因而可以用樣條函數分別描述載體主動段在3個方向上的運動,相對能比較真實地反映實際運動。圖8~圖10所示為基于樣條濾波算法的x、y、z這3個方向濾波位置誤差變化曲線圖。

從圖8~圖10可以看出,主動段樣條濾波算法精度較高,這是由于樣條濾波算法是一種基于動態多項式模型的算法,多項式隨軌跡特性變化,且軌跡平穩時服從嚴格的均加加速運動,模型的保持特性提高了隨機誤差的抑制能力。此外,這種算法拋棄了參數回歸思想,通過遞推濾波解算軌跡,極大地降低了計算負擔,消除了解算的滯后性,從而使軌跡解算具有較高的實時、動態性能。

5.3.2 單站測距信息修正捷聯慣導結果與分析

在陸基工作區域內,根據修正SINS模塊提供的位置精度閾值選擇切換可用陸基站,并利用其提供的站點位置信息和一維測距數據,分別得到單站測距信息修正捷聯慣導的導航參數結果。圖11~圖13所示為修正方案x、y、z這3個方向的位置誤差結果,圖14~圖16所示為修正方案位置誤差與純捷聯慣導系統輸出位置誤差的對比曲線。

從圖11~圖13可以看出,修正方案計算的位置誤差均在14m以內,而由于慣導系統自身的誤差發散特性,SINS單獨工作時位置誤差有隨時間逐漸積累的趨勢。其中x方向和y方向誤差較大,x方向誤差發散至500m,y方向誤差發散至800m,z方向誤差發散至300m。而圖11~圖13中修正方案輸出的位置誤差相對較小,x、y方向誤差最大值分別為6.0m、6.0m,z方向的誤差最大值為13.0m,且在修正后期由于慣導誤差積累以及信息權值無法準確預測,所以修正系統的誤差仍存在較大波動,但整體較純慣導系統有了較大的提高。

慣導/單站測距修正系統要求陸基系統提供偽距等直接量測數據以及由信號解碼得到的地面站點信息,涉及搭載接收機內部搜索回路的參數設置和輸出,技術實現比較復雜。由于陸基導航系統誤差源較多,而建模過程中難以對偽距誤差實現精確補償。對于高動態載體,其加速度、姿態等參數在短時間內變化劇烈,導航系統需要很高的數據更新率,而由于搭載接收機數據輸出速率的限制,使得導航系統數據更新率尚難以滿足實時性的要求。因此,提高搭載接收機的更新速率,完善組合系統對搜索回路的輔助是慣導/單站測距修正系統在工程應用的關鍵。

5.3.3 實時偏差求解與補償結果與分析

測距系統的測距誤差取高斯正態分布模型,其誤差的方差為2m(1σ),用蒙特卡羅方法分別對修正方案和純慣導兩種情況進行50次仿真,利用陸基單站測量地面站至載體的距離ρ,并融合慣導信息ωe、pg、Vg進行多觀測器的狀態估計后,進而得到徑向速度估計誤差如圖17所示;進一步,結合相對精確的徑向速度VL,從而得出等射程線法向速度Vd,并依據探測解算遞推的軌跡偏差完成縱向速度誤差和方向的修正過程,得到如表2所示的修正方案與純慣導估計精度對比。

表2 修正方案與純慣導精度對比

從圖17可以看出,采用慣導信息和測距信息融合后,可以得到較高精度的徑向速度估計,精度在±0.20m/s以內。而從表2可以看出,融合慣導信息和測距信息后,等射程線法線方向的速度精度可提高約3.5倍,則對應的落點縱向誤差將減小至原來的2/7。由此可見,提出的基于陸基單站的測距信息,通過測量發射點至載體的距離,并結合慣導輸出導航信息,可以有效提高徑向速度估計精度,從而減小縱向誤差,提高載體的縱向精度。此外,這種方法對地面站布局無特殊要求,甚至可以將單個基站安裝于載車上隨車機動。

6 結論

本文根據一段飛行軌跡測量數據,通過融合單個陸基無線電測距信息和樣條濾波軌跡外推預測信息,對一定區域內的載體進行慣導誤差修正,并融合慣導修正后信息進行多狀態觀測估計,解算出影響載體縱向誤差的等射程線法向速度,再度結合外推預測軌跡偏差完成縱向速度誤差和方向的修正過程,在實現慣導修正的同時,進一步提高了載體的縱向精度。該修正方案的優點總結如下:

1)通過融合樣條濾波預測信息,減小了載體高動態特性引起的短時方位系數矩陣奇異,提高了單站測距信息修正慣導數據的實時性。

2)利用修正后的SINS導航信息參與狀態估計器求解法向速度,可以實時解算得到更為精確的軌跡參數待修正量,通過執行載體速度和方向修正補償后,載體的縱向誤差減小至純慣導時的2/7。

3)利用SINS修正后的姿態信息調整陸基窄波束天線指向,保證了陸基單站的測距精度,陸基地面站無需特殊布站,甚至可將基站安置于載車上,隨車機動,且不影響修正精度。

猜你喜歡
信息
訂閱信息
中華手工(2017年2期)2017-06-06 23:00:31
展會信息
中外會展(2014年4期)2014-11-27 07:46:46
信息超市
大眾創業(2009年10期)2009-10-08 04:52:00
展會信息
展會信息
展會信息
展會信息
展會信息
信息
建筑創作(2001年3期)2001-08-22 18:48:14
健康信息
祝您健康(1987年3期)1987-12-30 09:52:32
主站蜘蛛池模板: 成人欧美在线观看| 日韩欧美国产三级| 激情视频综合网| 国产人碰人摸人爱免费视频| 91视频区| 亚洲国产清纯| 2021国产精品自产拍在线观看| 无码AV高清毛片中国一级毛片| 国产99精品久久| 色网在线视频| 在线欧美一区| 国内精品久久九九国产精品| 国产永久在线观看| 狠狠久久综合伊人不卡| 欧美综合中文字幕久久| 97精品伊人久久大香线蕉| 久久香蕉国产线看观看精品蕉| 久久毛片网| 国产欧美日韩另类精彩视频| 色综合久久88色综合天天提莫| 国产精品亚欧美一区二区三区| 国内精品视频| 色综合中文| 91 九色视频丝袜| 制服丝袜无码每日更新| 91久久青青草原精品国产| 国产真实二区一区在线亚洲| 亚洲欧美在线综合图区| 欧美在线观看不卡| 日本精品αv中文字幕| 成人日韩欧美| 一区二区日韩国产精久久| 久久香蕉欧美精品| 中文字幕无码中文字幕有码在线| 日韩毛片在线播放| 国产精品私拍在线爆乳| 成人在线观看一区| 午夜视频免费试看| 欧美成人午夜影院| 2020国产精品视频| 一级片免费网站| 国产成人精品免费视频大全五级| 香蕉视频在线精品| 中文字幕2区| 国产欧美日韩视频怡春院| 亚洲成人在线网| 久久人与动人物A级毛片| 人妻少妇乱子伦精品无码专区毛片| 欧美日韩国产精品综合 | www.youjizz.com久久| 无码精品福利一区二区三区| 国产精品久线在线观看| 香蕉eeww99国产在线观看| www中文字幕在线观看| 国产超碰一区二区三区| 亚洲国产中文在线二区三区免| 日本久久网站| 无码高清专区| 激情网址在线观看| 国产网友愉拍精品视频| 国产性爱网站| 91网红精品在线观看| 亚洲视频在线观看免费视频| 亚洲黄色高清| 谁有在线观看日韩亚洲最新视频| 精品成人一区二区三区电影 | 免费在线视频a| 成年女人a毛片免费视频| 国产精品久久国产精麻豆99网站| 视频一本大道香蕉久在线播放| 中文字幕第4页| 国产一级裸网站| 免费播放毛片| 国产精选小视频在线观看| 永久成人无码激情视频免费| 无码aaa视频| 色呦呦手机在线精品| 一本久道热中字伊人| 亚洲人人视频| 欧美va亚洲va香蕉在线| 国产人妖视频一区在线观看| 午夜性刺激在线观看免费|