      INCLUDE 'RIETCOM.FOR'
      IOTTI =  5
      IOTTO =  6
      IODCF = 10
      IOPTN = 11
      IOBKG = 12
      IOALT =  7
      IOPLT = 13
      CALL INPTR
      CALL REWRIT
      END
C---  PROGRAM: RIETVELD  (RIET10)
C---  VERSION: 881001  (YYMMDD)
C---
      SUBROUTINE INPTR
C-
      INCLUDE 'RIETCOM.FOR'
      CHARACTER QDATE*9,QTIME*8,QSWTCH(-1:1)*4,QSPG(20)*1
      DIMENSION XRYZ(5),MLTT(MREFS)
      DIMENSION IMATX(0:MMATX)
      DIMENSION CELDIR(6),CELREC(3,3)
      COMMON XDEL(16),IDEL(16)
      EQUIVALENCE ( MLTT     , F2TOT    )
      DATA XRYZ(1), XRYZ(2), XRYZ(3), XRYZ(4), XRYZ(5)
     1   / 2.28962, 1.93597, 1.54060, 0.70926, 0.556363 /
      DATA QSWTCH(-1),QSWTCH(0),QSWTCH(1) /'On  ','Off ','On  '/
C-
C--  SAY HELLO
C-
1000  CALL DATE( QDATE )
      CALL TIME( QTIME )
      WRITE(IOTTO,1005) ' ',IVER,QDATE,QTIME(1:5)
1005  FORMAT(A1,
     1'*** RIETVELD analysis program',10X,'V ',I6,9X,'Run ',A9,2X,A5)
C-
C--  INITIALIZE PARAMETERS
C-
      DO 1010 I=0,MMATX
1010  IMATX(I) = 0
      IOV = 0
      DO 1020 ITMP=1,MPHAS
1020  IPREFS(ITMP) = 0
      NSCTOT = 0
C-
C--  SETUP I/O FILE NAMES
C-
      CALL SETIO
C-
C--  SAY HELLO AGAIN- THIS TIME FOR THE RECORD
C-
2000  WRITE(IOTTO,1005) '1',IVER,QDATE,QTIME(1:5)
C-
C-  GIVE PROGRAM LIMITATIONS
C-
      WRITE(IOTTO,2072) MREFS,MATMS,MMATX,MNPTS,MOVRL,MPHAS
2072  FORMAT(/,' Program limitations:',/,
     1     '  Total number of reflections   ',I5,
     1 T41,' Total number of atoms         ',I5,
     1 T82,' Number of variable parameters ',I5,/,
     1     '  Points in observed pattern    ',I5,
     1 T41,' Overlapping reflections       ',I5,
     1 T82,' Number of phases              ',I5)
C-
C--  READ THE PROBLEM TITLE AND DESCRIPTIVE INFORMATION
C-
      READ(IODCF,2010) QTITLE
2010  FORMAT(A80)
      WRITE(IOTTO,2080) QTITLE(1:77)
2080  FORMAT(/,' Problem ID:',/,2X,A77)
C-
      READ(IODCF,2020) JOBXN,NPHASE,NBCKGD,NEXCRG,NSCAT,KINTSF,KAMOR,
     1                 KWSGON,MAGNIN
2020  FORMAT(9I4)
C-
      IF(JOBXN.EQ.1) THEN
        WRITE(IOTTO,2040)
2040    FORMAT(2X,'Refinement of X-Ray data')
       ELSE IF(JOBXN.EQ.2) THEN
        WRITE(IOTTO,2050)
2050    FORMAT(2X,'Refinement of neutron data, nuclear intensities')
       ELSE IF(JOBXN.EQ.3) THEN
        WRITE(IOTTO,2060)
2060    FORMAT(2X,' Pattern calculation, X-Ray')
       ELSE IF(JOBXN.EQ.4) THEN
        WRITE(IOTTO,2070)
2070    FORMAT(2X,' Pattern calculation, Neutron')
       ELSE
        STOP '--INPTR: UNKNOWN JOB TYPE SPECIFIED'
      ENDIF
C-
      WRITE(IOTTO,2090) NPHASE,NEXCRG,NSCAT
2090  FORMAT(/,' Parameters for current refinement:',/,
     1     '  Number of phases in specimen  ',I5,
     1 T41,' Number of excluded regions    ',I5,
     1 T82,' Number of scattering sets     ',I5)
      IF(NPHASE.GT.MPHAS) THEN
        WRITE(IOTTO,2092)
2092    FORMAT(/,
     1' --THE NUMBER OF PHASES TO BE INPUT EXCEEDS THE CAPABILITY OF ',
     1'THIS PROGRAM.',/,
     1'   INCREASE THE VALUE OF *MPHAS* IN ALL PARAMETER STATEMENTS TO',
     1' SUIT YOUR NEEDS.')
        STOP '--INPTR: NUMBER OF INPUT PHASES EXCEEDS PROGRAM STORAGE'
      ENDIF
      IF(NPHASE.LT.0) STOP '--INPTR: NUMBER OF PHASES MUST BE > 0'
C-
      IF(NBCKGD.LT.1) THEN
        NBCKGD = 0
        WRITE(IOTTO,2100)
2100    FORMAT(/,' Background to be refined')
       ELSE IF(NBCKGD.EQ.1) THEN
        WRITE(IOTTO,2110)
2110    FORMAT(/,' Background data to be read from file')
       ELSE
        WRITE(IOTTO,2120)
2120    FORMAT(/,' Background correction by interpolation ',
     1  'between supplied data points')
      ENDIF
C-
      IF(KAMOR.EQ.1) THEN
        WRITE(IOTTO,2122)
2122    FORMAT(/,
     1  ' A Gaussian profile will be used to fit an amorphous line')
       ELSE
        KAMOR = 0
      ENDIF
C-
      IF(KWSGON.EQ.1) WRITE(IOTTO,2140)
2140  FORMAT(/,' A convolution based PSF will be used.')
C-
      IF(MAGNIN.EQ.1) WRITE(IOTTO,2142)
2142  FORMAT(/,' Magnetic scattering vectors entered for atoms.')
C-
C--  READ AND PRINT OUTPUT FLAGS
C-
2200  READ(IODCF,2210) IOUT,IPLT,IPCH,MAT,NXT,LST1,LST2,LST3,IPLT1,
     1  IPLT2,KCLRCW,KCONCW,KSCALZ
2210  FORMAT(13I4)
      IF(IOUT  .NE.0) IOUT   = 1
      IF(IPLT  .NE.0) IPLT   = 1
      IF(IPCH  .NE.0) IPCH   = 1
      IF(MAT   .NE.0) MAT    = 1
      IF(NXT   .NE.0) NXT    = 1
      IF(LST1  .NE.0) LST1   = 1
      IF(LST2  .NE.0) LST2   = 1
      IF(LST3  .NE.0) LST3   = 1
      IF(IPLT1 .NE.0) IPLT1  = 1
      IF(IPLT2 .NE.0) IPLT2  = 1
      IF(KCLRCW.NE.0) KCLRCW = 1
      IF(KCONCW.NE.0) KCONCW = 1
      IF(KSCALZ.NE.0) KSCALZ = 1
      IF(JOBXN .GT.2) THEN
        MAT    = 0
        KCLRCW = 0
      ENDIF
      IF(NPHASE  .EQ.1) LST3   = 0
      WRITE(IOTTO,2220) QSWTCH(IOUT),QSWTCH(IPLT),QSWTCH(IPCH)
2220  FORMAT(/,' Program flag settings:',/,
     1         '  Output intensities on last cycle ',A4,
     1     T41,' Generate line printer plot       ',A4,
     1     T82,' Output structure factors         ',A4)
      WRITE(IOTTO,2230) QSWTCH(MAT),QSWTCH(NXT),QSWTCH(LST1)
2230  FORMAT(
     1         '  Output correlation matrix        ',A4,
     1     T41,' Write new input file             ',A4,
     1     T82,' Print reflection list for phases ',A4)
      WRITE(IOTTO,2240) QSWTCH(LST2),QSWTCH(LST3),QSWTCH(IPLT1)
2240  FORMAT(
     1         '  Print corrected data             ',A4,
     1     T41,' Print merged reflection list     ',A4,
     1     T82,' Print symmetry operators         ',A4)
      WRITE(IOTTO,2250) QSWTCH(IPLT2),QSWTCH(KCLRCW),QSWTCH(KCONCW)
