knitr::opts_chunk$set(echo = TRUE)
library(aceecostats)
library(sp)
data(aes_region_ll)

Previously we classified all pixels based on their region polygon.

library(raadtools)
topo <- readtopo("etopo2")
library(tabularaster)
cn <- cellnumbers(topo, aes_region_ll)
saveRDS(cn, file = "cache/cn_ETOPO2_cellnumbers_aes_region_ll.rds")
#xy <- xyFromCell(topo, cn$cell_)
knitr::kable(as.data.frame(aes_region))

Bathymetry

library(raadtools)
topo <- readAll(readtopo("etopo2"))
library(dplyr)
cn <- readRDS("cn_ETOPO2_cellnumbers_aes_region_ll.rds")
cn <- cn %>% mutate(topoval = extract(topo, cn$cell_))
cn <- cn %>% mutate(area = extract(area(topo), cell_))
saveRDS(cn, file = "cache/cn_ETOPO2_cellnumbers_area_topo_aes_region_ll.rds")
library(dplyr)
cn <- readRDS("../cache/cn_ETOPO2_cellnumbers_area_topo_aes_region_ll.rds")
depths <- cn %>% group_by(object_) %>% summarize(topo = mean(topoval)) %>% arrange(as.integer(object_))

aes_region$depths <- depths$topo

library(ggpolypath)

ggplot(fortify(aes_region, region = "Name") %>% inner_join(aes_region@data, c("id" = "Name"))) + aes(x = long, y = lat, group = group, fill = depths) + geom_polypath()

Day lengths

library(raadtools)
light <- readderivaadc("light_budget")
#     units: W.day/m^2"
#[1] "            long_name: annual solar radiation incident on the water surface, adjusted for sea ice cover"

aes_region$light_budget <- extract(light, aes_region_ll, fun = sum, na.rm = TRUE)

ggplot(fortify(aes_region, region = "Name") %>% inner_join(aes_region@data, c("id" = "Name"))) + aes(x = long, y = lat, group = group, fill = light_budget) + geom_polypath()


AustralianAntarcticDivision/aceecostats documentation built on May 5, 2019, 8:14 a.m.