Skip to content

Commit

Permalink
New function for calculating great circle distance
Browse files Browse the repository at this point in the history
Calculated using the Haversine formula.

I also removed some of the existing functions in pandana.utils
that depended on GeoPandas and other geo libraries.
  • Loading branch information
jiffyclub committed Dec 19, 2014
1 parent 1e1875c commit 55675f9
Show file tree
Hide file tree
Showing 2 changed files with 54 additions and 51 deletions.
16 changes: 16 additions & 0 deletions pandana/tests/test_great_circle_dist.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,16 @@
import numpy.testing as npt

from pandana.utils import great_circle_dist as gcd


def test_gcd():
# tested against geopy
# https://geopy.readthedocs.org/en/latest/#module-geopy.distance
lat1 = 41.49008
lon1 = -71.312796
lat2 = 41.499498
lon2 = -81.695391

expected = 864456.76162966

npt.assert_allclose(gcd(lat1, lon1, lat2, lon2), expected)
89 changes: 38 additions & 51 deletions pandana/utils.py
Original file line number Diff line number Diff line change
@@ -1,51 +1,38 @@
from shapely.geometry import Point
from fiona import crs
import geopandas as gpd
import pandas as pd
import numpy as np


def bbox_convert(bbox, from_epsg, to_epsg):
bbox = gpd.GeoSeries([Point(bbox[0], bbox[1]),
Point(bbox[2], bbox[3])],
crs=crs.from_epsg(from_epsg))
bbox = bbox.to_crs(epsg=to_epsg)
bbox = [bbox[0].x, bbox[0].y, bbox[1].x, bbox[1].y]
return bbox


def get_nodes_from_osm(bbox, query, to_epsg=3740):
gdf = gpd.io.osm.query_osm('node',
bbox=bbox,
tags=query)
gdf = gdf[gdf.type == 'Point'].to_crs(epsg=to_epsg)
print "Found %d nodes" % len(gdf)
x, y = zip(*[(p.x, p.y) for (i, p)
in gdf.geometry.iteritems()])
x = pd.Series(x)
y = pd.Series(y)
return x, y


def anything_score(net, config, max_distance, decay, bbox):
score = pd.Series(np.zeros(len(net.node_ids)), index=net.node_ids)
for query, weights in config.iteritems():
print "Computing for query: %s" % query
print "Fetching nodes from OSM"
x, y = get_nodes_from_osm(bbox, query)
print "Done"
net.set_pois(query, x, y)
print "Computing nearest"
df = net.nearest_pois(max_distance, query, num_pois=len(weights))
print "Done"

for idx, weight in enumerate(weights):
# want the 1st not the 0th
idx += 1
print "Adding contribution %f for number %d nearest" % \
(weight, idx)
score += decay(df[idx])*weight
# print score.describe()

assert score.min() > 0
return score/score.max()*100
from __future__ import division

import math


def great_circle_dist(lat1, lon1, lat2, lon2):
"""
Get the distance (in meters) between two lat/lon points
via the Haversine formula.
Parameters
----------
lat1, lon1, lat2, lon2 : float
Latitude and longitude in degrees.
Returns
-------
dist : float
Distance in meters.
"""
radius = 6372795 # meters

lat1 = math.radians(lat1)
lon1 = math.radians(lon1)
lat2 = math.radians(lat2)
lon2 = math.radians(lon2)

dlat = lat2 - lat1
dlon = lon2 - lon1

# formula from:
# http://en.wikipedia.org/wiki/Haversine_formula#The_haversine_formula
a = math.pow(math.sin(dlat / 2), 2)
b = math.cos(lat1) * math.cos(lat2) * math.pow(math.sin(dlon / 2), 2)
d = 2 * radius * math.asin(math.sqrt(a + b))

return d

0 comments on commit 55675f9

Please sign in to comment.