From 40379563357b94afd19cdc0c48c59fce195b6a05 Mon Sep 17 00:00:00 2001 From: Rose Pearson Date: Fri, 3 Nov 2023 13:19:25 +1300 Subject: [PATCH] Tidy-up to remove the now unused coarse DEM chunking code --- src/geofabrics/dem.py | 157 ------------------------------------------ 1 file changed, 157 deletions(-) diff --git a/src/geofabrics/dem.py b/src/geofabrics/dem.py index 106ccd68..1c908a06 100644 --- a/src/geofabrics/dem.py +++ b/src/geofabrics/dem.py @@ -1836,163 +1836,6 @@ def add_coarse_dems( 0, ) - """if chunk_size is None: - self._add_coarse_dem_no_chunking( - coarse_dem_path=coarse_dem_path, - mask=no_value_mask, - ) - else: - self._add_coarse_dem_chunked( - coarse_dem_path=coarse_dem_path, - chunk_size=chunk_size, - mask=no_value_mask, - radius=coarse_dem_resolution * numpy.sqrt(2), - )""" - - def _add_coarse_dem_no_chunking( - self, - coarse_dem_path: pathlib.Path, - mask: numpy.ndarray, - ): - """Fill gaps in dense DEM from areas with no LiDAR with the coarse DEM. - Perform linear interpolation. - - Parameters - ---------- - - coarse_dem - The coarse DEM to use for rasterising - mask - the pixel mask of where to set these values - - """ - - # Load in the coarse DEM - extents = { - "total": self.catchment_geometry.catchment, - "land": self.catchment_geometry.full_land, - "foreshore": self.catchment_geometry.foreshore, - } - coarse_dem = CoarseDem( - dem_file=coarse_dem_path, - extents=extents, # by reference, so "total" trimmed to coarse DEM bounds - set_foreshore=self.drop_offshore_lidar, - ) - - # create dictionary defining raster options - # Set search radius to the diagonal cell length to ensure corners covered - raster_options = { - "raster_type": geometry.RASTER_TYPE, - "radius": coarse_dem.resolution * numpy.sqrt(2), - "method": "linear", - } - # Get the grid locations overwhich to perform averaging - grid_x, grid_y = numpy.meshgrid(self._dem.x, self._dem.y) - - # Mask to only rasterise where there aren't already LiDAR derived values - xy_out = numpy.empty((mask.values.sum(), 2)) - xy_out[:, 0] = grid_x[mask] - xy_out[:, 1] = grid_y[mask] - - # Perform specified averaging from the coarse DEM where there is no data - z_flat = elevation_from_points( - point_cloud=coarse_dem.points, - xy_out=xy_out, - options=raster_options, - ) - # Update the DEM - self._dem.z.data[mask] = z_flat - # Update the data source layer - self._dem["data_source"] = self._dem.data_source.where( - ~(mask & self._dem.z.notnull().values), - self.SOURCE_CLASSIFICATION["coarse DEM"], - ) - - def _add_coarse_dem_chunked( - self, - coarse_dem_path: pathlib.Path, - chunk_size: int, - mask: numpy.ndarray, - radius: float, - ): - """Fill gaps in dense DEM from areas with no LiDAR with the coarse DEM. - Perform linear interpolation. - - Parameters - ---------- - - coarse_dem - The coarse DEM to use for rasterising - extents - a dictionary of area to rasterise using the coarse DEM - chunk_size - the size of each chunk - - """ - - # Define raster options - raster_options = { - "raster_type": geometry.RASTER_TYPE, - "elevation_range": self.elevation_range, - "radius": radius, - "method": "linear", - } - extents = { - "total": self.catchment_geometry.catchment, - "land": self.catchment_geometry.full_land, - "foreshore": self.catchment_geometry.foreshore, - } - - # get chunking information - chunked_dim_x, chunked_dim_y = self._chunks_from_dem(chunk_size, self._dem) - elevations = {} - - # cycle through index chunks - and collect in a delayed array - delayed_chunked_matrix = [] - for i, dim_y in enumerate(chunked_dim_y): - delayed_chunked_x = [] - for j, dim_x in enumerate(chunked_dim_x): - logging.info(f"\tCoarse chunk {[i, j]}") - - # Define the region of the chunk to rasterise - chunk_region_to_tile = self._define_chunk_region( - region_to_rasterise=extents["total"], - dim_x=dim_x, - dim_y=dim_y, - radius=raster_options["radius"], - ) - chunk_extents = { # reset to unbuffered land and foreshore - "total": chunk_region_to_tile, - "land": self.catchment_geometry.full_land, - "foreshore": self.catchment_geometry.foreshore, - } - - # Send through the clipped Coarse DEM points - coarse_points = delayed_chunk_coarse_dem( - dem_file=coarse_dem_path, - extents=chunk_extents, - set_foreshore=self.drop_offshore_lidar, - ) - - # Use the coase DEM to sample any missing values - delayed_chunked_x.append( - dask.array.from_delayed( - delayed_elevation_over_chunk( - dim_x=dim_x, - dim_y=dim_y, - tile_points=coarse_points, - options=raster_options, - ), - shape=(len(dim_y), len(dim_x)), - dtype=geometry.RASTER_TYPE, - ) - ) - delayed_chunked_matrix.append(delayed_chunked_x) - - # Combine chunks into a array and replace missing values in dataset - elevations = dask.array.block(delayed_chunked_matrix) - self._dem["z"] = self._dem.z.where(~mask, elevations) - # Update the data source layer - self._dem["data_source"] = self._dem.data_source.where( - ~(mask & self._dem.z.notnull()), - self.SOURCE_CLASSIFICATION["coarse DEM"], - ) - class RoughnessDem(LidarBase): """A class to add a roughness (zo) layer to a hydrologically conditioned DEM.