C===========================================================================
	SUBROUTINE DGEN(STATUS)
	INTEGER STATUS
C
C	D-space generation program based on code from the 
C     NATIONAL BUREAU OF STANDARDS D-SPACING PROGRAM (MODIFIED 7/7/79)
C     AS MODIFIED BY D. K. SMITH AND G. G. J. JOHNSON
C
C      B.H. Toby 1/92
C===========================================================================
	include 'cmpr.common'
	include '[-.shared]flabel.common'

	INTEGER NRMAX
	PARAMETER (NRMAX=1000)
	INTEGER*2 H1(NRMAX),K1(NRMAX),L1(NRMAX),ISEQ(NRMAX),I2
	REAL	D(NRMAX)
	real RCP(6),PRCP(6)
	real  A1, B1, C1, AL, BE, GA, RCV, VOL

	REAL DMIN,DMAX
	INTEGER LAUE, IEXT(13)/13*0/

	CHARACTER*1 NQUES,ANS
	INTEGER nsys, ncr, NREF, LINNUM, I, J, nfil, NCOL, NLINES, NR
	INTEGER j1, j2
	REAL XLIMIT, x1, VOLREAL, TTMIN, TTMAX
	REAL ST, SST, CT, TTH, XLP
	INTEGER RUNNUM/0/
	SAVE RUNNUM

	IF (MFILE .eq. M$FIL) THEN
	   WRITE (lo,'(A)') 'No more room for files'
	   goto 9911
	ENDIF

	nfil = mfile + 1

	CALL CELLPARM(lo,li,LAUE,A1,B1,C1,AL,BE,GA)

	if (laue .le. 0 .or. MIN(A1,B1,C1) .LE. 0.0) then
	  CALL GETLAUE(LAUE,NSYS,lo,li)
	  if (laue .le. 0 .or. nsys .le. 0) goto 9911

	  CALL GETCELL(NSYS,lo,li,A1,B1,C1,AL,BE,GA)
	ENDIF
	IF (MIN(A1,B1,C1) .LE. 0.0) goto 9911

	CALL CELCON(A1,B1,C1,AL,BE,GA,RCV,RCP,VOL,PRCP)
	IF (VOL .LE. 0.0) goto 9911

	MONO(nfil) = .false.
	if (WAVE0 .gt. 0) then 
	  WAVE(nfil) = WAVE0
	  WKA1(nfil) = WAVE1
	  WKA2(nfil) = WAVE2
	  MONO(nfil) = .true.
	ELSE
133	  CALL WAVEINQ(LO,LI,IRAD(nfil),WAVE(nfil),WKA1(nfil),WKA2(nfil)
	1	,*133)
	  MONO(nfil) = .true.
	  WAVE(nfil) = WKA1(nfil)
	  WKA2(nfil) = WAVE(nfil)
	ENDIF
	IF (WAVE(nfil) .le. 0) goto 133

74	WRITE (lo,76)
76	FORMAT (/' ENTER THE 2THETA LIMIT FOR D VALUES [90.00]: ',$)
	READ (li,'(Q,f10.0)',ERR=74,END=9911) NCR,XLIMIT
	IF (XLIMIT .LE. 0.0 .or. NCR .LE. 0) XLIMIT=90.00
	IF (XLIMIT .GT. 179.00) GO TO 74

	CALL TTTOD(DMIN,XLIMIT,wave(nfil))
	CALL TTTOD(DMAX,0.5,wave(nfil))

80	WRITE (lo,82)
82	FORMAT (/' ENTER REFLECTION CONDITIONS FOR LATTICES ',
	1	'AND SYMMETRY ELEMENTS.')
	CALL GETEXTN(IEXT,lo,li)

100	CALL DGENSUB(lo,li,LAUE,DMIN,DMAX,NREF,IEXT,
	1	A1,B1,C1,rcp,prcp,
	1	NRMAX,D,H1,K1,L1)

	IF (LAUE .EQ. 7 .or. LAUE .eq. 8) THEN
	  WRITE (lo,1001) NREF,'CUBIC'
1001	  FORMAT (/,1x,I5,' D-spacings calculated for ',a,
	1	' system with:')
	ELSEIF (LAUE .EQ. 5 .or. LAUE .eq. 6) THEN
	  WRITE (lo,1001) NREF,'TETRAGONAL'
	ELSEIF (LAUE .EQ. 4) THEN
	  WRITE (lo,1001) NREF, 'ORTHORHOMBIC'
	ELSEIF (LAUE .ge. 11 .and. LAUE .le. 15) THEN
	  WRITE (lo,1001) NREF, 'HEXAGONAL'
	ELSEIF (LAUE .EQ. 2 .or. LAUE .eq. 3) THEN
	  WRITE (lo,1001) NREF, 'MONOCLINIC'
	ELSEIF (LAUE .EQ. 1) THEN
	  WRITE (lo,1001) NREF, 'TRICLINIC'
	ELSEIF (LAUE .EQ. 9 .or. LAUE .eq. 10) THEN
	  WRITE (lo,1001) NREF, 'RHOMBOHEDRAL'
	ELSE
	  WRITE (lo,'(A)') ' Illegal value for LAUE'
	ENDIF

	VOLREAL=1.0/VOL
	WRITE (lo,1007) A1,B1,C1,AL,BE,GA,WAVE(nfil),DMIN,VOLREAL
1007	FORMAT (/T10,'A=        ',F8.4,'   B=     ',F11.4,'   C='
     *,1X,'       ',F8.4,
     */T10,'ALPHA= ',F10.2,4X,'BETA=  ',F10.2,4X,'GAMMA= ',F10.2,/
     *T10,'LAMBDA= ',F11.5,2X,'D-MIN= ',F11.4,2X,' VOLUME=  'F11.3)

	IF (NREF .LE. 0) THEN
	  WRITE (lo,76)
	  READ (li,'(Q,f10.0)',ERR=74,END=9911) NCR,XLIMIT
	  IF (XLIMIT .LE. 0.0 .or. NCR .LE. 0) XLIMIT=90.00
	  IF (XLIMIT .GT. 179.00) GO TO 74
	  goto 100
	ENDIF

C     ORDER ON D
	DO I=1,NREF
	  ISEQ(I) = I
	ENDDO
	i2 = NREF
	CALL SORT (ISEQ, D, I2, -1)

	WRITE(lo,1145)
1145	FORMAT('$Type the generated d-spaces? (Y/[N]): ')
	READ(li,1146,END=9911) NCR,ANS
1146	FORMAT(Q,1A1)
	IF(ANS.EQ.'Y'.OR.ANS.EQ.'y') THEN
        LINNUM=60
	DO J=1,NREF
	  I = ISEQ(J)
	  CALL DTOTT(D(J),TTH,wave(nfil))
	  ST = SIN(TTH * 3.14159 / 180.)
	  SST=ST**2
	  CT = SQRT(1.0 -SST)
	  XLP=SST*CT/(1.0 +(1.0 -2.0 *SST)**2)
	  IF (LINNUM .GE. 22) THEN
	    IF (LINNUM .EQ. 22) CALL PS_PAUSE(PROMPT)
	    WRITE (lo,1415)
1415	    FORMAT ('1      H      K      L       D-SPACING    SINE SQ.'
	1	,'     2THETA      LP')
	    LINNUM=0
	  ENDIF

	  WRITE (lo,120) H1(I),K1(I),L1(I),D(J),SST,TTH,XLP
120	  FORMAT (1X,3I7,F16.6,F12.5,F11.2,F12.4)
	  LINNUM=LINNUM+1
	ENDDO
	ENDIF

	WRITE(lo,1144)
1144	FORMAT('$Print a copy of the generated d-spaces? (Y/[N]): ')
	READ(li,1146,END=9911) NCR,ANS

	IF(ANS.EQ.'Y'.OR.ANS.EQ.'y') THEN
	  call opnpr(lo,li,out,'DGEN.LIS',ncol,nlines)
	  IF (LAUE .EQ. 7 .or. LAUE .eq. 8) THEN
	    WRITE (out,1001) NREF,'CUBIC'
	  ELSEIF (LAUE .EQ. 5 .or. LAUE .eq. 6) THEN
	    WRITE (out,1001) NREF,'TETRAGONAL'
	  ELSEIF (LAUE .EQ. 4) THEN
	    WRITE (out,1001) NREF, 'ORTHORHOMBIC'
	  ELSEIF (LAUE .ge. 11 .and. LAUE .le. 15) THEN
	    WRITE (out,1001) NREF, 'HEXAGONAL'
	  ELSEIF (LAUE .EQ. 2 .or. LAUE .eq. 3) THEN
	    WRITE (out,1001) NREF, 'MONOCLINIC'
	  ELSEIF (LAUE .EQ. 1) THEN
	    WRITE (out,1001) NREF, 'TRICLINIC'
	  ELSEIF (LAUE .EQ. 9 .or. LAUE .eq. 10) THEN
	    WRITE (out,1001) NREF, 'RHOMBOHEDRAL'
	  ENDIF
	  VOLREAL=1.0/VOL
	  WRITE (out,1007) A1,B1,C1,AL,BE,GA,WAVE(nfil),DMIN,VOLREAL
152	  FORMAT ('0Extinction codes are ',13I2,'  Laue Code is',i3/)
	  WRITE (out,152) (IEXT(I),I=1,13),laue
	  LINNUM=60
	  DO J=1,NREF
	    I = ISEQ(J)
	    CALL DTOTT(D(J),TTH,wave(nfil))
	    ST = SIN(TTH * 3.14159 / 180.)
	    SST=ST**2
	    CT= SQRT(1.0 -SST)
	    XLP=SST*CT/(1.0 +(1.0 - 2.0*SST)**2)
	    IF (LINNUM .GE. 55) THEN 
	      WRITE (out,1415)
	      LINNUM=0
	    ENDIF

	    WRITE (out,120) H1(I),K1(I),L1(I),D(J),SST,TTH,XLP
	    LINNUM=LINNUM+1
	  ENDDO
	  CALL CLPR(lo,li,out)
	ENDIF

	CALL DTOTT(D(1),TTMIN,wave(nfil))
	CALL DTOTT(D(NREF),TTMAX,wave(nfil))
	J1 = NREF
200	DO I=NREF,1,-1
	  CALL DTOTT(D(I),TTH,wave(nfil))
	  IF (TTH .GE. TTMIN)  J1 = I
	ENDDO
	J2 = 1
	DO I=1,NREF
	  CALL DTOTT(D(I),TTH,wave(nfil))
	  IF (TTH .LE. TTMAX)  J2 = I
	ENDDO
	NR = MAX(0,J2-j1+1)
	IF (NR .GT. M$PKS) THEN
	  WRITE (lo,'(A,i4,A)') ' There are',NR-M$PKS,' too many peaks.'
	  WRITE (lo,'(A,i4,A)') ' For the program dimensions (',M$PKS,
	1	'). Change the 2-theta limits.'
	  ans = 'N'
	ELSEIF (NR .Le. 0) THEN
	  WRITE (lo,'(2A)') ' There are no peaks in the 2theta ',
	1	'range. Change the 2-theta limits.'
	  ans = 'N'
	ELSE
	  WRITE (lo,'(A,i4,A,f6.2,a,f6.2)') ' There are',NR,
	1	' peaks between',TTMin,' and',TTMax
	  WRITE (lo,'(A)') '$Is this OK? ([Y]/N): '
	  READ (li,'(a)',end=9911) ans
	ENDIF

	IF (ans .eq. 'N' .or. ans .eq. 'n') THEN
	  WRITE (lo,'(A,f6.2,a)') '$Enter 2Theta minimum [default=',
	1	ttmin,']: '
	  READ (li,'(q,F10.0)',err=200,end=9911) I,x1
	  IF (i .gt. 0) ttmin = x1
	  WRITE (lo,'(A,f6.2,a)') '$Enter 2Theta maximum [default=',
	1	ttmax,']: '
	  READ (li,'(q,F10.0)',err=200,end=9911) I,x1
	  IF (i .gt. 0) ttmax = x1
	  goto 200
	ENDIF

	NPKS(nfil) = J2 - J1 + 1
	I = 0
	DO J = J1,J2 
	  I = I + 1
	  FINTEG(I,nfil) = 1.
	  PKHT(I,nfil) = 1.
	  DSPACE(I,nfil) = D(J)
	  RELINT(I,nfil) = 100. 
	  call dTOtt(DSPACE(I,nfil),TT(I,nfil),wave(nfil))
	  IH(I,nfil)= H1(ISEQ(J))
	  IK(I,nfil)= K1(ISEQ(J))
	  IL(I,nfil)= L1(ISEQ(J))
	ENDDO
C scale maximum relative intensity to 100.0
	START(nfil) = TTmin
	FINISH(nfil) = TTmax
	DELTA(nfil) = 0.0

	RUNNUM = RUNNUM + 1

C save the name
	WRITE (ATITLE(nfil),2007) A1,B1,C1,AL,BE,GA
2007	FORMAT ('A=',F6.2,' B=',F6.2,' C=',F6.2,
     	1 ' \ga=',f5.1,' \gb=',F5.1,' \gc=',F5.1)

	write(ofilnam(nfil),'(A,i3.3))') 'DGEN',RUNNUM
	cfilnam(nfil) = ofilnam(nfil)
	FYTYP(nfil) = YTYP_norm
	FXTYP(nfil) = XTYP_D
	IF (wave(nfil) .gt. 0) FXTYP(nfil) = XTYP_TT

	write (LO,'(I5,A)') NPKS(nfil),' Peaks Created by D-GEN'
	mfile = mfile + 1
	STATUS = .TRUE.
	RETURN

9911	STATUS = .FALSE.
	RETURN
	END
C===========================================================================
	SUBROUTINE GETCELL(NSYS,no,ni,A,B,C,AL,BE,GA)
	INTEGER NSYS,NO,NI
	REAL A,B,C,AL,BE,GA		
C===========================================================================
	REAL AMIN
	INTEGER NCR
   25 WRITE (NO,26)
   26 FORMAT (/' ENTER THE CELL CONSTANTS:'/T25,'A    :  ',$)
      READ (NI,28,ERR=25,END=9911) A
      IF (A.LT.2.0 .OR. A.GT.99.0) GO TO 25
   28 FORMAT (F10.0)
      IF (NSYS.EQ.6) THEN
      B=A
      C=A
      AL=90.0
      BE=90.0
      GA=90.0
      GO TO 69
      END IF
      IF (NSYS.EQ.7) THEN
      B=A
      C=A
      GO TO 42
      END IF
      IF (NSYS.EQ.5 .OR. NSYS.EQ.8) THEN
      B=A
      GO TO 38
      END IF
   30 WRITE (NO,32)
   32 FORMAT (T25,'B[=A]:  ',$)
   34 READ (NI,36,ERR=30,END=9911) NCR,B
   36 FORMAT (Q,F10.0)
      IF (NCR.EQ.0) B=A
      IF (B.LT.2.0 .OR. B.GT.99.0) GO TO 30
   38 WRITE (NO,40)
   40 FORMAT (T25,'C[=A]:  ',$)
      READ (NI,36,ERR=38,END=9911) NCR,C
      IF (NCR.EQ.0) C=A
      IF (C.LT.2.0 .OR. C.GT.99.0) GO TO 38
      IF (NSYS.EQ.4 .OR. NSYS.EQ.5 .OR. NSYS.EQ.8) THEN
      AL=90.0
      BE=90.0
      GA=90.0
      IF (NSYS.EQ.8) GA=120.0
      GO TO 69
      END IF
      IF (NSYS.EQ.2) THEN
      AL=90.0
      GA=90.0
      GO TO 52
      END IF
      IF (NSYS.EQ.3) THEN
      AL=90.0
      BE=90.0
      GO TO 62
      END IF
   42 WRITE (NO,44)
   44 FORMAT (T10,'ALPHA(DEGREES)[90.00]:  ',$)
      READ (NI,36,ERR=42,END=9911) NCR,AL
      IF (NCR.EQ.0) AL=90.00
      IF (AL.LT.30.0 .OR. AL.GT.150.0) GO TO 42
      IF (NCR.NE.0) THEN
      IF (MOD(AL,1.0).NE.0.0) GO TO 51
   46 WRITE (NO,48)
   48 FORMAT (T11,'ALPHA(MINUTES)[0.00]:  ',$)
      READ (NI,36,ERR=46,END=9911) NCR,AMIN
      IF (NCR.EQ.0) AMIN=0.0
      IF (AMIN.LT.0.0 .OR. AMIN.GT.59.9) GO TO 46
      AL=AL+AMIN/60.
      ELSE
      END IF
   51 IF (NSYS.EQ.7) THEN
      BE=AL
      GA=AL
      GO TO 69
      END IF
   52 WRITE (NO,54)
   54 FORMAT (T10,'BETA (DEGREES)[90.00]:  ',$)
      READ (NI,36,ERR=52,END=9911) NCR,BE
      IF (NCR.EQ.0) BE=90.00
      IF (BE.LT.30.0 .OR. BE.GT.150.0) GO TO 52
      IF (NCR.NE.0) THEN
      IF (MOD(BE,1.0).NE.0.0) GO TO 61
   56 WRITE (NO,58)
   58 FORMAT (T11,'BETA (MINUTES)[0.00]:  ',$)
      READ (NI,36,ERR=56,END=9911) NCR,AMIN
      IF (NCR.EQ.0) AMIN=0.0
      IF (AMIN.LT.0.0 .OR. AMIN.GT.59.9) GO TO 56
      BE=BE+AMIN/60.
      ELSE
      END IF
   61 IF (NSYS.EQ.2) GO TO 69
   62 WRITE (NO,64)
   64 FORMAT (T10,'GAMMA(DEGREES)[90.00]:  ',$)
      READ (NI,36,ERR=62,END=9911) NCR,GA
      IF (NCR.EQ.0) GA=90.00
      IF (GA.LT.30.0 .OR. GA.GT.150.0) GO TO 62
      IF (NCR.NE.0) THEN
      IF (MOD(GA,1.0).NE.0.0) GO TO 69
   66 WRITE (NO,68)
   68 FORMAT (T11,'GAMMA(MINUTES)[0.00]:  ',$)
      READ (NI,36,ERR=66,END=9911) NCR,AMIN
      IF (NCR.EQ.0) AMIN=0.0
      IF (AMIN.LT.0.0 .OR. AMIN.GT.59.9) GO TO 66
      GA=GA+AMIN/60.
      ELSE
      END IF
   69 CONTINUE
	RETURN
9911	A = 0.0
	RETURN
	END
C===========================================================================
	SUBROUTINE GETLAUE(LAUE,NSYS,no,ni)
	INTEGER LAUE,NSYS,NO,NI
C===========================================================================
100	WRITE (NO,21)
   21 FORMAT (/' ENTER THE CRYSTAL SYSTEM AND LAUE GROUP USING THE LAUE 
     *CODE'/
     *T10,'1  TRICLINIC       LAUE GROUP -1                  '/
     *T10,'2  MONOCLINIC      LAUE GROUP 2/M (Y AXIS UNIQUE) '/
     *T10,'3  MONOCLINIC      LAUE GROUP 2/M (Z AXIS UNIQUE) '/
     *T10,'4  ORTHORHOMBIC    LAUE GROUP MMM                 '/
     *T10,'5  TETRAGONAL      LAUE GROUP 4/M                 '/
     *T10,'6  TETRAGONAL      LAUE GROUP 4/MMM               '/
     *T10,'7  CUBIC           LAUE GROUP M3                  '/
     *T10,'8  CUBIC           LAUE GROUP M3M                 ')
      WRITE (NO,22)
   22 FORMAT (
     *T10,'9  RHOMBOHEDRAL    LAUE GROUP -3     RHOMBOHEDRAL AXES '/
     *T10,'10 RHOMBOHEDRAL    LAUE GROUP -3M    RHOMBOHEDRAL AXES '/
     *T10,'11 TRIGONAL        LAUE GROUP -3     HEXAGONAL AXES    '/
     *T10,'12 TRIGONAL        LAUE GROUP -3M1   HEXAGONAL AXES    '/
     *T10,'13 HEXAGONAL       LAUE GROUP 6/M                      '/
     *T10,'14 HEXAGONAL       LAUE GROUP 6/MMM                    '/
     *T10,'15 TRIGONAL        LAUE GROUP -31M   ROTATED           '/
     *T5,'Enter a code number [HELP]: ',$)
      READ (NI,23,ERR=20,END=9911) LAUE
   23 FORMAT (I)

      NSYS=0
      IF (LAUE.EQ.1) NSYS=1
      IF (LAUE.EQ.2) NSYS=2
      IF (LAUE.EQ.3) NSYS=3
      IF (LAUE.EQ.4) NSYS=4
      IF (LAUE.EQ.5 .OR. LAUE.EQ.6) NSYS=5
      IF (LAUE.EQ.7 .OR. LAUE.EQ.8) NSYS=6
      IF (LAUE.EQ.9 .OR. LAUE.EQ.10) NSYS=7
      IF (LAUE.GE.11 .AND. LAUE.LE.15) NSYS=8
	IF (NSYS .GT. 0) RETURN
