Skip to content
1 change: 1 addition & 0 deletions doc/api/index.rst
Original file line number Diff line number Diff line change
Expand Up @@ -90,6 +90,7 @@ Operations on grids:

grdclip
grdcut
grdfill
grdfilter
grdtrack

Expand Down
1 change: 1 addition & 0 deletions pygmt/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -35,6 +35,7 @@
grd2cpt,
grdclip,
grdcut,
grdfill,
grdfilter,
grdinfo,
grdtrack,
Expand Down
1 change: 1 addition & 0 deletions pygmt/src/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand Down
77 changes: 77 additions & 0 deletions pygmt/src/grdfill.py
Original file line number Diff line number Diff line change
@@ -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
Comment thread
willschlitzer marked this conversation as resolved.
Outdated
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 `r^2 = sqrt(X^2 + Y^2)`,

@core-man core-man May 19, 2021

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

LaTeX formula?

Suggested change
pixels [default radius is `r^2 = sqrt(X^2 + Y^2)`,
pixels [default radius is :math:`r^2 = \sqrt{X^2 + Y^2}`,

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

If you want to use the math directive, it should be:

:math:`r^2 = \sqrt{X^2 + Y^2}`

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

If you want to use the math directive, it should be:

:math:`r^2 = \sqrt{X^2 + Y^2}`

updated~ Thanks

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@willschlitzer previously had this in there - it did not work because fmt_docstring tries to insert common text into the brackets. So, it either needs to stay as non-sphinx math or fmt_docstring needs to be modified.

@seisman seisman May 19, 2021

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

If we want LaTeX math, we may need to use double curly braces

:math:`r^2 = \sqrt{{ X^2 + Y^2 }}`

See https://stackoverflow.com/questions/5466451

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

If we want LaTeX math, we may need to use double curly braces

:math:`r^2 = \sqrt{{ X^2 + Y^2 }}`

See https://stackoverflow.com/questions/5466451

This appears to be working!

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
44 changes: 44 additions & 0 deletions pygmt/tests/test_grdfill.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,44 @@
"""
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")
grid[10, 10] = np.nan
Comment thread
willschlitzer marked this conversation as resolved.
Outdated
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[10, 10] == 20
Comment thread
willschlitzer marked this conversation as resolved.
Outdated


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 == "-180 180 -90 90 -8182 5651.5 1 1 360 180 1 1"
Comment thread
willschlitzer marked this conversation as resolved.
Outdated