2250  FORMAT(
     1         '  Write data to plot file          ',A4,
     1     T41,' Clear unused code words          ',A4,
     1     T82,' Apply pre-selected constraints   ',A4)
      WRITE(IOTTO,2260) QSWTCH(KSCALZ-1)
2260  FORMAT(
     1         '  Scale normal equation matrix     ',A4)
C-
C-  READ AND PRINT AUTOPILOT FLAGS
C-
      READ(IODCF,2180) KAUTOV,KAUTOR
2180  FORMAT(20I4,/,20I4)
      DO 2182 I=1,20
2182  IF( KAUTOV(I).EQ.0 .OR. KAUTOR(I).EQ.0 ) GO TO 2184
      I = 21
2184  KAUTON = I - 1
      IF(KAUTON.EQ.0) STOP '--INPTR: NO AUTOPILOT DATA, ABORT'
      WRITE(IOTTO,2186) (I,KAUTOR(I),KAUTOV(I),I=1,KAUTON)
2186  FORMAT(/,' Autopilot refinement will follow:',/,
     1  (2X,'Round ',I2,':  Number of iterations ',I2,
     1  '  Number of variables ',I2))
      KMXVAR = KAUTOV(KAUTON)
      IF(KMXVAR.GT.MMATX) THEN
        WRITE(IOTTO,2922)
2922    FORMAT(/,
     1' --THE NUMBER OF VARIABLES EXCEEDS THE CAPABILITIES OF THIS ',
     1'PROGRAM.',/,
     1'   INCREASE THE VALUE OF *MMATX* ON ALL PARAMETER STATEMENTS TO',
     1' SUIT YOUR NEEDS.')
        STOP '--INPTR: NUMBER OF VARIABLES EXCEEDS PROGRAM STORAGE'
      ENDIF
C-
2400  READ(IODCF,2410) WAVELN,RATIO(2),WDTTHK,CTHM,THMIN,STEP,THMAX
2410  FORMAT(8F8.0)
      WRITE(IOTTO,2420) WAVELN,RATIO(2)
2420  FORMAT(/,' Wavelengths = ',2F10.6,'  Alpha2:Alpha1 ratio = ',F8.4)
      RATIO(1) = 1.0
      WRITE(IOTTO,2450) WDTTHK
2450  FORMAT(' Base of peak extends to = FWHM * ',F8.2)
      WRITE(IOTTO,2460) CTHM
2460  FORMAT(' Monochromator correction = ',F8.4)
C-
2500  READ(IODCF,2510) EPS,RELAX
2510  FORMAT(5F8.0)
      WRITE(IOTTO,2530) EPS
2530  FORMAT(/,
     1' Convergence achieved when parameter adjustmenst are <',F5.2,
     1' (EPS) * their ESDs')
      WRITE(IOTTO,2540) RELAX
2540  FORMAT(/,' Relaxation factors:',/,
     1         '  Atom xyz, iso temp and occup.   ',F5.2,2X,
     1         ' Atom anisotropic temp factor    ',F5.2,/,
     1         '  Phase related parameters        ',F5.2,2X,
     1         ' Global parameters               ',F5.2)
C-
C--  GET THE INSTRUMENT CALIBRATION DATA
C-
CZZZ  CALL GETWSG
C-
C-  SET IXRAY TO REFLECT RADIATION WAVELENGTH
      DO 2430 IXRAY=1,5
2430  IF( (WAVELN(1)+0.1) .GT. XRYZ(IXRAY) ) GO TO 3100
      IXRAY = 5
C-
C--  CHECK INPUT/OUTPUT PATTERN PARAMETERS AND READ IF NECESSARY
C-
C-  READ THE INPUT PATTERN IF NECESSARY
3100  IF(JOBXN.LT.3) THEN
CZZZ    CALL GETPAT
        GO TO 3002
      ENDIF
C-
C-  SET THE PARAMETERS FOR CALCULATED OUTPUT PATTERN
3200  IF(THMIN.EQ.0.0) THEN
        THMIN = 10.0
        WRITE(IOTTO,3210) THMIN
3210    FORMAT(2X,'2Theta minimum defaulting to ',F8.2,' degrees')
      ENDIF
      IF(THMAX.EQ.0.0) THEN
        THMAX = 90.0
        WRITE(IOTTO,3220) THMAX
3220    FORMAT(2X,'2Theta maximum defaulting to ',F8.2,' degrees')
      ENDIF
      IF(STEP.EQ.0.0) THEN
        STEP = 0.02
        WRITE(IOTTO,3230) STEP
3230    FORMAT(2X,'Angle increment defaulting to ',F8.2,' degrees')
      ENDIF
C-
C-  INSURE THAT OUTPUT PATTERN DOES NOT EXCEED STORAGE AREAS
      NPTS = NINT( (THMAX-THMIN) / STEP ) + 1
      WRITE(IOTTO,3240) THMIN,STEP,THMAX,NPTS
3240  FORMAT(/,' Calculated pattern:  THMIN ',F8.4,'  STEP ',F8.4,
     1 '  THMAX ',F8.4,'  NPTS ',I5)
      IF(NPTS.GT.MNPTS) THEN
        WRITE(IOTTO,3250) MNPTS
3250    FORMAT('  -LIMITING POINTS IN PATTERN TO ',I5)
        NPTS = MNPTS
        THMAX = THMIN + (MNPTS-1) * STEP
        WRITE(IOTTO,3240) THMIN,STEP,THMAX,NPTS
      ENDIF
C-
C--  READ AND PRINT TWO-THETA ERROR CORRECTION PARAMETERS
C-
3002  READ(IODCF,3010) PGLB(1),AGLB(1),PGLB(10),AGLB(10)
3010  FORMAT(6F8.0)
      WRITE(IOTTO,3020) PGLB(1),AGLB(1),PGLB(10),AGLB(10)
3020  FORMAT(/,' * GLOBAL PARAMETERS AND CODEWORDS:',//,
     1  ' Zeropoint = ',F8.4,2X,F8.2,6X,'Specimen displacement = ',
     1  F8.4,2X,F8.2)
C-
C--  HANDLING OF BACKGROUND FOR THE PATTERN
C-
C-  READ LINEAR INTERPOLATION POINTS FOR BACKGROUND IF REQUIRED
2600  IF(NBCKGD.EQ.0) THEN
        READ(IODCF,3030) (PGLB(I),I=2,7),(AGLB(I),I=2,7)
3030    FORMAT(6G12.5)
        WRITE(IOTTO,3040) (PGLB(I),I=2,7),(AGLB(I),I=2,7)
3040    FORMAT(/,' Polynomial background coefficients and codewords:',/,
     1  6X,'C0',11X,'C1',11X,'C2',11X,'C3',11X,'C4',11X,'C-1',/,
     1  1P,1X,6(1X,E12.5),/,0P,2X,6(1X,F8.2,4X))
       ELSE IF(NBCKGD.EQ.1) THEN
        OPEN(UNIT=IOBKG,FILE=QBKGFN,STATUS='OLD',ERR=2603)
        READ(IOBKG,2602) (BK(I),I=1,NPTS)
2602    FORMAT(10F7.0)
        CLOSE(IOBKG)
        GO TO 2700
2603    STOP'--THE BACKGROUND DATA FILE WAS NOT FOUND'
       ELSE IF(NBCKGD.GT.1) THEN
        NTOGET = NBCKGD
        NBCKGD = 0
2604    READ(IODCF,2606) TEMP1,TEMP2
2606    FORMAT(2F8.0)
        IF(TEMP1.NE.-100.0) THEN
          IF(NBCKGD.EQ.30) THEN
            STOP'--Max of 30 background pts for interpolation allowed.'
          ENDIF
          NBCKGD = NBCKGD + 1
          POS(NBCKGD) = TEMP1
          BCK(NBCKGD) = TEMP2
          IF(NBCKGD.EQ.NTOGET) GO TO 2607
          GO TO 2604
        ENDIF
2607    WRITE(IOTTO,2608) NBCKGD,(POS(I),BCK(I),I=1,NBCKGD)
2608    FORMAT(/,
     1  1X,'The ',I2,' points for background interpolation:',//,
     1  3X,'Position',4X,'Intensity',/,
     1  (5X,F8.4,4X,F9.4))
      ENDIF
C-
C-  DETERMINE BACKGROUND BY INTERPOLATION
      IF(NBCKGD.LT.2) GO TO 2700
      DIFB = POS(1) - THMIN
      IF(DIFB) 2610,2620,2620
2620  NBX = DIFB / STEP + 2.5
      DO 2630 I=1,NBX
2630  BK(I) = BCK(1)
      NXX = 1
      GO TO 2640
2610  NBX = (POS(2)-THMIN)/STEP + 1.5
      BK(1) = BCK(1) - DIFB / (POS(2)-POS(1)) * (BCK(2)-BCK(1))
      BSTEP = STEP / (POS(2)-POS(1)) * (BCK(2)-BCK(1))
      DO 2650 I=2,NBX
2650  BK(I) = BK(I-1) + BSTEP
      NXX = 2
2640  NBC = NBCKGD - 1
      IF(POS(NBCKGD).GE.THMAX) GO TO 2660
      POS(NBCKGD+1) = THMAX
      BCK(NBCKGD+1) = BCK(NBCKGD)
      NBC = NBC + 1
2660  DO 2670 J=NXX,NBC
      BSTEP = STEP * (BCK(J+1)-BCK(J)) / (POS(J+1)-POS(J))
      NINC = (POS(J+1)-POS(J))/STEP + 1.5
      N2X = MIN0( NPTS, NBX+NINC )
      DO 2680 I=NBX,N2X
2680  BK(I) = BK(I-1) + BSTEP
      NBX = N2X
      IF(NBX.EQ.NPTS) GO TO 2700
2670  CONTINUE
C-
C--  READ EXCLUDED REGIONS IN PATTERN
C-
2700  IF(NEXCRG.GT.0) THEN
        READ(IODCF,2606)  (ALOW(I),AHIGH(I),I=1,NEXCRG)
        WRITE(IOTTO,2710) (ALOW(I),AHIGH(I),I=1,NEXCRG)
2710    FORMAT(/,' Angular regions excluded from refinement:',/,
     1  (2X,'From ',F8.4,' to ',F8.4))
      ENDIF
C-
C--  READ SCATTERING SET ENTRIES
C-
2800  IF(NSCAT.LT.1) GO TO 2900
      IF( JOBXN.EQ.1 .OR. JOBXN.EQ.3 ) THEN
        WRITE(IOTTO,2802)
2802    FORMAT(/,' Form factors and analytical curve coefficients:')
       ELSE
        WRITE(IOTTO,2804)
2804    FORMAT(/,
     1  ' Scattering lengths and analytical curve coefficients:')
      ENDIF
C-
C-  REPEAT FOR INDICATED NUMBER OF SETS
      DO 2810 I=1,NSCAT
C-
      READ(IODCF,2820) QSCNAM(I),DFP(I),DFPP(I)
2820  FORMAT(A4,4X,2F8.0)
      IF( JOBXN.EQ.1 .OR. JOBXN.EQ.3 ) THEN
        WRITE(IOTTO,2892) QSCNAM(I),DFP(I),DFPP(I)
2892    FORMAT(/,3X,'Scatterer ',A4,':  DFP = ',F10.6,'  DFPP = ',F10.6)
       ELSE
        WRITE(IOTTO,2898) QSCNAM(I),DFP(I)
2898    FORMAT(/,3X,'Scatterer ',A4,':  Length = ',F10.6)
      ENDIF
      NSCOEF = 0
2830  READ(IODCF,2840) (SCATAC(J,I),J=1,9)
2840  FORMAT(9F8.0)
C-
C-  IF END-OF-LIST INDICATED THEN ASSUME SCATTERING CURVES WERE ENTERED
C-    AND CALL COEF TO GENERATE THE INTERNAL CURVE SET
      IF(SCATAC(1,I).EQ.-100.0) THEN
CZZZ    CALL COEF( I )
        DO 2850 J=1,NSCOEF
        TEMP = 0.0
        DO 2860 JJ=1,9,2
2860    TEMP = TEMP + SCATAC(JJ,I)*EXP( -SCATAC(JJ+1,I) * SINLAM(J)**2 )
        WRITE(IOTTO,2870) SINLAM(J),SCAT(J),TEMP
2870    FORMAT(
     1  5X,'Sinth/Lambda  ',F7.4,'  f = ',F7.4,'  f calc = ',F7.4)
2850    CONTINUE
        GO TO 2890
      ENDIF
C-
C-  IF THIRD COEFFICIENT IS 0.0 THEN ASSUME SCATTERING CURVE ENTERED
      IF(SCATAC(3,I).EQ.0.0) THEN
        NSCOEF = NSCOEF + 1
        SINLAM(NSCOEF) = SCATAC(1,I)
        SCAT(NSCOEF)   = SCATAC(2,I)
        IF(NSCOEF.LE.29) GO TO 2830
        STOP '--INPTR: TOO MANY SCATTERING CURVE ENTRIES'
      ENDIF
C-
2890  WRITE(IOTTO,2896) (SCATAC(J,I),J=1,9)
2896  FORMAT(5X,9F8.4)
C-
2810  CONTINUE
C-
C-  SET TOTAL NUMBER OF SCATTERING SETS TO THE NUMBER ENTERED
      NSCTOT = NSCAT
C-
C--  READ AND PRINT INTERNAL STANDARD DATA
C-
2900  IF(KINTSF.EQ.0) GO TO 3000
      READ(IODCF,2910) KINTSN,RINTSQ,KINTSA(1),RINTSA(1),
     1                 KINTSA(2),RINTSA(2),KINTSA(3),RINTSA(3)
2910  FORMAT(4(I8,F8.0))
C-
      IF(KINTSN.NE.0) THEN
        KINTSF = 1
        WRITE(IOTTO,2924) KINTSN,RINTSQ
2924    FORMAT(/,' Phase ',I2,' is an internal standard.  Weight',
     1  ' fraction ',F8.4)
        IF( KINTSN.LT.1 .OR. KINTSN.GT.NPHASE ) THEN
          WRITE(IOTTO,2926)
2926      FORMAT(' --INVALID PHASE:  NUMBER OUT OF RANGE.')
          KINTSF = 0
        ENDIF
      ENDIF
      IF(KINTSA(1).LT.0) THEN
        WRITE(IOTTO,2927) -KINTSA(1),RINTSA(1),-KINTSA(2),RINTSA(2),
     1  -KINTSA(3),RINTSA(3)
2927    FORMAT(/,
     1  ' Weight fraction data for scale factor corrections:',/,
     1  (2X,'Phase ',I2,': Weight fraction ',F8.4))
      ENDIF
      IF(KINTSA(1).GT.0) THEN
        WRITE(IOTTO,2928) KINTSA(1),RINTSA(1),KINTSA(2),RINTSA(2),
     1  KINTSA(3),RINTSA(3)
2928    FORMAT(/,' Data for scale factor corrections:',/,
     1  (2X,'Phase ',I2,': Correction factor ',F8.4))
      ENDIF
C-
C--  READ AND PRINT AMORPHOUS PROFILE PARAMETERS
C-
3000  IF(KAMOR.EQ.1) THEN
        READ(IODCF,3050) (PGLB(I),I=11,13), (AGLB(I),I=11,13)
3050    FORMAT(3F8.0,/,3F8.0)
        WRITE(IOTTO,3060) (PGLB(I),I=11,13), (AGLB(I),I=11,13)
3060    FORMAT(/,
     1  ' Amorphous line parameters: 2theta = ',F6.2,'  Peak hgt = ',
     1  F8.2,'  FWHM = ',F5.2,/,
     1  16X,'Code words:',10X,F6.2,13X,F8.2,8X,F6.2)
      ENDIF
C-
C-----------  REPEAT FOR EVERY PHASE IN THE PATTERN
C-
3500  DO 3510 IPHASE=1,NPHASE
C-
C-  OFFSET INTO ATOM PARAMETER LISTS
      IOF = 0
      IF(IPHASE.GT.1) IOF = NATM(2,IPHASE-1) + NATM(1,IPHASE-1) - 1
C-
      READ(IODCF,3520) QPHSNM(IPHASE)
3520  FORMAT(A32)
C-
      READ(IODCF,3560) (QSSYMB(I,IPHASE),I=1,20),QXTALS(IPHASE)
3560  FORMAT(20A1,A20)
C-
      READ(IODCF,3530) NATM(1,IPHASE)
3530  FORMAT(I4)
      NATOM = NATM(1,IPHASE)
      NATM(2,IPHASE) = IOF + 1
      WRITE(IOTTO,3540) IPHASE,QPHSNM(IPHASE),NATOM,
     1                  (QSSYMB(I,IPHASE),I=1,20),QXTALS(IPHASE)
