Perpendicular on a line segment from a given point

Question

I want to calculate a point on a given line that is perpendicular from a given point.

I have a line segment AB and have a point C outside line segment. I want to calculate a point D on AB such that CD is perpendicular to AB.

Find point D

I have to find point D.

It quite similar to this, but I want to consider to Z coordinate also as it does not show up correctly in 3D space.

In the absence of a *programming* question, this would be better on [math.se] (where it is almost certainly already duplicated) — AakashM, Apr 24 '12 at 15:36
http://stackoverflow.com/questions/1811549/perpendicular-on-a-line-from-a-given-point — jdbertron, Mar 03 '13 at 15:43

score 37 · Accepted Answer · answered Mar 03 '13 at 15:56

37

Proof: Point D is on a line CD perpendicular to AB, and of course D belongs to AB. Write down the Dot product of the two vectors CD.AB = 0, and express the fact D belongs to AB as D=A+t(B-A).

We end up with 3 equations:

 Dx=Ax+t(Bx-Ax)
 Dy=Ay+t(By-Ay)
(Dx-Cx)(Bx-Ax)+(Dy-Cy)(By-Ay)=0

Subtitute the first two equations in the third one gives:

(Ax+t(Bx-Ax)-Cx)(Bx-Ax)+(Ay+t(By-Ay)-Cy)(By-Ay)=0

Distributing to solve for t gives:

(Ax-Cx)(Bx-Ax)+t(Bx-Ax)(Bx-Ax)+(Ay-Cy)(By-Ay)+t(By-Ay)(By-Ay)=0

which gives:

t= -[(Ax-Cx)(Bx-Ax)+(Ay-Cy)(By-Ay)]/[(Bx-Ax)^2+(By-Ay)^2]

getting rid of the negative signs:

t=[(Cx-Ax)(Bx-Ax)+(Cy-Ay)(By-Ay)]/[(Bx-Ax)^2+(By-Ay)^2]

Once you have t, you can figure out the coordinates for D from the first two equations.

 Dx=Ax+t(Bx-Ax)
 Dy=Ay+t(By-Ay)

answered Mar 03 '13 at 15:56

jdbertron

733
7
11

2

This ignores the Z component of the original question. – Stephen Aug 10 '15 at 05:08
2

If you want to limit the found point D to the line between A and B limit t for 0 >= t <= 1. – Oliver Becker Oct 05 '19 at 16:12
1

If this not work for u.. Try this solution https://stackoverflow.com/questions/1811549/perpendicular-on-a-line-from-a-given-point – Zumbarlal Saindane Apr 03 '20 at 17:52

cuixiping · Answer 2 · 2015-06-20T02:11:46.927

27

function getSpPoint(A,B,C){
    var x1=A.x, y1=A.y, x2=B.x, y2=B.y, x3=C.x, y3=C.y;
    var px = x2-x1, py = y2-y1, dAB = px*px + py*py;
    var u = ((x3 - x1) * px + (y3 - y1) * py) / dAB;
    var x = x1 + u * px, y = y1 + u * py;
    return {x:x, y:y}; //this is D
}

question

edited Jun 20 '15 at 02:11

answered Sep 19 '12 at 17:01

cuixiping

24,167
8
82
93

8

Would be nice to see a little explanation of what you did. – Patrick Murphy Jul 21 '15 at 15:35
Hi. I have point a,b and d and i want to calculate x and y coordinates of point c. Can i use your function? If no please tell me the modifications. – Mean Coder Jun 11 '17 at 17:09
@MeanCoder point c has countless possibilities if only know point a,b,d. – cuixiping Jul 20 '17 at 02:59
What if A,B, C are GPS coordinates (lat,lon)? – TkrA Apr 12 '22 at 12:44
It doesn't matter if they're GPS or not. They're still coordinates. You will have to account for cases where A and B lie on different sides of the 0 coordinate. – Vlad Vyatkin Aug 12 '22 at 10:57

score 8 · Answer 3 · answered Apr 24 '12 at 16:01

There is a simple closed form solution for this (requiring no loops or approximations) using the vector dot product.

Imagine your points as vectors where point A is at the origin (0,0) and all other points are referenced from it (you can easily transform your points to this reference frame by subtracting point A from every point).

