From 6cbb18b3d7bc6b947ecd41d5152167b90c2c503f Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Fri, 10 Sep 2021 21:45:57 -0400 Subject: [PATCH 01/22] add dimfilter --- doc/api/index.rst | 1 + pygmt/__init__.py | 1 + pygmt/src/__init__.py | 1 + pygmt/src/dimfilter.py | 83 ++++++++++++++++++++++++++++++++++++++++++ 4 files changed, 86 insertions(+) create mode 100644 pygmt/src/dimfilter.py diff --git a/doc/api/index.rst b/doc/api/index.rst index 0426f79385d..5b8609dd3a6 100644 --- a/doc/api/index.rst +++ b/doc/api/index.rst @@ -88,6 +88,7 @@ Operations on grids: .. autosummary:: :toctree: generated + dimfilter grdclip grdcut grdfill diff --git a/pygmt/__init__.py b/pygmt/__init__.py index 4963c21a163..eec6e9ac2e1 100644 --- a/pygmt/__init__.py +++ b/pygmt/__init__.py @@ -33,6 +33,7 @@ blockmean, blockmedian, config, + dimfilter, grd2cpt, grdclip, grdcut, diff --git a/pygmt/src/__init__.py b/pygmt/src/__init__.py index e70e74367b6..76de8e8ee88 100644 --- a/pygmt/src/__init__.py +++ b/pygmt/src/__init__.py @@ -9,6 +9,7 @@ from pygmt.src.colorbar import colorbar from pygmt.src.config import config from pygmt.src.contour import contour +from pygmt.src.dimfilter import dimfilter from pygmt.src.grd2cpt import grd2cpt from pygmt.src.grdclip import grdclip from pygmt.src.grdcontour import grdcontour diff --git a/pygmt/src/dimfilter.py b/pygmt/src/dimfilter.py new file mode 100644 index 00000000000..53efee333bd --- /dev/null +++ b/pygmt/src/dimfilter.py @@ -0,0 +1,83 @@ +""" +dimfilter - Filter a grid file by dividing the filter circle. +""" + +from pygmt.clib import Session +from pygmt.helpers import ( + GMTTempFile, + build_arg_string, + fmt_docstring, + kwargs_to_strings, + use_alias, +) +from pygmt.io import load_dataarray + + +@fmt_docstring +@use_alias( + D="distance", + F="filter", + G="outgrid", + I="spacing", + N="sectors", + R="region", + V="verbose", +) +@kwargs_to_strings(R="sequence") +def dimfilter(grid, **kwargs): + r""" + Filter a grid file by dividing the filter circle. + + Filter a grid file in the space (or time) domain by + dividing the given filter circle into *n\_sectors*, applying one of the + selected primary convolution or non-convolution filters to each sector, + and choosing the final outcome according to the selected secondary + filter. It computes distances using Cartesian or Spherical geometries. + The output *.nc* file can optionally be generated as a subregion of the + input and/or with a new **-I**\ ncrement. In this way, one may have + "extra space" in the input data so that there will be no edge effects + for the output grid. If the filter is low-pass, then the output may be + less frequently sampled than the input. The **-Q** option is for the + error analysis mode and expects the input file to contains the filtered + depths. Finally, one should know that **dimfilter** will not produce a + smooth output as other spatial filters + do because it returns a minimum median out of *N* medians of *N* + sectors. The output can be rough unless the input data is noise-free. + Thus, an additional filtering (e.g., Gaussian via :doc:`grdfilter`) of the + DiM-filtered data is generally recommended. + + Full option list at :gmt-docs:`dimfilter.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. + {I} + {R} + {V} + + 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("dimfilter", arg_str) + + return load_dataarray(outgrid) if outgrid == tmpfile.name else None From c9affe3cf658963b2de7088d87ad7519be82ec8f Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Sat, 11 Sep 2021 08:03:57 -0400 Subject: [PATCH 02/22] add docstring for sectors, filters, and distance --- pygmt/src/dimfilter.py | 76 ++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 76 insertions(+) diff --git a/pygmt/src/dimfilter.py b/pygmt/src/dimfilter.py index 53efee333bd..fec8bb3139d 100644 --- a/pygmt/src/dimfilter.py +++ b/pygmt/src/dimfilter.py @@ -57,6 +57,82 @@ def dimfilter(grid, **kwargs): outgrid : str or None The name of the output netCDF file with extension .nc to store the grid in. + Distance *flag* determines how grid (x,y) relates to filter *width*, as follows: + distance : int or str + Distance flag tells how grid (x,y) relates to filter width, as follows: + + **0**\ : grid (x,y) in same units as *width*, Cartesian distances. + + **1**\ : grid (x,y) in degrees, *width* in kilometers, Cartesian distances. + + **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 two options are slower because they + recompute weights for each latitude. + + **3**\ : grid (x,y) in degrees, *width* in km, dx scaled by + cosine(y), Cartesian distance calculation. + + **4**\ : grid (x,y) in degrees, *width* in km, Spherical distance + calculation. + filter : str + **x**\ *width*\ [**+l**\|\ **u**] + Sets the primary filter type. Choose among convolution and + non-convolution filters. Append the filter code **x** followed by the full + diameter *width*. Available convolution filters are: + + (**b**) Boxcar: All weights are equal. + + (**c**) Cosine Arch: Weights follow a cosine arch curve. + + (**g**) Gaussian: Weights are given by the Gaussian function. + + Non-convolution filters are: + + (**m**) Median: Returns median value. + + (**p**) Maximum likelihood probability (a mode estimator): Return + modal value. If more than one mode is found we return their average + value. Append **+l** or **+h** to the filter width if you rather want to + return the smallest or largest of each sector's modal values.**x**\ *sectors*\ [**+l**\|\ **u**] + sectors : str + **x**\ *sectors*\ [**+l**\|\ **u**] + Sets the secondary filter type **x** and the number of bow-tie sectors. + *sectors* must be integer and larger than 0. When *sectors* is + set to 1, the secondary filter is not effective. Available secondary + filters **x** are: + + (**l**) Lower: Return the minimum of all filtered values. + + (**u**) Upper: Return the maximum of all filtered values. + + (**a**) Average: Return the mean of all filtered values. + + (**m**) Median: Return the median of all filtered values. + + (**p**) Mode: Return the mode of all filtered values: + If more than one mode is found we return their average + value. Append **+l** or **+h** to the sectors if you rather want to + return the smallest or largest of the modal values.**x**\ *sectors*\ [**+l**\|\ **u**] + Sets the secondary filter type **x** and the number of bow-tie sectors. + *sectors* must be integer and larger than 0. When *sectors* is + set to 1, the secondary filter is not effective. Available secondary + filters **x** are: + + (**l**) Lower: Return the minimum of all filtered values. + + (**u**) Upper: Return the maximum of all filtered values. + + (**a**) Average: Return the mean of all filtered values. + + (**m**) Median: Return the median of all filtered values. + + (**p**) Mode: Return the mode of all filtered values: + If more than one mode is found we return their average + value. Append **+l** or **+h** to the sectors if you rather want to + return the smallest or largest of the modal values. {I} {R} {V} From 17b3c54ac3f37e974657692d475ed74ed86ae06c Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Sat, 11 Sep 2021 08:13:27 -0400 Subject: [PATCH 03/22] run make format --- pygmt/src/dimfilter.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pygmt/src/dimfilter.py b/pygmt/src/dimfilter.py index fec8bb3139d..7d379c17c3b 100644 --- a/pygmt/src/dimfilter.py +++ b/pygmt/src/dimfilter.py @@ -77,7 +77,7 @@ def dimfilter(grid, **kwargs): **4**\ : grid (x,y) in degrees, *width* in km, Spherical distance calculation. - filter : str + filter : str **x**\ *width*\ [**+l**\|\ **u**] Sets the primary filter type. Choose among convolution and non-convolution filters. Append the filter code **x** followed by the full From 3373066be9b90a8d0c06c33f5a434d086127f513 Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Sat, 11 Sep 2021 17:24:34 -0400 Subject: [PATCH 04/22] Add check for required arguments --- pygmt/src/dimfilter.py | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/pygmt/src/dimfilter.py b/pygmt/src/dimfilter.py index 7d379c17c3b..1dc7fc3b6ca 100644 --- a/pygmt/src/dimfilter.py +++ b/pygmt/src/dimfilter.py @@ -3,6 +3,7 @@ """ from pygmt.clib import Session +from pygmt.exceptions import GMTInvalidInput from pygmt.helpers import ( GMTTempFile, build_arg_string, @@ -146,6 +147,11 @@ def dimfilter(grid, **kwargs): - None if ``outgrid`` is set (grid output will be stored in file set by ``outgrid``) """ + if ("D" not in kwargs) or ("F" not in kwargs) or ("N" not in kwargs): + raise GMTInvalidInput( + """The following parameters are required: distance, filters, sectors""" + ) + with GMTTempFile(suffix=".nc") as tmpfile: with Session() as lib: file_context = lib.virtualfile_from_data(check_kind="raster", data=grid) From e3a4f7a04119dab0f04eb22d0edf6e4a7cc9628d Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Sat, 11 Sep 2021 17:24:43 -0400 Subject: [PATCH 05/22] Add test_dimfilter.py --- pygmt/tests/test_dimfilter.py | 51 +++++++++++++++++++++++++++++++++++ 1 file changed, 51 insertions(+) create mode 100644 pygmt/tests/test_dimfilter.py diff --git a/pygmt/tests/test_dimfilter.py b/pygmt/tests/test_dimfilter.py new file mode 100644 index 00000000000..eefa0f598d6 --- /dev/null +++ b/pygmt/tests/test_dimfilter.py @@ -0,0 +1,51 @@ +""" +Tests for dimfilter. +""" +import os + +import numpy.testing as npt +import pytest +from pygmt import dimfilter +from pygmt.datasets import load_earth_relief +from pygmt.exceptions import GMTInvalidInput +from pygmt.helpers import GMTTempFile + + +@pytest.fixture(scope="module", name="grid") +def fixture_grid(): + """ + Load the grid data from the sample earth_relief file. + """ + return load_earth_relief(resolution="30m", region=[-5, 5, -5, 5]) + + +def test_dimfilter_outgrid(grid): + """ + Test the required parameters for dimfilter with a set outgrid. + """ + with GMTTempFile(suffix=".nc") as tmpfile: + result = dimfilter(grid=grid, outgrid=tmpfile.name, filter="m600", distance=4, sectors="l6") + assert result is None # return value is None + assert os.path.exists(path=tmpfile.name) # check that outgrid exists + + +def test_grdgradient_no_outgrid(grid): + """ + Test the required parameters for dimfilter with no set outgrid. + """ + temp_grid = dimfilter(grid=grid, filter="m600", distance=4, sectors="l6") + assert temp_grid.dims == ("lat", "lon") + assert temp_grid.gmt.gtype == 1 # Geographic grid + assert temp_grid.gmt.registration == 1 # Pixel registration + npt.assert_allclose(temp_grid.min(), -5125) + npt.assert_allclose(temp_grid.max(), -3750.5) + npt.assert_allclose(temp_grid.median(), -4826.5) + npt.assert_allclose(temp_grid.mean(), -4789.791) + +def test_dimfilter_fails(grid): + """ + Check that dimfilter fails correctly when sector, + filters, and distance are not specified.. + """ + with pytest.raises(GMTInvalidInput): + dimfilter(grid=grid) \ No newline at end of file From 4fbced3e6ffbdb4c92af7c0905b9486514469235 Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Sat, 11 Sep 2021 17:26:17 -0400 Subject: [PATCH 06/22] run make format --- pygmt/src/dimfilter.py | 2 +- pygmt/tests/test_dimfilter.py | 11 +++++++---- 2 files changed, 8 insertions(+), 5 deletions(-) diff --git a/pygmt/src/dimfilter.py b/pygmt/src/dimfilter.py index 1dc7fc3b6ca..8b1cd9cf6e7 100644 --- a/pygmt/src/dimfilter.py +++ b/pygmt/src/dimfilter.py @@ -150,7 +150,7 @@ def dimfilter(grid, **kwargs): if ("D" not in kwargs) or ("F" not in kwargs) or ("N" not in kwargs): raise GMTInvalidInput( """The following parameters are required: distance, filters, sectors""" - ) + ) with GMTTempFile(suffix=".nc") as tmpfile: with Session() as lib: diff --git a/pygmt/tests/test_dimfilter.py b/pygmt/tests/test_dimfilter.py index eefa0f598d6..61418cd5478 100644 --- a/pygmt/tests/test_dimfilter.py +++ b/pygmt/tests/test_dimfilter.py @@ -24,7 +24,9 @@ def test_dimfilter_outgrid(grid): Test the required parameters for dimfilter with a set outgrid. """ with GMTTempFile(suffix=".nc") as tmpfile: - result = dimfilter(grid=grid, outgrid=tmpfile.name, filter="m600", distance=4, sectors="l6") + result = dimfilter( + grid=grid, outgrid=tmpfile.name, filter="m600", distance=4, sectors="l6" + ) assert result is None # return value is None assert os.path.exists(path=tmpfile.name) # check that outgrid exists @@ -42,10 +44,11 @@ def test_grdgradient_no_outgrid(grid): npt.assert_allclose(temp_grid.median(), -4826.5) npt.assert_allclose(temp_grid.mean(), -4789.791) + def test_dimfilter_fails(grid): """ - Check that dimfilter fails correctly when sector, - filters, and distance are not specified.. + Check that dimfilter fails correctly when sector, filters, and distance are + not specified.. """ with pytest.raises(GMTInvalidInput): - dimfilter(grid=grid) \ No newline at end of file + dimfilter(grid=grid) From 3d55789026212324ffcf122b7b83ff92bc47f59c Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Sat, 11 Sep 2021 17:29:43 -0400 Subject: [PATCH 07/22] run make check --- pygmt/src/dimfilter.py | 29 ++++++----------------------- 1 file changed, 6 insertions(+), 23 deletions(-) diff --git a/pygmt/src/dimfilter.py b/pygmt/src/dimfilter.py index 8b1cd9cf6e7..e2c08bba2a1 100644 --- a/pygmt/src/dimfilter.py +++ b/pygmt/src/dimfilter.py @@ -58,13 +58,13 @@ def dimfilter(grid, **kwargs): outgrid : str or None The name of the output netCDF file with extension .nc to store the grid in. - Distance *flag* determines how grid (x,y) relates to filter *width*, as follows: distance : int or str Distance flag tells how grid (x,y) relates to filter width, as follows: **0**\ : grid (x,y) in same units as *width*, Cartesian distances. - **1**\ : grid (x,y) in degrees, *width* in kilometers, Cartesian distances. + **1**\ : grid (x,y) in degrees, *width* in kilometers, Cartesian + distances. **2**\ : grid (x,y) in degrees, *width* in km, dx scaled by cos(middle y), Cartesian distances. @@ -81,8 +81,8 @@ def dimfilter(grid, **kwargs): filter : str **x**\ *width*\ [**+l**\|\ **u**] Sets the primary filter type. Choose among convolution and - non-convolution filters. Append the filter code **x** followed by the full - diameter *width*. Available convolution filters are: + non-convolution filters. Append the filter code **x** followed by + the full diameter *width*. Available convolution filters are: (**b**) Boxcar: All weights are equal. @@ -96,8 +96,8 @@ def dimfilter(grid, **kwargs): (**p**) Maximum likelihood probability (a mode estimator): Return modal value. If more than one mode is found we return their average - value. Append **+l** or **+h** to the filter width if you rather want to - return the smallest or largest of each sector's modal values.**x**\ *sectors*\ [**+l**\|\ **u**] + value. Append **+l** or **+h** to the filter width if you want + to return the smallest or largest of each sector's modal values. sectors : str **x**\ *sectors*\ [**+l**\|\ **u**] Sets the secondary filter type **x** and the number of bow-tie sectors. @@ -113,23 +113,6 @@ def dimfilter(grid, **kwargs): (**m**) Median: Return the median of all filtered values. - (**p**) Mode: Return the mode of all filtered values: - If more than one mode is found we return their average - value. Append **+l** or **+h** to the sectors if you rather want to - return the smallest or largest of the modal values.**x**\ *sectors*\ [**+l**\|\ **u**] - Sets the secondary filter type **x** and the number of bow-tie sectors. - *sectors* must be integer and larger than 0. When *sectors* is - set to 1, the secondary filter is not effective. Available secondary - filters **x** are: - - (**l**) Lower: Return the minimum of all filtered values. - - (**u**) Upper: Return the maximum of all filtered values. - - (**a**) Average: Return the mean of all filtered values. - - (**m**) Median: Return the median of all filtered values. - (**p**) Mode: Return the mode of all filtered values: If more than one mode is found we return their average value. Append **+l** or **+h** to the sectors if you rather want to From f0e88f41e44aa6ef3bf283a8636ef5a6e7855f5a Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Thu, 23 Dec 2021 08:38:23 +0000 Subject: [PATCH 08/22] Apply suggestions from code review Co-authored-by: Meghan Jones --- pygmt/src/dimfilter.py | 70 ++++++++++++++++++------------------------ 1 file changed, 30 insertions(+), 40 deletions(-) diff --git a/pygmt/src/dimfilter.py b/pygmt/src/dimfilter.py index e2c08bba2a1..f9621a5b794 100644 --- a/pygmt/src/dimfilter.py +++ b/pygmt/src/dimfilter.py @@ -27,14 +27,14 @@ @kwargs_to_strings(R="sequence") def dimfilter(grid, **kwargs): r""" - Filter a grid file by dividing the filter circle. + Filter a grid by dividing the filter circle. - Filter a grid file in the space (or time) domain by - dividing the given filter circle into *n\_sectors*, applying one of the + Filter a grid in the space (or time) domain by + dividing the given filter circle into the given number of sectors, applying one of the selected primary convolution or non-convolution filters to each sector, and choosing the final outcome according to the selected secondary filter. It computes distances using Cartesian or Spherical geometries. - The output *.nc* file can optionally be generated as a subregion of the + The output grid can optionally be generated as a subregion of the input and/or with a new **-I**\ ncrement. In this way, one may have "extra space" in the input data so that there will be no edge effects for the output grid. If the filter is low-pass, then the output may be @@ -44,7 +44,7 @@ def dimfilter(grid, **kwargs): smooth output as other spatial filters do because it returns a minimum median out of *N* medians of *N* sectors. The output can be rough unless the input data is noise-free. - Thus, an additional filtering (e.g., Gaussian via :doc:`grdfilter`) of the + Thus, an additional filtering (e.g., Gaussian via :func:`pygmt.grdfilter`) of the DiM-filtered data is generally recommended. Full option list at :gmt-docs:`dimfilter.html` @@ -61,43 +61,37 @@ def dimfilter(grid, **kwargs): distance : int or str Distance flag tells how grid (x,y) relates to filter width, as follows: - **0**\ : grid (x,y) in same units as *width*, Cartesian distances. - - **1**\ : grid (x,y) in degrees, *width* in kilometers, Cartesian - distances. - - **2**\ : grid (x,y) in degrees, *width* in km, dx scaled by - cos(middle y), Cartesian distances. + - **0**\ : grid (x,y) in same units as *width*, Cartesian distances. + - **1**\ : grid (x,y) in degrees, *width* in kilometers, Cartesian + distances. + - **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 two options are slower because they recompute weights for each latitude. - **3**\ : grid (x,y) in degrees, *width* in km, dx scaled by - cosine(y), Cartesian distance calculation. - - **4**\ : grid (x,y) in degrees, *width* in km, Spherical distance - calculation. + - **3**\ : grid (x,y) in degrees, *width* in km, dx scaled by + cosine(y), Cartesian distance calculation. + - **4**\ : grid (x,y) in degrees, *width* in km, Spherical distance + calculation. filter : str **x**\ *width*\ [**+l**\|\ **u**] Sets the primary filter type. Choose among convolution and non-convolution filters. Append the filter code **x** followed by the full diameter *width*. Available convolution filters are: - (**b**) Boxcar: All weights are equal. - - (**c**) Cosine Arch: Weights follow a cosine arch curve. - - (**g**) Gaussian: Weights are given by the Gaussian function. + - (**b**) Boxcar: All weights are equal. + - (**c**) Cosine Arch: Weights follow a cosine arch curve. + - (**g**) Gaussian: Weights are given by the Gaussian function. Non-convolution filters are: - (**m**) Median: Returns median value. - - (**p**) Maximum likelihood probability (a mode estimator): Return - modal value. If more than one mode is found we return their average - value. Append **+l** or **+h** to the filter width if you want - to return the smallest or largest of each sector's modal values. + - (**m**) Median: Returns median value. + - (**p**) Maximum likelihood probability (a mode estimator): Return + modal value. If more than one mode is found we return their average + value. Append **+l** or **+h** to the filter width if you want + to return the smallest or largest of each sector's modal values. sectors : str **x**\ *sectors*\ [**+l**\|\ **u**] Sets the secondary filter type **x** and the number of bow-tie sectors. @@ -105,18 +99,14 @@ def dimfilter(grid, **kwargs): set to 1, the secondary filter is not effective. Available secondary filters **x** are: - (**l**) Lower: Return the minimum of all filtered values. - - (**u**) Upper: Return the maximum of all filtered values. - - (**a**) Average: Return the mean of all filtered values. - - (**m**) Median: Return the median of all filtered values. - - (**p**) Mode: Return the mode of all filtered values: - If more than one mode is found we return their average - value. Append **+l** or **+h** to the sectors if you rather want to - return the smallest or largest of the modal values. + - (**l**) Lower: Return the minimum of all filtered values. + - (**u**) Upper: Return the maximum of all filtered values. + - (**a**) Average: Return the mean of all filtered values. + - (**m**) Median: Return the median of all filtered values. + - (**p**) Mode: Return the mode of all filtered values: + If more than one mode is found we return their average + value. Append **+l** or **+h** to the sectors if you rather want to + return the smallest or largest of the modal values. {I} {R} {V} From 6866e2180d6b151d04bfa6c7e83cb91149c59f85 Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Thu, 23 Dec 2021 08:45:36 +0000 Subject: [PATCH 09/22] line length fixes --- pygmt/src/dimfilter.py | 28 ++++++++++++++-------------- 1 file changed, 14 insertions(+), 14 deletions(-) diff --git a/pygmt/src/dimfilter.py b/pygmt/src/dimfilter.py index f9621a5b794..e69f8594741 100644 --- a/pygmt/src/dimfilter.py +++ b/pygmt/src/dimfilter.py @@ -30,22 +30,22 @@ def dimfilter(grid, **kwargs): Filter a grid by dividing the filter circle. Filter a grid in the space (or time) domain by - dividing the given filter circle into the given number of sectors, applying one of the - selected primary convolution or non-convolution filters to each sector, - and choosing the final outcome according to the selected secondary - filter. It computes distances using Cartesian or Spherical geometries. - The output grid can optionally be generated as a subregion of the - input and/or with a new **-I**\ ncrement. In this way, one may have - "extra space" in the input data so that there will be no edge effects - for the output grid. If the filter is low-pass, then the output may be - less frequently sampled than the input. The **-Q** option is for the - error analysis mode and expects the input file to contains the filtered - depths. Finally, one should know that **dimfilter** will not produce a - smooth output as other spatial filters + dividing the given filter circle into the given number of sectors, + applying one of the selected primary convolution or non-convolution + filters to each sector, and choosing the final outcome according to the + selected secondary filter. It computes distances using Cartesian or + Spherical geometries. The output grid can optionally be generated as a + subregion of the input and/or with a new **-I**\ ncrement. In this way, + one may have "extra space" in the input data so that there will be no edge + effects for the output grid. If the filter is low-pass, then the output + may be less frequently sampled than the input. The **-Q** option is for + the error analysis mode and expects the input file to contains the + filtered depths. Finally, one should know that **dimfilter** will not + produce a smooth output as other spatial filters do because it returns a minimum median out of *N* medians of *N* sectors. The output can be rough unless the input data is noise-free. - Thus, an additional filtering (e.g., Gaussian via :func:`pygmt.grdfilter`) of the - DiM-filtered data is generally recommended. + Thus, an additional filtering (e.g., Gaussian via :func:`pygmt.grdfilter`) + of the DiM-filtered data is generally recommended. Full option list at :gmt-docs:`dimfilter.html` From 7d045a40674badae340f1cb0ea07cb325c7a8ac3 Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Thu, 23 Dec 2021 08:47:07 +0000 Subject: [PATCH 10/22] remove mention of Q alias --- pygmt/src/dimfilter.py | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/pygmt/src/dimfilter.py b/pygmt/src/dimfilter.py index e69f8594741..9bad2676e83 100644 --- a/pygmt/src/dimfilter.py +++ b/pygmt/src/dimfilter.py @@ -38,9 +38,7 @@ def dimfilter(grid, **kwargs): subregion of the input and/or with a new **-I**\ ncrement. In this way, one may have "extra space" in the input data so that there will be no edge effects for the output grid. If the filter is low-pass, then the output - may be less frequently sampled than the input. The **-Q** option is for - the error analysis mode and expects the input file to contains the - filtered depths. Finally, one should know that **dimfilter** will not + may be less frequently sampled than the input. **dimfilter** will not produce a smooth output as other spatial filters do because it returns a minimum median out of *N* medians of *N* sectors. The output can be rough unless the input data is noise-free. From 1703238956d965e924f428b2b1f40cec9f5809d8 Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Thu, 23 Dec 2021 08:58:36 +0000 Subject: [PATCH 11/22] updating spacing docstring --- pygmt/src/dimfilter.py | 7 ++++++- 1 file changed, 6 insertions(+), 1 deletion(-) diff --git a/pygmt/src/dimfilter.py b/pygmt/src/dimfilter.py index 9bad2676e83..587706114bb 100644 --- a/pygmt/src/dimfilter.py +++ b/pygmt/src/dimfilter.py @@ -105,7 +105,12 @@ def dimfilter(grid, **kwargs): If more than one mode is found we return their average value. Append **+l** or **+h** to the sectors if you rather want to return the smallest or largest of the modal values. - {I} + spacing : str + *x_inc* [and optionally *y_inc*] is the output Increment. Append + **m** to indicate minutes, or **c** to indicate seconds. If the new + *x_inc*, *y_inc* are NOT integer multiples of the old ones (in the + input data), filtering will be considerably slower. [Default: Same + as input.] {R} {V} From b000b5577046002294c0a6dc7760867985c332a7 Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Thu, 23 Dec 2021 09:06:06 +0000 Subject: [PATCH 12/22] update region docstring --- pygmt/src/dimfilter.py | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/pygmt/src/dimfilter.py b/pygmt/src/dimfilter.py index 587706114bb..153dbe88f42 100644 --- a/pygmt/src/dimfilter.py +++ b/pygmt/src/dimfilter.py @@ -111,7 +111,9 @@ def dimfilter(grid, **kwargs): *x_inc*, *y_inc* are NOT integer multiples of the old ones (in the input data), filtering will be considerably slower. [Default: Same as input.] - {R} + region : str or list + [*west*, *east*, *south*, *north*]. + Defines the region of the output points. [Default: Same as input.] {V} Returns From f1788497810bcdb2663143779b92fe8e353f6cb5 Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Thu, 23 Dec 2021 09:22:14 +0000 Subject: [PATCH 13/22] update test_dimfilter.py to use xarray testing --- pygmt/tests/test_dimfilter.py | 49 +++++++++++++++++++++++++---------- 1 file changed, 36 insertions(+), 13 deletions(-) diff --git a/pygmt/tests/test_dimfilter.py b/pygmt/tests/test_dimfilter.py index 61418cd5478..c758e0b54b5 100644 --- a/pygmt/tests/test_dimfilter.py +++ b/pygmt/tests/test_dimfilter.py @@ -5,7 +5,8 @@ import numpy.testing as npt import pytest -from pygmt import dimfilter +import xarray as xr +from pygmt import dimfilter, load_dataarray from pygmt.datasets import load_earth_relief from pygmt.exceptions import GMTInvalidInput from pygmt.helpers import GMTTempFile @@ -16,10 +17,33 @@ def fixture_grid(): """ Load the grid data from the sample earth_relief file. """ - return load_earth_relief(resolution="30m", region=[-5, 5, -5, 5]) + return load_earth_relief( + resolution="01d", registration="pixel", region=[124, 130, -25, -20] + ) -def test_dimfilter_outgrid(grid): +@pytest.fixture(scope="module", name="expected_grid") +def fixture_grid_result(): + """ + Load the expected dimfilter grid result. + """ + return xr.DataArray( + data=[ + [397.0, 393.0, 377.5, 382.75, 419.0, 435.5], + [361.75, 359.0, 364.5, 415.5, 411.0, 373.5], + [321.5, 361.75, 373.5, 379.5, 369.0, 379.5], + [292.5, 324.0, 356.0, 361.5, 337.0, 361.5], + [306.0, 282.5, 310.25, 324.0, 353.0, 372.5], + ], + coords=dict( + lon=[124.5, 125.5, 126.5, 127.5, 128.5, 129.5], + lat=[-24.5, -23.5, -22.5, -21.5, -20.5], + ), + dims=["lat", "lon"], + ) + + +def test_dimfilter_outgrid(grid, expected_grid): """ Test the required parameters for dimfilter with a set outgrid. """ @@ -29,26 +53,25 @@ def test_dimfilter_outgrid(grid): ) assert result is None # return value is None assert os.path.exists(path=tmpfile.name) # check that outgrid exists + temp_grid = load_dataarray(tmpfile.name) + xr.testing.assert_allclose(a=temp_grid, b=expected_grid) -def test_grdgradient_no_outgrid(grid): +def test_grdgradient_no_outgrid(grid, expected_grid): """ Test the required parameters for dimfilter with no set outgrid. """ - temp_grid = dimfilter(grid=grid, filter="m600", distance=4, sectors="l6") - assert temp_grid.dims == ("lat", "lon") - assert temp_grid.gmt.gtype == 1 # Geographic grid - assert temp_grid.gmt.registration == 1 # Pixel registration - npt.assert_allclose(temp_grid.min(), -5125) - npt.assert_allclose(temp_grid.max(), -3750.5) - npt.assert_allclose(temp_grid.median(), -4826.5) - npt.assert_allclose(temp_grid.mean(), -4789.791) + result = dimfilter(grid=grid, filter="m600", distance=4, sectors="l6") + assert result.dims == ("lat", "lon") + assert result.gmt.gtype == 1 # Geographic grid + assert result.gmt.registration == 1 # Pixel registration + xr.testing.assert_allclose(a=result, b=expected_grid) def test_dimfilter_fails(grid): """ Check that dimfilter fails correctly when sector, filters, and distance are - not specified.. + not specified. """ with pytest.raises(GMTInvalidInput): dimfilter(grid=grid) From 4e7165ffcd7d937f2bd04bb31bb841fd3e154b8f Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Thu, 23 Dec 2021 09:25:37 +0000 Subject: [PATCH 14/22] docstring change --- pygmt/src/dimfilter.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/pygmt/src/dimfilter.py b/pygmt/src/dimfilter.py index 153dbe88f42..1923f0007fb 100644 --- a/pygmt/src/dimfilter.py +++ b/pygmt/src/dimfilter.py @@ -35,8 +35,8 @@ def dimfilter(grid, **kwargs): filters to each sector, and choosing the final outcome according to the selected secondary filter. It computes distances using Cartesian or Spherical geometries. The output grid can optionally be generated as a - subregion of the input and/or with a new **-I**\ ncrement. In this way, - one may have "extra space" in the input data so that there will be no edge + subregion of the input and/or with a new increment using ``spacing``, + which may add an "extra space" in the input data to prevent edge effects for the output grid. If the filter is low-pass, then the output may be less frequently sampled than the input. **dimfilter** will not produce a smooth output as other spatial filters From caf6daee6cf29520e5dfb227450433dfccdacea3 Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Thu, 23 Dec 2021 09:26:19 +0000 Subject: [PATCH 15/22] remove unused import in test_dimfilter.py --- pygmt/tests/test_dimfilter.py | 1 - 1 file changed, 1 deletion(-) diff --git a/pygmt/tests/test_dimfilter.py b/pygmt/tests/test_dimfilter.py index c758e0b54b5..06d8a0c7fdf 100644 --- a/pygmt/tests/test_dimfilter.py +++ b/pygmt/tests/test_dimfilter.py @@ -3,7 +3,6 @@ """ import os -import numpy.testing as npt import pytest import xarray as xr from pygmt import dimfilter, load_dataarray From 87c7cb60d47180cb08a46663b632ebea3c82dea3 Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Thu, 23 Dec 2021 22:07:34 +0000 Subject: [PATCH 16/22] change docstring --- pygmt/src/dimfilter.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pygmt/src/dimfilter.py b/pygmt/src/dimfilter.py index 1923f0007fb..2b4ba485601 100644 --- a/pygmt/src/dimfilter.py +++ b/pygmt/src/dimfilter.py @@ -1,5 +1,5 @@ """ -dimfilter - Filter a grid file by dividing the filter circle. +dimfilter - Directional filtering of grids in the space domain. """ from pygmt.clib import Session From 8cb34e062989d2d35dba22ec535bc45a191432a5 Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Sun, 27 Feb 2022 08:24:56 +0000 Subject: [PATCH 17/22] modify test_dimfilter.py to use static earth relief --- pygmt/tests/test_dimfilter.py | 33 +++++++++++++++++++-------------- 1 file changed, 19 insertions(+), 14 deletions(-) diff --git a/pygmt/tests/test_dimfilter.py b/pygmt/tests/test_dimfilter.py index 06d8a0c7fdf..8de90abbe26 100644 --- a/pygmt/tests/test_dimfilter.py +++ b/pygmt/tests/test_dimfilter.py @@ -6,19 +6,17 @@ import pytest import xarray as xr from pygmt import dimfilter, load_dataarray -from pygmt.datasets import load_earth_relief from pygmt.exceptions import GMTInvalidInput from pygmt.helpers import GMTTempFile +from pygmt.helpers.testing import load_static_earth_relief @pytest.fixture(scope="module", name="grid") def fixture_grid(): """ - Load the grid data from the sample earth_relief file. + Load the grid data from the static_earth_relief file. """ - return load_earth_relief( - resolution="01d", registration="pixel", region=[124, 130, -25, -20] - ) + return load_static_earth_relief() @pytest.fixture(scope="module", name="expected_grid") @@ -28,15 +26,15 @@ def fixture_grid_result(): """ return xr.DataArray( data=[ - [397.0, 393.0, 377.5, 382.75, 419.0, 435.5], - [361.75, 359.0, 364.5, 415.5, 411.0, 373.5], - [321.5, 361.75, 373.5, 379.5, 369.0, 379.5], - [292.5, 324.0, 356.0, 361.5, 337.0, 361.5], - [306.0, 282.5, 310.25, 324.0, 353.0, 372.5], + [346.0, 344.5, 349.0, 349.0], + [344.5, 318.5, 344.5, 394.0], + [344.5, 356.5, 345.5, 352.5], + [367.5, 349.0, 385.5, 349.0], + [435.0, 385.5, 413.5, 481.5], ], coords=dict( - lon=[124.5, 125.5, 126.5, 127.5, 128.5, 129.5], - lat=[-24.5, -23.5, -22.5, -21.5, -20.5], + lon=[-54.5, -53.5, -52.5, -51.5], + lat=[-23.5, -22.5, -21.5, -20.5, -19.5], ), dims=["lat", "lon"], ) @@ -48,7 +46,12 @@ def test_dimfilter_outgrid(grid, expected_grid): """ with GMTTempFile(suffix=".nc") as tmpfile: result = dimfilter( - grid=grid, outgrid=tmpfile.name, filter="m600", distance=4, sectors="l6" + grid=grid, + outgrid=tmpfile.name, + filter="m600", + distance=4, + sectors="l6", + region=[-55, -51, -24, -19], ) assert result is None # return value is None assert os.path.exists(path=tmpfile.name) # check that outgrid exists @@ -60,7 +63,9 @@ def test_grdgradient_no_outgrid(grid, expected_grid): """ Test the required parameters for dimfilter with no set outgrid. """ - result = dimfilter(grid=grid, filter="m600", distance=4, sectors="l6") + result = dimfilter( + grid=grid, filter="m600", distance=4, sectors="l6", region=[-55, -51, -24, -19] + ) assert result.dims == ("lat", "lon") assert result.gmt.gtype == 1 # Geographic grid assert result.gmt.registration == 1 # Pixel registration From 26faa9d94909593e9951e11013c2f71ee2405baa Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Wed, 9 Mar 2022 08:25:08 +0000 Subject: [PATCH 18/22] Apply suggestions from code review Co-authored-by: Dongdong Tian Co-authored-by: Wei Ji <23487320+weiji14@users.noreply.github.com> --- pygmt/src/dimfilter.py | 7 ++++--- pygmt/tests/test_dimfilter.py | 2 +- 2 files changed, 5 insertions(+), 4 deletions(-) diff --git a/pygmt/src/dimfilter.py b/pygmt/src/dimfilter.py index 2b4ba485601..160050a34a2 100644 --- a/pygmt/src/dimfilter.py +++ b/pygmt/src/dimfilter.py @@ -112,7 +112,7 @@ def dimfilter(grid, **kwargs): input data), filtering will be considerably slower. [Default: Same as input.] region : str or list - [*west*, *east*, *south*, *north*]. + [*xmin*, *xmax*, *ymin*, *ymax*]. Defines the region of the output points. [Default: Same as input.] {V} @@ -127,14 +127,15 @@ def dimfilter(grid, **kwargs): """ if ("D" not in kwargs) or ("F" not in kwargs) or ("N" not in kwargs): raise GMTInvalidInput( - """The following parameters are required: distance, filters, sectors""" + """At least one of the following parameters must be specified: + distance, filters, or sectors.""" ) 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 + if "G" not in kwargs: # if outgrid is unset, output to tempfile kwargs.update({"G": tmpfile.name}) outgrid = kwargs["G"] arg_str = " ".join([infile, build_arg_string(kwargs)]) diff --git a/pygmt/tests/test_dimfilter.py b/pygmt/tests/test_dimfilter.py index 8de90abbe26..72a2074ba14 100644 --- a/pygmt/tests/test_dimfilter.py +++ b/pygmt/tests/test_dimfilter.py @@ -59,7 +59,7 @@ def test_dimfilter_outgrid(grid, expected_grid): xr.testing.assert_allclose(a=temp_grid, b=expected_grid) -def test_grdgradient_no_outgrid(grid, expected_grid): +def test_dimfilter_no_outgrid(grid, expected_grid): """ Test the required parameters for dimfilter with no set outgrid. """ From fd04ddf0f45e1970df294cc9c316ceb628c67691 Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Fri, 11 Mar 2022 08:36:29 +0000 Subject: [PATCH 19/22] Apply suggestions from code review Co-authored-by: Wei Ji <23487320+weiji14@users.noreply.github.com> --- pygmt/src/dimfilter.py | 2 +- pygmt/tests/test_dimfilter.py | 6 +++--- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/pygmt/src/dimfilter.py b/pygmt/src/dimfilter.py index 160050a34a2..b9669866c19 100644 --- a/pygmt/src/dimfilter.py +++ b/pygmt/src/dimfilter.py @@ -133,7 +133,7 @@ def dimfilter(grid, **kwargs): with GMTTempFile(suffix=".nc") as tmpfile: with Session() as lib: - file_context = lib.virtualfile_from_data(check_kind="raster", data=grid) + file_context = lib.virtualfile_from_data(check_kind=None, data=grid) with file_context as infile: if "G" not in kwargs: # if outgrid is unset, output to tempfile kwargs.update({"G": tmpfile.name}) diff --git a/pygmt/tests/test_dimfilter.py b/pygmt/tests/test_dimfilter.py index 72a2074ba14..879a81d6b6c 100644 --- a/pygmt/tests/test_dimfilter.py +++ b/pygmt/tests/test_dimfilter.py @@ -74,8 +74,8 @@ def test_dimfilter_no_outgrid(grid, expected_grid): def test_dimfilter_fails(grid): """ - Check that dimfilter fails correctly when sector, filters, and distance are - not specified. + Check that dimfilter fails correctly when not all of sectors, filters, and + distance are specified. """ with pytest.raises(GMTInvalidInput): - dimfilter(grid=grid) + dimfilter(grid=grid, sectors="l6", distance=4) From bb33afd28002377fbb980af3e8ff4c4eed674ba3 Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Wed, 16 Mar 2022 08:13:19 +0000 Subject: [PATCH 20/22] Update pygmt/src/dimfilter.py Co-authored-by: Wei Ji <23487320+weiji14@users.noreply.github.com> --- pygmt/src/dimfilter.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pygmt/src/dimfilter.py b/pygmt/src/dimfilter.py index b9669866c19..2c412483697 100644 --- a/pygmt/src/dimfilter.py +++ b/pygmt/src/dimfilter.py @@ -125,7 +125,7 @@ def dimfilter(grid, **kwargs): - None if ``outgrid`` is set (grid output will be stored in file set by ``outgrid``) """ - if ("D" not in kwargs) or ("F" not in kwargs) or ("N" not in kwargs): + if not all(arg in kwargs for arg in ["D", "F", "N"]) and "Q" not in kwargs: raise GMTInvalidInput( """At least one of the following parameters must be specified: distance, filters, or sectors.""" From 0efbd83d0d7fbde20726d74094ab1f672578a6a3 Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Thu, 17 Mar 2022 09:46:44 +0000 Subject: [PATCH 21/22] Update pygmt/src/dimfilter.py Co-authored-by: Dongdong Tian --- pygmt/src/dimfilter.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pygmt/src/dimfilter.py b/pygmt/src/dimfilter.py index 2c412483697..7753abf8348 100644 --- a/pygmt/src/dimfilter.py +++ b/pygmt/src/dimfilter.py @@ -24,7 +24,7 @@ R="region", V="verbose", ) -@kwargs_to_strings(R="sequence") +@kwargs_to_strings(I="sequence", R="sequence") def dimfilter(grid, **kwargs): r""" Filter a grid by dividing the filter circle. From 388b521c9677313a4aa395a61ccd48f6e34378eb Mon Sep 17 00:00:00 2001 From: Will Schlitzer Date: Sat, 19 Mar 2022 08:21:22 -0700 Subject: [PATCH 22/22] Apply suggestions from code review Co-authored-by: Wei Ji <23487320+weiji14@users.noreply.github.com> --- pygmt/src/dimfilter.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/pygmt/src/dimfilter.py b/pygmt/src/dimfilter.py index 7753abf8348..72e7cd695fe 100644 --- a/pygmt/src/dimfilter.py +++ b/pygmt/src/dimfilter.py @@ -74,9 +74,9 @@ def dimfilter(grid, **kwargs): - **4**\ : grid (x,y) in degrees, *width* in km, Spherical distance calculation. filter : str - **x**\ *width*\ [**+l**\|\ **u**] + **x**\ *width*\ [**+l**\|\ **u**]. Sets the primary filter type. Choose among convolution and - non-convolution filters. Append the filter code **x** followed by + non-convolution filters. Use the filter code **x** followed by the full diameter *width*. Available convolution filters are: - (**b**) Boxcar: All weights are equal. @@ -105,7 +105,7 @@ def dimfilter(grid, **kwargs): If more than one mode is found we return their average value. Append **+l** or **+h** to the sectors if you rather want to return the smallest or largest of the modal values. - spacing : str + spacing : str or list *x_inc* [and optionally *y_inc*] is the output Increment. Append **m** to indicate minutes, or **c** to indicate seconds. If the new *x_inc*, *y_inc* are NOT integer multiples of the old ones (in the