Google
 

Trailing-Edge - PDP-10 Archives - decus_20tap5_198111 - decus/20-0137/bmd/bmd02m.for
There is 1 other file named bmd02m.for in the archive. Click here to see a list.
C        REGRESSION ON PRINCIPLE COMPONENTS             JUNE 11,1967
C     THIS IS A SYSTEM/360 FORTRAN IV (LEVEL H) VERSION OF BMD02M
C     ORIGINALLY WRITTEN IN FORTRAN II.  THE PROGRAM HAS BEEN SIFTED
C     AND SLIGHTLY MODIFIED TO MAKE IT OPERABLE ON THE 7094.
C     IT HAS BEEN FURTHER MODIFIED TO ALLOW ITS EXECUTION ON THE 360.
C
      DIMENSION X(200,30),XMEAN(45),COV(200,15),VALU(25),SCALE(45),
     1C(200,15),Z(25,25),INDEX(20),YMEAN(20),SQ(20),RES(20),FMT(180)
     2,NAMES(45),CCC(25,25),DUMY1(6000),DUMY2(45)
      DOUBLE PRECISION TODE,B123,     NPROB,A123,NAMES
      COMMON  X      , XMEAN  , COV
      COMMON  VALU   , SCALE  , C      , Z      , INDEX  , YMEAN
      COMMON  SQ     , RES    , NK     , NY     , NV     , N
      COMMON  NT
      EQUIVALENCE (DUMY1,X),(DUMY2,XMEAN),(COV,CCC)
      NTAPE=5
	CALL USAGEB('BMD02M')
      DATA A123/6HPROBLM/,B123/6HFINISH/,D123/3HYES/
  204 FORMAT (' EXPECTED PROBLM CARD, BUT READ ',A6)
  213 FORMAT(55H1BMD02M-REGRESSION ON PRINCIPAL COMPONENTS - VERSION OF,
     118H  JUNE 11, 1967   ,/
     240H HEALTH SCIENCES COMPUTING FACILITY,UCLA/
     314H PROBLEM CODE A6,/
     427H NO. OF ORIGINAL VARIABLES I3,/
     520H NO. OF X VARIABLES I3,/
     620H NO. OF Y VARIABLES I3,//)
 10   READ (5,921)TODE,NPROB,NT,NV,NK,N,PNCR,GCK,NVG,NLV,MTAPE,KVR
      IERROR=0
      IF(TODE.NE.B123)GO TO 200
      GO TO 201
  202 WRITE (6,204) TODE
 201  IF(NTAPE-5)7,7,6
    6 REWIND NTAPE
    7 STOP
  200 IF(TODE.NE.A123)GO TO 202
  203 CALL TPWD(MTAPE,NTAPE)
 9    IF(NV*(NV-26))205,944,944
 205  IF(NK*(NK-20)) 206,206,946
  206 IF(NV+NK-45) 207,207,948
  207 IF((N-2)*(N-201)) 208,950,950
  208 IF(KVR.GT.0.AND.KVR.LE.10)GO TO 209
      KVR=1
      WRITE (6,4000)
  209 NT1=NV+NK
      WRITE (6,213)NPROB,NT,NV,NK
      CALL RDLBL(NLV,NT1,NAMES)
   17 KVR=KVR*18
      READ (5,942)(FMT(I),I=1,KVR)
      WRITE (6,943) (FMT(I),I=1,KVR)
  943 FORMAT (' VARIABLE FORMAT IS '/(5X,18A4))
      DO 13 I=1,N
   13 READ (NTAPE,FMT)(X(I,J),J=1,NT)
      IF(NVG) 952,19,1000
 1000 CALL TRANS(X,NT,N,IERROR,NVG)
      IF(IERROR)10,19,19
   19 IF(NT) 954,954,212
  212 NT=NT1
      ON=N
      DO 21 J=1,NT
      XMEAN(J)=0.0
      DO 20 I=1,N
  20  XMEAN(J)=XMEAN(J)+X(I,J)
  21  XMEAN(J)=XMEAN(J)/ON
      DO 90 I=1,NK
      IJ=NV+I
 90   YMEAN(I)=XMEAN(IJ)
      DO 22 I=1,NV
      DO 22 J=1,NV
      COV(I,J)=0.0
      DO 22 K=1,N
  22  COV(I,J)=COV(I,J)+(X(K,I)-XMEAN(I))*(X(K,J)-XMEAN(J))
      DO 23 I=1,NV
  23  SCALE(I)=SQRT(COV(I,I))
      DO 24 I=1,NV
      DO 24 J=1,NV
  24  Z(I,J)=COV(I,J)/(SCALE(I)*SCALE(J))
      WRITE (6,936)
      WRITE (6,904)
      CALL PATTY2(Z,NV,NAMES,1)
      CALL EIGEN(VALU,NV,NV)
      WRITE (6,936)
      WRITE (6,907)
      WRITE (6,906)(VALU(I),I=1,NV)
      RANK=0.0
      DO 26 I=1,NV
 26   RANK=RANK+VALU(I)
      SMALL=0.0
      DO 18 I=1,NV
      SMALL=SMALL+VALU(I)
   18 C(I,1)=SMALL/RANK
      WRITE (6,937)
      WRITE (6,936)
      WRITE (6,938)(C(I,1),I=1,NV)
      WRITE (6,908)
      CALL PATTY2(Z,NV,NAMES,0)
      DO 29 J=1,NV
      DO 29 I=1,N
  29  X(I,J)=(X(I,J)-XMEAN(J))/SCALE(J)
      HIP= N-1
      SQT=SQRT(HIP)
      DO 43 I=1,N
      DO 43 J=1,NV
      C(I,J)=0.0
      DO 43 K=1,NV
  43  C(I,J)=C(I,J)+X(I,K)*Z(K,J)
      IF(D123.NE.GCK)GO TO 57
  41  DO 51 J=1,NV
      XMEAN(J)=0.0
      DO 50 I=1,N
  50  XMEAN(J)=XMEAN(J)+C(I,J)
  51  XMEAN(J)=XMEAN(J)/ON
      DO 52 I=1,NV
      DO 52 J=1,NV
      CCC(I,J)=0.0
      DO 52 K=1,N
   52 CCC(I,J)=CCC(I,J)+(C(K,I)-XMEAN(I))*(C(K,J)-XMEAN(J))
      WRITE (6,936)
      WRITE (6,922)
      CALL PATTY2(CCC,NV,NAMES,1)
  57  DO 58 I=1,N
      DO 58 J=1,NV
  58  COV(I,J)=C(I,J)
      IF(D123.NE.PNCR)GO TO 40
  30  WRITE (6,936)
      WRITE (6,909)
      WRITE (6,910)
      SMALL=-(10.0**36.0)
      DO 39 II=1,NV
      WRITE (6,912)II
      DO 32 I=1,N
      C(I,1)=0.0
      C(I,2)=0.0
      DO 31 K=1,NV
  31  C(I,1)=C(I,1)+X(I,K)*Z(K,II)
 32   C(I,1)=C(I,1)*SQT
      DO 39 I=1,N
      RANK=SMALL
      DO 38 J=1,N
      IF(C(J,1)-RANK)38,38,36
  36  IF(C(J,2)-999.0)37,38,38
  37  RANK=C(J,1)
      NJ=J
  38  CONTINUE
      C(NJ,2)=999.0
      WRITE (6,911)RANK,NJ
  39  CONTINUE
      IF(NK.EQ.0) GO TO 10
  40  CALL REGRES
      GO TO 10
 900  FORMAT(12F6.0)
 903  FORMAT(12H PROBLEM NO.I4)
 904  FORMAT(31H0CORRELATION COEFFICIENT MATRIX)
  906 FORMAT(8F12.4)
  907 FORMAT(12H0EIGENVALUES/)
  908 FORMAT(13H0EIGENVECTORS/)
 909  FORMAT(48H0RANK ORDER OF EACH STANDARDIZED CASE ORDERED BY)
 910  FORMAT(44H SIZE OF EACH PRINCIPAL COMPONENT SEPARATELY)
  911 FORMAT(F18.6,I10)
 912  FORMAT(16H0  COMPONENT NO.I3,12H    CASE NO.)
 920  FORMAT(43H0REGRESSION ON PRIMARY PRINCIPAL COMPONENTS)
  921 FORMAT(2A6,3I2,I3,2A3,2I2,37X,2I2)
 922  FORMAT(25H0EIGEN VALUE CHECK MATRIX)
 936  FORMAT(1H0)
 937  FORMAT(40H0CUMULATIVE PROPORTION OF TOTAL VARIANCE)
 938  FORMAT(1H F11.2,7F15.2)
 939  FORMAT(7H VECTORI3)
  942 FORMAT (18A4)
  944  WRITE (6,945) NV
  945 FORMAT (1X,I3,' INDEPENDENT VARIABLES IS ILLEGAL ')
      GO TO 201
  946 WRITE (6,947) NK
  947 FORMAT (1X,I3,' DEPENDENT VARIABLES IS ILLEGAL')
      GO TO 201
  948  NVNK=NV+NK
      WRITE (6,949) NVNK
  949 FORMAT (' TOTAL OF ',I4,' VARIABLES IS ILLEGAL')
      GO TO 201
  950 WRITE (6,951) N
  951  FORMAT (1X,I4,' CASES IS ILLEGAL')
      GO TO 201
  952 WRITE (6,953) NVG
  953 FORMAT (1X,I3,' TRANSGENERATION CARDS IS ILLEGAL')
       GO TO 201
  954  WRITE (6,955) NT
  955 FORMAT (1X,I3,' ORIGINAL VARIABLES IS ILLEGAL')
      GO TO 201
 4000 FORMAT(1H023X71HNUMBER OF VARIABLE FORMAT CARDS INCORRECTLY SPECIF
     1IED, ASSUMED TO BE 1.)
      END
      SUBROUTINE EIGEN(VALU,N,M)
C
C        SUBROUTINE EIGEN FOR BMD02M                  AUGUST 17, 1966
C     EIGENVALUES AND EIGENVECTORS OF A REAL SYMMETRIC MATRIX
C
      DOUBLE PRECISION  ANORM2
      DIMENSION A(25,25), B(25,25), VALU(25), DIAG(25), SUPERD(24),
     1          Q(24), VALL(25), S(24), C(24), D(25), IND(25), U(25)
      DIMENSION DUMY4(3070)
      DIMENSION X123(200,30),DUMY1(45),XMEAN(45),DUMY2(3000),DUMY3(6000)
      COMMON  X123   , XMEAN  , DUMY2
      COMMON  DUMY4  , A
      EQUIVALENCE (X123,DUMY3),(XMEAN,DUMY1),(DUMY2,B),(SUPERD,DUMY2(675
     1)),(DIAG,DUMY2(700)),(VALL,D,DUMY2(725)),(Q,S,DUMY2(750)),(IND,U,D
     2UMY2(775)),(C,DUMY2(800)),(TAU,BETA),(ANORM,ANORM2),(P,PRODS),(T,S
     3MALLD),(II,MATCH)
C
C     CALCULATE NORM OF MATRIX
C
    3 ANORM2=0.0
    4 DO 6 I=1,N
    5 DO 6 J=1,N
    6 ANORM2=ANORM2+A(I,J)**2
    7 ANORM=ANORM2
      ANORM=SQRT(ANORM)
C
C     GENERATE IDENTITY MATRIX
C
    9 IF (M) 10, 45, 10
   10 DO 40 I=1,N
   12 DO 40 J=1,N
   20 IF(I-J) 35, 25, 35
   25 B(I,J)=1.0
   30 GO TO 40
   35 B(I,J)=0.0
   40 CONTINUE
C
C     PERFORM ROTATIONS TO REDUCE MATRIX TO JACOBI FORM
C
   45 IEXIT=1
   50 NN=N-2
   52 IF (NN) 890, 170, 55
   55 DO 160 I=1,NN
   60 II=I+2
   65 DO 160 J=II,N
   70 T1=A(I,I+1)
   75 T2=A(I,J)
   80 GO TO 900
   90 DO 105 K=I,N
   95 T2=COS*A(K,I+1)+SIN*A(K,J)
  100 A(K,J)=COS*A(K,J)-SIN*A(K,I+1)
  105 A(K,I+1)=T2
  110 DO 125 K=I,N
  115 T2=COS*A(I+1,K)+SIN*A(J,K)
  120 A(J,K)=COS*A(J,K)-SIN*A(I+1,K)
  125 A(I+1,K)=T2
  128 IF (M) 130, 160, 130
  130 DO 150 K=1,N
  135 T2=COS*B(K,I+1)+SIN*B(K,J)
  140 B(K,J)=COS*B(K,J)-SIN*B(K,I+1)
  150 B(K,I+1)=T2
  160 CONTINUE
C
C     MOVE JACOBI FORM ELEMENTS AND INITIALIZE EIGENVALUE BOUNDS
C
  170 DO 200 I=1,N
  180 DIAG(I)=A(I,I)
  190 VALU(I)=ANORM
  200 VALL(I)=-ANORM
  210 DO 230 I=2,N
  220 SUPERD(I-1)=A(I-1,I)
  230 Q(I-1)=(SUPERD(I-1))**2
C
C     DETERMINE SIGNS OF PRINCIPAL MINORS
C
  235 TAU=0.0
  240 I=1
  260 MATCH=0
  270 T2=0.0
  275 T1=1.0
  277 DO 450 J=1,N
  280 P=DIAG(J)-TAU
  290 IF(T2) 300, 330, 300
  300 IF(T1) 310, 370, 310
  310 T=P*T1-Q(J-1)*T2
  320 GO TO 410
  330 IF(T1) 335, 350, 350
  335 T1=-1.0
  340 T=-P
  345 GO TO 410
  350 T1=1.0
  355 T=P
  360 GO TO 410
  370 IF(Q(J-1)) 380, 350, 380
  380 IF(T2) 400, 390, 390
  390 T=-1.0
  395 GO TO 410
  400 T=1.0
C
C     COUNT AGREEMENTS IN SIGN
C
  410 IF(T1) 425, 420, 420
  420 IF(T) 440, 430, 430
  425 IF(T) 430, 440, 440
  430 MATCH=MATCH+1
  440 T2=T1
  450 T1=T
C
C     ESTABLISH TIGHTER BOUNDS ON EIGENVALUES
C
  460 DO 530 K=1,N
  465 IF (K-MATCH) 470, 470, 520
  470 IF(TAU-VALL(K)) 530, 530, 480
  480 VALL(K)=TAU
  490 GO TO 530
  520 IF(TAU-VALU(K)) 525, 530, 530
  525 VALU(K)=TAU
  530 CONTINUE
 540  IF(ABS(VALU(I)-VALL(I)).GT.ABS(VALU(I)+VALL(I))/1.E+6) GO TO 580
  570 I=I+1
  575 IF(I-N) 540, 540, 590
  580 TAU=(VALL(I)+VALU(I))/2.0
  585 GO TO 260
C
C     JACOBI EIGENVECTORS BY ROTATIONAL TRIANGULARIZATION
C
  590 IF (M) 593, 890, 593
  593 IEXIT=2
  595 DO 610 I=1,N
  600 DO 610 J=1,N
  610 A(I,J)=0.0
  615 DO 850 I=1,N
  620 IF (I-1) 625, 625, 621
 621  IF(ABS(VALU(I-1)-VALU(I)).LT.ABS(VALU(I-1)+VALU(I))/1.E+6) GOTO730
  625 COS=1.0
  628 SIN=0.0
  630 DO 700 J=1,N
  635 IF(J-1) 680, 680, 640
  640 GO TO 900
  650 S(J-1)=SIN
  660 C(J-1)=COS
  670 D(J-1)=T1*COS+T2*SIN
  680 T1=(DIAG(J)-VALU(I))*COS-BETA*SIN
  690 T2=SUPERD(J)
  700 BETA=SUPERD(J)*COS
  710 D(N)=T1
  720 DO 725 J=1,N
  725 IND(J)=0
  730 SMALLD=ANORM
  735 DO 780 J=1,N
  740 IF (IND(J)-1) 750, 780, 780
  750 IF (ABS(SMALLD)-ABS(D(J)))780, 780, 760
  760 SMALLD=D(J)
  770 NN=J
  780 CONTINUE
  790 IND(NN)=1
  800 PRODS=1.0
  805 IF (NN-1) 810, 850, 810
  810 DO 840 K=2,NN
  820 II=NN+1-K
  830 A(II+1,I)=C(II)*PRODS
  840 PRODS=-PRODS*S(II)
  850 A(1,I)=PRODS
C
C     FORM MATRIX PRODUCT OF ROTATION MATRIX WITH JACOBI VECTOR MATRIX
C
  855 DO 885 J=1,N
  860 DO 865 K=1,N
  865 U(K)=A(K,J)
  870 DO 885 I=1,N
  875 A(I,J)=0.0
  880 DO 885 K=1,N
  885 A(I,J)=B(I,K)*U(K)+A(I,J)
  890 GO TO 941
C
C     CALCULATE SINE AND COSINE OF ANGLE OF ROTATION
C
  900 IF (T2) 910, 940, 910
  910 T=SQRT(T1**2+T2**2)
  920 COS=T1/T
  925 SIN=T2/T
  930 GO TO (90,650), IEXIT
  940 GO TO (160,910), IEXIT
 941  CONTINUE
      RETURN
      END
      SUBROUTINE PATTY2(A,N,NAMES,JK)
C             SUBROUTINE PATTY2 FOR BMD02M            AUGUST 17, 1966
      DOUBLE PRECISION NAMES,NN
      DIMENSION A(25,25),NAMES(45),NN(8)
      IT=1
      KK=0
      K1=IT
      K2=MIN0(8,N)
    5 KK=KK+8
      IF(N-KK)3,3,4
    4 IT=IT+1
      GO TO 5
    3 DO 50 JX=1,IT
      LLL=K2-K1+1
      LL=0
      IF(JK)35,35,37
   35 WRITE (6,350)(IG,IG=1,LLL)
      GO TO 45
   37 DO 40 JJ=K1,K2
      LL=LL+1
   40 NN(LL)=NAMES(JJ)
      WRITE (6,300)(NN(II),II=1,LLL)
   45 DO 10 I=1,N
   10 WRITE (6,20)NAMES(I),(A(I,J),J=K1,K2)
      K1=K2+1
      K2=K1+7
      K2=MIN0(K2,N)
   50 CONTINUE
      RETURN
  300 FORMAT(1H013X,A6,7(8X,A6)/)
   20 FORMAT(1H A6,1X,8F14.4)
  350 FORMAT(1H017X,I2,7(12X,I2)/)
      END
      DOUBLE PRECISION FUNCTION ANUMB(I)
      DOUBLE PRECISION X
      DIMENSION B(2)
      ENCODE(8,10,B)I
10    FORMAT(4X,I4)
      DECODE(8,11,B)X
11    FORMAT(A8)
      ANUMB=X
      RETURN
      END
      SUBROUTINE RDLBL(NLBVAR,NVAR,ARRAY)
C
C        SUBROUTINE RDLBL FOR BMD02M                  AUGUST 17, 1966
C     SUBROUTINE TO READ IN LABELS CARDS, STORE THEM IN ARRAY,
C     AND SUBSTITUTE NUMBERS FOR UNLABELED VARIABLES
C     NVAR IS TOTAL NUMBER OF VARIABLES
C     NLBVAR IS NUMBER OF LABELED VARIABLES EXPECTED
      DOUBLE PRECISION ARRAY,ANUMB,DUMY
      DIMENSION ARRAY (45),IDUM(7),DUMY(7)
      DATA ALABEL/3HLAB/
C     NUMBER VARIABLES
      DO 1 I=1,NVAR
   1  ARRAY(I)=ANUMB(I)
C     IF NO LABELS, RETURN
      IF(NLBVAR) 9,9,2
   2  N=0
C     READ 1 LABELS CARD
  20  READ (5,3) TEST,(IDUM(J),DUMY(J),J=1,7)
   3  FORMAT(A3,3X,7(I4,A6))
C     TEST FOR 'LAB' IN FIRST 3 COLS.
      IF(TEST.EQ.ALABEL)GO TO 6
C     ERROR--PRINT MESSAGE AND QUIT
   4  WRITE (6,5)
   5  FORMAT(36H0LABELS CARD NOT FOUND WHEN EXPECTED)
      STOP
C     EXAMINE 7 FIELDS
   6  DO 8 J=1,7
      K=IDUM(J)
C     TEST INDEX.  IF 0, IGNORE.  IF ILLEGAL, PRINT MESSAGE AND
C     IGNORE EXCEPT TO COUNT
      IF(K) 11,8,10
  10  IF(K-NVAR) 7,7,11
  11  WRITE (6,12)K,DUMY(J)
  12  FORMAT(18H0LABELS CARD INDEX,I7,18H INCORRECT. LABEL ,A6,9H IGNORE
     1D.)
      GO TO 13
C     MOVE LABEL TO ARRAY
   7  ARRAY(K)=DUMY(J)
C     STEP NUMBER OF VARIABLES
  13  N=N+1
C     TEST FOR END. IF END, RETURN. IF NOT, SCAN OTHER FIELDS.
      IF(N-NLBVAR) 8,9,9
   8  CONTINUE
      WRITE (6,7000)(ARRAY (I),I=1,15)
 7000 FORMAT (9X,A6,5X,A6,5X,A6,5X,A6,5X,A6)
      GO TO 20
   9  RETURN
      END
      SUBROUTINE REGRES
C          SUBROUTINE REGRES FOR BMD02M               AUGUST 17,1966
      DIMENSION X(200,30),XMEAN(45),COV(200,15),VALU(25),SCALE(45),
     1C(200,15),Z(25,25),INDEX(20),YMEAN(20),SQ(20),RES(20)
      COMMON  X      , XMEAN  , COV    , VALU   , SCALE  , C
      COMMON  Z      , INDEX  , YMEAN  , SQ     , RES    , NK
      COMMON  NY     , NV     , N      , NT
      DO 61 I=1,NK
      NY=NV+I
      DO 61 J=1,NV
      C(I,J)=0.0
      DO 61 K=1,N
  61  C(I,J)=C(I,J)+X(K,NY)*COV(K,J)
      DO 63 J=1,NV
      DO 63 I=1,NK
 63   C(I,J)=C(I,J)/VALU(J)
      WRITE (6,936)
      WRITE (6,923)
      WRITE (6,924)
      DO 98 I=1,NK
 98   INDEX(I)=I
      ASSIGN 106 TO NNN
 200  NNK=NK
      KF=1
      KL=0
 99   IF(NNK-6)  97, 101, 101
 97   KL=KL+NNK
      GO TO 102
 101  KL=KL+6
 102  WRITE (6,925)(INDEX(I), I=KF,KL)
      WRITE (6,940)(YMEAN(I), I=KF,KL)
      WRITE (6,941)
      DO 103 J=1,NV
 103  WRITE (6,942)J,(C(I,J), I=KF,KL)
      WRITE (6,936)
      NNK=NNK-6
      IF(NNK) 105, 105, 104
 104  KF=KF+6
      GO TO 99
 105  GO TO NNN, (106,135,145,155,10)
 106  DO 161 J=1,NK
      NY=NV+J
      SQ(J)=0.0
      DO 161 I=1,N
 161  SQ(J)=SQ(J)+(X(I,NY)-YMEAN(J))**2
      DO 165 I=1,NK
      DO 165 J=1,NV
 165  COV(I,J)=C(I,J)
      DO 66 J=1,NV
      DO 66 I=1,NK
  66  C(I,J)=(COV(I,J)*COV(I,J))*VALU(J)
      WRITE (6,927)
      DO 164 I=1,NK
      SUM=0.0
      DO 163 J=1,NV
 163  SUM=SUM+C(I,J)
 164  RES(I)=SQ(I)-SUM
      NNK=NK
      KF=1
      KL=0
 119  IF(NNK-6) 120, 120, 121
 120  KL=KL+NNK
      GO TO 122
 121  KL=KL+6
 122  WRITE (6,925)(INDEX(I), I=KF,KL)
      WRITE (6,928)(RES(I), I=KF,KL)
      WRITE (6,943)(SQ(I), I=KF,KL)
      WRITE (6,944)
      DO 123 J=1,NV
 123  WRITE (6,942)J,(C(I,J), I=KF,KL)
      WRITE (6,936)
      NNK=NNK-6
      IF(NNK) 125, 125, 124
 124  KF=KF+6
      GO TO 119
 125  WRITE (6,930)
      WRITE (6,931)
      WRITE (6,932)
C
C
C
      DO 68 I=1,NK
      DO 68 J=1,NV
  68  C(I,J)=COV(I,1)*Z(J,1)
      ASSIGN 135 TO NNN
      GO TO 200
 135  IF(NV-2)10,70,70
  70  WRITE (6,933)
      DO 71 I=1,NK
      DO 71 J=1,NV
  71  C(I,J)=C(I,J)+COV(I,2)*Z(J,2)
      ASSIGN 145 TO NNN
      GO TO 200
 145  IF(NV-3)10,73,73
  73  WRITE (6,934)
      DO 74 I=1,NK
      DO 74 J=1,NV
  74  C(I,J)=C(I,J)+COV(I,3)*Z(J,3)
      ASSIGN 155 TO NNN
      GO TO 200
 155  IF(NV-4)10,76,76
  76  WRITE (6,935)
      DO 77 II=4,NV
      DO 77 I=1,NK
      DO 77 J=1,NV
  77  C(I,J)=C(I,J)+COV(I,II)*Z(J,II)
      ASSIGN 10 TO NNN
      GO TO 200
 923  FORMAT(54H0COEFFICIENTS OF REGRESSION EQUATIONS USING ORTHOGONAL)
 924  FORMAT(50H COMPONENTS FOR STANDARDIZED INDEPENDENT VARIABLES)
 925  FORMAT(15H0      Y       6I15)
 927  FORMAT(80H0REDUCTION IN SUM OF SQUARES OF THE RESIDUALS DUE TO USI
     1NG ORTHOGONAL COMPONENTS)
 928  FORMAT(20H0RESIDUAL SUM SQ.   6F15.7)
 930  FORMAT(99H0COEFFICIENTS OF REGRESSION EQUATION WHEN ONLY FIRST, FI
     1RST TWO, FIRST THREE AND ALL COMPONENTS ARE)
 931  FORMAT(99H USED AS INDEPENDENT VARIABLES, WHERE EACH COMPONENT IS   
     1EXPRESSED IN TERMS OF THE STANDARDIZED DATA)
 932  FORMAT(26H0FIRST PRINCIPAL COMPONENT)
 933  FORMAT(38H0FIRST AND SECOND PRINCIPAL COMPONENTS)
 934  FORMAT(45H0FIRST, SECOND AND THIRD PRINCIPAL COMPONENTS)
 935  FORMAT(25H0ALL PRINCIPAL COMPONENTS)
 936  FORMAT(1H0)
 940  FORMAT(20H0MEAN (INTERCEPT)   6F15.7)
 941  FORMAT(15H0COEFF. OF VAR.)
 942  FORMAT(1H I7,12H            6F15.7)
 943  FORMAT(20H0TOTAL SUM SQ.      6F15.7)
 944  FORMAT(27H0REDUCTION DUE TO COMPONENT)
  10  RETURN
      END
      SUBROUTINE TPWD(NT1,NT2)
C        SUBROUTINE TPWD FOR BMD02M                   JUNE 20, 1966
      IF(NT1)40,10,12
 10   NT1=5
 12   IF(NT1-NT2)14,19,14
   14 IF(NT2.EQ.5)GO TO 18
   17 REWIND NT2
   19 IF(NT1-5)18,24,18
 18   IF(NT1-6)22,40,22
 22   REWIND NT1
 24   NT2=NT1
 28   RETURN
 40   WRITE (6,49)
      STOP
 49   FORMAT(25H ERROR ON TAPE ASSIGNMENT)
      END
C        SUBROUTINE TRANS FOR BMD02M                  AUGUST 17, 1966
      SUBROUTINE TRANS(DATA,NVAR,NSAM,IERROR,NVG)
      DOUBLE PRECISION C123,TODE
       DIMENSION DATA (200,30)
      ASN(XX)=ATAN(XX/SQRT(1.0-XX**2))
      DATA C123/6HTRNGEN/
      MARY=0
      FN=NSAM
      WRITE (6,1403)
      WRITE (6,1400)
      DO1000 J=1,NVG
      READ (5,1100)TODE,NEWA,LCODE,LVA,BNEW
      IF(C123.EQ.TODE)GO TO 201
  200 NVAR=-NVAR
  201 WRITE (6,1402)J,NEWA,LCODE,LVA,BNEW
      IF(LCODE*(15-LCODE)) 710,710,714
  710 WRITE (6,712)
  712 FORMAT(29HERROR ON TRANSGENERATION CODE)
      GO TO 1000
  714 IF(LCODE-10)4,5,5
    5 NEWB=BNEW
    4 DO 300 I=1,NSAM
      D=DATA(I,LVA)
  202 GOTO(10,20,30,40,50,60,70,80,90,100,110,120,130,140),LCODE
   10 IF(D)99,7,8
    7 D2=0.0
      GO TO 3
    8 D2=SQRT(D)
      GO TO 3
   20 IF(D)99,11,12
   11 D2=1.0
      GO TO 3
   12 D2=SQRT(D)+SQRT(D+1.0)
      GO TO 3
   30 IF(-D)14,99,99
   14 D2=ALOG10(D)
      GO TO 3
   40 D2=EXP(D)
      GO TO 3
   50 IF(-D)17,7,99
   17 IF(D-1.0)18,19,99
   19 D2=3.14159265/2.0
      GO TO 3
   18 D2=ASN(SQRT(D))
      GO TO 3
   60 A=D/(FN+1.0)
      B=A+1.0/(FN+1.0)
      IF(A)99,23,24
   23 IF(-B)27,7,99
   27 D2=ASN(SQRT(B))
      GO TO 3
   24 IF(B)99,28,29
   28 D2=ASN(SQRT(A))
      GO TO 3
   29 D2=ASN(SQRT(A))+ASN(SQRT(B))
      GO TO 3
   70 IF(D)31,99,31
   31 D2=1.0/D
      GO TO 3
   80 D2=D+BNEW
      GO TO 3
   90 D2=D*BNEW
      GO TO 3
  100 IF(-D)33,7,99
   33 D2=D**NEWB
      GO TO 3
  110 D2=D+DATA(I,NEWB)
      GO TO 3
  120 D2=D-DATA(I,NEWB)
      GO TO 3
  130 D2=D*DATA(I,NEWB)
      GO TO 3
  140 IF(DATA(I,NEWB))34,99,34
   34 D2=D/DATA(I,NEWB)
      GO TO 3
   99 IF(MARY)43,44,44
   44 MARY=-999
      IERROR=-999
      WRITE (6,1404)J
   43 WRITE (6,1405)I
      GO TO 300
    3 DATA(I,NEWA)=D2
  300 CONTINUE
      IF(IERROR)42,1000,1000
 1000 CONTINUE
 1100 FORMAT(A6,I3,I2,I3,F6.0)
 1402 FORMAT(2H  I2,I8,2I9,F15.5)
 1403 FORMAT(1H06X,23HTRANS GENERATOR CARD(S))
 1405 FORMAT(10H ITEM NO. I3)
 1400 FORMAT(46H0CARD    NEW     TRANS    ORIG.   ORIG. VAR(B)/45H  NO.   
     1VARIABLE   CODE    VAR(A)   OR CONSTANT)
 1404 FORMAT(30H0THE INSTRUCTIONS INDICATED ON/25H TRANS GENERATOR CARD   
     1NO.I2,4H RE-/29H SULTED IN THE VIOLATION OF A/31H RESTRICTION FOR   
     2THIS TRANSFOR-/31H MATION. THE VIOLATION OCCURRED/28H FOR THE ITEM
     3S LISTED BELOW./)
 1401 FORMAT(41H0PROGRAM CANNOT CONTINUE FOR THIS PROBLEM)
      IF(IERROR)42,1111,1111
   42 WRITE (6,1401)
 1111 RETURN
      END