Merge branch 'v0.4'

docs/source/index.rst

@@ -14,6 +14,7 @@ Its primary application field is magnetotelluric (MT) and other electromagnetic
    citation
    overview
    tutorials/index
+   script
 
 
 .. toctree::

docs/source/script.rst

@@ -0,0 +1,103 @@
+Command line interface
+======================
+
+The :code:`rzb` command
+-----------------------
+
+Razorback provides some command line tools under the :code:`rzb` command.
+Typing the command shows the help message with the list of available commands :
+
+.. code-block:: bash
+
+  $ rzb
+  Usage: rzb [OPTIONS] COMMAND [ARGS]...
+  
+    The razorback command line interface.
+  
+  Options:
+    -h, --help  Show this message and exit.
+  
+  Commands:
+    path       Manipulate data path.
+    version    Razorback installed version.
+
+Getting the version number with :code:`rzb version`
+~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+To get the number of the installed version :
+
+.. code-block:: bash
+
+  $ rzb version
+  razorback 0.4.0
+
+Manipulating the data path with :code:`rzb path`
+~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+:code:`rzb path` provides several subcommands for manipulating the data path.
+Just typing :code:`rzb path` performs no action but shows the help and the list of available commands :
+
+.. code-block:: bash
+
+  $ rzb path
+  Usage: rzb path [OPTIONS] COMMAND [ARGS]...
+  
+    Manipulate data path.
+  
+    Data path is either global or local. If the local path is not available,
+    the global path is used instead.
+  
+    The path commands depend on the current directory where they are executed.
+  
+  Options:
+    -c, --create  Create the local path if missing.
+    -h, --help    Show this message and exit.
+  
+  Commands:
+    base      Current base data path.
+    metronix  Current path for Metronix calibration files.
+
+
+Without options, :code:`rzb path base` and :code:`rzb path metronix` just show the path where data, like calibration files, will be searched. This path depends on the working directory.
+
+The option :code:`--create` will create the corresponding local path.
+
+
+Extending the command line interface
+------------------------------------
+
+You can add commands to the :code:`rzb` command line interface by using `Click <https://click.palletsprojects.com/>`_ and `setuptools entry points <https://setuptools.readthedocs.io/en/latest/userguide/entry_point.html>`_.
+
+Let us look at an example.
+We have a simple `Click <https://click.palletsprojects.com/>`_ program in our package:
+
+.. code-block:: python
+
+  # mypkg/rzb_cli.py
+  import click
+
+  @click.command('say-hello')
+  def cli():
+      click.echo('Hello world !')
+
+We also have a `setup.py` for installing our package.
+To extend the :code:`rzb` command, we need to informs the `setup()` function in the following way:
+
+.. code-block:: python
+
+  # setup.py
+  setup(
+    # ...
+    entry_points={
+        'rzb.commands': [
+            'say-hello=mypkg.rzb_cli:cli',
+        ]
+    },
+  )
+
+Once `mypkg` is installed (:code:`python setup.py install` or :code:`pip install .`), the :code:`rzb` command can now expose our new subcommand:
+
+.. code-block:: bash
+
+  $ rzb say-hello
+  Hello world !

docs/source/tutorials/survey-study.ipynb

@@ -104,7 +104,7 @@
     "inv = rb.Inventory()\n",
     "for fname, [tag] in rb.utils.tags_from_path(files, pattern, tag_template):\n",
     "    calib = calibration(fname, name_converter)  # getting calibration for data file\n",
-    "    signal = rb.io.ats.load_ats([fname], [calib], lazy=True)  # loading data file\n",
+    "    signal = rb.io.ats.load_ats([fname], [calib])  # loading data file\n",
     "    inv.append(rb.SignalSet({tag:0}, signal))  # tagging and storing the signal"
    ]
   },

razorback/__init__.py

@@ -1,7 +1,7 @@
 """ razorback: tools for robust estimations of transfer functions.
 """
 
-__version__ = "0.3.1"
+__version__ = "0.4.0"
 
 try:
     import razorback_plus as plus

razorback/cli.py

@@ -0,0 +1,80 @@
+
+import pathlib
+import os
+from pkg_resources import iter_entry_points
+import click
+import razorback as rzb
+
+
+class PluginGroup(click.Group):
+
+    def __init__(self, *args, **kwds):
+        self.extra_commands = {
+            e.name: e.load() for e in iter_entry_points('rzb.commands')
+        }
+        super().__init__(*args, **kwds)
+
+    def list_commands(self, ctx):
+        return sorted(super().list_commands(ctx) + list(self.extra_commands))
+
+    def get_command(self, ctx, name):
+        return self.extra_commands.get(name) or super().get_command(ctx, name)
+
+
+@click.group(cls=PluginGroup, context_settings={'help_option_names': ('-h', '--help')})
+def cli():
+    """ The razorback command line interface.
+    """
+    pass
+
+
+@cli.command()
+def version():
+    """ Razorback installed version.
+    """
+    click.echo(f"razorback{'+' if rzb.plus else ''} {rzb.__version__}")
+
+
+@cli.group()
+@click.option('-c', '--create', is_flag=True,
+    help="Create the local path if missing."
+)
+@click.pass_context
+def path(ctx, create):
+    """ Manipulate data path.
+
+    Data path is either global or local.
+    If the local path is not available, the global path is used instead.
+
+    The path commands depend on the current directory where they are executed.
+
+    """
+    ctx.obj = {'create': create}
+
+
+def path_action(ctx, target):
+    if not ctx.obj['create']:
+        click.echo(rzb.data.data_path(target))
+    else:
+        p = rzb.data.local_data_path(target)
+        if p.exists():
+            click.echo(f"{p} already exists")
+        else:
+            p.mkdir(parents=True)
+            click.echo(f"{p} created")
+
+
+@path.command('base')
+@click.pass_context
+def path_base(ctx):
+    """ Current base data path.
+    """
+    path_action(ctx, '')
+
+
+@path.command('metronix')
+@click.pass_context
+def path_metronix(ctx):
+    """ Current path for Metronix calibration files.
+    """
+    path_action(ctx, rzb.calibrations.METRONIX_DATA_PATH)

