1596 lines
83 KiB
C++
1596 lines
83 KiB
C++
// clang-format off
|
|
|
|
#include "cafHexGridIntersectionTools.h"
|
|
|
|
#include "../LibCore/cvfPlane.h"
|
|
#include <cmath>
|
|
|
|
#include <algorithm>
|
|
#include <array>
|
|
|
|
namespace external
|
|
{
|
|
namespace caf
|
|
{
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
// HexGridIntersectionTools::ClipVx::ClipVx()
|
|
// : vx(cvf::Vec3d::ZERO),
|
|
// normDistFromEdgeVx1(HUGE_VAL),
|
|
// clippedEdgeVx1Id(-1),
|
|
// clippedEdgeVx2Id(-1),
|
|
// isVxIdsNative(true),
|
|
// derivedVxLevel(-1)
|
|
//{
|
|
//}
|
|
//
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
/// Find intersection between a line segment and a plane
|
|
///
|
|
/// \param a Start of line segment
|
|
/// \param b End of line segment
|
|
/// \param intersection Returns intersection point along the infinite line defined by a-b
|
|
/// \param normalizedDistFromA Returns the normalized (0..1) position from a to b of the intersection point.
|
|
/// Will return values along the infinite line defined by the a-b direcion,
|
|
/// and HUGE_VAL if plane and line are parallel.
|
|
/// \param epsilon Tolerance margin for accepting the position being within (0..1)
|
|
///
|
|
/// \return True if line segment intersects the plane
|
|
//--------------------------------------------------------------------------------------------------
|
|
bool HexGridIntersectionTools::planeLineIntersect( const cvf::Plane& plane,
|
|
const cvf::Vec3d& a,
|
|
const cvf::Vec3d& b,
|
|
cvf::Vec3d* intersection,
|
|
double* normalizedDistFromA,
|
|
double epsilon )
|
|
{
|
|
// From Real-Time Collision Detection by Christer Eriscon, published by Morgen Kaufmann Publishers, (c) 2005
|
|
// Elsevier Inc
|
|
|
|
cvf::Vec3d ab = b - a;
|
|
cvf::Vec3d normal = plane.normal();
|
|
|
|
double normDotAB = normal * ab;
|
|
if ( normDotAB == 0 )
|
|
{
|
|
( *normalizedDistFromA ) = HUGE_VAL;
|
|
return false;
|
|
}
|
|
|
|
double interpolationParameter = ( -plane.D() - ( normal * a ) ) / normDotAB;
|
|
|
|
( *intersection ) = a + interpolationParameter * ab;
|
|
( *normalizedDistFromA ) = interpolationParameter;
|
|
|
|
return ( interpolationParameter >= -epsilon && interpolationParameter <= 1.0 + epsilon );
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
/// Returns whether the triangle was hit by the plane.
|
|
/// isMostVxesOnPositiveSide returns true if all or two of the vxes is on the positive side of the plane.
|
|
/// newVx1/2.vx1ClippedEdge returns the index of the single vx that is alone on one side of the plane.
|
|
/// Going newVx1 to newVx2 will make the top triangle same winding as the original triangle,
|
|
/// and the quad opposite winding
|
|
|
|
// The permutations except for the trivial cases where all vertices are in front or behind plane:
|
|
//
|
|
//
|
|
// 1. Single vertex on positive side of plane => isMostVxesOnPositiveSide = false
|
|
//
|
|
// +\ /\3 /\3 /+ /\3 .
|
|
// \ / \ / \ / + / \ + .
|
|
// \2 \ / \/1 __1/____\2__ .
|
|
// / \ \ / /\ / \ .
|
|
// 1/___\1___\2 1/____2/__\2 1/________\2 .
|
|
// +\ /+
|
|
//
|
|
//
|
|
// 2. Two vertices vertex on positive side of plane => isMostVxesOnPositiveSide = true
|
|
//
|
|
// \+ /\3 /\3 +/ /\3 .
|
|
// \ / \ / \ / / \ .
|
|
// \2 \ / \/1 __1/____\2__ .
|
|
// / \ \ / /\ + / \ + .
|
|
// 1/___\1___\2 1/____2/__\2 1/________\2 .
|
|
// \+ +/
|
|
//
|
|
// 3. The special cases of touching one vertex, either exactly or "close enough"
|
|
// in finite precision. These occur for both 2. and 3 and in any rotation.
|
|
//
|
|
// a) Should not be counted b) May need a tolerance margin to intersect
|
|
// as intersecting: both 1->3 and 2->3 as it is theoretically required to:
|
|
// 3
|
|
// \ /\ /|\ .
|
|
// \ / \ / | \ .
|
|
// \ / \ / | \ .
|
|
// \ / \ / | \ .
|
|
// \/________\ /____|____\ .
|
|
// \ 1 | 2
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
|
|
bool HexGridIntersectionTools::planeTriangleIntersection( const cvf::Plane& plane,
|
|
const cvf::Vec3d& p1,
|
|
size_t p1Id,
|
|
const cvf::Vec3d& p2,
|
|
size_t p2Id,
|
|
const cvf::Vec3d& p3,
|
|
size_t p3Id,
|
|
ClipVx* newVx1,
|
|
ClipVx* newVx2,
|
|
bool* isMostVxesOnPositiveSide )
|
|
{
|
|
const double nonDimensionalTolerance = 1.0e-8;
|
|
|
|
double sqrSignedDistances[3];
|
|
sqrSignedDistances[0] = plane.distanceSquared( p1 );
|
|
sqrSignedDistances[1] = plane.distanceSquared( p2 );
|
|
sqrSignedDistances[2] = plane.distanceSquared( p3 );
|
|
|
|
double maxSqrAbsDistance = std::max( std::abs( sqrSignedDistances[0] ),
|
|
std::max( std::abs( sqrSignedDistances[1] ), std::abs( sqrSignedDistances[2] ) ) );
|
|
|
|
const double sqrDistanceTolerance = nonDimensionalTolerance * maxSqrAbsDistance;
|
|
|
|
int onPosSide[3];
|
|
onPosSide[0] = sqrSignedDistances[0] >= 0;
|
|
onPosSide[1] = sqrSignedDistances[1] >= 0;
|
|
onPosSide[2] = sqrSignedDistances[2] >= 0;
|
|
|
|
const int numPositiveVertices = onPosSide[0] + onPosSide[1] + onPosSide[2];
|
|
|
|
// The entire triangle is on the negative side
|
|
// Clip everything
|
|
if ( numPositiveVertices == 0 )
|
|
{
|
|
( *isMostVxesOnPositiveSide ) = false;
|
|
return false;
|
|
}
|
|
|
|
// All triangle vertices are on the positive side
|
|
if ( numPositiveVertices == 3 )
|
|
{
|
|
( *isMostVxesOnPositiveSide ) = true;
|
|
return false;
|
|
}
|
|
|
|
( *isMostVxesOnPositiveSide ) = ( numPositiveVertices == 2 );
|
|
|
|
int topVx = 0;
|
|
if ( numPositiveVertices == 1 )
|
|
{
|
|
if ( onPosSide[0] ) topVx = 1;
|
|
if ( onPosSide[1] ) topVx = 2;
|
|
if ( onPosSide[2] ) topVx = 3;
|
|
|
|
// Case 3a: Two negative distances and the last is within tolerance of zero.
|
|
if ( topVx > 0 && sqrSignedDistances[topVx - 1] < sqrDistanceTolerance )
|
|
{
|
|
return false;
|
|
}
|
|
}
|
|
else if ( numPositiveVertices == 2 )
|
|
{
|
|
if ( !onPosSide[0] ) topVx = 1;
|
|
if ( !onPosSide[1] ) topVx = 2;
|
|
if ( !onPosSide[2] ) topVx = 3;
|
|
|
|
// Case 3a: Two positive distances and the last is within tolerance of zero.
|
|
if ( topVx > 0 && sqrSignedDistances[topVx - 1] > -sqrDistanceTolerance )
|
|
{
|
|
return false;
|
|
}
|
|
}
|
|
else
|
|
{
|
|
CVF_ASSERT( false );
|
|
}
|
|
|
|
[[maybe_unused]] bool ok1 = false;
|
|
[[maybe_unused]] bool ok2 = false;
|
|
|
|
if ( topVx == 1 )
|
|
{
|
|
ok1 = planeLineIntersect( plane, p1, p2, &( ( *newVx1 ).vx ), &( ( *newVx1 ).normDistFromEdgeVx1 ), nonDimensionalTolerance );
|
|
( *newVx1 ).clippedEdgeVx1Id = p1Id;
|
|
( *newVx1 ).clippedEdgeVx2Id = p2Id;
|
|
ok2 = planeLineIntersect( plane, p1, p3, &( ( *newVx2 ).vx ), &( ( *newVx2 ).normDistFromEdgeVx1 ), nonDimensionalTolerance );
|
|
( *newVx2 ).clippedEdgeVx1Id = p1Id;
|
|
( *newVx2 ).clippedEdgeVx2Id = p3Id;
|
|
}
|
|
else if ( topVx == 2 )
|
|
{
|
|
ok1 = planeLineIntersect( plane, p2, p3, &( ( *newVx1 ).vx ), &( ( *newVx1 ).normDistFromEdgeVx1 ), nonDimensionalTolerance );
|
|
( *newVx1 ).clippedEdgeVx1Id = p2Id;
|
|
( *newVx1 ).clippedEdgeVx2Id = p3Id;
|
|
ok2 = planeLineIntersect( plane, p2, p1, &( ( *newVx2 ).vx ), &( ( *newVx2 ).normDistFromEdgeVx1 ), nonDimensionalTolerance );
|
|
( *newVx2 ).clippedEdgeVx1Id = p2Id;
|
|
( *newVx2 ).clippedEdgeVx2Id = p1Id;
|
|
}
|
|
else if ( topVx == 3 )
|
|
{
|
|
ok1 = planeLineIntersect( plane, p3, p1, &( ( *newVx1 ).vx ), &( ( *newVx1 ).normDistFromEdgeVx1 ), nonDimensionalTolerance );
|
|
( *newVx1 ).clippedEdgeVx1Id = p3Id;
|
|
( *newVx1 ).clippedEdgeVx2Id = p1Id;
|
|
ok2 = planeLineIntersect( plane, p3, p2, &( ( *newVx2 ).vx ), &( ( *newVx2 ).normDistFromEdgeVx1 ), nonDimensionalTolerance );
|
|
( *newVx2 ).clippedEdgeVx1Id = p3Id;
|
|
( *newVx2 ).clippedEdgeVx2Id = p2Id;
|
|
}
|
|
else
|
|
{
|
|
CVF_ASSERT( false );
|
|
}
|
|
|
|
// CVF_TIGHT_ASSERT(ok1 && ok2);
|
|
|
|
return true;
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
//
|
|
//
|
|
// P2 P2 P2 P2
|
|
// Keep Keep Keep Keep
|
|
// None Top 3 Quad All
|
|
// | | + | |
|
|
// | | / \ | |
|
|
// | | | | / \ | | | |
|
|
// | | | | / \ | | | |
|
|
// | | | | / \| | | |
|
|
// | | | |/ 1+ | | |
|
|
// | | | +2 |\ | | |
|
|
// | | | /| | \ | | |
|
|
// | | | / | | \ | _ | |
|
|
// | | | / | | \| |\Dir | |
|
|
// | | |/ | | 1+ \ | |
|
|
// | | +2 | | |\ \ | |
|
|
// | | /| | | | \ | |
|
|
// | | / |1 |1 2| 2| \ | |
|
|
// | | +--+----+----------+----+---+ | |
|
|
// | |1 | | | | 2 | |
|
|
// P1 P1 P1 P1
|
|
// Keep Keep Keep Keep
|
|
// All Quad Top None
|
|
//
|
|
//
|
|
// Clips the supplied triangles into new triangles returned in clippedTriangleVxes.
|
|
// New vertices have set isVxIdsNative = false and their vxIds is indices into triangleVxes
|
|
// The cellFaceForEachTriangleEdge refer to the edge after the corresponding triangle vertex.
|
|
// This method will keep the faces provided, while added edges is marked with no face = 6
|
|
//--------------------------------------------------------------------------------------------------
|
|
|
|
void HexGridIntersectionTools::clipTrianglesBetweenTwoParallelPlanes( const std::vector<ClipVx>& triangleVxes,
|
|
const std::vector<int>& cellFaceForEachTriangleEdge,
|
|
const cvf::Plane& p1Plane,
|
|
const cvf::Plane& p2Plane,
|
|
std::vector<ClipVx>* clippedTriangleVxes,
|
|
std::vector<int>* cellFaceForEachClippedTriangleEdge )
|
|
{
|
|
#define HT_NO_FACE 6
|
|
|
|
size_t triangleCount = triangleVxes.size() / 3;
|
|
|
|
for ( size_t tIdx = 0; tIdx < triangleCount; ++tIdx )
|
|
{
|
|
size_t triVxIdx = tIdx * 3;
|
|
|
|
ClipVx newVx1OnP1;
|
|
newVx1OnP1.isVxIdsNative = false;
|
|
ClipVx newVx2OnP1;
|
|
newVx2OnP1.isVxIdsNative = false;
|
|
|
|
bool isMostVxesOnPositiveSideOfP1 = false;
|
|
|
|
bool isIntersectingP1 = planeTriangleIntersection( p1Plane,
|
|
triangleVxes[triVxIdx + 0].vx,
|
|
triVxIdx + 0,
|
|
triangleVxes[triVxIdx + 1].vx,
|
|
triVxIdx + 1,
|
|
triangleVxes[triVxIdx + 2].vx,
|
|
triVxIdx + 2,
|
|
&newVx1OnP1,
|
|
&newVx2OnP1,
|
|
&isMostVxesOnPositiveSideOfP1 );
|
|
|
|
if ( !isIntersectingP1 && !isMostVxesOnPositiveSideOfP1 )
|
|
{
|
|
continue; // Discard triangle
|
|
}
|
|
|
|
ClipVx newVx1OnP2;
|
|
newVx1OnP2.isVxIdsNative = false;
|
|
ClipVx newVx2OnP2;
|
|
newVx2OnP2.isVxIdsNative = false;
|
|
|
|
bool isMostVxesOnPositiveSideOfP2 = false;
|
|
|
|
bool isIntersectingP2 = planeTriangleIntersection( p2Plane,
|
|
triangleVxes[triVxIdx + 0].vx,
|
|
triVxIdx + 0,
|
|
triangleVxes[triVxIdx + 1].vx,
|
|
triVxIdx + 1,
|
|
triangleVxes[triVxIdx + 2].vx,
|
|
triVxIdx + 2,
|
|
&newVx1OnP2,
|
|
&newVx2OnP2,
|
|
&isMostVxesOnPositiveSideOfP2 );
|
|
|
|
if ( !isIntersectingP2 && !isMostVxesOnPositiveSideOfP2 )
|
|
{
|
|
continue; // Discard triangle
|
|
}
|
|
|
|
bool p1KeepAll = ( !isIntersectingP1 && isMostVxesOnPositiveSideOfP1 );
|
|
bool p2KeepAll = ( !isIntersectingP2 && isMostVxesOnPositiveSideOfP2 );
|
|
bool p1KeepQuad = ( isIntersectingP1 && isMostVxesOnPositiveSideOfP1 );
|
|
bool p2KeepQuad = ( isIntersectingP2 && isMostVxesOnPositiveSideOfP2 );
|
|
bool p1KeepTop = ( isIntersectingP1 && !isMostVxesOnPositiveSideOfP1 );
|
|
bool p2KeepTop = ( isIntersectingP2 && !isMostVxesOnPositiveSideOfP2 );
|
|
|
|
if ( p1KeepAll && p2KeepAll )
|
|
{
|
|
// Keep the triangle
|
|
clippedTriangleVxes->push_back( triangleVxes[triVxIdx + 0] );
|
|
clippedTriangleVxes->push_back( triangleVxes[triVxIdx + 1] );
|
|
clippedTriangleVxes->push_back( triangleVxes[triVxIdx + 2] );
|
|
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[triVxIdx + 0] );
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[triVxIdx + 1] );
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[triVxIdx + 2] );
|
|
|
|
continue;
|
|
}
|
|
|
|
if ( p1KeepQuad && p2KeepAll )
|
|
{
|
|
// Split the resulting quad and add the two triangles
|
|
clippedTriangleVxes->push_back( newVx2OnP1 );
|
|
clippedTriangleVxes->push_back( newVx1OnP1 );
|
|
clippedTriangleVxes->push_back( triangleVxes[newVx1OnP1.clippedEdgeVx2Id] );
|
|
|
|
clippedTriangleVxes->push_back( triangleVxes[newVx1OnP1.clippedEdgeVx2Id] );
|
|
clippedTriangleVxes->push_back( triangleVxes[newVx2OnP1.clippedEdgeVx2Id] );
|
|
clippedTriangleVxes->push_back( newVx2OnP1 );
|
|
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[newVx1OnP1.clippedEdgeVx1Id] );
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[newVx1OnP1.clippedEdgeVx2Id] );
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[newVx2OnP1.clippedEdgeVx2Id] );
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
|
|
continue;
|
|
}
|
|
|
|
if ( p2KeepQuad && p1KeepAll )
|
|
{
|
|
// Split the resulting quad and add the two triangles
|
|
clippedTriangleVxes->push_back( newVx2OnP2 );
|
|
clippedTriangleVxes->push_back( newVx1OnP2 );
|
|
clippedTriangleVxes->push_back( triangleVxes[newVx2OnP2.clippedEdgeVx2Id] );
|
|
|
|
clippedTriangleVxes->push_back( newVx1OnP2 );
|
|
clippedTriangleVxes->push_back( triangleVxes[newVx1OnP2.clippedEdgeVx2Id] );
|
|
clippedTriangleVxes->push_back( triangleVxes[newVx2OnP2.clippedEdgeVx2Id] );
|
|
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[newVx2OnP2.clippedEdgeVx2Id] );
|
|
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[newVx1OnP2.clippedEdgeVx1Id] );
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[newVx1OnP2.clippedEdgeVx2Id] );
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
|
|
continue;
|
|
}
|
|
|
|
if ( p1KeepTop && p2KeepAll )
|
|
{
|
|
// Add the top triangle
|
|
clippedTriangleVxes->push_back( newVx1OnP1 );
|
|
clippedTriangleVxes->push_back( newVx2OnP1 );
|
|
clippedTriangleVxes->push_back( triangleVxes[newVx1OnP1.clippedEdgeVx1Id] );
|
|
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[newVx2OnP1.clippedEdgeVx2Id] );
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[newVx1OnP1.clippedEdgeVx1Id] );
|
|
|
|
continue;
|
|
}
|
|
|
|
if ( p2KeepTop && p1KeepAll )
|
|
{
|
|
// Add the top triangle
|
|
clippedTriangleVxes->push_back( newVx1OnP2 );
|
|
clippedTriangleVxes->push_back( newVx2OnP2 );
|
|
clippedTriangleVxes->push_back( triangleVxes[newVx1OnP2.clippedEdgeVx1Id] );
|
|
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[newVx2OnP2.clippedEdgeVx2Id] );
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[newVx1OnP2.clippedEdgeVx1Id] );
|
|
|
|
continue;
|
|
}
|
|
|
|
if ( p1KeepQuad && p2KeepQuad )
|
|
{
|
|
// We end up with a pentagon.
|
|
clippedTriangleVxes->push_back( newVx2OnP1 );
|
|
clippedTriangleVxes->push_back( newVx1OnP1 );
|
|
clippedTriangleVxes->push_back( newVx2OnP2 );
|
|
|
|
clippedTriangleVxes->push_back( newVx2OnP2 );
|
|
clippedTriangleVxes->push_back( newVx1OnP2 );
|
|
clippedTriangleVxes->push_back( newVx2OnP1 );
|
|
|
|
// Two variants. The original point might be along newVx1OnP1 to newVx2OnP2 or along newVx2OnP1 to newVx1OnP2
|
|
if ( newVx1OnP1.clippedEdgeVx2Id == newVx2OnP2.clippedEdgeVx1Id )
|
|
{
|
|
clippedTriangleVxes->push_back( newVx2OnP1 );
|
|
clippedTriangleVxes->push_back( newVx1OnP2 );
|
|
clippedTriangleVxes->push_back( triangleVxes[newVx2OnP1.clippedEdgeVx2Id] );
|
|
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[newVx1OnP1.clippedEdgeVx1Id] );
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[newVx1OnP2.clippedEdgeVx1Id] );
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[newVx2OnP1.clippedEdgeVx2Id] );
|
|
}
|
|
else
|
|
{
|
|
clippedTriangleVxes->push_back( newVx2OnP2 );
|
|
clippedTriangleVxes->push_back( newVx1OnP1 );
|
|
clippedTriangleVxes->push_back( triangleVxes[newVx2OnP2.clippedEdgeVx2Id] );
|
|
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[newVx1OnP2.clippedEdgeVx1Id] );
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[newVx1OnP1.clippedEdgeVx1Id] );
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[newVx2OnP2.clippedEdgeVx2Id] );
|
|
}
|
|
|
|
continue;
|
|
}
|
|
|
|
if ( p1KeepQuad && p2KeepTop )
|
|
{
|
|
// We end up with a quad.
|
|
clippedTriangleVxes->push_back( newVx1OnP1 );
|
|
clippedTriangleVxes->push_back( newVx1OnP2 );
|
|
clippedTriangleVxes->push_back( newVx2OnP1 );
|
|
|
|
clippedTriangleVxes->push_back( newVx1OnP2 );
|
|
clippedTriangleVxes->push_back( newVx2OnP2 );
|
|
clippedTriangleVxes->push_back( newVx2OnP1 );
|
|
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[newVx1OnP1.clippedEdgeVx1Id] );
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[newVx2OnP2.clippedEdgeVx2Id] );
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
|
|
continue;
|
|
}
|
|
|
|
if ( p2KeepQuad && p1KeepTop )
|
|
{
|
|
// We end up with a quad.
|
|
clippedTriangleVxes->push_back( newVx2OnP1 );
|
|
clippedTriangleVxes->push_back( newVx2OnP2 );
|
|
clippedTriangleVxes->push_back( newVx1OnP2 );
|
|
|
|
clippedTriangleVxes->push_back( newVx2OnP1 );
|
|
clippedTriangleVxes->push_back( newVx1OnP2 );
|
|
clippedTriangleVxes->push_back( newVx1OnP1 );
|
|
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[newVx2OnP1.clippedEdgeVx2Id] );
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
cellFaceForEachClippedTriangleEdge->push_back( cellFaceForEachTriangleEdge[newVx1OnP2.clippedEdgeVx1Id] );
|
|
cellFaceForEachClippedTriangleEdge->push_back( HT_NO_FACE );
|
|
|
|
continue;
|
|
}
|
|
|
|
CVF_ASSERT( false );
|
|
}
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
// Creates a plane with normal perpendicular to the edge, pointing in the direction of the pointInNormalDirection
|
|
//--------------------------------------------------------------------------------------------------
|
|
cvf::Plane createPlaneFromEdgeAndPointInNormalDirection( cvf::Vec3d ep1, cvf::Vec3d ep2, cvf::Vec3d pointInNormalDirection )
|
|
{
|
|
cvf::Vec3d ep1ep2 = ep2 - ep1;
|
|
cvf::Vec3d ep1pointforNorm = pointInNormalDirection - ep1;
|
|
cvf::Vec3d triNormal = ep1ep2 ^ ep1pointforNorm;
|
|
cvf::Vec3d pointInPlane = ep1 + triNormal;
|
|
|
|
cvf::Plane plane;
|
|
plane.setFromPoints( ep1, pointInPlane, ep2 );
|
|
|
|
return plane;
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
// Clips the supplied triangles into new triangles returned in clippedTriangleVxes.
|
|
// New vertices have set isVxIdsNative = false and their vxIds is indices into triangleVxes
|
|
// The cellFaceForEachTriangleEdge refer to the edge after the corresponding triangle vertex.
|
|
// This method will keep the faces provided, while added edges is marked with no face = 6
|
|
//--------------------------------------------------------------------------------------------------
|
|
void HexGridIntersectionTools::clipPlanarTrianglesWithInPlaneTriangle( const std::vector<cvf::Vec3d>& triangleVxes,
|
|
const std::vector<int>& cellFaceForEachTriangleEdge,
|
|
const cvf::Vec3d& tp1,
|
|
const cvf::Vec3d& tp2,
|
|
const cvf::Vec3d& tp3,
|
|
std::vector<cvf::Vec3d>* clippedTriangleVxes,
|
|
std::vector<int>* cellFaceForEachClippedTriangleEdge )
|
|
{
|
|
#define HT_NO_FACE 6
|
|
|
|
size_t triangleCount = triangleVxes.size() / 3;
|
|
|
|
// Creating a plane for each of the edges of the clipping triangle
|
|
std::array<cvf::Plane, 3> clipTrianglePlanes;
|
|
|
|
clipTrianglePlanes[0] = createPlaneFromEdgeAndPointInNormalDirection( tp1, tp2, tp3 );
|
|
clipTrianglePlanes[1] = createPlaneFromEdgeAndPointInNormalDirection( tp2, tp3, tp1 );
|
|
clipTrianglePlanes[2] = createPlaneFromEdgeAndPointInNormalDirection( tp3, tp1, tp2 );
|
|
|
|
#define reserveSize 100
|
|
|
|
std::vector<cvf::Vec3d> currentInputTriangleVxes;
|
|
currentInputTriangleVxes.reserve( reserveSize );
|
|
std::vector<int> currentInputCellFaceForEachTriangleEdge;
|
|
currentInputCellFaceForEachTriangleEdge.reserve( reserveSize );
|
|
std::vector<cvf::Vec3d> currentOutputTriangleVxes;
|
|
currentOutputTriangleVxes.reserve( reserveSize );
|
|
std::vector<int> currentOutputCellFaceForEachTriangleEdge;
|
|
currentOutputCellFaceForEachTriangleEdge.reserve( reserveSize );
|
|
|
|
for ( size_t tIdx = 0; tIdx < triangleCount; ++tIdx )
|
|
{
|
|
size_t triVxIdx = tIdx * 3;
|
|
|
|
currentInputTriangleVxes.clear();
|
|
currentInputCellFaceForEachTriangleEdge.clear();
|
|
currentOutputTriangleVxes.clear();
|
|
currentOutputCellFaceForEachTriangleEdge.clear();
|
|
|
|
currentOutputTriangleVxes.push_back( triangleVxes[triVxIdx + 0] );
|
|
currentOutputTriangleVxes.push_back( triangleVxes[triVxIdx + 1] );
|
|
currentOutputTriangleVxes.push_back( triangleVxes[triVxIdx + 2] );
|
|
|
|
currentOutputCellFaceForEachTriangleEdge.push_back( cellFaceForEachTriangleEdge[triVxIdx + 0] );
|
|
currentOutputCellFaceForEachTriangleEdge.push_back( cellFaceForEachTriangleEdge[triVxIdx + 1] );
|
|
currentOutputCellFaceForEachTriangleEdge.push_back( cellFaceForEachTriangleEdge[triVxIdx + 2] );
|
|
|
|
ClipVx newVx1;
|
|
newVx1.isVxIdsNative = false;
|
|
ClipVx newVx2;
|
|
newVx2.isVxIdsNative = false;
|
|
|
|
for ( int planeIdx = 0; planeIdx < 3; ++planeIdx )
|
|
{
|
|
currentInputTriangleVxes.swap( currentOutputTriangleVxes );
|
|
currentInputCellFaceForEachTriangleEdge.swap( currentOutputCellFaceForEachTriangleEdge );
|
|
|
|
currentOutputTriangleVxes.clear();
|
|
currentOutputCellFaceForEachTriangleEdge.clear();
|
|
|
|
size_t inTriangleCount = currentInputTriangleVxes.size() / 3;
|
|
|
|
for ( size_t inTrIdx = 0; inTrIdx < inTriangleCount; ++inTrIdx )
|
|
{
|
|
size_t inTriVxIdx = inTrIdx * 3;
|
|
|
|
bool isMostVxesOnPositiveSide = false;
|
|
|
|
bool isIntersectingPlane = planeTriangleIntersection( clipTrianglePlanes[planeIdx],
|
|
currentInputTriangleVxes[inTriVxIdx + 0],
|
|
inTriVxIdx + 0,
|
|
currentInputTriangleVxes[inTriVxIdx + 1],
|
|
inTriVxIdx + 1,
|
|
currentInputTriangleVxes[inTriVxIdx + 2],
|
|
inTriVxIdx + 2,
|
|
&newVx1,
|
|
&newVx2,
|
|
&isMostVxesOnPositiveSide );
|
|
|
|
if ( !isIntersectingPlane )
|
|
{
|
|
// All on negative side: Discard triangle
|
|
|
|
if ( !isMostVxesOnPositiveSide )
|
|
{
|
|
continue;
|
|
}
|
|
else // All on positive side: keep all
|
|
{
|
|
currentOutputTriangleVxes.push_back( currentInputTriangleVxes[inTriVxIdx + 0] );
|
|
currentOutputTriangleVxes.push_back( currentInputTriangleVxes[inTriVxIdx + 1] );
|
|
currentOutputTriangleVxes.push_back( currentInputTriangleVxes[inTriVxIdx + 2] );
|
|
|
|
currentOutputCellFaceForEachTriangleEdge.push_back(
|
|
currentInputCellFaceForEachTriangleEdge[inTriVxIdx + 0] );
|
|
currentOutputCellFaceForEachTriangleEdge.push_back(
|
|
currentInputCellFaceForEachTriangleEdge[inTriVxIdx + 1] );
|
|
currentOutputCellFaceForEachTriangleEdge.push_back(
|
|
currentInputCellFaceForEachTriangleEdge[inTriVxIdx + 2] );
|
|
}
|
|
}
|
|
else // intersecting
|
|
{
|
|
if ( isMostVxesOnPositiveSide )
|
|
{
|
|
// We need the Quad
|
|
|
|
currentOutputTriangleVxes.push_back( newVx2.vx );
|
|
currentOutputTriangleVxes.push_back( newVx1.vx );
|
|
currentOutputTriangleVxes.push_back( currentInputTriangleVxes[newVx1.clippedEdgeVx2Id] );
|
|
|
|
currentOutputCellFaceForEachTriangleEdge.push_back( HT_NO_FACE );
|
|
currentOutputCellFaceForEachTriangleEdge.push_back(
|
|
currentInputCellFaceForEachTriangleEdge[newVx1.clippedEdgeVx1Id] );
|
|
currentOutputCellFaceForEachTriangleEdge.push_back( HT_NO_FACE );
|
|
|
|
currentOutputTriangleVxes.push_back( currentInputTriangleVxes[newVx1.clippedEdgeVx2Id] );
|
|
currentOutputTriangleVxes.push_back( currentInputTriangleVxes[newVx2.clippedEdgeVx2Id] );
|
|
currentOutputTriangleVxes.push_back( newVx2.vx );
|
|
|
|
currentOutputCellFaceForEachTriangleEdge.push_back(
|
|
currentInputCellFaceForEachTriangleEdge[newVx1.clippedEdgeVx2Id] );
|
|
currentOutputCellFaceForEachTriangleEdge.push_back(
|
|
currentInputCellFaceForEachTriangleEdge[newVx2.clippedEdgeVx2Id] );
|
|
currentOutputCellFaceForEachTriangleEdge.push_back( HT_NO_FACE );
|
|
}
|
|
else
|
|
{
|
|
currentOutputTriangleVxes.push_back( newVx1.vx );
|
|
currentOutputTriangleVxes.push_back( newVx2.vx );
|
|
currentOutputTriangleVxes.push_back( currentInputTriangleVxes[newVx1.clippedEdgeVx1Id] );
|
|
|
|
currentOutputCellFaceForEachTriangleEdge.push_back( HT_NO_FACE );
|
|
currentOutputCellFaceForEachTriangleEdge.push_back(
|
|
currentInputCellFaceForEachTriangleEdge[newVx2.clippedEdgeVx2Id] );
|
|
currentOutputCellFaceForEachTriangleEdge.push_back(
|
|
currentInputCellFaceForEachTriangleEdge[newVx1.clippedEdgeVx1Id] );
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
// Append the result of the completely clipped triangle to the output
|
|
|
|
clippedTriangleVxes->insert( clippedTriangleVxes->end(),
|
|
currentOutputTriangleVxes.begin(),
|
|
currentOutputTriangleVxes.end() );
|
|
cellFaceForEachClippedTriangleEdge->insert( cellFaceForEachClippedTriangleEdge->end(),
|
|
currentOutputCellFaceForEachTriangleEdge.begin(),
|
|
currentOutputCellFaceForEachTriangleEdge.end() );
|
|
}
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
/// Will return the intersection point. If the plane is outside the line, it returns the closest line endpoint
|
|
//--------------------------------------------------------------------------------------------------
|
|
cvf::Vec3d HexGridIntersectionTools::planeLineIntersectionForMC( const cvf::Plane& plane,
|
|
const cvf::Vec3d& p1,
|
|
const cvf::Vec3d& p2,
|
|
double* normalizedDistFromP1 )
|
|
{
|
|
// From http://local.wasp.uwa.edu.au/~pbourke/geometry/planeline/
|
|
//
|
|
// P1 (x1,y1,z1) and P2 (x2,y2,z2)
|
|
//
|
|
// P = P1 + u (P2 - P1)
|
|
//
|
|
// A*x1 + B*y1 + C*z1 + D
|
|
// u = ---------------------------------
|
|
// A*(x1-x2) + B*(y1-y2) + C*(z1-z2)
|
|
|
|
CVF_TIGHT_ASSERT( normalizedDistFromP1 );
|
|
|
|
const cvf::Vec3d v = p2 - p1;
|
|
|
|
( *normalizedDistFromP1 ) = 0.0;
|
|
|
|
double denominator = -( plane.A() * v.x() + plane.B() * v.y() + plane.C() * v.z() );
|
|
if ( denominator != 0 )
|
|
{
|
|
double u = ( plane.A() * p1.x() + plane.B() * p1.y() + plane.C() * p1.z() + plane.D() ) / denominator;
|
|
( *normalizedDistFromP1 ) = u;
|
|
if ( u > 0.0 && u < 1.0 )
|
|
{
|
|
return ( p1 + u * v );
|
|
}
|
|
else
|
|
{
|
|
if ( u >= 1.0 )
|
|
{
|
|
return p2;
|
|
}
|
|
else
|
|
{
|
|
return p1;
|
|
}
|
|
}
|
|
}
|
|
else
|
|
{
|
|
return p1;
|
|
}
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
/// Based on description and implementation from Paul Bourke:
|
|
///
|
|
/// http://paulbourke.net/geometry/polygonise/
|
|
///
|
|
/// Note that the element is turned inside-out compared to what we use elsewhere in caf/ResInsight
|
|
/// So the winding of all the sides are opposite.
|
|
/// 4-----4------5
|
|
/// /| /| k POS_I = 0
|
|
/// 7 8 5 9 | NEG_I = 1
|
|
/// / | / | | POS_J = 2
|
|
/// 7------6-----6 | | NEG_J = 3
|
|
/// | 0-----0--|---1 *------i POS_K = 4
|
|
/// 11 / 10 / / NEG_K = 5
|
|
/// | 3 | 1 / NO_FACE = 6
|
|
/// |/ |/ j
|
|
/// 3------2-----2
|
|
///
|
|
// The cellFaceForEachTriangleEdge refer to the edge after the corresponding triangle vertex.
|
|
//--------------------------------------------------------------------------------------------------
|
|
int HexGridIntersectionTools::planeHexIntersectionMC( const cvf::Plane& plane,
|
|
const cvf::Vec3d cell[8],
|
|
const size_t hexCornersIds[8],
|
|
std::vector<ClipVx>* triangleVxes,
|
|
std::vector<int>* cellFaceForEachTriangleEdge )
|
|
{
|
|
// clang-format off
|
|
|
|
static const cvf::uint cubeIdxToCutEdgeBitfield[256] =
|
|
{
|
|
0x0, 0x109, 0x203, 0x30a, 0x406, 0x50f, 0x605, 0x70c,
|
|
0x80c, 0x905, 0xa0f, 0xb06, 0xc0a, 0xd03, 0xe09, 0xf00,
|
|
0x190, 0x99, 0x393, 0x29a, 0x596, 0x49f, 0x795, 0x69c,
|
|
0x99c, 0x895, 0xb9f, 0xa96, 0xd9a, 0xc93, 0xf99, 0xe90,
|
|
0x230, 0x339, 0x33, 0x13a, 0x636, 0x73f, 0x435, 0x53c,
|
|
0xa3c, 0xb35, 0x83f, 0x936, 0xe3a, 0xf33, 0xc39, 0xd30,
|
|
0x3a0, 0x2a9, 0x1a3, 0xaa, 0x7a6, 0x6af, 0x5a5, 0x4ac,
|
|
0xbac, 0xaa5, 0x9af, 0x8a6, 0xfaa, 0xea3, 0xda9, 0xca0,
|
|
0x460, 0x569, 0x663, 0x76a, 0x66, 0x16f, 0x265, 0x36c,
|
|
0xc6c, 0xd65, 0xe6f, 0xf66, 0x86a, 0x963, 0xa69, 0xb60,
|
|
0x5f0, 0x4f9, 0x7f3, 0x6fa, 0x1f6, 0xff, 0x3f5, 0x2fc,
|
|
0xdfc, 0xcf5, 0xfff, 0xef6, 0x9fa, 0x8f3, 0xbf9, 0xaf0,
|
|
0x650, 0x759, 0x453, 0x55a, 0x256, 0x35f, 0x55, 0x15c,
|
|
0xe5c, 0xf55, 0xc5f, 0xd56, 0xa5a, 0xb53, 0x859, 0x950,
|
|
0x7c0, 0x6c9, 0x5c3, 0x4ca, 0x3c6, 0x2cf, 0x1c5, 0xcc,
|
|
0xfcc, 0xec5, 0xdcf, 0xcc6, 0xbca, 0xac3, 0x9c9, 0x8c0,
|
|
0x8c0, 0x9c9, 0xac3, 0xbca, 0xcc6, 0xdcf, 0xec5, 0xfcc,
|
|
0xcc, 0x1c5, 0x2cf, 0x3c6, 0x4ca, 0x5c3, 0x6c9, 0x7c0,
|
|
0x950, 0x859, 0xb53, 0xa5a, 0xd56, 0xc5f, 0xf55, 0xe5c,
|
|
0x15c, 0x55, 0x35f, 0x256, 0x55a, 0x453, 0x759, 0x650,
|
|
0xaf0, 0xbf9, 0x8f3, 0x9fa, 0xef6, 0xfff, 0xcf5, 0xdfc,
|
|
0x2fc, 0x3f5, 0xff, 0x1f6, 0x6fa, 0x7f3, 0x4f9, 0x5f0,
|
|
0xb60, 0xa69, 0x963, 0x86a, 0xf66, 0xe6f, 0xd65, 0xc6c,
|
|
0x36c, 0x265, 0x16f, 0x66, 0x76a, 0x663, 0x569, 0x460,
|
|
0xca0, 0xda9, 0xea3, 0xfaa, 0x8a6, 0x9af, 0xaa5, 0xbac,
|
|
0x4ac, 0x5a5, 0x6af, 0x7a6, 0xaa, 0x1a3, 0x2a9, 0x3a0,
|
|
0xd30, 0xc39, 0xf33, 0xe3a, 0x936, 0x83f, 0xb35, 0xa3c,
|
|
0x53c, 0x435, 0x73f, 0x636, 0x13a, 0x33, 0x339, 0x230,
|
|
0xe90, 0xf99, 0xc93, 0xd9a, 0xa96, 0xb9f, 0x895, 0x99c,
|
|
0x69c, 0x795, 0x49f, 0x596, 0x29a, 0x393, 0x99, 0x190,
|
|
0xf00, 0xe09, 0xd03, 0xc0a, 0xb06, 0xa0f, 0x905, 0x80c,
|
|
0x70c, 0x605, 0x50f, 0x406, 0x30a, 0x203, 0x109, 0x0
|
|
};
|
|
|
|
// clang-format on
|
|
|
|
static const int cubeIdxToTriangleIndices[256][16] =
|
|
{ { -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 8, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 1, 9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 8, 3, 9, 8, 1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 2, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 8, 3, 1, 2, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 9, 2, 10, 0, 2, 9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 2, 8, 3, 2, 10, 8, 10, 9, 8, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 3, 11, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 11, 2, 8, 11, 0, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 9, 0, 2, 3, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 11, 2, 1, 9, 11, 9, 8, 11, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 3, 10, 1, 11, 10, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 10, 1, 0, 8, 10, 8, 11, 10, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 3, 9, 0, 3, 11, 9, 11, 10, 9, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 9, 8, 10, 10, 8, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 4, 7, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 4, 3, 0, 7, 3, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 1, 9, 8, 4, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 4, 1, 9, 4, 7, 1, 7, 3, 1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 2, 10, 8, 4, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 3, 4, 7, 3, 0, 4, 1, 2, 10, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 9, 2, 10, 9, 0, 2, 8, 4, 7, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 2, 10, 9, 2, 9, 7, 2, 7, 3, 7, 9, 4, -1, -1, -1, -1 },
|
|
{ 8, 4, 7, 3, 11, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 11, 4, 7, 11, 2, 4, 2, 0, 4, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 9, 0, 1, 8, 4, 7, 2, 3, 11, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 4, 7, 11, 9, 4, 11, 9, 11, 2, 9, 2, 1, -1, -1, -1, -1 },
|
|
{ 3, 10, 1, 3, 11, 10, 7, 8, 4, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 11, 10, 1, 4, 11, 1, 0, 4, 7, 11, 4, -1, -1, -1, -1 },
|
|
{ 4, 7, 8, 9, 0, 11, 9, 11, 10, 11, 0, 3, -1, -1, -1, -1 },
|
|
{ 4, 7, 11, 4, 11, 9, 9, 11, 10, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 9, 5, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 9, 5, 4, 0, 8, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 5, 4, 1, 5, 0, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 8, 5, 4, 8, 3, 5, 3, 1, 5, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 2, 10, 9, 5, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 3, 0, 8, 1, 2, 10, 4, 9, 5, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 5, 2, 10, 5, 4, 2, 4, 0, 2, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 2, 10, 5, 3, 2, 5, 3, 5, 4, 3, 4, 8, -1, -1, -1, -1 },
|
|
{ 9, 5, 4, 2, 3, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 11, 2, 0, 8, 11, 4, 9, 5, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 5, 4, 0, 1, 5, 2, 3, 11, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 2, 1, 5, 2, 5, 8, 2, 8, 11, 4, 8, 5, -1, -1, -1, -1 },
|
|
{ 10, 3, 11, 10, 1, 3, 9, 5, 4, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 4, 9, 5, 0, 8, 1, 8, 10, 1, 8, 11, 10, -1, -1, -1, -1 },
|
|
{ 5, 4, 0, 5, 0, 11, 5, 11, 10, 11, 0, 3, -1, -1, -1, -1 },
|
|
{ 5, 4, 8, 5, 8, 10, 10, 8, 11, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 9, 7, 8, 5, 7, 9, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 9, 3, 0, 9, 5, 3, 5, 7, 3, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 7, 8, 0, 1, 7, 1, 5, 7, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 5, 3, 3, 5, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 9, 7, 8, 9, 5, 7, 10, 1, 2, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 10, 1, 2, 9, 5, 0, 5, 3, 0, 5, 7, 3, -1, -1, -1, -1 },
|
|
{ 8, 0, 2, 8, 2, 5, 8, 5, 7, 10, 5, 2, -1, -1, -1, -1 },
|
|
{ 2, 10, 5, 2, 5, 3, 3, 5, 7, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 7, 9, 5, 7, 8, 9, 3, 11, 2, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 9, 5, 7, 9, 7, 2, 9, 2, 0, 2, 7, 11, -1, -1, -1, -1 },
|
|
{ 2, 3, 11, 0, 1, 8, 1, 7, 8, 1, 5, 7, -1, -1, -1, -1 },
|
|
{ 11, 2, 1, 11, 1, 7, 7, 1, 5, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 9, 5, 8, 8, 5, 7, 10, 1, 3, 10, 3, 11, -1, -1, -1, -1 },
|
|
{ 5, 7, 0, 5, 0, 9, 7, 11, 0, 1, 0, 10, 11, 10, 0, -1 },
|
|
{ 11, 10, 0, 11, 0, 3, 10, 5, 0, 8, 0, 7, 5, 7, 0, -1 },
|
|
{ 11, 10, 5, 7, 11, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 10, 6, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 8, 3, 5, 10, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 9, 0, 1, 5, 10, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 8, 3, 1, 9, 8, 5, 10, 6, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 6, 5, 2, 6, 1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 6, 5, 1, 2, 6, 3, 0, 8, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 9, 6, 5, 9, 0, 6, 0, 2, 6, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 5, 9, 8, 5, 8, 2, 5, 2, 6, 3, 2, 8, -1, -1, -1, -1 },
|
|
{ 2, 3, 11, 10, 6, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 11, 0, 8, 11, 2, 0, 10, 6, 5, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 1, 9, 2, 3, 11, 5, 10, 6, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 5, 10, 6, 1, 9, 2, 9, 11, 2, 9, 8, 11, -1, -1, -1, -1 },
|
|
{ 6, 3, 11, 6, 5, 3, 5, 1, 3, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 8, 11, 0, 11, 5, 0, 5, 1, 5, 11, 6, -1, -1, -1, -1 },
|
|
{ 3, 11, 6, 0, 3, 6, 0, 6, 5, 0, 5, 9, -1, -1, -1, -1 },
|
|
{ 6, 5, 9, 6, 9, 11, 11, 9, 8, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 5, 10, 6, 4, 7, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 4, 3, 0, 4, 7, 3, 6, 5, 10, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 9, 0, 5, 10, 6, 8, 4, 7, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 10, 6, 5, 1, 9, 7, 1, 7, 3, 7, 9, 4, -1, -1, -1, -1 },
|
|
{ 6, 1, 2, 6, 5, 1, 4, 7, 8, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 2, 5, 5, 2, 6, 3, 0, 4, 3, 4, 7, -1, -1, -1, -1 },
|
|
{ 8, 4, 7, 9, 0, 5, 0, 6, 5, 0, 2, 6, -1, -1, -1, -1 },
|
|
{ 7, 3, 9, 7, 9, 4, 3, 2, 9, 5, 9, 6, 2, 6, 9, -1 },
|
|
{ 3, 11, 2, 7, 8, 4, 10, 6, 5, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 5, 10, 6, 4, 7, 2, 4, 2, 0, 2, 7, 11, -1, -1, -1, -1 },
|
|
{ 0, 1, 9, 4, 7, 8, 2, 3, 11, 5, 10, 6, -1, -1, -1, -1 },
|
|
{ 9, 2, 1, 9, 11, 2, 9, 4, 11, 7, 11, 4, 5, 10, 6, -1 },
|
|
{ 8, 4, 7, 3, 11, 5, 3, 5, 1, 5, 11, 6, -1, -1, -1, -1 },
|
|
{ 5, 1, 11, 5, 11, 6, 1, 0, 11, 7, 11, 4, 0, 4, 11, -1 },
|
|
{ 0, 5, 9, 0, 6, 5, 0, 3, 6, 11, 6, 3, 8, 4, 7, -1 },
|
|
{ 6, 5, 9, 6, 9, 11, 4, 7, 9, 7, 11, 9, -1, -1, -1, -1 },
|
|
{ 10, 4, 9, 6, 4, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 4, 10, 6, 4, 9, 10, 0, 8, 3, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 10, 0, 1, 10, 6, 0, 6, 4, 0, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 8, 3, 1, 8, 1, 6, 8, 6, 4, 6, 1, 10, -1, -1, -1, -1 },
|
|
{ 1, 4, 9, 1, 2, 4, 2, 6, 4, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 3, 0, 8, 1, 2, 9, 2, 4, 9, 2, 6, 4, -1, -1, -1, -1 },
|
|
{ 0, 2, 4, 4, 2, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 8, 3, 2, 8, 2, 4, 4, 2, 6, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 10, 4, 9, 10, 6, 4, 11, 2, 3, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 8, 2, 2, 8, 11, 4, 9, 10, 4, 10, 6, -1, -1, -1, -1 },
|
|
{ 3, 11, 2, 0, 1, 6, 0, 6, 4, 6, 1, 10, -1, -1, -1, -1 },
|
|
{ 6, 4, 1, 6, 1, 10, 4, 8, 1, 2, 1, 11, 8, 11, 1, -1 },
|
|
{ 9, 6, 4, 9, 3, 6, 9, 1, 3, 11, 6, 3, -1, -1, -1, -1 },
|
|
{ 8, 11, 1, 8, 1, 0, 11, 6, 1, 9, 1, 4, 6, 4, 1, -1 },
|
|
{ 3, 11, 6, 3, 6, 0, 0, 6, 4, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 6, 4, 8, 11, 6, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 7, 10, 6, 7, 8, 10, 8, 9, 10, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 7, 3, 0, 10, 7, 0, 9, 10, 6, 7, 10, -1, -1, -1, -1 },
|
|
{ 10, 6, 7, 1, 10, 7, 1, 7, 8, 1, 8, 0, -1, -1, -1, -1 },
|
|
{ 10, 6, 7, 10, 7, 1, 1, 7, 3, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 2, 6, 1, 6, 8, 1, 8, 9, 8, 6, 7, -1, -1, -1, -1 },
|
|
{ 2, 6, 9, 2, 9, 1, 6, 7, 9, 0, 9, 3, 7, 3, 9, -1 },
|
|
{ 7, 8, 0, 7, 0, 6, 6, 0, 2, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 7, 3, 2, 6, 7, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 2, 3, 11, 10, 6, 8, 10, 8, 9, 8, 6, 7, -1, -1, -1, -1 },
|
|
{ 2, 0, 7, 2, 7, 11, 0, 9, 7, 6, 7, 10, 9, 10, 7, -1 },
|
|
{ 1, 8, 0, 1, 7, 8, 1, 10, 7, 6, 7, 10, 2, 3, 11, -1 },
|
|
{ 11, 2, 1, 11, 1, 7, 10, 6, 1, 6, 7, 1, -1, -1, -1, -1 },
|
|
{ 8, 9, 6, 8, 6, 7, 9, 1, 6, 11, 6, 3, 1, 3, 6, -1 },
|
|
{ 0, 9, 1, 11, 6, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 7, 8, 0, 7, 0, 6, 3, 11, 0, 11, 6, 0, -1, -1, -1, -1 },
|
|
{ 7, 11, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 7, 6, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 3, 0, 8, 11, 7, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 1, 9, 11, 7, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 8, 1, 9, 8, 3, 1, 11, 7, 6, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 10, 1, 2, 6, 11, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 2, 10, 3, 0, 8, 6, 11, 7, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 2, 9, 0, 2, 10, 9, 6, 11, 7, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 6, 11, 7, 2, 10, 3, 10, 8, 3, 10, 9, 8, -1, -1, -1, -1 },
|
|
{ 7, 2, 3, 6, 2, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 7, 0, 8, 7, 6, 0, 6, 2, 0, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 2, 7, 6, 2, 3, 7, 0, 1, 9, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 6, 2, 1, 8, 6, 1, 9, 8, 8, 7, 6, -1, -1, -1, -1 },
|
|
{ 10, 7, 6, 10, 1, 7, 1, 3, 7, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 10, 7, 6, 1, 7, 10, 1, 8, 7, 1, 0, 8, -1, -1, -1, -1 },
|
|
{ 0, 3, 7, 0, 7, 10, 0, 10, 9, 6, 10, 7, -1, -1, -1, -1 },
|
|
{ 7, 6, 10, 7, 10, 8, 8, 10, 9, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 6, 8, 4, 11, 8, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 3, 6, 11, 3, 0, 6, 0, 4, 6, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 8, 6, 11, 8, 4, 6, 9, 0, 1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 9, 4, 6, 9, 6, 3, 9, 3, 1, 11, 3, 6, -1, -1, -1, -1 },
|
|
{ 6, 8, 4, 6, 11, 8, 2, 10, 1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 2, 10, 3, 0, 11, 0, 6, 11, 0, 4, 6, -1, -1, -1, -1 },
|
|
{ 4, 11, 8, 4, 6, 11, 0, 2, 9, 2, 10, 9, -1, -1, -1, -1 },
|
|
{ 10, 9, 3, 10, 3, 2, 9, 4, 3, 11, 3, 6, 4, 6, 3, -1 },
|
|
{ 8, 2, 3, 8, 4, 2, 4, 6, 2, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 4, 2, 4, 6, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 9, 0, 2, 3, 4, 2, 4, 6, 4, 3, 8, -1, -1, -1, -1 },
|
|
{ 1, 9, 4, 1, 4, 2, 2, 4, 6, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 8, 1, 3, 8, 6, 1, 8, 4, 6, 6, 10, 1, -1, -1, -1, -1 },
|
|
{ 10, 1, 0, 10, 0, 6, 6, 0, 4, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 4, 6, 3, 4, 3, 8, 6, 10, 3, 0, 3, 9, 10, 9, 3, -1 },
|
|
{ 10, 9, 4, 6, 10, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 4, 9, 5, 7, 6, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 8, 3, 4, 9, 5, 11, 7, 6, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 5, 0, 1, 5, 4, 0, 7, 6, 11, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 11, 7, 6, 8, 3, 4, 3, 5, 4, 3, 1, 5, -1, -1, -1, -1 },
|
|
{ 9, 5, 4, 10, 1, 2, 7, 6, 11, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 6, 11, 7, 1, 2, 10, 0, 8, 3, 4, 9, 5, -1, -1, -1, -1 },
|
|
{ 7, 6, 11, 5, 4, 10, 4, 2, 10, 4, 0, 2, -1, -1, -1, -1 },
|
|
{ 3, 4, 8, 3, 5, 4, 3, 2, 5, 10, 5, 2, 11, 7, 6, -1 },
|
|
{ 7, 2, 3, 7, 6, 2, 5, 4, 9, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 9, 5, 4, 0, 8, 6, 0, 6, 2, 6, 8, 7, -1, -1, -1, -1 },
|
|
{ 3, 6, 2, 3, 7, 6, 1, 5, 0, 5, 4, 0, -1, -1, -1, -1 },
|
|
{ 6, 2, 8, 6, 8, 7, 2, 1, 8, 4, 8, 5, 1, 5, 8, -1 },
|
|
{ 9, 5, 4, 10, 1, 6, 1, 7, 6, 1, 3, 7, -1, -1, -1, -1 },
|
|
{ 1, 6, 10, 1, 7, 6, 1, 0, 7, 8, 7, 0, 9, 5, 4, -1 },
|
|
{ 4, 0, 10, 4, 10, 5, 0, 3, 10, 6, 10, 7, 3, 7, 10, -1 },
|
|
{ 7, 6, 10, 7, 10, 8, 5, 4, 10, 4, 8, 10, -1, -1, -1, -1 },
|
|
{ 6, 9, 5, 6, 11, 9, 11, 8, 9, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 3, 6, 11, 0, 6, 3, 0, 5, 6, 0, 9, 5, -1, -1, -1, -1 },
|
|
{ 0, 11, 8, 0, 5, 11, 0, 1, 5, 5, 6, 11, -1, -1, -1, -1 },
|
|
{ 6, 11, 3, 6, 3, 5, 5, 3, 1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 2, 10, 9, 5, 11, 9, 11, 8, 11, 5, 6, -1, -1, -1, -1 },
|
|
{ 0, 11, 3, 0, 6, 11, 0, 9, 6, 5, 6, 9, 1, 2, 10, -1 },
|
|
{ 11, 8, 5, 11, 5, 6, 8, 0, 5, 10, 5, 2, 0, 2, 5, -1 },
|
|
{ 6, 11, 3, 6, 3, 5, 2, 10, 3, 10, 5, 3, -1, -1, -1, -1 },
|
|
{ 5, 8, 9, 5, 2, 8, 5, 6, 2, 3, 8, 2, -1, -1, -1, -1 },
|
|
{ 9, 5, 6, 9, 6, 0, 0, 6, 2, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 5, 8, 1, 8, 0, 5, 6, 8, 3, 8, 2, 6, 2, 8, -1 },
|
|
{ 1, 5, 6, 2, 1, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 3, 6, 1, 6, 10, 3, 8, 6, 5, 6, 9, 8, 9, 6, -1 },
|
|
{ 10, 1, 0, 10, 0, 6, 9, 5, 0, 5, 6, 0, -1, -1, -1, -1 },
|
|
{ 0, 3, 8, 5, 6, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 10, 5, 6, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 11, 5, 10, 7, 5, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 11, 5, 10, 11, 7, 5, 8, 3, 0, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 5, 11, 7, 5, 10, 11, 1, 9, 0, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 10, 7, 5, 10, 11, 7, 9, 8, 1, 8, 3, 1, -1, -1, -1, -1 },
|
|
{ 11, 1, 2, 11, 7, 1, 7, 5, 1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 8, 3, 1, 2, 7, 1, 7, 5, 7, 2, 11, -1, -1, -1, -1 },
|
|
{ 9, 7, 5, 9, 2, 7, 9, 0, 2, 2, 11, 7, -1, -1, -1, -1 },
|
|
{ 7, 5, 2, 7, 2, 11, 5, 9, 2, 3, 2, 8, 9, 8, 2, -1 },
|
|
{ 2, 5, 10, 2, 3, 5, 3, 7, 5, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 8, 2, 0, 8, 5, 2, 8, 7, 5, 10, 2, 5, -1, -1, -1, -1 },
|
|
{ 9, 0, 1, 5, 10, 3, 5, 3, 7, 3, 10, 2, -1, -1, -1, -1 },
|
|
{ 9, 8, 2, 9, 2, 1, 8, 7, 2, 10, 2, 5, 7, 5, 2, -1 },
|
|
{ 1, 3, 5, 3, 7, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 8, 7, 0, 7, 1, 1, 7, 5, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 9, 0, 3, 9, 3, 5, 5, 3, 7, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 9, 8, 7, 5, 9, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 5, 8, 4, 5, 10, 8, 10, 11, 8, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 5, 0, 4, 5, 11, 0, 5, 10, 11, 11, 3, 0, -1, -1, -1, -1 },
|
|
{ 0, 1, 9, 8, 4, 10, 8, 10, 11, 10, 4, 5, -1, -1, -1, -1 },
|
|
{ 10, 11, 4, 10, 4, 5, 11, 3, 4, 9, 4, 1, 3, 1, 4, -1 },
|
|
{ 2, 5, 1, 2, 8, 5, 2, 11, 8, 4, 5, 8, -1, -1, -1, -1 },
|
|
{ 0, 4, 11, 0, 11, 3, 4, 5, 11, 2, 11, 1, 5, 1, 11, -1 },
|
|
{ 0, 2, 5, 0, 5, 9, 2, 11, 5, 4, 5, 8, 11, 8, 5, -1 },
|
|
{ 9, 4, 5, 2, 11, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 2, 5, 10, 3, 5, 2, 3, 4, 5, 3, 8, 4, -1, -1, -1, -1 },
|
|
{ 5, 10, 2, 5, 2, 4, 4, 2, 0, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 3, 10, 2, 3, 5, 10, 3, 8, 5, 4, 5, 8, 0, 1, 9, -1 },
|
|
{ 5, 10, 2, 5, 2, 4, 1, 9, 2, 9, 4, 2, -1, -1, -1, -1 },
|
|
{ 8, 4, 5, 8, 5, 3, 3, 5, 1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 4, 5, 1, 0, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 8, 4, 5, 8, 5, 3, 9, 0, 5, 0, 3, 5, -1, -1, -1, -1 },
|
|
{ 9, 4, 5, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 4, 11, 7, 4, 9, 11, 9, 10, 11, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 8, 3, 4, 9, 7, 9, 11, 7, 9, 10, 11, -1, -1, -1, -1 },
|
|
{ 1, 10, 11, 1, 11, 4, 1, 4, 0, 7, 4, 11, -1, -1, -1, -1 },
|
|
{ 3, 1, 4, 3, 4, 8, 1, 10, 4, 7, 4, 11, 10, 11, 4, -1 },
|
|
{ 4, 11, 7, 9, 11, 4, 9, 2, 11, 9, 1, 2, -1, -1, -1, -1 },
|
|
{ 9, 7, 4, 9, 11, 7, 9, 1, 11, 2, 11, 1, 0, 8, 3, -1 },
|
|
{ 11, 7, 4, 11, 4, 2, 2, 4, 0, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 11, 7, 4, 11, 4, 2, 8, 3, 4, 3, 2, 4, -1, -1, -1, -1 },
|
|
{ 2, 9, 10, 2, 7, 9, 2, 3, 7, 7, 4, 9, -1, -1, -1, -1 },
|
|
{ 9, 10, 7, 9, 7, 4, 10, 2, 7, 8, 7, 0, 2, 0, 7, -1 },
|
|
{ 3, 7, 10, 3, 10, 2, 7, 4, 10, 1, 10, 0, 4, 0, 10, -1 },
|
|
{ 1, 10, 2, 8, 7, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 4, 9, 1, 4, 1, 7, 7, 1, 3, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 4, 9, 1, 4, 1, 7, 0, 8, 1, 8, 7, 1, -1, -1, -1, -1 },
|
|
{ 4, 0, 3, 7, 4, 3, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 4, 8, 7, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 9, 10, 8, 10, 11, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 3, 0, 9, 3, 9, 11, 11, 9, 10, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 1, 10, 0, 10, 8, 8, 10, 11, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 3, 1, 10, 11, 3, 10, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 2, 11, 1, 11, 9, 9, 11, 8, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 3, 0, 9, 3, 9, 11, 1, 2, 9, 2, 11, 9, -1, -1, -1, -1 },
|
|
{ 0, 2, 11, 8, 0, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 3, 2, 11, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 2, 3, 8, 2, 8, 10, 10, 8, 9, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 9, 10, 2, 0, 9, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 2, 3, 8, 2, 8, 10, 0, 1, 8, 1, 10, 8, -1, -1, -1, -1 },
|
|
{ 1, 10, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 1, 3, 8, 9, 1, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 9, 1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ 0, 3, 8, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 },
|
|
{ -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1 } };
|
|
|
|
static const int edgeTable[12][2] =
|
|
{ { 0, 1 }, { 1, 2 }, { 2, 3 }, { 3, 0 }, { 4, 5 }, { 5, 6 }, { 6, 7 }, { 7, 4 }, { 0, 4 }, { 1, 5 }, { 2, 6 }, { 3, 7 } };
|
|
|
|
int cubeIndex = 0;
|
|
if ( plane.distanceSquared( cell[0] ) < 0 ) cubeIndex |= 1;
|
|
if ( plane.distanceSquared( cell[1] ) < 0 ) cubeIndex |= 2;
|
|
if ( plane.distanceSquared( cell[2] ) < 0 ) cubeIndex |= 4;
|
|
if ( plane.distanceSquared( cell[3] ) < 0 ) cubeIndex |= 8;
|
|
if ( plane.distanceSquared( cell[4] ) < 0 ) cubeIndex |= 16;
|
|
if ( plane.distanceSquared( cell[5] ) < 0 ) cubeIndex |= 32;
|
|
if ( plane.distanceSquared( cell[6] ) < 0 ) cubeIndex |= 64;
|
|
if ( plane.distanceSquared( cell[7] ) < 0 ) cubeIndex |= 128;
|
|
|
|
if ( cubeIdxToCutEdgeBitfield[cubeIndex] == 0 )
|
|
{
|
|
return 0;
|
|
}
|
|
|
|
cvf::Vec3d edgeIntersections[12];
|
|
double normDistAlongEdge[12]{};
|
|
|
|
// Compute vertex coordinates on the edges where we have intersections
|
|
if ( cubeIdxToCutEdgeBitfield[cubeIndex] & 1 )
|
|
edgeIntersections[0] = planeLineIntersectionForMC( plane, cell[0], cell[1], &normDistAlongEdge[0] );
|
|
if ( cubeIdxToCutEdgeBitfield[cubeIndex] & 2 )
|
|
edgeIntersections[1] = planeLineIntersectionForMC( plane, cell[1], cell[2], &normDistAlongEdge[1] );
|
|
if ( cubeIdxToCutEdgeBitfield[cubeIndex] & 4 )
|
|
edgeIntersections[2] = planeLineIntersectionForMC( plane, cell[2], cell[3], &normDistAlongEdge[2] );
|
|
if ( cubeIdxToCutEdgeBitfield[cubeIndex] & 8 )
|
|
edgeIntersections[3] = planeLineIntersectionForMC( plane, cell[3], cell[0], &normDistAlongEdge[3] );
|
|
if ( cubeIdxToCutEdgeBitfield[cubeIndex] & 16 )
|
|
edgeIntersections[4] = planeLineIntersectionForMC( plane, cell[4], cell[5], &normDistAlongEdge[4] );
|
|
if ( cubeIdxToCutEdgeBitfield[cubeIndex] & 32 )
|
|
edgeIntersections[5] = planeLineIntersectionForMC( plane, cell[5], cell[6], &normDistAlongEdge[5] );
|
|
if ( cubeIdxToCutEdgeBitfield[cubeIndex] & 64 )
|
|
edgeIntersections[6] = planeLineIntersectionForMC( plane, cell[6], cell[7], &normDistAlongEdge[6] );
|
|
if ( cubeIdxToCutEdgeBitfield[cubeIndex] & 128 )
|
|
edgeIntersections[7] = planeLineIntersectionForMC( plane, cell[7], cell[4], &normDistAlongEdge[7] );
|
|
if ( cubeIdxToCutEdgeBitfield[cubeIndex] & 256 )
|
|
edgeIntersections[8] = planeLineIntersectionForMC( plane, cell[0], cell[4], &normDistAlongEdge[8] );
|
|
if ( cubeIdxToCutEdgeBitfield[cubeIndex] & 512 )
|
|
edgeIntersections[9] = planeLineIntersectionForMC( plane, cell[1], cell[5], &normDistAlongEdge[9] );
|
|
if ( cubeIdxToCutEdgeBitfield[cubeIndex] & 1024 )
|
|
edgeIntersections[10] = planeLineIntersectionForMC( plane, cell[2], cell[6], &normDistAlongEdge[10] );
|
|
if ( cubeIdxToCutEdgeBitfield[cubeIndex] & 2048 )
|
|
edgeIntersections[11] = planeLineIntersectionForMC( plane, cell[3], cell[7], &normDistAlongEdge[11] );
|
|
|
|
// Create the triangles
|
|
|
|
const int* triangleIndicesToCubeEdges = cubeIdxToTriangleIndices[cubeIndex];
|
|
cvf::uint triangleVxIdx = 0;
|
|
|
|
int cubeEdgeIdx = triangleIndicesToCubeEdges[triangleVxIdx];
|
|
while ( cubeEdgeIdx != -1 )
|
|
{
|
|
ClipVx cvx;
|
|
cvx.vx = edgeIntersections[cubeEdgeIdx];
|
|
cvx.normDistFromEdgeVx1 = normDistAlongEdge[cubeEdgeIdx];
|
|
cvx.clippedEdgeVx1Id = hexCornersIds[edgeTable[cubeEdgeIdx][0]];
|
|
cvx.clippedEdgeVx2Id = hexCornersIds[edgeTable[cubeEdgeIdx][1]];
|
|
|
|
( *triangleVxes ).push_back( cvx );
|
|
++triangleVxIdx;
|
|
|
|
cubeEdgeIdx = triangleIndicesToCubeEdges[triangleVxIdx];
|
|
}
|
|
|
|
cvf::uint triangleCount = triangleVxIdx / 3;
|
|
|
|
static const int edgeEdgeCutsToCellFace[12][12] = {
|
|
// 0 1 2 3 4 5 6 7 8 9 10 11
|
|
{ 6, 5, 5, 5, 3, 6, 6, 6, 3, 3, 6, 6 }, // 0
|
|
{ 5, 6, 5, 5, 6, 0, 6, 6, 6, 0, 0, 6 }, // 1 POS_I = 0
|
|
{ 5, 5, 6, 5, 6, 6, 2, 6, 6, 6, 2, 2 }, // 2 NEG_I = 1
|
|
{ 5, 5, 5, 6, 6, 6, 6, 1, 1, 6, 6, 1 }, // 3 POS_J = 2
|
|
{ 3, 6, 6, 6, 6, 4, 4, 4, 3, 3, 6, 6 }, // 4 NEG_J = 3
|
|
{ 6, 0, 6, 6, 4, 6, 4, 4, 6, 0, 0, 6 }, // 5 POS_K = 4
|
|
{ 6, 6, 2, 6, 4, 4, 6, 4, 6, 6, 2, 2 }, // 6 NEG_K = 5
|
|
{ 6, 6, 6, 1, 4, 4, 4, 6, 1, 6, 6, 1 }, // 7 NO_FACE = 6
|
|
{ 3, 6, 6, 1, 3, 6, 6, 1, 6, 3, 6, 1 }, // 8
|
|
{ 3, 0, 6, 6, 3, 0, 6, 6, 3, 6, 0, 6 }, // 9
|
|
{ 6, 0, 2, 6, 6, 0, 2, 6, 6, 0, 6, 2 }, // 10
|
|
{ 6, 6, 2, 1, 6, 6, 2, 1, 1, 6, 2, 6 } // 11
|
|
};
|
|
|
|
( *cellFaceForEachTriangleEdge ).clear();
|
|
( *cellFaceForEachTriangleEdge ).resize( triangleVxIdx, 6 );
|
|
|
|
for ( cvf::uint tIdx = 0; tIdx < triangleCount; ++tIdx )
|
|
{
|
|
cvf::uint triVxIdx = 3 * tIdx;
|
|
|
|
int cubeEdgeIdx1 = triangleIndicesToCubeEdges[triVxIdx];
|
|
int cubeEdgeIdx2 = triangleIndicesToCubeEdges[triVxIdx + 1];
|
|
int cubeEdgeIdx3 = triangleIndicesToCubeEdges[triVxIdx + 2];
|
|
|
|
( *cellFaceForEachTriangleEdge )[triVxIdx + 0] = edgeEdgeCutsToCellFace[cubeEdgeIdx1][cubeEdgeIdx2];
|
|
( *cellFaceForEachTriangleEdge )[triVxIdx + 1] = edgeEdgeCutsToCellFace[cubeEdgeIdx2][cubeEdgeIdx3];
|
|
( *cellFaceForEachTriangleEdge )[triVxIdx + 2] = edgeEdgeCutsToCellFace[cubeEdgeIdx3][cubeEdgeIdx1];
|
|
}
|
|
|
|
#if 0
|
|
// Calculate what triangle edges are representing the cut of a cell face
|
|
// Do this by counting the times two specific cube edges are used for a triangle edge.
|
|
// Internal edges will have a count of 2, while external edges only 1
|
|
|
|
(*isTriEdgeCellContour).clear();
|
|
(*isTriEdgeCellContour).resize(triangleVxIdx);
|
|
|
|
|
|
int triangleEdgeCount[12][12] = {
|
|
{ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 },
|
|
{ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 },
|
|
{ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 },
|
|
{ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 },
|
|
{ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 },
|
|
{ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 },
|
|
{ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 },
|
|
{ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 },
|
|
{ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 },
|
|
{ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 },
|
|
{ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 },
|
|
{ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 }
|
|
};
|
|
|
|
(*isTriEdgeCellContour).clear();
|
|
(*isTriEdgeCellContour).resize(triangleVxIdx, false);
|
|
|
|
for (cvf::uint tIdx = 0; tIdx < triangleCount; ++tIdx)
|
|
{
|
|
cvf::uint triVxIdx = 3 * tIdx;
|
|
int cubeEdgeIdx1 = triangleIndicesToCubeEdges[triVxIdx];
|
|
int cubeEdgeIdx2 = triangleIndicesToCubeEdges[triVxIdx + 1];
|
|
int cubeEdgeIdx3 = triangleIndicesToCubeEdges[triVxIdx + 2];
|
|
|
|
cubeEdgeIdx1 < cubeEdgeIdx2 ? ++triangleEdgeCount[cubeEdgeIdx1][cubeEdgeIdx2] : ++triangleEdgeCount[cubeEdgeIdx2][cubeEdgeIdx1];
|
|
cubeEdgeIdx2 < cubeEdgeIdx3 ? ++triangleEdgeCount[cubeEdgeIdx2][cubeEdgeIdx3] : ++triangleEdgeCount[cubeEdgeIdx3][cubeEdgeIdx2];
|
|
cubeEdgeIdx3 < cubeEdgeIdx1 ? ++triangleEdgeCount[cubeEdgeIdx3][cubeEdgeIdx1] : ++triangleEdgeCount[cubeEdgeIdx1][cubeEdgeIdx3];
|
|
}
|
|
|
|
for (cvf::uint tIdx = 0; tIdx < triangleCount; ++tIdx)
|
|
{
|
|
cvf::uint triVxIdx = 3 * tIdx;
|
|
|
|
int cubeEdgeIdx1 = triangleIndicesToCubeEdges[triVxIdx];
|
|
int cubeEdgeIdx2 = triangleIndicesToCubeEdges[triVxIdx + 1];
|
|
int cubeEdgeIdx3 = triangleIndicesToCubeEdges[triVxIdx + 2];
|
|
|
|
// We have a contour if the count is exactly 1.
|
|
|
|
(*isTriEdgeCellContour)[triVxIdx + 0] = (1 == (cubeEdgeIdx1 < cubeEdgeIdx2 ? triangleEdgeCount[cubeEdgeIdx1][cubeEdgeIdx2] : triangleEdgeCount[cubeEdgeIdx2][cubeEdgeIdx1]));
|
|
(*isTriEdgeCellContour)[triVxIdx + 1] = (1 == (cubeEdgeIdx2 < cubeEdgeIdx3 ? triangleEdgeCount[cubeEdgeIdx2][cubeEdgeIdx3] : triangleEdgeCount[cubeEdgeIdx3][cubeEdgeIdx2]));
|
|
(*isTriEdgeCellContour)[triVxIdx + 2] = (1 == (cubeEdgeIdx3 < cubeEdgeIdx1 ? triangleEdgeCount[cubeEdgeIdx3][cubeEdgeIdx1] : triangleEdgeCount[cubeEdgeIdx1][cubeEdgeIdx3]));
|
|
}
|
|
|
|
#endif
|
|
return triangleCount;
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
/// Based on description and implementation from Paul Bourke:
|
|
///
|
|
/// http://paulbourke.net/geometry/polygonise/
|
|
///
|
|
/// Note that the element is turned inside-out compared to what we use elsewhere in caf/ResInsight
|
|
/// So the winding of all the sides are opposite.
|
|
/// 4-----4------5
|
|
/// /| /| k POS_I = 0
|
|
/// 7 8 5 9 | NEG_I = 1
|
|
/// / | / | | POS_J = 2
|
|
/// 7------6-----6 | | NEG_J = 3
|
|
/// | 0-----0--|---1 *------i POS_K = 4
|
|
/// 11 / 10 / / NEG_K = 5
|
|
/// | 3 | 1 / NO_FACE = 6
|
|
/// |/ |/ j
|
|
/// 3------2-----2
|
|
///
|
|
// The cellFaceForEachTriangleEdge refer to the edge after the corresponding triangle vertex.
|
|
|
|
/*
|
|
Based on description and implementation from Paul Bourke:
|
|
http://local.wasp.uwa.edu.au/~pbourke/geometry/polygonise/
|
|
|
|
Polygonise a tetrahedron given its vertices within a cube
|
|
This is an alternative algorithm to polygonisegrid.
|
|
It results in a smoother surface but more triangular facets.
|
|
|
|
+ 0 + 0
|
|
/|\ /|\
|
|
/ | \ / | \
|
|
/ | \ / | \
|
|
/ | \ / | \
|
|
/ | \ / 2 | \
|
|
/ | \ / __--+_ \
|
|
+-------------+ 1 3 +--__ --_ \
|
|
3 \ | / ---__ --_\
|
|
\ | / ---__-\
|
|
\ | / --+ 1
|
|
\ | /
|
|
\ | /
|
|
\|/ 2 is behind 1 and 3
|
|
+ 2
|
|
|
|
|
|
Build six tets from a cube to make sure the split direction is equal for opposite sides.
|
|
Surface normals are pointing outward.
|
|
|
|
See following comment is taken from http://www.iue.tuwien.ac.at/phd/wessner/node32.html
|
|
The decompositions of a cube into five tetrahedra yields an orientation switch of two opposite diagonal face edges of
|
|
the cube. Due to this fact, the tessellation of one cube, as part of a larger cubic grid, forces a particular
|
|
tessellation of all neighboring cubes to guarantee a conformal mesh. This means that, if such five-decompositions cubes
|
|
are stacked together to a chain, the mesh of each cube must be rotated by an angle of 90 deg
|
|
|
|
The tessellation makes sure opposite faces are divided along the same line
|
|
See figure http://www.ics.uci.edu/~eppstein/projects/tetra/
|
|
|
|
4, 5, 6, 0
|
|
0, 1, 5, 6
|
|
0, 2, 1, 6
|
|
4, 6, 7, 0
|
|
0, 7, 3, 6
|
|
0, 3, 2, 6
|
|
|
|
Introduces the additional diagonal edges in the Hex from 12 up to and including 18:
|
|
0 2 // 12 NEG_K
|
|
0 5 // 13 NEG_J
|
|
1 6 // 14 POS_I
|
|
3 6 // 15 POS_J
|
|
0 7 // 16 NEG_I
|
|
4 6 // 17 POS_K
|
|
0 6 // 18 Internal Diagonal
|
|
|
|
/// 4-----4------5
|
|
/// /| /| k POS_I = 0
|
|
/// 7 8 17 5 9 | NEG_I = 1
|
|
/// / | 13 / | | POS_J = 2
|
|
/// 7------6-----6 14| | NEG_J = 3
|
|
/// |16 0-----0--|---1 *------i POS_K = 4
|
|
/// 11 / 15 10 / / NEG_K = 5
|
|
/// | 3 12 | 1 / NO_FACE = 6
|
|
/// |/ |/ j
|
|
/// 3------2-----2
|
|
|
|
|
|
*/
|
|
//--------------------------------------------------------------------------------------------------
|
|
int HexGridIntersectionTools::planeHexIntersectionMCTet( const cvf::Plane& plane,
|
|
const cvf::Vec3d cell[8],
|
|
const size_t hexCornersIds[8],
|
|
std::vector<ClipVx>* triangleVxes,
|
|
std::vector<int>* cellFaceForEachTriangleEdge )
|
|
{
|
|
std::array<double, 8> cellCornerSqDistToPlane = {
|
|
plane.distanceSquared( cell[0] ),
|
|
plane.distanceSquared( cell[1] ),
|
|
plane.distanceSquared( cell[2] ),
|
|
plane.distanceSquared( cell[3] ),
|
|
plane.distanceSquared( cell[4] ),
|
|
plane.distanceSquared( cell[5] ),
|
|
plane.distanceSquared( cell[6] ),
|
|
plane.distanceSquared( cell[7] ),
|
|
};
|
|
|
|
int cubeIndex = 0;
|
|
if ( cellCornerSqDistToPlane[0] < 0 ) cubeIndex |= 1;
|
|
if ( cellCornerSqDistToPlane[1] < 0 ) cubeIndex |= 2;
|
|
if ( cellCornerSqDistToPlane[2] < 0 ) cubeIndex |= 4;
|
|
if ( cellCornerSqDistToPlane[3] < 0 ) cubeIndex |= 8;
|
|
if ( cellCornerSqDistToPlane[4] < 0 ) cubeIndex |= 16;
|
|
if ( cellCornerSqDistToPlane[5] < 0 ) cubeIndex |= 32;
|
|
if ( cellCornerSqDistToPlane[6] < 0 ) cubeIndex |= 64;
|
|
if ( cellCornerSqDistToPlane[7] < 0 ) cubeIndex |= 128;
|
|
|
|
if ( cubeIndex == 0 || cubeIndex == 255 ) return 0;
|
|
|
|
int tetCount = 0;
|
|
tetCount += planeMcTetIntersection( plane,
|
|
cell,
|
|
hexCornersIds,
|
|
cellCornerSqDistToPlane.data(),
|
|
{ 4, 5, 6, 0 },
|
|
triangleVxes,
|
|
cellFaceForEachTriangleEdge );
|
|
tetCount += planeMcTetIntersection( plane,
|
|
cell,
|
|
hexCornersIds,
|
|
cellCornerSqDistToPlane.data(),
|
|
{ 0, 1, 5, 6 },
|
|
triangleVxes,
|
|
cellFaceForEachTriangleEdge );
|
|
tetCount += planeMcTetIntersection( plane,
|
|
cell,
|
|
hexCornersIds,
|
|
cellCornerSqDistToPlane.data(),
|
|
{ 0, 2, 1, 6 },
|
|
triangleVxes,
|
|
cellFaceForEachTriangleEdge );
|
|
tetCount += planeMcTetIntersection( plane,
|
|
cell,
|
|
hexCornersIds,
|
|
cellCornerSqDistToPlane.data(),
|
|
{ 4, 6, 7, 0 },
|
|
triangleVxes,
|
|
cellFaceForEachTriangleEdge );
|
|
tetCount += planeMcTetIntersection( plane,
|
|
cell,
|
|
hexCornersIds,
|
|
cellCornerSqDistToPlane.data(),
|
|
{ 0, 7, 3, 6 },
|
|
triangleVxes,
|
|
cellFaceForEachTriangleEdge );
|
|
tetCount += planeMcTetIntersection( plane,
|
|
cell,
|
|
hexCornersIds,
|
|
cellCornerSqDistToPlane.data(),
|
|
{ 0, 3, 2, 6 },
|
|
triangleVxes,
|
|
cellFaceForEachTriangleEdge );
|
|
|
|
return tetCount;
|
|
}
|
|
|
|
//--------------------------------------------------------------------------------------------------
|
|
///
|
|
//--------------------------------------------------------------------------------------------------
|
|
cvf::uint HexGridIntersectionTools::planeMcTetIntersection( const cvf::Plane& plane,
|
|
const cvf::Vec3d hexCell[8],
|
|
const size_t hexCornersIds[8],
|
|
const double cornerDistToPlane[8],
|
|
const std::array<int, 4>& tetCell,
|
|
std::vector<ClipVx>* triangleVxes,
|
|
std::vector<int>* cellFaceForEachTriangleEdge )
|
|
{
|
|
// clang-format off
|
|
|
|
static const int edgeEdgeCutsToCellFace[19][19] = {
|
|
// 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 4--------------4---------------5
|
|
{ 6, 5, 5, 5, 3, 6, 6, 6, 3, 3, 6, 6, 5, 3, 6, 6, 6, 6, 6 }, // 0 /|\__ __//|
|
|
{ 5, 6, 5, 5, 6, 0, 6, 6, 6, 0, 0, 6, 5, 6, 0, 6, 6, 6, 6 }, // 1 POS_I = 0 / | \___ _/ / |
|
|
{ 5, 5, 6, 5, 6, 6, 2, 6, 6, 6, 2, 2, 5, 6, 6, 2, 6, 6, 6 }, // 2 NEG_I = 1 7 | \__ __/ / | k
|
|
{ 5, 5, 5, 6, 6, 6, 6, 1, 1, 6, 6, 1, 5, 6, 6, 6, 1, 6, 6 }, // 3 POS_J = 2 / | 17___ _/ 5 | |
|
|
{ 3, 6, 6, 6, 6, 4, 4, 4, 3, 3, 6, 6, 6, 3, 6, 6, 6, 4, 6 }, // 4 NEG_J = 3 / 8 \___ / 9 |
|
|
{ 6, 0, 6, 6, 4, 6, 4, 4, 6, 0, 0, 6, 6, 6, 0, 6, 6, 4, 6 }, // 5 POS_K = 4 / | __/ \___ / | |
|
|
{ 6, 6, 2, 6, 4, 4, 6, 4, 6, 6, 2, 2, 6, 6, 6, 2, 6, 4, 6 }, // 6 NEG_K = 5 7---------------6--------------6 | *------i
|
|
{ 6, 6, 6, 1, 4, 4, 4, 6, 1, 6, 6, 1, 6, 6, 6, 6, 1, 4, 6 }, // 7 NO_FACE = 6 |\_ | __13 ____/_/|\_ | /
|
|
{ 3, 6, 6, 1, 3, 6, 6, 1, 6, 3, 6, 1, 6, 3, 6, 6, 1, 6, 6 }, // 8 | 16 | __/ _18_/ __/ | 14 | /
|
|
{ 3, 0, 6, 6, 3, 0, 6, 6, 3, 6, 0, 6, 6, 3, 0, 6, 6, 6, 6 }, // 9 | \_ | __/____/ _/ | \_ | j
|
|
{ 6, 0, 2, 6, 6, 0, 2, 6, 6, 0, 6, 2, 6, 6, 0, 2, 6, 6, 6 }, // 10 | \|__/__/ __/ | \|
|
|
{ 6, 6, 2, 1, 6, 6, 2, 1, 1, 6, 2, 6, 6, 6, 6, 2, 1, 6, 6 }, // 11 | 0-----------_/----0-----|------1
|
|
{ 5, 5, 5, 5, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6 }, // 12 11 / \__ __15 10 /
|
|
{ 3, 6, 6, 6, 3, 6, 6, 6, 3, 3, 6, 6, 6, 6, 6, 6, 6, 6, 6 }, // 13 | / \_ _/ | /
|
|
{ 6, 0, 6, 6, 6, 0, 6, 6, 6, 0, 0, 6, 6, 6, 6, 6, 6, 6, 6 }, // 14 | 3 __/ \___ | 1
|
|
{ 6, 6, 2, 6, 6, 6, 2, 6, 6, 6, 2, 2, 6, 6, 6, 6, 6, 6, 6 }, // 15 | / __/ 12__ | /
|
|
{ 6, 6, 6, 1, 6, 6, 6, 1, 1, 6, 6, 1, 6, 6, 6, 6, 6, 6, 6 }, // 16 | / __/ \___ | /
|
|
{ 6, 6, 6, 6, 4, 4, 4, 4, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6 }, // 17 |/__/ \___|/
|
|
{ 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6 }, // 18 3---------------2--------------2
|
|
};
|
|
|
|
// clang-format on
|
|
|
|
static const int cellCornerCellCornerToEdge[8][8] = {
|
|
// 0 1 2 3 4 5 6 7
|
|
{ -1, 0, 12, 3, 8, 13, 18, 16 }, // 0
|
|
{ 0, -1, 1, -1, -1, 9, 14, -1 }, // 1
|
|
{ 12, 1, -1, 2, -1, -1, 10, -1 }, // 2
|
|
{ 3, -1, 2, -1, -1, -1, 15, 11 }, // 3
|
|
{ 8, -1, -1, -1, -1, 4, 17, 7 }, // 4
|
|
{ 13, 9, -1, -1, 4, -1, 5, -1 }, // 5
|
|
{ 18, 14, 10, 15, 17, 5, -1, 6 }, // 6
|
|
{ 16, -1, -1, 11, 7, -1, 6, -1 }, // 7
|
|
};
|
|
|
|
cvf::uint ntri = 0;
|
|
|
|
int triindex = 0;
|
|
if ( cornerDistToPlane[tetCell[0]] < 0 ) triindex |= 1;
|
|
if ( cornerDistToPlane[tetCell[1]] < 0 ) triindex |= 2;
|
|
if ( cornerDistToPlane[tetCell[2]] < 0 ) triindex |= 4;
|
|
if ( cornerDistToPlane[tetCell[3]] < 0 ) triindex |= 8;
|
|
|
|
auto clipEdgeFunc = [&]( int hexCornerIdx0, int hexCornerIdx1 ) {
|
|
ClipVx cvx;
|
|
cvx.vx =
|
|
planeLineIntersectionForMC( plane, hexCell[hexCornerIdx0], hexCell[hexCornerIdx1], &cvx.normDistFromEdgeVx1 );
|
|
cvx.clippedEdgeVx1Id = hexCornersIds[hexCornerIdx0];
|
|
cvx.clippedEdgeVx2Id = hexCornersIds[hexCornerIdx1];
|
|
return cvx;
|
|
};
|
|
|
|
auto addCellFaceStatusForTriangleEdges = [&]( int e11, int e12, int e21, int e22, int e31, int e32 ) {
|
|
int cutEdge1 = cellCornerCellCornerToEdge[e11][e12];
|
|
int cutEdge2 = cellCornerCellCornerToEdge[e21][e22];
|
|
int cutEdge3 = cellCornerCellCornerToEdge[e31][e32];
|
|
|
|
CVF_ASSERT( cutEdge1 >= 0 );
|
|
CVF_ASSERT( cutEdge2 >= 0 );
|
|
CVF_ASSERT( cutEdge3 >= 0 );
|
|
|
|
cellFaceForEachTriangleEdge->emplace_back( edgeEdgeCutsToCellFace[cutEdge1][cutEdge2] );
|
|
cellFaceForEachTriangleEdge->emplace_back( edgeEdgeCutsToCellFace[cutEdge2][cutEdge3] );
|
|
cellFaceForEachTriangleEdge->emplace_back( edgeEdgeCutsToCellFace[cutEdge3][cutEdge1] );
|
|
};
|
|
|
|
switch ( triindex )
|
|
{
|
|
case 0x00:
|
|
case 0x0F:
|
|
break;
|
|
case 0x0E:
|
|
case 0x01:
|
|
{
|
|
triangleVxes->emplace_back( clipEdgeFunc( tetCell[0], tetCell[1] ) );
|
|
triangleVxes->emplace_back( clipEdgeFunc( tetCell[0], tetCell[2] ) );
|
|
triangleVxes->emplace_back( clipEdgeFunc( tetCell[0], tetCell[3] ) );
|
|
|
|
addCellFaceStatusForTriangleEdges( tetCell[0], tetCell[1], tetCell[0], tetCell[2], tetCell[0], tetCell[3] );
|
|
ntri++;
|
|
}
|
|
break;
|
|
case 0x0D:
|
|
case 0x02:
|
|
{
|
|
triangleVxes->emplace_back( clipEdgeFunc( tetCell[1], tetCell[0] ) );
|
|
triangleVxes->emplace_back( clipEdgeFunc( tetCell[1], tetCell[3] ) );
|
|
triangleVxes->emplace_back( clipEdgeFunc( tetCell[1], tetCell[2] ) );
|
|
|
|
addCellFaceStatusForTriangleEdges( tetCell[1], tetCell[0], tetCell[1], tetCell[3], tetCell[1], tetCell[2] );
|
|
ntri++;
|
|
}
|
|
break;
|
|
case 0x0C:
|
|
case 0x03:
|
|
{
|
|
triangleVxes->emplace_back( clipEdgeFunc( tetCell[0], tetCell[3] ) );
|
|
ClipVx vx1 = clipEdgeFunc( tetCell[0], tetCell[2] );
|
|
triangleVxes->push_back( vx1 );
|
|
ClipVx vx2 = clipEdgeFunc( tetCell[1], tetCell[3] );
|
|
triangleVxes->push_back( vx2 );
|
|
|
|
addCellFaceStatusForTriangleEdges( tetCell[0], tetCell[3], tetCell[0], tetCell[2], tetCell[1], tetCell[3] );
|
|
|
|
ntri++;
|
|
|
|
triangleVxes->push_back( vx2 );
|
|
triangleVxes->emplace_back( clipEdgeFunc( tetCell[1], tetCell[2] ) );
|
|
triangleVxes->push_back( vx1 );
|
|
|
|
addCellFaceStatusForTriangleEdges( tetCell[1], tetCell[3], tetCell[1], tetCell[2], tetCell[0], tetCell[2] );
|
|
ntri++;
|
|
}
|
|
break;
|
|
case 0x0B:
|
|
case 0x04:
|
|
{
|
|
triangleVxes->emplace_back( clipEdgeFunc( tetCell[2], tetCell[0] ) );
|
|
triangleVxes->emplace_back( clipEdgeFunc( tetCell[2], tetCell[1] ) );
|
|
triangleVxes->emplace_back( clipEdgeFunc( tetCell[2], tetCell[3] ) );
|
|
|
|
addCellFaceStatusForTriangleEdges( tetCell[2], tetCell[0], tetCell[2], tetCell[1], tetCell[2], tetCell[3] );
|
|
ntri++;
|
|
}
|
|
break;
|
|
case 0x0A:
|
|
case 0x05:
|
|
{
|
|
ClipVx vx0 = clipEdgeFunc( tetCell[0], tetCell[1] );
|
|
triangleVxes->push_back( vx0 );
|
|
ClipVx vx1 = clipEdgeFunc( tetCell[2], tetCell[3] );
|
|
triangleVxes->push_back( vx1 );
|
|
ClipVx vx2 = clipEdgeFunc( tetCell[0], tetCell[3] );
|
|
triangleVxes->push_back( vx2 );
|
|
|
|
addCellFaceStatusForTriangleEdges( tetCell[0], tetCell[1], tetCell[2], tetCell[3], tetCell[0], tetCell[3] );
|
|
ntri++;
|
|
|
|
triangleVxes->push_back( vx0 );
|
|
triangleVxes->emplace_back( clipEdgeFunc( tetCell[1], tetCell[2] ) );
|
|
triangleVxes->push_back( vx1 );
|
|
|
|
addCellFaceStatusForTriangleEdges( tetCell[0], tetCell[1], tetCell[1], tetCell[2], tetCell[2], tetCell[3] );
|
|
|
|
ntri++;
|
|
}
|
|
break;
|
|
case 0x09:
|
|
case 0x06:
|
|
{
|
|
ClipVx vx0 = clipEdgeFunc( tetCell[0], tetCell[1] );
|
|
triangleVxes->push_back( vx0 );
|
|
ClipVx vx1 = clipEdgeFunc( tetCell[1], tetCell[3] );
|
|
triangleVxes->push_back( vx1 );
|
|
ClipVx vx2 = clipEdgeFunc( tetCell[2], tetCell[3] );
|
|
triangleVxes->push_back( vx2 );
|
|
|
|
addCellFaceStatusForTriangleEdges( tetCell[0], tetCell[1], tetCell[1], tetCell[3], tetCell[2], tetCell[3] );
|
|
|
|
ntri++;
|
|
|
|
triangleVxes->push_back( vx0 );
|
|
triangleVxes->emplace_back( clipEdgeFunc( tetCell[0], tetCell[2] ) );
|
|
triangleVxes->push_back( vx2 );
|
|
|
|
addCellFaceStatusForTriangleEdges( tetCell[0], tetCell[1], tetCell[0], tetCell[2], tetCell[2], tetCell[3] );
|
|
|
|
ntri++;
|
|
}
|
|
break;
|
|
case 0x07:
|
|
case 0x08:
|
|
{
|
|
triangleVxes->emplace_back( clipEdgeFunc( tetCell[3], tetCell[0] ) );
|
|
triangleVxes->emplace_back( clipEdgeFunc( tetCell[3], tetCell[2] ) );
|
|
triangleVxes->emplace_back( clipEdgeFunc( tetCell[3], tetCell[1] ) );
|
|
|
|
addCellFaceStatusForTriangleEdges( tetCell[3], tetCell[0], tetCell[3], tetCell[2], tetCell[3], tetCell[1] );
|
|
|
|
ntri++;
|
|
}
|
|
break;
|
|
}
|
|
|
|
return ntri;
|
|
}
|
|
|
|
} // namespace caf
|
|
|
|
// clang-format on
|
|
} //namespace external
|