| 16 | constexpr float NoAngleChangeLimit = 2 * PI_F; |
| 17 | |
| 18 | bool checkDeloneQuadrangle( const Vector3d& a, const Vector3d& b, const Vector3d& c, const Vector3d& d, double maxAngleChange ) |
| 19 | { |
| 20 | static constexpr double criticalDot = -0.9; |
| 21 | const auto nABC = normal( a, b, c ); |
| 22 | const auto nACD = normal( a, c, d ); |
| 23 | // true if the current triangles ABC and ACD are almost oppositely oriented |
| 24 | const bool oldPocket = dot( nABC, nACD ) < criticalDot; |
| 25 | |
| 26 | const auto nABD = normal( a, b, d ); |
| 27 | const auto nDBC = normal( d, b, c ); |
| 28 | // true if after flip the triangles ABD and DBC are almost oppositely oriented |
| 29 | const bool newPocket = dot( nABD, nDBC ) < criticalDot; |
| 30 | |
| 31 | if ( oldPocket != newPocket ) |
| 32 | return newPocket; // prefer the configuration without pocket |
| 33 | |
| 34 | // there should be significant difference in metrics (above floating point error) to return false |
| 35 | static constexpr double eps = 1e-7; // when we computed in floats then even 1e-5f was too small here and did not prevent infinite loop during resolveMeshDegenerations |
| 36 | if ( oldPocket ) |
| 37 | { |
| 38 | // before and after flip there are pockets, select the configuration with the smallest triangles |
| 39 | const auto metricAC = std::max( mincircleDiameterSq( a, c, d ), mincircleDiameterSq( c, a, b ) ); |
| 40 | const auto metricBD = std::max( mincircleDiameterSq( b, d, a ), mincircleDiameterSq( d, b, c ) ); |
| 41 | return metricAC <= metricBD + eps * ( metricAC + metricBD ); |
| 42 | } |
| 43 | |
| 44 | // before and after flip there are no pockets |
| 45 | |
| 46 | if ( maxAngleChange < NoAngleChangeLimit ) |
| 47 | { |
| 48 | const auto oldAngle = dihedralAngle( nABD, nDBC, d - b ); |
| 49 | const auto newAngle = dihedralAngle( nABC, nACD, a - c ); |
| 50 | const auto angleChange = std::abs( oldAngle - newAngle ); |
| 51 | if ( angleChange > maxAngleChange ) |
| 52 | return true; |
| 53 | } |
| 54 | |
| 55 | const auto metricAC = std::max( circumcircleDiameterSq( a, c, d ), circumcircleDiameterSq( c, a, b ) ); |
| 56 | const auto metricBD = std::max( circumcircleDiameterSq( b, d, a ), circumcircleDiameterSq( d, b, c ) ); |
| 57 | |
| 58 | if ( !std::isfinite( metricAC ) ) |
| 59 | { |
| 60 | if ( !std::isfinite( metricBD ) ) |
| 61 | { |
| 62 | // we are here if both configurations include a zero area obtuse triangle with all 3 vertices distinct; |
| 63 | // select the configuration with shorter diagonal |
| 64 | return distanceSq( a, c ) <= distanceSq( b, d ); |
| 65 | } |
| 66 | return metricAC <= metricBD; // (metricAC <= metricBD + eps * ( metricAC + metricBD ))==true if metricAC is +infinity and metricBD is finite |
| 67 | } |
| 68 | |
| 69 | return metricAC <= metricBD + eps * ( metricAC + metricBD ); // this shall work even if metricAC and metricBD are infinities, unlike ( metricAC - metricBD ), which becomes NaN |
| 70 | } |
| 71 | |
| 72 | bool checkDeloneQuadrangle( const Vector3f& a, const Vector3f& b, const Vector3f& c, const Vector3f& d, float maxAngleChange ) |
| 73 | { |