Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Fix #10656 - PolygonClipping: CLIPLINE produces NaN when line to clip is parallel to the maxX of the clip window #10657

Merged
merged 4 commits into from
Aug 14, 2024
Merged

Fix #10656 - PolygonClipping: CLIPLINE produces NaN when line to clip is parallel to the maxX of the clip window #10657

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
21527bb
Expose CLIPLINE and write a simple unit test that shows the NaN occur…
jmarrec Aug 12, 2024
55edc57
Complete unit test for CLIPLINE
jmarrec Aug 12, 2024
251ea56
Fix #10656 - Rewrite the CLIPLINE algorithm
jmarrec Aug 12, 2024
955463a
Forgot to handle the rev case when I removed it as a param
jmarrec Aug 12, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
233 changes: 81 additions & 152 deletions src/EnergyPlus/SolarShading.cc
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -46,6 +46,7 @@
// POSSIBILITY OF SUCH DAMAGE.

// C++ Headers
#include <algorithm>
#include <cassert>
#include <cmath>
#include <memory>
Expand Down Expand Up @@ -3868,162 +3869,99 @@ inline bool d_eq(Real64 a, Real64 b)
return std::abs(a - b) < 2.0;
}

void CLIPLINE(Real64 &x1, Real64 &x2, Real64 &y1, Real64 &y2, Real64 maxX, Real64 minX, Real64 maxY, Real64 minY, bool &visible, bool &rev)
void CLIPLINE(Real64 &x0, Real64 &x1, Real64 &y0, Real64 &y1, Real64 maxX, Real64 minX, Real64 maxY, Real64 minY, bool &visible)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is fantastic. So much cleaner!

{

// Line segment clipping
// Reference:
// Slater, M., Barsky, B.A.
// Liang, Y.D., Barsky, B.A., Slater, M.
// 2D line and polygon clipping based on space subdivision.
// The Visual Computer 10, 407–422 (1994).
Real64 dx, dy, e, xinc, yinc, tempVar;
bool needX = true, needY = true;
int c1, c2;

if (x1 > x2) { // reverse for efficiency
tempVar = x1;
x1 = x2;
x2 = tempVar;
tempVar = y1;
y1 = y2;
y2 = tempVar;

// Tweaked via microbenchmarking to improve efficiency

bool rev = false;
if (x0 > x1) { // reverse for efficiency
std::swap(x0, x1);
std::swap(y0, y1);
rev = true;
}
if (x1 > maxX || x2 < minX) return; // x is positive
if (x1 < minX) {
if (y1 < minY) {
if (y2 < minY) return;
c1 = 0;
dx = x2 - x1;
dy = y2 - y1;
e = dy * (minX - x1) + dx * (y1 - minY);
} else if (y1 > maxY) {
if (y2 > maxY) return;
c1 = 6;
dx = x2 - x1;
dy = y2 - y1;
e = dy * (minX - x1) + dx * (y1 - maxY);
} else {
c1 = 3;
dx = x2 - x1;
dy = y2 - y1;
if (dy > 0) {
e = dy * (minX - x1) + dx * (y1 - maxY);
} else {
e = dy * (minX - x1) + dx * (y1 - minY);
}

if (x0 > maxX || x1 < minX) {
// Both points are outside the clip window, so they can't cross it
return;
}

// defining variables
Real64 const dx = x1 - x0; // >= 0
Real64 const dy = y1 - y0;

Real64 const q1 = x0 - minX;
Real64 const q2 = maxX - x0;
Real64 const q3 = y0 - minY;
Real64 const q4 = maxY - y0;

Real64 u1 = 0;
Real64 u2 = 1;

if ((dx == 0 && (q1 < 0 || q2 < 0)) || (dy == 0 && (q3 < 0 || q4 < 0))) {
// Line is parallel to clipping window
return;
}
if (dx != 0) {
Real64 const r1 = q1 / -dx;
if (r1 > u1) {
u1 = r1;
}
} else {
if (y1 < minY) {
if (y2 < minY) return;
c1 = 1;
dx = x2 - x1;
dy = y2 - y1;
e = dy * (maxX - x1) + dx * (y1 - minY);
} else if (y1 > maxY) {
if (y2 > maxY) return;
c1 = 7;
dx = x2 - x1;
dy = y2 - y1;
e = dy * (maxX - x1) + dx * (y1 - maxY);
} else {
visible = true;
if (x2 <= maxX && (y2 >= minY && y2 <= maxY)) return;
c1 = 4;
dx = x2 - x1;
dy = y2 - y1;
if (dy > 0) {
e = dy * (maxX - x1) + dx * (y1 - maxY);
} else {
e = dy * (maxX - x1) + dx * (y1 - minY);
}
Real64 const r2 = q2 / dx;
if (r2 < u2) {
u2 = r2;
}
}
c2 = c1;
if (dy > 0) {
while (true) {
if (e < 0.0) {
if (c2 == 1)
return;
else if (c2 == 3) {
visible = true;
x1 = minX;
y1 = maxY + e / dx;
if (x2 <= maxX && y2 <= maxY) return;
} else if (c2 == 4) {
x2 = maxX;
y2 = maxY + e / dx;
return;
}
if (needX) {
xinc = dy * (maxX - minX);
needX = false;
}
e += xinc;
c2 += 1;
} else {
if (c2 == 3)
return;
else if (c2 == 1) {
visible = true;
x1 = maxX - e / dy;
y1 = minY;
if (x2 <= maxX && y2 <= maxY) return;
} else if (c2 == 4) {
x2 = maxX - e / dy;
y2 = maxY;
return;
}
if (needY) {
yinc = dx * (maxY - minY);
needY = false;
}
e -= yinc;
c2 += 3;
if (dy != 0) {
Real64 const r3 = q3 / -dy;
Real64 const r4 = q4 / dy;
if (dy > 0) {
if (r3 > u1) {
u1 = r3;
}
}
} else {
while (true) {
if (e >= 0.0) {
if (c2 == 7)
return;
else if (c2 == 3) {
visible = true;
x1 = minX;
y1 = minY + e / dx;
if (x2 <= maxX && y2 >= minY) return;
} else if (c2 == 4) {
x2 = maxX;
y2 = minY + e / dx;
return;
}
if (needX) {
xinc = dy * (maxX - minX);
needX = false;
}
e += xinc;
c2 += 1;
} else {
if (c2 == 3)
return;
else if (c2 == 7) {
visible = true;
x1 = maxX - e / dy;
y1 = maxY;
if (x2 <= maxX && y2 >= minY) return;
} else if (c2 == 4) {
x2 = maxX - e / dy;
y2 = minY;
return;
}
if (needY) {
yinc = dx * (maxY - minY);
needY = false;
}
e += yinc;
c2 -= 3;
if (r4 < u2) {
u2 = r4;
}
} else {
if (r4 > u1) {
u1 = r4;
}
if (r3 < u2) {
u2 = r3;
}
}
}

if (u1 > u2) { // reject
// Line is outside the clipping window
return;
}

visible = true;

Real64 const xn0 = x0 + dx * u1;
Real64 const yn0 = y0 + dy * u1;

Real64 const xn1 = x0 + dx * u2;
Real64 const yn1 = y0 + dy * u2;

if (rev) {
x0 = xn1;
y0 = yn1;
x1 = xn0;
y1 = yn0;
} else {
x0 = xn0;
y0 = yn0;
x1 = xn1;
y1 = yn1;
}
}

void CLIPRECT(EnergyPlusData &state, int const NS2, int const NV1, int &NV3)
Expand Down Expand Up @@ -4063,20 +4001,11 @@ void CLIPRECT(EnergyPlusData &state, int const NS2, int const NV1, int &NV3)
Real64 x1 = x_1, x2 = x_2, y1 = y_1, y2 = y_2;

bool visible = false;
bool rev = false;
CLIPLINE(x_1, x_2, y_1, y_2, maxX, minX, maxY, minY, visible, rev);
CLIPLINE(x_1, x_2, y_1, y_2, maxX, minX, maxY, minY, visible);
if (visible) {
if ((x_1 != x1 || y_1 != y1) || (x_2 != x2 || y_2 != y2)) {
INTFLAG = true;
}
if (rev) { // undo reverse
Real64 tempVar = x_1;
x_1 = x_2;
x_2 = tempVar;
tempVar = y_1;
y_1 = y_2;
y_2 = tempVar;
}
// if line on edge, or inside, add both points
if (arrc == 0 || ((neq(arrx[arrc - 1], x_1) || neq(arry[arrc - 1], y_1)) && (neq(arrx[0], x_1) || neq(arry[0], y_1)))) {
arrx[arrc] = x_1;
Expand Down
2 changes: 2 additions & 0 deletions src/EnergyPlus/SolarShading.hh
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -136,6 +136,8 @@ namespace SolarShading {

void CLIP(EnergyPlusData &state, int const NVT, Array1D<Real64> &XVT, Array1D<Real64> &YVT, Array1D<Real64> &ZVT);

void CLIPLINE(Real64 &x0, Real64 &x1, Real64 &y0, Real64 &y1, Real64 maxX, Real64 minX, Real64 maxY, Real64 minY, bool &visible);
Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I expose it so I can test it. I'm also handling the rev inside the routine, makes no sense to leave the caller to deal with it


void CTRANS(EnergyPlusData &state,
int const NS, // Surface number whose vertex coordinates are being transformed
int const NGRS, // Base surface number for surface NS
Expand Down
Loading
Loading