Add examples for upcoming ulab guide

ulab_Crunch_Numbers_Fast/benchmark.py

@@ -0,0 +1,40 @@
+import time
+import math
+import ulab
+import ulab.numerical
+
+def mean(values):
+    return sum(values) / len(values)
+
+def normalized_rms(values):
+    minbuf = int(mean(values))
+    samples_sum = sum(
+        float(sample - minbuf) * (sample - minbuf)
+        for sample in values
+    )
+
+    return math.sqrt(samples_sum / len(values))
+
+def normalized_rms_ulab(values):
+    minbuf = ulab.numerical.mean(values)
+    values = values - minbuf
+    samples_sum = ulab.numerical.sum(values * values)
+    return math.sqrt(samples_sum / len(values))
+
+
+# Instead of using sensor data, we generate some data
+# The amplitude is 5000 so the rms should be around 5000/1.414 = 3536
+nums_list = [int(8000 + math.sin(i) * 5000) for i in range(100)]
+nums_array = ulab.array(nums_list)
+
+def timeit(s, f, n=100):
+    t0 = time.monotonic_ns()
+    for _ in range(n):
+        x = f()
+    t1 = time.monotonic_ns()
+    r = (t1 - t0) * 1e-6 / n
+    print("%-20s : %8.3fms [result=%f]" % (s, r, x))
+
+print("Computing the RMS value of 1000 numbers")
+timeit("traditional", lambda: normalized_rms(nums_list))
+timeit("ulab", lambda: normalized_rms_ulab(nums_array))

ulab_Crunch_Numbers_Fast/cluebarometer.py

