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

捷聯慣導系統全狀態量現場標定方法

2023-08-09 14:49:10奔粵陽劉利強
系統工程與電子技術 2023年8期
關鍵詞:系統

張 睿, 奔粵陽, 劉利強, 王 坤, 侯 靚, 邱 天

(哈爾濱工程大學智能科學與工程學院, 黑龍江 哈爾濱 150001)

0 引 言

捷聯慣導系統標定技術能確定慣性器件誤差模型的相關參數,在提高導航精度方面有重要作用[1]。根據標定場所的不同,標定技術分為實驗室標定和現場標定。實驗室標定是較早被開展研究的標定技術,其相關技術已經發展得較為成熟;而現場標定方法的研究起步較晚,目前仍然在不斷發展的過程中[2]。

標定方法分為分立式標定和系統級標定,目前常用系統級標定方法進行現場標定,進行系統級標定前需要確定標定的誤差參數是可觀測的,一般使用可觀測性分析的相關方法對狀態量是否可觀測進行衡量。常用的可觀測性分析方法主要是文獻[3]和文獻[4]分別論述的分段線性定常系統(piece-wise constant system, PWCS)和奇異值分解(singular value decomposition, SVD)。PWCS雖然可以直觀地得到可觀測狀態量的數目,但是不能求出每個狀態量具體的可觀測度。而SVD雖然能量化計算狀態量的可觀測度,但是目前在進行狀態估計時經常把噪聲干擾假設為已知協方差的高斯白噪聲。盡管這樣處理便于進行狀態估計,但是這個嚴苛的假設在實際工程中很難滿足。文獻[5]提出了一種基于狀態量正交基的可觀測度計算方法,該方法無需限制噪聲干擾,就可以計算出每個狀態量的可觀測度,具有一定的工程意義。

現場標定的目的是在復雜的外場環境中提高標定的精度,此外還需要保證觀測誤差模型是全狀態量模型時,還能實現快速觀測。系統級標定方法盡管通過使用導航誤差作為觀測量便于進行現場標定,但是既不能提高狀態量的可觀測度,也不能提高標定效果。文獻[6]提出的設計轉位方案解決了以上兩個問題,不僅可以提高捷聯慣導系統的可觀測度,還能使標定效果更好。文獻[7]方法的缺陷是只能對常用的誤差參數進行標定,隨著對系統級標定精度的要求越來越高,需要激勵出更多的系統誤差。文獻[8]提出基于低成本轉臺的高精度光纖陀螺混合式標定方法,該方法標定了27維誤差參數,但只實現了較高精度的實驗室標定,沒有實現場外標定。

在全狀態量現場標定的過程中,選取的誤差參數包括陀螺漂移、加速度計零位偏移、陀螺儀標度因數誤差、加速度計標度因數誤差、陀螺儀安裝誤差、加速度計安裝誤差,以及加速度計二次非線性誤差。本文先是根據捷聯慣導系統的理論建立慣性器件誤差模型和標定誤差方程。之后在線性時變系統中,提出一種基于狀態量正交基的可觀測度計算方法。最后在建立的捷聯慣導誤差模型基礎上,根據計算可觀測度數值的變化設計出一種現場標定的轉位方案。理論分析和實驗驗證表明,該轉位方案不僅提高了系統狀態量的可觀測度,還能充分激勵標定誤差參數,獲得了更好的現場標定效果,是一種實用的現場標定方法。

1 捷聯慣導系統誤差方程

根據可觀測度方面的文獻綜述[9]和捷聯慣導系統現場標定的相關理論知識[1],建立陀螺儀和加速度計這兩個慣性器件全狀態量誤差模型和對應的標定誤差方程。

1.1 慣性器件誤差模型建立

(1)

類似地,在b系中建立加速度計的全狀態量誤差模型[11]為

(2)

1.2 標定誤差方程

一般地,捷聯慣導系統的速度誤差方程和姿態誤差方程可以寫為[13]如下形式:

(3)

(4)

(5)

(6)

(7)

(8)

(9)

式中:g表示當地地球重力加速度;ωie是地球自轉的角速度;L是當地的地理緯度。

在式(5)~式(9)的條件下,忽略牽連加速度,將式(3)、式(4)展開整理得到姿態誤差方程和速度誤差方程[14]:

(10)

(11)

2 線性時變系統狀態量可觀測度

在線性時變系統中,系統狀態的估計精度由3個因素決定:噪聲干擾、狀態估計方法和可觀測度。噪聲干擾是一個不可控的隨機過程,狀態估計采用卡爾曼濾波方法。為了提高狀態估計精度,就要研究可觀測度計算方法。

以往關于可觀測度的研究只是定性分析狀態量是否可觀,有些方法盡管求出了可觀測度,但是還要假設噪聲是已知協方差的高斯白噪聲。文獻[5]提出了一種基于狀態量正交基的可觀測度計算方法,不僅能夠定量地計算每一個狀態量的可觀測度,還沒有噪聲干擾假設的限制。

2.1 基于狀態量正交基的可觀測度

一般地,線性時變系統離散化后得到的狀態空間模型為

(12)

式中:Xk∈Rn×1是一個n維的狀態向量;Φk/k-1∈Rn×n是更新狀態量的狀態一步轉移矩陣;Zk∈Rm×1是量測向量;Hk∈Rm×n是量測矩陣;ΔZk∈Rm×1是由模型誤差、系統噪聲和數值計算共同造成的干擾,這個干擾是不可預測也不可控的,在基于狀態量正交基可觀測度計算方法中不對這個干擾作任何限制。

對式(12)進行進一步分析,測量序列、干擾和初始狀態的關系為

(13)

將式 (13)展開處理得:

(14)

以計算狀態量x1的可觀測度為例,將式(14)進行處理得:

(15)

(16)

又因為rank[q2,q3,…,qn]

量測的理想值不存在噪聲,則式(15)的量測值是理想值時該式可改寫為

(17)

將式(17)處理可得:

(18)

(19)

使用歐幾里得范數算子進行處理,得:

(20)

(21)

同時需要滿足ξqi=0,i=2,3,…,n成立。

2.2 線性時變系統可觀測度的計算方法

第2.1節介紹了基于狀態量正交基的可觀測度計算方法,下面使用此方法計算狀態量的可觀測度,仍以狀態量x1為例。

根據線性代數的相關基礎知識,可以求出ξqi=0,i=2,3,…,n的通解為

(22)

(23)

式中:b1,b2,…,bσ是組合系數,它們的平方和為1。

為了求出目標函數maxJ1和對應的組合系數,使用拉格朗日乘數法求解。拉格朗日函數為

(24)

其中,λ是拉格朗日算子。實現目標函數的必要條件是

(25)

將式(25)帶入式(24)計算得到:

(26)

此時,可以得到可觀測度的計算公式為

(27)

接著根據式(27)計算出轉位過程中各誤差參數的可觀測度數值,通過分析在不同轉動過程中數值的變化,設計出增大可觀測度的標定轉位方案,使得該系統可以通過系統級標定激勵誤差參數[15]。

3 標定轉位方案的分析與設計

在第1節建立的慣性器件誤差模型和慣性導航理論的基礎上,結合第2節提出的可觀測度計算方法,確定設計轉位方案的方法并且設計出新的標定方案,對比經過該方案轉動后誤差參數的可觀測度和未經轉動的誤差參數可觀測度,確定轉位方案的有效性[16]。

3.1 標定原則的確定

若在靜止狀態進行系統級標定,那么既不能提高狀態量的可觀測度,又無法提高標定效果,而通過轉動慣性器件可以改變式(13)中的Q矩陣,進而改變可觀測度,所以可以通過設計轉位方案來提高各狀態量的可觀測度[17]。

由第1節建立的慣性器件誤差模型可知,需要標定的誤差參數包括陀螺漂移、加速度計零偏、陀螺儀和加速度計的標度因數和安裝誤差、加速度計的二次非線性誤差。為了實現所有誤差參數的系統級標定,需要通過轉動提高各誤差參數的可觀測度[18]。在原有轉位方案的基礎上,使用第2節推導出的可觀測度計算公式計算誤差參數的可觀測度變化情況,確定標定原則,根據標定原則設計轉位方案。

根據第1節建立的慣性器件誤差模型確定所有狀態量,定義式(12)的狀態空間方程中對應的兩個向量分別為

(28)

Z=[δvT]T

(29)

在式(28)和式(29)中,確定的狀態量依次為:失準角、速度誤差、陀螺漂移、加速度計零偏、陀螺儀標度因數和安裝誤差組成矩陣的列向量、加速度計標度因數和安裝誤差組成矩陣的列向量、加速度計二次非線性誤差。量測量則只有速度誤差。

根據對應狀態量之間的關系,可以分別寫出處理后的一步轉移矩陣和量測矩陣為

(30)

(31)

得到一步轉移矩陣和量測矩陣后,可以根據第2.2節的計算方法計算各狀態量的可觀測度,進一步設計轉動方案。

根據文獻[19]的分析,可以把轉動方式分為兩類,分別是繞同一個坐標軸轉動和依次繞3個坐標軸轉動。需要驗證的是:繞同一個坐標軸轉動幾次可觀測度數值足夠大,繞3個軸轉動次序如何確定,才能夠使可觀測度最大[20]。下面使用第2節推導的式(27)分別計算兩種轉動情況下各狀態量可觀測度變化情況,分析并總結可觀測度變化規律。其中表1是依次繞x軸轉動3次過程中每個狀態量的可觀測度。

表1 繞x軸轉動3次各狀態量的可觀測度

根據表1數據分析可知:繞同一坐標軸連續轉動3次后各狀態量可觀測度的數值足夠大,此外深入分析可知,相比加速度計相關狀態量的可觀測度,陀螺相關的可觀測度受轉軸的轉動影響更大,通過這種轉動方式可以把加速度計的誤差和陀螺儀的誤差分離[21]。

接下來是對第2種情況的分析,繞3個坐標軸轉動順序可以分為6種情況,由于轉動情況較多,分別在表2和表3中對6種轉動情況下誤差參數的可觀測度進行了計算,之后根據計算結果確定轉動方式。

表2 3種轉動次序誤差參數的可觀測度

表3 其余3種轉動次序誤差參數的可觀測度

由于之前分析的結論是:繞同一坐標軸連續轉動對加速度計誤差參數的影響較小,因此三軸轉動次序需要對加速度計誤差參數的可觀測度有較大的提升[22]。由表2和表3計算的不同轉動方式下的誤差參數可觀測度數值中,顯然轉動方式為zxy和zyx時和加速度計相關的誤差參數可觀測度增加很小,再結合以往的標定轉位方法以及所有誤差參數的可觀測度進行綜合考慮,最終選擇xyz轉動方式[23]。

隨著以上兩種轉動方式各狀態量的可觀測度計算,可以總結出轉位方案設計的規律:繞同一個坐標軸連續轉動只會較大程度地改善陀螺相關誤差參數的可觀測度,繞三軸轉動則會同時對所有誤差參數的可觀測度產生影響。通過這個規律,可以把陀螺和加速度計的誤差參數解耦[24]。

由以上計算結果可以確定標定過程有如下原則:通過轉位方案的設計,建立預標定的單一器件參數誤差和系統狀態參數誤差之間的關系,使得式(13)中的Q矩陣滿秩[25],可以通過式(27)計算狀態參數的可觀測度并進行評估,通過提高誤差參數的可觀測度來提升系統級標定的效果。轉位方案的編排是為了提高可觀測度,以進一步提高標定效果[26]。

3.2 轉位方案的設計

根據標定誤差方程式(10)和式(11),若在系統級標定過程中忽略不可控的慣性器件隨機誤差,那么系統的速度誤差δvn主要由慣性器件的測量誤差產生,而這個測量誤差則受慣性器件的誤差模型式(1)、式(2)中的各項誤差參數影響[27]。當慣性器件從一個位置翻轉到另一個位置,這段時間內靜基座導航的速度誤差就是由慣性器件的一部分誤差參數造成的,可認為慣性器件的位置轉動激勵出了它的部分誤差參數[28]。如果設計出完整的轉位方案,就能充分激勵出慣性器件的各項誤差參數,同時使各誤差參數之間解耦,通過提高誤差參數的濾波估計效果提高系統級標定的精度[29]。

