search

r-spatial / stars

Spatiotemporal Arrays, Raster and Vector Data Cubes
https://r-spatial.github.io/stars/
Apache License 2.0
553 stars 94 forks source link

st_rasterize Error in rep_len #671

Closed manomuthus closed 6 months ago

manomuthus commented 6 months ago

I am trying to convert some small grid sf objects into rasters using st_rasterize(sf_object). It works most of the time, but sometimes I get the following error. 

Error in rep_len(x, prod(dim)) : invalid 'length.out' value
 >R_ReplConsole(): before "for(;;)" {main.c

It's weird because there is hardly any difference between the sf objects I am using except some minor changes in the values e.g. No issues with this sf 

structure(list(haz = c(2.485, 1.625, 1.325, 0.98, 1.125, 2.765, 
2.105, 1.385, 1.405, 1.3050001, 3.17, 1.465, 0.76, 3.9699998, 
1.125, 0.08, 5.8450003, 3.92, 3.02, 2.245, 0.52), cell = structure(list(
    structure(list(structure(c(-121.8889815, -121.8889815, -121.8890741, 
    -121.8890741, -121.8889815, 37.6368519, 37.6367593, 37.6367593, 
    37.6368519, 37.6368519), dim = c(5L, 2L))), class = c("XY", 
    "POLYGON", "sfg")), structure(list(structure(c(-121.8889815, 
    -121.8889815, -121.8890741, -121.8890741, -121.8889815, 37.6369444, 
    37.6368519, 37.6368519, 37.6369444, 37.6369444), dim = c(5L, 
    2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
        structure(c(-121.8889815, -121.8889815, -121.8890741, 
        -121.8890741, -121.8889815, 37.637037, 37.6369444, 37.6369444, 
        37.637037, 37.637037), dim = c(5L, 2L))), class = c("XY", 
    "POLYGON", "sfg")), structure(list(structure(c(-121.8889815, 
    -121.8889815, -121.8890741, -121.8890741, -121.8889815, 37.6371296, 
    37.637037, 37.637037, 37.6371296, 37.6371296), dim = c(5L, 
    2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
        structure(c(-121.8889815, -121.8889815, -121.8890741, 
        -121.8890741, -121.8889815, 37.6372222, 37.6371296, 37.6371296, 
        37.6372222, 37.6372222), dim = c(5L, 2L))), class = c("XY", 
    "POLYGON", "sfg")), structure(list(structure(c(-121.8890741, 
    -121.8890741, -121.8891667, -121.8891667, -121.8890741, 37.6368519, 
    37.6367593, 37.6367593, 37.6368519, 37.6368519), dim = c(5L, 
    2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
        structure(c(-121.8890741, -121.8890741, -121.8891667, 
        -121.8891667, -121.8890741, 37.6369444, 37.6368519, 37.6368519, 
        37.6369444, 37.6369444), dim = c(5L, 2L))), class = c("XY", 
    "POLYGON", "sfg")), structure(list(structure(c(-121.8890741, 
    -121.8890741, -121.8891667, -121.8891667, -121.8890741, 37.637037, 
    37.6369444, 37.6369444, 37.637037, 37.637037), dim = c(5L, 
    2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
        structure(c(-121.8890741, -121.8890741, -121.8891667, 
        -121.8891667, -121.8890741, 37.6371296, 37.637037, 37.637037, 
        37.6371296, 37.6371296), dim = c(5L, 2L))), class = c("XY", 
    "POLYGON", "sfg")), structure(list(structure(c(-121.8890741, 
    -121.8890741, -121.8891667, -121.8891667, -121.8890741, 37.6372222, 
    37.6371296, 37.6371296, 37.6372222, 37.6372222), dim = c(5L, 
    2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
        structure(c(-121.8891667, -121.8891667, -121.8892593, 
        -121.8892593, -121.8891667, 37.6368519, 37.6367593, 37.6367593, 
        37.6368519, 37.6368519), dim = c(5L, 2L))), class = c("XY", 
    "POLYGON", "sfg")), structure(list(structure(c(-121.8891667, 
    -121.8891667, -121.8892593, -121.8892593, -121.8891667, 37.6371296, 
    37.637037, 37.637037, 37.6371296, 37.6371296), dim = c(5L, 
    2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
        structure(c(-121.8891667, -121.8891667, -121.8892593, 
        -121.8892593, -121.8891667, 37.6372222, 37.6371296, 37.6371296, 
        37.6372222, 37.6372222), dim = c(5L, 2L))), class = c("XY", 
    "POLYGON", "sfg")), structure(list(structure(c(-121.8892593, 
    -121.8892593, -121.8893519, -121.8893519, -121.8892593, 37.6368519, 
    37.6367593, 37.6367593, 37.6368519, 37.6368519), dim = c(5L, 
    2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
        structure(c(-121.8892593, -121.8892593, -121.8893519, 
        -121.8893519, -121.8892593, 37.6371296, 37.637037, 37.637037, 
        37.6371296, 37.6371296), dim = c(5L, 2L))), class = c("XY", 
    "POLYGON", "sfg")), structure(list(structure(c(-121.8892593, 
    -121.8892593, -121.8893519, -121.8893519, -121.8892593, 37.6372222, 
    37.6371296, 37.6371296, 37.6372222, 37.6372222), dim = c(5L, 
    2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
        structure(c(-121.8893519, -121.8893519, -121.8894444, 
        -121.8894444, -121.8893519, 37.6368519, 37.6367593, 37.6367593, 
        37.6368519, 37.6368519), dim = c(5L, 2L))), class = c("XY", 
    "POLYGON", "sfg")), structure(list(structure(c(-121.8893519, 
    -121.8893519, -121.8894444, -121.8894444, -121.8893519, 37.6369444, 
    37.6368519, 37.6368519, 37.6369444, 37.6369444), dim = c(5L, 
    2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
        structure(c(-121.8893519, -121.8893519, -121.8894444, 
        -121.8894444, -121.8893519, 37.637037, 37.6369444, 37.6369444, 
        37.637037, 37.637037), dim = c(5L, 2L))), class = c("XY", 
    "POLYGON", "sfg")), structure(list(structure(c(-121.8893519, 
    -121.8893519, -121.8894444, -121.8894444, -121.8893519, 37.6371296, 
    37.637037, 37.637037, 37.6371296, 37.6371296), dim = c(5L, 
    2L))), class = c("XY", "POLYGON", "sfg")), structure(list(
        structure(c(-121.8893519, -121.8893519, -121.8894444, 
        -121.8894444, -121.8893519, 37.6372222, 37.6371296, 37.6371296, 
        37.6372222, 37.6372222), dim = c(5L, 2L))), class = c("XY", 
    "POLYGON", "sfg"))), class = c("sfc_POLYGON", "sfc"), precision = 0, bbox = structure(c(xmin = -121.8894444, 
ymin = 37.6367593, xmax = -121.8889815, ymax = 37.6372222), class = "bbox"), crs = structure(list(
    input = "EPSG:4326", wkt = "GEOGCRS[\"WGS 84\",\n    ENSEMBLE[\"World Geodetic System 1984 ensemble\",\n        MEMBER[\"World Geodetic System 1984 (Transit)\"],\n        MEMBER[\"World Geodetic System 1984 (G730)\"],\n        MEMBER[\"World Geodetic System 1984 (G873)\"],\n        MEMBER[\"World Geodetic System 1984 (G1150)\"],\n        MEMBER[\"World Geodetic System 1984 (G1674)\"],\n        MEMBER[\"World Geodetic System 1984 (G1762)\"],\n        MEMBER[\"World Geodetic System 1984 (G2139)\"],\n        ELLIPSOID[\"WGS 84\",6378137,298.257223563,\n            LENGTHUNIT[\"metre\",1]],\n        ENSEMBLEACCURACY[2.0]],\n    PRIMEM[\"Greenwich\",0,\n        ANGLEUNIT[\"degree\",0.0174532925199433]],\n    CS[ellipsoidal,2],\n        AXIS[\"geodetic latitude (Lat)\",north,\n            ORDER[1],\n            ANGLEUNIT[\"degree\",0.0174532925199433]],\n        AXIS[\"geodetic longitude (Lon)\",east,\n            ORDER[2],\n            ANGLEUNIT[\"degree\",0.0174532925199433]],\n    USAGE[\n        SCOPE[\"Horizontal component of 3D system.\"],\n        AREA[\"World.\"],\n        BBOX[-90,-180,90,180]],\n    ID[\"EPSG\",4326]]"), class = "crs"), n_empty = 0L)), row.names = c(NA, 
-21L), sf_column = "cell", agr = structure(c(haz = NA_integer_), levels = c("constant", 
"aggregate", "identity"), class = "factor"), class = c("sf", 
"tbl_df", "tbl", "data.frame"))

this one returns the above error

structure(list(haz = c(2.485, 1.625, 1.325, 0.98, 1.125, 2.765, 
2.105, 1.385, 1.405, 1.3050001, 3.17, 1.465, 0.76, 3.9699998, 
1.125, 0.08, 5.8450003, 3.92, 3.02, 2.245, 0.52), cell = structure(list(
    structure(list(structure(c(-121.9, -121.9, -121.9, -121.9, 
    -121.9, 37.64, 37.64, 37.64, 37.64, 37.64), dim = c(5L, 2L
    ))), class = c("XY", "POLYGON", "sfg")), structure(list(structure(c(-121.9, 
    -121.9, -121.9, -121.9, -121.9, 37.64, 37.64, 37.64, 37.64, 
    37.64), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg"
    )), structure(list(structure(c(-121.9, -121.9, -121.9, -121.9, 
    -121.9, 37.64, 37.64, 37.64, 37.64, 37.64), dim = c(5L, 2L
    ))), class = c("XY", "POLYGON", "sfg")), structure(list(structure(c(-121.9, 
    -121.9, -121.9, -121.9, -121.9, 37.64, 37.64, 37.64, 37.64, 
    37.64), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg"
    )), structure(list(structure(c(-121.9, -121.9, -121.9, -121.9, 
    -121.9, 37.64, 37.64, 37.64, 37.64, 37.64), dim = c(5L, 2L
    ))), class = c("XY", "POLYGON", "sfg")), structure(list(structure(c(-121.9, 
    -121.9, -121.9, -121.9, -121.9, 37.64, 37.64, 37.64, 37.64, 
    37.64), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg"
    )), structure(list(structure(c(-121.9, -121.9, -121.9, -121.9, 
    -121.9, 37.64, 37.64, 37.64, 37.64, 37.64), dim = c(5L, 2L
    ))), class = c("XY", "POLYGON", "sfg")), structure(list(structure(c(-121.9, 
    -121.9, -121.9, -121.9, -121.9, 37.64, 37.64, 37.64, 37.64, 
    37.64), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg"
    )), structure(list(structure(c(-121.9, -121.9, -121.9, -121.9, 
    -121.9, 37.64, 37.64, 37.64, 37.64, 37.64), dim = c(5L, 2L
    ))), class = c("XY", "POLYGON", "sfg")), structure(list(structure(c(-121.9, 
    -121.9, -121.9, -121.9, -121.9, 37.64, 37.64, 37.64, 37.64, 
    37.64), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg"
    )), structure(list(structure(c(-121.9, -121.9, -121.9, -121.9, 
    -121.9, 37.64, 37.64, 37.64, 37.64, 37.64), dim = c(5L, 2L
    ))), class = c("XY", "POLYGON", "sfg")), structure(list(structure(c(-121.9, 
    -121.9, -121.9, -121.9, -121.9, 37.64, 37.64, 37.64, 37.64, 
    37.64), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg"
    )), structure(list(structure(c(-121.9, -121.9, -121.9, -121.9, 
    -121.9, 37.64, 37.64, 37.64, 37.64, 37.64), dim = c(5L, 2L
    ))), class = c("XY", "POLYGON", "sfg")), structure(list(structure(c(-121.9, 
    -121.9, -121.9, -121.9, -121.9, 37.64, 37.64, 37.64, 37.64, 
    37.64), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg"
    )), structure(list(structure(c(-121.9, -121.9, -121.9, -121.9, 
    -121.9, 37.64, 37.64, 37.64, 37.64, 37.64), dim = c(5L, 2L
    ))), class = c("XY", "POLYGON", "sfg")), structure(list(structure(c(-121.9, 
    -121.9, -121.9, -121.9, -121.9, 37.64, 37.64, 37.64, 37.64, 
    37.64), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg"
    )), structure(list(structure(c(-121.9, -121.9, -121.9, -121.9, 
    -121.9, 37.64, 37.64, 37.64, 37.64, 37.64), dim = c(5L, 2L
    ))), class = c("XY", "POLYGON", "sfg")), structure(list(structure(c(-121.9, 
    -121.9, -121.9, -121.9, -121.9, 37.64, 37.64, 37.64, 37.64, 
    37.64), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg"
    )), structure(list(structure(c(-121.9, -121.9, -121.9, -121.9, 
    -121.9, 37.64, 37.64, 37.64, 37.64, 37.64), dim = c(5L, 2L
    ))), class = c("XY", "POLYGON", "sfg")), structure(list(structure(c(-121.9, 
    -121.9, -121.9, -121.9, -121.9, 37.64, 37.64, 37.64, 37.64, 
    37.64), dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg"
    )), structure(list(structure(c(-121.9, -121.9, -121.9, -121.9, 
    -121.9, 37.64, 37.64, 37.64, 37.64, 37.64), dim = c(5L, 2L
    ))), class = c("XY", "POLYGON", "sfg"))), class = c("sfc_POLYGON", 
"sfc"), precision = 0, bbox = structure(c(xmin = -121.9, ymin = 37.64, 
xmax = -121.9, ymax = 37.64), class = "bbox"), crs = structure(list(
    input = "EPSG:4326", wkt = "GEOGCRS[\"WGS 84\",\n    ENSEMBLE[\"World Geodetic System 1984 ensemble\",\n        MEMBER[\"World Geodetic System 1984 (Transit)\"],\n        MEMBER[\"World Geodetic System 1984 (G730)\"],\n        MEMBER[\"World Geodetic System 1984 (G873)\"],\n        MEMBER[\"World Geodetic System 1984 (G1150)\"],\n        MEMBER[\"World Geodetic System 1984 (G1674)\"],\n        MEMBER[\"World Geodetic System 1984 (G1762)\"],\n        MEMBER[\"World Geodetic System 1984 (G2139)\"],\n        ELLIPSOID[\"WGS 84\",6378137,298.257223563,\n            LENGTHUNIT[\"metre\",1]],\n        ENSEMBLEACCURACY[2.0]],\n    PRIMEM[\"Greenwich\",0,\n        ANGLEUNIT[\"degree\",0.0174532925199433]],\n    CS[ellipsoidal,2],\n        AXIS[\"geodetic latitude (Lat)\",north,\n            ORDER[1],\n            ANGLEUNIT[\"degree\",0.0174532925199433]],\n        AXIS[\"geodetic longitude (Lon)\",east,\n            ORDER[2],\n            ANGLEUNIT[\"degree\",0.0174532925199433]],\n    USAGE[\n        SCOPE[\"Horizontal component of 3D system.\"],\n        AREA[\"World.\"],\n        BBOX[-90,-180,90,180]],\n    ID[\"EPSG\",4326]]"), class = "crs"), n_empty = 0L)), row.names = c(NA, 
-21L), sf_column = "cell", agr = structure(c(haz = NA_integer_), levels = c("constant", 
"aggregate", "identity"), class = "factor"), class = c("sf", 
"tbl_df", "tbl", "data.frame"))

I am using version 0.6-4

edzer commented 6 months ago

Copy & paste for some reason doesn't work; could you provide the data e.g. as a zipped .rda file?

manomuthus commented 6 months ago

here you go sf_files.zip

edzer commented 6 months ago 
> st_is_valid(sf_with_issue)
 [1] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
[13] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE

don't expect too many things to work with invalid geometries.

manomuthus commented 6 months ago

ahh..not sure how I missed that. Sorry. The error message is not useful either. Closing this now