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

Use the Manaus data in the Composites tutorial #52

Merged
merged 1 commit into from
Sep 21, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
7 changes: 4 additions & 3 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -25,9 +25,10 @@ manipulation, processing, and visualization.

## Example

Here's a quick example of loading this
Here's a quick example of loading and plotting this
[Landsat 9 scene from the city of Manaus, Brazil](https://doi.org/10.6084/m9.figshare.24167235.v1),
which is where the Amazonas and Negro river merge:
which is where the Solimões (brown water) and Negro (black water) rivers merge
to form the Amazon river:

```python
import xlandsat as xls
Expand All @@ -51,7 +52,7 @@ plt.axis("scaled")
plt.show()
```

<img src="https://raw.githubusercontent.com/compgeolab/xlandsat/main/doc/_static/readme-example.jpg" alt="RGB image showing the city on the left and the black waters of the Negro river merging with the brown waters of the Amazonas river" width="100%">
<img src="https://raw.githubusercontent.com/compgeolab/xlandsat/main/doc/_static/readme-example.jpg" alt="RGB image showing the city on the left and the black waters of the Negro river merging with the brown waters of the Solimões river" width="100%">

## Project goals

Expand Down
126 changes: 89 additions & 37 deletions doc/composites.rst
Original file line number Diff line number Diff line change
@@ -1,24 +1,28 @@
.. _composites:

Composites
==========
Making composites
=================

Plotting individual bands is good but we usually want to make some composite
images to visualize information from multiple bands at once.
For that, we have to create **composites**.
We provide the :func:`xlandsat.composite` function to make this process easier.
Plotting individual bands is nice but we usually want to make some composite
images, both RGB and false-color, to visualize information from multiple bands
at once.
We provide the :func:`xlandsat.composite` function to make this process easier
and produce composites that are compatible with both :mod:`xarray` and
:mod:`matplotlib`.

As an example, let's load an example scene from the
`Brumadinho tailings dam disaster <https://en.wikipedia.org/wiki/Brumadinho_dam_disaster>`__:
As an example, let's load our example scene from Manaus, Brazil, which is where
the Solimões (brown water) and Negro (black water) rivers merge to form the
Amazon river:

.. jupyter-execute::

import xlandsat as xls
import matplotlib.pyplot as plt

path = xls.datasets.fetch_brumadinho_after()

path = xls.datasets.fetch_manaus()
scene = xls.load_scene(path)
scene


RGB composites
Expand All @@ -34,10 +38,12 @@ bands aren't specified.
rgb = xls.composite(scene)
rgb

The composite has a similar layout as the bands of a scene but with an extra
``"channel"`` dimension corresponding to red, green, blue, and
alpha/transparency. The values are scaled to the [0, 255] range and the
composite is an array of unsigned 8-bit integers.
The composite is also an :class:`xarray.DataArray` and is similar to the bands.
It has the same easting and northing dimensions but also an extra ``"channel"``
dimension, which corresponds to red, green, blue, and alpha/transparency. This
extra dimension is what combines the three bands into a single color image. The
values are scaled to the [0, 255] range and the composite is an array of
unsigned 8-bit integers.

.. admonition:: Transparency
:class: note
Expand All @@ -48,11 +54,21 @@ composite is an array of unsigned 8-bit integers.
blue).


Plotting a composite
--------------------
Plotting composites
-------------------

Composites can be plotted using :meth:`xarray.DataArray.plot.imshow` (using
:meth:`~xarray.DataArray.plot` won't work and will display histograms instead).
Composites can be plotted using the :meth:`xarray.DataArray.plot.imshow` method:

.. jupyter-execute::

rgb.plot.imshow()


With no arguments, xarray takes care of creating the new figure and adding a
lot of the different plot elements, like axis labels and units.
If we want more control over the plot, we must create a matplotlib figure and
axes separately and tell :meth:`~xarray.DataArray.plot.imshow` to plot on those
instead:

.. jupyter-execute::

Expand All @@ -61,17 +77,23 @@ Composites can be plotted using :meth:`xarray.DataArray.plot.imshow` (using
rgb.plot.imshow(ax=ax)

# The "long_name" of the composite is the band combination
ax.set_title(f"Composite: {rgb.attrs['long_name']}")
ax.set_title(rgb.attrs["long_name"].title())

# Make sure pixels are square and don't have any distortions from plotting
# Make sure pixels are square when plotting to avoid distortions
ax.set_aspect("equal")

plt.show()

Well, this looks bad because that bright cloud is making it so the ground
pixels have only a small share of the full range of available values. This can
be mitigated by rescaling the intensity of the image to a smaller range of
reflectance values.
Well, this looks bad because some very bright pixels in the city are making the
majority of the other pixels have only a small share of the full range of
available values. This can be mitigated by rescaling the intensity of the image
to a smaller range of reflectance values.

.. note::

Using :meth:`xarray.DataArray.plot` instead of
:meth:`xarray.DataArray.plot.imshow` won't work and will display histograms
instead.


Rescaling intensity (AKA contrast stretching)
Expand All @@ -84,16 +106,21 @@ the following values by trial and error until it looks nice:

.. jupyter-execute::

rgb = xls.composite(scene, rescale_to=[0.03, 0.2])
rgb = xls.composite(scene, rescale_to=[0.01, 0.2])

# Pretty much the same plotting code
fig, ax = plt.subplots(1, 1, figsize=(10, 6))
rgb.plot.imshow(ax=ax)
ax.set_title(f"Rescaled composite: {rgb.attrs['long_name']}")
ax.set_title(f"Rescaled {rgb.attrs['long_name'].title()}")
ax.set_aspect("equal")
plt.show()

Notice that we can more clearly see the colors of the ground but we lose a lot
of detail in the clouds.
Notice that we can more clearly see the different colors of the forest and the
rivers.
However, it can still be a bit hard to distinguish between some of the water
bodies and the forest in the right side of the scene.
Other band combinations can generate composites that better highlight these
features.

.. note::

Expand All @@ -106,23 +133,48 @@ of detail in the clouds.
Color infrared composites
-------------------------

Another common type of composite is the color infrared (CIR) composites. These
change the bands used to NIR, red, and green and serve primarily to distinguish
healthy vegetation from other objects in the scene. Let's make one by specifying
this band combination to :func:`xlandsat.composite` to see if we can more clearly spot the dam flood.
Another common type of composite is the color infrared (CIR) composites,
which uses the NIR, red, and green bands as the red, green, blue channels.
The added information of the NIR band helps highlight vegetation, which can
help us distinguish between the water and forest on the right.
Let's make one by specifying
this band combination to :func:`xlandsat.composite` to see what happens:

.. jupyter-execute::

cir = xls.composite(scene, bands=("nir", "red", "green"), rescale_to=[0, 0.4])
cir = xls.composite(
scene, rescale_to=[0, 0.4], bands=("nir", "red", "green"),
)

fig, ax = plt.subplots(1, 1, figsize=(10, 6))
cir.plot.imshow(ax=ax)
ax.set_title(f"Composite: {rgb.attrs['long_name']}")
ax.set_title(cir.attrs["long_name"].title())
ax.set_aspect("equal")
plt.show()

The flood region can be clearly spotted in the image above as the brown/gray
blog in the center.
In this composite, the contrast between the forest and water bodies on the
right are much clearer!


Other composites
----------------

You can make pretty much any composite you'd like by passing the correct band
combination to :func:`xlandsat.composite`.
For example, let's make one with NIR as red, SWIR1 as green, and red as blue:

.. jupyter-execute::

custom = xls.composite(
scene, rescale_to=[0, 0.4], bands=("nir", "swir1", "red"),
)

fig, ax = plt.subplots(1, 1, figsize=(10, 6))
custom.plot.imshow(ax=ax)
ax.set_title(custom.attrs["long_name"].title())
ax.set_aspect("equal")
plt.show()

**With this, you can now make composites using any other band combination you
may want!**
This particular composite is great at distinguishing between built structures
in the city and along the canals (light green), the water ways (dark blue and
purple), and the forest (orange).
3 changes: 2 additions & 1 deletion xlandsat/datasets.py
Original file line number Diff line number Diff line change
Expand Up @@ -259,7 +259,8 @@ def fetch_manaus(untar=False):
Download a sample scene from Manaus, Brazil

Manaus is located in the Brazilian Amazon. The scene shows a part of the
city and the meeting of the Amazon and Negro rivers.
city and the meeting of the Solimões and Negro rivers to form the Amazon
river.

This is a cropped version of a Landsat 9 scene from 2023/07/23, during
the annual Amazon river floods.
Expand Down
Loading