(part of the AQMEII-NA_N2O family of projects)

table of contents



In EPIC, the emittivity per gridcell per timestep per crop is the amount of N2O the crop would produce over the given timestep in the given gridcell if the full area of the gridcell was devoted to that crop over the full timestep.


Uses bash to drive NCL and R code to

  1. "reunit" an EDGAR-4.2 N2O emission inventory (in netCDF format) for agricultural soils from (per-gridcell) flux rate (kg/m^2/s) to molar-mass rate (mol/s, which CMAQ wants anyway) prior to regridding.
    1. ... and visualize it.
  2. 2D-regrid the EDGAR inventory from global/unprojected to a projected subdomain (AQMEII-NA).
    1. ... and visualize it.
  3. Process BELD crop-coverage data needed to process EPIC data.
    1. Read a BELD CSV file with structure {row=GRIDID, col=crop}, where
      • each crop maps one-to-one to a layer in the EPIC netCDF
      • each GRIDID maps one-to-one to an individual gridcell in the EPIC netCDF
      • each (GRIDID,crop) tuple is the percentage (0 .. 100) of the gridcell area devoted to that crop (i.e., the crop's coverage of that gridcell).
    2. Write output (with format=, e.g., RDS) with structure {col=AQMEII gridcell col, row=AQMEII gridcell row, layer=EPIC crop layer} such that
      • the index of each (col,row) tuple matches the corresponding gridcell in the EPIC data.
      • the index of each layer matches the corresponding crop layer in the EPIC data.
      • each (col,row,layer) tuple is the coverage of that gridcell=(col,row) by that crop=layer as a proportion (0 .. 1).
  4. Obtain EPIC N2O emittivities.
    1. Read per-crop N2O emittivity from the raw EPIC netCDF. Unfortunately the latter is too large to provide in this project, since it contains many species other than N2O, so, first, we remove the data variables for the other species.
    2. "Demonotonicize" the N2O emittivities, if necessary. Consider the 2D spatial grid to be "flattened" into a 1D vector, and to be "desparsified" such that all missing values are removed. The current dataset has the unfortunate problem that, for each crop/layer, each estimation in the gridcell/vector was summed with all its predecessors, instead of being recorded separately. E.g., if the true estimations were a[i]==u, a[j]==v, and a[k]==w (for vector indices i < j < k), the data values given are a[i]==u, a[j]==u+v, and a[k]==u+v+w.
    3. "Reunit" the emittivities from a mass flux (explicitly kg/ha, implicitly kg/ha/yr) to a molar-mass rate (mol/s) consumable by CMAQ.
  5. (completed) Obtain EPIC N2O emissions:
    1. Obtain the crop N2O-emittivity vector from each gridcell's EPIC data: a vector of the N2O emittivities of all the crops over that gridcell during that timestep.
    2. Obtain the crop coverage vector from each gridcell's BELD data: a vector of the coverages of all the crops over that gridcell during that timestep.
    3. Write to output, for each gridcell, the scalar of EPIC N2O emissions for that gridcell (i.e., the scalar product of the crop N2O-emittivity vector and the crop coverage vector). I.e., for output dimensions=(TSTEP, LAY, ROW, COL),
      1. TSTEP=0 (annual value)
      2. LAY=0 (corresponding to a surface/vertical layer--semantics of 'layer' revert back to "the usual")
      3. (ROW, COL) is the AQMEII-NA grid
      4. value(TSTEP, LAY, ROW, COL) is EPIC's estimate of the total emission of N2O (in mol/s) from soils planted in all crops modeled by EPIC in that gridcell over that timestep.
    4. Visualize the EPIC emissions.
  6. Combine EDGAR and EPIC emissions over AQMEII-NA.
    1. (algorithm may change) For each gridcell in AQMEII-NA: if the EPIC emissions exceed the EDGAR emissions, use the EPIC value.
    2. Visualize the combined AQMEII-NA emissions.
  7. Create CMAQ-style emissions files (e.g., this, when gunziped) containing emissions for every hour in 2008 (the model year).
  8. Check extent to which mass is conserved from global input to regional output.

Currently does not provide a clean or general-purpose (much less packaged) solution! but merely shows how to do these tasks using

  • bash (tested with version=3.2.25)
  • NCL (tested with version=6.1.2)
  • R (tested with version=3.0.0) and packages including

Visualizations of datasets are provided

  • statistically (summary statistics of datasets are generated and presented in console)
  • by plots (saved as PDF, displayed by launching the user-configured PDF viewer)
  • by VERDI (which also serves as a check on the IOAPI-compatibility of the generated netCDF)


To run this code,

  1. git clone this repo. (For more details, see Cloning a git repository.)
  2. cd to its working directory (where you cloned it to).
  3. Setup your applications and paths.
    1. Clone regrid_utils to a subdirectory of the working directory, copy regrid_utils/bash_utilities.sh up to the working directory, and open ./bash_utilities.sh in an editor! You will probably need to edit its functions setup_paths and setup_apps to make it work on your platform. Notably you will want to point it to your PDF viewer and NCL and R executables.
    2. Open the driver (bash) script ./AQMEII_driver.sh in an editor and ensure paths to your inputs are (or will be) correct.
  4. Download inputs that cannot currently be automated by this project, saving them to the paths you denoted in the driver. Most required data can be downloaded by the driver (provided the URIs pointing to the data remain unchanged), but some cannot as of the time of this writing:
    1. CONUS EPIC emittivities. These unfortunately are currently much too large to save on bitbucket; we hope to find a home for them online before too long.
  5. Run the driver: $ ./AQMEII_driver.sh This will download input, then run ...
    1. reunit_vis_EDGAR_global.r, which will prepare the raw EDGAR input (unfortunately not currently directly downloadable from EDGAR) for regridding, and plot it.
    2. regrid_vis.r, which will 2D-regrid the reunit-ed global EDGAR (from the previous step) to AQMEII-NA, and plot it.
    3. BELD_CSV_to_RDS.r, which will read AQMEII-NA BELD coverage data and convert it to RDS for later use.
    4. distill_EPIC_emittivities.ncl, which will read the previously-saved CONUS EPIC emittivities and create more consumable output by
      • removing species and processes other than N2O production
      • demonotonicizing the N2O emittivities, if necessary (see above)
      • reunit-ing the emittivities (see above)
    5. compute_EPIC_emissions.r, which calculates the scalar product of the BELD crop coverages and the EPIC emittivities, creating and visualizing CONUS-only, EPIC-only N2O emissions. Note: there appears to be a spurious maximum in the current version of EPIC, which is removed: see comments in the code, as well as visualization statistical output, for details.
    6. combine_EDGAR_and_EPIC_emissions.r, which (predicably) combines the EPIC emissions from the previous step with the reunit-ed, regridded EDGAR from step 2 into one dataset, and visualizes it.
    7. make_hourlies.ncl, which outputs the combined emissions in CMAQ format: i.e., 25 hourly (and in this case identical, since we have no temporality) grids.
    8. check_conservation.ncl to investigate conservation of mass from inputs to outputs. Current output includes

          EDGAR global emissions (mol/s) from ./v42_N2O_2008_IPCC_4C_4D.0.1x0.1_reunit.nc
          |obs|  =6.48e+06
          min    =0
          q1     =0
          med    =0
          mean   =5.177e-04
          q3     =0
          max    =1.003e+00
          sum    =3.352e+03
          EDGAR emissions (mol/s) over AQMEII-NA from ./v42_N2O_2008_IPCC_4C_4D.0.1x0.1_reunit_regrid.nc
          |obs|  =137241
          min    =0
          q1     =0
          med    =1.575e-04
          mean   =1.890e-03
          q3     =3.121e-03
          max    =8.609e-02
          sum    =2.593e+02
          EPIC emissions (mol/s) over CONUS from ./5yravg_20111219_pure_emis.nc
          |obs|  =40763
          min    =0
          q1     =1.021e-05
          med    =1.714e-04
          mean   =1.373e-03
          q3     =1.112e-03
          max    =7.571e-02
          sum    =5.596e+01
          EDGAR+EPIC emissions (mol/s) over AQMEII-NA from ./emis_mole_N2O_2008_12US1_cmaq_cb05_soa_2008ab_08c.nc
          |obs|  =137241
          min    =0
          q1     =0
          med    =1.600e-04
          mean   =1.986e-03
          q3     =3.198e-03
          max    =8.609e-02
          sum    =2.725e+02
          Comparison of EDGAR emissions
              EDGAR     global      EDGAR  AQMEII-NA  AQMEII-NA/
             global        NAs  AQMEII-NA        NAs     global
           3.35e+03          0   2.59e+02          0   7.74e-02
          EPIC vs EDGAR emissions over AQMEII-NA
               EPIC       EPIC      EDGAR      EDGAR      EPIC/
            (CONUS)        NAs  AQMEII-NA        NAs      EDGAR
           5.60e+01      96478   2.59e+02          0   2.16e-01
          EDGAR vs EPIC+EDGAR emissions (both over AQMEII-NA)
              EDGAR      EDGAR   combined   combined     EDGAR/
          AQMEII-NA        NAs  AQMEII-NA        NAs   combined
           2.59e+02          0   2.72e+02          0   9.52e-01
          EPIC+EDGAR AQMEII-NA vs EDGAR global emissions
           combined   combined      EDGAR      EDGAR  combined/
          AQMEII-NA        NAs     global        NAs      EDGAR
           2.72e+02          0   3.35e+03          0   8.13e-02


  1. Retest with newest regrid_utils! Currently, repo_diff.sh shows the following local diffs:
    • IOAPI.ncl
    • string.ncl
    • summarize.ncl
  2. Move all these TODOs to issue tracker.
  3. *.sh: use bash booleans à la N2O_integration_driver.sh.
  4. Handle {name, path to} VERDI in bash_utilities.sh.
  5. Generate plot of EPIC-EDGAR (i.e., show where EDGAR < EPIC, EDGAR > EPIC).
  6. Add to README, wiki: links to overviews of BELD, EPIC.
  7. Add to wiki: Simulation-of-N2O-Production-and-Transport-in-the-US-Cornbelt-Compared-to-Tower-Measurements#wiki-input-processing-EPIC-0509 (I have input processing instructions for other EIs, but not EPIC).
  8. Move explicit constants BELD_CSV_to_RDS.r -> AQMEII_driver.sh.
  9. Create common project for regrid_resources à la regrid_utils, so I don't hafta hunt down which resource is in which project.
  10. (this project only?) Use the EPIC input as template.
  11. All regrids: how to nudge off/onshore as required? e.g., soil or burning emissions should never be offshore, marine emissions should never be onshore.
  12. All regrid maps: add Caribbean islands (esp Bahamas! for offshore burning), Canadian provinces, Mexican states.
  13. Complain to ncl-talk about NCL "unsupported extensions," e.g., .ncf and <null/> (e.g., MCIP output).
  14. Determine why <- assignment is occasionally required in calls to visualize.*(...).
  15. Fully document platform versions (e.g., linux, compilers, bash, NCL, R).
  16. Test on
    • tlrPanP5 (which now has R package=ncdf4, but readAsciiTable of input .txt's is very slow compared to terrae)
    • HPCC (once problem with ncdf4 on amad1 is debugged: in process with JOB and KMF)