OpenStreetMap vs Global Building Atlas Explained -
Python Tutorial
import osmnx as ox
import math
import geopandas as gpd
import urllib.request
import math
import pandas as pd
from shapely import wkb
from shapely.geometry import box
import folium
from folium.plugins import GroupedLayerControl
import json
import warnings
warnings.filterwarnings(”ignore”)
1. Buildings from OSM
admin = ox.geocode_to_gdf(’5th District, Budapest’)
#admin = ox.geocode_to_gdf(’Siklos, Hungary’)
geom = admin.geometry.to_list()[0]
buildings_osm = ox.features_from_polygon(geom, tags={”building”: True})
print(len(buildings_osm))
buildings_osm.plot()2. Global Building Atlas
Full credits: https://mediatum.ub.tum.de/1782307
Full credits: https://huggingface.co/datasets/zhu-xlab/GBA.ODbLPolygon
Full technical guide: https://github.com/zhu-xlab/GlobalBuildingAtlas
Mirror source: https://source.coop/tge-labs/globalbuildingatlas-lod1
def get_tiles_from_geometry(geom):
“”“
Given a shapely geometry, returns all GBA tile names that intersect it.
Works for any geometry type (Point, Polygon, MultiPolygon, etc.)
“”“
minx, miny, maxx, maxy = geom.bounds
# Find all 5° tiles that overlap the bounding box
lon_start = math.floor(minx / 5) * 5
lat_start = math.floor(miny / 5) * 5
tiles = []
lon = lon_start
while lon < maxx:
lat = lat_start
while lat < maxy:
lon_min = lon
lon_max = lon + 5
lat_min = lat
lat_max = lat + 5
ew_min = “e” if lon_min >= 0 else “w”
ew_max = “e” if lon_max >= 0 else “w”
ns_max = “n” if lat_max >= 0 else “s”
ns_min = “n” if lat_min >= 0 else “s”
tile = (
f”{ew_min}{abs(lon_min):03d}_{ns_max}{abs(lat_max):02d}”
f”_{ew_max}{abs(lon_max):03d}_{ns_min}{abs(lat_min):02d}.parquet”
)
tiles.append(tile)
lat += 5
lon += 5
return tiles
# The tiles are named by 5°x5° bounding boxes in WGS84
# Budapest is around lon=19, lat=47.5 → tile: e015_n50_e020_n45
tiles = get_tiles_from_geometry(geom)
print(tiles)def check_file_size(url):
req = urllib.request.Request(url, method=”HEAD”)
with urllib.request.urlopen(req) as response:
size_bytes = int(response.headers.get(”Content-Length”, 0))
size_mb = size_bytes / (1024 ** 2)
size_gb = size_bytes / (1024 ** 3)
print(f”URL: {url}”)
print(f”Size: {size_bytes:,} bytes ({size_mb:.1f} MB / {size_gb:.2f} GB)”)
return size_bytes
tile_name = tiles[0]
url = f”https://data.source.coop/tge-labs/globalbuildingatlas-lod1/{tile_name}”
check_file_size(url)def download_gba_tile(tile_name):
url = f”https://data.source.coop/tge-labs/globalbuildingatlas-lod1/{tile_name}”
# check size first
req = urllib.request.Request(url, headers={”Range”: “bytes=0-0”})
with urllib.request.urlopen(req) as r:
size_bytes = int(r.headers.get(”Content-Range”, “0/0”).split(”/”)[-1])
print(f”Size: {size_bytes/1024**3:.2f} GB”)
confirm = input(”Proceed? (y/n): “)
if confirm.lower() != “y”:
print(”Cancelled.”)
return None
urllib.request.urlretrieve(url, tile_name)
print(f”Saved: {tile_name}”)
return tile_name
# usage
download_gba_tile(tiles[0])# Load it
df = pd.read_parquet(tile_name)
df.head(3)df[”geometry”] = df[”geometry”].apply(lambda x: wkb.loads(x))
gdf = gpd.GeoDataFrame(df, geometry=”geometry”, crs=”EPSG:4326”)
gdf.head()buildings_gba = gdf[gdf.geometry.intersects(geom)]
len(buildings_gba)buildings_gba.plot()3. Comparison
print(”OSM buildings:”, len(buildings_osm))
print(”GBA buildings:”, len(buildings_gba))
osm_total_area = buildings_osm.to_crs(”EPSG:3857”).geometry.area.sum()
gba_total_area = buildings_gba.to_crs(”EPSG:3857”).geometry.area.sum()
print(f”OSM total footprint: {osm_total_area/1e6:.4f} km²”)
print(f”GBA total footprint: {gba_total_area/1e6:.4f} km²”)
print(f”GBA covers {100*gba_total_area/osm_total_area:.1f}% of OSM footprint area”)print(buildings_gba[”source”].value_counts())
# ms = Microsoft
# goo = Google Open Buildings
# osm = OpenStreetMap
# ours2 = TUM satellite-derived
# OSM typically stores height as a string tag — convert first
buildings_osm[”height_m”] = pd.to_numeric(buildings_osm.get(”height”), errors=”coerce”)
# GBA height is already numeric
buildings_gba[”height_m”] = pd.to_numeric(buildings_gba[”height”], errors=”coerce”)
# Compare
osm_with_height = buildings_osm[”height_m”].notna().sum()
gba_with_height = buildings_gba[”height_m”].notna() & (buildings_gba[”height_m”] > 0)
gba_with_height = gba_with_height.sum()
print(f”OSM buildings with height: {osm_with_height} / {len(buildings_osm)} ({100*osm_with_height/len(buildings_osm):.1f}%)”)
print(f”GBA buildings with height: {gba_with_height} / {len(buildings_gba)} ({100*gba_with_height/len(buildings_gba):.1f}%)”)# prep
buildings_osm = buildings_osm.to_crs(”EPSG:4326”)
buildings_gba = buildings_gba.to_crs(”EPSG:4326”)
# center
center_point = buildings_gba.to_crs(”EPSG:3857”).geometry.centroid.to_crs(”EPSG:4326”)
center = [center_point.y.mean(), center_point.x.mean()]
# map
m = folium.Map(
location=center,
zoom_start=15,
tiles=”https://server.arcgisonline.com/ArcGIS/rest/services/World_Imagery/MapServer/tile/{z}/{y}/{x}”,
attr=”Esri World Imagery”
)
# OSM layer — blue with height, grey without
osm_layer = folium.FeatureGroup(name=”🔵 OSM Buildings”, show=True)
for _, row in buildings_osm.iterrows():
has_height = pd.notna(row[”height_m”]) and row[”height_m”] > 0
color = “#4A90D9” if has_height else “#888888”
opacity = 0.7 if has_height else 0.2
tooltip = f”Height: {row[’height_m’]}m” if has_height else “No height data”
folium.GeoJson(
row.geometry.__geo_interface__,
style_function=lambda x, c=color, o=opacity: {
“fillColor”: c, “color”: c, “weight”: 1, “fillOpacity”: o
},
tooltip=tooltip
).add_to(osm_layer)
osm_layer.add_to(m)
# GBA layer — orange with height, grey without
gba_layer = folium.FeatureGroup(name=”🟠 GBA Buildings”, show=True)
for _, row in buildings_gba.iterrows():
has_height = pd.notna(row[”height_m”]) and row[”height_m”] > 0
color = “#E8813A” if has_height else “#888888”
opacity = 0.7 if has_height else 0.2
tooltip = f”Height: {row[’height_m’]}m | Source: {row.get(’source’, ‘N/A’)}” if has_height else f”No height | Source: {row.get(’source’, ‘N/A’)}”
folium.GeoJson(
row.geometry.__geo_interface__,
style_function=lambda x, c=color, o=opacity: {
“fillColor”: c, “color”: c, “weight”: 1, “fillOpacity”: o
},
tooltip=tooltip
).add_to(gba_layer)
gba_layer.add_to(m)
# controls + legend
folium.LayerControl(collapsed=False, show_base_layers=False).add_to(m)
m.get_root().html.add_child(folium.Element(”“”
<div style=”position:fixed; bottom:30px; left:30px; z-index:1000;
background:rgba(20,20,20,0.85); padding:12px 18px; border-radius:8px;
color:white; font-family:monospace; font-size:13px; border:1px solid #444;”>
<b>Height Data Coverage</b><br><br>
<span style=”color:#4A90D9”>■</span> OSM — has height<br>
<span style=”color:#E8813A”>■</span> GBA — has height<br>
<span style=”color:#888888”>■</span> No height data<br><br>
<span style=”font-size:11px; color:#aaa”>Toggle layers top-right</span>
</div>
“”“))
m