Topographic Wetness Index

Topographic Analysis · Wetness, Slope and Extreme Values · geoprocessing tool · by Jeff Jenness
Works at every ArcGIS Pro license level

Summary

Computes the Topographic Wetness Index (TWI, also called the Compound Topographic Index or CTI) for every cell of a DEM: the natural log of the specific catchment area divided by the local slope gradient (Beven and Kirkby 1979). The index combines the two things about the shape of the land that decide where water collects: how much ground drains to a cell, and how steep the cell is, and consequently how quickly water drains off it. High values mark cells with a large upslope area and a gentle slope, the valley bottoms, hollows and flats where water arrives and lingers; low values mark ridges and steep, diverging slopes that shed it. It is a purely topographic estimate of relative wetness, built from the shape of the land alone, with no knowledge of rainfall, soils or subsurface flow. The sink filling, flow direction and flow accumulation it needs are computed internally, with a choice of D8, MFD or D-Infinity flow direction methods, so no Spatial Analyst extension is required; you may instead supply a filled DEM and a flow-accumulation raster of your own, and you may save the ones the tool computes.

Learn more This page follows the wetness-index slides of the Topographic Roughness lecture from my GIS course (slides; see the training page for the course), and the bonus exercise of the course's Watershed lab, which builds the same index by hand with Esri's hydrology tools and the Raster Calculator. About Topographic Roughness places this index among the other terrain measures in the toolbox.
Project the DEM first, and project it bilinearly TWI is built entirely from flow paths, and flow paths are the first casualty of a DEM projected with nearest-neighbor resampling. The regular lattice of repeated and skipped cells that nearest neighbor leaves behind acts like a network of shallow channels across the surface; flow direction follows them, modeled streams run in straight lines and right-angle turns along the artifact grid, and the wetness index inherits every one of those false channels. If your DEM has to change coordinate systems, project it with bilinear interpolation before running this tool, never with nearest neighbor. The Projecting Rasters page shows what the artifacts look like and why.

The idea

Beven and Kirkby (1979) proposed the index as part of TOPMODEL, a rainfall-runoff model they describe as “a model for humid temperate areas.” The question it answers is simple: where on a landscape does water tend to collect? Two things decide that. Water arrives at a place in proportion to how much ground drains to it, and it leaves in proportion to how steep the place is. A cell at the bottom of a broad basin receives water from every cell above it and, being nearly level, passes it on slowly. A cell on a ridge crest receives nothing from anywhere and sheds what falls on it at once. The wetness index puts the two together as a ratio and takes its logarithm:

TWI= ln( atanβ )

Here a is the specific catchment area, explained in the next section, and tan β is the slope gradient at the cell, its rise over run. That is the tangent of the slope angle β, which is percent slope divided by 100, and not the slope in degrees or in percent. The logarithm is there because the ratio spans an enormous range, from values near 0 on a steep ridge to millions or more at the mouth of a large basin; the log compresses that into a scale of roughly 0 to 25 that can be mapped and compared. The index rises as catchment area grows, because larger areas collect more water in general, and rises as slope falls, because water moves more slowly across gentle ground. The lecture calls it a sort of measure of water pressure or concentration at a location, and that is a fair way to picture it.

The same quantity goes by several names. Beven and Kirkby wrote it as ln(a/tan β) and called it a topographic index; the USGS distributed it with the HYDRO1k global elevation derivatives as the Compound Topographic Index, “commonly referred to as the Wetness Index,” following Moore et al. (1991); Topographic Wetness Index is the name most ecological papers use now, and the name this tool uses.

Specific catchment area and the unit contour length

Catchment area is easy: it is the area of ground that drains to a point, the watershed above it. Specific catchment area is the part that confuses people, and it took me years to get a handle on it, because watersheds do not always drain to a point. Picture a spring pouring water over a convex surface such as a cone. As the water runs down, it spreads out. There is no single point where the flow from the spring collects; the outlet of that little watershed is a line along the cone, and the watershed that feeds the line is a fan-shaped strip between two flow paths.

