通量及其不確定性對農(nóng)業(yè)區(qū)高塔CO2濃度模擬的影響*

胡 誠，張 彌**，肖 薇，王詠薇，王 偉，Tim Griffis，劉壽東，李旭輝

?

通量及其不確定性對農(nóng)業(yè)區(qū)高塔CO2濃度模擬的影響*

胡 誠1,2，張 彌1,2**，肖 薇1,2，王詠薇1,2，王 偉1,2，Tim Griffis3，劉壽東1,2，李旭輝1

（1.南京信息工程大學(xué)大氣環(huán)境中心，南京 210044；2.南京信息工程大學(xué)大氣環(huán)境與裝備技術(shù)協(xié)同創(chuàng)新中心，南京210044；3.明尼蘇達(dá)大學(xué)，美國圣保羅市，55108）

利用WRF-STILT模型模擬玉米種植區(qū)生長季（6-9月）小時CO2濃度，并基于美國最大農(nóng)業(yè)種植區(qū)‘玉米帶’100m高塔CO2濃度觀測數(shù)據(jù)，對WRF-STILT模型的模擬能力及CO2通量的不確定性對模擬結(jié)果的影響進(jìn)行分析。結(jié)果表明：（1）WRF-STILT能夠模擬高塔觀測的CO2濃度日變化特征，模擬值與觀測值的均方根誤差為13.70mmol×mol-1，模擬結(jié)果偏高7.26mmol×mol-1。（2）EDGAR和Carbon Tracker兩種典型化石燃料的CO2通量，其區(qū)域平均值相差＜6%，但兩者對CO2濃度增加值的模擬結(jié)果相差約10%；（3）CO2通量空間分辨率的差異會導(dǎo)致模擬結(jié)果產(chǎn)生偏差，使用區(qū)域邊長為1o的EDGAR化石燃料CO2通量模擬的濃度貢獻(xiàn)值僅為0.1o的0.4倍，且空間分辨率越低，模擬誤差越大；（4）白天和夜晚Carbon Tracker模擬的植被生態(tài)系統(tǒng)凈交換數(shù)據(jù)是高塔渦度相關(guān)方法觀測結(jié)果的2.26和1.56倍，下墊面分類的誤差以及相應(yīng)的通量模擬誤差使模擬的CO2濃度貢獻(xiàn)出現(xiàn)12mmol×mol-1的差異，這是模擬結(jié)果偏高7.26mmol×mol-1的潛在誤差來源。研究認(rèn)為，WRF-STILT模型和高空間及時間分辨率的CO2通量能夠較好模擬出農(nóng)業(yè)區(qū)生長季的CO2強日變化特征，CO2通量的誤差是模擬結(jié)果誤差的主要來源，研究結(jié)果表明該方法具有評估和優(yōu)化通量的巨大潛力。

WRF-STILT模型；渦度相關(guān)；化石燃料；通量不確定性

對陸地生態(tài)系統(tǒng)CO2通量的準(zhǔn)確估算是預(yù)測未來氣候變化的基礎(chǔ)[1-5]。傳統(tǒng)的計算和觀測方法，如植被生物量清單調(diào)查法、渦度相關(guān)方法、陸面過程模型或IPCC算法等由于自身的局限性，對陸地生態(tài)系統(tǒng)在區(qū)域尺度（102～106km2）碳交換估算上存在較大不確定性，這主要來自空間異質(zhì)性大以及不同方法觀測和模擬能力的限制[1,6-8]?；诖髿鈧鬏斈Ｐ秃拖闰灒ǔ跏技僭O(shè)）溫室氣體通量，并結(jié)合高精度的大氣濃度實際觀測，反演區(qū)域甚至全球尺度的后驗（真實）通量，已經(jīng)被越來越多地應(yīng)用于對溫室氣體（CO2、CH4、N2O）通量的估算上[9-12]。而大氣傳輸模型對CO2濃度的模擬能力則是其反演陸地生態(tài)系統(tǒng)CO2通量結(jié)果準(zhǔn)確性的基礎(chǔ)。

基于拉格朗日原理的WRF-STILT（Stochastic Time Inverted Lagrangian Transport model）大氣傳輸模型[13]，相對于歐拉傳輸模型，主要有以下兩點優(yōu)勢：（1）在已知觀測和模擬站點的情況下，對靠近站點的區(qū)域進(jìn)行氣象驅(qū)動場的多重嵌套模擬，使模擬的濃度貢獻(xiàn)源區(qū)（足跡權(quán)重）更準(zhǔn)確；（2）通過釋放大量空氣粒子的形式來模擬大氣的湍流運動比參數(shù)化方案更接近真實空氣隨機運動的情況[14]。此外，WRF-STILT模型更以其計算效率高，數(shù)值模擬穩(wěn)定等優(yōu)勢，逐漸被用于區(qū)域尺度的植被生態(tài)系統(tǒng)凈交換（NEE）的反演。不過，由于反演過程均假設(shè)化石燃料燃燒CO2通量無誤差，而僅調(diào)整陸地生態(tài)系統(tǒng)CO2通量，所以模型中所使用的先驗化石燃料燃燒CO2通量的不確定性可能會對最終的濃度模擬和生態(tài)系統(tǒng)凈交換的反演帶來很大誤差，這包括不同化石燃料燃燒CO2通量的選取，及其與真實通量大小和空間分布的差異；除此之外生態(tài)系統(tǒng)凈交換CO2通量在空間上的分布也會影響最后的通量優(yōu)化結(jié)果，因此，CO2通量及其不確定性對WRF-STILT模型模擬結(jié)果的影響是反演生態(tài)系統(tǒng)凈交換的基礎(chǔ)。農(nóng)田生態(tài)系統(tǒng)作為陸地生態(tài)系統(tǒng)的重要組成部分，占全球陸地面積的12%[15]，且其在所有陸地生態(tài)系統(tǒng)中受人為活動干擾最強，碳存儲量短期變化最大[16-17]，但研究農(nóng)田生態(tài)系統(tǒng)CO2通量及其不確定性對區(qū)域CO2濃度模擬的影響卻少有報道[18-19]。

