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

weiji14 · 2023-08-25T06:38:53Z

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

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 · 2023-08-25T06:44:02Z

examples/gallery/images/rgb_image.py

+fig.grdimage(
+    grid=image,
+    projection="x1:100000",
+    frame=[r"WSne+tL@!a\225hain@!a\225, Hawai\047i on 9 Aug 2023", "af"],


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

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?

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.

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

weiji14 · 2023-08-25T06:58:27Z

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.

examples/gallery/images/rgb_image.py

Co-Authored-By: Michael Grund <[email protected]>

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.

weiji14 · 2023-08-26T02:17:00Z

examples/gallery/images/rgb_image.py

+    fig.grdimage(
+        grid=image,
+        projection="x1:100000",
+        frame=[r"WSne+tL@!a\225hain@!a\225, Hawai\140i on 9 Aug 2023", "af"],


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:

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:

The 1:100000 scale means 1 centimetre on the map is equivalent to 1 kilometre on the ground.

examples/gallery/images/rgb_image.py

seisman · 2023-08-26T05:30:59Z

examples/gallery/images/rgb_image.py

+###############################################################################
+# 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",


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

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).

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.

weiji14 · 2023-08-30T03:45:19Z

examples/gallery/images/rgb_image.py

+) 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


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.

Yep, it works! See https://pygmt-dev--2641.org.readthedocs.build/en/2641/gallery/images/rgb_image.html

seisman

Looks good.

examples/gallery/images/rgb_image.py

yvonnefroehlich · 2023-08-30T17:19:34Z

examples/gallery/images/rgb_image.py

+---------
+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


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.

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

Ah, great. Thanks!

Co-authored-by: Yvonne Fröhlich <[email protected]>

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

9d38316

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 added the documentation Improvements or additions to documentation label Aug 25, 2023

weiji14 added this to the 0.10.0 milestone Aug 25, 2023

weiji14 self-assigned this Aug 25, 2023

weiji14 commented Aug 25, 2023

View reviewed changes

weiji14 mentioned this pull request Aug 25, 2023

Support non-ASCII characters in PyGMT arguments and text in Figure.text #2204

Closed

michaelgrund reviewed Aug 25, 2023

View reviewed changes

examples/gallery/images/rgb_image.py Outdated Show resolved Hide resolved

yvonnefroehlich reviewed Aug 25, 2023

View reviewed changes

examples/gallery/images/rgb_image.py Outdated Show resolved Hide resolved

yvonnefroehlich mentioned this pull request Aug 25, 2023

Expand documentation for Cartesian projections #2643

Merged

7 tasks

weiji14 and others added 3 commits August 26, 2023 13:15

Add abbreviation for Cloud-Optimized GeoTIFF

60d1be3

Co-Authored-By: Michael Grund <[email protected]>

Set title font to Times-Roman

b516c19

Makes the unequal spacing before and after the ā less obvious.

Change apostrophe to Okina symbol ʻ

e7b66e6

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.

weiji14 commented Aug 26, 2023

View reviewed changes

Inline comment on map scale projection

9085f65

The 1:100000 scale means 1 centimetre on the map is equivalent to 1 kilometre on the ground.

seisman reviewed Aug 26, 2023

View reviewed changes

examples/gallery/images/rgb_image.py Outdated Show resolved Hide resolved

seisman reviewed Aug 26, 2023

View reviewed changes

weiji14 and others added 3 commits August 26, 2023 20:11

Use ¯ and ` instead of \225 and \140

4fd90be

The non-octal code versions are much easier to see. Co-Authored-By: Dongdong Tian <[email protected]>

Merge branch 'main' into gallery/rgb_image

7586ef7

Force loading dataarray into memory

3c0712d

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.

weiji14 commented Aug 30, 2023

View reviewed changes

weiji14 added the needs review This PR has higher priority and needs review. label Aug 30, 2023

weiji14 marked this pull request as ready for review August 30, 2023 03:53

weiji14 mentioned this pull request Aug 28, 2023

Release PyGMT v0.10.0 #2640

Closed

35 tasks

seisman approved these changes Aug 30, 2023

View reviewed changes

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

michaelgrund approved these changes Aug 30, 2023

View reviewed changes

yvonnefroehlich reviewed Aug 30, 2023

View reviewed changes

Apply suggestions from code review

00b427b

Co-authored-by: Yvonne Fröhlich <[email protected]>

Lint

2553cce

weiji14 merged commit 26d4fed into main Aug 30, 2023

weiji14 deleted the gallery/rgb_image branch August 30, 2023 22:33

weiji14 removed the final review call This PR requires final review and approval from a second reviewer label Aug 30, 2023

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

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

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

weiji14 commented Aug 25, 2023 •

edited

Loading

weiji14 Aug 25, 2023

michaelgrund Aug 25, 2023

yvonnefroehlich Aug 25, 2023

weiji14 Aug 26, 2023

weiji14 commented Aug 25, 2023

weiji14 Aug 26, 2023

seisman Aug 26, 2023

weiji14 Aug 26, 2023

weiji14 Aug 30, 2023

weiji14 Aug 30, 2023

seisman left a comment

yvonnefroehlich Aug 30, 2023

weiji14 Aug 30, 2023

yvonnefroehlich Aug 31, 2023

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

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

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

Conversation

weiji14 commented Aug 25, 2023 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

weiji14 commented Aug 25, 2023

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

seisman left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

weiji14 commented Aug 25, 2023 •

edited

Loading