"""Clusters of wind farms and the placement of their boundaries."""
from __future__ import annotations
import itertools
import json
from collections.abc import Iterator, Mapping, Sequence
from pathlib import Path
from typing import TYPE_CHECKING, Any, overload
import numpy as np
import shapely
from shapely import affinity
from shapely.ops import unary_union
from ._random import as_rng, derive_seed
from ._sampling_specs import _json_safe, draw_scalar, is_number
from .boundaries.base import Boundary, _cells, _GeometryBoundary, _translation
from .layouts.base import Farm
if TYPE_CHECKING: # pragma: no cover
import matplotlib.axes
__all__ = ["Cluster"]
def _geometry_radius(geometry) -> float:
"""Largest distance from the centroid of ``geometry`` to its vertices."""
center = np.asarray(geometry.centroid.coords[0], dtype=float)
coordinates = shapely.get_coordinates(geometry)
return float(np.linalg.norm(coordinates - center, axis=1).max())
def _circle_intersections(first, second, first_radius, second_radius):
"""The intersection points of two circles, none when they do not meet."""
delta = second - first
distance = float(np.linalg.norm(delta))
if distance == 0 or distance > first_radius + second_radius:
return []
if distance < abs(first_radius - second_radius):
return []
along = (first_radius**2 - second_radius**2 + distance**2) / (2 * distance)
height_sq = first_radius**2 - along**2
if height_sq < -1e-9 * max(first_radius**2, 1.0):
return []
height = np.sqrt(max(height_sq, 0.0))
direction = delta / distance
foot = first + along * direction
normal = np.asarray([-direction[1], direction[0]])
return [foot + height * normal, foot - height * normal]
def _translated_geometry(geometry, center, centroid=None):
"""``geometry`` moved so its centroid sits at ``center``.
Callers that move one geometry many times pass its ``centroid`` once.
"""
if centroid is None:
centroid = np.asarray(geometry.centroid.coords[0], dtype=float)
offset = np.asarray(center, dtype=float) - centroid
return affinity.translate(geometry, xoff=float(offset[0]), yoff=float(offset[1]))
def _has_clearance(geometry, others, gap: float, area_tolerance: float) -> bool:
"""Whether ``geometry`` overlaps none of ``others`` and stays ``gap`` away."""
for other in others:
distance = geometry.distance(other)
if distance < gap:
return False
# Geometries any distance apart share no area.
if distance == 0 and _shared_area(geometry, other) > area_tolerance:
return False
return True
#: Share of a pair's span the overlay of :func:`_shared_area` snaps to.
_OVERLAY_GRID = 1e-12
def _shared_area(first, second) -> float:
"""The area ``first`` and ``second`` share, measured robustly.
A floating overlay took two touching Channels cells, one's corner on the
other's edge, for overlapping ones and returned most of a cell, as GEOS
3.14.1 did. Snapped to a grid far finer than any tolerance, the overlay
is robust and finds them apart.
"""
minx, miny, maxx, maxy = shapely.total_bounds([first, second])
grid = _OVERLAY_GRID * max(maxx - minx, maxy - miny)
return float(shapely.intersection(first, second, grid_size=grid).area)
class _PlacedGeometries:
"""Placed geometries with their bounding boxes for quick clearance tests."""
def __init__(self, geometries):
self.geometries = list(geometries)
self._bounds = np.empty((0, 4), dtype=float)
if self.geometries:
self._bounds = shapely.bounds(np.asarray(self.geometries, dtype=object))
def append(self, geometry) -> None:
"""Add ``geometry`` to the placed ones."""
self.geometries.append(geometry)
self._bounds = np.vstack([self._bounds, shapely.bounds(geometry)])
def near(self, geometry, gap: float) -> list:
"""The placed geometries whose bounding boxes come within ``gap`` of it.
Every other placed geometry lies farther than ``gap`` from ``geometry``
and cannot overlap it, so it passes any clearance test. They are
returned in placement order.
"""
bounds = shapely.bounds(geometry)
dx = np.maximum(self._bounds[:, 0] - bounds[2], bounds[0] - self._bounds[:, 2])
dy = np.maximum(self._bounds[:, 1] - bounds[3], bounds[1] - self._bounds[:, 3])
separation = np.hypot(np.maximum(dx, 0.0), np.maximum(dy, 0.0))
# The margin absorbs rounding here and in the exact distances.
margin = 1e-9 * (gap + float(np.abs(bounds).max()) + 1.0)
close = np.flatnonzero(separation <= gap + margin)
return [self.geometries[index] for index in close]
def _pull_towards_cluster(
geometry, center, target, placed, gap: float, reach: float = np.inf
):
"""Slide ``geometry`` from ``center`` towards ``target`` while it keeps clear.
Bisects the straight path and returns the moved geometry and its center,
or None as soon as the center is sure to stop farther than ``reach`` from
``target``.
"""
center = np.asarray(center, dtype=float)
travel = float(np.linalg.norm(target - center))
scale = max(travel, _geometry_radius(geometry), 1.0)
area_tolerance = 1e-12 * scale**2
centroid = np.asarray(geometry.centroid.coords[0], dtype=float)
# The margin absorbs rounding in the distance the center stops at.
margin = 1e-9 * (travel + np.abs(center).max() + np.abs(target).max() + 1.0)
def moved(position):
return _translated_geometry(geometry, position, centroid)
def clear(candidate):
others = placed.near(candidate, gap)
return _has_clearance(candidate, others, gap, area_tolerance)
target_geometry = moved(target)
if clear(target_geometry):
return target_geometry, target
low = 0.0
high = 1.0
while (high - low) * travel > 1e-9 * scale:
fraction = (low + high) / 2.0
candidate = moved(center + fraction * (target - center))
if clear(candidate):
low = fraction
else:
high = fraction
# The center stops short of ``high``, so it ends farther away.
if (1.0 - high) * travel - margin > reach:
return None
resolved_center = center + low * (target - center)
return moved(resolved_center), resolved_center
def _candidate_centers(disks, radius: float, gap: float, angular_offset: float):
"""Centers where a disk of ``radius`` touches the placed disks ``gap`` apart.
They lie in six directions around each placed disk, rotated by
``angular_offset``, and where the disk touches two placed disks at once.
"""
candidates = []
angles = angular_offset + np.deg2rad(np.arange(0.0, 360.0, 60.0))
directions = np.column_stack([np.cos(angles), np.sin(angles)])
for center, other_radius in disks:
distance = radius + other_radius + gap
candidates.extend(center + distance * directions)
if len(disks) < 2:
return candidates
# Only pairs of circles that come near each other can meet. The slack
# leaves borderline pairs to the exact test.
centers = np.asarray([center for center, _ in disks], dtype=float)
reach = radius + np.asarray([other for _, other in disks], dtype=float) + gap
delta = centers[:, None, :] - centers[None, :, :]
distances = np.hypot(delta[..., 0], delta[..., 1])
total = reach[:, None] + reach[None, :]
slack = 1e-6 * total
meet = (distances <= total + slack) & (
distances >= np.abs(reach[:, None] - reach[None, :]) - slack
)
for index, other_index in zip(*np.nonzero(np.triu(meet, k=1)), strict=True):
first, first_radius = disks[index]
second, second_radius = disks[other_index]
candidates.extend(
_circle_intersections(
first,
second,
radius + first_radius + gap,
radius + second_radius + gap,
)
)
return candidates
def _clear_of_disks(candidates, disks, radius: float, gap: float):
"""Whether each candidate center keeps its disk clear of the placed disks."""
centers = np.asarray([center for center, _ in disks], dtype=float)
limits = radius + np.asarray([other for _, other in disks], dtype=float) + gap
limits = limits - 1e-9
points = np.asarray(candidates, dtype=float).reshape(-1, 2)
delta = points[:, None, :] - centers[None, :, :]
distances = np.hypot(delta[..., 0], delta[..., 1])
# Distances this close to a limit are decided by the exact norm.
slack = 1e-9 * np.maximum(np.abs(limits), 1.0)
clear = ~np.any(distances < limits - slack, axis=1)
borderline = np.abs(distances - limits) <= slack
for row, column in zip(*np.nonzero(borderline & clear[:, None]), strict=True):
if not clear[row]:
continue
other_center, other_radius = disks[column]
clear[row] = (
np.linalg.norm(candidates[row] - other_center)
>= radius + other_radius + gap - 1e-9
)
return clear
def _compact_position(geometry, placed, disks, gap: float, angular_offset: float):
"""Place ``geometry`` beside the ``placed`` ones, as close to them as it fits.
Candidates are ranked by distance to the cluster center, then by the area
of the convex hull they leave, then by candidate order.
"""
radius = _geometry_radius(geometry)
centroid = np.asarray(geometry.centroid.coords[0], dtype=float)
candidates = _candidate_centers(disks, radius, gap, angular_offset)
tolerance = 1e-10 * max(geometry.area, 1.0)
placed_union = unary_union(placed.geometries)
cluster_center = np.asarray(placed_union.centroid.coords[0], dtype=float)
# The hull of the placed ones carries every vertex the cluster hull needs.
placed_hull = placed_union.convex_hull
# Quantized tie-breaks avoid platform-dependent coordinate noise.
unit = max(radius, 1.0) ** 2
starts = []
clear_of_disks = _clear_of_disks(candidates, disks, radius, gap)
for rank, center in enumerate(candidates):
if not clear_of_disks[rank]:
continue
candidate = _translated_geometry(geometry, center, centroid)
others = placed.near(candidate, gap)
if _has_clearance(candidate, others, gap, tolerance):
travel = float(np.linalg.norm(cluster_center - center))
starts.append((travel, rank, center, candidate))
# Candidates that start close tend to end close. Pulling them first lets
# the others stop as soon as they cannot rank first.
starts.sort(key=lambda item: item[:2])
valid = []
best = None
for _, rank, center, candidate in starts:
reach = np.inf
if best is not None:
# Any center farther than this has a larger quantized distance.
reach = np.sqrt((best + 1e-9) * unit) * (1.0 + 1e-9)
pulled = _pull_towards_cluster(
candidate, center, cluster_center, placed, gap, reach
)
if pulled is None:
continue
candidate, center = pulled
offset = center - cluster_center
center_distance = round(float(np.dot(offset, offset)) / unit, 9)
if best is not None and center_distance > best:
continue
best = center_distance
hull_area = shapely.convex_hull(
shapely.geometrycollections([placed_hull, candidate])
).area
valid.append((center_distance, round(hull_area / unit, 9), rank, candidate))
if not valid:
raise ValueError("could not arrange cluster boundaries without overlap")
*_, selected = min(valid, key=lambda item: item[:-1])
return selected
#: Most sites an uncapped domain fill places.
_SITE_LIMIT = 512
#: Candidate positions offered to one site before it is abandoned.
_SITE_ATTEMPTS = 96
#: Consecutive unplaceable sites after which the domain is treated as full.
_SITE_GIVE_UP = 12
#: Share of a site's inradius the domain is eroded by, so the rounded corners
#: of the eroded outline, drawn inside the exact one, keep every center the
#: site fits at.
_CENTER_EROSION = 0.98
def _center_region(developable, geometry, centroid):
"""Where the centroid of ``geometry`` may lie for it to fit in ``developable``.
A site holds the disk about its centroid out to its nearest edge, so its
centroid lies at least that far inside the domain wherever it fits.
Drawing centers there finds a large site the whole domain rarely offers
room for.
"""
point = shapely.Point(centroid)
if not geometry.contains(point):
return developable
inradius = float(geometry.boundary.distance(point))
return developable.buffer(-_CENTER_EROSION * inradius)
def _scatter_sites(make_site, domain, gap, max_sites, seed):
"""Scatter sites drawn from ``make_site(index)`` at random inside ``domain``.
Sites are added ``gap`` apart until ``max_sites`` are placed, or
:data:`_SITE_LIMIT` without a cap, or :data:`_SITE_GIVE_UP` sites in a
row find no room. A site that finds none does not count, so another is
drawn in its place. Returns ``(parts, groups, records)`` with the member
boundaries, the site group of each part and JSON-compatible placement
records.
"""
developable = domain.geometry
if developable.is_empty or developable.area <= 0:
raise ValueError("domain must enclose a positive developable area")
limit = _SITE_LIMIT if max_sites is None else int(max_sites)
tolerance = 1e-9 * developable.area
placed: list[tuple] = []
# Only sites whose boxes come within the gap can fail the clearance test,
# so the rest are not measured. Every position tried was measured against
# every placed site, four fifths of the slowest regions' time.
geometries = _PlacedGeometries([])
failures = 0
tried: list[float] = []
# Each draw places a site or adds to a run of failures that ends the fill,
# so the loop ends. A site that found no room used up one of the
# ``n_sites``, and 6 asked for placed 5 though smaller ones still fit.
for index in itertools.count():
if len(placed) >= limit:
break
site = make_site(index)
geometry = site.geometry
tried.append(float(geometry.area))
if geometry.area > developable.area:
raise ValueError("a site is larger than the domain's developable area")
site_seed = derive_seed(seed, "site", index)
rng = as_rng(site_seed)
clearance = draw_scalar(
gap, rng, "gap", index=index, choice_seed=seed, minimum=0.0
)
arranged = None
centroid = np.asarray(geometry.centroid.coords[0], dtype=float)
region = _center_region(developable, geometry, centroid)
# A site too wide for anywhere in the domain is offered no position.
attempts = 0 if region.is_empty else _SITE_ATTEMPTS
lowx, lowy, highx, highy = region.bounds
for _ in range(attempts):
center = (
float(rng.uniform(lowx, highx)),
float(rng.uniform(lowy, highy)),
)
# Checking the center first rejects most samples cheaply.
if not shapely.contains_xy(region, center[0], center[1]):
continue
candidate = _translated_geometry(geometry, center, centroid)
if not developable.covers(candidate):
continue
if not _has_clearance(
candidate, geometries.near(candidate, clearance), clearance, tolerance
):
continue
arranged = candidate
break
if arranged is None:
failures += 1
if failures >= _SITE_GIVE_UP:
break
continue
failures = 0
placed.append((index, site, arranged, clearance, site_seed))
geometries.append(arranged)
parts: list = []
groups: list[int] = []
records: list[dict] = []
for group, (index, site, arranged, clearance, site_seed) in enumerate(placed):
offset = np.asarray(arranged.centroid.coords[0], dtype=float) - np.asarray(
site.geometry.centroid.coords[0], dtype=float
)
site_parts = _boundary_parts(site)
# Partition cells move separately so they are not merged.
for part in site_parts:
moved = affinity.translate(
part.geometry, xoff=float(offset[0]), yoff=float(offset[1])
)
parts.append(_GeometryBoundary(moved, positioned=True))
groups.append(group)
records.append(
{
"site": index,
"group": group,
"boundary": type(site).__name__,
"seed": int(site_seed),
"clearance": clearance,
"center": [float(value) for value in arranged.centroid.coords[0]],
"area": float(arranged.area),
"n_parts": len(site_parts),
}
)
if not parts:
smallest = min(tried) / 1e6
raise ValueError(
f"no site fits inside the domain of {developable.area / 1e6:.4g} km2, "
f"where {len(tried)} sites of {smallest:.4g} km2 or more found no "
"room, so use smaller sites or a larger domain"
)
return tuple(parts), tuple(groups), records
def _boundary_parts(boundary: Boundary) -> tuple:
"""One boundary per partition cell or polygon of ``boundary``, else itself."""
parts = boundary.parts
if len(parts) != 1 or parts[0] is not boundary:
return tuple(parts)
geometry = boundary.geometry
if isinstance(geometry, shapely.MultiPolygon):
return tuple(_GeometryBoundary(polygon) for polygon in geometry.geoms)
return (boundary,)
def _manifest_entry(farm: Farm, summary: dict, filename: str) -> dict:
"""The manifest record of one saved farm, with any partition cells."""
entry = {
"file": filename,
"summary": summary,
"metadata": _json_safe(farm.metadata),
}
cells = None if farm.boundary is None else _cells(farm.boundary)
if cells is not None:
# Hex WKB keeps every coordinate bit, so the cells tile the site again.
entry["cells"] = [shapely.to_wkb(cell, hex=True) for cell in cells]
return entry
def _fixed_too_close(first: int, second: int, geometries, gap: float) -> ValueError:
"""The error for two members fixed in place that overlap or break the gap."""
distance = float(geometries[first].distance(geometries[second]))
where = (
"overlap"
if distance == 0.0
else f"are {distance:.6g} m apart, closer than the gap of {gap:.6g} m"
)
return ValueError(
f"cluster members {first} and {second} are fixed in place and {where}. "
"A partition's cells and a boundary given a center keep their place, so "
"space them apart, or set the gap within the partition instead"
)
def _arrange_boundaries(boundaries, gap, seed):
"""Keep positioned boundaries and pack unpositioned ones around them.
Returns the arranged boundaries and JSON-compatible placement records.
"""
# Positioned boundaries must keep a fixed gap between themselves.
fixed_floor = float(gap) if is_number(gap) else 0.0
boundaries = tuple(boundaries)
positioned = [getattr(boundary, "_positioned", True) for boundary in boundaries]
geometries = [boundary.geometry for boundary in boundaries]
fixed_indices = [index for index, value in enumerate(positioned) if value]
unplaced_indices = [index for index, value in enumerate(positioned) if not value]
tolerance = 1e-10 * max(sum(geometry.area for geometry in geometries), 1.0)
for offset, first_index in enumerate(fixed_indices):
for second_index in fixed_indices[offset + 1 :]:
if not _has_clearance(
geometries[first_index],
[geometries[second_index]],
fixed_floor,
tolerance,
):
raise _fixed_too_close(
first_index, second_index, geometries, fixed_floor
)
placed = _PlacedGeometries(geometries[index] for index in fixed_indices)
disks = [
(
np.asarray(geometry.centroid.coords[0], dtype=float),
_geometry_radius(geometry),
)
for geometry in placed.geometries
]
translations = [np.zeros(2, dtype=float) for _ in boundaries]
clearances = [None for _ in boundaries]
angular_offsets = [None for _ in boundaries]
# Quantize equal-area ordering to avoid platform-dependent placement.
scale = max((geometry.area for geometry in geometries), default=1.0) or 1.0
order = sorted(
unplaced_indices,
key=lambda index: (-round(geometries[index].area / scale, 9), index),
)
for index in order:
geometry = geometries[index]
original_center = np.asarray(geometry.centroid.coords[0], dtype=float)
if not placed.geometries:
arranged = _translated_geometry(geometry, (0.0, 0.0))
else:
member_rng = as_rng(derive_seed(seed, "member", index))
clearance = draw_scalar(
gap, member_rng, "gap", index=index, choice_seed=seed, minimum=0.0
)
angular_offset = float(member_rng.uniform(0.0, np.pi / 3.0))
arranged = _compact_position(
geometry,
placed,
disks,
clearance,
angular_offset,
)
clearances[index] = clearance
angular_offsets[index] = float(np.rad2deg(angular_offset))
resolved_center = np.asarray(arranged.centroid.coords[0], dtype=float)
translations[index] = resolved_center - original_center
geometries[index] = arranged
placed.append(arranged)
disks.append((resolved_center, _geometry_radius(arranged)))
arranged_boundaries = tuple(
boundary
if positioned[index]
else _GeometryBoundary(
geometries[index],
positioned=True,
# Partition cells move with the shape, so a farm keeps its parts.
cells=_cells(boundary, _translation(translations[index])),
)
for index, boundary in enumerate(boundaries)
)
records = [
{
"index": index,
"source": type(boundary).__name__,
"positioned": bool(positioned[index]),
"translation": translations[index].astype(float).tolist(),
"center": list(map(float, geometries[index].centroid.coords[0])),
"clearance": clearances[index],
"angular_offset": angular_offsets[index],
}
for index, boundary in enumerate(boundaries)
]
return arranged_boundaries, records
def _member_metadata(metadata: dict, kept: list[int], farms) -> dict:
"""Cut the per-member records of generated cluster ``metadata`` to ``kept``.
Member records follow ``kept`` in order. Site records of a scattered
cluster keep the sites that still hold a kept farm, and the totals derived
from the records are recomputed. Cluster-level keys stay as they are.
"""
if metadata.get("kind") != "generated_cluster":
return metadata
groups = {farm.metadata.get("group") for farm in farms}
def subset(records: list) -> list:
# Site records carry the group that the member farms record.
if records and all(isinstance(r, Mapping) and "group" in r for r in records):
return [record for record in records if record["group"] in groups]
return [records[index] for index in kept]
result = dict(metadata)
if isinstance(result.get("counts"), list):
result["counts"] = subset(result["counts"])
counts = result["counts"]
result["requested_n"] = None if None in counts else sum(counts)
if "achieved_n" in result:
result["achieved_n"] = sum(farm.n_turbines for farm in farms)
if isinstance(result.get("boundaries"), list):
result["boundaries"] = subset(result["boundaries"])
placement = result.get("placement")
if isinstance(placement, Mapping) and isinstance(placement.get("members"), list):
placement = dict(placement)
placement["members"] = subset(placement["members"])
if "coverage" in placement:
area = sum(record["area"] for record in placement["members"])
placement["coverage"] = area / placement["domain_area"]
if "n_sites" in placement:
placement["n_sites"] = len(placement["members"])
result["placement"] = placement
turbines = result.get("turbines")
if isinstance(turbines, Mapping) and isinstance(turbines.get("members"), list):
result["turbines"] = {**turbines, "members": subset(turbines["members"])}
return result
[docs]
class Cluster(Sequence[Farm]):
"""A sequence of :class:`~oowl.Farm` members with metadata for the whole set.
Cluster metadata holds the root seed, the placement record and any shared
wind resource and surface. A slice is a cluster of the kept members whose
per-member records are cut to match.
"""
def __init__(self, farms: Sequence[Farm], *, metadata: Mapping | None = None):
self.farms = tuple(farms)
if not all(isinstance(farm, Farm) for farm in self.farms):
raise TypeError("farms must contain only Farm objects")
self.metadata = _json_safe(dict(metadata or {}))
def __len__(self) -> int:
return len(self.farms)
def __repr__(self) -> str:
return f"Cluster(n_farms={self.n_farms}, n_turbines={self.n_turbines})"
def __eq__(self, other: object) -> bool:
"""Whether ``other`` is a cluster of equal farms and the same metadata."""
if other.__class__ is not self.__class__:
return NotImplemented
return self.farms == other.farms and self.metadata == other.metadata
@overload
def __getitem__(self, index: int) -> Farm: ...
@overload
def __getitem__(self, index: slice) -> Cluster: ...
def __getitem__(self, index: int | slice) -> Farm | Cluster:
if isinstance(index, slice):
kept = list(range(len(self.farms)))[index]
farms = [self.farms[position] for position in kept]
return Cluster(farms, metadata=_member_metadata(self.metadata, kept, farms))
return self.farms[index]
def __iter__(self) -> Iterator[Farm]:
return iter(self.farms)
[docs]
def summary(self, metrics: Sequence[str] | None = None) -> list[dict]:
"""Return one summary row per farm, with its sampling record prefixed.
Each row holds ``index``, ``algorithm`` and ``seed`` from the farm's
``sample`` metadata, which are ``None`` for a farm not drawn by
:meth:`~oowl.Sampler.sample`, plus ``Farm.summary(metrics)``.
"""
rows = []
for farm in self.farms:
record = farm.metadata.get("sample", {})
row = {
"index": record.get("index"),
"algorithm": record.get("layout", {}).get("algorithm"),
"seed": record.get("seed"),
}
row.update(farm.summary(metrics))
rows.append(row)
return rows
[docs]
def save(self, directory: str | Path, **kwargs: Any) -> Path:
"""Write the cluster into ``directory`` and return that directory.
Writes ``manifest.json`` holding ``schema_version``, the cluster
metadata and a ``{file, summary, metadata}`` entry per farm. A farm on
partition cells adds their hex WKB as ``cells``, so its parts reload
in the same order. Each farm is saved as windIO YAML
(``farm_%04d.yaml``), a system with its site where it has a boundary
and a ``wind_farm`` document of its positions and turbines where it
has none, as :meth:`Farm.save` writes them. Extra keyword arguments go
to ``Farm.save``, only ``name`` and ``crs`` to a farm without a
boundary. A resource or surface that several farms share is written
once, as ``resource_%04d.yaml`` or ``surface_%04d.yaml``, which their
files ``!include``. windIO records no guide, and a farm without a
boundary has no site to hold a resource or surface, so
:meth:`Cluster.load` does not recover those.
"""
from ._windio import _save, _SharedComponents
directory = Path(directory)
directory.mkdir(parents=True, exist_ok=True)
shared = _SharedComponents(
directory,
[farm for farm in self.farms if farm.boundary is not None],
kwargs.get("energy_resource"),
kwargs.get("surface_format", "yaml"),
)
layout_only = {key: kwargs[key] for key in ("name", "crs") if key in kwargs}
files = []
for index, farm in enumerate(self.farms):
filename = f"farm_{index:04d}.yaml"
if farm.boundary is None:
_save(farm, directory / filename, document="wind_farm", **layout_only)
else:
_save(farm, directory / filename, shared=shared, **kwargs)
files.append(filename)
manifest = {
"schema_version": 1,
"collection": self.metadata,
"farms": [
_manifest_entry(farm, summary, filename)
for farm, summary, filename in zip(self.farms, self.summary(), files)
],
}
(directory / "manifest.json").write_text(
json.dumps(_json_safe(manifest), indent=2), encoding="utf-8"
)
return directory
[docs]
@classmethod
def load(cls, directory: str | Path) -> "Cluster":
"""Rebuild a cluster from a directory written by :meth:`Cluster.save`.
Loads each farm through windIO, then restores its metadata and
partition cells from the manifest. A farm keeps the name and CRS its
file records, and farms that include one file share the object it
holds. Raises :class:`FileNotFoundError` for a missing manifest and
:class:`ValueError` for an unsupported schema or a malformed entry.
"""
directory = Path(directory)
manifest_path = directory / "manifest.json"
if not manifest_path.is_file():
raise FileNotFoundError(f"collection manifest not found: {manifest_path}")
manifest = json.loads(manifest_path.read_text(encoding="utf-8"))
if manifest.get("schema_version") != 1:
raise ValueError("unsupported Cluster schema")
from ._windio import _Documents, _farm_name, _load
documents = _Documents()
entries = manifest.get("farms", [])
farms = []
# The metadata each windIO file loads with.
recorded = []
for entry in entries:
filename = entry.get("file") if isinstance(entry, Mapping) else None
if not isinstance(filename, str):
raise ValueError("collection manifest entry has no file")
path = directory / filename
member = _load(path, documents)
if not isinstance(member, Farm):
raise ValueError(f"collection member is not a farm: {path}")
farms.append(member)
recorded.append(member.metadata)
farms = tuple(farms)
if len(farms) != len(entries):
raise ValueError("collection manifest and farm files disagree")
for farm, entry, written in zip(farms, entries, recorded):
if isinstance(entry, Mapping) and isinstance(
entry.get("metadata"), Mapping
):
metadata = dict(entry["metadata"])
# A name or CRS passed to save is only in the file, so it is
# kept for the next save to write again.
if written["name"] != _farm_name(metadata):
metadata["name"] = written["name"]
if written.get("crs") is not None:
metadata["crs"] = written["crs"]
farm.metadata = metadata
if isinstance(entry, Mapping) and entry.get("cells") is not None:
if farm.boundary is None:
raise ValueError("collection manifest gives cells without a site")
farm.boundary = _GeometryBoundary(
farm.boundary.geometry,
cells=tuple(shapely.from_wkb(cell) for cell in entry["cells"]),
)
return cls(farms, metadata=manifest.get("collection", {}))
@property
def n_farms(self) -> int:
"""Number of member farms."""
return len(self.farms)
@property
def n_turbines(self) -> int:
"""Total number of turbines across every farm."""
return sum(farm.n_turbines for farm in self.farms)
@property
def positions(self) -> np.ndarray:
"""Every farm's positions concatenated, in order."""
arrays = [farm.positions for farm in self.farms if farm.n_turbines]
return np.concatenate(arrays) if arrays else np.empty((0, 2), dtype=float)
@property
def boundaries(self) -> tuple[Boundary | None, ...]:
"""Boundary of each farm, in order."""
return tuple(farm.boundary for farm in self.farms)
@property
def bounds(self) -> tuple[float, float, float, float]:
"""Turbine-position bounding box ``(min_x, min_y, max_x, max_y)`` in meters."""
positions = self.positions
if not len(positions):
raise ValueError("an empty collection has no turbine bounds")
return (
float(np.min(positions[:, 0])),
float(np.min(positions[:, 1])),
float(np.max(positions[:, 0])),
float(np.max(positions[:, 1])),
)
@property
def centroid(self) -> tuple[float, float]:
"""Mean turbine position ``(x, y)`` across every farm, in meters."""
positions = self.positions
if not len(positions):
raise ValueError("an empty collection has no turbine centroid")
center = np.mean(positions, axis=0)
return (float(center[0]), float(center[1]))
@property
def extent(self) -> tuple[float, float]:
"""Turbine-position extent ``(width, height)`` in meters."""
min_x, min_y, max_x, max_y = self.bounds
return (max_x - min_x, max_y - min_y)
[docs]
def plot(
self, ax: "matplotlib.axes.Axes | None" = None, **kwargs: Any
) -> "matplotlib.axes.Axes":
"""Plot every member farm on one map and return its Matplotlib axes.
Each member is drawn in its own color of the Matplotlib cycle, in
member order. For a cluster scattered across a ``domain``,
``regions`` outlines the domain and shades its exclusion zones, and
``regions=False`` omits them. Other keywords go to
:meth:`~oowl.Farm.plot`, so turbine size is ``s``. Raises
:class:`ValueError` when the cluster holds no farms.
"""
from .plotting import _plot_cluster
if not self.farms:
raise ValueError("no farms to plot")
return _plot_cluster(self, ax=ax, **kwargs)
