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

基于分數階傅里葉變換的線性調頻引信定距方法

2015-02-28 10:46:08岳凱郝新紅栗蘋陶艷李永亮
兵工學報 2015年5期
關鍵詞:信號方法

岳凱,郝新紅,栗蘋,陶艷,李永亮

(北京理工大學 機電學院,北京100081)

0 引言

傳統雙通道諧波定距調頻(FM)引信的定距精度主要受調制頻偏影響。調制頻偏越大,定距精度越高。大的FM 頻偏必然會帶來大的寄生調幅和工程實現難題。因此,本文探索了一種新的FM 引信定距方法,以期在易于工程實現的調制頻偏內獲得較高的定距精度。參考現代數學成果在工程中的應用以及現代信號處理技術發展[1-2],分數階傅里葉變換(FRFT)作為一種基于Chirp 基的正交分解算法,在處理連續波線性FM 信號方面具有天生的優勢[3]。國內外很多學者展開了FRFT 在處理Chirp信號方面的研究[4-7],并取得了顯著成果[8-13]。但將FRFT 研究成果應用于FM 引信實現其定距功能還未見報道。為此,本文提出了基于分數階域瞬時初始頻率估計和基于分數階相關的兩種線性FM 引信定距方法。

1 FRFT 及其離散算法

假設信號f(t)的角度α 的FRFT 表示為Fα(u),則FRFT 的定義式[13]為

(1)式的變換核函數為

式中:α 為信號的Wigner 分布在時間-頻率平面上投影的順時針旋轉角度,如圖1所示。如果角度α或(α+π)為2π 的整數倍,那么核函數Kα(t,u)將會退化為δ(t-u)或δ(t+u).

圖1中,坐標軸u 為相應的分數階域坐標軸;fif(t)為線性FM 信號的瞬時頻率曲線;IIF 為fif(t)在零時刻的值,即線性FM 信號初始頻率;up為線性FM 信號在分數階域的最優FRFT 階次下能量聚集的坐標位置。根據FRFT 的性質[3],假設需處理的線性FM 信號與某一Chirp 基吻合,即在圖1中變換角度其FRFT 的結果必定為u 軸上的一個有限沖激,且其分數階域坐標為up.

圖1 定義在時頻平面上的FRFTFig.1 The definition of FRFT on phase plane

為了便于工程上實現,FRFT 必須具有快速的離散算法。各國學者對離散FRFT 算法的研究已經取得了較多成果,提出了很多可行的算法,尤其是Ozaktas 的采樣型算法,以速度快、精度高、具有與快速傅里葉變換(FFT)比擬的計算復雜度而受到廣泛采用[14],其算法可表示為

2 基于分數階域瞬時初始頻率估計的FM引信定距方法

2.1 定距實現原理

基于分數階域瞬時初始頻率估計的FM 引信定距方法實現原理框圖如圖2所示。該方法的定距原理是利用線性FM 信號在特定角度(FM 率角度)的FRFT 為一沖激函數的性質,求取最優FRFT 階次下發射與接收信號的FRFT 譜峰值在u 域的位置差,從而獲得發射與接收信號的初始頻率(如圖1所示的IIF1 與IIF2),并通過初始頻率差與時延估計的對應關系獲得彈目間距離信息。由于該方法通過估計發射與接收信號的初始頻率差得到時延估計值,且初始頻率時變的,因此稱該方法為瞬時初始頻率估計。

圖2 基于分數階域瞬時初始頻率估計的定距方法原理框圖Fig.2 Block diagram of the ranging method based on the initial frequency estimation in FRFT domain

首先,對FM 引信的發射信號和接收信號在射頻端進行下變頻,濾除載波頻率,獲得發射支路和接收支路的中頻信號;而后,對兩路中頻信號進行固定α 角度的FRFT,因為對引信而言,中頻信號的FM 率是已知的,所以變換角度其 中,K 為FM 率;再次,在固定α 角度的分數階域搜索最大值點位置up,根據up求得瞬時初始頻率;最后,根據瞬時初始頻率與延遲時間的對應關系求得延遲時間,進而求得彈目距離。當彈目距離的估計值與預設的引爆距離相等時,給出引信點火信號。

假設FM 引信的調制信號為鋸齒波,如圖3所示,那么FM 引信發射信號可表示為

式中:a 為發射信號幅值;f0為信號中心頻率;T 為調制信號周期;K =B/T 為調制信號的FM 率,B 為調制頻偏。

圖3 鋸齒波FM 信號頻率隨時間變化圖Fig.3 Instantaneous frequency of saw tooth wave

回波信號經目標反射后,相對于發射信號只是延遲了時間τ(t). 在實際引信應用中,每個運算周期內τ(t)的變換十分微小,可近似為常數,所以可將τ(t)在每個運算周期內看成是與時間t 無關的獨立變量τ,因此單個周期內(后面出現的時域信號表達式均指單個周期之內的信號)的回波信號可表示為

式中:Φ0n=πn(n+1)KT2+Φ0.

對兩路中頻信號的相位表達式求導,可得信號在零時刻的頻率,即初始頻率為

由于(7)式的初始頻率是時變的,所以稱之為瞬時初始頻率。

結合(7)式、(8)式和cot α = -K,可獲得延遲時間的估計

因此,根據關系式τ=2R/c,可獲得彈目距離的估計值為

由(10)式可知,該方法可用于線性FM 引信的彈目距離估計。如果回波的延遲時間為固定值,則可進行預定距離的估計;如果回波延遲時間是動態變化的,則可實現連續測距功能。

2.2 定距方法性能分析

2.2.1 計算量分析

根據圖2可知,要獲得延遲時間估計,需同時估計兩個中頻信號的初始頻率。該過程中只涉及到兩路信號的FRFT,同時,中頻信號的FM 率相對發射信號是不變的(即FM 率已知),所以只需進行固定角度的FRFT,無需在整個分數階域掃描,大大減少了算法的計算量。由于FRFT 采用Ozaktas 的采樣型離散算法[14],所以離散FRFT所需的乘法計算量為

式中:N 為離散FRFT 的數據點數。其他部分所需的乘法運算次數相對(11)式可忽略。所以基于分數階域瞬時初始頻率估計的FM 引信定距方法所需的乘法計算量即如(11)式所示。

2.2.2 距離分辨力分析

由(10)式可知,該方法的固有距離分辨力決定于信號的時延分辨力,即決定于在α 階分數階傅里葉域up估計值的主瓣寬度。理想情況下,線性FM信號最優階次的FRFT 的幅度譜是一個δ 函數,具有無限窄的主瓣寬度。文獻[15]經由推導線性FM信號的FRFT 表達式已證明此點。因此,固有距離分辨力是無窮小,即對具有任意延遲時間差異的兩個信號,該方法可以分辨。

實際工程應用中,采用離散方法求取彈目間距離估值時,該方法的距離分辨力還受到信號重復周期數及采樣點數的影響。因為在有限觀測時間Td內,線性FM 信號分數階域傅里葉變換的積分區間為觀測區間,可證明其在最優階次FRFT 的幅度譜由沖激函數演化為sinc 函數,且分數階域up估計值的-3 dB 主瓣寬度[16]。

考慮距離估計值來自離散采樣集合{us,2us,…,Nus},其中,us為最優階次分數階域離散步長。設定觀測區間Td=mT,將(12)式帶入(10)式可得該方法的距離分辨力為

設定采樣頻率為fs(fs≥2B,滿足采樣定理),采樣點數為N,則距離分辨力可表示為

(13)式與(14)式為距離分辨力的兩種不同表述方式。由(13)式、(14)式可知,該方法的距離分辨力主要受信號帶寬B、運算周期數m、調制周期T、運算點數N 等參數影響。為了提高距離分辨力,在不增加調制帶寬的條件下,只需增加運算點數或減小調制周期即可。相對于諧波定距FM 引信增大調制頻偏來提高距離分辨力的方法,這種方法更容易實現,且可提高定距引信的精度。

此外,由(13)式可知,當運算周期m =1 時,該方法的時延分辨力約為1/B,即距離分辨力約為c2B,等效于FM 諧波定距方法的定距精度。

3 基于分數階相關的FM 引信定距方法

3.1 定距實現原理

基于分數階相關的FM 引信定距方法實現原理如圖4所示。分別對參考信號(發射信號)和回波信號進行下變頻,獲得中頻信號;參考信號的中頻信號延遲特定時間(該時間與引信預定作用距離相對應);對延遲后的參考信號的中頻信號與回波信號的中頻信號進行離散FRFT;變換后的兩路信號相乘并累加,獲得兩路信號的相關信號;對相關信號進行閾值檢測,滿足條件時即可獲得預定距離處的點火信號。

圖4 基于分數階相關的FM 引信定距方法原理框圖Fig.4 Block diagram of the ranging method based on the fractional correlation

假設FM 引信以鋸齒波為調制信號,那么以(4)式發射信號作為參考信號,(5)式為回波信號,參考信號和回波信號通過與單載波信號x(t)=cos(2πf0t)混頻進行下變頻,并對下變頻后的發射信號做固定延遲τref,獲得兩路中頻信號為

分數階相關等價于對信號做FRFT,再在分數階域進行一般意義的相關運算,因此其數學表達式[3]可表示為

對(16)式取模,可得相關峰函數為

設x(t)為要檢測的目標信號,輸入信號s(t)=x(t-τ0),參考信號h(t)=x(t -τref),對(17)式求導并令其等于0,可得到分數階相關峰的位置為

本方法所用的分數階相關與如上所述的分數階相關并不完全相同,而是采用先延遲,而后進行FRFT,再求取內積的計算方法。如圖4所示,本方法中相關輸出可定義為

式中:up為分數階s(t-τref)的FRFT 所得的有限沖激函數在分數階域的位置。

最后,對(20)式所示的相關函數做閾值檢測,當滿足閾值條件時,給出點火信號。

該方法是一種基于發射信號支路固定延遲的定距方案,故可用于預定距離估計,不能用于連續估計瞬時距離。

3.2 定距方法性能分析

3.2.1 計算量分析

根據(15)式的中頻信號和(4)式的發射信號表達式,可知FM 率是保持不變的,即中頻信號的FM 率是已知的。因此,只需對中頻信號進行固定角度的FRFT,相對于階數搜素的算法可大大減少運算復雜度。

為了獲得一個相關值,首先需對兩個中頻信號分別進行FRFT,由于方法中的FRFT 采用的是Ozaktas的采樣型離散算法[14],故進行FRFT 的乘法計算量為完成FRFT 后,需對N 點的變換結果進行乘累加(即相關運算),其乘法計算量為N. 故計算一個相關值所需的乘法計算量為

3.2.2 檢測概率分析

假設系統輸入噪聲為高斯噪聲,由于FRFT 是一種線性變換,所以輸出的噪聲仍然服從高斯分布。設其均值和方差分別為un和σ2n,設檢測門限為λ.根據檢測理論,在貝葉斯準則下的檢測概率PD和虛警概率PF分別為

3.2.3 距離分辨力分析

由(19)式和(20)式可知,該方法的定距原理基于信號在最優階次下分數階傅里葉域up的估值,獲得彈目距離估計。本文定義的分數階相關函數在分數階域也是一個δ 函數,具有無限窄的主瓣寬度,所以延遲時間估計的理論分辨力無窮小,即對具有任意延遲時間差異的兩個信號,該算法是可分的。

因此,本質上該方法的距離分辨力等同于第一種方法。當采用離散方法求取彈目間距離估值時,該方法的距離分辨力同樣會受到信號重復周期數及采樣點數的影響。這里不再贅述。

與瞬態初始頻率估計方法不同的是,其利用分數階域幅度譜的模平方求取譜峰值(如(20)式所示),而前者是利用幅度譜的模求取譜峰值。

4 仿真與討論

本文分別對兩種方法進行了Matlab 仿真來驗證其可行性與測距性能。

仿真參數如下:FM 帶寬B =30 MHz,調制周期T=166 μs,彈目交會距離取[25 m,50 m]. 當運算周期為一個調制周期時,則距離分辨力約5 m.

分別設定引信的預定起爆距離為35 m、40 m、45 m,基于分數階相關方法和基于瞬時初始頻率估計方法的預定起爆距離估計仿真結果分別如圖5和圖6所示。

圖5 不同預定距離條件下中頻信號的分數階相關函數Fig.5 Fractional correlation of the intermediate frequency signal under the conditions of different fixed distances

圖6 不同距離處中頻信號的FRFTFig.6 Fractional Fourier transform of the intermediate frequency signal at different distances

由圖5與圖6可知:引信3 個預定的起爆距離分別對應于最優階次下的分數階域相關函數峰值與中頻信號的FRFT 譜峰值,證明了兩種定距方法的可行性;當運算周期為一個調制周期時,兩種方法譜峰值下降-3 dB 時對應的主瓣寬度為5 m,即圖5相關函數峰值與圖6FRFT 幅度譜峰值下降1/2 時所對應的主瓣寬度,與理論距離分辨力相吻合。

為了分析兩種方法的抗噪聲性能,分別對兩種方法進行了不同信噪比(-20 ~20 dB)條件下的1 000點蒙特卡洛仿真實驗。設定引信預定起爆距離為40 m,加載的噪聲類型為高斯白噪聲,其他仿真參數同上,不同信噪比下兩種方法蒙特卡洛實驗方差與均值的仿真結果分別如圖7和圖8所示。

圖7 不同信噪比條件下兩種算法的1 000 點蒙特卡洛仿真實驗結果(方差)Fig.7 Results of 1 000 points Monte-Carlo simulation at different SNRs(variance)

圖8 不同信噪比條件下兩種算法的1 000 點蒙特卡洛仿真實驗結果(均值)Fig.8 Results of 1 000 points Monte-Carlo simulation at different SNRs

對比分析圖7、圖8:信噪比大于0 dB 時,兩種方法均具有良好的抗噪聲性能;當信噪比小于0 dB時,初始頻率估計方法抗噪性能急劇惡化。仿真結果說明,分數階相關方法抗噪聲性能優于初始頻率估計方法,更適合于低信噪比環境下(如強背景噪聲干擾)信號檢測。

分析調制頻偏、運算點數等參數對距離分辨力的影響,其仿真結果如圖9所示,其中采樣率為10 倍調制頻偏,其他參數:調制頻偏10 ~50 MHz,運算點數100 ~5 000,調制信號頻率600 kHz.

圖9 調制頻偏和運算點數對距離分辨力影響Fig.9 Influences of the modulation frequency derivation and operation points on range resolution

圖9仿真結果表明:分數階域瞬時初始頻率估計的定距方法在有限增加的調制頻偏內(受限于探測器調制帶寬與體積約束),增加運算點數可顯著提高引信的距離分辨力。

將cot α= -K,K=B/T 代入(14)式,可得

(23)式表明:增加調制頻偏與運算點數可顯著提高引信的距離分辨力;若取采樣頻率為10 倍調制頻偏,則引信距離分辨力與運算點數呈反比。該結論與圖9的仿真結果相吻合。該結論同時表明,引信在不增大調制頻偏的前提下,通過稍微加大運算點數的方法可顯著提高其距離分辨力。

固定調制頻偏,取仿真參數為:調制頻偏30 MHz,運算點數250、700、1 200,采樣率10 倍調制頻偏,調制信號頻率600 kHz.

根據(14)式,在采樣率為10 倍調制頻偏的條件下可計算出距離分辨力,如表1所示。

表1 距離分辨力對比仿真結果Tab.1 Comparison of range resolutions of different ranging methods

表1結果表明:在相同的調制信號參數條件下,運算點數無需太大,兩種方法的定距性能優于傳統的雙諧波定距方法。

圖10 為仿真了運算點數為1 200 點情況下,不同距離處目標的FRFT. 由圖10 可知,30 m 處與30.5 m 處目標不可分辨,而與32 m、35 m 處目標峰值可區分,證明該運算點數下距離分辨力約2 m.

圖10 運算點數為1 200 點時不同距離處目標回波中頻信號的FRFTFig.10 FRFT of the intermediate frequency signals at different distances with calculation points of 1 200

圖11 為仿真了700 點運算點數的情況下,不同距離處目標的FRFT. 從中可見,30 m 與30.5 m、32 m 處的目標不能區分,而32 m 與35 m 處目標可分辨,證明距離分辨力約3 m. 仿真結果與理論計算距離分辨力相吻合。

圖11 運算點數為700 點時不同距離處目標回波中頻信號的FRFTFig.11 FRFT of the intermediate frequency signals at different distances with calculation points of 700

最后,仿真了0 dB 高斯白噪聲背景下,分數階域瞬時初始頻率估計方法的動態測距性能,仿真結果如圖12 所示,估計值與理想值吻合較好,證明了該定距方法的可行性。

5 結論

圖12 0 dB 信噪比時目標距離估計仿真結果Fig.12 Simulation results of target range estimation for SNR=0 dB

本文提出了基于分數階域瞬時初始頻率估計和基于分數階相關的FM 引信定距方法。理論推導了兩種方法的定距原理,仿真證明了兩種方法的可行性。在定距性能方面,兩種方法的距離分辨力相同,在相同的調制頻偏條件下只需不大的運算點數,其定距精度就可超過傳統的雙諧波定距方法;在算法計算量上,由于兩種方法均利用了回波信號FM 率,已知這一先驗條件,避免了對最優FRFT 階數的搜索,使計算量控制在與FFT 相當的水平,滿足工程應用的要求。此外,基于分數階域瞬時頻率估計的FM 引信定距方法可實現連續測距功能,可用于引信高度表或具有炸高精確分檔的引信。

References)

[1]Namias V. The fractional order Fourier transform and its application to quantum mechanics[J]. Ima Journal of Applied Mathematics,1980,25(3):241 -265.

[2]張賢達,保錚. 非平穩信號分析與處理[M]. 北京:國防工業出版社,1998.ZHANG Xian-da,BAO Zheng. Non-stationary signal analysis and processing[M]. Beijing:National Defense Industry Press,1998.(in Chinese)

[3]陶然,齊林,王越. 分數階Fourier 變換的原理與應用[M]. 北京:清華大學出版社,2009.TAO Ran,QI Lin,WANG Yue. Theory and application of fractional Fourier transform[M]. Beijing:Tsinghua University Press,2009. (in Chinese)

[4]溫景陽,張煥宇,王越. 線性調頻脈沖壓縮雷達信號參數估計方法[J]. 北京理工大學學報,2012,32(7):746 -750.WEN Jing-yang,ZHANG Huan-yu,WANG Yue. Parameters estimation algorithm of LFM pulse compression radar signal [J].Transactions of Beijing Institute of Technology,2012,32(7):746 -750. (in Chinese)

[5]馬艷,羅美玲. 基于分數階傅里葉變換水下目標距離及速度的聯合估計[J]. 兵工學報,2011,32(8):1030 -1035.MA Yan,LUO Mei-ling. FRFT-based joint range and radial velocity estimation of underwater target[J]. Acta Armamentarii,2011,32(8):1030 -1035. (in Chinese)

[6]鄧兵,王旭,陶然,等. 基于分數階傅里葉變換的線性調頻脈沖時延估計特性分析[J]. 兵工學報,2012,33(6):764 -768.DENG Bing,WANG Xu,TAO Ran,et al. Performance analysis of time delay estimation for linear frequency-modulated pulse based on fractional Fourier transform [J]. Acta Armamentarii,2012,33(6):764 -768. (in Chinese)

[7]鄧兵,陶然,董云龍. 基于分數階傅里葉變換的標量脫靶量測量新方法[J]. 兵工學報,2010,31(12):1627 -1631.DENG Bing,TAO Ran,DONG Yun-long. New method for scalar miss distance measurement based on the fractional Fourier transform[J]. Acta Armamentarii,2010,31(12):1627 -1631. (in Chinese)

[8]Ozaktas H M,Barshan B,Mendlovic D,et al. Convolution,filtering,and multiplexing in fractional Fourier domains and their relationship to chirp and wavelet transform[J]. Journal of the Optical Society of America A,1994,11(2):547 -559.

[9]Ozaktas H M,Aytur O. Fractional Fourier domains[J]. Signal Process,1995,46(1):119 -124.

[10]Ozaktas H M,Zalevsky Z,Kutay M A. The fractional Fourier transform with applications in optics and signal processing[M].New York:Wiley,2000:1 -513.

[11]Olcay A. Fractional convolution and correlation via operator methods and an application to detection of linear FM signals[J].IEEE Transactions on Signal Processing,2001,49(5):979 -993.

[12]Akay O. Unitary and Hermitian fractional operators and their extension:fractional Mellin transform,joint fractional representations and fractional correlations [D]. Kingston:University of Rhode Island,2000.

[13]Almeida L B. The fractional Fourier transform and time-frequency representations [J]. IEEE Transactions on Signal Process,1994,42(11):3084 -3091.

