How to find the watershed for each outlet point and sum values within it?

Question

I have the following files to work with:

1. Stream network of a basin;
2. Outlet points for each segment* of the network and
3. Points which represent water withdrawals along of the stream, which are not located in the stream, but near of segments.

*segment = a line with a start point and an outlet point.

My goal is to calculate the sum of withdrawals located upstream of each outlet points.

My first attempt was using model builder, iterating over each outlet point. For each of them, I delineated the watershed and calculated the sum using Zonal Statistics (summing withdrawals contained by the watershed). But it didn't work well because ArcGIS aren't calculating watershed properly. I do not know why...

I'm wondering if I can figure out how to select directly upstream watersheds each iteration, assuming I've already obtained a feature containing the watersheds for each segment.

NOTICE: I don't want to find the watershed for each segment. I want to find the watershed for each outlet point. That's more comprehensive.

I really need a solution..=( .... I'm exhausted...

Thanks.

Answer 1

I did! =)

Watershed command are working right now. I regenerate outlet points and I could get the correct pour points.

In other hand, I followed another method. Instead using Zonal Statistics, I did iteratively for each outlet point as follows:

1. Generate watershed [watershed];
2. Select withdrawals within the watershed [select by location];
3. If the count of selection is greater than zero copy selected [get count and copy features]; Sum field values of the feature selected obtaining as result a table called "sum.table" [summary statistics];
4. Make a layer of the outlet file [make raster layer];
5. Get field value of outlet raster layer [get raster properties];
6. Add field in "sum.table" which name is the same of the outlet raster layer, usually "value" [add field]
7. Insert the field value of outlet raster layer into the new field added in "sum.table" [calculate field];
8. Join the outlet raster layer table with "sum.table" based on the new field added in "sum.table" [add join];
9. Save outlet raster layer to file [save to layer file] - maybe this step is unnecessary; Convert layer file into raster file [copy raster];
10. Actually, each outlet file in these steps do not have the amount of withdrawals in the "value" field, but in another [called 'sum'], which was joined from "sum.table". So, I created a new raster by looking up values found in field which contains the sum [lookup];

So, I have a fold containing all the outlet raster points whose values are the amounts of withdrawals upstream. Next step is using [cell statistics] to sum all these rasters, ignoring NO DATA values.

Thank you guys for the support.

Here is the code generated from Model Builder. Some parts are written in Portuguese.

---------------------------------------------------------------------------

Estatisticas de montante.py

Created on: 2012-09-14 22:14:50.00000

(generated by ArcGIS/ModelBuilder)

Usage: Estatisticas de montante > > > >

Description:

--------------------------------------------------------------------------->

Import arcpy module

import arcpy

Check out any necessary licenses

arcpy.CheckOutExtension("spatial")

Load required toolboxes

arcpy.ImportToolbox("Model Functions") arcpy.ImportToolbox("C:/data/SOMA_OUT/1-Werisk.tbx")

Set Geoprocessing environments

arcpy.env.snapRaster = "" arcpy.env.extent = "" arcpy.env.cellSize = ""

Script arguments

Fozes_raster = arcpy.GetParameterAsText(0)

Campo_com_valores_de_0_ao_número_de_linhas__raster_fozes_ = arcpy.GetParameterAsText(1)

Shape_de_outorgas__exige_campo_com_nome_vazao_ = arcpy.GetParameterAsText(2)

Selecione_o_campo_e_o_tipo_de_estatística = arcpy.GetParameterAsText(3)

Direção_do_Escoamento = arcpy.GetParameterAsText(4)

Pasta_de_saída_de_arquivos = arcpy.GetParameterAsText(5)

Local variables:

Estat_n_ = Fozes_raster Watersh_n_ = Estat_n_ Bacia_shp_n_ = Watersh_n_ outorgas_selecionadas = Bacia_shp_n_ Row_Count = outorgas_selecionadas outorgas_na_bacia_n_ = Row_Count Sum_table = outorgas_na_bacia_n_ Sum_table_com_o_campo = Sum_table Sum_table_com_o_campo_e_com_o_valor = Sum_table_com_o_campo Foz_Extraida_com_valor_da_soma_da_vazao = Sum_table_com_o_campo_e_com_o_valor Est_n_lyr_lyr = Foz_Extraida_com_valor_da_soma_da_vazao Copy_n_ras = Est_n_lyr_lyr Est_n = Copy_n_ras Ext_lyr_n_ = Sum_table Valor_do_campo_da_foz_extraida = Ext_lyr_n_ Num_de_linhas = Fozes_raster Value = Num_de_linhas Campo = "VALUE" Field_Name = "%Campo%"

Process: Get Count

arcpy.GetCount_management(Fozes_raster)

Process: For

arcpy.IterateCount_mb("1", Num_de_linhas, "1")

Process: Extract by Attributes

tempEnvironment0 = arcpy.env.scratchWorkspace arcpy.env.scratchWorkspace = "C:\Users\PC\Documents\ArcGIS\Default.gdb" tempEnvironment1 = arcpy.env.workspace arcpy.env.workspace = "C:\Users\PC\Documents\ArcGIS\Default.gdb" arcpy.gp.ExtractByAttributes_sa(Fozes_raster, "\"%Campo com valores de 0 ao número de >linhas (raster fozes)%\" = %n%", Estat_n_) arcpy.env.scratchWorkspace = tempEnvironment0 arcpy.env.workspace = tempEnvironment1

Process: Watershed

arcpy.gp.toolbox = "C:/data/SOMA_OUT/1-Werisk.tbx";

