It looks like the problem is related to the choice of coordinate reference system.
EPSG:4326 is an ellipsoidal 2D system with units of degrees used for the whole world. See https://epsg.io/4326
EPSG:3348, on the other hand, is a 2D Cartesian system with units of metres used specifically for Canada. See https://epsg.io/3348
require( sf )
DATA_PATH <- '/path/to/your/data'
cfg <- list(
regions = file.path( DATA_PATH, 'HR_000a18a_e/HR_000a18a_e.shp' )
, postcodes = file.path( DATA_PATH, 'lfsa000b16a_e/lfsa000b16a_e.shp' )
, crs = 3348
)
regions <- sf::st_read(
cfg$regions, options = "ENCODING=WINDOWS-1252"
) %>% st_transform( crs = cfg$crs )
postcodes <- sf::st_read(
cfg$postcodes, options = "ENCODING=WINDOWS-1252"
) %>% st_transform( crs = cfg$crs )
joined_data <- st_join( regions, postcodes, join = st_contains, left = TRUE )
joined_data
# Simple feature collection with 1386 features and 8 fields
#Geometry type: MULTIPOLYGON
#Dimension: XY
#Bounding box: xmin: 3658000 ymin: 658900 xmax: 9019000 ymax: 6083000
#Projected CRS: NAD83(CSRS) / Statistics Canada Lambert
#First 10 features:
# HR_UID ENGNAME #FRENAME
#1 3537 City of Hamilton Health Unit Circonscription sanitaire de #la cité de Hamilton
#1.1 3537 City of Hamilton Health Unit Circonscription sanitaire de #la cité de Hamilton
#1.2 3537 City of Hamilton Health Unit Circonscription sanitaire de #la cité de Hamilton
#1.3 3537 City of Hamilton Health Unit Circonscription sanitaire de #la cité de Hamilton
#1.4 3537 City of Hamilton Health Unit Circonscription sanitaire de #la cité de Hamilton
#1.5 3537 City of Hamilton Health Unit Circonscription sanitaire de #la cité de Hamilton
#1.6 3537 City of Hamilton Health Unit Circonscription sanitaire de #la cité de Hamilton
#1.7 3537 City of Hamilton Health Unit Circonscription sanitaire de #la cité de Hamilton
#1.8 3537 City of Hamilton Health Unit Circonscription sanitaire de #la cité de Hamilton
#1.9 3537 City of Hamilton Health Unit Circonscription sanitaire de #la cité de Hamilton
# SHAPE_AREA SHAPE_LEN CFSAUID PRUID PRNAME #geometry
#1 1.213e+09 173964 L8G 35 Ontario MULTIPOLYGON (((7192649 886...
#1.1 1.213e+09 173964 L8L 35 Ontario MULTIPOLYGON (((7192649 886...
#1.2 1.213e+09 173964 L8M 35 Ontario MULTIPOLYGON (((7192649 886...
#1.3 1.213e+09 173964 L8N 35 Ontario MULTIPOLYGON (((7192649 886...
#1.4 1.213e+09 173964 L8S 35 Ontario MULTIPOLYGON (((7192649 886...
#1.5 1.213e+09 173964 L8T 35 Ontario MULTIPOLYGON (((7192649 886...
#1.6 1.213e+09 173964 L8V 35 Ontario MULTIPOLYGON (((7192649 886...
#1.7 1.213e+09 173964 L9G 35 Ontario MULTIPOLYGON (((7192649 886...
#1.8 1.213e+09 173964 L9H 35 Ontario MULTIPOLYGON (((7192649 886...
#1.9 1.213e+09 173964 L9K 35 Ontario MULTIPOLYGON (((7192649 886...
# ---------------------------------------------------------------------
# Using EPSG:4326, I see the same kind of error you reported:
joined_data <- st_join( regions, postcodes, join = st_contains, left = TRUE )
## Error in wk_handle.wk_wkb(wkb, s2_geography_writer(oriented = oriented, :
## Loop 0 is not valid: Edge 6509 has duplicate vertex with edge 6521