Wrap grdfill (#1276)

Wrap the gmt module grdfill and create test_grdfill.py.

Co-authored-by: Yao Jiayuan <coreman.seism@gmail.com>
Co-authored-by: Wei Ji <23487320+weiji14@users.noreply.github.com>

Author: 3 people

doc/api/index.rst

@@ -90,6 +90,7 @@ Operations on grids:
 
     grdclip
     grdcut
+    grdfill
     grdfilter
     grdtrack
 

pygmt/__init__.py

@@ -35,6 +35,7 @@
     grd2cpt,
     grdclip,
     grdcut,
+    grdfill,
     grdfilter,
     grdinfo,
     grdtrack,

pygmt/src/__init__.py

@@ -13,6 +13,7 @@
 from pygmt.src.grdclip import grdclip
 from pygmt.src.grdcontour import grdcontour
 from pygmt.src.grdcut import grdcut
+from pygmt.src.grdfill import grdfill
 from pygmt.src.grdfilter import grdfilter
 from pygmt.src.grdimage import grdimage
 from pygmt.src.grdinfo import grdinfo

pygmt/src/grdfill.py

@@ -0,0 +1,77 @@
+"""
+grdfill - Fill blank areas from a grid.
+"""
+
+import xarray as xr
+from pygmt.clib import Session
+from pygmt.helpers import (
+    GMTTempFile,
+    build_arg_string,
+    fmt_docstring,
+    kwargs_to_strings,
+    use_alias,
+)
+
+
+@fmt_docstring
+@use_alias(
+    A="mode",
+    G="outgrid",
+    R="region",
+)
+@kwargs_to_strings(R="sequence")
+def grdfill(grid, **kwargs):
+    r"""
+    Fill blank areas from a grid file.
+
+    Read a grid that presumably has unfilled holes that the user
+    wants to fill in some fashion. Holes are identified by NaN values but
+    this criteria can be changed. There are several different algorithms that
+    can be used to replace the hole values.
+
+    Full option list at :gmt-docs:`grdfill.html`
+
+    {aliases}
+
+    Parameters
+    ----------
+    grid : str or xarray.DataArray
+        The file name of the input grid or the grid loaded as a DataArray.
+    outgrid : str or None
+        The name of the output netCDF file with extension .nc to store the grid
+        in.
+    mode : str
+        Specify the hole-filling algorithm to use.  Choose from **c** for
+        constant fill and append the constant value, **n** for nearest
+        neighbor (and optionally append a search radius in
+        pixels [default radius is :math:`r^2 = \sqrt{{ X^2 + Y^2 }}`,
+        where (*X,Y*) are the node dimensions of the 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``)
+    """
+    with GMTTempFile(suffix=".nc") as tmpfile:
+        with Session() as lib:
+            file_context = lib.virtualfile_from_data(check_kind="raster", data=grid)
+            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 = " ".join([infile, build_arg_string(kwargs)])
+                lib.call_module("grdfill", 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

pygmt/tests/test_grdfill.py

@@ -0,0 +1,49 @@
+"""
+Tests for grdclip.
+"""
+import os
+
+import numpy as np
+import pytest
+import xarray as xr
+from pygmt import grdfill, grdinfo
+from pygmt.datasets import load_earth_relief
+from pygmt.helpers import GMTTempFile
+
+
+@pytest.fixture(scope="module", name="grid")
+def fixture_grid():
+    """
+    Load the grid data from the sample earth_relief file and set value(s) to
+    NaN.
+    """
+    grid = load_earth_relief(registration="pixel", region=[-5, 5, -5, 5])
+    grid[3:5, 3:5] = np.nan
+    grid[5:7, 5:7] = np.inf
+    return grid
+
+
+def test_grdfill_dataarray_out(grid):
+    """
+    grdfill with a DataArray output.
+    """
+    result = grdfill(grid=grid, mode="c20")
+    # check information of the output grid
+    assert isinstance(result, xr.DataArray)
+    assert result[4, 4] == 20
+    assert result[5, 5] == np.inf
+    assert not result.isnull().all()  # check that no NaN values exists
+    assert result.gmt.gtype == 1  # Geographic grid
+    assert result.gmt.registration == 1  # Pixel registration
+
+
+def test_grdfill_file_out(grid):
+    """
+    grdfill with an outgrid set.
+    """
+    with GMTTempFile(suffix=".nc") as tmpfile:
+        result = grdfill(grid=grid, mode="c20", outgrid=tmpfile.name)
+        assert result is None  # return value is None
+        assert os.path.exists(path=tmpfile.name)  # check that outgrid exists
+        result = grdinfo(tmpfile.name, per_column=True).strip()
+        assert result == "-5 5 -5 5 -5130.5 inf 1 1 10 10 1 1"