Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

grdfilter wrapper #612

Closed
wants to merge 4 commits into from
Closed
Show file tree
Hide file tree
Changes from 2 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion pygmt/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -19,7 +19,7 @@
from .sampling import grdtrack
from .mathops import makecpt
from .modules import GMTDataArrayAccessor, config, info, grdinfo, which
from .gridops import grdcut
from .gridops import grdcut, gridfilter
Copy link
Member

Choose a reason for hiding this comment

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

Suggested change
from .gridops import grdcut, gridfilter
from .gridops import grdcut, grdfilter

Just a small typo, that's why it's not importing 😃

from .x2sys import x2sys_init, x2sys_cross
from . import datasets

Expand Down
92 changes: 85 additions & 7 deletions pygmt/gridops.py
Original file line number Diff line number Diff line change
Expand Up @@ -20,13 +20,15 @@

@fmt_docstring
@use_alias(
Copy link
Member

Choose a reason for hiding this comment

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

You're changing the alias settings for grdcut. That's incorrect. You should revert the changes here and add decorators @fmt_docstring, @use_alias and @kwargs_to_strings just before the grdfilter function, i.e., each function has its own decorators.

Copy link
Contributor Author

Choose a reason for hiding this comment

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

Ok.. I thought they could be the same. Thanks!

G="outgrid",
R="region",
J="projection",
N="extend",
S="circ_subregion",
Z="z_subregion",
)
G="outgrid",
R="region",
J="projection",
N="extend",
S="circ_subregion",
Z="z_subregion",
F="filter",
D="distance"
)
@kwargs_to_strings(R="sequence")
def grdcut(grid, **kwargs):
"""
Expand Down Expand Up @@ -113,3 +115,79 @@ def grdcut(grid, **kwargs):
result = None # if user sets an outgrid, return None

return result


def grdfilter(grid, **kwargs):
"""

filter a grid file in the time domain using one of the selected convolution
or non-convolution isotropic or rectangular filters and compute distances
using Cartesian or Spherical geometries. The output grid file can optionally
be generated as a sub-region of the input (via -R) and/or with new increment
(via -I) or registration (via -T). In this way, one may have “extra space” in
the input data so that the edges will not be used and the output can be within
one half-width of the input edges. If the filter is low-pass, then the output
may be less frequently sampled than the input.

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.
{F} : str
Name of filter type you which to apply, followed by the width
b: Box Car; c: Cosine Arch; g: Gaussian; o: Operator; m: Median; p: Maximum Likelihood probability; h: histogram
{D}: str
Distance flag, that tells how grid (x,y) rrlated to the filter width as follows:
flag = p: grid (px,py) with width an odd number of pixels; Cartesian distances.

flag = 0: grid (x,y) same units as width, Cartesian distances.

flag = 1: grid (x,y) in degrees, width in kilometers, Cartesian distances.

flag = 2: grid (x,y) in degrees, width in km, dx scaled by cos(middle y), Cartesian distances.

The above options are fastest because they allow weight matrix to be computed only once. The next three options are slower because they recompute weights for each latitude.

flag = 3: grid (x,y) in degrees, width in km, dx scaled by cosine(y), Cartesian distance calculation.

flag = 4: grid (x,y) in degrees, width in km, Spherical distance calculation.

flag = 5: grid (x,y) in Mercator -Jm1 img units, width in km, Spherical distance calculation.

Returns
-------
ret: xarray.DataArray or None
Return type depends on whether the *outgrid* parameter is set:
- xarray.DataArray if *outgrid* is not set
- None if *outgrid* is set (grid output will be stored in *outgrid*)
"""
Copy link
Member

Choose a reason for hiding this comment

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

The indentation of the docstrings looks wrong.

kind = data_kind(grid)

with GMTTempFile(suffix=".nc") as tmpfile:
with Session() as lib:
if kind == "file":
file_context = dummy_context(grid)
elif kind == "grid":
file_context = lib.virtualfile_from_grid(grid)
else:
raise GMTInvalidInput("Unrecognized data type: {}".format(type(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("grdfilter", 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