2

I have a ton of consecutive points and I need to know the bearing between each of them (a to b, b to c, c to d,etc). Ultimately, I will convert the numbers (bearings, azimuths?) to directions: N, NE, NNE, NNW,etc.

I have found a python script for use in the field calculator that comes close to doing this but it returns angles: https://geonet.esri.com/blogs/dan_patterson/2016/09/01/distance-calculations-using-the-field-calculator

Another option would be a script that can convert the angles returned by the above script to bearings/azimuths. For instance, the above script returns -90 rather than 270 for due west.

4
  • 1
    What GIS software are you using? Have you attempted to modify the script to output what you require? Please edit your question to include a snippet of what you've tried, and where you're getting stuck.
    – Midavalo
    Commented Oct 7, 2016 at 18:51
  • I am using ArcGIS 10.4. No, I haven't tried modifying anything because I do not know Python or scripting in general. Dan Patterson has a script for returning directions from azimuths but he doesn't have a script for calculating azimuths between consecutive points.
    – Bee Cee
    Commented Oct 7, 2016 at 19:56
  • You can always add/subtract 360 to change the sign. From there, it's simple to convert to cardinal directions.
    – Paul
    Commented Oct 7, 2016 at 23:40
  • Try tool Generate near table
    – Bera
    Commented Nov 14, 2017 at 14:36

1 Answer 1

1

This will do it on a plane surface

import math
from itertools import cycle
def angle_wrt_x(B,A):
    ax, ay = A
    bx, by = B
    b = math.degrees(math.atan2(ax-bx,ay-by))
    if b < 0:
        b+=360
    elif b >= 360:
        b-=360
    return b

firstPoint = (0,0)
secondPoint = (1,1)

coord_pairs = [(0,0), (0,1), (1,1), (1,0)]

listcycle = cycle(coord_pairs)

for firstPoint in coord_pairs:
    secondPoint = listcycle.next()
    print firstPoint,secondPoint, angle_wrt_x(firstPoint,secondPoint)

Returns

(0, 0) (0, 1) 0.0
(0, 1) (1, 1) 90.0
(1, 1) (1, 0) 180.0
(1, 0) (0, 0) 270.0

Add this def, you'll need to finish it off for all other regions on the compass

def compassPoints(b):
    if b > 348.75 or b <= 11.25:
        compB = 'N'
    elif 11.25 < b <= 33.75:
        compB ='NNE'
    elif 33.75 < b <= 56.25:
        compB = 'NE'
    #continue this pattern for the rest of them....

   return compB

then call it with

for firstPoint in coord_pairs:
    secondPoint = listcycle.next()
    b = angle_wrt_x(firstPoint,secondPoint)
    print compassPoints(b)
1
  • To use it I need it to return NE, SE, NW, SW, etc and also NNW, NNE, SSE, SSW etc. Thanks for your help
    – Bee Cee
    Commented Oct 8, 2016 at 11:46

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.