How to Flat-Fix a DEM for Flow Analysis

635
1
12-16-2022 06:11 AM
ZachBodenner
MVP Regular Contributor

I'm working on writing a tutorial for some of my users to practice flow modeling. I grabbed the mouth of the St. Louis river in Duluth, thinking it would be a topographically interesting area with lots of rivers and streams. So far so good. I downloaded a 1m DEM from the state's database which looks pretty good. I start to run into problems actually running the flow model though.

The river flows Northeast, but the algorithm breaks down here because the DEM has identical values for a good chunk of the river mouth, basically suggesting it's a flat surface. As a result, the classified raster suggests that the river itself is flowing in multiple directions (dark green = east, light green = north, yellow = west, red=south (interesting that I didn't get any of the other directions for the "flat" area but maybe that's a different issue)). That might not be inherently bad, but the red sections of southerly flow essentially mean that the river never reaches lake superior and begins to flow back up a nearby river.
 
So my question is, how do I apply some kind of artificial sloping of the raster to trick it into thinking that it flows northeast continually? I'm sure there are some tools to help out, but I admit I'm kind of stumped.
0 Kudos
1 Reply
DanPatterson
MVP Esteemed Contributor

Through trickery and numpy/arcpy wizardary.

Now... you need to know the size of your raster in terms of rows and columns and the direction in which you want it to point.

As an example, consider a small area, say 5x5 units (so we can see the results).

# -- the delta e0
e0 =np.linspace(0., 10., num=25).reshape(5, 5)  # -- make some fake elevation data

e0  # -- so far so good, but it points north-west 
array([[ 0.        ,  0.41666667,  0.83333333,  1.25      ,  1.66666667],
       [ 2.08333333,  2.5       ,  2.91666667,  3.33333333,  3.75      ],
       [ 4.16666667,  4.58333333,  5.        ,  5.41666667,  5.83333333],
       [ 6.25      ,  6.66666667,  7.08333333,  7.5       ,  7.91666667],
       [ 8.33333333,  8.75      ,  9.16666667,  9.58333333, 10.        ]])

e1 = np.flip(e0, axis=1)  # -- cleverly flip it so it points north-east
e1
array([[ 1.66666667,  1.25      ,  0.83333333,  0.41666667,  0.        ],
       [ 3.75      ,  3.33333333,  2.91666667,  2.5       ,  2.08333333],
       [ 5.83333333,  5.41666667,  5.        ,  4.58333333,  4.16666667],
       [ 7.91666667,  7.5       ,  7.08333333,  6.66666667,  6.25      ],
       [10.        ,  9.58333333,  9.16666667,  8.75      ,  8.33333333]])

Now "simply" convert this to a raster using 

NumPyArrayToRaster—ArcGIS Pro | Documentation

which gives you the opportunity to define a spatial reference etc (see the help topic's code example.

Add your resultant raster to your DEM which effectively provides a NE underlying slope to all your existing elevation values.

Calculate you stuff.

Now be caution on the elevation difference that you want.  In my example I produced a 10 unit elevation difference in a 5x5 unit raster.  Be realistic or be dramatic if you want to help nature along or cause torrents of water flow


... sort of retired...
0 Kudos