From 0685078ab373d587924d0dc8e3d0928dea141523 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?J=C3=B8rgen=20Herje?= <82032112+jorgenherje@users.noreply.github.com> Date: Tue, 20 Jun 2023 10:08:10 +0200 Subject: [PATCH] Fix automatic part id detection for Fault Reactivation Result, and resampling bug in RigWellLogCurveData * Fix resampling bug and refactor code - Fix bug for resampling, prevent index increment. - Refactor functions into static functions. - Fix issue in interpolateSegment not using correct indices for depthType != resamplingDepthType - Add unit tests * Change WellAllocationPlot to use step left Remove dummy point and utilize step left for WellAllocationPont * Fix bug in creating resampled values and depths for RigWellLogCurveData * Fix automatic part detection for Fault Reactivation Result - Fix incorrect automatic part detection - Set default distance to intersection to 1.0 [m] --- .../GeoMech/GeoMechDataModel/RigFemPart.cpp | 31 +++ .../GeoMech/GeoMechDataModel/RigFemPart.h | 3 + .../Flow/RimWellAllocationPlot.cpp | 6 +- .../RimGeoMechFaultReactivationResult.cpp | 39 +++- .../RimGeoMechFaultReactivationResult.h | 2 + .../RigWellLogCurveData.cpp | 166 ++++++++++----- .../ReservoirDataModel/RigWellLogCurveData.h | 24 ++- .../UnitTests/CMakeLists_files.cmake | 2 + .../UnitTests/RigWellLogCurveData-Test.cpp | 197 ++++++++++++++++++ 9 files changed, 395 insertions(+), 75 deletions(-) create mode 100644 ApplicationLibCode/UnitTests/RigWellLogCurveData-Test.cpp diff --git a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPart.cpp b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPart.cpp index 8f7e1bdf7f..02d09fadaa 100644 --- a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPart.cpp +++ b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPart.cpp @@ -124,6 +124,37 @@ const int* RigFemPart::connectivities( size_t elementIdx ) const return &m_allElementConnectivities[m_elementConnectivityStartIndices[elementIdx]]; } +//-------------------------------------------------------------------------------------------------- +/// Returns state of success for fill element coordinates +//-------------------------------------------------------------------------------------------------- +bool RigFemPart::fillElementCoordinates( size_t elementIdx, std::array& coordinates ) const +{ + const auto elemType = elementType( elementIdx ); + const int* elemConnectivity = connectivities( elementIdx ); + const auto nodeCount = RigFemTypes::elementNodeCount( elemType ); + + // Only handle element of hex for now - false success + if ( nodeCount != 8 ) return false; + + // Retrieve the node indices + std::vector nodeIndices; + for ( int i = 0; i < nodeCount; ++i ) + { + const int nodeIdx = elemConnectivity[i]; + nodeIndices.push_back( nodeIdx ); + } + + // Fill coordinates for each node + const auto& partNodes = nodes(); + for ( int i = 0; i < nodeIndices.size(); ++i ) + { + coordinates[i].set( partNodes.coordinates[nodeIndices[i]] ); + } + + // Return true success + return true; +} + //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- diff --git a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPart.h b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPart.h index 34c9da23a6..91b7ffa189 100644 --- a/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPart.h +++ b/ApplicationLibCode/GeoMech/GeoMechDataModel/RigFemPart.h @@ -27,6 +27,7 @@ #include "cvfBoundingBox.h" #include "cvfObject.h" #include "cvfVector3.h" +#include #include #include @@ -64,6 +65,8 @@ public: bool isHexahedron( size_t elementIdx ) const; const int* connectivities( size_t elementIdx ) const; + bool fillElementCoordinates( size_t elementIdx, std::array& coordinates ) const; + size_t elementNodeResultIdx( int elementIdx, int elmLocalNodeIdx ) const; size_t elementNodeResultCount() const; int nodeIdxFromElementNodeResultIdx( size_t elmNodeResultIdx ) const; diff --git a/ApplicationLibCode/ProjectDataModel/Flow/RimWellAllocationPlot.cpp b/ApplicationLibCode/ProjectDataModel/Flow/RimWellAllocationPlot.cpp index d20ae95b77..67c8db9122 100644 --- a/ApplicationLibCode/ProjectDataModel/Flow/RimWellAllocationPlot.cpp +++ b/ApplicationLibCode/ProjectDataModel/Flow/RimWellAllocationPlot.cpp @@ -337,10 +337,6 @@ void RimWellAllocationPlot::updateFromWell() accFlow = ( m_flowType == ACCUMULATED ? wfCalculator->accumulatedTracerFlowPrConnection( tracerName, brIdx ) : wfCalculator->tracerFlowPrConnection( tracerName, brIdx ) ); - // Insert the first depth position again, to add a value pair - curveDepthValues.insert( curveDepthValues.begin(), curveDepthValues[0] ); - accFlow.insert( accFlow.begin(), 0.0 ); - if ( m_flowType == ACCUMULATED && brIdx == 0 && !accFlow.empty() ) // Add fictitious point to -1 // for first branch { @@ -474,6 +470,8 @@ void RimWellAllocationPlot::addStackedCurve( const QString& tracerNa curve->setColor( getTracerColor( tracerName ) ); } + curve->setInterpolation( RiuQwtPlotCurveDefines::CurveInterpolationEnum::INTERPOLATION_STEP_LEFT ); + plotTrack->addCurve( curve ); curve->loadDataAndUpdate( true ); diff --git a/ApplicationLibCode/ProjectDataModel/GeoMech/RimGeoMechFaultReactivationResult.cpp b/ApplicationLibCode/ProjectDataModel/GeoMech/RimGeoMechFaultReactivationResult.cpp index f138bf0fed..6ea2bf307a 100644 --- a/ApplicationLibCode/ProjectDataModel/GeoMech/RimGeoMechFaultReactivationResult.cpp +++ b/ApplicationLibCode/ProjectDataModel/GeoMech/RimGeoMechFaultReactivationResult.cpp @@ -27,6 +27,7 @@ #include "RigFemPartCollection.h" #include "RigGeoMechCaseData.h" +#include "RigHexIntersectionTools.h" #include "RigReservoirGridTools.h" #include "RimGeoMechCase.h" @@ -52,6 +53,8 @@ #include "cvfBoundingBox.h" +#include + CAF_PDM_SOURCE_INIT( RimGeoMechFaultReactivationResult, "RimGeoMechFaultReactivationResult" ); //-------------------------------------------------------------------------------------------------- @@ -64,7 +67,10 @@ RimGeoMechFaultReactivationResult::RimGeoMechFaultReactivationResult() CAF_PDM_InitFieldNoDefault( &m_intersection, "Intersection", "Intersection" ); - CAF_PDM_InitField( &m_distanceFromIntersection, "FaceDistanceFromIntersection", 0.0, "Face Distance From Intersection" ); + CAF_PDM_InitField( &m_distanceFromIntersection, + "FaceDistanceFromIntersection", + m_defaultDistanceFromIntersection, + "Face Distance From Intersection" ); CAF_PDM_InitField( &m_widthOutsideIntersection, "FaceWidthOutsideIntersection", 0.0, "Face Width Outside Intersection" ); CAF_PDM_InitFieldNoDefault( &m_createFaultReactivationResult, "CreateReactivationResult", "" ); @@ -156,7 +162,8 @@ void RimGeoMechFaultReactivationResult::fieldChangedByUi( const caf::PdmFieldHan } if ( changedField == &m_createFaultReactivationResult && m_intersection() ) { - onLoadDataAndUpdate(); + createWellGeometry(); + createWellLogCurves(); } } @@ -172,7 +179,7 @@ void RimGeoMechFaultReactivationResult::defineEditorAttribute( const caf::PdmFie caf::PdmUiPushButtonEditorAttribute* attrib = dynamic_cast( attribute ); if ( attrib ) { - attrib->m_buttonText = "Create"; + attrib->m_buttonText = "Create Plot"; } } } @@ -332,19 +339,31 @@ void RimGeoMechFaultReactivationResult::createWellLogCurves() //-------------------------------------------------------------------------------------------------- int RimGeoMechFaultReactivationResult::getPartIndexFromPoint( const RigFemPartCollection* const partCollection, const cvf::Vec3d& point ) const { - int idx = 0; + const int idx = 0; if ( !partCollection ) return idx; + // Find candidates for intersected global elements const cvf::BoundingBox intersectingBb( point, point ); - std::vector intersectedGlobalElementIndices; - partCollection->findIntersectingGlobalElementIndices( intersectingBb, &intersectedGlobalElementIndices ); + std::vector intersectedGlobalElementIndexCandidates; + partCollection->findIntersectingGlobalElementIndices( intersectingBb, &intersectedGlobalElementIndexCandidates ); - if ( intersectedGlobalElementIndices.empty() ) return idx; + if ( intersectedGlobalElementIndexCandidates.empty() ) return idx; - // Utilize first intersected element to detect part for point - const auto [partId, elementIndex] = partCollection->partIdAndElementIndex( intersectedGlobalElementIndices.front() ); - idx = partId; + // Iterate through global element candidates and check if point is in hexCorners + for ( const auto& globalElementIndex : intersectedGlobalElementIndexCandidates ) + { + const auto [part, elementIndex] = partCollection->partAndElementIndex( globalElementIndex ); + // Find nodes from element + std::array coordinates; + const bool isSuccess = part->fillElementCoordinates( elementIndex, coordinates ); + if ( !isSuccess ) continue; + + const bool isPointInCell = RigHexIntersectionTools::isPointInCell( point, coordinates.data() ); + if ( isPointInCell ) return part->elementPartId(); + } + + // Utilize first part to have an id return idx; } diff --git a/ApplicationLibCode/ProjectDataModel/GeoMech/RimGeoMechFaultReactivationResult.h b/ApplicationLibCode/ProjectDataModel/GeoMech/RimGeoMechFaultReactivationResult.h index a82bde7534..c751276def 100644 --- a/ApplicationLibCode/ProjectDataModel/GeoMech/RimGeoMechFaultReactivationResult.h +++ b/ApplicationLibCode/ProjectDataModel/GeoMech/RimGeoMechFaultReactivationResult.h @@ -77,4 +77,6 @@ private: caf::PdmField m_faceAWellPathPartIndex; caf::PdmField m_faceBWellPathPartIndex; + + const double m_defaultDistanceFromIntersection = 1.0; // [m] Default value from intersection and into each part }; diff --git a/ApplicationLibCode/ReservoirDataModel/RigWellLogCurveData.cpp b/ApplicationLibCode/ReservoirDataModel/RigWellLogCurveData.cpp index dcfcc2de10..c090d07001 100644 --- a/ApplicationLibCode/ReservoirDataModel/RigWellLogCurveData.cpp +++ b/ApplicationLibCode/ReservoirDataModel/RigWellLogCurveData.cpp @@ -345,41 +345,57 @@ cvf::ref RigWellLogCurveData::calculateResampledCurveData( return reSampledData; } -void RigWellLogCurveData::interpolateSegment( RiaDefines::DepthTypeEnum resamplingDepthType, - double depthValue, - size_t firstIndex, - std::vector& xValues, - std::map>& resampledDepths, - const double eps ) const +//-------------------------------------------------------------------------------------------------- +/// The function creates and adds interpolated values for property and depths. The interpolated +/// values are added to the end of the resampledValues. Target depth is added to the end of the +/// resampledDepths vector for the resampling type. The depth values for remaining depth types +/// are created by linear interpolation between first and second depth value of the resampling type. +//-------------------------------------------------------------------------------------------------- +void RigWellLogCurveData::createAndAddInterpolatedSegmentValueAndDepths( + std::vector& resampledValues, + std::map>& resampledDepths, + RiaDefines::DepthTypeEnum resamplingDepthType, + double targetDepthValue, + size_t firstIndex, + size_t secondIndex, + const std::map>& originalDepths, + const std::vector& propertyValues, + double eps ) { - auto depthIt = m_depths.find( resamplingDepthType ); + if ( !originalDepths.contains( resamplingDepthType ) ) return; - size_t secondIndex = firstIndex + 1; + const auto& depthValues = originalDepths.find( resamplingDepthType )->second; + if ( firstIndex >= depthValues.size() || secondIndex >= depthValues.size() ) return; - double depth0 = depthIt->second[firstIndex]; - double depth1 = depthIt->second[secondIndex]; - double x0 = m_propertyValues[firstIndex]; - double x1 = m_propertyValues[secondIndex]; + double depth0 = depthValues[firstIndex]; + double depth1 = depthValues[secondIndex]; + double x0 = propertyValues[firstIndex]; + double x1 = propertyValues[secondIndex]; double slope = 0.0; if ( std::fabs( depth1 - depth0 ) > eps ) { slope = ( x1 - x0 ) / ( depth1 - depth0 ); } - double xValue = slope * ( depthValue - depth0 ) + x0; - xValues.push_back( xValue ); + double resampledValue = slope * ( targetDepthValue - depth0 ) + x0; + resampledValues.push_back( resampledValue ); + resampledDepths[resamplingDepthType].push_back( targetDepthValue ); - for ( auto depthTypeValuesPair : m_depths ) + // Add depth values for remaining depth types + for ( const auto& [depthType, depthTypeValues] : originalDepths ) { - if ( depthTypeValuesPair.first != resamplingDepthType ) - { - double otherDepth0 = depthTypeValuesPair.second[0]; - double otherDepth1 = depthTypeValuesPair.second[1]; - double otherSlope = ( otherDepth1 - otherDepth0 ) / ( depth1 - depth0 ); - resampledDepths[depthTypeValuesPair.first].push_back( otherSlope * ( depthValue - depth0 ) + otherDepth0 ); - } + if ( depthType == resamplingDepthType ) continue; + if ( firstIndex >= depthTypeValues.size() || secondIndex >= depthTypeValues.size() ) continue; + + double otherDepth0 = depthTypeValues[firstIndex]; + double otherDepth1 = depthTypeValues[secondIndex]; + double otherSlope = ( otherDepth1 - otherDepth0 ) / ( depth1 - depth0 ); + resampledDepths[depthType].push_back( otherSlope * ( targetDepthValue - depth0 ) + otherDepth0 ); } } +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- bool isLeftOf( double x1, double x2, bool reverseOrder, double eps ) { if ( reverseOrder ) @@ -389,6 +405,9 @@ bool isLeftOf( double x1, double x2, bool reverseOrder, double eps ) return x2 - x1 > eps; } +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- bool isRightOf( double x1, double x2, bool reverseOrder, double eps ) { return isLeftOf( x2, x1, reverseOrder, eps ); @@ -397,53 +416,59 @@ bool isRightOf( double x1, double x2, bool reverseOrder, double eps ) //-------------------------------------------------------------------------------------------------- /// //-------------------------------------------------------------------------------------------------- -cvf::ref RigWellLogCurveData::calculateResampledCurveData( RiaDefines::DepthTypeEnum resamplingDepthType, - const std::vector& depths ) const +std::pair, std::map>> + RigWellLogCurveData::createResampledValuesAndDepths( RiaDefines::DepthTypeEnum resamplingDepthType, + const std::vector& targetDepths, + const std::map>& originalDepths, + const std::vector& propertyValues ) { const double eps = 1.0e-8; - std::vector xValues; - - std::map> resampledDepths; - resampledDepths.insert( std::make_pair( resamplingDepthType, depths ) ); - - auto depthIt = m_depths.find( resamplingDepthType ); - - cvf::ref reSampledData = new RigWellLogCurveData; - - if ( depthIt == m_depths.end() || depthIt->second.empty() ) return reSampledData; + auto depthIt = originalDepths.find( resamplingDepthType ); + if ( depthIt == originalDepths.end() || depthIt->second.empty() ) return {}; + const auto& depthValues = depthIt->second; bool reverseOrder = resamplingDepthType == RiaDefines::DepthTypeEnum::CONNECTION_NUMBER; + std::vector resampledValues; + std::map> resampledDepths; + size_t segmentSearchStartIdx = 0; - for ( auto depth : depths ) + for ( const auto& targetDepth : targetDepths ) { bool foundPoint = false; - for ( size_t segmentStartIdx = segmentSearchStartIdx; segmentStartIdx < depthIt->second.size(); ++segmentStartIdx ) + for ( size_t segmentStartIdx = segmentSearchStartIdx; segmentStartIdx < depthValues.size(); ++segmentStartIdx ) { - if ( std::fabs( depthIt->second[segmentStartIdx] - depth ) < eps ) // already have this depth point, - // reuse it + if ( std::fabs( depthValues[segmentStartIdx] - targetDepth ) < eps ) // already have this depth point, reuse it { - xValues.push_back( m_propertyValues[segmentStartIdx] ); - // Copy all depth types for this segment - for ( auto depthTypeValuesPair : m_depths ) + // Copy all depth data for this segment + resampledValues.push_back( propertyValues[segmentStartIdx] ); + for ( const auto& [depthType, depthValues] : originalDepths ) { - if ( depthTypeValuesPair.first != resamplingDepthType ) - { - resampledDepths[depthTypeValuesPair.first].push_back( depthTypeValuesPair.second[segmentStartIdx] ); - } + resampledDepths[depthType].push_back( depthValues[segmentStartIdx] ); } segmentSearchStartIdx = segmentStartIdx + 1; foundPoint = true; break; } - else if ( segmentStartIdx < depthIt->second.size() - 1 ) + else if ( segmentStartIdx > 0 && segmentStartIdx < depthValues.size() ) { - double minDepthSegment = std::min( depthIt->second[segmentStartIdx], depthIt->second[segmentStartIdx + 1] ); - double maxDepthSegment = std::max( depthIt->second[segmentStartIdx], depthIt->second[segmentStartIdx + 1] ); - if ( cvf::Math::valueInRange( depth, minDepthSegment, maxDepthSegment ) ) + // Interpolate between current and previous depth point + const size_t firstIndex = segmentStartIdx - 1; + const size_t secondIndex = segmentStartIdx; + double minDepthSegment = std::min( depthValues[firstIndex], depthValues[secondIndex] ); + double maxDepthSegment = std::max( depthValues[firstIndex], depthValues[secondIndex] ); + if ( cvf::Math::valueInRange( targetDepth, minDepthSegment, maxDepthSegment ) ) { - interpolateSegment( resamplingDepthType, depth, segmentStartIdx, xValues, resampledDepths, eps ); + createAndAddInterpolatedSegmentValueAndDepths( resampledValues, + resampledDepths, + resamplingDepthType, + targetDepth, + firstIndex, + secondIndex, + originalDepths, + propertyValues, + eps ); segmentSearchStartIdx = segmentStartIdx; foundPoint = true; break; @@ -452,17 +477,36 @@ cvf::ref RigWellLogCurveData::calculateResampledCurveData( } if ( !foundPoint ) { - if ( isLeftOf( depth, depthIt->second.front(), reverseOrder, eps ) ) + if ( isLeftOf( targetDepth, depthValues.front(), reverseOrder, eps ) ) { // Extrapolate from front two - interpolateSegment( resamplingDepthType, depth, 0, xValues, resampledDepths, eps ); + const size_t firstIndex = 0; + const size_t secondIndex = 1; + createAndAddInterpolatedSegmentValueAndDepths( resampledValues, + resampledDepths, + resamplingDepthType, + targetDepth, + firstIndex, + secondIndex, + originalDepths, + propertyValues, + eps ); foundPoint = true; } - else if ( isRightOf( depth, depthIt->second.back(), reverseOrder, eps ) ) + else if ( isRightOf( targetDepth, depthValues.back(), reverseOrder, eps ) ) { // Extrapolate from end two - const size_t N = depthIt->second.size() - 1; - interpolateSegment( resamplingDepthType, depth, N - 1, xValues, resampledDepths, eps ); + const size_t N = depthValues.size() - 1; + const size_t N_1 = N - 1; + createAndAddInterpolatedSegmentValueAndDepths( resampledValues, + resampledDepths, + resamplingDepthType, + targetDepth, + N_1, + N, + originalDepths, + propertyValues, + eps ); foundPoint = true; } } @@ -470,6 +514,18 @@ cvf::ref RigWellLogCurveData::calculateResampledCurveData( CAF_ASSERT( foundPoint ); } + return std::make_pair( resampledValues, resampledDepths ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +cvf::ref RigWellLogCurveData::calculateResampledCurveData( RiaDefines::DepthTypeEnum resamplingDepthType, + const std::vector& depths ) const +{ + const auto [xValues, resampledDepths] = createResampledValuesAndDepths( resamplingDepthType, depths, m_depths, m_propertyValues ); + + cvf::ref reSampledData = new RigWellLogCurveData; reSampledData->setValuesAndDepths( xValues, resampledDepths, m_rkbDiff, m_depthUnit, true, m_useLogarithmicScale ); return reSampledData; } diff --git a/ApplicationLibCode/ReservoirDataModel/RigWellLogCurveData.h b/ApplicationLibCode/ReservoirDataModel/RigWellLogCurveData.h index 1e1dbbc923..8f6a243a12 100644 --- a/ApplicationLibCode/ReservoirDataModel/RigWellLogCurveData.h +++ b/ApplicationLibCode/ReservoirDataModel/RigWellLogCurveData.h @@ -82,12 +82,24 @@ public: cvf::ref calculateResampledCurveData( double newMeasuredDepthStepSize ) const; cvf::ref calculateResampledCurveData( RiaDefines::DepthTypeEnum resamplingDepthType, const std::vector& depths ) const; - void interpolateSegment( RiaDefines::DepthTypeEnum resamplingDepthType, - double depthValue, - size_t firstIndex, - std::vector& xValues, - std::map>& resampledDepths, - const double eps ) const; + + // Made static due to unit testing + static void createAndAddInterpolatedSegmentValueAndDepths( std::vector& resampledValues, + std::map>& resampledDepths, + RiaDefines::DepthTypeEnum resamplingDepthType, + double targetDepthValue, + size_t firstIndex, + size_t secondIndex, + const std::map>& originalDepths, + const std::vector& propertyValues, + double eps ); + + // Made static due to unit testing + static std::pair, std::map>> + createResampledValuesAndDepths( RiaDefines::DepthTypeEnum resamplingDepthType, + const std::vector& targetDepths, + const std::map>& originalDepths, + const std::vector& propertyValues ); private: void calculateIntervalsOfContinousValidValues(); diff --git a/ApplicationLibCode/UnitTests/CMakeLists_files.cmake b/ApplicationLibCode/UnitTests/CMakeLists_files.cmake index 4af078d47a..7ec80e7d4d 100644 --- a/ApplicationLibCode/UnitTests/CMakeLists_files.cmake +++ b/ApplicationLibCode/UnitTests/CMakeLists_files.cmake @@ -93,6 +93,8 @@ set(SOURCE_GROUP_SOURCE_FILES ${CMAKE_CURRENT_LIST_DIR}/RigDeclineCurveCalculator-Test.cpp ${CMAKE_CURRENT_LIST_DIR}/RifReaderFmuRft-Test.cpp ${CMAKE_CURRENT_LIST_DIR}/RimSummaryRegressionAnalysisCurve-Test.cpp + ${CMAKE_CURRENT_LIST_DIR}/RimWellLogCalculatedCurve-Test.cpp + ${CMAKE_CURRENT_LIST_DIR}/RigWellLogCurveData-Test.cpp ) if(RESINSIGHT_ENABLE_GRPC) diff --git a/ApplicationLibCode/UnitTests/RigWellLogCurveData-Test.cpp b/ApplicationLibCode/UnitTests/RigWellLogCurveData-Test.cpp new file mode 100644 index 0000000000..68cad3f77b --- /dev/null +++ b/ApplicationLibCode/UnitTests/RigWellLogCurveData-Test.cpp @@ -0,0 +1,197 @@ +///////////////////////////////////////////////////////////////////////////////// +// +// Copyright (C) 2023- Equinor ASA +// +// ResInsight is free software: you can redistribute it and/or modify +// it under the terms of the GNU General Public License as published by +// the Free Software Foundation, either version 3 of the License, or +// (at your option) any later version. +// +// ResInsight is distributed in the hope that it will be useful, but WITHOUT ANY +// WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. +// +// See the GNU General Public License at +// for more details. +// +///////////////////////////////////////////////////////////////////////////////// + +#include "gtest/gtest.h" + +#include "RiaDefines.h" + +#include "RigWellLogCurveData.h" + +#include "cvfVector3.h" + +#include + +#include + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +TEST( RigWellLogCurveData, createAndAddInterpolatedSegmentValueAndDepths_first ) +{ + // Input data + const std::map> originalDepths = + { { RiaDefines::DepthTypeEnum::MEASURED_DEPTH, { 0.0, 20.0, 40.0 } }, + { RiaDefines::DepthTypeEnum::TRUE_VERTICAL_DEPTH, { 0.0, 30.0, 60.0 } }, + { RiaDefines::DepthTypeEnum::PSEUDO_LENGTH, { 0.0, 40.0, 80.0 } } }; + const std::vector propertyValues = { 0.0, 100.0, 150.0 }; + const double eps = 1e-6; + + // Output data + const auto resamplingDepthType = RiaDefines::DepthTypeEnum::MEASURED_DEPTH; + std::vector resampledValues; + std::map> resampledDepths; + + // Target data (resampling with MEASURED_DEPTH) + const double targetDepthValue = 10.0; // Halfway between index 0 and 1 for MEASURED_DEPTH in originalDepths + const size_t firstIndex = 0; + const size_t secondIndex = firstIndex + 1; + + // Call the function under test + RigWellLogCurveData::createAndAddInterpolatedSegmentValueAndDepths( resampledValues, + resampledDepths, + resamplingDepthType, + targetDepthValue, + firstIndex, + secondIndex, + originalDepths, + propertyValues, + eps ); + + // Check the results + ASSERT_EQ( resampledValues.size(), size_t( 1 ) ); + ASSERT_DOUBLE_EQ( resampledValues[0], 50.0 ); + + ASSERT_EQ( resampledDepths.size(), size_t( 3 ) ); + + ASSERT_EQ( resampledDepths[RiaDefines::DepthTypeEnum::MEASURED_DEPTH].size(), size_t( 1 ) ); + ASSERT_DOUBLE_EQ( resampledDepths[RiaDefines::DepthTypeEnum::MEASURED_DEPTH][0], 10.0 ); + + ASSERT_EQ( resampledDepths[RiaDefines::DepthTypeEnum::TRUE_VERTICAL_DEPTH].size(), size_t( 1 ) ); + ASSERT_DOUBLE_EQ( resampledDepths[RiaDefines::DepthTypeEnum::TRUE_VERTICAL_DEPTH][0], 15.0 ); + + ASSERT_EQ( resampledDepths[RiaDefines::DepthTypeEnum::PSEUDO_LENGTH].size(), size_t( 1 ) ); + ASSERT_DOUBLE_EQ( resampledDepths[RiaDefines::DepthTypeEnum::PSEUDO_LENGTH][0], 20.0 ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +TEST( RigWellLogCurveData, createAndAddInterpolatedSegmentValueAndDepths_second ) +{ + // Input data + const std::map> originalDepths = + { { RiaDefines::DepthTypeEnum::MEASURED_DEPTH, { 0.0, 20.0, 40.0 } }, + { RiaDefines::DepthTypeEnum::TRUE_VERTICAL_DEPTH, { 0.0, 30.0, 60.0 } }, + { RiaDefines::DepthTypeEnum::PSEUDO_LENGTH, { 0.0, 40.0, 80.0 } } }; + const std::vector propertyValues = { 0.0, 100.0, 150.0 }; + const double eps = 1e-6; + + // Output data + const auto resamplingDepthType = RiaDefines::DepthTypeEnum::MEASURED_DEPTH; + std::vector resampledValues; + std::map> resampledDepths; + + // Target data (resampling with MEASURED_DEPTH) + const double firstTargetDepthValue = 10.0; // Halfway between first and second index for MEASURED_DEPTH in originalDepths + const double secondTargetDepthValue = 30.0; // Halfway between second and third index for MEASURED_DEPTH in originalDepths + const size_t firstIndex = 0; + const size_t secondIndex = 1; + const size_t thirdIndex = 2; + + // Call the function under test with interpolating between first and second index + RigWellLogCurveData::createAndAddInterpolatedSegmentValueAndDepths( resampledValues, + resampledDepths, + resamplingDepthType, + firstTargetDepthValue, + firstIndex, + secondIndex, + originalDepths, + propertyValues, + eps ); + + // Call the function under test with interpolating between second and third index + RigWellLogCurveData::createAndAddInterpolatedSegmentValueAndDepths( resampledValues, + resampledDepths, + resamplingDepthType, + secondTargetDepthValue, + secondIndex, + thirdIndex, + originalDepths, + propertyValues, + eps ); + + // Check the results + ASSERT_EQ( resampledValues.size(), size_t( 2 ) ); + ASSERT_DOUBLE_EQ( resampledValues[0], 50.0 ); + ASSERT_DOUBLE_EQ( resampledValues[1], 125.0 ); + + ASSERT_EQ( resampledDepths.size(), size_t( 3 ) ); + + ASSERT_EQ( resampledDepths[RiaDefines::DepthTypeEnum::MEASURED_DEPTH].size(), size_t( 2 ) ); + ASSERT_DOUBLE_EQ( resampledDepths[RiaDefines::DepthTypeEnum::MEASURED_DEPTH][0], 10.0 ); + ASSERT_DOUBLE_EQ( resampledDepths[RiaDefines::DepthTypeEnum::MEASURED_DEPTH][1], 30.0 ); + + ASSERT_EQ( resampledDepths[RiaDefines::DepthTypeEnum::TRUE_VERTICAL_DEPTH].size(), size_t( 2 ) ); + ASSERT_DOUBLE_EQ( resampledDepths[RiaDefines::DepthTypeEnum::TRUE_VERTICAL_DEPTH][0], 15.0 ); + ASSERT_DOUBLE_EQ( resampledDepths[RiaDefines::DepthTypeEnum::TRUE_VERTICAL_DEPTH][1], 45.0 ); + + ASSERT_EQ( resampledDepths[RiaDefines::DepthTypeEnum::PSEUDO_LENGTH].size(), size_t( 2 ) ); + ASSERT_DOUBLE_EQ( resampledDepths[RiaDefines::DepthTypeEnum::PSEUDO_LENGTH][0], 20.0 ); + ASSERT_DOUBLE_EQ( resampledDepths[RiaDefines::DepthTypeEnum::PSEUDO_LENGTH][1], 60.0 ); +} + +//-------------------------------------------------------------------------------------------------- +/// +//-------------------------------------------------------------------------------------------------- +TEST( RigWellLogCurveData, CreateResampledValuesAndDepthsTest ) +{ + // Input data + RiaDefines::DepthTypeEnum resamplingDepthType = RiaDefines::DepthTypeEnum::MEASURED_DEPTH; + const std::vector targetDepths = { 0.0, 5.0, 10.0, 15.0 }; + const std::map> originalDepths = + { { RiaDefines::DepthTypeEnum::MEASURED_DEPTH, { 0.0, 10.0, 20.0 } }, + { RiaDefines::DepthTypeEnum::TRUE_VERTICAL_DEPTH, { 0.0, 20.0, 40.0 } }, + { RiaDefines::DepthTypeEnum::PSEUDO_LENGTH, { 0.0, 30.0, 60.0 } } }; + const std::vector propertyValues = { 0.0, 100.0, 200.0 }; + + // Call the function under test + auto result = RigWellLogCurveData::createResampledValuesAndDepths( resamplingDepthType, targetDepths, originalDepths, propertyValues ); + + // Check the results + std::vector& resampledPropertyValues = result.first; + std::map>& resampledDepths = result.second; + const auto expectedSize = targetDepths.size(); + + ASSERT_EQ( resampledDepths.size(), originalDepths.size() ); + ASSERT_EQ( resampledDepths[RiaDefines::DepthTypeEnum::MEASURED_DEPTH].size(), expectedSize ); + ASSERT_EQ( resampledDepths[RiaDefines::DepthTypeEnum::TRUE_VERTICAL_DEPTH].size(), expectedSize ); + ASSERT_EQ( resampledDepths[RiaDefines::DepthTypeEnum::PSEUDO_LENGTH].size(), expectedSize ); + + ASSERT_EQ( resampledPropertyValues.size(), expectedSize ); + + // Example assertions for the specific values + ASSERT_DOUBLE_EQ( resampledPropertyValues[0], 0.0 ); + ASSERT_DOUBLE_EQ( resampledPropertyValues[1], 50.0 ); + ASSERT_DOUBLE_EQ( resampledPropertyValues[2], 100.0 ); + ASSERT_DOUBLE_EQ( resampledPropertyValues[3], 150.0 ); + + ASSERT_DOUBLE_EQ( resampledDepths[RiaDefines::DepthTypeEnum::MEASURED_DEPTH][0], 0.0 ); + ASSERT_DOUBLE_EQ( resampledDepths[RiaDefines::DepthTypeEnum::MEASURED_DEPTH][1], 5.0 ); + ASSERT_DOUBLE_EQ( resampledDepths[RiaDefines::DepthTypeEnum::MEASURED_DEPTH][2], 10.0 ); + ASSERT_DOUBLE_EQ( resampledDepths[RiaDefines::DepthTypeEnum::MEASURED_DEPTH][3], 15.0 ); + + ASSERT_DOUBLE_EQ( resampledDepths[RiaDefines::DepthTypeEnum::TRUE_VERTICAL_DEPTH][0], 0.0 ); + ASSERT_DOUBLE_EQ( resampledDepths[RiaDefines::DepthTypeEnum::TRUE_VERTICAL_DEPTH][1], 10.0 ); + ASSERT_DOUBLE_EQ( resampledDepths[RiaDefines::DepthTypeEnum::TRUE_VERTICAL_DEPTH][2], 20.0 ); + ASSERT_DOUBLE_EQ( resampledDepths[RiaDefines::DepthTypeEnum::TRUE_VERTICAL_DEPTH][3], 30.0 ); + + ASSERT_DOUBLE_EQ( resampledDepths[RiaDefines::DepthTypeEnum::PSEUDO_LENGTH][0], 0.0 ); + ASSERT_DOUBLE_EQ( resampledDepths[RiaDefines::DepthTypeEnum::PSEUDO_LENGTH][1], 15.0 ); + ASSERT_DOUBLE_EQ( resampledDepths[RiaDefines::DepthTypeEnum::PSEUDO_LENGTH][2], 30.0 ); + ASSERT_DOUBLE_EQ( resampledDepths[RiaDefines::DepthTypeEnum::PSEUDO_LENGTH][3], 45.0 ); +}