shadow
is an R package for geometric shadow calculations in an urban environment. A detailed description can be found in the R Journal paper (2019):
CRAN version:
install.packages("shadow")
GitHub version:
install.packages("remotes") remotes::install_github("michaeldorman/shadow")
The complete documentation can be found at https://michaeldorman.github.io/shadow/.
library(shadow) #> Loading required package: sp library(raster) # Point location = rgeos::gCentroid(build) # Time time = as.POSIXct( "2004-12-24 13:30:00", tz = "Asia/Jerusalem" ) # Location in geographical coordinates proj4string(location) = CRS("+init=epsg:32636") location_geo = sp::spTransform( location, CRS("+init=epsg:4326") ) crs(location) = NA # Solar position solar_pos = maptools::solarpos( crds = location_geo, dateTime = time ) solar_pos #> [,1] [,2] #> [1,] 208.7333 28.79944 # Shadow height at a single point h = shadowHeight( location = location, obstacles = build, obstacles_height_field = "BLDG_HT", solar_pos = solar_pos ) # Result h #> [,1] #> [1,] 19.86451 # Visualization sun = shadow:::.sunLocation( location = location, sun_az = solar_pos[1, 1], sun_elev = solar_pos[1, 2] ) sun_ray = ray(from = location, to = sun) build_outline = as(build, "SpatialLinesDataFrame") inter = rgeos::gIntersection(build_outline, sun_ray) plot(build) plot(sun_ray, add = TRUE, col = "yellow") plot(location, add = TRUE) text(location, paste(round(h, 2), "m"), pos = 3) plot(inter, add = TRUE, col = "red")
# Raster template ext = as(extent(build)+50, "SpatialPolygons") r = raster(ext, res = 2) # Shadow height surface height_surface = shadowHeight( location = r, obstacles = build, obstacles_height_field = "BLDG_HT", solar_pos = solar_pos, parallel = 2 ) # Visualization plot(height_surface, col = grey(seq(0.9, 0.2, -0.01))) contour(height_surface, add = TRUE) plot(build, add = TRUE, border = "red") text(rgeos::gCentroid(build, byid = TRUE), build$BLDG_HT) text(location, paste(round(h, 2), "m"), pos = 3, col = "red", font = 2) plot(sun_ray, add = TRUE, col = "yellow") plot(inter, add = TRUE, col = "red") plot(location, add = TRUE)