[wrridgeway]

This PR creates a new raw script that ingests OSM data for "main" roads so that we don't have to do it again any time we rebuild corner lots. As @jeancochrane stated in #1100 (comment) we might want to consider renaming the raw script to signal the data's broader usefulness, but for this PR I'm leaving the naming scheme for the sake of a more easily understandable lineage.

I rebuilt corners for 2023 in s3://ccao-data-warehouse-dev-us-east-1/spatial/ccao/corner/ using this new warehouse script to compare its output to what's currently in spatial.corner. Below I have code that shows a comparison for Barrington. There are differences, but they are quite minor, and the disagreements do seem to be largely for parcels where defining a corner status leaves room for argument.

library(arrow)
library(DBI)
library(dplyr)
library(noctua)
library(purrr)
library(sf)
library(mapview)

noctua_options(unload = TRUE)
AWS_ATHENA_CONN_NOCTUA <- dbConnect(noctua::athena(), rstudio_conn_tab = FALSE)

# Grab old corner lots
old <- dbGetQuery(conn = AWS_ATHENA_CONN_NOCTUA, "
select pin10, geometry from spatial.corner where year = '2023' and township_code = '10'
") %>%
  st_as_sf(crs = 4326)

# Grab new corner lots
new <- dbGetQuery(conn = AWS_ATHENA_CONN_NOCTUA, "
select pin10, geometry from z_dev_wridgeway_spatial.corner where year = '2023' and township_code = '10'
") %>%
  st_as_sf(crs = 4326)

# Grab polygons for corners that appear in either only new or old data
diff <- full_join(
  new %>% st_drop_geometry() %>% mutate(new = TRUE),
  old %>% st_drop_geometry() %>% mutate(old = TRUE),
  by = "pin10"
) %>%
  filter(is.na(old) | is.na(new))

# Grab parcels for mapping
parcels <- dbGetQuery(conn = AWS_ATHENA_CONN_NOCTUA, "
select pin10, geometry_3435 as geometry from spatial.parcel where year = '2023' and town_code = '10'
") %>%
  st_as_sf(crs = 3435)

new_corners_only <- parcels %>% inner_join(diff %>% filter(new))
old_corners_only <- parcels %>% inner_join(diff %>% filter(old))

# Create map for comparison
mapview(old, color = "red") +
  mapview(new, color = "blue") +
  mapview(parcels, col.regions = "grey") +
  mapview(new_corners_only, col.regions = "purple") +
  mapview(old_corners_only, col.regions = "yellow")

More broadly, we also see that the total number of corner lots for 2023 is, again, very similar to what it is currently in our data:

select count() as corners,
	'new' as source
from z_dev_wridgeway_spatial.corner
where year = '2023'
union
select count() as corners,
	'old' as source
from spatial.corner
where year = '2023'

Corners	Source
177265	old
177251	new

[wrridgeway]

This is the only substantive change to this script. It was throwing an error when we joined cross_final to town_parcels because the resulting dataframe was no longer an sf object. This produces the same exact output, but retains the sf class for the resulting dataframe.

[wrridgeway]

I suggest reviewing this diff with whitespace hidden.

[jeancochrane]

Looks good to me, thanks for the cleanup and the thorough QC!