Intersection of an area in Python

Question

I have a code that accepts condition C as input, and computes the solution to my problem as an “allowed area” A in space (x, y). This area consists of several “tubes” that are defined by two lines that can never intersect.

The final result that I am looking for must satisfy k conditions {C1, .., Ck} and, therefore, is the intersection of S between k regions {A1, .., Ak}.

Here is an example with two conditions (A1: green, 3 handsets. A2: purple, 1 handset); S solution is red.

How can I find S when I deal with 4 sections of 10 tubes each? (The last plot is terrible!)

I would need to build it and find the average coordinate and the variance of the points in S (the variance of each coordinate). [If there is an effective way to find out if point P belongs to S or not, I just use the Monte Carlo method].

Ideally, Id would also like to implement “forbidden tubes” that I would remove from S [this can be a bit more complicated than crossing S with the limits of my restricted area, since two tubes from the same area can (even if the lines defining the tube, never intersect).

---

Note:

• The code also stores the length of the arc lines.

• Lines are saved as arrays of points (about 1000 points per line). The two lines defining the tube do not necessarily have the same number of points, but Python can interpolate ALL of them depending on the length of their arc in 1 second.

• Lines are parametric functions (i.e., we cannot write y = f (x), since lines are allowed vertically).

• The plot was edited with paint to get the result on the right ... Not very effective!

---

Edit:

• I don’t know how I can use plt.fill_between for multiple intersections (I can do it here for 2 conditions, but I need code to do this automatically when there are too many lines to evaluate eyes).

• For now, I'm just creating strings. I did not write anything to find the final solution, since I absolutely do not know which structure is most suitable for this. [However, the previous version of the code was able to find the intersection points between the lines from two different tubes, and I planned to transfer them as polygons to be beautiful, but this implied several other problems.]

• I don’t think I can do this with sets : scanning the entire area (x, y) with the required accuracy is about 6e8 points ... [Lines have only 1e3 points due to variable step size (adapts to curvature), but the whole problem is pretty great]

Answer 1

Problem solved with Shapely!

I defined each tube as a Polygon , and region A is a MultiPolygon object built as a union of its tubes.

Then, the intersection method calculates the solution I was looking for (overlap between all areas).

All this is almost instant. I did not know that slim was so good with large objects (about 2000 points per tube, 10 tubes per square, 4 areas).

Thanks for the help!:)

Edit:

Working example.

 import matplotlib.pyplot as plt import shapely from shapely.geometry import Polygon from descartes import PolygonPatch import numpy as np def create_tube(a,height): x_tube_up = np.linspace(-4,4,300) y_tube_up = a*x_tube_up**2 + height x_tube_down = np.flipud(x_tube_up) #flip for correct definition of polygon y_tube_down = np.flipud(y_tube_up - 2) points_x = list(x_tube_up) + list(x_tube_down) points_y = list(y_tube_up) + list(y_tube_down) return Polygon([(points_x[i], points_y[i]) for i in range(600)]) def plot_coords(ax, ob): x, y = ob.xy ax.plot(x, y, '+', color='grey') area_1 = Polygon() #First area, a MultiPolygon object for h in [-5, 0, 5]: area_1 = area_1.union(create_tube(2, h)) area_2 = Polygon() for h in [8, 13, 18]: area_2 = area_2.union(create_tube(-1, h)) solution = area_1.intersection(area_2) #What I was looking for ########## PLOT ########## fig = plt.figure() ax = fig.add_subplot(111) for tube in area_1: plot_coords(ax, tube.exterior) patch = PolygonPatch(tube, facecolor='g', edgecolor='g', alpha=0.25) ax.add_patch(patch) for tube in area_2: plot_coords(ax, tube.exterior) patch = PolygonPatch(tube, facecolor='m', edgecolor='m', alpha=0.25) ax.add_patch(patch) for sol in solution: plot_coords(ax, sol.exterior) patch = PolygonPatch(sol, facecolor='r', edgecolor='r') ax.add_patch(patch) plt.show() 

And the plot: