;nyquist plug-in
;version 4
;type generate
;categories "http://lv2plug.in/ns/lv2core#GeneratorPlugin"
;preview linear
;name "Additive Tone..."
;action "Generating Tones..."
;author "dm"
;copyright "Released under terms of the GNU General Public License version 2"
;release 1.1

;control coefs-string "Decay Coefficients" string "" "(2.25 3 2.0 5 3) (3.3 2.7) (5.4 2 7.5 (3.2 5) (4 2)) (7.8 9.2) (13.3)"
;control pch "Pitch (Steps)" float-text "" 60 1 96
;control vib-rate "Vibrato Speed (Hz)" float-text "" 36 0 1102.4
;control vib-depth "Vibrato Depth (Hz)" float-text "" 3 0 nil
;control dur "Duration" time "" 3 nil nil

;;;WHAT IT DO

;Generates a tone with pitch pch and duration dur, with vibrato. 
;coefs-string determines the timbre, by specifying how quickly harmonics should ;attenuate according to their frequency.
;coefs-string consists of terms delimited by parentheses (), whose elements may be ;either positive numbers or further sub-terms delimited by parentheses.
;If there are n terms, the first term gives the timbre at the start time, the second at ;dur/n, the third at 2*dur/n, and so on. Timbre is smoothly interpolated between these ;points.
;In a term with k elements, the rth element determines the amplitude of 
;the r-1 mod k harmonics.
;If the rth element is a number x, then the ith harmonic (where i is r-1 modulo k) has ;amplitude i^(-x). That is, higher numbers make the harmonics attenuate faster.
;If the rth element is a subterm with l elements, then the sth element of that subterm
;determines the amplitude of the (r+k*(s-1)) mod k*l harmonics.

;As an example, here are the harmonics and relative amplitudes when
;coefs-string is "((5 (4 2)) 3)" and pch is 69:
;Harmonic    Frequency    Amplitude    Amplitude
;1           440          1^-5         1.000000
;2           880          2^-3         0.125000
;3           1320         3^-4         0.012345
;4           1760         4^-3         0.015625
;5           2200         5^-5         0.000320
;6           2640         6^-3         0.004630
;7           3080         7^-2         0.020408
;8           3520         8^-3         0.001953

;;;ACKNOWLEDGEMENT

;Shout out to steve for a very thorough code review and loads of coding tips, and the ;rewritten gen-tables which I nicked verbatim from a post of his.
;Thank you David R. Sky for the handy string-to-list function.

;;;ADDITIONAL GLOBAL VARIABLES

;;Don't compute harmonics if their frequencies are too high.
(setq hi-freq (min (/ *sound-srate* 2) 17000.0))
(setq top-overtone (min
  1000    ;magic number to avoid argument stack overflow when calling simrep.
  (+ 1 (truncate (/ hi-freq (step-to-hz pch))))))

(defun get-coefficient (overtone lst divis)
  ;;;recursively search coefficient list for exponent
  (let* ((len (length lst))
         (branches (mult divis len))
         (term (nth (rem (/ overtone divis) len) lst)))
    (if (numberp term)
        term
        (get-coefficient overtone term branches))))

(defun get-amplitude (overtone coef-list)
  (expt (+ 1.0 overtone) 
    (- (get-coefficient overtone coef-list 1))))

;2022-07-06, I fucked up my backup and started from this old version
;which still had stack overflow errors from simrep.
;(defun gen-table (coef-list)
;  (simrep (current-overtone top-overtone)
;    (scale (get-amplitude current-overtone coef-list)
;      (build-harmonic (+ 1 current-overtone) 2048))))
;2022-07-06, this one should work better.
(defun gen-table (coef-list)
  (do ((table 0)
       (current-overtone 0))
   ((> current-overtone top-overtone) table)
   (setq current-overtone (+ 1 current-overtone))
   (setq table
      (sum table (scale (get-amplitude current-overtone coef-list)
        (build-harmonic (+ 1 current-overtone) 2048))))))

(defun normalize (sig)
  ;;;assume signal has 2048 samples
  (scale (/ 0.75 (peak sig 2048)) sig))

(defun gen-tables (coef-lists pitch interval)
  ;;; Push wavetable and breakpoint onto a list 'tables'.
  ;;; Discard the final breakpoint, then return 'tables' for siosc.
  (let ((count 1) tables)
    (dolist (coef-list coef-lists (reverse (cdr tables)))
      (push (normalize (gen-table coef-list)) tables)
      (push (* interval count) tables)
      (incf count))))

(defun modulation (vib-rate vib-depth duration)
  (if (= vib-rate 0)
      (s-rest duration)
      (scale vib-depth (lfo vib-rate duration))))

(defun gen-tone (pitch duration coef-lists)
  ;;;get wavetables and mod. env., pass to built-in Nyquist oscillators.
  (validate-coefs coef-lists)
  (let* ((mod (modulation vib-rate vib-depth duration))
         (interval (/ duration (length coef-lists)))
         (tables (gen-tables coef-lists pitch interval)))
    (if (= (length coef-lists) 1)
        (fmosc pitch mod (maketable (car tables)) 0)
        (siosc pitch mod tables))))

(defun validate-coefs (lst)
  ;;; All items in top level list must be lists.
  (dolist (coef lst)
    (when (not (listp coef))
      (throw 'err (format nil "Coefficient ~s is invalid (not a list)." coef)))
    (validate-items lst)))

(defun validate-items (item)
  ;;; Items must be lists or numbers. Test recursively.
  (cond
    ((and (numberp item) (>= item 0)))  ;valid
    ((listp item) (mapc 'validate-items item))  ;recurse
    (t (throw 'err (format nil "Coefficient ~s is invalid." item)))))

;;;Took this one-liner from David R. Sky's Sequencer 2, released under GPLv2.
;;;Replace this with eval-string when 2.3.1 releases.
(defun string-to-list (string)
  (read (make-string-input-stream (format nil "(~a)" string))))

;;read string, then generate tone (or throw error)
(setf coefs (string-to-list coefs-string))
(catch 'err (gen-tone pch dur coefs))