混合標(biāo)氣標(biāo)校溫室氣體數(shù)據(jù)的計(jì)算方法與實(shí)現(xiàn)

藍(lán)維嘉，韋漢勇，詹棠，邱國(guó)強(qiáng)

（珠海市國(guó)家氣候觀象臺(tái)，廣東珠海 519000）

精確的溫室氣體觀測(cè)數(shù)據(jù)能客觀反映相關(guān)區(qū)域碳排放狀況，是各類碳研究模式得以有效運(yùn)行、各級(jí)政府及職能部門應(yīng)對(duì)氣候變化、實(shí)施節(jié)能減排提供決策的依據(jù)［1－2］。目前國(guó)內(nèi)使用較多的Picarro系列溫室氣體在線分析儀（簡(jiǎn)稱“分析儀”）是基于波長(zhǎng)掃描光腔衰蕩光譜技術(shù)（Wavelength Scan Cavity Ring Down Spectroscopy，WS-CRDS）的高精度觀測(cè)系統(tǒng)［3－4］。但是分析儀在連續(xù)觀測(cè)期間會(huì)發(fā)生微小的漂移，因此對(duì)于痕量級(jí)的CO2（單位量級(jí)parts per million by volume，ppmv）、CH4（單位量級(jí)parts per billion by volume，ppbv）等觀測(cè)要素，原始數(shù)據(jù)的質(zhì)量控制和標(biāo)校訂正必不可少。張東啟等［5］通過保證觀測(cè)系統(tǒng)穩(wěn)定性和對(duì)標(biāo)氣進(jìn)行逐級(jí)嚴(yán)格標(biāo)定來實(shí)現(xiàn)CO2、CH4觀測(cè)數(shù)據(jù)的質(zhì)量控制；劉立新等［6］建立了科學(xué)規(guī)范的本底大氣CO2采樣觀測(cè)過程的質(zhì)量保證與質(zhì)量控制方法；鄧雪嬌等［7］對(duì)廣東省溫室氣體觀測(cè)數(shù)據(jù)進(jìn)行了對(duì)比分析，對(duì)溫室氣體觀測(cè)質(zhì)控方法進(jìn)行研究，掌握了觀測(cè)標(biāo)準(zhǔn)傳遞、標(biāo)校方法流程，建立了溫室氣體可視化應(yīng)用平臺(tái)；臧昆鵬等［8］組裝新型CO2和CH4混合標(biāo)氣標(biāo)校系統(tǒng)，研究建立了簡(jiǎn)便高效的標(biāo)校流程和方法。

基層監(jiān)測(cè)站往往不具備國(guó)家級(jí)標(biāo)校實(shí)驗(yàn)的條件，無法實(shí)時(shí)對(duì)設(shè)備進(jìn)行專業(yè)的校準(zhǔn)，混合標(biāo)氣標(biāo)校法便是各站對(duì)數(shù)據(jù)進(jìn)行后期標(biāo)校訂正的一種常用方法。對(duì)此，各站基本依托外部第三方開發(fā)的應(yīng)用系統(tǒng)，但由于其源代碼未公開，各個(gè)流程不可見、不可控，這不僅使應(yīng)用系統(tǒng)維護(hù)不便，還制約著臺(tái)站人員對(duì)數(shù)據(jù)的質(zhì)量控制和應(yīng)用。因此，開發(fā)出技術(shù)可控、運(yùn)維方便的溫室氣體數(shù)據(jù)混合標(biāo)氣標(biāo)校軟件或工具，就顯得十分必要。本研究應(yīng)用JAVA語(yǔ)言編譯了混合標(biāo)氣標(biāo)校溫室氣體數(shù)據(jù)的計(jì)算機(jī)方法［9］，對(duì)各個(gè)計(jì)算流程進(jìn)行詳細(xì)介紹（圖1）。

圖1 溫室氣體數(shù)據(jù)混合標(biāo)氣、標(biāo)校法流程設(shè)計(jì)示意圖

1 混合標(biāo)氣、標(biāo)校溫室氣體數(shù)據(jù)的流程

1.1 原始數(shù)據(jù)的質(zhì)量控制及平均值計(jì)算

1）原始數(shù)據(jù)文件及關(guān)鍵字段。

樣氣（空氣）經(jīng)過二級(jí)過濾、二級(jí)除水、流量控制等處理，在分析儀內(nèi)進(jìn)行WS-CRDS分析后，生成名為CFADSaaaa-yyyyMMdd-HHmmssZ-Data-Log＿User＿Sync（aaaa為設(shè)備序列號(hào)）的原始數(shù)據(jù)文件（每5 s一行數(shù)據(jù)），關(guān)鍵字段如下：

（1）DATE、TIME：日期、時(shí)間（世界時(shí)）；

（2）ALARM＿STATUS：報(bào)警狀態(tài)，0表示分析儀狀態(tài)正常，否則為異常提示；

（3）MPVPosition：多口閥代碼。以珠海為例，1為樣氣、4為高質(zhì)量濃度標(biāo)氣WH、5為中質(zhì)量濃度目標(biāo)氣T、6為低質(zhì)量濃度標(biāo)氣WL；WH（WL）的標(biāo)稱值一般應(yīng)高于（低于）站點(diǎn)樣氣質(zhì)量濃度的最高（最低）水平。

（4）CO2＿dry＿sync、CH4＿dry＿sync：經(jīng)水汽校正后，樣氣的CO2、CH4質(zhì)量濃度實(shí)測(cè)值［10］。

2）氣路切換與相關(guān)計(jì)算方法。

按照生成時(shí)間，遍歷所有原始數(shù)據(jù)文件，逐行解析。首先需要判斷本行的多口閥代碼與上一行是否一致，條件判別語(yǔ)句為：

if（（mpv＿1－mpv＿0）?。?）｛／／mpv＿1、mpv＿0：分別為本行、上一行數(shù)據(jù)的多口閥代碼

多口閥代碼不一致的，表明標(biāo)氣與樣氣的氣路已通過多口閥切換。氣路切換后2 min內(nèi)的數(shù)據(jù)應(yīng)剔除，因?yàn)闅饴分械臍怏w是標(biāo)氣和樣氣的混合態(tài)（對(duì)已運(yùn)行時(shí)間較長(zhǎng)或者響應(yīng)速率變慢的分析儀，以3 min為宜）；此后，氣路中的氣體質(zhì)量濃度基本達(dá)到平衡（圖2）。氣路切換與否，是計(jì)算方法變換的重要依據(jù)。

圖2 氣路切換前后的CO2質(zhì)量濃度變化曲線

本研究歸納了5種氣路切換情形，其對(duì)應(yīng)的計(jì)算步驟和計(jì)算機(jī)代碼如下：

（1）情形1。氣路不變，測(cè)標(biāo)氣。

①判斷多口閥位置：mpv＿0＝＝mpv＿1、mpv＿0＞1、mpv＿1＞1。

②確保距上次氣路切換時(shí)間已經(jīng)超過2～3 min（120～180 s）：

（sec.ss（t＿str＿1，t＿mpv＿changed））＞＝180；其中，sec.ss（）是判斷兩個(gè)時(shí)間相差秒數(shù)的方法，t＿str＿1為本行數(shù)據(jù)時(shí)間，t＿mpv＿changed為最近一次氣路切換時(shí)間。

③判斷本行數(shù)據(jù)時(shí)間與上一行數(shù)據(jù)時(shí)間處在同一個(gè)標(biāo)氣觀測(cè)時(shí)段（間隔小于10 s）。

④將本行標(biāo)氣值寫入標(biāo)氣原始觀測(cè)樣本數(shù)組（每個(gè)數(shù)組包含一個(gè)標(biāo)氣觀測(cè)時(shí)段的后2 min的數(shù)據(jù)，每5 s一個(gè)觀測(cè)數(shù)據(jù)，共計(jì)24個(gè)數(shù)據(jù)）：

co2＿mpv［k＿5min－1］＝Double.valueOf（zz.substring（416，428））*100；

