100 lines
3.6 KiB
Python
100 lines
3.6 KiB
Python
"""
|
|
HELPER SCRIPT: generate_keys.py
|
|
-------------------------------
|
|
This script generates a list of 'Quadkeys' (map tiles) that intersect a specific
|
|
geographic area (WKT Geometry). This list is used to seed the scraper in 'power2.py'
|
|
so it knows exactly where to look for outages without scanning the whole world.
|
|
|
|
PREREQUISITES:
|
|
pip install shapely mercantile
|
|
|
|
USAGE INSTRUCTIONS:
|
|
1. Obtain the WKT (Well-Known Text) geometry for your target area from your database.
|
|
|
|
Example SQL to get WKT for a specific county:
|
|
SELECT ST_AsText(geom) FROM county WHERE countyname = 'Kanawha' AND state = 'WV';
|
|
|
|
Example SQL to get WKT for a whole NWS Warning Area (Union of counties):
|
|
SELECT ST_AsText(ST_Union(geom)) FROM county WHERE cwa = 'RLX';
|
|
|
|
2. Run this script:
|
|
python generate_keys.py
|
|
|
|
3. Paste the huge text string returned by your SQL query into the prompt and hit Enter.
|
|
|
|
4. Copy the resulting 'KEY_LIST = [...]' output and paste it into the configuration
|
|
section of 'power2.py' (e.g., replace AEP_WV_KEYS).
|
|
"""
|
|
|
|
import mercantile
|
|
from shapely import wkt
|
|
from shapely.geometry import box
|
|
import sys
|
|
|
|
# --- CONFIGURATION ---
|
|
# The zoom level to generate keys for.
|
|
# Level 7 is standard for the "entry" lists in your scraper (e.g. '0320001').
|
|
# If the utility uses a different zoom level for its top-level clusters, adjust this.
|
|
TARGET_ZOOM = 7
|
|
|
|
# Increase the CSV field size limit to handle massive WKT strings
|
|
import csv
|
|
csv.field_size_limit(sys.maxsize)
|
|
|
|
def generate_keys_from_wkt(wkt_string):
|
|
"""
|
|
Takes a WKT geometry string, finds the bounding box, generates all tiles
|
|
within that box, and filters them to keep only those that actually
|
|
intersect the geometry.
|
|
"""
|
|
try:
|
|
print("Parsing WKT...")
|
|
# 1. Parse the WKT
|
|
service_area = wkt.loads(wkt_string)
|
|
|
|
# 2. Get the Bounding Box (minx, miny, maxx, maxy)
|
|
bounds = service_area.bounds
|
|
|
|
# 3. Generate all tiles covering the bounding box at the target zoom
|
|
# mercantile.tiles(west, south, east, north, zooms)
|
|
candidate_tiles = list(mercantile.tiles(bounds[0], bounds[1], bounds[2], bounds[3], zooms=[TARGET_ZOOM]))
|
|
|
|
valid_keys = []
|
|
|
|
print(f"Scanning {len(candidate_tiles)} candidate tiles...")
|
|
|
|
for tile in candidate_tiles:
|
|
# 4. Create a geometry for the tile itself
|
|
# mercantile.bounds returns (west, south, east, north)
|
|
t_bounds = mercantile.bounds(tile)
|
|
tile_geom = box(t_bounds.west, t_bounds.south, t_bounds.east, t_bounds.north)
|
|
|
|
# 5. Intersection Check: Does this tile actually touch the service area?
|
|
# We use intersects() because contains() might miss tiles that only partially overlap.
|
|
if service_area.intersects(tile_geom):
|
|
valid_keys.append(mercantile.quadkey(tile))
|
|
|
|
# 6. Output formatted for Python list
|
|
valid_keys.sort()
|
|
print(f"\nFound {len(valid_keys)} intersecting tiles.")
|
|
print("-" * 30)
|
|
print(f"KEY_LIST = {valid_keys}")
|
|
print("-" * 30)
|
|
|
|
except Exception as e:
|
|
print(f"Error processing geometry: {e}")
|
|
|
|
if __name__ == "__main__":
|
|
print("Paste your PostGIS WKT (Well-Known Text) string below and press Enter:")
|
|
|
|
# Use standard input reading to handle potentially large inputs better than input()
|
|
try:
|
|
user_wkt = input().strip()
|
|
if user_wkt:
|
|
generate_keys_from_wkt(user_wkt)
|
|
else:
|
|
print("No input provided.")
|
|
except EOFError:
|
|
print("Error reading input.")
|
|
except KeyboardInterrupt:
|
|
print("\nCancelled.") |