diff --git a/swmmanywhere/geospatial_utilities.py b/swmmanywhere/geospatial_utilities.py index 712e775a..51fce25d 100644 --- a/swmmanywhere/geospatial_utilities.py +++ b/swmmanywhere/geospatial_utilities.py @@ -614,8 +614,8 @@ def remove_(mp): return remove_zero_area_subareas(mp, removed_subareas) return polys_gdf def derive_rc(polys_gdf: gpd.GeoDataFrame, - G: nx.Graph, - building_footprints: gpd.GeoDataFrame) -> gpd.GeoDataFrame: + building_footprints: gpd.GeoDataFrame, + streetcover: gpd.GeoDataFrame) -> gpd.GeoDataFrame: """Derive the Runoff Coefficient (RC) of each subcatchment. The runoff coefficient is the ratio of impervious area to total area. The @@ -626,10 +626,10 @@ def derive_rc(polys_gdf: gpd.GeoDataFrame, Args: polys_gdf (gpd.GeoDataFrame): A GeoDataFrame containing polygons that represent subcatchments with columns: 'geometry', 'area', and 'id'. - G (nx.Graph): The input graph, with node 'ids' that match polys_gdf and - edges with the 'id', 'width' and 'geometry' property. building_footprints (gpd.GeoDataFrame): A GeoDataFrame containing building footprints with a 'geometry' column. + streetcover (gpd.GeoDataFrame): A GeoDataFrame containing street cover + with a 'geometry' column. Returns: gpd.GeoDataFrame: A GeoDataFrame containing polygons with columns: @@ -638,23 +638,7 @@ def derive_rc(polys_gdf: gpd.GeoDataFrame, polys_gdf = polys_gdf.copy() ## Format as swmm type catchments - - # TODO think harder about lane widths (am I double counting here?) - lines = [ - { - 'geometry': x['geometry'].buffer(x['width'], - cap_style=2, - join_style=2), - 'id': x['id'] - } - for u, v, x in G.edges(data=True) - ] - lines_df = pd.DataFrame(lines) - lines_gdf = gpd.GeoDataFrame(lines_df, - geometry=lines_df.geometry, - crs = polys_gdf.crs) - - result = gpd.overlay(lines_gdf[['geometry']], + result = gpd.overlay(streetcover[['geometry']], building_footprints[['geometry']], how='union') result = gpd.overlay(polys_gdf, result) diff --git a/swmmanywhere/graph_utilities.py b/swmmanywhere/graph_utilities.py index 9a121d0f..75c8cb30 100644 --- a/swmmanywhere/graph_utilities.py +++ b/swmmanywhere/graph_utilities.py @@ -266,36 +266,55 @@ def __call__(self, return G @register_graphfcn -class format_osmnx_lanes(BaseGraphFunction, - required_edge_attributes = ['lanes'], - adds_edge_attributes = ['width']): - """format_osmnx_lanes class.""" +class calculate_streetcover(BaseGraphFunction, + required_edge_attributes = ['lanes'] + ): + """calculate_streetcover class.""" # i.e., in osmnx format, i.e., empty for single lane, an int for a # number of lanes or a list if the edge has multiple carriageways def __call__(self, G: nx.Graph, - subcatchment_derivation: parameters.SubcatchmentDerivation, - **kwargs) -> nx.Graph: + subcatchment_derivation: parameters.SubcatchmentDerivation, + addresses: parameters.FilePaths, + **kwargs) -> nx.Graph: """Format the lanes attribute of each edge and calculates width. Args: G (nx.Graph): A graph subcatchment_derivation (parameters.SubcatchmentDerivation): A SubcatchmentDerivation parameter object + addresses (parameters.FilePaths): A FilePaths parameter object **kwargs: Additional keyword arguments are ignored. Returns: G (nx.Graph): A graph """ G = G.copy() + lines = [] for u, v, data in G.edges(data=True): lanes = data.get('lanes',1) if isinstance(lanes, list): lanes = sum([float(x) for x in lanes]) else: lanes = float(lanes) - data['width'] = lanes * subcatchment_derivation.lane_width + lines.append({'geometry' : data['geometry'].buffer(lanes * + subcatchment_derivation.lane_width, + cap_style=2, + join_style=2), + 'u' : u, + 'v' : v + } + ) + lines_df = pd.DataFrame(lines) + lines_gdf = gpd.GeoDataFrame(lines_df, + geometry=lines_df.geometry, + crs = G.graph['crs']) + if addresses.streetcover.suffix in ('.geoparquet','.parquet'): + lines_gdf.to_parquet(addresses.streetcover) + else: + lines_gdf.to_file(addresses.streetcover, driver='GeoJSON') + return G @register_graphfcn @@ -494,7 +513,7 @@ def __call__(self, G: nx.Graph, **kwargs) -> nx.Graph: @register_graphfcn class calculate_contributing_area(BaseGraphFunction, - required_edge_attributes = ['id', 'geometry', 'width'], + required_edge_attributes = ['id', 'geometry'], adds_edge_attributes = ['contributing_area'], adds_node_attributes = ['contributing_area']): """calculate_contributing_area class.""" @@ -507,7 +526,10 @@ def __call__(self, G: nx.Graph, This function calculates the contributing area for each edge. The contributing area is the area of the subcatchment that drains to the - edge. The contributing area is calculated from the elevation data. + edge. The contributing area is calculated from the elevation data. + Runoff coefficient (RC) for each contributing area is also calculated, + the RC is calculated using `addresses.buildings` and + `addresses.streetcover`. Also writes the file 'subcatchments.geojson' to addresses.subcatchments. @@ -543,7 +565,12 @@ def __call__(self, G: nx.Graph, buildings = gpd.read_parquet(addresses.building) else: buildings = gpd.read_file(addresses.building) - subs_rc = go.derive_rc(subs_gdf, G, buildings) + if addresses.streetcover.suffix in ('.geoparquet','.parquet'): + streetcover = gpd.read_parquet(addresses.streetcover) + else: + streetcover = gpd.read_file(addresses.streetcover) + + subs_rc = go.derive_rc(subs_gdf, buildings, streetcover) # Write subs # TODO - could just attach subs to nodes where each node has a list of subs diff --git a/swmmanywhere/parameters.py b/swmmanywhere/parameters.py index 05bbf187..b8b6c8b0 100644 --- a/swmmanywhere/parameters.py +++ b/swmmanywhere/parameters.py @@ -292,6 +292,9 @@ def _generate_elevation(self): def _generate_building(self): return self._generate_property(f'building.geo{self.extension}', 'download') + def _generate_streetcover(self): + return self._generate_property(f'streetcover.geo{self.extension}', + 'model') def _generate_precipitation(self): return self._generate_property(f'precipitation.{self.extension}', 'download') diff --git a/tests/test_data/demo_config.yml b/tests/test_data/demo_config.yml index ecbb15aa..ca4d8db8 100644 --- a/tests/test_data/demo_config.yml +++ b/tests/test_data/demo_config.yml @@ -14,8 +14,8 @@ real: starting_graph: null graphfcn_list: - assign_id - - format_osmnx_lanes - remove_non_pipe_allowable_links + - calculate_streetcover - double_directed - fix_geometries - split_long_edges diff --git a/tests/test_geospatial_utilities.py b/tests/test_geospatial_utilities.py index b158a367..34d952e5 100644 --- a/tests/test_geospatial_utilities.py +++ b/tests/test_geospatial_utilities.py @@ -279,20 +279,16 @@ def test_derive_rc(): crs = crs) subs['area'] = subs.geometry.area - subs_rc = go.derive_rc(subs, G, buildings).set_index('id') + subs_rc = go.derive_rc(subs, buildings, buildings).set_index('id') assert subs_rc.loc[6277683849,'impervious_area'] == 0 assert subs_rc.loc[107733,'impervious_area'] > 0 - for u,v,d in G.edges(data=True): - d['width'] = 10 - subs_rc = go.derive_rc(subs, G, buildings).set_index('id') + buildings.geometry = buildings.buffer(50) + subs_rc = go.derive_rc(subs, buildings, buildings).set_index('id') assert subs_rc.loc[6277683849,'impervious_area'] > 0 assert subs_rc.loc[6277683849,'rc'] > 0 assert subs_rc.rc.max() <= 100 - for u,v,d in G.edges(data=True): - d['width'] = 0 - def test_calculate_angle(): """Test the calculate_angle function.""" # Test with points forming a right angle diff --git a/tests/test_graph_utilities.py b/tests/test_graph_utilities.py index fe597056..28436f9c 100644 --- a/tests/test_graph_utilities.py +++ b/tests/test_graph_utilities.py @@ -52,14 +52,23 @@ def test_double_directed(): for u, v in G.edges(): assert (v,u) in G.edges -def test_format_osmnx_lanes(): - """Test the format_osmnx_lanes function.""" +def test_calculate_streetcover(): + """Test the calculate_streetcover function.""" G, _ = load_street_network() params = parameters.SubcatchmentDerivation() - G = gu.format_osmnx_lanes(G, params) - for u, v, data in G.edges(data=True): - assert 'width' in data.keys() - assert isinstance(data['width'], float) + addresses = parameters.FilePaths(base_dir = None, + project_name = None, + bbox_number = None, + model_number = None, + extension = 'json') + with tempfile.TemporaryDirectory() as temp_dir: + addresses.streetcover = Path(temp_dir) / 'streetcover.geojson' + _ = gu.calculate_streetcover(G, params, addresses) + # TODO test that G hasn't changed? or is that a waste of time? + assert addresses.streetcover.exists() + gdf = gpd.read_file(addresses.streetcover) + assert len(gdf) == len(G.edges) + assert gdf.geometry.area.sum() > 0 def test_split_long_edges(): """Test the split_long_edges function.""" @@ -82,6 +91,7 @@ def test_derive_subcatchments(): model_number = 1) addresses.elevation = Path(__file__).parent / 'test_data' / 'elevation.tif' addresses.building = temp_path / 'building.geojson' + addresses.streetcover = temp_path / 'building.geojson' addresses.subcatchments = temp_path / 'subcatchments.geojson' params = parameters.SubcatchmentDerivation() G, bbox = load_street_network() @@ -275,6 +285,7 @@ def test_iterate_graphfcns(): """Test the iterate_graphfcns function.""" G = load_graph(Path(__file__).parent / 'test_data' / 'graph_topo_derived.json') params = parameters.get_full_parameters() + params['topology_derivation'].omit_edges = ['primary', 'bridge'] with tempfile.TemporaryDirectory() as temp_dir: temp_path = Path(temp_dir) addresses = parameters.FilePaths(base_dir = None, @@ -286,12 +297,13 @@ def test_iterate_graphfcns(): addresses.model = temp_path G = iterate_graphfcns(G, ['assign_id', - 'format_osmnx_lanes'], + 'remove_non_pipe_allowable_links'], params, addresses) for u, v, d in G.edges(data=True): assert 'id' in d.keys() - assert 'width' in d.keys() + assert 'primary' not in get_edge_types(G) + assert len(set([d.get('bridge',None) for u,v,d in G.edges(data=True)])) == 1 def test_fix_geometries(): """Test the fix_geometries function."""