C----------------------------------------------------------------------
	Subroutine getindex(dindex,aindex,LUNC)
	integer dindex,aindex,LUNC
C Open the requested user index file. If the file cannot be opened
C find the next lowest file. Set aindex as the actual file opened or
C set aindex to 101 if no file is opened.
C  DINDEX is the threshold in percent or -1 for the 3strongest file
C----------------------------------------------------------------------
	CHARACTER*50 file
	INTEGER  I
C assume failure
	AINDEX = 101
	IF (dindex .le. -1) THEN
	  OPEN(unit=LUNC,file='PS_bitmap:DSP3S.BITMAP',status='old',
	1	shared,readonly,access='direct',form='unformatted',
	2	err=3)
	  AINDEX = -1
	  RETURN
	ENDIF
C assume success
	aindex = dindex
1	write(file,'(a,i2.2,a)') 'PS_bitmap:DSP',aindex,'.BITMAP'
	OPEN(unit=LUNC,file=file,status='old',shared,
	1	readonly,access='direct',form='unformatted',
	2	err=10)
	RETURN
C round down to next lowest 10%
10	i = 10 * int(aindex/10.)
	if (i .eq. aindex) i = Max(0,i-10)
	if (i .eq. aindex) THEN
	  aindex = 101
	  RETURN
	ENDIF
	aindex = i
	goto 1
3	RETURN
	END
C---------------------------------------------------------------------------
	LOGICAL FUNCTION DTOTT(DSPACE,TWOTHETA,WAVE)
C--
C  Converts a D space to a 2Theta value
C	the function returns .TRUE. if the conversion was successful
C--
C	DSPACE    is the D value 
C	TWOTHETA  is the (returned) TWOTHETA value
C	WAVE	  is the wavelength for 2Theta's
C
	REAL DSPACE,TWOTHETA,WAVE

	REAL X
	
	DTOTT = .false.

	IF (DSPACE .LE. 0.0 .or. WAVE .le. 0.0) RETURN

	X = 0.5 * WAVE / DSPACE
C	Is the position within the Ewald sphere?
	IF (X .LE. 1.0) THEN
	  twotheta = 360.0 * ASIN(X) / 3.14159265
	  DTOTT = .true.
	ENDIF
	RETURN
	END
C---------------------------------------------------------------------------
	LOGICAL FUNCTION TTTOD(DSPACE,TWOTHETA,WAVE)
C--
C  Converts a 2Theta value to a D space
C	the function returns .TRUE. if the conversion was successful
C--
C	DSPACE    is the (returned ) D value 
C	TWOTHETA  is the TWOTHETA value
C	WAVE	  is the wavelength for 2Theta's
C
	REAL DSPACE,TWOTHETA,WAVE

	REAL X

	TTTOD = .false.

	IF (WAVE .le. 0.0) RETURN

	IF (twotheta .le. 0.0 .or. twotheta .gt. 180.) RETURN
C convert to theta (radians) 
	x = 3.14159265 * twotheta / 360.0 
	X = SIN(X)
	IF (X .le. 0) RETURN
	DSPACE = 0.5 * WAVE / X
	TTTOD = .true.
	RETURN
	END
C---------------------------------------------------------------------------
	LOGICAL FUNCTION PDFTOD(PDF2T,DSPACE)
C--
C  Converts a PDF (D/I) encoding to a D space
C	the function returns .TRUE. if the conversion was successful
C--
C	DSPACE    is the (returned) D value 
C	PDF2T     is the encoded 2Theta value

C
	REAL DSPACE
	INTEGER*2 PDF2T
	REAL X

C     WAVE2	  is the wavelength/2 (1.540598A)
	REAL WAVE2
	PARAMETER (WAVE2 = 0.770299)

	PDFTOD = .false.

C    convert PDF2T to theta (in radians) at 1.540598A 
	X = 0.25*(0.0001*PDF2T + 3.14159265)
	IF (X .le. 0) RETURN
C    convert to D space 
	DSPACE = WAVE2/SIN(X)
	PDFTOD = .true.
	RETURN
	END
C---------------------------------------------------------------------------
	LOGICAL FUNCTION DTOPDF(PDF2T,DSPACE)
C--
C  Converts a D space to a PDF (D/I) encoded value
C	the function returns .TRUE. if the conversion was successful
C--
C	DSPACE    is the D value 
C	PDF2T     is the (returned) encoded 2Theta value

	REAL DSPACE
	INTEGER*2 PDF2T

C     WAVE2	  is the wavelength/2 (1.540598A)
	REAL WAVE2
	PARAMETER (WAVE2 = 0.770299)

	REAL X

	DTOPDF = .false.

	IF (DSPACE .LE. 0.0) RETURN
C    convert D to theta (in radians) at 1.540598A 
	X = WAVE2 / DSPACE
C	Is the position within the Ewald sphere?
	IF (X .LE. 1.0) THEN
	  X = ASIN(X)
	  X = 10000.*(4.0*X - 3.14159265)
C   now round to nearest integer 
	  IF (X .lt. 0) THEN
	    PDF2T = X - 0.5
	  ELSE
	    PDF2T = X + 0.5
	  ENDIF
	  DTOPDF = .true.
	ENDIF
	RETURN
	END
C---------------------------------------------------------------------------
	LOGICAL FUNCTION DTOIND(DSPACE,INDEX,NREGIONS)
C--
C  Converts a D space to an INDEX number for the ``D-space Indices''
C	the function returns .TRUE. if the conversion was successful
C
C  this approach is a little more complex than desired, but it is completely
C  compatible with the computation in Harlow's PDFENCODE program (subroutine
C  ENDSPINT) 12/87 BHT
C--
C	DSPACE    is the D value 
C	INDEX     is the (returned) record number 1<=INDEX<=NREGIONS
C	NREGIONS  is the number of records on the files

	REAL DSPACE
	INTEGER*4 NREGIONS,INDEX

C     WAVE2	  is the wavelength/2 (1.540598A)
	REAL WAVE2
	PARAMETER (WAVE2 = 0.770299)
	INTEGER*2 PDF2T

	REAL X

	DTOIND = .false.
	INDEX = NREGIONS

	IF (DSPACE .LE. 0.0) RETURN
C    convert D to theta (in radians) at 1.540598A 
	X = WAVE2 / DSPACE
C	Is the position within the Ewald sphere?
	IF (X .GT. 1.0) THEN
	  INDEX = NREGIONS
	ELSE
	  X = ASIN(X)
	  X = 10000.*(4.0*X - 3.14159265)
C   now round to nearest integer 
	  IF (X .lt. 0) THEN
	    PDF2T = X - 0.5
	  ELSE
	    PDF2T = X + 0.5
	  ENDIF
C subtract offset
	  INDEX = PDF2T/100. + 314.0 + 1.
	  IF (INDEX .GE. 1 .and. INDEX .LE. NREGIONS) DTOIND = .true.
	ENDIF
	INDEX = MAX(1,INDEX)
	INDEX = MIN(NREGIONS,INDEX)
	RETURN
	END
C---------------------------------------------------------------------------
	LOGICAL FUNCTION INDTOD(INDEX1,INDEX2,DMAX,DMIN,NREGIONS)
C--
C  Converts a range of INDEX numbers to the MIN and MAX d-space covered 
C  in the ``D-space Indices.''
C  The d-space range minimum d-space covered in INDEX1 to the maximum
C  in INDEX2 (or reverse, if INDEX2 < INDEX1). Values outside the legal 
C  range (1...NREGIONS) are set to the appropriate value
C	The function returns .TRUE. if the INDEX values are legal
C--
C  this approach is a little more complex than desired, but it is completely
C  compatible with the computation in Harlow's PDFENCODE program (subroutine
C  ENDSPINT) 12/87 BHT
C	INDEX1    is a record number
C	INDEX2    is a record number 
C	DMAX      is the maximum D value 
C	DMIN      is the minimum D value 

	REAL DMAX,DMIN
	INTEGER*4 INDEX1,INDEX2,NREGIONS

