Noskosteve wrote:

On Tuesday, January 27, 2015 at 3:01:56 PM UTC-6, GianLuigi Piacentini
wrote:

Hi all,

so far, I solved this kind of problem as a 2 by 2 system of linear
(parametric) equation, using a routine from a math library (SLATEC).
This
solver makes acrobacies to ensure numerical stability, minimize
round-off errors, and so on. That's fine, but seem to me an overkill for
a 2 by 2. Mantyla in his 'Introduction to Solid Modeling', solves
directly (I think it's Cramer's rule, but equations are so simple...)
without too much worry
about numerical aspects. He only checks for a non-zero denominator.
So far I do not have speed or memory constraint, I can go well with
SLATEC, but I would like to hear experts opinion on the subject.
If you will point me to alternative solvers, consider that it's an hobby
project and only public domain stuff is eligible (possibly Fortran,
alternatively Ada, Basic, c or Pascal - languages that I hope to be able
to understand) .

Thanks in advance.

Gigi

I know this is quite late and I do not know if this is what you are after, but.... From the C.G.A. FAQ:

Subject 1.03: How do I find intersections of 2 2D line segments?
This problem can be extremely easy or extremely difficult depends on
your applications. If all you want is the intersection point, the
following should work:

You can consider the word "vector" to mean "point" or "location".

Let A,B,C,D be 2-space position vectors. This means that A is located
at Ax, Ay ; B is located at Bx, By; etc. Then the directed line
segments AB & CD are given by:
AB=A+r(B-A), r in [0,1]
CD=C+s(D-C), s in [0,1]

This means that A & B are two points on line AB and C & D are two points
on line CD.

If AB & CD intersect, then
A+r(B-A)=C+s(D-C), or

Ax+r(Bx-Ax)=Cx+s(Dx-Cx)
Ay+r(By-Ay)=Cy+s(Dy-Cy) for some r,s in [0,1]

Solving the above for r and s yields
(Ay-Cy)(Dx-Cx)-(Ax-Cx)(Dy-Cy)
r = ----------------------------- (eqn 1)
(Bx-Ax)(Dy-Cy)-(By-Ay)(Dx-Cx)

(Ay-Cy)(Bx-Ax)-(Ax-Cx)(By-Ay)
s = ----------------------------- (eqn 2)
(Bx-Ax)(Dy-Cy)-(By-Ay)(Dx-Cx)

Let P be the position/location (vector) of the intersection point,
then
P=A+r(B-A) or

Px=Ax+r(Bx-Ax)
Py=Ay+r(By-Ay)

Note that "s" is not needed for position of intersection.

By examining the values of r & s, you can also determine some
other limiting conditions:

If 0<=r<=1 & 0<=s<=1, intersection exists
r<0 or r>1 or s<0 or s>1 line segments do not intersect

If the denominator in eqn 1 is zero, AB & CD are parallel
If the numerator in eqn 1 is also zero, AB & CD are collinear.

If they are collinear, then the segments may be projected to the x-
or y-axis, and overlap of the projected intervals checked.

If the intersection point of the 2 lines are needed (lines in this
context mean infinite lines) regardless whether the two line
segments intersect, then

If r>1, P is located on extension of AB
If r<0, P is located on extension of BA
If s>1, P is located on extension of CD
If s<0, P is located on extension of DC

Also note that the denominators of eqn 1 & 2 are identical and thus,
only need to be calculated once.

References:

