Google
 

Trailing-Edge - PDP-10 Archives - decus_20tap5_198111 - decus/20-0137/bmd/bmd01s.for
There is 1 other file named bmd01s.for in the archive. Click here to see a list.
C             LIFE TABLE AND SURVIVAL RATES            MAY 20, 1966
C        THIS IS A SIFTED VERSION OF BMD01S ORIGINALLY WRITTEN IN
C        FORTRAN II. SOME MODIFICATIONS WERE MADE TO MAKE IT OPERABLE
C        AND SLIGHTLY MORE EFFICIENT THAN THE SIFTED VERSION.
C        IT WAS THEN CONVERTED TO 360 FORTRAN IV (H-LEVEL)
      DIMENSION DATA(31,31,5),C(31,12),N(5),KSR(31),IDATE(12),DCODE(10),
     1FCODE(10),IORD(25),TD(25),         C1(31),C3(31),COH(31),PD(31),
     2X(8),PROB(15,31),Y(15)
      COMMON  C      , K      , SUM    , STD    , ITABLE , N
      COMMON  SRATE  , DATA   , KSR    , LK     , NTOT   , NV
      COMMON  NJ     , FORM   , JCO    , IPD    , UNIT   , NTAPE
      COMMON  MAT    , C1     , C3     , DCODE  , FCODE  , IORD
      COMMON  IDATE  , TD     , X      , PROB   , MN     , NKI
      DATA BLANK/' '/
      INTEGER SYM(15)
      DOUBLE PRECISION A1,A2,A3,A4
      DATA A2/'PROBLM'/,A3/'SUBPRO'/,A4/'FINISH'/
      COMMON  MTAPE  , IFLAG
C     SCRATCH TAPE 1 USED IN PTDATA TO CONVERT INTEGER TO ALPHA
      REWIND 1
      MTAPE=5
	CALL USAGEB('BMD01S')
C
  902 FORMAT(52H1BMD01S - LIFE TABLE AND SURVIVAL RATE - VERSION OF 
     118HMAY 20, 1966      ,/
     241H HEALTH SCIENCES COMPUTING FACILITY, UCLA//)
C
      IDATE( 1)=11
      IDATE( 2)=21
      IDATE( 3)=31
      IDATE( 4)=41
      IDATE( 5)=12
      IDATE( 6)=22
      IDATE( 7)=32
      IDATE( 8)=42
      IDATE( 9)=13
      IDATE(10)=23
      IDATE(11)=33
      IDATE(12)=43
   10 READ (5,910)A1,PRON,NSUB,NPLOT,LTAPE
      IF(A1.EQ.A2) GO TO 12
    5 IF(A1.EQ.A4)GO TO 40
    7 WRITE (6,915) A1
      GO TO 40
   12 CALL TPWD(LTAPE,MTAPE)
      IFLAG=0
      IF(NSUB*(NSUB-11))3,90,90
   90 WRITE(6,300)NSUB
  300 FORMAT('0NO. OF SUBPROBLEM CARDS SPECIFIED WAS'
     9 ,I3,'.  ILLEGAL.'  )
      GO TO 40
 3    DO 100 LP=1,NSUB
      MN=LP
      WRITE (6,902)
      WRITE (6,903)PRON,MN
      READ (5,900)A1,LK,NKI,KCO,FORM,JCO,NTOT,NV,IPD,UNIT, SYM(MN),IPRIN
     1T,MAT
      IF(A1.EQ.A3) GO TO 11
    8 WRITE (6,916) A1
      GO TO 40
   11 ITABLE=1
      IF(LK*(LK-32))110,91,91
  110 IF(NKI-LK)111,92,92
   91 WRITE(6,301) LK
      GO TO 40
   92 WRITE(6,302) NKI,LK
      GO TO 40
  301 FORMAT('0MAXIMUM NO. OF PERIODS AFTER DIAGNOSIS IS',I3,', WHICH IS
     9 GREATER THAN 31.')
  302 FORMAT( '0MAXIMUM NO. OF PERIODS FOR THE LIFE TABLE COMPUTATIONS W
     9AS',I3/   ' WHICH IS NOT LESS THAN THE MAXIMUM NO. OF PERIODS AFTE
     9R DIAGNOSIS,',I3)
  111 LNG=0
      IF(FORM.EQ.BLANK) GO TO 14
   13 WRITE (6,918)FORM
      IF(MAT.GT.0.AND.MAT.LE.10) GO TO 200
      WRITE(6,4000)
      MAT=1
  200 DO 16 I=1,31
      DO 16 J=1,31
      DO 16 KK=1,5
   16 DATA(I,J,KK)=0.0
      CALL READ
      IF(NV*(NV-10))41,93,93
   93 WRITE(6,303)NV
      GO TO 40
  303 FORMAT('0NO. OF VARIABLES,',I3,', INCORRECT')
   14 WRITE (6,917)
      IF(JCO)81,81,80
   80 MX=LK
      GO TO 82
   81 MX=1
   82 DO 15 II=1,MX
      KK=LK-(II-1)
      DO 15 I=1,KK
   15 READ (5,901)(DATA(II,I,J),J=1,5)
   41 IF(IPRINT)75,75,70
   70 WRITE (6,904)ITABLE
      ITABLE=ITABLE+1
      WRITE (6,905)
      WRITE (6,906)
      IF(JCO)72,72,71
   71 MX=LK
      GO TO 73
   72 MX=1
   73 DO 18 II=1,MX
      KK=LK-(II-1)
      WRITE (6,907)II
      N1=0
      N2=1
      DO 18 I=1,KK
      DO 17 J=2,5
   17 N(J)=DATA(II,I,J)
      WRITE (6,908)N1,N2,(N(J),J=2,5)
      N1=N1+1
   18 N2=N2+1
   75 CONTINUE
      IF(KCO) 30, 30, 31
   30 NK=1
      KSR(1)=NKI
      GO TO 19
   31 DO 32 I=1,NKI
      KB=NKI-(I-1)
   32 KSR(I)=KB
      NK=NKI
   19 DO 23 L=1,NK
      K=KSR(L)
      DO 20 I=1,LK
      DO 20 J=1,12
   20 C(I,J)=0.0
      IF(JCO)76,76,77
   76 MX=1
      GO TO 78
   77 MX=LK
   78 DO 22 II=1,MX
      KK=LK-(II-1)
      DO 22 I=1,KK
      DO 22 J=2,5
   22 C(I,J)=C(I,J)+DATA(II,I,J)
      CALL RATE
      CALL REPORT
      IF(JCO)23,23,42
   42 CALL COHORT
   23 CONTINUE
  100 CONTINUE
      IF(NPLOT)46,46,45
   45 WRITE (6,911)
      XMIN=0.0
      XMAX=NKI
      DO 60 I=1,NSUB
   60 Y(I)= 1.0
      X(1)=0.0
      CALL PLOTR(X,XMIN,XMAX,Y,SYM,0.0,1.0,NSUB,-1)
      XMIN=1.0
      DO 55 I=1,NKI
      X(1)=I
      DO 50 J=1,NSUB
 50   Y(J)=-PROB(J,I)
   55 CALL PLOTR(X,XMIN,XMAX,Y,SYM,0.0,1.0,NSUB,-1)
      CALL PLOTR(X,XMIN,XMAX,Y,SYM,0.0,1.0, -1 ,-1)
   46 GO TO 10
   40 IF(MTAPE-5)741,742,741
  741 REWIND MTAPE
  742 STOP
 900  FORMAT(A6,3I2,A1,I2,I4,I2,I4,A1,A2,4X,I2,36X,I2)
  901 FORMAT(5F4.0)
  903 FORMAT(12H PROBLEM NO.6X,A2/16H SUB-PROBLEM NO.I4)
  904 FORMAT(1H028X,5HTABLEI2,16H  SURVIVAL  DATA)
  905 FORMAT(85H0PERIODS AFTER    ALIVE AT BEGIN-   DIED DURING   LOST T
     1O FOLLOW-UP   WITHDRAWN ALIVE)
  906 FORMAT(85H   DIAGNOSIS     NING OF INTERVAL    INTERVAL      DURIN
     1G INTERVAL    DURING INTERVAL)
  907 FORMAT(6H0GROUPI3,20H  PATIENTS DIAGNOSED)
  908 FORMAT(I6,2H -I2,I16,I17,I16,I19)
  910 FORMAT(A6,A2,2I2,56X,I2)
  911 FORMAT(1H146X,19HSURVIVAL PROPORTION//)
  915 FORMAT (' PROBLM OR FINISH CARD EXPECTED,BUT READ ',A6)
  916 FORMAT (' SUBPRO CARD EXPECTED, BUT READ ',A6)
  917 FORMAT(21H0SURVIVAL TABLE INPUT//)
  918 FORMAT(7H0FORM  A1,7H  INPUT//)
 4000 FORMAT(1H0,23X,71HNUMBER OF VARIABLE FORMAT CARDS INCORRECTLY SPEC
     XIFIED, ASSUMED TO BE 1.)
      END
      SUBROUTINE AFORM(II,DIV,COH)
C             SUBROUTINE AFORM FOR BMD01S              MAY 20, 1966
      DIMENSION DATA(31,31,5),C(31,12),N(5),KSR(31),IDATE(12),DCODE(10),
     1FCODE(10),IORD(25),TD(25),         C1(31),C3(31),COH(31),
     2X(8),PROB(15,31)
      COMMON  C      , K      , SUM    , STD    , ITABLE , N
      COMMON  SRATE  , DATA   , KSR    , LK     , NTOT   , NV
      COMMON  NJ     , FORM   , JCO    , IPD    , UNIT   , NTAPE
      COMMON  MAT    , C1     , C3     , DCODE  , FCODE  , IORD
      COMMON  IDATE  , TD     , X      , PROB   , MN     , NKI
      IX=II-1
      DATA A/1HA/
      IF(FORM.NE.A)GO TO 45
    8 IF(X(2))15,15,10
   10 KK=2
      GO TO 20
   15 IF(X( 3))80,80,16
   16 KK=3
   20 IF(JCO)25,25,30
   25 X(KK)=(X(KK)-X(1))/DIV
   30 DO 35 I=1,LK
      L=I
      IF(X(KK)-COH(I))40,40,35
   35 CONTINUE
   40 KK=KK+1
      L=L-IX
      GO TO 75
   45 IF(X(6))55,55,50
   50 KK=6
      GO TO 60
   55 IF(X(7))80,80,56
   56 KK=7
   60 DO 65 I=1,LK
      L=I
      IF(X(KK)-COH(I))70,70,65
   65 CONTINUE
   70 KK=KK-3
   75 DATA(II,L,KK)=DATA(II,L,KK)+1.0
   80 RETURN
      END
      SUBROUTINE BFORM (II,DIV,COH)
C             SUBROUTINE BFORM FOR BMD01S              MAY 20, 1966
      DIMENSION DATA(31,31,5),C(31,12),N(5),KSR(31),IDATE(12),DCODE(10),
     1FCODE(10),IORD(25),TD(25),         C1(31),C3(31),COH(31),
     2X(8),PROB(15,31)
      COMMON  C      , K      , SUM    , STD    , ITABLE , N
      COMMON  SRATE  , DATA   , KSR    , LK     , NTOT   , NV
      COMMON  NJ     , FORM   , JCO    , IPD    , UNIT   , NTAPE
      COMMON  MAT    , C1     , C3     , DCODE  , FCODE  , IORD
      COMMON  IDATE  , TD     , X      , PROB   , MN     , NKI
      IX=II-1
      DATA B/1HB/
      IF(FORM.NE.B)GO TO 15
   10 MM=4
      IF(JCO)11,11,20
   11 X(4)=(X(4)-X(1))/DIV
      GO TO 20
   15 MM=8
   20 DO 25 I=1,LK
      JJ=I
      IF(X(MM)-COH(I))30,30,25
   25 CONTINUE
   30 IF(MM-8) 32, 34, 32
   32 JJ=JJ-IX
   34 DO 40 I=1,10
      IF(X(5)-DCODE(I))35,45,35
   35 IF(X(5)-FCODE(I))40,50,40
   40 CONTINUE
   45 DATA(II,JJ,3)=DATA(II,JJ,3)+1.0
      GO TO 55
   50 DATA(II,JJ,4)=DATA(II,JJ,4)+1.0
   55 RETURN
      END
      SUBROUTINE COHORT
C             SUBROUTINE COHORT FOR BMD01S             MAY 20, 1966
      DIMENSION DATA(31,31,5),C(31,12),N(5),KSR(31),IDATE(12),DCODE(10),
     1FCODE(10),IORD(25),TD(25),FMT(120),C1(31),C3(31),COH(31),PD(31),
     2X(8),PROB(15,31)
      COMMON  C      , K      , SUM    , STD    , ITABLE , N
      COMMON  SRATE  , DATA   , KSR    , LK     , NTOT   , NV
      COMMON  NJ     , FORM   , JCO    , IPD    , UNIT   , NTAPE
      COMMON  MAT    , C1     , C3     , DCODE  , FCODE  , IORD
      COMMON  IDATE  , TD     , X      , PROB   , MN     , NKI
      PRED=0.0
      PINC=0.0
      DO 10 I=1,LK
      DO 10 J=1,12
 10   C(I,J)=0.0
      K1=DATA(1,1,1)
      WRITE (6,920)
      WRITE (6,921)ITABLE,K
      ITABLE=ITABLE+1
      WRITE (6,922)K
      WRITE (6,923)K
      WRITE (6,924)
      WRITE (6,926)
      DO 21 I=1,LK
      KK=LK-(I-1)
      DO 15 II=1,KK
      DO 15 J=2,5
 15   C(II,J)=C(II,J)+DATA(I,II,J)
      CALL RATE
      IF(PRED)18,19,18
 18   PINC=((PRED-STD)/PRED)*100.0
      GO TO 20
 19   PRED=STD
 20   K2=DATA(1,I,1)
      NCASE=C(1,2)
 21   WRITE (6,925)K1,K2,NCASE,SRATE,STD,PINC
 920  FORMAT(1H0)
 921  FORMAT(6H0TABLEI3,I5,76H-PERIOD SURVIVAL RATES AND THEIR STANDARD 
     1ERROR WITH INCREASE IN COHORT SIZE)
 922  FORMAT(23H0             NUMBER OFI8, 47H-PERIOD      STANDARD ERRO
     1R     PERCENT RED. IN)
 923  FORMAT(47H   COHORT       CASES         SURVIVAL       OFI3,27H-PE
     1RIOD       STD. ERROR OF)
 924  FORMAT(77H              DIAGNOSED         RATE        SURVIVAL RAT
     1E       SURVIVAL RATE)
 925  FORMAT(I6,1H-I6, I7,F16.2,F17.2,F20.1)
 926  FORMAT(1H )
      RETURN
      END
C             SUBROUTINE PLOTR (IBM 360)              AUGUST 13, 1966
C
C     'PLOTR' IS A UTILITY SUBPROGRAM FOR THE BMD... PROGRAMS WHICH
C     PLOTS EITHER SINGLE-LINE OR WHOLE-PAGE GRAPHS AND SETS UP
C     APPROPRIATE SCALING.  THE CALLING PARAMETERS ARE AS FOLLOWS -
C
C     X - THE VALUE OF THE INDEPENDENT VARIABLE
C     ZMIN - THE MINIMUM VALUE OF X FOR THIS PLOT
C     ZMAX - THE MAXIMUM VALUE OF X FOR THIS PLOT
C     Y - THE ARRAY CONTAINING THE VALUES OF UP TO 15 DEPENDENT VAR.'S
C     SYM - THE ARRAY CONTAINING THE SYMBOLS TO BE PLOTTED
C     WMIN - THE MINIMUM VALUE OF ALL Y'S FOR THIS PLOT
C     WMAX - THE MAXIMUM VALUE OF ALL Y'S FOR THIS PLOT
C     NC - THE NUMBER OF DEPENDENT VARIABLES
C               NC=-1 CLOSES A SINGLE-LINE PLOT
C               NC= 0 PRINTS AND CLOSES A WHOLE-PAGE PLOT
C     NP - THE CONTROL VARIABLE
C               NP=-1 PRINTS A SINGLE LINE
C               NP=0 OR NP=1 SETS UP A WHOLE-PAGE PLOT
C
C     THE PLOTTING ROUTINE MUST BE CALLED ONCE FOR EACH VALUE OF THE
C     INDEPENDENT VARIABLE THAT IS TO BE PLOTTED NO MATTER WHETHER IN
C     THE SINGLE-LINE OR WHOLE-PAGE MODE
C
      SUBROUTINE PLOTR(X,ZMIN,ZMAX,Y,SYM,WMIN,WMAX,NC,NP)
C
      DIMENSION Y(15),CLAB(12),GF(10),FMT(12),XY(51,101),SYM(15)
      DATA TC,TP/1H.,1H+/
      DATA GF/            4H 1X,,4H 2X,,4H 3X,,4H 4X,,4H 5X,,4H 6X,,
     14H 7X,,4H 8X,,4H 9X,,4H 10X/
      DATA FMT/'(17X',' ','5(F1','2.3,','8X)/','7X, ',' ','4(F1','2.3,',
     1'8X),','F12.','3)  '/
      INTEGER XY,BLANKS
      DATA BLANKS/'  '/
C
 100   FORMAT(1H 6X5(F12.3,8X),F12.3/17X,5(F12.3,8X))
 101  FORMAT(1H F12.3,1X,103A1,F12.3)
  102 FORMAT (1H 13X,103A1)
 1000 FORMAT(1H  14X,101A1)
 1001 FORMAT(15X,20(5H+....),1H+)
C
      DATA NCC/2/
C    'NCC' ON THE INITIAL ENTRY TO PLOTR IS NON-ZERO BECAUSE OF THE DATA
C    STATEMENT ABOVE.
C
C    'NCC' IS 0 WHILE A PLOT IS BEING MADE.  IT IS 1 OR 2 AT OTHER TIMES
C
      IF(NCC) 50,48,50
C
C    THE VARIABLE 'KL' CONTROLS THE FUNCTIONING OF THE OPENING AND
C    CLOSING  SECTIONS OF PLOTR.  KL=0 INDICATES OPENING OF THE GRAPH,
C    KL=1 INDICATES CLOSING.
C
   50 KL=0
      CALL SCALE(WMIN,WMAX,100.0,JY,YMIN,YMAX,YIJ)
      YR=YMAX-YMIN
  230 J=JY
      IF(J*(J-10))204,201,201
C
C     THE FOLLOWING SECTION OPENS OR CLOSES A PLOT IN FIXED FORMAT
C     UNDER CONTROL OF KL
C
  201 IF(KL)220,220,231
C
  231 WRITE (6,1001)
      IF(KL)250,250,220
C
  220 CLAB(1)= YMIN
      DO 222 I=2,11
  222 CLAB(I)=CLAB(I-1)+YIJ
      WRITE (6,100)(CLAB(I),I=1,11,2),(CLAB(J),J=2,10,2)
      IF(KL)231,231,14
C
C     THE FOLLOWING SECTION OPENS OR CLOSES A PLOT IN A VARIABLE
C     FORMAT UNDER CONTROL OF KL AND JY FROM 'SCALE'
C
  204 IF(J-5)205,221,207
  207 J=J-5
  205 JYT=5-J
  221 CONTINUE
  226 FMT(2)=GF(JY)
      IF (KL) 225,225,227
C
  225 FMT(7)=GF(JY)
      TT=JY
      TT=TT*YIJ/10.0
      CLAB(1)= YMIN+TT
      DO 223 I=2,10
  223 CLAB(I)=CLAB(I-1) +YIJ
      WRITE (6,FMT)(CLAB(I),I=2,10,2),(CLAB(I),I=1,9 ,2)
      IF(KL)227,227,14
C
  227 IF(JY-5)208,209,208
  208 WRITE(6,1000)(TC,I=1,J    ),(TP,(TC,I=1,4),K=1,19),TP,(TC,I=1,JYT)
      IF (KL) 250,250,225
C
  209 WRITE (6,1001)
      IF (KL) 250,250,225
C
  250 CONTINUE
      NCC=0
      IC=0
      IF(NP)80,11,11
C
C    THIS SECTION PREPARES FOR A FULL PAGE PLOT BY FILLING IN XY WITH
C     BLANKS AND SETTING UP SCALING FOR THE INDEPENDENT VARIABLE 'X'
C
   11 DO 1 I=1,51
       DO 1 J=1,101
   1  XY(I,J)=BLANKS
      CALL SCALE (ZMIN,ZMAX,50.0,JX,XMIN,XMAX,XIJ)
      XR=XMAX-XMIN
      GO TO 48
C
C
C     ENTRY TO PLOTS CAN BE USED ONLY AFTER THE CALLING PARAMETERS
C     HAVE BEEN TRANSFERRED BY A CALL ON PLOTR.  THE CALL ON PLOTS
C     IS IDENTICAL WITH ENTRY TO PLOTR BUT IT ALLOWS THE PROGRAMMER TO
C     CALL THE PLOTTING ROUTINE WITHOUT HAVING TO INCLUDE THE PARAMETERS
C
C      ENTRY PLOTS
   48 IF(NC)52,13,49
   49 IF(NP)80,10,10
C    THE FOLLOWING SECTION SETS UP A FULL PAGE BUT DOES NO PRINTING.
C    THIS SECTION IS REACHED BY SPECIFYING NC POSITIVE AND NP POSITIVE.
C
   10 DO 9 N=1,NC
      SYMB=SYM(N)
      XDIFFR=XMAX-X
      IF(XDIFFR)105,106,106
  105 XDIFFR=0.0
  106 YDIFFR=YMAX-Y(N)
      IF(YDIFFR)107,108,108
  107 YDIFFR=0.0
  108 L=51.0-(50.0*XDIFFR)/XR+.5
      K=101.0-(100.0*YDIFFR)/YR+.5
      CALL FORM2(SYMB,XY(L,K))
 9     CONTINUE
      GO TO 15
C
C     THE FOLLOWING SECTION CONSTRUCTS AND PLOTS ONE LINE OF A MULTILINE
C     GRAPH.  LOCATION ALONG THE AXIS OF THE PAPER IS PRINTED AT EVERY
C     STEP.  THIS SECTION IS ACCESSIBLE BY SPECIFYING NC POSITIVE AND
C     NP NEGATIVE.
C
   80 DO 86 I=1,101
   86 XY(1,I)=BLANKS
       L=1
      DO 95 N=1,NC
      SYMB=SYM(N)
      YDIFFR=YMAX-Y(N)
      IF(YDIFFR)860,865,865
  860 YDIFFR=0.0
  865 K=101.0-(100.0*YDIFFR)/YR+.5
   95 CALL FORM2(SYMB,XY(L,K))
      IF(MOD(IC,5))97,96,97
   96 W=TP
      GO TO 98
   97 W=TC
   98 WRITE (6,101)X,W,(XY(1,N),N=1,101),W,X
      IC=IC+1
      GO TO 15
C
C    THIS SECTION PLOTS OUT THE PREVIOUSLY PREPARED WHOLE PAGE GRAPH.
C    IT PRINTS LOCATION ALONG THE PAPER'S AXIS EVERY FIFTH STEP.  THIS
C    SECTION IS ACCESSED BY SPECIFYING NC=0.
C
   13 M=6-JX
      LL=50+M
      T=JX
      IF(5-JX)131,131,135
  131 T=0.0
  135 RLAB=XMAX-(T*XIJ)/5.0
      W=TC
      K=52
      DO 31 L=M,LL
      K=K-1
      I=MOD(L,5)
      IF(I-1)2,3,2
    3 W=TP
      WRITE (6,101)RLAB,W,(XY(K,N),N=1,101),W,RLAB
      RLAB=RLAB-XIJ
      W=TC
      GO TO 31
    2 WRITE (6,102)W,(XY(K,N),N=1,101),W
   31 CONTINUE
C
   52 KL=1
      GO TO 230
C
   14 NCC=1
   15 RETURN
      END
C             SUBROUTINE SCALE FOR PLOTR              JUNE 21, 1966
C
C     SUBROUTINE 'SCALE' CALCULATES THE SCALING FOR 'PLOTR'
C
      SUBROUTINE SCALE(YMIN,YMAX,YINT,JY,TYMIN,TYMAX,YIJ)
      DIMENSION C(10)
      DATA C           /1.0,1.5,2.0,3.0,4.0,5.0,7.5,10.0,15.0,20.0/
      TEST=1.0/(2.0**17)
   50 YR=YMAX-YMIN
      TT=YR/YINT
      J = ALOG10(TT)+TEST
      E=10.0**J
      TT=TT/E
      I=0
      IF(TT-1.0+TEST)205,201,201
  205 TT=TT*10.0
      E=E/10.0
 201  I=I+1
      IF(9-I)1,2,2
    1 E=E*10.0
      I=1
    2 IF(TT-C(I))233,202,201
  233 YIJ=C(I)*E
      GO TO 203
  202 Y=YMIN/C(I)
      J=Y
      T=J
      IF(0.0001-ABS(T-Y))204,233,233
  204 YIJ=C(I+1)*E
  203 X=((YMAX+YMIN)/YIJ-YINT )/2.0+.00001
      K=X
      IF(K)235,240,240
  235 Y=K
      IF(X-Y)236,240,236
  236 K=K-1
  240 TYMIN=K
      TYMIN=YIJ*TYMIN
      TYMAX=TYMIN+YINT*YIJ
      IF (YMAX-TYMAX-TEST)11,11,201
   11 YIJJ=C(I)*E
      XT=((YMAX+YMIN)/YIJJ-YINT)/2.0+.00001
      KT=XT
      IF (KT) 1235,1240,1240
 1235 YT=KT
      IF (XT.NE.YT) KT=KT-1
 1240 TYMINT=KT
      TYMINT=YIJJ*TYMINT
      TYMAXT=TYMINT+YINT*YIJJ
      IF (YMAX-TYMAXT.GT.TEST) GO TO 10
      TYMIN=TYMINT
      TYMAX=TYMAXT
      YIJ=YIJJ
      K=KT
   10 TT=YINT/10.0
      JY=TT+.000001
      YIJ=YINT*(YIJ/10.0)
      J=TYMIN/ YIJ
      IF (K)242,241,241
  242 J=J-1
  241 J=J*JY+JY-K
      JY=J
      RETURN
      END
      SUBROUTINE FORM2(SYMB,XY)
C        SUBROUTINE FORM2 FOR PLOTR (IBM 360)         JUNE 21, 1966
      INTEGER XY,SYMB,BLANK,TEST(18)
      DATA BLANK/'  '/,TEST/'2 ','3 ','4 ','5 ','6 ',
     1'7 ','8 ','9 ','A ','B ','C ','D ','E ','F ','G ','H ','I ','/ '/
      IF(XY.EQ.BLANK)GO TO 50
      DO 30 I=1,17
      IF(XY.NE.TEST(I))GO TO 30
C        PUT IN NEXT SYMBOL OF ARRAY FOR MULTIPLE POINTS
      XY=TEST(I+1)
      GO TO 100
   30 CONTINUE
      IF(XY.EQ.TEST(18))GO TO 100
C        IF OTHER THAN CHARACTERS IN ARRAY TEST PUT IN CHARACTER 2.
      XY=TEST(1)
      GO TO 100
C        IF BLANK, PUT IN SYMBOL
   50 XY=SYMB
  100 RETURN
      END
       SUBROUTINE PTDATA
C            SUBROUTINE PTDATA FOR BMD01S              MAY 20, 1966
C     THE VARIABLE NAMED 'N' WAS CHANGED TO 'M' IN THE DIMENSION AND
C     COMMON STATEMENTS SINCE  A DIMENSIONED VARIABLE CANNOT BE USED
C     AS DO INDEX.
      DIMENSION DATA(31,31,5),C(31,12),M(5),KSR(31),IDATE(12),DCODE(10),
     1FCODE(10),IORD(25),TD(25),FMT(120),C1(31),C3(31),COH(31),PD(31),
     2X(8),PROB(15,31),IZ(9),AZ(9)
      COMMON  C      , K      , SUM    , STD    , ITABLE , M
      COMMON  SRATE  , DATA   , KSR    , LK     , NTOT   , NV
      COMMON  NJ     , FORM   , JCO    , IPD    , UNIT   , NTAPE
      COMMON  MAT    , C1     , C3     , DCODE  , FCODE  , IORD
      COMMON  IDATE  , TD     , X      , PROB   , MN     , NKI
      COMMON  MTAPE  , IFLAG
      DATA AMZ,BLANK,A,B,C2/2H-0,1H ,1HA,1HB,1HC/
      EQUIVALENCE(AZ,IZ)
      IF(FORM.NE.A)GO TO 39
 10   IF(IFLAG)9,8,9
 8    WRITE (6,900)
      WRITE (6,901)
      WRITE (6,907)
      IFLAG=1
 9    DO 38 N=1,NV
      IF(IORD(N)-12)11,12,14
 11   IZ(1)=TD(N)
      GO TO 38
 12   IZ(2)=TD(N)
      GO TO 38
 13   IZ(3)=TD(N)
      GO TO 38
 14   IF(IORD(N)-21)13,21,24
 21   IZ(4)=TD(N)
      GO TO 38
 22   IZ(5)=TD(N)
      GO TO 38
 23   IZ(6)=TD(N)
      GO TO 38
 24   IF(IORD(N)-23)22,23,34
 31   IZ(7)=TD(N)
      GO TO 38
 32   IZ(8)=TD(N)
      GO TO 38
 33   IZ(9)=TD(N)
      GO TO 38
 34   IF(IORD(N)-32)31,32,33
 38   CONTINUE
      WRITE(1,910) (IZ(N),N=1,NV)
      REWIND 1
      READ (1,911) (IZ(N),N=1,NV)
      REWIND 1
      DO 35 N=1,NV
      IF(AZ(N).EQ.AMZ) GO TO 135
      GO TO 35
 135  AZ(N)=BLANK
 35   CONTINUE
      WRITE (6,908)(IZ(N),N=1,NV)
      RETURN
   39 IF(FORM.NE.B)GO TO 58
 40   IF(IFLAG)37,36,37
 36   WRITE (6,902)
      WRITE (6,903)
      WRITE (6,907)
      IFLAG=1
 37   DO 57 N=1,NV
      IF(IORD(N)-12)45,46,47
 41   IZ(4)=TD(N)
      GO TO 57
 42   IZ(5)=TD(N)
      GO TO 57
 43   IZ(6)=TD(N)
      GO TO 57
 44   IZ(3)=TD(N)
      GO TO 57
 45   IZ(1)=TD(N)
      GO TO 57
 46   IZ(2)=TD(N)
      GO TO 57
 47   IF(IORD(N)-41)44,41,48
 48   IF(IORD(N)-43)42,43,50
 50   IZ(7)=TD(N)
 57   CONTINUE
      WRITE(1,910) (IZ(N),N=1,NV)
      REWIND 1
      READ (1,911) (IZ(N),N=1,NV)
      REWIND 1
      DO 55 N=1,NV
      IF(AZ(N).EQ.AMZ) GO TO 155
      GO TO 55
 155  AZ(N)=BLANK
 55   CONTINUE
      WRITE (6,908)(IZ(N),N=1,NV)
      RETURN
   58 IF(FORM.NE.C2)GO TO 75
 59   IF(IFLAG)66,67,66
 67   WRITE (6,904)
      WRITE (6,905)
      WRITE (6,907)
      IFLAG=1
 66   DO 74 N=1,NV
      IF(IORD(N)-12)60,61,62
 60   IZ(1)=TD(N)
      GO TO 74
 61   IZ(2)=TD(N)
      GO TO 74
 62   IF(IORD(N)-60)63,64,65
 63   IZ(3)=TD(N)
      GO TO 74
 64   IZ(4)=TD(N)
      GO TO 74
 65   IZ(5)=TD(N)
 74   CONTINUE
      WRITE(1,910) (IZ(N),N=1,NV)
      REWIND 1
      READ (1,911) (IZ(N),N=1,NV)
      REWIND 1
      DO 84 N=1,NV
      IF(AZ(N).EQ.AMZ) GO TO 184
      GO TO 84
 184  AZ(N)=BLANK
 84   CONTINUE
      WRITE (6,908)(IZ(N),N=1,NV)
      RETURN
 75   IF(IFLAG)83,82,83
 82   WRITE (6,906)
      WRITE (6,905)
      WRITE (6,907)
      IFLAG=1
 83   DO 95 N=1,NV
      IF(IORD(N)-12)76,77,78
 76   IZ(1)=TD(N)
      GO TO 95
 77   IZ(2)=TD(N)
      GO TO 95
 78   IF(IORD(N)-50)79,80,81
 79   IZ(3)=TD(N)
      GO TO 95
 80   IZ(4)=TD(N)
      GO TO 95
 81   IZ(5)=TD(N)
 95   CONTINUE
      WRITE(1,910) (IZ(N),N=1,NV)
      REWIND 1
      READ (1,911) (IZ(N),N=1,NV)
      REWIND 1
      DO 96 N=1,NV
      IF(AZ(N).EQ.AMZ) GO TO 196
      GO TO 96
 196  AZ(N)=BLANK
 96   CONTINUE
      WRITE (6,908)(IZ(N),N=1,NV)
      RETURN
 900  FORMAT(1H012X17HDATE OF DIAGNOSIS12X13HDATE OF DEATH13X18HDATE LAS
     1T OBSERVED)
 901  FORMAT(1H03(12X16HMONTH  DAY  YEAR))
 902  FORMAT(1H012X17HDATE OF DIAGNOSIS9X18HDATE LAST OBSERVED18X6HSTATU
     1S/40X15H(ALIVE OR DEAD)16X15H(ALIVE OR DEAD))
  903 FORMAT(1H02(12X16HMONTH  DAY  YEAR),19X4HCODE)
 904  FORMAT(1H012X17HDATE OF DIAGNOSIS9X19HLENGTH OF TIME FROM11X19HLEN
     1GTH OF TIME FROM/39X18HDIAGNOSIS TO DEATH12X17HDIAGNOSIS TO LAST/5
     22X33HOBSERVATION FOR THOSE STILL ALIVE)
 905  FORMAT(1H012X16HMONTH  DAY  YEAR2(19X4HCODE))
 906  FORMAT(1H012X17HDATE OF DIAGNOSIS16X6HSTATUS17X19HLENGTH OF TIME F
     1ROM/40X15H(ALIVE OR DEAD)15X17HDIAGNOSIS TO LAST/56X27HOBSERVATION
     2 (ALIVE OR DEAD))
 907  FORMAT(1H )
 908  FORMAT(3(15XA2,4XA2,3XA2))
 910  FORMAT(9I2)
 911  FORMAT(9A2)
      END
      SUBROUTINE RATE
C             SUBROUTINE RATE FOR BMD01S               MAY 20, 1966
      DIMENSION DATA(31,31,5),C(31,12),N(5),KSR(31),IDATE(12),DCODE(10),
     1FCODE(10),IORD(25),TD(25),FMT(120),C1(31),C3(31),COH(31),PD(31),
     2X(8),PROB(15,31)
      COMMON  C      , K      , SUM    , STD    , ITABLE , N
      COMMON  SRATE  , DATA   , KSR    , LK     , NTOT   , NV
      COMMON  NJ     , FORM   , JCO    , IPD    , UNIT   , NTAPE
      COMMON  MAT    , C1     , C3     , DCODE  , FCODE  , IORD
      COMMON  IDATE  , TD     , X      , PROB   , MN     , NKI
      DO 23 I=1,K
 23   C(I,6)=C(I,2)-0.5*(C(I,4)+C(I,5))
      DO 24 I=1,K
      C(I,7)=C(I,3)/C(I,6)
 24   C(I,8)=1.0-C(I,7)
      DO 25 I=1,K
      C(I,9)=1.0
      DO 25 J=1,I
 25   C(I,9)=C(I,9)*C(J,8)
      SUM=0.0
      DO 26 I=1,K
      C(I,10)=C(I,6)-C(I,3)
      C(I,11)=C(I,7)/C(I,10)
      SUM=SUM+C(I,11)
   26 C(I,12)=C(I,9)*SQRT(SUM)
      SRATE=C(K,9)
      STD=C(K,12)
      RETURN
      END
      SUBROUTINE READ
C             SUBROUTINE READ FOR BMD01S               MAY 20, 1966
      DIMENSION DATA(31,31,5),C(31,12),N(5),KSR(31),IDATE(12),DCODE(10),
     1FCODE(10),IORD(25),TD(25),FMT(180),C1(31),C3(31),COH(31),PD(31),
     2X(8),PROB(15,31)
      COMMON  C      , K      , SUM    , STD    , ITABLE , N
      COMMON  SRATE  , DATA   , KSR    , LK     , NTOT   , NV
      COMMON  NJ     , FORM   , JCO    , IPD    , UNIT   , NTAPE
      COMMON  MAT    , C1     , C3     , DCODE  , FCODE  , IORD
      COMMON  IDATE  , TD     , X      , PROB   , MN     , NKI
      COMMON  MTAPE
	DOUBLE PRECISION A1,A2,A3,A4
      DATA A2/'COHORT'/,A3/'ORDER '/,A4/'STATUS'/
      DATA B/'B'/,D/'D'/,W/'W'/,XM/'M'/,Y/'Y'/,A/'A'/,C2/'C'/
   21 IF(JCO)30,30,25
   25 READ (5,911)A1,(C1(I),COH(I),C3(I), I=1,LK)
      IF(A1.EQ.A2) GO TO 18
   15 WRITE (6,900) A1
      STOP
   18 DO 26 I=1,LK
      C4   =C3(I)*10000.0
      C5   =C1(I)*100.0
      C6=COH(I)*100.0
      C7=C1(I)*10000.0
      DATA(1,I,1)=C6+C7+C3(I)
   26 COH(I)=C4+C5+COH(I)
   30 IF(FORM.EQ.B)GO TO 32
   31 IF(FORM.NE.D)GO TO 35
   32 READ (5,912)A1,(DCODE(I),I=1,10),(FCODE(I),I=1,10)
      IF(A1.EQ.A4) GO TO 35
   33 WRITE (6,901) A1
      STOP
   35 READ (5,913)A1,(IORD(I),I=1,NV)
      IF(A1.EQ.A3) GO TO 37
   36 WRITE (6,902) A1
      STOP
   37 MAT=MAT*18
      READ (5,914)(FMT(I),I=1,MAT)
      WRITE (6,916) (FMT(I),I=1,MAT)
  916 FORMAT (' VARIABLE FORMAT IS'/(1X,18A4))
      IF(IPD)55,55,41
   41 DO 42 I=1,LK
   42 PD(I)=IPD*I
   43 IF(UNIT.EQ.D)GO TO 50
   44 IF(UNIT.EQ.W)GO TO 51
   45 IF(UNIT.EQ.XM)GO TO 52
   46 IF(UNIT.EQ.Y)GO TO 53
      GO TO 55
   50 DIV=1.0
      GO TO 55
   51 DIV=7.0
      GO TO 55
   52 DIV=30.4
      GO TO 55
   53 DIV=365.0
 55   DO 125 MM=1,NTOT
      READ (MTAPE,FMT)(TD(I),I=1,NV)
      CALL PTDATA
      CALL RESET (MM)
      IF(JCO)75,75,60
   60 DO 65 I=1,LK
      III=I
      IF(X(1)-COH(I))70,70,65
   65 CONTINUE
   70 II=III
      GO TO 80
   75 II=1
   80 DATA(II,1,2)=DATA(II,1,2)+1.0
      IF(FORM.EQ.A)GO TO 100
   85 IF(FORM.EQ.B)GO TO 115
   90 IF(FORM.EQ.C2)GO TO 110
   95 CALL BFORM(II,DIV,PD)
      GO TO 125
  100 IF(JCO)110,110,105
  105 CALL AFORM(II,DIV,COH)
      GO TO 125
  110 CALL AFORM(II,DIV,PD)
      GO TO 125
  115 IF(JCO)95,95,120
  120 CALL BFORM (II,DIV,COH)
  125 CONTINUE
      WRITE (6,915)
      IF(JCO)127,127,126
  126 LM=LK
      GO TO 121
  127 LM=1
  121 DO 135 I=1,LM
      LL=I-1
      KK=LK-LL
      DO 135 JJ=2,KK
      J=JJ-1
  135 DATA(I,JJ,2)=DATA(I,J,2)-DATA(I,J,3)-DATA(I,J,4)
  900 FORMAT (' COHORT CARD EXPECTED,BUT READ ',A6)
  901 FORMAT (' STATUS CARD EXPECTED,BUT READ ',A6)
  902 FORMAT (' ORDER CARD EXOECTED,BUT READ ',A6)
  911 FORMAT(A6,33F2.0/6X,(33F2.0))
  912 FORMAT(A6,20F3.0)
  913 FORMAT(A6,24(1X,I2))
  914 FORMAT(18A4)
  915 FORMAT(1H1)
      RETURN
      END
      SUBROUTINE REPORT
C             SUBROUTINE REPORT FOR BMD01S             MAY 20, 1966
      DIMENSION DATA(31,31,5),C(31,12),N(5),KSR(31),IDATE(12),DCODE(10),
     1FCODE(10),IORD(25),TD(25),FMT(120),C1(31),C3(31),COH(31),PD(31),
     2X(8),PROB(15,31)
      COMMON  C      , K      , SUM    , STD    , ITABLE , N
      COMMON  SRATE  , DATA   , KSR    , LK     , NTOT   , NV
      COMMON  NJ     , FORM   , JCO    , IPD    , UNIT   , NTAPE
      COMMON  MAT    , C1     , C3     , DCODE  , FCODE  , IORD
      COMMON  IDATE  , TD     , X      , PROB   , MN     , NKI
      WRITE (6,910)ITABLE,K
      WRITE (6,911)
      ITABLE=ITABLE+1
      WRITE (6,912)
      WRITE (6,913)
      WRITE (6,914)
      WRITE (6,919)
      IF(K-NKI)22,20,22
   20 DO 21 I=1,K
   21 PROB(MN,I)=C(I,9)*(-1.0)
   22 N1=0
      N2=1
      DO 30 I=1,K
      DO 29 J=2,5
 29   N(J)=C(I,J)
      WRITE (6,915)N1,N2,N(2),N(3),N(4),N(5),(C(I,M),M=6,9), C(I,11),C(I
     1,12)
      N1=N1+1
 30   N2=N2+1
      WRITE (6,919)
      WRITE (6,916)K,SRATE
      WRITE (6,917)STD
      SUM=(SRATE*(1.0-SRATE))/STD**2
      KK=SUM
      WRITE (6,918)KK
 910  FORMAT(15H1         TABLEI3,20H  COMPUTATION OF THEI3,44H-PERIOD S
     1URVIVAL RATE AND ITS STANDARD ERROR)
 911  FORMAT(120H0PERIODS     ALIVE AT     DIED     LOST TO   WITHDRAWN 
     1  EFFECTIVE NO.   PROPOR-  PROPOR-   CUMULATIVE   COL.7/   STAND.)
 912  FORMAT(120H  AFTER       BEG. OF    DURING    FOLLOW-   ALIVE DUR-
     1   EXPOSED TO      TION   TION SUR-  PROPORTION   COL.6-   ERROR )
 913  FORMAT(111H DIAGNOSIS   INTERVAL   INTERVAL     UP    ING INTERVAL
     1  RISK OF DYING    DYING   VIVING    SURVIVING    COL.3  )
 914  FORMAT(119H   (1)          (2)        (3)       (4)        (5)    
     1       (6)          (7)      (8)        (9)        (10)     (11))
 915  FORMAT(I3,2H -I2,I12,I11,I10,I11,F15.1,F12.2,F9.2,F11.2,F12.4,F10.
     14)
 916  FORMAT(1H0I2,24H-PERIOD SURVIVAL RATE...F7.2)
 917  FORMAT(27H0STANDARD ERROR............F7.2)
 918  FORMAT(27H0EFFECTIVE SAMPLE SIZE.....I4)
 919  FORMAT(1H )
      RETURN
      END
      SUBROUTINE RESET (NVAR)
C             SUBROUTINE RESET FOR BMD01S              MAY 20, 1966
      DIMENSION DATA(31,31,5),C(31,12),N(5),KSR(31),IDATE(12),DCODE(10),
     1FCODE(10),IORD(25),TD(25),FMT(120),C1(31),C3(31),COH(31),PD(31),
     2X(8),VAR(20),PROB(15,31)
      COMMON  C      , K      , SUM    , STD    , ITABLE , N
      COMMON  SRATE  , DATA   , KSR    , LK     , NTOT   , NV
      COMMON  NJ     , FORM   , JCO    , IPD    , UNIT   , NTAPE
      COMMON  MAT    , C1     , C3     , DCODE  , FCODE  , IORD
      COMMON  IDATE  , TD     , X      , PROB   , MN     , NKI
      DO 10 I=1,8
   10 X(I)=0.0
      IF(JCO)15,15,20
   15 DIV1=30.4
      DIV2=365.0
      GO TO 25
   20 DIV1=100.0
      DIV2=10000.0
   25 II=0
      DO 80 I=1,NV
      DO 45 J=1,4
      IF(IORD(I)-IDATE(J))45,40,45
   40 IJ=IDATE(J)/10
      X(IJ)=X(IJ)+TD(I)*DIV1
      GO TO 80
   45 CONTINUE
      DO 65 J=5,8
      J10=J*10
      IF(IORD(I)-J10)55,50,55
   50 X(J)=TD(I)
      GO TO 80
   55 IF(IORD(I)-IDATE(J))65,60,65
   60 IJ=IDATE(J)/10
      X(IJ)=X(IJ)+TD(I)
      GO TO 80
   65 CONTINUE
      DO 75 J=9,12
      IF(IORD(I)-IDATE(J))75,70,75
   70 IJ=IDATE(J)/10
      X(IJ)=X(IJ)+TD(I)*DIV2
      GO TO 80
   75 CONTINUE
   80 CONTINUE
      RETURN
      END
      SUBROUTINE TPWD(NT1,NT2)
C        SUBROUTINE TPWD FOR BMD01S                    MAY 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
   15 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