本文基于美國玉米帶（U.S. Corn Belt）農(nóng)業(yè)區(qū)CO2濃度觀測和WRF-STILT 模型開展研究，觀測站點高塔位于美國最大的農(nóng)業(yè)種植區(qū)，其100m高度處觀測的CO2濃度表現(xiàn)出強季節(jié)變化和夏季強日變化特征，生長季高塔濃度受下墊面內(nèi)C4農(nóng)作物和牧草強CO2吸收的影響，同時其區(qū)域內(nèi)還擁有多種類型的化石燃料燃燒排放源（石油提煉、道路交通、居民區(qū)、能源工業(yè)），其CO2濃度和通量綜合觀測為WRF-STILT模型模擬CO2濃度提供了驗證基礎(chǔ)。本研究的主要目的包括（1）評估WRF-STILT模型對農(nóng)業(yè)區(qū)觀測的CO2濃度強日變化特征的模擬能力，（2）評估CO2通量及其不確定性對模型模擬結(jié)果的影響。以期為中國建立高塔CO2濃度觀測，以及反演和評估區(qū)域尺度CO2通量提供技術(shù)和理論支撐。

1 材料與方法

1.1 觀測站點及資料來源

觀測站點位于美國（圖1a中紅色線條區(qū)域）玉米帶北部[11]，包括9個州，作為世界面積最大、產(chǎn)量最高的玉米種植區(qū)，玉米帶貢獻(xiàn)了全美80%和全球40%的玉米產(chǎn)量[20]。觀測塔位于明尼蘇達(dá)首府“明尼阿波利斯-圣保羅市”東南25km，塔高244m（44o41'19''N, 93o4'22''W；海拔高度290m，黑色‘+’所在位置），對CO2濃度進(jìn)行連續(xù)觀測的進(jìn)氣口位于離地面100m高度處，塔南為玉米種植區(qū)（圖1b中深黃色區(qū)域），塔北是牧草種植區(qū)（如圖1b淺綠色區(qū)域）。每小時對儀器進(jìn)行標(biāo)定[21]，標(biāo)定后的CO2濃度誤差小于0.03mmol×mol-1[22]，2007-2016年數(shù)據(jù)由明尼蘇達(dá)大學(xué)提供，鑒于2008年CO2濃度觀測數(shù)據(jù)相對其它年份的生長季缺測最少，且該年同時包含其它渦度相關(guān)等輔助觀測，所以本文選取2008年6-9月小時CO2濃度資料進(jìn)行模擬和對比研究。

在100m處安裝的渦度相關(guān)觀測系統(tǒng)，包括三維超聲風(fēng)速儀（型號CSAT3，Campbell Scientific）和閉路式濃度測量儀（型號TGA 100A，Campbell Scientific），對生態(tài)系統(tǒng)CO2通量進(jìn)行直接觀測，通量數(shù)據(jù)均經(jīng)過數(shù)據(jù)校正和質(zhì)量控制，詳見Griffis 等[22]；在32m和56m處同時進(jìn)行CO2濃度觀測（型號 TGA 100A，Campbell Scientific），用于計算從地面至100m高度的CO2儲存項，本文所使用的高塔生態(tài)系統(tǒng)凈交換數(shù)據(jù)為渦度相關(guān)直接觀測和儲存項之和[8,23]，觀測時段為與濃度觀測相同的2008年6-9月。

注：圖a中藍(lán)色、淺黃色和深黃色區(qū)域分別表示W(wǎng)RF三層嵌套氣象場，紅色線條區(qū)域代表美國各州邊界；圖b中深黃色區(qū)域為玉米作物，淺綠色為牧草，灰色為城市用地，測點高度100m

Note：The color of blue, light yellow and rich yellow represent 3 domains in WRF meterological setup, and the red line is states boundary for U.S.A in subfigure a; with yellow, light green and gray color indicating corn, pasture, and ruban land use categories in subfigure b

1.2 CO2濃度模擬

1.2.1 模擬大氣CO2濃度組成

在本研究中，模型中模擬的CO2濃度（CO2,m）由CO2初始場濃度（CO2,bg）、植被生態(tài)系統(tǒng)凈交換貢獻(xiàn)值（ΔCO2,NEE）與燃料燃燒貢獻(xiàn)值（ΔCO2,comb）組成[24-25]；其中燃料燃燒貢獻(xiàn)值又包含生物質(zhì)燃燒（ΔCO2,bb）和化石燃料燃燒（ΔCO2,ff），如簡化公式（1）-（3）所示，濃度單位為mmol×mol-1。

（2）

（3）

1.2.2 WRF-STILT模型

由1.2.1可知，WRF-STILT模型對足跡權(quán)重（footprint）的準(zhǔn)確計算是模擬CO2濃度的關(guān)鍵，而足跡權(quán)重是在WRF模型模擬和輸出的氣象場驅(qū)動下運行STILT模型得出。WRF模型是由美國國家大氣研究中心于20世紀(jì)90年代研發(fā)并經(jīng)不斷改進(jìn)的中尺度天氣預(yù)報模型；STILT（http://stilt-model.org/ pmwiki/pmwiki.php）是基于拉格朗日原理的粒子隨機游走模型。模型通過釋放大量的空氣粒子來模擬大氣的運動過程以計算足跡權(quán)重，即通量的倒數(shù)，代表單位CO2通量產(chǎn)生的CO2濃度貢獻(xiàn)值，它與下墊面所有格點的CO2通量相乘可得到CO2濃度貢獻(xiàn)值，如式（3）[13]。

