Skip to content

Instantly share code, notes, and snippets.

@tonyelhabr
Last active October 10, 2019 12:49
Show Gist options
  • Save tonyelhabr/525a62b7eb6cf8b5d7e6eee315b55a5a to your computer and use it in GitHub Desktop.
Save tonyelhabr/525a62b7eb6cf8b5d7e6eee315b55a5a to your computer and use it in GitHub Desktop.
R functions for creating a voronoi map for ERCOT.
download_spdf_bound_ercot <- function(..., path = config$path_spdf_bound, backup = FALSE) {
url <- 'http://www.ercot.com/content/cdr/static/ercot_boundary.kml'
download.file(url = url, destfile = config$path_kml_bound, quiet = TRUE)
bound <-
config$path_kml_bound %>%
xml2::read_xml() %>%
xml2::xml_children() %>%
xml2::xml_children() %>%
magrittr::extract(4) %>%
xml2::xml_children() %>%
magrittr::extract(3) %>%
xml2::xml_children() %>%
magrittr::extract(4) %>%
xml2::xml_children() %>%
magrittr::extract(4) %>%
xml2::xml_text() %>%
# readLines()
str_split(',0\\s+', simplify = TRUE) %>%
str_replace_all('[^(0-9|.|,|\\-)]','') %>%
str_split(',', simplify = TRUE) %>%
as_tibble() %>%
purrr::set_names(c('lon', 'lat')) %>%
mutate_all(list(as.numeric)) %>%
filter(!is.na(lon) & !is.na(lat)) %>%
bind_rows(
.,
slice(., 1L)
)
bound_nontop <-
bound %>%
arrange(desc(row_number())) %>%
slice(c(24:n()))
bound_state <-
teplot::get_map_data_state(state = 'tx') %>%
as_tibble()
bound_top <-
bound_state %>%
slice(c(840:870)) %>%
select(lon = long, lat)
bound_fix <-
bind_rows(
bound_top,
bound_nontop
) %>%
bind_rows(
.,
slice(., 1L)
)
# # bound_fix %>% count(lon, lat)
# bound_fix %>%
# ggplot(aes(lon, lat)) +
# geom_path()
bound_poly <-
bound_fix %>%
sp::Polygon() %>%
list() %>%
sp::Polygons(ID = 'dummy') %>%
list() %>%
sp::SpatialPolygons()
spdf_bound <-
sp::SpatialPolygonsDataFrame(
bound_poly,
data.frame(state = c('texas'), row.names = c('dummy'))
)
# sp::proj4string(spdf_bound) <- sp::CRS('+init=epsg:4326')
# # NOTE: This is actually the same coordinate system. By transforming here,
# # removing some of the attributes to match with `spdf_pnt`.
# spdf_bound <- sp::spTransform(spdf_bound, sp::CRS('+proj=longlat +datum=WGS84'))
# spdf_bound
write_rds(spdf_bound, path = path)
}
get_spdf_bound_ercot <- function(..., path = config$path_spdf_bound, overwrite = FALSE, backup = FALSE) {
exists <- fs::file_exists(path)
if(!overwrite & !exists) {
.display_warning(
'Needed file does not already exist. ',
'Ignoring `overwrite = FALSE` and trying to download now.',
...
)
overwrite <- TRUE
}
if(!exists | overwrite) {
f_download_possibly <-
purrr::possibly(
.f = ~download_sdpf_bound_ercot(..., path = path),
otherwise = NULL
)
res <- f_download_possibly()
} else {
res <- .import_data_from_path(..., path = path)
}
if (!is.null(res)) {
return(res)
}
if(!exists) {
.display_warning(
'Could not download data.',
...
)
} else if (overwrite) {
.display_warning(
'Could not download data. Nonetheless, the file may already exist. ',
'(Try `overwrite = FALSE`)',
...
)
}
invisible(NULL)
}
create_spdf <- function(data) {
stopifnot(is.data.frame(data))
stopifnot(sum(c('lat', 'lon') %in% colnames(data)) == 2)
stopifnot(any('value' %in% colnames(data)))
data <-
data %>%
filter(!is.na(lat), !is.na(lon))
spdf_pnt <-
sp::SpatialPointsDataFrame(
cbind(data$lon, data$lat),
data,
match.ID = TRUE
)
# spdf_pnt@bbox <- spdf_base_map@bbox
# sp::proj4string(spdf_pnt) <- sp::proj4string(spdf_base_map)
spdf_pnt@bbox <- spdf_bound@bbox
sp::proj4string(spdf_pnt) <- sp::proj4string(spdf_bound)
spdf_pnt
}
convert_spdf_pnt_to_sf_poly <- function(spdf_pnt) {
# stopifnot(class(spdf_pnt) == '')
th <-
spdf_pnt %>%
maptools::as.ppp.SpatialPointsDataFrame() %>%
spatstat::dirichlet() %>%
methods::as('SpatialPolygons')
sp::proj4string(th) <- sp::proj4string(spdf_pnt)
th_z <- sp::over(th, spdf_pnt)
spdf_pnt_poly <- sp::SpatialPolygonsDataFrame(th, th_z)
# spdf_pnt_poly_trim <- raster::intersect(spdf_base_map, spdf_pnt_poly)
spdf_pnt_poly_trim <- raster::intersect(spdf_bound, spdf_pnt_poly)
spdf_pnt_poly_trim %>% sf::st_as_sf()
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment