Skip to content

Commit 19601b9

Browse files
committed
pygmt.grdfilter: Add parameters filter_type/filter_width and more to set the filter
1 parent 4610d88 commit 19601b9

2 files changed

Lines changed: 196 additions & 18 deletions

File tree

pygmt/src/grdfilter.py

Lines changed: 188 additions & 16 deletions
Original file line numberDiff line numberDiff line change
@@ -9,16 +9,126 @@
99
from pygmt._typing import PathLike
1010
from pygmt.alias import Alias, AliasSystem
1111
from pygmt.clib import Session
12+
from pygmt.exceptions import GMTParameterError, GMTValueError
1213
from pygmt.helpers import build_arg_list, fmt_docstring, use_alias
1314

1415
__doctest_skip__ = ["grdfilter"]
1516

1617

18+
def _alias_option_F( # noqa: N802
19+
filter_type=None,
20+
filter_width=None,
21+
hist_bin_width=None,
22+
highpass=None,
23+
median_quantile=None,
24+
hist_center_bins=None,
25+
mode_extreme=None,
26+
filter=None, # noqa: A002
27+
):
28+
"""
29+
Helper function to create the alias list for the -F option.
30+
"""
31+
if filter is not None:
32+
kwdict = {
33+
"filter_type": filter_type,
34+
"filter_width": filter_width,
35+
"hist_bin_width": hist_bin_width,
36+
"highpass": highpass,
37+
"median_quantile": median_quantile,
38+
"hist_center_bins": hist_center_bins,
39+
"mode_extreme": mode_extreme,
40+
}
41+
if any(v is not None and v is not False for v in kwdict.values()):
42+
raise GMTParameterError(
43+
conflicts_with=("filter", kwdict.keys()),
44+
reason="'filter' is specified using the unrecommended GMT command string syntax.",
45+
)
46+
return Alias(filter, name="filter") # Deprecated raw GMT string.
47+
48+
if median_quantile is not None and filter_type != "median":
49+
raise GMTValueError(
50+
filter_type,
51+
description="filter type",
52+
reason="'filter_type' must be 'median' when 'median_quantile' is set.",
53+
)
54+
if hist_bin_width is not None and filter_type != "histogram":
55+
raise GMTValueError(
56+
filter_type,
57+
description="filter type",
58+
reason="'filter_type' must be 'histogram' when 'hist_bin_width' is set.",
59+
)
60+
if hist_center_bins is not None and filter_type != "histogram":
61+
raise GMTValueError(
62+
filter_type,
63+
description="filter type",
64+
reason="'filter_type' must be 'histogram' when 'hist_center_bins' is set.",
65+
)
66+
if mode_extreme is not None and filter_type not in {"mlprob", "histogram"}:
67+
raise GMTValueError(
68+
filter_type,
69+
description="filter type",
70+
reason="'filter_type' must be 'mlprob' or 'histogram' when 'mode_extreme' is set.",
71+
)
72+
73+
return [
74+
Alias(
75+
filter_type,
76+
name="filter_type",
77+
mapping={
78+
"boxcar": "b",
79+
"cosine_arch": "c",
80+
"gaussian": "g",
81+
"custom": "f",
82+
"operator": "o",
83+
"median": "m",
84+
"mlprob": "p",
85+
"histogram": "h",
86+
"minall": "l",
87+
"minpos": "L",
88+
"maxall": "u",
89+
"maxneg": "U",
90+
},
91+
),
92+
Alias(filter_width, name="filter_width", sep="/"),
93+
Alias(hist_bin_width, name="hist_bin_width", prefix="/"),
94+
Alias(hist_center_bins, name="hist_center_bins", prefix="+c"),
95+
Alias(highpass, name="highpass", prefix="+h"),
96+
Alias(median_quantile, name="median_quantile", prefix="+q"),
97+
Alias(
98+
mode_extreme,
99+
name="mode_extreme",
100+
mapping={"min": "+l", "max": "+u"},
101+
),
102+
]
103+
104+
17105
@fmt_docstring
18-
@use_alias(D="distance", F="filter", f="coltypes")
19-
def grdfilter(
106+
@use_alias(D="distance", f="coltypes")
107+
def grdfilter( # noqa: PLR0913
20108
grid: PathLike | xr.DataArray,
21109
outgrid: PathLike | None = None,
110+
filter_type: Literal[
111+
"boxcar",
112+
"cosine_arch",
113+
"gaussian",
114+
"custom",
115+
"operator",
116+
"median",
117+
"mlprob",
118+
"histogram",
119+
"minall",
120+
"minpos",
121+
"maxall",
122+
"maxneg",
123+
]
124+
| None = None,
125+
filter_width: Sequence[float] | None = None,
126+
highpass: bool = False,
127+
median_quantile: float | None = None,
128+
hist_bin_width: float | None = None,
129+
hist_center_bins: bool = False,
130+
mode_extreme: Literal["min", "max"] | None = None,
131+
filter: str | None = None, # noqa: A002
22132
spacing: Sequence[float | str] | None = None,
23133
nans: Literal["ignore", "replace", "preserve"] | None = None,
24134
toggle: bool = False,
@@ -58,19 +168,68 @@ def grdfilter(
58168
----------
59169
$grid
60170
$outgrid
61-
filter : str
62-
**b**\|\ **c**\|\ **g**\|\ **o**\|\ **m**\|\ **p**\|\ **h**\ *width*\
63-
[/*width2*\][*modifiers*].
64-
Name of the filter type you wish to apply, followed by the *width*:
65-
66-
- **b**: Box Car
67-
- **c**: Cosine Arch
68-
- **g**: Gaussian
69-
- **o**: Operator
70-
- **m**: Median
71-
- **p**: Maximum Likelihood probability
72-
- **h**: Histogram
171+
filter_type
172+
The filter type. Choose among convolution and non-convolution filters.
73173
174+
Convolution filters are:
175+
176+
- ``"boxcar"``: All weights are equal.
177+
- ``"cosine_arch"``: Weights follow a cosine arch curve.
178+
- ``"gaussian"``: Weights are given by the Gaussian function, where filter width
179+
is 6 times the conventional Gaussian sigma.
180+
- ``"custom"``: Weights are given by the precomputed values in the filter weight
181+
grid file *weight*, which must have odd dimensions; also requires ``distance=0``
182+
and output spacing must match input spacing or be integer multiples.
183+
- ``"operator"``: Weights are given by the precomputed values in the filter weight
184+
grid file *weight*, which must have odd dimensions; also requires ``distance=0``
185+
and output spacing must match input spacing or be integer multiples. Weights
186+
are assumed to sum to zero so no accumulation of weight sums and normalization
187+
will be done.
188+
189+
Non-convolution filters are:
190+
191+
- ``"median"``: Returns median value. To select another quantile, use the
192+
parameter ``median_quantile`` in the 0-1 range [Default is 0.5, i.e., median].
193+
- ``"mlprob"``: Maximum likelihood probability (a mode estimator). Return modal
194+
value. If more than one mode is found we return their average value. Set
195+
``mode_extreme`` to ``"min"`` or ``"max"`` to return the lowermost or uppermost
196+
of the modal values.
197+
- ``"histogram"``: Histogram mode (another mode estimator). Return the modal value
198+
as the center of the dominant peak in a histogram. Use parameter
199+
``histogram_center_bins`` to center the bins on multiples of bin width [Default
200+
has bin edges that are multiples of bin width]. Use parameter
201+
``histogram_bin_width`` to set the bin width. If more than one mode is found we
202+
return their average value. Set ``mode_extreme`` to ``"min"`` or ``"max"`` to
203+
return the lowermost or uppermost of the modal values.
204+
205+
- ``"minall"``: Return minimum of all values.
206+
- ``"minpos"``: Return minimum of all positive values only.
207+
- ``"maxall"``: Return maximum of all values.
208+
- ``"maxneg"``: Return maximum of all negative values only.
209+
median_quantile
210+
Quantile to use when ``filter_type="median"``. Must be a float in the range 0-1.
211+
[Default is 0.5 (median)].
212+
hist_bin_width
213+
Bin width to use when ``filter_type="histogram"``.
214+
hist_center_bins
215+
Center the histogram bins on multiples of *histogram_bin_width* when
216+
``filter_type="histogram"``. By default, the bins are aligned such that
217+
their edges are on multiples of *hist_bin_width*.
218+
mode_extreme
219+
Choose which extreme to return when ``filter_type="mlprob"`` or
220+
``filter_type="histogram"`` and multiple modes are found. Options are: ``"min"``
221+
to return the lowermost mode, or ``"max"`` to return the uppermost mode. By
222+
default, the average of all modes is returned.
223+
filter
224+
Set the filter type.
225+
226+
.. deprecated:: v0.19.0
227+
228+
This parameter is deprecated. Use the parameters ``filter_type``,
229+
``filter_width``, ``hist_bin_width``, ``highpass``, ``median_quantile``,
230+
``hist_center_bins``, and ``mode_extreme`` instead. This parameter still
231+
accepts raw GMT CLI strings for the ``-F`` option of the ``grdfilter``
232+
module for backward compatibility.
74233
distance : str
75234
State how the grid (x,y) relates to the filter *width*:
76235
@@ -130,7 +289,8 @@ def grdfilter(
130289
>>> # and return a filtered grid (saved as netCDF file).
131290
>>> pygmt.grdfilter(
132291
... grid="@earth_relief_30m_g",
133-
... filter="m600",
292+
... filter_type="median",
293+
... filter_width=600,
134294
... distance="4",
135295
... region=[150, 250, 10, 40],
136296
... spacing=0.5,
@@ -140,9 +300,21 @@ def grdfilter(
140300
>>> # Apply a Gaussian smoothing filter of 600 km to the input DataArray and return
141301
>>> # a filtered DataArray with the smoothed grid.
142302
>>> grid = pygmt.datasets.load_earth_relief()
143-
>>> smooth_field = pygmt.grdfilter(grid=grid, filter="g600", distance="4")
303+
>>> smooth_field = pygmt.grdfilter(
304+
... grid=grid, filter_type="gaussian", filter_width=600, distance="4"
305+
... )
144306
"""
145307
aliasdict = AliasSystem(
308+
F=_alias_option_F(
309+
fliter_type=filter_type,
310+
filter_width=filter_width,
311+
hist_bin_width=hist_bin_width,
312+
highpass=highpass,
313+
median_quantile=median_quantile,
314+
hist_center_bins=hist_center_bins,
315+
mode_extreme=mode_extreme,
316+
filter=filter,
317+
),
146318
I=Alias(spacing, name="spacing", sep="/", size=2),
147319
N=Alias(
148320
nans, name="nans", mapping={"ignore": "i", "replace": "r", "preserve": "p"}

pygmt/tests/test_grdfilter.py

Lines changed: 8 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -47,7 +47,12 @@ def test_grdfilter_dataarray_in_dataarray_out(grid, expected_grid):
4747
Test grdfilter with an input DataArray, and output as DataArray.
4848
"""
4949
result = grdfilter(
50-
grid=grid, filter="g600", distance="4", region=[-53, -49, -20, -17], cores=2
50+
grid=grid,
51+
filter_type="gaussian",
52+
filter_width=600,
53+
distance="4",
54+
region=[-53, -49, -20, -17],
55+
cores=2,
5156
)
5257
# check information of the output grid
5358
assert isinstance(result, xr.DataArray)
@@ -65,7 +70,8 @@ def test_grdfilter_dataarray_in_file_out(grid, expected_grid):
6570
result = grdfilter(
6671
grid,
6772
outgrid=tmpfile.name,
68-
filter="g600",
73+
filter_type="gaussian",
74+
filter_width=600,
6975
distance="4",
7076
region=[-53, -49, -20, -17],
7177
)

0 commit comments

Comments
 (0)