Geographic Information Systems Asked on June 11, 2021
I have a set of raster images (GeoTIFF, Landsat 1984-2018) which I cropped with my AOI using a shapefile.
Now I want to stack them and calculate the mean and standard deviation of every pixel of the stacked images (in space and time). How can I stack them and calculate the mean and standard deviation? For the outputfile I need one image for the mean and one for the standard deviation.
library(raster)
library(rgdal)
#load data
setwd("C:/Users/cathe/Documents/GEOTIFF_offsetcorrected")
dbase = "C:/Users/cathe/Documents/GEOTIFF_offsetcorrected"
#polygon with crop-extend
shape_data <- readOGR("C:/Users/cathe/OneDrive/Documents/ArcGIS/Projects/Cologne/WGS_1984_UTM_Zone_32N/Cologne.shp", stringsAsFactors=FALSE)
#load tif files
ALL_FILES <- list.files(path = dbase, pattern = "*.tif$|*.TIF$")
#output
outfiles = file.path("C:/Users/cathe/Documents/Output",
paste0(basename(tools::file_path_sans_ext(ALL_FILES)),
".tif"))
#crop
for (i in seq_along(ALL_FILES)) {
r = crop(stack(ALL_FILES[i]), shape_data)
writeRaster(r, filename=outfiles[i],
bylayer=FALSE,
format="GTiff",
datatype= "INT1U",
options="COMPRESS=ZIP",
overwrite=TRUE)
}
mean <- mean(r)
plot(mean)
writeRaster(x = mean, filename = "mean1.tif", driver = "GeoTiff", overwrite=TRUE)
I get a result, but I am not sure if it is correct?
I think you are looking for calc
function from library raster
. If you want a mean raster from your raster stack you can put fun = mean
, or if you want standard deviation you can put fun = sd
, etc.
For example, your rasters are r1
, r2
, r3
:
library(raster)
stacked <- stack(r1, r1, r3) # make a raster stack
# calculate a mean raster
meanR <- calc(stacked, fun = mean)
# make a new raster as standard deviation of stacked rasters
sdR <- calc(stacked, fun = sd)
# to calculate raster range
raR <- calc(stacked, fun = range) # this will give you a list of two rasters, one for lower range and one for upper range. To get the final range raster run the following code.
rangeR <- abs(raR[[1]] - raR[[2]])
Now you can save the output using writeRaster
function.
Answered by Tiny_hopper on June 11, 2021
Get help from others!
Recent Questions
Recent Answers
© 2024 TransWikia.com. All rights reserved. Sites we Love: PCI Database, UKBizDB, Menu Kuliner, Sharing RPP