Update data

This commit is contained in:
Andras Schmelczer 2026-05-14 08:17:10 +01:00
parent a4103b0896
commit 273d7a83ee
15 changed files with 716 additions and 316 deletions

View file

@ -7,6 +7,8 @@ from scipy.spatial import cKDTree
from .haversine import haversine_km
EARTH_RADIUS_KM = 6371.0088
KM_PER_DEGREE_LAT = 111.32
DEFAULT_GRID_SIZE_DEGREES = 0.02
def _build_poi_grid(
@ -34,16 +36,29 @@ def _build_poi_grid(
def _get_nearby_indices(
pc_lat: float, pc_lon: float, poi_grid: dict, grid_size: float = 0.05
pc_lat: float,
pc_lon: float,
poi_grid: dict,
radius_km: float,
grid_size: float = DEFAULT_GRID_SIZE_DEGREES,
) -> np.ndarray | None:
"""Get POI indices from grid cells near the given coordinate."""
grid_lat = int(np.floor(pc_lat / grid_size))
grid_lng = int(np.floor(pc_lon / grid_size))
"""Get POI indices from all grid cells intersecting the radius bounding box."""
if not np.isfinite(pc_lat) or not np.isfinite(pc_lon):
return None
lat_delta = radius_km / KM_PER_DEGREE_LAT
cos_lat = abs(np.cos(np.radians(pc_lat)))
lng_delta = 180.0 if cos_lat < 1e-12 else radius_km / (KM_PER_DEGREE_LAT * cos_lat)
min_grid_lat = int(np.floor((pc_lat - lat_delta) / grid_size))
max_grid_lat = int(np.floor((pc_lat + lat_delta) / grid_size))
min_grid_lng = int(np.floor((pc_lon - lng_delta) / grid_size))
max_grid_lng = int(np.floor((pc_lon + lng_delta) / grid_size))
nearby_indices = []
for dlat in [-1, 0, 1]:
for dlng in [-1, 0, 1]:
cell_key = (grid_lat + dlat, grid_lng + dlng)
for grid_lat in range(min_grid_lat, max_grid_lat + 1):
for grid_lng in range(min_grid_lng, max_grid_lng + 1):
cell_key = (grid_lat, grid_lng)
if cell_key in poi_grid:
nearby_indices.append(poi_grid[cell_key])
@ -83,7 +98,7 @@ def count_pois_per_postcode(
n_pois = len(pois)
print(f" {n_postcodes:,} postcodes, {n_pois:,} POIs")
grid_size = 0.05
grid_size = DEFAULT_GRID_SIZE_DEGREES
print(" Building POI spatial grid...")
poi_lats, poi_lngs, poi_cats, poi_grid = _build_poi_grid(pois, grid_size)
print(f" POI grid has {len(poi_grid):,} occupied cells")
@ -120,7 +135,9 @@ def count_pois_per_postcode(
# Process batch
for i in range(start_idx, end_idx):
nearby = _get_nearby_indices(pc_lats[i], pc_lons[i], poi_grid, grid_size)
nearby = _get_nearby_indices(
pc_lats[i], pc_lons[i], poi_grid, radius_km, grid_size
)
if nearby is None:
continue