@@ -0,0 +1,175 @@
+import time
+
+import adafruit_bmp280
+import board
+import displayio
+import ulab
+import ulab.filter
+
+# Blank the screen.  Scrolling text causes unwanted delays.
+d = displayio.Group()
+board.DISPLAY.show(d)
+
+# Sampling rate: 16Hz
+# Cutoff frequency: 0.16Hz
+# Transition bandwidth: 0.16Hz
+# Window type: Hamming
+# Filter has 311 coefficients
+taps = ulab.array([
+    -0.000050679794726066, -0.000041099278318167, -0.000031279920668665,
+    -0.000021183486597150, -0.000010770285292045, +0.000000000000000000,
+    +0.000011167446754809, +0.000022770999889941, +0.000034847558259864,
+    +0.000047431049079466, +0.000060551498686721, +0.000074234108254511,
+    +0.000088498343199344, +0.000103357045109305, +0.000118815575023601,
+    +0.000134870996840645, +0.000151511309510219, +0.000168714736477861,
+    +0.000186449080596567, +0.000204671152403140, +0.000223326279275218,
+    +0.000242347902542027, +0.000261657269119383, +0.000281163223679941,
+    +0.000300762106756334, +0.000320337763510928, +0.000339761667195315,
+    +0.000358893160569588, +0.000377579817760222, +0.000395657928211038,
+    +0.000412953103529159, +0.000429281007152519, +0.000444448205872873,
+    +0.000458253141344113, +0.000470487218795955, +0.000480936009263626,
+    +0.000489380560741255, +0.000495598812776238, +0.000499367108150093,
+    +0.000500461794444300, +0.000498660907473236, +0.000493745927786584,
+    +0.000485503600706003, +0.000473727809671115, +0.000458221492033063,
+    +0.000438798585855176, +0.000415285995764155, +0.000387525565446236,
+    +0.000355376044004699, +0.000318715033091691, +0.000277440901501588,
+    +0.000231474653767861, +0.000180761739242710, +0.000125273788160487,
+    +0.000065010261293197, +0.000000000000000000, -0.000069697336247377,
+    -0.000143989957415198, -0.000222752767634882, -0.000305826338672358,
+    -0.000393016043088374, -0.000484091357342654, -0.000578785344322494,
+    -0.000676794323931742, -0.000777777739462615, -0.000881358226495441,
+    -0.000987121890034750, -0.001094618794499868, -0.001203363670049808,
+    -0.001312836837542114, -0.001422485353209744, -0.001531724372895900,
+    -0.001639938734420840, -0.001746484755374530, -0.001850692242341569,
+    -0.001951866706278179, -0.002049291777482158, -0.002142231812333790,
+    -0.002229934682745978, -0.002311634738053158, -0.002386555927898205,
+    -0.002453915073551964, -0.002512925274028313, -0.002562799432345805,
+    -0.002602753886341418, -0.002632012127569287, -0.002649808591023194,
+    -0.002655392497711921, -0.002648031731496151, -0.002627016731069257,
+    -0.002591664377536210, -0.002541321857718479, -0.002475370483091317,
+    -0.002393229444145817, -0.002294359479963247, -0.002178266442894981,
+    -0.002044504738458277, -0.001892680620886388, -0.001722455325210333,
+    -0.001533548017297868, -0.001325738543930948, -0.001098869965763655,
+    -0.000852850856865069, -0.000587657355512251, -0.000303334951952833,
+    +0.000000000000000001, +0.000322159059450752, +0.000662880773589522,
+    +0.001021830060775982, +0.001398597909331569, +0.001792701398335994,
+    +0.002203584045179127, +0.002630616483032971, +0.003073097469789485,
+    +0.003530255228366684, +0.004001249116626688, +0.004485171623483914,
+    +0.004981050686118591, +0.005487852321559077, +0.006004483564265146,
+    +0.006529795699742466, +0.007062587782654920, +0.007601610426384373,
+    +0.008145569849526276, +0.008693132163411565, +0.009242927883419039,
+    +0.009793556645595150, +0.010343592108937170, +0.010891587022627668,
+    +0.011436078436539264, +0.011975593032464911, +0.012508652552774892,
+    +0.013033779302562583, +0.013549501700820601, +0.014054359855790191,
+    +0.014546911139352909, +0.015025735735186426, +0.015489442135386880,
+    +0.015936672560369614, +0.016366108277098043, +0.016776474791055797,
+    +0.017166546887869318, +0.017535153501103896, +0.017881182383493146,
+    +0.018203584559716979, +0.018501378539810983, +0.018773654273367416,
+    +0.019019576825867947, +0.019238389759765797, +0.019429418204303113,
+    +0.019592071599501125, +0.019725846101288819, +0.019830326636332028,
+    +0.019905188596781104, +0.019950199166862841, +0.019965218274992248,
+    +0.019950199166862841, +0.019905188596781104, +0.019830326636332028,
+    +0.019725846101288819, +0.019592071599501125, +0.019429418204303113,
+    +0.019238389759765800, +0.019019576825867947, +0.018773654273367420,
+    +0.018501378539810983, +0.018203584559716979, +0.017881182383493149,
+    +0.017535153501103892, +0.017166546887869318, +0.016776474791055797,
+    +0.016366108277098043, +0.015936672560369614, +0.015489442135386881,
+    +0.015025735735186426, +0.014546911139352912, +0.014054359855790193,
+    +0.013549501700820601, +0.013033779302562583, +0.012508652552774890,
+    +0.011975593032464912, +0.011436078436539264, +0.010891587022627668,
+    +0.010343592108937174, +0.009793556645595150, +0.009242927883419041,
+    +0.008693132163411567, +0.008145569849526276, +0.007601610426384373,
+    +0.007062587782654920, +0.006529795699742466, +0.006004483564265146,
+    +0.005487852321559078, +0.004981050686118592, +0.004485171623483914,
+    +0.004001249116626688, +0.003530255228366684, +0.003073097469789486,
+    +0.002630616483032971, +0.002203584045179127, +0.001792701398335993,
+    +0.001398597909331569, +0.001021830060775982, +0.000662880773589522,
+    +0.000322159059450752, +0.000000000000000001, -0.000303334951952833,
+    -0.000587657355512251, -0.000852850856865070, -0.001098869965763655,
+    -0.001325738543930948, -0.001533548017297868, -0.001722455325210333,
+    -0.001892680620886389, -0.002044504738458278, -0.002178266442894981,
+    -0.002294359479963247, -0.002393229444145818, -0.002475370483091317,
+    -0.002541321857718479, -0.002591664377536210, -0.002627016731069256,
+    -0.002648031731496151, -0.002655392497711923, -0.002649808591023195,
+    -0.002632012127569288, -0.002602753886341418, -0.002562799432345805,
+    -0.002512925274028314, -0.002453915073551965, -0.002386555927898205,
+    -0.002311634738053157, -0.002229934682745978, -0.002142231812333790,
+    -0.002049291777482159, -0.001951866706278179, -0.001850692242341569,
+    -0.001746484755374531, -0.001639938734420841, -0.001531724372895900,
+    -0.001422485353209745, -0.001312836837542114, -0.001203363670049807,
+    -0.001094618794499867, -0.000987121890034750, -0.000881358226495441,
+    -0.000777777739462615, -0.000676794323931742, -0.000578785344322494,
+    -0.000484091357342654, -0.000393016043088375, -0.000305826338672358,
+    -0.000222752767634882, -0.000143989957415198, -0.000069697336247377,
+    +0.000000000000000000, +0.000065010261293198, +0.000125273788160487,
+    +0.000180761739242710, +0.000231474653767861, +0.000277440901501588,
+    +0.000318715033091691, +0.000355376044004699, +0.000387525565446236,
+    +0.000415285995764155, +0.000438798585855176, +0.000458221492033064,
+    +0.000473727809671115, +0.000485503600706004, +0.000493745927786584,
+    +0.000498660907473236, +0.000500461794444300, +0.000499367108150093,
+    +0.000495598812776237, +0.000489380560741255, +0.000480936009263626,
+    +0.000470487218795956, +0.000458253141344114, +0.000444448205872873,
+    +0.000429281007152519, +0.000412953103529159, +0.000395657928211038,
+    +0.000377579817760223, +0.000358893160569588, +0.000339761667195315,
+    +0.000320337763510927, +0.000300762106756335, +0.000281163223679941,
+    +0.000261657269119383, +0.000242347902542027, +0.000223326279275218,
+    +0.000204671152403140, +0.000186449080596567, +0.000168714736477861,
+    +0.000151511309510219, +0.000134870996840645, +0.000118815575023601,
+    +0.000103357045109305, +0.000088498343199344, +0.000074234108254511,
+    +0.000060551498686721, +0.000047431049079466, +0.000034847558259864,
+    +0.000022770999889941, +0.000011167446754809, +0.000000000000000000,
+    -0.000010770285292045, -0.000021183486597150, -0.000031279920668665,
+    -0.000041099278318167, -0.000050679794726066,
+])
+
+# How often we are going to poll the sensor (If you change this, you need
+# to change the filter above and the integration time below)
+dt = 62500000 # 16Hz, 62.5ms
+
+# Wait until after deadline_ns has passed
+def sleep_deadline(deadline_ns):
+    while time.monotonic_ns() < deadline_ns:
+        pass
+
+
+# Initialize our sensor
+i2c = board.I2C()
+sensor = adafruit_bmp280.Adafruit_BMP280_I2C(i2c)
+sensor.standby_period = adafruit_bmp280.STANDBY_TC_1000
+# Disable in-sensor filtering, because we want to show how it's done in
+# CircuitPython
+sensor.iir_filter = adafruit_bmp280.IIR_FILTER_DISABLE
+sensor.overscan_pressure = adafruit_bmp280.OVERSCAN_X1
+
+# And our data structures
+# The most recent data samples, equal in number to the filter taps
+data = ulab.zeros(len(taps))
+t0 = deadline = time.monotonic_ns()
+n = 0
+# Take an initial reading to subtract off later, so that the graph in mu
+# accentuates the small short term changes in pressure rather than the large
+# DC offset of around 980
+offset = sensor.pressure
+
+while True:
+    deadline += dt
+    sleep_deadline(deadline)
+    # Move the trace near the origin so small differences can be seen in the mu
+    # plot window .. you wouldn't do this subtraction step if you are really
+    # interested in absolute barometric pressure.
+    value = sensor.pressure - offset
+    if n == 0:
+        # The first time, fill the filter with the initial value
+        data = data + value
+    else:
+        # Otherwise, add it as the next sample
+        ulab.numerical.roll(data, 1)
+        data[-1] = value
+    filtered = ulab.numerical.sum(data * taps)
+    # Actually print every 10th value.  This prints about 1.6 values per
+    # second.  You can print values more quickly by removing the 'if' and
+    # making the print unconditional, or change the frequency of prints up
+    # or down by changing the number '10'.
+    if n % 10 == 0:
+        print((filtered, value))
+    n += 1

ulab_Crunch_Numbers_Fast/cluepulse.py

@@ -0,0 +1,148 @@
+import time
+
+import adafruit_apds9960.apds9960
+import board
+import digitalio
+import ulab
+import ulab.filter
+
+# Blank the screen.  Scrolling text causes unwanted delays.
+import displayio
+d = displayio.Group()
+board.DISPLAY.show(d)
+
+# Filter computed at https://fiiir.com/
+# Sampling rate: 8Hz
+# Cutoff freqency: 0.5Hz
+# Transition bandwidth 0.25Hz
+# Window type: Regular
+# Number of coefficients: 31
+# Manually trimmed to 16 coefficients
+taps = ulab.array([
+    +0.861745279666917052/2,
+    -0.134728583242092248,
+    -0.124472980501612152,
+    -0.108421190967457198,
+    -0.088015688587190874,
+    -0.065052714580474319,
+    -0.041490993500537393,
+    -0.019246940463156042,
+    -0.000000000000000005,
+    +0.014969842582454691,
+    +0.024894596100322432,
+    +0.029569415718397409,
+    +0.029338562862396955,
+    +0.025020274838643962,
+    +0.017781854357373172,
+    +0.008981905549472832,
+])
+
+# How much reflected light is required before pulse sensor activates
+# These values are triggered when I bring my finger within a half inch.
+# The sensor works when the finger is pressed lightly against the sensor.
+PROXIMITY_THRESHOLD_HI = 225
+PROXIMITY_THRESHOLD_LO = 215
+
+# These constants control how much the sensor amplifies received light
+APDS9660_AGAIN_1X = 0
+APDS9660_AGAIN_4X = 1
+APDS9660_AGAIN_16X = 2
+APDS9660_AGAIN_64X = 3
+
+# How often we are going to poll the sensor (If you change this, you need
+# to change the filter above and the integration time below)
+dt = 125000000 # 8Hz, 125ms
+
+# Wait until after deadline_ns has passed
+def sleep_deadline(deadline_ns):
+    while time.monotonic_ns() < deadline_ns:
+        pass
+
+# Compute a high resolution crossing-time estimate for the sample, using a
+# linear model
+def estimated_cross_time(y0, y1, t0):
+    m = (y1 - y0) / dt
+    return t0 + round(-y1 / m)
+
+i2c = board.I2C()
+sensor = adafruit_apds9960.apds9960.APDS9960(i2c)
+white_leds = digitalio.DigitalInOut(board.WHITE_LEDS)
+white_leds.switch_to_output(False)
+
+def main():
+    sensor.enable_proximity = True
+    while True:
+        # Wait for user to put finger over sensor
+        while sensor.proximity() < PROXIMITY_THRESHOLD_HI:
+            time.sleep(.01)
+
+        # After the finger is sensed, set up the color sensor
+        sensor.enable_color = True
+        # This sensor integration time is just a little bit shorter than 125ms,
+        # so we should always have a fresh value when we ask for it, without
+        # checking if a value is available.
+        sensor.integration_time = 220
+        # In my testing, 64X gain saturated the sensor, so this is the biggest
+        # gain value that works properly.
+        sensor.color_gain = APDS9660_AGAIN_4X
+        white_leds.value = True
+
+        # And our data structures
+        # The most recent data samples, equal in number to the filter taps
+        data = ulab.zeros(len(taps))
+        # The filtered value on the previous iteration
+        old_value = 1
+        # The times of the most recent pulses registered.  Increasing this number
+        # makes the estimation more accurate, but at the expense of taking longer
+        # before a pulse number can be computed
+        pulse_times = []
+        # The estimated heart rate based on the recent pulse times
+        rate = None
+        # the number of samples taken
+        n = 0
+
+        # Rather than sleeping for a fixed duration, we compute a deadline
+        # in nanoseconds and wait for the new deadline time to arrive.  This
+        # helps the long term frequency of measurements better match the desired
+        # frequency.
+        t0 = deadline = time.monotonic_ns()
+        # As long as their finger is over the sensor, capture data
+        while sensor.proximity() >= PROXIMITY_THRESHOLD_LO:
+            deadline += dt
+            sleep_deadline(deadline)
+            value = sum(sensor.color_data) # Combination of all channels
+            ulab.numerical.roll(data, 1)
+            data[-1] = value
+            # Compute the new filtered variable by applying the filter to the
+            # recent data samples
+            filtered = ulab.numerical.sum(data * taps)
+
+            # We gathered enough data to fill the filters, and
+            # the light value crossed the zero line in the positive direction
+            # Therefore we need to record a pulse
+            if n > len(taps) and old_value < 0 and filtered >= 0:
+                # This crossing time is estimated, but it increases the pulse
+                # estimate resolution quite a bit.  If only the nearest 1/8s
+                # was used for pulse estimation, the smallest pulse increment
+                # that can be measured is 7.5bpm.
+                cross = estimated_cross_time(old_value, filtered, deadline)
+                # store this pulse time (in seconds since sensor-touch)
+                pulse_times.append((cross - t0) * 1e-9)
+                # and maybe delete an old pulse time
+                del pulse_times[:-10]
+                # And compute a rate based on the last recorded pulse times
+                if len(pulse_times) > 1:
+                    rate = 60/(pulse_times[-1]-pulse_times[0])*(len(pulse_times)-1)
+            old_value = filtered
+
+            # We gathered enough data to fill the filters, so report the light
+            # value and possibly the estimated pulse rate
+            if n > len(taps):
+                print((filtered, rate))
+            n += 1
+
+        # Turn off the sensor and the LED and go back to the top for another run
+        sensor.enable_color = False
+        white_leds.value = False
+        print()
+main()

ulab_Crunch_Numbers_Fast/ledwave.py

@@ -0,0 +1,76 @@
+import random
+
+import board
+import neopixel
+from _pixelbuf import wheel
+import ulab
+import ulab.filter
+
+# Customize your neopixel configuration here...
+pixel_pin = board.D5
+num_pixels = 96
+pixels = neopixel.NeoPixel(pixel_pin, num_pixels, brightness=0.1,
+                           auto_write=False, pixel_order=neopixel.RGB)
+
+ddt = ulab.array([1.,-2.,1.])
+def step(u, um, f, n, dx, dt, c):
+    dt2 = dt*dt
+    C2 = (c*dt/dx)**2
+    deriv = ulab.filter.convolve(u, ddt)[1:-1] * C2
+    up = -um + u * 2 + deriv + f * dt2
+    up[0] = 0
+    up[n-1] = 0
+
+    return up
+
+def main():
+    # This precomputes the color palette for maximum speed
+    # You could change it to compute the color palette of your choice
+    w = [wheel(i) for i in range(256)]
+
+    # This sets up the initial wave as a smooth gradient
+    u = ulab.zeros(num_pixels)
+    um = ulab.zeros(num_pixels)
+    f = ulab.zeros(num_pixels)
+
+    slope = ulab.linspace(0, 256, num=num_pixels)
+    th = 1
+
+    # the first time is always random (is that a contradiction?)
+    r = 0
+
+    while True:
+
+        # Some of the time, add a random new wave to the mix
+        # increase .15 to add waves more often
+        # decrease it to add waves less often
+        if r < .01:
+            ii = random.randrange(1, num_pixels-1)
+            # increase 2 to make bigger waves
+            f[ii] = (random.random() - .5) * 2
+
+        # Here's where to change dx, dt, and c
+        # try .2, .02, 2 for relaxed
+        # try 1., .7, .2 for very busy / almost random
+        u, um = step(u, um, f, num_pixels, .1, .02, 1), u
+
+        v = u * 200000 + slope + th
+        for i, vi in enumerate(v):
+            # Scale up by an empirical value, rotate by th, and look up the color
+            pixels[i] = w[round(vi) % 256]
+
+        # Take away a portion of the energy of the waves so they don't get out
+        # of control
+        u = u * .99
+
+        # incrementing th causes the wheel to slowly cycle even if nothing else is happening
+        th = (th + .25) % 256
+        pixels.show()
+
+        # Clear out the old random value, if any
+        f[ii] = 0
+
+        # and get a new random value
+        r = random.random()
+
+main()

