Google
 

Trailing-Edge - PDP-10 Archives - decus_20tap5_198111 - decus/20-0137/manova/manov1.for
There is 1 other file named manov1.for in the archive. Click here to see a list.
00100	C	WMU IMPLEMENTATION OF THE WAYNE STATE VERSION OF
00200	C
00300	C	THE MIAMI MULTIVARIATE ANALYSIS OF VARIANCE(MANOVA)
00400	C
00500	C	MODIFIED FOR WMU COMPUTER CENTER BY:	SAM ANEMA
00600	C
00700	C	FILES: MANOVA.F4,MANOV1.F4,MANOV2.F4,MANOV3.F4,MANOV4.F4
00800	C
00900	C	ADDITIONAL SUBROUTINES: USAGE - RUSS BARR(IN NGLIB).MAC
01000	C
01100	C	LOADING PROCEDURE:
01200	C
01300	C	R LOADER
01400	C	MANOV=MANOVA,USAGE/<MANOV1/>MANOV2/>MANOV3/>MANOV4/>/G
01500	C	SAV MANOVA
01600	C
01700	C	MANOV1.F4
01800	      SUBROUTINE RDCONT
01900	      DIMENSION XKONE(14,14,5), NCDHLD(100,8),VECHLD(100),
02000	     1LIST(720),     OPTION(75), INLET(75), RECODE(8,20),SSRES(25,26)
02100	      DIMENSION XNUM(50)
02200	      COMMON ORTHES(100,20),DUMMY(26,27),SSHYP(25,26),SSEAD(40,41),
02300	     1ESTIM(50,50),        NUMERR(8),NERRS(100),NDFCUM(100),NDFTST(
02400	     2100),VARNAM(2,50),HEAD(3,100),LEVEL(8),LEVSUB(8,10),LEVCUM(8,10),
02500	     3NCELL(100), NTABLE(27),ITABLE(9,9),OBS(100)
02600	     4      ,NVAR,NCOVAR,NERWIT,NERRES,HNUM(100),ERROR,NTESTS,RVARC,
02700	     5           FIRST,IORD(50),IPOSV(50),IVPOS(50),NCELLS,NVART,NDFTOT,
02800	     6SPECOR,VLIST,PRINTR,NFACT , READK,  PRINTK,CONTR ,TESTR, MFIRST,
02900	     7TRUTH,BLANK,MAXFAC,MAXCEL,MAXPAR,MAXLEV,MAXVAR,ATITLE,APROB,
03000	     8AANALY,AFINIS,WITHIN,SPACE(10)
03100	      EQUIVALENCE (ESTIM,XKONE), (RECODE,DUMMY),(TAGP,SPACE(3)),
03200	     1(SSEAD,OPTION,INLET),   (LIST,SSEAD(876) ),(SSRES,SSHYP),
03300	     2(TAG,SPACE(1) ), (NCDHLD,SSEAD(76) ), (EONLY,SPACE(2) )
03400	     3,(PRTSSH,SPACE(5)),(NPCT,SPACE(6)),(NSIG,SPACE(7))
03500	      DATA INP,IOUT,IAUX1,IAUX2/5,30,22,23/
03600	      DATA IAUX4/24/
03700	      LOGICAL TRUTH, FIRST, ERROR,       RDNAM, TESTR, TRANS, PRINTK,
03800	     1PRINTR, SPECOR, CONTR, VLIST, RVARC ,IN4,  FIRSTT,READK ,MFIRST,
03900	     2PRINTM, EONLY,PRTSSH
04000	      DOUBLE PRECISION ABLANK,VARNAM,VECHLD
04100	      INTEGER RECODE,TRANCD,                  TWO24,OPTION,BLANK
04200	      DATA XNUM/' 1  ',' 2  ',' 3  ',' 4  ',' 5  ',' 6  ',' 7  ',' 8  ',
04300	     C   ' 9  ','10  ','11  ','12  ','13  ','14  ','15  ','16  ','17  ',
04400	     C   '18  ','19  ','20  ','21  ','22  ','23  ','24  ','25  ','26  ',
04500	     C   '27  ','28  ','29  ','30  ','31  ','32  ','33  ','34  ','35  ',
04600	     C   '36  ','37  ','38  ','39  ','40  ','41  ','42  ','43  ','44  ',
04700	     C   '45  ','46  ','47  ','48  ','49  ','50  '/
04800	      DATA  LR/4HR   /, ABLANK/8H         /
04900	      FIRSTT = .TRUE.
05000	      IF (FIRST) GO TO 100
05100	      IF (ERROR) GO TO 120
05200	      IF (.NOT. TRUTH) GO TO 150
05300	      GO TO 310
05400	  100 CALL ALPHAN(19,NTABLE(9),76H0   1   2   3   4   5   6   7   8   9
05500	     1   ,       +   .   =   -   W   G   M        )
05600	      CALL ALPHAN(4,ATITLE,16HTITLPROBANALFINI  )
05700	      CALL ALPHAN(1,BLANK,4H         )
05800	      DO 105 J=1,50
05900	  105 HNUM(J)=XNUM(J)
06000	  110 ITABLE(1,1) = 1
06100	      ITABLE(1,2) = 2
06200	      ITABLE(1,3) = 1
06300	      ITABLE(1,4) = 5
06400	      ITABLE(1,5) = 1
06500	      ITABLE(1,6) = 1
06600	      ITABLE(1,7) = 1
06700	      ITABLE(1,8) = 1
06800	      ITABLE(1,9) = 1
06900	      ITABLE(2,1) = 4
07000	      ITABLE(2,2) = 10
07100	      ITABLE(2,3) = 4
07200	      ITABLE(2,4) = 5
07300	      ITABLE(2,5) = 4
07400	      ITABLE(2,6) = 4
07500	      ITABLE(2,7) = 4
07600	      ITABLE(2,8) = 4
07700	      ITABLE(2,9) = 4
07800	      ITABLE(3,1) = 3
07900	      ITABLE(3,2) = 5
08000	      ITABLE(3,3) = 5
08100	      ITABLE(3,4) = 3
08200	      ITABLE(3,5) = 5
08300	      ITABLE(3,6) = 5
08400	      ITABLE(3,7) = 5
08500	      ITABLE(3,8) = 5
08600	      ITABLE(3,9) = 3
08700	      ITABLE(4,1) = 4
08800	      ITABLE(4,2) = 5
08900	      ITABLE(4,3) = 5
09000	      ITABLE(4,4) = 4
09100	      ITABLE(4,5) = 5
09200	      ITABLE(4,6) = 5
09300	      ITABLE(4,7) = 11
09400	      ITABLE(4,8) = 11
09500	      ITABLE(4,9) = 4
09600	      ITABLE(5,1) = 6
09700	      ITABLE(5,2) = 5
09800	      ITABLE(5,3) = 5
09900	      ITABLE(5,4) = 5
10000	      ITABLE(5,5) = 5
10100	      ITABLE(5,6) = 5
10200	      ITABLE(5,7) = 5
10300	      ITABLE(5,8) = 5
10400	      ITABLE(5,9) = 5
10500	      ITABLE(7,1) = 5
10600	      ITABLE(7,2) = 7
10700	      ITABLE(7,3) = 5
10800	      ITABLE(7,4) = 5
10900	      ITABLE(7,5) = 5
11000	      ITABLE(7,6) = 5
11100	      ITABLE(7,7) = 5
11200	      ITABLE(7,8) = 5
11300	      ITABLE(7,9) = 5
11400	      ITABLE(8,1) = 5
11500	      ITABLE(8,2) = 8
11600	      ITABLE(8,3) = 5
11700	      ITABLE(8,4) = 5
11800	      ITABLE(8,5) = 5
11900	      ITABLE(8,6) = 5
12000	      ITABLE(8,7) = 5
12100	      ITABLE(8,8) = 5
12200	      ITABLE(8,9) = 5
12300	      ITABLE(9,1) = 9
12400	      ITABLE(9,2) = 5
12500	      ITABLE(9,3) = 9
12600	      ITABLE(9,4) = 5
12700	      ITABLE(9,5) = 5
12800	      ITABLE(9,6) = 9
12900	      ITABLE(9,7) = 5
13000	      ITABLE(9,8) = 9
13100	      ITABLE(9,9) = 5
13200	      NPCT = 0
13300	  120 FIRSTT = .TRUE.
13400	      IF (.NOT. ERROR) WRITE         (IOUT,125)
13500	  125 FORMAT(1H1//' WESTERN MICHIGAN UNIVERSITY'//,
13600	     1	 ' MULTI-VARIATE ANALYSIS OF VARIANCE'/)
13700	  130 READ(INP,140)TAG,TAGP,OPTION
13800	  140 FORMAT (A4,76A1)
13900	  150 IF (TAG .NE. ATITLE .AND. TAG .NE. APROB) GO TO 560
14000	      IF       (ERROR) WRITE         (IOUT,125)
14100	      FIRST = .TRUE.
14200	      EONLY = .FALSE.
14300	      NSIG = 1
14400	  170 IF (TAG .EQ. APROB) GO TO 190
14500	      IF (TAG .NE. ATITLE) GO TO 560
14600	      WRITE(IOUT,180)OPTION
14700	  180 FORMAT (1H 80A1)
14800	      READ(INP,140) TAG,TAGP,OPTION
14900	      GO TO 170
15000	  190 ERROR = .FALSE.
15100	      DO 200 I = 1,8
15200	  200 NTABLE(I) = BLANK + 1
15300	      TWO24 = 16777216
15400	      DO 210 J = 1,4
15500	      IF (OPTION(J) .NE. BLANK) GO TO 205
15600	      OPTION(J)=0
15700	      GO TO 210
15800	  205 OPTION(J)=(OPTION(J)/"4000000000)-"60
15900	  206 FORMAT(I1)
16000	  210 CONTINUE
16100	      NVART = OPTION(1)*10 + OPTION(2)
16200	      NFACT = OPTION(3)
16300	      NFMT = OPTION(4)
16400	      RDNAM  = OPTION( 5) .NE. BLANK
16500	      VLIST  = OPTION( 6) .NE. BLANK
16600	      SPECOR = OPTION( 7) .NE. BLANK
16700	      TESTR  = OPTION( 8) .NE. BLANK
16800	      RVARC  = OPTION( 9) .NE. BLANK
16900	C     READ DATA FROM TAPE 4
17000	      IN4    = OPTION(10) .NE. BLANK
17100	      IF (OPTION(10) .EQ. LR) REWIND IAUX4
17200	      PRINTM = OPTION(11) .EQ. BLANK
17300	      TRANS  = OPTION(12) .NE. BLANK
17400	      PRINTK = OPTION(13) .NE. BLANK
17500	      PRINTR = OPTION(14) .NE. BLANK
17600	      READK  = OPTION(15) .NE. BLANK
17700	      MFIRST = OPTION(16) .EQ. BLANK
17800	      PRTSSH = OPTION(17) .NE. BLANK
17900	      IF (OPTION(18) .EQ. BLANK) GO TO 215
18000	      REWIND IAUX4
18100	  212 READ(INP,140)TAG,TAGP,OPTION
18200	      WRITE(IAUX4,140) TAG,TAGP,OPTION
18300	      IF (TAG .NE. AFINIS) GO TO 212
18400	      REWIND IAUX4
18500	  215 NPCT = (NPCT/1000)+1
18600	      WRITE(IOUT,220)  NPCT,NVART,NFACT
18700	  220 FORMAT (8H0PROBLEM I4,I4, 10H VARIABLES I3, 8H FACTORS /)
18800	      NPCT  = NPCT*1000
18900	      IF (NVART .GT. MAXVAR .OR. NFACT .GT. MAXFAC) GO TO 540
19000	      DO 230 I = 1,NVART
19100	      IORD(I) = I
19200	      IPOSV(I) = I
19300	  230 IVPOS(I) = I
19400	      IF (RDNAM) GO TO 250
19500	      DO 240 J = 1,NVART
19600	      VARNAM(1,J) = ABLANK
19700	  240 VARNAM(2,J) = HNUM(J)
19800	      GO TO 270
19900	  250 READ(INP,260)((VARNAM(J,I),J = 1,2),I = 1,NVART)
20000	  260 FORMAT (8(A6,A4))
20100	  270 NVAR = NVART
20200	      IF(RDNAM)WRITE(IOUT,300)(VARNAM(1,J), VARNAM(2,J),J = 1,NVART)
20300	      NCOVAR = 0
20400	      IF (VLIST)READ(INP,280)NVAR,NCOVAR,(IORD(J),J = 1,38)
20500	  280 FORMAT (40I2)
20600	      NVCVAR = NVAR+NCOVAR
20700	      IF(VLIST.AND.NVCVAR.GT. 38) READ(INP,280)(IORD(J), J = 39,NVCVAR)
20800	      IF (NCOVAR .EQ. 0) NCOVAR = 0
20900	      WRITE(IOUT,290)NVAR,NCOVAR
21000	  290 FORMAT (1H0I2,9H CRITERIA I4, 43H COVARIATES  WITH THE FOLLOWING V
21100	     1ARIABLES      )
21200	      IF (NVCVAR .GT. NVART) GO TO 540
21300	      N2 = 2*NVCVAR
21400	      DO 295 J = 1,NVCVAR
21500	      J1 = IORD(J)
21600	      VECHLD(2*J-1) = VARNAM(1,J1)
21700	  295 VECHLD(2*J) = VARNAM(2,J1)
21800	      WRITE(IOUT,300)(VECHLD(J), J = 1,N2)
21900	  300 FORMAT (10(3X A6,A4) )
22000	  310 NVARTP = NVART+1
22100	      REWIND IAUX1
22200	      IF (FIRST) GO TO 330
22300	      DO 320 I = 1,NVART
22400	      IPOSV(I) = I
22500	  320 IVPOS(I) = I
22600	      IF (.NOT. CONTR) GO TO 430
22700	      READ(IAUX1)(((XKONE(I,J,K),I = 1,MAXLEV),J = 1,MAXLEV),K= 1,NFACT)
22800	      READ(IAUX1)   ((NCDHLD(J,K),J = 1,NCELLS),K = 1,NFACT)
22900	      READ(IAUX1)((ORTHES(J,K),J = 1,NCELLS), K = 1,NVART),((SSRES(J,K),
23000	     1J = 1,NVART),K = 1,NVARTP), ((VARNAM(J,K),J = 1,2),K = 1,NVART)
23100	      REWIND IAUX1
23200	  330 DO 420 I = 1,NFACT
23300	      READ(INP,340) NTABLE(I), ICONT,    M, (LEVSUB(I,J),J = 2,9),
23400	     1(RECODE(I,J), J = 1,20    ), (INLET(J),J = 1,5)
23500	  340 FORMAT (A1,I1,29I2,5A4)
23600	      NOREC=' '
23700	      IF (NOREC .NE. BLANK )  GO TO 360
23800	      DO 350 J = 1,M
23900	  350 RECODE(I,J) = J
24000	  360 LEVEL(I) = M
24100	      WRITE(IOUT,370)NTABLE(I),LEVEL(I), (INLET(J),J = 1,5)
24200	  370 FORMAT(8H0FACTOR A1,I5,7H LEVELS  10X,5A4)
24300	      IF (LEVEL(I) .EQ. 0 .OR. LEVEL(I) .GT. MAXLEV) GO TO 540
24400	      LEVSUB(I,1)    = LEVEL(I) - 1
24500	      LEVCUM(I,1)    = 1
24600	      IF (LEVSUB(I,2) .EQ. 0) GO TO 410
24700	      LEVCUM(I,2) = 1
24800	      DO 390 J = 2,9
24900	      IF (LEVSUB(I,J) .EQ. 0) GO TO 400
25000	      J1 = J-1
25100	      WRITE(IOUT,380)NTABLE(I),J1,LEVSUB(I,J )
25200	  380  FORMAT (1H 9X, A1,I1,4H HAS I3,3H DF)
25300	  390 LEVCUM(I,J+1 ) = LEVCUM(I,J)  + LEVSUB(I,J)
25400	      J = 10
25500	  400 IF (LEVEL(I) .NE. LEVCUM(I,J))GO TO 520
25600	  410 CALL GETK1(ICONT,I,LEVEL(I),XKONE,LEVSUB,NTABLE,ERROR)
25700	      IF (ERROR) GO TO 120
25800	  420 CONTINUE
25900	      WRITE(IAUX1)(((XKONE(I,J,K),I=1,MAXLEV),J= 1,MAXLEV), K = 1,NFACT)
26000	  430 IF (.NOT. SPECOR) GO TO 460
26100	      DO 450 N1 = 1,641,80
26200	      N2 = N1 + 79
26300	      READ(INP,440) (LIST(J),J = N1,N2)
26400	  440 FORMAT (80A1)
26500	      DO 450 J = N1,N2
26600	      IF (LIST(J) .EQ.'.' ) GO TO 460
26700	  450 CONTINUE
26800	      GO TO 510
26900	  460 IF (.NOT. FIRST .AND. .NOT. CONTR) GO TO 480
27000	      IF (.NOT. FIRST) GO TO 470
27100	      CALL RDDAT(NOBST,NFMT,TRANS,IN4,PRINTM)
27200	      NERWIT = NOBST-NCELLS
27300	      IF (ERROR) GO TO 120
27400	  470 WRITE(IAUX1)  ((NCDHLD(J,K),J = 1,NCELLS),K = 1,NFACT)
27500	      WRITE(IAUX1)((ORTHES(J,K),J = 1,NCELLS),K = 1,NVART),((SSRES(J,K),
27600	     1J = 1,NVART),K = 1,NVARTP), ((VARNAM(J,K),J = 1,2),K = 1,NVART)
27700	      REWIND IAUX1
27800	  480 DO 500 J = 1,MAXLEV
27900	      DO 490 K = 1,MAXLEV
28000	  490 XKONE(J,K,MAXFAC+1) = 0.0
28100	  500 XKONE(J,J,MAXFAC+1) = 1.0
28200	      READ(IAUX1)(((XKONE(I,J,K),I=1,MAXLEV),J = 1,MAXLEV), K = 1,NFACT)
28300	510   RETURN
28400	  520 WRITE(IOUT,530)
28500	  530 FORMAT (25H SUBEFFECTS ABOVE WRONG    )
28600	      GO TO 120
28700	  540 WRITE(IOUT,550)
28800	  550 FORMAT (13H ERROR ABOVE   )
28900	      GO TO 120
29000	  560 IF (TAG .EQ. AFINIS) GO TO 580
29100	      IF (.NOT. FIRSTT)GO TO 130
29200	      WRITE(IOUT,570)TAG,OPTION
29300	  570 FORMAT (1H A4,76A1/40H ABOVE NOT VALID HERE, SKIPPING PROBLEM   )
29400	      FIRSTT = .FALSE.
29500	      GO TO 130
29600	  580 WRITE(IOUT,590)
29700	  590 FORMAT  (12H1END OF DATA  )
29800	      STOP
29900	      END
30000	      SUBROUTINE GETK1(ICONT,N1,NLEV,XKONE,LEVSUB,NTABLE,ERROR)
30100	      LOGICAL ERROR
30200	      DIMENSION  XKONE(14,14,5), LEVSUB(8,10), NTABLE(27)
30300	      DATA INP,IOUT,IAUX1,IAUX2/5,30,22,23/
30400	      DATA IAUX4/24/
30500	      M = ICONT+1
30600	      NBLANK = NTABLE(20)
30700	      NSUM = 2
30800	      K=0
30900	      IF (M .GT. 4) GO TO 350
31000	      GO TO (100,150,220,280),M
31100	  100 IF (NLEV .EQ. 1) GO TO 140
31200	      WRITE (IOUT,110)
31300	  110 FORMAT (23H   DEVIATION CONTRASTS   )
31400	      DO 130 I = 2,NLEV
31500	      XKONE(I-1,1,N1) = 1.0
31600	      XKONE(NLEV,I,N1) = -1.0
31700	      DO 130 J = 2,NLEV
31800	      IF (I .NE. J) GO TO 120
31900	      XKONE(I-1,J,N1) = 1.0
32000	      GO TO 130
32100	  120 XKONE(I-1,J,N1) = 0.0
32200	  130 CONTINUE
32300	  140 XKONE(NLEV,1,N1) = 1.0
32400	      GO TO 340
32500	  150 SUM = 0.0
32600	      WRITE (IOUT,160)
32700	  160 FORMAT (20H   SPECIAL CONTRASTS )
32800	      I = 0
32900	  165 I = I+1
33000	      READ (INP,170)(XKONE(I,J,N1),J=1,NLEV)
33100	  170 FORMAT (16F5.0)
33200	      GO TO 330
33300	  180 WRITE (IOUT,190)LETTER,NUMB,(XKONE(I,J,N1),J=1,NLEV)
33400	  190 FORMAT (3X,2A1,10F12.5/(5X,10F12.5) )
33500	  200 SUM = SUM + XKONE(1,I,N1)
33600	      IF (I .LT. NLEV) GO TO 165
33700	      DO 210 I = 1,NLEV
33800	  210 XKONE(1,I,N1) = XKONE(1,I,N1)/SUM
33900	      CALL INVERT(NLEV,XKONE(1,1,N1),E    ,10HCONTRAST    )
34000	      GO TO 340
34100	C
34200	  220 READ (INP,170)(XKONE(I,2,N1),I=1,NLEV)
34300	      WRITE (IOUT,230)(XKONE(I,2,N1),I=1,NLEV)
34400	  230 FORMAT (45H   ORTHOGONAL POLYNOMIAL CONTRASTS WITH SCALE /
34500	     1(5X 10F12.5) )
34600	      WRITE (IOUT,190)
34700	      A = .5*(XKONE(NLEV,2,N1) + XKONE(1,2,N1) )
34800	      B = 2./(XKONE(NLEV,2,N1) - XKONE(1,2,N1) )
34900	      DO 240 I = 1,NLEV
35000	      XKONE(I,2,N1) = B*(XKONE(I,2,N1)-A)
35100	  240 XKONE(I,1,N1) = 1.0
35200	      DO 250 I = 3,NLEV
35300	      DO 250 J = 1,NLEV
35400	  250 XKONE(J,I,N1) = 2.0*XKONE(J,I-1,N1)*XKONE(J,2,N1)- XKONE(J,I-2,N1)
35500	      CALL GRAM(XKONE(1,1,N1),NLEV,NLEV)
35600	      DO 260 I = 1,NLEV
35700	  260 XKONE(I,1,N1) = 1.0
35800	      I = 0
35900	  265 I = I+1
36000	      GO TO 330
36100	  270 WRITE (IOUT,190)LETTER,NUMB,(XKONE(J,I,N1),J=1,NLEV)
36200	      IF (I .LT. NLEV) GO TO 265
36300	      GO TO 340
36400	  280 XKONE(1,1,N1) = 1.0
36500	      WRITE (IOUT,290)
36600	  290 FORMAT (25H   DIFFERENCE CONTRASTS   )
36700	      DO 320 K = 2,NLEV
36800	      RK = -1.0/FLOAT(K)
36900	      XKONE(K,1,N1) = 1.0
37000	      DO 300 J = K,NLEV
37100	  300 XKONE(J,K,N1) = 0.0
37200	      DO 310 J = 1,K
37300	  310 XKONE(J,K,N1) = RK
37400	  320 XKONE(K,K,N1) = 1.0 + RK
37500	      GO TO 340
37600	  330 LETTER = NBLANK
37700	      NUMB = NBLANK
37800	      IF (I .NE. NSUM) GO TO (180,270),ICONT
37900	      LETTER = NTABLE(N1)
38000	      NUMB = NTABLE(K+10)
38100	      K = K+1
38200	      NSUM = NSUM + LEVSUB(N1,K+1)
38300	      GO TO (180,270), ICONT
38400	340   RETURN
38500	  350 ERROR = .TRUE.
38600	      GO TO 340
38700	      END
38800	      SUBROUTINE RDDAT(NOBST,NFMT,TRANS,IN4,PRINTM)
38900	      DIMENSION TRANCD(50), RECODE(8,20),
39000	     1FMT(180) ,X(50),NCODE(8) ,SSWITH(25,25),CELLM(100,20),
39100	     2CLLSS(50,20),NCDHLD(100,8),SSRES(25,26)
39200	      COMMON ORTHES(100,20),DUMMY(26,27),SSHYP(25,26),SSEAD(40,41),
39300	     1ESTIM(50,50),            NUMERR(8),NERRS(100),NDFCUM(100),NDFTST(
39400	     2100),VARNAM(2,50),HEAD(3,100),LEVEL(8),LEVSUB(8,10),LEVCUM(8,10),
39500	     3NCELL(100), NTABLE(27),ITABLE(9,9),OBS(100)
39600	     4      ,NVAR,NCOVAR,NERWIT,NERRES,HNUM(100),ERROR,NTESTS,RVARC,
39700	     5           FIRST,IORD(50),IPOSV(50),IVPOS(50),NCELLS,NVART,NDFTOT,
39800	     6SPECOR,VLIST,PRINTR,NFACT , READK,  PRINTK,CONTR ,TESTR, MFIRST,
39900	     7TRUTH,BLANK,MAXFAC,MAXCEL,MAXPAR,MAXLEV,MAXVAR,ATITLE,AJOBCD,
40000	     8AANALY,AFINIS,WITHIN,SPACE(10)
40100	      EQUIVALENCE (ESTIM,CLLSS), (ORTHES,CELLM,DESIGN),
40200	     1(SSWITH,SSRES(1,2) ),  (RECODE,DUMMY),(NCDHLD,SSEAD(76) ),
40300	     2(TRANCD,SSEAD), (SSRES,SSHYP)
40400	      DATA INP,IOUT,IAUX1,IAUX2/5,30,22,23/
40500	      DATA IAUX4/24/
40600	      INTEGER RECODE,TRANCD
40700	      DOUBLE PRECISION VARNAM
40800	      LOGICAL ERROR,END,FIRSTC,TRANS,IN4,PRINTM
40900	C
41000	      FIRSTC = .TRUE.
41100	      NCELLS = 0
41200	      KOBS = 0
41300	      DO 100 J = 1,NVART
41400	      DO 100 K = J,NVART
41500	      SSRES(K,J) = 0.0
41600	  100 SSWITH(J,K) = 0.0
41700	      IF (NFMT .NE. 0) GO TO 120
41800	      N = 10-NFACT
41900	      ENCODE(37,601,FMT) NFACT,N
42000	601   FORMAT('(2X,',I4,'I1,',I4,'X,10F6.0/(12X,10F6.0))')
42100	      GO TO 160
42200	  120 NFMT20 = NFMT*20
42300	      WRITE(IOUT,130)
42400	  130 FORMAT   ( 21H0FORMAT OF DATA CARDS )
42500	      READ(INP,140) (FMT(J), J = 1,NFMT20)
42600	  140 FORMAT (20A4)
42700	      WRITE(IOUT,150)(FMT(J), J = 1,NFMT20)
42800	  150 FORMAT (4X 20A4)
42900	  160 IF (.NOT. TRANS) GO TO 280
43000	      READ(INP,170) (TRANCD(I),I = 1,NVART)
43100	  170 FORMAT (80I1)
43200	      WRITE(IOUT,180)
43300	  180 FORMAT (17H0TRANSFORMATIONS    )
43400	      DO 270 I = 1,NVART
43500	      M = TRANCD(I) + 1
43600	      IF (M .GT. 4) GO TO 250
43700	      GO TO (270,190,210,230), M
43800	  190 WRITE(IOUT,200)VARNAM(1,I), VARNAM(2,I)
43900	  200 FORMAT (16H SQRT  TRANS VAR   2X A6,A4)
44000	      GO TO 270
44100	  210 WRITE(IOUT,220) VARNAM(1,I), VARNAM(2,I)
44200	  220 FORMAT (16H LOGE  TRANS VAR   2X,A6,A4)
44300	      GO TO 270
44400	  230 WRITE(IOUT,240) VARNAM(1,I), VARNAM(2,I)
44500	  240 FORMAT (16H ASIN  TRANS VAR   2X A6,A4)
44600	      GO TO 270
44700	  250 WRITE(IOUT,260) TRANCD(  I), VARNAM(1,I), VARNAM(2,I)
44800	  260 FORMAT(5H USER I1,10H TRANS VAR  2X A6,A4)
44900	  270 CONTINUE
45000	  280 I = 0
45100	  290 IF (IN4)   GO TO 300
45200	      READ (INP,FMT)   (NCODE(K),K = 1,NFACT), (X(K), K = 1,NVART)
45300	      GO TO 305
45400	  300 READ (IAUX4,FMT) (NCODE(K),K = 1,NFACT), (X(K), K = 1,NVART)
45500	  305 DO 307 K = 1,NVART
45600	      IF (X(K) .NE. 0.) GO TO 311
45700	  307 CONTINUE
45800	      DO 310 K = 1,NFACT
45900	      IF (NCODE(K) .NE. 0)  GO TO 311
46000	  310 CONTINUE
46100	      GO TO 430
46200	  311 KOBS = KOBS + 1
46300	      ISUM = 0
46400	      DO 320 J = 1,NFACT
46500	      KLEV = LEVEL(J)
46600	      DO 318 K = 1,KLEV
46700	      IF (NCODE(J) .EQ. RECODE(J,K) ) GO TO 321
46800	  318 CONTINUE
46900	  322 WRITE(IOUT,319)KOBS,  (NCODE(NUMFAC), NUMFAC=1,NFACT)
47000	  319 FORMAT  (20X,11HOBSERVATION I10,47H OMITTED FROM ANALYSIS BECAUSE
47100	     1OF BAD CELL CODE        8I3)
47200	      GO TO 290
47300	  321 ISUM0 = ISUM
47400	      ISUM = ISUM*LEVEL(J) + K - 1
47500	      IF (ISUM .LT. ISUM0) GO TO 322
47600	  320 CONTINUE
47700	      IF (FIRSTC) GO TO 325
47800	      DO 323 I = 1,NCELLS
47900	      IF (ISUM .LT. NCELL(I)) GO TO 326
48000	      IF (ISUM .EQ. NCELL(I) ) GO TO 350
48100	  323 CONTINUE
48200	  325 I = NCELLS + 1
48300	      GO TO 329
48400	  326 IF (NCELLS .EQ. MAXCEL) GO TO 380
48500	      DO 328 J1 = I,NCELLS
48600	      J = NCELLS - J1 + I
48700	      JP = J+1
48800	      OBS(JP) = OBS(J)
48900	      NCELL(JP) = NCELL(J)
49000	      DO 327 K = 1,NVART
49100	      CELLM(JP,K) = CELLM(J,K)
49200	  327 CLLSS(JP,K) = CLLSS(J,K)
49300	      DO 328 K = 1,NFACT
49400	  328 NCDHLD(JP,K) = NCDHLD(J,K)
49500	  329 FIRSTC = .FALSE.
49600	      IF (NCELLS .EQ. MAXCEL) GO TO 380
49700	      NCELLS = NCELLS + 1
49800	      NCELL(I) = ISUM
49900	      DO 330 J = 1,NVART
50000	      CLLSS(I,J) = 0.0
50100	  330 CELLM(I,J) = 0.0
50200	      OBS(I) = 0.0
50300	      DO 340 J = 1,NFACT
50400	  340 NCDHLD(I,J) = NCODE(J)
50500	  350 OBS(I) = OBS(I)+1.
50600	      IF (TRANS) CALL TRANSF(NVART,X,TRANCD)
50700	      DO 360 K = 1,NVART
50800	      CELLM(I,K)=CELLM(I,K)+X(K)
50900	      IF (OBS(I) .EQ. 1.0) GO TO 360
51000	      DIFF1 = OBS(I)*X(K)-CELLM(I,K)
51100	      DIFF2 = DIFF1/(OBS(I)*(OBS(I)-1.))
51200	      CLLSS(I,K) = CLLSS(I,K) + DIFF1*DIFF2
51300	      DO 355 L = 1,K
51400	  355 SSWITH(L,K) =SSWITH(L,K) + (OBS(I)*X(L)-CELLM(I,L))*DIFF2
51500	  360 CONTINUE
51600	      GO TO 290
51700	  380 WRITE(IOUT,400) KOBS
51800	  400 FORMAT (46H0ERROR IN DATA OR TOO MANY CELLS  OBSERVATION  I10)
51900	      ERROR = .TRUE.
52000	  430 WRITE(IOUT,432)NCELLS
52100	  432 FORMAT (1H0,I3,6H CELLS  )
52200	      IF (NCELLS .NE. 0) GO TO 434
52300	      ERROR = .TRUE.
52400	      GO TO 530
52500	  434 OBST = 0.0
52600	      DO 450 I = 1,NCELLS
52700	      OBST = OBST + OBS(I)
52800	      DO 450 J = 1,NVART
52900	      IF (OBS(I) .EQ. 1.0) GO TO 450
53000	      CLLSS(I,J) = SQRT(CLLSS(I,J)/(OBS(I)-1.))
53100	  450 CELLM(I,J) = CELLM(I,J)/OBS(I)
53200	      NOBST = OBST
53300	      IF (.NOT. PRINTM) GO TO 530
53400	      ENCODE(23,540,FMT) NFACT
53500	540   FORMAT('(1H ,',I4,'I3,I10,4H OBS)')
53600	      N2 = 0
53700	      IF (NCELLS .LE. 8) WRITE (IOUT,470)
53800	  470 FORMAT (1H1 50X 30HMEANS AND STANDARD DEVIATIONS )
53900	  480 N1 = N2+1
54000	      N2 = N1+7
54100	      IF (N2 .GT. NVART) N2 = NVART
54200	      IF (NCELLS .GT. 8) WRITE (IOUT,470)
54300	      WRITE(IOUT,490) (NTABLE(J), J = 1,MAXFAC),(VARNAM(1,J),
54400	     1VARNAM(2,J),J = N1,N2)
54500	  490 FORMAT(1H0/ 7H/FACTOR 48X 8HVARIABLE/ 1H 4(2X A1), 23X 8(2XA6,A4))
54600	      DO 510 I = 1,NCELLS
54700	      NOBS = OBS(I)
54800	      WRITE(IOUT,FMT)  (NCDHLD(I,J), J = 1,NFACT),NOBS
54900	      WRITE(IOUT,500)          (CELLM(I,J),J = N1,N2)
55000	  500 FORMAT (1H 28X,2H M 5X 8F12.3)
55100	  510 IF (NOBS .GT. 1) WRITE(IOUT,520)(CLLSS(I,J), J = N1,N2)
55200	  520 FORMAT (1H 28X 2HSD 5X 8F12.3)
55300	      IF (N2 .NE. NVART) GO TO 480
55400	530   RETURN
55500	      END
55600	      SUBROUTINE INVERT(N,A,ERROR,HEAD)
55700	      DIMENSION A(14,14), NPERM(19),NSEQ(20),HEAD(3)
55800	      DATA INP,IOUT,IAUX1,IAUX2/5,30,22,23/
55900	      DATA IAUX4/24/
56000	      LOGICAL ERROR
56100	      ERROR = .FALSE.
56200	      DO 100 K=1,N
56300	  100 NSEQ(K)=K
56400	      AMAX=0.
56500	      K1=1
56600	  110 K=K1
56700	      K1=K1+1
56800	      IF(N-K) 210,210,120
56900	  120 TEMP=ABS (A(K,K))
57000	      IPIV=K
57100	      DO 140 I=K1,N
57200	      IF(TEMP-ABS (A(I,K))) 130,140,140
57300	  130 IPIV=I
57400	      TEMP=ABS (A(I,K))
57500	  140 CONTINUE
57600	      IF(IPIV-K) 170,170,150
57700	  150 DO 160 J=1,N
57800	      TEMP=A(IPIV,J)
57900	      A(IPIV,J)=A(K,J)
58000	  160 A(K,J)=TEMP
58100	      ITEMP=NSEQ(IPIV)
58200	      NSEQ(IPIV)=NSEQ(K)
58300	      NSEQ(K)=ITEMP
58400	  170 IF(AMAX-ABS (A(IPIV,IPIV))) 180,190,190
58500	  180 AMAX=ABS (A(IPIV,IPIV))
58600	  190 NPERM(K)=IPIV
58700	      DO 200 I=K1,N
58800	      TEMP=A(I,K)/A(K,K)
58900	      A(I,K)=TEMP
59000	      DO 200 J=K1,N
59100	  200 A(I,J)=A(I,J)-TEMP*A(K,J)
59200	      K=K+1
59300	      GO TO 110
59400	  210 DO 290 K=1,N
59500	      J=0
59600	  220 J=J+1
59700	      IF(K-J) 260,260,230
59800	  230 TEMP=A(K,J)
59900	      I=K
60000	  240 I=I+1
60100	      IF(N-I) 220,250,250
60200	  250 A(I,J)=A(I,J)-TEMP*A(I,K)
60300	      GO TO 240
60400	  260 I=K
60500	  270 I=I+1
60600	      IF(N-I) 290,280,280
60700	  280 A(I,K)=-A(I,K)
60800	      GO TO 270
60900	  290 CONTINUE
61000	      K=N
61100	  300 AKK=A(K,K)
61200	      IF(ABS (AKK/AMAX)-.01) 520,520,310
61300	  310 AKK=1./AKK
61400	      A(K,K)=AKK
61500	      J=K
61600	  320 J=J+1
61700	      IF(N-J) 360,330,330
61800	  330 TEMP=AKK*A(K,J)
61900	      A(K,J)=TEMP
62000	      I=0
62100	  340 I=I+1
62200	      IF(K-I) 320,320,350
62300	  350 A(I,J)=A(I,J)-TEMP*A(I,K)
62400	      GO TO 340
62500	  360 I=0
62600	  370 I=I+1
62700	      IF(K-I) 390,390,380
62800	  380 A(I,K)=-AKK*A(I,K)
62900	      GO TO 370
63000	  390 K=K-1
63100	      IF(K) 400,400,300
63200	  400 K1=1
63300	  410 K=K1
63400	      K1=K1+1
63500	      IF(N-K) 470,470,420
63600	  420 AKK=A(K,K)
63700	      DO 460 J=K1,N
63800	      AKK=AKK+A(K,J)*A(J,K)
63900	      SUMR=A(K,J)
64000	      SUMC=A(J,J)*A(J,K)
64100	      J1=J+1
64200	      IF(N-J1) 450,430,430
64300	  430 DO 440 I=J1,N
64400	      SUMR=SUMR+A(K,I)*A(I,J)
64500	  440 SUMC=SUMC+A(J,I)*A(I,K)
64600	  450 A(K,J)=SUMR
64700	  460 A(J,K)=SUMC
64800	      A(K,K)=AKK
64900	      GO TO 410
65000	  470 J=N
65100	  480 J=J-1
65200	      IF(J) 540,540,490
65300	  490 JPIV=NPERM(J)
65400	      IF(J-JPIV) 500,480,500
65500	  500 DO 510 I=1,N
65600	      TEMP=A(I,J)
65700	      A(I,J)=A(I,JPIV)
65800	  510 A(I,JPIV)=TEMP
65900	      GO TO 480
66000	  520 ERROR = .TRUE.
66100	      WRITE (IOUT,530)HEAD,(NSEQ(L),L=1,K)
66200	  530 FORMAT (14H0IN THE ABOVE 2A4,A2,13H MATRIX, ROW I3, 48H IS PROBABL
66300	     1Y A LINEAR COMBINATION OF ROWS        24I3)
66400	540   RETURN
66500	      END
66600	      SUBROUTINE GRAM(XMAT,NROW,NCOLX)
66700	      DIMENSION XMAT(14,14),CON(20)
66800	      N1=NCOLX-1
66900	      DO 150 K=1,NCOLX
67000	      SUM=0.
67100	      DO 100 I=K,NROW
67200	  100 SUM=SUM+XMAT(I,K)**2
67300	      TEMP=SUM
67400	      SUM=SQRT(SUM)
67500	      AKK=XMAT(K,K)
67600	      IF (AKK .LT. 0.0) SUM = -SUM
67700	  110 XMAT(K,K)=AKK+SUM
67800	      CONS=TEMP+SUM*AKK
67900	      CON(K)=CONS
68000	      IF (N1 .LT. K)  GO TO 150
68100	  120 DO 140 J=K,N1
68200	      SUM=0.
68300	      DO 130 I=K,NROW
68400	  130 SUM=SUM+XMAT(I,K)*XMAT(I,J+1)
68500	      TEMP=SUM/CONS
68600	      DO 140 I=K,NROW
68700	  140 XMAT(I,J+1)=XMAT(I,J+1)-TEMP*XMAT(I,K)
68800	  150 CONTINUE
68900	      DO 210 K=1,NCOLX
69000	      KK=NCOLX-K+1
69100	      TEMP=XMAT(KK,KK)/CON(KK)
69200	      XMAT(KK,KK)=1.-TEMP*XMAT(KK,KK)
69300	      IF (NROW .EQ. KK) GO TO 170
69400	      KK1=KK+1
69500	      DO 160 I=KK1,NROW
69600	  160 XMAT(I,KK)=-TEMP*XMAT(I,KK)
69700	  170 KM=KK-1
69800	      IF (KM .EQ. 0) GO TO 220
69900	      DO 180 I = 1,KM
70000	  180 XMAT(I,KK) = 0.0
70100	      DO 210 L=1,KM
70200	      LL=KK-L
70300	      SUM=0.
70400	      DO 190 I=LL,NROW
70500	  190 SUM=SUM+XMAT(I,LL)*XMAT(I,KK)
70600	      TEMP=SUM/CON(LL)
70700	      DO 200 I=LL,NROW
70800	  200 XMAT(I,KK)=XMAT(I,KK)-TEMP*XMAT(I,LL)
70900	  210 CONTINUE
71000	220   RETURN
71100	      END
71200	      SUBROUTINE TRANSF(NVART,X,TRANCD)
71300	      DIMENSION X(50), TRANCD(50),SAVE(50)
71400	      INTEGER TRANCD
71500	      DO 90 I = 1,NVART
71600	   90 SAVE(I) = X(I)
71700	      DO 140 I = 1,NVART
71800	      IF (TRANCD(I) .GT. 3) GO TO 130
71900	      M = TRANCD(I) + 1
72000	      GO TO (140,100,110,120), M
72100	  100 X(I) = SQRT(X(I))
72200	      GO TO 140
72300	  110 X(I) = ALOG(X(I))
72400	      GO TO 140
72500	  120 X(I) = ASIN(X(I))
72600	      GO TO 140
72700	  130 CALL SPECTR(I,X,TRANCD,SAVE)
72800	  140 CONTINUE
72900	      RETURN
73000	      END
73100	      SUBROUTINE SPECTR(I,X,TRANCD,SAVE)
73200	      DIMENSION X(50),TRANCD(50),SAVE(50)
73300	      INTEGER  TRANCD
73400	C     LIST TO INDICATE ONE OF SEVERAL TRANSFORMATIONS.  THIS AND THE
73500	      XI = SAVE(I)
73600	      IND = TRANCD(I)
73700	      GO TO (1,1,1,4,5,6,7,8,9), IND
73800	    4 X(I) = 1.0/XI
73900	      GO TO 1
74000	    5 X(I) = ALOG(XI  +.5)
74100	      GO TO 1
74200	    6 X(I) = XI  /(1.0-XI   )
74300	      GO TO 1
74400	    7 X(I) = .5*ALOG( (1.0+XI  )/(1.0 - XI  ) )
74500	      GO TO 1
74600	    8 X(I) = ALOG(XI  /(1.0-XI   ) )
74700	      GO TO 1
74800	    9 X(I) = ASIN(SQRT(XI   ) )
74900	  1   RETURN
75000	      END
75100	      FUNCTION MSIGN(N)
75200	      MSIGN=N
75300	      RETURN
75400	      END
75500	      SUBROUTINE ALPHAN(N,NTABLE,IHOLV)
75600	      DIMENSION NTABLE (1),IHOLV(20)
75700	      NO=N*5
75800	      DECODE(NO,3,IHOLV)(NTABLE(I),I=1,N)
75900	3     FORMAT(20A4)
76000	      RETURN
76100	      END