Skip to content

Commit f25ca7a

Browse files
blaylockbkblaylockbk
authored
Merge pull request #459 from nshankar/458-rotated-ll-proj
Add rotated_ll to supported CF CRS grid types
2 parents 6c2321b + 269b885 commit f25ca7a

4 files changed

Lines changed: 114 additions & 27 deletions

File tree

‎check_pygrib_vs_herbie_crs_extraction.ipynb‎

Lines changed: 26 additions & 24 deletions
Original file line numberDiff line numberDiff line change
@@ -15,19 +15,20 @@
1515
},
1616
{
1717
"cell_type": "code",
18-
"execution_count": 1,
18+
"execution_count": 2,
1919
"metadata": {},
2020
"outputs": [],
2121
"source": [
2222
"from herbie import Herbie\n",
2323
"import pygrib\n",
24+
"import pandas as pd\n",
2425
"from pyproj import CRS\n",
2526
"from herbie.crs import get_cf_crs"
2627
]
2728
},
2829
{
2930
"cell_type": "code",
30-
"execution_count": 6,
31+
"execution_count": 2,
3132
"metadata": {},
3233
"outputs": [
3334
{
@@ -49,7 +50,6 @@
4950
" pygrib {'a': 6371229, 'b': 6371229, 'proj': 'longlat'}\n",
5051
" Herbie {'a': 6371229, 'b': 6371229, 'proj': 'longlat'}\n",
5152
" equal= True\n",
52-
"👨🏻‍🏭 Created directory: [/home/blaylock/data/graphcast/20250101]\n",
5353
"\n",
5454
"Model: GRAPHCAST\n",
5555
" pygrib {'a': 4326.0, 'b': 4326.0, 'proj': 'longlat'}\n",
@@ -81,16 +81,6 @@
8181
" Herbie {'a': 6371229, 'b': 6371229, 'proj': 'longlat'}\n",
8282
" equal= True\n",
8383
"\n",
84-
"Model: HAFSA\n",
85-
" pygrib {'a': 6371229, 'b': 6371229, 'proj': 'longlat'}\n",
86-
" Herbie {'a': 6371229, 'b': 6371229, 'proj': 'longlat'}\n",
87-
" equal= True\n",
88-
"\n",
89-
"Model: HREF\n",
90-
" pygrib {'a': 6371229, 'b': 6371229, 'lat_0': 25.0, 'lat_1': 25.0, 'lat_2': 25.0, 'lon_0': 265.0, 'proj': 'lcc'}\n",
91-
" Herbie {'a': 6371229, 'b': 6371229, 'lat_0': 25.0, 'lat_1': 25.0, 'lat_2': 25.0, 'lon_0': 265.0, 'proj': 'lcc'}\n",
92-
" equal= True\n",
93-
"\n",
9484
"Model: NAM\n",
9585
" pygrib {'a': 6371229, 'b': 6371229, 'lat_0': 38.5, 'lat_1': 38.5, 'lat_2': 38.5, 'lon_0': 262.5, 'proj': 'lcc'}\n",
9686
" Herbie {'a': 6371229, 'b': 6371229, 'lat_0': 38.5, 'lat_1': 38.5, 'lat_2': 38.5, 'lon_0': 262.5, 'proj': 'lcc'}\n",
@@ -104,6 +94,11 @@
10494
"Model: RTMA\n",
10595
" pygrib {'a': 6371200.0, 'b': 6371200.0, 'lat_0': 25.0, 'lat_1': 25.0, 'lat_2': 25.0, 'lon_0': 265.0, 'proj': 'lcc'}\n",
10696
" Herbie {'a': 6371200.0, 'b': 6371200.0, 'lat_0': 25.0, 'lat_1': 25.0, 'lat_2': 25.0, 'lon_0': 265.0, 'proj': 'lcc'}\n",
97+
" equal= True\n",
98+
"\n",
99+
"Model: HRDPS\n",
100+
" pygrib {'a': 6371229, 'b': 6371229, 'lon_0': 245.305142, 'o_lat_p': 36.08852, 'o_lon_p': 0.0, 'o_proj': 'longlat', 'proj': 'ob_tran'}\n",
101+
" Herbie {'a': 6371229, 'b': 6371229, 'lon_0': 245.305142, 'o_lat_p': 36.08852, 'o_lon_p': 0.0, 'o_proj': 'longlat', 'proj': 'ob_tran'}\n",
107102
" equal= True\n"
108103
]
109104
}
@@ -119,24 +114,31 @@
119114
" (dict(model=\"aifs\"), \":2t:\"),\n",
120115
" (dict(model=\"gefs\", member=1), \":TMP:2 m above\"),\n",
121116
" (dict(model=\"gefs\", product=\"wave\", member=1), \":WIND:surf\"),\n",
122-
" (dict(model=\"hafsa\", product=\"storm.atm\", storm=\"07s\"), \":TMP:2 m above\"),\n",
123-
" (dict(model=\"href\", products=\"mean\", domain=\"conus\", fxx=1), \":TMP:2 m above\"),\n",
117+
" # FAILING\n",
118+
" # (dict(model=\"hafsa\", product=\"storm.atm\", storm=\"07s\"), \":TMP:2 m above\"),\n",
119+
" # FAILING\n",
120+
" # (dict(model=\"href\", products=\"mean\", domain=\"conus\", fxx=1), \":TMP:2 m above\"),\n",
124121
" (dict(model=\"nam\"), \":TMP:2 m above\"),\n",
125122
" (dict(model=\"urma\"), \":TMP:2 m above\"),\n",
126123
" (dict(model=\"rtma\"), \":TMP:2 m above\"),\n",
124+
" (dict(model=\"hrdps\", fxx=0, product=\"continental/2.5km\", variable=\"TMP\", level=\"AGL-2m\"), None),\n",
125+
"\n",
127126
"]:\n",
128127
" if model[\"model\"] in (\"hafsa\", \"href\"): # only on nomads\n",
129128
" date = \"2025-01-16 06:00\"\n",
129+
" elif model[\"model\"] in (\"hrdps\",): # no historical data\n",
130+
" date = pd.Timestamp(\"now\").floor(\"6h\") - pd.Timedelta(\"6h\")\n",
130131
" else:\n",
131-
" date = \"2025-01-01\"\n",
132+
" date = \"2025-07-05T06:00:00\"\n",
132133
"\n",
133-
" ds = Herbie(date, verbose=False, **model).xarray(\n",
134-
" search, remove_grib=False, _use_pygrib_for_crs=False\n",
135-
" )\n",
134+
" # ds = Herbie(date, verbose=False, **model).xarray(\n",
135+
" # search, remove_grib=False, _use_pygrib_for_crs=False\n",
136+
" # )\n",
137+
" ds = Herbie(date, verbose=False, **model).download()\n",
136138
"\n",
137139
" print()\n",
138140
" print(f\"Model: {model['model'].upper()}\")\n",
139-
" with pygrib.open(str(ds.local_grib)) as grb:\n",
141+
" with pygrib.open(str(ds)) as grb:\n",
140142
" msg = grb.message(1)\n",
141143
" projparams_dict_pygrib = dict(sorted(msg.projparams.items()))\n",
142144
" print(\" pygrib\", projparams_dict_pygrib)\n",
@@ -145,7 +147,7 @@
145147
" sorted(get_cf_crs(ds, _return_projparams=True).items())\n",
146148
" )\n",
147149
" print(\" Herbie\", projparams_dict_herbie)\n",
148-
" print(\" equal=\", projparams_dict_herbie == projparams_dict_pygrib)\n"
150+
" print(\" equal=\", projparams_dict_herbie == projparams_dict_pygrib)"
149151
]
150152
},
151153
{
@@ -5787,7 +5789,7 @@
57875789
],
57885790
"metadata": {
57895791
"kernelspec": {
5790-
"display_name": "herbie-dev",
5792+
"display_name": "Python 3 (ipykernel)",
57915793
"language": "python",
57925794
"name": "python3"
57935795
},
@@ -5801,9 +5803,9 @@
58015803
"name": "python",
58025804
"nbconvert_exporter": "python",
58035805
"pygments_lexer": "ipython3",
5804-
"version": "3.13.1"
5806+
"version": "3.12.9"
58055807
}
58065808
},
58075809
"nbformat": 4,
5808-
"nbformat_minor": 2
5810+
"nbformat_minor": 4
58095811
}

