;;; -*- syntax: common-lisp; base: 10; mode: lisp -*-

;;; Some probability distributions

;;; Uniform

(defun uniform ()
  (random 1.0))

;;; Linear

(defun linear ()
  (let* ((ua (random 1.0))
         (ub (random 1.0)))
    (if (< ua ub) ua ub)))

;;; Triangular

(defun triangular ()
  (let* ((ua (random 1.0))
         (ub (random 1.0)))
    (* 0.5 (+ ua ub))))

;;; Exponential
;;; delta controls the horizontal scale of the density function

(defun exponential (delta)
  (/ (- (log (random 1.0))) delta))

;;; Bilateral exponential

(defun biexp (delta)
  (let* ((val (random 2.0)))
    (if (> val 1.0)
        (/ (- (log val)) delta)
      (- (/ (- (log (- 2.0 val))) delta)))))
  
;;; Gaussian 
;;; mu = center of the distribution (between 0 and 1)
;;; sigma = standard deviation (% of span) 68.26% of results fall into
;;;         mu+-sigma

(defun gaussian (sigma mu &optional (n 12))
  (let* ((novr (/ n 2))
         (scale (/ (sqrt (/ n 12))))
         (sum (loop repeat n sum (random 1.0))))
    (+ (* sigma scale (- sum novr)) mu)))
    
;;; a function we can use to plot the distribution of the
;;; probability function, first argument is number of slots
;;; in which we split the 0..1 range, second argument is the
;;; number of trials (the bigger, the more we average the 
;;; random variations) and the last is the distribution
;;; we want to plot with an optional argument list. A key
;;; argument lets you select the range over which values 
;;; are counted (compile and load "plot.lisp" before
;;; trying to run this example)

(defun plot-distribution (slots times dist &optional (dist-args nil)
                          &key (range '(0.0 1.0)))
  (let* ((limits (make-array (+ slots 1) :initial-contents
                             (loop 
                               for l from (first range) 
                               by (* (- (second range)(first range))(/ slots))
                               repeat (+ slots 1) collect l)))
         (bins (make-array slots :initial-element 0)))
    (loop repeat times
      for each = (apply dist dist-args)
      do 
      (loop 
        for i from 0 below slots
        for low = (aref limits i)
        for high = (aref limits (+ i 1))
        do 
        (if (< low each high)
            (progn
              (incf (aref bins i))
              (loop-finish)))))
    (loop
      for i from 0
      for bin across bins
      collect (/ bin times 1.0) into data
      finally (plot-data data))))


;;; Silly example

(with-sound()
  (loop
    repeat 60
    for time from 0 by 0.1
    for freq = 440.0 then (+ 440 (- (* (triangular) 100) 50))
    do
    (fm-violin time 0.17 freq 0.1)))

(with-sound()
  (loop
    repeat 60
    for time from 0 by 0.1
    for freq = 440.0 then (+ 440 (- (* (exponential 10) 100) 50))
    do
    (fm-violin time 0.17 freq 0.1)))



(defparameter distribs (new cycle :of '(uniform triangular linear)))

(with-sound()
  (let* ((time 0))
    (loop
      repeat 4
      for dist = (next distribs)
      do 
      (format t "~s...~%" dist)
      (fm-violin time 0.1 3000 0.1)
      (loop
        repeat 70
        for freq = 440.0 then (+ 440 (- (* (funcall dist) 400) 200))
        do
        (fm-violin time 0.1 freq 0.1)
        (incf time 0.02)))))