Beyond
Mapping III

Map
Analysis book with companion CDROM
for handson exercises and further reading 
Use
Data to Characterize MicroTerrain Features — describes
techniques to identify convex and concave features
Characterizing
Local Terrain Conditions — discusses
the use of "roving windows" to distinguish localized variations
Characterizing
Terrain Slope and Roughness — discusses
techniques for determining terrain inclination and coarseness
The Long and Short of Slope
— investigates
longitudinal and transverse slope calculation
Generating
Mountains and Molehills from Field Sampled Data — creating
an elevation surface from field sampled data
Segmenting Our World — discusses
techniques for segmenting linear routes based on terrain inflection
Confluence Maps Further Characterize Microterrain
Features —
describes the use of optimal path density analysis for mapping surface flows
Modeling Erosion and
Sediment Loading — illustrates
a
Identify Valley Bottoms in Mountainous Terrain
— illustrates a
technique for identifying flat areas connected to streams
Use Surface Area for Realistic Calculations
— describes a
technique for adjusting planimetric area to surface area considering terrain
slope
Beware of Slope’s Slippery Slope
— describes
various slope calculations and compares results
Shedding Light on Terrain Analysis — discusses how terrain orientation is used to generate Hillshade maps
Identifying Upland Ridges — describes a procedure for locating extended upland ridges
Note: The processing and figures discussed in this topic were derived using MapCalc^{TM}
software. See www.innovativegis.com to download a
free MapCalc Learner version with tutorial materials for classroom and
selflearning map analysis concepts and procedures.
<Click here>
rightclick to download a
printerfriendly version of this topic (.pdf).
(Back to the Table of
Contents)
______________________________
Use Data
to Characterize MicroTerrain Features
(GeoWorld January, 2000; “” pg.
22)
The past several columns investigated surface modeling and
analysis. The data surfaces derived in
these instances weren't the familiar terrain surfaces you walk the dog, bike
and hike on. Nonetheless they form
surfaces that contain all of the recognizable gullies, hummocks, peaks and
depressions we see on most hillsides.
The "wrinkledcarpet" look in the real world is directly
transferable to the cognitive realm of the data world.
Figure 311. Identifying
Convex and Concave features by deviation from the trend of the terrain.
One of the most
frequently used procedures involves comparing the trend of the surface to the
actual elevation values. Figure 311
shows a terrain profile extending across a small gully. The dotted line identifies a smoothed curve
fitted to the data, similar to a draftsman's alignment of a French curve. It "splitsthedifference" in the
succession of elevation values—half above and half below. Locations above the trend line identify
convex features while locations below identify the concave ones. The further above or below determines how pronounced
the feature is.
In a
In generating the smoothed surface, elevation values were averaged for a 4by4
window moved throughout the area. Note
the subtle differences between the surfaces—the tendency to pulldown the
hilltops and pushup the gullies.
While you see (imagine?) these differences in the surfaces, the computer
quantifies them by subtracting. The
difference surface on the right contains values from 84 (prominent concave
feature) to +94 (prominent convex feature).
The big bump on the difference surface corresponds to the smaller
hilltop on elevation surface. Its actual
elevation is 2016 while the smoothed elevation is 1922 resulting in 2016  1922
= +94 difference. In microterrain
terms, these areas are likely drier than their surroundings as water flows
away.
The other arrows on the surface indicate other interesting locations. The "pockmark" in the foreground is
a small depression (764  796 = 32 difference) that is likely wetter as water
flows into it. The "deep cut"
at the opposite end of the difference surface (539  623 = 84) suggests a
prominent concavity. However
representing the water body as fixed elevation isn't a true account of terra
firma configuration and misrepresents the true microterrain pattern.
Figure 312. Example of a microterrain
deviation surface.
In fact the entire
concave feature in the upper left portion of 2D representation of the difference
surface is suspect due to its treatment of the water body as a constant
elevation value. While a fixed value for
water on a topographic map works in traditional mapping contexts it's
insufficient in most analytical applications.
Advanced
The 2D map of differences identifies areas that are concave (dark red), convex
(light blue) and transition (white portion having only 20 to +20 feet
difference between actual and smoothed elevation values). If it were a map of a farmer's field, the
groupings would likely match a lot of the farmer's recollection of crop
production—more water in the concave areas, less in the convex areas.
A
The idea of variable rate response to spatial conditions has been around for
thousands of years as indigenous peoples adjusted the spacing of holes they
poked in the ground to drop in a seed and a piece of fish. While the mechanical and green revolutions
enable farmers to cultivate much larger areas they do so in part by applying
broad generalizations of microterrain and other spatial variables over large
areas. The ability to continuously
adjust management actions to unique spatial conditions on expansive tracks of
land foretells the next revolution.
Investigate the effects of microterrain conditions goes
well beyond the farm. For example, the
Universal Soil Loss Equation uses "average" conditions, such as
stream channel length and slope, dominant soil types and existing land use
classes, to predict water runoff and sediment transport from entire
watersheds. These nonspatial models are
routinely used to determine the feasibility of spatial activities, such as
logging, mining, road building and housing development. While the procedures might be applicable to
typical conditions, they less accurately track unusual conditions clumped
throughout an area and provide no spatial guidance within the boundaries of the
modeled area.
Characterizing Localized
Terrain Conditions
(GeoWorld, February 2000, pg. 26)
Last month's column
described a technique for characterizing microterrain features involving the
difference between the actual elevation values and those on a smoothed
elevation surface (trend). Positive
values on the difference map indicate areas that "bumpup" while
negative values indicate areas that "dipdown" from the general trend
in the data.
A related technique
to identify the bumps and dips of the terrain involves moving a "roving
window" (termed a spatial filter) throughout an elevation surface. The profile of a gully can have
microfeatures that dip below its surroundings (termed concave) as shown on the
right side of figure 313.
Figure 313. Localized deviation uses a spatial filter to
compare a location to its surroundings.
The localized
deviation within a roving window is calculated by subtracting the average of
the surrounding elevations from the center location's elevation. As depicted in the example calculations for
the concave feature, the average elevation of the surroundings is 106, that computes to a 6.00 deviation when subtracted from the
center's value of 100. The negative sign
denotes a concavity while the magnitude of 6 indicates it's fairly significant
dip (a 6/100= .06). The protrusion above
its surroundings (termed a convex feature) shown on the right of the figure has
a localized deviation of +4.25 indicating a somewhat pronounced bump (4.25/114=
.04).
Figure 314. Applying Deviation and Coefficient of
Variation filters to an elevation surface.
The result of moving
a deviation filter throughout an elevation surface is shown in the top right
inset in figure 314. Its effect is
nearly identical to the trend analysis described last month comparison of
each location's actual elevation to the typical elevation (trend) in its
vicinity. Interpretation of a Deviation
Surface is the same as that for the difference surface discussed last
month—protrusions (large positive values) locate drier convex areas;
depressions (large negative values) locate wetter concave areas.
The implication of the "Localized Deviation" approach goes far beyond
simply an alternative procedure for calculating terrain irregularities. The use of "roving windows"
provides a host of new metrics and map surfaces for assessing microterrain
characteristics. For example, consider
the Coefficient of Variation (Coffvar) Surface shown
in the bottomright portion of figure 314.
In this instance, the standard deviation of the window is compared to
its average elevation—small "coffvar"
values indicate areas with minimal differences in elevation; large values
indicate areas with lots of different elevations. The large ridge in the coffvar
surface in the figure occurs along the shoreline of a lake. Note that the ridge is largest for the
steeplyrising terrain with big changes in elevation. The other bumps of surface variability noted
in the figure indicate areas of less terrain variation.
While a statistical summary of elevation values is useful as a general
indicator of surface variation or "roughness," it doesn't consider
the pattern of the differences. A
checkerboard pattern of alternating higher and lower values (very rough) cannot
be distinguished from one in which all of the higher values are in one portion
of the window and lower values in another.
Figure 315. Calculation of slope considers the
arrangement and magnitude of elevation differences within a roving window.
There are several
roving window operations that track the spatial arrangement of the elevation
values as well as aggregate statistics.
A frequently used one is terrain slope that calculates the
"slant" of a surface. In
mathematical terms, slope equals the difference in elevation (termed the
"rise") divided by the horizontal distance (termed the
"run").
As shown in figure 315, there are eight surrounding elevation values in a 3x3
roving window. An individual slope from
the center cell can be calculated for each one.
For example, the percent slope to the north (top of the window) is
((2332  2262) / 328) * 100 = 21.3%. The
numerator computes the rise while the denominator of 328 feet is the distance
between the centers of the two cells.
The calculations for the northeast slope is ((2420  2262) / 464) * 100
= 34.1%, where the run is increased to account for the diagonal distance (328 *
1.414 = 464).
The eight slope values can be used to identify the Maximum, the Minimum and the
Average slope as reported in the figure.
Note that the large difference between the maximum and minimum slope (53
 7 = 46) suggests that the overall slope is fairly variable. Also note that the sign of the slope value
