CODEWITHSUNDEEP
c#mathgeometry
VOX
VOX

Reputation: 2923

get closest point to a line

I'd like to have a straight forward C# function to get a closest point (from a point P) to a line-segment, AB. An abstract function may look like this. I've search through SO but not found a usable (by me) solution.

public Point getClosestPointFromLine(Point A, Point B, Point P);

Upvotes: 40

Views: 63821

Answers (13)

abishek arvind
abishek arvind

Reputation: 1

In case anybody looking for python implementation, here is the code:

p1 and p2 is line, p3 is the point

def p4(p1, p2, p3):
     x1, y1 = p1
     x2, y2 = p2
     x3, y3 = p3
     dx, dy = x2-x1, y2-y1
     det = dx*dx + dy*dy
     a = (dy*(y3-y1)+dx*(x3-x1))/det
     x= x1+a*dx, y1+a*dy
     # print(x)
     if x[0]<x1 or x[1]<y1:
         return p1
     elif x[0]>x2 or x[1]>y2:
         return p2
     else:
         return x

This was taken from another thread and modified a little. Python: point on a line closest to third point

Upvotes: 0

Bexhet Islamaj
Bexhet Islamaj

Reputation: 1

This is the right algorythm to get nearest point of a segment from a point(Tested)(vb.net)

s2 = ClosestPointToSegment(point_x, Point_y, Segment_start_x, Segment_start_y, Segment_end_X, Segment_end_Y)

Public Shared Function DistanceTo(x1 As Double, y1 As Double, x2 As Double, y2 As Double) As Double
    Return Math.Sqrt(Math.Pow(x1 - x2, 2) + Math.Pow(y1 - y2, 2))
End Function


Public Shared Function DistanceTo(point_x As Double, point_y As Double, lineStart_x As Double, lineStart_y As Double, lineEnd_x As Double, lineEnd_y As Double) As Double
    Dim tI As Double = ((lineEnd_x - lineStart_x) * (point_x - lineStart_x) + (lineEnd_y - lineStart_y) * (point_y - lineStart_x)) / Math.Pow(DistanceTo(lineStart_x, lineStart_y, lineEnd_x, lineEnd_y), 2)
    Dim dP As Double = ((lineEnd_x - lineStart_x) * (point_y - lineStart_y) - (lineEnd_y - lineStart_y) * (point_x - lineStart_x)) / DistanceTo(lineStart_x, lineStart_y, lineEnd_x, lineEnd_y)

    If tI >= 0R AndAlso tI <= 1.0R Then
        Return Math.Abs(dP)
    Else
        Return Math.Min(DistanceTo(point_x, point_y, lineStart_x, lineStart_y), DistanceTo(point_x, point_y, lineEnd_x, lineEnd_y))
    End If
End Function
Private Shared Function ClosestPointToSegment(P_x As Double, p_y As Double, A_x As Double, a_y As Double, B_x As Double, b_y As Double) As Double()
    Dim a_to_p As PointF = New PointF(), a_to_b As PointF = New PointF()
    Dim rikthex As Double, rikthey As Double
    Dim s1(1) As Double
    Dim p1_v1_X As Double, p1_v1_y As Double, distanca1 As Double, distanca2 As Double
    a_to_p.X = P_x - A_x
    a_to_p.Y = p_y - a_y
    a_to_b.X = B_x - A_x
    a_to_b.Y = b_y - a_y
    Dim atb2 As Single = a_to_b.X * a_to_b.X + a_to_b.Y * a_to_b.Y
    Dim atp_dot_atb As Single = a_to_p.X * a_to_b.X + a_to_p.Y * a_to_b.Y
    Dim t As Single = atp_dot_atb / atb2
    rikthex = A_x + a_to_b.X * t
    rikthey = a_y + a_to_b.Y * t
    If A_x > B_x Then
        If rikthex < A_x And rikthex > B_x Then 'pika duhet ne rregulll
            If a_y > b_y Then
                If rikthey < a_y And rikthey > b_y Then 'pika duhet ne rregulll

                Else
                    distanca1 = DistanceTo(P_x, p_y, A_x, a_y)
                    distanca2 = DistanceTo(P_x, p_y, B_x, b_y)
                    If distanca1 < distanca2 Then
                        rikthex = A_x
                        rikthey = a_y
                    Else
                        rikthex = B_x
                        rikthey = b_y
                    End If

                End If
            Else
                If rikthey > a_y And rikthey < b_y Then 'pika duhet ne rregulll

                Else
                    distanca1 = DistanceTo(P_x, p_y, A_x, a_y)
                    distanca2 = DistanceTo(P_x, p_y, B_x, b_y)
                    If distanca1 < distanca2 Then
                        rikthex = A_x
                        rikthey = a_y
                    Else
                        rikthex = B_x
                        rikthey = b_y
                    End If

                End If

            End If
        Else
            distanca1 = DistanceTo(P_x, p_y, A_x, a_y)
            distanca2 = DistanceTo(P_x, p_y, B_x, b_y)
            If distanca1 < distanca2 Then
                rikthex = A_x
                rikthey = a_y
            Else
                rikthex = B_x
                rikthey = b_y
            End If
        End If
    Else
        If rikthex > A_x And rikthex < B_x Then 'pika duhet ne rregulll
            If a_y > b_y Then
                If rikthey < a_y And rikthey > b_y Then 'pika duhet ne rregulll

                Else
                    distanca1 = DistanceTo(P_x, p_y, A_x, a_y)
                    distanca2 = DistanceTo(P_x, p_y, B_x, b_y)
                    If distanca1 < distanca2 Then
                        rikthex = A_x
                        rikthey = a_y
                    Else
                        rikthex = B_x
                        rikthey = b_y
                    End If

                End If
            Else
                If rikthey > a_y And rikthey < b_y Then 'pika duhet ne rregulll

                Else
                    distanca1 = DistanceTo(P_x, p_y, A_x, a_y)
                    distanca2 = DistanceTo(P_x, p_y, B_x, b_y)
                    If distanca1 < distanca2 Then
                        rikthex = A_x
                        rikthey = a_y
                    Else
                        rikthex = B_x
                        rikthey = b_y
                    End If

                End If

            End If
        Else
            distanca1 = DistanceTo(P_x, p_y, A_x, a_y)
            distanca2 = DistanceTo(P_x, p_y, B_x, b_y)
            If distanca1 < distanca2 Then
                rikthex = A_x
                rikthey = a_y
            Else
                rikthex = B_x
                rikthey = b_y
            End If
        End If
    End If
    s1(0) = rikthex
    s1(1) = rikthey
    Return s1

End Function

Upvotes: 0

Dave_cz
Dave_cz

Reputation: 1210

Here are extension methods that should do the trick:

public static double DistanceTo(this Point from, Point to)
    {
        return Math.Sqrt(Math.Pow(from.X - to.X, 2) + Math.Pow(from.Y - to.Y, 2));
    }

public static double DistanceTo(this Point point, Point lineStart, Point lineEnd)
    {
        double tI = ((lineEnd.X - lineStart.X) * (point.X - lineStart.X) + (lineEnd.Y - lineStart.Y) * (point.Y - lineStart.Y)) / Math.Pow(lineStart.DistanceTo(lineEnd), 2);
        double dP = ((lineEnd.X - lineStart.X) * (point.Y - lineStart.Y) - (lineEnd.Y - lineStart.Y) * (point.X - lineStart.X)) / lineStart.DistanceTo(lineEnd);

        if (tI >= 0d && tI <= 1d)
            return Math.Abs(dP);
        else
            return Math.Min(point.DistanceTo(lineStart), point.DistanceTo(lineEnd));
    }

Then just call:

P.DistanceTo(A, B);

To get distance of point "P" from line |AB|. It should be easy to modify this for PointF.

Finding the closest point then is just a matter of searching minimal distance. LINQ has methods for this.

Upvotes: 2

Justin L.
Justin L.

Reputation: 13600

Here's Ruby disguised as Pseudo-Code, assuming Point objects each have a x and y field.

def GetClosestPoint(A, B, P)

  a_to_p = [P.x - A.x, P.y - A.y]     # Storing vector A->P
  a_to_b = [B.x - A.x, B.y - A.y]     # Storing vector A->B

  atb2 = a_to_b[0]**2 + a_to_b[1]**2  # **2 means "squared"
                                      #   Basically finding the squared magnitude
                                      #   of a_to_b

  atp_dot_atb = a_to_p[0]*a_to_b[0] + a_to_p[1]*a_to_b[1]
                                      # The dot product of a_to_p and a_to_b

  t = atp_dot_atb / atb2              # The normalized "distance" from a to
                                      #   your closest point

  return Point.new( :x => A.x + a_to_b[0]*t,
                    :y => A.y + a_to_b[1]*t )
                                      # Add the distance to A, moving
                                      #   towards B

end

Alternatively:

From Line-Line Intersection, at Wikipedia. First, find Q, which is a second point that is to be had from taking a step from P in the "right direction". This gives us four points.

def getClosestPointFromLine(A, B, P)

  a_to_b = [B.x - A.x, B.y - A.y]   # Finding the vector from A to B
                                        This step can be combined with the next
  perpendicular = [ -a_to_b[1], a_to_b[0] ]
                                    # The vector perpendicular to a_to_b;
                                        This step can also be combined with the next

  Q = Point.new(:x => P.x + perpendicular[0], :y => P.y + perpendicular[1])
                                    # Finding Q, the point "in the right direction"
                                    # If you want a mess, you can also combine this
                                    # with the next step.

  return Point.new (:x => ((A.x*B.y - A.y*B.x)*(P.x - Q.x) - (A.x-B.x)*(P.x*Q.y - P.y*Q.x)) / ((A.x - B.x)*(P.y-Q.y) - (A.y - B.y)*(P.y-Q.y)),
                    :y => ((A.x*B.y - A.y*B.x)*(P.y - Q.y) - (A.y-B.y)*(P.x*Q.y - P.y*Q.x)) / ((A.x - B.x)*(P.y-Q.y) - (A.y - B.y)*(P.y-Q.y)) )

end

Caching, Skipping steps, etc. is possible, for performance reasons.

Upvotes: 55

Darien Pardinas
Darien Pardinas

Reputation: 6186

I wrote this a long time ago, it's not much different to what others have said, but it's a copy/paste solution in C# if you have a class (or struct) named PointF with members X and Y:

private static PointF ClosestPointToSegment(PointF P, PointF A, PointF B)
{
    PointF a_to_p = new PointF(), a_to_b = new PointF();
    a_to_p.X = P.X - A.X;
    a_to_p.Y = P.Y - A.Y; //     # Storing vector A->P  
    a_to_b.X = B.X - A.X;
    a_to_b.Y = B.Y - A.Y; //     # Storing vector A->B

    float atb2 = a_to_b.X * a_to_b.X + a_to_b.Y * a_to_b.Y;
    float atp_dot_atb = a_to_p.X * a_to_b.X + a_to_p.Y * a_to_b.Y; // The dot product of a_to_p and a_to_b
    float t = atp_dot_atb / atb2;  //  # The normalized "distance" from a to the closest point
    return new PointF(A.X + a_to_b.X * t, A.Y + a_to_b.Y * t);
}

Update: Looking at the comments it looks like I adapted it to C# from the same source code mentioned in the accepted answer.

Upvotes: 5

Nick Bilyk
Nick Bilyk

Reputation: 371

In case somebody is looking for a way to do this with Java + LibGdx:

Intersector.nearestSegmentPoint

Upvotes: 0

MvG
MvG

Reputation: 60858

This answer is based on ideas from projective geometry.

Compute the cross product (Ax,Ay,1)×(Bx,By,1)=(u,v,w). The resulting vector describes the line connecting A and B: it has the equation ux+vy+w=0. But you can also interpret (u,v,0) as a point infinitely far away in a direction perpendicular to that line. Doing another cross product you get the line joining hat point to P: (u,v,0)×(Px,Py,1). And to intersect that line with the line AB, you do another cross product: ((u,v,0)×(Px,Py,1))×(u,v,w). The result will be a homogenous coordinate vector (x,y,z) from which you can read the coordinates of this closest point as (x/z,y/z).

Take everything together and you get the following formula:

{\scriptsize\begin{pmatrix}x\y\z\end{pmatrix}}=\Bigl(\bigl({\scriptsize\begin{pmatrix}1&0&0\0&1&0\0&0&0\end{pmatrix}}(A\times B)\bigr)\times P\Bigr)\times(A\times B)

Using a computer algebra system, you can find the resulting coordinates to be the following:

x = ((Ax - Bx)*Px + (Ay - By)*Py)*(Ax - Bx) + (Ay*Bx - Ax*By)*(Ay - By)
y = -(Ay*Bx - Ax*By)*(Ax - Bx) + ((Ax - Bx)*Px + (Ay - By)*Py)*(Ay - By)
z = (Ax - Bx)^2 + (Ay - By)^2

As you notice, there are a lot of recurring terms. Inventing (pretty much arbitrary) names for these, you can get the following final result, written in pseudocode:

dx = A.x - B.x
dy = A.y - B.y
det = A.y*B.x - A.x*B.y
dot = dx*P.x + dy*P.y
x = dot*dx + det*dy
y = dot*dy - det*dx
z = dx*dx + dy*dy
zinv = 1/z
return new Point(x*zinv, y*zinv)

Benefits of this approach:

  • No case distinctions
  • No square roots
  • Only a single division

Upvotes: 4

N.Schilke
N.Schilke

Reputation: 521

if anyone is interested in a C# XNA function based on the above:

    public static Vector2 GetClosestPointOnLineSegment(Vector2 A, Vector2 B, Vector2 P)
    {
        Vector2 AP = P - A;       //Vector from A to P   
        Vector2 AB = B - A;       //Vector from A to B  

        float magnitudeAB = AB.LengthSquared();     //Magnitude of AB vector (it's length squared)     
        float ABAPproduct = Vector2.Dot(AP, AB);    //The DOT product of a_to_p and a_to_b     
        float distance = ABAPproduct / magnitudeAB; //The normalized "distance" from a to your closest point  

        if (distance < 0)     //Check if P projection is over vectorAB     
        {
            return A;

        }
        else if (distance > 1)             {
            return B;
        }
        else
        {
            return A + AB * distance;
        }
    }

