Creating a Landmask File¶

Pre-Made North American Landmask¶

Landmask files created by NOAA for North America are included in the distribution. If you want to use these files, then you do not need to follow the steps below for creating a landmask file. A global regions file created from CarbonTracker data is included, called ‘ctregions.nc’. A landmask file created by NOAA for North America is included called ‘landmask_na.npy’.

Purpose¶

A landmask file is needed to determine which grid cells in the domain under study are designated as land cells. This mask is based on the Carbontracker bio fluxes. A land mask is created where it has a value of 1 where there is a Carbontracker bio flux, and a value of 0 otherwise. There are 2 steps involved, the creation of a ‘regions.nc’ file which has the 0/1 values for a global grid, and the creation of a mask file which can be used for extracting land values from a 2d gridded array.

Regions Creation¶

The first step is to run the mklandmask.py script, which uses a carbontracker flux file and creates a netcdf file which has 0/1 values defining ocean/land grid cells. The flux file that is used is the 1x1 degree flux file for the same date as the start of the inversion, such as ‘fluxes/three-hourly/CT2013B.flux1x1.20120601.nc’ from the carbontracker ftp distribution at /aftp/products/carbontracker/co2/.

The format of the netcdf file that is created has the name ‘ctlpdmregions.nc’ and looks like:

netcdf ctlpdmregions {
dimensions:
        nlat = 180 ;
        nlon = 360 ;
variables:
        double lat(nlat) ;
                lat:units = "degrees_north" ;
                lat:long_name = "latitude" ;
        double lon(nlon) ;
                lon:units = "degrees_east" ;
                lon:long_name = "longitude" ;
        double regions(nlat, nlon) ;
                regions:_FillValue = -1.e+34 ;
                regions:units = "NA" ;
                regions:long_name = "land regions" ;
}

Landmask Creation¶

The script mklm.py is used to create the landmask file. The domain must be specified in the configuration file.

Prerequisites¶

The domain must be specfied in the configuration file.

Example:

# domain latitude and longitude boundaries
south = 10
north = 80
west = -170
east = -50

# landmask file
# numpy binary file containing land/ocean designation for domain
landmask_file = 'landmask_na.npy'

Usage¶

mklm.py is run with:

python mklm.py [-c configfile] region_file

where
        region_file : A netcdf file with a 'regions' variable for the globe,
                      which has 1/0 values for land/ocean grid cells.
        configfile  : configuration file. If not specified, use 'config.ini'.

python mklm.py /data/work/regions/ctregions.nc

Output¶

Three files are produced.

• region.npy - Contains a 2d array of nlatgrids x nlongrids with 0 or 1 designation for non land/land. a[0][0] is the south west corner of the domain.
• regions.txt - Text file with nlatgrids lines with nlongrids columns with 0 or 1 designation for non land/land
• landmask.npy - Contains a 2d array of 2 x nlandgridcells, where row 0 is array of grid cell indices for land cells by latitude, and row 1 is array of grid cell indices for land cells by longitude. For example, if row[0][0] = 0 and row[1][0] = 84, then the southernmost latitude band (0) and the 84th cell from the west edge is a land cell. The actual file name is specified in the configuration file with the key ‘landmask_file’.

Description¶

The landmask file is used for extracting values from a gridded data set for land only grid cells. For example, let’s say we have a domain defined as:

south = 10
north = 80
west = -170
east = -50
lat_resolution = 1
lon_resolution = 1

This would give us 70 latitude cells by 120 longitude cells, for a total of 8400 cells. If we have a 2d array of values over this grid a, and we wanted only the values that were in land cells, we can apply the landmask to this array. In python, we would do this like:

>>> # create an array 70x120 of random numbers
>>> a = numpy.random.rand(70,120)
>>> # extract only the land cells
>>> landonly = a[landmask[0], landmask[1]]
>>> print landonly.shape
(3470,)

which gives us a 1d array with 3470 elements, which is the number of land cells in our land mask. This is the same length as each of the rows in the landmask array. The indices in the landonly array can be used to get the actual latitude and longitude. So index 0 of landonly corresponds to latitude cell number landmask[0][0], longitude cell number landmask[1][0]. The actual latitude would be south + landmask[0][0] + 0.5 (the 0.5 centers the location in the grid cell). The longitude would be west + landmask[1][0] + 0.5. For the NOAA landmask, landmask[0][0] = 0, landmask[1][0] = 84, so latitude is 10.5, and longitude is -85.5 (Costa Rica).

These calculations are made when reading the configuration file, and are available as config[“lats”] and config[“lons”], for example:

>>> from utils import *
>>> config = getConfig("config.ini")
>>> config["lats"]
array([ 10.5,  10.5,  10.5, ...,  79.5,  79.5,  79.5])
>>> config["lons"]
array([-85.5, -84.5, -83.5, ..., -81.5, -80.5, -75.5])

So, given a land cell number, ranging from 0 to nlandcells-1, we can convert to latitude and longitude with:

latitude = south + landmask[0][cellnum] + 0.5
longitude = west + landmask[1][cellnum] + 0.5

or in python with:

>>> landmask = numpy.load("landmask_na.npy")
>>> cellnum = 0
>>> latitude = config["south"] + landmask[0][cellnum] + 0.5
>>> longitude = config["west"] + landmask[1][cellnum] + 0.5
>>> print latitude, longitude
10.5 -85.5
>>> # use precomputed values
>>> latitude = config["lats"][cellnum]
>>> longitude = config["lons"][cellnum]
>>> print latitude, longitude
10.5 -85.5
>>> # or use convenience function in utils.py
>>> getCellLatLon(config, 0)
(10.5, -85.5)