20      WRITE (NO,16) 
   16 FORMAT (' THE SYMMETRY AND CRYSTAL SYSTEM DEFINE THE D-VALUES AND
     * THE SET OF MILLER'/1X,'INDICES.  FOR MANY SYMMETRY GROUPS A GIVEN 
     * D-VALUE MAY HAVE MORE THAN ONE'/1X,'MILLER INDEX WHICH IS NOT EQV
     *IVALENT TO THE OTHER (FOR EXAMPLE IN THE'/1X,'TETRAGONAL SYSTEM H
     *KL MAY NOT BE EQUIVALENT TO KHL).  IT IS THEN NECESSARY'/1X,'TO GE
     *NERATE ALL THE NONEQUIVALENT MILLER INDICES.  THIS IS DONE BY USIN
     *G'/1X,'WHAT IS KNOWN AS THE LAUE GROUP.  THE LAUE GROUP IS DETERMI
     *NED FROM THE')
      WRITE (NO,14)
   14 FORMAT (' SPACE GROUP BY IGNORING THE LATTICE, REPLACING 
     *ALL TRANSLATIONAL SYMMETRY'/1X,'ELEMENTS WITH NON-TRANSLATIONAL EL
     *EMENTS AND ADDING A CENTER OF SYMMETRY.'/1X,'IF YOU DO NOT KNOW TH
     *E LAUE GROUP, USE THE HIGHEST LAUE NUMBER FOR THE'/1X,'CRYSTAL SY
     *STEM DESIRED (EXCEPT LAUE=15).  IF YOU WANT TO USE THE EXTINCTION' 
     */1X,'CONDITIONS, YOU MUST USE THE PROPER LAUE GROUP FOR THE SPACE 
     *GROUP.'///)
	GOTO 100
9911	NSYS = 0
	LAUE = 0
	RETURN
	END
C===========================================================================
	SUBROUTINE GETEXTN(EXT,lo,li)
C allow the user to input the extinction conditions
	INTEGER EXT(*), lo, li
C The EXT values determine Extinctions (set in GETEXTN, used in TSTEXTN)
C  EXT(1) IS CODE  FOR LATTICE TYPE
C     1  A CENTERED  K+L=2N
C     2  B CENTERED  H+L=2N
C     3  C CENTERED  H+K=2N
C     4  F CENTERED  H+K,K+L,L+H=2N
C     5  I CENTERED  H+K+L=2N
C     6  R OBVERSE   -H+K+L=3N
C     7  R REVERSE   H-K+L=3N
C     8  H HEXAGONAL H-K=3N
C     9  P PRIMITIVE NO CONDITIONS
C  EXT(2) IS CODE WORD FOR (HK0)
C     0  NO CONDITIONS
C     1  A GLIDE     H=2N
C     2  B GLIDE     K=2N
C     3  N GLIDE     H+K=2N
C     4  D GLIDE     H+K=4N (H=2N,K=2N)
C  EXT(3) IS CODE WORD FOR (H0L)
C     0  NO CONDITIONS
C     1  A GLIDE     H=2N
C     2  C GLIDE     L=2N
C     3  N GLIDE     H+L=2N
C     4  D GLIDE     H+L=4N (H=2N,L=2N)
C  EXT(4) IS CODE WORD FOR (0KL)
C     0  NO CONDITIONS
C     1  B GLIDE     K=2N
C     2  C GLIDE     L=2N
C     3  N GLIDE     K+L=2N
C     4  D GLIDE     K+L=4N (K=2N,L=2N)
C  EXT(5) IS CODE WORD FOR (HHL)
C     0  NO CONDITIONS
C     1  C(N) GLIDE  L=2N
C     2  D GLIDE     2H+L=4N
C  EXT(6) IS CODE WORD FOR (H-HL)
C     0  NO CONDITIONS
C     1  C(N) GLIDE  L=2N
C     2  D GLIDE     2H+L=4N
C  EXT(7) IS CODE WORD FOR (HKH)
C     0  NO CONDITIONS
C     1  B(N) GLIDE  K=2N
C     2  D GLIDE     2H+K=4N
C  EXT(8) IS CODE WORD FOR (HK-H)
C     0  NO CONDITIONS
C     1  B(N) GLIDE  K=2N
C     2  D GLIDE     2H+K=4N
C  EXT(9) IS CODE WORD FOR (HKK)
C     0  NO CONDITIONS
C     1  A(N) GLIDE  H=2N
C     2  D GLIDE     H+2K=4N
C  EXT(10) IS CODE WORD FOR (HK-K)
C     0  NO CONDITIONS
C     1  A(N) GLIDE  H=2N
C     2  D GLIDE  2K+H=4N
C  EXT(11) IS CODE WORD FOR (H00)
C     0  NO CONDITIONS
C     1  2(1) OR 4(2) SCREW  H=2N
C     2  4(1) OR 4(3) SCREW  H=4N
C  EXT(12) IS CODE WORD FOR (0K0)
C     0  NO CONDITIONS
C     1  2(1) OR 4(2) SCREW  K=2N
C     2  4(1) OR 4(3) SCREW  K=4N
C  EXT(13) IS CODE WORD FOR (00L)
C     0  NO CONDITIONS
C     1  2(1) OR 4(2) SCREW  L=2N
C     2  3(1),3(2),6(2) OR 6(4) SCREW  L=3N
C     3  4(1) OR 4(3) SCREW  L=4N
C     4  6(1) OR 6(5) SCREW  L=6N
C===========================================================================
	INTEGER NCR, I
	CHARACTER*1 ians


	IF (EXT(1) .eq. 0) THEN
	  WRITE (lo,'(A)') 
	1	'$Is there lattice centering? (Y/[N]): '
	  read(li,'(A1)',err=999) IANS
	ELSE
	  ians = 'Y'
	ENDIF
	IF (ians .eq. 'y' .or. ians .eq. 'Y') THEN
85	  WRITE (lo,86) EXT(1)
86	  FORMAT (/' EXTINCTION DUE TO LATTICE TYPE:'/1X,
     *'        0  P PRIMITIVE (NO CONDITIONS)'/1X,
     *'        1  A CENTERED  K+L=2N'/1X,
     *'        2  B CENTERED  H+L=2N'/1X,
     *'        3  C CENTERED  H+K=2N'/1X,
     *'        4  F CENTERED  H+K,K+L,L+H=2N'/1X,
     *'        5  I CENTERED  H+K+L=2N'/1X,
     *'        6  R OBVERSE   -H+K+L=3N'/1X,
     *'        7  R REVERSE   H-K+L=3N'/1X,
     *'        8  H HEXAGONAL H-K=3N'/1X,
     *' Enter a number for the lattice type [default=',i2,']: ',$)
          READ (li,87,ERR=85,END=999) NCR,I
	  IF (NCR .gt. 0) EXT(1) = I
87	  FORMAT (Q,I)
          IF (EXT(1).LT.0 .OR. EXT(1).GT.9) GO TO 85
	ENDIF

	IF (EXT(2)+EXT(3)+EXT(4) .eq. 0) THEN
	  WRITE (lo,'(2A)') 
	1	'$Are there glide planes perpendicular to ',
	1	'A, B or C? (Y/[N]): '
	  read(li,'(A1)',err=999) IANS
	ELSE
	  ians = 'Y'
	ENDIF
	IF (ians .eq. 'y' .or. ians .eq. 'Y') THEN
88	  WRITE (lo,93) 'C','HK0','A','H','B','K'
	1	,'H+K','H+K','H','K',ext(2)
          READ (li,87,ERR=88,END=999) NCR,I
	  IF (NCR .gt. 0) EXT(2) = I
          IF (EXT(2).LT.0 .OR. EXT(2).GT.4) GO TO 88

90	  WRITE (lo,93) 'B','H0L','A','H','C','L'
	1	,'H+L','H+L','H','L',ext(3)
          READ (li,87,ERR=90,END=999) NCR,I
	  IF (NCR .gt. 0) EXT(3) = I
          IF (EXT(3).LT.0 .OR. EXT(3).GT.4) GO TO 90

92	  WRITE (lo,93) 'A','0KL','B','K','C','L'
	1	,'K+L','K+L','K','L',ext(4)
          READ (li,87,ERR=92,END=999) NCR,I
	  IF (NCR .gt. 0) EXT(4) = I
          IF (EXT(4).LT.0 .OR. EXT(4).GT.4) GO TO 92

93	  FORMAT (/' Glide plane perpendicular to ',A,' [(',A,
	1	') Conditions]:'/1x,
     *'        0  None'/1X,
     *'        1  ',A,' GLIDE     ',A,'=2N',/,1X,
     *'        2  ',A,' GLIDE     ',A,'=2N',/,1X,
     *'        3  N GLIDE     ',a,'=2N',/,1X,
     *'        4  D GLIDE     ',a,'=4N (',a,'=2N,',a,'=2N)',/,1X,
     *' Enter a number for the glide plane [default=',i2,']: ',$)

	ENDIF

	IF (EXT(5)+EXT(6)+EXT(7)+EXT(8)+EXT(9)+EXT(10) .eq. 0) THEN
	  WRITE (lo,'(2A)') 
	1	'$Are there glide planes perpendicular to ',
	1	'A+B, B+C or A+C? (Y/[N]): '
	  read(li,'(A1)',err=999) IANS
	ELSE
	  ians = 'Y'
	ENDIF

	IF (ians .eq. 'y' .or. ians .eq. 'Y') THEN
95	  WRITE (lo,94) 'A+B [(HHL','C(N)','L','2H+L',EXT(5)
94	  FORMAT (/' Glide plane perpendicular to ',a,
	1	') Conditions]:'/1x,
     *'        0  None'/1X,
     *'        1  ',A,' GLIDE  ',A,'=2N'/1X,
     *'        2  D GLIDE     ',a,'=4N'/1X,
     *' Enter a number for the glide plane [default=',i2,']: ',$)
          READ (li,87,ERR=95,END=999) NCR,I
	  IF (NCR .gt. 0) EXT(5) = I
          IF (EXT(5).LT.0 .OR. EXT(5).GT.2) GO TO 95

96	  WRITE (lo,94) 'A+B [(H-HL','C(N)','L','2H+L',EXT(6)
          READ (li,87,ERR=96,END=999) NCR,I
	  IF (NCR .gt. 0) EXT(6) = I
          IF (EXT(6).LT.0 .OR. EXT(6).GT.2) GO TO 96

97	  WRITE (lo,94) 'A+C [(HKH','B(N)','K','2H+K',EXT(7)
          READ (li,87,ERR=97,END=999) NCR,I
	  IF (NCR .gt. 0) EXT(7) = I
          IF (EXT(7).LT.0 .OR. EXT(7).GT.2) GO TO 97

98	  WRITE (lo,94) 'A+C [(HK-H','B(N)','K','2H+K',EXT(8)
          READ (li,87,ERR=98,END=999) NCR,I
	  IF (NCR .gt. 0) EXT(8) = I
          IF (EXT(8).LT.0 .OR. EXT(8).GT.2) GO TO 98

99	  WRITE (lo,94) 'B+C [(HKK','A(N)','H','2K+H',EXT(9)
          READ (li,87,ERR=99,END=999) NCR,I
	  IF (NCR .gt. 0) EXT(9) = I
          IF (EXT(9).LT.0 .OR. EXT(9).GT.2) GO TO 99

910	  WRITE (lo,94) 'B+C [(HK-K','A(N)','H','2K+H',EXT(10)
          READ (li,87,ERR=910,END=999) NCR,I
	  IF (NCR .gt. 0) EXT(10) = I
          IF (EXT(10).LT.0 .OR. EXT(10).GT.2) GO TO 910
	ENDIF

	IF (EXT(11)+EXT(12)+EXT(13) .eq. 0) THEN
	  WRITE (lo,'(A)') 
	1	'$Are there screw axis || to A, B or C? (Y/[N]): '
	  read(li,'(A1)',err=999) IANS
	ELSE
	  ians = 'Y'
	ENDIF
	IF (ians .eq. 'y' .or. ians .eq. 'Y') THEN
811	  WRITE (lo,810) 'A [(H00','H','H',ext(11)
810	  FORMAT (/' Screw axis || to ',A,') Conditions]:'/1x,
     *'        0  None'/1X,
     *'        1  2_1 OR 4_2 SCREW  ',a,'=2N'/1X,
     *'        2  4_1 OR 4_3 SCREW  ',a,'=4N'/1X,
     *' Enter a number for the screw axis [default=',i2,']: ',$)
          READ (li,87,ERR=811,END=999) NCR,I
	  IF (NCR .gt. 0) EXT(11) = I
          IF (EXT(11).LT.0 .OR. EXT(11).GT.2) GO TO 811

812	  WRITE (lo,810) 'B [(0K0','K','K',ext(12)
          READ (li,87,ERR=812,END=999) NCR,I
	  IF (NCR .gt. 0) EXT(12) = I
          IF (EXT(12).LT.0 .OR. EXT(12).GT.2) GO TO 812

813	  WRITE (lo,814) ext(13)
814	  FORMAT (/' Screw axis || to C [(00L) Conditions]:'/1x,
     *'        0  None'/1X,
     *'        1  2_1, 4_2 or 6_3 SCREW       L=2N'/1X,
     *'        2  3_1, 3_2, 6_2 or 6_4 SCREW  L=3N'/1X,
     *'        3  4_1 or 4_3 SCREW            L=4N'/1X,
     *'        4  6_1 or 6_5 SCREW            L=6N'/1X,
     *' Enter a number for the screw axis [default=',i2,']: ',$)
          READ (li,87,ERR=813,END=999) NCR,I
	  IF (NCR .gt. 0) EXT(13) = I
          IF (EXT(13).LT.0 .OR. EXT(13).GT.4) GO TO 813
	ENDIF

999	RETURN
	END
C===========================================================================
	LOGICAL FUNCTION TSTEXTN(NH,NK,NL,EXT)
C     
C     REJECT SYSTEMATICALLY EXTINCT REFLECTIONS
C
C RETURNS TRUE IF REFLECTION (NH,NK,NL) IS EXTINCT BY THE SYMMETRY CODES
C IN EXT().
C
      INTEGER EXT(*)
      INTEGER NH, NK, NL
C
C TAKEN FROM NATIONAL BUREAU OF STANDARDS D-SPACING PROGRAM (MODIFIED 7/7/79)
C AS MODIFIED 20 MAR 1988 BY D K SMITH
C   (rewritten using structured FORTRAN -- B.H. Toby 1/92)
C===========================================================================
C
C     EXT(1) - LATTICE TYPE (ALL HKL)
C     
	IF (EXT(1) .eq. 1) THEN
C     A (K+L=2N)
	  IF (MOD(NK+NL,2) .ne. 0) goto 2000
	ELSEIF (EXT(1) .eq. 2) THEN
C     B (H+L=2N)
	  IF (MOD(NH+NL,2) .ne. 0) goto 2000
	ELSEIF (EXT(1) .eq. 3) THEN
C     C (H+K=2N)
	  IF (MOD(NH+NK,2) .ne. 0) goto 2000
	ELSEIF (EXT(1) .eq. 4) THEN
C     F (H+K=2N,H+L=2N,K+L=2N)
	  IF (MOD(NK*NH*NL,2).EQ.0 .and. 
	1	(MOD(NH,2).NE.0 .or. MOD(NK,2).NE.0 .or. 
	1	MOD(NL,2).NE. 0)) goto 2000
	ELSEIF (EXT(1) .eq. 5) THEN
C     I (H+K+L=2N)
	  IF (MOD(NH+NK+NL,2) .NE. 0) goto 2000
	ELSEIF (EXT(1) .eq. 7) THEN
C     R(R) (H-K+L=3N)
	  IF (MOD(NH-NK+NL,3).NE.0) GO TO 2000
	ELSEIF (EXT(1) .eq. 8) THEN
C     H (H-K=3N)
	  IF (MOD(NH-NK,3).NE.0) GO TO 2000
	ELSEIF (EXT(1) .eq. 6) THEN
C     R(O) (-H+K+L=3N)
	  IF (MOD(NK-NH+NL,3).NE.0) GO TO 2000
	ENDIF
C----------------------------------------------------------------------
C     EXT(2) (HK0 ZONE)
C     
	IF (NL.EQ.0 .and. ext(2) .eq. 1) THEN
C     A-GLIDE (H=2N)
	  IF (MOD(nh,2).NE.0) GO TO 2000
	ELSEIF (NL.EQ.0 .and. ext(2) .eq. 2) THEN
C     B-GLIDE (K=2N)
	  IF (MOD(nk,2).NE.0) GO TO 2000
	ELSEIF (NL.EQ.0 .and. ext(2) .eq. 3) THEN
C     N-GLIDE (H+K=2N)
	  IF (MOD(NH+NK,2).NE.0) GO TO 2000
	ELSEIF (NL.EQ.0 .and. ext(2) .eq. 4) THEN
C     D-GLIDE (H+K=4N,H=2N,K=2N)
	  IF (MOD(NH+NK,4).NE.0 .OR. MOD(NK,2).NE.0) GO TO 2000
	ENDIF
C----------------------------------------------------------------------
C     
C     EXT(3) (H0L ZONE)
C     
	IF (NK.EQ.0 .and. ext(3) .eq. 1) THEN
C     A-GLIDE (H=2N)
	  IF (MOD(NH,2).NE.0) GO TO 2000
	ELSEIF (NK.EQ.0 .and. ext(3) .eq. 2) THEN
C     C-GLIDE (L=2N)
	  IF (MOD(NL,2).NE.0) GO TO 2000
	ELSEIF (NK.EQ.0 .and. ext(3) .eq. 3) THEN
C     N-GLIDE (H+L=2N)
	  IF (MOD(NH+NL,2).NE.0) GO TO 2000
	ELSEIF (NK.EQ.0 .and. ext(3) .eq. 4) THEN
C     D-GLIDE (H+L=4N,H=2N,L=2N)
	  IF (MOD(NH+NL,4).NE.0 .OR. MOD(NL,2).NE.0) GO TO 2000
	ENDIF
C----------------------------------------------------------------------
C     
C     EXT(4) (0KL ZONE)
C     
	IF (NH.EQ.0 .and. ext(4) .eq. 1) THEN
C     B-GLIDE (K=2N)
	  IF (MOD(NK,2).NE.0) GO TO 2000
	ELSEIF (NH.EQ.0 .and. ext(4) .eq. 2) THEN
C     C-GLIDE (L=2N)
	  IF (MOD(NL,2).NE.0) GO TO 2000
	ELSEIF (NH.EQ.0 .and. ext(4) .eq. 3) THEN
C     N-GLIDE (K+L=2N)
	  IF (MOD(NK+NL,2).NE.0) GO TO 2000
	ELSEIF (NH.EQ.0 .and. ext(4) .eq. 4) THEN
C     D-GLIDE (K+L=4N,K=2N,L=2N)
	  IF (MOD(NK+NL,4).NE.0 .OR. MOD(NL,2).NE.0) GO TO 2000
	ENDIF
C----------------------------------------------------------------------
C     
C     EXT(5) (HHL ZONE)
C     
	IF (NH-NK.EQ.0 .and. EXT(5) .eq. 1) THEN
C     C(N)-GLIDE (L=2N)
	  IF (MOD(NL,2).NE.0) GO TO 2000
	ELSEIF (NH-NK.EQ.0 .and. EXT(5) .eq. 2) THEN
C     D-GLIDE (2H+L=4N)
	  IF (MOD(2*NH+NL,4).NE.0) GO TO 2000
	ENDIF
C----------------------------------------------------------------------
C     
C     EXT(6) (-HHL ZONE)
C     
	IF (NH+NK.EQ.0 .and. EXT(6) .eq. 1) THEN
C     C(N)-GLIDE (L=2N)
	  IF (MOD(NL,2).NE.0) GO TO 2000
	ELSEIF (NH+NK.EQ.0 .and. EXT(6) .eq. 2) THEN
C     D-GLIDE (2H+L=4N)
	  IF (MOD(2*NH+NL,4).NE.0) GO TO 2000
	ENDIF
C----------------------------------------------------------------------
C     
C     EXT(7) (HKH ZONE)
C     
	IF (NH-NL.EQ.0 .and. EXT(7) .eq. 1) THEN
C     B(N)-GLIDE (K=2N)
	  IF (MOD(NK,2).NE.0) GO TO 2000
	ELSEIF (NH-NL.EQ.0 .and. EXT(7) .eq. 2) THEN
C     D-GLIDE (2H+K=4N)
	  IF (MOD(2*NH+NK,4).NE.0) GO TO 2000
	ENDIF
C----------------------------------------------------------------------
C     
C     EXT(8) (-HKH ZONE)
C     
	IF (NH+NL.EQ.0 .and. EXT(8) .eq. 1) THEN
C     B(N)-GLIDE (K=2N)
	  IF (MOD(NK,2).NE.0) GO TO 2000
	ELSEIF (NH+NL.EQ.0 .and. EXT(8) .eq. 2) THEN
C     D-GLIDE (2H+K=4N)
	  IF (MOD(2*NH+NK,4).NE.0) GO TO 2000
	ENDIF
C----------------------------------------------------------------------
C     
C     EXT(9) (HKK ZONE)
C     
	IF (NL-NK.EQ.0 .and. EXT(9) .eq. 1) THEN
C     A(N)-GLIDE (H=2N)
	  IF (MOD(NH,2).NE.0) GO TO 2000
	ELSEIF (NL-NK.EQ.0 .and. EXT(9) .eq. 2) THEN
C     D-GLIDE (2K+H=4N)
	  IF (MOD(2*NK+NH,4).NE.0) GO TO 2000
	ENDIF
C----------------------------------------------------------------------
C     
C     EXT(10) (H-KK ZONE)
C     
	IF (NL+NK.EQ.0 .and. EXT(10) .eq. 1) THEN
C     A(N)-GLIDE (H=2N)
	  IF (MOD(NH,2).NE.0) GO TO 2000
	ELSEIF (NL+NK.EQ.0 .and. EXT(10) .eq. 2) THEN
C     D-GLIDE (2K+H=4N)
	  IF (MOD(2*NK+NH,4).NE.0) GO TO 2000
	ENDIF
C----------------------------------------------------------------------
C     
C     EXT(11) (H00 ROW)
C     
	IF (NK.EQ.0 .and. NL.EQ.0 .and. EXT(11) .eq. 1) THEN
C     2(1) OR 4(2) SCREW H=2N)
	  IF (MOD(NH,2).NE.0) GO TO 2000
	ELSEIF (NK.EQ.0 .and. NL.EQ.0 .and. EXT(11) .eq. 2) THEN
C     4(1) OR 4(3) SCREW (H=4N)
	  IF (MOD(NH,4).NE.0) GO TO 2000
	ENDIF
C----------------------------------------------------------------------
C     
C     EXT(12) (0K0 ROW)
C     
	IF (NH.EQ.0 .and. NL.EQ.0 .and. EXT(12) .eq. 1) THEN
C     2(1) OR 4(2) SCREW K=2N)
	  IF (MOD(NK,2).NE.0) GO TO 2000
	ELSEIF (NH.EQ.0 .and. NL.EQ.0 .and. EXT(12) .eq. 2) THEN
C     4(1) OR 4(3) SCREW (K=4N)
	  IF (MOD(NK,4).NE.0) GO TO 2000
	ENDIF
C----------------------------------------------------------------------
C     
C     EXT(13) (00L ROW)
C     
	IF (NH.EQ.0 .and. NK.EQ.0 .and. EXT(13) .eq. 1) THEN
C     2(1) OR 4(2) SCREW L=2N)
	  IF (MOD(NL,2).NE.0) GO TO 2000
	ELSEIF (NH.EQ.0 .and. NK.EQ.0 .and. EXT(13) .eq. 2) THEN
C     3(1),3(2),6(2) OR 6(4) SCREW (L=3N)
	  IF (MOD(NL,3).NE.0) GO TO 2000
	ELSEIF (NH.EQ.0 .and. NK.EQ.0 .and. EXT(13) .eq. 3) THEN
C     4(1) OR 4(3) SCREW (L=4N)
	  IF (MOD(NL,4).NE.0) GO TO 2000
	ELSEIF (NH.EQ.0 .and. NK.EQ.0 .and. EXT(13) .eq. 4) THEN
C     6(1) OR 6(5) SCREW (L=6N)
	  IF (MOD(NL,6).NE.0) GO TO 2000
	ENDIF

C     ELIMINATE 000 REFLECTION
	IF (NH.EQ.0 .and. NK.EQ.0 .and. NL .eq. 0) GO TO 2000

C  NOT EXTINCT
	TSTEXTN = .FALSE.
	RETURN

C  EXTINCT
2000	TSTEXTN = .TRUE.
	RETURN 

	END
C===========================================================================
      SUBROUTINE CELCON(A,B,C,AL,BE,GA,RCV,RCP,VOL,PRCP)
C     
C     CALCULATE RECIPROCAL CELL CONSTANTS
C     
      REAL A,B,C,AL,BE,GA,RCV,RCP(6),VOL,PRCP(6)
C===========================================================================

	REAL RAD
	PARAMETER (RAD = 3.1415926/180.0)
	REAL COSAL,COSBE,COSGA,SINAL,SINBE,SINGA

      COSAL = COS(AL*RAD)
      COSBE = COS(BE*RAD)
      COSGA = COS(GA*RAD)
      SINAL = SQRT(1.0-COSAL**2)
      SINBE = SQRT(1.0-COSBE**2)
      SINGA = SQRT(1.0-COSGA**2)
      RCV = A*B*C*SQRT(
	1	1.0-COSAL**2-COSBE**2-COSGA**2+2.0*COSAL*COSBE*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	1.0-RCP(4)**2-RCP(5)**2-RCP(6)**2+2.0*RCP(4)*RCP(5)*RCP(6))
C     
C     CALCULATE CROSS PRODUCTS OF RECIPROCAL CELL PARAMETERS
C     
      PRCP(1) = RCP(1)**2*0.25
      PRCP(2) = RCP(2)**2*0.25
      PRCP(3) = RCP(3)**2*0.25
      PRCP(4) = RCP(1)*RCP(2)*RCP(6)*0.25
      PRCP(5) = RCP(1)*RCP(3)*RCP(5)*0.25
      PRCP(6) = RCP(2)*RCP(3)*RCP(4)*0.25
      RETURN
      END
C===========================================================================
      SUBROUTINE DGENSUB(LO,LI,LAUE,DMIN,DMAX,NREF,IEXT,
	1	A,B,C,rcp,prcp,
	1	NRMAX,D,IH,IK,IL)
C     
C     GENERATE MILLER INDICES ACCORDING TO LAUE GROUP WITHOUT REDUNDANCIES
C
	INTEGER LO, LI, LAUE, IEXT(13), NRMAX, NREF
	REAL D(NRMAX)
	REAL DMIN, DMAX
	INTEGER*2 IH(NRMAX),IK(NRMAX),IL(NRMAX)
	real A,B,C,RCP(6),PRCP(6)
