You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.

530 lines
14KB

  1. // Mutable Instruments Shelves emulation for VCV Rack
  2. // Copyright (C) 2020 Tyler Coy
  3. //
  4. // This program is free software: you can redistribute it and/or modify
  5. // it under the terms of the GNU General Public License as published by
  6. // the Free Software Foundation, either version 3 of the License, or
  7. // (at your option) any later version.
  8. //
  9. // This program is distributed in the hope that it will be useful,
  10. // but WITHOUT ANY WARRANTY; without even the implied warranty of
  11. // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
  12. // GNU General Public License for more details.
  13. //
  14. // You should have received a copy of the GNU General Public License
  15. // along with this program. If not, see <https://www.gnu.org/licenses/>.
  16. #include <cmath>
  17. #include <algorithm>
  18. #include <rack.hpp>
  19. #include "aafilter.hpp"
  20. namespace shelves
  21. {
  22. using namespace rack;
  23. // Knob ranges
  24. static const float kFreqKnobMin = 20.f;
  25. static const float kFreqKnobMax = 20e3f;
  26. static const float kGainKnobRange = 18.f;
  27. static const float kQKnobMin = 0.5f;
  28. static const float kQKnobMax = 40.f;
  29. // Opamp saturation voltage
  30. static const float kClampVoltage = 10.5f;
  31. // Filter core
  32. static const float kFilterMaxCutoff = kFreqKnobMax;
  33. static const float kFilterR = 100e3f;
  34. static const float kFilterRC = 1.f / (2.f * M_PI * kFilterMaxCutoff);
  35. static const float kFilterC = kFilterRC / kFilterR;
  36. // Frequency CV amplifier
  37. static const float kFreqAmpR = 18e3f;
  38. static const float kFreqAmpC = 560e-12f;
  39. static const float kFreqAmpInputR = 100e3f;
  40. static const float kMinVOct = -kClampVoltage * kFreqAmpInputR / kFreqAmpR;
  41. static const float kFreqKnobVoltage = std::log2f(kFreqKnobMax / kFreqKnobMin);
  42. // The 2164's gain constant is -33mV/dB
  43. static const float kVCAGainConstant = -33e-3f;
  44. inline float QFactorToVoltage(float q_factor)
  45. {
  46. // Calculate the voltage at the VCA control port for a given Q factor
  47. return 20.f * -kVCAGainConstant * std::log10(2.f * q_factor);
  48. }
  49. // Q CV amplifier
  50. static const float kQAmpR = 22e3f;
  51. static const float kQAmpC = 560e-12f;
  52. static const float kQAmpGain = -kQAmpR / 150e3f;
  53. static const float kQKnobMinVoltage = QFactorToVoltage(kQKnobMin) / kQAmpGain;
  54. static const float kQKnobMaxVoltage = QFactorToVoltage(kQKnobMax) / kQAmpGain;
  55. // Gain CV amplifier
  56. static const float kGainPerVolt = 10.f * std::log10(2.f); // 3dB per volt
  57. static const float kMaximumGain = 24.f;
  58. // Clipping indicator
  59. static const float kClipLEDThreshold = 7.86f;
  60. static const float kClipInputR = 150e3f;
  61. static const float kClipInputC = 100e-9f;
  62. static const float kClipLEDRiseTime = 2e-3f;
  63. static const float kClipLEDFallTime = 10e-3f;
  64. // Solves an ODE system using the 2nd order Runge-Kutta method
  65. template <typename T, typename F>
  66. inline T StepRK2(float dt, T y, F f)
  67. {
  68. T k1 = f(y);
  69. T k2 = f(y + k1 * (dt / 2.f));
  70. return y + dt * k2;
  71. }
  72. template <int len, typename T, typename F>
  73. inline void StepRK2(float dt, T y[], F f)
  74. {
  75. T k1[len];
  76. T k2[len];
  77. T yi[len];
  78. f(y, k1);
  79. for (int i = 0; i < len; i++)
  80. {
  81. yi[i] = y[i] + k1[i] * (dt / 2.f);
  82. }
  83. f(yi, k2);
  84. for (int i = 0; i < len; i++)
  85. {
  86. y[i] += k2[i] * dt;
  87. }
  88. }
  89. template <typename T>
  90. class LPFilter
  91. {
  92. public:
  93. void Init(void)
  94. {
  95. voltage_ = 0.f;
  96. }
  97. T Process(float timestep, T v_in, T vca_level)
  98. {
  99. T rad_per_s = -vca_level / kFilterRC;
  100. voltage_ = StepRK2(timestep, voltage_, [&](T v_state)
  101. {
  102. return rad_per_s * (v_in + v_state);
  103. });
  104. voltage_ = simd::clamp(voltage_, -kClampVoltage, kClampVoltage);
  105. return voltage_;
  106. }
  107. T voltage(void)
  108. {
  109. return voltage_;
  110. }
  111. protected:
  112. T voltage_;
  113. };
  114. template <typename T>
  115. class SVFilter
  116. {
  117. public:
  118. void Init(void)
  119. {
  120. bp_ = 0.f;
  121. lp_ = 0.f;
  122. hp_ = 0.f;
  123. }
  124. T Process(float timestep, T in, T vca_level, T q_level)
  125. {
  126. // Thanks Émilie!
  127. // https://mutable-instruments.net/archive/documents/svf_analysis.pdf
  128. //
  129. // For the bandpass integrator,
  130. // dv/dt = -A/RC * hp
  131. // For the lowpass integrator,
  132. // dv/dt = -A/RC * bp
  133. // with
  134. // hp = -(in + lp - 2*Q*bp)
  135. T rad_per_s = -vca_level / kFilterRC;
  136. StepRK2<2>(timestep, voltage_, [&](const T v_state[], T v_stepped[])
  137. {
  138. T lp = v_state[0];
  139. T bp = v_state[1];
  140. T hp = -(in + lp - 2.f * q_level * bp);
  141. v_stepped[0] = rad_per_s * bp;
  142. v_stepped[1] = rad_per_s * hp;
  143. });
  144. lp_ = simd::clamp(lp_, -kClampVoltage, kClampVoltage);
  145. bp_ = simd::clamp(bp_, -kClampVoltage, kClampVoltage);
  146. T out = bp() * -2.f * q_level;
  147. hp_ = -(in + lp() + out);
  148. return out;
  149. }
  150. T hp(void)
  151. {
  152. return hp_;
  153. }
  154. T bp(void)
  155. {
  156. return bp_;
  157. }
  158. T lp(void)
  159. {
  160. return lp_;
  161. }
  162. protected:
  163. union
  164. {
  165. T voltage_[2];
  166. struct
  167. {
  168. T lp_;
  169. T bp_;
  170. };
  171. };
  172. T hp_;
  173. };
  174. class ShelvesEngine
  175. {
  176. public:
  177. struct Frame
  178. {
  179. // Parameters
  180. float hs_freq_knob; // 0 to 1 linear
  181. float hs_gain_knob; // -1 to 1 linear
  182. float p1_freq_knob; // 0 to 1 linear
  183. float p1_gain_knob; // -1 to 1 linear
  184. float p1_q_knob; // 0 to 1 linear
  185. float p2_freq_knob; // 0 to 1 linear
  186. float p2_gain_knob; // -1 to 1 linear
  187. float p2_q_knob; // 0 to 1 linear
  188. float ls_freq_knob; // 0 to 1 linear
  189. float ls_gain_knob; // -1 to 1 linear
  190. // Inputs
  191. float hs_freq_cv;
  192. float hs_gain_cv;
  193. float p1_freq_cv;
  194. float p1_gain_cv;
  195. float p1_q_cv;
  196. float p2_freq_cv;
  197. float p2_gain_cv;
  198. float p2_q_cv;
  199. float ls_freq_cv;
  200. float ls_gain_cv;
  201. float global_freq_cv;
  202. float global_gain_cv;
  203. float main_in;
  204. bool hs_freq_cv_connected;
  205. bool hs_gain_cv_connected;
  206. bool p1_freq_cv_connected;
  207. bool p1_gain_cv_connected;
  208. bool p1_q_cv_connected;
  209. bool p2_freq_cv_connected;
  210. bool p2_gain_cv_connected;
  211. bool p2_q_cv_connected;
  212. bool ls_freq_cv_connected;
  213. bool ls_gain_cv_connected;
  214. bool global_freq_cv_connected;
  215. bool global_gain_cv_connected;
  216. // Outputs
  217. float p1_hp_out;
  218. float p1_bp_out;
  219. float p1_lp_out;
  220. float p2_hp_out;
  221. float p2_bp_out;
  222. float p2_lp_out;
  223. float main_out;
  224. bool p1_hp_out_connected;
  225. bool p1_bp_out_connected;
  226. bool p1_lp_out_connected;
  227. bool p2_hp_out_connected;
  228. bool p2_bp_out_connected;
  229. bool p2_lp_out_connected;
  230. // Lights
  231. float clip;
  232. // Options
  233. bool pre_gain; // True = -6dB, False = 0dB
  234. };
  235. ShelvesEngine()
  236. {
  237. setSampleRate(1.f);
  238. }
  239. void setSampleRate(float sample_rate)
  240. {
  241. sample_time_ = 1.f / sample_rate;
  242. oversampling_ = OversamplingFactor(sample_rate);
  243. up_filter_[0].Init(sample_rate);
  244. up_filter_[1].Init(sample_rate);
  245. up_filter_[2].Init(sample_rate);
  246. down_filter_[0].Init(sample_rate);
  247. down_filter_[1].Init(sample_rate);
  248. low_high_.Init();
  249. mid_.Init();
  250. float freq_cut = 1.f / (2.f * M_PI * kFreqAmpR * kFreqAmpC);
  251. freq_lpf_.reset();
  252. freq_lpf_.setCutoffFreq(freq_cut / sample_rate);
  253. float q_cut = 1.f / (2.f * M_PI * kQAmpR * kQAmpC);
  254. q_lpf_.reset();
  255. q_lpf_.setCutoffFreq(q_cut / sample_rate);
  256. float clip_in_cut = 1.f / (2.f * M_PI * kClipInputR * kClipInputC);
  257. clip_hpf_.reset();
  258. clip_hpf_.setCutoffFreq(clip_in_cut / sample_rate);
  259. float rise = 1.f / kClipLEDRiseTime;
  260. float fall = 1.f / kClipLEDFallTime;
  261. clip_slew_.reset();
  262. clip_slew_.setRiseFall(rise, fall);
  263. }
  264. void process(Frame& frame)
  265. {
  266. auto f_knob = simd::float_4(
  267. frame.ls_freq_knob,
  268. frame.p1_freq_knob,
  269. frame.p2_freq_knob,
  270. frame.hs_freq_knob);
  271. auto f_cv = simd::float_4(
  272. frame.ls_freq_cv,
  273. frame.p1_freq_cv,
  274. frame.p2_freq_cv,
  275. frame.hs_freq_cv);
  276. bool f_cv_exists =
  277. frame.hs_freq_cv_connected ||
  278. frame.p1_freq_cv_connected ||
  279. frame.p2_freq_cv_connected ||
  280. frame.ls_freq_cv_connected ||
  281. frame.global_freq_cv_connected;
  282. auto q_knob = simd::float_4(
  283. 0.f,
  284. frame.p1_q_knob,
  285. frame.p2_q_knob,
  286. 0.f);
  287. auto q_cv = simd::float_4(
  288. 0.f,
  289. frame.p1_q_cv,
  290. frame.p2_q_cv,
  291. 0.f);
  292. bool q_cv_exists = frame.p1_q_cv_connected || frame.p2_q_cv_connected;
  293. auto gain_knob = simd::float_4(
  294. frame.ls_gain_knob,
  295. frame.p1_gain_knob,
  296. frame.p2_gain_knob,
  297. frame.hs_gain_knob);
  298. auto gain_cv = simd::float_4(
  299. frame.ls_gain_cv,
  300. frame.p1_gain_cv,
  301. frame.p2_gain_cv,
  302. frame.hs_gain_cv);
  303. bool gain_cv_exists =
  304. frame.hs_gain_cv_connected ||
  305. frame.p1_gain_cv_connected ||
  306. frame.p2_gain_cv_connected ||
  307. frame.ls_gain_cv_connected ||
  308. frame.global_gain_cv_connected;
  309. f_cv += frame.global_freq_cv;
  310. gain_cv += frame.global_gain_cv;
  311. // V/oct
  312. simd::float_4 v_oct = f_cv + kFreqKnobVoltage * (f_knob - 1.f);
  313. freq_lpf_.process(v_oct);
  314. v_oct = freq_lpf_.lowpass();
  315. // Q CV
  316. q_cv -= simd::rescale(q_knob,
  317. 0.f, 1.f, kQKnobMinVoltage, kQKnobMaxVoltage);
  318. q_cv *= -kQAmpGain;
  319. q_lpf_.process(q_cv);
  320. q_cv = q_lpf_.lowpass();
  321. // Gain CV
  322. simd::float_4 gain_db =
  323. gain_knob * kGainKnobRange + gain_cv * kGainPerVolt;
  324. // Stuff input into unused element of Q CV vector
  325. q_cv[0] = frame.main_in * (frame.pre_gain ? 0.25f : 0.5f);
  326. float timestep = sample_time_ / oversampling_;
  327. // If a CV input is not connected, we can perform the expensive
  328. // exponential calculation here only once instead of inside the loop,
  329. // since we needn't apply oversampling and anti-aliasing to a low-rate
  330. // UI control.
  331. simd::float_4 f_level;
  332. simd::float_4 q_level;
  333. simd::float_4 gain_level;
  334. if (!f_cv_exists)
  335. {
  336. f_level = FreqVCALevel(v_oct);
  337. }
  338. if (!q_cv_exists)
  339. {
  340. q_level = QVCALevel(q_cv);
  341. }
  342. if (!gain_cv_exists)
  343. {
  344. gain_level = GainVCALevel(gain_db);
  345. }
  346. // Outputs
  347. simd::float_4 out1;
  348. simd::float_4 out2;
  349. bool out2_connected =
  350. frame.p2_hp_out_connected ||
  351. frame.p2_bp_out_connected ||
  352. frame.p2_lp_out_connected;
  353. for (int i = 0; i < oversampling_; i++)
  354. {
  355. // Upsample and apply anti-aliasing filters if needed
  356. if (f_cv_exists)
  357. {
  358. v_oct = up_filter_[0].Process(
  359. (i == 0) ? (v_oct * oversampling_) : 0.f);
  360. f_level = FreqVCALevel(v_oct);
  361. }
  362. // We can't skip this one since it contains the input signal
  363. q_cv = up_filter_[1].Process(
  364. (i == 0) ? (q_cv * oversampling_) : 0.f);
  365. if (q_cv_exists)
  366. {
  367. q_level = QVCALevel(q_cv);
  368. }
  369. if (gain_cv_exists)
  370. {
  371. gain_db = up_filter_[2].Process(
  372. (i == 0) ? (gain_db * oversampling_) : 0.f);
  373. gain_level = GainVCALevel(gain_db);
  374. }
  375. // Unpack input from Q CV vector
  376. simd::float_4 in =
  377. _mm_shuffle_ps(q_cv.v, q_cv.v, _MM_SHUFFLE(0, 0, 0, 0));
  378. // Process VCFs
  379. low_high_.Process(timestep, in, f_level);
  380. float low = low_high_.voltage()[0];
  381. float high = low_high_.voltage()[3];
  382. simd::float_4 mid = mid_.Process(timestep, in, f_level, q_level);
  383. // Calculate output
  384. low *= 1.f - gain_level[0];
  385. mid *= 1.f - gain_level;
  386. high = -high + (high + in[0]) * gain_level[3];
  387. float sum = 2.f * (low + mid[1] + mid[2] + high);
  388. out1 = simd::float_4(sum, mid_.lp()[1], mid_.bp()[1], mid_.hp()[1]);
  389. out1 = simd::clamp(out1, -kClampVoltage, kClampVoltage);
  390. // Pre-downsample anti-alias filtering
  391. out1 = down_filter_[0].Process(out1);
  392. if (out2_connected)
  393. {
  394. out2 = simd::float_4(0.f, mid_.lp()[2], mid_.bp()[2], mid_.hp()[2]);
  395. out2 = simd::clamp(out2, -kClampVoltage, kClampVoltage);
  396. out2 = down_filter_[1].Process(out2);
  397. }
  398. }
  399. frame.main_out = out1[0];
  400. clip_hpf_.process(out1[0]);
  401. float clip = 1.f * (std::abs(clip_hpf_.highpass()) > kClipLEDThreshold);
  402. frame.clip = clip_slew_.process(sample_time_, clip);
  403. frame.p1_lp_out = out1[1];
  404. frame.p1_bp_out = out1[2];
  405. frame.p1_hp_out = out1[3];
  406. if (out2_connected)
  407. {
  408. frame.p2_lp_out = out2[1];
  409. frame.p2_bp_out = out2[2];
  410. frame.p2_hp_out = out2[3];
  411. }
  412. }
  413. protected:
  414. float sample_time_;
  415. int oversampling_;
  416. UpsamplingAAFilter<simd::float_4> up_filter_[3];
  417. DownsamplingAAFilter<simd::float_4> down_filter_[2];
  418. LPFilter<simd::float_4> low_high_;
  419. SVFilter<simd::float_4> mid_;
  420. dsp::TRCFilter<simd::float_4> freq_lpf_;
  421. dsp::TRCFilter<simd::float_4> q_lpf_;
  422. dsp::TRCFilter<float> clip_hpf_;
  423. dsp::SlewLimiter clip_slew_;
  424. template <typename T>
  425. T FreqVCALevel(T v_oct)
  426. {
  427. v_oct = simd::clamp(v_oct, kMinVOct, 0.f);
  428. return simd::pow(2.f, v_oct);
  429. }
  430. template <typename T>
  431. T QVCALevel(T q_cv)
  432. {
  433. q_cv = simd::clamp(q_cv, 0.f, kClampVoltage);
  434. return simd::pow(10.f, q_cv / kVCAGainConstant / 20.f);
  435. }
  436. template <typename T>
  437. T GainVCALevel(T gain_db)
  438. {
  439. gain_db = simd::fmin(gain_db, kMaximumGain);
  440. return simd::pow(10.f, gain_db / 20.f);
  441. }
  442. };
  443. }