CMAQ點源常數檔案之準備

Table of contents

緣起
- 點源nc檔案煙囪參數之版本差異
思路策略
程式設計
- 分段說明
- CAMx 6版與7版的差異
後處理
- REAS點源解讀結果(CMAQ格式)之比較

顧名思義，這支程式就是為了產生CMAQ所需的點源排放量中常數、不隨時間變化的部分。

因為是常數，沒有必要切割成（複製）成逐日檔。但至少在執行批次範圍內所有日期的變數檔可以共用此常數檔。簡言之，就是執行批次範圍內的煙囪個數、順序、必須保持一致，這使得很多未操作的污染源及煙道，必須仍然保留在排放檔中。

執行變數檔案轉換之前，必須先執行常數檔案之轉換，主要是REAS點源排放量的處理是在這一階段執行。

緣起

在CAMx(與傳統的UAM)中，點源檔案與其他檔案最大的不同就是有常數部分的檔頭，其他檔案可能只有4筆header，點源因有煙囪位置、排放條件等，這些不會隨著時間有所改變，因此放在檔頭。
- 即使CAMx 6/7等最新版nc格式，仍然採取此一策略。煙囪參數是1維變數，是沒有時間維度的。
CMAQ基本上也照此劃分原則。同時為要保持檔案可以符合ioapi的協定，CMAQ乃將常數項、變數項予以分開成2個檔案，檔案結構較為單純。

點源nc檔案煙囪參數之版本差異

項目	CAMx6	CAMx7	CMAQ tim_const	說明
煙囪參數	`[v+'STK' for v in 'XYHDTV']`	`xcoor, ycoor, stkheight, stkdiam, stktemp, stkspeed`	`XLOCA, YLOCA, STKHT, STKDM, STKTK, STKVE`	-
煙囪參數維度	`NSTK`	`COL`	`TSTEP, LAY, ROW, COL`	CAMx皆為1維。CMAQ為符合ioapi協定，NTSTEP, NLAYS, NCOLS皆為1
特殊參數	`(FLOW, PLMHT, IONE, ITWO, KCELL)`¹, `(NSTK)`²	`pigflag, saoverride`	`LATITUDE, LONGITUDE, STKFLW, COL, ROW, ISTACK, IFIP, LMAJOR, LPING, STKCNT`	CAMx6前2項與CMAQ前3項為實數，其餘皆整數
特殊參數維度	¹`(TSTEP, NSTK)`, ²`(TSTEP)`	`COL`	`TSTEP, LAY, ROW, COL`	CMAQ符合ioapi協定
模擬期間煙囪個數、順序一致性	不需要	必須	必須	CAMx6也可以使用隨時間改變的氣流量與煙流高

在新版的CMAQ中，也取消了隨時間變化的煙囪上升高度項，另外在模式內計算。這一點算是引進CAMx的優勢強項。
CMAQ/CAMx排放量檔案的演進、格式與單位的差異、轉換、等等，可以參考CMAQ/CAMx排放量檔案之轉換。

思路策略

如前述此二模式檔案有此背景，因此在準備時就不需猶豫該採取什麼策略，先將原始數據先整理成CAMx檔案的格式，再將其發展程常數及變數2個檔案，就會是最精簡、省工的方向。考慮因素：

CAMx檔案同時有常數及變數項，不會出現檔案配對的困擾(CMAQ必須在run_cctm.csh中才能確認其關係)。
就常數項而言，CAMx檔案會在全月範圍內保持一致，這對該月內所有的批次而言，都可以適用、不必擔心煙囪筆數順序的錯誤。
REAS、TEDS如果直接轉成這個常數檔就不容易處理了，因為每一天CEMS展開之後，煙囪的筆數與順序應該是都不一樣，還是要定義一段較長的時間來匯總（聯集）所有煙囪，這就是CAMx檔案的邏輯。

程式設計

分段說明

詳見pt_const程式說明的內容。

CAMx 6版與7版的差異

既然是從CAMx點源檔案開始讀起，CAMx改版就會對本程式造成衝擊。

早期ptsource格式，必須以pnc來讀取，6版以後可以接受nc格式，但內容還是有差異，見CMAQ/CAMx排放量檔案之轉換詳細說明。

  kuang@114-32-164-198 /Users/cmaqruns/2019base/data/ptse/REAS
$ diff /Users/cmaqruns/2016base/data/ptse/pt_const.py pt_const.py
3d2
< import PseudoNetCDF as pnc
9,10c8,9
<
< mon=int(sys.argv[1][-2:])
---
> #arvg=../twn/fortBE.413_teds11.ptse02.nc
> mon=int(sys.argv[1][-5:-3])
12c11
< pth='.' #~/mac/cmaqruns/2016base/data/ptse'
---
> pth='/nas1/cmaqruns/2016base/data/ptse'
56c55
< df.set_index('lon').to_csv('point_reas16'+sys.argv[1][-2:]+'.csv')
---
> df.set_index('lon').to_csv('point_reas16'+sys.argv[1][-5:-3]+'.csv')
59c58
< pt=pnc.pncopen(fname1,format='point_source')
---
> pt = netCDF4.Dataset(fname1,'r')
63c62
< nhr,nvar,dt=pt.variables[v3[0]].shape
---
> nhr,nvar,dt=pt.variables['TFLAG'].shape
66,67c65,66
< tb=pt.STIME[0]-8 #UTC
< fname1=fname1.replace('fortBE.14.','').replace('base','16')
---
> tb=pt.STIME - 80000
> fname1=fname1.split('/')[-1].replace('fortBE.413_','').replace('ptse','19').replace('.nc','')
78,79c77,78
< x=list(pt.variables['XSTK'][:nopts])+list(df.x_m)
< y=list(pt.variables['YSTK'][:nopts])+list(df.y_m)
---
> x=list(pt.variables['xcoord'][:nopts])+list(df.x_m)
> y=list(pt.variables['ycoord'][:nopts])+list(df.y_m)
95c94
< nc.SDATE=2000000+pt.SDATE[0]
---
> nc.SDATE=2000000+pt.SDATE
97c96
< mp={'STKDM':'DSTK','STKHT':'HSTK','STKTK':'TSTK','STKVE':'VSTK','XLOCA':'XSTK', 'YLOCA':'YSTK',}
---
> mp={'STKDM':'stkdiam','STKHT':'stkheight','STKTK':'stktemp','STKVE':'stkspeed','XLOCA':'xcoord', 'YLOCA':'ycoord',}
113c112
<     HSTK=pt.variables['HSTK'][i]
---
>     HSTK=pt.variables['stkheight'][i]

後處理

CMAQ點源檔案無法檢視單一檢視，需有常數及變數2個檔案同時進行。
為此設計pt2em_d01.py，將離散點之點源排放量，歸併在面源網格中，以便以VERDI進行檢視，如附圖所示。

REAS點源解讀結果(CMAQ格式)之比較


3月結果	1月結果

Resources

Barron Henderson, pseudonetcdf tutorial, http://www.barronh.com/software/tutorials/pseudonetcdf-tutorial
verdi usage https://www.airqualitymodeling.org/index.php/VERDI_1.5_User_Manual#3.1_Installation_Instructions_for_Linux_and_Mac
VERDI使用說明 : http://www.evernote.com/l/AH3leuVQTuBEF7Vrs0D1C8Q-Iff5CpHl7eU
pt2emNest:https://github.com/sinotec2/CAMx_utility/blob/master/pt2emNest.f
ncks