Standard library · Encoding and data

audio/dsp

A Fourier transform, windows, filters and resampling for sound.

import audio.dsp

On this page

Types
Constants and variables
Functions

This is the arithmetic behind sound: a Fourier transform to see what frequencies are in a block, windows to taper the block first, filters to shape it and resampling to move it between rates. The codecs in this family lean on the same code, which is why it lives here rather than inside any one of them.

Everything works on a clip or on plain arrays of f64. The transform is radix-2, so a length has to be a power of two: use NextPow2 to pad up to one.

Things to know

FFT and IFFT work in place on two arrays, the real parts and the imaginary parts. IFFT already divides by the length, so it undoes FFT exactly.
A real signal shows the same strength at bin k and at bin n - k, so only the first half of a spectrum tells you anything new.
Taper a block with ApplyWindow before transforming it, or a tone that does not fit a whole number of cycles will smear across every bin.
The filters follow the well known cookbook formulas. A Biquad remembers two samples, so BiquadReset it before a different clip, and prefer BiquadProcessAll, which gives each channel its own copy.
BiquadResponse tells you what a filter does at a frequency without running any sound through it, which makes filters easy to test.
ResampleSinc is the one to use on music; ToRate in audio is the quick linear version.

Types

`Spectrum`

struct Spectrum

A spectrogram: one magnitude row per analysis frame.

Name	Type	Description
`bins`	`int = 0`	How many frequency bins each row holds, which is half the window.
`frames`	`int = 0`	How many rows there are.
`rate`	`int = 0`	The sample rate the clip had, needed to turn a bin into a frequency.
`hop`	`int = 0`	How many samples apart the rows are.
`data`	`void*`	The magnitudes, row after row. Release them with FreeSpectrum.

dsp-Analyze.salam

import audio
import audio.dsp
import audio.synth
import math

func main:
    low := synth.Tone(110.0, 0.1, 22050, synth.WaveSine)
    high := synth.Tone(1760.0, 0.1, 22050, synth.WaveSine)
    println dsp.ZeroCrossingRate(high) > dsp.ZeroCrossingRate(low) * 4.0
    println dsp.ZeroCrossingRate(audio.Clip {})

    each f in [220.0, 440.0]:
        t := synth.Tone(f, 0.1, 22050, synth.WaveSine)
        got := dsp.DetectPitch(t, 0, 2048)
        println math.Abs(got - f) / f < 0.05
        mut d := t
        audio.Free(d)
    end
    println dsp.DetectPitch(audio.Clip {}, 0, 2048)

    rich := synth.Additive(220.0, 8, 0.7, 0.12, 22050)
    mut s := dsp.STFT(rich, 1024, 256, dsp.WindowHann)
    println s.bins, s.rate, s.hop, s.frames > 1
    println dsp.SpectrumAt(s, 0, 9999)
    cen := dsp.SpectralCentroid(s, 1)
    println cen > 100.0, cen < 5000.0
    println dsp.SpectralRolloff(s, 1, 0.85) > 0.0
    flat := dsp.SpectralFlatness(s, 1)
    println flat >= 0.0, flat <= 1.0
    println dsp.SpectralFlux(s, 0)

    noise := synth.NoiseWhite(0.12, 22050)
    mut ns := dsp.STFT(noise, 1024, 256, dsp.WindowHann)
    println dsp.SpectralFlatness(ns, 1) > flat

    loud := synth.Tone(440.0, 0.1, 22050, synth.WaveSine)
    lu := dsp.LoudnessLUFS(loud)
    println lu < 0.0, lu > -40.0
    mut soft := audio.Clone(loud)
    audio.Gain(soft, -20.0)
    println math.Abs((dsp.LoudnessLUFS(soft) - lu) + 20.0) < 1.0
    mut norm := audio.Clone(loud)
    dsp.NormalizeLUFS(norm, -16.0)
    println math.Abs(dsp.LoudnessLUFS(norm) + 16.0) < 1.0
    println dsp.TruePeakDb(loud) >= audio.PeakDb(loud) - 0.01

    dsp.FreeSpectrum(s)
    dsp.FreeSpectrum(ns)
    audio.Free(low)
    audio.Free(high)
    audio.Free(rich)
    audio.Free(noise)
    audio.Free(loud)
    audio.Free(soft)
    audio.Free(norm)
end

Outputtrue 0 true true 0 512 22050 256 true 0 true true true true true 0 true true true true true true

`Complex`

struct Complex

A complex number, used by the Fourier transform.

Name	Type	Description
`re`	`f64 = 0`	The real part.
`im`	`f64 = 0`	The imaginary part.

dsp-Complex.salam

import audio.dsp
import math

func main:
    a := dsp.Complex { re = 3.0, im = 4.0 }
    b := dsp.Complex { re = 1.0, im = -2.0 }
    println a.re, a.im
    println dsp.CAbs(a)
    println dsp.CAdd(a, b).re, dsp.CAdd(a, b).im
    println dsp.CSub(a, b).re, dsp.CSub(a, b).im
    println dsp.CMul(a, b).re, dsp.CMul(a, b).im
    println dsp.CScale(a, 2.0).re, dsp.CScale(a, 2.0).im
    println dsp.CConj(a).re, dsp.CConj(a).im
    up := dsp.Complex { re = 0.0, im = 1.0 }
    println math.Abs(dsp.CArg(up) - math.PI / 2.0) < 0.000000001
    p := dsp.CFromPolar(2.0, 0.0)
    println math.Abs(p.re - 2.0) < 0.000000001, math.Abs(p.im) < 0.000000001
    println dsp.CAbs(dsp.CFromPolar(5.0, 1.0))
end

Output3 4 5 4 2 2 6 11 -2 6 8 3 -4 true true true 5

`Compressor`

struct Compressor

How hard to squeeze a signal and how fast to react.

Name	Type	Description
`threshold`	`f64 = -20`	The level in decibels above which the squeezing starts.
`ratio`	`f64 = 4`	How much of the excess to keep. A ratio of 4 lets a quarter of it through.
`attack`	`f64 = 0.005`	How many seconds to take hold once the signal goes over.
`release`	`f64 = 0.1`	How many seconds to let go once it comes back under.
`makeup`	`f64 = 0`	Decibels of gain added afterwards, to win back what the squeezing cost.

dsp-Effect.salam

import audio
import audio.dsp
import audio.synth
import math

func main:
    src := synth.Tone(440.0, 0.1, 22050, synth.WaveSine)
    println src.frames, src.rate

    e := dsp.Echo(src, 0.02, 2, 0.5)
    println e.frames, e.frames == src.frames + 2 * 441
    println math.Abs(audio.Sample(e, 451, 0) - (audio.Sample(src, 451, 0) + audio.Sample(src, 10, 0) * 0.5)) < 0.001

    rv := dsp.Reverb(src, 0.5, 0.3)
    println rv.frames, audio.Peak(rv) <= 1.0
    dry := dsp.Reverb(src, 0.5, 0.0)
    println math.Abs(audio.Sample(dry, 50, 0) - audio.Sample(src, 50, 0)) < 0.0001

    ch := dsp.Chorus(src, 0.5, 2.0, 0.5)
    fl := dsp.Flanger(src, 0.5, 1.0, 0.5)
    println ch.frames, fl.frames
    println dsp.Correlation(src, ch) < 0.999, dsp.Correlation(src, fl) < 0.9999

    mut tr := audio.Clone(src)
    dsp.Tremolo(tr, 10.0, 1.0)
    println audio.RMS(tr) < audio.RMS(src), audio.Peak(tr) <= audio.Peak(src) + 0.0001

    mut dl := audio.Clone(src)
    dsp.Delay(dl, 0.01, 0.5, 0.5)
    println dl.frames == src.frames

    k := dsp.Compressor { threshold = -20.0, ratio = 8.0, attack = 0.001, release = 0.05 }
    mut hard := audio.Clone(src)
    dsp.CompressorProcess(k, hard)
    println audio.RMSDb(hard) < audio.RMSDb(src) - 6.0

    mut hot := audio.Clone(src)
    audio.Scale(hot, 4.0)
    dsp.LimiterProcess(hot, -1.0)
    println audio.PeakDb(hot) <= -1.0 + 0.001

    mut gated := audio.Clone(src)
    audio.Gain(gated, -60.0)
    dsp.GateProcess(gated, -40.0, 0.001)
    println audio.Peak(gated) < 0.0001

    longer := dsp.TimeStretch(src, 2.0)
    println longer.frames == src.frames * 2, longer.rate
    up := dsp.PitchShift(src, 12.0)
    println up.frames, up.rate
    println dsp.DetectPitch(up, 0, 2048) > 700.0

    audio.Free(src)
    audio.Free(e)
    audio.Free(rv)
    audio.Free(dry)
    audio.Free(ch)
    audio.Free(fl)
    audio.Free(tr)
    audio.Free(dl)
    audio.Free(hard)
    audio.Free(hot)
    audio.Free(gated)
    audio.Free(longer)
    audio.Free(up)
end

Output2205 22050 3087 true true 2205 true true 2205 2205 true true true true true true true true true 22050 2205 22050 true

