# NOT RUN {
#' # Define file names
ex_sel <- system.file(
"extdata/out/S2A2A_20170703_022_Barbellino_RGB432B_10.tif",
package = "sen2r"
)
ex_ref <- system.file(
"extdata/out/S2A2A_20170703_022_Barbellino_SCL_10.tif",
package = "sen2r"
)
crop_poly <- system.file("extdata/vector/dam.geojson", package = "sen2r")
crop_line <- sf::st_cast(sf::read_sf(crop_poly), "LINESTRING")
# Simple clip
test1 <- tempfile(fileext = "_test1.tif")
gdal_warp(ex_sel, test1, mask = crop_line)
# Clip and mask
test2 <- tempfile(fileext = "_test2.tif")
gdal_warp(ex_sel, test2, mask = crop_poly)
# Show output
crop_bbox <- sf::st_as_sfc(sf::st_bbox(crop_line))
oldpar <- par(mfrow = c(1,3), mar = rep(0,4))
image(stars::read_stars(ex_sel), rgb = 1:3)
plot(crop_line, add = TRUE, col = "blue", lwd = 2)
plot(crop_bbox, add = TRUE, border = "red", lwd = 2)
image(stars::read_stars(test1), rgb = 1:3)
plot(crop_bbox, add = TRUE, border = "red", lwd = 2)
image(stars::read_stars(test2), rgb = 1:3)
plot(crop_line, add = TRUE, col = "blue", lwd = 2)
# Warp on a reference raster
test3 <- tempfile(fileext = "_test3.tif")
gdal_warp(ex_sel, test3, ref = ex_ref)
# Show output
par(mfrow = c(1,3))
par(mar = rep(0,4)); image(stars::read_stars(ex_sel), rgb = 1:3)
par(mar = rep(2/3,4)); image(stars::read_stars(ex_ref))
par(mar = rep(0,4)); image(stars::read_stars(test3), rgb = 1:3)
# Reproject all the input file
test4 <- tempfile(fileext = "_test4.tif")
gdal_warp(ex_sel, test4, t_srs = "+init=epsg:32631")
# Reproject and clip on a bounding box
test5 <- tempfile(fileext = "_test5.tif")
gdal_warp(ex_sel, test5, t_srs = "+init=epsg:32631", mask = stars::read_stars(test1))
# Reproject and clip on polygon (masking outside)
test6 <- tempfile(fileext = "_test6.tif")
gdal_warp(ex_sel, test6, t_srs = "+init=epsg:32631", mask = crop_poly)
# Show output
crop_line_31N <- sf::st_transform(crop_line, 32631)
test1_bbox <- sf::st_as_sfc(sf::st_bbox(stars::read_stars(test1)))
test1_bbox_31N <- sf::st_transform(test1_bbox, 32631)
par(mfrow = c(1,4), mar = rep(0,4))
image(stars::read_stars(ex_sel), rgb = 1:3)
plot(crop_line, add = TRUE, col = "blue", lwd = 2)
plot(test1_bbox, add = TRUE, border = "red", lwd = 2)
image(stars::read_stars(test4), rgb = 1:3)
image(stars::read_stars(test5), rgb = 1:3)
plot(test1_bbox_31N, add = TRUE, border = "red", lwd = 2)
image(stars::read_stars(test6), rgb = 1:3)
plot(crop_line_31N, add = TRUE, col = "blue", lwd = 2)
# Use a reference raster with a different projection
test7 <- tempfile(fileext = "_test7.tif")
gdal_warp(ex_sel, test7, ref = test6)
# Use a reference raster with a different projection
# and specify a different bounding box
test8 <- tempfile(fileext = "_test8.tif")
gdal_warp(ex_sel, test8, mask = stars::read_stars(test1), ref = test6)
# Use a reference raster with a different projection and a mask
test9 <- tempfile(fileext = "_test9.tif")
gdal_warp(ex_sel, test9, mask = crop_poly, ref = test6)
# Show output
par(mfrow = c(1,4), mar = rep(0,4))
image(stars::read_stars(ex_sel), rgb = 1:3)
plot(crop_line, add = TRUE, col = "blue", lwd = 2)
image(stars::read_stars(test7), rgb = 1:3)
plot(crop_line_31N, add = TRUE, col = "blue", lwd = 2)
image(stars::read_stars(test8), rgb = 1:3)
plot(test1_bbox_31N, add = TRUE, border = "red", lwd = 2)
image(stars::read_stars(test9), rgb = 1:3)
plot(crop_line_31N, add = TRUE, col = "blue", lwd = 2)
par(oldpar)
# }
Run the code above in your browser using DataCamp Workspace