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

?

不同云微物理方案對(duì)青藏高原一次強(qiáng)降水的模擬影響分析

2022-04-22 02:43:06毛智朱志鵬張如翼周立旻
熱帶氣象學(xué)報(bào) 2022年1期
關(guān)鍵詞:云水水汽降水

毛智,朱志鵬,張如翼,周立旻

(華東師范大學(xué)地理科學(xué)學(xué)院/地理信息科學(xué)教育部重點(diǎn)實(shí)驗(yàn)室,上海 200241)

1 引 言

青藏高原,位于亞洲內(nèi)陸,平均海拔為4 000~5 000 m,是中國最大、世界海拔最高的高原,被稱為“世界屋脊”。其特殊的熱力和動(dòng)力作用不僅對(duì)高原本身及其周邊地區(qū)的氣候有重要影響,甚至對(duì)北半球乃至全球的大氣環(huán)流產(chǎn)生影響[1]。在全球變暖的背景下,高原地區(qū)的許多氣象要素也在發(fā)生著顯著變化。研究都表明高原地區(qū)是氣候變化的敏感區(qū),那么對(duì)高原地區(qū)氣候的研究就顯得尤為重要。Duan等[2]分析了1980—2013年高原近地表氣溫觀測(cè)數(shù)據(jù),指出高原仍在顯著變暖并且云-輻射反饋機(jī)制在調(diào)節(jié)變暖趨勢(shì)中發(fā)揮重要作用。Liu 等[3]基于 1980—2018年高原地區(qū) 86個(gè)站點(diǎn)近地表氣溫和地表溫度觀測(cè)數(shù)據(jù),研究高原感熱通量的變化以及對(duì)南亞大氣環(huán)流的可能影響,結(jié)果顯示感熱通量在2001—2018年間逐漸增加并且這一變化使得印度副熱帶高壓南移收縮。

基于模式的氣候模擬結(jié)果表明地形對(duì)氣候模擬結(jié)果具有顯著的影響。Zhu 等[4]將CMIP5 中30個(gè)氣候模式和CMIP6 中12個(gè)氣候模式對(duì)1961—2005年中國區(qū)域10個(gè)氣候指標(biāo)的模擬結(jié)果與地表觀測(cè)值進(jìn)行比較,結(jié)果顯示CMIP6 多模式集合對(duì)整個(gè)中國區(qū)域年平均氣溫等氣候指標(biāo)空間分布的模擬效果很好,但對(duì)高原地區(qū)的模擬卻存在冷偏差,總體來說CMIP6 的表現(xiàn)優(yōu)于CMIP5。Li等[5]將中國氣象科學(xué)研究院氣候系統(tǒng)模式(CAMSCSM)對(duì)2000—2015年東亞夏季逐時(shí)降水的模擬結(jié)果與地表觀測(cè)值比較,結(jié)果表明模擬偏差與地形密切相關(guān),具體表現(xiàn)為在我國西部高原(東部平原)地區(qū),降水量和降水頻次被高估(低估),降水持續(xù)時(shí)間變長(zhǎng)(縮短)。

區(qū)域尺度天氣過程的模擬表明,云微物理特征對(duì)降水模擬結(jié)果具有顯著影響。Yun 等[6]利用WRF 模式模擬2008—2017年中國東部地區(qū)暖季降水,以了解模式模擬能力并探究對(duì)模式缺點(diǎn)的改進(jìn)方法,結(jié)果表明模式能夠很好地再現(xiàn)季節(jié)與季節(jié)內(nèi)降水的空間分布以及與夏季風(fēng)相關(guān)的降水帶南北遷移等特征,并且改善對(duì)次網(wǎng)格云量和氣溶膠效應(yīng)的處理方式將有助于抑制模式高估降水量這一缺點(diǎn)。Orr 等[7]利用WRF 模式并選擇不同云微物理方案,研究云微物理方案對(duì)喜馬拉雅朗塘山谷夏季季風(fēng)性降水的敏感性,對(duì)所有方案水凝物的分析結(jié)果表明冷云過程是主要的降水形成機(jī)制。Gao 等[8]利用WRF 模式并結(jié)合中國氣象科學(xué)研究院微物理方案(WRF-CAMS)對(duì)青藏高原一次對(duì)流性降水過程進(jìn)行了模擬,旨在探討高原降水的微物理特性,結(jié)果表明該方案能夠合理地再現(xiàn)累積降水量的空間分布但會(huì)延長(zhǎng)降水時(shí)間1~3 h,并且應(yīng)該使用可變雨滴粒徑形狀參數(shù)來處理因降水速率改變而引起的雨滴粒徑改變。

為了解區(qū)域尺度上高原區(qū)域天氣過程模擬中不同云微物理參數(shù)化方案對(duì)高原強(qiáng)降水過程的模擬影響及可能原因,本研究利用WRF 模式中三種不同云微物理參數(shù)化方案對(duì)高原一次強(qiáng)降水天氣過程進(jìn)行模擬分析。通過與實(shí)測(cè)的對(duì)比和對(duì)模擬結(jié)果的分析,進(jìn)一步探索對(duì)高原降水過程有較好模擬效果的云微物理參數(shù)化方案,為提升高原區(qū)域的數(shù)值模擬和數(shù)值預(yù)報(bào)精度提供參考。

2 資料與模式設(shè)置

2.1 資料選取

本研究使用的實(shí)測(cè)降水?dāng)?shù)據(jù)是從中國氣象數(shù)據(jù)網(wǎng)獲取的中國自動(dòng)氣象站與CMORPH 降水產(chǎn)品融合的逐時(shí)降水量網(wǎng)格數(shù)據(jù)(經(jīng)度范圍為70.05~139.95 °E,緯度范圍為15.05~58.95 °N,水平分辨率為0.1 °×0.1 °)。前人的研究表明這一數(shù)據(jù)能夠很好地反映中國區(qū)域降水的主要特點(diǎn)[9-10],將這一數(shù)據(jù)作為評(píng)估模式模擬降水的真值。

本研究使用的初始場(chǎng)及邊界場(chǎng)資料來自NCEP 的FNL 資料,其時(shí)間間隔為6 h,水平空間范圍為全球,水平分辨率為1 °×1 °,垂直方向分為32層,最高層為1 hPa。

2.2 參數(shù)化方案選擇

邊界層參數(shù)化方案:不同于許魯君等[11]對(duì)下墊面類型相對(duì)均一的高原那曲地區(qū)的研究,此次研究的降水過程主要發(fā)生在雅魯藏布江中游地區(qū),而這一地區(qū)地形復(fù)雜多樣,高差起伏大[12]。鑒于前人[13-14]利用WRF模式中不同邊界層參數(shù)化方案對(duì)復(fù)雜地形模擬的研究表明非局地邊界層方案較局地邊界層方案能更好地模擬復(fù)雜地形的邊界層特征,本研究選擇非局地邊界層方案中的代表性方案(Yonsei University(YSU)方案)進(jìn)行模擬。

積云對(duì)流參數(shù)化方案:由于積云對(duì)流這一物理過程的復(fù)雜性,利用區(qū)域氣候模式精確模擬降水仍然存在巨大的難度。鑒于Huang 等[15]利用WRF 模式在低水平分辨率(30 km)下對(duì)中國區(qū)域降水的研究結(jié)果和Wang 等[16]利用WRF 模式在高水平分辨率(1~5 km)下對(duì)高原地區(qū)夏季降水的研究結(jié)果,本研究選擇Kain-Fritsch(KF)方案進(jìn)行模擬。

云微物理參數(shù)化方案:云微物理過程在影響水汽凝結(jié)潛熱加熱和雨水蒸發(fā)潛熱冷卻方面具有重要作用[17]。Lin 方案和WSM6 方案都是單參數(shù)體積方案[18-19],主要通過假設(shè)粒徑分布來預(yù)測(cè)水汽和水凝物(冰、雪、霰、云水和雨水)的混合比。Lin方案是基于文獻(xiàn)[20-22]的研究,并且采用了Tao等[23]對(duì)飽和調(diào)整和冰沉降的修改。在該方案中,當(dāng)溫度低于結(jié)冰點(diǎn)(-40 ℃)時(shí),云水和雨水分別處理為冰和雪;當(dāng)溫度高于0 ℃時(shí),冰全部處理為云水;當(dāng)溫度介于二者之間時(shí),冰、雪、霰、云水和雨水五種水凝物共存[24]。WSM6 方案則分別采用Dudhia[25]和Hong等[26]的研究結(jié)果來處理冰和云水的飽和調(diào)整過程,并且采用小步長(zhǎng)時(shí)間積分來處理凍結(jié)和融化過程以提高垂直加熱廓線的精度,通過雪和霰的混合比例加權(quán)來指定它們的單粒子下落速度。Eta(Ferrier)方案預(yù)測(cè)以云水、雨水和降水冰的形式出現(xiàn)的水汽和水凝物的變化,根據(jù)降水冰的混合比估算其密度,并以此來確定降水冰的形式為雪、霰或雨夾雪,利用Eta 方案來表示云中的微物理過程與診斷混合相態(tài)過程可使得在較大的時(shí)間步長(zhǎng)下有穩(wěn)定的結(jié)果[24]。

2.3 試驗(yàn)設(shè)計(jì)

本研究使用的是WRF 模式4.0 版本。模擬時(shí)間設(shè)置為2016年6月15日08時(shí)—16日08時(shí)(北京時(shí)間,下同),其中前6 h為spin-up時(shí)間。模擬中心經(jīng)緯度為90 °E,29.5 °N,模擬區(qū)域采用三重嵌套,從外層至內(nèi)層水平分辨率分別為27 km×27 km、9 km×9 km 和3 km×3 km,格點(diǎn)數(shù)分別為147×77、367×157 和 226×100,積分時(shí)間步長(zhǎng)分別為 90 s、30 s 和 10 s,垂直方向采用 eta 坐標(biāo),共有 33 層,最高層為50 hPa,每隔1 h輸出一次結(jié)果。

模式選用RRTM 長(zhǎng)波輻射方案,Dudhia 短波輻射方案,Noah 陸面過程方案,YSU 邊界層參數(shù)化方案,KF 積云對(duì)流參數(shù)化方案,Lin、Eta 和WSM6 云微物理參數(shù)化方案。外層只使用積云對(duì)流參數(shù)化方案,中層既使用積云對(duì)流參數(shù)化方案也使用云微物理參數(shù)化方案,內(nèi)層只使用云微物理參數(shù)化方案。

3 模擬結(jié)果及分析

3.1 降水量與雷達(dá)反射率

為研究不同云微物理方案模擬的總體差異,本文首先考察了累計(jì)降水量與雷達(dá)反射率。從累計(jì)降水量場(chǎng)(圖1)可看出,三種方案模擬的降水范圍均較實(shí)測(cè)偏大。整個(gè)降水過程表現(xiàn)出兩個(gè)主要降水區(qū)域,區(qū)域A 中最大累計(jì)降水量為28.1 mm,區(qū)域B 中最大累計(jì)降水量為37.7 mm。對(duì)于區(qū)域A,三種方案模擬的位置均較實(shí)測(cè)偏東,Lin 和WSM6 方案模擬的個(gè)別點(diǎn)最大累計(jì)降水量均較實(shí)測(cè)偏大(超過了40 mm),Eta 方案模擬的總體最大累計(jì)降水量與實(shí)測(cè)相當(dāng);對(duì)于區(qū)域B,三種方案模擬的位置與實(shí)測(cè)相當(dāng),Lin 方案模擬的個(gè)別點(diǎn)最大累計(jì)降水量較實(shí)測(cè)偏大(超過了40 mm),Eta 和WSM6 方案模擬的總體最大累計(jì)降水量與實(shí)測(cè)相當(dāng)。除區(qū)域A 和B 之外,三種方案均模擬出了一些虛假的大降水區(qū)域。

圖1 6月15日14時(shí)—16日08時(shí)16 h累計(jì)降水分布

三種方案模擬得到的15日19—20 時(shí)的雷達(dá)反射率(圖2a1~2c2)在強(qiáng)度和覆蓋范圍上基本一致,部分區(qū)域雷達(dá)反射率大于35 dBZ;Eta 方案模擬得到的16日00—01 時(shí)的雷達(dá)反射率覆蓋范圍比 Lin 和 WSM6 方案稍大(圖 2a3~2c4),部分區(qū)域雷達(dá)反射率大于40 dBZ??傮w而言,強(qiáng)雷達(dá)反射率區(qū)域與前文提到的主要降水區(qū)域有較好的對(duì)應(yīng),但形態(tài)各異。

圖2 Lin方案(a)、Eta方案(b)、WSM6方案(c)模擬的雷達(dá)反射率

3.2 水汽輸送情況

水汽通量與水汽通量散度可定量描述水汽輸送的大小和方向以及水汽在何處聚集[27]。三種方案模擬的水汽通量及其散度(圖3)大小和范圍分布差異不大,最大水汽通量在5 kg/(m·s·hPa)左右,最大水汽通量散度在3×10-6kg/(m2·s·hPa)(輻合)左右,均表明在強(qiáng)降水過程發(fā)生之前南北水汽的經(jīng)向輸送和聚集。水汽的輸送為強(qiáng)降水提供了有利的水汽條件。

圖3 Lin方案模擬(a)的6月15日14時(shí)500 hPa水汽通量(矢量,單位:kg/(m·s·hPa))和水汽通量散度合成(陰影,單位:10-6 kg/(m2·s·hPa));與對(duì)應(yīng)時(shí)刻Eta方案模擬(b)、WSM6方案(c)模擬的水汽通量和水汽通量散度

3.3 動(dòng)力熱力差異

為探究三種方案模擬的強(qiáng)雷達(dá)反射率區(qū)域差異的內(nèi)在原因,下面將從溫度場(chǎng)和垂直運(yùn)動(dòng)場(chǎng)兩方面進(jìn)行分析。相當(dāng)位溫能同時(shí)反映大氣的溫度和濕度狀況,其垂直梯度大,則有利于低層暖濕空氣向上抬升,與高層干冷空氣相遇形成相當(dāng)位溫線密集帶[28]。雖然三種方案模擬的相當(dāng)位溫線密集帶位置差異較大,但無論是在降水過程的前一階段(15日 20 時(shí))還是后一階段(16日 00 時(shí)),三種方案模擬的相當(dāng)位溫線密集帶位置與各自的強(qiáng)雷達(dá)反射率區(qū)域均有較好的對(duì)應(yīng)。同時(shí),絕對(duì)渦度大值區(qū)(>50×10-5s-1)配合垂直風(fēng)速大值區(qū)(>3 m/s)形成的垂直運(yùn)動(dòng)大值區(qū)和相當(dāng)位溫線密集帶位置也有較好的對(duì)應(yīng),亦即與強(qiáng)雷達(dá)反射率區(qū)域有較好的對(duì)應(yīng)。

圖4 Lin方案模擬(a、d),Eta方案(b、e),WSM6方案(c、f)的相當(dāng)位溫(等值線,單位:K)、絕對(duì)渦度(陰影,單位:10-5 s-1)和垂直風(fēng)速(矢量,單位:m/s)沿圖2所示線段的剖面

3.4 微物理差異

3.4.1 水凝物垂直分布

前文以宏觀特征為落腳點(diǎn)對(duì)三種方案模擬差異進(jìn)行了分析,為了從微觀過程進(jìn)行解釋,圖5 和圖6 給出了三種方案模擬的水凝物垂直含量隨時(shí)間的演變特征。Lin 方案模擬的云水含量最多(最大值約為0.44×10-4kg/kg)且絕大部分是分布在0 ℃層之上,其云系發(fā)展最高,WSM6 方案模擬的云水含量次之,Eta 方案模擬的云水含量最少;三種方案模擬的雨水均分布在0 ℃層之下(圖5a2~5c2),且含量差異不大;固態(tài)水凝物冰(圖5a3~5c3)、雪(圖 6a、6b)和霰(圖 6c、6d)主要分布在 0 ℃層之上。三種方案模擬的冰的高值中心雖在同一高度,但Eta方案和WSM6方案模擬的冰含量均高于Lin方案;Lin方案模擬的雪含量低于WSM6方案,但高值中心位置卻高于WSM6 方案;Lin 方案和WSM6方案模擬的霰含量和分布差異不大。

圖5 15日14時(shí)—16日08時(shí)Lin方案(a1~a3)、Eta方案(b1~b3)、WSM6方案(c1~c3)模擬的區(qū)域平均水凝物混合比(單位:10-4 kg/kg)的氣壓-時(shí)間剖面

圖6 同圖5,但為L(zhǎng)in方案(a、c)模擬、WSM6方案(b、d)模擬

