Quarter sine wave look up table synth with various oscillators in Nim

_lutsynth.nim

import math

const
  SampleRate {.intdefine.} = 44100
  SRate* = SampleRate.float
  LutSize = 256  # Quarter sine wave
  Lut = static:
    var arr: array[LutSize, float]
    for i in 0..<LutSize:
      let angle = (i.float / LutSize.float) * (PI / 2.0)
      arr[i] = sin(angle)
    arr

type
  Ticker* = object
    tick: uint

  WaveGenerator* = object
    #lut: seq[float]
    phase: float
    freq: float
    sampleRate: float


proc initWaveGenerator*(sampleRate: float = SRate): WaveGenerator =
  result.sampleRate = sampleRate
  result.phase = 0.0
  result.freq = 440.0  
  # Create quarter sine LUT, 0 to π/2 
  #result.lut = newSeq[float](LutSize)
  #for i in 0..<LutSize:
  #  let angle = (i.float / LutSize.float) * (PI / 2.0)
  #  result.lut[i] = sin(angle)


proc setFrequency*(gen: var WaveGenerator, freq: float) =
  gen.freq = freq


proc lookupSine*(gen: WaveGenerator, phase: float): float =
  # Convert phase (0-2π) to LUT index using quarter-wave symmetry
  let
    normalizedPhase = phase / (2.0 * PI)  # 0-1 range
    scaledPhase = normalizedPhase * 4.0   # 0-4 range for quarters
    quarter = int(scaledPhase)
    fractionalPart = scaledPhase - quarter.float
    index = int(fractionalPart * LutSize.float)
    clampedIndex = min(index, LutSize - 1) 
  case quarter:
  of 0: # 0 to π/2: direct lookup
    return Lut[clampedIndex]
  of 1: # π/2 to π: reverse lookup
    return Lut[LutSize - 1 - clampedIndex]
  of 2: # π to 3π/2: negative direct lookup
    return -Lut[clampedIndex]
  of 3: # 3π/2 to 2π: negative reverse lookup
    return -Lut[LutSize - 1 - clampedIndex]
  else:
    return 0.0


template advancePhase*()=
  let phaseIncrement = 2.0 * PI * gen.freq / gen.sampleRate
  gen.phase += phaseIncrement
  if gen.phase >= 2.0 * PI:
    gen.phase -= 2.0 * PI  


template advanceModPhase*()=
  gen.phase += modulatedIncrement
  if gen.phase >= 2.0 * PI:
    gen.phase -= 2.0 * PI


proc nextSample*(gen: var WaveGenerator): float =
  # Generate sine wave using quarter LUT
  result = gen.lookupSine(gen.phase) 
  advancePhase() 


proc nextSawtooth*(gen: var WaveGenerator): float =
  # "Sawtooth" by stepping through LUT linearly (ignoring symmetry)
  advancePhase()
  let
    normalizedPhase = gen.phase / (2.0 * PI)
    index = int(normalizedPhase * LutSize.float)
    clampedIndex = min(index, LutSize - 1)
  return Lut[clampedIndex] * 2.0 - 1.0  # Center around 0


proc nextTriangle*(gen: var WaveGenerator): float =
  # Triangle wave using bidirectional LUT traversal
  advancePhase()
  let normalizedPhase = gen.phase / (2.0 * PI)  
  if normalizedPhase < 0.5:
    # Forward through LUT
    let
      index = int(normalizedPhase * 2.0 * LutSize.float)
      clampedIndex = min(index, LutSize - 1)
    return Lut[clampedIndex]
  else:
    # Backward through LUT
    let
      index = int((1.0 - normalizedPhase) * 2.0 * LutSize.float)
      clampedIndex = min(index, LutSize - 1)
    return Lut[clampedIndex]


proc nextSquare*(gen: var WaveGenerator): float =
  # Square wave using LUT threshold
  advancePhase()
  let
    normalizedPhase = gen.phase / (2.0 * PI)
    index = int(normalizedPhase * LutSize.float)
    clampedIndex = min(index, LutSize - 1)
  return if Lut[clampedIndex] > 0.5: 1.0 else: -1.0


proc nextRingMod*(gen: var WaveGenerator, modFreq: float): float =
  ## Ring Modulation, Multiply two oscillators
  let
    carrier = gen.lookupSine(gen.phase)  
    # Modulator, different frequency
    modPhase = (gen.phase * modFreq) / gen.freq
    modulator = gen.lookupSine(modPhase)
  result = carrier * modulator
  advancePhase()


proc nextAM*(
  gen: var WaveGenerator, modFreq: float, modDepth: float = 0.5
): float =
  ## Amplitude Modulation
  let
    carrier = gen.lookupSine(gen.phase)
    # Amplitude modulator, slower than carrier
    modPhase = (gen.phase * modFreq) / gen.freq
    modulator = gen.lookupSine(modPhase)
  result = carrier * (1.0 + modDepth * modulator)
  advancePhase()


proc nextWaveshaped*(gen: var WaveGenerator, drive: float = 2.0): float =
  ## Waveshaping, use LUT as a transfer function
  let
    input = gen.lookupSine(gen.phase)
    # Use the input to index the LUT as a waveshaper, drive controls the input 
    # signal stretch
    drivenInput = input * drive
    clampedInput = max(-1.0, min(1.0, drivenInput))
    # Map -1..1 input to 0..1 for LUT indexing
    normalizedInput = (clampedInput + 1.0) * 0.5
    lutIndex = int(normalizedInput * (LutSize - 1).float)
    clampedIndex = min(max(lutIndex, 0), LutSize - 1)
  result = Lut[clampedIndex] * 2.0 - 1.0
  advancePhase()


proc nextDistorted*(gen: var WaveGenerator, distType: int = 0): float =
  ## Waveshaping: Use sine LUT to create different transfer curves
  let input = gen.lookupSine(gen.phase)
  case distType:
  of 0: # Soft clipping using tanh-like curve from LUT
    let absInput = abs(input)
    let sign = if input >= 0.0: 1.0 else: -1.0
    let lutIndex = int(absInput * (LutSize - 1).float)
    let clampedIndex = min(lutIndex, LutSize - 1)
    result = sign * Lut[clampedIndex]  
  of 1: # Asymmetric distortion
    if input >= 0.0:
      let lutIndex = int(input * (LutSize - 1).float)
      let clampedIndex = min(lutIndex, LutSize - 1)
      result = Lut[clampedIndex]
    else:
      # Different curve for negative values
      let lutIndex = int(abs(input) * (LutSize - 1).float)
      let clampedIndex = min(lutIndex, LutSize - 1)
      result = -Lut[clampedIndex] * 0.7  # Asymmetric scaling
  else: # Fold-back distortion
    var foldedInput = input
    while abs(foldedInput) > 1.0:
      if foldedInput > 1.0:
        foldedInput = 2.0 - foldedInput
      elif foldedInput < -1.0:
        foldedInput = -2.0 - foldedInput
    let
      absInput = abs(foldedInput)
      sign = if foldedInput >= 0.0: 1.0 else: -1.0
      lutIndex = int(absInput * (LutSize - 1).float)
      clampedIndex = min(lutIndex, LutSize - 1)
    result = sign * Lut[clampedIndex]
  advancePhase()


proc nextPhaseDistortion*(gen: var WaveGenerator, distAmount: float = 0.5): float =
  ## Phase Distorion, Casio CZ: Manipulate phase increment,
  ## varying the speed of wavetable traversel
  let
    basePhaseIncrement = 2.0 * PI * gen.freq / gen.sampleRate  
    # Phase distortion curve using LUT
    distPhase = gen.phase * 2.0  # Double frequency for distortion curve
    distortionCurve = gen.lookupSine(distPhase)  
    # Modulate phase increment, move faster/slower through the main waveform
    modulatedIncrement = basePhaseIncrement * (1.0 + distAmount * distortionCurve)
  result = gen.lookupSine(gen.phase)
  advanceModPhase()


