Building upon our previous explorations of Python, Shapely, and GeoPandas, this lesson introduces the broader Python geospatial ecosystem. We will focus on two particularly useful libraries: osmnx for acquiring urban data from OpenStreetMap, and momepy for conducting urban morphological analysis.
Ensure you have osmnx and momepy installed in your Python environment. If you are following along with a notebook, you can install these using pip.
# Use an exclamation mark to run in notebooks
!pip install osmnx momepyPackages such as matplotlib and geopandas are often already installed by default but, if necessary, you can install them in the same way.
Let’s import the necessary packages:
import geopandas as gpd
# use an alias for convenience
import osmnx as ox
import momepy
import matplotlib.pyplot as plt/home/runner/work/cityseer-examples/cityseer-examples/.venv/lib/python3.13/site-packages/tqdm/auto.py:21: TqdmWarning: IProgress not found. Please update jupyter and ipywidgets. See https://ipywidgets.readthedocs.io/en/stable/user_install.html
from .autonotebook import tqdm as notebook_tqdmosmnx is a Python package that simplifies downloading and working with geospatial data from OpenStreetMap (OSM), such as building footprints or points of interest.
This example demonstrates how to download building footprints within a 1km radius of a specified location in Nicosia, Cyprus, defined by its coordinates. osmnx provides several methods for downloading data; here, we will use features_from_point.
# Define the point of interest (latitude, longitude) and distance
center_point = (35.17526, 33.36402)
distance_m = 1000
# Download building footprints
gdf_buildings = ox.features_from_point(
center_point, tags={"building": True}, dist=distance_m
)osmnx returns GeoDataFrames which, as shown in the previous lesson, are ideal for spatial analysis in Python. Note the tags={"building": True} argument, which instructs osmnx to fetch all features tagged as buildings in OSM. By changing these tags, you can also download other types of features, such as roads or parks.
It is good practice to inspect the data you have downloaded. The head() function displays the first few rows of the GeoDataFrame, allowing you to quickly check the data structure and attributes.
# Display the first few rows of the buildings GeoDataFrame
gdf_buildings.head()| geometry | amenity | building | bus | name | public_transport | addr:city | addr:housenumber | addr:postcode | addr:street | ... | indoor_seating | tomb | name:hy | wikipedia:hy | construction | theatre:type | association | house | type | inscription | ||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| element | id | |||||||||||||||||||||
| node | 4338033483 | POINT (33.36625 35.18368) | bus_station | yes | yes | İtimat (Lefkoşa Mağusa Servis) | station | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN |
| 5206675859 | POINT (33.35992 35.17344) | NaN | yes | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |
| 12027040593 | POINT (33.35773 35.17435) | NaN | yes | NaN | Apostolic Nunciature to Cyprus | NaN | Nicosia | 1 | 1010 | Pafou | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | |
| relation | 2580980 | POLYGON ((33.36224 35.17641, 33.3628 35.17652,... | NaN | yes | NaN | Büyük Han | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | multipolygon | OTTOMAN |
| 2785751 | POLYGON ((33.36028 35.17739, 33.36019 35.17701... | NaN | yes | NaN | NaN | NaN | NaN | NaN | NaN | Müftü Raci Sokak | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | multipolygon | NaN |
5 rows × 170 columns
You can also plot the downloaded data:
# Set up a plot
fig, ax = plt.subplots(figsize=(10, 10))
# Plot the buildings
gdf_buildings.plot(ax=ax)
# Set the title and remove the axis for a cleaner look
ax.set_title(f"Buildings around {center_point} ({distance_m}m radius)")
ax.axis("off")(np.float64(33.351150465),
np.float64(33.376499635),
np.float64(35.164918945000004),
np.float64(35.185606754999995))
For more detailed information on different ways to query OSM features data, refer to the osmnx features documentation.
To streamline the subsequent analysis, it is advisable to first filter for the types of geometry you intend to work with. In this instance, we are interested in polygon or multi-polygon geometries and will discard other types, such as points or linestrings.
# Filter out non-polygon geometries
gdf_buildings = gdf_buildings[
gdf_buildings.geometry.type.isin(['Polygon', 'MultiPolygon'])
]
gdf_buildings.head()| geometry | amenity | building | bus | name | public_transport | addr:city | addr:housenumber | addr:postcode | addr:street | ... | indoor_seating | tomb | name:hy | wikipedia:hy | construction | theatre:type | association | house | type | inscription | ||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| element | id | |||||||||||||||||||||
| relation | 2580980 | POLYGON ((33.36224 35.17641, 33.3628 35.17652,... | NaN | yes | NaN | Büyük Han | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | multipolygon | OTTOMAN |
| 2785751 | POLYGON ((33.36028 35.17739, 33.36019 35.17701... | NaN | yes | NaN | NaN | NaN | NaN | NaN | NaN | Müftü Raci Sokak | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | multipolygon | NaN | |
| 3403727 | POLYGON ((33.36277 35.17718, 33.36274 35.17688... | NaN | yes | NaN | Kumarcılar Hanı | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | multipolygon | NaN | |
| 8756098 | POLYGON ((33.36483 35.17161, 33.36483 35.17161... | NaN | yes | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | multipolygon | NaN | |
| 8762332 | POLYGON ((33.36003 35.17852, 33.36002 35.17849... | NaN | house | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | multipolygon | NaN |
5 rows × 170 columns
Secondly, we will reset the index so that all features are neatly indexed from zero upwards, without duplicates.
# Reset the index
gdf_buildings = gdf_buildings.reset_index(drop=True)
gdf_buildings.head()| geometry | amenity | building | bus | name | public_transport | addr:city | addr:housenumber | addr:postcode | addr:street | ... | indoor_seating | tomb | name:hy | wikipedia:hy | construction | theatre:type | association | house | type | inscription | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | POLYGON ((33.36224 35.17641, 33.3628 35.17652,... | NaN | yes | NaN | Büyük Han | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | multipolygon | OTTOMAN |
| 1 | POLYGON ((33.36028 35.17739, 33.36019 35.17701... | NaN | yes | NaN | NaN | NaN | NaN | NaN | NaN | Müftü Raci Sokak | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | multipolygon | NaN |
| 2 | POLYGON ((33.36277 35.17718, 33.36274 35.17688... | NaN | yes | NaN | Kumarcılar Hanı | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | multipolygon | NaN |
| 3 | POLYGON ((33.36483 35.17161, 33.36483 35.17161... | NaN | yes | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | multipolygon | NaN |
| 4 | POLYGON ((33.36003 35.17852, 33.36002 35.17849... | NaN | house | NaN | NaN | NaN | NaN | NaN | NaN | NaN | ... | NaN | NaN | NaN | NaN | NaN | NaN | NaN | NaN | multipolygon | NaN |
5 rows × 170 columns
Thirdly, we will drop any columns not relevant to our analysis. In this case, we will retain the geometry column and the building column.
gdf_buildings = gdf_buildings[['geometry', 'building']]
gdf_buildings.head()| geometry | building | |
|---|---|---|
| 0 | POLYGON ((33.36224 35.17641, 33.3628 35.17652,... | yes |
| 1 | POLYGON ((33.36028 35.17739, 33.36019 35.17701... | yes |
| 2 | POLYGON ((33.36277 35.17718, 33.36274 35.17688... | yes |
| 3 | POLYGON ((33.36483 35.17161, 33.36483 35.17161... | yes |
| 4 | POLYGON ((33.36003 35.17852, 33.36002 35.17849... | house |
Before performing morphological analysis, it is necessary to ensure your data is in a projected Coordinate Reference System (CRS). Morphological metrics often involve measurements of distance and area, which are only accurate in a projected CRS. For Nicosia, we will use EPSG:3035 (ETRS89 / LAEA Europe), a European projection.
# Set the target CRS
gdf_buildings_proj = gdf_buildings.to_crs(3035)Let’s replot the building data to ensure everything is in order.
# Set up a plot
fig, ax = plt.subplots(figsize=(10, 10))
# Plot the buildings
gdf_buildings_proj.plot(ax=ax)
# Set the title and remove axis
ax.set_title(f"Buildings around {center_point} ({distance_m}m radius)")
ax.axis("off")(np.float64(6433775.66431334),
np.float64(6436461.588935598),
np.float64(1669215.3095963916),
np.float64(1672018.7854921184))
momepy is a library for the quantitative analysis of urban form – also known as urban morphology. It operates primarily on GeoDataFrames and provides a range of functions for calculating various morphological metrics.
By way of example, we will explore two of these functions.
momepy can calculate the orientation of building footprints using the orientation function.
# Calculate building orientation
gdf_buildings_proj['orientation'] = momepy.orientation(gdf_buildings_proj)
fig, ax = plt.subplots(figsize=(10, 10))
gdf_buildings_proj.plot(ax=ax, column='orientation', cmap='Spectral')
ax.set_title(f"Building orientations")
ax.axis("off")(np.float64(6433775.66431334),
np.float64(6436461.588935598),
np.float64(1669215.3095963916),
np.float64(1672018.7854921184))
momepy can calculate shared wall lengths between buildings using the shared_walls function.
# Calculate shared wall lengths
gdf_buildings_proj['shared_wall_length'] = momepy.shared_walls(gdf_buildings_proj)
fig, ax = plt.subplots(figsize=(10, 10))
gdf_buildings_proj.plot(ax=ax, column='shared_wall_length', cmap='hot')
ax.set_title(f"Building shared wall lengths")
ax.axis("off")(np.float64(6433775.66431334),
np.float64(6436461.588935598),
np.float64(1669215.3095963916),
np.float64(1672018.7854921184))
Many other functions are available in momepy for calculating various morphological metrics. For a comprehensive list, please refer to the momepy documentation.
This has been a brief exploration of the broader Python geospatial ecosystem, focusing on osmnx for downloading urban data from OpenStreetMap and momepy for conducting urban morphology analysis. This is just a small sample of what these and other tools can achieve.
Remember that with GeoPandas, you can export your data. So, after downloading data and running any variety of metrics using available Python packages, you can then export the file, which can be further visualised or manipulated with packages such as QGIS. For example, you can export your data to a GeoPackage file using the to_file method:
gdf_buildings_proj.to_file("nicosia_buildings_metrics.gpkg", driver="GPKG")