I have multipolygons that span across the globe. I was wondering if it was possible to re-project each row, i.e. each country/ ISO3 to a local projection, using a central grouping mechanism.
For example, my polygons within the a particular region say SE Asia is then re-projected to a local projection for that area. Preferably to a LAEA type projection as I want to create buffers around these polygons for further analyses.
Here is the structure of my data and the projection I have used for both multipolygons and multipoints.
Simple feature collection with 6 features and 2 fields
geometry type: MULTIPOLYGON
dimension: XY
bbox: xmin: -11959630 ymin: -4449011 xmax: 5370704 ymax: 10196680
epsg (SRID): NA
proj4string: +proj=laea +lat_0=45.5 +lon_0=-114.125 +x_0=0 +y_0=0 +a=6371007.181 +b=6371007.181 +units=m +no_defs
PARENT_ISO ISO3 geom
1 ABNJ ABNJ MULTIPOLYGON (((-7616720 -4...
2 ARE ARE MULTIPOLYGON (((2426453 101...
3 ATG ATG MULTIPOLYGON (((5364022 -14...
4 AUS AUS MULTIPOLYGON (((-10014124 -...
5 AUS CCK MULTIPOLYGON (((-8801631 79...
6 AUS CXR MULTIPOLYGON (((-9722602 63...
> head(c_pt)
Simple feature collection with 4 features and 2 fields
geometry type: MULTIPOINT
dimension: XY
bbox: xmin: -409654.4 ymin: -3100315 xmax: 3235020 ymax: -1443240
epsg (SRID): NA
proj4string: +proj=laea +lat_0=45.5 +lon_0=-114.125 +x_0=0 +y_0=0 +a=6371007.181 +b=6371007.181 +units=m +no_defs
PARENT_ISO ISO3 geom
1 BLZ BLZ MULTIPOINT (2757788 -286421...
2 GTM GTM MULTIPOINT (2773602 -287236...
3 HND HND MULTIPOINT (2883494 -283266...
4 MEX MEX MULTIPOINT (-409654.4 -1831...
Creating buffers from this projection seemed to create a lot of distortion, as can be seen by Madagascar.
I do realise I could split the data by country and reproject each new country file, but was wondering if there was a way that was possible whilst keeping all the data together?
Or is this just not the way to handle this kind of situation?
Example of my workflow
#Reproject coral reefs
c_pt <- st_transform(c_pt, crs = "+proj=laea +lat_0=45.5 +lon_0=-114.125 +x_0=0 +y_0=0 +a=6371007.181 +b=6371007.181 +units=m +no_defs")
c_py <- st_transform(c_py, crs = "+proj=laea +lat_0=45.5 +lon_0=-114.125 +x_0=0 +y_0=0 +a=6371007.181 +b=6371007.181 +units=m +no_defs")
#Create buffer of 100km around coral reefs
buffer_100_py <- st_buffer(c_py, dist = 100000)
buffer_100_pt <- st_buffer(c_pt, dist = 100000)
#this is gives overlapping portion of polygons
overlap <- st_intersection(buffer_100_py, buffer_100_pt)
#this is to delete the overlapping polygons
diffPoly <- st_difference(buffer_100_pt, st_union(overlap))
#joining the cleaned point buffer data with the polygon data
buffer_100_all<-rbind(diffPoly, buffer_100_py)
#Now to clean up the polygons to remove overlap by union by feature - there are still overlaps of polygons,
#But this should be fine for the population extraction
buffer_100_all_one<-ms_dissolve(buffer_100_all, field = "ISO3")
Extra:
Dealing with buffers crossing the dateline
I have attempted to correct this using:
#deal with buffers that cross the dateline
buffer_100_all_one <- buffer_100_all %>% st_wrap_dateline(options = c("WRAPDATELINE=YES", "DATELINEOFFSET=180"), quiet = FALSE)
#Note: I cannot seem to workout how to correct for buffers crossing the dateline, this does not correct the plotting of this.
#Lon_wrap in crs projections???
But think maybe this could be resolved with the projection?