‎pyproject.toml‎

Lines changed: 2 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -36,7 +36,7 @@ dependencies = [
3636
"pandas>=2.1",
3737
"pyproj>=3.7.0",
3838
"requests>=2.23.3",
39-
"toml>=0.10.2", # TODO: Drop in favor of tomllib when Python >=3.11 is required.
39+
"toml>=0.10.2", # TODO: Drop in favor of tomllib when Python >=3.11 is required.
4040
"xarray>=2025.1.1",
4141
]
4242
dynamic = ["version"]
@@ -55,6 +55,7 @@ build-backend = "hatchling.build"
5555
[dependency-groups]
5656
dev = [
5757
"ipykernel>=6.29.5",
58+
"jupyter>=1.1.1",
5859
"pytest>=8.3.5",
5960
"pytest-cov>=6.1.1",
6061
"ruff>=0.11.6",

‎src/herbie/crs.py‎

Lines changed: 10 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -6,9 +6,8 @@
66

77
from typing import Any
88

9-
from pyproj import CRS
10-
119
import xarray as xr
10+
from pyproj import CRS
1211

1312

1413
def get_cf_crs(
@@ -82,6 +81,15 @@ def get_cf_crs(
8281
projparams["a"] = a
8382
projparams["b"] = b
8483

84+
elif da.GRIB_gridType == "rotated_ll":
85+
projparams = {"proj": "ob_tran"}
86+
projparams["o_proj"] = "longlat"
87+
projparams["a"] = a
88+
projparams["b"] = b
89+
projparams["lon_0"] = da.GRIB_longitudeOfSouthernPoleInDegrees
90+
projparams["o_lon_p"] = 0.0
91+
projparams["o_lat_p"] = -da.GRIB_latitudeOfSouthernPoleInDegrees
92+
8593
elif da.GRIB_gridType == "regular_gg":
8694
projparams = {"proj": "longlat"}
8795
projparams["a"] = a

‎tests/test_hrdps.py‎

Lines changed: 76 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,76 @@
1+
## Nikhil Shankar
2+
## August 10, 2025
3+
4+
"""
5+
Tests for downloading HRDPS model.
6+
7+
DRAGONS:
8+
9+
At ~00Z the HRDPS servers are cleared of the previous day's forecasts.
10+
The next forecast arrives at ~06Z. During the interim period these tests will fail.
11+
"""
12+
13+
from datetime import datetime
14+
15+
import pandas as pd
16+
17+
from herbie import Herbie, config
18+
19+
now = datetime.now()
20+
latest = pd.Timestamp("now").floor("6h") - pd.Timedelta("6h")
21+
22+
save_dir = config["default"]["save_dir"] / "Herbie-Tests-Data/"
23+
24+
# Remove all previous test data
25+
for i in (save_dir / "hrdps").rglob("*"):
26+
if i.is_file():
27+
i.unlink()
28+
29+
30+
def test_hrdps_download():
31+
H = Herbie(
32+
latest,
33+
model="hrdps",
34+
product="continental/2.5km",
35+
variable="TMP",
36+
level="AGL-2m",
37+
overwrite=True,
38+
save_dir=save_dir,
39+
)
40+
f = H.download()
41+
assert H.get_localFilePath().exists()
42+
f.unlink()
43+
44+
45+
def test_hrdps_xarray():
46+
H = Herbie(
47+
latest,
48+
model="hrdps",
49+
product="continental/2.5km",
50+
variable="TMP",
51+
level="AGL-2m",
52+
overwrite=True,
53+
save_dir=save_dir,
54+
)
55+
H.xarray(remove_grib=False)
56+
assert H.get_localFilePath().exists()
57+
H.get_localFilePath().unlink()
58+
59+
60+
def test_hrdps_to_netcdf():
61+
"""Check that a xarray Dataset can be written to a NetCDF file.
62+
63+
It is important that I have haven't put any python objects in the
64+
xarray Dataset attributes.
65+
"""
66+
H = Herbie(
67+
latest,
68+
model="hrdps",
69+
product="continental/2.5km",
70+
variable="TMP",
71+
level="AGL-2m",
72+
overwrite=True,
73+
save_dir=save_dir,
74+
)
75+
ds = H.xarray(remove_grib=False)
76+
ds.to_netcdf(save_dir / "test_hrdps_to_netcdf.nc")

0 commit comments

Comments
 (0)