C TAKEN FROM NATIONAL BUREAU OF STANDARDS D-SPACING PROGRAM (MODIFIED 7/7/79)
C AS MODIFIED 20 MAR 1988 BY D K SMITH
C   (rewritten using structured FORTRAN -- B.H. Toby 1/92)
C===========================================================================
C     
C  THE CODE WORD, LAUE, HAS THE FOLLOWING INTERPRETATION
C     LAUE = 1   LAUE GROUP -1
C     LAUE = 2   LAUE GROUP 2/M (Y AXIS UNIQUE)
C     LAUE = 3   LAUE GROUP 2/M (Z AXIS UNIQUE)
C     LAUE = 4   LAUE GROUP MMM
C     LAUE = 5   LAUE GROUP 4/M
C     LAUE = 6   LAUE GROUP 4/MMM
C     LAUE = 7   LAUE GROUP M3
C     LAUE = 8   LAUE GROUP M3M
C     LAUE = 9   LAUE GROUP -3     RHOMBOHEDRAL AXES
C     LAUE =10   LAUE GROUP -3M    RHOMBOHEDRAL AXES
C     LAUE =11   LAUE GROUP -3     HEXAGONAL AXES
C     LAUE =12   LAUE GROUP -3M    HEXAGONAL AXES
C     LAUE =13   LAUE GROUP 6/M
C     LAUE =14   LAUE GROUP 6/MMM
C     LAUE =15   LAUE GROUP -3M    ROTATED
C
        INTEGER NH,NK,NL
	LOGICAL TSTEXTN

	INTEGER K,I,J,JJ,II,KK,MN,KQ,NINDX,IQ,JQ,M,N
	REAL HH(100),HK(100),HL(100)
	REAL RDMAX,RECD
	INTEGER HMAX, KMAX, LMAX


	RDMAX = 0.0
!	IF (DMAX .gt. 0) RDMAX = 1./DMAX
	DO I=1,100
          HH(I)=0.0
          HK(I)=0.0
 	  HL(I)=0.0
	ENDDO
	DO I=1,NRMAX
	  D(I)=0.0
	  IH(I)=0
	  IK(I)=0
	  IL(I)=0
	ENDDO
C     
C     DETERMINE MAXIMUM INDICES
C     
C     SET UP RANGE OF MILLER INDICES
C     
	LMAX = C/DMIN+1.0
	DO K=1,LMAX
	  HL(K+1)=HL(K)+1.0
	ENDDO
C     
	KMAX = B/DMIN+1.0
	IF (LAUE.NE.1) THEN
	  DO J=1,KMAX
	    HK(J+1)=HK(J)+1.0
	  ENDDO
	ELSE
	  KMAX = 2*KMAX-1
	  N = B/DMIN
	  HK(1)=N
	  DO J=1,KMAX
	    HK(J+1)=HK(J)-1.0
	  ENDDO
	ENDIF
C     
	HMAX = A/DMIN+1.0
	IF (LAUE.GT.3) THEN
	  DO I=1,HMAX
	    HH(I+1)=HH(I)+1.0
	  ENDDO
	ELSE
	  HMAX = 2*HMAX-1
	  M = A/DMIN
	  HH(1) = M
	  DO I=1,HMAX
	    HH(I+1)=HH(I)-1.0
	  ENDDO
	ENDIF
C     
  200 nref=0
      GO TO (220,260,310,350,350,820,460,460,460,460,820,820,820,820,820
     1), LAUE
C----------------------------------------------------------------------
C     
C     TRICLINIC SYSTEM - LAUE GROUP -1
C     
C     GENERATE MILLER INDICES (H00)
C     
  220 NK=0
      NL=0
      DO 250  I=1,M
	NH=HH(I)
	IF (nref+1.GT.NRMAX) GOTO 210
C Test for extinction
	IF (.not. TSTEXTN(NH,NK,NL,IEXT)) THEN
C  Compute D-space
	  RECD=HH(I)*RCP(1)
C if in range, add to the arrays
          IF (RECD .GT. RDMAX) THEN
	    D(nref+1) = 1.0 / RECD
	    IH(nref+1) = NH
	    IK(nref+1) = NK
	    IL(nref+1) = NL
	    IF (D(nref+1) .ge. DMIN) nref = nref + 1
	  ENDIF
	ENDIF
  250 CONTINUE
C     
C     GENERATE MILLER INDICES (HK0) (K NOT ZERO)
C     
      JQ=N
      IQ=HMAX
      NINDX = 400 
      GO TO 270
C----------------------------------------------------------------------
C     
C     MONOCLINIC SYSTEM (Y UNIQUE)
C     
  260 JQ=KMAX
      IQ=M+1
      NINDX = 410 
      NL=0
C     
C     GENERATE MILLER INDICES (HK0)
C     
  270 DO 300  J=1,JQ
      DO 300  I=1,IQ
      NH=HH(I)
      NK=HK(J)
	IF (nref+1.GT.NRMAX) GOTO 210
C Test for extinction
	IF (.not. TSTEXTN(NH,NK,NL,IEXT)) THEN
C  Compute D-space
	  IF (LAUE.EQ.1) THEN
	    RECD=HH(I)*HK(J)*PRCP(4)*2.0
	    RECD=2.0*SQRT(HH(I)**2*PRCP(1)+HK(J)**2*PRCP(2)+RECD)
	  ELSE
	    RECD=2.0*SQRT(HH(I)**2*PRCP(1)+HK(J)**2*PRCP(2))
	  ENDIF
C if in range, add to the arrays
          IF (RECD .GT. RDMAX) THEN
	    D(nref+1) = 1.0 / RECD
	    IH(NREF+1) = NH
	    IK(NREF+1) = NK
	    IL(NREF+1) = NL
	    IF (D(NREF+1) .ge. DMIN) NREF = NREF + 1
	  ENDIF
	ENDIF
  300 CONTINUE
C     
C     GENERATE MILLER INDICES (HKL) (L NOT ZERO)
C     
      JQ=1
      KQ=2
      GO TO 370
C----------------------------------------------------------------------
C     
C     MONOCLINIC SYSTEM (Z UNIQUE)
C     
  310 MN=M+1
      NK=0
C     
C     GENERATE MILLER INDICES (H0L)
C     
      DO 340  K=1,LMAX
      DO 340  I=1,MN
	NH=HH(I)
	NL=HL(K)
	IF (NREF+1.GT.NRMAX) GOTO 210
C Test for extinction
	IF (.not. TSTEXTN(NH,NK,NL,IEXT)) THEN
C  Compute D-space
	  RECD=2.0*SQRT(HH(I)**2*PRCP(1)+HL(K)**2*PRCP(3))
C if in range, add to the arrays
          IF (RECD .GT. RDMAX) THEN
	    D(NREF+1) = 1.0 / RECD
	    IH(NREF+1) = NH
	    IK(NREF+1) = NK
	    IL(NREF+1) = NL
	    IF (D(NREF+1) .ge. DMIN) NREF = NREF + 1
	  ENDIF
	ENDIF
  340 CONTINUE
C     
C     GENERATE MILLER INDICES (HKL) (K NOT ZERO)
C     
      JQ=2
      NINDX = 420 
      GO TO 360
C----------------------------------------------------------------------
C     
C     ORTHORHOMBIC SYSTEM
C     
C     TETRAGONAL SYSTEM - 4/M
C     
C     
C     GENERATE MILLER INDICES (HKL)
C     
  350 NINDX = 430 
      JQ=1
  360 KQ=1
  370 DO 450  K=KQ,LMAX
      DO 450  J=JQ,KMAX
      DO 450  I=1,HMAX
	IF (LAUE.EQ.5 .and. HH(I).EQ.0.AND.HK(J).NE.0.0) GO TO 450
  	NH=HH(I)
	NK=HK(J)
	NL=HL(K)
	IF (NREF+1.GT.NRMAX) GOTO 210
C Test for extinction
	IF (.not. TSTEXTN(NH,NK,NL,IEXT)) THEN
C  Compute D-space
	  RECD=HH(I)**2*PRCP(1)+HK(J)**2*PRCP(2)+HL(K)**2*PRCP(3)
	  IF (NINDX .eq. 400) THEN
	    RECD=RECD+2.0*(HH(I)*HK(J)*PRCP(4)+HH(I)*HL(K)*PRCP(5)+
	1	HK(J)*HL(K)*PRCP(6))
	  ELSEIF (NINDX .eq. 410) THEN
	    RECD=RECD+2.0*HH(I)*HL(K)*PRCP(5)
	  ELSEIF (NINDX .eq. 420) THEN
	    RECD=RECD+2.0*HH(I)*HK(J)*PRCP(4)
	  ENDIF
	  RECD=2.0*SQRT(RECD)
C if in range, add to the arrays
          IF (RECD .GT. RDMAX) THEN
	    D(NREF+1) = 1.0 / RECD
	    IH(NREF+1) = NH
	    IK(NREF+1) = NK
	    IL(NREF+1) = NL
	    IF (D(NREF+1) .ge. DMIN) NREF = NREF + 1
	  ENDIF
	ENDIF
  450 CONTINUE
      GO TO 1100
C----------------------------------------------------------------------
C     
C     RHOMBOHEDRAL SYSTEM
C     RHOMBOHEDRAL LAUE GROUPS WITH RHOMBOHEDRAL LATTICE
C     
C     CUBIC SYSTEM
C     
C     
C     LAUE GROUPS -3 AND -3M            3M AND M3M
C     
C     
C     GENERATE MILLER INDICES (HKL) (H GREATER OR EQUAL K GREATER OR
C     EQUAL +L)
C     
  460 DO 490  K=1,HMAX
      DO 490  J=K,HMAX
      DO 490  I=J,HMAX
      NH=HH(I)
      NK=HK(J)
      NL=HL(K)
	IF (NREF+1.GT.NRMAX) GOTO 210
C Test for extinction
	IF (.not. TSTEXTN(NH,NK,NL,IEXT)) THEN
C  Compute D-space
	IF (LAUE.GE.9) THEN
	  RECD=2.0*RCP(4)*(HH(I)*HK(J)+HH(I)*HL(K)+HK(J)*HL(K))
	  RECD=2.0*SQRT((HH(I)**2+HK(J)**2+HL(K)**2+RECD)*PRCP(1))
	ELSE
	  RECD=2.0*SQRT((HH(I)**2+HK(J)**2+HL(K)**2)*PRCP(1))
	ENDIF
