geom_imshow() renders NaN values as transparent pixels¶import numpy as np
import rasterio as ro
from lets_plot import *
LetsPlot.setup_html()
Data are available at
Japan DEM (c) USGS. Processed for Japan by Skinner Regional Systems Analysis. Cambridge: Center for Geographic Analysis, Harvard University, 2015.
# Load data with Rasterio
jp_tif = ro.open('https://raw.githubusercontent.com/JetBrains/lets-plot-docs/master/data/japan_dem_wgs84/japan_dem_wgs84.tif')
# Metadata
jp_crs=jp_tif.crs
jp_bounds=jp_tif.bounds
print("{}\n{}".format(jp_crs, jp_bounds))
EPSG:4326 BoundingBox(left=127.9421598662918, bottom=29.201489591057303, right=155.88170436191123, top=50.02079185141284)
# create the "extent" to use in geom_imshow()
xmin=jp_bounds[0]
xmax=jp_bounds[2]
ymin=jp_bounds[1]
ymax=jp_bounds[3]
jp_ext=[xmin,xmax, ymin, ymax]
# Get the first band as a 2D numpy array
jp_arr=jp_tif.read(1)
type(jp_arr), np.shape(jp_arr)
(numpy.ndarray, (2038, 2735))
ggplot() + geom_imshow(jp_arr, extent=jp_ext)
Easy to see (chart below) that values in jp_arr form 2 clusters: large negative values and values close to 0.
After "normalization", those 2 clasters fall into opposite sides of the luminosity spectrum, and as result we can only see blak&wite picture.
(ggplot()
+ geom_histogram(aes(x=jp_arr.flatten()))
+ geom_vline(xintercept=-32768, color="red")
+ geom_label(label="-32768", x=-30000, y=5.2E6)
+ ggtitle("Distribution of elevations in Japan DEM")
+ xlab("elevation (m)")
)
NaN values¶The sea level elevation is traditionally encoded as minimal 16-bit value: -32768.
jp_arr_nan=jp_arr.copy().astype(np.float32)
jp_arr_nan[jp_arr<-30000]=np.nan
(ggplot()
+ geom_histogram(aes(x=jp_arr_nan.flatten()))
+ ggtitle("Distribution of elevations in Japan DEM", subtitle="excluding 'sea level'")
+ xlab("elevation (m)")
)
ggplot() + geom_imshow(jp_arr_nan, extent=jp_ext)