3540  FORMAT('1Phase: ',I2,'  Name: ',A32,'  Number of atoms: ',I3,/,
     1     12X,'Space group: ',20A1,T53,'Crystal system: ',A20)
C-
C-  GENERATE SPACE GROUP ASSOCIATED PARAMETERS
      DO 3570 I=1,20
3570  QSPG(I) = QSSYMB(I,IPHASE)
CZZZ  CALL SPGP( QSPG )
C-
C--  READ AND PRINT ATOM PARAMETERS
C-
3600  WRITE(IOTTO,3610)
3610  FORMAT(/,' *  ATOM PARAMETERS:',//,
     1 ' Label   Atom',11X,'X',9X,'Y',9X,'Z',9X,'B',9X,'N',/,
     1 23X,'b11',7X,'b22',7X,'b33',7X,'b12',7X,'b13',7X,'b23')
      IF(MAGNIN.EQ.1) WRITE(IOTTO,3612)
3612  FORMAT(23x,'KX',8X,'KY',8X,'KZ')
      DO 3620 I=1,NATOM
      IPIOF = I + IOF
      IF(IPIOF.GT.MATMS) THEN
        WRITE(IOTTO,3622)
3622    FORMAT(/,
     1' --The number of atoms for this problem exceeds the alloted',
     1' storage space.',/,
     1'   Increase the value of *MATMS* in all parameter statements.',/)
        STOP '--INPTR:  NUMBER OF ATOMS EXCEEDS PROGRAM STORAGE'
      ENDIF
C-
      READ(IODCF,3630) QATLBL(IPIOF),QATTYP(IPIOF),
     1                 (PATOM(J,IPIOF),J=1,11),QASPG(IPIOF)
3630  FORMAT(2A4,8X,5F8.5,/,6F8.5,1X,A20)
      WRITE(IOTTO,3640) QATLBL(IPIOF),QATTYP(IPIOF),
     1 (PATOM(J,IPIOF),J=1,5),QASPG(IPIOF)(1:17),(PATOM(J,IPIOF),J=6,11)
3640  FORMAT(1X,A4,4X,A4,5X,5F10.5,/,1X,A17,6F10.5)
      IF(MAGNIN.EQ.1) THEN
        READ(IODCF,3642) (VMMAG(J,IPIOF),J=1,3)
3642    FORMAT(3F8.0)
        WRITE(IOTTO,3644) (VMMAG(J,IPIOF),J=1,3)
3644    FORMAT(18X,3F10.5)
      ENDIF
C-
C-  CHECK ANISOTROPIC THERMAL PARAMETERS TO SEE IF THEY ARE USED
      KANISU(IPIOF) = 0
      DO 3650 J=6,11
3650  IF( PATOM(J,IPIOF).NE.0.0 ) KANISU(IPIOF) = 1
C-
3620  CONTINUE
C-
C--  READ AND PRINT GLOBAL PARAMETERS
C-
C-  SCALE AND OVERALL TEMPERATURE FACTOR
3700  READ(IODCF,3720) PPHSCA(IPHASE),PPHOIT(IPHASE)
3720  FORMAT(2F8.0)
      WRITE(IOTTO,3722) PPHSCA(IPHASE), PPHOIT(IPHASE)
3722  FORMAT(/,' *  PHASE RELATED PARAMETERS:',//,
     1 ' Overall scale factor = ',F9.2,/,
     1 ' Overall temp. factor = ',F9.4)
C-
C-  PROFILE-SHAPE-FUNCTION PARAMETERS
      READ(IODCF,3730) (PPHUVW(I,1,IPHASE),I=1,3),NPSF(1,IPHASE),
     1  (PPHPSF(I,1,IPHASE),I=1,2)
3730  FORMAT(3F8.2,I8,2F8.4)
      WRITE(IOTTO,3732) NPSF(1,IPHASE)
3732  FORMAT(/,' Profile-shape-function =',7X,I1)
      NPSF(1,IPHASE) = NPSF(1,IPHASE)
      IF(NPSF(1,IPHASE).EQ.5) THEN
        WRITE(IOTTO,3740) PPHPSF(1,1,IPHASE),PPHPSF(2,1,IPHASE)
3740    FORMAT(' Skewness, Excess =',6X,F8.4,1X,F8.4)
       ELSE IF(NPSF(1,IPHASE).EQ.6) THEN
        WRITE(IOTTO,3742) PPHPSF(1,1,IPHASE),PPHPSF(2,1,IPHASE)
3742    FORMAT(' Pseudo-Voigt  Eta-0,1 =',2(1X,F8.4))
       ELSE IF( NPSF(1,IPHASE).EQ.7 .OR. NPSF(1,IPHASE).EQ.8 ) THEN
        WRITE(IOTTO,3744) PPHPSF(1,1,IPHASE)
3744    FORMAT(' Pearson exponent (m) =',2X,F8.4)
        IF(NPSF(1,IPHASE).EQ.8) WRITE(IOTTO,3746) PPHPSF(2,1,IPHASE)
3746    FORMAT('   m= m0 + m1*2theta, m1=',F8.4)
      ENDIF
      IF(KWSGON.EQ.1) THEN
        WRITE(IOTTO,3734) PPHUVW(1,1,IPHASE)
3734    FORMAT(' Crystallite size =',4X,F10.2,' (nm)')
        PPHUVW(1,1,IPHASE) = PPHUVW(1,1,IPHASE) * 10.0
        PPHUVW(2,1,IPHASE) = PPHUVW(2,1,IPHASE) * 1.0E-3
        WRITE(IOTTO,3736) PPHUVW(2,1,IPHASE)
3736    FORMAT(' Crystallite strain =',2X,F10.6)
       ELSE
        WRITE(IOTTO,3738) (PPHUVW(I,1,IPHASE),I=1,3)
3738    FORMAT(' FWHM parameters U,V,W =',3(1X,F8.4))
      ENDIF
C-
C-  LATTICE PARAMETERS
      READ(IODCF,3760) CELDIR
3760  FORMAT(6F8.2)
CZZZ  CALL CELL2( CELDIR, CELREC )
      DO 3762 I=1,6
      SAVE(I,IPHASE)  = CELDIR(I)
3762  CONTINUE
C     PPHLAT(1,IPHASE) = SQRT( CELREC(1,1) )
C     PPHLAT(2,IPHASE) = SQRT( CELREC(2,2) )
C     PPHLAT(3,IPHASE) = SQRT( CELREC(3,3) )
C     PPHLAT(4,IPHASE) = RTOD * ACOS(CELREC(2,3) /
C    1                  (2.0 * PPHLAT(2,IPHASE) * PPHLAT(3,IPHASE)) )
C     PPHLAT(5,IPHASE) = RTOD * ACOS(CELREC(1,3) /
C    1                  (2.0 * PPHLAT(3,IPHASE) * PPHLAT(1,IPHASE)) )
C     PPHLAT(6,IPHASE) = RTOD * ACOS(CELREC(1,2) /
C    1                  (2.0 * PPHLAT(2,IPHASE) * PPHLAT(1,IPHASE)) )
      WRITE(IOTTO,3770) CELDIR,(PPHLAT(I,IPHASE),I=1,6)
3770  FORMAT(/,' Direct and reciprocal lattice parameters',//,
     1 8X,'A',11X,'B',11X,'C',8X,'Alpha',7X,'Beta',8X,'Gamma',/,
     1 (1X,3(1X,F11.6),3(1X,F9.4,2X)))
      PPHLAT(1,IPHASE) = CELREC(1,1)
      PPHLAT(2,IPHASE) = CELREC(2,2)
      PPHLAT(3,IPHASE) = CELREC(3,3)
      PPHLAT(4,IPHASE) = CELREC(2,3)
      PPHLAT(5,IPHASE) = CELREC(1,3)
      PPHLAT(6,IPHASE) = CELREC(1,2)
C-
C-  PREFERRED ORIENTATION AND ASYMMETRY PARAMETERS
      READ(IODCF,3800) (PPHPOP(I,1,IPHASE),I=1,2),
     1  (POVEC(I,1,IPHASE),I=1,3),
     1  PPHASY(IPHASE),RLIM(IPHASE)
3800  FORMAT(7F8.4)
      WRITE(IOTTO,3550) (POVEC(I,1,IPHASE),I=1,3),
     1 (PPHPOP(I,1,IPHASE),I=1,2)
