diff --git a/doc/api/index.rst b/doc/api/index.rst index 1e74629c5af..40ab49ac430 100644 --- a/doc/api/index.rst +++ b/doc/api/index.rst @@ -92,6 +92,7 @@ Operations on grids: grdcut grdfill grdfilter + grdlandmask grdgradient grdtrack diff --git a/pygmt/__init__.py b/pygmt/__init__.py index 57e223d5c14..972f6474286 100644 --- a/pygmt/__init__.py +++ b/pygmt/__init__.py @@ -39,6 +39,7 @@ grdfilter, grdgradient, grdinfo, + grdlandmask, grdtrack, info, makecpt, diff --git a/pygmt/src/__init__.py b/pygmt/src/__init__.py index 4e79d0d96f7..f8bae667f02 100644 --- a/pygmt/src/__init__.py +++ b/pygmt/src/__init__.py @@ -18,6 +18,7 @@ from pygmt.src.grdgradient import grdgradient from pygmt.src.grdimage import grdimage from pygmt.src.grdinfo import grdinfo +from pygmt.src.grdlandmask import grdlandmask from pygmt.src.grdtrack import grdtrack from pygmt.src.grdview import grdview from pygmt.src.histogram import histogram diff --git a/pygmt/src/grdlandmask.py b/pygmt/src/grdlandmask.py new file mode 100644 index 00000000000..ad6fb693fdc --- /dev/null +++ b/pygmt/src/grdlandmask.py @@ -0,0 +1,76 @@ +""" +grdlandmask - Create a "wet-dry" mask grid from shoreline data base +""" + +import xarray as xr +from pygmt.clib import Session +from pygmt.exceptions import GMTInvalidInput +from pygmt.helpers import ( + GMTTempFile, + build_arg_string, + fmt_docstring, + kwargs_to_strings, + use_alias, +) + + +@fmt_docstring +@use_alias( + G="outgrid", + I="spacing", + R="region", + r="registration", +) +@kwargs_to_strings(R="sequence") +def grdlandmask(**kwargs): + r""" + Create a grid file with set values for land and water. + + Read the selected shoreline database and create a grid to specify which + nodes in the specified grid are over land or over water. The nodes defined + by the selected region and lattice spacing + will be set according to one of two criteria: (1) land vs water, or + (2) the more detailed (hierarchical) ocean vs land vs lake + vs island vs pond. + + Full option list at :gmt-docs:`grdlandmask.html` + + {aliases} + + Parameters + ---------- + outgrid : str or None + The name of the output netCDF file with extension .nc to store the grid + in. + {I} + {R} + {r} + + Returns + ------- + ret: xarray.DataArray or None + Return type depends on whether the ``outgrid`` parameter is set: + + - :class:`xarray.DataArray` if ``outgrid`` is not set + - None if ``outgrid`` is set (grid output will be stored in file set by + ``outgrid``) + """ + if "I" not in kwargs.keys() or "R" not in kwargs.keys(): + raise GMTInvalidInput("Both 'region' and 'spacing' must be specified.") + + with GMTTempFile(suffix=".nc") as tmpfile: + with Session() as lib: + if "G" not in kwargs.keys(): # if outgrid is unset, output to tempfile + kwargs.update({"G": tmpfile.name}) + outgrid = kwargs["G"] + arg_str = build_arg_string(kwargs) + lib.call_module("grdlandmask", arg_str) + + if outgrid == tmpfile.name: # if user did not set outgrid, return DataArray + with xr.open_dataarray(outgrid) as dataarray: + result = dataarray.load() + _ = result.gmt # load GMTDataArray accessor information + else: + result = None # if user sets an outgrid, return None + + return result diff --git a/pygmt/tests/test_grdlandmask.py b/pygmt/tests/test_grdlandmask.py new file mode 100644 index 00000000000..4b019525ea3 --- /dev/null +++ b/pygmt/tests/test_grdlandmask.py @@ -0,0 +1,44 @@ +""" +Tests for grdlandmask. +""" +import os + +import pytest +from pygmt import grdinfo, grdlandmask +from pygmt.exceptions import GMTInvalidInput +from pygmt.helpers import GMTTempFile + + +def test_grdlandmask_outgrid(): + """ + Creates a grid land mask with an outgrid argument. + """ + with GMTTempFile(suffix=".nc") as tmpfile: + result = grdlandmask(outgrid=tmpfile.name, spacing=1, region=[-5, 5, -5, 5]) + assert result is None # return value is None + assert os.path.exists(path=tmpfile.name) # check that outgrid exists + result = ( + grdinfo(grid=tmpfile.name, force_scan=0, per_column="n").strip().split() + ) + assert result == ["-5", "5", "-5", "5", "0", "1", "1", "1", "11", "11", "0", "1"] + + +def test_grdlandmask_no_outgrid(): + """ + Test grdlandmask with no set outgrid. + """ + temp_grid = grdlandmask(spacing=1, region=[-5, 5, -5, 5]) + assert temp_grid.dims == ("lat", "lon") + assert temp_grid.gmt.gtype == 1 # Geographic grid + assert temp_grid.gmt.registration == 0 # Pixel registration + assert temp_grid.min() == 0 + assert temp_grid.max() == 1 + + +def test_grdlandmask_fails(): + """ + Check that grdlandmask fails correctly when region and spacing are not + given. + """ + with pytest.raises(GMTInvalidInput): + grdlandmask()