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

高溫結構激光掃描模態測試的降噪與參數識別

2020-06-16 03:15:28朱天煦臧朝平
強度與環境 2020年2期
關鍵詞:模態信號方法

朱天煦 臧朝平

高溫結構激光掃描模態測試的降噪與參數識別

朱天煦 臧朝平

(南京航空航天大學能源與動力學院,南京,210016)

針對高溫結構激光掃描模態測試中存在的噪聲問題,提出了一種基于模態峰值漢克爾奇異值分解的降噪處理與模態參數識別方法。首先,將測試得到的頻率響應函數(Frequency Response Function(FRF))經過FFT逆變換得到對應的時間域脈沖響應函數(Impulse Response Function(IRF)),并通過漢克爾奇異值分解(Hankel Singular Value Decomposition (HSVD)),進一步將脈沖響應函數轉化為按能量從大到小順序排列的一系列分量信號組合;其次,以恢復所有關心的模態峰值為基準,將分量信號從前到后累加,并在所關心的模態峰值完全恢復后,將剩余分量信號當作噪聲舍棄掉;該步驟會去除掉信號中包含的大部分噪聲,但仍會有一些殘余噪聲不可避免地被恢復;再次,對步驟二中提取得到的分量信號,采用基于模態峰值頻率通帶的迭代選取進行二次濾波,以分離出屬于模態峰值的分量信號,進而將它們累加為降噪后的IRF信號,并轉換至頻域以獲取降噪后FRF信號;最后,對降噪后的頻響函數進行模態辨識以提取模態參數。本方法應用于600度高溫環境下一個直板葉片的激光掃描模態測試數據的處理,結果表明頻響函數中的噪聲被有效濾去,模態參數可準確地提取,顯示了方法的有效性和優越性。

高溫;非接觸測量;漢克爾奇異值分解;模態峰值;降噪處理

0 引言

在航空、航天工程中,如飛行器等結構受到高溫的作用,其剛度特性會發生變化,從而會改變結構的動力學特性。因此,開展高溫結構的模態測試與參數識別,以獲取高溫條件下的結構動力學特性,具有重要的意義。在高溫環境下,接觸式傳感器的使用受到限制,采用非接觸式的激光傳感器進行高溫環境的模態測試,具有很高的優越性。但是,高溫環境會改變結構的表面狀態,會不可避免地使結構表面反射而來的激光帶有一定的噪聲,進而使測試得到的頻響函數包含噪聲,并干擾模態參數的準確提取。因此,需要研究一種相應的去噪方法,實現強噪聲影響的測試頻響函數的模態參數識別。

在過去的幾十年內,信號降噪手段得到了很大的發展,常用的信號處理方法有奇異值分解(Singular Value Decomposition(SVD))降噪[1],小波降噪[2],經驗模態分解降噪(Empirical Modal Decomposition(EMD))[3]等。在這些方法中,由于奇異值分解具有簡單,非參數的特性,其應用最為廣泛。近年來,基于脈沖響應信號(IRF)的漢克爾奇異值分解(HSVD)方法在動力學測試數據降噪中得到了廣泛應用[4]。其工作原理是將含有噪聲的IRF信號分解為一系列的奇異值和對應的分量信號,通過選取真實信號的奇異值和分量信號進行重構,舍棄由噪聲信號產生的奇異值和分量信號以實現降噪處理。然而,如何選取有用的分量信號或奇異值至關重要。Zhao[5]提出了一種基于奇異值差分譜(Difference Spectrum of Singular Values(DSSV))的降噪方法,該方法提取了DSSV最大值及之前所對應的奇異值以完成降噪處理。Bao[6]提出了一種基于模型階數指標(Model Order Indicator(MOI))的奇異值選取方法,該方法選取了MOI最大值所對應的奇異值及之前的數據進行重構以完成降噪。以上兩種方法主要利用了真實信號和噪聲所對應的奇異值之間會有一個巨大的差異這一特點來選取有效的奇異值。但是,當信號中噪聲能量較大而待提取的信號特征能量較小時,上述方法就不再適用。近年來,學者們進一步提出了一些基于分量信號特性的選取方法來提升HSVD的降噪效果。文獻[7]提出了一種基于周期調制強度(Periodic Modulation Intensity(PMI))的分量信號選取方法,有效地提升了旋轉機械故障信號的故障特征。文獻[8]提出了一種基于信號頻率分量的奇異值選取方法,通過分析分量信號、奇異值和頻率分量之間的關系,有效地從多頻率分量信號中提取出了單一頻率分量。上述方法均在旋轉機械故障信號處理中取得了較好的效果。基于模態峰值的漢克爾奇異值分解的降噪方法可用于單個FRF信號的降噪處理,其采用累加迭代和二次選取迭代來確定IRF信號中與模態信息相關的分量信號,進而將它們與噪聲產生的分量信號分離出來,以合成降噪后IRF信號,再將其轉至頻域的FRF信號。該方法無需提取與頻率分量相關的奇異值,并以其作為中間參數提取相應的分量信號來完成頻率分量的提取和降噪,而是采用待提取頻率分量的頻率和振幅,直接選取對應的分量信號來完成頻率成分的提取和降噪處理,因此具有更高的可操縱性。本文將該方法應用于600度高溫環境下直板葉片的激光掃描模態測試數據的降噪處理和模態參數提取,以驗證該方法對高溫結構激光掃描模態測試數據的可適用性。結果表明測試數據中的噪聲被有效地消除,同時降噪后的頻響函數可用于模態參數的準確提取。

1 基于模態峰值漢克爾奇異值分解的降噪方法與參數識別

基于模態峰值漢克爾奇異值分解(Modal Peak-based Hankel Singular Value Decomposition (MPHSVD))的降噪方法是采用了模態峰值分量信號選取方法的漢克爾奇異值分解濾波器,可用于單通道脈沖響應函數或頻率響應函數信號的降噪處理,降噪處理后的頻率響應函數可以直接通過模態參數識別方法,獲得準確的模態參數。

1.1 脈沖響應函數的漢克爾奇異值分解

它所對應的漢克爾(Hankel)矩陣為

式(3)中,為Hankel矩陣,為Hankel矩陣的行數,為Hankel矩陣的列數,同時+-1=。為了使Hankel矩陣近似為一個方陣,的大小設置為(/2),‘’的意思為向下取整。設置Hankel矩陣接近方陣的目的是為了達到HSVD最大的降噪效果。在形成Hankel矩陣后,矩陣將被SVD分解為分量矩陣

1.2 分量信號累加迭代

對包含了一系列模態信息及部分噪聲的IRF信號來說,漢克爾奇異值分解可將其分解為一系列按奇異值大小從大到小排列的分量信號,IRF信號中的某個模態信息可分解為一對或多對相鄰的單頻率成分的分量信號。通常,某階模態信息分解得到的分量信號所對應的奇異值的大小要大于大部分噪聲信號所對應的奇異值,因此,若按照奇異值的大小順序迭代累加分量信號,即可在一個未知的迭代次數恢復IRF信號中所有關心的模態信息,并且通過IRF的傅里葉變換在FRF中恢復所關心的各階模態。迭代累加過程中舍棄掉的分量信號,為大部分低能量的噪聲。

式中,E代表了所對應的能量譜。之后,計算模態峰值幅度和能量譜E在峰值頻率處的幅度差異

累加迭代將會停止,否則,迭代將繼續進行。在累加結束后,所關心的模態特性被恢復,并記錄選取的分量信號個數為。

該步驟會去除大部分低能量的高斯噪聲和少部分諧波噪聲,然而,一些高能量噪聲所對應的分量信號在分量信號序列中的位置要高于某些模態峰值頻率成分所對應的分量信號,因此,在累加迭代中,這些高能量噪聲的恢復要早于某些模態峰值頻率成分,進而會產生一定的殘余噪聲,因此,需要進行二次濾波去除殘余的高能量噪聲。

1.3 基于模態峰值頻率通帶的分量信號迭代選取

對累加迭代所恢復的殘余高能量噪聲,其存在形式大多形如諧波分量,所對應的分量信號也近似為單頻率成分信號。進而,由于高能量噪聲和模態峰值頻率成分所對應的分量信號為單頻率成分信號,故均可通過其在頻域能量譜中峰值振幅所對應的峰值頻率指代該分量信號。同時,模態峰值頻率成分所對應的分量信號的峰值頻率會在模態峰值頻率附近,而高能量噪聲對應的分量信號的峰值頻率則會在一定程度上遠離模態峰值頻率,因此,若在模態峰值頻率附近設置一系列的頻率通帶,即可通過選取峰值頻率在通帶內部的分量信號用于重構,并舍棄不在通帶內部的分量信號以去除殘余高能量噪聲。

