Initialize a GMTDataArrayAccessor

.github/workflows/ci_tests.yaml

@@ -67,7 +67,7 @@ jobs:
           path: |
             ~/.gmt/cache
             ~/.gmt/server
-          key: cache-gmt-${{ runner.os }}-${{ github.ref }}-20200710-2
+          key: cache-gmt-${{ runner.os }}-${{ github.ref }}-20200711
           restore-keys: cache-gmt-${{ runner.os }}-refs/heads/master-
 
       # Workaround for the timeouts of 'gmt which' on Linux and Windows
@@ -80,7 +80,7 @@ jobs:
           for data in earth_relief_01d_p.grd earth_relief_01d_g.grd earth_relief_30m_p.grd earth_relief_30m_g.grd earth_relief_10m_p.grd earth_relief_10m_g.grd; do
             wget --no-check-certificate https://oceania.generic-mapping-tools.org/server/earth/earth_relief/${data} -P ~/.gmt/server/earth/earth_relief/
           done
-          for data in ridge.txt Table_5_11.txt tut_bathy.nc tut_quakes.ngdc tut_ship.xyz usgs_quakes_22.txt; do
+          for data in ridge.txt Table_5_11.txt test.dat.nc tut_bathy.nc tut_quakes.ngdc tut_ship.xyz usgs_quakes_22.txt; do
             wget --no-check-certificate https://oceania.generic-mapping-tools.org/cache/${data} -P ~/.gmt/cache/
           done
         if: steps.cache.outputs.cache-hit != 'true' && runner.os != 'macOS'
@@ -90,7 +90,7 @@ jobs:
         shell: bash -l {0}
         run: |
           gmt which -Ga @earth_relief_10m_p @earth_relief_10m_g @earth_relief_30m_p @earth_relief_30m_g @earth_relief_01d_p @earth_relief_01d_g
-          gmt which -Ga @ridge.txt @Table_5_11.txt @tut_bathy.nc @tut_quakes.ngdc @tut_ship.xyz @usgs_quakes_22.txt
+          gmt which -Ga @ridge.txt @Table_5_11.txt @test.dat.nc @tut_bathy.nc @tut_quakes.ngdc @tut_ship.xyz @usgs_quakes_22.txt
         if: steps.cache.outputs.cache-hit != 'true' && runner.os == 'macOS'
 
       # Install the package that we want to test

doc/api/index.rst

@@ -62,7 +62,6 @@ Operations on tabular data:
     :toctree: generated
 
     blockmedian
-    info
     surface
 
 Operations on grids:
@@ -71,7 +70,6 @@ Operations on grids:
     :toctree: generated
 
     grdcut
-    grdinfo
     grdtrack
 
 GMT Defaults
@@ -84,6 +82,19 @@ Operations on GMT defaults:
 
     config
 
+Metadata
+--------
+
+Getting metadata from tabular or grid data:
+
+.. autosummary::
+    :toctree: generated
+
+    GMTDataArrayAccessor
+    info
+    grdinfo
+
+
 Miscellaneous
 -------------
 

pygmt/__init__.py

@@ -18,7 +18,7 @@
 from .gridding import surface
 from .sampling import grdtrack
 from .mathops import makecpt
-from .modules import config, info, grdinfo, which
+from .modules import GMTDataArrayAccessor, config, info, grdinfo, which
 from .gridops import grdcut
 from . import datasets
 

pygmt/modules.py

@@ -1,6 +1,8 @@
 """
 Non-plot GMT modules.
 """
+import xarray as xr
+
 from .clib import Session
 from .helpers import (
     build_arg_string,
@@ -214,3 +216,71 @@ def __exit__(self, exc_type, exc_value, traceback):
         )
         with Session() as lib:
             lib.call_module("set", arg_str)
