;; An implementation of parametric curves.
;; 
;; This code is inspired by the Haskell rsagl library.
;;
;; Copyright 2012 Ivan Raikov and the Okinawa Institute of
;; Science and Technology.
;;
;; This program is free software: you can redistribute it and/or
;; modify it under the terms of the GNU General Public License as
;; published by the Free Software Foundation, either version 3 of the
;; License, or (at your option) any later version.
;;
;; This program is distributed in the hope that it will be useful, but
;; WITHOUT ANY WARRANTY; without even the implied warranty of
;; MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
;; General Public License for more details.
;;
;; A full copy of the GPL license can be found at
;; <http://www.gnu.org/licenses/>.
;;

(module parametric-curve

	(simple-curve linear-curve 
	 translate-curve scale-curve
         sample-curve inverse-sample-curve
	 sample-curve* inverse-sample-curve*
	 )
	
	(import scheme chicken data-structures)
  
	(require-library srfi-1 srfi-4 extras bvsp-spline)

	(import (only srfi-1 fold every list-tabulate)
		(only srfi-4 f64vector list->f64vector f64vector-ref)
		(only extras fprintf pp)
		(prefix bvsp-spline bvsp-spline:)
		)

;; A boundary-value preserving spline 
(define-record-type spline (make-spline n k x y d d2)
  spline? 
  (n spline-n)
  (k spline-k)
  (x spline-x)
  (y spline-y)
  (d spline-d)
  (d2 spline-d2)
  )


;; A parametric equation with an associated interval, interpolating
;; spline, inverse interpolating spline.
(define-record-type peq (make-peq fn xmin xmax ymin ymax spl invspl)
  peq? 
  (fn      peq-fn )
  (xmin    peq-xmin)
  (xmax    peq-xmax)
  (ymin    peq-ymin)
  (ymax    peq-ymax)
  (spl     peq-spl)
  (invspl  peq-invspl)
  )


;; A N-dimensional parametric curve is a set of N parametric equations.
(define (pcurve? x) (every pcurve? x))


;; Defines a simple parametric equation.
(define (simple-peq n k fn xmin xmax)
  (if (> xmin xmax)
      (simple-peq n k fn xmax xmin)
      (let ((dx (/ (- xmax xmin) n)))
	(let* ((x (list-tabulate (+ 1 n) (lambda (i) (+ xmin (* i dx)))))
	       (y (map fn x))
	       (ymin (fold min +inf.0 y))
	       (ymax (fold max -inf.0 y))
	       (xv (list->f64vector x))
	       (yv (list->f64vector y)))
	  (let-values (((d d2 constr errc diagn)
			(bvsp-spline:compute n k xv yv)))
	    (let-values (((inv-d inv-d2 constr errc diagn)
			  (bvsp-spline:compute n k yv xv)))
	      (make-peq fn xmin xmax ymin ymax
			(make-spline n k xv yv d d2) 
			(make-spline n k yv xv inv-d inv-d2) )
	      ))
	  ))
      ))


;; Samples a parametric equation at the given point of interest.
(define (sample-peq c)
  (let ((s (peq-spl c))
	(min (peq-xmin c))
	(max (peq-xmax c)))
    (lambda (xp)
      (and (>= xp min) (<= xp max)
	   (let ((v (bvsp-spline:evaluate (spline-n s) (spline-k s) 
					  (spline-x s) (spline-y s) 
					  (spline-d s) (spline-d2 s) 
					  (f64vector xp) 0)))
	     (and v (f64vector-ref v 0)))
	   ))
    ))


;; Inverse sampling: y -> x
(define (inverse-sample-peq c)
  (let ((s (peq-invspl c))
	(min (peq-ymin c))
	(max (peq-ymax c)))
    (lambda (yp)
      (and (>= yp min) (<= yp max)
	   (let ((v (bvsp-spline:evaluate (spline-n s) (spline-k s) 
					  (spline-x s) (spline-y s) 
					  (spline-d s) (spline-d2 s) 
					  (f64vector yp) 0)))
	     (and v (f64vector-ref v 0)))
	   ))
    ))


;; Samples a parametric equation at the given points of interest.
(define (sample-peq* c)
  (let ((s (peq-spl c)))
    (lambda (xps)
      (let ((xpsv (if (list? xps) (list->f64vector xps) xps)))
	(bvsp-spline:evaluate (spline-n s) (spline-k s) 
			      (spline-x s) (spline-y s) 
			      (spline-d s) (spline-d2 s) 
			      xpsv 0))
      )))


;; Inverse sampling: [y] -> [x]
(define (inverse-sample-peq* c)
  (let ((s (peq-invspl c)))
    (lambda (yps)
      (let ((ypsv (if (list? yps) (list->f64vector yps) yps)))
	(bvsp-spline:evaluate (spline-n s) (spline-k s) 
			      (spline-x s) (spline-y s) 
			      (spline-d s) (spline-d2 s) 
			      ypsv 0))
      )))


;; Samples a parametric equation at regular intervals in the range 0..1 inclusive.
(define (iterate-peq c n)
  (let ((f  (sample-peq c))
	(dx (/ 1. n)))
    (list-tabulate (+ n 1) (lambda (i) (f (* dx i))))
    ))
  

;; Transforms an equation using the given function.
(define (map-peq f x)
  (let ((fx (sample-peq x)) (splx (peq-spl x)))
    (let ((f1 (lambda (u) (f (fx u)))))
      (simple-peq (spline-n splx) (spline-k splx) f1 (peq-xmin x) (peq-xmax x))
      )))


;; Combines two equations using the given function.
(define (zip-peq f x y)
  (let ((fx   (sample-peq x)) 
	(splx (peq-spl x)) 
	(fy   (sample-peq y)))
    (let ((f1 (lambda (u) (f (fx u) (fy u)))))
      (simple-peq (spline-n splx) (spline-k splx) f1 (peq-xmin x) (peq-xmax x))
      )))


;; Translates a parametric equation.
(define (translate-peq x c) (map-peq (lambda (v) (+ x v)) c))


;; Scales a parametric equation.  
(define (scale-peq s c) (map-peq (lambda (v) (* s v)) c))


;; Defines a simple parametric curve.
(define (simple-curve n k fs tmin tmax)
  (if (null? fs) '()
      (let ((c (simple-peq n k (car fs) tmin tmax)))
	(cons c (simple-curve n k (cdr fs) tmin tmax))
	)))
  
  
;; Samples a curve at the given point
(define (sample-curve s)
  (let ((scs (map sample-peq s)))
    (lambda (t) (map (lambda (sc) (sc t)) scs))
    ))


;; Inverse curve sample
(define (inverse-sample-curve s)
  (let ((scs (map inverse-sample-peq s)))
    (lambda (t) (map (lambda (sc) (sc t)) scs))
    ))


;; Samples a curve at the given points
(define (sample-curve* s)
  (let ((scs (map sample-peq* s)))
    (lambda (ts) (map (lambda (sc) (sc ts)) scs))
    ))


;; Inverse curve sample
(define (inverse-sample-curve* s)
  (let ((scs (map sample-peq* s)))
    (lambda (ts) (map (lambda (sc) (sc ts)) scs))
    ))


;; Linear curve of the form c1 * x + c2
;; Argument coeffs supplies c1 and c2 for the different dimensions
(define (linear-curve n k coeffs tmin tmax )
  (simple-curve n k (map (lambda (s) (lambda (x) (+ (* (car s) x) (cadr s)))) coeffs) tmin tmax))


;; Translates a parametric curve.
(define (translate-curve xs c)
  (map (lambda (x p) (translate-peq x p)) xs c))

			 
;; Scales a parametric curve.
(define (scale-curve xs c)
  (map (lambda (x p) (scale-peq x p)) xs c))
			 
		
  

)