Fault distance from point (#13367)

* Add python interface for distance to closest fault
* Split up some older code into multiple files.
* Add python test

Author: jonjenssen

ApplicationLibCode/ReservoirDataModel/RigFault.cpp

@@ -231,3 +231,28 @@ void RigFaultsPrCellAccumulator::setFaultIdx( size_t reservoirCellIndex, cvf::St
 {
     m_faultIdxForCellFace[reservoirCellIndex][face] = faultIdx;
 }
+
+//--------------------------------------------------------------------------------------------------
+///
+//--------------------------------------------------------------------------------------------------
+std::pair<double, cvf::StructGridInterface::FaceType> RigFault::minimumDistanceToPoint( const cvf::Vec3d& point, const RigMainGrid* mainGrid ) const
+{
+    double                             minDistance = std::numeric_limits<double>::max();
+    cvf::StructGridInterface::FaceType minFace     = cvf::StructGridInterface::NO_FACE;
+    for ( const FaultFace& ff : m_faultFaces )
+    {
+        const RigCell& cell = mainGrid->cell( ff.m_nativeReservoirCellIndex );
+        if ( cell.isInvalid() ) continue;
+
+        for ( auto& c : cell.faceCorners( ff.m_nativeFace ) )
+        {
+            double distance = ( c - point ).length();
+            if ( distance < minDistance )
+            {
+                minDistance = distance;
+                minFace     = ff.m_nativeFace;
+            }
+        }
+    }
+    return std::make_pair( minDistance, minFace );
+}

ApplicationLibCode/ReservoirDataModel/RigFault.h

@@ -97,6 +97,8 @@ class RigFault : public cvf::Object
 
     static bool ordering( CellAndFace first, CellAndFace second );
 
+    std::pair<double, cvf::StructGridInterface::FaceType> minimumDistanceToPoint( const cvf::Vec3d& point, const RigMainGrid* mainGrid ) const;
+
 private:
     QString m_name;
 

ApplicationLibCode/ReservoirDataModel/RigMainGrid.cpp

@@ -809,6 +809,28 @@ const RigFault* RigMainGrid::findFaultFromCellIndexAndCellFace( size_t reservoir
     return nullptr;
 }
 
+//--------------------------------------------------------------------------------------------------
+/// Returns the name of the closest fault, the distance to the fault, and the face type of the closest face
+//--------------------------------------------------------------------------------------------------
+std::tuple<QString, double, cvf::StructGridInterface::FaceType> RigMainGrid::minimumDistanceFaultToPoint( const cvf::Vec3d& point ) const
+{
+    double                             minDistance = std::numeric_limits<double>::max();
+    cvf::StructGridInterface::FaceType minFace     = cvf::StructGridInterface::FaceType::NO_FACE;
+    QString                            minFaultName;
+    for ( const auto& fault : m_faults )
+    {
+        auto [faultMinDistance, faultMinFace] = fault->minimumDistanceToPoint( point, this );
+
+        if ( faultMinDistance < minDistance )
+        {
+            minDistance  = faultMinDistance;
+            minFace      = faultMinFace;
+            minFaultName = fault->name();
+        }
+    }
+    return { minFaultName, minDistance, minFace };
+}
+
 //--------------------------------------------------------------------------------------------------
 ///
 //--------------------------------------------------------------------------------------------------

ApplicationLibCode/ReservoirDataModel/RigMainGrid.h

@@ -27,6 +27,8 @@
 #include "cvfBoundingBox.h"
 #include "cvfCollection.h"
 
+#include <QString>
+
 #include <vector>
 
 class RigActiveCellInfo;
@@ -116,6 +118,8 @@ class RigMainGrid : public RigGridBase
     // invalidate all cells with I > iLimit (0 based index)
     void invalidateCellsAboveI( size_t iLimit );
 
+    std::tuple<QString, double, cvf::StructGridInterface::FaceType> minimumDistanceFaultToPoint( const cvf::Vec3d& point ) const;
+
 protected: // only for use by file readers and internal services. TODO: replace with a better API
     friend class RigGridBase;
     friend class RigReservoirBuilder;
@@ -124,7 +128,7 @@ class RigMainGrid : public RigGridBase
     friend class RifReaderEclipseOutput;
     friend class RifReaderOpmCommon;
     friend class RiaGrpcCaseService;
-    friend class RiaActiveCellInfoStateHandler;
+    friend class RiaGrpcActiveCellInfoStateHandler;
     friend class RicCreateTemporaryLgrFeature;
     friend class RimCornerPointCase;
     std::vector<RigCell>&       reservoirCells();

GrpcInterface/CMakeLists.txt

@@ -39,6 +39,8 @@ set(SOURCE_GROUP_HEADER_FILES
     RiaGrpcApplicationInterface.h
     RiaGrpcWellPathService.h
     RiaWellPathDataToGrpcConverter.h
+    RiaGrpcActiveCellInfoStateHandler.h
+    RiaGrpcSelectedCellsStateHandler.h
     )
 
 set(SOURCE_GROUP_SOURCE_FILES
@@ -58,6 +60,8 @@ set(SOURCE_GROUP_SOURCE_FILES
     RiaGrpcApplicationInterface.cpp
     RiaGrpcWellPathService.cpp
     RiaWellPathDataToGrpcConverter.cpp
+    RiaGrpcActiveCellInfoStateHandler.cpp
+    RiaGrpcSelectedCellsStateHandler.cpp
 )
 
 # Find Protobuf installation Looks for protobuf-config.cmake file installed by

GrpcInterface/GrpcProtos/Case.proto

@@ -21,6 +21,7 @@ service Case {
   rpc GetCaseInfo(CaseRequest) returns (CaseInfo) {}
   rpc GetPdmObject(CaseRequest) returns (PdmObject) {}
   rpc GetReservoirBoundingBox(CaseRequest) returns (BoundingBox) {}
+  rpc GetDistanceToClosestFault(ClosestFaultRequest) returns (ClosestFault) {}
 }
 
 message CaseRequest { int32 id = 1; }
@@ -103,3 +104,14 @@ message SelectedCell {
 }
 
 message SelectedCells { repeated SelectedCell cells = 1; }
+
+message ClosestFaultRequest {
+  CaseRequest case_request = 1;
+  Vec3d point = 2;
+}
+
+message ClosestFault {
+  string fault_name = 1;
+  double distance = 2;
+  string face_name = 3;
+}

GrpcInterface/Python/rips/PythonExamples/case_and_grid_operations/fault_distance.py

@@ -0,0 +1,35 @@
+###################################################################################
+# This example prints the distance to and the name of the fault closest to a point
+###################################################################################
+
+import rips
+
+resinsight = rips.Instance.find()
+if resinsight is None:
+    exit(1)
+
+cases = resinsight.project.cases()
+if len(cases) == 0:
+    exit(1)
+
+case = cases[0]
+print("Using case: " + case.name)
+
+# random test point (positive Z for depth)
+point_x = 5039.84
+point_y = 6303.76
+point_z = 4144.21
+
+print("Looking for closest fault to point %f, %f, %f:" % (point_x, point_y, point_z))
+
+faultname, distance, facename = case.distance_to_closest_fault(
+    point_x, point_y, point_z
+)
+
+if facename == "":
+    print("- No fault found!")
+else:
+    print(
+        "- Distance to closest fault %s is %f, closest face direction is %s"
+        % (faultname, distance, facename)
+    )

GrpcInterface/Python/rips/case.py

@@ -43,6 +43,7 @@
 import Case_pb2_grpc
 import Commands_pb2 as Cmd
 import PdmObject_pb2 as PdmObject_pb2
+import Definitions_pb2
 
 import Properties_pb2
 import Properties_pb2_grpc
@@ -300,6 +301,17 @@ def reservoir_boundingbox(self):
     return self.__case_stub.GetReservoirBoundingBox(self.__request())
 
 
+@add_method(Case)
+def distance_to_closest_fault(self, x: float, y: float, z: float):
+    """Find the closest fault to the given point and return the distance, fault name and fault face"""
+    request = Case_pb2.ClosestFaultRequest(
+        case_request=self.__request(), point=Definitions_pb2.Vec3d(x=x, y=y, z=z)
+    )
+    reply = self.__case_stub.GetDistanceToClosestFault(request)
+
+    return (reply.fault_name, reply.distance, reply.face_name)
+
+
 @add_method(Case)
 def reservoir_depth_range(self) -> Tuple[float, float]:
     """Get the reservoir depth range

GrpcInterface/Python/rips/tests/test_faults.py

@@ -0,0 +1,30 @@
+import sys
+import os
+import math
+import pytest
+
+sys.path.insert(1, os.path.join(sys.path[0], "../../"))
+import rips
+
+import dataroot
+
+
+def test_faultDistance(rips_instance, initialize_test):
+    case = rips_instance.project.load_case(
+        dataroot.PATH + "/TEST10K_FLT_LGR_NNC/TEST10K_FLT_LGR_NNC.EGRID"
+    )
+
+    # a test point
+    point_x = 4817.84
+    point_y = 5204.76
+    point_z = 4137.21
+
+    faultname, distance, facename = case.distance_to_closest_fault(
+        point_x, point_y, point_z
+    )
+
+    # Fault name is unstable between grid readers, so we skip this check
+    # assert faultname == "Undefined Grid Faults"
+
+    assert facename == "I+"
+    assert math.isclose(distance, 53.84, abs_tol=0.1)