Phase modulation synth in LuaJIT and C

portaudio.lua

local ffi = require 'ffi'

local pa = ffi.load('portaudio')

ffi.cdef [[
int Pa_GetVersion( void );
const char* Pa_GetVersionText( void );
typedef int PaError;
typedef enum PaErrorCode
{
    paNoError = 0,

    paNotInitialized = -10000,
    paUnanticipatedHostError,
    paInvalidChannelCount,
    paInvalidSampleRate,
    paInvalidDevice,
    paInvalidFlag,
    paSampleFormatNotSupported,
    paBadIODeviceCombination,
    paInsufficientMemory,
    paBufferTooBig,
    paBufferTooSmall,
    paNullCallback,
    paBadStreamPtr,
    paTimedOut,
    paInternalError,
    paDeviceUnavailable,
    paIncompatibleHostApiSpecificStreamInfo,
    paStreamIsStopped,
    paStreamIsNotStopped,
    paInputOverflowed,
    paOutputUnderflowed,
    paHostApiNotFound,
    paInvalidHostApi,
    paCanNotReadFromACallbackStream,
    paCanNotWriteToACallbackStream,
    paCanNotReadFromAnOutputOnlyStream,
    paCanNotWriteToAnInputOnlyStream,
    paIncompatibleStreamHostApi,
    paBadBufferPtr
} PaErrorCode;
const char *Pa_GetErrorText( PaError errorCode );
PaError Pa_Initialize( void );
PaError Pa_Terminate( void );

typedef int PaDeviceIndex;
enum
{
    paNoDevice=-1,
    paUseHostApiSpecificDeviceSpecification=-2
};
typedef int PaHostApiIndex;
PaHostApiIndex Pa_GetHostApiCount( void );
PaHostApiIndex Pa_GetDefaultHostApi( void );
typedef enum PaHostApiTypeId
{
    paInDevelopment=0, /* use while developing support for a new host API */
    paDirectSound=1,
    paMME=2,
    paASIO=3,
    paSoundManager=4,
    paCoreAudio=5,
    paOSS=7,
    paALSA=8,
    paAL=9,
    paBeOS=10,
    paWDMKS=11,
    paJACK=12,
    paWASAPI=13,
    paAudioScienceHPI=14
} PaHostApiTypeId;
typedef struct PaHostApiInfo
{
    int structVersion;
    PaHostApiTypeId type;
    const char *name;
    int deviceCount;
    PaDeviceIndex defaultInputDevice;
    PaDeviceIndex defaultOutputDevice;
} PaHostApiInfo;
const PaHostApiInfo * Pa_GetHostApiInfo( PaHostApiIndex hostApi );
PaHostApiIndex Pa_HostApiTypeIdToHostApiIndex( PaHostApiTypeId type );
PaDeviceIndex Pa_HostApiDeviceIndexToDeviceIndex( PaHostApiIndex hostApi,
        int hostApiDeviceIndex );
typedef struct PaHostErrorInfo{
    PaHostApiTypeId hostApiType;    /**< the host API which returned the error code */
    long errorCode;                 /**< the error code returned */
    const char *errorText;          /**< a textual description of the error if available, otherwise a zero-length string */
}PaHostErrorInfo;
const PaHostErrorInfo* Pa_GetLastHostErrorInfo( void );
PaDeviceIndex Pa_GetDeviceCount( void );
PaDeviceIndex Pa_GetDefaultInputDevice( void );
PaDeviceIndex Pa_GetDefaultOutputDevice( void );
typedef double PaTime;
typedef unsigned long PaSampleFormat;

enum
{
    paFloat32        = 0x00000001,
    paInt32          = 0x00000002,
    paInt24          = 0x00000004,
    paInt16          = 0x00000008,
    paInt8           = 0x00000010,
    paUInt8          = 0x00000020,
    paCustomFormat   = 0x00010000,

    paNonInterleaved = 0x80000000
};

typedef struct PaDeviceInfo
{
    int structVersion;
    const char *name;
    PaHostApiIndex hostApi;
    
    int maxInputChannels;
    int maxOutputChannels;

    PaTime defaultLowInputLatency;
    PaTime defaultLowOutputLatency;
    PaTime defaultHighInputLatency;
    PaTime defaultHighOutputLatency;

    double defaultSampleRate;
} PaDeviceInfo;

const PaDeviceInfo* Pa_GetDeviceInfo( PaDeviceIndex device );

typedef struct PaStreamParameters
{
    PaDeviceIndex device;
    int channelCount;
    PaSampleFormat sampleFormat;
    PaTime suggestedLatency;
    void *hostApiSpecificStreamInfo;
} PaStreamParameters;
enum
{
    paFormatIsSupported=0
};
PaError Pa_IsFormatSupported( const PaStreamParameters *inputParameters,
                              const PaStreamParameters *outputParameters,
                              double sampleRate );
typedef void PaStream;
enum
{
    paFramesPerBufferUnspecified=0
};
typedef unsigned long PaStreamFlags;

enum
{
    paNoFlag          = 0,
    paClipOff         = 0x00000001,
    paDitherOff       = 0x00000002,
    paNeverDropInput  = 0x00000004,
    paPrimeOutputBuffersUsingStreamCallback = 0x00000008,
    paPlatformSpecificFlags = 0xFFFF0000
};

typedef struct PaStreamCallbackTimeInfo {
    PaTime inputBufferAdcTime;
    PaTime currentTime;
    PaTime outputBufferDacTime;
} PaStreamCallbackTimeInfo;


typedef unsigned long PaStreamCallbackFlags;
enum
{
    paInputUnderflow  = 0x00000001,
    paInputOverflow   = 0x00000002,
    paOutputUnderflow = 0x00000004,
    paOutputOverflow  = 0x00000008,
    paPrimingOutput   = 0x00000010
};
enum PaStreamCallbackResult
{
    paContinue=0,
    paComplete=1,
    paAbort=2
} PaStreamCallbackResult;
typedef int PaStreamCallback(
    const void *input, void *output,
    unsigned long frameCount,
    const PaStreamCallbackTimeInfo* timeInfo,
    PaStreamCallbackFlags statusFlags,
    void *userData );
PaError Pa_OpenStream( PaStream** stream,
                       const PaStreamParameters *inputParameters,
                       const PaStreamParameters *outputParameters,
                       double sampleRate,
                       unsigned long framesPerBuffer,
                       PaStreamFlags streamFlags,
                       PaStreamCallback *streamCallback,
                       void *userData );
PaError Pa_OpenDefaultStream( PaStream** stream,
                              int numInputChannels,
                              int numOutputChannels,
                              PaSampleFormat sampleFormat,
                              double sampleRate,
                              unsigned long framesPerBuffer,
                              PaStreamCallback *streamCallback,
                              void *userData );
PaError Pa_CloseStream( PaStream *stream );
typedef void PaStreamFinishedCallback( void *userData );
PaError Pa_SetStreamFinishedCallback( PaStream *stream, PaStreamFinishedCallback* streamFinishedCallback ); 
PaError Pa_StartStream( PaStream *stream );
PaError Pa_StopStream( PaStream *stream );
PaError Pa_AbortStream( PaStream *stream );
PaError Pa_IsStreamStopped( PaStream *stream );
PaError Pa_IsStreamActive( PaStream *stream );
typedef struct PaStreamInfo
{
    int structVersion;
    PaTime inputLatency;
    PaTime outputLatency;
    double sampleRate;
} PaStreamInfo;
const PaStreamInfo* Pa_GetStreamInfo( PaStream *stream );
PaTime Pa_GetStreamTime( PaStream *stream );
double Pa_GetStreamCpuLoad( PaStream* stream );
PaError Pa_ReadStream( PaStream* stream,
                       void *buffer,
                       unsigned long frames );
PaError Pa_WriteStream( PaStream* stream,
                        const void *buffer,
                        unsigned long frames );
signed long Pa_GetStreamReadAvailable( PaStream* stream );
signed long Pa_GetStreamWriteAvailable( PaStream* stream );
PaError Pa_GetSampleSize( PaSampleFormat format );
void Pa_Sleep( long msec );
]]

local buffer_size = 1024

local function check_error(err)
  if err == pa.paUnanticipatedHostError then
    assert(false, ffi.string(pa.Pa_GetLastHostErrorInfo().errorText))
  elseif err ~= pa.paNoError then
    assert(false, ffi.string(pa.Pa_GetErrorText(err)))
  end
end

local stream_ptr = ffi.new('PaStream*[1]')

local function init()
  check_error(pa.Pa_Initialize())

  local outputParams = ffi.new('struct PaStreamParameters')
  outputParams.device = pa.Pa_GetDefaultOutputDevice()
  outputParams.channelCount = 1
  outputParams.sampleFormat = pa.paFloat32
  outputParams.suggestedLatency =
    pa.Pa_GetDeviceInfo(outputParams.device).defaultHighOutputLatency
  outputParams.hostApiSpecificStreamInfo = nil

  check_error(pa.Pa_OpenStream(
    stream_ptr, nil, outputParams, 44100, buffer_size, 0, nil, nil))

  check_error(pa.Pa_StartStream(stream_ptr[0]))
end

local buffer = ffi.new('float[?]', buffer_size)
local index = 0

local function put_sample(sample)
  buffer[index] = sample
  index = index + 1
  if index == buffer_size then
    index = 0
    --check_error(pa.Pa_WriteStream(stream_ptr[0], buffer, buffer_size))
    pa.Pa_WriteStream(stream_ptr[0], buffer, buffer_size)
  end
