Punkt in Polygon

Zur Vorstellung von Komponenten und Units für Lazarus
Antworten
Ekkehard
Beiträge: 12
Registriert: So 12. Feb 2023, 12:42

Punkt in Polygon

Beitrag von Ekkehard »

Hallo,
ich hatte gesucht aber nicht gefunden, deshalb etwas vorhandenes aus C in FPC übersetzt, welches ich gerne teile.
Die Quelle
https://web.archive.org/web/20130126163 ... usion.html
zeigt die Implementierung zweier Ansätze, des "crossing number algorithm" und des "winding number algorithm".
Ersterer ist die bekanntere Methode. Sie basiert auf der Überlegung, dass auf einer Geraden, die vom Testpunkt zu einem beliebigen, aber sicher außerhalb des Polygons liegenden Punkt führt, die Anzahl der Kreuzungen 0 oder gerade ist, wenn der Testpunkt außerhalb liegt, 1 oder ungerade ist, wenn er innerhalb liegt.
Die Methode hat ein Problem, wenn das Polygon deformiert ist und sich selber überlagert, dann liefert sie ggf. ein falsches Ergebnis.
Die zweitere Methode ist relativ frisch (2001), sie zählt die Windungszahl an den Grenzen des Polygons, allerdings nur deren Ergebnis, nämlich ob die Kante, die überschritten wird, nach "Oben" oder nach "Unten" führt. Heben sich die "Nach oben"-Kanten mit denen "Nach unten" auf, ist der Punkt außerhalb. Ist der Punkt innerhalb, so ist die Summe aus "Nach oben" und "Nach unten" ungleich null, kann aber bei deformierten Polygonen im Betrag auch größer als 1 werden.
Beide Algorithmen kommen ohne aufwendige Winkelberechnungen aus und sind entsprechend schnell.
Eine Erläuterung der Funktion findet sich in der englischsprachigen Wikipedia: https://en.wikipedia.org/wiki/Point_in_polygon
Es ist aber sehr sinnvoll zunächst an anderer Stelle, das das Polygon umgebende Rechteck zu bestimmen, damit lassen sich die Prüfungen für alle Punkte die außerhalb des Rechteckes liegen sehr stark verkürzen. Die ist insbesondere für komplexe und aus vielen Teilstücken bestehenden Polygone wichtig.
Die Lizenz für das Original ist Public Domain, dies gilt natürlich auch für die Übersetzung.

// This code may be freely used and modified for any purpose
// providing that this copyright notice is included with it.
// SoftSurfer makes no warranty for this code, and cannot be held
// liable for any real or imagined damage resulting from its use.
// Users of this code must verify correctness for their application.

Eine Anwendung sieht in etwa so aus:

Code: Alles auswählen

uses
(...)
  uWnPnPoly;  
(...)
procedure TForm1.Button1Click(Sender: TObject);
var
  poly : array of Double;
  rc : Integer;
  b : Boolean;
  ax,ay : Double;
  s : String;
begin
  SetLength(poly,5*2);
  //           X               Y
  poly[0] :=  0.0; poly[1] :=   0.0;
  poly[2] := 10.0; poly[3] :=   0.0;
  poly[4] := 10.0; poly[5] :=  10.0;
  poly[6] :=  0.0; poly[7] :=  10.0;
  poly[8] :=  0.0; poly[9] :=   0.0;
  ax := 5.0;
  ay := 5.0;
  rc := WindingPointInPolygon(ax,ay,
                              poly[0],
                              Length(poly) div 2,
                              False);
  b := CrossingNumberPointInPolygon(ax,ay,
                              poly[0],
                              Length(poly) div 2,
                              False);
  if b then
    s := 'Crossing result: Inside'
  else
    s := 'Crossing result: Outside';
  Label1.Caption := Format('%s , WindingResult: %d',[s,rc]);
end;
Dateianhänge
uwnpnpoly.pas
Diese Unit enthält eine fpc Implementation zweier "Punkt in Polygon"-Funktionen.
(10.29 KiB) 57-mal heruntergeladen

siro
Beiträge: 732
Registriert: Di 23. Aug 2016, 14:25
OS, Lazarus, FPC: Windows 11
CPU-Target: 64Bit
Wohnort: Berlin

Re: Punkt in Polygon

Beitrag von siro »

Moin Ekkehard,
gute Arbeit es in FPC zu codieren und zu veröffentlichen. Danke.

Bernd (Wennerer) hatte mal eine Unit "ptin"
mit diversen Funktionen zum Prüfen ob ein Punkt sich in einer Figur befindet.
Dort befindet sich unter anderem auch die Funktion PointInPoly.
Welche nochmals verbessert (beschleunigt) wurde.

Die Unit hatte er für seinen "Multibutton" erstellt.

viewtopic.php?p=118099#p118099
Grüße von Siro
Bevor ich "C" ertragen muß, nehm ich lieber Lazarus...

Ekkehard
Beiträge: 12
Registriert: So 12. Feb 2023, 12:42

Re: Punkt in Polygon

Beitrag von Ekkehard »

Moin siro,
Dort befindet sich unter anderem auch die Funktion PointInPoly.
Welche nochmals verbessert (beschleunigt) wurde.
Vielen Dank für den Link!
Jetzt war ich doch in Sorge mich durch nachlässige Suche blamiert zu haben, aber der Kelch ist glücklicherweise an mir vorrüber gegangen :)

Der verlinkte Code benutzt nämlich ArcTan2 um den Winkel des entsprechenden Teilstücks zu berechnen.
Trigonometrische Funktionen sind aber recht "teuer" und deshalb wo immer möglich zu vermeiden.
Und das macht der von mir gefundene und übersetzte Code, er kommt ohne diese Funktionen aus und benutzt nur Grundrechenarten.

So sollte "meine" Variante deutlich schneller unterwegs sein.
LG

Ekkehard
Beiträge: 12
Registriert: So 12. Feb 2023, 12:42

Re: Punkt in Polygon

Beitrag von Ekkehard »

Keine Software ohne Fehler :-/

Neben der Fehlerkorrektur in der Übersetzung des Algorithmus, habe ich auch die Parameter der Funktionen ergänzt, weil es für den Algorithmus einen Unterschied macht, ob der Wert für Y nach "oben" fällt (wie bei Anwendungen in TCanvas) oder steigt (wie bei Landkarten).

Jetzt sollte es aber passen.
Gruß Ekkehard
Dateianhänge
uwnpnpoly.pas
(12.09 KiB) 61-mal heruntergeladen

wp_xyz
Beiträge: 4889
Registriert: Fr 8. Apr 2011, 09:01

Re: Punkt in Polygon

Beitrag von wp_xyz »

Ekkehard hat geschrieben:
Do 8. Jun 2023, 12:04
weil es für den Algorithmus einen Unterschied macht, ob der Wert für Y nach "oben" fällt (wie bei Anwendungen in TCanvas) oder steigt (wie bei Landkarten).
Wirklich? Bei der "Crossing Point" Methode sollte die Zahl der Schnittpunkte doch egal sein, ob die horizontal Testgerade von einem noch oben oder nach unten laufenden Polygon-Abschnitt geschnitten wird. Und auch bei der "Non-zero Winding Rule" sollte es egal sein, ob einem nach oben laufenden Kurvenstück ein positiver oder ein negativer Beitrag zur Gesamt-Windungszahl zugeschrieben wird; es kommt nur darauf an, ob die Summe Null ist oder nicht.

Ekkehard
Beiträge: 12
Registriert: So 12. Feb 2023, 12:42

Re: Punkt in Polygon

Beitrag von Ekkehard »

Wirklich? Bei der "Crossing Point" Methode sollte die Zahl der Schnittpunkte doch egal sein, ob die horizontal Testgerade von einem noch oben oder nach unten laufenden Polygon-Abschnitt geschnitten wird. Und auch bei der "Non-zero Winding Rule" sollte es egal sein, ob einem nach oben laufenden Kurvenstück ein positiver oder ein negativer Beitrag zur Gesamt-Windungszahl zugeschrieben wird; es kommt nur darauf an, ob die Summe Null ist oder nicht.
Bedauerlicherweise sind meine mathematischen Fähigkeiten limitiert, so dass ich die Aussage nicht mathematisch widerlegen kann.
Aber praktischerweise, kann man den Unterschied ja austesten und anders erläutern.

Der Unterschied kommt zum Tragen, wenn die Frage gestellt wird, welcher der Punkte die das Polygon aufspannen, wenn er mit einem Testpunkt geprüft wird, als "inerhalb" oder "außerhalb" gesehen wird. Wenn man ein Rechteck definiert, so ist idR der linke und obere Rand Teil der Fläche, der linke und untere Rand dagegen nicht. Das erschließt sich unmittelbar, denn man möchte ja, dass sich zwei nebeneinander liegende Rechtecke
A ((0;0) (1;1))
B ((1;0) (2;1)),
die sich z.B. den Punkt (1;0) teilen, weder überlappen, noch dass es eine Lücke dazwischen gibt.
Diese Bedingung sollte der Algorithmus für Polygone auch beherrschen, mithin den linken und oberen Rand als zugehörig, den rechte und unteren Rand nicht. Wobei die Definition was oben, unten, rechts udn links ist, bei Polygonen sich aus den Winkeln ergibt.
(Natürlich kommt auf nicht-vertikalen bzw nicht-horizontalen Kanten des Polygons die Rechenungenauigkeit der ließkommaberechnung zum Tragen, nicht aber bei den Punkten, bzw. den vertikalen oder horizontalen Kanten. Dies ändert aber an den grundsätzlichen Überlegungen nichts.)

Zu meiner (nur kleinen) Überraschung machten beide Implementierungen das genau richtig, außer wenn die Y-Achse ein anderes Vorzeichen hat, dann ist das Ergebnis falsch.
Da der Einsatz der Algorithmen offensichtlich auf Polygone auf dem Bildschirm (Bitmaps, Canvas, etc.) abgestellt war, bei denen die Y-Werte nach "unten" zunehmen, ist diese Information intrinsisch vorhanden. Für meine Anwendung (Geokoordinaten, Grenzen) ist dagegen das Vorzeichen andersherum und relevant.

Lange Rede kurzer Sinn, ohne die Berücksichtigung der Achsenrichtung für die Y-Achse vertauschen sich die Zugehörigkeiten der entsprechenden Punkte zum Polygon, wenn die Y-Ache nach "oben" ansteigt.

Antworten