      PROGRAM SCORE( A2,TAPE2=A2, A3,TAPE3=A3 )
      CODE ANALYSIS
      CLICHE COMBLR
C
      PARAMETER (NLIM=111)
C
      COMMON /ICOM/ QERR, WL1, WLBAR, TITLE(8), CC(6), SG,
     , D(NLIM), Q(NLIM), SINSQ(NLIM), TWOTH(NLIM)
      COMMON /CCOM/ NASP, KS, QMAX, QMIN, NPOS, NINPUT, X(6),
     , QC(5*NLIM), QM(NLIM)
C
      DIMENSION BLR(NLIM, 4), QO(NLIM)
C
      EQUIVALENCE (BLR,D), (Q,QO), (A,CC(1)), (B,CC(2)), (C,CC(3)),
     , (ALPHA,CC(4)), (BETA,CC(5)), (GAMMA,CC(6)), (QC,MI)
C
      INTEGER H, MI(5*NLIM), SG
C
      DATA (IT=2), (NT=3), (ZERO=0.), (MASK=77777777777777000000B),
     , (DTOR=1.74532925199E-2), (RTOD=57.2957795131)
C
      NAMELIST /INPUTS/ QERR, WL1, WLBAR, CC, SG, D, Q, SINSQ, TWOTH
C
      ENDCLICHE
C
C
C
C
      USE COMBLR
C
      COMMON /SNYDER/ IOUT, RCP(6), RCV, SDQ, SDTT, SDT, VOL, WAVE
      INTEGER OLDSUM, OLDL
C
C.... CREATE DROP FILE
      CALL CHANGE( 6H+SCORE )
      L = 32768  $  CALL CREATE( 2HA3, L, IOC3 )
C.... CLEAR STORAGE
   10 CONTINUE
      CALL QACOPY( QERR, ZERO, 18+4*NLIM, 1, 0 )
C.... READ INPUT
C.... DEFAULT            .001      1.5406     B=A, C=A, ANGLES=90.
C.... INPUT? ALLOWABLE Q ERROR, WAVE LENGTH, CELL CONSTANTS, SPACE GROUP
C....        AND ONE OF THE FOLLOWING- D, Q, SINE SQUARED OR TWO THETA
      QERR = .001
      WL1 = 1.5406
      CALL RDINPT
C.... IDENTIFY SPACE GROUP BY INDEX
      CALL SGID( KS, SG, NASP )
C.... CHECK IS ASPECT HAS BEEN FOUND
      IF ( NASP .EQ. 0 ) GO TO 10
      IF ( KS - 5 ) 100, 500, 600
C.... CALCULATE CUBIC OR TETRAGONAL OR HEXAGONAL OR ORTHORHOMBIC INDICES
C....              FOR ASPECT = NASP
  100 CONTINUE
      CALL CTHO
      GO TO 700
C.... CALCULATE MONOCLINIC INDICES FOR ASPECT = NASP
  500 CONTINUE
      CALL MONO
      GO TO 700
C.... ANORTHIC NOT IMPLEMENTED
  600 CONTINUE
      CALL EXIT( 1 )
C.... ATTEMPT TO MATCH CALCULATED AND OBSERVED Q
  700 CONTINUE
      NQC = NPOS
      NMAX = NINPUT
      CALL QMATCH( NQC, NMAX )
      IMAX = NPOS + NINPUT
      WRITE (NT,701) NASP, TITLE
  701 FORMAT (10H ASPECT IS, I4, ////, 1X8A10, //)
      WAVE = WL1
      IOUT = NT
      SDQ = SDT = SDTT = 0.
C
CALCULATE RECIPROCAL LATTICE PARAMETERS AND VOLUME
      CALL CELCON( CC )
      WRITE (NT,702) CC, RCV, WAVE, RCP, VOL
  702 FORMAT(4X1HA, 8X1HB, 8X1HC, 6X5HCOS A, 2X5HCOS B, 2X5HCOS G,
     1 4X3HVOL, 5X4HWAVE / 3(1XF8.5), 3(1XF6.4), 1XF8.2, 1XF7.5 /
     2 3(1XF8.6), 3(1XF6.4), 1XF8.6///
     3 2X,~H K L NO NP  2T OBS 2T CALC  AV2T RMS2T  FN AV FN RMS   M(N)
     4  RATIO~)
C
      N = 1
      IMAX = IMAX + 1
      NPOS = 0
      NMINUS = 0
      DO 770 I=2,IMAX
        IF ( QC(I) .EQ. QM(N) ) GO TO 740
        IF ( QC(I) .EQ. QO(N) ) GO TO 750
C.... CHECK IF CALCULATED Q EXCEEDS MAXIMUM Q
        IF ( QC(I).NE.QM(NINPUT) .AND. QC(I).GT.QO(NINPUT) ) GO TO 780
C.... A CALCULATED Q IS NOT MATCHED
        DC = SQRT( 1.0 / QC(I) )
        H = ( QC(I) .INT. 770000B ) / 10000B
        K = ( QC(I) .INT. 7700B ) / 100B
        L = QC(I) .INT. 37B
        IF ( QC(I).INT.40B .NE. 0 ) H = -H
        IF ( KS .GE. 4 ) GO TO 720
C.... CHECK FOR EQUIVALENT REFLECTIONS IF CUBIC, TETRAGONAL OR HEXAGONAL
        NEWSUM = H*H + K*K
        IF ( KS .EQ. 1 ) NEWSUM = NEWSUM + L*L
        IF ( KS .EQ. 3 ) NEWSUM = NEWSUM + H*K
        NEWL = L
        IF ( KS .EQ. 1 ) NEWL = 0
        IF ( NEWSUM.EQ.OLDSUM .AND. NEWL.EQ.OLDL ) GO TO 730
        OLDL = NEWL
        OLDSUM = NEWSUM
C.... A CALCULATED Q IS UNMATCHED AND COUNTABLE
  720   CONTINUE
        NPOS = NPOS + 1
C.... EDIT UNMATCHED CALCULATED Q
  730   CONTINUE
        TWOTHC = 2.0 * RTOD * ASIN( .5 * WL1 / DC )
        WRITE (NT,771) H, K, L, NPOS, TWOTHC
        GO TO 770
C.... A CALCULATED Q MATCHES AN OBSERVED Q
  740   CONTINUE
        IF ( QC(I+1) .NE. QO(N) ) N = N + 1
        GO TO 770
C.... AN OBSERVED Q IS MATCHED
  750   CONTINUE
        IF ( QM(N) .EQ. 0. ) GO TO 760
        DC = SQRT( 1.0 / QM(N) )
        H = ( QM(N) .INT. 770000B ) / 10000B
        K = ( QM(N) .INT. 7700B ) / 100B
        L = QM(N) .INT. 37B
        IF ( QM(N).INT.40B .NE. 0 ) H = -H
        IF ( KS .GE. 4 ) GO TO 752
C.... CHECK FOR EQUIVALENT REFLECTIONS IF CUBIC, TETRAGONAL OR HEXAGONAL
        NEWSUM = H*H + K*K
        IF ( KS .EQ. 1 ) NEWSUM = NEWSUM + L*L
        IF ( KS .EQ. 3 ) NEWSUM = NEWSUM + H*K
        NEWL = L
        IF ( KS .EQ. 1 ) NEWL = 0
        IF ( NEWSUM.EQ.OLDSUM .AND. NEWL.EQ.OLDL ) GO TO 754
        OLDL = NEWL
        OLDSUM = NEWSUM
C.... A CALCULATED Q IS MATCHED AND COUNTABLE
  752   CONTINUE
        NPOS = NPOS + 1
  754   CONTINUE
C.... CALCULATE SNYDER FIGURE OF MERIT
        CALL FIG( H, K, L, D(N), NPOS, N-NMINUS )
        IF ( QC(I+1) .NE. QM(N) ) N = N + 1
        GO TO 770
C.... EDIT UNMATCHED OBSERVED Q
  760   CONTINUE
        TWOTHO = 2.0 * RTOD * ASIN( .5 * WL1 / D(N) )
        WRITE (NT,772) TWOTHO
        IF ( N .LT. NINPUT ) N = N + 1
        IF ( N .LT. NINPUT ) NMINUS = NMINUS + 1
  770 CONTINUE
  771 FORMAT (1X 3I2, 3X I3, 9X F7.3)
  772 FORMAT (13X, F8.3)
  780 CONTINUE
  781 FORMAT (///, 7H NPOS =, I4)
      WRITE (NT,781) NPOS
      GO TO 10
 
 
      END
      SUBROUTINE ASPKS1( H, K, L, NASP, KEEP )
C.... CUBIC
      CODE ANALYSIS
C
      PARAMETER (NN=17)
      INTEGER ASPECT(NN), H, HKL
C
C
C....   7654321, 7654321, 7654321, 7654321, 7654321, 7654321, 7654321
C....    POSITION 1 IS CRITERION FOR 00L TYPE
C....    POSITION 2                  0K0 TYPE BECOMES 00L TYPE
C....    POSITION 3 IS CRITERION FOR 0KL TYPE
C....    POSITION 4                  H00 TYPE BECOMES 00L TYPE
C....    POSITION 5                  H0L TYPE BECOMES 0KL TYPE
C....    POSITION 6                  HK0 TYPE BECOMES 0KL TYPE
C....    POSITION 7 IS CRITERION FOR HKL TYPE, H = K
C....    CRITERION 1 IS P - NO RESTRICTION
C....    CRITERION 2 IS I - SUM EVEN
C....    CRITERION 3 IS F - ALL ODD OR ALL EVEN
C....    CRITERION 4 IS J - HALF-SUM EVEN
C....    CRITERION 5 IS Q - ALL EVEN
C....    CRITERION 6 IS G - ALL EVEN AND HALF-SUM EVEN
C....    CRITERION 7 IS Z - TYPE 0KL WITH K EVEN
      DATA ( ASPECT =                       1000101, 1000102, 1000104,
     ,  1000202, 2000102, 2000202, 1000702, 2000202, 2000204, 2000302,
     ,  4000204, 4000304, 3000302, 3000304, 3000504, 5000302, 5000604 )
C
C
C
      I = NASP
C.... INITIALIZE FOR HKL TYPE WITH ALL INDICES DIFFERENT.
C.... I
      KODE = 2
C.... P
      IF ( NASP .LE. 7 ) KODE = 1
C.... F
      IF ( NASP .GE. 13 ) KODE = 3
  110 CONTINUE
      HKL = 4 * MIN0( H, 1 )  +  2 * MIN0( K, 1 )  +  MIN0( L, 1 )
C             00L, 0K0, 0KL, H00, H0L, HK0, HKL
      GO TO ( 124, 116, 118, 114, 112, 114, 120 ), HKL
C.... EXCHANGE H AND K IF H0L
  112 J = H  $  H = K  $  K = J  $  GO TO 110
C.... EXCHANGE H AND L IF H00, HK0 OR HKL (L .EQ. K .NE. H)
  114 J = H  $  H = L  $  L = J  $  GO TO 110
C.... EXCHANGE K AND L IF 0K0 OR HKL (H .EQ. L .NE. K)
  116 J = K  $  K = L  $  L = J  $  GO TO 110
C.... 0KL TYPE
  118 CONTINUE
      IF ( MOD( K, 2 ) .EQ. 0 ) GO TO 122
      IF ( MOD( L, 2 ) ) 122, 116, 122
C.... HKL TYPE - CHECK IF HHL
  120 CONTINUE
      IF ( H .EQ. K ) GO TO 122
      IF ( H .EQ. L ) GO TO 116
      IF ( K - L ) 124, 114, 124
C.... EXTRACT CRITERION
  122 CONTINUE
      KODE = ASPECT(I) / 10**(HKL-1)
      KODE = MOD( KODE, 10 )
  124 CONTINUE
      NODD = MOD( H, 2 )  + MOD( K, 2 )  +  MOD( L, 2 )
      NEVEN = 3 - NODD
C....           P    I    F    J    Q    G    Z
      GO TO ( 136, 126, 128, 126, 130, 132, 134 ), KODE
C.... I OR J
  126 CONTINUE
      IF ( MOD( H+K+L, KODE ) ) 138, 136, 138
C.... F  ( ALL ODD OR ALL EVEN )
  128 CONTINUE
      IF ( NODD * NEVEN ) 138, 136, 138
C.... Q  ( ALL EVEN )
  130 CONTINUE
      IF ( NODD ) 138, 136, 138
C.... G  ( ALL EVEN, HALF-SUM EVEN )
  132 CONTINUE
      IF ( MOD( H+K+L, 4 ) ) 138, 130, 138
C.... Z  ( 0KL, K EVEN )
  134 CONTINUE
      IF ( MOD( K, 2 ) ) 138, 136, 138
  136 CONTINUE
      KEEP = 1
  138 CONTINUE
      RETURN
      END
      SUBROUTINE ASPKS2( H, K, L, NASP, KEEP )
C.... TETRAGONAL
      CODE ANALYSIS
C
      PARAMETER (NN=31)
      INTEGER ASPECT(NN), H, HKL
C
C
C....   7654321, 7654321, 7654321, 7654321, 7654321, 7654321, 7654321
C....    POSITION 1 IS CRITERION FOR 00L TYPE
C....    POSITION 2 IS CRITERION FOR 0K0 TYPE
C....    POSITION 3 IS CRITERION FOR 0KL TYPE
C....    POSITION 4                  H00 TYPE BECOMES 0K0 TYPE
C....    POSITION 5                  H0L TYPE BECOMES 0KL TYPE
C....    POSITION 6 IS CRITERION FOR HK0 TYPE, H = K
C....    POSITION 7 IS CRITERION FOR HKL TYPE, H = K
C....    CRITERION 1 IS P - NO RESTRICTIONS
C....    CRITERION 2 IS I - SUM EVEN
C....    CRITERION 3 IS C - (H+K) EVEN
C....    CRITERION 4 IS J - HALF-SUM EVEN
C....    CRITERION 5 IS B - (H+L) EVEN
C....    CRITERION 6 IS F - ALL ODD OR ALL EVEN
      DATA ( ASPECT =                       1100111, 1100112, 1100114,
     ,  1100121, 1100122, 1100124, 2100112, 2100122, 1100321, 2100322,
     ,  1100512, 2100512, 1100222, 2100222, 1200121, 1200122, 2200122,
     ,  1200321, 2200322, 1200522, 2200522, 1200222, 2200222, 2200222,
     ,  2200224, 2200622, 4400224, 4400624, 2600224, 4600224, 4600624 )
C
C
C
      I = NASP
C.... EXCHANGE H AND K IF H00 OR H0L
      IF ( K .EQ. 0 ) GO TO 110
      J = H  $  H = K  $  K = J
  110 CONTINUE
      HKL = 4 * MIN0( H, 1 )  +  2 * MIN0( K, 1 )  +  MIN0( L, 1 )
      KODE = ASPECT(I) / 10**(HKL-1)
      KODE = MOD( KODE, 10 )
      IF ( HKL .LE. 3 ) GO TO 114
C.... IF H = K HK0 AND HKL TYPES SHOULD BE ALL RIGHT
      IF ( H .EQ. K ) GO TO 114
C.... CORRECT HK0 TYPE IF NASP = 27 OR 28
      IF ( KODE .EQ. 4 ) KODE = 2
      IF ( HKL .EQ. 6 ) GO TO 114
C.... HKL TYPE, H .NE. K
      KODE = 1
      IF ( NASP .GE. 24 ) KODE = 2
  114 CONTINUE
C....           P    I    C    J    B    F
      GO TO ( 124, 116, 118, 116, 120, 122 ), KODE
C.... I OR J
  116 CONTINUE
      IF ( MOD( H+K+L, KODE ) ) 126, 124, 126
C.... C  ( H+K EVEN )
  118 CONTINUE
      IF ( MOD( H+K, 2 ) ) 126, 124, 126
C.... B  ( L+H EVEN )
  120 CONTINUE
      IF ( MOD( L+H, 2 ) ) 126, 124, 126
C.... F  ( ALL ODD OR ALL EVEN )
  122 CONTINUE
      NODD = MOD( H, 2 )  +  MOD( K, 2 )  +  MOD( L, 2 )
      NEVEN = 3 - NODD
      IF ( NODD * NEVEN ) 126, 124, 126
  124 CONTINUE
      KEEP = 1
  126 CONTINUE
      RETURN
      END
      SUBROUTINE ASPKS3( H, K, L, NASP, KEEP )
C.... HEXAGONAL
      CODE ANALYSIS
C
      PARAMETER (NN=9)
      INTEGER ASPECT(NN), H, HKL
C
C
C....   7654321, 7654321, 7654321, 7654321, 7654321, 7654321, 7654321
C        POSITION 1 IS CRITERION FOR 00L TYPE
C....    POSITION 2 IS CRITERION FOR 0K0 TYPE ( H0L TYPE )
C....    POSITION 3 IS CRITERION FOR 0KL TYPE ( H0L TYPE )
C....    POSITION 4 IS CRITERION FOR H00 TYPE ( H0L TYPE )
C....    POSITION 5 IS CRITERION FOR H0L TYPE
C....    POSITION 6 IS CRITERION FOR HK0 TYPE ( HKL TYPE )
C....    POSITION 7 IS CRITERION FOR HKL TYPE, H = K
C....    CRITERION 1 IS P - NO RESTRICTIONS
C....    CRITERION 2 IS L - L EVEN
C....    CRITERION 3 FOR L DIVISIBLE BY 3
C....    CRITERION 4 IS R - (-H+K+L) OR (H-K+L) DIVISIBLE BY 3
C....    CRITERION 5 IS S - L AND R RESTRICTIONS
C....    CRITERION 6 FOR L DIVISIBLE BY 6
      DATA       ( ASPECT =                          1111111, 1111112,
     ,  1111113, 1111116, 1122222, 2211112, 2222222, 4444443, 4455556 )
C
C
C
      I = NASP
      IF ( HKL .LT. 7 ) GO TO 112
      IF ( H .EQ. K ) GO TO 112
C.... CORRECT ASPECTS 6 AND 7 FOR HKL, H .NE. K
      IF ( KODE .EQ. 2 ) KODE = 1
  112 CONTINUE
      HKL = 4 * MIN0( H, 1 )  +  2 * MIN0( K, 1 )  +  MIN0( L, 1 )
      KODE = ASPECT(I) / 10**(HKL-1)
      KODE = MOD( KODE, 10 )
C               P,   L, L/3,   R,   S, L/6
      GO TO ( 122, 114, 114, 116, 118, 114 ), KODE
  114 CONTINUE
C.... L DIVISIBLE BY 2, 4 OR 6
      IF ( MOD( L, KODE ) ) 124, 122, 124
C.... R - (-H+K+L) OR (H-K+L) DIVISIBLE BY 3
  116 CONTINUE
      IF ( MOD( -H+K+L, 3 )  .EQ.  0 ) GO TO 122
      IF ( MOD( H-K+L, 3 ) ) 124, 122, 124
C.... S - R AND L
  118 CONTINUE
      IF ( MOD( L, 2 ) ) 124, 118, 124
  122 CONTINUE
      KEEP = 1
  124 CONTINUE
      RETURN
      END
      SUBROUTINE ASPKS4( H, K, L, NASP, KEEP )
C.... ORTHORHOMBIC
      CODE ANALYSIS
C
      INTEGER A(198), ASPECT(33,6), H, HKL, HKSUM, HLSUM, HKLSUM
      EQUIVALENCE (A,ASPECT)
C
C....  7654321, 7654321, 7654321, 7654321, 7654321, 7654321, 7654321
C....    POSITION 1 IS CRITERION FOR 00L TYPE
C....    POSITION 2 IS CRITERION FOR 0K0 TYPE
C....    POSITION 3 IS CRITERION FOR 0KL TYPE
C....    POSITION 4 IS CRITERION FOR H00 TYPE
C....    POSITION 5 IS CRITERION FOR H0L TYPE
C....    POSITION 6 IS CRITERION FOR HK0 TYPE
C....    POSITION 7 IS CRITERION FOR HKL TYPE
C....    CRITERION 1 IS P - NO RESTRICTIONS
C....    CRITERION 2 IS I - SUM EVEN
C....    CRITERION 3 IS A - (K+L) EVEN
C....    CRITERION 4 IS G - ALL EVEN AND HALF-SUM EVEN
C....    CRITERION 5 IS B - (L+H) EVEN
C....    CRITERION 6 IS C - (H+K) EVEN
C....    CRITERION 7 IS F - ALL ODD OR ALL EVEN
C.... ORIENTATION A  -  (ASPECT(I,1), I=1,33)
      DATA ( (A(I),I=1,33) =
     1                   1111111, 1111112, 1112121, 1112122, 1512111,
     2 1122112, 1131512, 1162512, 1162621, 1131222, 1162222, 1122222,
     3 1532512, 1532622, 1232522, 1262621, 1232622, 1522222, 1222222,
     4 6262621, 6262622, 6272622, 3331722, 3372722, 6272722, 6772722,
     5 2222222, 2272222, 2272722, 2772722, 7772722, 7744444, 7444444 )
C.... ORIENTATION B  -  (ASPECT(I,2), I=1,33)
      DATA ( (A(I),I=34,66) =
     1                         1,       1,       1,       1, 1311121,
     2 1111222,       1, 1131622,       1, 1122512, 1122622,       1,
     3 1331522, 1362522,       1,       1, 1262522, 1322222,       1,
     4       1,       1, 6262722, 5572512, 5572722,       1,       1,
     5       1, 2222722,       1,       1,       1,       1,       1 )
C.... ORIENTATION C  -  (ASPECT(I,3), I=1,33)
      DATA ( (A(I),I=67,99) =
     1                         1, 1111121, 1112112,       1, 1131112,
     2       2, 1311621, 1321622, 1512512, 1212621, 1212522, 1212222,
     3 1331622,       1, 1322622, 1522512, 1522622, 1232222,       1,
     4 5522512, 5522522, 5522722, 6762621, 6762722, 5722722, 5772722,
     5       1,       2, 2722722,       1,       1, 7474444,       1 )
C.... ORIENTATION D  -  (ASPECT(I,4), I=1,33)
      DATA ( (A(I),I=100,132) =
     1                         1, 1112111, 1111122,       1, 1111512,
     2       1, 1562111, 1532112, 1331122, 1262121, 1232122, 1222122,
     3 1562512,       2, 1562222, 1331222, 1362222, 1222522,       1,
     4 3331222, 3332222, 3372222,       3, 6772622, 3772222, 3772722,
     5       1,       1, 2772222,       1,       1, 7444744,       1 )
C.... ORIENTATION E  -  (ASPECT(I,5), I=1,33)
      DATA ( (A(I),I=133,165) =
     1                         1,       4,       4,       1, 1111621,
     2 1212121,       4, 1362121,       4, 1522112, 1322122,       4,
     3 1562621,       1,       4,       4, 1532222, 1222622,       1,
     4       4,       4, 3732222,       2, 5772522,       4,       4,
     5       1, 2722222,       4,       1,       1,       4,       1 )
C.... ORIENTATION F  -  (ASPECT(I,6), I=1,33)
      DATA ( (A(I),I=166,198) =
     1                         1,       3,       3,       1, 1162111,
     2       5,       3, 1512621,       3, 1311222, 1512222,       3,
     3 1362621,       2,       3,       3, 1322522, 1262222,       1,
     4       3,       3, 5722522,       1, 3732722,       3,       3,
     5       1,       5,       3,       1,       1,       3,       1 )
C
C
C
C.... ORIENTATION IS IN HUNDRED POSITION
      JO = NASP/100
C.... ASPECT NUMBER IS IN LAST TWO DIGITS
      I = NASP - 100*JO
      HKL = 4 * MIN0( H, 1 )  +  2 * MIN0( K, 1 )  +  MIN0( L, 1 )
      NODD = MOD( H, 2 )  +  MOD( K, 2 )  +  MOD( L, 2 )
        KODE = ASPECT(I,JO) / 10**(HKL-1)
        IF ( ASPECT(I,JO) .LE. 5 ) GO TO 130
        KODE = MOD( KODE, 10 )
C....             P    I    A    G    B    C    F
        GO TO ( 119, 112, 113, 114, 115, 116, 117 ), KODE
C....   I  ( H+K+L EVEN )
  112   CONTINUE
        HKLSUM = H + K + L
        IF ( MOD( HKLSUM, KODE ) ) 120, 119, 120
C....   A  ( K+L EVEN )
  113   CONTINUE
        KLSUM = K + L
        IF ( MOD( KLSUM, 2 ) ) 120, 119, 120
C....   G  ( ALL EVEN, H+K+L DIVISIBLE BY 4 )
  114   CONTINUE
        IF ( NODD ) 120, 112, 120
C....   B  ( L+H EVEN )
  115   CONTINUE
        HLSUM = H + L
        IF ( MOD( HLSUM, 2 ) ) 120, 119, 120
C....   C  ( H+K EVEN )
  116   CONTINUE
        HKSUM = H + K
        IF ( MOD( HKSUM, 2 ) ) 120, 119, 120
C....   F  ( ALL ODD OR ALL EVEN )
  117   CONTINUE
        NEVEN = 3 - NODD
        IF ( NEVEN*NODD .NE. 0 ) GO TO 120
  119   CONTINUE
      KEEP = 1
  120 CONTINUE
      RETURN
  130 CONTINUE
      WRITE (59,131) NASP
  131 FORMAT (24H CANNOT PROVIDE ASPECT =, I4, 21H FOR CRYSTAL SYSTEM 4)
      END
      SUBROUTINE ASPKS5( H, K, L, I, KEEP )
C.... MONOCLINIC
      CODE ANALYSIS
C
      PARAMETER (NN=14)
      INTEGER ASPECT(NN), H, HKL, HKSUM, HLSUM, HKLSUM
C....   7654321, 7654321, 7654321, 7654321, 7654321, 7654321
C....    POSITION 1 IS CRITERION FOR 00L TYPE ( H0L TYPE )
C....    POSITION 2 IS CRITERION FOR 0K0 TYPE
C....    POSITION 3 IS CRITERION FOR 0KL TYPE ( HKL TYPE )
C....    POSITION 4 IS CRITERION FOR H00 TYPE ( H0L TYPE )
C....    POSITION 5 IS CRITERION FOR H0L TYPE
C....    POSITION 6 IS CRITERION FOR HK0 TYPE ( HKL TYPE )
C....    POSITION 7 IS CRITERION FOR HKL TYPE
      DATA ( ASPECT =                       1111111, 1111121,
     ,  1133113, 1155115, 1122112, 1133123, 1155125, 1122122,
     ,  5555525, 3333333, 2222222, 5544524, 3344324, 2244224 )
C....    CRITERION 1 IS P - NO RESTRICTUINS
C....    CRITERION 2 IS I - SUM EVEN
C....    CRITERION 3 IS A - (K+L) EVEN
C....    CRITERION 4 IS F - ALL ODD OR ALL EVEN
C....    CRITERION 5 IS C - (H+K) EVEN
C
C
C
        HKL = 4*MIN0( IABS(H), 1 ) + 2*MIN0( K, 1 ) + MIN0( L, 1 )
        KODE = ASPECT(I) / 10**(HKL-1)
        KODE = MOD( KODE, 10 )
C....             P    I    A    F    C
        GO TO ( 201, 112, 113, 114, 115 ), KODE
C....   I  ( SUM EVEN )
  112   CONTINUE
        HKLSUM = H + K + L
        IF ( MOD( HKLSUM, 2 ) ) 300, 201, 300
C....   A  ( K+L EVEN )
  113   CONTINUE
        KLSUM = K + L
        IF ( MOD( KLSUM, 2 ) ) 300, 201, 300
C....   F  ( ALL ODD OR ALL EVEN )
  114   CONTINUE
        NODD = MOD( IABS( H ), 2 )  +  MOD( K, 2 ) + MOD( L, 2 )
        NEVEN = 3 - NODD
      IF ( NODD*NEVEN ) 300, 201, 300
C....   C  ( H+K EVEN )
  115   CONTINUE
        HKSUM = H + K
        IF ( MOD( HKSUM, 2 ) ) 300, 201, 300
  201 CONTINUE
      KEEP = 1
  300 CONTINUE
      RETURN
      END
      SUBROUTINE CELCON( CC )
      CODE ANALYSIS
C
C      CALCULATE RECIPROCAL CELL CONSTANTS
C
      COMMON /SNYDER/ IOUT, RCP(6), RCV, SDQ, SDTT, SDT, VOL, WAVE
      DIMENSION CC(6)
      DATA (DTR=.01745329)
C
C
      A = CC(1)
      B = CC(2)
      C = CC(3)
      COSAL = CC(4)
      COSBE = CC(5)
      COSGA = CC(6)
C
      COSAL = COS(COSAL * DTR)
      COSBE = COS(COSBE * DTR)
      COSGA = COS(COSGA * DTR)
      SINAL=SQRT(1.0-COSAL**2)
      SINBE=SQRT(1.0-COSBE**2)
      SINGA=SQRT(1.0-COSGA**2)
      RCV=A*B*C*SQRT(1.0-COSAL**2-COSBE**2-COSGA**2+2.0*COSAL*COSBE
     1 *COSGA)
      RCP(1)=B*C*SINAL/RCV
      RCP(2)=A*C*SINBE/RCV
      RCP(3)=A*B*SINGA/RCV
      RCP(4)=(COSBE*COSGA-COSAL)/(SINBE*SINGA)
      RCP(5)=(COSAL*COSGA-COSBE)/(SINAL*SINGA)
      RCP(6)=(COSAL*COSBE-COSGA)/(SINAL*SINBE)
      VOL=RCP(1)*RCP(2)*RCP(3)*SQRT(1.0-RCP(4)**2-RCP(5)**2-RCP(6)**2+
     1 2.0*RCP(4)*RCP(5)*RCP(6))
C
      CC(4) = COSAL
      CC(5) = COSBE
      CC(6) = COSGA
C
      RETURN
      END
      SUBROUTINE CTHO
      CODE ANALYSIS
C.... CALCULATE ALL POSSIBLE Q VALUES WITHIN DATA RANGE FOR
C.... CUBIC? X(1)=X(2)=X(3), C2=0., KS=1
C.... TET?   X(1)=X(2), C2=0., KS=2
C.... HEX?   X(1)=X(2), C2=1., KS=3
C.... ORTHO? C2=0., KS=4
C.... QC(I) = X(1)*H*H + X(2)*K*(K+C2*H) + X(3)*L*L
C
C
      USE COMBLR
      INTEGER HMAX
C
      HMAX = SQRT( QO(NINPUT) / X(1) )  +  1.0
      KMAX = SQRT( QO(NINPUT) / X(2) )  +  1.0
      LMAX = SQRT( QO(NINPUT) / X(3) )  +  1.0
      QMIN = ( 2.0 * SIN( DTOR ) / WL1 )**2
      QMAX = 10000.
      C2 = 0.  $  IF ( KS .EQ. 3 ) C2 = 1.0
      I = 0
C
      DO 130 H=0,HMAX
        HHF = X(1) * FLOAT( H * H )
        KMIN = H
        IF ( KS .GT. 3 ) KMIN = 0
C
        DO 120 K=KMIN,KMAX
          HKF = HHF + X(2)*FLOAT( K )*( FLOAT( K ) + C2*FLOAT( H ) )
          IF ( HKF .GT. QMAX ) GO TO 130
          LMIN = 0
          IF ( KS .EQ. 1 ) LMIN = K
C
          DO 110 L=LMIN,LMAX
            QC(I+1) = HKF + X(3) * FLOAT( L * L )
            IF ( QC(I+1) .EQ. 0. ) GO TO 110
            IF ( QC(I+1) .GT. QMAX ) GO TO 120
            IF ( QC(I+1) .LT. QMIN ) GO TO 110
            IH = H
            IK = K
            IL = L
            KEEP = 0
C
C....       DETERMINE IF THIS REFLECTION ACCEPTABLE FOR ASPECT = NASP
C....                 C,   T,   H     O
            GO TO ( 101, 102, 103, 104 ), KS
C....         KS=1, CUBIC
  101         CALL ASPKS1( IH, IK, IL, NASP, KEEP )  $  GO TO 105
C....        KS=2, TETRAGONAL
  102         CALL ASPKS2( IH, IK, IL, NASP, KEEP )  $  GO TO 105
C....        KS=3, HEXAGONAL
  103         CALL ASPKS3( IH, IK, IL, NASP, KEEP )  $  GO TO 105
C....         KS=4, ORTHORHOMBIC
  104         CALL ASPKS4( IH, IK, IL, NASP, KEEP )  $  GO TO 105
  105       CONTINUE
            IF ( KEEP .EQ. 0 ) GO TO 110
C....       REFLECTION IS ACCEPTABLE
            I = I + 1
            IF ( QC(I) .GT. QO(NINPUT) ) QMAX = QC(I)
C....       STORE INDICES IN LAST 18 BITS
            QC(I) = ( QC(I) .INT. MASK ) .UN. ( 4096*IH + 64*IK + IL )
  110     CONTINUE
  120   CONTINUE
  130 CONTINUE
C.... NUMBER OF POSSIBLE Q-VALUES
      NPOS = I
      RETURN
      END
      SUBROUTINE FIG( J, K, L, DO, NP, NO )
      COMMON /SNYDER/ IOUT, RCP(6), RCV, SDQ, SDTT, SDT, VOL, WAVE
C
C  PROGRAM TO COMPUTE  FIGURES OF MERIT FOR POWDER PATTERNS OF ANY
C  CRYSTAL SYSTEM                                R. SNYDER  6/8/77
C
      CODE ANALYSIS
      DATA (DTR=.01745329)
C
C  COMPUTE CALCULATED D AND 2 THETA
C
      QC = J**2*RCP(1)**2 + K**2*RCP(2)**2 + L**2*RCP(3)**2
     1 + 2.0*J*K*RCP(1)*RCP(2)*RCP(6)
     2 + 2.0*J*L*RCP(1)*RCP(3)*RCP(5)
     3 + 2.0*K*L*RCP(2)*RCP(3)*RCP(4)
      DC = SQRT(1.0 / QC )
      TC = 2.0 * ASIN(WAVE / (2.0*DC)) / DTR
C
C
C  COMPUTE OBS 2 THETA
      TO = 2.0 * ASIN(WAVE / (2.0*DO)) / DTR
      QO = 1.0 / DO**2
C
C  SAVE SUMS
C
      DT = ABS(TO - TC)
      SDT = SDT + DT
      SDTT = SDTT + DT**2
      SDQ = SDQ + ABS(QO-QC)
      AT = SDT / NO
      RT = SQRT(SDTT / NO)
      IF(NP .GT. 0) GO TO 250
      WRITE(IOUT,40)J,K,L,NO,TO,TC,AT,RT
   40 FORMAT(1X,3I2,I3,3X,2(1XF7.3),2(1XF5.4))
      RETURN
  250 FNA = NO / (AT * NP)
      FNR = NO / (RT * NP)
      DQ = SDQ / NO
      DEWOLF = QO * NO
      DEWOLF = NO * QO / (2.0 * DQ * NP)
      RAT = FLOAT(NO) / FLOAT(NP)
      WRITE(IOUT,45)J,K,L,NO,NP,TO,TC,AT,RT,FNA,FNR,DEWOLF,RAT
   45 FORMAT(1X,3I2,2I3, 2(1XF7.3), 2(1XF5.4), 2(1XF6.2), 1XF6.0,1XF6.4)
      RETURN
      END
      SUBROUTINE GETASP( NASP, SG, XSG, LSGX, KS )
      CODE ANALYSIS
C.... ON ENTRY NASP IS TOTAL NUMBER OF ASPECTS IN THE CRYSTAL SYSTEM
C.... ON EXIT NASP IS THE ASPECT OF THE INPUT SPACE GROUP
C.... SG IS THE SPACE GROUP
C.... XSG IS THE TABLE OF SPACE GROUPS
C.... LSGX(I) = INDEX OF LAST SPACE GROUP IN ASPECT(I)
C.... KS IS CRYSTAL SYSTEM  ( C=1, T=2, H=3, O=4, M=5 )
      DIMENSION LSGX(2), XSG(2)
      INTEGER XSG, SG
C
C.... NSG IS TOTAL NUMBER OF SPACEGROUPS FOR THE SPECIFIC CRYSTAL SYSTEM
      NSG = LSGX(NASP)
C.... FIND INDEX OF SPACE GROUP
      DO 110 N = 1,NSG
        IF ( SG .EQ. XSG(N) ) GO TO 120
  110 CONTINUE
      WRITE (59,113) SG, KS
  113 FORMAT (29H CANNOT IDENTIFY SPACE GROUP , A10,
     , 19H FOR CRYSTAL SYSTEM, I2)
      NASP = 0
      RETURN
  120 CONTINUE
C.... GET ASPECT NUMBER OF SPACE GROUP INDEX
      DO 130 I=1,NASP
      IF ( N .LE. LSGX(I) ) GO TO 140
  130 CONTINUE
      PAUSE
  140 CONTINUE
      NASP = I
      IF ( KS .NE. 4 ) RETURN
C.... DETERMINE ORTHORHOMBIC ORIENTATION
      J = MOD( N, 6 )
      IF ( J .EQ. 0 ) J = 6
      NASP = 100*J + NASP
      RETURN
      END
      SUBROUTINE MONO
      CODE ANALYSIS
C.... CALCULATE ALL POSSIBLE Q WITHIN DATA RANGE FOR MONOCLINIC  ( KS = 5 <
C.... X(1) = 1. / ( A*SIN(BETA) )**2
C.... X(2) = 1.0 / ( B*B )
C.... X(3) = 1. / ( C*SIN(BETA) )**2
C.... X(4) = 2.0*COS(BETA) / (A*C*(SIN(BETA)**2)
C.... QC(I) = X(1)*H*H + X(2)*K*K + X(3)*L*L - X(4)*L*H
C
      USE COMBLR
      INTEGER HMIN, HMAX
      F1 = X(1) / ( X(1)*X(3) - .25*X(4)*X(4) )
      IF ( F1 .LT. 0. ) GO TO 1313
      F2 = 4.0 * X(1)
      F3 = 0.5 / X(1)
      KMAX = SQRT( QO(NINPUT) / X(2) )  +  1.0
      IF ( MOD( KMAX, 2 ) .NE. 0 ) KMAX = KMAX + 1
      QMAX = X(2) * FLOAT( KMAX * KMAX )
      QMIN = 0.
      I = 0
      DO 130 K=0,KMAX
        X2TERM = X(2) * FLOAT( K*K )
        IF ( X2TERM .GT. QMAX ) GO TO 140
        LMAX = SQRT( F1 * ( QMAX - X2TERM ) )
        DO 120 L=0,LMAX
          X4TERM = X(4) * FLOAT( L )
          SUMKL = X2TERM + X(3) * FLOAT( L*L )
          HPART = X4TERM*X4TERM - F2*( SUMKL - QMAX )
          IF ( HPART .LT. 0. ) GO TO 120
          HDEL = SQRT( HPART )
          HMIN = F3 * ( X4TERM - HDEL )
          HMAX = F3 * ( X4TERM + HDEL ) + 1.0
          DO 110 H=HMIN,HMAX
            IF ( L .GT. 0 ) GO TO 105
            JMAX = MAX0( H, K, L )
            IF ( JMAX .LE. 0 ) GO TO 110
            IF ( H .GE. 0 ) GO TO 105
            IF ( IABS( H ) .LE. HMAX ) GO TO 110
  105       CONTINUE
            QC(I+1) = SUMKL + FLOAT( H ) * ( X(1)*FLOAT( H ) - X4TERM )
            IF ( QC(I+1) .GT. QMAX ) GO TO 110
            IF ( QC(I+1) .LT. QMIN ) GO TO 110
            KEEP = 0
C....       DETERMINE IF THIS REFLECTION ACCEPTABLE FOR ASPECT = NASP
            CALL ASPKS5( H, K, L, NASP, KEEP )
            IF ( KEEP .EQ. 0 ) GO TO 110
            I = I + 1
C....       REFLECTION IS ACCEPTABLE
            IF ( QC(I) .GT. QO(NINPUT) ) QMAX = QC(I)
C....       PACK INDICES IN LAST 18 BITS
            QC(I) = (QC(I).INT.MASK) .UN. (4096*IABS( H ) + 64*K + L)
C....       RETAIN SIGN AND HOPE H _ 32
            IF ( H .LT. 0 ) QC(I) = QC(I) .UN. 40B
  110     CONTINUE
  120   CONTINUE
  130 CONTINUE
  140 CONTINUE
C... TOTAL NUMBER OF REFLECTIONS POSSIBLE FOR ASPECT = NASP
      NPOS = I
      RETURN
 1313 CONTINUE
      WRITE (59,1314)
 1314 FORMAT (15H F1 IS NEGATIVE)
      CALL EXIT( 1 )
      END
      SUBROUTINE QMATCH( IARG1, IARG2 )
C.... SUBROUTINE QMATCH MATCHES THE CALCULATED Q-VALUES TO THE INPUT Q-VALUES
      CODE ANALYSIS
      USE COMBLR
      DIMENSION DELTAQ((-1,1))
C...............................................................................
          IF ( IARG1 .EQ. 0 ) RETURN
C.... IARG1 CONTAINS NUMBER OF CALCULATED Q-VALUES
      NQC = IARG1
C.... IARG2 CONTAINS NUMBER OF INPUT Q-VALUES TO ATTEMPT TO MATCH
      NMAX = IARG2
      NQM = 0
      IL = 2
C
C.... ADD THE INPUT Q-VALUES TO THE CALCULATED Q-VALUE ARRAY
      DO 110 N=1,NMAX
        QC(NQC+N) = QO(N)
  110 CONTINUE
C
      IMAX = NQC + NMAX
C
C.... SORT IN ASCENDING ORDER
      CALL QSORT( QC, QC, IMAX )
      IMAX = IMAX + 1
C.... SUBROUTINE QSORT HAS PROVIDED BOUNDARY Q-VALUES AT Q(1) AND Q(IMAX&1<.
C.... PROVIDE AN ADDITIONAL BOUNDARY VALUE AT (IMAX+2)
      QC(IMAX+2) = QC(IMAX+1)
C
      DO 130 N=1,NMAX
C....   INITIALIZE QM(N)
        QM(N) = 0.
        DO 120 I=IL,IMAX
C....     FIND THE NEXT INPUT Q-VALUE
          IF ( QO(N) .NE. QC(I) ) GO TO 120
C
C....     INITIALIZE ARRAY DELTAQ TO UNACCEPTABLE SOLUTION
          DELTAQ(-1) = DELTAQ(0) = DELTAQ(1) = QMAX
C....     QC ARRAY IS IN ASCENDING ORDER.  START WHERE PRECEDING MATCH WAS MADE.
          IL = I
C....     DO NOT MATCH ONE INPUT Q TO ANOTHER INPUT Q
          IF ( QC(I-1) .NE. QO(N-1) ) DELTAQ(-1) = QO(N) - QC(I-1)
          IF ( QC(I+1) .NE. QO(N+1) ) DELTAQ(1) = QC(I+1) - QO(N)
C....     FIND MINIMUM OF THE THREE DIFFERENCES
          DQMIN = AMINAF( DELTAQ(1), -1, 1, 1, J)
C....     J = -1 IF CLOSEST MATCH FOR QO(N) = QC(I) IS QC(I-1)
C....     J = +1 IF CLOSEST MATCH FOR QO(N) = QC(I) IS QC(I+1)
          IF ( DQMIN .LE. QERR ) QM(N) = ABS( QC(I+J) * FLOAT( J ) )
C....     TALLY THE MATCHES
          IF ( QM(N) .EQ. 0. ) GO TO 130
          NQM = NQM + 1
          IF ( N .EQ. 1 ) GO TO 130
C
C....     ALLOW COMPUTED Q TO MATCH ONLY ONE OBSERVED Q
          IF ( QM(N-1) .NE. QM(N) ) GO TO 130
            NQM = NQM - 1
            DQR = QC(I-1) - QC(I-2)
            IF ( DQR .LT. DQMIN ) QM(N) = 0.
            IF ( DQR .GT. DQMIN ) QM(N-1) = 0.
            GO TO 130
C
  120   CONTINUE
C
  130 CONTINUE
C
C.... IARG1 RETURNS TOTAL NUMBER OF ACCEPTABLE MATCHES.
      IARG1 = NQM
      RETURN
C...............................................................................
      END
      SUBROUTINE QSORT(KEYS, REC, N)
      CODE ANALYSIS
C
C   QUICKSORT (EXCHANGE-PARTITION SORT)
C     ALGORITHM DEVELOPED BY C.A.R. HOARE  (COMPUT. J.15 (1962), PP. 10-15<
C
C     IN TEMPORARY STORAGE AN AUXILLIARY STACK OF _ LOG2(N) ENTRIES IS NECESSARY
C     TWO ARTIFICIAL KEYS ARE NECESSARY
C        AT KEYS(1) A KEY SMALLER THAN ANY OTHER
C     SUBFILES OF M OR FEWER ELEMENTS WILL BE SORTED BY STRAIGHT INSERTION
C        TYPICALLY A VALUE JUST BELOW 10 IS A GOOD CHOICE FOR M
C     THIS SORT IS UNSTABLE (A SPEED CONSIDERATION)
C    DECLARE VARIABLES
      INTEGER KEYS(N), REC(N), N
      INTEGER LSTACK(20), RSTACK(20), STACKP, I, J, L, KCAP, RCAP, R, M
C
      DATA (BOUND = 1.0E+10)
      EQUIVALENCE (BOUND,NOTIN)
C
C
C    SET M = 9
      DATA M/9/
C
C    SET UP THE DUMMY KEYS
      KEYS(N+1) = KEYS(1)
      REC(N+1) = REC(1)
      KEYS(1) = -NOTIN
      KEYS(N+2) = NOTIN
C
C  INITIALIZE
   10 CONTINUE
        STACKP = 0
        L = 2
        R = N + 1
C
C  BEGIN SORT OF THE SUBFILE REC(L)...REC(R)
C    R.GE.(L-1),  KEYS(L-1).LE.KEYS(I),  KEYS(I).LE.KEYS(R+1),  I IN RANGE %L,R<
   20 CONTINUE
C    USE STRAIGHT INSERTION  FOR FILE OF M OR LESS ENTRIES
          IF ((R-L).LT.M) GO TO 80
          I = L
          J = R
          KCAP = KEYS(L)
          RCAP = REC(L)
C
C  COMPARE KCAP AND KEYS(J)
   30 CONTINUE
        IF (KCAP.GE.KEYS(J)) GO TO 40
          J = J - 1
          GO TO 30
C
C  MOVE
C    KEYS(I).GE.KCAP,  KEYS(J).LE.KCAP.   KEYS(I) IS IRRELEVANT.
   40 CONTINUE
        IF (J.GT.I) GO TO 45
          REC(I) = RCAP
          GO TO 70
   45   CONTINUE
          REC(I) = REC(J)
          I = I + 1
C
C  COMPARE KEYS(I) AND KCAP
   50 CONTINUE
        IF (KEYS(I).GE.KCAP) GO TO 60
          I = I + 1
          GO TO 50
C
C  MOVE
C    KEYS(J).LE.KCAP,  KCAP.LE.KEYS(I).  KEYS(J) IS IRRELEVANT.
   60 CONTINUE
        IF (J.LE.I) GO TO 65
          REC(J) = REC(I)
          J = J - 1
          GO TO 30
   65   CONTINUE
          REC(J) = RCAP
          I = J
C
C  MODIFY STACK
C    THE SUBFILE REC(L)...REC(R) HAS BEEN PARTITIONED SO THAT KEYS(K).LE.KEYS%I<
C      FOR K IN THE RANGE (L,I).  KEYS(I).LE.KEYS(K) FOR K IN THE RANGE %I,R<.
   70 CONTINUE
        STACKP = STACKP + 1
          IF ((R-I).GE.(I-L)) GO TO 75
C    L IS FURTHER FROM I THAN R IS
            LSTACK(STACKP) = L
            RSTACK(STACKP) = I - 1
            L = I + 1
            GO TO 20
   75     CONTINUE
C    R IS FURTHER FROM I THAN L IS
            LSTACK(STACKP) = I + 1
            RSTACK(STACKP) = R
            R = I - 1
            GO TO 20
C
C  STRAIGHT INSERTION SORT OF SUBFILE R(L+1)...R(R)
   80 CONTINUE
        J = L + 1
   82     CONTINUE
            IF (J.GT.R) GO TO 90
            KCAP = KEYS(J)
            RCAP = REC(J)
            I = J - 1
   84         CONTINUE
                IF (KEYS(I).LE.KCAP) GO TO 86
                  REC(I+1) = REC(I)
                  I = I - 1
                  IF(I.GE.L) GO TO 84
   86         CONTINUE
                REC(I+1) = RCAP
                J = J + 1
                GO TO 82
C
C  TAKE OFF STACK
   90 CONTINUE
        IF (STACKP.EQ.0) RETURN
          L = LSTACK(STACKP)
          R = RSTACK(STACKP)
          STACKP = STACKP - 1
          GO TO 20
      END
      SUBROUTINE RDINPT
      CODE ANALYSIS
C.... READ INPUT, DETERMINE CRYSTAL SYSTEM, CALCULATE X(1), X(2), X(3), X%4<
      DIMENSION SYSTEM(6,5), CS1D(12), CS2D(2,6)
      EQUIVALENCE (CS1D,CS2D)
      DATA ( CS1D = 5HCUBIC, 1H , 10HTETRAGONAL, 1H , 9HHEXAGONAL, 1H ,
     , 10HORTHORHOMB, 2HIC, 10HMONOCLINIC, 1H , 8HANORTHIC, 1H  )
      USE COMBLR
  100 FORMAT (8A10)
      READ (IT,100) TITLE
C... EXIT ON BLANK TITLE CARD
      DO 101 I=1,8
        IF ( TITLE(I) .NE. 1H  ) GO TO 102
  101 CONTINUE
      CALL EXIT( 1 )
  102 CONTINUE
      WRITE (NT,104) TITLE
  104 FORMAT (1H1, //, 1X 8A10)
C.... USE NAMELIST TO READ INPUT
      INPUT DATA INPUTS, IT, NT
      CALL EMPTY( NT )
C.... DETERMINE INPUT TYPE - D, Q, SINSQ, TWOTH
      DO 105 J=1,4
        IF ( BLR(1,J) .NE. 0. ) GO TO 110
  105 CONTINUE
      WRITE (59,106)
  106 FORMAT (28H AN OBSERVED DATA POINT = 0.)
      CALL EXIT( 1 )
  110 CONTINUE
      IF ( J .LT. 4 ) GO TO 115
        FACTOR = 0.5 * DTOR
        IF ( BLR(1,J) .LT. 1. ) FACTOR = 0.5
  115 CONTINUE
      DO 120 N=1,NLIM
        IF ( BLR(N,J) .EQ. 0. ) GO TO 125
C....             D    Q SINSQ TWOTH
        GO TO ( 119, 120, 116, 117 ), J
C....   SINSQ INPUT
  116   CONTINUE
          SINTH = SQRT( SINSQ(N) )
          GO TO 118
C....   TWO THETA INPUT
  117   CONTINUE
          SINTH = SIN( FACTOR * TWOTH(N) )
  118     CONTINUE
          D(N) = .5 * WL1 / SINTH
C....   D INPUT
  119   CONTINUE
          QO(N) = 1.0 / ( D(N) * D(N) )
  120 CONTINUE
  125 CONTINUE
      NINPUT = N - 1
C.... DETERMINE CRYSTAL SYSTEM
C.... DEFAULT FOR B AND C IS A.  DEFAULT FOR ALPHA, BETA AND GAMMA IS 90 DEGREES
      DO 130 I=1,3
        IF ( CC(I) .EQ. 0. ) CC(I) = CC(1)
        IF ( CC(I+3) .EQ. 0. ) CC(I+3) = 90.
  130 CONTINUE
      DO 135 J=1,5
        DO 135 I=1,3
          SYSTEM(I,J) = CC(I)
          SYSTEM(I+3,J) = 90.
  135 CONTINUE
C.... CUBIC B       CUBIC C       TET B         HEX B
      SYSTEM(2,1) = SYSTEM(3,1) = SYSTEM(2,2) = SYSTEM(2,3) = CC(1)
C.... MONOCLINIC BETA
      SYSTEM(5,5) = CC(5)
C... HEXAGONAL GAMMA
      SYSTEM(6,3) = 120.
      DO 145 J=1,5
        DO 140 I=1,6
          IF ( CC(I) .NE. SYSTEM(I,J) ) GO TO 145
  140   CONTINUE
          GO TO 150
  145 CONTINUE
      WRITE (59,146)
  146 FORMAT (25H ANORTHIC NOT IMPLEMENTED)
      CALL EXIT( 1 )
  150 CONTINUE
  151 FORMAT (//, 19H CRYSTAL SYSTEM IS , 2A10)
      WRITE (NT,151) CS2D(1,J), CS2D(2,J)
      CALL EMPTY( NT )
C.... KS=1 FOR CUBIC, KS=2 FOR TET, KS=3 FOR HEX, KS=4 FOR ORTHO, KS=5 FOR MONO
      KS = J
C.... CALCULATE X(I), I=1,4
C.... X(1) = C1 / ( A*SIN(BETA) )**2
C.... X(2) = 1.0 / ( B*B )
C.... X(3) = C1 / ( C*SIN(BETA) )**2
C.... X(4) = 2.0*COS(BETA) / (A*C*(SIN(BETA)**2)
C.... FOR CUBIC, TETRAGONAL AND ORTHORHOMBIC C1 = 1..  FOR HEXAGONAL C1 # 4./3.
      C1 = 1.0
      IF ( KS .EQ. 3 ) C1 = 1.333333333333
      X(1) = C1 / ( A * A )
      X(2) = C1 / ( B * B )
      X(3) = 1.0 / ( C * C )
      X(4) = 1.0 / SIN( BETA * DTOR )
      X(1) = X(1) * X(4) * X(4)
      X(3) = X(3) * X(4) * X(4)
      X(4) = 2.0 * COS( BETA * DTOR ) * X(4) * X(4) / ( A * C )
      RETURN
      END
      SUBROUTINE SGID( KS, SG, NASP )
      CODE ANALYSIS
C.... THIS SUBROUTINE SUPPLIES THE TABLES WHICH DETERMINE THE ASPECT NUMBER
C.... SEE ~CRYSTAL DATA~  -  DONNAY AND ONDIK
C
C.... BLOCKS XSG (X=C,T,H,O,M) ARE DATA LOADED WITH THE HOLLERITH REPRESENTATION
C....        OF THE SPACE GROUP.  SUBSCRIPTS ARE FOLLOWED BY A BLANK.
C.... BLOCKS LSGX (X=C,T,H,O,M) ARE DATA LOADED WITH THE INDEX OF THE LAST SPACE
C....        GROUP OF THE ASPECT.
C
      INTEGER OSG1D(384), NAS(6)
      INTEGER CSG(36), TSG(84), HSG(56), OSG2D(6,64), MSG(29), ASG(2),
     , NSG(6), LSGC(17), LSGT(31), LSGH(52), LSGO(33), LSGM(14), LSGA
      EQUIVALENCE (OSG2D,OSG1D)
C
C
C.... NUMBER OF ASPECTS ACCORDING TO CRYSTAL SYSTEM
C....               C,  T, H,  O,  M, A
      DATA ( NAS = 17, 31, 9, 33, 14, 1 )
C.... NUMBER OF SPACE GROUPS ACCORDING TO CRYSTAL SYSTEM
C....               C,  T,  H,   O,  M, A
      DATA ( NSG = 36, 84, 56, 384, 29, 2 )
C.... CUBIC
      DATA ( LSGC = 5, 7, 9, 11, 13, 14, 15, 21, 22, 23, 24, 25, 30, 31,
     , 33, 35, 36 )
      DATA ( CSG =
     , 4HPM3M, 4HP432, 5HP-43M, 3HPM3, 3HP23, 6HP42 32, 5HP21 3,
     , 6HP41 32, 6HP43 32, 4HPN3M, 3HPN3, 4HPM3N, 5HP-43N, 4HPN3N,
     , 3HPA3, 4HIM3M, 4HI432, 5HI-43M, 3HIM3, 3HI23, 5HI21 3, 6HI41 32,
     , 3HIA3, 5HI-43D, 4HIA3D,
     , 4HFM3M, 4HF432, 5HF-43M, 3HFM3, 3HF23, 6HF41 32, 4HFD3M, 3HFD3,
     , 4HFM3C, 5HF-43C, 4HFD3C )
C.... TETRAGONAL
      DATA ( LSGT =
     ,       8, 11, 15, 17, 18, 20, 25, 26, 29, 33, 38, 40, 45, 47, 49,
     ,  50, 52, 53, 55, 57, 58, 60, 61, 69, 71, 74, 77, 79, 80, 82, 84 )
      DATA ( (TSG(I), I=1,57) =
     , 6HP4/MMM, 4HP422, 4HP4MM, 5HP-42M, 5HP-4M2, 4HP4/M, 2HP4, 3HP-4,
     , 6HP42 22, 6HP42 /M, 4HP42 , 6HP41 22, 6HP43 22, 4HP41 , 4HP43 ,
     , 6HP421 2, 7HP-421 M, 8HP42 21 2, 8HP41 21 2, 8HP43 21 2,
     , 6HP4/MMC, 8HP42 /MMC, 4HP4MC, 6HP42 MC, 5HP-42C, 6HP421 C,
     , 6HP4/MBM, 4HP4BM, 5HP-4B2, 6HP4/MBC, 8HP42 /MBC, 4HP4BC,
     , 6HP42 BC, 6HP4/MCM, 8HP42 /MCM, 4HP4CM, 6HP42 CM, 5HP-4C2,
     , 6HP4/MCC, 4HP4CC, 6HP4/MNM, 8HP42 /MNM, 4HP4NM, 6HP42 NM,
     , 5HP-4N2, 6HP4/MNC, 4HP4NC, 6HP4/NMM, 4HP4/N, 6HP42 /N, 6HP4/NMC,
     , 8HP42 /NMC, 6HP4/NBM, 6HP4/NBC, 8HP42 /NBC, 6HP4/NCM, 8HP42 /NCM)
      DATA ( (TSG(I), I=58,84) =
     , 6HP4/NCC, 6HP4/NNM, 8HP42 /NNM, 6HP4/NNC,
     , 6HI4/MMM, 4HI422, 4HI4MM, 5HI-42M, 5HI-4M2, 4HI4/M, 2HI4, 3HI-4,
     , 6HI41 22, 4HI41 , 6HI4/MCM, 4HI4CM, 5HI-4C2, 4HI4MD, 6HI41 MD,
     , 5HI-42D, 4HI4CD, 6HI41 CD, 6HI41 /A, 6HI4/AMD, 8HI41 /AMD,
     , 6HI4/ACD, 8HI41 /ACD )
C.... HEXAGONAL
      DATA ( LSGH = 16, 19, 29, 33, 40, 47, 49, 54, 56 )
      DATA ( HSG =
     , 6HP6/MMM, 4HP622, 4HP6MM, 5HP-6M2, 5HP-62M, 4HP6/M, 2HP6, 3HP-6,
     , 5HP-3M1, 5HP-31M, 4HP321, 4HP312, 4HP3M1, 4HP31M, 3HP-3, 2HP3,
     , 6HP63 22, 6HP63 /M, 4HP63 , 6HP62 22, 6HP64 22, 4HP62 , 4HP64 ,
     , 6HP31 21, 6HP32 21, 6HP31 12, 6HP32 12, 4HP31 , 4HP32 ,
     , 6HP61 22, 6HP65 22, 4HP61 , 4HP65 ,
     , 6HP6/MCM, 8HP63 /MCM, 4HP6C3, 6HP63 CM, 5HP-6C2, 5HP-3C1, 4HP3C1,
     , 6HP6/MMC, 8HP63 /MMC, 4HP6MC, 6HP63 MC, 5HP-62C, 5HP-31C, 4HP31C,
     , 6HP6/MCC, 4HP6CC,
     , 4HR-3M, 3HR32, 3HR3M, 3HR-3, 2HR3, 4HR-3C, 3HR3C )
C.... ORTHORHOMBIC
      DATA ( LSGO =                   18,  24,  30,  36,  60,  78,  90,
     , 108, 120, 132, 150, 162, 168, 174, 180, 186, 192, 198, 204, 228,
     , 234, 258, 270, 282, 294, 300, 324, 336, 348, 354, 372, 378, 384 )
      DATA ( (OSG1D(I), I=1,66) =
     , 6(4HPMMM), 6(4HP222), 2(4HPMM2), 4HPM2M, 2(4HP2MM), 4HPM2M,
     , 2(5HP2221), 6HP221 2, 2(6HP21 22), 6HP221 2,
     , 2(8HP21 21 2), 8HP21 221 , 2(8HP221 21 ), 8HP21 221 ,
     , 6(10HP21 21 21 ),
     , 4HPMMA, 4HPMMB, 4HPMCM, 4HPCMM, 4HPBMM, 4HPMAM,
     , 4HPM2A, 4HP2MB, 4HP2CM, 4HPC2M, 4HPBM2, 4HPMA2,
     , 4HP2MA, 4HPM2B, 4HPMC2, 4HPCM2, 4HPB2M, 4HP2AM,
     , 6HP21 MA, 6HPM21 B, 6HPMC21 , 6HPCM21 , 6HPB21 M, 6HP21 AM,
     , 4HPMNM, 2(4HPNMM)     , 4HPMNM, 2(4HPMMN) )
      DATA ( (OSG1D(I), I=67,120) =
     , 4HPMN2, 4HPNM2, 4HPN2M, 4HP2NM, 4HP2MN, 4HPM2N,
     , 6HPMN21 , 6HPNM21 , 6HPN21 M, 6HP21 NM, 6HP21 MN, 6HPM21 N,
     , 2(4HPCCM), 4HPBMB, 2(4HPMAA), 4HPBMB,
     , 2(4HPCC2), 4HPB2B, 2(4HP2AA), 4HPB2B,
     , 4HPCAM, 4HPBCM, 4HPCMB, 4HPMCA, 4HPMAB, 4HPBMA,
     , 4HPCA2, 4HPBC2, 4HPC2B, 4HP2CA, 4HP2AB, 4HPB2A,
     , 6HPCA21 , 6HPBC21 , 6HPC21 B, 6HP21 CA, 6HP21 AB, 6HPB21 A,
     , 2(4HPBAM), 4HPCMA, 2(4HPMCB), 4HPCMA,
     , 2(4HPBA2), 4HPC2A, 2(4HP2CB), 4HPC2A )
      DATA ( (OSG1D(I), I=121,174) =
     , 4HPNCM, 4HPCNM, 4HPBMN, 4HPMAN, 4HPMNA, 4HPNMB,
     , 4HPNC2, 4HPCN2, 4HPB2N, 4HP2AN, 4HP2NA, 4HPN2B,
     , 4HPNAM, 4HPBNM, 4HPCMN, 4HPMCN, 4HPMNB, 4HPNMA,
     , 6HPNA21 , 6HPBN21 , 6HPC21 N, 6HP21 CN, 6HP21 NB, 6HPN21 A,
     , 4HPNA2, 4HPBN2, 4HPC2N, 4HP2CN, 4HP2NB, 4HPN2A,
     , 2(4HPNNM), 4HPNMN, 2(4HPMNN), 4HPNMN,
     , 2(4HPNN2), 4HPN2N, 2(4HP2NN), 4HPN2N,
     , 4HPCCA, 4HPCCB, 4HPBCB, 4HPCAA, 4HPBAA, 4HPBAB,
     , 4HPBCA, 4HPCAB, 4HPBCA, 4HPCAB, 4HPBCA, 4HPCAB )
      DATA ( (OSG1D(I), I=175,234) =
     , 2(4HPCCN), 4HPBNB, 2(4HPNAA), 4HPBNB,
     , 2(4HPBAN), 4HPCNA, 2(4HPNCB), 4HPCNA,
     , 4HPBCN, 4HPCAN, 4HPBNA, 4HPNAB, 4HPNCA, 4HPCNB,
     , 4HPNNA, 4HPNNB, 4HPNCN, 4HPCNN, 4HPBNN, 4HPNAN,        6(4HPNNN),
     , 2(4HCMMM), 4HBMMM, 2(4HAMMM), 4HBMMM,
     , 2(4HC222), 4HB222, 2(4HA222), 4HB222,
     , 2(4HCMM2), 4HBM2M, 2(4HA2MM), 4HBM2M,
     , 4HCM2M, 4HC2MM, 4HB2MM, 4HAM2M, 4HAMM2, 4HBMM2,
     , 2(6HC2221 ), 6HB221 2, 2(6HA21 22), 6HB221 2 )
      DATA ( (OSG1D(I), I=235,288) =
     , 4HCMCM, 4HCCMM, 4HBBMM, 4HAMAM, 4HAMMA, 4HBMMB,
     , 4HCMC2, 4HCCM2, 4HBB2M, 4HA2AM, 4HA2MA, 4HBM2B,
     , 6HCMC21 , 6HCCM21 , 6HBB21 M, 6HA21 AM, 6HA21 MA, 6HBM21 B,
     , 4HC2CM, 4HCC2M, 4HBBM2, 4HAMA2, 4HAM2A, 4HB2MB,
     , 4HABMM, 4HBMAM, 2(4HCMMA), 4HBMAM, 4HABMM,
     , 4HABM2, 4HBMA2, 4HCM2A, 4HC2MA, 4HB2AM, 4HAB2M,
     , 4HABAM, 4HBBAM, 4HCCMA, 4HCMCA, 4HBMAB, 4HABMA,
     , 4HABA2, 4HBBA2, 4HCC2A, 4HC2CA, 4HB2AB, 4HAB2A,
     , 2(4HCCCM), 4HBBMB, 2(4HAMAA), 4HBBMB )
      DATA ( (OSG1D(I), I=289,372) =
     , 2(4HCCC2), 4HBB2B, 2(4HA2AA), 4HBB2B,
     , 2(4HCCCA), 4HBBAB, 2(4HABAA), 4HBBAB,
     , 6(4HIMMM), 6(4HI222), 6(10HI21 21 21 ),
     , 2(4HIMM2), 4HIM2M, 2(4HI2MM), 4HIM2M,
     , 4HIMAM, 2(4HIBMM), 4HIMAM, 2(4HIMMA),
     , 4HIMA2, 4HIBM2, 4HIB2M, 4HI2AM, 4HI2MA, 4HIM2A,
     , 2(4HIBAM), 4HIBMA, 2(4HIMAA), 4HIBMA,
     , 2(4HIBA2), 4HIB2A, 2(4HI2AA), 4HIB2A,                  6(4HIBCA),
     , 6(4HFMMM), 6(4HF222), 2(4HFMM2), 4HFM2M, 2(4HF2MM), 4HFM2M )
      DATA ( (OSG1D(I), I=373,384) =
     , 2(4HFDD2), 4HFD2D, 2(4HF2DD), 4HFD2D,               6(4HFDDD)   )
C.... MONOCLINIC
C.... SOME MONO CLINIC ASPECT NUMBERS NOT SAME AS ~CRYSTAL DATA~ ENTRIES
      DATA ( LSGM = 3, 5, 7, 9, 11, 12, 13, 14, 17, 20, 23, 25, 27, 29 )
      DATA (MSG = 4HP2/M, 2HP2, 2HPM, 6HP21 /M, 4HP21 , 4HP2/C, 2HPC,
     , 4HP2/A, 2HPA, 4HP2/N, 2HPN, 6HP21 /C, 6HP21 /A, 6HP21 /N, 4HC2/M,
     , 2HC2, 2HCM, 4HA2/M, 2HA2, 2HAM, 4HI2/M, 2HI2, 2HIM, 4HC2/C, 2HCC,
     , 4HA2/A, 2HAA, 4HI2/A, 2HIA )
C.... ANORTHIC
      DATA (LSGA=2), (ASG = 3HP-1, 2HP1)
C
C
C
C.... NUMBER OF ASPECTS IN CRYSTAL SYSTEM
      NASP = NAS(KS)
C.... DETERMINE ASPECT NUMBER
      GO TO ( 100, 200, 300, 400, 500, 600 ), KS
C.... CUBIC
  100 CALL GETASP( NASP, SG, CSG, LSGC, KS )  $  GO TO 700
C.... TETRAGONAL
  200 CALL GETASP( NASP, SG, TSG, LSGT, KS )  $  GO TO 700
C.... HEXAGONAL
  300 CALL GETASP( NASP, SG, HSG, LSGH, KS )  $  GO TO 700
C.... ORTHORHOMBIC
  400 CALL GETASP( NASP, SG, OSG1D, LSGO, KS )  $  GO TO 700
C.... MONOCLINIC
  500 CALL GETASP( NASP, SG, MSG, LSGM, KS )  $  GO TO 700
C.... ANORTHIC
  600 CALL GETASP( NASP, SG, ASG, LSGA, KS )
      WRITE (59,601)
  601 FORMAT (25H ANORTHIC NOT IMPLEMENTED)
      CALL EXIT( 1 )
  700 CONTINUE
      WRITE (59,701) KS, NASP, SG
  701 FORMAT (5H KS =, I2, 11H,  ASPECT =, I4, 8H,  SG = , A10)
      RETURN
      END
