Make CMAQ Emissions from a CSV File

Make a CMAQ emissions input for NOx from a CSV file with latitude, longitude, and hour-of-day coordinates. The code below will make a file from a GRIDDESC file. Then, add a variable for NO and NO2. Finally, populate the file and save it as an IOAPI-like file.

import PseudoNetCDF as pnc
import pandas as pd
import numpy as np

# create a dummy input file
lat = (np.sin(np.linspace(0, 2*np.pi, 100)) * 7 + 40).round(2)
lon = np.linspace(-120, -80, 100)
hour = np.floor(np.linspace(0, 23.99, 100)).astype('i')
NO = np.random.lognormal(0.09, sigma=0.1, size=100)
NO2 = np.random.lognormal(0.01, sigma=0.01, size=100)
data = pd.DataFrame.from_dict(dict(
    lat=lat, lon=lon, hour=hour, NO=NO, NO2=NO2
)).to_csv('emis.csv', index=False)
# create a GRIDDESC
with open('GRIDDESC', 'w') as gf:
    gf.write(
       "' '\n'LamCon_40N_97W'\n 2 33.000 45.000 -97.000 -97.000 40.000\n" +
       "' '\n'12US1'\n'LamCon_40N_97W' " +
       "-2556000.0 -1728000.0 12000.0 12000.0 459 299 1\n' '"
    )

# Create a Template
gf = pnc.pncopen('GRIDDESC', format='griddesc')
gf.SDATE = 2016001
gf.TSTEP = 10000
gf.updatetflag(overwrite=True)
tmpf = gf.sliceDimensions(TSTEP=[0]*25)

#  Read in the emissions
data = pd.read_csv('emis.csv')

# Add I,J and time coordinate
data['I'], data['J'] = tmpf.ll2ij(data.lon.values, data.lat.values)
cellhourtotal = data.groupby(['I', 'J', 'hour'], as_index=False).sum()
# Track emission keys
metakeys = ('lat', 'lon', 'hour')
emiskeys = [k for k in data.columns if k not in metakeys]

h = cellhourtotal.hour
i = cellhourtotal.I
j = cellhourtotal.J

for emiskey in emiskeys:
    evar = tmpf.createVariable(emiskey, 'f', ('TSTEP', 'LAY', 'ROW', 'COL'))
    evar.setncatts(dict(
      units='moles/s', long_name=emiskey, var_desc=emiskey
    ))
    evar[h, h*0, j, i] = cellhourtotal[emiskey].values

# Get rid of initial DUMMY variable
del tmpf.variables['DUMMY']

# Update TFLAG to be consistent with variables
tmpf.updatetflag(tstep=10000, overwrite=True)

# Remove VAR-LIST so that it can be inferred
delattr(tmpf, 'VAR-LIST')
tmpf.updatemeta()

# Save out
tmpf.save('emis.nc', format='NETCDF3_CLASSIC')