A square, both in the correct order and clockwise:
correct_json = '{
"type": "Polygon",
"coordinates": [
[
[30.12370823171228, 51.44399790755677],
[30.12370823171228, 50.42971299122449],
[31.72930089180182, 50.42971299122449],
[31.72930089180182, 51.44399790755677],
[30.12370823171228, 51.44399790755677]
]
]
}'
cw_json = '{
"type": "Polygon",
"coordinates": [
[
[30.12370823171228, 51.44399790755677],
[31.72930089180182, 51.44399790755677],
[31.72930089180182, 50.42971299122449],
[30.12370823171228, 50.42971299122449],
[30.12370823171228, 51.44399790755677]
]
]
}'
GeoJSON spanning the antimeridian
anti_json = '{
"type": "Polygon",
"coordinates": [
[
[-177.1, 48.6],
[179.6, 49.0],
[179.3, 47.5],
[-178.6, 47.1],
[-177.1, 48.6]
]
]
}'
We’ll use the lower-level geojsonR
package, because it
allows us to load in coordinates as a matrix, and see whether they are
reordered later.
If we only wanted to load geojson files into spatial objects, the
geojsonsf
package would be easier to use.
library(geojsonR)
gj_anti = FROM_GeoJson(anti_json)
gj_correct = FROM_GeoJson(correct_json)
gj_cw = FROM_GeoJson(cw_json)
gj_anti
## $type
## [1] "Polygon"
##
## $coordinates
## [,1] [,2]
## [1,] -177.1 48.6
## [2,] 179.6 49.0
## [3,] 179.3 47.5
## [4,] -178.6 47.1
## [5,] -177.1 48.6
Use the SF package to turn the matrices into polygons
library(sf)
## Linking to GEOS 3.9.0, GDAL 3.2.2, PROJ 7.2.1; sf_use_s2() is TRUE
sf_anti = st_sfc(st_polygon(list(gj_anti$coordinates)), crs = st_crs(4326))
Check if SF reordered the coonrdinates
st_coordinates(sf_anti)
## X Y L1 L2
## [1,] -177.1 48.6 1 1
## [2,] 179.6 49.0 1 1
## [3,] 179.3 47.5 1 1
## [4,] -178.6 47.1 1 1
## [5,] -177.1 48.6 1 1
Throw the polygons on a web map
library(mapview)
mapview(sf_anti)
Apparently GeoJSON should
be split in two when crossing the antimeridian. sf
provides a function for that:
sf_split = st_wrap_dateline(sf_anti)
st_as_text(sf_split)
## [1] "MULTIPOLYGON (((179.3 47.5, 179.6 49, 180 48.95152, 180 47.36667, 179.3 47.5)), ((-180 48.95152, -177.1 48.6, -178.6 47.1, -180 47.36667, -180 48.95152)))"
mapview(sf_split)
It looks like the SF package does not reorder the coordinates:
sf_corr = st_sfc(st_polygon(list(gj_correct$coordinates)), crs = st_crs(4326))
sf_cw = st_sfc(st_polygon(list(gj_cw$coordinates)), crs = st_crs(4326))
st_coordinates(sf_corr)
## X Y L1 L2
## [1,] 30.12371 51.44400 1 1
## [2,] 30.12371 50.42971 1 1
## [3,] 31.72930 50.42971 1 1
## [4,] 31.72930 51.44400 1 1
## [5,] 30.12371 51.44400 1 1
st_coordinates(sf_cw)
## X Y L1 L2
## [1,] 30.12371 51.44400 1 1
## [2,] 31.72930 51.44400 1 1
## [3,] 31.72930 50.42971 1 1
## [4,] 30.12371 50.42971 1 1
## [5,] 30.12371 51.44400 1 1
But they plot the same, either way
mapview(sf_corr, col.regions = 'blue')
mapview(sf_cw, col.regions = 'red')
plot(sf_corr, lwd = 10, border = 'blue')
plot(sf_cw, col = 'red', add = TRUE)
Are these geometries considered to be valid?
st_is_valid(sf_anti)
## [1] TRUE
st_is_valid(sf_corr)
## [1] TRUE
st_is_valid(sf_cw)
## [1] TRUE