#2029 Well Log Plot: Handle curves with different depth in 'Show Plot Data'

ApplicationLibCode/Application/Tools/CMakeLists_files.cmake

@@ -28,6 +28,7 @@ set(SOURCE_GROUP_HEADER_FILES
     ${CMAKE_CURRENT_LIST_DIR}/RiaCurveMerger.h
     ${CMAKE_CURRENT_LIST_DIR}/RiaCurveMerger.inl
     ${CMAKE_CURRENT_LIST_DIR}/RiaCurveDataTools.h
+    ${CMAKE_CURRENT_LIST_DIR}/RiaWellLogCurveMerger.h
     ${CMAKE_CURRENT_LIST_DIR}/RiaTimeHistoryCurveResampler.h
     ${CMAKE_CURRENT_LIST_DIR}/RiaStatisticsTools.h
     ${CMAKE_CURRENT_LIST_DIR}/RiaOffshoreSphericalCoords.h
@@ -77,6 +78,7 @@ set(SOURCE_GROUP_SOURCE_FILES
     ${CMAKE_CURRENT_LIST_DIR}/RiaFilePathTools.cpp
     ${CMAKE_CURRENT_LIST_DIR}/RiaCurveMerger.cpp
     ${CMAKE_CURRENT_LIST_DIR}/RiaCurveDataTools.cpp
+    ${CMAKE_CURRENT_LIST_DIR}/RiaWellLogCurveMerger.cpp
     ${CMAKE_CURRENT_LIST_DIR}/RiaTimeHistoryCurveResampler.cpp
     ${CMAKE_CURRENT_LIST_DIR}/RiaStatisticsTools.cpp
     ${CMAKE_CURRENT_LIST_DIR}/RiaWeightedGeometricMeanCalculator.cpp

ApplicationLibCode/Application/Tools/RiaCurveMerger.h

@@ -57,6 +57,9 @@ public:
     // See ExpressionParserImpl::assignVector()
     std::vector<double>& interpolatedYValuesForAllXValues( size_t curveIdx );
 
+    static void removeValuesForPartialCurves( std::set<XValueType, XComparator>&                    unionOfXValues,
+                                              const std::vector<std::pair<XValueType, XValueType>>& originalXBounds );
+
 public:
     // Helper methods, available as public to be able to access from unit tests
 

ApplicationLibCode/Application/Tools/RiaCurveMerger.inl

@@ -194,7 +194,17 @@ void RiaCurveMerger<XValueType>::computeUnionOfXValues( bool includeValuesForPar
         originalXBounds.push_back( std::make_pair( *( minmax_it.first ), *( minmax_it.second ) ) );
     }
 
-    if ( !includeValuesForPartialCurves )
+    if ( !includeValuesForPartialCurves ) removeValuesForPartialCurves( unionOfXValues, originalXBounds );
+
+    m_allXValues = std::vector<XValueType>( unionOfXValues.begin(), unionOfXValues.end() );
+}
+
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+template <typename XValueType>
+void RiaCurveMerger<XValueType>::removeValuesForPartialCurves( std::set<XValueType, XComparator>& unionOfXValues,
+                                                               const std::vector<std::pair<XValueType, XValueType>>& originalXBounds )
 {
     for ( auto it = unionOfXValues.begin(); it != unionOfXValues.end(); )
     {
@@ -218,9 +228,6 @@ void RiaCurveMerger<XValueType>::computeUnionOfXValues( bool includeValuesForPar
     }
 }
 
-    m_allXValues = std::vector<XValueType>( unionOfXValues.begin(), unionOfXValues.end() );
-}
-
 //--------------------------------------------------------------------------------------------------
 ///
 //--------------------------------------------------------------------------------------------------

ApplicationLibCode/Application/Tools/RiaWellLogCurveMerger.cpp

@@ -0,0 +1,195 @@
+/////////////////////////////////////////////////////////////////////////////////
+//
+//  Copyright (C) 2021     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 <http://www.gnu.org/licenses/gpl.html>
+//  for more details.
+//
+/////////////////////////////////////////////////////////////////////////////////
+
+#include "RiaWellLogCurveMerger.h"
+
+#include "RiaCurveMerger.h"
+
+#include <algorithm>
+#include <cmath>
+
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+RiaWellLogCurveMerger::RiaWellLogCurveMerger()
+{
+}
+
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+void RiaWellLogCurveMerger::addCurveData( const std::vector<double>& xValues, const std::vector<double>& yValues )
+{
+    CVF_ASSERT( xValues.size() == yValues.size() );
+
+    if ( !xValues.empty() )
+    {
+        m_originalValues.push_back( std::make_pair( xValues, yValues ) );
+    }
+}
+
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+size_t RiaWellLogCurveMerger::curveCount() const
+{
+    return m_lookupValuesForAllCurves.size();
+}
+
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+RiaCurveDataTools::CurveIntervals RiaWellLogCurveMerger::validIntervalsForAllXValues() const
+{
+    return m_validIntervalsForAllXValues;
+}
+
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+const std::vector<double>& RiaWellLogCurveMerger::allXValues() const
+{
+    return m_allXValues;
+}
+
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+const std::vector<double>& RiaWellLogCurveMerger::lookupYValuesForAllXValues( size_t curveIdx ) const
+{
+    CVF_ASSERT( curveIdx < m_lookupValuesForAllCurves.size() );
+
+    return m_lookupValuesForAllCurves[curveIdx];
+}
+
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+void RiaWellLogCurveMerger::computeLookupValues( bool includeValuesFromPartialCurves )
+{
+    m_validIntervalsForAllXValues.clear();
+    m_allXValues.clear();
+    m_lookupValuesForAllCurves.clear();
+
+    computeUnionOfXValues( includeValuesFromPartialCurves );
+
+    const size_t curveCount = m_originalValues.size();
+    if ( curveCount == 0 )
+    {
+        return;
+    }
+
+    const size_t dataValueCount = m_allXValues.size();
+    if ( dataValueCount == 0 )
+    {
+        return;
+    }
+
+    m_lookupValuesForAllCurves.resize( curveCount );
+
+    std::vector<double> accumulatedValidValues( dataValueCount, 1.0 );
+
+    for ( size_t curveIdx = 0; curveIdx < curveCount; curveIdx++ )
+    {
+        std::vector<double>& curveValues = m_lookupValuesForAllCurves[curveIdx];
+        curveValues.resize( dataValueCount );
+
+        for ( size_t valueIndex = 0; valueIndex < dataValueCount; valueIndex++ )
+        {
+            double interpolValue = lookupYValue( m_allXValues[valueIndex],
+                                                 m_originalValues[curveIdx].first,
+                                                 m_originalValues[curveIdx].second );
+            if ( !RiaCurveDataTools::isValidValue( interpolValue, false ) )
+            {
+                accumulatedValidValues[valueIndex] = HUGE_VAL;
+            }
+
+            curveValues[valueIndex] = interpolValue;
+        }
+    }
+
+    m_validIntervalsForAllXValues = RiaCurveDataTools::calculateIntervalsOfValidValues( accumulatedValidValues, false );
+}
+
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+
+void RiaWellLogCurveMerger::computeUnionOfXValues( bool includeValuesForPartialCurves )
+{
+    m_allXValues.clear();
+
+    std::set<double, XValueComparator<double>> unionOfXValues;
+
+    std::vector<std::pair<double, double>> originalXBounds;
+    for ( const auto& curveData : m_originalValues )
+    {
+        if ( curveData.first.empty() ) continue;
+
+        // Well log data has top and bottom depth for each zone
+        // Find the mid point in the zone to
+        const std::vector<double>& xValues = curveData.first;
+        for ( size_t i = 0; i < xValues.size(); i += 2 )
+        {
+            if ( i + 1 < xValues.size() )
+            {
+                double top    = xValues.at( i );
+                double bottom = xValues.at( i + 1 );
+                double mid    = ( top + bottom ) / 2.0;
+                unionOfXValues.insert( mid );
+            }
+        }
+        auto minmax_it = std::minmax_element( curveData.first.begin(), curveData.first.end() );
+        originalXBounds.push_back( std::make_pair( *( minmax_it.first ), *( minmax_it.second ) ) );
+    }
+
+    if ( !includeValuesForPartialCurves )
+        RiaCurveMerger<double>::removeValuesForPartialCurves( unionOfXValues, originalXBounds );
+
+    m_allXValues = std::vector<double>( unionOfXValues.begin(), unionOfXValues.end() );
+}
+
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+double RiaWellLogCurveMerger::lookupYValue( const double&              interpolationXValue,
+                                            const std::vector<double>& xValues,
+                                            const std::vector<double>& yValues )
+{
+    if ( yValues.size() != xValues.size() ) return HUGE_VAL;
+
+    const bool removeInterpolatedValues = false;
+
+    if ( interpolationXValue < xValues[0] ) return HUGE_VAL;
+
+    for ( size_t firstI = 0; firstI < xValues.size(); firstI++ )
+    {
+        if ( xValues.at( firstI ) >= interpolationXValue )
+        {
+            double firstValue = yValues.at( firstI );
+            if ( !RiaCurveDataTools::isValidValue( firstValue, removeInterpolatedValues ) )
+            {
+                return HUGE_VAL;
+            }
+
+            return firstValue;
+        }
+    }
+
+    return HUGE_VAL;
+}

ApplicationLibCode/Application/Tools/RiaWellLogCurveMerger.h

@@ -0,0 +1,55 @@
+/////////////////////////////////////////////////////////////////////////////////
+//
+//  Copyright (C) 2021     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 <http://www.gnu.org/licenses/gpl.html>
+//  for more details.
+//
+/////////////////////////////////////////////////////////////////////////////////
+
+#pragma once
+
+#include "RiaCurveDataTools.h"
+
+#include <ctime>
+
+//==================================================================================================
+///
+//==================================================================================================
+class RiaWellLogCurveMerger
+{
+public:
+    RiaWellLogCurveMerger();
+
+    void   addCurveData( const std::vector<double>& xValues, const std::vector<double>& yValues );
+    size_t curveCount() const;
+
+    void computeLookupValues( bool includeValuesFromPartialCurves = true );
+
+    RiaCurveDataTools::CurveIntervals validIntervalsForAllXValues() const;
+    const std::vector<double>&        allXValues() const;
+    const std::vector<double>&        lookupYValuesForAllXValues( size_t curveIdx ) const;
+
+private:
+    void computeUnionOfXValues( bool includeValuesFromPartialCurves );
+
+    static double lookupYValue( const double&              xValue,
+                                const std::vector<double>& curveXValues,
+                                const std::vector<double>& curveYValues );
+
+    std::vector<std::pair<std::vector<double>, std::vector<double>>> m_originalValues;
+
+    RiaCurveDataTools::CurveIntervals m_validIntervalsForAllXValues;
+
+    std::vector<double>              m_allXValues;
+    std::vector<std::vector<double>> m_lookupValuesForAllCurves;
+};

ApplicationLibCode/ProjectDataModel/WellLog/RimWellLogTrack.cpp

@@ -25,6 +25,7 @@
 #include "RiaLogging.h"
 #include "RiaPreferences.h"
 #include "RiaSimWellBranchTools.h"
+#include "RiaWellLogCurveMerger.h"
 
 #include "RigEclipseCaseData.h"
 #include "RigEclipseResultAddress.h"
@@ -844,6 +845,9 @@ QString RimWellLogTrack::asciiDataForPlotExport() const
         }
     }
 
+    RiaWellLogCurveMerger curveMerger;
+    bool                  foundNonMatchingDepths = false;
+
     for ( RimWellLogCurve* curve : m_curves() )
     {
         if ( !curve->isCurveVisible() ) continue;
@@ -858,7 +862,7 @@ QString RimWellLogTrack::asciiDataForPlotExport() const
         }
 
         std::vector<double> xPlotValues = curveData->xPlotValues();
-        if ( curveDepths.size() != xPlotValues.size() || xPlotValues.empty() )
+        if ( xPlotValues.empty() )
         {
             curveNames.pop_back();
 
@@ -868,6 +872,15 @@ QString RimWellLogTrack::asciiDataForPlotExport() const
             }
             continue;
         }
+
+        if ( curveDepths.size() != xPlotValues.size() )
+        {
+            foundNonMatchingDepths = true;
+        }
+
+        std::vector<double> depths = curveData->depthPlotValues( depthType, depthUnit );
+        curveMerger.addCurveData( depths, xPlotValues );
+
         curvesPlotXValues.push_back( xPlotValues );
     }
 
@@ -900,6 +913,26 @@ QString RimWellLogTrack::asciiDataForPlotExport() const
     }
     out += "\n";
 
+    // Resample when curves have different depth
+    if ( foundNonMatchingDepths )
+    {
+        curvesPlotXValues.clear();
+        curveDepths.clear();
+
+        curveMerger.computeLookupValues();
+
+        const std::vector<double>& allDepths = curveMerger.allXValues();
+        curveDepths                          = allDepths;
+        for ( size_t depthIdx = 0; depthIdx < allDepths.size(); depthIdx++ )
+        {
+            for ( size_t curveIdx = 0; curveIdx < curveMerger.curveCount(); ++curveIdx )
+            {
+                const std::vector<double>& curveValues = curveMerger.lookupYValuesForAllXValues( curveIdx );
+                curvesPlotXValues.push_back( curveValues );
+            }
+        }
+    }
+
     for ( size_t dIdx = 0; dIdx < curveDepths.size(); ++dIdx )
     {
         size_t i          = dIdx;