由于累加迭代共選取了個分量信號,迭代次數為的迭代選取被用于去除殘余噪聲,迭代中的第次迭代如下。帶寬的大小靠經驗選取,同時通過后續的試驗數據分析表明帶寬不能選取太大,也不能選取太小,5Hz或10Hz的帶寬便足夠濾波器保留模態特性并去除噪聲。

1.4 基于有理分式多項式擬合的模態辨識

2 高溫結構激光掃描測量數據降噪和模態參數辨識

2.1 試驗裝置和測試數據分析

本文所采用的試驗裝置為一端固支的直板葉片,材料為45鋼。直板葉片尺寸為200mm×60mm×3mm。左下角帶有一個3mm的通孔,用于與激振器相連,孔中心距自由端底部6mm,左邊12mm。固支端是一個60×60×7mm的鋼板,共有四個10mm的通孔分布在固支端的四周,用于固定,其孔中心距相鄰兩邊的距離均為12mm,固支端與自由端之間采用2mm的圓弧接觸連接在一起。試驗時,固支端被四個螺栓固定在一鋼板上,鋼板與試驗臺采用螺栓鎖緊以模擬剛性支撐。直板葉片后方固定有一連有電磁加熱器的加熱板,以用于加熱直板來模擬高溫環境。激振器采用懸掛支撐,激振器的頂針和直板葉片之間通過一碳桿透過加熱板相連,以保證在高溫環境下激振器可正常工作,并將高溫環境與激振器隔離開來。采用非接觸掃描多普勒激光測振儀(Scanning Laser Doppler Vibrometer(SLDV))進行測量,獲取激勵的響應信號。試驗的測點為通過激光對焦得到的9×5個測點,如圖1(a)所示,這些測點均勻在直板葉片的自由端四周。試驗激勵為電腦控制的偽隨機激勵,采樣頻率設置為6400Hz,并通過激勵點和響應點的測試得到測點的頻率響應函數。

圖1 測點及結構表面情況

測量時用電磁加熱器加熱,熱電偶測溫。當直板葉片加熱到600度時,鋼板出現了圖1(b)所示的明顯的紅熱現象,同時可發現由于加熱不均勻,鋼板左右兩側的溫度明顯高于中間及上下兩端。進而通過檢查600度下的頻響函數數據可發現,溫度較高的左右兩端及其邊緣的部分測點出現了嚴重的噪聲,其余部分測點則噪聲較小。推測噪聲產生的原因為較高的溫度改變了結構的表面反射情況,使測試得到的頻響函數帶有一定的散斑噪聲,因為散斑噪聲經常被認為是激光測振中出現的主要噪聲[10],而左右兩側極高的溫度同測試用到的紅色激光發生了干涉,進一步增強了測試中存在的散斑噪聲,使得該區域部分測點的噪聲嚴重。本文選取了圖1(a)所示包含極大噪聲的43號測點和包含較小噪聲的38號測點的600度頻響函數進行分析,頻響函數的對比見圖2。

圖2 測試頻響函數的對比

圖2 (a)和(b)給出了頻響函數的對比圖像。圖2 (a)中,淡藍色曲線是測試得到的43號測點頻響函數,黑色曲線為38號測點的頻響函數,可以明顯地觀察到43號測點的測試數據包含極大的噪聲,并且多數噪聲的振幅要高于模態頻率的振幅。圖2 (a)中紅色和藍色圓環所包含的峰值是所關心的模態峰值,紅色圓環內的模態峰值是在38號測點和43號測點均能觀察到的模態峰值,而藍色圓環是只能在38號測點中觀察到的模態峰值,43號測點相應的模態峰值已經完全湮沒在噪聲當中,故無法恢復。從圖2 (b)中可以發現,38號測點頻響函數相位譜包含一定程度得噪聲,同時,相比于38號測點的數據,可發現43號測點的大部分相位譜已湮沒在噪聲當中,無法分清各階模態的真實相位。

2.2 測量數據降噪處理

在累加迭代結束后,將進行基于模態峰值頻率通帶的迭代選取,選取峰值頻率在設置的頻率通帶內的分量信號做進一步的降噪處理。考慮到38號測點共有8個關心的模態峰值頻率而43號測點共有4個關心的模態峰值頻率。8個在模態峰值頻率附近帶寬為10Hz的頻率帶被設置為38號測點的頻率通帶。相應的,4個帶寬為10Hz的頻率帶被設置為43號測點的頻率通帶。下面,進行迭代選取以分離出屬于模態峰值的分量信號,以模態5為例,屬于38號測點頻響函數模態5的分量信號選取方法如圖4所示。

圖4 模態分量信號選取方法

如圖4所示,模態5的模態峰值頻率為1044Hz,根據該頻率設置的頻率通帶為[1039Hz,1049Hz],如圖4中黑色虛線所示。該模態主要分解為了圖4中紫色和綠色曲線所代表的分量信號57和59,同時,該模態還有小于百分之1的能量混疊在了一峰值頻率為531.5的高能量噪聲中。分量信號57和59的峰值頻率為1044Hz,在黑色虛線所圍成的通帶區間內部,進而被選為屬于模態5的分量信號以用于恢復降噪后IRF信號,而由于分量信號55所代表的高能量噪聲峰值頻率為531.5Hz,不在頻率通帶內,分量信號55即被舍棄。

對38號測點在累加迭代中選取得到的59個分量信號,按照圖4所示的方法共選取出了36個屬于模態特性的分量信號,對43號測點在累加迭代中選取得到的181個分量信號,則選取出了10個屬于模態特性的分量信號。這些選取得到的屬于模態特性的分量信號被進一步加合在一起以構成降噪后的IRF信號,進而被轉換至頻域以獲取相應的降噪后FRF信號。降噪后FRF信號與測量FRF信號的對比分析見圖5。圖5 (a)和(b)給出了38號測點降噪后頻響函數的對比圖像,可以發現,圖像上所能觀察到的峰值只有所關心的模態峰值,FRF中存在的噪聲和不關心的峰值被很好的濾除掉了。圖5 (c)給出了43號測點的降噪前后FRF的對比,可以發現,FRF信號中大部分噪聲都被去除掉了,但仍剩余一些能量遠低于模態特性的噪聲,這是由于噪聲較大,部分噪聲混疊在了模態特性所對應的分量信號中,進而隨著模態特性的恢復而恢復所造成的,但這些未能去除的噪聲相比于模態特性能量較小,已不影響進一步的模態分析。圖5 (d)給出了降噪前后FRF相位譜的對比。相比于降噪前的FRF相位譜,相位譜中大部分噪聲被去除掉了,同時屬于四個模態峰值的相位突變被很好的恢復,同理,由于存在小能量峰值,FRF信號的能量譜并不光滑,并存在部分從180至-180的相位突變。從降噪前后的FRF對比來看,降噪保留了所有關心的模態特性,并去除了大部分噪聲,這證明了方法的有效性和優越性。

圖5 降噪后信號對比分析

2.3 降噪后數據的模態辨識

在完成了降噪處理后,要對降噪后的數據進行模態辨識,本文選取的模態辨識方法為RFP法,該方法對降噪后得到的FRF信號進行擬合以獲取相應的模態參數,降噪后FRF信號的擬合是分段處理的。擬合所選取的頻段和擬合結果見圖6所示。圖6 (a)和(b)給出了38號測點降噪后FRF擬合結果,可以觀察到曲線與降噪后曲線吻合度很高。圖6 (c)和(d)給出了43號測點降噪后FRF信號的擬合結果,同樣可以觀察到降噪后曲線和擬合曲線有較高的吻合度,同時可以發現,擬合頻段內能量較小的噪聲并未對擬合結果造成影響。

在完成了曲線擬合后,將從擬合得到的曲線中提取模態參數,提取的參數為模態頻率和模態阻尼,并計算了兩測點模態參數之間的差異,提取和計算結果見表1所示。表中,模態頻率差異計算的是提取前后模態頻率的絕對差,模態阻尼比的差計算如下

式中,為模態阻尼比的差異,為由38號測點提取得到的降噪后FRF信號模態阻尼比,則是由43號測點提取得到的降噪后FRF信號的模態阻尼比。由表1知,由38號測點和43號測點降噪后頻響函數提取出來的模態頻率最大差異為1.6Hz,模態阻尼比的最大差異約為16%,表明模態頻率和模態阻尼比的提取是準確的,驗證了降噪手段的有效性和優越性。

表1 提取得到的模態參數

3 結論

針對高溫結構激光掃描模態測試中出現的噪聲問題,提出了一種基于模態峰值漢克爾奇異值分解的降噪和模態識別方法。該方法首先通過分量信號累加迭代恢復所關心的模態信息,同時會一并恢復一定的殘余噪聲,進而通過基于模態頻率通帶的迭代選取分離出屬于模態特性的分量信號,并將其累加實現2次降噪處理。本文采用600度高溫環境下一直板葉片結構的激光掃描模態測試數據進行降噪和模態參數識別,結果表明,本方法能有效地去除頻響函數中的強噪聲,并保留結構的模態特性,準確識別結構的模態參數。該方法對航空、航天工程中的高溫模態試驗數據的分析有指導和借鑒作用。

[1] Golafshan R, Yuce Sanliturk K. SVD and Hankel matrix based de-noising approach for ball bearing fault detection and its assessment using artificial faults[J]. Mechanical Systems and Signal Processing, 2016(70-71): 36-50.

[2] Fang H-T, Huang D-S. Noise reduction in lidar signal based on discrete wavelet transform [J]. Optics Communications, 2004, 233(1-3): 67-76.

[3] Rong H, Gao Y , Guan L, et al. GAM-Based Mooring Alignment for SINS Based on An Improved CEEMD Denoising Method[J]. Sensors, 2019, 19(16): 3564.

[4] Sanliturk K Y, Cakar O. Noise elimination from measured frequency response functions[J]. Mechanical systems and signal processing, 2005, 19(3): 615-631.

[5] Zhao X, Ye B. Selection of effective singular values using difference spectrum and its application to fault diagnosis of headstock[J]. mechanical systems & signal processing, 2011, 25(5): 1617-1631.

[6] Bao X X, Li C L, Xiong C B. Noise elimination algorithm for modal analysis [J]. Applied Physics Letters, 2015, 107(4).

[7] Zhao M, Jia X . A novel strategy for signal denoising using reweighted SVD and its applications to weak fault feature enhancement of rotating machinery[J]. Mechanical Systems & Signal Processing, 2017, 94: 129-147.

[8] Zhao X, Ye B. Separation of Single Frequency Component Using Singular Value Decomposition[J]. Circuits, Systems, and Signal Processing, 2018, 38(1): 191-217.

[9] Formenti, DaveRichardson, et.al. Parameter estimation from frequency response measurements using rational fraction polynomials (twenty years of progress)[C]. Proceedings of SPIE - The International Society for Optical Engineering, 1982.

[10] MartinP, RothbergS. Introducing speckle noise maps for Laser Vibrometry[J]. Optics and Lasers in Engineering, 2009, 47 (3-4): 431-442.

[11] 王智勇, 王則力, 宮文然, 等. 熱結構高溫應變光學測量技術發展探討[J]. 強度與環境,2019,46(6): 1-8. [WANG Zhi-yong, WANG Ze-li, GONG Wen-ran, et al. Discussion on the Development of High Temperature Strain Optical Measurement Technology for Thermal Structure[J]. Structure & Environment Engineering, 2019,46(6): 1-8.]

[12] 白志富, 張忠, 朱儀凡, 等.高溫動力學模型誤差定位與修正方法研究[J]. 強度與環境,2019,46(2):21-26.[BAI Zhi-fu, ZHANG Zhong, ZHU Yi-fan, et al. Discussion on the Development of High Temperature Strain Optical Measurement Technology for Thermal Structure[J]. Structure & Environment Engineering, 2019,46(2): 21-26.]

Noise Reduction and Parameter Identification of Modal Test on a High Temperature Structure Measured by Scanning Laser Doppler Vibrometer

ZHU Tian-xu ZANG Chao-ping

(College of Energy and Power Engineering, Nanjing University of Aeronautics and Astronautics, Nanjing 210016, China)

The noise reduction and modal identification method of the modal peak based Hanke lSingular Value Decomposition(HSVD) is proposed in this paper to deal with the noise contaminated in modal testing of structures in the high temperature environment using a non-contact Scanning Laser Doppler Vibrometer (SLDV). Firstly, the measured frequency response function (FRF) signal is transformed to the impulse response function (IRF) by the inverse Fast Fourier Transform (iFFT), and the IRF signal is subsequently decomposed to a set of component signals using the HSVD method. These component signals are ordered according to the singular values from large to small. Secondly, the component signals are sequentially accumulated together with the ordered amplitudes of singular values until the all concerned modal modes are recovered and the other rest abandoned component signals are treated to be the noise. This step can eliminate most parts of the noise contaminated in the signal, but some residue noise is still inevitably recovered. Thirdly, the further iterative noise reduction based on the narrow frequency pass bands of the modal mode is conducted to separate the component signals of the modal features from the noise and the final deposing IRF signal is transformed to the FRF signal. Fourthly, the modal analysis is conducted to the noise reduced FRF signal to extract the modal parameters. The method is applied to measured data from modal testing of a straight plate at 600 centigrade degrees. Results shows that the strong noise contaminated in FRFs is effectively removed and modal parameters are precisely extracted. It demonstrates the validity and superiority of the method.

high temperature; noncontact measurement; the modal peak; HankelSingularValueDecomposition; noise reduction

V414.1

A

1006-3919(2020)02-0046-10

10.19447/j.cnki.11-1773/v.2020.02.008

2019-11-20;

2020-01-17

國家自然科學基金委員會與中國工程物理研究院聯合基金(U1730129)

朱天煦(1998—),男,本科生,研究方向:結構動力學;(210016)南京航空航天大學能源與動力學院.

猜你喜歡
模態信號方法
信號
鴨綠江(2021年35期)2021-04-19 12:24:18
完形填空二則
基于FPGA的多功能信號發生器的設計
電子制作(2018年11期)2018-08-04 03:25:42
用對方法才能瘦
Coco薇(2016年2期)2016-03-22 02:42:52
基于LabVIEW的力加載信號采集與PID控制
國內多模態教學研究回顧與展望
四大方法 教你不再“坐以待病”!
Coco薇(2015年1期)2015-08-13 02:47:34
捕魚
基于HHT和Prony算法的電力系統低頻振蕩模態識別
由單個模態構造對稱簡支梁的抗彎剛度
計算物理(2014年2期)2014-03-11 17:01:39
主站蜘蛛池模板: 久久精品欧美一区二区| 国产极品美女在线| 91在线播放免费不卡无毒| 国产精品成人第一区| 巨熟乳波霸若妻中文观看免费 | 国产美女91呻吟求| 毛片免费视频| 日韩AV无码一区| 91亚洲视频下载| 亚洲精品无码抽插日韩| AⅤ色综合久久天堂AV色综合| 超碰色了色| 99久久精品免费观看国产| 亚洲国产精品日韩专区AV| 日韩在线成年视频人网站观看| 国产福利在线免费| 亚洲精品午夜无码电影网| 中国丰满人妻无码束缚啪啪| 国产精品亚洲一区二区三区z| 99久久国产精品无码| 国产成人精品一区二区不卡| 亚洲人在线| 国产成人亚洲无码淙合青草| av在线5g无码天天| 成人免费黄色小视频| 呦系列视频一区二区三区| 亚洲第一香蕉视频| 精品国产欧美精品v| 欧洲日本亚洲中文字幕| 亚洲中文精品人人永久免费| 青草国产在线视频| 日日碰狠狠添天天爽| 午夜欧美理论2019理论| 免费观看国产小粉嫩喷水| 国产成人一区免费观看| 国产制服丝袜无码视频| 日韩欧美国产中文| 日韩不卡高清视频| 国产午夜在线观看视频| 精品少妇人妻无码久久| 国产人碰人摸人爱免费视频| 99在线视频免费观看| 亚洲国产精品一区二区第一页免| 精品国产福利在线| 国产精品开放后亚洲| 欧美精品三级在线| 国产人人射| 91麻豆精品国产高清在线 | 国产在线自揄拍揄视频网站| 国产伦片中文免费观看| 亚洲av无码片一区二区三区| 国产真实乱子伦精品视手机观看| 日韩无码视频网站| 波多野结衣第一页| 国产精品成人免费视频99| 欧美日韩国产在线播放| 国产高清无码第一十页在线观看| 午夜不卡视频| 国产成人精品高清在线| 亚洲五月激情网| 在线国产三级| 老司机aⅴ在线精品导航| 99久久精品无码专区免费| 欧美午夜在线观看| 国产精品国产主播在线观看| 丁香五月亚洲综合在线| 欧美va亚洲va香蕉在线| 国产第一色| 免费看美女自慰的网站| 亚洲视频无码| 亚洲综合精品第一页| 亚洲欧美综合精品久久成人网| 精品色综合| 国产综合无码一区二区色蜜蜜| 亚洲手机在线| 青青青国产视频手机| 亚洲综合18p| 国产精品yjizz视频网一二区| 伊人查蕉在线观看国产精品| 成人午夜网址| 亚洲中文字幕无码爆乳| 国产成人亚洲无码淙合青草|