gempy icon indicating copy to clipboard operation
gempy copied to clipboard

KeyError when show_topography=True with layer not in lith_block (e.g., layer below resolution of blocks)

Open PeterBetlem opened this issue 2 years ago • 0 comments

Describe the bug Whenever a layer is too thin/substantially thinner than the resolution, its colour gets mapped but is not actually present in the lith_block data. Plotting this in 2D/3D is fine as long as the plot has Topography disabled, even if topography is applied to the generated surfaces. However, plotting in 3D with show_topography= True results in the following KeyError:

Log dump/KeyError
---------------------------------------------------------------------------
KeyError                                  Traceback (most recent call last)
Input In [63], in <cell line: 1>()
----> 1 gp.plot_3d(geo_model, plotter_type="basic", off_screen=False,
      2                       show_topography=True,
      3                       show_scalar=False,
      4                       show_lith=True, show_data=False)

File ~\.conda\envs\gemgis_deltaneset\lib\site-packages\gempy\plot\plot_api.py:336, in plot_3d(model, plotter_type, show_data, show_results, show_surfaces, show_lith, show_scalar, show_boundaries, show_topography, scalar_field, ve, kwargs_plot_structured_grid, kwargs_plot_topography, kwargs_plot_data, image, off_screen, **kwargs)
    333     gpv.plot_data(**kwargs_plot_data)
    335 if show_topography and model._grid.topography is not None:
--> 336     gpv.plot_topography(**kwargs_plot_topography)
    338 if ve is not None:
    339     gpv.p.set_scale(zscale=ve)

File ~\.conda\envs\gemgis_deltaneset\lib\site-packages\gempy\plot\vista.py:557, in GemPyToVista.plot_topography(self, topography, scalars, contours, clear, **kwargs)
    553 colors_rgb = pd.DataFrame(colors_rgb_.to_list(), index=colors_hex.index) * 255
    555 sel = np.round(self.model.solutions.geological_map[0]).astype(int)[0]
--> 557 scalars_val = numpy_to_vtk(colors_rgb.loc[sel], array_type=3)
    558 cm = mcolors.ListedColormap(list(self._get_color_lot(is_faults=True, is_basement=True)))
    559 rgb = True

File ~\.conda\envs\gemgis_deltaneset\lib\site-packages\pandas\core\indexing.py:931, in _LocationIndexer.__getitem__(self, key)
    928 axis = self.axis or 0
    930 maybe_callable = com.apply_if_callable(key, self.obj)
--> 931 return self._getitem_axis(maybe_callable, axis=axis)

File ~\.conda\envs\gemgis_deltaneset\lib\site-packages\pandas\core\indexing.py:1153, in _LocIndexer._getitem_axis(self, key, axis)
   1150     if hasattr(key, "ndim") and key.ndim > 1:
   1151         raise ValueError("Cannot index with multidimensional key")
-> 1153     return self._getitem_iterable(key, axis=axis)
   1155 # nested tuple slicing
   1156 if is_nested_tuple(key, labels):

File ~\.conda\envs\gemgis_deltaneset\lib\site-packages\pandas\core\indexing.py:1093, in _LocIndexer._getitem_iterable(self, key, axis)
   1090 self._validate_key(key, axis)
   1092 # A collection of keys
-> 1093 keyarr, indexer = self._get_listlike_indexer(key, axis)
   1094 return self.obj._reindex_with_indexers(
   1095     {axis: [keyarr, indexer]}, copy=True, allow_dups=True
   1096 )

File ~\.conda\envs\gemgis_deltaneset\lib\site-packages\pandas\core\indexing.py:1314, in _LocIndexer._get_listlike_indexer(self, key, axis)
   1311 else:
   1312     keyarr, indexer, new_indexer = ax._reindex_non_unique(keyarr)
-> 1314 self._validate_read_indexer(keyarr, indexer, axis)
   1316 if needs_i8_conversion(ax.dtype) or isinstance(
   1317     ax, (IntervalIndex, CategoricalIndex)
   1318 ):
   1319     # For CategoricalIndex take instead of reindex to preserve dtype.
   1320     #  For IntervalIndex this is to map integers to the Intervals they match to.
   1321     keyarr = ax.take(indexer)

File ~\.conda\envs\gemgis_deltaneset\lib\site-packages\pandas\core\indexing.py:1377, in _LocIndexer._validate_read_indexer(self, key, indexer, axis)
   1374     raise KeyError(f"None of [{key}] are in the [{axis_name}]")
   1376 not_found = list(ensure_index(key)[missing_mask.nonzero()[0]].unique())
-> 1377 raise KeyError(f"{not_found} not in index")

KeyError: '[21] not in index'

Indeed, there is no lithology with key 21 in the lith_block:

21 in geo_model.solutions.lith_block
> False

To Reproduce Provide detailed steps to reproduce the behavior:

See above.

Expected behavior

  1. Prior to calling numpy_to_vtk, make sure to check all lith_block IDs are actually in the lith_block
  2. Only apply numpy_to_vkt (and colouring) with the lith_block IDs that are present.
  3. Plot the topography with lith_block colours dynamically -> If a higher resolution topography would should more of the lith_blocks, than these should become visible after rerunning the models.

Desktop (please complete the following information):

  • OS: Windows 11
  • GemPy Version: '2.2.11'

PeterBetlem avatar May 19 '22 14:05 PeterBetlem