Re: Line intersection point is wrong

Emre Hasegeli  writes:
>> You haven't done anything to exclude the possibility that l1->B is zero,
>> so you could be getting zero-divide in the y calculation.

> Ah, yes.  We can use l2 when l1->B is zero.

After working out the algebra by hand, I think the attached is correct
(and it does pass the regression tests, yay).  I also made the
line_parallel and line_perp tests more symmetric and zero-divide-free.

			regards, tom lane

diff --git a/src/backend/utils/adt/geo_ops.c b/src/backend/utils/adt/geo_ops.c
index 657bcee..ce796c7 100644
*** a/src/backend/utils/adt/geo_ops.c
--- b/src/backend/utils/adt/geo_ops.c
*************** line_parallel(PG_FUNCTION_ARGS)
*** 1108,1117 ****
  	LINE	   *l1 = PG_GETARG_LINE_P(0);
  	LINE	   *l2 = PG_GETARG_LINE_P(1);
  
! 	if (FPzero(l1->B))
  		PG_RETURN_BOOL(FPzero(l2->B));
  
! 	PG_RETURN_BOOL(FPeq(l2->A, l1->A * (l2->B / l1->B)));
  }
  
  Datum
--- 1108,1119 ----
  	LINE	   *l1 = PG_GETARG_LINE_P(0);
  	LINE	   *l2 = PG_GETARG_LINE_P(1);
  
! 	if (FPzero(l1->A))			/* horizontal? */
! 		PG_RETURN_BOOL(FPzero(l2->A));
! 	if (FPzero(l1->B))			/* vertical? */
  		PG_RETURN_BOOL(FPzero(l2->B));
  
! 	PG_RETURN_BOOL(FPeq(l2->A * l1->B, l1->A * l2->B));
  }
  
  Datum
*************** line_perp(PG_FUNCTION_ARGS)
*** 1120,1131 ****
  	LINE	   *l1 = PG_GETARG_LINE_P(0);
  	LINE	   *l2 = PG_GETARG_LINE_P(1);
  
! 	if (FPzero(l1->A))
  		PG_RETURN_BOOL(FPzero(l2->B));
! 	else if (FPzero(l1->B))
  		PG_RETURN_BOOL(FPzero(l2->A));
  
! 	PG_RETURN_BOOL(FPeq(((l1->A * l2->B) / (l1->B * l2->A)), -1.0));
  }
  
  Datum
--- 1122,1133 ----
  	LINE	   *l1 = PG_GETARG_LINE_P(0);
  	LINE	   *l2 = PG_GETARG_LINE_P(1);
  
! 	if (FPzero(l1->A))			/* horizontal? */
  		PG_RETURN_BOOL(FPzero(l2->B));
! 	if (FPzero(l1->B))			/* vertical? */
  		PG_RETURN_BOOL(FPzero(l2->A));
  
! 	PG_RETURN_BOOL(FPeq(l2->A * l1->B, -(l1->A * l2->B)));
  }
  
  Datum
*************** line_interpt_internal(LINE *l1, LINE *l2
*** 1233,1250 ****
  
  	if (FPzero(l1->B))			/* l1 vertical? */
  	{
! 		x = l1->C;
! 		y = (l2->A * x + l2->C);
  	}
  	else if (FPzero(l2->B))		/* l2 vertical? */
  	{
! 		x = l2->C;
! 		y = (l1->A * x + l1->C);
  	}
  	else
  	{
! 		x = (l1->C - l2->C) / (l2->A - l1->A);
! 		y = (l1->A * x + l1->C);
  	}
  	result = point_construct(x, y);
  
--- 1235,1252 ----
  
  	if (FPzero(l1->B))			/* l1 vertical? */
  	{
! 		x = -(l1->C / l1->A);
! 		y = -(l2->A * x + l2->C) / l2->B;
  	}
  	else if (FPzero(l2->B))		/* l2 vertical? */
  	{
! 		x = -(l2->C / l2->A);
! 		y = -(l1->A * x + l1->C) / l1->B;
  	}
  	else
  	{
! 		x = (l1->B * l2->C - l2->B * l1->C) / (l1->A * l2->B - l2->A * l1->B);
! 		y = -(l1->A * x + l1->C) / l1->B;
  	}
  	result = point_construct(x, y);