Higher quality wavetable mipmaps

src/sfizz/Wavetables.cpp

@@ -283,65 +283,70 @@ const HarmonicProfile& HarmonicProfile::getSquare()
 }
 
 //------------------------------------------------------------------------------
-constexpr unsigned WavetableRange::countOctaves;
+constexpr unsigned MipmapRange::N;
-constexpr float WavetableRange::frequencyScaleFactor;
+constexpr float MipmapRange::F1;
+constexpr float MipmapRange::FN;
 
-unsigned WavetableRange::getOctaveForFrequency(float f)
+const float MipmapRange::K = 1.0 / F1;
+const float MipmapRange::LogB = std::log(FN / F1) / (N - 1);
+
+const std::array<float, 1024> MipmapRange::FrequencyToIndex = []()
 {
-    int oct = fp_exponent(frequencyScaleFactor * f);
+    std::array<float, 1024> table;
-    return clamp<int>(oct, 0, countOctaves - 1);
-}
 
-static const auto octaveForFrequencyTable = []()
+    for (unsigned i = 0; i < table.size() - 1; ++i) {
-{
+        double r = i * (1.0 / (table.size() - 1));
-    static constexpr unsigned N = 1024;
+        double f = F1 + r * (FN - F1);
-    std::array<float, N> table;
+        double t = std::log(K * f) / LogB;
-
+        table[i] = clamp<float>(t, 0, N - 1);
-    constexpr double fmin = 1 / WavetableRange::frequencyScaleFactor;
-    constexpr double fmax = (1 << (WavetableRange::countOctaves - 1)) / WavetableRange::frequencyScaleFactor;
-
-    for (unsigned i = 0; i < N; ++i) {
-        double f = fmin + (i * (1.0 / (N - 1))) * (fmax - fmin);
-        table[i] = std::log2(f * WavetableRange::frequencyScaleFactor);
     }
+    // ensure the last element to be exact
+    table[table.size() - 1] = N - 1;
 
     return table;
 }();
 
-float WavetableRange::getFractionalOctaveForFrequency(float f)
+float MipmapRange::getIndexForFrequency(float f)
 {
-    static constexpr unsigned N = octaveForFrequencyTable.size();
+    static constexpr unsigned tableSize = FrequencyToIndex.size();
 
-    constexpr double fmin = 1 / WavetableRange::frequencyScaleFactor;
+    float pos = (f - F1) * ((tableSize - 1) / static_cast<float>(FN - F1));
-    constexpr double fmax = (1 << (WavetableRange::countOctaves - 1)) / WavetableRange::frequencyScaleFactor;
+    pos = clamp<float>(pos, 0, tableSize - 1);
 
-    float pos = (f - fmin) * ((N - 1) / static_cast<float>(fmax - fmin));
     int index1 = static_cast<int>(pos);
+    int index2 = std::min<int>(index1 + 1, tableSize - 1);
     float frac = pos - index1;
-    index1 = clamp<int>(index1, 0, N - 1);
-    int index2 = std::min<int>(index1 + 1, N - 1);
 
-    return (1.0f - frac) * octaveForFrequencyTable[index1] +
+    return (1.0f - frac) * FrequencyToIndex[index1] +
-        frac * octaveForFrequencyTable[index2];
+        frac * FrequencyToIndex[index2];
 }
 
-WavetableRange WavetableRange::getRangeForOctave(int o)
+const std::array<float, MipmapRange::N + 1> MipmapRange::IndexToStartFrequency = []()
 {
-    WavetableRange range;
+    std::array<float, N + 1> table;
+    for (unsigned t = 0; t < N; ++t)
+        table[t] = std::exp(t * LogB) / K;
+    // end value for final table
+    table[N] = 22050.0;
 
-    Fraction<uint64_t> mant = fp_mantissa(0.0f);
+    return table;
-    float k = 1.0f / frequencyScaleFactor;
+}();
 
-    range.minFrequency = k * fp_from_parts<float>(0, o, 0);
+MipmapRange MipmapRange::getRangeForIndex(int o)
-    range.maxFrequency = k * fp_from_parts<float>(0, o, mant.den - 1);
+{
+    o = clamp<int>(o, 0, N - 1);
+
+    MipmapRange range;
+    range.minFrequency = IndexToStartFrequency[o];
+    range.maxFrequency = IndexToStartFrequency[o + 1];
 
     return range;
 }
 
-WavetableRange WavetableRange::getRangeForFrequency(float f)
+MipmapRange MipmapRange::getRangeForFrequency(float f)
 {
-    int oct = getOctaveForFrequency(f);
+    int index = static_cast<int>(getIndexForFrequency(f));
-    return getRangeForOctave(oct);
+    return getRangeForIndex(index);
 }
 
 //------------------------------------------------------------------------------
@@ -356,7 +361,7 @@ WavetableMulti WavetableMulti::createForHarmonicProfile(
     wm.allocateStorage(tableSize);
 
     for (unsigned m = 0; m < numTables; ++m) {
-        WavetableRange range = WavetableRange::getRangeForOctave(m);
+        MipmapRange range = MipmapRange::getRangeForIndex(m);
 
         double freq = range.maxFrequency;
 

src/sfizz/Wavetables.h

@@ -11,6 +11,7 @@
 #include "MathHelpers.h"
 #include <absl/types/span.h>
 #include <absl/container/flat_hash_map.h>
+#include <array>
 #include <memory>
 #include <complex>
 
@@ -57,6 +58,11 @@ public:
      */
     void setQuality(int q) { _quality = q; }
 
+    /**
+       Get the quality of this oscillator. (cf. `oscillator_quality`)
+    */
+    int quality() const { return _quality; }
+
     /**
        Compute a cycle of the oscillator, with constant frequency.
      */
@@ -117,36 +123,41 @@ public:
 };
 
 /**
-   A helper to select ranges of a multi-sampled oscillator, according to the
+   A helper to select ranges of a mip-mapped wave, according to the
    frequency of an oscillator.
 
    The ranges are identified by octave numbers; not octaves in a musical sense,
    but as logarithmic divisions of the frequency range.
  */
-class WavetableRange {
+class MipmapRange {
 public:
     float minFrequency = 0;
     float maxFrequency = 0;
 
-    static constexpr unsigned countOctaves = 10;
+    // number of tables in the mipmap
-    static constexpr float frequencyScaleFactor = 0.05;
+    static constexpr unsigned N = 24;
+    // start frequency of the first table in the mipmap
+    static constexpr float F1 = 20.0;
+    // start frequency of the last table in the mipmap
+    static constexpr float FN = 12000.0;
 
-    static unsigned getOctaveForFrequency(float f);
+    static float getIndexForFrequency(float f);
-    static float getFractionalOctaveForFrequency(float f);
+    static MipmapRange getRangeForIndex(int o);
-    static WavetableRange getRangeForOctave(int o);
+    static MipmapRange getRangeForFrequency(float f);
-    static WavetableRange getRangeForFrequency(float f);
 
-    // Note: using the frequency factor 0.05, octaves are as follows:
+    // the frequency mapping of the mipmap is defined by formula:
-    //     octave 0: 20 Hz - 40 Hz
+    //     T(f) = log(k*f)/log(b)
-    //     octave 1: 40 Hz - 80 Hz
+    // - T is the table number, converted to index by rounding down
-    //     octave 2: 80 Hz - 160 Hz
+    // - f is the oscillation frequency
-    //     octave 3: 160 Hz - 320 Hz
+    // - k and b are adjustment parameters according to constant parameters
-    //     octave 4: 320 Hz - 640 Hz
+    //     k = 1/F1
-    //     octave 5: 640 Hz - 1280 Hz
+    //     b = exp(log(FN/F1)/(N-1))
-    //     octave 6: 1280 Hz - 2560 Hz
+
-    //     octave 7: 2560 Hz - 5120 Hz
+    static const float K;
-    //     octave 8: 5120 Hz - 10240 Hz
+    static const float LogB;
-    //     octave 9: 10240 Hz - 20480 Hz
+
+    static const std::array<float, 1024> FrequencyToIndex;
+    static const std::array<float, N + 1> IndexToStartFrequency;
 };
 
 /**
@@ -159,7 +170,7 @@ public:
     unsigned tableSize() const { return _tableSize; }
 
     // number of tables in the multisample
-    static constexpr unsigned numTables() { return WavetableRange::countOctaves; }
+    static constexpr unsigned numTables() { return MipmapRange::N; }
 
     // get the N-th table in the multisample
     absl::Span<const float> getTable(unsigned index) const
@@ -170,7 +181,7 @@ public:
     // get the table which is adequate for a given playback frequency
     absl::Span<const float> getTableForFrequency(float freq) const
     {
-        return getTable(WavetableRange::getOctaveForFrequency(freq));
+        return getTable(MipmapRange::getIndexForFrequency(freq));
     }
 
     // adjacent tables with interpolation factor between them
@@ -186,15 +197,15 @@ public:
         DualTable dt;
         int index = static_cast<int>(position);
         dt.delta = position - index;
-        dt.table1 = getTablePointer(clamp<int>(index, 0, WavetableRange::countOctaves - 1));
+        dt.table1 = getTablePointer(clamp<int>(index, 0, MipmapRange::N - 1));
-        dt.table2 = getTablePointer(clamp<int>(index + 1, 0, WavetableRange::countOctaves - 1));
+        dt.table2 = getTablePointer(clamp<int>(index + 1, 0, MipmapRange::N - 1));
         return dt;
     }
 
     // get the pair of tables for the given playback frequency (range checked)
     DualTable getInterpolationPairForFrequency(float freq) const
     {
-        float position = WavetableRange::getFractionalOctaveForFrequency(freq);
+        float position = MipmapRange::getIndexForFrequency(freq);
         return getInterpolationPair(position);
     }
 

tests/DemoWavetables.cpp

@@ -36,6 +36,7 @@ private:
 
 private:
     void valueChangedWave(int value);
+    void valueChangedQuality(int value);
     void buttonClickedPlaySweep();
 
 private:
@@ -46,6 +47,7 @@ private:
     sfz::WavetableOscillator fOsc;
     unsigned fWavePlaying = 0;
     std::atomic<int> fNewWavePending { -1 };
+    std::atomic<int> fNewQualityPending { -1 };
     std::atomic<bool> fStartNewSweep { false };
 
     static constexpr float sweepMin = 0.0;
@@ -88,13 +90,13 @@ bool DemoApp::initSound()
     fTmpFrequency.reset(new float[bufferSize]);
 
     fMulti[0] = sfz::WavetableMulti::createForHarmonicProfile(
-        sfz::HarmonicProfile::getSine(), 1.0, 2048);
+        sfz::HarmonicProfile::getSine(), sfz::config::amplitudeSine, 2048);
     fMulti[1] = sfz::WavetableMulti::createForHarmonicProfile(
-        sfz::HarmonicProfile::getTriangle(), 1.0, 2048);
+        sfz::HarmonicProfile::getTriangle(), sfz::config::amplitudeTriangle, 2048);
     fMulti[2] = sfz::WavetableMulti::createForHarmonicProfile(
-        sfz::HarmonicProfile::getSaw(), 1.0, 2048);
+        sfz::HarmonicProfile::getSaw(), sfz::config::amplitudeSaw, 2048);
     fMulti[3] = sfz::WavetableMulti::createForHarmonicProfile(
-        sfz::HarmonicProfile::getSquare(), 1.0, 2048);
+        sfz::HarmonicProfile::getSquare(), sfz::config::amplitudeSquare, 2048);
 
     fClient.reset(client);
 
@@ -128,10 +130,21 @@ void DemoApp::initWindow()
     fUi.valWave->addItem(tr("3 - Saw"));
     fUi.valWave->addItem(tr("4 - Square"));
 
+    fUi.valQuality->addItem(tr("1 - Nearest"));
+    fUi.valQuality->addItem(tr("2 - Linear"));
+    fUi.valQuality->addItem(tr("3 - High"));
+    fUi.valQuality->addItem(tr("4 - Dual-High"));
+
+    fUi.valQuality->setCurrentIndex(fOsc.quality());
+
     connect(
         fUi.valWave, QOverload<int>::of(&QComboBox::currentIndexChanged),
         this, [this](int index) { valueChangedWave(index); });
 
+    connect(
+        fUi.valQuality, QOverload<int>::of(&QComboBox::currentIndexChanged),
+        this, [this](int index) { valueChangedQuality(index); });
+
     connect(
         fUi.btnPlaySweep, &QPushButton::clicked,
         this, [this]() { buttonClickedPlaySweep(); });
@@ -152,6 +165,10 @@ int DemoApp::processAudio(jack_nframes_t nframes, void* cbdata)
     if (newWave != -1)
         self->fWavePlaying = newWave;
 
+    int newQuality = self->fNewQualityPending.exchange(-1);
+    if (newQuality != -1)
+        osc.setQuality(newQuality);
+
     osc.setWavetable(&self->fMulti[self->fWavePlaying]);
 
     float* left = reinterpret_cast<float*>(
@@ -183,6 +200,11 @@ void DemoApp::valueChangedWave(int value)
     fNewWavePending.store(value);
 }
 
+void DemoApp::valueChangedQuality(int value)
+{
+    fNewQualityPending.store(value);
+}
+
 void DemoApp::buttonClickedPlaySweep()
 {
     fStartNewSweep.store(true);

tests/DemoWavetables.ui

@@ -6,7 +6,7 @@
    <rect>
     <x>0</x>
     <y>0</y>
-    <width>170</width>
+    <width>255</width>
     <height>103</height>
    </rect>
   </property>
@@ -20,6 +20,13 @@
      </widget>
     </item>
     <item row="0" column="1">
+     <widget class="QLabel" name="label_3">
+      <property name="text">
+       <string>Select quality</string>
+      </property>
+     </widget>
+    </item>
+    <item row="0" column="2">
      <widget class="QLabel" name="label_2">
       <property name="text">
        <string>Play sweep</string>
@@ -30,6 +37,9 @@
      <widget class="QComboBox" name="valWave"/>
     </item>
     <item row="1" column="1">
+     <widget class="QComboBox" name="valQuality"/>
+    </item>
+    <item row="1" column="2">
      <widget class="QPushButton" name="btnPlaySweep">
       <property name="minimumSize">
        <size>
@@ -38,7 +48,8 @@
        </size>
       </property>
       <property name="icon">
-       <iconset theme="media-playback-start"/>
+       <iconset theme="media-playback-start">
+        <normaloff>.</normaloff>.</iconset>
       </property>
      </widget>
     </item>