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

基于天波重構技術的eLORAN信號周期識別算法

2022-10-29 03:36:12劉時堯張首剛
電子與信息學報 2022年10期
關鍵詞:信號

劉時堯 華 宇 張首剛

①(中國科學院國家授時中心 西安 710600)

②(中國科學院精密導航定位與定時技術重點研究室 西安 710600)

③(中國科學院時間頻率基準重點研究室 西安 710600)

④(中國科學院大學 北京 100049)

1 引言

羅蘭(LOng RAnge Navigation, LORAN)C系統是美國海軍于20世紀50年代出于軍事目的而發展起來的陸基導航系統,是一種遠程雙曲線導航定位系統[1]。增強型羅蘭(enhanced LORAN, eLORAN)是由Loran-C系統演變而來,其采用Eurofix數據鏈技術進行三態脈位調制(Pulse Position Modulation,PPM)實現增強數據信息發播[2],具有作用距離遠、穩定性好、抗干擾能力強等優點。近年來,美國、韓國、英國、俄羅斯等多國均重新啟動了新型eLORAN系統的建設研究,而我國也于2017年啟動了“高精度地基授時系統”項目,旨在成為北斗衛星系統的可靠地基備份以完善我國的星地一體化高精度定位、導航、授時(Positioning, Navigation and Timing, PNT)系統。

周期識別是信號處理中的核心問題,其穩定性極易受到天波、連續波、交叉干擾以及各種噪聲的影響,可能引起跳周(1載波周期=10 μs)。另外,對于添加了數據調制的eLORAN信號,很容易因為將天波信號誤捕獲為地波而造成解碼誤碼率增大,最終導致授時錯誤及導航定位的極大誤差。天波干擾指的是經過電離層反射的天波輻射信號與地面波間接收時延差小于信號的有效時長時對地面波造成的重疊效應[3]。根據國際海事無線電技術委員會(Radio Technical Commission for Maritime, RTCM)發布的長波接收機最低性能標準,天波一般位于地波信號后的37.5~1500 μs[4]。為了保證PNT精度,eLORAN接收機要求利用地波信號進行周期識別,因此一般將地波脈沖組第1個脈沖的第3載波周期正向過0點(30 μs處)作為標準過0點。

近年來出現了不少周期識別算法的研究[5—9],包括優化包絡相關法、多徑估計延遲鎖定環(Multipath Estimation Delay Lock Loop, MEDLL)算法等,但均是單純針對尋找標準過0點的研究,而忽略了討論如何判斷天波干擾和地波識別問題。另外這些均是時域方法,容易受到各種強干擾而導致信號畸變,因此在惡劣條件下的適用性不高。針對以上難點,本文通過比較結合現有方法并提出聯合算法,在低信噪比、強天波干擾環境下解決了天地波的識別問題,并在此基礎上有效構造天波信號,利用去除天波后的偽地波成功實現周期識別,為后續解調解碼等過程提供了保障。

2 信號體制

Loran-C信號以脈沖組的形式發射,而eLORAN信號與傳統的Loran-C信號的區別主要是在第3~8個脈沖增加了數據調制,標準eLORAN脈沖時域波形用公式表示為[10]

實際中由于多跳天波離地波較遠,因此只考慮距離較近的一跳天波的影響。受到一跳天波干擾的信號時域波形可用公式表示為

其中,sground(t),ssky(t),τg,τs,A,B分別為地波、天波信號及其時延差和幅值;n(t)為噪聲干擾。

設定天地波時延差Δτ為100 μs,天地波幅度比(Sky-Ground Ratio, SGR)為6 dB并加以10 dB的白噪聲,根據式(1)、式(2)可以得到eLORAN標準脈沖及混合信號波形如圖1所示。

