PROGRAM d13r2(input,output);
(* driver for routine MDIAN2 *)
CONST
   npts=50;
TYPE
   glsarray = ARRAY [1..npts] OF real;
   narray = glsarray;
VAR
   glinext,glinextp : integer;
   glma : ARRAY [1..55] OF real;
   gliset : integer;
   glgset : real;
   i,idum : integer;
   xmed : real;
   data : narray;

(*$I MODFILE.PAS *)
(*$I RAN3.PAS *)

(*$I GASDEV.PAS *)

(*$I SORT.PAS *)

(*$I MDIAN1.PAS *)

(*$I MDIAN2.PAS *)

BEGIN
   gliset := 0;
   idum := -5;
   FOR i := 1 to npts DO BEGIN
      data[i] := gasdev(idum)
   END;
   mdian2(data,npts,xmed);
   writeln('Data drawn from a gaussian distribution');
   writeln('with zero mean, unit variance');
   writeln;
   writeln('median according to mdian2 is',xmed:9:6);
   mdian1(data,npts,xmed);
   writeln('median according to mdian1 is',xmed:9:6)
END.