ulab_Crunch_Numbers_Fast/waterfall.py

@@ -0,0 +1,96 @@
+"""Waterfall FFT demo adapted from
+https://teaandtechtime.com/fft-circuitpython-library/
+to work with ulab on Adafruit CLUE"""
+
+import array
+
+import board
+import audiobusio
+import displayio
+import ulab
+import ulab.fft
+import ulab.vector
+
+display = board.DISPLAY
+
+# Create a heatmap color palette
+palette = displayio.Palette(52)
+for i, pi in enumerate((0xff0000, 0xff0a00, 0xff1400, 0xff1e00,
+                        0xff2800, 0xff3200, 0xff3c00, 0xff4600,
+                        0xff5000, 0xff5a00, 0xff6400, 0xff6e00,
+                        0xff7800, 0xff8200, 0xff8c00, 0xff9600,
+                        0xffa000, 0xffaa00, 0xffb400, 0xffbe00,
+                        0xffc800, 0xffd200, 0xffdc00, 0xffe600,
+                        0xfff000, 0xfffa00, 0xfdff00, 0xd7ff00,
+                        0xb0ff00, 0x8aff00, 0x65ff00, 0x3eff00,
+                        0x17ff00, 0x00ff10, 0x00ff36, 0x00ff5c,
+                        0x00ff83, 0x00ffa8, 0x00ffd0, 0x00fff4,
+                        0x00a4ff, 0x0094ff, 0x0084ff, 0x0074ff,
+                        0x0064ff, 0x0054ff, 0x0044ff, 0x0032ff,
+                        0x0022ff, 0x0012ff, 0x0002ff, 0x0000ff)):
+    palette[51-i] = pi
+
+class RollingGraph(displayio.TileGrid):
+    def __init__(self, scale=2):
+        # Create a bitmap with heatmap colors
+        self.bitmap = displayio.Bitmap(display.width//scale,
+                                       display.height//scale, len(palette))
+        super().__init__(self.bitmap, pixel_shader=palette)
+
+        self.scroll_offset = 0
+
+    def show(self, data):
+        y = self.scroll_offset
+        bitmap = self.bitmap
+
+        board.DISPLAY.auto_refresh = False
+        offset = max(0, (bitmap.width-len(data))//2)
+        for x in range(min(bitmap.width, len(data))):
+            bitmap[x+offset, y] = int(data[x])
+
+        board.DISPLAY.auto_refresh = True
+
+        self.scroll_offset = (y + 1) % self.bitmap.height
+
+group = displayio.Group(scale=3)
+graph = RollingGraph(3)
+fft_size = 128
+
+# Add the TileGrid to the Group
+group.append(graph)
+
+# Add the Group to the Display
+display.show(group)
+
+# instantiate board mic
+mic = audiobusio.PDMIn(board.MICROPHONE_CLOCK, board.MICROPHONE_DATA,
+                       sample_rate=16000, bit_depth=16)
+
+#use some extra sample to account for the mic startup
+samples_bit = array.array('H', [0] * (fft_size+3))
+
+# Main Loop
+def main():
+    max_all = 1
+
+    while True:
+        mic.record(samples_bit, len(samples_bit))
+        samples = ulab.array(samples_bit[3:])
+        spectrogram1 = ulab.fft.spectrum(samples)
+        # spectrum() is always nonnegative, but add a tiny value
+        # to change any zeros to nonzero numbers
+        spectrogram1 = ulab.vector.log(spectrogram1 + 1e-7)
+        spectrogram1 = spectrogram1[1:(fft_size//2)-1]
+        min_curr = min(spectrogram1)
+        max_curr = max(spectrogram1)
+
+        if max_curr > max_all:
+            max_all = max_curr
+        else:
+            max_curr = max_curr-1
+
+        # Plot FFT
+        data = (spectrogram1 - min_curr) * (51. / (max_all - min_curr))
+        graph.show(data)
+
+main()