indicates the direction of surface flow—positive slopes indicate flows into the
center cell while negative ones indicate flows out. While the flow into the center cell depends
on the uphill conditions (we'll worry about that in a subsequent column), the
flow away from the cell will take the steepest downhill slope (southwest flow
in the example… you do the math).
In practice, the Average slope can be misleading. It is supposed to indicate the overall slope
within the window but fails to account for the spatial arrangement of the slope
values. An alternative technique fits a
“plane” to the nine individual elevation values. The procedure determines the best fitting
plane by minimizing the deviations from the plane to the elevation values. In the example, the Fitted
slope is 65%… more than the maximum individual slope.
At first this might seem a bit fishy—overall slope more than the maximum
slope—but believe me, determination of fitted slope is a different kettle of
fish than simply scrutinizing the individual slopes. Next time we'll look a bit deeper into this
fitted slope thing and its applications in microterrain analysis.
_______________________
Author's Note: An Excel worksheet
investigating Maximum, Minimum, and Average slope calculations is available
online at the "Column Supplements" page at http://www.innovativegis.com/basis.
(GeoWorld, July 2007, pg. 1819)
Recall from the past sections dealing with slope calculation that the nine elevation values in a 3x3 window are generally used. There are two basic approaches for characterizing relative steepness—surface fitting and individual slope statistics.
Surface fitting aligns a plane to the elevation values that minimizes the deviations from the plane to the values. If all of the values are the same, a horizontal plane is fitted and the slope is zero. However as the values on one side of the window increase and the values on the other decrease, the fitted plane tilts indicating increasing terrain steepness.
The surface fitting algorithm first establishes the 3x3 window, then identifies the nine elevation values, fits the plane and records the slope for the center cell in the window. The “roving window” shifts to the next adjacent location and repeats the process until all of the cells in a project area have been evaluated.
The algorithm for summarizing individual slopes also utilizes a 3x3 window, but instead of fitting a plane it calculates the slopes formed by the center location and each of its eight neighbors. The minimum, maximum or average of the individual slope values is recorded for the center cell and the window advanced.
Both slope techniques result in a continuous map surface indicating the relative steepness throughout a project area. Sometimes, however, slope along a linear feature, such as a pipeline, is sought (figure 1). In this case, only three cells are involved—the center cell and its input and output neighbors at each location along the pipeline.
Figure 1. Longitudinal slope calculates the steepness
along a linear feature, such as a pipeline.
Longitudinal Slope (LS) requires isolating just the route’s elevation values. In gridbased modeling this means first creating a “masking map” of the route where the cells defining the pipeline are assigned the value 1 and the nonpipeline cells are assigned a special value of Null (or No Data). This map is multiplied by the elevation surface having the effect of retaining the elevation values along the route but all other locations are assigned Null.
Most map analysis packages have been “trained” to recognize a Null value and simply skip over these locations when processing. The effect in this case is to calculate slope using just the three aligned elevation values in the 3x3 roving window—either fitted or maximum/minimum/average.
A simple extension to the algorithm, Directional Slope (DS), enables a user to enter a starting location and the direction of flow, and then it calculates individual slopes for each step along the route using just the input and center cell. This technique characterizes where gravity is helping or hindering flow and is useful in determining drawdown if the pipe is ruptured.
Figure 2. Modified transverse slope is calculated as
the average of offroute slopes surrounding a location.
A more radical extension considers the offroute cells in calculating Transverse Slope (TS). Traditionally, transverse slope was manually calculated by estimating the slope of a line perpendicular to the route. For example in the upper portion of figure 2, the slopes from the center cell (5) to the NE (53) and SW (57) cells identify the longitudinal slope along the route while the slopes to the NW (511) and SE (59) cells indicate the transverse slope.
The perpendicular approach works for any orthogonal and diagonal transect of the window. The elevation grid, however, can take an oblique bend involving diagonal and orthogonal directions as shown and there isn’t a true perpendicular. In this instance, GISderived transverse slope requires a new definition.
Modified Transverse Slope (MTS) calculates slope for the offroute cells, regardless of the linear shape. In the lowerleft portion of figure 2 all of green cells (1, 2, 4, 6, 8 and 9) are considered in calculating terrain steepness surrounding the pipeline. Ideally, routes with fairly gentle modified transverse slopes are preferred as subsidence risks are lower.
These fairly innocuous extensions for calculating slope point to a bigger issue in map analysis—mainly that our gridbased analytical toolbox and GIS modeling expressions are a long way from complete. Most of the development occurred in the 1970s and 80s and coding of new capabilities in flagship systems have languished. As a result, solution providers are encoding specialized tools for their clients and running them outside the GIS or hooking them in as objects or addins.
What I find interesting is that this is where GIS was in the 1970s and early 80s—a collection of disparate proprietary systems. With all that you hear about data standards, open systems and GIS community access it’s peculiar that map analysis seems to be moving in the opposite direction.
Generating Mountains
and Molehills from Field Sampled Data
(GeoWorld, August 2013)
Geotechnology has revolutionized the dogma, development and
dissemination of spatial information—duh.
Today’s youth is growing up in an increasingly cyberbased society with
instantaneous connection to the world about them, and like Peter Pan they can
be whiskedaway by Google Earth to distant lands in the blink of an eye.
But what if they are interested in detail beyond the maximum level of
the zoom slider bar? Or they want to
simulate localized surface flow paths, convergence and path density? Or identify subtle hummocks and
depressions? Or model subsurface
moisture regimes? Or relate any of these characteristics to the spatial
relationships with vegetative patterns?
That’s the dilemma of many aspiring scientists who want to go beyond simply visualizing a landscape to investigating system interactions that are spatially driven. But their data analysis quiver is primarily full of nonspatial analytic tools that rely on assessing and comparing the central tendencies of field collected data without consideration of the spatial distributions inherent in the data.
Two major obstacles stand in their way— appropriate spatial data and spatially aware analysis tools. Most readily downloadable data is too course for detailed scientific study, and even if it were available, most potential users are unacquainted with the quantitative analysis procedures involved.
Let’s consider the data limitation hurdle first and then see how its practical solution can spark spatial reasoning among nonGIS’ers. Suppose a budding grad student wants to develop a microterrain surface for a project area. Unless there is a Lidar over flight or total station survey of the area available, they will need to construct a detailed elevation surface from scratch.
First impulse is to use an inexpensive GPS unit (or Google’s GPS Surveyor app for Smartphones) to capture XYZ measurements throughout a project area. But alas, it is soon discovered that precision of the coordinates is insufficient for detailed mapping. The vertical precision in a GPS unit is characteristically less than half of that of the horizontal precision which in itself is plus or minus several meters—hardly the precision required for surface modeling.
Another possibility is to acquire sophisticated surveying equipment (transit, theodolite, laser rangefinder or GPS total station), but their costs are far too high both in terms of a tight budget and a demanding learning curve. A Dumpy level, tripod, graduated staff and a couple of meter tapes form a far more practical option (they are likely buried in an old storeroom on campus).
Dumpy levels work by precisely measuring vertical distances in a horizontal plane about a base location (figure 1). A bubble level is used to establish a horizontal plane in each quartercircle to make sure it is level throughout the entire 360 degree swing of the instrument. An operator looks through the eyepiece of the telescope, while an assistant slowly rocks a graduated staff at a sample location and the maximum reading on the staff is recorded. Subtracting the Staff Reading from the Base Height determines the relative elevation of any measured point to a fraction of centimeter.
Figure 1. A Dumpy level establishes a
horizontal plane to measure the relative elevation differences throughout a
project area.
To establish realworld positioning, an inexpensive GPS unit can be used. The average of several readings taken throughout the day provides a “reasonable estimate” of the Base Location’s earth position (topright portion of figure 2).
Figure 2. The XYZ coordinates of
sample points are easily derived from field data.
Before the staff is moved to another location, the distance along a center line is measured (EW as the X axis in the example) and the perpendicular deflection to the point is measured (NS as the Y axis). The result is a record of relative horizontal location and elevation difference that is easily evaluated for XYZ coordinates within a working coordinate system (bottomright portion of figure 2) or relative to the base station’s Lat/Lon realworld coordinates.
The final step involves surface modeling software to complete the
3dimentional surface by “spatially interpolating” the elevation values for all
of the grid cells within the project area.
Figure 3 shows the result using the powerful yet inexpensive Surfer
software package (by Golden Software).
The derived terrain surface can be exported to ArcGIS, MapCalc or other
map analysis package, as well as to Google Earth as contour lines or a
colorshaded raster image.
Figure 3. A detailed elevation
surface is generated using inexpensive surface modeling software.
I have presented this 3hour workshop in various forms to numerous GIS
and nonGIS students and professionals since the 1980s. The practical fieldwork coupled with
easytouse software encourages dialog about the inherent spatial nature
ingrained in most field collected data.
Discussion quickly extends to data other than elevation—the “Z values”
could easily be relative concentrations of an environmental variable; or
disease occurrences in an epidemiological study; or relative biomass or stem
density measurements of a plant species; or crime incidences for mapping
pockets of crime in a city; or point of sales values for a particular product;
or questionnaire response levels in a social survey; etc.
The peaks and valleys in any map surface characterize the spatial
distribution of that mapped variable— identifying where the relative highs and
lows occur in geographic space. And
since a map surface is a spatially organized set of numbers, mathematical and
statistical analysis tools can be applied (see Author’s Note).
The bottom line is that many (most?) fieldbased mapped data contain
useful information about the spatial distribution/pattern inherent in the
data. While the ability to surface map
this characteristic is wellknown within the GIS community, it is less familiar
and woefully underutilized in the applied sciences. A simple experience in surface mapping
terrain elevation can be useful in crossing this conceptual chasm and often
serves as a Rosetta stone in stimulating crossdiscipline discussion and
recognition of “maps as data” awaiting quantitative analysis— as opposed to
traditional “graphic images” for human viewing and interpretation.
_____________________________
Author’s Note:
for more discussion on math/stat analytics for mapped data, see the online book
Beyond Mapping III, Topic 30 “A
Math/Stat Framework for Map Analysis” posted at www.innovativegis.com/basis/MapAnalysis/.
(GeoWorld, June 2007, pg. 1819)
Dividing geographic space into meaningful groupings is the basis of all mapping. This holds for abstract maps like customer propensity clusters for purchasing a product or animal habitat suitability, as well as more traditional mapping applications such as vegetation maps, ownership parcels, highways and pipeline routes. Further divisions of the patterns in a base map often involve segmentation.
Figure 1. Basic approaches used in route segmentation.
For example, a pipeline is frequently divided into “stationing coordinates” that partitions a serpentining route into uniform linear distances along a pipeline (Uniform segmentation). Alternatively, the route can be segmented by identifying a point at each major change in direction as shown in the top right portion of figure 1. The result divides the route into a set of variable length segments responding to planimetric orientation of straight line sections (Deflection Angle segmentation).
A more complex form of segmentation is often referred to as “dynamic segmentation.” This approach uses changes in conditions on other map layers to subdivide a route. The Conditionsbased segmentation procedure first derives a “universal conditions” map by overlaying a set of map layers to establish individual coincident polygons containing the same combination of conditions throughout their interiors (left side of figure 1). A route is intersected with the universal conditions map and “broken” into a series of line segments at the entry and exit of each universal condition polygon.
The result is a series of variable length line segments with uniform conditions throughout their length. The line segment set is written to a table containing x, y and z coordinates followed by fields identifying the combination of conditions along each segment. In addition, the table can report “crossing counts” that note the number of intersections of the derived line segments and other linear features, such as river and road crossings.
Another advanced form of segmentation involves breaking a route into segments considering elevation changes along a terrain surface as shown in the bottom right portion of figure 1. Terrainbased segmentation divides the route into segments representing similar terrain configuration as determined by major inflection points and slope conditions.
The first step is to establish the elevation profile of the route by “masking” the terrain surface with the route. To eliminate subtle changes, the elevation values are smoothed to characterize the overall trend of the terrain. This is done by passing a movingaverage window along the elevation profile resulting in a smoothed profile.
The width of the smoothing window is critical because if it is too large the smoothing will eliminate potentially important inflection points; if it is too small there will be a multitude of insignificant segment breaks. The dilemma is similar to choosing an appropriate angular change factor in Deflectionangle segmentation and is as much art as it is science. Experience for liquid pipelines suggests an 11 to 15cell diameter window on a 10 meter digital elevation model (DEM) works fairly well.
Figure 2. Analyzing the slope differences along a proposed route.
Figure 2 illustrates the final step in the procedure. The tool evaluates the smoothed elevation values by subtracting the current value from the previous location’s value as the solution progresses lefttoright along the smoothed profile. If the sign is the same, a continuous “running slope” is indicated. A negative difference indicates an upward incline; a negative difference indicates a downward descent; zero difference indicates no change (flat portion). A location where the sign of the difference changes (negative to positive; positive to negative) identifies an inflection point. The coordinates and elevation value for the location that changed sign is written to a Segments Table.
For example, the calculations in the upper left portion of figure 2 show a continued upward incline until the elevation point of 124 is encountered. At this location the sign of the difference changes indicating the start of a downward descent (from a negative to a positive difference).
A refinement to the approach avoids storing and processing the entire elevation profile of actual values by calculating and comparing the moving average values (smoothed elevation) onthefly. This approach is much more efficient and greatly improves performance. In this application a 1x15 cell window is used to smooth the actual elevation profile and store the information in a temporary table with just two records—previous and current average.
The difference is calculated and tested for a change in sign; if different, the coordinates and elevation value for the current location is written to the Segments Table. The window is advanced and the old previous average is replaced and a new smoothed value for the current location is considered. The sequence of steps of 1) replace, 2) calculate, 3) compare and 4) write if sign change is repeated until the entire profile of a route has been evaluated.
The result is a set of map segments with consistent terrain configuration. Summarizing the average slope along the segmented file provides important information for assessing the hydraulics along a pipeline and pumpstation positioning and design required for anticipated product flows …a big step beyond simply mapping a pipeline route.
(GeoWorld, March 2000, pg. 2324)
The past few columns
discussed several techniques for generating maps that identify the bumps
(convex features), the dips (concave features) and the tilt (slope) of a
terrain surface. Although the procedures
have a wealth of practical applications, the hidden agenda of the discussions
was to get you to think of geographic space in a less traditional way—as an
organized set of numbers (numerical data), instead of points, lines and areas
represented by various colors and patterns (graphic map features).
A terrain surface is organized as a rectangular "analysis grid" with
each cell containing an elevation value.
Gridbased processing involves retrieving values from one or more of
these "input data layers" and performing a mathematical or
statistical operation on the subset of data to produce a new set numbers. While computer mapping or spatial database management
often operates with the numbers defining a map, these types of processing
simply repackage the existing information.
A spatial query to "identify all the locations above 8000'
elevation in
Map analysis operations, on the other hand, create entirely new spatial
information. For example, a map of
terrain slope can be derived from an elevation surface and then used to expand
the geoquery to "identify all the locations above 8000' elevation in
Now back to business. Last month's
column described several approaches for calculating terrain slope from an
elevation surface. Each of the
approaches used a "3x3 roving window" to retrieve a subset of data,
yet applied a different analysis function (maximum, minimum, average or
"fitted" summary of the data).
Figure 316. 2D, 3D and
draped displays of terrain slope.
Figure 316 shows
the slope map derived by "fitting a plane" to the nine elevation
values surrounding each map location.
The inset in the upper left corner of the figure shows a 2D display of
the slope map. Note that steeper
locations primarily occur in the upper central portion of the area, while
gentle slopes are concentrated along the left side.
The inset on the right shows the elevation data as a wireframe plot with the
slope map draped over the surface. Note
the alignment of the slope classes with the surface configuration—flat slopes
where it looks flat; steep slopes where it looks steep. So far, so good.
The 3D view of slope in the lower left, however, looks a bit strange. The big bumps on this surface identify steep
areas with large slope values. The
upsanddowns (undulations) are particularly interesting. If the area was perfectly flat, the slope
value would be zero everywhere and the 3D view would be flat as well. But what do you think the 3D view would look
like if the surface formed a steeply sloping plane?
Are you sure? The slope values at each
location would be the same, say 65% slope, therefore the 3D view would be a
flat plane "floating" at a height of 65. That suggests a useful principle—as a slope
map progresses from a constant plane (single value everywhere) to more
upsanddowns (bunches of different values), an increase in terrain roughness
is indicated.
Figure 317 outlines
this concept by diagramming the profiles of three different terrain
crosssections. An elevation surface's 2^{nd}
derivative (slope of a slope map) quantifies the amount of upsanddowns of the
terrain. For the straight line on the
left, the "rate of change in elevation per unit distance" is constant
with the same difference in elevation everywhere along the line—slope = 65%
everywhere. The resultant slope map
would have the value 65 stored at each grid cell, therefore the "rate of
change in slope" is zero everywhere along the line (no change in
slope)—slope2 = 0% everywhere.
Figure 317. Assessing terrain roughness through the 2^{nd} derivative
of an elevation surface.
A slope2 value of zero
is interpreted as a perfectly smooth condition, which in this case happens to
be steep. The other profiles on the
right have varying slopes along the line, therefore
the "rate of change in slope" will produce increasing larger slope2
values as the differences in slope become increasingly larger.
So who cares? Water drops for one, as
steep smooth areas are ideal for downhill racing, while "steep 'n rough
terrain" encourages more infiltration, with "gentle yet rough
terrain" the most.
Figure 318 shows a
roughness map based on the 2^{nd} derivative for the same terrain as
depicted in Figure 316. Note the
relationships between the two figures.
The areas with the most "upsanddowns" on the slope map in
figure 316 correspond to the areas of highest values on the roughness map in
figure 318. Now focus your attention on
the large steep area in the upper central portion of the map. Note the roughness differences for the same
area as indicated in figure 318—the favorite raindrop racing areas are the
smooth portions of the steep terrain.
Figure 318. 2D, 3D and draped displays of terrain
roughness.
Whew! That's a lot of mapematical
explanation for a couple of pretty simple concepts—steepness and
roughness. Next month we'll continue the
trek on steep part of the map analysis learning curve by considering
"confluence patterns" in micro terrain analysis.
Confluence Maps Further
Characterize Microterrain Features
(GeoWorld, April 2000, pg. 2223)
Last section's
discussion focused on terrain steepness and roughness. While the concepts are simple and straightforward,
the mechanics of computing them are a bit more challenging. As you hike in the mountains your legs sense
the steepness and your mind is constantly assessing terrain roughness. A smooth, steeplysloped area would have you
clinging to things, while a rough steeplysloped area would look more like
stair steps.
Water has a similar
vantage point of the slopes it encounters, except given its head,
water will take the steepest downhill path (sort of like an outofcontrol
skier). Figure 3191 shows a 3D grid
map of an elevation surface and the resulting flow confluence. It is based on the assumption that water will
follow a path that chooses the steepest downhill step at each point (grid cell
"step") along the terrain surface.
Figure 319. Map of surface flow confluence.
In effect, a drop of water is placed at each location and allowed to pick its
path down the terrain surface. Each grid
cell that is traversed gets the value of one added to it. As the paths from other locations are
considered the areas sharing common paths get increasing larger values (one +
one + one, etc.).
The inset on the
right shows the path taken by a couple of drops into a slight depression. The inset on the left shows the considerable
inflow for the depression as a high peak in the 3D display. The high value indicates that a lot of uphill
locations are connected to this feature. However, note that the pathways to the
depression are concentrated along the southern edge of the area.
Now turn your attention to figure 3110.
Ridges on the confluence density surface (lower left) identify areas of
high surface flow. Note how these areas
(darker) align with the creases in the terrain as shown on the draped elevation
surface on the right inset. The water
collection in the "saddle" between the two hills is obvious, as are
the two westerly facing confluences on the side of the hills. The 2D map in the upper left provides a more
familiar view of where not to unroll your sleeping bag if flash floods are a
concern.
The various spatial
analysis techniques for characterizing terrain surfaces introduced in this
series provide a wealth of different perspectives on surface
configuration. Deviation from Trend,
Difference Maps and Deviation Surfaces are used to identify areas that
"bumpup" (convex) or "dipdown" (concave). A Coefficient of Variation Surface looks at
the overall disparity in elevation values occurring within a small area. A Slope Map shares a similar algorithm
(roving window) but the summary of is different and reports the
"tilt" of the surface. An
Aspect Map extends the analysis to include the direction of the tilt as well as
the magnitude. The Slope of a Slope Map
(2^{nd} derivative) summarizes the frequency of the changes along an
incline and reports the roughness throughout an elevation surface. Finally, a Confluence Map takes an extended
view and characterizes the number of uphill locations connected to each
location.
Figure 3110. 2D, 3D and draped displays of surface flow
confluence.
The coincidence of
these varied perspectives can provide valuable input to decisionmaking. Areas that are smooth, steep and coincide
with high confluence are strong candidates for gullywashers that gouge the
landscape. On the other hand, areas that are rough, gentlysloped and with minimal confluence
are relatively stable. Concave features
in these areas tend to trap water and recharge soil moisture and the water
table. Convex features under erosive
conditions tend to become more prominent as the confluence of water flows
around it.
Similar interpretations can be made for hikers, who like raindrops react to
surface configuration in interesting ways.
While steep, smooth surfaces are avoided by all but the rockclimber,
too gentle surfaces tend too provide boring
hikes. Prominent convex features can
make interesting areas for viewing—from the top for hearty and from the bottom
for the aesthetically bent. Areas of
water confluence don't mix with hiking trail unless a considerable number of
waterbars are placed in the trail.
These "rulesofthumb" make sense in a lot of situations,
however, there are numerous exceptions that can undercut them. Two concerns in particular are important—
conditions and resolution. First,
conditions along the surface can alter the effect of terrain
characteristics. For example, soil
properties and the vegetation at a location greatly effects surface runoff and
sediment transport. The nature of accumulated
distance along the surface is also a determinant. If the uphill slopes are long steep, the
water flow has accumulated force and considerable erosion potential. A hiker that has been hiking up a steep slope
for a long time might collapse before reaching the summit. If that steep slope is southerly oriented and
without shade trees, then exhaustion is reached even sooner.
In addition, the resolution of the elevation grid can effect
the calculations. In the case of water
drops the gridding resolution and accurate "Z" values must be high to
capture the subtle twists and bends that direct water flow. A hiker on the other hand, is less sensitive
to subtle changes in elevation. The rub
is that collection of the appropriate elevation is prohibitively expensive in
most practical applications. The result
is that existing elevation data, such as the USGS Digital Terrain Models (DTM),
are used in most cases by default. Since
the
The recognition of the importance of spatial analysis and surface modeling is
imperative, both for today and into the future.
Its effective use requires informed and wary users. However, as with all technological things,
what appears to be a data barrier today, becomes
routine in the future. For example,
The more important limitation is intellectual.
For decades, manual measurement, photo interpretation and process
modeling approaches have served as input for decisionmaking involving terrain
conditions. Instead of using
_______________________
Author's
Note: The figures presented in this series on
"Characterizing MicroTerrain Features" plus several other
illustrative ones are available online as a set of annotated PowerPoint slides
at the "Column Supplements" page at http://www.innovativegis.com/basis.
(GeoWorld, May 2000, pg. 2223)
Previous discussions
suggested that combining derived maps often is necessary for a complete
expression of an application. A simple
erosion potential model, for example, can be developed by characterizing the
coincidence of a Slope Map and a Flow Map (see the previous two columns in this
series). The flowchart in the figure
3111 identifies the processing steps that form the model— generate slope and
flow maps, establish relative classes for both, then combine.
While a flowchart of
the processing might appear unfamiliar, the underlying assumptions are quite
straightforward. The slope map
characterizes the relative "energy" of water flow at a location and
the confluence map identifies the "volume" of flow. It's common sense that as energy and volume
increases, so does erosion potential.
Figure 3111. A simple erosion potential model combines information
on terrain steepness and water flow confluence.
The various
combinations of slope and flow span from high erosion potential to deposition
conditions. On the map in the lower
right, the category "33 Heavy Flow; Steep" (dark blue) identifies areas
that are steep and have a lot of uphill locations contributing water. Loosened dirtballs under these circumstances
are easily washed downhill. However,
category "12 Light Flow; Moderate" (light green) identifies locations
with minimal erosion potential. In fact,
deposition (the opposite of erosion) can occur in areas of gentle slope, such
as category "11 Light Flow; Gentle" (dark red).
Before we challenge the scientific merit of the simplified model, note the
basic elements of the
The remaining sentences in the macro and the corresponding boxes/lines in the
flowchart complete the model. The macro
enables entering, editing, executing, storing and retrieving the individual
operations that form the application.
The flowchart provides an effective means for communicating the
processing steps. Most "
Figure 3112. Extended erosion model
considering sediment loading potential considering intervening terrain
conditions.
It provides a
starting place for model refinement as well.
Suppose the user wants to extend the simple erosion model to address
sediment loading potential to open water.
The added logic is captured by the additional boxes/lines shown in
figure 3112. Note that the upper left
box (Erosion_classes) picks up where the flowchart in
figure 3111 left off.
A traditional
approach would generate a simple buffer of a couple of hundred feet around the
stream and restrict all dirt disturbing actives to outside the buffered
area. But are all bufferfeet the
same? Is a simple geographic reach on
either side sufficient to hold back sediment?
Do physical laws apply or merely civil ones that placate planners?
Common sense suggests that the intervening conditions play a role. In areas that are steep and have high water
volume, the setback ought to be a lot as erosion potential is high. Areas with minimal erosion potential require
less of a setback. In the schematic in
figure 3112, a dirt disturbing activity on the steep hillside, though 200 feet
away, would likely rain dirtballs into the stream. A similar activity on the other side of the
stream, however, could proceed almost adjacent to the stream.
The first step in extending the erosion potential model to sediment loading
involves "calibrating" the intervening conditions for dirtball
impedance. The friction map identified
in the flowchart ranges from 1 (very low friction for the 33 Heavy Flow: Steep
condition) to 10 (very high friction for 11 Light Flow; Gentle). A loose dirtball in an area with a high
friction factor isn't going anywhere, while one in an area of very low friction
almost has legs of its own.
The second processing step calculates the effective distance from open water
based on the relative friction. The
command, "SPREAD Water TO 50 THRU Friction OVER Elevation Uphill Across For Water_Prox," is
entered simply by completing a dialog box.
The result is a variablewidth buffer that reaches farther in areas of
high erosion potential and less into areas of low potential. The lighter red tones identify locations that
are effectively close to water from the perspective of a dirtball on the
move. The darker green tones indicate
areas effectively far away.
But notice the small dark green area in the lower left corner of the map of
sediment loading potential. How can it
be effectively far away though right near a stream? Actually it is a small depression that traps
dirtballs and can't contribute sediment to the stream— effectively infinitely
far away.
Several other realworld extensions are candidates to improve model. Shouldn't one consider the type of soils? The surface roughness? Or the time of year? The possibilities are numerous. In part, that's the trouble with
Identify Valley Bottoms in Mountainous Terrain
(GeoWorld, November 2002, pg. 2223)
On occasion I am asked how I come up with ideas for the
Beyond Mapping column every month for over a dozen years. It’s easy—just hang
around bright folks with interesting questions.
It’s also the same rationale I use for keeping academic ties in
For example, last section’s discussion on “Connecting Riders with Their Stops”
is an outgrowth of a grad student’s thesis tackling an automated procedure for
routing passengers through alternative public transportation (bike, bus and
light rail). This month’s topic picks up
on another grad student’s need to identify valley
bottoms in mountainous terrain.
The upperleft inset in figure 1 identifies a portion of a
30m digital elevation model (DEM). The “floor”
of the composite plot is a 50meter contour map of the terrain coregistered
with an exaggerated 3D lattice display of the elevation values. Note the steep gradients in the southwestern
portion of the area.
Figure 1. Gentle slopes surrounding streams are
visually apparent on a 3D plot of elevation.
The enlarged 3D grid display on the right graphically
overlays the stream channels onto the terrain surface. Note that the headwaters occur on fairly
gentle slopes then flow into the steep canyon.
While your eye can locate areas of gentle slope surrounding streams in
the plot, an automated technique is needed to delineate and characterize the
valley bottoms for millions of acres. In
addition, a precise rule set is needed to avoid inconsistencies among
subjective visual interpretations.
Like your eye, the computer needs to 1) locate areas with gentle slopes that are 2) connected to the stream channels. In evaluating the first step the computer doesn’t “see” a colorful 3D plot like you do. It “sees” relative elevation differences—rather precisely—that are stored in a matrix of elevation values. The slope algorithm retrieves the elevation value for a location and its eight surrounding values then mathematically compares the subset of values. If the nine elevation values were all the same, a slope of zero percent (flat) would be computed. As the elevation values in one portion of the 3x3 window get relatively large compared to another portion, increasing slope values are calculated. You see a steep area on the 3D plot; the computer sees a large computed slope value.
Figure 2 shows the slope calculation results for the project area. The small 3D grid plot at the top identifies terrain steepness from 0 to 105% slope. Gently sloped areas of 06% are shown in green tones with steeper slopes shown in red tones. The gently sloped areas are the ones of interest for deriving valley bottoms, and as your eye detected, the computer locates them primarily along the rim of the canyon. But these areas only meet half of the valley bottom definition—which of the flat areas are connected to stream channels?
Figure 2. Slope values for the terrain surface are
calculated and gently sloped areas identified.
Figure 3 identifies the procedure for linking the two criteria. First the flat areas are isolated for a binary map of 06% slope=1 and greater slopes=0. Then an effective distance operation is used to measure distance away from streams considering just the areas of gentle slopes. The enlarged inset shows the result with connected valley bottom as far away as 900 meters.
While the prototype bottoming technique appears to work,
several factors need to be considered before applying it to millions of
areas. First is the question whether the
30 meter DEM data supports the slope calculations. In figure 2 the horizontal stripping in the
slope values suggests that the elevation data contains some
inconsistencies. In addition what type
of slope calculation (average, fitted, maximum, etc) ought to be used? And is the assumption of 06% for defining
valley bottoms valid? Should it be more
or less? Does the definition change for
different regions?
Figure 3. Valley bottoms are identified as flat areas
connected to streams.
Another concern involves the guiding surface used in
determining effective proximity—distance measurement uphill/downhill only or
across slopes? Do intervening soil and vegetation
conditions need to be considered in establishing connectivity? A final consideration is designing a field
methodology for empirically evaluating model results. Is there a strong correlation with riparian
vegetation and/or soil maps?
That’s the fun of scholarly pursuit.
Take a simple idea and grind it to dust in the research crucible—only
the best ideas and solutions survive.
But keep in mind that the major ingredient is a continuous flow of
bright and energetic minds.
_________________
Author's
Note: Sincerest
thanks to grad students Jeff Gockley and Dennis
Staley with the
Use Surface Area for Realistic Calculations
(GeoWorld, December 2002, pg. 2223)
The earth is flat …or so most
Common sense suggests that if you walk up and down in steep terrain you will
walk a lot farther than the planimetric length of your path. While we have known this fact for thousands
of years, surface area and length calculations were practically impossible to
apply to paper maps. However, mapematical processing in a
In a vectorbased system, area calculations use plane geometry equations. In a gridbased system, planimetric area is calculated by multiplying the number of cells defining a map feature times the area of an individual cell. For example, if a forest cover type occupies 500 cells with a spatial resolution of 1 hectare each, then the total planimetric area for the cover type is 500 cells * 1ha/cell = 500ha.
However, the actual surface area of each grid cell is
dependent on its inclination (slope). Surface
area increases as a grid cell is tilted from a horizontal reference
plane (planimetric) to its threedimensional position in a digital elevation
surface (see figure 1).
Figure 1. Surface area increases with increasing
terrain slope.
Surface area can be calculated by dividing the planimetric area by the cosine of the slope angle (see Author’s Note). For example, a grid location with a 20 percent slope has a surface area of 1.019803903 hectares based on solving the following sequence of equations:
Tan_Slope_Angle = 20 %_ slope / 100
= .20 slope_ratio
Slope_Angle = arctan (.20 slope_ratio) =
11.30993247 degrees
Surface_Area_Factor = cos (11.30993247 degrees) = .980580676
Surface_Area = 1 ha / .980580676
= 1.019803903 ha
The table in figure 2 identifies the surface area calculations for a 1hectare gridding resolution under several terrain slope conditions. Note that the column Surface_Area_Factor is independent of the gridding resolution. Deriving the surface area for a cell on a 5hectare resolution map, simply changes the last step in the 20% slope example (above) to 5 ha / .9806 = 5.0989 ha Surface_Area.
Figure 2. Surface areas for selected
terrain slopes.
For an empirical test of the surface area conversion
procedure, I have students cut a stick of butter at two angles— one at 0
degrees (perpendicular) and the other at 45 degrees. They then stamp an imprint of each patty on a
piece of graph paper and count the number of grid spaces covered by the grease
spots to determine their areas.
Comparing the results to the calculations in the table above confirms
that the 45degree patty is about 1.414 larger.
What do you think the area difference would be for a 60degree patty?
In a gridbased
Step 1) SLOPE
ELEVATION FOR SLOPE_
Step 2) CALCULATE ( COS
( ( ARCTAN ( SLOPE_
Step 3) CALCULATE 1.0 / SURFACE_
Step 4) COMPOSITE HABITAT_DISTRICT_
FOR HABITAT_SURFACE_
The command macro is scaled to a 1hectare gridding resolution project area by
assigning the value 1.0 (ha) in the third step.
If a standard DEM (Digital Elevation Model with 30m resolution) surface
is used, the GRID_
A regionwide, or zonal overlay procedure (Composite), is used
to sum the surface areas for the cells defining any map feature—habitat
districts in this case. Figure 2 shows
the planimetric and surface area results for the districts considering the
terrain surface. Note that District 2
(light green) aligns with most of the steep slopes on the elevation
surface. As a result, the surface area
is considerably greater (9.62%) than its simple planimetric area.
Figure 3. Planimetric vs. surface
area differences for habitat districts.
Surface length of a line also is affected by terrain. In this instance, azimuth as well as slope of
the tilted plane needs to be considered as it relates to the direction of the
line. If the grid cell is flat or the
line is perpendicular to the slope of a tilted plane, there is no correction to
the planimetric length of a line—from orthogonal (1.0 grid space to diagonal
(1.414 grid space) length. If the line
is parallel to the slope (same direction as the azimuth of the cell) the full
cosine correction to its planimetric projection takes hold. The geometric calculations are a bit more
complex and reserved for online discussion (see Author’s Note below).
So who cares about surface area and length?
Suppose you needed to determine the amount of pesticide needed for weed
spraying along an updown powerline route or estimating the amount of seeds
for reforestation or are attempting to calculate surface infiltration in rough
terrain. Keep in mind that “simple
area/length calculations can significantly underestimate real world
conditions.” Just ask a pooped hiker
that walked five planimetricmiles in the Colorado Rockies.
_________________
Author's
Note: For background
theory and equations for calculating surface area and length… see www.innovativegis.com/basis, select
“Column Supplements.”
Beware Slope’s Slippery Slope
(GeoWorld,
January 2003, pg. 2223)
If you are a
skier it’s likely that you have a good feel for terrain slope. If it’s a steep doublediamond run, the timid
will zigzag across the slope while the boomers will follow the fallline
straight down the hill.
This suggests that there are myriad slopes—uphill, across and downhill; N, E, S
and West; 0360 azimuth—at any location.
So how can a slope map report just one slope value at each
location? Which slope is it? How is it calculated? How might it be used?
Figure 1. A terrain
surface is stored as a grid of elevation values.
Digital
Elevation Model (DEM) data is readily available at several spatial
resolutions. Figure 1 shows a portion of
a standard elevation map with a cell size of 30 meters (98.43 feet). While your eye easily detects locations of
steep and gentle slopes in the 3D display, the computer doesn‘t have the
advantage of a visual representation.
All it “sees” is an organized set of elevation values—10,000 numbers
organized as 100 columns by 100 rows in this case.
The enlarged portion of figure 1 illustrates the relative positioning of nine
elevation values and their corresponding grid storage locations. Your eye detects slope as relative vertical
alignment of the cells. However, the
computer calculates slope as the relative differences in elevation values.
The simplest approach to calculating slope focuses on the eight surrounding
elevations of the grid cell in the center of a 3by3 window. In the example, the individual percent slope
for de is change in elevation
(vertical rise) divided by change in
position (horizontal run) equals
[((8071 ft – 8136 ft) / 98.43 ft) * 100] = 66%. For diagonal positions, such as ae, the calculation changes to [((8071
ft – 8136 ft) / 139.02 ft) * 100] = 47% using an adjusted horizontal run of SQRT(98.43**2 + 98.43**2) = 139.02 ft.
Applying the calculations to all of the neighboring slopes results in eight
individual slope values yields (clockwise from ae) 47, 0, 38, 54, 36, 20, 45 and 66 percent slope. The minimum slope of 0% would be the
choice the timid skier while the boomer would go for the maximum slope of 66%
provided they stuck to one of the eight directions. The simplest calculation for an overall slope
is an arithmetic average slope of 38% based on the eight individual slopes.
Figure 2. Visual comparison of different slope maps.
Most
The equations and example calculations for the advanced techniques are beyond
the scope of this column but are available online (see author’s note). The bottom line is that gridbased slope
calculators use a roving window to collect neighboring elevations and relate
changes in the values to their relative positions. The result is a slope value assigned to each
grid cell.
Figure 2 shows several slope maps derived from the elevation data shown in
figure 1 and organized by increasing overall slope estimates. The techniques show somewhat similar spatial
patterns with the steepest slopes in the northeast quadrant. Their derived values, however, vary
widely. For example, the average slope estimates range from .4 to
53.9% whereas the maximum slope
estimates are from 3.3 to 120%. Slope
values for a selected grid location in the central portion are identified on each
map and vary from 0 to 59.4%.
So what’s
the takehome on slope calculations?
Different algorithms produce different results and a conscientious user
should be aware of the differences. A
corollary is that you should be hesitant to download a derived map of slope
without knowing its heritage. The chance
of edgematching slope maps from a variety of systems is unlikely. Even more insidious is the questionable
accuracy of utilizing a mapematical equation
calibrated using one slope technique then evaluated using another.
Figure 3. Comparison of Fitted versus Generalized slope maps.
The two advanced techniques result in very similar slope estimates. Figure 3 “mapematically” compares the two maps by simply subtracting them. Note that the slope estimates for about twothirds of the project area is within one percent agreement (grey). Disagreement between the two maps is predominantly where the fitted technique estimates a steeper slope than the generalized technique (bluetones). The locations where generalized slopes exceed the fitted slopes (redtones) are primarily isolated along the river valley.
Shedding Light on Terrain Analysis
(GeoWorld, April 2008,)
A lot of GIS is straightforward and mimics our boy/girl scout days wrestling with paper maps. We learned that North is at the top (at least for us in the northern latitudes) and red roads and blue streams wind their way around green globs of forested areas. However, the brown concentric rings presented a bit more of a conceptual challenge.
When the contour lines formed a fairly small circle it was likely a hilltop; or a depression, if the line sprouted whiskers. Sharp “Vshaped” contour lines pointed upstream when a blue line was present; and somewhat rounded V’s pointed downhill along a ridge. Once these and a few other subtle nuances were mastered and you survived a day or so in the woods, a merit badge and sense power over geography was attained.
Your computer is devoid of such a nostalgic experience. All it sees are organized sets of colorless numbers. In the case of vector systems, these number sets identify lines defining discrete spatial features in a manner fairly analogous to traditional mapping. On the surface, direction measurement fundamentals remain the same and your fond memory of the 0360^{o} markings around a compass ring holds for most GIS mapping applications.
Figure 1. Azimuth is a discontinuous
numeric scale as it wraps around on itself; Facing
Angle is a continuous gradient indicating relative alignment with a specified
direction.
The upperleft portion of figure 1 unfolds the compass ring into a continuum of azimuth from 0^{o} (north) to 360^{o} (north). Whoa …you mean north is at both ends of the continuum? As azimuth increases, it becomes less northlike for awhile and then more northlike until 0 = 360. Now that’s enough to really messup a numeric beancounter like a computer. Actually, azimuthal direction is termed a discontinuous number set as it wraps around on itself (spiral in the topright side of figure 1).
The lower scales in figure 1 transform direction into a more stable continuous number set that indicates the relative alignment with a specified direction. A Facing Angle utilizes the concept of a back azimuth to indicate an orientation as different as it gets—180^{o} is completely opposite; and 0^{o} is exactly the same direction.
For a south facing location the continuum starts at 0 (South) and progresses “less southlike” in both the East and West directions until the orientation “is 180 degrees off” at due North. The importance of the consistency of this continuous scale might be lost on humans, but it means the world to a computer attempting to statistically analyze a set of directional data. The procedure normalizes terrain data for any given facing angle into a consistent scale of 0 to 180 indicating the degree of orientation similarly throughout an entire project area.
Figure 2 summarizes the geometry between Azimuth and Facing
angles. Also it introduces the concept
of a Horizontal angle that takes direction from 2D planar to 3D solid
geometry. The “Hillshade”
map of
Figure 2. A Hillshade
map identifies the relative brightness at every grid location in a project
area.
Figure 3 expresses the relationship between the azimuth and horizontal angle gradients derived from aspect and slope maps. Each grid cell (30x30 meters in this example) can be thought of as a tilted plane in threedimensional space. In the case of the brightest location, it identifies a grid cell that is perpendicular to the sun—like a contented lizard or sunbather positioned to absorb the most rays. All other possible orientations of the grid cell facets are some combination of the facing and horizontal angle gradients.
However, the “superstar mapematicians” among us know that things are a bit more complex than an independent peek at each grid cell. A grid cell could be directly oriented toward the sun but on a small hill in the shadow of a much larger hill. To account for the surrounding terrain configuration one needs to solve the angles using solid geometry algorithms involving direction cosines to connect the grid cell to the sun and test if there is any intervening terrain blocking connection.
So what’s the takehome from all this discussion? Actually three points should stand out. First, that GIS technology encompasses our paper map thinking (e.g., discontinuous azimuthal direction), but the digital map enables us to go beyond mapping. Secondly, that a map is an organized set of numbers first, and a picture later; we must understand the nature of the data (e.g., continuous facing angle) to fully capitalize on the potential of map analysis.
Figure 3. Terrain orientation combines azimuthal and horizontal angles for each grid cell facet.
Finally, the new mapematical expression of maps supports radically new applications. For example, one can integrate a series of brightness maps over a day for a map of solar influx (insolation) that drives a wealth of spatial systems from wildlife habitat to global warming. In addition, these data can be statistically analyzed for insight into spatial patterns, relationships and dependencies that were beyond discovery a few years ago. In short, GIS is ushering in a new era of decisionmaking, not simply keeping track of where is what.
_____________________________
Author’s Note: related discussion on Terrain
Analysis is in Topic 6, Summarizing Neighbors in the book Map Analysis (
Correction:
the following letter (
Joseph—In the many years of Beyond Mapping I’ve never had a serious
disagreement. In the May column it seems as though your slant on angles is a
bit off course. I’m no mapematician but azimuth is not discontinuous; rather it is a
cyclic attribute type that continuously wraps around (and around). If we
thought of phase angles as discontinuous we would have a hard time using sines and cosines to compute mean wind directions slope
aspects or buoy longitudes east of
Best
wishes, Peter H. Dana, Ph.D.,
Research Fellow and Lecturer, Department of Geography,
Peter—right
on …excellent points well taken. Azimuth is not “discontinuous” and your
classification of “cyclical” is a good one. The bottom line is that
azimuth isn’t your normal numerical gradient and being a bit whacko one can’t
simply utilize raw azimuth values in spatial data mining. My use of the
term “brightest areas” was misleading …I was referring to surface illumination
intensity (amount of sunlight impacting a location as contained in the facing
angle) not the relative brightness in the image which is dependent on viewer
angle as well as solar angle as compounded by the unique lay of the
terrain. Your recounting of the imaging interpretation rule also holds
true for interpreting the graphic. My point however was less on the
interpretation of the map graphic as on the ability to track solar insolation …mapped data versus visualization. The
bottom line is that I am delighted that someone out there actually reads the
column, particularly with your level of interest and expertise. The
dialog is much appreciated and I stand corrected. Thank you. Joe
Identifying Upland
Ridges
(GeoWorld, May 2009)
The good news is that there is a tsunami of mapped data out there; the bad news is that there is a tsunami of mapped data out there. Rarely is it as simple as downloading and displaying the ideal map for decisionmaking. More often than not the base data that is available is just that—a base for further exploration and translation into meaningful information for solutions wellbeyond a basic wall decoration.
Digital Elevation Models (DEMs) are no exception. While there is an ever increasing wealth of elevation data at ever increasing detail, most applications require a translation of the ups and downs into decision contexts. For example, suppose you are interested in identifying upland ridges in mountainous terrain for wildfire considerations, wind power location or animal corridors. Downloading and displaying a set of DEMs provides a qualitative glimpse at the topography but most human interpretation becomes overwhelmed by the sheer magnitude of the possibilities and the inability to apply a consistent definition.
A topographic ridge can be defined as a long narrow upper section or raised crest of elevation formed by the juncture of two sloping planes. While the definition of a ridge is straightforward, its practical expression is a bit murkier. Figure 1 outlines a commonly used gridbased map analysis procedure for identifying ridges. The first step involves “smoothing” the surface to eliminate subtle and often artificial elevation changes. This is accomplished by moving a small “roving window” over the surface that averages the localized elevation values in the vicinity of each map location (1 Scan).
The next step simulates water falling across the entire project area and moving downhill by the steepest path (2 Drain). The computer keeps track of how many uphill locations contribute water to each grid location to generate the flow accumulation map draped over the smoothed elevation surface shown in the upperleft portion of the figure. Note that lower values (grey tones) indicate locations with minimal uphill contributors and align with the elevation crests.
The third step reclassifies the flow map to isolate the locations containing just one runoff path—candidate ridges (3 Renumber). However note the cluttered pattern of the map in the upperright that confuses a consistent and practical definition of a ridge. The final step in the figure uses another roving window operation to assign the majority value in the vicinity of each map location thereby eliminating the “salt and pepper” pattern of the raw ridges (4 Scan). If most of the values in the window are 0 (not a potential ridge) then the location is assigned 0. However if there is a majority of potential ridges in the window then a 1 is assigned. The draped display in the lowerleft confirms the alignment of the derived ridges with the terrain surface.
Figure 1. Accumulated surface flow is used to identify
candidate ridges.
(MapCalc^{TM}
software commands are indicated to illustrate processing options)
Figure 2 continues the processing to identify just the upland ridges. Contiguous ridge locations are individually identified (5 Clump) and then the average elevation for each ridge grouping is determined (6 Composite). Note the low average elevation of the three small ridges in the center portion of the project area. The final step in the figure (7 Renumber) eliminates candidate ridges that are below the average elevation of the project area thereby leaving only the high ridge locations (black).
Figure 3 illustrates the processing steps for expanding the ridges as a function of terrain steepness. A slope map is generated (8 Slope) and calibrated to isolate areas of gentle slopes of 0 to 5 percent (9 Renumber). In turn, these areas are processed to eliminate the “salt and pepper” pattern (10 Scan) leaving relatively large expanses of gentle slopes.
The lowerright
map identifies the upland ridges (black) superimposed on the gently sloped
areas (red). An effective distance
operation (11 Spread) is used to
extend the ridges considering the nongently sloped locations (grey) as
absolute barriers. The warmer tones in
the proximity map draped on the elevation surface indicate increasing distance
from a ridge into the surrounding relatively flat areas. Depending on the application the upland ridge
expansion could be reclassified from the full extent to just the immediate
surrounding areas (green).
Figure 2. Low ridges are identified and eliminated.
Figure 3. The ridges are extended into surrounding
gently sloped areas to identify effective upland ridges.
The simple model
may be extended to characterize conditions within the derived ridges. For example, an aspect map can be combined
with the extended ridges to indicate their overall terrain orientation or
precise azimuth each grid location defining the upland ridges to consider for
prevailing wind direction or solar exposure.
Another extension could generate an overland traveltime or a road
construction surface to consider relative access potential of the ridges.
The bottom line
is that it takes a bit of spatial reasoning to translate base maps into
decision contexts …it’s the analytic nature of GIS that moves it beyond
mapping.
_____________________________
Author’s Note: see the online book Beyond
Mapping III, Topic 11, Characterizing Micro Terrain Features at www.innovativegis.com/basis/MapAnalysis/ for additional terrain analysis techniques.
(Back to the Table of Contents)