STILT模型由WRF3.5模型輸出的高精度（小時分辨率率）氣象場驅(qū)動[11-12]，主要包括不同模擬層高度的三維風(fēng)速、氣壓、虛位溫度、相對濕度、空氣密度等，以及用于計算下墊面湍流交換參數(shù)的粗糙度長度、感熱、潛熱、摩擦風(fēng)速等，所有輸出參數(shù)見Nehrkorn等[23]。WRF模型所采用的參數(shù)化方案見文獻(xiàn)[11-12]，其所基于的氣象場初始和邊界條件采用NCEP（National Centers for Environmental Prediction）的FNL資料（http://rda.ucar.edu/datasets/ds083.2），氣象場的模擬時間段為與高塔CO2濃度相同的2008年的生長季6-9月，時間分辨率為小時。為給STILT模型提供更準(zhǔn)確的WRF3.5氣象驅(qū)動場，在WRF3.5模型中采用3層嵌套（Domain1、Domain2、Domain3）雙向反饋的設(shè)置，把3層嵌套的模擬區(qū)域的空間分辨率分別設(shè)置為27、9、3km，東西和南北格點數(shù)分別為250×180、385×409、670×532；3層模擬分別為圖1a中深藍(lán)、淺黃和深黃色區(qū)域。而在STILT模型中設(shè)置空氣粒子運動范圍為125-65oW，25-62oN，區(qū)域覆蓋美國本土。式（1）-式（3）可進(jìn)一步由式（4）表達(dá)，即

式中，CO2,m（xr, tr）為位于xr位置tr時刻模擬的CO2濃度，本研究中，模型設(shè)置xr為濃度觀測所在位置44o41'19''N，93o4'22''W和高度100m，tr則與高塔濃度觀測時間對應(yīng)，右邊第一項為CO2濃度貢獻(xiàn)項，表示氣體流經(jīng)上游所有模擬區(qū)域時，由于源匯項S(x，t)（即式（1）-式（3）中CO2通量）的影響產(chǎn)生的t0-tr時間段累積濃度貢獻(xiàn)值，因為模型顯示絕大多數(shù)粒子在7d前來源于相對干凈的背景場濃度區(qū)域（如東太平洋和加拿大北部），所以模型設(shè)置該累積過程時間為168h（即n=168）。第二項為初始場濃度項，其中為影響函數(shù)，代表氣流從流入模擬區(qū)域開始，經(jīng)過的所有區(qū)域CO2通量在累計時間段（t0-tr）對濃度模擬的影響，則是在初始t0時刻的影響函數(shù)。在STILT模型中通過設(shè)置每小時在100m高度處釋放500個粒子，計算在任意時間段和模擬區(qū)域內(nèi)任意格點所含粒子數(shù)占總粒子數(shù)的總停留時間比例，可得到最終需要的足跡權(quán)重，換算過程詳見Lin 等[13,26]。模型輸出的足跡權(quán)重的時間分辨率為小時，空間分辨率為0.1o。

1.2.3 CO2通量

如1.2.1所述，本研究所使用的CO2通量即源匯項S(x，t)由3部分組成，（1）化石燃料燃燒CO2通量（CO2,ff），來自EDGAR（Emission Database for Global Atmospheric Research，4.2 FT2010; European Commission，2009[28]）和Carbon Tracker[9]，兩者是被廣泛使用的化石燃料燃燒CO2排放源，EDGAR的空間分辨率為0.1°×0.1°，時間分辨率為年，為了得到與足跡權(quán)重和濃度觀測相同的小時CO2通量，基于“VULCAN”提供的小時變化系數(shù)得到小時尺度的化石燃料燃燒CO2通量[29]； Carbon Tracker的空間分辨率為1°×1°，時間分辨率為月，其化石燃料通量取自Miller（Boden 等[30]）和ODIAC（Carbon Dioxide Information and Analysis Center , Oda and Maksyutov[31]），取兩者的平均值，這兩種數(shù)據(jù)在國家尺度上差異很小，但是在區(qū)域尺度上卻有較大差異。為了定量研究EDGAR和Carbon Tracker兩種化石燃料燃燒CO2通量對濃度模擬結(jié)果的影響，本文將對比兩者濃度模擬的差異。

（2）生物質(zhì)燃燒（CO2,bb），其空間分辨率為1o×1o，時間分辨率為3h，由Carbon Tracker提供（ftp:// aftp.cmdl.noaa.gov/products/carbontracker/co2/fluxes/），是GFED4.1s（Global Fire Emissions Database）[32-33]和FINN（Fire Inventory from NCAR）[34]兩個數(shù)據(jù)集的平均。其中GFED4.1s（http://www.globalfiredata. org/）基于MODIS觀測到的過火面積與程度，CASA模擬的作物生物量以及釋放CO2的排放因子等計算得到；FINN（https://www2.acom.ucar.edu/modeling/ finn-fire-inventory-ncar）的火點信息和作物生物量數(shù)據(jù)均基于MODIS衛(wèi)星遙感觀測，結(jié)合排放因子，轉(zhuǎn)化為生物質(zhì)燃燒CO2通量；Carbon Tracker同化系統(tǒng)使用Mu 等[32]的方法降尺度至3h時間分辨率，由于其日變化很小，本文所使用的小時生物質(zhì)燃燒資料CO2通量來自其3h資料。

（3）植被生態(tài)系統(tǒng)凈交換，雖然渦度相關(guān)方法能進(jìn)行植被生態(tài)系統(tǒng)凈交換的直接觀測，但由于代表源區(qū)面積小，且代表下墊面類型單一，遠(yuǎn)不能滿足模型中所需要的模擬區(qū)域內(nèi)包含不同下墊面的CO2通量，所以本文所使用的模擬區(qū)域內(nèi)植被生態(tài)系統(tǒng)凈交換來自Carbon Tracker，它是一種基于大氣的濃度觀測，在CASA（Carnegie-Ames Stanford Approach）模型模擬的基礎(chǔ)上進(jìn)行優(yōu)化后的結(jié)果[9,35]，使植被生態(tài)系統(tǒng)凈交換更接近于真實值，其空間分辨率為1o×1o，時間分辨率為3h，本文用線性內(nèi)插的方法得到小時分辨率的通量，高塔渦度相關(guān)觀測數(shù)據(jù)將用于其所在格點的Carbon Tracker對比分析，以探討其不確定性的影響。

1.3 集合誤差

集合誤差（Aggregation error）[18]是指由于數(shù)據(jù)觀測和計算手段的限制，把非均質(zhì)下墊面的通量在格點里用區(qū)域平均值表示時所帶來的模擬大氣CO2濃度的誤差。采用Zhao等[10]的方法，本文基于已有高空間分辨率的0.1o×0.1oCO2通量數(shù)據(jù)，將其進(jìn)行平均得到不同空間分辨率的通量資料，再進(jìn)行CO2濃度模擬，利用CO2濃度貢獻(xiàn)值的差異計算相應(yīng)的集合誤差（S），即

式中，ΔCO2(0.10)和ΔCO2(x0)分別為使用空間分辨率0.1o和xoCO2通量數(shù)據(jù)時模擬的CO2濃度貢獻(xiàn)值（mmol×mol-1）。

2 結(jié)果與分析

2.1 生長季的足跡權(quán)重

在本研究中，WRF-STILT模型首先模擬了生長季的小時足跡權(quán)重。圖3為生長季4個月的平均足跡權(quán)重，由于足跡權(quán)重的空間變異性可達(dá)7個數(shù)量級，所以采用常用對數(shù)log10的形式在圖中展示。按照Chen 等[11-12]的定義，把足跡權(quán)重大于-4的區(qū)域定義為對高塔濃度觀測的強貢獻(xiàn)區(qū)，它代表對高塔100m觀測點CO2濃度模擬影響最大的范圍。由圖可知，由于湍流強度和平均風(fēng)速風(fēng)向的差異，不同月份的強貢獻(xiàn)區(qū)面積（9月＞8月＞6月＞7月）和形狀差異顯著。從其覆蓋范圍可知，不僅有位于塔南的農(nóng)業(yè)作物和塔北牧草植被生態(tài)系統(tǒng)凈交換影響，也有人為活動化石燃料燃燒的CO2源的貢獻(xiàn)。

2.2 CO2通量對濃度模擬的影響

2.2.1 不同化石燃料的影響

為評估Carbon Tracker 和EDGAR這兩種被廣泛使用的化石燃料燃燒排放源的差異，以觀測塔為中心，分別以4o、6o、10o、14o、20o為邊長所在區(qū)域（其中1o代表的長度約為90～100km），計算兩者2008年的區(qū)域平均化石燃料燃燒CO2通量，結(jié)果顯示，對應(yīng)區(qū)域的差異都在6%以內(nèi)，其中邊長為4o時兩者相對誤差為5.9%，隨著區(qū)域面積的增加，相對誤差呈現(xiàn)逐漸減小的趨勢，當(dāng)區(qū)域邊長為20o時相對誤差降至1.8%（圖4）。

為了單獨分析兩者空間分布差異對CO2濃度模擬的影響，先把0.1°×0.1°空間分辨率的EDGAR數(shù)據(jù)集合到與Carbon Tracker相同的1°×1°分辨率，再計算模擬得到的CO2濃度貢獻(xiàn)值的差異（如表1），結(jié)果顯示，除7月外（R2=0.65），兩者相關(guān)系數(shù)都大于0.97，表明相關(guān)性很高（P＜0.001），這主要是因為兩者的區(qū)域平均通量相近；比例系數(shù)的比較也表現(xiàn)出了相同的規(guī)律，除7月外（0.79），6、8、9月分別為0.87、0.91和0.91，說明在區(qū)域平均上兩者差異均小于6%，但是區(qū)域內(nèi)空間分布的差異會導(dǎo)致最后模擬CO2濃度貢獻(xiàn)值相差近10%。

表1 Carbon Tracker和EDGAR空間分布差異對CO2濃度增加值影響的比較

2.2.2 集和誤差對CO2濃度模擬的影響

研究區(qū)域位于下墊面均一的農(nóng)業(yè)區(qū)，故選擇0.1o高空間分辨率的EDGAR農(nóng)業(yè)土壤，其CO2通量可近似代表生態(tài)系統(tǒng)凈交換的空間分布狀況，從而定量探討化石燃料燃燒和植被生態(tài)系統(tǒng)凈交換集合誤差對模擬CO2濃度的影響。所以基于化石燃料和農(nóng)業(yè)土壤0.1o的CO2通量數(shù)據(jù)，集合到13組不同空間分辨率（0.2o，0.3o，0.4o，0.5o，0.6o，0.7o，0.8o，0.9o，1o，1.5o，2o，3o，4o），再模擬得到CO2濃度貢獻(xiàn)值。

結(jié)果表明（圖5），化石燃料燃燒CO2通量空間分辨率與模擬的生長季平均CO2濃度貢獻(xiàn)值為冪函數(shù)關(guān)系，隨著分辨率的降低，集合誤差增大（集合誤差率越?。?；當(dāng)空間分辨率為1o時，模擬的CO2濃度貢獻(xiàn)值僅0.1o的0.4倍，說明若使用空間分辨率為1o的Carbon Tracker，或者更低空間分辨率的化石燃料燃燒CO2源，模擬結(jié)果會嚴(yán)重偏低。而由于在農(nóng)業(yè)種植區(qū)，下墊面均一，不同空間分辨率農(nóng)業(yè)土壤釋放CO2通量模擬得到的CO2濃度貢獻(xiàn)值變化都在5%以內(nèi)，說明生態(tài)系統(tǒng)凈交換的集合誤差和分辨率并無明顯的關(guān)系。所以在模擬CO2通量過程中，對于下墊面空間異質(zhì)性高的人為化石源排放區(qū)域，應(yīng)使用更高空間分辨率的CO2通量數(shù)據(jù)。

注：誤差線為6-9月的標(biāo)準(zhǔn)差

Note：The error bar represents standard deviation in the 4 months

2.2.3 生態(tài)系統(tǒng)凈交換的不確定性

為了評估Carbon Tracker生態(tài)系統(tǒng)凈交換的誤差，對高塔100m處渦度相關(guān)的觀測數(shù)據(jù)與其所在格點的Carbon Tracker的生態(tài)系統(tǒng)凈交換通量進(jìn)行對比分析，選擇白天（10：00-15：00）和夜晚（1：00-6：00）兩個時間段進(jìn)行對比。如圖6所示，散點表示在生長季相對應(yīng)的3h平均值。雖然Carbon Tracker的格點代表的是邊長為1o（-90km）區(qū)域的平均，而100m高塔渦度相關(guān)的通量平均源區(qū)半徑小于10km[36-37]，但是由于該農(nóng)業(yè)種植區(qū)下墊面均一，可認(rèn)為兩者都能觀測到相同下墊面類型的CO2通量。由圖可見，對于夜晚和白天，其對應(yīng)的擬合斜率分別為1.56和2.26，說明所使用的生態(tài)系統(tǒng)凈交換偏高，而兩者呈極顯著相關(guān)（P<0.001）。觀測點西北邊的牧草區(qū)域由于沒有直接的渦度相關(guān)觀測，所以本研究不作分析。

2.2.4 生態(tài)系統(tǒng)凈交換與CO2濃度模擬的定量關(guān)系

由于觀測塔北部為牧草區(qū)域，南部為玉米種植區(qū)，兩者植被類型差異大，所以風(fēng)向或下墊面類型的誤差帶來差異顯著的生態(tài)系統(tǒng)凈交換信號。為了分析下墊面作物類型和其相應(yīng)的通量對CO2濃度貢獻(xiàn)值的影響，在模型中分別把所有下墊面換作牧草和玉米地，其通量分別對應(yīng)Carbon Tracker中牧草和玉米地生態(tài)系統(tǒng)凈交換的通量，結(jié)果如表2所示。由表可見，8月下墊面為牧草，9月下墊面為玉米時，兩者的月平均生態(tài)系統(tǒng)凈交換均為負(fù)，然而玉米下墊面的濃度貢獻(xiàn)值為負(fù)，牧草下墊面的濃度貢獻(xiàn)值為正，這反映了不同月份的足跡權(quán)重的差異。整個生長季草地的平均生態(tài)系統(tǒng)凈交換為-1.29μmol·m-2·s-1（碳匯），然而其CO2濃度貢獻(xiàn)值為正（2.46mmol×mol-1），說明對于CO2濃度，夜晚正生態(tài)系統(tǒng)凈交換的貢獻(xiàn)權(quán)重高于白天負(fù)通量，這是因為夜晚邊界層低于白天，所以導(dǎo)致即使是相同的CO2通量卻產(chǎn)生了差異顯著的濃度貢獻(xiàn)值。而當(dāng)下墊面全為玉米時，7、8月的濃度貢獻(xiàn)值分別為-9.09mmol×mol-1和-11.19mmol×mol-1，且其平均CO2濃度貢獻(xiàn)值與下墊面全為草地的差異達(dá)到12mmol×mol-1，說明下墊面分類產(chǎn)生的生態(tài)系統(tǒng)凈交換差異和模擬風(fēng)向誤差會導(dǎo)致CO2濃度的顯著模擬偏差。

表2 不同下墊面的生態(tài)系統(tǒng)凈交換（NEE）及其濃度貢獻(xiàn)值

2.3 CO2濃度模擬

基于CO2通量及其不確定性對CO2濃度模擬的影響分析，本研究選取0.1o空間分辨率的EDGAR化石燃料燃燒通量和Carbon Tracker生態(tài)系統(tǒng)凈交換CO2通量進(jìn)行濃度模擬研究，如圖7a所示，WRF-STILT能夠模擬出CO2強的日變化特征，結(jié)果顯示，整個生長季觀測和模擬的CO2濃度均方根誤差為13.70mmol×mol-1。生長季平均生物質(zhì)燃燒濃度貢獻(xiàn)為0.07mmol×mol-1，相對于化石燃料濃度貢獻(xiàn)值6.43mmol×mol-1，可以被忽略，但是在某些時段，如6月30日，生物質(zhì)燃燒濃度貢獻(xiàn)的最大值達(dá)到1.39mmol×mol-1，表明北方的火點燃燒在特定天氣條件下也會對高塔的CO2濃度產(chǎn)生明顯影響。植被生態(tài)系統(tǒng)凈交換的濃度貢獻(xiàn)值為-0.36mmol×mol-1，而CO2的背景值濃度為381.99mmol×mol-1。對比WRF-STILT模擬的濃度貢獻(xiàn)值，更能客觀反映WRF-STILT真實的模擬能力。從圖7b可以看出，模擬的濃度貢獻(xiàn)值與觀測值相關(guān)系數(shù)高（R=0.52，P<0.001），由模擬方程可知，其模擬值偏高7.26mmol×mol-1。

3 結(jié)論與討論

（1）本研究通過WRF-STILT模型，對位于美國玉米帶的高塔100m處CO2濃度進(jìn)行模擬，結(jié)果表明，模型能夠很好地模擬出生長季CO2濃度日變化特征，但模擬結(jié)果偏高7.26mmol×mol-1，這可能是模型氣象場（風(fēng)速和邊界層高度等）和所使用的CO2通量誤差導(dǎo)致的（石油化石燃料燃燒的高估和植被生態(tài)系統(tǒng)凈交換的低估）。刁一偉等[38]使用WRF-VPRM模型，模擬了位于長三角區(qū)域地面站點夏季6d的CO2濃度，發(fā)現(xiàn)模擬結(jié)果偏低5～15mmol×mol-1，認(rèn)為模擬誤差由氣象場和CO2通量偏差導(dǎo)致，其誤差大于本研究的方法，且其研究時段短。Mallia 等[25]使用WRF-STILT模型對鹽湖城的高塔CO2濃度模擬結(jié)果顯示，模型能很好地表示其日變化特征，而日變化的主要貢獻(xiàn)是由化石燃料燃燒導(dǎo)致的，而本研究由于在農(nóng)業(yè)區(qū)進(jìn)行，所以植被生態(tài)系統(tǒng)凈交換是CO2濃度的主要貢獻(xiàn)。

（2）通過定量評估CO2通量對模型模擬的影響，發(fā)現(xiàn)CO2通量大小和分辨率的差異是模擬結(jié)果的主要誤差來源，這是因為其空間異質(zhì)性大，當(dāng)觀測塔附近有明顯的強點源化石燃料燃燒釋放CO2時，低空間分辨率會降低其對高塔濃度的影響，進(jìn)而導(dǎo)致模擬的CO2濃度偏低，使用1o空間分辨率的化石燃料CO2通量模擬的CO2濃度增加值僅0.1o模擬結(jié)果的0.4倍，且分辨率越低，模擬誤差越大；因此集合誤差是模擬高塔CO2濃度和反演CO2通量的潛在誤差來源。Kaminski等[18-19]也強調(diào)了分析集合誤差的重要性，Zhao 等[10]發(fā)現(xiàn)不同分辨率的CH4通量同樣會導(dǎo)致模擬和反演結(jié)果的誤差，本研究更定量分析了誤差大小隨分辨率變化的關(guān)系，建議今后基于WRF-STILT模型對其它氣體進(jìn)行濃度模擬和通量反演時，首先需要對所使用通量進(jìn)行集合誤差分析。模型使用的白天和夜晚生態(tài)系統(tǒng)凈交換分別是高塔渦度相關(guān)觀測的2.26和1.56倍，而這是模擬CO2濃度偏高的原因，麥博儒等[39]分別使用代表珠江三角洲的3個不同類型下墊面的生態(tài)系統(tǒng)凈交換觀測值，與Carbon Tracker模擬值進(jìn)行對比，發(fā)現(xiàn)雖然兩者相關(guān)系數(shù)較高，但是后者同樣會高于觀測值。Carbon Tracker和EDGAR兩種重要的化石燃料通量，其區(qū)域平均值相差＜6%，而對模擬的CO2濃度增加值的平均差異約為10%，這主要是由于兩者區(qū)域平均和空間分布差異所導(dǎo)致，Miller等[40-41]同樣認(rèn)為EDGAR燃燒對模擬的CO2濃度增加值誤差在國家尺度上為5%～10%，但在區(qū)域尺度上，其不確定性甚至?xí)h(yuǎn)遠(yuǎn)超出該范圍，所以降低其CO2通量的誤差是未來準(zhǔn)確模擬CO2濃度的基礎(chǔ)。研究結(jié)果也表明下墊面類型的差異產(chǎn)生的植被生態(tài)系統(tǒng)凈交換模擬的誤差，同樣是導(dǎo)致模擬CO2偏差的原因。

研究結(jié)果表明WRF-STILT模型具有強的CO2濃度模擬能力，CO2通量的偏差及其在空間分布上的差異是主要誤差來源，該方法表明WRF-STILT模型在濃度模擬和區(qū)域通量反演上具有可行性，可為將來在中國農(nóng)業(yè)種植區(qū)建立高塔CO2濃度觀測網(wǎng)絡(luò)，反演區(qū)域尺度的生態(tài)系統(tǒng)凈交換和計算區(qū)域植被生產(chǎn)力提供理論和技術(shù)指導(dǎo)。由于觀測的原因，本研究并未分析WRF模型輸出的氣象驅(qū)動場與實際觀測的差異（主要包括風(fēng)速和邊界層高度），將在未來研究中加以開展。

[1]于貴瑞,王秋鳳,朱先進(jìn).區(qū)域尺度陸地生態(tài)系統(tǒng)碳收支評估方法及其不確定性[J].地理科學(xué)進(jìn)展,2011,30(1):103-113.

Yu G R,Wang Q F,Zhu X J.Method and uncertainties in evaluating the carbon budgets of regional terrestrial ecosystems[J].Progress in Geography,2011,30(1):103-113.(in Chinese)

[2]Maki T,Ikegami M,Fujita T,et al.New technique to analyse global distributions of CO2,concentrations and fluxes from non-processed observational data[J].Tellus Series B-chemical & Physical Meteorology, 2010, 62(5):797-809.

[3]Saeki T,Maksyutov S,Sasakawa M,et al.Carbon flux estimation for Siberia by inverse modeling constrained by aircraft and tower CO2,measurements[J].Journal of Geophysical Research Atmospheres, 2013,118(2):1100-1122.

[4]Houghton R A.Balancing the global carbon budget[J].Earth and Planetary Sciences,2007,35(35):313-347.

[5]Quéré C L,Raupach M R,Canadell J G,et al.Trends in the sources and sinks of carbon dioxide[J].Nature Geoscience, 2009,2(12):831-836.

[6]于成龍,劉丹.小興安嶺天然闊葉混交林生長季CO2通量特征分析[J].中國農(nóng)業(yè)氣象,2011,32(4):525-529.

Yu C L,Liu D.Analysis on CO2flux during growth season of natural broadleaved mixed forest in Xiaoxinganling Mountains[J].Chinese Journal of Agrometeorology,2011, 32(4):525-529.(in Chinese)

[7]李永秀,楊再強,張富存,等. 光合作用模型在長江下游冬麥區(qū)的適用性研究[J].中國農(nóng)業(yè)氣象,2011,32(4):588-592.

Li Y X,Yang Z Q,Zhang F C,et al.Applicability of different photosynthesis models for winter wheat in the lower Yangtze river[J].Chinese Journal of Agrometeorology,2011,32(4): 588-592.(in Chinese)

[8]Chen M,Griffis T J,Baker J,et al.Simulating crop phenology in the Community Land Model and its impact on energy and carbon fluxes[J].Journal of Geophysical Research- Biogeosciences,2015,120(2):310-325.

[9]Peters W,Jacobson A R,Sweeney C,et al.An atmospheric perspective on North American carbon dioxide exchange: CarbonTracker[J].Proceedings of the National Academy of Sciences of the United States of America,2007,104(48): 18925-18930.

[10]Zhao C,Andrews A E,Bianco L,et al.Atmospheric inverse estimates of methane emissions from Central California[J]. Journal of Geophysical Research Atmospheres, 2009, 114(D16):4723-4734.

[11]Chen Z,Griffis T J,Millet D B,et al.Partitioning N2O Emissions within the US Corn Belt using an Inverse Modeling Approach[J].Global Biogeochemical Cycles, 2016,30:1192-1205.

[12]Su Y K,Millet D B,Hu L,et al.Constraints on carbon monoxide emissions based on tall tower measurements in the U.S. upper midwest[J].Environmental Science & Technology, 2013,47(15):8316-24.

[13]Lin J C,Gerbig C,Wofsy S C,et al.A near-field tool for simulating the upstream influence of atmospheric observations:the Stochastic Time-Inverted Lagrangian Transport (STILT) model[J].Journal of Geophysical Research Atmospheres,2003,108(4493):1211-1222.

[14]Gerbig C,Lin J C,Wofsy S C,et al.Toward constraining regional-scale fluxes of CO2,with atmospheric observations over a continent:analysis of COBRA data using a receptor-oriented framework[J].Journal of Geophysical Research Atmospheres,2003,108(D24):4757-ACH6.

[15]Alexandratos N,Bruinsma A J.World agriculture towards 2030/2050: the 2012 revision[M].Rome:Food and Agriculture Organization of the United Nations,2008.

[16]潘根興,趙其國.我國農(nóng)田土壤碳庫演變研究:全球變化和國家糧食安全[J].地球科學(xué)進(jìn)展,2005,20(4):384-393.

Pan G X,Zhao Q G.Study on evolution of organic carbon stock in agricultural soils of China:facing the challenge of global change and food security[J].Advances in Earth Science,2005, 20(4):384-393.(in Chinese)

[17]陳紀(jì)波,胡慧,陳克垚,等.基于非線性PLSR模型的氣候變化對糧食產(chǎn)量的影響分析[J].中國農(nóng)業(yè)氣象,2016,37(6): 674-681.

Chen J B,Hu H,Chen K Y,et al.Effects of climate change on the grain yield based on nonlinear PLSR model[J].Chinese Journal of Agrometeorology,2016,37(6):674-681.(in Chinese)

[18]Kaminski T,Rayner P J,Heimann M,et al.On aggregation errors in atmospheric transport inversions[J].Journal of Geophysical Research Atmospheres,2001,106(D5):4703- 4715.

[19]Turner A J,Jacob D J.Balancing aggregation and smoothing errors in inverse models[J].Atmospheric Chemistry & Physics,2015,15(1):1001-1026.

[20]Wang E,Little B B,Williams J R,et al.Simulation of hail effects on crop yield losses for corn-belt states in USA[J].Transactions of the CSAE,2012,258(21):10012- 10016.

[21]高韻秋,劉壽東,胡凝,等.南京夏季城市冠層大氣CO2濃度時空分布規(guī)律的觀測[J].環(huán)境科學(xué),2015(7):2367-2373.

Gao Y Q,Liu S D,Hu N,et al.Direct observation on the temporal and spatial patterns of the CO2concentration in the atmospheric of Nanjing urban canyon in summer[J]. Environmental Science,2015(7):2367-2373.(in Chinese)

[22]Griffis T J,Baker J M,Sargent S D,et al.Influence of C4 vegetation on13CO2discrimination and isoforcing in the upper Midwest,United States[J].Global Biogeochemical Cycles,2010,24(4):16.

[23]Nehrkorn T,Eluszkiewicz J,Wofsy S C,et al.Coupled weather research and forecasting-stochastic time-inverted lagrangian transport (WRF–STILT) model[J].Meteorology & Atmospheric Physics,2010,107(1):564.

[24]Ahmadov R,Gerbig C,Kretschmer R,et al.Comparing high resolution WRF-VPRM simulations and two global CO2transport models with coastal tower measurements of CO2[J]. Biogeosciences,2009,6(5):807-817.

[25]Mallia D V,Lin J C,Urbanski S,et al.Impacts of upwind wildfire emissions on CO,CO2,and PM2.5, concentrations in Salt Lake City,Utah[J].Journal of Geophysical Research Atmospheres,2015, 120(1):147-166.

[26]Gerbig C,Lin J C,Wofsy S C,et al.Towards constraining regional scale fluxes of CO2with atmospheric observations over a continent 1:observed spatial variability from airborne platforms[J].Journal of Geophysical Research Atmospheres, 2003,108(D24):ACH 5-1.

[27]Pillai D,Gerbig C,Kretschmer R,et al.Comparing lagrangian and eulerian models for CO2transport – a step towards Bayesian inverse modeling using WRF/STILT-VPRM[J]. Atmospheric Chemistry & Physics,2012,12(1):1267-1298.

[28]European Commission,Joint Research Centre/Netherlands Environmental Assessment Agency. Emission Database for Global Atmospheric Research (EDGAR),release version 4.0[M].2009.

[29]Gurney K R,Mendoza D L,Zhou Y,et al.The vulcan project:high resolution fossil fuel combustion CO2emissions fluxes for the United States[J].Environmental Science & Technology,2009,43(14):5535-5541.

[30]Boden T A,Marland G,Andres R J.Global,regional,and national fossil-fuel CO2emissions[J].Carbon Dioxide Information Analysis Center, Oak Ridge National Laboratory, U.S.Department of Energy,Oak Ridge, Tenn.,U.S. A.,2013.

[31]Oda T,Maksyutov S.A very high-resolution (1km×1km) global fossil fuel CO2emission inventory derived using a point source database and satellite observations of nighttime lights[J].Atmospheric Chemistry & Physics,2010,10(7): 16307-16344.

