Point on Line Segment in 2D. Which code is faster ? Can you improve it ?

[pete]

Skybuck Flying wrote:

Using epsilon's and stuff like should probably be prevented as much as
possible since those are common forms of errors etc and limitations
etc...

.... otherwise you might wind up learning how to use
relevant features of the C standard library.

 

[nobody]

Using epsilon's and stuff like should probably be prevented as much as
possible since those are common forms of errors etc and limitations etc...

You have a gross and dangerous misunderstanding of principles of doing
numerical processing with the computer. I suggest you study the
fundmentals carefully first before wasting your time writing what will
undoubtedly be some ridicolously naive and rather useless code.

 

[John Bode]

Hi,

I needed a method to determine if a point was on a line segment in 2D. So I
googled for some help and so far I have evaluated two methods.

Here's another one to play with. It goes from the assumption that if
the slope from p1 to p2 is the same as the slope from p1 to p3, then p3
is colinear with p1 and p2. Then it checks x coordinates to see if p3
is on the segment between p1 and p2.

It expresses slope as rational numbers, so the arithmetic is integral.
It probably needs to be bulletproofed.

#include <stdio.h>

typedef struct {
long x;
long y;
} point_t;

typedef struct {
long num;
long denom;
} rational_t;

#define REDUCE(frac) \
do { \
long e1=frac.num, \
e2=frac.denom, \
quot, rem = -1; \
long tmp; \
int done = 0; \
while(!done) \
{ \
if (e1 < e2) \
{ \
tmp = e1; \
e1 = e2; \
e2 = tmp; \
} \
quot = e1/e2; \
rem = e1 % e2; \
if (rem) \
e1 = rem; \
else \
done = 1; \
} \
frac.num /= e2; \
frac.denom /= e2; \
} while(0)

/**
* We are assuming that p1 and p2 have been ordered
* so that p1.x <= p2.x
*/
int isOnLine(point_t p1, point_t p2, point_t p3)
{
rational_t m1 = {p3.y - p1.y, p3.x - p1.x};
rational_t m2 = {p2.y - p1.y, p2.x - p1.x};

/**
* special case -- p1.x == p2.x
*/
if (p1.x == p2.x)
{
return p3.x == p1.x &&
(p1.y <= p3.y && p3.y <= p2.y ||
p1.y >= p3.y && p3.y >= p2.y);
}

/**
* Reduce both fractions before comparing. This is where
* any performance issues would be.
*/
REDUCE(m1);
REDUCE(m2);

if (m1.num == m2.num && m1.denom == m2.denom)
{
return p1.x <= p3.x && p3.x <= p2.x;
}

return 0;
}

 

[glen herrmannsfeldt]

(snip)

You are more then welcome to prove this very soon.
I'll shall do or attempt to do two things:

I pick the points (0,0) and (512,997) slightly easier to see as
integers, but you can shift the binary point over and make them floating
point. Find any point on the segment between them.

-- glen

 

[websnarf]

You ll see quite soon the test program is done

Ok, you are pissing people off with statements like this.

Here's the thing, the "epsilon thing" is far and away the most
practical way of doing this; but just as a matter of correctness, you
actually should not choose a constant epsilon. The correct epsilon
will depend on the magnitude of the coordinates for the original
segment points. For example it would be easy to construct an example
where the best correct epsilon is a million.

Ok, perhaps a better way to ask you for more information, is this.

From what sort of *SOURCE* are your input points comming from? Are

they just chosen literally at random? (Using say, the Mersenne Twister
random number generator, which includes floating point random.) Or
(more likely) are they actually constructed to be *ON* the line with
some likelihood, but for some reason you loose track of this fact, and
wish to recalculate it?

This is important because, the *WAY* you construct the point (say, but
being given the x intercept, then computing the y that is supposed to
correspond to it) may actually give an exact matching algorithm without
the need for any epsilons. One could, for example, just try to
reproduce the procedure for finding the point, from one of the
coordinates, and see if the other matches exactly.

The problem with this is that starting with an x and computing a y, or
the other way around will not necessarily yield the same results.
Further more if you compute the point as (alpha * p_0 + (1-alpha)* p_1)
(which might be *WAY* harder to match -- intuitively it seems so to me,
but I don't have a 100% clear idea), you will again get different LSB
round offs.

So I hope you understand that these accuracy issues are pretty integral
to what it is that you want to do. Without more information from you,
I don't believe that anyone can suggest anything else more with any
expectation of it necessarily working better.

 

[Skybuck Flying]

Ok, you are pissing people off with statements like this.

Here's the thing, the "epsilon thing" is far and away the most
practical way of doing this; but just as a matter of correctness, you
actually should not choose a constant epsilon. The correct epsilon
will depend on the magnitude of the coordinates for the original
segment points. For example it would be easy to construct an example
where the best correct epsilon is a million.

Ok, perhaps a better way to ask you for more information, is this.
they just chosen literally at random? (Using say, the Mersenne Twister
random number generator, which includes floating point random.) Or
(more likely) are they actually constructed to be *ON* the line with
some likelihood, but for some reason you loose track of this fact, and
wish to recalculate it?

This is important because, the *WAY* you construct the point (say, but
being given the x intercept, then computing the y that is supposed to
correspond to it) may actually give an exact matching algorithm without
the need for any epsilons. One could, for example, just try to
reproduce the procedure for finding the point, from one of the
coordinates, and see if the other matches exactly.

The problem with this is that starting with an x and computing a y, or
the other way around will not necessarily yield the same results.
Further more if you compute the point as (alpha * p_0 + (1-alpha)* p_1)
(which might be *WAY* harder to match -- intuitively it seems so to me,
but I don't have a 100% clear idea), you will again get different LSB
round offs.

So I hope you understand that these accuracy issues are pretty integral
to what it is that you want to do. Without more information from you,
I don't believe that anyone can suggest anything else more with any
expectation of it necessarily working better.

The test program is "done" well not really but a first version is available
on the net.

VerificationProgramVersion009.zip

https://sourceforge.net/project/showfiles.php?group_id=53726

See the other usenet thread called "VerificationProgram009 etc" for more
details etc

Bye,
Skybuck.

 

[Skybuck Flying]

I pick the points (0,0) and (512,997) slightly easier to see as
integers, but you can shift the binary point over and make them floating
point. Find any point on the segment between them.

The test program is "done" well not really but a first version is available
on the net.

Which also includes your line. "My" second implementation has problems with
this line. The second method uses some sort of dot product or cross product.
Which is kinda surprising because I think Nicholas's implementation also
works like that... oh well

My first implementation however... which I do fully understand still works
perfectly =D

Can you find/think up any other possibly problems ?

See this link for the source code, binaries, data, etc, etc, etc.

VerificationProgramVersion009.zip

https://sourceforge.net/project/showfiles.php?group_id=53726

See the other usenet thread called "VerificationProgram009 etc" for more
details etc

Bye,
Skybuck.

 

[Skybuck Flying]

Here's another one to play with. It goes from the assumption that if
the slope from p1 to p2 is the same as the slope from p1 to p3, then p3
is colinear with p1 and p2. Then it checks x coordinates to see if p3
is on the segment between p1 and p2.

It expresses slope as rational numbers, so the arithmetic is integral.
It probably needs to be bulletproofed.

#include <stdio.h>

typedef struct {
long x;
long y;
} point_t;

typedef struct {
long num;
long denom;
} rational_t;

#define REDUCE(frac) \
do { \
long e1=frac.num, \
e2=frac.denom, \
quot, rem = -1; \
long tmp; \
int done = 0; \
while(!done) \
{ \
if (e1 < e2) \
{ \
tmp = e1; \
e1 = e2; \
e2 = tmp; \
} \
quot = e1/e2; \
rem = e1 % e2; \
if (rem) \
e1 = rem; \
else \
done = 1; \
} \
frac.num /= e2; \
frac.denom /= e2; \
} while(0)

/**
* We are assuming that p1 and p2 have been ordered
* so that p1.x <= p2.x
*/
int isOnLine(point_t p1, point_t p2, point_t p3)
{
rational_t m1 = {p3.y - p1.y, p3.x - p1.x};
rational_t m2 = {p2.y - p1.y, p2.x - p1.x};

/**
* special case -- p1.x == p2.x
*/
if (p1.x == p2.x)
{
return p3.x == p1.x &&
(p1.y <= p3.y && p3.y <= p2.y ||
p1.y >= p3.y && p3.y >= p2.y);
}

/**
* Reduce both fractions before comparing. This is where
* any performance issues would be.
*/
REDUCE(m1);
REDUCE(m2);

if (m1.num == m2.num && m1.denom == m2.denom)
{
return p1.x <= p3.x && p3.x <= p2.x;
}

return 0;
}

Hi, cool stuff. I haven't had the time yet to look into this.

It would help if someone could make it more suited for my test program and
maybe make a nice little dll ?

Anyway maybe somebody can convert this to delphi as well.

The test program is "done" well not really but a first version is available
on the net.

VerificationProgramVersion009.zip

https://sourceforge.net/project/showfiles.php?group_id=53726

See the other usenet thread called "VerificationProgram009 etc" for more
details etc

Bye,
Skybuck.

 

[alanglloyd]

Why don't you make a region slightly larger than the line, and use
PtInRegion to check if the click is on the line.

var
hndRgn : hRGN;

const
LineSt : TPoint = (X:0; Y:0);
LineEnd : TPoint = (X:997; Y:512);
LineWidth : integer = 2;

procedure TForm1.FormCreate(Sender: TObject);
var
PointArray : array[0..3] of TPoint;
LW : integer;
begin
LW := LineWidth;
PointArray[0] := Point(LineSt.X - LW, LineSt.Y - LW);
PointArray[1] := Point(LineEnd.X + LW, LineEnd.Y - LW);
PointArray[2] := Point(LineEnd.X + LW, LineEnd.Y + LW);
PointArray[3] := Point(LineSt.X - LW, LineSt.Y + LW);
hndRgn := CreatePolygonRgn(PointArray, 4, WINDING);
end;

procedure TForm1.FormMouseUp(Sender: TObject; Button: TMouseButton;
Shift: TShiftState; X, Y: Integer);
var
OnLine : boolean;
begin
StartPeriod; // timer start
OnLine := PtInRegion(hndRgn, X, Y);
EndPeriod; // timer stop
{indicate if on line}
if OnLine then
Shape1.Brush.Color := clLime
else
Shape1.Brush.Color := clRed;
{display time to check in region}
Label1.Caption := Format('%d microsecs',
[trunc(ScaleTime(ElapsedTime, ttMicroSecs))]);
end;

This gives about 13 microsecs for the first click (2 microsecs for
later ones) at the start of the line, and 17 microsecs for the first
click (7 microsecs for later ones) at the end of the line. This is with
a 2.5GHz PC.

I can't believe that this is too slow for you <g>.

Alan Lloyd

 

[Sander Vesik]

Skybuck Flying wrote:

(snip)



You really need to tell us what "this case" is.

You expect everyone to help you, but don't seem
interested in helping much yourself.

Its rather clear he is interested in generating lots of usenet traffic
and not much else...

 

[Skybuck Flying]

Here's another one to play with. It goes from the assumption that if
the slope from p1 to p2 is the same as the slope from p1 to p3, then p3
is colinear with p1 and p2. Then it checks x coordinates to see if p3
is on the segment between p1 and p2.

It expresses slope as rational numbers, so the arithmetic is integral.
It probably needs to be bulletproofed.

#include <stdio.h>

typedef struct {
long x;
long y;
} point_t;

typedef struct {
long num;
long denom;
} rational_t;

#define REDUCE(frac) \
do { \
long e1=frac.num, \
e2=frac.denom, \
quot, rem = -1; \
long tmp; \
int done = 0; \
while(!done) \
{ \
if (e1 < e2) \
{ \
tmp = e1; \
e1 = e2; \
e2 = tmp; \
} \
quot = e1/e2; \
rem = e1 % e2; \
if (rem) \
e1 = rem; \
else \
done = 1; \
} \
frac.num /= e2; \
frac.denom /= e2; \
} while(0)

/**
* We are assuming that p1 and p2 have been ordered
* so that p1.x <= p2.x
*/
int isOnLine(point_t p1, point_t p2, point_t p3)
{
rational_t m1 = {p3.y - p1.y, p3.x - p1.x};
rational_t m2 = {p2.y - p1.y, p2.x - p1.x};

/**
* special case -- p1.x == p2.x
*/
if (p1.x == p2.x)
{
return p3.x == p1.x &&
(p1.y <= p3.y && p3.y <= p2.y ||
p1.y >= p3.y && p3.y >= p2.y);
}

/**
* Reduce both fractions before comparing. This is where
* any performance issues would be.
*/
REDUCE(m1);
REDUCE(m2);

if (m1.num == m2.num && m1.denom == m2.denom)
{
return p1.x <= p3.x && p3.x <= p2.x;
}

return 0;
}

Ok, your method has been added to the verification program.

However during verification the algorithm hangs at this data:

Verifieing:
Type index: 1
Type description: Diagonal line with point inside line segment
Data index: 2
Data description: positive area, third crossed

Which is this one:

AddGeneratedVerificationData( 'positive area, third crossed', 200,2000,
3000,100, 0.3, 0.0 );

Which means:
X1 := 200;
Y1 := 2000;
X2 := 3000;
Y2 := 100;
PX := 1040; // generated/calculated
PY := 1430; // generated/calculated

I put your source code in a library/DLL which looks like this:

JOHNBODERATIONALMETHOD_API BOOL __stdcall PointOnLineSegmentIn2D(
double StartX, double StartY,
double EndX, double EndY,
double PointX, double PointY )
{
point_t p1;
point_t p2;
point_t p3;
if (StartX <= EndX)
{
p1.x = long(StartX); // typecast to prevent loss of data warning
p1.y = long(StartY); // typecast to prevent loss of data warning
p2.x = long(EndX); // typecast to prevent loss of data warning
p2.y = long(EndY); // typecast to prevent loss of data warning
} else
{
p1.x = long(EndX); // typecast to prevent loss of data warning
p1.y = long(EndY); // typecast to prevent loss of data warning
p2.x = long(StartX); // typecast to prevent loss of data warning
p2.y = long(StartY); // typecast to prevent loss of data warning
}
p3.x = long(PointX); // typecast to prevent loss of data warning
p3.y = long(PointY); // typecast to prevent loss of data warning
return isOnLine(p1,p2,p3);
}

It might be hanging because Y2 is smaller than Y1 ? or maybe there is
another problem ?

Bye,
Skybuck.

 

[Skybuck Flying]

Hi John, your posted method has now also been added to the
program/project/etc in
VerificationProgram0.10 for PointOnLineSegmentIn2D =D
(Though the dll has been disabled by simply renaming it until the hang/bug
is fixed )

Hello Peeps.

Here is what's new

version 0.10 created on 18 and 19 september 2005 by Skybuck Flying

+ Some code moved to other/new units etc.
+ New source and binary implementation added: Nils from www.simdesign.nl
+ New source and binary implementation added: Skybuck Flying method 3
+ New binary implementation added: JohnBodeRationalMethod (source in visual
cpp).
+ ProjectGroups added
+ Each verification is displayed in case a binary plug in hangs.. then at
least
it's possible to see at which verification it hangs.

Get the source and binaries from:

https://sourceforge.net/project/showfiles.php?group_id=53726

=D

Bye,
Skybuck.

 

[Skybuck Flying]

Why don't you make a region slightly larger than the line, and use
PtInRegion to check if the click is on the line.

Not too shabby

This would/could be pretty good for gui's where the user needs to pick a
line etc.

That's just one application of this method.

Another application is for massive calculations of all kinds of geometry
algorithms etc.

Your method is kinda funny since it uses windows to do it's thing. Tomorrow
or so I will add it as well to see how it performance, speed wise, compared
to all the other methods =D

That's going to be much fun

I also wonder how stable it would be... probably pretty stable though I
also wonder how large or small the values could be.

So this method works with integers etc... small values would get rounded to
zero or so or one like
0.05 and 0.0002 etc.. and point 0.06 whatever... then this method would
always return true or something like that could be fun and maybe even
handy as well =D

Bye,
Skybuck.

var
hndRgn : hRGN;

const
LineSt : TPoint = (X:0; Y:0);
LineEnd : TPoint = (X:997; Y:512);
LineWidth : integer = 2;

procedure TForm1.FormCreate(Sender: TObject);
var
PointArray : array[0..3] of TPoint;
LW : integer;
begin
LW := LineWidth;
PointArray[0] := Point(LineSt.X - LW, LineSt.Y - LW);
PointArray[1] := Point(LineEnd.X + LW, LineEnd.Y - LW);
PointArray[2] := Point(LineEnd.X + LW, LineEnd.Y + LW);
PointArray[3] := Point(LineSt.X - LW, LineSt.Y + LW);
hndRgn := CreatePolygonRgn(PointArray, 4, WINDING);
end;

procedure TForm1.FormMouseUp(Sender: TObject; Button: TMouseButton;
Shift: TShiftState; X, Y: Integer);
var
OnLine : boolean;
begin
StartPeriod; // timer start
OnLine := PtInRegion(hndRgn, X, Y);
EndPeriod; // timer stop
{indicate if on line}
if OnLine then
Shape1.Brush.Color := clLime
else
Shape1.Brush.Color := clRed;
{display time to check in region}
Label1.Caption := Format('%d microsecs',
[trunc(ScaleTime(ElapsedTime, ttMicroSecs))]);
end;

This gives about 13 microsecs for the first click (2 microsecs for
later ones) at the start of the line, and 17 microsecs for the first
click (7 microsecs for later ones) at the end of the line. This is with
a 2.5GHz PC.

Nice one

I can't believe that this is too slow for you <g>.

For one click no

For massive calculations I dont know lol =D

Bye,
Skybuck.

 

[Randy Howard]

(e-mail address removed) wrote
(in article

Ok, you are pissing people off with statements like this.

Shouldn't be a surprise. If you search newsgroup archives for
"Skybuck" you'll find his been doing that sort of thing in a lot
of different groups for a long time.

 

[Alan Balmer]

(e-mail address removed) wrote
(in article


Shouldn't be a surprise. If you search newsgroup archives for
"Skybuck" you'll find his been doing that sort of thing in a lot
of different groups for a long time.

Yep. I changed his filter from specific newsgroup to global a long
time ago.

 

[Skybuck Flying]

I pick the points (0,0) and (512,997) slightly easier to see as
integers, but you can shift the binary point over and make them floating
point. Find any point on the segment between them.

I wrote a little program which uses rational numbers and floating point
numbers and shows the difference.

The rational number version is able to perfectly calculate the point.
The floating point version fails miserably at it

The rational number version is about 20 times slower than the floating point
version

program started

(Rational) PointX in rational notation: 768/5
(Rational) PointY in rational notation: 2991/10

(Rational) PointX in real notation: 153.6000000000000000
(Rational) PointY in real notation: 299.1000000000000000
(Rational) ticks: 282
(Rational) time (in seconds): 0.0000787809623849

(Floating point) PointX: 153.5999999999999940
(Floating point) PointY: 299.0999999999999660
(Floating point) ticks: 14
(Floating point) time (in seconds): 0.0000039111116078

press enter to exit

Here is the source code. Thanks to delphi fundamentals for the Trational
class etc =D

Now for the final question

Can you still find a point which would in theory be on the line segment but
would still be mis-calculated by the rational number version ?

program TestRationalNumbers;

{$APPTYPE CONSOLE}

uses
SysUtils,
cUtils in 'Utils\cUtils.pas',
cStrings in 'Utils\cStrings.pas',
cRational in 'Maths\cRational.pas',
cMaths in 'Maths\cMaths.pas',
Windows;


// rational number versions:

procedure calc_perpendicular( var x, y : Trational );
var
temp_x, temp_y : Trational;
begin
temp_x := Trational.Create;
temp_y := Trational.Create;

// temp_x := -y;
temp_x.Assign( y );
temp_x.Negate;

// temp_y := x;
temp_y.Assign( x );

// x := temp_x;
// y := temp_y;
x.Assign( temp_x );
y.Assign( temp_y );

temp_x.Destroy;
temp_y.Destroy;
end;

procedure calc_normalize( var x, y : Trational );
var
temp_squared_x : Trational;
temp_squared_y : Trational;
temp_length : Trational;

begin
temp_squared_x := Trational.Create;
temp_squared_y := Trational.Create;
temp_length := Trational.Create;

// temp_length := sqrt( (x*x) + (y*y) );
temp_squared_x.Assign( x );
temp_squared_y.Assign( y );

temp_squared_x.Sqr;
temp_squared_y.Sqr;

temp_length.Assign( temp_squared_x );
temp_length.Add( temp_squared_y );

temp_length.Sqrt;

// if temp_length<>0 then
if not temp_length.IsZero then
begin
// x := x / temp_length;
// y := y / temp_length;

x.Divide( temp_length );
y.Divide( temp_length );
end;

temp_squared_x.Destroy;
temp_squared_y.Destroy;
temp_length.Destroy;
end;

procedure calc_normal( _x1,_y1, _x2,_y2 : Trational; var _normal_x,
_normal_y : Trational );
begin
// _normal_x := _x2-_x1;
// _normal_y := _y2-_y1;

_normal_x.Assign( _x2 );
_normal_y.Assign( _y2 );

_normal_x.Subtract( _x1 );
_normal_y.Subtract( _y1 );

// calculate perpendicular "vector".
calc_perpendicular( _normal_x, _normal_y );

// divide vector by it's own length there by creating a "normal"
calc_normalize( _normal_x, _normal_y );
end;

// original double versions
// using flawed floating point shit

// helper procedure to generate points on or near line segments
procedure GeneratePointForLineSegment( StartX, StartY,
EndX, EndY : Trational;
DistanceFromStart : Trational; // distance from start
// < 0 is before start
// 0 is on start
// 0.5 is on center
// 1 is one end
// > 1 is after end
DistanceFromLine : Trational; // closest distance from line

// follows the normal
// positive is left side of line
// zero is on line
// negative is right side of line
var px, py : Trational );
var
NX, NY : Trational;

begin
NX := Trational.Create;
NY := Trational.Create;

// theory of method1
// PX := StartX + DistanceFromStart * (EndX - StartX);
// PY := StartY + DistanceFromStart * (EndY - StartY);

PX.Assign( EndX );
PX.Subtract( StartX );
PX.Multiply( DistanceFromStart );
PX.Add( StartX );

PY.Assign( EndY );
PY.Subtract( StartY );
PY.Multiply( DistanceFromStart );
PY.Add( StartY );

// create a perpendicular point... we do this by calculating a normal and
multiplieing
// this with the distance from the line and adding this to px and py

// calculate normilized normal for line
calc_normal( StartX,StartY, EndX,EndY, NX,NY );

// offset the point along the normal with distance from line
// PX := PX + NX * DistanceFromLine;
// PY := PY + NY * DistanceFromLine;

NX.Multiply( DistanceFromLine );
NY.Multiply( DistanceFromLine );

PX.Add( NX );
PY.Add( NY );

NX.Destroy;
NY.Destroy;
end;

procedure org_calc_perpendicular( var x, y : double );
var
temp_x, temp_y : double;
begin
temp_x := -y;
temp_y := x;
x := temp_x;
y := temp_y;
end;

procedure org_calc_normalize( var x, y : double );
var
temp_length : double;
begin
temp_length := sqrt( (x*x) + (y*y) );
if temp_length<>0 then
begin
x := x / temp_length;
y := y / temp_length;
end;
end;

procedure org_calc_normal( _x1,_y1, _x2,_y2 : double; var _normal_x,
_normal_y : double );
begin
_normal_x := _x2-_x1;
_normal_y := _y2-_y1;

// calculate perpendicular "vector".
org_calc_perpendicular( _normal_x, _normal_y );

// divide vector by it's own length there by creating a "normal"
org_calc_normalize( _normal_x, _normal_y );
end;

// helper procedure to generate points on or near line segments
procedure org_GeneratePointForLineSegment( StartX, StartY,
EndX, EndY : double;
DistanceFromStart : double; // distance from start
// < 0 is before start
// 0 is on start
// 0.5 is on center
// 1 is one end
// > 1 is after end
DistanceFromLine : double; // closest distance from line

// follows the normal
// positive is left side of line
// zero is on line
// negative is right side of line
var px, py : double );
var
NX, NY : double;
begin
// theory of method1
PX := StartX + DistanceFromStart * (EndX - StartX);
PY := StartY + DistanceFromStart * (EndY - StartY);

// create a perpendicular point... we do this by calculating a normal and
multiplieing
// this with the distance from the line and adding this to px and py

// calculate normilized normal for line
org_calc_normal( StartX,StartY, EndX,EndY, NX,NY );

// offset the point along the normal with distance from line
PX := PX + NX * DistanceFromLine;
PY := PY + NY * DistanceFromLine;
end;

var
vHighPerformanceCounterFrequency : int64;
vCounter1 : int64;
vCounter2 : int64;

vStartX,
vStartY,
vEndX,
vEndY,
vPointX,
vPointY,
vDistanceFromStart,
vDistanceFromLine : Trational;

vOrgStartX,
vOrgStartY,
vOrgEndX,
vOrgEndY,
vOrgPointX,
vOrgPointY,
vOrgDistanceFromStart,
vOrgDistanceFromLine : double;

vMeasurementsEnabled : boolean;

vRationalTicks : int64;
vRationalTime : double;

vFloatingPointTicks : int64;
vFloatingPointTime : double;

begin
writeln('program started');

vMeasurementsEnabled := true;;

vRationalTicks := 0;
vRationalTime := 0;

vFloatingPointTicks := 0;
vFloatingPointTime := 0;

if vMeasurementsEnabled then
begin
if not QueryPerformanceFrequency( vHighPerformanceCounterFrequency ) then
begin
writeln('High performance counter not present, speed measurements
disabled');
vMeasurementsEnabled := false;
end;
end;

vStartX := Trational.Create;
vStartY := Trational.Create;
vEndX := Trational.Create;
vEndY := Trational.Create;
vPointX := Trational.Create;
vPointY := Trational.Create;
vDistanceFromStart := Trational.Create;
vDistanceFromLine := Trational.Create;

vStartX.Assign( 0, 1 );
vStartY.Assign( 0, 1 );

vEndX.Assign( 512, 1 );
vEndY.Assign( 997, 1 );

vDistanceFromStart.Assign( 3, 10 );
vDistanceFromLine.Assign( 0, 1 );

if vMeasurementsEnabled then
begin
QueryPerformanceCounter( vCounter1 );
end;

GeneratePointForLineSegment(
vStartX, vStartY,
vEndX, vEndY,
vDistanceFromStart,
vDistanceFromLine,
vPointX, vPointY );

if vMeasurementsEnabled then
begin
QueryPerformanceCounter( vCounter2 );

vRationalTicks := vCounter2-vCounter1;
vRationalTime := ( (vCounter2-vCounter1)/
vHighPerformanceCounterFrequency );
end;

writeln;
writeln( '(Rational) PointX in rational notation: ', vPointX.AsString );
writeln( '(Rational) PointY in rational notation: ', vPointY.AsString );

writeln;
writeln( '(Rational) PointX in real notation: ', vPointX.AsReal : 16 :
16 );
writeln( '(Rational) PointY in real notation: ', vPointY.AsReal : 16 :
16 );

if vMeasurementsEnabled then
begin
writeln('(Rational) ticks: ', vRationalTicks );
writeln('(Rational) time (in seconds): ', vRationalTime : 16 : 16 );
end;

vOrgStartX := 0.0;
vOrgStartY := 0.0;
vOrgEndX := 512.0;
vOrgEndY := 997.0;
vOrgDistanceFromStart := 0.3;
vOrgDistanceFromLine := 0.0;

if vMeasurementsEnabled then
begin
QueryPerformanceCounter( vCounter1 );
end;

org_GeneratePointForLineSegment(
vOrgStartX, vOrgStartY,
vOrgEndX, vOrgEndY,
vOrgDistanceFromStart,
vOrgDistanceFromLine,
vOrgPointX, vOrgPointY );

if vMeasurementsEnabled then
begin
QueryPerformanceCounter( vCounter2 );

vFloatingPointTicks := vCounter2-vCounter1;
vFloatingPointTime := ( (vCounter2-vCounter1)/
vHighPerformanceCounterFrequency );
end;

writeln;
writeln( '(Floating point) PointX: ', vOrgPointX : 16 : 16 );
writeln( '(Floating point) PointY: ', vOrgPointY : 16 : 16 );

if vMeasurementsEnabled then
begin
writeln('(Floating point) ticks: ', vFloatingPointTicks );
writeln('(Floating point) time (in seconds): ', vFloatingPointTime : 16 :
16 );
end;

vStartX.Destroy;
vStartY.Destroy;
vEndX.Destroy;
vEndY.Destroy;
vPointX.Destroy;
vPointY.Destroy;
vDistanceFromStart.Destroy;
vDistanceFromLine.Destroy;

writeln;
writeln('press enter to exit');
readln;

writeln('program finished');
end.

Bye,
Skybuck.

 

[Skybuck Flying]

[Massive crossposting reduced to single group. Please, Skybuck, try
to follow netiquette a bit better in the future. If you feel a need
to X-post, at be a bit more selective to where --- this had *no*
business being posted to comp.arch or sci.electronics. And *always*
select a single group for the follow-ups.]

No, I like to have input from all those newsgroups since this is about
optimizations and those newsgroup are related to software and hardware
optimizations. C/Delphi/Asm for software/cpu optimizations and
comp.arch/asm/eletronics about hardware optimizations and algorithms
possibly about algorithm optimization. Those newsgroups are likely to
contain people who know something about this.

I'm with Hans on this one. Skybuck, by your logic you should have
posted to ALL newsgroups, since there's probably someone in every
newsgroup who will be somehow "related to optimization". :-(

Comp.graphics.algorithms was the right place. Comp.lang.c is hard to
justify, since your question has nothing to do with C. But there's no
way to justify posting to sci.electronics.design, comp.arch, or
alt.comp.lang.borland-delphi. Nor comp.lang.asm, unless you want a
response written in assembly language.

That said, other appropriate forums might have included comp.programming
or *maybe* comp.theory.

...

A 2D line segment defined by (x1,y1) - (x2,y2)

Integer or floating point? In either case, how big are your points and
how fat is the line? Is it acceptable to return a false positive or
negative? If so, how often? If it's not acceptable, then you have a
problem. On a any machine that lacks infinite numerical precision, you
can never know whether an infinitely small point lies on an infinitely
thin line segment.

In the end, you're going to have to do something like this:

"Draw two lines, one from each point in the line segment to the
candidate point. If these lines' slopes are equal but differ in sign
(e.g. X and -X), then your point lies on the segment."

Of course, what does it mean for slope values to be "equal" on a digital
computer? Unless you can represent the problem entirely symbolically,
"equal" will be necessarily limited by the representation of the data
and the accuracy of the math. In the real world, epsilon is a necessary
evil.

As precise as possible. (Floating point precision, no rounding)

This is an oxymoron. Floating point is inescapably about rounding.
Digital representations of floating point are approximations, as are the
math operations that manipulate them. The programmer has to manage
numerical imprecision on computers or live with wrong answers.

For some perspective, read David Goldberg's "What Every Computer
Scientist Should Know About Floating-Point Arithmetic":

http://docs.sun.com/source/806-3568/ncg_goldberg.html

Yes I am sure you all right about this... But instead of wasting my time
reading through this highly theoretical shit I rather find practical
solutions to this potential floating point rounding problem. And using
epsilon's is probably not it since it still leaves room for errors etc.

So far I have found one nice solution... Trational class from delphi
fundamentals.

Another solution might be a fixed point library which might be faster =D

Bye,
Skybuck.

 

[Skybuck Flying]

Sorry, I missed some of this code from my app:

function PointOnLine2D(x1, y1, x2, y2, x3, y3: double): boolean;
begin
result :=
(((x3 >= x1) and (x3 <= x2)) or ((x3 >= x2) and (x3 <= x1))) and
(((y3 >= y1) and (y3 <= y2)) or ((y3 >= y2) and (y3 <= y1))) and
(abs(((x2 - x1) * (y3 - y1)) - ((x3 - x1) * (y2 - y1))) < 0.00001);
end;

Ok I first I didn't understand this code, but not I get it.

The other C example was like this example.

The C example compares the first multiplication result with the second
multiplication result.

If they are the same then it's true else false.

Because floating point rounding errors could occur the C example could fail.
Which was demonstrated by my method 2 which was the ported C version.

So your/the fastgeo version solves this rounding error by using an
epsilon...

I still think it's a shabby solution since it could horribly fail if the
points are small enough.

Maybe something like:
0.00000003245
0.000000003245
0.000000003456
0.000000000123

etc...

OUCH



Bye,
Skybuck =D

 

[Skybuck Flying]

Skybuck Flying wrote:

(snip)



You really need to tell us what "this case" is.

I don't see the relevance of the case.

The math is clear and simple and is case independant =D

You expect everyone to help you, but don't seem
interested in helping much yourself.

No on the contrary See above lol.

Bye,
Skybuck.

 

[Nicholas Sherlock]

The epsilon method is flawed in my mind since it assumes the points are
quite large but in fact could be quite small thereby completely failing.

Moron.



Nicholas Sherlock