Add EXPORT_RFT hook workflow

Author: HakonSohoel

docs/ert/reference/configuration/data_types.rst

@@ -519,6 +519,37 @@ For loading multiple RFT observations from a CSV file, see the
 :ref:`RFT_OBSERVATION <rft_observation>` documentation.
 
 
+Exporting RFT data for visualization
+~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+The `EXPORT_RFT <../workflows/added_workflow_jobs.html#EXPORT_RFT>`_ workflow job exports RFT observations and simulated
+responses to CSV files compatible with the webviz-subsurface RftPlotter plugin.
+This workflow job is a replacement for the `MERGE_RFT_ERTOBS` forward model when
+using :ref:`RFT_OBSERVATION <rft_observation>` instead of the legacy `GENDATA_RFT` method.
+
+To use it, create a workflow file e.g. ``export_rft_data`` with content::
+
+   EXPORT_RFT
+
+Then add the workflow to your ERT configuration::
+
+   LOAD_WORKFLOW export_rft_data
+   HOOK_WORKFLOW export_rft_data POST_SIMULATION
+
+By default, the output file is written to
+``share/results/tables/rft_ert.csv`` in each realization's runpath.
+A custom filename can be specified as a parameter::
+
+   EXPORT_RFT custom_rft.csv
+
+.. note::
+
+   To include SGAS, SWAT, and SOIL saturation values in the exported CSV,
+   SGAS and SWAT properties must be listed in the :ref:`RFT <rft>` response
+   configuration no make ERT extract these properties from the RFT files. For example::
+
+      RFT WELL:PROD DATE:2015-02-01 PROPERTIES:PRESSURE,SWAT,SGAS
+
 General data: ``GEN_DATA``
 --------------------------
 

src/ert/config/_create_observation_dataframes.py

@@ -179,6 +179,8 @@ def _handle_rft_observation(
                 f"{rft_observation.date}:"
                 f"{rft_observation.property}"
             ),
+            "well": rft_observation.well,
+            "date": rft_observation.date,
             "observation_key": rft_observation.name,
             "east": pl.Series([location[0]], dtype=pl.Float32),
             "north": pl.Series([location[1]], dtype=pl.Float32),

src/ert/config/rft_config.py