[32]Mu M,Randerson J T,Van d W G R,et al.Daily and 3‐hourly variability in global fire emissions and consequences for atmospheric model predictions of carbon monoxide[J]. Journal of Geophysical Research Atmospheres, 2011, 116(D24):191-200.

[33]Giglio L,Randerson J T,Werf G R V D.Analysis of daily, monthly,and annual burned area using the fourth-generation global fire emissions database (GFED4)[J].Journal of Geophysical Research Biogeosciences,2013,118(1):317-328.

[34]Wiedinmyer C,Akagi S K,Yokelson R J,et al.The Fire INventory from NCAR (FINN):a high resolution global model to estimate the emissions from open burning[J]. Geoscientific Model Development Discussions,2010, 3(4): 625-641.

[35]Potter C S,Randerson J T,Field C B,et al.Terrestrial ecosystem production:a process model based on global satellite and surface data[J].Global Biogeochemical Cycles,1993,7(4):811-841.

[36]Zhang X,Lee X,Griffis T J,et al.Estimating regional greenhouse gas fluxes: an uncertainty analysis of planetary boundary layer techniques and bottom-up inventories[J]. Atmospheric Chemistry & Physics,2014,14(19): 10705- 10719.

[37]Davis K J,Bakwin P S,Yi C,et al.The annual cycles of CO2,and H2O exchange over a northern mixed forest as observed from a very tall tower[J].Global Change Biology, 2003,9(9):1278-1293.

[38]刁一偉,黃建平,劉誠,等.長江三角洲地區(qū)凈生態(tài)系統(tǒng)二氧化碳通量及濃度的數(shù)值模擬[J].大氣科學(xué),2015,39(5): 849-860.

Diao Y W,Huang J P,Liu C,et al.A modeling study of CO2flux and concentrations over the Yangtze River Delta using the WRF-GHG model[J].Chinese Journal of Atmospheric Sciences,2015, 39(5):849-860.(in Chinese)

[39]麥博儒,安興琴,鄧雪嬌,等.珠江三角洲近地層CO2通量模擬分析與評估驗證[J].中國環(huán)境科學(xué),2014,34(8): 1960-1971.

Mai B R,An X Q,Deng X J,et al.Simulation analysis and verification of surface CO2flux over Pearl River Delta,China[J].China Enviromental Science,2014,34(8): 1960-1971.(in Chinese)

[40]Miller J B,Lehman S J,Montzka S A,et al.Linking emissions of fossil fuel CO2,and other anthropogenic trace gases using atmospheric 14CO2[J].Journal of Geophysical Research Atmospheres,2012,117(D8):8302.

[41]Nassar R,Napier-Linton L,Gurney K R,et al.Improving the temporal and spatial distribution of CO2emissions from global fossil fuel emission data sets[J].Journal of Geophysical Research Atmospheres,2013,118(2):917-933.

Effect of Flux and its Uncertainty on Tall Tower CO2Concentration Simulation in the Agricultural Domain

HU Cheng1,2, ZHANG Mi1,2, XIAO Wei1,2, WANG Yong-wei1,2,WANG Wei1,2, TIM Griffis3,LIU Shou-dong1,2, LI Xu-hui1

(1.Yale-NUIST Center on Atmospheric Environment, Nanjing University of Information Science and Technology, Nanjing 210044, China;2. Collaborative Innovation Center of Atmospheric Environment and Equipment Technology, Nanjing University of Information Science and Technology, Nanjing 210044, China;3. University of Minnesota-Twin Cities, Saint Paul 55108, U.S.A)

Based on the CO2concentration observations in U.S. corn belt, which was measured at 100m height of a tall tower, hourly CO2concentration was simulated for the growing season (June–September, 2008) with the WRF-STILT model. And the effect of flux uncertainty on modeled CO2concentration was also analyzed. The results showed as: (1) WRF-STILT model can simulate the observed strong diurnal variation in growing season, with RMSE be 13.70mmol×mol-1, and it was overestimated by 7.26mmol×mol-1, the shape and area of intense footprint zonesare different for different months (September＞August＞June＞July) .(2) The difference of regional average anthropogenic CO2flux for EDGAR and Carbon Tracker was within 6%, when both of them were at the same spatial resolution, the simulated CO2enhancement difference was close to 10%. (3) Spatial resolution can lead to large bias in the modeled CO2enhancement, when using 1oemissions, the simulated CO2enhancement was only 0.4 times of the results using 0.1oemissions, and with the decreases of spatial resolution, the modeled bias increases. (4) Daytime and nighttime NEE of Carbon Tracker is 2.26 and 1.56 times that of tall tower NEE observations, and the misrepresentatives of underlying land use categories can lead to about 12mmol×mol-1bias in the modeled results, which may be the potential reason of bias high for 7.26mmol×mol-1. Our study concludes that when combing WRF-STILT model with high quality CO2flux, the strong diurnal variation of CO2concentration can be well simulated, and the uncertainty of CO2flux is the main reason for modeled CO2concentration bias, it also indicates the potential of evaluating and retrieving prior CO2flux.

WRF-STILT model；Eddy covariance；Fossil emissions；Flux uncertainty

10.3969/j.issn.1000-6362.2017.08.001

2016-12-21

。E-mail: zhangm.80@nuist.edu.cn

國家自然科學(xué)基金項目（41575147；41475141；41505005）；江蘇省高校優(yōu)勢學(xué)科建設(shè)工程項目（PAPD）；教育部長江學(xué)者和創(chuàng)新團隊發(fā)展計劃項目（PCSIRT）；2016年度江蘇省高校研究生科技創(chuàng)新項目（1354051601006）；國家公派聯(lián)合培養(yǎng)博士研究生項目（201508320287）

胡誠（1989-），博士生，主要研究方向為基于高塔濃度觀測的區(qū)域尺度溫室氣體通量反演。E-mail: nihaohucheng@163.com

胡誠,張彌,肖薇,等.通量及其不確定性對農(nóng)業(yè)區(qū)高塔CO2濃度模擬的影響[J].中國農(nóng)業(yè)氣象,2017,38(8):469-480