diff --git a/code/scipy/signal/signal.c b/code/scipy/signal/signal.c index f930a943..53554e2d 100644 --- a/code/scipy/signal/signal.c +++ b/code/scipy/signal/signal.c @@ -10,6 +10,7 @@ * 2020 Scott Shawcroft for Adafruit Industries * 2020-2021 Zoltán Vörös * 2020 Taku Fukada + * 2024 Edouard Leroy (oaconvolve implementation) */ #include @@ -19,6 +20,8 @@ #include "../../ulab.h" #include "../../ndarray.h" #include "../../numpy/carray/carray_tools.h" +#include "../../numpy/fft/fft_tools.h" +#include "../../ulab_tools.h" #if ULAB_SCIPY_SIGNAL_HAS_SOSFILT & ULAB_MAX_DIMS > 1 static void signal_sosfilt_array(mp_float_t *x, const mp_float_t *coeffs, mp_float_t *zf, const size_t len) { @@ -120,11 +123,157 @@ mp_obj_t signal_sosfilt(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_ar MP_DEFINE_CONST_FUN_OBJ_KW(signal_sosfilt_obj, 2, signal_sosfilt); #endif /* ULAB_SCIPY_SIGNAL_HAS_SOSFILT */ +#if ULAB_SCIPY_SIGNAL_HAS_OACONVOLVE + +mp_obj_t signal_oaconvolve(size_t n_args, const mp_obj_t *pos_args, mp_map_t *kw_args) { + static const mp_arg_t allowed_args[] = { + {MP_QSTR_a, MP_ARG_REQUIRED | MP_ARG_OBJ, {.u_rom_obj = MP_ROM_NONE}}, + {MP_QSTR_v, MP_ARG_REQUIRED | MP_ARG_OBJ, {.u_rom_obj = MP_ROM_NONE}}, + }; + + mp_arg_val_t args[MP_ARRAY_SIZE(allowed_args)]; + mp_arg_parse_all(n_args, pos_args, kw_args, MP_ARRAY_SIZE(allowed_args), allowed_args, args); + + if (!mp_obj_is_type(args[0].u_obj, &ulab_ndarray_type) || !mp_obj_is_type(args[1].u_obj, &ulab_ndarray_type)) { + mp_raise_TypeError(MP_ERROR_TEXT("oa convolve arguments must be ndarrays")); + } + + ndarray_obj_t *a = MP_OBJ_TO_PTR(args[0].u_obj); + ndarray_obj_t *c = MP_OBJ_TO_PTR(args[1].u_obj); +// deal with linear arrays only +#if ULAB_MAX_DIMS > 1 + if ((a->ndim != 1) || (c->ndim != 1)) { + mp_raise_TypeError(MP_ERROR_TEXT("convolve arguments must be linear arrays")); + } +#endif + size_t signal_len = a->len; + size_t kernel_len = c->len; + if (signal_len == 0 || kernel_len == 0) { + mp_raise_TypeError(MP_ERROR_TEXT("convolve arguments must not be empty")); + } + + + + size_t segment_len = kernel_len; // Min segment size, at least size of kernel + size_t fft_size = 1; + while (fft_size < segment_len + kernel_len - 1) { + fft_size *= 2; // Find smallest power of 2 for fft size + } + segment_len = fft_size - kernel_len + 1; // Final segment size + size_t output_len = signal_len + kernel_len - 1; + uint8_t dtype = NDARRAY_FLOAT; +#if ULAB_SUPPORTS_COMPLEX + if ((a->dtype == NDARRAY_COMPLEX) || (c->dtype == NDARRAY_COMPLEX)) { + dtype = NDARRAY_COMPLEX; + } +#endif + + uint8_t sz = 2 * sizeof(mp_float_t); + + // Buffer allocation + mp_float_t *kernel_fft_array = m_new0(mp_float_t, fft_size * 2); + mp_float_t *segment_fft_array = m_new0(mp_float_t, fft_size * 2); + ndarray_obj_t *output = ndarray_new_linear_array(output_len, dtype); + mp_float_t *output_array = (mp_float_t *)output->array; + uint8_t *kernel_array = (uint8_t *)c->array; + uint8_t *signal_array = (uint8_t *)a->array; + + // Copy kernel data + if (c->dtype == NDARRAY_COMPLEX) { + for (size_t i = 0; i < kernel_len; i++){ + memcpy(kernel_fft_array, kernel_array, sz); + kernel_fft_array += 2; + kernel_array += c->strides[ULAB_MAX_DIMS - 1]; + } + } + else { + mp_float_t (*func)(void *) = ndarray_get_float_function(c->dtype); + for (size_t i = 0; i < kernel_len; i++) { + // real part; the imaginary part is 0, no need to assign + *kernel_fft_array = func(kernel_array); + kernel_fft_array += 2; + kernel_array += c->strides[ULAB_MAX_DIMS - 1]; + } + } + kernel_fft_array -= 2 * kernel_len; + + // Compute kernel fft in place + fft_kernel(kernel_fft_array, fft_size, 1); + + mp_float_t real, imag; // For complex multiplication + size_t current_segment_size = segment_len; + mp_float_t (*funca)(void *) = ndarray_get_float_function(a->dtype); + + for (size_t i = 0; i < signal_len; i += segment_len) { + // Check if remaining data is less than segment length + if (i + segment_len > signal_len) { + current_segment_size = signal_len - i; + } + // Load segment in buffer + memset(segment_fft_array, 0, fft_size * sz); // Reset buffer to zero + if (a->dtype == NDARRAY_COMPLEX) { + for (size_t k = 0; k < current_segment_size; k++){ + memcpy(segment_fft_array, signal_array, sz); + segment_fft_array += 2; + signal_array += a->strides[ULAB_MAX_DIMS - 1]; + } + } + else { + for (size_t k = 0; k < current_segment_size; k++) { + // real part; the imaginary part is 0, no need to assign + *segment_fft_array = funca(signal_array); + segment_fft_array += 2; + signal_array += a->strides[ULAB_MAX_DIMS - 1]; + } + } + segment_fft_array -= 2 * current_segment_size; + + // Compute segment fft in place + fft_kernel(segment_fft_array, fft_size, 1); + + // Product of kernel and segment fft (complex multiplication) + for (size_t j = 0; j < fft_size; j++) { + real = segment_fft_array[j * 2] * kernel_fft_array[j * 2] - segment_fft_array[j * 2 + 1] * kernel_fft_array[j * 2 + 1]; + imag = segment_fft_array[j * 2] * kernel_fft_array[j * 2 + 1] + segment_fft_array[j * 2 + 1] * kernel_fft_array[j * 2]; + segment_fft_array[j * 2] = real; + segment_fft_array[j * 2 + 1] = imag; + } + // Inverse FFT in place + fft_kernel(segment_fft_array, fft_size, -1); + + #if ULAB_SUPPORTS_COMPLEX + if (dtype == NDARRAY_COMPLEX) { + // Overlap, Add and normalized inverse fft + for (size_t j = 0; j < fft_size * 2; j++) { + if ((i * 2 + j) < (output_len * 2)) { + output_array[i * 2 + j] += (segment_fft_array[j] / fft_size); + } + } + continue; + } + #endif + + + for (size_t j = 0; j < fft_size; j++) { + if ((i + j) < (output_len)) { // adds only real part + output_array[i + j] += (segment_fft_array[j * 2] / fft_size); + } + } + } + return MP_OBJ_FROM_PTR(output); +} +MP_DEFINE_CONST_FUN_OBJ_KW(signal_oaconvolve_obj, 2, signal_oaconvolve); +#endif /* ULAB_SCIPY_SIGNAL_HAS_OACONVOLVE */ + + static const mp_rom_map_elem_t ulab_scipy_signal_globals_table[] = { { MP_ROM_QSTR(MP_QSTR___name__), MP_ROM_QSTR(MP_QSTR_signal) }, #if ULAB_SCIPY_SIGNAL_HAS_SOSFILT & ULAB_MAX_DIMS > 1 { MP_ROM_QSTR(MP_QSTR_sosfilt), MP_ROM_PTR(&signal_sosfilt_obj) }, #endif + #if ULAB_SCIPY_SIGNAL_HAS_OACONVOLVE + { MP_ROM_QSTR(MP_QSTR_oaconvolve), MP_ROM_PTR(&signal_oaconvolve_obj) }, + #endif }; static MP_DEFINE_CONST_DICT(mp_module_ulab_scipy_signal_globals, ulab_scipy_signal_globals_table); diff --git a/code/scipy/signal/signal.h b/code/scipy/signal/signal.h index 033f6e4c..3b713187 100644 --- a/code/scipy/signal/signal.h +++ b/code/scipy/signal/signal.h @@ -19,5 +19,6 @@ extern const mp_obj_module_t ulab_scipy_signal_module; MP_DECLARE_CONST_FUN_OBJ_KW(signal_sosfilt_obj); +MP_DECLARE_CONST_FUN_OBJ_KW(signal_oaconvolve_obj); #endif /* _SCIPY_SIGNAL_ */ diff --git a/code/ulab.c b/code/ulab.c index 12f37027..17058454 100644 --- a/code/ulab.c +++ b/code/ulab.c @@ -33,7 +33,7 @@ #include "user/user.h" #include "utils/utils.h" -#define ULAB_VERSION 6.6.1 +#define ULAB_VERSION 6.6.2 #define xstr(s) str(s) #define str(s) #s diff --git a/code/ulab.h b/code/ulab.h index 78b8a1ba..30c237b5 100644 --- a/code/ulab.h +++ b/code/ulab.h @@ -740,6 +740,10 @@ #define ULAB_SCIPY_SIGNAL_HAS_SOSFILT (1) #endif +#ifndef ULAB_SCIPY_SIGNAL_HAS_OACONVOLVE +#define ULAB_SCIPY_SIGNAL_HAS_OACONVOLVE (1) +#endif + #ifndef ULAB_SCIPY_HAS_OPTIMIZE_MODULE #define ULAB_SCIPY_HAS_OPTIMIZE_MODULE (1) #endif diff --git a/docs/manual/source/scipy-signal.rst b/docs/manual/source/scipy-signal.rst index d1f34818..0f078eba 100644 --- a/docs/manual/source/scipy-signal.rst +++ b/docs/manual/source/scipy-signal.rst @@ -2,9 +2,10 @@ scipy.signal ============ -This module defines the single function: +This module defines the functions: 1. `scipy.signal.sosfilt <#sosfilt>`__ +2. `scipy.signal.oaconvolve <#oaconvolve>`__ sosfilt ------- @@ -64,6 +65,39 @@ initial values are assumed to be 0. ======================================== zf: array([[37242.0, 74835.0], [1026187.0, 1936542.0]], dtype=float) + +oaconvolve +------- + +``scipy``: +https://docs.scipy.org/doc/scipy/reference/generated/scipy.signal.oaconvolve.html + +Convolve two N-dimensional arrays using the overlap-add method. + +Convolve in1 and in2 using the overlap-add method. Similarly to numpy.convolve, +this method works in full mode and other modes can be obtained by slicing the result? + +This is generally much faster than linear convolve for large arrays (n > ~500), +and generally much faster than fftconvolve (not implemented yet in ulab) when one array is much larger +than the other (e.g. for a matched filter algorithm), but can be slower when only a few output values +are needed or when the arrays are very similar in shape, and can only output float arrays (int or object array inputs will be cast to float). + +.. code:: + + # code to be run in micropython + + from ulab import numpy as np + from ulab import scipy as spy + + x = np.array((1,2,3)) + y = np.array((1,10,100,1000)) + result = spy.signal.oaconvolve(x, y) + print('result: ', result) + +.. parsed-literal:: + + result: array([1.0, 12.00024, 123.0001, 1230.0, 2300.0, 3000.0], dtype=float32) + diff --git a/docs/scipy-signal.ipynb b/docs/scipy-signal.ipynb index ec10d2e6..8abfcde5 100644 --- a/docs/scipy-signal.ipynb +++ b/docs/scipy-signal.ipynb @@ -2,7 +2,7 @@ "cells": [ { "cell_type": "code", - "execution_count": 2, + "execution_count": 1, "metadata": { "ExecuteTime": { "end_time": "2021-01-12T16:11:12.111639Z", @@ -11,10 +11,19 @@ }, "outputs": [ { - "name": "stdout", - "output_type": "stream", - "text": [ - "Populating the interactive namespace from numpy and matplotlib\n" + "ename": "ModuleNotFoundError", + "evalue": "No module named 'matplotlib'", + "output_type": "error", + "traceback": [ + "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m", + "\u001b[1;31mModuleNotFoundError\u001b[0m Traceback (most recent call last)", + "Cell \u001b[1;32mIn[1], line 1\u001b[0m\n\u001b[1;32m----> 1\u001b[0m \u001b[43mget_ipython\u001b[49m\u001b[43m(\u001b[49m\u001b[43m)\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mrun_line_magic\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;124;43m'\u001b[39;49m\u001b[38;5;124;43mpylab\u001b[39;49m\u001b[38;5;124;43m'\u001b[39;49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;124;43m'\u001b[39;49m\u001b[38;5;124;43minline\u001b[39;49m\u001b[38;5;124;43m'\u001b[39;49m\u001b[43m)\u001b[49m\n", + "File \u001b[1;32mc:\\Users\\EL221079\\Documents\\MICROPYTHON\\oaconvolve\\.venv\\Lib\\site-packages\\IPython\\core\\interactiveshell.py:2480\u001b[0m, in \u001b[0;36mInteractiveShell.run_line_magic\u001b[1;34m(self, magic_name, line, _stack_depth)\u001b[0m\n\u001b[0;32m 2478\u001b[0m kwargs[\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mlocal_ns\u001b[39m\u001b[38;5;124m'\u001b[39m] \u001b[38;5;241m=\u001b[39m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mget_local_scope(stack_depth)\n\u001b[0;32m 2479\u001b[0m \u001b[38;5;28;01mwith\u001b[39;00m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mbuiltin_trap:\n\u001b[1;32m-> 2480\u001b[0m result \u001b[38;5;241m=\u001b[39m \u001b[43mfn\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;241;43m*\u001b[39;49m\u001b[43margs\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;241;43m*\u001b[39;49m\u001b[38;5;241;43m*\u001b[39;49m\u001b[43mkwargs\u001b[49m\u001b[43m)\u001b[49m\n\u001b[0;32m 2482\u001b[0m \u001b[38;5;66;03m# The code below prevents the output from being displayed\u001b[39;00m\n\u001b[0;32m 2483\u001b[0m \u001b[38;5;66;03m# when using magics with decorator @output_can_be_silenced\u001b[39;00m\n\u001b[0;32m 2484\u001b[0m \u001b[38;5;66;03m# when the last Python token in the expression is a ';'.\u001b[39;00m\n\u001b[0;32m 2485\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;28mgetattr\u001b[39m(fn, magic\u001b[38;5;241m.\u001b[39mMAGIC_OUTPUT_CAN_BE_SILENCED, \u001b[38;5;28;01mFalse\u001b[39;00m):\n", + "File \u001b[1;32mc:\\Users\\EL221079\\Documents\\MICROPYTHON\\oaconvolve\\.venv\\Lib\\site-packages\\IPython\\core\\magics\\pylab.py:159\u001b[0m, in \u001b[0;36mPylabMagics.pylab\u001b[1;34m(self, line)\u001b[0m\n\u001b[0;32m 155\u001b[0m \u001b[38;5;28;01melse\u001b[39;00m:\n\u001b[0;32m 156\u001b[0m \u001b[38;5;66;03m# invert no-import flag\u001b[39;00m\n\u001b[0;32m 157\u001b[0m import_all \u001b[38;5;241m=\u001b[39m \u001b[38;5;129;01mnot\u001b[39;00m args\u001b[38;5;241m.\u001b[39mno_import_all\n\u001b[1;32m--> 159\u001b[0m gui, backend, clobbered \u001b[38;5;241m=\u001b[39m \u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mshell\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43menable_pylab\u001b[49m\u001b[43m(\u001b[49m\u001b[43margs\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mgui\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[43mimport_all\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mimport_all\u001b[49m\u001b[43m)\u001b[49m\n\u001b[0;32m 160\u001b[0m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39m_show_matplotlib_backend(args\u001b[38;5;241m.\u001b[39mgui, backend)\n\u001b[0;32m 161\u001b[0m \u001b[38;5;28mprint\u001b[39m(\n\u001b[0;32m 162\u001b[0m \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124m%\u001b[39m\u001b[38;5;124mpylab is deprecated, use \u001b[39m\u001b[38;5;124m%\u001b[39m\u001b[38;5;124mmatplotlib inline and import the required libraries.\u001b[39m\u001b[38;5;124m\"\u001b[39m\n\u001b[0;32m 163\u001b[0m )\n", + "File \u001b[1;32mc:\\Users\\EL221079\\Documents\\MICROPYTHON\\oaconvolve\\.venv\\Lib\\site-packages\\IPython\\core\\interactiveshell.py:3719\u001b[0m, in \u001b[0;36mInteractiveShell.enable_pylab\u001b[1;34m(self, gui, import_all, welcome_message)\u001b[0m\n\u001b[0;32m 3692\u001b[0m \u001b[38;5;250m\u001b[39m\u001b[38;5;124;03m\"\"\"Activate pylab support at runtime.\u001b[39;00m\n\u001b[0;32m 3693\u001b[0m \n\u001b[0;32m 3694\u001b[0m \u001b[38;5;124;03mThis turns on support for matplotlib, preloads into the interactive\u001b[39;00m\n\u001b[1;32m (...)\u001b[0m\n\u001b[0;32m 3715\u001b[0m \u001b[38;5;124;03m This argument is ignored, no welcome message will be displayed.\u001b[39;00m\n\u001b[0;32m 3716\u001b[0m \u001b[38;5;124;03m\"\"\"\u001b[39;00m\n\u001b[0;32m 3717\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mIPython\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01mcore\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01mpylabtools\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m import_pylab\n\u001b[1;32m-> 3719\u001b[0m gui, backend \u001b[38;5;241m=\u001b[39m \u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43menable_matplotlib\u001b[49m\u001b[43m(\u001b[49m\u001b[43mgui\u001b[49m\u001b[43m)\u001b[49m\n\u001b[0;32m 3721\u001b[0m \u001b[38;5;66;03m# We want to prevent the loading of pylab to pollute the user's\u001b[39;00m\n\u001b[0;32m 3722\u001b[0m \u001b[38;5;66;03m# namespace as shown by the %who* magics, so we execute the activation\u001b[39;00m\n\u001b[0;32m 3723\u001b[0m \u001b[38;5;66;03m# code in an empty namespace, and we update *both* user_ns and\u001b[39;00m\n\u001b[0;32m 3724\u001b[0m \u001b[38;5;66;03m# user_ns_hidden with this information.\u001b[39;00m\n\u001b[0;32m 3725\u001b[0m ns \u001b[38;5;241m=\u001b[39m {}\n", + "File \u001b[1;32mc:\\Users\\EL221079\\Documents\\MICROPYTHON\\oaconvolve\\.venv\\Lib\\site-packages\\IPython\\core\\interactiveshell.py:3665\u001b[0m, in \u001b[0;36mInteractiveShell.enable_matplotlib\u001b[1;34m(self, gui)\u001b[0m\n\u001b[0;32m 3662\u001b[0m \u001b[38;5;28;01mimport\u001b[39;00m \u001b[38;5;21;01mmatplotlib_inline\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01mbackend_inline\u001b[39;00m\n\u001b[0;32m 3664\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mIPython\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01mcore\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m pylabtools \u001b[38;5;28;01mas\u001b[39;00m pt\n\u001b[1;32m-> 3665\u001b[0m gui, backend \u001b[38;5;241m=\u001b[39m \u001b[43mpt\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mfind_gui_and_backend\u001b[49m\u001b[43m(\u001b[49m\u001b[43mgui\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;28;43mself\u001b[39;49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mpylab_gui_select\u001b[49m\u001b[43m)\u001b[49m\n\u001b[0;32m 3667\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m gui \u001b[38;5;241m!=\u001b[39m \u001b[38;5;28;01mNone\u001b[39;00m:\n\u001b[0;32m 3668\u001b[0m \u001b[38;5;66;03m# If we have our first gui selection, store it\u001b[39;00m\n\u001b[0;32m 3669\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m \u001b[38;5;28mself\u001b[39m\u001b[38;5;241m.\u001b[39mpylab_gui_select \u001b[38;5;129;01mis\u001b[39;00m \u001b[38;5;28;01mNone\u001b[39;00m:\n", + "File \u001b[1;32mc:\\Users\\EL221079\\Documents\\MICROPYTHON\\oaconvolve\\.venv\\Lib\\site-packages\\IPython\\core\\pylabtools.py:338\u001b[0m, in \u001b[0;36mfind_gui_and_backend\u001b[1;34m(gui, gui_select)\u001b[0m\n\u001b[0;32m 321\u001b[0m \u001b[38;5;28;01mdef\u001b[39;00m \u001b[38;5;21mfind_gui_and_backend\u001b[39m(gui\u001b[38;5;241m=\u001b[39m\u001b[38;5;28;01mNone\u001b[39;00m, gui_select\u001b[38;5;241m=\u001b[39m\u001b[38;5;28;01mNone\u001b[39;00m):\n\u001b[0;32m 322\u001b[0m \u001b[38;5;250m \u001b[39m\u001b[38;5;124;03m\"\"\"Given a gui string return the gui and mpl backend.\u001b[39;00m\n\u001b[0;32m 323\u001b[0m \n\u001b[0;32m 324\u001b[0m \u001b[38;5;124;03m Parameters\u001b[39;00m\n\u001b[1;32m (...)\u001b[0m\n\u001b[0;32m 335\u001b[0m \u001b[38;5;124;03m 'WXAgg','Qt4Agg','module://matplotlib_inline.backend_inline','agg').\u001b[39;00m\n\u001b[0;32m 336\u001b[0m \u001b[38;5;124;03m \"\"\"\u001b[39;00m\n\u001b[1;32m--> 338\u001b[0m \u001b[38;5;28;01mimport\u001b[39;00m \u001b[38;5;21;01mmatplotlib\u001b[39;00m\n\u001b[0;32m 340\u001b[0m \u001b[38;5;28;01mif\u001b[39;00m _matplotlib_manages_backends():\n\u001b[0;32m 341\u001b[0m backend_registry \u001b[38;5;241m=\u001b[39m matplotlib\u001b[38;5;241m.\u001b[39mbackends\u001b[38;5;241m.\u001b[39mregistry\u001b[38;5;241m.\u001b[39mbackend_registry\n", + "\u001b[1;31mModuleNotFoundError\u001b[0m: No module named 'matplotlib'" ] } ], @@ -22,6 +31,11 @@ "%pylab inline" ] }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [] + }, { "cell_type": "markdown", "metadata": {}, @@ -31,7 +45,7 @@ }, { "cell_type": "code", - "execution_count": 1, + "execution_count": null, "metadata": { "ExecuteTime": { "end_time": "2022-01-29T20:50:20.813162Z", @@ -49,7 +63,7 @@ }, { "cell_type": "code", - "execution_count": 2, + "execution_count": null, "metadata": { "ExecuteTime": { "end_time": "2022-01-29T20:50:21.613220Z", @@ -131,7 +145,7 @@ }, { "cell_type": "code", - "execution_count": 57, + "execution_count": null, "metadata": { "ExecuteTime": { "end_time": "2020-05-07T07:35:35.126401Z", @@ -147,7 +161,7 @@ }, { "cell_type": "code", - "execution_count": 9, + "execution_count": null, "metadata": { "ExecuteTime": { "end_time": "2020-05-19T19:11:18.145548Z", @@ -162,7 +176,7 @@ }, { "cell_type": "code", - "execution_count": 58, + "execution_count": null, "metadata": { "ExecuteTime": { "end_time": "2020-05-07T07:35:38.725924Z", @@ -225,9 +239,10 @@ "source": [ "# scipy.signal\n", "\n", - "This module defines the single function:\n", + "This module defines the functions:\n", "\n", - "1. [scipy.signal.sosfilt](#sosfilt)" + "1. [scipy.signal.sosfilt](#sosfilt)\n", + "2. [scipy.signal.oaconvolve](#oaconvolve)" ] }, { @@ -245,7 +260,7 @@ }, { "cell_type": "code", - "execution_count": 7, + "execution_count": null, "metadata": { "ExecuteTime": { "end_time": "2020-06-19T20:24:10.529668Z", @@ -277,7 +292,7 @@ }, { "cell_type": "code", - "execution_count": 11, + "execution_count": null, "metadata": { "ExecuteTime": { "end_time": "2020-06-19T20:27:39.508508Z", @@ -314,11 +329,60 @@ "print('y: ', y)\n", "print('\\n' + '='*40 + '\\nzf: ', zf)" ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## oaconvolve\n", + "\n", + "``scipy``:\n", + "https://docs.scipy.org/doc/scipy/reference/generated/scipy.signal.oaconvolve.html\n", + "\n", + "Convolve two N-dimensional arrays using the overlap-add method.\n", + "\n", + "Convolve in1 and in2 using the overlap-add method. Similarly to numpy.convolve, \n", + "this method works in full mode and other modes can be obtained by slicing the result?\n", + "\n", + "This is generally much faster than linear convolve for large arrays (n > ~500), \n", + "and generally much faster than fftconvolve (not implemented yet in ulab) when one array is much larger \n", + "than the other (e.g. for a matched filter algorithm), but can be slower when only a few output values \n", + "are needed or when the arrays are very similar in shape, and can only output float arrays (int or object array inputs will be cast to float).\n" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [ + { + "ename": "ModuleNotFoundError", + "evalue": "No module named 'ulab'", + "output_type": "error", + "traceback": [ + "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m", + "\u001b[1;31mModuleNotFoundError\u001b[0m Traceback (most recent call last)", + "Cell \u001b[1;32mIn[2], line 3\u001b[0m\n\u001b[0;32m 1\u001b[0m \u001b[38;5;66;03m# code to be run in micropython\u001b[39;00m\n\u001b[1;32m----> 3\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mulab\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m numpy \u001b[38;5;28;01mas\u001b[39;00m np\n\u001b[0;32m 4\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mulab\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m scipy \u001b[38;5;28;01mas\u001b[39;00m spy\n\u001b[0;32m 6\u001b[0m x \u001b[38;5;241m=\u001b[39m np\u001b[38;5;241m.\u001b[39marray((\u001b[38;5;241m1\u001b[39m,\u001b[38;5;241m2\u001b[39m,\u001b[38;5;241m3\u001b[39m))\n", + "\u001b[1;31mModuleNotFoundError\u001b[0m: No module named 'ulab'" + ] + } + ], + "source": [ + " # code to be run in micropython\n", + " \n", + "from ulab import numpy as np\n", + "from ulab import scipy as spy\n", + "\n", + "x = np.array((1,2,3))\n", + "y = np.array((1,10,100,1000))\n", + "result = spy.signal.oaconvolve(x, y)\n", + "print('result: ', result)" + ] } ], "metadata": { "kernelspec": { - "display_name": "Python 3", + "display_name": ".venv", "language": "python", "name": "python3" }, @@ -332,7 +396,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.8.5" + "version": "3.11.9" }, "toc": { "base_numbering": 1, diff --git a/tests/1d/scipy/oaconvolve.py b/tests/1d/scipy/oaconvolve.py new file mode 100644 index 00000000..48e8492a --- /dev/null +++ b/tests/1d/scipy/oaconvolve.py @@ -0,0 +1,16 @@ +import math + +try: + from ulab import scipy, numpy as np +except ImportError: + import scipy + import numpy as np + +x = np.array((1,2,3)) +y = np.array((1,10,100,1000)) +result = (scipy.signal.oaconvolve(x, y)) +ref_result = np.array([1, 12, 123, 1230, 2300, 3000],dtype=np.float) +cmp_result = [] +for p,q in zip(list(result), list(ref_result)): + cmp_result.append(math.isclose(p, q, rel_tol=1e-06, abs_tol=5e-04)) +print(cmp_result) diff --git a/tests/1d/scipy/oaconvolve.py.exp b/tests/1d/scipy/oaconvolve.py.exp new file mode 100644 index 00000000..63a3ac63 --- /dev/null +++ b/tests/1d/scipy/oaconvolve.py.exp @@ -0,0 +1 @@ +[True, True, True, True, True, True]