afsc-gap-products / akgfmaps

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

Invalid geometries in built-in shapefiles #105

Open sean-rohan-NOAA opened 7 months ago

sean-rohan-NOAA commented 7 months ago

Issue

Four shapefiles in the package have invalid geometries. Two are survey shapefiles (bs_grid.shp, bs_grid_w_corners.shp), one is the AKRO NMFS reporting area file (NMFS Reporting Areas.shp), and one is a land shapefile (alaska_canada_dcw.shp). The survey grid and land shapefiles should be corrected to avoid errors in spatial operations.

Errors:

Code

library(akgfmaps)

file_paths <- list.files(system.file("extdata", package = "akgfmaps"), full.names = TRUE, pattern = ".shp")

file_paths <- file_paths[!grepl(pattern = ".xml", x = file_paths)]

valid_geom <- logical()
valid_geom_wgs84 <- logical()
# file_crs <- character()

for(ii in 1:length(file_paths)) {

  dat <- sf::st_read(file_paths[ii])

  dat_wgs84 <- sf::st_transform(dat, crs = "WGS84")

  valid_geom <- c(valid_geom,  all(sf::st_is_valid(dat)))
  valid_geom_wgs84 <- c(valid_geom_wgs84,  all(sf::st_is_valid(dat_wgs84)))

}

invalid_paths <- file_paths[!valid_geom]

invalid_info <- vector(mode = "list", length = length(invalid_paths))
invalid_geoms <- vector(mode = "list", length = length(invalid_paths))

for(jj in 1:length(invalid_paths)) {

  dat <- sf::st_read(invalid_paths[jj])

  dat$reason <- sf::st_is_valid(dat, reason = TRUE)

  invalid_info[[jj]] <- dat

  invalid_geoms[[jj]] <- dplyr::filter(dat, reason != "Valid Geometry")

}

basename(file_paths[!valid_geom])
basename(file_paths[!valid_geom_wgs84])

invalid_geoms

invalid_geoms [[1]] Simple feature collection with 2 features and 3 fields Geometry type: MULTIPOLYGON Dimension: XY Bounding box: xmin: -136.4445 ymin: 51.1736 xmax: -61.13481 ymax: 83.11062 Geodetic CRS: NAD83 POPYCOUN POPYADMIN geometry reason 1 MULTIPOLYGON (((-110 55.833... Loop 8: Edge 4 crosses edge 7 2 CA NORTHWEST TERRITORIES MULTIPOLYGON (((-78.87537 5... Loop 3172: Edge 5 has duplicate vertex with edge 8

[[2]] Simple feature collection with 1 feature and 12 fields Geometry type: MULTIPOLYGON Dimension: XY Bounding box: xmin: -298595.7 ymin: 953450.7 xmax: -256172.7 ymax: 991824.3 Projected CRS: NAD83 / Alaska Albers OBJECTID FID_1 NM2 HECTARES NUM_ID STATION_ID STATION LON LAT Shape_Leng Shape_Area geometry 1 645 644 205.9528 70639.84 644 L-15 L15 -158.6515 58.66617 127981.7 706398406 MULTIPOLYGON (((-298457.2 9... reason 1 Ring Self-intersection[-276937.6875 961103.5]

[[3]] Simple feature collection with 1 feature and 2 fields Geometry type: MULTIPOLYGON Dimension: XY Bounding box: xmin: -298595.7 ymin: 953450.7 xmax: -256172.7 ymax: 991824.3 Projected CRS: NAD_1983_Albers STATIONID geometry reason 1 L-15 MULTIPOLYGON (((-261361.7 9... Ring Self-intersection[-276937.6875 961103.5]

[[4]] Simple feature collection with 2 features and 6 fields Geometry type: MULTIPOLYGON Dimension: XY Bounding box: xmin: -15062.46 ymin: 283196.6 xmax: 1360571 ymax: 1306651 Projected CRS: NAD83 / Alaska Albers REP_AREA Area AreaSqM Shape_Leng Shape_Area geometry reason 1 630 275164.9 275164939368 9480495 275164939388 MULTIPOLYGON (((52598.34 79... Ring Self-intersection[97554.5511000007 938299.317299999] 2 650 179263.3 179263268369 5173867 179263268369 MULTIPOLYGON (((999533.1 10... Ring Self-intersection[999533.125 1048954.25]

basename(file_paths[!valid_geom]) [1] "alaska_canada_dcw.shp" "bs_grid.shp" "bs_grid_w_corners.shp" "NMFS Reporting Areas.shp"