`Biquad`

struct Biquad

A two pole, two zero filter: five coefficients and the two samples it remembers.

Name	Type	Description
`b0`	`f64 = 1`	The first feed forward coefficient.
`b1`	`f64 = 0`	The second feed forward coefficient.
`b2`	`f64 = 0`	The third feed forward coefficient.
`a1`	`f64 = 0`	The first feedback coefficient, already divided through by a0.
`a2`	`f64 = 0`	The second feedback coefficient, already divided through by a0.
`z1`	`f64 = 0`	The first remembered sample. Reset it before filtering a new clip.
`z2`	`f64 = 0`	The second remembered sample.

dsp-Filter.salam

import audio
import audio.dsp
import math

func tone(rate: int, freq: f64, frames: int): audio.Clip:
    mut c := audio.New(rate, 1, frames)
    mut i := 0
    until i < frames:
        audio.SetSample(c, i, 0, math.Sin(2.0 * math.PI * freq * (i as f64) / (rate as f64)))
        i += 1
    end
    ret c
end

func main:
    lp := dsp.LowPass(1000.0, 0.707, 48000)
    println dsp.BiquadResponse(lp, 100.0, 48000) > 0.99
    println math.Abs(dsp.BiquadResponse(lp, 1000.0, 48000) - 0.7071) < 0.01
    println dsp.BiquadResponse(lp, 10000.0, 48000) < 0.02
    hp := dsp.HighPass(1000.0, 0.707, 48000)
    println dsp.BiquadResponse(hp, 10000.0, 48000) > 0.99, dsp.BiquadResponse(hp, 100.0, 48000) < 0.02
    bp := dsp.BandPass(1000.0, 1.0, 48000)
    println math.Abs(dsp.BiquadResponse(bp, 1000.0, 48000) - 1.0) < 0.01, dsp.BiquadResponse(bp, 100.0, 48000) < 0.2
    nt := dsp.Notch(1000.0, 1.0, 48000)
    println dsp.BiquadResponse(nt, 1000.0, 48000) < 0.01, dsp.BiquadResponse(nt, 100.0, 48000) > 0.9
    ap := dsp.AllPass(1000.0, 0.707, 48000)
    println math.Abs(dsp.BiquadResponse(ap, 500.0, 48000) - 1.0) < 0.000001
    pk := dsp.Peaking(1000.0, 1.0, 6.0, 48000)
    println dsp.BiquadResponse(pk, 1000.0, 48000) > 1.98, dsp.BiquadResponse(pk, 1000.0, 48000) < 2.01
    ls := dsp.LowShelf(1000.0, 0.707, 6.0, 48000)
    println dsp.BiquadResponse(ls, 50.0, 48000) > 1.9, math.Abs(dsp.BiquadResponse(ls, 20000.0, 48000) - 1.0) < 0.05
    hs := dsp.HighShelf(1000.0, 0.707, 6.0, 48000)
    println dsp.BiquadResponse(hs, 20000.0, 48000) > 1.9

    mut hi := tone(48000, 10000.0, 2048)
    before := audio.RMS(hi)
    dsp.BiquadProcessAll(lp, hi)
    println audio.GainToDb(audio.RMS(hi) / before) < -20.0
    mut lo := tone(48000, 100.0, 2048)
    lo_before := audio.RMS(lo)
    dsp.BiquadProcessAll(lp, lo)
    println audio.GainToDb(audio.RMS(lo) / lo_before) > -1.0

    mut one := dsp.LowPass(1000.0, 0.707, 48000)
    dsp.BiquadStep(one, 1.0)
    println one.z1 != 0.0
    dsp.BiquadReset(one)
    println one.z1, one.z2
    println math.Abs(one.b0 - lp.b0) < 0.000000001

    mut ch := tone(8000, 100.0, 400)
    audio.Scale(ch, 0.5)
    mut k := 0
    until k < ch.frames:
        audio.AddSample(ch, k, 0, 0.25)
        k += 1
    end
    println audio.DcOffset(ch) > 0.2
    dsp.DCRemove(ch)
    println math.Abs(audio.DcOffset(ch)) < 0.000001
    mut single := dsp.LowPass(1000.0, 0.707, 8000)
    dsp.BiquadProcess(single, ch, 0)
    println ch.frames
    audio.Free(hi)
    audio.Free(lo)
    audio.Free(ch)
end

Outputtrue true true true true true true true true true true true true true true true true true 0 0 true true true 400

Constants and variables

Resampling quality `Quality…`

How carefully to resample. Linear is quickest, sinc is best on music.

