Skip to content

Commit

Permalink
Added CLIPPER2_HI_PRECISION option in CMakeLists.txt
Browse files Browse the repository at this point in the history
  • Loading branch information
AngusJohnson committed Dec 13, 2023
1 parent 326f127 commit 928379c
Show file tree
Hide file tree
Showing 5 changed files with 124 additions and 75 deletions.
83 changes: 39 additions & 44 deletions CPP/BenchMark/GetIntersectPtBenchmark.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -38,6 +38,40 @@ struct SetConsoleTextColor
//////////////////////////////////////////////////////////////////////////////////////


/////////////////////////////////////////////////////////
// Miscellaneous functions
/////////////////////////////////////////////////////////

double GetSineAbc(const Point64& a, const Point64& b, const Point64& c)
{
double dpB = DotProduct(a, b, c);
double SqrCosB = dpB * dpB / (DistanceSqr(a, b) * DistanceSqr(b, c));
double cos2B = SqrCosB * 2 - 1; // trig. itentity
return std::sqrt(1 - cos2B); // sin(B) = Sqrt(1-cos(2B))
}

static inline Point64 ReflectPoint(const Point64& pt, const Point64& pivot)
{
return Point64(pivot.x + (pivot.x - pt.x), pivot.y + (pivot.y - pt.y));
}

static inline Point64 MidPoint(const Point64& p1, const Point64& p2)
{
Point64 result;
result.x = (p1.x + p2.x) / 2;
result.y = (p1.y + p2.y) / 2;
return result;
}

static inline Point64 MakeRandomPoint(int64_t min_val, int64_t max_val)
{
std::random_device rd;
std::mt19937 gen(rd());
std::uniform_int_distribution<int64_t> x(min_val, max_val);
std::uniform_int_distribution<int64_t> y(min_val, max_val);
return Point64(x(gen), y(gen));
}

typedef std::function<bool(const Point64&, const Point64&,
const Point64&, const Point64&, Point64&)> GipFunction;

