external/clem/src/typed-ops/defmatrix-transform.lisp

   1 ;;; -*- Mode: lisp; outline-regexp: ";;;;;*"; indent-tabs-mode: nil -*-;;;
   2 ;;;
   3 ;;; file: defmatrix-transform.cl
   4 ;;; author: cyrus harmon
   5 ;;;
   6
   7 (in-package :clem)
   8
   9 (defconstant +epsilon+ 0.00001d0)
  10
  11 (defmacro def-matrix-transform (type-1 type-2 transform-type)
  12   (let ((element-type-1 (element-type (find-class `,type-1)))
  13         (element-type-2 (element-type (find-class `,type-2)))
  14         (transform-element-type (element-type (find-class `,transform-type))))
  15     (let ((zero (coerce 0 `,element-type-2))
  16           (one (coerce 1 `,transform-element-type))
  17           (half (coerce 0.5d0 `,transform-element-type))
  18           (epsilon (coerce +epsilon+ `,transform-element-type)))
  19       `(progn
  20          (defmethod %transform-matrix ((m ,type-1) (n ,type-2) (xfrm ,transform-type)
  21                                       &key (background ,zero)
  22                                        (interpolation :nearest-neighbor))
  23            (let ((mr (rows m)) (mc (cols m))
  24                  (nr (rows n)) (nc (cols n)))
  25              (declare (type fixnum mr mc nr nc)
  26                       (type ,element-type-2 background)
  27                       (optimize (speed 3) (safety 0)))
  28              (let  ((inv-xfrm (clem::invert-matrix xfrm))
  29                     (coord1 (make-instance ',transform-type :rows 3 :cols 1))
  30                     (coord2 (make-instance ',transform-type :rows 3 :cols 1))
  31                     #+nil (mrf (coerce mr ',transform-element-type))
  32                     #+nil (mcf (coerce mc ',transform-element-type)))
  33                (let ((a (clem::matrix-vals m))
  34                      (b (clem::matrix-vals n))
  35                      (c (clem::matrix-vals coord1))
  36                      (d (clem::matrix-vals coord2))
  37                      (g (clem::matrix-vals inv-xfrm)))
  38                  (declare (type (simple-array ,element-type-1 *) a)
  39                           (type (simple-array ,element-type-2 *) b)
  40                           (type (simple-array ,transform-element-type *) c d g)
  41                           (optimize (speed 3) (safety 0)))
  42                  (setf (aref c 2 0) ,one)
  43                  (dotimes (i nr)
  44                    (declare (type fixnum i))
  45                    (setf (aref c 0 0) (coerce i ',transform-element-type))
  46                    (dotimes (j nc)
  47                      (declare (type fixnum j))
  48                      (setf (aref c 1 0) (coerce j ',transform-element-type))
  49 ;;; the slow(er) way to do this is:
  50 ;;;   (setf (aref d 0 0) ,zero
  51 ;;;         (aref d 1 0) ,zero
  52 ;;;         (aref d 2 0) ,zero)
  53 ;;; (clem::mat-mult3 inv-xfrm coord1 coord2)
  54 ;;;
  55 ;;; but since we don't need the full matrix multiply,
  56 ;;; based on what we know is in the affine
  57 ;;; transformation matrix, we can get away with
  58 ;;; fewer operations (Foley et al., 1996, p. 213)
  59                      (setf (aref d 1 0) (the ,transform-element-type
  60                                           (+ (* (aref c 1 0) (aref g 0 0))
  61                                              (* (aref c 0 0) (aref g 0 1))
  62                                              (aref g 0 2)))
  63                            (aref d 0 0) (the ,transform-element-type
  64                                           (+ (* (aref c 1 0) (aref g 1 0))
  65                                              (* (aref c 0 0) (aref g 1 1))
  66                                              (aref g 1 2))))
  67 ;;; this does nearest neighbor interpolation
  68 ;;; we should also offer a way to do bilinear,
  69 ;;; bicubic, etc... interpolation
  70                      (case interpolation
  71                        ((:quadratic)
  72                         (if (and
  73                              (< ,(coerce most-negative-fixnum `,transform-element-type)
  74                                 (aref d 0 0)
  75                                 ,(coerce most-positive-fixnum `,transform-element-type))
  76                              (< ,(coerce most-negative-fixnum `,transform-element-type)
  77                                 (aref d 1 0)
  78                                 ,(coerce most-positive-fixnum `,transform-element-type)))
  79                             (multiple-value-bind (l ry)
  80
  81                                 (truncate (the ,transform-element-type
  82                                             (+ (the (,transform-element-type
  83                                                      ,(coerce most-negative-fixnum `,transform-element-type)
  84                                                      ,(coerce most-positive-fixnum `,transform-element-type))
  85                                                  (aref d 0 0))
  86                                                ,epsilon)))
  87                               (declare (type fixnum l)
  88                                        (type ,transform-element-type ry))
  89                               (multiple-value-bind (k rx)
  90                                   (truncate (the ,transform-element-type
  91                                               (+
  92                                                (the (,transform-element-type
  93                                                      ,(coerce most-negative-fixnum `,transform-element-type)
  94                                                      ,(coerce most-positive-fixnum `,transform-element-type))
  95                                                  (aref d 1 0))
  96                                                ,epsilon)))
  97                                 (declare (type fixnum k)
  98                                          (type ,transform-element-type rx))
  99                                 (cond
 100                                   ((and (< -1 l mr)
 101                                         (< -1 k mc))
 102                                    (let ((l0 (max (the fixnum (1- l)) 0))
 103                                          (l2 (min (the fixnum (1+ l)) (the fixnum (1- mr))))
 104                                          (k0 (max (the fixnum (1- k)) 0))
 105                                          (k2 (min (the fixnum (1+ k)) (the fixnum (1- mc)))))
 106                                      (declare (type fixnum l0 l2 k0 k2))
 107                                      (setf (aref b i j)
 108                                            (maybe-truncate
 109                                             (quadratic-interpolate
 110                                              (aref a l0 k0)
 111                                              (aref a l0 k)
 112                                              (aref a l0 k2)
 113                                              (aref a l k0)
 114                                              (aref a l k)
 115                                              (aref a l k2)
 116                                              (aref a l2 k0)
 117                                              (aref a l2 k)
 118                                              (aref a l2 k2)
 119                                              ry rx
 120                                              ,transform-element-type)
 121                                             ,transform-element-type
 122                                             ,element-type-2
 123                                             ))))
 124                                   (t
 125                                    (setf (aref b i j) background)))))
 126                             (setf (aref b i j) background)))
 127                        ((:bilinear :bi-linear)
 128                         (if (and
 129                              (< ,(coerce most-negative-fixnum `,transform-element-type)
 130                                 (aref d 0 0)
 131                                 ,(coerce most-positive-fixnum `,transform-element-type))
 132                              (< ,(coerce most-negative-fixnum `,transform-element-type)
 133                                 (aref d 1 0)
 134                                 ,(coerce most-positive-fixnum `,transform-element-type)))
 135                             (multiple-value-bind (l ry)
 136                                 (floor (the ,transform-element-type
 137                                          (+ (the (,transform-element-type
 138                                                   ,(coerce most-negative-fixnum `,transform-element-type)
 139                                                   ,(coerce most-positive-fixnum `,transform-element-type))
 140                                               (aref d 0 0))
 141                                             ,epsilon)))
 142                               (declare (type fixnum l)
 143                                        (type ,transform-element-type ry))
 144                               (multiple-value-bind (k rx)
 145                                   (floor (the ,transform-element-type
 146                                            (+
 147                                             (the (,transform-element-type
 148                                                   ,(coerce most-negative-fixnum `,transform-element-type)
 149                                                   ,(coerce most-positive-fixnum `,transform-element-type))
 150                                               (aref d 1 0))
 151                                             ,epsilon)))
 152                                 (declare (type fixnum k)
 153                                          (type ,transform-element-type rx))
 154                                 (cond
 155                                   ((and (< -1 l mr)
 156                                         (< -1 k mc))
 157                                    (let ((l1 (the fixnum (min (the fixnum (1+ l)) (the fixnum (1- mr)))))
 158                                          (k1 (the fixnum (min (the fixnum (1+ k)) (the fixnum (1- mc))))))
 159                                      (declare (type fixnum l1 k1))
 160                                      (setf (aref b i j)
 161                                            (maybe-truncate
 162                                             (bilinear-interpolate
 163                                              (aref a l k)
 164                                              (aref a l k1)
 165                                              (aref a l1 k)
 166                                              (aref a l1 k1)
 167                                              ry rx)
 168                                             ,transform-element-type
 169                                             ,element-type-2
 170                                             ))))
 171                                   (t
 172                                    (setf (aref b i j) background)))))
 173                             (setf (aref b i j) background)))
 174                        ((nil :nearest-neighbor)
 175                         (let ((oldx (the fixnum (truncate (+ (aref d 0 0) ,half))))
 176                               (oldy (the fixnum (truncate (+ (aref d 1 0) ,half)))))
 177                           (declare (type fixnum oldx oldy))
 178                           (if (and (< -1 oldx mr)
 179                                    (< -1 oldy mc))
 180                               (setf (aref b i j) (aref a oldx oldy))
 181                               (setf (aref b i j) background))))
 182                        (t (error 'matrix-argument-error
 183                                  :format-control "Unknown interpolation type: ~A"
 184                                  :format-arguments (list  interpolation)))))))))
 185            n)))))
 186
 187 (macrolet ((frob (type-1 type-2 type-3)
 188              `(def-matrix-transform ,type-1 ,type-2 ,type-3)))
 189   (frob double-float-matrix double-float-matrix double-float-matrix)
 190   (frob single-float-matrix single-float-matrix double-float-matrix)
 191   (frob single-float-matrix single-float-matrix single-float-matrix)
 192   (frob ub8-matrix ub8-matrix double-float-matrix)
 193   (frob ub16-matrix ub16-matrix double-float-matrix)
 194   (frob ub32-matrix ub32-matrix double-float-matrix)
 195   (frob sb8-matrix sb8-matrix double-float-matrix)
 196   (frob sb16-matrix sb16-matrix double-float-matrix)
 197   (frob sb32-matrix sb32-matrix double-float-matrix)
 198   (frob fixnum-matrix fixnum-matrix double-float-matrix)
 199   (frob bit-matrix bit-matrix double-float-matrix))