(* Copyright (C) 1992, Xerox                                                 *)
(* All rights reserved.                                                      *)

(* Last modified on Tue Feb 11 15:18:43 PST 1992 by muller                   *)
(*      modified on Wed Sep 25 00:33:01 1991 by goldberg@xerox.parc.com      *)

UNSAFE MODULE RealFloat;

IMPORT SunOsIeee, Word;

(* SunOs only has IEEE routines for the C 'double' data type.  The
   first 3 routines can be computed by converting to LONGREAL and back.
   The others require more work *)

PROCEDURE Scalb(x: REAL; n: INTEGER): REAL =
BEGIN
	RETURN(FLOAT(SunOsIeee.scalbn(FLOAT(x, LONGREAL), n)));
END Scalb;

PROCEDURE ILogb(x: REAL): INTEGER =
BEGIN
	RETURN(SunOsIeee.ilogb(FLOAT(x, LONGREAL)));
END ILogb;

PROCEDURE Sqrt(x: REAL): REAL =
BEGIN
   RETURN(FLOAT(SunOsIeee.sqrt(FLOAT(x, LONGREAL))));
END Sqrt;

PROCEDURE Differs(x, y: T): BOOLEAN =
BEGIN
  RETURN (x < y OR y < x);
END Differs;
    
PROCEDURE Unordered(x, y: T): BOOLEAN =
BEGIN
  RETURN (NOT (x <= y OR y <= x));
END Unordered;

PROCEDURE Class(x: REAL): IEEEClass =
  VAR
    w: Word.T;
    y: REAL;
    exp, frac: INTEGER;
  BEGIN
    y := x; (* make sure that x is converted from C double to C float *)
    w := LOOPHOLE(y, Word.T);
    exp := Word.Extract(w, i := 23, n := 8);
    frac := Word.Extract(w, i := 0, n := 23);
    IF exp = 0 THEN
      IF frac = 0 THEN
	RETURN (IEEEClass.Zero)
      ELSE
	RETURN (IEEEClass.Denormal);
      END;
    ELSIF exp = 255 THEN
      IF frac = 0 THEN
	RETURN (IEEEClass.Infinity)
      ELSE
	(* XXX: what about Signaling NaN? *)
	RETURN (IEEEClass.QuietNaN); END;
    ELSE
      RETURN (IEEEClass.Normal);
    END;
  END Class;

PROCEDURE Sign(x: REAL): [0..1] =
  VAR
    y: REAL;
    w: Word.T;
  BEGIN
    y := x; (* make sure that x is converted from C double to C float *)
    w := LOOPHOLE(y, Word.T);
    RETURN(Word.Extract(w, i := 31, n := 1));
  END Sign;

PROCEDURE Finite(x: T): BOOLEAN =
  VAR
    y: REAL;
    w: Word.T;
  BEGIN
    y := x; (* make sure that x is converted from C double to C float *)
    w := LOOPHOLE(y, Word.T);
    RETURN(Word.Extract(w, i := 23, n := 8) # 255);
  END Finite;

PROCEDURE IsNaN(x: T): BOOLEAN =
  VAR
    y: REAL;
    w: Word.T;
    exp, frac: INTEGER;
  BEGIN
    y := x; (* make sure that x is converted from C double to C float *)
    w := LOOPHOLE(y, Word.T);
    exp := Word.Extract(w, i := 23, n := 8);
    frac := Word.Extract(w, i := 0, n := 23);
    RETURN(exp = 255 AND frac # 0);
  END IsNaN;

PROCEDURE Logb(x: T): T =
  VAR
    y: REAL;
    w: Word.T;
    exp: INTEGER;
  BEGIN
    IF x = 0.0 THEN
      RETURN(-1.0/ABS(x));
    ELSIF IsNaN(x) OR NOT Finite(x) THEN
      RETURN(x*x);
    ELSE
      y := x; (* make sure that x is converted from C double to C float *)
      w := LOOPHOLE(y, Word.T);
      exp := Word.Extract(w, i := 23, n := 8);
      IF exp = 0 THEN 
	RETURN(-126.0);
      ELSE
	RETURN(FLOAT(exp - 127));
      END;
    END;
  END Logb;

PROCEDURE CopySign(x, y: T): T =
  VAR
    x1, y1: REAL;
    xw, yw: Word.T;
  BEGIN
    x1 := x; (* make sure that x is converted from C double to C float *)
    xw := LOOPHOLE(x1, Word.T);

    y1 := y; (* make sure that y is converted from C double to C float *)
    yw := LOOPHOLE(y1, Word.T);

    xw := Word.Insert(xw, Word.Extract(yw, i := 31, n := 1), i := 31, n := 1);
    RETURN(LOOPHOLE(xw, REAL));
  END CopySign;

PROCEDURE NextAfter(x, y: T): T =
  VAR
    x1: REAL;
    xw: Word.T;
    exp: INTEGER;
  BEGIN
    IF x = y THEN RETURN(x); END;
    IF IsNaN(x) THEN
      RETURN x
    ELSIF IsNaN(y)
      THEN RETURN(y);
    END;

    IF x = 0.0 THEN 
      xw := 0;
      xw := Word.Insert(xw, 1, i := 0, n := 1);
      xw := Word.Insert(xw, Sign(y), i := 31, n := 1);
      RETURN(LOOPHOLE(xw, REAL));
    END;

    x1 := x; (* make sure that x is converted from C double to C float *)
    xw := LOOPHOLE(x1, Word.T);
    IF (x > 0.0 AND x > y) OR (x < 0.0 AND x < y) THEN 
	DEC(xw);
    ELSE
	INC(xw);
    END;

    exp := Word.Extract(xw, i := 23, n := 8);
    x1 := LOOPHOLE(xw, REAL);
    IF exp = 255 THEN
      RETURN(x + x);  (* generate overflow *)
    ELSIF exp = 0 THEN
      RETURN((2.0*x1)/2.0);  (* generate underflow *)
    ELSE
      RETURN(x1);
    END;
  END NextAfter;

BEGIN
(* start code *)
END RealFloat.