Expand Down Expand Up @@ -106,7 +140,9 @@ static bool GIP_Func_F(const Point64& ln1a, const Point64& ln1b,

/////////////////////////////////////////////////////////
// GIP_F_Mod: Modified GIP_Func_F (see above).
// Uses static_cast instead of nearbyint.
// Uses static_cast instead of nearbyint. Surprisingly,
// while faster here, this is much slower than GIP_Func_F
// when replacing GetIntersectPoint() in clipper.core.h.
/////////////////////////////////////////////////////////
#define CC_MIN(x,y) ((x)>(y)?(y):(x))
#define CC_MAX(x,y) ((x)<(y)?(y):(x))
Expand Down Expand Up @@ -160,46 +196,6 @@ typedef std::vector<TestRecord>::iterator test_iter;
// global data
std::vector<TestRecord> tests;

/////////////////////////////////////////////////////////
// Miscellaneous functions
/////////////////////////////////////////////////////////

double GetSine(const Point64& a, const Point64& b, const Point64& c)
{
Point64 ab = Point64(b.x - a.x, b.y - a.y);
Point64 bc = Point64(b.x - c.x, b.y - c.y);
double dp = DotProduct(ab,bc);
double distSqr1 = DistanceSqr(a, b);
double distSqr2 = DistanceSqr(b, c);
// cos(A)^2 = dp^2/DistSqr(ab)/DistSqr(bc)
double cosSqr = dp * dp / distSqr1 / distSqr2;
// cos2A = cos(A)^2 * 2 - 1
double cos2 = cosSqr * 2 - 1;
// sinA = Sqrt(1-cos2A)
return std::sqrt(1 - cos2);
}

static inline Point64 ReflectPoint(const Point64& pt, const Point64& pivot)
{
return Point64(pivot.x + (pivot.x - pt.x), pivot.y + (pivot.y - pt.y));
}

static inline Point64 MidPoint(const Point64& p1, const Point64& p2)
{
Point64 result;
result.x = (p1.x + p2.x) / 2;
result.y = (p1.y + p2.y) / 2;
return result;
}

static inline Point64 MakeRandomPoint(int64_t min_val, int64_t max_val)
{
std::random_device rd;
std::mt19937 gen(rd());
std::uniform_int_distribution<int64_t> x(min_val, max_val);
std::uniform_int_distribution<int64_t> y(min_val, max_val);
return Point64(x(gen), y(gen));
}

static inline GipFunction GetGipFunc(int index)
{
Expand Down Expand Up @@ -285,7 +281,6 @@ int main(int argc, char** argv)
BENCHMARK(BM_GIP_Current);
BENCHMARK(BM_GIP_Func_F);
BENCHMARK(BM_GIP_F_Mod);
benchmark::Counter();

bool first_pass = true;
for (int current_pow10 = 12; current_pow10 <= 18; ++current_pow10)
Expand All @@ -294,7 +289,7 @@ int main(int argc, char** argv)
// at their midpoints, while using random coordinates that are
// restricted to the specified power of 10 range
int64_t max_coord = static_cast<int64_t>(pow(10, current_pow10));
for (int64_t i = 0; i < 10000; ++i)
for (int64_t i = 0; i < 100000; ++i)
{
Point64 ip1 = MakeRandomPoint(-max_coord, max_coord);
Point64 ip2 = MakeRandomPoint(-max_coord, max_coord);
Expand All @@ -306,7 +301,7 @@ int main(int argc, char** argv)
// the closer segments are to collinear, the less accurate are
// calculations that determine intersection points.
// So exclude segments that are **almost** collinear. eg. sin(1deg) ~= 0.017
if (std::abs(GetSine(ip1, actual, ip3)) < 0.017) continue;
if (std::abs(GetSineAbc(ip1, actual, ip3)) < 0.017) continue;
// alternatively, only exclude segments that are collinear
//if (!CrossProduct(ip1, actual, ip3)) continue;
tests.push_back(TestRecord(participants, actual, ip1, ip2, ip3, ip4));
Expand Down
26 changes: 18 additions & 8 deletions CPP/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -9,21 +9,27 @@ set(CMAKE_CXX_STANDARD_REQUIRED ON)
set(CMAKE_CXX_EXTENSIONS OFF)
set_property(GLOBAL PROPERTY USE_FOLDERS ON)

# CLIPPER2_HI_PRECISION: See GetIntersectPoint() in clipper.core.h
option(CLIPPER2_HI_PRECISION "Caution: enabling this will compromise performance" OFF)

option(CLIPPER2_UTILS "Build utilities" ON)
option(CLIPPER2_EXAMPLES "Build examples" ON)
option(CLIPPER2_TESTS "Build tests" ON)
option(USE_EXTERNAL_GTEST "Use system-wide installed GoogleTest" OFF)
option(USE_EXTERNAL_GBENCHMARK "Use the googlebenchmark" OFF)
option(USE_EXTERNAL_GBENCHMARK "Use the googlebenchmark" OFF)
option(BUILD_SHARED_LIBS "Build shared libs" OFF)
set(CLIPPER2_USINGZ "ON" CACHE STRING "Build Clipper2Z, either \"ON\" or \"OFF\" or \"ONLY\"")
set(CLIPPER2_MAX_PRECISION 8 CACHE STRING "Maximum precision allowed for double to int64 scaling")

# CLIPPER2_MAX_DECIMAL_PRECISION: maximum decimal precision when scaling PathsD to Paths64.
# Caution: excessive scaling will increase the likelihood of integer overflow errors.
set(CLIPPER2_MAX_DECIMAL_PRECISION 8 CACHE STRING "Maximum decimal precision range")


if (APPLE)
set(CMAKE_SHARED_LIBRARY_SUFFIX ".dylib")
endif ()

include(GNUInstallDirs)

set(CLIPPER2_INC_FOLDER ${PROJECT_SOURCE_DIR}/Clipper2Lib/include/clipper2)
configure_file(clipper.version.in ${CLIPPER2_INC_FOLDER}/clipper.version.h)

Expand Down Expand Up @@ -52,10 +58,12 @@ if (NOT (CLIPPER2_USINGZ STREQUAL "ONLY"))
add_library(Clipper2 ${CLIPPER2_INC} ${CLIPPER2_SRC})

target_compile_definitions(
Clipper2
PUBLIC CLIPPER2_MAX_PRECISION=${CLIPPER2_MAX_PRECISION}
Clipper2 PUBLIC
CLIPPER2_MAX_DECIMAL_PRECISION=${CLIPPER2_MAX_DECIMAL_PRECISION}
$<$<BOOL:${CLIPPER2_HI_PRECISION}>:CLIPPER2_HI_PRECISION>
)


target_include_directories(Clipper2
PUBLIC Clipper2Lib/include
)
Expand All @@ -74,9 +82,11 @@ if (NOT (CLIPPER2_USINGZ STREQUAL "OFF"))
add_library(Clipper2Z ${CLIPPER2_INC} ${CLIPPER2_SRC})

target_compile_definitions(
Clipper2Z
PUBLIC USINGZ CLIPPER2_MAX_PRECISION=${CLIPPER2_MAX_PRECISION}
)
Clipper2Z PUBLIC
USINGZ
CLIPPER2_MAX_DECIMAL_PRECISION=${CLIPPER2_MAX_DECIMAL_PRECISION}
$<$<BOOL:${CLIPPER2_HI_PRECISION}>:CLIPPER2_HI_PRECISION>
)

target_include_directories(Clipper2Z
PUBLIC Clipper2Lib/include
Expand Down
70 changes: 57 additions & 13 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 : 9 December 2023 *
* Date : 13 December 2023 *
* Website : http://www.angusj.com *
* Copyright : Angus Johnson 2010-2023 *
* Purpose : Core Clipper Library structures and functions *
Expand Down Expand Up @@ -58,10 +58,10 @@ namespace Clipper2Lib
static const double PI = 3.141592653589793238;
#endif

#ifdef CLIPPER2_MAX_PRECISION
const int MAX_DECIMAL_PRECISION = CLIPPER2_MAX_PRECISION;
#ifdef CLIPPER2_MAX_DECIMAL_PRECISION
const int CLIPPER2_MAX_DEC_PRECISION = CLIPPER2_MAX_DECIMAL_PRECISION;
#else
const int MAX_DECIMAL_PRECISION = 8; // see Discussions #564
const int CLIPPER2_MAX_DEC_PRECISION = 8; // see Discussions #564
#endif

static const int64_t MAX_COORD = INT64_MAX >> 2;
Expand Down Expand Up @@ -242,6 +242,14 @@ namespace Clipper2Lib
(std::numeric_limits<double>::max)(),
(std::numeric_limits<double>::max)());

template<typename T>
static inline Point<T> MidPoint(const Point<T>& p1, const Point<T>& p2)
{
Point<T> result;
result.x = (p1.x + p2.x) / 2;
result.y = (p1.y + p2.y) / 2;
return result;
}

// Rect ------------------------------------------------------------------------

Expand Down Expand Up @@ -623,18 +631,19 @@ namespace Clipper2Lib

// Miscellaneous ------------------------------------------------------------

inline void CheckPrecision(int& precision, int& error_code)
inline void CheckPrecisionRange(int& precision, int& error_code)
{
if (precision >= -MAX_DECIMAL_PRECISION && precision <= MAX_DECIMAL_PRECISION) return;
if (precision >= -CLIPPER2_MAX_DEC_PRECISION &&
precision <= CLIPPER2_MAX_DEC_PRECISION) return;
error_code |= precision_error_i; // non-fatal error
DoError(precision_error_i); // does nothing unless exceptions enabled
precision = precision > 0 ? MAX_DECIMAL_PRECISION : -MAX_DECIMAL_PRECISION;
DoError(precision_error_i); // does nothing when exceptions are disabled
precision = precision > 0 ? CLIPPER2_MAX_DEC_PRECISION : -CLIPPER2_MAX_DEC_PRECISION;
}

inline void CheckPrecision(int& precision)
inline void CheckPrecisionRange(int& precision)
{
int error_code = 0;
CheckPrecision(precision, error_code);
CheckPrecisionRange(precision, error_code);
}

template <typename T>
Expand Down Expand Up @@ -721,6 +730,40 @@ namespace Clipper2Lib
return Area<T>(poly) >= 0;
}

#if CLIPPER2_HI_PRECISION

This comment has been minimized.

Copy link
@olologin

olologin Dec 15, 2023

@AngusJohnson
Do you have an approximate slowdown factor? Might be helpful to put somewhere in code comments.

This comment has been minimized.

Copy link
@AngusJohnson

AngusJohnson Dec 15, 2023

Author Owner

https://github.com/AngusJohnson/Clipper2/blob/main/CPP/BenchMark/GetIntersectPtBenchmark.cpp

However, for reasons I don't yet understand, the benchmarking above seems to grossly underestimate the performance slow down that I observe here:
https://github.com/AngusJohnson/Clipper2/blob/main/CPP/Examples/Benchmarks/Benchmarks.cpp
And the relative performace hit depends very much on the number of intersections in the clipping op.

// caution: this will compromise performance
// https://github.com/AngusJohnson/Clipper2/issues/317#issuecomment-1314023253
// See also CPP/BenchMark/GetIntersectPtBenchmark.cpp
#define CC_MIN(x,y) ((x)>(y)?(y):(x))
#define CC_MAX(x,y) ((x)<(y)?(y):(x))
inline bool GetIntersectPoint(const Point64& ln1a, const Point64& ln1b,
const Point64& ln2a, const Point64& ln2b, Point64& ip)
{
double ln1dy = (double)(ln1b.y - ln1a.y);
double ln1dx = (double)(ln1a.x - ln1b.x);
double ln2dy = (double)(ln2b.y - ln2a.y);
double ln2dx = (double)(ln2a.x - ln2b.x);
double det = (ln2dy * ln1dx) - (ln1dy * ln2dx);
if (det == 0.0) return false;
int64_t bb0minx = CC_MIN(ln1a.x, ln1b.x);
int64_t bb0miny = CC_MIN(ln1a.y, ln1b.y);
int64_t bb0maxx = CC_MAX(ln1a.x, ln1b.x);
int64_t bb0maxy = CC_MAX(ln1a.y, ln1b.y);
int64_t bb1minx = CC_MIN(ln2a.x, ln2b.x);
int64_t bb1miny = CC_MIN(ln2a.y, ln2b.y);
int64_t bb1maxx = CC_MAX(ln2a.x, ln2b.x);
int64_t bb1maxy = CC_MAX(ln2a.y, ln2b.y);
int64_t originx = (CC_MIN(bb0maxx, bb1maxx) + CC_MAX(bb0minx, bb1minx)) >> 1;
int64_t originy = (CC_MIN(bb0maxy, bb1maxy) + CC_MAX(bb0miny, bb1miny)) >> 1;
double ln0c = (ln1dy * (double)(ln1a.x - originx)) + (ln1dx * (double)(ln1a.y - originy));
double ln1c = (ln2dy * (double)(ln2a.x - originx)) + (ln2dx * (double)(ln2a.y - originy));
double hitx = ((ln1dx * ln1c) - (ln2dx * ln0c)) / det;
double hity = ((ln2dy * ln0c) - (ln1dy * ln1c)) / det;
ip.x = originx + (int64_t)nearbyint(hitx);
ip.y = originy + (int64_t)nearbyint(hity);
return true;
}
#else
inline bool GetIntersectPoint(const Point64& ln1a, const Point64& ln1b,
const Point64& ln2a, const Point64& ln2b, Point64& ip)
{
Expand All @@ -733,15 +776,16 @@ namespace Clipper2Lib
double det = dy1 * dx2 - dy2 * dx1;
if (det == 0.0) return false;
double t = ((ln1a.x - ln2a.x) * dy2 - (ln1a.y - ln2a.y) * dx2) / det;
if (t <= 0.0) ip = ln1a; // ?? check further (see also #568)
else if (t >= 1.0) ip = ln1b; // ?? check further
if (t <= 0.0) ip = ln1a;
else if (t >= 1.0) ip = ln1b;
else
{
ip.x = static_cast<int64_t>(ln1a.x + t * dx1);
ip.y = static_cast<int64_t>(ln1a.y + t * dy1);
}
}
return true;
}
#endif

inline bool SegmentsIntersect(const Point64& seg1a, const Point64& seg1b,
const Point64& seg2a, const Point64& seg2b, bool inclusive = false)
Expand Down
4 changes: 2 additions & 2 deletions CPP/Clipper2Lib/include/clipper2/clipper.engine.h
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
/*******************************************************************************
* Author : Angus Johnson *
* Date : 22 November 2023 *
* Date : 13 December 2023 *
* Website : http://www.angusj.com *
* Copyright : Angus Johnson 2010-2023 *
* Purpose : This is the main polygon clipping module *
Expand Down Expand Up @@ -530,7 +530,7 @@ namespace Clipper2Lib {
public:
explicit ClipperD(int precision = 2) : ClipperBase()
{
CheckPrecision(precision, error_code_);
CheckPrecisionRange(precision, error_code_);
// to optimize scaling / descaling precision
// set the scale to a power of double's radix (2) (#25)
scale_ = std::pow(std::numeric_limits<double>::radix,
Expand Down
16 changes: 8 additions & 8 deletions CPP/Clipper2Lib/include/clipper2/clipper.h
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
/*******************************************************************************
* Author : Angus Johnson *
* Date : 18 November 2023 *
* Date : 13 December 2023 *
* Website : http://www.angusj.com *
* Copyright : Angus Johnson 2010-2023 *
* Purpose : This module provides a simple interface to the Clipper Library *
Expand Down Expand Up @@ -47,7 +47,7 @@ namespace Clipper2Lib {
const PathsD& subjects, const PathsD& clips, int precision = 2)
{
int error_code = 0;
CheckPrecision(precision, error_code);
CheckPrecisionRange(precision, error_code);
PathsD result;
if (error_code) return result;
ClipperD clipper(precision);
Expand All @@ -63,7 +63,7 @@ namespace Clipper2Lib {
{
polytree.Clear();
int error_code = 0;
CheckPrecision(precision, error_code);
CheckPrecisionRange(precision, error_code);
if (error_code) return;
ClipperD clipper(precision);
clipper.AddSubject(subjects);
Expand Down Expand Up @@ -104,7 +104,7 @@ namespace Clipper2Lib {
{
PathsD result;
int error_code = 0;
CheckPrecision(precision, error_code);
CheckPrecisionRange(precision, error_code);
if (error_code) return result;
ClipperD clipper(precision);
clipper.AddSubject(subjects);
Expand Down Expand Up @@ -149,7 +149,7 @@ namespace Clipper2Lib {
int precision = 2, double arc_tolerance = 0.0)
{
int error_code = 0;
CheckPrecision(precision, error_code);
CheckPrecisionRange(precision, error_code);
if (!delta) return paths;
if (error_code) return PathsD();
const double scale = std::pow(10, precision);
Expand Down Expand Up @@ -219,7 +219,7 @@ namespace Clipper2Lib {
{
if (rect.IsEmpty() || paths.empty()) return PathsD();
int error_code = 0;
CheckPrecision(precision, error_code);
CheckPrecisionRange(precision, error_code);
if (error_code) return PathsD();
const double scale = std::pow(10, precision);
Rect64 r = ScaleRect<int64_t, double>(rect, scale);
Expand Down Expand Up @@ -251,7 +251,7 @@ namespace Clipper2Lib {
{
if (rect.IsEmpty() || lines.empty()) return PathsD();
int error_code = 0;
CheckPrecision(precision, error_code);
CheckPrecisionRange(precision, error_code);
if (error_code) return PathsD();
const double scale = std::pow(10, precision);
Rect64 r = ScaleRect<int64_t, double>(rect, scale);
Expand Down Expand Up @@ -545,7 +545,7 @@ namespace Clipper2Lib {
inline PathD TrimCollinear(const PathD& path, int precision, bool is_open_path = false)
{
int error_code = 0;
CheckPrecision(precision, error_code);
CheckPrecisionRange(precision, error_code);
if (error_code) return PathD();
const double scale = std::pow(10, precision);
Path64 p = ScalePath<int64_t, double>(path, scale, error_code);
Expand Down

0 comments on commit 928379c

Please sign in to comment.