panosz / alpha_shapes

A Python library for working with alpha shapes
MIT License
44 stars 8 forks source link

vertices missing from triangulation #3

Closed panosz closed 1 year ago

panosz commented 1 year ago
          Hi!

Thanks for attempting to fix it! Not sure why but it works for some sets of points but not others for me. I have a large data frame of species occurrences with their latitude and longitude coordinates and basically want to plot the alpha shapes of all the coordinates belonging to different decades.

for i, d in reversed(list(enumerate(decades))):
    decade_points = cane_toads[["decimalLongitude", "decimalLatitude"]] [cane_toads["decade"] == d]
    decade_points.drop_duplicates(subset=['decimalLongitude', 'decimalLatitude'], inplace=True)
    shaper = Alpha_Shaper(decade_points)
    alpha_opt, alpha_shape = shaper.optimize()
    print(alpha_opt)

Depending on the decade, some sets of points appear to work fine and produce an alpha shape while others result in the error below. I am being careful to remove duplicate coordinates from my dataframe but still get this error.

---------------------------------------------------------------------------
OptimizationFailure                       Traceback (most recent call last)
Cell In[11], line 5
      3 decade_points.drop_duplicates(subset=['decimalLongitude', 'decimalLatitude'], inplace=True)
      4 shaper = Alpha_Shaper(decade_points)
----> 5 alpha_opt, alpha_shape = shaper.optimize()
      6 print(alpha_opt)

File ~/Documents/ALA/galah-env/lib/python3.10/site-packages/alpha_shapes/alpha_shapes.py:153, in Alpha_Shaper.optimize(self)
    151 def optimize(self):
    152     # At least N//3 triangles are needed to connect N points.
--> 153     n_min = self._get_minimum_fully_covering_index_of_simplices()
    154     alpha_opt = 1 / np.sqrt(self._sorted_circumradii_sw()[n_min]) - 1e-10
    155     simplices = self._sorted_simplices()

File ~/Documents/ALA/galah-env/lib/python3.10/site-packages/alpha_shapes/alpha_shapes.py:149, in Alpha_Shaper._get_minimum_fully_covering_index_of_simplices(self)
    146         return n
    148 if uncovered_vertices:
--> 149     raise OptimizationFailure("Maybe there are duplicate points?")

OptimizationFailure: Maybe there are duplicate points?

_Originally posted by @due23 in https://github.com/panosz/alpha_shapes/issues/1#issuecomment-1411423427_

panosz commented 1 year ago

the point set creating the problem can be found here https://github.com/panosz/alpha_shapes/issues/1#issuecomment-1418333993

panosz commented 1 year ago

It looks like there are some vertices that are not included in the original triangulation

original triplot with the missing points highlighted: triplot_with_missing_points

and zooming in the uppermost missing point:

tri_plot_with_missing_points_zoomed Note the scaling in the axes. This vertex is very close to the one included, seen on the left.

There is a related comment in the documentation for scipy.spatial.Delaunay :

Unless you pass in the Qhull option “QJ”, Qhull does not guarantee that each input point appears as a vertex in the Delaunay triangulation. Omitted points are listed in the coplanar attribute.

I will attempt a temporary fix, but I think that I will eventually switch to scipy for creating the triangulation with the option QJ on.

panosz commented 1 year ago

https://github.com/panosz/alpha_shapes/commit/99b9a3c1e3c9287c8ceaa05d10bc8d3e973024d4 is a temporary fix, but as I said I will probably need a better solution. One problem I can see with the current solution is that, while it will compute the alpha shape, the shaper will probably not be able to be used as a triangulation for plotting or interpolating.