TaipanRex / pyvisgraph

Given a list of simple obstacle polygons, build the visibility graph and find the shortest path between two points
MIT License
221 stars 45 forks source link

Simplifying polygons #7

Open TaipanRex opened 8 years ago

TaipanRex commented 8 years ago

Implement algorithm to reduce number of points in a polygon:

https://bost.ocks.org/mike/simplify/

http://stackoverflow.com/questions/10558299/visvalingam-whyatt-polyline-simplification-algorithm-clarification

lucasloisp commented 8 years ago

I began implementing this but as I thought more about it... were should this be implemented?

It might be possible to simplify the whole graph, for example (self.graph), but perhaps we should simplify every polygon in self.polygon (or choose which one to).

Wouldn't it be useful to add a separate class Polygon to manage each one independent of the graph?

TaipanRex commented 8 years ago

Thanks for the interest!

The reason I did not implement a separate Polygon class is that the visible vertices algorithm doesn't care about polygons, only Point's and Edge's. So from a performance standpoint (I use fairly large polygon obstacles), its faster to use a dictionary where the key is a Point and value is a set of Edge's. This is in a similar spirit to the canonical way of implementing graphs in Python, see Guido's essay. I do store the polygons in self.polygon, as its needed for point_in_polygon checks.

My thinking for implementing simplify_polygon would be as follows: I would not do it as a Graph method that runs on self.graph, because think of the way the visible vertices algorithm uses the obstacle graph - it uses it once to build the visibility graph. Once the visibility graph is built, you would not update or change/simplify the obstacle graph, because that would break the visibility graph.

You will always run simplify_polygon before building the visibility graph with VisGraph.build. Therefore, I consider simplify_polygon more of a helper function you access through VisGraph. You could implement it as a function in graph.py, _simplify_polygon(polys, resolution), where polys is a list of polygons (where each polygon is a list of Point's). A user would call it through VisGraph.simplify_polygon.

You can use it stand alone to test out and visualize your simplification:

import pyvisgraph as vg
w = vg.VisGraph()
polys = ... # Obstacle polygons from some shapefile or other input
simple = w.simplify_polygon(polys, 0.5)
# draw the simplified polygon with matplotlib or similar to verify your simplification

Or you could use it directly when building the visibility graph:

import pyvisgraph as vg
w = vg.VisGraph()
polys = ... # Obstacle polygons from some shapefile or other input
w.build(w.simplify_polygon(polys, 0.5))

Thats how I am thinking, but you might have better ideas!

lucasloisp commented 8 years ago

Sure, I thought about the performance cost but I guessed the cost wouldn''t be much and the benefit of abstraction could facilitate the implementation of the simplify algorithm.

However, it is true that changing the polygon would break the visibility graph, and so I implemented the function (simplify_polygon) on VisGraph. And, given that a polygon is simply a list of Point, instead of a list of Edge, it is much simpler to manipulate them (removing the need for a Polygon class abstraction).

My pull request will be coming soon.

One silly question. I assume we define "resolution" to mean len(new_poly)/len(old_poly), right?

TaipanRex commented 8 years ago

I haven't looked into Visvalingam's algorithm yet, so I am not sure if you specify a resolution like you defined or if you specify a triangle threshold. I just wrote resolution in my above post to illustrate that some value is given to specify how large of a simplification you want.

Btw, the shapefile I use for testing is available here. I use GSHHS_c_L1.shp, which is the lowest resolution. I use pyshp to read the shapefile and get the obstacle graph in the following way:

import shapefile
import pyvisgraph as vg
sf = shapefile.Reader("GSHHS_c_L1.dbf")
shapes = sf.shapes()
polys = []
polys.append([vg.Point(p[0], p[1]) for p in shapes[0].points])

This is for getting just one of the shapes, if you want all you just loop on shape in shapes. Note that in that shapefile, shapes[0] is the Eurasian continent, with 1003 points. shapes[3] is North America, 667 points.