Skip to content

Instantly share code, notes, and snippets.

@amitkgupta
Created February 23, 2013 09:55
Show Gist options
  • Save amitkgupta/5019163 to your computer and use it in GitHub Desktop.
Save amitkgupta/5019163 to your computer and use it in GitHub Desktop.
Find a point on the sphere minimizing the average distance to a given set of points, within a given degree of precision.
class Point
attr_reader :lat, :long
def initialize(lat, long)
@lat = lat.to_f
@long = long.to_f
if (-Math::PI/2 > @lat || @lat > Math::PI/2)
raise "lat must between -PI/2 and PI/2, inclusive"
end
if (0 > @long || @long >= Math::PI*2)
raise "long must be between 0 and 2*PI (including 0 but not 2*PI)"
end
end
def distance(point)
Math.acos(Math.sin(lat)*Math.sin(point.lat) + Math.cos(lat)*Math.cos(point.lat)*Math.cos((long-point.long).abs))
end
end
###########################
class Algorithm
def initialize(points, precision)
@points = points
@precision = precision
end
def compute
distance = Float::INFINITY
lat = -Math::PI/2
long = 0
point = nil
while (lat <= Math::PI/2) do
while (long < 2*Math::PI) do
new_point = Point.new(lat, long)
new_distance = (@points.map{|pt| pt.distance(new_point)}).inject(:+)
if new_distance < distance
point = new_point
distance = new_distance
end
long = long + @precision
end
long = 0
lat = lat + @precision
end
point
end
end
###########################
point_1 = Point.new(-Math::PI/2, 0)
point_2 = Point.new(0, 0)
point_3 = Point.new(0, Math::PI/2)
point_4 = Point.new(-1.5, 1)
point_5 = Point.new(0.1, 4)
point_6 = Point.new(1.6, 6)
point_7 = Point.new(-0.4, 2)
points = [point_1, point_2, point_3, point_4, point_5, point_6, point_7]
precision = 1.0/360
# This will give you the latitude and longitude of "the" point minimizing
# the average distance to the 7 points above, good to within one 360th
# of a radian. This means it tests about 2.5 million points.
p Algorithm.new(points, precision).compute
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment