C     REVISION NO. 0                                         PROGRAM NO.EXP  10
C                                                                       EXP  20
C      PROGRAM INDEX--2      REVISED JULY 1,1965                        EXP  30
      COMMON INT(1254),EL(60),TWOTH(150),NAME(20),S2SQ(150),ES(600),S(15EXP  40
     10),NCUT,XLAMDA,YLAMDA,NI,ETEST,TETBLK,DENSTY,FORMWT,NOCUBE,NOHEX,NEXP  50
     2OTET,NONO                                                         EXP  60
      COMMON/IOUNIT/ IREAD,  IPRNT
1     FORMAT (5I10)                                                     EXP  70
      IREAD = 5
      IPRNT = 6
C *
C *  READ MATRIX FROM PERMENANT SEQUENTIAL FILE ON LOGICAL 1
C *
      OPEN(UNIT=1,FILE='POWD:[PUBLIC.EXE.INDEX]GOBWIL.DAT',
	1STATUS='OLD',READONLY)
      READ(1,1) (INT(I),I=1,1254)
      CLOSE(UNIT=1)
C *
C     GENERATE THE L LIST OF ALL LL                                     EXP 100
      DO 2 I=1,30                                                       EXP 110
2     EL(I)=(I-1)*(I-1)                                                 EXP 120
    3 CALL INPUT
      IF (NOCUBE) 5,4,5                                                 EXP 160
    4 CALL CUBIC
    5 CALL HEXTET
      IF (NONO) 3,6,3                                                   EXP 210
    6 CALL ORTHO
      GO TO 3                                                           EXP 240
      END                                                               EXP 250
      SUBROUTINE INPUT                                                       10
      COMMON INT(1254),EL(60),TWOTH(150),NAME(20),S2SQ(150),ES(600),S(15     20
     10),NCUT,XLAMDA,YLAMDA,NI,ETEST,TETBLK,DENSTY,FORMWT,NOCUBE,NOHEX,N     30
     2OTET,NONO                                                              40
      COMMON/IOUNIT/ IREAD,  IPRNT
1     FORMAT (20A4)                                                          50
2     FORMAT (7F10.3)                                                        60
3     FORMAT (I5,6F10.5,4I1)                                                 70
      READ(IREAD,1,END=20)NAME
      READ(IREAD,3)NCUT,XLAMDA,YLAMDA,ETEST,TETBLK,FORMWT,DENSTY,NOCUBE,     90
     1NOHEX,NOTET,NONO                                                      100
      IF (ETEST) 5,4,5                                                      110
4     ETEST=0.0005                                                          120
5     CONTINUE                                                              130
C                                                                           140
C                                                                           150
C     ETEST IS THE PROBABLE EXPERIMENTAL ERROR IN THE DATA.  IF UNKNOWN,    160
C     LEAVE ETEST BLANK AND PROGRAM WILL SET ETEST=.0005.     IF UNCER-     170
C     TAIN, GUESS LOW.  E.G. IF ERROR IS ABOUT .00025, SET ETEST=.0002.     180
C     NO SOLUTION WILL BE PRINTED FOR E LESS THAN ETEST.                    190
C                                                                           200
C                                                                           210
C     TETBLK WILL PREVENT THE PROGRAM FROM PRODUCING A TETRAGONAL           220
C     INDEXING UNLESS NO HEXAGONAL INDEXING HAS BEEN FOUND.  IN THE         230
C     LATTER CASE THE TETRAGONAL CASE WILL BE TRIED IN SPITE OF THE         240
C     BLOCK.  IF A TETRAGONAL INDEXING IS NOT DESIRED, PUT A POSITIVE       250
C     DECIMAL NUMBER, E.G. .5 IN COLUMNS 36 TO 45 OF THE CUT NUMBER         260
C     CARD. OTHERWISE LEAVE THAT SPACE BLANK.                               270
C                                                                           280
C                                                                           290
C     XLAMDA IS THE UNRESOLVED WAVELENGTH AND YLAMDA IS THE WAVELENGTH      300
C     OF THE K-ALPHA-1 RADIATION.                                           310
C                                                                           320
C                                                                           330
      IF (XLAMDA) 6,6,7                                                     340
6     XLAMDA=1.5418                                                         350
      YLAMDA=1.54051                                                        360
7     CONTINUE                                                              370
      RATIO=YLAMDA/XLAMDA                                                   380
C                                                                           390
C                                                                           400
C     IF NO WAVELENGTHS ARE READ IN THE PROGRAM USES THOSE OF COPPER.       410
C                                                                           420
      DO 8 I=1,22000,7                                                      430
      J=I+6                                                                 440
      READ(IREAD,2)(TWOTH(K),K=I,J)                                         450
      IF (ABS(TWOTH(I))+ABS(TWOTH(I+1))+ABS(TWOTH(I+2))+ABS(TWOTH(I+3))+    460
     1ABS(TWOTH(I+4))+ABS(TWOTH(I+5))+ABS(TWOTH(J))) 8,9,8                  470
8     CONTINUE                                                              480
9     DO 10 J=1,7                                                           490
      N=I-J                                                                 500
      IF (TWOTH(N)) 11,10,11                                                510
10    CONTINUE                                                              520
11    NI=MIN0(N,150)                                                        530
      IF (NCUT) 14,14,12                                                    540
12    DO 13 I=1,NCUT                                                        550
      S(I)=(RATIO**2)*SIN(TWOTH(I)*8.726646E-03)**2                         560
      S2SQ(I)=SIN(TWOTH(I)*1.745329E-03)**2                                 570
13    CONTINUE                                                              580
14    NKUT=NCUT+1                                                           590
      DO 15 I=NKUT,NI                                                       600
      S(I)=SIN(TWOTH(I)*8.726646E-03)**2                                    610
      S2SQ(I)=SIN(TWOTH(I)*1.745329E-02)**2                                 620
15    CONTINUE                                                              630
      MARK=8                                                                640
      GO TO 25
   20 CALL EXIT
   25 RETURN
      END                                                                   660
      SUBROUTINE CUBIC                                                       10
      COMMON INT(1254),EL(60),TWOTH(150),NAME(20),S2SQ(150),ES(600),S(15     20
     10),NCUT,XLAMDA,YLAMDA,NI,ETEST,TETBLK,DENSTY,FORMWT,NOCUBE,NOHEX,N     30
     2OTET,NONO                                                              40
      COMMON/IOUNIT/ IREAD,  IPRNT
      DIMENSION YHAT(150), R(150), C(150), B(150,2), A(999), OMIT(175),      50
     1BE(175)                                                                60
      WRITE(IPRNT,1)                                                         70
1     FORMAT (1H0,'ENTERING SUBROUTINE CUBIC')                               80
2     FORMAT (I5,F15.2,F20.5,3F25.5)                                         90
3     FORMAT (1X,46H THE NUMBER OF FORMULA UNITS PER UNIT CELL IS F10.2)    100
4     FORMAT (52X8HOBSERVED16X10HCALCULATED/4X1HI6X9HTWO THETA8X12HCALCU    110
     1LATED N7X44HSINE SQUARED THETA       SINE SQUARED THETA 15X10HDIFF    120
     2ERENCE)                                                               130