[14]Ozaktas H M,Arikan O,Kutay M A,et al. Digital computation of the fractional Fourier transform [J]. IEEE Transactions on Signal Process,1996,44(9):2141 -2150.

[15]張淑寧,趙惠昌,吳兵. 基于分數階傅里葉變換的偽碼體制引信線性調頻干擾抑制技術[J]. 兵工學報,2006,27(1):32 -36.ZHANG Shu-ning,ZHAO Hui-chang,WU Bing. LFM interference excision technique in pseudo random code fuse based on fractional Fourier transform[J]. Acta Armamentarii,2006,27(1):32 -36.(in Chinese)

[16]張南,陶然,單濤,等. 基于分數階傅里葉變換的線性調頻信號分辨率分析[J]. 電子學報,2007,35(12A):8 -13.ZHANG Nan,TAO ran,SHAN Tao,et al. Analysis of the resolution of the linear frequency modulated signal based on the fractional Fourier transform[J]. Acta Electronica Sinica,2007,35(12A):8 -13. (in Chinese)

猜你喜歡
信號方法
信號
鴨綠江(2021年35期)2021-04-19 12:24:18
完形填空二則
學習方法
孩子停止長個的信號
用對方法才能瘦
Coco薇(2016年2期)2016-03-22 02:42:52
基于LabVIEW的力加載信號采集與PID控制
一種基于極大似然估計的信號盲抽取算法
四大方法 教你不再“坐以待病”!
Coco薇(2015年1期)2015-08-13 02:47:34
賺錢方法
捕魚
主站蜘蛛池模板: 在线国产你懂的| 无码专区在线观看| 麻豆精品在线播放| 日韩天堂在线观看| 亚洲精品国产日韩无码AV永久免费网| 一级毛片不卡片免费观看| 亚洲一区二区三区国产精华液| 人妻无码中文字幕第一区| 精品国产www| 日韩无码视频网站| 久久大香香蕉国产免费网站| 91毛片网| 亚洲一区二区视频在线观看| 亚洲国产成人久久77| 操国产美女| 亚洲精品无码成人片在线观看| 国产视频一二三区| 国产精品手机在线播放| 国产成人亚洲日韩欧美电影| 国产免费久久精品44| 中日无码在线观看| 日韩精品无码免费一区二区三区 | 午夜一区二区三区| 免费无遮挡AV| 高清欧美性猛交XXXX黑人猛交| 日韩av资源在线| 国产免费福利网站| 91成人在线免费视频| 欧美精品伊人久久| 久久久久九九精品影院| 人禽伦免费交视频网页播放| 国产女人水多毛片18| 成人久久精品一区二区三区| 日本亚洲国产一区二区三区| 免费一极毛片| 四虎亚洲国产成人久久精品| 亚洲综合片| 一本大道香蕉高清久久| 四虎永久免费在线| 国产精品短篇二区| 蜜桃臀无码内射一区二区三区| 女人18毛片水真多国产| 一级毛片免费播放视频| 午夜成人在线视频| 91国内在线视频| 欧美视频在线播放观看免费福利资源| 永久在线播放| 亚洲第一成年网| 人人妻人人澡人人爽欧美一区 | 国产av一码二码三码无码 | 欧美成人影院亚洲综合图| 日韩欧美国产三级| 毛片一区二区在线看| 欧美一级视频免费| 国产特级毛片| 人妻21p大胆| 成人福利在线观看| 国产99免费视频| 一本大道香蕉久中文在线播放| 欧美日韩一区二区三| 伊人精品成人久久综合| 国产福利免费视频| 国产精品性| 青青久在线视频免费观看| 91色在线观看| 91国内外精品自在线播放| 欧美无专区| 国产aⅴ无码专区亚洲av综合网 | 亚洲人成网站18禁动漫无码| 2048国产精品原创综合在线| 色综合久久久久8天国| 亚洲丝袜中文字幕| 91日本在线观看亚洲精品| 亚洲综合专区| 国产成人精品免费视频大全五级 | 91免费片| 8090成人午夜精品| 久久无码高潮喷水| 992tv国产人成在线观看| 欧美国产精品不卡在线观看| 日韩亚洲高清一区二区| 人人91人人澡人人妻人人爽|