import math


def haversine_bearing(lat1, lon1, lat2, lon2):
    """
    Calculate the bearing from 1 point to 1 other
    on the earth (specified in decimal degrees)
    """
    # convert decimal degrees to radians 
    lat1, lon1, lat2, lon2 = map(math.radians, [lat1, lon1, lat2, lon2])

    dlon = lon2 - lon1
 
    b = math.atan2(math.sin(dlon) * math.cos(lat2),
        math.cos(lat1) * math.sin(lat2)
        - math.sin(lat1) * math.cos(lat2) * math.cos(dlon)) # bearing calc

    bd = math.degrees(b)

    br, bn = divmod(bd + 360, 360) # the bearing remainder and final bearing
    
    return bn