Warning: the toolbox C:/data/SOMA_OUT/1-Werisk.tbx DOES NOT have an alias.

Please assign this toolbox an alias to avoid tool name collisions

And replace arcpy.gp.Watershed(...) with arcpy.Watershed_ALIAS(...)

arcpy.gp.Watershed(Direção_do_Escoamento, Estat_n_, Watersh_n_, Campo)

Process: Raster to Polygon

arcpy.RasterToPolygon_conversion(Watersh_n_, Bacia_shp_n_, "SIMPLIFY", "VALUE")

Process: Select Layer By Location

arcpy.SelectLayerByLocation_management(Shape_de_outorgas__exige_campo_com_nome_vazao_, >"WITHIN", Bacia_shp_n_, "", "NEW_SELECTION")

Process: Get Count (2)

arcpy.GetCount_management(Shape_de_outorgas__exige_campo_com_nome_vazao_)

Process: Copy Features

arcpy.CopyFeatures_management(Shape_de_outorgas__exige_campo_com_nome_vazao_, >outorgas_na_bacia_n_, "", "0", "0", "0")

Process: Summary Statistics

tempEnvironment0 = arcpy.env.scratchWorkspace arcpy.env.scratchWorkspace = "C:\Users\PC\Documents\ArcGIS\Default.gdb" arcpy.Statistics_analysis(outorgas_na_bacia_n_, Sum_table, >Selecione_o_campo_e_o_tipo_de_estatística, "") arcpy.env.scratchWorkspace = tempEnvironment0

Process: Make Raster Layer

tempEnvironment0 = arcpy.env.scratchWorkspace arcpy.env.scratchWorkspace = "C:\Users\PC\Documents\ArcGIS\Default.gdb" arcpy.MakeRasterLayer_management(Estat_n_, Ext_lyr_n_, "", "302545,967240849 >8506968,76613879 725280,246664134 8894612,97717014", "") arcpy.env.scratchWorkspace = tempEnvironment0

Process: Get Raster Properties

arcpy.GetRasterProperties_management(Ext_lyr_n_, "MAXIMUM")

Process: Add Field

arcpy.AddField_management(Sum_table, Field_Name, "LONG", "", "", "", "", "NULLABLE", >"NON_REQUIRED", "")

Process: Calculate Field

arcpy.CalculateField_management(Sum_table_com_o_campo, "%Campo%", "%Valor do campo da foz >extraida%", "VB", "")

Process: Add Join

arcpy.AddJoin_management(Ext_lyr_n_, Campo, Sum_table, Campo, "KEEP_COMMON")

Process: Save To Layer File

arcpy.SaveToLayerFile_management(Ext_lyr_n_, Est_n_lyr_lyr, "")

Process: Copy Raster

arcpy.CopyRaster_management(Est_n_lyr_lyr, Copy_n_ras, "", "", "", "NONE", "NONE", "")

Process: Lookup

arcpy.gp.Lookup_sa(Copy_n_ras, "SUM_VAZAO", Est_n)

Answer 2

Assuming that each withdrawal point is related to the closest stream segment, we can do the following:

1. Use the split line at point tool to split the stream line into segments based on outlet points (depending on your data you may already have this). Make sure that each point has an ID, which should relate to stream order (see one method below for an example in Python).
2. Use the generate near table tool to generate a table that will link your withdrawal points to the nearest stream segment. Check the only get closest checkbox so you only get one result for each point.
3. Join the withdrawal amout to the near table to associate with the segment IDs.
4. Use the summary statistics tool to sum the withdrawal amount grouped by segment ID.
5. Join the output withdrawal sums to your stream segment, select those segments upstream from the point of interest, and sum the total withdrawals (again, see below for one way to do this).

And example of IDing the segments using Python.

import arcpy

#Written assuming arcGIS 10
#Note this will require modification to run...

start_pts = {}
end_pts = {}

for segment in arcpy.SearchCursor("segmented_line"):
    start_pts[segment.SEGMENT_ID] = segment.SHAPE.getObject(0) #SHAPE returns an arcpy.Array
                                                               #object containing points for
                                                               #geom type line
    end_pts[segment.SEGMENT_ID] = segment.SHAPE.getObject(segment.SHAPE.count)

try:
    #cleanup as you can only have one cursor open at once, otherwise wackiness
    del segment
except:
    pass

rows = arcpy.InsertCursor("start_pts")
for start in start_pts:
    row = rows.newRow()
    row.SHAPE = start_pts[start]
    row.SEGMENT_ID = start
    rows.insertRow(row)

try:
    #more cleanup
    del row
    del rows
except:
    pass

arcpy.SpatialJoin_analysis(start_pts, outlet_pts, "start_outlet_pts",
                           "JOIN_ONE_TO_ONE", match_option="WITHIN_A_DISTANCE",
                           search_radius="some metres")

arcpy.JoinField_management("segmented_line", "SEGMENT_ID", "start_outlet_pts", "SEGMENT_ID", "OUTLET_ID")

################################################################################
##Repeat for end_pts
################################################################################

Once I get a bit of time later I'll come back and update this answer with a bit of script to build a network using networkx to help sum the values...

Hope this helps!

Answer 3

The solution to your problem may be to use Infrastructure Network Editing, specifically the Flow Accumulation tool on the Infrastructure Reporting Toolbar.

You'll need at least a Standard (Editor) license to build a geometric network. You'll also need to specify a weight. At each point or line, Calculate Accumulation will add up all your upstream weights and store them in a field.

Google "Infrastructure Network Editing 10.0" or "Infrastructure Network Editing 10.1" to get the right version. This only works with ArcGIS 10.0 and above.