Gridded files (#264)

* stream from from get_gridded_data()

* stream for get_gridded_data

* add optional file_path arg to get_gridded_data_from_api

* implement get_gridded_file() running locally

* support remote get_gridded_file()

* update docstrings

* gridded_file() for remote

* deprecate get_gridded_files() with a warning

* fix unit tests

* change test to get_gridded_file() and pass from option when getting huc maps from internal function

* update Docstring comments

* fix error message when no network connection

* add warning for deprecation

Author: wh3248

pyproject.toml

@@ -1,6 +1,6 @@
 [tool.poetry]
 name = "hf_hydrodata"
-version = "1.4.6"
+version = "1.4.7"
 description = "hydroframe tools and utilities"
 authors = ["William M. Hasling", "Laura Condon", "Reed Maxwell",  "George Artavanis", "Amy M. Johnson", "Amy C. Defnet"]
 license = "MIT"

src/hf_hydrodata/__init__.py

@@ -23,6 +23,7 @@
 from hf_hydrodata.gridded import (
     get_gridded_data,
     get_gridded_files,
+    get_gridded_file,
     get_raw_file,
     get_date_range,
     get_huc_from_latlon,

src/hf_hydrodata/data_model_access.py

@@ -129,9 +129,13 @@ def _query_data_catalog(self, options: dict):
             response = requests.get(url, timeout=120, headers=headers)
             if response.status_code == 200:
                 response_json = json.loads(response.text)
+            elif response.status_code == 403:
+                raise ValueError(
+                    f"Unable to connect to '{HYDRODATA_URL}' to get data catalog information."
+                )
             else:
                 raise ValueError(
-                    f"Unable to connect to '{HYDRODATA_URL}' code = '{response.status_code}' to get data catalog information."
+                    f"Error ({response.status_code}) from server while trying to get data catalog information."
                 )
         return response_json
 
@@ -229,7 +233,7 @@ def _get_api_headers(required=True) -> dict:
             response = requests.get(url_security, timeout=(10, 100))
         except:
             raise ValueError(
-                f"Unable to authenticate with your email/pin with '{HYDRODATA_URL}' server."
+                f"Unable to connect to the server '{HYDRODATA_URL}' server."
             )
         if not response.status_code == 200:
             if not required:

src/hf_hydrodata/gridded.py

@@ -142,7 +142,7 @@ def generate_netcdf_file(
 
     # Add metadata attributes to the NC file
     for key, value in metadata.items():
-        setattr(nc, key, value)
+        setattr(nc, key, str(value))
 
     nc.close()
 
@@ -157,7 +157,9 @@ def create_time_variable(
 
     if temporal_resolution in ["daily", "hourly", "monthly", "weekly"]:
         date_start_str = get_date_start(query_parameters)
-        if ":" in date_start_str:
+        if isinstance(date_start_str, datetime.datetime):
+            date_start = date_start_str
+        elif ":" in date_start_str:
             date_start = datetime.datetime.strptime(date_start_str, "%Y-%m-%d %H:%M:%S")
         else:
             date_start = datetime.datetime.strptime(date_start_str, "%Y-%m-%d")
@@ -176,6 +178,8 @@ def create_time_variable(
             delta = datetime.timedelta(weeks=1)
         elif temporal_resolution == "monthly":
             delta = dateutil.relativedelta.relativedelta(months=1)
+        else:
+            delta = datetime.timedelta(days=1)
         dt = date_start
         for _ in range(0, time_steps):
             dates.append(dt)

src/hf_hydrodata/netcdf_util.py

@@ -1274,6 +1274,10 @@ def _validate_user():
         response = requests.get(url_security, headers=None, timeout=(60, 60))
         if response.status_code == 502:
             raise ValueError("Server is unavailable. Try again later.")
+        if response.status_code == 403:
+            raise ValueError(
+                f"Unable connect to '{HYDRODATA_URL}' server to validate user."
+            )
         if not response.status_code == 200:
             raise ValueError(
                 f"Unable to authenticate with your email/pin with '{HYDRODATA_URL}' server."

src/hf_hydrodata/point.py

@@ -1700,7 +1700,7 @@ def test_topographic_index(tmp_path):
         "grid_bounds": bounds,
         "variable": "topographic_index",
     }
-    hf.get_gridded_files(options, filename_template="foo.nc")
+    hf.get_gridded_file("foo.nc", options)
     ds = xr.open_dataset("foo.nc")
     da = ds["topographic_index"]
     assert da.shape == (5, 5)
@@ -2580,3 +2580,121 @@ def test_get_bounds_limits():
     with pytest.raises(ValueError) as info:
         hf.get_gridded_data(options)
     assert "is outside the grid" in str(info.value)
+
+
+def test_gridded_file_2wy():
+    """Test downloading subset in 2 waters years using get_gridded_file."""
+
+    with tempfile.TemporaryDirectory() as tempdirname:
+        bounds = [1399, 1784, 1447, 1803]
+        date_start = "1990-09-28"
+        date_end = "1990-10-04"
+        options = {
+            "dataset": "CW3E",
+            "variable": "air_temp",
+            "grid_bounds": bounds,
+            "temporal_resolution": "daily",
+            "date_start": date_start,
+            "date_end": date_end,
+        }
+
+        # Download to 2 different files for each water year
+        file_path = tempdirname + "/" + "{dataset}_{variable}_{wy}.nc"
+        hf.get_gridded_file(file_path, options)
+        path = f"{tempdirname}/CW3E_air_temp_1990.nc"
+        if not os.path.exists(path):
+            raise ValueError(f"File not created: '{path}'")
+        ds = xr.open_dataset(path)
+        assert ds["time"].shape[0] == 3
+        assert str(ds["time"].values[0]).startswith("1990-09-28")
+        assert ds["x"].shape[0] == 48
+        assert ds["y"].shape[0] == 19
+
+        # Download to one file for both water years
+        file_path = tempdirname + "/" + "{dataset}_{variable}.nc"
+        hf.get_gridded_file(file_path, options)
+        path = f"{tempdirname}/CW3E_air_temp.nc"
+        if not os.path.exists(path):
+            raise ValueError(f"File not created: '{path}'")
+        ds = xr.open_dataset(path)
+        assert ds["time"].shape[0] == 6
+        assert str(ds["time"].values[0]).startswith("1990-09-28")
+        assert ds["x"].shape[0] == 48
+        assert ds["y"].shape[0] == 19
+
+        # Download as one tiff (this should be an error)
+        with pytest.raises(ValueError):
+            file_path = tempdirname + "/" + "{dataset}_{variable}.tif"
+            hf.get_gridded_file(file_path, options)
+            path = f"{tempdirname}/CW3E_air_temp.nc"
+
+
+@pytest.mark.private_dataset
+def test_gridded_file_wtd():
+    """Test downloading subset of 30m WTD usign get_gridded_file."""
+
+    with tempfile.TemporaryDirectory() as tempdirname:
+        tempdirname = "."
+        bounds = [3000, 3000, 3010, 3020]
+        options = {
+            "dataset": "ma_2025",
+            "variable": "water_table_depth",
+            "grid_bounds": bounds,
+            "grid": "conus2_wtd",
+        }
+
+        # Download it as a tif
+        file_path = tempdirname + "/" + "{dataset}_{variable}.tif"
+        hf.get_gridded_file(file_path, options)
+        path = f"{tempdirname}/ma_2025_water_table_depth.tif"
+        if not os.path.exists(path):
+            raise ValueError(f"File not created: '{path}'")
+        ds = xr.open_dataset(path)
+        assert ds["x"].shape[0] == 10
+        assert ds["y"].shape[0] == 20
+
+        # Download it as netcdf
+        file_path = tempdirname + "/" + "{dataset}_{variable}.nc"
+        hf.get_gridded_file(file_path, options)
+        path = f"{tempdirname}/ma_2025_water_table_depth.nc"
+        if not os.path.exists(path):
+            raise ValueError(f"File not created: '{path}'")
+        ds = xr.open_dataset(path)
+        assert ds["x"].shape[0] == 10
+        assert ds["y"].shape[0] == 20
+
+
+def test_gridded_file_pfb():
+    """Test downloading subset in 2 waters years using get_gridded_file to pfb."""
+
+    with tempfile.TemporaryDirectory() as tempdirname:
+        bounds = [1399, 1784, 1447, 1803]
+        date_start = "1990-09-28"
+        date_end = "1990-10-04"
+        options = {
+            "dataset": "CW3E",
+            "variable": "air_temp",
+            "grid_bounds": bounds,
+            "temporal_resolution": "hourly",
+            "date_start": date_start,
+            "date_end": date_end,
+        }
+
+        # Download to 6 files for each day in separate water year directories
+        file_path = (
+            tempdirname
+            + "/"
+            + "WY{wy}/{dataset}_{variable}.{wy_start_24hr:06d}_to_{wy_end_24hr:06}.pfb"
+        )
+        hf.get_gridded_file(file_path, options)
+        path = f"{tempdirname}/WY1990/CW3E_air_temp.008689_to_008712.pfb"
+        if not os.path.exists(path):
+            raise ValueError(f"File not created: '{path}'")
+        data = hf.gridded.read_fast_pfb(path)
+        assert data.shape == (1, 24, 19, 48)
+
+        # Download to one file
+        file_path = tempdirname + "precip.pfb"
+        hf.get_gridded_file(file_path, options)
+        data = hf.gridded.read_fast_pfb(file_path)
+        assert data.shape == (1, 144, 19, 48)

tests/hf_hydrodata/test_gridded.py

@@ -97,7 +97,7 @@ def test_hourly():
         "start_time": "2006-01-01 00:00:00",
         "end_time": "2006-01-01 01:00:00",
         "file_type": "pfb",
-        "huc_id": "14010002",
+        "grid_bounds": [1000, 1000, 1020, 1010],
         "dataset": "CW3E",
     }
     with tempfile.TemporaryDirectory() as temp_dir:
@@ -109,13 +109,13 @@ def test_hourly():
         )
         ds = xr.open_dataset(file_name)
         assert os.path.exists(file_name)
-        assert len(ds["x"].values) == 52
-        assert len(ds["y"].values) == 76
+        assert len(ds["x"].values) == 20
+        assert len(ds["y"].values) == 10
         assert len(ds["time"].values) == 1
-        assert len(ds.coords["x"]) == 52
-        assert ds.coords["x"][0] == pytest.approx(-779500.3088117298)
+        assert len(ds.coords["x"]) == 20
+        assert ds.coords["x"].values[0] == pytest.approx(-1207500.3088117298)
         da = ds["air_temp"]
-        assert da.shape == (1, 76, 52)
-        assert da.values[0, 20, 20] == pytest.approx(273.6231)
+        assert da.shape == (1, 10, 20)
+        assert da.values[0, 5, 5] == pytest.approx(287.357452)
         assert ds.attrs["variable"] == "air_temp"
         assert ds.attrs["temporal_resolution"] == "hourly"