oce
oce copied to clipboard
possible solutions to the EOW and UHL problems
I have an idea that might help with the EOW (edge of world) problem, and the related UHL (ugly horizontal line) problem: just find the edge of the world, and use that to trim points and complete polygons by tracing the edge of the world.
This seems preferable to:
- looking for large xy jumps -- which is not very successful, partly because political boundaries could be straight
-
ad hoc additions like
coastlineCut-- which fails if the pole is in view, which means the cut point will not be at constant longitude, and which fails for Antarctic (well, that might be for another reason...)
The code below is a brute-force method that seems promising. The next step could be to fit a better function to the convex hull (I'm thinking a harmonic function, to get repeat when you get back to the first spot on the trace of the edge).
library(oce)
if (!interactive()) png("eoe.png")
par(mfrow=c(1, 1))
data(coastlineWorld)
LO <- seq(-360,360,10)
LA <- seq(-90,90,10)
g <- expand.grid(longitude=LO, latitude=LA)
zero <- list(lat=30, lon=-30)
projection <- sprintf("+proj=ortho +lat_0=30 +lon_0=-30", zero$lat, zero$lon)
xy <- lonlat2map(longitude=g$longitude, latitude=g$latitude, projection=projection)
ok <- is.finite(xy$x) & is.finite(xy$y)
par(mar=c(3, 3, 1, 1), mgp=c(2, 0.7, 0))
plot(xy$x, xy$y, col=ifelse(ok, 1, 2), cex=1/2)
## Show coastline, for reference
clxy <- lonlat2map(longitude=coastlineWorld[["longitude"]],
latitude=coastlineWorld[["latitude"]],
projection=projection)
lines(clxy$x, clxy$y, col="gray")
## find convex hull of finite points
ch <- chull(xy$x[ok], xy$y[ok])
hullx <- xy$x[ok][ch]
hullx <- c(hullx, hullx[1])
hully <- xy$y[ok][ch]
hully <- c(hully, hully[1])
points(hullx, hully, col=2, cex=1/2)
lines(hullx, hully, col=2, pch=20)
## Need a function to define the edge. Spline isn't
## much good. FIXME: think of a better function, perhaps
## harmonic, so it will repeat.
i <- seq_along(hullx)
sx <- smooth.spline(i, hullx)
sy <- smooth.spline(i, hully)
lines(predict(sx)$y, predict(sy)$y, col="blue", lwd=3)
if (!interactive()) dev.off()

The s2 package (and, eventually, an update sf package) may prove relevant to this issue and to #616 and #656; see https://www.r-spatial.org//r/2020/06/17/s2.html