Skip to content

Instantly share code, notes, and snippets.

@brfitzpatrick
Created September 13, 2018 13:47
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save brfitzpatrick/05bf7e518477a69a815b790864f49a1b to your computer and use it in GitHub Desktop.
Save brfitzpatrick/05bf7e518477a69a815b790864f49a1b to your computer and use it in GitHub Desktop.

Preparing Elevation Data for Visualisation with rayshader

library(magrittr)
library(viridis)
library(raster)
library(rayshader)

Using the DEM available here which the source states is using the LV03 (a.k.a. EPSG 21781) coordinate reference system:

dem.200m.rast <- raster('~/Downloads/data/DHM200.asc')

crs(dem.200m.rast) <- crs('+init=epsg:21781')

print(dem.200m.rast)
## class       : RasterLayer 
## dimensions  : 1201, 1926, 2313126  (nrow, ncol, ncell)
## resolution  : 200, 200  (x, y)
## extent      : 479900, 865100, 61900, 302100  (xmin, xmax, ymin, ymax)
## coord. ref. : +init=epsg:21781 +proj=somerc +lat_0=46.95240555555556 +lon_0=7.439583333333333 +k_0=1 +x_0=600000 +y_0=200000 +ellps=bessel +towgs84=674.4,15.1,405.3,0,0,0,0 +units=m +no_defs 
## data source : /home/ben/Downloads/data/DHM200.asc 
## names       : DHM200

Thus is appears that cells in this raster are squares.

Visualising the elevation raster with raster::plot

plot(dem.200m.rast, col = viridis(n = 1e3))

plot of chunk unnamed-chunk-4

My extraction of the data from the raster to a matrix results in a distorted image despite attempting to set the aspect ratio to 1:

dem.200m.mat <- raster::as.matrix(dem.200m.rast)

apply(X = dem.200m.mat, MARGIN = 2, FUN = rev) %>%
  t() %>%
    graphics::image(x = ., col = viridis(n = 1e3), asp = 1)

plot of chunk unnamed-chunk-5

dem.200m.mat.t <- t(dem.200m.mat)

system.time(
  rayshader::sphere_shade(dem.200m.mat.t, texture = 'bw') %>%
  rayshader::plot_map(asp = 1)
)

plot of chunk unnamed-chunk-6

##    user  system elapsed 
##   3.203   0.260   3.463
@mdsumner
Copy link

mdsumner commented Sep 14, 2018

Oh I see, the matrix is implicitly mapped to 0, 1 in image, I just never would do that bare matrix vis. rayshader also needs to keep the original axes but I don't yet understand why it's not 0, ncol; 0, nrow already

@brfitzpatrick
Copy link
Author

Thanks again for looking at this. In my first go at this I plotted the DEM with raster::plot then tried rayshader. The image plot was just a diagnostic to see where in the process the distortion was occurring. Could it be that raster::plot is actually using an aspect ratio different from one that is estimated from the CRS and that I need to supply this to image or rayshader?

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment