Scripts for generating basins

actions/make-basins/DESCRIPTION.md

@@ -0,0 +1,37 @@
+# Make Watershed Basins
+
+## Purpose
+Watersheds in British Columbia can be organized into "basins" - a basin being either an entire river, or "Coast" as sort of a catchall for small rivers that flow to the coast without joining a larger water system. These scripts generate shapefiles that describe basins, given a shapefile of all the watersheds and a list of which watershed is in which basin.
+
+## virtual environment
+
+You will need to set up a virtual environment with python GDAL and shapely installed.
+
+## concave_hull.py
+
+Run this script first. It converts all multipolygon watersheds (coastal areas with lots of islands) into single-polygon watersheds, leaving polygonal watersheds and all watershed attributes unchanged. If you don't do this, the next step will get bogged down in errors about invalid geometries, which are caused by the features with literally thousands of islands.
+
+### Selecting a ratio
+The `ratio` argument is passed directly to shapely's `concave_hull` function. Ratio must be between 0 and 1; higher ratios result in fewer points in the generated shape.
+
+You might think that more points -  a higher resolution geometry -  is ideal, but given the way the concave hull algorithm works, this is not always the case. You might imagine concave hull as attempting to shrinkwrap all the individual points with one sheet of shrinkwrap with the smallest possible amount of air inside the shrink wrap. Lower ratios corresponding to giving the algorithm longer sheets of shrinkwrap. 
+
+Imagine you have three points, and concave hull is attempting to shrinkwrap them. You give the algorithm a short piece of shrinkwrap, and it wraps around the outside of the three points, making a triangle. If, instead, you give the algorithm a longer piece of shrink wrap, it might wrap around the outside from point one to point two and from point two to point three, then turn and wrap arond the inside from point three back to point two, and along the inside back to point one. In this case, rather than a triangle, you get the hollow outline of one, a narrow wrapped strip just big enough to contain the three points.
+
+A ratio of .8 was good for the shapefiles included in this directory, but different ratios might be good for different files. You'll have to open your shapefile in QGIS or something to make sure the generated polygons make sense, and rerun with a different ratio if they don't.
+
+This repository includes an example of a poorly selected ratio. ratio-original.png shows the original shapefile (brown) for Haida Gwai. ratio-08.png shows the resulting geometry (yellow) run with a ratio of .8, which came out well. Ratio-01.png shows the resulting geometry (green) run with a ratio of .1. You can see that the large island at the top has become a thin strip around the outside of the region, instead of enclosing the entire region.
+
+### Selecting a name
+
+This argument is only used in logging. When the script reports issues with a region, this is the shapefile attribute it will refer to the region as so a human can investigate. Pick an attribute that is different for each region, like its name.
+
+## merge-basins.py
+
+This script accepts a shapefile full of watersheds and a CSV list of which watershed belongs to which river basin, and outputs a shapefile full of river basins. In order to build a correspondance between the CSV and the shapefile, one column of the CSV should be named to match an attribute of the shapefile, and the name of this attribute supplied as the "key" argument.
+
+This script removes all holes from the basins, even if individual watersheds contain holes.
+
+## Files Processed
+
+watershed-degrees.shp was processed to concave.shp by concave_hull.py. basins.shp was generated from concave.shp with MjrDrain_Formatted.csv as a list of which watershed is in which basin.

actions/make-basins/concave_hull.py

@@ -0,0 +1,95 @@
+# This script reads a shapefile and outputs a version where all multipolygon
+# features have been replaced by a single polygon representing the concave hull
+# of the original shape. All feature attributes are unchanged.
+
+# This is intended for use with the BC Freshwater Atlas's very detailed maps,
+# where some watersheds have hundreds of little islands.
+
+
+import argparse
+from osgeo import ogr
+import logging, traceback
+from shapely import from_wkt, to_wkt, union, concave_hull
+
+parser = argparse.ArgumentParser(
+    "Merge watershed from a shapefile to make larger basins"
+)
+parser.add_argument("input", help="an input shapefile")
+parser.add_argument("output", help="an output shapefile")
+parser.add_argument("ratio", help="number between 0 and 1, higher = less points")
+parser.add_argument(
+    "name", help="name of the attribute to use in troubleshooting reports"
+)
+
+args = parser.parse_args()
+
+# set up logging
+formatter = logging.Formatter(
+    "%(asctime)s %(levelname)s: %(message)s", "%Y-%m-%d %H:%M:%S"
+)
+handler = logging.StreamHandler()
+handler.setFormatter(formatter)
+
+logger = logging.getLogger(__name__)
+logger.addHandler(handler)
+logger.setLevel(logging.DEBUG)
+
+try:
+    # load input shapefile
+    logger.info("Loading input files")
+
+    driver = ogr.GetDriverByName("ESRI Shapefile")
+    input = driver.Open(args.input)
+
+    if input is None:
+        raise Exception("Could not open shapefile {}".format(args.input))
+
+    in_layer = input.GetLayer(0)  # shapefiles have only one layer
+    in_layer_definition = in_layer.GetLayerDefn()
+
+    # create output files
+    logger.info("Creating output files")
+    output = driver.CreateDataSource(args.output)
+    out_layer = output.CreateLayer(in_layer.GetName())
+
+    # copy all the fields present in the old file to the new one.
+    for i in range(0, in_layer_definition.GetFieldCount()):
+        out_layer.CreateField(in_layer_definition.GetFieldDefn(i))
+
+    hulled = 0
+    # copy each feature to the new file.
+    for in_feature in in_layer:
+        out_feature = ogr.Feature(out_layer.GetLayerDefn())
+
+        for i in range(0, out_layer.GetLayerDefn().GetFieldCount()):
+            out_field = out_layer.GetLayerDefn().GetFieldDefn(i)
+            field_name = out_field.GetName()
+            if field_name == args.name:
+                logger.info("Copying {}".format(in_feature.GetField(i)))
+
+            out_feature.SetField(
+                out_layer.GetLayerDefn().GetFieldDefn(i).GetNameRef(),
+                in_feature.GetField(i),
+            )
+
+        in_geom = in_feature.GetGeometryRef()
+        geom = in_geom.ExportToWkt()
+        if "MULTIPOLYGON" in geom:
+            # convert geometry to shapely, calculate concave hull, and convert it back
+            logger.info("Calculating concave hull")
+            con_hull = to_wkt(concave_hull(from_wkt(geom), args.ratio))
+            out_feature.SetGeometry(ogr.CreateGeometryFromWkt(con_hull))
+            hulled = hulled + 1
+        else:
+            logger.info("Feature is not a multipolygon, left unchanged")
+            out_feature.SetGeometry(in_geom.Clone())
+
+        out_layer.CreateFeature(out_feature)
+        out_feature = None
+
+
+except Exception as e:
+    logger.error(traceback.format_exc())
+
+finally:
+    logger.info("Finished - calculated concave hulls of {} features".format(hulled))

actions/make-basins/inputs/MjrDrain_Formatted.csv

@@ -0,0 +1,247 @@
+Basin,WTRSHDGRPC
+Coast,ALBN
+Coast,ATNA
+Coast,BELA
+Coast,BRKS
+Coast,CAMB
+Coast,CLAY
+Coast,COMX
+Coast,COWN
+Coast,GOLD
+Coast,GRAI
+Coast,HOLB
+Coast,HOMA
+Coast,JERV
+Coast,KEEC
+Coast,KHTZ
+Coast,KITL
+Coast,KITR
+Coast,KLIN
+Coast,KNIG
+Coast,KSHR
+Coast,KTSU
+Coast,KUMR
+Coast,LDEN
+Coast,LRDO
+Coast,MBNK
+Coast,MORI
+Coast,NASC
+Coast,NBNK
+Coast,NECL
+Coast,NEVI
+Coast,NIEL
+Coast,NIMP
+Coast,OWIK
+Coast,PARK
+Coast,PORI
+Coast,SALM
+Coast,SANJ
+Coast,SEYM
+Coast,SKGT
+Coast,SQAM
+Coast,TAHS
+Coast,TATR
+Coast,TOBA
+Coast,TSAY
+Coast,TSIT
+Coast,UDEN
+Coast,UNUR
+Coast,VICT
+Coast,WORC
+Columbia River,BULL
+Columbia River,CANO
+Columbia River,CLRH
+Columbia River,COLR
+Columbia River,DUNC
+Columbia River,ELKR
+Columbia River,KETL
+Columbia River,KHOR
+Columbia River,KOTL
+Columbia River,KOTR
+Columbia River,LARL
+Columbia River,OKAN
+Columbia River,REVL
+Columbia River,SIML
+Columbia River,SLOC
+Columbia River,SMAR
+Columbia River,UARL
+Fraser River,ADMS
+Fraser River,BBAR
+Fraser River,BIGC
+Fraser River,BLAR
+Fraser River,BONP
+Fraser River,BOWR
+Fraser River,BRID
+Fraser River,CARR
+Fraser River,CHES
+Fraser River,CHIL
+Fraser River,CHIR
+Fraser River,CHWK
+Fraser River,CLWR
+Fraser River,COTR
+Fraser River,DEAD
+Fraser River,DOGC
+Fraser River,DRIR
+Fraser River,EUCH
+Fraser River,EUCL
+Fraser River,FRAN
+Fraser River,FRCN
+Fraser River,GRNL
+Fraser River,GUIC
+Fraser River,HARR
+Fraser River,HERR
+Fraser River,HORS
+Fraser River,LCHL
+Fraser River,LCHR
+Fraser River,LEUT
+Fraser River,LFRA
+Fraser River,LILL
+Fraser River,LNIC
+Fraser River,LNRS
+Fraser River,LNTH
+Fraser River,LSAL
+Fraser River,LTRE
+Fraser River,MAHD
+Fraser River,MCGR
+Fraser River,MFRA
+Fraser River,MIDR
+Fraser River,MORK
+Fraser River,MURT
+Fraser River,MUSK
+Fraser River,NARC
+Fraser River,NAZR
+Fraser River,NECR
+Fraser River,NICL
+Fraser River,QUES
+Fraser River,SAJR
+Fraser River,SALR
+Fraser River,SETN
+Fraser River,SHUL
+Fraser River,STHM
+Fraser River,STUL
+Fraser River,STUR
+Fraser River,TABR
+Fraser River,TAKL
+Fraser River,TASR
+Fraser River,THOM
+Fraser River,TWAC
+Fraser River,UCHR
+Fraser River,UEUT
+Fraser River,UFRA
+Fraser River,UNRS
+Fraser River,UNTH
+Fraser River,USHU
+Fraser River,UTRE
+Fraser River,WILL
+Mackenzie River,BEAV
+Mackenzie River,BLUR
+Mackenzie River,CARP
+Mackenzie River,COAL
+Mackenzie River,CRKD
+Mackenzie River,CRYL
+Mackenzie River,DEAL
+Mackenzie River,DEAR
+Mackenzie River,DUNE
+Mackenzie River,FINA
+Mackenzie River,FINL
+Mackenzie River,FIRE
+Mackenzie River,FONT
+Mackenzie River,FOXR
+Mackenzie River,FROG
+Mackenzie River,GATA
+Mackenzie River,HAYR
+Mackenzie River,INGR
+Mackenzie River,KAHN
+Mackenzie River,KCHL
+Mackenzie River,KISK
+Mackenzie River,LBTN
+Mackenzie River,LFRT
+Mackenzie River,LHAF
+Mackenzie River,LIAR
+Mackenzie River,LKEC
+Mackenzie River,LMUS
+Mackenzie River,LOMI
+Mackenzie River,LPCE
+Mackenzie River,LPET
+Mackenzie River,LPRO
+Mackenzie River,LRAN
+Mackenzie River,LSIK
+Mackenzie River,MDEA
+Mackenzie River,MESI
+Mackenzie River,MFRT
+Mackenzie River,MILL
+Mackenzie River,MMUS
+Mackenzie River,MPRO
+Mackenzie River,MURR
+Mackenzie River,NATR
+Mackenzie River,OSPK
+Mackenzie River,PARA
+Mackenzie River,PARS
+Mackenzie River,PCEA
+Mackenzie River,PINE
+Mackenzie River,SAHD
+Mackenzie River,SAHT
+Mackenzie River,SHEK
+Mackenzie River,SMOK
+Mackenzie River,TOAD
+Mackenzie River,TOOD
+Mackenzie River,TSEA
+Mackenzie River,TURN
+Mackenzie River,UBTN
+Mackenzie River,UFRT
+Mackenzie River,UHAF
+Mackenzie River,UKEC
+Mackenzie River,ULRD
+Mackenzie River,UMUS
+Mackenzie River,UOMI
+Mackenzie River,UPCE
+Mackenzie River,UPET
+Mackenzie River,UPRO
+Mackenzie River,USIK
+Nass River,KINR
+Nass River,LBIR
+Nass River,LNAR
+Nass River,NASR
+Nass River,TAYR
+Nass River,UBIR
+Nass River,UNAR
+Skeena River,BABL
+Skeena River,BABR
+Skeena River,BULK
+Skeena River,KISP
+Skeena River,KLUM
+Skeena River,LKEL
+Skeena River,LSKE
+Skeena River,MORR
+Skeena River,MSKE
+Skeena River,SUST
+Skeena River,USKE
+Skeena River,ZYMO
+Stikine River,BARR
+Stikine River,CHUK
+Stikine River,ISKR
+Stikine River,KAKC
+Stikine River,KLAR
+Stikine River,LISR
+Stikine River,LSTR
+Stikine River,MESC
+Stikine River,MSTR
+Stikine River,PITR
+Stikine River,SPAT
+Stikine River,STIR
+Stikine River,TAHR
+Stikine River,TUYR
+Stikine River,UISR
+Stikine River,USTK
+Taku River,INKR
+Taku River,NAHR
+Taku River,NAKR
+Taku River,SHER
+Yukon River,ATLL
+Yukon River,GLAR
+Yukon River,JENR
+Yukon River,KUSR
+Yukon River,SWIR
+Yukon River,TESR
+Yukon River,TUTR
+Yukon River,UJER