(require 'asdf)
(mapc 'require '(#:alexandria #:bordeaux-fft))

(defun wave-1 (x)
(* (sin (* x 2 pi (/ 200)))
  (sin (* x 2 pi (/ 17)))))

(defun wave-2 (x)
(* (sin (* x 2 pi (/ 300)))
  (sin (* x 2 pi (/ 23)))))

(defun wave-3 (x)
(* (sin (* x 2 pi (/ 250)))
  (sin (* x 2 pi (/ 37)))))

(defvar *array*
(make-array 4096 :initial-contents
           (loop for x below 4096 collecting 0.0)))

(loop for x below 4096 do
(setf (aref *array* x)
     (+ (wave-3 x)
       (cond ((< 100 x 300) (wave-1 x))
            ((< 3700 x 4000) (wave-2 x))
            (t 0)))))

(defvar *kernel*
(make-array 4096 :initial-contents
           (loop for x below 4096 collecting ; just change to wave-2
            (if (< x 200) (wave-1 x) 0.0)))) ; for 2nd kernel

(defvar *convolution-1*
(make-array 4096 :initial-contents
           (loop for x below 4096 collecting #C(0.0 0.0))))

(let ((farray (bordeaux-fft:sfft *array*))
     (fkernel (bordeaux-fft:sfft *kernel*))
     (fconv   (make-array 4096))
     (auto    (make-array 4096)))
(loop for a across farray
     for k across fkernel
     for x from 0
     do (setf (aref fconv x) (* a k))
     do (setf (aref auto x) (* a a))
     finally (loop for n from 0
              for o across (bordeaux-fft:sifft auto)
              for p across (bordeaux-fft:sifft fconv) do
              (setf (aref *convolution-1* n) p))))

(with-open-file (out "to-graph.txt" :direction :output :if-exists :supersede)
(loop for c across *convolution-1*
 do (print (realpart c) out)
 finally (terpri)))

(ext:quit)