LTS, LTS impedance, slope impedance

CyclingImpedances.R

@@ -0,0 +1,74 @@
+# Using a network created by 'NetworkGenerator.R', add impedances for 
+# calculating cycling accessibility
+
+# Assumes input network is a one-way network that includes elevation, 
+# and that one-way daily traffic volumes are available
+
+addImpedances <- function() {
+
+  # Parameters -----------------------------------------------------------------
+  # Input network, to which cycling impedances are to be added, with layer names
+  # input.network <- "./output/test/network.sqlite"
+  input.network <- "./output/test/melbourne_network.sqlite"  #<<< OLD VERSION FOR TESTING
+  input.node.layer <- "nodes"
+  input.link.layer <- "links"
+
+  # Traffic file, with links layer name - file must match input.network, 
+  # and with a 'total_vol' column containing 1-way daily traffic
+  # traffic.file <- "./output/test/network_traffic.sqlite"
+  traffic.file <- "./output/test/links_with_traffic.sqlite" #<<< OLD VERSION FOR TESTING
+  traffic.link.layer <- "cars_aht" 
+
+  # Traffic multiplier - where volumes are for a sample of traffic only (eg
+  # multiplier of 20 if the volumes are a 5% sample; use 1 if full volumes)
+  traffic.multiplier <- 10
+
+  # Output location - same directory as input.network
+  output.location <- paste0(path_dir(input.network), "/networkWeighted.sqlite")
+
+
+  # Packages and functions -----------------------------------------------------
+  library(dplyr)
+  library(sf)
+  library(fs)
+
+  dir_walk(path="./functions/",source, recurse=T, type = "file")
+
+
+  # Load input network and traffic file ----------------------------------------
+  input.nodes <- st_read(input.network, layer = input.node.layer)  
+  input.links <- st_read(input.network, layer = input.link.layer)
+  traffic.links <- st_read(traffic.file, layer = traffic.link.layer)
+
+
+  # Add LTS and its impedance --------------------------------------------------
+  echo("Adding daily traffic volumes\n")
+  networkTraffic <- addTraffic(input.nodes,
+                               input.links, 
+                               traffic.links, 
+                               traffic.multiplier)
+  ## TO DO - maybe traffic can just be joined on link_id?  See whether traffic
+  ## file neatly uses the link_id's from the one-way input
+
+  echo("Adding LTS and its impedance\n")
+  networkLTS <- addLTS(networkTraffic[[1]], networkTraffic[[2]])
+
+
+  # Add slope impedance --------------------------------------------------------
+  echo("Adding slope impedance")
+  networkSlope <- addSlopeImped(networkLTS[[1]], networkLTS[[2]])
+
+
+  # Calculate total weight -----------------------------------------------------
+  echo("Calculating cycling weight")
+  networkWeighted <- 
+    list(networkSlope[[1]],
+         networkSlope[[2]] %>%
+           mutate(cycle.weight = length + LTS_imped + slope_imped))
+
+
+  # write output ---------------------------------------------------------------
+  st_write(networkWeighted[[1]], "./output/test/networkWeighted.sqlite", layer = "nodes")
+  st_write(networkWeighted[[2]], "./output/test/networkWeighted.sqlite", layer = "links")
+
+}

Destinations.R

