I want to clip a polygon shapefile using
sf::st_intersection, however the territory is not cut straight.
Here is a RE. I use a Canadian province shapefile that can be downloaded easily:
library(sf) library(dplyr) the_zip <- tempfile(fileext=".zip") download.file("http://www12.statcan.gc.ca/census-recensement/2011/geo/bound-limit/files-fichiers/2016/lpr_000b16a_e.zip",the_zip) unzip(the_zip, exdir = dirname(the_zip)) can <- st_read(paste0(dirname(the_zip), "/lpr_000b16a_e.shp")) %>% filter(PRNAME %in% c("Quebec / Québec", "Ontario")) %>% st_simplify %>% st_transform(st_crs(4326)) clip <- st_as_sfc(st_bbox(c(xmin = -100, xmax = -60, ymax = 50, ymin = 40), crs = st_crs(4326)))
To help visualize the problem here is the shape I want to cut and the clipping polygon in blue:
dev.off() plot(can[,1], reset=F, axes=T) plot(clip, add=T, border="blue", lwd=2)
Now, I can clip the polygon:
can_clipped <- st_intersection(can, clip)
However, when I plot it, I get:
dev.off() plot(can_clipped[,1], reset=F, axes=T) plot(clip, add=T, border="blue", lwd=2)
As you can see, the clip do not fit well with the desired output. Also, both polygons in the shapefile seems bent.
What causes this problem and can it be fixed?
Could it be linked with the fact I'm clipping on shapes in
I do want to cut my shape on a strait line at the 50th parallel, so I don't really want to play around another projection.
Or maybe it's related to the warning I get when doing the
Warning message: attribute variables are assumed to be spatially constant throughout all geometries