Two green cones with gridlines. On the left, a narrow blue triangle runs from a point partway up the cone down to its base, widening as it descends: water diverging as it flows down a convex surface. On the right, the same fan is drawn in purple above a short red line along a contour near the base: the red line is the watershed outlet, and the purple fan is the area draining to it.
From the lecture: water leaving a point on a convex surface diverges as it flows down (left). The outlet of that watershed is a line, drawn in red on the right, and the watershed feeding it is the purple fan. Catchment area alone cannot describe this; the area has to be standardized by the length of the outlet line.

That is what specific catchment area does. It is the catchment area divided by the length of the line across which the water leaves it, the length of the contour segment through which the flow passes. The literature calls that line the unit contour length (Beven and Kirkby 1979; Sørensen et al. 2006) and defines the specific catchment area as the upstream area per unit width of contour (Moore et al. 1991; Yang et al. 2011). Its units are square meters divided by meters, which reduces mathematically to meters, but the units make better sense left as they are: read them as square meters of watershed per meter of outlet line. On a diverging surface like the cone, the area per meter of outlet is small; where flow converges into a channel, all of the watershed passes through a short outlet line and the number is large. That is why the index needs the specific area and not the plain one: a point on the cone has a catchment area of essentially zero, and the index would have nothing to work with.

What a raster does with this

Raster flow analysis changes the picture in one important way. Flow direction and flow accumulation treat every cell as a watershed outlet in its own right, and the accumulation raster tells you how many cells drain into each cell. You cannot set the whole red line on the cone as an outlet; each cell along it gets its own count. So in a raster the outlet line of a cell's watershed is the edge of that cell, the unit contour length is the cell size, and the derivation runs in three lines:

QuantityFormulaUnits
Catchment area of a cell(flow accumulation + 1) × cell size × cell sizem²
Unit contour lengthcell sizem
Specific catchment area, a(flow accumulation + 1) × cell sizem² per m

One cell size cancels. The “+ 1” is there because Esri's Flow Accumulation, and this tool's own accumulation, count the cells upstream of a cell and not the cell itself. This means that the Flow Accumulation algorithms do not consider the cell itself to be contributing to its own accumulated flow, so we add “1” to the flow accumulation value in order to force the cell's catchment area to include the cell itself. Rain that falls on a cell is in that cell. Adding one also keeps the logarithm defined at the head of every drainage, where the count is zero and ln(0) would otherwise be the result. This tool computes the specific catchment area exactly as the table shows.

A note on Geographic DEMs: On a geographic (latitude/longitude) DEM a cell is not square on the ground: at 36° north a one-arc-second cell is about 25 m east to west and 31 m north to south, and the east-west span changes with every row. The tool measures both spans in meters on the spheroid for each row, takes their product as the cell area, and takes the square root of that area, the side of a square of the same size, as the unit contour length, about 28 m in this example. One length serves every cell whatever direction its water leaves by, for the same reason the diagonal is not used for diagonal flow: a divisor that changed with the direction of flow would tilt the index by compass direction for no physical reason.

One more raster question is which edge of the cell to divide by. Water leaves a cell in one of eight directions under the D8 rule, and when it leaves diagonally the line it crosses is, strictly, the cell's diagonal, 1.414 times the edge. Dividing by the diagonal for diagonal flow and the edge for the rest would make the index systematically higher wherever water happens to flow north, east, south or west, and nothing in the real world justifies that, so this tool divides every cell by the same constant, the cell size, whatever direction the water takes.

A worked example

The table shows the index for a 30 m DEM at a few combinations of flow accumulation and slope. Flow accumulation is Esri's count of upslope cells, so the first row is a cell at the head of a drainage with nothing above it and the last is a cell with 10,000 cells, nine square kilometers, draining through it. The constant added to the slope is the tool's default of 0.001, which matters only in the level column.

Upslope cellsCatchment area (m²)Specific catchment area a (m²/m)TWI, level groundTWI, 5° (tan 0.087)TWI, 20° (tan 0.364)TWI, 45° (tan 1.0)
09003010.35.84.43.4
99,00030012.68.16.75.7
9990,0003,00014.910.49.08.0
999900,00030,00017.212.711.310.3
9,9999,000,000300,00019.515.013.612.6

