;;;-*- Mode:LISP; Package:NEW-MATH; Base:10.; Readtable:CL -*-

;;Code overhaul performed by WKF from 5/20 to 5/26/88.

;**************************************************
;*           Floating point functions		  *
;**************************************************

(defmacro float-error-p (status)                   ;;WKF added 5/19/88 to work around defconstant limitations
  `(hw:32logbitp (hw:32logical-shift-up  1. ,status)
		 ,(hw:unboxed-constant  ;;+++ Turn this into a global register.
		    (lisp:logior (lisp:ash 1. (1+ hw:$$fpu-exact-infinity))
				 (lisp:ash 1. (1+ hw:$$fpu-overflow-inexact))
				 (lisp:ash 1. (1+ hw:$$fpu-x-denormalized))
				 (lisp:ash 1. (1+ hw:$$fpu-y-denormalized))
				 (lisp:ash 1. (1+ hw:$$fpu-xy-denormalized))
				 (lisp:ash 1. (1+ hw:$$fpu-zero-divide))
				 (lisp:ash 1. (1+ hw:$$fpu-x-nan))
				 (lisp:ash 1. (1+ hw:$$fpu-y-nan))
				 (lisp:ash 1. (1+ hw:$$fpu-xy-nan))
				 (lisp:ash 1. (1+ hw:$$fpu-invalid-operation))))))

(defun float-error-message (status)
  (li:error "Floating point operation error"
	    (dispatch (byte 4 0) status
		      (hw:$$fpu-exact-zero        "Result =  0, exact")
		      (hw:$$fpu-exact-infinity    "Result =  infinity, exact")
		      (hw:$$fpu-finite-exact      "Result finite and non-zero, exact")
		      (hw:$$fpu-finite-inexact    "Result finite and non-zero, inexact")
		      (hw:$$fpu-overflow-inexact  "Overflow, inexact")
		      (hw:$$fpu-underflow-exact   "Underflow, exact")
		      (hw:$$fpu-underflow-inexact "Underflow, inexact")
		      (hw:$$fpu-x-denormalized    "Denormalized operand X")
		      (hw:$$fpu-y-denormalized    "Denormalized operand Y")
		      (hw:$$fpu-xy-denormalized   "X & Y operands denormalized")
		      (hw:$$fpu-zero-divide       "Zero divide")
		      (hw:$$fpu-X-NAN             "X operand Not-A-Number")
		      (hw:$$fpu-Y-NAN             "Y operand Not-A-Number")
		      (hw:$$fpu-XY-NAN            "X & Y operands Not-A-Number")
		      (hw:$$fpu-invalid-operation "Invalid operation requested"))))
     
(defmacro fpu-status ()
#+unhandled-by-fleabit  '(the (type (integer 0 15))  ;;@@@ Compiler should be able to determine this itself from byte-spec.
	(hw:field-extract (hw:read-processor-status)
		     hw:%%processor-status-floating-point-status))
     '(hw:field-extract (hw:read-processor-status)
			hw:%%processor-status-floating-point-status))


;;+++  When does the weitek produce negative zero?  This should be looked into.  --wkf
;;+++  Currently, 7/3/88, Negative zero is setting both zero and negative status bits.

;**************************************************
;*     Weitek floating point chip initializer     *
;**************************************************

(defafun init-float ()
  (falu fpu-default-mode0 nop ignore ignore fpu-load-mode fpu-unload-high)
  (falu fpu-default-mode1 nop ignore ignore fpu-load-mode fpu-unload-high)
  (falu fpu-default-mode2 nop ignore ignore fpu-load-mode fpu-unload-high)
  (falu fpu-default-mode3 nop ignore ignore fpu-load-mode fpu-unload-high)
  (fmul fpu-default-mode0 nop ignore ignore fpu-load-mode fpu-unload-high)
  (fmul fpu-default-mode1 nop ignore ignore fpu-load-mode fpu-unload-high)
  (fmul fpu-default-mode2 nop ignore ignore fpu-load-mode fpu-unload-high)
  (fmul fpu-default-mode3 nop ignore ignore fpu-load-mode fpu-unload-high)
  (return gr:*nil* boxed-right))

;**************************************************
;*               Short precision                  *
;**************************************************

(defconstant %%dpb-short-in-single (byte 26. 6.))
(defconstant %%dpb-short-round-bit (byte  1  5.))

(defun add-short (xx yy)
  (let* ((x      (hw:dpb-unboxed xx %%dpb-short-in-single gr:*all-zero*))
	 (y      (hw:dpb-unboxed yy %%dpb-short-in-single gr:*all-zero*))
	 (tflag  (hw:trap-off))
	 (result (hw:float-add-single x y))
	 (status (hw:read-processor-status)))
    (trap:%trap-restore tflag)  ;;WARNING: No trapping or functional source frobbing for 4 cycles!!!  -wkf
    (let ((status (hw:field-extract status hw:%%processor-status-floating-point-status)))
      (when (float-error-p status)
	(float-error-message status)))
    (fixup-short result)))

(defun subtract-short (xx yy)
  (let* ((x      (hw:dpb-unboxed xx %%dpb-short-in-single gr:*all-zero*))
	 (y      (hw:dpb-unboxed yy %%dpb-short-in-single gr:*all-zero*))
	 (tflag  (hw:trap-off))
	 (result (hw:float-subtract-single x y))
	 (status (hw:read-processor-status)))
    (trap:%trap-restore tflag)  ;;WARNING: No trapping or functional source frobbing for 4 cycles!!!  -wkf
    (let ((status (hw:field-extract status hw:%%processor-status-floating-point-status)))
      (when (float-error-p status)
	(float-error-message status)))
    (fixup-short result)))

(defun multiply-short (xx yy)
  (let* ((x      (hw:dpb-unboxed xx %%dpb-short-in-single gr:*all-zero*))
	 (y      (hw:dpb-unboxed yy %%dpb-short-in-single gr:*all-zero*))
	 (tflag  (hw:trap-off))
	 (result (hw:float-multiply-single x y))
	 (status (hw:read-processor-status)))
    (trap:%trap-restore tflag)  ;;WARNING: No trapping or functional source frobbing for 4 cycles!!!  -wkf
    (let ((status (hw:field-extract status hw:%%processor-status-floating-point-status)))
      (when (float-error-p status)
	(float-error-message status)))
    (fixup-short result)))

(defun divide-short (xx yy)
  (li:error "Currently divide wedges the K unrecoverably!!!")
  (let* ((x      (hw:dpb-unboxed xx %%dpb-short-in-single gr:*all-zero*))
	 (y      (hw:dpb-unboxed yy %%dpb-short-in-single gr:*all-zero*))
	 (tflag  (hw:trap-off))
	 (result (hw:float-divide-single x y))
	 (status (hw:read-processor-status)))
    (trap:%trap-restore tflag)  ;;WARNING: No trapping or functional source frobbing for 4 cycles!!!  -wkf
    (let ((status (hw:field-extract status hw:%%processor-status-floating-point-status)))
      (when (float-error-p status)
	(float-error-message status)))
    (fixup-short result)))

(defun fixup-short (result)
  (setq result
	(if (hw:32logbitp (byte-position %%dpb-short-round-bit) result)
	    (hw:32-1+ (hw:ldb-boxed result %%dpb-short-in-single (hw:unboxed-constant (ash vinc:$$dtp-short-float 26.))))
	  (hw:ldb-boxed result %%dpb-short-in-single (hw:unboxed-constant (ash vinc:$$dtp-short-float 26.)))))
        ;; We can round up by adding one.  If this overflows the mantissa then the exponent needs the overflow.
        ;; Exponent won't overflow into the sign since legal exponents are 1 up to 254.
        ;;+++ Handle denormalized numbers correctly.  Carry into exponent needs floating point excess added back in.
  (values result (short-status result)))

(defun compare-short (xx yy)
  (let* ((x      (hw:dpb-unboxed xx %%dpb-short-in-single gr:*all-zero*))
	 (y      (hw:dpb-unboxed yy %%dpb-short-in-single gr:*all-zero*))
	 (tflag  (hw:trap-off))
	 (result (hw:float-compare-single x y))
	 (status (hw:read-processor-status)))
    (trap:%trap-restore tflag)
    (dispatch (byte 2. 4.) status ;;(byte 2. 0) of hw:%%processor-status-floating-point-status = (byte 2 4)
	                          ;; @@@ compiler should do this optimization
       (hw:$$fpu-equal        (values  0 hw:$$alu-status-equal))
       (hw:$$fpu-less-than    (values -1 hw:$$alu-status-less-than))
       (hw:$$fpu-greater-than (values  1 hw:$$alu-status-greater-than))
       (hw:$$fpu-unordered    (progn (li:error "Short float compare, Unordered") nil)))));;@@@Have compiler make not tail call

(defmacro short-status (short)
  `(hw:dpb (hw:ldb ,short hw:%%short-float-sign gr:*all-zero*) hw:%%alu-status-negative
	   (progn (setf (hw:r0) ,short)   ;;@@@ Figure out how to remove this instruction.
		  (hw:read-alu-status)))) ;;pretend we have a fixnum, test for zero (BW-24).

