	SUBROUTINE RIR(ISEL)
C*********************************************************************
C
C  SUBROUTINE TO DETERMINE THE REFERENCE INTENSITY RATIO (RIR).
C
C  IF CORUNDUM IS THE REFERENCE MATERIAL THEN RIR=I/Ic
C
C        WRITTEN BY
C  R. L. SNYDER NYS COLLEGE OF CERAMICS ALFERD UNIVERSITY
C       AND
C  C. R. HUBBARD NATIONAL BUREAU OF STANDARDS
C
C
C*********************************************************************
C  RELATIVE INTENSITY VALUES FOR THE INTERNAL STANDARD METHOD
C  IF ENTERED, WILL BE USED TO SCALE ALL LINES FROM EACH
C  STANDARD TO THE 100% LINE.  THIS ALLOWS COMPARISON
C  OF THE STANDARD INTENSITY VALUES WITH ALL ON THE SAME
C  SCALE.  IF RELATIVE INTENSITY VALUES WERE NOT INPUT
C  TO PROGRAM AUTO THEN THEY ARE
C  COMPUTED FOR EACH LINE OF EACH STANDARD.
C
C
C  THE FOLLOWING PARAMETER STATMENT DEFINES THE DIMENSIONS FOR
C  THE QUANTATATIVE ANALYSIS SUBROUTINES.  TO REDIMENSION
C  SIMPLY CHANGE THIS STATMENT HERE AND IN EACH SUBROUTINE.
C      NF = THE NUMBER OF DATA FILES ALLOWED FOR INPUT
C      NR = THE MAXIMUM NUMBER OF REFLECTIONS WHICH CAN BE PROCESSED
C      NP = THE MAXIMUM NUMBER OF PHASES ALLOWED
C      NT = THE MAXIMUM NUMBER OF DATA POINTS FOR ANY REFLECTION
C      NS = THE MAXIMUM NUMBER OF LINES FOR THE STANDARD PHASE
C      NU = THE MAXIMUM NUMBER OF LINES FOR ALL UNKNOWN PHASES
C
	PARAMETER (NF=30, NR=100, NP=10, NT=1250, NS=5, NU=20)
	COMMON/IO/ IRD,ITYP,IPRT,ISCRT,ICOBRG,LIMPRT,IDBUG
	COMMON /FILE/ LU,ITYPE,IEND
	COMMON /FILEC/ NAME
	COMMON /RECRD1/ IFLAG(10)
	COMMON /RECR1C/ TITLE,DATE,TIME,RNFLNM
	COMMON /RECRD2/ XLAMBD,TAKEOF,DIVANG,RADIUS,SOLLIN,
     1 SOLLOT,RECSLT,DEADTM,STABIL,TBASE
	COMMON /RECR2C/ UNIT,TARGET,LINE,FOCAL,DIVSLT,MONOCH,DETECT
	COMMON /RECRD6/ XWT(10),SIG(10)
	COMMON /RECR6C/ ELEMNT(10)
	COMMON /RECRD8/ FRELOV(2),SRELOV(2),WTPCPH(8)
	COMMON /RECR8C/ NAMEOV(2),JKLOV(2),ELMTPH(8)
	COMMON/RECRD7/BANG,PANG,EANG,DANG,RATE,MVCNT,TIMLIM,FIXERR,
     * LIMCNT,KODE,IOVLP,IPK,NPK,NBG,FREL,STDREL,STD,SIGSTD,
     * WTPRCT,PKTIME,BKGTIM,IOPTIM,ILOC,CONSTD
	COMMON /RECR7C/ NAM,JKL
	COMMON/SVDAT1/ AINT(NR,NF),VINT(NR,NF),NREFL,IFILE
	COMMON/LABELS/ PHASE(NR),HKL(NR)
	COMMON/SVDAT2/ FRELS(NR),SRELS(NR),WT(NF),KODES(NR),STDS(NR),
     * SSTDS(NR),IOVS(NR),IPKS(NR),RELOV(2,NR),
     * RELOVS(2,NR), WTP(8,NR), CONC(NF), XMURHO(NR)
	COMMON /SVDA2C/ MAN(2,NR),KHL(2,NR),ELM(8,NR)
	COMMON /RNFIL/ KIND(NR), NUNK, L,IWT
	COMMON /RNFILC/ AVL, BL, ALL, ESDHDG(2)
	COMMON/DATA/ PEAK,VPEAK,AREA,VAREA,XBKG,VBKG
	DIMENSION NOV(NU), SUM(NU), AVS(NU), AVSS(NU), NAVS(NU)
	CHARACTER*6 PHASE,HKL,BL,ALL
	CHARACTER*100 TITLE
	CHARACTER*6 DATE,TIME,MAN,JKL,NAMEOV,JKLOV,NAM,KHL,RNFLNM
	CHARACTER*20 NAME
	CHARACTER*2 ELMTPH, ELM, ELEMNT
	CHARACTER*4 UNIT
	CHARACTER*2 TARGET,LINE,FOCAL,DIVSLT,MONOCH,DETECT
	CHARACTER*7 AVL
	CHARACTER*4 LABEL(2)
	CHARACTER*9 ESDHDG
	DATA LABEL/' RIR','I/Io'/