In this reference frame point D is simply the vector projection of point C on the vector B which is expressed as:

// Per wikipedia this is more efficient than the standard (A . Bhat) * Bhat
Vector projection = Vector.DotProduct(A, B) / Vector.DotProduct(B, B) * B

The result vector can be transformed back to the original coordinate system by adding point A to it.

score 6 · Answer 4 · answered Apr 24 '12 at 18:01

A point on line AB can be parametrized by:

M(x)=A+x*(B-A), for x real.

You want D=M(x) such that DC and AB are orthogonal:

dot(B-A,C-M(x))=0.

That is: dot(B-A,C-A-x*(B-A))=0, or dot(B-A,C-A)=x*dot(B-A,B-A), giving:

x=dot(B-A,C-A)/dot(B-A,B-A) which is defined unless A=B.

score 2 · Answer 5 · answered Apr 24 '12 at 15:57

2

What you are trying to do is called vector projection

answered Apr 24 '12 at 15:57

Nicolas Repiquet

9,097
2
31
53

score 1 · Answer 6 · answered Mar 09 '16 at 11:50

Here i have converted answered code from "cuixiping" to matlab code.

function Pr=getSpPoint(Line,Point)
% getSpPoint(): find Perpendicular on a line segment from a given point
x1=Line(1,1);
y1=Line(1,2);
x2=Line(2,1);
y2=Line(2,1);
x3=Point(1,1);
y3=Point(1,2);

px = x2-x1;
py = y2-y1;
dAB = px*px + py*py;

u = ((x3 - x1) * px + (y3 - y1) * py) / dAB;
x = x1 + u * px;
y = y1 + u * py;

Pr=[x,y];

end

score 0 · Answer 7 · answered Dec 21 '16 at 05:12

I didn't see this answer offered, but Ron Warholic had a great suggestion with the Vector Projection. ACD is merely a right triangle.

Create the vector AC i.e (Cx - Ax, Cy - Ay)
Create the Vector AB i.e (Bx - Ax, By - Ay)
Dot product of AC and AB is equal to the cosine of the angle between the vectors. i.e cos(theta) = ACx*ABx + ACy*ABy.
Length of a vector is sqrt(x*x + y*y)
Length of AD = cos(theta)*length(AC)
Normalize AB i.e (ABx/length(AB), ABy/length(AB))
D = A + NAB*length(AD)

score 0 · Answer 8 · edited Mar 19 '21 at 07:43

0

For anyone who might need this in C# I'll save you some time:

double Ax = ;
double Ay = ;
double Az = ;
    
double Bx = ;
double By = ;
double Bz = ;
    
double Cx = ;
double Cy = ;
double Cz = ; 
    
double t = ((Cx - Ax) * (Bx - Ax) + (Cy - Ay) * (By - Ay)) / (Math.Pow(Bx - Ax, 2) + Math.Pow(By - Ay, 2));
    
double Dx = Ax + t*(Bx - Ax);
double Dy = Ay + t*(By - Ay);

edited Mar 19 '21 at 07:43

Dharman

30,962
25
85
135

answered Mar 19 '21 at 07:38

Michael684

9
3

MBV · Answer 9 · 2021-12-30T03:57:59.593

Here is another python implementation without using a for loop. It works for any number of points and any number of line segments. Given p_array as a set of points, and x_array , y_array as continues line segments or a polyline.

