9 Define Strata
In this target, strata definitions are assigned to the numerous logging coupes within the project area. Mosaic planting is associated with less steep regions and therefore those coupes which largely intersected the “flat” regions, defined in Chapter 8, were classified as Mosaic. Where mosaic regions had a slope between 15 \(\textdegree\) and 25 \(\textdegree\), they were classified as mosaic-unplanted, regions with a slope < 15 \(\textdegree\) were classified as mosaic-planted.
Any regions with a slope > 25 \(\textdegree\) were classified as no-go and excluded from the claimable region. Any regions within: 100 m of a main river, 30 m of a small river, 10 m of a road (excluding tracks), or 50 m from a protected area, were also classified as no-go and removed from the claimable area.
All remaining coupes were classified as being under Natural Forest Management (NFM)
These strata are shown in the following interactive map in addition to the “flat areas”, derived from the k-means slope analysis in Chapter 8.
The areas for each strata type are shown in the Table 9.1. Note that the area of the gazetted project area is 83381 ha, however measured area, based on project boundary is 84800 ha. In order to align with the official gazetted area, all credit generating calculations in Chapter 20 and Chapter 21 are based on measured areas adjusted by the percentage difference between the gazetted and measured areas. This adjustment is also applied to the strata areas shown in the Table 9.1 to illustrate the effect of this adjustment. The project started in 2015, but was extended in 2016 by 15511 ha was in 2016. Table 9.1 shows the area of each strata type for the pre-2016 and post-2016 project areas - all downstream calculations, relevant to the first year of the project, account for this difference in project area.
9.1 Code
Define Strata Code R/DefineStrata.R
#' Define The Strata
#'
#' Define the baseline scenario management strata for all logging coupes.
#'
#' @param slope_src file path for slope raster from `R/slopeAnalysis.R`
#' @param boundary_src file path for the aoi boundary
#' @param flat_zones_src file path for zones considered flat from
#' `R/slopeAnalysis.R`
#' @param logg_coups_src file path for logging coupe data
#' @param roads_path file path for road network
#' @param rivers_path file path for river network
#' @param protectedArea_path file path for protected area
#' @param KuamutRiver_path file path for the Kuamut river
#' @param roadsBuffer buffer distance for roads in meters; default is 10
#' @param riversBuffer buffer distance for rivers in meters; default is 30
#' @param BigRiverBuffer buffer distance for main rivers in meters;
#' default is 100
#' @param protectedBuffer buffer distance for protected areas in meters;
#' default is 50
#'
#' @return shapefile paths for each strata No-Go areas, The two Mosaic strata
#' and NFM Strata
#'
#' @details
#' This file takes Kuamut Historical logging coupe data file from Eva,
#' (tidied original from RBJ in 2021) and compares it to the cluster analysis
#' note Eva matched files to our boundary as they where inconsistent.
define_strata <- function(slope_src, boundary_src, flat_zones_src,
logg_coups_src, roads_path, rivers_path,
protectedArea_path, KuamutRiver_path,
roadsBuffer = 10, riversBuffer = 30,
BigRiverBuffer = 100, protectedBuffer = 50) {
slp <- terra::rast(slope_src)
Boundary <- sf::read_sf(boundary_src)
boundary_area <- sf::st_area(Boundary) |> units::set_units(ha)
# Flat contiguous areas identified suitable for mosaic
ClusterFlat <- sf::read_sf(flat_zones_src)
# Load in Kuamut logging coupes data.
Kcoupes2 <- sf::read_sf(logg_coups_src)
Kcoupes2_area <- sum(sf::st_area(Kcoupes2)) |> units::set_units(ha)
diff <- boundary_area - Kcoupes2_area
units::set_units(diff, m^2)
## subset overlapping coupes by hand as forestry do
Kcoupes3 <- select(Kcoupes2, -c(
"FMU_NO", "RBJ_REG", "RECNO", "LOGGER", "GENERAT",
"CoupeNm", "Cntrctr", "OprtnTy", "TtlPr_3", "Remarks",
"CttngLm", "LggngSy"
))
# These Coupes for Mosaic part of Mosaic (not including NFM within mosaic)
selected4Mosaic <- c(
"YSH7/06",
"YL3/02", "YL2/04", "YL2/01",
"YL3/03",
"YLH3/06(2)Area B", "YLH3/06(1)",
"YLH1/05", "YLH2/05", "YSH8/06(Area A) Partial",
"Steep area YLH406", "YLH4/06", "YS9/06", "YSH2/06", "YLH4/05",
"YL2/00", "YL1/01"
)
MosaicCoupes <- Kcoupes3 |>
dplyr::filter(BLOCK_N %in% selected4Mosaic)
c_flat_area <- sum(sf::st_area(ClusterFlat)) |> units::set_units(ha)
mos_coupes_area <- sum(sf::st_area(MosaicCoupes)) |> units::set_units(ha)
# define NFM zones.
NFMCoupes <- Kcoupes3 |>
dplyr::filter(!BLOCK_N %in% selected4Mosaic)
# ------------------- REMOVE NO-GO STRATA FROM BOTH AREAS -------------------
# rivers and buffer, roads and buffer, Protected area and buffer,
# above 25 degrees
#----------------------------------------------------------------------------#
# function to read and buffer an sf feature.
read_buff <- function(x, y) {
sf::read_sf(x) |>
sf::st_buffer(y) |>
sf::st_geometry()
}
# import roads and rivers layers for both main area and extension (all)
# From digitisation. Buffer on load.
roads_buf <- read_buff(roads_path, roadsBuffer)
rivers_buf <- read_buff(rivers_path, riversBuffer)
BigRiver_buf <- read_buff(KuamutRiver_path, BigRiverBuffer)
protect_buf <- read_buff(protectedArea_path, protectedBuffer)
# ----------- Separate steep slope Strata polygon -------------------------
# 1 ha slope (100m pixels). Steep areas above 25 degrees (will always be
# no-go). i.e not logged in the baseline
#----------------------------------------------------------------------------#
# Reclass table: if its in, gets NA, if its out gets 1;above 25 degrees,
# is always out
rclmat <- matrix(c(0, 25, NA, 25, 100, 1), ncol = 3, byrow = TRUE)
# reclassify slope bands, go and no-go
reclass_slp <- terra::classify(slp, rclmat)
# convert raster to polygon
abv25out <- terra::as.polygons(reclass_slp) |>
sf::st_as_sf() |>
sf::st_geometry() |>
sf::st_as_sf()
# combine all No-Go areas (make others st-objects also so can rbind)
NoGo <- rbind(
roads_buf, rivers_buf, BigRiver_buf,
protect_buf, abv25out
) |>
sf::st_union() |> # make a single feature
sf::st_intersection(Boundary) # crop to boundary
# Total area of NoGo region.
units::set_units(x = sum(sf::st_area(NoGo)), value = ha)
## Remove No-Go strata from NFM strata
NFMCoupesGo <- sf::st_difference(NFMCoupes, NoGo)
## Remove No-Go strata from mosaic
MosaicCoupesGo <- sf::st_difference(MosaicCoupes, NoGo)
# ------ Separate Mosaic area into planted and unplanted strata -------------#
# Reclass table: if its in, gets NA, if its out gets 1. 0-15 will be planted
rclmat2 <- matrix(c(0, 15, NA, 15, 25, 1), ncol = 3, byrow = TRUE)
# reclassify slope bands 4 grid
reclass_slp2 <- terra::classify(slp, rclmat2)
# convert raster to polygon
outOfPlant <- terra::as.polygons(reclass_slp2) |>
sf::st_as_sf() |>
sf::st_geometry() |>
sf::st_as_sf() |>
sf::st_union()
# take only the particular strata of the mosaic (i.e minus out)
MosaicPLANT <- sf::st_difference(MosaicCoupesGo, outOfPlant)
MosaicUNplanted <- sf::st_difference(
MosaicCoupesGo,
sf::st_union(MosaicPLANT)
)
# -------------------------- Save Shape files ------------------------------#
feature_list <- list(
NoGo, NFMCoupesGo,
MosaicCoupesGo, MosaicPLANT, MosaicUNplanted
)
file_list <- list(
"NoGo/NoGo.shp", "/NFMCoupes/NFMCoupesGo.shp",
"MosaicCoupes/MosaicCoupesGo.shp", "MosaicPLANT/MosaicPLANT.shp",
"MosaicUNplanted/MosaicUNplanted.shp"
)
file_list <- purrr::map2(feature_list, file_list, ~ save_feature(.x, .y))
names(file_list) <- c(
"NoGo", "NFMCoupesGo", "MosaicCoupesGo", "MosaicPLANT",
"MosaicUNplanted"
)
# merge all strata-classified coupes into a single sf object and dissolve
# then save output to fgb for use in ML pipeline/downstream.
coupes_merge <- purrr::map2(
feature_list, names(file_list),
\(.x, .y) {
v <- sf::st_as_sf(.x) |>
dplyr::mutate(mgmt_strata = .y)
st_geometry(v) <- "geometry"
return(v)
}
) |>
dplyr::bind_rows() |>
dplyr::group_by(mgmt_strata) |>
dplyr::summarise()
mgmt_strat_path <- save_feature(
coupes_merge,
"mgmt_strata_merge/coupes_merge.fgb"
)
return(c(file_list, mgmt_strata_union = mgmt_strat_path))
}Code to generate aggregated strata for pre and post 2016; R/split_strata_ext_region.R
#' @title Split strata into pre- and post-2016
#' @param coupes_src Path to kuamut coupes spat vector
#' @param strata_src Path to strata spat vector
#' @details
#' This function takes the aggregated strata for the whole project area and
#' removes the region which was not purchased in the first year of the project.
#' this enables the correct areal calculation of project sequestration for the
#' first year of the project.
split_strata <- function(coupes_src, strata_src) {
coupes <- read_sf(coupes_src)
strata <- read_sf(strata_src)
strat_union <- dplyr::filter(strata, mgmt_strata != "MosaicCoupesGo")
total <- sf::st_union(strat_union) |>
sf::st_area() |>
units::set_units(ha) |>
as.numeric()
GazettedArea <- 69454 + 13927
TakePercent <- GazettedArea / total
coupes_add <- coupes |>
dplyr::filter(BLOCK_N %in% c(
"YLH3/06(2)Area B", "YLH3/06(1)", "YLH2/05", "YLH1/05",
"Steep area YLH406", "YLH4/06", "YS4/05(1)", "YT1/06(I)"
)) |>
sf::st_union()
coupes15 <- sf::st_difference(
sf::st_union(coupes),
coupes_add
) |>
wk::wk_flatten() |>
st_as_sf()
strat_union_15 <- sf::st_intersection(strat_union, coupes15) |>
dplyr::select(mgmt_strata) |>
wk::wk_flatten()
strat_union_15 <- strat_union_15 |>
dplyr::filter(st_geometry_type(strat_union_15) == "POLYGON") |>
group_by(mgmt_strata) |>
dplyr::summarise()
add_attrs <- function(x) {
dplyr::mutate(x,
area_ha = units::set_units(sf::st_area(x), "ha"),
area_adj = area_ha * TakePercent,
mgmt_strata = dplyr::case_when(
mgmt_strata == "NFMCoupesGo" ~ "NFM",
mgmt_strata == "MosaicPLANT" ~ "Mosaic-planted",
mgmt_strata == "MosaicUNplanted" ~ "Mosaic-unplanted",
mgmt_strata == "NoGo" ~ "Inoperable (No Go)"
)
)
}
strat_union_15 <- strat_union_15 |>
add_attrs() |>
save_feature("mgmt_strata_merge/strata_union_pre16.fgb")
strat_union <- strat_union |>
add_attrs() |>
save_feature("mgmt_strata_merge/strata_union_post16.fgb")
return(list(
strata_15_16 = strat_union_15,
strata_post16 = strat_union
))
}