Skip to content
Snippets Groups Projects
Commit 61f97d3e authored by Alex Reisner's avatar Alex Reisner
Browse files

Add :linear method for bearing calculation.

Also clean up code comments.
parent b98c0a36
No related branches found
No related tags found
No related merge requests found
......@@ -48,23 +48,41 @@ module Geocoder
# Calculate bearing between two sets of coordinates.
# Returns a number of degrees from due north (clockwise).
#
# Also accepts an options hash:
#
# * <tt>:method</tt> - <tt>:linear</tt> (default) or <tt>:spherical</tt>;
# the spherical method is "correct" in that it returns the shortest path
# (one along a great circle) but the linear method is the default as it
# is less confusing (returns due east or west when given two points with
# the same latitude)
#
# Based on: http://www.movable-type.co.uk/scripts/latlong.html
#
def bearing_between(lat1, lon1, lat2, lon2)
def bearing_between(lat1, lon1, lat2, lon2, options = {})
options[:method] = :linear unless options[:method] == :spherical
# convert degrees to radians
lat1, lon1, lat2, lon2 = to_radians(lat1, lon1, lat2, lon2)
# compute delta
# compute deltas
dlat = lat2 - lat1
dlon = lon2 - lon1
y = Math.sin(dlon) * Math.cos(lat2)
x = Math.cos(lat1) * Math.sin(lat2) - Math.sin(lat1) *
Math.cos(lat2) * Math.cos(dlon)
brng = Math.atan2(x,y)
# brng is in radians counterclockwise from due east.
case options[:method]
when :linear
y = dlon
x = dlat
when :spherical
y = Math.sin(dlon) * Math.cos(lat2)
x = Math.cos(lat1) * Math.sin(lat2) -
Math.sin(lat1) * Math.cos(lat2) * Math.cos(dlon)
end
bearing = Math.atan2(x,y)
# Answer is in radians counterclockwise from due east.
# Convert to degrees clockwise from due north:
(90 - to_degrees(brng) + 360) % 360
(90 - to_degrees(bearing) + 360) % 360
end
##
......
......@@ -27,9 +27,11 @@ module Geocoder::Orm
"OR #{geocoder_options[:longitude]} IS NULL"}}
##
# Find all objects within a radius (in miles) of the given location
# (address string). Location (the first argument) may be either a string
# to geocode or an array of coordinates (<tt>[lat,long]</tt>).
# Find all objects within a radius of the given location.
# Location may be either a string to geocode or an array of
# coordinates (<tt>[lat,lon]</tt>). Also takes an options hash
# (see Geocoder::Orm::ActiveRecord::ClassMethods.near_scope_options
# for details).
#
scope :near, lambda{ |location, *args|
latitude, longitude = location.is_a?(Array) ?
......@@ -56,6 +58,10 @@ module Geocoder::Orm
# * +:units+ - <tt>:mi</tt> (default) or <tt>:km</tt>; to be used
# for interpreting radius as well as the +distance+ attribute which
# is added to each found nearby object
# * +:bearing+ - <tt>:linear</tt> (default) or <tt>:spherical</tt>;
# the method to be used for calculating the bearing (direction)
# between the given point and each found nearby point;
# set to false for no bearing calculation
# * +:select+ - string with the SELECT SQL fragment (e.g. “id, name”)
# * +:order+ - column(s) for ORDER BY SQL clause
# * +:limit+ - number of records to return (for LIMIT SQL clause)
......@@ -76,7 +82,8 @@ module Geocoder::Orm
##
# Scope options hash for use with a database that supports POWER(),
# SQRT(), PI(), and trigonometric functions (SIN(), COS(), and ASIN()).
# SQRT(), PI(), and trigonometric functions SIN(), COS(), ASIN(),
# ATAN2(), DEGREES(), and RADIANS().
#
# Distance calculations based on the excellent tutorial at:
# http://www.scribd.com/doc/2569355/Geo-Distance-Search-with-MySQL
......@@ -87,15 +94,28 @@ module Geocoder::Orm
def full_near_scope_options(latitude, longitude, radius, options)
lat_attr = geocoder_options[:latitude]
lon_attr = geocoder_options[:longitude]
bearing = "(DEGREES(ATAN2( " +
"SIN(RADIANS(#{lon_attr} - #{longitude})) * " +
"COS(RADIANS(#{lat_attr})), (" +
"COS(RADIANS(#{latitude})) * SIN(RADIANS(#{lat_attr}))" +
") - (" +
"SIN(RADIANS(#{latitude})) * COS(RADIANS(#{lat_attr})) * " +
"COS(RADIANS(#{lon_attr} - #{longitude}))" +
"))) + 360)"
bearing = "CAST(#{bearing} AS decimal) % 360"
options[:bearing] = :linear unless options.include?(:bearing)
bearing = case options[:bearing]
when :linear
"CAST(" +
"DEGREES(ATAN2( " +
"RADIANS(#{lon_attr} - #{longitude}), " +
"RADIANS(#{lat_attr} - #{latitude})" +
")) + 360 " +
"AS decimal) % 360"
when :spherical
"CAST(" +
"DEGREES(ATAN2( " +
"SIN(RADIANS(#{lon_attr} - #{longitude})) * " +
"COS(RADIANS(#{lat_attr})), (" +
"COS(RADIANS(#{latitude})) * SIN(RADIANS(#{lat_attr}))" +
") - (" +
"SIN(RADIANS(#{latitude})) * COS(RADIANS(#{lat_attr})) * " +
"COS(RADIANS(#{lon_attr} - #{longitude}))" +
")" +
")) + 360 " +
"AS decimal) % 360"
end
distance = "#{Geocoder::Calculations.earth_radius} * 2 * ASIN(SQRT(" +
"POWER(SIN((#{latitude} - #{lat_attr}) * PI() / 180 / 2), 2) + " +
"COS(#{latitude} * PI() / 180) * COS(#{lat_attr} * PI() / 180) * " +
......@@ -103,7 +123,8 @@ module Geocoder::Orm
options[:order] ||= "#{distance} ASC"
default_near_scope_options(latitude, longitude, radius, options).merge(
:select => "#{options[:select] || '*'}, " +
"#{distance} AS distance, #{bearing} AS bearing",
"#{distance} AS distance" +
(bearing ? ", #{bearing} AS bearing" : ""),
:having => "#{distance} <= #{radius}"
)
end
......
......@@ -199,17 +199,17 @@ class GeocoderTest < Test::Unit::TestCase
:w => [40, -76]
}
directions = [:n, :e, :s, :w]
types = [:spherical]
methods = [:linear, :spherical]
types.each do |t|
methods.each do |m|
directions.each_with_index do |d,i|
opp = directions[(i + 2) % 4] # opposite direction
p1 = points[d]
p2 = points[opp]
b = Geocoder::Calculations.bearing_between(*(p1 + p2))
b = Geocoder::Calculations.bearing_between(*(p1 + p2), :method => m)
assert (b - bearings[opp]).abs < 1,
"Bearing (#{t}) should be close to #{bearings[opp]} but was #{b}."
"Bearing (#{m}) should be close to #{bearings[opp]} but was #{b}."
end
end
end
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment