[R-sig-Geo] Not perfect match between the bounding box for a EPSG when projected (st_bbox from sf)

Roger Bivand Roger@B|v@nd @end|ng |rom nhh@no
Fri Mar 26 15:47:01 CET 2021


On Thu, 25 Mar 2021, Manuel Spínola wrote:

> Dear list members,
>
> I am working with the epsg 5367.
>
> The information for this epsg according to this website (
> https://epsg.io/5367) is:
>
> WGS84 bounds:
> -90.45 2.15
> -81.43 11.77
>
> Projected bounds:
> -218963.01 239235.14
> 780173.56 1302793.77
>
> When I create a box with st_bbox from the package sf and transform the
> object I don't get the same bounding box that I was expecting.  Is it a
> fair expectation?

No. Firstly, epsg.io is not an EPSG service. Secondly, two different 
definitions are offered. I think that this works (sf 0.9-8):

pts1 <- st_sfc(list(st_point(c(-90.45, 2.15)), st_point(c(-81.43,
  11.77))), crs="OGC:CRS84")
st_transform(pts1, "EPSG:5367")
# POINT (-218962 239233.5)
# POINT (780173.7 1302793)
o <- sf_proj_pipelines("OGC:CRS84", "EPSG:5367")
o$definition
st_transform(pts1, "EPSG:5367", pipeline = o$definition[1])
# POINT (-218962 239233.5)
# POINT (780173.7 1302793)
# Costa Rica - onshore and offshore, accuracy 1.0 m, code 8969 [grid] 
st_transform(pts1, "EPSG:5367", pipeline = o$definition[2])
# POINT (-218963 239235.1)
# POINT (780173.6 1302794)
# Costa Rica - onshore and offshore, accuracy 1.0 m, code 5376 (default) 
# [grid]

So, as in https://github.com/r-spatial/sf/issues/1634, treating epsg.io as 
authoritative (especially when there are different transformation 
pipelines in play) may be unhelpful.sf/PROJ chooses a different pipeline 
as default than epsg.io.

This roundtrips:

st_transform(st_transform(pts1, "EPSG:5367"), "OGC:CRS84")
# POINT (-90.45 2.15)
# POINT (-81.43 11.77)
# Costa Rica - onshore and offshore, accuracy 1.0 m, code 8969 [grid]

but for epsg.io's default, one needs to intervene manually to choose 
pipelines both ways that are not sf's default choices.

Of course, here I am not using the bounding box polygon, because its 
projected version is not a rectangle, so its bounding box expands the 
range, especially further from the equator.

Hope this clarifies,

Roger

>
> b <- st_bbox(c(xmin = -90.45, xmax = -81.43, ymax = 11.77, ymin = 2.15),
> crs = st_crs(4326))
>
> b <- st_as_sfc(b)
>
> b_5367 <- st_transform(b, crs = 5367)
>
> b_5367
>
> Geometry set for 1 feature
> Geometry type: POLYGON
> Dimension:     XY
> Bounding box:  xmin: -218962 ymin: 237951.9 xmax: 785959.7 ymax: 1309622
> Projected CRS: CR05 / CRTM05POLYGON ((-218962 239233.5, 785959.7 237951.9, ...
>
>
>
>

-- 
Roger Bivand
Department of Economics, Norwegian School of Economics,
Postboks 3490 Ytre Sandviken, 5045 Bergen, Norway.
e-mail: Roger.Bivand using nhh.no
https://orcid.org/0000-0003-2392-6140
https://scholar.google.no/citations?user=AWeghB0AAAAJ&hl=en


More information about the R-sig-Geo mailing list