Generate land sea mask

This notebook provides a method to generate a landmask for lat-lon grid using the PMP’s create_target_grid and create_land_sea_mask functions.

Jiwoo Lee (LLNL), 2025-1-15

from pcmdi_metrics.utils import create_target_grid
from pcmdi_metrics.utils import create_land_sea_mask

Generate your grid

Define your resolution in degree.

#resolution = '2.5x2.5'
resolution = '1x1'
#resolution = '0.5x0.5'

grid = create_target_grid(target_grid_resolution=resolution)

grid

<xarray.Dataset> Size: 13kB
Dimensions:   (lon: 360, bnds: 2, lat: 180)
Coordinates:
  * lon       (lon) float64 3kB 0.5 1.5 2.5 3.5 4.5 ... 356.5 357.5 358.5 359.5
  * lat       (lat) float64 1kB -89.5 -88.5 -87.5 -86.5 ... 86.5 87.5 88.5 89.5
Dimensions without coordinates: bnds
Data variables:
    lon_bnds  (lon, bnds) float64 6kB 0.0 1.0 1.0 2.0 ... 359.0 359.0 360.0
    lat_bnds  (lat, bnds) float64 3kB -90.0 -89.0 -89.0 -88.0 ... 89.0 89.0 90.0

xarray.Dataset

• Dimensions:
  • lon: 360
  • bnds: 2
  • lat: 180
• Coordinates: (2)
  • lon
    (lon)
    float64
    0.5 1.5 2.5 ... 357.5 358.5 359.5

    units :

    degrees_east

    axis :

    X

    coordinate :

    longitude

    bounds :

    lon_bnds

    array([  0.5,   1.5,   2.5, ..., 357.5, 358.5, 359.5])

  • lat
    (lat)
    float64
    -89.5 -88.5 -87.5 ... 88.5 89.5

    units :

    degrees_north

    axis :

    Y

    coordinate :

    latitude

    bounds :

    lat_bnds

    array([-89.5, -88.5, -87.5, -86.5, -85.5, -84.5, -83.5, -82.5, -81.5, -80.5,
           -79.5, -78.5, -77.5, -76.5, -75.5, -74.5, -73.5, -72.5, -71.5, -70.5,
           -69.5, -68.5, -67.5, -66.5, -65.5, -64.5, -63.5, -62.5, -61.5, -60.5,
           -59.5, -58.5, -57.5, -56.5, -55.5, -54.5, -53.5, -52.5, -51.5, -50.5,
           -49.5, -48.5, -47.5, -46.5, -45.5, -44.5, -43.5, -42.5, -41.5, -40.5,
           -39.5, -38.5, -37.5, -36.5, -35.5, -34.5, -33.5, -32.5, -31.5, -30.5,
           -29.5, -28.5, -27.5, -26.5, -25.5, -24.5, -23.5, -22.5, -21.5, -20.5,
           -19.5, -18.5, -17.5, -16.5, -15.5, -14.5, -13.5, -12.5, -11.5, -10.5,
            -9.5,  -8.5,  -7.5,  -6.5,  -5.5,  -4.5,  -3.5,  -2.5,  -1.5,  -0.5,
             0.5,   1.5,   2.5,   3.5,   4.5,   5.5,   6.5,   7.5,   8.5,   9.5,
            10.5,  11.5,  12.5,  13.5,  14.5,  15.5,  16.5,  17.5,  18.5,  19.5,
            20.5,  21.5,  22.5,  23.5,  24.5,  25.5,  26.5,  27.5,  28.5,  29.5,
            30.5,  31.5,  32.5,  33.5,  34.5,  35.5,  36.5,  37.5,  38.5,  39.5,
            40.5,  41.5,  42.5,  43.5,  44.5,  45.5,  46.5,  47.5,  48.5,  49.5,
            50.5,  51.5,  52.5,  53.5,  54.5,  55.5,  56.5,  57.5,  58.5,  59.5,
            60.5,  61.5,  62.5,  63.5,  64.5,  65.5,  66.5,  67.5,  68.5,  69.5,
            70.5,  71.5,  72.5,  73.5,  74.5,  75.5,  76.5,  77.5,  78.5,  79.5,
            80.5,  81.5,  82.5,  83.5,  84.5,  85.5,  86.5,  87.5,  88.5,  89.5])

• Data variables: (2)
  • lon_bnds
    (lon, bnds)
    float64
    0.0 1.0 1.0 ... 359.0 359.0 360.0

    xcdat_bounds :

    True

    array([[  0.,   1.],
           [  1.,   2.],
           [  2.,   3.],
           [  3.,   4.],
           [  4.,   5.],
           [  5.,   6.],
           [  6.,   7.],
           [  7.,   8.],
           [  8.,   9.],
           [  9.,  10.],
           [ 10.,  11.],
           [ 11.,  12.],
           [ 12.,  13.],
           [ 13.,  14.],
           [ 14.,  15.],
           [ 15.,  16.],
           [ 16.,  17.],
           [ 17.,  18.],
           [ 18.,  19.],
           [ 19.,  20.],
    ...
           [340., 341.],
           [341., 342.],
           [342., 343.],
           [343., 344.],
           [344., 345.],
           [345., 346.],
           [346., 347.],
           [347., 348.],
           [348., 349.],
           [349., 350.],
           [350., 351.],
           [351., 352.],
           [352., 353.],
           [353., 354.],
           [354., 355.],
           [355., 356.],
           [356., 357.],
           [357., 358.],
           [358., 359.],
           [359., 360.]])

  • lat_bnds
    (lat, bnds)
    float64
    -90.0 -89.0 -89.0 ... 89.0 90.0

    xcdat_bounds :

    True

    array([[-90., -89.],
           [-89., -88.],
           [-88., -87.],
           [-87., -86.],
           [-86., -85.],
           [-85., -84.],
           [-84., -83.],
           [-83., -82.],
           [-82., -81.],
           [-81., -80.],
           [-80., -79.],
           [-79., -78.],
           [-78., -77.],
           [-77., -76.],
           [-76., -75.],
           [-75., -74.],
           [-74., -73.],
           [-73., -72.],
           [-72., -71.],
           [-71., -70.],
    ...
           [ 70.,  71.],
           [ 71.,  72.],
           [ 72.,  73.],
           [ 73.,  74.],
           [ 74.,  75.],
           [ 75.,  76.],
           [ 76.,  77.],
           [ 77.,  78.],
           [ 78.,  79.],
           [ 79.,  80.],
           [ 80.,  81.],
           [ 81.,  82.],
           [ 82.,  83.],
           [ 83.,  84.],
           [ 84.,  85.],
           [ 85.,  86.],
           [ 86.,  87.],
           [ 87.,  88.],
           [ 88.,  89.],
           [ 89.,  90.]])

• Indexes: (2)
  • lon
    PandasIndex

    PandasIndex(Index([  0.5,   1.5,   2.5,   3.5,   4.5,   5.5,   6.5,   7.5,   8.5,   9.5,
           ...
           350.5, 351.5, 352.5, 353.5, 354.5, 355.5, 356.5, 357.5, 358.5, 359.5],
          dtype='float64', name='lon', length=360))

  • lat
    PandasIndex

    PandasIndex(Index([-89.5, -88.5, -87.5, -86.5, -85.5, -84.5, -83.5, -82.5, -81.5, -80.5,
           ...
            80.5,  81.5,  82.5,  83.5,  84.5,  85.5,  86.5,  87.5,  88.5,  89.5],
          dtype='float64', name='lat', length=180))

• Attributes: (0)

Generate land sea mask for the grid

mask = create_land_sea_mask(grid)

mask.plot()

<matplotlib.collections.QuadMesh at 0x151925630>

Optional: Apply different landsea mask generation method

This option uses the PCMDI method developed by Taylor and Doutriaux (2000) (https://pcmdi.llnl.gov/report/ab58.html)

mask2 = create_land_sea_mask(grid, method="pcmdi")

mask2.plot()

<matplotlib.collections.QuadMesh at 0x151b6f4f0>

Save to a netCDF file

ncfile = f"landmask_{resolution}_regionmask.nc"
ncfile2 = f"landmask_{resolution}_pcmdi.nc"

mask.to_netcdf(ncfile)
mask2.to_netcdf(ncfile2)

Quick check for the saved file

import xarray as xr
ds = xr.open_dataset(ncfile)

ds

<xarray.Dataset> Size: 523kB
Dimensions:  (lat: 180, lon: 360)
Coordinates:
  * lat      (lat) float64 1kB -89.5 -88.5 -87.5 -86.5 ... 86.5 87.5 88.5 89.5
  * lon      (lon) float64 3kB 0.5 1.5 2.5 3.5 4.5 ... 356.5 357.5 358.5 359.5
Data variables:
    mask     (lat, lon) int64 518kB ...

xarray.Dataset

• Dimensions:
  • lat: 180
  • lon: 360
• Coordinates: (2)
  • lat
    (lat)
    float64
    -89.5 -88.5 -87.5 ... 88.5 89.5

    units :

    degrees_north

    axis :

    Y

    coordinate :

    latitude

    bounds :

    lat_bnds

    array([-89.5, -88.5, -87.5, -86.5, -85.5, -84.5, -83.5, -82.5, -81.5, -80.5,
           -79.5, -78.5, -77.5, -76.5, -75.5, -74.5, -73.5, -72.5, -71.5, -70.5,
           -69.5, -68.5, -67.5, -66.5, -65.5, -64.5, -63.5, -62.5, -61.5, -60.5,
           -59.5, -58.5, -57.5, -56.5, -55.5, -54.5, -53.5, -52.5, -51.5, -50.5,
           -49.5, -48.5, -47.5, -46.5, -45.5, -44.5, -43.5, -42.5, -41.5, -40.5,
           -39.5, -38.5, -37.5, -36.5, -35.5, -34.5, -33.5, -32.5, -31.5, -30.5,
           -29.5, -28.5, -27.5, -26.5, -25.5, -24.5, -23.5, -22.5, -21.5, -20.5,
           -19.5, -18.5, -17.5, -16.5, -15.5, -14.5, -13.5, -12.5, -11.5, -10.5,
            -9.5,  -8.5,  -7.5,  -6.5,  -5.5,  -4.5,  -3.5,  -2.5,  -1.5,  -0.5,
             0.5,   1.5,   2.5,   3.5,   4.5,   5.5,   6.5,   7.5,   8.5,   9.5,
            10.5,  11.5,  12.5,  13.5,  14.5,  15.5,  16.5,  17.5,  18.5,  19.5,
            20.5,  21.5,  22.5,  23.5,  24.5,  25.5,  26.5,  27.5,  28.5,  29.5,
            30.5,  31.5,  32.5,  33.5,  34.5,  35.5,  36.5,  37.5,  38.5,  39.5,
            40.5,  41.5,  42.5,  43.5,  44.5,  45.5,  46.5,  47.5,  48.5,  49.5,
            50.5,  51.5,  52.5,  53.5,  54.5,  55.5,  56.5,  57.5,  58.5,  59.5,
            60.5,  61.5,  62.5,  63.5,  64.5,  65.5,  66.5,  67.5,  68.5,  69.5,
            70.5,  71.5,  72.5,  73.5,  74.5,  75.5,  76.5,  77.5,  78.5,  79.5,
            80.5,  81.5,  82.5,  83.5,  84.5,  85.5,  86.5,  87.5,  88.5,  89.5])

  • lon
    (lon)
    float64
    0.5 1.5 2.5 ... 357.5 358.5 359.5

    units :

    degrees_east

    axis :

    X

    coordinate :

    longitude

    bounds :

    lon_bnds

    array([  0.5,   1.5,   2.5, ..., 357.5, 358.5, 359.5])

• Data variables: (1)
  • mask
    (lat, lon)
    int64
    ...

    [64800 values with dtype=int64]

• Indexes: (2)
  • lat
    PandasIndex

    PandasIndex(Index([-89.5, -88.5, -87.5, -86.5, -85.5, -84.5, -83.5, -82.5, -81.5, -80.5,
           ...
            80.5,  81.5,  82.5,  83.5,  84.5,  85.5,  86.5,  87.5,  88.5,  89.5],
          dtype='float64', name='lat', length=180))

  • lon
    PandasIndex

    PandasIndex(Index([  0.5,   1.5,   2.5,   3.5,   4.5,   5.5,   6.5,   7.5,   8.5,   9.5,
           ...
           350.5, 351.5, 352.5, 353.5, 354.5, 355.5, 356.5, 357.5, 358.5, 359.5],
          dtype='float64', name='lon', length=360))

• Attributes: (0)

ds["mask"].plot()

<matplotlib.collections.QuadMesh at 0x151e08a30>