Upvotes: 52

stbn
stbn

Reputation: 425

The answer from Justin L. is almost fine, but it doesn't check if the normalized distance is less than 0, or higher than the AB vector magnitude. Then it won't work well when P vector proyection is out of bounds (from the line segment AB). Here's the corrected pseudocode:

    function GetClosestPoint(A, B, P)
{
  vectorAP = (p.x - a.x, p.y - a.y)     //Vector from A to P
  vectorAB = (b.x - a.x, b.y - a.y)     //Vector from A to B

  magnitudeAB = vectorAB[0]^2 + vectorAB[1]^2  
  //Magnitude of AB vector (it's length)


  ABAPproduct = vectorAB[0]*vectorAP[0] + vectorAB[1]*vectorAP[1] 
  //The product of a_to_p and a_to_b


  distance = ABAPproduct / magnitudeAB       
  //The normalized "distance" from a to your closest point

  if ( distance < 0)     //Check if P projection is over vectorAB
    {
        returnPoint.x = a.x
        returnPoint.y = a.y
    }   
  else if (distance > magnitudeAB)
    {
        returnPoint.x = b.x
        returnPoint.y = b.y
    }
  else
    {
        returnPoint.x = a.x + vectorAB[0]*distance
        returnPoint.y = a.y + vectorAB[1]*distance
    }

}

Upvotes: 7

Amadan
Amadan

Reputation: 198294

Find the slope a1 of AB by dividing the y-difference with the x-difference; then draw a perpendicular line (with slope a2 = -1/a1, you need to solve for the offset (b2) by putting P's coordinates into y = a2*x + b2); then you have two lines (i.e. two linear equations), and you need to solve the intersection. That will be your closest point.

Do the math right, and the function will be pretty trivial to write.

To elaborate a bit:

Original line:
y = a1 * x + b1
a1 = (By - Ay) / (Bx - Ax)   <--
b1 = Ay - a1 * Ax            <--

Perpendicular line:
y = a2 * x + b2
a2 = -1/a1                   <--
b2 = Py - a2 * Px            <--

Now you have P which lies on both lines:
y = a1 * x + b1
y = a2 * x + b2
--------------- subtract:
0 = (a1 - a2) * Px + (b1 - b2)
x = - (b1 - b2) / (a1 - a2)  <--
y = a1 * x + b1              <--

Hope I didn't mess up somewhere :) UPDATE Of course I did. Serve me right for not working things out on paper first. I deserved every downvote, but I'd've expected someone to correct me. Fixed (I hope).

Arrows point the way.

UPDATE Ah, the corner cases. Yeah, some languages don't handle infinities well. I did say the solution was language-free...

You can check the special cases, they're quite easy. The first one is when the x difference is 0. That means the line is vertical, and the closest point is on a horizontal perpendicular. Thus, x = Ax, y = Px.

The second one is when y difference is 0, and the opposite is true. Thus, x = Px, y = Ay

Upvotes: 4

comingstorm
comingstorm

Reputation: 26087

Your point (X) will be a linear combination of points A and B:

X = k A + (1-k) B

For X to be actually on the line segment, the parameter k must be between 0 and 1, inclusive. You can compute k as follows:

k_raw = (P-B).(A-B)  /  (A-B).(A-B)

(where the period denotes the dot product)

Then, to make sure the point is actually on the line segment:

if k_raw < 0:
    k= 0
elif k_raw > 1:
    k= 1
else:
    k= k_raw

Upvotes: 13

John Feminella
John Feminella

Reputation: 311436

The closest point C will be on a line whose slope is the reciprocal of AB and which intersects with P. This sounds like it might be homework, but I'll give some pretty strong hints, in order of increasing spoiler-alert level:

  • There can be only one such line.

  • This is a system of two line equations. Just solve for x and y.

  • Draw a line segment between A and B; call this L. The equation for L is y = mx + b, where m is the ratio of the y-coordinates to the x-coordinates. Solve for b using either A or B in the expression.

  • Do the same as above, but for CP. Now solve the simultaneous linear system of equations.

  • A Google search will give you a bevy of examples to choose from.

Upvotes: 1

Ryan Ternier
Ryan Ternier

Reputation: 8804

The algorithm would be quite easy:

you have 3 points - triangle. From there you should be able to find AB, AC, BC.

Theck this out: http://www.topcoder.com/tc?d1=tutorials&d2=geometry1&module=Static#line_point_distance

Upvotes: -3

Related Questions