Read down a column and each tenfold increase in catchment adds 2.3 to the index, which is ln(10): the log turns multiplication into addition. Read across a row and steepness pulls the index down, but more gently, because going from 5° to 45° multiplies the tangent by only 11.4, worth 2.4 on the log scale. Now look at the level column. A single cell with nothing draining into it scores 10.3 on level ground, about the same as a hundred-cell watershed at 5°. That is the added constant talking: on a perfectly level cell the tangent is zero, the ratio would be infinite, and the small constant added to the slope is what keeps it finite, but it keeps it very large. A lake surface, a playa, a flat roof of a mesa, or a level stretch of DEM produced by an integer elevation grid will all score near the top of the map. Whether that is right depends on what you are asking; a flat is a place where water lingers, but a flat on a summit is not wet. The slope section says more about that constant.

The slope term

The slope is computed from the original DEM, not from the filled one. The DEM is filled only so that flow routing does not get stuck in sinks; for every other purpose the true, unmodified surface is the one to use, and a filled DEM would report zero slope across the floor of every filled depression. The gradient is Horn's method, the same 3 × 3 operator as Esri's Slope tool, and its output is a rise over run, the tangent of the slope angle. Because it is a tangent, the DEM's elevation units must match its horizontal units: the tool reads the units from the DEM's vertical coordinate system when it has one and asks you only when the DEM does not say, and it converts internally so that a DEM in feet gives the same gradient as the same terrain in meters.

On a projected DEM the slope is Planar by default, measured in the projection's own plane with the raster's cell size, which matches Esri's Slope tool with the Planar method and is the right choice for nearly every projected DEM. The Geodesic choice measures the distances between cell centers on the spheroid instead, correcting for the projection's scale distortion; on a UTM or State Plane DEM it is essentially identical to Planar, and it earns its place only on large-extent or highly distorted projections such as Web Mercator across a wide range of latitudes. It is a true-ground-distance correction and does not attempt to reproduce Esri's geocentric Geodesic method. On a geographic DEM the choice is disabled: slope is always computed per row on the spheroid, because a degree of longitude is a different distance at every latitude.

The value added to slope is a small constant added to the gradient before dividing, so that a level cell does not divide by zero. The default 0.001 is the value the USGS used for flat cells in the HYDRO1k wetness index, and the value my Watershed lab uses; the tool adds it to every cell's slope, level or not. Its effect is confined to nearly level ground: at 5° the tangent is 0.0875, the constant changes the index by 0.01, and at any real slope it does nothing you could see. On a perfectly level cell it is the whole denominator. The HYDRO1k documentation substitutes 0.001 for the slope of a flat cell rather than adding it everywhere. The two conventions give the same value on flat cells, and elsewhere they differ by ln((tan β + 0.001) / tan β): about a hundredth at 5°, 0.06 at 1°, 0.11 at 0.5°, and 0.21 on the gentlest slope a 30 m DEM with whole-meter elevations can express, a rise of 1 m across the 240 m span of Horn's window. On ground that steep or steeper the difference is invisible; on nearly level ground it is a few tenths. If the level areas of your DEM are scoring higher than you believe they should, raising the constant to 0.01 lowers a level cell's index by 2.3 and moves a 5° cell by only 0.1. When the flats are lakes, mask the lakes out of the output rather than out of the DEM: a NoData hole in the DEM is an outlet, as the next section explains, and every cell downstream of a masked lake would lose the lake's catchment.

How the hydrology is computed

Everything is done inside the tool with NumPy, so it runs at any ArcGIS Pro license and needs no Spatial Analyst extension. Three steps produce the flow accumulation.

Fill. Real DEMs are full of small closed depressions, some real and most of them artifacts, and water routed across such a surface stops in every one of them. Sinks are removed with the priority-flood algorithm of Barnes, Lehman and Mulla (2014), which floods the surface inward from its edges, raising each depression to the level of its spill point. The result matches Esri's Fill tool cell for cell.

Flow direction. Each cell of the filled surface passes its water to one or more of its eight neighbors, by the routing method you choose; the next section describes the three. Cells on a flat, which have no downhill neighbor, drain across the flat to its nearest exit, so every cell drains somewhere and no water circles forever.

