# NOT RUN {
# create image collection from example Landsat data only
# if not already done in other examples
if (!file.exists(file.path(tempdir(), "L8.db"))) {
L8_files <- list.files(system.file("L8NY18", package = "gdalcubes"),
".TIF", recursive = TRUE, full.names = TRUE)
create_image_collection(L8_files, "L8_L1TP", file.path(tempdir(), "L8.db"))
}
L8.col = image_collection(file.path(tempdir(), "L8.db"))
v = cube_view(extent=list(left=388941.2, right=766552.4,
bottom=4345299, top=4744931, t0="2018-01", t1="2018-06"),
srs="EPSG:32618", nx = 497, ny=526, dt="P1M")
L8.cube = raster_cube(L8.col, v)
L8.rgb = select_bands(L8.cube, c("B02", "B03", "B04"))
L8.rgb.median = reduce_time(L8.rgb, "median(B02)", "median(B03)", "median(B04)")
L8.rgb.median
# }
# NOT RUN {
plot(L8.rgb.median, rgb=3:1)
# }
# NOT RUN {
# user defined reducer calculating interquartile ranges
L8.rgb.iqr = reduce_time(L8.rgb, names=c("iqr_R", "iqr_G","iqr_B"), FUN = function(x) {
c(diff(quantile(x["B04",],c(0.25,0.75), na.rm=TRUE)),
diff(quantile(x["B03",],c(0.25,0.75), na.rm=TRUE)),
diff(quantile(x["B02",],c(0.25,0.75), na.rm=TRUE)))
})
L8.rgb.iqr
# }
# NOT RUN {
plot(L8.rgb.iqr, key.pos=1)
# }
# NOT RUN {
# }
Run the code above in your browser using DataLab