Skip to content

Commit

Permalink
Merge branch 'development' of github.com:eliotmcintire/Boreal_LBMRDat…
Browse files Browse the repository at this point in the history 
…aPrep into development
  • Loading branch information
@eliotmcintire
eliotmcintire committed Nov 7, 2018
2 parents 242bb14 + b75e536 commit f5b4fa0
Show file tree
Hide file tree
Showing 4 changed files with 83 additions and 190 deletions.
99 changes: 60 additions & 39 deletions Boreal_LBMRDataPrep.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -135,23 +135,23 @@ checkSpeciesLayerNames <- function(specieslayers) {
}
names(specieslayers) <- speciesEquivalencyTable$Latin[matches]
}


if (FALSE) {
# not needed, but this is for calculating vegetation type maps, from a species abundances stack
VTM <- Cache(pemisc::makeVegTypeMap, sim$specieslayers, sim$vegLeadingProportion)
tableCC <- table(factorValues(CCvtm, CCvtm[], att = "Species"))
tablePP <- table(factorValues(VTM, VTM[], att = "Species"))

propCC <- round(tableCC/sum(tableCC),2)
propPP <- round(tablePP/sum(tablePP),2)
names(propPP) <- speciesEquivalencyTable$Common[pmatch(names(propPP), speciesEquivalencyTable$Latin)]

propCC[match(names(propPP), names(propCC))]
propPP
}
}

return(specieslayers)
return(specieslayers)
}

estimateParameters <- function(sim) {
Expand All @@ -178,14 +178,19 @@ estimateParameters <- function(sim) {
ecoregion = 1:1031)

message("ecoregionProducer: ", Sys.time())
# Note: this ecoregionMap is NOT the Canadian EcoRegion -- it is for LBMR, which uses "ecoregion"
ecoregionMap <- Cache( postProcess, sim$ecoDistrict,
studyArea = sim$shpStudyArea)
ecoregionFiles <- Cache(ecoregionProducer,
studyAreaRaster = initialCommFiles$initialCommunityMap,
ecoregionMapFull = sim$ecoDistrict,
#studyAreaRaster = initialCommFiles$initialCommunityMap,
#ecoregionMapFull = sim$ecoDistrict,
ecoregionMap = ecoregionMap,
ecoregionName = "ECODISTRIC",
ecoregionActiveStatus = ecoregionstatus,
studyArea = sim$studyArea,
rstStudyArea = rstStudyRegionBinary,
maskFn = fastMask,
rasterToMatch = initialCommFiles$initialCommunityMap, #sim$rasterToMatch,
#studyArea = sim$studyArea,
#rstStudyArea = rstStudyRegionBinary,
#maskFn = fastMask,
userTags = "stable")

message("3: ", Sys.time())
Expand Down Expand Up @@ -444,14 +449,50 @@ Save <- function(sim) {

polyCenter <- SpatialPoints(coords = data.frame(x = c(-1349980), y = c(6986895)),
proj4string = crsUsed)

seedToKeep <- .GlobalEnv$.Random.seed
set.seed(1234)
sim$shpStudyAreaLarge <- SpaDES.tools::randomPolygon(x = polyCenter, hectares = 10000)
.GlobalEnv$.Random.seed <- seedToKeep
}

needRstSR <- FALSE
if (!suppliedElsewhere("rasterToMatch", sim)) {
needRstSR <- TRUE
if (!suppliedElsewhere("shpStudyArea", sim)) {
message("'shpStudyArea' was not provided by user. Using the same as 'shpStudyAreaLarge'")
sim$shpStudyArea <- sim$shpStudyAreaLarge
}

needRTM <- FALSE
if (is.null(sim$rasterToMatch)) {
if (!suppliedElsewhere("rasterToMatch", sim)) {
needRTM <- TRUE
message("There is no rasterToMatch supplied; will attempt to use biomassMap")
} else {
stop("rasterToMatch is going to be supplied, but ", currentModule(sim), " requires it ",
"as part of its .inputObjects. Please make it accessible to ", currentModule(sim),
" in the .inputObjects by passing it in as an object in simInit(objects = list(rasterToMatch = aRaster)",
" or in a module that gets loaded prior to ", currentModule(sim))
}
}
if (needRstSR) {

if (!suppliedElsewhere("biomassMap", sim) || needRTM) {
sim$biomassMap <- Cache(prepInputs,
targetFile = asPath(basename(biomassMapFilename)),
archive = asPath(c("kNN-StructureBiomass.tar",
"NFI_MODIS250m_kNN_Structure_Biomass_TotalLiveAboveGround_v0.zip")),
#url = extractURL("biomassMap"),
destinationPath = dPath,
studyArea = sim$shpStudyArea,
rasterToMatch = sim$rasterToMatch,
useSAcrs = TRUE,
method = "bilinear",
datatype = "INT2U",
filename2 = TRUE, overwrite = TRUE,
userTags = c("stable", currentModule(sim)))
}

if (needRTM) {
# if we need rasterToMatch, that means a) we don't have it, but b) we will have biomassMap
sim$rasterToMatch <- sim$biomassMap
message(" Rasterizing the shpStudyAreaLarge polygon map")
if (!is(sim$shpStudyAreaLarge, "SpatialPolygonsDataFrame")) {
dfData <- if (is.null(rownames(sim$shpStudyAreaLarge))) {
Expand Down Expand Up @@ -481,11 +522,6 @@ Save <- function(sim) {
datatype = "INT2U", overwrite = TRUE)
}

if (!suppliedElsewhere("shpStudyArea", sim)) {
message("'shpStudyArea' was not provided by user. Using the same as 'shpStudyAreaLarge'")
sim$shpStudyArea <- sim$shpStudyAreaLarge
}

if (!identical(crsUsed, crs(sim$shpStudyAreaLarge))) {
sim$shpStudyAreaLarge <- spTransform(sim$shpStudyAreaLarge, crsUsed) #faster without Cache
}
Expand All @@ -496,21 +532,6 @@ Save <- function(sim) {

cacheTags = c(currentModule(sim), "function:.inputObjects", "function:spades")

if (!suppliedElsewhere("biomassMap", sim)) {
sim$biomassMap <- Cache(prepInputs,
targetFile = biomassMapFilename,
archive = asPath(c("kNN-StructureBiomass.tar",
"NFI_MODIS250m_kNN_Structure_Biomass_TotalLiveAboveGround_v0.zip")),
url = extractURL("biomassMap"),
destinationPath = dPath,
studyArea = sim$shpStudyArea,
useSAcrs = TRUE,
method = "bilinear",
datatype = "INT2U",
filename2 = TRUE,
userTags = c("stable", currentModule(sim)))
}

# LCC2005
if (!suppliedElsewhere("LCC2005", sim)) {
sim$LCC2005 <- Cache(prepInputs,
Expand All @@ -522,7 +543,7 @@ Save <- function(sim) {
rasterToMatch = sim$rasterToMatch,
method = "bilinear",
datatype = "INT2U",
filename2 = TRUE,
filename2 = TRUE, overwrite = TRUE,
userTags = currentModule(sim))

projection(sim$LCC2005) <- projection(sim$rasterToMatch)
Expand Down Expand Up @@ -575,17 +596,17 @@ Save <- function(sim) {
# stand age map
if (!suppliedElsewhere("standAgeMap", sim)) {
sim$standAgeMap <- Cache(prepInputs, #notOlderThan = Sys.time(),
targetFile = standAgeMapFilename,
targetFile = basename(standAgeMapFilename),
archive = asPath(c("kNN-StructureStandVolume.tar",
"NFI_MODIS250m_kNN_Structure_Stand_Age_v0.zip")),
destinationPath = dPath,
url = extractURL("standAgeMap"),
#url = extractURL("standAgeMap"),
fun = "raster::raster",
studyArea = sim$shpStudyAreaLarge,
rasterToMatch = sim$rasterToMatch,
method = "bilinear",
datatype = "INT2U",
filename2 = TRUE,
filename2 = TRUE, overwrite = TRUE,
userTags = c("stable", currentModule(sim)))
}

Expand Down Expand Up @@ -613,7 +634,7 @@ Save <- function(sim) {
writeRaster(specieslayersList$specieslayers, file.path(outputPath(sim), "speciesLayers.grd"), overwrite = TRUE)
sim$specieslayers <- specieslayersList$specieslayers
sim$speciesList <- specieslayersList$speciesList
}
}

# 3. species maps
sim$speciesTable <- Cache(prepInputs, "speciesTraits.csv",
Expand Down
16 changes: 8 additions & 8 deletions Boreal_LBMRDataPrep.Rmd
View file
Original file line number Diff line number Diff line change
Expand Up @@ -62,11 +62,10 @@ paths <- getPaths()

```{r get-study-area}
library(raster)
cachePath <- file.path("Boreal_LBMRDataPrep", "cache")
modulePath <- Cache(readline, paste0("Where is the module path? (e.g., ~/module, with no quotes).\n",
"Press Enter to accept the path in getPaths()$modulePath: "),
cacheRepo = cachePath)
setPaths(cachePath = cachePath, modulePath = modulePath)
# modulePath <- Cache(readline, paste0("Where is the module path? (e.g., ~/module, with no quotes).\n",
# "Press Enter to accept the path in getPaths()$modulePath: "),
# cacheRepo = cachePath)
# setPaths(cachePath = cachePath, modulePath = modulePath)
## do you want to hand-draw a map or use defaults?
# - note that large areas will take longer to compute
Expand All @@ -89,16 +88,17 @@ if (handDrawMap) {
shpStudyAreaLarge <- SpatialPolygons(list(Polygons(list(Polygon(severalrandompoints$coords)), ID = 1)),
proj4string = crs(canadaMap))
}
Plot(shpStudyAreaLarge, addTo = "canadaMap", col = "red")
}
Plot(shpStudyAreaLarge, addTo = "canadaMap", col = "red")
times <- list(start = 0, end = 10)
modules <- list("Boreal_LBMRDataPrep")
objects <- if (handDrawMap) list("shpStudyAreaLarge" = shpStudyAreaLarge,
"shpStudyArea" = shpStudyAreaLarge) else list()
mySim <- simInit(times = times, params = parameters, modules = modules,
objects = objects)
mySim <- simInit(times = times, #params = parameters,
modules = append(modules, "LBMR"),
objects = objects, paths = getPaths())
```

# Run `spades`
Expand Down
45 changes: 6 additions & 39 deletions R/ecoregionProducers.R
View file
Original file line number Diff line number Diff line change
@@ -1,48 +1,15 @@
ecoregionProducer <- function(studyAreaRaster, ecoregionMapFull, ecoregionName,
ecoregionActiveStatus, studyArea, rstStudyArea, maskFn) {
ecoregionProducer <- function(ecoregionMap, ecoregionName,
ecoregionActiveStatus, rasterToMatch) {
# change the coordinate reference for all spatialpolygons
message("ecoregionProducer 1: ", Sys.time())
ecoregionMapInStudy <- raster::intersect(ecoregionMapFull, fixErrors(aggregate(studyArea)))
# ecoregions <- ecoregionMapInStudy@data[,ecoregionName]
# ecoregionTable <- data.table(mapcode = numeric(),
# ecoregion = character())
# mapcode <- 1
# for(ecoregion in unique(ecoregions)){
# # for(ecoregion in ecoregions){
# singleecoMapPoly <- ecoregionMapInStudy[ecoregionMapInStudy@data[,ecoregionName]==ecoregion,]
# studyAreaRaster <- setValues(studyAreaRaster, mapcode)
# singleecoMapRaster <- crop(studyAreaRaster, singleecoMapPoly)
# singleecoMapRaster <- suppressWarnings(maskFn(singleecoMapRaster, singleecoMapPoly))
# if(length(unique(getValues(singleecoMapRaster)))==1){
# if(is.na(unique(getValues(singleecoMapRaster)))){
# ecoregionTable <- rbind(ecoregionTable,
# data.table(mapcode = NA,
# ecoregion = ecoregion))
# } else {
# ecoregionTable <- rbind(ecoregionTable,
# data.table(mapcode = mapcode,
# ecoregion = ecoregion))
# }
# } else {
# ecoregionTable <- rbind(ecoregionTable,
# data.table(mapcode = mapcode,
# ecoregion = ecoregion))
# }
#
# if(mapcode == 1){
# ecoregionMap <- singleecoMapRaster
# } else {
# ecoregionMap <- merge(ecoregionMap, singleecoMapRaster)
# }
# mapcode <- mapcode + 1
# }
#ecoregionMapInStudy <- raster::intersect(ecoregionMapFull, fixErrors(aggregate(studyArea)))

# Alternative
message("ecoregionProducer fastRasterize: ", Sys.time())
ecoregionMap <- fasterize::fasterize(sf::st_as_sf(ecoregionMapInStudy), studyAreaRaster, field = "ECODISTRIC")
ecoregionMap <- fasterize::fasterize(sf::st_as_sf(ecoregionMap), raster(rasterToMatch), field = "ECODISTRIC")
ecoregionMap[is.na(rasterToMatch[])] <- NA

#ecoregionMap1 <- rasterize(ecoregionMapInStudy, studyAreaRaster, field = "ECODISTRIC")
ecoregionFactorValues <- unique(ecoregionMap[])
ecoregionFactorValues <- na.omit(unique(ecoregionMap[]))

ecoregionTable <- data.table(mapcode = seq_along(ecoregionFactorValues[!is.na(ecoregionFactorValues)]),
ecoregion = as.numeric(ecoregionFactorValues[!is.na(ecoregionFactorValues)]))
Expand Down
Loading

0 comments on commit f5b4fa0

Please sign in to comment.