Created
August 6, 2019 16:57
-
-
Save dbaston/c65ad0e7f800c4693bc6508fa3996a8f to your computer and use it in GitHub Desktop.
Cropping a raster from s3, masking using exactextract
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
library(exactextractr) | |
library(raster) | |
library(stars) | |
library(sf) | |
# define quad bounds using a GeoJSON string | |
quad_bounds <- st_read("{\n \"coordinates\": [\n [\n [\n [\n -98.97075, \n 46.7021\n ], \n [\n -98.6657, \n 47.44555\n ], \n [\n -97.48445, \n 47.27325\n ], \n [\n -97.7895, \n 46.5298\n ], \n [\n -98.97075, \n 46.7021\n ]\n ]\n ]\n ], \n \"type\": \"MultiPolygon\"\n}\n") | |
# create a stars proxy for the s3 object | |
scene <- read_stars('/vsis3/bucket-name-here/dataGlobal/landsat/8/031/027/LC08_L1TP_031027_20161107_20170219_01_T1/LC08_L1TP_031027_20161107_20170219_01_T1_TCT.tif', | |
proxy=TRUE) | |
# crop the stars proxy to the bbbox of the quad (cropping to the quad itself takes forever) | |
quad <- st_as_stars(scene[st_bbox(quad_bounds)]) | |
# convert to raster | |
quad.r <- as(quad, 'Raster') | |
# plot to check | |
plot(quad.r[[1]]) | |
plot(st_geometry(quad_bounds), add=TRUE) | |
# compute coverage fractions (0-1) | |
cov_frac <- coverage_fraction(quad.r[[1]], quad_bounds)[[1]] | |
# generate a mask based on pct of cell that must be within quad | |
mask <- (cov_frac >= 0.5) | |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment