I have downloaded SRTM GDEM (~90 km resolution). I am using ArcGIS 10. I have tried to use spatial analyst to compute for slope. However, I can not compute for the slope. The output values has only two ranges 0 and 0.1-90. I am not really sure what is the problem?
I'm guessing the horizontal units for your raster are in either degrees or arcseconds. You need to reproject this raster to a spatial projection where your horizontal and vertical units are the same (i.e., if the vertical units are in metres, then I suggest using UTM, which has horizontal units of metres).
To reproject a raster with ArcCatalog/ArcGIS, look in:
Choose a projected spatial reference that covers your region of interest, e.g., try a UTM zone. There are many other options which are best documented in the manual.
Better answer, using GDAL with a scale
Now that SRTM data are globally available I can actually see and work with the files. The
This seems like a good place to describe a simple, fast, and more than reasonably accurate way to compute slopes for a globally extensive DEM.
Recall that the slope of a surface at a point is essentially the largest ratio of "rise" to "run" encountered at all possible bearings from that point. The issue is that when a projection has scale distortion, the values of "run" will be incorrectly computed. Even worse, when the scale distortion varies with bearing--which is the case with all projections that are not conformal--how the slope varies with bearing will be incorrectly estimated, preventing accurate identification of the maximum rise:run ratio (and skewing the calculation of the aspect).
We can solve this by using a conformal projection to ensure that the scale distortion does not vary with bearing, and then correcting the slope estimates to account for the scale distortion (which varies from point to point throughout the map). The trick is to use a global conformal projection that allows a simple expression for its scale distortion.
The Mercator projection fits the bill: assuming scale is correct at the Equator, its distortion equals the secant of the latitude. That is, distances on the map appear to be multiplied by the secant. This causes any slope calculation to compute rise:(sec(f)*run) (which is a ratio), where f is the latitude. To correct this, we need to multiply the computed slopes by sec(f); or, equivalently, divide them by cos(f). This gives us the simple recipe:
To do this with a grid given in decimal degrees (such as an SRTM DEM), perform the following steps:
The errors in the slope calculations will be up to 0.3% (because this procedure uses a spherical earth model rather than an ellipsoidal one, which is flattened by 0.3%). That error is substantially smaller than other errors that go into slope calculations and so can be neglected.
Fully global calculations
The Mercator projection cannot handle either pole. For work in polar regions, consider using a polar Stereographic projection with true scale at the pole. The scale distortion equals 2 / (1 + sin(f)). Use this expression in place of sec(f) in the workflow. Specifically, instead of computing a cosine(latitude) grid, compute a grid whose values are (1 + sin(latitude))/2. Then proceed exactly as before.
For a complete global solution, consider breaking the terrestrial grid into three parts--one around each pole and one around the equator--, performing a slope calculation separately in each part using a suitable projection, and mosaicing the results. A reasonable place to split the globe is along circles of latitude at latitudes of 2*ArcTan(1/3), which is about 37 degrees, because at these latitudes the Mercator and Stereographic correction factors are equal to each other (having a common value of 5/4) and it would be nice to minimize the sizes of the corrections made. As a check of the computations, the grids should be in very close agreement where they overlap (tiny amounts of floating point imprecision and differences due to resampling of the projected grids ought to be the only sources of discrepancies).
John P. Snyder, Map Projections--A Working Manual. USGS Professional Paper 1395, 1987.
Simply put, there isn't one. By definition a coordinate system based on degrees is un-projected. In common parlance we say WGS84 is a "geographic" projection, but that's untrue, just for convenience.
I think I remember reading about a software or process for accurately working with elevation models in un-projected geographic space but I can't locate it right now. In any case it would have been an experimental or build it yourself from code kind of process.
Ahhh, found it: Development of a Global Slope Dataset for Estimation of Landslide Occurrence Resulting from Earthquakes (USGS). Page 4 describes the problem well
The paper goes on to describe the specific calculations and software tools (gdal, python, numpy) they used to workaround this fundamental issue. The paper doesn't include the code, but if asked nicely they might share. In any case though I'd probably just ask where the results are, being the USGS it's probably already online somewhere. :)
Slope is rise / run. Compute rise and compute run and you have your answer. It is simple to compute the distance between geographic coordinates. This will introduce less resampling error as compared to conversion to UTM, etc.