faster lookup
This commit is contained in:
parent
a4073284ee
commit
be5938ffb0
1 changed files with 21 additions and 41 deletions
60
lookup.py
60
lookup.py
|
|
@ -2,9 +2,10 @@ import json
|
||||||
import click
|
import click
|
||||||
import pygeohash as pgh
|
import pygeohash as pgh
|
||||||
import pyarrow.parquet as pq
|
import pyarrow.parquet as pq
|
||||||
|
import geopandas as gpd
|
||||||
from shapely.geometry import Point, LineString, MultiLineString
|
from shapely.geometry import Point, LineString, MultiLineString
|
||||||
from shapely import wkb
|
|
||||||
from pyproj import Geod
|
from pyproj import Geod
|
||||||
|
import numpy as np
|
||||||
|
|
||||||
@click.command()
|
@click.command()
|
||||||
@click.argument("latitude", type=float)
|
@click.argument("latitude", type=float)
|
||||||
|
|
@ -15,15 +16,13 @@ def geohash_lookup(latitude, longitude, distance, parquet_dir):
|
||||||
"""
|
"""
|
||||||
Finds polygons containing the given point and lines within a threshold (in meters).
|
Finds polygons containing the given point and lines within a threshold (in meters).
|
||||||
"""
|
"""
|
||||||
|
|
||||||
# Compute geohash (precision can be adjusted)
|
# Compute geohash (precision can be adjusted)
|
||||||
geohash = pgh.encode(latitude, longitude, precision=4) # Adjust precision as needed
|
geohash = pgh.encode(latitude, longitude, precision=4) # Adjust precision as needed
|
||||||
parquet_file = f"{parquet_dir}/{geohash}.parquet"
|
parquet_file = f"{parquet_dir}/{geohash}.parquet"
|
||||||
|
|
||||||
try:
|
try:
|
||||||
# Load the Parquet file
|
# Load the Parquet file directly into a GeoDataFrame
|
||||||
table = pq.read_table(parquet_file)
|
df = gpd.read_parquet(parquet_file)
|
||||||
df = table.to_pandas()
|
|
||||||
except FileNotFoundError:
|
except FileNotFoundError:
|
||||||
click.echo(json.dumps({"error": f"No Parquet file found for geohash {geohash}."}))
|
click.echo(json.dumps({"error": f"No Parquet file found for geohash {geohash}."}))
|
||||||
return
|
return
|
||||||
|
|
@ -33,46 +32,27 @@ def geohash_lookup(latitude, longitude, distance, parquet_dir):
|
||||||
|
|
||||||
found_geometries = []
|
found_geometries = []
|
||||||
|
|
||||||
for _, row in df.iterrows():
|
# Vectorized distance calculation for LineStrings and MultiLineStrings
|
||||||
# Decode geometry from WKB
|
def calculate_min_distance(geom):
|
||||||
try:
|
if isinstance(geom, LineString):
|
||||||
geom = wkb.loads(bytes(row["geometry"])) # Ensure WKB is properly read
|
return min(geod.inv(point.x, point.y, p[0], p[1])[2] for p in geom.coords)
|
||||||
except Exception as e:
|
elif isinstance(geom, MultiLineString):
|
||||||
click.echo(json.dumps({"error": f"Failed to decode geometry: {str(e)}"}))
|
return min(min(geod.inv(point.x, point.y, p[0], p[1])[2] for p in line.coords) for line in geom.geoms)
|
||||||
continue
|
return np.inf
|
||||||
|
|
||||||
|
for row in df.itertuples(index=False):
|
||||||
|
geom = row.geometry # Directly use the Shapely geometry object
|
||||||
|
|
||||||
# ✅ Ensure `tags` is always a dictionary
|
|
||||||
metadata = {}
|
metadata = {}
|
||||||
if isinstance(row["tags"], dict):
|
if isinstance(row.tags, dict):
|
||||||
metadata = {k: v for k, v in row["tags"].items() if v is not None} # Remove empty values
|
metadata = {k: v for k, v in row.tags.items() if v is not None} # Remove empty values
|
||||||
|
|
||||||
# ✅ Check if the point is inside a polygon
|
|
||||||
if geom.contains(point):
|
if geom.contains(point):
|
||||||
found_geometries.append({"type": "polygon", **metadata})
|
found_geometries.append({"type": "polygon", **metadata})
|
||||||
|
elif isinstance(geom, (LineString, MultiLineString)):
|
||||||
# ✅ Handle LineString and MultiLineString
|
min_distance = calculate_min_distance(geom)
|
||||||
elif isinstance(geom, LineString):
|
if min_distance <= distance:
|
||||||
try:
|
found_geometries.append({"type": "line" if isinstance(geom, LineString) else "multi_line", **metadata})
|
||||||
min_distance = min(geod.inv(point.x, point.y, p[0], p[1])[2] for p in geom.coords)
|
|
||||||
if min_distance <= distance:
|
|
||||||
found_geometries.append({"type": "line", **metadata})
|
|
||||||
except Exception as e:
|
|
||||||
click.echo(json.dumps({"error": f"Failed distance calculation (LineString): {str(e)}"}))
|
|
||||||
continue
|
|
||||||
|
|
||||||
elif isinstance(geom, MultiLineString):
|
|
||||||
try:
|
|
||||||
# ✅ Iterate through all LineStrings inside the MultiLineString
|
|
||||||
min_distance = min(
|
|
||||||
min(geod.inv(point.x, point.y, p[0], p[1])[2] for p in line.coords)
|
|
||||||
for line in geom.geoms
|
|
||||||
)
|
|
||||||
|
|
||||||
if min_distance <= distance:
|
|
||||||
found_geometries.append({"type": "multi_line", **metadata})
|
|
||||||
except Exception as e:
|
|
||||||
click.echo(json.dumps({"error": f"Failed distance calculation (MultiLineString): {str(e)}"}))
|
|
||||||
continue
|
|
||||||
|
|
||||||
# Output JSON result
|
# Output JSON result
|
||||||
click.echo(json.dumps(found_geometries, indent=2, ensure_ascii=False)) # Ensure proper Unicode display
|
click.echo(json.dumps(found_geometries, indent=2, ensure_ascii=False)) # Ensure proper Unicode display
|
||||||
|
|
|
||||||
Loading…
Add table
Add a link
Reference in a new issue