此外,由于每個轉動過程中速度誤差δvn的數值較小,速度誤差和慣性器件的誤差參數可以近似認為是線性關系,并且這個線性關系隨時間不斷變化,因此把這個系統近似為線性時變系統[30]。為了提高線性時變系統的可觀測度,根據可觀測度計算公式(27)可知,狀態量x1的可觀測度大小主要和Q矩陣中x1對應第一列向量q1以及Q矩陣從第二列到最后一列的列向量標準化的零空間正交基p1,p2,…,pσ有關。由于p1,p2,…,pσ經過了標準化處理后對可觀測度數值影響較小,而p1,p2,…,pσ又是[q2,q3,…,qn]T的零空間正交基,改變p1,p2,…,pσ的數值又會影響狀態量x2,x3,…,xn的可觀測度,可以簡化問題為主要分析轉動過程對qi(i=1,2,…,n)向量的影響來設計轉位方案,將能夠使βi(i=1,2,…,n)變大的轉動過程進行合理的組合,最終設計出標定轉位方案[31]。根據理論分析流程,得到設計轉位方案的基本流程圖如圖1所示。

圖1 轉位方案設計流程圖Fig.1 Design flowchart of transposition scheme

根據第2節的計算結果和本節的設計轉位方案的方法,確定了轉動過程是繞x軸、y軸、z軸的轉動順序依次轉動3次。但是在實際的轉位過程中,需要進行適當的調整才有更好的標定效果。因此在原有19位置轉位方案的基礎上,根據文獻[17]的設計原理,以及轉動過程標定效果的不斷試驗,設計出一種新的23位置轉位方案(見表4)。該方案首先繞x軸正、反向各轉3次,之后繞y軸正向轉動一次后繞z軸正反向各轉3次,隨后再繞y軸正向轉動一次后繞x軸正反向各轉3次,最后繞y軸正向轉動一次再反向轉回原位。該方案是一種綜合考慮了轉動次序和轉動順序的可觀測度計算結果,結合了傳統的19位置轉動方法,同時結合了實際的標定測試結果而設計出的轉位方案。該方案充分激勵了慣性器件的各項誤差參數,還能夠提高系統各誤差參數的可觀測度。

表4 新設計的標定轉位方案

2) 以+oy為例,表示繞y軸正向旋轉90°。

得到理論適宜的轉位方案后,接下來需要通過真實的實驗結果驗證所設計的轉位方案是否合理。驗證的方法就是在實際的仿真環境中,對比傳統標定方案和新的標定方案的可觀測度數值以及標定效果,從兩個方面判斷新的轉位方案的合理性。

4 仿真實驗結果分析

實驗的目的是驗證轉位方案對狀態量可觀測度和標定的影響。若濾波后狀態量協方差快速收斂到零,說明這個狀態量波動小,可觀測度大。而驗證標定效果則是將狀態誤差估計值和設置的真實值相減,若標定結果收斂到零,則效果更好。為了直接比較轉位方案對系統級標定的影響,分別對無轉位和有轉位的標定結果進行對比。

實驗進行前,首先設置標定誤差參數的真實值如表5所示,接下來進行兩種轉位方案的對比實驗。

表5 標定真實值設定

4.1 原始轉位方案系統級標定仿真實驗結果

首先列舉原始使用的19位置系統級標定方案,如表6所示。

表6 標定轉位方案

傳統19位置轉動方案是先繞y軸正反向各轉3次,之后繞z軸正向轉動一次后繞x軸正反向各轉3次,之后繞z軸正向轉動兩次后反向轉動3次并轉回原位。相比新設計的轉位方案,傳統的轉位方案少了1次繞x軸的正反向轉動,轉動軸的轉動次序也不同,是主要的產生變化的因素。因此,接下來分別驗證兩種方案可觀測度的變化和標定效果的變化。

當標定轉動過程按照表6中的傳統轉位方案進行時,可觀測度的計算結果如表7所示。

表7 傳統轉位方案的可觀測度計算結果

為了驗證表7中的計算結果,進行對應傳統轉位方式的系統級標定實驗。標定結果和對應的協方差結果如圖2~圖11所示。

圖2 陀螺漂移標定Fig.2 Gyro drift calibration

圖3 加速度計零偏標定Fig.3 Accelerometer bias calibration

圖4 標度因數誤差標定Fig.4 Scale factor error calibration

圖5 安裝誤差標定Fig.5 Installation error calibration

圖2~圖6是傳統轉位方式下系統級標定各個誤差參數的實驗結果,圖7~圖11則是這些誤差參數的協方差估計結果。仿真結果結合表7的數據分析可知:可觀測度計算結果越大,協方差估計結果越能夠收斂到零,此時系統級標定的濾波估計結果和設置的真實值相減的結果越能夠收斂到零。傳統轉位方式下的多數狀態量可觀測度都較大,但是也有個別誤差參數的可觀測度較小,協方差矩陣也沒能完全收斂,同時仍然存在不能全狀態量可觀測的缺陷,因此傳統的轉位方案不能對全部狀態量的協方差取得較好的收斂效果。標定的結果相對無轉位方式的可觀測度有了較大提升,協方差矩陣也可更好的收斂,標定效果也有較好的改善,但是仍未能實現全狀態量的可觀測。

圖6 二次非線性誤差標定Fig.6 Quadratic nonlinear error calibration

圖7 陀螺漂移協方差Fig.7 Gyro drift covariance

圖8 加速度計零偏協方差Fig.8 Accelerometer bias covariance

圖10 安裝誤差協方差Fig.10 Installation error covariance

圖11 二次非線性誤差協方差Fig.11 Quadratic nonlinear error covariance

傳統轉位方法系統級標定的結果較無轉位的標定已經有了很大的提升,但是仍存在部分狀態量的可觀測度不大的問題,協方差矩陣也沒有能夠全部收斂,也沒有實現全狀態量的可觀測。接下來設置新的轉位方案并進行進一步的實驗,以驗證新的轉位方案是否能夠實現全狀態量的系統級標定。

4.2 新轉位方案系統級標定仿真實驗結果

當系統級標定過程中設置的是新設計的轉位方案時,計算各狀態量的可觀測度結果,如表8所示。

表8 轉位可觀測度計算結果

將表8的結果與表7進行對比,可以看到表8轉動方案的優勢是:對于傳統方案中可觀測度較高的誤差參數,例如陀螺漂移、陀螺的標度因數等狀態量,改進的新方法仍能夠保持較高的可觀測度,在保證較高的可觀測度的前提下甚至有所提升;對于傳統方案中可觀測度提高不多的誤差參數,例如加速度計零偏、加速度計的部分安裝誤差等狀態量,新方法對它們的可觀測度有明顯的提高,新方法對全狀態量可觀測度的數值有了很大改善。改進后的轉位方案進行系統級標定的實驗結果如圖12~圖21所示。

圖12 陀螺漂移標定Fig.12 Gyro drift calibration

圖14 標度因數誤差標定Fig.14 Scale factor error calibration

圖15 安裝誤差標定Fig.15 Installation error calibration

圖12~圖16是設置了新轉位方案的系統級標定的結果,由圖12~圖16可以得到的結論是:各誤差參數最終都穩定地收斂到零,這說明系統級標定的濾波估計結果能夠收斂到設置的真實值附近。更新轉位方案后的系統級標定結果對比傳統轉位方案的標定結果有了進一步的提高。

圖17~圖21是設置了新轉位方案的系統級標定的協方差估計結果。無論開始時間狀態量的協方差多么大,在結束系統級標定后,各狀態量的協方差都能夠收斂到零,說明最終狀態量濾波結果波動較小,標定結果穩定。新轉位方案的協方差估計效果也有明顯的提升。

圖17 陀螺漂移協方差Fig.17 Gyro drift covariance

圖18 加速度計零偏協方差Fig.18 Accelerometer bias covariance

圖19 標度因數誤差協方差Fig.19 Scale factor error covariance

圖20 安裝誤差協方差Fig.20 Installation error covariance

圖21 二次非線性誤差協方差Fig.21 Quadratic nonlinear error covariance

綜合對比分析傳統轉位方案和新轉位方案的系統級標定實驗的結果可知:改變轉動過程后,濾波估計狀態量的可觀測度整體能夠變得更大,對應的協方差也都能夠向零收斂,系統級標定結果和真實值相減的結果最終收斂到零。總結兩種轉位方案的系統級標定實驗的結果,驗證了根據捷聯慣導基本方程和基于狀態量正交基的可觀測度計算方法能夠設計出更加合理的轉位方案。該轉位方案不僅能夠改善各個誤差參數的可觀測度,使各狀態量協方差有更好的收斂效果,還能有效提高系統級標定的結果。

5 結 論

