Kürzeste Entfernung zwischen einem Punkt und einem Linienabschnitt

Question

Ich benötige eine Grundfunktion, um die kürzeste Entfernung zwischen einem Punkt und einem Linienabschnitt zu finden. Fühlen Sie sich frei, die Lösung in jeder beliebigen Sprache zu schreiben; ich kann sie in das übersetzen, was ich benutze (Javascript).

BEARBEITEN: Mein Linienabschnitt wird durch zwei Endpunkte definiert. Also ist mein Linienabschnitt AB definiert durch die beiden Punkte A (x1, y1) und B (x2, y2). Ich versuche die Entfernung zwischen diesem Linienabschnitt und einem Punkt C (x3, y3) zu finden. Meine geometrischen Fähigkeiten sind eingerostet, daher sind die Beispiele, die ich gesehen habe, verwirrend, das muss ich zugeben.

Answer 1

Es sieht so aus, als hätten fast alle anderen auf StackOverflow eine Antwort abgegeben (bisher 23 Antworten), also hier ist mein Beitrag für C#. Dies basiert größtenteils auf der Antwort von M. Katz, die wiederum auf der Antwort von Grumdrig basiert.

   public struct MyVector
   {
      private readonly double _x, _y;

      // Konstruktor
      public MyVector(double x, double y)
      {
         _x = x;
         _y = y;
      }

      // Entfernung von diesem Punkt zu einem anderen Punkt, quadriert
      private double DistanceSquared(MyVector otherPoint)
      {
         double dx = otherPoint._x - this._x;
         double dy = otherPoint._y - this._y;
         return dx * dx + dy * dy;
      }

      // Finde die Entfernung von diesem Punkt zu einem Linienabschnitt (der nicht dasselbe ist wie von diesem
      //  Punkt zu irgendwo auf einer unendlichen Linie). Gibt auch den nächsten Punkt zurück.
      public double DistanceToLineSegment(MyVector lineSegmentPoint1, MyVector lineSegmentPoint2,
                                          out MyVector closestPoint)
      {
         return Math.Sqrt(DistanceToLineSegmentSquared(lineSegmentPoint1, lineSegmentPoint2, 
                          out closestPoint));
      }

      // Das gleiche wie oben, aber vermeiden Sie die Verwendung von Sqrt(), spart ein paar Nanosekunden in Fällen, in denen Sie nur mehrere Entfernungen vergleichen möchten, um die kleinste oder größte zu finden, aber die Entfernung nicht benötigen
      public double DistanceToLineSegmentSquared(MyVector lineSegmentPoint1, 
                                              MyVector lineSegmentPoint2, out MyVector closestPoint)
      {
         // Länge des Linienabschnitts berechnen (quadriert) und den besonderen Fall von koinzidenten Punkten behandeln
         double segmentLengthSquared = lineSegmentPoint1.DistanceSquared(lineSegmentPoint2);
         if (segmentLengthSquared < 1E-7f)  // Beliebter "gut genug für Regierungsarbeit" Wert
         {
            closestPoint = lineSegmentPoint1;
            return this.DistanceSquared(closestPoint);
         }

         // Verwenden Sie die magische Formel, um die "Projektion" dieses Punktes auf die unendliche Linie zu berechnen
         MyVector lineSegment = lineSegmentPoint2 - lineSegmentPoint1;
         double t = (this - lineSegmentPoint1).DotProduct(lineSegment) / segmentLengthSquared;

         // Behandeln Sie die beiden Fälle, in denen die Projektion nicht auf dem Linienabschnitt liegt, und den Fall, in dem
         //  die Projektion auf dem Segment liegt
         if (t <= 0)
            closestPoint = lineSegmentPoint1;
         else if (t >= 1)
            closestPoint = lineSegmentPoint2;
         else 
            closestPoint = lineSegmentPoint1 + (lineSegment * t);
         return this.DistanceSquared(closestPoint);
      }

      public double DotProduct(MyVector otherVector)
      {
         return this._x * otherVector._x + this._y * otherVector._y;
      }

      public static MyVector operator +(MyVector leftVector, MyVector rightVector)
      {
         return new MyVector(leftVector._x + rightVector._x, leftVector._y + rightVector._y);
      }

      public static MyVector operator -(MyVector leftVector, MyVector rightVector)
      {
         return new MyVector(leftVector._x - rightVector._x, leftVector._y - rightVector._y);
      }

      public static MyVector operator *(MyVector aVector, double aScalar)
      {
         return new MyVector(aVector._x * aScalar, aVector._y * aScalar);
      }

      // Hinzugefügt mit ReSharper aufgrund von CodeAnalysis-Mahnungen

      public bool Equals(MyVector other)
      {
         return _x.Equals(other._x) && _y.Equals(other._y);
      }

      public override bool Equals(object obj)
      {
         if (ReferenceEquals(null, obj)) return false;
         return obj is MyVector && Equals((MyVector) obj);
      }

      public override int GetHashCode()
      {
         unchecked
         {
            return (_x.GetHashCode()*397) ^ _y.GetHashCode();
         }
      }

      public static bool operator ==(MyVector left, MyVector right)
      {
         return left.Equals(right);
      }

      public static bool operator !=(MyVector left, MyVector right)
      {
         return !left.Equals(right);
      }
   }

Und hier ist ein kleines Testprogramm.

   public static class JustTesting
   {
      public static void Main()
      {
         Stopwatch stopwatch = new Stopwatch();
         stopwatch.Start();

         for (int i = 0; i < 10000000; i++)
         {
            TestIt(1, 0, 0, 0, 1, 1, 0.70710678118654757);
            TestIt(5, 4, 0, 0, 20, 10, 1.3416407864998738);
            TestIt(30, 15, 0, 0, 20, 10, 11.180339887498949);
            TestIt(-30, 15, 0, 0, 20, 10, 33.541019662496844);
            TestIt(5, 1, 0, 0, 10, 0, 1.0);
            TestIt(1, 5, 0, 0, 0, 10, 1.0);
         }

         stopwatch.Stop();
         TimeSpan timeSpan = stopwatch.Elapsed;
      }

      private static void TestIt(float aPointX, float aPointY, 
                                 float lineSegmentPoint1X, float lineSegmentPoint1Y, 
                                 float lineSegmentPoint2X, float lineSegmentPoint2Y, 
                                 double expectedAnswer)
      {
         // Katz
         double d1 = DistanceFromPointToLineSegment(new MyVector(aPointX, aPointY), 
                                              new MyVector(lineSegmentPoint1X, lineSegmentPoint1Y), 
                                              new MyVector(lineSegmentPoint2X, lineSegmentPoint2Y));
         Debug.Assert(d1 == expectedAnswer);

         /*
         // Katz using squared distance
         double d2 = DistanceFromPointToLineSegmentSquared(new MyVector(aPointX, aPointY), 
                                              new MyVector(lineSegmentPoint1X, lineSegmentPoint1Y), 
                                              new MyVector(lineSegmentPoint2X, lineSegmentPoint2Y));
         Debug.Assert(Math.Abs(d2 - expectedAnswer * expectedAnswer) < 1E-7f);
          */

         /*
         // Matti (optimiert)
         double d3 = FloatVector.DistanceToLineSegment(new PointF(aPointX, aPointY), 
                                                new PointF(lineSegmentPoint1X, lineSegmentPoint1Y), 
                                                new PointF(lineSegmentPoint2X, lineSegmentPoint2Y));
         Debug.Assert(Math.Abs(d3 - expectedAnswer) < 1E-7f);
          */
      }

      private static double DistanceFromPointToLineSegment(MyVector aPoint, 
                                             MyVector lineSegmentPoint1, MyVector lineSegmentPoint2)
      {
         MyVector closestPoint;  // Nicht verwendet
         return aPoint.DistanceToLineSegment(lineSegmentPoint1, lineSegmentPoint2, 
                                             out closestPoint);
      }

      private static double DistanceFromPointToLineSegmentSquared(MyVector aPoint, 
                                             MyVector lineSegmentPoint1, MyVector lineSegmentPoint2)
      {
         MyVector closestPoint;  // Nicht verwendet
         return aPoint.DistanceToLineSegmentSquared(lineSegmentPoint1, lineSegmentPoint2, 
                                                    out closestPoint);
      }
   }

