Geometry - shortest distance between a point and line

Mark Oates

I'm having trouble grasping this simple concept for some reason. Basically I have a point (A) and a line (BC), and I want to find the shortest distance between the two.

591076

Knowing that the shortest distance will be a line (AX) that has a slope perpendicular to the slope of BC and goes through A. So the intersection of the two lines is X.

591077

For some reason I can't seem to get the math to work out right?! Basically I find the slope-intercept form of both lines and then work the two out to find the x and y coordinates of point X, but this code isn't giving me the correct answers. Here's what I have:

1#include "allegro.h"
2#include <math.h>
3 
4float d ( float x1, float y1, float x2, float y2)
5{
  return sqrt( ((x1 - x2) * (x1 - x2)) + ((y1 - y2) * (y1 - y2)) ) ;
7}
8 
9float shortest_distance(float point_x, float point_y, float line_x1, float line_y1, float line_x2, float line_y2)
10{
  // find the slope
  float slope_of_line = (line_y1 - line_y2) / (line_x1 - line_x2);
  
  // find the perpendicular slope
  float perpendicular_slope = (line_x1 - line_x2) / (line_y1 - line_y2) * -1;
  
  // find the y_intercept of line BC
  float y_intercept = slope_of_line * line_x2 - line_y2;
  
  // find the y_intercept of line AX
  float new_line_y_intercept = perpendicular_slope * -1 * point_x - point_y;
  
  // get the x_coordinate of point X
  // equation of BC is    y = slope_of_line * x + y_intercept;
  // equation of AX is    y = perpendicular_slope * x + new_line_y_intercept;
  //   perpendicular_slope * x + new_line_y_intercept == slope_of_line * x + y_intercept;
  //   perpendicular_slope * x == slope_of_line * x + y_intercept - new_line_y_intercept;
  //   (perpendicular_slope - slope_of_line) * x == (y_intercept - new_line_y_intercept);
  float intersect_x = (y_intercept - new_line_y_intercept) / (perpendicular_slope - slope_of_line);
  
  // get the y_coordinate of point X
  float intersect_y = slope_of_line * intersect_x + y_intercept;
  
  // measure the distance between A and X
  return d(point_x, point_y, intersect_x, intersect_y);
36}
37 
38 
39 
40void main()
41{
   allegro_init();
   
   float point_x = 0;
   float point_y = -10;
   float line_x1 = 0;
   float line_y1 = 0;
   float line_x2 = 10;
   float line_y2 = -10;
   
   float dist = shortest_distance(point_x, point_y, line_x1, line_y1, line_x2, line_y2);
   dist = d(line_x1, line_y1, line_x2, line_y2);
   
   allegro_message("%i", (int)dist);
55} END_OF_MAIN();

What am I doing wrong?

SiegeLord

