afsc-gap-products / akgfmaps

Make AFSC bottom trawl survey maps and retrieve map layers
Other
17 stars 6 forks source link

spatial boundaries for each area_id in the gap_products.area table #112

Open MattCallahan-NOAA opened 4 months ago

MattCallahan-NOAA commented 4 months ago

Data product requested: Would it be possible to add the spatial boundaries for each area_id in the gap_products.area table as built in data? This would probably be a big lift but area_id is the denominator for many of the other gap_product tables and I think it could really help clarify what users are getting when they query a particular area. Feel free to punt this as a gap_products issue if it fits better there.

Output type: Spatial data

Research team making the request: AKFIN with stock assessors in mind as users.

sean-rohan-NOAA commented 4 months ago

Hi Matt,

This would need to be a gap_products request because it falls outside of the scope of the akgfmaps package.

Before transferring the issue to the gap_products repo, can you clarify what you mean by having spatial boundaries for each area_id? Does that mean adding (for example) a well-known text/well-known binary representation of polygon points?

Currently, our program does not treat shapefiles as 'data products,' which seems like a legacy from when OFIS used to produce and maintain survey shapefiles. However, I feel there is a strong case treating our shapefiles as data products and we've had internal discussions about how we would include spatial information in Oracle tables.

Sean

MattCallahan-NOAA commented 4 months ago

I think the most useful product (for me at least) would be to have a multipolygon with each area_id as a polygon. It could be a filegeodatabase, SF object as a .rda, or some other format that one could read into R and visualize using ggplot/sf.

My understanding is that AKGFMAPS has stratum level polygons, and the area_ids represent either single or grouped strata so it should be possible to derive area_id boundaries from the stratum boundaries.

sean-rohan-NOAA commented 4 months ago

Currently, the akgfmaps package contains stratum shapefiles for each region and you could write the outputs of get_base_layers() to an .rda file. I recently added a function that writes all of the properly masked shapefiles to a .zip file along with a text file with metadata (since polygons are likely to change over time).

With that in mind, would it solve the problem if the shapefiles in akgfmaps simply had AREA_ID fields in addition to STRATUM? AREA_ID is a new field that was created for GAP_PRODUCTS so I haven't encountered use cases for it in my own mapping (yet). However, if there is an emergent need to include AREA_ID in the shapefiles, I could definitely add it.

On Wed, Feb 21, 2024 at 11:59 AM MattCallahan-NOAA @.***> wrote:

I think the most useful product (for me at least) would be to have a multipolygon with each area_id as a polygon. It could be a filegeodatabase, SF object as a .rda, or some other format that one could read into R and visualize using ggplot/sf.

My understanding is that AKGFMAPS has stratum level polygons, and the area_ids represent either single or grouped strata so it should be possible to derive area_id boundaries from the stratum boundaries.

— Reply to this email directly, view it on GitHub https://github.com/afsc-gap-products/akgfmaps/issues/112#issuecomment-1957802891, or unsubscribe https://github.com/notifications/unsubscribe-auth/APULYIMFXVJC3JYYJBARIVDYUZG3BAVCNFSM6AAAAABDTPGNWKVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMYTSNJXHAYDEOBZGE . You are receiving this because you commented.Message ID: @.***>

-- Sean K. Rohan, PhD (he/him) Research Fish Biologist Bering Sea Bottom Trawl Survey Group Resource Assessment and Conservation Engineering Division NOAA Alaska Fisheries Science Center 7600 Sand Point Way NE Seattle, WA 98115 @.***

MattCallahan-NOAA commented 4 months ago

I don't think that adding an area_id column would quite solve it since each stratum will be a component of multiple area_ids. I played around with this a little bit and I think I have a solution. @EmilyMarkowitz-NOAA @zoyafuso-NOAA Does this look correct for determining area_id boundaries?

# devtools::install_github("afsc-gap-products/akgfmaps", build_vignettes = TRUE)
# devtools::install_github("MattCallahan-NOAA/akfingapdata")
# devtools::install_github("MattCallahan-NOAA/akmarineareas2")
library(akgfmaps)
library(tidyverse)
library(sf)
library(akfingapdata)
library(httr)
library(jsonlite)
library(akmarineareas2)

# load akfin_areas and akfin_stratum_groups tables
# This uses the akfingapdata package but you could also use a SQL query
# get_gap_area() pulls "select * from gap_products.akfin_area"
token <- akfingapdata::create_token("web services/Callahan_token.txt")

area <- akfingapdata::get_gap_area()
stratum_groups <- akfingapdata::get_gap_stratum_groups()

# load stratum polygons from akgfmaps
bs_strata <- akgfmaps::get_base_layers(select.region="bs.all", set.crs = "auto")$survey.strata %>%
  rename_with(tolower)
ai_strata <- akgfmaps::get_base_layers(select.region="ai", set.crs = "auto")$survey.strata %>%
  rename_with(tolower) 
goa_strata <- akgfmaps::get_base_layers(select.region="goa", set.crs = "auto")$survey.strata %>%
  rename_with(tolower) 
ebsslope_strata <- akgfmaps::get_base_layers(select.region="ebs.slope", set.crs = "auto")$survey.strata %>%
  rename_with(tolower) 

# plot survey strata
ggplot()+
  geom_sf(data=bs_strata)+
  geom_sf(data=ai_strata, color="red")+
  geom_sf(data=goa_strata, color="blue")+
  geom_sf(data=ebsslope_strata, color="darkgreen")

# Generate area_id polygons by region
# Bering Sea
bs_area <- area %>%
  filter(survey_definition_id %in% c(98, 143) & design_year == 2022)

bs_stratum_groups <- stratum_groups %>%
  filter(survey_definition_id %in% c(98, 143) & design_year == 2022) %>%
  bind_rows(bs_area %>%
  filter(area_type == "STRATUM") %>%
    mutate(stratum = area_id) %>%
  dplyr::select(area_id, survey_definition_id, design_year, stratum, akfin_load_date))

bs_stratum_area <- bs_strata %>%
  mutate(stratum=as.integer(stratum)) %>%
  right_join(bs_stratum_groups, by="stratum")

bs_stratum_area2 <-bs_stratum_area %>%
  group_by(area_id, survey_definition_id) %>%
  summarize()

ggplot(data=bs_stratum_area2)+
  geom_sf(data=ak)+
  geom_sf( fill="red")+
  facet_wrap(~area_id)

# repeat for Other regions
sean-rohan-NOAA commented 4 months ago

Cross posted with gap_products issue 27