Tell me more ×
Geographic Information Systems Stack Exchange is a question and answer site for cartographers, geographers and GIS professionals. It's 100% free, no registration required.
up vote 0 down vote favorite

I am calculating Topographic Wetness Index for a study area using SRTM DEM (90 and 30m).

I understand there are many "flavors" of TWI (CTI). I have been using: TWI = ln(As +1 / tanB) where As = catchment area, B = slope (degrees/radians??)

When I run this function, the output contains negative values (using slope in degrees and radians). Are negative values acceptable? Or does this indicate an error in the analysis??

I understand the output should be scaled to 1-10, but negative values seems counter intuitive. Could someone explain what negative values could mean (due to flat areas, no data values?)

On a side not, the calculations are being run on on a very large area, which overaps several watersheds. Is it appropriate to run TWI over such a large extent?

share|improve this question
Where does this formula for TWI come from? As stated, there's nothing to prevent it being negative: all you need is to measure catchment area in (say) millions of square kilometers, so that As is tiny, and have a slope appreciably greater than 45 degrees (so that 1/tan(B) is less than 1). At a minimum, then, the formula has to assume some specific units of measurement for the catchment area and that the the area is not too small and the slope is not too steep. (BTW, if you use degrees where radians are intended, you are likely to get many undefined values for the TWI.) – whuber♦ Dec 5 '12 at 21:19
The original reference for the CTI is: Moore, ID., P.E. Gessler, G.A. Nielsen, and G.A. Petersen (1993) Terrain attributes: estimation methods and scale effects. In Modeling Change in Environmental Systems, edited by A.J. Jakeman M.B. Beck and M. McAleer Wiley , London, pp. 189 - 214. Another common one that is in a journal: Gessler, P.E., I.D. Moore, N.J. McKenzie, and P.J. Ryan. (1995). Soil-landscape modeling and spatial prediction of soil attributes. International Journal of GIS. 9(4):421-432. – Jeffrey Evans Dec 6 '12 at 18:59

1 Answer

up vote 1 down vote

You should not be seeing negative values in the CTI. Since you did not provide a reproducible example I cannot speculate as to why you are getting incorrect results. The expected range is not limited 1-10. The range will be defined by flow accumulation which is influenced by the size of the basins that are accumulating flow. The index does not rely on washed boundaries so the extent of the DEM should not matter.

The common implementation of CTI can be shown as: CTI = ln(a / (tan (beta)) where; a=Upland contributing area [(flow accumulation + 1 ) * (cellsize)]and beta=slope in radians.

Here are the raster algebra steps to calculate CTI (you will have to translate syntax to your version of ArcGIS):

cellsize=30
fd = flowdirection(dem)
sca = flowaccumulation(fd)
slope = ( slope(dem) * 1.570796 ) / 90
tan_slp = con( slope > 0, tan(slope), 0.001 )
sca_scaled = ( sca + 1 ) * cellsize
cti = ln ( sca_scaled / tan_slp )
share|improve this answer
This looks good overall, but I am puzzled by the formula you quote: if I measure cellsize in kilometers, say (or worse, in decimal degrees), then a / tan(beta) can easily be less than 1, leading to a negative CTI. Perhaps your formula assumes cellsize is given in some standard, fixed units such as meters? – whuber♦ Dec 5 '12 at 21:22
Sorry, I should have specified that the assumed projection units are meters or feet. – Jeffrey Evans Dec 5 '12 at 23:30
Thanks for the reply, everyone. Jeffrey, I found the Geomorphometry Arc toolbox online and some of your other blogs. I will investigate further today. Thanks again. – SoilSciGuy Dec 6 '12 at 18:18
One more question regarding the ArcToolbox! In the python script, flow accumulation is as follows: "corr_sca = ( sca +1 ) * 30 " but in the help file it states cell size should be in sq. meters: As = (flow accumulation + 1) * (pixel area in m2) So, should cell size be area, such as square meters (30x30 m), or spatial resolution (30m)? Also, will the function automatically recognize different cell sizes, or should the cell size be changed in the script before running? – SoilSciGuy Dec 6 '12 at 19:10
Thanks for catching this. My programmer was suppose to add this as an toolbox argument and then calculate area in the code. I will get a new version of the toolbox posted in the next day or so. – Jeffrey Evans Dec 6 '12 at 19:37
show 3 more comments

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.