Wie Sie sehen können, habe ich versucht, den Unterschied zwischen der Verwendung der Version, die die Methode Sqrt() vermeidet und der normalen Version zu messen. Meine Tests zeigen, dass Sie vielleicht etwa 2,5% sparen können, aber ich bin mir nicht einmal sicher darüber - die Variationen innerhalb der verschiedenen Testläufe waren in derselben Größenordnung. Ich habe auch versucht, die Version zu messen, die von Matti veröffentlicht wurde (plus einer offensichtlichen Optimierung), und diese Version scheint etwa 4% langsamer zu sein als die Version, die auf dem Katz/Grumdrig-Code basiert.

Bearbeitet: Übrigens habe ich auch versucht, eine Methode zu messen, die die Entfernung zu einer unendlichen Linie (nicht zu einem Linienabschnitt) mit einem Kreuzprodukt (und einem Sqrt()) findet, und sie ist etwa 32% schneller.

Answer 2

Eine 2D- und 3D-Lösung

Betrachten Sie eine Basisänderung, bei der das Liniensegment zu (0, 0, 0)-(d, 0, 0) und der Punkt zu (u, v, 0) wird. Der kürzeste Abstand liegt in dieser Ebene und wird gegeben durch

    u  0 -> d(A, C)
0  u  d -> |v|
d  u     -> d(B, C)

(der Abstand zu einem der Endpunkte oder zur Stützlinie, je nach der Projektion auf die Linie. Der Iso-Abstandsort besteht aus zwei Halbkreisen und zwei Linienabschnitten.)

In obigem Ausdruck ist d die Länge des Segments AB, und u, v sind jeweils das Skalarprodukt und das (Betrag des) Kreuzprodukts von AB/d (Einheitsvektor in Richtung AB) und AC. Folglich vektoriell,

AB.AC  0             -> |AC|
    0  AB.AC  AB²   -> |ABxAC|/|AB|
          AB²  AB.AC -> |BC|

Answer 3

Hier ist devnullicus' C++ Version in C# konvertiert. Für meine Implementierung musste ich den Schnittpunkt kennen und habe seine Lösung als gut funktionierend empfunden.