ch4＿mpv［k＿5min－1］＝Double.valueOf（zz.substring（468，480））*1000；

其中，co2＿mpv［］、ch4＿mpv［］為CO2、CH4的標(biāo)氣原始觀測(cè)樣本數(shù)組。

（2）情形2。氣路不變，測(cè)樣氣，且要求本行與上一行數(shù)據(jù)的時(shí)間在同一個(gè)5 min時(shí)段——“5 min時(shí)段”，指的是正5 min時(shí)刻（含）之前的300 s。例如，00：00：01—00：05：00的60個(gè)樣氣有效數(shù)據(jù)的平均值，表示為“00：05：00的平均值”；若因氣路切換、數(shù)據(jù)中斷等原因造成5 min時(shí)段的樣氣有效數(shù)據(jù)不滿60個(gè)，如00：00：01—00：04：35或00：02：30—00：05：00的樣氣有效數(shù)據(jù)的平均值，均表示為“00：05：00的平均值”。

①判斷多口閥位置：mpv＿0＝＝1、mpv＿1＝＝1；

②氣路切換時(shí)間已經(jīng)超過3 min（180 s）；

③判斷本行與上一行數(shù)據(jù)在同一個(gè)5 min時(shí)段：sec.ss（w5＿1，w5＿0）＝＝0；其中：w5＿1、w5＿0分別為本行、上一行數(shù)據(jù)所對(duì)應(yīng)的正5 min時(shí)刻。

④5 s樣氣實(shí)測(cè)值累加：

co2＿5min＿all＝co2＿5min＿all＋Double.valueOf（zz.substring（416，428））*100；

ch4＿5min＿all＝ch4＿5min＿all＋Double.valueOf（zz.substring（468，480））*1000；

k＿5min＋＋。

（3）情形3。氣路不變，測(cè)樣氣，但是本行與上一行數(shù)據(jù)的時(shí)間已經(jīng)不在同一個(gè)5 min時(shí)段。

①判斷多口閥位置：mpv＿0＝＝1、mpv＿1＝＝1；

②判斷本行與上一行數(shù)據(jù)不在同一個(gè)5 min時(shí)段：sec.ss（w5＿1，w5＿0）！＝0；

③計(jì)算樣氣5 min平均值：

data＿co2＿ave＝Double.valueOf（String.format（″%.4f″，co2＿5min＿all／k＿5min））；

data＿ch4＿ave＝Double.valueOf（String.format（″%.4f″，ch4＿5min＿all／k＿5min））；

④將上一行數(shù)據(jù)對(duì)應(yīng)的正5 min時(shí)刻作為樣氣5 min平均值的探測(cè)時(shí)間。

（4）情形4。氣路從樣氣切換至標(biāo)氣。

①判斷多口閥位置：mpv＿0＝＝1、mpv＿1＞1；

②計(jì)算樣氣5 min平均值；

③將上一行數(shù)據(jù)對(duì)應(yīng)的正5 min時(shí)刻作為樣氣5 min平均值的探測(cè)時(shí)間。

（5）情形5。氣路從標(biāo)氣切換至樣氣。

①判斷多口閥位置：mpv＿0＞1、mpv＿1＝＝1；

②某標(biāo)氣測(cè)定時(shí)段內(nèi)（通常為5 min時(shí)長(zhǎng)，每5 s 1個(gè)、共計(jì)60個(gè)實(shí)測(cè)值），標(biāo)氣測(cè)定值序列為｛c1、c2、……c59、c60｝。取后2 min的觀測(cè)序列，作為1.1節(jié)情形1中提出的標(biāo)氣原始觀測(cè)樣本數(shù)組co2＿mpv［］、ch4＿mpv［］，依次計(jì)算1 min時(shí)長(zhǎng)的數(shù)據(jù)序列｛c49、c50…c59、c60｝、｛c48、c49…c58、c59｝、……｛c38、c39…c48、c49｝、｛c37、c38…c47、c48｝的標(biāo)準(zhǔn)差σ；