C     WAVE2	  is the wavelength/2 (1.540598A)
	REAL WAVE2
	PARAMETER (WAVE2 = 0.770299)
	INTEGER*2 PDF2T
	INTEGER*4 Ilow,Ihigh

	REAL X,X1

	INDTOD = .true.
	Ilow = MIN(INDEX1,INDEX2)
	Ihigh = MAX(INDEX1,INDEX2)
	IF (Ilow .lt. 1) then 
	  ILOW = 1
	  INDTOD = .false.
	ENDIF
	IF (Ihigh .gt. NREGIONS) then 
	  Ihigh = NREGIONS
	  INDTOD = .false.
	ENDIF

	x = 100.*(Ilow-1.0049999999) - 31400
	x1 = 100.*(Ihigh-0.0049999999) - 31400
	x = (X/10000. + 3.14159265)/4.0
	x1 = (X1/10000. + 3.14159265)/4.0
	DMAX = WAVE2 / sin(x)
	DMIN = WAVE2 / sin(x1)
	RETURN
	END
C---------------------------------------------------------------------------
	LOGICAL FUNCTION QTOTT(Q,TWOTHETA,WAVE)
C--
C  Converts a Q value to a 2Theta value
C	the function returns .TRUE. if the conversion was successful
C--
C	Q	is the Q value 
C	TWOTHETA  is the (returned) TWOTHETA value
C	WAVE	  is the wavelength for 2Theta's
C
	REAL Q,TWOTHETA,WAVE

	REAL X
	
	QTOTT = .false.

	IF (WAVE .le. 0.0) RETURN

C    Is the position within the Ewald sphere?
	X = WAVE * Q / (4.* 3.14159265)
	IF (X .lt. 0.0 .or. X .gt. 1.0) RETURN

	twotheta = 360.0 * ASIN(X) / 3.14159265
	QTOTT = .true.

	RETURN
	END
C---------------------------------------------------------------------------
	LOGICAL FUNCTION TTTOQ(Q,TWOTHETA,WAVE)
C--
C  Converts a 2Theta value to a Q space
C	the function returns .TRUE. if the conversion was successful
C--
C	Q is the (returned) Q value 
C	TWOTHETA  is the TWOTHETA value
C	WAVE	  is the wavelength for 2Theta's
C
	REAL Q,TWOTHETA,WAVE

	REAL X

	TTTOQ = .false.

	IF (WAVE .le. 0.0) RETURN

	IF (twotheta .le. 0.0 .or. twotheta .gt. 180.) RETURN
C convert to theta (radians)
	x = 3.14159265 * twotheta / 360.0
	Q = SIN(X) * 4.* 3.14159265 / wave
	TTTOQ = .true.
	RETURN
	END
C---------------------------------------------------------------------------
	LOGICAL FUNCTION PDFTOCU(PDF2T,TWOTHETA)
C--
C  Converts a PDF (D/I) encoding to a 2Theta angle at Cu Ka1 (1.540598A)
C	the function returns .TRUE. if the conversion was successful
C--
C	TWOTHETA  is the (returned) TWOTHETA value
C	PDF2T     is the encoded 2Theta value
C
	REAL TWOTHETA
	INTEGER*2 PDF2T

	REAL X

	PDFTOCU = .false.

C    convert PDF2T to 2theta (in degrees) at 1.540598A 
	X = 180. * 0.5*(PDF2T/10000. + 3.14159265) / 3.14159265
	IF (X .lt. 0. .or. X .gt. 180.) RETURN
	TWOTHETA = x
	PDFTOCU = .true.
	RETURN
	END
C---------------------------------------------------------------------------
	LOGICAL FUNCTION CUTOPDF(PDF2T,TWOTHETA)
C--
C  Converts a 2Theta angle at Cu Ka1 (1.540598A) to a PDF (D/I) encoded value
C	the function returns .TRUE. if the conversion was successful
C--
C	TWOTHETA  is the TWOTHETA value
C	PDF2T     is the (returned) encoded 2Theta value
C
	REAL TWOTHETA
	INTEGER*2 PDF2T

	REAL X

	CUTOPDF = .false.

	IF (TWOTHETA .lt. 0. .or. TWOTHETA .gt. 180.) RETURN
C    convert 2theta angle to PDF value 
	X = (2. * 3.14159265 * TWOTHETA / 180. - 3.14159265) * 10000. 
	IF (X .lt. 0) THEN
	  PDF2T = X - 0.5
	ELSE
	  PDF2T = X + 0.5
	ENDIF
	CUTOPDF = .true.
	RETURN
	END
