diff --git a/doc/api/index.rst b/doc/api/index.rst
index fc69a95d1e2..fcda91042ee 100644
--- a/doc/api/index.rst
+++ b/doc/api/index.rst
@@ -93,6 +93,7 @@ Operations on grids:
     grdfill
     grdfilter
     grdtrack
+    xyz2grd
 
 Crossover analysis with x2sys:
 
diff --git a/pygmt/__init__.py b/pygmt/__init__.py
index 5cbea0fc853..3041bf3f21c 100644
--- a/pygmt/__init__.py
+++ b/pygmt/__init__.py
@@ -45,6 +45,7 @@
     which,
     x2sys_cross,
     x2sys_init,
+    xyz2grd,
 )
 
 # Get semantic version through setuptools-scm
diff --git a/pygmt/src/__init__.py b/pygmt/src/__init__.py
index b4855a37737..6330ecb8f93 100644
--- a/pygmt/src/__init__.py
+++ b/pygmt/src/__init__.py
@@ -39,3 +39,4 @@
 from pygmt.src.wiggle import wiggle
 from pygmt.src.x2sys_cross import x2sys_cross
 from pygmt.src.x2sys_init import x2sys_init
+from pygmt.src.xyz2grd import xyz2grd
diff --git a/pygmt/src/xyz2grd.py b/pygmt/src/xyz2grd.py
new file mode 100644
index 00000000000..34e09faff4f
--- /dev/null
+++ b/pygmt/src/xyz2grd.py
@@ -0,0 +1,97 @@
+"""
+xyz2grd - Create a grid file from an xyz table.
+"""
+import xarray as xr
+from pygmt.clib import Session
+from pygmt.exceptions import GMTInvalidInput
+from pygmt.helpers import (
+    GMTTempFile,
+    build_arg_string,
+    data_kind,
+    dummy_context,
+    fmt_docstring,
+    kwargs_to_strings,
+    use_alias,
+)
+
+
+@fmt_docstring
+@use_alias(
+    G="outgrid",
+    R="region",
+    I="increment",
+)
+@kwargs_to_strings(R="sequence")
+def xyz2grd(table, **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.
+    increment : str
+        *xinc*\ [**+e**\|\ **n**][/\ *yinc*\ [**+e**\|\ **n**]].
+        *x_inc* [and optionally *y_inc*] is the grid spacing. **Geographical
+        (degrees) coordinates**: Optionally, append a increment unit. Choose
+        among **m** to indicate arc minutes or **s** to indicate arc seconds.
+        If one of the units **e**, **f**, **k**, **M**, **n** or **u** is
+        appended instead, the increment is assumed to be given in meter, foot,
+        km, mile, nautical mile or US survey foot, respectively, and will be
+        converted to the equivalent degrees longitude at the middle latitude
+        of the region. If *y_inc* is given but set to 0 it will be reset equal
+        to *x_inc*; otherwise it will be converted to degrees latitude. **All
+        coordinates**: If **+e** is appended then the corresponding max
+        *x* (*east*) or *y* (*north*) may be slightly adjusted to fit exactly
+        the given increment [by default the increment may be adjusted slightly
+        to fit the given domain]. Finally, instead of giving an increment you
+        may specify the *number of nodes* desired by appending **+n** to the
+        supplied integer argument; the increment is then recalculated from the
+        number of nodes, the *registration*, and the domain. The resulting
+        increment value depends on whether you have selected a
+        gridline-registered or pixel-registered grid.
+    {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("Region and increment must be specified.")
+    kind = data_kind(table)
+    with GMTTempFile(suffix=".nc") as tmpfile:
+        with Session() as lib:
+            if kind == "file":
+                file_context = dummy_context(table)
+            elif kind == "matrix":
+                file_context = lib.virtualfile_from_matrix(matrix=table)
+            else:
+                raise GMTInvalidInput("Unrecognized data type: {}".format(type(table)))
+            with file_context as infile:
+                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)
+                arg_str = " ".join([infile, arg_str])
+                lib.call_module("xyz2grd", 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