end

local function put_buffer(b)
  for i = 1, 1024 do
    buffer[i-1] = b[i]
  end
  pa.Pa_WriteStream(stream_ptr[0], buffer, 1024)
end

local function uninit()
  check_error(pa.Pa_StopStream(stream_ptr[0]))
  check_error(pa.Pa_CloseStream(stream_ptr[0]))
  pa.Pa_Terminate()
end

return {init = init, put_buffer = put_buffer, put_sample = put_sample, uninit = uninit}

synth.c

// compiled with -O3 -lportaudio -lm -std=c99
#include <portaudio.h>
#include <stdlib.h>
#include <stdio.h>
#include <math.h>

typedef double number_t;

#define SAMPLE_RATE 44100

#define TABLE_MASK 0x3ff
#define TABLE_SIZE (TABLE_MASK + 1)

number_t sine_table[TABLE_SIZE];
void sine_table_init()
{
    for(int i = 0; i < TABLE_SIZE; i++)
        sine_table[i] = sin((number_t)i/TABLE_SIZE * 2 * 3.14159);
}

number_t sine_wave(number_t phase)
{
    return sine_table[(int)(phase * TABLE_SIZE) & TABLE_MASK];
}

typedef struct {
    number_t frequency;
    number_t phase;
} operator_s;

int operator_init(operator_s *operator)
{
    operator->frequency = 440;
    operator->phase = 0;
}

number_t operator_process_sample(operator_s *operator, number_t offset)
{
    operator->phase += operator->frequency/SAMPLE_RATE;
    while(operator->phase >= 1) operator->phase -= 1;
    return sine_wave(operator->phase + offset);
}

typedef struct {
    operator_s operators[6];
    number_t buffer[6];
    number_t back_buffer[6];
    number_t modulation[7*6]; // [a+7*b] is modulation of b on a
} synth_s;

void synth_init(synth_s *synth)
{
    for(int i = 0; i < 6; i++)
    {
        operator_init(&synth->operators[i]);
        synth->buffer[i] = 0;
        synth->back_buffer[i] = 0;
    }
    for(int i = 0; i < 7*6; i++)
        synth->modulation[i] = 0;
}

number_t synth_process_sample(synth_s *synth)
{
    for(int i = 0; i < 6; i++)
    {
        synth->back_buffer[i] = operator_process_sample(&synth->operators[i],
            synth->modulation[i+7*0] * synth->buffer[0] +
            synth->modulation[i+7*1] * synth->buffer[1] +
            synth->modulation[i+7*2] * synth->buffer[2] +
            synth->modulation[i+7*3] * synth->buffer[3] +
            synth->modulation[i+7*4] * synth->buffer[4] +
            synth->modulation[i+7*5] * synth->buffer[5]);
    }
    for(int i = 0; i < 6; i++)
    {
        synth->buffer[i] = synth->back_buffer[i];
    }
    return
        synth->modulation[6+7*0] * synth->buffer[0] +
        synth->modulation[6+7*1] * synth->buffer[1] +
        synth->modulation[6+7*2] * synth->buffer[2] +
        synth->modulation[6+7*3] * synth->buffer[3] +
        synth->modulation[6+7*4] * synth->buffer[4] +
        synth->modulation[6+7*5] * synth->buffer[5];
}

// portaudio

#define PORTAUDIO_BUFFER_SIZE 1024
int portaudio_index = 0;
float portaudio_buffer[PORTAUDIO_BUFFER_SIZE];

PaStream *portaudio_stream;

void portaudio_check_error(int err)
{
    if(err == paUnanticipatedHostError)
    {
        puts(Pa_GetLastHostErrorInfo()->errorText);
        abort();
    }
    else if(err != paNoError)
    {
        puts(Pa_GetErrorText(err));
        abort();
    }
}

void portaudio_init()
{
    portaudio_check_error(Pa_Initialize());
    struct PaStreamParameters outputParams;
    outputParams.device = Pa_GetDefaultOutputDevice();
    outputParams.channelCount = 1;
    outputParams.sampleFormat = paFloat32;
    outputParams.suggestedLatency =
        Pa_GetDeviceInfo(outputParams.device)->defaultHighOutputLatency;
    outputParams.hostApiSpecificStreamInfo = 0;
    
    portaudio_check_error(Pa_OpenStream(
        &portaudio_stream, 0, &outputParams, SAMPLE_RATE, PORTAUDIO_BUFFER_SIZE, 0, 0, 0));
    portaudio_check_error(Pa_StartStream(portaudio_stream));
}

void portaudio_put_sample(number_t sample)
{
    portaudio_buffer[portaudio_index] = sample;
    portaudio_index += 1;
    if(portaudio_index == PORTAUDIO_BUFFER_SIZE)
    {
        portaudio_index = 0;
        Pa_WriteStream(portaudio_stream, portaudio_buffer, PORTAUDIO_BUFFER_SIZE);
    }
}

void portaudio_uninit()
{
    portaudio_check_error(Pa_StopStream(portaudio_stream));
    portaudio_check_error(Pa_CloseStream(portaudio_stream));
    Pa_Terminate();
}

// test

#define VOICE_COUNT 20

int main()
{
    sine_table_init();

    synth_s voices[VOICE_COUNT];
    for(int v = 0; v < VOICE_COUNT; v++)
    {
        synth_s *s = &voices[v];
        synth_init(s);

        number_t base_frequency = 110;
        s->operators[0].frequency = base_frequency*3.5;
        s->operators[1].frequency = base_frequency*0.5001;
        s->operators[2].frequency = base_frequency*4.5;
        s->operators[3].frequency = base_frequency*0.501;
        s->operators[4].frequency = base_frequency*1.011;
        s->operators[5].frequency = base_frequency*1;

        s->modulation[0+7*0] = 0.08;
        s->modulation[1+7*0] = 0.16;
        s->modulation[1+7*1] = 0.04;

        s->modulation[2+7*2] = 0.11;
        s->modulation[3+7*2] = 0.19;
        s->modulation[3+7*3] = 0.05;

        s->modulation[4+7*4] = 0.20;
        s->modulation[5+7*4] = 0.10;
        s->modulation[5+7*5] = 0.04;

        s->modulation[6+7*1] = 0.30;
        s->modulation[6+7*3] = 0.31;
        s->modulation[6+7*5] = 0.30;
    }

    number_t buffer[1024];
    //portaudio_init();
    for(int x = 0; x < SAMPLE_RATE * 180 / 1024; x++)
    {
        for(int i = 0; i < 1024; i++)
        {
            buffer[i] = 0;
        }
        for(int v = 0; v < VOICE_COUNT; v++)
        {
            for(int i = 0; i < 1024; i++)
            {
                buffer[i] += synth_process_sample(&voices[v]);
            }
        }
        //for(int i = 0; i < 1024; i++)
        //{
        //    portaudio_put_sample(buffer[i]/VOICE_COUNT/4);
        //}
    }
    //portaudio_uninit();
}

synth.lua

local sample_rate = 44100

---- Sine Wave

local table_mask = 0x3ff
local table_size = table_mask + 1;
local sine_table = {}
for i = 0, table_size-1 do
  sine_table[i] = math.sin(i/table_size * 2 * math.pi)
end

local function sine_wave(phase)
  return sine_table[bit.band(phase * table_size, table_mask)]
end

---- Class stub

local function class()
  local c = {}
  c.__index = c
  return setmetatable(c, {
    __call = function (...)
      local obj = setmetatable({}, c)
      if obj._init then
        obj:_init(...)
      end
      return obj
    end
  })
end

---- Phase modulation operator

local operator = class()
operator._name = 'operator'

function operator:_init()
  self.frequency = 440
  self.phase = 0
end

function operator:process_sample(offset)
  self.phase = (self.phase + self.frequency/sample_rate) % 1
  return sine_wave(self.phase + offset)
end

---- 6-operator phase modulation synth

local synth = class()
synth.name = 'synth'

function synth:_init()
  self.operators = {}
  self.buffer = {}
  self.back_buffer = {}
  self.modulation = {}

  for i = 1, 6 do
    self.operators[i] = operator()
    self.buffer[i] = 0
    self.back_buffer[i] = 0
  end

  for i = 1, 7 do
    self.modulation[i] = {}
    for j = 1, 6 do
      self.modulation[i][j] = 0
    end
  end
end

function synth:process_buffer(output, first, last)
  for i = first, last do
    output[i] = output[i] + self:process_sample()
  end
end

function synth:process_sample()
  -- run the operators
  for i = 1, 6 do
    self.back_buffer[i] = self.operators[i]:process_sample(
      self.modulation[i][1] * self.buffer[1] +
      self.modulation[i][2] * self.buffer[2] +
      self.modulation[i][3] * self.buffer[3] +
      self.modulation[i][4] * self.buffer[4] +
      self.modulation[i][5] * self.buffer[5] +
      self.modulation[i][6] * self.buffer[6])
  end

  for i = 1, 6 do
    self.buffer[i] = self.back_buffer[i]
  end

  -- mix final output
  return
    self.modulation[7][1] * self.buffer[1] +
    self.modulation[7][2] * self.buffer[2] +
    self.modulation[7][3] * self.buffer[3] +
    self.modulation[7][4] * self.buffer[4] +
    self.modulation[7][5] * self.buffer[5] +
    self.modulation[7][6] * self.buffer[6]
end

---- test

local voices = {}
for i = 1, 20 do
  local s = synth()

  -- cool distorted bass-type preset
  local base_frequency = 110
  s.operators[1].frequency = base_frequency*3.5
  s.operators[2].frequency = base_frequency*0.5001
  s.operators[3].frequency = base_frequency*4.5
  s.operators[4].frequency = base_frequency*0.501
  s.operators[5].frequency = base_frequency*1.011
  s.operators[6].frequency = base_frequency*1

  s.modulation[1][1] = 0.08
  s.modulation[2][1] = 0.16
  s.modulation[2][2] = 0.04

  s.modulation[3][3] = 0.11
  s.modulation[4][3] = 0.19
  s.modulation[4][4] = 0.05

  s.modulation[5][5] = 0.20
  s.modulation[6][5] = 0.10
  s.modulation[6][6] = 0.04

  s.modulation[7][2] = 0.30
  s.modulation[7][4] = 0.31
  s.modulation[7][6] = 0.30

  table.insert(voices, s)
end

local play = (...) == '--play'

local portaudio = play and require 'portaudio'

if play then
  portaudio.init()
end
local buffer = {}
for x = 1, sample_rate * 180 / 1024 do
  for i = 1, 1024 do
    buffer[i] = 0
  end
  for i = 1, #voices do
    voices[i]:process_buffer(buffer, 1, #buffer)
  end
  if play then
    for i = 1, #buffer do
      portaudio.put_sample(buffer[i]/#voices/4)
    end
  end
end
if play then
  portaudio.uninit()
end