C if in range, add to the arrays
          IF (RECD .GT. RDMAX) THEN
	    D(NREF+1) = 1.0 / RECD
	    IH(NREF+1) = NH
	    IK(NREF+1) = NK
	    IL(NREF+1) = NL
	    IF (D(NREF+1) .ge. DMIN) NREF = NREF + 1
	  ENDIF
	ENDIF
  490 CONTINUE
      IF (LAUE-8) 630,1100,500
C     
C     GENERATE MILLER INDICES (HKL) (-H GREATER OR EQUAL K GREATER OR
C     EQUAL L, K NOT ZERO)
C     
  500 DO 540  K=1,HMAX
      DO 540  J=K,HMAX
      IF (HK(J).EQ.0) GO TO 540
      DO 530  I=J,HMAX
      NH=-HH(I)
      NK=HK(J)
      NL=HL(K)
	IF (NREF+1.GT.NRMAX) GOTO 210
C Test for extinction
	IF (.not. TSTEXTN(NH,NK,NL,IEXT)) THEN
C  Compute D-space
	  RECD=2.0*SQRT((HH(I)**2+HK(J)**2+HL(K)**2+2.0*RCP(4)*
	1	(-HH(I)*HK(J)-HH(I)*HL(K)+HK(J)*HL(K)))*PRCP(1))
C if in range, add to the arrays
          IF (RECD .GT. RDMAX) THEN
	    D(NREF+1) = 1.0 / RECD
	    IH(NREF+1) = NH
	    IK(NREF+1) = NK
	    IL(NREF+1) = NL
	    IF (D(NREF+1) .ge. DMIN) NREF = NREF + 1
	  ENDIF
	ENDIF
  530 CONTINUE
  540 CONTINUE
C     
C     GENERATE MILLER INDICES (HKL) (H GREATER OR EQUAL K GREATER OR
C     EQUAL -L, L NOT ZERO)
C     
      DO 570  K=2,HMAX
      DO 570  J=K,HMAX
      DO 570  I=J,HMAX
      IF (HH(I).EQ.HL(K)) GO TO 570
      NH=HH(I)
      NK=HK(J)
      NL=-HL(K)
	IF (NREF+1.GT.NRMAX) GOTO 210
C Test for extinction
	IF (.not. TSTEXTN(NH,NK,NL,IEXT)) THEN
C  Compute D-space
	  RECD=2.0*SQRT((HH(I)**2+HK(J)**2+HL(K)**2+2.0*RCP(4)*
	1	(HH(I)*HK(J)-HH(I)*HL(K)-HK(J)*HL(K)))*PRCP(1))
C if in range, add to the arrays
          IF (RECD .GT. RDMAX) THEN
	    D(NREF+1) = 1.0 / RECD
	    IH(NREF+1) = NH
	    IK(NREF+1) = NK
	    IL(NREF+1) = NL
	    IF (D(NREF+1) .ge. DMIN) NREF = NREF + 1
	  ENDIF
	ENDIF
  570 CONTINUE
C     
C     GENERATE MILLER INDICES (HKL) (H GREATER THAN -K GREATER THAN L)
C     
      DO 620  K=2,HMAX
      KK=K+1
      IF (KK.GT.HMAX) GO TO 620
      DO 610  J=KK,HMAX
      JJ=J+1
      IF (JJ.GT.HMAX) GO TO 610
      DO 600  I=JJ,HMAX
      NH=HH(I)
      NK=-HK(J)
      NL=HL(K)
	IF (NREF+1.GT.NRMAX) GOTO 210
C Test for extinction
	IF (.not. TSTEXTN(NH,NK,NL,IEXT)) THEN
C  Compute D-space
	  RECD=2.0*SQRT((HH(I)**2+HK(J)**2+HL(K)**2+2.0*RCP(4)*
	1	(-HH(I)*HK(J)+HH(I)*HL(K)-HK(J)*HL(K)))*PRCP(1))
C if in range, add to the arrays
          IF (RECD .GT. RDMAX) THEN
	    D(NREF+1) = 1.0 / RECD
	    IH(NREF+1) = NH
	    IK(NREF+1) = NK
	    IL(NREF+1) = NL
	    IF (D(NREF+1) .ge. DMIN) NREF = NREF + 1
	  ENDIF
	ENDIF
  600 CONTINUE
  610 CONTINUE
  620 CONTINUE
      IF (LAUE.NE.9) GO TO 1100
C     
C     GENERATE MILLER INDICES (HKL) (H LESS THAN K LESS THAN L)
C     
  630 DO 680  I=1,HMAX
      II=I+1
      IF (II.GT.HMAX) GO TO 680
      DO 670  J=II,HMAX
      JJ=J+1
      IF (JJ.GT.HMAX) GO TO 670
      DO 660  K=JJ,HMAX
      NH=HH(I)
      NK=HK(J)
      NL=HL(K)
	IF (NREF+1.GT.NRMAX) GOTO 210
C Test for extinction
	IF (.not. TSTEXTN(NH,NK,NL,IEXT)) THEN
C  Compute D-space
          IF (LAUE.NE.7) THEN
	    RECD=2.0*RCP(4)*(HH(I)*HH(J)+HH(I)*HL(K)+HK(J)*HL(K))
	    RECD=2.0*SQRT((HH(I)**2+HK(J)**2+HL(K)**2+RECD)*PRCP(1))
	  ELSE
            RECD=2.0*SQRT((HH(I)**2+HK(J)**2+HL(K)**2)*PRCP(1))
	  ENDIF
C if in range, add to the arrays
          IF (RECD .GT. RDMAX) THEN
	    D(NREF+1) = 1.0 / RECD
	    IH(NREF+1) = NH
	    IK(NREF+1) = NK
	    IL(NREF+1) = NL
	    IF (D(NREF+1) .ge. DMIN) NREF = NREF + 1
	  ENDIF
	ENDIF
  660 CONTINUE
  670 CONTINUE
  680 CONTINUE
      IF (LAUE.EQ.7) GO TO 1100
C     
C     GENERATE MILLER INDICES (HKL) (H LESS THAN -K LESS THAN L)
C     
      DO 730  I=1,HMAX
      II=I+1
      IF (II.GT.HMAX) GO TO 730
      DO 720  J=II,HMAX
      JJ=J+1
      IF (JJ.GT.HMAX) GO TO 720
      DO 710  K=JJ,HMAX
      NH=HH(I)
      NK=-HK(J)
      NL=HL(K)
	IF (NREF+1.GT.NRMAX) GOTO 210
C Test for extinction
	IF (.not. TSTEXTN(NH,NK,NL,IEXT)) THEN
C  Compute D-space
	  RECD=2.0*SQRT((HH(I)**2+HK(J)**2+HL(K)**2+2.0*RCP(4)*
	1	(-HH(I)*HK(J)+HH(I)*HL(K)-HK(J)*HL(K)))*PRCP(1))
C if in range, add to the arrays
          IF (RECD .GT. RDMAX) THEN
	    D(NREF+1) = 1.0 / RECD
	    IH(NREF+1) = NH
	    IK(NREF+1) = NK
	    IL(NREF+1) = NL
	    IF (D(NREF+1) .ge. DMIN) NREF = NREF + 1
	  ENDIF
	ENDIF
  710 CONTINUE
  720 CONTINUE
  730 CONTINUE
C     
C     GENERATE MILLER INDICES (HKL) (-H LESS OR EQUAL K LESS THAN L,
C     H NOT ZERO)
C     
      DO 770  I=2,HMAX
      DO 770  J=I,HMAX
      JJ=J+1
      IF (JJ.GT.HMAX) GO TO 770
      DO 760  K=JJ,HMAX
      NH=-HH(I)
      NK=HK(J)
      NL=HL(K)
	IF (NREF+1.GT.NRMAX) GOTO 210
C Test for extinction
	IF (.not. TSTEXTN(NH,NK,NL,IEXT)) THEN
C  Compute D-space
	  RECD=2.0*SQRT((HH(I)**2+HK(J)**2+HL(K)**2+2.0*RCP(4)*
	1	(-HH(I)*HK(J)-HH(I)*HL(K)+HK(J)*HL(K)))*PRCP(1))
C if in range, add to the arrays
          IF (RECD .GT. RDMAX) THEN
	    D(NREF+1) = 1.0 / RECD
	    IH(NREF+1) = NH
	    IK(NREF+1) = NK
	    IL(NREF+1) = NL
	    IF (D(NREF+1) .ge. DMIN) NREF = NREF + 1
	  ENDIF
	ENDIF
  760 CONTINUE
  770 CONTINUE