@@ -285,6 +285,9 @@ def _props_matcher(props: list[str]) -> str:
             return pl.DataFrame(
                 {
                     "response_key": [],
+                    "well": [],
+                    "date": [],
+                    "property": [],
                     "time": [],
                     "depth": [],
                     "values": [],
@@ -304,6 +307,9 @@ def _props_matcher(props: list[str]) -> str:
                     pl.DataFrame(
                         {
                             "response_key": [f"{well}:{time.isoformat()}:{prop}"],
+                            "well": [well],
+                            "date": [time.isoformat()],
+                            "property": [prop],
                             "time": [time],
                             "depth": [fetched[well, time]["DEPTH"]],
                             "values": [vals],

src/ert/dark_storage/endpoints/responses.py

@@ -182,7 +182,7 @@ def data_for_gradient(ensemble: Ensemble, key: str) -> pd.DataFrame:
 response_to_pandas_x_axis_fns: dict[str, Callable[[tuple[Any, ...]], Any]] = {
     "summary": lambda t: pd.Timestamp(t[2]).isoformat(),
     "gen_data": lambda t: str(t[3]),
-    "rft": lambda t: str(t[4]),
+    "rft": lambda t: str(t[6]),
 }
 
 

src/ert/plugins/hook_implementations/workflows/__init__.py

@@ -7,6 +7,7 @@
 from .csv_export import CSVExportJob
 from .disable_parameters import DisableParametersUpdate
 from .export_misfit_data import ExportMisfitDataJob
+from .export_rft import ExportRFTJob
 from .export_runpath import ExportRunpathJob
 from .gen_data_rft_export import GenDataRFTCSVExportJob
 from .misfit_preprocessor import MisfitPreprocessor
@@ -21,6 +22,7 @@ def ertscript_workflow(config: WorkflowConfigs) -> None:
         ExportMisfitDataJob, "EXPORT_MISFIT_DATA", category="observations.correlation"
     )
     config.add_workflow(ExportRunpathJob, "EXPORT_RUNPATH")
+    config.add_workflow(ExportRFTJob, "EXPORT_RFT")
     config.add_workflow(DisableParametersUpdate, "DISABLE_PARAMETERS")
     config.add_workflow(
         MisfitPreprocessor, "MISFIT_PREPROCESSOR", category="observations.correlation"

src/ert/plugins/hook_implementations/workflows/export_rft.py

@@ -0,0 +1,51 @@
+from __future__ import annotations
+
+from pathlib import Path
+from typing import TYPE_CHECKING, Any
+
+from ert import ErtScript
+
+if TYPE_CHECKING:
+    from ert.runpaths import Runpaths
+    from ert.storage import Ensemble
+
+
+class ExportRFTJob(ErtScript):
+    """
+    Export RFT observations and simulated responses to CSV files.
+
+    By default, the output file is "share/results/tables/rft_ert.csv" written to
+    each realization's runpath. The filename can be overridden by giving it as
+    the first parameter:
+
+        EXPORT_RFT custom_filename.csv
+
+    Each realization gets its own file containing observation data joined with
+    simulated response values for that realization. The output is compatible
+    with the webviz-subsurface RftPlotter.
+    """
+
+    def run(
+        self,
+        run_paths: Runpaths,
+        ensemble: Ensemble,
+        workflow_args: list[Any],
+    ) -> None:
+        filename = (
+            workflow_args[0] if workflow_args else "share/results/tables/rft_ert.csv"
+        )
+
+        observations_and_responses = ensemble.get_rft_observations_and_responses()
+
+        iteration = ensemble.iteration
+        realizations = ensemble.get_realization_list_with_responses()
+        paths = run_paths.get_paths(realizations, iteration)
+
+        for realization, runpath in zip(realizations, paths, strict=True):
+            realization_data = observations_and_responses.filter(
+                observations_and_responses["realization"] == realization
+            ).drop("realization")
+
+            target_file = Path(runpath) / filename
+            target_file.parent.mkdir(parents=True, exist_ok=True)
+            realization_data.write_csv(target_file)

src/ert/storage/local_ensemble.py

@@ -27,6 +27,7 @@
     ParameterConfig,
     SummaryConfig,
 )
+from ert.exceptions import StorageError
 from ert.substitutions import substitute_runpath_name
 
 from .load_status import LoadResult
@@ -1086,6 +1087,149 @@ def get_observations_and_responses(
                 pl.col("response_key").cast(pl.String).alias("response_key")
             )
 
+    def get_rft_observations_and_responses(
+        self,
+    ) -> pl.DataFrame:
+        """Fetches and aligns RFT observations with their corresponding
+        simulated responses from an ensemble.
+
+        Returns a DataFrame with observation/response data using
+        column names equal to the ones used by the subscript forward model
+        MERGE_RFT_ERTOBS, and compatible with the webviz-subsurface RftPlotter.
+        """
+        rft_observations = self.experiment.observations.get("rft")
+        if rft_observations is None or rft_observations.is_empty():
+            raise StorageError("No RFT observations found in experiment")
+
+        if "rft" not in self.experiment.response_configuration:
+            raise KeyError("No RFT response configuration found in experiment")
+
+        realizations = self.get_realization_list_with_responses()
+        if not realizations:
+            raise StorageError("No realizations with responses found")
+
+        # Build date-to-report_step mapping from summary responses if available
+        date_to_report_step: dict[str, int] = {}
+        if "summary" in self.experiment.response_configuration:
+            try:
+                summary_df = self.load_responses("summary", (realizations[0],))
+                times = summary_df["time"].unique().sort()
+                for report_step, time in enumerate(times):
+                    date_str = time.strftime("%Y-%m-%d")
+                    date_to_report_step[date_str] = report_step
+            except (KeyError, IndexError):
+                pass  # No summary data available, will use default
+
+        observations = rft_observations.with_columns(
+            pl.int_range(pl.len()).over("well").alias("order"),
+        )
+
+        join_keys = ["well", "date", "east", "north", "tvd"]
+        observed_values = {k: observations[k].unique() for k in join_keys}
+
+        pivot_index = [
+            "well",
+            "date",
+            "realization",
+            "response_zone",
+            "east",
+            "north",
+            "tvd",
+            "i",
+            "j",
+            "k",
+        ]
+
+        output_columns = [
+            "order",
+            "east",
+            "north",
+            "md",
+            "tvd",
+            "zone",
+            "pressure",
+            "swat",
+            "sgas",
+            "soil",
+            "valid_zone",
+            "is_active",
+            "i",
+            "j",
+            "k",
+            "well",
+            "date",
+            "realization",
+            "report_step",
+            "observations",
+            "std",
+        ]
+
+        result_frames: list[pl.DataFrame] = []
+
+        for real in sorted(realizations):
+            responses = (
+                self.load_responses("rft", (real,))
+                .with_columns(
+                    pl.col("property").str.to_lowercase(),
+                )
+                .rename({"zone": "response_zone"})
+            )
+
+            for col, values in observed_values.items():
+                responses = responses.filter(
+                    pl.col(col).is_in(values.implode(), nulls_equal=True)
+                )
+
+            pivoted = responses.pivot(
+                on="property",
+                index=pivot_index,
+                values="values",
+            )
+
+            for col in ["pressure", "sgas", "swat"]:
+                if col not in pivoted.columns:
+                    pivoted = pivoted.with_columns(
+                        pl.lit(None).cast(pl.Float32).alias(col)
+                    )
+
+            pivoted = pivoted.with_columns(
+                pl.col("pressure").is_not_null().alias("is_active")
+            )
+
+            joined = (
+                observations.join(
+                    pivoted,
+                    how="left",
+                    on=join_keys,
+                    nulls_equal=True,
+                )
+                .with_columns(
+                    pl.col("zone")
+                    .eq_missing(pl.col("response_zone"))
+                    .alias("valid_zone"),
+                    (1 - pl.col("sgas") - pl.col("swat")).alias("soil"),
+                    pl.col("is_active").fill_null(False),
+                    pl.col("date")
+                    .replace_strict(date_to_report_step, default=0)
+                    .alias("report_step"),
+                )
+                .select(output_columns)
+            )
+
+            result_frames.append(joined)
+
+        return pl.concat(result_frames, how="vertical").rename(
+            {
+                "east": "utm_x",
+                "north": "utm_y",
+                "md": "measured_depth",
+                "tvd": "true_vertical_depth",
+                "date": "time",
+                "observations": "observed",
+                "std": "error",
+            }
+        )
+
     @property
     def everest_realization_info(self) -> dict[int, EverestRealizationInfo] | None:
         return self._index.everest_realization_info

test-data/ert/rft_example/export_rft_data

@@ -0,0 +1 @@
+EXPORT_RFT

test-data/ert/rft_example/rft.ert

@@ -21,4 +21,7 @@ FORWARD_MODEL TEMPLATE_RENDER( \
 
 FORWARD_MODEL FLOW
 
-RFT WELL:PROD DATE:2015-02-01 PROPERTIES:PRESSURE
+RFT WELL:PROD DATE:2015-02-01 PROPERTIES:PRESSURE,SWAT,SGAS
+
+LOAD_WORKFLOW export_rft_data
+HOOK_WORKFLOW export_rft_data POST_SIMULATION

tests/ert/unit_tests/config/test_observations.py

@@ -328,6 +328,8 @@ def test_that_rft_config_is_created_from_observations():
         pl.DataFrame(
             {
                 "response_key": "well:2013-03-31:PRESSURE",
+                "well": "well",
+                "date": "2013-03-31",
                 "observation_key": "NAME",
                 "east": pl.Series([30.0], dtype=pl.Float32),
                 "north": pl.Series([71.0], dtype=pl.Float32),