根據《船用羅蘭C接收設備技術標準》[10],羅蘭C信號的信噪比(Signal Noise Ratio, SNR)指信號接收電平與大氣噪聲均方根電平值之比,可用無量綱的數或分貝數表示。其中信號電平定義為脈沖包絡起點之后25 μs處連續波有效值電平,本文之后的仿真SNR均是根據以上定義為標準。另外,信號進入接收機后,首先會經過梳狀濾波器(40次左右線性累加),在保證信號不畸變的同時提高信噪比(15 dB左右),故文中的仿真SNR如無特別說明均是在天線接收端信噪比基礎上增加15 dB后的結果。

3 算法描述

eLORAN接收機進行信號處理,實現高精度PNT功能是依賴地波信號的,而天波信號的混入會影響周期識別效果,造成以下兩類致命的錯誤:(1)將天波誤捕獲為地波信號造成的地波識別錯誤;(2)天波影響地波的位置估計而導致的跳周錯誤。為了解決以上問題,提高尋找地波第3周過0點的準確度,本文對天波特性進行了深入分析,并提出了一種針對天波干擾的周期識別聯合算法。該算法可以分為以下3個部分,如圖2所示。

(1) 天地波識別模塊。當天地波位置較近或天波信號較強時,加之各種噪聲干擾等因素的影響,會造成地波識別錯誤,導致數十微秒甚至上百微秒級別的時延估計(Time Delay Estimation, TDE)誤差并最終引起授時定位偏差。而該模塊的目的正是在于將識別信號首先定位在地波信號上。

(2) 天地波分離模塊。經研究發現,高強度的天波干擾會對地波的位置估計產生較大影響。因此,經天地波識別、決策后,需要進行天地波分離的信號進入此模塊,通過自適應匹配算法得到天波的準確位置及強度并消除天波,得到較為干凈的偽地波信號。

(3) 周期識別模塊。該模塊將對以上兩個模塊的輸出結果進行進一步處理,利用基本無天波干擾的偽地波信號尋找到標準過零點位置,實現周期識別。

3.1 天地波識別模塊

天地波識別模塊的核心功能是判斷天波的有無,而該過程的實現主要是利用信號時延差、幅度比估計,以及二者的相互關系進行分析判斷。目前已經存在一些較為成熟的天地波時延差估計方法,包括多重信號分類(MUltiple SIgnal Classification,MUSIC), (Auto Regressive Moving Average,ARMA)及ESPRIT(Estimation of Signal Parameters via Rotational Invariance Techniques)算法等[11—14],都可以在不是很極端的信噪比環境下達到TDE需求。其中MUSIC算法分辨率更高,峰值更尖銳,但信噪比要求也最高;ARMA模型及ESPRIT算法抗噪能力稍強,但分辨率及估計準確度較差。但以上提到的算法均存在一個致命問題,即均是在“天波存在”假設下提出的技術,研究重點主要為天、地波的時延差估計,不能體現信號幅度特征,無法在盲測時利用幅度信息輔助判斷是否存在天波。近年來的一些新的研究成果,如EXIP-WRELAX(EXtended Invariance Principle Weighted fourier transform and RELAXation)算法、神經網絡算法等[15—17]也都存在同樣的問題,因此在實際應用時容易受到天波或交叉干擾等影響而導致錯捕或誤識別且無法糾正,這也是目前各種接收終端中周期識別一直無法較好實現的關鍵所在。

經典頻譜相除算法的最大優勢在于可以較好地反映信號幅度特征,這是判斷信號存在與否,實現天地波信號的識別功能的關鍵。因此本文對此算法進行了深入研究以改進其在強干擾環境下的分辨率及抗干擾性能要求。其基本原理推導為[18,19]

eLORAN信號的能量主要集中在85~115 kHz,因此在頻譜相除算法中通常會使用窗函數以降低帶外噪聲及干擾。文獻[18,19]均使用了帶寬50 kHz的漢寧窗進行頻域截斷,但是此種窗在天地波時延差較小時由于分辨率不足而無法達到效果,因此本文對窗函數進行了改進。通過對各類窗函數的研究比較可以得到一致的結論:窗口越窄,去噪效果越好,IFFT變換后峰越平坦(分辨率低),估計越準確;反之,窗口越寬,峰越尖銳,分辨率越高,而去噪效果也越差。另外分析表明,窗函數的下降梯度越大,兩側越平坦,則干擾抑制能力越強。因此為了兼顧高分辨率和干擾抑制這兩方面需求,可以在更寬窗口的基礎上對其梯度進行優化。經過大量對比驗證,本文選擇了在式(4)中取140 kHz寬的Parzen窗及h=2.3,以提高窗口下降梯度并使帶外部分更加平緩。根據以上改進,通過式(3)、式(4)步驟提取到了更惡劣條件下(SNR=5 dB(天線接收端—10 dB),SGR=18 dB, Δτ=38 μs)的天地波信息。

通過圖3(a)可以看出,改進的Parzen窗在保證窗寬度的同時對梯度有了明顯改善,圖3(b)的峰信息對比也顯示出改進窗的高分辨率可以成功將天波和地波區分開來,并維持了較好的去噪效果。但超過強天波閾值( Δτ≤60 μs, SGR>5 dB)的高強天波會嚴重影響TDE的準確度,加之窗函數引起δ峰左右兩側的旁瓣效應(Δτ<35.5 μs,分別小于—22 dB,—33 dB)及噪聲干擾等會使TDE誤差遠大于5 μs。這對周期識別來說極其致命,同時也是以往許多研究中都無法克服的問題。另外,強噪聲及其他臺站的交叉干擾信號等造成的尖峰也可能會使真實信號峰被淹沒,導致信號識別錯誤,故僅利用1次頻譜相除的峰信息分辨天波和地波是缺乏可信度的。為此,本文基于大數法則的思想設計了天地波識別算法,其步驟如下:

(1) 利用捕獲框將信號控制在采樣數據中心附近,連續采集N組峰信息(建議不少于50組以保證樣本數充足,每組1 ms),累加并取極大值點Pm,則由于信號已捕獲,Pm必為天波或地波信號標識;在其左右各取1個次大峰Pl, Pr,計算并預存峰值比及時延差信息;

(2) 在每組峰信息中尋找最接近Pm位置的峰Qm(i),并在其左右各統計1個次大峰Ql(i)和Qr(i)共N組。若某側存在信號,則應有大量次大峰集中在Pl或Pr附近,因此可分別統計此兩點±10 μs范圍內Ql(i)和Qr(i)的個數,利用經驗個數閾值(可取N的50%~75%)分別判斷左右兩側是否存在eLORAN信號,并根據旁瓣的幅度、位置特性排除將旁瓣誤判為信號的可能性;

(3) 根據天波判斷結果進行信號輸出抉擇:

(a) 若存在天波則根據信號時延關系確定天波、地波,并保存預估幅度比D、時延差 Δτ及等信息。由于弱天波干擾對地波的位置判斷的影響不大,因此可根據前段提到的強天波閾值進行判決,將參數D和 Δτ在閾值范圍內的信號輸入天地波分離模塊。

(b) 若判斷無天波信號干擾或參數值在強天波閾值范圍之外,則直接將當前信號及預估的地波位置信息τg0輸入周期識別模塊。

該模塊是后續模塊的基礎,其準確率主要受交叉干擾、突發干擾及噪聲的強度及頻度的影響,其中窗函數及帶通濾波器(Band-Pass Filter, BPF)可以高效抑制大量的帶外干擾及噪聲,而梳狀濾波器可壓制絕大多數的交叉干擾信號。另外,在接收設備非常接近某干擾臺站時,交叉干擾信號強度也有可能高于遠處臺站的跟蹤信號。此種情況下可同時捕獲該干擾信號并利用其脈沖組重復周期關系對消大部分交叉干擾能量。

接收信號經過各種前端處理后,仍可能有少量無法預測的同頻干擾或突發干擾等影響信號峰,因此在工程應用時應根據不同地點測試數據綜合考慮并給出經驗個數閾值以供接收設備設定。