3550  FORMAT(/,' Preferred orientation vector = ',4X,3F8.4,/,
     1         ' Preferred orientation parameters = ',2F8.4)
      WRITE(IOTTO,3802) RLIM(IPHASE),PPHASY(IPHASE)
3802  FORMAT(/,
     1' Asymmetry correction applied for angles less than ',F6.2,
     1' degrees',/,
     1' Profile asymmetry parameter = ',F8.4)
C-
C--  READ AND ASSIGN CODEWORDS FOR THE VARIABLES
C-
3900  WRITE(IOTTO,3910)
3910  FORMAT(//,' *  CODING OF VARIABLES',//,
     1 ' LABEL',18X,'X',9X,'Y',9X,'Z',9X,'B',9X,'N',/,
     1 23X,'b11',7X,'b22',7X,'b33',7X,'b12',7X,'b13',7X,'b23')
      IF(MAGNIN.EQ.1) WRITE(IOTTO,3612)
C-
C--  READ AND PRINT CODEWORDS FOR ATOM PARAMETERS
C-
      DO 3920 I=1,NATOM
      IPIOF = I + IOF
C-
      READ(IODCF,3930) (AATOM(J,IPIOF),J=1,11)
3930  FORMAT(4X,4X,8X,5F8.2,/,6F8.2)
      WRITE(IOTTO,3940) QATLBL(IPIOF),(AATOM(J,IPIOF),J=1,11)
3940  FORMAT(1X,A4,13X,5F10.2,/,18X,6F10.2)
      IF(MAGNIN.EQ.1) THEN
        READ(IODCF,3642) (VCMAG(J,IPIOF),J=1,3)
        WRITE(IOTTO,3944) (VCMAG(J,IPIOF),J=1,3)
3944    FORMAT(18X,3F10.2)
      ENDIF
C-
C-  GENERATE MATRIX POINTERS AND FRACTIONAL SHIFTS FOR ATOM PARAMETERS
      DO 3950 J=1,11
      X = AATOM(J,IPIOF)
      CPYATM(J,IPIOF) = X
      IYY = INT( ABS(X)/10.0 )
      IF(IYY.GT.KMXVAR) THEN
        IYY = 0
        X = 0.0
        IOV = IOV + 1
      ENDIF
      IF( KCLRCW.EQ.1 .AND. IYY.GT.KMXVAR ) THEN
        AATOM(J,IPIOF) = 0.0
        LATOM(J,IPIOF) = 0
        CPYATM(J,IPIOF) = 0.0
       ELSE
        LATOM(J,IPIOF) = IYY
        AATOM(J,IPIOF) = ( ABS(X) - 10.0 * REAL(IYY) ) * SIGN(1.0,X)
        IMATX(IYY) = IMATX(IYY) + 1
      ENDIF
3950  CONTINUE
C-
C-  GENERATE MATRIX POINTERS AND FRACTIONAL SHIFTS FOR MAG PARAMETERS
      IF(MAGNIN.EQ.1) THEN
        DO 3952 J=1,3
        X = VCMAG(J,IPIOF)
        CWMAGC(J,IPIOF) = X
        IYY = INT( ABS(X)/10.0 )
        IF(IYY.GT.KMXVAR) THEN
          IYY = 0
          X = 0.0
          IOV = IOV + 1
        ENDIF
        IF( KCLRCW.EQ.1 .AND. IYY.GT.KMXVAR ) THEN
          VCMAG(J,IPIOF) = 0.0
          LCMAG(J,IPIOF) = 0
          CWMAGC(J,IPIOF) = 0.0
         ELSE
          LCMAG(J,IPIOF) = IYY
          VCMAG(J,IPIOF) = ( ABS(X) - 10.0 * REAL(IYY) ) * SIGN(1.0,X)
          IMATX(IYY) = IMATX(IYY) + 1
        ENDIF
3952    CONTINUE
      ENDIF
C-
C-  CHECK ANISO THERMAL PARAMETER CODEWORDS TO SEE IF THEY ARE USED
      DO 3960 J=6,11
3960  IF(LATOM(J,IPIOF).GT.0) KANISU(IPIOF) = 1
C-
3920  CONTINUE
C-
C------  READ AND PRINT CODEWORDS FOR PHASE AND GLOBAL PARAMETERS
C-
C-  OVERALL SCALE AND TEMPERATURE FACTOR
4000  READ(IODCF,4010) APHSCA(IPHASE),APHOIT(IPHASE)
4010  FORMAT(2F8.0)
      WRITE(IOTTO,4012) APHSCA(IPHASE),APHOIT(IPHASE)
4012  FORMAT(/,' Overall scale factor = ',F8.2,/,
     1         ' Overall temp. factor = ',F8.2,/)
C-
C-  PROFILE-SHAPE-FUNCTION
      READ(IODCF,4030) (APHUVW(I,1,IPHASE),I=1,3),
     1 (APHPSF(I,1,IPHASE),I=1,2)
4030  FORMAT(3F8.2,8X,2F8.2)
      IF(NPSF(1,IPHASE).EQ.5) THEN
        WRITE(IOTTO,4052) APHPSF(1,1,IPHASE),APHPSF(2,1,IPHASE)
4052    FORMAT(' Skewness, Excess =',5X,2F8.2)
       ELSE IF(NPSF(1,IPHASE).EQ.6) THEN
        WRITE(IOTTO,4054) APHPSF(1,1,IPHASE),APHPSF(2,1,IPHASE)
4054    FORMAT(' Pseudo-Voigt  Eta-0,1 =',2F8.2)
       ELSE IF( NPSF(1,IPHASE).EQ.7 .OR. NPSF(1,IPHASE).EQ.8 ) THEN
        WRITE(IOTTO,4056) APHPSF(1,1,IPHASE)
4056    FORMAT(' Pearson exponent (m) = ',F8.2)
        IF(NPSF(1,IPHASE).EQ.8) WRITE(IOTTO,4058) APHPSF(2,1,IPHASE)
4058    FORMAT(' m= m0 + m1*2theta, m1= ',F8.2)
      ENDIF
      IF(KWSGON.EQ.1) THEN
        WRITE(IOTTO,4032) (APHUVW(I,1,IPHASE),I=1,2)
4032    FORMAT(' Crystallite Size = ',F8.2,'  Strain = ',F8.2)
       ELSE
        WRITE(IOTTO,4034) (APHUVW(I,1,IPHASE),I=1,3)
4034    FORMAT(' FWHM parameters U,V,W =',3F8.2)
      ENDIF
C-
C-  CELL PARAMETERS
      READ(IODCF,4060) (APHLAT(I,IPHASE),I=1,6)
4060  FORMAT(6F8.2)
      WRITE(IOTTO,4062) (APHLAT(I,IPHASE),I=1,6)
4062  FORMAT(/,' Lattice parameters =',3X,6F8.2)
C-
C-  PREFERRED ORIENTATION AND PROFILE ASYMMETRY PARAMETERS
      READ(IODCF,4070) (APHPOP(I,1,IPHASE),I=1,2),APHASY(IPHASE)
4070  FORMAT(2F8.0,8X,8X,8X,F8.0)
      WRITE(IOTTO,4072) (APHPOP(I,1,IPHASE),I=1,2)
4072  FORMAT(/,' Preferred orientation parameters = ',2F8.2)
      WRITE(IOTTO,4074) APHASY(IPHASE)
4074  FORMAT(' Profile asymmetry parameter =',6X,F8.2)
C-
C-- MATRIX POINTERS AND FRACTIONAL SHIFTS FOR PHASE RELATED PARAMETERS
C-
C-  PHASE-RELATED PARAMETERS
4100  CPYPHS( 1,IPHASE) = APHSCA(IPHASE)
      CPYPHS( 2,IPHASE) = APHOIT(IPHASE)
      CPYPHS( 3,IPHASE) = APHUVW(1,1,IPHASE)
      CPYPHS( 4,IPHASE) = APHUVW(2,1,IPHASE)
      CPYPHS( 5,IPHASE) = APHUVW(3,1,IPHASE)
      CPYPHS( 6,IPHASE) = APHLAT(1,IPHASE)
      CPYPHS( 7,IPHASE) = APHLAT(2,IPHASE)
      CPYPHS( 8,IPHASE) = APHLAT(3,IPHASE)
      CPYPHS( 9,IPHASE) = APHLAT(4,IPHASE)
      CPYPHS(10,IPHASE) = APHLAT(5,IPHASE)
      CPYPHS(11,IPHASE) = APHLAT(6,IPHASE)
      CPYPHS(12,IPHASE) = APHPOP(1,1,IPHASE)
      CPYPHS(13,IPHASE) = APHPOP(2,1,IPHASE)
      CPYPHS(14,IPHASE) = APHASY(IPHASE)
      CPYPHS(15,IPHASE) = APHPSF(1,1,IPHASE)
      CPYPHS(16,IPHASE) = APHPSF(2,1,IPHASE)
      DO 4110 I=1,16
      X   = CPYPHS(I,IPHASE)
      IYY = NINT( ABS(X)/10.0 )
      X   = SIGN(1.0,X) * ( ABS(X) - 10.0*REAL(IYY) )
      IF(IYY.GT.KMXVAR) THEN
        IYY = 0
        IOV = IOV + 1
        IF(KCLRCW.EQ.1) THEN
          X = 0.0
          CPYPHS(I,IPHASE) = 0.0
        ENDIF
      ENDIF
      IMATX(IYY) = IMATX(IYY) + 1
      XDEL(I) = X
      IDEL(I) = IYY
4110  CONTINUE
      APHSCA(IPHASE)   = XDEL( 1)
      LPHSCA(IPHASE)   = IDEL( 1)
      APHOIT(IPHASE)   = XDEL( 2)
      LPHOIT(IPHASE)   = IDEL( 2)
      APHUVW(1,1,IPHASE) = XDEL( 3)
      LPHUVW(1,1,IPHASE) = IDEL( 3)
      APHUVW(2,1,IPHASE) = XDEL( 4)
      LPHUVW(2,1,IPHASE) = IDEL( 4)
      APHUVW(3,1,IPHASE) = XDEL( 5)
      LPHUVW(3,1,IPHASE) = IDEL( 5)
      APHLAT(1,IPHASE) = XDEL( 6)
      LPHLAT(1,IPHASE) = IDEL( 6)
      APHLAT(2,IPHASE) = XDEL( 7)
      LPHLAT(2,IPHASE) = IDEL( 7)
      APHLAT(3,IPHASE) = XDEL( 8)
      LPHLAT(3,IPHASE) = IDEL( 8)
      APHLAT(4,IPHASE) = XDEL( 9)
      LPHLAT(4,IPHASE) = IDEL( 9)
      APHLAT(5,IPHASE) = XDEL(10)
      LPHLAT(5,IPHASE) = IDEL(10)
      APHLAT(6,IPHASE) = XDEL(11)
      LPHLAT(6,IPHASE) = IDEL(11)
      APHPOP(1,1,IPHASE) = XDEL(12)
      LPHPOP(1,1,IPHASE) = IDEL(12)
      APHPOP(2,1,IPHASE) = XDEL(13)
      LPHPOP(2,1,IPHASE) = IDEL(13)
      APHASY(IPHASE)   = XDEL(14)
      LPHASY(IPHASE)   = IDEL(14)
      APHPSF(1,1,IPHASE) = XDEL(15)
      LPHPSF(1,1,IPHASE) = IDEL(15)
      APHPSF(2,1,IPHASE) = XDEL(16)
      LPHPSF(2,1,IPHASE) = IDEL(16)
C-
C--  CHECKS ON GLOBAL PARAMETERS
C-
4200  IF((PPHPOP(1,1,IPHASE).EQ.0.0.AND.PPHPOP(2,1,IPHASE).EQ.0.0).AND.
     1   (LPHPOP(1,1,IPHASE).NE.0  .OR. LPHPOP(2,1,IPHASE).NE.0 ) ) THEN
      WRITE(IOTTO,4210)
4210    FORMAT(' --PREFERRED ORIENTATION PARAMETERS MUST BE NON-ZERO',
     1  ' FOR REFINEMENT.')
      ENDIF
C-
C-- LOOKUP ATOM PARAMTERS AND SET REFLECTION DATA
C-
4400  CONTINUE
C4400  CALL LOOKUP( NSCTOT )
C      CALL REFGEN( IPHASE, CELDIR, CELREC )
C      CALL RTMT( IPLT1, IPHASE )
C-
C--  WRITE REFLECTION LIST FOR EACH PHASE
C-
C-  CALCULATE POINTERS INTO REFLECTION LIST FOR LISTING AND MULTS
4500  IEND = 0
      DO 4502 ITMP=1,IPHASE
4502  IEND = IEND + IPREFS(ITMP)
      ISTR = 1
      IF(IPHASE.GT.1) THEN
        DO 4504 ITMP=1,(IPHASE-1)
4504    ISTR = ISTR + IPREFS(ITMP)
      ENDIF
      IOFF = ISTR - 1
C-
C-  WRITE LIST
      IF(LST1.NE.1) GO TO 4600
      IF(LST1.EQ.1) GO TO 4600
      DO 4510 IX=ISTR,IEND
      IRL = MOD( IREFS(IX),             256 ) - 128
      IRK = MOD( IREFS(IX)/256,         256 ) - 128
      IRH = MOD( IREFS(IX)/256/256,     256 ) - 128
      IRC = MOD( IREFS(IX)/256/256/256,   4 )
      D = WAVELN(IRC) / ( 2.0 * SIN(RFTT(IX)*DTORD2) )
      IP = IX - IOFF
      IF( MOD(IP-1,60) .EQ. 0 ) WRITE(IOTTO,4520)
4520  FORMAT('1No. Alpha    H   K   L   Mult  FWHM  ',
     1 ' 2Theta    L/P      D spacing',/)
      WRITE(IOTTO,4530) IP,IRC,IRH,IRK,IRL,MLTT(IX),RFWD(IX),
     1 RFTT(IX),RFLP(IX),D
4530  FORMAT(2I4,3X,3I4,I6,2F8.3,F10.6,2X,F10.6)
4510  CONTINUE
C-
C-  MULTIPLY LP OF EACH REFLECTION BY THE ITS MULTIPLICITY
4600  DO 4610 IX=ISTR,IEND
4610  RFLP(IX) = RFLP(IX) * REAL( MLTT(IX) )
C-
C-  NEXT PHASE
3510  CONTINUE
C-
C--  READ AND ECHO PLOTTING PARAMETER CARDS
C-
C-  LINE PRINTER PLOTS
4800  IF(IPLT.EQ.0) GO TO 4850
      READ(IODCF,4804,END=4808) ISCALE,IDIF
4804  FORMAT(2I8)
4810  WRITE(IOTTO,4812) ISCALE,IDIF
4812  FORMAT(/,
     1' *  Line printer plot parameters: ISCALE = ',I8,'  IDIF = ',I8)
      GO TO 4850
4808  ISCALE = 1200
      IDIF   =  400
      WRITE(IOTTO,4816)
4816  FORMAT(/,
     1' --The printer plot parameter card is missing: defaults assumed')
      GO TO 4810
C-
C-  PLOT FILE PARAMETERS
4850  IF(IPLT2.EQ.0) GO TO 4700
      READ(IODCF,4854,END=4858) IFY,IFYC,IFM,IFD
4854  FORMAT(4I8)
4862  WRITE(IOTTO,4866) IFY,IFYC,IFM,IFD
4866  FORMAT(/,
     1' Plot file parameters:',/,
     1 2X,'Observed pattern...',I1,'  Calculated pattern...',I1,/,
     1 2X,'Peak Markers.......',I1,'  Difference pattern...',I1)
      GO TO 4874
4858  WRITE(IOTTO,4870)
4870  FORMAT(/,
     1' --The plot file parameter card is missing: defaults assumed')
      IFY  = 1
      IFYC = 1
      IFM  = 1
      IFD  = 1
      GO TO 4862
4874  IF(JOBXN.GT.2) THEN
        IFY = 0
        IFD = 0
      ENDIF
C-
C--  MATRIX, FRACTIONAL SHIFTS FOR GLOBAL PARAMTERS
C-
4700  DO 4710 I=1,13
      X = AGLB(I)
      CPYGLB(I) = X
      IYY = INT( ABS(X)/10.0 )
      IF(IYY.GT.KMXVAR) THEN
        IYY = 0
        X = 0.0
        IOV = IOV + 1
      ENDIF
      IF( KCLRCW.EQ.1 .AND. IYY.GT.KMXVAR ) THEN
        AGLB(I) = 0.0
        LGLB(I) = 0
        CPYGLB(I) = 0.0
       ELSE
        LGLB(I) = IYY
        AGLB(I) = ( ABS(X) - 10.0 * REAL(IYY) ) * SIGN(1.0,X)
        IMATX(IYY) = IMATX(IYY) + 1
      ENDIF
