Skip to content
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

Hg gdal rasterize fix #34

Open
wants to merge 10 commits into
base: master
Choose a base branch
from
2 changes: 1 addition & 1 deletion DESCRIPTION
Original file line number Diff line number Diff line change
Expand Up @@ -22,7 +22,7 @@ Depends:
raster,
ModelMap (>= 3.4.0.1),
SystemRequirements: SAGA-GIS (v2.3+)
RoxygenNote: 7.1.1
RoxygenNote: 7.1.2
Suggests:
knitr,
rmarkdown
Expand Down
73 changes: 56 additions & 17 deletions R/raster_stack.R
Original file line number Diff line number Diff line change
Expand Up @@ -29,6 +29,13 @@
#' with any of the following: c("SLOPE", "ASPECT", "DAH", "MRVBF", "TPI",
#' "CPLAN", "CPROF", "TOPOWET", "CAREA"). If not provided, this parameter
#' defaults to generating all products.
#' @param param_vect vector (character). Similar to the \code{products} parameter, this is a vector of
#' the SAGA module names (see Details section below) for which non-default
#' parameters will be passed. Defaults to NULL.
#' @param param_values This is a list of lists, and the order corresponding to \code{param_vect}. For each
#' parameter specified in \code{param_vect}, the additional parameters for the corresponding
#' SAGA module must be specified using a list of those additional SAGA parameters and their values.
#' Defaults to NULL.
#'
#' @return NULL
#'
Expand All @@ -38,7 +45,15 @@
#' products = c("SLOPE", "CAREA"))
#' }
#' @export
create_dem_products <- function(dem,stream_vec = NULL,burn_val=NULL,outdir, products = NULL) {
create_dem_products <- function(dem,stream_vec = NULL,burn_val=NULL,outdir, products = NULL,param_vect=NULL,param_values=NULL) {

if(is.null(param_vect))
{
param_vect<-c(NULL)
}else{
param_values<-append(list(NULL),param_values)
}

env <- RSAGA::rsaga.env()

#Get dem raster params
Expand All @@ -62,11 +77,17 @@ create_dem_products <- function(dem,stream_vec = NULL,burn_val=NULL,outdir, prod
{
streams_r <- file.path(tempdir(), "streams_r.tif")
lyr_name <- tools::file_path_sans_ext(basename(stream_vec))
cmd <- paste("gdal_rasterize -l ", lyr_name, " -burn 1 -tr ",
x_res, " ", y_res, " -te ", x_min, " ", y_min, " ",
x_max, " ", y_max, " -ot Byte ", stream_vec, " ",
streams_r, sep = "")
system(cmd)

gdalUtils::gdal_rasterize(src_datasource = stream_vec,
dst_filename = streams_r,
burn = 1,
l = lyr_name,
tr = c(x_res,y_res),
te = c(x_min,y_min,x_max,y_max),
ot='Byte',
output_Raster = FALSE)


streams_saga <- file.path(tempdir(), "streams_saga_r.sgrd")
streams_resamp <- file.path(tempdir(), "streams_resamp.sgrd")
RSAGA::rsaga.import.gdal(in.grid = streams_r, out.grid = streams_saga, env = env)
Expand Down Expand Up @@ -99,6 +120,24 @@ create_dem_products <- function(dem,stream_vec = NULL,burn_val=NULL,outdir, prod
products.out <- file.path(outdir, products)
raster::extension(products.out) <- "sgrd"


get_param_lst_idx <- function(product)
{
idx<-1
if(!is.null(param_vect[1])){
for(i in c(1:length(param_vect)))
{
if(product == param_vect[i])
{
idx<-i+1
}
i<-i+1
}
}
return(idx)
}


for (i in 1:length(products)) {
p <- products[i]

Expand All @@ -112,7 +151,7 @@ create_dem_products <- function(dem,stream_vec = NULL,burn_val=NULL,outdir, prod
out.aspect = products.out[i],
unit.aspect = "radians",
env = env)

RSAGA::rsaga.grid.calculus(file.path(outdir,'ASPECT.sgrd'),file.path(outdir,"NORTHNESS.sgrd"),~cos(a))


Expand All @@ -132,8 +171,8 @@ create_dem_products <- function(dem,stream_vec = NULL,burn_val=NULL,outdir, prod
USE_NODATA=1,
GRIDS=file.path(outdir,'EASTNESS.sgrd'),
RESULT=file.path(outdir,'EASTNESS.sgrd')))


} else if (p == "CPLAN") {
RSAGA::rsaga.slope.asp.curv(in.dem = dem.sgrd,
out.cplan = products.out[i],
Expand All @@ -145,20 +184,20 @@ create_dem_products <- function(dem,stream_vec = NULL,burn_val=NULL,outdir, prod
} else if (p == "DAH") {
RSAGA::rsaga.geoprocessor(lib = "ta_morphometry",
module = 12,
param = list(DEM = dem.sgrd,
DAH = products.out[i]),
param = append(list(DEM = dem.sgrd,
DAH = products.out[i]),param_values[[get_param_lst_idx("DAH")]]),
env = env)
} else if (p == "MRVBF") {
RSAGA::rsaga.geoprocessor(lib = "ta_morphometry",
module = 8,
param = list(DEM = dem.sgrd,
MRVBF = products.out[i]),
param = append(list(DEM = dem.sgrd,
MRVBF = products.out[i]),param_values[[get_param_lst_idx("MRVBF")]]),
env = env)
} else if (p == "TPI") {
RSAGA::rsaga.geoprocessor(lib = "ta_morphometry",
module = 18,
param = list(DEM = dem.sgrd,
TPI = products.out[i]),
param = append(list(DEM = dem.sgrd,
TPI = products.out[i]),param_values[[get_param_lst_idx("TPI")]]),
env = env)
} else if (p == "TOPOWET") {
RSAGA::rsaga.wetness.index(in.dem = dem.sgrd,
Expand All @@ -167,8 +206,8 @@ create_dem_products <- function(dem,stream_vec = NULL,burn_val=NULL,outdir, prod
} else if (p == "CAREA") {
RSAGA::rsaga.geoprocessor(lib = "ta_hydrology",
module = 0,
param = list(ELEVATION = dem.sgrd,
FLOW = products.out[i]),
param = append(list(ELEVATION = dem.sgrd,
FLOW = products.out[i]),param_values[[get_param_lst_idx("CAREA")]]),
env = env)
}
}
Expand Down
4 changes: 3 additions & 1 deletion README.md
Original file line number Diff line number Diff line change
Expand Up @@ -37,7 +37,7 @@ This function assumes `set_grass_env` has been called. This function attributes
## Installation
Get the latest version from GitHub with:
```r
devtools::install_github("bcgov/wetlandmapR")
remotes::install_github("bcgov/wetlandmapR",ref="master")
```

## Or with Docker
Expand All @@ -51,5 +51,7 @@ docker run -e PASSWORD=URPassword -p 8787:8787 --rm huntgdok/wetlandmapr:latest
```
Where `URPassword` is any password of your choice, and username `rstudio`. The running container can be viewed by passing [localhost](http://localhost:8787/) to your browser. Be sure to copy all outputs locally before exiting as all data will be lost.

PLEASE NOTE: A Docker subscription is required if using the *wetlandmapR* Docker image for commercial use, see service agreement [Docker Subscription Service Agreement](https://www.docker.com/legal/docker-subscription-service-agreement?utm_campaign=2021-08-31-business-tier-launch&utm_medium=email&utm_source=mailgun&utm_content=service-agreement) for more information.

## Examples
See the example vignette which describes [wetlandmapR_example](https://htmlpreview.github.io/?https://github.com/bcgov/wetlandmapR/blob/significant_update_HG/vignettes/wetlandmapR_example.html) how to use the functions in this package together for mapping wetlands.
13 changes: 12 additions & 1 deletion man/create_dem_products.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.