diff --git a/.github/workflows/test.yml b/.github/workflows/test.yml index 02032d5b..4cb38130 100644 --- a/.github/workflows/test.yml +++ b/.github/workflows/test.yml @@ -23,6 +23,7 @@ jobs: - name: test run: | + docker compose run -e QT_QPA_PLATFORM=offscreen qgis-desktop make check docker compose run -e QT_QPA_PLATFORM=offscreen qgis-desktop make test docker compose run -e QT_QPA_PLATFORM=offscreen qgis-desktop make flake8 docker compose run qgis-desktop make docstrings diff --git a/.gitignore b/.gitignore index 9dc00eac..c82f40c2 100644 --- a/.gitignore +++ b/.gitignore @@ -26,3 +26,4 @@ tests/data/small_2019_01_2entries_4tiff.sqlite .pytest_cache/v/cache/nodeids .h5py_marker constraints.txt +requirements.in diff --git a/CHANGES.rst b/CHANGES.rst index e10e7dfc..d65fadf7 100755 --- a/CHANGES.rst +++ b/CHANGES.rst @@ -1,12 +1,44 @@ 3Di Results Analysis changelog ======================== -3.9.4 (unreleased) ------------------- +3.13 (unreleased) +----------------- - Nothing changed yet. +3.12 (2024-12-10) +----------------- + +- Compatibility with Python 3.12 (#1061) +- Bumped h5py to 3.10.0 and scipy to 1.13.0 for python 3.12 compatibility (#1061) +- Added Processing Algorithm "Extract structure control actions" (#926) +- Fixed attributeError when loading a QGIS project (#1063) +- Fix in Rasters to NetCDF algorithm to properly convert the units Enum to string (#1067) + +3.11 (2024-11-12) +----------------- + +- Allow LinestringZ input by using WKB instead of WKT as conversion format (#1040) +- Make "Use selected features" behaviour explicit (#1057) +- Added animation settings (#1046) +- Bump hydxlib to 1.5.3 + + +3.10.0 (2024-09-12) +------------------- + +- Add pump support to Result Aggregation tool +- Add preset "Total pumped volume" to Result Aggregation tool +- Add preset "Pumps: % of time at max capacity" to Result Aggregation tool +- Cross-sectional discharge algorithm: allow cross-section lines to have different crs than 3Di results, automatically reproject +- Watershed tool: 2D flowlines intersecting obstacles are no longer shown as 1D flowlines (#1034) +- Model selection dialog: fixed order bug when sorting. +- Bump threedi-mi-utils to 0.1.4 +- Fixed attribute error when selecting substance in a model that contains pumps (#1044) +- Added Flow Summary Tool (#725) + + 3.9.3 (2024-08-14) ------------------ diff --git a/Dockerfile b/Dockerfile index 832cb5b0..d95ebd23 100644 --- a/Dockerfile +++ b/Dockerfile @@ -1,4 +1,4 @@ -from qgis/qgis:final-3_28_4 +FROM qgis/qgis:final-3_28_4 RUN apt-get update && apt-get install -y python3-pyqt5.qtwebsockets wget python3-scipy python3-h5py zip && apt-get clean # RUN mkdir -p /tests_directory COPY requirements-dev.txt /root @@ -9,6 +9,6 @@ RUN qgis_setup.sh # Copied the original PYTHONPATH and added the profile's python dir to # imitate qgis' behaviour. -ENV PYTHONPATH /usr/share/qgis/python/:/usr/share/qgis/python/plugins:/usr/lib/python3/dist-packages/qgis:/usr/share/qgis/python/qgis:/root/.local/share/QGIS/QGIS3/profiles/default/python +ENV PYTHONPATH=/usr/share/qgis/python/:/usr/share/qgis/python/plugins:/usr/lib/python3/dist-packages/qgis:/usr/share/qgis/python/qgis:/root/.local/share/QGIS/QGIS3/profiles/default/python # Note: we'll mount the current dir into this WORKDIR WORKDIR /root/.local/share/QGIS/QGIS3/profiles/default/python/plugins/threedi_results_analysis diff --git a/Makefile b/Makefile index 5adc9943..1a19f2a1 100755 --- a/Makefile +++ b/Makefile @@ -84,6 +84,12 @@ zip: compile find /tmp/$(PLUGINNAME) -iname "*.pyc" -delete cd /tmp; zip -9r $(CURDIR)/$(PLUGINNAME).zip $(PLUGINNAME) +check: constraints.txt + # Use pip-compile to check whether all dependencies version constraints are met. + cp constraints.txt requirements.in + pip-compile --dry-run + rm requirements.in + package: compile # Create a zip package of the plugin named $(PLUGINNAME).zip. # This requires use of git (your plugin development directory must be a diff --git a/datasource/result_constants.py b/datasource/result_constants.py index a5baf7e3..6bd1592e 100644 --- a/datasource/result_constants.py +++ b/datasource/result_constants.py @@ -135,3 +135,14 @@ # TODO: QH is also defined above. LAYER_OBJECT_TYPE_MAPPING = dict([(a[0], a[1]) for a in layer_information]) + +""" + Maps action type to affected structure. + This can be used when no control action is present for an object. +""" +ACTION_TYPE_ATTRIBUTE_MAP = { + "set_crest_level": {"variable": "dpumax", "unit": "m MSL", "applicable_structures": ["v2_weir", "v2_orifice"]}, + "set_pump_capacity": {"variable": "capacity", "unit": "", "applicable_structures": None}, + "set_discharge_coefficients": {"variable": "discharge_coefficient_positive", "unit": "", "applicable_structures": None}, + "set_gate_level": {"variable": None, "unit": "m MSL", "applicable_structures": None} +} diff --git a/datasource/test_datasources.py b/datasource/test_datasources.py index 7b1563b7..c790fd3a 100644 --- a/datasource/test_datasources.py +++ b/datasource/test_datasources.py @@ -1,9 +1,9 @@ -from threedigrid.admin import gridresultadmin -from threedigrid.admin.constants import NO_DATA_VALUE from threedi_results_analysis.datasource.threedi_results import find_aggregation_netcdf from threedi_results_analysis.datasource.threedi_results import normalized_object_type from threedi_results_analysis.datasource.threedi_results import ThreediResult from threedi_results_analysis.tests.utilities import TemporaryDirectory +from threedigrid.admin import gridresultadmin +from threedigrid.admin.constants import NO_DATA_VALUE import mock import numpy as np @@ -110,16 +110,6 @@ def test_get_timeseries_filter_node(threedi_result): np.testing.assert_equal(time_series[:, 1], data.return_value[:, 0]) -def test_get_timeseries_filter_content_pk(threedi_result): - with mock.patch( - "threedigrid.orm.base.models.Model.get_filtered_field_value" - ) as data: - data.return_value = np.ones((len(threedi_result.timestamps), 1)) - time_series = threedi_result.get_timeseries("s1", content_pk=5) - np.testing.assert_equal(time_series[:, 0], threedi_result.get_timestamps()) - np.testing.assert_equal(time_series[:, 1], data.return_value[:, 0]) - - def test_get_timeseries_filter_fill_value(threedi_result): with mock.patch( "threedigrid.orm.base.models.Model.get_filtered_field_value" diff --git a/datasource/threedi_results.py b/datasource/threedi_results.py index 1695f92c..7825e776 100644 --- a/datasource/threedi_results.py +++ b/datasource/threedi_results.py @@ -1,11 +1,19 @@ from functools import cached_property -from threedigrid.admin.constants import NO_DATA_VALUE -from threedi_results_analysis.datasource.result_constants import LAYER_OBJECT_TYPE_MAPPING +from threedi_results_analysis.datasource.result_constants import ( + ACTION_TYPE_ATTRIBUTE_MAP, +) +from threedi_results_analysis.datasource.result_constants import ( + LAYER_OBJECT_TYPE_MAPPING, +) from threedi_results_analysis.datasource.result_constants import SUBGRID_MAP_VARIABLES +from threedigrid.admin.constants import NO_DATA_VALUE from threedigrid.admin.gridadmin import GridH5Admin from threedigrid.admin.gridresultadmin import GridH5AggregateResultAdmin from threedigrid.admin.gridresultadmin import GridH5ResultAdmin +from threedigrid.admin.gridresultadmin import GridH5StructureControl from threedigrid.admin.gridresultadmin import GridH5WaterQualityResultAdmin +from threedigrid.admin.structure_controls.models import StructureControlSourceTypes +from threedigrid.admin.structure_controls.models import StructureControlTypes import glob import h5py @@ -37,6 +45,7 @@ class ThreediResult(): - GridH5ResultAdmin - GridH5AggregateResultAdmin - GridH5WaterQualityResultAdmin + - GridH5StructureControl For more information about threedigrid see https://threedigrid.readthedocs.io/en/latest/ @@ -119,7 +128,30 @@ def available_water_quality_vars(self): @property def available_vars(self): """Return a list of all available variables""" - return self.available_subgrid_map_vars + self.available_aggregation_vars + self.available_water_quality_vars + return self.available_subgrid_map_vars + self.available_aggregation_vars + self.available_water_quality_vars + self.available_structure_control_actions_vars + + @property + def available_structure_control_actions_vars(self): + """Return a list of all structure control actions variables""" + ga = self.structure_control_actions_result_admin + if not ga: + return [] + available_vars = [] + for control_type in StructureControlTypes.__members__.values(): + control_type_data = getattr(ga, control_type.name) + action_types = np.unique(control_type_data.action_type) + for action_type in action_types: + source_types = [cta.source_type.value for cta in control_type_data.group_by_action_type(action_type)] + var = { + "name": action_type[4:].replace("_", " ").capitalize(), # sanitize + "unit": ACTION_TYPE_ATTRIBUTE_MAP[action_type]["unit"], + "parameters": action_type, + "types": source_types + } + if var not in available_vars: + available_vars.append(var) + + return available_vars @cached_property def timestamps(self): @@ -180,13 +212,13 @@ def get_gridadmin(self, variable=None): """Return the gridadmin where the variable is stored. If no variable is given, a gridadmin without results is returned. - Results are either stored in the 'results_3di.nc', 'aggregate_results_3di.nc' - or 'water_quality_results_3di.nc'. These make use of the GridH5ResultAdmin, - GridH5AggregateResultAdmin or GridH5WaterQualityResultAdmin to query the data + Results are either stored in the 'results_3di.nc', 'aggregate_results_3di.nc', + 'water_quality_results_3di.nc' or 'structure_control_actions_3di.nc'. These make use of the GridH5ResultAdmin, + GridH5AggregateResultAdmin, GridH5WaterQualityResultAdmin or GridH5StructureControl to query the data respectively. :param variable: str of the variable name, e.g. 's1', 'q_pump' - :return: handle to GridAdminResult, AggregateGridAdminResult or GridH5WaterQualityResultAdmin + :return: handle to GridAdminResult, AggregateGridAdminResult, GridH5WaterQualityResultAdmin or GridH5StructureControl """ if variable is None: return self.gridadmin @@ -196,17 +228,19 @@ def get_gridadmin(self, variable=None): return self.aggregate_result_admin elif variable in [v["parameters"] for v in self.available_water_quality_vars]: return self.water_quality_result_admin + elif variable in [v["parameters"] for v in self.available_structure_control_actions_vars]: + return self.structure_control_actions_result_admin else: raise AttributeError(f"Unknown subgrid or aggregate or water quality variable: {variable}") def get_timeseries( - self, nc_variable, node_id=None, content_pk=None, fill_value=None + self, nc_variable, node_id=None, fill_value=None, selected_object_type=None ): """Return a time series array of the given variable A 2d array is given, with first column being the timestamps in seconds. The next columns are the values of the nodes of the given variable. - You can also filter on a specific node using node_id or content_pk, + You can also filter on a specific node using node_id, in which case only the timeseries of the given node is returned. If there is no values of the given variable, only the timestamps are @@ -214,27 +248,112 @@ def get_timeseries( :param nc_variable: :param node_id: - :param content_pk: :param fill_value: + :param selected_object_type: layer type of selected feature :return: 2D array, first column being the timestamps """ ga = self.get_gridadmin(nc_variable) - filtered_result = ga.get_model_instance_by_field_name(nc_variable).timeseries( - indexes=slice(None) - ) - if node_id: - filtered_result = filtered_result.filter(id=node_id) - elif content_pk: - filtered_result = filtered_result.filter(content_pk=content_pk) + if isinstance(ga, GridH5StructureControl): + # GridH5StructureControl has a different interface compared to the other GridAdmin structures + return self.get_structure_control_action_timeseries(ga, nc_variable, node_id, selected_object_type, fill_value) + else: + filtered_result = ga.get_model_instance_by_field_name(nc_variable).timeseries( + indexes=slice(None) + ) + if node_id: + filtered_result = filtered_result.filter(id=node_id) + values = self.get_timeseries_values(filtered_result, nc_variable) + if fill_value is not None: + values[values == NO_DATA_VALUE] = fill_value + + timestamps = self.get_timestamps(nc_variable) + timestamps = timestamps.reshape(-1, 1) # reshape (n,) to (n, 1) + + return np.hstack([timestamps, values]) + + def get_structure_control_action_timeseries(self, ga, nc_variable, node_id, selected_object_type, fill_value): + assert nc_variable + assert node_id + timestamps = [] + values = [] + + prev_value = None + for control_type in StructureControlTypes.__members__.values(): + control_type_data = getattr(ga, control_type.name) + structure_controls_for_id = control_type_data.group_by_grid_id(node_id) + structure_controls = [sc for sc in structure_controls_for_id if sc.action_type == nc_variable] + + # It could be that the same action is applied on nodes, lines and pumps, we need to find the right one. + desired_type = StructureControlSourceTypes.LINES if selected_object_type == "flowline" else StructureControlSourceTypes.PUMPS + structure_controls = [sc for sc in structure_controls_for_id if sc.source_type == desired_type] + + for structure_control in structure_controls: + assert len(structure_control.time) == len(structure_control.action_value_1) + for time_key, value in zip(structure_control.time, structure_control.action_value_1): + if prev_value is not None: + timestamps.append(time_key) + values.append(prev_value) + prev_value = value + timestamps.append(time_key) + values.append(value) + + # Retrieve gridadmin structure + if selected_object_type == "flowline": + structure = self.gridadmin.lines.filter(id=node_id) + elif selected_object_type == "pump": + structure = self.gridadmin.pumps.filter(id=node_id) + else: + raise NotImplementedError("Plotting node control actions is not yet implemented") + + # Check whether this object's content type is applicable for this control action + applicable_structures = ACTION_TYPE_ATTRIBUTE_MAP[nc_variable]["applicable_structures"] + if applicable_structures: + content_type = structure.content_type[0].decode() + if content_type not in applicable_structures: + # This action is not applicable to this object + logger.info(f"Parameter {nc_variable} not applicable for type {str(content_type)}") + return np.column_stack(([], [])) + + orig_timestamps = self.get_timestamps() + max_time_stamp_orig = max(orig_timestamps) + min_time_stamp_orig = min(orig_timestamps) + if timestamps: + assert values + max_time_stamp_structure = max(timestamps) + # Possibly append the structure control actions till the end + if values and max_time_stamp_orig > max_time_stamp_structure: + values.append(values[-1]) + timestamps.append(max_time_stamp_structure) + values.append(values[-1]) + timestamps.append(max_time_stamp_orig) + + affected_nc_variable = ACTION_TYPE_ATTRIBUTE_MAP[nc_variable]["variable"] + if affected_nc_variable: + # Check if we need to prepend the plot with non-controlled (static) values + orig_value = getattr(structure, affected_nc_variable)[0] + if not timestamps: + # No actions at all, take original value to plot + assert not values + values = [orig_value, orig_value] + timestamps = [min_time_stamp_orig, max_time_stamp_orig] + else: + min_time_stamp_structure = min(timestamps) + if min_time_stamp_orig < min_time_stamp_structure: + values.insert(0, orig_value) + timestamps.insert(0, min_time_stamp_structure) + values.insert(0, orig_value) + timestamps.insert(0, min_time_stamp_orig) + + if not values: + return np.column_stack(([], [])) - values = self.get_timeseries_values(filtered_result, nc_variable) if fill_value is not None: - values[values == NO_DATA_VALUE] = fill_value + for index in range(len(values)): + if values[index] == NO_DATA_VALUE: + values[index] = fill_value - timestamps = self.get_timestamps(nc_variable) - timestamps = timestamps.reshape(-1, 1) # reshape (n,) to (n, 1) - return np.hstack([timestamps, values]) + return np.column_stack((timestamps, values)) def get_values_by_timestep_nr(self, variable, timestamp_idx, node_ids): """Return an array of values of the given variable on the specified @@ -317,7 +436,7 @@ def result_admin(self): # Note: passing a file-like object due to an issue in threedigrid # https://github.com/nens/threedigrid/issues/183 file_like_object_h5 = open(h5, 'rb') - file_like_object_h5.startswith = lambda x: '' + file_like_object_h5.startswith = lambda x: False file_like_object_nc = open(self.file_path, 'rb') return GridH5ResultAdmin(file_like_object_h5, file_like_object_nc) @@ -353,6 +472,22 @@ def water_quality_result_admin(self): file_like_object_nc = open(wq_path, 'rb') return GridH5WaterQualityResultAdmin(file_like_object_h5, file_like_object_nc) + @cached_property + def structure_control_actions_result_admin(self): + try: + # Note: both of these might raise the FileNotFoundError + sca_path = find_structure_control_actions_netcdf(self.file_path) + h5 = self.h5_path + except FileNotFoundError: + logger.exception("Structure control actions result not found") + return None + # Note: passing a file-like object due to an issue in threedigrid + # https://github.com/nens/threedigrid/issues/183 + file_like_object_h5 = open(h5, 'rb') + file_like_object_h5.startswith = lambda x: False + file_like_object_nc = open(sca_path, 'rb') + return GridH5StructureControl(file_like_object_h5, file_like_object_nc) + @property def short_model_slug(self): model_slug = self.gridadmin.model_slug @@ -388,6 +523,28 @@ def find_aggregation_netcdf(netcdf_file_path): ) +def find_structure_control_actions_netcdf(netcdf_file_path): + """An ad-hoc way to find the structure control actions netcdf file + + Args: + netcdf_file_path: path to the result netcdf + + Returns: + the structure control actions netcdf path + + Raises: + FileNotFoundError if nothing is found + """ + pattern = "structure_control_actions_3di.nc" + result_dir = os.path.dirname(netcdf_file_path) + sca_result_files = glob.glob(os.path.join(result_dir, pattern)) + if sca_result_files: + return sca_result_files[0] + raise FileNotFoundError( + "'structure_control_actions_3di.nc' file not found relative to %s" % result_dir + ) + + def find_water_quality_netcdf(netcdf_file_path): """An ad-hoc way to find the water quality netcdf file diff --git a/dependencies.py b/dependencies.py index ac4443fa..ac12dad9 100644 --- a/dependencies.py +++ b/dependencies.py @@ -75,7 +75,7 @@ Dependency("threedigrid-builder", "threedigrid_builder", "==1.17.*", False), Dependency("h5netcdf", "h5netcdf", "", False), Dependency("greenlet", "greenlet", "!=0.4.17", False), - Dependency("threedi-mi-utils", "threedi_mi_utils", "==0.1.2", False), + Dependency("threedi-mi-utils", "threedi_mi_utils", "==0.1.4", False), ] # On Windows, the hdf5 binary and thus h5py version depends on the QGis version @@ -84,14 +84,24 @@ if QGIS_VERSION < 32806 and platform.system() == "Windows": SUPPORTED_HDF5_VERSIONS = ["1.10.7"] H5PY_DEPENDENCY = Dependency("h5py", "h5py", "==2.10.0", False) +elif QGIS_VERSION >= 34000 and platform.system() == "Windows": + SUPPORTED_HDF5_VERSIONS = ["1.14.0"] + H5PY_DEPENDENCY = Dependency("h5py", "h5py", "==3.10.0", False) else: SUPPORTED_HDF5_VERSIONS = ["1.14.0"] H5PY_DEPENDENCY = Dependency("h5py", "h5py", "==3.8.0", True) -WINDOWS_PLATFORM_DEPENDENCIES = [Dependency("scipy", "scipy", "==1.6.2", False)] -if QGIS_VERSION >= 32811 and platform.system() == "Windows": +if QGIS_VERSION < 32811 and platform.system() == "Windows": + WINDOWS_PLATFORM_DEPENDENCIES = [ + Dependency("scipy", "scipy", "==1.6.2", True), + ] +elif QGIS_VERSION >= 34000 and platform.system() == "Windows": + WINDOWS_PLATFORM_DEPENDENCIES = [ + Dependency("scipy", "scipy", "==1.13.0", True), + ] +else: WINDOWS_PLATFORM_DEPENDENCIES = [ - Dependency("scipy", "scipy", "==1.10.1", True), + Dependency("scipy", "scipy", "==1.10.1", False), ] # If you add a dependency, also adjust external-dependencies/populate.sh @@ -104,7 +114,7 @@ def create_progress_dialog(progress, text): dialog = QProgressDialog() - dialog.setWindowTitle("3Di Toolbox install progress") + dialog.setWindowTitle("3Di Results Analysis install progress") dialog.setLabelText(text) dialog.setWindowFlags(Qt.WindowStaysOnTopHint) bar = QProgressBar(dialog) @@ -523,7 +533,7 @@ def _install_dependencies(dependencies, target_dir): # sticking around. if bar: - bar.setValue((count / len(dependencies)) * 100) + bar.setValue(int((count / len(dependencies)) * 100)) bar.update() QApplication.processEvents() diff --git a/external-dependencies/populate.sh b/external-dependencies/populate.sh index 6a46c5b0..6de70624 100755 --- a/external-dependencies/populate.sh +++ b/external-dependencies/populate.sh @@ -40,17 +40,31 @@ mkdir build cd build # Download the custom compiled qgis version tar of h5py, create a tar from the distro subfolder +# Download h5py 3.8.0 for QGis versions before 3.40 wget http://download.osgeo.org/osgeo4w/v2/x86_64/release/python3/python3-h5py/python3-h5py-3.8.0-1.tar.bz2 tar -xvf python3-h5py-3.8.0-1.tar.bz2 tar -cf h5py-3.8.0.tar -C ./apps/Python39/Lib/site-packages/ . cp h5py-3.8.0.tar .. +# Download h5py 3.10.0 for QGis versions after 3.40 +wget http://download.osgeo.org/osgeo4w/v2/x86_64/release/python3/python3-h5py/python3-h5py-3.10.0-1.tar.bz2 +tar -xvf python3-h5py-3.10.0-1.tar.bz2 +tar -cf h5py-3.10.0.tar -C ./apps/Python312/Lib/site-packages/ . +cp h5py-3.10.0.tar .. + # as well as scipy +# Download scipy 1.6.2 for QGis versions before 3.40 wget http://download.osgeo.org/osgeo4w/v2/x86_64/release/python3/python3-scipy/python3-scipy-1.10.1-1.tar.bz2 tar -xvf python3-scipy-1.10.1-1.tar.bz2 tar -cf scipy-1.10.1.tar -C ./apps/Python39/Lib/site-packages/ . cp scipy-1.10.1.tar .. +# Download scipy 1.13.0 for QGis versions after 3.40 +wget http://download.osgeo.org/osgeo4w/v2/x86_64/release/python3/python3-scipy/python3-scipy-1.13.0-1.tar.bz2 +tar -xvf python3-scipy-1.13.0-1.tar.bz2 +tar -cf scipy-1.13.0.tar -C ./apps/Python312/Lib/site-packages/ . +cp scipy-1.13.0.tar .. + # Back up a level and clean up the build/ directory. cd .. rm -rf build @@ -64,15 +78,16 @@ cp h5py/h5py-2.10.0-cp39-cp39-win_amd64.whl . # Copy pure wheels to prevent pip in docker (or Windows) to select platform dependent version wget https://files.pythonhosted.org/packages/cd/84/66072ee12c3e79061f183c09a24be24f45bb1286600589640363d9d416b0/SQLAlchemy-2.0.6-py3-none-any.whl#sha256=c5d754665edea1ecdc79e3023659cb5594372e10776f3b3734d75c2c3ce95013 -# Download windows wheels (cp39, win, amd64) +# Download windows wheels (cp39, cp312, win, amd64) wget https://files.pythonhosted.org/packages/b2/8e/83d9e3bff5c0ff7a0ec7e850c785916e616ab20d8793943f9e1d2a987fab/shapely-2.0.0-cp39-cp39-win_amd64.whl wget https://files.pythonhosted.org/packages/fe/e3/15a630b4cfd787bba84fb43beac7a25ae29c2b417ad7c4cfadae129a0819/threedigrid_builder-1.17.1-cp39-cp39-win_amd64.whl#sha256=24e06c136e399e1d6299c0a9b3e4c869084aba5878c833487fb13390aec10af8 wget https://files.pythonhosted.org/packages/b3/89/1d3b78577a6b2762cb254f6ce5faec9b7c7b23052d1cdb7237273ff37d10/greenlet-2.0.2-cp39-cp39-win_amd64.whl#sha256=db1a39669102a1d8d12b57de2bb7e2ec9066a6f2b3da35ae511ff93b01b5d564 wget https://files.pythonhosted.org/packages/5f/d6/5f59a5e5570c4414d94c6da4c97731deab832cbd14eaf23189d54a92d1e1/cftime-1.6.2-cp39-cp39-win_amd64.whl#sha256=86fe550b94525c327578a90b2e13418ca5ba6c636d5efe3edec310e631757eea +wget https://files.pythonhosted.org/packages/17/98/ba5b4a2f37c6c88454b696dd5c7a4e76fc8bfd014364b47ddd7e2cec0fcd/cftime-1.6.4-cp312-cp312-win_amd64.whl#sha256=5b5ad7559a16bedadb66af8e417b6805f758acb57aa38d2730844dfc63a1e667 - -# Download linux wheels (both cp38 and cp310) +# Download linux wheels (cp38, cp310, cp312) +wget https://files.pythonhosted.org/packages/d5/7d/9a57e187cbf2fbbbdfd4044a4f9ce141c8d221f9963750d3b001f0ec080d/shapely-2.0.6-cp312-cp312-manylinux_2_17_x86_64.manylinux2014_x86_64.whl wget https://files.pythonhosted.org/packages/06/07/0700e5e33c44bc87e19953244c29f73669cfb6f19868899170f9c7e34554/shapely-2.0.0-cp310-cp310-manylinux_2_17_x86_64.manylinux2014_x86_64.whl wget https://files.pythonhosted.org/packages/4e/03/f3bcb7d96aef6d56b62e2f25996f161c05f92a45d452165be2007b756e0f/shapely-2.0.0-cp38-cp38-manylinux_2_17_x86_64.manylinux2014_x86_64.whl wget https://files.pythonhosted.org/packages/53/67/a3f7e0ec8936a2882c97eb1e595ea91c3a5de4869d023064bbc8886d6ab7/threedigrid_builder-1.17.1-cp38-cp38-manylinux_2_17_x86_64.manylinux2014_x86_64.whl @@ -80,5 +95,6 @@ wget https://files.pythonhosted.org/packages/d0/6f/19252e21244e19adf6816cbc6b9d5 wget https://files.pythonhosted.org/packages/6e/11/a1f1af20b6a1a8069bc75012569d030acb89fd7ef70f888b6af2f85accc6/greenlet-2.0.2-cp310-cp310-manylinux_2_17_x86_64.manylinux2014_x86_64.whl#sha256=d75209eed723105f9596807495d58d10b3470fa6732dd6756595e89925ce2470 wget https://files.pythonhosted.org/packages/e1/17/d8042d82f44c08549b535bf2e7d1e87aa1863df5ed6cf1cf773eb2dfdf67/cftime-1.6.2-cp310-cp310-manylinux_2_17_x86_64.manylinux2014_x86_64.whl#sha256=acb294fdb80e33545ae54b4421df35c4e578708a5ffce1c00408b2294e70ecef wget https://files.pythonhosted.org/packages/44/51/bc9d47beee47afda1d335f05efa848dc403bd183344f03d431281518e8ab/cftime-1.5.0-cp38-cp38-manylinux_2_5_x86_64.manylinux1_x86_64.whl#sha256=7a820e16357dbdc9723b2059f7178451de626a8b2e5f80b9d91a77e3dac42133 +wget https://files.pythonhosted.org/packages/04/56/233d817ef571d778281f3d639049b342f6ff0bb4de4c5ee630befbd55319/cftime-1.6.4-cp312-cp312-manylinux_2_17_x86_64.manylinux2014_x86_64.whl#sha256=f92f2e405eeda47b30ab6231d8b7d136a55f21034d394f93ade322d356948654 touch .generated.marker diff --git a/gui/threedi_plugin_dockwidget.py b/gui/threedi_plugin_dockwidget.py index de7c5e96..371ea490 100644 --- a/gui/threedi_plugin_dockwidget.py +++ b/gui/threedi_plugin_dockwidget.py @@ -1,20 +1,25 @@ from pathlib import Path - -from threedi_results_analysis.threedi_plugin_model import ThreeDiGridItem, ThreeDiResultItem -from threedi_results_analysis.utils.constants import TOOLBOX_QGIS_SETTINGS_GROUP -from threedi_results_analysis.gui.threedi_plugin_grid_result_dialog import ThreeDiPluginGridResultDialog +from qgis.core import QgsSettings from qgis.PyQt import uic +from qgis.PyQt.QtCore import pyqtSignal +from qgis.PyQt.QtCore import pyqtSlot +from qgis.PyQt.QtCore import QModelIndex from qgis.PyQt.QtCore import Qt +from qgis.PyQt.QtGui import QPixmap +from qgis.PyQt.QtWidgets import QAction from qgis.PyQt.QtWidgets import QDockWidget -from qgis.PyQt.QtCore import QModelIndex -from qgis.PyQt.QtCore import pyqtSignal, pyqtSlot from qgis.PyQt.QtWidgets import QMenu -from qgis.PyQt.QtWidgets import QAction -from qgis.core import QgsSettings -from qgis.PyQt.QtGui import QPixmap from threedi_results_analysis import PLUGIN_DIR +from threedi_results_analysis.gui.threedi_plugin_grid_result_dialog import ( + ThreeDiPluginGridResultDialog, +) +from threedi_results_analysis.threedi_plugin_model import ThreeDiGridItem +from threedi_results_analysis.threedi_plugin_model import ThreeDiResultItem +from threedi_results_analysis.utils.constants import TOOLBOX_QGIS_SETTINGS_GROUP + import logging + logger = logging.getLogger(__name__) FORM_CLASS, _ = uic.loadUiType( @@ -65,15 +70,57 @@ def __init__(self, parent, iface): self.dialog.grid_file_selected.connect(self.grid_file_selected) self.dialog.result_grid_file_selected.connect(self.result_grid_file_selected) + self.custom_actions = {} + + def add_custom_actions(self, actions): + self.custom_actions |= actions + def customMenuRequested(self, pos): index = self.treeView.indexAt(pos) menu = QMenu(self) - action_delete = QAction("Delete", self) + action_remove = QAction("Remove", self) + action_remove.triggered.connect(lambda _, sel_index=index: self._remove_current_index_clicked(sel_index)) + menu.addAction(action_remove) + + for custom_action in self.custom_actions: + if custom_action.isSeparator(): + menu.addSeparator() + continue + else: + menu.addAction(custom_action) + custom_action.triggered.disconnect() + custom_action.triggered.connect(lambda _, sel_index=index: self._current_index_clicked(sel_index)) - action_delete.triggered.connect(lambda _, sel_index=index: self._remove_current_index_clicked(sel_index)) - menu.addAction(action_delete) menu.popup(self.treeView.viewport().mapToGlobal(pos)) + @pyqtSlot(QModelIndex) + def _current_index_clicked(self, index=None): + # note that index is the "current", not the "selected" + if not index: + index = self.treeView.selectionModel().currentIndex() + if index is not None and index.isValid(): + item = self.treeView.model().itemFromIndex(index) + action = self.sender() + if isinstance(item, ThreeDiGridItem): + self.custom_actions[action][0](item) + elif isinstance(item, ThreeDiResultItem): + self.custom_actions[action][1](item) + else: + raise RuntimeError("Unknown model item type") + + def _remove_current_index_clicked(self, index=None): + # note that index is the "current", not the "selected" + if not index: + index = self.treeView.selectionModel().currentIndex() + if index is not None and index.isValid(): + item = self.treeView.model().itemFromIndex(index) + if isinstance(item, ThreeDiGridItem): + self.grid_removal_selected.emit(item) + elif isinstance(item, ThreeDiResultItem): + self.result_removal_selected.emit(item) + else: + raise RuntimeError("Unknown model item type") + def set_model(self, model): tree_view = self.treeView @@ -94,19 +141,6 @@ def _add_clicked(self): self.dialog.refresh() self.dialog.exec() - def _remove_current_index_clicked(self, index=None): - # note that index is the "current", not the "selected" - if not index: - index = self.treeView.selectionModel().currentIndex() - if index is not None and index.isValid(): - item = self.treeView.model().itemFromIndex(index) - if isinstance(item, ThreeDiGridItem): - self.grid_removal_selected.emit(item) - elif isinstance(item, ThreeDiResultItem): - self.result_removal_selected.emit(item) - else: - raise RuntimeError("Unknown model item type") - def _align_starts_clicked(self): self.align_starts_checked.emit(self.alignStartsCheckBox.isChecked()) diff --git a/gui/threedi_plugin_grid_result_dialog.py b/gui/threedi_plugin_grid_result_dialog.py index 2c465aa7..520cc611 100644 --- a/gui/threedi_plugin_grid_result_dialog.py +++ b/gui/threedi_plugin_grid_result_dialog.py @@ -1,15 +1,22 @@ from functools import wraps from pathlib import Path +from qgis.core import QgsSettings +from qgis.PyQt import QtWidgets +from qgis.PyQt import uic +from qgis.PyQt.QtCore import pyqtSignal +from qgis.PyQt.QtCore import pyqtSlot +from qgis.PyQt.QtCore import QModelIndex +from qgis.PyQt.QtCore import QSortFilterProxyModel +from qgis.PyQt.QtCore import Qt +from qgis.PyQt.QtGui import QStandardItem +from qgis.PyQt.QtGui import QStandardItemModel +from qgis.PyQt.QtWidgets import QAbstractItemView +from threedi_mi_utils import list_local_schematisations +from threedi_results_analysis.utils.constants import TOOLBOX_QGIS_SETTINGS_GROUP + import logging import os -from threedi_results_analysis.utils.constants import TOOLBOX_QGIS_SETTINGS_GROUP -from qgis.PyQt import QtWidgets, uic -from qgis.PyQt.QtCore import pyqtSignal, pyqtSlot, QModelIndex, Qt, QSortFilterProxyModel -from qgis.PyQt.QtWidgets import QAbstractItemView -from qgis.core import QgsSettings -from qgis.PyQt.QtGui import QStandardItemModel, QStandardItem -from threedi_mi_utils import list_local_schematisations logger = logging.getLogger(__name__) @@ -168,14 +175,16 @@ def _retrieve_selected_grid_folder(self, index: QModelIndex) -> str: @pyqtSlot() @disable_dialog def _add_grid_from_table(self) -> None: - grid_file = os.path.join(self._retrieve_selected_grid_folder(self.tableView.currentIndex()), "gridadmin.h5") + index = self.proxy_model.mapToSource(self.tableView.currentIndex()) + grid_file = os.path.join(self._retrieve_selected_grid_folder(index), "gridadmin.h5") self.grid_file_selected.emit(grid_file) @pyqtSlot() @disable_dialog def _add_result_from_table(self) -> None: - result_file = os.path.join(self._retrieve_selected_result_folder(self.tableView.currentIndex()), "results_3di.nc") - grid_file = os.path.join(self._retrieve_selected_grid_folder(self.tableView.currentIndex()), "gridadmin.h5") + index = self.proxy_model.mapToSource(self.tableView.currentIndex()) + result_file = os.path.join(self._retrieve_selected_result_folder(index), "results_3di.nc") + grid_file = os.path.join(self._retrieve_selected_grid_folder(index), "gridadmin.h5") # Also emit corresponding grid file, if existst if not os.path.isfile(grid_file): @@ -184,6 +193,7 @@ def _add_result_from_table(self) -> None: self.result_grid_file_selected.emit(result_file, grid_file) def _item_selected(self, index: QModelIndex): + index = self.proxy_model.mapToSource(index) self.loadGridPushButton.setEnabled(True) # Only activate result button when revision contain results if self.model.item(index.row(), 2): @@ -193,6 +203,8 @@ def _item_selected(self, index: QModelIndex): @disable_dialog def _item_double_clicked(self, index: QModelIndex): + index = self.proxy_model.mapToSource(index) + # The selection contains a result if self.model.item(index.row(), 2): result_file = os.path.join(self._retrieve_selected_result_folder(index), "results_3di.nc") diff --git a/icons/icon_custom_statistics.png b/icons/icon_custom_statistics.png index 517714fd..ad1d1450 100644 Binary files a/icons/icon_custom_statistics.png and b/icons/icon_custom_statistics.png differ diff --git a/icons/icon_summary.png b/icons/icon_summary.png new file mode 100644 index 00000000..daf949a3 Binary files /dev/null and b/icons/icon_summary.png differ diff --git a/icons/sliders.svg b/icons/sliders.svg new file mode 100644 index 00000000..6df245fe --- /dev/null +++ b/icons/sliders.svg @@ -0,0 +1,18 @@ + diff --git a/metadata.txt b/metadata.txt index 7ad6a6ff..dfad14b4 100644 --- a/metadata.txt +++ b/metadata.txt @@ -2,7 +2,7 @@ name=3Di Results Analysis qgisMinimumVersion=3.22 description=3Di Results Analysis -version=3.9.4 +version=3.12 author=3Di Water Management email=servicedesk@nelen-schuurmans.nl @@ -24,8 +24,8 @@ about=Analyse 3Di results and visualize computational grids in the 3Di Modeller For questions on 3Di or the use of this plugin please contact us via servicedesk@nelen-schuurmans.nl -tracker=https://github.com/nens/ThreeDiToolbox -repository=https://github.com/nens/ThreeDiToolbox +tracker=https://github.com/nens/threedi-results-analysis +repository=https://github.com/nens/threedi-results-analysis changelog=https://docs.3di.live/a_releasenotes_3di_mi.html diff --git a/processing/cross_sectional_discharge_algorithm.py b/processing/cross_sectional_discharge_algorithm.py index 7838e950..f508ad8e 100644 --- a/processing/cross_sectional_discharge_algorithm.py +++ b/processing/cross_sectional_discharge_algorithm.py @@ -33,6 +33,7 @@ from osgeo import ogr from qgis.core import QgsCoordinateReferenceSystem +from qgis.core import QgsCoordinateTransform from qgis.core import QgsFeatureSink from qgis.core import QgsField from qgis.core import QgsFields @@ -40,6 +41,7 @@ from qgis.core import QgsProcessing from qgis.core import QgsProcessingAlgorithm from qgis.core import QgsProcessingContext +from qgis.core import QgsProcessingParameterBoolean from qgis.core import QgsProcessingParameterNumber from qgis.core import QgsProcessingParameterEnum from qgis.core import QgsProcessingParameterFeatureSink @@ -55,7 +57,7 @@ from qgis.core import QgsGeometry from qgis.core import QgsFeature from qgis.PyQt.QtCore import QCoreApplication, QVariant -from shapely import wkt +from shapely import wkb from threedigrid.admin.gridresultadmin import GridH5ResultAdmin from threedigrid.admin.constants import TYPE_V2_CHANNEL from threedigrid.admin.constants import TYPE_V2_CULVERT @@ -108,6 +110,7 @@ class CrossSectionalDischargeAlgorithm(QgsProcessingAlgorithm): GRIDADMIN_GPKG = "GRIDADMIN_GPKG" RESULTS_3DI_INPUT = "RESULTS_3DI_INPUT" CROSS_SECTION_LINES_INPUT = "CROSS_SECTION_LINES_INPUT" + SELECTED_CROSS_SECTION_LINES = "SELECTED_CROSS_SECTION_LINES" START_TIME = "START_TIME" END_TIME = "END_TIME" SUBSET = "SUBSET" @@ -162,6 +165,13 @@ def initAlgorithm(self, config): ) ) + self.addParameter( + QgsProcessingParameterBoolean( + self.SELECTED_CROSS_SECTION_LINES, + self.tr("Selected features only"), + ) + ) + self.addParameter( QgsProcessingParameterNumber( self.START_TIME, @@ -239,6 +249,7 @@ def processAlgorithm(self, parameters, context, feedback): parameters, self.CROSS_SECTION_LINES_INPUT, context ) self.cross_section_lines_id = cross_section_lines.id() + selected_cross_section_lines_only = self.parameterAsBool(parameters, self.SELECTED_CROSS_SECTION_LINES, context) start_time = ( self.parameterAsInt(parameters, self.START_TIME, context) if parameters[self.START_TIME] is not None @@ -287,14 +298,17 @@ def processAlgorithm(self, parameters, context, feedback): QgsField(name="q_net_sum", type=QVariant.Double) ) - crs = QgsCoordinateReferenceSystem(f"EPSG:{gr.epsg_code}") + threedi_results_crs = QgsCoordinateReferenceSystem(f"EPSG:{gr.epsg_code}") + cross_section_lines_crs = cross_section_lines.sourceCrs() + coordinate_transform = QgsCoordinateTransform(cross_section_lines_crs, threedi_results_crs, context.project()) + (flowlines_sink, self.flowlines_sink_dest_id) = self.parameterAsSink( parameters, self.OUTPUT_FLOWLINES, context, fields=flowlines_sink_fields, geometryType=QgsWkbTypes.LineString, - crs=crs, + crs=threedi_results_crs, ) self.target_field_idx = ( @@ -314,7 +328,8 @@ def processAlgorithm(self, parameters, context, feedback): feedback.setProgress(0) self.total_discharges = dict() - if cross_section_lines.selectedFeatureCount() > 0: + + if selected_cross_section_lines_only: iterator = cross_section_lines.getSelectedFeatures() nr_features = cross_section_lines.selectedFeatureCount() else: @@ -326,7 +341,9 @@ def processAlgorithm(self, parameters, context, feedback): feedback.setProgressText( f"Processing cross-section line {gauge_line.id()}..." ) - shapely_linestring = wkt.loads(gauge_line.geometry().asWkt()) + transformed_geometry = gauge_line.geometry() + transformed_geometry.transform(coordinate_transform) + shapely_linestring = wkb.loads(bytes(transformed_geometry.asWkb())) tgt_ds = MEMORY_DRIVER.CreateDataSource("") ts_gauge_line, total_discharge = left_to_right_discharge_ogr( gr=gr, @@ -368,24 +385,25 @@ def processAlgorithm(self, parameters, context, feedback): ) feedback.setProgress(100 * i / nr_features) - np.savetxt( - self.csv_output_file_path, - ts_all_cross_section_lines, - delimiter=",", - header=",".join(column_names), - fmt=formatting, - comments="", - ) - layer = QgsVectorLayer(self.csv_output_file_path, "Time series output") - context.temporaryLayerStore().addMapLayer(layer) - layer_details = QgsProcessingContext.LayerDetails( - "Output: Time series", context.project(), "Output: Time series" - ) - context.addLayerToLoadOnCompletion(layer.id(), layer_details) + if nr_features > 0: + np.savetxt( + self.csv_output_file_path, + ts_all_cross_section_lines, + delimiter=",", + header=",".join(column_names), + fmt=formatting, + comments="", + ) + layer = QgsVectorLayer(self.csv_output_file_path, "Time series output") + context.temporaryLayerStore().addMapLayer(layer) + layer_details = QgsProcessingContext.LayerDetails( + "Output: Time series", context.project(), "Output: Time series" + ) + context.addLayerToLoadOnCompletion(layer.id(), layer_details) return { self.OUTPUT_FLOWLINES: self.flowlines_sink_dest_id, - self.OUTPUT_TIME_SERIES: self.csv_output_file_path, + self.OUTPUT_TIME_SERIES: self.csv_output_file_path if nr_features > 0 else None, } def postProcessAlgorithm(self, context, feedback): diff --git a/processing/deps/__init__.py b/processing/deps/__init__.py new file mode 100644 index 00000000..e69de29b diff --git a/processing/deps/discharge/cross_sectional_discharge.py b/processing/deps/discharge/cross_sectional_discharge.py index df657499..7fce6b84 100644 --- a/processing/deps/discharge/cross_sectional_discharge.py +++ b/processing/deps/discharge/cross_sectional_discharge.py @@ -154,7 +154,7 @@ def left_to_right_discharge( content_type__in=content_types ) ts, tintervals = prepare_timeseries( - nodes_or_lines=intersecting_lines, + threedigrid_object=intersecting_lines, start_time=start_time, end_time=end_time, aggregation=Q_NET_SUM, diff --git a/processing/deps/discharge/discharge_reduction.py b/processing/deps/discharge/discharge_reduction.py index 0b15b794..bf9cf58a 100644 --- a/processing/deps/discharge/discharge_reduction.py +++ b/processing/deps/discharge/discharge_reduction.py @@ -67,7 +67,7 @@ def __init__( feedback.setProgressText("Calculate cumulative discharges...") all_2d_open_water_flowlines = grid_result_admin.lines.subset('2D_OPEN_WATER').filter(id__in=flowline_ids) discharges, self.tintervals = prepare_timeseries( - nodes_or_lines=all_2d_open_water_flowlines, + threedigrid_object=all_2d_open_water_flowlines, aggregation=self.Q_NET_SUM ) q_net_sum = aggregate_prepared_timeseries( diff --git a/processing/deps/rasters_to_netcdf/__init__.py b/processing/deps/rasters_to_netcdf/__init__.py new file mode 100644 index 00000000..e69de29b diff --git a/processing/deps/rasters_to_netcdf/rain.tif b/processing/deps/rasters_to_netcdf/rain.tif new file mode 100644 index 00000000..8555821e Binary files /dev/null and b/processing/deps/rasters_to_netcdf/rain.tif differ diff --git a/processing/deps/rasters_to_netcdf/rasters_to_netcdf.py b/processing/deps/rasters_to_netcdf/rasters_to_netcdf.py new file mode 100644 index 00000000..90513b15 --- /dev/null +++ b/processing/deps/rasters_to_netcdf/rasters_to_netcdf.py @@ -0,0 +1,190 @@ +from datetime import datetime, timedelta +from pathlib import Path +from typing import List, Union, Dict + +from cftime import date2num +import numpy as np +# import netCDF4 +from osgeo import gdal, osr # IMPORT THIS BEFORE IMPORTING h5netcdf!!! +import h5netcdf.legacyapi as netCDF4 +from h5netcdf.legacyapi import Variable as h5netcdf_Variable +from pyproj import CRS + +gdal.UseExceptions() +osr.UseExceptions() + + +def get_datasets(filepaths: List[Union[str, Path]]) -> List[gdal.Dataset]: + result = list() + for filepath in filepaths: + dataset = gdal.Open(str(filepath)) + if not dataset: + raise Exception(f"Unable to open {filepath}") + result.append(dataset) + return result + + +def get_srs(dataset: gdal.Dataset) -> osr.SpatialReference: + projection = dataset.GetProjection() + srs = osr.SpatialReference() + srs.ImportFromWkt(projection) + return srs + + +def rasters_have_same_srs(datasets: List[gdal.Dataset]) -> bool: + if not datasets: + raise ValueError("The filepaths list is empty.") + + base_crs = get_srs(datasets[0]) + + for dataset in datasets[1:]: + current_crs = get_srs(dataset) + if not base_crs.IsSame(current_crs): + return False + return True + + +def rasters_have_same_geotransform(datasets: List[gdal.Dataset]) -> bool: + if not datasets: + raise ValueError("The filepaths list is empty.") + + base_geotransform = datasets[0].GetGeoTransform() + + for dataset in datasets[1:]: + current_geotransform = dataset.GetGeoTransform() + if base_geotransform != current_geotransform: + return False + return True + + +def rasters_have_same_dimensions(datasets: List[gdal.Dataset]): + if not datasets: + raise ValueError("The filepaths list is empty.") + + base_dimensions = datasets[0].RasterXSize, datasets[0].RasterYSize + + for dataset in datasets[1:]: + current_dimensions = dataset.RasterXSize, dataset.RasterYSize + if base_dimensions != current_dimensions: + return False + return True + + +def rasters_have_same_nodatavalue(datasets: List[gdal.Dataset]): + if not datasets: + raise ValueError("The filepaths list is empty.") + + base_ndv = datasets[0].GetRasterBand(1).GetNoDataValue() + + for dataset in datasets[1:]: + current_ndv = dataset.GetRasterBand(1).GetNoDataValue() + if base_ndv != current_ndv: + return False + return True + + +def setncatts(variable: h5netcdf_Variable, attributes: Dict): + """Reproduces the behaviour of netcdf4.Variable.setncatts(), which is not supported in h5netcdf lecagy API""" + for key, value in attributes.items(): + variable.setncattr(key, value) + + +def rasters_to_netcdf( + rasters: List[Union[str, Path]], + start_time: datetime, + interval: int, + units: str, + output_path: Union[str, Path], + time_units: str = 'seconds since 1970-01-01 00:00:00.0 +0000', + calendar: str = 'standard', + offset: int = 0 +) -> None: + """ + :param interval: interval in seconds + :param units: one of 'mm', 'm/s', 'mm/h', 'mm/hr'. Note: in case of `mm` the rate is determined by looking at the + next `time` value. + :param offset: offset in seconds + """ + datasets = get_datasets(rasters) + assert rasters_have_same_srs(datasets), "Not all input rasters have the same Spatial Reference System" + assert rasters_have_same_geotransform(datasets), "Not all input rasters have the same origin, pixel size, and skew" + assert rasters_have_same_dimensions(datasets), "Not all input rasters have the same width and height" + assert rasters_have_same_nodatavalue(datasets), "Not all input rasters have the same nodatavalue" + + srs = get_srs(datasets[0]) + crs = CRS.from_wkt(srs.ExportToWkt()) + geotransform = datasets[0].GetGeoTransform() + + # create netcdf + output_dataset = netCDF4.Dataset( + str(output_path), + mode='w', + # format="NETCDF4" + ) + + # dataset attributes + output_dataset.setncattr(name='OFFSET', value=offset) + crs_var = output_dataset.createVariable(varname='crs', datatype='int') + setncatts(crs_var, crs.to_cf()) + crs_var.setncattr(name="spatial_ref", value=crs.to_wkt()) + crs_var.setncattr(name="GeoTransform", value=" ".join([str(i) for i in geotransform])) + + # set dimensions + output_dataset.createDimension('lon', size=datasets[0].RasterXSize) + output_dataset.createDimension('lat', size=datasets[0].RasterYSize) + output_dataset.createDimension('time', size=len(rasters)) + + # x and y or lon and lat + x_attrs, y_attrs = crs.cs_to_cf() + + max_x_ordinate = geotransform[0] + datasets[0].RasterXSize * geotransform[1] + x_ordinates = np.arange( + start=geotransform[0], + stop=max_x_ordinate - 0.5 * geotransform[1], # to prevent unexpected behaviour due to rounding differences + step=geotransform[1] + ) + lon_var = output_dataset.createVariable(varname='lon', datatype='float32', dimensions=('lon',)) + lon_var[:] = x_ordinates + setncatts(lon_var, x_attrs) + + max_y_ordinate = geotransform[3] + datasets[0].RasterYSize * geotransform[5] + y_ordinates = np.arange( + start=geotransform[3], + stop=max_y_ordinate - 0.5 * geotransform[5], # to prevent unexpected behaviour due to rounding differences + step=geotransform[5] + ) + lat_var = output_dataset.createVariable(varname='lat', datatype='float32', dimensions=('lat',)) + lat_var[:] = y_ordinates + setncatts(lat_var, y_attrs) + + # time + time_var = output_dataset.createVariable(varname='time', datatype='float64', dimensions=('time',)) + time_attrs = { + 'standard_name': 'time', + 'long_name': 'Time', + 'units': time_units, + 'calendar': calendar, + 'axis': 'T' + } + setncatts(time_var, time_attrs) + time_delta = timedelta(seconds=interval) + end_time = start_time + len(rasters) * time_delta + time_steps_numpy = np.arange(start_time, end_time, time_delta, dtype='datetime64[s]') + time_steps_datetime = [datetime.utcfromtimestamp(dt.astype('int')) for dt in time_steps_numpy] + time_var[:] = date2num(time_steps_datetime, units=time_units, calendar=calendar) + + # rain data + rain_var = output_dataset.createVariable(varname='values', datatype='float', dimensions=('time', 'lat', 'lon')) + rain_attrs = { + 'long_name': 'rain', + 'grid_mapping': 'crs', + '_FillValue': datasets[0].GetRasterBand(1).GetNoDataValue(), + 'missing_value': datasets[0].GetRasterBand(1).GetNoDataValue(), + 'units': units + } + setncatts(rain_var, rain_attrs) + for i, dataset in enumerate(datasets): + rain_var[i] = dataset.GetRasterBand(1).ReadAsArray() + output_dataset.sync() + output_dataset.close() + diff --git a/processing/deps/rasters_to_netcdf/test_rasters_to_netcdf.py b/processing/deps/rasters_to_netcdf/test_rasters_to_netcdf.py new file mode 100644 index 00000000..741832c1 --- /dev/null +++ b/processing/deps/rasters_to_netcdf/test_rasters_to_netcdf.py @@ -0,0 +1,30 @@ +from datetime import datetime +from pathlib import Path +import tempfile + +import pytest + +from threedi_results_analysis.processing.deps.rasters_to_netcdf.rasters_to_netcdf import rasters_to_netcdf +from osgeo import gdal +gdal.UseExceptions() + +DATA_DIR = Path(__file__).parent + + +def test_rasters_to_netcdf(): + filepaths = [ + DATA_DIR / "rain.tif", + DATA_DIR / "rain.tif", + DATA_DIR / "rain.tif" + ] + with tempfile.TemporaryDirectory() as tmpdir: + output_path = Path(tmpdir) / "output4.nc" + rasters_to_netcdf( + rasters=filepaths, + start_time=datetime.strptime('2020-01-01T12:00:00', "%Y-%m-%dT%H:%M:%S"), + interval=3600, + units="mm", + output_path=output_path + ) + result_nc = gdal.Open(str(output_path)) + assert result_nc.RasterCount == 3 diff --git a/processing/dwf_calculation_algorithm.py b/processing/dwf_calculation_algorithm.py index 4a4443c5..0eec410c 100644 --- a/processing/dwf_calculation_algorithm.py +++ b/processing/dwf_calculation_algorithm.py @@ -295,21 +295,10 @@ def displayName(self): return self.tr("DWF Calculator") def group(self): - """ - Returns the name of the group this algorithm belongs to. This string - should be localised. - """ - return self.tr("Dry weather flow") + return "Pre-process simulation inputs" def groupId(self): - """ - Returns the unique ID of the group this algorithm belongs to. This - string should be fixed for the algorithm, and must not be localised. - The group id should be unique within each provider. Group id should - contain lowercase alphanumeric characters only and no spaces or other - formatting characters. - """ - return "dwf" + return "pre_process_sim_inputs" def shortHelpString(self): diff --git a/processing/providers.py b/processing/providers.py index 56f80d4d..6748bb5c 100644 --- a/processing/providers.py +++ b/processing/providers.py @@ -7,8 +7,13 @@ from threedi_results_analysis.processing.cross_sectional_discharge_algorithm import CrossSectionalDischargeAlgorithm from threedi_results_analysis.processing.leak_detector_algorithms import ( DetectLeakingObstaclesAlgorithm, +) +from threedi_results_analysis.processing.leak_detector_algorithms import ( DetectLeakingObstaclesWithDischargeThresholdAlgorithm, ) +from threedi_results_analysis.processing.rasters_to_netcdf_algorithm import ( + RastersToNetCDFAlgorithm, +) from threedi_results_analysis.processing.schematisation_algorithms import ( CheckSchematisationAlgorithm, MigrateAlgorithm, @@ -16,8 +21,21 @@ # GuessIndicatorAlgorithm, # ImportHydXAlgorithm, ) - -from threedi_results_analysis.processing.threedidepth_algorithms import ThreediDepthAlgorithm, ThreediMaxDepthAlgorithm +from threedi_results_analysis.processing.schematisation_algorithms import ( + ImportSufHydAlgorithm, +) +from threedi_results_analysis.processing.schematisation_algorithms import ( + MigrateAlgorithm, +) +from threedi_results_analysis.processing.structure_control_action_algorithms import ( + StructureControlActionAlgorithm, +) +from threedi_results_analysis.processing.threedidepth_algorithms import ( + ThreediDepthAlgorithm, +) +from threedi_results_analysis.processing.threedidepth_algorithms import ( + ThreediMaxDepthAlgorithm, +) import os @@ -39,6 +57,8 @@ def loadAlgorithms(self, *args, **kwargs): self.addAlgorithm(CrossSectionalDischargeAlgorithm()) self.addAlgorithm(DetectLeakingObstaclesAlgorithm()) self.addAlgorithm(DetectLeakingObstaclesWithDischargeThresholdAlgorithm()) + self.addAlgorithm(RastersToNetCDFAlgorithm()) + self.addAlgorithm(StructureControlActionAlgorithm()) def id(self, *args, **kwargs): """The ID of your plugin, used for identifying the provider. diff --git a/processing/rasters_to_netcdf_algorithm.py b/processing/rasters_to_netcdf_algorithm.py new file mode 100644 index 00000000..788dccab --- /dev/null +++ b/processing/rasters_to_netcdf_algorithm.py @@ -0,0 +1,162 @@ +from typing import List + +from threedi_results_analysis.processing.deps.rasters_to_netcdf.rasters_to_netcdf import rasters_to_netcdf +from qgis.core import ( + QgsMapLayer, + QgsProcessing, + QgsProcessingAlgorithm, + QgsProcessingContext, + QgsProcessingParameterDateTime, + QgsProcessingParameterEnum, + QgsProcessingParameterFileDestination, + QgsProcessingParameterMultipleLayers, + QgsProcessingParameterNumber, + QgsRasterLayer, +) + + +class RastersToNetCDFAlgorithm(QgsProcessingAlgorithm): + """ + Processing algorithm to create a NetCDF file with data of rain or other forcings that vary in space and time. + """ + + INPUT_RASTERS = "INPUT_RASTERS" + INPUT_START_TIME = "INPUT_START_TIME" + INPUT_INTERVAL = "INPUT_INTERVAL" + INPUT_UNITS = "INPUT_UNITS" + INPUT_OUTPUT_PATH = "INPUT_OUTPUT_PATH" + INPUT_OFFSET = "INPUT_OFFSET" + + OUTPUT = "OUTPUT" + + def initAlgorithm(self, config): + self.addParameter( + QgsProcessingParameterMultipleLayers( + name=self.INPUT_RASTERS, + description="Input rasters", + layerType=QgsProcessing.TypeRaster, + ) + ) + + self.addParameter( + QgsProcessingParameterDateTime( + name=self.INPUT_START_TIME, + description="Start", + ) + ) + + self.addParameter( + QgsProcessingParameterNumber( + name=self.INPUT_INTERVAL, + description="Interval [s]", + type=QgsProcessingParameterNumber.Integer + ) + ) + + self.addParameter( + QgsProcessingParameterNumber( + name=self.INPUT_OFFSET, + description="Offset [s]", + type=QgsProcessingParameterNumber.Integer, + defaultValue=0 + ) + ) + + self.addParameter( + QgsProcessingParameterEnum( + name=self.INPUT_UNITS, + description="Units", + options=["mm", "m/s", "mm/h"], + usesStaticStrings=True + ) + ) + + self.addParameter( + QgsProcessingParameterFileDestination( + name=self.INPUT_OUTPUT_PATH, + description="Output file", + fileFilter="NetCDF (*.nc *.NC)" + ) + ) + + def processAlgorithm(self, parameters, context, feedback): + input_raster_layers: List[QgsMapLayer] = self.parameterAsLayerList(parameters, self.INPUT_RASTERS, context) + rasters = [layer.dataProvider().dataSourceUri() for layer in input_raster_layers] + start_datetime = self.parameterAsDateTime(parameters, self.INPUT_START_TIME, context) \ + .toPyDateTime() + if start_datetime.microsecond != 0: + feedback.pushInfo("Setting start datetime milliseconds to 0...") + start_datetime.replace(microsecond=0) + if start_datetime.second != 0: + feedback.pushInfo("Setting start datetime seconds to 0...") + start_datetime.replace(second=0) + interval = self.parameterAsInt(parameters, self.INPUT_INTERVAL, context) + offset = self.parameterAsInt(parameters, self.INPUT_OFFSET, context) + units = self.parameterAsEnumString(parameters, self.INPUT_UNITS, context) + output_path = self.parameterAsFile(parameters, self.INPUT_OUTPUT_PATH, context) + rasters_to_netcdf( + rasters=rasters, + start_time=start_datetime, + interval=interval, + units=units, + output_path=output_path, + offset=offset + ) + + layer = QgsRasterLayer(output_path, "Spatiotemporal NetCDF") + context.temporaryLayerStore().addMapLayer(layer) + layer_details = QgsProcessingContext.LayerDetails( + "Spatiotemporal NetCDF", + context.project(), + self.OUTPUT + ) + context.addLayerToLoadOnCompletion( + layer.id(), + layer_details + ) + + return { + self.OUTPUT: output_path, + } + + def group(self): + return "Pre-process simulation inputs" + + def groupId(self): + return "pre_process_sim_inputs" + + def shortHelpString(self): + return """ +
Create a NetCDF file with data of rain or other forcings that vary in space and time.
+The algorithm takes a list of rasters and stacks them into a NetCDF, one raster for each time step.
+ⓘ Note that 3Di also offers services that seamlessly integrate historical and forecast rain with 3Di. For example, to set up flood early warning systems or operational water management systems. Get in touch via www.3diwatermanagement.com to learn the possibilities for your area.
+A list of rasters (e.g. GeoTIFF) to be stacked.
+Date and time of the first time step in the output NetCDF. Seconds and milliseconds are ignored (set to 0).
+Time in seconds between time steps.
+If greater than 0, the forcing will only be applied after offset seconds have passed in the simulation.
+The units of the forcing's data. Choose 'mm' to indicate that the values are total amounts per time interval.
+Name and location of the NetCDF output. Must have the .nc extension.
+ """ + + def name(self): + """ + Returns the algorithm name, used for identifying the algorithm + """ + return "rasters_to_netcdf" + + def displayName(self): + """ + Returns the algorithm name, which should be used for any + user-visible display of the algorithm name. + """ + return "Rasters to spatiotemporal NetCDF" + + def createInstance(self): + return RastersToNetCDFAlgorithm() diff --git a/processing/structure_control_action_algorithms.py b/processing/structure_control_action_algorithms.py new file mode 100644 index 00000000..d06d5fd4 --- /dev/null +++ b/processing/structure_control_action_algorithms.py @@ -0,0 +1,116 @@ +# -*- coding: utf-8 -*- + +""" +*************************************************************************** +* * +* This program is free software; you can redistribute it and/or modify * +* it under the terms of the GNU General Public License as published by * +* the Free Software Foundation; either version 2 of the License, or * +* (at your option) any later version. * +* * +*************************************************************************** +""" +from qgis.core import QgsProcessingAlgorithm +from qgis.core import QgsProcessingException +from qgis.core import QgsProcessingParameterBoolean +from qgis.core import QgsProcessingParameterFile +from qgis.core import QgsProcessingParameterFileDestination +from qgis.core import QgsProject +from qgis.core import QgsVectorLayer +from qgis.PyQt.QtCore import QCoreApplication +from threedigrid.admin.gridresultadmin import GridH5StructureControl +from threedigrid.admin.structure_controls.exporters import ( + structure_control_actions_to_csv, +) + +import os + + +class StructureControlActionAlgorithm(QgsProcessingAlgorithm): + """ + Converts a structure control actions NetCDF to CSV + """ + + SCA_INPUT = "SCA_INPUT" + OUTPUT_FILENAME = "OUTPUT_FILENAME" + GRIDADMIN_INPUT = "GRIDADMIN_INPUT" + ADD_TO_PROJECT = "ADD_TO_PROJECT" + + def tr(self, string): + return QCoreApplication.translate("Processing", string) + + def createInstance(self): + return StructureControlActionAlgorithm() + + def name(self): + return "structure_control_actions_csv" + + def displayName(self): + return self.tr("Convert structure control actions") + + def group(self): + return self.tr("Post-process results") + + def groupId(self): + return "postprocessing" + + def shortHelpString(self): + return self.tr("Convert a structure control actions NetCDF to CSV") + + def initAlgorithm(self, config=None): + # Input parameters + self.addParameter( + QgsProcessingParameterFile(self.GRIDADMIN_INPUT, self.tr("Gridadmin.h5 file"), extension="h5") + ) + self.addParameter( + QgsProcessingParameterFile(self.SCA_INPUT, self.tr("structure_control_actions_3di.nc file"), extension="nc") + ) + self.addParameter( + QgsProcessingParameterBoolean( + self.ADD_TO_PROJECT, self.tr("Add result to project"), defaultValue=True + ) + ) + # output parameters + self.addParameter( + QgsProcessingParameterFileDestination( + self.OUTPUT_FILENAME, + self.tr("Destination CSV file path"), + fileFilter="csv", + ) + ) + + def processAlgorithm(self, parameters, context, feedback): + """ + Create the water depth raster with the provided user inputs + """ + gridadmin_path = parameters[self.GRIDADMIN_INPUT] + results_3di_path = parameters[self.SCA_INPUT] + generated_output_file_path = self.parameterAsFileOutput( + parameters, self.OUTPUT_FILENAME, context + ) + self.csv_output_file = f"{os.path.splitext(generated_output_file_path)[0]}.csv" + + self.add_to_project = self.parameterAsBoolean( + parameters, self.ADD_TO_PROJECT, context + ) + + if not self.csv_output_file: + raise QgsProcessingException(self.invalidSourceError(parameters, self.OUTPUT_FILENAME)) + + # https://threedigrid.readthedocs.io/en/latest/structure_control.html + gst = GridH5StructureControl(gridadmin_path, results_3di_path) + try: + structure_control_actions_to_csv(gst, self.csv_output_file) + except Exception as e: + return {"result": False, "error": str(e)} + + return {"result": True} + + def postProcessAlgorithm(self, context, feedback): + if self.add_to_project: + if self.csv_output_file: + result_layer = QgsVectorLayer( + self.csv_output_file, "Structure control actions" + ) + QgsProject.instance().addMapLayer(result_layer) + return {self.OUTPUT_FILENAME: self.csv_output_file} diff --git a/requirements-dev.txt b/requirements-dev.txt index cfe41755..d3cedc0e 100644 --- a/requirements-dev.txt +++ b/requirements-dev.txt @@ -10,6 +10,7 @@ docstr-coverage >= 1.0.4 flake8 isort mock +pip-tools pytest < 8.0 pytest-cov pytest-flake8 diff --git a/threedi_plugin.py b/threedi_plugin.py index d5feab45..1bf0ea66 100644 --- a/threedi_plugin.py +++ b/threedi_plugin.py @@ -1,30 +1,47 @@ +from qgis.core import QgsApplication +from qgis.core import QgsMapLayer +from qgis.core import QgsPathResolver +from qgis.core import QgsProject +from qgis.core import QgsSettings from qgis.PyQt.QtCore import QObject from qgis.PyQt.QtGui import QIcon from qgis.PyQt.QtWidgets import QAction -from qgis.PyQt.QtXml import QDomDocument, QDomElement +from qgis.PyQt.QtXml import QDomDocument +from qgis.PyQt.QtXml import QDomElement from qgis.utils import iface -from qgis.core import QgsApplication, QgsProject, QgsPathResolver, QgsSettings, QgsMapLayer +from threedi_results_analysis.gui.threedi_plugin_dockwidget import ( + ThreeDiPluginDockWidget, +) from threedi_results_analysis.misc_tools import About from threedi_results_analysis.misc_tools import ShowLogfile from threedi_results_analysis.misc_tools import ToggleResultsManager from threedi_results_analysis.processing.providers import ThreediProvider -from threedi_results_analysis.gui.threedi_plugin_dockwidget import ThreeDiPluginDockWidget -from threedi_results_analysis.threedi_plugin_layer_manager import ThreeDiPluginLayerManager -from threedi_results_analysis.threedi_plugin_model import ThreeDiPluginModel -from threedi_results_analysis.threedi_plugin_model_validation import ThreeDiPluginModelValidator from threedi_results_analysis.temporal import TemporalManager -from threedi_results_analysis.threedi_plugin_model_serialization import ThreeDiPluginModelSerializer +from threedi_results_analysis.threedi_plugin_layer_manager import ( + ThreeDiPluginLayerManager, +) +from threedi_results_analysis.threedi_plugin_model import ThreeDiPluginModel +from threedi_results_analysis.threedi_plugin_model_serialization import ( + ThreeDiPluginModelSerializer, +) +from threedi_results_analysis.threedi_plugin_model_validation import ( + ThreeDiPluginModelValidator, +) from threedi_results_analysis.tool_animation.map_animator import MapAnimator +from threedi_results_analysis.tool_flow_summary.flow_summary import FlowSummaryTool from threedi_results_analysis.tool_graph.graph import ThreeDiGraph from threedi_results_analysis.tool_sideview.sideview import ThreeDiSideView from threedi_results_analysis.tool_statistics.statistics import StatisticsTool from threedi_results_analysis.tool_water_balance import WaterBalanceTool -from threedi_results_analysis.tool_watershed.watershed_analysis import ThreeDiWatershedAnalyst +from threedi_results_analysis.tool_watershed.watershed_analysis import ( + ThreeDiWatershedAnalyst, +) from threedi_results_analysis.utils import color from threedi_results_analysis.utils.qprojects import ProjectStateMixin import logging + logger = logging.getLogger(__name__) @@ -65,6 +82,9 @@ def initGui(self): self.actions = [] self.menu = "&3Di Results Analysis" + assert not hasattr(self, "dockwidget") # Should be destroyed on unload + self.dockwidget = ThreeDiPluginDockWidget(None, iface) + # Set toolbar and init a few toolbar widgets self.toolbar = self.iface.addToolBar("ThreeDiResultAnalysis") self.toolbar.setObjectName("ThreeDiResultAnalysisToolBar") @@ -79,16 +99,18 @@ def initGui(self): self.watershed_tool = ThreeDiWatershedAnalyst(iface, self.model) self.logfile_tool = ShowLogfile(iface) self.temporal_manager = TemporalManager(self.model) + self.flow_summary_tool = FlowSummaryTool(self.dockwidget.get_tools_widget(), iface, self.model) self.tools = [ # second item indicates enabled on startup (self.about_tool, True), (self.toggle_results_manager, True), + (self.flow_summary_tool, False), (self.graph_tool, False), (self.sideview_tool, False), (self.stats_tool, False), (self.water_balance_tool, False), (self.watershed_tool, False), - (self.logfile_tool, True), + (self.logfile_tool, True) ] # Styling (TODO: can this be moved to where it is used?) @@ -96,17 +118,15 @@ def initGui(self): color.add_color_ramp(color_ramp) for tool, enabled in self.tools: - self._add_action( - tool, - tool.icon_path, - text=tool.menu_text, - callback=tool.run, - parent=self.iface.mainWindow(), - enabled_flag=enabled - ) - - assert not hasattr(self, "dockwidget") # Should be destroyed on unload - self.dockwidget = ThreeDiPluginDockWidget(None, iface) + if tool.icon_path is not None: + self._add_action( + tool, + tool.icon_path, + text=tool.menu_text, + callback=tool.run, + parent=self.iface.mainWindow(), + enabled_flag=enabled + ) # Connect the signals @@ -147,8 +167,14 @@ def initGui(self): self.model.result_checked.connect(self.map_animator.results_changed) self.model.result_unchecked.connect(self.map_animator.results_changed) self.model.result_added.connect(self.map_animator.results_changed) + self.model.result_removed.connect(self.map_animator.results_changed) self.temporal_manager.updated.connect(self.map_animator.update_results) + # flow summary signals + self.model.result_removed.connect(self.flow_summary_tool.result_removed) + self.model.result_changed.connect(self.flow_summary_tool.result_changed) + self.model.result_added.connect(self.flow_summary_tool.result_added) + # graph signals self.model.result_added.connect(self.graph_tool.result_added) self.model.result_removed.connect(self.graph_tool.result_removed) @@ -183,6 +209,9 @@ def initGui(self): self.model.result_changed.connect(self.water_balance_tool.result_changed) self.model.grid_changed.connect(self.water_balance_tool.grid_changed) + for tool, _ in self.tools: + self.dockwidget.add_custom_actions(tool.get_custom_actions()) + # Further administrative signals that need to happens last: # https://doc.qt.io/qt-5/signalsandslots.html#signals # If several slots are connected to one signal, the slots will be executed one after the other, @@ -209,6 +238,10 @@ def write(self, doc: QDomDocument) -> bool: if not tool.write(doc, node): return False + # Also allow animator to persist settings + if not self.map_animator.write(doc, node): + return False + return True def write_map_layer(self, layer: QgsMapLayer, elem: QDomElement, _: QDomDocument): @@ -229,9 +262,15 @@ def read(self, doc: QDomDocument) -> bool: self.model.clear() return False - # Allow each tool to read additional info from the dedicated xml node - for tool, _ in self.tools: - if not tool.read(tool_node): + if tool_node: + # Allow each tool to read additional info from the dedicated xml node + for tool, _ in self.tools: + if not tool.read(tool_node): + self.model.clear() + return False + + # Also allow animator to read saved settings + if not self.map_animator.read(tool_node): self.model.clear() return False diff --git a/threedi_plugin_model.py b/threedi_plugin_model.py index 445d9c09..7e58517e 100644 --- a/threedi_plugin_model.py +++ b/threedi_plugin_model.py @@ -1,14 +1,18 @@ -from pathlib import Path -from typing import List from functools import cached_property -from qgis.PyQt.QtCore import pyqtSignal, pyqtSlot -from qgis.PyQt.QtCore import QModelIndex, Qt -from qgis.PyQt.QtGui import QStandardItem, QStandardItemModel +from pathlib import Path +from qgis.PyQt.QtCore import pyqtSignal +from qgis.PyQt.QtCore import pyqtSlot +from qgis.PyQt.QtCore import QModelIndex +from qgis.PyQt.QtCore import Qt +from qgis.PyQt.QtGui import QStandardItem +from qgis.PyQt.QtGui import QStandardItemModel from threedi_results_analysis.datasource.threedi_results import ThreediResult +from typing import List import logging import uuid + logger = logging.getLogger(__name__) already_used_ids = [] diff --git a/threedi_plugin_tool.py b/threedi_plugin_tool.py index 898e5614..2c520bd6 100644 --- a/threedi_plugin_tool.py +++ b/threedi_plugin_tool.py @@ -1,8 +1,12 @@ from qgis.PyQt.QtCore import QObject -from qgis.PyQt.QtXml import QDomElement, QDomDocument -import logging - -logger = logging.getLogger(__name__) +from qgis.PyQt.QtWidgets import QAction +from qgis.PyQt.QtXml import QDomDocument +from qgis.PyQt.QtXml import QDomElement +from threedi_results_analysis.threedi_plugin_model import ThreeDiGridItem +from threedi_results_analysis.threedi_plugin_model import ThreeDiResultItem +from typing import Callable +from typing import Dict +from typing import Tuple class ThreeDiPluginTool(QObject): @@ -27,3 +31,24 @@ def read(self, _: QDomElement) -> bool: def on_unload(self): """Called when the plugin is unloaded. Tool can cleanup necessary items""" pass + + def get_custom_actions(self) -> Dict[QAction, Tuple[Callable[[ThreeDiGridItem], None], Callable[[ThreeDiResultItem], None]]]: + """Called to retrieve the tool specific actions for the context-menu (right-button click) in Result Manager tree, including + optional separators. Tool needs to provide an implementation for both a grid item and for a result item, e.g: + + @pyqtSlot(ThreeDiGridItem) + def add_summary_grid(self, item:ThreeDiGridItem) -> None: + logger.info(f"grid {item.id}") + + @pyqtSlot(ThreeDiResultItem) + def add_summary_result(self, item:ThreeDiGridItem) -> None: + logger.info(f"result {item.id}") + + def get_custom_actions(self) -> Dict[QAction, Tuple[Callable[[ThreeDiGridItem], None], Callable[[ThreeDiResultItem], None]]]: + separator = QAction() + separator.setSeparator(True) + return {separator: (None, None), + QAction("Show flow summary"): (self.add_summary_grid, self.add_summary_result) + } + """ + return {} diff --git a/tool_animation/animation_styler.py b/tool_animation/animation_styler.py index c065edde..1a8b289f 100644 --- a/tool_animation/animation_styler.py +++ b/tool_animation/animation_styler.py @@ -1,25 +1,22 @@ from pathlib import Path -from typing import List -import logging - -import numpy as np - from qgis.core import QgsFeatureRequest from qgis.core import QgsMarkerSymbol from qgis.core import QgsProperty from qgis.core import QgsSymbolLayer from qgis.core import QgsVectorLayer from qgis.utils import iface - from threedi_results_analysis.datasource.result_constants import WET_CROSS_SECTION_AREA +from threedi_results_analysis.utils.color import color_ramp_from_data from threedi_results_analysis.utils.color import COLOR_RAMP_OCEAN_CURL from threedi_results_analysis.utils.color import COLOR_RAMP_OCEAN_DEEP from threedi_results_analysis.utils.color import COLOR_RAMP_OCEAN_HALINE -from threedi_results_analysis.utils.color import color_ramp_from_data +from typing import List + +import logging +import numpy as np + STYLES_ROOT = Path(__file__).parent / "layer_styles" -ANIMATION_LAYERS_NR_LEGEND_CLASSES = 24 -assert ANIMATION_LAYERS_NR_LEGEND_CLASSES % 2 == 0 DEFAULT_LOWER_THRESHOLD = 1e-6 logger = logging.getLogger(__name__) @@ -174,7 +171,7 @@ def style_animation_node_current( def style_animation_node_difference( - lyr: QgsVectorLayer, percentiles: List[float], variable: str, cells: bool, field_postfix="" + lyr: QgsVectorLayer, percentiles: List[float], variable: str, cells: bool, nr_of_classes, field_postfix="" ): """Applies styling to Animation Toolbar node layer in 'difference' mode""" @@ -201,7 +198,7 @@ def style_animation_node_difference( stop=abs_high, step=( (abs_high - abs_high * -1) - / (ANIMATION_LAYERS_NR_LEGEND_CLASSES - 2) + / (nr_of_classes - 2) ), ) ) diff --git a/tool_animation/map_animator.py b/tool_animation/map_animator.py index f22f0949..760a767e 100644 --- a/tool_animation/map_animator.py +++ b/tool_animation/map_animator.py @@ -1,37 +1,156 @@ +from bisect import bisect_left +from enum import Enum +from functools import lru_cache from qgis.core import NULL from qgis.core import QgsProject -from qgis.PyQt.QtCore import Qt, pyqtSlot +from qgis.PyQt.QtCore import pyqtSlot +from qgis.PyQt.QtCore import QSize +from qgis.PyQt.QtCore import Qt +from qgis.PyQt.QtGui import QDoubleValidator +from qgis.PyQt.QtGui import QIcon +from qgis.PyQt.QtGui import QIntValidator from qgis.PyQt.QtWidgets import QCheckBox from qgis.PyQt.QtWidgets import QComboBox -from qgis.PyQt.QtWidgets import QHBoxLayout, QGridLayout -from qgis.PyQt.QtWidgets import QWidget +from qgis.PyQt.QtWidgets import QDialog +from qgis.PyQt.QtWidgets import QDialogButtonBox +from qgis.PyQt.QtWidgets import QGridLayout from qgis.PyQt.QtWidgets import QGroupBox -from threedigrid.admin.constants import NO_DATA_VALUE +from qgis.PyQt.QtWidgets import QHBoxLayout +from qgis.PyQt.QtWidgets import QLabel +from qgis.PyQt.QtWidgets import QLineEdit +from qgis.PyQt.QtWidgets import QPushButton +from qgis.PyQt.QtWidgets import QVBoxLayout +from qgis.PyQt.QtWidgets import QWidget +from qgis.PyQt.QtXml import QDomDocument +from qgis.PyQt.QtXml import QDomElement +from threedi_results_analysis import PLUGIN_DIR +from threedi_results_analysis.datasource.result_constants import AGGREGATION_OPTIONS from threedi_results_analysis.datasource.result_constants import DISCHARGE from threedi_results_analysis.datasource.result_constants import H_TYPES from threedi_results_analysis.datasource.result_constants import NEGATIVE_POSSIBLE from threedi_results_analysis.datasource.result_constants import Q_TYPES from threedi_results_analysis.datasource.result_constants import WATERLEVEL -from threedi_results_analysis.datasource.result_constants import AGGREGATION_OPTIONS from threedi_results_analysis.datasource.threedi_results import ThreediResult -from threedi_results_analysis.threedi_plugin_model import ThreeDiResultItem, ThreeDiGridItem -from threedi_results_analysis.utils.user_messages import StatusProgressBar -from threedi_results_analysis.utils.utils import generate_parameter_config, is_substance_variable, pretty +from threedi_results_analysis.threedi_plugin_model import ThreeDiGridItem +from threedi_results_analysis.threedi_plugin_model import ThreeDiResultItem from threedi_results_analysis.utils.timing import timing +from threedi_results_analysis.utils.user_messages import pop_up_critical +from threedi_results_analysis.utils.user_messages import StatusProgressBar +from threedi_results_analysis.utils.utils import generate_parameter_config +from threedi_results_analysis.utils.utils import is_substance_variable +from threedi_results_analysis.utils.utils import pretty +from threedigrid.admin.constants import NO_DATA_VALUE from typing import List -import threedi_results_analysis.tool_animation.animation_styler as styler import copy import logging import math import numpy as np -from bisect import bisect_left -from functools import lru_cache +import threedi_results_analysis.tool_animation.animation_styler as styler logger = logging.getLogger(__name__) +class MethodEnum(str, Enum): + PRETTY = "Pretty Breaks" + PERCENTILE = "Equal Count (Quantile)" + + +class MapAnimatorSettings(object): + lower_cutoff_percentile: float = 2.0 + upper_cutoff_percentile: float = 98.0 + method: MethodEnum = MethodEnum.PRETTY + nr_classes: int = 24 # Must be EVEN! + + def __str__(self): + return f"{self.lower_cutoff_percentile} {self.upper_cutoff_percentile} {self.method.value} {self.nr_classes}" + + +class MapAnimatorSettingsdialog(QDialog): + def __init__(self, parent, title: str, default_settings: MapAnimatorSettings): + super().__init__(parent) + self.setWindowTitle(f"Visualisation settings for {title}") + + layout = QVBoxLayout(self) + self.setLayout(layout) + + settings_group = QGroupBox(self) + settings_group.setLayout(QGridLayout()) + + # Set up GUI and populate with settings + settings_group.layout().addWidget(QLabel("Lower cutoff percentile:"), 0, 0) + self.lower_cutoff_percentile_lineedit = QLineEdit(str(default_settings.lower_cutoff_percentile), settings_group) + lower_percentile_validator = QDoubleValidator(0.0, 100.0, 2, self.lower_cutoff_percentile_lineedit) + lower_percentile_validator.setNotation(QDoubleValidator.StandardNotation) + self.lower_cutoff_percentile_lineedit.setValidator(lower_percentile_validator) + settings_group.layout().addWidget(self.lower_cutoff_percentile_lineedit, 0, 1) + + settings_group.layout().addWidget(QLabel("Upper cutoff percentile:"), 1, 0) + self.upper_cutoff_percentile_lineedit = QLineEdit(str(default_settings.upper_cutoff_percentile), settings_group) + upper_percentile_validator = QDoubleValidator(0.0, 100.0, 2, self.upper_cutoff_percentile_lineedit) + upper_percentile_validator.setNotation(QDoubleValidator.StandardNotation) + self.upper_cutoff_percentile_lineedit.setValidator(upper_percentile_validator) + settings_group.layout().addWidget(self.upper_cutoff_percentile_lineedit, 1, 1) + + settings_group.layout().addWidget(QLabel("Method:"), 2, 0) + self.method_combo = QComboBox(settings_group) + self.method_combo.addItems([s.value for s in MethodEnum]) + for c in range(self.method_combo.count()): + if default_settings.method.value == self.method_combo.itemText(c): + self.method_combo.setCurrentIndex(c) + break + settings_group.layout().addWidget(self.method_combo, 2, 1) + + explanation_msg = "The number of classes used in the styling may differ slightly from the number of classes set here." + class_label = QLabel("Number of classes: 🛈") + class_label.setToolTip(explanation_msg) + settings_group.layout().addWidget(class_label, 3, 0) + self.nr_classes_lineedit = QLineEdit(str(default_settings.nr_classes), settings_group) + self.nr_classes_lineedit.setToolTip(explanation_msg) + self.nr_classes_lineedit.setValidator(QIntValidator(2, 42, self.nr_classes_lineedit)) + settings_group.layout().addWidget(self.nr_classes_lineedit, 3, 1) + + layout.addWidget(settings_group) + + buttonBox = QDialogButtonBox(QDialogButtonBox.Ok | QDialogButtonBox.Cancel) + buttonBox.accepted.connect(self.accept) + buttonBox.rejected.connect(self.reject) + + layout.addWidget(buttonBox) + + def accept(self) -> None: + # Check logic + if int(self.nr_classes_lineedit.text()) % 2 != 0: + pop_up_critical("Number of classes should be even.") + return + if int(self.nr_classes_lineedit.text()) <= 0 or int(self.nr_classes_lineedit.text()) > 42: + pop_up_critical("The maximum number of classes is 42") + return + upper_cutoff_percentile = float(self.upper_cutoff_percentile_lineedit.text()) + lower_cutoff_percentile = float(self.lower_cutoff_percentile_lineedit.text()) + if upper_cutoff_percentile < 0 or upper_cutoff_percentile > 100: + pop_up_critical("The upper cutoff percentile should be greater than 0 and less than 100.") + return + if lower_cutoff_percentile < 0 or lower_cutoff_percentile > 100: + pop_up_critical("The lower cutoff percentile should be greater than 0 and less than 100.") + return + + if upper_cutoff_percentile <= lower_cutoff_percentile: + pop_up_critical("The upper cutoff percentile should be greater than the lower cutoff percentile.") + return + + return super().accept() + + def get_settings(self) -> MapAnimatorSettings: + result = MapAnimatorSettings() + result.lower_cutoff_percentile = float(self.lower_cutoff_percentile_lineedit.text()) + result.upper_cutoff_percentile = float(self.upper_cutoff_percentile_lineedit.text()) + result.nr_classes = int(self.nr_classes_lineedit.text()) + result.method = MethodEnum(self.method_combo.currentText()) + return result + + def get_layer_by_id(layer_id): return QgsProject.instance().mapLayer(layer_id) @@ -54,9 +173,9 @@ def threedi_result_legend_class_bounds( lower_cutoff_percentile: float, upper_cutoff_percentile: float, relative_to_t0: bool, + nr_classes: int, simple=False, - method: str = "pretty", - nr_classes: int = styler.ANIMATION_LAYERS_NR_LEGEND_CLASSES + method: str = "Pretty Breaks", ) -> List[float]: """ Calculate percentile values given variable in a 3Di results netcdf @@ -74,7 +193,7 @@ def threedi_result_legend_class_bounds( :param lower_threshold: ignore values below this threshold :param relative_to_t0: calculate percentiles on difference w/ initial values (applied before absolute) :param nodatavalue: ignore these values - :param method: 'pretty' (pretty breaks) or 'percentile' (equal count) + :param method: 'Pretty Breaks' or 'Equal Count (Quantile)' """ class_bounds_empty = [0] * nr_classes @@ -158,25 +277,32 @@ def threedi_result_legend_class_bounds( ) ] - if method == "pretty": + if method == MethodEnum.PRETTY.value: try: result = pretty(values_cutoff, n=nr_classes) except ValueError: # All values are the same result = class_bounds_empty - elif method == "percentile": + elif method == MethodEnum.PERCENTILE.value: result = np.nanpercentile( values_cutoff, class_bounds_percentiles ).tolist() else: - raise ValueError("'method' must be one of 'pretty', 'percentile'") + raise ValueError(f"'method' must be one of '{MethodEnum.PRETTY.value}', '{MethodEnum.PERCENTILE.value}'") real_min = 0 if absolute else np.nanmin(values).item() real_max = np.nanmax(values).item() if result[0] == real_min: if lower_threshold > real_min: result = np.insert(result, 1, lower_threshold) # create a class for all values that can be regarded as 0 - else: + elif real_min < result[0]: result = np.insert(result, 0, real_min) + if absolute: + assert real_min == 0 + if lower_threshold > result[0] and lower_threshold < result[1]: + result = np.insert(result, 1, lower_threshold) # insert a 0-regarded class + else: # real_min > result[0]: + pass + if result[-1] != real_max: result = np.insert(result, len(result), real_max) @@ -188,11 +314,14 @@ class MapAnimator(QGroupBox): def __init__(self, parent, model): - super().__init__("Visualization settings", parent) + super().__init__("Visualisation settings", parent) self.model = model self.node_parameters = None self.line_parameters = None + self.node_parameter_setting = {} + self.line_parameter_setting = {} + self.current_datetime = None self.setup_ui(parent) @@ -208,6 +337,7 @@ def results_changed(self, item: ThreeDiResultItem): self._update_parameter_attributes() self._update_parameter_combo_boxes() + self._update_parameter_settings() if not active: return @@ -216,19 +346,128 @@ def results_changed(self, item: ThreeDiResultItem): self.update_results() # iface.mapCanvas().refresh() + def _update_parameter_settings(self): + # Update cached parameter settings, remove param if no longer present + param_settings_to_delete = [] + for param_key in self.node_parameter_setting: + found = False + for param in self.node_parameters.values(): + if param_key == f"{param['name']}-{param['unit']}-{param['parameters']}": + # This parameter is still present in results, so keep it. + found = True + break + if not found: + param_settings_to_delete.append(param_key) + for param_key in param_settings_to_delete: + logger.info(f"Removing settings for {param_key} as no longer present in all results") + del self.node_parameter_setting[param_key] + + param_settings_to_delete.clear() + for param_key in self.line_parameter_setting: + found = False + for param in self.line_parameters.values(): + if param_key == f"{param['name']}-{param['unit']}-{param['parameters']}": + found = True + break + if not found: + param_settings_to_delete.append(param_key) + for param_key in param_settings_to_delete: + logger.info(f"Removing settings for {param_key} as no longer present in all results") + del self.line_parameter_setting[param_key] + + def write(self, doc: QDomDocument, xml_elem: QDomElement) -> bool: + """Called when a QGS project is written, allowing each tool to presist + additional info int the dedicated xml tools node.""" + + animator_node = doc.createElement("map_animator") + xml_elem.appendChild(animator_node) + + node_parameter_setting_element = doc.createElement("node_parameter_setting") + for param_key, settings in self.node_parameter_setting.items(): + node_parameter_setting_element.appendChild(MapAnimator._setting_to_xml(doc, param_key, settings)) + + line_parameter_setting_element = doc.createElement("line_parameter_setting") + for param_key, settings in self.line_parameter_setting.items(): + line_parameter_setting_element.appendChild(MapAnimator._setting_to_xml(doc, param_key, settings)) + + animator_node.appendChild(node_parameter_setting_element) + animator_node.appendChild(line_parameter_setting_element) + + return True + + @staticmethod + def _setting_to_xml(doc: QDomDocument, param_key: str, settings: MapAnimatorSettings) -> QDomElement: + settings_element = doc.createElement("ParameterSetting") + # We write the parameter key as attribute as it might contain spaces + settings_element.setAttribute("parameter", param_key) + lower_cutoff_percentile_element = doc.createElement("lower_cutoff_percentile") + lower_cutoff_percentile_element.appendChild(doc.createTextNode(str(settings.lower_cutoff_percentile))) + upper_cutoff_percentile_element = doc.createElement("upper_cutoff_percentile") + upper_cutoff_percentile_element.appendChild(doc.createTextNode(str(settings.upper_cutoff_percentile))) + method_element = doc.createElement("method") + method_element.appendChild(doc.createTextNode(settings.method.value)) + nr_classes_element = doc.createElement("nr_classes") + nr_classes_element.appendChild(doc.createTextNode(str(settings.nr_classes))) + settings_element.appendChild(lower_cutoff_percentile_element) + settings_element.appendChild(upper_cutoff_percentile_element) + settings_element.appendChild(method_element) + settings_element.appendChild(nr_classes_element) + return settings_element + + def read(self, xml_elem: QDomElement) -> bool: + animator_node = xml_elem.firstChildElement("map_animator") + if not animator_node: + logger.info("No animation settings in project") + return True + + self.node_parameter_setting.clear() + nodes_parameter_settings = animator_node.elementsByTagName("node_parameter_setting") + assert nodes_parameter_settings.count() == 1 + nodes_parameter_settings = nodes_parameter_settings.item(0).toElement() + param_nodes = nodes_parameter_settings.childNodes() + for i in range(param_nodes.count()): + param_node = param_nodes.at(i).toElement() + param_key = param_node.attribute("parameter") + self.node_parameter_setting[param_key] = MapAnimator._setting_from_xml(param_node) + + self.line_parameter_setting.clear() + line_parameter_settings = animator_node.elementsByTagName("line_parameter_setting") + assert line_parameter_settings.count() == 1 + line_parameter_settings = line_parameter_settings.item(0).toElement() + param_lines = line_parameter_settings.childNodes() + for i in range(param_lines.count()): + param_line = param_lines.at(i).toElement() + param_key = param_line.attribute("parameter") + self.line_parameter_setting[param_key] = MapAnimator._setting_from_xml(param_line) + + return True + + @staticmethod + def _setting_from_xml(param_node: QDomElement) -> MapAnimatorSettings: + setting = MapAnimatorSettings() + assert param_node.elementsByTagName("lower_cutoff_percentile").count() == 1 + setting.lower_cutoff_percentile = float(param_node.elementsByTagName("lower_cutoff_percentile").item(0).toElement().text()) + assert param_node.elementsByTagName("upper_cutoff_percentile").count() == 1 + setting.upper_cutoff_percentile = float(param_node.elementsByTagName("upper_cutoff_percentile").item(0).toElement().text()) + assert param_node.elementsByTagName("method").count() == 1 + setting.method = MethodEnum(param_node.elementsByTagName("method").item(0).toElement().text()) + assert param_node.elementsByTagName("nr_classes").count() == 1 + setting.nr_classes = int(param_node.elementsByTagName("nr_classes").item(0).toElement().text()) + return setting + def _update_parameter_attributes(self): config = self._get_active_parameter_config() self.line_parameters = {r["name"]: r for r in config["q"]} self.node_parameters = {r["name"]: r for r in config["h"]} def _style_line_layers(self, result_item: ThreeDiResultItem, progress_bar): + current_line_settings = self.line_parameter_setting.get(self.current_line_parameter_key, MapAnimatorSettings()) threedi_result = result_item.threedi_result line_parameter_class_bounds, _ = self._get_class_bounds_line( - threedi_result, self.current_line_parameter["parameters"] + threedi_result, self.current_line_parameter["parameters"], current_line_settings ) grid_item = result_item.parent() assert isinstance(grid_item, ThreeDiGridItem) - logger.info("Styling flowline layer") layer_id = grid_item.layer_ids["flowline"] virtual_field_name = result_item._result_field_names[layer_id][0] postfix = virtual_field_name[6:] # remove "result" prefix @@ -243,16 +482,16 @@ def _style_line_layers(self, result_item: ThreeDiResultItem, progress_bar): def _style_node_layers(self, result_item: ThreeDiResultItem, progress_bar): """ Compute class bounds and apply style to node and cell layers. """ + current_node_settings = self.node_parameter_setting.get(self.current_node_parameter_key, MapAnimatorSettings()) threedi_result = result_item.threedi_result node_parameter_class_bounds, _ = self._get_class_bounds_node( - threedi_result, self.current_node_parameter["parameters"], + threedi_result, self.current_node_parameter["parameters"], current_node_settings ) # Adjust the styling of the grid layer based on the bounds and result field name grid_item = result_item.parent() assert isinstance(grid_item, ThreeDiGridItem) - logger.info("Styling node layer") layer_id = grid_item.layer_ids["node"] layer = get_layer_by_id(layer_id) virtual_field_name = result_item._result_field_names[layer_id][0] @@ -263,6 +502,7 @@ def _style_node_layers(self, result_item: ThreeDiResultItem, progress_bar): node_parameter_class_bounds, self.current_node_parameter["parameters"], False, + current_node_settings.nr_classes, postfix, ) else: @@ -288,6 +528,7 @@ def _style_node_layers(self, result_item: ThreeDiResultItem, progress_bar): node_parameter_class_bounds, self.current_node_parameter["parameters"], True, + current_node_settings.nr_classes, postfix, ) else: @@ -308,6 +549,14 @@ def current_line_parameter(self): def current_node_parameter(self): return self.node_parameters[self.node_parameter_combo_box.currentText()] + @property + def current_node_parameter_key(self): + return f"{self.current_node_parameter['name']}-{self.current_node_parameter['unit']}-{self.current_node_parameter['parameters']}" + + @property + def current_line_parameter_key(self): + return f"{self.current_line_parameter['name']}-{self.current_line_parameter['unit']}-{self.current_line_parameter['parameters']}" + def _restyle(self, lines, nodes): result_items = self.model.get_results(checked_only=True) total = (int(lines) + 2 * int(nodes)) * len(result_items) @@ -330,7 +579,7 @@ def _restyle_and_update_nodes(self): self._restyle(lines=False, nodes=True) self.update_results() - def _get_class_bounds_node(self, threedi_result, node_variable): + def _get_class_bounds_node(self, threedi_result, node_variable, settings: MapAnimatorSettings): base_nc_name = strip_agg_options(node_variable) if ( base_nc_name in NEGATIVE_POSSIBLE and NEGATIVE_POSSIBLE[base_nc_name] @@ -345,10 +594,11 @@ def _get_class_bounds_node(self, threedi_result, node_variable): variable=node_variable, absolute=False, lower_threshold=lower_threshold, - lower_cutoff_percentile=2, - upper_cutoff_percentile=98, + lower_cutoff_percentile=settings.lower_cutoff_percentile, + upper_cutoff_percentile=settings.upper_cutoff_percentile, relative_to_t0=self.difference_checkbox.isChecked(), - method="pretty", + nr_classes=settings.nr_classes, + method=settings.method, ) with timing('percentiles1'): surfacewater_bounds = threedi_result_legend_class_bounds( @@ -360,16 +610,17 @@ def _get_class_bounds_node(self, threedi_result, node_variable): ) return surfacewater_bounds, groundwater_bounds - def _get_class_bounds_line(self, threedi_result, line_variable): + def _get_class_bounds_line(self, threedi_result, line_variable, settings: MapAnimatorSettings): kwargs = dict( threedi_result=threedi_result, variable=line_variable, absolute=True, lower_threshold=styler.DEFAULT_LOWER_THRESHOLD, - lower_cutoff_percentile=2, - upper_cutoff_percentile=98, + lower_cutoff_percentile=settings.lower_cutoff_percentile, + upper_cutoff_percentile=settings.upper_cutoff_percentile, relative_to_t0=self.difference_checkbox.isChecked(), - method="pretty", + nr_classes=settings.nr_classes, + method=settings.method, ) with timing('percentiles3'): surfacewater_bounds = threedi_result_legend_class_bounds( @@ -437,7 +688,7 @@ def _get_active_parameter_config(self): available_wq_vars = threedi_result.available_water_quality_vars[:] # a copy parameter_config = generate_parameter_config( - available_subgrid_vars, agg_vars=agg_vars, wq_vars=available_wq_vars + available_subgrid_vars, agg_vars=agg_vars, wq_vars=available_wq_vars, sca_vars=None ) def _intersection(a: List, b: List): @@ -577,7 +828,13 @@ def setup_ui(self, parent_widget: QWidget): self.line_parameter_combo_box = QComboBox(line_group) self.line_parameter_combo_box.setSizeAdjustPolicy(QComboBox.AdjustToContents) self.line_parameter_combo_box.setToolTip("Choose flowline variable to display") - line_group.layout().addWidget(self.line_parameter_combo_box, 0, 0, Qt.AlignTop) + line_group.layout().addWidget(self.line_parameter_combo_box, 0, 0, 1, 2, Qt.AlignTop) + equalizer_icon = QIcon(str(PLUGIN_DIR / "icons" / "sliders.svg")) + + setting_button_line = QPushButton(equalizer_icon, "", line_group) + setting_button_line.setFixedSize(QSize(26, 26)) + setting_button_line.clicked.connect(self.show_line_settings) + line_group.layout().addWidget(setting_button_line, 1, 1) self.HLayout.addWidget(line_group) @@ -586,7 +843,7 @@ def setup_ui(self, parent_widget: QWidget): self.node_parameter_combo_box = QComboBox(node_group) self.node_parameter_combo_box.setSizeAdjustPolicy(QComboBox.AdjustToContents) self.node_parameter_combo_box.setToolTip("Choose node variable to display") - node_group.layout().addWidget(self.node_parameter_combo_box) + node_group.layout().addWidget(self.node_parameter_combo_box, 0, 0, 1, 2) self.difference_checkbox = QCheckBox("Relative", self) self.difference_checkbox.setToolTip( @@ -595,6 +852,11 @@ def setup_ui(self, parent_widget: QWidget): node_group.layout().addWidget(self.difference_checkbox, 1, 0) + setting_button_node = QPushButton(equalizer_icon, "", line_group) + setting_button_node.setFixedSize(QSize(26, 26)) + setting_button_node.clicked.connect(self.show_node_settings) + node_group.layout().addWidget(setting_button_node, 1, 1) + self.HLayout.addWidget(node_group) self.line_parameter_combo_box.activated.connect(self._restyle_and_update_lines) @@ -603,6 +865,28 @@ def setup_ui(self, parent_widget: QWidget): self.setEnabled(False) + @pyqtSlot(bool) + def show_node_settings(self, _: bool): + current_node_settings = MapAnimatorSettings() + if self.current_node_parameter_key in self.node_parameter_setting: + current_node_settings = self.node_parameter_setting[self.current_node_parameter_key] + + dialog = MapAnimatorSettingsdialog(self, self.current_node_parameter['name'], current_node_settings) + if dialog.exec(): + self.node_parameter_setting[self.current_node_parameter_key] = dialog.get_settings() + self._restyle(lines=False, nodes=True) + + @pyqtSlot(bool) + def show_line_settings(self, _: bool): + current_line_settings = MapAnimatorSettings() + if self.current_line_parameter_key in self.line_parameter_setting: + current_line_settings = self.line_parameter_setting[self.current_line_parameter_key] + + dialog = MapAnimatorSettingsdialog(self, self.current_line_parameter['name'], current_line_settings) + if dialog.exec(): + self.line_parameter_setting[self.current_line_parameter_key] = dialog.get_settings() + self._restyle(lines=True, nodes=False) + @staticmethod def index_to_duration(index, timestamps): """Return the duration between start of simulation and the selected time index diff --git a/tool_flow_summary/__init__.py b/tool_flow_summary/__init__.py new file mode 100644 index 00000000..e69de29b diff --git a/tool_flow_summary/flow_summary.py b/tool_flow_summary/flow_summary.py new file mode 100644 index 00000000..ab92471d --- /dev/null +++ b/tool_flow_summary/flow_summary.py @@ -0,0 +1,205 @@ +# -*- coding: utf-8 -*- +from qgis.PyQt.QtCore import pyqtSlot +from qgis.PyQt.QtCore import Qt +from qgis.PyQt.QtWidgets import QAction +from qgis.PyQt.QtWidgets import QDialog +from qgis.PyQt.QtWidgets import QGridLayout +from qgis.PyQt.QtWidgets import QGroupBox +from qgis.PyQt.QtWidgets import QHBoxLayout +from qgis.PyQt.QtWidgets import QPushButton +from qgis.PyQt.QtWidgets import QSizePolicy +from qgis.PyQt.QtWidgets import QSpacerItem +from qgis.PyQt.QtWidgets import QWidget +from threedi_results_analysis.threedi_plugin_model import ThreeDiGridItem +from threedi_results_analysis.threedi_plugin_model import ThreeDiResultItem +from threedi_results_analysis.threedi_plugin_tool import ThreeDiPluginTool +from threedi_results_analysis.tool_flow_summary.variable_table import VariableTable +from typing import Callable +from typing import Dict +from typing import List +from typing import Tuple + +import json +import logging +import os + + +logger = logging.getLogger(__name__) + +GROUP_NAMES = [("general_information", Qt.AlignLeft), ("volume_balance", Qt.AlignRight), ("volume_balance_of_0d_model", Qt.AlignRight)] + + +class FlowSummaryTool(ThreeDiPluginTool): + + def __init__(self, parent, iface, model): + super().__init__(parent) + self.iface = iface + self.model = model + self.first_time = True + + # The list of shown results in the summary, idx+1 corresponding to column idx in the table + self.result_ids : List[int] = [] + # The map of group name to table + self.tables: Dict[str, VariableTable] = {} + + self.setup_ui() + + def setup_ui(self) -> None: + self.icon_path = os.path.join(os.path.dirname(os.path.dirname(__file__)), "icons", "icon_summary.png") + self.menu_text = "Flow summary tool" + self.main_widget = QDialog(None) + self.main_widget.setWindowTitle("Flow summary") + self.main_widget.setLayout(QGridLayout()) + + for group_name, info_alignment in GROUP_NAMES: + group_title = (group_name[0].capitalize() + group_name[1:]).replace("_", " ") + group_title = group_title.replace("0d", "0D") + variable_group = QGroupBox(group_title, self.main_widget) + variable_group.setStyleSheet("QGroupBox { font-weight: bold; }") + variable_group.setLayout(QGridLayout()) + + self.tables[group_name] = VariableTable(info_alignment, variable_group) + variable_group.layout().addWidget(self.tables[group_name]) + + self.main_widget.layout().addWidget(variable_group) + + self.main_widget.setEnabled(True) + self.main_widget.hide() + self.main_widget.setWindowFlags(Qt.WindowStaysOnTopHint) + + # Add Ok button + button_widget = QWidget(self.main_widget) + button_widget.setLayout(QHBoxLayout(button_widget)) + reset_button = QPushButton("Reset column widths", button_widget) + button_widget.layout().addWidget(reset_button) + reset_button.clicked.connect(self._reset_column_widths) + spacer_item = QSpacerItem(40, 20, QSizePolicy.Expanding, QSizePolicy.Minimum) + button_widget.layout().addItem(spacer_item) + ok_button = QPushButton("OK", button_widget) + button_widget.layout().addWidget(ok_button, alignment=Qt.AlignRight) + self.main_widget.layout().addWidget(button_widget) + ok_button.clicked.connect(self.main_widget.hide) + + # set comfortable start width + self.main_widget.resize(800, self.main_widget.height()) + + def _reset_column_widths(self) -> None: + for table in self.tables.values(): + if table.columnCount() == 0: + return + + max_width_variables = 0 + for table in self.tables.values(): + group_width = table.get_preferred_variable_column_width() + if group_width > max_width_variables: + max_width_variables = group_width + + for table in self.tables.values(): + table.setColumnWidth(0, max_width_variables) + if table.columnCount() == 1: + continue + scrollbar_width = table.horizontalScrollBar().geometry().height() + # evenly spread remainder of columns + column_width = (table.width() - max_width_variables - scrollbar_width) // (table.columnCount() - 1) + column_width -= 3 # prevent scrollbar + for c in range(1, table.columnCount()): + table.setColumnWidth(c, column_width) + + def add_summary_grid(self, item: ThreeDiGridItem) -> None: + results = [] + self.model.get_results_from_item(item=item, checked_only=False, results=results) + for result in results: + self.add_summary_result(result) + + def add_summary_result(self, item: ThreeDiResultItem, show: bool = True) -> None: + + if item.id in self.result_ids: + logger.warning("Result already added to flow summary, ignoring...") + return + + self.result_ids.append(item.id) + + # find and parse the result files + flow_summary_path = item.path.parent / "flow_summary.json" + if not flow_summary_path.exists(): + logger.warning(f"Flow summary file from Result {item.text()} cannot be found.") + # TODO: ideally make red, but unclear how to style individual items in header + for group_name, _ in GROUP_NAMES: + self.tables[group_name].add_summary_results(item.text(), dict()) + else: + # retrieve all the entries in this file + with flow_summary_path.open() as file: + data = json.load(file) + + for group_name, _ in GROUP_NAMES: + assert group_name in self.tables + if group_name in data: + group_data = data[group_name] + else: + group_data = {} + self.tables[group_name].add_summary_results(item.text(), group_data) + + self._reset_column_widths() + if show: + self.main_widget.show() + + def remove_summary_grid(self, item: ThreeDiGridItem) -> None: + results = [] + self.model.get_results_from_item(item=item, checked_only=False, results=results) + for result in results: + self.remove_summary_result(result) + + def remove_summary_result(self, item: ThreeDiResultItem) -> None: + self.result_removed(item) + + def get_custom_actions(self) -> Dict[QAction, Tuple[Callable[[ThreeDiGridItem], None], Callable[[ThreeDiResultItem], None]]]: + separator = QAction() + separator.setSeparator(True) + return {separator: (None, None), + QAction("Add to flow summary"): (self.add_summary_grid, self.add_summary_result), + QAction("Remove from flow summary"): (self.remove_summary_grid, self.remove_summary_result) + } + + def on_unload(self) -> None: + del self.main_widget + self.main_widget = None + + @pyqtSlot(ThreeDiResultItem) + def result_removed(self, result_item: ThreeDiResultItem): + # Remove column if required, pop from self.result_ids + try: + idx = self.result_ids.index(result_item.id) + except ValueError: + return # result not in summary + + self.result_ids.pop(idx) + for table in self.tables: + self.tables[table].remove_result(idx+1) + + # if empty: fully clean tables + if len(self.result_ids) == 0: + for table in self.tables: + self.tables[table].clean_results() + + self._reset_column_widths() + + @pyqtSlot(ThreeDiResultItem) + def result_changed(self, result_item: ThreeDiResultItem): + try: + idx = self.result_ids.index(result_item.id) + except ValueError: + return # result not in summary + + for table in self.tables: + self.tables[table].change_result(idx+1, result_item.text()) + + @pyqtSlot(ThreeDiResultItem) + def result_added(self, item: ThreeDiResultItem): + self.action_icon.setEnabled(self.model.number_of_results() > 0) + self.add_summary_result(item, False) + + def run(self) -> None: + self.main_widget.show() + if self.first_time: + self._reset_column_widths() + self.first_time = False diff --git a/tool_flow_summary/test_tool_flow_summary.py b/tool_flow_summary/test_tool_flow_summary.py new file mode 100644 index 00000000..d384dc11 --- /dev/null +++ b/tool_flow_summary/test_tool_flow_summary.py @@ -0,0 +1,116 @@ +from qgis.PyQt.QtCore import Qt +from threedi_results_analysis.tool_flow_summary.flow_summary import FlowSummaryTool +from threedi_results_analysis.tool_flow_summary.flow_summary import VariableTable + +import mock +import unittest + + +TEST_DATA = { + "calculation_node_with_max_volume_error": 13973, + "default_timestep": { + "units": "s", + "value": 5.0146 + }, + "maximum_timestep": { + "units": "s", + "value": 5.0379 + }, + "minimum_timestep": { + "units": "s", + "value": 5.0 + }, + "model_id": 65899, + "model_name": "ilan2023-16", + "model_type": "1D/2D", + "revision_id": 58400, + "schematisation_id": 6745, + "simulation_id": 217430, + "simulation_start": "2024-09-03 12:00:00", + "simulation_time": { + "units": "s", + "value": 14400.1904 + } + } + +TEST_DATA2 = { + "calculation_node_with_max_volume_error": 13973, + "additional": 13973, + "default_timestep": { + "units": "s", + "value": 10000 + }, + "additional2": { + "units": "s", + "value": 5.0146 + }, + } + +TEST_DATA3 = { + "additional": 2, + "additional2": { + "units": "s", + "value": 3 + }, + "simulation_time": { + "units": "s", + "value": 3.14 + } + } + +# same param, but different unit +TEST_DATA4 = { + "additional2": { + "units": "h", + "value": 44444 + }, + } + + +class TestFlowSummaryTool(unittest.TestCase): + def setUp(self): + """test whether FlowSummaryTool can be instantiated""" + iface = mock.Mock() + self.flow_summary = FlowSummaryTool(None, iface, None) + + def test_icon_path_is_set(self): + self.assertEqual( + self.flow_summary.icon_path, "/root/.local/share/QGIS/QGIS3/profiles/default/python/plugins/threedi_results_analysis/icons/icon_summary.png" + ) + + def test_result_addition_removal(self): + table = VariableTable(Qt.AlignRight, None) + assert table.columnCount() == 1 + assert table.rowCount() == 0 + + table.add_summary_results("test", TEST_DATA) + assert table.columnCount() == 2 + assert table.rowCount() == 12 + + table.add_summary_results("test", TEST_DATA2) + assert table.columnCount() == 3 + assert table.rowCount() == 14 + + table.add_summary_results("test", TEST_DATA3) + assert table.columnCount() == 4 + assert table.rowCount() == 14 + + assert table.item(1, 2).text() == "10,000" + + # same data, but different unit + table.add_summary_results("test", TEST_DATA4) + assert table.columnCount() == 5 + assert table.rowCount() == 15 + + assert table.item(14, 4).text() == "44,444" + + # remove TEST_DATA + table.remove_result(1) + assert table.columnCount() == 4 + assert table.rowCount() == 15 + + assert table.item(14, 3).text() == "44,444" + + table.clean_results() + assert table.columnCount() == 1 + assert table.rowCount() == 0 diff --git a/tool_flow_summary/variable_table.py b/tool_flow_summary/variable_table.py new file mode 100644 index 00000000..34059a03 --- /dev/null +++ b/tool_flow_summary/variable_table.py @@ -0,0 +1,139 @@ + +from qgis.PyQt.QtCore import QLocale +from qgis.PyQt.QtCore import Qt +from qgis.PyQt.QtWidgets import QAbstractItemView +from qgis.PyQt.QtWidgets import QApplication +from qgis.PyQt.QtWidgets import QHeaderView +from qgis.PyQt.QtWidgets import QTableWidget +from qgis.PyQt.QtWidgets import QTableWidgetItem +from typing import List +from typing import Tuple + + +class VariableTable(QTableWidget): + def __init__(self, variable_alignment, parent): + super().__init__(0, 1, parent) + self.variable_alignment = variable_alignment + self.setHorizontalHeaderLabels([""]) + self.horizontalHeader().setSectionResizeMode(0, QHeaderView.Interactive) + + self.verticalHeader().hide() + self.setSortingEnabled(False) + self.setSelectionMode(QAbstractItemView.ContiguousSelection) + + # for proper aligning, we always need to reserve space for the scrollbar + self.setVerticalScrollBarPolicy(Qt.ScrollBarAlwaysOn) + self.horizontalHeader().setStretchLastSection(True) + + # The list of parameters shown in the summary, idx corresponding to row idx in the table + self.param_names : List[str] = [] + + def add_summary_results(self, result_text: str, group_data): + header_item = QTableWidgetItem(result_text) + self.insertColumn(self.columnCount()) + self.setHorizontalHeaderItem(self.columnCount()-1, header_item) + + for param in group_data: + param_name, param_value = self._format_variable(param, group_data[param]) + + # Check if we've added this parameter before, then use that row idx, + # otherwise append to bottom of table + try: + param_index = self.param_names.index(param_name) + except ValueError: + param_index = len(self.param_names) + self.param_names.append(param_name) + # Add a new row and set the parameter name + assert param_index == self.rowCount() + self.insertRow(param_index) + item = QTableWidgetItem(param_name) + item.setFlags(item.flags() ^ Qt.ItemIsEditable) + item.setTextAlignment(Qt.AlignLeft) + self.setItem(param_index, 0, item) + + item = QTableWidgetItem(param_value) + item.setTextAlignment(self.variable_alignment) + item.setFlags(item.flags() ^ Qt.ItemIsEditable) + self.setItem(param_index, self.columnCount()-1, item) + + for idx in range(0, self.columnCount()): + self.horizontalHeader().setSectionResizeMode(idx, QHeaderView.Interactive) + + def resizeEvent(self, event): + self.resizeRowsToContents() + super().resizeEvent(event) + + def keyPressEvent(self, event): + # https://stackoverflow.com/questions/1230222/selected-rows-in-qtableview-copy-to-qclipboard/24133289#24133289 + if event.key() == Qt.Key_C and event.modifiers() & Qt.ControlModifier: + + indexes = self.selectedIndexes() + current_text = "" + current_row = 0 # To determine when to insert newlines + + indexes = sorted(indexes) # Necessary, otherwise they are in column order + + for index in indexes: + if current_text == "": # first line + pass + elif index.row() != current_row: # new row + current_text += "\n" + else: + current_text += "\t" + + current_row = index.row() + current_item = self.item(index.row(), index.column()) + if current_item: + current_text += str(self.item(index.row(), index.column()).text()) + + QApplication.clipboard().setText(current_text) + return + + super().keyPressEvent(event) + + def clean_results(self) -> None: + self.clearContents() + self.setColumnCount(1) + self.setRowCount(0) + self.setHorizontalHeaderLabels([""]) + self.param_names.clear() + + def remove_result(self, idx: int) -> None: + self.removeColumn(idx) + + def get_preferred_variable_column_width(self) -> int: + # iterate over the texts and determine the max width when resized to contents + self.resizeColumnToContents(0) + max_width = 0 + for r in range(self.rowCount()): + width = self.columnWidth(0) + if width > max_width: + max_width = width + + return max_width + + def change_result(self, idx: int, text: str) -> None: + self.setHorizontalHeaderItem(idx, QTableWidgetItem(text)) + + def _format_variable(self, param_name: str, param_data: dict) -> Tuple[str, str]: + + param_name = param_name.replace("_", " ") + if type(param_data) is dict: + param_name = f'{param_name} [{param_data["units"]}]' + param_data = param_data["value"] + + locale = QLocale() + + # numbers in 4 decimals, or scientific notation when not possible + if isinstance(param_data, float): + if param_data < 0.0001: + param_data = locale.toString(param_data, "g", 5) + else: + param_data = locale.toString(param_data, "f", 4) + elif isinstance(param_data, int): + param_data = locale.toString(param_data) + else: + param_data = str(param_data) + + param_name = param_name[0].capitalize() + param_name[1:] + return param_name, param_data diff --git a/tool_graph/graph_model.py b/tool_graph/graph_model.py index 3113b4fc..638ff60a 100644 --- a/tool_graph/graph_model.py +++ b/tool_graph/graph_model.py @@ -1,15 +1,17 @@ from collections import OrderedDict +from qgis.PyQt.QtCore import Qt +from qgis.PyQt.QtGui import QColor from random import randint from threedi_results_analysis.models.base import BaseModel -from threedi_results_analysis.models.base_fields import CheckboxField, CHECKBOX_FIELD +from threedi_results_analysis.models.base_fields import CHECKBOX_FIELD +from threedi_results_analysis.models.base_fields import CheckboxField from threedi_results_analysis.models.base_fields import ValueField +from threedi_results_analysis.utils.color import COLOR_LIST import logging import numpy as np import pyqtgraph as pg -from qgis.PyQt.QtCore import Qt -from qgis.PyQt.QtGui import QColor -from threedi_results_analysis.utils.color import COLOR_LIST + logger = logging.getLogger(__name__) @@ -143,13 +145,20 @@ def timeseries_table(self, parameters, absolute, time_units): if (parameters not in threedi_result.available_subgrid_map_vars and parameters not in threedi_result.available_aggregation_vars and - parameters not in [v["parameters"] for v in threedi_result.available_water_quality_vars]): + parameters not in [v["parameters"] for v in threedi_result.available_water_quality_vars] and + parameters not in [v["parameters"] for v in threedi_result.available_structure_control_actions_vars]): logger.warning(f"Parameter {parameters} not available in result {self.result.value.text()}") return EMPTY_TIMESERIES ga = threedi_result.get_gridadmin(parameters) if ga.has_pumpstations: - pump_fields = set(list(ga.pumps.Meta.composite_fields.keys())) + # In some gridadmin types pumps do not have a Meta attribute... In + # such cases (e.g. water quality) the attribute does not have a meaning and + # the timeserie should be empty. + try: + pump_fields = set(list(ga.pumps.Meta.composite_fields.keys())) + except AttributeError: + pump_fields = {} else: pump_fields = {} if self.object_type.value == "pump_linestring" and parameters not in pump_fields: @@ -158,7 +167,7 @@ def timeseries_table(self, parameters, absolute, time_units): return EMPTY_TIMESERIES timeseries = threedi_result.get_timeseries( - parameters, node_id=self.object_id.value, fill_value=np.NaN + parameters, node_id=self.object_id.value, fill_value=np.NaN, selected_object_type=self.object_type.value ) if timeseries.shape[1] == 1: logger.info("1-element timeserie, plotting empty serie") diff --git a/tool_graph/graph_view.py b/tool_graph/graph_view.py index 06c8df75..49df0cca 100644 --- a/tool_graph/graph_view.py +++ b/tool_graph/graph_view.py @@ -1,45 +1,47 @@ -from qgis.core import QgsFeatureRequest from qgis.core import Qgis -from qgis.core import QgsWkbTypes -from qgis.core import QgsValueMapFieldFormatter from qgis.core import QgsFeature +from qgis.core import QgsFeatureRequest +from qgis.core import QgsProject +from qgis.core import QgsValueMapFieldFormatter +from qgis.core import QgsVectorLayer +from qgis.core import QgsWkbTypes from qgis.gui import QgsMapToolIdentify from qgis.gui import QgsRubberBand -from qgis.core import QgsProject from qgis.PyQt.QtCore import pyqtSignal from qgis.PyQt.QtCore import pyqtSlot from qgis.PyQt.QtCore import QEvent from qgis.PyQt.QtCore import QMetaObject from qgis.PyQt.QtCore import QSize from qgis.PyQt.QtCore import Qt +from qgis.PyQt.QtGui import QColor +from qgis.PyQt.QtWidgets import QAbstractItemView +from qgis.PyQt.QtWidgets import QAction +from qgis.PyQt.QtWidgets import QActionGroup from qgis.PyQt.QtWidgets import QCheckBox +from qgis.PyQt.QtWidgets import QColorDialog from qgis.PyQt.QtWidgets import QComboBox -from qgis.PyQt.QtWidgets import QSplitter from qgis.PyQt.QtWidgets import QDockWidget from qgis.PyQt.QtWidgets import QHBoxLayout +from qgis.PyQt.QtWidgets import QMenu from qgis.PyQt.QtWidgets import QMessageBox -from qgis.PyQt.QtWidgets import QActionGroup, QToolButton from qgis.PyQt.QtWidgets import QSizePolicy from qgis.PyQt.QtWidgets import QSpacerItem +from qgis.PyQt.QtWidgets import QSplitter from qgis.PyQt.QtWidgets import QTableView -from qgis.PyQt.QtWidgets import QAbstractItemView from qgis.PyQt.QtWidgets import QTabWidget -from qgis.PyQt.QtWidgets import QMenu -from qgis.PyQt.QtWidgets import QAction +from qgis.PyQt.QtWidgets import QToolButton from qgis.PyQt.QtWidgets import QVBoxLayout from qgis.PyQt.QtWidgets import QWidget -from qgis.PyQt.QtWidgets import QColorDialog -from qgis.PyQt.QtGui import QColor +from threedi_results_analysis.datasource.threedi_results import normalized_object_type +from threedi_results_analysis.threedi_plugin_model import ThreeDiGridItem +from threedi_results_analysis.threedi_plugin_model import ThreeDiPluginModel +from threedi_results_analysis.threedi_plugin_model import ThreeDiResultItem from threedi_results_analysis.tool_graph.graph_model import LocationTimeseriesModel +from threedi_results_analysis.utils.constants import TOOLBOX_MESSAGE_TITLE from threedi_results_analysis.utils.user_messages import messagebar_message from threedi_results_analysis.utils.user_messages import statusbar_message from threedi_results_analysis.utils.utils import generate_parameter_config from threedi_results_analysis.utils.widgets import PenStyleWidget -from threedi_results_analysis.utils.constants import TOOLBOX_MESSAGE_TITLE -from qgis.core import QgsVectorLayer -from threedi_results_analysis.datasource.threedi_results import normalized_object_type -from threedi_results_analysis.threedi_plugin_model import ThreeDiPluginModel, ThreeDiResultItem, ThreeDiGridItem - from typing import List import logging @@ -843,11 +845,12 @@ def _get_active_parameter_config(self, result_item_ignored: ThreeDiResultItem = available_subgrid_vars = threedi_result.available_subgrid_map_vars available_agg_vars = threedi_result.available_aggregation_vars[:] # a copy available_wq_vars = threedi_result.available_water_quality_vars[:] # a copy + available_sca_vars = threedi_result.available_structure_control_actions_vars[:] # a copy if not available_agg_vars: messagebar_message("Warning", "No aggregation netCDF was found.", level=1, duration=5) parameter_config = generate_parameter_config( - available_subgrid_vars, agg_vars=available_agg_vars, wq_vars=available_wq_vars + available_subgrid_vars, agg_vars=available_agg_vars, wq_vars=available_wq_vars, sca_vars=available_sca_vars ) def _union(a: List, b: List): diff --git a/tool_statistics/presets.py b/tool_statistics/presets.py index 35863781..e8fd6918 100644 --- a/tool_statistics/presets.py +++ b/tool_statistics/presets.py @@ -2,8 +2,6 @@ from threedi_results_analysis.utils.threedi_result_aggregation.constants import ( AGGREGATION_VARIABLES, AGGREGATION_METHODS, - THRESHOLD_DRAIN_LEVEL, - THRESHOLD_EXCHANGE_LEVEL, ) from .style import ( Style, @@ -18,6 +16,8 @@ STYLE_MANHOLE_WATER_DEPTH_1D2D_NODE, STYLE_MANHOLE_MIN_FREEBOARD_0D1D, STYLE_MANHOLE_MIN_FREEBOARD_1D2D, + STYLE_SINGLE_COLUMN_GRADUATED_PUMP, + STYLE_SINGLE_COLUMN_GRADUATED_PUMP_LINESTRING, ) @@ -28,16 +28,26 @@ def __init__( description: str = "", aggregations=None, resample_point_layer: bool = False, + flowlines_style: Style = None, cells_style: Style = None, nodes_style: Style = None, + pumps_style: Style = None, + pumps_linestring_style: Style = None, + flowlines_style_param_values: dict = None, cells_style_param_values: dict = None, nodes_style_param_values: dict = None, + pumps_style_param_values: dict = None, + pumps_linestring_style_param_values: dict = None, + flowlines_layer_name: str = None, cells_layer_name: str = None, nodes_layer_name: str = None, + pumps_layer_name: str = None, + pumps_linestring_layer_name: str = None, raster_layer_name: str = None, + only_manholes: bool = False, ): if aggregations is None: @@ -46,16 +56,26 @@ def __init__( self.description = description self.__aggregations = aggregations self.resample_point_layer = resample_point_layer + self.flowlines_style = flowlines_style self.cells_style = cells_style self.nodes_style = nodes_style + self.pumps_style = pumps_style + self.pumps_linestring_style = pumps_linestring_style + self.flowlines_style_param_values = flowlines_style_param_values self.cells_style_param_values = cells_style_param_values self.nodes_style_param_values = nodes_style_param_values + self.pumps_style_param_values = pumps_style_param_values + self.pumps_linestring_style_param_values = pumps_linestring_style_param_values + self.flowlines_layer_name = flowlines_layer_name self.cells_layer_name = cells_layer_name self.nodes_layer_name = nodes_layer_name + self.pumps_layer_name = pumps_layer_name + self.pumps_linestring_layer_name = pumps_linestring_layer_name self.raster_layer_name = raster_layer_name + self.only_manholes = only_manholes def add_aggregation(self, aggregation: Aggregation): @@ -260,7 +280,7 @@ def aggregations(self): Aggregation( variable=AGGREGATION_VARIABLES.get_by_short_name("s1"), method=AGGREGATION_METHODS.get_by_short_name("time_above_threshold"), - threshold=THRESHOLD_DRAIN_LEVEL, + threshold="drain_level", ), ] @@ -268,7 +288,7 @@ def aggregations(self): Aggregation( variable=AGGREGATION_VARIABLES.get_by_short_name("s1"), method=AGGREGATION_METHODS.get_by_short_name("time_above_threshold"), - threshold=THRESHOLD_EXCHANGE_LEVEL, + threshold="exchange_level_1d2d", ), ] @@ -301,7 +321,7 @@ def aggregations(self): "connection to the 2D domain, so the 'water on street duration' will be 0 for all manholes.", aggregations=water_on_street_aggregations_1d2d, nodes_style=STYLE_WATER_ON_STREET_DURATION_NODE, - nodes_style_param_values={"column": "s1_time_above_threshold_exchange_level"}, + nodes_style_param_values={"column": "s1_time_above_threshold_exchange_level_1d2d"}, nodes_layer_name="Manhole: Water on street duration (1D2D)", only_manholes=True, ) @@ -393,6 +413,51 @@ def aggregations(self): only_manholes=True ) +# Pump: Total pumped volume +total_pumped_volume_aggregations = [ + Aggregation( + variable=AGGREGATION_VARIABLES.get_by_short_name("q_pump"), + method=AGGREGATION_METHODS.get_by_short_name("sum") + ), +] + +TOTAL_PUMPED_VOLUME_PRESETS = Preset( + name="Pump: Total pumped volume", + description="Total volume pumped by each pump in the selected time period.", + aggregations=total_pumped_volume_aggregations, + + pumps_style=STYLE_SINGLE_COLUMN_GRADUATED_PUMP, + pumps_style_param_values={"column": "q_pump_sum"}, + pumps_layer_name="Pump (point): Total pumped volume [m3]", + + pumps_linestring_style=STYLE_SINGLE_COLUMN_GRADUATED_PUMP_LINESTRING, + pumps_linestring_style_param_values={"column": "q_pump_sum"}, + pumps_linestring_layer_name="Pump (line): Total pumped volume [m3]", +) + +# Pump: time at max capacity +pump_time_at_max_capacity_aggregations = [ + Aggregation( + variable=AGGREGATION_VARIABLES.get_by_short_name("q_pump"), + method=AGGREGATION_METHODS.get_by_short_name("on_thres"), + threshold="capacity" + ), +] + +PUMP_TIME_AT_MAX_CAPACITY_PRESETS = Preset( + name="Pump: % of time at max capacity", + description="Percentage of time that each pump is pumping at its maximum capacity in the selected time period.\n\n" + "Note that both the pump implicit factor and the output time step will affect the result.", + aggregations=pump_time_at_max_capacity_aggregations, + + pumps_style=STYLE_SINGLE_COLUMN_GRADUATED_PUMP, + pumps_style_param_values={"column": "q_pump_on_thres_capacity"}, + pumps_layer_name="Pump (point): % of time at max capacity", + + pumps_linestring_style=STYLE_SINGLE_COLUMN_GRADUATED_PUMP_LINESTRING, + pumps_linestring_style_param_values={"column": "q_pump_on_thres_capacity"}, + pumps_linestring_layer_name="Pump (line): % of time at max capacity", +) PRESETS = [ NO_PRESET, @@ -407,4 +472,6 @@ def aggregations(self): MIN_FREEBOARD_1D2D_PRESETS, WATER_ON_STREET_DURATION_0D1D_PRESET, WATER_ON_STREET_DURATION_1D2D_PRESET, + TOTAL_PUMPED_VOLUME_PRESETS, + PUMP_TIME_AT_MAX_CAPACITY_PRESETS, ] diff --git a/tool_statistics/statistics.py b/tool_statistics/statistics.py index 47af4ba7..9b2cecac 100644 --- a/tool_statistics/statistics.py +++ b/tool_statistics/statistics.py @@ -24,13 +24,6 @@ * * ***************************************************************************/ """ -from typing import List - -from osgeo.gdal import GetDriverByName -from qgis.core import Qgis, QgsApplication, QgsProject, QgsTask, QgsRasterLayer -from threedi_results_analysis.utils.threedi_result_aggregation.base import aggregate_threedi_results -from threedi_results_analysis.utils.ogr2qgis import as_qgis_memory_layer - # Import the code for the dialog from .threedi_custom_stats_dialog import ThreeDiCustomStatsDialog from threedi_results_analysis.threedi_plugin_tool import ThreeDiPluginTool @@ -46,236 +39,11 @@ # TODO: cfl strictness factors instelbaar maken # TODO: berekening van max timestep ook op basis van volume vs. debiet -# TODO: opties af laten hangen van wat er in het model aanwezig is; is wel tricky ivm presets -GROUP_NAME = "Result aggregation outputs" - - -class Aggregate3DiResults(QgsTask): - def __init__( - self, - description: str, - parent: ThreeDiCustomStatsDialog, - layer_groups, - result: ThreeDiResultItem, - demanded_aggregations: List, - bbox, - start_time: int, - end_time: int, - only_manholes: bool, - interpolation_method, - resample_point_layer: bool, - resolution, - output_flowlines: bool, - output_cells: bool, - output_nodes: bool, - output_rasters: bool, - ): - super().__init__(description, QgsTask.CanCancel) - self.exception = None - self.parent = parent - self.parent.setEnabled(False) - self.result = result - self.layer_groups = layer_groups - self.demanded_aggregations = demanded_aggregations - self.bbox = bbox - self.start_time = start_time - self.end_time = end_time - self.only_manholes = only_manholes - self.interpolation_method = interpolation_method - self.resample_point_layer = resample_point_layer - self.resolution = resolution - self.output_flowlines = output_flowlines - self.output_cells = output_cells - self.output_nodes = output_nodes - self.output_rasters = output_rasters - - self.parent.iface.messageBar().pushMessage( - "3Di Statistics", - "Started aggregating 3Di results", - level=Qgis.Info, - duration=3, - ) - self.parent.iface.mainWindow().repaint() # to show the message before the task starts - - def run(self): - grid_admin = str(self.result.parent().path.with_suffix('.h5')) - grid_admin_gpkg = str(self.result.parent().path.with_suffix('.gpkg')) - results_3di = str(self.result.path) - - try: - self.ogr_ds, self.mem_rasts = aggregate_threedi_results( - gridadmin=grid_admin, - gridadmin_gpkg=grid_admin_gpkg, - results_3di=results_3di, - demanded_aggregations=self.demanded_aggregations, - bbox=self.bbox, - start_time=self.start_time, - end_time=self.end_time, - only_manholes=self.only_manholes, - interpolation_method=self.interpolation_method, - resample_point_layer=self.resample_point_layer, - resolution=self.resolution, - output_flowlines=self.output_flowlines, - output_cells=self.output_cells, - output_nodes=self.output_nodes, - output_rasters=self.output_rasters, - ) - - return True - - except Exception as e: - self.exception = e - - return False - - def _get_or_create_result_group(self, result: ThreeDiResultItem, group_name: str): - # We'll place the result layers in a dedicated result group - grid_item = result.parent() - assert grid_item - tool_group = grid_item.layer_group.findGroup(group_name) - if not tool_group: - tool_group = grid_item.layer_group.insertGroup(0, group_name) - tool_group.willRemoveChildren.connect(lambda n, i1, i2: self._group_removed(n, i1, i2)) - - # Add result group - result_group = tool_group.findGroup(result.text()) - if not result_group: - result_group = tool_group.addGroup(result.text()) - self.layer_groups[result.id] = result_group - # Use to modify result name when QgsLayerTreeNode is renamed. Note that this does not cause a - # infinite signal loop because the model only emits the result_changed when the text has actually - # changed. - result_group.nameChanged.connect(lambda _, txt, result_item=result: result_item.setText(txt)) - - return result_group - - def _group_removed(self, n, idxFrom, idxTo): - for result_id in list(self.layer_groups): - group = self.layer_groups[result_id] - for i in range(idxFrom, idxTo+1): - if n.children()[i] is group: - del self.layer_groups[result_id] - - def finished(self, result): - if self.exception is not None: - self.parent.setEnabled(True) - self.parent.repaint() - raise self.exception - if result: - # Add layers to layer tree - # They are added in order so the raster is below the polygon is below the line is below the point layer - - # raster layer - if len(self.mem_rasts) > 0: - for rastname, rast in self.mem_rasts.items(): - raster_output_dir = ( - self.parent.mQgsFileWidgetRasterFolder.filePath() - ) - raster_output_fn = os.path.join( - raster_output_dir, rastname + ".tif" - ) - drv = GetDriverByName("GTiff") - drv.CreateCopy( - utf8_path=raster_output_fn, src=rast - ) - layer_name = self.parent.lineEditOutputRasterLayer.text() + f": {rastname}" - raster_layer = QgsRasterLayer(raster_output_fn, layer_name if layer_name else f"Aggregation results: raster {rastname}") - result_group = self._get_or_create_result_group(self.result, GROUP_NAME) - QgsProject.instance().addMapLayer(raster_layer, addToLegend=False) - result_group.insertLayer(0, raster_layer) - - # cell layer - ogr_lyr = self.ogr_ds.GetLayerByName("cell") - if ogr_lyr is not None: - if ogr_lyr.GetFeatureCount() > 0: - layer_name = self.parent.lineEditOutputCellLayer.text() - qgs_lyr = as_qgis_memory_layer(ogr_lyr, layer_name if layer_name else "Aggregation results: cells") - result_group = self._get_or_create_result_group(self.result, GROUP_NAME) - QgsProject.instance().addMapLayer(qgs_lyr, addToLegend=False) - result_group.insertLayer(0, qgs_lyr) - - style = self.parent.comboBoxCellsStyleType.currentData() - style_kwargs = self.parent.get_styling_parameters( - output_type=style.output_type - ) - style.apply(qgis_layer=qgs_lyr, style_kwargs=style_kwargs) - - # flowline layer - ogr_lyr = self.ogr_ds.GetLayerByName("flowline") - if ogr_lyr is not None: - if ogr_lyr.GetFeatureCount() > 0: - layer_name = self.parent.lineEditOutputFlowLayer.text() - qgs_lyr = as_qgis_memory_layer(ogr_lyr, layer_name if layer_name else "Aggregation results: flowlines") - result_group = self._get_or_create_result_group(self.result, GROUP_NAME) - QgsProject.instance().addMapLayer(qgs_lyr, addToLegend=False) - result_group.insertLayer(0, qgs_lyr) - style = ( - self.parent.comboBoxFlowlinesStyleType.currentData() - ) - style_kwargs = self.parent.get_styling_parameters( - output_type=style.output_type - ) - style.apply(qgis_layer=qgs_lyr, style_kwargs=style_kwargs) - - # node layer - ogr_lyr = self.ogr_ds.GetLayerByName("node") - if ogr_lyr is not None: - if ogr_lyr.GetFeatureCount() > 0: - layer_name = self.parent.lineEditOutputNodeLayer.text() - qgs_lyr = as_qgis_memory_layer(ogr_lyr, layer_name if layer_name else "Aggregation results: nodes") - result_group = self._get_or_create_result_group(self.result, GROUP_NAME) - QgsProject.instance().addMapLayer(qgs_lyr, addToLegend=False) - result_group.insertLayer(0, qgs_lyr) - style = self.parent.comboBoxNodesStyleType.currentData() - style_kwargs = self.parent.get_styling_parameters( - output_type=style.output_type - ) - style.apply(qgis_layer=qgs_lyr, style_kwargs=style_kwargs) - - # resampled point layer - ogr_lyr = self.ogr_ds.GetLayerByName("node_resampled") - if ogr_lyr is not None: - if ogr_lyr.GetFeatureCount() > 0: - layer_name = self.parent.lineEditOutputNodeLayer.text() - qgs_lyr = as_qgis_memory_layer(ogr_lyr, (layer_name + "_resampled_nodes") if layer_name else "Aggregation results: resampled nodes") - result_group = self._get_or_create_result_group(self.result, GROUP_NAME) - QgsProject.instance().addMapLayer(qgs_lyr, addToLegend=False) - result_group.insertLayer(0, qgs_lyr) - style = self.parent.comboBoxNodesStyleType.currentData() - style_kwargs = self.parent.get_styling_parameters( - output_type=style.output_type - ) - style.apply(qgis_layer=qgs_lyr, style_kwargs=style_kwargs) - - self.parent.setEnabled(True) - self.parent.iface.messageBar().pushMessage( - "3Di Result aggregation", - "Finished custom aggregation", - level=Qgis.Success, - duration=3, - ) - - else: - self.parent.setEnabled(True) - self.parent.iface.messageBar().pushMessage( - "3Di Result aggregation", - "Aggregating 3Di results returned no results", - level=Qgis.Warning, - duration=3, - ) - - def cancel(self): - self.parent.iface.messageBar().pushMessage( - "3Di Result aggregation", - "Pre-processing simulation results cancelled by user", - level=Qgis.Info, - duration=3, - ) - super().cancel() - class StatisticsTool(ThreeDiPluginTool): + group_name = "Result aggregation outputs" + def __init__(self, iface, model): super().__init__() @@ -291,8 +59,6 @@ def __init__(self, iface, model): # Check if plugin was started the first time in current QGIS session self.first_start = True - self.tm = QgsApplication.taskManager() - def read(self, _) -> bool: """A new project is loaded, see if we can fetch some precreated groups""" return self._collect_result_groups() @@ -305,7 +71,7 @@ def _collect_result_groups(self): for result in results: grid_item = result.parent() assert grid_item - tool_group = grid_item.layer_group.findGroup(GROUP_NAME) + tool_group = grid_item.layer_group.findGroup(self.group_name) if tool_group: tool_group.willRemoveChildren.connect(lambda n, i1, i2: self._group_removed(n, i1, i2)) result_group = tool_group.findGroup(result.text()) @@ -328,70 +94,14 @@ def run(self): if self.first_start: self._collect_result_groups() self.first_start = False - self.dlg = ThreeDiCustomStatsDialog(self.iface, self.model) + self.dlg = ThreeDiCustomStatsDialog( + iface=self.iface, + model=self.model, + owner=self, + ) # show the dialog self.dlg.show() - # Run the dialog event loop - result = self.dlg.exec_() - # See if OK was pressed - if result: - # 3Di results - result = self.model.get_result(self.dlg.result_id) - - # Filtering parameters - start_time = self.dlg.doubleSpinBoxStartTime.value() - end_time = self.dlg.doubleSpinBoxEndTime.value() - bbox_qgs_rectangle = ( - self.dlg.mExtentGroupBox.outputExtent() - ) # bbox is now a https://qgis.org/pyqgis/master/core/QgsRectangle.html#qgis.core.QgsRectangle - - bbox = None - if bbox_qgs_rectangle is not None: - if not bbox_qgs_rectangle.isEmpty(): - bbox = [ - bbox_qgs_rectangle.xMinimum(), - bbox_qgs_rectangle.yMinimum(), - bbox_qgs_rectangle.xMaximum(), - bbox_qgs_rectangle.yMaximum(), - ] - only_manholes = self.dlg.onlyManholeCheckBox.isChecked() - - # Resolution - resolution = self.dlg.doubleSpinBoxResolution.value() - - # Outputs - output_flowlines = self.dlg.groupBoxFlowlines.isChecked() - output_nodes = self.dlg.groupBoxNodes.isChecked() - output_cells = self.dlg.groupBoxCells.isChecked() - output_rasters = self.dlg.groupBoxRasters.isChecked() - - # Resample point layer - resample_point_layer = self.dlg.checkBoxResample.isChecked() - if resample_point_layer: - interpolation_method = "linear" - else: - interpolation_method = None - - aggregate_threedi_results_task = Aggregate3DiResults( - description="Aggregate 3Di Results", - parent=self.dlg, - layer_groups=self.layer_groups, - result=result, - demanded_aggregations=self.dlg.demanded_aggregations, - bbox=bbox, - start_time=start_time, - end_time=end_time, - only_manholes=only_manholes, - interpolation_method=interpolation_method, - resample_point_layer=resample_point_layer, - resolution=resolution, - output_flowlines=output_flowlines, - output_cells=output_cells, - output_nodes=output_nodes, - output_rasters=output_rasters, - ) - self.tm.addTask(aggregate_threedi_results_task) @pyqtSlot(ThreeDiResultItem) def result_added(self, result_item: ThreeDiResultItem) -> None: @@ -417,7 +127,7 @@ def result_removed(self, result_item: ThreeDiResultItem) -> None: tool_group.removeChildNode(result_group) # In case the tool ("statistics") group is now empty, we'll remove that too - tool_group = result_item.parent().layer_group.findGroup(GROUP_NAME) + tool_group = result_item.parent().layer_group.findGroup(self.group_name) if len(tool_group.children()) == 0: tool_group.parent().removeChildNode(tool_group) diff --git a/tool_statistics/style.py b/tool_statistics/style.py index aef825ee..b5aa9c50 100644 --- a/tool_statistics/style.py +++ b/tool_statistics/style.py @@ -24,7 +24,7 @@ def __init__( styling_method, ): self.name = name - assert output_type in ("flowline", "node", "cell", "raster") + assert output_type in ("flowline", "node", "cell", "pump", "pump_linestring", "raster") self.output_type = output_type self.params = params if os.path.isabs(qml): @@ -48,6 +48,13 @@ def style_on_single_column(layer, qml: str, column: str, update_classes: bool = mode=layer.renderer().mode(), nclasses=len(layer.renderer().ranges()), ) + + # Add a class for 0 if the lowest value is 0 + range_0 = layer.renderer().ranges()[0] + if range_0.lowerValue() == 0 and range_0.upperValue() > 0.000001: + layer.renderer().addBreak(breakValue=0.000001, updateSymbols=True) + layer.renderer().updateRangeLabel(rangeIndex=0, label="0") + layer.triggerRepaint() utils.iface.layerTreeView().refreshLayerSymbology(layer.id()) @@ -407,6 +414,23 @@ def style_ts_reduction_analysis( styling_method=style_balance, ) +STYLE_SINGLE_COLUMN_GRADUATED_PUMP = Style( + name="Single column graduated", + output_type="pump", + params={"column": "column"}, + qml="pump.qml", + styling_method=style_on_single_column, +) + +STYLE_SINGLE_COLUMN_GRADUATED_PUMP_LINESTRING = Style( + name="Single column graduated", + output_type="pump_linestring", + params={"column": "column"}, + qml="pump_linestring.qml", + styling_method=style_on_single_column, +) + + STYLES = [ STYLE_FLOW_DIRECTION, STYLE_GRADIENT, @@ -421,9 +445,12 @@ def style_ts_reduction_analysis( STYLE_MANHOLE_WATER_DEPTH_0D1D_NODE, STYLE_MANHOLE_WATER_DEPTH_1D2D_NODE, STYLE_MANHOLE_MIN_FREEBOARD_0D1D, - STYLE_MANHOLE_MIN_FREEBOARD_1D2D + STYLE_MANHOLE_MIN_FREEBOARD_1D2D, + STYLE_SINGLE_COLUMN_GRADUATED_PUMP, + STYLE_SINGLE_COLUMN_GRADUATED_PUMP_LINESTRING, ] + DEFAULT_STYLES = { # Flowlines "q": {"flowline": STYLE_FLOW_DIRECTION}, @@ -436,6 +463,12 @@ def style_ts_reduction_analysis( "bed_grad": {"flowline": STYLE_GRADIENT}, "wl_at_xsec": {"flowline": STYLE_SINGLE_COLUMN_GRADUATED_FLOWLINE}, + # Pumps + "q_pump": { + "pump": STYLE_SINGLE_COLUMN_GRADUATED_PUMP, + "pump_linestring": STYLE_SINGLE_COLUMN_GRADUATED_PUMP_LINESTRING, + }, + # Nodes "s1": { "node": STYLE_SINGLE_COLUMN_GRADUATED_NODE, diff --git a/tool_statistics/style/pump.qml b/tool_statistics/style/pump.qml new file mode 100644 index 00000000..80d15313 --- /dev/null +++ b/tool_statistics/style/pump.qml @@ -0,0 +1,584 @@ + +