③對(duì)CO2（CH4），選取σ最小、且小于0.1 ppmv（或2 ppbv）的1 min數(shù)據(jù)序列，其平均值為該時(shí)段的標(biāo)氣實(shí)測(cè)平均值。關(guān)于標(biāo)氣采集樣本的選取原則，不在于數(shù)量多，而在于質(zhì)量精，可保證最終所選的標(biāo)氣樣本集合的平衡態(tài)最佳［11］；

④若有效測(cè)定序列為｛ck－11、ck－10……ck－1、ck｝，則選取ck的測(cè)定時(shí)間作為標(biāo)氣實(shí)測(cè)平均值的對(duì)應(yīng)時(shí)間。

3）計(jì)算程序與生成文件。

利用Data＿5min＿ave.java來實(shí)現(xiàn)樣氣5 min平均值和標(biāo)氣實(shí)測(cè)平均值的計(jì)算，并形成“STA＿5min＿ave.txt”文件（表1）。；STA為站名，下同。

表1 STA＿5min＿ave.txt1）文本格式示例

1.2 標(biāo)校參數(shù)的獲?。ɑ诨旌蠘?biāo)氣標(biāo)校法）

1）標(biāo)校參數(shù)計(jì)算。

設(shè)溫室氣體濃度的實(shí)測(cè)值為X、標(biāo)稱值為Y，對(duì)符合線性條件的設(shè)備［12］，X、Y符合最小二乘法原則，即X、Y滿足線性回歸方程Y＝aX＋b（a、b為標(biāo)校參數(shù)）。以珠海為例，分析儀觀測(cè)一輪樣氣、標(biāo)氣的時(shí)長(zhǎng)分別為2 h和5 min，順序依次為：WH、樣氣、T、樣氣、WL、樣氣、WH……。WH、WL用于計(jì)算標(biāo)校參數(shù)；T作為目標(biāo)氣，用于檢驗(yàn)標(biāo)校參數(shù)。

舉例假設(shè)分析儀的觀測(cè)時(shí)間和要素如下：07：55—08：00觀測(cè)WH、08：00—10：00觀測(cè)樣氣、10：00—10：05觀測(cè)T、10：05—12：05觀測(cè)樣氣、12：05—12：10觀測(cè)WL……。聯(lián)立方程組：

其中，WH、WL為標(biāo)氣實(shí)測(cè)平均值；WH1、WL1為標(biāo)稱值；解得標(biāo)校參數(shù)a、b，用于標(biāo)校08：00—12：05的樣氣5 min平均值。

2）標(biāo)校參數(shù)檢驗(yàn)。

分析儀在連續(xù)觀測(cè)期間會(huì)發(fā)生微小的漂移，因此，不同時(shí)段內(nèi)的標(biāo)校參數(shù)不同。1.2節(jié)例子中，令目標(biāo)氣計(jì)算誤差k＝｜（aT＋b）－T1｜（T為目標(biāo)氣實(shí)測(cè)平均值，T1為目標(biāo)氣標(biāo)稱值）。如果對(duì)于CO2（CH4），k＜0.1 ppmv（或2 ppbv），則a、b值是有效的［10］；否則，應(yīng)采用上一個(gè)樣氣觀測(cè)時(shí)段的a、b值進(jìn)行替換。

3）計(jì)算程序與生成文件。

利用Data＿calibration.java計(jì)算標(biāo)校參數(shù)，并生成名為“STA＿a＿b.xlxs”的標(biāo)校參數(shù)文件。

（1）標(biāo)氣實(shí)測(cè)平均值：從“STA＿5 min＿ave.txt”中讀取的標(biāo)氣實(shí)測(cè)平均值；

（2）標(biāo)氣標(biāo)稱值：對(duì)應(yīng)時(shí)間點(diǎn)的標(biāo)氣標(biāo)稱值；

（3）差值：標(biāo)氣實(shí)測(cè)平均值與標(biāo)稱值的差；

（4）目標(biāo)氣計(jì)算誤差：將目標(biāo)氣實(shí)測(cè)平均值T、標(biāo)稱值T1代入公式｜（a×T＋b）－T1｜的結(jié)果；

（5）a＿CO2、b＿CO2、a＿CH4、b＿CH4：通過公式（1）和（2）計(jì)算所得、且經(jīng)過檢驗(yàn)或替換后的標(biāo)校參數(shù)。

標(biāo)校參數(shù)與時(shí)間點(diǎn)的對(duì)應(yīng)關(guān)系如下：

（1）t1、t3為相鄰的WH、WL標(biāo)氣觀測(cè)時(shí)間點(diǎn)，根據(jù)t1、t3的標(biāo)氣實(shí)測(cè)平均值、標(biāo)稱值計(jì)算出時(shí)間t1至t3的標(biāo)校參數(shù)a＿CO2、b＿CO2；之后再用t1、t3之間的目標(biāo)氣T觀測(cè)時(shí)間點(diǎn)（即t2）的實(shí)測(cè)平均值、標(biāo)稱值來驗(yàn)證a＿CO2、b＿CO2。

（2）時(shí)間t3至t4沒有目標(biāo)氣T觀測(cè)，則選取t4之后、距離t4最近的一次目標(biāo)氣T觀測(cè)時(shí)間點(diǎn)（即t5）的實(shí)測(cè)平均值、標(biāo)稱值來驗(yàn)證時(shí)間t3至t4的標(biāo)校參數(shù)。

（3）t1對(duì)應(yīng)行的a＿CO2、b＿CO2是時(shí)間t1至t3的標(biāo)校參數(shù)，t3對(duì)應(yīng)行的a＿CO2、b＿CO2是時(shí)間t3至t4的標(biāo)校參數(shù)，依此類推。

4）算法的可靠性評(píng)估。

利用混合標(biāo)氣標(biāo)校溫室氣體數(shù)據(jù)的計(jì)算工具，計(jì)算珠海市某溫室氣體監(jiān)測(cè)站2014年全年的CO2、CH4目標(biāo)氣計(jì)算誤差k。該監(jiān)測(cè)站在2013年12月底正式運(yùn)行，2014年業(yè)務(wù)運(yùn)行期間滿足以下條件：（1）新投入使用的分析儀性能穩(wěn)定、漂移較小；（2）WH、T、WL等3瓶標(biāo)氣均由中國(guó)氣象科學(xué)研究院溫室氣體及相關(guān)微量成分實(shí)驗(yàn)室標(biāo)定，定值結(jié)果可溯源至世界氣象組織CO2、CH4中心標(biāo)校實(shí)驗(yàn)室維持的國(guó)際一級(jí)標(biāo)氣，期間標(biāo)氣質(zhì)量濃度較穩(wěn)定且不更換。2014年，該站的絕大部分CO2目標(biāo)氣計(jì)算誤差k＜0.1 ppmv，所有CH4目標(biāo)氣計(jì)算誤差k＜2 ppbv（表2），這說明編譯出的混合標(biāo)氣標(biāo)校溫室氣體數(shù)據(jù)的計(jì)算方法，其計(jì)算精度滿足WMO的要求。

表2 珠海市某溫室氣體監(jiān)測(cè)站2014年CO2、CH4目標(biāo)氣計(jì)算誤差k值分布

1.3 樣氣5 min平均值的標(biāo)校

利用Data＿output.java計(jì)算標(biāo)校值，將“STA＿5min＿ave.txt”文件的樣氣5 min平均值，代入相應(yīng)時(shí)段的標(biāo)校公式Y(jié)＝aX＋b，可得出標(biāo)校值，生成“STA＿output.txt”文件。

2 算法集成與運(yùn)行

運(yùn)行算法集成模塊QC＿manual.java（圖3），一鍵執(zhí)行Data＿5 min＿ave.java、Data＿calibration.java、Data＿output.java，完成多年歷史溫室氣體原始數(shù)據(jù)的質(zhì)量控制和標(biāo)校訂正。該計(jì)算機(jī)工具使溫室氣體數(shù)據(jù)標(biāo)校更方便快捷，根據(jù)標(biāo)校法計(jì)算工具的測(cè)試情況，完成一個(gè)站點(diǎn)7年原始5 s數(shù)據(jù)（超過4 000萬行）的質(zhì)量控制、標(biāo)校訂正，并生成樣氣5 min平均值、標(biāo)校參數(shù)、樣氣標(biāo)校值等數(shù)據(jù)文件，時(shí)長(zhǎng)約為8 min。

圖3 算法集成模塊QC＿manual.java（手動(dòng)質(zhì)控）界面

該工具還提供更多排查系統(tǒng)故障的途徑，以標(biāo)校參數(shù)文件“STA＿a＿b.xlxs”為例。當(dāng)標(biāo)氣實(shí)測(cè)平均值與標(biāo)氣標(biāo)稱值差異過大時(shí)，提示標(biāo)氣質(zhì)量濃度發(fā)生改變、標(biāo)氣更換記錄未輸入電腦等可能（圖4），可為維護(hù)人員更換標(biāo)氣、完善標(biāo)氣記錄提供參考依據(jù)。

圖4 標(biāo)氣實(shí)測(cè)平均值與標(biāo)稱值差異過大（紅色方框處）界面（提示可能未及時(shí)錄入標(biāo)氣更換記錄）

當(dāng)標(biāo)氣實(shí)測(cè)平均值為0時(shí)，提示可能出現(xiàn)以下情況：（1）標(biāo)氣質(zhì)量濃度穩(wěn)定，但分析儀出現(xiàn)震蕩性漂移，從而導(dǎo)致對(duì)標(biāo)氣的觀測(cè)達(dá)不到規(guī)定標(biāo)準(zhǔn)；如果分析儀對(duì)所有標(biāo)氣的測(cè)量都無法達(dá)到平衡態(tài)，則應(yīng)對(duì)分析儀進(jìn)行人工校準(zhǔn)甚至檢查激光發(fā)生器、光腔等關(guān)鍵組件是否故障，以提高分析儀運(yùn)行效率和觀測(cè)數(shù)據(jù)精度。（2）分析儀不存在漂移，但是標(biāo)氣質(zhì)量濃度出現(xiàn)震蕩性變化；若分析儀只是對(duì)某種標(biāo)氣的測(cè)量無法達(dá)到平衡態(tài)，則需要對(duì)相應(yīng)的標(biāo)氣進(jìn)行重新標(biāo)定或更換標(biāo)氣。

3 結(jié)論

混合標(biāo)氣標(biāo)校溫室氣體數(shù)據(jù)的計(jì)算機(jī)工具，逐一實(shí)現(xiàn)了溫室氣體原始數(shù)據(jù)質(zhì)量控制、標(biāo)定參數(shù)計(jì)算、標(biāo)校值計(jì)算等流程，目前已應(yīng)用于珠海市溫室氣體監(jiān)測(cè)站網(wǎng)的數(shù)據(jù)標(biāo)校。

珠海市氣象工作者全程負(fù)責(zé)編譯該計(jì)算工具，在掌握計(jì)算機(jī)源代碼與核心編譯技術(shù)的基礎(chǔ)上，能自行承擔(dān)計(jì)算工具的日常維護(hù)、代碼修改和后期功能擴(kuò)展，可提升基層監(jiān)測(cè)站關(guān)于數(shù)據(jù)應(yīng)用的自主性；該計(jì)算工具基于Java編程語(yǔ)言，架構(gòu)簡(jiǎn)潔，可移植性強(qiáng)，無需新建數(shù)據(jù)庫(kù)，可推動(dòng)各基層溫室氣體監(jiān)測(cè)站之間的業(yè)務(wù)交流。

此外，該計(jì)算工具對(duì)于觀測(cè)主設(shè)備符合線性條件的大氣成分（如N2O、CO、NH3以及含氟溫室氣體等）觀測(cè)系統(tǒng)同樣適用，用戶只需針對(duì)不同的原始文件格式，更改Java代碼中的Substring函數(shù)參數(shù)，重新抓取各要素原始數(shù)據(jù)，并結(jié)合準(zhǔn)確的標(biāo)氣標(biāo)稱值，便可進(jìn)行數(shù)據(jù)標(biāo)校。