dft for chunked data along the axis of FT #51

roxyboy · 2018-08-28T18:38:23Z

@rabernat Something like this for Issue #49 ?

codecov-io · 2018-08-28T18:42:44Z

Codecov Report

❗ No coverage uploaded for pull request base (master@0c8833d). Click here to learn what that means.
The diff coverage is 95.55%.

@@            Coverage Diff            @@
##             master      #51   +/-   ##
=========================================
  Coverage          ?   88.92%           
=========================================
  Files             ?        1           
  Lines             ?      271           
  Branches          ?       90           
=========================================
  Hits              ?      241           
  Misses            ?       15           
  Partials          ?       15

Impacted Files	Coverage Δ
xrft/xrft.py	`88.92% <95.55%> (ø)`

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update 0c8833d...86a71da. Read the comment docs.

rabernat

Amazing start Takaya. You have become so fast!

I have a few comments on the api and the implementation.

rabernat · 2018-08-29T14:16:49Z

xrft/xrft.py

+        Whether the data is chunked along the axis to take FFT.
+    kwargs : dict
+        Need to provide the information of chunk length (`seglen`: int)
+        when `chunks_to_segment=True`.


a) we shouldn't hide arguments in kwargs unnecessarily

b) we shouldn't need to specify seglen at all! It should be determined by dask chunks themselves.

rabernat · 2018-08-29T14:19:05Z

xrft/xrft.py

+        if len(dim)>1:
+            raise NotImplementedError("Currently, only a one-dimensional "
+                                     "DFT is supported when the dimension "
+                                     "to FT is chunked.")


I don't see why we need this limitation. It should be possible to to 2d dft on 2d chunks (or any number of dimensions).

One needs to check each dimension in dim to see if reshaping is needed.

The current major issue with this is that the 2D, 3D detrending won't work if the axes are chunked. I'll see if I can change this easily.

rabernat · 2018-08-29T14:21:02Z

xrft/xrft.py

+                newdims.append(d + suffix)
+                newdims.append(d)
+                n = len(da[d])
+                seglen = kwargs['seglen']


Here we should instead figure out the chunks by examining the dask chunks.

An error should be raised if the chunks are not uniform.

rabernat · 2018-08-29T14:23:34Z

xrft/xrft.py

+                newdims.append(d)
+                newshape.append(len(da[d]))
+                newcoords[d] = da[d].data
+        da = xr.DataArray(data.reshape(newshape), dims=newdims, coords=newcoords)


I recommend putting the logic in the above ~20 lines into a standalone private function

def _stack_chunks(da, dim, suffix='_segment): """Reshape a DataArray so there is only one chunk along dimension `dim` ... """

rabernat · 2018-08-29T14:24:28Z

xrft/xrft.py

+        Whether the data is chunked along the axis to take FFT.
+    kwargs : dict
+        Need to provide the information of chunk length (`seglen`: int)
+        when `chunks_to_segment=True`.


Don't use kwargs this way.

You should review the zen of python. Line 2: "Explicit is better than implicit."

Still need to change this doc string to match the function signature.

rabernat · 2018-08-29T14:25:14Z

xrft/tests/test_xrft.py

+                 chunks_to_segments=True, seglen=8)
+    data = da.chunk({'time':8}).data.reshape((2,8,16,16,16))
+    npt.assert_almost_equal(ft.values, dsar.fft.fftn(data, axes=[1]).compute(),
+                           decimal=7)


Rather than absorbing this into an existing test, I would prefer to have a new test function for this feature.

roxyboy · 2018-08-29T15:48:15Z

@rabernat I incorporated everything except for the linear detrending for chunked axes. That's gonna take a bit longer so can we merge this for now and I'll make a new pull request for the chunked linear detrending?

roxyboy · 2018-08-29T16:08:19Z

Nevermind, got linear detrending for 2D cases to work.

…ys have a chunk length of one

rabernat · 2018-08-29T17:49:59Z

As far as I understand it, you should first transform the data to have reshaped dimensions, and then all the old windowing / detrending stuff should just work. It's really a pre-processing step.

roxyboy · 2018-08-29T17:51:15Z

Yea, so everything is working fine now will all your reviews incorporated :)

rabernat

Looks great. Almost there! Just a few more questions.

