Close polygons using linestrings and selection of polygons created by area

Question

I'd like to close all the polygons using linestrings and the selection of polygons by area. In my example:

# Packages
library(stars)
library(sf)
library(ggplot2)

#Vectorizing a raster object to an sf object
tif = system.file("tif/L7_ETMs.tif", package = "stars")
# Let's stay with only the first band, indicated in the final dimension
x = read_stars(tif)[, 1:50, 1:50, 1]
x = round(x/5)
# Calculate the min and max of the raster values
purrr::map(x, min)
# 10
purrr::map(x, max)
# 28
# Change values lower than 15 for 1 or 2
x[[1]] = ifelse(x[[1]]<15,1,2)

# Take a look at the image
plot(x)

# Obtain the contours
l =  st_contour(x, 
                # Remove NA
                na.rm = T,
                # Obtain contour lines instead of polygons
                contour_lines = TRUE,
                # raster values at which to draw contour levels
                breaks = 1.5)

# Plot the contours
raster_line <- st_as_sf(l, coords = c("x", "y"), crs = 32725) %>%
  st_cast("LINESTRING")
ggplot(raster_line) +
  geom_sf(aes(color =  L7_ETMs.tif), show.legend = "line")  +
  theme_bw() +
  theme(
   panel.grid.major = element_blank(),
   panel.grid.minor = element_blank(),
   legend.position = "none")
#

Now, I'd like to create new contours just only if the polygons are greater than 900 m^2, but in the first step I try to make the conversion on "LINESTRING" to "MULTIPOLYGON" (close polygons) using directly st_cast() and doesn't work. My final goal is to close the polygons, make the selections of polygons if it area is > 900m^2? Please, any help with it? Thanks in advance

Answer 1

Assuming your example demonstrates the workflow you have, you could also directly polygonize the stars raster:

library(stars)
#> Loading required package: abind
#> Loading required package: sf
#> Linking to GEOS 3.8.1, GDAL 3.2.1, PROJ 7.2.1; sf_use_s2() is TRUE
library(sf)
library(dplyr)
#> 
#> Attaching package: 'dplyr'
#> The following objects are masked from 'package:stats':
#> 
#>     filter, lag
#> The following objects are masked from 'package:base':
#> 
#>     intersect, setdiff, setequal, union

#Vectorizing a raster object to an sf object
tif = system.file("tif/L7_ETMs.tif", package = "stars")
# Let's stay with only the first band, indicated in the final dimension
x = read_stars(tif)[, 1:50, 1:50, 1]
x = round(x/5)

x %>% 
  cut(., breaks = c(14,Inf)) %>% 
  st_as_sf(., as_points = FALSE, merge = TRUE) %>% 
  mutate(area = st_area(.)) %>% 
  filter(area > units::set_units(900, "m^2", mode = "standard")) %>% 
  dplyr::select(-area) %>% 
  plot()

^{Created on 2022-01-13 by the reprex package (v2.0.1)}

Answer 2

# Convert to polygon and multipolygons
p <- st_cast(raster_line,"POLYGON")
mult.p <- st_cast(p,"MULTIPOLYGON")

# Calculate polygons area
mult.p$area <- st_area(mult.p)

# Selection of polygons on interest
mult.p.target <-  mult.p %>% 
        dplyr::filter(as.numeric(mult.p$area) > 900) 
#

# Plot the new contours
raster_poly <- st_as_sf(mult.p.target, coords = c("x", "y"), crs = 32725) %>%
  st_cast("MULTIPOLYGON")
ggplot(raster_poly) +
  geom_sf(aes(color =  L7_ETMs.tif))  +
  theme_bw() +
  theme(
   panel.grid.major = element_blank(),
   panel.grid.minor = element_blank(),
   legend.position = "none")
#