#### Based on libsamplerate
####
#### Copyright (c) 2002-2021, Erik de Castro Lopo <erikd@mega-nerd.com>
#### Copyright (c) 2024, Remilia Scarlet <remilia@posteo.jp>
#### All rights reserved.
####
#### This code is released under 2-clause BSD license. Please see the
#### file at : https://github.com/libsndfile/libsamplerate/blob/master/COPYING
require "./resampler"
@[Link(ldflags: "#{__DIR__}/libsinctables.a")]
lib LibRASincTables
$sinc_fast_coeffs_size : LibC::Int
$sinc_medium_coeffs_size : LibC::Int
$sinc_best_coeffs_size : LibC::Int
$sinc_fast_coeffs_index_inc : LibC::Int
$sinc_medium_coeffs_index_inc : LibC::Int
$sinc_best_coeffs_index_inc : LibC::Int
$sinc_fast_coeffs : StaticArray(Float32, 2464)
$sinc_medium_coeffs : StaticArray(Float32, 22438)
$sinc_best_coeffs : StaticArray(Float32, 340239)
end
module RemiAudio::Resampler
# A bandlimited interpolator derived from the mathematical sinc function.
# Based on the techniques of [Julius
# O. Smith](http://ccrma.stanford.edu/~jos/resample/), although this code was
# developed independently.
abstract class SincResampler < Resampler
MAX_CHANNELS = 128
# :nodoc:
SHIFT_BITS = 12
# :nodoc:
FP_ONE = (1 << SHIFT_BITS).to_f64!
# :nodoc
INV_FP_ONE = 1.0 / FP_ONE
# Defines the quality of a `SincResampler`.
enum Quality
# The fastest, but lowest quality, for a `SincResampler`. The resampler
# will have a signal-to-noise ratio of 97dB and a bandwidth of 80%.
Fast
# A fairly fast, medium quality, for a `SincResampler`. The resampler will
# have a signal-to-noise ratio of 97dB and a bandwidth of 90%
Medium
# A very slow, but very high quality, for a `SincResampler`. This
# provides a worst case signal-to-noise ratio of 97dB at a bandwidth of
# 97%
Best
end
# :nodoc:
alias Coefficients = Pointer(Float32)
@inCount : Int64 = 0i64
@inUsed : Int64 = 0i64
@outCount : Int64 = 0i64
@outGen : Int64 = 0i64
@coeffHalfLen : Int32 = 0
@indexInc : Int32 = 0
@srcRatio : Float64 = 0.0
@inputIndex : Float64 = 0.0
@coeffs : Coefficients
@bCurrent : Int32 = 0
@bEnd : Int32 = 0
@bRealEnd : Int32 = 0
@bLen : Int32 = 0
@leftCalc : Array(Float64) = Array(Float64).new(MAX_CHANNELS, 0.0)
@rightCalc : Array(Float64) = Array(Float64).new(MAX_CHANNELS, 0.0)
@buffer : Array(Float32)
# The `Quality` for this instance.
getter quality : Quality
# Appease the compiler deities
private def initialize(@channels, @quality, @coeffs, @buffer)
end
# :inherit:
def reset : Nil
super
@bCurrent = 0
@bEnd = 0
@bRealEnd = -1
@srcRatio = 0.0
@inputIndex = 0.0
@buffer.fill(0.0f32)
end
@[AlwaysInline]
protected def toFp(x : Float64) : Int32
(x.round * FP_ONE).to_i32!
end
@[AlwaysInline]
protected def toFp(x : Int) : Int32
(x << SHIFT_BITS).to_i32!
end
@[AlwaysInline]
protected def fpToInt(x : Int32) : Int32
x >> SHIFT_BITS
end
@[AlwaysInline]
protected def fpFracPart(x : Int32)
x & ((1 << SHIFT_BITS) &- 1)
end
@[AlwaysInline]
protected def fpToF64(x : Int32) : Float64
fpFracPart(x) * INV_FP_ONE
end
@[AlwaysInline]
protected def intDivCeil(divident : Int32, divisor : Int32) # == (divident.to_f32! / divisor).ceiling.to_i32!
RemiLib.assert(divident >= 0 && divisor > 0) # For positive numbers only
(divident &+ (divisor &- 1)).tdiv(divisor)
end
protected def prepareData(data : Data, halfFilterChanLen : Int32) : Nil
return if @bRealEnd >= 0 # Should be terminating, just return
return if data.dataIn.empty?
len : Int32 = 0
if @bCurrent == 0
# Initial state. Set up zeros at the start of the buffer and then load
# new data after that.
len = @bLen - 2 * halfFilterChanLen
@bCurrent = halfFilterChanLen
@bEnd = halfFilterChanLen
elsif (@bEnd + halfFilterChanLen + @channels) < @bLen
# Load data at current end position.
len = Math.max(@bLen - @bCurrent - halfFilterChanLen, 0)
else
# Move data at end of buffer back to the start of the buffer.
len = @bEnd - @bCurrent
@buffer.to_unsafe.move_from(@buffer.to_unsafe + @bCurrent - halfFilterChanLen,
halfFilterChanLen + len)
@bCurrent = halfFilterChanLen
@bEnd = @bCurrent + len
# Now load data at current end of buffer.
len = Math.max(@bLen - @bCurrent - halfFilterChanLen, 0)
end
len = Math.min((@inCount - @inUsed).to_i32!, len)
len -= (len % @channels)
if len < 0 || (@bEnd + len) > @bLen
raise Error.new("prepareData: bad length")
end
(@buffer.to_unsafe + @bEnd).copy_from(data.dataIn.to_unsafe + @inUsed, len)
@bEnd += len
@inUsed += len
if @inUsed == @inCount && (@bEnd - @bCurrent) < (2 * halfFilterChanLen) && data.endOfInput?
# Handle the case where all data in the current buffer has been consumed
# and this is the last buffer.
if (@bLen - @bEnd) < (halfFilterChanLen + 5)
# If necessary, move data down to the start of the buffer.
len = @bEnd - @bCurrent
@buffer.to_unsafe.move_from(@buffer.to_unsafe + @bCurrent - halfFilterChanLen,
halfFilterChanLen + len)
@bCurrent = halfFilterChanLen
@bEnd = @bCurrent + len
end
@bRealEnd = @bEnd
len = halfFilterChanLen + 5
if len < 0 || (@bEnd + len) > @bLen
len = @bLen - @bEnd
end
@buffer.fill(0.0, @bEnd, len)
@bEnd += len
end
end
end
# A `SincResampler` optimized for two channels.
class SincResamplerStereo < SincResampler
# Creates a new `SincResamplerStereo` instance that will operate at the
# given `Quality`.
def initialize(@quality : Quality)
@channels = 2
case @quality
in .fast?
@coeffs = LibRASincTables.sinc_fast_coeffs.to_unsafe
@coeffHalfLen = LibRASincTables.sinc_fast_coeffs_size.to_i32 - 2
@indexInc = LibRASincTables.sinc_fast_coeffs_index_inc.to_i32
in .medium?
@coeffs = LibRASincTables.sinc_medium_coeffs.to_unsafe
@coeffHalfLen = LibRASincTables.sinc_medium_coeffs_size.to_i32 - 2
@indexInc = LibRASincTables.sinc_medium_coeffs_index_inc.to_i32
in .best?
@coeffs = LibRASincTables.sinc_best_coeffs.to_unsafe
@coeffHalfLen = LibRASincTables.sinc_best_coeffs_size.to_i32 - 2
@indexInc = LibRASincTables.sinc_best_coeffs_index_inc.to_i32
end
@bLen = 3 * ((@coeffHalfLen + 2.0) / @indexInc * SRC_MAX_RATIO + 1).round.to_i32!
@bLen = Math.max(@bLen, 4096)
@bLen *= @channels
@bLen += 1
@buffer = Array(Float32).new(@bLen * @channels, 0.0f32)
reset
end
private def calcOutput(increment : Int32, startFilterIndex : Int32, scale : Float64, output : Slice(Float32)) : Nil
fraction : Float64 = 0.0
left : StaticArray(Float64, 2) = StaticArray(Float64, 2).new(0.0)
right : StaticArray(Float64, 2) = StaticArray(Float64, 2).new(0.0)
icoeff : Float64 = 0.0
idx : Int32 = 0
# Convert input parameters into fixed point.
maxFilterIndex : Int32 = toFp(@coeffHalfLen)
# First apply the left half of the filter.
filterIndex : Int32 = startFilterIndex
coeffCount : Int32 = (maxFilterIndex - filterIndex).tdiv(increment)
filterIndex = filterIndex + coeffCount * increment
dataIndex : Int32 = @bCurrent - @channels * coeffCount
if dataIndex < 0 # Avoid underflow access to filter->buffer.
steps : Int32 = intDivCeil(-dataIndex, 2)
# If the assert triggers we would have to take care not to underflow/overflow
RemiLib.assert(steps <= intDivCeil(filterIndex, increment))
filterIndex -= (increment * steps)
dataIndex += (steps * 2)
end
while filterIndex >= 0
fraction = fpToF64(filterIndex)
idx = fpToInt(filterIndex)
#RemiLib.assert(idx >= 0 && (idx + 1) < (@coeffHalfLen + 2))
icoeff = @coeffs[idx].to_f64! + fraction * (@coeffs[idx + 1] - @coeffs[idx])
#RemiLib.assert(dataIndex >= 0 && (dataIndex + 1) < @bLen)
#RemiLib.assert((dataIndex + 1) < @bEnd)
2.times do |ch|
left[ch] += (icoeff * @buffer[dataIndex + ch])
end
filterIndex -= increment
dataIndex = dataIndex + 2
end
# Now apply the right half of the filter.
filterIndex = increment - startFilterIndex
coeffCount = (maxFilterIndex - filterIndex).tdiv(increment)
filterIndex = filterIndex + coeffCount * increment
dataIndex = @bCurrent + @channels * (1 + coeffCount)
loop do
fraction = fpToF64(filterIndex)
idx = fpToInt(filterIndex)
#RemiLib.assert(idx >= 0 && (idx + 1) < (@coeffHalfLen + 2))
icoeff = @coeffs[idx].to_f64! + fraction * (@coeffs[idx + 1] - @coeffs[idx])
#RemiLib.assert(dataIndex >= 0 && (dataIndex + 1) < @bLen)
#RemiLib.assert("dataIndex = #{dataIndex}, bEnd: #{@bEnd}",
# (dataIndex + 1) < @bEnd)
2.times do |ch|
right[ch] += (icoeff * @buffer[dataIndex + ch])
end
filterIndex -= increment
dataIndex = dataIndex - 2
break unless filterIndex > 0
end
2.times do |ch|
output[ch] = (scale * (left[ch] + right[ch])).to_f32!
end
end
protected def variProcess(data : Data) : Nil
samplesInHand : Int32 = 0
floatIncrement : Float64 = 0.0
startFilterIndex : Int32 = 0
increment : Int32 = 0
@inCount = data.inputFrames * @channels
@outCount = data.outputFrames * @channels
@inUsed = 0
@outGen = 0
ratio : Float64 = @lastRatio
# Check the sample rate ratio wrt the buffer len.
count : Float64 = (@coeffHalfLen + 2.0) / @indexInc
if Math.min(@lastRatio, data.srcRatio) < 1.0
count /= Math.min(@lastRatio, data.srcRatio)
end
# Maximum coefficients on either side of center point.
halfFilterChanLen : Int32 = @channels * (count.round.to_i32! &+ 1)
inputIndex : Float64 = @lastPosition
rem : Float64 = fmodOne(inputIndex)
@bCurrent = ((@bCurrent.to_i64! + @channels * (inputIndex - rem).round.to_i64!) % @bLen).to_i32!
inputIndex = rem
terminate : Float64 = 1.0 / ratio + 1e-20
# Main processing loop.
while @outGen < @outCount
# Need to reload buffer?
samplesInHand = (@bEnd - @bCurrent + @bLen) % @bLen
if samplesInHand <= halfFilterChanLen
prepareData(data, halfFilterChanLen)
samplesInHand = (@bEnd - @bCurrent + @bLen) % @bLen
break if samplesInHand <= halfFilterChanLen
end
# This is the termination condition.
if @bRealEnd >= 0
break if (@bCurrent + inputIndex + terminate) >= @bRealEnd
end
if @outCount > 0 && (@lastRatio - data.srcRatio).abs > 1e-10
ratio = @lastRatio + @outGen * (data.srcRatio - @lastRatio) / @outCount
end
floatIncrement = @indexInc * (ratio < 1.0 ? ratio : 1.0)
increment = toFp(floatIncrement)
startFilterIndex = toFp(inputIndex * floatIncrement)
calcOutput(increment, startFilterIndex, floatIncrement / @indexInc, data.dataOut[@outGen..])
@outGen += 2
# Figure out the next index.
inputIndex += (1.0 / ratio)
rem = fmodOne(inputIndex)
@bCurrent = ((@bCurrent.to_i64! + @channels * (inputIndex - rem).round.to_i64!) % @bLen).to_i32!
inputIndex = rem
end
@lastPosition = inputIndex
# Save current ratio rather then target ratio.
@lastRatio = ratio
data.inputFramesUsed = @inUsed.tdiv(@channels)
data.outputFramesGen = @outGen.tdiv(@channels)
end
end
# # A callback version of the `SincResamplerStereo` class.
class SincResamplerStereoCb < SincResamplerStereo
include CallbackResampler
private def initialize(@channels : Int32)
@coeffs = [] of Float32
@buffer = [] of Float32
@quality = Quality::Fast
@callbackFunc = ->{ {Slice(Float32).empty, 0i64} }
raise NotImplementedError.new("Use initialize(Int32, Float64, CallbackProc)")
end
# Creates a new `SincResamplerStereoCb` that will operate at the given
# `Quality`. The *ratio* parameter sets the initial ratio of the resampler
# (`targetRate / sourceRate`), while *callbackFunc* is the method that will
# be called whenever this resampler needs more input data.
def initialize(@quality : Quality, ratio : Float64, @callbackFunc : CallbackProc)
@channels = 2
case @quality
in .fast?
@coeffs = LibRASincTables.sinc_fast_coeffs.to_unsafe
@coeffHalfLen = LibRASincTables.sinc_fast_coeffs_size.to_i32 - 2
@indexInc = LibRASincTables.sinc_fast_coeffs_index_inc.to_i32
in .medium?
@coeffs = LibRASincTables.sinc_medium_coeffs.to_unsafe
@coeffHalfLen = LibRASincTables.sinc_medium_coeffs_size.to_i32 - 2
@indexInc = LibRASincTables.sinc_medium_coeffs_index_inc.to_i32
in .best?
@coeffs = LibRASincTables.sinc_best_coeffs.to_unsafe
@coeffHalfLen = LibRASincTables.sinc_best_coeffs_size.to_i32 - 2
@indexInc = LibRASincTables.sinc_best_coeffs_index_inc.to_i32
end
@bLen = 3 * ((@coeffHalfLen + 2.0) / @indexInc * SRC_MAX_RATIO + 1).round.to_i32!
@bLen = Math.max(@bLen, 4096)
@bLen *= @channels
@bLen += 1
@buffer = Array(Float32).new(@bLen * @channels, 0.0f32)
reset
@data.srcRatio = ratio
self.ratio = ratio
end
# This is not supported by this class since `LinearResamplerCb` uses the
# `CallbackResampler` mixin. Calling this always raises a
# `NotImplementedError`.
def process(source : Array(Float32)|Slice(Float32), dest : Array(Float32)|Slice(Float32),
ratio : Float64) : Tuple(Int64, Int64, Bool)
raise NotImplementedError.new("Not supported by callback resamplers")
end
# :ditto:
def process(source : Array(Float64)|Slice(Float64), dest : Array(Float64)|Slice(Float64),
ratio : Float64) : Tuple(Int64, Int64, Bool)
raise NotImplementedError.new("Not supported by callback resamplers")
end
protected def processData : Nil
variProcess(@data)
end
# :inherit:
def reset : Nil
super
@data.reset
end
end
end