Creating valleys with DEM, centerline, and bedrock pts

This weekend I created Python classes for a Digital Elevation Model (DEM) file, for a Valley consisting of 3 files for elevation, centerline, and bedrock points, and for a plot to display these. In addition, Valleys can contain other geometric elements such as local minimal points and lines perpendicular to the centerline. It is expected that this will be the primary class for this project, and will shortly contain the bedrock surface model, and will calculate volumes between the elevation and bedrock data and other boundary surfaces.

I have tested both the Skate Creek and Lower Knowles Creek datasets with these new classes. Here is a plot of the LKC data:

(Click for full-size image)

Here it is in 3d:

(Click for full-size image)

To create a Valley montage, 3 files for elevation, centerline, and bedrock points are read as inputs. At present, in addition to these, local minimal points and vectors perpendicular to the centerline are also calculated. Specifically, the following steps are performed:

1. The elevation data file is read.
   
   
2. The elevation array is smoothed by a 5x5 m filter.
   
   
3. The centerline curve is read as a raster image. These points are then transformed to the coordinates of the elevation dataset.
   
   
4. The centerline is converted to vector format (parameterized by path length), is smoothed with an 11 m filter, and then is resampled and interpolated on circular arcs at 10 m intervals.
   
   
5. The bedrock points are read from another raster file and are transformed to the coordinates of the elevation dataset.
   
   
6. Local minima in the elevation array are calculated.
   
   
7. Vectors perpendicular to the centerline curve are calculated every 20 m. These vector are based on the circular arcs used to interpolate the centerline.
   
   
8. The elevation data, centerline, bedrock pts, minimal pts, and perpendicular lines are plotted in the same window.

Because the overlayed items are in color, it is usually easier to see them if the elevation data is plotted in black & white. Here is the result for the SC dataset:

(Click for full-size image)

In these plots:

• The centerline is red.
• The bedrock pts are large light blue crosses.
• The minimal pts are small yellow or white crosses.
• The perpendicular lines are green, and have been extended 100 m on either side of the centerline.

Note that:

• Neither the bedrock, nor the minimal, points are all on the centerline, nor does it pass directly through the middle of their clusters.
• Minimal pts also occur in the beds of other tributaries, at higher elevations.
• Perpendicular lines often cross one another, especially in regions where the centerline is curved.

I have plotted the minimal pts primarily because Stephen has suggested that they could be used as alternate bedrock pts. They might also allow the software to automatically construct a centerline which passes evenly through them. Note that smoothing the elevation data is necessary to eliminate many spurious minimal pts not in the valley bottom. Here is the result without smoothing the elevation array:

(Click for full-size image)

Here is the Valley for the LKC dataset:

(Click for full-size image)

(Click for full-size image)

One question that immediately arises is: exactly what is your centerline used for, and how did you get it? It might be better to have left and right valley boundaries ('walls', etc...), which are probably not symmetric and which significantly alter the valley width along its length. Once these two curves are defined, a more accurate centerline between them could be calculated (for an example, see the MNRR blog for May). The centerline your've drawn, especialy for the LKC dataset, is not very good. It is clearly too far to the left (south), and doesn't 'wiggle' enough to accurately follow the valley floor.

However, I don't think you really need a centerline, since the volume you want to calculate is enclosed between the elevation surface, the bedrock (yet-to-be-determined), and the left and right valley boundaries. The limits of the valley floor might be determined by some combination of elevation, distribution of local minima, and the gradient of the elevation (or, digitized by hand using the same program I wrote for the MNRR project). For example, here are plots of the magnitude of the gradient of the (unsmoothed) elevation datasets for SC and LKC:

(Click for full-size image)

(Click for full-size image)

Once the valley bounds are determined, the bedrock might be modeled as a kind of 'smooth curving ramp' underneath the surface, either from the bedrock pts you have, plus the distribution of other local minima, or by some other means (e.g. extrapolating the slopes of the mts until they meet). You might also take into consideration the surface shape of the material accumulated in the dihedral. For example, the surface of the SC valley is clearly not level (maybe it's a catenary?):

Finally, there is the question of the shape of the vertical cross-sections: whether they are generally triangular, or trapezoidal (with curved top sides), or something else. All of these considerations will affect the volume estimates.

I think I have a pretty good idea of what you're looking for now. Working with the LKC dataset has helped considerably. There are still some unanswered questions, and certainly more than one way to go about things. However, as I go along the rest of this month I'll try to keep these things in mind, and allow room for making changes to the methodology as we discover new things about the data.