rabernat · 2018-08-29T17:50:29Z

xrft/tests/test_xrft.py

+    #                    )
+    # da = da5d.chunk({'time':1})
+    # with pytest.raises(ValueError):
+    #     func(da.data).compute()


Why is this commented out?

The current code should be able to take data that has more than 5 dimensions now so I commented it out.

rabernat · 2018-08-29T17:50:55Z

xrft/tests/test_xrft.py

+    # with pytest.raises(NotImplementedError):
+    #     xrft.dft(da, detrend='linear')
+    # with pytest.raises(NotImplementedError):
+    #     xrft.dft(da, dim=['time','y','x'], detrend='linear')


why do we have commented-out code here?

I'll just delete this as it is implemented.

rabernat · 2018-08-29T17:51:11Z

xrft/tests/test_xrft.py

+
+    # with pytest.raises(NotImplementedError):
+    #     xrft.dft(da.chunk({'time':16}), dim=['time','y'], detrend='linear',
+    #             chunks_to_segments=True)


Also implemented so I'll just delete it.

xrft/tests/test_xrft.py

rabernat · 2018-08-29T17:53:42Z

xrft/xrft.py

-            raise ValueError("Data has too many dimensions "
-                            "and/or too many axes to detrend over.")
+        if len(axes) > 3:
+            raise ValueError("Data has too many axes to detrend over.")


Let's be more explicit with this error, e.g. "Detrending is only possible over 3 dimensions or less, but input had 5 dimensions."

It works now with data that has more than 4 dimensions. Detrending itself is only supported up to 3 dimensions at the moment.

Made the error message more explicit to: "Detrending is only supported up to 3 dimensions."

rabernat · 2018-08-29T17:56:38Z

xrft/xrft.py

+            newshape.append(len(da[d]))
+            newcoords[d] = da[d].data
+
+    return xr.DataArray(data.reshape(newshape), dims=newdims, coords=newcoords)


Really nice function! 👍

xrft/xrft.py

rabernat · 2018-08-29T17:58:00Z

xrft/xrft.py

@@ -234,6 +270,9 @@ def dft(da, spacing_tol=1e-3, dim=None, shift=True, detrend=None, window=False):
    if detrend == 'constant':
        da = da - da.mean(dim=dim)
    elif detrend == 'linear':
+        for d in da.dims:
+            if d not in dim:
+                da = da.chunk({d:1})


Can you help me understand why you need to do this?

The linear algebra for detrending currently can only take one chunk for each node/processor.

But why does dask have to be involved? Are we using map_blocks to apply the detrending?

I'm worried about this. It is generally not a good idea to silently rechunk the user's data. What if this accidentally creates a billion chunks. I would prefer to raise an error and make the user do it.

Yes, we use map_blocks.

Ok, I'll make it raise an error but modify it so that it'll still automatically chunk the _segment dimensions.

rabernat · 2018-08-29T17:58:38Z

xrft/xrft.py

+        Whether the data is chunked along the axis to take FFT.
+    kwargs : dict
+        Need to provide the information of chunk length (`seglen`: int)
+        when `chunks_to_segment=True`.


Still need to change this doc string to match the function signature.

xrft/xrft.py

rabernat · 2018-08-29T18:01:09Z

I guess we probably eventually want to add the chunks_to_segments options to cross_spectrum and isotropic_power_spectrum too, no?

roxyboy · 2018-08-29T18:25:51Z

Adding chunks_to_segments to isotropic_power_spectrum would fail currently due to the same reasonings as in issue #9 . I haven't been able to parallelize the algorithm for azimuthal averaging.

roxyboy · 2018-08-30T05:06:54Z

xrft/xrft.py

+                if d not in dim and da.chunks[a_n][0]>1:
+                    raise ValueError("Linear detrending utilizes the `dask.map_blocks` "
+                                    "API so the dimensions not being detrended "
+                                    "must have the chunk length of 1.")


@rabernat I added this to raise an error if the user didn't provide data chunked with lengths of one for dimensions that won't be detrended.

roxyboy · 2018-08-30T05:09:37Z

xrft/xrft.py

