Wielokąty w międzynarodowej linii danych [-180 .. + 180 długości geograficznej]

10

Próbuję generować wielokąty dla satelitarnych pokosów orbit. Do tej pory mam metodę generowania dwóch linii, które reprezentują krawędź każdego pokosu w [łac, długo]. Niektóre pokosy przekraczają międzynarodową linię danych i dlatego obejmują:

pokos owija się wokół

Udało mi się to rozwiązać za pomocą ogr2ogr -wrapdateline:

ogr2ogr -wrapdateline  -f "ESRI Shapefile" test.shp orbits.shp

Co prawdopodobnie dzieli linie

Chcę teraz móc generować wielokąty we wnętrzu obu linii. Na przykład w przypadku, gdy jedna krawędź pokosu przecina linię danych, wielokąt wypełnia się, gdy pojawia się po drugiej stronie, na przykład:

napełnić

Potrzebuję metody zautomatyzowanej, ponieważ muszę często powtarzać to zadanie. Najlepiej w pythonie, ponieważ tak wygenerowałem linie. Oto dwa pliki kształtów zawierające linie: wraparound ; naprawiono linię danych

python ogr wgs84 shapely antimeridian James
źródło
Dodatkowe pomysły można znaleźć w powiązanych wątkach na gis.stackexchange.com/questions/429 i gis.stackexchange.com/questions/18562 . Pomysły przedstawione na gis.stackexchange.com/questions/17788 mogą być również pomocne. Zastanawiam się jednak, co rozumiesz przez „wnętrze”: te wielokąty nie są dobrze zdefiniowane, więc przynajmniej musisz podać informacje, aby wskazać (a), która strona każdej polilinii jest uważana za „wnętrze” i (b) jak odciąć je w pobliżu biegunów.
whuber

Odpowiedzi:

3

Możesz zbudować niestandardową projekcję mercatora wyśrodkowaną w przybliżeniu na środku pokosu. Na przykład użyj dla pokosu 25:

+proj=merc +lon_0=-140 +k=1 +x_0=0 +y_0=0 +ellps=WGS84 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs

W tej projekcji pokos nie jest łamany przez linię danych. Możesz utworzyć wielokąt z linii.

wprowadź opis zdjęcia tutaj

Następnie utwórz wycięty wielokąt między -179,95 ° E a 179,95 ° E w EPSG: 4326:

Nr;WKT
1;POLYGON ((-179.95 89, 179.95 89, 179.95 -89, -179.95 -89, -179.95 89))

Odrzuć go również do własnego CRS i odejmij od wielokąta pokosu.

Po ponownym rzutowaniu z powrotem do EPSG: 4326 pokos jest poprawnie dzielony przez linię danych:

wprowadź opis zdjęcia tutaj

Kontynuuj z wszystkimi pokosami, które przecinają linię danych.

AndreJ
źródło
1

Przepisałbym proces generowania linii pokosu, aby rozpocząć i zakończyć w tej samej ciągłej przestrzeni podskórnej. tzn. jeśli linia zaczyna się przy 170 ° i kończy przy -170 °, przepisałbym proces, aby zakończyć przy 190 ° zamiast owijania przy -180,180

Następnie możesz tworzyć nieprzerwane wielokąty między liniami.

Następnie użyj procesu przycinania, aby podzielić wielokąty na linii 180, -180 i przesunąć dowolne części, które leżą poza przestrzenią -180,180, dodając lub odejmując odpowiednio 360 °.

Po prostu zrób to wszystko, zanim zapiszesz go z konkretną projekcją / punktem odniesienia

Mr Purple
źródło
1

Dzięki @AndreJ za ten pomysł, użycie Django GEOS API tutaj jest prostym rozwiązaniem, które pozwala uniknąć konieczności ponownego wyświetlania czegokolwiek:

1) Utwórz MultiPolygon, który graniczy z linią danych:

from django.contrib.gis.geos.collections import MultiPolygon, LinearRing, Polygon
box1 = ((180.0, 89), (179.95, 89), (179.95, -89), (180.0, -89), (180.0, 89))
box2 = ((-180.0, 89), (-179.95, 89), (-179.95, -89), (-180.0, -89), (-180.0, 89))
poly1 = Polygon(LinearRing(box1))
poly2 = Polygon(LinearRing(box2))
poly = MultiPolygon(poly1, poly2)

2) Jeśli przecinająca się geometria przecina się, zwróć różnicę: 

from django.contrib.gis.geos.geometry import GEOSGeometry
geometry = GEOSGeometry(WKT)  # WKT is your polygon in WKT string format
if geometry.intersects(poly):
    print("Geometry crosses dateline... splitting")
    geometry = geometry.difference(poly) # clip with dateline polygon

Wynik pokazano w następujący sposób:

Przed

Po

MarMat
źródło