Rasterio Visualizations¶

Now we'll do a brief overview of visualizations--although I will admit that integrating these techniques with Matplotlib can be a very deep dive... this is just a primer.

Imports:

In [1]:
import rasterio as rio
from rasterio.plot import show
from matplotlib import pyplot as plt
from rasterio.plot import adjust_band
import numpy as np
import os

Important: only run this once or you'll wind up in the wrong directory.

os.chdir('../')
os.getcwd()
In [2]:
os.chdir('../')
os.getcwd()
Out[2]:
'/home/jupyter-chain.saw'

Open the Flatirons DEM raster:

src = rio.open('workshopdata/Flatirons_DEM_1m.tif')
In [3]:
src = rio.open('workshopdata/Flatirons_DEM_1m.tif')

The matplotlib way... What we're doing here is reading it as an array and having matplotlib plot the array as an image:

plt.imshow(src.read(1))
In [4]:
plt.imshow(src.read(1))
Out[4]:
<matplotlib.image.AxesImage at 0x7f0c70d48460>

OR, as we've used already, the built in Rasterio method:

show(src)
In [5]:
show(src)
Out[5]:
<AxesSubplot: >

Note that this one is in the proper units...

If you use show() to plot the single band, you lose the units...

show(src.read(1))
In [6]:
show(src.read(1))
Out[6]:
<AxesSubplot: >

But, you can provide the transform and it will provide the proper units... this may seem pointless, but it is actually a useful method, as you'll soon see...

show(src.read(1), transform = src.transform)
In [7]:
show(src.read(1), transform = src.transform)
Out[7]:
<AxesSubplot: >

Colormaps!¶

Let's make it good ol' "DEM gray"...

show(src.read(1), transform = src.transform, cmap = 'gray')
In [8]:
show(src.read(1), transform = src.transform, cmap = 'gray')
Out[8]:
<AxesSubplot: >

Here are Matplotlib's colormaps: https://matplotlib.org/stable/tutorials/colors/colormaps.html

Let's have some fun with this:

show(src.read(1), transform = src.transform, cmap = 'Reds')
In [9]:
show(src.read(1), transform = src.transform, cmap = 'Reds')
Out[9]:
<AxesSubplot: >

Size¶

Tired of looking at these tiny plots? We can use the Matplotlib conventions...

fig, ax = plt.subplots(1, figsize=(12,12))
show(src.read(1), transform = src.transform, cmap = 'gray', ax=ax)
In [10]:
fig, ax = plt.subplots(1, figsize = (12,12,))
show(src.read(1), transform = src.transform, cmap = 'gray', ax=ax)
Out[10]:
<AxesSubplot: >

It's quite easy to create contours from a DEM:

fig, ax = plt.subplots(1, figsize=(12,12))
show(src.read(1), transform = src.transform,contour=True, ax=ax)
In [12]:
fig, ax = plt.subplots(1, figsize=(12,12))
show(src.read(1), transform = src.transform, contour=True, ax=ax)
Out[12]:
<AxesSubplot: >

Or do both Contours and a DEM together:

fig, ax = plt.subplots(1, figsize=(12,12))
show(src.read(1), transform = src.transform, cmap = 'gray', ax=ax)
show(src.read(1), transform = src.transform,contour=True, ax=ax)
In [13]:
fig, ax = plt.subplots(1, figsize=(12,12))
show(src.read(1), transform = src.transform, cmap='gray', ax=ax)
show(src.read(1), transform = src.transform, contour=True, ax=ax)
Out[13]:
<AxesSubplot: >

Save some memory...

In [14]:
del src

Okay, let's work on some multiband imagery...¶

We'll use our cropped NAIP image of campus:

naip = rio.open('NAIP_Campus_Clip.tif')
In [16]:
naip = rio.open('NAIP_Campus_Clip.tif')

show(naip)
In [17]:
show(naip)
Out[17]:
<AxesSubplot: >

You can also specify particular bands:

show((naip,2))
In [18]:
show((naip,2))
Out[18]:
<AxesSubplot: >

Near infrared:

show((naip,4))
In [19]:
show((naip,4))
Out[19]:
<AxesSubplot: >

Let's work on displaying the true color version...

To display a true color RGB image, you need to normalize the bands....

Start by reading the image into a Numpy array. This will create an n-dimensional array--essentially the four different bands stacked.

image = naip.read()
In [20]:
image = naip.read()

But, to create the true color visualization, we really just need the rgb bands--that is, the first 3 bands. So we'll use Numpy to grab just those three:

rgb = image[0:3]
In [23]:
rgb = image[0:3]

np.shape(rgb)
In [24]:
np.shape(rgb)
Out[24]:
(3, 3841, 3174)

We'll need to normalize the rgb array so that all the values range betwee 0 and 1. We'll use rastio.plot.adjust_band for this:

rgb_norm = adjust_band(rgb)
In [25]:
rgb_norm = adjust_band(rgb)

So what did we just do?

We normalized the rgb image to have values between 0 and 1, rather than rgb values...

In [26]:
print(rgb_norm.min())
print(rgb_norm.max())
print(image.min())
print(image.max())

0.0
1.0
0
255

Now, finally we can plot our normalized rgb image

Might crush memory!!!

fig, ax = plt.subplots(1, figsize=(12,12))
show(rgb_norm, transform = naip.transform, ax=ax, title = 'Aerial Image of CU Campus')
In [27]:
fig, ax = plt.subplots(1, figsize=(12,12))
show(rgb_norm, transform= naip.transform, ax=ax, title = 'Aerial Image of CU Campus')
Out[27]:
<AxesSubplot: title={'center': 'Aerial Image of CU Campus'}>

Sweet!¶