該算法的計算復雜度主要取決于每組數據采樣點數M,以及頻譜相除和峰值檢測等步驟。其中:傅里葉變換(Fast Fourier Transform, FFT)、IFFT的時間復雜度均為O(M ·log2M);除法及窗函數運算復雜度均為O(M);檢測最大峰、次大峰的復雜度為O(M)。故利用N組峰信息進行天地波識別的時間復雜度為:O(N ·(2·M ·log2M+M)+N ·M+2·M+N ·M)=O(N ·M ·log2M)。相比于天地波分離模塊所需的10 MHz采樣率,此模塊對采樣率要求不高,在實際應用中可使用現有接收機通常使用的2 MHz采樣率以節省80%的計算量。以上算法利用了大數法則的理論,利用多組峰信息并選擇合適的個數閾值進行聯合判決,極好地解決了交叉干擾、同頻干擾以及其他突發干擾的影響,在利用模擬源及新型接收終端進行測試時幾乎沒有出現過錯誤,較好地實現了不同天波干擾環境下的天地波識別,準確地將捕獲框鎖定在地波信號上。

3.2 天地波分離模塊

天地波分離的核心是準確地重構天波,而強天波干擾會極大影響地波的位置估計精度,使TDE誤差大于5 μs并導致周期識別跳周。針對此問題,本文設計了自適應去天波算法,包含了兩階搜索過程以及用于匹配的混合信號幅度2維網格,以保證準確估計天波時延τs及幅度B的同時盡可能節省計算量。其中網格中映射的數值為無干擾、無噪聲環境下不同幅度(范圍0.1~18倍,步長0.1)、時延差(范圍37.5~200 μs,步長0.1 μs)天波干擾下的混合信號幅值,由仿真軟件生成后存儲于DSP存儲器中。該算法其流程圖如圖4所示,并包括如圖4的步驟。

(1) 將預估的時延差Δτ及幅度比D與預存的混合信號幅度網格C進行匹配,根據其數值與經過帶通濾波器的實際信號幅度間的比例關系以及 Δτ的搜索范圍,通過式(5)預估出第1階搜索中時延循環對應的天波幅度搜索中心B0(k)。由于不同時延差對應的天波幅度估計初值差別可能較大,因此,在第1階段時延差循環中得到不同的B0(k)可以達到縮小幅度循環范圍,提高效率的目的

(2) 以τs0,B0(k)為中心進行第1階段搜索,選擇相對較寬的搜索范圍以及較大的步長循環構造并去除偽天波,得到偽地波組;利用偽地波組的峰值信息與基準信號組合并與偽地波組對消得到殘差組并計算均方根(Root Mean Square, RMS),利用最小均方原則得到粗略的天波時延差及幅度參數值,B′。

其中,信號輸入接收機后會進行電平調整,故本文后續仿真統一將采樣信號幅度標準化為10,并選擇10 MHz采樣率以滿足天波幅度及時延估計精度需求。另外,進入此模塊的信號中天波功率相對較大,因此利用頻譜相除技術時對天波參數的估計誤差相對較小,時延誤差范圍一般不超過5 μs,而幅度誤差一般不超過25%。據此,兩階段搜索范圍及步長推薦設置如下:

根據上文中天地波分離模塊的詳細步驟說明以及圖4的具體實現流程,式(6)—式(9)結合各過程變量給出了幾處關鍵節點的計算方法,其中式(6)、式(7)表示兩階搜索中偽天波對消、偽地波對消及殘差RMS的循環計算方法;式(8)、式(9)為構造最優偽天波、產生最優偽地波的計算方法

根據圖5(a)所示,去天波前的天、地波時延估計分別為138.9 μs, 95.0 μs,與真實值偏差分別為1.1 μs, 5.0 μs,天地波幅度比D估計為2.1562。而運用本節提出的天波分離算法后,根據仿真過程量可知天波時延估計誤差為0,幅度比估計為2.0662,偏差為3.31%,即基本上去除了天波,如圖5(b)、圖5(c)所示。

3.3 周期識別模塊

以上兩個模塊成功地實現了天地波的識別和分離,極大降低了天波干擾的影響并準確定位地波信號。而周期識別模塊的目的是準確尋找偽地波信號的標準過零點,其中的關鍵正是識別正確的載波周期。目前對此問題的方法研究主要分為以下兩類:

(1) 時域方法:周期累加值法、波形匹配法、優化包絡相關算法、MEDLL算法等。

時域算法對信噪比要求較高,其中周期累加值法、周期和值法、波形匹配法的信噪比下限為23~25 dB;峰值比算法有了一定的性能提升,在信噪比達到15 dB 以上時能夠有效檢測標準過0點[9],因此有一定的工程應用價值;包絡相關算法分辨率較低,容易受到波形畸變影響且計算量相對較大;MEDLL算法在包絡相關算法的基礎上結合了包絡求導技術[6],一定程度提高了分辨率及信噪比性能,但當信噪比低于14 dB或天波過強時甚至無法分辨。

(2) 變換域方法:MUSIC算法、ARMA模型、ESPRIT算法、IFFT頻譜相除算法等。

相對于時域算法,各種變換域方法復雜度更高,但在信噪比條件及穩定性方面有著明顯優勢,3.1節總結了幾種方法的性能特點。其中,ARMA模型及ESPRIT算法抗噪能力稍強,但分辨率及估計準確度較差,容易導致跳周;MUSIC算法及頻譜相除算法相對穩定,在天波干擾較小時估計精確度極高,但MUSIC算法的信噪比要求相對較高[14]且計算量更大。

經過對上述算法理論的深入分析可以看出,各種時域算法較為簡單實用,但對信號接收環境有較高的要求。如峰值比算法,在噪聲及干擾環境適宜時有較高的穩定性,但也容易受到天波殘差以及各種噪聲干擾引起的波形畸變的影響,導致跳周錯誤,在工程應用中可結合環境估計結果決策使用。而本文研究重點為惡劣環境下的周期識別,噪聲環境較差,因此更適合變換域方法。下文選擇了其中穩定性高,計算量相對較小且信噪比要求更低的頻譜相除算法進行周期識別過程(經仿真章節驗證,SNR下限接近—2 dB,遠低于上述方法)。

利用頻譜相除方法推算標準過零點位置的基本流程與其他各種算法大致相同,即首先識別正確的載波周期,再檢測標準過零點位置,其具體步驟描述如下:

(1) 不需要天地波分離的信號在天地波識別后直接輸入此模塊,并將預存的地波位置τg0記為PA;對天地波分離模塊輸出的偽地波則利用頻譜相除技術,預估信號的起始位置點PA。由于偽地波中的天波殘差較小,對IFFT計算的分辨率要求不高,故可以在計算中選擇較窄的窗函數以實現更大程度的噪聲抑制,突出信號尖峰。

(2) 根據估計所得地波起始位置點PA推后30 μs預估第3載波周期過0點PB。雖然基本無天波干擾,但PA點的估計誤差會導致PB點與準確的標準過0點PC之間存在一個小的偏差。

(3) 根據PB點位置尋找最接近PB點的載波正向近0點PC并作為推算的標準過0點。

根據以上方式繼續對3.2節示例中輸出的偽地波進行處理,可得信號位置信息如圖6所示。

從幅度歸一化位置估計對比圖6(a)中可以看出,去除天波干擾后,偽地波信號的起始位置估計誤差已控制在半個載波周期范圍內(本例為0.3 μs),隨后如圖6(b)所示,點PC準確地標識了標準過0點的位置。該方法利用了頻譜相除技術的高穩定性,后續仿真結果可以驗證,在去除大部分天波干擾或無天波干擾時,利用偽地波可以得到精確的標準過0點位置以及地波的起始位置。

4 仿真實驗及結果分析

第3節分解討論了該自適應周期識別算法的各個模塊。本節對整體算法進行了充分的仿真實驗及結果展示,并對主要過程量進行了對比分析。其中,由于天地波識別算法在仿真實驗及不同地點的拷機測試中均無出錯,故沒有給出實驗數據。

