MicroDexed is a compatible 6-operator-FM-synth based on the Teensy(-3.6/-4.0) Microcontroller. https://www.parasitstudio.de
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.

355 lines
11KB

  1. /*
  2. * Copyright (C) 2015-2017 Pascal Gauthier.
  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, write to the Free Software Foundation,
  16. * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
  17. *
  18. * The code is based on ppplay https://github.com/stohrendorf/ppplay and opl3
  19. * math documentation :
  20. * https://github.com/gtaylormb/opl3_fpga/blob/master/docs/opl3math/opl3math.pdf
  21. *
  22. */
  23. #include "EngineMkI.h"
  24. #define _USE_MATH_DEFINES
  25. #include <cmath>
  26. #include <cstdlib>
  27. #include "sin.h"
  28. #include "exp2.h"
  29. #ifdef DEBUG
  30. #include "time.h"
  31. //#define MKIDEBUG
  32. #endif
  33. #ifdef _WIN32
  34. #if _MSC_VER < 1800
  35. FRAC_NUM log2(FRAC_NUM n) {
  36. //return log(n) / log(2.0);
  37. return LOG_FUNC(n) / LOG_FUNC(2.0);
  38. }
  39. FRAC_NUM round(FRAC_NUM n) {
  40. return n < 0.0 ? ceil(n - 0.5) : floor(n + 0.5);
  41. }
  42. #endif
  43. __declspec(align(16)) const int zeros[N] = {0};
  44. #else
  45. const int32_t __attribute__ ((aligned(16))) zeros[_N_] = {0};
  46. #endif
  47. static const uint16_t NEGATIVE_BIT = 0x8000;
  48. static const uint16_t ENV_BITDEPTH = 14;
  49. static const uint16_t SINLOG_BITDEPTH = 10;
  50. static const uint16_t SINLOG_TABLESIZE = 1<<SINLOG_BITDEPTH;
  51. static uint16_t sinLogTable[SINLOG_TABLESIZE];
  52. static const uint16_t SINEXP_BITDEPTH = 10;
  53. static const uint16_t SINEXP_TABLESIZE = 1<<SINEXP_BITDEPTH;
  54. static uint16_t sinExpTable[SINEXP_TABLESIZE];
  55. const uint16_t ENV_MAX = 1<<ENV_BITDEPTH;
  56. static inline uint16_t sinLog(uint16_t phi) {
  57. const uint16_t SINLOG_TABLEFILTER = SINLOG_TABLESIZE-1;
  58. const uint16_t index = (phi & SINLOG_TABLEFILTER);
  59. switch( ( phi & (SINLOG_TABLESIZE * 3) ) ) {
  60. case 0:
  61. return sinLogTable[index];
  62. case SINLOG_TABLESIZE:
  63. return sinLogTable[index ^ SINLOG_TABLEFILTER];
  64. case SINLOG_TABLESIZE * 2 :
  65. return sinLogTable[index] | NEGATIVE_BIT;
  66. default:
  67. return sinLogTable[index ^ SINLOG_TABLEFILTER] | NEGATIVE_BIT;
  68. }
  69. }
  70. EngineMkI::EngineMkI() {
  71. float bitReso = SINLOG_TABLESIZE;
  72. for(int i=0;i<SINLOG_TABLESIZE;i++) {
  73. //float x1 = sin(((0.5+i)/bitReso) * M_PI/2.0);
  74. float x1 = SIN_FUNC(((0.5+i)/bitReso) * M_PI/2.0);
  75. sinLogTable[i] = round(-1024 * log2(x1));
  76. }
  77. bitReso = SINEXP_TABLESIZE;
  78. for(int i=0;i<SINEXP_TABLESIZE;i++) {
  79. float x1 = (pow(2, float(i)/bitReso)-1) * 4096;
  80. sinExpTable[i] = round(x1);
  81. }
  82. #ifdef MKIDEBUG
  83. char buffer[4096];
  84. int pos = 0;
  85. for(int i=0;i<SINLOG_TABLESIZE;i++) {
  86. pos += sprintf(buffer+pos, "%d ", sinLogTable[i]);
  87. if ( pos > 90 ) {
  88. buffer[0] = 0;
  89. pos = 0;
  90. }
  91. }
  92. buffer[0] = 0;
  93. pos = 0;
  94. for(int i=0;i<SINEXP_TABLESIZE;i++) {
  95. pos += sprintf(buffer+pos, "%d ", sinExpTable[i]);
  96. if ( pos > 90 ) {
  97. buffer[0] = 0;
  98. pos = 0;
  99. }
  100. }
  101. #endif
  102. }
  103. inline int32_t mkiSin(int32_t phase, uint16_t env) {
  104. uint16_t expVal = sinLog(phase >> (22 - SINLOG_BITDEPTH)) + (env);
  105. //int16_t expValShow = expVal;
  106. const bool isSigned = expVal & NEGATIVE_BIT;
  107. expVal &= ~NEGATIVE_BIT;
  108. const uint16_t SINEXP_FILTER = 0x3FF;
  109. uint16_t result = 4096 + sinExpTable[( expVal & SINEXP_FILTER ) ^ SINEXP_FILTER];
  110. //uint16_t resultB4 = result;
  111. result >>= ( expVal >> 10 ); // exp
  112. #ifdef MKIDEBUG
  113. if ( ( time(NULL) % 5 ) == 0 ) {
  114. if ( expValShow < 0 ) {
  115. expValShow = (expValShow + 0x7FFF) * -1;
  116. }
  117. }
  118. #endif
  119. if( isSigned )
  120. return (-result - 1) << 13;
  121. else
  122. return result << 13;
  123. }
  124. void EngineMkI::compute(int32_t *output, const int32_t *input,
  125. int32_t phase0, int32_t freq,
  126. int32_t gain1, int32_t gain2, bool add) {
  127. int32_t dgain = (gain2 - gain1 + (_N_ >> 1)) >> LG_N;
  128. int32_t gain = gain1;
  129. int32_t phase = phase0;
  130. const int32_t *adder = add ? output : zeros;
  131. for (int i = 0; i < _N_; i++) {
  132. gain += dgain;
  133. int32_t y = mkiSin((phase+input[i]), gain);
  134. output[i] = y + adder[i];
  135. phase += freq;
  136. }
  137. }
  138. void EngineMkI::compute_pure(int32_t *output, int32_t phase0, int32_t freq,
  139. int32_t gain1, int32_t gain2, bool add) {
  140. int32_t dgain = (gain2 - gain1 + (_N_ >> 1)) >> LG_N;
  141. int32_t gain = gain1;
  142. int32_t phase = phase0;
  143. const int32_t *adder = add ? output : zeros;
  144. for (int i = 0; i < _N_; i++) {
  145. gain += dgain;
  146. int32_t y = mkiSin(phase , gain);
  147. output[i] = y + adder[i];
  148. phase += freq;
  149. }
  150. }
  151. void EngineMkI::compute_fb(int32_t *output, int32_t phase0, int32_t freq,
  152. int32_t gain1, int32_t gain2,
  153. int32_t *fb_buf, int fb_shift, bool add) {
  154. int32_t dgain = (gain2 - gain1 + (_N_ >> 1)) >> LG_N;
  155. int32_t gain = gain1;
  156. int32_t phase = phase0;
  157. const int32_t *adder = add ? output : zeros;
  158. int32_t y0 = fb_buf[0];
  159. int32_t y = fb_buf[1];
  160. for (int i = 0; i < _N_; i++) {
  161. gain += dgain;
  162. int32_t scaled_fb = (y0 + y) >> (fb_shift + 1);
  163. y0 = y;
  164. y = mkiSin((phase+scaled_fb), gain);
  165. output[i] = y + adder[i];
  166. phase += freq;
  167. }
  168. fb_buf[0] = y0;
  169. fb_buf[1] = y;
  170. }
  171. // exclusively used for ALGO 6 with feedback
  172. void EngineMkI::compute_fb2(int32_t *output, FmOpParams *parms, int32_t gain01, int32_t gain02, int32_t *fb_buf, int fb_shift) {
  173. int32_t dgain[2];
  174. int32_t gain[2];
  175. int32_t phase[2];
  176. int32_t y0 = fb_buf[0];
  177. int32_t y = fb_buf[1];
  178. phase[0] = parms[0].phase;
  179. phase[1] = parms[1].phase;
  180. parms[1].gain_out = (ENV_MAX-(parms[1].level_in >> (28-ENV_BITDEPTH)));
  181. gain[0] = gain01;
  182. gain[1] = parms[1].gain_out == 0 ? (ENV_MAX-1) : parms[1].gain_out;
  183. dgain[0] = (gain02 - gain01 + (_N_ >> 1)) >> LG_N;
  184. dgain[1] = (parms[1].gain_out - (parms[1].gain_out == 0 ? (ENV_MAX-1) : parms[1].gain_out));
  185. for (int i = 0; i < _N_; i++) {
  186. int32_t scaled_fb = (y0 + y) >> (fb_shift + 1);
  187. // op 0
  188. gain[0] += dgain[0];
  189. y0 = y;
  190. y = mkiSin(phase[0]+scaled_fb, gain[0]);
  191. phase[0] += parms[0].freq;
  192. // op 1
  193. gain[1] += dgain[1];
  194. y = mkiSin(phase[1]+y, gain[1]);
  195. phase[1] += parms[1].freq;
  196. output[i] = y;
  197. }
  198. fb_buf[0] = y0;
  199. fb_buf[1] = y;
  200. }
  201. // exclusively used for ALGO 4 with feedback
  202. void EngineMkI::compute_fb3(int32_t *output, FmOpParams *parms, int32_t gain01, int32_t gain02, int32_t *fb_buf, int fb_shift) {
  203. int32_t dgain[3];
  204. int32_t gain[3];
  205. int32_t phase[3];
  206. int32_t y0 = fb_buf[0];
  207. int32_t y = fb_buf[1];
  208. phase[0] = parms[0].phase;
  209. phase[1] = parms[1].phase;
  210. phase[2] = parms[2].phase;
  211. parms[1].gain_out = (ENV_MAX-(parms[1].level_in >> (28-ENV_BITDEPTH)));
  212. parms[2].gain_out = (ENV_MAX-(parms[2].level_in >> (28-ENV_BITDEPTH)));
  213. gain[0] = gain01;
  214. gain[1] = parms[1].gain_out == 0 ? (ENV_MAX-1) : parms[1].gain_out;
  215. gain[2] = parms[2].gain_out == 0 ? (ENV_MAX-1) : parms[2].gain_out;
  216. dgain[0] = (gain02 - gain01 + (_N_ >> 1)) >> LG_N;
  217. dgain[1] = (parms[1].gain_out - (parms[1].gain_out == 0 ? (ENV_MAX-1) : parms[1].gain_out));
  218. dgain[2] = (parms[2].gain_out - (parms[2].gain_out == 0 ? (ENV_MAX-1) : parms[2].gain_out));
  219. for (int i = 0; i < _N_; i++) {
  220. int32_t scaled_fb = (y0 + y) >> (fb_shift + 1);
  221. // op 0
  222. gain[0] += dgain[0];
  223. y0 = y;
  224. y = mkiSin(phase[0]+scaled_fb, gain[0]);
  225. phase[0] += parms[0].freq;
  226. // op 1
  227. gain[1] += dgain[1];
  228. y = mkiSin(phase[1]+y, gain[1]);
  229. phase[1] += parms[1].freq;
  230. // op 2
  231. gain[2] += dgain[2];
  232. y = mkiSin(phase[2]+y, gain[2]);
  233. phase[2] += parms[2].freq;
  234. output[i] = y;
  235. }
  236. fb_buf[0] = y0;
  237. fb_buf[1] = y;
  238. }
  239. void EngineMkI::render(int32_t *output, FmOpParams *params, int algorithm, int32_t *fb_buf, int32_t feedback_shift) {
  240. const int kLevelThresh = ENV_MAX-100;
  241. FmAlgorithm alg = algorithms[algorithm];
  242. bool has_contents[3] = { true, false, false };
  243. bool fb_on = feedback_shift < 16;
  244. switch(algorithm) {
  245. case 3 : case 5 :
  246. if ( fb_on )
  247. alg.ops[0] = 0xc4;
  248. }
  249. for (int op = 0; op < 6; op++) {
  250. int flags = alg.ops[op];
  251. bool add = (flags & OUT_BUS_ADD) != 0;
  252. FmOpParams &param = params[op];
  253. int inbus = (flags >> 4) & 3;
  254. int outbus = flags & 3;
  255. int32_t *outptr = (outbus == 0) ? output : buf_[outbus - 1].get();
  256. int32_t gain1 = param.gain_out == 0 ? (ENV_MAX-1) : param.gain_out;
  257. int32_t gain2 = ENV_MAX-(param.level_in >> (28-ENV_BITDEPTH));
  258. param.gain_out = gain2;
  259. if (gain1 <= kLevelThresh || gain2 <= kLevelThresh) {
  260. if (!has_contents[outbus]) {
  261. add = false;
  262. }
  263. if (inbus == 0 || !has_contents[inbus]) {
  264. // PG: this is my 'dirty' implementation of FB for 2 and 3 operators...
  265. if ((flags & 0xc0) == 0xc0 && fb_on) {
  266. switch ( algorithm ) {
  267. // three operator feedback, process exception for ALGO 4
  268. case 3 :
  269. compute_fb3(outptr, params, gain1, gain2, fb_buf, min((feedback_shift+2), 16));
  270. params[1].phase += params[1].freq << LG_N; // hack, we already processed op-5 - op-4
  271. params[2].phase += params[2].freq << LG_N; // yuk yuk
  272. op += 2; // ignore the 2 other operators
  273. break;
  274. // two operator feedback, process exception for ALGO 6
  275. case 5 :
  276. compute_fb2(outptr, params, gain1, gain2, fb_buf, min((feedback_shift+2), 16));
  277. params[1].phase += params[1].freq << LG_N; // yuk, hack, we already processed op-5
  278. op++; // ignore next operator;
  279. break;
  280. default:
  281. // one operator feedback, normal proces
  282. compute_fb(outptr, param.phase, param.freq, gain1, gain2, fb_buf, feedback_shift, add);
  283. break;
  284. }
  285. } else {
  286. compute_pure(outptr, param.phase, param.freq, gain1, gain2, add);
  287. }
  288. } else {
  289. compute(outptr, buf_[inbus - 1].get(), param.phase, param.freq, gain1, gain2, add);
  290. }
  291. has_contents[outbus] = true;
  292. } else if (!add) {
  293. has_contents[outbus] = false;
  294. }
  295. param.phase += param.freq << LG_N;
  296. }
  297. }