objects.ss

(define-scheme-record object
  ((immutable color)
   (immutable opacity)
   (immutable surface)                  ; Surface shader, cannot be #f
   (immutable volume)                   ; Volume shader | #f
   (immutable displacement)             ; Displacement shader | #f
   (immutable center)
   (immutable M)                        ; Transform matrix
   (immutable Mi)                       ; Inverse transform
   ))

(define-scheme-record sphere object ())
(define-scheme-record plane object ())
(define-scheme-record quadric object ((immutable coefficients)))
(define-scheme-record polyhedron object ((immutable planes)))
(define-scheme-record torus object ((immutable radius2)))

(define-scheme-record csg object ((immutable object1) (immutable object2)))
(define-scheme-record csg-union csg ())
(define-scheme-record csg-intersect csg ())
(define-scheme-record csg-difference csg ())

(define (sort-intersections ls)
  (sort (lambda (x y) (< (<intersect> time x) (<intersect> time y)))
    ls))

(define (object-intersections object ray)
  ((cond
    [(sphere? object) sphere-intersections]
    [(plane? object) plane-intersections]
    [(polyhedron? object) polyhedron-intersections]
    [(quadric? object) quadric-intersections]
    [(torus? object) torus-intersections]
    [(csg-union? object) csg-union-intersections]
    [(csg-intersect? object) csg-intersect-intersections]
    [(csg-difference? object) csg-difference-intersections]
    [else (errorf 'object-intersections "unknown object type: ~s" object)])
   object
   (<ray> make
     [origin (mat-vec-mul (object-Mi object)
               (vec-sub (<ray> origin ray) (object-center object)))]
     [direction (mat-vec-mul (object-Mi object) (<ray> direction ray))])))

(define (object-normal object extra intersect-point)
  ((cond
    [(sphere? object) sphere-normal]
    [(plane? object) plane-normal]
    [(polyhedron? object) polyhedron-normal]
    [(quadric? object) quadric-normal]
    [(torus? object) torus-normal]
    [(csg? object) csg-normal]
    [else (errorf 'object-normal "unknown object type: ~s" object)])
   object extra intersect-point))

(define (object-displace obj intersect ip n)
  (cond
   [(object-displacement obj) =>
    (lambda (shader)
      (let ([st (delay (point->texture obj intersect (point->surface obj ip)))])
        (parameterize ([$object obj]
                       [$P ip]
                       [$Ng n]
                       [$N n]
                       [$s (delay (vec-i (force st)))]
                       [$t (delay (vec-j (force st)))])
          (call-with-values shader
            (case-lambda
             [(normal) (values ip normal)]
             [(intersect-point normal) (values intersect-point normal)])))))]
   [else
    (values ip n)]))

(define (object-shade obj intersect extra ip ng n i d)
  (cond
   [(object-surface obj) =>
    (lambda (shader)
      (let ([st (delay (point->texture obj intersect (point->surface obj ip)))])
        (parameterize ([$object obj]
                       [$depth d]
                       [$P ip]
                       [$Ng ng]
                       [$N n]
                       [$I i]
                       [$s (delay (vec-i (force st)))]
                       [$t (delay (vec-j (force st)))])
          (call-with-values shader
            (case-lambda
             [(color) (values color Os)]
             [(color opacity) (values color opacity)])))))]
   [(and (csg? obj) (<intersect> object extra))
    (object-shade (<intersect> object extra) extra (<intersect> extra extra)
      ip ng n i d)]
   [else
    (errorf 'object-shade "no object shader defined for ~a" obj)]))

(define (volume-shade obj ip i ci oi)
  (cond
   [(object-volume obj) =>
    (lambda (shader)
      (parameterize ([$object obj]
                     [$P ip]
                     [$I i]
                     [$Ci ci]
                     [$Oi oi])
        (call-with-values shader
          (case-lambda
           [(color) (values color (object-opacity obj))]
           [(color opacity) (values color opacity)]))))]
   [else #f]))

(define (point->surface object point)
  ;; Maps intersect point to surface shader coordinate space
  (cond
   [(plane? object) (plane-point->surface object point)]
   [else
    (mat-vec-mul (object-Mi object)
      (vec-sub point (object-center object)))]))

(define (point->texture object intersect point)
  ;; Maps point from surface shader coordinates to texture coordinates
  (cond
   [(sphere? object) (sphere-point->texture object intersect point)]
   [(plane? object) (plane-point->texture object intersect point)]
   [(torus? object) (torus-point->texture object intersect point)]
   [(cube? object) (cube-point->texture object intersect point)]
   [else (errorf 'point->texture "unknown object type: ~s" object)]))

(define (sphere-intersections object ray)
  (let ([origin (<ray> origin ray)]
        [direction (<ray> direction ray)])
    (let* ([a (vec-dot direction direction)]
           [b (vec-dot direction origin)]
           [c (- (vec-dot origin origin) 1)]
           [dis (- (* b b) (* a c))])
      (if (>= dis 0)
          (let ([d (sqrt dis)])
            (list
             (<intersect> make [time (/ (- (- b) d) a)]
               [object object] [extra #f])
             (<intersect> make [time (/ (+ (- b) d) a)]
               [object object] [extra #f])))
          '()))))

(define (sphere-normal object extra intersect-point)
  (vec-sub intersect-point (object-center object)))

(define (sphere-point->texture object intersect point)
  (let* ([N (make-vec 0 0 1)]
         [E (make-vec 1 0 0)]
         [P (vec-normalize point)]
         [phi (acos (- (vec-dot N P)))]
         [v (/ phi pi)]
         [theta (/ (acos (/ (vec-dot P E) (sin phi))) (* 2 pi))]
         [u (if (> (vec-dot (vec-cross N E) P) 0)
                theta
                (- 1 theta))])
    (make-vec u v 0)))

(define (plane-intersections object ray)
  (list (<intersect> make
          [time (- (/ (vec-k (<ray> origin ray))
                      (vec-k (<ray> direction ray))))]
          [object object]
          [extra #f])))

(define (plane-normal object extra intersect-point)
  (mat-vec-mul (object-M object) (make-vec 0 0 1)))

(define (plane-point->surface object point)
  (vec-num-mul
   (vec-num-add
    (mat-vec-mul (object-Mi object)
      (vec-sub point (object-center object)))
    1)
   0.5))

(define (plane-point->texture object intersect point)
  (make-vec (mod (vec-i point) 1) (mod (vec-j point) 1) 0))

(define (polyhedron-intersections obj ray)
  (let ([origin (<ray> origin ray)]
        [direction (<ray> direction ray)])
    (let lp ([planes (polyhedron-planes obj)] [t-in -inf.0] [t-out +inf.0]
             [plane-in #f] [plane-out #f])
      (if (null? planes)
          (if (< t-in t-out)
              (list                     ; found intersections
               (<intersect> make [time t-in] [object obj] [extra plane-in])
               (<intersect> make [time t-out] [object obj] [extra plane-out]))
              '())
          (let ([plane (car planes)])
            (if (>= t-in t-out)
                '()                     ; No hit
                (let ([numer (- 1 (vec-dot plane origin))]
                      [denom (vec-dot plane direction)])
                  (cond
                   [(= denom 0)         ; ray is parallel
                    (if (< numer EPSILON)
                        '()             ; done
                        (lp (cdr planes) t-in t-out plane-in plane-out))]
                   [else                ; ray is not parallel
                    (let ([t-hit (/ numer denom)])
                      (cond
                       [(> denom 0)
                        (if (< t-out t-hit)
                            (lp (cdr planes) t-in t-out plane-in plane-out)
                            (lp (cdr planes) t-in t-hit plane-in plane))]
                       [else
                        (if (> t-in t-hit)
                            (lp (cdr planes) t-in t-out plane-in plane-out)
                            (lp (cdr planes) t-hit t-out plane plane-out))]))]))))))))

(define (polyhedron-normal obj extra intersect-point)
  (mat-vec-mul (object-M obj) extra))

(define (cube-point->texture object intersect point)
  (define (st s t)
    (make-vec (mod s 1) (mod t 1) 0))
  (let ([x (/ (+ (vec-i point) 1) 2)]
        [y (/ (+ (vec-j point) 1) 2)]
        [z (/ (+ (vec-k point) 1) 2)]
        [plane (<intersect> extra intersect)])
    (match plane
      [`(<vec> [i 0] [j 0] [k 1])
       (st (+ (/ x 4) 1/4) (+ (/ (- 1 y) 3) 1/3))] ; forward
      [`(<vec> [i 0] [j 0] [k -1])
       (st (+ (/ (- 1 x) 4) 3/4) (+ (/ (- 1 y) 3) 1/3))] ; back
      [`(<vec> [i 0] [j 1] [k 0])
       (st (+ (/ x 4) 1/4) (+ (/ z 3) 0))] ; up
      [`(<vec> [i 0] [j -1] [k 0])
       (st (+ (/ x 4) 1/4) (+ (/ (- 1 z) 3) 2/3))] ; down
      [`(<vec> [i 1] [j 0] [k 0])
       (st (+ (/ (- 1 z) 4) 2/4) (- (+ (/ y 3) 1/3)))] ; right
      [`(<vec> [i -1] [j 0] [k 0])
       (st (+ (/ z 4) 0) (+ (/ (- 1 y) 3) 1/3))] ; left
      )))

(define (quadric-intersections object ray)
  ;; A x^2 + B y^2 + C z^2 + D xy + E xz + F yz + G x + H y + I z + J = 0
  (let ([origin (<ray> origin ray)]
        [direction (<ray> direction ray)])
    (let ([Xo (vec-i origin)] [Yo (vec-j origin)] [Zo (vec-k origin)]
          [Xd (vec-i direction)] [Yd (vec-j direction)] [Zd (vec-k direction)])
      (match (quadric-coefficients object)
        [#(,A ,E ,H ,B ,C ,F ,D ,G ,I ,J)
         (let ([a (+ (* Xd (+ (* A Xd) (* B Yd) (* C Zd)))
                     (* Yd (+ (* E Yd) (* F Zd)))
                     (* H Zd Zd))]
               [b (+ (* Xd (+ (* A Xo) (* 0.5 (+ (* B Yo) (* C Zo) D))))
                     (* Yd (+ (* E Yo) (* 0.5 (+ (* B Xo) (* F Zo) G))))
                     (* Zd (+ (* H Zo) (* 0.5 (+ (* C Xo) (* F Yo) I)))))]
               [c (+ (* Xo (+ (* A Xo) (* B Yo) (* C Zo) D))
                     (* Yo (+ (* E Yo) (* F Zo) G))
                     (* Zo (+ (* H Zo) I))
                     J)])
           (if (not (= a 0))            ; quadratic
               (let ([dis (- (* b b) (* a c))])
                 (if (>= dis 0)
                     (let ([d (sqrt dis)])
                       (list
                        (<intersect> make [time (/ (- (- b) d) a)]
                          [object object] [extra #f])
                        (<intersect> make [time (/ (+ (- b) d) a)]
                          [object object] [extra #f])))
                     '()))
               (if (= b 0)
                   '()
                   (list (<intersect> make [time (- 0.5 (/ c b))]
                           [object object] [extra #f])))))]))))

(define (quadric-normal object extra intersect-point)
  (let ([ip (mat-vec-mul (object-Mi object)
              (vec-sub intersect-point (object-center object)))])
    (let ([ipx (vec-i ip)] [ipy (vec-j ip)] [ipz (vec-k ip)])
      (match (quadric-coefficients object)
        [#(,A ,E ,H ,B ,C ,F ,D ,G ,I ,J)
         (mat-vec-mul (object-M object)
           (make-vec
            (+ (* 2 A ipx) (* B ipy) (* C ipz) D)
            (+ (* B ipx) (* 2 E ipy) (* F ipz) G)
            (+ (* C ipx) (* F ipy) (* 2 H ipz) I)))]))))

(define (torus-intersections object ray)
  (let ([origin (<ray> origin ray)]
        [direction (<ray> direction ray)])
    (let* ([len (vec-length direction)]
           [direction (vec-normalize direction)]
           [R2 (sqr 1)]                 ; outer radius is one
           [r2 (sqr (torus-radius2 object))]
           [oy2 (sqr (vec-j origin))]
           [dy2 (sqr (vec-j direction))]
           [ody2 (* (vec-j origin) (vec-j direction))]
           [k1 (+ (sqr (vec-i origin)) (sqr (vec-k origin)) oy2 (- R2) (- r2))]
           [k2 (+ (* (vec-i origin) (vec-i direction))
                  (* (vec-k origin) (vec-k direction))
                  ody2)])
      (map (lambda (t)
             (<intersect> make [time (/ t len)] [object object] [extra #f]))
        (solve-quartic
         (+ (sqr k1) (* 4 R2 (- oy2 r2)))
         (* 4 (+ (* k2 k1) (* 2 R2 ody2)))
         (* 2 (+ k1 (* 2 (+ (sqr k2) (* R2 dy2)))))
         (* 4 k2)
         1)))))

(define (torus-normal object extra intersect-point)
  (let ([ip (mat-vec-mul (object-Mi object)
              (vec-sub intersect-point (object-center object)))])
    (let ([d (sqrt (+ (sqr (vec-i ip)) (sqr (vec-k ip))))])
      (mat-vec-mul (object-M object)
        (if (> d EPSILON)
            (vec-sub ip (make-vec (/ (vec-i ip) d) 0 (/ (vec-k ip) d)))
            ip)))))

(define (torus-point->texture object intersect point)
  (let* ([x (vec-i point)]
         [y (vec-j point)]
         [z (vec-k point)]
         [u (- 1 (/ (+ (atan z (- x)) pi) (* 2 pi)))]
         [len (sqrt (+ (sqr x) (sqr z)))]
         [x (- len 1)]
         [v (/ (+ (atan y x) pi) (* 2 pi))])
    (make-vec u v 0)))

(define (csg-union-intersections object ray)
  (map
   (lambda (inter)
     (<intersect> copy inter [object object] [extra inter]))
   (let union ([a (sort-intersections
                   (object-intersections (csg-object1 object) ray))]
               [b (sort-intersections
                   (object-intersections (csg-object2 object) ray))])
     (cond
      [(null? a) b]
      [(null? b) a]
      [else
       (let ([a1 (car a)]
             [a2 (cadr a)]
             [b1 (car b)]
             [b2 (cadr b)])
         (let ([ca1 (<intersect> time a1)]
               [ca2 (<intersect> time a2)]
               [cb1 (<intersect> time b1)]
               [cb2 (<intersect> time b2)])
           (cond
            [(<= ca2 cb1)               ; A is completely less than B
             (append (list a1 a2) (union (cddr a) b))]
            [(<= cb2 ca1)               ; B is completely less than A
             (append (list b1 b2) (union a (cddr b)))]
            [(and (<= ca1 cb1) (>= ca2 cb2)) ; A is completely around B
             (union a (cddr b))]
            [(and (<= cb1 ca1) (>= cb2 ca2)) ; B is completely around A
             (union (cddr a) b)]
            [(and (< ca1 cb1) (> cb2 ca2)) ; A and B cross, A first
             (union (append (list a1 b2) (cddr a)) b)]
            [(and (< cb1 ca1) (> ca2 cb2)) ; A and B cross, B first
             (union a (append (list b1 a2) (cddr b)))]
            [else (errorf 'union "Unhandled condition A=~s, B=~s" a b)])))]))))

(define (csg-intersect-intersections object ray)
  (map
   (lambda (inter)
     (<intersect> copy inter [object object] [extra inter]))
   (let intersection ([a (sort-intersections
                          (object-intersections (csg-object1 object) ray))]
                      [b (sort-intersections
                          (object-intersections (csg-object2 object) ray))])
     (cond
      [(null? a) '()]
      [(null? b) '()]
      [else
       (let ([a1 (car a)]
             [a2 (cadr a)]
             [b1 (car b)]
             [b2 (cadr b)])
         (let ([ca1 (<intersect> time a1)]
               [ca2 (<intersect> time a2)]
               [cb1 (<intersect> time b1)]
               [cb2 (<intersect> time b2)])
           (cond
            [(and (= ca1 cb1) (= ca2 cb2)) ; A and B are the same segment
             (append (list a1 a2) (intersection (cddr a) (cddr b)))]
            [(<= ca2 cb1)               ; A is completely less than B
             (intersection (cddr a) b)]
            [(<= cb2 ca1)               ; B is completely less than A
             (intersection a (cddr b))]
            [(and (< ca1 cb1) (> ca2 cb2)) ; A is completely around B
             (intersection (append (list a1 b1 b1 b2 b2 a2) (cddr a)) b)]
            [(and (< cb1 ca1) (> cb2 ca2)) ; B is completely around A
             (intersection a (append (list b1 a1 a1 a2 a2 b2) (cddr b)))]
            [(and (< ca1 cb1) (> cb2 ca2)) ; A and B cross, A first
             (intersection (append (list a1 b1 b1 a2) (cddr a))
               (append (list b1 a2 a2 b2) (cddr b)))]
            [(and (< cb1 ca1) (> ca2 cb2)) ; A and B cross, B first
             (intersection (append (list a1 b2 b2 a2) (cddr a))
               (append (list b1 a1 a1 b2) (cddr b)))]
            [else (errorf 'intersection "Unhandled condition A=~s, B=~s" a b)])))]))))

(define (csg-difference-intersections object ray)
  (map
   (lambda (inter)
     (<intersect> copy inter [object object] [extra inter]))
   (let difference ([a (sort-intersections
                        (object-intersections (csg-object1 object) ray))]
                    [b (sort-intersections
                        (object-intersections (csg-object2 object) ray))])
     (cond
      [(null? a) '()]
      [(null? b) a]
      [else
       (let ([a1 (car a)]
             [a2 (cadr a)]
             [b1 (car b)]
             [b2 (cadr b)])
         (let ([ca1 (<intersect> time a1)]
               [ca2 (<intersect> time a2)]
               [cb1 (<intersect> time b1)]
               [cb2 (<intersect> time b2)])
           (cond
            [(<= ca2 cb1)               ; A is completely less than B
             (append (list a1 a2) (difference (cddr a) b))]
            [(<= cb2 ca1)               ; B is completely less than A
             (difference a (cddr b))]
            [(and (<= ca1 cb1) (> ca2 cb2)) ; A is completely around B
             (difference (append (list a1 b1 b2 a2) (cddr a)) b)]
            [(and (<= cb1 ca1) (> cb2 ca2)) ; B is completely around A
             (difference (cddr a) b)]
            [(and (< ca1 cb1) (>= cb2 ca2)) ; A and B cross, A first
             (difference (append (list a1 b1) (cddr a)) b)]
            [(and (< cb1 ca1) (>= ca2 cb2)) ; A and B cross, B first
             (difference (append (list b2 a2) (cddr a)) b)]
            [(and (= ca1 cb1) (= ca2 cb2)) ; A and B are exactly the same
             (difference (cddr a) (cddr b))]
            [else (errorf 'difference "Unhandled condition ca1=~s, ca2=~s, cb1=~s, cb2=~s" ca1 ca2 cb1 cb2)])))]))))

(define (csg-normal object extra intersect-point)
  (let ([ip (mat-vec-mul (object-Mi object)
              (vec-sub intersect-point (object-center object)))])
    (mat-vec-mul (object-M object)
      (object-normal (<intersect> object extra) (<intersect> extra extra) ip))))