From cc5127e3342565592ba8e85b7d04fa189cb3479f Mon Sep 17 00:00:00 2001 From: Brad Barnett <127794626+bdbarnett@users.noreply.github.com> Date: Thu, 29 Oct 2026 03:31:19 -0500 Subject: [PATCH 1/3] audiofeatures: speech features for a wake word model, from PCM at any rate TensorFlow Lite Micro's audio microfrontend with microWakeWord's settings, carried in third_party with the kissfft it calls and built as one C translation unit. 16 kHz audio gives the reference build's bytes exactly; other rates go through a windowed-sinc resampler to 16 kHz first. feed() from Python anywhere, attach() to an audiopump.Tap on a board. Not built into CircuitPython, like audiometer. --- MANIFEST.in | 2 + NOTICE | 7 + ROADMAP.md | 4 + docs/get-started/install.md | 6 +- docs/guides/speech-features.md | 63 +++ docs/how-it-works/additions.md | 8 +- docs/index.md | 1 + docs/reference/modules.md | 1 + examples/audiofeatures_simpletest.py | 67 +++ micropython.cmake | 16 + micropython.mk | 14 + pyproject.toml | 10 +- setup.py | 6 +- src/audiofeatures/module.c | 269 +++++++++++ src/cpython/_audiodsp.c | 120 +++++ src/cpython/audiofeatures.py | 93 ++++ src/shared/audiodsp_features.c | 434 ++++++++++++++++++ src/shared/audiodsp_features.h | 68 +++ src/shared/audiodsp_microfrontend.c | 167 +++++++ src/shared/audiodsp_microfrontend_names.h | 62 +++ tests/data/audiofeatures/signal_16k.s8 | Bin 0 -> 11920 bytes tests/data/audiofeatures/signal_16k.u16 | Bin 0 -> 23840 bytes tests/parity/features_probe.py | 93 ++++ tests/parity/verify_dsp.py | 6 + tests/test_cpython_audiofeatures.py | 272 +++++++++++ tests/test_cpython_free_threading.py | 2 +- third_party/kissfft/COPYING | 11 + third_party/kissfft/ORIGIN.txt | 16 + third_party/kissfft/_kiss_fft_guts.h | 168 +++++++ third_party/kissfft/kiss_fft.c | 408 ++++++++++++++++ third_party/kissfft/kiss_fft.h | 124 +++++ third_party/kissfft/tools/kiss_fftr.c | 159 +++++++ third_party/kissfft/tools/kiss_fftr.h | 46 ++ third_party/microfrontend/LICENSE | 202 ++++++++ third_party/microfrontend/ORIGIN.txt | 20 + .../experimental/microfrontend/lib/bits.h | 102 ++++ .../experimental/microfrontend/lib/fft.cc | 52 +++ .../lite/experimental/microfrontend/lib/fft.h | 50 ++ .../microfrontend/lib/fft_util.cc | 70 +++ .../experimental/microfrontend/lib/fft_util.h | 34 ++ .../microfrontend/lib/filterbank.c | 134 ++++++ .../microfrontend/lib/filterbank.h | 63 +++ .../microfrontend/lib/filterbank_util.c | 220 +++++++++ .../microfrontend/lib/filterbank_util.h | 50 ++ .../experimental/microfrontend/lib/frontend.c | 72 +++ .../experimental/microfrontend/lib/frontend.h | 64 +++ .../microfrontend/lib/frontend_util.c | 85 ++++ .../microfrontend/lib/frontend_util.h | 52 +++ .../microfrontend/lib/kiss_fft_common.h | 48 ++ .../microfrontend/lib/kiss_fft_int16.cc | 10 + .../microfrontend/lib/kiss_fft_int16.h | 33 ++ .../experimental/microfrontend/lib/log_lut.c | 30 ++ .../experimental/microfrontend/lib/log_lut.h | 40 ++ .../microfrontend/lib/log_scale.c | 83 ++++ .../microfrontend/lib/log_scale.h | 39 ++ .../microfrontend/lib/log_scale_util.c | 27 ++ .../microfrontend/lib/log_scale_util.h | 45 ++ .../microfrontend/lib/noise_reduction.c | 51 ++ .../microfrontend/lib/noise_reduction.h | 46 ++ .../microfrontend/lib/noise_reduction_util.c | 45 ++ .../microfrontend/lib/noise_reduction_util.h | 50 ++ .../microfrontend/lib/pcan_gain_control.c | 56 +++ .../microfrontend/lib/pcan_gain_control.h | 47 ++ .../lib/pcan_gain_control_util.c | 92 ++++ .../lib/pcan_gain_control_util.h | 57 +++ .../experimental/microfrontend/lib/window.c | 70 +++ .../experimental/microfrontend/lib/window.h | 49 ++ .../microfrontend/lib/window_util.c | 73 +++ .../microfrontend/lib/window_util.h | 45 ++ 69 files changed, 5016 insertions(+), 13 deletions(-) create mode 100644 docs/guides/speech-features.md create mode 100644 examples/audiofeatures_simpletest.py create mode 100644 src/audiofeatures/module.c create mode 100644 src/cpython/audiofeatures.py create mode 100644 src/shared/audiodsp_features.c create mode 100644 src/shared/audiodsp_features.h create mode 100644 src/shared/audiodsp_microfrontend.c create mode 100644 src/shared/audiodsp_microfrontend_names.h create mode 100644 tests/data/audiofeatures/signal_16k.s8 create mode 100644 tests/data/audiofeatures/signal_16k.u16 create mode 100644 tests/parity/features_probe.py create mode 100644 tests/test_cpython_audiofeatures.py create mode 100644 third_party/kissfft/COPYING create mode 100644 third_party/kissfft/ORIGIN.txt create mode 100644 third_party/kissfft/_kiss_fft_guts.h create mode 100644 third_party/kissfft/kiss_fft.c create mode 100644 third_party/kissfft/kiss_fft.h create mode 100644 third_party/kissfft/tools/kiss_fftr.c create mode 100644 third_party/kissfft/tools/kiss_fftr.h create mode 100644 third_party/microfrontend/LICENSE create mode 100644 third_party/microfrontend/ORIGIN.txt create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/bits.h create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/fft.cc create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/fft.h create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/fft_util.cc create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/fft_util.h create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/filterbank.c create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/filterbank.h create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/filterbank_util.c create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/filterbank_util.h create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/frontend.c create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/frontend.h create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/frontend_util.c create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/frontend_util.h create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/kiss_fft_common.h create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/kiss_fft_int16.cc create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/kiss_fft_int16.h create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_lut.c create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_lut.h create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_scale.c create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_scale.h create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_scale_util.c create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_scale_util.h create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/noise_reduction.c create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/noise_reduction.h create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/noise_reduction_util.c create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/noise_reduction_util.h create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control.c create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control.h create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control_util.c create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control_util.h create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/window.c create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/window.h create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/window_util.c create mode 100644 third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/window_util.h diff --git a/MANIFEST.in b/MANIFEST.in index 33c2abdb..bffc5145 100644 --- a/MANIFEST.in +++ b/MANIFEST.in @@ -3,3 +3,5 @@ recursive-include src/shared *.h recursive-include lib/audiorender *.py *.md recursive-include lib/audiodsp_util *.py include third_party/tremor/*.c third_party/tremor/*.h third_party/tremor/COPYING third_party/tremor/ORIGIN.txt third_party/tremor/README +recursive-include third_party/microfrontend *.c *.cc *.h LICENSE ORIGIN.txt +recursive-include third_party/kissfft *.c *.h COPYING ORIGIN.txt diff --git a/NOTICE b/NOTICE index 0cc0fc11..e5f471cc 100644 --- a/NOTICE +++ b/NOTICE @@ -17,3 +17,10 @@ third_party/tremor/ is Xiph.Org's Tremor, an integer Ogg Vorbis decoder, used under the BSD-3-Clause license in third_party/tremor/COPYING (Copyright (c) 2002, Xiph.org Foundation). third_party/tremor/ORIGIN.txt names the upstream commit and the changes made to it. + +third_party/microfrontend/ is TensorFlow Lite Micro's audio microfrontend, +used under the Apache-2.0 license in third_party/microfrontend/LICENSE +(Copyright The TensorFlow Authors). third_party/kissfft/ is KISS FFT, used +under the BSD-3-Clause license in third_party/kissfft/COPYING (Copyright (c) +2003-2010 Mark Borgerding). Each directory's ORIGIN.txt names the upstream +commit it came from and how it is built. diff --git a/ROADMAP.md b/ROADMAP.md index db51b919..078aeaff 100644 --- a/ROADMAP.md +++ b/ROADMAP.md @@ -15,6 +15,10 @@ audio pump holding up under live use such as a stompbox or a stage rig. - A level ramp on `audiomixer` voices: a level move ramps over a set number of frames instead of stepping, so an effect can fade its wet path in or out without a click when you switch it in mid-note. +- Voice activity in `audiofeatures`: a flag per slice saying whether it holds + speech, from the frontend's own noise estimate, so a wake word model need + only run while someone is talking. The frontend costs about 14 KB of flash + on an ESP32-P4 today; this should add little to it. Bugs, and things you need that don't work yet, go to [issues](https://github.com/PyDevices/audiodsp/issues). diff --git a/docs/get-started/install.md b/docs/get-started/install.md index bddda5f2..08aa4f09 100644 --- a/docs/get-started/install.md +++ b/docs/get-started/install.md @@ -18,7 +18,7 @@ This gets you `audiocore`, `synthio`, `audiomixer`, `audiofilters`, `audiodelays`, `audiofreeverb`, `audiospeed`, `audiodynamics`, `audioroute`, `audiomath`, `audioecho`, `audioshaper`, `audioladder`, `audiooscillator`, `audioconvolve`, `audiobiquad`, `audioverb`, `audiomodal`, `audiosampler`, `audiometer`, -`audiopump`, `audiovorbis`, `audiograph`, and the `audiorender` package. +`audiofeatures`, `audiopump`, `audiovorbis`, `audiograph`, and the `audiorender` package. `audiomp3` remains firmware-only. The distribution declares no *required* runtime dependencies and does not itself publish an `audiodsp` import; its version is the `VERSION` file, which is also what `_audiodsp.__version__` @@ -105,5 +105,5 @@ own `audiobusio`, `audioio` and `audiopwmio`; `audiocore`, `synthio`, `audiomixer` and the effects modules stay CircuitPython's. Two modules differ there: `audiopump` carries only `Ring` (see [Play computed audio](../guides/play-computed-audio.md)), and `audiometer` -isn't built, because its sources are a board's pump and a USB sound card, and -CircuitPython has neither. +and `audiofeatures` aren't built, because their sources are a board's pump and +a USB sound card, and CircuitPython has neither. diff --git a/docs/guides/speech-features.md b/docs/guides/speech-features.md new file mode 100644 index 00000000..dcc64e56 --- /dev/null +++ b/docs/guides/speech-features.md @@ -0,0 +1,63 @@ +# Speech features + +`audiofeatures` turns audio into what a microWakeWord model listens to: every +10 ms, a slice of 40 values, one per mel band from 125 Hz to 7.5 kHz. Hand +the slices to a TensorFlow Lite Micro model and you have a wake word +detector; the model and the interpreter that runs it aren't part of audiodsp. + +```python +fe = audiofeatures.Frontend() +fe.feed(pcm, 48000, 2) # interleaved s16, any rate from 16 kHz +buf = bytearray(audiofeatures.FEATURE_SIZE * 10) +n = fe.readinto(buf) # n whole slices of 40 int8 +``` + +[`examples/audiofeatures_simpletest.py`](../../examples/audiofeatures_simpletest.py) +prints the slices of a WAV file, one row of characters each. + +On a board it can listen without Python touching a sample. +`fe.attach(tap, 48000)` reads an `audiopump.Tap` in C each time you call +`readinto()`, so whatever the pump is playing or recording, an I2S +microphone or a USB sound card, becomes slices. A tap runs at the pump's +rate and channel count; give `attach()` the rate. + +## What a slice is + +The frontend is TensorFlow Lite Micro's audio microfrontend, carried in +[`third_party/microfrontend`](../../third_party/microfrontend/ORIGIN.txt), +with the settings ESPHome's micro_wake_word models were trained on: a 30 ms +window every 10 ms, 40 channels, noise reduction, PCAN gain control and a +log scale. `readinto()` gives each value as the int8 those models take, +`(v * 256 + 333) // 666 - 128`, and `audiofeatures.quantize()` is that sum. +`readinto_raw()` gives the frontend's own uint16 values instead, about 0 to +670. + +At 16 kHz the samples go to the frontend untouched, and the bytes are the +ones TensorFlow Lite Micro's own C++ build makes, in any block size. +[`tests/test_cpython_audiofeatures.py`](../../tests/test_cpython_audiofeatures.py) +holds goldens that build wrote. + +## Other rates + +Audio at any other rate up to 384 kHz is resampled to 16 kHz first, through +a windowed-sinc low-pass that is flat to 7.5 kHz, where the top mel band +ends, and about 70 dB down from 8.5 kHz up. Stereo and wider audio is mixed +to mono before it. The features can't be byte-identical to the 16 kHz ones, +but they are close: on 74 spoken clips fed at 48 kHz, 95% of the int8 values +are the same as at 16 kHz and the mean difference is 0.27 of a step, and +three microWakeWord models detected and missed exactly the same clips, with +peak probabilities within 0.03. + +## The queue + +Slices wait in a queue of `slices` (100 by default, a second of audio) until +you read them. When it is full, the oldest slice goes and `stats()["dropped"]` +counts it, so a reader that falls behind loses old audio, never new. +`reset()` forgets the audio and the queue, for the start of a new utterance. + +A frontend takes about 17 to 26 KB of heap with the default queue, depending +on the rate: 8 KB is the queue, the rest the frontend's state and the +resampler's filter. + +CircuitPython doesn't get `audiofeatures`; see +[additions](../how-it-works/additions.md). diff --git a/docs/how-it-works/additions.md b/docs/how-it-works/additions.md index 834bac35..5865b38d 100644 --- a/docs/how-it-works/additions.md +++ b/docs/how-it-works/additions.md @@ -40,10 +40,10 @@ These have no ancestor anywhere and are audiodsp's own: Every one of them builds into CircuitPython 11 too, as a user C module; see [Install](../get-started/install.md#circuitpython-11). -Three more are for cases CircuitPython doesn't cover. `audiometer` (band -levels for a spectrum meter) is the one CircuitPython doesn't get: its -sources are a board's pump and a USB sound card, and CircuitPython has -neither. `audiovorbis` plays Ogg Vorbis files, which CircuitPython has no +Four more are for cases CircuitPython doesn't cover. `audiometer` (band +levels for a spectrum meter) and `audiofeatures` (speech features for a wake +word model) are the two CircuitPython doesn't get: their sources are a board's +pump and a USB sound card, and CircuitPython has neither. `audiovorbis` plays Ogg Vorbis files, which CircuitPython has no module for. And `audiopump.Ring` lets CircuitPython play audio your own code computes. diff --git a/docs/index.md b/docs/index.md index f1613733..42fa4e88 100644 --- a/docs/index.md +++ b/docs/index.md @@ -21,6 +21,7 @@ have. Start with install, then your first program. - [Filter a sound](guides/filter-a-sound.md) with `audiobiquad`. - [Play recorded notes](guides/play-recorded-notes.md) with `audiosampler`. - [Meter levels](guides/meter-levels.md) for a spectrum display. +- [Speech features](guides/speech-features.md) for a wake word model. - [Render on several threads](guides/use-threads.md), including free-threaded Python. ## Reference diff --git a/docs/reference/modules.md b/docs/reference/modules.md index 67d0cfbf..f9b0cbd3 100644 --- a/docs/reference/modules.md +++ b/docs/reference/modules.md @@ -25,6 +25,7 @@ It adds modules CircuitPython doesn't have. Each has a page here or a guide: | `audiograph` | A whole patch rendered by one C call per block | [reference](sources-and-synthesis.md#audiograph), [how it works](../how-it-works/audiograph.md) | | `audiosampler` | Recorded notes played at any pitch, looped, under an envelope | [guide](../guides/play-recorded-notes.md) | | `audiometer` | Band levels for a spectrum meter | [guide](../guides/meter-levels.md) | +| `audiofeatures` | Speech features for a wake word model | [guide](../guides/speech-features.md) | | `audiovorbis` | Ogg Vorbis files | [guide](../guides/decode-vorbis.md), [C interface](vorbis.md) | | `audiopump` | The audio pull, off the interpreter thread | [guide](../guides/run-the-audio-pump.md), [reference](audiopump.md) | diff --git a/examples/audiofeatures_simpletest.py b/examples/audiofeatures_simpletest.py new file mode 100644 index 00000000..f6f33813 --- /dev/null +++ b/examples/audiofeatures_simpletest.py @@ -0,0 +1,67 @@ +# SPDX-FileCopyrightText: 2026 Brad Barnett +# +# SPDX-License-Identifier: MIT +"""Print the speech features a wake word model sees in a WAV file. + +An audiofeatures Frontend turns every 10 ms of audio into a slice of 40 +values, one per mel band from 125 Hz to 7.5 kHz. The program feeds it a WAV +file a block at a time and prints each slice as a row of characters, low +bands on the left, darker for more energy, so speech shows up as dense rows +and silence as blank ones. + +Needs: audiodsp, built into MicroPython or installed on CPython with +pip install pydevices-audiodsp, and a 16-bit PCM WAV file at 16 kHz or +higher, mono or stereo. + +Run it: + micropython audiofeatures_simpletest.py + +Tested on: the unix port of MicroPython and CPython 3.12 on Linux. +""" + +import struct + +import audiofeatures + +# Change these to match your setup. +FILENAME = "speech.wav" +BLOCK = 480 # frames read and fed at a time +SHADES = " .:-=+*#%@" # from no energy to the most + + +def main(): + with open(FILENAME, "rb") as f: + # Walk the WAV's chunks to its format and the start of its samples. + assert f.read(12)[8:] == b"WAVE", "not a WAV file" + while True: + tag, size = struct.unpack("<4sI", f.read(8)) + if tag == b"data": + break + body = f.read(size + (size & 1)) + if tag == b"fmt ": + kind, channels, rate = struct.unpack(" 0: + block = f.read(min(left, BLOCK * 2 * channels)) + left -= len(block) + frontend.feed(block, rate, channels) + count = frontend.readinto(features) + for i in range(count * audiofeatures.FEATURE_SIZE): + # Each value is a signed byte, -128 to 127, held unsigned here. + value = (features[i] + 128) & 0xFF + print(SHADES[value * len(SHADES) // 256], end="") + if i % audiofeatures.FEATURE_SIZE == audiofeatures.FEATURE_SIZE - 1: + print() + rows += count + + # Every 10 ms after the first 30 ms makes a slice. + print("%d slices from %d Hz audio, %d dropped" % (rows, rate, frontend.stats()["dropped"])) + + +if __name__ == "__main__": + main() diff --git a/micropython.cmake b/micropython.cmake index 0d0be950..2b3290cd 100644 --- a/micropython.cmake +++ b/micropython.cmake @@ -249,6 +249,22 @@ target_sources(usermod_mpaudio INTERFACE ${MPAUDIO_SRC_DIR}/audiovorbis/module.c ) +# --- audiofeatures (Frontend) --- +# +# See micropython.mk: TensorFlow Lite Micro's microfrontend and kissfft from +# third_party, compiled as one translation unit by +# src/shared/audiodsp_microfrontend.c, which renames their symbols and silences +# their warnings for their own code only. +target_sources(usermod_mpaudio INTERFACE + ${MPAUDIO_SRC_DIR}/shared/audiodsp_microfrontend.c + ${MPAUDIO_SRC_DIR}/shared/audiodsp_features.c + ${MPAUDIO_SRC_DIR}/audiofeatures/module.c +) +target_include_directories(usermod_mpaudio INTERFACE + ${MPAUDIO_MOD_DIR}/third_party/microfrontend + ${MPAUDIO_MOD_DIR}/third_party/kissfft +) + # CIRCUITPY_SYNTHIO_MAX_CHANNELS (phase 8e decision, made while actually # building the ESP32-P4 firmware -- see docs/porting-plan.md "8e"): the diff --git a/micropython.mk b/micropython.mk index 411e1e70..763a31f3 100644 --- a/micropython.mk +++ b/micropython.mk @@ -388,4 +388,18 @@ SRC_USERMOD_C += \ $(MPAUDIO_SRC_DIR)/audiovorbis/VorbisDecoder.c \ $(MPAUDIO_SRC_DIR)/audiovorbis/module.c +# --- audiofeatures (Frontend) --- +# +# TensorFlow Lite Micro's audio microfrontend and the kissfft it calls, carried +# in third_party/microfrontend (Apache-2.0) and third_party/kissfft +# (BSD-3-Clause), compiled as one translation unit by +# src/shared/audiodsp_microfrontend.c the way Tremor is. The library names its +# headers by their full tensorflow/... path, hence the two include directories. +# docs/guides/speech-features.md has what it does. +CFLAGS_USERMOD += -I$(MPAUDIO_MOD_DIR)/third_party/microfrontend -I$(MPAUDIO_MOD_DIR)/third_party/kissfft +SRC_USERMOD_LIB_C += \ + $(MPAUDIO_SRC_DIR)/shared/audiodsp_microfrontend.c \ + $(MPAUDIO_SRC_DIR)/shared/audiodsp_features.c +SRC_USERMOD_C += $(MPAUDIO_SRC_DIR)/audiofeatures/module.c + endif # not CircuitPython diff --git a/pyproject.toml b/pyproject.toml index 1f7a8636..58445e3c 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -9,9 +9,11 @@ description = "CircuitPython-compatible audio synthesis and effects for PyDevice readme = "README.md" requires-python = ">=3.11" # audiovorbis compiles Xiph's Tremor into the extension, under its own -# BSD-3-Clause licence, whose text travels with the wheel. -license = "MIT AND BSD-3-Clause" -license-files = ["LICENSE", "NOTICE", "third_party/tremor/COPYING"] +# BSD-3-Clause licence, and audiofeatures compiles TensorFlow Lite Micro's +# microfrontend (Apache-2.0) and kissfft (BSD-3-Clause). Their texts travel +# with the wheel. +license = "MIT AND BSD-3-Clause AND Apache-2.0" +license-files = ["LICENSE", "NOTICE", "third_party/tremor/COPYING", "third_party/microfrontend/LICENSE", "third_party/kissfft/COPYING"] authors = [{name = "PyDevices contributors"}] classifiers = [ "Programming Language :: Python :: 3", @@ -41,7 +43,7 @@ Repository = "https://github.com/PyDevices/audiodsp" Issues = "https://github.com/PyDevices/audiodsp/issues" [tool.setuptools] -py-modules = ["audiocore", "synthio", "audiomixer", "audiofilters", "audiodelays", "audiofreeverb", "audiospeed", "audiodynamics", "audioroute", "audiomath", "audioecho", "audioconvolve", "audiobiquad", "audioshaper", "audioladder", "audiooscillator", "audioverb", "audiomodal", "audiosampler", "audiometer", "audiopump", "audiovorbis", "audiograph"] +py-modules = ["audiocore", "synthio", "audiomixer", "audiofilters", "audiodelays", "audiofreeverb", "audiospeed", "audiodynamics", "audioroute", "audiomath", "audioecho", "audioconvolve", "audiobiquad", "audioshaper", "audioladder", "audiooscillator", "audioverb", "audiomodal", "audiosampler", "audiometer", "audiofeatures", "audiopump", "audiovorbis", "audiograph"] # Keep py-modules explicit: src/cpython is the root package directory now, and # an automatic scan there would also pick up anything else that lands in it. packages = ["audiorender", "audiodsp_util"] diff --git a/setup.py b/setup.py index 03c4730e..00a4d0ef 100644 --- a/setup.py +++ b/setup.py @@ -75,8 +75,12 @@ # Xiph's Tremor, included into this one file from # third_party/tremor; docs/reference/vorbis.md. "src/shared/audiodsp_vorbis.c", + # TensorFlow Lite Micro's microfrontend and kissfft, included + # into this one file from third_party; audiofeatures' engine. + "src/shared/audiodsp_microfrontend.c", + "src/shared/audiodsp_features.c", ], - include_dirs=["src"], + include_dirs=["src", "third_party/microfrontend", "third_party/kissfft"], define_macros=[("AUDIODSP_VERSION", '"%s"' % VERSION), ("AUDIODSP_REVISION", '"%s"' % REVISION)], extra_compile_args=MACOS_COMPILE_ARGS, diff --git a/src/audiofeatures/module.c b/src/audiofeatures/module.c new file mode 100644 index 00000000..4dfb6622 --- /dev/null +++ b/src/audiofeatures/module.c @@ -0,0 +1,269 @@ +// audiofeatures: speech features for a wake word model. audiodsp's own; the +// engine is shared/audiodsp_features.c. +// +// fe = audiofeatures.Frontend() +// fe.feed(pcm, 48000, channel_count=2) # from Python, anywhere +// fe.attach(tap, 48000) # an audiopump.Tap, read in C +// n = fe.readinto(buf) # whole 10 ms slices of 40 int8 +// +// As with audiometer.Meter, the two ways in differ in where the work runs. +// feed() runs it in the caller. A tap is drained in C each time readinto() +// is called, so the work lands on the reader's thread and the pump pays +// nothing. +// +// SPDX-License-Identifier: MIT + +#include + +#include "py/obj.h" +#include "py/runtime.h" + +#include "audiopump/audiopump_tap.h" +#include "cp_compat/audiodsp_build.h" +#include "shared/audiodsp_features.h" + +#define AUDIOFEATURES_DEFAULT_SLICES 100 + +typedef struct { + mp_obj_base_t base; + audiodsp_features_t *f; + mp_obj_t tap; // MP_OBJ_NULL when no tap is attached + uint32_t tap_cursor, tap_rate, tap_channels; + uint32_t lapped; +} audiofeatures_frontend_obj_t; + +static audiofeatures_frontend_obj_t *audiofeatures_get(mp_obj_t self_in) { + audiofeatures_frontend_obj_t *self = MP_OBJ_TO_PTR(self_in); + if (self->f == NULL) { + mp_raise_ValueError(MP_ERROR_TEXT("frontend is deinitialized")); + } + return self; +} + +// Frontend(*, slices=100) +static mp_obj_t audiofeatures_frontend_make_new(const mp_obj_type_t *type, size_t n_args, size_t n_kw, + const mp_obj_t *all_args) { + enum { ARG_slices }; + static const mp_arg_t allowed[] = { + { MP_QSTR_slices, MP_ARG_INT | MP_ARG_KW_ONLY, { .u_int = AUDIOFEATURES_DEFAULT_SLICES } }, + }; + mp_arg_val_t a[MP_ARRAY_SIZE(allowed)]; + mp_arg_parse_all_kw_array(n_args, n_kw, all_args, MP_ARRAY_SIZE(allowed), allowed, a); + mp_int_t slices = a[ARG_slices].u_int; + if (slices < 8 || slices > 100000) { + mp_raise_ValueError(MP_ERROR_TEXT("slices 8..100000")); + } + audiofeatures_frontend_obj_t *self = mp_obj_malloc_with_finaliser(audiofeatures_frontend_obj_t, type); + self->f = NULL; + self->tap = MP_OBJ_NULL; + self->lapped = 0; + self->f = audiodsp_features_new((uint32_t)slices); + if (self->f == NULL) { + mp_raise_msg(&mp_type_MemoryError, MP_ERROR_TEXT("audiofeatures: no memory for the frontend")); + } + return MP_OBJ_FROM_PTR(self); +} + +static void audiofeatures_check_rate(mp_int_t rate) { + if (rate < AUDIODSP_FEATURES_MIN_RATE || rate > AUDIODSP_FEATURES_MAX_RATE) { + mp_raise_ValueError(MP_ERROR_TEXT("sample_rate 16 kHz..384 kHz")); + } +} + +static void audiofeatures_run(audiofeatures_frontend_obj_t *self, const int16_t *pcm, size_t frames, + uint32_t channels, uint32_t rate) { + if (audiodsp_features_feed_s16(self->f, pcm, frames, channels, rate) != 0) { + mp_raise_msg(&mp_type_MemoryError, MP_ERROR_TEXT("audiofeatures: no memory for the resampler")); + } +} + +// feed(buffer, sample_rate, channel_count=1) -> slices waiting +static mp_obj_t audiofeatures_frontend_feed(size_t n_args, const mp_obj_t *args) { + audiofeatures_frontend_obj_t *self = audiofeatures_get(args[0]); + if (self->tap != MP_OBJ_NULL) { + mp_raise_ValueError(MP_ERROR_TEXT("frontend is attached to a source; detach() first")); + } + mp_buffer_info_t buf; + mp_get_buffer_raise(args[1], &buf, MP_BUFFER_READ); + mp_int_t rate = mp_obj_get_int(args[2]); + mp_int_t ch = n_args > 3 ? mp_obj_get_int(args[3]) : 1; + audiofeatures_check_rate(rate); + if (ch < 1 || ch > 8 || buf.len % (2 * ch)) { + mp_raise_ValueError(MP_ERROR_TEXT("whole s16 frames, 1..8 channels")); + } + audiofeatures_run(self, (const int16_t *)buf.buf, buf.len / (2 * ch), (uint32_t)ch, (uint32_t)rate); + return MP_OBJ_NEW_SMALL_INT(audiodsp_features_waiting(self->f)); +} +static MP_DEFINE_CONST_FUN_OBJ_VAR_BETWEEN(audiofeatures_frontend_feed_obj, 3, 4, audiofeatures_frontend_feed); + +// attach(tap, sample_rate) +static mp_obj_t audiofeatures_frontend_attach(mp_obj_t self_in, mp_obj_t src, mp_obj_t rate_in) { + audiofeatures_frontend_obj_t *self = audiofeatures_get(self_in); + self->tap = MP_OBJ_NULL; + if (!mp_obj_is_type(src, &audiopump_tap_type)) { + mp_raise_TypeError(MP_ERROR_TEXT("attach() takes an audiopump.Tap")); + } + mp_int_t rate = mp_obj_get_int(rate_in); + audiofeatures_check_rate(rate); + self->tap = src; + self->tap_rate = (uint32_t)rate; + self->tap_channels = audiopump_tap_frame_bytes(src) / 2; + self->tap_cursor = audiopump_tap_position(src); + return mp_const_none; +} +static MP_DEFINE_CONST_FUN_OBJ_3(audiofeatures_frontend_attach_obj, audiofeatures_frontend_attach); + +static mp_obj_t audiofeatures_frontend_detach(mp_obj_t self_in) { + audiofeatures_get(self_in)->tap = MP_OBJ_NULL; + return mp_const_none; +} +static MP_DEFINE_CONST_FUN_OBJ_1(audiofeatures_frontend_detach_obj, audiofeatures_frontend_detach); + +// Reads the tap until `want` slices wait or it has nothing new. A chunk is +// at most 1024 bytes, which at 16 kHz mono is 32 ms: a few slices at most, so +// the loop stops while the queue still has room for them and never drops one. +static void audiofeatures_drain(audiofeatures_frontend_obj_t *self, uint32_t want) { + if (self->tap == MP_OBJ_NULL) { + return; + } + uint8_t chunk[1024]; + while (audiodsp_features_waiting(self->f) < want && audiodsp_features_room(self->f) >= 5) { + bool lapped; + uint32_t got = audiopump_tap_read_since(self->tap, &self->tap_cursor, chunk, sizeof(chunk), &lapped); + if (lapped) { + self->lapped++; + } + if (got == 0) { + break; + } + audiofeatures_run(self, (const int16_t *)(void *)chunk, got / (2 * self->tap_channels), + self->tap_channels, self->tap_rate); + } +} + +static mp_obj_t audiofeatures_read(mp_obj_t self_in, mp_obj_t buf_in, bool raw) { + audiofeatures_frontend_obj_t *self = audiofeatures_get(self_in); + mp_buffer_info_t buf; + mp_get_buffer_raise(buf_in, &buf, MP_BUFFER_WRITE); + size_t width = raw ? 2 * AUDIODSP_FEATURES_SIZE : AUDIODSP_FEATURES_SIZE; + size_t max = buf.len / width; + audiofeatures_drain(self, max > UINT32_MAX ? UINT32_MAX : (uint32_t)max); + size_t n = raw ? audiodsp_features_read_u16(self->f, (uint16_t *)buf.buf, max) + : audiodsp_features_read_s8(self->f, (int8_t *)buf.buf, max); + return MP_OBJ_NEW_SMALL_INT(n); +} + +// readinto(buf) -> slices: 40 int8 each, ready for a model +static mp_obj_t audiofeatures_frontend_readinto(mp_obj_t self_in, mp_obj_t buf_in) { + return audiofeatures_read(self_in, buf_in, false); +} +static MP_DEFINE_CONST_FUN_OBJ_2(audiofeatures_frontend_readinto_obj, audiofeatures_frontend_readinto); + +// readinto_raw(buf) -> slices: 40 uint16 each, the frontend's own values +static mp_obj_t audiofeatures_frontend_readinto_raw(mp_obj_t self_in, mp_obj_t buf_in) { + return audiofeatures_read(self_in, buf_in, true); +} +static MP_DEFINE_CONST_FUN_OBJ_2(audiofeatures_frontend_readinto_raw_obj, audiofeatures_frontend_readinto_raw); + +static mp_obj_t audiofeatures_frontend_reset(mp_obj_t self_in) { + audiofeatures_frontend_obj_t *self = audiofeatures_get(self_in); + audiodsp_features_reset(self->f); + self->lapped = 0; + if (self->tap != MP_OBJ_NULL) { + self->tap_cursor = audiopump_tap_position(self->tap); + } + return mp_const_none; +} +static MP_DEFINE_CONST_FUN_OBJ_1(audiofeatures_frontend_reset_obj, audiofeatures_frontend_reset); + +static mp_obj_t audiofeatures_frontend_stats(mp_obj_t self_in) { + audiofeatures_frontend_obj_t *self = audiofeatures_get(self_in); + audiodsp_features_stats_t st; + audiodsp_features_stats(self->f, &st); + mp_obj_t d = mp_obj_new_dict(8); + #define PUT(k, v) mp_obj_dict_store(d, MP_OBJ_NEW_QSTR(k), (v)) + PUT(MP_QSTR_slices, mp_obj_new_int_from_uint(st.slices)); + PUT(MP_QSTR_dropped, mp_obj_new_int_from_uint(st.dropped)); + PUT(MP_QSTR_waiting, mp_obj_new_int_from_uint(st.waiting)); + PUT(MP_QSTR_capacity, mp_obj_new_int_from_uint(st.capacity)); + PUT(MP_QSTR_rate, mp_obj_new_int_from_uint(st.rate)); + PUT(MP_QSTR_resampling, mp_obj_new_bool(st.resampling)); + PUT(MP_QSTR_lapped, mp_obj_new_int_from_uint(self->lapped)); + PUT(MP_QSTR_source, self->tap != MP_OBJ_NULL ? MP_OBJ_NEW_QSTR(MP_QSTR_tap) : mp_const_none); + #undef PUT + return d; +} +static MP_DEFINE_CONST_FUN_OBJ_1(audiofeatures_frontend_stats_obj, audiofeatures_frontend_stats); + +static mp_obj_t audiofeatures_frontend_deinit(mp_obj_t self_in) { + audiofeatures_frontend_obj_t *self = MP_OBJ_TO_PTR(self_in); + if (self->f != NULL) { + self->tap = MP_OBJ_NULL; + audiodsp_features_free(self->f); + self->f = NULL; + } + return mp_const_none; +} +static MP_DEFINE_CONST_FUN_OBJ_1(audiofeatures_frontend_deinit_obj, audiofeatures_frontend_deinit); + +static mp_obj_t audiofeatures_frontend_exit(size_t n_args, const mp_obj_t *args) { + (void)n_args; + return audiofeatures_frontend_deinit(args[0]); +} +static MP_DEFINE_CONST_FUN_OBJ_VAR_BETWEEN(audiofeatures_frontend_exit_obj, 4, 4, audiofeatures_frontend_exit); + +static mp_obj_t audiofeatures_frontend_enter(mp_obj_t self_in) { + return self_in; +} +static MP_DEFINE_CONST_FUN_OBJ_1(audiofeatures_frontend_enter_obj, audiofeatures_frontend_enter); + +static const mp_rom_map_elem_t audiofeatures_frontend_locals_table[] = { + { MP_ROM_QSTR(MP_QSTR_feed), MP_ROM_PTR(&audiofeatures_frontend_feed_obj) }, + { MP_ROM_QSTR(MP_QSTR_attach), MP_ROM_PTR(&audiofeatures_frontend_attach_obj) }, + { MP_ROM_QSTR(MP_QSTR_detach), MP_ROM_PTR(&audiofeatures_frontend_detach_obj) }, + { MP_ROM_QSTR(MP_QSTR_readinto), MP_ROM_PTR(&audiofeatures_frontend_readinto_obj) }, + { MP_ROM_QSTR(MP_QSTR_readinto_raw), MP_ROM_PTR(&audiofeatures_frontend_readinto_raw_obj) }, + { MP_ROM_QSTR(MP_QSTR_reset), MP_ROM_PTR(&audiofeatures_frontend_reset_obj) }, + { MP_ROM_QSTR(MP_QSTR_stats), MP_ROM_PTR(&audiofeatures_frontend_stats_obj) }, + { MP_ROM_QSTR(MP_QSTR_deinit), MP_ROM_PTR(&audiofeatures_frontend_deinit_obj) }, + { MP_ROM_QSTR(MP_QSTR___del__), MP_ROM_PTR(&audiofeatures_frontend_deinit_obj) }, + { MP_ROM_QSTR(MP_QSTR___enter__), MP_ROM_PTR(&audiofeatures_frontend_enter_obj) }, + { MP_ROM_QSTR(MP_QSTR___exit__), MP_ROM_PTR(&audiofeatures_frontend_exit_obj) }, +}; +static MP_DEFINE_CONST_DICT(audiofeatures_frontend_locals, audiofeatures_frontend_locals_table); + +MP_DEFINE_CONST_OBJ_TYPE( + audiofeatures_frontend_type, + MP_QSTR_Frontend, + MP_TYPE_FLAG_NONE, + make_new, audiofeatures_frontend_make_new, + locals_dict, &audiofeatures_frontend_locals + ); + +// quantize(value): the int8 a model takes for one raw value +static mp_obj_t audiofeatures_quantize(mp_obj_t value) { + mp_int_t v = mp_obj_get_int(value); + if (v < 0 || v > 65535) { + mp_raise_ValueError(MP_ERROR_TEXT("value 0..65535")); + } + return MP_OBJ_NEW_SMALL_INT(audiodsp_features_quantize((uint16_t)v)); +} +static MP_DEFINE_CONST_FUN_OBJ_1(audiofeatures_quantize_obj, audiofeatures_quantize); + +static const mp_rom_map_elem_t audiofeatures_module_globals_table[] = { + { MP_ROM_QSTR(MP_QSTR___name__), MP_ROM_QSTR(MP_QSTR_audiofeatures) }, + AUDIODSP_BUILD_GLOBALS, + { MP_ROM_QSTR(MP_QSTR_Frontend), MP_ROM_PTR(&audiofeatures_frontend_type) }, + { MP_ROM_QSTR(MP_QSTR_FEATURE_SIZE), MP_ROM_INT(AUDIODSP_FEATURES_SIZE) }, + { MP_ROM_QSTR(MP_QSTR_SAMPLE_RATE), MP_ROM_INT(AUDIODSP_FEATURES_RATE) }, + { MP_ROM_QSTR(MP_QSTR_STEP_MS), MP_ROM_INT(AUDIODSP_FEATURES_STEP_MS) }, + { MP_ROM_QSTR(MP_QSTR_quantize), MP_ROM_PTR(&audiofeatures_quantize_obj) }, +}; +static MP_DEFINE_CONST_DICT(audiofeatures_module_globals, audiofeatures_module_globals_table); + +const mp_obj_module_t audiofeatures_module = { + .base = { &mp_type_module }, + .globals = (mp_obj_dict_t *)&audiofeatures_module_globals, +}; + +MP_REGISTER_MODULE(MP_QSTR_audiofeatures, audiofeatures_module); diff --git a/src/cpython/_audiodsp.c b/src/cpython/_audiodsp.c index c5d1bfa8..a50dda1c 100644 --- a/src/cpython/_audiodsp.c +++ b/src/cpython/_audiodsp.c @@ -34,6 +34,7 @@ #include "shared/audiodsp_modal.h" #include "shared/audiodsp_sampler.h" #include "shared/audiodsp_meter.h" +#include "shared/audiodsp_features.h" #include "shared/audiodsp_filter_f32.h" #include "shared/audiodsp_vorbis.h" #include "shared/audiodsp_graph.h" @@ -127,6 +128,7 @@ typedef struct { PyObject *suboctave_state_type; PyObject *convolver_state_type; PyObject *meter_state_type; + PyObject *features_state_type; PyObject *vorbis_stream_type; PyObject *tank_state_type; PyObject *modal_state_type; @@ -3402,6 +3404,114 @@ static PyType_Spec meter_state_spec = { .slots = meter_state_slots, }; +// --- FeaturesState: audiofeatures' engine (shared/audiodsp_features.c) ------ + +typedef struct { + PyObject_HEAD + audiodsp_features_t *f; +} audiodsp_features_object_t; + +static int features_state_init(audiodsp_features_object_t *self, PyObject *args, PyObject *kwargs) { + unsigned int slices; + if (!PyArg_ParseTuple(args, "I:FeaturesState", &slices)) return -1; + if (slices < 8 || slices > 100000) { + PyErr_SetString(PyExc_ValueError, "slices 8..100000"); + return -1; + } + audiodsp_features_free(self->f); + self->f = audiodsp_features_new(slices); + if (self->f == NULL) { + PyErr_NoMemory(); + return -1; + } + return 0; +} + +static void features_state_dealloc(audiodsp_features_object_t *self) { + PyTypeObject *type = Py_TYPE(self); + audiodsp_features_free(self->f); + self->f = NULL; + type->tp_free((PyObject *)self); + Py_DECREF(type); +} + +static PyObject *features_state_feed(audiodsp_features_object_t *self, PyObject *args) { + Py_buffer pcm = {0}; + unsigned int rate, channels = 1; + if (!PyArg_ParseTuple(args, "y*I|I:feed", &pcm, &rate, &channels)) return NULL; + if (rate < AUDIODSP_FEATURES_MIN_RATE || rate > AUDIODSP_FEATURES_MAX_RATE) { + PyBuffer_Release(&pcm); + PyErr_SetString(PyExc_ValueError, "sample_rate 16 kHz..384 kHz"); + return NULL; + } + if (channels < 1 || channels > 8 || pcm.len % (2 * (Py_ssize_t)channels)) { + PyBuffer_Release(&pcm); + PyErr_SetString(PyExc_ValueError, "whole s16 frames, 1..8 channels"); + return NULL; + } + int rc = audiodsp_features_feed_s16(self->f, (const int16_t *)pcm.buf, + (size_t)(pcm.len / (2 * (Py_ssize_t)channels)), channels, rate); + PyBuffer_Release(&pcm); + if (rc != 0) { + return PyErr_NoMemory(); + } + return PyLong_FromUnsignedLong(audiodsp_features_waiting(self->f)); +} + +// read(buf, raw): whole slices into a writable buffer, int8 or (raw) uint16 +static PyObject *features_state_read(audiodsp_features_object_t *self, PyObject *args) { + Py_buffer out = {0}; + int raw; + if (!PyArg_ParseTuple(args, "w*p:read", &out, &raw)) return NULL; + size_t width = raw ? 2 * AUDIODSP_FEATURES_SIZE : AUDIODSP_FEATURES_SIZE; + size_t max = (size_t)out.len / width; + size_t n = raw ? audiodsp_features_read_u16(self->f, (uint16_t *)out.buf, max) + : audiodsp_features_read_s8(self->f, (int8_t *)out.buf, max); + PyBuffer_Release(&out); + return PyLong_FromSize_t(n); +} + +static PyObject *features_state_reset(audiodsp_features_object_t *self, PyObject *unused) { + audiodsp_features_reset(self->f); + Py_RETURN_NONE; +} + +static PyObject *features_state_stats(audiodsp_features_object_t *self, PyObject *unused) { + audiodsp_features_stats_t st; + audiodsp_features_stats(self->f, &st); + return Py_BuildValue("(kkkkkO)", (unsigned long)st.slices, (unsigned long)st.dropped, + (unsigned long)st.waiting, (unsigned long)st.capacity, (unsigned long)st.rate, + st.resampling ? Py_True : Py_False); +} + +AUDIODSP_LOCKED_ARGS(features_state_feed, audiodsp_features_object_t) +AUDIODSP_LOCKED_ARGS(features_state_read, audiodsp_features_object_t) +AUDIODSP_LOCKED_ARGS(features_state_reset, audiodsp_features_object_t) +AUDIODSP_LOCKED_ARGS(features_state_stats, audiodsp_features_object_t) +static PyMethodDef features_state_methods[] = { + {"feed", (PyCFunction)features_state_feed_locked, METH_VARARGS, NULL}, + {"read", (PyCFunction)features_state_read_locked, METH_VARARGS, NULL}, + {"reset", (PyCFunction)features_state_reset_locked, METH_NOARGS, NULL}, + {"stats", (PyCFunction)features_state_stats_locked, METH_NOARGS, NULL}, + {NULL, NULL, 0, NULL}, +}; + +AUDIODSP_LOCKED_INIT(features_state_init, audiodsp_features_object_t) +static PyType_Slot features_state_slots[] = { + {Py_tp_new, PyType_GenericNew}, + {Py_tp_init, features_state_init_locked}, + {Py_tp_dealloc, features_state_dealloc}, + {Py_tp_methods, features_state_methods}, + {0, NULL}, +}; + +static PyType_Spec features_state_spec = { + .name = "_audiodsp.FeaturesState", + .basicsize = sizeof(audiodsp_features_object_t), + .flags = Py_TPFLAGS_DEFAULT | Py_TPFLAGS_HEAPTYPE, + .slots = features_state_slots, +}; + // --- Graph: audiograph's runner (shared/audiodsp_graph.c) ------------------ // // The same calls as the MicroPython binding (src/audiograph/Graph.c), with @@ -5596,6 +5706,14 @@ static int audiodsp_exec(PyObject *module) { state->vorbis_stream_type) < 0) return -1; if (PyModule_AddIntConstant(module, "METER_MAX_BANDS", AUDIODSP_METER_MAX_BANDS) < 0) return -1; + state->features_state_type = PyType_FromModuleAndSpec(module, + &features_state_spec, NULL); + if (state->features_state_type == NULL) return -1; + if (PyModule_AddObjectRef(module, "FeaturesState", + state->features_state_type) < 0) return -1; + if (PyModule_AddIntConstant(module, "FEATURES_SIZE", AUDIODSP_FEATURES_SIZE) < 0) return -1; + if (PyModule_AddIntConstant(module, "FEATURES_RATE", AUDIODSP_FEATURES_RATE) < 0) return -1; + if (PyModule_AddIntConstant(module, "FEATURES_STEP_MS", AUDIODSP_FEATURES_STEP_MS) < 0) return -1; state->tank_state_type = PyType_FromModuleAndSpec(module, &tank_state_spec, NULL); if (state->tank_state_type == NULL) return -1; @@ -5668,6 +5786,7 @@ static int audiodsp_traverse(PyObject *module, visitproc visit, void *arg) { Py_VISIT(state->suboctave_state_type); Py_VISIT(state->convolver_state_type); Py_VISIT(state->meter_state_type); + Py_VISIT(state->features_state_type); Py_VISIT(state->vorbis_stream_type); Py_VISIT(state->tank_state_type); Py_VISIT(state->modal_state_type); @@ -5695,6 +5814,7 @@ static int audiodsp_clear(PyObject *module) { Py_CLEAR(state->suboctave_state_type); Py_CLEAR(state->convolver_state_type); Py_CLEAR(state->meter_state_type); + Py_CLEAR(state->features_state_type); Py_CLEAR(state->vorbis_stream_type); Py_CLEAR(state->tank_state_type); Py_CLEAR(state->modal_state_type); diff --git a/src/cpython/audiofeatures.py b/src/cpython/audiofeatures.py new file mode 100644 index 00000000..583619bf --- /dev/null +++ b/src/cpython/audiofeatures.py @@ -0,0 +1,93 @@ +"""Speech features for a wake word model, from PCM at any rate. + +audiodsp's own, not a CircuitPython module. ``Frontend`` is TensorFlow Lite +Micro's audio microfrontend with the settings microWakeWord's models were +trained on (ESPHome's): a 30 ms window every 10 ms, 40 mel channels from +125 Hz to 7.5 kHz, noise reduction, PCAN gain control and a log scale. The +arithmetic is the same C on every interpreter (``src/shared/audiodsp_features.c``), +so a desktop and a board make the same bytes from the same samples. + + fe = audiofeatures.Frontend() + fe.feed(pcm, 48000, 2) # interleaved s16, any rate from 16 kHz + buf = bytearray(40 * 10) + n = fe.readinto(buf) # n whole slices of 40 int8 + +Audio that isn't 16 kHz is resampled to 16 kHz first, through a low-pass that +is flat to 7.5 kHz; 16 kHz audio goes to the frontend untouched. Each 10 ms +makes a slice, which waits in a queue of ``slices`` until it is read. When +the queue is full the oldest slice is dropped, and ``stats()`` counts it. + +On a board the frontend can also listen without Python touching a sample: +``attach(tap, sample_rate)`` reads an ``audiopump.Tap`` in C. There is no pump +to attach to on CPython, so ``attach()`` raises there. +""" + +import _audiodsp + +#: The audiodsp this was built from (src/cp_compat/audiodsp_build.h). +__version__ = _audiodsp.__version__ +__revision__ = _audiodsp.__revision__ + +#: Values in a slice. +FEATURE_SIZE = _audiodsp.FEATURES_SIZE +#: The rate the frontend itself runs at. +SAMPLE_RATE = _audiodsp.FEATURES_RATE +#: Audio per slice, in milliseconds. +STEP_MS = _audiodsp.FEATURES_STEP_MS + + +def quantize(value): + """The int8 a model takes for one raw value: value * 256 / 666 - 128, rounded.""" + if not 0 <= value <= 65535: + raise ValueError("value 0..65535") + v = (value * 256 + 333) // 666 - 128 + return max(-128, min(127, v)) + + +class Frontend: + """Feature slices from PCM, ten a second, held until read.""" + + def __init__(self, *, slices=100): + self._state = _audiodsp.FeaturesState(slices) + + def _live(self): + if self._state is None: + raise ValueError("frontend is deinitialized") + return self._state + + def feed(self, buffer, sample_rate, channel_count=1): + """Interleaved s16 frames at any rate from 16 to 384 kHz. Returns the slices waiting.""" + return self._live().feed(buffer, sample_rate, channel_count) + + def attach(self, source, sample_rate=None): + self._live() + raise NotImplementedError("attach() needs a board's pump: on CPython, feed() the samples") + + def detach(self): + self._live() + + def readinto(self, buf): + """Whole slices of 40 int8 into ``buf``, oldest first. Returns how many.""" + return self._live().read(buf, False) + + def readinto_raw(self, buf): + """Whole slices of 40 uint16 (the frontend's own values) into ``buf``. Returns how many.""" + return self._live().read(buf, True) + + def reset(self): + """Forget the audio so far and every queued slice.""" + self._live().reset() + + def stats(self): + slices, dropped, waiting, capacity, rate, resampling = self._live().stats() + return {"slices": slices, "dropped": dropped, "waiting": waiting, "capacity": capacity, + "rate": rate, "resampling": resampling, "lapped": 0, "source": None} + + def deinit(self): + self._state = None + + def __enter__(self): + return self + + def __exit__(self, *exc): + self.deinit() diff --git a/src/shared/audiodsp_features.c b/src/shared/audiodsp_features.c new file mode 100644 index 00000000..91aea553 --- /dev/null +++ b/src/shared/audiodsp_features.c @@ -0,0 +1,434 @@ +// audiofeatures' engine; audiodsp_features.h says what it does. +// SPDX-License-Identifier: MIT +// SPDX-FileCopyrightText: Copyright (c) 2026 PyDevices + +#include "shared/audiodsp_features.h" + +#include +#include +#include +#include + +#include "shared/audiodsp_microfrontend_names.h" +#include "tensorflow/lite/experimental/microfrontend/lib/frontend.h" +#include "tensorflow/lite/experimental/microfrontend/lib/frontend_util.h" + +// --- the resampler ----------------------------------------------------------- +// +// A polyphase windowed-sinc filter. The kernel is sinc(u) under a Kaiser +// window (beta 6) that is RS_ZEROS output samples wide on each side, where u +// counts 16 kHz output samples: its cutoff is 8 kHz, it is flat (within +// 0.01 dB) to 7.5 kHz, and from 8.5 kHz up it is at least 62 dB down. What it +// lets through between 8 and 8.5 kHz folds to between 7.5 and 8 kHz, above +// the frontend's top channel, so nothing it aliases reaches a feature. +// +// The ratio is held exactly as L/M (16000/rate in lowest terms), so the +// position of every output sample is exact. When the L rows of taps fit in +// RS_MAX_EXACT floats (every integer ratio, and 24 kHz) they are worked out +// once; otherwise (44.1 kHz and its family) each output's row comes from a +// table of the kernel, linearly interpolated. Each row is scaled to sum to 1, +// so a constant input comes out unchanged. + +#define RS_ZEROS 32 +#define RS_BETA 6.0 +#define RS_TABLE_STEPS 64 // kernel table points per output sample +#define RS_MAX_EXACT 4096 +#define RS_BLOCK 256 // input frames taken at a time +#define RS_PI 3.14159265358979323846 + +typedef struct { + uint32_t rate, L, M; + int J0; // taps before the centre sample + int T; // taps in a row + float *rows; // L rows of T, or NULL + float *table; // RS_ZEROS * RS_TABLE_STEPS + 2 points, or NULL + float *row; // the current row, when rows is NULL + float *buf; // mono input still needed + int cap, n; + int k0; // the next output's centre is buf[k0] + ph / L + uint32_t ph; +} rs_t; + +static double rs_bessel_i0(double x) { + double sum = 1.0, term = 1.0, q = x * x / 4.0; + for (int k = 1; k < 64; k++) { + term *= q / ((double)k * (double)k); + sum += term; + if (term < sum * 1e-17) { + break; + } + } + return sum; +} + +// The kernel at u output samples from the centre. +static double rs_kernel(double u) { + u = fabs(u); + if (u >= RS_ZEROS) { + return 0.0; + } + double r = u / RS_ZEROS; + double w = rs_bessel_i0(RS_BETA * sqrt(1.0 - r * r)) / rs_bessel_i0(RS_BETA); + double s = u < 1e-12 ? 1.0 : sin(RS_PI * u) / (RS_PI * u); + return s * w; +} + +static uint32_t rs_gcd(uint32_t a, uint32_t b) { + while (b) { + uint32_t t = a % b; + a = b; + b = t; + } + return a; +} + +static void rs_free(rs_t *rs) { + free(rs->rows); + free(rs->table); + free(rs->row); + free(rs->buf); + memset(rs, 0, sizeof(*rs)); +} + +// u for tap j of phase ph, in output samples. +static double rs_u(const rs_t *rs, uint32_t ph, int j) { + return ((double)ph + (double)(rs->J0 - j) * rs->L) / rs->M; +} + +static void rs_normalize(float *row, int T) { + double sum = 0.0; + for (int j = 0; j < T; j++) { + sum += (double)row[j]; + } + if (sum != 0.0) { + float g = (float)(1.0 / sum); + for (int j = 0; j < T; j++) { + row[j] *= g; + } + } +} + +static void rs_restart(rs_t *rs) { + // History before the first sample is silence, so the first output is + // centred on the first input sample. + memset(rs->buf, 0, (size_t)rs->cap * sizeof(float)); + rs->n = rs->J0; + rs->k0 = rs->J0; + rs->ph = 0; +} + +static bool rs_configure(rs_t *rs, uint32_t rate) { + rs_free(rs); + uint32_t g = rs_gcd(AUDIODSP_FEATURES_RATE, rate); + rs->rate = rate; + rs->L = AUDIODSP_FEATURES_RATE / g; + rs->M = rate / g; + double W = (double)RS_ZEROS * rs->M / rs->L; // the kernel's half-width in input samples + rs->J0 = (int)ceil(W); + rs->T = 2 * rs->J0 + 2; + rs->cap = rs->T + RS_BLOCK; + rs->buf = malloc((size_t)rs->cap * sizeof(float)); + if (rs->buf == NULL) { + goto fail; + } + if ((size_t)rs->L * (size_t)rs->T <= RS_MAX_EXACT) { + rs->rows = malloc((size_t)rs->L * (size_t)rs->T * sizeof(float)); + if (rs->rows == NULL) { + goto fail; + } + for (uint32_t ph = 0; ph < rs->L; ph++) { + float *row = rs->rows + (size_t)ph * rs->T; + for (int j = 0; j < rs->T; j++) { + row[j] = (float)rs_kernel(rs_u(rs, ph, j)); + } + rs_normalize(row, rs->T); + } + } else { + int pts = RS_ZEROS * RS_TABLE_STEPS + 2; + rs->table = malloc((size_t)pts * sizeof(float)); + rs->row = malloc((size_t)rs->T * sizeof(float)); + if (rs->table == NULL || rs->row == NULL) { + goto fail; + } + for (int i = 0; i < pts; i++) { + rs->table[i] = (float)rs_kernel((double)i / RS_TABLE_STEPS); + } + } + rs_restart(rs); + return true; +fail: + rs_free(rs); + return false; +} + +static const float *rs_row(rs_t *rs) { + if (rs->rows != NULL) { + return rs->rows + (size_t)rs->ph * rs->T; + } + float step = (float)rs->L / (float)rs->M * RS_TABLE_STEPS; + float x = ((float)rs->ph / (float)rs->M + (float)rs->J0 * (float)rs->L / (float)rs->M) * RS_TABLE_STEPS; + for (int j = 0; j < rs->T; j++, x -= step) { + float a = fabsf(x); + int i = (int)a; + float v = 0.0f; + if (i < RS_ZEROS * RS_TABLE_STEPS) { + float fr = a - (float)i; + v = rs->table[i] + (rs->table[i + 1] - rs->table[i]) * fr; + } + rs->row[j] = v; + } + rs_normalize(rs->row, rs->T); + return rs->row; +} + +// --- the engine -------------------------------------------------------------- + +struct audiodsp_features { + struct FrontendConfig config; + struct FrontendState fe; + bool fe_ready; + rs_t rs; + uint32_t rate; // last rate fed + uint16_t *queue; // capacity slices + uint32_t capacity, head, count; + uint32_t slices, dropped; +}; + +int8_t audiodsp_features_quantize(uint16_t value) { + int32_t v = (((int32_t)value * 256) + 333) / 666 - 128; + return (int8_t)(v < -128 ? -128 : v > 127 ? 127 : v); +} + +static void features_config(struct FrontendConfig *c) { + // ESPHome's micro_wake_word settings, which microWakeWord's models were + // trained with. + FrontendFillConfigWithDefaults(c); + c->window.size_ms = AUDIODSP_FEATURES_WINDOW_MS; + c->window.step_size_ms = AUDIODSP_FEATURES_STEP_MS; + c->filterbank.num_channels = AUDIODSP_FEATURES_SIZE; + c->filterbank.lower_band_limit = 125.0f; + c->filterbank.upper_band_limit = 7500.0f; + c->noise_reduction.smoothing_bits = 10; + c->noise_reduction.even_smoothing = 0.025f; + c->noise_reduction.odd_smoothing = 0.06f; + c->noise_reduction.min_signal_remaining = 0.05f; + c->pcan_gain_control.enable_pcan = 1; + c->pcan_gain_control.strength = 0.95f; + c->pcan_gain_control.offset = 80.0f; + c->pcan_gain_control.gain_bits = 21; + c->log_scale.enable_log = 1; + c->log_scale.scale_shift = 6; +} + +audiodsp_features_t *audiodsp_features_new(uint32_t capacity) { + if (capacity == 0) { + capacity = 1; + } + audiodsp_features_t *f = calloc(1, sizeof(*f)); + if (f == NULL) { + return NULL; + } + f->capacity = capacity; + f->queue = malloc((size_t)capacity * AUDIODSP_FEATURES_SIZE * sizeof(uint16_t)); + features_config(&f->config); + if (f->queue == NULL) { + audiodsp_features_free(f); + return NULL; + } + memset(&f->fe, 0, sizeof(f->fe)); + if (!FrontendPopulateState(&f->config, &f->fe, AUDIODSP_FEATURES_RATE)) { + // What it allocated before failing is freed with it. + f->fe_ready = true; + audiodsp_features_free(f); + return NULL; + } + f->fe_ready = true; + return f; +} + +void audiodsp_features_free(audiodsp_features_t *f) { + if (f == NULL) { + return; + } + if (f->fe_ready) { + FrontendFreeStateContents(&f->fe); + } + rs_free(&f->rs); + free(f->queue); + free(f); +} + +void audiodsp_features_reset(audiodsp_features_t *f) { + FrontendReset(&f->fe); + if (f->rs.buf != NULL) { + rs_restart(&f->rs); + } + f->head = f->count = 0; + f->slices = f->dropped = 0; +} + +static void features_push(audiodsp_features_t *f, const uint16_t *values) { + if (f->count == f->capacity) { + f->head = (f->head + 1) % f->capacity; + f->count--; + f->dropped++; + } + uint32_t tail = (f->head + f->count) % f->capacity; + memcpy(f->queue + (size_t)tail * AUDIODSP_FEATURES_SIZE, values, AUDIODSP_FEATURES_SIZE * sizeof(uint16_t)); + f->count++; + f->slices++; +} + +// 16 kHz mono samples through the frontend. +static void features_run(audiodsp_features_t *f, const int16_t *x, size_t n) { + while (n > 0) { + size_t read = 0; + struct FrontendOutput out = FrontendProcessSamples(&f->fe, x, n, &read); + if (read == 0) { + break; + } + x += read; + n -= read; + if (out.size == AUDIODSP_FEATURES_SIZE) { + features_push(f, out.values); + } + } +} + +static int16_t features_round(float y) { + if (y >= 32767.0f) { + return 32767; + } + if (y <= -32768.0f) { + return -32768; + } + return (int16_t)(y >= 0.0f ? (int32_t)(y + 0.5f) : -(int32_t)(0.5f - y)); +} + +int audiodsp_features_feed_s16(audiodsp_features_t *f, const int16_t *pcm, size_t frames, + uint32_t channels, uint32_t rate) { + if (rate < AUDIODSP_FEATURES_MIN_RATE || rate > AUDIODSP_FEATURES_MAX_RATE || channels == 0) { + return -1; + } + if (rate != f->rate) { + if (rate == AUDIODSP_FEATURES_RATE) { + rs_free(&f->rs); + } else if (!rs_configure(&f->rs, rate)) { + f->rate = 0; + return -1; + } + f->rate = rate; + } + int16_t mono[RS_BLOCK]; + if (rate == AUDIODSP_FEATURES_RATE) { + if (channels == 1) { + features_run(f, pcm, frames); + return 0; + } + while (frames > 0) { + size_t m = frames < RS_BLOCK ? frames : RS_BLOCK; + for (size_t i = 0; i < m; i++) { + int32_t s = 0; + for (uint32_t c = 0; c < channels; c++) { + s += pcm[i * channels + c]; + } + int32_t h = (int32_t)channels / 2; + mono[i] = (int16_t)((s >= 0 ? s + h : s - h) / (int32_t)channels); + } + features_run(f, mono, m); + pcm += m * channels; + frames -= m; + } + return 0; + } + rs_t *rs = &f->rs; + float scale = 1.0f / (float)channels; + while (frames > 0) { + int room = rs->cap - rs->n; + size_t m = frames < (size_t)room ? frames : (size_t)room; + float *dst = rs->buf + rs->n; + for (size_t i = 0; i < m; i++) { + int32_t s = 0; + for (uint32_t c = 0; c < channels; c++) { + s += pcm[i * channels + c]; + } + dst[i] = (float)s * scale; + } + rs->n += (int)m; + pcm += m * channels; + frames -= m; + // Every output whose taps are all here. + size_t out = 0; + while (rs->k0 - rs->J0 + rs->T <= rs->n) { + const float *row = rs_row(rs); + const float *x = rs->buf + rs->k0 - rs->J0; + float y = 0.0f; + for (int j = 0; j < rs->T; j++) { + y += row[j] * x[j]; + } + mono[out++] = features_round(y); + if (out == RS_BLOCK) { + features_run(f, mono, out); + out = 0; + } + rs->ph += rs->M; + rs->k0 += (int)(rs->ph / rs->L); + rs->ph %= rs->L; + } + features_run(f, mono, out); + // Keep only what the next output needs. + int drop = rs->k0 - rs->J0; + if (drop > 0) { + if (drop > rs->n) { + drop = rs->n; + } + memmove(rs->buf, rs->buf + drop, (size_t)(rs->n - drop) * sizeof(float)); + rs->n -= drop; + rs->k0 -= drop; + } + } + return 0; +} + +uint32_t audiodsp_features_room(const audiodsp_features_t *f) { + return f->capacity - f->count; +} + +uint32_t audiodsp_features_waiting(const audiodsp_features_t *f) { + return f->count; +} + +size_t audiodsp_features_read_u16(audiodsp_features_t *f, uint16_t *out, size_t max_slices) { + size_t n = 0; + while (n < max_slices && f->count > 0) { + memcpy(out + n * AUDIODSP_FEATURES_SIZE, f->queue + (size_t)f->head * AUDIODSP_FEATURES_SIZE, + AUDIODSP_FEATURES_SIZE * sizeof(uint16_t)); + f->head = (f->head + 1) % f->capacity; + f->count--; + n++; + } + return n; +} + +size_t audiodsp_features_read_s8(audiodsp_features_t *f, int8_t *out, size_t max_slices) { + size_t n = 0; + while (n < max_slices && f->count > 0) { + const uint16_t *v = f->queue + (size_t)f->head * AUDIODSP_FEATURES_SIZE; + for (int i = 0; i < AUDIODSP_FEATURES_SIZE; i++) { + out[n * AUDIODSP_FEATURES_SIZE + i] = audiodsp_features_quantize(v[i]); + } + f->head = (f->head + 1) % f->capacity; + f->count--; + n++; + } + return n; +} + +void audiodsp_features_stats(const audiodsp_features_t *f, audiodsp_features_stats_t *st) { + st->slices = f->slices; + st->dropped = f->dropped; + st->waiting = f->count; + st->capacity = f->capacity; + st->rate = f->rate; + st->resampling = f->rate != 0 && f->rate != AUDIODSP_FEATURES_RATE; +} diff --git a/src/shared/audiodsp_features.h b/src/shared/audiodsp_features.h new file mode 100644 index 00000000..083cdcd1 --- /dev/null +++ b/src/shared/audiodsp_features.h @@ -0,0 +1,68 @@ +// audiofeatures' engine: speech features for a wake word model, from PCM at +// any rate. +// SPDX-License-Identifier: MIT +// SPDX-FileCopyrightText: Copyright (c) 2026 PyDevices +// +// Samples go in as interleaved s16 frames at any rate from 16 kHz up. They are +// mixed to mono and, unless they are already 16 kHz, resampled to 16 kHz +// through a windowed-sinc low-pass that is flat to 7.5 kHz and at least 60 dB +// down from 8.5 kHz. 16 kHz input skips the resampler, so the samples reach the +// frontend untouched. +// +// The frontend is TensorFlow Lite Micro's microfrontend +// (audiodsp_microfrontend.c) with the settings microWakeWord's models were +// trained on: a 30 ms window every 10 ms, 40 mel channels from 125 Hz to +// 7.5 kHz, noise reduction, PCAN gain control and a log scale. Each 10 ms of +// audio makes one slice of 40 values, which wait in a queue until read, either +// raw (uint16, about 0 to 670) or as the int8 a model takes. + +#pragma once + +#include +#include + +#define AUDIODSP_FEATURES_SIZE 40 // values per slice +#define AUDIODSP_FEATURES_RATE 16000 // the frontend's own rate +#define AUDIODSP_FEATURES_STEP_MS 10 // audio per slice +#define AUDIODSP_FEATURES_WINDOW_MS 30 +#define AUDIODSP_FEATURES_MIN_RATE 16000 +#define AUDIODSP_FEATURES_MAX_RATE 384000 + +typedef struct audiodsp_features audiodsp_features_t; + +typedef struct { + uint32_t slices; // slices made since creation or reset + uint32_t dropped; // slices lost because the queue was full + uint32_t waiting; // slices in the queue now + uint32_t capacity; // the queue's size in slices + uint32_t rate; // the input rate last fed, 0 before any + uint8_t resampling; // 1 when that rate goes through the resampler +} audiodsp_features_stats_t; + +// A frontend whose queue holds up to `capacity` slices; NULL when out of memory. +audiodsp_features_t *audiodsp_features_new(uint32_t capacity); +void audiodsp_features_free(audiodsp_features_t *f); + +// Forget all audio and queued slices, as if new. +void audiodsp_features_reset(audiodsp_features_t *f); + +// Interleaved s16 frames. A rate other than the last one fed starts the +// resampler over. Returns 0, or -1 when the rate is out of range or the +// resampler couldn't get memory. When the queue is full the oldest slice goes. +int audiodsp_features_feed_s16(audiodsp_features_t *f, const int16_t *pcm, size_t frames, + uint32_t channels, uint32_t rate); + +// How many more slices the queue can take before it drops one. +uint32_t audiodsp_features_room(const audiodsp_features_t *f); +uint32_t audiodsp_features_waiting(const audiodsp_features_t *f); + +// Up to max_slices of the oldest queued slices, removed from the queue. +// Returns how many were written: AUDIODSP_FEATURES_SIZE values each. +size_t audiodsp_features_read_u16(audiodsp_features_t *f, uint16_t *out, size_t max_slices); +size_t audiodsp_features_read_s8(audiodsp_features_t *f, int8_t *out, size_t max_slices); + +// The int8 a microWakeWord model takes for one raw value: ESPHome's +// value * 256 / 666 - 128, rounded, in integers. +int8_t audiodsp_features_quantize(uint16_t value); + +void audiodsp_features_stats(const audiodsp_features_t *f, audiodsp_features_stats_t *st); diff --git a/src/shared/audiodsp_microfrontend.c b/src/shared/audiodsp_microfrontend.c new file mode 100644 index 00000000..2cd43b13 --- /dev/null +++ b/src/shared/audiodsp_microfrontend.c @@ -0,0 +1,167 @@ +// TensorFlow Lite Micro's audio microfrontend, compiled for every audiodsp +// target in this one translation unit. +// SPDX-License-Identifier: MIT +// SPDX-FileCopyrightText: Copyright (c) 2026 PyDevices +// +// The library lives unchanged in third_party/microfrontend, and the FFT it +// calls, kissfft, in third_party/kissfft (each with an ORIGIN.txt). They are +// included here rather than listed in each build, as Tremor is in +// audiodsp_vorbis.c, for the same reasons: one file for every build system to +// add, every external symbol renamed (audiodsp_microfrontend_names.h) so a +// firmware that links another copy still links, and the library's warnings +// silenced for its code only. +// +// Upstream builds its FFT from three C++ files (fft.cc, fft_util.cc and +// kiss_fft_int16.cc), which put kissfft in a namespace and compile it as +// 16-bit fixed point. They are carried in third_party but not compiled: this +// file is their C equivalent, statement for statement, with kissfft compiled +// as FIXED_POINT 16 the way kiss_fft_int16.cc does it, so the arithmetic is +// the same and no target needs a C++ compiler. The only difference is that +// the failure messages fft_util.cc prints are left out. +// +// The build adds third_party/microfrontend and third_party/kissfft to the +// include path, because the library names its headers by their full +// "tensorflow/lite/..." path. + +#include +#include +#include +#include +#include +#include +#include + +#include "shared/audiodsp_microfrontend_names.h" + +#if defined(_MSC_VER) +#pragma warning(push, 0) +#elif defined(__clang__) +#pragma clang diagnostic push +#pragma clang diagnostic ignored "-Weverything" +#elif defined(__GNUC__) +#pragma GCC diagnostic push +#pragma GCC diagnostic ignored "-Wpragmas" +#pragma GCC diagnostic ignored "-Wunused-variable" +#pragma GCC diagnostic ignored "-Wunused-but-set-variable" +#pragma GCC diagnostic ignored "-Wunused-parameter" +#pragma GCC diagnostic ignored "-Wunused-function" +#pragma GCC diagnostic ignored "-Wsign-compare" +#pragma GCC diagnostic ignored "-Wtype-limits" +#pragma GCC diagnostic ignored "-Wcast-align" +#pragma GCC diagnostic ignored "-Wmissing-field-initializers" +#pragma GCC diagnostic ignored "-Wold-style-definition" +#pragma GCC diagnostic ignored "-Wstrict-prototypes" +#pragma GCC diagnostic ignored "-Wmissing-prototypes" +#pragma GCC diagnostic ignored "-Wmissing-declarations" +#pragma GCC diagnostic ignored "-Wredundant-decls" +#pragma GCC diagnostic ignored "-Wshadow" +#pragma GCC diagnostic ignored "-Wconversion" +#pragma GCC diagnostic ignored "-Wsign-conversion" +#pragma GCC diagnostic ignored "-Wfloat-conversion" +#pragma GCC diagnostic ignored "-Wdouble-promotion" +#pragma GCC diagnostic ignored "-Wfloat-equal" +#pragma GCC diagnostic ignored "-Wpedantic" +#pragma GCC diagnostic ignored "-Wbad-function-cast" +#pragma GCC diagnostic ignored "-Wimplicit-fallthrough" +#endif + +// --- kissfft, 16-bit fixed point (kiss_fft_int16.cc) ----------------------- + +#define FIXED_POINT 16 +#include "kiss_fft.c" +#include "tools/kiss_fftr.c" + +// --- the FFT stage (fft.cc and fft_util.cc) --------------------------------- + +#include "tensorflow/lite/experimental/microfrontend/lib/fft.h" +#include "tensorflow/lite/experimental/microfrontend/lib/fft_util.h" + +void FftCompute(struct FftState *state, const int16_t *input, int input_scale_shift) { + const size_t input_size = state->input_size; + const size_t fft_size = state->fft_size; + int16_t *fft_input = state->input; + // First, scale the input by the given shift. + size_t i; + for (i = 0; i < input_size; ++i) { + fft_input[i] = (int16_t)((uint16_t)input[i] << input_scale_shift); + } + // Zero out whatever else remains in the top part of the input. + for (; i < fft_size; ++i) { + fft_input[i] = 0; + } + // Apply the FFT. + kiss_fftr((kiss_fftr_cfg)state->scratch, state->input, (kiss_fft_cpx *)state->output); +} + +void FftInit(struct FftState *state) { + // All the initialization is done in FftPopulateState(). + (void)state; +} + +void FftReset(struct FftState *state) { + memset(state->input, 0, state->fft_size * sizeof(*state->input)); + memset(state->output, 0, (state->fft_size / 2 + 1) * sizeof(*state->output)); +} + +int FftPopulateState(struct FftState *state, size_t input_size) { + state->input_size = input_size; + state->fft_size = 1; + while (state->fft_size < state->input_size) { + state->fft_size <<= 1; + } + state->input = (int16_t *)malloc(state->fft_size * sizeof(*state->input)); + if (state->input == NULL) { + return 0; + } + state->output = (struct complex_int16_t *)malloc((state->fft_size / 2 + 1) * sizeof(*state->output) * 2); + if (state->output == NULL) { + return 0; + } + // Ask kissfft how much memory it wants. + size_t scratch_size = 0; + kiss_fftr_cfg kfft_cfg = kiss_fftr_alloc((int)state->fft_size, 0, NULL, &scratch_size); + if (kfft_cfg != NULL) { + return 0; + } + state->scratch = malloc(scratch_size); + if (state->scratch == NULL) { + return 0; + } + state->scratch_size = scratch_size; + // Let kissfft configure the scratch space just allocated. + kfft_cfg = kiss_fftr_alloc((int)state->fft_size, 0, state->scratch, &scratch_size); + if (kfft_cfg != state->scratch) { + return 0; + } + return 1; +} + +void FftFreeStateContents(struct FftState *state) { + free(state->input); + free(state->output); + free(state->scratch); +} + +// --- the rest of the library, as it is --------------------------------------- + +#include "tensorflow/lite/experimental/microfrontend/lib/filterbank.c" +#include "tensorflow/lite/experimental/microfrontend/lib/filterbank_util.c" +#include "tensorflow/lite/experimental/microfrontend/lib/frontend.c" +#include "tensorflow/lite/experimental/microfrontend/lib/frontend_util.c" +#include "tensorflow/lite/experimental/microfrontend/lib/log_lut.c" +#include "tensorflow/lite/experimental/microfrontend/lib/log_scale.c" +#include "tensorflow/lite/experimental/microfrontend/lib/log_scale_util.c" +#include "tensorflow/lite/experimental/microfrontend/lib/noise_reduction.c" +#include "tensorflow/lite/experimental/microfrontend/lib/noise_reduction_util.c" +#include "tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control.c" +#include "tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control_util.c" +#include "tensorflow/lite/experimental/microfrontend/lib/window.c" +#include "tensorflow/lite/experimental/microfrontend/lib/window_util.c" + +#if defined(_MSC_VER) +#pragma warning(pop) +#elif defined(__clang__) +#pragma clang diagnostic pop +#elif defined(__GNUC__) +#pragma GCC diagnostic pop +#endif diff --git a/src/shared/audiodsp_microfrontend_names.h b/src/shared/audiodsp_microfrontend_names.h new file mode 100644 index 00000000..f08fa8a7 --- /dev/null +++ b/src/shared/audiodsp_microfrontend_names.h @@ -0,0 +1,62 @@ +// Every external symbol the microfrontend and kissfft define, renamed into +// audiodsp's namespace. +// SPDX-License-Identifier: MIT +// SPDX-FileCopyrightText: Copyright (c) 2026 PyDevices +// +// audiodsp_microfrontend.c includes this before any of the library, and so +// does every file that calls it, so declarations, definitions and calls all +// come out as audiodsp_mf_*. Without it, a firmware that also links another +// copy of the frontend or of kissfft fails to link on the first duplicate. +// +// The list is every global symbol in the compiled translation unit. +// tests/test_cpython_audiofeatures.py fails on any symbol that escaped, so an +// update that adds a function shows up there. + +#pragma once + +#define FftCompute audiodsp_mf_FftCompute +#define FftFreeStateContents audiodsp_mf_FftFreeStateContents +#define FftInit audiodsp_mf_FftInit +#define FftPopulateState audiodsp_mf_FftPopulateState +#define FftReset audiodsp_mf_FftReset +#define FilterbankAccumulateChannels audiodsp_mf_FilterbankAccumulateChannels +#define FilterbankConvertFftComplexToEnergy audiodsp_mf_FilterbankConvertFftComplexToEnergy +#define FilterbankFillConfigWithDefaults audiodsp_mf_FilterbankFillConfigWithDefaults +#define FilterbankFreeStateContents audiodsp_mf_FilterbankFreeStateContents +#define FilterbankPopulateState audiodsp_mf_FilterbankPopulateState +#define FilterbankReset audiodsp_mf_FilterbankReset +#define FilterbankSqrt audiodsp_mf_FilterbankSqrt +#define FrontendFillConfigWithDefaults audiodsp_mf_FrontendFillConfigWithDefaults +#define FrontendFreeStateContents audiodsp_mf_FrontendFreeStateContents +#define FrontendPopulateState audiodsp_mf_FrontendPopulateState +#define FrontendProcessSamples audiodsp_mf_FrontendProcessSamples +#define FrontendReset audiodsp_mf_FrontendReset +#define LogScaleApply audiodsp_mf_LogScaleApply +#define LogScaleFillConfigWithDefaults audiodsp_mf_LogScaleFillConfigWithDefaults +#define LogScalePopulateState audiodsp_mf_LogScalePopulateState +#define NoiseReductionApply audiodsp_mf_NoiseReductionApply +#define NoiseReductionFillConfigWithDefaults audiodsp_mf_NoiseReductionFillConfigWithDefaults +#define NoiseReductionFreeStateContents audiodsp_mf_NoiseReductionFreeStateContents +#define NoiseReductionPopulateState audiodsp_mf_NoiseReductionPopulateState +#define NoiseReductionReset audiodsp_mf_NoiseReductionReset +#define PcanGainControlApply audiodsp_mf_PcanGainControlApply +#define PcanGainControlFillConfigWithDefaults audiodsp_mf_PcanGainControlFillConfigWithDefaults +#define PcanGainControlFreeStateContents audiodsp_mf_PcanGainControlFreeStateContents +#define PcanGainControlPopulateState audiodsp_mf_PcanGainControlPopulateState +#define PcanGainLookupFunction audiodsp_mf_PcanGainLookupFunction +#define PcanShrink audiodsp_mf_PcanShrink +#define WideDynamicFunction audiodsp_mf_WideDynamicFunction +#define WindowFillConfigWithDefaults audiodsp_mf_WindowFillConfigWithDefaults +#define WindowFreeStateContents audiodsp_mf_WindowFreeStateContents +#define WindowPopulateState audiodsp_mf_WindowPopulateState +#define WindowProcessSamples audiodsp_mf_WindowProcessSamples +#define WindowReset audiodsp_mf_WindowReset +#define kLogLut audiodsp_mf_kLogLut +#define kiss_fft audiodsp_mf_kiss_fft +#define kiss_fft_alloc audiodsp_mf_kiss_fft_alloc +#define kiss_fft_cleanup audiodsp_mf_kiss_fft_cleanup +#define kiss_fft_next_fast_size audiodsp_mf_kiss_fft_next_fast_size +#define kiss_fft_stride audiodsp_mf_kiss_fft_stride +#define kiss_fftr audiodsp_mf_kiss_fftr +#define kiss_fftr_alloc audiodsp_mf_kiss_fftr_alloc +#define kiss_fftri audiodsp_mf_kiss_fftri diff --git a/tests/data/audiofeatures/signal_16k.s8 b/tests/data/audiofeatures/signal_16k.s8 new file mode 100644 index 0000000000000000000000000000000000000000..2b369e3861a11daae858b16459a944e18be809f5 GIT binary patch literal 11920 zcmeG?*={0PmMmVt7_%E=Fss>Dvlz3O&HIv>&dSQFtXw*)Yr3YprJk-H%@3S^y1#Ir z<1L=zAs&#BkPt!$AwVp1Zh*70O6p8M&Pa7E!V&kLd+yoq`tHr&-u>f;53m3Br@y}W z>HXh7{OkSC@Bi!lFCTyT_}@SN^H0Bg{1?9c^5g&f`PJ34)60{Ockj;MKL7Lc|G0m+ z|NE=A7w_-?^6>HYudn{|;jhnse)!A7PY*vme^RuA4>Oy?Af$#K6)0tiaaE4Qq{ou zMk9H3a#i@IaBHl_>ZyD8jp;P>kbdqdCvURPBRAHY*xAu-tf9ClzY1I>uN92zV}vp8U5*DDSL)AAJ(5m8_KftI(BY3HJ8`R8)vTbNYz#|)b}f? z3vIDM zpv;&m=87T{jEd!r`Xa3=8F$HHE6XWoHFt5~u~j@3dwy21y9`-tu2nGF)fvn2T2YbR zDoi~j8;Sr7YMMg`t`}!k-p^ack!FsM&C{KtIk!_?OKbcPn&#Y>G_0^vp?3r!V%sax zXw$ouhJ-WQ&UEwJC&pAa01?@S<&u^{0iocH2_f-xl!Vek-%Iy)P1&wvoL8M_W0S;| zPnKw>mEL~YT2+LG2XHK;W@VX+=8379Yel`$abXe|orry??&S0?$m{*kCfd#v0vlnH zZuvBctvGckQfm<*CgeL>Q+2WqW-w0T-~rRXE84{9B*=vfrxX9o&`&h|#M>&4ha{ra zK#sy7u$|tGPvUc5FC;jLdsKiq=!l7i!~ho?&!7zp&;x*dY@GPjCp6M`4#N=TwXi>r zFZ|TnvhjA1gRUJH=3~d~5Zr^!iC38-lcAlDb=BaUS=-V64j=A$I+|ug6N2&F&~j}$ zT2`3z2D(Y)B*vKs3y5$*suSE)82c0)%IPBBL1eC{VSbJDJC1=L{Gjcai3gNm>;-72 zZx;Y2!m#Io**rV5plsfum&&Oee8R5K1y3hUET9?e4yIrMAW;rL+w}nlO*m+K=O)k) zLrX^V9TiM~w(?@6c>p8SfF9*O61L#56QX;h4?K*)j@Wvlq5w9qnaoZ9SPnf1@Sc_l zvU7LSO-HbXyeVK#oS@^h#0J7>qB)j}G{9M@6~s=+fx3bP4Cn$pwC`>pj z>bVQR@2-W$nsJEOxjAwr@JgcRjigP6O%LCJCT)R;m;)H|3$mwj(cPg&RA7k+4$>h* zIz_9YUYb3wi?|$_pe}ik*EBuf2e_umDaNio1xn zfK6DSyaX|UMEEucf&zGI;s^*>3au$Ghzqb${dgqoAouOIlJ=2uuWQFPG8iIW56EP@ zON2HVV;Rpe?#W0E3d(_1G|?V(NDnayw5U-WQ8jv;8hd8i5DBxDAZQT_TM1VDA|nTM zi54_y5A+k%9m+}C*iW-(L@dgywgQN!Uhji~sE~%7(?XX=L709Kb|+D^)x3KE6nSDj zX=sw;vEMZ!f=(o$0vJue#Lb2()B@v}LrG9zur+zE4V$e#ia)~sNs9xHDn$-Du)s2k zcx(I`ECeF)Gz0^f#>h3NYj490T4;N8Ll|^9eji#Kp#*i`<@zuaMC7rJT=bAmyVXUB zMtlo~ZRDB=NO_9<3*bQWQf$K`6@YEzsSmwY8`~@KyF+-&BgX_aD#M;o60nb$e zq5=$CO?>NPxrn|l(hca7-C=|YDV0qm4?`ZwbR0$sBm6P{-vp7fGzr=Uno}NGfIgj2 zpfjvhz<{hH_ao$)c^?WOSR{=J_^L&ZU>q=-ebh<^d2aF!3^^Y0MYIm{XVh-1gHPtW zoFp&fkP_M^Oi9S)CicRlo*)JGUxrMe1@Ldie$GZdhP)9h-ET==hFmC-5|b6Jg;j0z zi|*Y&-@Rr90_qkiG2oCpIoflI??mc@Bki02e)ksRghLMiFfj8YFkOhi^8`6sw)Nfo z^ZO5d>=31)(4!bD66k<(6vSRj(nt%!SAYNJ1GWXBiH$`ND3v(tBZue%QyD%1{XBi~ z^RxFzf%z-zVQClvBOfCP`5d6Ms4J3h{B-*k`2m=T0`ksq*U%E!k~Kv>*y#SF?pF(9mr(9D?Q*h(-sq!XR(Im0%I9<%^^p9LLnY`NZ%m6&PpSQ zidMo$lM2`a#qg1M54dk5XRML|cG!_RMU8;aW}G1!nyyV9yBD6yiX}M{jm1QBFk9-8 ziNcA@-`(p~Q@%4~)DX}nBgQ%FtBXG*PFm~{-;SU1Dio!~e@h$u3Qk)%~g z1pHd+n4q=DD&~e{j7|T#Q=2 zOLQ1FOmNgaC6dx5Q?hCake(o*3$!p>ktV}zA_B$y&OZjjBi^zBdlIxD<|5=oqWDR% z_$3oLj^zU&A}^)w8RMAw6^nKW=%9%_~qX%=17* zhNP|-#n15I+czIE4pG8Y(`7vcHNLP0j}+g&{xg1c@Uy*SS*j4x#-Dv6#@m-aF3@^X z%LzUax-BVASUh3cT$aJx`;XlL;~(h^P!v?P1k4vvfsfGq;pT%!!&*@`TBF1a0msru zOFYth-S|OLorHuKl5a}>j(&}Otnj(#rEk#q$Qg5D{1bpjimT4i2CJZfo6J3PU+fu2 z6f^}dJNd4(~q-h_{bFN2dO8oT*#{ za*d0t>rA!w&E;+Nu71-vJvqO)EB(?y92{ui_<7=4y;83nU*xLCDc@24dEu^jo2#ET z@|T6{<9}SBcB_b2zB42<)#Gxy5v?Ab1C)0l&mEGhIkOiU-sM>v)&eu zk87by;xuv@itQgqYq6^NDpGQt9{!AY)P;BvpEcAvi&UzmXnm0jMv^6CVLzE{92P8v zgHrS-#@`zLRE63_xV&|dt;I8{f?MZEb6LZYy^%UGp1Ml0_sAPgg#e7zY}{BosM?FR zl)acvsS5-9EUGQ~a?Yf!=KL;;IDY(X1P#q@DVSUYG{tbjUDz%K^UCa2&XRJab~Eb_ zF~ngWFtZcwMot-=K(QhAtbGFcoxt}nX4F8bZE_OT`<9B}KLTJW|_?afYbw#)it0}zk7Pd$Ys z@5l9{a}fnM@9cOc-ic>q6s@?%MDVHc5oPX5BV#X^MxjLn*UYxb0c`a{qusU|*G~wW zn28q|svWd&b`xlM2IH+5#CB}r8gcNtj&~r6k{E2!*mvXHz}i%`_01;EU}tjgun$}O z4!_H95kX3TY(;_4jx6FW!*m~4xBv&bDGtvgve!4@N^wUbH)es$_hcqqh3#RyM_LQm zYi|wJWM`mg3q;8|9y1AmTu0f;HN_3VPIX;us~J zq_2<91)6wsRY~1=Zd13_=cPCG>*~YlRrH~9Q+R$_&wh^)-Y#T3F$$4;e=&4be3`vV zCrsDL%cJYueeN#*P3~C@_q~d^jz1|Ul3jbcZamLbGdG3XU^;Z;eUWU$Z_C$_>%>|9 zK3fT4Ix}A#yY!jki{l4h{-}|yZ)NkP!pnH!B(Mwk* za$aoMGNVkq;kyZ!`{|9O`|7CfD|<1>H`1xiwMnGBkqw-c&n@LbR(7DW_^P2=yyU4S z(&0k96hysnmT;8wh5Y0sb>|2{Y_lONY}Z2P-kL9eP`6~~CD)M?gdQ>B=)pA$HqXIn zFE#Wzj~55F3xC;|*)8oB%w<_kd&(RS_aOQ(ZFx3S+Eg#D4Eys&pE@R<*lVhiA-R*+ zRh(&DKoHVY?t*ffFvn)@c4|%Et2-ilB}H~8Z_C;fnxkfZFT^qsqz_s zHfBpe+Vu}{2gGq=iSt}Hvwdvy%6v#UI8Vrf@}pjKn!~!L zOAn4G$+dlw>V*ZrHahW3eyYI4W;F-q5u>YBG1|KBW*k$Zh4YdK__t#Ev?VJ=h7s3I zYMrf@lcgcv;Z@<0m!vccF0z?a>PFeMJ2aQ$K4zs;mHP< z>V>#M*V)wIkVkOOjaY+*1F#VVaZ()@;XVBn=jE-fBTOj<@b`S{=vQn;Km!;7EzYF~ zVH3BeL7X4#f~K1`xw$rP>pEbQET%phVv9PE;|xHB^BR9VqVnL3V!f9d8T(4u zZQ5~t-|Lwsd(6B;bXHug%!!}cnhriBhJ+dJgKcoOu3a#Ig(`SIF0Vy4GZ05S)_`u} z;qx5@?#}2TGH@WSx39x56N_`)LtsfdYSLiuoVdA~4{U&ipNu`KZJgPRe=|2(uH0G6 zFGKFj)sMCG^6~A<|0nmW1;v-qK>tGa6^v!I<3axuY-%S`S{mCkF(LIfDB|$LBkpI(>#R@y8!eYmKKM$bu^(?OgZp~}*NS17}Ek{PY& zj}HP+O+MuB+<9lE_Uz;~T#4KUUu5pm*Rh+_E!1L%h_^-0-6bd|j>^xDZhWOi<3+uA zpS(&HkLt0r$W^ogiOmG!WgYcO74n7r-HXOqxt=-8yN=GvHJ>+d8c!!1sRD^h2PXl2 zGMk1%H1%-jNEL4i;ZPtSjYW=p!BR2e%N-|SZeX)qN3@tRfh(hafm?R9awN8|btlV7 zcP!@f#Y5R*4$d;nPRVjLf%wYi!?i!A55@cmzdxK!W{<-DSSsXnhu}18L_D>9pZM)&#B1y{S@r;wPSn zGiLVS)V0g_bET#RO4dT&Y9y8n?Zr`Co!stbv}S?zp!?Xa?#sGH&SUj)`_|x|Q{&J%Cc!Zn zg3CRa)CL|FaB!qI=~eUThI216(0BDc!(m8n?zl$^o0aKhap2WzVK)SEnRD118a6k# z!$U$nqa8kCHMD)|o#iCFml*4(!0AX(9uV7o_1HSJwatovz5@o8YpiamFPvi|^3=t81_llc=7G2CIf(Z+2YXYM=xQ6*c1Bu$Dw_M8 zY=uXc8;cs^9rolZvfZ(^%tH^@C*6UO7>6pYwCoftdfPg10lX)|EqvSv+p9ZwbTuidh?sGu4uTSMFx)Q{G1WzNio z8x3rOUplm$W2Vt-R4^SVyUMnL7{wkts$Ky)1Mb-2_5tc0uLzCx!bqHQ9i+l_CU}1E z2mGOdKQ!?FSOcFT{mU{`U+gdE#65ynu-L9Ef?k>hV2VRvC3wvE&Fe$y6Ancq$XIS2 zvTF(KCJx(u!qY!|iqL)g>X;=D6TY*fcuX_+DPk?V1jH^;NS8PO-aq7*^x0LLM^Irh z5)g-791-xGPypfGeM7dYW|yuQERS~CMGki0$aXq)h<9CINVRMwpjctDI@SxrW4!i& zzB#}A{iKeKi_7fJE@?5Um=a%>Dt0VoxVR}RK?g66tzMBTj(Y6lE}dJpUmnFW8QFm@ z<5JeciQ#H|7OtYw?&WfwGLud`m+gRTNhh^S(xmQ=u%xV6Nr11PK)SMmUy}Aq<}eOX z2uq#di$q{4+maTePmYxiwuisOk~EDMUVeq84N&gnaF3OK_SD%u#*!_jy0S;-(M>Dk8 zWjBWM8Cnl{LD^XMkH(U&bFo8x29KzC*$H3qoLagh!BBtoz2Cn7D#@~(@%$Oz5|1oy x*u{oN>i8+?nO&z^!snj0as;gS#NsF0vEF`#N8L*{lZdm&vX{?zzyiD$^uO;fn4bUu literal 0 HcmV?d00001 diff --git a/tests/data/audiofeatures/signal_16k.u16 b/tests/data/audiofeatures/signal_16k.u16 new file mode 100644 index 0000000000000000000000000000000000000000..7167898a3ea3d16fa8f8c7d9421bfd3a8b4e9580 GIT binary patch literal 23840 zcmeI4d2m(NeaG)PiQ@!2POU6f65Gn6EF_K6sEyT-MoE;Yp;fEWA){(&h=yTQ28O6m z7&O2D!y*jEfZ4-p79k@W2oDf4ViB9!9w2tJNr25}Teh){jf|7HpYJ{I=sw9fGag?4 z)6pEAd(ZOQ&-y#}{n-A%{;v{+ ziT`2K?1RJ-`wjbuH9$?QLpd)Yi^s?14KZT@Hun$Os+<{>-6oHrZI{WfllU1aXH`_0>?&OBgWH0id{ zWZ4??Dcb;^lhARVSOvZ-rn5a@zCioC?M7hiGM#9D-K;TxGCj>Lvm5@Kpz#H0xM)6Y zhf`J$+_?FiT>_3#rrB(u-C=MnGcTDv?(5um@k{*@Uq98Gjq+CY_+**cF}G#`m6=la`0xFBc>dDPt)!Mb>rp@ z6N8@~(7(mBm?B`zgU5PkA4k10=9=k2>*XHkI`;`Sz=14JBs#G%mq_sUNb9rE;JXx^PZ`s{&L=% zgG;0h)VbSEhKA1OUGpe)2BMe8$X`#}B6HM?G_TX{Rb+D1tTS6pOCUPmHU-q_1%Jb7 zRS`Ume2$xVP-OCfz1UnbSzwwJF`1av-6{R?IUIi<*ItT}Cpp}9L{ zeW;gbmICvHsWF@2EeNuZdtdVm&ufv>C15>fnu3+Yv8D^svgQ=e+e}wzN~e4@?*ox} zG4FlNYTB&_=JQC#+E`Evp2tiqXbv*5$ROANuBm1ca)||7fVCVN4wJvj%!i&d_|Aa- z{@CmYG?`)kKv^B~?1}CK?+UVtAO)l2c$2W&1mv9cgI1KfS^eT{mHpz|y=kA&8f$Z<0* zHUskluwvkC0)J=PTt^0{z`ls*M$(m}X|yXtMtMB{0S=z!{V*_AAScnW6WFuC8%K5< z!MTmtNc}Qse-rpUXg8WTg6Ey2ONk@V^>T1dCT=Ew9d=M-hGJ8>$UYy9cLufve??FQ zuFKGx&GSb1?@Rkr;FONnAmhVGX(#QrqswY!wjTWRsQV;*U8HVT;vn8f(q;qr)_|)X zc=_a~5etarM1|?(4+qvNZ0tPsH$umJFh5KA7VKpJ=``@{Ku4=7p9pPJ$e%%eH{xLQ zu^s$ffinvXvcJXfS^~`F=xHBx_gNW%X2?)jwi~-#u0;{BG?*ID;hGOJqSd*CH+0J zsshg_p8L{%D0m7<%YoOMyj)@|=)_pt5^O?->I*qU*}xD|$>$(6^dWtQHuS{camr-> zP2lK4I+NHhtX~oIqFiALc~ikU8knoluk@Wsj3K*9Xc)=!ap+>iU^I(`IE%@X+`1FH z0&hO)VE7uq^JJczsg;Y4PXnz5s%4uKXe<4;AmtjgAv}YK($DxXzcYO+2nO;#le&_> zg7Apa+z>{X70Ov$ji+y}hvvS-^|YHo9rX*zs+sa4^2*_31gY#J8y+foPp901gP9C^ zT*%_N3uOh=n}ozQwmgcS`@(_vR8ZYHl&6QVVuAFkzNU~1>`d@yB;dx;wdfFD$zQgt z_c(teNTmWgkF!HE72|JOl(~ytSIGkofqhu&}$}0&| zr{4!b$`k$$dN3>kLoo+xtDzf6402@9hl;M0{loWYA;zU6-?3B ziTn|i=a8z6Y+tE(RuFxXo$puTbwC1^%0WBv6olPq029RnDR^L-w$w?A)`jq`R zexzmTPay^uy_qXwybE3_vw~#O9T@7fqN5w~E#o-@nL3R)8XbfH8_Q*kbp5I;a#FvJ z(N1#mJoTe|4{&4y8q?KhB)dV#K{638jRvwW=`-JChhv8H=VdIOT#u3r1Ny0jqIk$x z3H`h?FdY=fe+9MFCJr@4;L)@CjrxazcoDo}plCcOz&aAq(Ec;S*t*Y*JEicCQl2H6I);)Ddpskt{M3^rXRqbT3XoQBou0#M zh*Nb12a5^uX9HKVQJ(N92$ouC-0T_VOHSFW?Yt(#YGFNAIU$WIadNC+jOb!&X{UF= zQ7XLNo}|gKU=_XUH?myycP9nWp|MvIiG_MniN=a7TJ)fvN(5JO7TpTYBAUZF&)b-+ zK|!{km8d7Vdpc!rYA>6ToMnx(Lko5?6@1aE_6idTw~DY&i5Kyupgtpb{=JZA!4aL}KvDMU()sLJ&aqY{6N;6my3?S^Ye#vqF3mln z$?H|QlABtIZ|OaQ9Y%lqW82xPRG~Rxw~i*4=V+1aG)HKzkPRz{H&3@@kj74K5b*xd zy7@#rXbeyohV3Q?Zx^Dwi9f+oNS1B2xre>`_wCPY=LCmNH{v-H8xar5VKDC+!ByT8 z_L|5>q~mAVLH^c$o&B}bk@%H7gjsYn(^74{%%m?-t@0+gmTnxCZfT9%AP`P9$I8MIm6G5=Vg7W}amK zzlr^~Z|PIVVR|^uc$i-2X<@EvIBxY{wj5)8eiRKaW#6tgk(9a_uQfkUSXSBj;aIAg z3f>m9%Fwe~FJ^we&Wcw1xa@W`x_%-pJ1RU3*@>rFNAyR$O_cj_Y2vZ zr$g|x_UjtXs(KR~!E-S6-iKax7|QArv^h%Rq0Oq?^$|U}AUeEdF`t1(?Wf9bL%{Bz|R ze?8|;ZZ6;d_ho{szNIvIp34Se7ovh_l=VrcQD{$*@VUO`Ybi~(3ZE}LD$$DHJv;o# zecREees}xe=UkO)y{UDr_J@+U=grx@sERDz+w<+9I{LJVcB4Qi;7U&w;b`hwXzo@J zU#j2K+5^JuHvc@oZ8F1`iq^SG>mp@oZRaKG>5tmEvchm}8gE@6h%&7aqt$kPYmf8e zm1yvA9lRLrGZLeU=5E#OoZ!l>Pa?c--goZ~N?TYei3j&`M*UDYPB$@EIBO9PPIoN{ z2}N8b%jJXcg5cXWo3n2>A37bjt(DwfP)ieBdkk%B)0gdcIDr?f8ZTU_#@r}GEwm0( zU%h>)Yms9!Y*$+luHd|V3WEA#6f#21`uW7qPlB5aTlF8YpT%R4!MN__@94=&XvJA| z1lji!XMs_^d$RZd?qS_g8J@E_T(NLX>}z&vo%`Irw0)mwPv3X#Entee+l6@2xEYP_ zf~he>p)LG)kZ!MXUhjIolbIVORa;T~X;-e`O2S9}xqgs5rl_P7zKYLUKl(cP3XaYS zl*(Sboa>+!GiR>Z;hD0BFWy9}@3)=~LDLzBe^-6ICBG}L5w+|C(GlRMckqpz2fRm(eA$mQ{VdW3`S4z9G1iYw9q&2}1Oj*48USWtQs@eWY z`&aFMd?)bC~n<7=GFs`-ELpE`Kc66q{_U^m%+ohE&i zy~!H>#h3zp*#jSKY$*(`5ck&HEYOe`i|~a9&HKXYjl3 z!WVjzJ&xiTV4cEWcN)J_HD1`W_-HQz>k@w4lhk+xZ08tNg2^+8H{&k3;)# z&O0W8|1uEIz+nyTo`bKA(6QcZGf$Ae6<_Bz@QlGvT!$}o5tKe{evRK=aAJ679K1UG zdygUIPSCRvoY$$l72j?FKFPKCg?5L2)C<7e&ELDgxM0Sy$C^w2W>a9_O~BEpq|l49 z8_28{e;{+6Jpmm%;l4NR7Mnq6v4J_@7(TcI@ZEy{_&oeHL*Es+`iXrbgyT+y1aCMz z90$S(;LPQDx_Qw|rqxV*sfX}$l_I|^a`T~P0A(xC|9m7j1uEvV&p(!+OCCmro^Umc z9askRucclAzP*05dltzrfyVw}{y9$5k3s)QJWO$Di?i=NigiAW58lD?V_iPDbynSz zx|w5&izlrfShE724Xe*Y=@7)j&6nL$HlU+yeo7Kqs;{P%!8gR zRu8j~k~^>8%_{6^^75eX9Pl?Ghw)|<{&k(gw+;31x*NGA+Yal#!=Urz6MZlazn-7L=jG6NHxcNFx zJ@q^JV+Ui?D{0*gIFAEU=gFDC8;9fy!}5r}X!VsIVwnnz40xGNDp_e%6L$H@#vy<8 z`k82GApL6r&)R>FCc2PGe>oM}LV)!?Xm=yCzVz2zeD({u$J59-RL6W+%{VsQ>_yA_ zDcj5(azN5U7SrJTOxRPkA8iQd&iA1|7Cgtju)FN*^w@hiXS~nuXP&#xoO^|PVegnv z+Bd-W3a8dvNMEPU1AIwzT@&$3oHX+F>0v@bHZrrQ_cT>kj;%w@Y6%QIM!G=OV3ccV@+P8|jR zMQB{d8Q})j7FFiUb}#e$4rqQE_?v+J9P{cq?p_^37BwNf97c`1?HOj3a!!v=BdgQQ zE7eHtAZ?Fvf2j}qoU@Ga7tL4f4!CKe{&plXpRsfncd)iY&kXifPlfO%Gd>Q$g7;%# zOX2wpz24e|z^3DgP8y$@;O2fiE0jn*m~H{@B6H3SMwvg+`)>fZ5e#pGL3bUlAgu$TCCbJtf>A)$ z-Xl)WyufI*jJnO(YMg$s0y&*V7IVoD>`K;HBjA5L^ZKLEH5z^jp|_0nkRS2fJqY!K zJZy6r{Pw0i4!rHm>9tU-y{tyadiub5)*f#$@?2pJvxgpXjxpjjXgh>Wi(9RJE$!t0 zFJkoAg}&<6b$3a3uHHtE8@bnYl{P1U(TI+Y(^sCOA01>gYox3Oxd@)ZKGp-9 z!MUAr<``{mpg*msmoVN{q4P`FOe2(TpzmsYoQka;0M-TMP)pt>=-A3Q;h!EK9&~h0w zwR#;ye2RKHE9gmk+2?xF3(!9w`ul-f?RCFnE%_C2tF^A?>B+Ry*ghYb=w8r9;4Z_4 z_VK)&XKtCIhak`@WhikKu&Y_`3SJKMILPmjk;7 zt5@IG96g?Pf~U|OnTqZ$tQU3HuMhUL4wnGZJ^~6=#Y0yI@dju67UvKw;K9u z!upO=-4RjH`bbb#VjUa7zm#Xu-2?gTN0%dzq2}gCX*+{dBT6P^Be6fN30=_LH1{-- zRQI!1fo}?J=0JslmL03-$)a_yLU+h|5(j~E zDevr1k*;`l_c7||C*oCbmqV|7!YUh;(4`=0r13tEEl{uWUPPU|sw}d6hie>oa_Rq- zj7id&&TP6QWF~m+gT~fWhUZfXkF&YNxqE78%w@jUYByD(oOMNE3gOLWy)h((@D{Qb zQ?E=_7{k6;JtVjm^RFzt9&yI))hc&67ufj}D=I zi}fkghG1)tlFFdDB$YvPL@I+Ht5XFtS}COpVh5ZK33N!HLjoNV=#W5%1Ue+pA%Q!S zKq^tXqD z$dZ@ytH@gt-JsR|&`b0hdB}7YxEUYTpXlK)bCM^|&rY5XkjlyqGcr{XU*`Ffc%0Xf zD!ho#Ydfd?`W1-l_?sTIjd)p?;7gu|r|T^|+vo92*YdoH(RVGU%10jP01-!+bY&G7`7XH3uZ$kCw>>j+xaeU1O88x?s zH_1L>x8fDmFIAM_>2vg*z~AvUxz~9=%bDj{eC}Vi?{Z_{Wp11urQKP)J_kv!;o*D= zuSpCq?M^&5d%?9I4{-%@h`xY{nIUj_g7P4^h;J|s$DgsY@oHU$#tV4i*5ko!#Phrj zkKX}gyAED1@vK{7XMieijD7>c(Kj2JzliTmXW0v(rxgC~;}*basJVv!ryf6J8S>c6 z+ivdpY$K-RAp=RF<&@?hZIB<+dDv+&yVg1%DHXP|Ha zoYw$?#I-!M`f{)mx-WCL@DP5?VcfGSgj)IOCPPOFkW2CP z?BFIsC308*=VN$sr>?q3;^emw9p&MN)4i)=uvc*(U=#f8hQmqBSiO<|QtI?%wwVB5 zPg1%F%F2x%{GZWF%VHltn!wYH$kalq>VEsL>%8DQ-Uw;4z|9xP>^*8Qv}xOcY}InM*u zY~(x(-`d*H-#wIe{;WO|io}hRr$&|oXm%9*^um9cN6CC(Y(h6ffIkE|4?_zz_`F9$ zpME!_JMz^%Zx?iDp$ZEtgTC?j@^$O85Gm_dFLb}BgunjCG8-LAV>xi#gr9mMWogh| z%DaB~!G}!r-vw$!pGJ`i{P_Gr2K1-Fb17vLfHjr2W1wUWH6`g`P``+>-XB>1=-z;S zzrn%p3RUvkZUxd<@^tHJD7=p2euRA2Q_+(wYX&}C-R2#hkW)sYy{B?$=&KltO!QkV z8Q8&G^x7F-vw)um9l9%*29L$?m&08beR-<;NS$dZ&3Q;lV_10f-WU4x+aQ_n-HYe` zaQIm0SI-WuY&v??7cKfdis97Ly^m<7)c4cn+z-%r9mgl_eorKh$9M&pbI5%X8uax> ze`wO3F8$t0DZ1D18+7AYH`6+W9#YZN1P*16!kje!nLj_#U21HRtJ`x3hDt zs&!vacjlfXRXc9ja>qk=(%e^qv!S~Qi;~XT-qrVS*|eS>ViIqP`n>_&V-Tj^=(33C zZ1Sb2A=EA=)!gcPo8J5vv+l;aB>gboRGNT5Rk9TMn}K!*fAiUg8Z z2Op)~?dx7=kFDSFi%NBSSt|q0_Q|+(BIVBqlg}u89i_40YkVXAbP_cEvRq5Jsw^be zoyBN(;BNH`t5#b6?vwJhQc}=LC4=8k{f0fwn#6^Xlqp}hl)BZ7E7h~wiblV_)4oD= zJ>5a@6TYx-WNi}#JxAqHzE)OR840@JsE(q(hVc3PnjpV1+rl0wYNsdR^fdc-k5hE2 z%&j@y?M*%F)?;(5?_JjbuY<5r1Zu~#`CnE`AC57(_rU7r4nW_>NwK8RMy($gFSI^B653yFfW zBlqt0D0s?OP+G{kT50=lTAKr%>YPka7s*cTWMkgWlJNN{+HZ5Je5Yaf0dnVEir2IL3j}bUcRaw<$AkU-f33-v)PF{c*(FtgA48# zUi@2fnoOhMHwC-dkxXjUhcEYOTe(G$Bp4A@-Z^%4)z`DA>-8#|QgEsJ zwL_N^={uS4INuEDbyO-UG?qjfOy(%6Cz&WZy?H6N?^)5WbE&1bOFg%TDmaq2g8QCY{J3`yHL{Vg zTI;_L%lA~>4wo8xIs3E~VEEv%M>}`flXz$k&$sjTta3jN`m#v7KF6hPPt;tV=Qp{g hhwbw`Z~i?wZEyEc+lX@i 1 else "audiofeatures" +features_module = __import__(MODULE) + +SIZE = 40 + + +def signal(quarters=6): + out = array("h") + s = 2024 + seg = 4000 + plan = [(0, 0, 300), (7721, 9000, 300), (3861, 12000, 2000), (0, 0, 9000), + (25736, 6000, 600), (38604, 3000, 4000)] + for i in range(quarters): + step, level, noise = plan[i % len(plan)] + x, y = level << 8, 0 + for n in range(seg): + env = n * 4 if n < seg // 4 else seg - (n - seg // 4) * 4 // 3 + env = max(0, min(seg, env)) + x -= (y * step) >> 16 + y += (x * step) >> 16 + s = (s * 1103515245 + 12345) & 0x7FFFFFFF + r = ((s >> 15) & 0xFFFF) - 32768 + v = ((x >> 8) * env) // seg + r * noise // 32768 + out.append(max(-32768, min(32767, v))) + return out + + +def hold(pcm, rate, channels): + """pcm (16 kHz mono) at `rate`, each output sample the input one before it.""" + out = array("h") + for i in range(len(pcm) * rate // 16000): + v = pcm[i * 16000 // rate] + for _ in range(channels): + out.append(v) + return out + + +def run(pcm, rate, channels, block): + fe = features_module.Frontend(slices=64) + raw = array("H", bytes(2 * SIZE * 64)) + got = array("H") + mv = memoryview(pcm) + step = block * channels + for at in range(0, len(pcm), step): + fe.feed(mv[at:at + step], rate, channels) + n = fe.readinto_raw(raw) + got.extend(raw[:n * SIZE]) + dropped = fe.stats()["dropped"] + fe.deinit() + return got, dropped + + +def main(): + base = signal() + cases = ( + ("16k_mono", base, 16000, 1), + ("48k_stereo", hold(base, 48000, 2), 48000, 2), + ("44k1_mono", hold(base, 44100, 1), 44100, 1), + ) + for name, pcm, rate, channels in cases: + a, dropped = run(pcm, rate, channels, 480) + b, _ = run(pcm, rate, channels, 37) + h = 2166136261 + for v in a: + h = ((h ^ v) * 16777619) & 0xFFFFFFFF + slices = len(a) // SIZE + print("%s slices %d dropped %d fnv %08x" % (name, slices, dropped, h)) + for k in range(0, slices, 25): + print("%s %d %s" % (name, k, " ".join("%d" % v for v in a[k * SIZE:(k + 1) * SIZE]))) + print("%s blocks %s" % (name, "same" if a == b else "DIFFERENT")) + print("quantize", " ".join("%d" % features_module.quantize(v) for v in (0, 1, 332, 333, 334, 665, 666, 670, 65535))) + + +main() diff --git a/tests/parity/verify_dsp.py b/tests/parity/verify_dsp.py index 61456de5..b27bd804 100644 --- a/tests/parity/verify_dsp.py +++ b/tests/parity/verify_dsp.py @@ -146,6 +146,12 @@ {"circuitpython": "audiometer is not added to CircuitPython (the meter's " "sources are the sound card's and the pump's, neither of " "which CircuitPython has)"}, None), + # audiofeatures prints feature slices, not PCM. Not built into + # CircuitPython either, for the same reason. + ("features_probe.py", "audiofeatures", + {"circuitpython": "audiofeatures is not added to CircuitPython (its " + "sources are the sound card's and the pump's, neither of " + "which CircuitPython has)"}, None), # synthio.Biquad and audiofilters.Filter are CircuitPython's, so they are # held to CircuitPython's bytes and NOT to a stored digest -- which is the # whole of audiodsp#77: graded against its own capture, this probe reported diff --git a/tests/test_cpython_audiofeatures.py b/tests/test_cpython_audiofeatures.py new file mode 100644 index 00000000..25cb1626 --- /dev/null +++ b/tests/test_cpython_audiofeatures.py @@ -0,0 +1,272 @@ +"""audiofeatures: the frontend's bytes, the resampler, and the surface. + +The bytes. tests/data/audiofeatures/signal_16k.u16 and .s8 are the features +TensorFlow Lite Micro's own microfrontend makes from ``signal()`` at 16 kHz: +the library at tflite-micro commit 230c9d9531bd9b08d445b4f2ed0760756ffd5f78, +built from its own C++ FFT files (fft.cc, fft_util.cc, kiss_fft_int16.cc) +rather than this repository's C, set up with ESPHome's micro_wake_word +settings, and quantized with ESPHome's integer scaling. They were written by +that build, not by this module, so a match is a match with the reference, not +with our own last answer. 16 kHz input skips the resampler, so the match must +be exact, in every block size. + +The resampler. Audio at another rate goes through a low-pass to 16 kHz first, +so its features can't be byte-identical; the bar is that a band-limited +signal fed at 48 kHz (and 44.1 kHz, the interpolated path) makes nearly the +features the same signal makes at 16 kHz, and that a tone above 8.5 kHz, which +the resampler must remove, adds nothing. +""" + +import os +import shutil +import subprocess +import unittest +from array import array + +import numpy as np + +import _audiodsp +import audiofeatures + +HERE = os.path.dirname(os.path.abspath(__file__)) +DATA = os.path.join(HERE, "data", "audiofeatures") +SIZE = audiofeatures.FEATURE_SIZE + + +def signal(seconds=3): + """Integer-only 16 kHz test audio: a noise floor, then tone and noise bursts. + + Integers only, so it is the same audio on every platform. Changing it + invalidates the goldens. + """ + out = array("h") + s = 2024 + seg = 16000 // 4 + # (Minsky circle step, 2*pi*f/16000 in Q16; tone level; noise level) + plan = [(0, 0, 300), (0, 0, 300), (7721, 9000, 300), (15442, 4000, 300), + (3861, 12000, 2000), (0, 0, 9000), (25736, 6000, 600), (0, 0, 300), + (11581, 15000, 300), (38604, 3000, 4000), (0, 0, 150), (5148, 8000, 1200)] + for i in range(seconds * 4): + step, level, noise = plan[i % len(plan)] + x, y = level << 8, 0 + for n in range(seg): + env = n * 4 if n < seg // 4 else seg - (n - seg // 4) * 4 // 3 + if env > seg: + env = seg + if env < 0: + env = 0 + x -= (y * step) >> 16 + y += (x * step) >> 16 + s = (s * 1103515245 + 12345) & 0x7FFFFFFF + r = ((s >> 15) & 0xFFFF) - 32768 + v = ((x >> 8) * env) // seg + r * noise // 32768 + out.append(max(-32768, min(32767, v))) + return out + + +def golden(ext): + with open(os.path.join(DATA, "signal_16k." + ext), "rb") as f: + return f.read() + + +def run(pcm, rate, block, channels=1, raw=True): + """Every slice the frontend makes from pcm, fed `block` frames at a time.""" + fe = audiofeatures.Frontend(slices=64) + width = 2 * SIZE if raw else SIZE + buf = bytearray(width * 64) + out = bytearray() + mv = memoryview(pcm).cast("B") + step = block * 2 * channels + for at in range(0, len(mv), step): + fe.feed(mv[at:at + step], rate, channels) + n = fe.readinto_raw(buf) if raw else fe.readinto(buf) + out += buf[:n * width] + st = fe.stats() + fe.deinit() + assert st["dropped"] == 0, st + return bytes(out) + + +def upsample(x, factor_num, factor_den=1): + """Band-limited resampling by FFT: the same audio at another rate.""" + x = np.asarray(x, dtype=np.float64) + n = len(x) + m = n * factor_num // factor_den + spec = np.fft.rfft(x) + out = np.zeros(m // 2 + 1, dtype=complex) + k = min(len(spec), len(out)) + out[:k] = spec[:k] + y = np.fft.irfft(out, m) * (m / n) + return np.clip(np.round(y), -32768, 32767).astype(" is included. + +Nothing in this directory is compiled on its own. src/shared/audiodsp_microfrontend.c +defines FIXED_POINT as 16 and includes kiss_fft.c and tools/kiss_fftr.c, which +makes the int16 fixed-point FFT TensorFlow Lite Micro's kiss_fft_int16.cc +makes, and renames its symbols to audiodsp_mf_*. diff --git a/third_party/kissfft/_kiss_fft_guts.h b/third_party/kissfft/_kiss_fft_guts.h new file mode 100644 index 00000000..1a0f4c26 --- /dev/null +++ b/third_party/kissfft/_kiss_fft_guts.h @@ -0,0 +1,168 @@ +#ifndef _KISS_FFT_GUTS_H +#define _KISS_FFT_GUTS_H + +/* +Copyright (c) 2003-2010, Mark Borgerding + +All rights reserved. + +Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met: + + * Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer. + * Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution. + * Neither the author nor the names of any contributors may be used to endorse or promote products derived from this software without specific prior written permission. + +THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. +*/ + +/* kiss_fft.h + defines kiss_fft_scalar as either short or a float type + and defines + typedef struct { kiss_fft_scalar r; kiss_fft_scalar i; }kiss_fft_cpx; */ +#include "kiss_fft.h" +#include + +#define MAXFACTORS 32 +/* e.g. an fft of length 128 has 4 factors + as far as kissfft is concerned + 4*4*4*2 + */ + +struct kiss_fft_state{ + int nfft; + int inverse; + int factors[2*MAXFACTORS]; + kiss_fft_cpx twiddles[1]; +}; + +/* + Explanation of macros dealing with complex math: + + C_MUL(m,a,b) : m = a*b + C_FIXDIV( c , div ) : if a fixed point impl., c /= div. noop otherwise + C_SUB( res, a,b) : res = a - b + C_SUBFROM( res , a) : res -= a + C_ADDTO( res , a) : res += a + * */ +#ifdef FIXED_POINT +#if (FIXED_POINT==32) +# define FRACBITS 31 +# define SAMPPROD int64_t +#define SAMP_MAX 2147483647 +#else +# define FRACBITS 15 +# define SAMPPROD int32_t +#define SAMP_MAX 32767 +#endif + +#define SAMP_MIN -SAMP_MAX + +#if defined(CHECK_OVERFLOW) +# define CHECK_OVERFLOW_OP(a,op,b) \ + if ( (SAMPPROD)(a) op (SAMPPROD)(b) > SAMP_MAX || (SAMPPROD)(a) op (SAMPPROD)(b) < SAMP_MIN ) { \ + fprintf(stderr,"WARNING:overflow @ " __FILE__ "(%d): (%d " #op" %d) = %ld\n",__LINE__,(a),(b),(SAMPPROD)(a) op (SAMPPROD)(b) ); } +#endif + + +# define smul(a,b) ( (SAMPPROD)(a)*(b) ) +# define sround( x ) (kiss_fft_scalar)( ( (x) + (1<<(FRACBITS-1)) ) >> FRACBITS ) + +# define S_MUL(a,b) sround( smul(a,b) ) + +# define C_MUL(m,a,b) \ + do{ (m).r = sround( smul((a).r,(b).r) - smul((a).i,(b).i) ); \ + (m).i = sround( smul((a).r,(b).i) + smul((a).i,(b).r) ); }while(0) + +# define DIVSCALAR(x,k) \ + (x) = sround( smul( x, SAMP_MAX/k ) ) + +# define C_FIXDIV(c,div) \ + do { DIVSCALAR( (c).r , div); \ + DIVSCALAR( (c).i , div); }while (0) + +# define C_MULBYSCALAR( c, s ) \ + do{ (c).r = sround( smul( (c).r , s ) ) ;\ + (c).i = sround( smul( (c).i , s ) ) ; }while(0) + +#else /* not FIXED_POINT*/ + +# define S_MUL(a,b) ( (a)*(b) ) +#define C_MUL(m,a,b) \ + do{ (m).r = (a).r*(b).r - (a).i*(b).i;\ + (m).i = (a).r*(b).i + (a).i*(b).r; }while(0) +# define C_FIXDIV(c,div) /* NOOP */ +# define C_MULBYSCALAR( c, s ) \ + do{ (c).r *= (s);\ + (c).i *= (s); }while(0) +#endif + +#ifndef CHECK_OVERFLOW_OP +# define CHECK_OVERFLOW_OP(a,op,b) /* noop */ +#endif + +#define C_ADD( res, a,b)\ + do { \ + CHECK_OVERFLOW_OP((a).r,+,(b).r)\ + CHECK_OVERFLOW_OP((a).i,+,(b).i)\ + (res).r=(a).r+(b).r; (res).i=(a).i+(b).i; \ + }while(0) +#define C_SUB( res, a,b)\ + do { \ + CHECK_OVERFLOW_OP((a).r,-,(b).r)\ + CHECK_OVERFLOW_OP((a).i,-,(b).i)\ + (res).r=(a).r-(b).r; (res).i=(a).i-(b).i; \ + }while(0) +#define C_ADDTO( res , a)\ + do { \ + CHECK_OVERFLOW_OP((res).r,+,(a).r)\ + CHECK_OVERFLOW_OP((res).i,+,(a).i)\ + (res).r += (a).r; (res).i += (a).i;\ + }while(0) + +#define C_SUBFROM( res , a)\ + do {\ + CHECK_OVERFLOW_OP((res).r,-,(a).r)\ + CHECK_OVERFLOW_OP((res).i,-,(a).i)\ + (res).r -= (a).r; (res).i -= (a).i; \ + }while(0) + + +#ifdef FIXED_POINT +# define KISS_FFT_COS(phase) floor(.5+SAMP_MAX * cos (phase)) +# define KISS_FFT_SIN(phase) floor(.5+SAMP_MAX * sin (phase)) +# define HALF_OF(x) ((x)>>1) +#elif defined(USE_SIMD) +# define KISS_FFT_COS(phase) _mm_set1_ps( cos(phase) ) +# define KISS_FFT_SIN(phase) _mm_set1_ps( sin(phase) ) +# define HALF_OF(x) ((x)*_mm_set1_ps(.5)) +#else +# define KISS_FFT_COS(phase) (kiss_fft_scalar) cos(phase) +# define KISS_FFT_SIN(phase) (kiss_fft_scalar) sin(phase) +# define HALF_OF(x) ((x)*(kiss_fft_scalar).5) +#endif + +#define kf_cexp(x,phase) \ + do{ \ + (x)->r = KISS_FFT_COS(phase);\ + (x)->i = KISS_FFT_SIN(phase);\ + }while(0) + + +/* a debugging function */ +#define pcpx(c)\ + fprintf(stderr,"%g + %gi\n",(double)((c)->r),(double)((c)->i) ) + + +#ifdef KISS_FFT_USE_ALLOCA +// define this to allow use of alloca instead of malloc for temporary buffers +// Temporary buffers are used in two case: +// 1. FFT sizes that have "bad" factors. i.e. not 2,3 and 5 +// 2. "in-place" FFTs. Notice the quotes, since kissfft does not really do an in-place transform. +#include +#define KISS_FFT_TMP_ALLOC(nbytes) alloca(nbytes) +#define KISS_FFT_TMP_FREE(ptr) +#else +#define KISS_FFT_TMP_ALLOC(nbytes) KISS_FFT_MALLOC(nbytes) +#define KISS_FFT_TMP_FREE(ptr) KISS_FFT_FREE(ptr) +#endif +#endif // _KISS_FFT_GUTS_H diff --git a/third_party/kissfft/kiss_fft.c b/third_party/kissfft/kiss_fft.c new file mode 100644 index 00000000..9133a013 --- /dev/null +++ b/third_party/kissfft/kiss_fft.c @@ -0,0 +1,408 @@ +/* +Copyright (c) 2003-2010, Mark Borgerding + +All rights reserved. + +Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met: + + * Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer. + * Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution. + * Neither the author nor the names of any contributors may be used to endorse or promote products derived from this software without specific prior written permission. + +THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. +*/ + + +#include "_kiss_fft_guts.h" +/* The guts header contains all the multiplication and addition macros that are defined for + fixed or floating point complex numbers. It also delares the kf_ internal functions. + */ + +static void kf_bfly2( + kiss_fft_cpx * Fout, + const size_t fstride, + const kiss_fft_cfg st, + int m + ) +{ + kiss_fft_cpx * Fout2; + kiss_fft_cpx * tw1 = st->twiddles; + kiss_fft_cpx t; + Fout2 = Fout + m; + do{ + C_FIXDIV(*Fout,2); C_FIXDIV(*Fout2,2); + + C_MUL (t, *Fout2 , *tw1); + tw1 += fstride; + C_SUB( *Fout2 , *Fout , t ); + C_ADDTO( *Fout , t ); + ++Fout2; + ++Fout; + }while (--m); +} + +static void kf_bfly4( + kiss_fft_cpx * Fout, + const size_t fstride, + const kiss_fft_cfg st, + const size_t m + ) +{ + kiss_fft_cpx *tw1,*tw2,*tw3; + kiss_fft_cpx scratch[6]; + size_t k=m; + const size_t m2=2*m; + const size_t m3=3*m; + + + tw3 = tw2 = tw1 = st->twiddles; + + do { + C_FIXDIV(*Fout,4); C_FIXDIV(Fout[m],4); C_FIXDIV(Fout[m2],4); C_FIXDIV(Fout[m3],4); + + C_MUL(scratch[0],Fout[m] , *tw1 ); + C_MUL(scratch[1],Fout[m2] , *tw2 ); + C_MUL(scratch[2],Fout[m3] , *tw3 ); + + C_SUB( scratch[5] , *Fout, scratch[1] ); + C_ADDTO(*Fout, scratch[1]); + C_ADD( scratch[3] , scratch[0] , scratch[2] ); + C_SUB( scratch[4] , scratch[0] , scratch[2] ); + C_SUB( Fout[m2], *Fout, scratch[3] ); + tw1 += fstride; + tw2 += fstride*2; + tw3 += fstride*3; + C_ADDTO( *Fout , scratch[3] ); + + if(st->inverse) { + Fout[m].r = scratch[5].r - scratch[4].i; + Fout[m].i = scratch[5].i + scratch[4].r; + Fout[m3].r = scratch[5].r + scratch[4].i; + Fout[m3].i = scratch[5].i - scratch[4].r; + }else{ + Fout[m].r = scratch[5].r + scratch[4].i; + Fout[m].i = scratch[5].i - scratch[4].r; + Fout[m3].r = scratch[5].r - scratch[4].i; + Fout[m3].i = scratch[5].i + scratch[4].r; + } + ++Fout; + }while(--k); +} + +static void kf_bfly3( + kiss_fft_cpx * Fout, + const size_t fstride, + const kiss_fft_cfg st, + size_t m + ) +{ + size_t k=m; + const size_t m2 = 2*m; + kiss_fft_cpx *tw1,*tw2; + kiss_fft_cpx scratch[5]; + kiss_fft_cpx epi3; + epi3 = st->twiddles[fstride*m]; + + tw1=tw2=st->twiddles; + + do{ + C_FIXDIV(*Fout,3); C_FIXDIV(Fout[m],3); C_FIXDIV(Fout[m2],3); + + C_MUL(scratch[1],Fout[m] , *tw1); + C_MUL(scratch[2],Fout[m2] , *tw2); + + C_ADD(scratch[3],scratch[1],scratch[2]); + C_SUB(scratch[0],scratch[1],scratch[2]); + tw1 += fstride; + tw2 += fstride*2; + + Fout[m].r = Fout->r - HALF_OF(scratch[3].r); + Fout[m].i = Fout->i - HALF_OF(scratch[3].i); + + C_MULBYSCALAR( scratch[0] , epi3.i ); + + C_ADDTO(*Fout,scratch[3]); + + Fout[m2].r = Fout[m].r + scratch[0].i; + Fout[m2].i = Fout[m].i - scratch[0].r; + + Fout[m].r -= scratch[0].i; + Fout[m].i += scratch[0].r; + + ++Fout; + }while(--k); +} + +static void kf_bfly5( + kiss_fft_cpx * Fout, + const size_t fstride, + const kiss_fft_cfg st, + int m + ) +{ + kiss_fft_cpx *Fout0,*Fout1,*Fout2,*Fout3,*Fout4; + int u; + kiss_fft_cpx scratch[13]; + kiss_fft_cpx * twiddles = st->twiddles; + kiss_fft_cpx *tw; + kiss_fft_cpx ya,yb; + ya = twiddles[fstride*m]; + yb = twiddles[fstride*2*m]; + + Fout0=Fout; + Fout1=Fout0+m; + Fout2=Fout0+2*m; + Fout3=Fout0+3*m; + Fout4=Fout0+4*m; + + tw=st->twiddles; + for ( u=0; ur += scratch[7].r + scratch[8].r; + Fout0->i += scratch[7].i + scratch[8].i; + + scratch[5].r = scratch[0].r + S_MUL(scratch[7].r,ya.r) + S_MUL(scratch[8].r,yb.r); + scratch[5].i = scratch[0].i + S_MUL(scratch[7].i,ya.r) + S_MUL(scratch[8].i,yb.r); + + scratch[6].r = S_MUL(scratch[10].i,ya.i) + S_MUL(scratch[9].i,yb.i); + scratch[6].i = -S_MUL(scratch[10].r,ya.i) - S_MUL(scratch[9].r,yb.i); + + C_SUB(*Fout1,scratch[5],scratch[6]); + C_ADD(*Fout4,scratch[5],scratch[6]); + + scratch[11].r = scratch[0].r + S_MUL(scratch[7].r,yb.r) + S_MUL(scratch[8].r,ya.r); + scratch[11].i = scratch[0].i + S_MUL(scratch[7].i,yb.r) + S_MUL(scratch[8].i,ya.r); + scratch[12].r = - S_MUL(scratch[10].i,yb.i) + S_MUL(scratch[9].i,ya.i); + scratch[12].i = S_MUL(scratch[10].r,yb.i) - S_MUL(scratch[9].r,ya.i); + + C_ADD(*Fout2,scratch[11],scratch[12]); + C_SUB(*Fout3,scratch[11],scratch[12]); + + ++Fout0;++Fout1;++Fout2;++Fout3;++Fout4; + } +} + +/* perform the butterfly for one stage of a mixed radix FFT */ +static void kf_bfly_generic( + kiss_fft_cpx * Fout, + const size_t fstride, + const kiss_fft_cfg st, + int m, + int p + ) +{ + int u,k,q1,q; + kiss_fft_cpx * twiddles = st->twiddles; + kiss_fft_cpx t; + int Norig = st->nfft; + + kiss_fft_cpx * scratch = (kiss_fft_cpx*)KISS_FFT_TMP_ALLOC(sizeof(kiss_fft_cpx)*p); + + for ( u=0; u=Norig) twidx-=Norig; + C_MUL(t,scratch[q] , twiddles[twidx] ); + C_ADDTO( Fout[ k ] ,t); + } + k += m; + } + } + KISS_FFT_TMP_FREE(scratch); +} + +static +void kf_work( + kiss_fft_cpx * Fout, + const kiss_fft_cpx * f, + const size_t fstride, + int in_stride, + int * factors, + const kiss_fft_cfg st + ) +{ + kiss_fft_cpx * Fout_beg=Fout; + const int p=*factors++; /* the radix */ + const int m=*factors++; /* stage's fft length/p */ + const kiss_fft_cpx * Fout_end = Fout + p*m; + +#ifdef _OPENMP + // use openmp extensions at the + // top-level (not recursive) + if (fstride==1 && p<=5) + { + int k; + + // execute the p different work units in different threads +# pragma omp parallel for + for (k=0;k floor_sqrt) + p = n; /* no more factors, skip to end */ + } + n /= p; + *facbuf++ = p; + *facbuf++ = n; + } while (n > 1); +} + +/* + * + * User-callable function to allocate all necessary storage space for the fft. + * + * The return value is a contiguous block of memory, allocated with malloc. As such, + * It can be freed with free(), rather than a kiss_fft-specific function. + * */ +kiss_fft_cfg kiss_fft_alloc(int nfft,int inverse_fft,void * mem,size_t * lenmem ) +{ + kiss_fft_cfg st=NULL; + size_t memneeded = sizeof(struct kiss_fft_state) + + sizeof(kiss_fft_cpx)*(nfft-1); /* twiddle factors*/ + + if ( lenmem==NULL ) { + st = ( kiss_fft_cfg)KISS_FFT_MALLOC( memneeded ); + }else{ + if (mem != NULL && *lenmem >= memneeded) + st = (kiss_fft_cfg)mem; + *lenmem = memneeded; + } + if (st) { + int i; + st->nfft=nfft; + st->inverse = inverse_fft; + + for (i=0;iinverse) + phase *= -1; + kf_cexp(st->twiddles+i, phase ); + } + + kf_factor(nfft,st->factors); + } + return st; +} + + +void kiss_fft_stride(kiss_fft_cfg st,const kiss_fft_cpx *fin,kiss_fft_cpx *fout,int in_stride) +{ + if (fin == fout) { + //NOTE: this is not really an in-place FFT algorithm. + //It just performs an out-of-place FFT into a temp buffer + kiss_fft_cpx * tmpbuf = (kiss_fft_cpx*)KISS_FFT_TMP_ALLOC( sizeof(kiss_fft_cpx)*st->nfft); + kf_work(tmpbuf,fin,1,in_stride, st->factors,st); + /* memcpy(fout,tmpbuf,sizeof(kiss_fft_cpx)*st->nfft); */ + KISS_FFT_TMP_FREE(tmpbuf); + }else{ + kf_work( fout, fin, 1,in_stride, st->factors,st ); + } +} + +void kiss_fft(kiss_fft_cfg cfg,const kiss_fft_cpx *fin,kiss_fft_cpx *fout) +{ + kiss_fft_stride(cfg,fin,fout,1); +} + + +void kiss_fft_cleanup(void) +{ + // nothing needed any more +} + +int kiss_fft_next_fast_size(int n) +{ + while(1) { + int m=n; + while ( (m%2) == 0 ) m/=2; + while ( (m%3) == 0 ) m/=3; + while ( (m%5) == 0 ) m/=5; + if (m<=1) + break; /* n is completely factorable by twos, threes, and fives */ + n++; + } + return n; +} diff --git a/third_party/kissfft/kiss_fft.h b/third_party/kissfft/kiss_fft.h new file mode 100644 index 00000000..24e4d0c0 --- /dev/null +++ b/third_party/kissfft/kiss_fft.h @@ -0,0 +1,124 @@ +#ifndef KISS_FFT_H +#define KISS_FFT_H + +#include +#include +#include +#include + +#ifdef __cplusplus +extern "C++" { +#endif + +/* + ATTENTION! + If you would like a : + -- a utility that will handle the caching of fft objects + -- real-only (no imaginary time component ) FFT + -- a multi-dimensional FFT + -- a command-line utility to perform ffts + -- a command-line utility to perform fast-convolution filtering + + Then see kfc.h kiss_fftr.h kiss_fftnd.h fftutil.c kiss_fastfir.c + in the tools/ directory. +*/ + +#ifdef USE_SIMD +# include +# define kiss_fft_scalar __m128 +#define KISS_FFT_MALLOC(nbytes) _mm_malloc(nbytes,16) +#define KISS_FFT_FREE _mm_free +#else +#define KISS_FFT_MALLOC(X) (void*)(0x0) /* Patched. */ +#define KISS_FFT_FREE(X) /* Patched. */ +#endif + + +#ifdef FIXED_POINT +#include /* Patched. */ +# if (FIXED_POINT == 32) +# define kiss_fft_scalar int32_t +# else +# define kiss_fft_scalar int16_t +# endif +#else +# ifndef kiss_fft_scalar +/* default is float */ +# define kiss_fft_scalar float +# endif +#endif + +typedef struct { + kiss_fft_scalar r; + kiss_fft_scalar i; +}kiss_fft_cpx; + +typedef struct kiss_fft_state* kiss_fft_cfg; + +/* + * kiss_fft_alloc + * + * Initialize a FFT (or IFFT) algorithm's cfg/state buffer. + * + * typical usage: kiss_fft_cfg mycfg=kiss_fft_alloc(1024,0,NULL,NULL); + * + * The return value from fft_alloc is a cfg buffer used internally + * by the fft routine or NULL. + * + * If lenmem is NULL, then kiss_fft_alloc will allocate a cfg buffer using malloc. + * The returned value should be free()d when done to avoid memory leaks. + * + * The state can be placed in a user supplied buffer 'mem': + * If lenmem is not NULL and mem is not NULL and *lenmem is large enough, + * then the function places the cfg in mem and the size used in *lenmem + * and returns mem. + * + * If lenmem is not NULL and ( mem is NULL or *lenmem is not large enough), + * then the function returns NULL and places the minimum cfg + * buffer size in *lenmem. + * */ + +kiss_fft_cfg kiss_fft_alloc(int nfft,int inverse_fft,void * mem,size_t * lenmem); + +/* + * kiss_fft(cfg,in_out_buf) + * + * Perform an FFT on a complex input buffer. + * for a forward FFT, + * fin should be f[0] , f[1] , ... ,f[nfft-1] + * fout will be F[0] , F[1] , ... ,F[nfft-1] + * Note that each element is complex and can be accessed like + f[k].r and f[k].i + * */ +void kiss_fft(kiss_fft_cfg cfg,const kiss_fft_cpx *fin,kiss_fft_cpx *fout); + +/* + A more generic version of the above function. It reads its input from every Nth sample. + * */ +void kiss_fft_stride(kiss_fft_cfg cfg,const kiss_fft_cpx *fin,kiss_fft_cpx *fout,int fin_stride); + +/* If kiss_fft_alloc allocated a buffer, it is one contiguous + buffer and can be simply free()d when no longer needed*/ +#define kiss_fft_free free + +/* + Cleans up some memory that gets managed internally. Not necessary to call, but it might clean up + your compiler output to call this before you exit. +*/ +void kiss_fft_cleanup(void); + + +/* + * Returns the smallest integer k, such that k>=n and k has only "fast" factors (2,3,5) + */ +int kiss_fft_next_fast_size(int n); + +/* for real ffts, we need an even size */ +#define kiss_fftr_next_fast_size_real(n) \ + (kiss_fft_next_fast_size( ((n)+1)>>1)<<1) + +#ifdef __cplusplus +} +#endif + +#endif diff --git a/third_party/kissfft/tools/kiss_fftr.c b/third_party/kissfft/tools/kiss_fftr.c new file mode 100644 index 00000000..0d22a04d --- /dev/null +++ b/third_party/kissfft/tools/kiss_fftr.c @@ -0,0 +1,159 @@ +/* +Copyright (c) 2003-2004, Mark Borgerding + +All rights reserved. + +Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met: + + * Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer. + * Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution. + * Neither the author nor the names of any contributors may be used to endorse or promote products derived from this software without specific prior written permission. + +THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. +*/ + +#include "kiss_fftr.h" +#include "_kiss_fft_guts.h" + +struct kiss_fftr_state{ + kiss_fft_cfg substate; + kiss_fft_cpx * tmpbuf; + kiss_fft_cpx * super_twiddles; +#ifdef USE_SIMD + void * pad; +#endif +}; + +kiss_fftr_cfg kiss_fftr_alloc(int nfft,int inverse_fft,void * mem,size_t * lenmem) +{ + int i; + kiss_fftr_cfg st = NULL; + size_t subsize, memneeded; + + if (nfft & 1) { + /* fprintf(stderr,"Real FFT optimization must be even.\n"); */ + return NULL; + } + nfft >>= 1; + + kiss_fft_alloc (nfft, inverse_fft, NULL, &subsize); + memneeded = sizeof(struct kiss_fftr_state) + subsize + sizeof(kiss_fft_cpx) * ( nfft * 3 / 2); + + if (lenmem == NULL) { + st = (kiss_fftr_cfg) KISS_FFT_MALLOC (memneeded); + } else { + if (*lenmem >= memneeded) + st = (kiss_fftr_cfg) mem; + *lenmem = memneeded; + } + if (!st) + return NULL; + + st->substate = (kiss_fft_cfg) (st + 1); /*just beyond kiss_fftr_state struct */ + st->tmpbuf = (kiss_fft_cpx *) (((char *) st->substate) + subsize); + st->super_twiddles = st->tmpbuf + nfft; + kiss_fft_alloc(nfft, inverse_fft, st->substate, &subsize); + + for (i = 0; i < nfft/2; ++i) { + double phase = + -3.14159265358979323846264338327 * ((double) (i+1) / nfft + .5); + if (inverse_fft) + phase *= -1; + kf_cexp (st->super_twiddles+i,phase); + } + return st; +} + +void kiss_fftr(kiss_fftr_cfg st,const kiss_fft_scalar *timedata,kiss_fft_cpx *freqdata) +{ + /* input buffer timedata is stored row-wise */ + int k,ncfft; + kiss_fft_cpx fpnk,fpk,f1k,f2k,tw,tdc; + + if ( st->substate->inverse) { + /* fprintf(stderr,"kiss fft usage error: improper alloc\n"); */ + return; /* exit(1); */ + } + + ncfft = st->substate->nfft; + + /*perform the parallel fft of two real signals packed in real,imag*/ + kiss_fft( st->substate , (const kiss_fft_cpx*)timedata, st->tmpbuf ); + /* The real part of the DC element of the frequency spectrum in st->tmpbuf + * contains the sum of the even-numbered elements of the input time sequence + * The imag part is the sum of the odd-numbered elements + * + * The sum of tdc.r and tdc.i is the sum of the input time sequence. + * yielding DC of input time sequence + * The difference of tdc.r - tdc.i is the sum of the input (dot product) [1,-1,1,-1... + * yielding Nyquist bin of input time sequence + */ + + tdc.r = st->tmpbuf[0].r; + tdc.i = st->tmpbuf[0].i; + C_FIXDIV(tdc,2); + CHECK_OVERFLOW_OP(tdc.r ,+, tdc.i); + CHECK_OVERFLOW_OP(tdc.r ,-, tdc.i); + freqdata[0].r = tdc.r + tdc.i; + freqdata[ncfft].r = tdc.r - tdc.i; +#ifdef USE_SIMD + freqdata[ncfft].i = freqdata[0].i = _mm_set1_ps(0); +#else + freqdata[ncfft].i = freqdata[0].i = 0; +#endif + + for ( k=1;k <= ncfft/2 ; ++k ) { + fpk = st->tmpbuf[k]; + fpnk.r = st->tmpbuf[ncfft-k].r; + fpnk.i = - st->tmpbuf[ncfft-k].i; + C_FIXDIV(fpk,2); + C_FIXDIV(fpnk,2); + + C_ADD( f1k, fpk , fpnk ); + C_SUB( f2k, fpk , fpnk ); + C_MUL( tw , f2k , st->super_twiddles[k-1]); + + freqdata[k].r = HALF_OF(f1k.r + tw.r); + freqdata[k].i = HALF_OF(f1k.i + tw.i); + freqdata[ncfft-k].r = HALF_OF(f1k.r - tw.r); + freqdata[ncfft-k].i = HALF_OF(tw.i - f1k.i); + } +} + +void kiss_fftri(kiss_fftr_cfg st,const kiss_fft_cpx *freqdata,kiss_fft_scalar *timedata) +{ + /* input buffer timedata is stored row-wise */ + int k, ncfft; + + if (st->substate->inverse == 0) { + /* fprintf (stderr, "kiss fft usage error: improper alloc\n"); */ + return; /* exit (1); */ + } + + ncfft = st->substate->nfft; + + st->tmpbuf[0].r = freqdata[0].r + freqdata[ncfft].r; + st->tmpbuf[0].i = freqdata[0].r - freqdata[ncfft].r; + C_FIXDIV(st->tmpbuf[0],2); + + for (k = 1; k <= ncfft / 2; ++k) { + kiss_fft_cpx fk, fnkc, fek, fok, tmp; + fk = freqdata[k]; + fnkc.r = freqdata[ncfft - k].r; + fnkc.i = -freqdata[ncfft - k].i; + C_FIXDIV( fk , 2 ); + C_FIXDIV( fnkc , 2 ); + + C_ADD (fek, fk, fnkc); + C_SUB (tmp, fk, fnkc); + C_MUL (fok, tmp, st->super_twiddles[k-1]); + C_ADD (st->tmpbuf[k], fek, fok); + C_SUB (st->tmpbuf[ncfft - k], fek, fok); +#ifdef USE_SIMD + st->tmpbuf[ncfft - k].i *= _mm_set1_ps(-1.0); +#else + st->tmpbuf[ncfft - k].i *= -1; +#endif + } + kiss_fft (st->substate, st->tmpbuf, (kiss_fft_cpx *) timedata); +} diff --git a/third_party/kissfft/tools/kiss_fftr.h b/third_party/kissfft/tools/kiss_fftr.h new file mode 100644 index 00000000..b888a28b --- /dev/null +++ b/third_party/kissfft/tools/kiss_fftr.h @@ -0,0 +1,46 @@ +#ifndef KISS_FTR_H +#define KISS_FTR_H + +#include "kiss_fft.h" +#ifdef __cplusplus +extern "C++" { +#endif + + +/* + + Real optimized version can save about 45% cpu time vs. complex fft of a real seq. + + + + */ + +typedef struct kiss_fftr_state *kiss_fftr_cfg; + + +kiss_fftr_cfg kiss_fftr_alloc(int nfft,int inverse_fft,void * mem, size_t * lenmem); +/* + nfft must be even + + If you don't care to allocate space, use mem = lenmem = NULL +*/ + + +void kiss_fftr(kiss_fftr_cfg cfg,const kiss_fft_scalar *timedata,kiss_fft_cpx *freqdata); +/* + input timedata has nfft scalar points + output freqdata has nfft/2+1 complex points +*/ + +void kiss_fftri(kiss_fftr_cfg cfg,const kiss_fft_cpx *freqdata,kiss_fft_scalar *timedata); +/* + input freqdata has nfft/2+1 complex points + output timedata has nfft scalar points +*/ + +#define kiss_fftr_free free + +#ifdef __cplusplus +} +#endif +#endif diff --git a/third_party/microfrontend/LICENSE b/third_party/microfrontend/LICENSE new file mode 100644 index 00000000..d6456956 --- /dev/null +++ b/third_party/microfrontend/LICENSE @@ -0,0 +1,202 @@ + + Apache License + Version 2.0, January 2004 + http://www.apache.org/licenses/ + + TERMS AND CONDITIONS FOR USE, REPRODUCTION, AND DISTRIBUTION + + 1. Definitions. + + "License" shall mean the terms and conditions for use, reproduction, + and distribution as defined by Sections 1 through 9 of this document. + + "Licensor" shall mean the copyright owner or entity authorized by + the copyright owner that is granting the License. + + "Legal Entity" shall mean the union of the acting entity and all + other entities that control, are controlled by, or are under common + control with that entity. For the purposes of this definition, + "control" means (i) the power, direct or indirect, to cause the + direction or management of such entity, whether by contract or + otherwise, or (ii) ownership of fifty percent (50%) or more of the + outstanding shares, or (iii) beneficial ownership of such entity. + + "You" (or "Your") shall mean an individual or Legal Entity + exercising permissions granted by this License. + + "Source" form shall mean the preferred form for making modifications, + including but not limited to software source code, documentation + source, and configuration files. + + "Object" form shall mean any form resulting from mechanical + transformation or translation of a Source form, including but + not limited to compiled object code, generated documentation, + and conversions to other media types. + + "Work" shall mean the work of authorship, whether in Source or + Object form, made available under the License, as indicated by a + copyright notice that is included in or attached to the work + (an example is provided in the Appendix below). + + "Derivative Works" shall mean any work, whether in Source or Object + form, that is based on (or derived from) the Work and for which the + editorial revisions, annotations, elaborations, or other modifications + represent, as a whole, an original work of authorship. For the purposes + of this License, Derivative Works shall not include works that remain + separable from, or merely link (or bind by name) to the interfaces of, + the Work and Derivative Works thereof. + + "Contribution" shall mean any work of authorship, including + the original version of the Work and any modifications or additions + to that Work or Derivative Works thereof, that is intentionally + submitted to Licensor for inclusion in the Work by the copyright owner + or by an individual or Legal Entity authorized to submit on behalf of + the copyright owner. For the purposes of this definition, "submitted" + means any form of electronic, verbal, or written communication sent + to the Licensor or its representatives, including but not limited to + communication on electronic mailing lists, source code control systems, + and issue tracking systems that are managed by, or on behalf of, the + Licensor for the purpose of discussing and improving the Work, but + excluding communication that is conspicuously marked or otherwise + designated in writing by the copyright owner as "Not a Contribution." + + "Contributor" shall mean Licensor and any individual or Legal Entity + on behalf of whom a Contribution has been received by Licensor and + subsequently incorporated within the Work. + + 2. Grant of Copyright License. Subject to the terms and conditions of + this License, each Contributor hereby grants to You a perpetual, + worldwide, non-exclusive, no-charge, royalty-free, irrevocable + copyright license to reproduce, prepare Derivative Works of, + publicly display, publicly perform, sublicense, and distribute the + Work and such Derivative Works in Source or Object form. + + 3. Grant of Patent License. Subject to the terms and conditions of + this License, each Contributor hereby grants to You a perpetual, + worldwide, non-exclusive, no-charge, royalty-free, irrevocable + (except as stated in this section) patent license to make, have made, + use, offer to sell, sell, import, and otherwise transfer the Work, + where such license applies only to those patent claims licensable + by such Contributor that are necessarily infringed by their + Contribution(s) alone or by combination of their Contribution(s) + with the Work to which such Contribution(s) was submitted. If You + institute patent litigation against any entity (including a + cross-claim or counterclaim in a lawsuit) alleging that the Work + or a Contribution incorporated within the Work constitutes direct + or contributory patent infringement, then any patent licenses + granted to You under this License for that Work shall terminate + as of the date such litigation is filed. + + 4. Redistribution. You may reproduce and distribute copies of the + Work or Derivative Works thereof in any medium, with or without + modifications, and in Source or Object form, provided that You + meet the following conditions: + + (a) You must give any other recipients of the Work or + Derivative Works a copy of this License; and + + (b) You must cause any modified files to carry prominent notices + stating that You changed the files; and + + (c) You must retain, in the Source form of any Derivative Works + that You distribute, all copyright, patent, trademark, and + attribution notices from the Source form of the Work, + excluding those notices that do not pertain to any part of + the Derivative Works; and + + (d) If the Work includes a "NOTICE" text file as part of its + distribution, then any Derivative Works that You distribute must + include a readable copy of the attribution notices contained + within such NOTICE file, excluding those notices that do not + pertain to any part of the Derivative Works, in at least one + of the following places: within a NOTICE text file distributed + as part of the Derivative Works; within the Source form or + documentation, if provided along with the Derivative Works; or, + within a display generated by the Derivative Works, if and + wherever such third-party notices normally appear. The contents + of the NOTICE file are for informational purposes only and + do not modify the License. You may add Your own attribution + notices within Derivative Works that You distribute, alongside + or as an addendum to the NOTICE text from the Work, provided + that such additional attribution notices cannot be construed + as modifying the License. + + You may add Your own copyright statement to Your modifications and + may provide additional or different license terms and conditions + for use, reproduction, or distribution of Your modifications, or + for any such Derivative Works as a whole, provided Your use, + reproduction, and distribution of the Work otherwise complies with + the conditions stated in this License. + + 5. Submission of Contributions. Unless You explicitly state otherwise, + any Contribution intentionally submitted for inclusion in the Work + by You to the Licensor shall be under the terms and conditions of + this License, without any additional terms or conditions. + Notwithstanding the above, nothing herein shall supersede or modify + the terms of any separate license agreement you may have executed + with Licensor regarding such Contributions. + + 6. Trademarks. This License does not grant permission to use the trade + names, trademarks, service marks, or product names of the Licensor, + except as required for reasonable and customary use in describing the + origin of the Work and reproducing the content of the NOTICE file. + + 7. Disclaimer of Warranty. Unless required by applicable law or + agreed to in writing, Licensor provides the Work (and each + Contributor provides its Contributions) on an "AS IS" BASIS, + WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or + implied, including, without limitation, any warranties or conditions + of TITLE, NON-INFRINGEMENT, MERCHANTABILITY, or FITNESS FOR A + PARTICULAR PURPOSE. You are solely responsible for determining the + appropriateness of using or redistributing the Work and assume any + risks associated with Your exercise of permissions under this License. + + 8. Limitation of Liability. In no event and under no legal theory, + whether in tort (including negligence), contract, or otherwise, + unless required by applicable law (such as deliberate and grossly + negligent acts) or agreed to in writing, shall any Contributor be + liable to You for damages, including any direct, indirect, special, + incidental, or consequential damages of any character arising as a + result of this License or out of the use or inability to use the + Work (including but not limited to damages for loss of goodwill, + work stoppage, computer failure or malfunction, or any and all + other commercial damages or losses), even if such Contributor + has been advised of the possibility of such damages. + + 9. Accepting Warranty or Additional Liability. While redistributing + the Work or Derivative Works thereof, You may choose to offer, + and charge a fee for, acceptance of support, warranty, indemnity, + or other liability obligations and/or rights consistent with this + License. However, in accepting such obligations, You may act only + on Your own behalf and on Your sole responsibility, not on behalf + of any other Contributor, and only if You agree to indemnify, + defend, and hold each Contributor harmless for any liability + incurred by, or claims asserted against, such Contributor by reason + of your accepting any such warranty or additional liability. + + END OF TERMS AND CONDITIONS + + APPENDIX: How to apply the Apache License to your work. + + To apply the Apache License to your work, attach the following + boilerplate notice, with the fields enclosed by brackets "[]" + replaced with your own identifying information. (Don't include + the brackets!) The text should be enclosed in the appropriate + comment syntax for the file format. We also recommend that a + file or class name and description of purpose be included on the + same "printed page" as the copyright notice for easier + identification within third-party archives. + + Copyright [yyyy] [name of copyright owner] + + Licensed under the Apache License, Version 2.0 (the "License"); + you may not use this file except in compliance with the License. + You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + + Unless required by applicable law or agreed to in writing, software + distributed under the License is distributed on an "AS IS" BASIS, + WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + See the License for the specific language governing permissions and + limitations under the License. diff --git a/third_party/microfrontend/ORIGIN.txt b/third_party/microfrontend/ORIGIN.txt new file mode 100644 index 00000000..25faa34b --- /dev/null +++ b/third_party/microfrontend/ORIGIN.txt @@ -0,0 +1,20 @@ +TensorFlow Lite Micro's audio microfrontend: the feature extractor microWakeWord's +models were trained on (window, FFT, mel filterbank, noise reduction, PCAN gain +control, log scale). + +Source: https://github.com/tensorflow/tflite-micro + (tensorflow/lite/experimental/microfrontend/lib) +Commit: 230c9d9531bd9b08d445b4f2ed0760756ffd5f78 (2026-10-09) +License: Apache-2.0, see LICENSE + +The files here are that commit's library files, unchanged. The tests, the +_io dumpers and the main programs are left out. + +Nothing in this directory is compiled on its own. src/shared/audiodsp_microfrontend.c +includes the .c files into one translation unit and renames every external +symbol to audiodsp_mf_*, so a firmware that links a second copy of the +frontend (TensorFlow Lite Micro's own, say) still links. The three C++ files +(fft.cc, fft_util.cc and kiss_fft_int16.cc) are not compiled either: audiodsp_microfrontend.c holds their C equivalent, so no +build needs a C++ compiler. Its output is byte-identical to the C++ build's, +which tests/test_cpython_audiofeatures.py checks against goldens that build +wrote. diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/bits.h b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/bits.h new file mode 100644 index 00000000..04b3ba6f --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/bits.h @@ -0,0 +1,102 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#ifndef TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_BITS_H_ +#define TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_BITS_H_ + +#ifdef __cplusplus +#include + +extern "C" { +#endif + +static inline int CountLeadingZeros32Slow(uint64_t n) { + int zeroes = 28; + if (n >> 16) zeroes -= 16, n >>= 16; + if (n >> 8) zeroes -= 8, n >>= 8; + if (n >> 4) zeroes -= 4, n >>= 4; + return "\4\3\2\2\1\1\1\1\0\0\0\0\0\0\0"[n] + zeroes; +} + +static inline int CountLeadingZeros32(uint32_t n) { +#if defined(_MSC_VER) + unsigned long result = 0; // NOLINT(runtime/int) + if (_BitScanReverse(&result, n)) { + return 31 - result; + } + return 32; +#elif defined(__GNUC__) + + // Handle 0 as a special case because __builtin_clz(0) is undefined. + if (n == 0) { + return 32; + } + return __builtin_clz(n); +#else + return CountLeadingZeros32Slow(n); +#endif +} + +static inline int MostSignificantBit32(uint32_t n) { + return 32 - CountLeadingZeros32(n); +} + +static inline int CountLeadingZeros64Slow(uint64_t n) { + int zeroes = 60; + if (n >> 32) zeroes -= 32, n >>= 32; + if (n >> 16) zeroes -= 16, n >>= 16; + if (n >> 8) zeroes -= 8, n >>= 8; + if (n >> 4) zeroes -= 4, n >>= 4; + return "\4\3\2\2\1\1\1\1\0\0\0\0\0\0\0"[n] + zeroes; +} + +static inline int CountLeadingZeros64(uint64_t n) { +#if defined(_MSC_VER) && defined(_M_X64) + // MSVC does not have __builtin_clzll. Use _BitScanReverse64. + unsigned long result = 0; // NOLINT(runtime/int) + if (_BitScanReverse64(&result, n)) { + return 63 - result; + } + return 64; +#elif defined(_MSC_VER) + // MSVC does not have __builtin_clzll. Compose two calls to _BitScanReverse + unsigned long result = 0; // NOLINT(runtime/int) + if ((n >> 32) && _BitScanReverse(&result, n >> 32)) { + return 31 - result; + } + if (_BitScanReverse(&result, n)) { + return 63 - result; + } + return 64; +#elif defined(__GNUC__) + + // Handle 0 as a special case because __builtin_clzll(0) is undefined. + if (n == 0) { + return 64; + } + return __builtin_clzll(n); +#else + return CountLeadingZeros64Slow(n); +#endif +} + +static inline int MostSignificantBit64(uint64_t n) { + return 64 - CountLeadingZeros64(n); +} + +#ifdef __cplusplus +} // extern "C" +#endif + +#endif // TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_BITS_H_ diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/fft.cc b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/fft.cc new file mode 100644 index 00000000..bcdd9cc0 --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/fft.cc @@ -0,0 +1,52 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#include "tensorflow/lite/experimental/microfrontend/lib/fft.h" + +#include + +#include "tensorflow/lite/experimental/microfrontend/lib/kiss_fft_int16.h" + +void FftCompute(struct FftState* state, const int16_t* input, + int input_scale_shift) { + const size_t input_size = state->input_size; + const size_t fft_size = state->fft_size; + + int16_t* fft_input = state->input; + // First, scale the input by the given shift. + size_t i; + for (i = 0; i < input_size; ++i) { + fft_input[i] = static_cast(static_cast(input[i]) + << input_scale_shift); + } + // Zero out whatever else remains in the top part of the input. + for (; i < fft_size; ++i) { + fft_input[i] = 0; + } + + // Apply the FFT. + kissfft_fixed16::kiss_fftr( + reinterpret_cast(state->scratch), + state->input, + reinterpret_cast(state->output)); +} + +void FftInit(struct FftState* state) { + // All the initialization is done in FftPopulateState() +} + +void FftReset(struct FftState* state) { + memset(state->input, 0, state->fft_size * sizeof(*state->input)); + memset(state->output, 0, (state->fft_size / 2 + 1) * sizeof(*state->output)); +} diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/fft.h b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/fft.h new file mode 100644 index 00000000..aaffa69d --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/fft.h @@ -0,0 +1,50 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#ifndef TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_FFT_H_ +#define TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_FFT_H_ + +#include +#include + +#ifdef __cplusplus +extern "C" { +#endif + +struct complex_int16_t { + int16_t real; + int16_t imag; +}; + +struct FftState { + int16_t* input; + struct complex_int16_t* output; + size_t fft_size; + size_t input_size; + void* scratch; + size_t scratch_size; +}; + +void FftCompute(struct FftState* state, const int16_t* input, + int input_scale_shift); + +void FftInit(struct FftState* state); + +void FftReset(struct FftState* state); + +#ifdef __cplusplus +} // extern "C" +#endif + +#endif // TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_FFT_H_ diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/fft_util.cc b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/fft_util.cc new file mode 100644 index 00000000..ed3dc8fb --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/fft_util.cc @@ -0,0 +1,70 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#include "tensorflow/lite/experimental/microfrontend/lib/fft_util.h" + +#include + +#include "tensorflow/lite/experimental/microfrontend/lib/kiss_fft_int16.h" + +int FftPopulateState(struct FftState* state, size_t input_size) { + state->input_size = input_size; + state->fft_size = 1; + while (state->fft_size < state->input_size) { + state->fft_size <<= 1; + } + + state->input = reinterpret_cast( + malloc(state->fft_size * sizeof(*state->input))); + if (state->input == nullptr) { + fprintf(stderr, "Failed to alloc fft input buffer\n"); + return 0; + } + + state->output = reinterpret_cast( + malloc((state->fft_size / 2 + 1) * sizeof(*state->output) * 2)); + if (state->output == nullptr) { + fprintf(stderr, "Failed to alloc fft output buffer\n"); + return 0; + } + + // Ask kissfft how much memory it wants. + size_t scratch_size = 0; + kissfft_fixed16::kiss_fftr_cfg kfft_cfg = kissfft_fixed16::kiss_fftr_alloc( + state->fft_size, 0, nullptr, &scratch_size); + if (kfft_cfg != nullptr) { + fprintf(stderr, "Kiss memory sizing failed.\n"); + return 0; + } + state->scratch = malloc(scratch_size); + if (state->scratch == nullptr) { + fprintf(stderr, "Failed to alloc fft scratch buffer\n"); + return 0; + } + state->scratch_size = scratch_size; + // Let kissfft configure the scratch space we just allocated + kfft_cfg = kissfft_fixed16::kiss_fftr_alloc(state->fft_size, 0, + state->scratch, &scratch_size); + if (kfft_cfg != state->scratch) { + fprintf(stderr, "Kiss memory preallocation strategy failed.\n"); + return 0; + } + return 1; +} + +void FftFreeStateContents(struct FftState* state) { + free(state->input); + free(state->output); + free(state->scratch); +} diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/fft_util.h b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/fft_util.h new file mode 100644 index 00000000..6a471301 --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/fft_util.h @@ -0,0 +1,34 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#ifndef TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_FFT_UTIL_H_ +#define TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_FFT_UTIL_H_ + +#include "tensorflow/lite/experimental/microfrontend/lib/fft.h" + +#ifdef __cplusplus +extern "C" { +#endif + +// Prepares and FFT for the given input size. +int FftPopulateState(struct FftState* state, size_t input_size); + +// Frees any allocated buffers. +void FftFreeStateContents(struct FftState* state); + +#ifdef __cplusplus +} // extern "C" +#endif + +#endif // TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_FFT_UTIL_H_ diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/filterbank.c b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/filterbank.c new file mode 100644 index 00000000..80f8738f --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/filterbank.c @@ -0,0 +1,134 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#include "tensorflow/lite/experimental/microfrontend/lib/filterbank.h" + +#include + +#include "tensorflow/lite/experimental/microfrontend/lib/bits.h" + +void FilterbankConvertFftComplexToEnergy(struct FilterbankState* state, + struct complex_int16_t* fft_output, + int32_t* energy) { + const int end_index = state->end_index; + int i; + energy += state->start_index; + fft_output += state->start_index; + for (i = state->start_index; i < end_index; ++i) { + const int32_t real = fft_output->real; + const int32_t imag = fft_output->imag; + fft_output++; + const uint32_t mag_squared = (real * real) + (imag * imag); + *energy++ = mag_squared; + } +} + +void FilterbankAccumulateChannels(struct FilterbankState* state, + const int32_t* energy) { + uint64_t* work = state->work; + uint64_t weight_accumulator = 0; + uint64_t unweight_accumulator = 0; + + const int16_t* channel_frequency_starts = state->channel_frequency_starts; + const int16_t* channel_weight_starts = state->channel_weight_starts; + const int16_t* channel_widths = state->channel_widths; + + int num_channels_plus_1 = state->num_channels + 1; + int i; + for (i = 0; i < num_channels_plus_1; ++i) { + const int32_t* magnitudes = energy + *channel_frequency_starts++; + const int16_t* weights = state->weights + *channel_weight_starts; + const int16_t* unweights = state->unweights + *channel_weight_starts++; + const int width = *channel_widths++; + int j; + for (j = 0; j < width; ++j) { + weight_accumulator += *weights++ * ((uint64_t)*magnitudes); + unweight_accumulator += *unweights++ * ((uint64_t)*magnitudes); + ++magnitudes; + } + *work++ = weight_accumulator; + weight_accumulator = unweight_accumulator; + unweight_accumulator = 0; + } +} + +static uint16_t Sqrt32(uint32_t num) { + if (num == 0) { + return 0; + } + uint32_t res = 0; + int max_bit_number = 32 - MostSignificantBit32(num); + max_bit_number |= 1; + uint32_t bit = 1U << (31 - max_bit_number); + int iterations = (31 - max_bit_number) / 2 + 1; + while (iterations--) { + if (num >= res + bit) { + num -= res + bit; + res = (res >> 1U) + bit; + } else { + res >>= 1U; + } + bit >>= 2U; + } + // Do rounding - if we have the bits. + if (num > res && res != 0xFFFF) { + ++res; + } + return res; +} + +static uint32_t Sqrt64(uint64_t num) { + // Take a shortcut and just use 32 bit operations if the upper word is all + // clear. This will cause a slight off by one issue for numbers close to 2^32, + // but it probably isn't going to matter (and gives us a big performance win). + if ((num >> 32) == 0) { + return Sqrt32((uint32_t)num); + } + uint64_t res = 0; + int max_bit_number = 64 - MostSignificantBit64(num); + max_bit_number |= 1; + uint64_t bit = 1ULL << (63 - max_bit_number); + int iterations = (63 - max_bit_number) / 2 + 1; + while (iterations--) { + if (num >= res + bit) { + num -= res + bit; + res = (res >> 1U) + bit; + } else { + res >>= 1U; + } + bit >>= 2U; + } + // Do rounding - if we have the bits. + if (num > res && res != 0xFFFFFFFFLL) { + ++res; + } + return res; +} + +uint32_t* FilterbankSqrt(struct FilterbankState* state, int scale_down_shift) { + const int num_channels = state->num_channels; + const uint64_t* work = state->work + 1; + // Reuse the work buffer since we're fine clobbering it at this point to hold + // the output. + uint32_t* output = (uint32_t*)state->work; + int i; + for (i = 0; i < num_channels; ++i) { + *output++ = Sqrt64(*work++) >> scale_down_shift; + } + return (uint32_t*)state->work; +} + +void FilterbankReset(struct FilterbankState* state) { + memset(state->work, 0, (state->num_channels + 1) * sizeof(*state->work)); +} diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/filterbank.h b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/filterbank.h new file mode 100644 index 00000000..1e6d3885 --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/filterbank.h @@ -0,0 +1,63 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#ifndef TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_FILTERBANK_H_ +#define TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_FILTERBANK_H_ + +#include +#include + +#include "tensorflow/lite/experimental/microfrontend/lib/fft.h" + +#define kFilterbankBits 12 + +#ifdef __cplusplus +extern "C" { +#endif + +struct FilterbankState { + int num_channels; + int start_index; + int end_index; + int16_t* channel_frequency_starts; + int16_t* channel_weight_starts; + int16_t* channel_widths; + int16_t* weights; + int16_t* unweights; + uint64_t* work; +}; + +// Converts the relevant complex values of an FFT output into energy (the +// square magnitude). +void FilterbankConvertFftComplexToEnergy(struct FilterbankState* state, + struct complex_int16_t* fft_output, + int32_t* energy); + +// Computes the mel-scale filterbank on the given energy array. Output is cached +// internally - to fetch it, you need to call FilterbankSqrt. +void FilterbankAccumulateChannels(struct FilterbankState* state, + const int32_t* energy); + +// Applies an integer square root to the 64 bit intermediate values of the +// filterbank, and returns a pointer to them. Memory will be invalidated the +// next time FilterbankAccumulateChannels is called. +uint32_t* FilterbankSqrt(struct FilterbankState* state, int scale_down_shift); + +void FilterbankReset(struct FilterbankState* state); + +#ifdef __cplusplus +} // extern "C" +#endif + +#endif // TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_FILTERBANK_H_ diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/filterbank_util.c b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/filterbank_util.c new file mode 100644 index 00000000..430fc6e3 --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/filterbank_util.c @@ -0,0 +1,220 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#include "tensorflow/lite/experimental/microfrontend/lib/filterbank_util.h" + +#include +#include +#include + +#define kFilterbankIndexAlignment 4 +#define kFilterbankChannelBlockSize 4 + +void FilterbankFillConfigWithDefaults(struct FilterbankConfig* config) { + config->num_channels = 32; + config->lower_band_limit = 125.0f; + config->upper_band_limit = 7500.0f; + config->output_scale_shift = 7; +} + +static float FreqToMel(float freq) { return 1127.0f * log1pf(freq / 700.0f); } + +static void CalculateCenterFrequencies(const int num_channels, + const float lower_frequency_limit, + const float upper_frequency_limit, + float* center_frequencies) { + assert(lower_frequency_limit >= 0.0f); + assert(upper_frequency_limit > lower_frequency_limit); + + const float mel_low = FreqToMel(lower_frequency_limit); + const float mel_hi = FreqToMel(upper_frequency_limit); + const float mel_span = mel_hi - mel_low; + const float mel_spacing = mel_span / ((float)num_channels); + int i; + for (i = 0; i < num_channels; ++i) { + center_frequencies[i] = mel_low + (mel_spacing * (i + 1)); + } +} + +static void QuantizeFilterbankWeights(const float float_weight, int16_t* weight, + int16_t* unweight) { + *weight = floorf(float_weight * (1 << kFilterbankBits) + 0.5f); + *unweight = floorf((1.0f - float_weight) * (1 << kFilterbankBits) + 0.5f); +} + +int FilterbankPopulateState(const struct FilterbankConfig* config, + struct FilterbankState* state, int sample_rate, + int spectrum_size) { + state->num_channels = config->num_channels; + const int num_channels_plus_1 = config->num_channels + 1; + + // How should we align things to index counts given the byte alignment? + const int index_alignment = + (kFilterbankIndexAlignment < sizeof(int16_t) + ? 1 + : kFilterbankIndexAlignment / sizeof(int16_t)); + + state->channel_frequency_starts = + malloc(num_channels_plus_1 * sizeof(*state->channel_frequency_starts)); + state->channel_weight_starts = + malloc(num_channels_plus_1 * sizeof(*state->channel_weight_starts)); + state->channel_widths = + malloc(num_channels_plus_1 * sizeof(*state->channel_widths)); + state->work = malloc(num_channels_plus_1 * sizeof(*state->work)); + + float* center_mel_freqs = + malloc(num_channels_plus_1 * sizeof(*center_mel_freqs)); + int16_t* actual_channel_starts = + malloc(num_channels_plus_1 * sizeof(*actual_channel_starts)); + int16_t* actual_channel_widths = + malloc(num_channels_plus_1 * sizeof(*actual_channel_widths)); + + if (state->channel_frequency_starts == NULL || + state->channel_weight_starts == NULL || state->channel_widths == NULL || + center_mel_freqs == NULL || actual_channel_starts == NULL || + actual_channel_widths == NULL) { + free(center_mel_freqs); + free(actual_channel_starts); + free(actual_channel_widths); + fprintf(stderr, "Failed to allocate channel buffers\n"); + return 0; + } + + CalculateCenterFrequencies(num_channels_plus_1, config->lower_band_limit, + config->upper_band_limit, center_mel_freqs); + + // Always exclude DC. + const float hz_per_sbin = 0.5f * sample_rate / ((float)spectrum_size - 1); + state->start_index = 1.5f + config->lower_band_limit / hz_per_sbin; + state->end_index = 0; // Initialized to zero here, but actually set below. + + // For each channel, we need to figure out what frequencies belong to it, and + // how much padding we need to add so that we can efficiently multiply the + // weights and unweights for accumulation. To simplify the multiplication + // logic, all channels will have some multiplication to do (even if there are + // no frequencies that accumulate to that channel) - they will be directed to + // a set of zero weights. + int chan_freq_index_start = state->start_index; + int weight_index_start = 0; + int needs_zeros = 0; + + int chan; + for (chan = 0; chan < num_channels_plus_1; ++chan) { + // Keep jumping frequencies until we overshoot the bound on this channel. + int freq_index = chan_freq_index_start; + while (FreqToMel((freq_index)*hz_per_sbin) <= center_mel_freqs[chan]) { + ++freq_index; + } + + const int width = freq_index - chan_freq_index_start; + actual_channel_starts[chan] = chan_freq_index_start; + actual_channel_widths[chan] = width; + + if (width == 0) { + // This channel doesn't actually get anything from the frequencies, it's + // always zero. We need then to insert some 'zero' weights into the + // output, and just redirect this channel to do a single multiplication at + // this point. For simplicity, the zeros are placed at the beginning of + // the weights arrays, so we have to go and update all the other + // weight_starts to reflect this shift (but only once). + state->channel_frequency_starts[chan] = 0; + state->channel_weight_starts[chan] = 0; + state->channel_widths[chan] = kFilterbankChannelBlockSize; + if (!needs_zeros) { + needs_zeros = 1; + int j; + for (j = 0; j < chan; ++j) { + state->channel_weight_starts[j] += kFilterbankChannelBlockSize; + } + weight_index_start += kFilterbankChannelBlockSize; + } + } else { + // How far back do we need to go to ensure that we have the proper + // alignment? + const int aligned_start = + (chan_freq_index_start / index_alignment) * index_alignment; + const int aligned_width = (chan_freq_index_start - aligned_start + width); + const int padded_width = + (((aligned_width - 1) / kFilterbankChannelBlockSize) + 1) * + kFilterbankChannelBlockSize; + + state->channel_frequency_starts[chan] = aligned_start; + state->channel_weight_starts[chan] = weight_index_start; + state->channel_widths[chan] = padded_width; + weight_index_start += padded_width; + } + chan_freq_index_start = freq_index; + } + + // Allocate the two arrays to store the weights - weight_index_start contains + // the index of what would be the next set of weights that we would need to + // add, so that's how many weights we need to allocate. + state->weights = calloc(weight_index_start, sizeof(*state->weights)); + state->unweights = calloc(weight_index_start, sizeof(*state->unweights)); + + // If the alloc failed, we also need to nuke the arrays. + if (state->weights == NULL || state->unweights == NULL) { + free(center_mel_freqs); + free(actual_channel_starts); + free(actual_channel_widths); + fprintf(stderr, "Failed to allocate weights or unweights\n"); + return 0; + } + + // Next pass, compute all the weights. Since everything has been memset to + // zero, we only need to fill in the weights that correspond to some frequency + // for a channel. + const float mel_low = FreqToMel(config->lower_band_limit); + for (chan = 0; chan < num_channels_plus_1; ++chan) { + int frequency = actual_channel_starts[chan]; + const int num_frequencies = actual_channel_widths[chan]; + const int frequency_offset = + frequency - state->channel_frequency_starts[chan]; + const int weight_start = state->channel_weight_starts[chan]; + const float denom_val = (chan == 0) ? mel_low : center_mel_freqs[chan - 1]; + + int j; + for (j = 0; j < num_frequencies; ++j, ++frequency) { + const float weight = + (center_mel_freqs[chan] - FreqToMel(frequency * hz_per_sbin)) / + (center_mel_freqs[chan] - denom_val); + + // Make the float into an integer for the weights (and unweights). + const int weight_index = weight_start + frequency_offset + j; + QuantizeFilterbankWeights(weight, state->weights + weight_index, + state->unweights + weight_index); + } + if (frequency > state->end_index) { + state->end_index = frequency; + } + } + + free(center_mel_freqs); + free(actual_channel_starts); + free(actual_channel_widths); + if (state->end_index >= spectrum_size) { + fprintf(stderr, "Filterbank end_index is above spectrum size.\n"); + return 0; + } + return 1; +} + +void FilterbankFreeStateContents(struct FilterbankState* state) { + free(state->channel_frequency_starts); + free(state->channel_weight_starts); + free(state->channel_widths); + free(state->weights); + free(state->unweights); + free(state->work); +} diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/filterbank_util.h b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/filterbank_util.h new file mode 100644 index 00000000..781d1024 --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/filterbank_util.h @@ -0,0 +1,50 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#ifndef TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_FILTERBANK_UTIL_H_ +#define TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_FILTERBANK_UTIL_H_ + +#include "tensorflow/lite/experimental/microfrontend/lib/filterbank.h" + +#ifdef __cplusplus +extern "C" { +#endif + +struct FilterbankConfig { + // number of frequency channel buckets for filterbank + int num_channels; + // maximum frequency to include + float upper_band_limit; + // minimum frequency to include + float lower_band_limit; + // unused + int output_scale_shift; +}; + +// Fills the frontendConfig with "sane" defaults. +void FilterbankFillConfigWithDefaults(struct FilterbankConfig* config); + +// Allocates any buffers. +int FilterbankPopulateState(const struct FilterbankConfig* config, + struct FilterbankState* state, int sample_rate, + int spectrum_size); + +// Frees any allocated buffers. +void FilterbankFreeStateContents(struct FilterbankState* state); + +#ifdef __cplusplus +} // extern "C" +#endif + +#endif // TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_FILTERBANK_UTIL_H_ diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/frontend.c b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/frontend.c new file mode 100644 index 00000000..9de2a879 --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/frontend.c @@ -0,0 +1,72 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#include "tensorflow/lite/experimental/microfrontend/lib/frontend.h" + +#include "tensorflow/lite/experimental/microfrontend/lib/bits.h" + +struct FrontendOutput FrontendProcessSamples(struct FrontendState* state, + const int16_t* samples, + size_t num_samples, + size_t* num_samples_read) { + struct FrontendOutput output; + output.values = NULL; + output.size = 0; + + // Try to apply the window - if it fails, return and wait for more data. + if (!WindowProcessSamples(&state->window, samples, num_samples, + num_samples_read)) { + return output; + } + + // Apply the FFT to the window's output (and scale it so that the fixed point + // FFT can have as much resolution as possible). + int input_shift = + 15 - MostSignificantBit32(state->window.max_abs_output_value); + FftCompute(&state->fft, state->window.output, input_shift); + + // We can re-ruse the fft's output buffer to hold the energy. + int32_t* energy = (int32_t*)state->fft.output; + + FilterbankConvertFftComplexToEnergy(&state->filterbank, state->fft.output, + energy); + + FilterbankAccumulateChannels(&state->filterbank, energy); + uint32_t* scaled_filterbank = FilterbankSqrt(&state->filterbank, input_shift); + + // Apply noise reduction. + NoiseReductionApply(&state->noise_reduction, scaled_filterbank); + + if (state->pcan_gain_control.enable_pcan) { + PcanGainControlApply(&state->pcan_gain_control, scaled_filterbank); + } + + // Apply the log and scale. + int correction_bits = + MostSignificantBit32(state->fft.fft_size) - 1 - (kFilterbankBits / 2); + uint16_t* logged_filterbank = + LogScaleApply(&state->log_scale, scaled_filterbank, + state->filterbank.num_channels, correction_bits); + + output.size = state->filterbank.num_channels; + output.values = logged_filterbank; + return output; +} + +void FrontendReset(struct FrontendState* state) { + WindowReset(&state->window); + FftReset(&state->fft); + FilterbankReset(&state->filterbank); + NoiseReductionReset(&state->noise_reduction); +} diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/frontend.h b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/frontend.h new file mode 100644 index 00000000..883df5fd --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/frontend.h @@ -0,0 +1,64 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#ifndef TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_FRONTEND_H_ +#define TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_FRONTEND_H_ + +#include +#include + +#include "tensorflow/lite/experimental/microfrontend/lib/fft.h" +#include "tensorflow/lite/experimental/microfrontend/lib/filterbank.h" +#include "tensorflow/lite/experimental/microfrontend/lib/log_scale.h" +#include "tensorflow/lite/experimental/microfrontend/lib/noise_reduction.h" +#include "tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control.h" +#include "tensorflow/lite/experimental/microfrontend/lib/window.h" + +#ifdef __cplusplus +extern "C" { +#endif + +struct FrontendState { + struct WindowState window; + struct FftState fft; + struct FilterbankState filterbank; + struct NoiseReductionState noise_reduction; + struct PcanGainControlState pcan_gain_control; + struct LogScaleState log_scale; +}; + +struct FrontendOutput { + const uint16_t* values; + size_t size; +}; + +// Main entry point to processing frontend samples. Updates num_samples_read to +// contain the number of samples that have been consumed from the input array. +// Returns a struct containing the generated output. If not enough samples were +// added to generate a feature vector, the returned size will be 0 and the +// values pointer will be NULL. Note that the output pointer will be invalidated +// as soon as FrontendProcessSamples is called again, so copy the contents +// elsewhere if you need to use them later. +struct FrontendOutput FrontendProcessSamples(struct FrontendState* state, + const int16_t* samples, + size_t num_samples, + size_t* num_samples_read); + +void FrontendReset(struct FrontendState* state); + +#ifdef __cplusplus +} // extern "C" +#endif + +#endif // TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_FRONTEND_H_ diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/frontend_util.c b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/frontend_util.c new file mode 100644 index 00000000..27224f6d --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/frontend_util.c @@ -0,0 +1,85 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#include "tensorflow/lite/experimental/microfrontend/lib/frontend_util.h" + +#include +#include + +#include "tensorflow/lite/experimental/microfrontend/lib/bits.h" + +void FrontendFillConfigWithDefaults(struct FrontendConfig* config) { + WindowFillConfigWithDefaults(&config->window); + FilterbankFillConfigWithDefaults(&config->filterbank); + NoiseReductionFillConfigWithDefaults(&config->noise_reduction); + PcanGainControlFillConfigWithDefaults(&config->pcan_gain_control); + LogScaleFillConfigWithDefaults(&config->log_scale); +} + +int FrontendPopulateState(const struct FrontendConfig* config, + struct FrontendState* state, int sample_rate) { + memset(state, 0, sizeof(*state)); + + if (!WindowPopulateState(&config->window, &state->window, sample_rate)) { + fprintf(stderr, "Failed to populate window state\n"); + return 0; + } + + if (!FftPopulateState(&state->fft, state->window.size)) { + fprintf(stderr, "Failed to populate fft state\n"); + return 0; + } + FftInit(&state->fft); + + if (!FilterbankPopulateState(&config->filterbank, &state->filterbank, + sample_rate, state->fft.fft_size / 2 + 1)) { + fprintf(stderr, "Failed to populate filterbank state\n"); + return 0; + } + + if (!NoiseReductionPopulateState(&config->noise_reduction, + &state->noise_reduction, + state->filterbank.num_channels)) { + fprintf(stderr, "Failed to populate noise reduction state\n"); + return 0; + } + + int input_correction_bits = + MostSignificantBit32(state->fft.fft_size) - 1 - (kFilterbankBits / 2); + if (!PcanGainControlPopulateState( + &config->pcan_gain_control, &state->pcan_gain_control, + state->noise_reduction.estimate, state->filterbank.num_channels, + state->noise_reduction.smoothing_bits, input_correction_bits)) { + fprintf(stderr, "Failed to populate pcan gain control state\n"); + return 0; + } + + if (!LogScalePopulateState(&config->log_scale, &state->log_scale)) { + fprintf(stderr, "Failed to populate log scale state\n"); + return 0; + } + + FrontendReset(state); + + // All good, return a true value. + return 1; +} + +void FrontendFreeStateContents(struct FrontendState* state) { + WindowFreeStateContents(&state->window); + FftFreeStateContents(&state->fft); + FilterbankFreeStateContents(&state->filterbank); + NoiseReductionFreeStateContents(&state->noise_reduction); + PcanGainControlFreeStateContents(&state->pcan_gain_control); +} diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/frontend_util.h b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/frontend_util.h new file mode 100644 index 00000000..895ce6cd --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/frontend_util.h @@ -0,0 +1,52 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#ifndef TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_FRONTEND_UTIL_H_ +#define TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_FRONTEND_UTIL_H_ + +#include "tensorflow/lite/experimental/microfrontend/lib/fft_util.h" +#include "tensorflow/lite/experimental/microfrontend/lib/filterbank_util.h" +#include "tensorflow/lite/experimental/microfrontend/lib/frontend.h" +#include "tensorflow/lite/experimental/microfrontend/lib/log_scale_util.h" +#include "tensorflow/lite/experimental/microfrontend/lib/noise_reduction_util.h" +#include "tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control_util.h" +#include "tensorflow/lite/experimental/microfrontend/lib/window_util.h" + +#ifdef __cplusplus +extern "C" { +#endif + +struct FrontendConfig { + struct WindowConfig window; + struct FilterbankConfig filterbank; + struct NoiseReductionConfig noise_reduction; + struct PcanGainControlConfig pcan_gain_control; + struct LogScaleConfig log_scale; +}; + +// Fills the frontendConfig with "sane" defaults. +void FrontendFillConfigWithDefaults(struct FrontendConfig* config); + +// Allocates any buffers. +int FrontendPopulateState(const struct FrontendConfig* config, + struct FrontendState* state, int sample_rate); + +// Frees any allocated buffers. +void FrontendFreeStateContents(struct FrontendState* state); + +#ifdef __cplusplus +} // extern "C" +#endif + +#endif // TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_FRONTEND_UTIL_H_ diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/kiss_fft_common.h b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/kiss_fft_common.h new file mode 100644 index 00000000..33556dab --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/kiss_fft_common.h @@ -0,0 +1,48 @@ +/* Copyright 2019 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ + +#ifndef TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_KISS_FFT_COMMON_H_ +#define TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_KISS_FFT_COMMON_H_ + +// This header file should be included in all variants of kiss_fft_$type.{h,cc} +// so that their sub-included source files do not mistakenly wrap libc header +// files within their kissfft_$type namespaces. +// E.g, This header avoids kissfft_int16.h containing: +// namespace kiss_fft_int16 { +// #include "kiss_fft.h" +// } +// where kiss_fft_.h contains: +// #include +// +// TRICK: By including the following header files here, their preprocessor +// header guards prevent them being re-defined inside of the kiss_fft_$type +// namespaces declared within the kiss_fft_$type.{h,cc} sources. +// Note that the original kiss_fft*.h files are untouched since they +// may be used in libraries that include them directly. + +#include +#include +#include +#include +#include + +#ifdef FIXED_POINT +#include +#endif + +#ifdef USE_SIMD +#include +#endif +#endif // TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_KISS_FFT_COMMON_H_ diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/kiss_fft_int16.cc b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/kiss_fft_int16.cc new file mode 100644 index 00000000..f1e781b2 --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/kiss_fft_int16.cc @@ -0,0 +1,10 @@ +#include + +#include "tensorflow/lite/experimental/microfrontend/lib/kiss_fft_common.h" + +#define FIXED_POINT 16 +namespace kissfft_fixed16 { +#include "kiss_fft.c" +#include "tools/kiss_fftr.c" +} // namespace kissfft_fixed16 +#undef FIXED_POINT diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/kiss_fft_int16.h b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/kiss_fft_int16.h new file mode 100644 index 00000000..fbb3f3a0 --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/kiss_fft_int16.h @@ -0,0 +1,33 @@ +/* Copyright 2019 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ + +#ifndef TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_KISS_FFT_INT16_H_ +#define TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_KISS_FFT_INT16_H_ + +#include "tensorflow/lite/experimental/microfrontend/lib/kiss_fft_common.h" + +// Wrap 16-bit kiss fft in its own namespace. Enables us to link an application +// with different kiss fft resolutions (16/32 bit integer, float, double) +// without getting a linker error. +#define FIXED_POINT 16 +namespace kissfft_fixed16 { +#include "kiss_fft.h" +#include "tools/kiss_fftr.h" +} // namespace kissfft_fixed16 +#undef FIXED_POINT +#undef kiss_fft_scalar +#undef KISS_FFT_H + +#endif // TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_KISS_FFT_INT16_H_ diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_lut.c b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_lut.c new file mode 100644 index 00000000..f59618e0 --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_lut.c @@ -0,0 +1,30 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#include "tensorflow/lite/experimental/microfrontend/lib/log_lut.h" +const uint16_t kLogLut[] +#ifndef _MSC_VER + __attribute__((aligned(4))) +#endif // _MSV_VER + = {0, 224, 442, 654, 861, 1063, 1259, 1450, 1636, 1817, 1992, 2163, + 2329, 2490, 2646, 2797, 2944, 3087, 3224, 3358, 3487, 3611, 3732, 3848, + 3960, 4068, 4172, 4272, 4368, 4460, 4549, 4633, 4714, 4791, 4864, 4934, + 5001, 5063, 5123, 5178, 5231, 5280, 5326, 5368, 5408, 5444, 5477, 5507, + 5533, 5557, 5578, 5595, 5610, 5622, 5631, 5637, 5640, 5641, 5638, 5633, + 5626, 5615, 5602, 5586, 5568, 5547, 5524, 5498, 5470, 5439, 5406, 5370, + 5332, 5291, 5249, 5203, 5156, 5106, 5054, 5000, 4944, 4885, 4825, 4762, + 4697, 4630, 4561, 4490, 4416, 4341, 4264, 4184, 4103, 4020, 3935, 3848, + 3759, 3668, 3575, 3481, 3384, 3286, 3186, 3084, 2981, 2875, 2768, 2659, + 2549, 2437, 2323, 2207, 2090, 1971, 1851, 1729, 1605, 1480, 1353, 1224, + 1094, 963, 830, 695, 559, 421, 282, 142, 0, 0}; diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_lut.h b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_lut.h new file mode 100644 index 00000000..b2448a32 --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_lut.h @@ -0,0 +1,40 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#ifndef TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_LOG_LUT_H_ +#define TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_LOG_LUT_H_ + +#include + +#ifdef __cplusplus +extern "C" { +#endif + +// Number of segments in the log lookup table. The table will be kLogSegments+1 +// in length (with some padding). +#define kLogSegments 128 +#define kLogSegmentsLog2 7 + +// Scale used by lookup table. +#define kLogScale 65536 +#define kLogScaleLog2 16 +#define kLogCoeff 45426 + +extern const uint16_t kLogLut[]; + +#ifdef __cplusplus +} // extern "C" +#endif + +#endif // TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_LOG_LUT_H_ diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_scale.c b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_scale.c new file mode 100644 index 00000000..c27a50a6 --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_scale.c @@ -0,0 +1,83 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#include "tensorflow/lite/experimental/microfrontend/lib/log_scale.h" + +#include "tensorflow/lite/experimental/microfrontend/lib/bits.h" +#include "tensorflow/lite/experimental/microfrontend/lib/log_lut.h" + +#define kuint16max 0x0000FFFF + +// The following functions implement integer logarithms of various sizes. The +// approximation is calculated according to method described in +// www.inti.gob.ar/electronicaeinformatica/instrumentacion/utic/ +// publicaciones/SPL2007/Log10-spl07.pdf +// It first calculates log2 of the input and then converts it to natural +// logarithm. + +static uint32_t Log2FractionPart(const uint32_t x, const uint32_t log2x) { + // Part 1 + int32_t frac = x - (1LL << log2x); + if (log2x < kLogScaleLog2) { + frac <<= kLogScaleLog2 - log2x; + } else { + frac >>= log2x - kLogScaleLog2; + } + // Part 2 + const uint32_t base_seg = frac >> (kLogScaleLog2 - kLogSegmentsLog2); + const uint32_t seg_unit = + (((uint32_t)1) << kLogScaleLog2) >> kLogSegmentsLog2; + + const int32_t c0 = kLogLut[base_seg]; + const int32_t c1 = kLogLut[base_seg + 1]; + const int32_t seg_base = seg_unit * base_seg; + const int32_t rel_pos = ((c1 - c0) * (frac - seg_base)) >> kLogScaleLog2; + return frac + c0 + rel_pos; +} + +static uint32_t Log(const uint32_t x, const uint32_t scale_shift) { + const uint32_t integer = MostSignificantBit32(x) - 1; + const uint32_t fraction = Log2FractionPart(x, integer); + const uint32_t log2 = (integer << kLogScaleLog2) + fraction; + const uint32_t round = kLogScale / 2; + const uint32_t loge = (((uint64_t)kLogCoeff) * log2 + round) >> kLogScaleLog2; + // Finally scale to our output scale + const uint32_t loge_scaled = ((loge << scale_shift) + round) >> kLogScaleLog2; + return loge_scaled; +} + +uint16_t* LogScaleApply(struct LogScaleState* state, uint32_t* signal, + int signal_size, int correction_bits) { + const int scale_shift = state->scale_shift; + uint16_t* output = (uint16_t*)signal; + uint16_t* ret = output; + int i; + for (i = 0; i < signal_size; ++i) { + uint32_t value = *signal++; + if (state->enable_log) { + if (correction_bits < 0) { + value >>= -correction_bits; + } else { + value <<= correction_bits; + } + if (value > 1) { + value = Log(value, scale_shift); + } else { + value = 0; + } + } + *output++ = (value < kuint16max) ? value : kuint16max; + } + return ret; +} diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_scale.h b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_scale.h new file mode 100644 index 00000000..a383f32f --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_scale.h @@ -0,0 +1,39 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#ifndef TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_LOG_SCALE_H_ +#define TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_LOG_SCALE_H_ + +#include +#include + +#ifdef __cplusplus +extern "C" { +#endif + +struct LogScaleState { + int enable_log; + int scale_shift; +}; + +// Applies a fixed point logarithm to the signal and converts it to 16 bit. Note +// that the signal array will be modified. +uint16_t* LogScaleApply(struct LogScaleState* state, uint32_t* signal, + int signal_size, int correction_bits); + +#ifdef __cplusplus +} // extern "C" +#endif + +#endif // TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_LOG_SCALE_H_ diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_scale_util.c b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_scale_util.c new file mode 100644 index 00000000..0e3dd1d1 --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_scale_util.c @@ -0,0 +1,27 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#include "tensorflow/lite/experimental/microfrontend/lib/log_scale_util.h" + +void LogScaleFillConfigWithDefaults(struct LogScaleConfig* config) { + config->enable_log = 1; + config->scale_shift = 6; +} + +int LogScalePopulateState(const struct LogScaleConfig* config, + struct LogScaleState* state) { + state->enable_log = config->enable_log; + state->scale_shift = config->scale_shift; + return 1; +} diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_scale_util.h b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_scale_util.h new file mode 100644 index 00000000..11f7d9ee --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/log_scale_util.h @@ -0,0 +1,45 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#ifndef TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_LOG_SCALE_UTIL_H_ +#define TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_LOG_SCALE_UTIL_H_ + +#include +#include + +#include "tensorflow/lite/experimental/microfrontend/lib/log_scale.h" + +#ifdef __cplusplus +extern "C" { +#endif + +struct LogScaleConfig { + // set to false (0) to disable this module + int enable_log; + // scale results by 2^(scale_shift) + int scale_shift; +}; + +// Populates the LogScaleConfig with "sane" default values. +void LogScaleFillConfigWithDefaults(struct LogScaleConfig* config); + +// Allocates any buffers. +int LogScalePopulateState(const struct LogScaleConfig* config, + struct LogScaleState* state); + +#ifdef __cplusplus +} // extern "C" +#endif + +#endif // TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_LOG_SCALE_UTIL_H_ diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/noise_reduction.c b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/noise_reduction.c new file mode 100644 index 00000000..16b30e66 --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/noise_reduction.c @@ -0,0 +1,51 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#include "tensorflow/lite/experimental/microfrontend/lib/noise_reduction.h" + +#include + +void NoiseReductionApply(struct NoiseReductionState* state, uint32_t* signal) { + int i; + for (i = 0; i < state->num_channels; ++i) { + const uint32_t smoothing = + ((i & 1) == 0) ? state->even_smoothing : state->odd_smoothing; + const uint32_t one_minus_smoothing = (1 << kNoiseReductionBits) - smoothing; + + // Update the estimate of the noise. + const uint32_t signal_scaled_up = signal[i] << state->smoothing_bits; + uint32_t estimate = + (((uint64_t)signal_scaled_up * smoothing) + + ((uint64_t)state->estimate[i] * one_minus_smoothing)) >> + kNoiseReductionBits; + state->estimate[i] = estimate; + + // Make sure that we can't get a negative value for the signal - estimate. + if (estimate > signal_scaled_up) { + estimate = signal_scaled_up; + } + + const uint32_t floor = + ((uint64_t)signal[i] * state->min_signal_remaining) >> + kNoiseReductionBits; + const uint32_t subtracted = + (signal_scaled_up - estimate) >> state->smoothing_bits; + const uint32_t output = subtracted > floor ? subtracted : floor; + signal[i] = output; + } +} + +void NoiseReductionReset(struct NoiseReductionState* state) { + memset(state->estimate, 0, sizeof(*state->estimate) * state->num_channels); +} diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/noise_reduction.h b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/noise_reduction.h new file mode 100644 index 00000000..46d3f52e --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/noise_reduction.h @@ -0,0 +1,46 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#ifndef TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_NOISE_REDUCTION_H_ +#define TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_NOISE_REDUCTION_H_ + +#define kNoiseReductionBits 14 + +#include +#include + +#ifdef __cplusplus +extern "C" { +#endif + +struct NoiseReductionState { + int smoothing_bits; + uint16_t even_smoothing; + uint16_t odd_smoothing; + uint16_t min_signal_remaining; + int num_channels; + uint32_t* estimate; +}; + +// Removes stationary noise from each channel of the signal using a low pass +// filter. +void NoiseReductionApply(struct NoiseReductionState* state, uint32_t* signal); + +void NoiseReductionReset(struct NoiseReductionState* state); + +#ifdef __cplusplus +} // extern "C" +#endif + +#endif // TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_NOISE_REDUCTION_H_ diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/noise_reduction_util.c b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/noise_reduction_util.c new file mode 100644 index 00000000..a6c9234e --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/noise_reduction_util.c @@ -0,0 +1,45 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#include "tensorflow/lite/experimental/microfrontend/lib/noise_reduction_util.h" + +#include + +void NoiseReductionFillConfigWithDefaults(struct NoiseReductionConfig* config) { + config->smoothing_bits = 10; + config->even_smoothing = 0.025; + config->odd_smoothing = 0.06; + config->min_signal_remaining = 0.05; +} + +int NoiseReductionPopulateState(const struct NoiseReductionConfig* config, + struct NoiseReductionState* state, + int num_channels) { + state->smoothing_bits = config->smoothing_bits; + state->odd_smoothing = config->odd_smoothing * (1 << kNoiseReductionBits); + state->even_smoothing = config->even_smoothing * (1 << kNoiseReductionBits); + state->min_signal_remaining = + config->min_signal_remaining * (1 << kNoiseReductionBits); + state->num_channels = num_channels; + state->estimate = calloc(state->num_channels, sizeof(*state->estimate)); + if (state->estimate == NULL) { + fprintf(stderr, "Failed to alloc estimate buffer\n"); + return 0; + } + return 1; +} + +void NoiseReductionFreeStateContents(struct NoiseReductionState* state) { + free(state->estimate); +} diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/noise_reduction_util.h b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/noise_reduction_util.h new file mode 100644 index 00000000..fa555391 --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/noise_reduction_util.h @@ -0,0 +1,50 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#ifndef TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_NOISE_REDUCTION_UTIL_H_ +#define TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_NOISE_REDUCTION_UTIL_H_ + +#include "tensorflow/lite/experimental/microfrontend/lib/noise_reduction.h" + +#ifdef __cplusplus +extern "C" { +#endif + +struct NoiseReductionConfig { + // scale the signal up by 2^(smoothing_bits) before reduction + int smoothing_bits; + // smoothing coefficient for even-numbered channels + float even_smoothing; + // smoothing coefficient for odd-numbered channels + float odd_smoothing; + // fraction of signal to preserve (1.0 disables this module) + float min_signal_remaining; +}; + +// Populates the NoiseReductionConfig with "sane" default values. +void NoiseReductionFillConfigWithDefaults(struct NoiseReductionConfig* config); + +// Allocates any buffers. +int NoiseReductionPopulateState(const struct NoiseReductionConfig* config, + struct NoiseReductionState* state, + int num_channels); + +// Frees any allocated buffers. +void NoiseReductionFreeStateContents(struct NoiseReductionState* state); + +#ifdef __cplusplus +} // extern "C" +#endif + +#endif // TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_NOISE_REDUCTION_UTIL_H_ diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control.c b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control.c new file mode 100644 index 00000000..22d58767 --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control.c @@ -0,0 +1,56 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#include "tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control.h" + +#include "tensorflow/lite/experimental/microfrontend/lib/bits.h" + +int16_t WideDynamicFunction(const uint32_t x, const int16_t* lut) { + if (x <= 2) { + return lut[x]; + } + + const int16_t interval = MostSignificantBit32(x); + lut += 4 * interval - 6; + + const int16_t frac = + ((interval < 11) ? (x << (11 - interval)) : (x >> (interval - 11))) & + 0x3FF; + + int32_t result = ((int32_t)lut[2] * frac) >> 5; + result += (int32_t)((uint32_t)lut[1] << 5); + result *= frac; + result = (result + (1 << 14)) >> 15; + result += lut[0]; + return (int16_t)result; +} + +uint32_t PcanShrink(const uint32_t x) { + if (x < (2 << kPcanSnrBits)) { + return (x * x) >> (2 + 2 * kPcanSnrBits - kPcanOutputBits); + } else { + return (x >> (kPcanSnrBits - kPcanOutputBits)) - (1 << kPcanOutputBits); + } +} + +void PcanGainControlApply(struct PcanGainControlState* state, + uint32_t* signal) { + int i; + for (i = 0; i < state->num_channels; ++i) { + const uint32_t gain = + WideDynamicFunction(state->noise_estimate[i], state->gain_lut); + const uint32_t snr = ((uint64_t)signal[i] * gain) >> state->snr_shift; + signal[i] = PcanShrink(snr); + } +} diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control.h b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control.h new file mode 100644 index 00000000..3f6222be --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control.h @@ -0,0 +1,47 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#ifndef TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_PCAN_GAIN_CONTROL_H_ +#define TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_PCAN_GAIN_CONTROL_H_ + +#include +#include + +#define kPcanSnrBits 12 +#define kPcanOutputBits 6 + +#ifdef __cplusplus +extern "C" { +#endif + +// Details at https://research.google/pubs/pub45911.pdf +struct PcanGainControlState { + int enable_pcan; + uint32_t* noise_estimate; + int num_channels; + int16_t* gain_lut; + int32_t snr_shift; +}; + +int16_t WideDynamicFunction(const uint32_t x, const int16_t* lut); + +uint32_t PcanShrink(const uint32_t x); + +void PcanGainControlApply(struct PcanGainControlState* state, uint32_t* signal); + +#ifdef __cplusplus +} // extern "C" +#endif + +#endif // TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_PCAN_GAIN_CONTROL_H_ diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control_util.c b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control_util.c new file mode 100644 index 00000000..e850d439 --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control_util.c @@ -0,0 +1,92 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#include "tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control_util.h" + +#include +#include + +#define kint16max 0x00007FFF + +void PcanGainControlFillConfigWithDefaults( + struct PcanGainControlConfig* config) { + config->enable_pcan = 0; + config->strength = 0.95; + config->offset = 80.0; + config->gain_bits = 21; +} + +int16_t PcanGainLookupFunction(const struct PcanGainControlConfig* config, + int32_t input_bits, uint32_t x) { + const float x_as_float = ((float)x) / ((uint32_t)1 << input_bits); + const float gain_as_float = + ((uint32_t)1 << config->gain_bits) * + powf(x_as_float + config->offset, -config->strength); + + if (gain_as_float > kint16max) { + return kint16max; + } + return (int16_t)(gain_as_float + 0.5f); +} + +int PcanGainControlPopulateState(const struct PcanGainControlConfig* config, + struct PcanGainControlState* state, + uint32_t* noise_estimate, + const int num_channels, + const uint16_t smoothing_bits, + const int32_t input_correction_bits) { + state->enable_pcan = config->enable_pcan; + if (!state->enable_pcan) { + return 1; + } + state->noise_estimate = noise_estimate; + state->num_channels = num_channels; + state->gain_lut = malloc(kWideDynamicFunctionLUTSize * sizeof(int16_t)); + if (state->gain_lut == NULL) { + fprintf(stderr, "Failed to allocate gain LUT\n"); + return 0; + } + state->snr_shift = config->gain_bits - input_correction_bits - kPcanSnrBits; + + const int32_t input_bits = smoothing_bits - input_correction_bits; + state->gain_lut[0] = PcanGainLookupFunction(config, input_bits, 0); + state->gain_lut[1] = PcanGainLookupFunction(config, input_bits, 1); + state->gain_lut -= 6; + int interval; + for (interval = 2; interval <= kWideDynamicFunctionBits; ++interval) { + const uint32_t x0 = (uint32_t)1 << (interval - 1); + const uint32_t x1 = x0 + (x0 >> 1); + const uint32_t x2 = + (interval == kWideDynamicFunctionBits) ? x0 + (x0 - 1) : 2 * x0; + + const int16_t y0 = PcanGainLookupFunction(config, input_bits, x0); + const int16_t y1 = PcanGainLookupFunction(config, input_bits, x1); + const int16_t y2 = PcanGainLookupFunction(config, input_bits, x2); + + const int32_t diff1 = (int32_t)y1 - y0; + const int32_t diff2 = (int32_t)y2 - y0; + const int32_t a1 = 4 * diff1 - diff2; + const int32_t a2 = diff2 - a1; + + state->gain_lut[4 * interval] = y0; + state->gain_lut[4 * interval + 1] = (int16_t)a1; + state->gain_lut[4 * interval + 2] = (int16_t)a2; + } + state->gain_lut += 6; + return 1; +} + +void PcanGainControlFreeStateContents(struct PcanGainControlState* state) { + free(state->gain_lut); +} diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control_util.h b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control_util.h new file mode 100644 index 00000000..d4bfaa2e --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control_util.h @@ -0,0 +1,57 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#ifndef TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_PCAN_GAIN_CONTROL_UTIL_H_ +#define TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_PCAN_GAIN_CONTROL_UTIL_H_ + +#include "tensorflow/lite/experimental/microfrontend/lib/pcan_gain_control.h" + +#define kWideDynamicFunctionBits 32 +#define kWideDynamicFunctionLUTSize (4 * kWideDynamicFunctionBits - 3) + +#ifdef __cplusplus +extern "C" { +#endif + +struct PcanGainControlConfig { + // set to false (0) to disable this module + int enable_pcan; + // gain normalization exponent (0.0 disables, 1.0 full strength) + float strength; + // positive value added in the normalization denominator + float offset; + // number of fractional bits in the gain + int gain_bits; +}; + +void PcanGainControlFillConfigWithDefaults( + struct PcanGainControlConfig* config); + +int16_t PcanGainLookupFunction(const struct PcanGainControlConfig* config, + int32_t input_bits, uint32_t x); + +int PcanGainControlPopulateState(const struct PcanGainControlConfig* config, + struct PcanGainControlState* state, + uint32_t* noise_estimate, + const int num_channels, + const uint16_t smoothing_bits, + const int32_t input_correction_bits); + +void PcanGainControlFreeStateContents(struct PcanGainControlState* state); + +#ifdef __cplusplus +} // extern "C" +#endif + +#endif // TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_PCAN_GAIN_CONTROL_UTIL_H_ diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/window.c b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/window.c new file mode 100644 index 00000000..10da6762 --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/window.c @@ -0,0 +1,70 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#include "tensorflow/lite/experimental/microfrontend/lib/window.h" + +#include + +int WindowProcessSamples(struct WindowState* state, const int16_t* samples, + size_t num_samples, size_t* num_samples_read) { + const int size = state->size; + + // Copy samples from the samples buffer over to our local input. + size_t max_samples_to_copy = state->size - state->input_used; + if (max_samples_to_copy > num_samples) { + max_samples_to_copy = num_samples; + } + memcpy(state->input + state->input_used, samples, + max_samples_to_copy * sizeof(*samples)); + *num_samples_read = max_samples_to_copy; + state->input_used += max_samples_to_copy; + + if (state->input_used < state->size) { + // We don't have enough samples to compute a window. + return 0; + } + + // Apply the window to the input. + const int16_t* coefficients = state->coefficients; + const int16_t* input = state->input; + int16_t* output = state->output; + int i; + int16_t max_abs_output_value = 0; + for (i = 0; i < size; ++i) { + int16_t new_value = + (((int32_t)*input++) * *coefficients++) >> kFrontendWindowBits; + *output++ = new_value; + if (new_value < 0) { + new_value = -new_value; + } + if (new_value > max_abs_output_value) { + max_abs_output_value = new_value; + } + } + // Shuffle the input down by the step size, and update how much we have used. + memmove(state->input, state->input + state->step, + sizeof(*state->input) * (state->size - state->step)); + state->input_used -= state->step; + state->max_abs_output_value = max_abs_output_value; + + // Indicate that the output buffer is valid for the next stage. + return 1; +} + +void WindowReset(struct WindowState* state) { + memset(state->input, 0, state->size * sizeof(*state->input)); + memset(state->output, 0, state->size * sizeof(*state->output)); + state->input_used = 0; + state->max_abs_output_value = 0; +} diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/window.h b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/window.h new file mode 100644 index 00000000..bad81514 --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/window.h @@ -0,0 +1,49 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#ifndef TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_WINDOW_H_ +#define TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_WINDOW_H_ + +#include +#include + +#define kFrontendWindowBits 12 + +#ifdef __cplusplus +extern "C" { +#endif + +struct WindowState { + size_t size; + int16_t* coefficients; + size_t step; + + int16_t* input; + size_t input_used; + int16_t* output; + int16_t max_abs_output_value; +}; + +// Applies a window to the samples coming in, stepping forward at the given +// rate. +int WindowProcessSamples(struct WindowState* state, const int16_t* samples, + size_t num_samples, size_t* num_samples_read); + +void WindowReset(struct WindowState* state); + +#ifdef __cplusplus +} // extern "C" +#endif + +#endif // TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_WINDOW_H_ diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/window_util.c b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/window_util.c new file mode 100644 index 00000000..6fdffdc4 --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/window_util.c @@ -0,0 +1,73 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#include "tensorflow/lite/experimental/microfrontend/lib/window_util.h" + +#include +#include +#include +#include + +// Some platforms don't have M_PI +#ifndef M_PI +#define M_PI 3.14159265358979323846 +#endif + +void WindowFillConfigWithDefaults(struct WindowConfig* config) { + config->size_ms = 25; + config->step_size_ms = 10; +} + +int WindowPopulateState(const struct WindowConfig* config, + struct WindowState* state, int sample_rate) { + state->size = config->size_ms * sample_rate / 1000; + state->step = config->step_size_ms * sample_rate / 1000; + + state->coefficients = malloc(state->size * sizeof(*state->coefficients)); + if (state->coefficients == NULL) { + fprintf(stderr, "Failed to allocate window coefficients\n"); + return 0; + } + + // Populate the window values. + const float arg = (float)M_PI * 2.0f / ((float)state->size); + size_t i; + for (i = 0; i < state->size; ++i) { + float float_value = 0.5f - (0.5f * cosf(arg * (i + 0.5f))); + // Scale it to fixed point and round it. + state->coefficients[i] = + floorf(float_value * (1 << kFrontendWindowBits) + 0.5f); + } + + state->input_used = 0; + state->input = malloc(state->size * sizeof(*state->input)); + if (state->input == NULL) { + fprintf(stderr, "Failed to allocate window input\n"); + return 0; + } + + state->output = malloc(state->size * sizeof(*state->output)); + if (state->output == NULL) { + fprintf(stderr, "Failed to allocate window output\n"); + return 0; + } + + return 1; +} + +void WindowFreeStateContents(struct WindowState* state) { + free(state->coefficients); + free(state->input); + free(state->output); +} diff --git a/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/window_util.h b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/window_util.h new file mode 100644 index 00000000..68e4de9e --- /dev/null +++ b/third_party/microfrontend/tensorflow/lite/experimental/microfrontend/lib/window_util.h @@ -0,0 +1,45 @@ +/* Copyright 2018 The TensorFlow Authors. All Rights Reserved. + +Licensed under the Apache License, Version 2.0 (the "License"); +you may not use this file except in compliance with the License. +You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + +Unless required by applicable law or agreed to in writing, software +distributed under the License is distributed on an "AS IS" BASIS, +WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +See the License for the specific language governing permissions and +limitations under the License. +==============================================================================*/ +#ifndef TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_WINDOW_UTIL_H_ +#define TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_WINDOW_UTIL_H_ + +#include "tensorflow/lite/experimental/microfrontend/lib/window.h" + +#ifdef __cplusplus +extern "C" { +#endif + +struct WindowConfig { + // length of window frame in milliseconds + size_t size_ms; + // length of step for next frame in milliseconds + size_t step_size_ms; +}; + +// Populates the WindowConfig with "sane" default values. +void WindowFillConfigWithDefaults(struct WindowConfig* config); + +// Allocates any buffers. +int WindowPopulateState(const struct WindowConfig* config, + struct WindowState* state, int sample_rate); + +// Frees any allocated buffers. +void WindowFreeStateContents(struct WindowState* state); + +#ifdef __cplusplus +} // extern "C" +#endif + +#endif // TENSORFLOW_LITE_EXPERIMENTAL_MICROFRONTEND_LIB_WINDOW_UTIL_H_ From d022a875f7182aa8644d608099f4bb29e9e81431 Mon Sep 17 00:00:00 2001 From: Brad Barnett <127794626+bdbarnett@users.noreply.github.com> Date: Thu, 29 Oct 2026 03:43:10 -0500 Subject: [PATCH 2/3] audiofeatures: taps without libm's sin, no fused multiply-adds The resampler's taps came from sin(), which may differ in the last bit between C libraries, and the unix MicroPython build in CI made different 48 kHz features from the CPython wheel. sin(pi x) is now a series in exact IEEE arithmetic, and both files take audiodsp_fp_contract.h like the rest of src/shared. The symbol test needs GNU nm, so it runs on Linux only. --- src/shared/audiodsp_features.c | 31 ++++++++++++++++++++++++++++- src/shared/audiodsp_microfrontend.c | 2 ++ tests/test_cpython_audiofeatures.py | 3 ++- 3 files changed, 34 insertions(+), 2 deletions(-) diff --git a/src/shared/audiodsp_features.c b/src/shared/audiodsp_features.c index 91aea553..8bca11c1 100644 --- a/src/shared/audiodsp_features.c +++ b/src/shared/audiodsp_features.c @@ -2,6 +2,8 @@ // SPDX-License-Identifier: MIT // SPDX-FileCopyrightText: Copyright (c) 2026 PyDevices +#include "shared/audiodsp_fp_contract.h" + #include "shared/audiodsp_features.h" #include @@ -61,6 +63,33 @@ static double rs_bessel_i0(double x) { return sum; } +// sin(pi * x) from + - * / alone. libm's sin() is allowed to differ in the +// last bit between libraries, and one bit in a tap can move a sample, so a +// board's newlib and a desktop's glibc would make different taps. These +// operations are exact IEEE ones everywhere. +static double rs_sinpi(double x) { + double sign = 1.0; + if (x < 0.0) { + x = -x; + sign = -1.0; + } + x -= 2.0 * (double)(long)(x / 2.0); // now 0 <= x < 2, exactly + if (x >= 1.0) { + x -= 1.0; + sign = -sign; + } + if (x > 0.5) { + x = 1.0 - x; // sin(pi x) = sin(pi (1 - x)) + } + // Taylor series of sin(t), t = pi x <= pi / 2: 13 terms reach 1e-17. + double t = RS_PI * x, t2 = t * t, term = t, sum = t; + for (int k = 1; k < 14; k++) { + term *= -t2 / ((double)(2 * k) * (double)(2 * k + 1)); + sum += term; + } + return sign * sum; +} + // The kernel at u output samples from the centre. static double rs_kernel(double u) { u = fabs(u); @@ -69,7 +98,7 @@ static double rs_kernel(double u) { } double r = u / RS_ZEROS; double w = rs_bessel_i0(RS_BETA * sqrt(1.0 - r * r)) / rs_bessel_i0(RS_BETA); - double s = u < 1e-12 ? 1.0 : sin(RS_PI * u) / (RS_PI * u); + double s = u < 1e-12 ? 1.0 : rs_sinpi(u) / (RS_PI * u); return s * w; } diff --git a/src/shared/audiodsp_microfrontend.c b/src/shared/audiodsp_microfrontend.c index 2cd43b13..78cb1389 100644 --- a/src/shared/audiodsp_microfrontend.c +++ b/src/shared/audiodsp_microfrontend.c @@ -23,6 +23,8 @@ // include path, because the library names its headers by their full // "tensorflow/lite/..." path. +#include "shared/audiodsp_fp_contract.h" + #include #include #include diff --git a/tests/test_cpython_audiofeatures.py b/tests/test_cpython_audiofeatures.py index 25cb1626..158e2cd1 100644 --- a/tests/test_cpython_audiofeatures.py +++ b/tests/test_cpython_audiofeatures.py @@ -20,6 +20,7 @@ import os import shutil import subprocess +import sys import unittest from array import array @@ -253,7 +254,7 @@ def test_context_manager(self): fe.stats() -@unittest.skipUnless(shutil.which("nm") and _audiodsp.__file__.endswith(".so"), "needs nm and an ELF build") +@unittest.skipUnless(sys.platform.startswith("linux") and shutil.which("nm"), "needs GNU nm and an ELF build") class SymbolTest(unittest.TestCase): """Every symbol of the microfrontend and kissfft is renamed, so another copy links.""" From ebf7fe1f4b3a8747729d1d2ab30e273e76545787 Mon Sep 17 00:00:00 2001 From: Brad Barnett <127794626+bdbarnett@users.noreply.github.com> Date: Thu, 29 Oct 2026 03:47:21 -0500 Subject: [PATCH 3/3] ROADMAP: the frontend's flash cost, measured at -O2 on the P4's toolchain --- ROADMAP.md | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/ROADMAP.md b/ROADMAP.md index 078aeaff..5d2dda40 100644 --- a/ROADMAP.md +++ b/ROADMAP.md @@ -17,8 +17,8 @@ audio pump holding up under live use such as a stompbox or a stage rig. without a click when you switch it in mid-note. - Voice activity in `audiofeatures`: a flag per slice saying whether it holds speech, from the frontend's own noise estimate, so a wake word model need - only run while someone is talking. The frontend costs about 14 KB of flash - on an ESP32-P4 today; this should add little to it. + only run while someone is talking. The frontend and resampler take about + 16 KB of flash on an ESP32-P4 today; this should add little to it. Bugs, and things you need that don't work yet, go to [issues](https://github.com/PyDevices/audiodsp/issues).