"""Show the same XYZ basemaps in Python with folium, leafmap and contextily. Writes two HTML maps and one PNG figure (built from QGIS renders plus a contextily render) to ./output. """ from pathlib import Path import contextily as cx import folium import leafmap.foliumap as leafmap import matplotlib.pyplot as plt import matplotlib.image as mpimg import xyzservices.providers as xyz from pyproj import Transformer OUT = Path(__file__).parent / "output" OUT.mkdir(exist_ok=True) OSM_URL = "https://tile.openstreetmap.org/{z}/{x}/{y}.png" OSM_ATTR = '© OpenStreetMap contributors' EOX_URL = "https://tiles.maps.eox.at/wmts/1.0.0/s2cloudless_3857/default/g/{z}/{y}/{x}.jpg" EOX_ATTR = ('Sentinel-2 cloudless - https://s2maps.eu ' "by EOX IT Services GmbH (contains modified Copernicus Sentinel data 2016)") CENTRE = [53.75, 0.01] # Withernsea, Holderness coast # folium: OpenStreetMap is the default; extra layers need url + attribution m = folium.Map(location=CENTRE, zoom_start=11, tiles=None) folium.TileLayer(OSM_URL, attr=OSM_ATTR, name="OpenStreetMap", max_zoom=19).add_to(m) folium.TileLayer(EOX_URL, attr=EOX_ATTR, name="Sentinel-2 cloudless 2016", max_zoom=14).add_to(m) folium.TileLayer(xyz.Esri.WorldImagery.build_url(), attr=xyz.Esri.WorldImagery.html_attribution, name="Esri World Imagery").add_to(m) folium.LayerControl().add_to(m) m.save(OUT / "folium_basemaps.html") print("folium", folium.__version__, "-> folium_basemaps.html") # leafmap: named basemaps or any XYZ URL lm = leafmap.Map(center=CENTRE, zoom=11) lm.add_tile_layer(url=EOX_URL, name="Sentinel-2 cloudless 2016", attribution=EOX_ATTR) lm.to_html(str(OUT / "leafmap_basemaps.html")) print("leafmap -> leafmap_basemaps.html") # contextily: static basemap under a matplotlib axis (EPSG:3857 metres) west, south, east, north = -0.10, 53.70, 0.12, 53.80 fig, axes = plt.subplots(1, 3, figsize=(15, 3.9), constrained_layout=True) for ax, fname, title in [ (axes[0], "qgis_osm.png", "QGIS 3.40 XYZ: OpenStreetMap"), (axes[1], "qgis_eox_2016.png", "QGIS 3.40 XYZ: Sentinel-2 cloudless 2016"), ]: ax.imshow(mpimg.imread(OUT / fname)) ax.set_title(title, fontsize=11) ax.set_axis_off() axes[0].text(0.99, 0.01, "© OpenStreetMap contributors", transform=axes[0].transAxes, ha="right", va="bottom", fontsize=8, bbox=dict(fc="white", ec="none", alpha=0.8)) axes[1].text(0.99, 0.01, "Sentinel-2 cloudless 2016 by EOX (CC BY 4.0)\nContains modified Copernicus Sentinel data 2016", transform=axes[1].transAxes, ha="right", va="bottom", fontsize=7, bbox=dict(fc="white", ec="none", alpha=0.8)) to_merc = Transformer.from_crs("EPSG:4326", "EPSG:3857", always_xy=True) x0, y0 = to_merc.transform(west, south) x1, y1 = to_merc.transform(east, north) axes[2].set_xlim(x0, x1) axes[2].set_ylim(y0, y1) cx.add_basemap(axes[2], crs="EPSG:3857", source=xyz.OpenStreetMap.Mapnik, zoom=11, attribution=False, headers={"User-Agent": "TechnicalGIS-tutorial/1.0 (+https://technicalgis.com)"}) print("contextily", cx.__version__, "extent (EPSG:3857 m)", [round(v) for v in (x0, y0, x1, y1)]) axes[2].set_title("contextily: OpenStreetMap.Mapnik, zoom 11", fontsize=11) axes[2].set_axis_off() axes[2].text(0.99, 0.01, "© OpenStreetMap contributors", transform=axes[2].transAxes, ha="right", va="bottom", fontsize=8, bbox=dict(fc="white", ec="none", alpha=0.8)) fig.savefig(OUT / "basemaps_qgis_and_python.png", dpi=110) print("saved basemaps_qgis_and_python.png")