Name	Value
`QualityLinear`	`0`
`QualityCubic`	`1`
`QualitySinc`	`2`

dsp-Resample.salam

import audio
import audio.dsp
import math

func tone(rate: int, freq: f64, frames: int): audio.Clip:
    mut c := audio.New(rate, 1, frames)
    mut i := 0
    until i < frames:
        audio.SetSample(c, i, 0, math.Sin(2.0 * math.PI * freq * (i as f64) / (rate as f64)))
        i += 1
    end
    ret c
end

func main:
    src := tone(44100, 440.0, 4410)
    up := dsp.ResampleSinc(src, 48000, 16)
    println up.rate, up.frames
    back := dsp.ResampleSinc(up, 44100, 16)
    println back.rate, back.frames
    inner := audio.Slice(src, 64, src.frames - 128)
    inner2 := audio.Slice(back, 64, back.frames - 128)
    println dsp.Correlation(inner, inner2) > 0.999
    cub := dsp.ResampleCubic(src, 48000)
    cback := dsp.ResampleCubic(cub, 44100)
    inner3 := audio.Slice(cback, 64, cback.frames - 128)
    println cub.rate, dsp.Correlation(inner, inner3) > 0.999
    println dsp.ResampleTo(src, 48000, dsp.QualityLinear).rate
    println dsp.ResampleTo(src, 48000, dsp.QualityCubic).rate
    println dsp.ResampleTo(src, 48000, dsp.QualitySinc).rate
    same := dsp.ResampleSinc(src, 44100, 16)
    println math.Abs(dsp.Correlation(src, same) - 1.0) < 0.000001
    println dsp.Correlation(src, audio.Clip {})
    println dsp.KaiserBeta
    down := dsp.ResampleSinc(src, 22050, 16)
    println down.rate, down.frames
    audio.Free(src)
    audio.Free(up)
    audio.Free(back)
    audio.Free(inner)
    audio.Free(inner2)
    audio.Free(cub)
    audio.Free(cback)
    audio.Free(inner3)
    audio.Free(same)
    audio.Free(down)
end

Output48000 4800 44100 4410 true 48000 true 48000 48000 48000 true 0 8.6 22050 2205

`KaiserBeta`

const KaiserBeta: f64 = 8.6

The beta the resampler uses for its Kaiser window.

The example for this item is shown with the ones above: dsp-Resample.salam

Window shapes `Window…`

Which window to taper a block with before a transform. Rectangular is no window at all.

Name	Value
`WindowRectangular`	`0`
`WindowHann`	`1`
`WindowHamming`	`2`
`WindowBlackman`	`3`
`WindowBlackmanHarris`	`4`
`WindowSine`	`5`
`WindowVorbis`	`6`
`WindowKaiser`	`7`

dsp-Window.salam

import audio.dsp
import math

func main:
    println dsp.BesselI0(0.0)
    hann := dsp.MakeWindow(dsp.WindowHann, 9, 0.0)
    println hann.len(), hann.get(0) < 0.000000001, math.Abs(hann.get(4) - 1.0) < 0.000000001
    println hann.get(8) < 0.000000001
    println dsp.WindowValue(dsp.WindowRectangular, 0, 9, 0.0)
    println math.Abs(dsp.WindowValue(dsp.WindowHamming, 0, 9, 0.0) - 0.08) < 0.000000001
    println math.Abs(dsp.MakeWindow(dsp.WindowBlackman, 9, 0.0).get(4) - 1.0) < 0.000000001
    println math.Abs(dsp.MakeWindow(dsp.WindowBlackmanHarris, 9, 0.0).get(4) - 1.0) < 0.000000001
    println math.Abs(dsp.MakeWindow(dsp.WindowSine, 9, 0.0).get(4) - 1.0) < 0.000000001
    println dsp.MakeWindow(dsp.WindowVorbis, 9, 0.0).get(0) < 0.000000001
    println math.Abs(dsp.MakeWindow(dsp.WindowKaiser, 9, 8.6).get(4) - 1.0) < 0.000000001
    println dsp.WindowName(dsp.WindowHann), dsp.WindowName(dsp.WindowKaiser), dsp.WindowName(99)
    println dsp.WindowSum(dsp.WindowRectangular, 8, 0.0)
    mut buf := Vector {} as Vector<f64>
    mut i := 0
    until i < 8:
        buf.push(1.0)
        i += 1
    end
    dsp.ApplyWindow(buf.ref(0), 8, dsp.WindowHann, 0.0)
    println buf.get(0) < 0.000000001, buf.get(4) > 0.9
end

Output1 9 true true true 1 true true true true true true hann kaiser unknown 8 true true