I’ve been writing some classes to calculate coverage and revisit in python, which uses a PolygonalFieldOfView coupled with long/lat points on the ground. I have considered using FieldOfViewDetectors however the implementation is slow.
Instead I have used the long/lat points output by FieldOfViewDetector. These appear to be wrapped around at -180/180 degrees of latitude, however I cannot find if this behaviour is documented anywhere, and it results in self-intersections of the polygon:
Is it possible to either have the option of getting the unwrapped polygon or separate polygons either side of the seam?
The classical way to handle this is to consider the first longitude l_0 in the boundary as a reference,
and then for each other longitude l_i to normalize it using MathUtils.normalizeAngle(l_i, l_0); rather than l_i itself. This works as long as the points are all confined within one hemisphere.
You have to define each contiguous regions first and then to use RegionFactory.union() from the Hipparchus library to combine them. Hipparchus knows how to handle regions that are not path-connected, or that have holes, or that contain the pole.
Beware about the order in which you define the boundary, it is important as the part that is on the left hand side as you travel along the boundary is the interior and the part on the right hand side is the exterior.
Latitude is limited to ±π/2, so there is no need for normalizing.
I’ve created some functions in which deal with the issue. Considering the polygon as a list of long/lat coordinate tuples/lists:
This first one corrects out the seam wrapping
def normalise_polygon(poly_list):
long_ref = poly_list[0][0]
long_norm = [math.degrees(math.atan(math.tan(
math.radians(x[0] - long_ref)))) + long_ref for x in poly_list]
return [(x[0], x[1][1]) for x in zip(long_norm, poly_list)]
This second one uses commands from the shapely library to split up an intersecting polygon into several separate polygons, and returns a list of polygons, which are made up of lists of long/lat pairs again.
from shapely.geometry import LineString
from shapely.ops import polygonize, unary_union
def de_intersect_polygon(poly):
# original data
ls = LineString(np.c_[poly])
# closed, non-simple
lr = LineString(ls.coords[:] + ls.coords[0:1])
mls = unary_union(lr)
mls.geom_type # MultiLineString’
polys = []
for polygon in polygonize(mls):
polys.append([(x[0], x[1]) for x in
zip(polygon.exterior.xy[0], polygon.exterior.xy[1])])
return polys
