Get geographical coordinates in a given square

Question

My task is as follows: knowing the center (starting point), for example - [{'lat': -7.7940023, 'lng': 110.3656535}] and knowing the radius 5km I need to get all the points included in this square in 1 km increments. How do I achieve this?

P.S Using the Haversine formula I can check if a point is in a given square according to the radius

Image

This could be useful https://stackoverflow.com/questions/1253499/simple-calculations-for-working-with-lat-lon-and-km-distance — Riccardo Bucco, Aug 01 '22 at 14:52
I think it's actually a spherical circle, not a square, unless you mean 1 km in either direction, not 1 km total — Barry Carter, Aug 01 '22 at 14:54
you need to clarify if you need the points inside a circle of a given radius or inside a square with a given side. in addition, what are your original points you need to check? do you have a grid? do you have and increment in radius and angle? — mauro, Aug 01 '22 at 15:08
@mauro I added an image. I need to get inside this area all the points in 1 km increments — Lite, Aug 01 '22 at 15:15
OK, the edit shows it's a square, which is a lot easier. Just add 1km of latitude (which is fixed if you assume spherical Earth) to the latitude and cosine of the latitude times that to the longitude in each direction. — Barry Carter, Aug 01 '22 at 15:16

mauro · Accepted Answer · 2022-08-01T15:54:04.930

1

if you consider a spherical Earth with radius R, the associated angle of a segment with length L (5km in your case) is:

import numpy as np
R = 6378.0 # km
L = 5.0 # km
angle = np.degrees(L/R)

so now, you can easily check if a point is inside your square:

center = {'lat': -7.7940023, 'lng': 110.3656535}
point = {'lat': 'point latitude', 'lng': 'point longitude'} # insert here your values

if (point['lat']-center['lat'] < angle) and (point['lng']-center['lng'] < angle):
    print('Point is inside')
else:
    print('Point is outside')

EDIT: check the one below.

import numpy as np

R = 6378.0  # km
L = 5  # side of the square
center = {'lat': -7.7940023, 'lng': 110.3656535}

square_side = np.linspace(-L/2, L/2, L+1)
angle = np.degrees(square_side/R)
latitude, longitude = np.meshgrid(angle+center['lat'], angle+center['lng'])
points = []
for lat, lng in zip(latitude.flatten(), longitude.flatten()):
    points.append({'lat': lat, 'lng': lng})

edited Aug 01 '22 at 15:54

answered Aug 01 '22 at 15:29

mauro

504
3
14

that's right, thank you. but my question was how to get me a list of all the points that enter this area in 1km increments – Lite Aug 01 '22 at 15:32
do you need a list of dictionaries with ```lat``` and ```lng```? – mauro Aug 01 '22 at 15:33
yes, mate. I need a list of dictionaries of geographical coordinates (lat, lng) from this area in 1 km increments – Lite Aug 01 '22 at 15:35

score 0 · Answer 2 · edited Aug 11 '22 at 19:47

This example should illustrate the point:

import pandas as pd 
import numpy as np

#for_map = pd.read_csv('campaign_contributions_for_map.tsv', sep='\t')
df_airports = pd.read_csv('C:\\airports.csv')
print(df_airports.head(3))

df_cities = pd.read_csv('C:\\worldcities.csv')
print(df_cities.head(3))

# join the two dataframes - must be the same length
df = pd.concat([df_cities, df_airports], axis=1)

# cast latitudes and longitudes to numeric
cols = ["lat", "lng", "latitude_deg", "longitude_deg"]
df[cols] = df[cols].apply(pd.to_numeric, errors='coerce', axis=1)

# create a mask where our conditions are met (difference between lat fuze and lat air < 5 and difference between long fuze and long air < 0.1)
mask = ((abs(df["lat"] - df["latitude_deg"]) < .5) & (abs(df["lng"] - df["longitude_deg"]) < .5))

# fill the type column
df.loc[mask, 'Type'] = "Airport"
df.shape

pd.set_option('display.max_columns', None)
pd.set_option('display.max_rows', None)
df.head()

# The haversine formula determines the great-circle distance between two points on a sphere given 
# their longitudes and latitudes. Important in navigation, it is a special case of a more general 
# formula in spherical trigonometry, the law of haversines, that relates the sides and angles of 
# spherical triangles. 

lat1 = df['lat']
lon1 = df['lng']
lat2 = df['latitude_deg']
lon2 = df['longitude_deg']

from math import radians, cos, sin, asin, sqrt
def haversine(lat1, lon1, lat2, lon2):
    """
    Calculate the great circle distance between two points 
    on the earth (specified in decimal degrees)
    """
    # convert decimal degrees to radians 
    lon1, lat1, lon2, lat2 = map(radians, [lon1, lat1, lon2, lat2])
    # haversine formula 
    dlon = lon2 - lon1 
    dlat = lat2 - lat1 
    a = sin(dlat/2)**2 + cos(lat1) * cos(lat2) * sin(dlon/2)**2
    c = 2 * asin(sqrt(a)) 
    # Radius of earth in kilometers is 6371
    km = 6371* c
    return km

# Creating a new column to generate the output by passing lat long information to Haversine Equation
df['distance'] = [haversine(df.lat[i],df.lng[i],df.latitude_deg[i],df.longitude_deg[i]) for i in range(len(df))]
df['distance'] = df['distance'].round(decimals=3)

# Printing the data table 
df.sort_values(by=['distance'], inplace=True)
df.head()

# Let's sort by our newly created field, which identifies airport lat/lonn coordinates within .5 places of 
# a city's lat/long coordinates

# Create a mask where our conditions are met (difference between lat and latitude_deg < 0.1 and 
# difference between lng and longitude_deg < 0.1)
mask = ((abs(df["lat"] - df["latitude_deg"]) < 0.1) & (abs(df["lng"] - df["longitude_deg"]) < 0.1))

# Fill the type column
df.loc[mask, 'Type'] = "Airport"

df.sort_values(by=['Type'], inplace=True)
df.head()

More details here.

https://github.com/ASH-WICUS/Notebooks/blob/master/Haversine%20Distance%20-%20Airport%20or%20Not.ipynb

Get geographical coordinates in a given square

2 Answers2