捷聯慣導系統的標定技術是慣導系統研究和應用的一項關鍵技術,但是靜基座無轉動條件下的系統級標定不能充分激勵各項誤差參數,需要設計出合理的轉位方案。此外,各狀態量并不完全可觀測,因此在系統級標定時無法通過濾波估計提高標定效果。

本文分析狀態量可觀測度的計算方法,提出了基于狀態量正交基的可觀測度計算方法以設計系統級標定的轉位方案。此方法根據可觀測度的計算思路,分析每次轉動過程中可觀測度如何變化,設計出最大化可觀測度同時還能激勵誤差參數的轉位方案。經過傳統轉位方案和新轉位方案的系統級標定實驗對比,所設計的轉位方案不僅能夠改善各狀態量的可觀測度,而且充分激勵了全狀態量誤差參數,得到了更好的標定結果。該方法適用于激勵出更多誤差參數、對結果要求較高的系統級標定具有一定的工程應用參考價值。

猜你喜歡
系統
Smartflower POP 一體式光伏系統
工業設計(2022年8期)2022-09-09 07:43:20
WJ-700無人機系統
ZC系列無人機遙感系統
北京測繪(2020年12期)2020-12-29 01:33:58
基于PowerPC+FPGA顯示系統
基于UG的發射箱自動化虛擬裝配系統開發
半沸制皂系統(下)
FAO系統特有功能分析及互聯互通探討
連通與提升系統的最后一塊拼圖 Audiolab 傲立 M-DAC mini
一德系統 德行天下
PLC在多段調速系統中的應用
主站蜘蛛池模板: 国产一区二区三区在线无码| 91亚洲国产视频| 亚洲香蕉伊综合在人在线| 伊人天堂网| 国产精品永久久久久| 人妻免费无码不卡视频| 欧美特级AAAAAA视频免费观看| 午夜色综合| 波多野结衣亚洲一区| 日本精品一在线观看视频| 国产福利一区二区在线观看| 亚洲人网站| 亚洲无码高清一区二区| 国产精品区视频中文字幕| 亚洲精品中文字幕无乱码| 无码综合天天久久综合网| 成人一级黄色毛片| 久久精品中文字幕少妇| 国产精品视频系列专区| 伊人久久大线影院首页| 免费jjzz在在线播放国产| Jizz国产色系免费| 青草视频在线观看国产| 在线欧美日韩| 久草网视频在线| 国产va视频| 欧美啪啪精品| 亚洲国产成人超福利久久精品| 国产精品无码久久久久AV| 四虎成人在线视频| 国产一区二区福利| 极品私人尤物在线精品首页| 国产精品网曝门免费视频| 精品免费在线视频| 91国内外精品自在线播放| 国产精品欧美日本韩免费一区二区三区不卡 | 九色在线观看视频| 免费a级毛片视频| AV老司机AV天堂| 亚洲侵犯无码网址在线观看| 99这里只有精品免费视频| 性做久久久久久久免费看| 国产精品一线天| 亚洲AV无码精品无码久久蜜桃| 国产香蕉国产精品偷在线观看| 一本色道久久88| 国产另类乱子伦精品免费女| 欧美a在线视频| 久久精品国产999大香线焦| 国产精品亚洲五月天高清| 免费在线国产一区二区三区精品 | 中文字幕人成人乱码亚洲电影| vvvv98国产成人综合青青| 国产91精品调教在线播放| 日韩精品亚洲人旧成在线| 日韩在线播放欧美字幕| 色噜噜狠狠狠综合曰曰曰| 丁香婷婷在线视频| 亚洲天堂啪啪| 婷婷六月综合网| 国产三级成人| 久久亚洲综合伊人| 国产一区二区三区精品欧美日韩| 国产自在线播放| 国产精品成人免费视频99| 欧美19综合中文字幕| 波多野结衣在线se| 狠狠色噜噜狠狠狠狠奇米777| 99精品国产高清一区二区| 国产成人综合日韩精品无码首页 | 国产屁屁影院| 久久精品国产亚洲麻豆| 国产麻豆精品久久一二三| 日韩专区欧美| 国产成a人片在线播放| 国产欧美专区在线观看| 亚洲天堂久久久| 国产熟睡乱子伦视频网站| 在线看免费无码av天堂的| 一本大道AV人久久综合| 午夜综合网| 日本五区在线不卡精品|