Przecięcie obszaru w Pythonie

Question

Mam kod, który pobiera warunek C jako dane wejściowe i oblicza rozwiązanie mojego problemu jako "dozwolony obszar" A na przestrzeni (x, y). Obszar ten składa się z kilku "rur", które są zdefiniowane przez 2 linie, których nigdy nie można przekroczyć.

Rezultat, którego szukam, musi spełniać warunki k {C1, .., Ck}, a zatem jest przecięciem S między obszarami k {A1, ..., Ak}.

Oto przykład z 2 warunkami (A1: zielony, 3 rurki: A2: fioletowy, 1 tuba); rozwiązanie S jest w kolorze czerwonym.

Jak znajdę S kiedy mam do czynienia z 4 obszarów wokół 10 tubki? (Ostateczna fabuła jest okropna!)

Potrzebowałbym móc ją wykreślić i znaleźć średnią współrzędną i wariancję punktów w S (wariancja każdej współrzędnej). [Jeśli istnieje skuteczny sposób sprawdzenia, czy punkt P należy do S, czy nie, użyję tylko metody Monte Carlo].

Idealnie, chciałbym też móc wprowadzić "zakazane rury", które usunęłbym z S [to może być nieco bardziej skomplikowane niż przecięcie S z zewnętrzną częścią mojego zabronionego obszaru, ponieważ dwie rurki z ten sam obszar może się przeciąć (nawet jeśli linie definiujące rurę nigdy nie przecinają się)].

---

Uwaga:

• Kod przechowuje również długość łuku linii.

• Linie są przechowywane jako tablice punktów (około 1000 punktów na linię). Dwie linie definiujące rurę niekoniecznie mają tę samą liczbę punktów, ale Python może interpolować WSZYSTKIE z nich jako funkcję ich długości łuku w ciągu 1 sekundy.

• Linie są funkcjami parametrycznymi (tzn. Nie można pisać y = f (x), ponieważ dozwolone są pionowe linie).

• Działka została zmontowana farbą, aby uzyskać wynik po prawej ... Niezbyt wydajna!

---

Edit:

• Nie wiem, w jaki sposób mogę korzystać plt.fill_between na stwardnienie skrzyżowaniu (można to zrobić tutaj 2 warunki, ale potrzebny jest kod zrobić to automatycznie, gdy jest zbyt wiele linii do oceny wzroku).

• Na razie po prostu generuję linie. Nie napisałem nic, aby znaleźć ostateczne rozwiązanie, ponieważ absolutnie nie wiem, która struktura jest najbardziej do tego przystosowana. [Jednak poprzednia wersja kodu była w stanie znaleźć punkty przecięcia między liniami 2 różnych rur, a ja planowałem przekazać je jako wielokąty, aby uzyskać kształt, ale to implikowało kilka innych problemów ..]

• Nie sądzę, żebym mógł to zrobić z sets: skanowanie całego obszaru (x, y) z wymaganą precyzją reprezentuje około 6e8 punktów ...[Linie mają tylko 1e3 punktów dzięki zmiennej wielkości kroku (dostosowuje się do krzywizny), ale cały problem jest dość duży]

Answer 1

Problem rozwiązany z Shapely!

Zdefiniowałem każdą rurkę jako Polygon, a obszar A to obiekt MultiPolygon zbudowany jako połączenie jego rur.

Następnie metoda oblicza szukane rozwiązanie (nakładanie się wszystkich obszarów).

Cała sprawa jest niemal natychmiastowa. Nie wiedziałem, że jest tak dobrze z dużymi obiektami [około 2000 punktów na probówkę, 10 rur na obszar, 4 obszary].

Dziękuję za pomoc! :)

Edytuj:

Przykład działania.

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() 

a fabuła: