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
C     DMIN IS THE SMALLEST D-SPACING TO BE CALCULATED, AND 
C     XLAMDA IS THE WAVELENGTH. 
C     IEXT values determine Extinctions (set in GETEXTN, used in TSTEXTN)
C     
	INTEGER NRMAX
	PARAMETER (NRMAX=1000)
	INTEGER*2 IH(NRMAX),IK(NRMAX),IL(NRMAX),ISEQ(NRMAX),I2
	REAL	D(NRMAX)
	real RCP(6),PRCP(6)
	real  A,B,C,AL,BE,GA,RCV,VOL

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

      CHARACTER*1 NQUES,ANS
      CHARACTER*12 USER
	INTEGER NI,NO,NOUT,nsys,ncr, NREF, LINNUM, I, J
	REAL XLAMDA, XLIMIT
	REAL ST, SST, CT, TTH, XLP

      NI=5
      NO=6
	nout = 25

	CALL GETLAUE(LAUE,NSYS,no,ni)
	if (laue .le. 0 .or. nsys .le. 0) goto 9911

	CALL GETCELL(NSYS,no,ni,A,B,C,AL,BE,GA)
	IF (MIN(A,B,C) .LE. 0.0) goto 9911

	CALL CELCON(A,B,C,AL,BE,GA,RCV,RCP,VOL,PRCP)
	IF (VOL .LE. 0.0) goto 9911

   70 WRITE (NO,72)
   72 FORMAT (/' ENTER THE WAVELENGTH TO BE USED [1.540598]:  ',$)
      READ (NI,36,ERR=70,END=9911) NCR,XLAMDA
   36 FORMAT (Q,F10.0)
      IF (NCR.EQ.0) XLAMDA=1.540598
      IF (XLAMDA.LT.0.7 .OR. XLAMDA.GT.2.5) GO TO 70
   74 WRITE (NO,76)
   76 FORMAT (/' ENTER THE 2THETA LIMIT FOR D VALUES [90.00]:  ',$)
      READ (NI,36,ERR=74,END=9911) NCR,XLIMIT
      IF (NCR.EQ.0) XLIMIT=90.00
      IF (XLIMIT .GT. 179.00) GO TO 74
      DMIN = XLAMDA/(2.0*SIN(XLIMIT/114.59156))
      DMAX = XLAMDA/(2.0*SIN(0.5/114.59156))


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

	CALL DGENSUB(no,ni,LAUE,DMIN,DMAX,NREF,IEXT,
	1	A,B,C,rcp,prcp,
	1	NRMAX,D,IH,IK,IL)

	IF (LAUE .EQ. 7 .or. LAUE .eq. 8) THEN
	  WRITE (NO,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 (NO,1001) NREF,'TETRAGONAL'
	ELSEIF (LAUE .EQ. 4) THEN
	  WRITE (NO,1001) NREF, 'ORTHORHOMBIC'
	ELSEIF (LAUE .ge. 11 .and. LAUE .le. 15) THEN
	  WRITE (NO,1001) NREF, 'HEXAGONAL'
	ELSEIF (LAUE .EQ. 2 .or. LAUE .eq. 3) THEN
	  WRITE (NO,1001) NREF, 'MONOCLINIC'
	ELSEIF (LAUE .EQ. 1) THEN
	  WRITE (NO,1001) NREF, 'TRICLINIC'
	ELSEIF (LAUE .EQ. 9 .or. LAUE .eq. 10) THEN
	  WRITE (NO,1001) NREF, 'RHOMBOHEDRAL'
	ELSE
	  WRITE (NO,'(A)') ' Illegal value for LAUE'
	ENDIF

      VOLREAL=1.0/VOL
      WRITE (NO,1007) A,B,C,AL,BE,GA,XLAMDA,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.4)

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

	WRITE(NO,1145)
1145	FORMAT('$Type the generated d-spaces? (Y/[N]): ')
	READ(NI,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)
	  ST=XLAMDA/(2.0D0*D(J))
	  SST=ST**2
	  CT=DSQRT(1.0D0-SST)
	  TTH=114.59156D0*ATAN(ST/CT)
	  XLP=SST*CT/(1.0D0+(1.0D0-2.0D0*SST)**2)
	  IF (LINNUM .GE. 35) THEN 
	    WRITE (NO,1415)
1415	    FORMAT ('1      H      K      L       D-SPACING    SINE SQ.'
	1	,'     2THETA      LP')
	    LINNUM=0
	  ENDIF

	  WRITE (NO,120) IH(I),IK(I),IL(I),D(J),SST,TTH,XLP
120	  FORMAT (1X,3I7,F16.6,F12.5,F11.2,F12.4)
	  LINNUM=LINNUM+1
	ENDDO
	ENDIF

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

	IF(ANS.EQ.'Y'.OR.ANS.EQ.'y') THEN
          OPEN(UNIT=nout,FILE='DGEN.OUT',FORM='FORMATTED',STATUS='NEW',
     1     RECL=80)
	IF (LAUE .EQ. 7 .or. LAUE .eq. 8) THEN
	  WRITE (nout,1001) NREF,'CUBIC'
	ELSEIF (LAUE .EQ. 5 .or. LAUE .eq. 6) THEN
	  WRITE (nout,1001) NREF,'TETRAGONAL'
	ELSEIF (LAUE .EQ. 4) THEN
	  WRITE (nout,1001) NREF, 'ORTHORHOMBIC'
	ELSEIF (LAUE .ge. 11 .and. LAUE .le. 15) THEN
	  WRITE (nout,1001) NREF, 'HEXAGONAL'
	ELSEIF (LAUE .EQ. 2 .or. LAUE .eq. 3) THEN
	  WRITE (nout,1001) NREF, 'MONOCLINIC'
	ELSEIF (LAUE .EQ. 1) THEN
	  WRITE (nout,1001) NREF, 'TRICLINIC'
	ELSEIF (LAUE .EQ. 9 .or. LAUE .eq. 10) THEN
	  WRITE (nout,1001) NREF, 'RHOMBOHEDRAL'
	ENDIF
	VOLREAL=1.0/VOL
	WRITE (nout,1007) A,B,C,AL,BE,GA,XLAMDA,DMIN,VOLREAL
152	FORMAT ('0Extinction codes are ',13I2,'  Laue Code is',i3/)
	WRITE (nout,152) (IEXT(I),I=1,13),laue
	LINNUM=60
	DO J=1,NREF
	  I = ISEQ(J)
	  ST=XLAMDA/(2.0D0*D(J))
	  SST=ST**2
	  CT=DSQRT(1.0D0-SST)
	  TTH=114.59156D0*ATAN(ST/CT)
	  XLP=SST*CT/(1.0D0+(1.0D0-2.0D0*SST)**2)
	  IF (LINNUM .GE. 55) THEN 
	    WRITE (nout,1415)
	    LINNUM=0
	  ENDIF

	  WRITE (nout,120) IH(I),IK(I),IL(I),D(J),SST,TTH,XLP
	  LINNUM=LINNUM+1
	ENDDO
	CLOSE(UNIT=nout,STATUS='PRINT/DELETE')
	ENDIF

      STOP 'END OF D-GEN PROGRAM'

9911  STOP 'ABNORMAL ENDING TO CODE'
      END
	SUBROUTINE SORT (NSEQ,ITEM,KDIM,NORDER)
C	SORT ARRAY ITEM AND STORE SEQUENCE IN ARRAY NSEQ.
C	DIMENSION OF ITEM AND NSEQ IS KDIM.
C	NORDER = +1 IMPLIES INCREASING ORDER.
C	NORDER = -1 IMPLIES DECREASING ORDER.
C
	INTEGER*2 NSEQ(*),IT,KDIM
	REAL ITEM(*),ITEMP
	N=KDIM-1
	DO 50 I=1,N
	  M=I
	  K=I+1
	  IF(NORDER.EQ.1) GO TO 30
	  DO 20 J=K,KDIM
	    IF(ITEM(M).LT.ITEM(J)) M=J
20	  CONTINUE
	  GO TO 45

30	  DO 40 J=K,KDIM
	    IF(ITEM(M).GT.ITEM(J)) M=J
40	  CONTINUE

45	  IF (M .NE. I) THEN
	    ITEMP=ITEM(M)
	    ITEM(M)=ITEM(I)
	    ITEM(I)=ITEMP
	    IT=NSEQ(M)
	    NSEQ(M)=NSEQ(I)
	    NSEQ(I)=IT
	  ENDIF
50	CONTINUE
	RETURN
	END
