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

Allow passing in Red, Green, Blue grids into grdimage #370

Closed
wants to merge 2 commits into from

Conversation

weiji14
Copy link
Member

@weiji14 weiji14 commented Nov 6, 2019

Description of proposed changes

Allow for multiple grid inputs into grdimage via a Python list. Using the sequence_space converter (created in #325) to do that. Edit: Using if-then parsing instead to handle a list of xarray.DataArray grids better. Useful when users want to make an RGB image plot e.g. from Landsat images (Bands 4,3,2).

Below is a bit of a convoluted example based on this tutorial (in Chinese) which has preprocessed grd files scaled to 0-255 range, all courtesy of @whyjz (Hi 👋!).

import pygmt

for band in [4, 3, 2]:
    pygmt.which(
        f"https://gmt-tutorials.org/en/_downloads/LC08_L1TP_160045_20170408_20170414_01_T1_B{band}_s.grd",
        download="l",
    )

fig = pygmt.Figure()
fig.grdimage(
    grid=[
        "LC08_L1TP_160045_20170408_20170414_01_T1_B4_s.grd",
        "LC08_L1TP_160045_20170408_20170414_01_T1_B3_s.grd",
        "LC08_L1TP_160045_20170408_20170414_01_T1_B2_s.grd",
    ],
    frame=["af", "wSnE"],
)
fig.show()

Rub’ al Khali desert image

Fixes #

Reminders

  • Run make format and make check to make sure the code follows the style guide.
  • Add tests for new features or tests that would have caught the bug that you're fixing.
  • Add new public functions/methods/classes to doc/api/index.rst.
  • Write detailed docstrings for all functions/methods.
  • If adding new functionality, add an example to docstrings or tutorials.

Copy link
Member Author

@weiji14 weiji14 left a comment

Choose a reason for hiding this comment

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

The implementation here only allows for file inputs at the moment. Wondering if there's a demand out there from our users to pass in xarray.DataArray RGB grids as well? It'll be a lot more work to implement (and might actually fit in a separate Pull Request), but just throwing the bone out there.
Edit: Updated implementation to work with list of xarray.DataArray grids as well.

def test_grdimage_rgb_files():
"Plot an image using Red, Green, and Blue file inputs"
fig = Figure()
fig.grdimage(grid=["@earth_relief_60m", "@earth_relief_60m", "@earth_relief_60m"])
Copy link
Member Author

@weiji14 weiji14 Nov 6, 2019

Choose a reason for hiding this comment

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

I don't quite like this unit test example, since they're clearly not RGB grids. Does anyone have a better sample dataset recommendation? Preferably a small one that has a fairly stable URL.

Copy link
Member

Choose a reason for hiding this comment

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

@weiji14 ideally, we should use a dataset that is open-access and published somewhere (figshare, zenodo, etc) or in the public domain. Maybe @djhoese can help us out here 🙂

Copy link

Choose a reason for hiding this comment

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

Hm I can only think of raw data at the moment. What format do you need as input (I'm not super familiar with the limitations of gmt)?

Copy link
Member Author

Choose a reason for hiding this comment

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

GeoTiffs would be nice. Or anything GDAL can read I think. I quite like the zenodo idea, something with a persistent doi. But happy with any other good public domain ones too!

Copy link

Choose a reason for hiding this comment

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

There has to be some landsat data or something. @snowman2 does a lot with rioxarray and may have an idea of a persistent/versioned RGB geotiff. I'll ask the Pytroll core developers too.

Choose a reason for hiding this comment

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

If I understand what you are asking, there is an RGB tif in the example here: https://corteva.github.io/rioxarray/html/examples/COG.html

Copy link
Member Author

Choose a reason for hiding this comment

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

I did look at the rioxarray example a while ago actually, but the test here requires tiffs of individual bands, e.g r.tif, g.tif, b.tif. Landsat on AWS does provide individual bands but they're quite big...

Copy link

Choose a reason for hiding this comment

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

Oh if you need the bands to be in separate tif files then that is going to be harder to track down.

Copy link

@snowman2 snowman2 Nov 25, 2019

Choose a reason for hiding this comment

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

You could just take a regular RGB file and write each band to its own tif file.

import rioxarray
rds = rioxarray,open_rasterio(....)
rdssub = rds.sel(band=1)
rdssub.rio.to_raster(....)

Copy link
Member Author

@weiji14 weiji14 Nov 25, 2019

Choose a reason for hiding this comment

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

Yeah, I thought I might need to create some dummy test files. @seisman just mentioned a similar idea with splitting a raster into R,G,B bands. Just seems like quite a bit of work for a unit test (requires set-up and tear down). What about daily/weather satellites e.g. MODIS/Himawari with a coarse resolution, anyone know where to get a sample?

Using the sequence_space converter to allow for multiple grid inputs into `grdimage`. Useful when users want to make an RGB image plot e.g. from Landsat images (Bands 4,3,2).
Extends `grdimage` to accept a list of individual red, green and blue xarray.DataArrays for plotting. Had to revert using the 'sequence_space' decorator as it would not work nicely for multiple grids... Implementation uses the contextlib.ExitStack method to create multiple virtual grid files and enter multiple contexts (similar to 3448972). Also added one test case with dummy data, but could do with more to capture edge cases.
@seisman
Copy link
Member

seisman commented Nov 25, 2019

FYI, there was a module grd2rgb in GMT5, but was removed from GMT6. Not sure if it helps with this issue.

@PaulWessel Is there an equivalent module or command in GMT6, which can convert a grid/image to rgb?

@weiji14 weiji14 added the enhancement Improving an existing feature label Jun 21, 2020
@weiji14 weiji14 mentioned this pull request Sep 3, 2020
@seisman
Copy link
Member

seisman commented Sep 4, 2020

FYI, it seems GMT will "disable" the feature to pass R, G, B grids into grdimage and grdview. Instead, users should use grdmix to convert RGB grids to a single image, and then pass it to grdimage/grdview.

See related discussions at GenericMappingTools/gmt#4044 (comment), GenericMappingTools/gmt#4044 (comment), GenericMappingTools/gmt#4163.

@weiji14
Copy link
Member Author

weiji14 commented Sep 4, 2020

Ok, I'll will put this on hold first until GMT 6.2.0. Probably should just close this, but I'll leave it open so that we can cherry-pick the tests from this branch into the grdmix wrapper later for solving #578.

@weiji14 weiji14 marked this pull request as draft September 4, 2020 21:19
@weiji14 weiji14 self-assigned this Mar 24, 2021
@weiji14 weiji14 added this to the 0.5.0 milestone Jun 3, 2021
@weiji14 weiji14 mentioned this pull request Aug 11, 2021
18 tasks
@weiji14 weiji14 removed this from the 0.5.0 milestone Oct 3, 2021
@seisman
Copy link
Member

seisman commented Apr 12, 2022

@weiji14 Do you think we should close this PR?

@weiji14
Copy link
Member Author

weiji14 commented Apr 12, 2022

@weiji14 Do you think we should close this PR?

Prefer to keep this open until #1437 is merged. Sorry that I haven't got around to this lately, been quite busy at uni...

@weiji14
Copy link
Member Author

weiji14 commented Jul 3, 2023

Since GMT doesn't support Red, Green and Blue grids separately into grdimage anymore as mentioned at #370 (comment), let's close this down then and only allow passing a single file/grid/dataarray into grdimage instead of a list of 3 things.

For those following, I've started a new implementation at #2590 to allow grdimage to plot xarray.DataArray objects of shape (3, Y, X) using a similar logic to what we did with fig.tilemap at #2394, and this also matches what other libraries like rioxarray are doing. The implementation of grdmix at #1437 should also be able to tie into that later.

@weiji14 weiji14 closed this Jul 3, 2023
@weiji14 weiji14 deleted the grdimage/rgb_grids branch July 3, 2023 06:17
@PaulWessel
Copy link
Member

Hi @weiji14 and @seisman and @maxrjones: Could you please send me a list of major activities for the pyGMT team since last July? I am writing the annual NSF report and would at least like to mention any workshops or presentations you may have had. Given time-zones, please send me something this morning (if you can).

@seisman
Copy link
Member

seisman commented Jul 3, 2023

Hi @weiji14 and @seisman and @maxrjones: Could you please send me a list of major activities for the pyGMT team since last July? I am writing the annual NSF report and would at least like to mention any workshops or presentations you may have had. Given time-zones, please send me something this morning (if you can).

@maxrjones gave a SciPy Talk on July, 13, 2022. The related material is available at https://github.com/maxrjones/scipy2022.

@weiji14
Copy link
Member Author

weiji14 commented Jul 3, 2023

Hi @weiji14 and @seisman and @maxrjones: Could you please send me a list of major activities for the pyGMT team since last July? I am writing the annual NSF report and would at least like to mention any workshops or presentations you may have had. Given time-zones, please send me something this morning (if you can).

Hmm, don't think we had any workshops for PyGMT since Jul 2022, the last one was the EGU 2022 workshop in May 2022. And yes, Max's SciPy talk is at https://www.youtube.com/watch?v=nCktihu9bWg

@PaulWessel
Copy link
Member

Great, thanks!

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
enhancement Improving an existing feature
Projects
None yet
Development

Successfully merging this pull request may close these issues.

6 participants