diff --git a/pyproject.toml b/pyproject.toml index 2fd48af..26c39a5 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -11,7 +11,8 @@ dependencies = [ "pygef", "cems-nuclei==0.4.0.b2", "geopandas", - "pyogrio" + "pyogrio", + "scipy" ] license = { file = "LICENSE" } readme = "README.md" @@ -92,7 +93,8 @@ module = [ "shapely.*", "pygef.*", "mpl_toolkits.*", - "nuclei.*" + "nuclei.*", + "numpy.*" ] ignore_missing_imports = true diff --git a/requirements.txt b/requirements.txt index b7ca1ad..6d3b4ce 100644 --- a/requirements.txt +++ b/requirements.txt @@ -86,7 +86,7 @@ gef-file-to-map==0.1.0 # via pygef geographiclib==2.0 # via geopy -geopandas==0.14.0 +geopandas==0.14.1 # via py-vibracore (pyproject.toml) geopy==2.4.0 # via contextily @@ -98,17 +98,24 @@ importlib-metadata==6.8.0 # via # fiona # sphinx -importlib-resources==6.1.0 +importlib-resources==6.1.1 # via matplotlib iniconfig==2.0.0 # via pytest +ipympl==0.9.3 + # via py-vibracore (pyproject.toml) ipython==8.11.0 # via # black + # ipympl # ipywidgets # py-vibracore (pyproject.toml) +ipython-genutils==0.2.0 + # via ipympl ipywidgets==8.1.1 - # via tqdm + # via + # ipympl + # tqdm isort==5.12.0 # via py-vibracore (pyproject.toml) jedi==0.19.1 @@ -128,6 +135,7 @@ markupsafe==2.1.3 matplotlib==3.8.1 # via # contextily + # ipympl # pygef matplotlib-inline==0.1.6 # via ipython @@ -144,15 +152,17 @@ mypy-extensions==1.0.0 # black # mypy # py-vibracore (pyproject.toml) -numpy==1.26.1 +numpy==1.26.2 # via # contourpy + # ipympl # matplotlib # pandas # pandas-stubs # pyarrow # pyogrio # rasterio + # scipy # shapely # snuggs packaging==23.2 @@ -164,7 +174,7 @@ packaging==23.2 # pyogrio # pytest # sphinx -pandas==2.1.2 +pandas==2.1.3 # via geopandas pandas-stubs==2.1.1.230928 # via py-vibracore (pyproject.toml) @@ -181,6 +191,7 @@ pickleshare==0.7.5 pillow==10.1.0 # via # contextily + # ipympl # matplotlib platformdirs==3.5.1 # via @@ -188,17 +199,17 @@ platformdirs==3.5.1 # py-vibracore (pyproject.toml) pluggy==1.3.0 # via pytest -polars[pyarrow]==0.19.12 +polars[pyarrow]==0.19.13 # via # polars # pygef -prompt-toolkit==3.0.39 +prompt-toolkit==3.0.40 # via ipython ptyprocess==0.7.0 # via pexpect pure-eval==0.2.2 # via stack-data -pyarrow==14.0.0 +pyarrow==14.0.1 # via polars pycodestyle==2.10.0 # via @@ -239,7 +250,12 @@ requests==2.31.0 # cems-nuclei # contextily # coveralls + # requests-mock # sphinx +requests-mock==1.11.0 + # via py-vibracore (pyproject.toml) +scipy==1.11.3 + # via py-vibracore (pyproject.toml) shapely==2.0.2 # via geopandas six==1.16.0 @@ -247,6 +263,7 @@ six==1.16.0 # asttokens # fiona # python-dateutil + # requests-mock snowballstemmer==2.2.0 # via sphinx snuggs==1.4.7 @@ -296,12 +313,13 @@ tqdm[notebook]==4.66.1 traitlets==5.13.0 # via # comm + # ipympl # ipython # ipywidgets # matplotlib-inline types-pytz==2023.3.1.1 # via pandas-stubs -types-tqdm==4.66.0.3 +types-tqdm==4.66.0.4 # via py-vibracore (pyproject.toml) typing-extensions==4.7.1 # via diff --git a/src/pyvibracore/__init__.py b/src/pyvibracore/__init__.py index 5967dff..a98c6c1 100644 --- a/src/pyvibracore/__init__.py +++ b/src/pyvibracore/__init__.py @@ -1,11 +1,5 @@ from ._version import __version__ -from .input import impact_force_properties, vibration_properties -from .results import impact_force_result, vibration_result __all__ = [ "__version__", - "impact_force_properties", - "impact_force_result", - "vibration_result", - "vibration_properties", ] diff --git a/src/pyvibracore/results/plot_utils.py b/src/pyvibracore/results/plot_utils.py new file mode 100644 index 0000000..eda5343 --- /dev/null +++ b/src/pyvibracore/results/plot_utils.py @@ -0,0 +1,35 @@ +from __future__ import annotations + +from matplotlib import pyplot as plt +from mpl_toolkits.axes_grid1.anchored_artists import AnchoredSizeBar + + +def _north_arrow(axes: plt.Axes) -> None: + """Add north arrow to axes""" + x, y, arrow_length = 0.05, 0.98, 0.1 + axes.annotate( + "N", + xy=(x, y), + xytext=(x, y - arrow_length), + arrowprops=dict(facecolor="black", width=5, headwidth=15), + ha="center", + va="center", + fontsize=20, + xycoords=axes.transAxes, + ) + + +def _scalebar(axes: plt.Axes) -> None: + """Add size bar to axes""" + scalebar = AnchoredSizeBar( + axes.transData, + 20, + "20 m", + "lower left", + pad=1, + color="black", + frameon=True, + size_vertical=2, + ) + + axes.add_artist(scalebar) diff --git a/src/pyvibracore/results/sound_result.py b/src/pyvibracore/results/sound_result.py new file mode 100644 index 0000000..6a40247 --- /dev/null +++ b/src/pyvibracore/results/sound_result.py @@ -0,0 +1,238 @@ +from __future__ import annotations + +import logging +from typing import Any, List, Tuple + +import geopandas as gpd +import matplotlib.patches as patches +import matplotlib.pyplot as plt +import numpy as np +from numpy.typing import NDArray +from scipy import interpolate +from shapely.geometry import LineString, Point, Polygon + +from pyvibracore.results.plot_utils import _north_arrow, _scalebar + + +def _sound_prediction( + power: float, k2: float, period: float, levels: List[float] +) -> NDArray: + distance = np.arange(1e-5, 150, step=0.2) + noise = ( + power + - (-10 * np.log(period / 12)) + - (20 * np.log(distance) + 0.005 * distance + 9.1) + + k2 + ) + + f = interpolate.interp1d(noise, distance, kind="cubic") + return f(levels) + + +def map_sound( + buildings: gpd.GeoDataFrame, + source_location: Point | LineString | Polygon, + building_name: str, + power: float, + k2: float, + period: float, + title: str = "Legend:", + figsize: Tuple[float, float] = (10.0, 12.0), + settings: dict | None = None, + **kwargs: Any, +) -> plt.Figure: + """ + Create map of the input building settings. + + Parameters + ---------- + buildings: + GeoDataFrame of the input buildings + source_location: + location of the vibration source + building_name: + name of the building + power: + source power [dB] + k2: + Correction term [dB] + period: + Operating period of the building code [hours] + title: + Legend title + figsize: + Size of the activate figure, as the `plt.figure()` argument. + settings: + Plot settings used in plot: default settings are: + + { + "source_location": {"label": "Trillingsbron", "color": "blue"}, + "levels": [ + { + "label": ">80 db [0 dagen]", + "level": 80, + "color": "darkred", + }, + { + "label": ">75 db [5 dagen]", + "level": 75, + "color": "red", + }, + { + "label": ">70 db [15 dagen]", + "level": 70, + "color": "orange", + }, + { + "label": ">65 db [30 dagen]", + "level": 65, + "color": "darkgreen", + }, + { + "label": ">60 db [50 dagen]", + "level": 60, + "color": "green", + }, + ], + } + **kwargs: + All additional keyword arguments are passed to the `pyplot.subplots()` call. + + Returns + ------- + Figure + """ + if settings is None: + settings = { + "source_location": {"label": "Trillingsbron", "color": "blue"}, + "levels": [ + { + "label": ">80 db [0 dagen]", + "level": 80, + "color": "darkred", + }, + { + "label": ">75 db [5 dagen]", + "level": 75, + "color": "red", + }, + { + "label": ">70 db [15 dagen]", + "level": 70, + "color": "orange", + }, + { + "label": ">65 db [30 dagen]", + "level": 65, + "color": "darkgreen", + }, + { + "label": ">60 db [50 dagen]", + "level": 60, + "color": "lightgreen", + }, + ], + } + + kwargs_subplot = { + "figsize": figsize, + "tight_layout": True, + } + + kwargs_subplot.update(kwargs) + + fig, axes = plt.subplots(**kwargs_subplot) + + gpd.GeoSeries(source_location).plot( + ax=axes, color=settings["source_location"]["color"], alpha=1, zorder=1, aspect=1 + ) + + building = buildings.get(buildings["name"] == building_name) + if building.empty: + raise ValueError(f"No buildings with name {building_name}.") + + building.plot(ax=axes, zorder=2, color="gray", aspect=1) + buildings.where(buildings["name"] != building_name).plot( + ax=axes, zorder=2, color="lightgray", aspect=1 + ) + + # plot contour + levels = [values["level"] for values in settings["levels"]] + distances = _sound_prediction(power, k2, period, levels=levels) + colors = [values["color"] for values in settings["levels"]] + for distance, color in zip(distances, colors): + gpd.GeoSeries(building.buffer(distance).exterior).plot( + ax=axes, zorder=3, color=color, aspect=1 + ) + + # plot name + for idx, row in buildings.iterrows(): + x = row.geometry.centroid.xy[0][0] + y = row.geometry.centroid.xy[1][0] + + axes.annotate( + idx, + xy=(x, y), + horizontalalignment="center", + ) + + # add legend + axes.legend( + title=title, + title_fontsize=18, + fontsize=15, + loc="lower right", + handles=[ + patches.Patch( + facecolor=value["color"], + label=value["label"], + alpha=0.9, + linewidth=2, + edgecolor="black", + ) + for value in settings["levels"] + ], + ) + + _north_arrow(axes) + _scalebar(axes) + + return fig + + +def get_normative_building( + buildings: gpd.GeoDataFrame, + location: Polygon | LineString | Point, +) -> str | None: + """ + Get the name of the closest building with one of the follwing category: + + - "woonfunctie", + - "gezondheidsfunctie", + - "onderwijsfunctie" + + Parameters + ---------- + buildings: + GeoDataFrame that holds the building information + location: + Geometry of the source location + + Returns + ------- + name: str + """ + + category = ["woonfunctie", "gezondheidsfunctie", "onderwijsfunctie"] + + gdf = buildings.get( + [ + any(item in category for item in row.split(",")) if row else False + for row in buildings["gebruiksdoel"].to_list() + ] + ) + if gdf.empty: + logging.error(f"ValueError: No buildings with category {category}.") + return None + gdf["distance"] = gdf.distance(location) + return gdf.sort_values("distance", na_position="last").iloc[0].get("name") diff --git a/src/pyvibracore/results/vibration_result.py b/src/pyvibracore/results/vibration_result.py index 48dd1fa..580508e 100644 --- a/src/pyvibracore/results/vibration_result.py +++ b/src/pyvibracore/results/vibration_result.py @@ -8,39 +8,9 @@ import matplotlib.patches as patches import matplotlib.pyplot as plt import numpy as np -from mpl_toolkits.axes_grid1.anchored_artists import AnchoredSizeBar from shapely.geometry import LineString, Point, Polygon - -def _north_arrow(axes: plt.Axes) -> None: - """Add north arrow to axes""" - x, y, arrow_length = 0.05, 0.98, 0.1 - axes.annotate( - "N", - xy=(x, y), - xytext=(x, y - arrow_length), - arrowprops=dict(facecolor="black", width=5, headwidth=15), - ha="center", - va="center", - fontsize=20, - xycoords=axes.transAxes, - ) - - -def _scalebar(axes: plt.Axes) -> None: - """Add size bar to axes""" - scalebar = AnchoredSizeBar( - axes.transData, - 20, - "20 m", - "lower left", - pad=1, - color="black", - frameon=True, - size_vertical=2, - ) - - axes.add_artist(scalebar) +from pyvibracore.results.plot_utils import _north_arrow, _scalebar @dataclass(frozen=True) @@ -219,12 +189,16 @@ def map_payload( { "source_location": {"label": "Trillingsbron", "color": "black"}, - "insufficient_cat1": { - "label": "Voldoet Niet - Cat.1", - "color": "orange", + "sensitive_cat1": { + "label": "Monumentaal/ gevoelig - Cat.1", + "color": "blue", + }, + "sensitive_cat2": { + "label": "Monumentaal/ gevoelig - Cat.2", + "color": "cyan", }, - "insufficient_cat2": {"label": "Voldoet Niet - Cat.2", "color": "red"}, - "sufficient": {"label": "Voldoet", "color": "green"}, + "normal_cat1": {"label": "Normaal - Cat.1", "color": "orange"}, + "normal_cat2": {"label": "Normaal - Cat.2", "color": "olive"}, } **kwargs: All additional keyword arguments are passed to the `pyplot.subplots()` call. diff --git a/tests/test_sound.py b/tests/test_sound.py new file mode 100644 index 0000000..5081227 --- /dev/null +++ b/tests/test_sound.py @@ -0,0 +1,38 @@ +import geopandas as gpd +import matplotlib.pyplot as plt + +from pyvibracore.input.vibration_properties import ( + BAG_WFS_URL, + get_buildings_geodataframe, +) +from pyvibracore.results.sound_result import get_normative_building, map_sound + + +def test_map_sound(requests_mock, mock_bag_response, mock_source_location): + requests_mock.get(BAG_WFS_URL, json=mock_bag_response) + gdf = get_buildings_geodataframe(1, 2, 3, 4) + + fig = map_sound( + gdf, + source_location=mock_source_location, + building_name="0", + power=140, + k2=5, + period=5, + ) + + plt.savefig("tpm.png") + + assert isinstance(fig, plt.Figure) + + +def test_get_normative_building(requests_mock, mock_bag_response, mock_source_location): + requests_mock.get(BAG_WFS_URL, json=mock_bag_response) + gdf = get_buildings_geodataframe(1, 2, 3, 4) + + name = get_normative_building( + gdf, + location=mock_source_location, + ) + + assert isinstance(name, str)