:py:mod:`wolfhece.eikonal` ========================== .. py:module:: wolfhece.eikonal Module Contents --------------- .. py:function:: _evaluate_distance_and_data_first_order_iso(i, j, fixed: numpy.ndarray, where_compute: numpy.ndarray, dx_mesh: float, dy_mesh: float, times: numpy.ndarray, base_data: numpy.ndarray, test_data: numpy.ndarray, speed: float = 1.0) -> tuple[float, float] Evaluate the time and data using a first order isotropic method. :param i: (int): The row index. :param j: (int): The column index. :param fixed: (2D numpy array): A boolean array where True indicates fixed points. :param where_compute: (2D numpy array): A boolean array where True indicates cells to be included in computation. :param dx_mesh: (float): The mesh size in x direction. :param dy_mesh: (float): The mesh size in y direction. :param times: (2D numpy array): The time function. :param base_data: (2D numpy array): The base data that will propagate. :param test_data: (2D numpy array): The test data that will be used to validate propagation. :param speed: (float): The isotropic propagation speed. :return: tuple[float, float]: The time and data. .. py:function:: __solve_eikonal_with_data(sources: list[list[int, int]], where_compute: numpy.ndarray, base_data: numpy.ndarray, test_data: numpy.ndarray, speed: numpy.ndarray, dx_mesh: float, dy_mesh: float) -> numpy.ndarray Solve the Eikonal equation using the Fast Marching Method (FMM). Jit version of the function. The next one is the non-jit version which uses this. .. py:function:: _solve_eikonal_with_data(sources: list[tuple[int, int]], where_compute: numpy.ndarray = None, base_data: numpy.ndarray = None, test_data: numpy.ndarray = None, speed: numpy.ndarray = None, dx_mesh: float = 1.0, dy_mesh: float = 1.0) -> numpy.ndarray Solve the Eikonal equation using the Fast Marching Method (FMM). :param sources: (list of tuples): The coordinates of the source points. :param where_compute: (2D numpy array, optional): A boolean array where True indicates cells to be included in computation. :param base_data: (2D numpy array, optional): The base data that will propagate. :param test_data: (2D numpy array, optional): The test data that will be used to validate propagation. :param speed: (2D numpy array): The speed function. :param dx_mesh: (float, optional): The mesh size in x direction. :param dy_mesh: (float, optional): The mesh size in y direction. :return: 2D numpy array The solution to the Eikonal equation. .. py:function:: _extract_patch_slices_with_margin(shape, labels, margin=2) Extract the slices of the patches with a margin around the labels. :param shape: (tuple): The shape of the array. :param labels: (2D numpy array): The labels of the patches. :param margin: (int, optional): The margin around the labels. .. py:function:: _process_submatrix(args) Function to process a submatrix in a multiprocess. .. py:function:: _solve_eikonal_with_value_on_submatrices(where_compute: numpy.ndarray, base_data: numpy.ndarray, test_data: numpy.ndarray = None, speed: numpy.ndarray = None, multiprocess: bool = False, dx: float = 1.0, dy: float = 1.0, ignore_last_patches: int = 1, NoDataValue: float = 0.0) -> numpy.ndarray Propagate data inside the mask using the Fast Marching Method (FMM). "base_data" will be updated with the new values. :param where_compute: (2D numpy array): A boolean array where True indicates cells to be included in computation. :param base_data: (2D numpy array): The base data that will propagate. :param test_data: (2D numpy array, optional): The test data that will be used to validate propagation (we used upstream data only if base_data > test_data). :param speed: (2D numpy array, optional): The isotropic propagation speed. :param multiprocess: (bool, optional): Whether to use multiprocessing. :param dx: (float, optional): The mesh size in x direction. :param dy: (float, optional): The mesh size in y direction. :return: 2D numpy array The solution to the Eikonal equation. .. py:function:: count_holes(mask: numpy.ndarray = None) Count the number of holes in the mask. .. py:function:: inpaint_array(data: numpy.ndarray | numpy.ma.MaskedArray, where_compute: numpy.ndarray, test: numpy.ndarray, ignore_last_patches: int = 1, inplace: bool = True, dx: float = 1.0, dy: float = 1.0, NoDataValue: float = 0.0, multiprocess: bool = True) -> tuple[numpy.ndarray, numpy.ndarray, numpy.ndarray] Inpaint the array using the Fast Marching Method (FMM). Main idea: - We have a 2D array "data" that we want to inpaint. - We have a 2D array "test" that will be used to validate the inpainting. - We have a 2D array "where_compute" that indicates where to compute the inpainting. - We will use the FMM to propagate the data inside the mask. - We will update the data only if the new value is greater than the test data. - (We can ignore n greatest patches to avoid computation in some large areas.) :param data: (2D numpy array): The water level to inpaint. :param mask: (2D numpy array): The simulation's Digital Elevation Model (DEM). :param test: (2D numpy array, optional): The digital terrain model (DTM). :param ignore_last_patches: (int, optional): The number of last patches to ignore. :param inplace: (bool, optional): Whether to update the water_level in place. :param dx: (float, optional): The mesh size in x direction. :param dy: (float, optional): The mesh size in y direction. :param NoDataValue: (float, optional): The NoDataValue, used if mask is not explicitly provided (mask atribute or water_level as a Numpy MaskedArray). Default is 0. :param multiprocess: (bool, optional): Whether to use multiprocessing. .. py:function:: inpaint_waterlevel(water_level: numpy.ndarray | numpy.ma.MaskedArray, dem: numpy.ndarray, dtm: numpy.ndarray, ignore_last_patches: int = 1, inplace: bool = True, dx: float = 1.0, dy: float = 1.0, NoDataValue: float = 0.0, multiprocess: bool = True, epsilon: float = 0.001) -> tuple[numpy.ndarray, numpy.ndarray, numpy.ndarray] Inpaint the water level using the Fast Marching Method (FMM). Similar to the HOLES.EXE Fortran program. Assumptions: - The simulations are in a steady state. - The flood extent is limited by: - natural topography (where DEM == DTM) - buildings or blocks of buildings - protective structures The goal is to propagate the free surface elevations into the buildings. We calculate the difference between the DEM (including buildings and walls) and the DTM to identify where the buildings are. Specifically: - if it is natural topography, the differences will be zero or almost zero - if it is a building or anthropogenic structure, the differences are significant We set the elevations of the cells where the difference is zero to a value unreachable by the flood. Thus, only buildings in contact with the flood will be affected and filled. HOLES.EXE vs Python code: - In "holes.exe", we must provide "in", "mask", and "dem" files: - "in" is the water level - "dem" is the digital terrain model (DTM) associated with the simulation's topo-bathymetry (not the topo-bathymetry itself) - "mask" is the DTM where the buildings are identified and a large value (above the maximum water level) if the cell is not inside a building. - The patches are identified by the water level "0.0" in the "in" file. - FMM is used and new water levels are calculated and retained if greater than the "mask" value. - The DTM is only used in the final part when evaluating the new water depths (Z-DTM). - In Python, we must provide the water level, the DEM, and the DTM: - The patches will be identified by the buildings array (filtered DEM - DTM). - FMM is used and new water levels are calculated and retained if greater than the DTM value. - FMM is only propagated in the cells where "buildings = True". - Final results must be the same even if the algorithm is a little bit different. - Fortran code demands the "mask" file to be provided (pre-computed/modified by the user), but in Python, we calculate it automatically from the DEM and DTM. - We can use "inpaint_array" to be more flexible... by manually providing a "where_compute" and a "test" arrays :param water_level: (2D numpy array): The water level to inpaint. :param dem: (2D numpy array): The simulation's Digital Elevation Model (DEM). :param dtm: (2D numpy array, optional): The digital terrain model (DTM). :param ignore_last_patches: (int, optional): The number of last patches to ignore. :param inplace: (bool, optional): Whether to update the water_level in place. :param dx: (float, optional): The mesh size in x direction. :param dy: (float, optional): The mesh size in y direction. :param NoDataValue: (float, optional): The NoDataValue, used if mask is not explicitly provided (mask attribute or water_level as a Numpy MaskedArray). Default is 0. :param multiprocess: (bool, optional): Whether to use multiprocessing.