4710  CONTINUE
C-
C--  PRINT MATRIX USAGE STATISTICS
C-
5000  WRITE(IOTTO,5010)
5010  FORMAT(/,' Matrix usage statistics (position/references):')
      ISTRT = 1
5020  IFINI = MIN( ISTRT+29, MMATX )
      WRITE(IOTTO,5030) (I,I=ISTRT,IFINI)
      WRITE(IOTTO,5040) (IMATX(I),I=ISTRT,IFINI)
5030  FORMAT(/,1X,30I4)
5040  FORMAT(  1X,30I4)
      ISTRT = ISTRT + 30
      IF(IFINI.LT.MMATX) GO TO 5020
      IMT = 0
      DO 5050 I=1,MMATX
5050  IF( I.LE.KMXVAR .AND. IMATX(I).EQ.0 ) IMT = IMT + 1
      WRITE(IOTTO,5060) IMT,IOV
5060  FORMAT(/,' Number of holes in matrix ',I3,
     1 '  References beyond max # variables ',I3)
      IF(IMT.GT.0) THEN
        WRITE(IOTTO,5070)
5070    FORMAT(/,
     1  ' --THE CURRENT CODEWORDS LEAVE A HOLE IN THE NORMAL EQUATION',
     1  ' MATRIX.',/,
     1  '   MATRIX USAGE MUST BE CONTIGUOUS IF THE PROGRAM IS TO',
     1  ' EXECUTE CORRECTLY.')
        STOP '--HOLES WERE FOUND IN THE NORMAL EQUATION MATRIX'
      ENDIF
C-
      RETURN
      END
C--  SETIO - SETUP FILE I/O NAMES
C-
      SUBROUTINE SETIO
C-
      INCLUDE 'RIETCOM.FOR'
      CHARACTER QFMT1S*36,QFMT2S*39,QFMT3S*33,QFMT4S*35,QFMT5S*37
      CHARACTER QFMT1N*2,QFMT2N*2,QFMT3N*2,QFMT4N*2,QFMT5N*2,QFILE1*40
      EQUIVALENCE ( QFMT1S(26:27), QFMT1N(1:2) )
      EQUIVALENCE ( QFMT2S(29:30), QFMT2N(1:2) )
      EQUIVALENCE ( QFMT3S(23:24), QFMT3N(1:2) )
      EQUIVALENCE ( QFMT4S(26:27), QFMT4N(1:2) )
      EQUIVALENCE ( QFMT5S(28:29), QFMT5N(1:2) )
      DATA QFMT1S/'('' Pattern file name <'',A##,''> ? '')'/
      DATA QFMT2S/'('' Background file name <'',A##,''> ? '')'/
      DATA QFMT3S/'('' Plot file name <'',A##,''> ? '')'/
      DATA QFMT4S/'('' Listing file name <'',A##,''> ? '')'/
      DATA QFMT5S/'('' New input file name <'',A##,''> ? '')'/
C-
C-  GET THE NAME OF THE FILE CONTAINING THE CONTROL PARAMETERS
C-
1000  WRITE(IOTTO,1020)
1020  FORMAT(/,' Control file name <EXIT> ? ')
      READ(IOTTI,1030) QDCFFN
1030  FORMAT(A40)
      IF( QDCFFN.EQ.' ' .OR. QDCFFN.EQ.'/' ) THEN
        WRITE(IOTTO,1035)
1035    FORMAT(/,' -NO CONTROL FILE NAME; EXECUTION TERMINATED')
        STOP '-NO CONTROL FILE NAME; EXECUTION TERMINATED'
      ENDIF
C-
C-  DETERMINE FILE NAME POINTERS
1100  DO 1110 NEND=40,1,-1
1110  IF(QDCFFN(NEND:NEND).EQ.'.') GO TO 1120
      DO 1130 NEND=40,1,-1
1130  IF(QDCFFN(NEND:NEND).NE.' ') GO TO 1140
      STOP '--SETIO:  BLANK CONTROL FILE NAME; ABORT'
1140  NEND = NEND + 1
1120  NENDP3 = NEND + 3
C-
C-  OPEN THE CONTROL DATA FILE NOW
CVAX
C      OPEN(UNIT=IODCF,FILE=QDCFFN,STATUS='OLD',ERR=1040,READONLY)
CVAX
CIBM
      OPEN(UNIT=IODCF,FILE=QDCFFN,STATUS='OLD',ERR=1040)
CIBM
      GO TO 6000
1040  WRITE(IOTTO,1045)
1045  FORMAT(' --SETIO: CONTROL FILE IS BUSY OR NOT THERE; ABORT.')
      STOP '--SETIO: CONTROL FILE IS BUSY OR NOT THERE; ABORT'
C-
6000  QNEWFN = QDCFFN
      QNEWFN(NEND:) = '.NEW'
      WRITE(QFMT5N,2010) NENDP3
2010  FORMAT(I2)
      WRITE(IOTTO,QFMT5S) QNEWFN(1:NENDP3)
      READ(IOTTI,1030) QFILE1
      IF( QFILE1(1:1).NE.' ' .AND. QFILE1(1:1).NE.'/' ) QNEWFN = QFILE1
C-
9000  RETURN
      END
C---  PROGRAM: RIETVELD  (RIET41)
C---  VERSION: 881001  (YYMMDD)
C---
      SUBROUTINE REWRIT
C-
      INCLUDE 'RIETCOM.FOR'
      CHARACTER QFMT*11
      DATA QFMT/'(F8.M,F8.2)'/
C-
C---- OPEN THE NEW OUTPUT FILE -----------------------------------------
C-
CVAX
C1000  OPEN(UNIT=IONEW,FILE=QNEWFN,STATUS='NEW',CARRIAGECONTROL='LIST',
C     1 ERR=1002)
CVAX
CIBM
1000  OPEN(UNIT=IONEW,FILE=QNEWFN,STATUS='UNKNOWN',ERR=1002)
CIBM
      WRITE(IOTTO,1004) QNEWFN
1004  FORMAT(/,' Writing new control data on file: ',A40)
      GO TO 1100
1002  WRITE(IOTTO,1003)
1003  FORMAT(/,' --ERROR WRITING NEW CONTROL DATA FILE:  SKIP IT.')
      GO TO 9000
C-
C---- RECODE THE PARAMETER CODEWORDS -----------------------------------
C-
C-  RECODE THE ATOM PARAMETER CODEWORDS
1100  NATOM = NATM(2,NPHASE) + NATM(1,NPHASE) - 1
      DO 1110 I=1,NATOM
      DO 1112 J=1,11
1112  AATOM(J,I)= CPYATM(J,I)
      DO 1114 J=1,11
1114  VCMAG(J,I)= CWMAGC(J,I)
1110  CONTINUE
C-
C-  RECODE THE SPECIMEN PARAMETER CODEWORDS
      DO 1120 IPHASE=1,NPHASE
      DO 1120 J=1,16
      IF(J.EQ. 1) APHSCA(IPHASE) = CPYPHS(J,IPHASE)
      IF(J.EQ. 2) APHOIT(IPHASE) = CPYPHS(J,IPHASE)
      IF(J.EQ. 3) APHUVW(1,1,IPHASE) = CPYPHS(J,IPHASE)
      IF(J.EQ. 4) APHUVW(2,1,IPHASE) = CPYPHS(J,IPHASE)
      IF(J.EQ. 5) APHUVW(3,1,IPHASE) = CPYPHS(J,IPHASE)
      IF(J.EQ. 6) APHLAT(1,IPHASE) = CPYPHS(J,IPHASE)
      IF(J.EQ. 7) APHLAT(2,IPHASE) = CPYPHS(J,IPHASE)
      IF(J.EQ. 8) APHLAT(3,IPHASE) = CPYPHS(J,IPHASE)
      IF(J.EQ. 9) APHLAT(4,IPHASE) = CPYPHS(J,IPHASE)
      IF(J.EQ.10) APHLAT(5,IPHASE) = CPYPHS(J,IPHASE)
      IF(J.EQ.11) APHLAT(6,IPHASE) = CPYPHS(J,IPHASE)
      IF(J.EQ.12) APHPOP(1,1,IPHASE) = CPYPHS(J,IPHASE)
      IF(J.EQ.13) APHPOP(2,1,IPHASE) = CPYPHS(J,IPHASE)
      IF(J.EQ.14) APHASY(IPHASE) = CPYPHS(J,IPHASE)
      IF(J.EQ.15) APHPSF(1,1,IPHASE) = CPYPHS(J,IPHASE)
      IF(J.EQ.16) APHPSF(2,1,IPHASE) = CPYPHS(J,IPHASE)
1120  CONTINUE
C-
C-  RECODE THE GLOBAL PARAMETER CODEWORDS
      DO 1135 I=1,13
      AGLB(I) = CPYGLB(I)
1135  CONTINUE
C-
C---- WRITE THE PARAMETERS TO THE NEW OUTPUT FILE ----------------------
C-
      WRITE(IONEW,1010) QTITLE(1:79)
1010  FORMAT(A79)
C-
      WRITE(IONEW,1020) JOBXN,NPHASE,NBCKGD,NEXCRG,NSCAT,KINTSF,KAMOR,
     1                  KWSGON
1020  FORMAT(8I4)
C-
      WRITE(IONEW,1030) IOUT,IPLT,IPCH,MAT,NXT,LST1,LST2,LST3,IPLT1,
     1                  IPLT2,KCLRCW,KCONCW
1030  FORMAT(12I4)
C-
      WRITE(IONEW,1022) (KAUTOV(I),I=1,KAUTON)
      WRITE(IONEW,1022) (KAUTOR(I),I=1,KAUTON)
1022  FORMAT(20I4)
C-
      WRITE(IONEW,1040) WAVELN,RATIO(2),WDTTHK,CTHM,THMIN,STEP,THMAX
1040  FORMAT(2F8.5,6F8.4)
C-
      WRITE(IONEW,1050) EPS,RELAX
1050  FORMAT(5F8.4)
C-
      WRITE(IONEW,1140) PGLB(1),PGLB(10),AGLB(1),AGLB(10)
1140  FORMAT('   0',4X,2F8.5,/,8X,2F8.2)
C-
      IF(NBCKGD.EQ.0) THEN
        WRITE(IONEW,1060) (PGLB(J),J=2,7),(AGLB(J),J=2,7)
1060    FORMAT(1P,6E12.5,/,0P,6(F8.2,4X))
       ELSE IF(NBCKGD.GT.1) THEN
        WRITE(IONEW,1062) (POS(I),BCK(I),I=1,NBCKGD)
1062    FORMAT(2F8.2)
      ENDIF
C-
      IF(NEXCRG.GT.0) THEN
        WRITE(IONEW,1062) (ALOW(I),AHIGH(I),I=1,NEXCRG)
      ENDIF
C-
      IF(NSCAT.GT.0) THEN
        DO 1070 I=1,NSCAT
        WRITE(IONEW,1080) QSCNAM(I),DFP(I),DFPP(I)
1080    FORMAT(A4,4X,2F8.4)
C       IF(JOBXN.EQ.1.OR.JOBXN.EQ.3) THEN
          WRITE(IONEW,1090) (SCATAC(J,I),J=1,9)
1090      FORMAT(9F8.5)
C       ENDIF
1070    CONTINUE
      ENDIF
C-
      IF(KINTSF.EQ.1) THEN
        WRITE(IONEW,1150) KINTSN,RINTSQ,KINTSA(1),RINTSA(1),
     1                    KINTSA(2),RINTSA(2),KINTSA(3),RINTSA(3)
1150    FORMAT(4(I8,F8.4))
      ENDIF
C-
      IF(KAMOR.EQ.1) THEN
        WRITE(IONEW,1160) (PGLB(J),J=11,13),(AGLB(J),J=11,13)
1160    FORMAT(F8.4,F8.2,F8.4,/,3F8.2)
      ENDIF
C-
C--  LINE PRINTER AND PLOT FILE PARAMETERS
C-
      IF(IPLT.NE.0) WRITE(IONEW,191) ISCALE,IDIF
191   FORMAT(2I4)
      IF(IPLT2.NE.0) WRITE(IONEW,171) IFY,IFYC,IFM,IFD
171   FORMAT(4I4)
C-
C-- WRITE THE PHASE INFORMATION
C-
2000  DO 2010 IPHASE=1,NPHASE
      NATOM = NATM(1,IPHASE)
      IOF   = NATM(2,IPHASE) - 1
C-
      WRITE(IONEW,2020) QPHSNM(IPHASE)
2020  FORMAT(A32)
      WRITE(IONEW,2030) (QSSYMB(I,IPHASE),I=1,20),QXTALS(IPHASE)
2030  FORMAT(20A1,A20)
      WRITE(IONEW,2040) NATM(1,IPHASE),KATMPF(IPHASE),KMAGSF(IPHASE)
2040  FORMAT(3I4)
C-
C-  ATOM PARAMETERS
      DO 2050 I=1,NATOM
      IPIOF = I + IOF
      WRITE(IONEW,2060) QATLBL(IPIOF),QATTYP(IPIOF),KMAGSI(IPIOF),
     1 (PATOM(J,IPIOF),J=1,5),(AATOM(J,IPIOF),J=1,5)
2060  FORMAT(A8,A4,I4,8X,5F8.5,/,24X,5F8.2)
      WRITE(IONEW,2062) QASPG(IPIOF),(PATOM(J,IPIOF),J=6,11),
     1 (AATOM(J,IPIOF),J=6,11)
2062  FORMAT(A20,4X,6F8.5,/,24X,6F8.2)
      IF(KMAGSI(IPIOF).EQ.1) THEN
        WRITE(IONEW,2064) (VMMAG(J,IPIOF),J=1,3),(VCMAG(J,IPIOF),J=1,3)
2064    FORMAT(24X,3F8.2,/,24X,3F8.5)
      ENDIF
2050  CONTINUE
C-
C-  OVERALL SCALE AND ISO TEMP FACTORS
      N = 7 - ( INT( ALOG10( AMAX1( PPHSCA(IPHASE), 1.0 ) ) ) + 1 )
      WRITE(QFMT(5:5),'(I1)') N
      WRITE(IONEW,QFMT) PPHSCA(IPHASE),PPHOIT(IPHASE)
      WRITE(IONEW,2160) APHSCA(IPHASE),APHOIT(IPHASE)
2160  FORMAT(2F8.2)
C-
C-  PROFILE-SHAPE-FUNCTION PARAMETERS
      WRITE(IONEW,2080) NPSF(1,IPHASE),RLIM(IPHASE),PPHASY(IPHASE),
     1 (PPHPSF(I,1,IPHASE),I=1,2),APHASY(IPHASE),
     1 (APHPSF(I,1,IPHASE),I=1,2)
2080  FORMAT(I4,4X,4F8.4,/,16X,3F8.2)
C-
C-  PROFILE BROADENING COEFFICIENTS
      IF(KWSGON.EQ.0) THEN
        WRITE(IONEW,2082) (PPHUVW(I,1,IPHASE),I=1,3),
     1  (APHUVW(I,1,IPHASE),I=1,3)
2082    FORMAT(8X,3F8.5,/,8X,3F8.2)
       ELSE
        PPHUVW(2,1,IPHASE) = PPHUVW(2,1,IPHASE) * 1.0E+3
        WRITE(IONEW,2084) PPHUVW(1,1,IPHASE)/10.0,PPHUVW(2,1,IPHASE),
     1  (APHUVW(I,1,IPHASE),I=1,2)
2084    FORMAT(F8.2,F8.6,/,2F8.2)
      ENDIF
C-
C-  LATTICE PARAMETERS
      WRITE(IONEW,2100) (SAVE(J,IPHASE),J=1,6),(APHLAT(I,IPHASE),I=1,6)
2100  FORMAT(3F8.5,3F8.4,/,6F8.2)
C-
C-  PREFERRED ORIENTATION PARAMETERS
      WRITE(IONEW,2110) (POVEC(I,1,IPHASE),I=1,3),
     1 (PPHPOP(I,1,IPHASE),I=1,2),(APHPOP(I,1,IPHASE),I=1,2)
2110  FORMAT('   0',4X,5F8.4,/,32X,2F8.2)
C-
C-  NEXT PHASE
2010  CONTINUE
C-
C-  CLOSE THE OUTPUT FILE
      CLOSE(UNIT=IONEW)
C-
C-  RETURN TO CALLER
9000  RETURN
      END
      SUBROUTINE DATE( QDATE )
      CHARACTER QDATE*9
      QDATE = 'DD-MMM-YY'
      RETURN
      END
      SUBROUTINE TIME( QDATE )
      CHARACTER QDATE*8
      QDATE = 'HH:MM:SS'
      RETURN
      END

