source: project/release/4/parametric-curve/trunk/parametric-curve.scm @ 33041

Last change on this file since 33041 was 33041, checked in by Ivan Raikov, 4 years ago

parametric-curve: added matchable as a test dependency [thanks to mario for reporting]

File size: 9.6 KB
Line 
1
2;; An implementation of parametric curves.
3;;
4;; This code is inspired by the Haskell rsagl library.
5;;
6;; Copyright 2012-2015 Ivan Raikov
7;;
8;; This program is free software: you can redistribute it and/or
9;; modify it under the terms of the GNU General Public License as
10;; published by the Free Software Foundation, either version 3 of the
11;; License, or (at your option) any later version.
12;;
13;; This program is distributed in the hope that it will be useful, but
14;; WITHOUT ANY WARRANTY; without even the implied warranty of
15;; MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
16;; General Public License for more details.
17;;
18;; A full copy of the GPL license can be found at
19;; <http://www.gnu.org/licenses/>.
20;;
21
22(module parametric-curve
23
24        (
25         parametric-curve?
26         simple-curve linear-curve line-segment
27         map-curve translate-curve scale-curve
28         compose-curve
29         sample-curve sample-curve*
30         iterate-curve fold-curve foldi-curve
31         bbox-curve arc-length
32         )
33       
34        (import scheme chicken data-structures)
35 
36        (require-library srfi-1 srfi-4 extras bvsp-spline)
37
38        (import (only srfi-1 fold every list-tabulate zip concatenate)
39                (only srfi-4 f64vector make-f64vector list->f64vector f64vector->list f64vector-length f64vector-ref f64vector-set!
40                      list->u32vector u32vector-ref)
41                (only extras fprintf pp)
42                (prefix bvsp-spline bvsp-spline:)
43                )
44
45;; A boundary-value preserving spline
46(define-record-type spline (make-spline n k x y d d2)
47  spline? 
48  (n spline-n)
49  (k spline-k)
50  (x spline-x)
51  (y spline-y)
52  (d spline-d)
53  (d2 spline-d2)
54  )
55
56
57;; A parametric equation with an associated interval and interpolating spline
58(define-record-type peq (make-peq fn xmin xmax ymin ymax spl)
59  peq? 
60  (fn      peq-fn )
61  (xmin    peq-xmin)
62  (xmax    peq-xmax)
63  (ymin    peq-ymin)
64  (ymax    peq-ymax)
65  (spl     peq-spl)
66  )
67
68
69(define-record-printer (peq x out)
70  (fprintf out "#(peq fn=~A xmin=~A xmax=~A ymin=~A ymax=~A spl=~A)"
71           (peq-fn x)
72           (peq-xmin x)
73           (peq-xmax x)
74           (peq-ymin x)
75           (peq-ymax x)
76           (peq-spl x)
77           ))
78
79
80
81;; A N-dimensional parametric curve is a set of N parametric equations.
82(define (parametric-curve? x) (every peq? x))
83
84
85;; Defines a simple parametric equation.
86(define (simple-peq n k fn xmin xmax)
87
88  (if (> xmin xmax)
89
90      (simple-peq n k fn xmax xmin)
91
92      (let ((dx (/ (- xmax xmin) n)))
93
94        (let* ((x (list-tabulate (+ 1 n) (lambda (i) (+ xmin (* i dx)))))
95               (y (map fn x))
96               (ymin (fold min +inf.0 y))
97               (ymax (fold max -inf.0 y))
98               (xv (list->f64vector x))
99               (yv (list->f64vector y))
100               )
101
102          (let ((spl
103                 (let-values (((d d2 constr errc diagn)
104                               (bvsp-spline:compute n k xv yv)))
105                   (if (zero? errc)
106                       (make-spline n k xv yv d d2) 
107                       (error 'simple-peq "unable to compute interpolating spline"))
108                   ))
109                )
110
111              (make-peq fn xmin xmax ymin ymax spl)
112              ))
113          ))
114      )
115
116
117;; Samples a parametric equation at the given point of interest.
118(define (sample-peq c)
119  (let ((s (peq-spl c))
120        (min (peq-xmin c))
121        (max (peq-xmax c)))
122    (lambda (xp)
123      (and (>= xp min) (<= xp max)
124           (let ((v (bvsp-spline:evaluate (spline-n s) (spline-k s) 
125                                          (spline-x s) (spline-y s) 
126                                          (spline-d s) (spline-d2 s) 
127                                          (f64vector xp) 0)))
128             (and v (f64vector-ref v 0)))
129           ))
130    ))
131
132
133;; Samples a parametric equation at the given points of interest.
134(define (sample-peq* c)
135  (let ((s (peq-spl c)))
136    (lambda (xps)
137      (let ((xpsv (if (list? xps) (list->f64vector xps) xps)))
138        (let ((res (bvsp-spline:evaluate (spline-n s) (spline-k s) 
139                                         (spline-x s) (spline-y s) 
140                                         (spline-d s) (spline-d2 s) 
141                                         xpsv 0)))
142          res
143          ))
144      ))
145  )
146
147
148
149;; Samples a parametric equation at regular intervals in the range xmin..xmax inclusive.
150(define (iterate-peq c n)
151  (let* (
152         (f     (sample-peq* c))
153         (xmin  (peq-xmin c))
154         (xmax  (peq-xmax c))
155         (delta (- xmax xmin))
156         (dx    (if (zero? delta) 0
157                    (if (< n 2)
158                        (error 'iterate-peq "number of iterations must be >= 2")
159                        (/ (- xmax xmin) (- n 1)))))
160         )
161    (f (list-tabulate n (lambda (i) (+ xmin (* dx i)))))
162    ))
163
164
165
166;; Transforms an equation using the given function.
167(define (map-peq f x)
168  (let ((fx (sample-peq x)) 
169        (splx (peq-spl x)))
170    (let ((f1 (lambda (u) (f (fx u)))))
171      (simple-peq (spline-n splx) (spline-k splx) f1 (peq-xmin x) (peq-xmax x))
172      )))
173
174
175
176
177;; Combines two equations using the given function.
178(define (compose-peq f x y)
179  (let ((fx   (sample-peq x)) 
180        (splx (peq-spl x)) 
181        (fy   (sample-peq y)))
182    (let ((f1 (lambda (u) (f (fx u) (fy u)))))
183      (simple-peq (spline-n splx) (spline-k splx) f1 (peq-xmin x) (peq-xmax x))
184      )))
185
186
187;; Translates a parametric equation.
188(define (translate-peq x p) (map-peq (lambda (v) (+ x v)) p))
189
190
191;; Scales a parametric equation. 
192(define (scale-peq s p) (map-peq (lambda (v) (* s v)) p))
193
194
195;; Defines a simple parametric curve.
196(define (simple-curve n k fs tmin tmax)
197  (if (null? fs) '()
198      (let ((c (simple-peq n k (car fs) tmin tmax)))
199        (cons c (simple-curve n k (cdr fs) tmin tmax))
200        )))
201 
202 
203;; Samples a curve at the given point
204(define (sample-curve c)
205  (let ((scs (map sample-peq c)))
206    (lambda (t) (map (lambda (sc) (sc t)) scs))
207    ))
208
209;; Samples a curve at the given points
210(define (sample-curve* c)
211  (let ((scs (map sample-peq* c)))
212    (lambda (ts) (map (lambda (sc) (sc ts)) scs))
213    ))
214
215
216;; Linear curve of the form c1 * x + c2
217;; Argument coeffs supplies c1 and c2 for the different dimensions
218(define (linear-curve n coeffs tmin tmax )
219  (simple-curve n 1 (map (lambda (s) 
220                           (let ((c1 (car s)) (c2 (cadr s)))
221                             (cond ((and (zero? c2) (zero? c1)) (lambda (x) 0.0))
222                                   ((zero? c1) (lambda (x) c2))
223                                   (else (lambda (x) (+ (* c1 x) c2))))))
224                           coeffs)
225                tmin tmax))
226
227
228;; Line segment curve of the form (x1,xn) defined on the parameter range 0.0 .. 1.0
229(define (line-segment n coeffs )
230  (simple-curve n 1 (map (lambda (s) (lambda (x) (* x s))) coeffs) 0.0 1.0))
231
232
233;; Maps the given functions to the parametric curve.
234(define (map-curve fs c)
235  (map (lambda (f p) (map-peq f p)) fs c))
236
237;; Composes the parametric curves using the given functions.
238(define (compose-curve fs c1 c2)
239  (map (lambda (f p1 p2) (compose-peq f p1 p2)) fs c1 c2))
240
241
242;; Translates a parametric curve.
243(define (translate-curve xs c)
244  (map (lambda (x p) (translate-peq x p)) xs c))
245
246                         
247;; Scales a parametric curve.
248(define (scale-curve xs c)
249  (map (lambda (x p) (scale-peq x p)) xs c))
250
251                         
252;; Obtain the bounding box of a curve
253(define (bbox-curve c) (append (map peq-ymin c) (map peq-ymax c)))
254
255
256;; Samples a parametric curve at regular intervals in the range xmin..xmax inclusive.
257(define (iterate-curve c n) (map (lambda (p) (iterate-peq p n)) c))
258
259
260;; Folds a parametric curve at regular intervals in the range xmin..xmax inclusive.
261(define (fold-curve c n f init) 
262  (let* (
263         (gs     (map sample-peq* c))
264         (xmins  (map peq-xmin c))
265         (xmaxs  (map peq-xmax c))
266         (deltas (map - xmaxs xmins))
267         (dxs    (map (lambda (xmin xmax delta)
268                        (if (zero? delta) 0
269                            (if (< n 2)
270                                (error 'fold-curve "number of iterations must be >= 2")
271                                (/ (- xmax xmin) (- n 1)))))
272                      xmins xmaxs deltas))
273         (inds   (map (lambda (xmin dx) 
274                        (list-tabulate n (lambda (i) (+ xmin (* dx i)))))
275                      xmins dxs))
276         (vs     (map (lambda (g ind) (g ind)) gs inds))
277         )
278    (let recur ((i 0) (init init))
279      (if (< i n)
280          (let* ((vsi   (map (lambda (vect) (f64vector-ref vect i)) vs))
281                 (init1 (f vsi init)))
282            (recur (+ 1 i) init1))
283          init))
284    ))
285         
286
287;; Like fold-curve, but F is of the form F(I,V,INIT)
288
289;; Folds a parametric curve at regular intervals in the range xmin..xmax inclusive.
290(define (foldi-curve c n f init) 
291  (let* (
292         (gs     (map sample-peq* c))
293         (xmins  (map peq-xmin c))
294         (xmaxs  (map peq-xmax c))
295         (deltas (map - xmaxs xmins))
296         (dxs    (map (lambda (xmin xmax delta)
297                        (if (zero? delta) 0
298                            (if (< n 2)
299                                (error 'fold-curve "number of iterations must be >= 2")
300                                (/ (- xmax xmin) (- n 1)))))
301                      xmins xmaxs deltas))
302         (inds   (map (lambda (xmin dx) 
303                        (list-tabulate n (lambda (i) (+ xmin (* dx i)))))
304                      xmins dxs))
305         (vs     (map (lambda (g ind) (g ind)) gs inds))
306         )
307    (let recur ((i 0) (init init))
308      (if (< i n)
309          (let* ((vsi   (map (lambda (vect) (f64vector-ref vect i)) vs))
310                 (init1 (f i vsi init)))
311            (recur (+ 1 i) init1))
312          init))
313    ))
314         
315
316
317;; Computes the arc length of the parametric curve given step dx
318(define (arc-length c dx)
319  (let* ((n (inexact->exact (round (/ 1.0 dx))))
320         (v (iterate-curve c n)))
321    (let recur ((i 1) (l 0.0) 
322                (s (map (lambda (x) (f64vector-ref x 0)) v)))
323      (if (< i n)
324          (let ((s1 (map (lambda (x) (f64vector-ref x i)) v)))
325            (recur (+ 1 i)
326                   (+ l (sqrt (fold (lambda (x x1 l) (+ l (expt (- x1 x) 2))) 0.0 s s1)))
327                   s1))
328          l))
329    ))
330           
331   
332 
333
334)
Note: See TracBrowser for help on using the repository browser.