[O'Rourke (C)] pp. 249-51
[Gems III] pp. 199-202 "Faster Line Segment Intersection,"

Thanks.

Seems to me that we have the usual problems: possible cancellation in the floating point subtraction (I think we cannot overcome this) and the
equalities to zero.
We are living in a floating point world, so probably the equalities are to
be transformed into ABS(value) < epsilon, where epsilon is a small number.
And now the problem is in properly choosing epsilon.

Any guidance ?

Gigi

--- SoupGate-Win32 v1.05
* Origin: fsxNet Usenet Gateway (21:1/5)

On Tuesday, January 27, 2015 at 3:01:56 PM UTC-6, GianLuigi Piacentini wrote:

Hi all,

so far, I solved this kind of problem as a 2 by 2 system of linear (parametric) equation, using a routine from a math library (SLATEC). This solver makes acrobacies to ensure numerical stability, minimize round-off errors, and so on. That's fine, but seem to me an overkill for a 2 by 2. Mantyla in his 'Introduction to Solid Modeling', solves directly (I think it's Cramer's rule, but equations are so simple...) without too much worry about numerical aspects. He only checks for a non-zero denominator.
So far I do not have speed or memory constraint, I can go well with SLATEC, but I would like to hear experts opinion on the subject.
If you will point me to alternative solvers, consider that it's an hobby project and only public domain stuff is eligible (possibly Fortran, alternatively Ada, Basic, c or Pascal - languages that I hope to be able to understand) .

Thanks in advance.

Gigi

I know this is quite late and I do not know if this is what you are after, but....
From the C.G.A. FAQ:

Subject 1.03: How do I find intersections of 2 2D line segments?
This problem can be extremely easy or extremely difficult depends on your applications. If all you want is the intersection point, the following should work:

You can consider the word "vector" to mean "point" or "location".

Let A,B,C,D be 2-space position vectors. This means that A is located at Ax, Ay ; B is located at Bx, By; etc. Then the directed line segments AB & CD are given by:
AB=A+r(B-A), r in [0,1]
CD=C+s(D-C), s in [0,1]

This means that A & B are two points on line AB and C & D are two points on line CD.

If AB & CD intersect, then
A+r(B-A)=C+s(D-C), or

Ax+r(Bx-Ax)=Cx+s(Dx-Cx)
Ay+r(By-Ay)=Cy+s(Dy-Cy) for some r,s in [0,1]

Solving the above for r and s yields
(Ay-Cy)(Dx-Cx)-(Ax-Cx)(Dy-Cy)
r = ----------------------------- (eqn 1)
(Bx-Ax)(Dy-Cy)-(By-Ay)(Dx-Cx)

(Ay-Cy)(Bx-Ax)-(Ax-Cx)(By-Ay)
s = ----------------------------- (eqn 2)
(Bx-Ax)(Dy-Cy)-(By-Ay)(Dx-Cx)

Let P be the position/location (vector) of the intersection point, then
P=A+r(B-A) or

Px=Ax+r(Bx-Ax)
Py=Ay+r(By-Ay)

Note that "s" is not needed for position of intersection.

By examining the values of r & s, you can also determine some
other limiting conditions:

If 0<=r<=1 & 0<=s<=1, intersection exists
r<0 or r>1 or s<0 or s>1 line segments do not intersect

If the denominator in eqn 1 is zero, AB & CD are parallel
If the numerator in eqn 1 is also zero, AB & CD are collinear.

If they are collinear, then the segments may be projected to the x-
or y-axis, and overlap of the projected intervals checked.

If the intersection point of the 2 lines are needed (lines in this
context mean infinite lines) regardless whether the two line
segments intersect, then

If r>1, P is located on extension of AB
If r<0, P is located on extension of BA
If s>1, P is located on extension of CD
If s<0, P is located on extension of DC

Also note that the denominators of eqn 1 & 2 are identical and thus, only need to be calculated once.

References:

[O'Rourke (C)] pp. 249-51
[Gems III] pp. 199-202 "Faster Line Segment Intersection,"

--- SoupGate-Win32 v1.05
* Origin: fsxNet Usenet Gateway (21:1/5)

On Monday, November 2, 2015 at 12:35:54 PM UTC-6, GianLuigi Piacentini wrote:

Noskosteve wrote:

On Tuesday, January 27, 2015 at 3:01:56 PM UTC-6, GianLuigi Piacentini wrote:

Hi all,

so far, I solved this kind of problem as a 2 by 2 system of linear
(parametric) equation, using a routine from a math library (SLATEC).
This
solver makes acrobacies to ensure numerical stability, minimize
round-off errors, and so on. That's fine, but seem to me an overkill for >> a 2 by 2. Mantyla in his 'Introduction to Solid Modeling', solves
directly (I think it's Cramer's rule, but equations are so simple...)
without too much worry
about numerical aspects. He only checks for a non-zero denominator.
So far I do not have speed or memory constraint, I can go well with
SLATEC, but I would like to hear experts opinion on the subject.
If you will point me to alternative solvers, consider that it's an hobby >> project and only public domain stuff is eligible (possibly Fortran,
alternatively Ada, Basic, c or Pascal - languages that I hope to be able >> to understand) .

Thanks in advance.

Gigi

I know this is quite late and I do not know if this is what you are after, but.... From the C.G.A. FAQ:

Subject 1.03: How do I find intersections of 2 2D line segments?
This problem can be extremely easy or extremely difficult depends on
your applications. If all you want is the intersection point, the
following should work:

You can consider the word "vector" to mean "point" or "location".

Let A,B,C,D be 2-space position vectors. This means that A is located
at Ax, Ay ; B is located at Bx, By; etc. Then the directed line
segments AB & CD are given by:
AB=A+r(B-A), r in [0,1]
CD=C+s(D-C), s in [0,1]

This means that A & B are two points on line AB and C & D are two points
on line CD.

If AB & CD intersect, then
A+r(B-A)=C+s(D-C), or

Ax+r(Bx-Ax)=Cx+s(Dx-Cx)
Ay+r(By-Ay)=Cy+s(Dy-Cy) for some r,s in [0,1]

Solving the above for r and s yields
(Ay-Cy)(Dx-Cx)-(Ax-Cx)(Dy-Cy)
r = ----------------------------- (eqn 1)
(Bx-Ax)(Dy-Cy)-(By-Ay)(Dx-Cx)

(Ay-Cy)(Bx-Ax)-(Ax-Cx)(By-Ay)
s = ----------------------------- (eqn 2)
(Bx-Ax)(Dy-Cy)-(By-Ay)(Dx-Cx)

Let P be the position/location (vector) of the intersection point,
then
P=A+r(B-A) or

Px=Ax+r(Bx-Ax)
Py=Ay+r(By-Ay)

Note that "s" is not needed for position of intersection.

By examining the values of r & s, you can also determine some
other limiting conditions:

If 0<=r<=1 & 0<=s<=1, intersection exists
r<0 or r>1 or s<0 or s>1 line segments do not intersect

If the denominator in eqn 1 is zero, AB & CD are parallel
If the numerator in eqn 1 is also zero, AB & CD are collinear.

If they are collinear, then the segments may be projected to the x-
or y-axis, and overlap of the projected intervals checked.

If the intersection point of the 2 lines are needed (lines in this
context mean infinite lines) regardless whether the two line
segments intersect, then

If r>1, P is located on extension of AB
If r<0, P is located on extension of BA
If s>1, P is located on extension of CD
If s<0, P is located on extension of DC

Also note that the denominators of eqn 1 & 2 are identical and thus,
only need to be calculated once.

References:

[O'Rourke (C)] pp. 249-51
[Gems III] pp. 199-202 "Faster Line Segment Intersection,"

Thanks.

Seems to me that we have the usual problems: possible cancellation in the floating point subtraction (I think we cannot overcome this) and the equalities to zero.
We are living in a floating point world, so probably the equalities are to
be transformed into ABS(value) < epsilon, where epsilon is a small number. And now the problem is in properly choosing epsilon.

Any guidance ?

Gigi

I can not add anything. I used this method on a floating point platform several years ago and had no problems.
Ciao, Steve

--- SoupGate-Win32 v1.05
* Origin: fsxNet Usenet Gateway (21:1/5)