根據eLORAN接收機的最低性能標準要求[4],當天線接收端信號SNR≥—10 dB, SGR≤12~26 dB,Δτ在37.5~60 μs范圍時接收設備應能進入鎖定狀態,而新型eLORAN接收機更要求在此范圍內能夠成功解碼。本文算法與3.1節提到的各種類似的變換域方法均有其信噪比極限,在噪聲功率過大時噪聲峰會淹沒信號峰,算法有效性也會急劇下降。其中在利用不同的窗及BPF處理且保證信號不失真時,頻譜相除算法信噪比下限不低于—2 dB(天線接收端—17 dB左右)。而在實際應用中,一般認為天線接收端信噪比小于0 dB時為低信噪比環境。基于以上分析,本文針對各種惡劣環境,對仿真中的可變參數作如下設計:

(1) 輸入信號:標準eLORAN信號+天波信號+隨機白噪聲,采樣率為10 MHz;

(2) 天地波時延差Δτ:38~60 μs,步長為1 μs;天地波幅度比D:SGR=0~24 dB,步長為6 dB;

(3) 信噪比SNR:0~10 dB(天線接收端—15~—5 dB),步長為1 dB。

本文在MATLAB2017環境下,獨立調整各參數,并在不同情況下分別運行了1000次。由于篇幅所限,文中只展示了SGR=0 dB, 18 dB和24 dB幾種情況下的部分仿真結果及數據統計,如表1所示。

表1體現了在SGR=0 dB時該算法極高的性能。由于算法運行過程中幅度比估計值超過強天波門限范圍,因此跳過了天地波分離而直接進入周期識別模塊。結果表明,此方式在節省了效率的同時非常有效地實現了周期識別,驗證了算法的合理性。

表1 周期識別錯誤次數及準確率(SGR=0 dB)

圖7、圖8分別展示了SGR=18 dB和24 dB時的天波位置及幅度的估計誤差。可以看出,位置和幅度的估計精度隨著SNR的增加而提升。經統計,位置誤差RMS最大值分別為0.0126 μs和0.0105 μs,幅值誤差率RMS最大值分別為1.09%和0.53%。通過對兩個參數非常準確的估計,可以有效地重構并去除天波。

通過3.1節討論可知,在強天波干擾環境下,改進前的頻譜相除技術甚至無法區分天地波,即使改進了窗函數提高分辨率也仍存在極大的位置誤差,導致周期識別跳周頻發;而如表2、表3所示,經過了天地波分離過程,極大地提高了在不同條件下周期識別的正確率,基本高于99%。

表2 天地波分離后周期識別錯誤次數及準確率(SGR=18 dB)

表3 天地波分離后周期識別錯誤次數及準確率(SGR=24 dB)

5 結論

本文所提聯合算法主要是針對各種惡劣環境所設計,這也是目前各種周期識別算法研究的難點所在。根據上文的算法討論和部分仿真結果展示以及更多的仿真數據統計,可以得出以下結論:

(1) 對窗函數的改進效果明顯,有效地將高強度天波干擾下的eLORAN天、地波信號區分開來,并可適配提出的天地波識別算法,成功實現對天、地波信號的準確定位。另外通過決策過程準確判斷是否需要進行天波分離,大幅節省了后續過程的計算量;

(2) 對天波的參數估計準確度指標均隨SNR及天地波時延差的增加而提高,在極端惡劣的噪聲及天波干擾環境下仍然能夠達到極高的水平,且本文所提2階搜索方式在大幅提高效率的同時可以有效地重構并去除天波信號,輸出較為穩定的偽地波信號;

(3) 本聯合算法可以保證在接收天線輸入端SNR=—15 dB的極低信噪比下也能夠保持極高的水準,周期識別準確率基本超過99%,完全可以滿足新型eLORAN接收機的需求。

