Trailing-Edge
-
PDP-10 Archives
-
decus_20tap5_198111
-
decus/20-0137/twoaov/twoaov.for
There is 1 other file named twoaov.for in the archive. Click here to see a list.
C WESTERN MICHIGAN UNIVERSITY
C TWOAOV.F4 (FILENAME ON LIBRARY DECTAPE)
C TWOAOV, 1.9.2 (CALLING NAME, SUBLST #)
C TWO WAY ANALYSIS OF VARIANCE
C PROGRAMMED BY EVA GAINES AT WMU, LATER MODIFIED BY SAM ANEMA,
C RUSS BARR
C STATISTICAL CONSULTANT - DR. M. STOLINE
C FORWMU PROGRAMS USED: TTYPTY, EXISTS, ALLCOR, PRINTS, DEVCHG
C APLIB PROGRAMS USED: GETFOR, IO, FISHER
C LIBRARY DECTAPE PROGRAMS USED: USAGE
C INTERNAL SUBR. USED: MNL, GAINE
C ABOVE COMMENTS AND RIGHT ADJUSTED COMMENTS PUT IN BY WG
C
C---------------
C---------------MAX NO. OF OBS. IN ANY ONE CELL. I.E., 400 IN METHOD
C--------------- 1 (SEE ST. 42 AND 42-1 IN SUBR. MNL)
C---------------ID(16) LIMITS OUTPUT IDENTIFICATION TO 80 CHARACTERS
C---------------JEFFMT LIMITS OBJECTIVE FORMAT TO 80 CHARACTERS
C---------------S(2) IS USED IN CALL MNL
DIMENSION JEFFMT(16),Y(400),ID(16),S(1)
C---------------ASSIGNMENT OF I/O CHANNELS (B. GRANET)
IUT=3
INP=2
NOUTD=-1
IND=-4
C---------------PRINTS 5 LEVELS PER LINE ON FACTOR 2 FOR TERMINAL
C--------------- OR 12 FOR BATCH
KS=5
C---------------TTYPTY RETURNS ZERO - TTY JOB, MINUS ONE - BATCH JOB
CALL TTYPTY(ICOD)
IF(ICOD.EQ.-1)KS=12
WRITE(NOUTD,9872)
9872 FORMAT(20X,'WMU TWO-WAY ANALYSIS OF VARIANCE'/)
C CALL USAGE('TWOAOV')
C---------------IO DEALS WITH INPUT? AND OUTPUT? DIALOGUE
CALL IO(IUT,NDV1,NOUTD,IND,1,ICOD)
1000 CALL IO(INP,NDV2,NOUTD,IND,0,ICOD)
C---------------DEALS WITH FORMAT DIALOGUE
CALL GETFOR(NOUTD,IND,JEFFMT,ISTD,16,4)
WRITE(NOUTD,8291)
8291 FORMAT(' ENTER OUTPUT IDENTIFICATION IF DESIRED, OTHERWISE RETURN.
1'/)
READ(IND,8292)ID
8292 FORMAT(16A5)
WRITE(NOUTD,450)
450 FORMAT(' ENTER INPUT METHOD DESIRED.(1,2,OR 3)'/)
READ(IND,1002)METHOD
GO TO (1004,456,456),METHOD
1004 WRITE(NOUTD,1001)
1001 FORMAT(' NUMBER OF LEVELS ON FACTOR 1 = ',$)
READ(IND,1002)N1
IF(N1.LT.1)GO TO 7882
1002 FORMAT(12I)
WRITE(NOUTD,1003)
1003 FORMAT(' NUMBER OF LEVELS ON FACTOR 2 = ',$)
READ(IND,1002)N2
IF(N2.LT.1)GO TO 7882
C DYNAMIC ALLOCATION OF MEMORY
460 NM=MAX0(N1,N2)
NMAX=4*N1+5*N2+2*N1*N2+N2*N2+NM*N2
CALL ALLCOR(NMAX,IER,K1,S)
IF(IER.LT.0)GO TO 7782
IYS=K1
IY1=IYS+NM*N2
IY2=IY1+N1
IS1=IY2+N2
IS2=IS1+N1
IB1=IS2+N2
IYM1=IB1+N2
IT=IYM1+N1
IC=IT+N1*N2
IAM=IC+N2
IYT=IAM+N1
IYM2=IYT+N1*N2
IYR=IYM2+N2
IXXX=IYR+N2*N2
CALL MNL(N1,N2,NM,METHOD,IUT,INP,NOUTD,IND,KS,ICOD,NDV1,NDV2,
1S(IYS),S(IY1),S(IY2),S(IS1),S(IS2),S(IB1),S(IYM1),S(IT),S(IC),
2S(IAM),S(IYT),S(IYM2),S(IYR),JEFFMT,Y,ID,ISTD)
GO TO 1000
7882 WRITE(NOUTD,7883)
7883 FORMAT(' IMPROPER VALUE',/)
GO TO 1004
7782 WRITE(NOUTD,7783)
7783 FORMAT(' PROBLEM TOO LARGE',/)
GO TO 1000
456 N1=12
N2=12
GO TO 460
END
C---------------N1, N2, NM, METHOD, IUT, INP, NOUTD, INDP, KS,
C--------------- ICOD, NDV1, NDV2 ARE INPUT. OTHER ARGS. ARE SPACES
C--------------- RESERVED BY DYN. ALLOC.
SUBROUTINE MNL(N1,N2,NM,METHOD,IUT,INP,NOUTD,IND,KS,ICOD,
1NDV1,NDV2,YS,Y1,Y2,S1,S2,B1,YM1,T,C,AM,YT,YM2,YR,
2JEFFMT,Y,ID,ISTD)
DIMENSION YS(NM,N2),Y1(1),Y2(1),S1(1),S2(1),B1(1),YM1(1)
DIMENSION T(N1,N2),C(1),AM(1),YT(N1,N2),JEFFMT(16),Y(400),ID(16)
DIMENSION YM2(1),YR(N2,N2)
N1S=N1
N2S=N2
NMS=NM
DO 66 I=1,N1
AM(I)=0.
Y1(I)=0.
S1(I)=0.
DO 66 J=1,N2
T(I,J)=0.
YT(I,J)=0.
66 CONTINUE
DO 68 J=1,N2
Y2(J)=0.
C(J)=0.
B1(J)=0.
S2(J)=0.
DO 68 I=1,NM
68 YS(I,J)=0.
RAB=0.
SST=0.
SS=0.
G=0.
SSC=0.
SSAB=0.
AN=0.
SSBA=0.
C---------------INITIALIZE TO PRINT PRELIMINARY ANOVA (SEE LINE
C--------------- FOLLOWING "PRINT MEANS ETC." LOOP)
ITYP=1
TW=0.
TC=0.
GO TO(1004,456,456),METHOD
1004 WRITE(NOUTD,1020)
1020 FORMAT(' ENTER CELL SIZES'/)
DO 5 I=1,N1
5 READ(IND,1008) (T(I,J),J=1,N2)
1008 FORMAT(12F)
456 IF(ISTD.EQ.1)JEFFMT(1)='(10F)'
457 GO TO (500,703,471),METHOD
7782 WRITE(NOUTD,7783)
7783 FORMAT(' PROBLEM IS TOO LARGE'/)
GO TO 1004
C---------------METHOD 2 OF INPUT
703 IF(NDV2.EQ.'TTY')WRITE(NOUTD,700)
700 FORMAT(' ENTER DATA'/)
IF(NDV2.NE.'TTY')WRITE(NOUTD,830)
830 FORMAT(' YOUR DATA IS BEING READ.'/)
N1=1
N2=1
I=1
J=1
701 READ(INP,702,END=1500)IAST
702 FORMAT(A1)
IF(IAST.EQ.'*')GO TO 750
REREAD JEFFMT,YIN
YT(I,J)=YT(I,J)+YIN
YS(I,J)=YS(I,J)+YIN*YIN
T(I,J)=T(I,J)+1.0
GO TO 701
750 REREAD 751,I,J
751 FORMAT(1X,2I)
IF(I.LT.0.OR.I.GT.12.OR.J.LT.0.OR.J.GT.12)GO TO 752
IF(I.EQ.0.OR.J.EQ.0)GO TO 1500
IF(I.GT.N1)N1=I
IF(J.GT.N2)N2=J
GO TO 701
752 WRITE(NOUTD,753)
753 FORMAT(' INDEX OUT OF RANGE ENTRY IGNORED'//)
GO TO 701
C---------------METHOD 3 OF INPUT
471 IF(NDV2.EQ.'TTY')WRITE(NOUTD,700)
IF(NDV2.NE.'TTY')WRITE(NOUTD,830)
IF(ISTD.EQ.1)JEFFMT(1)='(2I,F'
IF(ISTD.EQ.1)JEFFMT(2)=') '
N1=0
N2=0
474 READ(INP,JEFFMT,END=1500)I,J,YIN
IF(I.LT.0.OR.I.GT.12.OR.J.LT.0.OR.J.GT.12)GO TO 472
IF(I.EQ.0.OR.J.EQ.0)GO TO 1500
IF(I.GT.N1)N1=I
IF(J.GT.N2)N2=J
YT(I,J)=YT(I,J)+YIN
YS(I,J)=YS(I,J)+YIN*YIN
T(I,J)=T(I,J)+1.0
GO TO 474
472 WRITE(NOUTD,753)
GO TO 474
C---------------METHOD 1 OF INPUT
500 IF(NDV2.NE.'TTY')WRITE(NOUTD,830)
DO 30 I=1,N1
DO 25 J=1,N2
IF(T(I,J))40,25,40
40 JK=T(I,J)
IF(NDV2.EQ.'TTY')WRITE(NOUTD,1025)I,J
1025 FORMAT(' ENTER DATA FOR CELL (',I2,',',I2,')'/)
READ (INP,JEFFMT)(Y(K),K=1,JK)
DO 42 K=1,JK
YT(I,J)=YT(I,J)+Y(K)
42 YS(I,J)=YS(I,J)+Y(K)*Y(K)
25 CONTINUE
30 CONTINUE
1500 IF(METHOD.NE.1.AND.(N1.GT.12.OR.N2.GT.12))GO TO 1005
C---------------PROCEED IF METHOD 2 OR 3 AND N1 AND N2 ARE NOT
C--------------- TOO LARGE
DO 4250 I=1,N1
DO 4251 J=1,N2
IF(T(I,J).NE.0.0)GO TO 4250
4251 CONTINUE
DO 4260 II=I,N1-1
DO 4261 J=1,N2
T(II,J)=T(II+1,J)
YT(II,J)=YT(II+1,J)
4261 YS(II,J)=YS(II+1,J)
4260 CONTINUE
N1=N1-1
4250 CONTINUE
DO 4252 J=1,N2
DO 4253 I=1,N1
IF(T(I,J).NE.0.0)GO TO 4252
4253 CONTINUE
DO 4262JJ=J,N2-1
DO 4263 I=1,N1
T(I,JJ)=T(I,JJ+1)
YT(I,JJ)=YT(I,JJ+1)
4263 YS(I,JJ)=YS(I,JJ+1)
4262 CONTINUE
N2=N2-1
4252 CONTINUE
DO 31 I=1,N1
DO 32 J=1,N2
IF(T(I,J))41,32,41
41 TC=TC+1.0
TW=TW+T(I,J)-1.0
Y1(I)=Y1(I)+YT(I,J)
S1(I)=S1(I)+T(I,J)
SSC=SSC+(YT(I,J)*YT(I,J))/T(I,J)
SST=SST+YS(I,J)
32 CONTINUE
SSAB=SSAB+(Y1(I)*Y1(I))/S1(I)
SS=SS+S1(I)
31 G=G+Y1(I)
IWM=0
ITYP=0
C---------------CHECK FOR EMPTY CELLS AND CALC. MEANS AND ST. DEV.
C--------------- FOR CELLS WITH MORE THAN ONE OBS. SEE LINE #250
C--------------- AND 825.
DO 400 J=1,N2
DO 300 I=1,N1
IF(T(I,J))250,825,250
250 IF(T(I,J).NE.1.0)ITYP=1
Y2(J)=Y2(J)+YT(I,J)
S2(J)=S2(J)+T(I,J)
YS(I,J)=YS(I,J)-(YT(I,J)*YT(I,J))/T(I,J)
IF(T(I,J)-1.)15,300,15
15 YS(I,J)=YS(I,J)/(T(I,J)-1.)
YS(I,J)=SQRT(YS(I,J))
YT(I,J)=YT(I,J)/T(I,J)
GO TO 300
C---------------WE HAVE AN EMPTY CELL
825 IWM=1
300 CONTINUE
400 SSBA=SSBA+(Y2(J)*Y2(J))/S2(J)
DO 312 I=1,N1
312 YM1(I)=Y1(I)/S1(I)
DO 313 I=1,N2
313 YM2(I)=Y2(I)/S2(I)
G1=(G*G)/SS
SST=SST-G1
SSC=SSC-G1
SSAB=SSAB-G1
SSBA=SSBA-G1
SSW=SST-SSC
NP=0
NP=NP+1
IF(TC-1.)83,997,83
83 SMC=SSC/(TC-1.)
NP=NP+1
IF(TW)93,6780,93
93 SMW=SSW/TW
NP=NP+1
6780 IF(SMW)53,6791,53
53 F=SMC/SMW
6791 KLL=FLOAT(N2)/FLOAT(KS)+.9
IF(ID(1).NE.' ')WRITE(IUT,8293)ID
8293 FORMAT('1',16A5)
K2=0
C---------------PRINT MEANS, SIZES, AND ST. DEV. OF CELLS
C--------------- AND MEANS OF LEVELS
DO 760 KL=1,KLL
K1=K2+1
K2=MIN0(N2,K1+KS-1)
WRITE(IUT,900)
900 FORMAT(//' ',35X,'FACTOR 2')
WRITE(IUT,90)(I,I=K1,K2)
90 FORMAT(/6X,'LEVEL - - - - - - -',12(I4,5X))
WRITE(IUT,91)
91 FORMAT(8X,'.')
WRITE(IUT,92)(YM2(I),I=K1,K2)
92 FORMAT(8X,'.',4X,'MEAN - - -',12F9.2)
WRITE(IUT,937)
937 FORMAT(8X,'.',5X,'.')
WRITE(IUT,94)(T(1,J),J=K1,K2)
94 FORMAT(19X,'SIZE',12F9.2)
WRITE(IUT,95)YM1(1),(YT(1,J),J=K1,K2)
95 FORMAT(8X,'1',1X,F8.2,' MEAN',12F9.2)
WRITE(IUT,96)(YS(1,J),J=K1,K2)
96 FORMAT(18X,'ST DEV',12(F8.2,1X))
IF(N1.EQ.1)GO TO 760
WRITE(IUT,977)(T(2,J),J=K1,K2)
977 FORMAT(/1X,'FACTOR',12X,'SIZE',12F9.2)
WRITE(IUT,98)YM1(2),(YT(2,J),J=K1,K2)
98 FORMAT(4X,'1',3X,'2',F9.2,' MEAN',12F9.2)
WRITE(IUT,96)(YS(2,J),J=K1,K2)
IF(N1.LT.3)GO TO 760
DO 897 I=3,N1
WRITE(IUT,898)(T(I,J),J=K1,K2)
898 FORMAT(/19X,'SIZE',12F9.2)
WRITE(IUT,899)I,YM1(I),(YT(I,J),J=K1,K2)
899 FORMAT(7X,I2,F9.2,' MEAN',12F9.2)
897 WRITE(IUT,96)(YS(I,J),J=K1,K2)
760 CONTINUE
IF(ITYP.EQ.0)GO TO 6781
C---------------PREL. ANOVA ONLY WHEN AT LEAST ONE CELL HAS
C--------------- MORE THAN ONE OBS.
80 WRITE (IUT,39)
39 FORMAT(//25X17HPRELIMINARY ANOVA)
WRITE (IUT,49)
49 FORMAT(6X6HSOURCE,17X2HDF,8X2HSS,8X2HMS,8X1HF,6X,4HPROB)
TP=TC-1.
PROB=FISHER(IFIX(TP),IFIX(TW),F)
WRITE (IUT,59)TP,SSC,SMC,F,PROB
59 FORMAT(6X5HCELLS,10X4F10.2,F9.3)
TP=N1-1
I=1
J=2
WRITE (IUT,69)I,J,TP,SSAB
69 FORMAT(4XI1,9H IGNORING,I2,5X2F10.2)
TP=N2-1
WRITE (IUT,69)J,I,TP,SSBA
WRITE (IUT,79)TW,SSW,SMW
79 FORMAT(6X6HWITHIN,9X3F10.2)
TP=SS-1.
WRITE (IUT,89)TP,SST
89 FORMAT(6X5HTOTAL,10X2F10.2)
IF(N1.EQ.1.OR.N2.EQ.1)GO TO 1000
C---------------START WORK ON LEAST SQ. ANOVA
C---------------PUT IDENTITY MATRIX IN YS FOR CALL GAINE BELOW,
C--------------- INIT. YR TO 0
6781 DO 33 I=1,N2
DO 33 J=1,N2
YR(I,J)=0.
IF(I-J)133,37,133
37 YS(I,J)=1.
GO TO 33
133 YS(I,J)=0.
33 CONTINUE
C---------------CALC. NORMAL (LEAST SQ.) EQS. COEFFS. IN YR AND C
DO 180 I=1,N2
DO 180 J=1,N2
DO 180 K=1,N1
IF(S1(K))180,999,180
999 WRITE(NOUTD,409)K,S1(K)
409 FORMAT(1X,3HS1(,I3,2H)=,F6.0)
180 YR(I,J)=YR(I,J)+(T(K,I)*T(K,J))/S1(K)
DO 102 I=1,N2
102 YR(I,I)=YR(I,I)-S2(I)
DO 120 J=1,N2
171 DO 110 K=1,N1
IF(S1(K))110,992,110
992 WRITE(NOUTD,409)K,S1(K)
110 C(J)=C(J)+(T(K,J)*Y1(K))/S1(K)
YR(N2,J)=1.
120 C(J)=C(J)-Y2(J)
C(N2)=0.
C---------------CALC. INVERSE OF YR, RETURN RESULT IN YS
CALL GAINE(YR,YS,N2,N2S,NMS,IUT)
C---------------SOLVE NORMAL (LEAST SQ.) EQS.
DO 280 J=1,N2
DO 280 K=1,N2
280 B1(J)=B1(J)+YS(J,K)*C(K)
DO 210 I=1,N1
DO 220 J=1,N2
220 AM(I)=AM(I)+T(I,J)*B1(J)
AM(I)=Y1(I)-AM(I)
IF(S1(I))13,993,13
993 WRITE(NOUTD,409)I,S1(I)
13 AM(I)=AM(I)/S1(I)
210 AN=AN+AM(I)
P1=N1
AN=AN/P1
DO 230 I=1,N1
A=AM(I)-AN
230 RAB=RAB+A*Y1(I)
DO 260 J=1,N2
260 RAB=RAB+B1(J)*Y2(J)
RAB=RAB+AN*G-G1
SI=SSC-RAB
P2=N2
TI=TC-P1-P2+1.
MP=0
MP=MP+1
IF(TI)77,998,77
77 SMI=SI/TI
TAB=RAB-SSBA
MP=MP+1
IF(P1-1.)87,998,87
87 TMAB=TAB/(P1-1.)
TBA=RAB-SSAB
MP=MP+1
IF(P2-1.)97,998,97
97 TMBA=TBA/(P2-1.)
MP=MP+1
IDF=IFIX(TW)
IF(ITYP.EQ.0)SMW=SMI
IF(ITYP.EQ.0)IDF=IFIX(TI)
IF(SMW)107,998,107
107 FI=SMI/SMW
FA=TMAB/SMW
FB=TMBA/SMW
WRITE (IUT,109)
109 FORMAT(//25X,'LEAST SQUARES ANOVA')
WRITE (IUT,49)
TP=TC-1.
WRITE (IUT,59)TP,SSC
TP=P1-1.
I=1
J=2
PROB=FISHER(IFIX(TP),IDF,FA)
WRITE (IUT,204)I,J,TP,TAB,TMAB,FA,PROB
204 FORMAT(3XI1,12H ELIMINATING,I2,3X4F10.2,F9.3)
TP=P2-1.
PROB=FISHER(IFIX(TP),IDF,FB)
WRITE (IUT,304)J,I,TP,TBA,TMBA,FB,PROB
304 FORMAT(3XI1,12H ELIMINATING,I2,3X4F10.2,F9.3)
PROB=FISHER(IFIX(TI),IDF,FI)
C---------------ITYP=0 MEANS EVERY CELL HAS ONE OBS., ITYP=1
C--------------- MEANS AT LEAST ONE CELL HAS MORE THAN ONE OBS.
IF(ITYP.EQ.0)WRITE(IUT,404)I,J,TI,SI,SMI
IF(ITYP.NE.0)WRITE (IUT,404)I,J,TI,SI,SMI,FI,PROB
404 FORMAT(6XI1,3H BY,I2,9X4F10.2,F9.3)
IF(ITYP.EQ.0)GO TO 6783
WRITE (IUT,79)TW,SSW,SMW
6783 TP=SS-1.
WRITE (IUT,89)TP,SST
C
C CHECK FOR BALANCED CASE
C
T11=T(1,1)
DO 821 I=1,N1
DO 821 J=1,N2
IF(T11.NE.T(I,J))GO TO 823
821 CONTINUE
WRITE(IUT,822)
822 FORMAT(/,' THE DESIGN IS BALANCED. IN THIS CASE THE WEIGHTED',/,
#' MEANS AND LEAST SQUARES ANOVA''S ARE IDENTICAL.',//)
GO TO 1000
C
C CHECK FOR PROPORTIONALITY
C
823 DO 824 I=2,N1
TI1=T(I,1)
DO 824 J=2,N2
IF(T11*T(I,J).NE.T(1,J)*TI1)GO TO 826
824 CONTINUE
WRITE(IUT,831)
831 FORMAT(/,' THE DESIGN IS PROPORTIONAL, BUT NOT BALANCED.',/,
#' IN THIS CASE THE LEAST SQUARES ANOVA GIVES THE',/,
#' CORRECT SUM OF SQUARES AND MEAN SQUARES, BUT NOT',/,
#' CORRECT F-TESTS. SEE T.A. BANCROFT "TOPICS IN',/,
#' INTERMEDIATE STATISTICAL METHODS(P. 14-15)" FOR',/,
#' THE PROPER F-TEST PROCEDURES FOR THE PROPOR-',/,
#' TIONAL CASE.',//)
GO TO 1000
C
C WEIGHTED MEANS ANALYSIS
C
C---------------WE HAVE NON-PROPORTIONALITY. IWM=1 MEANS AT LEAST
C--------------- ONE EMPTY CELL, SEE LINE ABOVE #250
826 IF(IWM.EQ.1)GO TO 870
DO 811 I=1,N1
YM1(I)=0.0
811 S1(I)=0.0
DO 812 J=1,N2
YM2(J)=0.0
812 S2(J)=0.0
DO 813 I=1,N1
DO 814 J=1,N2
YM1(I)=YM1(I)+YT(I,J)
S1(I)=S1(I)+1.0/T(I,J)
YM2(J)=YM2(J)+YT(I,J)
814 S2(J)=S2(J)+1.0/T(I,J)
YM1(I)=YM1(I)/FLOAT(N2)
813 S1(I)=FLOAT(N2**2)/S1(I)
DO 815 J=1,N2
YM2(J)=YM2(J)/FLOAT(N1)
815 S2(J)=FLOAT(N1**2)/S2(J)
A=0.0
B=0.0
D=0.0
DO 816 I=1,N1
A=A+S1(I)*YM1(I)**2
B=B+S1(I)*YM1(I)
816 D=D+S1(I)
SSR=A-B**2/D
A=0.0
B=0.0
D=0.0
DO 817 J=1,N2
A=A+S2(J)*YM2(J)**2
B=B+S2(J)*YM2(J)
817 D=D+S2(J)
SSCO=A-B**2/D
SMR=SSR/FLOAT(N1-1)
SMCO=SSCO/FLOAT(N2-1)
C
C OUTPUT ANOVA TABLE FOR WEIGHTED MEANS
C
WRITE(IUT,818)
818 FORMAT(//25X,'WEIGHTED MEANS ANOVA')
WRITE(IUT,49)
TP=TC-1
WRITE(IUT,59)TP,SSC,SMC
TP=P1-1.0
F=SMR/SMW
PROB=FISHER(IFIX(TP),IDF,F)
I=1
J=2
WRITE(IUT,819)I,TP,SSR,SMR,F,PROB
819 FORMAT(8X,I1,12X,4F10.2,F9.3)
TP=P2-1.0
F=SMCO/SMW
PROB=FISHER(IFIX(TP),IDF,F)
WRITE(IUT,819)J,TP,SSCO,SMCO,F,PROB
PROB=FISHER(IFIX(TI),IDF,FI)
IF(ITYP.EQ.0)WRITE(IUT,404)I,J,TI,SI,SMI
IF(ITYP.NE.0)WRITE(IUT,404)I,J,TI,SI,SMI,FI,PROB
IF(ITYP.EQ.0)GO TO 820
WRITE(IUT,79)TW,SSW,SMW
820 TP=SS-1.0
WRITE(IUT,89)TP,SST
GO TO 1000
870 WRITE(IUT,871)
871 FORMAT(//' NO WEIGHTED MEANS ANALYSIS CAN BE PERFORMED WITH'
1,' EMPTY CELLS'//)
GO TO 1000
1005 WRITE(NOUTD,1007)
1007 FORMAT(' PROBLEM TOO LARGE'/)
IF(ICODE.EQ.-1)CALL EXIT
GO TO 1000
C VALUES OF NP MAY BE 1,2,3
C 1 IMPLIES TC-1=0,2 IMPLIES TW=0,3 IMPLIES SMW=0
997 WRITE(NOUTD,410)NP
410 FORMAT(1X,3HNP=,I2)
GO TO(83,93,53),NP
C VALUES OF MP MAY BE 1,2,3,4
C 1 IMPLIES TI=0,2 IMPLIES P1-1=0,3 IMPLIES P2-1=0,4 IMPLIES SMW=0
998 WRITE(NOUTD,1410)MP
1410 FORMAT(1X,3HMP=,I2)
GO TO(77,87,97,107),MP
1000 RETURN
END
C---------------A, N, IUT, N2 (USED IN DIM. ST.), NM (USED IN
C--------------- DIM. ST.) ARE INPUT. B RETURNED
C--------------- CALCULATE INVERSE OF A, STORE INVERSE IN B
SUBROUTINE GAINE(A,B,N,N2,NM,IUT)
DIMENSION A(N2,N2),B(NM,N2)
I=0
100 I=I+1
IF(A(I,I)-100.+100.)19,8,19
19 DIA=A(I,I)
C MAKE DIAGONAL ELEMENT ONE
DO 39 J=1,N
A(I,J)=A(I,J)/DIA
39 B(I,J)=B(I,J)/DIA
IF(I-N)29,18,29
C MAKE ALL ELEMENTS BELOW THE DIAGONAL ELEMENT ZERO
29 M=I+1
DO 11 L=M,N
XA=A(L,I)
DO 11 K=1,N
A(L,K)=A(L,K)-XA*A(I,K)
11 B(L,K)=B(L,K)-XA*B(I,K)
IF(I-1)18,100,18
C MAKE ALL ELEMENTS ABOVE THE DIAGONAL ELEMENT ZERO
18 KP=0
MN=I-1
DO 40 J=1,MN
KP=KP+1
YA=A(KP,I)
DO 40 L=1,N
A(J,L)=A(J,L)-YA*A(I,L)
40 B(J,L)=B(J,L)-YA*B(I,L)
IF(I-N)100,99,100
8 MN=I
48 MN=MN+1
IF(MN-N)9,9,90
C MAKE DIAGONAL ELEMENT NON-ZERO
9 DO 6 J=1,N
A(I,J)=A(MN,J)+A(I,J)
6 B(I,J)=B(MN,J)+B(I,J)
IF(A(I,I)-100.+100.)19,48,19
90 WRITE (IUT,80)
80 FORMAT(1X,26HTHE INVERSE DOES NOT EXIST)
99 RETURN
END