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

基于諧波小波的電力系統(tǒng)諧波分析

2009-05-12 03:14:34張?zhí)扈?/span>
現(xiàn)代電子技術(shù) 2009年5期
關(guān)鍵詞:電力系統(tǒng)

摘 要:電力系統(tǒng)中的諧波對(duì)電網(wǎng)危害巨大,對(duì)其進(jìn)行監(jiān)測(cè)和分析就顯得非常重要。在諧波小波以及諧波小波包的基礎(chǔ)上,提出諧波小波變換的表達(dá)式以及諧波小波算法,給出電力系統(tǒng)諧波分析的仿真示例。仿真結(jié)果表明,利用諧波小波變換分解,并通過(guò)最小二乘法擬合出的各次諧波頻率和幅度的誤差率完全符合諧波分析的精度要求。在電力系統(tǒng)諧波的分析中,諧波小波算法具有其他算法無(wú)可比擬的優(yōu)越性。

關(guān)鍵詞:諧波小波;諧波分析;電力系統(tǒng);諧波;間諧波;最小二乘法擬合

中圖分類號(hào):TM711文獻(xiàn)標(biāo)識(shí)碼:A

文章編號(hào):1004-373X(2009)05-159-04

Harmonic Analysis in Power System Based on Harmonic Wavelet

ZHANG Tianyu

(Wuxi Radio & Television University,Wuxi,214011,China)

Abstract:Harmonics in power system are harmful to power network.It is significant to detect and analyze them.On the basis of harmonic wavelet and harmonic wavelet packet,the expression of harmonic wavelet transform and the algorithm of harmonic wavelet are proposed and then a simulation example of harmonic analysis in power system is given.The simulation results show that the error rate of each harmonic frequency and amplitude corresponds to the accuracy requirement of harmonic analysis completely with the decomposition of harmonic wavelet transform and least square method fitting.The algorithm of harmonic wavelet is superior to the others in harmonic analysis of power system.

Keywords:harmonic wavelet;harmonic analysis;power system;harmonics;inter-harmonics;least square method fitting

0 引 言

由于電力系統(tǒng)中大量非線性設(shè)備的存在,導(dǎo)致它們?cè)诠ぷ鲿r(shí)不僅會(huì)產(chǎn)生基波頻率的整數(shù)次諧波,還可能產(chǎn)生基波頻率的非整次諧波,即間諧波,這會(huì)對(duì)電能造成嚴(yán)重的污染,增加能量損失,威脅電力設(shè)備的安全運(yùn)行[1-4]。因此,諧波和間諧波的分析對(duì)于電力系統(tǒng)的監(jiān)控與保護(hù)都具有十分重要的意義。

傳統(tǒng)的正交小波包變換在電力系統(tǒng)諧波分析與檢測(cè)中有著廣泛的應(yīng)用。但是由于小波包變換固有的性質(zhì),如小波包變換的混疊現(xiàn)象比小波變換的混疊現(xiàn)象更為直觀形象,其影響也比小波變換嚴(yán)重,這主要是由于分解濾波器之間存在頻帶混疊現(xiàn)象,小波頻譜的起始頻率和截止頻率之間存在過(guò)渡帶[5]。諧波小波變換是一種基于快速傅里葉變換(Fast Fourier Transform,FFT)及其逆變換(Inverse Fast Fourier Transform,IFFT)的快速算法,在數(shù)值上容易實(shí)現(xiàn),其算法快,精度高,具有很好的工程實(shí)用價(jià)值[6-8]。通常的小波算法(如Mallat算法,Daubechies小波)在分解信號(hào)時(shí)要隔二取一,從而使得在小波分解時(shí)各層的數(shù)據(jù)點(diǎn)數(shù)和采樣頻率隨尺度的增加逐漸減小。諧波小波相對(duì)于傳統(tǒng)的小波函數(shù)而言,具有更普遍意義上的正交性以及優(yōu)異的視頻分解能力,其明顯優(yōu)勢(shì)就是信號(hào)任意頻段的“細(xì)化”能力,雖然它在時(shí)域中的局部化能力一般,但在頻域分析中對(duì)精度有特殊要求的場(chǎng)合,這種優(yōu)勢(shì)就非常符合需求[9,10]。

1 諧波小波分析

1.1 經(jīng)典諧波小波

設(shè)時(shí)域函數(shù)h璭(t)和h璷(t)的傅里葉變換所對(duì)應(yīng)的頻域函數(shù)為璭(ω)和璷(ω),它們的表達(dá)式見(jiàn)式(1):

璭(ω)=1/(4π), ω∈[-4π,-2π]∪[2π,4π]

0,其他

璷(ω)=i/(4π), ω∈[-4π,-2π]

-i/(4π),ω∈[2π,4π]

0,其他

(1)

式中:下標(biāo)e和o分別表示該函數(shù)是變量ω的偶函數(shù)和奇函數(shù)。

將頻域函數(shù)璭(ω)和璷(ω)組成復(fù)合函數(shù)(ω),可得:

(ω)=璭(ω)+i璷(ω)=1/(2π),ω∈[2π,4π]

0,其他

(2)

(ω)具有良好的緊支撐特性和盒形特征。對(duì)式(1)作廣義的傅里葉逆變換(忽略系數(shù)1/(2π)),可得:

h璭(t)=∫∞-∞璭(ω)exp(iωt)dω=sin (4πt)-sin (2πt)2πt

h璷(t)=∫∞-∞璷(ω)exp(iωt)dω=cos (2πt)-cos (4πt)2πt

(3)

將時(shí)域函數(shù)h璭(t)和h璷(t)組成復(fù)合函數(shù)h(t),可得:

h(t)=h璭(t)+ih璷(t)=exp(i4πt)-exp(i2πt)i2πt

(4)

由此定義的復(fù)合函數(shù)h(t)稱為諧波小波函數(shù),亦稱為經(jīng)典諧波小波或二進(jìn)諧波小波,其實(shí)部h璭(t)和虛部h璷(t)的波形如圖1所示。

圖1 諧波小波h(t)實(shí)部與虛部的波形

由圖1可以看出,諧波小波h(t)是由相差90°的實(shí)部偶小波和虛部奇小波構(gòu)成。虛部奇小波所構(gòu)成的濾波器都是零相移濾波器,具有鎖定信號(hào)相位的功能。它在時(shí)域上的衰減速度較慢(與時(shí)間t成反比),導(dǎo)致其時(shí)域局部化特性較弱。

為了獲得諧波小波h(t)的二進(jìn)伸縮平移系,令:

t=2jt-k

(5)

式中:j為非負(fù)整數(shù);k為整數(shù)。

把式(5)代入式(4),可得:

h(2jt-k)=exp[i4π(2jt-k)]-exp[i2π(2jt-k)]i2π(2jt-k)

(6)

在式(6)中,小波的形狀沒(méi)有改變,只是在水平尺度上被壓縮了2j,并且位置在新的尺度上被平移了k個(gè)單位,這與二進(jìn)小波變換的形式是一致的。其j值決定諧波小波的尺度或?qū)訑?shù)。例如當(dāng)j=0時(shí),諧波小波的傅里葉變換位于[2π,4π]頻帶中;若在第j層時(shí),則諧波小波的傅里葉變換位于[2j+1π,2j+2π]頻帶之間。即隨著j值的增大,其頻譜的帶寬以二進(jìn)方式逐漸加大。諧波小波對(duì)信號(hào)的分解從低頻到高頻是以2倍的關(guān)系逐漸增加的,它對(duì)信號(hào)的低頻部分劃分比較細(xì),而對(duì)信號(hào)的高頻部分劃分比較粗,這說(shuō)明經(jīng)典諧波小波分解也屬于二進(jìn)小波分解的范疇。

1.2 諧波小波的改進(jìn)

為了使分析頻帶的選取更為靈活,不受二進(jìn)方式的限制,對(duì)經(jīng)典諧波小波加以改進(jìn),拓寬諧波小波的概念及應(yīng)用范圍。引入正整數(shù)m=2j,n=2j+1(m<n),把m,n代入式(4),并通過(guò)伸縮平移生成的諧波小波族為:

Ψ璵,n(t)=exp(in2πt)-exp(im2πt)i2π(n-m)t

(7)

其頻域表達(dá)式為:

璵,n(ω)=F[Ψ璵,n(t)]=12π(n-m),

ω∈[2πm,2πn]0,其他

(8)