C
C  ZERO THE AVERAGE     ARRAYS
C
	WAV = 0
	WAVS = 0
	K = 0
	DO 35 I=1,NU
	SUM(I) = 0.0
	AVS(I) = 0.0
	AVSS(I) = 0.0
35	NAVS(I) = 0
C
C  IF THIS IS THE INTENSITY RATIO METHOD JUMP TO 1000
C
	IF(IFLAG(3).GT.1)GO TO 1000
C******************************************************************
C
C  ALL CODE FROM HERE TO 1000 IS FOR INTERNAL STANDARD METHOD
C
C*******************************************************************
	WRITE(IPRT,1900)
1900	FORMAT(/' RIR Calibration for quantitative analysis by the',
     & ' internal standard method')
	WRITE(IPRT,2000) LABEL(1)
2000	FORMAT(/' ***** Line by line ', A4, ' analysis *****')
C
C-----FIND THE NEXT REFERENCE REFLECTION
C
38	NSTD=0
C
C  SET THE NOV FLAG IF ANY OVERLAPPED LINES ARE TO BE USED
C
	DO 45 I=1,NREFL
	IF(IOVS(I) - 1)45,41,44
C
C  CHECK IF THE OVERLAPPED PHASE BELONGS TO THE STANDARD
C
41	DO 42 J=1,NREFL
	IF(KODES(J).EQ.2)GO TO 43
42	CONTINUE
43	IF(MAN(1,I).NE.PHASE(J))GO TO 44
C
C  THIS LINE WILL BE USED IN THE COMPUTATION
C
	NOV(I) = J
	GO TO 45
44	NOV(I) = -1
45	CONTINUE
	JR = 0
	NRFL=0
	DO 300 I=1,NREFL
	IF(KODES(I).NE.2) GOTO 300
C
C  SAVE AN INDEX OF A REFERENCE LINE FOR OUTPUT
C
	IF(JR.GT.0)GO TO 46
	JR = I
C
C-----FOUND THE NEXT REFERENCE REFLECTION.  USE IT TO CALCULATE THE
C-----RIR FOR ALL COMBINATIONS OF THIS REFERENCE REFLECTION WITH UP TO
C-----NU REFLECTIONS FROM THE "UNKNOWN"
C
46	NSTD=NSTD+1
	IF(NSTD.LE.NS) GOTO 50
	WRITE(ITYP,40)NS
	IF(ITYP.NE.IPRT) WRITE(IPRT,40)NS
40	FORMAT(/' # of reference lines exceeds array limits NS=',I3,
     *' skip the extra lines')
	NSTD=NSTD-1
	GOTO 300
C
C  THIS IS A STANDARD LINE AND MAY NOT BE LABELED AS OVERLAPPED!
C
50	  IF(IOVS(I).EQ.0)GO TO 55
	WRITE(ITYP,2010)PHASE(I),HKL(I)
	IF(ITYP.NE.IPRT) WRITE(IPRT,2010)PHASE(I),HKL(I)
2010	FORMAT(/' The standard line ',2A6,' is labeled as overlapped'
     1,' the phase to be determined must be used as the principal',
     2' label'/' Change the data file and start over')
	CALL EXIT
C
C-----BEGIN LOOP OVER SAMPLE LINES
C
55	NRFL = 0
	DO 200 J=1,NREFL
	IF(KODES(J).NE.1) GOTO 200
C
C  SKIP ALL OVERLAPPED LINES EXCEPT THE CASE WHERE A SINGLE LINE
C  OVERLAPS A STANDARD LINE (THIS ONE CAN BE COMPUTED!)
C
	IF(NOV(J).LT.0)GO TO 200
C
C-----FOUND A SAMPLE REFLECTION.  NOW CALCULATE THE RIR  FOR ALL OBSERVA
C-----OF THE PAIR (I,J). (IE FOR ALL FILES)
C
	NRFL=NRFL+1
	IF(NRFL.LE.NU) GOTO 70
	WRITE(ITYP,60)NU
	IF(IPRT.NE.ITYP) WRITE(IPRT,60)NU
60	FORMAT(/' # of sample lines exceeds array limits NU=',I3,
     *' skip the extra lines')
	NRFL=NRFL-1
	GOTO 200
C
C  CALCULATE THE RIR FROM EACH MEASURED LINE
C
70	WRITE(IPRT,2020) LABEL(1), ESDHDG
2020	FORMAT(/' ','File Analyte  Line  Standard Line',6X,A4,1X,A9,
     1'  Root Mean  Number of  Degrees of'/
     1' ','  #   Phase  h k l   Name',
     1'   h k l  Constant ', A9,  ' square dev. measurments  freedom')
	S1 = 0.0
	S2 = 0.0
	S3 = 0.0
	NFREE = IFILE -1
	DO 100 K=1,IFILE
	RIOBS=AINT(J,K)/AINT(I,K)
	VRIOBS=RIOBS*RIOBS*(VINT(J,K)/(AINT(J,K)*AINT(J,K))+
     * VINT(I,K)/(AINT(I,K)*AINT(I,K)))
	RWTS = WT(K) / ((1.0 - WT(K))*CONC(K))
C
C  ADJUST RATIO IF LINE IS OVERLAPPED
C
	IF(NOV(J) .LE.0)GO TO 90
	RIOBS = RIOBS - RELOV(1,J) / FRELS(I)
	VRIOBS=RIOBS*RIOBS*(VINT(J,K)/(AINT(J,K)*AINT(J,K))+
     * VINT(I,K)/(AINT(I,K)*AINT(I,K)))
85	VRIOBS = VRIOBS + RIOBS*RIOBS*(RELOVS(1,J)*RELOVS(1,J)
     1 / (RELOV(1,J)*RELOV(1,J)) + SRELS(I)*SRELS(I) /
     2 (FRELS(I)*FRELS(I)))
90	COMTMP = RWTS * (FRELS(I)/FRELS(J)) * RIOBS
	VTMP = 1.0
	IF(IWT.NE.1) VTMP=COMTMP*COMTMP*(SRELS(I)*SRELS(I)/(FRELS(I)*
     *FRELS(I)) + SRELS(J)*SRELS(J)/(FRELS(J)*FRELS(J))+VRIOBS/
     *(RIOBS*RIOBS))
C
C  SAVE THE SUMS TO TAKE WEIGHTED AVERAGES AND RMS SIGMAS
C
C  WE USE 1/SIGMA**2 AS WEIGHTS (WHERE SIGMA IS FROM COUNTING
C  STATISTICS OR WE USE UNIT WEIGHTS)
C
	COM = COMTMP / VTMP
	S2 = S2 + COM
	AVS(NRFL) = AVS(NRFL) + COM
	COM = COMTMP*COM
	S1 = S1 + COM
	AVSS(NRFL) = AVSS(NRFL) + COM
	S3 = S3 + 1.0/VTMP
	SUM(NRFL) = SUM(NRFL) + 1.0/VTMP
	NAVS(NRFL) = NAVS(NRFL) + 1
	VV = SQRT(VTMP)
	WRITE(IPRT,2200)K,PHASE(J),HKL(J),PHASE(I),HKL(I),
     1COMTMP,VV
2200	FORMAT(1X,I3,2(2X,A6,1X,A6),F10.4,F10.4,1X,F10.4,I7)
100	CONTINUE
	WAV = S2 / S3
	WAVS = 0.0
	IF(IWT.NE.1) WAVS = SQRT(1.0 / S3)
	RMS = 0.0
	IF(IFILE.EQ.1)GO TO 200
	T = (S1/S3 - WAV*WAV)/(IFILE-1)
	IF(T.LE.0.0)GO TO 195
	RMS = SQRT(T)
	IF(NOV(J).GT.0)WRITE(IPRT,2300)
2300	FORMAT(' *************The following RIR is  computed using ',
     1'the input relative intensities************')
195	WRITE(IPRT,110)PHASE(J),HKL(J),AVL,HKL(I),WAV,WAVS,RMS,IFILE,NFREE
110	FORMAT(6X,A6,1X,A6,1X,A7,1X,A6,F10.4,F10.4,1X,F10.4,I7,I10)
200	CONTINUE
300	CONTINUE
C
C-----CHECK THAT AT LEAST 1 STANDARD AND AT LEAST 1 SAMPLE LINE WAS INPU
C
	IF(NSTD.NE.0 .AND. NRFL.NE.0) GOTO 400
	WRITE(ITYP,320) NSTD,NRFL
	IF(ITYP.NE.IPRT) WRITE(IPRT,320) NSTD,NRFL
