diff --git a/src/radclss/util/column_utils.py b/src/radclss/util/column_utils.py index f6e5f27..6b22bac 100644 --- a/src/radclss/util/column_utils.py +++ b/src/radclss/util/column_utils.py @@ -669,7 +669,15 @@ def _accumulate_to_grid(grd_ds, column_time, default_step): """ target = np.asarray(column_time.values).ravel() step = _column_time_step(column_time, default_step) - edges = np.concatenate([[target[0] - step.to_timedelta64()], target]) + # pd.Timedelta.to_timedelta64 carries sub-second precision, so subtracting + # it promotes the leading edge -- and with it the whole array -- off the + # column's own datetime unit. interp hands that unit to the result, and + # _apply_match then assigns the result into a column whose time coordinate + # is still datetime64[s]; xarray compares the two coordinates by dtype as + # well as by value and rejects the write. Stay on the column's unit. + edges = np.concatenate([[target[0] - step.to_timedelta64()], target]).astype( + target.dtype + ) # A zero anchor ahead of the record lets the first output interval # difference against "nothing accumulated yet" rather than against a NaN. diff --git a/tests/test_column_utils.py b/tests/test_column_utils.py index 114bd0e..56c37c0 100644 --- a/tests/test_column_utils.py +++ b/tests/test_column_utils.py @@ -7,6 +7,7 @@ from radclss.util.column_utils import ( _accumulate_to_grid, + _apply_match, _column_time_step, get_nexrad_column, subset_points, @@ -244,6 +245,41 @@ def test_accumulate_to_grid_keeps_gaps_missing(): assert np.isnan(regridded["accum_nrt"].values).any() +def test_accumulate_to_grid_keeps_the_column_datetime_unit(): + """ + A real column time coordinate is built from base_time and so carries + second resolution, not the nanoseconds pd.date_range hands the tests + above. Subtracting a pandas Timedelta promoted the interpolation edges to + nanoseconds, and interp passed that unit on to the result. _apply_match + could then no longer write the result back into the column: xarray + compares coordinates by dtype as well as by value, so the two disagreed + despite holding identical timestamps. + """ + gauge = _synthetic_gauge() + seconds = gauge.time.values.astype("datetime64[s]")[::5] + column_time = xr.DataArray(seconds, dims="time", coords={"time": seconds}) + + regridded = _accumulate_to_grid(gauge, column_time, "5Min") + + assert regridded["time"].dtype == column_time.dtype + + # The write _apply_match performs in the pipeline, which is where the + # mismatched unit actually surfaced. + matched = regridded.assign_coords(station="M1").expand_dims("station") + column = xr.Dataset( + {"accum_nrt": (("station", "time"), np.zeros((1, seconds.size)))}, + coords={"station": ["M1"], "time": seconds}, + ) + + _apply_match(column, "M1", matched) + + np.testing.assert_allclose( + column["accum_nrt"].values[0], + regridded["accum_nrt"].values, + atol=1e-9, + ) + + def test_column_time_step_falls_back_when_unmeasurable(): """Degenerate grids fall back to the supplied default rather than raising.""" single = xr.DataArray(pd.to_datetime(["2025-06-19T00:00"]), dims="time")