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

Migrate xarray examples #676

Open
dcherian opened this issue Jun 6, 2023 · 0 comments
Open

Migrate xarray examples #676

dcherian opened this issue Jun 6, 2023 · 0 comments
Labels
documentation Documentation related issue proposal Idea for a new feature.

Comments

@dcherian
Copy link

dcherian commented Jun 6, 2023

From the visualization gallery

The data is here


imshow() and rasterio map projections

Using rasterio's projection information for more accurate plots.

This example extends recipes.rasterio and plots the image in the
original map projection instead of relying on pcolormesh and a map
transformation.

# opening raster files requires the rioxarray package
da = xr.open_rasterio("RGB.byte")

# The data is in UTM projection. We have to set it manually until
# https://github.com/SciTools/cartopy/issues/813 is implemented
crs = ccrs.UTM("18")

# Plot on a map
ax = plt.subplot(projection=crs)
da.plot.imshow(ax=ax, rgb="band", transform=crs)
ax.coastlines("10m", color="r")

Parsing rasterio geocoordinates

Converting a projection's cartesian coordinates into 2D longitudes and
latitudes.

These new coordinates might be handy for plotting and indexing, but it should
be kept in mind that a grid which is regular in projection coordinates will
likely be irregular in lon/lat. It is often recommended to work in the data's
original map projection (see recipes.rasterio_rgb).

from pyproj import Transformer
import numpy as np

# opening raster files requires the rioxarray package
da = xr.tutorial.open_rasterio("RGB.byte")

x, y = np.meshgrid(da["x"], da["y"])
transformer = Transformer.from_crs(da.crs, "EPSG:4326", always_xy=True)
lon, lat = transformer.transform(x, y)
da.coords["lon"] = (("y", "x"), lon)
da.coords["lat"] = (("y", "x"), lat)

# Compute a greyscale out of the rgb image
greyscale = da.mean(dim="band")

# Plot on a map
ax = plt.subplot(projection=ccrs.PlateCarree())
greyscale.plot(
    ax=ax,
    x="lon",
    y="lat",
    transform=ccrs.PlateCarree(),
    cmap="Greys_r",
    shading="auto",
    add_colorbar=False,
)
ax.coastlines("10m", color="r")
@dcherian dcherian added the proposal Idea for a new feature. label Jun 6, 2023
@snowman2 snowman2 added the documentation Documentation related issue label Jun 6, 2023
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
documentation Documentation related issue proposal Idea for a new feature.
Projects
None yet
Development

No branches or pull requests

2 participants