#+never
(defun short-status (short)
  (cond ((hw:field= short gr:*all-zero* %%pointer)
	 hw:$$alu-status-positive-zero)
	((zerop (hw:ldb-boxed short hw:%%short-float-sign 0))
	 hw:$$alu-status-positive)
	((hw:field= short gr:*all-zero* hw:%%short-float-exponent-and-mantissa)
	 hw:$$alu-status-negative-zero)
	(t hw:$$alu-status-negative)))

(defun test-short (xx)
  (let ((status (short-status xx)))
    (values status status)))

#+never
(defun test-short (xx)
  (let ((status (hw:dpb (hw:ldb xx hw:%%short-float-sign 0) hw:%%alu-status-negative
			(if (hw:24= xx gr:*all-zero*)
			    hw:$$alu-status-zero
			  0))))
    (values status status)))

#+never
(defun test-short (xx)
  (if (zerop (hw:ldb xx hw:%%short-float-sign 0))
      (if (hw:24= xx gr:*all-zero*)
	  (values hw:$$alu-status-positive-zero hw:$$alu-status-positive-zero)
	(values hw:$$alu-status-positive hw:$$alu-status-positive))
    (if (hw:24= xx gr:*all-zero*)
	(values hw:$$alu-status-negative-zero hw:$$alu-status-negative-zero)
      (values hw:$$alu-status-negative hw:$$alu-status-negative))))

(defun negate-short (xx)
  (let ((result (hw:dpb-xor 1 hw:%%short-float-sign xx)))
    (values result (short-status result))))

;**************************************************
;*               Single precision                 *
;**************************************************

(defun add-single (xx yy)
  (let* ((x      (array:%vm-read32 xx 1))
	 (y      (array:%vm-read32 yy 1))
	 (tflag  (hw:trap-off))
	 (result (hw:float-add-single x y))
	 (status (fpu-status)))
    (trap:trap-restore tflag)
    (fixup-single result status)))

(defun subtract-single (xx yy)
  (let* ((x      (array:%vm-read32 xx 1))
	 (y      (array:%vm-read32 yy 1))
	 (tflag  (hw:trap-off))
	 (result (hw:float-subtract-single x y))
	 (status (fpu-status)))
    (trap:trap-restore tflag)
    (fixup-single result status)))

(defun multiply-single (xx yy)
  (let* ((x      (array:%vm-read32 xx 1))
	 (y      (array:%vm-read32 yy 1))
	 (tflag  (hw:trap-off))
	 (result (hw:float-multiply-single x y))
	 (status (fpu-status)))
    (trap:trap-restore tflag)
    (fixup-single result status)))

(defun divide-single (xx yy)
  (let* ((x      (array:%vm-read32 xx 1))
	 (y      (array:%vm-read32 yy 1))
	 (tflag  (hw:trap-off))
	 (result (hw:float-divide-single x y))
	 (status (fpu-status)))
    (trap:trap-restore tflag)
    (fixup-single result status)))

(defun fixup-single (result status)
  (when (float-error-p status)
    (li:error "Single float operation" (float-error-message status)))
  (values
    (array:make-single-float result)
    (single-status           result)))

(defun compare-single (xx yy)
  (let* ((x      (array:%vm-read32 xx 1))
	 (y      (array:%vm-read32 yy 1))
	 (tflag  (hw:trap-off))
	 (result (hw:float-compare-single x y))
	 (status (fpu-status)))
    (trap:trap-restore tflag)
    (dispatch (byte 2 0) status
       (hw:$$fpu-equal        (values  0 hw:$$alu-status-equal))
       (hw:$$fpu-less-than    (values -1 hw:$$alu-status-less-than))
       (hw:$$fpu-greater-than (values  1 hw:$$alu-status-greater-than))
       (hw:$$fpu-unordered    (li:tail-error "Single float compare, Unordered" status)))))

;;@@@ Turn single-status into a macro for speed?  --wkf
(defun single-status (single-word)
  (cond ((hw:32zerop single-word)
	 hw:$$alu-status-positive-zero)
	((zerop (hw:ldb-boxed single-word hw:%%single-float-sign 0))
	 hw:$$alu-status-positive)
	((hw:32zerop (hw:ldb single-word hw:%%single-float-exponent-and-mantissa gr:*all-zero*))
	 hw:$$alu-status-negative-zero)
	(t hw:$$alu-status-negative)))

(defun test-single (xx)
  (let ((status (single-status (array:%vm-read32 xx 1))))
    (values status status)))

(defun negate-single (xx)
  (let ((result (hw:dpb-xor 1 hw:%%single-float-sign (array:%vm-read-32 xx 1))))
    (values (array:make-single-float result) (single-status result))))

;**************************************************
;*               Double precision                 *
;**************************************************

(defun add-double (xx yy)
  (let ((x-lo  (array:%vm-read32 xx 1))
	(x-hi  (array:%vm-read32 xx 2))
	(y-lo  (array:%vm-read32 yy 1))
	(y-hi  (array:%vm-read32 yy 2))
	(tflag (hw:trap-off)))
    (multiple-value-bind (result-hi result-lo)
	(hw:float-add-double x-hi x-lo y-hi y-lo)
      (let ((status (fpu-status)))
	(trap:trap-restore tflag)
	(fixup-double result-hi result-lo status)))))

(defun subtract-double (xx yy)
  (let ((x-lo  (array:%vm-read32 xx 1))
	(x-hi  (array:%vm-read32 xx 2))
	(y-lo  (array:%vm-read32 yy 1))
	(y-hi  (array:%vm-read32 yy 2))
	(tflag (hw:trap-off)))
    (multiple-value-bind (result-hi result-lo)
	(hw:float-subtract-double x-hi x-lo y-hi y-lo)
      (let ((status (fpu-status)))
	(trap:trap-restore tflag)
	(fixup-double result-hi result-lo status)))))

(defun multiply-double (xx yy)
  (let ((x-lo  (array:%vm-read32 xx 1))
	(x-hi  (array:%vm-read32 xx 2))
	(y-lo  (array:%vm-read32 yy 1))
	(y-hi  (array:%vm-read32 yy 2))
	(tflag (hw:trap-off)))
    (multiple-value-bind (result-hi result-lo)
	(hw:float-multiply-double x-hi x-lo y-hi y-lo)
      (let ((status (fpu-status)))
	(trap:trap-restore tflag)
	(fixup-double result-hi result-lo status)))))

(defun divide-double (xx yy)
  (let ((x-lo  (array:%vm-read32 xx 1))
	(x-hi  (array:%vm-read32 xx 2))
	(y-lo  (array:%vm-read32 yy 1))
	(y-hi  (array:%vm-read32 yy 2))
	(tflag (hw:trap-off)))
    (multiple-value-bind (result-hi result-lo)
	(hw:float-divide-double x-hi x-lo y-hi y-lo)
      (let ((status (fpu-status)))
	(trap:trap-restore tflag)
	(fixup-double result-hi result-lo status)))))

(defun fixup-double (result-hi result-lo status)
  (when (float-error-p status)
    (li:error "Double float operation" (float-error-message status)))
  (values
    (array:make-double-float result-hi result-lo)
    (double-status           result-hi result-lo)))

(defun compare-double (xx yy)
  (let ((x-lo  (array:%vm-read32 xx 1))
	(x-hi  (array:%vm-read32 xx 2))
	(y-lo  (array:%vm-read32 yy 1))
	(y-hi  (array:%vm-read32 yy 2))
	(tflag (hw:trap-off)))
    (let ((ignore (hw:float-compare-double x-hi x-lo y-hi y-lo))
	  (status (fpu-status)))
      (trap:trap-restore tflag)
      (dispatch (byte 2 0) status
	 (hw:$$fpu-equal        (values  0 hw:$$alu-status-equal))
	 (hw:$$fpu-less-than    (values -1 hw:$$alu-status-less-than))
	 (hw:$$fpu-greater-than (values  1 hw:$$alu-status-greater-than))
	 (hw:$$fpu-unordered    (li:tail-error "Double float compare, Unordered" status))))))

(defun double-status (double-word-high double-word-low)
  (cond ((hw:32zerop (hw:32logior double-word-high double-word-low))
	 hw:$$alu-status-positive-zero)
	((zerop (hw:ldb-boxed double-word-high hw:%%double-float-sign-word2 0))
	 hw:$$alu-status-positive)
	((hw:32zerop (logior double-word-low
			     (hw:ldb double-word-high hw:%%double-float-exponent-and-mantissa-word2 gr:*all-zero*)))
	 hw:$$alu-status-negative-zero)
	(t hw:$$alu-status-negative)))

(defun test-double (xx)
  (let* ((x-hi (array:%vm-read32 xx 1))
	 (x-lo (array:%vm-read32 xx 2))
	 (status (double-status x-hi x-lo)))
    (values status status)))

(defun negate-double (xx)
  (let ((x-lo   (array:%vm-read32 xx 1))
	(x-hi   (hw:dpb-xor 1 hw:%%double-float-sign-word2 (array:%vm-read32 xx 2))))
    (values (array:make-double-float x-hi x-lo)
	    (double-status           x-hi x-lo))))

;**************************************************
;*            Convert fix to short                *
;**************************************************

(defafun convert-fixnum-to-short (i)
  (alu sex-r a1 a0 a0 bw-24 unboxed dt-both-fixnum)
  (move a2 trap-off unboxed)
  (falu single-float nop ignore a1 fpu-load-x   fpu-unload-high unboxed)
  (falu single-float nop ignore a1 fpu-load-nop fpu-unload-high unboxed)
  (falu single-float nop ignore a1 fpu-load-nop fpu-unload-high unboxed)
  (falu single-float nop ignore a1 fpu-load-nop fpu-unload-high unboxed)
  (falu single-float nop ignore a1 fpu-load-nop fpu-unload-high unboxed)
  (falu single-float a3  ignore a1 fpu-load-nop fpu-unload-high unboxed)
  (alu-field field-pass memory-control a2 memory-control hw:%%single-float-sign unboxed)
  (movei a4 #.(ash vinc:$$dtp-short-float 26.) unboxed)
  (alu-field field-pass return a3 a4 (byte 26. -6.) boxed ch-return next-pc-return)
 )

;**************************************************
;*            Convert fix to single               *
;**************************************************

(defafun convert-fixnum-to-single (i)
  (alu sex-r a1 ignore a0 bw-24 unboxed dt-both-fixnum)
  (move a2 trap-off unboxed)
  (falu single-float nop ignore a1 fpu-load-x   fpu-unload-high unboxed)
  (falu single-float nop ignore a1 fpu-load-nop fpu-unload-high unboxed)
  (falu single-float nop ignore a1 fpu-load-nop fpu-unload-high unboxed)
  (falu single-float nop ignore a1 fpu-load-nop fpu-unload-high unboxed)
  (falu single-float nop ignore a1 fpu-load-nop fpu-unload-high unboxed)
  (falu single-float a3  ignore a1 fpu-load-nop fpu-unload-high unboxed)
  (alu-field field-pass memory-control a2 memory-control hw:%%single-float-sign unboxed)
  (tail-open-call (array:make-single-float 1) (o0 a3))
 )

;**************************************************
;*             Convert fix to double              *
;**************************************************

(defafun convert-fixnum-to-double (i)
  (alu sex-r a1 ignore a0 bw-24 unboxed)
  (move a2 trap-off unboxed)
  (falu double-float nop ignore a1 fpu-load-x   fpu-unload-high unboxed)
  (falu double-float nop ignore a1 fpu-load-nop fpu-unload-high unboxed)
  (falu double-float nop ignore a1 fpu-load-nop fpu-unload-high unboxed)
  (falu double-float nop ignore a1 fpu-load-nop fpu-unload-high unboxed)
  (falu double-float nop ignore a1 fpu-load-nop fpu-unload-low  unboxed)
  (falu double-float o0  ignore a1 fpu-load-nop fpu-unload-low  unboxed ch-tail-open)
  (falu double-float o1  ignore a1 fpu-load-nop fpu-unload-low  unboxed)
  (alu-field field-pass memory-control a2 memory-control hw:%%double-float-sign-word2 unboxed)
  (tail-call (array:make-double-float 2) ())
 )

;**************************************************
;*           Convert short to single              *
;**************************************************

(defun convert-short-to-single (x)
  (array:make-single-float (hw:dpb-unboxed x %%dpb-short-in-single gr:*all-zero*)))

;**************************************************
;*           Convert short to double              *
;**************************************************

(defafun convert-short-to-double (x)
  (alu-field field-extract-r a1 ignore a0 %%dpb-short-in-single unboxed-md)
  (move a2 trap-off unboxed)
  (falu single-to-double nop a1 ignore fpu-load-x   fpu-unload-high unboxed)
  (falu single-to-double nop a1 ignore fpu-load-nop fpu-unload-high unboxed)
  (falu single-to-double nop a1 ignore fpu-load-nop fpu-unload-high unboxed)
  (falu single-to-double nop a1 ignore fpu-load-nop fpu-unload-high unboxed)
  (falu single-to-double nop a1 ignore fpu-load-nop fpu-unload-low  unboxed)
  (falu single-to-double o0  a1 ignore fpu-load-nop fpu-unload-low  unboxed ch-tail-open)
  (falu single-to-double o1  a1 ignore fpu-load-nop fpu-unload-low  unboxed)
  (alu-field field-pass memory-control a3 memory-control hw:%%double-float-sign-word2 unboxed)
  (tail-call (array:make-double-float 2) ())
 )

;**************************************************
;*           Convert single to double             *
;**************************************************

(defafun convert-single-to-double (x)
  (alu r+1 vma-start-read-no-transport ignore a0 unboxed-vma unboxed-md)
  (nop)
  (move a1 md unboxed)
  (move a2 trap-off unboxed)
  (falu single-to-double nop a1 ignore fpu-load-x   fpu-unload-high unboxed)
  (falu single-to-double nop a1 ignore fpu-load-nop fpu-unload-high unboxed)
  (falu single-to-double nop a1 ignore fpu-load-nop fpu-unload-high unboxed)
  (falu single-to-double nop a1 ignore fpu-load-nop fpu-unload-high unboxed)
  (falu single-to-double nop a1 ignore fpu-load-nop fpu-unload-low  unboxed)
  (falu single-to-double o0  a1 ignore fpu-load-nop fpu-unload-low  unboxed ch-tail-open)
  (falu single-to-double o1  a1 ignore fpu-load-nop fpu-unload-low  unboxed)
  (alu-field field-pass memory-control a3 memory-control hw:%%double-float-sign-word2 unboxed)
  (tail-call (array:make-double-float 2) ())
 )

;**************************************************
;*               SIGN-VALUE                       *
;**************************************************
;;; SIGN-VALUE dispatches to these depending on the type of float.
;; Returns one value:
;;    - a floating point 1.0 of the same format with the same sign as the argument supplied.
;;
;; These functions should later be made into macros for speed. @@@

(defun sign-value-short (short)
  (if (hw:field= short gr:*all-zero* hw:%%short-float-sign)
      1.0s0
    -1.0s0))

(defun sign-value-single (single-word)
  (if (hw:field= single-word gr:*all-zero* hw:%%single-float-sign)
      1.0
    -1.0))

(defun sign-value-double (double-word2)
  (if (hw:field= double-word2 gr:*all-zero* hw:%%double-float-sign-word2)
      1.0d0
    -1.0d0))

(defun SIGN-VALUE (float)
  (prims:dispatch vinc::%%data-type float
     ($$dtp-short-float     (sign-value-short  float))
     ($$dtp-single-float    (sign-value-single (array::%vm-read32 float 1)))
     ($$dtp-double-float    (sign-value-double (array::%vm-read32 float 2)))
     (t (li:tail-error "~d is not a floating point number" float))))

;**************************************************
;*               ZERO-FLOATP                      *
;**************************************************
;;; ZERO-FLOATP dispatches to these depending on the type of float suppled.
;; Returns one value:
;;    - T if argument is a floating point zero, else nil.
;;
;; These functions should later be made into macros for speed. @@@

(defun zero-floatp-short (short)
  (hw:field= short gr:*all-zero* hw:%%short-float-exponent-and-mantissa))

(defun zero-floatp-single (single-word)
  (hw:field= single-word gr:*all-zero*  hw:%%single-float-exponent-and-mantissa))

(defun zero-floatp-double (double-word1 double-word2)
  (and (hw:field= double-word1 gr:*all-zero* hw:%%double-float-exponent-and-mantissa-word1)
       (hw:field= double-word2 gr:*all-zero* hw:%%double-float-exponent-and-mantissa-word2)))

(defun ZERO-FLOATP (float)
  (prims:dispatch vinc::%%data-type float
     ($$dtp-short-float     (zero-floatp-short  float))
     ($$dtp-single-float    (zero-floatp-single (array::%vm-read32 float 1)))
     ($$dtp-double-float    (zero-floatp-double (array::%vm-read32 float 1) (array::%vm-read32 float 2)))
     (t (li:tail-error "~d is not a floating point number" float))))

;**************************************************
;*              SCALE-MANTISSA                    *
;**************************************************
;;; SCALE-MANTISSA dispatches to these depending on the type of float suppled.
;; Returns one value:
;;    - a floating point number of the same format as the argument but
;;      between zero inclusive and one exclusive [0,1) (a negative one exponent)
;;
;; These functions should later be made into macros for speed.  @@@

(defun scale-mantissa-short (short)
  (hw:dpb-aligned short hw:%%short-float-mantissa 1.0s-1))

#+This-breaks-fleabit	 (hw:dpb (1- hw:$$short-float-exponent-excess)
				 hw:%%short-float-exponent 0.0s0)

(defun scale-mantissa-single (single-word)
  (array:make-single-float (hw:dpb-aligned single-word hw:%%single-float-mantissa
					   (hw:dpb-unboxed (1- hw:$$single-float-exponent-excess)
							   hw:%%single-float-exponent gr:*all-zero*))))

(defun scale-mantissa-double (double-word1 double-word2)
  (array:make-double-float (hw:dpb-aligned double-word2 hw:%%double-float-mantissa-word2
					   (hw:dpb-unboxed (1- hw:$$double-float-exponent-excess)
							   hw:%%double-float-exponent-word2
							   gr:*all-zero*))
			   double-word1))

(defun SCALE-MANTISSA (float)
  (prims:dispatch vinc::%%data-type float
     ($$dtp-short-float     (scale-mantissa-short  float))
     ($$dtp-single-float    (scale-mantissa-single (array::%vm-read32 float 1)))
     ($$dtp-double-float    (scale-mantissa-double (array::%vm-read32 float 1) (array::%vm-read32 float 2)))
     (t (li:tail-error "~d is not a floating point number" float))))

;**************************************************
;*              FIND-MANTISSA                     *
;**************************************************
;;; FIND-MANTISSA dispatches to these depending on the type of float suppled.
;; Returns one value:
;;    - an integer representing the mantissa
;;
;; These functions should later be made into macros for speed. @@@

(defun find-mantissa-short (short)
  (hw:ldb short hw:%%short-float-mantissa
	  (hw:dpb 1 (byte 1 (byte-size hw:%%short-float-mantissa))
		  0)))

(defun find-mantissa-single (single-word)
  (let ((bignum (allocate-bignum 1))) ;;This is at least one bigger than most-positive-fixnum
    (array::%vm-write32 bignum 1
			(hw:ldb single-word hw:%%single-float-mantissa
				(hw:dpb-unboxed 1 (byte 1 (byte-size hw:%%single-float-mantissa))
						gr:*all-zero*)))
    bignum))

(defun find-mantissa-single-and-shift-with-sign (single-word shift sign) ;;shift must be positive
  (let* ((mant-size  (1+ (byte-size hw:%%single-float-mantissa)))
	 (total-size (+ shift mant-size))
	 (words      (ceiling (1+ total-size) 32.))
	 (bignum     (allocate-bignum words))  ;;This code assumes bignum is initialized to all zeros.
	 (mantissa   (hw:ldb single-word hw:%%single-float-mantissa
			     (hw:dpb-unboxed 1 (byte 1 (byte-size hw:%%single-float-mantissa))
					     gr:*all-zero*)))
	 (plus-sign  (plusp sign))
	 (signed-mant (if plus-sign
			  mantissa
			(hw:32+ (hw:32logxor mantissa gr:*all-ones*) (hw:unboxed-constant 1)))))
    (cond (( words (ceiling total-size 32.))
	   (array::%vm-write32 bignum (1- words)
			       (hw:dpb signed-mant (byte mant-size (- 32. mant-size)) gr:*all-zero*))
	   (unless plus-sign
	     (array::%vm-write32 bignum words gr:*all-ones*)))
	  (t (let* ((high-bits (- total-size (* 32. (1- words))))
		    (low-bits  (- mant-size high-bits)))
	       (cond ((not (plusp low-bits))
		      (let ((high-word (hw:dpb-unboxed signed-mant (byte mant-size (- low-bits)) gr:*all-zero*)))
			(unless plus-sign
			  (setq high-word (hw:dpb-unboxed gr:*all-ones* (byte (- 32. mant-size) mant-size) high-word)))
			(array::%vm-write32 bignum words high-word)))
		     (t (let ((high-word (hw:ldb signed-mant (byte high-bits low-bits) gr:*all-zero*)))
			  (unless plus-sign
			    (setq high-word (hw:dpb-unboxed gr:*all-ones* (byte (- 32. high-bits) high-bits) high-word)))
			  (array::%vm-write32 bignum words high-word)
			  (array::%vm-write32 bignum (1- words)
					      (hw:dpb-unboxed signed-mant (byte low-bits (- 32. low-bits))
							      gr:*all-zero*))))))))))

(defun find-mantissa-double (double-word1 double-word2)
  (let ((bignum (allocate-bignum 2)))
    (array::%vm-write32 bignum 1 double-word1)
    (array::%vm-write32 bignum 2
			(hw:ldb double-word2 hw:%%double-float-mantissa-word2
				(hw:dpb-unboxed 1 (byte 1 (byte-size hw:%%double-float-mantissa-word2))
						gr:*all-zero*)))
    bignum))

(defconstant mant-size-const (1+ 
#+This-breaks-fleabit	         (byte-size hw:%%double-float-mantissa)
                        52.))
(defconstant high-size-const (1+
#+This-breaks-fleabit   	 (byte-size hw:%%double-float-mantissa-word2)
                        21.))

(defun find-mantissa-double-and-shift-with-sign (word1 word2 shift sign)
  (let* ((total-size (+ shift mant-size-const))
	 (words      (ceiling (1+ total-size) 32.))               ;;one for sign
	 (bignum     (allocate-bignum words))                     ;;This code assumes bignum is initialized to all zeros.
	 (plus-sign  (plusp sign)))
    (setq word1 (hw:ldb word1 hw:%%double-float-mantissa-word1 gr:*all-zero*)
	  word2 (hw:ldb word2 hw:%%double-float-mantissa-word2
			     (hw:dpb-unboxed 1 (byte 1 (byte-size hw:%%double-float-mantissa-word2))
					     gr:*all-zero*)))
    (unless plus-sign
      (setq word2 (if (hw:32zerop word1)                   ;;carry bit?
			     (hw:32logxor word2 gr:*all-ones*)
			   (hw:32+ (hw:32logxor word2 gr:*all-ones*) (hw:unboxed-constant 1)))
	    word1 (hw:32+ (hw:32logxor word1 gr:*all-ones*) (hw:unboxed-constant 1))))
    (cond (( words (ceiling total-size 32.))
	   (array::%vm-write32 bignum (1- words)
			       (hw:dpb-unboxed word2
					       (byte high-size-const (- 32. high-size-const))
					       (hw:ldb word1
						       (byte (- 32. high-size-const) (- 32. high-size-const))
						       gr:*all-zero*)))
	   (array::%vm-write32 bignum (- words 2)
			       (hw:dpb-unboxed word1
					       (byte high-size-const (- 32. high-size-const))
					       gr:*all-zero*))
	   (unless plus-sign
	     (array::%vm-write32 bignum words gr:*all-ones*)))
	  (t (find-mantissa-double-and-shift-with-sign-1 word1 word2 total-size words bignum plus-sign)))))

(defun find-mantissa-double-and-shift-with-sign-1 (word1 word2 total-size words bignum plus-sign)
  (let* ((high-bits   (- total-size (* 32. (1- words))))
	 (hi-hi-bits  (- high-bits high-size-const))
	 (remain-bits (- mant-size-const high-bits))
	 (low-bits    (- remain-bits 32.)))
    (cond ((not (plusp low-bits))
	   (let ((high-word (hw:dpb-unboxed word2
					    (byte high-size-const hi-hi-bits)
					    (hw:ldb word1
						    (byte hi-hi-bits (- 32. hi-hi-bits))
						    gr:*all-zero*))))
	     (unless plus-sign
	       (setq high-word (hw:dpb-unboxed gr:*all-ones* (byte (- 32. high-bits) high-bits) high-word)))
	     (array::%vm-write32 bignum words high-word)
	     (array::%vm-write32 bignum (1- words)
				 (hw:dpb-unboxed word1
						 (byte remain-bits (- 32. remain-bits))
						 gr:*all-zero*))))
	  (t (find-mantissa-double-and-shift-with-sign-2 word1 word2 high-bits plus-sign bignum words
							 hi-hi-bits remain-bits low-bits)))))

(defun find-mantissa-double-and-shift-with-sign-2 (word1 word2 high-bits plus-sign bignum words
						   hi-hi-bits remain-bits low-bits)
  (let* ((hi-lo-bits (- hi-hi-bits))
	 (high-word  (hw:ldb word2 (byte high-bits hi-lo-bits) gr:*all-zero*)))
    (unless plus-sign
      (setq high-word (hw:dpb-unboxed gr:*all-ones* (byte (- 32. high-bits) high-bits) high-word)))
    (array::%vm-write32 bignum words high-word)
    (array::%vm-write32 bignum (1- words)
			(hw:dpb-unboxed word2
					(byte hi-lo-bits (- 32. hi-lo-bits))
					(hw:ldb word1
						(byte (- 32. hi-lo-bits) hi-lo-bits)
						gr:*all-zero*)))
    (array::%vm-write32 bignum (- words 2)
			(hw:dpb-unboxed word1 (byte low-bits (- 32. low-bits))
					gr:*all-zero*))))

(defun find-mantissa-double-and-neg-shift-with-sign (double-word1 double-word2 neg-shift sign)
  (let* ((mant-size  (1+ (byte-size hw:%%double-float-mantissa))) ;;one for hidden bit
	 (ans-size   (+ mant-size neg-shift)))
    (cond
      ((< ans-size (byte-size vinc:%%fixnum-field))
       (* sign
	  (if (= ans-size (1- (byte-size vinc:%%fixnum-field)))
	      (hw:dpb double-word2 (byte (byte-size hw:%%double-float-mantissa-word2) 1)
		      (hw:dpb 1 (byte 1 (1+ (byte-size hw:%%double-float-mantissa-word2)))
			      (hw:ldb-boxed double-word1 (byte 1 (1- (byte-size hw:%%double-float-mantissa-word1)))
					    0)))
	    (hw:dpb double-word2 (byte (1- ans-size) 0)
		    (hw:dpb 1 (byte 1 (1- ans-size)) 0)))))
      (t (let* ((high-size  (1+ (byte-size hw:%%double-float-mantissa-word2)))
		(total-size (+ neg-shift mant-size))
		(words      (ceiling (1+ total-size) 32.))         ;;one for sign
		(bignum     (allocate-bignum words))               ;;This code assumes bignum is initialized to all zeros.
		(mant-low   (hw:ldb double-word1 hw:%%double-float-mantissa-word1 gr:*all-zero*))
		(mant-high  (hw:ldb double-word2 hw:%%double-float-mantissa-word2
				    (hw:dpb-unboxed 1 (byte 1 (byte-size hw:%%double-float-mantissa-word2))
						    gr:*all-zero*)))
		(plus-sign  (plusp sign))
		(sign-mant-low (if plus-sign
				   mant-low
				 (hw:32+ (hw:32logxor mant-low gr:*all-ones*) (hw:unboxed-constant 1))))
		(sign-mant-high (if plus-sign
				    mant-high
				  (if (hw:32zerop mant-low)                   ;;carry bit?
				      (hw:32logxor mant-high gr:*all-ones*)
				    (hw:32+ (hw:32logxor mant-high gr:*all-ones*) (hw:unboxed-constant 1))))))
	   (find-mant-dbl-and-neg-hard-case high-size total-size words bignum plus-sign sign-mant-low sign-mant-high))))))

(defun find-mant-dbl-and-neg-hard-case (high-size total-size words bignum plus-sign sign-mant-low sign-mant-high)
  (cond (( words (ceiling total-size 32.))
	 (array::%vm-write32 bignum (1- words)
			     (hw:dpb-unboxed sign-mant-high
					     (byte high-size (- 32. high-size))
					     (hw:ldb sign-mant-low
						     (byte (- 32. high-size) (- 32. high-size))
						     gr:*all-zero*)))
	 (unless plus-sign
	   (array::%vm-write32 bignum words gr:*all-ones*)))
	(t (let* ((high-bits   (- total-size (* 32. (1- words))))
		  (hi-hi-bits  (- high-bits high-size)))
	     (cond ((= total-size high-bits)
		    (let ((hi-hi-bits (max 0 hi-hi-bits))
			  (high-word (hw:dpb-unboxed sign-mant-high
						     (byte high-size hi-hi-bits)
						     (hw:ldb sign-mant-low
							     (byte hi-hi-bits (- 32. hi-hi-bits))
							     gr:*all-zero*))))
		      (unless plus-sign
			(setq high-word (hw:dpb-unboxed gr:*all-ones* (byte (- 32. high-bits) high-bits) high-word)))
		      (array::%vm-write32 bignum words high-word)))
		   (t (let* ((hi-lo-bits (- hi-hi-bits))
			     (high-word  (hw:ldb sign-mant-high (byte high-bits hi-lo-bits) gr:*all-zero*)))
			(unless plus-sign
			  (setq high-word (hw:dpb-unboxed gr:*all-ones* (byte (- 32. high-bits) high-bits) high-word)))
			(array::%vm-write32 bignum words high-word)
			(array::%vm-write32 bignum (1- words)
					    (hw:dpb-unboxed sign-mant-high
							    (byte hi-lo-bits (- 32. hi-lo-bits))
							    (hw:ldb sign-mant-low
								    (byte (- 32. hi-lo-bits) hi-lo-bits)
								    gr:*all-zero*))))))))))
(defun FIND-MANTISSA (float)
  (prims:dispatch vinc::%%data-type float
     ($$dtp-short-float     (find-mantissa-short  float))
     ($$dtp-single-float    (find-mantissa-single (array::%vm-read32 float 1)))
     ($$dtp-double-float    (find-mantissa-double (array::%vm-read32 float 1) (array::%vm-read32 float 2)))
     (t (li:tail-error "~d is not a floating point number" float))))

;**************************************************
;*              FIND-EXPONENT                     *
;**************************************************
;;; FIND-EXPONENT dispatches to these depending on the type of float suppled.
;; Returns one value:
;;    - an integer representing the exponent.
;;
;; These functions should later be made into macros for speed.

(defun find-exponent-short (short)
  (- (hw:ldb short        hw:%%short-float-exponent        0) hw:$$short-float-exponent-excess))

(defun find-exponent-single (single-word)
  (- (hw:ldb single-word  hw:%%single-float-exponent       0) hw:$$single-float-exponent-excess))

(defun find-exponent-double (double-word2)
  (- (hw:ldb double-word2 hw:%%double-float-exponent-word2 0) hw:$$double-float-exponent-excess))

(defun FIND-EXPONENT (float)
  (prims:dispatch vinc::%%data-type float
     ($$dtp-short-float     (find-exponent-short  float))
     ($$dtp-single-float    (find-exponent-single (array::%vm-read32 float 1)))
     ($$dtp-double-float    (find-exponent-double (array::%vm-read32 float 2)))
     (t (li:tail-error "~d is not a floating point number" float))))

;**************************************************
;*              SCALE-EXPONENT                     *
;**************************************************
;;; SCALE-EXPONENT dispatches to these depending on the type of float suppled.
;; Returns one value:
;;    - an integer representing the exponent.
;;
;; These functions should later be made into macros for speed.

(defun scale-exponent-short (short)
  (- (find-exponent-short  short)        (byte-size hw:%%short-float-mantissa)))

(defun scale-exponent-single (single-word)
  (- (find-exponent-single single-word)  (byte-size hw:%%single-float-mantissa)))

(defun scale-exponent-double (double-word2)
  (- (find-exponent-double double-word2) (byte-size hw:%%double-float-mantissa)))

(defun SCALE-EXPONENT (float)
  (prims:dispatch vinc::%%data-type float
     ($$dtp-short-float     (scale-exponent-short  float))
     ($$dtp-single-float    (scale-exponent-single (array::%vm-read32 float 1)))
     ($$dtp-double-float    (scale-exponent-double (array::%vm-read32 float 2)))
     (t (li:tail-error "~d is not a floating point number" float))))

;**************************************************
;*              DECODE-FLOAT                      *
;**************************************************
;;; DECODE-FLOAT dispatches to these depending on the type of float.
;;; We should return three values:
;;;    - a floating point number of the same format as our argument but
;;;      between zero inclusive and one exclusive [0,1) (a negative one exponent).
;;;    - an integer representing the exponent.
;;;    - a floating point 1.0 of the same format with the same sign as the argument.

;;; short floats are the top 26 bits of the single float representation
(defun decode-float-short (float)
  (let ((sign-value (sign-value-short float)))
    (if (zero-floatp-short float)
	(values 0.0s0 0 sign-value)
      (values (scale-mantissa-short float)
	      (find-exponent-short float)
	      sign-value))))

(defun decode-float-single (float)
  (let* ((word       (array::%vm-read32 float 1))
	 (sign-value (sign-value-single word)))
    (if (zero-floatp-single word)
	(values 0.0 0 sign-value)
      (values (scale-mantissa-single word)
	      (find-exponent-single word)
	      sign-value))))

(defun decode-float-double (float)
  (let* ((word1      (array::%vm-read32 float 1))
	 (word2      (array::%vm-read32 float 2))
	 (sign-value (sign-value-double word2)))
    (if (zero-floatp-double word1 word2)
	(values 0.0d0 0 sign-value)
      (values (scale-mantissa-double word1 word2)
	      (find-exponent-double word2)
	      sign-value))))

(defun DECODE-FLOAT (float)
  (prims:dispatch vinc::%%data-type float
     ($$dtp-short-float     (decode-float-short float))
     ($$dtp-single-float    (decode-float-single float))
     ($$dtp-double-float    (decode-float-double float))
     (t (li:tail-error "~d is not a floating point number" float))))

;**************************************************
;*  INTEGER-DECODE-FLOAT
;**************************************************
;;; These are similar to the DECODE-FLOAT functions except
;;; that the first value is an integer which represents the
;;; scaled mantissa.
(defun integer-decode-float-short (float)
  (let ((sign-value (sign-value-short float)))
    (if (zero-floatp-short float)
	(values 0 0 sign-value)
      (values (find-mantissa-short float)
	      (scale-exponent-short float)
	      sign-value))))

(defun integer-decode-float-single (float)
  (let* ((word       (array::%vm-read32 float 1))
	 (sign-value (sign-value-single word)))
    (if (zero-floatp-single word)
	(values 0 0 sign-value)
	(values	(find-mantissa-single word)
  		(scale-exponent-single word)
		sign-value))))

(defun integer-decode-float-double (float)
  (let* ((word1      (array::%vm-read32 float 1))
	 (word2      (array::%vm-read32 float 2))
	 (sign-value (sign-value-double word2)))
    (if (zero-floatp-double word1 word2)
	(values 0 0 sign-value)
	(values (find-mamtissa-double word1 word2)
	        (scale-exponent-double word2)
		sign-value))))

(defun INTEGER-DECODE-FLOAT (float)
  (prims:dispatch vinc::%%data-type float
     ($$dtp-short-float     (integer-decode-float-short float))
     ($$dtp-single-float    (integer-decode-float-single float))
     ($$dtp-double-float    (integer-decode-float-double float))
     (t (li:tail-error "~d is not a floating point number" float))))


;**************************************************
;*  things to TRUNCATE floats                     *
;**************************************************

(defun significant-bits (n)
  (if (zerop n) 0 (- 32. (hw:32-prioritize n)))) ;;Is this zero test neccessary?

(defun truncate-short (float)
  (let* ((mantissa (find-mantissa-short  float))
	 (exponent (scale-exponent-short float))
	 (sign-bit (hw:ldb-boxed float hw:%%short-float-sign 0))
	 (sign     (if (zerop sign-bit) 1 -1)))
    (cond ((not (minusp exponent))              		  ;all bits are integer, none are fraction:
	   (values (lisp:ash (* sign mantissa) exponent) 0 sign))
	  ((> exponent (- (1+ (byte-size hw:%%short-float-mantissa))))
	   ;;; remember that with the hidden bit, the mantissa is one bit wider than the field it is stored in.
	   ;;; Some bits are integer, some are fraction.  Exponent is enough to shift the binary point
	   ;;; somewhere into the mantissa but not so big as to move it beyond the left edge.
	   (let* ((remain-mant	          (hw:ldb mantissa (byte (- exponent) 0) gr:*all-zero*))
		  (mantissa-byte-width	  (1- (significant-bits remain-mant))))
	     (if (minusp mantissa-byte-width)
		 (values (* sign (hw:ldb mantissa (byte (+ (1+ (byte-size hw:%%short-float-mantissa)) exponent)
							(- exponent))
					 0))
			 (hw:dpb sign-bit hw:%%short-float-sign 0.0s0))
	       (truncate-short-1 remain-mant mantissa mantissa-byte-width
				 exponent sign sign-bit))))
	  (t ;;;all the bits are fraction so the argument is the remainder
	   (values 0 float)))))

(defun truncate-short-1 (remain-mant mantissa mantissa-byte-width
			 exponent sign sign-bit)
  (let* ((mantissa-byte-position (- (byte-size hw:%%short-float-mantissa) mantissa-byte-width))
	 (new-expt               (+ exponent mantissa-byte-width))
	 (new-remain-mant        (hw:dpb remain-mant
					 (byte mantissa-byte-width mantissa-byte-position)
					 0.0s0))
	 (quotient-mant          (hw:ldb mantissa (byte (+ (1+ (byte-size hw:%%short-float-mantissa)) exponent)
							(- exponent))
					 0)))
    (values (* sign quotient-mant)
	    (hw:dpb (+ new-expt hw:$$short-float-exponent-excess)
		    hw:%%short-float-exponent
		    (hw:dpb sign-bit hw:%%short-float-sign new-remain-mant)))))

(defun truncate-single (float)
  (let* ((word     (array::%vm-read32 float 1))
	 (exponent (scale-exponent-single word))
	 (sign-bit (hw:ldb-boxed word hw:%%single-float-sign 0))
	 (sign     (if (zerop sign-bit) 1 -1)))
    (cond ((not (minusp exponent))
	   ;;; all bits are integer, none are fraction:
	   (values (find-mantissa-single-and-shift-with-sign word exponent sign) 0))
	  ((> exponent (- (1+ (byte-size hw:%%single-float-mantissa))))
	   ;;; remember that with the hidden bit, the mantissa is one bit wider than the field it is stored in.
	   ;;; Some bits are integer, some are fraction.  Exponent is enough to shift the binary point
	   ;;; somewhere into the mantissa but not so big as to move it beyond the left edge.
	   (let* ((remain-mant            (hw:ldb word (byte (- exponent) 0) gr:*all-zero*))
		  (mantissa-byte-width	  (1- (significant-bits remain-mant)))
		  (quotient               (* sign (hw:ldb-boxed word
								(byte (+ (byte-size hw:%%single-float-mantissa) exponent)
								      (- exponent))
								(hw:dpb 1 (byte 1 (byte-size hw:%%single-float-mantissa))
									0)))))
	     (if (minusp mantissa-byte-width)
		 (values quotient (if (plusp sign) 0.0 -0.0))
	       (let* ((mantissa-byte-position (- (byte-size hw:%%single-float-mantissa)
						 mantissa-byte-width))
		      (new-expt               (+ exponent mantissa-byte-width))
		      (new-remain-mant        (hw:dpb-unboxed remain-mant
							      (byte mantissa-byte-width mantissa-byte-position)
							      gr:*all-zero*)))
		 (values quotient
			 (array:make-single-float
			   (hw:dpb-unboxed (+ new-expt hw:$$single-float-exponent-excess)
					   hw:%%single-float-exponent
					   (hw:dpb-unboxed sign-bit hw:%%single-float-sign new-remain-mant))))))))
	  (t ;;;all the bits are fraction so the argument is the remainder
	   (values 0 float)))))

(defun truncate-double (float)
  (let* ((word1    (array::%vm-read32 float 1))
	 (word2    (array::%vm-read32 float 2))
	 (exponent (scale-exponent-double word2))
	 (sign-bit (hw:ldb-boxed word2 hw:%%double-float-sign-word2 0))
	 (sign     (if (zerop sign-bit) 1 -1)))
    (cond ((not (minusp exponent))   	       ;;; all bits are integer, none are fraction:
	   (values (find-mantissa-double-and-shift-with-sign word1 word2 exponent sign) 0))
	  ((plusp (+ (1+ (byte-size hw:%%double-float-mantissa)) exponent))
	   ;;; remember that with the hidden bit, the mantissa is one bit wider than the field it is stored in.
	   ;;; Some bits are integer, some are fraction.  Exponent is enough to shift the binary point
	   ;;; somewhere into the mantissa but not so big as to move it beyond the left edge.
	   (if ( exponent -32.)
	       (values (find-mantissa-double-and-neg-shift-with-sign word1 word2 exponent sign)
		       (trunc-dbl-remain1 (hw:ldb word1 (byte (- exponent) 0) gr:*all-zero*)
					 exponent sign-bit))
	     (let* ((remain-word1 word1)
		    (remain-word2 (hw:ldb word2 (byte (- -32. exponent) 0) gr:*all-zero*))
		    (remain-width (1- (significant-bits remain-word2))))
	       (values (find-mantissa-double-and-neg-shift-with-sign word1 word2 exponent sign)
		       (if (minusp remain-width)
			   (trunc-dbl-remain1 remain-word1 exponent sign-bit)
			 (trunc-dbl-remain2 remain-word1 remain-word2 (+ remain-width 32.) exponent sign-bit))))))
	  (t ;;;all the bits are fraction so the argument is the remainder
	   (values 0 float)))))

(defun trunc-dbl-remain2 (remain-word1 remain-word2 width exponent sign-bit)
  (let* ((position (- (byte-size hw:%%double-float-mantissa) width))
	 (new-expt (+ exponent width))
	 (high     (hw:dpb-unboxed remain-word2
				   (byte (- (byte-size hw:%%double-float-mantissa-word2) position)  position)
				   (hw:ldb remain-word1
					   (byte position (- (byte-size hw:%%double-float-mantissa-word1) position))
					   gr:*all-zero*)))
	 (low      (hw:dpb-unboxed remain-word1
				   (byte (- (byte-size hw:%%double-float-mantissa-word1) position) position)
				   gr:*all-zero*)))
    (array:make-double-float
      (hw:dpb-unboxed (+ new-expt hw:$$double-float-exponent-excess)
		      hw:%%double-float-exponent-word2
		      (hw:dpb-unboxed sign-bit hw:%%double-float-sign-word2 high))
      low)))

(defun trunc-dbl-remain1 (remain-word1 exponent sign-bit)
  (let* ((width    (1- (significant-bits remain-word1)))
	 (one-word ( width (byte-size hw:%%double-float-mantissa-word2))))
    (if (minusp width)
	(if (zerop sign-bit) 0.0d0 -0.0d0)         ;;+++ The reader may not work on 0.0 --wkf
      (let ((position (- (byte-size hw:%%double-float-mantissa) width))
	    (new-expt (+ exponent width))
	    (high     (if one-word
			  (hw:dpb-unboxed remain-word1
					  (byte width
						(- (byte-size hw:%%double-float-mantissa-word2)
						   width))
					  gr:*all-zero*)
			(hw:ldb remain-word1 (byte (byte-size hw:%%double-float-mantissa-word2)
						   (- width
						      (byte-size hw:%%double-float-mantissa-word2)))
				gr:*all-zero*)))
	    (low      (if one-word
			  gr:*all-zero*
			(hw:dpb remain-word1
				(byte (- width
					 (byte-size hw:%%double-float-mantissa-word2))
				      (- 32. (- width
						(byte-size hw:%%double-float-mantissa-word2))))
				gr:*all-zero*))))
	(array:make-double-float
	  (hw:dpb-unboxed (+ new-expt hw:$$double-float-exponent-excess)
			  hw:%%double-float-exponent-word2
			  (hw:dpb-unboxed sign-bit hw:%%double-float-sign-word2 high))
	  low)))))