eLORAN接收終端的PNT精度很大程度上取決于周期識別的準確度,但高強度的天波干擾一直是信號處理中的一個關鍵問題,嚴重限制了周期識別算法的完整性及有效性。針對這一問題,本文提出一種自適應天波識別、分離和周期識別的聯合算法。經過理論分析和各種不同接收環境下的仿真實驗,驗證了算法的可行性。本算法首先克服了經典IFFT頻譜相除算法無法同時兼顧高分辨率和時延估計準確度需求的弊端,并結合其能夠體現信號幅度特征的特點,有效解決了以往研究中常被忽略的天地波信號識別問題。其次,能夠在惡劣的條件下成功將天地波的識別、估計、分離及尋找標準過零點等過程有機結合,在提高效率的同時保證了周期識別的正確率,具有良好的魯棒性和可實現性。同時,精確的天波估計,也為接收機處理交叉干擾,提高解調解碼的成功率提出了新的思路和理論參考。

猜你喜歡
信號
信號
鴨綠江(2021年35期)2021-04-19 12:24:18
完形填空二則
7個信號,警惕寶寶要感冒
媽媽寶寶(2019年10期)2019-10-26 02:45:34
孩子停止長個的信號
《鐵道通信信號》訂閱單
基于FPGA的多功能信號發生器的設計
電子制作(2018年11期)2018-08-04 03:25:42
基于Arduino的聯鎖信號控制接口研究
《鐵道通信信號》訂閱單
基于LabVIEW的力加載信號采集與PID控制
Kisspeptin/GPR54信號通路促使性早熟形成的作用觀察
主站蜘蛛池模板: 亚洲成人精品| 波多野结衣在线se| AV在线天堂进入| 欧美日韩专区| 久久特级毛片| 精品欧美一区二区三区久久久| 国产精品欧美亚洲韩国日本不卡| 丁香六月激情综合| 日韩AV无码免费一二三区 | 美女内射视频WWW网站午夜| 99视频在线精品免费观看6| 日本91视频| 久夜色精品国产噜噜| 国产一二视频| 久久网欧美| 青青青视频蜜桃一区二区| 毛片最新网址| 四虎免费视频网站| 精品国产Av电影无码久久久| 亚洲中久无码永久在线观看软件| 国产成人精品高清不卡在线| 国产亚洲精品无码专| 在线视频亚洲欧美| 最新亚洲人成网站在线观看| 福利视频一区| 国产乱人伦精品一区二区| 国产欧美日韩精品综合在线| a毛片在线| 日韩欧美国产成人| 91视频99| 人妻熟妇日韩AV在线播放| 国产成人精品免费视频大全五级 | 免费在线一区| 男女猛烈无遮挡午夜视频| 亚洲国产精品成人久久综合影院| 日韩精品一区二区深田咏美| 精品人妻无码中字系列| 亚洲色无码专线精品观看| 日韩欧美中文字幕在线韩免费| 免费看美女自慰的网站| 国外欧美一区另类中文字幕| 91在线播放国产| V一区无码内射国产| 无码内射在线| 国产综合日韩另类一区二区| 国产精品真实对白精彩久久| 亚洲最大看欧美片网站地址| 99久久国产综合精品2023| 色哟哟国产精品一区二区| 欧美一级黄色影院| 国产精品亚欧美一区二区| 精品国产黑色丝袜高跟鞋| 婷婷中文在线| 天堂在线www网亚洲| 欧美性色综合网| 无码日韩视频| 欧美五月婷婷| 91精品啪在线观看国产91九色| 免费人成网站在线观看欧美| 亚洲va视频| 久久午夜夜伦鲁鲁片无码免费| 亚洲国产成人自拍| 欧美日韩国产精品va| 久久免费观看视频| 久久亚洲国产视频| 亚欧成人无码AV在线播放| 大学生久久香蕉国产线观看| 婷婷综合亚洲| 综合色天天| 99久久精品免费看国产电影| 99热亚洲精品6码| 波多野吉衣一区二区三区av| 精品久久久无码专区中文字幕| yy6080理论大片一级久久| 亚洲人成人伊人成综合网无码| 国产精品深爱在线| 亚洲中文字幕在线一区播放| 国产网站免费观看| 亚洲乱码精品久久久久..| 欧美在线视频不卡第一页| 国产亚洲欧美在线人成aaaa| 亚洲国产天堂在线观看|