;; -*- Mode: LISP; Package: (JPEG :use (common-lisp)) -*-
;;;
;;; copyright © 1998, 1999 Christopher J. Vogt  All Rights Reserved
;;;                  vogt@computer.org
;;;
;;; This code is free for individual use.  You may not distribute this
;;; code or any derivite work, without prior written consent from author.  If 
;;; you want to use it in a commercial product, email to negotiate for 
;;; reasonable pricing.
;;;
;;; The API is:
;;;       (GET-IMAGE-SIZE NAME)
;;;       (READ-FILE NAME &KEY PROGRESS-FUNCTION IMAGE)
;;;       (WRITE-FILE NAME IMAGE WIDTH HEIGHT &KEY PROGRESS-FUNCTION SAMPLING 
;;;                                                QUANTS AMOUNT COMPONENTS)
;;;       *COMPONENT-ORDER*
;;; GET-IMAGE-SIZE - returns 2 values: WIDTH and HEIGHT of image to be read
;;;    NAME is a file pathname to a jpeg file
;;;
;;; READ-FILE - returns 4 values: IMAGE WIDTH HEIGHT COMPONENTS
;;;    NAME is a file pathname to a jpeg file
;;;    PROGRESS-FUNCTION takes 2 arguments: step-number & text, where text is 
;;;        to be displayed describing what the program is doing, and step-number
;;;        is an integer between 0 and 100, 0 representing nothing done, and 100
;;;        represent 100% done.  This function is to be used to keep track of
;;;        the progress of the read-file operation.
;;;    IMAGE is the array where the image pixels are to be stored.  It should 
;;;        contain numbers, preferably 0.  If image is not passed in, one will
;;;        be created and returned.  GET-IMAGE-SIZE can be used to get width
;;;        and height to create an appropiately sized image array.  The pixels
;;;        are stored as <BGR> w/ Red in the LSB by default, but this is 
;;;        configurable via the *component-order* variable.  Must be of type
;;;        simple-vector (for performance hacking reasons).
;;;
;;; WRITE-FILE - returns no useful value
;;;    NAME is a file pathname to a jpeg file
;;;    IMAGE is the array where the image pixels are to be read.  The pixels
;;;        are stored as <BGR> w/ Red in the LSB by default, but this is 
;;;        configurable via the *component-order* variable.  Must be of type
;;;        simple-vector (for performance hacking reasons).
;;;    WIDTH is the width of the image
;;;    HEIGHT is the height of the image
;;;    PROGRESS-FUNCTION takes 2 arguments: step-number & text, where text is
;;;        to be displayed describing what the program is doing, and step-number
;;;        is an integer between 0 and 100, 0 representing nothing done, and 100
;;;        represent 100% done.  This function is to be used to keep track of
;;;        the progress of the write-file operation.  (this is the same as in
;;;        READ-FILE.
;;;    SAMPLING is a list of lists that specify the amount of sampling to do
;;;        on each component in both dimensions.  there is one list for
;;;        each component, and each list has two values specifying the sampling
;;;        do be done in the horizontal and vertical dimension.  the spec
;;;        limits the amount of sampling that can be done: sum the products
;;;        of the horizontal and vertical sampling for each component, and it
;;;        may not exceed 10.  An example with no sampling of 3 components:
;;;        ((1 1) (1 1) (1 1)).  An example of maximum allowable sampling:
;;;        ((8 1) (1 1) (1 1)).
;;;    QUANTS allows for the specification of a quantization tables.  This is a
;;;        64 element table that is used to quantize the coefficients before
;;;        huffman encoding.  Values of all 1's would cause no quantization
;;;        and minimal loss of quality.  There may be 1,2,3, or 4 tables listed.
;;;        If one table, than it is used for all components.  If 3 or 4 tables
;;;        are listed, then each one is used for each component (if 3 tables
;;;        are listed it is expected that there will be 3 components, and if
;;;        4 tables are listed it is expected that there are 4 components).  If
;;;        2 tables are listed, then the first table is for the first
;;;        component (luminence) and the second table is for the lsat two
;;;        components (chromanance).
;;;    AMOUNT is a value between 0 and 4 specifying how much to compress
;;;        the image.  A value of 0 should result in a nearly lossless image
;;;        a value of 1 should be vitually unrecognizable from 0, but
;;;        should generate a smaller file.  Values of 2 to 4 are more
;;;        aggressive in reducing the file size, with commensurate reduction
;;;        in image quality.  Amount defaults to 0.
;;;    COMPONENTS is the number of components in the image, should be 1,2,3,or 4
;;;
;;; *COMPONENT-ORDER* - one of :RGB :RBG :BGR :BRG :GRB :GBR specifies the
;;;        component order of the final image written.  e.g. :BGR puts Red
;;;        in the LS byte, then Green, then Blue in the MS byte.
;;;
;;; This file contains code that reads most JPEG files.  It supports the
;;; reading of the baseline JPEG files i.e.:
;;;   DCT based 
;;;   8-bit samples
;;;   sequential
;;;   huffman, 4 AC tables and 4 DC tables
;;;   1,2,3, and 4 component scans
;;;   interleaved and non-interleaved scans
;;; It does not support lossless, hierarchical, progressive, or 12-bit files
;;;
;;; Much of the documentation for this file is really in the specification
;;; I used the CCITT Recommendation T.81 (09/92) "Digital Compression and
;;; coding of continuous-tone still images - Requirements and Guidelines"
;;;
;;; The DCT is based an the algorithm in the book "JPEG Still Image Data
;;; Compression Standard" by William B. Pennebaker and Joan L. Mitchell
;;;
;;; at the bottom of the file is an example (commented out) of useage that
;;; was written under LWW 4.0.1
;;;
;;; Design Notes:
;;;    Originally the math was done in floating point, but the performance 
;;;    was horrendous (the development implementation is LWW 4.0.1, which
;;;    promotes all floats to double-float, with lots of boxing, and no
;;;    support for performance improvement via. declareing types.  This
;;;    might not be the case with Alegro or other implementation, but I don't 
;;;    have an implementation to test, and it might be the case that the fixnum 
;;;    version would still be faster). 
;;;    
;;;    The performanc bottleneck seems to be read-byte/write-byte, which in
;;;    LWW takes about 50% of the total time.  I didn't bother trying to
;;;    perform futher optimization as the best I could hope for is a 2x
;;;    improvement, which doesn't seem worth the effort.  read-sequence
;;;    was tried, but I found it to be no faster than read-byte w/ declartions
;;; 
;;;    The bottom line on performance is that this is about 1/5th to 1/15th
;;;    as fast as photoshop in reading/writing jpeg images (call it an order
;;;    of magnitude slower.  Not particularly stellar, but not much I can do 
;;;    about it in a CL portable fasion.
;;;
;;;
;;; CHANGES:
;;; 11/27/98 CJ Vogt - fixed non-compliance with ANSI CL in using macros that 
;;;                    expand to declarations.  Replaced with variables:
;;;                    *optimize* and *fixnum-or-int* to serve the same 
;;;                    functionality, and modified callers of declare-opt.
;;; 11/28/98 CJ Vogt - fixed type declaration compatibility problems, removed
;;;                    *fixnum-or-int*.  changed definition of defsubst to use
;;;                    declaim to make it comform to ANSI.  Also wrapped export
;;;                    in an eval-when
;;;  3/ 1/99 CJ Vogt - make all make-array specify an initial value.  The code
;;;                    was relying on the initial element of an array to be NIL.
;;;                    However, the ANSI CL spec says that the value of the
;;;                    elements are undefined.  Change all calls to make-array
;;;                    that didn't initialize the elements to a value, to now
;;;                    initialize to NIL.
;;;  3/ 4/99 CJ Vogt - another compatability problem, convert all 'unsigned-byte
;;;                    to '(unsigned-byte 8)
;;;  3/16/99 CJ Vogt - when writing a jpeg file, APP0 header should be 16 bytes, 
;;;                    not 14, with the added 2 bytes representing the x,y size 
;;;                    of the thumbnail, so I'll set them both to 0.
;;; 12/16/99 CJ Vogt - fixed bug in macro write-components-ycbcr where it was
;;;                    reading off the end of the array in the case where the
;;;                    height was not an even multiple of the block size.  The
;;;                    code had a min test in the wrong place.
;;; 6/23/05  CJ Vogt - fix bug in copy-image that was manifesting itself when
;;;                    writing one-component files.  Also document that the
;;;                    image must be of type simple-vector.
;;;          CJ Vogt - Finally adding fix from reti in decode-image when 
;;;                    sample factor is > 1.  Untested by me, but I trust Kalman
;;;          CJ Vogt - Added components to be a return value of read-file as
;;;                    suggested by Kalman
;;;

(defpackage "JPEG" (:use "COMMON-LISP"))


(in-package "JPEG")


(eval-when (compile)
  (export '(read-file write-file get-image-size *component-order*)))


;;;
;;; for inlining
;;;
(defmacro defsubst (name args &rest body)
  `(progn
     (declaim (inline ,name))
     (defun ,name ,args ,@body)))


(eval-when (compile load eval)
  (defvar *optimize* '(optimize (safety 0) (space 0) (debug 0) (speed 3))
    "Used to easily switch between debug and delivery of performance sensitive
   code."))


;;;
;;; These macros are for fixnum performance improvement
;;;    
(defmacro i+ (one two)
  (if (>= most-positive-fixnum #x7FFFFF)
      `(the fixnum (+ (the fixnum ,one) (the fixnum ,two)))
    `(+ ,one ,two)))


(defmacro i- (one two)
  (if (>= most-positive-fixnum #x7FFFFF)
      `(the fixnum (- (the fixnum ,one) (the fixnum ,two)))
    `(- ,one ,two)))


(defmacro i* (one two)
  (if (>= most-positive-fixnum #x7FFFFF)
      `(the fixnum (* (the fixnum ,one) (the fixnum ,two)))
    `(* ,one ,two)))


;;;
;;; loops executing body while a restart is thrown;  This is to support
;;; the handling of restart markers in reading files.  If a restart
;;; marker is read, it throws 'restart, which is caught here, and just
;;; initializes *bit-pointer* and dc-delta, and continues executing
;;; body until a restart is not thrown
;;;
(defmacro with-restart (&body body)
  `(loop for result = (catch 'restart ,@body)
         finally (return result)
         while (eq 'restart result) do
         (loop for i fixnum from 0 below components do
               (setf (svref dc-delta i) 0))
         (setq *bit-pointer* 0)))

;;;
;;; converts from YCbCr to RGB and stores the RGB values in the order specified
;;; by *component-order* (which is coded here as MS middle & LS.  This is just
;;; for performanc optimization so that the componnet order testing is done
;;; at compile time generating unique executable code for each possible 
;;; component ordering, rather than performing the test on a per pixel basis.
;;; This is intended to be used only by convert-image-ycbcr-to-rgb.
;;;

(defmacro write-components-rgb (MS middle LS &optional alpha)
  `(loop for y fixnum from 0 below height
         for ybase fixnum from 0 by width do
         (loop for x fixnum from 0 below width
               for addr fixnum = (+ x ybase)
               for data ,(if (not alpha) 'fixnum t) = (svref image addr)
               for y fixnum = (ash (ldb (byte 8 0) data) 14) ; scale y
               for cb fixnum = (i- (ash (ldb (byte 8 8) data) 1) 255)
               for cr fixnum = (i- (ash (ldb (byte 8 16) data) 1) 255)
               ,@(if alpha `(for alpha = (ash (i- (ldb (byte 8 24) data) 128) 24)))
               for red fixnum = (int-clamp (ash (i+ y (i* *1.402* cr)) -14))
               for green fixnum = (int-clamp (ash (i- (i- y (i* *.34414* cb))
                                                      (i* *.71414* cr)) -14))
               for blue fixnum = (int-clamp (ash (i+ y (i* *1.772* cb)) -14)) do
               ,(if alpha
                    `(setf (svref image addr) (+ alpha ,(cond ((eq ls :r)
                                                               'red)
                                                              ((eq ls :g)
                                                               'green)
                                                              (t 'blue))
                                                 (ash ,(cond ((eq middle :r)
                                                              'red)
                                                             ((eq middle :g)
                                                              'green)
                                                             (t 'blue)) 8)
                                                 (ash ,(cond ((eq ms :r)
                                                              'red)
                                                             ((eq ms :g)
                                                              'green)
                                                             (t 'blue)) 16)))
                  `(setf (svref image addr) (i+ (i+ ,(cond ((eq ls :r)
                                                            'red)
                                                           ((eq ls :g)
                                                            'green)
                                                           (t 'blue))
                                                    (ash ,(cond ((eq middle :r)
                                                                 'red)
                                                                ((eq middle :g)
                                                                 'green)
                                                                (t 'blue)) 8))
                                                (ash ,(cond ((eq ms :r)
                                                             'red)
                                                            ((eq ms :g)
                                                             'green)
                                                            (t 'blue)) 16)))))))


;;;
;;; converts YCbCr to RGB and stores the values according to the component 
;;; ordering specified in *component-order*  This is very similar in intent
;;; and rationale as write-components-rgb (see above)  This is intended to
;;; be used only by convert-image-rgb-to-ycbcr
;;;
(defmacro write-components-ycbcr (MS middle LS &optional alpha)
  ms
  `(loop with maxy = (* (1- iheight) ispan)
         for y fixnum from iy below (+ iy height)
         for dy fixnum from 0
         for dbase fixnum from 0 by width
         for ybase fixnum from (* y ispan) by ispan do
	 (setq ybase (min ybase maxy))
         (loop for x fixnum from ix below (+ ix width)
               for dx fixnum from 0
               for daddr fixnum = (+ dx dbase)
               for addr fixnum = (+ (min x (1- ispan)) ybase)
               for data ,(if (not alpha) 'fixnum t) = (svref image addr)
               for r fixnum = (ldb (byte 8 ,(cond ((eq ls :r) 0)
                                                  ((eq middle :r) 8)
                                                  (t 16))) data)
               for g fixnum = (ldb (byte 8 ,(cond ((eq ls :g) 0)
                                                  ((eq middle :g) 8)
                                                  (t 16))) data)
               for b fixnum = (ldb (byte 8 ,(cond ((eq ls :b) 0)
                                                  ((eq middle :b) 8)
                                                  (t 16))) data)
               ,@(if alpha `(for alpha = (ash (i- (ldb (byte 8 24) data) 128)
                                              24)))
               for yy fixnum = (int-clamp (ash (i+ (i+ (i+ (i* *.299* R)
                                                           (i* *.587* G))
                                                       (i* *.114* B))
                                                   *round*) -15))
               for cb fixnum = (int-clamp (ash (i- (i- (i+ (ash b 14)
                                                           (+ *127.5* *round*))
                                                       (i* *.1687* R))
                                                   (i* *.3313* G)) -15))
               for cr fixnum = (int-clamp (ash (i- (i- (i+ (ash r 14)
                                                           (+ *127.5* *round*))
                                                       (i* *.4187* G))
                                                   (i* *.0813* B)) -15)) do
               ,(if alpha
                    `(setf (svref dest daddr) (+ alpha yy
                                                 (ash cb 8)
                                                 (ash cr 16)))
                  `(setf (svref dest daddr)
                         (i+ (i+ yy (ash cb 8)) (ash cr 16)))))))


;;;
;;; This macro is intended to be used only by sample-and-shift,
;;; and generates unique code for each component.  This allows the:
;;; (ldb (byte 8 component-positioni) data) to be compiled efficiently.
;;;
(defmacro sample-and-shift-body (component)
  `(loop with index fixnum = 0
         with samples fixnum = (i* h-sample v-sample)
         for y fixnum from iy by v-sample
         for ybase = (i* y iw)
         repeat 8 do
         (loop for sum fixnum = 0
               for x fixnum from ix by h-sample
               repeat 8 do
               (loop for sy fixnum from 0 repeat v-sample
                     for sybase fixnum = (i* sy iw) do
                     (loop for sx fixnum from 0 repeat h-sample do
                           (setf sum (i+ sum (ldb (byte 8 ,(ash component 3))
                                                  (svref temp (i+ (i+ (i+ x sx)
                                                                      ybase)
                                                                  sybase)))))))
               (setf (svref zz index)
                     (i- (ash sum (svref *sample-shift* samples)) 128))
               (setf index (i+ index 1)))))


(defconstant *sof0* #xC0 "Start Of Frame Baseline DCT Huffman")
(defconstant *sof1* #xC1 "Start Of Frame Extended DCT Huffman")
(defconstant *sof2* #xC2 "Start Of Frame Progressive DCT Huffman")
(defconstant *sof3* #xC3 "Start Of Frame Lossless Huffman")
(defconstant *sof5* #xC5 "Start Of Frame Hierarchical Sequential DCT Huffman")
(defconstant *sof6* #xC6 "Start Of Frame Hierarchical Progressive DCT Huffman")
(defconstant *sof7* #xC7 "Start Of Frame Hierarchical Lossless DCT Huffman")
(defconstant *sof9* #xC9 "Start Of Frame Extended DCT Arithmetic")
(defconstant *sofA* #xCA "Start Of Frame Progressive DCT Arithmetic")
(defconstant *sofB* #xCB "Start Of Frame Lossless Arithmetic")
(defconstant *sofD* #xCD "Start Of Frame Hierarchical DCT Arithmetic")
(defconstant *sofE* #xCE "Start Of Frame Hierarchical Progressive DCT 
                          Arithmetic")
(defconstant *sofF* #xCF "Start Of Frame Hierarchical Lossless DCT Arithmetic")

(defconstant *soi*  #xD8 "Start Of Image")
(defconstant *eoi*  #xd9 "End Of Image")
(defconstant *dht* #xC4 "Huffman Coding Table")
(defconstant *dqt* #xDb "Quantization table")
(defconstant *sos* #xDA "Start of Scan - this begins the compressed data")
(defconstant *com* #xFE "COM marker, Comment")
(defconstant *app0* #xE0 "Application specific marker #0")
(defconstant *appC* #xEC "Application specific marker #C") ; treat like COM
(defconstant *rst* #xD0 "Restart with modulo 8 count m (D0 to D7)")
(defconstant *dnl* #xDC "Define number of lines")

(defconstant *ratios* #(0 0 1 2 1 4 2 0 1 8  4  2  0  0  0  0  1)
  "Hack to quickly figure out how many times a component has been sampled,
   and therefore how many times we need to duplicate it.")

(defconstant *zigzag* #(0 1  8 16  9  2  3 10
                          17 24 32 25 18 11  4  5
                          12 19 26 33 40 48 41 34
                          27 20 13  6  7 14 21 28
                          35 42 49 56 57 50 43 36
                          29 22 15 23 30 37 44 51
                          58 59 52 45 38 31 39 46
                          53 60 61 54 47 55 62 63) 
  "The AC coefficients & quantization values are encoded in the JPEG file in
   zigzag order.  We use this to do/undo this and get them in lineal order.")

(defconstant *dct-fixed-constants*
  (make-array 64 :initial-contents
              (loop for float in 
                    (loop with line =
                          (cons 1 (loop for k from 1 to 7
                                        collect (* (cos (* k (/ PI 16)))
                                                   (sqrt 2))))
                          for y in line appending
                          (loop for x in line
                                collect (* x y)))
                    collect (round (* 2048 float))))
  "Scaled (by 2048 or 11-bits) fixed point versions of DCT coefficients")

;;;
;;; the following constants support the encoding of jpeg images
;;; by providing some default quantization tables and huffman tables.
;;;
(defconstant *ones* (make-array 64 :initial-element 1) "Quantization table")
(defconstant *base-luminance* #(16 11 10 16 24 40 51 61
                                   12 12 14 19 26 58 60 55
                                   14 13 16 24 40 57 69 56
                                   14 17 22 29 51 87 80 62
                                   18 22 37 56 68 109 103 77
                                   24 35 55 64 81 104 113 92
                                   49 64 78 87 103 121 120 101
                                   72 92 95 98 112 100 103 99)
  "Quantization table")
(defconstant *base-chrominance* #(17 18 24 47 99 99 99 99
                                     18 21 26 66 99 99 99 99
                                     24 26 56 99 99 99 99 99
                                     47 66 99 99 99 99 99 99
                                     99 99 99 99 99 99 99 99
                                     99 99 99 99 99 99 99 99
                                     99 99 99 99 99 99 99 99
                                     99 99 99 99 99 99 99 99)
  "Quantization table")
(defconstant *luminance* #(4  3  2  4  6 10 13 15
			      3  3  4  5  6 14 15 14
 			      4  3  4  6 10 14 17 14
			      4  4  6  7 13 22 20 16
                              4  6  9 14 17 27 26 19
                              6  9 14 16 20 26 28 23
                              12 16 20 22 26 30 30 25
	                      18 23 24 24 28 25 26 25)
  "Quantization table")
(defconstant *chrominance* #(4  4  6 12 25 25 25 25
			        4  5  6 16 25 25 25 25
	 		        6  6 14 25 25 25 25 25
	                        12 16 25 25 25 25 25 25
		                25 25 25 25 25 25 25 25
	                        25 25 25 25 25 25 25 25
                                25 25 25 25 25 25 25 25
	                        25 25 25 25 25 25 25 25)
  "Quantization table")
(defconstant *max-luminance* #(32 24 16 32 48 80 104 120
                                  24 24 32 40 48 112 120 112
                                  32 24 32 48 80 112 136 112
                                  32 32 48 56 104 176 160 128
                                  32 48 72 112 136 216 208 152
                                  48 72 112 128 160 208 224 184
                                  96 128 160 176 208 240 240 200
                                  144 184 192 192 224 200 208 200))
(defconstant *max-chrominance* #(32 32 48 96 200 200 200 200
                                    32 40 48 128 200 200 200 200
                                    48 48 112 200 200 200 200 200
                                    96 128 200 200 200 200 200 200
                                    200 200 200 200 200 200 200 200
                                    200 200 200 200 200 200 200 200
                                    200 200 200 200 200 200 200 200
                                    200 200 200 200 200 200 200 200))
(defconstant *dc-huffman-luminance-bits* '(0 0 1 5 1 1 1 1 1 1 0 0 0 0 0 0 0))
(defconstant *dc-huffman-chrominance-bits* '(0 0 3 1 1 1 1 1 1 1 1 1 0 0 0 0 0))
(defconstant *ac-huffman-luminance-bits* '(0 0 2 1 3 3 2 4 3 5 5 4 4 0 0 1 #x7d))
(defconstant *ac-huffman-chrominance-bits* '(0 0 2 1 2 4 4 3 4 7 5 4 4 0 1 2 #x77))
(defconstant *ac-huffman-luminance-values*
  #(#x01 #x02 #x03 #x00 #x04 #x11 #x05 #x12 #x21 #x31 #x41 #x06 #x13 #x51 #x61 #x07
         #x22 #x71 #x14 #x32 #x81 #x91 #xa1 #x08 #x23 #x42 #xb1 #xc1 #x15 #x52 #xd1 #xf0
         #x24 #x33 #x62 #x72 #x82 #x09 #x0a #x16 #x17 #x18 #x19 #x1a #x25 #x26 #x27 #x28
         #x29 #x2a #x34 #x35 #x36 #x37 #x38 #x39 #x3a #x43 #x44 #x45 #x46 #x47 #x48 #x49
         #x4a #x53 #x54 #x55 #x56 #x57 #x58 #x59 #x5a #x63 #x64 #x65 #x66 #x67 #x68 #x69
         #x6a #x73 #x74 #x75 #x76 #x77 #x78 #x79 #x7a #x83 #x84 #x85 #x86 #x87 #x88 #x89
         #x8a #x92 #x93 #x94 #x95 #x96 #x97 #x98 #x99 #x9a #xa2 #xa3 #xa4 #xa5 #xa6 #xa7
         #xa8 #xa9 #xaa #xb2 #xb3 #xb4 #xb5 #xb6 #xb7 #xb8 #xb9 #xba #xc2 #xc3 #xc4 #xc5
         #xc6 #xc7 #xc8 #xc9 #xca #xd2 #xd3 #xd4 #xd5 #xd6 #xd7 #xd8 #xd9 #xda #xe1 #xe2
         #xe3 #xe4 #xe5 #xe6 #xe7 #xe8 #xe9 #xea #xf1 #xf2 #xf3 #xf4 #xf5 #xf6 #xf7 #xf8
         #xf9 #xfa))
(defconstant *ac-huffman-chrominance-values*
  #(0 1 2 3 17 4 5 33 49 6 18 65 81 7 97 113
      19 34 50 129 8 20 66 145 161 177 193 9 35 51 82 240
      21 98 114 209 10 22 36 52 225 37 241 23 24 25 26 38
      39 40 41 42 53 54 55 56 57 58 67 68 69 70 71 72
      73 74 83 84 85 86 87 88 89 90 99 100 101 102 103 104
      105 106 115 116 117 118 119 120 121 122 130 131 132 133 134 135
      136 137 138 146 147 148 149 150 151 152 153 154 162 163 164 165
      166 167 168 169 170 178 179 180 181 182 183 184 185 186 194 195
      196 197 198 199 200 201 202 210 211 212 213 214 215 216 217 218
      226 227 228 229 230 231 232 233 234 242 243 244 245 246 247 248 249 250))
        
(defconstant *integer-length*
  (make-array 2047 :initial-contents (loop for i from -1023 to 1023
                                           collect (integer-length (abs i))))
  "A table to find the length of integers i.e. the lisp function integer-length")

(defconstant *sample-shift* (make-array 9 :initial-contents '(0 0 -1 0 -2 0 0 0 -3))
  "Provides a mapping between the number of samples, and the amount need to shift
   to effect a divide")

(defconstant *dc0-code* (make-array 4096 :initial-element nil))
(defconstant *dc0-size* (make-array 4096 :initial-element nil))
(defconstant *dc1-code* (make-array 4096 :initial-element nil))
(defconstant *dc1-size* (make-array 4096 :initial-element nil))

;;;
;;; these constants are used in conversion from/to RGB/YCbCr
;;;

(defconstant *.299* (round (* 0.299 (ash 1 15))))
(defconstant *.587* (round (* 0.587 (ash 1 15))))
(defconstant *.114* (round (* 0.114 (ash 1 15))))
(defconstant *.1687* (round (* 0.1687 (ash 1 15))))
(defconstant *.3313* (round (* 0.3313 (ash 1 15))))
(defconstant *.4187* (round (* 0.4187 (ash 1 15))))
(defconstant *.0813* (round (* 0.0813 (ash 1 15))))
(defconstant *127.5* (round (* 127.5 (ash 1 15))))
(defconstant *round* (ash 1 14))
(defconstant *1.402* (round (* 1.402 (ash 1 13))))
(defconstant *.34414* (round (* 0.34414 (ash 1 13))))
(defconstant *.71414* (round (* 0.71414 (ash 1 13))))
(defconstant *1.772* (round (* 1.772 (ash 1 13))))


(defvar *file* nil "For debugging, gets set to the file object")

(defvar *current-byte* "Was a slot in File. Due to performance issues, I wanted 
                        the methods that use this slot to be in-lined.
                        Due to compiler limitations, I found this untenable.
                        What I want is something like defsubst-method, but all
                        I have is defsubst, so these became variables.
                        (see huffman-nextbit)")

(defvar *bit-pointer* "same comment as current-byte")

(defvar *component-order* :BGR "Defines the order in which color coomponents are
                                are stored in a memory word.  The default is
                                :BGR with R in the LSB.  Should be one of:
                                :RGB :RBG :BGR :BRG :GRB :GBR.  If there is an
                                Alpha component, it is always the MSB")


;;;
;;; read 2 bytes from stream and create a halfword MSB first
;;; return the value or NIL if EOF reached
;;;
(defsubst read-halfword (stream)
  (let ((one (read-byte stream nil nil))
        (two (read-byte stream nil nil)))
    (and one two (+ (ash one 8) two))))


;;;
;;; all markers are of the form #xFFyz.  #xFFFF is padding
;;; return NIL if a marker is not found, otherwise return the LSB of the marker
;;;
(defsubst read-marker (stream)
  (and (equal (read-byte stream nil nil) #xFF) ; start of marker
       (loop for marker? = (read-byte stream nil nil)
             when (or (null marker?)
                      (not (= marker? #xff))) do
             (return marker?))))


;;;
;;; Test that the next 2 bytes read are a marker, and that the LSB is = to the 
;;; specified marker
;;; Return T or NIL
;;;
(defsubst test-marker (stream marker)
  (let ((byte (read-marker stream)))
    (and byte (= marker byte))))

  
;;;
;;; Clamps value to a value between 0 and 255 inclusive
;;;
(defsubst int-clamp (value)
  (declare #.*optimize*
           (type fixnum value))
  (if (> value 254)
      255
    (if (< value 1)
        0
      value)))

;;;
;;; The component-table slot in the file object has 5 values packed into 1
;;; value, the routines below unpack/pack the specified information
;;; see also the describe-component-table method below
;;;
(defsubst horizontal-sample-factor (table index)
  (ldb (byte 3 0) (svref table index)))


(defsubst vertical-sample-factor (table index)
  (ldb (byte 3 3) (svref table index)))


(defsubst quantization-index (table index)
  (ldb (byte 2 6) (svref table index)))


(defsubst ac-huffman-index (table index)
  (ldb (byte 2 8) (svref table index)))


(defsubst dc-huffman-index (table index)
  (ldb (byte 2 10) (svref table index)))

(defsubst make-component-table-value (horizontal vertical quant)
  (dpb quant (byte 2 6) (dpb vertical (byte 3 3) (ldb (byte 3 0) horizontal))))


;;;
;;; Gets and returns the next bit from the stream.  #xFF is a marker indicator,
;;; but it is also used as data if it is followed by #x0 (the #x0 in this case
;;; is ignored) If it is a marker, and it is one of the restart markers (e.g. 
;;; *rst*) then we throw out to the restart catcher, so we can start reading
;;; data again.
;;;
;;; *bit-pointer* & *current-byte* were originally slots in the File object, but 
;;; i wanted this routine in-lined for improved performance.  What I really
;;; need is something like: defsubst-method, but I'm not sure how to write it
;;;
(defsubst huffman-nextbit (byte-stream)
  (when (zerop *bit-pointer*) ; we have used all the bits, get another byte
    (setq *current-byte* (read-byte byte-stream))
    (setq *bit-pointer* 128)
    (when (= #xff *current-byte*) ; is it a marker?
      (let ((next-byte (read-byte byte-stream)))
        (when (not (= 0 next-byte)) ; if the next byte is 0, ignore it
          (if (or (= *eoi* next-byte) (= *dnl* next-byte))
              (error "We need more bytes, but we reached the end of input!")
            (if (or (< next-byte *rst*) (> next-byte (i+ 8 *rst*)))
                (error "Error processing byte-stream, we need more bits,
                        but we found a marker?")
              (throw 'restart 'restart))))))) ; reached restart marker
  (prog1 (if (zerop (logand *bit-pointer* *current-byte*)) 0 1)
    (setq *bit-pointer* (ash *bit-pointer* -1))))


;;;
;;; Using the huffman tables maxcode, mincode, valptr, & huffval, decodes the
;;; input stream, returning the next decoded byte of data
;;;
(defsubst huffman-decode (maxcode mincode valptr huffval byte-stream)
  (declare #.*optimize*
           (type (simple-vector *) maxcode mincode valptr huffval))
  (loop with code fixnum = (huffman-nextbit byte-stream)
        for i fixnum from 1
        while (> code (svref maxcode i))
        finally (return (svref huffval (i- (i+ (svref valptr i) code)
                                           (svref mincode i)))) do
        (setq code (i+ (ash code 1) (huffman-nextbit byte-stream)))))


;;;
;;; get nbits from the input stream
;;;
(defsubst huffman-receive (nbits byte-stream)
  (declare #.*optimize*
           (type fixnum nbits))
  (loop with v fixnum = 0
        for i fixnum from 0 below nbits
        finally (return v) do
        (setq v (i+ (ash v 1) (huffman-nextbit byte-stream)))))


;;;
;;; decodes DC value, adds it to the contents of dc-delta for this component
;;;
(defsubst dc-huffman-decode (component maxcode mincode valptr huffval
                                       byte-stream dc-delta)
  (declare #.*optimize*
           (type fixnum component)
           (type (simple-vector *)
                 maxcode mincode valptr huffval dc-delta))
  (let* ((value (huffman-decode maxcode mincode valptr huffval byte-stream))
         (diff (huffman-receive value byte-stream)))
    (setf (svref dc-delta component) (i+ (svref dc-delta component)
                                         (huffman-extend diff value)))))


;;;
;;; decodes the 63 AC values and stores them in zz
;;;
(defsubst ac-huffman-decode (zz maxcode mincode valptr huffval byte-stream)
  (declare #.*optimize*
           (type (simple-vector *)
                 zz *zigzag* maxcode mincode valptr huffval))
  (loop for i fixnum from 1 to 63 do
        (setf (svref zz i) 0))
  (loop with k fixnum = 1
        while (< k 64)
        for rs fixnum = (huffman-decode maxcode mincode valptr huffval byte-stream)
        for ssss fixnum = (ldb (byte 4 0) rs)
        for rrrr fixnum = (ash rs -4) do
        (if (= 0 ssss)
            (if (= 15 rrrr)
                (setq k (i+ k 16))  ; 16 zeros
              (return)) ; EOB rest are 0
          (progn
            (setq k (i+ k rrrr))
            (if (> k 63)
                (return))
            (setf (svref zz (svref *zigzag* k))
                  (huffman-extend (huffman-receive ssss byte-stream) ssss))
            (setq k (i+ 1 k))))))
  

;;;
;;; returns the number of tables in the tables array, used when writing files
;;;
(defsubst count-quantization-tables (tables)
  (loop for i from 0
        for table = (svref tables i)
        while table
        finally (return i)))


;;;
;;; gets the index to the quantization table based on the number of tables and
;;; the component, this is for writing files
;;;
(defsubst get-quantization-index (component tables)
  (if (= 1 tables) ; if 1 table, index 0 for all componnents
      0
    (if (= 2 tables) ; if 2 tables, index 0 for first
        (if (zerop component)
            0
          1)
      component))) ; if there is a table for each component


;;;
;;; Writes size bits of data to stream.  Since #xff is a special marker
;;; if we have #xff to be written out, it is stuffed with a #x0 following
;;; it to indicate that it is not a marker, but rather to be interpreted
;;; literally.
;;;
(defsubst write-bits (stream size data)
  (declare #.*optimize*
           (type fixnum size data))
  (loop with pointer fixnum = (ash 1 (i- size 1))
        with bit-pointer fixnum = *bit-pointer*
        with current-byte fixnum = *current-byte*
        while (> pointer 0)
        finally (progn (setq *bit-pointer* bit-pointer)
                  (setq *current-byte* current-byte)) do
        (setq current-byte (i+ (ash (if (zerop (logand pointer data)) 0 1)
                                    bit-pointer) current-byte))
        (setq bit-pointer (i- bit-pointer 1))
        (setq pointer (ash pointer -1))
        (when (< bit-pointer 0)
          (write-byte current-byte stream)
          (if (= #xff current-byte)
              (write-byte 0 stream))
          (setq bit-pointer 7)
          (setq current-byte 0))))


;;;
;;; Each non-zero AC coefficient in zz is described by a composite 8-bit value
;;; of the form RRRRSSSS where r is the number of 0s we have seen so far (and
;;; this is where they are going to get encoded) and s is the length of the
;;; data.  We huffman-encode RRRRSSSS and write it to the stream, then we
;;; write out ssss bits of data to the stream
;;;
(defsubst encode (stream r data huffman-size huffman-code)
  (declare #.*optimize*
           (type (simple-vector *)
                 huffman-size huffman-code *integer-length*)
           (type fixnum r data))
  (let* ((ssss (svref *integer-length* (i+ 1023 data)))
         (rs (i+ ssss (ash r 4))))
    (declare (type fixnum ssss rs))
    (write-bits stream (svref huffman-size rs) (svref huffman-code rs))
    (if (< data 0)
        (setq data (i- data 1)))
    (write-bits stream ssss data)))


;;;
;;; Huffman encodes the 63 (1...63) AC coefficients of zz to stream
;;; 15 zeros in a row get special encodeing (rather than calling 
;;; encode) and if the last "n" coefficients are zero, they get a 
;;; special encoding as well.
;;;
(defsubst encode-ac-coefficients (stream zz huffman-size huffman-code)
  (declare #.*optimize*
           (type (simple-vector *)
                 huffman-size huffman-code zz *zigzag*))
  (loop with r fixnum = 0
        finally (if (plusp r) (write-bits stream (svref huffman-size 0)
                                          (svref huffman-code 0)))
        for k fixnum from 1 to 63
        for data fixnum = (svref zz (svref *zigzag* k)) do
        (if (zerop data)
            (setq r (i+ r 1))
          (progn
            (if (> data 1023)
                (setq data 1023)
              (if (< data -1023)
                  (setq data -1023)))
            (loop while (> r 15) do
                  (write-bits stream (svref huffman-size #xf0)
                              (svref huffman-code #xf0))
                  (setq r (i- r 16)))
            (encode stream r data huffman-size huffman-code)
            (setq r 0)))))


;;;
;;; The dc coefficient is encoded as a differential value.  dc-delta
;;; holds the previous value, and gets updated with the new value.
;;; we encode the difference between what is stored in dc-delta and
;;; value.
;;;
(defsubst encode-dc-coefficient (byte-stream dc-delta value component)
  (declare #.*optimize*
           (type (simple-vector *) dc-delta)
           (type fixnum value component))
  (let* ((pred (svref dc-delta component))
         (diff (i- value pred)))
    (declare (type fixnum pred diff))
    (if (> diff 2047)
        (setq diff 2047)
      (if (< diff -2047)
          (setq diff -2047)))
    (if (zerop component)
        (write-bits byte-stream
                    (svref *dc0-size* (i+ 2047 diff))
                    (svref *dc0-code* (i+ 2047 diff)))
      (write-bits byte-stream
                  (svref *dc1-size* (i+ 2047 diff))
                  (svref *dc1-code* (i+ 2047 diff))))
    (setf (svref dc-delta component) value)))


;;;
;;; the file object is mostly used for encapsulating data, and there is a
;;; lot of data to be encapsulated
;;;
(defclass file ()
  ((byte-stream :initform nil
                :initarg :byte-stream
                :accessor file-byte-stream)
   (width :initform nil
          :initarg :width
          :accessor file-width)
   (height :initform nil
           :initarg :height
           :accessor file-height)
   (components :initform nil
               :initarg :components
               :accessor file-components)
   (bits-per-component :initform nil
                       :initarg :bits-per-component
                       :accessor file-bits-per-component)
   (version :initform ""
            :initarg :version
            :accessor file-version)
   (units :initform nil
          :initarg :units
          :accessor file-units)
   (x-density :initform nil
              :initarg :x-density
              :accessor file-x-density)
   (y-density :initform nil
              :initarg :y-density
              :accessor file-y-density)
   (quantization-tables :initform (make-array 4 :initial-element nil)
                        :accessor file-quantization-tables)
   (huffman-values :initform (make-array 8 :initial-element nil)                 ; DC: 0-3;  AC: 4-7
                   :accessor file-huffman-values)
   (huffman-bits :initform (make-array 8 :initial-element nil)                   ; DC: 0-3;  AC: 4-7
                 :accessor file-huffman-bits)
   (huffman-size :initform (make-array 8 :initial-element nil)
                 :accessor file-huffman-size)
   (huffman-code :initform (make-array 8 :initial-element nil)
                 :accessor file-huffman-code)
   (huffman-ordered-code :initform (make-array 256 :initial-element nil)
                         :accessor file-huffman-ordered-code)
   (huffman-ordered-size :initform (make-array 256 :initial-element nil)
                         :accessor file-huffman-ordered-size)
   (huffman-mincode :initform (make-array 8 :initial-element nil)
                    :accessor file-huffman-mincode)
   (huffman-maxcode :initform (make-array 8 :initial-element nil)
                    :accessor file-huffman-maxcode)
   (huffman-valptr :initform (make-array 8 :initial-element nil)
                   :accessor file-huffman-valptr)
   (component-tables :initform (make-array 4 :initial-element nil)
                     :accessor file-component-tables)  ; contains packed data
   (dc-delta :initform (make-array 4 :initial-element nil)
             :accessor file-dc-delta)
   (ss :initform nil
       :accessor file-ss)
   (se :initform nil
       :accessor file-se)
   (ah :initform nil
       :accessor file-ah)
   (al :initform nil
       :accessor file-al)
   (zz :initform (make-array 64 :initial-element nil)
       :accessor file-zz)
   (image :initform nil
          :accessor file-image)
   (progress :initform nil
             :accessor file-progress)))


;;;
;;; just prints the name of the file and the width and height
;;;
(defmethod print-object ((self file) stream)
  (with-slots (width height) self
    (format stream "#<FILE: ~a ~d x ~d>"
            (truename (file-byte-stream self)) width height)))


;;;
;;; Read and optionally print out header information
;;; Return T if header read correctly, else return nil
;;;
(defmethod read-header ((self file) &optional verbose)
  (with-slots (byte-stream version units x-density y-density) self
    (let* (length
           (valid-jpeg-file (and (test-marker byte-stream *soi*)
                                 (test-marker byte-stream *app0*)
                                 (setq length (read-halfword byte-stream))
                                 (equal (read-byte byte-stream nil nil)
                                        (char-code #\J))
                                 (equal (read-byte byte-stream nil nil)
                                        (char-code #\F))
                                 (equal (read-byte byte-stream nil nil)
                                        (char-code #\I))
                                 (equal (read-byte byte-stream nil nil)
                                        (char-code #\F))
                                 (equal (read-byte byte-stream nil nil) 0)
                                 (setq version (read-halfword byte-stream))
                                 (setq units (read-byte byte-stream nil nil))
                                 (setq x-density (read-halfword byte-stream))
                                 (setq y-density
                                       (read-halfword byte-stream)))))
      (when valid-jpeg-file
        (setq version (format nil "~d.~d"
                              (ldb (byte 8 8) version)
                              (ldb (byte 8 0) version)))
        ;;
        ;; read & ignore remaining bytes in header
        ;;
        (loop repeat (max 0 (- length 14)) do
              (read-byte byte-stream nil nil))
        (when verbose
          (format t "~%Version: ~a" version)
          (when (> units 0)
            (format t "~%Horizontal Density: ~a ~a"
                    x-density (if (= units 1) "dots/inch" "dots/cm"))
            (format t "~%Vertical Density: ~a ~a"
                    y-density (if (= units 1) "dots/inch" "dots/cm"))))
        t))))
      

;;;
;;; Returns 2 values: marker length, where marker is the last marker read
;;; and length is the length of the marker.  It reads and sets up the
;;; file object until a *EOI* or *SOS* marker is reached (or some 
;;; error condition occurs such as EOF.
;;;
(defmethod process-markers ((self file))
  (with-slots (byte-stream width height components bits-per-component) self
    (loop for marker = (read-marker byte-stream)
          for length = (read-halfword byte-stream)
          until (or (null marker) (null length) (= marker *eoi*)
                    (= marker *sos*))
          finally (return (values marker length)) do 
          (cond ((= marker *dht*) ; huffman table
                 (get-dht self length))
                ((= marker *dqt*) ; quantization table
                 (get-dqt self length))
                ((or (= marker *sof0*))
                 (setq bits-per-component (read-byte byte-stream nil nil))
                 (setq height (read-halfword byte-stream))
                 (setq width (read-halfword byte-stream))
                 (setq components (read-byte byte-stream nil nil))
                 (if (or (null height) (zerop height) (null width)
                         (zerop width) (null components) (< components 1)
                         (> components 4))
                     (error "~%SOF marker doesn't specify a legal width,
                               height, or # of components")
                   (get-sof self components)))
                ((or (= marker *sof3*) (= marker *sof5*) (= marker *sof6*)
                     (= marker *sof7*) (= marker *sofB*) (= marker *sofD*)
                     (= marker *sofE*) (= marker *sof1*) (= marker *sof2*)
                     (= marker *sof9*) (= marker *sofA*) (= marker *sofF*))
                 (error "~%Contains an unsupported encoding"))
                (t ; lots of markers fall in here
                 (loop repeat (max 0 (- length 2))
                       for byte = (read-byte byte-stream nil nil)
                       while byte))))))


;;;
;;; If a quantization table does not already exist at the index, create one
;;;
(defmethod make-quantization-table ((self file) index)
  (with-slots (quantization-tables) self
    (if (null (svref quantization-tables index))
        (setf (svref quantization-tables index) (make-array 64 :initial-element nil)))))


;;;
;;; If huffman tables for the index do not yet exist, create them
;;;
(defmethod make-huffman-table ((self file) index)
  (with-slots (huffman-values huffman-mincode huffman-maxcode huffman-valptr
                              huffman-bits huffman-size huffman-code) self
    (when (null (svref huffman-values index))
      (setf (svref huffman-values index) (make-array 256 :initial-element nil))
      (setf (svref huffman-size index) (make-array 256 :initial-element 0))
      (setf (svref huffman-code index) (make-array 256 :initial-element 0))
      (setf (svref huffman-mincode index) (make-array 17 :initial-element 0))
      (setf (svref huffman-maxcode index) (make-array 17 :initial-element 0))
      (setf (svref huffman-valptr index) (make-array 17 :initial-element 0))
      (setf (svref huffman-bits index) (make-array 17 :initial-element nil)))))


;;;
;;; Reads jpeg files saving the image data in image
;;;

(defun read-file (name &key progress-function image)
  "Reads the jpeg file name and returns 4 values: image, width, height and components.
   Where image is an array of pixels, width is the width of the read image
   and height is the height of the image.  If image is supplied as a keyword
   argument, it is used, otherwise an array is created and returned.
   NAME is a file pathname to a jpeg file
   PROGRESS-FUNCTION takes 2 arguments: step-number & text, where text is 
        to be displayed describing what the program is doing, and step-number
        is an integer between 0 and 100, 0 representing nothing done, and 100
        represent 100% done.  This function is to be used to keep track of
        the progress of the read-file operation.
   IMAGE is the array where the image pixels are to be stored.  It should 
        contain numbers, preferably 0.  If image is not passed in, one will
        be created and returned.  GET-IMAGE-SIZE can be used to get width
        and height to create an appropiately sized image array.  The pixels
        are stored as <BGR> w/ Red in the LSB by default, but this is 
        configurable via the *component-order* variable.  Must be of type 
        simple-vector for performance hacking reasons."
  (if (< most-positive-fixnum #x7FFFFF)
      (error "Sorry, due to the heavily optimized code, and required precision
              fixnums need to be at least 24 bits to work")
    (with-open-file (stream name :element-type '(unsigned-byte 8))
      (let ((*bit-pointer* nil)
            (*current-byte* nil)
            (file (make-instance 'file :byte-stream stream)))
        (setq *file* file)
        (setf (file-progress file) progress-function)
        (file-read file image)))))


;;;
;;; Performs necessary initialization, reads the file header information
;;; and calls decode-image and convert-image-ycbcr-to-rgb which are the
;;; two real workhorses in reading jpeg files.
;;;
(defmethod file-read ((self file) &optional image-array)
  (with-slots (byte-stream width height components image dc-delta
                           quantization-tables) self
    (if (null (read-header self))
        (error " does not appear to be a valid JPEG file")
      (multiple-value-bind (marker length) (process-markers self)
        (when (= marker *sos*)
          (huffman-generate-size-tables self) ; create *real* huffman tables
          (huffman-generate-code-tables self)
          (huffman-make-decode-tables self)
          (get-sos self length)
          (cond ((null image-array)
                 (if (> (* width height) array-total-size-limit)
                     (error "Image too large to display"))
                 (setq image (make-array (* width height) :initial-element 0)))
                (t
                 (if (< (length image-array) (* width height))
                     (error "Image-array is too small for file"))
                 (setq image image-array)
                 (if (null (svref image 0))
                     (loop for i from 0 below (* width height) do
                           (setf (svref image i) 0)))))
          (loop for i from 0 to 3
                for table = (svref quantization-tables i)
                when table do
                (loop for i from 0 to 63 do
                      (setf (svref table i)
                            (i* (svref table i)
                                (svref *dct-fixed-constants* i)))))
          (setq *bit-pointer* 0) ; we haven't read anything yet
          (loop for i from 0 below components do 
                (setf (svref dc-delta i) 0)) 
          ;;
          ;; This is the meat of the program
          ;;
          (decode-image self)
          ;;
          ;; If there is 1 or 2 components, do nothing more.
          ;; If there are 3 components, they are in YCbCr format, convert
          ;; them to RGB format using *component-order* to put each byte
          ;; in the appropriate place.
          ;;
          (if (> components 2)
              (convert-image-ycbcr-to-rgb image width height (= components 4))))))
    (values image width height components)))


;;;
;;; reads the bytes that are huffman encoded, and decodes them and writes
;;; them out to the image.  This image is not RGB, but rather YCbCr
;;; assuming it has 3 components.
;;;

(defmethod decode-image ((self file))
  (declare #.*optimize*)
  (with-slots (width height components component-tables image dc-delta zz
                     quantization-tables progress byte-stream huffman-maxcode
                     huffman-mincode huffman-valptr huffman-values) self
    (declare (type fixnum width height components)
             (type (simple-vector *) huffman-maxcode huffman-mincode
                   huffman-valptr huffman-values quantization-tables)
             (type (simple-vector *) zz image dc-delta component-tables))
    (loop with max-width fixnum =
          (loop for i fixnum from 0 below components
                for table = (svref component-tables i)
                for factor = (and table (horizontal-sample-factor
                                         component-tables i))
                when factor
                maximize factor)
          with max-height fixnum =
          (loop for i fixnum from 0 below components
                for table = (svref component-tables i)
                for factor = (and table
                                  (vertical-sample-factor component-tables i))
                when factor
                maximize factor)
          with maxy = (ceiling (ceiling height max-height) 8)
          with maxx = (ceiling (ceiling width max-width) 8)
          for y fixnum from 0 below maxy do
          (if progress
              (funcall progress (round (* 100 y) maxy) "Decoding Image"))
          (loop for x fixnum from 0 below maxx do
                ;;
                ;; this is an MCU (there has got to be a better way than this)
                ;; maybe unroll the loop and case on how many components so
                ;;  this all gets done outside the loop
                ;;
                (loop for i fixnum from 0 below components
                      for hs fixnum = (horizontal-sample-factor
                                       component-tables i)
                      for vs fixnum = (vertical-sample-factor
                                       component-tables i)
                      for h-ratio fixnum = (svref *ratios* (i+ max-width hs))
                      for v-ratio fixnum = (svref *ratios* (i+ max-height vs))
                      for ac-index fixnum =
                      (logior 4 (the fixnum (ac-huffman-index component-tables
                                                              i)))
                      for dc-index fixnum = (dc-huffman-index component-tables
                                                              i)
                      for dmaxcode = (svref huffman-maxcode dc-index)
                      for dmincode = (svref huffman-mincode dc-index)
                      for dvalptr = (svref huffman-valptr dc-index)
                      for dhuffval = (svref huffman-values dc-index)
                      for amaxcode = (svref huffman-maxcode ac-index)
                      for amincode = (svref huffman-mincode ac-index)
                      for avalptr = (svref huffman-valptr ac-index)
                      for ahuffval = (svref huffman-values ac-index)
                      for quantization-table =
                      (svref quantization-tables
                             (quantization-index component-tables i)) do
                      (loop for v fixnum from 0 below vs do
                            (loop for h fixnum from 0 below hs 
				  for xpos fixnum = (ash (i+ h (i* x hs)) 3)
				  for ypos fixnum = (ash (i+ v (i* y vs)) 3)
				  unless (or (>= xpos width)
					     (>= ypos height))
				  do
				  (with-restart
                                    (dc-huffman-decode i dmaxcode dmincode
                                                       dvalptr dhuffval
                                                       byte-stream dc-delta))
                                  (with-restart
                                    (ac-huffman-decode zz amaxcode amincode
                                                       avalptr ahuffval
                                                       byte-stream))
                                  (setf (svref zz 0) (svref dc-delta i))
                                  (dequantize zz quantization-table)
                                  (inverse-dct zz)
                                  (write-inverse-dct zz image
                                                     xpos
                                                     ypos
                                                     i width height h-ratio
                                                     v-ratio))))))))


;;; 
;;; converts the image array from YCbCr to RGB using scaled fixnums
;;;
(defun convert-image-ycbcr-to-rgb (image width height alpha)
  (declare #.*optimize*
           (type (simple-vector *) image)
           (type fixnum width height))
  (if alpha
      (cond ((eq *component-order* :BGR)
             (write-components-rgb :B :G :R :A))
            ((eq *component-order* :BRG)
             (write-components-rgb :B :R :G :A))
            ((eq *component-order* :RGB)
             (write-components-rgb :R :G :B :A))
            ((eq *component-order* :RBG)
             (write-components-rgb :R :B :G :A))
            ((eq *component-order* :GRB)
             (write-components-rgb :G :R :B :A))
            ((eq *component-order* :GBR)
             (write-components-rgb :G :B :R :A))
            (t (error "Illegal component specified, must be one of: :BGR :BRG :RGB
                   :RBG :GRB :GBR")))
    (cond ((eq *component-order* :BGR)
           (write-components-rgb :B :G :R))
          ((eq *component-order* :BRG)
           (write-components-rgb :B :R :G))
          ((eq *component-order* :RGB)
           (write-components-rgb :R :G :B))
          ((eq *component-order* :RBG)
           (write-components-rgb :R :B :G))
          ((eq *component-order* :GRB)
           (write-components-rgb :G :R :B))
          ((eq *component-order* :GBR)
           (write-components-rgb :G :B :R))
          (t (error "Illegal component specified, must be one of: :BGR :BRG :RGB
                   :RBG :GRB :GBR")))))


;;;
;;; get the 64 element quantization table
;;; there might be more than one table per marker
;;; there might be more than one marker per file
;;;
(defmethod get-dqt ((self file) length)
  (with-slots (byte-stream quantization-tables) self
    (decf length 2)
    (loop while (> length 0)
          for info = (read-byte byte-stream)
          for one-byte-entries = (and info (zerop (ldb (byte 4 4) info)))
          for index = (and info (ldb (byte 4 0) info))
          when info do
          (make-quantization-table self index)
          (loop for i from 0 to 63
                for element = (if one-byte-entries
                                  (read-byte byte-stream)
                                (read-halfword byte-stream)) do
                (setf (svref (svref quantization-tables index)
                             (svref *zigzag* i)) element))
          (if one-byte-entries
              (decf length 65)
            (decf length 129)))))


;;;
;;; when processing the file header, this is called to get process
;;; the sof header.
;;;
(defmethod get-sof ((self file) components)
  (with-slots (byte-stream component-tables) self
    (loop for i from 0 below components
          for index = (1- (read-byte byte-stream))
          for sample-factor = (read-byte byte-stream) do
          (setf (svref component-tables index)
                (make-component-table-value (ldb (byte 4 4) sample-factor)
                                            (ldb (byte 4 0) sample-factor)
                                            (read-byte byte-stream))))))


;;;
;;; when processing the file header, this is called to get process
;;; the sos header.
;;;
(defmethod get-sos ((self file) length)
  (with-slots (byte-stream component-tables ss se ah al) self
    (decf length 2)
    (loop while (> length 0)
          with components = (read-byte byte-stream)
          for i from 0 below components
          for index = (1- (read-byte byte-stream))
          for table-selector = (read-byte byte-stream) do
          (setf (svref component-tables index)
                (dpb (ldb (byte 2 4) table-selector) (byte 2 10)
                     (dpb (ldb (byte 2 0) table-selector) (byte 2 8)
                          (svref component-tables index)))))
    (setq ss (read-byte byte-stream))
    (setq se (read-byte byte-stream))
    (setq ah (read-byte byte-stream))
    (setq al (ldb (byte 4 0) ah))
    (setq ah (ldb (byte 4 4) ah))))


;;;
;;; when processing the file header, this is called to get process
;;; the dht header.
;;;
(defmethod get-dht ((self file) length)
  (with-slots (byte-stream huffman-values huffman-bits) self
    (decf length 2)
    (loop with huffman-table
          with bits-table
          while (> length 16)
          for count = 0
          for index = (ldb (byte 5 0) (read-byte byte-stream)) do
          (setq index (ldb (byte 3 0)
                           (dpb (ldb (byte 1 4) index) (byte 1 2) index)))
          (make-huffman-table self index)
          (setq huffman-table (svref huffman-values index))
          (setq bits-table (svref huffman-bits index))
          (setf (svref bits-table  0) 0)
          (loop for i from 1 to 16 do
                (setf (svref bits-table i) (read-byte byte-stream))
                (incf count (svref bits-table i)))
          (loop for i from 0 to (1- count) do
                (setf (svref huffman-table i) (read-byte byte-stream)))
          (decf length (+ 17 count)))))


;;;
;;; returns list of last-k for each table
;;;
(defmethod huffman-generate-size-tables ((self file))
  (with-slots (huffman-bits huffman-size) self
    (loop for i from 0 to 7
          for bits = (svref huffman-bits i)
          for size = (svref huffman-size i)
          when bits
          collect
          (loop with k = 0
                finally (progn (setf (svref size k) 0) (return k))
                for i from 1 to 16 do
                (loop for j from 1 to (svref bits i) do
                      (setf (svref size k) i)
                      (incf k))))))


(defmethod huffman-generate-code-tables ((self file))
  (with-slots (huffman-code huffman-size) self
    (loop for i from 0 to 7
          for size = (svref huffman-size i)
          for code-table = (svref huffman-code i)
          when size do
          (loop with k = 0
                with code = 0
                with si = (svref size 0) do
                (loop do 
                      (setf (svref code-table k) code)
                      (incf code)
                      (incf k)
                      while (= si (svref size k)))
                (if (= 0 (svref size k))
                    (return)
                  (loop do 
                        (setq code (ash code 1))
                        (incf si)
                        while (not (= si (svref size k)))))))))


(defmethod huffman-make-decode-tables ((self file))
  (with-slots (huffman-bits huffman-code huffman-maxcode huffman-mincode
                            huffman-valptr) self
    (loop for i from 0 to 7
          for code = (svref huffman-code i)
          for maxcode = (svref huffman-maxcode i)
          for mincode = (svref huffman-mincode i)
          for valptr = (svref huffman-valptr i)
          for bits = (svref huffman-bits i)
          when code do
          (loop with j = 0
                for i from 1 to 16 do
                (if (= 0 (svref bits i))
                    (setf (svref maxcode i) -1)
                  (progn
                    (setf (svref valptr i) j)
                    (setf (svref mincode i) (svref code j))
                    (incf j (1- (svref bits i)))
                    (setf (svref maxcode i) (svref code j))
                    (incf j)))))))


(defun huffman-extend (value precision)
  (declare #.*optimize*
           (type fixnum value precision))
  (let ((temp (ash 1 (i- precision 1))))
    (declare (type fixnum temp))
    (when (< value temp)
      (setq temp (i+ 1 (ash -1 precision)))
      (setq value (i+ value temp)))
    value))


;;;
;;; mutliples the two arrays element for element
;;; stores the result in zz
;;;
(defun dequantize (zz quant)
  (loop for i from 0 to 63 do
        (setf (svref zz i) (* (svref zz i) (svref quant i)))))


;;;
;;; performs the inverse-dct on the zz array which contains the
;;; values that have already been quantized.  Does this in situ
;;; in two passes.
;;;
(defun inverse-dct (zz)
  (declare #.*optimize*
           (type (simple-vector *) zz))
  ;;
  ;; first pass, columns
  ;;
  (loop for i from 0 to 7
        for tmp0 = (svref zz i)
        for tmp1 = (svref zz (i+ i 16))
        for tmp2 = (svref zz (i+ i 32))
        for tmp3 = (svref zz (i+ i 48))
        for tmp4 = (svref zz (i+ i 8))
        for tmp5 = (svref zz (i+ i 24))
        for tmp6 = (svref zz (i+ i 40))
        for tmp7 = (svref zz (i+ i 56))
        for tmp10 = (i+ tmp0 tmp2)
        for tmp11 = (i- tmp0 tmp2)
        for tmp13 = (i+ tmp1 tmp3)
        for tmp12 = (i- (i* (ash (i- tmp1 tmp3) -8)
                            (round (* 1.414213562 256))) tmp13)
        for z13 = (i+ tmp6 tmp5)
        for z10 = (i- tmp6 tmp5)
        for z11 = (i+ tmp4 tmp7)
        for z12 = (i- tmp4 tmp7)
        for z5 = (i* (ash (i+ z10 z12) -8) (round (* 1.847759065 256))) do
        (setq tmp0 (i+ tmp10 tmp13))
        (setq tmp1 (i+ tmp11 tmp12))
        (setq tmp2 (i- tmp11 tmp12))
        (setq tmp3 (i- tmp10 tmp13))
        (setq tmp7 (i+ z11 z13))
        (setq tmp11 (i* (ash (i- z11 z13) -8) (round (* 1.414213562 256))))
        (setq tmp10 (i- (i* (ash z12 -8) (round (* 1.082392200 256))) z5))
        (setq tmp12 (i+ (i* (ash z10 -8)(round (* -2.613125930 256))) z5))
        (setq tmp6 (i- tmp12 tmp7))
        (setq tmp5 (i- tmp11 tmp6))
        (setq tmp4 (i+ tmp10 tmp5))
        (setf (svref zz i) (i+ tmp0 tmp7))
        (setf (svref zz (i+ i 8)) (i+ tmp1 tmp6))
        (setf (svref zz (i+ i 16)) (i+ tmp2 tmp5))
        (setf (svref zz (i+ i 24)) (i- tmp3 tmp4))
        (setf (svref zz (i+ i 32)) (i+ tmp3 tmp4))
        (setf (svref zz (i+ i 40)) (i- tmp2 tmp5))
        (setf (svref zz (i+ i 48)) (i- tmp1 tmp6))
        (setf (svref zz (i+ i 56)) (i- tmp0 tmp7)))
  ;;
  ;; 2nd pass, rows
  ;;
  (loop for i from 0 to 7
        for j from 0 by 8
        with tmp4
        with tmp5
        with tmp6
        for tmp10 = (i+ (svref zz j) (svref zz (i+ j 4)))
        for tmp11 = (i- (svref zz j) (svref zz (i+ j 4)))
        for tmp13 = (i+ (svref zz (i+ j 2)) (svref zz (i+ j 6)))
        for tmp12 = (i- (i* (ash (i- (svref zz (i+ j 2)) (svref zz (i+ j 6)))
                                 -8)
                            (round (* 1.414213562 256))) tmp13)
        for tmp0 = (i+ tmp10 tmp13)
        for tmp3 = (i- tmp10 tmp13)
        for tmp1 = (i+ tmp11 tmp12)
        for tmp2 = (i- tmp11 tmp12)
        for z13 = (i+ (svref zz (i+ j 5)) (svref zz (i+ j 3)))
        for z10 = (i- (svref zz (i+ j 5)) (svref zz (i+ j 3)))
        for z11 = (i+ (svref zz (i+ j 1)) (svref zz (i+ j 7)))
        for z12 = (i- (svref zz (i+ j 1)) (svref zz (i+ j 7)))
        for tmp7 = (i+ z11 z13)
        for z5 = (i* (ash (i+ z10 z12) -8) (round (* 1.847759065 256))) do
        (setq tmp11 (i* (ash (i- z11 z13) -8) (round (* 1.414213562 256))))
        (setq tmp10 (i- (i* (ash z12 -8) (round (* 1.082392200 256))) z5))
        (setq tmp12 (i+ (i* (ash z10 -8)(round (* -2.613125930 256))) z5))
        (setq tmp6 (i- tmp12 tmp7))
        (setq tmp5 (i- tmp11 tmp6))
        (setq tmp4 (i+ tmp10 tmp5))
        (setf (svref zz j) (int-clamp (i+ 128 (ash (i+ tmp0 tmp7) -14))))
        (setf (svref zz (i+ j 1))
              (int-clamp (i+ 128 (ash (i+ tmp1 tmp6) -14))))
        (setf (svref zz (i+ j 2))
              (int-clamp (i+ 128 (ash (i+ tmp2 tmp5) -14))))
        (setf (svref zz (i+ j 3))
              (int-clamp (i+ 128 (ash (i- tmp3 tmp4) -14))))
        (setf (svref zz (i+ j 4))
              (int-clamp (i+ 128 (ash (i+ tmp3 tmp4) -14))))
        (setf (svref zz (i+ j 5))
              (int-clamp (i+ 128 (ash (i- tmp2 tmp5) -14))))
        (setf (svref zz (i+ j 6))
              (int-clamp (i+ 128 (ash (i- tmp1 tmp6) -14))))
        (setf (svref zz (i+ j 7))
              (int-clamp (i+ 128 (ash (i- tmp0 tmp7) -14))))))

           
;;;
;;; this gets called once for each dct (each 8x8 block for each component)
;;; it writes the dct to destination, and if the component has been sampled
;;; we duplicate the data when writing to the destination
;;;
(defun write-inverse-dct (dct destination startx starty component span height hrep
                              vrep)
  (declare #.*optimize*
           (type fixnum startx starty component span hrep vrep height)
           (type (simple-vector *) dct)
           (type (simple-vector *) destination))
  (setq startx (i* startx hrep))
  (setq starty (i* starty vrep))
  (loop for i fixnum from 0 to 63
        with index fixnum = (i+ startx (i* starty span))
        with shift-amount fixnum = (ash component 3)
        for xpart fixnum = (ldb (byte 3 0) i)
        for ypart fixnum = (ldb (byte 3 3) i)
        for data fixnum = (svref dct i) do
        ;;
        ;; If the component was sampled before encoding, we will need
        ;; to write it out multiple times to the destination.  Also
        ;; since the width and height are filled to gaurantee full
        ;; 8x8 sections, it may have been packed with bits we need
        ;; to ignore, so we have to test to be sure we don't write
        ;; off the end of the image array
        ;;
        (loop for y fixnum from 0 below vrep
              for ybase fixnum from (i* span (i* ypart vrep)) by span
              when (< (i+ (i+ starty (i* ypart vrep)) y) height) do
              (loop for x fixnum from 0 below hrep
                    for xbase fixnum from (+ index (i* xpart hrep))
                    for next-loc fixnum = (i+ ybase xbase)
                    when (< (i+ (i+ startx (i* xpart hrep)) x) span) do
                    (setf (svref destination next-loc)
                          (i+ (ash data shift-amount)
                              (svref destination next-loc)))))))



;;;
;;; parses the image file, verifying it is a JPEG file, and returns 2
;;; values: width height of the image contained in the file
;;;
(defun get-image-size (name)
  "Returns 2 values: opens and reads enough information from the file name so
   that the width and  height of the image in the file can be determined, and
   then returns 2 values: width height"
  (with-open-file (stream name :element-type '(unsigned-byte 8))
    (let ((file (make-instance 'file :byte-stream stream)))
      (if (null (read-header file))
          (error " does not appear to be a valid JPEG file")
        (multiple-value-bind (marker length) (process-markers file)
          (when (= marker *sos*)
            (get-sos file length)
            (values (file-width file) (file-height file))))))))


;;;
;;; Top level function
;;;
(defun write-file (name image width height &key progress-function sampling
                        quants (amount 0) (components 3))
  "The JPEG algorithm can compress files from 2x to 100x and more, with
   very little loss in image quality.  This compression is achieved in
   three different ways: 1. Huffman coding of image data; 2. Optional sub-
   sampling of image components; 3. Optional quantization of image data;  The 
   huffman coding is lossless, and achieves about 2x compression by itself.  
   Sub sampling is inherently lossy, but the sampling is done in the Y,Cb,Cr
   domain, not in the R,G,B domain, and CbCr can be subsampled with very little
   loss in observable image fidelity, with obvious gains in reduction of 
   amount of data to be stored.  Quantization can be performed on the coefficients
   of the DCT (discrete cosine transform).  If the quantization tables are 
   all 1's, then no quntization is done.  A JPEG file can be saved virtually
   lossless using no subsampling, and using quantization tables of all 1s.
   The reason it is only 'virtually' lossless is because of rounding errors
   and color space conversion loss of accurracy.  If the quantization tables are
   not all 1s, then the DCT coefficients are divided by these values, pushing
   the resulting values close to 0.  The huffman encoding is *very* efficient
   at encoding 0s, so this is the other way in which image data can get compressed
   and if you choose the quantization values csvrefully, you can get increased
   image compression with limited effects on perceptible image quality.

   This function returns no useful value.

    NAME is a file pathname to a jpeg file
    IMAGE is the array where the image pixels are to be read.  The pixels
        are stored as <BGR> w/ Red in the LSB by default, but this is 
        configurable via the *component-order* variable.  Must be of type 
        simple-vector for performance hacking reasons.
    WIDTH is the width of the image
    HEIGHT is the height of the image
    PROGRESS-FUNCTION takes 2 arguments: step-number & text, where text is
        to be displayed describing what the program is doing, and step-number
        is an integer between 0 and 100, 0 representing nothing done, and 100
        represent 100% done.  This function is to be used to keep track of
        the progress of the read-file operation.  (this is the same as in
        READ-FILE.
    SAMPLING is a list of lists that specify the amount of sampling to do
        on each component in both dimensions.  there is one list for
        each component, and each list has two values specifying the sampling
        do be done in the horizontal and vertical dimension.  the spec
        limits the amount of sampling that can be done: sum the products
        of the horizontal and vertical sampling for each component, and it
        may not exceed 10.  An example with no sampling of 3 components:
        ((1 1) (1 1) (1 1)).  An example of maximum allowable sampling:
        ((8 1) (1 1) (1 1)).
    QUANTS allows for the specification of a quantization tables.  This is a
        64 element table that is used to quantize the coefficients before
        huffman encoding.  Values of all 1's would cause no quantization
        and minimal loss of quality.  There may be 1,2,3, or 4 tables listed.
        If one table, than it is used for all components.  If 3 or 4 tables
        are listed, then each one is used for each component (if 3 tables
        are listed it is expected that there will be 3 components, and if
        4 tables are listed it is expected that there are 4 components).  If
        2 tables are listed, then the first table is for the first
        component (luminence) and the second table is for the lsat two
        components (chromanance).
    AMOUNT is a value between 0 and 4 specifying how much to compress
        the image.  A value of 0 should result in a nearly lossless image
        a value of 1 should be vitually unrecognizable from 0, but
        should generate a smaller file.  Values of 2 to 4 are more
        aggressive in reducing the file size, with commensurate reduction
        in image quality.  Amount defaults to 0.
    COMPONENTS is the number of components in the image, should be 1,2,3,or 4"
  (if (< most-positive-fixnum #x7FFFFF)
      (error "Sorry, due to the heavily optimized code, fixnums need to be
              at least 24 bits to work")
    (if (not (= (length image) (* width height)))
        (error "Image doesn't match the width and height specification given")
      (with-open-file (stream name :element-type '(unsigned-byte 8) :direction :output)
        (let ((file (make-instance 'file :byte-stream stream)))
          (setq *file* file)
          (setf (file-progress file) progress-function)
          (initialize-dc-huffman-encoding-tables)
          (if (null sampling)
              (setq sampling (case amount
                               (0 '((1 1) (1 1) (1 1) (1 1)))
                               (1 '((1 1) (2 1) (2 1) (1 1)))
                               (2 '((1 1) (2 1) (2 1) (1 1)))
                               (3 '((1 1) (1 1) (1 1) (1 1)))
                               (t '((1 1) (2 2) (2 2) (1 1))))))
          (if (null quants)
              (setq quants (case amount
                             (0 (list *ones* *ones*)) ; about 2x compression
                             (1 (list *ones* *ones*)) ; about 3x compression
                             (2 (list *luminance* *chrominance*)) ; about 7x compression
                             (3 (list *base-luminance* *base-chrominance*)) ; about 15x compression
                             (t (list *max-luminance* *max-chrominance*))))) ; about 30x compression
          (loop for i in quants
                for j from 0 do
                (setf (svref (file-quantization-tables file) j) (copy-seq i)))
          (setf (file-components file) components)
          (setf (file-image file) image)
          (setf (file-width file) width)
          (setf (file-height file) height)
          (file-write file sampling))))))


(defmethod initialize-huffman-tables ((self file))
  (with-slots (huffman-bits huffman-values huffman-size huffman-code) self
    (setf (svref huffman-bits 0)
          (make-array 17 :initial-contents *dc-huffman-luminance-bits*))
    (setf (svref huffman-bits 1)
          (make-array 17 :initial-contents *dc-huffman-chrominance-bits*))
    (setf (svref huffman-bits 4)
          (make-array 17 :initial-contents *ac-huffman-luminance-bits*))
    (setf (svref huffman-bits 5)
          (make-array 17 :initial-contents *ac-huffman-chrominance-bits*))
    (setf (svref huffman-values 0)
          (make-array 12 :initial-contents '(0 1 2 3 4 5 6 7 8 9 10 11)))
    (setf (svref huffman-values 1)
          (make-array 12 :initial-contents '(0 1 2 3 4 5 6 7 8 9 10 11)))
    (setf (svref huffman-values 4) *ac-huffman-luminance-values*)
    (setf (svref huffman-values 5) *ac-huffman-chrominance-values*)
    (setf (svref huffman-size 0) (make-array 256 :initial-element 0))
    (setf (svref huffman-size 1) (make-array 256 :initial-element 0))
    (setf (svref huffman-size 4) (make-array 256 :initial-element 0))
    (setf (svref huffman-size 5) (make-array 256 :initial-element 0))
    (setf (svref huffman-code 0) (make-array 256 :initial-element 0))
    (setf (svref huffman-code 1) (make-array 256 :initial-element 0))
    (setf (svref huffman-code 4) (make-array 256 :initial-element 0))
    (setf (svref huffman-code 5) (make-array 256 :initial-element 0))
    (let ((last-ks (huffman-generate-size-tables self)))
      (huffman-generate-code-tables self)
      (huffman-order-code-tables self last-ks))))


;;;
;;; Performs necessary initialization and writes the file header information
;;; converts the image from RGB to YCbCr, perfforms sampling as needed,
;;; perfoms the DCT, quantizes the coeeficients, and encodes and writes out
;;; the data.
;;;
(defmethod file-write ((self file) sampling)
  (with-slots (byte-stream image components quantization-tables dc-delta width
                           height zz progress) self
    (declare #.*optimize*
             (type fixnum components width height)
             (type (simple-vector *) image zz dc-delta))
    (let ((*bit-pointer* 7)
          (*current-byte* 0))
      (initialize-huffman-tables self)
      (write-header self sampling) ; includes "JFIF" and *soi*
      (loop for i from 0 to 3
            for table = (svref quantization-tables i)
            when table do
            (loop for i from 0 to 63 do
                  (setf (svref table i)
                        (* (svref *dct-fixed-constants* i)
                           (svref table i)))))
      (loop for i from 0 below components do 
            (setf (svref dc-delta i) 0))
      (let* ((max-h (loop for i in sampling maximize (first i)))
             (max-v (loop for i in sampling maximize (second i)))
             (vertical-count (ceiling (ceiling height max-v) 8))
             (horizontal-count (ceiling (ceiling width max-h) 8))
             (iw (* 8 max-h))
             (ih (* 8 max-v))
             (size (* iw ih))
             (temp (make-array size :initial-element 0))
             (tables (count-quantization-tables quantization-tables)))
        (declare (type fixnum max-h max-v vertical-count horizontal-count iw ih
                       size)
                 (type (simple-vector *) temp))
        (loop for iy fixnum from 0 by ih repeat vertical-count do
              (if progress
                  (funcall progress (round (* 100 iy) (* vertical-count ih)) "Encoding Image"))
              (loop for ix fixnum from 0 by iw repeat horizontal-count do
                    (if (> components 2)
                        (convert-image-rgb-to-ycbcr image ix iy width height
                                                    temp iw ih (= components 4))
                      (copy-image image ix iy width height temp iw ih))
                    (loop for component fixnum from 0 below components
                          for quant = (svref quantization-tables
                                             (get-quantization-index
                                              component tables))
                          for samples in sampling
                          for h-sample fixnum = (first samples)
                          for v-sample fixnum = (second samples) do
                          (loop for iy fixnum from 0 by 8
                                repeat (ash max-v
                                            (svref *sample-shift* v-sample)) do
                                (loop for ix fixnum from 0 by 8
                                      repeat (ash max-h (svref *sample-shift*
                                                               h-sample)) do
                                      (sample-and-shift temp ix iy iw h-sample
                                                        v-sample zz component)
                                      (forward-dct zz)
                                      (quantize zz quant)
                                      (encode-image self component)))))))
      (if (< *bit-pointer* 7) ; we need to pack in final bits
          (write-bits byte-stream (1+ *bit-pointer*) #xff))
      (write-marker *eoi* byte-stream))))


;;;
;;; copys pixels from the specified portion of image array into the temp
;;; array, filling in extra rows/columns via duplication
;;;
(defun copy-image (image ix iy width height temp iw ih)
  (loop with index = 0
        for y from iy
        for ybase = (* (min y (1- height)) width)
        repeat ih do
        (loop for x from ix
              repeat iw do
              (setf (svref temp index)
                    (svref image (+ (min x (1- width)) ybase)))
              (incf index))))


;;;
;;; performs sampling as specified, and level shifts and stores in zz
;;;
(defun sample-and-shift (temp ix iy iw h-sample v-sample zz component)
  (declare #.*optimize*
           (type fixnum ix iy iw h-sample v-sample component)
           (type (simple-vector *) zz)
           (type (simple-vector *) temp))
  (cond ((= component 0)
         (sample-and-shift-body 0))
        ((= component 1)
         (sample-and-shift-body 1))
        ((= component 2)
         (sample-and-shift-body 2))
        ((= component 3)
         (sample-and-shift-body 3))))
  

;;;
;;; converts huffman-code and huffman-size into ordered tables for encodeing
;;;
(defmethod huffman-order-code-tables ((self file) last-ks)
  (with-slots (huffman-values huffman-code huffman-size huffman-ordered-code
                              huffman-ordered-size) self
    (loop for i from 0 to 7
          for value = (svref huffman-values i)
          for code = (svref huffman-code i)
          for size = (svref huffman-size i)
          with last-k
          when value do
          (setq last-k (pop last-ks))
          (loop for i from 0 to 255 do ;initialize tables to 0
                (setf (svref huffman-ordered-code i) 0)
                (setf (svref huffman-ordered-size i) 0))
          (loop for k from 0 below last-k
                for i = (svref value k) do
                (setf (svref huffman-ordered-code i) (svref code k))
                (setf (svref huffman-ordered-size i) (svref size k)))
          (loop for i from 0 to 255 do ;copy tables
                (setf (svref code i) (svref huffman-ordered-code i))
                (setf (svref size i) (svref huffman-ordered-size i))))))


;;;
;;; creates lookup tables for the size and code of DC coefficients
;;;
(defun initialize-dc-huffman-encoding-tables ()
  (when (null (svref *dc0-code* 0))
    (setf (svref *dc0-code* 2047) 0)
    (setf (svref *dc1-code* 2047) 0)
    (setf (svref *dc0-size* 2047) 2)
    (setf (svref *dc1-size* 2047) 2)
    (loop for size in '(3 3 3 3 3 4 5 6 7 8 9)
          for code in '(2 3 4 5 6 14 30 62 126 254 510)
          for ssss from 1 do
          (loop for i from (expt 2 (1- ssss)) below (expt 2 ssss) do
                (setf (svref *dc0-size* (+ 2047 i)) (+ size ssss))
                (setf (svref *dc0-size* (- 2047 i)) (+ size ssss))
                (setf (svref *dc0-code* (+ 2047 i)) (dpb code (byte 11 ssss) i))
                (setf (svref *dc0-code* (- 2047 i)) (dpb code (byte 11 ssss)
                                                         (ldb (byte ssss 0)
                                                              (1- (- i)))))))
    (loop for size in '(2 2 3 4 5 6 7 8 9 10 11)
          for code in '(1 2 6 14 30 62 126 254 510 1022 2046)
          for ssss from 1 do
          (loop for i from (expt 2 (1- ssss)) below (expt 2 ssss) do
                (setf (svref *dc1-size* (+ 2047 i)) (+ size ssss))
                (setf (svref *dc1-size* (- 2047 i)) (+ size ssss))
                (setf (svref *dc1-code* (+ 2047 i)) (dpb code (byte 11 ssss) i))
                (setf (svref *dc1-code* (- 2047 i)) (dpb code (byte 11 ssss)
                                                         (ldb (byte ssss 0)
                                                              (1- (- i)))))))))


;;;
;;; perform huffman encoding of the block
;;;
(defmethod encode-image ((self file) component)
  (with-slots (dc-delta zz byte-stream huffman-size huffman-code) self
    (encode-dc-coefficient byte-stream dc-delta (svref zz 0) component) 
    (encode-ac-coefficients byte-stream zz (svref huffman-size
                                                  (if (zerop component) 4 5))
                            (svref huffman-code (if (zerop component) 4 5)))))

;;;
;;; writes the marker, which entails writing #xff followed by marker to the 
;;; stream
;;;
(defun write-marker (marker stream)
  (write-byte #xff stream)
  (write-byte marker stream))


;;;
;;; writes 2 bytes MSB first
;;;
(defun write-halfword (halfword stream)
  (write-byte (ldb (byte 8 8) halfword) stream)
  (write-byte (ldb (byte 8 0) halfword) stream))


;;;
;;; bytes is an array, and start and end are indexes into that array that
;;; specify the inclusive indicies that reference which bytes are to be
;;; written to the stram
;;;
(defun write-bytes (bytes start end stream)
  (loop for i from start to end
        for byte = (svref bytes i) do
        (write-byte byte stream)))


;;;
;;; At some point may want to add restart markers, especially for large files
;;; this is not requried by the spec, but could allow for faster reading of 
;;; these images by some programs if they have multiple processors working
;;;
(defmethod write-header ((self file) sampling)
  (with-slots (byte-stream) self
    ;;
    ;; Write the initial 22 bytes that define this to be a JFIF (JPEG) file
    ;;
    (write-marker *soi* byte-stream)
    (write-marker *app0* byte-stream)
    (write-halfword 16 byte-stream) ; length of header
    (write-byte (char-code #\J) byte-stream) ; write "JFIF"
    (write-byte (char-code #\F) byte-stream)
    (write-byte (char-code #\I) byte-stream)
    (write-byte (char-code #\F) byte-stream)
    (write-byte 0 byte-stream)
    (write-halfword (dpb 1 (byte 8 8) 2) byte-stream) ; version 1.2
    (write-byte 0 byte-stream) ;units
    (write-halfword 0 byte-stream) ; x-density
    (write-halfword 0 byte-stream) ; y-density
    (write-halfword 0 byte-stream) ; x,y thumbnail size
    (write-quantization-tables self)
    (write-huffman-tables self)
    (write-frame-header self sampling)
    (write-scan-header self)))

;;;
;;; writes out the frame header to the byte-stream
;;;
(defmethod write-frame-header ((self file) sampling)
  (with-slots (byte-stream components width height quantization-tables) self
    (let ((tables (count-quantization-tables quantization-tables))
          (max-h (loop for i in sampling maximize (first i)))
          (max-v (loop for i in sampling maximize (second i))))
      (write-marker *sof0* byte-stream)
      (write-halfword (+ 8 (* 3 components)) byte-stream) ; length
      (write-byte 8 byte-stream) ; precision i.e. 8-bit components
      (write-halfword height byte-stream)
      (write-halfword width byte-stream)
      (write-byte components byte-stream)
      (loop for i from 0 below components
            for sample in sampling do
            (write-byte (1+ i) byte-stream) ; component identifier
            (write-byte (dpb (/ max-h (first sample))
                             (byte 4 4)
                             (/ max-v (second sample))) byte-stream)
            (write-byte (get-quantization-index i tables) byte-stream)))))


;;;
;;; writes out the start of scan marker and header
;;;
(defmethod write-scan-header ((self file))
  (with-slots (byte-stream components) self
    (write-marker *sos* byte-stream)
    (write-halfword (+ 6 (* 2 components)) byte-stream) ; length
    (write-byte components byte-stream)
    ;;
    ;; write component specification parameters
    ;;
    (loop for i from 0 below components do
          (write-byte (1+ i) byte-stream) ; scan component selector
          (write-byte (if (zerop i)       ; DC/AC huffman-table selector
                          0
                        (dpb 1 (byte 4 4) 1))
                      byte-stream))
    (write-byte 0 byte-stream) ; SS
    (write-byte 63 byte-stream) ; SE
    (write-byte 0 byte-stream))) ; AH/AL


;;;
;;; writes out the quntization tables
;;;
(defmethod write-quantization-tables ((self file))
  (with-slots (byte-stream quantization-tables) self
    (loop for i from 0
          for table = (svref quantization-tables i)
          while table do
          (write-marker *dqt* byte-stream)
          (write-halfword 67 byte-stream) ; length
          (write-byte i byte-stream) ; identifier
          (loop for j from 0 to 63 do
                (write-byte (svref table (svref *zigzag* j)) byte-stream)))))
    

;;;
;;; for simplicity sake, and performance, we are going to use 4 fixed
;;; huffman tables, DC/AC luminance/chrominance
;;; 
(defmethod write-huffman-tables ((self file))
  (with-slots (byte-stream huffman-bits huffman-values) self
    ;;
    ;; write the DC luminance huffman table
    ;;
    (write-marker *dht* byte-stream)
    (write-halfword 31 byte-stream) ; length
    (write-byte 0 byte-stream) ; DC table index 0
    (write-bytes (svref huffman-bits 0) 1 16 byte-stream)
    (write-bytes (svref huffman-values 0) 0 11 byte-stream) ;huffval-table
    ;;
    ;; write the DC Chrominance huffman table
    ;;
    (write-marker *dht* byte-stream)
    (write-halfword 31 byte-stream) ; length
    (write-byte 1 byte-stream) ; DC table index 1
    (write-bytes (svref huffman-bits 1) 1 16 byte-stream)
    (write-bytes (svref huffman-values 1) 0 11 byte-stream) ;huffval-table
    ;;
    ;; write the AC luminance
    ;;
    (write-marker *dht* byte-stream)
    (write-halfword 181 byte-stream) ; length
    (write-byte (dpb 1 (byte 4 4) 0) byte-stream) ; AC table index 0
    (write-bytes (svref huffman-bits 4) 1 16 byte-stream)
    (write-bytes (svref huffman-values 4) 0 161 byte-stream) ; huffvals
    ;;
    ;; write the AC chrominance
    ;;
    (write-marker *dht* byte-stream)
    (write-halfword 181 byte-stream) ; length
    (write-byte (dpb 1 (byte 4 4) 1) byte-stream) ; AC table index 1
    (write-bytes (svref huffman-bits 5) 1 16 byte-stream)
    (write-bytes (svref huffman-values 5) 0 161 byte-stream)))


;;;
;;; zz and quant are both scaled by 10 on input, result is not scaled!
;;;
(defun quantize (zz quant)
  (declare #.*optimize*
           (type (simple-vector *) zz quant))
  (loop for i fixnum from 0 to 63 do
        (setf (svref zz i) (the fixnum (round (svref zz i) (svref quant i))))))


;;;
;;; converts the rgb image to YCbCr storing into dest
;;;
(defun convert-image-rgb-to-ycbcr (image ix iy ispan iheight dest width height
                                         alpha)
  (declare #.*optimize*
           (type (simple-vector *) image dest)
           (type fixnum width height ix iy ispan iheight))
  (if alpha
      (cond ((eq *component-order* :BGR)
             (write-components-ycbcr :B :G :R :A))
            ((eq *component-order* :BRG)
             (write-components-ycbcr :B :R :G :A))
            ((eq *component-order* :RGB)
             (write-components-ycbcr :R :G :B :A))
            ((eq *component-order* :RBG)
             (write-components-ycbcr :R :B :G :A))
            ((eq *component-order* :GRB)
             (write-components-ycbcr :G :R :B :A))
            ((eq *component-order* :GBR)
             (write-components-ycbcr :G :B :R :A))
            (t (error "Illegal component specified, must be one of: :BGR :BRG
                       :RGB :RBG :GRB :GBR")))
    (cond ((eq *component-order* :BGR)
           (write-components-ycbcr :B :G :R))
          ((eq *component-order* :BRG)
           (write-components-ycbcr :B :R :G))
          ((eq *component-order* :RGB)
           (write-components-ycbcr :R :G :B))
          ((eq *component-order* :RBG)
           (write-components-ycbcr :R :B :G))
          ((eq *component-order* :GRB)
           (write-components-ycbcr :G :R :B))
          ((eq *component-order* :GBR)
           (write-components-ycbcr :G :B :R))
          (t (error "Illegal component specified, must be one of: :BGR :BRG
                     :RGB :RBG :GRB :GBR")))))


;;;
;;;  Forward DCT - output is scaled by 10-bits
;;;
(defun forward-dct (zz)
  (declare #.*optimize*
           (type (simple-vector *) zz))
  ;;
  ;; Rows
  ;;
  (loop for i from 0 to 56 by 8
        for tmp0 fixnum = (i+ (svref zz (i+ i 0)) (svref zz (i+ i 7))) ; 9-bits
        for tmp1 fixnum = (i+ (svref zz (i+ i 1)) (svref zz (i+ i 6)))
        for tmp2 fixnum = (i+ (svref zz (i+ i 2)) (svref zz (i+ i 5)))
        for tmp3 fixnum = (i+ (svref zz (i+ i 3)) (svref zz (i+ i 4)))
        for tmp4 fixnum = (i- (svref zz (i+ i 3)) (svref zz (i+ i 4)))
        for tmp5 fixnum = (i- (svref zz (i+ i 2)) (svref zz (i+ i 5)))
        for tmp6 fixnum = (i- (svref zz (i+ i 1)) (svref zz (i+ i 6)))
        for tmp7 fixnum = (i- (svref zz (i+ i 0)) (svref zz (i+ i 7)))
        for tmp10 fixnum = (i+ tmp0 tmp3) ; 10-bits
        for tmp11 fixnum = (i+ tmp1 tmp2)
        for tmp12 fixnum = (i- tmp1 tmp2)
        for tmp13 fixnum = (i- tmp0 tmp3)
        for z1 fixnum = (i* (i+ tmp12 tmp13) (round (* 0.707106781 256)))
        with z2 with z3 with z4 with z5 with z6 with z7 do
        (setf (svref zz (i+ i 0)) (ash (i+ tmp10 tmp11) 8))
        (setf (svref zz (i+ i 2)) (i+ (ash tmp13 8) z1))
        (setf (svref zz (i+ i 4)) (ash (i- tmp10 tmp11) 8))
        (setf (svref zz (i+ i 6)) (i- (ash tmp13 8) z1)) 
        (setq tmp10 (i+ tmp4 tmp5))
        (setq tmp11 (i+ tmp5 tmp6))
        (setq tmp12 (i+ tmp6 tmp7))
        (setq z5 (i* (i- tmp10 tmp12) (round (* 0.382683433 256))))
        (setq z2 (i+ (ash (i* (round (* 0.541196100 512)) tmp10) -1) z5))
        (setq z4 (i+ (ash (i* (round (* 1.306562965 512)) tmp12) -1) z5))
        (setq z3 (i* tmp11 (round (* 0.707106781 256))))
        (setq z6 (i+ (ash tmp7 8) z3))
        (setq z7 (i- (ash tmp7 8) z3))
        (setf (svref zz (i+ i 1)) (i+ z6 z4))
        (setf (svref zz (i+ i 3)) (i- z7 z2))
        (setf (svref zz (i+ i 5)) (i+ z7 z2))
        (setf (svref zz (i+ i 7)) (i- z6 z4))) ; scaled by 8
  ;;
  ;; Columns - generally values are 20-bits scaled by 8
  ;;
  (loop for i from 0 to 7
        for tmp0 = (i+ (svref zz (i+ i (* 8 0))) (svref zz (i+ i (* 8 7))))
        for tmp1 = (i+ (svref zz (i+ i (* 8 1))) (svref zz (i+ i (* 8 6))))
        for tmp2 = (i+ (svref zz (i+ i (* 8 2))) (svref zz (i+ i (* 8 5))))
        for tmp3 = (i+ (svref zz (i+ i (* 8 3))) (svref zz (i+ i (* 8 4))))
        for tmp4 = (i- (svref zz (i+ i (* 8 3))) (svref zz (i+ i (* 8 4))))
        for tmp5 = (i- (svref zz (i+ i (* 8 2))) (svref zz (i+ i (* 8 5))))
        for tmp6 = (i- (svref zz (i+ i (* 8 1))) (svref zz (i+ i (* 8 6))))
        for tmp7 = (i- (svref zz (i+ i (* 8 0))) (svref zz (i+ i (* 8 7))))
        for tmp10 = (i+ tmp0 tmp3)
        for tmp11 = (i+ tmp1 tmp2)
        for tmp12 = (i- tmp1 tmp2)
        for tmp13 = (i- tmp0 tmp3) ; 22-bits
        for z1 = (i* (ash (i+ tmp12 tmp13) -8) (round (* 0.707106781 256)))
        with z2 with z3 with z4 with z5 with z6 with z7 do
        (setf (svref zz (i+ i (* 8 0))) (i+ tmp10 tmp11)) ; scaled by 8, but /
        (setf (svref zz (i+ i (* 8 2))) (i+ tmp13 z1))    ;  final result by 4
        (setf (svref zz (i+ i (* 8 4))) (i- tmp10 tmp11)) ;   means *really*
        (setf (svref zz (i+ i (* 8 6))) (i- tmp13 z1))    ;    scaled by 10 now!
        (setq tmp10 (i+ tmp4 tmp5))
        (setq tmp11 (i+ tmp5 tmp6))
        (setq tmp12 (i+ tmp6 tmp7)) ; 22-bits
        (setq z5 (i* (ash (i- tmp10 tmp12) -7) (round (* 0.382683433 128))))
        (setq z2 (i+ (i* (round (* 0.541196100 512)) (ash tmp10 -9)) z5))
        (setq z4 (i+ (i* (round (* 1.306562965 512)) (ash tmp12 -9)) z5))
        (setq z3 (i* (ash tmp11 -8) (round (* 0.707106781 256))))
        (setq z6 (i+ tmp7 z3))
        (setq z7 (i- tmp7 z3))
        (setf (svref zz (i+ i (* 8 1))) (i+ z6 z4))
        (setf (svref zz (i+ i (* 8 3))) (i- z7 z2))
        (setf (svref zz (i+ i (* 8 5))) (i+ z7 z2))
        (setf (svref zz (i+ i (* 8 7))) (i- z6 z4))))


#||

;;;
;;; Some debug support routines
;;;


(defun find-restart-marker (name)
  (with-open-file (stream name :element-type '(unsigned-byte 8))
    (loop with first = (read-byte stream nil nil)
          for second = (read-byte stream nil nil)
          while (and first second) do
          (when (and (= first #xff) (>= second *rst*) (<= second (+ 7 *rst*)))
            (format t " Found Restart")
            (return))
          (setq first second))))


(defun print-file (name)
  (with-open-file (stream name :element-type '(unsigned-byte 8))
    (format t "~%")
    (loop for byte = (read-byte stream nil nil)
          with count = 0
          while byte do
          (format t "~3x" byte)
          (when (= 32 (incf count))
            (setq count 0)
            (format t "~%")))))


(defmethod show-tables ((self file) tables)
  (loop for i from 0 below (length tables)
        for table = (svref tables i)
        when table do
        (format t "~%Table ~d:" i)
        (loop for j from 0 to 7 do
              (format t "~%")
              (loop for k from 0 to 7 do
                    (format t "~3d " (svref table (+ k (* j 8))))))))


(defmethod show-table ((self file) table)
  (loop with length = (length table)
        while (> length 0) do
        (loop for j from 0 to 7
              while (> length 0) do
              (format t "~%")
              (loop for k from 0 to 7
                    while (> length 0) do
                    (decf length)
                    (format t "~3d " (svref table (+ k (* j 8))))))))


(defmethod show-huffman-tables ((self file))
  (with-slots (huffman-maxcode huffman-mincode huffman-valptr huffman-size
                               huffman-code) self
    (loop for i from 0 to 7
          for maxcode = (svref huffman-maxcode i)
          for mincode = (svref huffman-mincode i)
          for valptr = (svref huffman-valptr i)
          for size = (svref huffman-size i)
          for code = (svref huffman-code i)
          when maxcode do
          (format t "~%~%Table ~d:~%Maxcode:" i)
          (loop for j from 1 to 16 do
                (format t "~5d " (svref maxcode j)))
          (format t "~%Mincode:")
          (loop for j from 1 to 16 do
                (format t "~5d " (svref mincode j)))
          (format t "~%Valptr: ")
          (loop for j from 1 to 16 do
                (format t "~5d " (svref valptr j)))
          (format t "~%Index Size Code")
          (loop for j from 0 to 255
                for code-size = (svref size j)
                when (> code-size 0) do
                (format t "~%~3d   ~3d  ~b" j code-size (svref code j))))))


(defun describe-file (name)
  (with-open-file (stream name :element-type '(unsigned-byte 8))
    (let ((file (make-instance 'file :byte-stream stream)))
      (setq *file* file)
      (if (null (read-header file))
          (format t " does not appear to be a valid JPEG file"
                  (file-byte-stream file)))
      (process-markers file)
      (format t "~a: ~dx~d"
              (pathname-name (truename (file-byte-stream *file*)))
              (file-width file) (file-height file))
      (describe-component-table file))))
  

(defmethod describe-component-table ((self file))
  (with-slots (component-tables) self
    (loop for i from 0 to 3
          for entry = (svref component-tables i)
          when entry do
          (format t "~%Index: ~a;  DC: ~a; AC: ~a; Quant: ~a; VS: ~a; HS: ~a"
                  i (ldb (byte 2 10) entry) (ldb (byte 2 8) entry)
                  (ldb (byte 2 6) entry)
                  (ldb (byte 3 3) entry) (ldb (byte 3 0) entry)))))


(defmethod show-image ((self file) x y w h)
  (with-slots (image width) self
    (loop for yy from y repeat h do
          (format t "~%")
          (loop for xx from x repeat w do
                (format t "~8x" (svref image (+ xx (* yy width))))))))


(defun show-8x8 (array)
  (loop for y from 0 to 7 do
        (format t "~%")
        (loop for x from 0 to 7 do
              (format t "~4x" (svref array (+ x (* y 8)))))))


(defun show-quant (index)
  (let ((table (svref (file-quantization-tables *file*) index)))
    (loop for y from 0 to 7 do
          (format t "~%")
          (loop for x from 0 to 7 do
                (format t "~3x" (svref table (+ x (* y 8))))))))

(defun ycbcr-rgb (y cb cr)
  (list (int-clamp (round (+ y (* 1.402 (- cr 127.5)))))
        (int-clamp (round (- (- y (* .34414 (- cb 127.5)))
                             (* .71414 (- cr 127.5)))))
        (int-clamp (round (+ y (* 1.772 (- cb 127.5)))))))


(defun rgb-ycbcr (r g b)
  (list (int-clamp (round (+ (* .299 r) (* .587 g) (* .114 b))))
        (int-clamp (round (- (- (+ 127.5 (* .5 b)) (* .3313 g)) (* .1687 r))))
        (int-clamp (round (- (- (+ 127.5 (* .5 r)) (* .4187 g)) (* .0813 b))))))



;;; 
;;; The following code is LWW 4.0.1 implementation specific code for displaying
;;; the JPEG image read in the above implementation independent code.
;;;
;;; (display-jpeg-file file)
;;;


(export '(display-jpeg-file))


(defvar *progress* nil "To be used in progress-step")


(defclass color-pane (capi:output-pane capi:color-screen)
  ((screen :initform nil
           :initarg :screen
           :accessor color-screen)
   (width :initform nil
          :initarg :width
          :accessor color-pane-width)
   (height :initform nil
           :initarg :height
           :accessor color-pane-height))
  (:default-initargs
   :display-callback 'display-screen
   :enabled nil
   :horizontal-scroll t
   :vertical-scroll t))


(defmethod display-screen ((self color-pane) x y w h)
  (with-slots (screen width height) self
    (when (null screen)
      (setq screen (gp:create-pixmap-port self width height)))
    (gp:copy-pixels self screen x y w h x y)))


(capi:define-interface jpeg-reader-progress ()
  ()
  (:panes
   (progress-bar capi:progress-bar
                 :reader progress-bar-pane
                 :height 26
                 :width 550
                 :x 10
                 :y 20)
   (percent-done-output capi:title-pane
                        :reader percent-done-output-pane
                        :height 26
                        :width 20
                        :x 570
                        :y 20)
   (progress-output capi:title-pane
                    :reader progress-output-pane
                    :height 26
                    :width 550
                    :x 10
                    :y 50))
  (:layouts
   (layout-with-progress
    capi:pinboard-layout
    '(progress-bar percent-done-output progress-output)
    :height 90
    :min-height 90
    :min-width 600
    :max-height 90
    :max-width 600
    :width 600
    :x 10
    :y 10))
  (:default-initargs
   :auto-menus nil
   :best-height 90
   :best-width 600
   :x 10
   :y 10
   :layout 'layout-with-progress
   :title "Read JPEG File"))


(defun display-jpeg-file (name)
  (let ((*progress* (capi:display (make-instance 'jpeg-reader-progress))))
    (unwind-protect
        (multiple-value-bind (width height) (get-image-size name)
          (let ((image-array (make-array (* width height) :initial-element 0)))
            (multiple-value-bind (ignore width height)
                (read-file name
                           :progress-function 'progress-step
                           :image image-array)
              ignore
              (let* ((pane (capi:contain (make-instance 'color-pane
                                                        :width width
                                                        :height height)
                                         :title (file-namestring name)
                                         :width (min width 640)
                                         :height (min height 480)
                                         :auto-menus nil
                                         :x 10
                                         :y 100))
                     (screen (loop until (color-screen pane) do (sleep 2)
                                   finally (return (color-screen pane)))))
                (loop for y from 0 below height
                      for ybase from 0 by width do
                      (progress-step (round (* 100 y) height) "Displaying Image")
                      (loop for x from 0 below width do
                            (gp:draw-point screen x y :foreground
                                           (svref image-array (+ x ybase)))))
                (capi:set-scroll-range pane width height)
                (display-screen pane 0 0 (min width 640) (min height 480))))
            (mp:process-kill
             (capi-internals:interface-process *progress*)))))))


(defun progress-step (step-number text)
  (setf (capi:range-slug-start (progress-bar-pane *progress*)) step-number)
  (setf (capi:title-pane-text (progress-output-pane *progress*)) text)
  (setf (capi:title-pane-text (percent-done-output-pane *progress*))
        (format nil "~d%" step-number)))

||#