320	FORMAT(/' **Input error** #reference lines=',I3,3X,
     * '# of sample lines=',I3)
	GOTO 999
C
C-----CALCULATE THE AVERAGE RIR FOR EACH ANALYTE LINE
C-----I.E. AVERAGE OVER ALL STANDARD LINES
C
400	WRITE(IPRT,1020) AVL
1020	FORMAT(///14X,A,  ' Over the standard lines used:'/)
	WRITE(IPRT,2020) LABEL(1), ESDHDG
C
C-----DETERMINE THE NUMBER OF DEGREES OF FREEDOM REMOVED BY DETERMINING
C-----RELATIVE INTENSITIES FROM THE SAME DATA SETS (ISEL.NE.3)
C-----EACH RELATIVE INTENSITY DETERMINED REDUCES THE NUMBER OF
C-----DEGREES OF FREEDOM BY ONE.
C     NSFIX IS FOR STANDARD LINES
C     NAFIX IS FOR ANALYTE LINES
C
	IF(ISEL.NE.3) THEN
	   NSFIX = NSTD-1
	   NAFIX = NRFL-1
	ELSE
	   NSFIX = 0
	   NAFIX = 0
	ENDIF
	NRFL = 0
	NTFREE = 0
	DO 420 I=1,NREFL
	IF(KODES(I).NE.1)GO TO 420
	IF(NOV(I).LT.0)GO TO 420
	NRFL = NRFL + 1
C
C  WRITE OMIT MESSAGE FOR OVERLAPPED LINES
C
	IF(NOV(I).GT.0)WRITE(IPRT,1025)
1025	FORMAT(' The following line will be excluded from the final',
     1' average due to overlap')
	WAV = AVS(NRFL) / SUM(NRFL)
	WAVS = 0.0
	IF(IWT.NE.1) WAVS = SQRT(1.0 / SUM(NRFL))
	RMS = 0.0
C
C-----NUMBER OF DEGREES OF FREEDOM IS EQUAL TO
C       # OF OBSERVATIONS - 1 FOR THE AVERAGE - NSFIX (# OF I-REL DETERM
C
	NFREE = NAVS(NRFL) - 1 - NSFIX
C-----IF THIS IS NOT THE I(REL)=100 LINE THEN DECREASE DEGREES OF FREEDO
C-----BY ONE SINCE THESE DATA ARE SCALED TO BE IN AGREEMENT WITH
C-----THE I(REL)=100 LINE ON AVERAGE.
	IF(FRELS(I).LT.99.999) NFREE = NFREE -1
C-----SUM UP TOTAL DEGREES OF FREEDOM FOR OVERALL AVERAGE
	NTFREE = NTFREE + NFREE
	IF(NFREE.LE.0) GOTO 410
	T = (AVSS(NRFL)/SUM(NRFL) - WAV*WAV)/NFREE
	IF(T.LE.0.0)GO TO 410
	RMS = SQRT(T)
410	WRITE(IPRT,1030)PHASE(I),HKL(I),PHASE(JR),ALL,WAV,WAVS,
     1 RMS,NAVS(NRFL) , NFREE
1030	FORMAT(6X,A6,1X,A6,2X,A6,1X,A6,F10.4,F10.4,1X,F10.4,I7,I10)
420	CONTINUE
C
C  LASTLY AVERAGE OVER THE LINES OF EACH PHASE
C
	WRITE(IPRT,1040) AVL
1040	FORMAT(///9X,' Overall ', A,  ' over all lines for each phase:'/)
	WRITE(IPRT,2020) LABEL(1), ESDHDG
	DO 450 I=1,NUNK
	S1 = 0.0
	S2 = 0.0
	S3 = 0.0
	NRFL = 0
	N = 0
	DO 440 J=1,NREFL
	IF(KODES(J).NE.1)GO TO 440
	IF(NOV(J).LT.0)GO TO 440
	NRFL = NRFL + 1
C
C  SKIP OVERLAPPED LINES IN THE FINAL AVERAGE
C
	IF(KIND(J).NE.I)GO TO 440
	IF(NOV(J).GT.0)GO TO 440
	S1 = S1 + AVSS(NRFL)
	S2 = S2 + AVS(NRFL)
	S3 = S3 + SUM(NRFL)
	N = N + NAVS(NRFL)
	K = J
440	CONTINUE
	WAV = S2 / S3
	WAVS = 0.0
	IF(IWT.NE.1) WAVS = SQRT(1.0 / S3)
	RMS = 0.0
C
C-----NUMBER OF DEGREES OF FREEDOM IS EQUAL TO
C       # OF OBSERVATIONS - 1 FOR THE AVERAGE - NSFIX - NAFIX
C         (NSFIX AND NAFIX CORRECT FOR NUMBER OF RELATIVE INTENSITIES
C          SCALED FROM THIS DATA)
C
	IF(NTFREE.LE.0) GOTO 445
	T = (S1/S3 - WAV*WAV)/NTFREE
	IF(T.LE.0.0)GO TO 445
	RMS = SQRT(T)
445	WRITE(IPRT,1030)PHASE(K),ALL,PHASE(JR),ALL,WAV,WAVS,RMS,N,NTFREE
450	CONTINUE
C
C  SAVE THE OVERALL AVERAGE FOR OUTPUT TO THE NEW RUNFILE
C
	DO 455 I=1,NREFL
	IF(KODES(I).NE.1)GO TO 455
	STDS(I) = WAV
	SSTDS(I) = WAVS
455	CONTINUE
	GO TO 999
C**********************************************************************
C
C  INTENSITY RATIO METHOD CALCULATION
C
C**********************************************************************
1000	WRITE(IPRT,1340)
1340	FORMAT(//' ',5X,'Calibration for the intensity ratio analysis'/)
C
C  CALCULATE STANDARD I  FOR EACH MEASURED LINE IN EACH FILE
C
	WRITE(IPRT,2000) LABEL(2)
	NRFL = 0
	LL = L
	IF(L.GT.0)GO TO 1050
	LL = NR
	PHASE(LL) = BL
	HKL(LL) = BL
1050	DO 1500 I=1,NREFL
	IF(I.EQ.L)GO TO 1500
	NRFL = NRFL + 1
	WRITE(IPRT,2020) LABEL(2), ESDHDG
	IF(NRFL.LE.NU)GO TO 1060
	WRITE(ITYP,60)NU
	IF(ITYP.NE.IPRT) WRITE(IPRT,60)NU
1060	DO 1400 J=1,IFILE
C
C  THE STANDARD INTENSITY IS SIMPLY THE I OF THE PURE SAMPLE
C
	RIOBS = AINT(I,J)
	VRIOBS = VINT(I,J)
	IF(L.EQ.0)GO TO 1350
C
C  A SPECIAL REFERENCE REFLECTION WAS USED
C
	RIOBS = RIOBS / AINT(L,J)
	VRIOBS=RIOBS*RIOBS*(VINT(I,J)/(AINT(I,J)*AINT(I,J))+
     * VINT(L,J)/(AINT(L,J)*AINT(L,J)))
C
C  NOTE AUTO STORES THE MU/RO OF THE STANDARD AS THE LAST
C  FIELD OF RECORD 7
C
1350	COMTMP = RIOBS
	VTMP = 1.0
	IF(IWT.NE.1) VTMP=COMTMP*COMTMP*(VRIOBS /
     *(RIOBS*RIOBS))
C
C  SAVE THE SUMS TO TAKE WEIGHTED AVERAGES AND RMS SIGMAS
C
C  WE USE 1/SIGMA**2 AS WEIGHTS (WHERE SIGMA IS FROM COUNTING STATS)
C
	COM = COMTMP / VTMP
	AVS(NRFL) = AVS(NRFL) + COM
	COM = COMTMP*COM
	AVSS(NRFL) = AVSS(NRFL) + COM
	SUM(NRFL) = SUM(NRFL) + 1.0/VTMP
	NAVS(NRFL) = NAVS(NRFL) + 1
	VV = SQRT(VTMP)
	WRITE(IPRT,3200)J,PHASE(I),HKL(I),PHASE(LL),HKL(LL),COMTMP,VV
3200	FORMAT(1X,I3,2(2X,A6,1X,A6),F10.1,F10.1,1X,F10.1,I7)
1400	CONTINUE
	WAV = AVS(NRFL) / SUM(NRFL)
	WAVS = 0.0
	IF(IWT.NE.1) WAVS = SQRT(1.0 / SUM(NRFL))
	STDS(I) = WAV
	SSTDS(I) = WAVS
	RMS = 0.0
	IF(NAVS(NRFL).EQ.1)GO TO 1500
	T = (AVSS(NRFL)/SUM(NRFL) - WAV*WAV)/(NAVS(NRFL)-1)
	IF(T.LE.0.0)GO TO 1450
	RMS = SQRT(T)
1450	WRITE(IPRT,3300)PHASE(I),HKL(I),AVL,HKL(LL),WAV,WAVS,RMS,IFILE
3300	FORMAT(6X,A6,1X,A6,1X,A7,1X,A6,F10.1,F10.1,1X,F10.1,I7)
1500	CONTINUE
999	RETURN
	END
