Adding the output of sessionInfo() and sf_extSoftVersion() might also help. The difficulty Zivan points out is that the structure as read is valid, but before conversion to text by dput() may have actually had non-identical first and last coordinates, because of the number of digits used. Also try putting the file created by saveRDS(aoi, ...) on a website with a link here. It would be useful to see the code in the workflow at the point of failure (you point to Error in MtrxSet, maybe also run traceback() after the failure to give more precision in where this is occurring, maybe in stars?
Hop this helps,
Roger
--
Roger Bivand
Emeritus Professor
Norwegian School of Economics
Postboks 3490 Ytre Sandviken, 5045 Bergen, Norway
roger.biv...@nhh.no
________________________________________
From: R-sig-Geo <r-sig-geo-boun...@r-project.org> on behalf of Zivan Karaman <zivan.kara...@gmail.com>
Sent: 19 August 2024 18:40
To: Micha Silver
Cc: R-sig-Geo@r-project.org
Subject: Re: [R-sig-Geo] Error in MtrxSet ... polygons not (all) closed
Could you please provide an example of code that actually raises this error?
Best regards,
Zivan
On 19 Aug 2024, at 16:23, Micha Silver <tsvi...@gmail.com> wrote:
Searching for the error in the subject line returns some discussions from a few years ago. But I'm not able to overcome this error.
Here's the polygon in question:
dput(aoi)
structure(list(x =
structure(list(structure(list(structure(c(11.0127043, 11.0127061, 11.0134526, 11.0134554, 11.0134626, 11.0134692, 11.0134776, 11.0134813, 11.0157343, 11.0157373, 11.0157437, 11.0157495, 11.0157546, 11.0157585, 11.015761, 11.0157617, 11.0157616, 11.0157613, 11.0153992, 11.0153983, 11.0153956, 11.0153922, 11.0153869, 11.0153804, 11.0153729, 11.0153694, 11.0127478, 11.0127447, 11.0127379, 11.012732, 11.0127269, 11.0127232, 11.0127212, 11.01272, 11.0127197, 11.0126956, 11.0126958, 11.0126965, 11.0126978, 11.0127004, 11.0127043, 46.8484069, 46.8484048, 46.8476584, 46.847656, 46.8476511, 46.8476481, 46.8476461, 46.8476456, 46.8475213, 46.8475214, 46.847522, 46.8475235, 46.8475266, 46.8475311, 46.8475365, 46.8475425, 46.8475489, 46.8475519, 46.8495623, 46.8495656, 46.8495729, 46.8495788, 46.8495833, 46.8495856, 46.849587, 46.8495874, 46.8496081, 46.8496079, 46.8496069, 46.849605, 46.8496014, 46.8495965, 46.8495906, 46.8495839, 46.8495808, 46.8484309, 46.8484282, 46.8484218, 46.8484166, 46.848412, 46.8484069), dim = c(41L, 2L), dimnames = list( NULL, c("X", "Y")))), class = c("XY", "POLYGON", "sfg"))), class = c("sfc_POLYGON", "sfc"), precision = 0, bbox = structure(c(xmin = 11.0126956, ymin = 46.8475213, xmax = 11.0157617, ymax = 46.8496081), 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 = 1L, class = c("sf", "data.frame"), sf_column = "x", agr = structure(integer(0), class = "factor", levels = c("constant", "aggregate", "identity"), names = character(0)))
This polygon is used as area of interest to crop Copernicus imagery that I am downloading using the {CDSE} package. (I have run this workflow successfully many times with several other aoi polygons)
I have tried:
1- sf::st_make_valid()
2- transforming to a UTM CRS
3- extracting coordinates and recreating the polygon (i.e.:
aoi_p <- st_polygon(list(st_coordinates(aoi)[,1:2]))
aoi3 <- aoi_p |>
st_sfc() |>
st_as_sf()
st_crs(aoi3) <- "EPSG:4326"
4- buffering by a small amount
The above error recurs in all cases (only with this problem polygon). Any suggestions?
Thanks
Micha Silver
Ben Gurion Univ.
Sde Boker, Remote Sensing Lab
cell: +972-523-665918
https://orcid.org/0000-0002-1128-1325
--
Micha Silver
Ben Gurion Univ.
Sde Boker, Remote Sensing Lab
cell: +972-523-665918
_______________________________________________
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
_______________________________________________
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