+
+
+@xr.register_dataarray_accessor("gmt")
+class GMTDataArrayAccessor:
+    """
+    This is the GMT extension for :class:`xarray.DataArray`.
+
+    You can access various GMT specific metadata about your grid as follows:
+
+    >>> from pygmt.datasets import load_earth_relief
+    >>> # Use the global Earth relief grid with 1 degree spacing
+    >>> grid = load_earth_relief(resolution="01d")
+
+    >>> # See if grid uses Gridline (0) or Pixel (1) registration
+    >>> grid.gmt.registration
+    1
+    >>> # See if grid uses Cartesian (0) or Geographic (1) coordinate system
+    >>> grid.gmt.gtype
+    1
+    """
+
+    def __init__(self, xarray_obj):
+        self._obj = xarray_obj
+        try:
+            self._source = self._obj.encoding["source"]  # filepath to NetCDF source
+            # From the shortened summary information of `grdinfo`,
+            # get grid registration in column 10, and grid type in column 11
+            self._registration, self._gtype = map(
+                int, grdinfo(self._source, C="n", o="10,11").split()
+            )
+        except KeyError:
+            self._registration = 0  # Default to Gridline registration
+            self._gtype = 0  # Default to Cartesian grid type
+
+    @property
+    def registration(self):
+        """
+        Registration type of the grid, either Gridline (0) or Pixel (1).
+        """
+        return self._registration
+
+    @registration.setter
+    def registration(self, value):
+        if value in (0, 1):
+            self._registration = value
+        else:
+            raise GMTInvalidInput(
+                f"Invalid grid registration value: {value}, should be a boolean of "
+                "either 0 for Gridline registration or 1 for Pixel registration"
+            )
+
+    @property
+    def gtype(self):
+        """
+        Coordinate system type of the grid, either Cartesian (0) or Geographic
+        (1).
+        """
+        return self._gtype
+
+    @gtype.setter
+    def gtype(self, value):
+        if value in (0, 1):
+            self._gtype = value
+        else:
+            raise GMTInvalidInput(
+                f"Invalid coordinate system type: {value}, should be a boolean of "
+                "either 0 for Cartesian or 1 for Geographic"
+            )

pygmt/tests/test_accessor.py

@@ -0,0 +1,68 @@
+"""
+Test the behaviour of the GMTDataArrayAccessor class
+"""
+import pytest
+import xarray as xr
+
+from .. import which
+from ..exceptions import GMTInvalidInput
+
+
+def test_accessor_gridline_cartesian():
+    """
+    Check that a grid returns a registration value of 0 when Gridline
+    registered, and a gtype value of 1 when using Geographic coordinates.
+    """
+    fname = which(fname="@test.dat.nc", download="a")
+    grid = xr.open_dataarray(fname)
+    assert grid.gmt.registration == 0  # gridline registration
+    assert grid.gmt.gtype == 0  # cartesian coordinate type
+
+
+def test_accessor_pixel_geographic():
+    """
+    Check that a grid returns a registration value of 1 when Pixel registered,
+    and a gtype value of 0 when using Cartesian coordinates.
+    """
+    fname = which(fname="@earth_relief_01d_p", download="a")
+    grid = xr.open_dataarray(fname)
+    assert grid.gmt.registration == 1  # pixel registration
+    assert grid.gmt.gtype == 1  # geographic coordinate type
+
+
+def test_accessor_set_pixel_registration():
+    """
+    Check that we can set a grid to be Pixel registered with a registration
+    value of 1.
+    """
+    grid = xr.DataArray(data=[[0.1, 0.2], [0.3, 0.4]])
+    assert grid.gmt.registration == 0  # default to gridline registration
+    grid.gmt.registration = 1  # set to pixel registration
+    assert grid.gmt.registration == 1  # ensure changed to pixel registration
+
+
+def test_accessor_set_geographic_cartesian_roundtrip():
+    """
+    Check that we can set a grid to switch between the default Cartesian
+    coordinate type using a gtype of 1, set it to Geographic 0, and then back
+    to Cartesian again 1.
+    """
+    grid = xr.DataArray(data=[[0.1, 0.2], [0.3, 0.4]])
+    assert grid.gmt.gtype == 0  # default to cartesian coordinate type
+    grid.gmt.gtype = 1  # set to geographic type
+    assert grid.gmt.gtype == 1  # ensure changed to geographic coordinate type
+    grid.gmt.gtype = 0  # set back to cartesian type
+    assert grid.gmt.gtype == 0  # ensure changed to cartesian coordinate type
+
+
+def test_accessor_set_non_boolean():
+    """
+    Check that setting non boolean values on registration and gtype do not work
+    """
+    grid = xr.DataArray(data=[[0.1, 0.2], [0.3, 0.4]])
+
+    with pytest.raises(GMTInvalidInput):
+        grid.gmt.registration = "2"
+
+    with pytest.raises(GMTInvalidInput):
+        grid.gmt.gtype = 2