C---------------------------------------------------------------------------
	LOGICAL FUNCTION PDFTOIND(PDF2T,INDEX,NREGIONS)
C--
C  Converts a PDF (D/I) encoded value to an INDEX number for 
C  the ``D-space Indices''
C	the function returns .TRUE. if the conversion was successful
C
C  this approach is a little more complex than desired, but it is completely
C  compatible with the computation in Harlow's PDFENCODE program (subroutine
C  ENDSPINT) 12/87 BHT
C--
C	angle     is the input two-theta value (in degrees)
C	INDEX     is the (returned) record number 1<=INDEX<=NREGIONS
C	NREGIONS  is the number of records on the files

	INTEGER*2 PDF2T
	INTEGER*4 NREGIONS,INDEX
	INTEGER*2 I

	PDFTOIND = .false.

	I = PDF2T/100. + 314.0 + 1.

	INDEX = 1
	IF (I .lt. 1) RETURN
	INDEX = NREGIONS
	IF (I .gt. NREGIONS) return

	INDEX = I
	PDFTOIND = .true.
	RETURN
	END
C---------------------------------------------------------------------------
	SUBROUTINE DATCONV(XIN,INTYP,npin,wavein,
	1		   XOUT,OUTYP,npout,WAVEout,*)
C---------------------------------------------------------------------------
C converts array XIN to XOUT 
C---------------------------------------------------------------------------
C Input:
C    XIN	array of input X values
C    INTYP	data type of XIN
C    npin	number of points in array XIN
C    wavein	wavelength (in Angstrom) for INTYP == 2theta
C    INTYP	data type for XOUT
C    waveout	wavelength (in Angstrom) for OUTYP == 2theta
C Output: 
C    XOUT	array of output X values
C    npout 	number of points in array XOUT
C Note that XIN and XOUT may be the same array 
C    npo is a temporary copy of NPOUT in case NPOUT and NPIN are the same
C	variable
C Units for INTYP and OUTYP are defined in FLABEL.COMMON
C---------------------------------------------------------------------------
	REAL XIN(*),XOUT(*),WAVEin,waveout
	INTEGER*2 INTYP,OUTYP
	INTEGER   NPOUT,NPIN
	include 'flabel.common'
	LOGICAL FLAG,DTOTT,TTTOD

	REAL X

	INTEGER  NPO,I,ITEMP
	NPO = 0 

D	type *,'DATA IN:'
D	type *,' type = ',intyp,' WAVE = ',wavein
D	type *,npin,' points from ',XIN(1),' to ',XIN(npin)
D	type *,'DATA OUT:'
D	type *,' type = ',outyp,' WAVE = ',waveout

C test for data types which cannot be converted at this time
	IF (INTYP .eq. XTYP_EV  .or. OUTYP .eq. XTYP_EV) goto  1
	IF (INTYP .eq. XTYP_TOF .or. OUTYP .eq. XTYP_TOF) goto  1
	IF (INTYP .eq. XTYP_R   .or. OUTYP .eq. XTYP_R) goto  1
	IF (INTYP .eq. XTYP_CH  .or. OUTYP .eq. XTYP_CH) goto  1
	IF (INTYP .eq. XTYP_LAM .or. OUTYP .eq. XTYP_LAM) goto  1
	IF (INTYP .eq. XTYP_TT .and. wavein .eq. 0) goto 2
	IF (OUTYP .eq. XTYP_TT .and. waveout .eq. 0) goto 2

C copy the in array to the out array; set a temporary variable with the data
C type
	npo = npin
	DO i=1,npin
	  XOUT(i) = XIN(i)
	ENDDO
	ITEMP = INTYP

C simple case -- no conversion
	IF (INTYP .eq. OUTYP .and. INTYP .ne. XTYP_TT) goto 100
	IF (INTYP .eq. OUTYP .and. wavein .eq. waveout) goto 100

