0

I have the two following files in R:

I am trying to now perform a merge on these files:

library(sf)

Bigger_Units  <- sf::st_read("C:/Users/ME/OneDrive/Documents/hr_shape/HR_000a18a_e.shp", options = "ENCODING=WINDOWS-1252") %>% 
  st_transform(crs = 4326) 


Smaller_Units <- sf::st_read("C:/Users/ME/OneDrive/Documents/shape7/lfsa000b16a_e.shp", options = "ENCODING=WINDOWS-1252") %>%  st_transform(crs = 4326) 


## Joined data
Dat <- st_join(Bigger_Units, Smaller_Units, join = st_contains, left = TRUE)

Unfortunately, this is returning the following error:

Error in s2_geography_from_wkb(x, oriented = oriented, check = check) : 
  Evaluation error: Found 3 features with invalid spherical geometry.
[21] Loop 0 is not valid: Edge 6509 has duplicate vertex with edge 6521
[27] Loop 0 is not valid: Edge 49 has duplicate vertex with edge 5931
[49] Loop 9 is not valid: Edge 1402 has duplicate vertex with edge 1574.

I tried researching this error online and came across the following links:

But I am unsure as to how I can apply the answers provided in these links to resolve my question.

Can someone please give me a hand with this?

Thank you!

stats_noob
  • 5,401
  • 4
  • 27
  • 83

1 Answers1

1

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
Karl Edwards
  • 305
  • 4
  • @ Karl Edwards: thank you so much for your answer! I will try your code when I get home! Just a question - what happens if a smaller unit is situated within two bigger units? Will this code output two rows and show that the smaller unit appears in two bigger units? Thank you so much! – stats_noob Nov 12 '22 at 18:56