From ec5a9e1682ed7efbd4aae469e537a2ef969e21ba Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Jacob=20St=C3=B8ren?= Date: Fri, 9 Mar 2018 18:01:13 +0100 Subject: [PATCH] #2552 Generalized the flatting transforms to be able to reuse for simwells and wellpaths --- .../RivIntersectionGeometryGenerator.cpp | 182 ++++-------------- .../RivIntersectionGeometryGenerator.h | 2 +- .../Intersections/RivIntersectionPartMgr.cpp | 4 +- .../Intersections/RivIntersectionPartMgr.h | 2 +- .../Intersections/RivSectionFlattner.h | 148 ++++++++++++++ 5 files changed, 187 insertions(+), 151 deletions(-) create mode 100644 ApplicationCode/ModelVisualization/Intersections/RivSectionFlattner.h diff --git a/ApplicationCode/ModelVisualization/Intersections/RivIntersectionGeometryGenerator.cpp b/ApplicationCode/ModelVisualization/Intersections/RivIntersectionGeometryGenerator.cpp index d134b593f1..85bc42b96a 100644 --- a/ApplicationCode/ModelVisualization/Intersections/RivIntersectionGeometryGenerator.cpp +++ b/ApplicationCode/ModelVisualization/Intersections/RivIntersectionGeometryGenerator.cpp @@ -36,7 +36,8 @@ #include "cvfPrimitiveSetIndexedUInt.h" #include "cvfScalarMapper.h" #include "cvfRay.h" -//#include "cvfTrace.h" + +#include "RivSectionFlattner.h" //-------------------------------------------------------------------------------------------------- @@ -69,127 +70,42 @@ RivIntersectionGeometryGenerator::~RivIntersectionGeometryGenerator() } -//-------------------------------------------------------------------------------------------------- -/// Origo in the intersection of the ray P1-ExtrDir with the XY plane -/// Ez in upwards extrusionDir -/// Ey normal tio the section pplane -/// Ex in plane along p1-p2 -//-------------------------------------------------------------------------------------------------- -cvf::Mat4d calculateSectionLocalFlatteningCS(const cvf::Vec3d& p1, const cvf::Vec3d& p2, const cvf::Vec3d& extrusionDir) -{ - using namespace cvf; - - Vec3d Ez = extrusionDir.z() > 0.0 ? extrusionDir: -extrusionDir; - - Vec3d sectionLineDir = p2 - p1; - sectionLineDir.normalize(); - - Vec3d Ey = Ez ^ sectionLineDir; - Ey.normalize(); - Vec3d Ex = Ey ^ Ez; - Ex.normalize(); - - Ray extrusionRay; - extrusionRay.setOrigin(p1); - - if (p1.z() > 0) extrusionRay.setDirection(-Ez); - else extrusionRay.setDirection(Ez); - - - Vec3d tr(Vec3d::ZERO); - extrusionRay.planeIntersect(Plane(0.0, 0.0 , 1.0, 0.0), &tr); - - return Mat4d(Ex[0], Ey[0], Ez[0], tr[0], - Ex[1], Ey[1], Ez[1], tr[1], - Ex[2], Ey[2], Ez[2], tr[2], - 0.0, 0.0, 0.0, 1.0); -} - //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- void RivIntersectionGeometryGenerator::calculateSegementTransformPrLinePoint() { - cvf::Vec3d displayOffset = m_hexGrid->displayOffset(); - cvf::Mat4d invSectionCS = cvf::Mat4d::fromTranslation(-displayOffset); - - m_segementTransformPrLinePoint.clear(); - - double previousSectionFlattenedEndPosX = m_horizontalLengthAlongWellToPolylineStart; - cvf::Vec3d previousSectionOrigo(cvf::Vec3d::ZERO); - - - for ( size_t pLineIdx = 0; pLineIdx < m_polyLines.size(); ++pLineIdx ) + if ( m_isFlattened ) { - m_segementTransformPrLinePoint.emplace_back(); - const std::vector& polyLine = m_polyLines[pLineIdx]; + if ( !(m_polyLines.size() && m_polyLines.back().size()) ) return; - size_t pointCount = polyLine.size(); + cvf::Vec3d startOffset ={ m_horizontalLengthAlongWellToPolylineStart, 0.0, m_polyLines[0][0].z() }; - size_t lIdx = 0; - while ( lIdx < pointCount ) + for ( size_t pLineIdx = 0; pLineIdx < m_polyLines.size(); ++pLineIdx ) { - size_t idxToNextP = indexToNextValidPoint(polyLine, m_extrusionDirection, lIdx); - if (idxToNextP == size_t(-1)) - { - size_t inc = 0; - while ((lIdx + inc) < pointCount) - { - m_segementTransformPrLinePoint.back().push_back(invSectionCS); - ++inc; - } - break; - } - - if (m_isFlattened) - { - cvf::Vec3d p1 = polyLine[lIdx]; - cvf::Vec3d p2 = polyLine[idxToNextP]; - - cvf::Mat4d sectionLocalCS = calculateSectionLocalFlatteningCS(p1, p2, m_extrusionDirection); - - if ( pLineIdx == 0 && lIdx == 0 ) previousSectionOrigo = sectionLocalCS.translation(); - - previousSectionFlattenedEndPosX += (sectionLocalCS.translation() - previousSectionOrigo).length(); - previousSectionOrigo = sectionLocalCS.translation(); - - invSectionCS = sectionLocalCS.getInverted(); - cvf::Vec3d flattenedTranslation(previousSectionFlattenedEndPosX, 0.0, 0.0); - - invSectionCS.setTranslation(invSectionCS.translation() + flattenedTranslation ); - } - - size_t inc = 0; - while ((lIdx + inc) < idxToNextP) - { - m_segementTransformPrLinePoint.back().push_back(invSectionCS); - inc++; - } - - lIdx = idxToNextP; + const std::vector& polyLine = m_polyLines[pLineIdx]; + m_segementTransformPrLinePoint.emplace_back(RivSectionFlattner::calculateFlatteningCSsForPolyline(polyLine, + m_extrusionDirection, + startOffset, + &startOffset)); } } + else + { + m_segementTransformPrLinePoint.clear(); - // for (auto mx: m_segementTransformPrLinePoint[0]) - // { - // cvf::String text; - // - // for (int r = 0; r < 4; ++r) - // { - // for (int c = 0; c < 4; ++c) - // { - // text += cvf::String::number(mx(r, c)); - // - // if (r * c < 9) - // { - // text += " "; - // } - // } - // text += "\n"; - // } - // - // cvf::Trace::show( text ); - // } + cvf::Mat4d invSectionCS = cvf::Mat4d::fromTranslation(-m_hexGrid->displayOffset()); + + for ( const auto & polyLine : m_polyLines ) + { + m_segementTransformPrLinePoint.emplace_back(); + std::vector& segmentTransforms = m_segementTransformPrLinePoint.back(); + for ( size_t lIdx = 0; lIdx < polyLine.size(); ++lIdx ) + { + segmentTransforms.push_back(invSectionCS); + } + } + } } //-------------------------------------------------------------------------------------------------- @@ -239,7 +155,7 @@ void RivIntersectionGeometryGenerator::calculateArrays() size_t lIdx = 0; while ( lIdx < lineCount - 1) { - size_t idxToNextP = indexToNextValidPoint(polyLine, m_extrusionDirection, lIdx); + size_t idxToNextP = RivSectionFlattner::indexToNextValidPoint(polyLine, m_extrusionDirection, lIdx); if (idxToNextP == size_t(-1)) break; @@ -569,34 +485,6 @@ cvf::ref RivIntersectionGeometryGenerator::createPointsFromPol } -//-------------------------------------------------------------------------------------------------- -/// Find the next point in the polyline that avoids making the line nearly parallel to the extrusion direction -/// Returns size_t(-1) if no point is found -//-------------------------------------------------------------------------------------------------- -size_t RivIntersectionGeometryGenerator::indexToNextValidPoint(const std::vector& polyLine, - const cvf::Vec3d extrDir, - size_t idxToStartOfLineSegment) -{ - size_t lineCount = polyLine.size(); - if ( !(idxToStartOfLineSegment + 1 < lineCount) ) return -1; - - - cvf::Vec3d p1 = polyLine[idxToStartOfLineSegment]; - - for ( size_t lIdx = idxToStartOfLineSegment+1; lIdx < lineCount; ++lIdx ) - { - cvf::Vec3d p2 = polyLine[lIdx]; - cvf::Vec3d p1p2 = p2 - p1; - - if ( (p1p2 - (p1p2 * extrDir)*extrDir).length() > 0.1 ) - { - return lIdx; - } - } - - return -1; -} - //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- @@ -635,23 +523,23 @@ RimIntersection* RivIntersectionGeometryGenerator::crossSection() const //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -cvf::Mat4d RivIntersectionGeometryGenerator::unflattenTransformMatrix(const cvf::Vec3d& intersectionPointUtm) +cvf::Mat4d RivIntersectionGeometryGenerator::unflattenTransformMatrix(const cvf::Vec3d& intersectionPointFlat) { cvf::Vec3d startPt = cvf::Vec3d::ZERO; int polyLineIdx = -1; int segIdx = -1; - for (size_t i = 0; i < m_flattenedOrOffsettedPolyLines.size(); i++) + for (size_t pLineIdx = 0; pLineIdx < m_flattenedOrOffsettedPolyLines.size(); pLineIdx++) { - std::vector pts = m_flattenedOrOffsettedPolyLines[i]; - for(size_t j = 0; j < pts.size(); j++) + std::vector polyLine = m_flattenedOrOffsettedPolyLines[pLineIdx]; + for(size_t pIdx = 0; pIdx < polyLine.size(); pIdx++) { // Assumes ascending sorted list - if (j > 0 && intersectionPointUtm.x() < pts[j].x()) + if (pIdx > 0 && intersectionPointFlat.x() < polyLine[pIdx].x()) { - polyLineIdx = static_cast(i); - segIdx = static_cast(j) - 1; - startPt = pts[segIdx]; + polyLineIdx = static_cast(pLineIdx); + segIdx = static_cast(pIdx) - 1; + startPt = polyLine[segIdx]; break; } } diff --git a/ApplicationCode/ModelVisualization/Intersections/RivIntersectionGeometryGenerator.h b/ApplicationCode/ModelVisualization/Intersections/RivIntersectionGeometryGenerator.h index e92cb0b41d..8338e073a8 100644 --- a/ApplicationCode/ModelVisualization/Intersections/RivIntersectionGeometryGenerator.h +++ b/ApplicationCode/ModelVisualization/Intersections/RivIntersectionGeometryGenerator.h @@ -75,7 +75,7 @@ public: RimIntersection* crossSection() const; - cvf::Mat4d unflattenTransformMatrix(const cvf::Vec3d& intersectionPointUtm); + cvf::Mat4d unflattenTransformMatrix(const cvf::Vec3d& intersectionPointFlat); private: void calculateArrays(); diff --git a/ApplicationCode/ModelVisualization/Intersections/RivIntersectionPartMgr.cpp b/ApplicationCode/ModelVisualization/Intersections/RivIntersectionPartMgr.cpp index 40cc3a2384..33a7653e30 100644 --- a/ApplicationCode/ModelVisualization/Intersections/RivIntersectionPartMgr.cpp +++ b/ApplicationCode/ModelVisualization/Intersections/RivIntersectionPartMgr.cpp @@ -974,9 +974,9 @@ const RimIntersection* RivIntersectionPartMgr::intersection() const //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -cvf::Mat4d RivIntersectionPartMgr::unflattenTransformMatrix(const cvf::Vec3d& intersectionPointUtm) +cvf::Mat4d RivIntersectionPartMgr::unflattenTransformMatrix(const cvf::Vec3d& intersectionPointFlat) { - return m_crossSectionGenerator->unflattenTransformMatrix(intersectionPointUtm); + return m_crossSectionGenerator->unflattenTransformMatrix(intersectionPointFlat); } //-------------------------------------------------------------------------------------------------- diff --git a/ApplicationCode/ModelVisualization/Intersections/RivIntersectionPartMgr.h b/ApplicationCode/ModelVisualization/Intersections/RivIntersectionPartMgr.h index 933532b5b3..5009e81451 100644 --- a/ApplicationCode/ModelVisualization/Intersections/RivIntersectionPartMgr.h +++ b/ApplicationCode/ModelVisualization/Intersections/RivIntersectionPartMgr.h @@ -76,7 +76,7 @@ public: const RimIntersection* intersection() const; - cvf::Mat4d unflattenTransformMatrix(const cvf::Vec3d& intersectionPointUtm); + cvf::Mat4d unflattenTransformMatrix(const cvf::Vec3d& intersectionPointFlat); public: static void calculateEclipseTextureCoordinates(cvf::Vec2fArray* textureCoords, diff --git a/ApplicationCode/ModelVisualization/Intersections/RivSectionFlattner.h b/ApplicationCode/ModelVisualization/Intersections/RivSectionFlattner.h new file mode 100644 index 0000000000..c0973123b7 --- /dev/null +++ b/ApplicationCode/ModelVisualization/Intersections/RivSectionFlattner.h @@ -0,0 +1,148 @@ +#pragma once + + +class RivSectionFlattner +{ +public: + //-------------------------------------------------------------------------------------------------- + /// Returns the next index higher than idxToStartOfLineSegment that makes the line + // polyline[idxToStartOfLineSegment] .. polyline[nextIdx] not parallel to extrDir + /// + /// Returns size_t(-1) if no point is found + //-------------------------------------------------------------------------------------------------- + static size_t indexToNextValidPoint(const std::vector& polyLine, + const cvf::Vec3d extrDir, + size_t idxToStartOfLineSegment) + { + size_t lineCount = polyLine.size(); + if ( !(idxToStartOfLineSegment + 1 < lineCount) ) return -1; + + + cvf::Vec3d p1 = polyLine[idxToStartOfLineSegment]; + + for ( size_t lIdx = idxToStartOfLineSegment+1; lIdx < lineCount; ++lIdx ) + { + cvf::Vec3d p2 = polyLine[lIdx]; + cvf::Vec3d p1p2 = p2 - p1; + + if ( (p1p2 - (p1p2 * extrDir)*extrDir).length() > 0.1 ) + { + return lIdx; + } + } + + return -1; + } + + //-------------------------------------------------------------------------------------------------- + /// Returns one CS pr point, valid for the next segment + //-------------------------------------------------------------------------------------------------- + static std::vector calculateFlatteningCSsForPolyline(const std::vector & polyLine, + const cvf::Vec3d& extrusionDir, + const cvf::Vec3d& startOffset, + cvf::Vec3d* endOffset) + { + CVF_ASSERT(endOffset); + size_t pointCount = polyLine.size(); + CVF_ASSERT(pointCount > 1); + + std::vector segmentTransforms; + segmentTransforms.reserve(pointCount); + + // Find initial transform, used if all is vertical + + cvf::Mat4d invSectionCS; + { + cvf::Vec3d p1 = polyLine[0]; + cvf::Vec3d p2 = polyLine[1]; + + cvf::Mat4d sectionLocalCS = calculateSectionLocalFlatteningCS(p1, p2, extrusionDir); + cvf::Mat4d invSectionCS = sectionLocalCS.getInverted(); + invSectionCS.setTranslation(invSectionCS.translation() + startOffset); + } + + cvf::Vec3d previousFlattenedSectionEndPoint = startOffset; + + size_t lIdx = 0; + while ( lIdx < pointCount ) + { + size_t idxToNextP = indexToNextValidPoint(polyLine, extrusionDir, lIdx); + + // If the rest is nearly parallel to extrusionDir, use the current inverse matrix for the rest of the points + + if ( idxToNextP == size_t(-1) ) + { + size_t inc = 0; + while ( (lIdx + inc) < pointCount ) + { + segmentTransforms.push_back(invSectionCS); + ++inc; + } + break; + } + + cvf::Vec3d p1 = polyLine[lIdx]; + cvf::Vec3d p2 = polyLine[idxToNextP]; + + cvf::Mat4d sectionLocalCS = calculateSectionLocalFlatteningCS(p1, p2, extrusionDir); + invSectionCS = sectionLocalCS.getInverted(); + cvf::Vec3d flattenedSectionEndPoint = p2.getTransformedPoint(invSectionCS); + + invSectionCS.setTranslation(invSectionCS.translation() + previousFlattenedSectionEndPoint ); + + previousFlattenedSectionEndPoint += flattenedSectionEndPoint; + + // Assign the matrix to the points in between + + size_t inc = 0; + while ( (lIdx + inc) < idxToNextP ) + { + segmentTransforms.push_back(invSectionCS); + inc++; + } + + lIdx = idxToNextP; + } + + *endOffset = previousFlattenedSectionEndPoint; + return segmentTransforms; + } + +private: + + //-------------------------------------------------------------------------------------------------- + /// Origo in P1 + /// Ez in upwards extrusionDir + /// Ey normal to the section plane + /// Ex in plane along p1-p2 + //-------------------------------------------------------------------------------------------------- + static cvf::Mat4d calculateSectionLocalFlatteningCS(const cvf::Vec3d& p1, + const cvf::Vec3d& p2, + const cvf::Vec3d& extrusionDir) + { + using namespace cvf; + + Vec3d Ez = extrusionDir.z() > 0.0 ? extrusionDir: -extrusionDir; + + Vec3d sectionLineDir = p2 - p1; + + if ( cvf::GeometryTools::getAngle(sectionLineDir, extrusionDir) < 0.01 ) + { + sectionLineDir = Ez.perpendicularVector(); + } + + Vec3d Ey = Ez ^ sectionLineDir; + Ey.normalize(); + Vec3d Ex = Ey ^ Ez; + Ex.normalize(); + + return Mat4d(Ex[0], Ey[0], Ez[0], p1[0], + Ex[1], Ey[1], Ez[1], p1[1], + Ex[2], Ey[2], Ez[2], p1[2], + 0.0, 0.0, 0.0, 1.0); + } + + +}; + +