進(jìn)一步分析水凝物在整個(gè)降水過程的含量差異情況,圖7給出了三種方案模擬的水凝物的垂直廓線。對(duì)于冰廓線,Eta 方案的模擬值遠(yuǎn)遠(yuǎn)高于其他兩種方案,最大值約為2.8×10-4kg /kg;對(duì)于雪廓線,WSM 方案模擬值高于Lin方案,最大值約為1.1×10-4kg/kg,但模擬高度較Lin 方案稍低;對(duì)于霰廓線,Lin方案和WSM6方案模擬的最大值約為0.5×10-4kg/kg,且高度分布差異不大。

圖7 15日14時(shí)—16日08時(shí)平均的Lin方案(a)、Eta方案(b)、WSM6方案(c)模擬的水凝物混合比(單位:10-4 kg/kg)隨氣壓的變化

3.4.2 云微物理轉(zhuǎn)化項(xiàng)

為深入探討固態(tài)水凝物含量差異的原因,并且考慮到Lin和WSM6方案微物理轉(zhuǎn)化過程較Eta方案復(fù)雜得多,圖8給出了這兩種方案云微物理過程的轉(zhuǎn)化特征,過程詳細(xì)意義見表1、表2。

表1 Lin方案云微物理轉(zhuǎn)化過程

表2 WSM6方案云微物理轉(zhuǎn)化過程

圖8 15日14時(shí)—16日08時(shí)Lin方案(a)、WSM6方案(b)模擬的水凝物(單位:10-11 kg/kg)的轉(zhuǎn)化過程

首先就水汽收支而言,兩種方案表現(xiàn)出不同的過程:Lin 方案中,水汽主要存在水汽向霰的凝華Pgdep(霰向水汽的升華Pgsub)、水汽向雪的凝華Psdep(霰向水汽的升華Pssub)以及雨水向水汽的蒸發(fā)Prevp等過程,這表明水汽主要向固態(tài)水凝物轉(zhuǎn)化;而在WSM6方案中,則主要是水汽與云水的相互轉(zhuǎn)化Pcond、水汽與雨水的相互轉(zhuǎn)化Prevp、水汽與冰的相互轉(zhuǎn)化Pigen 和Pidep、水汽與雪的相互轉(zhuǎn)化Psevp 和Psdep、水汽與霰的相互轉(zhuǎn)化Pgevp 和Pgdep 等過程,這表明水汽主要向液態(tài)水凝物轉(zhuǎn)化。

就云水而言,Lin 方案中主要表現(xiàn)為支出,即被霰碰并Pgacw、被雪碰并Psacw、被雨水碰并Pracw、同質(zhì)凍結(jié)成冰 Pihom 等匯項(xiàng);WSM6 方案中存在前文提到的主要源匯項(xiàng)Pcond、被雪碰并轉(zhuǎn)化為霰或者雨水Psacw、被霰碰并轉(zhuǎn)化為霰或者雨水Pgacw、被雨水碰并Pracw。就冰而言,Lin 方案中主要存在一個(gè)大的匯項(xiàng),即通過貝吉龍過程形成雪Psfi,此過程比其它過程大一個(gè)量級(jí);WSM6方案中除了前文提到的源匯項(xiàng)Pigen 和Pidep、還存在向雪的轉(zhuǎn)化Psaci 和Psaut 以及向霰的轉(zhuǎn)化Pgaci。就雪而言,Lin 和WSM6 方案中都存在來自冰的轉(zhuǎn)化項(xiàng)以及向霰的匯項(xiàng),同時(shí)WSM6 方案還存在另外一個(gè)大的匯項(xiàng),即冰融化生成雨水Psmlt。

對(duì)雨水收支過程的分析表明,Lin和WSM6方案都表現(xiàn)出霰向雨水轉(zhuǎn)化的特點(diǎn):Lin 方案中,水汽、云水和雪均大量向霰轉(zhuǎn)化,霰轉(zhuǎn)化為雨水;WSM6 方案中,云水和雪均大量向霰轉(zhuǎn)化,霰轉(zhuǎn)化為雨水。

綜上所述,對(duì)于此次降水,Lin 方案模擬的微物理過程表現(xiàn)為:通過水汽向霰粒子凝華的過程Pgdep,通過霰碰并水汽凝華生成的雪粒子而造成霰粒子增長(zhǎng)的過程Pgacs,以及霰碰并云水造成霰粒子增長(zhǎng)的過程Pgacw,此三種過程生成的霰粒子最終融化為雨水(Pgmlt)。WSM6 方案模擬的微物理過程表現(xiàn)為:一方面水汽凝結(jié)成云水,云水被雪和霰碰并收集轉(zhuǎn)化為霰,之后霰融化為雨水;另一方面水汽凝華為冰,一部分冰轉(zhuǎn)化為雪,雪直接融化為雨水或轉(zhuǎn)化為霰融化為雨水,另一部分冰轉(zhuǎn)化為霰,霰融化為雨水。

4 結(jié) 論

本文利用 WRF 模式對(duì) 2016年 6月 15—16日發(fā)生在青藏高原雅魯藏布江中游地區(qū)一次強(qiáng)降水天氣過程進(jìn)行了模擬研究,探討了三種云微物理參數(shù)化方案對(duì)高原降水的影響,分析了三種方案水凝物的分布并著重分析了Lin 和WSM6 這兩種復(fù)雜云微物理方案的微物理轉(zhuǎn)化過程,得到以下結(jié)論。

(1) 整個(gè)降水過程表現(xiàn)出兩個(gè)降水區(qū)域。對(duì)于區(qū)域A,三種方案模擬的位置均偏東,Lin 和WSM6 方案模擬的降水量偏多,Eta 方案模擬的降水量與實(shí)測(cè)相當(dāng);對(duì)于區(qū)域B,三種方案模擬的位置與實(shí)測(cè)相當(dāng),Lin 方案模擬的降水量偏多,Eta 和WSM6方案模擬的降水量與實(shí)測(cè)相當(dāng)。

(2) 對(duì)比分析三種方案模擬的液態(tài)水凝物(云水、雨水)和固態(tài)水凝物(冰、雪、霰)垂直分布和垂直廓線可知:三種方案模擬的冰粒子差異較明顯;Lin和WSM6方案模擬的雪粒子差異較大,但霰粒子無明顯差異。

(3) 對(duì)比分析了Lin 和WSM6 方案的微物理轉(zhuǎn)化過程,結(jié)果表明:這兩種方案都表現(xiàn)出了霰向雨水轉(zhuǎn)化的特點(diǎn)。在Lin方案中,通過水汽向霰粒子凝華的過程,通過霰碰并水汽凝華生成的雪粒子而造成霰粒子增長(zhǎng)的過程,以及霰碰并云水造成霰粒子增長(zhǎng)的過程,這三種過程生成的霰粒子最終融化為雨水。而在WSM6 方案中,一方面水汽凝結(jié)成云水,云水被雪和霰碰并收集轉(zhuǎn)化為霰,之后霰融化為雨水;另一方面水汽凝華為冰,一部分冰轉(zhuǎn)化為雪,雪直接融化為雨水或轉(zhuǎn)化為霰融化為雨水,另一部分冰轉(zhuǎn)化為霰,霰融化為雨水。

猜你喜歡
云水水汽降水
云水禪心
青藏高原上空平流層水汽的時(shí)空演變特征
黑龍江省玉米生長(zhǎng)季自然降水與有效降水對(duì)比分析
黑龍江氣象(2021年2期)2021-11-05 07:07:00
云水謠
幸福家庭(2019年14期)2019-01-06 09:14:52
為什么南極降水很少卻有很厚的冰層?
家教世界(2018年16期)2018-06-20 02:22:00
悠然云水(七絕)
寶藏(2017年4期)2017-05-17 03:33:48
1979~2011年間平流層溫度及平流層水汽的演變趨勢(shì)
降水現(xiàn)象儀模擬軟件設(shè)計(jì)與實(shí)現(xiàn)
深圳“5·11”特大暴雨過程的水汽輸送特征分析
Fast Scheme for Projective Geometric Correction and Edge Blending Based on High Dynamic Range Images
手机| 兴宁市| 柞水县| 武陟县| 霞浦县| 乐清市| 札达县| 寻甸| 宜兴市| 青神县| 名山县| 安平县| 浑源县| 阿城市| 五指山市| 巨鹿县| 宿迁市| 舒兰市| 江山市| 永定县| 长春市| 延安市| 巴东县| 通州市| 水富县| 龙海市| 和顺县| 博白县| 阜南县| 东源县| 南木林县| 福建省| 高陵县| 方山县| 嘉定区| 曲阜市| 梧州市| 开阳县| 永德县| 监利县| 松江区|