按距离过滤Python地理编码

Question

我需要将地理编码过滤到某个位置。例如，我想筛选餐馆地理编码列表以识别距我目前位置10英里内的餐馆。

有人可以指向我的函数，将距离转换成纬度&经度德尔塔？例如：

class GeoCode(object): 
    """Simple class to store geocode as lat, lng attributes.""" 
    def __init__(self, lat=0, lng=0, tag=None): 
     self.lat = lat 
     self.lng = lng 
     self.tag = None 

def distance_to_deltas(geocode, max_distance): 
    """Given a geocode and a distance, provides dlat, dlng 
     such that 

     |geocode.lat - dlat| <= max_distance 
     |geocode.lng - dlng| <= max_distance 
    """ 
    # implementation 
    # uses inverse Haversine, or other function? 
    return dlat, dlng 

注意：我正在使用距离的上等规范。

Answer 1

似乎不是都得到了良好的Python实现。幸运的是，“相关文章”边栏是我们的朋友。 This SO article指向给数学和Java实现的excellent article。您需要的实际功能相当短，并嵌入在下面的Python代码中。测试程度显示。在评论中阅读警告。

from math import sin, cos, asin, sqrt, degrees, radians 

Earth_radius_km = 6371.0 
RADIUS = Earth_radius_km 

def haversine(angle_radians): 
    return sin(angle_radians/2.0) ** 2 

def inverse_haversine(h): 
    return 2 * asin(sqrt(h)) # radians 

def distance_between_points(lat1, lon1, lat2, lon2): 
    # all args are in degrees 
    # WARNING: loss of absolute precision when points are near-antipodal 
    lat1 = radians(lat1) 
    lat2 = radians(lat2) 
    dlat = lat2 - lat1 
    dlon = radians(lon2 - lon1) 
    h = haversine(dlat) + cos(lat1) * cos(lat2) * haversine(dlon) 
    return RADIUS * inverse_haversine(h) 

def bounding_box(lat, lon, distance): 
    # Input and output lats/longs are in degrees. 
    # Distance arg must be in same units as RADIUS. 
    # Returns (dlat, dlon) such that 
    # no points outside lat +/- dlat or outside lon +/- dlon 
    # are <= "distance" from the (lat, lon) point. 
    # Derived from: http://janmatuschek.de/LatitudeLongitudeBoundingCoordinates 
    # WARNING: problems if North/South Pole is in circle of interest 
    # WARNING: problems if longitude meridian +/-180 degrees intersects circle of interest 
    # See quoted article for how to detect and overcome the above problems. 
    # Note: the result is independent of the longitude of the central point, so the 
    # "lon" arg is not used. 
    dlat = distance/RADIUS 
    dlon = asin(sin(dlat)/cos(radians(lat))) 
    return degrees(dlat), degrees(dlon) 

if __name__ == "__main__": 

    # Examples from Jan Matuschek's article 

    def test(lat, lon, dist): 
     print "test bounding box", lat, lon, dist 
     dlat, dlon = bounding_box(lat, lon, dist) 
     print "dlat, dlon degrees", dlat, dlon 
     print "lat min/max rads", map(radians, (lat - dlat, lat + dlat)) 
     print "lon min/max rads", map(radians, (lon - dlon, lon + dlon)) 

    print "liberty to eiffel" 
    print distance_between_points(40.6892, -74.0444, 48.8583, 2.2945) # about 5837 km 
    print 
    print "calc min/max lat/lon" 
    degs = map(degrees, (1.3963, -0.6981)) 
    test(*degs, dist=1000) 
    print 
    degs = map(degrees, (1.3963, -0.6981, 1.4618, -1.6021)) 
    print degs, "distance", distance_between_points(*degs) # 872 km 

Answer 2

这是你如何计算使用haversine公式经/纬对之间的距离：

import math 

R = 6371 # km 
dLat = (lat2-lat1) # Make sure it's in radians, not degrees 
dLon = (lon2-lon1) # Idem 
a = math.sin(dLat/2) * math.sin(dLat/2) + 
    math.cos(lat1) * math.cos(lat2) * 
    math.sin(dLon/2) * math.sin(dLon/2) 
c = 2 * math.atan2(math.sqrt(a), math.sqrt(1-a)) 
d = R * c; 

现在是微不足道的测试“d”（也公里）对你的阈值。如果你想要的东西不是公里，请调整半径。

对不起，我不能给你一个插件解决方案，但我不明白你的代码框架（见评论）。

另请注意，这些天你可能想要使用余弦的球形定律而不是Haversine。数值稳定性的优点不再值得，而且它的理解，编码和使用都非常简单。

Answer 3

如果您在MongoDB中存储数据，它会为您编制好的地理定位搜索索引，并且优于上面的纯Python解决方案，因为它会为您处理优化。

http://www.mongodb.org/display/DOCS/Geospatial+Indexing

Answer 4

John Machin的回答对我很有帮助。只有一个小小的错误：经纬度被换成了boundigbox：

dlon = distance/RADIUS 
dlat = asin(sin(dlon)/cos(radians(lon))) 
return degrees(dlat), degrees(dlon) 

这就解决了这个问题。原因在于经度不会改变它们每度的距离 - 但纬度确实如此。它们的距离取决于经度。