MAINT: Convert hysteresis_relay to C++.

Author: WarrenWeckesser

README.md

@@ -1703,7 +1703,7 @@ It generates the following plot:
 
 #### `backlash`
 
-`backlash(x, deadband, initial)`, a gufunc with signature `(i),(),()->(i)`,
+`backlash(x, deadband, initial)`, a gufunc with signature `(n),(),()->(n)`,
 computes the "backlash" response of a signal; see the Wikipedia article
 [Backlash (engineering)](https://en.wikipedia.org/wiki/Backlash_(engineering)).
 The function emulates the
@@ -1733,7 +1733,7 @@ the plot
 #### `hysteresis_relay`
 
 `hysteresis_relay(x, low_threshold, high_threshold, low_value, high_value, init)`,
-a gufunc with signature `(i),(),(),(),(),()->(i)`, passes `x` through a relay with
+a gufunc with signature `(n),(),(),(),(),()->(n)`, passes `x` through a relay with
 hysteresis (like a [Schmitt trigger](https://en.wikipedia.org/wiki/Schmitt_trigger)).
 The function is similar to the
 [relay block](https://www.mathworks.com/help/simulink/slref/relay.html)

pyproject.toml

@@ -8,7 +8,7 @@ requires = [
 
 [project]
 name = 'ufunclab'
-version = '0.0.8.dev11'
+version = '0.0.8.dev12'
 description = 'NumPy ufuncs and utilities.'
 readme = 'README.md'
 requires-python = '>=3.9'

src/hysteresis_relay/define_cxx_gufunc_extmod.py

@@ -0,0 +1,77 @@
+
+
+from ufunc_config_types import UFuncExtMod, UFunc, UFuncSource
+
+
+HYSTERESIS_RELAY_DOCSTRING = """\
+hysteresis_relay(x, low_threshold, high_threshold, low_value, high_value, init, /, ...)
+
+Pass x through a 'relay' with hysteresis.
+
+Parameters
+----------
+x : array_like
+    Input signal
+low_threshold : scalar
+    Low end of hysteresis interval.
+high_threshold : scalar
+    High end of the hysteresis interval.
+low_value : scalar
+    Output value for x < low_threshold.
+high_value : scalar
+    Outout value for x > high_threshold.
+init : scalar
+    Initial output value if the initial value of x is
+    between low_threshold and high_threshold.  Normally
+    this would be either low_value or high_value, but
+    the function does not require it.
+
+Returns
+-------
+out : ndarray
+    Output of the relay.
+
+Notes
+-----
+The function expects ``low_threshold <= high_threshold``.
+It does not check that this condition is satisfied.
+
+Examples
+--------
+>>> import numpy as np
+>>> from ufunclab import hysteresis_relay
+>>> x = np.array([-0.2, -0.6, -2, 0.2, 1.2, 2, 0.5, -0.7, -0.2, 0.7])
+
+`x` is the input signal.  The lower and upper thresholds
+are -0.5 and 0.5, respectively. The low and high output
+values are -1 and 1 (except for the initial output, which
+is 0, as determined by the last argument of `hysteresis_relay`).
+
+>>> hysteresis_relay(x, -0.5, 0.5, -1, 1, 0)
+array([ 0., -1., -1., -1.,  1.,  1.,  1., -1., -1.,  1.])
+"""
+
+
+hysteresis_relay_core_source = UFuncSource(
+    funcname='hysteresis_relay_core',
+    typesignatures=['ffffff->f', 'dddddd->d', 'gggggg->g'],
+)
+
+hysteresis_relay_gufunc = UFunc(
+    name='hysteresis_relay',
+    docstring=HYSTERESIS_RELAY_DOCSTRING,
+    header='hysteresis_relay_gufunc.h',
+    signature='(n),(),(),(),(),() -> (n)',
+    sources=[hysteresis_relay_core_source],
+)
+
+
+MODULE_DOCSTRING = """\
+This module defines the hysteresis_relay function.
+"""
+
+extmod = UFuncExtMod(
+    module='_hysteresis_relay',
+    docstring=MODULE_DOCSTRING,
+    ufuncs=[hysteresis_relay_gufunc],
+)

src/hysteresis_relay/hysteresis_relay_gufunc.c.src

@@ -0,0 +1,53 @@
+#ifndef HYSTERESIS_RELAY_GUFUNC_H
+#define HYSTERESIS_RELAY_GUFUNC_H
+
+#define PY_SSIZE_T_CLEAN
+#include "Python.h"
+
+#define NPY_NO_DEPRECATED_API NPY_API_VERSION
+#include "numpy/ndarraytypes.h"
+
+#include "../src/util/strided.hpp"
+
+
+template<typename T>
+static void
+hysteresis_relay_core(
+        npy_intp n,                     // core dimension n
+        T *p_x,                         // pointer to first element of x, a strided 1-d array with shape (n,)
+        const npy_intp x_stride,        // stride (in bytes) of x
+        T *p_low_threshold,
+        T *p_high_threshold,
+        T *p_low_value,
+        T *p_high_value,
+        T *p_init,
+        T *p_out,                       // pointer to out, a strided 1-d array with shape (n,)
+        const npy_intp out_stride       // stride (in bytes) of out
+)
+{
+    T low_threshold = *p_low_threshold;
+    T high_threshold = *p_high_threshold;
+    T low_value = *p_low_value;
+    T high_value = *p_high_value;
+    T init = *p_init;
+
+    for (npy_intp k = 0; k < n; ++k) {
+        T x = get(p_x, x_stride, k);
+        if (x < low_threshold) {
+            set(p_out, out_stride, k, low_value);
+        }
+        else if (x > high_threshold) {
+            set(p_out, out_stride, k, high_value);
+        }
+        else {
+            if (k == 0) {
+                set(p_out, out_stride, k, init);
+            }
+            else {
+                set(p_out, out_stride, k, get(p_out, out_stride, k - 1));
+            }
+        }
+    }
+}
+
+#endif

src/hysteresis_relay/hysteresis_relay_gufunc.h

@@ -49,6 +49,7 @@ gufunc_cxx_src_dirs = [
     'all_same',
     'backlash',
     'corr',
+    'hysteresis_relay',
     'mad',
     'meanvar',
     'minmax',
@@ -192,7 +193,6 @@ numpy_templated_c_gufunc_dirs = [
     'cross',
     'fillnan1d',
     'first',
-    'hysteresis_relay',
     'linear_interp1d',
     'means',
     'peaktopeak',