Created
March 29, 2025 05:45
-
-
Save leighklotz/8400a496cea83d9590e691161ee98728 to your computer and use it in GitHub Desktop.
Extract tiles from a Zip of OSM tiles, for Meshastic, APRS, etc.
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| #!/usr/bin/env python3 | |
| # $ ./extract-bbox.py --zipfile OxedOSM_USA_0-13.zip --depth 13 --center-lat 37.123456 --center-lon -122.876543 --radius 10 | |
| import os | |
| import zipfile | |
| import argparse | |
| from math import floor, pi, tan, cos, log | |
| def main(): | |
| parser = argparse.ArgumentParser(description="Extract tiles from a zip file based on a bounding box and zoom depth.") | |
| parser.add_argument("--zipfile", required=True, help="Path to the zip file containing the tiles.") | |
| parser.add_argument("--min-lon", type=float, help="Minimum longitude of the bounding box.") | |
| parser.add_argument("--min-lat", type=float, help="Minimum latitude of the bounding box.") | |
| parser.add_argument("--max-lon", type=float, help="Maximum longitude of the bounding box.") | |
| parser.add_argument("--max-lat", type=float, help="Maximum latitude of the bounding box.") | |
| parser.add_argument("--depth", type=int, required=True, help="Maximum zoom depth to extract tiles for.") | |
| parser.add_argument("--center-lat", type=float, help="Center latitude.") | |
| parser.add_argument("--center-lon", type=float, help="Center longitude.") | |
| parser.add_argument("--radius", type=float, help="Radius in units around the center point.") | |
| parser.add_argument("--units", type=str, default="miles", choices=["miles", "kilometers"], help="Units for the radius (miles or kilometers).") | |
| args = parser.parse_args() | |
| # Calculate bounding box from center and radius if provided | |
| if args.center_lat and args.center_lon and args.radius: | |
| # Approximate conversion factor from degrees to miles/kilometers at the given center latitude | |
| # This is a simplification; actual conversion varies with latitude. | |
| if args.units == "miles": | |
| deg_to_miles = 69.0 * cos(args.center_lat * pi / 180.0) # Approximate miles per degree at center latitude | |
| else: # kilometers | |
| deg_to_miles = 111.0 * cos(args.center_lat * pi / 180.0) # Approximate km per degree at center latitude | |
| min_lat = args.center_lat - (args.radius / deg_to_miles) | |
| max_lat = args.center_lat + (args.radius / deg_to_miles) | |
| min_lon = args.center_lon - (args.radius / (deg_to_miles * cos(args.center_lat * pi / 180.0))) | |
| max_lon = args.center_lon + (args.radius / (deg_to_miles * cos(args.center_lat * pi / 180.0))) | |
| else: | |
| min_lat = args.min_lat | |
| max_lat = args.max_lat | |
| min_lon = args.min_lon | |
| max_lon = args.max_lon | |
| # Calculate tile boundaries | |
| depth = args.depth | |
| min_x = long_to_tile_x(min_lon, args.depth) | |
| max_x = long_to_tile_x(max_lon, args.depth) | |
| min_y = lat_to_tile_y(min_lat, args.depth) | |
| max_y = lat_to_tile_y(max_lat, args.depth) | |
| min_x, max_x = min(min_x, max_x), max(min_x, max_x) | |
| min_y, max_y = min(min_y, max_y), max(min_y, max_y) | |
| print(f"/{depth}/[{min_x}-{max_x}]/[{min_y}-{max_y}].png") | |
| try: | |
| with zipfile.ZipFile(args.zipfile, 'r') as zf: | |
| for filename in zf.namelist(): | |
| if not filename.endswith(".png"): | |
| continue # Skip non-png files | |
| # Extract zoom level, x, and y from the filename | |
| parts = filename.split('/') # Assuming the path within the zip matches the directory structure | |
| try: | |
| zoom_level = int(parts[1]) | |
| x = int(parts[2]) | |
| y = int(parts[3].split('.')[0]) #handle the .png extension | |
| except (IndexError, ValueError): | |
| print(f"Skipping file {filename}: Invalid filename format.") | |
| continue | |
| if zoom_level == depth: | |
| if is_tile_in_area(x, y, zoom_level, min_x, max_x, min_y, max_y): | |
| # Construct the output directory path | |
| output_dir = os.path.join(os.getcwd(), "OxedOSM-USA", str(zoom_level), str(x)) | |
| os.makedirs(output_dir, exist_ok=True) # Create directories if they don't exist | |
| output_path = os.path.join(output_dir, f"{y}.png") | |
| # Extract the tile to the output path | |
| with zf.open(filename) as source, open(output_path, "wb") as dest: | |
| dest.write(source.read()) | |
| print(f"Extracted: {filename} to {output_path}") | |
| else: | |
| # print(f"Skipping: {filename} (outside area)") | |
| pass | |
| except FileNotFoundError: | |
| print(f"Error: Zip file not found at {args.zipfile}") | |
| except zipfile.BadZipFile: | |
| print(f"Error: Invalid zip file: {args.zipfile}") | |
| except Exception as e: | |
| print(f"An unexpected error occurred: {e}", exc_info=True) | |
| def long_to_tile_x(lon, zoom): | |
| return int(floor(((lon + 180.0) / 360.0) * (2**zoom))) | |
| def lat_to_tile_y(lat, zoom): | |
| """ | |
| Converts latitude to Y tile coordinate at the given zoom level. | |
| """ | |
| n = 2.0 ** zoom | |
| lat_rad = lat * pi / 180.0 # Convert latitude to radians | |
| return int(floor((1 - log(tan(lat_rad) + (1 / cos(lat_rad))) / pi) / 2 * n)) | |
| def is_tile_in_area(x, y, zoom, min_x, max_x, min_y, max_y): | |
| return min_x <= x <= max_x and min_y <= y <= max_y | |
| if __name__ == "__main__": | |
| main() |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment