How to regridding Netcdf emission to model domain without use SMOKE?

SMOKE Emissions Inventory
1

Hi all,
i have trouble when i regridding MEIC emission to CMAQ model domain(i use the fortran language) ,it can’t directly used by the SMOKE model，i used the bilinear interpolation method， but the results is low. is the IOAPI tools or other methods can do this? Thank you!!

2

the MEIC emission inventory is 0.25°，my CMAQ domain is 36km.

3

AFAIK, you’re probably going to have to “roll your own code” a bit. As you do so, there are several issues:

  • What type of files does MEIC have? Are they “raw” netCDF? If so, the I/O API module MODULE MODNCFIO has high-level routines for reading netCDF that may be useful.
  • For the grid-to-grid transforms, MODULE MODGCTP has lots of map-transformation utilities, including particularly routines GRID2INDX, PNTS2INDX, and INDXMULT that give you a highly-optimized, OpenMP-parallel bilinear interpolation package.
  • Be wary of the emissions units: they usually have an implicit “per grid cell” that isn’t mentioned, so you’ll have to re-scale by the ratio of grid-cell areas – cell-areas are constant for Lambert grids, but variable for Lat-Lon grids; where LAT is the latitude at the cell-centers, a good approximation for these areas is

DLAT x DLON x COS(LAT) x R_Earth^2

4

The regridding and speciation could be done similar to what we just did for MIX. Of course, but the details depend on the format of the raw data which we don’t have

MIX is a new anthropogenic emission inventory for Asia for the years 2008 and 2010. It is … a mosaic of up-to-date regional and national emission inventories including REAS inventory version 2.1 …, the Multi-resolution Emission Inventory for China (MEIC) …, a high-resolution NH3 emission inventory …, an Indian emission inventory …, and the official Korean emission inventory…

http://meicmodel.org/?page_id=87&lang=en

See this file from Barron Henderson: https://gist.github.com/barronh/97633470997ca77b4bb949df20d9a0a3

1 Like
5

Thank you for above professors:
The code partlity in the below:
Take MIX emission inventory(MICS_Asia_CO_2010_0.25x0.25.nc) as example
(it equal to MEIC emission inventory,both in NETCDF formats)
the resolution is 0.25°,the unit is Mg/month/grid.the (lon,lat,month) is (560,441,12), my question
is how to regridding MIX_3D(560,441,12) to my model domain REGRID_MIX_3D(65,99,12)

! MIX emission inventory domain
INTEGER,PARAMETER :: MIXLON=560,MIXLAT=441,MIXMONTH=12
! my model dpmain
INTEGER,PARAMETER :: NCOLS=65,NROWS=99,NLAYS=1
REAL :: MIX_3D(MIXLON,MIXLAT,MIXMONTH)
!!!
REAL :: REGRID_MIX_3D(NCOLS,NROWS,MIXMONTH)
REAL :: CO_3D(NCOLS,NROWS,NLAYS)

PRINT *, ’ — Reading MIX emissions from NetCDF files —’
i will use the professor mentioned method ‘MODULE MODNCFIO’ have a try.
FILENAME=‘MICS_Asia_CO_2010_0.25x0.25.nc’
STATUS = NF90_OPEN(TRIM(FILENAME), NF90_NOWRITE, NCID)
IF(STATUS /= NF90_NOERR) CALL HANDLE_ERR(STATUS)

STATUS = NF90_INQ_VARID(NCID, “CO_RESIDENTIAL”, VAR_ID)
STATUS = NF90_GET_VAR(NCID, VAR_ID, MIX_3D)
IF(STATUS /= NF90_NOERR) CALL HANDLE_ERR(STATUS)

!!!
The key questions is here:
i want to regridding MIX_3D into REGRID_MIX_3D, is the IOAPI tools
has the magic functions to do this? i will have a try in the above the professor
mentioned method in IOAPI

LOGDEV = INIT3()
WRITE(,) " THE INIT3 NUMBER IS ", LOGDEV
…

IOAPI_FILE=‘CO.nc’
IF ( .NOT. OPEN3( IOAPI_FILE, FSNEW3, ‘OUTPUT’ ) ) THEN
CALL M3EXIT( ‘OUTPUT’, 0, 0, ‘COULD NOT OPEN IOAPI_FILE’, 2 )
END IF
WRITE(,) " SUCCESS OPENED IOAPI_FILE!"

MXREC3D=12
DO T0=1,MXREC3D

DO R=1,NROWS
DO C=1,NCOLS
!!!
The question may be also in this:
i did not resclae the ratio of grid-cell areas.
ETEMP=REGRID_MIX_3D(C,R,2)*1000000/24./3600./28. ! only test in here
CO_3D(C,R,1)=C0_3D(C,R,1) + ETEMP
END DO
END DO
IF ( .NOT. WRITE3 (IOAPI_FILE,‘CO’,JDATE,JTIME,CO_3D)) THEN
MESG = ‘ERROR WRITTING "’//VNAME//’" to IOAPI_FILE :’
CALL M3EXIT(IOAPI_FILE,0,0,MESG,2)
END IF
CALL NEXTIME( JDATE, JTIME, TSTEP3D )
END DO

END PROGRAM

SUBROUTINE HANDLE_ERR(STATUS)
USE NETCDF
IMPLICIT NONE

INTEGER,INTENT(IN)  :: STATUS

IF (STATUS /= NF90_NOERR) THEN
 WRITE(*,*) NF90_STRERROR(STATUS)
STOP "STOPPED"
END IF
END SUBROUTINE HANDLE_ERR
6

Thank you professors.

7

Here’s an outline of the algorithm (as best as I can present it given this Forum software):

  • Compute RATIO( MIXLON,MIXLAT) as area(output-cell)/area(lat-lon-cell) using the area-formula above (noting that R_Earth needs to be in meters).
  • Use `GRID2INDX` to compute `IX2, PX2, PY2`

  • For each time step: 

    • read `MIX` from the input-file and multiply it by `RATIO`.

    • Use `INDXMULT` to compute `REGRID` from `MIX x RATIO` and  `IX2, PX2, PY2`

    • Write REGRID to your output file
8

Thank you, professors, i will have a try.

9

Hi professors:
i have another questions, whrn i use GRID2INDX, i find the MIX emission inventory did not have GDTYP1, P_ALP1, P_BET1, P_GAM1 variables, how to do this? thank you.

10

For normal Lat-Lon, use 0.0d0 for map-projection parameters.
XORIG_I,YORIG_I are the lon and lat of the SW corner of the (1,1)-cell.

11

thanks you!Professors.I have solved my question!