X-Git-Url: http://git.kpe.io/?p=kmrcl.git;a=blobdiff_plain;f=math.lisp;h=fcbab3c4baa2e922b0eba049af0bd04079a19739;hp=682cd74027776f7dfc5a3d4346c7a0ab2be55f50;hb=5470e786c59285e2f09497598505654672be847d;hpb=4a5b626f01db51b02f969adb33ddad6aa9ee303a diff --git a/math.lisp b/math.lisp index 682cd74..fcbab3c 100644 --- a/math.lisp +++ b/math.lisp @@ -7,7 +7,7 @@ ;;;; Programmer: Kevin M. Rosenberg ;;;; Date Started: Nov 2002 ;;;; -;;;; $Id: math.lisp,v 1.4 2003/06/06 21:59:29 kevin Exp $ +;;;; $Id$ ;;;; ;;;; This file, part of KMRCL, is Copyright (c) 2002 by Kevin M. Rosenberg ;;;; @@ -34,3 +34,65 @@ `(if (stringp ,obj) (parse-integer ,obj) ,obj)) + +(defun histogram (v n-bins &key min max) + (declare (fixnum n-bins)) + (when (listp v) + (setq v (coerce v 'vector))) + (when (zerop (length v)) + (return-from histogram (values nil nil nil)) ) + (let ((n (length v)) + (bins (make-array n-bins :element-type 'integer :initial-element 0)) + found-min found-max) + (declare (fixnum n)) + (unless (and min max) + (setq found-min (aref v 0) + found-max (aref v 0)) + (loop for i fixnum from 1 to (1- n) + do + (let ((x (aref v i))) + (cond + ((> x found-max) + (setq found-max x)) + ((< x found-min) + (setq found-min x))))) + (unless min + (setq min found-min)) + (unless max + (setq max found-max))) + (let ((width (/ (- max min) n-bins))) + (setq width (+ width (* double-float-epsilon width))) + (dotimes (i n) + (let ((bin (nth-value 0 (truncate (- (aref v i) min) width)))) + (declare (fixnum bin)) + (when (and (not (minusp bin)) + (< bin n-bins)) + (incf (aref bins bin)))))) + (values bins min max))) + + +(defun fixnum-width () + (nth-value 0 (truncate (+ (/ (log (1+ most-positive-fixnum)) (log 2)) .5)))) + +(defun scaled-epsilon (float &optional (operation '+)) + "Return the smallest number that would return a value different from + FLOAT if OPERATION were applied to FLOAT and this number. OPERATION + should be either + or -, and defauls to +." + (multiple-value-bind (significand exponent) + (decode-float float) + (multiple-value-bind (1.0-significand 1.0-exponent) + (decode-float (float 1.0 float)) + (if (and (eq operation '-) + (= significand 1.0-significand)) + (scale-float (typecase float + (short-float short-float-negative-epsilon) + (single-float single-float-negative-epsilon) + (double-float double-float-negative-epsilon) + (long-float long-float-negative-epsilon)) + (- exponent 1.0-exponent)) + (scale-float (typecase float + (short-float short-float-epsilon) + (single-float single-float-epsilon) + (double-float double-float-epsilon) + (long-float long-float-epsilon)) + (- exponent 1.0-exponent))))))