search

GenericMappingTools / pygmt

A Python interface for the Generic Mapping Tools.
https://www.pygmt.org
BSD 3-Clause "New" or "Revised" License
758 stars 220 forks source link

Better return values for grdinfo #593

Open weiji14 opened 4 years ago

weiji14 commented 4 years ago

Description of the desired feature

Certain pygmt.grdinfo return values can be in an array form rather than a string. This would be a good first issue to take up, and would be similar to #575 but for grids instead of tables.

Originally posted by @liamtoney in https://github.com/GenericMappingTools/pygmt/issues/147#issuecomment-450018956

tl;dr

After reading the above discussion, I would like to work on more fully wrapping grdinfo. It seems like a potentially good "first issue," and I can see a benefit in wrapping it even with Python containers also available. See below, for example.

Motivation

In my use case, I want to easily see the min/max z-values of a DEM so that I can choose my colorbar limits.

1. No specified region

With GRIDFILE being the full path to my grid file, I could use

import gmt
print(gmt.grdinfo(GRIDFILE).strip())
...
../srtm15_plus_w180n90.nc: z_min: -7840.125 z_max: 6057.64355469 name: z
...

or (with xarray, for example)

import xarray
print(xarray.open_dataset(GRIDFILE).z.actual_range)
[-7840.125       6057.64355469]

These are about equal effort (and xarray returns things in a more usable output).

2. Specified region of interest

For viewing the min/maz z-values within a specified REGION = (lonmin, lonmax, latmin, latmax), I think that grdinfo makes this easier:

import gmt
print(gmt.grdinfo(GRIDFILE, R='{}/{}/{}/{}'.format(*REGION)).strip())
...
../srtm15_plus_w180n90.nc: z_min: -4996.11816406 z_max: 6057.64355469 name: z
...

vs.

import xarray
z = xarray.open_dataset(GRIDFILE).sel(lon=slice(*REGION[0:2]), lat=slice(*REGION[2:4])).z.values
print(z.min(), z.max())
-4996.118 6057.6436

Plus grdinfo is more format agnostic. By minimally including the @use_alias and @kwargs_to_strings decorators, the grdinfo command above could be just

print(gmt.grdinfo(GRIDFILE, region=REGION).strip())

and (with more work) could perhaps output a dictionary? Please let me know if this is a worthwhile effort and I'll make a PR for modules.py.

Are you willing to help implement and maintain this feature? Yes/No

maxrjones commented 3 years ago

This issue was mentioned in https://github.com/GenericMappingTools/pygmt/issues/1120, but it is not clear to me how this would be changed in a backwards compatible way. Does anyone have ideas for how to follow the deprecation policy here? Or should it be considered a bug where the deprecation policy may not apply?

willschlitzer commented 3 years ago

This issue was mentioned in #1120, but it is not clear to me how this would be changed in a backwards compatible way. Does anyone have ideas for how to follow the deprecation policy here? Or should it be considered a bug where the deprecation policy may not apply?

My guess for why this was brought up is that if the default option is return a dictionary, this will break functions that expect a string returned. Truth be told, I think the people most affected by this deprecation would be pygmt maintainers because of the number of tests that use grdinfo to check the answers.

My thought for the implementation of this is a parameter along the lines of info_type=dictionary that would return a dictionary by default, but passing string instead would get the original format returned.

weiji14 commented 3 years ago

Just wanted to surface a point I made in https://github.com/GenericMappingTools/pygmt/pull/1418#discussion_r700746488. Please try and refrain from using grdinfo while writing unit tests for new wrapper functions (otherwise we'll have lots of code to change if the pygmt.grdinfo implementation changes). Instead, use numpy.testing.assert_allclose or xarray.testing.assert_allclose (or similar) to test the output data array directly. E.g.

https://github.com/GenericMappingTools/pygmt/blob/2edaa63362feb58224b9e585ec3958c1e4df8b41/pygmt/tests/test_sphdistance.py#L43-L46

Try to keep this in mind when reviewing PRs or wrapping new modules. We should also maybe start migrating tests using grdinfo to check output values to use npt.assert_allclose instead before refactoring pygmt.grdinfo itself.

maxrjones commented 3 years ago

I would be pretty excited to see this in the next release, so I will get started on refactoring the tests. Here is a list of modules in pygmt/tests/ that contain grdinfo. If anyone wants to help, please edit this comment to add your username next to the module or post a comment below to avoid redundant work.

weiji14 commented 3 years ago

I would be pretty excited to see this in the next release, so I will get started on refactoring the tests. Here is a list of modules in pygmt/tests/ that contain grdinfo. If anyone wants to help, please edit this comment to add your username next to the module or post a comment below to avoid redundant work.

Thanks Meghan for making the list! I'm just wondering though, since the SRTM+v2.3 grids will be coming soon-ish (https://github.com/GenericMappingTools/gmtserver-admin/issues/105), if we should refactor the grd* tests now, or wait until the new @earth_relief grids are available. Just want to try to cut down on the amount of maintenance work we need to do :slightly_smiling_face:

Oh, and by next release, do you mean PyGMT v0.5.0 or v0.6.0?

maxrjones commented 3 years ago

I would be pretty excited to see this in the next release, so I will get started on refactoring the tests. Here is a list of modules in pygmt/tests/ that contain grdinfo. If anyone wants to help, please edit this comment to add your username next to the module or post a comment below to avoid redundant work.

Thanks Meghan for making the list! I'm just wondering though, since the SRTM+v2.3 grids will be coming soon-ish (GenericMappingTools/gmtserver-admin#105), if we should refactor the grd* tests now, or wait until the new @earth_relief grids are available. Just want to try to cut down on the amount of maintenance work we need to do 🙂

Oh, and by next release, do you mean PyGMT v0.5.0 or v0.6.0?

I was (probably naively) hoping for v0.5.0. But since we have a bit of a backlog in PRs and not much time until the release, we could hold off until v0.6.0 and after the remaking of earth_relief grids.