C     
C     GENERATE MILLER INDICES (HKL) (H LESS THAN K LESS OR EQUAL -L,
C     H NOT ZERO
C     
      DO 810  I=2,HMAX
      II=I+1
      IF (II.GT.HMAX) GO TO 810
      DO 800  J=II,HMAX
      DO 800  K=J,HMAX
      NH=HH(I)
      NK=HK(J)
      NL=-HL(K)
	IF (NREF+1.GT.NRMAX) GOTO 210
C Test for extinction
	IF (.not. TSTEXTN(NH,NK,NL,IEXT)) THEN
C  Compute D-space
	  RECD=2.0*SQRT((HH(I)**2+HK(J)**2+HL(K)**2+2.0*RCP(4)*
	1	(HH(I)*HK(J)-HH(I)*HL(K)-HK(J)*HL(K)))*PRCP(1))
C if in range, add to the arrays
          IF (RECD .GT. RDMAX) THEN
	    D(NREF+1) = 1.0 / RECD
	    IH(NREF+1) = NH
	    IK(NREF+1) = NK
	    IL(NREF+1) = NL
	    IF (D(NREF+1) .ge. DMIN) NREF = NREF + 1
	  ENDIF
	ENDIF
  800 CONTINUE
  810 CONTINUE
      GO TO 1100
C----------------------------------------------------------------------
C     
C     HEXAGONAL SYSTEM
C     
C     TETRAGONAL SYSTEM
C     
C     LAUE GROUPS -3, -3M, 6/M, AND 6/MMM        4/MMM
C     
C     GENERATE MILLER INDICES FOR (HKL) (H GREATER OR EQUAL K ALL +L)
C     
  820 DO 850  J=1,HMAX
      DO 850  I=J,HMAX
      DO 850  K=1,LMAX
      NH=HH(I)
      NK=HK(J)
      NL=HL(K)
	IF (NREF+1.GT.NRMAX) GOTO 210
C Test for extinction
	IF (.not. TSTEXTN(NH,NK,NL,IEXT)) THEN
C  Compute D-space
	  IF (LAUE.NE.6) THEN
	    RECD=2.0*SQRT((HH(I)**2+HK(J)**2+HH(I)*HK(J))*PRCP(1)+
	1	HL(K)**2*PRCP(3))
	  ELSE
	    RECD=2.0*SQRT((HH(I)**2+HK(J)**2)*PRCP(1)+HL(K)**2*PRCP(3))
	  ENDIF
C if in range, add to the arrays
          IF (RECD .GT. RDMAX) THEN
	    D(NREF+1) = 1.0 / RECD
	    IH(NREF+1) = NH
	    IK(NREF+1) = NK
	    IL(NREF+1) = NL
	    IF (D(NREF+1) .ge. DMIN) NREF = NREF + 1
	  ENDIF
	ENDIF
  850 CONTINUE
      IF (LAUE.EQ.6) GO TO 1100
      IF (LAUE-13) 910,860,970
C     
C     LAUE GROUP 6/M
C     
C     GENERATE MILLER INDICES (HKL) (H LESS THAN K ALL +L, H NOT ZERO)
C     
  860 DO 900  I=2,HMAX
      II=I+1
      IF (II.GT.HMAX) GO TO 900
      DO 890  J=II,HMAX
      DO 890  K=1,LMAX
      NH=HH(I)
      NK=HK(J)
      NL=HL(K)
	IF (NREF+1.GT.NRMAX) GOTO 210
C Test for extinction
	IF (.not. TSTEXTN(NH,NK,NL,IEXT)) THEN
C  Compute D-space
	  RECD=2.0*SQRT((HH(I)**2+HK(J)**2+HH(I)*HK(J))*PRCP(1)+HL(K)**2
	1	*PRCP(3))
C if in range, add to the arrays
          IF (RECD .GT. RDMAX) THEN
	    D(NREF+1) = 1.0 / RECD
	    IH(NREF+1) = NH
	    IK(NREF+1) = NK
	    IL(NREF+1) = NL
	    IF (D(NREF+1) .ge. DMIN) NREF = NREF + 1
	  ENDIF
	ENDIF
  890 CONTINUE
  900 CONTINUE
C     
C     LAUE GROUP -3M
C     
C     GENERATE MILLER INDICES (HKL) (H LESS THAN K, L NOT ZERO)
C     
  910 IF (LAUE-12) 1010,920,1100
  920 DO 960  I=1,HMAX
      II=I+1
      IF (II.GT.HMAX) GO TO 960
      DO 950  J=II,HMAX
      DO 950  K=2,LMAX
      NH=HH(I)
      NK=HK(J)
      NL=HL(K)
	IF (NREF+1.GT.NRMAX) GOTO 210
C Test for extinction
	IF (.not. TSTEXTN(NH,NK,NL,IEXT)) THEN
C  Compute D-space
	  RECD=2.0*SQRT((HH(I)**2+HK(J)**2+HH(I)*HK(J))*PRCP(1)+HL(K)**2
	1	*PRCP(3))
C if in range, add to the arrays
          IF (RECD .GT. RDMAX) THEN
	    D(NREF+1) = 1.0 / RECD
	    IH(NREF+1) = NH
	    IK(NREF+1) = NK
	    IL(NREF+1) = NL
	    IF (D(NREF+1) .ge. DMIN) NREF = NREF + 1
	  ENDIF
	ENDIF
  950 CONTINUE
  960 CONTINUE
      GO TO 1100
C     
C     LAUE GROUP -3M ROTATED
C     
C     GENERATE MILLER INDICES (HKL)  (H GREATER OR EQUAL K, K NOT EQUAL
C     ZERO, ALL +L)
C     
  970 IF (LAUE.NE.15) GO TO 1100
      DO 1000  J=2,HMAX
      DO 1000  I=J,HMAX
      DO 1000  K=2,LMAX
      NH=-HH(I)
      NK=-HK(J)
      NL=HL(K)
	IF (NREF+1.GT.NRMAX) GOTO 210
C Test for extinction
	IF (.not. TSTEXTN(NH,NK,NL,IEXT)) THEN
C  Compute D-space
	  RECD=2.0*SQRT((HH(I)**2+HK(J)**2-HH(I)*(-HK(J)))*PRCP(1)+HL(K)**2
	1	*PRCP(3))
C if in range, add to the arrays
          IF (RECD .GT. RDMAX) THEN
	    D(NREF+1) = 1.0 / RECD
	    IH(NREF+1) = NH
	    IK(NREF+1) = NK
	    IL(NREF+1) = NL
	    IF (D(NREF+1) .ge. DMIN) NREF = NREF + 1
	  ENDIF
	ENDIF
 1000 CONTINUE
      GO TO 1100
C     
C     LAUE GROUP -3
C     
C     GENERATE MILLER INDICES (HKL) (H LESS THAN K, ALL +L)
C     
 1010 DO 1060  I=1,HMAX
      II=I+1
      IF (II.GT.HMAX) GO TO 1060
      DO 1050  J=II,HMAX
      DO 1050  K=1,LMAX
      IF (HL(K).NE.0) GO TO 1020
      IF (HH(I).EQ.0) GO TO 1050
 1020 NH=HH(I)
      NK=HK(J)
      NL=HL(K)
	IF (NREF+1.GT.NRMAX) GOTO 210
C Test for extinction
	IF (.not. TSTEXTN(NH,NK,NL,IEXT)) THEN
C  Compute D-space
	  RECD=2.0*SQRT((HH(I)**2+HK(J)**2+HH(I)*HK(J))*PRCP(1)+HL(K)**2
	1	*PRCP(3))
C if in range, add to the arrays
          IF (RECD .GT. RDMAX) THEN
	    D(NREF+1) = 1.0 / RECD
	    IH(NREF+1) = NH
	    IK(NREF+1) = NK
	    IL(NREF+1) = NL
	    IF (D(NREF+1) .ge. DMIN) NREF = NREF + 1
	  ENDIF
	ENDIF
 1050 CONTINUE
 1060 CONTINUE
C     
C     GENERATE MILLER INDICES (HKL) (H AND K NOT ZERO, ALL -L)
C     
      DO 1090  I=2,HMAX
      DO 1090  J=2,HMAX
      DO 1090  K=2,LMAX
      NH=HH(I)
      NK=HK(J)
      NL=-HL(K)
	IF (NREF+1.GT.NRMAX) GOTO 210
C Test for extinction
	IF (.not. TSTEXTN(NH,NK,NL,IEXT)) THEN
C  Compute D-space
	  RECD=2.0*SQRT((HH(I)**2+HK(J)**2+HH(I)*HK(J))*PRCP(1)+HL(K)**2
	1	*PRCP(3))
C if in range, add to the arrays
          IF (RECD .GT. RDMAX) THEN
	    D(NREF+1) = 1.0 / RECD
	    IH(NREF+1) = NH
	    IK(NREF+1) = NK
	    IL(NREF+1) = NL
	    IF (D(NREF+1) .ge. DMIN) NREF = NREF + 1
	  ENDIF
	ENDIF
 1090 CONTINUE
C----------------------------------------------------------------------
 1100 RETURN
210	WRITE(lo,'(A)') 
	1	' There are too many reflections -- increasing DMIN'
	  DMIN=DMIN*1.1
	GOTO 200
      END
C----------------------------------------------------------------------
	SUBROUTINE CELLPARM(lo,li,LAUE,A1,B1,C1,AL,BE,GA)
	INTEGER lo,li,LAUE
	REAL A1,B1,C1,AL,BE,GA
C----------------------------------------------------------------------
	LOGICAL GET_NUMBER,GET_FNUMBER,flag

	LAUE = 0

	FLAG = GET_FNUMBER('A',A1)
	IF (FLAG) RETURN

	B1 = A1
	C1 = A1
	AL = 90.
	BE = 90.
	GA = 90.

	FLAG = GET_FNUMBER('B',B1)
	FLAG = GET_FNUMBER('C',C1)
	FLAG = GET_FNUMBER('ALPHA',AL)
	FLAG = GET_FNUMBER('BETA',BE)
	FLAG = GET_FNUMBER('GAMMA',GA)

C now classify the lattice symmetry
	IF (AL .EQ. 90. .and. BE .EQ. 90. .and. GA .EQ. 90.) THEN
	  IF (a1 .eq. a2 .and. a1 .eq. a3) THEN
	    LAUE = 8
	  ELSEIF (a1 .eq. a2) THEN
	    LAUE = 5
	  ELSE
	    LAUE = 4
	  ENDIF
	ELSEIF (AL .EQ. 90. .and. BE .EQ. 90. .and. GA .EQ. 120.
	1	.and. a1 .eq. a2) THEN
	  LAUE = 13
	ELSEIF (AL .EQ. BE .and. BE .EQ. GA 
	1	.and. a1 .eq. a2 .and. a1 .eq. a3) THEN
	  LAUE = 9
	ELSEIF (AL .EQ. 90. .and. GA .EQ. 90.) THEN
	  LAUE = 2
	ELSEIF (AL .EQ. 90. .and. BE .EQ. 90.) THEN
	  LAUE = 3
	ELSE
	  LAUE = 1
	ENDIF
	RETURN
	END
