From 11da11b17eb6f7749378bbc31afdf37d6499194b Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?J=C3=BCrgen=20Hock?= Date: Tue, 10 Jan 2023 18:47:10 +0100 Subject: [PATCH] Add more chroma features #2 --- examples/chromatest.py | 2 +- src/python/qdft/chroma.py | 80 +++++++++------------------------ src/python/qdft/fafe.py | 95 +++++++++++++++++++++++++++++++++++++++ src/python/qdft/qdft.py | 10 +---- 4 files changed, 117 insertions(+), 70 deletions(-) create mode 100644 src/python/qdft/fafe.py diff --git a/examples/chromatest.py b/examples/chromatest.py index f32c4c8..6ec7657 100644 --- a/examples/chromatest.py +++ b/examples/chromatest.py @@ -43,7 +43,7 @@ def main(): # 2) analyze previously built test signal - chroma = Chroma(sr, cp) + chroma = Chroma(sr, cp, feature='cent') chromagram = chroma.chroma(x)[-1] # pick last chromagram entry diff --git a/src/python/qdft/chroma.py b/src/python/qdft/chroma.py index d384d37..81fa17a 100644 --- a/src/python/qdft/chroma.py +++ b/src/python/qdft/chroma.py @@ -9,13 +9,14 @@ import numpy +from .fafe import QFAFE from .qdft import QDFT from .scale import Scale class Chroma: - def __init__(self, samplerate, concertpitch=440, bandwidth=('A0', 'C#8'), decibel=True): + def __init__(self, samplerate, concertpitch=440, bandwidth=('A0', 'C#8'), decibel=True, feature=None): scale = Scale(concertpitch) @@ -39,6 +40,7 @@ def __init__(self, samplerate, concertpitch=440, bandwidth=('A0', 'C#8'), decibe self.concertpitch = concertpitch self.bandwidth = bandwidth self.decibel = decibel + self.feature = feature self.semitones = semitones self.frequencies = frequencies self.notes = notes @@ -49,77 +51,35 @@ def __init__(self, samplerate, concertpitch=440, bandwidth=('A0', 'C#8'), decibe def chroma(self, samples): - stash = { 'cents': None } - - def analysis(dfts, mode=None): - stash['cents'] = self.analyze(dfts, mode) - - # TODO: analyze raw dfts - # dfts = self.qdft.qdft(samples, analysis) - - # TODO: analyze windowed dfts dfts = self.qdft.qdft(samples) - stash['cents'] = self.analyze(dfts, 'p') - magnis = numpy.abs(dfts) - cents = stash['cents'] + magnitudes = numpy.abs(dfts) + features = None if self.decibel: with numpy.errstate(all='ignore'): - magnis = 20 * numpy.log10(magnis) - - chromagram = magnis + 1j * cents - - chromagram = chromagram[..., ::2] - assert chromagram.shape[-1] == self.frequencies.shape[-1] - - return chromagram + magnitudes = 20 * numpy.log10(magnitudes) - def analyze(self, dfts, mode=None): + if str(self.feature).lower() in 'phase': - l = numpy.roll(dfts, +1, axis=-1) - m = dfts - r = numpy.roll(dfts, -1, axis=-1) + features = numpy.angle(dfts) - if mode is None: + if str(self.feature).lower() in 'hz': - with numpy.errstate(all='ignore'): - drifts = -numpy.real((r - l) / (2 * m - r - l)) - - elif str(mode).lower() == 'p': + fafe = QFAFE(self.qdft) + features = fafe.hz(dfts) - p = 1.36 + if str(self.feature).lower() in 'cent': - l = numpy.abs(l) - m = numpy.abs(m) - r = numpy.abs(r) + fafe = QFAFE(self.qdft) + features = fafe.cent(dfts) - with numpy.errstate(all='ignore'): - drifts = p * (r - l) / (m + r + l) + chromagram = (magnitudes + 1j * features) \ + if features is not None \ + else magnitudes - elif str(mode).lower() == 'q': - - q = 0.55 - - with numpy.errstate(all='ignore'): - drifts = -numpy.real(q * (r - l) / (2 * m + r + l)) - - else: - - drifts = numpy.zeros(dfts.shape) - - drifts[..., 0] = 0 - drifts[..., -1] = 0 - - oldfreqs = self.qdft.frequencies - oldbins = numpy.arange(oldfreqs.size) - newbins = oldbins + drifts - # TODO: is approximation possible? https://en.wikipedia.org/wiki/Cent_(music) - newfreqs = self.qdft.bandwidth[0] * numpy.power(2, newbins / self.qdft.resolution) - # TODO: does interp make sense? - # newfreqs = numpy.interp(newbins, oldbins, oldfreqs) - - cents = 1200 * numpy.log2(newfreqs / oldfreqs) + chromagram = chromagram[..., ::2] + assert chromagram.shape[-1] == self.frequencies.shape[-1] - return cents + return chromagram diff --git a/src/python/qdft/fafe.py b/src/python/qdft/fafe.py new file mode 100644 index 0000000..414c7f7 --- /dev/null +++ b/src/python/qdft/fafe.py @@ -0,0 +1,95 @@ +""" +Copyright (c) 2023 Juergen Hock + +SPDX-License-Identifier: MIT + +Fast, Accurate Frequency Estimator according to [1]. + +[1] Eric Jacobsen and Peter Kootsookos + Fast, Accurate Frequency Estimators + IEEE Signal Processing Magazine (2007) + https://ieeexplore.ieee.org/document/4205098 + +Source: https://github.com/jurihock/qdft +""" + + +import numpy + + +class FAFE: + + def __init__(self, mode=None): + + self.mode = mode + + def __call__(self, dfts): + + dfts = numpy.atleast_2d(dfts) + assert dfts.ndim == 2 + + l = numpy.roll(dfts, +1, axis=-1) + m = dfts + r = numpy.roll(dfts, -1, axis=-1) + + if self.mode is None: + + with numpy.errstate(all='ignore'): + drifts = -numpy.real((r - l) / (2 * m - r - l)) + + elif str(self.mode).lower() == 'p': + + p = 1.36 # TODO: hann? + + l = numpy.abs(l) + m = numpy.abs(m) + r = numpy.abs(r) + + with numpy.errstate(all='ignore'): + drifts = p * (r - l) / (m + r + l) + + elif str(self.mode).lower() == 'q': + + q = 0.55 # TODO: hann? + + with numpy.errstate(all='ignore'): + drifts = -numpy.real(q * (r - l) / (2 * m + r + l)) + + else: + + drifts = numpy.zeros(dfts.shape) + + drifts[..., 0] = 0 + drifts[..., -1] = 0 + + return drifts + + +class QFAFE: + + def __init__(self, qdft): + + self.qdft = qdft + self.fafe = FAFE('p' if qdft.window is not None else None) + + def hz(self, dfts): + + oldfreqs = self.qdft.frequencies + + oldbins = numpy.arange(oldfreqs.size) + newbins = oldbins + self.fafe(dfts) + + # TODO: is approximation possible? https://en.wikipedia.org/wiki/Cent_(music) + newfreqs = self.qdft.bandwidth[0] * numpy.power(2, newbins / self.qdft.resolution) + + # TODO: does interp make sense? + # newfreqs = numpy.interp(newbins, oldbins, oldfreqs) + + return newfreqs + + def cent(self, dfts): + + newfreqs = self.hz(dfts) + oldfreqs = self.qdft.frequencies + + return 1200 * numpy.log2(newfreqs / oldfreqs) diff --git a/src/python/qdft/qdft.py b/src/python/qdft/qdft.py index 94143fd..b9ea299 100644 --- a/src/python/qdft/qdft.py +++ b/src/python/qdft/qdft.py @@ -67,7 +67,7 @@ def __init__(self, samplerate, bandwidth, resolution=24, latency=0, window=(+0.5 self.outputs = outputs self.kernels = kernels - def qdft(self, samples, callback=None): + def qdft(self, samples): samples = numpy.atleast_1d(samples).astype(float) assert samples.ndim == 1 @@ -103,14 +103,6 @@ def qdft(self, samples, callback=None): outputs[k] = dfts[k][-1] - if callback is not None: - - assert callable(callback) - - dfts[0].setflags(write=False) - try: callback(dfts[0]) - finally: dfts[0].setflags(write=True) - if window is not None: a, b = window[0], window[1] / 2