Testing geo bounds http://janmatuschek.de/LatitudeLongitudeBoundingCoordinates and http://gis.stackexchange.com/questions/19221/find-tangent-point-on-circle-furthest-east-or-west ; Detailed diagram for the same here http://gis.stackexchange.com/a/195148/61719
from __future__ import division | |
import numpy as np | |
RADIUS_OF_EARTH_IN_KM = 6371.01 | |
def haversine(lat1, lon1, lat2, lon2): | |
""" | |
Utility to calcutlate distance between two pointtodo explain regarding height | |
coverting from geodisc co-ordinate to cartestian gives errors when distances are further apart | |
""" | |
dLat = np.radians(lat2 - lat1) | |
dLon = np.radians(lon2 - lon1) | |
lat1 = np.radians(lat1) | |
lat2 = np.radians(lat2) | |
a = np.sin(dLat/2)**2 + np.cos(lat1)*np.cos(lat2)*np.sin(dLon/2)**2 | |
c = 2*np.arcsin(np.sqrt(a)) | |
return RADIUS_OF_EARTH_IN_KM * c | |
from __future__ import division | |
import numpy as np | |
def get_bounding_box(distance, latittude, longitude): | |
""" | |
Given a Geographic co-ordintate in degress, aim is to get a bounding box of a given point, | |
basically a maximum latitude, minmun latitude and max longitude and min longitude based on explanation | |
and formulae by these two sites | |
http://janmatuschek.de/LatitudeLongitudeBoundingCoordinates | |
http://gis.stackexchange.com/questions/19221/find-tangent-point-on-circle-furthest-east-or-west | |
:param distance: bounding box distance in kilometer | |
:param latittude: of the point of the center of the box | |
:param longitude: of the point of the center of the box | |
:return: top left co-ordinates, and top-right co-ordinates of the bounding box | |
""" | |
# Convert to radians | |
latittude = np.radians(latittude) | |
longitude = np.radians(longitude) | |
if not MIN_LATITUDE <= latittude <= MAX_LATITUDE: | |
raise ValueError("Invalid latitude passed in") | |
if not MIN_LONGITUDE <= longitude <= MAX_LONGITUDE: | |
raise ValueError("Invalid longitude passed in") | |
distance_from_point_km = distance | |
angular_distance = distance_from_point_km / RADIUS_OF_EARTH_IN_KM | |
# Moving along a latitude north or south, keeping the longitude fixed, is travelling along | |
# a great circle; which means basically that doing so,the Radius from center of the earth , | |
# to the circle traced by such a path is constant; So to derive the top latitude and bottom latitude | |
# without any error you just need to add or suntract the angular distance | |
lat_min = latittude - angular_distance | |
lat_max = latittude + angular_distance | |
# Handle poles here | |
if lat_min < MIN_LATITUDE or lat_max > MAX_LATITUDE: | |
lat_min = max(lat_min, MIN_LATITUDE) | |
lat_max = min(lat_max, MAX_LATITUDE) | |
lon_max = MAX_LONGITUDE | |
lon_min = MIN_LONGITUDE | |
lat_max = np.degrees(lat_max) | |
lat_min = np.degrees(lat_min) | |
lon_max = np.degrees(lon_max) | |
lon_min = np.degrees(lon_min) | |
return (lat_min, lon_min), (lat_max, lon_max), lat_max | |
# the latitude which intersects T1 and T2 in [1] and [2] is based on speherical trignomerty | |
# this is just calculated for verifying this method | |
latitude_of_intersection = np.arcsin(np.sin(latittude) / np.cos(angular_distance)) | |
latitude_of_intersection = np.degrees(latitude_of_intersection) | |
# if we draw a circle with distance around the lat and longitude, then the point where the | |
# logitude touches the circle will be at a latitude a bit apart from the original latitiude as the | |
# circle is on the sphere , the earths sphere | |
delta_longitude = np.arcsin(np.sin(angular_distance) / np.cos(latittude)) | |
lon_min = longitude - delta_longitude | |
lon_max = longitude + delta_longitude | |
lon_min = np.degrees(lon_min) | |
lat_max = np.degrees(lat_max) | |
lon_max = np.degrees(lon_max) | |
lat_min = np.degrees(lat_min) | |
return (lat_min, lon_min), (lat_max, lon_max), latitude_of_intersection |
This comment has been minimized.
This comment has been minimized.
Poles are not handled currently |
This comment has been minimized.
This comment has been minimized.
poles also handled |
This comment has been minimized.
This comment has been minimized.
Looks like you forgot to convert lat_min to degrees. |
This comment has been minimized.
This comment has been minimized.
good catch; will change; However the (latitude_of_intersection, min_longitude) is what defines the top left corner of the bounding box
|
This comment has been minimized.
This comment has been minimized.
Nice! Appreciate you sharing this. Thought you might want to know latitude is misspelled in the get_bounding_box function in several places stemming from the typo in the parameters for the function. Cheers |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
This comment has been minimized.
Bounding Box
lat,long,angular_distance 13.08955 78.0335055556 0.000784584484057
delta_longitude 0.000805513923644 0.0461525481638
lat_min,lat_max 13.0445966204 13.1345033796
lon_min,lon_max 77.9873530074 78.0796581038
Distance from Point to T1 4.99859559377
Distnace from Point to T2 4.99859559377
Total time = 0.00556182861328