国产日韩欧美一区二区三区三州_亚洲少妇熟女av_久久久久亚洲av国产精品_波多野结衣网站一区二区_亚洲欧美色片在线91_国产亚洲精品精品国产优播av_日本一区二区三区波多野结衣 _久久国产av不卡

?

基于自適應層數(shù)分解的小波變換濾除沖擊信號

2015-11-30 18:27:57王鍇劉志國劉剛易龍濤陳曼彭詩棋
湖南大學學報·自然科學版 2015年10期
關鍵詞:小波變換沖擊振動

王鍇 劉志國 劉剛 易龍濤 陳曼 彭詩棋

摘要:針對摻雜在沖擊信號中的微弱振動信息被淹沒的問題,提出了一種基于小波變換濾除沖擊信號提取微弱振動信號的方法。該方法采用自適應層數(shù)分解的小波變換分析原理,根據(jù)沖擊信號特征,合理選擇小波基,對沖擊信號進行分析濾除。為獲得較好的濾除效果,針對分解層數(shù)的自適應性,提出了一種分解層數(shù)的自適應確定方法。對含有沖擊干擾的微弱標準信號進行數(shù)值模擬,對含有沖擊干擾的微弱金剛砂振動信號進行實測分析。結果表明該方法能自適應性地確定小波分解層數(shù),并有效濾除沖擊提取被淹沒的微弱振動信號。

關鍵詞:沖擊;振動;小波變換;頻率分析

中圖分類號:TP395。02 文獻標識碼:A

沖擊信號在機械振動中是一種常見的特征信號\[1-2\],具有持續(xù)時間短,衰減快,能量在廣闊的頻率范圍內(nèi)發(fā)散等特點,發(fā)生沖擊振動的同時往往還包含其他重要微弱振動信息,微弱的振動信號容易被沖擊信號淹沒,不易識別\[3\]。因此為了有效提取淹沒在沖擊信號中的微弱振動信號,濾除沖擊信號顯得尤為重要。

時頻分析是信號分析的傳統(tǒng)方法,頻譜分析是以傅里葉變換為基礎的經(jīng)典譜分析和非傅里葉變換分析為基礎的現(xiàn)代譜分析,時頻分析方法已經(jīng)廣泛應用于非平穩(wěn)信號分析。針對摻雜沖擊的典型非平穩(wěn)振動信號\[4-5\],小波變換具有優(yōu)越的時頻局部化分析和細化能力,改變了傳統(tǒng)傅里葉變換不能在時域局部化的缺陷,可以迅速聚焦到信號任意頻段細節(jié),準確提取特征信號進行分析。

湖南大學學報(自然科學版)2015年

第10期王鍇等:基于自適應層數(shù)分解的小波變換濾除沖擊信號

利用長時間進行低頻局部化處理,短時間進行高頻局部化處理,這種自適應窗口分析的能力改進了短時傅里葉分析方法中時頻分辨率固定的不足。根據(jù)時域卷積,頻域相乘原理和自適應的多分辨率特征,小波變換已經(jīng)成為研究非線性非平穩(wěn)信號的重要工具\[6-7\]。

譚繼勇等\[8\]針對沖擊信號的檢測問題,提出了一種沖擊信號的自適應隨機共振檢測算法,該算法將加權鞘度指標與變尺度隨機共振相結合,可以有效對微弱沖擊信號進行檢測,對提取沖擊信號特征的自適應算法有一定參考價值。滕軍等\[9\]針對常見的噪聲信號提出了一種基于中值濾波與小波閾值降噪相結合的分解層數(shù)自適應算法,通過確定各層閾值的取值方法,有效提取淹沒在噪聲中的信號。本文基于濾除沖擊提取有效微弱信號實驗為背景,針對小波在信噪分離,提取微弱信號和突變點信號優(yōu)勢\[10\],結合沖擊信號特征,合理選擇小波基函數(shù),提出了一種濾除沖擊有效提取淹沒信號的自適應層數(shù)分解的小波變換分析方法。

1小波分析濾除沖擊信號原理

1。1小波變換定義

信號或函數(shù)f(t)∈L(R),與小波基進行卷積計算,可得連續(xù)小波變換\[11\],其表達式為:

WTfa,b=〈f(t),ψa,b(t)〉=

1a∫Rf(t)ψt-badt=

a2π∫Rf(ω)Ψ(ω)ejωtdt。(1)

式中:a為尺度因子(a>0);b為時移因子,可正可負;WTf(a,b)為小波變換系數(shù);Ψ a,b(t)為小波母函數(shù);f(ω),Ψ(ω)分別為f(t),Ψ(t)的傅立葉變換。母函數(shù)在t=b處通過尺度a變換可得Ψt-τ/a,其變換原理如圖1所示。鏡頭相當于Ψ(t)起的作用。a相當于鏡頭向f(t)推近或遠離,推近時為高頻分析,遠離時為低頻分析。b相當于鏡頭相對f(t)平行移動。結合這2個參數(shù)的調(diào)整,就可以在不同的頻率尺度上進行不同的時間分析。

圖1小波變換作用原理

Fig。1Wavelet transform theory

利用尺度因子a和時移因子b,對連續(xù)小波函數(shù)Ψ a, b(t)進行離散化處理,可得離散小波函數(shù),并推導出相應的離散小波變換表達式,再利用時移因子對時間軸歸一化處理后,即得到相應的二進制離散小波變換。在具有時移不變和雙尺度特性的函數(shù)空間里,由多分辨率分析可得低通平滑尺度函數(shù)Φ(t)和高通細節(jié)小波函數(shù)Ψ(t)的雙尺度差分方程為:

Φ(t)=2∑k∈Zh0,kφ(2t-k); (2)

Ψ(t)=2∑k∈Zh1,kφ(2t-k)。 (3)

式中:h0,k為尺度系數(shù);h1,k為小波系數(shù);h0,k,h1,k為共軛正交濾波器系數(shù)。先通過尺度函數(shù)Φ(t)的傅里葉變換得出H(ω),再通過無窮乘積定義Φ(ω),使Φ(t-k)具有標準正交條件,同時利用數(shù)學家Daubechies推導的濾波器數(shù)組,即可構造Daubechies緊支撐小波函數(shù)。根據(jù)振動信號特性,合理選擇小波

基,使得與小波基函數(shù)相似的能量集中在高幅值的時間尺度上,與小波基函數(shù)不相似的能量發(fā)散到零幅值的時間尺度平面上,實現(xiàn)小波變換分析。

1。2小波分解層數(shù)自適應性確定

經(jīng)過限幅處理的原始信號,可以看成微弱振動信息與沖擊信號的疊加。在小波域中,隨著小波分解層數(shù)的增多,分解細節(jié)系數(shù)的模極大值即微弱振動信號模極大值及稠密程度增大,沖擊和噪聲信號的模極大值幅值及稠密程度減小,且個數(shù)基本保持不變\[9\]。通過模極大值的幅值即可確定小波分解后信號所占比重,從而確定小波分解的自適應層數(shù)。

對原始信號進行m(一般m≤6)層分解,獲取每層小波細節(jié)系數(shù)的模極大值序列W2,…,WN-1,如果|Wi|≥|Wi±1|,則i為信號f(t)該層分解的一個模極大值點。在對信號分解過程中,當滿足最大模極大值條件式(4)時,則認為原始信號中以微弱振動信號為主,否則原始信號中沖擊信號為主。

Wn≥Wn-1+k×σ。 (4)

式中:σ為隨機干擾的小波變換系數(shù)模的標準差,由隨機干擾信號的第2層小波分解的系數(shù)值確定,通過它判定局部極大值是否由隨機干擾引起,k為優(yōu)化參數(shù)(一般k=2),主要是調(diào)節(jié)檢測微弱振動信號的能力值大小,其值越小,微弱信號被提取出來的能力越強,同時模極大值受隨機干擾的影響越大。因此,k×σ項主要消除隨機干擾導致的局部極大值的增加,避免微弱振動信號對應的局部極大值點淹沒在隨機干擾產(chǎn)生的局部極大值中。小波分解層數(shù)自適應步驟如圖2所示。

1)對原始信號進行一層小波分解;

2)保留一層小波分解的尺度系數(shù),對細節(jié)系數(shù)重構,得到新的信號,同時判斷新信號的模極大值序列是否滿足式(4)的模極大值條件,若滿足則確定分解層數(shù)為1,結束自適應性層數(shù)分解,若不滿足則對原始信號進一步分解。

3)重復步驟2),直至分解出的細節(jié)系數(shù)重構信號的模極大值滿足式(4)。

4)若最后一次分解層數(shù)為m,則最終確定分解層數(shù)為m-1。

1。3小波分析濾除沖擊信號流程

在摻雜了沖擊的振動信號中,沖擊信號的幅值一般遠大于淹沒在沖擊信號中的微弱信號幅值,為了快速精準地獲取濾除沖擊后的微弱信號,先根據(jù)實際情況確定是否需要限定沖擊信號幅度,并根據(jù)實際經(jīng)驗限定沖擊信號幅度值,即確定濾波過程中的第1層分解的閾值。限幅后的信號通過自適應層數(shù)分解的小波分析后,可獲得濾除沖擊后的所需信號,其簡單信號處理如圖3所示。

圖2小波分解層數(shù)自適應步驟

Fig。2The step of adaptive layer decomposition

of wavelet transform

圖3濾除沖擊信號流程圖

Fig。3Filtering shock signal flow diagram

2小波濾除模擬沖擊信號

2。1小波基函數(shù)選取

小波基函數(shù)Ψj,k(t)具有多樣性,不同的小波函數(shù)Ψ(t)和相應的尺度函數(shù)Φ(t)可以構成不同的小波基,其作用相當于濾波器,采用不同的小波基函數(shù),會產(chǎn)生不同的信號分析結果。

沖擊信號具有非明顯的高斯分布特性,為了獲得較好的數(shù)字信號濾除效果,選取的小波基函數(shù)需滿足有限的支撐長度,較好的頻率分辨率。本文根據(jù)振動信號特征選取Daubechies小波對沖擊信號進行小波分析,由濾波系數(shù)可得N分別為2,4,12,20時的尺度函數(shù)和小波函數(shù),如圖4和圖5所示。

由圖4和圖5可以看出,隨著濾波器長度增加,即N值的增加,小波高頻系數(shù)減小,分解能量越來越集中,尺度函數(shù)和小波函數(shù)的光滑程度越來越高,頻率特性越來越好,即小波重構信號的光滑程度和頻率特性提高。但是濾波器長度增加,其支撐長度增加,導致邊緣信息和奇異點定位不準確,同時小波變換的計算量增大。在實際應用中綜合沖擊信號特征以及時域分辨率等要求,在濾除沖擊信號的分析中選取DB4小波對信號進行小波分析。DB4小波分解重構對應的濾波系數(shù)見表1。

圖4幾種DB小波尺度函數(shù)圖

Fig。4Some DB scaling function

圖5幾種DB小波函數(shù)圖

Fig。5 Some DB wavelet function

表1DB4小波濾波系數(shù)

Tab。1DB4 wavelet filter coefficients

k

分解低

通濾波系數(shù)

分解高

通濾波系數(shù)

重構低

通濾波系數(shù)

重構高

通濾波系數(shù)

0

-0。010 6

-0。230 4

0。230 4

-0。010 6

1

0。032 9

0。714 8

0。714 8

-0。032 9

2

0。030 8

-0。630 9

0。630 9

0。030 8

3

-0。187 0

-0。028 0

-0。028 0

0。187 0

4

-0。028 0

0。187 0

-0。187 0

-0。028 0

5

0。630 9

0。030 8

0。030 8

-0。630 9

6

0。714 8

-0。032 9

0。032 9

0。714 8

7

0。230 4

-0。010 6

-0。010 6

-0。230 4

2。2模擬沖擊信號時域分析

2。2。1濾除理想沖擊信號

選取的模擬信號為理想沖擊信號和不同頻率正弦信號的疊加f1(t),并在f1(t)信號中加入信號采集系統(tǒng)產(chǎn)生的白噪聲n(t),且噪聲信號的σ值為0。342 4,f1(t)的表達式為:

f1(t)=sin2πωt+sin4πωt+sin6πωt+

A+n(t)。(5)

式中:ω為角頻率。理想沖擊信號幅值趨于無窮即A→∞,模擬沖擊信號經(jīng)過經(jīng)驗限幅后其幅值為遠小于A。圖6(a)為限幅后的f1(t)信號,其中局部放大部分顯示f1(t)信號加噪聲后效果。f1(t)信號經(jīng)過分解后,其1~6層的沖擊和噪聲信號的模極大值序列中的最大值分別為9。4,4。6,2。3,1。4,0。73,0。63。其中第4與第5層的模極大值差值為0。67,小于k×σ值0。684 8。因此自適應分解的層數(shù)確定為5層。圖6(b)為經(jīng)過小波分析處理得到不同頻率正弦的疊加信號。

由圖6可以看出,采用本文的自適應小波分析方法可以有效濾除模擬信號中的沖擊信號,得到所需微弱振動信息。

圖6分析模擬理想沖擊信號

Fig。6Analyze analog ideal shock signal

2。2。2濾除非理想沖擊信號

選取模擬沖擊信號x(t)=sintt,微弱振動信號為不同頻率的正弦疊加信號。記非理想沖擊信號為 f2(t),并在f2(t)信號中加入信號采集系統(tǒng)產(chǎn)生的白噪聲n(t),且噪聲信號的σ值為0。342 4,f2(t)的表達式為:

f2(t)=sin2πωt+sin4πωt+sin6πωt+Ax(t)+n(t)。 (6)

式中:ω為角頻率;A為振幅。圖7(a)為f2(t)信號,其中局部放大部分顯示f2(t)信號加噪聲后效果。f2(t)信號經(jīng)過分解后,其1-5層的沖擊和噪聲信號的模極大值序列中的最大值分別為20。0,7。5,6。6,6。1,5。6。其中第3與第4層的模極大值差值為0。5小于k×σ值0。684 8,因此自適應分解的層數(shù)確定為4層。圖7(b)為經(jīng)過小波濾除沖擊信號后得到的微弱振動信號。

由圖7可以看出,對于非理想振動信號,采用本文自適應小波分析方法可以有效濾除沖擊,得到淹沒在沖擊信號中幅值較小且不同頻率的正弦疊加信號。

圖7分析模擬非理想沖擊信號

Fig。7Analyze analog nonideal shock signal

3小波濾除實測沖擊信號

3。1實驗系統(tǒng)構成

本實驗系統(tǒng)主要由硬件和軟件2部分組成,其中硬件部分包括:壓電式加速度傳感器、高速信號采集卡、電荷放大器、振動信號導向管、力錘和金剛砂。軟件部分包括信號采集和處理程序。由力錘和金剛砂產(chǎn)生的振動信號通過安裝在振動信號導向管上的壓電式加速度傳感器送至電荷放大器,將電荷量轉換成模擬電壓量,然后通過多通道高速采集卡將模擬信號轉換為數(shù)字信號傳輸給計算機,最后通過計算機采集與處理程序?qū)崿F(xiàn)信號的處理與顯示,完成對沖擊信號的濾除。實驗系統(tǒng)如圖8所示。

實驗系統(tǒng)選用的壓電式加速度傳感器為PCB357B03電荷型傳感器,靈敏度為10 pc/g。采用的高速采集卡主要參數(shù)為12位分辨率,連續(xù)實時傳輸速率32 MB/s,傳輸接口為USB2。0高速屏蔽線。軟件設置采樣率為80 000 Hz,采樣長度8 192個點,系統(tǒng)噪聲信號的σ值為0。001 25。

圖8實驗系統(tǒng)

Fig。8Experimental system

3。2實驗結果分析

實驗信號由2部分組成,模擬油井生產(chǎn)出砂監(jiān)測過程中對出砂信號的提取\[12\]。第1部分為力錘敲擊振動信號導向管后產(chǎn)生的沖擊信號,以模擬旋轉機械、螺桿泵等對生產(chǎn)管道造成的沖擊信號;第2部分為在導向管研磨金剛砂產(chǎn)生的微弱振動信號,以模擬砂礫沖擊輸油管壁產(chǎn)生的振動信號。

3。2。1對金剛砂振動信號分析

在導向管研磨金剛砂后產(chǎn)生微弱振動信號,壓電式加速度傳感器接收到受系統(tǒng)噪聲干擾的微弱振動信號,信號通過放大采集后,經(jīng)過信號采集小波分析處理程序后,得到所需振動信號,如圖9所示。

圖9金剛砂振動信號

Fig。9Emery vibration signal

圖9(a)為實驗系統(tǒng)采集的原始振動信號,振動信號經(jīng)過分解后,其1~3層的沖擊和噪聲信號的模極大值序列中的最大值分別為0。006 2,0。005 6,0。005 2。其中第1與第2層的模極大值差值為0。000 6小于k×σ值0。002 5,因此,自適應分解的層數(shù)確定為2層。圖9(b)為經(jīng)過小波濾波處理后振動信號,橫坐標為振動時間,縱坐標為振動加速度。通過時域波形對比,經(jīng)過小波濾波后信號幅值基本不變,信號基本保持原有信號特征。

圖10為對金剛砂振動信號進行頻域分析的結果,通過對比發(fā)現(xiàn)頻譜特征基本相符。因此利用本文自適應小波分析方法濾除沖擊信號時不會對混雜在沖擊信號中的重要振動信號造成畸變和干擾。

圖10金剛砂振動信號FFT

Fig。10Emery vibration signal FFT

3。2。2對混有沖擊信號的金剛砂振動信號分析

在導向管上研磨金剛砂的同時,使用力錘敲擊導向管產(chǎn)生振動信號,經(jīng)過數(shù)據(jù)采集處理后如圖11所示。

圖11混合振動信號

Fig。11Mixed vibration signal

圖11(a)為本實驗系統(tǒng)采集的原始振動信號時域波形。振動信號經(jīng)過分解后,其1~6層的沖擊和噪聲信號的模極大值序列中的最大值分別為0。073 2,0。052 5,0。027 6,0。018 7,0。010 2,0。008 3。其中第5與第6層的模極大值差值為0。001 9小于k×σ值0。002 5,因此自適應分解的層數(shù)確定為6層。圖11(b)為經(jīng)過濾波處理后振動信號時域波形,橫坐標為振動時間,縱坐標為振動加速度值。沖擊信號加速度最高幅值可達1 m/s2,經(jīng)過濾波處理后提取的微弱振動信號加速度幅值為0。005 m/s2,其幅值與圖9(a)中只存在單一金剛砂振動信號加速度幅值一致。從時間角度分析沖擊信號持續(xù)時間僅為0。024 s,而經(jīng)過濾波分析后信號時間長度為一個完整采樣周期時間0。01 s。即混雜在沖擊信號中的微弱振動信號被提取出來,結果表明一個振動周期內(nèi)都存在研磨金剛砂產(chǎn)生的振動信號,其中在0。036 s和0。096 s附近金剛砂振動信號幅值最高。

圖12為混合信號頻譜分析,從頻域角度分析可知,原始混合信號頻域范圍大概在0~20 000 Hz,經(jīng)過濾波處理后頻域范圍集中在12 000~18 000 Hz附近,由文獻\[13\]可知此頻段為金剛砂振動頻段,因此淹沒在沖擊信號中的微弱金剛砂振動信號通過本文方法被成功提取出來。

圖12混合信號頻譜分析

Fig。12Mixed vibration signal FFT

4結論

本文采用自適應層數(shù)分解小波分析方法,針對微弱振動淹沒在沖擊信號中的問題,根據(jù)信號特征合理選擇小波基,對模擬理想和非理想沖擊信號進行自適應層數(shù)DB4小波分析,結果表明該方法可以有效濾除沖擊信號并提取微弱振動信息。同時,對本實驗系統(tǒng)采集的混合振動信號進行自適應DB4小波分析,成功濾除沖擊信號。通過模擬計算和實驗結果的時頻分析可知,自適應小波濾除沖擊信號算法正確、結果可靠,自適應層數(shù)分解小波分析方法可以濾除沖擊提取微弱振動信號,解決工程實際問題。

參考文獻

[1]霍靜思,何遠明,肖莉平,等。 高溫后鋼管混凝土抗多次沖擊力學性能試驗研究\[J\]。 湖南大學學報:自然科學版, 2012, 39(9): 6-10。

HUO Jingsi, HE Yuanming, XIAO Liping, et al。 Experimental study on the dynamic behavior of concretefilled steel tube after exposure to high temperatures under multiple impact loadings\[J\]。 Journal of Hunan University: Natural Sciences, 2012, 39(9): 6-10。(In Chinese)

\[2\]鄭雪蓮,李顯生,任園園,等。 非滿載罐體液體沖擊等效機械模型參數(shù)確定\[J\]。 湖南大學學報: 自然科學版, 2013,40(6): 53-58。

ZHENG Xuelian, LI Xiansheng, REN Yuanyuan, et al。 Parameter values of equivalent mechanical model for liquid sloshing in partiallyfilled tanks \[J\]。 Journal of Hunan University: Natural Sciences, 2013, 40(6): 53-58。 (In Chinese)

\[3\]LANGLEY R S。 The analysis of impact forces in randomly vibrating elastic systems\[J\]。 Journal of Sound and Vibration, 2012, 331(16): 3738-3750。