razorback/data.py

@@ -7,7 +7,14 @@
 
 def default_data_path(target):
     " returns the default (internal) path to the data directory "
-    return pathlib.Path(__file__).parent / 'data' / target
+    p = pathlib.Path(__file__).parent / 'data' / target
+    return p.resolve().absolute()
+
+
+def local_data_path(target):
+    " returns the local path to the data directory "
+    p = pathlib.Path() / 'data' / target
+    return p.resolve().absolute()
 
 
 def data_path(target):

razorback/fourier_transform.py

@@ -4,6 +4,7 @@
 
 import numpy as np
 from scipy import signal
+import dask.array as da
 
 
 __all__ = ['time_to_freq', 'slepian_window']
@@ -20,7 +21,7 @@ def _slepian(N):
     return _slepian
 
 
-def time_to_freq(data, sampling_freq, freq, Nper, overlap, window=None):
+def time_to_freq(data, sampling_freq, freq, Nper, overlap, window=None, compute=True):
     """ Compute the fourier coefficients on sliding windows.
 
     Parameters
@@ -63,10 +64,22 @@ def time_to_freq(data, sampling_freq, freq, Nper, overlap, window=None):
         window = window(Lw)
     assert window is None or window.shape == (Lw,)
 
-    result = [
-        ft(sliding_window(y, Nw, Lw, shift).T, freq, sampling_freq, window)
+    data = (da.asarray(y) for y in data)
+    sw_views = [
+        da.overlap.sliding_window_view(y, Lw)[::shift]
         for y in data
     ]
+    pulsation = 2 * np.pi * nf
+    x = np.exp((-1j * pulsation) * da.arange(Lw)) * (2/Lw)
+    if window is not None:
+        x *= window
+    result = [da.dot(x, y.T) for y in sw_views]
+    if compute:
+        result = da.compute(*result)
+
+    length = set(len(e) for e in result)
+    assert len(length) == 1
+    Nw = length.pop()
 
     return result, (Nw, Lw, shift)
 
@@ -95,29 +108,6 @@ def discrete_window(size, normalized_freq, Nper, overlap):
     return int(Nw), int(Lw), shift
 
 
-def sliding_window(arr, nb_window, size_window, shift):
-    """ sliding window over arr
-
-    arr must be 1d
-    result is 2d: result[i] is the i-th window
-
-    >>> sliding_window(np.arange(15), 4, 3, 2)
-    array([[0, 1, 2],
-           [2, 3, 4],
-           [4, 5, 6],
-           [6, 7, 8]])
-
-    """
-    assert arr.ndim == 1, "array must be 1D"
-    min_ = (nb_window - 1) * shift + size_window
-    assert len(arr) >= min_, "array is too small (min=%d)" % min_
-    size = arr.itemsize
-    return np.lib.stride_tricks.as_strided(
-        arr, shape=(nb_window, size_window),
-        strides=(shift*size, size)
-    )
-
-
 def slepian(N, tau, N_MAX=1000):
     """ return the slepian window of size N with main lobe end at +/- tau
 
@@ -140,6 +130,52 @@ def slepian(N, tau, N_MAX=1000):
         res = np.interp(np.arange(N), a * np.arange(N_MAX), ref)
     return res
 
+#############################
+
+def time_to_freq_old(data, sampling_freq, freq, Nper, overlap, window=None):
+    length = set(len(e) for e in data)
+    assert len(length) == 1, "all data must have the same length"
+    Ntot = length.pop()
+
+    nf = np.true_divide(freq, sampling_freq)
+    Nw, Lw, shift = discrete_window(Ntot, nf, Nper, overlap)
+
+    if callable(window):
+        window = window(Lw)
+    assert window is None or window.shape == (Lw,)
+
+    result = [
+        ft(sliding_window(y, Nw, Lw, shift).T, freq, sampling_freq, window)
+        for y in data
+    ]
+
+    return result, (Nw, Lw, shift)
+
+
+def sliding_window(arr, nb_window, size_window, shift):
+    """ sliding window over arr
+
+    arr must be 1d
+    result is 2d: result[i] is the i-th window
+
+    >>> sliding_window(np.arange(15), 4, 3, 2)
+    array([[0, 1, 2],
+           [2, 3, 4],
+           [4, 5, 6],
+           [6, 7, 8]])
+
+    """
+    assert arr.ndim == 1, "array must be 1D"
+    min_ = (nb_window - 1) * shift + size_window
+    assert len(arr) >= min_, "array is too small (min=%d)" % min_
+    size = arr.itemsize
+    return np.lib.stride_tricks.as_strided(
+        # arr[:min_+1],
+        arr,
+        shape=(nb_window, size_window),
+        strides=(shift*size, size)
+    )
+
 
 def ft(y, freq, sampling_freq=1., window=None):
     """ return the complex Fourier coefficient of y at freq
@@ -155,3 +191,4 @@ def ft(y, freq, sampling_freq=1., window=None):
     if window is not None:
         x *= window
     return x.dot(y) / N * 2
+