[GRASS5] [mshapiro@ncsa.uiuc.edu: bug in libes/gis/line_dist.c]

[Radim Blazek]

On Thursday 21 February 2002 13:45, Markus Neteler wrote:
> Hi developers,
>
> find attached a bug report from Michael Shapiro.
> Can anyone confirm this bug?
>
> Markus
>
> Hi Markus,
>
> I was looking through the grass src code and found a bug in the
> G_distance2_point_to_line() function in libes/gis/line_dist.c
>
> lines 47,48 should be changed from
>
> 	dx = dx * t + x1;
> 	dy = dy * t + y1;
> to
> 	dx = dx * t;
> 	dy = dy * t;
>
> This bug has been there ever since I wrote this piece of code,
> but only now just discovered it. You should have someone double
> check me on this, though, just to be sure.
>
> Michael Shapiro

I think that (dx * t + x1, dy * t + y1) is point x,y projected to line 
x1,y1->x2,y2 (comment in code:  /* go t from x1,y1 towards x2,y2 */)
and we want distance between this point and and x,y. That is:
dx = dx * t + x1 - x;
dy = dy * t + y1 - y;

Example: line: 1,0 -> 3,0, point 2,3 (distance^2 = 9)
dx = x2 - x1 = 2
dy = y2 - y1 = 0
dx' = x - x1 = 1
dy' = y - y1 = 3 
t = (2*1 + 0*3) / (2*2 + 0*0) = 0.5
original:
dx = 2*0.5 + 1 = 2
dy = 0*0.5 + 0 = 0
dx*dx + dy*dy = 2*2 + 0*0 = 4
Shapiro's fix:
dx = 2*0.5 = 1
dy = 0*0.5 = 0
dx*dx + dy*dy = 1*1 + 0*0 = 1
my proposal:
dx = 2*0.5 + 1 - 2 = 0 
dy = 0*0.5 + 0 - 3 = -3
dx*dx + dy*dy = (-3)*(-3) + 0*0 = 9

BTW, version from dig_distance2_point_to_line():
  if (t < 0.0)                  /* go to x1,y1 */
    t = 0.0;
  else if (t > 1.0)             /* go to x2,y2 */
    t = 1.0;

  dx = dx * t + x1 - x;
  dy = dy * t + y1 - y;

You should have someone double check me on this, though, just to be sure :)

Radim

-- 
Radim Blazek
Praga, Bohemia