Flow accumulation. Following the directions downstream, each cell's value is the area whose water passes through it, in cells, not counting itself, which is the convention of Esri's Flow Accumulation tool. Under D8 that is a whole number of cells; under the two dispersive methods (Multiple Flow Direction [MFD] and D-∞ [D-Infinity or DINF]), it is fractional, because upstream cells send only shares of their area. If you need a result that reproduces an Esri workflow exactly, run Esri's Fill, Flow Direction and Flow Accumulation, supply the accumulation raster in the optional group of the dialog, and the tool will use it directly and skip its own hydrology; it still needs the original DEM for the slope.

Three ways to route the water: D8, D-Infinity and MFD

These are the same three algorithms that Esri's Flow Direction and Flow Accumulation tools offer, but not Esri's implementations of them. This tool is designed to run without Spatial Analyst, so it carries its own versions of all three, written from the published descriptions. How closely their results match Esri's is reported below.

D8 (O'Callaghan and Mark 1984), the default, sends all of a cell's water to the one neighbor that lies most steeply below it. That is the right thing in a channel and the wrong thing on the cone of the lecture, where water genuinely diverges: under D8 a cell on a convex hillside receives water only from the one thread of cells that happens to point at it, and its catchment stays small no matter how much ground lies above it. As the lecture puts it, calculating specific catchment area is difficult on a convex surface, and fortunately we are usually more interested in the index in concave regions, where flow converges and D8 does well. D8 is also Esri's default method for calculating flow accumulation, and the simplest method to understand.

D-Infinity (Tarboton 1997) fits eight triangular facets around the cell, one between each pair of adjacent neighbors, and takes the steepest facet's downslope direction as a continuous angle rather than one of eight compass points. The cell's water is then split between the two neighbors that bracket that angle, in proportion to how close the angle is to each: a direction 30° from the east neighbor and 15° from the northeast one sends two thirds northeast and one third east. Each receiving cell does the same in its turn, so on the cone the water from a point spreads into a widening fan, which is the purple strip of the lecture slide, and in a channel the two shares collapse toward a single thread.

The steepest facet direction is essentially the cell's aspect, the direction the ground faces, computed on the facets rather than on Horn's 3 × 3 window. Tarboton reports it in the mathematical convention, as a polar angle with zero at east and increasing counter-clockwise, rather than as a compass bearing with zero at north and increasing clockwise, so the angle in a D-Infinity flow-direction raster does not read like an aspect raster even though it points the same way. The aspect tools accept such a raster through their Input aspect convention setting, described on the Aspect Zonal Statistics as Table page.

MFD, the multiple-flow-direction method of Quinn et al. (1991), sends a share to every lower neighbor, weighted by the slope toward it and by the contour length that neighbor represents: half a cell width for a cardinal neighbor and 0.354 of a cell width for a diagonal one, from the geometric construction in their Figure 1, which this tool uses. The tool uses the adaptive form of Qin et al. (2007), which is also what Esri's MFD option implements: the slope is raised to an exponent that grows with the local steepness, so that flow spreads widely on gentle ground and concentrates on steep ground. MFD disperses more than D-Infinity, and on a plane it spreads water sideways where none should go; Tarboton's tests on cones and planes found D-Infinity closer to the true contributing area, while Sørensen et al. (2006), who tested many routing choices against measured groundwater levels, soil moisture, soil pH and plant species richness at two boreal forest sites, found that no single method was best for every variable.

How do our outputs compare with Esri's? On a projected floating-point DEM of the Grand Canyon, the tool's D8 accumulation agrees with Esri's Flow Accumulation on 99.9 percent of cells and its D-Infinity on 99.8 percent, the remainder being flats, where the exit a flat drains to is a matter of convention. Integer DEMs are harder to match, because two or three neighbors often tie for steepest and something has to break the tie; this tool breaks it in the same order Esri's Flow Direction does, and on a 30 m integer DEM its D8 accumulation matched Esri's exactly on 85 percent of cells and within 5 percent on 91 percent, with the rest again on flats. Its MFD agrees with Esri's within 5 percent on 83 percent of cells, with a correlation of 0.97 between the two on the log scale; Esri does not publish every detail of its MFD implementation, so the two are close rather than identical. If you need an Esri result exactly, then just calculate that Flow Accumulation manually using Esri's Flow Accumulation tool, then give that raster to this TWI tool, as described above.

Whichever method routes the water, the specific catchment area is the same formula: the accumulation plus one, times the cell size. It is tempting to think that a cell sending water to two neighbors should be divided by two cell widths, but the divisor is the width of the cell whose value is being computed, not the number of cells it sends water to; Tarboton (1997) defines the specific catchment area as A/L with L the pixel width for every routing method, and his cone test uses it. The spreading is already expressed in the smaller shares each downslope cell receives. Divide again by the number of receivers and the values on a cone come out about half of the truth. Quinn et al. (1991) did set it up differently in their own paper: they divided the upslope area by the sum of the contour lengths of all the downhill directions, and used a weighted mean of the downhill gradients as the slope, so their divisor grows with the number of directions a cell drains in. This tool follows Tarboton's convention instead, which is also what Esri's Flow Accumulation delivers, so that one specific-catchment-area formula serves all three methods. The table shows the tool's D8 and D-Infinity results on a cone with unit cells, where the true specific catchment area on the ring at radius r is r/2 for every cell of the ring:

Ring radiusTrue value r/2D8: min, median, maxD-Infinity: min, median, max
1051, 5.0, 103.8, 4.9, 6.8
20101, 8.5, 207.3, 9.2, 12.9
40201, 18, 4014.1, 17.9, 25.1
60301, 25, 6021.0, 26.7, 37.3

D8 shows the thread problem. Around any ring its values scatter from 1, on the cells the threads bypass, up to double the truth on the eight cells where the cardinal and diagonal rays run; the median sits near the true value, but no single cell can be trusted. D-Infinity, divided by one cell width, stays near the true value all the way around.

The raster's edge and every NoData cell act as outlets. The fill starts from the cells that touch NoData or the edge, and a flat beside NoData drains into it, so water that reaches either one leaves the analysis. Two consequences follow. A DEM clipped across a drainage loses the catchment of every stream that enters from outside the clip, and the wetness index along those streams comes out too low; use a DEM that covers the full upstream watershed of the area of interest, and clip the output instead. And a lake or any other area set to NoData inside the DEM swallows the water that flows into it, cutting off every cell downstream, which is why the section above recommends masking lakes out of the output, not the DEM.

Reading the map

The output is drawn with a tan-to-blue ramp, dry to wet. The high, blue values trace the drainage network and pool in every basin, meadow and flat; the low, tan values sit on the ridges, the noses of spurs and the steep faces. In the Flagstaff DEM the Watershed lab works on, the high values occur in both the meadows and the drainage bottoms, which is what the index was built to find: the places on a landscape where water arrives and stays. If those features matter to a species you manage, this raster is a cheap and repeatable way to map them from elevation alone. The lecture suggests two other uses: finding locations vulnerable to flooding, and finding places where water may collect, tinajas and natural catchments in the desert.

The wetness index draped over a hillshade of the area north of Flagstaff, Arizona, around Schultz Pass and Doney Park: tan hills and ridges, a broad pale valley floor, and a blue network of washes and a braided main drainage scoring as the wettest ground on the map
The index over the topography just north of Flagstaff, Arizona, showing the Schultz Pass and Doney Park areas. The broad valley floor and its washes are the wettest ground on the map, because the terrain would collect water there. Whether any water reaches them is a different question, answered below.

Keep in mind what the index leaves out. It is built entirely from the surface topography, so it probably describes actual soil wetness best in mesic areas that receive steady rain, where the water the terrain would collect is actually there to be collected. Beven and Kirkby framed the model it comes from as one for humid temperate areas. The index pays no attention to aquifers or subsurface structure, and it is normal for xeric or desert areas to have enormously large drainages, be flat as a table, and still be dry as a bone: even substantial rainfall upstream can drain entirely into the ground before it reaches some downstream areas, and those areas will score near the top of the map regardless. In the desert Southwest the water that matters ecologically is more often associated with springs and exposed aquifers than with the wetness index, and a spring on a hillside will never show up in it. Use the index as a map of where the terrain would put water, and let soils and hydrology say whether the water is there.

Comparing wetness indices from different sources

Published TWI values are not comparable across studies unless the recipes match, and the recipes rarely do. The table lists the constant offsets that different formulas introduce; each shifts every cell by the same amount, so the maps look the same and the rankings agree, but the numbers differ.

Recipe differenceEffect on every cellSize, for a 30 m DEM
Catchment area used in place of specific catchment area (flow accumulation multiplied by the cell area instead of the cell size)adds ln(cell size)+3.4
Percent slope used in place of the gradient (tan β)subtracts ln(100)−4.6
Both of the above togetheradds ln(cell size) − ln(100)−1.2
Distances in kilometers instead of meterssubtracts ln(1000)−6.9

Two other differences are not constant offsets and change the pattern itself. The flow-routing method, D8 against MFD or D-Infinity, redistributes the upslope area on every hillside, as the hydrology section describes. And the DEM's resolution changes the catchment counts: Yang et al. (2011) show that the upslope contributing area at a given location grows as the DEM becomes coarser, most strongly in upslope positions and least near the bottoms of drainages, so a wetness index from a 10 m DEM and one from a 90 m DEM are different surfaces, not the same surface at two levels of detail. Within one study, use one DEM, one routing method and one formula, and compare values only with each other.

A tour of the dialog

The examples use a DEM of the country just north of Flagstaff, Arizona, with 28.08 m cells in NAD 1983 UTM Zone 12N: the eastern flank of the San Francisco Peaks on the left, US 89 running north through the middle, the cinder cones and the broad, gently sloping plain of Doney Park on the right, and the washes that gather it all toward the San Francisco Wash in the southeast corner. It is the Schultz Pass and Doney Park country of the figure above.

The Topographic Analysis Tools gallery open on the ribbon, with the Topographic Wetness Index button, in the Wetness, Slope and Extreme Values row, outlined in blue
Where to find it: Topographic Wetness Index is in the Wetness, Slope and Extreme Values row of the Topographic Analysis Tools gallery, in the Topographic Analysis group of the Wildlife and Forestry tab.
A reference map of the Flagstaff area: national forest in green over hillshaded terrain on the left, the tan Doney Park area and the town of Flagstaff on the right and bottom, US 89 and Interstate 40 in orange and blue, mapped stream paths in light blue, and a 5 kilometer scale bar
The study area, with the mapped stream paths in light blue and a 5 km scale bar. The forested slopes of the Peaks are on the left; Doney Park and its cinder cones on the right; Flagstaff along the bottom.
The Topographic Wetness Index geoprocessing pane with both optional groups expanded: input DEM_Flagstaff_Area, output Flagstaff_Area_TWI_DINF, Slope method Geodesic, Elevation units Meters, Value added to slope 0.001, Flow routing method DINF (D-Infinity; Tarboton 1997), empty Existing filled DEM and Existing flow-accumulation raster rows, and empty Output filled DEM and Output flow-accumulation raster rows
The pane for the D-Infinity run, with both optional groups opened to show every row. The D8 and MFD runs below differ only in the flow routing method.

The dialog asks for the DEM, an output name, the slope method for projected DEMs, the value added to slope, and the flow routing method, which is greyed out when you supply an accumulation raster of your own. An Elevation units row is filled in and locked when the DEM carries a vertical coordinate system; when it does not, the row is editable and shows the units from the last run, so check it. The slope method, the constant and the routing method are also remembered from the last run. Below those sit two collapsed groups. Use existing intermediate rasters takes a filled DEM, a flow-accumulation raster, or both; either must share the input DEM's grid exactly, the same extent, cell size and alignment, and supplying the accumulation raster makes the filled DEM unnecessary. Save computed intermediate rasters writes the filled DEM and the flow accumulation the tool computes, for reuse or inspection. Each save option is offered only when the tool actually computes that raster; supply a flow-accumulation raster and both save rows are disabled, because there is nothing of the tool's own to save.

The same ground under the three flow direction methods

The Topographic Wetness Index using the D8 method over a hillshade of the Flagstaff area, 2.81 to 23.18: thin blue channel lines threading the hills, and closely spaced parallel streaks running across the gently sloping plain of Doney Park
D8, the default: the index runs from 2.81 to 23.18. Water is confined to narrow threads, and on the gently sloping plain of Doney Park the threads run parallel across the ground like the teeth of a comb.
The Topographic Wetness Index using the D-Infinity method over the same area, 3.06 to 23.14: the same channels, but the ground around them grades smoothly from tan to blue, with broad soft aprons of wetter ground below the slopes and across the plain
D-Infinity: 3.06 to 23.14. The channels are in the same places, but the hillsides now grade into them, and the plain reads as a broad, soft apron of wetter ground rather than a set of parallel lines.
The Topographic Wetness Index using the MFD method over the same area, 3.10 to 23.11: the smoothest of the three, with the washes widening into broad blue fans and the gradation from hillside to channel softer still
MFD: 3.10 to 23.11. The smoothest of the three: the washes below Doney Park widen into broad fans, and the gradation from hillside into channel is softer than under D-Infinity.

Look first at the plain east of the Peaks, the pale ground between the forest and Doney Park. Under D8 it is covered with closely spaced parallel streaks. That is the thread problem of the cone on a real landscape: every cell sends all of its water to one neighbor, so on a smooth slope the flow lines never merge and never spread, and the index draws each one as a separate hairline of wetness with dry ground between. Under D-Infinity and MFD the same ground grades naturally from the drier hillsides down into the washes, with wetter aprons below the slopes and broad fans where the washes open out onto the plain, which is how water actually behaves on ground like this. D8 is much more constrained into narrow channels; the two dispersive methods spread the contributing area sideways across the slope, and the difference is what the eye picks up first. The main washes, where flow converges, look nearly the same in all three, and the maxima agree to within a tenth, 23.18, 23.14 and 23.11, because the total area draining through a channel does not depend on how the hillsides delivered it. The minima rise a little from D8 to the dispersive methods, 2.81 to 3.06 and 3.10, because under MFD and D-Infinity even a ridge cell receives some share from its neighbors, so no cell is left with only its own area.

The tool reports its progress through the fill, the flow direction, the accumulation and the slope. The messages record the formula and the constant added to the slope, and the output raster's metadata carries the same provenance. The output is a floating-point raster with statistics and a histogram already computed, NoData wherever the DEM is NoData, and the dry-to-wet ramp applied by default. Outputs may go into a geodatabase or a folder; in a folder, a name without an extension becomes a GeoTIFF.

The tool holds the whole DEM in memory at once, so the memory it needs grows with the number of cells. On most DEMs that is no concern. On a very large one the tool may need more memory than your computer has free, and then one of two things happens: Windows starts using the disk as overflow memory and the tool slows to a crawl, or the tool stops with an out-of-memory error. There is no fixed limit; it depends on how much memory your computer has free. If a DEM is too large, clip it to the area you need first.

ModelBuilder

A ModelBuilder diagram: DEM_Flagstaff_Area feeding Topographic Wetness Index, producing Flagstaff_Area_TWI_DINF, Output filled DEM and Output flow-accumulation raster
The DEM in; the wetness index out, here from the D-Infinity run, with the two optional intermediate rasters saved as well.

Parameters

LabelExplanationData type
Input elevation raster (single band)Required · in_raster The DEM. Must be single band; projected or geographic. Unless you supply a flow-accumulation raster below, this DEM is filled and flow-routed internally. The slope term always comes from this original, unfilled DEM. Raster Layer
Output TWI rasterRequired · out_raster The wetness index, ln(a / (tan β + the added constant)), as a floating-point raster. Higher is wetter. NoData outside the valid DEM area. A name is suggested from the input and a dry-to-wet color ramp is applied. Raster Dataset
Slope method (projected DEMs)Optional · slope_method Planar (the default) measures slope in the projection's plane with the raster's cell size, matching Esri's Slope tool; Geodesic uses true ground distances on the spheroid. Disabled for geographic DEMs, which are always computed per row on the spheroid. String
Elevation unitsOptional · elev_units Meters or Feet. Filled in and locked from the DEM's vertical coordinate system when it has one; you are asked only when the DEM does not say. Keeps the gradient dimensionless whatever units the DEM uses. String
Value added to slope (avoids divide-by-zero)Required · slope_epsilon A small positive constant added to every cell's slope gradient before dividing. Default 0.001, the HYDRO1k value. Larger values pull level cells toward lower TWI; smaller values push them higher; sloping cells are unaffected. Double
Flow routing method (when the tool computes the flow accumulation)Optional · flow_method D8 (single flow direction), the default; MFD (multiple flow direction; Qin et al. 2007); or DINF (D-Infinity; Tarboton 1997), as described above. Disabled when you supply your own flow-accumulation raster, because the routing is then already done. String
Existing filled DEM (optional)Optional · in_filled_dem A depression-filled version of the input DEM to route flow on, so the tool skips its own fill. Must share the input DEM's grid. Ignored if a flow-accumulation raster is also supplied. In the group Use existing intermediate rasters. Raster Layer
Existing flow-accumulation raster (optional)Optional · in_flow_accumulation A pre-computed flow accumulation, upstream cell count, from Esri's Flow Accumulation (any flow-direction type) or from this tool. The tool then skips filling, flow direction and accumulation and uses it directly. Must share the input DEM's grid. Same group. Raster Layer
Output filled DEMOptional · out_filled_dem Saves the filled DEM the tool computes. Available only when the tool computes the fill. In the group Save computed intermediate rasters. Raster Dataset
Output flow-accumulation rasterOptional · out_flow_accumulation Saves the flow accumulation the tool computes, upslope area in cells excluding the cell itself, fractional under MFD and D-Infinity. Available only when the tool computes the accumulation. Same group. Raster Dataset

Python

import arcpy
arcpy.ImportToolbox(r"C:\path\to\JennessEnterprisesTools.pyt")  # your install path
# TWI from a DEM, letting the tool fill and flow-route internally,
# and keeping the flow accumulation it computes.
arcpy.jenness.TopographicWetnessIndex(
    in_raster=r"C:\Project\Elev.gdb\DEM",
    out_raster=r"C:\Project\Elev.gdb\DEM_TWI",
    slope_method="Planar", elev_units="Meters", slope_epsilon=0.001,
    flow_method="D8 (single flow direction)",
    out_flow_accumulation=r"C:\Project\Elev.gdb\DEM_FlowAcc")
# The same with D-Infinity routing. The other string is
# "MFD (multiple flow direction; Qin et al. 2007)".
arcpy.jenness.TopographicWetnessIndex(
    in_raster=r"C:\Project\Elev.gdb\DEM",
    out_raster=r"C:\Project\Elev.gdb\DEM_TWI_Dinf",
    slope_method="Planar", elev_units="Meters", slope_epsilon=0.001,
    flow_method="DINF (D-Infinity; Tarboton 1997)")
# TWI on an Esri flow accumulation (for example one computed with MFD).
arcpy.jenness.TopographicWetnessIndex(
    in_raster=r"C:\Project\Elev.gdb\DEM",
    out_raster=r"C:\Project\Elev.gdb\DEM_TWI_MFD",
    slope_method="Planar", elev_units="Meters", slope_epsilon=0.001,
    in_flow_accumulation=r"C:\Project\Elev.gdb\DEM_FlowAcc_MFD")

Recommended citation

Jenness, J. 2026. Topographic Wetness Index. Wildlife and Forestry Tools add-in for ArcGIS Pro, v. 1.99 (September 2026). Jenness Enterprises. Available at: https://github.com/JeffJenness/Wildlife_Tools. Please also cite Beven and Kirkby (1979) for the index.

Credits and references

By Jeff Jenness, Jenness Enterprises (www.jennessent.com). The index is Beven and Kirkby's; the sink filling is the priority-flood algorithm of Barnes, Lehman and Mulla; the flow routing follows O'Callaghan and Mark (D8), Tarboton (D-Infinity) and Quinn et al. with the exponent of Qin et al. (MFD). The explanation of specific catchment area follows my Topographic Roughness lecture and the Watershed lab exercises of my GIS course.

Licensing information

Works at every ArcGIS Pro license level (Basic, Standard, Advanced). No extension licenses are required; the sink filling, flow direction, flow accumulation and slope are all computed internally, without Spatial Analyst.