FTN,L,B 
C 
C HP92403A STATISTICS PACKAGE 
C 
C SOURCE TAPE 92403-80001 REV. A
C RELOC. TAPE 92403-60001 REV. A
C 
C AUTHOR - T.A. SAPONAS 
C 
C VERSION OF DECEMBER 1973
C 
C 
      SUBROUTINE HISTI(IDATA,NPTS,ISTRT,IDLTA,NHIST,I,IERR) 
C  THIS SUBROUTINE PRODUCES A HISTOGRAM IN ARRAY "NHIST" FROM THE 
C   DATA IN ARRAY "IDATA" 
      DIMENSION IDATA(1),NHIST(1) 
      NBARS = IABS(I) 
C 
C  INITIALIZE ALL BARS OF THE HISTOGRAM TO ZERO 
      DO 1 J = 1,NBARS
    1 NHIST(J) = 0
      IERR = 0
C 
C   LOOP THROUGH ALL OF THE DATA
      DO 10 J = 1,NPTS
C 
C  COMPUTE BAR NUMBER  (M = BAR NUMBER - 1) 
      M = (IDATA(J)-ISTRT)/IDLTA
C 
C CHECK TO SEE DATA IS IN BOUNDS OF HISTOGRAM 
      IF(M)400,301
  301 IF(M-NBARS)302,450
C 
C  DATA GREATER THAN UPPER BOUND OF HISTOGRAM, INCREMENT IERR 
  450 IERR = IERR-1 
C  IF I>0 INCLUDE POINT IN LAST BAR COUNT, ELSE IGNORE IT 
      M = NBARS-1 
      IF(I)10,402 
C 
C  DATA LESS THAN LOWER BOUND OF HISTOGRAM, INCREMENT IERR
  400 IERR = IERR-1 
C  IF I>0 INCLUDE POINT IN FIRST BAR COUNT, ELSE IGNORE IT
      IF(I)10,401 
  401 M = 0 
  402 IERR = IERR+2 
C 
C  INCREMENT PROPER BAR OF HISTOGRAM
  302 NHIST(M+1) = NHIST(M+1)+1 
   10 CONTINUE
      RETURN
      END 
      SUBROUTINE HISTF(DATA,NPTS,START,DELTA,NHIST,I,IERR)
C  THIS SUBROUTINE PRODUCES A HISTOGRAM IN ARRAY "NHIST" FROM THE 
C   DATA IN ARRAY "DATA"
      DIMENSION DATA(1),NHIST(1)
      NBARS = IABS(I) 
C 
C  INITIALIZE ALL BARS OF THE HISTOGRAM TO ZERO 
      DO 1 J = 1,NBARS
    1 NHIST(J) = 0
      IERR = 0
C 
C   LOOP THROUGH ALL OF THE DATA
      DO 10 J = 1,NPTS
C 
C  COMPUTE BAR NUMBER  (M = BAR NUMBER - 1) 
      M = (DATA(J)-START)/DELTA 
C 
C CHECK TO SEE DATA IS IN BOUNDS OF HISTOGRAM 
      IF(M)400,301
  301 IF(M-NBARS)302,450
C 
C  DATA GREATER THAN UPPER BOUND OF HISTOGRAM, INCREMENT IERR 
  450 IERR = IERR-1 
C  IF I>0 INCLUDE POINT IN LAST BAR COUNT, ELSE IGNORE IT 
      M = NBARS-1 
      IF(I)10,402 
C 
C  DATA LESS THAN LOWER BOUND OF HISTOGRAM, INCREMENT IERR
  400 IERR = IERR-1 
C  IF I>0 INCLUDE POINT IN FIRST BAR COUNT, ELSE IGNORE IT
      IF(I)10,401 
  401 M = 0 
  402 IERR = IERR +2
C 
C  INCREMENT PROPER BAR OF HISTOGRAM
  302 NHIST(M+1) = NHIST(M+1)+1 
   10 CONTINUE
      RETURN
      END 
      SUBROUTINE HISTB(DATA,NPTS,START,DELTA,RHIST,I,IERR)
C  THIS SUBROUTINE PRODUCES A HISTOGRAM IN ARRAY "RHIST" FROM THE 
C   DATA IN ARRAY "DATA".  RHIST IS A FLOATING POINT ARRAY SO THAT
C   HISTB IS "BASIC CALLABLE".
      DIMENSION DATA(1),RHIST(1)
      NBARS = IABS(I) 
C 
C  INITIALIZE ALL BARS OF THE HISTOGRAM TO ZERO 
      DO 1 J = 1,NBARS
    1 RHIST(J) = 0. 
      IERR = 0
C 
C   LOOP THROUGH ALL OF THE DATA
      DO 10 J = 1,NPTS
C 
C  COMPUTE BAR NUMBER  (M = BAR NUMBER - 1) 
      M = (DATA(J)-START)/DELTA 
C 
C CHECK TO SEE DATA IS IN BOUNDS OF HISTOGRAM 
      IF(M)400,301
  301 IF(M-NBARS)302,450
C 
C  DATA GREATER THAN UPPER BOUND OF HISTOGRAM, INCREMENT IERR 
  450 IERR = IERR-1 
C  IF I>0 INCLUDE POINT IN LAST BAR COUNT, ELSE IGNORE IT 
      M = NBARS-1 
      IF(I)10,402 
C 
C  DATA LESS THAN LOWER BOUND OF HISTOGRAM, INCREMENT IERR
  400 IERR = IERR-1 
C  IF I>0 INCLUDE POINT IN FIRST BAR COUNT, ELSE IGNORE IT
      IF(I)10,401 
  401 M = 0 
  402 IERR = IERR +2
C 
C  INCREMENT PROPER BAR OF HISTOGRAM
  302 RHIST(M+1) = RHIST(M+1)+1.
   10 CONTINUE
      RETURN
      END 
      SUBROUTINE STATI(IDATA,NPTS,RMEAN,STDEV)
C "STATI" COMPUTES THE MEAN AND STANDARD DEVIATION FROM THE MEAN
C  OF THE DATA IN ARRAY "IDATA".
      DIMENSION IDATA(1)
C 
C  CLEAR THE ACCUMULATORS FOR THE SUM AND SUM OF SQUARES
      SUM = 0.
      SUMSQ = 0.
C 
C COMPUTE SUM AND SUM OF SQUARES OF THE DATA
      DO 10 J = 1,NPTS
      DATA = IDATA(J) 
      SUM = SUM+DATA
   10 SUMSQ = SUMSQ+DATA*DATA 
C 
C  COMPUTE MEAN AND STANDARD DEVIATION
      RNPTS = NPTS
      RMEAN = SUM/RNPTS 
      STDEV = SQRT((SUMSQ-SUM*RMEAN)/(RNPTS-1.))
      END 
      SUBROUTINE STATF(DATA,NPTS,RMEAN,STDEV) 
C "STATF" COMPUTES THE MEAN AND STANDARD DEVIATION FROM THE MEAN
C  OF THE DATA IN ARRAY "DATA". 
      DIMENSION DATA(1) 
C 
C  CLEAR THE ACCUMULATORS FOR THE SUM AND SUM OF SQUARES
      SUM = 0.
      SUMSQ = 0.
C 
C COMPUTE SUM AND SUM OF SQUARES OF THE DATA
      DO 10 J = 1,NPTS
      SUM = SUM+DATA(J) 
   10 SUMSQ = SUMSQ+DATA(J)*DATA(J) 
C 
C  COMPUTE MEAN AND STANDARD DEVIATION
      RNPTS = NPTS
      RMEAN = SUM/RNPTS 
      STDEV = SQRT((SUMSQ-SUM*RMEAN)/(RNPTS-1.))
      END 
      END$
FTN,L,B 
      SUBROUTINE INTLI(ITYPE,A,ISTRT,IDLTA,NHIST,I) 
C   "INTLI" INITIALIZES THE HISTOGRAM AND OTHER STATISTICAL 
C     ACCUMULATORS FOR RUNNING STATISTICS ON INTEGER DATA.
C  ARRAY "A" MUST BE PROVIDED FOR EACH VARIABLE FOR WHICH STATISTICS
C  ARE MAINTAINED.  "A" MUST BE A REAL ARRAY DIMENSIONED BY 5.
C  THE CONTENTS OF "A" ARE AS FOLLOWS:
C       A(1) - RUNNING SUM OF DATA
C       A(2) - RUNNING SUM OF SQUARES OF DATA 
C       A(3) - WORD 1 - NUMBER OF POINTS
C              WORD 2 - MODE OF STATISTICS
C       A(4) - WORD 1 - LOWER BOUND OF HISTOGRAM
C              WORD 2 - WIDTH OF EACH BAR IN HISTOGRAM
C       A(5) - WORD 1 - NUMBER OF BARS IN HISTOGRAM 
C              WORD 2 - MODE OF HISTOGRAM 
      DIMENSION A(5),N(2),NHIST(1)
C 
C THE FOLLOWING EQUIVALENCE STATEMENT GIVES ACCESS TO THE TWO 
C WORDS OF THE FLOATING POINT VARIABLE "RN", WHERE IWRD1 IS THE 
C FIRST WORD AND IWRD2 IS THE SECOND WORD.
      EQUIVALENCE (RN,N,IWRD1),(N(2),IWRD2) 
C 
C INITIALIZE SUM AND SUM OF SQUARES TO 0.0
      A = 0.
      A(2) = 0. 
C 
C  INITIALIZE NUMBER OF POINTS TO 0 AND MODE OF STATISTICS TO ITYPE 
      IWRD1 = 0 
      IWRD2 = ITYPE 
      A(3) = RN 
C 
C  IF ITYPE < OR = 0 THEN INITIALIZE HISTOGRAM OTHERWISE RETURN 
      IF(ITYPE)10,10,150
   10 IWRD1 = IABS(I) 
      DO 20 J = 1,IWRD1 
   20 NHIST(J) = 0
C 
C A(5) IS SET TO NUMBER OF BARS IN HISTOGRAM AND MODE OF HISTOGRAM
      IWRD2 = I 
      A(5) = RN 
C 
C A(4) IS SET TO LOWER BOUND AND WIDTH OF HISTOGRAM 
      IWRD1 = ISTRT 
      IWRD2 = IDLTA 
      A(4) = RN 
150   RETURN
      END 
      SUBROUTINE INTLF(ITYPE,A,START,DELTA,NHIST,I) 
C   "INTLF" INITIALIZES THE HISTOGRAM AND OTHER STATISTICAL 
C     ACCUMULATORS FOR RUNNING STATISTICS ON FLOATING POINT DATA. 
C  ARRAY "A" MUST BE PROVIDED FOR EACH VARIABLE FOR WHICH STATISTICS
C  ARE MAINTAINED.  "A" MUST BE A REAL ARRAY DIMENSIONED BY 6.
C  THE CONTENTS OF "A" ARE AS FOLLOWS:
C       A(1) - RUNNING SUM OF DATA
C       A(2) - RUNNING SUM OF SQUARES OF DATA 
C       A(3) - WORD 1 - NUMBER OF POINTS
C              WORD 2 - MODE OF STATISTICS
C       A(4) - LOWER BOUND OF HISTOGRAM 
C       A(5) - WIDTH OF EACH BAR IN HISTOGRAM 
C       A(6) - WORD 1 - NUMBER OF BARS IN HISTOGRAM 
C              WORD 2 - MODE OF HISTOGRAM 
      DIMENSION A(6),N(2),NHIST(1)
