diff --git a/download_files_from_LGL_BW.py b/download_files_from_LGL_BW.py index 2c57c00664b0323fb6c577ca77147c7f6ee301e1..0d25e47372cd8d0835c6dffaaefee47403839451 100644 --- a/download_files_from_LGL_BW.py +++ b/download_files_from_LGL_BW.py @@ -41,9 +41,11 @@ EXTRACT_REGIONS = True CITYGML_SERVER = "https://opengeodata.lgl-bw.de/data/lod2" RASTER = 2 # [km] BUNDESLAND = 'bw' -WGS84 = CRS.from_epsg(4326) +WGS84 = CRS.from_epsg(4326) LOCAL_CRS = CRS.from_epsg(32632) # UTM32N, used in BW. https://epsg.io/32632 +TO_LOCAL_CRS = Transformer.from_crs(WGS84, LOCAL_CRS, always_xy=True) + UTM = 32 SCRIPT_DIR = Path(__file__).parent @@ -62,10 +64,9 @@ if EXTRACT_REGIONS: GML_GLOB = "LoD2_*/LoD2_*.gml" -def coordinates_to_grid(longitude: float, latitude: float, local_utm: CRS) -> tuple[int, int]: +def coordinates_to_grid(longitude: float, latitude: float) -> tuple[int, int]: """Returns (x, y) of the tile on CITYGML_SERVER containing a given point.""" - transformer = Transformer.from_crs(WGS84, local_utm, always_xy=True) - x, y = transformer.transform(longitude, latitude) + x, y = TO_LOCAL_CRS.transform(longitude, latitude) x = floor(x / 1000) - 1 # Odd x y = floor(y / 1000) # Even y x -= x % RASTER @@ -88,8 +89,8 @@ def wkt_polygon_to_grid_coords(location_name: str, wkt: str) -> tuple[int, int, print("%s (%.3f°N %.3f°E -> %.3f°N %.3f°E)" % (location_name, max_lat, min_lon, min_lat, max_lon)) - x1, y1 = coordinates_to_grid(min_lon, min_lat, LOCAL_CRS) - x2, y2 = coordinates_to_grid(max_lon, max_lat, LOCAL_CRS) + x1, y1 = coordinates_to_grid(min_lon, min_lat) + x2, y2 = coordinates_to_grid(max_lon, max_lat) return (x1, x2, y1, y2) @@ -181,8 +182,7 @@ def main(polygons: dict[str, str]) -> None: def show_coordinates(match): longitude, latitude = match.groups() - transformer = Transformer.from_crs(WGS84, LOCAL_CRS, always_xy=True) - x, y = transformer.transform(longitude, latitude) + x, y = TO_LOCAL_CRS.transform(longitude, latitude) return f"{x} {y}" def convert_wkt_to_local(wkt):