"""windIO serialization for farms, turbines and wind resources.
:func:`save` and :func:`load` read and write the windIO subset that OOWL
produces. The YAML reader and writer format numeric lists directly, because
PyYAML's per-scalar path dominates the cost of large layouts.
"""
from __future__ import annotations
import json
import math
import re
import typing
import uuid
from collections import Counter
from collections.abc import Mapping
from functools import partial
from pathlib import Path
from typing import TYPE_CHECKING
import numpy as np
from shapely.errors import GEOSException
from shapely.geometry import MultiPolygon, Point, Polygon
if TYPE_CHECKING: # pragma: no cover
from .layouts.base import Farm
from .resources import Dimitrov, WindResource, _Resource
from .turbines import Turbine
__all__ = ["save", "load"]
def _require_json_compatible(document: dict) -> None:
"""Reject values on which JSON and YAML would disagree.
Leaves must be plain ``str``, ``int``, ``float``, ``bool`` or ``None``, or
an :class:`_Include`, which only YAML can write. Safe YAML cannot write
NumPy scalars or arrays, which JSON partly accepts.
"""
stack = [document]
while stack:
node = stack.pop()
if isinstance(node, dict):
for key, value in node.items():
if not isinstance(key, str):
raise TypeError(
f"windIO document keys must be strings; got {key!r} "
f"of type {type(key).__name__}"
)
stack.append(value)
elif isinstance(node, (list, tuple)):
# Whole float lists, the bulk of a document, are checked in C.
if set(map(type, node)) == {float}:
if not all(map(math.isfinite, node)):
raise ValueError(
"windIO documents cannot contain non-finite numbers"
)
else:
stack.extend(node)
elif type(node) is float:
if not math.isfinite(node):
raise ValueError("windIO documents cannot contain non-finite numbers")
elif type(node) not in (str, int, bool, type(None), _Include):
raise TypeError(
"windIO document values must be str, int, float, bool or None; "
f"got {type(node).__name__}, so convert NumPy values with tolist()"
)
# Decimal scalars that PyYAML's YAML 1.1 resolver reads as exactly int() or
# float() of their text. Anything else in a flow sequence is left to PyYAML.
_YAML_FLOAT = r"[-+]?[0-9]+\.[0-9]*(?:[eE][-+][0-9]+)?"
_YAML_NUMBER = rf"(?>{_YAML_FLOAT}|[-+]?(?:0|[1-9][0-9]*))"
# YAML whitespace only. Python's \s also matches Unicode spaces YAML keeps.
_SPACE = "[ \\t\\r\\n]"
# The bracket comes first so the search skips to each "[", and possessive
# repeats never backtrack, since no shorter number can be followed by a comma.
_NUMERIC_FLOW = re.compile(
rf"\[(?<=: \[|- \[){_SPACE}*+{_YAML_NUMBER}"
rf"(?:{_SPACE}*+,{_SPACE}*+{_YAML_NUMBER})*+{_SPACE}*+\]"
)
_LINE_WIDTH = 88
# Items converted per batch before the conversion cache is judged.
_CACHE_BATCH = 4096
# A float written by repr as an exponent whose mantissa has no dot, as
# 1e-05, at the start of a flow sequence item.
_UNDOTTED_EXPONENT = re.compile(r"(?<=[\[ ])(-?[0-9]+)e")
# Numeric values beyond which a list is written on one line. Only a surface is
# that long, and its values are not read line by line.
_WRAPPED_MOST = 10_000
# Characters a flow sequence between its brackets reaches beyond which its
# floats are read without ``_NUMERIC_FLOW``, see ``_long_floats``.
_LONG_FLOW = 100_000
# Deletes the characters of decimal numbers and the space between them.
_NUMBER_CHARACTERS = str.maketrans("", "", "0123456789.eE+-, \t\r\n")
# An item's dot that no digit leads, which YAML 1.1 reads as a string.
_BARE_DOTS = tuple(f"{lead}." for lead in " ,-+\t\r\n")
# An exponent without a sign, which YAML 1.1 reads as a string.
_UNSIGNED_EXPONENTS = tuple(f"{e}{digit}" for e in "eE" for digit in "0123456789")
class _Include:
"""A YAML ``!include`` of ``file``, relative to the including document."""
__slots__ = ("file",)
def __init__(self, file: str):
self.file = file
def _dotted_mantissa(text: str) -> str:
"""Give an exponent's mantissa a dot, as YAML 1.1 reads 1e-05 as a string."""
if "e" in text and "." not in text:
return text.replace("e", ".0e", 1)
return text
def _yaml_number(value: int | float) -> str:
if type(value) is int:
return str(value)
return _dotted_mantissa(repr(value))
class _Converted(dict):
"""``convert`` of each item, so that a repeated item converts once.
A zero float is not stored, because 0.0 and -0.0 are one key with two texts.
"""
def __init__(self, convert):
super().__init__()
self._convert = convert
def __missing__(self, key):
value = self._convert(key)
if key:
self[key] = value
return value
def _convert_all(convert, items: list | tuple) -> list:
"""Return ``list(map(convert, items))``, cached while most items repeat.
The point coordinates of a surface repeat each grid line's x or y, so each
distinct value converts once. Items that are mostly new skip the cache.
"""
cache = _Converted(convert)
values: list = []
for start in range(0, len(items), _CACHE_BATCH):
batch = items[start : start + _CACHE_BATCH]
known = len(cache)
values.extend(map(cache.__getitem__, batch))
if 2 * (len(cache) - known) > len(batch):
values.extend(map(convert, items[start + _CACHE_BATCH :]))
break
return values
def _number_sequence(values: list | tuple, kinds: set, column: int) -> str:
"""``values`` of the types ``kinds`` as a flow sequence starting at ``column``.
YAML 1.1 reads it back exactly, a float written as ``repr``. A list of more
than ``_WRAPPED_MOST`` floats, a surface's, is written on one line, by the
list's own ``repr`` where its values are mostly distinct, as depths, and
from cached strings where they repeat, as grid coordinates.
"""
if kinds == {float} and len(values) > _WRAPPED_MOST:
sample = values[:_CACHE_BATCH]
if 2 * len(set(sample)) > len(sample):
text = repr(values if isinstance(values, list) else list(values))
else:
text = "[" + ", ".join(_convert_all(float.__repr__, values)) + "]"
if "e" in text:
# YAML 1.1 reads 1e-05 as a string, so its mantissa gets a dot.
text = _UNDOTTED_EXPONENT.sub(r"\1.0e", text)
return text
if kinds == {int}:
strings = list(map(str, values))
elif kinds == {float}:
strings = _convert_all(float.__repr__, values)
else:
strings = list(map(_yaml_number, values))
return _flow_sequence(strings, column)
def _flow_sequence(strings: list[str], column: int) -> str:
"""Join ``strings`` as a flow sequence whose values start at ``column``.
Each line holds as many values as the widest fits within ``_LINE_WIDTH``,
at least one, and continuation lines align under the first value, deeper
than any key. A fixed count per line places the breaks by slicing, so a
surface of millions of values is not walked line by line.
"""
text = _joined_lines(strings, column)
if "e" in text and _UNDOTTED_EXPONENT.search(text):
# YAML 1.1 reads 1e-05 as a string, so its mantissa gets a dot.
text = _joined_lines(list(map(_dotted_mantissa, strings)), column)
return text
def _joined_lines(strings: list[str], column: int) -> str:
"""``strings`` in brackets, a fixed number per line, see :func:`_flow_sequence`."""
widest = max(map(len, strings)) + 2
per_line = max(1, (_LINE_WIDTH - column) // widest)
pieces = [", "] * (2 * len(strings) - 1)
pieces[::2] = strings
breaks = range(2 * per_line - 1, len(pieces), 2 * per_line)
pieces[2 * per_line - 1 :: 2 * per_line] = [",\n" + " " * column] * len(breaks)
return "[" + "".join(pieces) + "]"
def _dump_yaml(document: dict, yaml) -> str:
"""Return ``document`` as block YAML with numeric lists as flow sequences.
PyYAML represents every scalar in Python, which dominates the cost of large
arrays, so plain int and float lists are formatted here instead. An
:class:`_Include` is written as an ``!include`` tag.
"""
token = f"oowl-array-{uuid.uuid4().hex}-"
arrays: list = []
includes: list[str] = []
def swap(node):
if isinstance(node, dict):
return {key: swap(value) for key, value in node.items()}
if isinstance(node, _Include):
includes.append(node.file)
return f"{token}i{len(includes) - 1}"
if isinstance(node, (list, tuple)):
kinds = set(map(type, node))
if node and kinds <= {int, float}:
arrays.append((node, kinds))
return f"{token}{len(arrays) - 1}"
return [swap(item) for item in node]
return node
# libyaml emits the same document about four times faster.
dumper = getattr(yaml, "CSafeDumper", yaml.SafeDumper)
text = yaml.dump(swap(document), Dumper=dumper, sort_keys=False)
def render(match) -> str:
if match.group(1):
return f"!include {includes[int(match.group(2))]}"
column = match.start() - text.rfind("\n", 0, match.start())
values, kinds = arrays[int(match.group(2))]
return _number_sequence(values, kinds, column)
return re.sub(re.escape(token) + r"(i?)([0-9]+)", render, text)
def _parse_flow(sequence: str) -> list:
"""Read a flow sequence that matched ``_NUMERIC_FLOW`` as PyYAML does.
An item is a float exactly when it has a dot, and int() and float() ignore
the YAML spaces around it.
"""
items = sequence[1:-1].split(",")
dots = sequence.count(".")
if dots == len(items):
return _convert_all(float, items)
if dots == 0:
return list(map(int, items))
return [float(item) if "." in item else int(item) for item in items]
def _long_floats(inner: str) -> list | None:
"""The floats between a flow sequence's brackets, or ``None``.
``None`` unless PyYAML surely reads every item as a float of the same
value. Only number characters may appear and each item has its dot. A
dot no digit leads or an unsigned exponent makes YAML 1.1 read a string,
which plain substring searches find. ``float`` then reads every item as
PyYAML does, and raises for a malformed one.
"""
if inner.translate(_NUMBER_CHARACTERS):
return None
items = inner.split(",")
if inner.count(".") != len(items):
return None
if inner.lstrip().startswith(".") or any(pair in inner for pair in _BARE_DOTS):
return None
if ("e" in inner or "E" in inner) and any(
pair in inner for pair in _UNSIGNED_EXPONENTS
):
return None
try:
return _convert_all(float, items)
except ValueError:
return None
def _swap_long_floats(text: str, token: str, arrays: list) -> str:
"""``text`` with each long float flow sequence swapped for a placeholder.
A surface's values fill most of a large document, and the general
pattern reads them about ten times slower. Each sequence swapped adds its
values to ``arrays``. A sequence nested in another, or in a string, is
left to :data:`_NUMERIC_FLOW`, or found in a string later and read again.
"""
pieces = []
position = 0
start = text.find("[")
while start != -1:
end = text.find("]", start)
if end == -1:
break
if end - start > _LONG_FLOW and text[start - 2 : start] in (": ", "- "):
values = _long_floats(text[start + 1 : end])
if values is not None:
pieces.append(text[position:start])
pieces.append(f"{token}{len(arrays)}")
arrays.append(values)
position = end + 1
start = text.find("[", end + 1)
pieces.append(text[position:])
return "".join(pieces)
def _load_yaml_text(text: str, loader, yaml, included: Mapping[int, object]):
"""Load YAML, reading numeric flow sequences without PyYAML's scalar path.
The result equals a plain ``yaml.load``. Each matched sequence is swapped for
a unique placeholder, and any placeholder that does not come back as a whole
value, such as one inside a string or a key, falls back to a plain load.
Documents in ``included``, keyed by their ``id``, were read from other
files and hold no placeholder, so they are not searched.
"""
token = f"oowl-array-{uuid.uuid4().hex}-"
arrays: list[list] = []
def swap(match) -> str:
arrays.append(_parse_flow(match.group(0)))
return f"{token}{len(arrays) - 1}"
swapped = _NUMERIC_FLOW.sub(swap, _swap_long_floats(text, token, arrays))
if not arrays:
return yaml.load(text, Loader=loader)
document = yaml.load(swapped, Loader=loader)
restored = 0
seen: set[int] = set()
stack = [document]
while stack:
node = stack.pop()
if not isinstance(node, (dict, list)) or id(node) in seen:
continue
seen.add(id(node))
if included.get(id(node)) is node:
continue
if isinstance(node, dict):
if any(isinstance(key, str) and token in key for key in node):
return yaml.load(text, Loader=loader)
slots = node.items()
else:
slots = enumerate(node)
for slot, value in list(slots):
if isinstance(value, str) and token in value:
index = value[len(token) :]
if not (value.startswith(token) and index.isdigit()):
return yaml.load(text, Loader=loader)
node[slot] = arrays[int(index)]
restored += 1
else:
stack.append(value)
if restored != len(arrays):
return yaml.load(text, Loader=loader)
return document
def _write_document(document: dict, path: Path) -> None:
"""Write a finite, JSON-compatible document as safe YAML.
Parent directories are created only after the document validates.
"""
_require_json_compatible(document)
path.parent.mkdir(parents=True, exist_ok=True)
try:
import yaml
except ImportError:
path.write_text(json.dumps(document, indent=2, allow_nan=False) + "\n")
else:
path.write_text(_dump_yaml(document, yaml))
class _Documents:
"""Included documents read once, and the objects built from them.
Every file that includes the same document gets the same parsed mapping,
and :meth:`build` makes one object from it, shared by all of them.
"""
def __init__(self):
self._paths: dict[Path, object] = {}
# Included mappings and lists by id, kept alive so ids stay unique.
self.included: dict[int, object] = {}
self._built: dict[tuple[int, str], object] = {}
def read(self, resolved: Path, load):
"""Return the document at ``resolved``, calling ``load`` the first time."""
if resolved not in self._paths:
document = load()
self._paths[resolved] = document
if isinstance(document, (dict, list)):
self.included[id(document)] = document
return self._paths[resolved]
def build(self, entry, kind: str, make):
"""Return ``make(entry)``, made once for an included ``entry``."""
if self.included.get(id(entry)) is not entry:
return make(entry)
key = (id(entry), kind)
if key not in self._built:
self._built[key] = make(entry)
return self._built[key]
def _read_document(path: Path, documents: _Documents | None = None) -> dict:
"""Read ``path``, sharing included documents through ``documents``."""
try:
import yaml # noqa: F401
except ImportError:
return json.loads(path.read_text())
if documents is None:
documents = _Documents()
return _load_yaml(path, frozenset(), documents)
def _load_yaml(path: Path, including: frozenset[Path], documents: _Documents):
"""Safely load ``path``, resolving ``!include`` relative to the including file.
YAML includes and NetCDF surfaces are read, each file once per
``documents``. Other includes raise ``ValueError``.
"""
import yaml
class Loader(getattr(yaml, "CSafeLoader", yaml.SafeLoader)):
pass
def include(loader, node):
target = path.parent / loader.construct_scalar(node)
suffix = target.suffix.lower()
resolved = target.resolve()
if suffix == ".nc":
return documents.read(resolved, lambda: _read_surface_netcdf(target))
if suffix not in (".yaml", ".yml"):
raise ValueError(
f"{path}: cannot read !include {target.name}; only YAML and "
"NetCDF surface includes are supported"
)
if resolved in including:
raise ValueError(f"{path}: !include {target.name} includes itself")
return documents.read(
resolved, lambda: _load_yaml(target, including | {resolved}, documents)
)
Loader.add_constructor("!include", include)
return _load_yaml_text(path.read_text(), Loader, yaml, documents.included)
_DOCUMENT_KINDS = (
"wind_energy_system",
"wind_farm",
"site",
"turbine",
"energy_resource",
)
def _number(value, key: str) -> float:
"""Return a windIO number as a float, naming ``key`` if it is not one."""
try:
return float(value)
except (TypeError, ValueError):
raise ValueError(f"windIO {key} must be a number") from None
def _number_list(value, key: str) -> np.ndarray:
"""Return a windIO list of numbers as a 1-D float array, naming ``key``."""
message = f"windIO {key} must be a list of numbers"
try:
array = np.asarray(value, dtype=float)
except (TypeError, ValueError):
raise ValueError(message) from None
if array.ndim != 1:
raise ValueError(message)
return array
def _validate_layouts(wind_farm: dict) -> None:
layouts = wind_farm.get("layouts")
if not isinstance(layouts, (dict, list)):
raise ValueError("windIO wind_farm is missing layouts")
layout = _first_layout(layouts)
if not isinstance(layout, dict) or not isinstance(layout.get("coordinates"), dict):
raise ValueError("windIO layout is missing its coordinates mapping")
coordinates = layout["coordinates"]
if "x" not in coordinates or "y" not in coordinates:
raise ValueError("windIO layout coordinates must include x and y")
x = _number_list(coordinates["x"], "layout coordinates x")
y = _number_list(coordinates["y"], "layout coordinates y")
if len(x) != len(y):
raise ValueError("windIO layout coordinate x/y lengths differ")
if "z" in coordinates and len(
_number_list(coordinates["z"], "layout coordinates z")
) != len(x):
raise ValueError("windIO layout coordinate z length differs")
layout_types = layout.get("turbine_types")
if layout_types is not None and not isinstance(layout_types, list):
raise ValueError("windIO layout turbine_types must be a list")
if layout_types is not None and len(layout_types) != len(x):
raise ValueError("windIO layout turbine_types length must match coordinates")
def _validate_site_mapping(site: dict) -> None:
"""Check the site keys OOWL reads. A missing ``energy_resource`` is allowed."""
boundaries = site.get("boundaries")
if not isinstance(boundaries, dict) or not (
"polygons" in boundaries or "circle" in boundaries
):
raise ValueError("windIO site is missing boundary polygons or a circle")
if "exclusions" in site:
exclusions = site["exclusions"]
if not isinstance(exclusions, dict):
raise ValueError("windIO site exclusions must be a mapping")
if not isinstance(exclusions.get("polygons") or [], list):
raise ValueError("windIO exclusion polygons must be a list")
def _validate_oowl_windio_document(document: dict, kind: str) -> None:
"""Validate the small windIO subset OOWL reads and writes."""
if "name" not in document:
raise ValueError("windIO document is missing 'name'")
if kind == "turbine":
for key in ("performance", "hub_height", "rotor_diameter"):
if key not in document:
raise ValueError(f"windIO turbine document is missing {key!r}")
elif kind == "wind_farm":
_validate_layouts(document)
elif kind == "site":
_validate_site_mapping(document)
elif kind == "energy_resource":
if not isinstance(document["wind_resource"], dict):
raise ValueError("windIO 'wind_resource' must be a mapping")
elif kind == "wind_energy_system":
for key in ("site", "wind_farm"):
if not isinstance(document[key], dict):
raise ValueError(f"windIO {key!r} must be a mapping")
_validate_site_mapping(document["site"])
_validate_layouts(document["wind_farm"])
def _document_kind(document: dict) -> str:
"""Name the windIO document kind from the keys it carries."""
if "site" in document and "wind_farm" in document:
return "wind_energy_system"
if "layouts" in document:
return "wind_farm"
if "boundaries" in document:
return "site"
if "performance" in document:
return "turbine"
if "wind_resource" in document:
return "energy_resource"
raise ValueError(
"windIO document is none of a wind energy system, wind farm, site, "
"turbine or energy resource"
)
def _placeholder_resource() -> dict:
"""Return a minimal wind resource marked as a placeholder."""
return {
"name": "PLACEHOLDER - replace with a real wind resource",
"wind_resource": {"probability": {"data": 1.0, "dims": []}},
}
def _ring_coordinates(ring, crs: str | None) -> dict:
coords = np.asarray(ring.coords, dtype=float)
if len(coords) > 1 and np.allclose(coords[0], coords[-1]):
coords = coords[:-1]
entry = {"x": coords[:, 0].tolist(), "y": coords[:, 1].tolist()}
if crs:
entry["crs"] = crs
return entry
class _SharedComponents:
"""Writes each resource or surface that several farms share to one file.
The farm documents saved to ``directory`` then ``!include`` that file, so
the component is formatted once and read once. Without PyYAML documents
are written as JSON, which cannot include, so every farm keeps a copy.
"""
def __init__(
self,
directory: Path,
farms,
energy_resource: dict | None = None,
surface_format: str = "yaml",
):
try:
import yaml # noqa: F401
except ImportError:
farms = ()
counts = Counter(
id(component)
for farm in farms
for component in (
farm.resource if energy_resource is None else energy_resource,
farm.surface,
)
if component is not None
)
self._directory = directory
self._shared = {key for key, count in counts.items() if count > 1}
self._includes: dict[tuple[int, str | None], _Include] = {}
self._files: Counter[str] = Counter()
self._surface_format = surface_format
def _include(self, component, crs: str | None, prefix: str, write, suffix):
"""The include of ``component``, calling ``write(path)`` the first time."""
if id(component) not in self._shared:
return None
key = (id(component), crs)
if key not in self._includes:
filename = f"{prefix}_{self._files[prefix]:04d}{suffix}"
write(self._directory / filename)
self._files[prefix] += 1
self._includes[key] = _Include(filename)
return self._includes[key]
def energy_resource(self, component) -> _Include | None:
"""The include of a shared resource or ``energy_resource`` mapping."""
def write(path: Path) -> None:
document = (
component
if isinstance(component, dict)
else component._to_windio_energy_resource()
)
_write_document(document, path)
return self._include(component, None, "resource", write, ".yaml")
def surface(self, surface, crs: str | None) -> _Include | None:
"""The include of a shared surface with its coordinates in ``crs``."""
if self._surface_format == "netcdf":
return self._include(
surface,
crs,
"surface",
lambda path: _write_surface_netcdf(surface, path, crs),
".nc",
)
return self._include(
surface,
crs,
"surface",
lambda path: _write_document(_surface_mapping(surface, crs), path),
".yaml",
)
def _resolve_energy_resource(
farm, name: str, energy_resource, shared: _SharedComponents | None
):
component = energy_resource
if component is None:
component = getattr(farm, "resource", None)
if component is None:
return _placeholder_resource()
include = None if shared is None else shared.energy_resource(component)
if include is not None:
return include
if component is energy_resource:
return energy_resource
return component._to_windio_energy_resource(name=f"{name} wind resource")
def _surface_mapping(surface, crs: str | None) -> dict:
"""Map a regular surface to windIO point form, ``x`` varying fastest."""
if surface.kind not in ("bathymetry", "terrain"):
raise TypeError(f"unsupported surface kind: {surface.kind!r}")
x, y = np.meshgrid(
np.asarray(surface.x, dtype=float), np.asarray(surface.y, dtype=float)
)
coordinates = {"x": x.ravel().tolist(), "y": y.ravel().tolist()}
if crs:
coordinates["crs"] = crs
field = "depth" if surface.kind == "bathymetry" else "elevation"
values = np.asarray(surface.values, dtype=float).ravel().tolist()
return {"coordinates": coordinates, field: values}
def _write_surface_netcdf(
surface, path: Path, crs: str | None, metadata: str | None = None
) -> None:
"""Write ``surface`` to ``path`` as a NetCDF grid, its field indexed ``(y, x)``.
The file is NetCDF 3 with 64-bit offsets, which SciPy writes and xarray
and windIO read, and holds ``x``, ``y`` and ``depth`` or ``elevation`` in
meters, with the CRS as a global attribute when known, and ``metadata``,
JSON text, as the global attribute ``oowl_metadata`` when given.
"""
from scipy.io import netcdf_file
if surface.kind not in ("bathymetry", "terrain"):
raise TypeError(f"unsupported surface kind: {surface.kind!r}")
field = "depth" if surface.kind == "bathymetry" else "elevation"
path.parent.mkdir(parents=True, exist_ok=True)
with netcdf_file(path, "w", version=2) as handle:
if crs:
handle.crs = crs
if metadata is not None:
handle.oowl_metadata = metadata
for name, values in (("x", surface.x), ("y", surface.y)):
handle.createDimension(name, len(values))
variable = handle.createVariable(name, "d", (name,))
variable[:] = np.asarray(values, dtype=float)
variable.units = "m"
variable = handle.createVariable(field, "d", ("y", "x"))
variable[:] = np.asarray(surface.values, dtype=float)
variable.units = "m"
def _read_surface_netcdf(path: Path) -> dict:
"""The surface in the NetCDF file at ``path``, as windIO grid form."""
return _read_surface_netcdf_file(path)[0]
def _read_surface_netcdf_file(path: Path) -> tuple[dict, str | None]:
"""The surface in the NetCDF file at ``path`` and its ``oowl_metadata``.
windIO reads the same file as ``x``, ``y`` and the field with its
``dims``, so the mapping is the one its loader makes, holding arrays, and
``crs`` when the file records one. The metadata is ``None`` when the file
records none.
"""
from scipy.io import netcdf_file
try:
handle = netcdf_file(path, "r", mmap=False)
except (OSError, TypeError, ValueError) as error:
raise ValueError(f"cannot read {path.name} as a NetCDF 3 file") from error
with handle:
fields = set(handle.variables) & {"depth", "elevation"}
if len(fields) != 1 or not {"x", "y"} <= set(handle.variables):
raise ValueError(
f"{path.name} must hold x, y and one of depth or elevation, as "
"OOWL reads NetCDF includes of surfaces only"
)
field = fields.pop()
variable = handle.variables[field]
entry: dict = {
"x": np.array(handle.variables["x"].data, dtype=float),
"y": np.array(handle.variables["y"].data, dtype=float),
field: {
"dims": [
name.decode() if isinstance(name, bytes) else name
for name in variable.dimensions
],
"data": np.array(variable.data, dtype=float),
},
}
crs = handle._attributes.get("crs")
metadata = handle._attributes.get("oowl_metadata")
if crs is not None:
entry["crs"] = crs.decode() if isinstance(crs, bytes) else str(crs)
if isinstance(metadata, bytes):
metadata = metadata.decode()
return entry, None if metadata is None else str(metadata)
def _build_site(
boundary,
name: str,
crs: str | None,
energy_resource,
surface,
shared: _SharedComponents | None,
own_surface: Path | None = None,
) -> dict:
"""The windIO site of ``boundary`` and its components.
A surface is written in point form, or included from the file ``shared``
holds for several farms, or else from the NetCDF file ``own_surface``.
"""
geometry = boundary.geometry
polygons = (
list(geometry.geoms) if isinstance(geometry, MultiPolygon) else [geometry]
)
exteriors = [_ring_coordinates(polygon.exterior, crs) for polygon in polygons]
holes = [
_ring_coordinates(hole, crs)
for polygon in polygons
for hole in polygon.interiors
]
site = {
"name": name,
"boundaries": {"polygons": exteriors},
"energy_resource": energy_resource,
}
if holes:
site["exclusions"] = {"polygons": holes}
if surface is not None:
key = "bathymetry" if surface.kind == "bathymetry" else "elevation"
include = None if shared is None else shared.surface(surface, crs)
if include is None and own_surface is not None:
_write_surface_netcdf(surface, own_surface, crs)
include = _Include(own_surface.name)
site[key] = _surface_mapping(surface, crs) if include is None else include
return site
def _build_wind_farm(farm, name: str, crs: str | None) -> dict:
"""Map a farm's positions and turbine types to a windIO ``wind_farm``."""
positions = np.asarray(farm.positions, dtype=float).reshape(-1, 2)
coordinates = {"x": positions[:, 0].tolist(), "y": positions[:, 1].tolist()}
surface = getattr(farm, "surface", None)
if surface is not None and surface.kind == "terrain":
coordinates["z"] = farm.coordinates[:, 2].tolist()
if crs:
coordinates["crs"] = crs
layout = {"coordinates": coordinates}
wind_farm = {"name": name, "layouts": layout}
assignment = getattr(farm, "turbines", None)
if assignment is not None:
from .turbines import _compress_turbines
types, type_indices = _compress_turbines(assignment, len(positions))
layout["turbine_types"] = type_indices.astype(int).tolist()
wind_farm["turbine_types"] = {
str(index): turbine._to_windio_turbine(
name=turbine.name or f"turbine type {index}"
)
for index, turbine in enumerate(types)
}
return wind_farm
def _resolved_kind(obj, document: str | None) -> str:
"""Name the document kind to write for *obj*, honoring an explicit choice."""
from .layouts.base import Farm
from .resources import _Resource
from .turbines import Turbine
for cls, label, kind in (
(Turbine, "turbine", "turbine"),
(_Resource, "resource", "energy_resource"),
):
if isinstance(obj, cls):
if document not in (None, kind):
raise ValueError(
f"a {label} is only written as a windIO {kind} document"
)
return kind
if not isinstance(obj, Farm):
raise TypeError(
"only farms, turbines and wind resources can be written to windIO"
)
if document is None:
# A windIO site, alone or in a system, requires its boundaries.
return "wind_farm" if obj.boundary is None else "wind_energy_system"
if document in ("turbine", "energy_resource"):
raise ValueError(f"a farm is not a windIO {document} document")
if document not in _DOCUMENT_KINDS:
raise ValueError(f"document must be one of {_DOCUMENT_KINDS}")
return document
def _save_document(document: dict, kind: str, path: Path) -> Path:
_validate_oowl_windio_document(document, kind)
_write_document(document, path)
return path
def _save_resource(resource: "_Resource", path: Path, name: str | None) -> Path:
"""Write a resource as a windIO ``energy_resource`` document."""
document = resource._to_windio_energy_resource(name=name or "oowl wind resource")
return _save_document(document, "energy_resource", path)
def _farm_name(metadata: Mapping) -> str:
"""The name a farm with ``metadata`` is saved under by default."""
name = metadata.get("name")
if name is not None and not isinstance(name, str):
raise TypeError("farm metadata['name'] must be a string or None")
return "oowl farm" if name is None else name
[docs]
def save(
obj: "Farm | Turbine | WindResource | Dimitrov",
path: str | Path,
*,
name: str | None = None,
crs: str | None = None,
energy_resource: dict | None = None,
document: typing.Literal[
"wind_energy_system", "wind_farm", "site", "turbine", "energy_resource"
]
| None = None,
surface_format: typing.Literal["yaml", "netcdf"] = "yaml",
) -> Path:
"""Write *obj* to *path* as a windIO document and return the path.
Parameters
----------
obj
A :class:`~oowl.layouts.base.Farm`, a :class:`~oowl.turbines.Turbine`,
a :class:`~oowl.resources.WindResource` or a
:class:`~oowl.resources.Dimitrov` design.
path
Destination ``.yaml`` file. Missing parent directories are created.
name
Name recorded in the document. Defaults to the farm's
``metadata["name"]``, as kept by :func:`load`, then to ``"oowl farm"``,
to the turbine's name for a turbine and to ``"oowl wind resource"`` for
a resource.
crs
Optional PROJ/EPSG string stored with every coordinate set. Defaults to
a CRS retained by :func:`load`. A turbine or resource document has no
coordinates, so passing ``crs`` for one raises :class:`ValueError`.
energy_resource
Optional windIO ``energy_resource`` mapping. Defaults to the farm's
attached resource. Only a site or system document carries one, so
passing it for a turbine, a resource or ``document="wind_farm"``
raises :class:`ValueError`.
document
Document kind to write. A turbine is always ``"turbine"`` and a
resource ``"energy_resource"``. A farm defaults to
``"wind_energy_system"``, and to ``"wind_farm"`` when it has no
boundary, since windIO sites need one. ``"wind_farm"`` writes the
layout alone and ``"site"`` the boundary, resource and surface alone.
A site or system for a farm without a resource gets a labeled uniform
placeholder.
surface_format
``"yaml"`` writes a surface into the document in windIO point form,
which validates against the windIO schema. ``"netcdf"`` writes it
beside the document as ``<name>_surface.nc``, a grid of ``x``, ``y``
and ``depth`` or ``elevation``, and includes it with ``!include``.
That file is about five times smaller and reads and writes in a
fraction of the time, and both :func:`load` and windIO load it, but
windIO's schema does not accept the included grid, so the site no
longer validates. It needs PyYAML, and a document without a site
raises :class:`ValueError` for it.
"""
return _save(
obj,
path,
name=name,
crs=crs,
energy_resource=energy_resource,
document=document,
surface_format=surface_format,
)
def _validated_surface_format(surface_format) -> str:
"""``surface_format`` checked to be ``"yaml"`` or ``"netcdf"``."""
if not isinstance(surface_format, str):
raise TypeError("surface_format must be 'yaml' or 'netcdf'")
if surface_format not in ("yaml", "netcdf"):
raise ValueError(
f"surface_format must be 'yaml' or 'netcdf', not {surface_format!r}"
)
return surface_format
def _save(
obj: "Farm | Turbine | WindResource | Dimitrov",
path: str | Path,
*,
shared: _SharedComponents | None = None,
name: str | None = None,
crs: str | None = None,
energy_resource: dict | None = None,
document: str | None = None,
surface_format: str = "yaml",
) -> Path:
""":func:`save`, writing components in ``shared`` once to include them."""
kind = _resolved_kind(obj, document)
path = Path(path)
from .layouts.base import Farm
from .turbines import Turbine
surface_format = _validated_surface_format(surface_format)
netcdf = surface_format if surface_format != "yaml" else None
unused = {"crs": crs, "energy_resource": energy_resource, "surface_format": netcdf}
if isinstance(obj, Farm):
unused = (
{"energy_resource": energy_resource, "surface_format": netcdf}
if kind == "wind_farm"
else {}
)
for option, value in unused.items():
if value is not None:
raise ValueError(f"{option} does not apply to a windIO {kind} document")
if netcdf is not None:
try:
import yaml # noqa: F401
except ImportError:
raise ValueError(
"surface_format='netcdf' needs PyYAML to write the !include of "
"the surface file"
) from None
if isinstance(obj, Turbine):
return _save_document(obj._to_windio_turbine(name=name), kind, path)
if not isinstance(obj, Farm):
return _save_resource(obj, path, name)
farm = obj
if name is None:
name = _farm_name(farm.metadata)
if farm.boundary is None and kind != "wind_farm":
raise ValueError(
f"a windIO {kind} document needs a boundary, and this farm has "
"none, so write it with document='wind_farm'"
)
resolved_crs = crs if crs is not None else farm.metadata.get("crs")
if resolved_crs is not None and not isinstance(resolved_crs, str):
raise TypeError("crs must be a string or None")
if kind == "wind_farm":
result = _build_wind_farm(farm, name, resolved_crs)
else:
own_surface = (
path.with_name(f"{path.stem}_surface.nc")
if surface_format == "netcdf"
else None
)
site = _build_site(
farm.boundary,
name,
resolved_crs,
_resolve_energy_resource(farm, name, energy_resource, shared),
getattr(farm, "surface", None),
shared,
own_surface,
)
result = (
site
if kind == "site"
else {
"name": name,
"site": site,
"wind_farm": _build_wind_farm(farm, name, resolved_crs),
}
)
return _save_document(result, kind, path)
def _first_layout(layouts) -> dict:
"""Return the first layout from single or multi-layout forms."""
if isinstance(layouts, list):
if not layouts:
raise ValueError("windIO wind_farm.layouts is empty")
return layouts[0]
return layouts
def _polygon_from_coordinates(entry) -> Polygon:
if not isinstance(entry, dict) or "x" not in entry or "y" not in entry:
raise ValueError("windIO boundary coordinates must contain x and y")
x = _number_list(entry["x"], "boundary coordinates x")
y = _number_list(entry["y"], "boundary coordinates y")
if len(x) != len(y):
raise ValueError("windIO boundary coordinate x/y lengths differ")
if not np.all(np.isfinite(x)) or not np.all(np.isfinite(y)):
raise ValueError("windIO boundary coordinates must be finite")
try:
return Polygon(np.column_stack([x, y]))
except (GEOSException, TypeError, ValueError) as error:
raise ValueError("windIO boundary coordinates do not form a polygon") from error
def _validate_boundary_geometry(geometry) -> None:
if not isinstance(geometry, (Polygon, MultiPolygon)):
raise ValueError("windIO boundary must be a Polygon or MultiPolygon")
if (
geometry.is_empty
or not geometry.is_valid
or not np.isfinite(geometry.area)
or geometry.area <= 0.0
):
raise ValueError(
"windIO boundary must be finite, nonempty, valid, and have positive area"
)
def _resolve_crs(site: dict | None, coordinates: dict | None) -> str | None:
entries = []
if site is not None:
entries.extend(site["boundaries"].get("polygons") or ())
entries.extend(site.get("exclusions", {}).get("polygons") or ())
for key in ("bathymetry", "elevation"):
surface = site.get(key)
if isinstance(surface, dict):
# The grid form, as a NetCDF surface reads, keeps its own CRS.
entries.append(surface.get("coordinates", surface))
if coordinates is not None:
entries.append(coordinates)
values = []
for entry in entries:
if isinstance(entry, dict) and entry.get("crs") is not None:
value = entry["crs"]
if not isinstance(value, str) or not value:
raise ValueError("windIO coordinate CRS must be a nonempty string")
values.append(value)
if len(set(values)) > 1:
raise ValueError("windIO boundary and layout coordinates have conflicting CRS")
return values[0] if values else None
def _surface_grid_from_windio(entry, kind, field_name):
"""Read the windIO grid form, with a 2-D field indexed ``(y, x)``."""
x = _number_list(entry["x"], f"{kind} x")
y = _number_list(entry["y"], f"{kind} y")
field = entry[field_name]
values = field.get("data") if isinstance(field, dict) else field
shape_message = f"windIO {kind} {field_name} must have shape (len(y), len(x))"
try:
values = np.asarray(values, dtype=float)
except (TypeError, ValueError):
raise ValueError(shape_message) from None
dims = field.get("dims") if isinstance(field, dict) else None
if dims is not None and list(dims) == ["x", "y"]:
values = values.T
if values.shape != (len(y), len(x)):
raise ValueError(shape_message)
return x, y, values
def _surface_from_windio(entry, kind):
"""Rebuild a regular surface from windIO point form or grid form."""
field_name = "depth" if kind == "bathymetry" else "elevation"
if not isinstance(entry, dict) or field_name not in entry:
raise ValueError(f"windIO {kind} is missing {field_name}")
from .surfaces import Bathymetry, Terrain
surface_type = Bathymetry if kind == "bathymetry" else Terrain
metadata = {"kind": kind, "source": "windio"}
if "coordinates" not in entry and {"x", "y"} <= set(entry):
return surface_type(
*_surface_grid_from_windio(entry, kind, field_name), metadata=metadata
)
coordinates = entry.get("coordinates")
if not isinstance(coordinates, dict) or not {"x", "y"} <= set(coordinates):
raise ValueError(f"windIO {kind} must contain coordinates with x and y")
x = _number_list(coordinates["x"], f"{kind} coordinates x")
y = _number_list(coordinates["y"], f"{kind} coordinates y")
values = _number_list(entry[field_name], f"{kind} {field_name}")
if not len(x) == len(y) == len(values):
raise ValueError(f"windIO {kind} needs one {field_name} per coordinate")
# OOWL surfaces are regular grids, so the points must fill one exactly.
xs, ys = np.unique(x), np.unique(y)
cells = np.searchsorted(ys, y) * len(xs) + np.searchsorted(xs, x)
if len(x) != len(xs) * len(ys) or np.any(np.bincount(cells, minlength=len(x)) != 1):
raise ValueError(
f"windIO {kind} points must fill a regular x-y grid, once each"
)
grid = np.empty(len(x), dtype=float)
grid[cells] = values
return surface_type(xs, ys, grid.reshape(len(ys), len(xs)), metadata=metadata)
def _circle_from_windio(circle) -> Polygon:
"""Return ``boundaries.circle`` as a circumscribed 64-gon containing the disc."""
try:
center = circle["center"]
x, y, radius = float(center["x"]), float(center["y"]), float(circle["radius"])
except (KeyError, TypeError, ValueError):
raise ValueError(
"windIO boundary circle needs a center with x and y and a radius"
) from None
if not all(np.isfinite((x, y, radius))) or radius <= 0:
raise ValueError("windIO boundary circle must be finite with a positive radius")
# A regular 64-gon with vertex radius R has inradius R cos(pi / 64).
return Point(x, y).buffer(radius / np.cos(np.pi / 64), quad_segs=16)
def _with_holes(polygons: list[Polygon], holes: list[Polygon]):
"""Return ``polygons`` with each hole as an interior ring, or ``None``.
``None`` means an exclusion is not strictly inside one polygon and must be
cut out instead.
"""
interiors: list[list] = [[] for _ in polygons]
for hole in holes:
owners = [i for i, polygon in enumerate(polygons) if polygon.contains(hole)]
if len(owners) != 1:
return None
interiors[owners[0]].append(hole.exterior.coords)
rebuilt = [
Polygon(polygon.exterior.coords, rings)
for polygon, rings in zip(polygons, interiors)
]
geometry = rebuilt[0] if len(rebuilt) == 1 else MultiPolygon(rebuilt)
return geometry if geometry.is_valid else None
def _geometry_from_site(site: dict):
"""Return the site boundary geometry with its exclusions removed."""
boundaries = site["boundaries"]
if "polygons" not in boundaries:
polygons = [_circle_from_windio(boundaries["circle"])]
else:
boundary_entries = boundaries["polygons"]
if not isinstance(boundary_entries, list) or not boundary_entries:
raise ValueError("windIO boundary polygons must be a nonempty list")
polygons = [_polygon_from_coordinates(entry) for entry in boundary_entries]
geometry = polygons[0] if len(polygons) == 1 else MultiPolygon(polygons)
exclusions = site.get("exclusions", {}).get("polygons")
if exclusions:
holes = [_polygon_from_coordinates(entry) for entry in exclusions]
rebuilt = _with_holes(polygons, holes)
if rebuilt is not None:
geometry = rebuilt
else:
try:
for hole in holes:
geometry = geometry.difference(hole)
except GEOSException as error:
raise ValueError(
"windIO exclusions do not define a valid boundary"
) from error
_validate_boundary_geometry(geometry)
return geometry
def _surface_from_site(site: dict, documents: _Documents):
"""Return the site's bathymetry or terrain surface, if it carries one."""
surface = None
for kind in ("bathymetry", "elevation"):
if kind in site:
if surface is not None:
raise ValueError(
"windIO site cannot contain both bathymetry and elevation"
)
surface = documents.build(
site[kind], kind, partial(_surface_from_windio, kind=kind)
)
return surface
def _turbines_from_wind_farm(wind_farm: dict, n_positions: int):
"""Read per-position turbine types or one turbine shared by the whole farm."""
from .turbines import _from_windio_turbines
layout = _first_layout(wind_farm["layouts"])
type_indices = layout.get("turbine_types")
turbine_types = wind_farm.get("turbine_types")
if turbine_types is not None and type_indices is not None:
return _from_windio_turbines(turbine_types, type_indices)
shared = wind_farm.get("turbines")
if shared is not None and n_positions:
# The shared turbine is type 0, so indices to any other type are errors.
indices = [0] * n_positions if type_indices is None else type_indices
return _from_windio_turbines({"0": shared}, indices)
if type_indices is not None:
raise ValueError(
"windIO layout turbine_types needs wind_farm.turbine_types or "
"wind_farm.turbines to define them"
)
return None
[docs]
def load(path: str | Path) -> "Farm | Turbine | WindResource":
"""Read a windIO document into the OOWL object it describes.
A ``wind_energy_system`` becomes a :class:`~oowl.layouts.base.Farm` with
boundary, resource, surface and turbines. A ``wind_farm`` document gives a
farm with only positions and turbines, and a ``site`` document a farm with
only boundary, resource and surface. A ``turbine`` document becomes a
:class:`~oowl.turbines.Turbine` and an ``energy_resource`` document a
:class:`~oowl.resources.WindResource`.
The boundary is the site boundaries minus exclusions, and positions come
from the first layout. A circle boundary becomes a circumscribed 64-gon, so
turbines on its rim stay inside. YAML ``!include`` files are read relative
to the including file, and a file included twice is read and built once,
so the objects that include it share it. A NetCDF include is read as a
surface grid, as ``save(..., surface_format="netcdf")`` writes it, and any
other include is refused. A CRS found in the coordinates is kept in the
farm metadata. Layout ``z`` coordinates are not read, because a farm takes
turbine heights from its terrain, so a ``wind_farm`` document saved from a
terrain farm reloads without them.
A resource loads from per-sector Weibull parameters, from a probability
table over direction and speed, from a rose with one wind speed for every
direction, which every case then blows at, or from a time series, counted
into 12 sectors and 1 m/s speed bins. A resource mapped over the site,
given per turbine or varying with height raises :class:`ValueError`.
windIO has no field for the stability climate, so a loaded
resource uses the default :class:`~oowl.resources.KellyGryning` and
ignores any stored ``LMO``. windIO also
cannot store :attr:`Farm.parts <oowl.Farm.parts>`, so a multipart farm
reloads as one farm. :meth:`Cluster.save <oowl.Cluster.save>` keeps farms
apart.
"""
return _load(path, _Documents())
def _load(path: str | Path, documents: _Documents) -> "Farm | Turbine | WindResource":
""":func:`load`, sharing included documents and objects through ``documents``."""
from .boundaries.base import _GeometryBoundary
from .layouts.base import Farm
path = Path(path)
document = _read_document(path, documents)
if not isinstance(document, dict):
raise ValueError(f"{path} did not contain a mapping document")
kind = _document_kind(document)
_validate_oowl_windio_document(document, kind)
from .resources import _from_windio_energy_resource
from .turbines import _from_windio_turbines
if kind == "turbine":
return _from_windio_turbines({"0": document}, [0])
if kind == "energy_resource":
resource = _from_windio_energy_resource(document)
if resource is None:
raise ValueError("windIO energy resource is a placeholder without wind")
return resource
if kind == "wind_energy_system":
site, wind_farm = document["site"], document["wind_farm"]
else:
site = document if kind == "site" else None
wind_farm = document if kind == "wind_farm" else None
boundary = None
resource = None
surface = None
if site is not None:
boundary = _GeometryBoundary(_geometry_from_site(site))
resource = documents.build(
site.get("energy_resource"),
"energy_resource",
_from_windio_energy_resource,
)
surface = _surface_from_site(site, documents)
coordinates = None
positions = np.empty((0, 2), dtype=float)
turbine_assignment = None
if wind_farm is not None:
coordinates = _first_layout(wind_farm["layouts"])["coordinates"]
positions = np.asarray(
np.column_stack([coordinates["x"], coordinates["y"]]), dtype=float
).reshape(-1, 2)
turbine_assignment = _turbines_from_wind_farm(wind_farm, len(positions))
metadata = {
"kind": "windio",
"document": kind,
"source": str(path),
"name": document.get("name"),
"modifiers": [],
}
resolved_crs = _resolve_crs(site, coordinates)
if resolved_crs is not None:
metadata["crs"] = resolved_crs
return Farm(
positions,
boundary,
metadata=metadata,
resource=resource,
turbines=turbine_assignment,
surface=surface,
)