由式(7)可以看出,實(shí)際上m,n既可以取正整數(shù),也可以取負(fù)整數(shù),這樣它們之間就不必滿足n=2m這一條件的限制(二進(jìn)限制),只要保證m<n即可,這就使得諧波小波在分析頻帶的選取上具有更大的靈活性。這就是改進(jìn)的諧波小波相對(duì)于經(jīng)典諧波小波的一個(gè)明顯優(yōu)勢(shì)。

若給定諧波小波的位移步長(zhǎng)為k/(n-m),k為整數(shù),對(duì)式(7)進(jìn)行平移變換可得:

Ψ璵,nt-kn-m=

expin2πt-kn-m-expim2πt-kn-mi2π(n-m)t-kn-m

(9)

其頻域表達(dá)式為:

璵,n[(n-m)ω]=FΨ璵,nt-kn-m=

12π(n-m)exp-jωkn-m,

ω∈[2πm,2πn]0,其他

(10)

由此可見(jiàn),式(10)是分析頻率帶寬為(n-m)2π,分析時(shí)間中心在t=k/(n-m)處的諧波小波一般表達(dá)式。文獻(xiàn)[12]證明了諧波小波族Ψ璵,n(t)是一個(gè)正交的解析信號(hào),它構(gòu)成了空間L2(R)的一組正交基[13]。

1.3 諧波小波包

由式(9)可知,諧波小波的關(guān)鍵在于尺度參數(shù)m,n的選取。令信號(hào)的奈奎斯特頻率為f璼,則第j(j為非負(fù)整數(shù))層各小波的分析頻率帶寬為:

B=f璼/2j+1

(11)

這樣可以設(shè)定分析頻帶的上、下限頻率分別為:

m=rB,n=(r+1)B,

r=0,1,2,…,2j-1 (12)

隨著分解層數(shù)j的逐漸增大,可以體現(xiàn)出諧波小波包對(duì)信號(hào)任意頻段的“細(xì)化”能力。如果要對(duì)信號(hào)的某一頻段進(jìn)行重點(diǎn)分析,則先由式(11)確定信號(hào)的分解層數(shù)j,再由式(12)確定所要分析頻帶的上、下限頻率,也就是定義諧波小波的尺度參數(shù)m,n。

由于諧波小波沒(méi)有尺度函數(shù),因此諧波小波包的思想與傳統(tǒng)的小波包理論有所不同,不能采用正交濾波器組對(duì)信號(hào)進(jìn)行頻帶分解[14]。由式(9)可知,諧波小波具有可調(diào)的尺度參數(shù)m,n,對(duì)在不同頻帶的信號(hào)進(jìn)行分解時(shí)采用不同的m,n,這樣就可以將諧波小波良好的濾波效果應(yīng)用到諧波小波包的分析中。信號(hào)經(jīng)過(guò)小波包分解后,在各個(gè)頻帶中的信號(hào)仍具有與原始信號(hào)相同的頻率分辨率,而且分解后信號(hào)的數(shù)據(jù)長(zhǎng)度并沒(méi)有減少,這克服了Mallat算法的小波包分解帶來(lái)數(shù)據(jù)長(zhǎng)度減少的問(wèn)題。由于小波濾波器不具有理想“盒形”的頻譜特性,起始頻率和截止頻率之間存在過(guò)渡帶,這導(dǎo)致在信號(hào)的分解過(guò)程中往往會(huì)發(fā)生頻帶間的能量冗余,造成誤差,而諧波小波包濾波器則完全可以克服以上問(wèn)題。具體方法是首先得到待分析信號(hào)的頻譜,確定譜線的頻點(diǎn)數(shù)值,然后根據(jù)預(yù)設(shè)的窗寬來(lái)確定尺度參數(shù)。

2 諧波小波變換及算法

2.1 諧波小波變換

根據(jù)小波變換的定義,對(duì)某一尺度的小波函數(shù)Ψ璵,n(t),信號(hào)x(t)∈L2(R)的小波變換可表示為[15]:

W(t)=∫∞-∞x(τ)Ψ*璵,n(t-τ)dτ

(13)

信號(hào)x(t)的諧波小波變換為:

W(m,n,k)=(n-m)∫∞-∞x(t)Ψ*璵,nt-kn-mdt

(14)

對(duì)式(14)進(jìn)行Fourier變換,可得:

(m,n,ω) =(ω)*璵,n [(n-m)ω]

(15)

式(14)和式(15)分別稱作信號(hào)x(t)在m,n尺度下的時(shí)域和頻域的諧波小波變換表達(dá)式。

對(duì)于離散信號(hào)序列x(r),r=0,1,2,…,N-1,其諧波小波變換為:

W(m,n,k)=n-mN∑N-1r=0x(r)Ψ*璵,nr-kn-m

(16)

由式(13)~式(16)可以看出,信號(hào)的諧波小波變換非常簡(jiǎn)潔,容易實(shí)現(xiàn)。同時(shí),由于諧波小波對(duì)信號(hào)各次諧波分量的相位有保持功能,所以對(duì)信號(hào)進(jìn)行諧波小波分解后,也可以對(duì)信號(hào)進(jìn)行重構(gòu),從而實(shí)現(xiàn)信號(hào)的濾波和降噪。

2.2 諧波小波算法

首先對(duì)諧波源信號(hào)x(t)進(jìn)行FFT運(yùn)算,對(duì)變換得到的結(jié)果(ω)進(jìn)行頻率搜索,以確定諧波小波的尺度參數(shù)m璲,n璲,進(jìn)而確定諧波小波函數(shù)h璵璲,n璲(t),然后將諧波小波函數(shù)h璵璲,n璲(t)進(jìn)行FFT運(yùn)算的結(jié)果璵璲,n璲(ω)與(ω)相乘,再對(duì)其相乘的結(jié)果(m璲,n璲,ω)進(jìn)行IFFT運(yùn)算,通過(guò)對(duì)時(shí)域的小波系數(shù)W(m璲,n璲,t)進(jìn)行重構(gòu),得到各次諧波和間諧波的瞬時(shí)值,最后利用最小二乘法對(duì)各頻率分量進(jìn)行擬合,得到諧波小波分析的結(jié)果,其流程圖如圖2所示。

圖2 諧波小波算法的流程圖

3 仿真實(shí)驗(yàn)與結(jié)果分析

為了更好地驗(yàn)證諧波小波算法在電力系統(tǒng)諧波與間諧波分析中的有效性,進(jìn)行如下的仿真實(shí)驗(yàn)。

設(shè)電網(wǎng)中的諧波源信號(hào)為:

u(t)=∑6i=1A璱sin 2πf璱+e(t)

(17)

式中:基波頻率為50 Hz,并且含有3,5,7,9次諧波和頻率為75 Hz(基波頻率的1.5倍)的間諧波共6個(gè)頻率分量以及隨機(jī)噪聲e(t),具體的參數(shù)設(shè)置如表1所示。

表1 諧波源參數(shù)的設(shè)置

諧波源分量基波u1(t)間諧波u1.5(t)3次諧波u3(t)5次諧波u5(t)7次諧波u7(t)9次諧波u9(t)

頻率 f璱 /Hz5075150250350450

幅度 A璱 /V22010151085

設(shè)采樣頻率f璼=1 250 Hz,采樣點(diǎn)數(shù)N=1 024。利用諧波小波變換(Harmonic Wavelet Transform,HWT)對(duì)諧波源信號(hào)u(t)進(jìn)行分解,通過(guò)Matlab仿真得到分解后各頻率分量的波形如圖3所示。

圖3 HWT分解后諧波與間諧波的波形

由圖3可以看出,諧波源中的各次諧波和間諧波分量被分解到了不同的頻帶中,這表明利用諧波小波算法來(lái)實(shí)現(xiàn)電力系統(tǒng)諧波和間諧波信號(hào)的分離是完全有效的。下一步需要對(duì)分解出的各個(gè)頻帶分量進(jìn)行參數(shù)提取,以計(jì)算出各次諧波的頻率和幅值。

最小二乘法擬合是一個(gè)基于全局觀念的擬合方法,針對(duì)某一樣本數(shù)據(jù)集合,利用該方法可以求得該集合中的主流趨勢(shì)。利用最小二乘法對(duì)6個(gè)頻帶內(nèi)的諧波和間諧波分量進(jìn)行擬合,并且定義頻率和幅度的誤差率分別為:

δ璮=(f′璱-f璱)/f璱

(18)

δ瑼=(A′璱-A璱/A璱

(19)

其計(jì)算結(jié)果如表2所示。

由表2可以看出,利用HWT法分解并擬合出的各次諧波頻率的誤差率在10-4數(shù)量級(jí),幅度的誤差率在10-2數(shù)量級(jí),完全符合諧波分析的精度要求。由此可見(jiàn),HWT法在諧波頻率和幅值的檢測(cè)中具有非常明顯的優(yōu)勢(shì)。

表2 HWT法的計(jì)算結(jié)果

真實(shí)值HWT法

頻率f璱 /Hz幅度A璱 /V頻率

f′璱 /Hz頻率誤差率 δ璮 /%幅度

A′璱/V幅度誤差率 δ瑼 /%

5022049.999 7-6.000 0e-4219.963 7-1.650 0e-2

751075.000 56.666 7e-410.003 53.500 0e-2

15015149.999 7-2.000 0e-414.996 5-2.333 3e-2

25010249.999 3-2.800 0e-49.996 5-3.500 0e-2

3508349.999 1-2.571 4e-48.003 74.625 0e-2

4505450.000 92.000 0e-45.001 83.600 0e-2

4 結(jié) 語(yǔ)

將諧波小波引入電力系統(tǒng)的諧波分析中,首先闡述了經(jīng)典諧波小波及其改進(jìn)及諧波小波包的概念,接著利用推導(dǎo)出的諧波小波算法對(duì)電網(wǎng)中的諧波源信號(hào)進(jìn)行諧波參數(shù)提取。仿真結(jié)果表明,諧波小波變換可以快速有效地對(duì)電力系統(tǒng)中的電壓諧波以及間諧波進(jìn)行檢測(cè),并能準(zhǔn)確地分解出各次諧波分量。可以預(yù)計(jì),隨著諧波小波理論的不斷發(fā)展和完善,諧波小波變換必將在電力系統(tǒng)間的諧波分析中發(fā)揮更大作用。

參考文獻(xiàn)

[1]Kumar A,Das B,Sharma J.Robust Dynamic State Estimation of Power System Harmonics [J].International Journal of Electrical Power & Energy Systems,2006,28(1):65-74.

[2]Madtharad C,Premrudeepreechacharn S,Watson N R,et al.An Optimal Measurement Placement Method for Power System Harmonic State Estimation [J].IEEE Trans.on Power Delivery,2005,20(2):1 514-1 521.

[3]Lin H C.Remote Power System Harmonics Measurement and Monitoring via the Internet [J].International Journal of Electrical Engineering Education,2005,42(3):267-279.

[4]Kanao N,Yamashita M,Yanagida H,et al.Power System Harmonic Analysis Using State-Estimation Method for Japanese Field Data [J].IEEE Trans.on Power Delivery,2005,20(2):970-977.

[5]梁學(xué)章,何甲興,王新民,等.小波分析 [M].北京:國(guó)防工業(yè)出版社,2005.

[6]Chang G W,Chen C I,Liu Y J,et al.Measuring Power System Harmonics and Interharmonics by an Improved Fast Fourier Transform-Based Algorithm [J].IET Generation,Transmission & Distribution,2008,2(2):192-201.

[7]Lu S L.Application of DFT Filter Bank to Power Frequency Harmonic Measurement [J].IEE Proceedings-Generation,Transmission and Distribution,2005,152(1):132-136.

[8]Testa,A,Gallo D,Langella,R.On the Processing of Harmonics and Interharmonics:Using Hanning Window in Standard Framework [J].IEEE Trans.on Power Delivery,2004,19(1):28-34.

[9]Musha T,Kumazawa T.Instantaneous Structural Intensity by the Harmonic Wavelet Transform [J].Journal of Sound and Vibration,2007,306(1-2):377-388.

[10]Chancey V C,Flowers G T,Howard C L.A Harmonic Wavelets Approach for Extracting Transient Patterns from Measured Rotor Vibration Data [J].Journal of Engineering for Gas Turbines and Power,2003,125(1):81-89.

[11]Liu B.Adaptive Harmonic Wavelet Transform with Applications in Vibration Analysis [J].Journal of Sound and Vibration,2003,262(1):45-64.

[12]飛思科技產(chǎn)品研究中心.小波分析理論與Matlab7實(shí)現(xiàn) [M].北京:電子工業(yè)出版社,2005.

[13]Pham V L,Wong K P.Antidistortion Method for Wavelet Transform Filter Banks and Nonstationary Power System Waveform Harmonic Analysis [J].IEEE Proceedings-Generation,Transmission and Distribution,2001,148(2):117-122.

[14]Barros J,Diego R I.Application of the Wavelet-Packet Transform to the Estimation of Harmonic Groups in Current and Voltage Waveforms [J].IEEE Trans.on Power Delivery,2006,21(1):533-535.

[15]劉濤,曾祥利,曾軍.實(shí)用小波分析入門 [M].北京:國(guó)防工業(yè)出版社,2006.

作者簡(jiǎn)介 張?zhí)扈?男,1980年出生,江蘇無(wú)錫人,無(wú)錫市廣播電視大學(xué)機(jī)電工程系教師,講師,江南大學(xué)05控制工程碩士研究生。主要從事小波理論和電力系統(tǒng)智能監(jiān)測(cè)方面的研究工作。

猜你喜歡
電力系統(tǒng)
淺析繼電保護(hù)自動(dòng)化技術(shù)在電力系統(tǒng)中的應(yīng)用
電力系統(tǒng)變電運(yùn)行設(shè)備維護(hù)技術(shù)分析
并聯(lián)電抗器在高壓電網(wǎng)中的應(yīng)用
試論電力繼電保護(hù)運(yùn)行及可靠性問(wèn)題
信息時(shí)代的電力系統(tǒng)信息化經(jīng)驗(yàn)談
關(guān)于智能變電站繼電保護(hù)技術(shù)優(yōu)化探討
探討電力系統(tǒng)中配網(wǎng)自動(dòng)化技術(shù)
科技視界(2016年21期)2016-10-17 20:00:58
配網(wǎng)自動(dòng)化技術(shù)的應(yīng)用探討
科技視界(2016年21期)2016-10-17 19:54:47
一種新型高壓電纜頭制作方法
電力系統(tǒng)繼電保護(hù)技術(shù)常見(jiàn)問(wèn)題及對(duì)策
主站蜘蛛池模板: 亚洲国产成人久久精品软件| 精品亚洲国产成人AV| 理论片一区| 国产精品无码久久久久久| 亚洲精品无码不卡在线播放| 97国产在线播放| 一本无码在线观看| 中国毛片网| 免费国产高清精品一区在线| 精品黑人一区二区三区| 欧美色图久久| 免费一看一级毛片| 激情乱人伦| 在线国产毛片| 一级一毛片a级毛片| 性色生活片在线观看| 成人午夜天| 无码福利视频| 中文字幕 日韩 欧美| 国产一区二区影院| 精品久久久久久中文字幕女| 91久久偷偷做嫩草影院电| 国产偷国产偷在线高清| 国产一级在线观看www色 | 日本人妻一区二区三区不卡影院| 五月婷婷综合色| 国产成人一区二区| 国产一区二区三区在线观看视频| 精品视频一区二区观看| 国产成人三级| 国产另类乱子伦精品免费女| 国产在线高清一级毛片| 26uuu国产精品视频| 精品一区二区久久久久网站| 国产毛片高清一级国语| 欧美激情视频二区三区| 999国产精品永久免费视频精品久久 | 久久综合丝袜日本网| 午夜色综合| 青青草原国产av福利网站| 在线看AV天堂| 韩日无码在线不卡| 青青国产视频| 国产在线无码av完整版在线观看| 伊人色天堂| 中文字幕欧美成人免费| 欧美中文字幕在线视频| 91成人在线免费观看| 四虎亚洲精品| 亚洲人成影视在线观看| 91久久偷偷做嫩草影院免费看 | 试看120秒男女啪啪免费| 日韩精品一区二区三区大桥未久 | 免费一级毛片在线播放傲雪网| 亚洲黄色片免费看| 久久综合色天堂av| 欧美怡红院视频一区二区三区| 四虎亚洲国产成人久久精品| 在线精品亚洲一区二区古装| av在线无码浏览| 国产精品视频系列专区| 欧美在线精品怡红院| 国产va在线观看免费| 国产精品久久久久久久伊一| 91外围女在线观看| 日韩国产综合精选| 亚洲第一黄色网址| 五月天福利视频| 美女高潮全身流白浆福利区| 欧美曰批视频免费播放免费| 欧美中文字幕第一页线路一| 97国产精品视频自在拍| 99久久性生片| 色综合成人| 国产精品毛片一区视频播| 国产欧美日韩一区二区视频在线| 欧美午夜精品| 精品视频在线观看你懂的一区| 亚洲无码高清视频在线观看| 亚洲一区二区三区麻豆| 亚洲一级毛片在线观播放| 日韩a级片视频|