4

I need to derive a raster representing aspect for North America (from STRM data, originally in geographic lat/long). In most projections, true north is only map north along eg the central meridian so aspect values are distorted. I attempted calculating aspect in Mercator and then reprojectiong to eg Lambert Conformal Conic, but the aspects values are very pixelated and distorted.

The following answer by whuber suggests how slope could be calculated accurately across a continent: "Compute the slope (as rise:run or a percent) using a Mercator projection, then divide the result by the cosine of the latitude". Using SRTM Global DEM for Slope calculation?

How would this be extended to calculate aspect?

Would one use flow direction based on the corrected slope?

Or can other methods be used?

I see in Compute global DEM topographic indices and compromise projections that @whuber states "Having done this [corrected slope calculation], perform the followup calculation of flow directions, aspects, etc. using the corrected slopes."

I'm just not sure how to derive aspect from a slope raster (vs from a DEM) in ArcGIS or otherwise.

Improve this question
1
  • You can't get an aspect from a slope, you need to use your DEM... but if you're happy with how the aspect is created your slope can be done the same way (except use Aspect instead of Slope). Note that slope isn't affected by true north issues in the same way as Aspect, there will be slight differences between calculating in geographic and projected (Lambers, Albers etc..) however that is mostly due to the resampling that is unavoidable when projecting. The SRTM DEM isn't fantastic (<=16 m height accuracy) so slight changes are likely to be less than the accuracy of the product using bilinear.
    – Michael Stimson
    Commented May 4, 2015 at 23:00

1 Answer 1

Reset to default
3

There are many possible solutions, but a convenient and accessible way will be to use an appropriate projection. What properties must it have?

  1. It should be conformal (or very nearly so throughout the region of interest). This means there is no relative scale distortion as the bearing is varied around any fixed point but--necessarily--there will be distortion of scale from one point to the next. However, this latter distortion is going to be so small over the neighborhoods in which aspect is computed that it will have negligible influence on the accuracy. The lack of relative local distortion means any aspect computed in the projected coordinates will not be affected by the projection. (The slope, however, might be seriously wrong, which is why for geographically extensive regions it often needs to be computed separately as described elsewhere.)

  2. It should be cylindrical. This means that true north always points in the same direction across the map. When using this common direction to designate the zero angle, aspects (which are computed in the projected coordinate system) will always be correct--or at least will have a constant offset from north, which is easily corrected afterwards.

The family of Mercator projections is conformal and cylindrical. (In fact, this more or less completely determines the Mercator projection, up to scale.) The solution for North America, then, which is covered by the Mercator, is to reproject the original DEM into Mercator and compute its aspect.

For those interested in working near poles, where the Mercator breaks down, it might be most convenient to use a Stereographic projection centered at the pole. Because this is not cylindrical, the aspect will have to be adjusted by the apparent direction of each meridian, which varies systematically from point to point in a simple way: the projected meridian always points directly towards where the pole is projected. For instance, suppose you were to use a Northern Stereographic projection with the Prime Meridian directed to the right. Then a calculated aspect of 270 degrees, which ordinarily would be interpreted as due West, actually corresponds to true North. To make this adjustment, 90 degrees has to be added to the aspect. At another location on the map, corresponding (say) to -90 degrees longitude, the meridian would point straight up to the pole: now, no adjustment is necessary. In general, in this situation, 90 degrees plus the longitude has to be added to the computed aspect (and then reduced modulo 360). But that's a simple sequence of grid operations that can be easily performed after an aspect computation.

For worldwide calculations, perform computations on overlapping areas covering the two poles and the near-Equatorial regions. Reproject those results to a common coordinate system and mosaic them. As a double-check, the aspects computed in the overlapping regions should be in close agreement throughout.

Because the projection used to compute the aspect might not be the one desired for further analysis, you might need to reproject the aspect grid. This must be done by computing a cosine and sine grid of the aspects, reprojecting them, and recombining those via an arctangent calculation. See https://gis.stackexchange.com/a/133709 for further discussion.

Improve this answer
6
  • Thanks for this info. The key point was in the linked question: "One of the simpler approaches is to compute two grids from the aspect: its cosine and its sine. Reproject those, then recompute the angle from the reprojected values (using an arctangent function). Watch the process carefully, because much GIS software is schizophrenic concerning trig functions: it tends to use degrees to report angular values (like aspects) but to insist on radians as input to trig functions. Spot-check the cosine and sine grids to make sure the calculations are correct."
    – ccarroll
    Commented May 5, 2015 at 19:36
  • I had a followup question: 1) Mercator projections are offered in ArcGIS based on both spheres and spheroids. STRM data is based on the WGS84 datum, so I would think a Mercator projection based on WGS84 would be preferable to a sphere-based one. But since the latitudinal correction to slope that you suggest is based on a sphere model, would it be better (less error) to use a sphere-based Mercator projection?
    – ccarroll
    Commented May 5, 2015 at 19:46
  • Good question. The correction I suggest is independent of the datum: it requires only that you be able to compute the longitude with reasonable accuracy. That can be done in numerous ways (such as by creating an x-coordinate grid in geographic coordinates and projecting it). In any event, the errors created by using (say) a spherical datum instead of WGS84 would be so small (typically less than 3 parts per thousand) that it is of little practical importance which datum is used. For work to be published it's still a good idea to use the correct datum :-).
    – whuber
    Commented May 5, 2015 at 19:50
  • One more question: I see the point in resampling sine and cosine values, but I cannot recreate the original angles from them. I had thought that since sin(x)/cos(x) = tan(x), then ATan(Sin(Aspect[in radians])/Cos(Aspect)) would do so. However, the tangent values repeat such that tan(x) = tan(x+180). Maybe I am misinterpreting the suggestion to use the arctangent?
    – ccarroll
    Commented May 6, 2015 at 20:32
  • 1
    A couple of minor notes to help others: One needs to add 2*pi to the output of atan2(sin(x), cos(x)) where it is negative: Con(z >= 0, z, z + (math.pi * 2)) Also I had to use bilinear interpolation rather than cubic, as the latter produced values outside the valid range of sine and cosine values (-1,1). I wanted to thank whuber for his assistance.
    – ccarroll
    Commented May 6, 2015 at 22:52

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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