Skip to content

Commit

Permalink
Updated IsCollinear functions in C# and Delphi (#834)
Browse files Browse the repository at this point in the history
Otherwise a minor code tidy and updated code comments
  • Loading branch information
AngusJohnson committed May 12, 2024
1 parent fdf3700 commit 235879b
Show file tree
Hide file tree
Showing 3 changed files with 148 additions and 34 deletions.
45 changes: 24 additions & 21 deletions CPP/Clipper2Lib/include/clipper2/clipper.core.h
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
/*******************************************************************************
* Author : Angus Johnson *
* Date : 11 May 2024 *
* Date : 12 May 2024 *
* Website : http://www.angusj.com *
* Copyright : Angus Johnson 2010-2024 *
* Purpose : Core Clipper Library structures and functions *
Expand Down Expand Up @@ -659,42 +659,43 @@ namespace Clipper2Lib
CheckPrecisionRange(precision, error_code);
}

inline int Sign(int64_t x)
inline int TriSign(int64_t x) // returns 0, 1 or -1
{
return (x > 0) - (x < 0);
return (x > 0) - (x < 0);
}

inline uint64_t CalcOverflowCarry(uint64_t a, uint64_t b)
inline uint64_t CalcOverflowCarry(uint64_t a, uint64_t b) // #834
{
const uint64_t aLo = a & 0xFFFFFFFF;
//const uint64_t aHi = a & 0xFFFFFFFF00000000;
const uint64_t aHiShr = a >> 32;

const uint64_t bLo = b & 0xFFFFFFFF;
//const uint64_t bHi = b & 0xFFFFFFFF00000000;
const uint64_t bHiShr = b >> 32;

// given aLo = (a & 0xFFFFFFFF) and
// aHi = (a & 0xFFFFFFFF00000000) and similarly with b, then
// a * b == (aHi + aLo) * (bHi + bLo)
// a * b == (aHi * bHi) + (aHi * bLo) + (aLo * bHi) + (aLo * bLo)
//int overflow of 64bit a * 64bit b ==>
return aHiShr * bHiShr + ((aHiShr * bLo) >> 32) + ((bHiShr * aLo) >> 32);

const uint64_t aLo = a & 0xFFFFFFFF;
const uint64_t aHi = a >> 32; // this avoids repeating shifts
const uint64_t bLo = b & 0xFFFFFFFF;
const uint64_t bHi = b >> 32;
// integer overflow of multiplying the unsigned 64bits a and b ==>
return aHi * bHi + ((aHi * bLo) >> 32) + ((bHi * aLo) >> 32);
}

// returns true if (and only if) a * b == c * d
inline bool ProductIsEqual(int64_t a, int64_t b, int64_t c, int64_t d)
inline bool ProductsAreEqual(int64_t a, int64_t b, int64_t c, int64_t d)
{
// nb: unsigned values will be needed for CalcOverflowCarry()
const auto abs_a = static_cast<uint64_t>(std::abs(a));
const auto abs_b = static_cast<uint64_t>(std::abs(b));
const auto abs_c = static_cast<uint64_t>(std::abs(c));
const auto abs_d = static_cast<uint64_t>(std::abs(d));

// the multiplications here can potentially overflow
// but any overflows will be compared further below.
// the multiplications here can potentially overflow, but
// any overflows will be compared using CalcOverflowCarry()
const auto abs_ab = abs_a * abs_b;
const auto abs_cd = abs_c * abs_d;

const auto sign_ab = Sign(a) * Sign(b);
const auto sign_cd = Sign(c) * Sign(d);
// nb: it's important to differentiate 0 values here from other values
const auto sign_ab = TriSign(a) * TriSign(b);
const auto sign_cd = TriSign(c) * TriSign(d);

const auto carry_ab = CalcOverflowCarry(abs_a, abs_b);
const auto carry_cd = CalcOverflowCarry(abs_c, abs_d);
Expand All @@ -709,10 +710,12 @@ namespace Clipper2Lib
const auto b = pt2.y - sharedPt.y;
const auto c = sharedPt.y - pt1.y;
const auto d = pt2.x - sharedPt.x;

return ProductIsEqual(a, b, c, d);
// When checking for collinearity with very large coordinate values
// then ProductsAreEqual is more accurate than using CrossProduct.
return ProductsAreEqual(a, b, c, d);
}


template <typename T>
inline double CrossProduct(const Point<T>& pt1, const Point<T>& pt2, const Point<T>& pt3) {
return (static_cast<double>(pt2.x - pt1.x) * static_cast<double>(pt3.y -
Expand Down
63 changes: 58 additions & 5 deletions CSharp/Clipper2Lib/Clipper.Core.cs
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
/*******************************************************************************
* Author : Angus Johnson *
* Date : 10 May 2024 *
* Date : 12 May 2024 *
* Website : http://www.angusj.com *
* Copyright : Angus Johnson 2010-2024 *
* Purpose : Core structures and functions for the Clipper Library *
Expand Down Expand Up @@ -585,7 +585,7 @@ public static double CrossProduct(Point64 pt1, Point64 pt2, Point64 pt3)
return ((double) (pt2.X - pt1.X) * (pt3.Y - pt2.Y) -
(double) (pt2.Y - pt1.Y) * (pt3.X - pt2.X));
}

#if USINGZ
public static Path64 SetZ(Path64 path, long Z)
{
Expand All @@ -595,22 +595,75 @@ public static Path64 SetZ(Path64 path, long Z)
}
#endif

[MethodImpl(MethodImplOptions.AggressiveInlining)]
internal static void CheckPrecision(int precision)
{
if (precision < -8 || precision > 8)
throw new Exception(precision_range_error);
}

[MethodImpl(MethodImplOptions.AggressiveInlining)]
internal static bool IsAlmostZero(double value)
{
return (Math.Abs(value) <= floatingPointTolerance);
}

[MethodImpl(MethodImplOptions.AggressiveInlining)]
internal static bool IsCollinear(Point64 pt1, Point64 pt2, Point64 pt3)
internal static int TriSign(long x) // returns 0, 1 or -1
{
return (pt2.X - pt1.X) * (pt3.Y - pt2.Y) ==
(pt2.Y - pt1.Y) * (pt3.X - pt2.X);
if (x < 0) return -1;
else if (x > 1) return 1;
else return 0;
}

[MethodImpl(MethodImplOptions.AggressiveInlining)]
internal static ulong CalcOverflowCarry(ulong a, ulong b) // #834
{
// given aLo = (a & 0xFFFFFFFF) and
// aHi = (a & 0xFFFFFFFF00000000) and similarly with b, then
// a * b == (aHi + aLo) * (bHi + bLo)
// a * b == (aHi * bHi) + (aHi * bLo) + (aLo * bHi) + (aLo * bLo)

ulong aLo = a & 0xFFFFFFFF;
ulong aHi = a >> 32;
ulong bLo = b & 0xFFFFFFFF;
ulong bHi = b >> 32;
// integer overflow of multiplying the unsigned 64bits a and b ==>
return aHi * bHi + ((aHi * bLo) >> 32) + ((bHi * aLo) >> 32);
}

// returns true if (and only if) a * b == c * d
internal static bool ProductsAreEqual(long a, long b, long c, long d)
{
// nb: unsigned values will be needed for CalcOverflowCarry()
ulong absA = (ulong) Math.Abs(a);
ulong absB = (ulong) Math.Abs(b);
ulong absC = (ulong) Math.Abs(c);
ulong absD = (ulong) Math.Abs(d);
// the multiplications here can potentially overflow, but
// any overflows will be compared using CalcOverflowCarry()
ulong abs_ab = absA * absB;
ulong abs_cd = absC * absD;

// nb: it's important to differentiate 0 values here from other values
int sign_ab = TriSign(a) * TriSign(b);
int sign_cd = TriSign(c) * TriSign(d);

ulong carry_ab = CalcOverflowCarry(absA, absB);
ulong carry_cd = CalcOverflowCarry(absC, absD);
return abs_ab == abs_cd && sign_ab == sign_cd && carry_ab == carry_cd;
}

[MethodImpl(MethodImplOptions.AggressiveInlining)]
internal static bool IsCollinear(Point64 pt1, Point64 sharedPt, Point64 pt2)
{
long a = sharedPt.X - pt1.X;
long b = pt2.Y - sharedPt.Y;
long c = sharedPt.Y - pt1.Y;
long d = pt2.X - sharedPt.X;
// When checking for collinearity with very large coordinate values
// then ProductsAreEqual is more accurate than using CrossProduct.
return ProductsAreEqual(a, b, c, d);
}

[MethodImpl(MethodImplOptions.AggressiveInlining)]
Expand Down
74 changes: 66 additions & 8 deletions Delphi/Clipper2Lib/Clipper.Core.pas
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@

(*******************************************************************************
* Author : Angus Johnson *
* Date : 3 May 2024 *
* Date : 12 May 2024 *
* Website : http://www.angusj.com *
* Copyright : Angus Johnson 2010-2024 *
* Purpose : Core Clipper Library module *
Expand Down Expand Up @@ -154,8 +154,7 @@ function IsPositive(const path: TPath64): Boolean; overload;
function IsPositive(const path: TPathD): Boolean; overload;
{$IFDEF INLINING} inline; {$ENDIF}

function IsCollinear(const pt1, pt2, pt3: TPoint64): Boolean;
{$IFDEF INLINING} inline; {$ENDIF}
function IsCollinear(const pt1, sharedPt, pt2: TPoint64): Boolean;

function CrossProduct(const pt1, pt2, pt3: TPoint64): double; overload;
{$IFDEF INLINING} inline; {$ENDIF}
Expand Down Expand Up @@ -1864,18 +1863,77 @@ function IsPositive(const path: TPathD): Boolean;
end;
//------------------------------------------------------------------------------

function TriSign(val: Int64): integer; // returns 0, 1 or -1
{$IFDEF INLINING} inline; {$ENDIF}
begin
if (val < 0) then Result := -1
else if (val > 1) then Result := 1
else Result := 0;
end;
//------------------------------------------------------------------------------

function CalcOverflowCarry(a, b: UInt64): UInt64; // #834
{$IFDEF INLINING} inline; {$ENDIF}
var
aLo, aHi, bLo, bHi: UInt64;
begin
// given aLo = (a and $FFFFFFFF) and
// aHi = (a and $FFFFFFFF00000000) and similarly with b, then
// a * b == (aHi + aLo) * (bHi + bLo)
// a * b == (aHi * bHi) + (aHi * bLo) + (aLo * bHi) + (aLo * bLo)
aLo := a and $FFFFFFFF;
aHi := a shr 32; // this avoids multiple shifts
bLo := b and $FFFFFFFF;
bHi := b shr 32;
//integer overflow of multiplying the unsigned 64bits a and b ==>
Result := (aHi * bHi) + ((aHi * bLo) shr 32) + ((bHi * aLo) shr 32);
end;
//------------------------------------------------------------------------------

{$OVERFLOWCHECKS OFF}
function IsCollinear(const pt1, pt2, pt3: TPoint64): Boolean;
function ProductsAreEqual(a, b, c, d: Int64): Boolean;
var
a,b: Int64;
absA,absB,absC,absD: UInt64;
absAB, absCD : UInt64;
carryAB, carryCD : UInt64;
signAB, signCD : integer;
begin
a := (pt2.X - pt1.X) * (pt3.Y - pt2.Y);
b := (pt2.Y - pt1.Y) * (pt3.X - pt2.X);
result := a = b;
// nb: unsigned values will be needed for CalcOverflowCarry()
absA := UInt64(Abs(a));
absB := UInt64(Abs(b));
absC := UInt64(Abs(c));
absD := UInt64(Abs(d));

// the multiplications here can potentially overflow, but
// any overflows will be compared using CalcOverflowCarry()
absAB := absA * absB;
absCD := absC * absD;

// nb: it's important to differentiate 0 values here from other values
signAB := TriSign(a) * TriSign(b);
signCD := TriSign(c) * TriSign(d);

carryAB := CalcOverflowCarry(absA, absB);
carryCD := CalcOverflowCarry(absC, absD);
Result := (absAB = absCD) and (signAB = signCD) and (carryAB = carryCD);
end;
{$OVERFLOWCHECKS ON}
//------------------------------------------------------------------------------

function IsCollinear(const pt1, sharedPt, pt2: TPoint64): Boolean;
var
a,b,c,d: Int64;
begin
a := sharedPt.X - pt1.X;
b := pt2.Y - sharedPt.Y;
c := sharedPt.Y - pt1.Y;
d := pt2.X - sharedPt.X;
// When checking for collinearity with very large coordinate values
// then ProductsAreEqual is more accurate than using CrossProduct.
Result := ProductsAreEqual(a, b, c, d);
end;
//------------------------------------------------------------------------------

function CrossProduct(const pt1, pt2, pt3: TPoint64): double;
begin
result := CrossProduct(
Expand Down

0 comments on commit 235879b

Please sign in to comment.