C    interconversion of 2theta's between wavelengths
	IF (INTYP .eq. OUTYP .and. INTYP .eq. XTYP_TT) then
	  npo = 0
	  DO i=1,npin
	    FLAG = TTTOD(temp,XIN(i),wavein)
	    IF (FLAG) then
	      FLAG = DTOTT(temp,XOUT(npo+1),waveout)
	      IF (FLAG) npo = npo + 1
	    ENDIF
	  ENDDO
	  goto 100
	ENDIF

C for input type=Q convert to d-space
	IF (INTYP .eq. XTYP_Q) then
	  npo = 0
	  DO i=1,npin
	    IF (XIN(i) .gt. 0) then
	      npo = npo + 1
	      XOUT(npo) = 2*3.141592654/XIN(i)
	    ENDIF
	  ENDDO
	  ITEMP = XTYP_D
C for input type=2theta convert to d-space
	ELSEIF (INTYP .eq. XTYP_TT) then
	  npo = 0            
	  DO i=1,npin
	    IF (XIN(i) .gt. 0) then
	      FLAG = TTTOD(XOUT(npo+1),XIN(i),wavein)
	      IF (FLAG) npo = npo + 1
	    ENDIF
	  ENDDO
	  ITEMP = XTYP_D
	ENDIF
	IF (npo .le. 0) goto 3

C Now our XOUT array should be in d-space units
	IF (ITEMP .ne. XTYP_D) then
	  type *,'DATCONV ERROR -- UNKNOWN data type: ',ITEMP
	  RETURN 1
	ENDIF

	IF (OUTYP .eq. XTYP_D) then
	  goto 100
C convert from d-space to Q
	ELSEIF (OUTYP .eq. XTYP_Q) then
	  nptemp = npo
	  npo = 0
	  DO i=1,nptemp
	    IF (XOUT(i) .gt. 0) then
	      npo = npo + 1
	      XOUT(npo) = 2*3.141592654/XOUT(i)
	    ENDIF
	  ENDDO
	  ITEMP = XTYP_Q
C convert from d-space to 2Theta
	ELSEIF (OUTYP .eq. XTYP_TT) then
	  nptemp = npo
	  npo = 0
	  DO i=1,nptemp
	    IF (XIN(i) .gt. 0) then
	      FLAG = DTOTT(XOUT(npo+1),XIN(i),waveout)
	      IF (FLAG) npo = npo + 1
	    ENDIF
	  ENDDO
	  ITEMP = XTYP_TT
	ENDIF
	IF (npo .le. 0) goto 3
	IF (ITEMP .eq. OUTYP) goto 100
	type *,'DATCONV ERROR -- UNKNOWN data type: ',OUTYP
	RETURN 1
C----------------------------------------------------------------------
1	type *,'DATCONV ERROR -- Unable to convert types: ',INTYP,OUTYP
	return 1
2	type *,'DATCONV ERROR -- Invalid wavelength: ',
	1	INTYP,wavein,OUTYP,waveout
	return 1
3	type *,'DATCONV ERROR -- No valid points in input/output'
	return 1
C----------------------------------------------------------------------
100	NPOUT = NPO
D	type *,npout,' points from ',Xout(1),' to ',Xout(npout)
	return
	end
C---------------------------------------------------------------------------
	SUBROUTINE XYCONV(XIN,YIN,INTYP,npin,wavein,
	1		   XOUT,YOUT,OUTYP,npout,WAVEout,warn,*)
C---------------------------------------------------------------------------
C converts array XIN to XOUT and copies YIN to YOUT 
C---------------------------------------------------------------------------
C Input:
C    XIN	array of input X values
C    YIN	array of input Y values
C    INTYP	data type of XIN
C    npin	number of points in array XIN/YIN
C    wavein	wavelength (in Angstrom) for INTYP == 2theta
C    INTYP	data type for XOUT
C    waveout	wavelength (in Angstrom) for OUTYP == 2theta
C    warn	if TRUE generate a warning message
C Output: 
C    XOUT	array of output X values
C    YOUT	array of output Y values
C    npout 	number of points in array XOUT/YOUT
C Note that XIN and XOUT may be the same array 
C    npo is a temporary copy of NPOUT in case NPOUT and NPIN are the same
C	variable
C Units for INTYP and OUTYP are defined in FLABEL.COMMON
C---------------------------------------------------------------------------
	REAL XIN(*),XOUT(*),YIN(*),YOUT(*),WAVEin,waveout
	INTEGER*2 INTYP,OUTYP
	INTEGER   NPOUT,NPIN
	LOGICAL   WARN

	include 'flabel.common'
	LOGICAL FLAG,DTOTT,TTTOD

	INTEGER  NPO,I

D	type *,'DATA IN:'
D	type *,' type = ',intyp,' WAVE = ',wavein
D	type *,npin,' points from ',XIN(1),' to ',XIN(npin)
D	type *,'DATA OUT:'
D	type *,' type = ',outyp,' WAVE = ',waveout

	NPO = 0 
C test for data types which cannot be converted at this time
	IF (INTYP .eq. XTYP_EV  .or. OUTYP .eq. XTYP_EV) goto  1
	IF (INTYP .eq. XTYP_TOF .or. OUTYP .eq. XTYP_TOF) goto  1
	IF (INTYP .eq. XTYP_R   .or. OUTYP .eq. XTYP_R) goto  1
	IF (INTYP .eq. XTYP_CH  .or. OUTYP .eq. XTYP_CH) goto  1
	IF (INTYP .eq. XTYP_LAM .or. OUTYP .eq. XTYP_LAM) goto  1
	IF (INTYP .eq. XTYP_TT .and. wavein .eq. 0) goto 2
	IF (OUTYP .eq. XTYP_TT .and. waveout .eq. 0) goto 2

C copy the in array to the out array; set a temporary variable with the data
C type
	DO 200 i=1,npin

	  XOUT(npo+1) = XIN(i)
	  YOUT(npo+1) = YIN(i)

C simple case -- no conversion needed
	  IF (INTYP .eq. OUTYP .and. (INTYP .ne. XTYP_TT .or. 
	1	wavein .eq. waveout)) goto 150

C for input type=Q convert to d-space
	  IF (INTYP .eq. XTYP_Q) then
	    IF (XIN(i) .gt. 0) then
	      XOUT(npo+1) = 2*3.141592654/XIN(i)
	    ELSE
	      GOTO 200
	    ENDIF

C for input type=2theta convert to d-space
	  ELSEIF (INTYP .eq. XTYP_TT) then
	    IF (XIN(i) .gt. 0) then
	      FLAG = TTTOD(XOUT(npo+1),XIN(i),wavein)
	      IF (.not. FLAG) GOTO 200
	    ELSE
	      GOTO 200
	    ENDIF

C not d-space
	  ELSEIF (INTYP .ne. XTYP_D) then
	    goto 1
	  ENDIF

C Now our XOUT array element should be in d-space units

C no further convert needed
	  IF (OUTYP .eq. XTYP_D) then
	    goto 150

C convert from d-space to Q
	  ELSEIF (OUTYP .eq. XTYP_Q) then
	    IF (XOUT(npo+1) .gt. 0) then
	      XOUT(npo+1) = 2*3.141592654/XOUT(npo+1)
	      goto 150
	    ELSE
	      GOTO 200
	    ENDIF

C convert from d-space to 2Theta
	  ELSEIF (OUTYP .eq. XTYP_TT) then
	    IF (XOUT(npo+1) .gt. 0) then
	      FLAG = DTOTT(XOUT(npo+1),XOUT(npo+1),waveout)
	      IF (.not. FLAG) GOTO 200
	    ELSE
	      GOTO 200
	    ENDIF

	  ELSE
	    goto 1
	  ENDIF

150	  NPO = NPO + 1
200	CONTINUE
	NPOUT = NPO
D	type *,npout,' points from ',Xout(1),' to ',Xout(npout)
	return
C----------------------------------------------------------------------
1	IF (WARN) type 
	1	*,'XYCONV ERROR -- Unable to convert types: ',INTYP,OUTYP
	return 1
2	IF (WARN) type *,'XYCONV ERROR -- Invalid wavelength: ',
	1	INTYP,wavein,OUTYP,waveout
	return 1
3	IF (WARN) type 
	1	*,'XYCONV ERROR -- No valid points in input/output'
	return 1
C----------------------------------------------------------------------
	end
