Merge pull request #61 from fema-ffrd/60-observed-flow-and-stage-data

60 observed flow and stage data

.devcontainer/Dockerfile

@@ -0,0 +1 @@
+FROM mcr.microsoft.com/devcontainers/base:jammy

.devcontainer/devcontainer.json

@@ -0,0 +1,27 @@
+{
+	"name": "rashdf-devcontainer",
+	"build": {
+		"dockerfile": "Dockerfile",
+		"context": ".."
+	},
+	"features": {
+		"ghcr.io/devcontainers/features/git:1": {},
+		"ghcr.io/devcontainers/features/python:1":{"version":"3.12"}
+	},
+  "mounts": [],
+	"customizations": {
+		"vscode": {
+			"extensions": [
+				"ms-python.python",
+				"charliermarsh.ruff",
+				"GitHub.copilot"
+			],
+			"settings": {
+				"python.defaultInterpreterPath": "/opt/conda/envs/rashdf/bin/python"
+			}
+		}
+	},
+	// avoid dubious ownership of the workspace folder https://www.kenmuse.com/blog/avoiding-dubious-ownership-in-dev-containers/
+	"postStartCommand": "git config --global --add safe.directory ${containerWorkspaceFolder}",
+	"postCreateCommand": "sudo chown -R vscode:vscode ${containerWorkspaceFolder}"
+}

README.md

@@ -110,8 +110,12 @@ $ python -m venv venv-rashdf
 
 Activate the virtual environment:
 ```
-$ source ./venv/bin/activate
+# For macOS/Linux
+$ source ./venv-rashdf/bin/activate
 (venv-rashdf) $
+
+# For Windows
+> ./venv-rashdf/Scripts/activate
 ```
 
 Install dev dependencies:

pyproject.toml

@@ -16,7 +16,7 @@ version = "0.6.0"
 dependencies = ["h5py", "geopandas>=1.0,<2.0", "pyarrow", "xarray"]
 
 [project.optional-dependencies]
-dev = ["pre-commit", "ruff", "pytest", "pytest-cov", "fiona", "kerchunk", "zarr", "dask", "fsspec", "s3fs"]
+dev = ["pre-commit", "ruff", "pytest", "pytest-cov", "kerchunk", "zarr", "dask", "fsspec", "s3fs", "fiona==1.9.6"]
 docs = ["sphinx", "numpydoc", "sphinx_rtd_theme"]
 
 [project.urls]

src/rashdf/plan.py

@@ -6,6 +6,7 @@
     ras_timesteps_to_datetimes,
     parse_ras_datetime_ms,
     deprecated,
+    convert_ras_hdf_value,
 )
 
 from geopandas import GeoDataFrame
@@ -156,6 +157,7 @@ class RasPlanHdf(RasGeomHdf):
     PLAN_INFO_PATH = "Plan Data/Plan Information"
     PLAN_PARAMS_PATH = "Plan Data/Plan Parameters"
     PRECIP_PATH = "Event Conditions/Meteorology/Precipitation"
+    OBS_DATA_PATH = "Event Conditions/Observed Data"
     RESULTS_UNSTEADY_PATH = "Results/Unsteady"
     RESULTS_UNSTEADY_SUMMARY_PATH = f"{RESULTS_UNSTEADY_PATH}/Summary"
     VOLUME_ACCOUNTING_PATH = f"{RESULTS_UNSTEADY_PATH}/Volume Accounting"
@@ -166,6 +168,8 @@ class RasPlanHdf(RasGeomHdf):
     UNSTEADY_TIME_SERIES_PATH = f"{BASE_OUTPUT_PATH}/Unsteady Time Series"
     REFERENCE_LINES_OUTPUT_PATH = f"{UNSTEADY_TIME_SERIES_PATH}/Reference Lines"
     REFERENCE_POINTS_OUTPUT_PATH = f"{UNSTEADY_TIME_SERIES_PATH}/Reference Points"
+    OBS_FLOW_OUTPUT_PATH = f"{OBS_DATA_PATH}/Flow"
+    OBS_STAGE_OUTPUT_PATH = f"{OBS_DATA_PATH}/Stage"
 
     RESULTS_STEADY_PATH = "Results/Steady"
     BASE_STEADY_PATH = f"{RESULTS_STEADY_PATH}/Output/Output Blocks/Base Output"
@@ -1117,6 +1121,74 @@ def reference_lines_timeseries_output(self) -> xr.Dataset:
         """
         return self.reference_timeseries_output(reftype="lines")
 
+    def observed_timeseries_input(self, vartype: str = "Flow") -> dict:
+        """Return observed timeseries input data for reference lines and points from a HEC-RAS HDF plan file.
+
+        Parameters
+        ----------
+        vartype : str, optional
+            The type of observed data to retrieve. Must be either "Flow" or "Stage".
+            (default: "Flow")
+
+        Returns
+        -------
+        xr.Dataset
+            An xarray Dataset with observed timeseries input data for both reference lines and reference points.
+        """
+        if vartype == "Flow":
+            output_path = self.OBS_FLOW_OUTPUT_PATH
+        elif vartype == "Stage":
+            output_path = self.OBS_STAGE_OUTPUT_PATH
+        else:
+            raise ValueError('vartype must be either "Flow" or "Stage".')
+
+        observed_group = self.get(output_path)
+        if observed_group is None:
+            raise RasPlanHdfError(
+                f"Could not find HDF group at path '{output_path}'."
+                f" Does the Plan HDF file contain reference {vartype} output data?"
+            )
+        if "Attributes" in observed_group.keys():
+            attr_path = observed_group["Attributes"]
+            attrs_df = pd.DataFrame(attr_path[:]).map(convert_ras_hdf_value)
+
+        das = {}
+        for idx, site in enumerate(observed_group.keys()):
+            if site != "Attributes":
+                # Site Ex: 'Ref Point: Grapevine_Lake_RP'
+                site_path = observed_group[site]
+                site_name = site.split(":")[1][1:]  # Grapevine_Lake_RP
+                ref_type = site.split(":")[0]  # Ref Point
+                if ref_type == "Ref Line":
+                    ref_type = "refln"
+                else:
+                    ref_type = "refpt"
+                df = pd.DataFrame(site_path[:]).map(convert_ras_hdf_value)
+                # rename Date to time
+                df = df.rename(columns={"Date": "time"})
+                # Ensure the Date index is unique
+                df = df.drop_duplicates(subset="time")
+                # Package into an 1D xarray DataArray
+                values = df["Value"].values
+                times = df["time"].values
+                da = xr.DataArray(
+                    values,
+                    name=vartype,
+                    dims=["time"],
+                    coords={
+                        "time": times,
+                    },
+                    attrs={
+                        "hdf_path": f"{output_path}/{site}",
+                    },
+                )
+                # Expand dimensions to add additional coordinates
+                da = da.expand_dims({f"{ref_type}_id": [idx - 1]})
+                da = da.expand_dims({f"{ref_type}_name": [site_name]})
+                das[site_name] = da
+        das = xr.concat([das[site] for site in das.keys()], dim="time")
+        return das
+
     def reference_points_timeseries_output(self) -> xr.Dataset:
         """Return timeseries output data for reference points from a HEC-RAS HDF plan file.