@@ -0,0 +1,169 @@
+# Make file of destinations required for accessibilty routing
+
+library(dplyr)
+library(sf)
+library(osmextract)
+
+# 1 Download OSM extract ----
+# -----------------------------------------------------------------------------#
+# Download from https://www.interline.io/osm/extracts/ 
+
+## Downloaded for Melbourne - melbourne_australia.osm.pbf
+
+
+# 2 Converting to .gpkg format ----
+# -----------------------------------------------------------------------------#
+# input file name and project CRS
+INPUTFILE <- "./data/melbourne_australia.osm.pbf"
+
+PROJECT.CRS = 28355
+
+# check layers
+st_layers(INPUTFILE)
+
+# check keys
+options(max.print = 2000)
+polygon.tags <- oe_get_keys(INPUTFILE, layer = "multipolygons") %>% sort()
+point.tags <- oe_get_keys(INPUTFILE, layer = "points") %>% sort()
+line.tags <- oe_get_keys(INPUTFILE, layer = "lines") %>% sort()
+
+# create gpkg file in same directory as INPUTFILE, using the 'extra_tags' 
+# argument for specific extra tags required for various destination types
+oe_vectortranslate(INPUTFILE, layer = "multipolygons", 
+                   extra_tags = c("access", "building", "grades", "healthcare", 
+                                  "healthcare:speciality","isced:level", 
+                                  "network", "operator",
+                                  "operator:type", "public_transport", "railway",
+                                  "school", "social_facility", "sport", 
+                                  "tourism", "train"))
+oe_vectortranslate(INPUTFILE, layer = "points",
+                   extra_tags = c("access", "amenity", "building", "grades", 
+                                  "healthcare", "healthcare:speciality", 
+                                  "isced:level", "landuse", "leisure", 
+                                  "network", "operator",
+                                  "operator:type", "public_transport", "railway",
+                                  "school", "shop", "social_facility",  "sport",
+                                  "tourism", "train"))
+oe_vectortranslate(INPUTFILE, layer = "lines",
+                   extra_tags = c("access", "amenity", "building", "grades", 
+                                  "healthcare", "healthcare:speciality", 
+                                  "isced:level", "landuse", "leisure", 
+                                  "network", "operator",
+                                  "operator:type", "public_transport", "railway",
+                                  "school", "shop", "social_facility",  "sport",
+                                  "tourism", "train",
+                                  "smoothness", "surface"))
+oe_vectortranslate(INPUTFILE, layer = "multilinestrings")
+oe_vectortranslate(INPUTFILE, layer = "other_relations")
+
+
+# 3 Read in the .gpkg file ----
+# -----------------------------------------------------------------------------#
+GPKG <- "./data/melbourne_australia.gpkg"
+
+polygons <- st_read(GPKG, layer = "multipolygons") %>% st_transform(PROJECT.CRS)
+points <- st_read(GPKG, layer = "points") %>% st_transform(PROJECT.CRS)
+lines <- st_read(GPKG, layer = "lines") %>% st_transform(PROJECT.CRS)
+multilines <- st_read(GPKG, layer = "multilinestrings") %>% st_transform(PROJECT.CRS)
+other_relations <- st_read(GPKG, layer = "other_relations") %>% st_transform(PROJECT.CRS)
+
+
+# 4 Extract required destinations ----
+# -----------------------------------------------------------------------------#
+
+## 4.1 Tag combinations for feature types  and network ----
+## ----------------------------------------------------------------------------#
+# load functions for locating specific feature types
+source("./functions/getDestinationTypes.R")
+
+# load network
+NETWORK <- "./output/test/network.sqlite"  #<<< CHECK FINAL NAME
+NODE.LAYER <- "nodes"
+LINK.LAYER <- "links"
+
+network.nodes <- st_read(NETWORK, layer = NODE.LAYER)  
+network.links <- st_read(NETWORK, layer = LINK.LAYER)
+
+
+## 4.2 Compile point and polygon destinations ----
+## ----------------------------------------------------------------------------#
+destination.layer <- function(layer) {
+  return(
+    bind_rows(
+      getPlayground(layer) %>% mutate(dest_type = "playground"),
+      getPark(layer) %>% mutate(dest_type = "park"),
+      getSport(layer) %>% mutate(dest_type = "sport"),
+      getKindergarten(layer) %>% mutate(dest_type = "kindergarten"),
+      getCommunity(layer) %>% mutate(dest_type = "community_centre"),
+      getLibrary(layer) %>% mutate(dest_type = "library"),
+      getPrimary(layer) %>% mutate(dest_type = "primary_school"),
+      getSecondary(layer) %>% mutate(dest_type = "secondary_school"),
+      getClinic(layer) %>% mutate(dest_type = "health_clinic"),
+      getDentist(layer) %>% mutate(dest_type = "dentist"),
+      getPharmacy(layer) %>% mutate(dest_type = "pharmacy"),
+      getConvenience(layer) %>% mutate(dest_type = "convenience_store"),
+      getSupermarket(layer) %>% mutate(dest_type = "supermarket"),
+      getShop(layer) %>% mutate(dest_type = "shop"),
+      getPost(layer) %>% mutate(dest_type = "post_office"),
+      getBank(layer) %>% mutate(dest_type = "bank"),
+      getRestaurant(layer) %>% mutate(dest_type = "restaurant"),
+      getCafe(layer) %>% mutate(dest_type = "cafe")
+    ))
+}
+
+# create tables of destinations, and allocate unique id's (so features with
+# multiple nodes can be grouped by the id where required)
+destination.pt <- 
+  bind_rows(destination.layer(points),
+            getStation() %>% mutate(dest_type = "railway_station")) %>%
+  mutate(dest_id = row_number())
+
+destination.poly <- 
+  destination.layer(polygons) %>%
+  mutate(dest_id = max(destination.pt$dest_id) + row_number())
+
+
+## 4.3 Find relevant nodes ----
+## ----------------------------------------------------------------------------#
+# TO CONFIRM:-
+# For all destinations except parks and schools ('small features'), relevant 
+# node is nearest node to point or to polygon centroid
+# For parks and schools ('large features'):
+# - points are buffered to 50m to create a polygon feature,
+# - for buffered points and polygons, relevant nodes are all nodes within the 
+#   feature and terminal nodes of links within 30m of boundary, or if none, 
+#   then nearest node to boundary
+
+# Maybe this should be all nodes within 30m of buffered feature, and if link is within
+# 30m of boundary but doesn't have a node within the buffer, then also its closest terminal
+# node ???
+
+dest.small <- bind_rows(destination.pt,
+                        destination.poly %>% st_centroid()) %>%
+  filter(!(dest_type %in% c("park", "primary_school", "secondary_school")))
+near_node <- network.nodes$id[st_nearest_feature(dest.small, network.nodes)]
+dest.small.with.nodes <- cbind(dest.small %>% st_drop_geometry(), near_node) 
+
+
+## NOTE - the code below is a simplified version which just finds nodes within
+## feature or its 30m buffer, or nearest node if none - doesn't extend to terminal
+## nodes of nearby features
+dest.large <- bind_rows(destination.pt %>% st_buffer(50),
+                        destination.poly) %>%
+  filter(dest_type %in% c("park", "primary_school", "secondary_school"))
+
+dest.large.found.nodes <- dest.large %>%
+  st_buffer(30) %>%
+  st_intersection(., network.nodes %>% dplyr::select(near_node = id))
+
+dest.large.need.nodes <- dest.large %>%
+  filter(!(dest_id %in% dest.large.found.nodes$dest_id))
+near_node <- network.nodes$id[st_nearest_feature(dest.large.need.nodes, network.nodes)]
+
+dest.large.with.nodes <- bind_rows(dest.large.found.nodes %>% st_drop_geometry(),
+                                   cbind(dest.large.need.nodes %>% st_drop_geometry(),
+                                         near_node))
+
+dest.with.nodes <- bind_rows(dest.small.with.nodes,
+                             dest.large.with.nodes)
+

NetworkGenerator.R

@@ -270,7 +270,13 @@ makeNetwork<-function(outputFileName="test"){
                                                outputCrs=outputCrs)) 
   }
 
-  networkFinal <- networkRestructured
+  # Make network oneway (required because cycling impedances such as level of 
+  # traffic stress and slope may be different in each direction)
+  echo("Making all links one way\n")
+  networkOneway <- makeEdgesOneway(networkRestructured[[1]], 
+                                   networkRestructured[[2]])
+
+  networkFinal <- networkOneway
 
   # writing outputs
   echo("========================================================\n")