module Geocoder::Calculations

def bearing_between(point1, point2, options = {})
  # set default options
  options[:method] ||= Geocoder::Configuration.distances
  options[:method] = :linear unless options[:method] == :spherical
  # convert to coordinate arrays
  point1 = extract_coordinates(point1)
  point2 = extract_coordinates(point2)
  # convert degrees to radians
  point1 = to_radians(point1)
  point2 = to_radians(point2)
  # compute deltas
  dlat = point2[0] - point1[0]
  dlon = point2[1] - point1[1]
  case options[:method]
  when :linear
    y = dlon
    x = dlat
  when :spherical
    y = Math.sin(dlon) * Math.cos(point2[0])
    x = Math.cos(point1[0]) * Math.sin(point2[0]) -
        Math.sin(point1[0]) * Math.cos(point2[0]) * 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(bearing) + 360) % 360
end

def bounding_box(point, radius, options = {})
  lat,lon = extract_coordinates(point)
  radius  = radius.to_f
  units   = options[:units] || Geocoder::Configuration.units
  [
    lat - (radius / latitude_degree_distance(units)),
    lon - (radius / longitude_degree_distance(lat, units)),
    lat + (radius / latitude_degree_distance(units)),
    lon + (radius / longitude_degree_distance(lat, units))
  ]
end

def compass_point(bearing, points = COMPASS_POINTS)
  seg_size = 360 / points.size
  points[((bearing + (seg_size / 2)) % 360) / seg_size]
end

def distance_between(point1, point2, options = {})
  # set default options
  options[:units] ||= Geocoder::Configuration.units
  # convert to coordinate arrays
  point1 = extract_coordinates(point1)
  point2 = extract_coordinates(point2)
  # convert degrees to radians
  point1 = to_radians(point1)
  point2 = to_radians(point2)
  # compute deltas
  dlat = point2[0] - point1[0]
  dlon = point2[1] - point1[1]
  a = (Math.sin(dlat / 2))**2 + Math.cos(point1[0]) *
      (Math.sin(dlon / 2))**2 * Math.cos(point2[0])
  c = 2 * Math.atan2( Math.sqrt(a), Math.sqrt(1-a))
  c * earth_radius(options[:units])
end

def distance_to_radians(distance, units = nil)
  units ||= Geocoder::Configuration.units
  distance.to_f / earth_radius(units)
end

def earth_radius(units = nil)
  units ||= Geocoder::Configuration.units
  units == :km ? EARTH_RADIUS : to_miles(EARTH_RADIUS)
end

def extract_coordinates(point)
  case point
  when Array
    if point.size == 2
      lat, lon = point
      if !lat.nil? && lat.respond_to?(:to_f) and
        !lon.nil? && lon.respond_to?(:to_f)
      then
        return [ lat.to_f, lon.to_f ]
      end
    end
  when String
    point = Geocoder.coordinates(point) and return point
  else
    if point.respond_to?(:to_coordinates)
      if Array === array = point.to_coordinates
        return extract_coordinates(array)
      end
    end
  end
  [ NAN, NAN ]
end

def geographic_center(points)
  # convert objects to [lat,lon] arrays and convert degrees to radians
  coords = points.map{ |p| to_radians(extract_coordinates(p)) }
  # convert to Cartesian coordinates
  x = []; y = []; z = []
  coords.each do |p|
    x << Math.cos(p[0]) * Math.cos(p[1])
    y << Math.cos(p[0]) * Math.sin(p[1])
    z << Math.sin(p[0])
  end
  # compute average coordinate values
  xa, ya, za = [x,y,z].map do |c|
    c.inject(0){ |tot,i| tot += i } / c.size.to_f
  end
  # convert back to latitude/longitude
  lon = Math.atan2(ya, xa)
  hyp = Math.sqrt(xa**2 + ya**2)
  lat = Math.atan2(za, hyp)
  # return answer in degrees
  to_degrees [lat, lon]
end

def km_in_mi
  KM_IN_MI
end

def latitude_degree_distance(units = nil)
  units ||= Geocoder::Configuration.units
  2 * Math::PI * earth_radius(units) / 360
end

def longitude_degree_distance(latitude, units = nil)
  units ||= Geocoder::Configuration.units
  latitude_degree_distance(units) * Math.cos(to_radians(latitude))
end

def mi_in_km
  1.0 / KM_IN_MI
end

def radians_to_distance(radians, units = nil)
  units ||= Geocoder::Configuration.units
  radians * earth_radius(units)
end

def to_degrees(*args)
  args = args.first if args.first.is_a?(Array)
  if args.size == 1
    (args.first * 180.0) / Math::PI
  else
    args.map{ |i| to_degrees(i) }
  end
end

def to_kilometers(mi)
  mi * mi_in_km
end

def to_miles(km)
  km * km_in_mi
end

def to_radians(*args)
  args = args.first if args.first.is_a?(Array)
  if args.size == 1
    args.first * (Math::PI / 180)
  else
    args.map{ |i| to_radians(i) }
  end
end