C 
C THE FOLLOWING EQUIVALENCE STATEMENT GIVES ACCESS TO THE TWO 
C WORDS OF THE FLOATING POINT VARIABLE "RN", WHERE IWRD1 IS THE 
C FIRST WORD AND IWRD2 IS THE SECOND WORD.
      EQUIVALENCE (RN,N,IWRD1),(N(2),IWRD2) 
C 
C INITIALIZE SUM AND SUM OF SQUARES TO 0.0
      A = 0.
      A(2) = 0. 
C 
C  INITIALIZE NUMBER OF POINTS TO 0 AND MODE OF STATISTICS TO ITYPE 
      IWRD1 = 0 
      IWRD2 = ITYPE 
      A(3) = RN 
C 
C  IF ITYPE < OR = 0 THEN INITIALIZE HISTOGRAM OTHERWISE RETURN 
      IF(ITYPE)10,10,150
   10 IWRD1 = IABS(I) 
      DO 20 J = 1,IWRD1 
   20 NHIST(J) = 0
C 
C A(6) IS SET TO NUMBER OF BARS IN HISTOGRAM AND MODE OF HISTOGRAM
      IWRD2 = I 
      A(6) = RN 
C 
C A(4) AND A(5) ARE SET TO LOWER BOUND AND WIDTH OF HISTOGRAM 
      A(4) = START
      A(5) = DELTA
150   RETURN
      END 
      SUBROUTINE RCRDI(IDATA,A,IERR,NHIST)
C  "RCRDI" UPDATES THE RUNNING STATISTICS BEING MAINTAINED IN 
C   IN ARRAYS "NHIST" AND "A" WITH THE VALUE OF "IDATA".  THE 
C   FORMAT OF "A" IS DESCRIBED IN SUBROUTINE INTLI
C 
      DIMENSION A(5),N(2),NHIST(1)
      EQUIVALENCE (RN,N,ISTRT),(N(2),MODE,IDLTA)
C 
C   ASSUME NO ERRORS
      IERR = 0
C 
C   IF THE MODE OF THE STATISTICS > =  0 THEN UPDATE THE SUM
C   AND SUM OF THE SQUARES OF "IDATA".
      RN = A(3) 
      IF(MODE)100,200 
C 
C   UPDATE NUMBER OF POINTS, SUM AND SUM OF SQUARES 
200   N = N+1 
      A(3) = RN 
      DATA = IDATA
      A = A+DATA
      A(2) = A(2)+DATA*DATA 
C 
C   IF THE MODE < =  0 UPDATE HISTOGRAM, OTHERWISE RETURN 
      IF(MODE)100,100,600 
C 
C 
C  UPDATE HISTOGRAM 
C 
C   COMPUTE BAR NUMBER (M = BAR NUMBER - 1 )
100   RN = A(4) 
      M = (IDATA-ISTRT)/IDLTA 
      RN = A(5) 
C 
C  CHECK TO SEE IDATA IS IN BOUNDS OF HISTOGRAM 
      IF(M)400,301
301   IF(M-N)302,450
C 
C  IDATA GREATER THAN UPPER BOUND OF HISTOGRAM
C  IF HISTOGRAM MODE > 0 INCLUDE POINT IN LAST BAR COUNT
450   M = N-1 
      IF(MODE)601,402 
C 
C   IDATA LESS THAN LOWER BOUND OF HISTOGRAM
C  IF HISTOGRAM MODE > 0 INCLUDE POINT IN FIRST BAR COUNT 
400   IF(MODE)601,401 
401   M = 0 
402   IERR = 1
C 
C  INCREMENT PROPER BAR OF HISTOGRAM
302   NHIST(M+1) = NHIST(M+1)+1 
600   RETURN
601   IERR = -1 
      RETURN
      END 
      SUBROUTINE RCRDF(DATA,A,IERR,NHIST) 
C  "RCRDF" UPDATES THE RUNNING STATISTICS BEING MAINTAINED IN 
C   IN ARRAYS "NHIST" AND "A" WITH THE VALUE OF "DATA".  THE
C   FORMAT OF "A" IS DESCRIBED IN SUBROUTINE INTLB
C 
      DIMENSION A(6),N(2),NHIST(1)
      EQUIVALENCE (RN,N),(N(2),MODE)
C 
C   ASSUME NO ERRORS
      IERR = 0
C 
C   IF THE MODE OF THE STATISTICS > =  0 THEN UPDATE THE SUM
C   AND SUM OF THE SQUARES OF "DATA". 
      RN = A(3) 
      IF(MODE)100,200 
C 
C   UPDATE NUMBER OF POINTS, SUM AND SUM OF SQUARES 
200   N = N+1 
      A(3) = RN 
      A = A+DATA
      A(2) = A(2)+DATA*DATA 
C 
C   IF THE MODE < =  0 UPDATE HISTOGRAM, OTHERWISE RETURN 
      IF(MODE)100,100,600 
C 
C 
C  UPDATE HISTOGRAM 
C 
C   COMPUTE BAR NUMBER (M = BAR NUMBER - 1 )
100   M = (DATA-A(4))/A(5)
      RN = A(6) 
C 
C  CHECK TO SEE DATA IS IN BOUNDS OF HISTOGRAM
      IF(M)400,301
301   IF(M-N)302,450
C 
C  DATA GREATER THAN UPPER BOUND OF HISTOGRAM 
C  IF HISTOGRAM MODE > 0 INCLUDE POINT IN LAST BAR COUNT
450   M = N-1 
      IF(MODE)601,402 
C 
C   DATA LESS THAN LOWER BOUND OF HISTOGRAM 
C  IF HISTOGRAM MODE > 0 INCLUDE POINT IN FIRST BAR COUNT 
400   IF(MODE)601,401 
401   M = 0 
402   IERR = 1
C 
C  INCREMENT PROPER BAR OF HISTOGRAM
302   NHIST(M+1) = NHIST(M+1)+1 
600   RETURN
601   IERR = -1 
      RETURN
      END 
      SUBROUTINE REPRT(A,RMEAN,STDEV,NUM) 
C 
C  "REPRT" COMPUTES MEAN AND STANDARD DEVIATION FROM THE
C   SUM, SUM OF SQUARES AND NUMBER OF POINTS CONTAINED
C   IN ARRAY "A". 
      DIMENSION A(3)
C 
C  THIS EQUIVALENCE ALLOWS ACCESS TO THE FIRST WORD OF A FLOATING 
C   POINT NUMBER. 
      EQUIVALENCE (RN,N)
C 
C  THE FIRST WORD OF A(3) IS THE NUMBER OF POINTS 
      RN = A(3) 
      NUM = N 
      RN = N
      RMEAN = A/RN
      STDEV = SQRT((A(2)-A*RMEAN)/(RN-1.))
      RETURN
      END 
      END$
FTN,L,B 
      SUBROUTINE INTLB(ITYPE,A,START,DELTA,RHIST,I) 
C   "INTLB" INITIALIZES THE HISTOGRAM AND OTHER STATISTICAL 
C     ACCUMULATORS FOR RUNNING STATISTICS ON FLOATING POINT DATA. 
C  ARRAY "A" MUST BE PROVIDED FOR EACH VARIABLE FOR WHICH STATISTICS
C  ARE MAINTAINED.  "A" MUST BE A REAL ARRAY DIMENSIONED BY 6.
C  THE CONTENTS OF "A" ARE AS FOLLOWS:
C       A(1) - RUNNING SUM OF DATA
C       A(2) - RUNNING SUM OF SQUARES OF DATA 
C       A(3) - WORD 1 - NUMBER OF POINTS
C              WORD 2 - MODE OF STATISTICS
C       A(4) - LOWER BOUND OF HISTOGRAM 
C       A(5) - WIDTH OF EACH BAR IN HISTOGRAM 
C       A(6) - WORD 1 - NUMBER OF BARS IN HISTOGRAM 
C              WORD 2 - MODE OF HISTOGRAM 
      DIMENSION A(6),N(2),RHIST(1)
C 
C THE FOLLOWING EQUIVALENCE STATEMENT GIVES ACCESS TO THE TWO 
C WORDS OF THE FLOATING POINT VARIABLE "RN", WHERE IWRD1 IS THE 
C FIRST WORD AND IWRD2 IS THE SECOND WORD.
      EQUIVALENCE (RN,N,IWRD1),(N(2),IWRD2) 
C 
C INITIALIZE SUM AND SUM OF SQUARES TO 0.0
      A = 0.
      A(2) = 0. 
C 
C  INITIALIZE NUMBER OF POINTS TO 0 AND MODE OF STATISTICS TO ITYPE 
      IWRD1 = 0 
      IWRD2 = ITYPE 
      A(3) = RN 
C 
C  IF ITYPE < OR = 0 THEN INITIALIZE HISTOGRAM OTHERWISE RETURN 
      IF(ITYPE)10,10,150
   10 IWRD1 = IABS(I) 
      DO 20 J = 1,IWRD1 
   20 RHIST(J) = 0
C 
C A(6) IS SET TO NUMBER OF BARS IN HISTOGRAM AND MODE OF HISTOGRAM
      IWRD2 = I 
      A(6) = RN 
C 
C A(4) AND A(5) ARE SET TO LOWER BOUND AND WIDTH OF HISTOGRAM 
      A(4) = START
      A(5) = DELTA
150   RETURN
      END 
      SUBROUTINE RCRDB(DATA,A,IERR,RHIST) 
C  "RCRDB" UPDATES THE RUNNING STATISTICS BEING MAINTAINED IN 
C   IN ARRAYS "RHIST" AND "A" WITH THE VALUE OF "DATA".  THE
C   FORMAT OF "A" IS DESCRIBED IN SUBROUTINE INTLF
C 
      DIMENSION A(6),N(2),RHIST(1)
      EQUIVALENCE (RN,N),(N(2),MODE)
C 
C   ASSUME NO ERRORS
      IERR = 0
C 
C   IF THE MODE OF THE STATISTICS > =  0 THEN UPDATE THE SUM
C   AND SUM OF THE SQUARES OF "DATA". 
      RN = A(3) 
      IF(MODE)100,200 
C 
C   UPDATE NUMBER OF POINTS, SUM AND SUM OF SQUARES 
200   N = N+1 
      A(3) = RN 
      A = A+DATA
      A(2) = A(2)+DATA*DATA 
C 
C   IF THE MODE < =  0 UPDATE HISTOGRAM, OTHERWISE RETURN 
      IF(MODE)100,100,600 
C 
C 
C  UPDATE HISTOGRAM 
C 
C   COMPUTE BAR NUMBER (M = BAR NUMBER - 1 )
100   M = (DATA-A(4))/A(5)
      RN = A(6) 
C 
C  CHECK TO SEE DATA IS IN BOUNDS OF HISTOGRAM
      IF(M)400,301
301   IF(M-N)302,450
C 
C  DATA GREATER THAN UPPER BOUND OF HISTOGRAM 
C  IF HISTOGRAM MODE > 0 INCLUDE POINT IN LAST BAR COUNT
450   M = N-1 
      IF(MODE)601,402 
C 
C   DATA LESS THAN LOWER BOUND OF HISTOGRAM 
C  IF HISTOGRAM MODE > 0 INCLUDE POINT IN FIRST BAR COUNT 
400   IF(MODE)601,401 
401   M = 0 
402   IERR = 1
C 
C  INCREMENT PROPER BAR OF HISTOGRAM
302   RHIST(M+1) = RHIST(M+1)+1.
600   RETURN
601   IERR = -1 
      RETURN
      END 
      END$