proc nextResonantSweep*(gen: var WaveGenerator, sweepRate: float = 0.1): float =
  ## PD: Resonant sweep, slow sweep to modulate phase incr.
  let
    basePhaseIncrement = 2.0 * PI * gen.freq / gen.sampleRate
    sweepPhase = gen.phase * sweepRate
    sweep = gen.lookupSine(sweepPhase)
    modulatedIncrement = basePhaseIncrement * (1.0 + 0.8 * sweep)
  result = gen.lookupSine(gen.phase)
  advanceModPhase()


proc nextCZStyle*(gen: var WaveGenerator, segments: seq[float]): float =
  ## PD: Multi-segment more like real CZ
  ## CZ synths divided the waveform into segments with different rates
  ## segments should be 4 values representing rate multipliers for each quarter
  let
    basePhaseIncrement = 2.0 * PI * gen.freq / gen.sampleRate
    normalizedPhase = gen.phase / (2.0 * PI)  # 0-1
    # Which segment? (4 segments)
    #segmentIndex = int(normalizedPhase * 4.0)
    #clampedSegment = min(segmentIndex, 3)
    segmentIndex = int(normalizedPhase * segments.len.float)
    clampedSegment = min(segmentIndex, segments.len - 1)
    # Get rate multiplier for current segment
    segmentRate = if clampedSegment < segments.len: 
                    segments[clampedSegment] 
                  else: 1.0
    modulatedIncrement = basePhaseIncrement * segmentRate
  result = gen.lookupSine(gen.phase)
  advanceModPhase()


proc nextPDSaw*(gen: var WaveGenerator): float =
  ## PD: Sawtooth-like, using variable rate by making first half of cycle fast,
  ## second half slow
  let
    basePhaseIncrement = 2.0 * PI * gen.freq / gen.sampleRate
    normalizedPhase = gen.phase / (2.0 * PI)
    rateMultiplier = if normalizedPhase < 0.5: 3.0 else: 0.2
    modulatedIncrement = basePhaseIncrement * rateMultiplier
  result = gen.lookupSine(gen.phase)
  advanceModPhase()


proc nextPDSquare*(gen: var WaveGenerator): float =
  ## PD: Square-like, using phase jumps
  let
    basePhaseIncrement = 2.0 * PI * gen.freq / gen.sampleRate
    normalizedPhase = gen.phase / (2.0 * PI)  
    # Slow down near 0.25 and 0.75
    distanceFromPeak1 = abs(normalizedPhase - 0.25)
    distanceFromPeak2 = abs(normalizedPhase - 0.75)
    nearPeak = min(distanceFromPeak1, distanceFromPeak2)  
    # Speed up elsewhere
    rateMultiplier = if nearPeak < 0.1: 0.1 else: 2.0
    modulatedIncrement = basePhaseIncrement * rateMultiplier
  result = gen.lookupSine(gen.phase)
  gen.phase += modulatedIncrement
  if gen.phase >= 2.0 * PI:
    gen.phase -= 2.0 * PI


when isMainModule:
  var gen = initWaveGenerator()
  gen.setFrequency(440.0)  # A4
  #var tickerTape = Ticker(tick: 0)
  var sample = 0.0

  proc tick*(): float =
    sample =  gen.nextSample()
    #sample =  gen.nextSawtooth()
    #sample =  gen.nextTriangle()
    #sample =  gen.nextSquare()
    #sample =  gen.nextPhaseDistortion(0.5)
    #sample =  gen.nextResonantSweep(0.8)
    #let czSegments = @[2.0, 0.1, 0.5, 0.3, 1.5, 0.7]
    #sample =  gen.nextCZStyle(czSegments)
    #sample =  gen.nextPDSaw()
    #sample =  gen.nextPDSquare()
    #sample =  gen.nextRingMod(100.0)
    #sample =  gen.nextAM(110.0, 0.9)
    #sample =  gen.nextWaveshaped(3.0)
    #sample =  gen.nextDistorted(0)
    #sample =  gen.nextDistorted(1)
    #sample =  gen.nextDistorted(2)
    #inc tickerTape.tick
    return sample

  include io    #libSoundIO
  var ss = newSoundSystem()
  ss.outstream.sampleRate = SampleRate
  let outstream = ss.outstream
  let sampleRate = outstream.sampleRate.toFloat
  
  echo "Format:\t\t", outstream.format
  echo "Sample Rate:\t", sampleRate
  echo "Latency:\t", outstream.softwareLatency
  
  while true:
    ss.sio.flushEvents
    let s = stdin.readLine
    if s == "q":
      break

io.nim

#import std/[math]
import soundio

type
  SoundSystem* = object
    sio*: ptr SoundIo
    indevice*: ptr SoundIoDevice
    instream*: ptr SoundIoInStream
    outdevice*: ptr SoundIoDevice
    outstream*: ptr SoundIoOutStream
    #outsource*: proc()


proc `=destroy`(s: SoundSystem) =
  if not isNil s.outstream:
    echo "destroy outstream"
    s.outstream.destroy
    dealloc(s.outstream.userdata)
  if not isNil s.outdevice:
    echo "destroy outdevice"
    s.outdevice.unref
  echo "destroy SoundSystem"
  s.sio.destroy
  echo "Quit"


proc writeCallback(outStream: ptr SoundIoOutStream, frameCountMin: cint, frameCountMax: cint) {.cdecl.} =
  let csz = sizeof SoundIoChannelArea
  var areas: ptr SoundIoChannelArea
  var framesLeft = frameCountMax
  var err: cint
  
  while true:
    var frameCount = framesLeft
    err = outStream.beginWrite(areas.addr, frameCount.addr)
    if err > 0:
      quit "Unrecoverable stream error: " & $err.strerror
    if frameCount <= 0:
      break
    let layout = outstream.layout
    let ptrAreas = cast[int](areas)
    for frame in 0..<frameCount:
      
      let sample = tick()

      for channel in 0..<layout.channelCount:
        let ptrArea = cast[ptr SoundIoChannelArea](ptrAreas + channel*csz)        
        var ptrSample = cast[ptr float32](cast[int](ptrArea.pointer) + frame*ptrArea.step)        
        ptrSample[] = sample
        
    err = outstream.endWrite
    if err > 0 and err != cint(SoundIoError.Underflow):
      quit "Unrecoverable stream error: " & $err.strerror
    framesLeft -= frameCount
    if framesLeft <= 0:
      break

proc newSoundSystem*(): SoundSystem =
  echo "SoundIO version : ", version_string()
  result.sio = soundioCreate()
  if isNil result.sio: quit "out of mem"
  var err = result.sio.connect
  if err > 0: quit "Unable to connect to backend: " & $err.strerror
  echo "Backend: \t", result.sio.currentBackend.name
  result.sio.flushEvents

  echo "Output:"
  let outdevID = result.sio.defaultOutputDeviceIndex
  echo "  Device Index: \t", outdevID
  if outdevID < 0: quit "Output device is not found"
  result.outdevice = result.sio.getOutputDevice(outdevID)
  if isNil result.outdevice: quit "out of mem"
  if result.outdevice.probeError > 0: quit "Cannot probe device"
  echo "  Device Name:\t", result.outdevice.name

  result.outstream = result.outdevice.outStreamCreate
  result.outstream.write_callback = writeCallback

  err = result.outstream.open
  if err > 0: quit "Unable to open output device: " & $err.strerror 
  if result.outstream.layoutError > 0:
    quit "Unable to set channel layout: " & $result.outstream.layoutError.strerror  
  err = result.outstream.start
  if err > 0: quit "Unable to start stream: " & $err.strerror