Link Search Menu Expand Document
Focus on Air Quality

地面排放檔案之轉換

Table of contents

背景

  • CMAQ模式的地面排放量為一逐日拆解的檔案（每檔案需25個小時）。其開始時間為UTC 0時（LST -8）。
  • 可以按照排放類別拆分或合併。
  • 空間網格系統之定義乃按照mcip之定義方式。
  • 由CAMx模式排放量檔案而來，其產生過程詳如下表：
範圍格式類別
東亞範圍nc面排放源含生物源、
  航空排放
  工業排放
  陸上交通排放源
  船舶排放
臺灣範圍nc面排放源
  植物源
  陸上交通排放源
  船舶排放
  地面點源排放

污染物項目之對照camx2cmaq.py

執行方式

  • 引數：0～4之序號，依序為area biog line ptse ship
  • I/O檔案
    • 輸入：home+'/'+p+'/'+'fortBE.413_teds10.'+ext[p]+mm+'.nc'
    • 模版：template.nc
    • 結果：p+'_TWN_3X3.16'+mm+'.nc'

重要轉換

  • 時間
    • CAMx模式為LST、CMAQ為UTC：需要減8
    • Datetime to/from Julian day：dtconvertor
  • 排放量
    • CAMx為g/hour、CMAQ為g/s：需要除3600
  • 污染項目
    • 在$EMISSCTRL_NML檔案中定義。此處不作轉換

camx2cmaq.py程式碼

kuang@master /nas1/cmaqruns/2019base/data/emis/TEDS
$ cat camx2cmaq.py
import os,sys
import numpy as np
import netCDF4
import datetime
from dtconvertor import dt2jul, jul2dt

ipth=int(sys.argv[1])
pth='area biog line ptse ship'.split()
ext={i:i for i in pth}
ext.update({'ship':'51Ab','ptse':'ptsG'})
home='/nas1/TEDS/teds11'
for p in pth[ipth:ipth+1]:
  for m in range(1,13):
    mm='{:02d}'.format(m)
    fname=home+'/'+p+'/'+'fortBE.413_teds10.'+ext[p]+mm+'.nc'
    fnameO=p+'_TWN_3X3.16'+mm+'.nc'
    os.system('cp template.nc '+fnameO)
    nc0= netCDF4.Dataset(fname, 'r')
    nc = netCDF4.Dataset(fnameO, 'r+')
    V0=[list(filter(lambda x:nc0.variables[x].ndim==j, [i for i in nc0.variables])) for j in [1,2,3,4]]
    V=[list(filter(lambda x:nc.variables[x].ndim==j, [i for i in nc.variables])) for j in [1,2,3,4]]
    nv=len(V[3])
    nt,nlay,nrow,ncol=nc0.variables[V0[3][0]].shape
    sdatetime=[jul2dt(nc0.variables['TFLAG'][t,0,:]) for t in range(nt)]
    jul2=[dt2jul(i+datetime.timedelta(hours=-8)) for i in sdatetime]
    for t in range(nt):
      nc.variables['TFLAG'][t,:,0]=[jul2[t][0] for i in range(nv)]
      nc.variables['TFLAG'][t,:,1]=[jul2[t][1] for i in range(nv)]
      for iv in range(nv):
        v=V[3][iv]
        if v in V0[3]:
          nc.variables[v][t,:,:,:]=nc0.variables[v][t,:,:,:]/3600.
        else:
          nc.variables[v][t,:,:,:]=0
    nc.SDATE, nc.STIME=(jul2[0][0],jul2[0][1])
    nc.close()
    print(fnameO)

按日拆分m3.nc檔案

for yy in 18;do 
  for mm in 0{1..9} {10..12};do
    for path in area biog line ptse ship;do 
      brk_day2.cs ${path}_TWN_3X3.$yy$mm.nc
    done
  done
done

Reference