Видимость 2D: полигоны видимости, метрики изовист и охват

get_visibility вычисляет точный полигон видимости (изовисту) каждого наблюдателя в пределах дальности обзора; препятствия — стены бесконечной высоты. Этот ноутбук показывает

  1. полигон видимости одного наблюдателя и нескольких наблюдателей с индивидуальными радиусами;

  2. метрики изовист (isovist_metrics=True) и поле изовист по территории (get_isovist_grid, rasterize);

  3. охват видимости — сколько наблюдателей видят каждую ячейку (get_visibility_catchment).

Области видимости 2.5D с учётом высот зданий разобраны в ноутбуке областей видимости.

Последний раздел, Своя территория, показывает, что заменить, чтобы получить такие же картинки для другой территории.

# !pip install objectnat pyarrow matplotlib
import geopandas as gpd
import shapely
from shapely.geometry import Point

from objectnat import get_visibility, get_visibility_catchment

obstacles = gpd.read_parquet("examples_data/buildings.parquet")  # EPSG:32636, metres
import matplotlib.pyplot as plt
import numpy as np
from mpl_toolkits.axes_grid1.anchored_artists import AnchoredSizeBar


def map_extent(ax, bounds, bar=100, outline=True):
    """Show `bounds` (x0, y0, x1, y1 in metres) without coordinate axes, with a scale bar and the area outline."""
    x0, y0, x1, y1 = bounds
    pad = 0.012 * max(x1 - x0, y1 - y0)
    ax.set_xlim(x0 - pad, x1 + pad)
    ax.set_ylim(y0 - pad, y1 + pad)
    ax.set_aspect("equal")
    ax.set_axis_off()
    if outline:
        ax.plot([x0, x1, x1, x0, x0], [y0, y0, y1, y1, y0], color="#3d3d3d", linewidth=1, zorder=8)
    scale = AnchoredSizeBar(ax.transData, bar, f"{bar:g} m", loc="lower right", pad=0.4, sep=4, size_vertical=bar / 40)
    scale.patch.set(facecolor="white", edgecolor="none", alpha=0.85)
    ax.add_artist(scale)
from matplotlib.patches import Patch

1. Один наблюдатель

Наблюдатели могут быть в любой системе координат; расчёт идёт в оценённой зоне UTM, а результат возвращается в системе координат наблюдателей.

observer = gpd.GeoDataFrame(geometry=[Point(30.2312112, 59.9482336)], crs=4326)
visible = get_visibility(observer, obstacles, view_distance=800)
local = visible.to_crs(obstacles.crs)
print(f"visible area {local.area.iloc[0] / 1e4:.1f} ha of {3.14159 * 800 ** 2 / 1e4:.1f} ha")

c = observer.to_crs(obstacles.crs).geometry.iloc[0]
bounds = (c.x - 800, c.y - 800, c.x + 800, c.y + 800)
fig, ax = plt.subplots(figsize=(7, 7))
obstacles.plot(ax=ax, color="#bdbdbd", edgecolor="#8c8c8c", linewidth=0.3)
local.plot(ax=ax, color="#2b7bba", alpha=0.55, linewidth=0)
observer.to_crs(obstacles.crs).plot(ax=ax, color="#d7301f", markersize=30, edgecolor="white")
map_extent(ax, bounds, bar=200, outline=False)
ax.set_title("visibility polygon within 800 m")
plt.show()
visible area 39.0 ha of 201.1 ha
../../_images/3a899ca7b7b318ea5902a0d61e6f521d7d17bde87fc692247dc1e78d6ccc8f13.png

2. Несколько наблюдателей

Столбец visibility_distance задаёт радиус каждого наблюдателя (метры). Наблюдатели обрабатываются параллельно расширением на Rust; наблюдатели внутри препятствия получают пустой полигон.

points = gpd.GeoDataFrame(
    {"visibility_distance": [300, 500, 800]},
    geometry=[Point(30.2312112, 59.9482336), Point(30.2390, 59.9465), Point(30.2250, 59.9500)],
    crs=4326,
)
many = get_visibility(points, obstacles)
many.to_crs(obstacles.crs).assign(area_ha=lambda g: g.area / 1e4)[["visibility_distance", "area_ha"]].round(2)
visibility_distance area_ha
0 300 15.54
1 500 21.70
2 800 1.87

3. Метрики изовист и поля изовист

isovist_metrics=True описывает каждый полигон видимости числами; длины — в метрах, площади — в м². Страница видимости объясняет каждую метрику на примерах. Кратко:

  • area и visible_share: сколько видит наблюдатель, в м² и как долю круга обзора;

  • perimeter и compactness (4πA/P²): длина границы и показатель от 1 для полного круга почти до 0 для звездообразной формы с длинными узкими взглядами;

  • occlusivity и n_occlusions: суммарная длина и число рёбер, на которых взгляд уходит за угол здания, то есть сколько пространства скрыто прямо за углами;

  • min_radial, max_radial, mean_radial, std_radial и skew_radial: статистики расстояния от наблюдателя до границы обзора по всем направлениям;

  • drift: расстояние от наблюдателя до центроида видимой области, то есть насколько обзор смещён в одну сторону;

  • open_share: доля границы, лежащая на границе дальности обзора, то есть не на стенах.

get_isovist_grid строит карту метрик по территории: функция ставит наблюдателя в центр каждой ячейки размером cell_size на открытой территории и возвращает наблюдателей с их метриками. Здания за пределами территории, но в пределах view_distance, тоже закрывают обзор. rasterize превращает метрику в Raster, который plot рисует в matplotlib; to_geotiff и to_folium записывают его в файл или на веб-карту (см. ноутбук инсоляции).

import time

from objectnat import get_isovist_grid, rasterize

centre = (346050, 6647450)
square = shapely.box(centre[0] - 250, centre[1] - 250, centre[0] + 250, centre[1] + 250)
scene = obstacles[obstacles.intersects(square.buffer(200))]

start = time.perf_counter()
field = get_isovist_grid(obstacles, square, view_distance=200, cell_size=5)
print(f"{len(field)} observers in {time.perf_counter() - start:.1f} s")
field.drop(columns="geometry").describe().T[["mean", "min", "50%", "max"]].round(2)
7889 observers in 0.6 s
mean min 50% max
area 26285.21 743.01 19527.09 70462.60
perimeter 1468.42 128.57 1490.88 2729.06
compactness 0.15 0.03 0.13 0.66
occlusivity 810.70 73.71 784.69 1823.30
n_occlusions 17.21 2.00 17.00 36.00
min_radial 14.74 0.00 11.42 78.83
max_radial 196.79 46.41 200.00 200.00
mean_radial 68.15 6.71 60.02 141.53
std_radial 51.42 12.46 50.24 92.14
skew_radial 1.16 -0.52 1.20 6.33
drift 39.68 0.22 37.42 116.97
visible_share 0.21 0.01 0.16 0.56
open_share 0.07 0.00 0.04 0.40
fig, axs = plt.subplots(2, 2, figsize=(11, 10), constrained_layout=True)
for ax, (col, label) in zip(axs.ravel(), [
    ("visible_share", "visible share of the 200 m disc"),
    ("compactness", "4πA/P²"),
    ("occlusivity", "length of occluding edges, m"),
    ("drift", "distance to the centroid, m"),
]):
    raster = rasterize(field, col)
    lo, hi = np.nanpercentile(raster.values, [2, 98])
    gpd.GeoSeries([square]).plot(ax=ax, color="#bdbdbd", zorder=0)  # cells without an observer lie under buildings
    im = raster.plot(ax=ax, cmap="viridis", vmin=lo, vmax=hi)
    scene.clip(square).plot(ax=ax, color="#bdbdbd", edgecolor="#8c8c8c", linewidth=0.3)
    map_extent(ax, square.bounds)
    ax.set_title(col)
    fig.colorbar(im, ax=ax, shrink=0.8, label=label)
plt.show()
../../_images/e2f5e95e983125f7612a6c39f01389c7d57d2c58a34d5dfc1e3cb40096987075.png

Каждая ячейка показывает, что видит наблюдатель, стоящий в её центре. У открытых площадей и широких улиц высокий visible_share, у узких дворов — низкий. compactness высока там, где видимая область простой, почти выпуклой формы, — на открытой площади или внутри закрытого двора, и низка там, где длинные узкие взгляды между зданиями делают её звездообразной. occlusivity высока там, где много углов зданий скрывают пространство рядом с наблюдателем, — на пересечениях проходов между зданиями. drift велик у стен, обращённых к открытому пространству: наблюдатель видит в основном одну сторону.

4. Охват видимости

get_visibility_catchment подсчитывает для каждой ячейки размером cell_size наблюдателей, полигон видимости которых содержит центр ячейки, и возвращает видимые ячейки с их count; rasterize их рисует. Наблюдателей вдоль улиц обычно расставляют функцией get_observer_points(urban_graph, spacing=15); например, с пешеходным графом из IduEdu:

import iduedu
from objectnat import get_observer_points

with iduedu.use_pbf(osm_id=1114252, region="northwestern-fed-district"):
    walk = iduedu.get_walk_graph(osm_id=1114252)
observers = get_observer_points(walk, spacing=15)

Здесь, чтобы не зависеть от загрузок, наблюдатели размещены через 15 м на открытой территории не дальше 12 м от зданий, по территории, расширенной на дальность обзора (250 м): ячейку у края территории видят и наблюдатели за её пределами, и без них её счёт получился бы заниженным.

# observers over the area extended by the view distance; buildings far enough around them
near = obstacles[obstacles.intersects(square.buffer(500))]
footprints = shapely.union_all(near.geometry.values)
x0, y0, x1, y1 = square.buffer(250, join_style="mitre").bounds
xx, yy = np.meshgrid(np.arange(x0, x1, 15), np.arange(y0, y1, 15))
pts = shapely.points(xx.ravel(), yy.ravel())
pts = pts[~shapely.intersects(footprints.buffer(2), pts) & shapely.intersects(footprints.buffer(12), pts)]
street_observers = gpd.GeoDataFrame(geometry=pts, crs=obstacles.crs)

catchment = get_visibility_catchment(street_observers, near, view_distance=250, cell_size=4)
print(len(street_observers), "observers,", len(catchment), "cells seen")
catchment["count"].describe().round(1)
1547 observers, 83357 cells seen
count    83357.0
mean        25.0
std         18.8
min          1.0
25%          9.0
50%         21.0
75%         39.0
max         91.0
Name: count, dtype: float64
seen = rasterize(catchment, "count")
fig, ax = plt.subplots(figsize=(7.5, 6.5), constrained_layout=True)
im = seen.plot(ax=ax, cmap="viridis", vmin=1, vmax=np.nanpercentile(seen.values, 98))
scene.clip(square).plot(ax=ax, color="#f0f0f0", edgecolor="#8c8c8c", linewidth=0.3)
street_observers.plot(ax=ax, color="#d7301f", markersize=1.5)
map_extent(ax, square.bounds)
fig.colorbar(im, ax=ax, shrink=0.8, extend="max", label="observers that see the cell")
ax.set_title("visibility catchment")
plt.show()
../../_images/33c9dfe38d2a56ee75f888a0c2adea21cd7a63a2348ae1a3b4651176efc78736.png

Значения выражены в числе наблюдателей, поэтому растры разных территорий или сценариев с одинаковым шагом наблюдателей можно рисовать в одной цветовой шкале. weight_column заменяет число наблюдателей суммой их весов, например пешеходного потока.

Метрики изовист для своих наблюдателей

get_isovist_grid ставит наблюдателей по сетке. Для наблюдателей, расставленных иначе, например вдоль улиц, get_visibility(observers, obstacles, view_distance, isovist_metrics=True) возвращает те же метрики, по строке на наблюдателя; рисуйте их точками. В отличие от счёта охвата, у этих значений нет краевого эффекта: каждое описывает своего наблюдателя.

street_metrics = get_visibility(street_observers, obstacles, view_distance=200, isovist_metrics=True)
street_points = street_metrics.set_geometry(street_observers.geometry)

fig, ax = plt.subplots(figsize=(7, 6))
scene.clip(square).plot(ax=ax, color="#bdbdbd", linewidth=0)
street_points.plot(ax=ax, column="visible_share", cmap="viridis", markersize=8, legend=True,
                   legend_kwds={"label": "visible share of the 200 m disc", "shrink": 0.8})
map_extent(ax, square.bounds)
plt.show()
../../_images/ca44a314a007cb768cd99451c0e59c295173c1fd4dafebd4640267b1760ac190.png

5. Своя территория

Наблюдатели через каждые 20 м вдоль пешеходного графа района из IduEdu (get_observer_points); здания служат только препятствиями и наблюдателей не добавляют. Около 44 тысяч наблюдателей рассчитываются за несколько секунд. Растр покрывает всё, что видят наблюдатели, в том числе за границей района. Значения зависят от шага наблюдателей, поэтому сравнивайте растры, построенные с одинаковым шагом. Raster.plot рисует растр, Raster.explore показывает его на веб-карте с цветовой шкалой, а Raster.to_geotiff записывает GeoTIFF для ГИС (нужен rasterio).

Охват и visible_share отвечают на разные вопросы. visible_share (метрика изовисты, раздел 3) смотрит из места наружу: какую долю круга обзора видит наблюдатель, стоящий в нём. Охват смотрит на место со стороны всех наблюдателей: сколько из них его видят, то есть насколько оно заметно с улиц. У середины большой площади visible_share высокий, а охват низкий, если к ней ведёт мало улиц. Поле изовист из раздела 3 строится для любой территории так же, например rasterize(get_isovist_grid(my_obstacles, my_area, view_distance=200, cell_size=5), "visible_share").

import iduedu

from objectnat import get_observer_points

OSM_ID, REGION = 1114252, "northwestern-fed-district"  # your territory
with iduedu.use_pbf(osm_id=OSM_ID, region=REGION):  # downloaded once, then cached
    my_boundary = iduedu.get_4326_boundary(osm_id=OSM_ID)
    my_walk = iduedu.get_walk_graph(territory=my_boundary, clip_by_territory=True)
my_observers = get_observer_points(my_walk, spacing=20)  # every 20 m along the pedestrian graph
my_obstacles = obstacles  # your building footprints, any CRS

start = time.perf_counter()
cells = get_visibility_catchment(my_observers, my_obstacles, view_distance=250, cell_size=5)
print(f"{len(my_observers)} observers, {len(cells)} cells seen, {time.perf_counter() - start:.1f} s")

seen = rasterize(cells, "count")
fig, ax = plt.subplots(figsize=(11, 7.5))
im = seen.plot(ax=ax, cmap="viridis", vmin=1, vmax=300)
my_obstacles.to_crs(seen.crs).plot(ax=ax, color="#bdbdbd", linewidth=0)
map_extent(ax, seen.bounds, bar=500, outline=False)
fig.colorbar(im, ax=ax, shrink=0.7, extend="max", label="observers that see the cell")
plt.show()
# seen.to_geotiff("my_visibility_count.tif")  # for GIS
2026-11-06 20:43:45.799 | WARNING  | Removing 336 nodes outside the largest connected component. Retaining 19802 of 20138 nodes.
43799 observers, 782044 cells seen, 1.2 s
../../_images/cb440107d3879a3327b9b1d8d98417c1461352ed48211b3c0b7962501327644f.png
seen.explore(cmap="viridis", vmin=1, vmax=300, name="observers that see the cell")
Make this Notebook Trusted to load map: File -> Trust Notebook