\[4\]VIJAYAN K, WOODHOUSE J。Shock amplification, curve veering and the role of damping\[J\]。 Journal of Sound and Vibration, 2013, 333(5):1379-1389。

\[5\]LI Hui,ZHANG Yuping,ZHENG Haiqi。 Application of hermitian wavelet to crack fault detection in gearbox\[J\]。 Mechanical Systems and Signal Processing, 2013, 25(4): 1353-1363。

\[6\]SHARK L K, YU C。Denoising by optimal fuzzy thresholding in wavelet domain\[J\]。 Electronics Letters, 2000, 36(6): 581-582。

\[7\]WANG Zuocai,REN Weixin, LIU Jingliang。A synchrosqueezed wavelet transform enhanced by extended analytical mode decomposition method for dynamic signal reconstruction\[J\]。 Journal of Sound and Vibration, 2013, 332(22): 6016-6028。

\[8\]譚繼勇,陳雪峰,何正嘉。 沖擊信號的隨機共振自適應檢測方法\[J\]。 機械工程學報, 2010, 46(23): 61—67。

TAN Jiyong, CHEN Xuefeng, HE Zhengjia。 Detection method with adaptive stochastic resonance\[J\]。 Journal of Mechanical Engineering, 2010,46(23): 61-67。(In Chinese)

\[9\]騰軍, 朱焰煌, 周峰, 等。 自適應分解層數(shù)的小波域中值濾波振動信號降噪法\[J\]。 振動與沖擊,2009,28(12): 58-62。

TENG Jun, ZHU Yanhuang,ZHOU Feng,et al。 Vibration signal denoising method based on median filter in wavelet domain with selfadaptive level decomposition\[J\]。 Journal of Vibration and Shock, 2009, 28(12): 58-62。 (In Chinese)

\[10\]LEE Y S, LEE S K, SHIN K。Piezoceramic cantilever sensor design for weakimpact detection on plates\[J\]。 Shock and Vibration, 2012, 19(5): 879-890。

\[11\]李媛。小波變換及其工程應用\[M\]。北京: 北京郵電大學出版社, 2010:18-24。

LI Yuan。Wavelet transform and engineering application\[M\]。 Beijing: Beijing University of Posts and Telecommunications Publishing House, 2010:18-24。(In Chinese)

\[12\]劉剛,劉彭濤,韓金良,等。油井出砂監(jiān)測技術現(xiàn)狀及發(fā)展趨勢\[J\]??萍紝螅?013,31(25):75-79。

LIU Gang, LIU Pengtao, HAN Jinliang, et al。 Status and development trends of sand monitoring technology of oil wells\[J\]。 Science & Technology Review,2013,31 (25): 75-79。(In Chinese)

\[13\]秦敏, 劉志國, 劉剛, 等?;诩铀俣葌鞲衅鞯母哳l振動信號檢測分析\[J\]。測控技術, 2013,32(5): 1-4。

QING Min, LIU Zhiguo, LIU Gang, et al。 Research on high frequency vibration signal detection technology based on accelerometer\[J\]。Journal of Measurement and Control Technology,2013, 32(5): 1-4。(In Chinese)

猜你喜歡
小波變換沖擊振動
振動的思考
科學大眾(2023年17期)2023-10-26 07:39:14
振動與頻率
天天愛科學(2020年6期)2020-09-10 07:22:44
中立型Emden-Fowler微分方程的振動性
MATLAB在《數(shù)字圖像處理》課程中的輔助教學
基于互信息和小波變換的圖像配準的研究
保持細節(jié)的Retinex紅外圖像增強算法
基于小波包絡分析的滾動軸承典型故障診斷技術研究
商情(2016年11期)2016-04-15 07:39:54
奧迪Q5換擋沖擊
奧迪A8L換擋沖擊
一汽奔騰CA7165AT4尊貴型車換擋沖擊
武安市| 宜都市| 甘南县| 商水县| 黄浦区| 靖远县| 察雅县| 新野县| 伊宁市| 额敏县| 洮南市| 长沙县| 丽江市| 和林格尔县| 两当县| 乳山市| 胶州市| 纳雍县| 康定县| 永福县| 山东省| 泰顺县| 阳城县| 峨边| 昭通市| 双峰县| 东辽县| 三江| 惠东县| 威宁| 武宁县| 视频| 开原市| 原阳县| 高碑店市| 陆丰市| 富川| 晋州市| 辛集市| 麻栗坡县| 台东市|