Google
 

Trailing-Edge - PDP-10 Archives - decus_20tap5_198111 - decus/20-0137/stp/stp5.for
There is 1 other file named stp5.for in the archive. Click here to see a list.
C                                            *** STAT PACK ***
C     SUBROUTINE FOR MANN-WHITNEY U'S
C     CALLING SEQUENCE: CALL MANN(NV,NC,MV,MC,DATA,AV,BV,NAMES)
C     WHERE NV - NUMBER OF VARIABLES
C           NC - NUMBER OF OBSERVATIONS
C           MV - MAXIMUM NUMBER OF VARIABLES (AS SPECIFIFED IN MAIN)
C           MC - MAXIMUM NUMBER OF OBSERVATIONS(AS SPECIFIED IN MAIN)
C           DATA - STORAGE FOR DATA DIMENSIONED FOR MAXIMUM
C           AV - VECTOR DIMENSIONED AT LEAST NV LONG
C           BV - VECTOR DIMENSIONED AT LEAST NV LONG
C           NAMES - VECTOR CONTAINING NAMES OF VARIABLES
C     
C     MANN-WHITNEY U TEST AS SPECIFIED IN BASIC STATISTICAL METHODS
C     BY DOWNE AND HEATH. INDEPENDENT SAMPLES ARE DETERMINED BY TWO
C     DIFFERENT METHODS: 1) EACH VARIABLE FORMS A SAMPLE OF 2) ONE
C     VARIABLE IS BROKEN INTO THE SAMPLE BY MEANS OF THE VALUE OF 
C     A SECOND VARIABLE.  IN THE FIRST CASE THE VARIABLE 1 IS CHOSEN
C     AND SORTED THEN 2 IS CHOSEN AND SORTED, THEY ARE RANKED
C     WITHOUT REPLACEING THEE VALUES WITH THE RANKS; THE RANKS
C     ARE SIMPLYCALCULATED AND ADDED INTO THE CORRECT TOTAL.  THE 
C     MEANS STANDARD DEVIATIONS AND Z'S ARE CALCULATED AND OUTPUT, THEN
C     VARIABLE 3 IS SORTED RANKED, REPLACEING 2, WHEN ALL VARIABLES
C     HAVE BEEN PROCESSED, VARIABLE 2 REPLACES 1 AND THE SAME
C     METHOD OF ANALYSIS IS EMPLOYED.  IN THE SECOND METHOD OF 
C     ANALYSIS, WHERE SAMPLES ARE DETERMINED BY A BREAK A TWO
C     VECTOR SORT IS USED IN A MANNER SIMILIAR TO THE XTAB PROGRAM
C     FOR STP.  THESE ARE THEN TREATED IN THE SAME FASHION AS THE
C     FIRST CASE.  SORTS USED ARE THE ACM PARTITIONING.  DOWNE AND 
C     HEATH BASIC STATISTICAL METHODS PAGES 213,214, ALSO 
C     NON PARAMETERIC STATISTICS BY SIEGEL PAGES 119-124.
C
      SUBROUTINE MANN(NV,NC,MV,MC,DATA,AV,BV,NAMES)
      COMMON /DEV/ICC,IDATA,IOUT,IDLG,IDSK
      COMMON /PRNT/ LINPP,ICOPS,RUNPRG
      COMMON /EXTRA/ HEDR(70),NSZ
      DIMENSION DATA(MC,MV),AV(1),BV(1),OPT(6),NAMES(1)
      DIMENSION IA(20),RNG(100,2)
C     FOLLOWING DIMENSION FOR SORTS ONLY
      DIMENSION IU(16),IL(16)
      NAUTO=0
      ALL=0
1     IF(ICC.NE.2) WRITE(IDLG,2)
2     FORMAT('0LIST OPTIONS SEPARATED BY COMMAS'/)
      DISCR=0
      BREAK=0
      AUTO=0
      RANGE=0
      READ(ICC,3) OPT
3     FORMAT(10(A5,1X))
      IF(OPT(1).EQ.'!') RETURN
      DO 4 I=1,6
      IF(OPT(I).EQ.' ') GO TO 20
      IF(OPT(I).NE.'HELP') GO TO 7
      WRITE(IDLG,5)
5     FORMAT('0MANN-WHITNEY U-TEST ASSUMES EACH VARIABLE TO',
     1' BE AN INDEPENDENT'/' SAMPLE, AND WILL CALCULATE THE MEANS,',
     2' STANDARD DEVIATIONS AND'/' U''S FOR ALL COMBINATIONS.  OPTIONS',
     3' EXIST TO CHANGE THIS BY CREATING'/' INDEPENDENT',
     4' SAMPLES FROM ONE VARIABLE BASED ON VALUES OF A SECOND'/
     5' VARIABLE.  IN THIS CASE THE PROGRAM WILL ACCEPT RANGES OF',
     6' VALUES'/' FOR THE SECOND VARIABLE TO BREAK ON.  IF YOU DESIRE',
     7' THE VARIABLE'/' TO BE BROKE INTO INDEPENDENT SAMPLES BASED ON',
     8' ALL DIFFERENT VALUES'/' OF THE SECOND VARIABLE YOU MAY',
     9' SPECIFY "DISCR" AS AN OPTION OR'/' "AUTO" WHEN ASKED FOR THE',
     1' RANGES.  IF "DISCR" OR "AUTO" IS USED THE'/' PROGRAM',
     2' WILL DO THE BREAKDOWN WITHOUT ANY ASSISTANCE, IF YOU'/' DESIRE',
     3' A LIST OF THE RANGES THE PROGRAM USES, TYPE',
     4' "RANGE" AS AN'/' OPTION')
      WRITE(IDLG,6)
6     FORMAT('0OPTIONS ARE:'/' "BREAK" - CREATE INDEPENDENT',
     1' SAMPLES ON VALUES OF A SECOND VARIABLE'/' "DISCR" - ALLOW FOR',
     2' BREAKDOWNS BASED ON INDIVIDUAL VALUES'/11X,' OF THE BREAKDOWN',
     3' VARIABLE'/' "RANGE" - LIST RANGES WHEN AUTOMATIC BREAKDOWN IS',
     4' USED'/' "AUTO" - AUTOAMATIC BREAKDOWN (SPECIFIED WHEN ASKED',
     5' FOR RANGES)'/
     6'0IF NO OPTIONS ARE DESIRED TYPE A RETURN')
      GO TO 1
7     IF(OPT(I).NE.'BREAK') GO TO 8
      BREAK=1
      GO TO 4
8     IF(OPT(I).NE.'RANGE') GO TO 9
      RANGE=1
      GO TO 4
9     IF(OPT(I).NE.'DISCR') GO TO 10
      DISCR=1
      GO TO 4
10    IF(OPT(I).NE.'AUTO') GO TO 13
11    WRITE(IDLG,12)
12    FORMAT('0"AUTO" MAY ONLY BE USED WHEN ASKED FOR RANGES')
      GO TO 1
13    IF(OPT(I).EQ.'AUTO,') GO TO 11
      WRITE(IDLG,14)OPT(I)
14    FORMAT('0OPTION "',A5,'" DOES NOT EXIST')
      GO TO 1
4     CONTINUE
20    IF(BREAK.EQ.1) GO TO 80
C
C     EACH VARIABLE IS TREATED AS AN INDEPENDENT SAMPLE
C
      IF(NV.LT.2) RETURN
      IF(IOUT.NE.21) WRITE(IOUT,5566) (HEDR(I),I=1,NSZ)
5566  FORMAT('1',80A1)
      IF(IOUT.EQ.21) CALL PRNTHD
      WRITE(IOUT,22)
22    FORMAT('0',10X,'***** MANN-WHITNEY U-TEST *****')
      IF(IOUT.NE.21)WRITE(IOUT,25)
25    FORMAT('0VAR. VS VAR.',7X,'MEAN',8X,'STANDARD DEVIATION',6X,'Z'/
     117X,'N1',4X,'N2',7X,'U1',16X,'U2'/1X,65(1H-))
      IF(IOUT.EQ.21) WRITE(IOUT,225)
225   FORMAT('0VAR. VS VAR.',7X,'MEAN',8X,'STANDARD DEVIATION',6X,'Z',
     113X,'N1',4X,'N2',8X,'U1',15X,'U2'/1X,113(1H-))
      LINES=7
      DO 23 I=2,NV
C     
C     SORT FIRST VECTOR
C
      DO 24 J=1,NC
24    AV(J)=DATA(J,I)
      M=1
      II=1
      J=NC
31    IF(II.GE.J) GO TO 38
32    K=II
      IJ=(J+II)/2
      T=AV(IJ)
      IF(AV(II).LE.T) GO TO 33
      AV(IJ)=AV(II)
      AV(II)=T
      T=AV(IJ)
33    LL=J
      IF(AV(J).GE.T) GO TO 35
      AV(IJ)=AV(J)
      AV(J)=T
      T=AV(IJ)
      IF(AV(II).LE.T) GO TO 35
      AV(IJ)=AV(II)
      AV(II)=T
      T=AV(IJ)
      GO TO 35
34    AV(LL)=AV(K)
      AV(K)=TT
35    LL=LL-1
      IF(AV(LL).GT.T) GO TO 35
      TT=AV(LL)
36    K=K+1
      IF(AV(K).LT.T) GO TO 36
      IF(K.LE.LL) GO TO 34
      IF((LL-II).LE.(J-K)) GO TO 37
      IL(M)=II
      IU(M)=LL
      II=K
      M=M+1
      GO TO 39
37    IL(M)=K
      IU(M)=J
      J=LL
      M=M+1
      GO TO 39
38    M=M-1
      IF(M.EQ.0) GO TO 43
      II=IL(M)
      J=IU(M)
39    IF((J-II).GE.11) GO TO 32
      IF(II.EQ.1) GO TO 31
      II=II-1
40    II=II+1
      IF(II.EQ.J) GO TO 38
      T=AV(II+1)
      IF(AV(II).LE.T) GO TO 40
      K=II
41    AV(K+1)=AV(K)
      K=K-1
      IF(T.LT.AV(K)) GO TO 41
      AV(K+1)=T
      GO TO 40
C
C     END SORT ON FIRST VARIABLE
C
43    DO 23 N=1,I-1
      DO 42 J=1,NC
42    BV(J)=DATA(J,N)
C
C     SORT SECOND VARIABLE
C
      M=1
      II=1
      J=NC
51    IF(II.GE.J) GO TO 58
52    K=II
      IJ=(J+II)/2
      T=BV(IJ)
      IF(BV(II).LE.T) GO TO 53
      BV(IJ)=BV(II)
      BV(II)=T
      T=BV(IJ)
53    LL=J
      IF(BV(J).GE.T) GO TO 55
      BV(IJ)=BV(J)
      BV(J)=T
      T=BV(IJ)
      IF(BV(II).LE.T) GO TO 55
      BV(IJ)=BV(II)
      BV(II)=T
      T=BV(IJ)
      GO TO 55
54    BV(LL)=BV(K)
      BV(K)=TT
55    LL=LL-1
      IF(BV(LL).GT.T) GO TO 55
      TT=BV(LL)
56    K=K+1
      IF(BV(K).LT.T) GO TO 56
      IF(K.LE.LL) GO TO 54
      IF((LL-II).LE.(J-K)) GO TO 57
      IL(M)=II
      IU(M)=LL
      II=K
      M=M+1
      GO TO 59
57    IL(M)=K
      IU(M)=J
      J=LL
      M=M+1
      GO TO 59
58    M=M-1
      IF(M.EQ.0) GOTO 62
      II=IL(M)
      J=IU(M)
59    IF((J-II).GE.11) GOTO 52
      IF(II.EQ.1) GOTO 51
      II=II-1
60    II=II+1
      IF(II.EQ.J) GO TO 58
      T=BV(II+1)
      IF(BV(II).LE.T) GOTO 60
      K=II
61    BV(K+1)=BV(K)
      K=K-1
      IF(T.LT.BV(K)) GO TO 61
      BV(K+1)=T
      GOTO 60
C
C     END SORT ON SECOND VARIABLE
C
C
C     BEGIN RANKING AND WHILE RANKING SUM UP RANKS OF U1
C     (SUM OF FIRST VAR RANKS)
C
62    K=1
      L=1
      SUMU2=0
      IRNK=1
      NUM=0
      NUMAV=0
      TOTAL=0
      T=AV(K)
      IF(T.GT.BV(L)) T=BV(L)
63    IF(K.GT.NC) GO TO 64
      IF(T.NE.AV(K)) GO TO 64
      NUM=NUM+1
      TOTAL=TOTAL+IRNK
      IRNK=IRNK+1
      K=K+1
      NUMAV=NUMAV+1
      GO TO 63
64    IF(L.GT.NC) GO TO 65
      IF(T.NE.BV(L)) GO TO 65
      NUM=NUM+1
      TOTAL=TOTAL+IRNK
      IRNK=IRNK+1
      L=L+1
      GO TO 64
65    RANK=TOTAL/NUM
      SUMU2=SUMU2+RANK*NUMAV
      TOTAL=0
      NUM=0
      NUMAV=0
      IF(K.GT.NC) GO TO 66
      T=AV(K)
      IF(L.GT.NC) GO TO 63
      IF(T.GT.BV(L)) T=BV(L)
      GO TO 63
66    XMN=(NC**2)/2.
      STD=SQRT(((NC**2)*(2.*NC+1.))/12.)
      SUMU1=(2.*NC*(2.*NC+1.))/2.-SUMU2
      U1=NC**2+(NC*(NC+1.))/2.-SUMU1
      U2=NC**2+(NC*(NC+1.))/2.-SUMU2
      Z=(U1-XMN)/STD
      IF(IOUT.NE.21) GO TO 68
      LINES=LINES+1
      IF(LINES.LE.LINPP) GO TO 68
      CALL PRNTHD
      WRITE(IOUT,225)
      LINES=6
68    IF(IOUT.NE.21)WRITE(IOUT,67) NAMES(N),NAMES(I),XMN,STD,Z,
     1NC,NC,U1,U2
67    FORMAT(1X,A5,2X,A5,3X,G15.7,2X,G15.7,3X,G15.7/15X,I4,2X,
     1I4,3X,G15.7,2X,G15.7)
      IF(IOUT.EQ.21) WRITE(IOUT,267) NAMES(N),NAMES(I),
     1XMN,STD,Z,NC,NC,U1,U2
267   FORMAT(1X,A5,2X,A5,3X,G15.7,2X,G15.7,3X,G15.7,2X,I4,2X,I4,3X,
     1G15.7,2X,G15.7)
23    CONTINUE
      RETURN
C
C     END OF ANALYSIS TREATING EACH VARIABLE AS AN INDEPENDENT SAMPLE
C     @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
C
C     BEGIN INDEPENDENT SAMPLES DUE TO BREAKDOWN
C
80    IF(ICC.NE.2) WRITE(IDLG,81)
81    FORMAT('0WHICH VARIABLES ARE TO BE ANALYSED? ',$)
      CALL ALPHA(IA,20,NANALY,IRET,IHELP,IERR,NAMES,NV)
      IF(IRET.EQ.1) RETURN
      IF(IERR.EQ.1) GO TO 80
      IF(IHELP.EQ.1) GO TO 80
      DO 82 I=1,NANALY
      IF(IA(I).GT.0) GO TO 82
      NANALY=NV
      ALL=1
      GO TO 87
82    CONTINUE
87    IF(ICC.NE.2) WRITE(IDLG,88)
88    FORMAT('0WHICH VARIABLE IS THE BREAKDOWN VARIABLE? ',$)
      CALL ALPHA(IB,1,I,IRET,IHELP,IERR,NAMES,NV)
      IF(IRET.EQ.1) RETURN
      IF(IERR.EQ.1) GO TO 87
      IF(IHELP.EQ.1) GO TO 87
      IF(IB.GT.0) GO TO 91
      WRITE(IDLG,89)
89    FORMAT(' ?, *, ALL MAY NOT BE USED HERE')
      GO TO 87
91    IF(DISCR.EQ.1) GO TO 103
C
C     USER TO ENTER RANGES
C
      I=1
96    IF(ICC.NE.2) WRITE(IDLG,92) NAMES(IB)
92    FORMAT('0ENTER RANGES FOR BREAKDOWN VARIABLE: ',A5/)
93    IF(ICC.NE.2) WRITE(IDLG,94)
94    FORMAT('+? ',$)
      READ(ICC,3,END=102) ALANS
      IF(ALANS.EQ.'!') RETURN
      IF(ALANS.EQ.'AUTO') GO TO 103
      IF(ALANS.EQ.'ALL') GO TO 103
      IF(ALANS.NE.'HELP') GO TO 97
      WRITE(IDLG,95)
95    FORMAT('0IF AUTOAMATIC BREAKDOWN IS DESIRED TYPE "AUTO"',
     1' OTHERWISE ENTER'/' THE RANGE FOR EACH BREAKDOWN SMALLER',
     2' FIRST ONE PER LINE.'/' EXAMPLE:'/' 17,35.2'/' CONTROL Z',
     3' (^Z), "STOP", OR "END" WILL TERMINATE ENTRY'//)
      GO TO 93
97    IF(ALANS.EQ.'STOP') GO TO 102
      IF(ALANS.EQ.'END') GO TO 102
      REREAD 98,(RNG(I,J),J=1,2)
98    FORMAT(2F)
      IF(RNG(I,1).LE.RNG(I,2)) GO TO 100
      WRITE(IDLG,99)
99    FORMAT(' SMALLER FIRST PLEASE REENTER LAST RANGE'/)
      GO TO 93
100   I=I+1
      IF(I.LE.100) GO TO 93
      WRITE(IDLG,101) 
101   FORMAT(' NO MORE RANGES ACCEPTED')
102   NAUTO=1
      NR=I-1
C
C     BEGIN ANALYSIS
C
103   DO 104 I=1,NANALY
      IF(ALL.EQ.1) GO TO 105
      N=IA(I)
      GO TO 106
105   IF(I.EQ.IB) GO TO 104
      N=I
106   K=1
      DO 108 J=1,NC
      BV(K)=DATA(J,IB)
      AV(K)=DATA(J,N)
      IF(NAUTO.NE.1) GO TO 110
      DO 109 L=1,NR
      IF(BV(K).LT.RNG(L,1)) GO TO 109
      IF(BV(K).GT.RNG(L,2)) GO TO 109
      BV(K)=L
      GO TO 110
109   CONTINUE
      GO TO 108
110   K=K+1
108   CONTINUE
C
C     AT THIS POINT RATHER THE USER ENTERED THE TANGES OR ALLOWED
C     THE MACHINE TO DO IT FOR HIM THE RESULTS ARE THE SAME.
C
C     SORT NOW MAJOR ON BV(BREAKDOWN VARIABLE) MINOR ON AV (ANALYSIS
C     VARIABLE)
C
      NCOO=K-1
      M=1
      II=1
      J=NCOO
111   IF(II.GE.J) GO TO 118
112   K=II
      IJ=(J+II)/2
      T=BV(IJ)
      S=AV(IJ)
      IF(BV(II).GT.T) GO TO 123
      IF(BV(II).LT.T) GO TO 113
      IF(AV(II).LE.S) GO TO 113
123   BV(IJ)=BV(II)
      BV(II)=T
      T=BV(IJ)
      AV(IJ)=AV(II)
      AV(II)=S
      S=AV(IJ)
113   LL=J
      IF(BV(J).LT.T) GO TO 124
      IF(BV(J).GT.T) GO TO 115
      IF(AV(J).GE.S) GO TO 115
124   BV(IJ)=BV(J)
      BV(J)=T
      T=BV(IJ)
      AV(IJ)=AV(J)
      AV(J)=S
      S=AV(IJ)
      IF(BV(II).GT.T) GO TO 125
      IF(BV(II).LT.T) GO TO 115
      IF(AV(II).LE.S) GO TO 115
125   BV(IJ)=BV(II)
      BV(II)=T
      T=BV(IJ)
      AV(IJ)=AV(II)
      AV(II)=S
      S=AV(IJ)
      GO TO 115
114   BV(LL)=BV(K)
      BV(K)=TT
      AV(LL)=AV(K)
      AV(K)=SS
115   LL=LL-1
      IF(BV(LL).GT.T) GO TO 115
      IF((BV(LL).EQ.T).AND.(AV(LL).GT.S)) GO TO 115
      TT=BV(LL)
      SS=AV(LL)
116   K=K+1
      IF(BV(K).LT.T) GO TO 116
      IF((BV(K).EQ.T).AND.(AV(K).LT.S)) GO TO 116
      IF(K.LE.LL) GO TO 114
      IF((LL-II).LE.(J-K)) GO TO 117
      IL(M)=II
      IU(M)=LL
      II=K
      M=M+1
      GO TO 119
117   IL(M)=K
      IU(M)=J
      J=LL
      M=M+1
      GO TO119
118   M=M-1
      IF(M.EQ.0) GO TO 130
      II=IL(M)
      J=IU(M)
119   IF((J-II).GE.11) GO TO 112
      IF(II.EQ.1) GOTO 111
      II=II-1
120   II=II+1
      IF(II.EQ.J) GO TO 118
      S=AV(II+1)
      T=BV(II+1)
      IF(BV(II).LT.T) GO TO 120
      IF((BV(II).EQ.T).AND.(AV(II).LE.S)) GO TO 120
      K=II
121   BV(K+1)=BV(K)
      AV(K+1)=AV(K)
      K=K-1
      IF(T.LT.BV(K)) GO TO 121
      IF((T.EQ.BV(K)).AND.(S.LT.AV(K))) GO TO 121
      BV(K+1)=T
      AV(K+1)=S
      GO TO 120
C
C     END OF SORT
C     RANGE WAS USED?
C     
130   IF((NAUTO.EQ.1).OR.(RANGE.NE.1)) GO TO 140
      T=BV(1)
      WRITE(IDLG,131) NAMES(IB)
131   FORMAT('0BREAKDOWN RANGES FOR VARIABLE: ',A5)
      WRITE(IDLG,132) T,T
132   FORMAT(1X,G10.4,',',G10.4)
      DO 133 J=2,NCOO
      IF(T.EQ.BV(J)) GO TO 133
      T=BV(J)
      WRITE(IDLG,132) T,T
133   CONTINUE
C
C     FORM U'S
C
140   IF(IOUT.NE.21) WRITE(IOUT,5566)(HEDR(J),J=1,NSZ)
      IF(IOUT.EQ.21) CALL PRNTHD
      WRITE(IOUT,22)
      WRITE(IOUT,107) NAMES(N),NAMES(IB)
107   FORMAT(' ANALYSIS ON VARIABLE: ',A5,' WITH INDEPENDENT SAMPLES',
     1' DETERMINED '/' BY A BREAKDOWN ON VARIABLE: ',A5)
      IF(IOUT.NE.21)WRITE(IOUT,25)
      IF(IOUT.EQ.21) WRITE(IOUT,225)
      LINES=9
      N1=1
      IV1=0
141   T=BV(N1)
      TT=T
      IV1=IV1+1
      IV2=IV1
      N2=N1
143   J=N2
      DO 142 N2=J,NCOO
      IF(TT.NE.BV(N2)) GO TO 151
142   CONTINUE
      J=N1
      DO 144 N1=J,NCOO
      IF(T.NE.BV(N1)) GO TO 141
144   CONTINUE
      GO TO 104
151   TT=BV(N2)
      IV2=IV2+1
      K=N1
      L=N2
C
C     RANK THE TWO NOW KEEPING ONLY THOSE RANKS FOUND BY INDEX N1
C
      SUMU1=0
      IRNK=1
      NUM=0
      NUMAV=0
      NUMU1=0
      NUMU2=0
      TOTAL=0
      S=AV(K)
      IF(S.GT.AV(L)) S=AV(L)
145   IF(T.NE.BV(K)) GO TO 146
      IF(S.NE.AV(K)) GO TO 146
      NUM=NUM+1
      TOTAL=TOTAL+IRNK
      IRNK=IRNK+1
      K=K+1
      NUMAV=NUMAV+1
      NUMU1=NUMU1+1
      GO TO 145
146   IF(L.GT.NCOO) GO TO 147
      IF(TT.NE.BV(L)) GO TO 147
      IF(S.NE.AV(L)) GO TO 147
      NUM=NUM+1
      TOTAL=TOTAL+IRNK
      IRNK=IRNK+1
      L=L+1
      NUMU2=NUMU2+1
      GO TO 146
147   RANK=TOTAL/NUM
      SUMU1=SUMU1+RANK*NUMAV
      TOTAL=0
      NUM=0
      NUMAV=0
      IF(BV(K).EQ.T) GO TO 148
      IF(BV(L).NE.TT) GO TO 150
      IF(L.GT.NCOO) GO TO 150
      S=AV(L)
      GO TO 145
148   S=AV(K)
      IF(L.GT.NCOO) GO TO 145
      IF(BV(L).NE.TT) GO TO 145
      IF(S.GT.AV(L)) S=AV(L)
      GO TO 145
150   XMN=(NUMU1*NUMU2)/2.
      STD=SQRT(((NUMU1*NUMU2)*(NUMU1+NUMU2+1.))/12.)
      SUMU2=((NUMU1+NUMU2)*(NUMU1+NUMU2+1))/2.-SUMU1
      U1=NUMU1*NUMU2+(NUMU1*(NUMU1+1.))/2.-SUMU1
      U2=NUMU1*NUMU2+(NUMU2*(NUMU2+1.))/2.-SUMU2
      Z=(U1-XMN)/STD
      ENCODE(5,152,NIV1) IV1
      ENCODE(5,152,NIV2) IV2
152   FORMAT(I4,1X)
      IF(IOUT.NE.21) GO TO 154
      LINES=LINES+1
      IF(LINES.LE.LINPP) GO TO 154
      CALL PRNTHD
      WRITE(IOUT,225)
      LINES=6
154   IF(IOUT.NE.21)WRITE(IOUT,67) NIV1,NIV2,XMN,STD,Z,NUMU1,NUMU2,U1,U2
      IF(IOUT.EQ.21) WRITE(IOUT,267) NIV1,NIV2,XMN,STD,Z,NUMU1,NUMU2
     1,U1,U2
      GO TO 143
104   CONTINUE
      RETURN
      END
C                                     *** STAT PACK ***
C     SUBROUTINE FOR SPEARMAN RANK ORDER CORRELATIONS
C     CALLING SEQUENCE: CALL STSRNK(NV,NC,MV,MC,DATA,LM,R,NAMES)
C     WHERE NV - IS THE NUMBER OF VARIABLES.
C            NC - IS THE NUMBER OF OBSERVATIONS.
C            MV - IS MAXIMUM VARIABLES, AS SPECIFIED IN MAIN.
C            MC - IS THE MAXIMUM OBSERVATIONS, AS SPECIFIED IN MAIN.
C            DATA- IS STORAGE FOR DATA, DIMENSIONED FOR MAXIMUM
C            LM - IS A VECTOR AT LEAST NV LONG.
C            R - IS A VECTOR AT LEAST NV LONG.
C           NAMES - IS A VECTOR CONTAINING VARIABLE NAMES
C
C     THE CORRELATIONS ARE GIVEN FOR ALL VARIABLES IN THE 
C     FORM OF A MATRIX.  RANKS ARE COMPLETED BY USING THE ACM
C     SORT.  FIRST VECTOR SORTED IS RANKED AND SAVED. THE REMAINDER
C     OF THE ROW IS RANKED, BUT ONLY SAVED LONG ENOUGH TO FIND THE
C     DEVIATIONS SQUARED, VARIABLE BY VARIABLE.  SOURCE FOR THE
C     SPEARMAN RANK CORRELATION METHOD  IS FOUND IN BASIC 
C     STATISTICAL METHODS BY DOWNE AND HEATH PAGES 178 AND 179.
C     HARPER & BROTHERS, PUBLISHERS 1959.
C
C     ORIGINAL SOURCE DOES NOT ADJUST FOR TIES IN RANKS, AS POINTED
C     OUT BY JANET C. WEBER PH. D., STATISTICAL CONSULTANT
C     INDIANA UNIVERSITY, PURDUE UNIVERSITY AT INDIANAPOLIS.  AS A RESULT
C     A SWITCH HAS BEEN MADE TO "NON-PARAMETRIC STATISTICS" 
C     BY SIEGEL AS THE REFERENCE.  PROGRAM HAS BEEN ALTERED
C     TO ADJUST FOR TIES IN RANKING.
C
      SUBROUTINE STSRNK(NV,NC,MV,MC,DATA,LM,R,NAMES)
      DIMENSION DATA(MC,MV),R(1),LM(1),C(15),NAMES(1)
C     FOLLOWING DIMENSION FOR SORTS ONLY
      DIMENSION IL(16),IU(16)
      COMMON /DEV/ICC,IDATA,IOUT,IDLG,IDSK
      COMMON /PRNT/ LINPP,ICOPS,RUNPRG
      COMMON /EXTRA/HEDR(70),NSZ
      ISQ=7
      IF(IOUT.EQ.21) ISQ=15
      TC=NC
      IF(IOUT.NE.21) WRITE(IOUT,5566)(HEDR(I),I=1,NSZ)
5566  FORMAT('1',70A1)
      IF(IOUT.EQ.21) CALL PRNTHD
      WRITE(IOUT,1)
1     FORMAT('0',10X,'***** SPEARMAN RANK-ORDER CORR. MATRIX *****')
      LINES=4
      DO 2 I=1,NV
      DO 4 J=1,NC
4     LM(J)=J
C
C     SORT INDIRECTLY THE FIRST VARIABLE (I) ACM SORT PARTITIONING
C
      M=1
      II=1
      J=NC
11    IF(II.GE.J) GO TO 18
12    K=II
      IJ=(J+II)/2
      T=DATA(LM(IJ),I)
      IF(DATA(LM(II),I).LE.T) GO TO 13
      ISAV=LM(IJ)
      LM(IJ)=LM(II)
      LM(II)=ISAV
      T=DATA(LM(IJ),I)
13    LL=J
      IF(DATA(LM(J),I).GE.T) GO TO 15
      ISAV=LM(IJ)
      LM(IJ)=LM(J)
      LM(J)=ISAV
      T=DATA(LM(IJ),I)
      IF(DATA(LM(II),I).LE.T) GO TO 15
      ISAV=LM(IJ)
      LM(IJ)=LM(II)
      LM(II)=ISAV
      T=DATA(LM(IJ),I)
      GO TO 15
14    ISAV=LM(LL)
      LM(LL)=LM(K)
      LM(K)=ISAV
15    LL=LL-1
      IF(DATA(LM(LL),I).GT.T) GO TO 15
      TT=DATA(LM(LL),I)
16    K=K+1
      IF(DATA(LM(K),I).LT.T) GO TO 16
      IF(K.LE.LL) GO TO 14
      IF((LL-II).LE.(J-K)) GO TO 17
      IL(M)=II
      IU(M)=LL
      II=K
      M=M+1
      GO TO 19
17    IL(M)=K
      IU(M)=J
      J=LL
      M=M+1
      GO TO 19
18    M=M-1
      IF(M.EQ.0) GO TO 23
      II=IL(M)
      J=IU(M)
19    IF((J-II).GE.11) GO TO 12
      IF(II.EQ.1) GO TO 11
      II=II-1
20    II=II+1
      IF(II.EQ.J) GOTO 18
      NEXTRA=LM(II+1)
      T=DATA(LM(II+1),I)
      IF(DATA(LM(II),I).LE.T) GO TO 20
      K=II
21    LM(K+1)=LM(K)
      K=K-1
      IF(T.LT.DATA(LM(K),I)) GO TO 21
      LM(K+1)=NEXTRA
      GO TO 20
C
C     END OF FIRST SORT, RANK NOW
C
23    K=1
      T=DATA(LM(1),I)
      TOTAL=NC
      NUM=1
      TOTX=0
      DO 24 J=2,NC
      IRNK=NC-J+1
      IF(T.EQ.DATA(LM(J),I)) GO TO 26
      RANK=TOTAL/NUM
      IF(NUM.GT.1) TOTX=TOTX+(NUM**3-NUM)/12.
      DO 25 L=K,J-1
25    R(LM(L))= RANK
      NUM=1
      K=J
      TOTAL=IRNK
      T=DATA(LM(J),I)
      GO TO 24
26    NUM=NUM+1
      TOTAL=TOTAL+IRNK
24    CONTINUE
      RANK=TOTAL/NUM
      IF(NUM.GT.1) TOTX=TOTX+(NUM**3-NUM)/12.
      DO 27 L=K,NC
27    R(LM(L))=RANK
      SUMX2=(TC**3-TC)/12.-TOTX
      IF(IOUT.NE.21) GO TO 28
      LL=(I+ISQ-1)/ISQ+1
      LINES=LINES+LL
      IF(LINES.LE.(LINPP-LL-1)) GO TO 28
      WRITE(IDLG,62)
      DO 29 NN=1,I-1,ISQ
      NEND=NN+ISQ-1
      IF(NEND.GT.(I-1)) NEND=I-1
29    WRITE(IOUT,61)(NAMES(J),J=NN,NEND)
      CALL PRNTHD
      LINES=2+LL
28    DO 2 NN=1,I,ISQ
      NEND=NN+ISQ-1
      IF(NEND.GT.I) NEND=I
      NCOR=NEND-NN+1
      DO 3 N=NN,NEND
      IF(N.EQ.I) GO TO 55
C
C     SORT ON THE VARIABLE TO BE COMPARED AGAINST I, RANK IT AND 
C     FIND THE CORRELATION, IN THIS MANNER THE FIRST VECTOR NEED
C     BE RANKED ONLY ONCE FOR THE TOTAL ROW TO BE CORRELATED.
C
C     SORT FIRST SAME AS LAST TIME ACM SORT
C
      DO 30 J=1,NC
30    LM(J)=J
      M=1
      II=1
      J=NC
31    IF(II.GE.J) GO TO 38
32    K=II
      IJ=(J+II)/2
      T=DATA(LM(IJ),N)
      IF(DATA(LM(II),N).LE.T) GO TO 33
      ISAV=LM(IJ)
      LM(IJ)=LM(II)
      LM(II)=ISAV
      T=DATA(LM(IJ),N)
33    LL=J
      IF(DATA(LM(J),N).GE.T) GO TO 35
      ISAV=LM(IJ)
      LM(IJ)=LM(J)
      LM(J)=ISAV
      T=DATA(LM(IJ),N)
      IF(DATA(LM(II),N).LE.T) GOTO 35
      ISAV=LM(IJ)
      LM(IJ)=LM(II)
      LM(II)=ISAV
      T=DATA(LM(IJ),N)
      GO TO 35
34    ISAV=LM(LL)
      LM(LL)=LM(K)
      LM(K)=ISAV
35    LL=LL-1
      IF(DATA(LM(LL),N).GT.T) GO TO 35
      TT=DATA(LM(LL),N)
36    K=K+1
      IF(DATA(LM(K),N).LT.T) GOTO 36
      IF(K.LE.LL) GOTO 34
      IF((LL-II).LE.(J-K)) GO TO 37
      IL(M)=II
      IU(M)=LL
      II=K
      M=M+1
      GOTO 39
37    IL(M)=K
      IU(M)=J
      J=LL
      M=M+1
      GO TO 39
38    M=M-1
      IF(M.EQ.0) GO TO 45
      II=IL(M)
      J=IU(M)
39    IF((J-II).GE.11) GOTO 32
      IF(II.EQ.1) GO TO 31
      II=II-1
40    II=II+1
      IF(II.EQ.J) GO TO 38
      NEXTRA=LM(II+1)
      T=DATA(LM(II+1),N)
      IF(DATA(LM(II),N).LE.T) GO TO 40
      K=II
41    LM(K+1)=LM(K)
      K=K-1
      IF(T.LT.DATA(LM(K),N)) GO TO 41
      LM(K+1)=NEXTRA
      GO TO 40
C
C     END SECOND SORT - RANK AND FIND DEVIATIONS SQUARE
C
45    NUM=1
      K=1
      D2=0
      T=DATA(LM(1),N)
      TOTAL=NC
      TOTY=0
      DO 46 J=2,NC
      IRNK=NC-J+1
      IF(T.EQ.DATA(LM(J),N)) GO TO 48
      RANK=TOTAL/NUM
      IF(NUM.GT.1) TOTY=TOTY+(NUM**3-NUM)/12.
      DO 47 L=K,J-1
47    D2=D2+(R(LM(L))-RANK)**2
      NUM=1
      K=J
      TOTAL=IRNK
      T=DATA(LM(J),N)
      GO TO 46
48    NUM=NUM+1
      TOTAL=TOTAL+IRNK
46    CONTINUE
      RANK=TOTAL/NUM
      IF(NUM.GT.1) TOTY=TOTY+(NUM**3-NUM)/12.
      DO 49 L=K,NC
49    D2=D2+(R(LM(L))-RANK)**2
      SUMY2=(TC**3-TC)/12.-TOTY
      GO TO 3
55    D2=0
      SUMY2=SUMX2
3     C(N-NN+1)=(SUMX2+SUMY2-D2)/(2.*SQRT(SUMX2*SUMY2))
      IF(NN.EQ.1) WRITE(IOUT,50) NAMES(I),(C(J),J=1,NCOR)
50    FORMAT('0',A5,2X,15F8.4)
      IF(NN.NE.1) WRITE(IOUT,51)(C(J),J=1,NCOR)
51    FORMAT(8X,15F8.4)
2     CONTINUE
      WRITE(IOUT,62)
62    FORMAT(1X)
      DO 60 I=1,NV,ISQ
      NN=I+ISQ-1
      IF(NN.GT.NV) NN=NV
60    WRITE(IOUT,61) (NAMES(J),J=I,NN)
61    FORMAT(8X,15(2X,A5,1X))
      RETURN
      END
C                                             *** STAT PACK ***
C     SUBROUTINE FOR NAMEING OF VARIABLES
C     CALLING SEQUENCE: CALL STPNAM(NV,NAMES)
C     WHERE NV - IS NUMBER OF VARIABLES
C           NAMES - IS A VECTOR CONTAINING VARIABLE NAMES
C
C     SUBROUTINE ALLOWS FOR THE  NAMEING OF VARIABLES.
C     NAMES MUST BE 5 CHARACTERS OR LESS AND BEGIN WITH A LETTER
C     NAMES "HELP" AND "ALL" ARE NOT ALLOWED.
C
      SUBROUTINE STPNAM(NV,NAMES)
      COMMON /DEV/ICC,IDATA,IOUT,IDLG,IDSK
      COMMON /PRNT/ LINPP,ICOPS,RUNPRG
      DIMENSION NAMES(1),NAME(5)
      DO 1 I=1,NV
5     IF(ICC.NE.2) WRITE(IDLG,2) I
2     FORMAT('+VAR',I3,'? ',$)
      READ(ICC,3,END=99) NAME
3     FORMAT(5A1)
      IF(NAME(1).EQ.'!') RETURN
      IF((NAME(1).NE.'H').OR.(NAME(2).NE.'E').OR.(NAME(3).NE.'L')
     1.OR.(NAME(4).NE.'P').OR.(NAME(5).NE.' ')) GO TO 6
      WRITE(IDLG,4)
4     FORMAT('0NAME IS USED TO ASSIGN UP TO FIVE CHARACTER VARIABLE',
     1' NAMES'/' WHEN ASKED FOR THE VARIABLE NAME "VAR.---?", ',
     2'SIMPLY TYPE IN'/' A 5 CHARACTER OR LESS NAME SUBJECT TO THE',
     3' FOLLOWING RESTRICTIONS:'/' THE FIRST CHARACTER MUST ',
     4'BE A LETTER AND THE NAMES "ALL" AND'/' "HELP" MAY NOT BE USED.',
     5'  IF A BLANK IS USED THE MACHINE WILL'/' RETAIN ',
     6'THE PREVIOUS NAME'/)
      GO TO 5
6     DO 11 J=5,1,-1
      IF(NAME(J).NE.' ') GO TO 12
11    CONTINUE
      GO TO 1
12    IF((NAME(1).GE.'A').AND.(NAME(1).LE.'Z')) GOTO 8
      WRITE(IDLG,7)
7     FORMAT(' NAME MUST BEGIN WITH A LETTER'/)
      GO TO 5
8     DO 13 K=1,J
      IF(NAME(K).EQ.',') GO TO 14
      IF(NAME(K).EQ.'-') GO TO 14
      IF(NAME(K).EQ.';') GO TO 14
      IF(NAME(K).EQ.' ') GO TO 14
13    CONTINUE
      GO TO 10
14    WRITE(IDLG,15)NAME(K)
15    FORMAT('0THE CHARACTER "',A1,'" MAY NOT BE USED IN A NAME'/)
      GO TO 5
10    ENCODE(5,3,NAME1) NAME
      IF(NAME1.EQ.'ALL') GO TO 16
      IF(NAME1.EQ.'STOP') GO TO 16
      IF(NAME1.EQ.'EMPTY') GO TO 16
      DO 20 J=1,NV
      IF(NAME1.NE.NAMES(J)) GO TO 20
      WRITE(IDLG,21) J
21    FORMAT(' VARIABLE NUMBER ',I3,' ALREADY HAS THE NAME'/)
      GO TO 5
20    CONTINUE
      NAMES(I)=NAME1
      GO TO 1
16    WRITE(IDLG,17) NAME1
17    FORMAT(' THE NAME "',A5,'" IS AN ILLEGAL NAME'/)
      GO TO 5
1     CONTINUE
99    RETURN
      END
C                                                         **** STAT PACK ****
C
C     ROUTINE FOR SORTING DATA INTO ASCENDING ORDER.
C     CALLING SEQUENCE: CALL SORTCR(NV,NC,MV,MC,DATA,IV,SP,NAMES)
C     WHERE: NV - NUMBER OF VARIABLES USED.
C            NC - NUMBER OF OBSERVATIONS USED.
C            MV - MAXIMUM NUMBER OF VARIABLES POSSIBLE.
C            MC - MAXIMUM NUMBER OF OBSERVATIONS POSSIBLE.
C            DATA - MATRIX CONTAINING DATA (MC X MV)
C            IV - EXTRA VECTOR AT LEAST NC LONG
C            SP - EXTRA VECTOR AT LEAST NC LONG.
C            NAMES - VECTOR CONTAINING VARIABLE NAMES.
C
C     SORT MADE SEPARATE FROM MANIP TO DECREASE OVERALL SIZE
C     OF STP
C
      SUBROUTINE SORTCR(NV,NC,MV,MC,DATA,IV,SP,NAMES)
      DIMENSION DATA(MC,MV),IV(1),SP(1),IX(20),NAMES(1)
      COMMON /DEV/ICC,IDATA,IOUT,IDLG,IDSK
      COMMON /PRNT/ LINPP,ICOPS,RUNPRG
2     IF(ICC.NE.2) WRITE(IDLG,1)
1     FORMAT('0LIST SORT VARIABLES MAJOR TO MINOR? ',$)
      CALL ALPHA (IX,20,NN,IRET,IHELP,IERR,NAMES,NV)
      IF(IRET.EQ.1) RETURN
      IF(IERR.EQ.1) GO TO 2
      IF(IHELP.NE.1) GO TO 5
      WRITE(IDLG,4)
4     FORMAT(' THE STP SORT IS AN INCORE SORT USED TO SORT DATA IN'/
     1' ASCENDING ORDER.  UP TO 20 VARIABLSE MAY BE SPECIFIED BY THE'/
     2' USER AS SORT FIELDS.  OBSERVATIONS WILL REMAIN UNCHANGED AT'/
     3' THE END OF THE SORT, ONLY THE ORDER IN WHICH THEY OCCUR WILL'/
     4' CHANGE.  VARIABLES TO BE USED AS SORT FIELDS SHOULD BE LISTED'/
     5' FROM MAJOR TO MINOR.  THE MAJOR SORT FIELD WILL BE IN ORDER'/
     6' WHEN FINISHED, HOWEVER, THOSE CASES WHERE TWO OBSERVATIONS'/
     7' HAVE THE SAME VALUE FOR THE MAJOR SORT VARIABLE WILL BE'/
     8' ORDERED BY THE VALUES OF THE SECOND MAJOR SORT VARIABLE.'/
     9' IN THIS MANNER THE SORT PROCEEDS ON.   THE ONLY'/
     1' CASE WHERE THE MINOR SORT FIELD HAS ANY AFFECT, IS WHEN 2'/
     2' OBSERVATIONS AHAVE EQUAL VALUES FOR ALL OTHER SORT VARIABLES'/
     3' IF AFTER ALL SORT FIELDS HAVE BEEN TESTED, NO DECISION CAN '/
     4' BE MADE, THE OBSERVATION NUMBERS WILL BE USED.')
      GO TO 2
C
C     IF "ALL" IS SPECIFIED IN OTHER THAN FIRST POSITION, TREAT
C     IT AS "ALL THE REST", BUT ALLOW THOSE VARIABLES SPECIFICALLY
C     MENTIONED PRIOR TO THE "ALL" TO BE IN ORDER THEY WERE
C     SPECIFIED
C
5     DO 6 M=1,NN
      IF(IX(M).GT.0) GO TO 6
      IF(NV.LE.20) GO TO 8
      WRITE(IDLG,7)
7     FORMAT(' CANNOT SPECIFY MORE THAN 20 VARIABLES')
      GO TO 2
8     L=M-1
      K=M
      DO 9 I=1,NV
      DO 10 J=1,L
      IF(I.EQ.IX(J)) GO TO 9
10    CONTINUE
      IX(K)=I
      K=K+1
9     CONTINUE
      NN=NV
      GO TO 20
6     CONTINUE
C
C     CALL SORT
C
20    CALL SORT (NV,NC,MV,MC,DATA,IX,NN,IV,SP)
      IF(ICC.NE.2) WRITE(IDLG,21)(NAMES(IX(I)),I=1,NN)
21    FORMAT(' DATA SORTED BY VARIABLES: ',7(A5,',')/10(A5,','))
      RETURN
      END
C                               *** STAT PACK ***
C     SUBROUTINE PART OF SORTCR USED TO SORT DATA
C     CALLING SEQUENCE: CALL SORT(NV,NC,MV,MC,DATA,IS,KKL,IV,SP)
C     WHERE: NV - NUMBER OF VARIABLES
C            NC - NUMBER OF OBSERVATIONS
C            MV - MAXIMUM NUMBER OF VARIABLES POSSIBLE
C            MC - MAXIMUM NUMBER OF OBSERVATIONS POSSIBLE
C            DATA - MATRIX CONTAINING DATA
C            IS - IS A VECTOR CONTAINING SORT KEYS MUST BE AT LEAST
C                 KKL LONG
C            KKL - NUMBER OF SORT KEYS
C            IV - EXTRA VECTOR AT LEAST NV LONG
C            SP - EXTRA VECTOR AT LEAST NC LONG
C
C     SORTING METHOD USED IS PARTITIONING SORT ORIGINALLY OBTAINED
C     FROM ACM COMMUNICATIONS
C
      SUBROUTINE SORT(NV,NC,MV,MC,DATA,IS,KKL,IV,SP)
      DIMENSION DATA(MC,MV),IV(1),IS(1),IU(16),IL(16),SP(1)
      DIMENSION GIP(25)
      DO 1 I=1,NC
1     IV(I)=I
      M=1
      II=1
      J=NC
11    IF(II.GE.J) GO TO 18
12    K=II
      IJ=(J+II)/2
      I=0
31    I=I+1
      IF(I.GT.KKL) GO TO 33
      T1=DATA(IV(IJ),IS(I))
      T2=DATA(IV(II),IS(I))
      IF(T2.EQ.T1) GO TO 31
      IF(T2.LT.T1) GO TO 13
      GO TO 32
33    IF(IV(II).LE.IV(IJ)) GO TO 13
32    ISAV=IV(IJ)
      IV(IJ)=IV(II)
      IV(II)=ISAV
13    LL=J
      I=0
34    I=I+1
      IF(I.GT.KKL) GO TO 36
      T1=DATA(IV(IJ),IS(I))
      T2=DATA(IV(J),IS(I))
      IF(T2.EQ.T1) GO TO 34
      IF(T2.GT.T1) GO TO 5
      GO TO 35
36    IF(IV(J).GE.IV(IJ)) GO TO 5
35    ISAV=IV(IJ)
      IV(IJ)=IV(J)
      IV(J)=ISAV
      I=0
37    I=I+1
      IF(I.GT.KKL) GO TO 39
      T1=DATA(IV(IJ),IS(I))
      T2=DATA(IV(II),IS(I))
      IF(T2.EQ.T1) GOTO 37
      IF(T2.LT.T1) GO TO 5
      GO TO 38
39    IF(IV(II).LE.IV(IJ)) GO TO 5
38    ISAV=IV(IJ)
      IV(IJ)=IV(II)
      IV(II)=ISAV
      GO TO 5
5     DO 6 L=1,KKL
6     GIP(L)=DATA(IV(IJ),IS(L))
      NEXTRA=IV(IJ)
      GO TO 15
14    ISAV=IV(LL)
      IV(LL)=IV(K)
      IV(K)=ISAV
15    LL=LL-1
      I=0
40    I=I+1
      IF(I.GT.KKL) GO TO 41
      T1=GIP(I)
      T2=DATA(IV(LL),IS(I))
      IF(T2.EQ.T1) GO TO 40
      IF(T2.GT.T1) GO TO 15
      GO TO 16
41    IF(IV(LL).GT.NEXTRA) GO TO 15
16    K=K+1
      I=0
42    I=I+1
      IF(I.GT.KKL) GO TO 44
      T1=GIP(I)
      T2=DATA(IV(K),IS(I))
      IF(T2.EQ.T1) GO TO 42
      IF(T2.LT.T1) GO TO 16
      GO TO 43
44    IF(IV(K).LT.NEXTRA) GO TO 16
43    IF(K.LE.LL) GO TO 14
      IF((LL-II).LE.(J-K)) GO TO 17
      IL(M)=II
      IU(M)=LL
      II=K
      M=M+1
      GOTO 19
17    IL(M)=K
      IU(M)=J
      J=LL
      M=M+1
      GOTO 19
18    M=M-1
      IF(M.EQ.0) GO TO 90
      II=IL(M)
      J=IU(M)
19    IF((J-II).GE.11) GO TO 12
      IF(II.EQ.1) GO TO 11
C
C     BUBLE SORT PORTION (FASTER THAN PARTITION ONLY IF SUBSET
C     BEING LOOKED AT IS 11 OBSERVATIONS OR LESS)
C
106   FORMAT(1X,I4)
      II=II-1
20    II=II+1
      IF(II.EQ.J) GO TO 18
      I=0
      NEXTRA=IV(II+1)
45    I=I+1
      IF(I.GT.KKL) GO TO 47
      T1=DATA(NEXTRA,IS(I))
      T2=DATA(IV(II),IS(I))
      IF(T2.EQ.T1) GO TO 45
      IF(T2.LT.T1) GO TO 20
      GO TO 46
47    IF(IV(II).LE.NEXTRA) GO T O 20
46    K=II
21    IV(K+1)=IV(K)
      K=K-1
      I=0
48    I=I+1
      IF(I.GT.KKL) GO TO 50
      T1=DATA(NEXTRA,IS(I))
      T2=DATA(IV(K),IS(I))
      IF(T2.EQ.T1) GOTO 48
      IF(T1.LT.T2) GO TO 21
      GO TO 49
50    IF(NEXTRA.LT.IV(K)) GO TO 21
49    IV(K+1)=NEXTRA
      GO TO 20
C
C     END SORT NOW PLACE TAGS IN CORRECT ORDER
C
90    DO 91 J=1,NC
      IF(IV(J).EQ.0) GOTO 91
      IF(IV(J).EQ.J) GO TO 91
      DO 92 K=1,NV
92    SP(K)=DATA(J,K)
      M=J
      L=J
93    DO 94 K=1,NV
94    DATA(M,K)=DATA(IV(L),K)
      M=IV(L)
      IV(L)=0
      L=M
      IF(IV(L).NE.J) GO TO 93
      DO 96 K=1,NV
96    DATA(L,K)=SP(K)
      IV(L)=0
91    CONTINUE
      RETURN
      END
C                                              *** STAT PACK ***
C     SUBROUTINE TO OUTPUT MORE THAN 1 COPY TO THE LINE PRINTER.
C     CALLING SEQUENCE: CALL STCOPY
C
C     SUBROUTINE SIMPLY SETS A CONSTANT FOR OUTPUT WHEN PRINTS
C     IS CALLED.
C
      SUBROUTINE STCOPY
      COMMON/DEV/ICC,IDATA,IOUT,IDLG,IDSK
      COMMON /PRNT/ LINPP,ICOPS,RUNPRG
1     WRITE (IDLG,2)
2     FORMAT('0HOW MANY OUTPUT COPIES? ',$)
      READ(ICC,3)ICOPS
3     FORMAT(I)
      IF(ICOPS.LE.63) RETURN
      WRITE(IDLG,4)
4     FORMAT(' NOT MORE THAN 63 COPIES')
      GO TO 1
      END