MODULE Maths;
(* =========================================================================  
   Additional Maths Functions 
 
   Author : Chris Burrows  
   Created: May 2007   
   (c) 2007-2008 CFB Software   
   http://www.cfbsoftware.com/gpcp  
   
   ========================================================================= *)

IMPORT
  Sys := "[mscorlib]System";

TYPE
  Random* = POINTER TO RECORD
    z: INTEGER 
  END;
  
(* ==================================================================== *)

PROCEDURE Sum*(IN r: ARRAY OF REAL): REAL;
VAR
  i: INTEGER;
  sum: REAL;
BEGIN
  sum := 0.0;
  FOR i := 0 TO LEN(r) - 1 DO sum := sum + r[i] END;
  RETURN sum
END Sum;

(* ==================================================================== *)

PROCEDURE SumOfSquares*(IN r: ARRAY OF REAL): REAL;
VAR
  i: INTEGER;
  sum: REAL;
BEGIN
  sum := 0.0;
  FOR i := 0 TO LEN(r) - 1 DO sum := sum + (r[i] * r[i]) END;
  RETURN sum
END SumOfSquares;

(* ==================================================================== *)

PROCEDURE Mean*(IN r: ARRAY OF REAL): REAL;
VAR
  i: INTEGER;
BEGIN
  RETURN Sum(r) / LEN(r)
END Mean;

(* ==================================================================== *)

PROCEDURE Variance*(IN r: ARRAY OF REAL): REAL;
VAR
  mean, sum, diff: REAL;
  i: INTEGER;
BEGIN
  sum := 0.0;
  mean := Mean(r);
  FOR i := 0 TO LEN(r) - 1 DO 
    diff := r[i] - mean;
    sum := sum + (diff * diff)
  END;
  RETURN sum / (LEN(r) - 1)
END Variance;

(* ==================================================================== *)

PROCEDURE StdDev*(IN r: ARRAY OF REAL): REAL;
BEGIN
  RETURN Sys.Math.Sqrt(Variance(r))
END StdDev;

(* ==================================================================== *)

PROCEDURE (rand: Random) Next*(range: INTEGER): INTEGER, NEW;
(* REF: Programming in Oberon - Reiser & Wirth, ACM Press 1992 *)
CONST 
  a = 16807; m = 2147483647; q = m DIV a; r = m MOD a;
VAR 
  gamma: INTEGER;
BEGIN
  gamma := a * (rand.z MOD q) - r * (rand.z DIV q);
  IF gamma > 0 THEN rand.z := gamma ELSE rand.z := gamma + m END;
  RETURN SHORT((LONG(rand.z) * range) DIV m) 
END Next;

(* ==================================================================== *)

PROCEDURE (rand: Random) InitSeed* (seed: INTEGER), NEW;
BEGIN
  rand.z := seed
END InitSeed;

(* ==================================================================== *)

PROCEDURE InitRandom*(): Random;
VAR
  rand: Random;
BEGIN
  NEW(rand);
  rand.InitSeed(314159);
  RETURN rand
END InitRandom;

END Maths.