(*-------------------------------------------------------------------------*)
(*          Ninv -- Find Percentage Point of Normal Distribution           *)
(*-------------------------------------------------------------------------*)

FUNCTION Ninv( P : REAL ) : REAL;

(*-------------------------------------------------------------------------*)
(*                                                                         *)
(*       Function:  Ninv                                                   *)
(*                                                                         *)
(*       Purpose:   Finds percentage point of normal distribution          *)
(*                  (low accuracy)                                         *)
(*                                                                         *)
(*       Calling Sequence:                                                 *)
(*                                                                         *)
(*            Perc  := Ninv( P : REAL ) : REAL;                            *)
(*                                                                         *)
(*                 P      --- input probability                            *)
(*                 Perc   --- resultant percentage point                   *)
(*                                                                         *)
(*       Calls:  None                                                      *)
(*                                                                         *)
(*       Remark:  This method provides about 6.5 decimal digits of         *)
(*                accuracy.                                                *)
(*                                                                         *)
(*-------------------------------------------------------------------------*)

CONST
   Lim = 1.0E-20;

VAR
   Y:  REAL;
   Pr: REAL;
   Nv: REAL;

(* Structured *) CONST

   PN:  ARRAY[1..5] OF REAL =
        ( -0.322232431088  ,
          -1.0             ,
          -0.342242088547  ,
          -0.0204231210245 ,
          -0.453642210148E-4 );

   QN:  ARRAY[1..5] OF REAL =
        ( 0.0993484626060 ,
          0.588581570495  ,
          0.531103462366  ,
          0.103537752850  ,
          0.38560700634E-2  );

BEGIN (* Ninv *)

   Ninv   := 0.0;

   IF( P > 0.5 ) THEN
      Pr := 1.0 - P
   ELSE
      Pr := P;

   IF( ( Pr >= Lim ) AND ( Pr <> 0.5 ) ) THEN

      BEGIN

         Y      := SQRT ( LN( 1.0 / Pr / Pr ) );

         Nv     := Y + ((((Y * PN[ 5 ] + PN[ 4 ]) * Y + PN[ 3 ] ) * Y
                                       + PN[ 2 ]) * Y + PN[ 1 ] ) /
                       ((((Y * QN[ 5 ] + QN[ 4 ]) * Y + QN[ 3 ] ) * Y
                                       + QN[ 2 ]) * Y + QN[ 1 ] );

         IF( P < 0.5 ) THEN
            Ninv := -Nv
         ELSE
            Ninv := Nv;

      END;

END   (* Ninv *);