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

Error in st_sample() with type = "hexagonal" for unprojected polygons #1220

Closed
Tracked by #15
dieghernan opened this issue Dec 19, 2019 · 1 comment
Closed
Tracked by #15

Comments

@dieghernan
Copy link

dieghernan commented Dec 19, 2019

It seems that there is a bug on st_sample when type="hexagonal" is performed over an unprojected polygon:

library(sf)
nc <- st_read(system.file("shape/nc.shp", package = "sf"))[1, ]
nc = st_transform(nc, 4326)
hexsample_4326 = st_sample(nc, size = 50, type = "hexagonal")

This returns this error:

Error in seq_len(nrow(xy)) : 
  argument must be coercible to non-negative integer
In addition: Warning message:
In seq_len(nrow(xy)) : first element used of 'length.out' argument

Whereas this works fine:

hexsample_other = st_sample(st_transform(nc, 3857), size = 50, type = "hexagonal")

I think the issue comes from this line

sf/R/sample.R

Line 101 in b96ffe8

a0 = as.numeric(st_area(st_make_grid(x, n = c(1,1))))

as it returns the area on m2, and then

sf/R/sample.R

Line 113 in b96ffe8

dx = sqrt(a0 / size / (sqrt(3)/2))

comes in meters as well, and finally

sf/R/sample.R

Lines 180 to 186 in b96ffe8

bb = st_bbox(obj)
dy = sqrt(3) * dx / 2
xlim = bb[c("xmin", "xmax")]
ylim = bb[c("ymin", "ymax")]
offset = c(x = (pt[1] - xlim[1]) %% dx, y = (pt[2] - ylim[1]) %% (2 * dy))
x = seq(xlim[1] - dx, xlim[2] + dx, dx) + offset[1]
y = seq(ylim[1] - 2 * dy, ylim[2] + 2 * dy, dy) + offset[2]

performs the operations mixing longitude and latitude with meters.

@edzer
Copy link
Member

edzer commented Jan 9, 2020

Yes, I can see this, but wonder whether this is better than giving an error message telling that users should first convert data to an equidistant rectangular (+proj=eqc) projection to get meaningful regular and hexagonal samples.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants