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

Add gallery example for plotting an RGB image from an xarray.DataArray #2641

Merged
merged 10 commits into from
Aug 30, 2023

Conversation

weiji14
Copy link
Member

@weiji14 weiji14 commented Aug 25, 2023

Description of proposed changes

Gallery example using pygmt.Figure.grdimage to plot an RGB image from a 3-band GeoTIFF loaded into an xarray.DataArray via rioxarray.open_rasterio. Example is over Lāhainā, Hawai'i on 9 Aug 2023.

Preview at https://pygmt-dev--2641.org.readthedocs.build/en/2641/gallery/images/rgb_image.html

Worldview 2 satellite image over Lāhainā, Hawai'i on 9 Aug 2023.

Image is sourced from https://map.openaerialmap.org/#/-156.64478302001953,20.89666397351934,11/square/022300003203/64d6ae6619cb3a000147a65f?resolution=high&_k=nh7cia

Followup from #2590.

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 wrapping a new module, open a 'Wrap new GMT module' issue and submit reasonably-sized PRs.
  • If adding new functionality, add an example to docstrings or tutorials.
  • Use underscores (not hyphens) in names of Python files and directories.

Slash Commands

You can write slash commands (/command) in the first line of a comment to perform
specific operations. Supported slash commands are:

  • /format: automatically format and lint the code
  • /test-gmt-dev: run full tests on the latest GMT development version

Gallery example using pygmt.Figure.grdimage to plot an RGB image from a 3-band GeoTIFF loaded into an xarray.DataArray via rioxarray.open_rasterio. Example is over Lāhainā, Hawai'i on 9 Aug 2023.
@weiji14 weiji14 added the documentation Improvements or additions to documentation label Aug 25, 2023
@weiji14 weiji14 added this to the 0.10.0 milestone Aug 25, 2023
@weiji14 weiji14 self-assigned this Aug 25, 2023
fig.grdimage(
grid=image,
projection="x1:100000",
frame=[r"WSne+tL@!a\225hain@!a\225, Hawai\047i on 9 Aug 2023", "af"],
Copy link
Member Author

Choose a reason for hiding this comment

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

Had to use the composite character @!a\225 to plot ā following https://docs.generic-mapping-tools.org/6.4/tutorial/session-2.html#plotting-text-strings, because the character ā is in ISO-8859-4, see https://en.wikipedia.org/wiki/%C4%80

Copy link
Member

Choose a reason for hiding this comment

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

I think that's fine. Maybe it's possible to avoid that the letters are too close to each other by using another font style?

Copy link
Member

Choose a reason for hiding this comment

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

Agree with @michaelgrund, these unequally spaced letters are a bit strange. Looking at the example at https://docs.generic-mapping-tools.org/6.4/tutorial/session-2.html#gmt-tut-10, using the font Times-Roman does not lead to unequally spaced letters.

Copy link
Member Author

Choose a reason for hiding this comment

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

Agree that the spacing looks weird. I've changed the font to Times-Roman at b516c19 and it looks a bit better now:

Lāhainā in Times-Roman font

@weiji14
Copy link
Member Author

weiji14 commented Aug 25, 2023

Note that the docs build for this example segfaults on both Readthedocs and Ubuntu, see e.g. https://github.com/GenericMappingTools/pygmt/actions/runs/5972759386/job/16203804119?pr=2641#step:7:82. Running this locally, I actually get a segfault too on the first run, but running cd doc && make all a second time works. Not quite sure what the problem is.

weiji14 and others added 3 commits August 26, 2023 13:15
Makes the unequal spacing before and after the ā less obvious.
See https://en.wikipedia.org/wiki/%CA%BBOkina. Using octal code 140 instead of 047 on the plot title for Hawaiʻi, so that it looks like a 6 instead of a 9, but unsure if this is still the correct Okina symbol.
fig.grdimage(
grid=image,
projection="x1:100000",
frame=[r"WSne+tL@!a\225hain@!a\225, Hawai\140i on 9 Aug 2023", "af"],
Copy link
Member Author

Choose a reason for hiding this comment

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

Changed to using \140 instead of \047 as it looks a bit closer to the Okina ʻ letter (U+02BB) (see also #1474 (comment)). Now it looks like this:

Lāhainā, Hawaiʻi in Times-Roman font.

If anyone has ideas on how to use the actual Okina letter, let me know! I did try using ʻ directly, but it shows up like this:

Hawaiʻi with the Okina letter plotted incorrectly.

The 1:100000 scale means 1 centimetre on the map is equivalent to 1 kilometre on the ground.
###############################################################################
# Read 3-band data from GeoTIFF into an xarray.DataArray object
image = rioxarray.open_rasterio(
filename="https://oin-hotosm.s3.us-east-1.amazonaws.com/64d6a49a19cb3a000147a65b/0/64d6a49a19cb3a000147a65c.tif",
Copy link
Member

Choose a reason for hiding this comment

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

This GeoTiff file is almost 1 GB, which is too big to download.

Copy link
Member Author

Choose a reason for hiding this comment

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

The next line overview_level=5 gets a reduced resolution version. Original image has shape (y: 115712, x: 99328), but at overview_level=5, the shape is (y: 1808, x: 1552).

weiji14 and others added 3 commits August 26, 2023 20:11
The non-octal code versions are much easier to see.

Co-Authored-By: Dongdong Tian <[email protected]>
Opening the GeoTIFF in a with statement, and loading the RGB image data fully into memory using `image.load()` as suggested in corteva/rioxarray#550 (comment), which should fix the segfault hopefully.
) as img:
# Subset to area of Lāhainā in EPSG:32604 coordinates
image = img.rio.clip_box(minx=738000, maxx=755000, miny=2300000, maxy=2318000)
image = image.load() # force loading dataarray into memory
Copy link
Member Author

Choose a reason for hiding this comment

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

Note that the docs build for this example segfaults on both Readthedocs and Ubuntu, see e.g. https://github.com/GenericMappingTools/pygmt/actions/runs/5972759386/job/16203804119?pr=2641#step:7:82. Running this locally, I actually get a segfault too on the first run, but running cd doc && make all a second time works. Not quite sure what the problem is.

Running .load() should fix the segfault? See corteva/rioxarray#550 (comment), and we actually use dataarray.load() in the pygmt.load_dataarray() function in #1439.

Copy link
Member Author

Choose a reason for hiding this comment

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

@weiji14 weiji14 added the needs review This PR has higher priority and needs review. label Aug 30, 2023
@weiji14 weiji14 marked this pull request as ready for review August 30, 2023 03:53
@weiji14 weiji14 mentioned this pull request Aug 28, 2023
35 tasks
Copy link
Member

@seisman seisman left a comment

Choose a reason for hiding this comment

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

Looks good.

@seisman seisman added final review call This PR requires final review and approval from a second reviewer and removed needs review This PR has higher priority and needs review. labels Aug 30, 2023
examples/gallery/images/rgb_image.py Outdated Show resolved Hide resolved
examples/gallery/images/rgb_image.py Outdated Show resolved Hide resolved
examples/gallery/images/rgb_image.py Outdated Show resolved Hide resolved
examples/gallery/images/rgb_image.py Outdated Show resolved Hide resolved
---------
The :meth:`pygmt.Figure.grdimage` method can be used to plot Red, Green, Blue
(RGB) images, or any 3-band false color combination. Here, we'll use
:meth:`rioxarray.open_rasterio` to read a GeoTIFF file into an
Copy link
Member

Choose a reason for hiding this comment

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

Is it possible to make this work as a link to https://corteva.github.io/rioxarray/stable/rioxarray.html#rioxarray.open_rasterio?

If this is not easy to do, then we leave it as is. At least in the code example below, the link works.

Copy link
Member Author

Choose a reason for hiding this comment

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

Ah yes, good catch. I used the wrong intersphinx directive, should be :py:func:

Suggested change
:meth:`rioxarray.open_rasterio` to read a GeoTIFF file into an
:py:func:`rioxarray.open_rasterio` to read a GeoTIFF file into an

Copy link
Member

Choose a reason for hiding this comment

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

Ah, great. Thanks!

Co-authored-by: Yvonne Fröhlich <[email protected]>
@weiji14 weiji14 merged commit 26d4fed into main Aug 30, 2023
8 checks passed
@weiji14 weiji14 deleted the gallery/rgb_image branch August 30, 2023 22:33
@weiji14 weiji14 removed the final review call This PR requires final review and approval from a second reviewer label Aug 30, 2023
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
documentation Improvements or additions to documentation
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

4 participants