Karplus Strong pluck generation I want to use a simple implementation of the Karplus-Strong algorithm to generate pluck sounds and chords.

My existing code seems correct but is slow when more than a few plucks are being generated. Is there a way to improve the performance of this. In particular is there a better way to use numpy here, or would something like a collections.deque be faster.

import numpy as np

sample_rate = 44100
damping = 0.999

def generate(f, vol, nsamples):
    """Generate a Karplus-Strong pluck.

    Arguments:
    f -- the frequency, as an integer
    vol -- volume, a float in [0.0, 1.0]
    nsamples -- the number of samples to generate. To generate a pluck t
    seconds long, one needs t * sample_rate samples. nsamples is an int.

    Return value:
    A numpy array of floats with length nsamples.
    """

    N = sample_rate // f
    buf = np.random.rand(N) * 2 - 1
    samples = np.empty(nsamples, dtype=float)

    for i in range(nsamples):
        samples[i] = buf[i % N]
        avg = damping * 0.5 * (buf[i % N] + buf[(1 + i) % N])
        buf[i % N] = avg

    return samples * vol


def generate_chord(f, nsamples):
    """Generate a chord of several plucks

    Arguments:
    f -- a (usually short) list of integer frequencies
    nsamples -- the length of the chord, a chord of length t seconds needs
    t * sample_rate, an integer.

    Return value:
    A numpy array
    """
    samples = np.zeros(nsamples, dtype=float)
    for note in f:
        samples += generate(note, 0.5, nsamples)
    return samples


if __name__ == "__main__":
    import matplotlib.pyplot as plt
    from scipy.io.wavfile import write

    strum = generate_chord(  # A Major
        [220, 220 * 5 // 4, 220 * 6 // 4, 220 * 2], sample_rate * 5)
    plt.plot(strum)
    plt.show()

    write("pluck.wav", sample_rate, strum)


This produces a sound file in pluck.wav with a waveform like this:
Congratulations

You mostly follows standards which is good. You also have good documentation which is great. All in all the code reads well and is easily understandable.

Style-wise, I just have some nitpicks:


from PEP8, you should name your constants using all caps (SAMPLE_RATE, DAMPING);
you should not abreviate that much your variable names frequency and volume are to be prefered to f and vol;
you may want to name things to avoid comments: A_major = [220, 220 * 5 // 4, 220 * 6 // 4, 220 * 2].


Summation

In pure Python world the following code:

value = 0
for element in array:
    value += compute(element)
return value


is prefered written as:

return sum(compute(element) for element in array)


As it is both cleaner and faster. In numpy world, it is pretty much the same. You can write:

def generate_chord(frequencies, samples):
    plucks = (generate(note, 0.5, samples) for note in frequencies)
    return np.sum(plucks, axis=0)


and get the same result.

Interface

When reading your docstrings, I wondered why you would ask your users (or yourself) to pass in a number of samples if what is more natural to them is to use a duration in seconds. Do not expose implementation details such as your sample rate to the users and let them use what is more natural to them.

I also wonder why you hardcode the volume in generate_chord and let it variable in generate. If the aim is to not have to choose it for each test, you can still use a parameter with default value instead of hardcoding it. It will make things more obvious.

Lastly, I don't really understand the need for integral frequencies as the computations can be performed just as well with floating point values. Just make sure to truncate the computation of N.

I would use the following code as a base for improvements:

import numpy as np


SAMPLE_RATE = 44100
DAMPING = 0.999


def generate(frequency, volume, duration):
    """Generate a Karplus-Strong pluck.

    Arguments:
    frequency -- the frequency, as a float
    volume -- volume, a float in [0.0, 1.0]
    duration -- the length of the generated pluck, in seconds.

    Return value:
    A numpy array of floats with length nsamples.
    """

    samples_count = duration * SAMPLE_RATE
    N = int(SAMPLE_RATE / frequency)
    buf = np.random.rand(N) * 2 - 1
    samples = np.empty(samples_count, dtype=float)

    for i in range(samples_count):
        samples[i] = buf[i % N]
        avg = DAMPING * 0.5 * (buf[i % N] + buf[(1 + i) % N])
        buf[i % N] = avg

    return samples * volume


def generate_chord(frequencies, duration, volume=0.5):
    """Generate a chord of several plucks

    Arguments:
    frequencies -- a (usually short) list of frequencies
    describing a chord.
    duration -- the length of the generated pluck, in seconds.

    Return value:
    A numpy array
    """

    plucks = (generate(note, volume, duration) for note in frequencies)
    return np.sum(plucks, axis=0)


if __name__ == "__main__":
    import scipy.io.wavfile
    import matplotlib.pyplot as plt

    A_major = [220, 220 * 5 // 4, 220 * 6 // 4, 220 * 2]
    strum = generate_chord(A_major, duration=5)

    plt.plot(strum)
    plt.show()

    scipy.io.wavfile.write(filename, SAMPLE_RATE, strum)


Beyond frequencies

Now, I’m not musician. And even if I were, I would need to look up frequencies of required note as I don't think your average user would know them by heart. It would be best to provide named constants to remember that for us.

Reading at https://en.wikipedia.org/wiki/Musical_note one can easily define something like:

C = 2 ** (-9/12) * 440
C_sharp = 2 ** (-8/12) * 440
D_flat = 2 ** (-8/12) * 440
D = 2 ** (-7/12) * 440
D_sharp = 2 ** (-6/12) * 440
E_flat = 2 ** (-6/12) * 440
E = 2 ** (-5/12) * 440
F = 2 ** (-4/12) * 440
F_sharp = 2 ** (-3/12) * 440
G_flat = 2 ** (-3/12) * 440
G = 2 ** (-2/12) * 440
G_sharp = 2 ** (-1/12) * 440
A_flat = 2 ** (-1/12) * 440
A = 2 ** (0/12) * 440
A_sharp = 2 ** (1/12) * 440
B_flat = 2 ** (1/12) * 440
B = 2 ** (2/12) * 440


But this does not account for every octave. Let's fix that and simplify the reading using an enum:

import enum


def octave(order):
    frequency_of_A = 440 * 2 ** (order - 4)

    class Octave(enum.Enum):
        C = -9
        C_sharp = -8
        D_flat = -8
        D = -7
        D_sharp = -6
        E_flat = -6
        E = -5
        F = -4
        F_sharp = -3
        G_flat = -3
        G = -2
        G_sharp = -1
        A_flat = -1
        A = 0
        A_sharp = 1
        B_flat = 1
        B = 2

        def __float__(self):
            return 2 ** (self.value/12) * frequency_of_A

    return Octave


SubSubContra = octave(-1)
SubContra = octave(0)
Contra = octave(1)
Great = octave(2)
Small = octave(3)
OneLined = octave(4)
TwoLined = octave(5)
ThreeLined = octave(6)
FourLined = octave(7)
FiveLined = octave(8)
SixLined = octave(9)


Now I can easily define A Major using:

A_major = [Small.A, OneLined.C_sharp, OneLined.E, OneLined.A]


Proposed improvements

import enum
import numpy as np


SAMPLE_RATE = 44100
DAMPING = 0.999


def octave(order):
    """Create an enum of all the notes in the octave
    of the given order.
    """
    frequency_of_A = 440 * 2 ** (order - 4)

    class Octave(enum.Enum):
        C = -9
        C_sharp = -8
        D_flat = -8
        D = -7
        D_sharp = -6
        E_flat = -6
        E = -5
        F = -4
        F_sharp = -3
        G_flat = -3
        G = -2
        G_sharp = -1
        A_flat = -1
        A = 0
        A_sharp = 1
        B_flat = 1
        B = 2

        def __float__(self):
            return 2 ** (self.value/12) * frequency_of_A

    return Octave


SubSubContra = octave(-1)
SubContra = octave(0)
Contra = octave(1)
Great = octave(2)
Small = octave(3)
OneLined = octave(4)
TwoLined = octave(5)
ThreeLined = octave(6)
FourLined = octave(7)
FiveLined = octave(8)
SixLined = octave(9)


def generate(frequency, volume, duration):
    """Generate a Karplus-Strong pluck.

    Arguments:
    frequency -- the frequency, as a float
    volume -- volume, a float in [0.0, 1.0]
    duration -- the length of the generated pluck, in seconds.

    Return value:
    A numpy array of floats with length nsamples.
    """

    samples_count = duration * SAMPLE_RATE
    N = int(SAMPLE_RATE / frequency)
    buf = np.random.rand(N) * 2 - 1
    samples = np.empty(samples_count, dtype=float)

    for i in range(samples_count):
        samples[i] = buf[i % N]
        avg = DAMPING * 0.5 * (buf[i % N] + buf[(1 + i) % N])
        buf[i % N] = avg

    return samples * volume


def generate_chord(frequencies, duration, volume=0.5):
    """Generate a chord of several plucks

    Arguments:
    frequencies -- a (usually short) list of frequencies
    describing a chord.
    duration -- the length of the generated pluck, in seconds.

    Return value:
    A numpy array
    """

    notes = map(float, frequencies)
    plucks = (generate(note, volume, duration) for note in notes)
    return np.sum(plucks, axis=0)


if __name__ == "__main__":
    import scipy.io.wavfile
    import matplotlib.pyplot as plt

    A_major = [Small.A, OneLined.C_sharp, OneLined.E, OneLined.A]
    strum = generate_chord(A_major, duration=5)

    plt.plot(strum)
    plt.show()

    scipy.io.wavfile.write(filename, SAMPLE_RATE, strum)