class GeoDistance::Vincenty
Public Class Methods
distance(*args)
click to toggle source
Calculate the distance between two Locations using the Vincenty
formula
# File lib/geo-distance/formula/vincenty.rb, line 9 def self.distance *args begin from, to, units = get_points(args) lat1, lon1, lat2, lon2 = [from.lat, from.lng, to.lat, to.lng] from_longitude = lon1.to_radians from_latitude = lat1.to_radians to_longitude = lon2.to_radians to_latitude = lat2.to_radians earth_major_axis_radius = earth_major_axis_radius_map[:kilometers] earth_minor_axis_radius = earth_minor_axis_radius_map[:kilometers] f = (earth_major_axis_radius - earth_minor_axis_radius) / earth_major_axis_radius l = to_longitude - from_longitude u1 = atan((1-f) * tan(from_latitude)) u2 = atan((1-f) * tan(to_latitude)) sin_u1 = sin(u1) cos_u1 = cos(u1) sin_u2 = sin(u2) cos_u2 = cos(u2) lambda = l lambda_p = 2 * Math::PI iteration_limit = 20 while (lambda-lambda_p).abs > 1e-12 && (iteration_limit -= 1) > 0 sin_lambda = sin(lambda) cos_lambda = cos(lambda) sin_sigma = sqrt((cos_u2*sin_lambda) * (cos_u2*sin_lambda) + (cos_u1*sin_u2-sin_u1*cos_u2*cos_lambda) * (cos_u1*sin_u2-sin_u1*cos_u2*cos_lambda)) return 0 if sin_sigma == 0 # co-incident points cos_sigma = sin_u1*sin_u2 + cos_u1*cos_u2*cos_lambda sigma = atan2(sin_sigma, cos_sigma) sin_alpha = cos_u1 * cos_u2 * sin_lambda / sin_sigma cosSqAlpha = 1 - sin_alpha*sin_alpha cos2SigmaM = cos_sigma - 2*sin_u1*sin_u2/cosSqAlpha cos2SigmaM = 0 if cos2SigmaM.nan? # equatorial line: cosSqAlpha=0 (ยง6) c = f/16*cosSqAlpha*(4+f*(4-3*cosSqAlpha)) lambda_p = lambda lambda = l + (1-c) * f * sin_alpha * (sigma + c*sin_sigma*(cos2SigmaM+c*cos_sigma*(-1+2*cos2SigmaM*cos2SigmaM))) end # formula failed to converge (happens on antipodal points) # We'll call Haversine formula instead. return Haversine.distance(from, to, units) if iteration_limit == 0 uSq = cosSqAlpha * (earth_major_axis_radius**2 - earth_minor_axis_radius**2) / (earth_minor_axis_radius**2) a = 1 + uSq/16384*(4096+uSq*(-768+uSq*(320-175*uSq))) b = uSq/1024 * (256+uSq*(-128+uSq*(74-47*uSq))) delta_sigma = b*sin_sigma*(cos2SigmaM+b/4*(cos_sigma*(-1+2*cos2SigmaM*cos2SigmaM)- b/6*cos2SigmaM*(-3+4*sin_sigma*sin_sigma)*(-3+4*cos2SigmaM*cos2SigmaM))) c = earth_minor_axis_radius * a * (sigma-delta_sigma) c = (c / unkilometer).to_deg units ? c.radians_to(units) : c rescue Errno::EDOM 0.0 end end
Private Class Methods
unkilometer()
click to toggle source
# File lib/geo-distance/formula/vincenty.rb, line 76 def self.unkilometer 6378.135 end