"""LINCOM Wind Point Data
A dataset that contains the results from a LINCOM simulation for the requested points.
These include the wind speed and direction, flow inclination angle, friction velocity,
6 horizontal wind derivatives, 3 vertical wind derivatives, and the shear exponent as
calculated by LINCOM. They are interpolated from the wind level, and are typically the
"Final" result from a LINCOM simulation.
"""
import xarray as xr
from windkit.spatial import create_dataset, get_crs, spatial_stack, spatial_unstack
from windkit.metadata import _LINCOM_WIND_POINT_ATTRS, update_var_attrs, update_history
from ..core import Rvea0305
FORT_CALC_WINDPTS = Rvea0305.calculate_points
[docs]
def get_wind_points(wind_level, out_ds):
"""Calculate point results from wind level for requested domain
.. warning::
This function is experimental and its signature may change.
Parameters
----------
wind_level: WindLevel
WindLevel object
out_ds: xarray.Dataset
Dataset containing the spatial structure to calculate
Result
------
Xarray.Dataset
xr.dataset containing results at the requested locations
"""
run_locs = spatial_stack(out_ds, wind_level.crs)
result = FORT_CALC_WINDPTS(
wind_level.nx_proj,
wind_level.ny_proj,
wind_level.orography,
wind_level.u,
wind_level.du_terrain,
wind_level.du_roughness,
wind_level.ustarp,
wind_level.iarr,
wind_level.rarr,
run_locs.west_east,
run_locs.south_north,
run_locs.height,
)
##############################
# Build the weibwc output file
##############################
result_names = (
"elevation",
"height",
"WS",
"WD",
"flow_inclination",
"USTAR",
"DU_DX",
"DU_DY",
"DU_DZ",
"DV_DX",
"DV_DY",
"DV_DZ",
"DTilt_DX",
"DTilt_DY",
"DTilt_DZ",
"ALPHA",
)
wind_pt = xr.Dataset()
for i, name in enumerate(result_names):
# First two variables are x and y coordinates so we skip them
wind_pt[name] = xr.DataArray(result[i + 2, :], dims=["point"])
wind_pt[name].attrs["_pwio_data_is_2d"] = False
# Add global attributes from run_locs
wind_pt.attrs = {**wind_pt.attrs, **run_locs.attrs}
# Assign the GIS coordinates
wind_pt = wind_pt.assign_coords(
west_east=run_locs.west_east, south_north=run_locs.south_north, crs=run_locs.crs
)
# Convert to raster if input was a raster
wind = spatial_unstack(wind_pt)
# Add global metadata
wind.attrs["title"] = "LINCOM wind point"
res = update_var_attrs(wind, _LINCOM_WIND_POINT_ATTRS)
return update_history(res)