-
Notifications
You must be signed in to change notification settings - Fork 228
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Add downscale_function and covariates scripts #3272
Changes from all commits
6011bf6
1a450cc
9ffd4df
f147849
84d1147
668540c
eb371cf
6b57164
3f20e8c
ed8394f
71040be
95fb57f
File filter
Filter by extension
Conversations
Jump to
Diff view
Diff view
There are no files selected for viewing
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,109 @@ | ||
##' @title Covariate Data Prep | ||
##' @author Joshua Ploshay | ||
##' | ||
##' @description This script lists the code needed to download, resample, and | ||
##' stack the covariate data used in the downscale_function.R sript. | ||
##' | ||
##' @return A spatraster object consisting of a stack of maps with labeled layers | ||
|
||
#### WorldClim #### | ||
# WorldClim v2.1 (1970-2000) historical climate data | ||
# Data can be found at https://www.worldclim.org/data/worldclim21.html | ||
# Product is a zip file containing 12 .tif files, one for each month | ||
# Use code below to average the 12 files to obtain one map | ||
# 2023 REU project used 10 minute spatial resolution | ||
|
||
## Solar Radiation (kJ m-2 day-1) | ||
srad <- terra::rast(list.files(path = "/projectnb/dietzelab/jploshay/pecan_copy/jploshay/10m_srad", | ||
pattern='.tif$', | ||
all.files= T, | ||
full.names= T)) | ||
srad <- terra::app(srad, mean) | ||
|
||
|
||
## Vapor Pressure (kPa) | ||
vapr <- terra::rast(list.files(path = "/projectnb/dietzelab/jploshay/pecan_copy/jploshay/10m_vapr", | ||
pattern='.tif$', | ||
all.files= T, | ||
full.names= T)) | ||
vapr <- terra::app(vapr, mean) | ||
|
||
|
||
## Average Temperature (*C) | ||
tavg <- terra::rast(list.files(path = "/projectnb/dietzelab/jploshay/pecan_copy/jploshay/avg_temp_prep/WorldClim", | ||
pattern='.tif$', | ||
all.files= T, | ||
full.names= T)) | ||
tavg <- terra::app(tavg, mean) | ||
|
||
|
||
## Total Precipitation (mm) | ||
prec <- terra::rast(list.files(path = "/projectnb/dietzelab/jploshay/pecan_copy/jploshay/total_prec", | ||
pattern='.tif$', | ||
all.files= T, | ||
full.names= T)) | ||
prec <- terra::app(prec, mean) | ||
|
||
|
||
#### SoilGrids #### | ||
# geodata::soil_world() pulls data from the SoilGRIDS database | ||
# More information on pulling different soil data can be found in the geodata | ||
# manual https://cran.r-project.org/web/packages/geodata/geodata.pdf | ||
|
||
## Soil pH | ||
phh2o <- geodata::soil_world(var = "phh2o", depth = 5, stat = "mean", path = tempdir()) | ||
|
||
## Soil Nitrogen (g kg-1) | ||
nitrogen <- geodata::soil_world(var = "nitrogen", depth = 5, stat = "mean", path = tempdir()) | ||
|
||
## Soil Organic Carbon (g kg-1) | ||
soc <- geodata::soil_world(var = "soc", depth = 5, stat = "mean", path = tempdir()) | ||
|
||
## Soil Sand (%) | ||
sand <- geodata::soil_world(var = "sand", depth = 5, stat = "mean", path = tempdir()) | ||
|
||
|
||
#### Land Cover #### | ||
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. For files that are manually downloaded, your earlier examples provided some indication of where the data is located and how to download it. Equivalent info is missing here. |
||
GLanCE_extract <- function(pattern, path) { | ||
files <- list.files(path = "/projectnb/dietzelab/dietze/glance2012/e4ftl01.cr.usgs.gov/MEASURES/GLanCE30.001/2012.07.01", #make this path default | ||
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. I'm pretty sure you don't want to hard code a file path within a function. Much better to make it an argument. indeed, this function already has a |
||
all.files = T, | ||
full.names = T, | ||
pattern) | ||
# empty .vrt file path | ||
vrtfile <- paste0(tempfile(), ".vrt") | ||
# "connects" tiles together that returns SpatRaster | ||
GLanCE_product <- terra::vrt(files, vrtfile, overwrite = T) | ||
return(GLanCE_product) | ||
} | ||
|
||
# Integer identifier for class in the current year | ||
land_cover <- GLanCE_extract(pattern = "NA_LC.tif$") | ||
|
||
|
||
#### Resample and Stack Maps #### | ||
# Define the extent to crop the covariates to North America | ||
NA_extent <- terra::ext(-178.19453125, -10, 7.22006835937502, 83.5996093750001) | ||
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. This would be good to have at the top of the script in case a user wants/needs to set a different extent. |
||
|
||
# Stack WorldClim maps | ||
WorldClim <- c(tavg, srad, prec, vapr) | ||
|
||
# Crop WolrdClim stack to North America using North America extent | ||
NA_WorldClim <- terra::crop(WorldClim, NA_extent) | ||
names(NA_WorldClim) <- c("tavg", "srad", "prec", "vapr") | ||
|
||
# Stack SoilGrids maps | ||
SoilGrids <- c(phh2o, nitrogen, soc, sand) | ||
|
||
# Crop SoilGrids stack to North America using North America extent | ||
NA_SoilGrids <- terra::crop(SoilGrids, NA_extent) | ||
names(NA_SoilGrids) <- c("phh2o", "nitrogen", "soc", "sand") | ||
|
||
# Resample SoilGrids to match WorldClim maps | ||
NA_SoilGrids <- terra::resample(NA_WorldClim, NA_SoilGrids, method = 'bilinear') # made change | ||
|
||
# Resample land cover to match WorldClim maps (~25 min) | ||
land_cover <- terra::resample(land_cover, NA_WorldClim, method = 'near') | ||
names(land_cover) <- "land_cover" | ||
|
||
# Stack all maps | ||
covariates <- c(NA_WorldClim, NA_SoilGrids, land_cover) |
Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
For anyone to be able to reuse this script they will have to redefine a LOT of local paths. Better would be to define ONE folder where all the data is located at the top of the script and then express all the other paths relative to that (e.g. using
file.path()
)