Skip to content

can use pl33 filter with dask-backed DataArrays #88

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 3 commits into from
May 19, 2023
Merged

can use pl33 filter with dask-backed DataArrays #88

Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
19fcdb3
can use pl33 filter with dask-backed DataArrays
kthyng May 19, 2023
ab97de4
added attributes in too
kthyng May 19, 2023
cf5001b
changed attr desc to filtered from subtidal
kthyng May 19, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
38 changes: 33 additions & 5 deletions oceans/filters.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -415,10 +415,7 @@ def medfilt1(x, L=3):
>>> L = 103
>>> xout = medfilt1(x=x, L=L)
>>> ax = plt.subplot(212)
>>> (
... l1,
... l2,
... ) = ax.plot(
>>> (l1, l2,) = ax.plot(
... x
... ), ax.plot(xout)
>>> ax.grid(True)
Expand Down Expand Up @@ -582,6 +579,10 @@ def pl33tn(x, dt=1.0, T=33.0, mode="valid"):
and cosine tapered at each end to return a filtered time series
xf of the same length. Assumes length of x greater than 67.

Can input a DataArray and use dask-supported for lazy execution. In that
use case, dt is ignored and calculated from the input DataArray.
cf-xarray is also required.

Examples
--------
>>> from oceans.filters import pl33tn
Expand All @@ -603,6 +604,16 @@ def pl33tn(x, dt=1.0, T=33.0, mode="valid"):

"""

import cf_xarray # noqa: F401
import xarray as xr

if isinstance(x, xr.Dataset):
raise TypeError("Input a DataArray not a Dataset.")

if isinstance(x, xr.DataArray):
# find dt in units of hours
dt = (x.cf["T"][1] - x.cf["T"][0]) * 1e-9 / 3600
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This will get clobbered below, no?

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Those other lines use an input dt though, they just assume 1 hour, then modify the dt according to T.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

(though I did test with an hourly dataset so...)


pl33 = np.array(
[
-0.00027,
Expand Down Expand Up @@ -686,5 +697,22 @@ def pl33tn(x, dt=1.0, T=33.0, mode="valid"):
pl33 = np.interp(filter_time, _dt, pl33)
pl33 /= pl33.sum()

xf = np.convolve(x, pl33, mode=mode)
if isinstance(x, xr.DataArray):
weight = xr.DataArray(pl33, dims=["window"])
xf = (
x.rolling({x.cf["T"].name: len(pl33)}, center=True)
.construct({x.cf["T"].name: "window"})
.dot(weight)
)
# update attrs
attrs = {
key: f"{value}, filtered"
for key, value in x.attrs.items()
if key != "units"
}
xf.attrs = attrs

else: # use numpy
xf = np.convolve(x, pl33, mode=mode)

return xf