Instead of trying to fix your algorithm I shall suggest you a better one. It is faster, and it can handle vertical/horizontal lines(which yours can't as of now).

What you do is you find a vector from the first point on the line to the second point. Then you find the absolute value of the cross product between the two vectors(which find the area of the parallelogram formed by the two vectors) and then divide it by the magnitude of the first vector to find the height(which is what you want).

Here's the code for that:

int A = x - x1;
int B = y - y1;
int C = x2 - x1;
int D = y2 - y1;

int dist = abs(A * D - C * B) / sqrt(C * C + D * D);

Mark Oates

seems to work, but I'm guessing this assumes the line is of infinate length? How should I accomidate for this?

Matthew Leverton

Then it would be the end point.

Mark Oates

Hmm, I still don't understand. How do I know when it's going to be the end point or a point along the line?

SiegeLord

(x1,y1) and (x2,y2) are any two(non-coinciding) on the line. There is no order requirement.

This method does indeed assume an infinite line. It shall return the distance from the line to the point regardless whether the perpendicular projection of it onto the line(point X in your diagram) is actually between the two points.

If you want the distance of a point from a segment(which is what I assume you are referring to) you need a different approach. It is obviously slower.

First you get the above vectors the same way. Then you calculate the dot product between them. You divide that by the length of the line segment to get the length of the projection. Then you divide it again by the same value to get the parameter value of that point. If that point parameter is less than 0, then it is behind the start of the segment, if it is more than 1 then it is in front of the end of the segment. If it is between those two then it is on the segment. You use the parameter to calculate the point(or choose one of the endpoints) and then use the distance formula to find the distance from the two points.

1float A = x - x1;
2float B = y - y1;
3float C = x2 - x1;
4float D = y2 - y1;
5 
6float dot = A * C + B * D;
7float len_sq = C * C + D * D;
8float param = dot / len_sq;
9 
10float xx,yy;
11 
12if(param < 0)
13{
14    xx = x1;
15    yy = y1;
16}
17else if(param > 1)
18{
19    xx = x2;
20    yy = y2;
21}
22else
23{
24    xx = x1 + param * C;
25    yy = y1 + param * D;
26}
27 
28float dist = d(x,y,xx,yy);//your distance function

Mark Oates

Thank you SO much! It's working perfectly now.

Fladimir da Gorf

SiegeLord, could I by any chance add this function to OpenLayer (Line::GetShortestDistanceTo( const Vec2D &point ))? I think it'd be a nice addition

SiegeLord

Sure.

It should also work for 3D lines/points as well.

Fladimir da Gorf

OK, there it is now, thanks!

Thread #589720. Printed from Allegro.cc

1	#include "allegro.h"
2	#include <math.h>
3
4	float d ( float x1, float y1, float x2, float y2)
5	{
6	return sqrt( ((x1 - x2) * (x1 - x2)) + ((y1 - y2) * (y1 - y2)) ) ;
7	}
8
9	float shortest_distance(float point_x, float point_y, float line_x1, float line_y1, float line_x2, float line_y2)
10	{
11	// find the slope
12	float slope_of_line = (line_y1 - line_y2) / (line_x1 - line_x2);
13
14	// find the perpendicular slope
15	float perpendicular_slope = (line_x1 - line_x2) / (line_y1 - line_y2) * -1;
16
17	// find the y_intercept of line BC
18	float y_intercept = slope_of_line * line_x2 - line_y2;
19
20	// find the y_intercept of line AX
21	float new_line_y_intercept = perpendicular_slope * -1 * point_x - point_y;
22
23	// get the x_coordinate of point X
24	// equation of BC is y = slope_of_line * x + y_intercept;
25	// equation of AX is y = perpendicular_slope * x + new_line_y_intercept;
26	// perpendicular_slope * x + new_line_y_intercept == slope_of_line * x + y_intercept;
27	// perpendicular_slope * x == slope_of_line * x + y_intercept - new_line_y_intercept;
28	// (perpendicular_slope - slope_of_line) * x == (y_intercept - new_line_y_intercept);
29	float intersect_x = (y_intercept - new_line_y_intercept) / (perpendicular_slope - slope_of_line);
30
31	// get the y_coordinate of point X
32	float intersect_y = slope_of_line * intersect_x + y_intercept;
33
34	// measure the distance between A and X
35	return d(point_x, point_y, intersect_x, intersect_y);
36	}
37
38
39
40	void main()
41	{
42	allegro_init();
43
44	float point_x = 0;
45	float point_y = -10;
46	float line_x1 = 0;
47	float line_y1 = 0;
48	float line_x2 = 10;
49	float line_y2 = -10;
50
51	float dist = shortest_distance(point_x, point_y, line_x1, line_y1, line_x2, line_y2);
52	dist = d(line_x1, line_y1, line_x2, line_y2);
53
54	allegro_message("%i", (int)dist);
55	} END_OF_MAIN();

1	float A = x - x1;
2	float B = y - y1;
3	float C = x2 - x1;
4	float D = y2 - y1;
5
6	float dot = A * C + B * D;
7	float len_sq = C * C + D * D;
8	float param = dot / len_sq;
9
10	float xx,yy;
11
12	if(param < 0)
13	{
14	xx = x1;
15	yy = y1;
16	}
17	else if(param > 1)
18	{
19	xx = x2;
20	yy = y2;
21	}
22	else
23	{
24	xx = x1 + param * C;
25	yy = y1 + param * D;
26	}
27
28	float dist = d(x,y,xx,yy);//your distance function