diff --git a/R/get_component_from_GDB.R b/R/get_component_from_GDB.R index 5fd10979..be1478b7 100644 --- a/R/get_component_from_GDB.R +++ b/R/get_component_from_GDB.R @@ -3,11 +3,11 @@ get_component_from_GDB <- function(dsn = "gNATSGO_CONUS.gdb", WHERE = NULL, childs = FALSE, droplevels = TRUE, stringsAsFactors = NULL) { # check - co_vars <- "comppct_l|comppct_r|comppct_h|compname|compkind|majcompflag|otherph|localphase|slope_l|slope_r|slope_h|slopelenusle_l|slopelenusle_r|slopelenusle_h|runoff|tfact|wei|weg|erocl|earthcovkind1|earthcovkind2|hydricon|hydricrating|drainagecl|elev_l|elev_r|elev_h|aspectccwise|aspectrep|aspectcwise|geomdesc|albedodry_l|albedodry_r|albedodry_h|airtempa_l|airtempa_r|airtempa_h|map_l|map_r|map_h|reannualprecip_l|reannualprecip_r|reannualprecip_h|ffd_l|ffd_r|ffd_h|nirrcapcl|nirrcapscl|nirrcapunit|irrcapcl|irrcapscl|irrcapunit|cropprodindex|constreeshrubgrp|wndbrksuitgrp|rsprod_l|rsprod_r|rsprod_h|foragesuitgrpid|wlgrain|wlgrass|wlherbaceous|wlshrub|wlconiferous|wlhardwood|wlwetplant|wlshallowwat|wlrangeland|wlopenland|wlwoodland|wlwetland|soilslippot|frostact|initsub_l|initsub_r|initsub_h|totalsub_l|totalsub_r|totalsub_h|hydgrp|corcon|corsteel|taxclname|taxorder|taxsuborder|taxgrtgroup|taxsubgrp|taxpartsize|taxpartsizemod|taxceactcl|taxreaction|taxtempcl|taxmoistscl|taxtempregime|soiltaxedition|castorieindex|flecolcomnum|flhe|flphe|flsoilleachpot|flsoirunoffpot|fltemik2use|fltriumph2use|indraingrp|innitrateleachi|misoimgmtgrp|vasoimgtgrp|cokey|mukey" - co_idx <- grepl(co_vars, WHERE, ignore.case = TRUE) - if (!length(co_idx) > 0) co_idx <- FALSE + vars_co <- "comppct_l|comppct_r|comppct_h|compname|compkind|majcompflag|otherph|localphase|slope_l|slope_r|slope_h|slopelenusle_l|slopelenusle_r|slopelenusle_h|runoff|tfact|wei|weg|erocl|earthcovkind1|earthcovkind2|hydricon|hydricrating|drainagecl|elev_l|elev_r|elev_h|aspectccwise|aspectrep|aspectcwise|geomdesc|albedodry_l|albedodry_r|albedodry_h|airtempa_l|airtempa_r|airtempa_h|mast_l|mast_r|mast_h|map_l|map_r|map_h|reannualprecip_l|reannualprecip_r|reannualprecip_h|ffd_l|ffd_r|ffd_h|nirrcapcl|nirrcapscl|nirrcapunit|irrcapcl|irrcapscl|irrcapunit|cropprodindex|constreeshrubgrp|wndbrksuitgrp|rsprod_l|rsprod_r|rsprod_h|foragesuitgrpid|wlgrain|wlgrass|wlherbaceous|wlshrub|wlconiferous|wlhardwood|wlwetplant|wlshallowwat|wlrangeland|wlopenland|wlwoodland|wlwetland|soilslippot|frostact|initsub_l|initsub_r|initsub_h|totalsub_l|totalsub_r|totalsub_h|hydgrp|corcon|corsteel|taxclname|taxorder|taxsuborder|taxgrtgroup|taxsubgrp|taxpartsize|taxpartsizemod|taxceactcl|taxreaction|taxtempcl|taxmoistscl|taxtempregime|soiltaxedition|castorieindex|flecolcomnum|flhe|flphe|flsoilleachpot|flsoirunoffpot|fltemik2use|fltriumph2use|indraingrp|innitrateleachi|misoimgmtgrp|vasoimgtgrp|cokey|mukey" + idx_co <- grepl(vars_co, WHERE, ignore.case = TRUE) + if (!length(idx_co) > 0) idx_co <- FALSE - if (!co_idx & !is.null(WHERE)) { + if (!idx_co & !is.null(WHERE)) { stop("the WHERE argument is not targeting the component table") } @@ -28,7 +28,8 @@ get_component_from_GDB <- function(dsn = "gNATSGO_CONUS.gdb", WHERE = NULL, chil # recode metadata domains co <- uncode( co, - droplevels = droplevels + droplevels = droplevels, + stringsAsFactors = stringsAsFactors ) @@ -39,10 +40,10 @@ get_component_from_GDB <- function(dsn = "gNATSGO_CONUS.gdb", WHERE = NULL, chil # exec sub query ## pm - idx <- seq(0, 2e8, 1000) - co$idx <- as.character(cut(1:nrow(co), breaks = idx)) - if (!is.null(WHERE)) { + idx <- seq(0, 2e8, 1000) + co$idx <- as.character(cut(1:nrow(co), breaks = idx)) + copm <- by(co, co$idx, function(x) { temp <- .get_copmgrp_from_GDB(dsn = dsn, x) }) @@ -53,27 +54,30 @@ get_component_from_GDB <- function(dsn = "gNATSGO_CONUS.gdb", WHERE = NULL, chil ## gmd - idx <- seq(0, 2e8, 1000) - co$idx <- as.character(cut(1:nrow(co), breaks = idx)) - if (!is.null(WHERE)) { - cogmd <- by(co, co$idx, function(x) { + idx <- seq(0, 2e8, 1000) + co$idx <- as.character(cut(1:nrow(co), breaks = idx)) + + cogmd <- by(co, co$idx, function(x) { temp <- .get_cogeomordesc_from_GDB(dsn = dsn, x) }) cogmd <- do.call("rbind", cogmd) } else { - cogmd <- .get_cogeomordesc_from_GDB(dsn = dsn) + cogmd <- .get_cogeomordesc_from_GDB(dsn) } # prep - copm <- .copm_prep(copm, db = "SDA") - cogmd <- .cogmd_prep(cogmd, db = "SDA") + # copm <- soilDB:::.copm_prep(copm, db = "SDA") + # cogmd <- soilDB:::.cogmd_prep(cogmd, db = "SDA") + copm2 <- .copm_prep2(copm, "cokey") + cogmd2 <- .cogmd_prep2(cogmd, "cokey") + # merge co$idx <- NULL - co <- merge(co, copm, by = "cokey", all.x = TRUE, sort = FALSE) - co <- merge(co, cogmd, by = "cokey", all.x = TRUE, sort = FALSE) + co <- merge(co, copm2, by = "cokey", all.x = TRUE, sort = FALSE) + co <- merge(co, cogmd2, by = "cokey", all.x = TRUE, sort = FALSE) } @@ -88,10 +92,10 @@ get_legend_from_GDB <- function(dsn = "gNATSGO_CONUS.gdb", WHERE = NULL, droplev if (!is.null(WHERE)) { # check - le_vars <- "mlraoffice|areasymbol|areaname|areatypename|areaacres|ssastatus|projectscale|cordate|lkey" - le_idx <- grepl(le_vars, WHERE, ignore.case = TRUE) + vars_le <- "mlraoffice|areasymbol|areaname|areatypename|areaacres|ssastatus|projectscale|cordate|lkey" + idx_le <- grepl(vars_le, WHERE, ignore.case = TRUE) - if (! le_idx) { + if (! idx_le) { stop("the WHERE argument is not targeting the legend table") } @@ -133,17 +137,17 @@ get_mapunit_from_GDB <- function(dsn = "gNATSGO_CONUS.gdb", WHERE = NULL, drople # tests if (!is.null(WHERE)) { - le_vars <- "mlraoffice|areasymbol|areaname|areatypename|areaacres|ssastatus|projectscale|cordate" - le_idx <- grepl(le_vars, WHERE, ignore.case = TRUE) + vars_le <- "mlraoffice|areasymbol|areaname|areatypename|areaacres|ssastatus|projectscale|cordate" + idx_le <- grepl(vars_le, WHERE, ignore.case = TRUE) - mu_vars <- "lkey|mukey|musym|muname|mukind|mustatus|invesintens|muacres|farmlndcl" - mu_idx <- grepl(mu_vars, WHERE, ignore.case = TRUE) + vars_mu <- "lkey|mukey|musym|muname|mukind|mustatus|invesintens|muacres|farmlndcl" + idx_mu <- grepl(vars_mu, WHERE, ignore.case = TRUE) - if (le_idx & mu_idx) { + if (idx_le & idx_mu) { stop("the WHERE argument can not target both the legend and mapunit table at the same time") } - if (! le_idx & ! mu_idx) { + if (! idx_le & ! idx_mu) { stop("the WHERE argument is not targeting either the legend or mapunit table") } @@ -151,7 +155,7 @@ get_mapunit_from_GDB <- function(dsn = "gNATSGO_CONUS.gdb", WHERE = NULL, drople message("getting mapunits from WHERE ", WHERE) - if (mu_idx) { + if (idx_mu) { qry <- paste0("SELECT * FROM mapunit WHERE ", WHERE) mu <- sf::read_sf(dsn = dsn, query = qry, as_tibble = FALSE) @@ -159,7 +163,7 @@ get_mapunit_from_GDB <- function(dsn = "gNATSGO_CONUS.gdb", WHERE = NULL, drople le <- sf::read_sf(dsn = dsn, query = qry, as_tibble = FALSE) } - if (le_idx) { + if (idx_le) { qry <- paste0("SELECT * FROM legend WHERE ", WHERE) le <- sf::read_sf(dsn = dsn, query = qry, as_tibble = FALSE) @@ -185,50 +189,72 @@ get_mapunit_from_GDB <- function(dsn = "gNATSGO_CONUS.gdb", WHERE = NULL, drople if (stats == TRUE) { - idx <- c(0, rep(3000, 1000) * 1:1000) - mu$idx <- as.character(cut(1:nrow(mu), breaks = idx)) - - co <- by(mu, mu$idx, function(x) { - - qry <- paste0( - "SELECT mukey, cokey, comppct_r, majcompflag, hydricrating - - FROM component - - WHERE mukey IN ('", paste0(x$mukey, collapse = "', '"), "')" - ) - message("getting components from subset ", x$idx[1]) - sf::read_sf(dsn = dsn, query = qry, as_tibble = FALSE) - }) - co <- do.call("rbind", co) - mu$idx <- NULL + if (!is.null(WHERE)) { + idx <- c(0, rep(3000, 1000) * 1:1000) + mu$idx <- as.character(cut(1:nrow(mu), breaks = idx)) + + co <- by(mu, mu$idx, function(x) { + qry <- paste0( + "SELECT mukey, cokey, comppct_r, majcompflag, hydricrating + + FROM component + + WHERE mukey IN ('", paste0(x$mukey, collapse = "', '"), "')" + ) + message("getting components from subset ", x$idx[1]) + co <- sf::read_sf(dsn = dsn, query = qry, as_tibble = FALSE) + }) + co <- do.call("rbind", co) + mu$idx <- NULL + } else { + qry <- "SELECT mukey, cokey, comppct_r, majcompflag, hydricrating + FROM component" + co <- sf::read_sf(dsn, query = qry, as_tibble = FALSE) + } if (nrow(co) > 0) { - co2 <- { - temp <- data.frame( - mukey = as.character(0), - pct_component = as.integer(0), - pct_hydric = as.integer(0), - n_component = as.integer(0), - n_majcompflag = as.integer(0), - stringsAsFactors = FALSE - ) + co2 <- data.table::as.data.table(co) + . = NULL + comppct_r = NULL + hydricrating = NULL + cokey = NULL + majcompflag = NULL + mukey = NULL + co2 <- co2[, .( + pct_component = sum(comppct_r, na.rm = TRUE), + pct_hydric = sum((hydricrating == "Yes") * comppct_r, na.rm = TRUE), + n_component = length(cokey), + n_majcompflag = sum(majcompflag == "Yes", na.rm = TRUE) + ), + by = mukey + ] |> + as.data.frame() - co$df <- list(temp)[rep(1, times = nrow(co))] - split(co, co$mukey, drop = TRUE) ->.; - lapply(., function(x) { - df = x$df[[1]] - df$mukey = x$mukey[1] - df$pct_component = sum(x$comppct_r, na.rm = TRUE) - df$pct_hydric = sum((x$hydricrating == "Yes") * x$comppct_r, na.rm = TRUE) - df$n_component = length(x$cokey) - df$n_majcompflag = sum(x$majcompflag == "Yes", na.rm = TRUE) - return(df) - }) ->.; - do.call("rbind", .) ->.; - } + # co2 <- { + # temp <- data.frame( + # mukey = as.character(0), + # pct_component = as.integer(0), + # pct_hydric = as.integer(0), + # n_component = as.integer(0), + # n_majcompflag = as.integer(0), + # stringsAsFactors = FALSE + # ) + # + # co$df <- list(temp)[rep(1, times = nrow(co))] + # split(co, co$mukey, drop = TRUE) ->.; + # lapply(., function(x) { + # df = x$df[[1]] + # df$mukey = x$mukey[1] + # df$pct_component = sum(x$comppct_r, na.rm = TRUE) + # df$pct_hydric = sum((x$hydricrating == "Yes") * x$comppct_r, na.rm = TRUE) + # df$n_component = length(x$cokey) + # df$n_majcompflag = sum(x$majcompflag == "Yes", na.rm = TRUE) + # return(df) + # }) ->.; + # do.call("rbind", .) ->.; + # } mu <- merge(mu, co2, by = "mukey", all.x = TRUE, sort = FALSE) } else { mu = cbind(mu, pct_component = NA_integer_, pct_hydric = NA_integer_, n_component = NA_integer_, n_majcompflag = NA_integer_) @@ -258,92 +284,134 @@ get_mapunit_from_GDB <- function(dsn = "gNATSGO_CONUS.gdb", WHERE = NULL, drople .get_copmgrp_from_GDB <- function(dsn = dsn, co = NULL) { - # message("getting ", unique(substr(x$areasymbol, 1, 2))) + # query copmgrp table qry <- paste0( - "SELECT cokey, copmgrpkey, pmgroupname + "SELECT * - FROM copmgrp" + FROM copmgrp + + WHERE rvindicator = 'Yes' + " , if (!is.null(co)) { - paste0("WHERE rvindicator = 'Yes' AND cokey IN ", soilDB::format_SQL_in_statement(co$cokey)) + paste0("AND cokey IN ", soilDB::format_SQL_in_statement(co$cokey)) } ) pmg <- sf::read_sf(dsn = dsn, query = qry, as_tibble = FALSE) - - qry <- paste0( - "SELECT copmgrpkey, pmorder, pmkind, pmorigin + + + # remove duplicate rvindicators + dat <- table(pmg$cokey, pmg$rvindicator) |> + as.data.frame.matrix() + n_rvindicator = NULL + dat <- data.frame(cokey = row.names(dat), n_rvindicator = dat$Yes) |> + subset(n_rvindicator > 1) + assign('dup.compmgrp.cokeyrvindictor', value=dat, envir=get_soilDB_env()) + message("-> QC: ", formatC(nrow(dat), format = "fg", big.mark = ","), " duplicate 'representative' rvindicators in the copmgrp table. \n\tUse `get('dup.compmgrp.cokeyrvindictor', envir=get_soilDB_env())` for offending component record IDs (cokey)") + + pmg$rvindicator <- NULL + pmg2 <- .flatten_gmd(pmg[c("cokey", "pmgroupname")], "cokey", "copmgrp", sep = "; ") + pmg$pmgroupname <- NULL + pmg <- merge(pmg, pmg2, by = "cokey", all.x = TRUE, sort = FALSE) + + + # query copm table + qry <- paste( + "SELECT * FROM copm" , if (!is.null(co)) { - paste0("WHERE copmgrpkey IN ", soilDB::format_SQL_in_statement(co$copmgrpkey)) + paste0("WHERE copmgrpkey IN ", soilDB::format_SQL_in_statement(pmg$copmgrpkey)) } ) pm <- sf::read_sf(dsn = dsn, query = qry, as_tibble = FALSE) pm <- merge(pmg, pm, by = "copmgrpkey", all.x = TRUE, sort = FALSE) - + + pm <- data.table::as.data.table(pm) pm <- pm[with(pm, order(cokey, copmgrpkey, pmorder)), ] + pm <- as.data.frame(pm) return(pm) } -.get_cogeomordesc_from_GDB <- function(dsn = dsn, co) { - +.get_cogeomordesc_from_GDB <- function(dsn = dsn, co = NULL) { + + + # geomorphic description ---- qry <- paste0( "SELECT cokey, geomftname, geomfname, geomfeatid, existsonfeat, cogeomdkey FROM cogeomordesc - WHERE rvindicator = 'Yes' AND - cokey IN ", format_SQL_in_statement(co$cokey) - ) + WHERE rvindicator = 'Yes'", + if (!is.null(co)) { + paste0("AND cokey IN ", format_SQL_in_statement(co$cokey)) + }) cogmd <- sf::read_sf(dsn, query = qry, as_tibble = FALSE) cogmd_ls <- cogmd[cogmd$geomftname == "Landscape", ] cogmd_lf <- cogmd[cogmd$geomftname == "Landform", ] cogmd_lf$geomftname <- NULL - + names(cogmd_ls)[names(cogmd_ls) == "geomfname"] <- "landscape" names(cogmd_lf)[names(cogmd_lf) == "geomfname"] <- "landform" - - # geomorphic components + + cogmd_ls <- .flatten_gmd(cogmd_ls, "cokey", "cogeomordesc") + + + + # geomorphic components ---- qry <- paste0( - "SELECT geomposmntn, geomposhill, geompostrce, geomposflats, cogeomdkey - - FROM cosurfmorphgc - - WHERE cogeomdkey IN ", format_SQL_in_statement(cogmd_lf$cogeomdkey) - ) + "SELECT cogeomdkey, geomposmntn, geomposhill, geompostrce, geomposflats + + FROM cosurfmorphgc ", + + if (!is.null(co)) { + paste0("WHERE cogeomdkey IN ", format_SQL_in_statement(cogmd_lf$cogeomdkey)) + }) lf_3d <- sf::read_sf(dsn, query = qry, as_tibble = FALSE) - names(lf_3d)[1:4] <- gsub("geompos", "", names(lf_3d)[1:4]) - - - # slope shape + names(lf_3d)[1:5] <- gsub("geompos", "", names(lf_3d)[1:5]) + + lf_3d <- .flatten_gmd(lf_3d, key = "cogeomdkey", table = "cosurfmorphgc") + + + + # slope shape ---- qry <- paste0( - "SELECT shapeacross, shapedown, cogeomdkey + "SELECT cogeomdkey, shapeacross, shapedown - FROM cosurfmorphss + FROM cosurfmorphss ", - WHERE cogeomdkey IN ", format_SQL_in_statement(cogmd_lf$cogeomdkey) - ) + if (!is.null(co)) { + paste0("WHERE cogeomdkey IN ", format_SQL_in_statement(cogmd_lf$cogeomdkey)) + }) lf_ss <- sf::read_sf(dsn, query = qry, as_tibble = FALSE) + + lf_ss <- .format_slopeshape(lf_ss) + lf_ss <- .flatten_gmd(lf_ss, "cogeomdkey", "cosurfmorphss") + + - - # hillslope position + # hillslope position ---- qry <- paste0( - "SELECT hillslopeprof, cogeomdkey - - FROM cosurfmorphhpp + "SELECT cogeomdkey, hillslopeprof - WHERE cogeomdkey IN ", format_SQL_in_statement(cogmd_lf$cogeomdkey) - ) + FROM cosurfmorphhpp ", + + if (!is.null(co)) { + paste0("WHERE cogeomdkey IN ", format_SQL_in_statement(cogmd_lf$cogeomdkey)) + }) lf_2d <- sf::read_sf(dsn, query = qry, as_tibble = FALSE) - - + + lf_2d <- .flatten_gmd(lf_2d, "cogeomdkey", "cosurfmorphhpp") + + + # merge results - cogmd <- merge(cogmd_ls[c("cokey", "landscape")], cogmd_lf, by = "cokey", all.x = TRUE, sort = FALSE) + cogmd <- merge(cogmd_ls[c("cokey", "landscape")], cogmd_lf, by = "cokey", all = TRUE, sort = FALSE) cogmd <- merge(cogmd, lf_3d, by = "cogeomdkey", all.x = TRUE, sort = FALSE) cogmd <- merge(cogmd, lf_ss, by = "cogeomdkey", all.x = TRUE, sort = FALSE) cogmd <- merge(cogmd, lf_2d, by = "cogeomdkey", all.x = TRUE, sort = FALSE) @@ -366,18 +434,16 @@ get_mapunit_from_GDB <- function(dsn = "gNATSGO_CONUS.gdb", WHERE = NULL, drople if (!is.null(cokey)) { idx2 <- as.character(cut(1:length(cokey), breaks = idx)) co <- data.frame(cokey, idx2) - - ch <- by(co, co$idx, function(x) { - qry <- paste0( - "SELECT * - - FROM chorizon - - WHERE cokey IN ", format_SQL_in_statement(x$cokey) - ) - ch <- sf::read_sf(dsn = dsn, query = qry, as_tibble = FALSE) - }) - ch <- do.call("rbind", ch) + + ch <- by(co, co$idx, function(x) { + qry <- paste0( + "SELECT * + FROM chorizon + WHERE cokey IN ", format_SQL_in_statement(x$cokey) + ) + ch <- sf::read_sf(dsn = dsn, query = qry, as_tibble = FALSE) + }) + ch <- do.call("rbind", ch) } else ch <- sf::read_sf(dsn = dsn, query = "SELECT * FROM chorizon", as_tibble = FALSE) @@ -563,27 +629,27 @@ fetchGDB <- function(dsn = "gNATSGO_CONUS.gdb", ) { # checks - le_vars <- "mlraoffice|areasymbol|areaname|areatypename|areaacres|ssastatus|projectscale|cordate|lkey" - mu_vars <- "mukey|musym|muname|mukind|mustatus|invesintens|muacres|farmlndcl" - co_vars <- "comppct_l|comppct_r|comppct_h|compname|compkind|majcompflag|otherph|localphase|slope_l|slope_r|slope_h|slopelenusle_l|slopelenusle_r|slopelenusle_h|runoff|tfact|wei|weg|erocl|earthcovkind1|earthcovkind2|hydricon|hydricrating|drainagecl|elev_l|elev_r|elev_h|aspectccwise|aspectrep|aspectcwise|geomdesc|albedodry_l|albedodry_r|albedodry_h|airtempa_l|airtempa_r|airtempa_h|map_l|map_r|map_h|reannualprecip_l|reannualprecip_r|reannualprecip_h|ffd_l|ffd_r|ffd_h|nirrcapcl|nirrcapscl|nirrcapunit|irrcapcl|irrcapscl|irrcapunit|cropprodindex|constreeshrubgrp|wndbrksuitgrp|rsprod_l|rsprod_r|rsprod_h|foragesuitgrpid|wlgrain|wlgrass|wlherbaceous|wlshrub|wlconiferous|wlhardwood|wlwetplant|wlshallowwat|wlrangeland|wlopenland|wlwoodland|wlwetland|soilslippot|frostact|initsub_l|initsub_r|initsub_h|totalsub_l|totalsub_r|totalsub_h|hydgrp|corcon|corsteel|taxclname|taxorder|taxsuborder|taxgrtgroup|taxsubgrp|taxpartsize|taxpartsizemod|taxceactcl|taxreaction|taxtempcl|taxmoistscl|taxtempregime|soiltaxedition|castorieindex|flecolcomnum|flhe|flphe|flsoilleachpot|flsoirunoffpot|fltemik2use|fltriumph2use|indraingrp|innitrateleachi|misoimgmtgrp|vasoimgtgrp|cokey" + vars_le <- "mlraoffice|areasymbol|areaname|areatypename|areaacres|ssastatus|projectscale|cordate|lkey" + vars_mu <- "mukey|musym|muname|mukind|mustatus|invesintens|muacres|farmlndcl" + vars_co <- "comppct_l|comppct_r|comppct_h|compname|compkind|majcompflag|otherph|localphase|slope_l|slope_r|slope_h|slopelenusle_l|slopelenusle_r|slopelenusle_h|runoff|tfact|wei|weg|erocl|earthcovkind1|earthcovkind2|hydricon|hydricrating|drainagecl|elev_l|elev_r|elev_h|aspectccwise|aspectrep|aspectcwise|geomdesc|albedodry_l|albedodry_r|albedodry_h|airtempa_l|airtempa_r|airtempa_h|map_l|map_r|map_h|reannualprecip_l|reannualprecip_r|reannualprecip_h|ffd_l|ffd_r|ffd_h|nirrcapcl|nirrcapscl|nirrcapunit|irrcapcl|irrcapscl|irrcapunit|cropprodindex|constreeshrubgrp|wndbrksuitgrp|rsprod_l|rsprod_r|rsprod_h|foragesuitgrpid|wlgrain|wlgrass|wlherbaceous|wlshrub|wlconiferous|wlhardwood|wlwetplant|wlshallowwat|wlrangeland|wlopenland|wlwoodland|wlwetland|soilslippot|frostact|initsub_l|initsub_r|initsub_h|totalsub_l|totalsub_r|totalsub_h|hydgrp|corcon|corsteel|taxclname|taxorder|taxsuborder|taxgrtgroup|taxsubgrp|taxpartsize|taxpartsizemod|taxceactcl|taxreaction|taxtempcl|taxmoistscl|taxtempregime|soiltaxedition|castorieindex|flecolcomnum|flhe|flphe|flsoilleachpot|flsoirunoffpot|fltemik2use|fltriumph2use|indraingrp|innitrateleachi|misoimgmtgrp|vasoimgtgrp|cokey" if (!is.null(WHERE)) { - le_idx <- grepl(le_vars, WHERE, ignore.case = TRUE) - mu_idx <- grepl(mu_vars, WHERE, ignore.case = TRUE) - co_idx <- grepl(co_vars, WHERE, ignore.case = TRUE) + idx_le <- grepl(vars_le, WHERE, ignore.case = TRUE) + idx_mu <- grepl(vars_mu, WHERE, ignore.case = TRUE) + idx_co <- grepl(vars_co, WHERE, ignore.case = TRUE) - if (sum(c(le_idx, mu_idx, co_idx)) > 1) { + if (sum(c(idx_le, idx_mu, idx_co)) > 1) { stop("WHERE can only target 1 table at a time") } } else { - le_idx <- FALSE - mu_idx <- FALSE - co_idx <- FALSE + idx_le <- FALSE + idx_mu <- FALSE + idx_co <- FALSE } # target legend table - if (le_idx) { + if (idx_le) { le <- get_legend_from_GDB(dsn = dsn, WHERE = WHERE, stats = FALSE) qry <- paste0("lkey IN ('", paste(le$lkey, collapse = "', '"), "')") @@ -592,21 +658,21 @@ fetchGDB <- function(dsn = "gNATSGO_CONUS.gdb", } # target mapunit table - if (mu_idx) { + if (idx_mu) { mu <- suppressMessages(get_mapunit_from_GDB(dsn = dsn, WHERE = WHERE)) mu <- mu[order(mu$areasymbol), ] } # get components - if (le_idx | mu_idx) { + if (idx_le | idx_mu) { - if (le_idx) mu$idx <- mu$areasymbol - if (mu_idx) mu$idx <- "1" + if (idx_le) mu$idx <- mu$areasymbol + if (idx_mu) mu$idx <- "1" temp <- by(mu, mu$idx, function(x) { - if (le_idx) { + if (idx_le) { message("getting components and horizons from areasymbol = '", unique(x$idx), "'") } else message("getting components and horizons from ", WHERE) @@ -636,7 +702,7 @@ fetchGDB <- function(dsn = "gNATSGO_CONUS.gdb", # horizons tryCatch({ - h <- .get_chorizon_from_GDB(dsn = dsn, co) + h <- .get_chorizon_from_GDB(dsn = dsn, co$cokey) }, error = function(err) { print(paste("Error occured: ", err)) @@ -650,7 +716,9 @@ fetchGDB <- function(dsn = "gNATSGO_CONUS.gdb", co <- do.call("rbind", lapply(temp, function(x) x$co)) h <- do.call("rbind", lapply(temp, function(x) x$h)) - } else { + } + + if (idx_co | is.null(WHERE)) { message("getting components and horizons from ", WHERE) diff --git a/R/utils.R b/R/utils.R index 77b41638..50530581 100644 --- a/R/utils.R +++ b/R/utils.R @@ -989,3 +989,383 @@ return(df) } + + + + +.copm_prep2 <- function(x, key = NULL) { + + idx_key <- grep(key, names(x)) + idx_pmkey <- grep("pmgrpkey", names(x)) + nm_pmkey <- names(x)[idx_pmkey] + names(x)[c(idx_key, idx_pmkey)] <- c("key", "pmgrpkey") + + vars <- c("key", "pmgrpkey", "pmorder", "pmkind") + + # pmk1 <- pmk + # pmk1 <- pmk1[with(pmk1, order(key, pmgrpkey, pmorder)), ] + # pmk1 <- aggregate(pmkind ~ key, data = pmk1, FUN = function(x) paste0(x, collapse = " over ")) + # + # test <- { + # strsplit(pmk1$pmkind, " over ") ->.; + # lapply(., function(x) { + # rle(x)$lengths |> + # cumsum() |> + # x[i = _] ->.; + # paste(., collapse = " over ") + # }) ->.; + # unlist(.) + # } + # pmk1$pmkind <- test + + + pm <- data.table::as.data.table(x); rm(x) + pm <- pm[order(pm$key, pm$pmgrpkey, pm$pmorder)] + + pm$id_k <- paste(pm$key, pm$pmkind) + pm$id_o <- paste(pm$key, pm$pmorigin) + + + # pmkind + # remove duplicate pmkind by cokey + pm_k <- { + vars <- c("key", "pmkind", "id_k") + # ..vars = NULL + pm_k <- pm[!is.na(pm$pmkind), vars, with = FALSE] + idx <- rle(pm_k$id_k)$lengths |> + cumsum() + pm_k[idx, ] + } + . = NULL + pmkind = NULL + pm_k <- pm_k[, .(pmkind = paste0(pmkind, collapse = " over ")), by = .(key)] + + + # pmorigin + pmorigin = NULL + pm_o <- { + vars <- c("key", "pmorigin", "id_o") + # ..vars = NULL + pm_o <- pm[!is.na(pm$pmorigin), vars, with = FALSE] + idx <- rle(pm_o$id_o)$lengths |> + cumsum() + pm_o[idx, ] + } + pm_o <- pm_o[, .(pmorigin = paste0(pmorigin, collapse = " over ")), by = .(key)] + + + # merge + pm <- merge(pm_k, pm_o, by = "key", all = TRUE, sort = FALSE) |> + as.data.frame() + + + names(pm)[1] <- c(key) + + return(pm) +} + + + +.cogmd_prep2 <- function(data, key = "cokey") { + + idx_key <- grep(key, names(data)) + names(data)[c(idx_key)] <- c("key") + + + # find sites with overlapping landforms ---- + n_bot = NULL + n_mis_geomfeatid = NULL + geomfeatid = NULL + .N = NULL + N = NULL + + test <- data.table::as.data.table(data)[ + , .( + # n_bot = sum(! existsonfeat %in% geomfeatid, na.rm = TRUE), + n_bot = sum(! match(existsonfeat, geomfeatid, nomatch = 0, incomparables = NA_integer_) > 0, na.rm = TRUE), + n_geomfeatid = sum(!is.na(geomfeatid)), + n_existonfeat = sum(!is.na(existsonfeat)), + .N, + n_mis_geomfeatid = sum(is.na(geomfeatid)) + ), + by = key + ] + data <- merge(test, data, by = "key", all.y = TRUE) + + + # determine row direction ---- + # ordered + data <- within(data, { + existsonfeat = ifelse(geomfeatid == existsonfeat, NA, existsonfeat) + existsonfeat = ifelse(n_bot == N, NA, existsonfeat) + + row_dir = ifelse(geomfeatid < existsonfeat, "top2bot", "bot2top") + row_dir = ifelse( + geomfeatid == existsonfeat + 1 | geomfeatid == existsonfeat - 1, + row_dir, + "chaos" + ) + row_dir = ifelse(n_bot == N | is.na(existsonfeat), "missing", row_dir) + row_dir = factor(row_dir, levels = c("top2bot", "bot2top", "chaos", "missing")) + }) + + + # find chaos within a component ---- + tb <- with(data, table(key, row_dir)) |> + as.data.frame.matrix() + chaos <- tb |> + within({ + tot = rowSums(cbind(top2bot > 0, bot2top > 0, chaos > 0)) + co_dir = ifelse(tot > 1, "chaos", "ordered") + }) |> + cbind(key = row.names(tb)) + data <- merge(data, chaos, by = "key", all.x = TRUE, sort = FALSE) + data <- within(data, { + co_dir = ifelse(N == n_bot | N == n_mis_geomfeatid | N == missing, "missing", co_dir) + }) + + + # replace NA row direction, where the component direction == "ordered" ---- + vars <- c("top2bot", "bot2top", "chaos") + # ..vars = NULL + idx <- which(data$row_dir == "missing" & data$co_dir == "ordered") + ordered_mis <- names(data[, vars, with = FALSE])[max.col(data[idx, vars, with = FALSE])] + data[idx, "row_dir"] <- ordered_mis + + # subset(data, key == "22230267") + + + # subset and sort different ordering conventions ---- + # top2bot & NA + top2bot <- { + subset(data, row_dir == "top2bot" & co_dir == "ordered") ->.; + .[order(.$key, .$geomfeatid, .$existsonfeat), ] + } + # bot2top + bot2top <- { + subset(data, row_dir == "bot2top" & co_dir == "ordered") ->.; + .[order(.$key, - .$geomfeatid, - .$existsonfeat), ] + } + # chaos and missing ordered + chaos2 <- { + subset(data, co_dir %in% c("chaos", "missing") | row_dir == "chaos") ->.; + .[order(.$key, .$geomfeatid, .$existsonfeat), ] + } + + + # recombine + data <- rbind(top2bot, bot2top, chaos2) + rm(top2bot, bot2top, chaos2) + + + # find N tops ---- + # test <- data.table::as.data.table(data)[ + # , .( + # # n_bot = sum(! existsonfeat %in% geomfeatid, na.rm = TRUE), + # n_bot = sum(! match(existsonfeat, geomfeatid, nomatch = 0, incomparables = NA_integer_) > 0, na.rm = TRUE), + # + # .N, + # n_mis_geomfeatid = sum(is.na(geomfeatid)) + # ), + # by = key + # ] + # data <- merge(test, data, by = "key", all.y = TRUE) + + + # flatten duplicated ids ---- + # create unique key + data$key2 <- with(data, paste(key, geomfeatid, existsonfeat)) + + ## find duplicates ---- + idx <- which(duplicated(data$key2)) + + if (length(idx) > 0) { + tb <- table(data$key2) + dups <- which(data$key2 %in% names(tb)[tb > 1]) + + if (length(dups) > 0) { + nodups <- { + len <- {data$key2[dups] |> rle() ->.; .$lengths} + len <- cumsum(len) - len + 1 + } + nodups <- dups[nodups] + } else nodups <- 1:nrow(data) + + # subset duplicates + vars <- c("key2", "landform", "mntn", "hill","trce", "flats", "shapeacross", "shapedown", "slopeshape", "hillslopeprof") + # ..vars = NULL + data_sub <- data[dups, vars, with = FALSE] + + # flatten duplicates + data_sub <- .flatten_gmd(as.data.frame(data_sub), key = "key2") |> data.table::as.data.table() + + # replace duplicates + data[nodups, vars] <- data_sub + + # remove duplicates + data <- data[- idx, ] + } + data$key2 <- NULL + + # subset different conventions ---- + data_comb <- subset(data, n_bot > 1 & co_dir != "missing") + data_mis <- subset(data, n_bot > 1 & co_dir == "missing") + data_simp <- subset(data, n_bot < 2) + + + vars <- c("key", "landform", "mntn", "hill","trce", "flats", "shapeacross", "shapedown", "slopeshape", "hillslopeprof") + # ..vars = NULL + data_mis <- data_mis[, vars, with = FALSE] |> + as.data.frame() |> + .flatten_gmd(sep = " and ") |> + data.table::as.data.table() + data_simp <- data_simp[, vars, with = FALSE] |> + as.data.frame() |> + .flatten_gmd(sep = " on ") |> + data.table::as.data.table() + + + # iterate over sites with unsorted overlapping landforms ---- + if (nrow(data_comb) > 0) { + data_comb_l <- split(data_comb, data_comb$key) + data_comb_l <- lapply(data_comb_l, function(x) { + + # replace landscape existsonfeat with NA + x$existsonfeat <- sapply(x$existsonfeat, function(y) { + ifelse(any(y == x$geomfeatid), y, NA) + }) + + + # find bottom landform + bot <- subset(x, is.na(existsonfeat)) + + + # iterate over bottoms + len <- nrow(x) + sep <- ifelse(x$co_dir == "ordered", " on ", " and ") + + bot <- split(bot, bot$geomfeatid) + x_sorted <- lapply(bot, function(y) { + rank <- integer(len) + rank[1] <- y$geomfeatid[1] + for (i in 1:len) { + idx <- x$geomfeatid[which(x$existsonfeat == rank[i])] + if (length(idx) > 0) rank[i + 1] = idx + } + rank <- rank[which(rank > 0)] + rank <- which(x$geomfeatid %in% rank) + vars <- c("key", "landform", "mntn", "hill","trce", "flats", "shapeacross", "shapedown", "slopeshape", "hillslopeprof") + # ..vars = NULL + x2 <- x[rank, vars, with = FALSE] + suppressMessages(y <- .flatten_gmd(x2, sep = sep, SORT = FALSE)) + return(y) + }) + x_sorted <- do.call("rbind", x_sorted) + + return(x_sorted) + }) + data_comb_l3 <- do.call("rbind", data_comb_l) + data_comb <- .flatten_gmd(as.data.frame(data_comb_l3), key = "key") |> data.table::as.data.table() + } else data_comb <- data_comb[, vars, with = FALSE] + + data <- rbind(data_simp, data_mis, data_comb) |> + as.data.frame() + names(data)[names(data) == "key"] <- key + + + # # uncode + # data("metadata", package = "soilDB") + # vars <- c("mntn", "hill", "trce", "flats") + # idx <- names(data) %in% vars + # names(data)[idx] <- paste0("geompos", vars) + # idx <- names(data) %in% metadata$ColumnPhysicalName + # data + + + return(data) +} + + +# vars <- c("geomfeatid", "existsonfeat") +# idx <- sapply(1:nrow(test), function(i) { +# test[i, vars, drop = TRUE] |> unlist() |> unname() +# }, +# simplify = FALSE +# ) |> unlist() +# idx <- idx[!duplicated(idx) & !is.na(idx)] + + +.format_slopeshape <- function(dat) { + + shapeacross = NA + shapedown = NA + + dat <- within(dat, { + ssa = NA # slope shape across + ssd = NA # slope shape down + slopeshape = NA + + ssa = gsub("Concave", "C", shapeacross) + ssa = gsub("Linear", "L", ssa) + ssa = gsub("Convex", "V", ssa) + + ssd = gsub("Concave", "C", shapedown) + ssd = gsub("Linear", "L", ssd) + ssd = gsub("Convex", "V", ssd) + + slopeshape = paste0(ssd, ssa, sep = "") + slopeshape[slopeshape %in% c("NANA", "")] = NA + }) + dat[c("ssa", "ssd")] <- NULL + + # ss_vars <- c("CC", "CV", "CL", "LC", "LL", "LV", "VL", "VC", "VV") + # if (all(dat$slopeshape[!is.na(dat$slopeshape)] %in% ss_vars)) { + # dat$slopeshape <- factor(dat$slopeshape, levels = ss_vars) + # dat$slopeshape <- droplevels(dat$slopeshape) + # } + return(dat) +} + + +.flatten_gmd <- function(data, key = "key", table = NULL, sep = " and ", SORT = TRUE) { + + idx_key <- grep(key, names(data)) + if (length(idx_key) != 1L) stop("the key/id argument does not match any of the column names in the data.frame") + names(data)[c(idx_key)] <- c("key") + + + tb <- table(data$key) + idx <- names(tb[tb > 1]) + idx <- data$key %in% idx + + test <- sum(idx, na.rm = TRUE) + if (test > 0) { + message(test, " ", key, " values were found in the ", table, " table that contain multiple entries, the resulting values will be flattened/combined into 1 record per ", key, " and separated with 'and'") + + data_sub <- data[idx, ] |> data.table::as.data.table() + .SD = NULL + data_sub <- data_sub[ + , + lapply(.SD, function(x) { + if (SORT) {paste0(sort(unique(x[!is.na(x)])), collapse = sep) + } else {paste0( unique(x[!is.na(x)]), collapse = sep)} + }), + by = key + ] |> + as.data.frame() + + data <- rbind(data[!idx, ], data_sub) + + } + + # replace "" values with NA + idx <- 1:ncol(data) + data[idx] <- lapply(data, function(x) ifelse(x == "", NA, x)) + + + names(data)[idx_key] <- key + + return(data) +} + + diff --git a/misc/gdb_examples.R b/misc/gdb_examples.R index 14eec910..b76ed580 100644 --- a/misc/gdb_examples.R +++ b/misc/gdb_examples.R @@ -1,34 +1,150 @@ +library(aqp) library(soilDB) library(sf) +source("C:/workspace2/github/ncss-tech/soilDB/R/get_component_from_GDB.R") +source("C:/workspace2/github/ncss-tech/soilDB/R/utils.R") + + # new GDB examples dsn <- "D:/geodata/soils/gNATSGO_CONUS_Oct2021/gNATSGO_CONUS.gdb" +dsn <- "D:/geodata/soils/gSSURGO_CONUS_202210.gdb" + + +# legend ---- +le1 <- get_legend_from_GDB(dsn = dsn) +le2 <- get_legend_from_GDB(dsn = dsn, stats = TRUE) +le3 <- get_legend_from_GDB(dsn = dsn, WHERE = "areasymbol LIKE 'CA%'") +le4 <- get_legend_from_GDB(dsn = dsn, WHERE = "areasymbol LIKE 'CA%'", stats = TRUE) + +le12 <- soilDB::get_legend_from_GDB(dsn = dsn) +le22 <- soilDB::get_legend_from_GDB(dsn = dsn, stats = TRUE) +le32 <- soilDB::get_legend_from_GDB(dsn = dsn, WHERE = "areasymbol LIKE 'CA%'") +le42 <- soilDB::get_legend_from_GDB(dsn = dsn, WHERE = "areasymbol LIKE 'CA%'", stats = TRUE) + + +daff::diff_data(le1, le12) |> daff::render_diff() +daff::diff_data(le2, le22) |> daff::render_diff() +daff::diff_data(le3, le32) |> daff::render_diff() +daff::diff_data(le4, le42) |> daff::render_diff() + + +le12 <- get_legend_from_SDA(WHERE = "areasymbol != 'US'") +le32 <- get_legend_from_SDA(WHERE = "areasymbol LIKE 'CA%'") +le42 <- get_legend_from_SDA(WHERE = "areasymbol LIKE 'CA%'") + + + +# mapunit ---- +mu1 <- get_mapunit_from_GDB(dsn = dsn) +mu2 <- get_mapunit_from_GDB(dsn = dsn, stats = TRUE) +mu3 <- get_mapunit_from_GDB(dsn = dsn, WHERE = "muname LIKE 'Miami%'") +mu4 <- get_mapunit_from_GDB(dsn = dsn, WHERE = "muname LIKE 'Miami%'", stats = TRUE) +mu5 <- get_mapunit_from_GDB(dsn = dsn, WHERE = "areasymbol LIKE 'IN%'") +mu6 <- get_mapunit_from_GDB(dsn = dsn, WHERE = "areasymbol LIKE 'IN%'", stats = TRUE) +mu7 <- get_mapunit_from_GDB(dsn = dsn, WHERE = "compname = 'Miami'") +mu8 <- get_mapunit_from_GDB(dsn = dsn, WHERE = "muname LIKE 'Miami%' AND areasymbol = 'IN001'") + +mu12 <- soilDB::get_mapunit_from_GDB(dsn = dsn) +mu22 <- soilDB::get_mapunit_from_GDB(dsn = dsn, stats = TRUE) +mu32 <- soilDB::get_mapunit_from_GDB(dsn = dsn, WHERE = "muname LIKE 'Miami%'") +mu42 <- soilDB::get_mapunit_from_GDB(dsn = dsn, WHERE = "muname LIKE 'Miami%'", stats = TRUE) +mu52 <- soilDB::get_mapunit_from_GDB(dsn = dsn, WHERE = "areasymbol LIKE 'IN%'") +mu62 <- soilDB::get_mapunit_from_GDB(dsn = dsn, WHERE = "areasymbol LIKE 'IN%'", stats = TRUE) +mu72 <- soilDB::get_mapunit_from_GDB(dsn = dsn, WHERE = "compname = 'Miami'") +mu82 <- soilDB::get_mapunit_from_GDB(dsn = dsn, WHERE = "muname LIKE 'Miami%' AND areasymbol = 'IN001'") + +daff::diff_data(mu1, mu12) |> daff::render_diff() +daff::diff_data(mu2, mu22) |> daff::render_diff() +daff::diff_data(mu3, mu32) |> daff::render_diff() +daff::diff_data(mu4, mu42) |> daff::render_diff() +daff::diff_data(mu5, mu52) |> daff::render_diff() +daff::diff_data(mu6, mu62) |> daff::render_diff() + + + +# component ---- +co1 <- get_component_from_GDB(dsn) +co2 <- get_component_from_GDB(dsn, childs = TRUE) +co3 <- get_component_from_GDB(dsn, WHERE = "compname = 'Miami' AND majcompflag = 'Yes'", childs = FALSE) + +co12 <- soilDB::get_component_from_GDB(dsn) +co22 <- soilDB::get_component_from_GDB(dsn, childs = TRUE) +co32 <- soilDB::get_component_from_GDB(dsn, WHERE = "compname = 'Miami' AND majcompflag = 'Yes'", childs = FALSE) -le <- get_legend_from_GDB(dsn = dsn) -le <- get_legend_from_GDB(dsn = dsn, stats = TRUE) -le <- get_legend_from_GDB(dsn = dsn, WHERE = "areasymbol LIKE 'CA%'") -le <- get_legend_from_GDB(dsn = dsn, WHERE = "areasymbol LIKE 'CA%'", stats = TRUE) +vars <- names(co12)[names(co12) %in% names(co2)] +summary(d1 <- daff::diff_data(co1[vars], co12)) +summary(d3 <- daff::diff_data(co3[vars], co32)) +co2 <- co2[order(co2$mukey, co2$cokey), ] +co22 <- co22[order(co22$mukey, co2$cokey), ] +vars2 <- names(co22)[names(co22) %in% names(co2)] +idx <- sample(1:length(co2$mukey), 1000) +all(co2$mukey == co22$mukey) +all(co2$cokey == co22$cpkey) +test1 <- co2[idx, vars2] +test2 <- co22[idx, vars2] +wtf <- daff::diff_data(test1, test2) -# mapunit -mu <- get_mapunit_from_GDB(dsn = dsn) -mu <- get_mapunit_from_GDB(dsn = dsn, stats = TRUE) -mu <- get_mapunit_from_GDB(dsn = dsn, WHERE = "muname LIKE 'Miami%'") -mu <- get_mapunit_from_GDB(dsn = dsn, WHERE = "muname LIKE 'Miami%'", stats = TRUE) -mu <- get_mapunit_from_GDB(dsn = dsn, WHERE = "areasymbol LIKE 'IN%'") -mu <- get_mapunit_from_GDB(dsn = dsn, WHERE = "areasymbol LIKE 'IN%'", stats = TRUE) -mu <- get_mapunit_from_GDB(dsn = dsn, WHERE = "compname = 'Miami'") -mu <- get_mapunit_from_GDB(dsn = dsn, WHERE = "muname LIKE 'Miami%' AND areasymbol = 'IN001'") +test0 <- merge(co2[c("cokey", "pmkind", "pmorigin")], co22[c("cokey", "pmkind", "pmorigin")], by = "cokey", all.x = TRUE) +test0[test0[3] != test0[5] & !is.na(test0[5]), ] |> View() -co <- get_component_from_GDB(dsn, WHERE = "compname = 'Miami' AND majcompflag = 'Yes'", childs = FALSE) +# fetchGDB ---- +f0 <- fetchGDB(dsn = "D:/geodata/soils/gSSURGO_RI_202110/gSSURGO_RI.gdb") +f00 <- fetchGDB(dsn = "D:/geodata/soils/gSSURGO_RI_202110/gSSURGO_RI.gdb", childs = TRUE) +f1 <- fetchGDB(dsn, WHERE = "areasymbol = 'IN001'") +f2 <- fetchGDB(dsn, WHERE = "areasymbol = 'IN001'", childs = TRUE) +# f3 <- fetchGDB(dsn, WHERE = "areasymbol LIKE 'CA%'", childs = TRUE) +f4 <- fetchGDB(dsn, WHERE = "muname LIKE 'Miami%'") +f5 <- fetchGDB(dsn, WHERE = "muname LIKE 'Miami%' AND areasymbol = 'IN001'") +f6 <- fetchGDB(dsn, WHERE = "Compname LIKE 'Miami%'") + + +f02 <- soilDB::fetchGDB(dsn = "D:/geodata/soils/gSSURGO_RI_202110/gSSURGO_RI.gdb") +f002 <- soilDB::fetchGDB(dsn = "D:/geodata/soils/gSSURGO_RI_202110/gSSURGO_RI.gdb", childs = TRUE) +f12 <- soilDB::fetchGDB(dsn, WHERE = "areasymbol = 'IN001'") +f22 <- soilDB::fetchGDB(dsn, WHERE = "areasymbol = 'IN001'", childs = TRUE) +# f32 <- soilDB::fetchGDB(dsn, WHERE = "areasymbol LIKE 'CA%'", childs = TRUE) +f42 <- soilDB::fetchGDB(dsn, WHERE = "muname LIKE 'Miami%'") +f52 <- soilDB::fetchGDB(dsn, WHERE = "muname LIKE 'Miami%' AND areasymbol = 'IN001'") +f62 <- soilDB::fetchGDB(dsn, WHERE = "Compname LIKE 'Miami%'") + +# save(f0, f00, f1, f2, f4, f6, f02, f002, f12, f22, f42, f62, file = "fetchGDA_test.RData") +load(file = "fetchGDA_test.RData") + + +test <- list( + list(f0, f02), + list(f00, f002), + list(f1, f12), + list(f2, f22), + # list(f3, f32), + list(f4, f42), + list(f6, f62) +) + + +test <- lapply(test, function(x) { + + s1 = x[[1]]@site + s2 = x[[2]]@site + + h1 = x[[1]]@horizons + h2 = x[[2]]@horizons + + s <- daff::diff_data(s1, s2) + h <- daff::diff_data(h1, h2) + return(list(s = s, h = h)) +}) + + -f_in001 <- fetchGDB(dsn, WHERE = "areasymbol = 'IN001'") f_sda <- fetchSDA(WHERE = "areasymbol = 'IN001'") -f_cokey <- fetchGDB(dsn, WHERE = "cokey = '18224902'") +f_cokey <- fetchGDB(dsn, WHERE = "cokey = '22617647'") f_mi <- fetchGDB(dsn, WHERE = "areasymbol LIKE 'MI%'") @@ -38,3 +154,28 @@ save(f_us, file = "gnatsgo.RData") WHERE <- paste0("areasymbol IN ('", paste0(le2, collapse = "', '"), "')") system.time(test <- fetchGDB(dsn, WHERE = WHERE)) + + +# formatLandformString vs .cogmd_prep2 ---- +vars <- c("cokey", "landscape", "landform", "mntn", "hill","trce", "flats", "shapeacross", "shapedown", "slopeshape", "hillslopeprof") +test <- cogmd[vars] +names(test)[c(3, 6:8)] <- c("geomfname", paste0("geompos", vars[4:6])) +test <- within(test, { + geomicrorelief = NA + geommicelev = NA + geomfmod = NA + geomfiidref = NA +}) +test <- as.data.table(test) +test2 <- test[, soilDB:::.formatLandformString(.SD, uid = .BY$cokey, name.sep = '|'), + by = list(cokey = test$cokey)] + +ed.lf <- data.table::as.data.table(test$geomorph) +lf <- ed.lf[, soilDB:::.formatLandformString(.SD, uid = .BY$peiid, name.sep = '|'), + by = list(peiid = ed.lf$peiid)] + +system.time(.cogmd_prep2(cogmd)) + + + +