Skip to content

Instantly share code, notes, and snippets.

@soegaard
Created July 28, 2020 16:30
Show Gist options
  • Save soegaard/d105bdc36abda5c92e44e510974392d5 to your computer and use it in GitHub Desktop.
Save soegaard/d105bdc36abda5c92e44e510974392d5 to your computer and use it in GitHub Desktop.
#lang racket/base
(require racket/match racket/list racket/port racket/string
"flomat.rkt")
;;;
;;; Homogeneous coordinates for 3 dimensions
;;;
; Follow the tutorial
; http://www.opengl-tutorial.org/beginners-tutorials/tutorial-3-matrices/
(define (point x y z) (column x y z 1))
(define (dir x y z) (column x y z 0))
(define (normalize p) (column (xc p) (yc p) (zc p)))
(define (xc p)
(define x (ref p 0 0))
(define w (ref p 3 0))
(if (= w 1.0) x (/ x w)))
(define (yc p)
(define y (ref p 1 0))
(define w (ref p 3 0))
(if (= w 1.0) y (/ y w)))
(define (zc p)
(define z (ref p 2 0))
(define w (ref p 3 0))
(if (= w 1.0) z (/ z w)))
;;;
;;; Translation
;;;
(define (translation/xyz x y z)
(flomat: [[1 0 0 x]
[0 1 0 y]
[0 0 1 z]
[0 0 0 1]]))
(define (translation . args)
(match args
[(list p) (translation/xyz (xc p) (yc p) (zc p))]
[(list x y z) (translation/xyz x y z)]
[_ (error 'translation)]))
;;;
;;; Scaling
;;;
(define (scaling/xyz sx sy sz)
(flomat: [[sx 0 0 0]
[0 sy 0 0]
[0 0 sz 0]
[0 0 0 1]]))
(define (scaling . args)
(match args
[(list s) (scaling/xyz s s s)]
[(list sx sy sz) (scaling/xyz sx sy sz)]
[_ (error 'scaling)]))
;;;
;;; Rotation
;;;
(define (rotation-z γ)
(define cγ (cos γ)) (define sγ (sin γ))
(flomat:
[[cγ (- sγ) 0 0]
[sγ cγ 0 0]
[0 0 1 0]
[0 0 0 1]]))
(define (rotation-y β)
(define cβ (cos β)) (define sβ (sin β))
(flomat:
[[ cβ 0 sβ 0]
[ 0 1 0 0]
[(- sβ) 0 cβ 0]
[ 0 0 0 1]]))
(define (rotation-x α)
(define cα (cos α)) (define sα (sin α))
(flomat:
[[1 0 0 0]
[0 cα (- sα) 0]
[0 sα cα 0]
[0 0 0 1]]))
; x y z
(define (rotation/angles γ β α) ; unfortunate naming ...
(define cα (cos α)) (define sα (sin α))
(define cβ (cos β)) (define sβ (sin β))
(define cγ (cos γ)) (define sγ (sin γ))
(flomat:
[[(* cα cβ) (- (* cα sβ sγ) (* sα cγ)) (+ (* cα sβ cγ) (* sα sγ)) 0]
[(* sα cβ) (+ (* sα sβ sγ) (* cα cγ)) (- (* sα sβ cγ) (* cα sγ)) 0]
[ (- sβ) (* cβ sγ) (* cβ cγ) 0]
[0 0 0 1]]))
(define (angle v w)
(acos (/ (dot v w) (* (norm v) (norm w)))))
(define (rotation . args)
(match args
[(list (? number?) (? number?) (? number?))
(apply rotation/angles args)]
[(list (? flomat? direction)) ; a vector
(define α (angle (dir 1 0 0) direction))
(define β (angle (dir 0 1 0) direction))
(define γ (angle (dir 0 0 1) direction))
(rotation/angles α β γ)]))
;;;
;;; Cube
;;;
(define-values (cube-vertices cube-faces)
; Model space coordinates (i.e. (0,0,0) is the center.
(let* ([T (λ (p) ; put center in origo
(times (translation -0.5 -0.5 -0.5) p))]
[A (T (point 0 0 0))]
[B (T (point 1 0 0))]
[C (T (point 1 0 1))]
[D (T (point 0 0 1))]
[E (T (point 0 1 0))]
[F (T (point 1 1 0))]
[G (T (point 1 1 1))]
[H (T (point 0 1 1))])
(define vertices (list A B C D E F G H))
(define-values (a b c d e f g h) (values 0 1 2 3 4 5 6 7))
(define face1 (list a b c d))
(define face2 (list e f g h))
(define face3 (list b f g c))
(define face4 (list a e h d))
(define face5 (list a b f e))
(define face6 (list d c g h))
(define face7 (list a a c c))
(define faces (list face1 face2 face3 face4 face5 face6 #;face7))
(values vertices faces)))
(define cube (apply augment cube-vertices))
;;;
;;; Model Space
;;;
; The model center is in local coordinates in origo.
; To put the model in world coordinates with the proper scaling and rotation,
; we get:
(define (model position direction scale)
(define S (scaling scale))
(define R (rotation direction))
(displayln (list 'position position))
(define T (translation position))
(displayln (list 'S S))
(displayln (list 'R R))
(displayln (list 'T T))
(times T R S))
;;;
;;; View Space (aka Camera Space)
;;;
; (define (view position direction)
;;;
;;; Projections
;;;
; The camera is pointed in direction (0,0,-1).
; It sees a cuboid width x height x (far - near).
(define (ortographic-projection w h near far)
; From View Space into Ortho Projected Space
; The result is with in the [-1;1]³ cube
(define Δz (- far near))
(define Σz (+ far near))
(flomat:
[[(/ w) 0 0 0]
[0 (/ h) 0 0]
[0 0 (/ -2 Δz) (- (/ Σz Δz))]
[0 0 0 1]]))
(define (perspective-projection left right top bottom near far)
(define-values (l r t b n f) (values left right top bottom near far))
(define w (- r l))
(define h (- t b))
(flomat:
[[(/ (* 2 n) w) 0 (/ (+ r l) w) 0]
[0 (/ (* 2 n) h) (/ (+ t b) h) 0]
[0 0 (/ (+ f n) (- n f)) (/ (* 2 f n) (- n f))]
[0 0 -1 0]]))
;;;
;;; Face Sorting
;;;
(define (col! A j)
(define r (nrows A))
(sub! A 0 j r 1))
(define (sort-faces S faces)
; the columns in S are points
; todo: the norm below assumes w=1, which might not be the case...
(define indices (for/list ([j (ncols S)]) j))
(define norms (for/vector ([j (ncols S)]) (norm (col! S j))))
(define (dist f) (apply + (map (λ (i) (vector-ref norms i)) f)))
(sort faces > #:key dist))
;;;
;;; OBJ Parser
;;;
(define (read-obj [port (current-input-port)])
(define lines (for/list ([line (in-lines port)])
(with-input-from-string line
(λ ()
(if (string-prefix? line "#")
(list 'comment line)
(for/list ([x (in-port read)]) x))))))
(define vertice-lines (filter (λ (l) (and (pair? l) (eq? (first l) 'v))) lines))
(define face-lines (filter (λ (l) (and (pair? l) (eq? (first l) 'f))) lines))
(define nverts (length vertice-lines))
(define nfaces (length face-lines))
(define S (make-flomat 4 nverts))
(for ([v vertice-lines] [j (in-naturals)])
(match v
[(list 'v x y z)
(mset! S 0 j x)
(mset! S 1 j y)
(mset! S 2 j z)
(mset! S 3 j 1)]))
(define (convert-triple t)
; input: 1/2/3
; we only use the vertex number (and ignore the texture vertex and
; vertex normal indices.)
(first (map sub1 (map string->number (string-split (symbol->string t) "/")))))
(define faces
(for/list ([fl face-lines])
(match fl
[(list 'f (? number? a) ...) (map sub1 a)]
[(list 'f (? symbol? a) ...) (map convert-triple a)])))
(values S faces))
;;;
;;; Example
;;;
;;;
;;; MetaPict Rendering
;;;
;; (require (except-in metapict shape norm mat? dot dir angle))
;; (define pi (* 4 (atan 1)))
;; (define (polygon pts)
;; (define (line p q) (bez p p q q))
;; (define bezs (let loop ([pts pts])
;; (match pts
;; ['() '()]
;; [(list a) '()]
;; [(list* a b d) (cons (line a b) (loop (cons b d)))])))
;; (curve: #t bezs))
;(set-curve-pict-size 100 100)
#;(apply
beside
(for/list ([v (in-range 0 (+ pi (/ pi 8)) (/ pi 8))])
(with-window (window -1.1 1.1 -1.1 1.1)
(define d (dir (cos v) (sin v) 0.5))
; pos, dir, scale
(define M (model (point 0 0 -10) d 4)) ; model -> world
; (define P (ortographic-projection 3 3 -2 2)) ; world -> device
(define P (perspective-projection -1.5 1.5 -1.5 1.5 -4 4 )) ; world -> device
(define PM (times P M))
(define S (times PM cube)) ; columns in curve are points
(define vs (for/vector ([j (ncols S)]) (define p (col S j)) (pt (xc p) (yc p))))
(define faces (sort-faces S cube-faces))
(for/draw ([f (in-list faces)])
(define pts (for/list ([i f]) (vector-ref vs i)))
(define c (curve* (append (add-between pts --) (list -- 'cycle))))
(draw (color "white" (filldraw c))
c)))))
(define-values (teapot teapot-faces) (read-obj (open-input-file "assets/teapot.obj")))
;(define-values (teapot teapot-faces) (read-obj (open-input-file "assets/teapot-low.obj")))
(define xmax (for/fold ([m -inf.0]) ([j (ncols teapot)]) (max m (ref teapot 0 j))))
(define ymax (for/fold ([m -inf.0]) ([j (ncols teapot)]) (max m (ref teapot 1 j))))
(define zmax (for/fold ([m -inf.0]) ([j (ncols teapot)]) (max m (ref teapot 2 j))))
(define xmin (for/fold ([m +inf.0]) ([j (ncols teapot)]) (min m (ref teapot 0 j))))
(define ymin (for/fold ([m +inf.0]) ([j (ncols teapot)]) (min m (ref teapot 1 j))))
(define zmin (for/fold ([m +inf.0]) ([j (ncols teapot)]) (min m (ref teapot 2 j))))
(list xmin xmax)
(list ymin ymax)
(list zmin zmax)
(define width (* 2 (max (abs xmin) (abs xmax))))
(define height (* 2 (max (abs ymin) (abs ymax))))
(define depth (* 2 (max (abs zmin) (abs zmax))))
; (set-curve-pict-size 600 600)
#;(with-window (window -1.1 1.1 -1.1 1.1)
; (define d (dir (cos v) (sin v) 0.5))
(define d (dir 0 1 0))
; pos, dir, scale
(define M (model (point 0 0 (- -10 7)) d 1)) ; model -> world
(define P (ortographic-projection width height -10 0.5)) ; world -> device
; (define P (perspective-projection -3.5 3.5 -3.5 3.5 -4 4 )) ; world -> device
(define PM (times P M))
(define S (times PM teapot)) ; columns in curve are points
(define vs (for/vector ([j (ncols S)]) (define p (col S j)) (pt (xc p) (yc p))))
(define faces (sort-faces S teapot-faces))
(for/draw ([f (in-list faces)])
(define pts (for/list ([i f]) (vector-ref vs i)))
(define c (polygon pts))
(draw (color "white" (filldraw c))
c)))
(require (except-in pict3d trace dir))
; (define d (dir (cos v) (sin v) 0.5))
(define d (dir -1 -1 0))
; pos, dir, scale
(define M (model (point 0 0 -4) d 1)) ; model -> world
(displayln (list 'M M))
(define P (ortographic-projection width height -10 0.5)) ; world -> device
; (define P (perspective-projection -3.5 3.5 -3.5 3.5 -4 4 )) ; world -> device
(displayln (list 'P P))
(define PM (times P M))
(displayln (list 'PM PM))
(define S (times PM teapot)) ; columns in curve are points
(define vs (for/vector ([j (ncols S)])
(define p (col S j))
(define x (xc p))
(define y (yc p))
(define z (zc p))
(if (member +nan.0 (list x y z))
(begin (display ".") (pos 0 0 0))
(pos x y z))))
(define faces (sort-faces S teapot-faces)) ; no need to sort fact in pict3d
(define p
(combine
(for/list ([f (in-list teapot-faces)])
(define pts (for/list ([i f]) (vector-ref vs i)))
(case (length pts)
[(3) (apply triangle pts #:back? #t)]
[(4) (apply quad pts #:back? #f)]
[else (error)]))
(light (pos 0 2 2) (emitted "white" 2))))
(current-pict3d-fov 80)
(current-pict3d-add-grid? #t)
(current-pict3d-add-wireframe 'color)
(pict3d->bitmap p 600 600)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment