Subroutine WRITEPHASE for program GSAS2CIF
This subroutine is used to write parameters and results for each phase to the
output CIF file.
See the gsas2cif documentation for an
explanation of this code.
SUBROUTINE WRITEPHASE(IUCIF,IUEXP,IUTRM,IPHAS,NPHAS,DAYTIME,
1 ONEBLOCK,MATRX,NUMPAR,MBW,IUDIS)
!PURPOSE: write information about the phase
INCLUDE '../INCLDS/COPYRIGT.FOR'
!PSEUDOCODE:
!CALLING ARGUMENTS:
INTEGER*4 IUCIF
INTEGER*4 IUEXP
INTEGER*4 IUTRM
INTEGER*4 IPHAS
INTEGER*4 NPHAS(9) !Phase existence flags
CHARACTER*20 DAYTIME
LOGICAL*4 ONEBLOCK ! true if the CIF will have one block
REAL*4 MATRX(1)
INTEGER*4 NUMPAR !Number of refined parameters
INTEGER*4 MBW !Matrix bandwidth
!INCLUDE STATEMENTS:
INCLUDE '../INCLDS/ARRAYSZE.FOR'
INCLUDE '../INCLDS/SPGCOMI.FOR'
INCLUDE '../INCLDS/HEADSCOM.FOR'
INCLUDE '../INCLDS/DISAGLCM.FOR'
INCLUDE '../INCLDS/CELLCOM.FOR'
!LOCAL VARIABLES:
INTEGER*4 IOPRTNS(50)
INTEGER*4 ISAM
INTEGER*4 JMLT(MAXATM) !Atom site multiplicities
INTEGER*4 NSYS(14)
1 /1,2,3,4,4,5,5,6,6,6,7,7,8,8/
REAL*4 ANGLES(3)
REAL*4 ANGSIG(3)
REAL*4 RM(6) !Recip. metr. tensor
REAL*4 VOLUME !Unit cell volume (=0 for error)
REAL*4 UB(3,3) !UB-matrix
LOGICAL*4 ANIFLAG
CHARACTER*1 CLBL(3)
1 /'a','b','c'/
CHARACTER*5 ALBL(3)
1 /'alpha','beta ','gamma'/
CHARACTER*20 STRING(10)
CHARACTER*80 TEXT !ISAM file read write buffer
CHARACTER*12 SYST(8)
1 /'triclinic ','monoclinic ','orthorhombic',
1 'tetragonal ','trigonal ','trigonal ',
1 'hexagonal ','cubic '/
CHARACTER*4 XYZLBL(9)
1 /'-z ','-y ','-x ','x-y ','ERR ','y-x ',
1 '+x ','+y ','+z '/
CHARACTER*4 TRA(13)
1 /' ','ERR ','+1/6','+1/4','+1/3','ERR ',
1 '+1/2','ERR ','+2/3','+3/4','+5/6','ERR ',' '/
CHARACTER*4 OUTL(6,2)
REAL*4 CONC(MAXELEM)
LOGICAL*4 NOTDONE
INTEGER*4 NOFFSET ! number of symmetry positions
INTEGER*4 OFFSYMID(192) ! symmetry ID needing offset correction
INTEGER*4 OFFSET(192) ! offset to be added to 100*x+10*y+z
CHARACTER*1 MSG
REAL*4 COFF(MAXODF) ! Spherical harmonic coefficients
INTEGER*4 INDX(3,MAXODF) ! Sph. harmonic index
INTEGER*4 ISAMSYM ! Sample symmetry number (1-4)
CHARACTER*12 KEYVAL ! ISAM key
CHARACTER*2 ELEMTBL(MAXELEM) ! Table of unique elements
CHARACTER*2 ELEM
REAL*4 COMPTBL(MAXELEM) ! Number of atoms for each unique element/cell -- full occupied sites
REAL*4 FRACTBL(MAXELEM) ! Number of atoms for each unique element/cell -- partially occupied sites
REAL*4 MASSTBL(MAXELEM) ! Mass for each type of atom
INTEGER*2 SEQTBL(MAXELEM) ! Sequence to show elements
INTEGER*4 Z
LOGICAL FLAG
CHARACTER*1 PUBFLG ! Y is distances/angles will be published
!SUBROUTINES CALLED:
!FUNCTION DEFINITIONS:
INTEGER*4 READEXP !ISAM file read function
CHARACTER*6 CRSKEY !ISAM key building routine
!Code:
NOFFSET = 0
TEXT = ' '
IULST = 6
CALL DSGREAD(IUEXP,IULST,IPHAS,NPHAS(IPHAS),'NOFA',NUMPAR,
1 MBW,MATRX) !Read unit cell and atom data
ISAM = READEXP(IUEXP,CRSKEY(IPHAS)//' PNAM',TEXT)
CALL WRVAL(IUCIF, '_pd_phase_name', text)
ISAM = READEXP(IUEXP,CRSKEY(IPHAS)//'ANGLES',TEXT)
READ (TEXT,'(3F10.0)') ANGLES
ISAM = READEXP(IUEXP,CRSKEY(IPHAS)//'ANGSIG',TEXT)
IF ( ISAM.EQ.0 ) THEN
READ (TEXT,'(3F10.0)') ANGSIG
ELSE
DO I=1,3
ANGSIG(I) = 0.0
END DO
END IF
CALL BMATRX(ABC(1,IPHAS),ANGLES,UB,
1 ABCST(1,IPHAS),CANGST(1,IPHAS))
DO I=1,3
IF ( (LAUE.GT.3 .AND. I.GT.1) .AND. !Laue symmetry above mmm
1 (I.EQ.2 .OR.
1 (((LAUE.EQ.6 .OR. LAUE.EQ.7) .OR. !Rhombohedral symmetry
1 LAUE.GT.11) .AND. I.EQ.3)) ) THEN !Cubic symmetry
CALL FESD(ABC(I,IPHAS), -CELSIG(I), text, ln)
CALL WRVAL(IUCIF, '_cell_length_'//clbl(I),text)
ELSE IF (CELSIG(I) .le. 0.0) THEN
CALL FESD(ABC(I,IPHAS), 0.0, text, ln)
CALL WRVAL(IUCIF, '_cell_length_'//clbl(I),text)
ELSE
CALL FESD(ABC(I,IPHAS), CELSIG(I), text, ln)
CALL WRVAL(IUCIF,'_cell_length_'//clbl(I),text)
END IF
END DO
DO I=1,3
NOTDONE = .TRUE.
IF ( LAUE.GT.1 ) THEN
IF ( (LAUE.EQ.6 .OR. LAUE.EQ.7) ) THEN !Rhombohedral setting
IF (I .gt. 1) THEN
NOTDONE = .FALSE.
END IF
ELSE IF ( (LAUE.GT.7 .AND. LAUE.LT.12) .AND. I.EQ.3 ) !Hexagonal cell, Gamma angle
1 THEN
NOTDONE = .FALSE.
ELSE IF ( (LAUE.EQ.2 .AND. NAXIS.NE.I) .OR. !Monoclinic, not the unique axis
1 LAUE.GT.2 ) THEN !Anything else
NOTDONE = .FALSE.
END IF
END IF
IF ( NOTDONE .and. ANGSIG(I) .gt. 0) THEN
CALL FESD(ANGLES(I), ANGSIG(I), text, ln)
CALL WRVAL(IUCIF, '_cell_angle_'//albl(I),text)
ELSE IF (ANGSIG(I) .gt. 0) THEN
CALL FESD(ANGLES(I), -ANGSIG(I), text, ln)
CALL WRVAL(IUCIF, '_cell_angle_'//albl(I),text)
ELSE
CALL FESD(ANGLES(I), 0.0, text, ln)
CALL WRVAL(IUCIF, '_cell_angle_'//albl(I),text)
END IF
END DO
CALL CELVOL(ABC(1,IPHAS),ANGLES,RM,VOLUME)
CALL FESD(VOLUME, 0.0, text, ln)
CALL WRVAL(IUCIF, '_cell_volume',text)
CALL WRVAL(IUCIF,'_symmetry_cell_setting',SYST(NSYS(LAUE)))
WRITE(text,'(20a1)') SPG
ln = LENCH(text)
C a R suffix is a GSAS code for a rhombohedral setting & should be removed
if (text(ln:ln) .eq. 'R') text(ln:ln) = ' '
CALL WRVAL(IUCIF,'_symmetry_space_group_name_H-M',
1 text(1:LENCH(TEXT)))
WRITE (IUCIF,'(2A)') 'loop_ _symmetry_equiv_pos_site_id',
1 ' _symmetry_equiv_pos_as_xyz'
DO ICV=1,NCV !Loop over the lattice centering
DO JCEN=0,ICEN !Loop over the identity and inversion
DO I=1,NSYM !Loop through the matrices
IM = 100
IOFF = 0
DO J=1,3
IJ = 2*JRT(J,1,I)+3*JRT(J,2,I)+4*JRT(J,3,I)+5
IK1 = JRT(J,4,I)+NINT(CEN(J,ICV)*12.0)+1
IK = MOD(JRT(J,4,I)+NINT(CEN(J,ICV)*12.0),12)+1
C has a offset been applied to the symmetry operator?
IF (IK .NE. IK1) THEN
IF (JCEN .EQ. 1) THEN
IOFF = IOFF - IM*(IK1-IK)/12
ELSE
IOFF = IOFF + IM*(IK1-IK)/12
END IF
END IF
IM = IM/10
IF ( JCEN.EQ.0 ) THEN
OUTL(J,1) = XYZLBL(IJ)
OUTL(J,2) = TRA(IK)
ELSE
IJ = 10-IJ
OUTL(J,1) = XYZLBL(IJ)
IK = 14-IK
OUTL(J,2) = TRA(IK)
END IF
END DO
I1MX = 3
TEXT = ' '
LN = 1
DO I1=1,I1MX
DO I2=1,2
LNX = LENCH(outl(i1,i2))
IF ( LNX.GT.0 ) THEN
TEXT(LN:) = outl(i1,i2)(:LNX)
LN = LN+LNX
END IF
END DO
IF ( MOD(I1,3).NE.0 ) THEN
TEXT(LN:LN) = ','
LN = LN+1
ELSE
K = 100*(ICV-1) + I
IF (JCEN .EQ. 1) K = -K
WRITE (IUCIF,'(3X,I5,1x,A)') K, TEXT(:LN)
IF (IOFF .NE. 0) THEN
NOFFSET = NOFFSET +1
IF (NOFFSET .GT. 192) THEN
PRINT '(A)','More than 192'//
1 'Offset symmetry ops -- how did this happen!'
STOP 'NOFFSET > 192'
END IF
OFFSYMID(NOFFSET) = K
OFFSET(NOFFSET) = IOFF
END IF
LN = 1
END IF
END DO
END DO
END DO
END DO
C initialize chemical formula arrays
DO I=1,MAXELEM
ELEMTBL(I) = ' '
COMPTBL(I) = 0.
FRACTBL(I) = 0.
MASSTBL(I) = 0.
SEQTBL(I) = 0
END DO
Z = 1
WRITE(IUCIF,'(/A/)') '# ATOMIC COORDINATES'//
1 ' AND DISPLACEMENT PARAMETERS'
WRITE(IUCIF,'(/a5)') 'loop_'
WRITE(IUCIF,'(6x,A)') '_atom_site_type_symbol',
1 '_atom_site_label',
1 '_atom_site_fract_x',
1 '_atom_site_fract_y','_atom_site_fract_z',
1 '_atom_site_occupancy',
1 '_atom_site_thermal_displace_type',
1 '_atom_site_U_iso_or_equiv',
1 '_atom_site_symmetry_multiplicity'
ANIFLAG = .false.
C Warn on duplicate labels
DO J=2,NATOM
DO I=1,J-1
IF (ATMNAM(I)(:LENCH(ATMNAM(I))) .EQ.
1 ATMNAM(J)(:LENCH(ATMNAM(J)))
1 .AND. FRACT(I) .NE. 0 .AND. FRACT(J) .NE. 0) THEN
PRINT '(A,i5,A,I5,2A)','Atoms',I,' and ',J,
1 ' are both labeled',ATMNAM(I)
CALL REDTRML('This is not allowed in a CIF'//
1 'Continue anyway? (,N)',MSG)
CALL UPCASE(MSG)
IF (MSG .EQ. 'N') STOP
END IF
END DO
END DO
NA = 0
DO I=1,NATOM
IF (FRACT(I) .NE. 0 .OR. SGFRAC(I) .NE. 0) THEN
NA = NA + 1
CALL SYTSYM (XYZ(1,I),ICEN,NSYM,JRT,NCV,CEN,LAUE, !Get site symmetries and multiplicities
1 JMLT(I),JSYM,IOPRTNS,STSYM(I))
C _atom_site_label
C Note: atom labels must be unique -- if need be, concatenate (_xxx)
C where xxx is the atom number -- not implemented at this time, instead warn (above)
string(1) = ATMNAM(I)(:LENCH(ATMNAM(I)))
DO K=1,3
CALL FESD(XYZ(K,I), SXYZ(K,I), string(K+1), ln)
END DO
CALL FESD(FRACT(I), SGFRAC(I), string(5), ln)
IF ( REFCODE(I)(1:1).EQ.'I' ) THEN
STRING(6) = 'Uiso'
ELSE
ANIFLAG = .true.
STRING(6) = 'Uani'
END IF
IF ( REFCODE(I)(1:1).EQ.'I' ) THEN
CALL FESD(BIJ(1,I), SBIJ(1,I), string(7), ln)
ELSE
UEQV = 0.0
DO IJ=1,3
IJ1 = IJ+1
IF ( IJ1.EQ.4 ) IJ1=3
IJ2 = 1
IF ( IJ.EQ.3 ) IJ2=2
UEQV = UEQV
1 +BIJ(IJ,I)*(ABC(IJ,IPHAS)*ABCST(IJ,IPHAS))**2
1 +2.0*BIJ(IJ+3,I)*ABC(IJ2,IPHAS)*ABC(IJ1,IPHAS)
1 *ABCST(IJ2,IPHAS)*ABCST(IJ1,IPHAS)
1 *CANG(4-IJ,IPHAS)
END DO
CALL FESD(UEQV/3.0, -0.0001, string(7), ln)
END IF
WRITE(string(8),'(i4)') JMLT(I)
CALL VSTRNG(ATMTYP(I),LENCH(ATMTYP(I)),.true.,.false.)
write(IUCIF,'(A)') ATMTYP(I)
write(IUCIF,'(A6,4a13,a5,a13,a4,A)') (string(J),J=1,8)
C enter the atom into the composition table
ELEM = ATMTYP(I)(1:2)
C is this a one-letter or two-letter element?
JCH = ICHAR(ELEM(2:2))
IF (JCH .LT. ICHAR('A') .OR. JCH .GT. ICHAR('Z')) THEN
KEYVAL = ' AFAC '//ELEM(1:1)//'_'
ELEM(2:2) = ' '
ELSE
KEYVAL = ' AFAC '//ELEM(1:2)//'_'
ELEM(2:2) = CHAR(JCH+32)
ENDIF
J = 1
DO WHILE (ELEMTBL(J) .NE. ' ' .AND. ELEMTBL(J) .NE. ELEM)
J = J + 1
END DO
ELEMTBL(J) = ELEM
TEXT = ' '
ISAM = READEXP(IUEXP,KEYVAL,TEXT)
READ (TEXT,'(F7.0)') MASSTBL(J)
IF (FRACT(I) .EQ. 1.0) THEN
COMPTBL(J) = COMPTBL(J) + JMLT(I)
ELSE
FRACTBL(J) = FRACTBL(J) + JMLT(I)*FRACT(I)
ENDIF
END IF
END DO
IF (K .EQ. 0) THEN
WRITE(IUCIF,'(A)') ' ? ? ? ? ? ? ? ? ?'
ELSE IF (ANIFLAG) THEN
WRITE(IUCIF,'(/a5,1x,A)') 'loop_', '_atom_site_aniso_label'
WRITE(IUCIF,'(6x,A)') '_atom_site_aniso_U_11',
1 '_atom_site_aniso_U_12','_atom_site_aniso_U_13',
1 '_atom_site_aniso_U_22','_atom_site_aniso_U_23',
1 '_atom_site_aniso_U_33'
DO I=1,NATOM
IF (FRACT(I) .NE. 0 .OR. SGFRAC(I) .NE. 0) THEN
IF ( REFCODE(I)(1:1).EQ.'A' ) THEN
string(1) = ATMNAM(I)(:LENCH(ATMNAM(I)))
DO K=1,6
CALL FESD(BIJ(K,I), SBIJ(K,I), string(K+1), ln)
END DO
write(IUCIF,'(A6,6a13)') (string(J),J=1,7)
END IF
END IF
END DO
END IF
C ======================================================================
C Loop over element types -- but only if the histogram info goes in a
C separate block. In a single-block histogram, this info is included
C with the scattering factor information (WRPOWDHI)
C
IF (.NOT. ONEBLOCK) THEN
C determine unit cell contents
DO I=1,MAXELEM
conc(i) = 0.0
END DO
DO I=1,NATOM
J = ID(I) !Get the atom type count flag
conc(J) = conc(j) + JMLT(I)*FRACT(I)
END DO
WRITE(IUCIF,'(/a5,1x,A)') 'loop_', '_atom_type_symbol'
WRITE(IUCIF,'(6x,A)') '_atom_type_number_in_cell'
DO j=1,MAXELEM
if (conc(j) .gt. 0.0) then
string(1) = ATYPE(j)(1:2)
CALL VSTRNG(string(1),2,.false.,.false.)
CALL FESD(conc(j), -0.01, string(2), ln)
write(IUCIF,'(t20,A2,a13)') (string(I),I=1,2)
END IF
END DO
END IF
C process the chemical formula: pick a Z value & generate molecular weight
C find the maximum possible Z value
N = 0
DO I=1,MAXELEM
IF (ELEMTBL(I) .NE. ' ') N = I
END DO
C factors of 2
FLAG = .TRUE.
DO WHILE(FLAG)
DO I=1,N
IF (Z*2.*INT(COMPTBL(I)/(Z*2.)) .NE. COMPTBL(I))
1 FLAG = .FALSE.
END DO
IF (FLAG) Z = Z * 2
END DO
C factors of 3
FLAG = .TRUE.
DO WHILE(FLAG)
DO I=1,N
IF (Z*3.*INT(COMPTBL(I)/(Z*3.)) .NE. COMPTBL(I))
1 FLAG = .FALSE.
END DO
IF (FLAG) Z = Z * 3
END DO
C order the elements in "Hill" order: C,H & alphabetical or alphabetical
C is C present?
FLAG = .FALSE.
J = 1
DO I=1,N
IF (ELEMTBL(I) .EQ. 'C ') THEN
FLAG = .TRUE.
SEQTBL(I) = J
J = J + 1
END IF
END DO
C if yes, get H
IF (FLAG) THEN
DO I=1,N
IF (ELEMTBL(I) .EQ. 'H ') THEN
SEQTBL(I) = J
J = J + 1
END IF
END DO
DO I=1,N
IF (ELEMTBL(I) .EQ. 'D ') THEN
SEQTBL(I) = J
J = J + 1
END IF
END DO
END IF
DO K=1,N
ELEM = 'Zz'
NUMELEM = 100*ICHAR(ELEM(1:1)) + ICHAR(ELEM(2:2))
NN = 0
DO I=1,N
NUMELEM1 = 100*ICHAR(ELEMTBL(I)(1:1)) + ICHAR(ELEMTBL(I)(2:2))
IF (NUMELEM1 .LT. NUMELEM .AND. SEQTBL(I) .EQ. 0) THEN
NN = I
NUMELEM = NUMELEM1
END IF
END DO
IF (NN .NE. 0) THEN
SEQTBL(NN) = J
J = J + 1
END IF
END DO
K = 1
ATMASS = 0
DO J=1,N
DO I=1,N
IF (SEQTBL(I) .EQ. J) THEN
TEXT(K:) = ELEMTBL(I)
IF (ELEMTBL(I)(2:2) .EQ. ' ') THEN
K = K + 1
ELSE
K = K + 2
ENDIF
IF (FRACTBL(I) .NE. 0) THEN
WRITE(KEYVAL,'(F12.2)') (COMPTBL(I) + FRACTBL(I))/Z
ELSE
WRITE(KEYVAL,'(I12)') NINT(COMPTBL(I)/Z)
ENDIF
ATMASS = ATMASS + MASSTBL(I) * (COMPTBL(I) + FRACTBL(I))/Z
NN = 1
DO WHILE (KEYVAL(NN:NN) .EQ. ' ')
NN = NN + 1
END DO
IF (KEYVAL(NN:) .NE. '1') THEN ! values of 1 are assumed
TEXT(K:) = KEYVAL(NN:)
K = K + 13 - NN
END IF
END IF
END DO
TEXT(K:K) = ' ' ! leave a blank space
K = K + 1
END DO
WRITE(IUCIF,'(A)') ' ',
1 '# If you change Z, be sure to change all 3 of the following'
CALL WRVAL(IUCIF, '_chemical_formula_sum',text)
WRITE(TEXT,'(F15.2)') ATMASS
CALL WRVAL(IUCIF, '_chemical_formula_weight',text)
WRITE(TEXT,'(I4)') Z
CALL WRVAL(IUCIF, '_cell_formula_units_Z',text)
C Spherical harmonic ODF
IODF = 0
IF (.NOT. ONEBLOCK) THEN
KEYVAL = CRSKEY(IPHAS)//'ODF '
ISAM = READEXP(IUEXP,KEYVAL,TEXT)
READ (TEXT,'(3I5)') NORD,NODFCOF,ISAMSYM
IF (NODFCOF .GT. 0) THEN
WRITE(IUCIF,'(/A)') '_pd_proc_ls_pref_orient_corr'
WRITE(IUCIF,'(2A)') ';',' Spherical Harmonic ODF'
WRITE(IUCIF,'(A,I2,A,I3)')
1 ' PHASE',I,' spherical harmonic order=',NORD
IF ( ISAMSYM.EQ.1 ) THEN
WRITE(IUCIF,'(A)') ' No sample symmetry'
ELSE IF ( ISAMSYM.EQ.2 ) THEN
WRITE(IUCIF,'(2A)') ' The sample symmetry is:',
1 ' 2/m (shear texture)'
ELSE IF ( ISAMSYM.EQ.3 ) THEN
WRITE(IUCIF,'(2A)') ' The sample symmetry is:',
1 ' mmm (rolling texture)'
ELSE IF ( ISAMSYM.EQ.0 ) THEN
WRITE(IUCIF,'(2A)') ' The sample symmetry is:',
1 ' cylindrical (fiber texture)'
END IF
NREC = 0
IF ( NODFCOF.GT.0 ) NREC = (NODFCOF-1)/6+1
IBEG = 1
DO IREC=1,NREC
WRITE(KEYVAL(10:12),'(I2,A)')IREC,'A'
IFIN = MIN(IBEG+5,NODFCOF)
ISAM = READEXP(IUEXP,KEYVAL,TEXT)
READ (TEXT,'(6(I4,2I3))') ((INDX(K,J),K=1,3),
1 J=IBEG,IFIN)
KEYVAL(12:12) = 'B'
ISAM = READEXP(IUEXP,KEYVAL,TEXT)
READ (TEXT,'(6(F10.0))') (COFF(K),K=IBEG,IFIN)
IBEG = IBEG+6
END DO
DO J=1,NODFCOF
WRITE(IUCIF,'(A,3I3,3x,A,F10.4)')
1 ' Index =',(INDX(K,J),K=1,3),
1 'Coeff=',COFF(J)
END DO
WRITE(IUCIF,'(a/)') ';'
END IF
END IF
C now loop over interatomic distances for this phase
WRITE(IUCIF,'(/a)') '# MOLECULAR GEOMETRY'
WRITE(IUCIF,'(/a5)') 'loop_'
WRITE(IUCIF,'(6x,A)') '_geom_bond_atom_site_label_1'
WRITE(IUCIF,'(6x,A)') '_geom_bond_atom_site_label_2'
WRITE(IUCIF,'(6x,A)') '_geom_bond_distance'
WRITE(IUCIF,'(6x,A)') '_geom_bond_site_symmetry_1'
WRITE(IUCIF,'(6x,A)') '_geom_bond_site_symmetry_2'
WRITE(IUCIF,'(6x,A)') '_geom_bond_publ_flag'
IDIS = 0
IF (IUDIS .NE. 0) THEN
REWIND(IUDIS)
READ (IUDIS,'(A)') ! skip the first record
KPHAS = 0
DO WHILE(KPHAS .LE. IPHAS)
1 READ (IUDIS,'(A1,2I2,2F10.4,7I5)',ERR=1,END=2)
1 PUBFLG,KPHAS,ITYP,D,STD,I,J,ISYM,IOFF
IF (KPHAS .EQ. IPHAS .AND. ITYP .EQ. 0 .AND.
1 (FRACT(I) .NE. 0 .OR. SGFRAC(I) .NE. 0) .AND.
1 (FRACT(J) .NE. 0 .OR. SGFRAC(J) .NE. 0)) THEN
IF (STD .LE. 0) STD = -0.0001
CALL FESD(D,STD, text, ln)
DO K=1,NOFFSET
IF (ISYM .EQ. OFFSYMID(K)) THEN
IOFF = OFFSET(K) + IOFF
GOTO 3
END IF
END DO
3 CONTINUE
CALL UPCASE(PUBFLG)
IF (PUBFLG .NE. 'Y') THEN
PUBFLG = 'N'
END IF
WRITE (IUCIF,'(2(2x,A),2x,A16,2x,A1,2X,I5,A1,I3,2x,A1)')
1 ATMNAM(I),ATMNAM(J),text(:ln),'.',ISYM,'_',IOFF,PUBFLG
IDIS = IDIS + 1
END IF
END DO
2 CONTINUE
END IF
IF (IDIS .EQ. 0) WRITE(IUCIF,'(A)') ' ? ? ? ? ? ?'
C now loop over interatomic angles for this phase
WRITE(IUCIF,'(/a5)') 'loop_'
WRITE(IUCIF,'(6x,A)') '_geom_angle_atom_site_label_1'
WRITE(IUCIF,'(6x,A)') '_geom_angle_atom_site_label_2'
WRITE(IUCIF,'(6x,A)') '_geom_angle_atom_site_label_3'
WRITE(IUCIF,'(6x,A)') '_geom_angle'
WRITE(IUCIF,'(6x,A)') '_geom_angle_site_symmetry_1'
WRITE(IUCIF,'(6x,A)') '_geom_angle_site_symmetry_2'
WRITE(IUCIF,'(6x,A)') '_geom_angle_site_symmetry_3'
WRITE(IUCIF,'(6x,A)') '_geom_angle_publ_flag'
IANG = 0
IF (IUDIS .NE. 0) THEN
REWIND(IUDIS)
READ (IUDIS,'(A)') ! skip the first record
KPHAS = 0
DO WHILE(KPHAS .LE. IPHAS)
11 READ (IUDIS,'(A1,2I2,2F10.4,7I5)',ERR=11,END=12)
1 PUBFLG,KPHAS,ITYP,D,STD,I,J,K,ISYM1,IOFF1,ISYM3,IOFF3
IF (KPHAS .EQ. IPHAS .AND. ITYP .EQ. 1 .AND.
1 (FRACT(I) .NE. 0 .OR. SGFRAC(I) .NE. 0) .AND.
1 (FRACT(K) .NE. 0 .OR. SGFRAC(K) .NE. 0) .AND.
1 (FRACT(J) .NE. 0 .OR. SGFRAC(J) .NE. 0)) THEN
IF (STD .LE. 0) STD = -0.001
CALL FESD(D,STD, text, ln)
DO K1=1,NOFFSET
IF (ISYM1 .EQ. OFFSYMID(K1)) THEN
IOFF1 = OFFSET(K) + IOFF1
END IF
IF (ISYM3 .EQ. OFFSYMID(K1)) THEN
IOFF3 = OFFSET(K) + IOFF3
END IF
END DO
CALL UPCASE(PUBFLG)
IF (PUBFLG .NE. 'Y') THEN
PUBFLG = 'N'
END IF
WRITE (IUCIF,'(3(2x,A),2x,A16,2x,
1 I5,A1,I3,2x,A1,2X,I5,A1,I3,2X,A1)')
1 ATMNAM(I),ATMNAM(J),ATMNAM(K),text(:ln),
1 ISYM1,'_',IOFF1,'.',ISYM3,'_',IOFF3,PUBFLG
IANG = IANG + 1
END IF
END DO
12 CONTINUE
END IF
IF (IANG .EQ. 0)
1 WRITE(IUCIF,'(A)') ' ? ? ? ? ? ? ? ?'
RETURN
END