This uses the equation Y = mX + n and considering that the m factor for a perpendicular line segment is -1/m.


      import numpy as np

      def ortoSegmentPoint(self, p_array, x_array, y_array):
            """
    
            :param p_array: np.array([[ 718898.941  9677612.901 ], [ 718888.8227 9677718.305 ], [ 719033.0528 9677770.692 ]])
            :param y_array: np.array([9677656.39934991 9677720.27550726 9677754.79])
            :param x_array: np.array([718895.88881594 718938.61392781 718961.46])
            :return: [POINT, LINE] indexes where point is orthogonal to line segment
            """
            # PENDIENTE "m" de la recta, y = mx + n
            m_array = np.divide(y_array[1:] - y_array[:-1], x_array[1:] - x_array[:-1])
            # PENDIENTE INVERTIDA, 1/m
            inv_m_array = np.divide(1, m_array)
            # VALOR "n", y = mx + n
            n_array = y_array[:-1] - x_array[:-1] * m_array
            # VALOR "n_orto" PARA LA RECTA PERPENDICULAR
            n_orto_array = np.array(p_array[:, 1]).reshape(len(p_array), 1) + inv_m_array * np.array(p_array[:, 0]).reshape(len(p_array), 1)
            # PUNTOS DONDE SE INTERSECTAN DE FORMA PERPENDICULAR
            x_intersec_array = np.divide(n_orto_array - n_array, m_array + inv_m_array)
            y_intersec_array = m_array * x_intersec_array + n_array
            # LISTAR COORDENADAS EN PARES
            x_coord = np.array([x_array[:-1], x_array[1:]]).T
            y_coord = np.array([y_array[:-1], y_array[1:]]).T
            # FILAS: NUMERO DE PUNTOS, COLUMNAS: NUMERO DE TRAMOS
            maskX = np.where(np.logical_and(x_intersec_array < np.max(x_coord, axis=1), x_intersec_array > np.min(x_coord, axis=1)), True, False)
            maskY = np.where(np.logical_and(y_intersec_array < np.max(y_coord, axis=1), y_intersec_array > np.min(y_coord, axis=1)), True, False)
            mask = maskY * maskX
            return np.argwhere(mask == True)

score 0 · Answer 10 · answered Aug 03 '22 at 16:41

As Ron Warholic and Nicolas Repiquet answered, this can be solved using vector projection. For completeness I'll add a python/numpy implementation of this here in case it saves anyone else some time:

import numpy as np

# Define some test data that you can solve for directly.
first_point = np.array([4, 4])
second_point = np.array([8, 4])
target_point = np.array([6, 6])

# Expected answer
expected_point = np.array([6, 4])

# Create vector for first point on line to perpendicular point.
point_vector = target_point - first_point
# Create vector for first point and second point on line.
line_vector = second_point - first_point

# Create the projection vector that will define the position of the resultant point with respect to the first point.
projection_vector = (np.dot(point_vector, line_vector) / np.dot(line_vector, line_vector)) * line_vector

# Alternative method proposed in another answer if for whatever reason you prefer to use this.
_projection_vector = (np.dot(point_vector, line_vector) / np.linalg.norm(line_vector)**2) * line_vector

# Add the projection vector to the first point
projected_point = first_point + projection_vector

# Test
(projected_point == expected_point).all()

garlimidon · Answer 11 · 2012-04-24T15:38:18.903

Since you're not stating which language you're using, I'll give you a generic answer:

Just have a loop passing through all the points in your AB segment, "draw a segment" to C from them, get the distance from C to D and from A to D, and apply pithagoras theorem. If AD^2 + CD^2 = AC^2, then you've found your point.

Also, you can optimize your code by starting the loop by the shortest side (considering AD and BD sides), since you'll find that point earlier.

score -1 · Answer 12 · edited May 23 '17 at 12:17

Here is a python implementation based on Corey Ogburn's answer from this thread.
It projects the point q onto the line segment defined by p1 and p2 resulting in the point r.
It will return null if r falls outside of line segment:

def is_point_on_line(p1, p2, q):

    if (p1[0] == p2[0]) and (p1[1] == p2[1]):
        p1[0] -= 0.00001

    U = ((q[0] - p1[0]) * (p2[0] - p1[0])) + ((q[1] - p1[1]) * (p2[1] - p1[1]))
    Udenom = math.pow(p2[0] - p1[0], 2) + math.pow(p2[1] - p1[1], 2)
    U /= Udenom

    r = [0, 0]
    r[0] = p1[0] + (U * (p2[0] - p1[0]))
    r[1] = p1[1] + (U * (p2[1] - p1[1]))

    minx = min(p1[0], p2[0])
    maxx = max(p1[0], p2[0])
    miny = min(p1[1], p2[1])
    maxy = max(p1[1], p2[1])

    is_valid = (minx <= r[0] <= maxx) and (miny <= r[1] <= maxy)

    if is_valid:
        return r
    else:
        return None

Perpendicular on a line segment from a given point

12 Answers12

Linked