Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

possible solutions to the EOW and UHL problems #1307

Open
dankelley opened this issue Oct 1, 2017 · 1 comment
Open

possible solutions to the EOW and UHL problems #1307

dankelley opened this issue Oct 1, 2017 · 1 comment
Assignees
Labels
development low priority mapping pinned Prevent automatic stale designation

Comments

@dankelley
Copy link
Owner

dankelley commented Oct 1, 2017

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:

  1. looking for large xy jumps -- which is not very successful, partly because political boundaries could be straight
  2. 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()

eoe

@dankelley dankelley self-assigned this Oct 1, 2017
@dankelley dankelley added the pinned Prevent automatic stale designation label Sep 22, 2018
@dankelley
Copy link
Owner Author

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

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
development low priority mapping pinned Prevent automatic stale designation
Projects
None yet
Development

No branches or pull requests

1 participant