From 94898aca13cab8a160b37d3f6620fa0c1b18f683 Mon Sep 17 00:00:00 2001 From: Matthieu Bessat Date: Wed, 29 May 2024 15:04:49 +0200 Subject: [PATCH] feat: coordinates bound computation to geo tiles --- .../plot_osm_map_with_coordinates.py | 72 +++++++++++++++++++ ...sm_map.py => plot_osm_map_with_tile_nb.py} | 2 + 2 files changed, 74 insertions(+) create mode 100644 geo_tiles_assembler/plot_osm_map_with_coordinates.py rename geo_tiles_assembler/{plot_osm_map.py => plot_osm_map_with_tile_nb.py} (97%) diff --git a/geo_tiles_assembler/plot_osm_map_with_coordinates.py b/geo_tiles_assembler/plot_osm_map_with_coordinates.py new file mode 100644 index 0000000..62e5cd3 --- /dev/null +++ b/geo_tiles_assembler/plot_osm_map_with_coordinates.py @@ -0,0 +1,72 @@ +import matplotlib.pyplot as plt +import matplotlib.image as img +import requests +import numpy as np +import os +import math + +def deg2num(lat_deg, lon_deg, zoom): + lat_rad = math.radians(lat_deg) + n = 1 << zoom + xtile = int((lon_deg + 180.0) / 360.0 * n) + ytile = int((1.0 - math.asinh(math.tan(lat_rad)) / math.pi) / 2.0 * n) + return xtile, ytile + +BOUNDS = ((49.17, 1.31), (49.18, 1.35)) + +DEFAULT_HEADERS = {"User-Agent": "Mozilla/5.0 (Windows NT 10.0; rv:127.0) Gecko/20100101 Firefox/127.0"} + +iz = 0 +for z in range(16, 17): + # https://tile.openstreetmap.org/11/1029/700.png + # 1029, 1030 going to the est, increasing the longitude angle (x) + # 700, 701, going to the south, deacreasing the flatitude latitude angle (y) + + TILES_BOUNDS = ( + deg2num(BOUNDS[0][0], BOUNDS[0][1], z), + deg2num(BOUNDS[1][0], BOUNDS[1][1], z) + ) + print(TILES_BOUNDS) + + fig, ax = plt.subplots() + + maxx = 0 + maxy = 0 + xrange = range(TILES_BOUNDS[0][0], TILES_BOUNDS[1][0]) + yrange = range(TILES_BOUNDS[1][1], TILES_BOUNDS[0][1]) + print(xrange, yrange) + + iy = 0 + for y in reversed(yrange): + ix = 0 + for x in xrange: + # https://tile.openstreetmap.org/17/65956/44733.png + url = f"https://tile.openstreetmap.org/{z}/{x}/{y}.png" + dir_path = f'tiles/{z}/{x}' + tile_path = dir_path + f'/{y}.png' + img_data = 0 + print(tile_path) + if not os.path.exists(tile_path): + res = requests.get(url, headers=DEFAULT_HEADERS) + print(res) + if res.status_code != 200: + print(res.content) + exit() + img_data = res.content + os.makedirs(dir_path, exist_ok=True) + with open(tile_path, 'wb') as handler: + handler.write(img_data) + im = img.imread(tile_path) + print("read") + print(ix, iy) + ax.imshow(im, aspect='auto', extent=(ix, ix+1, iy, iy+1), zorder=-1) + print("imshow") + ix += 1 + iy += 1 + ax.set_xlim(0, len(xrange)) + ax.set_ylim(0, len(yrange)) + ax.set_aspect('equal', adjustable='box') + plt.show() + break + iz += 1 + diff --git a/geo_tiles_assembler/plot_osm_map.py b/geo_tiles_assembler/plot_osm_map_with_tile_nb.py similarity index 97% rename from geo_tiles_assembler/plot_osm_map.py rename to geo_tiles_assembler/plot_osm_map_with_tile_nb.py index f46b95f..66449fe 100644 --- a/geo_tiles_assembler/plot_osm_map.py +++ b/geo_tiles_assembler/plot_osm_map_with_tile_nb.py @@ -6,6 +6,8 @@ import os DEFAULT_HEADERS = {"User-Agent": "Mozilla/5.0 (Windows NT 10.0; rv:127.0) Gecko/20100101 Firefox/127.0"} +BOUNDS = ((49, 1.1), (49.3, 1.4)) + fig, ax = plt.subplots() maxx = 0