Geographic Information Systems Asked on October 8, 2020
I have a point data set with terrain aspect at the location of the points, measured in degrees. I want to compute the similarity of aspect between the points, e.g. 70° is similar to 90°. This is straight forward, but the problem arises for aspects near the north direction. A aspect of 350° should be similar to 10°. How can I detect these two aspects to be similar?
I am looking for a transformation of the values that takes into account that 358° and 2° is almost the same, even though it is far apart in my 0-360 metric.
Compute the difference, and then take 360-the absolute difference if the absolute difference is greater than 180.
Here's one way of doing that in an R function:
d180 = function(a,b){
d = a-b
d = ifelse(abs(d)<180 , d , 360-abs(d))
d
}
and some examples:
> d180(0,0)
[1] 0
> d180(0,1)
[1] -1
> d180(0,180)
[1] 180
> d180(0,181)
[1] 179
> d180(358,2)
[1] 4
> d180(2,358)
[1] 4
Whether the distance is positive or negative is a bit odd, return the abs(d)
if you don't care...
Answered by Spacedman on October 8, 2020
A bit of a different perspective, you could just linearize aspect. Here is pseudo-code pulled from the linear aspect tool in our ArcGIS Geomorphometry & Gradient Metrics toolbox. Note; we use (mod (..) * 100)
and (36000(360*100) / 100)
allow for two decimal points since Fmod
is no longer available in ArcGIS.
a = aspect(elev) # degrees
t1 = focal(sin(a/57.296), window=3,"sum")
t2 = focal(cos(a/57.296), window=3,"sum")
linear.aspect = mod(((450-(atan2(t1, t2) * 57.296)) * 100), 36000) / 100
Another nice trick is to look at the Stage (1978) slope-aspect transformation that provides an interaction. Assuming slope in degrees (Stage paper assumes percent) and aspect in degrees, you can calculate it thus:
s = (slope(elev) / 0.572957795786) * 0.01
a = aspect(elev) * (pi / 180)
scosa = s * cos(a)
scosa = s * sin(a)
Example transformed values for 50% slope across 10 aspects.
Aspect cosine sine
N 0.500 0.000
N30E 0.433 0.250
N45E 0.345 0.345
N60E 0.250 0.433
E 0.000 0.500
ESE -0.354 0.354
S -0.500 0.000
SSW -0.354 -0.354
W 0.000 -0.500
Or, linearize with a recentered origin. This is done in vegetation modeling to provide a metric that is optimized for site prodictivity (eg. Roberts and Cooper 1989 TRASP transformation) Assuming degrees, 1 - cos( (pi / 180) (a - 30)) / 2
a = aspect(elev)
trasp = (1 - cos( (3.142 / 180) * (a - 30)) ) / 2
These transformations are available in ArcGIS in the Geomorphometry & Gradient Metrics toolbox and in the R package spatialEco (on CRAN).
Answered by Jeffrey Evans on October 8, 2020
Retaining sign (pseudo):
a = modulus((targetangle-sourceangle),6.283185307179586)
if (a>3.141592653589793)
return a-6.283185307179586
else
return a
(Talking degrees:) If modulus(-45, 360) doesn't return 315 in your language, then use modulus((-45 + 360), 360). Your modulus fn must support decimal values.
Answered by Stormwind on October 8, 2020
Get help from others!
Recent Questions
Recent Answers
© 2024 TransWikia.com. All rights reserved. Sites we Love: PCI Database, UKBizDB, Menu Kuliner, Sharing RPP