5     FORMAT (77H THE CRYSTAL HAS BEEN INDEXED IN THE CUBIC SYSTEM WITH     140
     1A LATTICE CONSTANT OF F7.4,35H ANGSTROMS,BUT IN VIEW OF THE LARGE/    150
     241H ERROR THIS INDEXING IS PROBABLY INVALID./66H AN ESTIMATE OF TH    160
     3E STANDARD DEVIATION OF THE LATTICE CONSTANT IS F8.6///)              170
6     FORMAT (10H0 ERROR ISF9.4//)                                          180
7     FORMAT (1H1//24X20A4///)                                              190
      DO 8 I=1,125                                                          200
      OMIT(I)=8*I-1                                                         210
8     CONTINUE                                                              220
      DO 9 I=1,31                                                           230
      OMIT(I+125)=4*(8*I-1)                                                 240
9     CONTINUE                                                              250
      DO 10 I=1,7                                                           260
      OMIT(I+156)=16*(8*I-1)                                                270
10    CONTINUE                                                              280
      DO 11 I=1,2                                                           290
      OMIT(I+163)=64*(8*I-1)                                                300
11    CONTINUE                                                              310
      DO 16 I=1,164                                                         320
      IF (OMIT(I)-OMIT(I+1)) 16,16,12                                       330
12    KB=OMIT(I)                                                            340
      OMIT(I)=OMIT(I+1)                                                     350
      OMIT(I+1)=KB                                                          360
      JB=I                                                                  370
13    IF (JB-1) 16,16,14                                                    380
14    IF (OMIT(JB)-OMIT(JB-1)) 15,16,16                                     390
15    KB=OMIT(JB)                                                           400
      OMIT(JB)=OMIT(JB-1)                                                   410
      OMIT(JB-1)=KB                                                         420
      JB=JB-1                                                               430
      GO TO 13                                                              440
16    CONTINUE                                                              450
      I=0                                                                   460
      DO 20 J=1,165                                                         470
      N=I+1                                                                 480
      DO 19 I=N,999                                                         490
      Z=FLOAT(I)                                                            500
      IF (Z-OMIT(J)) 17,18,17                                               510
17    A(I)=I                                                                520
      GO TO 19                                                              530
18    A(I)=0.0000                                                           540
      GO TO 20                                                              550
19    CONTINUE                                                              560
20    CONTINUE                                                              570
      E=0.2                                                                 580
      U=0.0                                                                 590
      V=0.0                                                                 600
      YLE=E                                                                 610
      NB=500                                                                620
      LGO=1                                                                 630
      L=1                                                                   640
      M=1                                                                   650
      MSW=0                                                                 660
      KAY=NI/2                                                              670
21    K=(NI/2)-1                                                            680
22    K=K+1                                                                 690
      N=A(K)                                                                700
      I=(NI/2)-1                                                            710
      J=(NI/2)-1                                                            720
      IF (N) 23,22,23                                                       730
23    IF (N-NB) 24,32,32                                                    740
24    X=S(KAY)/FLOAT(N)                                                     750
25    I=I+1                                                                 760
      IF (I-NI) 27,27,26                                                    770
26    GO TO (37,47), LGO                                                    780
27    C(I)=S(I)/X                                                           790
28    J=J+1                                                                 800
      IF (J-999) 29,22,29                                                   810
29    IF (C(I)-A(J)-E) 30,30,28                                             820
30    IF (C(I)-A(J)+E) 22,31,31                                             830
31    B(I,1)=S(I)                                                           840
      B(I,2)=A(J)                                                           850
      GO TO 25                                                              860
32    GO TO (33,42), M                                                      870
33    IF (MSW+1) 34,35,34                                                   880
34    IF (MSW) 36,35,36                                                     890
35    YLE=E                                                                 900
      E=E*2.                                                                910
      MSW=-1                                                                920
      GO TO 21                                                              930
36    M=2                                                                   940
      GO TO 41                                                              950
37    IF (MSW) 38,46,38                                                     960
38    GO TO (39,44), L                                                      970
39    IF (MSW-1) 40,46,40                                                   980
40    L=2                                                                   990
      MSW=1                                                                1000
41    U=AMIN1(E,YLE)                                                       1010
      V=AMAX1(E,YLE)                                                       1020
      E=(U+V)/2.                                                           1030
      YLE=E                                                                1040
      M=2                                                                  1050
      L=2                                                                  1060
      GO TO 21                                                             1070
42    U=E                                                                  1080
      YLE=E                                                                1090
      E=(E+V)/2.                                                           1100
      YLE=E                                                                1110
      IF (V-U-0.05) 43,43,21                                               1120
43    E=V                                                                  1130
      LGO=2                                                                1140
      GO TO 21                                                             1150
44    V=E                                                                  1160
      IF (V-U-.05) 47,47,45                                                1170
45    E=(E+U)/2.                                                           1180
      GO TO 21                                                             1190
46    YLE=E                                                                1200
      E=E/2.                                                               1210
      MSW=1                                                                1220
      GO TO 21                                                             1230
47    SUMSQ=0.0                                                            1240
      DO 48 I=KAY,NI                                                       1250
48    SUMSQ=SUMSQ+(B(I,2)*B(I,2))                                          1260
      DNUM=0.0                                                             1270
      DO 49 I=KAY,NI                                                       1280
49    DNUM=DNUM+(B(I,1)*B(I,2))                                            1290
      AHAT=DNUM/SUMSQ                                                      1300
      SUM=0.0                                                              1310
      SUMSQ=0.                                                             1320
      DO 50 I=1,NI                                                         1330
      B(I,1)=S(I)                                                          1340
      B(I,2)=AINT((S(I)/AHAT)+.5)                                          1350
      YHAT(I)=AHAT*B(I,2)                                                  1360
      R(I)=B(I,1)-YHAT(I)                                                  1370
      SUM=SUM+R(I)**2                                                      1380
      SUMSQ=SUMSQ+B(I,2)**2                                                1390
50    CONTINUE                                                             1400
      DEGFRD=NI-1                                                          1410
      SIGMAA=SQRT(SUM/(DEGFRD*SUMSQ))                                      1420
      XLC=YLAMDA/(2.*SQRT(AHAT))                                           1430
      WRITE(IPRNT,7)NAME                                                   1440
      WRITE(IPRNT,6)E                                                      1450
      IF (FORMWT) 52,51,52                                                 1460
51    ZNMBER=XLC**3*DENSTY*0.6023                                          1470
      GO TO 53                                                             1480
52    ZNMBER=XLC**3*DENSTY*0.6023/FORMWT                                   1490
53    WRITE(IPRNT,3)ZNMBER                                                 1500
      IF (E-.3) 54,54,55                                                   1510
54    WRITE(IPRNT,5)XLC,SIGMAA                                             1520
      GO TO 56                                                             1530
55    WRITE(IPRNT,5)XLC,SIGMAA                                             1540
56    WRITE(IPRNT,4)                                                       1550
      DO 59 I=1,NI                                                         1560
      LTMP=I+1                                                             1570
      IF (MOD(LTMP,40)) 58,57,58                                           1580
57    WRITE(IPRNT,7)NAME                                                   1590
      WRITE(IPRNT,4)                                                       1600
58    WRITE(IPRNT,2)(I,TWOTH(I),B(I,2),S(I),YHAT(I),R(I))                  1610
59    CONTINUE                                                             1620
      WRITE(IPRNT,60)                                                      1630
60    FORMAT (1H0,'EXITING FROM SUBROUTINE CUBIC')                         1640
      RETURN                                                               1650
      END                                                                  1660
      SUBROUTINE HEXTET                                                      10
C     PROGRAM TO DETERMINE CRYSTAL INDEXES, HEXAGONAL CASE                   20
      COMMON INT(1254),EL(60),TWOTH(150),NAME(20),S2SQ(150),ES(600),S(15     30
     10),NCUT,XLAMDA,YLAMDA,NI,ETEST,TETBLK,DENSTY,FORMWT,NOCUBE,NOHEX,N     40
     2OTET,NONO                                                              50
      COMMON/IOUNIT/ IREAD,  IPRNT
      DIMENSION O(150), D(150), KA(150,2), ASTOR(25), CSTOR(25), FACTOR(     60
     125), C(150), DUM1(150,2), A(150,2), DUM2(155,2)                        70
1     FORMAT (44X9H OBSERVED15X10HCALCULATED/2X1HI4X9HTWO THETA4X8HHH+KK     80
     1+HK5X2HLL5X18HSINE SQUARED THETA5X18HSINE SQUARED THETA8X10HDIFFER     90
     2ENCE//)                                                               100
2     FORMAT (44X9H OBSERVED15X10HCALCULATED/2X1HI4X9HTWO THETA4X5HHH+KK    110
     18X2HLL5X18HSINE SQUARED THETA5X18HSINE SQUARED THETA8X10HDIFFERENC    120
     2E//)                                                                  130
3     FORMAT (1XI3,F11.3,2I10,F20.5,F23.5,F20.5)                            140
4     FORMAT (90H0THE FOLLOWING INDEXING IS ON THE BASIS OF THE TETRAGON    150
     1AL SYSTEM WITH THE ERROR LESS THAN F8.6)                              160
5     FORMAT (89H0THE FOLLOWING INDEXING IS ON THE BASIS OF THE HEXAGONA    170
     1L SYSTEM WITH THE ERROR LESS THAN F8.6)                               180
6     FORMAT (64H AN ESTIMATE OF THE STANDARD DEVIATION OF SINE SQUARED     190
     1THETA IS F10.4/51H AN ESTIMATE OF THE STANDARD DEVIATION OF THETA     200
     2IS F10.4/51H AN ESTIMATE OF THE STANDARD DEVIATION OF AZERO IS F10    210
     3.4/51H AN ESTIMATE OF THE STANDARD DEVIATION OF CZERO IS F10.4/9H     220
     4AHAT IS F10.6,2X9H BHAT IS F10.6/31H0THE LATTICE CONSTANT AZERO IS    230
     5 F10.5,11H ANGSTROMS./31H0THE LATTICE CONSTANT CZERO IS F10.5,11H     240
     6ANGSTROMS.//)                                                         250
7     FORMAT (36H0THE UNRESOLVED WAVELENGTH USED WAS F8.5,41H AND THE K-    260
     1ALPHA-ONE WAVELENGTH USED WAS F8.5/16H THE DENSITY IS F6.3,1H,,22H    270
     2 THE FORMULA WEIGHT ISF8.3,21H, AND THE NUMBER IS  F6.3,1H.)          280
8     FORMAT (11H0THE FIRST I3,53H LINES ARE UNRESOLVED. THE REST ARE K-    290
     1ALPHA-1 LINES. )                                                      300
9     FORMAT (92H0THE PROGRAM IS UNABLE TO INDEX THE CRYSTAL ON A HEXAGO    310
     1NAL BASIS WITHIN THE MAXIMUM ERROR OFF10.5)                           320
10    FORMAT (93H0THE PROGRAM IS UNABLE TO INDEX THE CRYSTAL ON A TETRAG    330
     1ONAL BASIS WITHIN THE MAXIMUM ERROR OFF10.5)                          340
11    FORMAT (64H0THE PROGRAM IS UNABLE TO INDEX THE CRYSTAL WITHIN THE     350
     1ERROR OF F8.6)                                                        360
12    FORMAT (1H1/24X20A4//)                                                370
      NSTOR=0                                                               380
      NSTORE=0                                                              390
      KKKK=0                                                                400
C     GET SET FOR THE HEXAGONAL CASE.                                       410
      R=1.0                                                                 420
      KT=1                                                                  430
C     GENERATE THE S LIST OF ALL HH+KK+HK                                   440
      TETBLK=NOTET                                                          450
      IF (NOHEX) 121,13,121                                                 460
13    I=0                                                                   470
      NAY=1                                                                 480
      KAY=1                                                                 490
      DO 15 M=1,32                                                          500
      DO 14 K=1,M                                                           510
      I=I+1                                                                 520
      ES(I)=(M-1)*(M-1)+(K-1)*(K-1)+KT*(M-1)*(K-1)                          530
14    CONTINUE                                                              540
15    CONTINUE                                                              550
      DO 16 K=1,25                                                          560
      ASTOR(K)=0.                                                           570
      CSTOR(K)=0.                                                           580
16    CONTINUE                                                              590
      KSKIP=2                                                               600
C     ORDER THE S LIST                                                      610
      CALL SORTAC(ES,528)                                                   620
      LP=0                                                                  630
      DO 20 I=2,277                                                         640
17    IF (ES(I-1)-ES(I)) 20,18,165                                          650
18    LP=LP+1                                                               660
      M=529-LP                                                              670
      DO 19 J=I,M                                                           680
      ES(J-1)=ES(J)                                                         690
19    CONTINUE                                                              700
      GO TO 17                                                              710
20    CONTINUE                                                              720
C     OMMANDS 594 TO 3 DELETE DUPLICATE INTEGERS.                           730
      NES=277                                                               740
      NEL=30                                                                750
      B=1.0                                                                 760
      DO 21 I=2,NI                                                          770
      D(I)=ABS(S(I)-S(I-1))                                                 780
      B=AMIN1(B,D(I))                                                       790
21    CONTINUE                                                              800
C      B IS THE MINIMUM DIFFERENCE BETWEEN S(I) AND S(I+1).                 810
22    E=.49*B                                                               820
C     OMMANDS   TO HERE FIX THE MAXIMUM ERROR E                             830
      YLE=4.*E                                                              840
      KAT=0                                                                 850
C     KAT IS INCREASED BY 1 WHEN WE FAIL TO INDEX THE CRYSTAL.              860
      NOD=0                                                                 870
      IGO=1                                                                 880
C     THE VALUE OF IGO DETERMINES WHETHER THE ERROR IS BEING DECREASED      890
C     OR IS FIXED.                                                          900
23    KAP=0                                                                 910
C     KAP IS INCREASED BY 1 EACH ITERATION.                                 920
      IF (KAT) 24,24,28                                                     930
C     KAT=1 MEANS THE MINIMUM ERROR HAS BEEN DETERMINED AND THE             940
C     PROGRAM IS NOW FINDING ALL SOLUTIONS FOR THAT FIXED ERROR.            950
24    GO TO (25,26), IGO                                                    960
25    E=YLE/2.0                                                             970
      IF (E-ETEST) 90,90,26                                                 980
26    MI=1                                                                  990
27    MI=MI+1                                                              1000
      IF (MI-5) 29,29,90                                                   1010
28    NOD=NOD+1                                                            1020
      IF (NOD-5) 36,36,90                                                  1030
29    NUT=1                                                                1040
      MUT=0                                                                1050
      KUP=1                                                                1060
      LUP=1                                                                1070
      KQ=MI-1                                                              1080
30    KU=MI                                                                1090
      KV=0                                                                 1100
31    KV=KV+1                                                              1110
      IF (KV-KQ) 32,32,54                                                  1120
32    KX=0                                                                 1130
33    KX=KX+1                                                              1140
      IF (KX-KQ) 35,35,34                                                  1150
34    GO TO (31,54), KUP                                                   1160
35    KY=0                                                                 1170
36    KY=KY+1                                                              1180
      IF (KY-KQ) 38,38,37                                                  1190
37    GO TO (33,54), LUP                                                   1200
38    GO TO (39,40,41,42,43,44,45,46,47,48,49,50,51,52,53), NUT            1210
39    IE=KU                                                                1220
      JE=KV                                                                1230
      NE=KX                                                                1240
      KE=KY                                                                1250
      GO TO 58                                                             1260
40    IE=KV                                                                1270
      JE=KU                                                                1280
      NE=KX                                                                1290
      KE=KY                                                                1300
      GO TO 58                                                             1310
41    IE=KX                                                                1320
      JE=KV                                                                1330
      NE=KU                                                                1340
      KE=KY                                                                1350
      GO TO 58                                                             1360
42    IE=KV                                                                1370
      JE=KX                                                                1380
      NE=KY                                                                1390
      KE=KU                                                                1400
      GO TO 58                                                             1410
43    IE=KU                                                                1420
      JE=KV                                                                1430
      NE=KX                                                                1440
      KE=KY                                                                1450
      GO TO 58                                                             1460
44    IE=KX                                                                1470
      JE=KU                                                                1480
      NE=KY                                                                1490
      KE=KV                                                                1500
      GO TO 58                                                             1510
45    IE=KU                                                                1520
      JE=KX                                                                1530
      NE=KV                                                                1540
      KE=KY                                                                1550
      GO TO 58                                                             1560
46    IE=KX                                                                1570
      JE=KY                                                                1580
      NE=KU                                                                1590
      KE=KV                                                                1600
      GO TO 58                                                             1610
47    IE=KX                                                                1620
      JE=KU                                                                1630
      NE=KV                                                                1640
      KE=KY                                                                1650
      GO TO 58                                                             1660
48    IE=KU                                                                1670
      JE=KX                                                                1680
      NE=KY                                                                1690
      KE=KV                                                                1700
      GO TO 58                                                             1710
49    IE=KU                                                                1720
      JE=KV                                                                1730
      NE=KX                                                                1740
      KE=KY                                                                1750
      GO TO 58                                                             1760
50    IE=KU                                                                1770
      JE=KV                                                                1780
      NE=KY                                                                1790
      KE=KX                                                                1800
      GO TO 58                                                             1810
51    IE=KU                                                                1820
      JE=KY                                                                1830
      NE=KV                                                                1840
      KE=KX                                                                1850
      GO TO 58                                                             1860
52    IE=KY                                                                1870
      JE=KU                                                                1880
      NE=KV                                                                1890
      KE=KX                                                                1900
      GO TO 58                                                             1910
53    IE=KU                                                                1920
      JE=KV                                                                1930
      NE=KX                                                                1940
      KE=KY                                                                1950
      GO TO 58                                                             1960
54    NUT=NUT+1                                                            1970
      IF (NUT-14) 55,55,27                                                 1980
55    MUT=MUT+1                                                            1990
      GO TO (30,30,30,56,32,32,32,32,32,57,35,35,35,53), MUT               2000
56    KUP=2                                                                2010
      GO TO 32                                                             2020
57    LUP=2                                                                2030
      GO TO 35                                                             2040
58    T=ES(IE)*EL(KE)-ES(NE)*EL(JE)                                        2050
C                                                                          2060
C     OMMANDS 9  THROUGH 76 COUNT THROUGH ALL 2X2 MATRICES WITH ENTRIES    2070
C     THAN MI+2                                                            2080
C                                                                          2090
C     T IS THE COEFFICIENT DETERMINANT                                     2100
C     IF T=0. WE HAVE DEPENDENCE                                           2110
      KAN=KAY+NAY                                                          2120
C     X AND Y ARE DETERMINED FROM LINES KAY AND KAN.                       2130
      IF (T) 70,59,70                                                      2140
C     OMMANDS 101 THROUGH 109 CHECK FOR DEPENDENCE                         2150
59    IF (ES(IE)) 61,60,61                                                 2160
60    IF (EL(JE)) 64,36,64                                                 2170
61    IF (ES(NE)) 62,63,62                                                 2180
62    SA=ES(IE)/ES(NE)                                                     2190
      GO TO 67                                                             2200
63    SA=0.0                                                               2210
      GO TO 67                                                             2220
64    IF (EL(KE)) 65,66,65                                                 2230
65    SA=EL(JE)/EL(KE)                                                     2240
      GO TO 67                                                             2250
66    SA=0.0                                                               2260
67    IF (SA-(S(KAY)/S(KAN))-0.001) 68,36,36                               2270
68    IF (SA-(S(KAY)/S(KAN))+0.001) 36,36,69                               2280
69    NAY=NAY+1                                                            2290
      IF (KAY+NAY-NI) 22,22,95                                             2300
70    X=(S(KAY)*EL(KE)-S(KAN)*EL(JE))/T                                    2310
      IF (X-.001) 36,36,71                                                 2320
71    Y=(S(KAN)*ES(IE)-S(KAY)*ES(NE))/T                                    2330
      IF (Y-.001) 36,36,72                                                 2340
72    KP=0                                                                 2350
73    NND=0                                                                2360
C     NND=0 UNTIL WE FIND AN ACCEPTABLE SET OF INDICES AT COMMANDS 244     2370
C     AND FOLLOWING.                                                       2380
      OLDDE=1.1                                                            2390
      KP=KP+1                                                              2400
C     KP IS THE NUMBER OF THE REFLECTION CURRENTLY BEING INDEXED           2410
      YLE=E                                                                2420
      IF (KP-4) 75,75,74                                                   2430
74    KP=4                                                                 2440
      GO TO 98                                                             2450
75    NR=0                                                                 2460
76    NR=NR+1                                                              2470
      IF (NR-NES) 77,77,85                                                 2480
77    IF (X*ES(NR)-S(KP)-E) 78,78,85                                       2490
78    NS=0                                                                 2500
79    NS=NS+1                                                              2510
      IF (NS-NEL) 80,80,76                                                 2520
80    C(KP)=X*ES(NR)+Y*EL(NS)                                              2530
C     %KP< IS THE CALCULATED SINE SQUARED THETA TO COMPARE WITH THE        2540
C     OBSERVED SINE SQUARED THETA.                                         2550
      IF (C(KP)-S(KP)-E) 81,81,76                                          2560
81    IF (C(KP)-S(KP)+E) 79,82,82                                          2570
82    DE=ABS(S(KP)-C(KP))                                                  2580
      IF (OLDDE-DE) 79,79,83                                               2590
83    IF (KP-1) 86,86,87                                                   2600
84    A(KP,1)=B1                                                           2610
      A(KP,2)=B2                                                           2620
C     WE ACCEPT THE TWO INDICES ES(NR) AND EL(NS) SO GO INCREASE KP        2630
C     BY 1 AND TRY TO INDEX THE NEXT REFLECTION.                           2640
      GO TO 73                                                             2650
85    IF (NND) 36,36,84                                                    2660
86    B1=ES(NR)                                                            2670
      B2=EL(NS)                                                            2680
      OLDDE=DE                                                             2690
      NND=1                                                                2700
      GO TO 79                                                             2710
C     OMMAND 244 STORES THE CURRENT BEST SELECTION OF INDICES FOR LINE K   2720
87    IF (A(KP-1,1)-ES(NR)) 86,88,86                                       2730
88    IF (A(KP-1,2)-EL(NS)) 86,89,86                                       2740
C     OMMANDS 246 AND 247 AVOID DUPLICATE SETS OF INDICES.                 2750
89    GO TO 79                                                             2760
90    KAT=KAT+1                                                            2770
      IF (KAT-1) 91,91,92                                                  2780
91    E=2.*E                                                               2790
      KSKIP=1                                                              2800
C     OMMANDS 250 AND 260 ACKNOWLEDGE THE FIRST UNSUCCESSFUL E AND SET T   2810
C     PROGRAM TO FIND AND PRINT ALL SOLUTIONS OBTAINABLE FOR THE LAST      2820
C     SUCCESSFUL E.                                                        2830
      IGO=2                                                                2840
      GO TO 26                                                             2850
92    IF (KT) 93,93,121                                                    2860
C     KT=0 IS THE TETRAGONAL CASE, KT=1 IS THE HEXAGONAL CASE.             2870
93    IF (KAT-1) 91,91,94                                                  2880
94    IF (KKKK) 163,163,164                                                2890
95    XKT=KT                                                               2900
      IF (XKT-0.5) 97,97,96                                                2910
96    WRITE(IPRNT,12)NAME                                                  2920
      WRITE(IPRNT,9)E                                                      2930
      GO TO 121                                                            2940
C                                                                          2950
C     OMPLETE FAILURE FOR HEXAGONAL CASE.  GO TO TETRAGONAL CASE.          2960
C                                                                          2970
97    WRITE(IPRNT,12)NAME                                                  2980
      WRITE(IPRNT,10)E                                                     2990
      GO TO 164                                                            3000
C                                                                          3010
C          COMPLETE FAILURE FOR BOTH CASES.                                3020
C                                                                          3030
C                                                                          3040
C     OMMANDS 600 TO 628 PREPARE FOR LEAST SQUARES CALCULATION OF          3050
C     AHAT AND BHAT.                                                       3060
98    ESIX=S(1)*A(1,1)+S(2)*A(2,1)+S(3)*A(3,1)                             3070
      ESIY=S(1)*A(1,2)+S(2)*A(2,2)+S(3)*A(3,2)                             3080
      YX=A(1,1)*A(1,2)+A(2,1)*A(2,2)+A(3,1)*A(3,2)                         3090
      YI=A(1,2)**2+A(2,2)**2+A(3,2)**2                                     3100
      XI=A(1,1)**2+A(2,1)**2+A(3,1)**2                                     3110
99    ESIX=ESIX+S(KP)*A(KP,1)                                              3120
      ESIY=ESIY+S(KP)*A(KP,2)                                              3130
      YX=YX+A(KP,1)*A(KP,2)                                                3140
      YI=YI+A(KP,2)**2                                                     3150
      XI=XI+A(KP,1)**2                                                     3160
100   AHAT=(ESIX*YI-ESIY*YX)/(XI*YI-YX**2)                                 3170
      BHAT=(ESIY-YX*AHAT)/YI                                               3180
C     AHAT AND BHAT ARE THE LEAST SQUARES ESTIMATES OF A AND C IN THE      3190
C     EQUATION S(I)=A*ES(I)+C*EL(I)                                        3200
      IF (KP-NI) 101,125,125                                               3210
C                                                                          3220
C     OMMAND 628 ASKS IF WE ARE DONE.  IF NOT, COMMANDS 629 THROUGH 690    3230
C     THE REST OF THE REFLECTIONS.                                         3240
C                                                                          3250
101   KP=KP+1                                                              3260
      ND=0                                                                 3270
      OLDDE=1.1                                                            3280
      IR=0                                                                 3290
102   IR=IR+1                                                              3300
C     OMMAND 657 SAYS THIS IS THE FIRST ATTEMPT TO INDEX THIS REFLECTION   3310
      IF (IR-NES) 104,104,103                                              3320
103   IF (ND) 165,111,117                                                  3330
C     OMMAND 667 ASKS IF WE FAILED TO INDEX THIS REFLECTION YES,YES,NO     3340
104   IS=0                                                                 3350
105   IS=IS+1                                                              3360
      IF (IS-NEL) 106,106,102                                              3370
106   C(KP)=AHAT*ES(IR)+BHAT*EL(IS)                                        3380
      IF (C(KP)-S(KP)-E) 109,109,107                                       3390
107   IF (EL(IS)) 165,108,102                                              3400
108   IF (ND) 165,111,117                                                  3410
109   IF (C(KP)-S(KP)+E) 105,110,110                                       3420
110   DE=ABS(S(KP)-C(KP))                                                  3430
      IF (OLDDE-DE) 105,105,113                                            3440
C                                                                          3450
C     WE ASK IF THE NEW SET OF INDICES IS BETTER THAN THE PREVIOUS SET.    3460
C                                                                          3470
111   IF (KAP) 165,36,36                                                   3480
C     IS THIS THE FIRST ITERATION                                          3490
112   B1=ES(IR)                                                            3500
      B2=EL(IS)                                                            3510
      OLDDE=DE                                                             3520
      ND=1                                                                 3530
      GO TO 105                                                            3540
C     OMMAND 686 STORES PRESENT BEST INDICES AS B1 AND B2.                 3550
113   IF (KP-1) 112,112,114                                                3560
114   IF (A(KP-1,1)-ES(IR)) 112,115,112                                    3570
115   IF (A(KP-1,2)-EL(IS)) 112,116,112                                    3580
C     OMMANDS 687 AND 688 CHECK FOR DUPLICATE SETS OF INDICES.             3590
116   GO TO 105                                                            3600
117   A(KP,1)=B1                                                           3610
      A(KP,2)=B2                                                           3620
C     OMMAND 690 ACCEPTS B1 AND B2 AS THE BEST SET OF INDICES FOR LINE K   3630
      IF (KAP) 99,99,118                                                   3640
C     IS THIS THE FIRST ITERATION YES,YES,NO                               3650
118   IF (KP-NI) 101,119,119                                               3660
119   ESIX=0.                                                              3670
      ESIY=0.                                                              3680
      YX=0.                                                                3690
      YI=0.                                                                3700
      XI=0.                                                                3710
      DO 120 L=1,NI                                                        3720
      ESIX=ESIX+S(L)*A(L,1)                                                3730
      ESIY=ESIY+S(L)*A(L,2)                                                3740
      YX=YX+A(L,1)*A(L,2)                                                  3750
      YI=YI+A(L,2)**2                                                      3760
      XI=XI+A(L,1)**2                                                      3770
120   CONTINUE                                                             3780
      X=AHAT                                                               3790
      Y=BHAT                                                               3800
      GO TO 100                                                            3810
C     OMMANDS 692 TO 693 PERFORM ALL LEAST SQUARES FITS EXCEPT THE FIRST   3820
121   KT=0                                                                 3830
      NSTOR=0                                                              3840
      KSKIP=2                                                              3850
      R=0.8660255                                                          3860
      IF (TETBLK) 13,13,164                                                3870
C     OMMAND 702 BRANCHES THE PROGRAM INTO THE TETRAGONAL CASE.            3880
C                                                                          3890
122   IF (KAP-5) 123,123,126                                               3900
C     IS THIS THE SIXTH ITERATION NO,NO,YES                                3910
123   U=ABS(X-AHAT)                                                        3920
      V=ABS(Y-BHAT)                                                        3930
      IF (U+V-.002) 126,126,124                                            3940
124   X=AHAT                                                               3950
      Y=BHAT                                                               3960
      GO TO 101                                                            3970
125   KAP=KAP+1                                                            3980
C     THE PROGRAM ADVANCES TO 706 UPON SUCCESSFULLY INDEXING THE CRYSTAL   3990
      KP=0                                                                 4000
      IF (KAP-1) 101,101,122                                               4010
126   DO 127 K=1,NI                                                        4020
      O(K)=AHAT*A(K,1)+BHAT*A(K,2)                                         4030
127   CONTINUE                                                             4040
C     O(K) IS THE CALCULATED SINE SQUARED THETA.                           4050
      DO 128 K=1,NI                                                        4060
      D(K)=S(K)-O(K)                                                       4070
128   CONTINUE                                                             4080
      SIGMA=0.                                                             4090
      SIGTH=0.                                                             4100
      DO 129 I=1,NI                                                        4110
      SIGMA=SIGMA+D(I)**2                                                  4120
      SIGTH=SIGTH+((D(I)**2)/S2SQ(I))                                      4130
129   CONTINUE                                                             4140
      SIGMA=SIGMA/FLOAT(NI-2)                                              4150
      SIGMA=SQRT(SIGMA)                                                    4160
      SIGTH=SIGTH/FLOAT(NI-2)                                              4170
      SIGTH=SQRT(SIGTH)*57.2957795                                         4180
      SIGHAT=YI*(SIGTH**2)/(XI*YI-YX**2)                                   4190
      SIGHBT=XI*SIGHAT/YI                                                  4200
      SIGHAT=SQRT(SIGHAT)                                                  4210
      SIGHBT=SQRT(SIGHBT)                                                  4220
C     OMMANDS FOLLOWING 731 CALCULATE THE VARIOUS STANDARD DEVIATIONS.     4230
      AZERO=YLAMDA*R/SQRT(3.0*AHAT)                                        4240
      CZERO=YLAMDA/(2.*SQRT(BHAT))                                         4250
      GO TO (130,144), KSKIP                                               4260
130   ASTORE=AZERO                                                         4270
      CSTORE=CZERO                                                         4280
      DO 132 I=1,25                                                        4290
      DO 131 K=1,20                                                        4300
      FACTOR(K)=SQRT(ES(K))                                                4310
      IF (ABS(ASTORE-FACTOR(K)*ASTOR(I))-0.001) 133,133,131                4320
131   CONTINUE                                                             4330
132   CONTINUE                                                             4340
133   CONTINUE                                                             4350
      DO 135 I=1,25                                                        4360
      DO 134 K=1,20                                                        4370
      FACTOR(K)=SQRT(EL(K))                                                4380
      IF (ABS(CSTORE-FACTOR(K)*CSTOR(I))-0.001) 162,162,134                4390
134   CONTINUE                                                             4400
135   CONTINUE                                                             4410
      IF (KT) 165,136,142                                                  4420
136   CONTINUE                                                             4430
      DO 138 I=1,25                                                        4440
      DO 137 K=1,20                                                        4450
      FACTOR(K)=K                                                          4460
      IF (ABS(ASTORE-FACTOR(K)*CSTOR(I))-.001) 139,139,137                 4470
137   CONTINUE                                                             4480
138   CONTINUE                                                             4490
      GO TO 142                                                            4500
139   CONTINUE                                                             4510
      DO 141 I=1,25                                                        4520
      DO 140 K=1,20                                                        4530
      FACTOR(K)=K                                                          4540
      IF (ABS(CSTORE-FACTOR(K)*ASTOR(I))-.001) 162,162,140                 4550
140   CONTINUE                                                             4560
141   CONTINUE                                                             4570
142   CONTINUE                                                             4580
      NSTOR=NSTOR+1                                                        4590
      ASTOR(NSTOR)=ASTORE                                                  4600
      CSTOR(NSTOR)=CSTORE                                                  4610
      IF (NSTOR-25) 143,143,121                                            4620
143   CONTINUE                                                             4630
144   CONTINUE                                                             4640
C                                                                          4650
C     ONCE THE MINIMUM ERROR HAS BEEN DETERMINED, THE COMMANDS FROM        4660
C     721 TO 727 WILL ELIMINATE SOLUTIONS WHICH ARE MULTIPLES OF           4670
C     PREVIOUSLY ACCEPTED SOLUTIONS FOR THAT ERROR.                        4680
C                                                                          4690
      SIGAZO=AZERO*SIGHAT/(2.*AHAT)                                        4700
      SIGCZO=CZERO*SIGHBT/(2.*BHAT)                                        4710
      DO 145 I=1,NI                                                        4720
      KA(I,1)=A(I,1)                                                       4730
      KA(I,2)=A(I,2)                                                       4740
145   CONTINUE                                                             4750
      WRITE(IPRNT,12)NAME                                                  4760
      IF (KT-1) 146,147,147                                                4770
C     WAS THIS CASE HEXAGONAL OR TETRAGONAL                                4780
146   WRITE(IPRNT,4)E                                                      4790
      KK=1                                                                 4800
      KKK=1                                                                4810
      VOLUME=AZERO*CZERO*AZERO                                             4820
      GO TO 148                                                            4830
147   WRITE(IPRNT,5)E                                                      4840
      KK=2                                                                 4850
      KKK=2                                                                4860
      VOLUME=AZERO*AZERO*CZERO*0.866026                                    4870
148   CONTINUE                                                             4880
      IF (FORMWT) 149,149,150                                              4890
149   ZNMBER=0.0                                                           4900
      GO TO 153                                                            4910
150   CONTINUE                                                             4920
      IF (FORMWT) 152,151,152                                              4930
151   ZNMBER=VOLUME*DENSTY*0.6023                                          4940
      GO TO 153                                                            4950
152   ZNMBER=VOLUME*DENSTY*0.6023/FORMWT                                   4960
153   CONTINUE                                                             4970
      WRITE(IPRNT,7)XLAMDA,YLAMDA,DENSTY,FORMWT,ZNMBER                     4980
      WRITE(IPRNT,8)NCUT                                                   4990
      WRITE(IPRNT,6)SIGMA,SIGTH,SIGAZO,SIGCZO,AHAT,BHAT,AZERO,CZERO        5000
      GO TO (154,155), KK                                                  5010
154   WRITE(IPRNT,2)                                                       5020
      GO TO 156                                                            5030
155   WRITE(IPRNT,1)                                                       5040
156   DO 161 I=1,NI                                                        5050
      ICHK=I+6                                                             5060
      IF (MOD(ICHK,36)) 160,157,160                                        5070
157   WRITE(IPRNT,12)NAME                                                  5080
C     OMMANDS 770 TO 771 COUNT THE LINES TO A PAGE SO YOU DO               5090
C     NOT WRITE OVER THE FOLD.                                             5100
      GO TO (158,159), KKK                                                 5110
158   WRITE(IPRNT,2)                                                       5120
      GO TO 160                                                            5130
159   WRITE(IPRNT,1)                                                       5140
160   WRITE(IPRNT,3)(I,TWOTH(I),KA(I,1),KA(I,2),S(I),O(I),D(I))            5150
161   CONTINUE                                                             5160
162   CONTINUE                                                             5170
      KKKK=KKKK+1                                                          5180
      GO TO 23                                                             5190
163   WRITE(IPRNT,12)NAME                                                  5200
      WRITE(IPRNT,11)E                                                     5210
164   RETURN                                                               5220
165   WRITE(IPRNT,166)                                                     5230
166   FORMAT (1H0,'SOMETHING IS WRONG IN HEXTET. GO BACK AND GET MORE DA   5240
     1TA AND TRY ANOTHER CASE')                                            5250
      RETURN                                                               5260
      END                                                                  5270
      SUBROUTINE SORTAC (A,N)                                                10
      DIMENSION A(528)                                                       20
      COMMON INT(1254),EL(60),TWOTH(150),NAME(20),S2SQ(150),ES(600),S(15     30
     10),NCUT,XLAMDA,YLAMDA,NI,ETEST,TETBLK,DENSTY,FORMWT,NOCUBE,NOHEX,N     40
     2OTET,NONO                                                              50
      COMMON/IOUNIT/ IREAD,  IPRNT
      NN=N-1                                                                 60
      DO 3 I=1,NN                                                            70
      K=I+1                                                                  80
      DO 2 J=K,N                                                             90
      IF (A(J)-A(I)) 1,2,2                                                  100
1     TEMP=A(I)                                                             110
      A(I)=A(J)                                                             120
      A(J)=TEMP                                                             130
2     CONTINUE                                                              140
3     CONTINUE                                                              150
      RETURN                                                                160
      END                                                                   170
      SUBROUTINE ORTHO                                                       10
C     PROGRAM TO DETERMINE CRYSTAL INDEXES, ORTHORHOMBIC CASE.               20
      COMMON INT(1254),SQ(60),TWOTH(150),NAME(20),S2SQ(150),ES(600),S(15     30
     10),NCUT,XLAMDA,YLAMDA,NI,ETEST,TETBLK,DENSTY,FORMWT,NOCUBE,NOHEX,N     40
     2OTET,NONO                                                              50
      COMMON/IOUNIT/ IREAD,  IPRNT
      DIMENSION AH(150), AK(150), AL(150), ESS(150), C(150), DE(150), KA     60
     1H(150), KAL(150), K9(9), STOREA(30), STOREB(30), STOREC(30), KAK(1     70
     250)                                                                    80
1     FORMAT (40X8HOBSERVED11X10HCALCULATED/2X1HI3X9HTWO THETA3X2HHH3X2H     90
     1KK3X2HLL5X38HSINE SQUARED THETA  SINE SQUARED THETA5X10HDIFFERENCE    100
     2//)                                                                   110
2     FORMAT (1H1/1X20A4/93H THE FOLLOWING INDEXING IS ON THE BASIS OF      120
     1THE ORTHORHOMBIC SYSTEM WITH THE ERROR LESS THAN F8.6/31H THE LATT    130
     2ICE CONSTANT AZERO IS F10.5,11H ANGSTROMS./31H THE LATTICE CONSTAN    140
     3T BZERO IS F10.5,11H ANGSTROMS./31H THE LATTICE CONSTANT CZERO IS     150
     4F10.5,11H ANGSTROMS.)                                                 160
3     FORMAT (57H0THE PROGRAM FINDS NO SOLUTION FOR VALUES OF E LESS THA    170
     1N E15.6///)                                                           180
4     FORMAT (1XI3,F12.5,3I5,F17.5,F20.5,F18.5)                             190
5     FORMAT (36H0THE UNRESOLVED WAVELENGTH USED WAS F8.5,41H AND THE K-    200
     1ALPHA-ONE WAVELENGTH USED WAS F8.5/16H THE DENSITY IS F6.3,22H THE    210
     2 FORMULA WEIGHT ISF8.3,21H, AND THE NUMBER IS  F6.2,1H.)              220
6     FORMAT (11H0THE FIRST I3,53H LINES ARE UNRESOLVED. THE REST ARE K-    230
     1ALPHA-1 LINES. )                                                      240
7     FORMAT (52H AN ESTIMATE OF THE STANDARD DEVIATION OF THETA IS  F10    250
     1.5/65H AN ESTIMATE OF THE STANDARD DEVIATION OF SINE SQUARED THETA    260
     2 IS  F10.5/51H AN ESTIMATE OF THE STANDARD DEVIATION OF AZERO IS F    270
     310.5/51H AN ESTIMATE OF THE STANDARD DEVIATION OF BZERO IS F10.5/5    280
     41H AN ESTIMATE OF THE STANDARD DEVIATION OF CZERO IS F10.5//)         290
8     FORMAT (31HOERROR. PROBABLY CUBIC CRYSTAL.//)                         300
9     FORMAT (1H1//24X20A4//)                                               310
      KKKK=0                                                                320
      I=40                                                                  330
      J=41                                                                  340
      K=42                                                                  350
      L=43                                                                  360
      M=44                                                                  370
      N=45                                                                  380
      LE=46                                                                 390
      ME=47                                                                 400
      NE=48                                                                 410
      KIT=0                                                                 420
      KAT=0                                                                 430
      IGO=1                                                                 440
      NOD=0                                                                 450
      MAY=NI                                                                460
      B=1.0                                                                 470
      DO 10 II=2,NI                                                         480
      B=AMIN1(B,ABS(S(II)-S(II-1)))                                         490
10    CONTINUE                                                              500
      E=.99*B                                                               510
      EMAX=E/2.                                                             520
      E=32.*EMAX                                                            530
C     OMMAND 230 AND FOLLOWING COMMANDS FIX THE MAXIMUM ERROR               540
      KAN=1                                                                 550
      KS=1                                                                  560
C     KS=1 WILL START THE VALUE OF K AT 1 SO THE PROGRAM DOES NOT LOOK      570
C     FOR THE 0,0,0 LINE.                                                   580
      JAN=KAN                                                               590
      NAN=JAN+1                                                             600
C     X,Y, AND Z WILL BE DETERMINED FROM REFLECTIONS KAN,JAN, AND NAN.      610
      JAN=JAN+1                                                             620
      NAN=NAN+1                                                             630
11    IF (KIT) 12,12,16                                                     640
C     IF KIT=1 WE SEEK ALL SOLUTIONS FOR A FIXED ERROR                      650
12    GO TO (13,14), IGO                                                    660
13    E=E/2.0                                                               670
      KSKIP=1                                                               680
      NSTORE=0                                                              690
      IF (E-ETEST) 19,19,14                                                 700
14    DO 62 IT=1,1255                                                       710
      CALL FETCH(IT,K9)                                                     720
      DO 15 ICK=1,9                                                         730
      SQ(ICK+39)=K9(ICK)**2                                                 740
15    CONTINUE                                                              750
      GO TO 21                                                              760
16    NOD=NOD+1                                                             770
      IF (NOD-20) 62,62,99                                                  780
17    KIT=KIT+1                                                             790
      IF (KIT-1) 18,18,62                                                   800
18    IGO=2                                                                 810
      KSKIP=2                                                               820
      GO TO 14                                                              830
19    IF (KIT) 20,20,17                                                     840
20    E=E*2.                                                                850
      GO TO 17                                                              860
C     THIS FIXES THE ERROR AND SEEKS ALL SOLUTIONS                          870
21    D=SQ(I)*(SQ(M)*SQ(NE)-SQ(N)*SQ(ME))-SQ(J)*(SQ(L)*SQ(NE)-SQ(N)*SQ(L    880
     1E))+SQ(K)*(SQ(L)*SQ(ME)-SQ(M)*SQ(LE))                                 890
C     D IS THE COEFFICIENT DETERMINANT                                      900
      IF (D) 39,22,39                                                       910
22    RANK=SQ(I)*SQ(L)+SQ(J)*SQ(M)+SQ(K)*SQ(N)                              920
      IF (RANK) 62,62,23                                                    930
23    RANK=SQ(I)*SQ(LE)+SQ(J)*SQ(ME)+SQ(K)*SQ(NE)                           940
      IF (RANK) 62,62,24                                                    950
C                                                                           960
C     OMMANDS 501 AND 502 REJECT ROWS WITH ALL ENTRIES ZERO.                970
C                                                                           980
24    U=S(KAN)*(SQ(M)*SQ(NE)-SQ(N)*SQ(ME))-S(JAN)*(SQ(J)*SQ(NE)-SQ(K)*SQ    990
     1(ME))+S(NAN)*(SQ(J)*SQ(N)-SQ(K)*SQ(M))                               1000
      IF (U) 62,25,62                                                      1010
25    V=SQ(I)*(S(JAN)*SQ(NE)-S(NAN)*SQ(N))-S(KAN)*(SQ(L)*SQ(NE)-SQ(N)*SQ   1020
     1(LE))+SQ(K)*(S(NAN)*SQ(L)-S(JAN)*SQ(LE))                             1030
      IF (V) 62,26,62                                                      1040
26    W=SQ(I)*(S(NAN)*SQ(M)-S(JAN)*SQ(ME))-SQ(J)*(S(NAN)*SQ(L)-S(JAN)*SQ   1050
     1(LE))+S(KAN)*(SQ(L)*SQ(ME)-SQ(M)*SQ(LE))                             1060
      IF (W) 62,29,62                                                      1070
C                                                                          1080
C     OMMANDS 505 THROUGH 507 REJECT INCONSISTENT EQUATIONS.               1090
C                                                                          1100
27    NAN=NAN+1                                                            1110
      IF (NAN-5) 14,14,28                                                  1120
28    JAN=JAN+1                                                            1130
      NAN=JAN+1                                                            1140
      IF (NAN-5) 14,14,96                                                  1150
29    IF (SQ(I)) 30,36,30                                                  1160
30    A22=SQ(M)-SQ(J)*SQ(L)/SQ(I)                                          1170
      A32=SQ(ME)-SQ(LE)*SQ(J)/SQ(I)                                        1180
      A23=SQ(N)-SQ(K)*SQ(L)/SQ(I)                                          1190
      A33=SQ(NE)-SQ(K)*SQ(LE)/SQ(I)                                        1200
      IF (A32**2+A33**2) 31,34,31                                          1210
31    IF (A22**2+A23**2) 32,35,32                                          1220
32    CONTINUE                                                             1230
      B=-A23/A22                                                           1240
      A=-(SQ(K)-B*SQ(J))/SQ(I)                                             1250
      GO TO 33                                                             1260
33    IF (ABS(A*S(1)+B*S(2)+1.*S(3))-0.001) 27,27,62                       1270
34    IF (ABS(SQ(LE)/SQ(I)-S(3)/S(1))-.001) 27,27,62                       1280
35    IF (ABS(SQ(L)/SQ(I)-S(2)/S(1))-.001) 27,27,62                        1290
36    IF (SQ(J)) 37,38,37                                                  1300
37    KOLDI=I                                                              1310
      KOLDL=L                                                              1320
      KOLDLE=LE                                                            1330
      I=J                                                                  1340
      L=M                                                                  1350
      LE=ME                                                                1360
      J=KOLDI                                                              1370
      M=KOLDL                                                              1380
      ME=KOLDLE                                                            1390
      GO TO 30                                                             1400
38    KOLDI=I                                                              1410
      KOLDL=L                                                              1420
      KOLDLE=LE                                                            1430
      I=K                                                                  1440
      L=N                                                                  1450
      LE=NE                                                                1460
      K=KOLDI                                                              1470
      N=KOLDL                                                              1480
      NE=KOLDLE                                                            1490
      GO TO 30                                                             1500
C                                                                          1510
C     OMMANDS 510 TO 526 PERFORM AN ANALYSIS OF ALL DEPENDENT CASES.       1520
C                                                                          1530
39    X=S(KAN)*(SQ(M)*SQ(NE)-SQ(N)*SQ(ME))-S(JAN)*(SQ(J)*SQ(NE)-SQ(K)*SQ   1540
     1(ME))+S(NAN)*(SQ(J)*SQ(N)-SQ(K)*SQ(M))                               1550
      X=X/D                                                                1560
      IF (X-.001) 62,62,40                                                 1570
40    Y=SQ(I)*(S(JAN)*SQ(NE)-S(NAN)*SQ(N))-S(KAN)*(SQ(L)*SQ(NE)-SQ(N)*SQ   1580
     1(LE))+SQ(K)*(S(NAN)*SQ(L)-S(JAN)*SQ(LE))                             1590
      Y=Y/D                                                                1600
      IF (Y-.001) 62,62,41                                                 1610
41    Z=SQ(I)*(S(NAN)*SQ(M)-S(JAN)*SQ(ME))-SQ(J)*(S(NAN)*SQ(L)-S(JAN)*SQ   1620
     1(LE))+S(KAN)*(SQ(L)*SQ(ME)-SQ(M)*SQ(LE))                             1630
      Z=Z/D                                                                1640
      IF (Z-.001) 62,62,42                                                 1650
42    JAY=0                                                                1660
43    JAY=JAY+1                                                            1670
C     JAY IS THE NUMBER OF THE REFLECTION BEING INDEXED                    1680
      NR=0                                                                 1690
      MD=0                                                                 1700
      OLDDE=1.1                                                            1710
44    NR=NR+1                                                              1720
      IF (NR-20) 45,45,49                                                  1730
45    NS=0                                                                 1740
46    NS=NS+1                                                              1750
      IF (NS-20) 47,47,44                                                  1760
47    NT=0                                                                 1770
48    NT=NT+1                                                              1780
      IF (NT-20) 50,50,46                                                  1790
49    IF (MD) 62,62,61                                                     1800
C                                                                          1810
C         MD NON-ZERO MEANS AT LEAST ONE ACCEPTABLE SET OF INDICES HAVE    1820
C         BEEN FOUND FOR REFLECTION JAY.                                   1830
50    C(JAY)=X*SQ(NR)+Y*SQ(NS)+Z*SQ(NT)                                    1840
C     %JAY< IS THE CALCULATED SINE SQUARED THETA                           1850
      IF (C(JAY)-S(JAY)-E) 54,54,51                                        1860
51    IF (X*SQ(NR)+Y*SQ(NS)-S(JAY)-E) 46,52,52                             1870
52    IF (X*SQ(NR)-S(JAY)-E) 44,53,53                                      1880
53    IF (MD) 62,62,61                                                     1890
C     WERE THERE SOLUTIONS FOR THIS X,Y, AND Z  NO,NO,YES                  1900
54    IF (C(JAY)-S(JAY)+E) 48,55,55                                        1910
C     OMMAND 557 ASKS IF THIS IS A SOLUTION NO,YES,YES.                    1920
55    DEE=ABS(C(JAY)-S(JAY))                                               1930
      IF (OLDDE-DEE) 48,48,56                                              1940
C     THIS COMMAND ASKS IF THIS SET OF INDICES IS BETTER THAN ALL          1950
C     PREVIOUS SOLUTIONS NO, NO, YES.                                      1960
56    IF (JAY-1) 57,57,58                                                  1970
57    OLDDE=DEE                                                            1980
      MD=1                                                                 1990
      B1=SQ(NR)                                                            2000
      B2=SQ(NS)                                                            2010
      B3=SQ(NT)                                                            2020
C     THE BI ARE THE CURRENT BEST SOLUTION FOR THE S(JAY) REFLECTION.      2030
      GO TO 48                                                             2040
58    IF (AH(JAY-1)-SQ(NR)) 57,59,57                                       2050
59    IF (AK(JAY-1)-SQ(NS)) 57,60,57                                       2060
60    IF (AL(JAY-1)-SQ(NT)) 57,48,57                                       2070
C                                                                          2080
C     OMMANDS 567,568, AND 569 REJECT DUPLICATE SETS OF INDICES.           2090
C                                                                          2100
61    AH(JAY)=B1                                                           2110
      AK(JAY)=B2                                                           2120
      AL(JAY)=B3                                                           2130
      MD=0                                                                 2140
      OLDDE=1.1                                                            2150
      C(JAY)=X*B1+Y*B2+Z*B3                                                2160
      IF (JAY-4) 43,64,65                                                  2170
62    CONTINUE                                                             2180
      IF (KKKK) 98,98,63                                                   2190
63    RETURN                                                               2200
64    AO=AH(1)*S(1)+AH(2)*S(2)+AH(3)*S(3)                                  2210
      BO=AH(1)**2+AH(2)**2+AH(3)**2                                        2220
      CO=AH(1)*AK(1)+AH(2)*AK(2)+AH(3)*AK(3)                               2230
      DO=AH(1)*AL(1)+AH(2)*AL(2)+AH(3)*AL(3)                               2240
      AA=AK(1)*S(1)+AK(2)*S(2)+AK(3)*S(3)                                  2250
      BB=CO                                                                2260
      CC=AK(1)**2+AK(2)**2+AK(3)**2                                        2270
      DD=AK(1)*AL(1)+AK(2)*AL(2)+AK(3)*AL(3)                               2280
      AT=AL(1)*S(1)+AL(2)*S(2)+AL(3)*S(3)                                  2290
      BT=DO                                                                2300
      CT=DD                                                                2310
      DT=AL(1)**2+AL(2)**2+AL(3)**2                                        2320
C                                                                          2330
C     OMMANDS FOLLOWING 575 DO THE LEAST SQUARES FIT FOR THE FIRST 3 LIN   2340
C                                                                          2350
65    AO=AO+AH(JAY)*S(JAY)                                                 2360
      BO=BO+AH(JAY)**2                                                     2370
      CO=CO+AH(JAY)*AK(JAY)                                                2380
      DO=DO+AH(JAY)*AL(JAY)                                                2390
      AA=AA+AK(JAY)*S(JAY)                                                 2400
      BB=CO                                                                2410
      CC=CC+AK(JAY)**2                                                     2420
      DD=DD+AK(JAY)*AL(JAY)                                                2430
      AT=AT+AL(JAY)*S(JAY)                                                 2440
      BT=DO                                                                2450
      CT=DD                                                                2460
      DT=DT+AL(JAY)**2                                                     2470
      DEN=BO*(CC*DT-DD*CT)-CO*(BB*DT-DD*BT)+DO*(BB*CT-CC*BT)               2480
      AZERO=(AO*(CC*DT-CT*DD)-CO*(AA*DT-AT*DD)+DO*(AA*CT-AT*CC))/DEN       2490
      BZERO=(BO*(AA*DT-AT*DD)-AO*(BB*DT-BT*DD)+DO*(BB*AT-BT*AA))/DEN       2500
      CZERO=(BO*(CC*AT-CT*AA)-CO*(BB*AT-BT*AA)+AO*(BB*CT-BT*CC))/DEN       2510
      IF (JAY-NI) 66,66,67                                                 2520
66    X=AZERO                                                              2530
      Y=BZERO                                                              2540
      Z=CZERO                                                              2550
      GO TO 43                                                             2560
67    IF (AZERO) 69,68,69                                                  2570
68    AZERO=YLAMDA/2.                                                      2580
      GO TO 70                                                             2590
69    AZERO=(YLAMDA/2.)/(SQRT(AZERO))                                      2600
70    IF (BZERO) 72,71,72                                                  2610
71    BZERO=YLAMDA/2.                                                      2620
      GO TO 73                                                             2630
72    BZERO=(YLAMDA/2.)/SQRT(BZERO)                                        2640
73    IF (CZERO) 75,74,75                                                  2650
74    CZERO=YLAMDA/2.                                                      2660
      GO TO 76                                                             2670
75    CZERO=(YLAMDA/2.)/SQRT(CZERO)                                        2680
C     OMMANDS 575 THROUGH 586 CALCULATE AZERO,BZERO, AND CZERO.            2690
76    GO TO (87,77), KSKIP                                                 2700
77    EX=AZERO                                                             2710
      KEX=0                                                                2720
78    DO 82 IC=1,20                                                        2730
      DO 81 IB=1,10                                                        2740
      ZB=IB                                                                2750
      IF (EX-ZB*STOREA(IC)-0.001) 83,83,79                                 2760
79    IF (EX-ZB*STOREB(IC)-0.001) 83,83,80                                 2770
80    IF (EX-ZB*STOREC(IC)-0.001) 83,83,81                                 2780
81    CONTINUE                                                             2790
82    CONTINUE                                                             2800
      GO TO 86                                                             2810
83    KEX=KEX+1                                                            2820
      IF (KEX-2) 84,85,62                                                  2830
84    EX=BZERO                                                             2840
      GO TO 78                                                             2850
85    EX=CZERO                                                             2860
      GO TO 78                                                             2870
86    NSTORE=NSTORE+1                                                      2880
      STOREA(NSTORE)=AZERO                                                 2890
      STOREB(NSTORE)=BZERO                                                 2900
      STOREC(NSTORE)=CZERO                                                 2910
87    CONTINUE                                                             2920
      Q=(YLAMDA/2.)**2                                                     2930
      DO 88 IC=1,NI                                                        2940
      ESS(IC)=Q*(AH(IC)/AZERO**2+AK(IC)/BZERO**2+AL(IC)/CZERO**2)          2950
88    CONTINUE                                                             2960
C     OMMAND 590 FINDS CALCULATED SINE SQUARED THETA                       2970
      SUMA=0.                                                              2980
      SUM=0.                                                               2990
      DO 89 IC=1,NI                                                        3000
      DE(IC)=S(IC)-ESS(IC)                                                 3010
      KAH(IC)=AH(IC)                                                       3020
      KAK(IC)=AK(IC)                                                       3030
      KAL(IC)=AL(IC)                                                       3040
      SUM=SUM+DE(IC)**2                                                    3050
      SUMA=SUMA+(DE(IC)**2)/S2SQ(IC)                                       3060
89    CONTINUE                                                             3070
      SIGMA=SQRT(SUM/FLOAT(NI-3))                                          3080
      SIGTH=SQRT(SUMA/FLOAT(NI-3))*57.2957795                              3090
      SIGMAX=SIGMA*SQRT((CC*DT-CT*DD)/DEN)                                 3100
      SIGMAY=SIGMA*SQRT((BO*DT-BT*DO)/DEN)                                 3110
      SIGMAZ=SIGMA*SQRT((BO*CC-BB*CO)/DEN)                                 3120
      SIGMAA=SIGMAX*AZERO/(2.*X)                                           3130
      SIGMAB=SIGMAY*BZERO/(2.*Y)                                           3140
      SIGMAC=SIGMAZ*CZERO/(2.*Z)                                           3150
      WRITE(IPRNT,2)NAME,E,AZERO,BZERO,CZERO                               3160
      VOLUME=AZERO*BZERO*CZERO                                             3170
      IF (FORMWT) 91,91,90                                                 3180
90    ZNMBER=DENSTY*VOLUME*.6023/FORMWT                                    3190
      GO TO 92                                                             3200
91    ZNMBER=0.                                                            3210
92    CONTINUE                                                             3220
      WRITE(IPRNT,5)XLAMDA,YLAMDA,DENSTY,FORMWT,ZNMBER                     3230
      WRITE(IPRNT,6)NCUT                                                   3240
      WRITE(IPRNT,7)SIGTH,SIGMA,SIGMAA,SIGMAB,SIGMAC                       3250
      WRITE(IPRNT,1)                                                       3260
      DO 95 IC=1,NI                                                        3270
      IF (MOD(IC+6,32)) 94,93,94                                           3280
93    WRITE(IPRNT,9)NAME                                                   3290
      WRITE(IPRNT,1)                                                       3300
94    WRITE(IPRNT,4)(IC,TWOTH(IC),KAH(IC),KAK(IC),KAL(IC),S(IC),ESS(IC),   3310
     1DE(IC))                                                              3320
95    CONTINUE                                                             3330
      KKKK=KKKK+1                                                          3340
      GO TO 11                                                             3350
96    IF (KKKK) 97,97,99                                                   3360
97    WRITE(IPRNT,9)NAME                                                   3370
      WRITE(IPRNT,8)                                                       3380
      RETURN                                                               3390
98    WRITE(IPRNT,9)NAME                                                   3400
      WRITE(IPRNT,3)E                                                      3410
99    RETURN                                                               3420
      END                                                                  3430
      SUBROUTINE FETCH (N,L)                                                 10
      COMMON INT(1254),EL(60),TWOTH(150),NAME(20),S2SQ(150),ES(600),S(15     20
     10),NCUT,XLAMDA,YLAMDA,NI,ETEST,TETBLK,DENSTY,FORMWT,NOCUBE,NOHEX,N     30
     2OTET,NONO                                                              40
      COMMON/IOUNIT/ IREAD,  IPRNT
      DIMENSION L(9)                                                         50
      IF (N-1254) 3,3,1                                                      60
1     DO 2 I=1,9                                                             70
2     L(I)=0                                                                 80
      RETURN                                                                 90
3     IL=INT(N)                                                             100
      IL1=IL/262144                                                         110
      IL2=IL/512-IL1*512                                                    120
      IL3=IL-IL1*262144-IL2*512                                             130
      L(1)=IL1/64                                                           140
      L(2)=IL1/8-L(1)*8                                                     150
      L(3)=IL1-L(1)*64-L(2)*8                                               160
      L(4)=IL2/64                                                           170
      L(5)=IL2/8-L(4)*8                                                     180
      L(6)=IL2-L(4)*64-L(5)*8                                               190
      L(7)=IL3/64                                                           200
      L(8)=IL3/8-L(7)*8                                                     210
      L(9)=IL3-L(7)*64-L(8)*8                                               220
      DO 4 I=1,9                                                            230
4     L(I)=L(I)-1                                                           240
      RETURN                                                                250
      END                                                                   260
