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 [ ]:
os.chdir('../')
os.getcwd()

Open the Flatirons DEM raster:

src = rio.open('workshopdata/dem.tif')
In [29]:
src = rio.open('workshopdata/dem.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 [30]:
plt.imshow(src.read(1))
Out[30]:
<matplotlib.image.AxesImage at 0x28d012e1b08>

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

show(src)
In [31]:
show(src)
Out[31]:
<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 [32]:
show(src.read(1))
Out[32]:
<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 [33]:
show(src.read(1),transform=src.transform)
Out[33]:
<AxesSubplot:>

Colormaps!¶

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

show(src.read(1), transform = src.transform, cmap = 'gray')
In [34]:
show(src.read(1), transform=src.transform, cmap='gray')
Out[34]:
<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 [35]:
show(src.read(1), transform = src.transform, cmap = 'Reds')
Out[35]:
<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 [36]:
fig, ax = plt.subplots(1,figsize=(12,12))
show(src.read(1), transform = src.transform, cmap='gray', ax=ax)
Out[36]:
<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 [37]:
fig, ax = plt.subplots(figsize=(12,12))
show(src.read(1), transform = src.transform, contour=True, ax=ax)
Out[37]:
<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 [38]:
fig, ax = plt.subplots(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[38]:
<AxesSubplot:>

Save some memory...

In [39]:
del src

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

We'll use our cropped NAIP image of campus:

naip = rio.open('NAIP_swQuad.tif')
In [42]:
naip = rio.open('NAIP_swQuad.tif')
show(naip)
In [53]:
show(naip)
Out[53]:
<AxesSubplot:>

You can also specify particular bands:

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

Near infrared:

show((naip,4))
In [45]:
show((naip,4))
Out[45]:
<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 [46]:
image = naip.read()

Take a quick look...

show(image)
In [47]:
show(image)
Out[47]:
<AxesSubplot:>

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 [48]:
rgb = image[0:3]
show(rgb)
In [49]:
show(rgb)
Out[49]:
<AxesSubplot:>

Now, finally we can plot our rgb image properly...

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 [50]:
fig, ax = plt.subplots(figsize=(12,12))
show(rgb, transform = naip.transform, ax=ax, title='Aerial Image of CU Campus')
Out[50]:
<AxesSubplot:title={'center':'Aerial Image of CU Campus'}>

Sweet!¶