@@ -283,17 +273,10 @@ def dft(da, spacing_tol=1e-3, dim=None, shift=True, detrend=None, window=False,
    if detrend == 'constant':
        da = da - da.mean(dim=dim)
    elif detrend == 'linear':
+        for d in da.dims:
+            if d not in dim:
+                da = da.chunk({d:1})


Although this may be a bit redundant, I kept it to make sure the chunks have unit length.

For the future, I recommend we use xarray.apply_ufunc to apply detrending, rather than map blocks. But that can wait for a future PR.

so... we can merge this?

rabernat

Final review I promise! 😉

rabernat · 2018-08-30T18:41:18Z

xrft/tests/test_xrft.py

    da_prime = xrft.detrendn(da[:,0].values, [0,1,2]) # cubic detrend over time, y, and x
    npt.assert_almost_equal(xrft.dft(da[:,0].drop('z'),
                                    dim=['time','y','x'],
                                    shift=False, detrend='linear'
                                    ).values,
                            np.fft.fftn(da_prime))

+def test_bartlett():


can we just give this test a more informative name. Like test_chunks_to_segments()?

rabernat · 2018-08-30T18:44:03Z

xrft/xrft.py

        Whether to shift the fft output.
    detrend : str (optional)
        If `constant`, the mean across the transform dimensions will be
        subtracted before calculating the Fourier transform (FT).
        If `linear`, the linear least-square fit will be subtracted before
        the FT.
-    window : bool (optional)
+    window : bool (default)


why did you change these from optional to default? That doesn't match the numpydoc spec. In fact, looking closely at the spec, it should be

window: bool, optional

rabernat · 2018-08-30T18:45:27Z

xrft/tests/test_xrft.py

+    npt.assert_almost_equal(ps.values,
+                           (ft*np.conj(ft)).real.values,
+                           )
+


There is no test for cross_spectrum.

roxyboy · 2018-08-30T19:26:31Z

Ok, the tests are passing so I'm gonna merge this.

dft for chunked data along the axis of FT

f6b98f3

roxyboy requested a review from rabernat August 28, 2018 18:39

update readthedocs for Bartlett

747e468

rabernat requested changes Aug 29, 2018

View reviewed changes

Takaya Uchida added 2 commits August 29, 2018 11:24

removed kwargs and refactored stacking

cacbeae

incorporated everything except linear detrending for chunked axes

0d19d85

linear detrending is robust now, I think

305ce6c

Takaya Uchida added 2 commits August 29, 2018 12:10

making sure when detrend is linear, the axes not being detrended alwa…

cb28483

…ys have a chunk length of one

refactoring chunking

7d025f4

rabernat requested changes Aug 29, 2018

View reviewed changes

Takaya Uchida added 5 commits August 29, 2018 14:32

added stuff suggested by rabernat

35c405b

made the error message more explicit

58a9eeb

keep attributes

aceca7e

no automatic rechunking

a94aa27

keeping the implicit functions simple

b1386af

roxyboy commented Aug 30, 2018

View reviewed changes

Takaya Uchida added 2 commits August 30, 2018 12:27

edited math

b65b4cd

before creating azimuth wrapper

d103c97

rabernat requested changes Aug 30, 2018

View reviewed changes

final push

86a71da

roxyboy merged commit 1ca9de9 into xgcm:master Aug 30, 2018

dft for chunked data along the axis of FT #51

dft for chunked data along the axis of FT #51

Conversation

roxyboy commented Aug 28, 2018

codecov-io commented Aug 28, 2018 • edited Loading

Codecov Report

rabernat left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

roxyboy commented Aug 29, 2018 • edited Loading

roxyboy commented Aug 29, 2018

rabernat commented Aug 29, 2018

roxyboy commented Aug 29, 2018

rabernat left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

roxyboy Aug 29, 2018 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

roxyboy Aug 29, 2018 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

roxyboy Aug 29, 2018 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

rabernat commented Aug 29, 2018

roxyboy commented Aug 29, 2018 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

rabernat left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

roxyboy commented Aug 30, 2018

codecov-io commented Aug 28, 2018 •

edited

Loading

roxyboy commented Aug 29, 2018 •

edited

Loading

roxyboy Aug 29, 2018 •

edited

Loading

roxyboy Aug 29, 2018 •

edited

Loading

roxyboy Aug 29, 2018 •

edited

Loading

roxyboy commented Aug 29, 2018 •

edited

Loading