use std::alloc::LayoutError;
use no_denormals::*;
use crate::buffer::*;
#[inline]
pub fn ratio_to_db(ratio : f64) -> f64 { 20.0 * ratio.log10() }
#[inline]
pub fn db_to_ratio(db : f64) -> f64 { 10.0f64.powf(db / 20.0) }
pub struct Convolution
{
buffer : PushBuffer<f64>,
kernel : Box<[f64]>
}
impl Convolution
{
pub fn new(kernel : &[f64]) -> Result<Self, LayoutError>
{
let conv = Self
{
buffer : PushBuffer::<f64>::new(kernel.len())?,
kernel : kernel.to_vec().into_boxed_slice()
};
conv.buffer.set_index(conv.buffer.len());
Ok(conv)
}
pub fn kernel_len(&self) -> usize { self.kernel.len() }
#[inline]
pub fn process(&self, input : f64) -> f64
{
let mut guard = self.buffer.write();
guard.push(input);
let mut sum = 0.0;
for i in 0..self.kernel.len()
{
sum += guard[i] * self.kernel[i];
}
sum
}
pub fn run(&self, input : &Buffer<f64>, output : &Buffer<f64>)
{
let input_guard = input.read();
let mut output_guard = output.write();
let mut buffer_guard = self.buffer.write();
no_denormals(||
{
for index in 0..input_guard.len().min(output_guard.len())
{
buffer_guard.push(input_guard[index]);
let mut sum = 0.0;
for i in 0..self.kernel.len()
{
sum += buffer_guard[i] * self.kernel[i];
}
output_guard[index] = sum;
}
});
}
}
pub struct Saturation
{
drive_alpha_plus : f64,
drive_alpha_minus : f64,
compression_beta_plus : f64,
compression_beta_minus : f64,
bias_delta : f64,
flip_polarity : bool,
norm_factor_plus : f64,
norm_factor_minus : f64
}
impl Saturation
{
pub fn new(alpha_plus : f64, alpha_minus : f64,
beta_plus : f64, beta_minus : f64,
delta_bias : f64, flip : bool) -> Self
{
let drive_alpha_plus = alpha_plus.max(1e-4);
let drive_alpha_minus = alpha_minus.max(1e-4);
let norm_factor_plus = 1.0 / (1.0 + drive_alpha_plus).log2();
let norm_factor_minus = 1.0 / (1.0 + drive_alpha_minus).log2();
Self
{
drive_alpha_plus,
drive_alpha_minus,
compression_beta_plus : beta_plus,
compression_beta_minus : beta_minus,
bias_delta : delta_bias,
flip_polarity : flip,
norm_factor_plus,
norm_factor_minus
}
}
#[inline]
pub fn process(&self, input_sample : f64) -> f64
{
let output_value = if input_sample >= self.bias_delta
{
let relative_input = input_sample - self.bias_delta;
let log_out = (1.0 + self.drive_alpha_plus * relative_input).log2();
self.compression_beta_plus * (log_out * self.norm_factor_plus)
}
else
{
let relative_input = self.bias_delta - input_sample;
let log_out = (1.0 + self.drive_alpha_minus * relative_input).log2();
-self.compression_beta_minus * (log_out * self.norm_factor_minus)
};
if self.flip_polarity { -output_value } else { output_value }
}
pub fn run(&self, input : &Buffer<f64>, output : &Buffer<f64>)
{
let input_guard = input.read();
let mut output_guard = output.write();
no_denormals(||
{
for index in 0..input_guard.len().min(output_guard.len())
{
output_guard[index] = self.process(input_guard[index]);
}
});
}
}
pub trait Component
{
fn nodes(&self) -> (i32, i32);
fn get_conductance(&self, dt : f64) -> f64;
fn get_current_source(&self, dt : f64) -> f64;
fn update_state(&mut self, v_a : f64, v_b : f64, dt : f64);
}
pub struct Resistor
{
node_a : i32,
node_b : i32,
conductance : f64
}
impl Resistor
{
pub fn new(n1 : i32, n2 : i32, resistance : f64) -> Self
{
Self { node_a : n1, node_b : n2, conductance : 1.0 / resistance }
}
}
impl Component for Resistor
{
fn nodes(&self) -> (i32, i32) { (self.node_a, self.node_b) }
fn get_conductance(&self, _dt : f64) -> f64 { self.conductance }
fn get_current_source(&self, _dt : f64) -> f64 { 0.0 }
fn update_state(&mut self, _v_a : f64, _v_b : f64, _dt : f64) {}
}
pub struct Capacitor
{
node_a : i32,
node_b : i32,
capacitance : f64,
prev_voltage : f64
}
impl Capacitor
{
pub fn new(n1 : i32, n2 : i32, capacitance : f64) -> Self
{
Self { node_a : n1, node_b : n2, capacitance, prev_voltage : 0.0 }
}
}
impl Component for Capacitor
{
fn nodes(&self) -> (i32, i32) { (self.node_a, self.node_b) }
fn get_conductance(&self, dt : f64) -> f64 { self.capacitance / dt }
fn get_current_source(&self, dt : f64) -> f64 { (self.capacitance / dt) * self.prev_voltage }
fn update_state(&mut self, v_a : f64, v_b : f64, _dt : f64) { self.prev_voltage = v_a - v_b; }
}
pub struct Inductor
{
node_a : i32,
node_b : i32,
inductance : f64,
prev_current : f64
}
impl Inductor
{
pub fn new(n1 : i32, n2 : i32, inductance : f64) -> Self
{
Self { node_a : n1, node_b : n2, inductance, prev_current : 0.0 }
}
}
impl Component for Inductor
{
fn nodes(&self) -> (i32, i32) { (self.node_a, self.node_b) }
fn get_conductance(&self, dt : f64) -> f64 { dt / self.inductance }
fn get_current_source(&self, _dt : f64) -> f64 { -self.prev_current }
fn update_state(&mut self, v_a : f64, v_b : f64, dt : f64)
{
let voltage = v_a - v_b;
self.prev_current += (voltage * dt) / self.inductance;
}
}
pub struct Circuit
{
components : Vec<Box<dyn Component + Send + Sync>>,
num_nodes : usize,
y_static : Box<[f64]>,
y_work : Box<[f64]>,
j : Box<[f64]>,
nodes : Box<[f64]>,
dt : f64
}
impl Circuit
{
pub fn new(sample_rate : f64, num_nodes : usize) -> Self
{
let matrix_size = num_nodes * num_nodes;
Self
{
components : Vec::new(),
num_nodes,
y_static : vec![0.0; matrix_size].into_boxed_slice(),
y_work : vec![0.0; matrix_size].into_boxed_slice(),
j : vec![0.0; num_nodes].into_boxed_slice(),
nodes : vec![0.0; num_nodes].into_boxed_slice(),
dt : 1.0 / sample_rate
}
}
pub fn get_nodes(&self) -> usize { self.num_nodes }
pub fn get_devices(&self) -> usize { self.components.len() }
pub fn add_component(&mut self, component : Box<dyn Component + Send + Sync>)
{
self.components.push(component);
}
pub fn preprocess(&mut self, impedance : f64)
{
let n = self.num_nodes;
self.y_static.fill(0.0);
for comp in &self.components
{
let (n1, n2) = comp.nodes();
let g = comp.get_conductance(self.dt);
if n1 > 0 { self.y_static[(n1 as usize - 1) * n + (n1 as usize - 1)] += g; }
if n2 > 0 { self.y_static[(n2 as usize - 1) * n + (n2 as usize - 1)] += g; }
if n1 > 0 && n2 > 0
{
self.y_static[(n1 as usize - 1) * n + (n2 as usize - 1)] -= g;
self.y_static[(n2 as usize - 1) * n + (n1 as usize - 1)] -= g;
}
}
if n >= 1 { self.y_static[0] += impedance; }
}
fn solve_linear_system(&mut self)
{
let n = self.num_nodes;
self.y_work.copy_from_slice(&self.y_static);
for i in 0..n
{
let mut pivot = i;
let mut max_val = self.y_work[i * n + i].abs();
for k in (i + 1)..n
{
let val = self.y_work[k * n + i].abs();
if val > max_val { max_val = val; pivot = k; }
}
if pivot != i
{
for col in i..n
{
self.y_work.swap(i * n + col, pivot * n + col);
}
self.j.swap(i, pivot);
}
let pivot_val = self.y_work[i * n + i];
if pivot_val.abs() < 1e-9 { continue; }
for k in (i + 1)..n
{
let factor = self.y_work[k * n + i] / pivot_val;
for j in i..n
{
self.y_work[k * n + j] -= factor * self.y_work[i * n + j];
}
self.j[k] -= factor * self.j[i];
}
}
for i in (0..n).rev()
{
let mut sum = 0.0;
for j in (i + 1)..n { sum += self.y_work[i * n + j] * self.nodes[j]; }
self.nodes[i] = (self.j[i] - sum) / self.y_work[i * n + i];
}
}
pub fn process(&mut self, input_voltage : f64, probe_node : usize) -> f64
{
let n = self.num_nodes;
self.j.fill(0.0);
let g_source = 1.0 / 0.1;
self.j[0] += input_voltage * g_source;
for comp in &self.components
{
let is = comp.get_current_source(self.dt);
if is == 0.0 { continue; }
let (n1, n2) = comp.nodes();
if n1 > 0 { self.j[n1 as usize - 1] -= is; }
if n2 > 0 { self.j[n2 as usize - 1] += is; }
}
self.solve_linear_system();
for comp in &mut self.components
{
let (n1, n2) = comp.nodes();
let v1 = if n1 == 0 { 0.0 } else { self.nodes[n1 as usize - 1] };
let v2 = if n2 == 0 { 0.0 } else { self.nodes[n2 as usize - 1] };
comp.update_state(v1, v2, self.dt);
}
if probe_node == 0 || probe_node > n { return 0.0; }
self.nodes[probe_node - 1]
}
}
pub struct Compression
{
pub threshold : f64,
pub ratio : f64,
pub attack : f64,
pub release : f64,
pub makeup : f64,
pub knee : f64,
envelope : f64,
attack_coeff : f64,
release_coeff : f64
}
impl Compression
{
pub fn new(sample_rate : f64) -> Self
{
let mut comp = Self
{
threshold : -20.0,
ratio : 4.0,
attack : 10.0,
release : 100.0,
makeup : 0.0,
knee : 0.0,
envelope : 0.0,
attack_coeff : 0.0,
release_coeff : 0.0
};
comp.update_coefficients(sample_rate);
comp
}
pub fn update_coefficients(&mut self, sample_rate : f64)
{
self.attack_coeff = (-1.0 / (self.attack * 0.001 * sample_rate)).exp();
self.release_coeff = (-1.0 / (self.release * 0.001 * sample_rate)).exp();
}
#[inline]
fn compute_gain(&self, input_db : f64) -> f64
{
if self.knee > 0.0
{
let half_knee = self.knee * 0.5;
let lower = self.threshold - half_knee;
let upper = self.threshold + half_knee;
if input_db <= lower { 0.0 }
else if input_db >= upper
{
(self.threshold + (input_db - self.threshold) / self.ratio) - input_db
}
else
{
let x = input_db - lower;
let slope = 1.0 / self.ratio - 1.0;
slope * x * x / (2.0 * self.knee)
}
}
else
{
if input_db <= self.threshold { 0.0 }
else { (self.threshold + (input_db - self.threshold) / self.ratio) - input_db }
}
}
pub fn run(&mut self, input : &Buffer<f64>, output : &Buffer<f64>)
{
let input_guard = input.read();
let mut output_guard = output.write();
let makeup_linear = db_to_ratio(self.makeup);
no_denormals(||
{
for index in 0..input_guard.len().min(output_guard.len())
{
let input_abs = input_guard[index].abs();
let input_db = if input_abs > 1e-10 { 20.0 * input_abs.log10() } else { -200.0 };
let target_gr = self.compute_gain(input_db);
let coeff = if target_gr < self.envelope { self.attack_coeff } else { self.release_coeff };
self.envelope = target_gr + coeff * (self.envelope - target_gr);
let gain = db_to_ratio(self.envelope) * makeup_linear;
output_guard[index] = input_guard[index] * gain;
}
});
}
}
pub struct Limit
{
pub gain : f64,
pub ceiling : f64,
pub release : f64,
envelope : f64,
release_coeff : f64
}
impl Limit
{
pub fn new(sample_rate : f64) -> Self
{
let mut lim = Self
{
gain : 0.0,
ceiling : 0.0,
release : 100.0,
envelope : 1.0,
release_coeff : 0.0
};
lim.update_coefficients(sample_rate);
lim
}
pub fn update_coefficients(&mut self, sample_rate : f64)
{
self.release_coeff = (-1.0 / (self.release * 0.001 * sample_rate)).exp();
}
#[inline]
pub fn process(&mut self, input : f64) -> f64
{
let gain_linear = db_to_ratio(self.gain);
let ceiling_linear = db_to_ratio(self.ceiling);
let amplified = input * gain_linear;
let abs_sample = amplified.abs();
let target = if abs_sample > ceiling_linear
{
ceiling_linear / abs_sample
}
else { 1.0 };
if target < self.envelope
{
self.envelope = target;
}
else
{
self.envelope = target + self.release_coeff * (self.envelope - target);
}
amplified * self.envelope
}
pub fn run(&mut self, input : &Buffer<f64>, output : &Buffer<f64>)
{
let input_guard = input.read();
let mut output_guard = output.write();
no_denormals(||
{
for index in 0..input_guard.len().min(output_guard.len())
{
output_guard[index] = self.process(input_guard[index]);
}
});
}
}
pub struct Delay
{
time : f64,
sample_rate : f64,
pub feedback : f64,
pub mix : f64,
buffer : CircularBuffer<f64>
}
impl Delay
{
pub fn new(time : f64, sample_rate : f64) -> Self
{
let delay_samples = ((time * 0.001 * sample_rate) as usize).max(1);
Self
{
time,
sample_rate,
feedback : 0.5,
mix : 0.5,
buffer : CircularBuffer::new(delay_samples).unwrap()
}
}
pub fn get_time(&self) -> f64 { self.time }
pub fn set_time(&mut self, time : f64)
{
self.time = time;
let delay_samples = ((time * 0.001 * self.sample_rate) as usize).max(1);
self.buffer.resize(delay_samples).unwrap();
}
pub fn set_sample_rate(&mut self, sample_rate : f64)
{
self.sample_rate = sample_rate;
let delay_samples = ((self.time * 0.001 * sample_rate) as usize).max(1);
self.buffer.resize(delay_samples).unwrap();
}
#[inline]
pub fn process(&self, input : f64) -> f64
{
let mut guard = self.buffer.write();
let delayed = guard.next();
guard.push(input + delayed * self.feedback);
input * (1.0 - self.mix) + delayed * self.mix
}
pub fn run(&self, input : &Buffer<f64>, output : &Buffer<f64>)
{
let input_guard = input.read();
let mut output_guard = output.write();
let mut buffer_guard = self.buffer.write();
no_denormals(||
{
for index in 0..input_guard.len().min(output_guard.len())
{
let delayed = buffer_guard.next();
buffer_guard.push(input_guard[index] + delayed * self.feedback);
output_guard[index] = input_guard[index] * (1.0 - self.mix) + delayed * self.mix;
}
});
}
}