Thanks for the update. Indeed, I have noticed that this polygon is on the edge of Sentinel-2 tiles. As a result, the AOI gets covered with two different tiles on every overpass, but not completely (bad luck for a 4.5 ha field!). Even the images that display 100% area coverage actually miss a tiny, tiny fraction of it (approximately 4e-10). catalog <- SearchCatalog(aoi = aoi, from = "2023-06-01", to = "2023-08-31", collection = "sentinel-2-l2a”, with_geometry = TRUE, client = OAuthClient) 100 - subset(catalog, areaCoverage > 99)$areaCoverage ## [1] 4.085194e-10 4.085194e-10 4.085194e-10 4.085194e-10 4.085194e-10 ... I'm not sure if this could be the source of the problem, but it's good that you found a workaround. Best regards, Zivan
On 21 Aug 2024, at 18:59, Micha Silver <tsvi...@gmail.com> wrote: Follow up: After some gnawing of my teeth, I found the problem: Nothing to do with R, sf, nor the CDSE package. My date range included one date where the Sentinel-2 tile from Copernicus data space was damaged/clipped/mangled. And my study area was "hanging over the edge". I found the problem image using the Copernicus Dataspace browser, and going thru date by date... I'm not sure why exactly that caused the error, but once I changed the date range to exclude that problem image, the procedure worked as expected. Thanks for your help. On 19/08/2024 21:29, Roger Bivand wrote: > 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 <mailto:roger.biv...@nhh.no> > > ________________________________________ > From: R-sig-Geo <r-sig-geo-boun...@r-project.org> > <mailto:r-sig-geo-boun...@r-project.org> on behalf of Zivan Karaman > <zivan.kara...@gmail.com> <mailto:zivan.kara...@gmail.com> > Sent: 19 August 2024 18:40 > To: Micha Silver > Cc: R-sig-Geo@r-project.org <mailto: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> > <mailto: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 <mailto: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 <mailto:R-sig-Geo@r-project.org> > https://stat.ethz.ch/mailman/listinfo/r-sig-geo -- Micha Silver Ben Gurion Univ. Sde Boker, Remote Sensing Lab cell: +972-523-665918 [[alternative HTML version deleted]] _______________________________________________ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo
