Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

v.stream.network: Fix tostream name, fix for v8 #936

Merged
merged 1 commit into from
Jul 29, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
16 changes: 15 additions & 1 deletion src/vector/v.stream.network/v.stream.network.html
Original file line number Diff line number Diff line change
Expand Up @@ -7,6 +7,10 @@ <h2>DESCRIPTION</h2>
<li><em>tostream</em>: category number of the next stream; is equal to 0 if the stream flows off of the map</li>
</ul>

Any stream-like network will work, but the lines need to be connected.
In terms of graph theory, a tree and forest are supported. Behavior for
a cyclic graph is undefined.

<h2>NOTES</h2>

<b>streams</b> is a set of vector lines that is generated by r.stream.extract. It is recommended to be built as follows (Python code):
Expand All @@ -33,7 +37,17 @@ <h2>NOTES</h2>

<h2>REFERENCES</h2>

None (yet)
<ul>
<li>
Ng, G-H. Crystal, Andrew D. Wickert, Lauren D. Somers, Leila Saberi,
Collin Cronkite-Ratcliff, Richard G. Niswonger,
and Jeffrey M. McKenzie.
"GSFLOW–GRASS v1. 0.0: GIS-enabled hydrologic modeling of coupled
groundwater–surface-water systems."
<em>Geoscientific Model Development</em> 11 (2018): 4755-4777.
<a href="https://doi.org/10.5194/gmd-11-4755-2018">DOI 10.5194/gmd-11-4755-2018</a>
</li>
</ul>

<h2>SEE ALSO</h2>

Expand Down
164 changes: 47 additions & 117 deletions src/vector/v.stream.network/v.stream.network.py
Original file line number Diff line number Diff line change
Expand Up @@ -4,12 +4,13 @@
# MODULE: v.stream.network
#
# AUTHOR(S): Andrew Wickert
# Vaclav Petras (v8 fixes and interface improvements)
#
# PURPOSE: Attach IDs of upstream and downstream nodes as well as the
# category value of the next downstream stream segment
# (0 if the stream exits the map)
#
# COPYRIGHT: (c) 2016-2017 Andrew Wickert
# COPYRIGHT: (c) 2016-2023 Andrew Wickert and the GRASS Development Team
#
# This program is free software under the GNU General Public
# License (>=v2). Read the file COPYING that comes with GRASS
Expand All @@ -20,9 +21,6 @@
# REQUIREMENTS:
# - uses inputs from r.stream.extract

# More information
# Started 14 October 2016

# %module
# % description: Build a linked stream network: each link knows its downstream link
# % keyword: vector
Expand Down Expand Up @@ -73,150 +71,82 @@
# %option
# % key: tostream_cat_column
# % type: string
# % description: Adjacent downstream segment cat (0 = offmap flow)
# % label: Adjacent downstream segment category
# % description: Zero (0) indicates off-map flow
# % answer: tostream
# % required : no
# %end

##################
# IMPORT MODULES #
##################
# PYTHON
"""Add stream network links to attribute table"""

import numpy as np

# GRASS
from grass.pygrass.modules.shortcuts import general as g
from grass.pygrass.modules.shortcuts import raster as r
from grass import script as gs

from grass.pygrass.modules.shortcuts import vector as v
from grass.pygrass.gis import region
from grass.pygrass import vector # Change to "v"?
from grass.pygrass.vector import VectorTopo
from grass.script import vector_db_select
from grass.pygrass.vector import Vector, VectorTopo
from grass.pygrass.raster import RasterRow
from grass.pygrass import utils
from grass import script as gscript

###############
# MAIN MODULE #
###############


def main():
"""
"""Links each segment to its downstream segment

Links each river segment to the next downstream segment in a tributary
network by referencing its category (cat) number in a new column. "0"
means that the river exits the map.
"""

options, flags = gscript.parser()
options, unused_flags = gs.parser()
streams = options["map"]
x1 = options["upstream_easting_column"]
y1 = options["upstream_northing_column"]
x2 = options["downstream_easting_column"]
y2 = options["downstream_northing_column"]
to_stream = options["tostream_cat_column"]

streamsTopo = VectorTopo(streams)
# streamsTopo.build()

# 1. Get vectorTopo
streamsTopo.open(mode="rw")
"""
points_in_streams = []
cat_of_line_segment = []

# 2. Get coordinates
for row in streamsTopo:
cat_of_line_segment.append(row.cat)
if type(row) == vector.geometry.Line:
points_in_streams.append(row)
"""

# 3. Coordinates of points: 1 = start, 2 = end
try:
streamsTopo.table.columns.add(x1, "double precision")
except:
pass
try:
streamsTopo.table.columns.add(y1, "double precision")
except:
pass
try:
streamsTopo.table.columns.add(x2, "double precision")
except:
pass
try:
streamsTopo.table.columns.add(y2, "double precision")
except:
pass
try:
streamsTopo.table.columns.add("tostream", "int")
except:
pass
streamsTopo.table.conn.commit()

# Is this faster than v.to.db?
"""
cur = streamsTopo.table.conn.cursor()
for i in range(len(points_in_streams)):
cur.execute("update streams set x1="+str(points_in_streams[i][0].x)+" where cat="+str(cat_of_line_segment[i]))
cur.execute("update streams set y1="+str(points_in_streams[i][0].y)+" where cat="+str(cat_of_line_segment[i]))
cur.execute("update streams set x2="+str(points_in_streams[i][-1].x)+" where cat="+str(cat_of_line_segment[i]))
cur.execute("update streams set y2="+str(points_in_streams[i][-1].y)+" where cat="+str(cat_of_line_segment[i]))
streamsTopo.table.conn.commit()
streamsTopo.build()
"""
# v.to.db Works more consistently, at least
streamsTopo.close()
# Get coordinates of points: 1 = start, 2 = end
v.to_db(map=streams, option="start", columns=x1 + "," + y1)
v.to_db(map=streams, option="end", columns=x2 + "," + y2)

# 4. Read in and save the start and end coordinate points
colNames = np.array(vector_db_select(streams)["columns"])
colValues = np.array(vector_db_select(streams)["values"].values())
cats = colValues[:, colNames == "cat"].astype(int).squeeze() # river number
xy1 = colValues[:, (colNames == "x1") + (colNames == "y1")].astype(
float
) # upstream
xy2 = colValues[:, (colNames == "x2") + (colNames == "y2")].astype(
float
) # downstream

# 5. Build river network
tocat = []
# Read in and save the start and end coordinate points
col_names = np.array(vector_db_select(streams)["columns"])
col_values = np.array(list(vector_db_select(streams)["values"].values()))
# river number
cats = col_values[:, col_names == "cat"].astype(int).squeeze()
# upstream
xy1 = col_values[:, (col_names == "x1") + (col_names == "y1")].astype(float)
# downstream
xy2 = col_values[:, (col_names == "x2") + (col_names == "y2")].astype(float)

# Build river network
to_cats = []
for i in range(len(cats)):
tosegment_mask = np.prod(xy1 == xy2[i], axis=1)
if np.sum(tosegment_mask) == 0:
tocat.append(0)
to_segment_mask = np.prod(xy1 == xy2[i], axis=1)
if np.sum(to_segment_mask) == 0:
to_cats.append(0)
else:
tocat.append(cats[tosegment_mask.nonzero()[0][0]])
tocat = np.asarray(tocat).astype(int)
to_cats.append(cats[to_segment_mask.nonzero()[0][0]])
to_cats = np.asarray(to_cats).astype(int)

streams_vector = VectorTopo(streams)
streams_vector.open("rw")
streams_vector.table.columns.add(f"{to_stream}", "int")
streams_vector.table.conn.commit()
# This gives us a set of downstream-facing adjacencies.
# We will update the database with it.
streamsTopo.build()
streamsTopo.open("rw")
cur = streamsTopo.table.conn.cursor()
cur = streams_vector.table.conn.cursor()
# Default to 0 if no stream flows to it
cur.execute("update " + streams + " set tostream=0")
for i in range(len(tocat)):
cur.execute(
"update "
+ streams
+ " set tostream="
+ str(tocat[i])
+ " where cat="
+ str(cats[i])
)
streamsTopo.table.conn.commit()
# streamsTopo.build()
streamsTopo.close()

gscript.message("")
gscript.message(
'Drainage topology built. Check "tostream" column for the downstream cat.'
cur.execute(f"update {streams} set {to_stream}=0")
for to_cat, where_cat in zip(to_cats, cats):
cur.execute(f"update {streams} set {to_stream}={to_cat} where cat={where_cat}")
streams_vector.table.conn.commit()
streams_vector.close()

gs.message(
_(
'Drainage topology built. See "{to_stream}" column for the downstream cat.'
).format(to_stream=to_stream)
)
gscript.message("A cat value of 0 indicates the downstream-most segment.")
gscript.message("")
gs.message(_("A cat value of 0 indicates the downstream-most segment."))


if __name__ == "__main__":
Expand Down
Loading