From 08f136a186143a6810cc4b1a666b0db577948773 Mon Sep 17 00:00:00 2001 From: Golnesa Date: Wed, 20 Dec 2023 19:22:53 +0100 Subject: [PATCH 01/12] adding yz_coordinate --- .../grid/cross_section_definitions.py | 39 ++++++++++++++----- threedigrid_builder/interface/gridadmin.py | 1 + .../tests/test_cross_section_definitions.py | 3 +- 3 files changed, 32 insertions(+), 11 deletions(-) diff --git a/threedigrid_builder/grid/cross_section_definitions.py b/threedigrid_builder/grid/cross_section_definitions.py index 397ea9c4..611f8fd1 100644 --- a/threedigrid_builder/grid/cross_section_definitions.py +++ b/threedigrid_builder/grid/cross_section_definitions.py @@ -43,12 +43,14 @@ def convert(self, ids): content_pk=ids, code=self.code[idx], count=0, + count_yz = 0, ) if len(result) == 0: return result - offset = 0 tables = [] + yz_coordinate = [] + # Numpy array views on width/height based on idx width_idx = self.width[idx] @@ -56,23 +58,32 @@ def convert(self, ids): for i, shape in enumerate(self.shape[idx]): tabulator = tabulators[shape] - result.shape[i], result.width_1d[i], result.height_1d[i], table = tabulator( + result.shape[i], result.width_1d[i], result.height_1d[i], table, yz = tabulator( shape, width_idx[i], height_idx[i] ) if table is not None: result.count[i] = len(table) - result.offset[i] = offset - offset += len(table) tables.append(table) + if yz is not None: + result.count_yz[i] = len(yz) + yz_coordinate.append(yz) + result.offset[:] = np.roll(np.cumsum(result.count), 1) result.offset[0] = 0 + result.offset_yz[:] = np.roll(np.cumsum(result.count_yz), 1) + result.offset_yz[0] = 0 if len(tables) > 0: result.tables = np.concatenate(tables, axis=0) else: result.tables = np.empty((0, 2)) + if len(yz_coordinate) > 0: + result.yz_coordinate = np.concatenate(yz_coordinate, axis=0) + else: + result.yz_coordinate = np.empty((0, 4)) + return result @@ -85,11 +96,14 @@ class CrossSection: height_1d: float offset: int count: int + offset_yz: int + count_yz: int # tables: Tuple[float, float] has different length so is specified on CrossSections class CrossSections(Array[CrossSection]): tables = None + yz_coordinate = None def tabulate_builtin(shape, width, height): @@ -112,7 +126,7 @@ def tabulate_builtin(shape, width, height): f"Unable to parse cross section definition width (got: '{width}')." ) - return shape, width, None, None + return shape, width, None, None, None def tabulate_egg(shape, width, height): @@ -151,7 +165,7 @@ def tabulate_egg(shape, width, height): widths = np.sqrt(p / q) * 2 table = np.array([heights, widths]).T - return CrossSectionShape.TABULATED_TRAPEZIUM, width, height, table + return CrossSectionShape.TABULATED_TRAPEZIUM, width, height, table, None def tabulate_inverted_egg(shape, width, height): @@ -161,7 +175,7 @@ def tabulate_inverted_egg(shape, width, height): """ type_, width, height, table = tabulate_egg(shape, width, height) table[:, 1] = table[::-1, 1] - return type_, width, height, table + return type_, width, height, table, None def tabulate_closed_rectangle(shape, width, height): @@ -185,7 +199,7 @@ def tabulate_closed_rectangle(shape, width, height): f"(got: '{width}', '{height}')." ) table = np.array([[0.0, width], [height, 0.0]], order="F") - return CrossSectionShape.TABULATED_RECTANGLE, width, height, table + return CrossSectionShape.TABULATED_RECTANGLE, width, height, table, None def _parse_tabulated(width, height): @@ -229,7 +243,7 @@ def tabulate_tabulated(shape, width, height): f"(got: {height})." ) - return shape, np.max(widths), np.max(heights), np.array([heights, widths]).T + return shape, np.max(widths), np.max(heights), np.array([heights, widths]).T, None def tabulate_yz(shape, width, height): @@ -275,6 +289,11 @@ def tabulate_yz(shape, width, height): if is_closed: ys = ys[:-1] zs = zs[:-1] + yz = None + else: + yz = np.empty((len(ys), 4), dtype=float) + yz[:, 0] = ys + yz[:, 1] = zs # Adapt non-unique height coordinates. Why? # Because if a segment of the profile is exactly horizontal, we need 2 widths @@ -318,7 +337,7 @@ def tabulate_yz(shape, width, height): table = table[1:] height_1d, width_1d = table.max(axis=0).tolist() - return CrossSectionShape.TABULATED_TRAPEZIUM, width_1d, height_1d, table + return CrossSectionShape.TABULATED_TRAPEZIUM, width_1d, height_1d, table, yz tabulators = { diff --git a/threedigrid_builder/interface/gridadmin.py b/threedigrid_builder/interface/gridadmin.py index cfd5d12e..1889abe6 100644 --- a/threedigrid_builder/interface/gridadmin.py +++ b/threedigrid_builder/interface/gridadmin.py @@ -690,6 +690,7 @@ def write_cross_sections(self, cross_sections: CrossSections): self.write_dataset(group, "offset", cross_sections.offset) self.write_dataset(group, "count", cross_sections.count) self.write_dataset(group, "tables", cross_sections.tables.T, insert_dummy=False) + self.write_dataset(group, "yz", cross_sections.yz_coordinate.T, insert_dummy=False) def write_obstacles(self, obstacles: Obstacles): """For backwards compat, the group is named 'levees'""" diff --git a/threedigrid_builder/tests/test_cross_section_definitions.py b/threedigrid_builder/tests/test_cross_section_definitions.py index 9f11fd9c..2cfceed3 100644 --- a/threedigrid_builder/tests/test_cross_section_definitions.py +++ b/threedigrid_builder/tests/test_cross_section_definitions.py @@ -156,13 +156,14 @@ def test_tabulate_tabulated_err(width, height, match): def test_tabulate_inverted_egg(): - shape, width_1d, height_1d, table = tabulate_inverted_egg( + shape, width_1d, height_1d, table, yz = tabulate_inverted_egg( "my-shape", "1.52", "ignored" ) assert shape == CrossSectionShape.TABULATED_TRAPEZIUM assert width_1d == 1.52 assert height_1d == 1.52 * 1.5 + assert yz is None # the expected table is exactly what inpy returns for a width of 1.52 expected_table = np.array( From b9fdf4b881c0b3726d4b57fb2c1233e2666cf6d4 Mon Sep 17 00:00:00 2001 From: Golnesa Date: Mon, 8 Jan 2024 14:15:49 +0100 Subject: [PATCH 02/12] Added vegetation to loc and def, tests updated --- threedigrid_builder/base/lines.py | 2 + .../grid/cross_section_definitions.py | 32 ++++++++--- .../grid/cross_section_locations.py | 6 ++ threedigrid_builder/interface/db.py | 9 +++ threedigrid_builder/interface/gridadmin.py | 7 ++- .../tests/test_cross_section_definitions.py | 55 +++++++++++++------ .../tests/test_cross_section_locations.py | 6 ++ threedigrid_builder/tests/test_gridadmin.py | 1 + 8 files changed, 92 insertions(+), 26 deletions(-) diff --git a/threedigrid_builder/base/lines.py b/threedigrid_builder/base/lines.py index e4b90988..d52aa25a 100644 --- a/threedigrid_builder/base/lines.py +++ b/threedigrid_builder/base/lines.py @@ -38,6 +38,8 @@ class Line: frict_type2: int frict_value1: float # For 1D-2Dgw: resistance factor-out (hydr. cond / thickness) frict_value2: float # For 1D-2Dgw: resistance factor-in (hydr. cond / thickness) + veg_coef1: float + veg_coef2: float cross_weight: float # For 1D-2Dgw: the exchange length (tables: cross_width) discharge_coefficient_positive: float discharge_coefficient_negative: float diff --git a/threedigrid_builder/grid/cross_section_definitions.py b/threedigrid_builder/grid/cross_section_definitions.py index 611f8fd1..bf276c0e 100644 --- a/threedigrid_builder/grid/cross_section_definitions.py +++ b/threedigrid_builder/grid/cross_section_definitions.py @@ -17,6 +17,12 @@ class CrossSectionDefinition: shape: CrossSectionShape height: str # space-separated list of floats width: str # space-separated list of floats + friction_values: str # space-separated list of floats + vegetation_stem_densities: str # space-separated list of floats + vegetation_stem_diameters: str # space-separated list of floats + vegetation_heights: str # space-separated list of floats + vegetation_drag_coefficients: str # space-separated list of floats + class CrossSectionDefinitions(Array[CrossSectionDefinition]): @@ -49,7 +55,7 @@ def convert(self, ids): return result tables = [] - yz_coordinate = [] + tables_yz = [] # Numpy array views on width/height based on idx @@ -66,7 +72,7 @@ def convert(self, ids): tables.append(table) if yz is not None: result.count_yz[i] = len(yz) - yz_coordinate.append(yz) + tables_yz.append(yz) result.offset[:] = np.roll(np.cumsum(result.count), 1) @@ -79,10 +85,10 @@ def convert(self, ids): else: result.tables = np.empty((0, 2)) - if len(yz_coordinate) > 0: - result.yz_coordinate = np.concatenate(yz_coordinate, axis=0) + if len(tables_yz) > 0: + result.tables_yz = np.concatenate(tables_yz, axis=0) else: - result.yz_coordinate = np.empty((0, 4)) + result.tables_yz = np.empty((0, 4)) return result @@ -103,7 +109,7 @@ class CrossSection: class CrossSections(Array[CrossSection]): tables = None - yz_coordinate = None + tables_yz = None def tabulate_builtin(shape, width, height): @@ -173,7 +179,7 @@ def tabulate_inverted_egg(shape, width, height): See tabulate_egg. """ - type_, width, height, table = tabulate_egg(shape, width, height) + type_, width, height, table, yz = tabulate_egg(shape, width, height) table[:, 1] = table[::-1, 1] return type_, width, height, table, None @@ -246,7 +252,7 @@ def tabulate_tabulated(shape, width, height): return shape, np.max(widths), np.max(heights), np.array([heights, widths]).T, None -def tabulate_yz(shape, width, height): +def tabulate_yz(shape, width, height, friction_values, vegetation_stem_densities, vegetation_stem_diameters, vegetation_heights, vegetation_drag_coefficients): """Tabulate an (open or closed) YZ profile Args: @@ -291,9 +297,17 @@ def tabulate_yz(shape, width, height): zs = zs[:-1] yz = None else: - yz = np.empty((len(ys), 4), dtype=float) + yz = np.zeros((len(ys), 4), dtype=float) yz[:, 0] = ys yz[:, 1] = zs + fric = np.array([float(x) for x in friction_values.split(" ")]) + veg_stemden = np.array([float(x) for x in vegetation_stem_densities.split(" ")]) + veg_stemdia = np.array([float(x) for x in vegetation_stem_diameters.split(" ")]) + veg_hght = np.array([float(x) for x in vegetation_heights.split(" ")]) + veg_drag = np.array([float(x) for x in vegetation_drag_coefficients.split(" ")]) + yz[:-1, 2] = fric + yz[:-1, 3] = veg_stemden * veg_stemdia * veg_hght * veg_drag + # Adapt non-unique height coordinates. Why? # Because if a segment of the profile is exactly horizontal, we need 2 widths diff --git a/threedigrid_builder/grid/cross_section_locations.py b/threedigrid_builder/grid/cross_section_locations.py index 7c9bd97d..ba9c895a 100644 --- a/threedigrid_builder/grid/cross_section_locations.py +++ b/threedigrid_builder/grid/cross_section_locations.py @@ -17,6 +17,10 @@ class CrossSectionLocation: bank_level: float friction_type: FrictionType friction_value: float + vegetation_stem_density: float + vegetation_stem_diameter: float + vegetation_height: float + vegetation_drag_coefficient: float class CrossSectionLocations(Array[CrossSectionLocation]): @@ -64,6 +68,8 @@ def apply_to_lines(self, lines, channels, extrapolate=False): lines.frict_type2 = self.friction_type[idx2] lines.frict_value1 = self.friction_value[idx1] lines.frict_value2 = self.friction_value[idx2] + lines.veg_coef1 = self.vegetation_stem_density[idx1] * self.vegetation_stem_diameter[idx1] * self.vegetation_height[idx1] * self.vegetation_drag_coefficient[idx1] + lines.veg_coef2 = self.vegetation_stem_density[idx2] * self.vegetation_stem_diameter[idx2] * self.vegetation_height[idx2] * self.vegetation_drag_coefficient[idx2] # Compute invert levels and start and end lines.invert_level_start_point = compute_bottom_level( diff --git a/threedigrid_builder/interface/db.py b/threedigrid_builder/interface/db.py index 9329530c..2d56f5be 100644 --- a/threedigrid_builder/interface/db.py +++ b/threedigrid_builder/interface/db.py @@ -468,6 +468,11 @@ def get_cross_section_definitions(self) -> CrossSectionDefinitions: models.CrossSectionDefinition.shape, models.CrossSectionDefinition.width, models.CrossSectionDefinition.height, + models.CrossSectionDefinition.friction_values, + models.CrossSectionDefinition.vegetation_stem_densities, + models.CrossSectionDefinition.vegetation_stem_diameters, + models.CrossSectionDefinition.vegetation_heights, + models.CrossSectionDefinition.vegetation_drag_coefficients, ) .order_by(models.CrossSectionDefinition.id) .as_structarray() @@ -497,6 +502,10 @@ def get_cross_section_locations(self) -> CrossSectionLocations: models.CrossSectionLocation.bank_level, models.CrossSectionLocation.friction_type, models.CrossSectionLocation.friction_value, + models.CrossSectionLocation.vegetation_stem_density, + models.CrossSectionLocation.vegetation_stem_diameter, + models.CrossSectionLocation.vegetation_height, + models.CrossSectionLocation.vegetation_drag_coefficient, ) .order_by(models.CrossSectionLocation.id) .as_structarray() diff --git a/threedigrid_builder/interface/gridadmin.py b/threedigrid_builder/interface/gridadmin.py index 1889abe6..250a849b 100644 --- a/threedigrid_builder/interface/gridadmin.py +++ b/threedigrid_builder/interface/gridadmin.py @@ -442,6 +442,8 @@ def write_lines(self, lines: Lines, cross_sections): self.write_dataset(group, "frict_type2", lines.frict_type2) self.write_dataset(group, "frict_value1", lines.frict_value1) self.write_dataset(group, "frict_value2", lines.frict_value2) + self.write_dataset(group, "veg_coef1", lines.veg_coef1) + self.write_dataset(group, "veg_coef2", lines.veg_coef2) self.write_dataset(group, "cross_weight", lines.cross_weight) self.write_dataset( group, "hydraulic_resistance_in", lines.hydraulic_resistance_in @@ -690,7 +692,10 @@ def write_cross_sections(self, cross_sections: CrossSections): self.write_dataset(group, "offset", cross_sections.offset) self.write_dataset(group, "count", cross_sections.count) self.write_dataset(group, "tables", cross_sections.tables.T, insert_dummy=False) - self.write_dataset(group, "yz", cross_sections.yz_coordinate.T, insert_dummy=False) + if cross_sections.tables_yz is not None: + self.write_dataset(group, "offset_yz", cross_sections.offset_yz) + self.write_dataset(group, "count_yz", cross_sections.count_yz) + self.write_dataset(group, "tables_yz", cross_sections.tables_yz.T, insert_dummy=False) def write_obstacles(self, obstacles: Obstacles): """For backwards compat, the group is named 'levees'""" diff --git a/threedigrid_builder/tests/test_cross_section_definitions.py b/threedigrid_builder/tests/test_cross_section_definitions.py index 2cfceed3..d2a9be79 100644 --- a/threedigrid_builder/tests/test_cross_section_definitions.py +++ b/threedigrid_builder/tests/test_cross_section_definitions.py @@ -1,3 +1,4 @@ +from cmath import nan from unittest import mock import numpy as np @@ -36,9 +37,9 @@ def test_convert_multiple(cross_section_definitions): with mock.patch.dict( "threedigrid_builder.grid.cross_section_definitions.tabulators", { - SHP.CIRCLE: mock.Mock(return_value=(1, 0.1, None, None)), - SHP.TABULATED_TRAPEZIUM: mock.Mock(return_value=(5, 15.0, 2.0, table_1)), - SHP.TABULATED_RECTANGLE: mock.Mock(return_value=(6, 11.0, 2.0, table_2)), + SHP.CIRCLE: mock.Mock(return_value=(1, 0.1, None, None, None)), + SHP.TABULATED_TRAPEZIUM: mock.Mock(return_value=(5, 15.0, 2.0, table_1, None)), + SHP.TABULATED_RECTANGLE: mock.Mock(return_value=(6, 11.0, 2.0, table_2, None)), }, ): actual = cross_section_definitions.convert([1, 3, 9]) @@ -59,7 +60,7 @@ def test_convert_multiple_filtered(cross_section_definitions): with mock.patch.dict( "threedigrid_builder.grid.cross_section_definitions.tabulators", { - SHP.CIRCLE: mock.Mock(return_value=(1, 0.1, 0.1, None)), + SHP.CIRCLE: mock.Mock(return_value=(1, 0.1, 0.1, None, None)), }, ): actual = cross_section_definitions.convert([1]) @@ -85,11 +86,11 @@ def test_convert_empty(cross_section_definitions): def test_tabulate_builtin(): actual = tabulate_builtin("my-shape", "1.52", "1.33") - assert actual == ("my-shape", 1.52, None, None) + assert actual == ("my-shape", 1.52, None, None, None) def test_tabulate_closed_rectangle(): - shape, width_1d, height_1d, table = tabulate_closed_rectangle( + shape, width_1d, height_1d, table, yz = tabulate_closed_rectangle( "my-shape", "1.52", "5.2" ) @@ -97,14 +98,16 @@ def test_tabulate_closed_rectangle(): assert width_1d == 1.52 assert height_1d == 5.2 assert_almost_equal(table, np.array([[0.0, 1.52], [5.2, 0.0]])) + assert yz is None def test_tabulate_egg(): - shape, width_1d, height_1d, table = tabulate_egg("my-shape", "1.52", "ignored") + shape, width_1d, height_1d, table, yz = tabulate_egg("my-shape", "1.52", "ignored") assert shape == CrossSectionShape.TABULATED_TRAPEZIUM assert width_1d == 1.52 assert height_1d == 1.52 * 1.5 + assert yz is None # the expected table is exactly what inpy returns for a width of 1.52 expected_table = np.array( @@ -131,13 +134,13 @@ def test_tabulate_egg(): def test_tabulate_tabulated(): - shape, width_1d, height_1d, table = tabulate_tabulated("my-shape", "1 2 3", "0 1 2") + shape, width_1d, height_1d, table, yz = tabulate_tabulated("my-shape", "1 2 3", "0 1 2") assert shape == "my-shape" assert width_1d == 3.0 # the max assert height_1d == 2.0 assert_almost_equal(table, np.array([[0, 1], [1, 2], [2, 3]], dtype=float)) - + assert yz is None @pytest.mark.parametrize( "width,height,match", @@ -190,41 +193,61 @@ def test_tabulate_inverted_egg(): @pytest.mark.parametrize( - "width,height,exp_width,exp_height,exp_table", + "width,height,friction_values,vegetation_stem_densities,vegetation_stem_diameters,vegetation_heights,vegetation_drag_coefficients,exp_width,exp_height,exp_table,exp_yz", [ - ("0 0.5 1 1.5", "0.5 0 0 0.5", 1.5, 0.5, [[0, 0.5], [0.5, 1.5]]), + ("0 0.5 1 1.5", "0.5 0 0 0.5", "1 1 1", "1 1 1", "1 1 1", "1 1 1", "1 1 1", 1.5, 0.5, [[0, 0.5], [0.5, 1.5]], [[0, 0.5, 1, 1], [0.5, 0, 1, 1], [1, 0, 1, 1], [1.5, 0.5, 0, 0]]), ( "0 0.5 1 1.5", "0.5 0 0 0.25", + "1 1 1", + "0.5 1 1", + "1 1 1", + "3 1 1", + "2 1 1", 1.5, 0.5, [[0, 0.5], [0.25, 1.25], [0.5, 1.5]], + [[0, 0.5, 1, 3], [0.5, 0, 1, 1], [1, 0, 1, 1], [1.5, 0.25, 0, 0]] ), ( "0 1 2 3 4 5", "1 0 0.5 0.5 0 1", + "1 1 1 1 1", + "1 1 1 1 1", + "1 1 1 1 1", + "1 1 1 1 1", + "1 1 1 1 1", 5, 1, [[0, 0], [0.5, 3], [0.5, 4], [1, 5]], + [[0, 1, 1, 1], [1, 0, 1, 1], [2, 0.5, 1, 1], [3, 0.5, 1, 1], [4, 0, 1, 1], [5, 1, 0, 0]] ), ( "0 1 2 2 0 0", "0.5 0 0.5 1.5 1.5 0.5", + "1 1 1 1", + "1 1 1 1", + "1 1 1 1", + "1 1 1 1", + "1 1 1 1", 2.0, 1.5, [[0, 0], [0.5, 2.0], [1.5, 2.0], [1.5, 0.0]], + None ), - ("0 0.5 0.75 1.0 1.5", "0.5 0 0 0 0.5", 1.5, 0.5, [[0, 0.5], [0.5, 1.5]]), - ("0 1 0 1 0", "0 1 1 0 0", 1, 1, [[0, 1], [0.5, 0], [1, 1], [1, 0]]), + ("0 0.5 0.75 1.0 1.5", "0.5 0 0 0 0.5", "1 1 1 1", "1 1 1 1", "1 1 1 1", "1 1 1 1", "1 1 1 1", 1.5, 0.5, [[0, 0.5], [0.5, 1.5]], [[0, 0.5, 1, 1], [0.5, 0, 1, 1], [0.75, 0, 1, 1], [1, 0, 1, 1], [1.5, 0.5, 0, 0]]), + ("0 1 0 1 0", "0 1 1 0 0", "1 1 1 1", "1 1 1 1", "1 1 1 1", "1 1 1 1", "1 1 1 1", 1, 1, [[0, 1], [0.5, 0], [1, 1], [1, 0]], None), ], ) -def test_tabulate_yz(width, height, exp_width, exp_height, exp_table): - shape, width_1d, height_1d, table = tabulate_yz("my-shape", width, height) +def test_tabulate_yz(width, height, friction_values,vegetation_stem_densities,vegetation_stem_diameters,vegetation_heights,vegetation_drag_coefficients,exp_width, exp_height, exp_table, exp_yz): + shape, width_1d, height_1d, table, yz = tabulate_yz("my-shape", width, height,friction_values,vegetation_stem_densities,vegetation_stem_diameters,vegetation_heights,vegetation_drag_coefficients) assert shape == CrossSectionShape.TABULATED_TRAPEZIUM assert width_1d == exp_width assert height_1d == exp_height assert_almost_equal(table, np.array(exp_table, dtype=float)) + if yz is not None: + assert_almost_equal(yz, np.array(exp_yz, dtype=float)) @pytest.mark.parametrize( @@ -248,4 +271,4 @@ def test_tabulate_yz(width, height, exp_width, exp_height, exp_table): ) def test_tabulate_yz_err(width, height, match): with pytest.raises(SchematisationError, match=match): - tabulate_yz("my-shape", width, height) + tabulate_yz("my-shape", width, height, None, None, None, None, None) diff --git a/threedigrid_builder/tests/test_cross_section_locations.py b/threedigrid_builder/tests/test_cross_section_locations.py index 32ec223a..30a583b7 100644 --- a/threedigrid_builder/tests/test_cross_section_locations.py +++ b/threedigrid_builder/tests/test_cross_section_locations.py @@ -62,6 +62,10 @@ def locations(): channel_id=[51, 53, 53, 52, 54], friction_type=[1, 1, 1, 2, 2], friction_value=[30, 35, 40, 0.02, 0.03], + vegetation_stem_density=[1, 0.3, 0.4, 0.02, 0.03], + vegetation_stem_diameter=[0.5, 0.5, 0.5, 0.3, 0.4], + vegetation_height=[0.2, 1, 1, 1, 3], + vegetation_drag_coefficient=[1, 1, 1, 1, 1], reference_level=[1.0, 2.0, 3.0, 5.0, 6.0], definition_id=[3, 3, 4, 4, 5], ) @@ -165,6 +169,8 @@ def test_apply_to_lines(channels, channel_lines, locations): assert_equal(channel_lines.frict_type2, [1, 2, 1, 1, 1, 2, 2]) assert_equal(channel_lines.frict_value1, [30, 0.02, 40, 40, 40, 0.03, 0.03]) assert_equal(channel_lines.frict_value2, [30, 0.02, 35, 35, 35, 0.03, 0.03]) + assert_almost_equal(channel_lines.veg_coef1, [0.1, 0.006, 0.2, 0.2, 0.2, 0.036, 0.036]) + assert_almost_equal(channel_lines.veg_coef2, [0.1, 0.006, 0.15, 0.15, 0.15, 0.036, 0.036]) assert_almost_equal( channel_lines.cross_weight, [1.0, 1.0, 1.0, 0.65, 0.0, 1.0, 1.0] ) diff --git a/threedigrid_builder/tests/test_gridadmin.py b/threedigrid_builder/tests/test_gridadmin.py index 052d828f..2f1dc334 100644 --- a/threedigrid_builder/tests/test_gridadmin.py +++ b/threedigrid_builder/tests/test_gridadmin.py @@ -140,6 +140,7 @@ def test_write_nodes_embedded(h5_out, dataset, shape, dtype): ("dist_calc_points", (6,), "float64"), ("friction_type", (6,), "int32"), ("friction_value", (6,), "float64"), + # ("vegetation_stem_density", (6,), "float64"), ("id", (6,), "int32"), ("invert_level_end_point", (6,), "float64"), ("invert_level_start_point", (6,), "float64"), From 322728fa5439ca2575cc3b2299612880052ebd5f Mon Sep 17 00:00:00 2001 From: Golnesa Date: Tue, 9 Jan 2024 10:27:51 +0100 Subject: [PATCH 03/12] Linters --- .../grid/cross_section_definitions.py | 29 +++-- .../grid/cross_section_locations.py | 14 ++- threedigrid_builder/interface/gridadmin.py | 4 +- .../tests/test_cross_section_definitions.py | 100 +++++++++++++++--- .../tests/test_cross_section_locations.py | 8 +- 5 files changed, 128 insertions(+), 27 deletions(-) diff --git a/threedigrid_builder/grid/cross_section_definitions.py b/threedigrid_builder/grid/cross_section_definitions.py index bf276c0e..79e1db94 100644 --- a/threedigrid_builder/grid/cross_section_definitions.py +++ b/threedigrid_builder/grid/cross_section_definitions.py @@ -24,7 +24,6 @@ class CrossSectionDefinition: vegetation_drag_coefficients: str # space-separated list of floats - class CrossSectionDefinitions(Array[CrossSectionDefinition]): def convert(self, ids): """Convert to CrossSections. @@ -49,7 +48,7 @@ def convert(self, ids): content_pk=ids, code=self.code[idx], count=0, - count_yz = 0, + count_yz=0, ) if len(result) == 0: return result @@ -57,16 +56,19 @@ def convert(self, ids): tables = [] tables_yz = [] - # Numpy array views on width/height based on idx width_idx = self.width[idx] height_idx = self.height[idx] for i, shape in enumerate(self.shape[idx]): tabulator = tabulators[shape] - result.shape[i], result.width_1d[i], result.height_1d[i], table, yz = tabulator( - shape, width_idx[i], height_idx[i] - ) + ( + result.shape[i], + result.width_1d[i], + result.height_1d[i], + table, + yz, + ) = tabulator(shape, width_idx[i], height_idx[i]) if table is not None: result.count[i] = len(table) tables.append(table) @@ -74,7 +76,6 @@ def convert(self, ids): result.count_yz[i] = len(yz) tables_yz.append(yz) - result.offset[:] = np.roll(np.cumsum(result.count), 1) result.offset[0] = 0 result.offset_yz[:] = np.roll(np.cumsum(result.count_yz), 1) @@ -88,7 +89,7 @@ def convert(self, ids): if len(tables_yz) > 0: result.tables_yz = np.concatenate(tables_yz, axis=0) else: - result.tables_yz = np.empty((0, 4)) + result.tables_yz = np.empty((0, 4)) return result @@ -252,7 +253,16 @@ def tabulate_tabulated(shape, width, height): return shape, np.max(widths), np.max(heights), np.array([heights, widths]).T, None -def tabulate_yz(shape, width, height, friction_values, vegetation_stem_densities, vegetation_stem_diameters, vegetation_heights, vegetation_drag_coefficients): +def tabulate_yz( + shape, + width, + height, + friction_values, + vegetation_stem_densities, + vegetation_stem_diameters, + vegetation_heights, + vegetation_drag_coefficients, +): """Tabulate an (open or closed) YZ profile Args: @@ -308,7 +318,6 @@ def tabulate_yz(shape, width, height, friction_values, vegetation_stem_densities yz[:-1, 2] = fric yz[:-1, 3] = veg_stemden * veg_stemdia * veg_hght * veg_drag - # Adapt non-unique height coordinates. Why? # Because if a segment of the profile is exactly horizontal, we need 2 widths seen = set() diff --git a/threedigrid_builder/grid/cross_section_locations.py b/threedigrid_builder/grid/cross_section_locations.py index ba9c895a..8313ad34 100644 --- a/threedigrid_builder/grid/cross_section_locations.py +++ b/threedigrid_builder/grid/cross_section_locations.py @@ -68,8 +68,18 @@ def apply_to_lines(self, lines, channels, extrapolate=False): lines.frict_type2 = self.friction_type[idx2] lines.frict_value1 = self.friction_value[idx1] lines.frict_value2 = self.friction_value[idx2] - lines.veg_coef1 = self.vegetation_stem_density[idx1] * self.vegetation_stem_diameter[idx1] * self.vegetation_height[idx1] * self.vegetation_drag_coefficient[idx1] - lines.veg_coef2 = self.vegetation_stem_density[idx2] * self.vegetation_stem_diameter[idx2] * self.vegetation_height[idx2] * self.vegetation_drag_coefficient[idx2] + lines.veg_coef1 = ( + self.vegetation_stem_density[idx1] + * self.vegetation_stem_diameter[idx1] + * self.vegetation_height[idx1] + * self.vegetation_drag_coefficient[idx1] + ) + lines.veg_coef2 = ( + self.vegetation_stem_density[idx2] + * self.vegetation_stem_diameter[idx2] + * self.vegetation_height[idx2] + * self.vegetation_drag_coefficient[idx2] + ) # Compute invert levels and start and end lines.invert_level_start_point = compute_bottom_level( diff --git a/threedigrid_builder/interface/gridadmin.py b/threedigrid_builder/interface/gridadmin.py index 250a849b..c4a43286 100644 --- a/threedigrid_builder/interface/gridadmin.py +++ b/threedigrid_builder/interface/gridadmin.py @@ -695,7 +695,9 @@ def write_cross_sections(self, cross_sections: CrossSections): if cross_sections.tables_yz is not None: self.write_dataset(group, "offset_yz", cross_sections.offset_yz) self.write_dataset(group, "count_yz", cross_sections.count_yz) - self.write_dataset(group, "tables_yz", cross_sections.tables_yz.T, insert_dummy=False) + self.write_dataset( + group, "tables_yz", cross_sections.tables_yz.T, insert_dummy=False + ) def write_obstacles(self, obstacles: Obstacles): """For backwards compat, the group is named 'levees'""" diff --git a/threedigrid_builder/tests/test_cross_section_definitions.py b/threedigrid_builder/tests/test_cross_section_definitions.py index d2a9be79..99efe7fa 100644 --- a/threedigrid_builder/tests/test_cross_section_definitions.py +++ b/threedigrid_builder/tests/test_cross_section_definitions.py @@ -1,4 +1,3 @@ -from cmath import nan from unittest import mock import numpy as np @@ -38,8 +37,12 @@ def test_convert_multiple(cross_section_definitions): "threedigrid_builder.grid.cross_section_definitions.tabulators", { SHP.CIRCLE: mock.Mock(return_value=(1, 0.1, None, None, None)), - SHP.TABULATED_TRAPEZIUM: mock.Mock(return_value=(5, 15.0, 2.0, table_1, None)), - SHP.TABULATED_RECTANGLE: mock.Mock(return_value=(6, 11.0, 2.0, table_2, None)), + SHP.TABULATED_TRAPEZIUM: mock.Mock( + return_value=(5, 15.0, 2.0, table_1, None) + ), + SHP.TABULATED_RECTANGLE: mock.Mock( + return_value=(6, 11.0, 2.0, table_2, None) + ), }, ): actual = cross_section_definitions.convert([1, 3, 9]) @@ -134,7 +137,9 @@ def test_tabulate_egg(): def test_tabulate_tabulated(): - shape, width_1d, height_1d, table, yz = tabulate_tabulated("my-shape", "1 2 3", "0 1 2") + shape, width_1d, height_1d, table, yz = tabulate_tabulated( + "my-shape", "1 2 3", "0 1 2" + ) assert shape == "my-shape" assert width_1d == 3.0 # the max @@ -142,6 +147,7 @@ def test_tabulate_tabulated(): assert_almost_equal(table, np.array([[0, 1], [1, 2], [2, 3]], dtype=float)) assert yz is None + @pytest.mark.parametrize( "width,height,match", [ @@ -195,7 +201,19 @@ def test_tabulate_inverted_egg(): @pytest.mark.parametrize( "width,height,friction_values,vegetation_stem_densities,vegetation_stem_diameters,vegetation_heights,vegetation_drag_coefficients,exp_width,exp_height,exp_table,exp_yz", [ - ("0 0.5 1 1.5", "0.5 0 0 0.5", "1 1 1", "1 1 1", "1 1 1", "1 1 1", "1 1 1", 1.5, 0.5, [[0, 0.5], [0.5, 1.5]], [[0, 0.5, 1, 1], [0.5, 0, 1, 1], [1, 0, 1, 1], [1.5, 0.5, 0, 0]]), + ( + "0 0.5 1 1.5", + "0.5 0 0 0.5", + "1 1 1", + "1 1 1", + "1 1 1", + "1 1 1", + "1 1 1", + 1.5, + 0.5, + [[0, 0.5], [0.5, 1.5]], + [[0, 0.5, 1, 1], [0.5, 0, 1, 1], [1, 0, 1, 1], [1.5, 0.5, 0, 0]], + ), ( "0 0.5 1 1.5", "0.5 0 0 0.25", @@ -207,7 +225,7 @@ def test_tabulate_inverted_egg(): 1.5, 0.5, [[0, 0.5], [0.25, 1.25], [0.5, 1.5]], - [[0, 0.5, 1, 3], [0.5, 0, 1, 1], [1, 0, 1, 1], [1.5, 0.25, 0, 0]] + [[0, 0.5, 1, 3], [0.5, 0, 1, 1], [1, 0, 1, 1], [1.5, 0.25, 0, 0]], ), ( "0 1 2 3 4 5", @@ -220,7 +238,14 @@ def test_tabulate_inverted_egg(): 5, 1, [[0, 0], [0.5, 3], [0.5, 4], [1, 5]], - [[0, 1, 1, 1], [1, 0, 1, 1], [2, 0.5, 1, 1], [3, 0.5, 1, 1], [4, 0, 1, 1], [5, 1, 0, 0]] + [ + [0, 1, 1, 1], + [1, 0, 1, 1], + [2, 0.5, 1, 1], + [3, 0.5, 1, 1], + [4, 0, 1, 1], + [5, 1, 0, 0], + ], ), ( "0 1 2 2 0 0", @@ -233,14 +258,65 @@ def test_tabulate_inverted_egg(): 2.0, 1.5, [[0, 0], [0.5, 2.0], [1.5, 2.0], [1.5, 0.0]], - None + None, + ), + ( + "0 0.5 0.75 1.0 1.5", + "0.5 0 0 0 0.5", + "1 1 1 1", + "1 1 1 1", + "1 1 1 1", + "1 1 1 1", + "1 1 1 1", + 1.5, + 0.5, + [[0, 0.5], [0.5, 1.5]], + [ + [0, 0.5, 1, 1], + [0.5, 0, 1, 1], + [0.75, 0, 1, 1], + [1, 0, 1, 1], + [1.5, 0.5, 0, 0], + ], + ), + ( + "0 1 0 1 0", + "0 1 1 0 0", + "1 1 1 1", + "1 1 1 1", + "1 1 1 1", + "1 1 1 1", + "1 1 1 1", + 1, + 1, + [[0, 1], [0.5, 0], [1, 1], [1, 0]], + None, ), - ("0 0.5 0.75 1.0 1.5", "0.5 0 0 0 0.5", "1 1 1 1", "1 1 1 1", "1 1 1 1", "1 1 1 1", "1 1 1 1", 1.5, 0.5, [[0, 0.5], [0.5, 1.5]], [[0, 0.5, 1, 1], [0.5, 0, 1, 1], [0.75, 0, 1, 1], [1, 0, 1, 1], [1.5, 0.5, 0, 0]]), - ("0 1 0 1 0", "0 1 1 0 0", "1 1 1 1", "1 1 1 1", "1 1 1 1", "1 1 1 1", "1 1 1 1", 1, 1, [[0, 1], [0.5, 0], [1, 1], [1, 0]], None), ], ) -def test_tabulate_yz(width, height, friction_values,vegetation_stem_densities,vegetation_stem_diameters,vegetation_heights,vegetation_drag_coefficients,exp_width, exp_height, exp_table, exp_yz): - shape, width_1d, height_1d, table, yz = tabulate_yz("my-shape", width, height,friction_values,vegetation_stem_densities,vegetation_stem_diameters,vegetation_heights,vegetation_drag_coefficients) +def test_tabulate_yz( + width, + height, + friction_values, + vegetation_stem_densities, + vegetation_stem_diameters, + vegetation_heights, + vegetation_drag_coefficients, + exp_width, + exp_height, + exp_table, + exp_yz, +): + shape, width_1d, height_1d, table, yz = tabulate_yz( + "my-shape", + width, + height, + friction_values, + vegetation_stem_densities, + vegetation_stem_diameters, + vegetation_heights, + vegetation_drag_coefficients, + ) assert shape == CrossSectionShape.TABULATED_TRAPEZIUM assert width_1d == exp_width diff --git a/threedigrid_builder/tests/test_cross_section_locations.py b/threedigrid_builder/tests/test_cross_section_locations.py index 30a583b7..ba29c562 100644 --- a/threedigrid_builder/tests/test_cross_section_locations.py +++ b/threedigrid_builder/tests/test_cross_section_locations.py @@ -169,8 +169,12 @@ def test_apply_to_lines(channels, channel_lines, locations): assert_equal(channel_lines.frict_type2, [1, 2, 1, 1, 1, 2, 2]) assert_equal(channel_lines.frict_value1, [30, 0.02, 40, 40, 40, 0.03, 0.03]) assert_equal(channel_lines.frict_value2, [30, 0.02, 35, 35, 35, 0.03, 0.03]) - assert_almost_equal(channel_lines.veg_coef1, [0.1, 0.006, 0.2, 0.2, 0.2, 0.036, 0.036]) - assert_almost_equal(channel_lines.veg_coef2, [0.1, 0.006, 0.15, 0.15, 0.15, 0.036, 0.036]) + assert_almost_equal( + channel_lines.veg_coef1, [0.1, 0.006, 0.2, 0.2, 0.2, 0.036, 0.036] + ) + assert_almost_equal( + channel_lines.veg_coef2, [0.1, 0.006, 0.15, 0.15, 0.15, 0.036, 0.036] + ) assert_almost_equal( channel_lines.cross_weight, [1.0, 1.0, 1.0, 0.65, 0.0, 1.0, 1.0] ) From f48fe07eb9eff5ae6406b31b5b36f2c2e2bcae57 Mon Sep 17 00:00:00 2001 From: Golnesa Date: Tue, 30 Jan 2024 12:10:09 +0100 Subject: [PATCH 04/12] Update shape in yz & unit test, fric/veg function --- setup.py | 2 +- threedigrid_builder/base/lines.py | 4 +- .../grid/cross_section_definitions.py | 67 ++++++++++------- .../tests/test_cross_section_definitions.py | 71 ++++++++++++++----- threedigrid_builder/tests/test_db.py | 2 +- threedigrid_builder/tests/test_gridadmin.py | 3 +- 6 files changed, 100 insertions(+), 49 deletions(-) diff --git a/setup.py b/setup.py index d301f23e..575e9999 100644 --- a/setup.py +++ b/setup.py @@ -52,7 +52,7 @@ def get_version(): install_requires = [ "numpy>=1.15,<1.25.0", - "threedi-schema>=0.217.0", + "threedi-schema==0.218.*", "shapely>=2", "pyproj>=3", "condenser[geo]>=0.1.1", diff --git a/threedigrid_builder/base/lines.py b/threedigrid_builder/base/lines.py index d52aa25a..8e0123c2 100644 --- a/threedigrid_builder/base/lines.py +++ b/threedigrid_builder/base/lines.py @@ -38,8 +38,8 @@ class Line: frict_type2: int frict_value1: float # For 1D-2Dgw: resistance factor-out (hydr. cond / thickness) frict_value2: float # For 1D-2Dgw: resistance factor-in (hydr. cond / thickness) - veg_coef1: float - veg_coef2: float + veg_coef1: float # the product of vegetation properties + veg_coef2: float # the product of vegetation properties cross_weight: float # For 1D-2Dgw: the exchange length (tables: cross_width) discharge_coefficient_positive: float discharge_coefficient_negative: float diff --git a/threedigrid_builder/grid/cross_section_definitions.py b/threedigrid_builder/grid/cross_section_definitions.py index 79e1db94..8fa63640 100644 --- a/threedigrid_builder/grid/cross_section_definitions.py +++ b/threedigrid_builder/grid/cross_section_definitions.py @@ -56,11 +56,8 @@ def convert(self, ids): tables = [] tables_yz = [] - # Numpy array views on width/height based on idx - width_idx = self.width[idx] - height_idx = self.height[idx] - - for i, shape in enumerate(self.shape[idx]): + for i, self_i in enumerate(idx): + shape = self.shape[self_i] tabulator = tabulators[shape] ( result.shape[i], @@ -68,12 +65,20 @@ def convert(self, ids): result.height_1d[i], table, yz, - ) = tabulator(shape, width_idx[i], height_idx[i]) + ) = tabulator(shape, self.width[self_i], self.height[self_i]) if table is not None: result.count[i] = len(table) tables.append(table) if yz is not None: result.count_yz[i] = len(yz) + yz = set_friction_vegetation_values( + yz, + self.friction_values[self_i], + self.vegetation_stem_densities[self_i], + self.vegetation_stem_diameters[self_i], + self.vegetation_heights[self_i], + self.vegetation_drag_coefficients[self_i], + ) tables_yz.append(yz) result.offset[:] = np.roll(np.cumsum(result.count), 1) @@ -180,7 +185,7 @@ def tabulate_inverted_egg(shape, width, height): See tabulate_egg. """ - type_, width, height, table, yz = tabulate_egg(shape, width, height) + type_, width, height, table, _ = tabulate_egg(shape, width, height) table[:, 1] = table[::-1, 1] return type_, width, height, table, None @@ -253,16 +258,7 @@ def tabulate_tabulated(shape, width, height): return shape, np.max(widths), np.max(heights), np.array([heights, widths]).T, None -def tabulate_yz( - shape, - width, - height, - friction_values, - vegetation_stem_densities, - vegetation_stem_diameters, - vegetation_heights, - vegetation_drag_coefficients, -): +def tabulate_yz(shape, width, height): """Tabulate an (open or closed) YZ profile Args: @@ -272,10 +268,11 @@ def tabulate_yz( Returns: tuple: shape, width_1d (float), - height_1d (float), table (ndarray of shape (M, 2)) + height_1d (float), table (ndarray of shape (M, 2)), yz (ndarray of shape (M, 4)) """ ys, zs = _parse_tabulated(width, height) is_closed = ys[0] == ys[-1] and zs[0] == zs[-1] + shape_return = CrossSectionShape.TABULATED_TRAPEZIUM if is_closed and len(zs) < 4: raise SchematisationError( f"Cross section definitions of closed profiles must have at least " @@ -310,13 +307,7 @@ def tabulate_yz( yz = np.zeros((len(ys), 4), dtype=float) yz[:, 0] = ys yz[:, 1] = zs - fric = np.array([float(x) for x in friction_values.split(" ")]) - veg_stemden = np.array([float(x) for x in vegetation_stem_densities.split(" ")]) - veg_stemdia = np.array([float(x) for x in vegetation_stem_diameters.split(" ")]) - veg_hght = np.array([float(x) for x in vegetation_heights.split(" ")]) - veg_drag = np.array([float(x) for x in vegetation_drag_coefficients.split(" ")]) - yz[:-1, 2] = fric - yz[:-1, 3] = veg_stemden * veg_stemdia * veg_hght * veg_drag + shape_return = CrossSectionShape.TABULATED_YZ # Adapt non-unique height coordinates. Why? # Because if a segment of the profile is exactly horizontal, we need 2 widths @@ -360,7 +351,7 @@ def tabulate_yz( table = table[1:] height_1d, width_1d = table.max(axis=0).tolist() - return CrossSectionShape.TABULATED_TRAPEZIUM, width_1d, height_1d, table, yz + return shape_return, width_1d, height_1d, table, yz tabulators = { @@ -373,3 +364,27 @@ def tabulate_yz( CrossSectionShape.TABULATED_YZ: tabulate_yz, CrossSectionShape.INVERTED_EGG: tabulate_inverted_egg, } + + +def set_friction_vegetation_values( + yz, + friction_values, + vegetation_stem_densities, + vegetation_stem_diameters, + vegetation_heights, + vegetation_drag_coefficients, +): + """Convert friction and vegetation properties from list into arrays, if available, + and add to yz""" + if len(friction_values) > 0: + fric = np.array([float(x) for x in friction_values.split(" ")]) + yz[:-1, 2] = fric + + if len(vegetation_drag_coefficients) > 0: + veg_stemden = np.array([float(x) for x in vegetation_stem_densities.split(" ")]) + veg_stemdia = np.array([float(x) for x in vegetation_stem_diameters.split(" ")]) + veg_hght = np.array([float(x) for x in vegetation_heights.split(" ")]) + veg_drag = np.array([float(x) for x in vegetation_drag_coefficients.split(" ")]) + yz[:-1, 3] = veg_stemden * veg_stemdia * veg_hght * veg_drag + + return yz diff --git a/threedigrid_builder/tests/test_cross_section_definitions.py b/threedigrid_builder/tests/test_cross_section_definitions.py index 99efe7fa..8b8ac2c0 100644 --- a/threedigrid_builder/tests/test_cross_section_definitions.py +++ b/threedigrid_builder/tests/test_cross_section_definitions.py @@ -8,6 +8,7 @@ from threedigrid_builder.exceptions import SchematisationError from threedigrid_builder.grid import CrossSectionDefinitions, CrossSections from threedigrid_builder.grid.cross_section_definitions import ( + set_friction_vegetation_values, tabulate_builtin, tabulate_closed_rectangle, tabulate_egg, @@ -204,15 +205,15 @@ def test_tabulate_inverted_egg(): ( "0 0.5 1 1.5", "0.5 0 0 0.5", - "1 1 1", - "1 1 1", - "1 1 1", - "1 1 1", - "1 1 1", + "", + "", + "", + "", + "", 1.5, 0.5, [[0, 0.5], [0.5, 1.5]], - [[0, 0.5, 1, 1], [0.5, 0, 1, 1], [1, 0, 1, 1], [1.5, 0.5, 0, 0]], + [[0, 0.5, 0, 0], [0.5, 0, 0, 0], [1, 0, 0, 0], [1.5, 0.5, 0, 0]], ), ( "0 0.5 1 1.5", @@ -264,16 +265,16 @@ def test_tabulate_inverted_egg(): "0 0.5 0.75 1.0 1.5", "0.5 0 0 0 0.5", "1 1 1 1", - "1 1 1 1", - "1 1 1 1", + "1 0.1 1 1", + "2 1 1 1", "1 1 1 1", "1 1 1 1", 1.5, 0.5, [[0, 0.5], [0.5, 1.5]], [ - [0, 0.5, 1, 1], - [0.5, 0, 1, 1], + [0, 0.5, 1, 2], + [0.5, 0, 1, 0.1], [0.75, 0, 1, 1], [1, 0, 1, 1], [1.5, 0.5, 0, 0], @@ -311,19 +312,24 @@ def test_tabulate_yz( "my-shape", width, height, - friction_values, - vegetation_stem_densities, - vegetation_stem_diameters, - vegetation_heights, - vegetation_drag_coefficients, ) - assert shape == CrossSectionShape.TABULATED_TRAPEZIUM assert width_1d == exp_width assert height_1d == exp_height assert_almost_equal(table, np.array(exp_table, dtype=float)) if yz is not None: - assert_almost_equal(yz, np.array(exp_yz, dtype=float)) + yz_1d = set_friction_vegetation_values( + yz, + friction_values, + vegetation_stem_densities, + vegetation_stem_diameters, + vegetation_heights, + vegetation_drag_coefficients, + ) + assert_almost_equal(yz_1d, np.array(exp_yz, dtype=float)) + assert shape == CrossSectionShape.TABULATED_YZ + else: + assert shape == CrossSectionShape.TABULATED_TRAPEZIUM @pytest.mark.parametrize( @@ -347,4 +353,33 @@ def test_tabulate_yz( ) def test_tabulate_yz_err(width, height, match): with pytest.raises(SchematisationError, match=match): - tabulate_yz("my-shape", width, height, None, None, None, None, None) + tabulate_yz("my-shape", width, height) + + +# @pytest.mark.parametrize( +# "yz,friction_values,vegetation_stem_densities,vegetation_stem_diameters,vegetation_heights,vegetation_drag_coefficients,exp_yz", +# [ +# ( +# [[0, 0, 0, 0], [0, 0, 0, 0], [0, 0, 0, 0], [0, 0, 0, 0]], +# "0.2 0.3 0.1", +# "1 2 3", +# "0.5 0.1 0.3", +# "2 3 1", +# "1 1 1", +# [ +# [0, 0, 0.2, 1], +# [0, 0, 0.3, 0.6], +# [0, 0, 0.1, 0.9], +# [0, 0, 0, 0], +# ], +# ), +# ], +# ) +# def test_set_friction_vegetation_values(yz,friction_values,vegetation_stem_densities,vegetation_stem_diameters,vegetation_heights,vegetation_drag_coefficients,exp_yz): +# yz_1d = set_friction_vegetation_values(yz, +# friction_values, +# vegetation_stem_densities, +# vegetation_stem_diameters, +# vegetation_heights, +# vegetation_drag_coefficients) +# assert_almost_equal(yz_1d, np.array(exp_yz, dtype=float)) diff --git a/threedigrid_builder/tests/test_db.py b/threedigrid_builder/tests/test_db.py index 7f329c9d..435d2202 100644 --- a/threedigrid_builder/tests/test_db.py +++ b/threedigrid_builder/tests/test_db.py @@ -65,7 +65,7 @@ def test_init_bad_version(tmp_path): def test_get_version(db): - assert db.get_version() == 217 + assert db.get_version() == 218 def test_get_boundary_conditions_1d(db): diff --git a/threedigrid_builder/tests/test_gridadmin.py b/threedigrid_builder/tests/test_gridadmin.py index 2f1dc334..ff8bbc45 100644 --- a/threedigrid_builder/tests/test_gridadmin.py +++ b/threedigrid_builder/tests/test_gridadmin.py @@ -140,7 +140,8 @@ def test_write_nodes_embedded(h5_out, dataset, shape, dtype): ("dist_calc_points", (6,), "float64"), ("friction_type", (6,), "int32"), ("friction_value", (6,), "float64"), - # ("vegetation_stem_density", (6,), "float64"), + ("veg_coef1", (6,), "float64"), + ("veg_coef2", (6,), "float64"), ("id", (6,), "int32"), ("invert_level_end_point", (6,), "float64"), ("invert_level_start_point", (6,), "float64"), From 1d4ba2d7bf085cc8f759a3dae5381a0558d4bd71 Mon Sep 17 00:00:00 2001 From: Golnesa Date: Wed, 31 Jan 2024 11:56:03 +0100 Subject: [PATCH 05/12] Clean up & bergemeer sql update --- .../tests/data/v2_bergermeer.sqlite | 2 +- .../tests/test_cross_section_definitions.py | 29 ------------------- 2 files changed, 1 insertion(+), 30 deletions(-) diff --git a/threedigrid_builder/tests/data/v2_bergermeer.sqlite b/threedigrid_builder/tests/data/v2_bergermeer.sqlite index aed1a37a..f1a67ec5 100755 --- a/threedigrid_builder/tests/data/v2_bergermeer.sqlite +++ b/threedigrid_builder/tests/data/v2_bergermeer.sqlite @@ -1,3 +1,3 @@ version https://git-lfs.github.com/spec/v1 -oid sha256:656f80f442a7daf84a7e69e55f4737769bcb9734531747aec081bd133202c6b8 +oid sha256:629fa87b841f65f613f45a41238ca903f1b49a271f61117bddda1538107273af size 6317056 diff --git a/threedigrid_builder/tests/test_cross_section_definitions.py b/threedigrid_builder/tests/test_cross_section_definitions.py index 8b8ac2c0..dc26fc13 100644 --- a/threedigrid_builder/tests/test_cross_section_definitions.py +++ b/threedigrid_builder/tests/test_cross_section_definitions.py @@ -354,32 +354,3 @@ def test_tabulate_yz( def test_tabulate_yz_err(width, height, match): with pytest.raises(SchematisationError, match=match): tabulate_yz("my-shape", width, height) - - -# @pytest.mark.parametrize( -# "yz,friction_values,vegetation_stem_densities,vegetation_stem_diameters,vegetation_heights,vegetation_drag_coefficients,exp_yz", -# [ -# ( -# [[0, 0, 0, 0], [0, 0, 0, 0], [0, 0, 0, 0], [0, 0, 0, 0]], -# "0.2 0.3 0.1", -# "1 2 3", -# "0.5 0.1 0.3", -# "2 3 1", -# "1 1 1", -# [ -# [0, 0, 0.2, 1], -# [0, 0, 0.3, 0.6], -# [0, 0, 0.1, 0.9], -# [0, 0, 0, 0], -# ], -# ), -# ], -# ) -# def test_set_friction_vegetation_values(yz,friction_values,vegetation_stem_densities,vegetation_stem_diameters,vegetation_heights,vegetation_drag_coefficients,exp_yz): -# yz_1d = set_friction_vegetation_values(yz, -# friction_values, -# vegetation_stem_densities, -# vegetation_stem_diameters, -# vegetation_heights, -# vegetation_drag_coefficients) -# assert_almost_equal(yz_1d, np.array(exp_yz, dtype=float)) From 71caea3b8c03f9b15b9e986f054fee0a7f678828 Mon Sep 17 00:00:00 2001 From: Golnesa Date: Thu, 1 Feb 2024 10:53:21 +0100 Subject: [PATCH 06/12] comments --- threedigrid_builder/grid/cross_section_locations.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/threedigrid_builder/grid/cross_section_locations.py b/threedigrid_builder/grid/cross_section_locations.py index 8313ad34..67dccb45 100644 --- a/threedigrid_builder/grid/cross_section_locations.py +++ b/threedigrid_builder/grid/cross_section_locations.py @@ -45,6 +45,8 @@ def apply_to_lines(self, lines, channels, extrapolate=False): - frict_type2: the friction type of the second cross section location - frict_value1: the friction value of the first cross section location - frict_value2: the friction value of the second cross section location + - veg_coef1: the product of vegetation properties of the first cross section location + - veg_coef2: the product of vegetation properties of the second cross section location - invert_level_start_point: 'reference_level' interpolated at the line end - invert_level_end_point: 'reference_level' interpolated at the line start - dpumax: the largest of the two invert levels From 4061d139c73d475ef375d1b1655ebd618f4d16c8 Mon Sep 17 00:00:00 2001 From: Golnesa Date: Thu, 1 Feb 2024 10:54:48 +0100 Subject: [PATCH 07/12] test.yml --- .github/workflows/test.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.github/workflows/test.yml b/.github/workflows/test.yml index ad6ca8f5..54e4511e 100644 --- a/.github/workflows/test.yml +++ b/.github/workflows/test.yml @@ -95,4 +95,4 @@ jobs: - name: Run integration tests shell: bash run: | - pytest -c /dev/null integration_tests + pytest integration_tests From 4a259bb99a23fc8be0b93cbfdd244d9e842ab113 Mon Sep 17 00:00:00 2001 From: Golnesa Date: Thu, 1 Feb 2024 14:23:03 +0100 Subject: [PATCH 08/12] update no 'single' vegetation --- .../grid/cross_section_locations.py | 4 ++++ .../tests/test_cross_section_locations.py | 20 ++++++++----------- 2 files changed, 12 insertions(+), 12 deletions(-) diff --git a/threedigrid_builder/grid/cross_section_locations.py b/threedigrid_builder/grid/cross_section_locations.py index 67dccb45..8c6c0ee2 100644 --- a/threedigrid_builder/grid/cross_section_locations.py +++ b/threedigrid_builder/grid/cross_section_locations.py @@ -70,18 +70,22 @@ def apply_to_lines(self, lines, channels, extrapolate=False): lines.frict_type2 = self.friction_type[idx2] lines.frict_value1 = self.friction_value[idx1] lines.frict_value2 = self.friction_value[idx2] + lines.veg_coef1 = ( self.vegetation_stem_density[idx1] * self.vegetation_stem_diameter[idx1] * self.vegetation_height[idx1] * self.vegetation_drag_coefficient[idx1] ) + lines.veg_coef1[np.isnan(lines.veg_coef1)] = 0.0 + lines.veg_coef2 = ( self.vegetation_stem_density[idx2] * self.vegetation_stem_diameter[idx2] * self.vegetation_height[idx2] * self.vegetation_drag_coefficient[idx2] ) + lines.veg_coef2[np.isnan(lines.veg_coef2)] = 0.0 # Compute invert levels and start and end lines.invert_level_start_point = compute_bottom_level( diff --git a/threedigrid_builder/tests/test_cross_section_locations.py b/threedigrid_builder/tests/test_cross_section_locations.py index ba29c562..de0d380d 100644 --- a/threedigrid_builder/tests/test_cross_section_locations.py +++ b/threedigrid_builder/tests/test_cross_section_locations.py @@ -62,10 +62,10 @@ def locations(): channel_id=[51, 53, 53, 52, 54], friction_type=[1, 1, 1, 2, 2], friction_value=[30, 35, 40, 0.02, 0.03], - vegetation_stem_density=[1, 0.3, 0.4, 0.02, 0.03], - vegetation_stem_diameter=[0.5, 0.5, 0.5, 0.3, 0.4], - vegetation_height=[0.2, 1, 1, 1, 3], - vegetation_drag_coefficient=[1, 1, 1, 1, 1], + vegetation_stem_density=[1, 0.3, 0.4, np.nan, np.nan], + vegetation_stem_diameter=[0.5, 0.5, 0.5, np.nan, np.nan], + vegetation_height=[0.2, 1, 1, np.nan, np.nan], + vegetation_drag_coefficient=[1, 1, 1, np.nan, np.nan], reference_level=[1.0, 2.0, 3.0, 5.0, 6.0], definition_id=[3, 3, 4, 4, 5], ) @@ -167,14 +167,10 @@ def test_apply_to_lines(channels, channel_lines, locations): assert_equal(channel_lines.cross_id2, [3, 4, 3, 3, 3, 5, 5]) assert_equal(channel_lines.frict_type1, [1, 2, 1, 1, 1, 2, 2]) assert_equal(channel_lines.frict_type2, [1, 2, 1, 1, 1, 2, 2]) - assert_equal(channel_lines.frict_value1, [30, 0.02, 40, 40, 40, 0.03, 0.03]) - assert_equal(channel_lines.frict_value2, [30, 0.02, 35, 35, 35, 0.03, 0.03]) - assert_almost_equal( - channel_lines.veg_coef1, [0.1, 0.006, 0.2, 0.2, 0.2, 0.036, 0.036] - ) - assert_almost_equal( - channel_lines.veg_coef2, [0.1, 0.006, 0.15, 0.15, 0.15, 0.036, 0.036] - ) + assert_equal(channel_lines.frict_value1, [30.0, 0.02, 40.0, 40.0, 40.0, 0.03, 0.03]) + assert_equal(channel_lines.frict_value2, [30.0, 0.02, 35.0, 35.0, 35.0, 0.03, 0.03]) + assert_almost_equal(channel_lines.veg_coef1, [0.1, 0.0, 0.2, 0.2, 0.2, 0.0, 0.0]) + assert_almost_equal(channel_lines.veg_coef2, [0.1, 0.0, 0.15, 0.15, 0.15, 0.0, 0.0]) assert_almost_equal( channel_lines.cross_weight, [1.0, 1.0, 1.0, 0.65, 0.0, 1.0, 1.0] ) From 572616bb8665195b267cec62589c9bb55d75b4af Mon Sep 17 00:00:00 2001 From: Golnesa Date: Thu, 8 Feb 2024 09:16:56 +0100 Subject: [PATCH 09/12] bergemeer update + corrections --- integration_tests/test_bergermeer.py | 6 +++--- threedigrid_builder/grid/cross_section_definitions.py | 4 ++-- threedigrid_builder/tests/data/v2_bergermeer.sqlite | 2 +- .../tests/test_cross_section_definitions.py | 10 +++++----- threedigrid_builder/tests/test_db.py | 2 +- 5 files changed, 12 insertions(+), 12 deletions(-) mode change 100755 => 100644 threedigrid_builder/tests/data/v2_bergermeer.sqlite diff --git a/integration_tests/test_bergermeer.py b/integration_tests/test_bergermeer.py index fe329cfb..e323fbc2 100644 --- a/integration_tests/test_bergermeer.py +++ b/integration_tests/test_bergermeer.py @@ -20,9 +20,7 @@ def count_unique(arr): return dict(zip(*np.unique(arr, return_counts=True))) -@pytest.mark.parametrize( - "filename", ["v2_bergermeer.sqlite", "v2_bergermeer_spatialite4.sqlite"] -) +@pytest.mark.parametrize("filename", ["v2_bergermeer.sqlite"]) def test_integration(tmp_path, filename): shutil.copyfile(unittests_data_path / filename, tmp_path / filename) progress_callback = Mock() @@ -108,6 +106,8 @@ def test_integration(tmp_path, filename): -9999: 1, CrossSectionShape.CIRCLE: 8, CrossSectionShape.RECTANGLE: 3, + CrossSectionShape.TABULATED_TRAPEZIUM: 1, + CrossSectionShape.TABULATED_YZ: 1, } ## EMBEDDED NODES diff --git a/threedigrid_builder/grid/cross_section_definitions.py b/threedigrid_builder/grid/cross_section_definitions.py index 8fa63640..1081e166 100644 --- a/threedigrid_builder/grid/cross_section_definitions.py +++ b/threedigrid_builder/grid/cross_section_definitions.py @@ -376,11 +376,11 @@ def set_friction_vegetation_values( ): """Convert friction and vegetation properties from list into arrays, if available, and add to yz""" - if len(friction_values) > 0: + if friction_values is not None: fric = np.array([float(x) for x in friction_values.split(" ")]) yz[:-1, 2] = fric - if len(vegetation_drag_coefficients) > 0: + if vegetation_drag_coefficients is not None: veg_stemden = np.array([float(x) for x in vegetation_stem_densities.split(" ")]) veg_stemdia = np.array([float(x) for x in vegetation_stem_diameters.split(" ")]) veg_hght = np.array([float(x) for x in vegetation_heights.split(" ")]) diff --git a/threedigrid_builder/tests/data/v2_bergermeer.sqlite b/threedigrid_builder/tests/data/v2_bergermeer.sqlite old mode 100755 new mode 100644 index f1a67ec5..c5d762dc --- a/threedigrid_builder/tests/data/v2_bergermeer.sqlite +++ b/threedigrid_builder/tests/data/v2_bergermeer.sqlite @@ -1,3 +1,3 @@ version https://git-lfs.github.com/spec/v1 -oid sha256:629fa87b841f65f613f45a41238ca903f1b49a271f61117bddda1538107273af +oid sha256:cbf6018bef282854cdf591d4b568d67390b4a66a5bd8cead0951934865530ab0 size 6317056 diff --git a/threedigrid_builder/tests/test_cross_section_definitions.py b/threedigrid_builder/tests/test_cross_section_definitions.py index dc26fc13..12290479 100644 --- a/threedigrid_builder/tests/test_cross_section_definitions.py +++ b/threedigrid_builder/tests/test_cross_section_definitions.py @@ -205,11 +205,11 @@ def test_tabulate_inverted_egg(): ( "0 0.5 1 1.5", "0.5 0 0 0.5", - "", - "", - "", - "", - "", + None, + None, + None, + None, + None, 1.5, 0.5, [[0, 0.5], [0.5, 1.5]], diff --git a/threedigrid_builder/tests/test_db.py b/threedigrid_builder/tests/test_db.py index 435d2202..c0f31f45 100644 --- a/threedigrid_builder/tests/test_db.py +++ b/threedigrid_builder/tests/test_db.py @@ -143,7 +143,7 @@ def test_get_cross_section_definitions(db): assert isinstance(definitions, CrossSectionDefinitions) # some test samples - assert len(definitions.id) == 11 + assert len(definitions.id) == 13 assert definitions.id[8] == 97 assert definitions.shape[7] == CrossSectionShape.RECTANGLE assert definitions.height[10] == "0.4" From bf063ea4900808291540a3cc1e086559af0f5324 Mon Sep 17 00:00:00 2001 From: Golnesa Date: Thu, 8 Feb 2024 11:48:23 +0100 Subject: [PATCH 10/12] update docstrings + minor corrections + changelog --- setup.py | 2 +- threedigrid_builder/grid/cross_section_definitions.py | 10 +++++----- .../tests/data/v2_bergermeer_spatialite4.sqlite | 3 --- 3 files changed, 6 insertions(+), 9 deletions(-) delete mode 100644 threedigrid_builder/tests/data/v2_bergermeer_spatialite4.sqlite diff --git a/setup.py b/setup.py index 575e9999..d6ecf44a 100644 --- a/setup.py +++ b/setup.py @@ -52,7 +52,7 @@ def get_version(): install_requires = [ "numpy>=1.15,<1.25.0", - "threedi-schema==0.218.*", + "threedi-schema==0.219.1*", "shapely>=2", "pyproj>=3", "condenser[geo]>=0.1.1", diff --git a/threedigrid_builder/grid/cross_section_definitions.py b/threedigrid_builder/grid/cross_section_definitions.py index 1081e166..9bf564d2 100644 --- a/threedigrid_builder/grid/cross_section_definitions.py +++ b/threedigrid_builder/grid/cross_section_definitions.py @@ -129,7 +129,7 @@ def tabulate_builtin(shape, width, height): height (str): ignored Returns: - tuple: shape, width_1d (float), None, None + tuple: shape, width_1d (float), None, None, None """ try: width = float(width) @@ -151,7 +151,7 @@ def tabulate_egg(shape, width, height): Returns: tuple: TABULATED_TRAPEZIUM, width_1d (float), - height_1d (float), table (ndarray of shape (M, 2)) + height_1d (float), table (ndarray of shape (M, 2)), None """ NUM_INCREMENTS = 16 @@ -200,7 +200,7 @@ def tabulate_closed_rectangle(shape, width, height): Returns: tuple: TABULATED_RECTANGLE, width_1d (float), - height (float), table (ndarray of shape (M, 2)) + height (float), table (ndarray of shape (M, 2)), None """ try: width = float(width) @@ -246,7 +246,7 @@ def tabulate_tabulated(shape, width, height): Returns: tuple: shape, width_1d (float), - height_1d (float), table (ndarray of shape (M, 2)) + height_1d (float), table (ndarray of shape (M, 2)), None """ widths, heights = _parse_tabulated(width, height) if len(heights) > 1 and np.any(np.diff(heights) < 0.0): @@ -272,7 +272,6 @@ def tabulate_yz(shape, width, height): """ ys, zs = _parse_tabulated(width, height) is_closed = ys[0] == ys[-1] and zs[0] == zs[-1] - shape_return = CrossSectionShape.TABULATED_TRAPEZIUM if is_closed and len(zs) < 4: raise SchematisationError( f"Cross section definitions of closed profiles must have at least " @@ -303,6 +302,7 @@ def tabulate_yz(shape, width, height): ys = ys[:-1] zs = zs[:-1] yz = None + shape_return = CrossSectionShape.TABULATED_TRAPEZIUM else: yz = np.zeros((len(ys), 4), dtype=float) yz[:, 0] = ys diff --git a/threedigrid_builder/tests/data/v2_bergermeer_spatialite4.sqlite b/threedigrid_builder/tests/data/v2_bergermeer_spatialite4.sqlite deleted file mode 100644 index 4ba53f45..00000000 --- a/threedigrid_builder/tests/data/v2_bergermeer_spatialite4.sqlite +++ /dev/null @@ -1,3 +0,0 @@ -version https://git-lfs.github.com/spec/v1 -oid sha256:ca42495e19176ad75b0f2db63a057e64fc454ab50a19d8e89a9980e62d9e4808 -size 7467008 From 895ab2bc034da975ada9309796b4b5afac68831b Mon Sep 17 00:00:00 2001 From: Golnesa Date: Thu, 8 Feb 2024 12:04:14 +0100 Subject: [PATCH 11/12] schema version update --- CHANGES.rst | 8 +++++++- 1 file changed, 7 insertions(+), 1 deletion(-) diff --git a/CHANGES.rst b/CHANGES.rst index f5137843..43982c1a 100644 --- a/CHANGES.rst +++ b/CHANGES.rst @@ -4,7 +4,13 @@ Changelog of threedigrid-builder 1.12.3 (unreleased) ------------------- -- Nothing changed yet. +- Add single vegetation and variable friction/vegetation for 1D elements. + +- Add a table for open YZ profile to support 1D variable friction and vegetation. + +- Update unit tests and the integration test (bergemeer test). + +- Remove the spatialie4 version of the integration test (bergemeer test). 1.12.2 (2023-12-07) From 224a37d0cdb784bc811ba7602c1578ea15dfbd69 Mon Sep 17 00:00:00 2001 From: Golnesa Date: Thu, 8 Feb 2024 12:06:59 +0100 Subject: [PATCH 12/12] schema version update --- setup.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/setup.py b/setup.py index d6ecf44a..af4a7fd0 100644 --- a/setup.py +++ b/setup.py @@ -52,7 +52,7 @@ def get_version(): install_requires = [ "numpy>=1.15,<1.25.0", - "threedi-schema==0.219.1*", + "threedi-schema==0.219.*", "shapely>=2", "pyproj>=3", "condenser[geo]>=0.1.1",