public static bool PointSegmentDistanceSquared(PointF point, PointF lineStart, PointF lineEnd, out double distance, out PointF intersectPoint)
{
    const double kMinSegmentLenSquared = 0.00000001; // anpassen. Wenn Sie float verwenden, möchten Sie wahrscheinlich etwas wie 0.000001f
    const double kEpsilon = 1.0E-14; // anpassen. Wenn Sie float verwenden, möchten Sie wahrscheinlich etwas wie 1E-7f
    double dX = lineEnd.X - lineStart.X;
    double dY = lineEnd.Y - lineStart.Y;
    double dp1X = point.X - lineStart.X;
    double dp1Y = point.Y - lineStart.Y;
    double segLenSquared = (dX * dX) + (dY * dY);
    double t = 0.0;

    if (segLenSquared >= -kMinSegmentLenSquared && segLenSquared <= kMinSegmentLenSquared)
    {
        // Das Segment ist ein Punkt.
        intersectPoint = lineStart;
        t = 0.0;
        distance = ((dp1X * dp1X) + (dp1Y * dp1Y));
    }
    else
    {
        // Projektion einer Linie von p zum Segment [p1,p2]. Durch Betrachten der Linie
        // die das Segment verlängert, parametrisiert als p1 + (t * (p2 - p1)),
        // finden wir die Projektion des Punktes p auf die Linie. 
        // Es fällt dort, wo t = [(p - p1) . (p2 - p1)] / |p2 - p1|^2
        t = ((dp1X * dX) + (dp1Y * dY)) / segLenSquared;
        if (t < kEpsilon)
        {
            // schneidet sich bei oder zu "links" des ersten Segmentpunkts (lineStart.X, lineStart.Y). Wenn t ungefähr 0,0 ist, dann
            // der Schnittpunkt ist bei p1. Wenn t kleiner ist als das, gibt es keinen Schnittpunkt (d.h. p liegt nicht innerhalb
            // der 'Grenzen' des Segments)
            if (t > -kEpsilon)
            {
                // schnittpunkt bei 1. Segmentpunkt
                t = 0.0;
            }
            // setzen Sie unseren 'Schnittpunkt' auf p1.
            intersectPoint = lineStart;
            // Hinweis: Wenn Sie den tatsächlichen Schnittpunkt dort haben wollten, wo sich die projizierten Linien schneiden würden, wenn
            // wir PointLineDistanceSquared machen würden, dann wäre intersectPoint.X (lineStart.X + (t * dx)) und intersectPoint.Y (lineStart.Y + (t * dy)).
        }
        else if (t > (1.0 - kEpsilon))
        {
            // schneidet sich bei oder zu "rechts" des zweiten Segmentpunkts (lineEnd.X, lineEnd.Y). Wenn t ungefähr 1,0 ist, dann
            // der Schnittpunkt ist bei p2. Wenn t größer ist als das, gibt es keinen Schnittpunkt (d.h. p liegt nicht innerhalb
            // der 'Grenzen' des Segments)
            if (t < (1.0 + kEpsilon))
            {
                // schnittpunkt bei 2. Segmentpunkt
                t = 1.0;
            }
            // setzen Sie unseren 'Schnittpunkt' auf p2.
            intersectPoint = lineEnd;
            // Hinweis: Wenn Sie den tatsächlichen Schnittpunkt dort haben wollten, wo sich die projizierten Linien schneiden würden, wenn
            // wir PointLineDistanceSquared machen würden, dann wäre intersectPoint.X (lineStart.X + (t * dx)) und intersectPoint.Y (lineStart.Y + (t * dy)).
        }
        else
        {
            // Die Projektion des Punktes auf den Punkt auf dem Segment, der senkrecht steht, war erfolgreich und der Punkt
            // befindet sich 'innerhalb' der Grenzen des Segments. Setzen Sie den Schnittpunkt als diesen projizierten Punkt.
            intersectPoint = new PointF((float)(lineStart.X + (t * dX)), (float)(lineStart.Y + (t * dY)));
        }
        // geben Sie den quadratischen Abstand von p zum Schnittpunkt zurück. Beachten Sie, dass wir den quadrierten Abstand zurückgeben
        // als Optimierung, weil Sie in vielen Fällen nur relative Abstände vergleichen müssen und die quadrierten Werte
        // dafür gut funktionieren. Wenn Sie den tatsächlichen Abstand möchten, nehmen Sie einfach die Quadratwurzel dieses Wertes.
        double dpqX = point.X - intersectPoint.X;
        double dpqY = point.Y - intersectPoint.Y;

        distance = ((dpqX * dpqX) + (dpqY * dpqY));
    }

    return true;
}

Answer 4

Sehen Sie sich das Matlab GEOMETRY-Toolbox auf der folgenden Website an: http://people.sc.fsu.edu/~jburkardt/m_src/geometry/geometry.html

Strg + F und tippen Sie "segment", um Funktionen im Zusammenhang mit Liniensegmenten zu finden. Die Funktionen "segment_point_dist_2d.m" und "segment_point_dist_3d.m" sind das, was Sie benötigen.

Die GEOMETRY-Codes sind in einer C-Version und einer C++-Version und einer FORTRAN77-Version und einer FORTRAN90-Version und einer MATLAB-Version verfügbar.

Answer 5

Ich habe einen interaktiven Desmos-Graphen erstellt, um zu zeigen, wie das erreicht werden kann:

https://www.desmos.com/calculator/kswrm8ddum

Der rote Punkt ist A, der grüne Punkt ist B, und der Punkt C ist blau. Du kannst die Punkte im Graphen verschieben, um die Werte zu ändern. Links ist der Wert 's' der Parameter des Streckenabschnitts (d.h. s = 0 bedeutet der Punkt A und s = 1 bedeutet der Punkt B). Der Wert 'd' ist der Abstand vom dritten Punkt zur Linie durch A und B.

BEARBEITEN:

Kleiner Spaßfakt: Die Koordinate (s, d) ist die Koordinate des dritten Punktes C im Koordinatensystem, in dem AB die Einheits-x-Achse ist und die Einheits-y-Achse senkrecht zu AB verläuft.