      PROGRAM FOUR
*****************************  File: U625002 FOUR FORTRAN N  =  FOUR.FOR
****  Source: MULTAN80 / SGG 87
****  The second part: 'FOUR2' is also used in TRACOR
*****************************************  Last update:     11 Nov. 1999
      COMMON /SYSTA/ IFILE(20), IFSTAT(20), ITIME(4),   KEYS(28),
     *               NFNUM,     NLIT,       NCOLN(32),  NCOLL(32),
     *               NFDOT(32), NFDOL(32),  NLUSER(32), FNUM(32),
     *               SWITCH(28)
      EQUIVALENCE (LIS1, IFILE(7)), (LIS2, IFILE(8))
      EQUIVALENCE (IRUN, KEYS(13))
      PARAMETER (KUSER2=30000, KUSER1=KUSER2/3)
      CALL KEPROG('FOUR')
      WRITE (LIS2, FMT = '('' Last FOUR update: 11 Nov. 1999'')')
      IF (IRUN .GT. 1) WRITE (LIS1, FMT = '(66X, ''RUN'', I3)') IRUN
      CALL FFTIN(KUSER1)
      CALL PP1
      CALL SEARCH
      IDDOKA = KEYS(10)
      CALL KEPROX
      WRITE (LIS2, FMT='(/'' Test: DDOKA exit MAIN SUBPROGRAM ''/)')
      STOP 99
      END
      SUBROUTINE FFTIN (KUSER1)
      COMMON /SYSTA/ IFILE(20), IFSTAT(20), ITIME(4),   KEYS(28),
     *               NFNUM,     NLIT,       NCOLN(32),  NCOLL(32),
     *               NFDOT(32), NFDOL(32),  NLUSER(32), FNUM(32),
     *               SWITCH(28)
      LOGICAL        SWITCH
      COMMON /SYSTB/ PROGNM,    PROSNM,     CCODE,      TITLE,
     *               CHIN,      LIT(32),    CHOUT
      CHARACTER      PROGNM *8, PROSNM *6,  CCODE *6,   TITLE *64,
     *               CHIN *80,  LIT *6,     CHOUT *72
      EQUIVALENCE (ICRYS, IFILE(3)), (ICOND, IFILE(4))
      EQUIVALENCE (LIS1, IFILE(7)), (LIS2, IFILE(8))
      EQUIVALENCE (IBINFF, IFILE(16)), (IFMAP,  IFILE(17))
      EQUIVALENCE (ISCRA,  IFILE(18))
      EQUIVALENCE (KEYS(25), KEYDS)
      EQUIVALENCE (KEYS(27), IMAP), (KEYS(28), IHALF)
      LOGICAL      SWPRI, PRIMAP, SWNOCY
      EQUIVALENCE (SWNOCY, SWITCH(18))
      EQUIVALENCE (SWPRI, SWITCH(10)), (PRIMAP, SWITCH(11))
      COMMON /CRYSA/ CELL(6),  CELLSD(6),   RCELL(6),   VOLUM,
     *               WAVE,     CELALL(10),  AMOLW,      ZET,
     *               NELEC,    F000,        ABSMU,      ICENT,
     *               ILATT,    ISYST,       ILAUE,      IMULT,
     *               IUNIQ,    IPOLA,       NTYPE,      NSYMM,
     *         IRSYMM(3,3,24), TSYMM(3,24), NLATT,      TLATT(3,4),
     *         FRAC2C(3,3),    CART2F(3,3), RRMAT(3,3), SSMAT(3,3)
      COMMON /CRYSB/ SPGR,     WAVEAT,      CELATY(10)
      CHARACTER      SPGR *16, WAVEAT *2,   CELATY *2
      COMMON /FFTDA/ SCALE, MH(3), NPP(3), XLMIN(3), XLMAX(3)
      COMMON /SEARDA/ D2R, DMPIC, DMAXB, DMOUT, DMINB, ANGM(2), MCON,
     *        SEARDX, NPIC, NATIN, NAT, NATS, NATX, NATSN, BOV, IPRY,
     *        NRECY, NRECYR, NRECYS, PSQ, NATREC, KPROG, SCALEX, R2X
      PARAMETER (MAXBUF = 198)
      DIMENSION FITFFT(5), BUFFFT(MAXBUF), IGM3(3)
      DIMENSION MAXHKL(3)
      PARAMETER (LCMAX = 16)
      CHARACTER * 6 LCONDA(LCMAX)
      DATA LCONDA / 'FOUR',   'PRIMAP', 'GRID' , 'MAXXYZ', 'MINXYZ',
     *              'MAXHKL', 'GRIDMO', 'PEAKS', 'PROJEC', 'DOUT',
     *              'DMIN',   'DMAX',   'AMIN',  'AMAX',   'dummy$',
     *              'NORECY'/
      R2X = 999.
      IACTOR  = 0
      FACTOR  = 0.25
      CALL KERNZA (0.0, XLMIN, 3)
      CALL KERNZA (1.0, XLMAX, 3)
      IGM3(1) = 1
      IGM3(2) = 1
      IGM3(3) = 2
      SCALEX = 0.
      CALL KERNZI (999, MH, 3)
      NPIC = 0
      NATIN = 0
      DMOUT = -1.
      DMINB = -1.
      DMAXB = -1.
      ANGM(1) = -1.
      ANGM(2) = -1.
      NAT = 0
      NATS = 0
      NRECY = 0
      NATREC = 0
      NATX = 0
      NATSN = 0
  85  CALL RDCOND (ICOND, LCONDA, LCMAX , KEND)
      IF (KEND.LE.0) GOTO 101
      GOTO (1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 85, 16), KEND
  1   IF (LIT(3) .EQ. 'NRECY') NRECY = NINT(FNUM(1))
      KEYS(26) = NRECY
      GOTO 85
  2   PRIMAP = .TRUE.
      GOTO 85
  3   FACTOR = FNUM(1)
      IACTOR  = 1
      GOTO 85
  4   XLMAX(1) = AMIN1 (FNUM(1), 1.0)
      XLMAX(2) = AMIN1 (FNUM(2), 1.0)
      XLMAX(3) = AMIN1 (FNUM(3), 1.0)
      GOTO 85
  5   XLMIN(1) = AMAX1 (FNUM(1), 0.0)
      XLMIN(2) = AMAX1 (FNUM(2), 0.0)
      XLMIN(3) = AMAX1 (FNUM(3), 0.0)
      GOTO 85
  6   DO 706 I=1,3
  706 MH(I) = NINT (FNUM(I))
      GOTO 85
  7   DO 707 I=1,3
      IF (FNUM(I).GE.1.0) IGM3(I) = NINT (FNUM(I))
      IF (IGM3(I).EQ.5 .OR. IGM3(I).GE.7)
     *    CALL KERROR ('BAD GRIDMO CARD', 6, 'FFTIN')
  707 CONTINUE
      GOTO 85
  8   NPIC = IFIX(FNUM(1))
      GOTO 85
  9   NATIN = IFIX(FNUM(1))
      GOTO 85
  10  DMOUT = FNUM(1)
      GOTO 85
  11  DMINB = FNUM(1)
      GOTO 85
  12  DMAXB = FNUM(1)
      GOTO 85
  13  ANGM(1) = FNUM(1)
      GOTO 85
  14  ANGM(2) = FNUM(1)
      GOTO 85
  16  SWNOCY = .TRUE.
      GOTO 85
  101 CALL FILCLO (ICOND, 'KEEP')
      NRECYR = NRECY / 10
      NRECYS = MOD(NRECY,10)
      CALL BINIFF (1, IBINFF, 'BINFFT', FITFFT, NITFFT, BUFFFT, NEND)
      PSQ = BUFFFT(27)
      IMAP = NINT(BUFFFT(28))
      KEYDS = NINT(BUFFFT(29))
      IF (SWPRI) WRITE (LIS1, FMT='('' Option IMAP = '', I2)') IMAP
      R2X = BUFFFT(30)
      IF (R2X .GT. 0.001)
     *   WRITE (LIS2, FMT='(48X, '' input atoms: R2 ='', F6.3)') R2X
      SCALEX = BUFFFT(31)
      WRITE (LIS2, FMT='(48X, '' BINFFT SCALE= '', F9.5)') SCALEX
      CALL FILINQ (ISCRA, 'BINBIG', 'UNFORMATTED', 'OUTPUT', KINQ)
      IF (IMAP .NE. 5) THEN
         CALL FILINQ (IFMAP, 'FMAP',  'UNFORMATTED', 'OUTPUT', KINQ)
      ELSE
         CALL FILINQ (IFMAP, 'FMAPT', 'UNFORMATTED', 'OUTPUT', KINQ)
         ENDIF
      CALL KERF2I (BUFFFT(22), MAXHKL, 3)
      IF (MAXHKL(1) .EQ. 0) CALL KERF2I (BUFFFT(7), MAXHKL, 3)
      IF (IMAP .EQ. 5) CALL KERNAI (MH, MAXHKL, 3)
      DO 105 I=1,3
      MH(I) = MIN0 (MH(I), MAXHKL(I))
  105 IF (MH(I) .LE. 0) MH(I) = 1
      CALL RDCRYS (ICRYS)
      CALL FILCLO (ICRYS, 'KEEP')
      IHALF = 0
      GOTO (110, 120, 130, 130, 215, 120), IMAP
      CALL KERNER(-3, 'FFTIN')
  110 WRITE (CHOUT,111) CCODE
  111 FORMAT (' Fourier in space group P1 for compound ', A6)
      CALL SHOUT2
      IF (ILATT.NE.1) WRITE (CHOUT, 112)
  112 FORMAT ('+', 47X, 'in non-promitive setting')
      CALL SHOUT2
      NSYMM = 1
      IMULT = NLATT
      ICENT = 1
      ILAUE = 1
      SWNOCY = .TRUE.
      GOTO 190
  120 ICENT = 2
      IMULT = NSYMM * ICENT * NLATT
      CALL KERNZA(0.0, TSYMM, 72)
      IF (IMAP .EQ. 2) WRITE (CHOUT, 121)
  121 FORMAT (' PATOR: sharpened Patterson for program ORIENT')
      IF (IMAP .EQ. 6) WRITE (CHOUT, 122)
  122 FORMAT (' Sharpened Patterson for program PATTY')
      CALL SHOUT2
      IF (ILAUE.EQ.2 .AND. IUNIQ.EQ.3) GOTO 220
      IF (ILAUE.EQ.1 .OR. ILAUE.EQ.4) GOTO 220
      IF (ILAUE.GE.6 .AND. ILAUE.LE.12) GOTO 220
      IHALF = -1
      GOTO 220
  130 IF (ICENT.EQ.2) THEN
         IHALF = -1
         GOTO 220
         ENDIF
  190 GOTO (192,191,192,191,191,191,191), ILATT
  191 IHALF = 1
      GOTO 220
  192 IF (IMAP.EQ.1) GOTO 220
      DO 197 II=1,NSYMM
      IF (IRSYMM(2,2,II).EQ.-1 .AND. TSYMM(2,II).LT.0.01) IHALF = -1
  197 IF (IRSYMM(2,2,II).EQ.1 .AND. (ABS(TSYMM(2,II)-0.5)).LT.0.01)
     *   IHALF = 1
      GOTO 220
  215 CONTINUE
      CALL KERROR (' IMAP=5 not aacepted!', 215, 'FFTIN')
  220 IF (IACTOR .EQ.1) GOTO 230
      IF (IMAP .NE. 2 .AND. IMAP .NE. 6) GOTO 223
      IF (VOLUM .LE. 4000.) GOTO 223
      RESOL = ALOG(VOLUM/4000.) / 3.
      RESOL = 0.3 * EXP(RESOL)
      IF (RESOL.LT.FACTOR) GOTO 223
      FACTOR = RESOL
      WRITE (LIS2, 222) FACTOR
  222 FORMAT (' The GRID spacing is approximately', F6.3, ' Angstrom')
  223 RESOL = AMIN1 (CELL(1) / FLOAT(MH(1)),
     +               CELL(2) / FLOAT(MH(2)),
     +               CELL(3) / FLOAT(MH(3))) * 0.25
      IF (RESOL .GE. FACTOR) THEN
         IF (RESOL .LT. FACTOR + 0.06) RESOL = FACTOR + 0.01
         IF (RESOL .GE. FACTOR + 0.06) RESOL = RESOL - 0.05
         FACTOR = RESOL
         WRITE (LIS1, 222) FACTOR
         ENDIF
  230 IF (IHALF.NE.0 .AND. IGM3(3).EQ.1) IGM3(2) = 2
  240 DO 241 I=1,3
  241 NPP(I) = CELL(I) / FACTOR + 0.5
      DO 280 I=1,3
      ISGG = MOD (NPP(I), IGM3(I))
      IF (ISGG.NE.0) NPP(I) = NPP(I) + IGM3(I) - ISGG
  250 NTEST = NPP(I)
      DO 270 J=2,5
  260 IF (NTEST.NE.(NTEST/J)*J) GOTO 270
      NTEST = NTEST / J
      IF (NTEST.EQ.1) GOTO 280
      GOTO 260
  270 CONTINUE
      NPP(I) = NPP(I) + IGM3(I)
      GOTO 250
  280 CONTINUE
      WRITE (LIS2, 222) FACTOR
      IF (NPP(1) .LE. 250) GOTO 400
      WRITE (LIS2, 320)
  320 FORMAT (' NX GREATER THAN 250 (SEE SUBR. -OUTPUT-). RESET.'/)
      FACTOR = FACTOR * FLOAT(NPP(1)) / 245.
      GOTO 240
  400 I = (NPP(1)+2) * (NPP(3)+2)
      IF (I.LT.KUSER1) GOTO 406
      FACTOR = FACTOR * 1.02 * SQRT(FLOAT(I)/FLOAT(KUSER1))
      WRITE (LIS1, 405)
      WRITE (LIS2, 405)
  405 FORMAT (' TOO MANY GRID POINTS FOR PEAK SEARCH. RESET.'/)
      GOTO 240
  406 WRITE (LIS2, 407) MH, NPP
  407 FORMAT (
     + ' Maximum indices allowed   h:', I4, '    k:', I4, '    l:', I4/
     + ' Number of grid points    Nx:', I4, '   Ny:', I4, '   Nz:', I4)
      DO 408 I = 1,3
      IF (NPP(I) .GE. 2 * MH(I) + 2) GOTO 408
      MH(I) = NPP(I) / 2 -1
      WRITE (LIS2, FMT='('' Reset MAXHKL:'')')
      GOTO 406
  408 CONTINUE
      IF (SWPRI .AND. PRIMAP) WRITE (LIS2, 410) XLMAX
  410 FORMAT (' FOURIER MAP TO BE PRINTED FROM:'/' X =  0.0  TO',
     +         F7.3,',  Y =  0.0  TO',F7.3,',  Z =  0.0  TO',F7.3)
      RETURN
      END
      SUBROUTINE SEARCH
      COMMON /SYSTA/ IFILE(20), IFSTAT(20), ITIME(4),   KEYS(28),
     *               NFNUM,     NLIT,       NCOLN(32),  NCOLL(32),
     *               NFDOT(32), NFDOL(32),  NLUSER(32), FNUM(32),
     *               SWITCH(28)
      LOGICAL        SWITCH
      COMMON /SYSTB/ PROGNM,    PROSNM,     CCODE,      TITLE,
     *               CHIN,      LIT(32),    CHOUT
      CHARACTER      PROGNM *8, PROSNM *6,  CCODE *6,   TITLE *64,
     *               CHIN *80,  LIT *6,     CHOUT *72
      EQUIVALENCE (IATOMS, IFILE(1))
      EQUIVALENCE (IPR1,   IFILE(6)), (LIS1, IFILE(7)), (LIS2, IFILE(8))
      EQUIVALENCE (IFMAP, IFILE(17))
      EQUIVALENCE (KEYS(27), IMAP), (KEYS(28), IHALF)
      LOGICAL SWRECY, DMAXCH
      EQUIVALENCE (SWITCH(17), SWRECY), (SWITCH(28), DMAXCH)
      CHARACTER SYMB(10) *2
      EQUIVALENCE (CHIN, SYMB(1))
      COMMON /CRYSA/ CELL(6),  CELLSD(6),   RCELL(6),   VOLUM,
     *               WAVE,     CELALL(10),  AMOLW,      ZET,
     *               NELEC,    F000,        ABSMU,      ICENT,
     *               ILATT,    ISYST,       ILAUE,      IMULT,
     *               IUNIQ,    IPOLA,       NTYPE,      NSYMM,
     *         IRSYMM(3,3,24), TSYMM(3,24), NLATT,      TLATT(3,4),
     *         FRAC2C(3,3),    CART2F(3,3), RRMAT(3,3), SSMAT(3,3)
      COMMON /CRYSB/ SPGR,     WAVEAT,      CELATY(10)
      CHARACTER      SPGR *16, WAVEAT *2,   CELATY *2
      PARAMETER (MAXAT=993)
      COMMON /XATXYZ/ X(4,MAXAT), ATXYZ(10,MAXAT), IZAT(MAXAT)
      COMMON /ATNAMA/ ATNAME(MAXAT)
      CHARACTER ATNAME *6
      COMMON / / IFRAG(MAXAT), ISYM(MAXAT), IDUM(MAXAT),
     +           DUM(MAXAT),  IBOND(MAXAT*10), JBOND(MAXAT*10),
     +           XGEO
      DIMENSION BLACOM(42000)
      EQUIVALENCE (BLACOM(1), IFRAG(1))
      COMMON /SEARDA/ D2R, DMPIC, DMAXB, DMOUT, DMINB, ANGM(2), MCON,
     *        SEARDX, NPIC, NATIN, NAT, NATS, NATX, NATSN, BOV, IPRY,
     *        NRECY, NRECYR, NRECYS, PSQ, NATREC, KPROG, SCALEX, R2X
      COMMON /SIZEX/ KFRAG(20), NFRAG, LFRAG(20), NOFRAG, NNA
      COMMON /DIRBFA/ NCELTY(10), NCELLZ(10), NCELIN(10), NCELIX(10)
      COMMON /DIRBFB/ ACELTY(10)
      CHARACTER ACELTY *2
      DIMENSION ITLE(20), XLOCK(3), CELPAR(10), IZTYPE(10), BUFFOX(10)
      IPRY = LIS1
      SCALAT = 0.
      NPROJ  = 2
      IF (NATIN .GT. 0) THEN
         IF (NATIN.EQ.1 .OR. NATIN.EQ.3) NPROJ = NATIN
         NATIN = 0
         ENDIF
      NPC = MIN0 (NPIC, MAXAT-50)
      IF (DMAXB .LT. 0.) THEN
         DMAXB  = 1.95
      ELSE
         IF (DMAXB .GT. 1.95) DMAXCH = .TRUE.
         ENDIF
      IF (ANGM(1) .LT. 0.) ANGM(1) = 80.0
      IF (ANGM(2) .LT. 0.) ANGM(2) = 145.0
      IF (DMINB .LT. 0.) DMINB = 0.90
      IF (DMOUT .LT. 0.) DMOUT = 2.40
      DMPIC  = 0.85
      REWIND IFMAP
      READ (IFMAP) ITLE, IMAP, IHALF
      WRITE (LIS2, FMT='('' OPTION IMAP = '', I2)') IMAP
      IF (IMAP .LE. 0 .OR. IMAP .EQ. 5 .OR. IMAP.GE.7) CALL KERROR
     *  ('Error reading output Fourier map for search', 0, 'SEARCH')
      IF (IMAP .EQ. 2) GOTO 120
      IF (IMAP .NE. 6) THEN
         CALL FILINQ (IATOMS, 'ATOMS', 'FORMATTED', 'INPUT', KINQ)
         IF (KINQ.NE.0) CALL KERROR ('No atoms file found', 0, 'SEARCH')
      ELSE
         CALL FILINQ (IATOMS, 'ATOMS', 'FORMATTED', 'OUTPUT', KINQ)
         ENDIF
  120 IF (IMAP.EQ.2 .OR. IMAP.EQ.6) GOTO 190
      CALL ATOMIN (IATOMS, ATXYZ, ATNAME, IZAT, MAXAT, NATIN, KEYT)
      REWIND IATOMS
      IF (NFNUM .GT. 0) THEN
         IF (LIT(NLIT). EQ. 'SC=' .AND.
     *   FNUM(NFNUM) .GT. 0.0001) SCALAT = FNUM(NFNUM)
         ENDIF
      IF (SCALEX .GT. 0.0001) SCALAT = SCALEX
      NATQ = NATIN
      NATH = 0
      N = 1
  131 CONTINUE
      IF (NRECYR.NE.0 .AND.ATNAME(N)(1:1).EQ.'H'.AND.IZAT(N).EQ.1) THEN
         NATH = NATH + 1
         ATNAME(N)(1:1) = 'Q'
         ENDIF
      IF (ATNAME(N)(1:1) .EQ. 'Q') THEN
         IF (N .EQ. NATIN) GOTO 135
         DO 133 N1 = N, NATIN - 1
         CALL KERNAB (ATXYZ(1,N1+1), ATXYZ(1,N1), 10)
         ATNAME(N1) = ATNAME(N1+1)
  133    IZAT(N1) = IZAT(N1+1)
  135    NATIN = NATIN - 1
         N = N - 1
         ENDIF
      N = N + 1
      IF (N .LE. NATIN) GOTO 131
      IF (NATIN.LT.NATQ) THEN
         NATQ = NATQ - NATIN - NATH
         IF (NATH .GT. 0) WRITE (LIS2, FMT=
     *      '('' Nr of H atoms rejected:'', I3)') NATH
         IF (NATH .GT. 0) WRITE (LIS2, FMT=
     *      '('' Nr of Q-atoms (= peaks) rejected:'', I3)') NATQ
         IF (NATIN.LE.0) CALL KERROR ('No atoms left!', 135, 'SEARCH')
         WRITE (LIS2, FMT='(1X)')
         ENDIF
      IF (NRECYR .LE. 1) GOTO 7167
      CALL KERNZI (0, IZTYPE, 10)
      DO 7157 J=1,NTYPE
      CALL ATOMIZ (CELATY(J), NLET, IZ)
      IZTYPE(J) = IZ
 7157 CONTINUE
      CALL KERNZA (0.0, CELPAR, 10)
      AAMULT = FLOAT(IMULT)
      DO 7161 I=1,NATIN
      DO 7160 J=1,NTYPE
      IF (IZAT(I).NE.IZTYPE(J)) GOTO 7160
      CELPAR(J) = CELPAR(J) + ATXYZ(4,I) * AAMULT
 7160 CONTINUE
 7161 CONTINUE
      IIII = 0
      DO 7162 J=1,NTYPE
      IF ( ( IZTYPE(J) .GE. 10 .AND. CELPAR(J) .GT. CELALL(J) ) .OR.
     *     ( IZTYPE(J) .GE. 2 .AND. CELPAR(J) .GT. CELALL(J) .AND.
     *      NRECYR .GE. 7) ) THEN
         CELALL(J) = CELPAR(J)
         IIII = 1
         ENDIF
 7162 CONTINUE
      IF (IIII .EQ. 0) GOTO 7167
      DO 7165 I=1,NTYPE
 7165 BUFFOX(I) = CELALL(I) / ZET
      I = NTYPE
      J = NINT(ZET)
      WRITE (LIS1, 7166) J, (CELATY(K), BUFFOX(K), K=1,I)
 7166 FORMAT (/' NOTE: Cell Contents reset [ output DDMAIN! ] :'/
     *  ' Z:', I3 / ' FORMUL:', 6(2X,A2,F6.1) /
     *                           ( 8X, 6(2X,A2,F6.1))/)
 7167 CONTINUE
      CALL CELZAT (ACELTY, NCELTY, NCELLZ)
      CALL CELZIN (ATXYZ, IZAT, NATIN, NCELLZ, NCELIN)
      WRITE (LIS2, FMT='('' = symmetry included ='')')
      CALL KERNZI (0, NCELIX, 10)
      DO 137 I = 1, NTYPE
      NCELIX(I) = NCELTY(I)
  137 CONTINUE
      WRITE (LIS1, 139) NATIN
  139 FORMAT (/' NUMBER OF ATOMS INPUT FROM FILE IS',I6)
      NATTR = 0
      NATTH = 0
      NATTHA= 0
      NATR = 0
      NATH = 0
      DO 140 I=1,NTYPE
      IF (NCELLZ(I) .NE. 1) THEN
         NATTR = NATTR + NCELTY(I)
         NATR = NATR + NCELIN(I)
         IF (NCELLZ(I) .GT. 40) NATTHA = NATTHA + NCELTY(I)
      ELSE
         NATTH = NATTH + NCELTY(I)
         NATH = NATH + NCELIN(I)
         ENDIF
  140 CONTINUE
      IF (NATTHA .GE. 1) NATTHA = NATTR / NATTHA
      IF (NATH .GT. 0) THEN
         NATR = NATR + NATH
         NATTR = NATTR + NATTH
         ENDIF
      NATX = MAX0(0, (NATTR - NATR)/IMULT)
      IF (NRECY .GT. 0) THEN
         NRECY = NRECY + 11
         NRECYR = NRECYR + 1
         NRECYS = NRECYS + 1
         WRITE (IPR1, 144) NRECYR
         WRITE (LIS1, 144) NRECYR
         WRITE (LIS2, 144) NRECYR
  144    FORMAT (56X, 'Fourier cycle ', I2)
         IF (NRECYR .GE. 3) THEN
            DO 147 I = 1, NTYPE
            NCELIX(I) = MAX0 (NCELTY(I), NCELIN(I))
  147       CONTINUE
            ENDIF
         ENDIF
      WRITE (LIS2, 154) NATIN, NATX
  154 FORMAT (/' NUMBER OF ATOMS INPUT FROM FILE IS',I6/
     + ' NUMBER OF NEW ATOMS TO BE FOUND IS',I6, ' (if on genl posn)')
      NATXX = NATX + NATIN
      IF (DMAXB .GT. DMOUT) DMOUT = DMAXB
      IF (DMAXB .LT. DMINB) THEN
         DMINB   = 0.
         ANGM(1) = 0.
         ANGM(2) = 180.
         ENDIF
      IF (NPC .GT. 0) THEN
         NPIC = NPC
         NATX = NPC
         WRITE (LIS2, 163) NATX
  163    FORMAT (' NUMBER OF ATOMS (PEAKS) TO BE CONSIDERED IS', I5)
      ELSE
         NATX = MIN0((20*NATXX + 5)/19, MAXAT-50) + 2
         NPIC = MIN0((15*NATXX + 6)/14, MAXAT-50) + 1
         NPIC = MAX0(NPIC, MIN0(20, NATTR) +3)
         IF (NATTHA .GT. 0)
     *      NPIC = NPIC + MIN0 (7 * (NATIN + NATX) / NATTHA, NATX/2)
         IF (NPIC .EQ. NATX) NPIC = NPIC + 1
         WRITE (LIS2, 172) NATX, NPIC
  172    FORMAT (' NUMBER OF ATOMS TO BE CONSIDERED IS', I5/
     +           ' NUMBER OF PEAKS to be searched is ',I6)
            III = MAX0(150, NATIN * 3/2)
            IF (IMAP .EQ. 1 .AND. NATX .GT. III) THEN
            WRITE (LIS2, 174) III
  174       FORMAT (' because of DIRP1 this is reduced to:', I4)
            NPIC = III + 1
            NATX = III
            ENDIF
         ENDIF
      WRITE (LIS2, 184) DMINB, DMAXB, ANGM
  184 FORMAT (' STEREOCHEMICAL CRITERIA'/ ' For molecular clusters:',
     *   15X, 'MINIMUM BONDING DISTANCE  =', F6.2 /
     *  ' defaults, or',  26X, 'MAXIMUM BONDING DISTANCE  =', F6.2 /
     *  ' user-defined;', 31X, 'MINIMUM BOND ANGLE  =', F6.1 /
     *  ' may be modified later.', 22X, 'MAXIMUM BOND ANGLE  =', F6.1 )
      IF (NPROJ .NE. 2) WRITE (LIS2, 186) NPROJ
  186 FORMAT (9X,'NUMBER OF PROJECTIONS OF EACH CLUSTER TO BE OUTPUT',
     +' IS',I3)
      GOTO 200
  190 IF (NPC .LE. 0) NPIC = 60
      NATX = NPIC
      NATIN = 0
  200 D2R = ATAN(1.0) / 45.0
      CALL PKSRCH (X, MAXAT, IHALF, IFMAP)
      DMAXLI = 0.5
      IF (IMAP .NE. 3) DMAXLI = 0.1
      DO 225 I = 1, NPIC
      CALL LOCKIN ( X(1,I), DMAXLI, XLOCK, DISTLI, NPOSLI )
      IF (DISTLI .GT. 0.15)
     *   WRITE (LIS1,223) I, DISTLI, (X(J,I),J=1,3), XLOCK
  223 FORMAT (' Atom ',I3,' locked in: peak shifted over', F6.2,
     *   ' Angstrom' / ' Peak xyz:',3F9.5, '  LOCKED xyz:',3F9.5)
  225 IF (NPOSLI .GT. 1) CALL KERNAB (XLOCK, X(1,I), 3)
      IF (IMAP .EQ. 2 .OR. IMAP .EQ. 6) THEN
         CALL FILCLO (IFMAP, 'KEEP')
      ELSE
         CALL FILCLO (IFMAP, 'DELETE')
         ENDIF
      CALL DIRBON (IMAP, ZSCAL)
      IDDOKA = KEYS(10)
      IF (IDDOKA .EQ. 17) RETURN
      IF (NRECYR .GE. 3) SCALAT = SCALEX
      NNA = 0
      CALL CLSTRS (LIS2)
      IF (MCON.GT.0) GOTO 350
      WRITE (LIS1, 227) ((X(J,I),J=1,4),I=1,NAT)
      WRITE (LIS2, 227) ((X(J,I),J=1,4),I=1,NAT)
  227 FORMAT (' NO BONDS FOUND'//
     + ' ATOMIC POSITIONS'//46X,'X',9X,'Y',9X,'Z',3X,
     + 'HEIGHT'/(41X,3F10.4,F8.0))
      CALL FILINQ (IATOMS, 'ATOMS', 'FORMATTED', 'OUTPUT', KINQ)
      IF (SCALAT .GT. 0.01) FNUM(32) = SCALAT
      CHOUT = ' Uninterpreted peaks, denoted C-atoms'
      CALL ATOMWA (IATOMS)
      DO 344 I=1,NAT
  344 WRITE (IATOMS, 341) I, (X(J,I),J=1,3)
  341 FORMAT ('ATOM   C', I3, 2X, 3F8.5)
      WRITE (IATOMS, FMT = '(''END'')')
      CALL FILCLO(IATOMS, 'KEEP')
      CHOUT =
     *' Warning: No atoms interpretation; No recycling: check file!'
      CALL SHOUT
      CALL KERROR
     * ('Savety stop: peaks output to ATOMS file.', 0, 'SEARCH')
  350 CONTINUE
      DO 600 NOFRG=1,NFRAG
      NOFRAG = NOFRG
      IF (NOFRAG .EQ. 20) GOTO 600
      IF (KFRAG(NOFRAG) .LT. 4) GOTO 600
      CALL PICTUR (NPROJ)
  600 CONTINUE
      CALL SCHOUT (KEYT, SCALAT, ZSCAL)
      IF (NRECYR .LE. 0) RETURN
      WRITE (CHOUT, FMT='('' Cycle'',I2,'' is finished'')')  NRECYR
      CALL SHOUT2
      IF (.NOT.SWRECY) THEN
         CHOUT = ' Fourier recycling procedure completed'
         CALL SHOUT2
         IF (NRECYR .GT. 1) THEN
            CHOUT =
     *    ' ATOMS from the each cycle are appended to the ATOLD file.'
            CALL SHOUT2
            ENDIF
         CHOUT = ' Final atomic parameters written to the ATOMS file,'
         CALL SHOUT
         CHOUT = ' and appended also to the ATOLD (= back-up) file.'
         CALL SHOUT
         ENDIF
      RETURN
      END
      SUBROUTINE PKSRCH (X, MAXAT, IHALF, LIN)
      DIMENSION  X(4, MAXAT)
      COMMON /CRYSA/ CELL(6),  CELLSD(6),   RCELL(6),   VOLUM,
     *               WAVE,     CELALL(10),  AMOLW,      ZET,
     *               NELEC,    F000,        ABSMU,      ICENT,
     *               ILATT,    ISYST,       ILAUE,      IMULT,
     *               IUNIQ,    IPOLA,       NTYPE,      NSYMM,
     *         IRSYMM(3,3,24), TSYMM(3,24), NLATT,      TLATT(3,4),
     *         FRAC2C(3,3),    CART2F(3,3), RRMAT(3,3), SSMAT(3,3)
      COMMON /SEARDA/ D2R, DMPIC, DMAXB, DMOUT, DMINB, ANGM(2), MCON,
     *        SEARDX, NPIC, NATIN, NAT, NATS, NATX, NATSN, BOV, IPRY,
     *        NRECY, NRECYR, NRECYS, PSQ, NATREC, KPROG, SCALEX, R2X
      PARAMETER (KUSER2=30000)
      COMMON / / NR3D
      INTEGER*2 NR3D(KUSER2)
      INTEGER*2 BLACOM(84000)
      EQUIVALENCE (BLACOM(1), NR3D(1))
      DIMENSION XS(3), X1(3), IDIFF(19), B(19)
      DIMENSION DXYZM(3)
      DIMENSION ITLE(20)
      EQUIVALENCE (SCALE, ITLE(18))
      DATA MAX, E, D / 0, 0.0, 0.0 /
      DO 101 I = 1, 3
  101 DXYZM(I) = DMPIC * RCELL(I)
      READ (LIN) NNX, NNZ, NNY, NNYT
      NNYOLD=NNY
      IF (IHALF.NE.0) NNY=NNY-3
      NNXP2 = NNX + 2
      NXZ = NNXP2 * (NNZ + 2)
      NXZ3 = 3 * NXZ
      IF (NXZ3 .GT. KUSER2) CALL KERNER (940, 'PKSRCH')
      DX = 1.0 / FLOAT(NNX)
      DY = 1.0 / FLOAT(NNYT)
      DZ = 1.0 / FLOAT(NNZ)
      LEVEL = 0
      LIMIT = MIN0(MAXAT, 2*NATX)
 1100 IDIFF(1) = -NXZ - 1
      IDIFF(2) = -NXZ - NNXP2
      IDIFF(3) = -NXZ
      IDIFF(4) = -NXZ + NNXP2
      IDIFF(5) = -NXZ + 1
      IDIFF(6) = -NNXP2 - 1
      IDIFF(7) = -1
      IDIFF(8) = NNXP2 - 1
      IDIFF(9) = -NNXP2
      IDIFF(10) = 0
      DO 1120 I=1,9
      J=20-I
      IDIFF(J) = -IDIFF(I)
 1120 CONTINUE
      NO = 0
      IY = -1
      NY = 0
 1200 REWIND LIN
      READ(LIN) ITLE, IMAP
      READ(LIN)
      IF (IMAP.EQ.1 .OR. IMAP.EQ.3 .OR. IMAP.EQ.4) THEN
         RESCAL = 4000. / VOLUM / SCALE
      ELSE
         RESCAL = 1.
         IMAP = 0
         ENDIF
      IF (IY+2.EQ.NNYOLD) GOTO 1400
      MAX=NXZ
      ISKIP = (NNYOLD-1) * NNZ
      DO 1305 I=1,ISKIP
 1305 READ (LIN)
      CALL RDSECT (MAX, NNXP2, NNZ, NXZ3, LIN)
      REWIND LIN
      READ(LIN)
      READ(LIN)
      CALL RDSECT (MAX, NNXP2, NNZ, NXZ3, LIN)
 1400 MX = MAX - NXZ + NNX + 1
      CALL RDSECT (MAX, NNXP2, NNZ, NXZ3, LIN)
      IY = IY + 1
      NY = MOD(NY+2, 3) - 1
      KK = NXZ3
      IF (NY) 1440, 1460, 1500
 1440 KK = -NXZ3
 1460 DO 1480 I=1,5
      IDIFF(I) = IDIFF(I) - KK
 1480 CONTINUE
      IF (NY .EQ. 0) GO TO 1540
 1500 DO 1520 I=15,19
      IDIFF(I) = IDIFF(I) - KK
 1520 CONTINUE
 1540 DO 2000 IZ=1,NNZ
      MN = MX + 3
      MX = MX + NNXP2
      DO 1980 IX=MN,MX
      IF (NR3D(IX) .LT. LEVEL) GO TO 1980
      DO 1560 I=1,9
      J = IDIFF(I) + IX
      IF (NR3D(IX) .LE. NR3D(J)) GO TO 1980
 1560 CONTINUE
      DO 1580 I=11,19
      J = IDIFF(I) + IX
      IF (NR3D(IX) .LT. NR3D(J)) GO TO 1980
 1580 CONTINUE
      DO 1600 I=1,19
      J = IDIFF(I) + IX
      B(I) = NR3D(J)
 1600 CONTINUE
      B1 = B(3) + B(7) + B(9) + B(11) + B(13) + B(17)
      B2 = B(1) + B(2) + B(4) + B(5) + B(6) + B(8) + B(12) + B(14) +
     +  B(15) + B(16) + B(18) + B(19)
      F = (30.0 * B(10) + 11.0 * B1 - 8.0 * B2) / 63.0
      C = (B(5)+B(12)+B(13)+B(14)+B(19)-B(1)-B(6)-B(7)-B(8)-B(15))/10.0
      DELTAX = C / F
      IF (ABS(DELTAX) .GT. 1.0) GO TO 1620
      D = (B(15)+B(16)+B(17)+B(18)+B(19)-B(1)-B(2)-B(3)-B(4)-B(5))/10.0
      DELTAY = D / F
      IF (ABS(DELTAY) .GT. 1.0) GO TO 1620
      E = (B(4)+B(8)+B(11)+B(14)+B(18)-B(2)-B(6)-B(9)-B(12)-B(16))/10.0
      DELTAZ = E / F
      IF (ABS(DELTAZ) .LE. 1.0) GO TO 1640
 1620 DELTAX = 0.0
      DELTAY = 0.0
      DELTAZ = 0.0
 1640 XX = (FLOAT(IX-MN+1) + DELTAX) * DX
      YY = (FLOAT(IY) + DELTAY) * DY
      ZZ = (FLOAT(IZ) + DELTAZ) * DZ
      A = (9.0 * B(10) + 4.0 * B1 - B2) / 21.0
      BINT = A + 0.5 * (C * DELTAX + D * DELTAY + E * DELTAZ)
      IF (BINT .GT. 1.05 * B(10)) BINT = 1.05 * B(10)
      B10 = B(10)
      B(10) =  AMAX1(B(10), BINT)
      IF (IMAP .EQ. 0) GOTO 1660
         IF (B(10) .LT. 1.) B(10) = 1.
         B(10) = SQRT (B(10) * RESCAL)
      B30 = 0.
      BX2 = B(3)  + B(2)  + B(4)  + B(1)  + B(5)
      BX3 = B10   + B(9)  + B(11) + B(7)  + B(13)
      BX4 = B(17) + B(16) + B(18) + B(15) + B(19)
      CALL EXTPOL (BX2, BX3, BX4, BX1, BX5)
      B30 = B30 + BX1 + BX5
      BX2 = B(9)  + B(6)  + B(12) + B(2)  + B(16)
      BX3 = B10   + B(7)  + B(13) + B(3)  + B(17)
      BX4 = B(11) + B(8)  + B(14) + B(4)  + B(18)
      CALL EXTPOL (BX2, BX3, BX4, BX1, BX5)
      B30 = B30 + BX1 + BX5
      BX2 = B(13) + B(12) + B(14) + B(5)  + B(19)
      BX3 = B10   + B(9)  + B(11) + B(3)  + B(17)
      BX4 = B(7)  + B(6)  + B(8)  + B(1)  + B(15)
      CALL EXTPOL (BX2, BX3, BX4, BX1, BX5)
      B30 = B30 + BX1 + BX5
      B8 = 0.
      BX2 = B(1)  + B(2)
      BX3 = B(3) * 2.
      BX4 = B(4)  + B(5)
      CALL EXTPOL (BX2, BX3, BX4, BX1, BX5)
      B8 = B8 + BX1 + BX5
      BX2 = B(1)  + B(4)
      BX4 = B(2)  + B(5)
      CALL EXTPOL (BX2, BX3, BX4, BX1, BX5)
      B8 = B8 + BX1 + BX5
      BX2 = B(15)  + B(16)
      BX3 = B(17) * 2.
      BX4 = B(18)  + B(19)
      CALL EXTPOL (BX2, BX3, BX4, BX1, BX5)
      B8 = B8 + BX1 + BX5
      BX2 = B(15)  + B(18)
      BX4 = B(16)  + B(19)
      CALL EXTPOL (BX2, BX3, BX4, BX1, BX5)
      B8 = B8 + BX1 + BX5
      BX2 = B(1)   + B(6)
      BX3 = B(7) * 2.
      BX4 = B(8)  + B(15)
      CALL EXTPOL (BX2, BX3, BX4, BX1, BX5)
      B8 = B8 + BX1 + BX5
      BX2 = B(1)  + B(8)
      BX4 = B(6)  + B(15)
      CALL EXTPOL (BX2, BX3, BX4, BX1, BX5)
      B8 = B8 + BX1 + BX5
      BX2 = B(5)  + B(12)
      BX3 = B(13) * 2.
      BX4 = B(14) + B(19)
      CALL EXTPOL (BX2, BX3, BX4, BX1, BX5)
      B8 = B8 + BX1 + BX5
      BX2 = B(5)  + B(14)
      BX4 = B(12) + B(19)
      CALL EXTPOL (BX2, BX3, BX4, BX1, BX5)
      B8 = B8 + BX1 + BX5
      BX2 = B(8)  + B(18)
      BX3 = B(11) * 2.
      BX4 = B(14) + B(4)
      CALL EXTPOL (BX2, BX3, BX4, BX1, BX5)
      B8 = B8 + BX1 + BX5
      BX2 = B(18) + B(14)
      BX4 = B(8) + B(4)
      CALL EXTPOL (BX2, BX3, BX4, BX1, BX5)
      B8 = B8 + BX1 + BX5
      BX2 = B(6)  + B(16)
      BX3 = B(9) * 2.
      BX4 = B(12) + B(2)
      CALL EXTPOL (BX2, BX3, BX4, BX1, BX5)
      B8 = B8 + BX1 + BX5
      BX2 = B(12) + B(16)
      BX4 = B(6) + B(2)
      CALL EXTPOL (BX2, BX3, BX4, BX1, BX5)
      B8 = B8 + BX1 + BX5
      B8 = B8 / 6.
      PX2 = B10
      PX3 = (B1+B2) / 18.
      PX4 = (B8+B30)/ 38.
      CALL EXTPOL (PX2, PX3, PX4, PX1, PX5)
      PX51  = PX5
      PX5 = AMIN1(PX51, 0.8 * PX4)
      PX5 = AMAX1(PX51, 0.2 * PX4)
      B58  = 58. * PX5
      BSUM = B10 + B1 + B2 + B8 + B30 + B58
      IF (BSUM .LT. 1.) BSUM = 1.
      BSUM = SQRT(BSUM)
      IF (BSUM .GT. 999.) BSUM = 999.
      I10 = NINT(B(10))
      B(10) = FLOAT(I10) + BSUM / 1000.
 1660 NOP1 = NO + 1
      IF(NOP1.GT.MAXAT) GOTO 1821
      X(1,NOP1) = XX
      X(2,NOP1) = YY
      X(3,NOP1) = ZZ
      X(4,NOP1) = B(10)
      IF (NO .EQ. 0) GO TO 1820
      IR=0
      DO 1800 K=1, IMULT
      CALL OPER1 (K, XS, X(1,NOP1))
      DO 1780 I=1,NO
      DO 1720 L=1,3
      X1(L) = X(L,I) - XS(L)
 1680 IF (ABS(X1(L)) .LE. 0.5) GO TO 1700
      X1(L) = X1(L) - SIGN(1.0, X1(L))
      GO TO 1680
 1700 IF (ABS(X1(L)) .GT. DXYZM(L)) GO TO 1780
 1720 CONTINUE
      IF (QUAD2 (X1, X1) .GT. DMPIC) GOTO 1780
      IF (IR.GT.0) X(4,IR)=0.0
      IR=0
      IF (B(10) .LE. X(4,I)) GOTO 1980
      X(1,I) = XX
      X(2,I) = YY
      X(3,I) = ZZ
      X(4,I) = B(10)
      IR=I
 1780 CONTINUE
 1800 CONTINUE
      IF(IR.GT.0) GO TO 1980
 1820 NO = NOP1
 1821 IF (NO .LT. LIMIT) GO TO 1980
      CALL SORT (X, MAXAT, NO, 4)
      NO = LIMIT
      NPIC = NO
      IF (IMAP .EQ. 0) THEN
         LEVEL = X(4,NPIC) + 0.5
      ELSE
         LEVEL = X(4,NPIC)**2 / RESCAL + 0.5
         ENDIF
 1980 CONTINUE
 2000 CONTINUE
      IF (IY .GE. NNY) GO TO 2100
      IF (IY - NNYOLD + 2) 1400, 1200, 1400
 2100 CONTINUE
      IF (IMAP .NE. 1 .AND. IMAP .NE. 3) GOTO 2200
      DO 2155 I = 1, NO
      I10 = X(4,I)
      ITOT = X(4,I) * 10000.
      ISUM = ITOT - 10000 * I10
      IF (I10 .GT. 999) I10 = 999
      X(4,I) = FLOAT(ISUM) + FLOAT(I10)/1000.
 2155 CONTINUE
 2200 CALL SORT (X, MAXAT, NO, 4)
      NNN = MIN0 (NO, NPIC)
      IF (NNN .EQ. NPIC) GOTO 3000
      LEVEL = LEVEL - 100
      IF(LEVEL.GE.(-200)) GO TO 1100
      NPIC=NO
 3000 CONTINUE
      RETURN
      END
      SUBROUTINE EXTPOL (BX2, BX3, BX4, BX1, BX5)
      A = BX3
      B = (BX4 - BX2) / 2.0
      C = BX4 - A - B
      BX1 = A - 2.0 * B + 4.0 * C
      BX5 = A + 2.0 * B + 4.0 * C
      IF (BX1 .GT. BX2) BX1 = BX2
      IF (BX1 .LT. 0.0) BX1 = 0.0
      IF (BX5 .GT. BX4) BX5 = BX4
      IF (BX5 .LT. 0.0) BX5 = 0.0
      RETURN
      END
      SUBROUTINE DIRBON (IMAP, ZSCAL)
      GOTO (1,   2,   1,   4,   5,   2), IMAP
  1   CALL DIRBF (ZSCAL)
      GOTO 100
  4   CALL DIRBD
  100 CALL DIRBB
      RETURN
  2   CALL DIRBP
  5   RETURN
      END
      SUBROUTINE DIRBF (ZSCAL)
      COMMON /SYSTA/ IFILE(20), IFSTAT(20), ITIME(4),   KEYS(28),
     *               NFNUM,     NLIT,       NCOLN(32),  NCOLL(32),
     *               NFDOT(32), NFDOL(32),  NLUSER(32), FNUM(32),
     *               SWITCH(28)
      LOGICAL        SWITCH
      COMMON /SYSTB/ PROGNM,    PROSNM,     CCODE,      TITLE,
     *               CHIN,      LIT(32),    CHOUT
      CHARACTER      PROGNM *8, PROSNM *6,  CCODE *6,   TITLE *64,
     *               CHIN *80,  LIT *6,     CHOUT *72
      EQUIVALENCE (IDDL, IFILE(1)), (IDDS, IFILE(2)), (ICRYS,IFILE(3))
      EQUIVALENCE (ICOND,IFILE(4))
      EQUIVALENCE (IPR1, IFILE(6)), (LIS1, IFILE(7)), (LIS2, IFILE(8))
      EQUIVALENCE (KEYS(25), KEYDS)
      LOGICAL SWRECY, SWNOCY, REN98, FREE98, NOFREE
      EQUIVALENCE (SWITCH(14), NOFREE)
      EQUIVALENCE (SWITCH(16), REN98), (SWITCH(15), FREE98)
      EQUIVALENCE (SWITCH(17), SWRECY), (SWITCH(18), SWNOCY)
      COMMON /CRYSA/ CELL(6),  CELLSD(6),   RCELL(6),   VOLUM,
     +               WAVE,     CELALL(10),  AMOLW,      ZET,
     +               NELEC,    F000,        ABSMU,      ICENT,
     +               ILATT,    ISYST,       ILAUE,      IMULT,
     +               IUNIQ,    IPOLA,       NTYPE,      NSYMM,
     +         IRSYMM(3,3,24), TSYMM(3,24), NLATT,      TLATT(3,4),
     +         FRAC2C(3,3),    CART2F(3,3), RRMAT(3,3), SSMAT(3,3)
      COMMON /CRYSB/ SPGR,     WAVEAT,      CELATY(10)
      CHARACTER      SPGR *16, WAVEAT *2,   CELATY *2
      PARAMETER (MAXAT=993)
      COMMON /XATXYZ/ X(4,MAXAT), ATXYZ(10,MAXAT), IZAT(MAXAT)
      COMMON /ATNAMA/ ATNAME(MAXAT)
      CHARACTER ATNAME *6
      COMMON / / IFRAG(MAXAT), ISYM(MAXAT), IDUM(MAXAT),
     +           DUM(MAXAT),  IBOND(MAXAT*10), JBOND(MAXAT*10),
     +           XGEO
      DIMENSION BLACOM(42000)
      EQUIVALENCE (BLACOM(1), IFRAG(1))
      DIMENSION LW(MAXAT), JCON(MAXAT), XDUM(MAXAT)
      EQUIVALENCE (LW(1),IFRAG(1)), (JCON(1),ISYM(1)),(XDUM(1),IDUM(1))
      COMMON /SEARDA/ D2R, DMPIC, DMAXB, DMOUT, DMINB, ANGM(2), MCON,
     *        SEARDX, NPIC, NATIN, NAT, NATS, NATX, NATSN, BOV, IPRY,
     *        NRECY, NRECYR, NRECYS, PSQ, NATREC, KPROG, SCALEX, R2X
      COMMON /DIRBFA/ NCELTY(10), NCELLZ(10), NCELIN(10), NCELIX(10)
      COMMON /DIRBFB/ ACELTY(10)
      CHARACTER ACELTY *2
      DIMENSION XS(3), XST(3), XSHIFT(3)
      CHARACTER ATNAM *6 , ATNAMX(MAXAT) *6
      DIMENSION RCELTY(10), WCELTY(10), TCELTY(10), LCELTY(10)
      DIMENSION SATLD(133), ISATL(133), JSATL(133)
      DIMENSION XBTEMP(MAXAT), IZTEMP(MAXAT)
      DIMENSION IOCC(1), IRECY(8)
      LOGICAL SWMATE
      DATA DPAR2, II, NRECYO / 0.0 , 0 , 0 /
      NTYPEZ = NTYPE
  101 IF (ACELTY(NTYPEZ) .EQ. ' ') THEN
         NTYPEZ = NTYPEZ - 1
         GOTO 101
         ENDIF
      DAV = 0.48
      DAV2 = DAV * DAV
      NQQQ = 0
      CALL KERNZI (0, LW, MAXAT)
      KPROG = 0
      IRUN = 0
      FREE98 = .FALSE.
      NOFREE = .FALSE.
      CALL LOGRD(IDDL, 'CALL', KLOG)
      IF (KLOG .GT. 0) THEN
         WRITE (LIS2, FMT='('' $TE Call?? LIT(NLIT):'', A6)') LIT(NLIT)
         IF (LIT(NLIT) .EQ. 'FREE') FREE98 = .TRUE.
         IF (LIT(NLIT) .EQ. 'NOFREE' .AND.
     *      LIT(NLIT-1) .EQ. 'NORECY') NOFREE = .TRUE.
         ENDIF
      CALL LOGRD(IDDL, 'NAT=', KLOG)
      IF (KLOG .GT. 0) IRUN = NINT(FNUM(2))
      IF (IRUN .NE. KEYS(13)) KLOG = 0
      IF (KLOG .GT. 0) THEN
         NATS = NINT(FNUM(3))
         KPROG= NINT(FNUM(4))
         IF (NATS .GT. NATIN) NATS = NATIN
         KEYS(23) = NATS
         ENDIF
      IF (KPROG .GE. 4 .AND. KPROG .LE. 6) KPROG = 0
      IF (KPROG .EQ. 8 .OR. KPROG .GT. 10) KPROG = 0
      IF (KLOG .LE. 0 .OR. (KPROG .GT. 0 .AND. NATS .LE. 0)) THEN
         NATS = NATIN
         IF (KPROG. LE. 0) KPROG = 9
         WRITE(LIS1, FMT='(/'' Old CONDA FILE USED ? No good.. ''/)')
         ENDIF
      IF (KPROG .NE. 0) SWRECY = .TRUE.
      IF (SWNOCY) SWRECY = .FALSE.
      IF (KEYDS.GT. 30 .AND. KEYDS .LT. 40) SWRECY = .FALSE.
      IF (SWRECY .AND. NRECY .EQ. 0) THEN
         NRECY = 11
         NRECYR = 1
         NRECYS = 1
         WRITE (IPR1, 111)
         WRITE (LIS1, 111)
         WRITE (LIS2, 111)
  111    FORMAT (/25X,'PHASEX or Fourier recycling procedure: cycle  1')
         ENDIF
      NRECYO = NRECY
      REN98 = .FALSE.
      IF (FREE98 .AND. NRECYS .GE. 7) THEN
         REN98 = .TRUE.
         NATS = 0
         IF (NRECYS .GT. 7) GOTO 112
         CHOUT = ' From now on all input atoms may be renamed!'
         CALL SHOUT
         CHOUT = ' -------------------------------------------'
         CALL SHOUT
         ENDIF
  112 IF (NOFREE) THEN
         CHOUT = ' All input atoms are renamed following CRYSIN data!'
         CALL SHOUT
         NATS = 0
         ENDIF
      IF (.NOT.SWRECY) THEN
         IPRY=LIS1
      ELSEIF (NRECY .LE. 11) THEN
         IPRY=LIS2
      ELSE
         IPRY=0
         ENDIF
      NATSN = NATS
      CALL LOGRD (IDDL, 'BOV', KLOG)
      CALL FILCLO (IDDL, 'KEEP')
      BOV = 3.
      IF (KLOG .GT. 0) BOV = FNUM(3)
      IF (.NOT. REN98) GOTO 134
      SUMBZ = 0.
      SUMZ = 0.
      DO 121 I = 1, NATIN
      Z2 = FLOAT ( IZAT(I) * IZAT(I) )
      SUMZ = SUMZ + Z2
      SUMBZ = SUMBZ + Z2 * ATXYZ(5,I)
  121 CONTINUE
      BAV = SUMBZ / SUMZ
      WRITE (LIS2, 122) NRECY, BOV, BAV
  122 FORMAT (' $TEMP NRECY, BOV, BAV = ', I4, 2F6.3)
      IF (BAV .LT. 0.1) BAV = BOV
      BOVBAV = BOV / BAV
      DO 125 I = 1, NATIN
      ATXYZ(5,I) = ATXYZ(5,I) * BOVBAV
  125 CONTINUE
  134 CONTINUE
      CALL FILINQ (ICRYS, 'CRYSDA', 'FORMATTED', 'INPUT', KINQ)
      DO 159 N = 1, NTYPE
      CALL RDCRYB (ICRYS, 'ELEM' , KEND)
      IF (KEND.LE.0) THEN
         WRITE (CHOUT, 142) N
  142    FORMAT (' CRYSDA file: ELEM for atom No ', I2, ' not found')
         CALL KERROR (CHOUT, 142, 'DIRBF')
         ENDIF
      CALL KERINB (LIT, 1)
      IZ = NINT(FNUM(1))
      DO 153 I = 1, NTYPE
      IF (NCELLZ(I) .EQ. IZ) GOTO 154
  153 CONTINUE
      CALL KERROR ('Atom-type conflict', 153, 'DIRBF')
  154 RCELTY(I) = AMIN1 (FNUM(4), FNUM(5), FNUM(6))
      IF (NCELLZ(I) .GT. 18) RCELTY(I) = RCELTY(I) + 0.35
      IF (NCELLZ(I) .LT. 10) RCELTY(I) = 0.01
  159 CONTINUE
      CALL FILCLO (ICRYS, 'KEEP')
      NSATL = 0
      NNOMAT = 0
      NNOMAS = 0
      NNOMAX = 0
      NNOMAW = 0
      CALL KERNZI (0, JCON, NPIC)
      MSGT = 0
      MSGS = 0
      ZSCAL = 100. * NCELLZ(1) / X(4,1) ** 2
      ZSCALA = 0.
      ZSCALX = 0.
      CALL EXPEAK(0, 0, 0, IPPQQ)
      QPPQQ = FLOAT(IPPQQ) / 1000.
      IF (QPPQQ .GT. 1.0 .AND. SCALEX .GT. 0.0001) THEN
         WRITE (LIS2, FMT='('' $TE scale /QPPQQ , old SCALE'',
     *      2F9.4)') QPPQQ, SCALEX
         IF (QPPQQ .GT. 1.05) QPPQQ = 1.05
         IF (QPPQQ .LT. 0.95) QPPQQ = 0.95
         SCALEX = SCALEX / QPPQQ
         WRITE (LIS2, FMT='('' $TE scale /QPPQQ: new SCALE='',
     *      2F9.4)') QPPQQ, SCALEX
         WRITE (CHOUT,FMT='(''RUN '',I3, '' CY'', I3,
     *     '' SCALEX '', F10.5)') IRUN, NRECYR, SCALEX
         CALL LOGWR (IDDL)
         CALL FILCLO (IDDL, 'KEEP')
         ENDIF
      CALL KERNZI (0, IRECY, 8)
      IRECYD = 0
      CALL KERNZA (0.0, XBTEMP, MAXAT)
      CALL KERNZI (0  , IZTEMP, MAXAT)
      DO 420 I= 1, NATIN
      SWMATE = .FALSE.
      LW(I) = 0
      IBOND(I) = 0
      ATNAMX(I) = ATNAME(I)
      DUM(I) = 0.
      XDUM(I) = 0.
      CALL ATOMOC (0, ATXYZ(1,I), IOCC, 1)
      ATXYZ(9,I) = IMULT/IOCC(1)
      DO 310 N = 1, NTYPE
      IF (NCELLZ(N) .NE. IZAT(I)) GOTO 310
      ATXYZ(10,I) = FLOAT(N)
      DPAR = RCELTY(N)
      DPAR2 = DPAR * DPAR
      DTEST = AMAX1 (DPAR, DAV)
      GOTO 312
  310 CONTINUE
      CALL KERROR (' Impossible...', 310, 'DIRBF')
  312 DO 410 J = 1, NPIC
      IF (X(4,J) .LE. 0.1) GOTO 410
      DO 400 IS = 1, NSYMM
      CALL SYMOP1 (IS, X(1,J), XS)
      DO 400 IC = 1, ICENT
      DO 400 IL = 1, NLATT
      CALL SYMOP2 (IC, IL, XS, XST)
      CALL DISTSQ (ATXYZ(1,I), XST, DTEST, XSHIFT, DIST2)
      IF (DIST2 .GT. DAV2) GOTO 380
      DIST = SQRT(DIST2)
      IF (DIST .GT. 0.20) GOTO 330
      WAV = 1.
      IF (DIST .GT. 0.05 .AND. DIST .LT. 0.10) WAV = 0.1 / DIST
      IF (DIST .LE. 0.05) WAV = 2.
      GOTO 340
  330 WAV = 5. * DIST
      WRED = AMIN1(6.0, FLOAT(J-NATIN-1)/2.0)
      IF (WRED.GT.1.0) WAV = AMIN1 (WRED*WAV - WRED + 1., DIST/0.08)
  340 DO 350 L = 1, 3
  350 ATXYZ(L,I) = ATXYZ(L,I) + XSHIFT(L) / WAV
      SHIFT = DIST / WAV
      IF (DIST.GT.0.2298 .AND. DIST.LT.0.2302) DIST = 0.2298
      SWMATE = .TRUE.
      DUM(I) = DIST
      XDUM(I) = SHIFT
      IF ((NRECYR.LE.6 .AND. DIST .GT. 0.05) .OR.
     *                       DIST .GT. 0.10) IRECYD = IRECYD + 1
      ATXYZ(6,I) = ABS (X(4,J))
      ATXYZ(8,I) = X(4,J) ** 2 * ZSCAL
      X(4,J) = - ABS(X(4,J))
      LW(I) = J
      IF (JCON(J) .LT. 0) WRITE (LIS1, 357) ATNAME(I), I
  357 FORMAT (' Warning: suspect input atom ', A6, ' (number', I3,
     *        ') : check bonds' )
      JCON(J) = I
      IF (DIST .GT. 0.23) THEN
         MSGT = MSGT + 1
         IF (I .LE. NATS) MSGS = MSGS + 1
         ENDIF
      CALL EXPEAK(1, I, NATSN, IZOLD)
      IF (IZOLD .NE. 0) THEN
         ATNAMX(I) = ATNAME(I)
         IRECY(1) = IRECY(1) + 1
         ENDIF
      ZSCALA = ZSCALA + FLOAT(IZAT(I))
      ZSCALX = ZSCALX + X(4,J) ** 2
      IF (IZAT(I) .GT. 18) GOTO 410
      IF (IZAT(I) .GT. 10 .AND. I .LE. NATS) GOTO 410
      IF (IZAT(I) .GT. 10 .AND. REN98) GOTO 410
      GOTO 420
  380 IF (IZAT(I) .LE. 10) GOTO 400
      IF (IZAT(I) .LE. 18 .AND. I .GT. NATS .AND. .NOT. REN98) GOTO 400
      IF (DIST2 .GT. DPAR2 .OR. X(4,J) .LE. 0.1) GOTO 400
      IF (ABS(X(4,J)) .GT. 0.25 * ATXYZ(6,I)) THEN
         IRECY(2) = IRECY(2) + 1
         IF (ABS(X(4,J)) .GT. 0.50 * ATXYZ(6,I)) THEN
            WRITE(LIS1, 383) J, I, ATNAME(I)
  383       FORMAT (' WARNING: strong peak',I3, '  close to heavy atom',
     *      ' (nr',I3, ' = ', A6,')  retained' )
            GOTO 400
         ELSE
            WRITE(LIS1, 384) J, I, ATNAME(I)
  384       FORMAT (' WARNING: strong peak',I3, '  close to heavy atom',
     *      ' (nr',I3, ' = ', A6,')  deleted' )
            ENDIF
         ENDIF
      NSATL = MIN0 (NSATL+1, 133)
      SATLD(NSATL) = SQRT(DIST2)
      ISATL(NSATL) = I
      JSATL(NSATL) = J
      X(4,J) = -ABS(X(4,J))
      JCON(J) = -I
  400 CONTINUE
  410 CONTINUE
      IF (SWMATE) GOTO 420
      ATXYZ(8,I) = -0.0001
      NNOMAT = NNOMAT + 1
      IRECY(3) = IRECY(3) + 1
      IF (LW(I) .NE. 0) CALL KERROR ('Kanniet', 410, 'DIRBF')
      IF (I .LE. NATS) THEN
         NNOMAS = NNOMAS + 1
         IF (NRECYR .GT. 1) THEN
            IF ( NATS  .GT. 1) THEN
               NNOMAX = NNOMAX + 1
               LW(I) = -1
            ELSE
               LW(I) = -4
               ATXYZ(8,I) = 0.0001
               ENDIF
            ENDIF
         ENDIF
      ATXYZ(6,I) = ATXYZ(8,I)
  420 CONTINUE
      IF (NSATL .EQ. 0) GOTO 460
      IF (NATIN .GT. NATS .AND. .NOT. REN98) WRITE (LIS2, 431)
  431 FORMAT (/' (May be based on the tentative chemical assignement:)')
      WRITE (LIS2, 433) NSATL
  433 FORMAT (' We have removed ', I4,
     *   ' peaks from the Fourier-peaklist,' /
     *   '    because they are too close to the input heavy atom.  ',
     *   '   They are:' /
     *   '  No. peakhght     x       y       z     close to atom:  No.')
      WRITE (LIS1, 435) NSATL
  435 FORMAT (' We have removed ', I4, ' peaks from the peaklist,' /
     * '    because they are too close to the heavy atom: see LIS2' )
      DO 440 N = 1, NSATL
      I = ISATL(N)
      J = JSATL(N)
      DDSATL = ABS ( X(4,J)) ** 2 * ZSCAL
      WRITE (LIS2, 439) J, DDSATL, (X(NN,J), NN=1,3),
     *   ATNAME(I), I, SATLD(N)
  439 FORMAT (I5, F8.0, 3F8.4, '  close to ', A6, I5, '   Dist=', F5.2)
  440 CONTINUE
  460 IF (NRECYR .LE. 1) GOTO 478
      IF (REN98 .OR. NOFREE) GOTO 492
      CALL KERNZA (0.0,  WCELTY, 10)
      CALL KERNZA (0.01, TCELTY, 10)
      DO 471 I = 1, NATIN
      IF (ATXYZ(8,I) .LE. 0.1) GOTO 471
      N = NINT(ATXYZ(10,I))
      IF (TCELTY(N) .GT. 8.5) THEN
         WCELTY(N) = 0.9 * (WCELTY(N) + ATXYZ(8,I))
      ELSE
         WCELTY(N) = WCELTY(N) + ATXYZ(8,I)
         TCELTY(N) = TCELTY(N) + 1.
         ENDIF
  471 CONTINUE
      DO 472 N = 1, NTYPE
  472 WCELTY(N) = WCELTY(N) / TCELTY(N)
      DO 477 I = 1, NATS
      IF (ATXYZ(8,I) .LT. 0.1) GOTO 477
      IF (NRECYR.EQ.3 .AND. IZAT(I).GE.5 .AND. IZAT(I).LE.9
     *    .AND. ATXYZ(8,I) .LT. 200.) NQQQ = 1
      N = NINT(ATXYZ(10,I))
      IF (ATXYZ(8,I) .LT. 0.3 * WCELTY(N)) THEN
         ATXYZ(8,I) = - ABS(ATXYZ(8,I))
         ATXYZ(6,I) = - ABS(ATXYZ(6,I))
         J = LW(I)
         IF (J .LE. 0) GOTO 477
         X(4,J) = - ABS(X(4,J))
         LW(I) = -2
         NNOMAW = NNOMAW + 1
         ENDIF
  477 CONTINUE
      NQQQ = NQQQ + NNOMAW + NNOMAX
  478 CONTINUE
      IF (NATS .EQ. 1 .AND. NATIN .LE. 10) THEN
         IBOND(1) = 1
         GOTO 489
         ENDIF
      II = 0
      DO 487 I= 1, NATS
      J = LW(I)
      IF (J .LT. 0 .AND. J .NE. -4) GOTO 487
      II = II + 1
      IBOND(I) = II
      IF (I .GT. II) THEN
         CALL KERNAB (ATXYZ(1,I), ATXYZ(1,II), 10)
         IZAT(II) = IZAT(I)
         ATNAME(II) = ATNAME(I)
         JCON(J) = II
         ENDIF
  487 CONTINUE
      NATSN  = II
      IF (II .NE. NATS - NNOMAX - NNOMAW) WRITE (LIS1,
     *   FMT= '(/'' Error in count of messages.. ''/)')
  489 IF (NATS .GE. NATIN) GOTO 502
  492 CONTINUE
      DO 501 I = NATS + 1 , NATIN
      J = LW(I)
      IF (J .LE. 0) GOTO 501
      IF (JCON(J) .LT. 0) GOTO 501
      X(1,J) = ATXYZ(1,I)
      X(2,J) = ATXYZ(2,I)
      X(3,J) = ATXYZ(3,I)
      X(4,J) = ATXYZ(6,I)
      XBTEMP(J) = ATXYZ(5,I)
      IZTEMP(J) = IZAT(I)
  501 CONTINUE
  502 CONTINUE
      IF (REN98 .OR. NOFREE) GOTO 8505
      WRITE (LIS2, 503) (ACELTY(I), I=1, NTYPE)
  503 FORMAT (/' Cell contents:  atoms:  ', 10(3X, A2))
      WRITE (LIS2, FMT='('' Count  original input,'')')
      CALL CELZIN (ATXYZ, IZAT, NATSN, NCELLZ, NCELIN)
      WRITE (LIS2, FMT='('' = symmetry included ='')')
      IF (NCELLZ(1) .GE. 10) WRITE (LIS2, 504) (RCELTY(I), I=1, NTYPEZ)
  504 FORMAT (' Satellite skip-distance ', 10F5.2)
 8505 CONTINUE
      DO 505 N = 1, NTYPE
      RCELTY(N) = 1.25
      IF (NCELLZ(N) .LE. 12) RCELTY(N) = 1.55
      IF (NCELLZ(N) .LT. 10) RCELTY(N) = 0.85
      IF (NCELLZ(N) .GT. 18) RCELTY(N) = 1.60
  505 CONTINUE
      WRITE (LIS2, 507) (RCELTY(I), I=1, NTYPE)
  507 FORMAT (' Approx. atomic radius:  ', 10F5.2)
      IF (NNOMAW .GT. 0) WRITE (LIS2, 511) (WCELTY(I), I=1, NTYPE)
  511 FORMAT (' Averaged low peak level:', 10F5.0)
      WRITE (LIS2, FMT='('' '')')
      IF (NNOMAT .EQ. 0) THEN
         CHOUT = ' Results: all input atoms recognized.'
         CALL SHOUT2
      ELSE
         WRITE (CHOUT, FMT='('' Nr of input atoms'',
     *     '' rejected because no mate was found:'',I4)')  NNOMAT
         CALL SHOUT
         IF (NNOMAS .GT. 0 .AND. NNOMAT .GT. NNOMAS) THEN
            WRITE (CHOUT, FMT='(I6, '' of these belong to the'',
     *           '' original list of input atoms'')')  NNOMAS
            CALL SHOUT
            ENDIF
         ENDIF
      IF (NNOMAW .GT. 0) THEN
         WRITE (CHOUT, 515) NNOMAW
  515    FORMAT (' Nr of original input atoms rejected: ',
     *       'too low Peakheight:', I3)
         CALL SHOUT
         IRECY(4) = IRECY(4) + 1
         ENDIF
      I = NNOMAX + NNOMAW
      IF (I .GT. 0) THEN
         WRITE (LIS1, 517) I
         WRITE (LIS2, 517) I
  517 FORMAT (' Because', I3,
     *       ' of the original input atoms will be rejected,' /
     *       ' the original input atoms have been renumbered.'/
     *       ' The atom-names are not changed  '/
     *       '    except when the atom type is changed !  .' )
         ENDIF
      ZSCAL1 = ZSCAL
      ZSCAL = 100. * ZSCALA / ZSCALX
      WRITE (LIS2, 1531) ZSCAL1, ZSCAL
 1531 FORMAT (/' $TEMP :   old and new ZSCAL =', 2F15.7)
      WRITE (LIS2, 531)
  531 FORMAT (/' INPUT ATOMS AND THEIR AVERAGED MATES'//
     +   '   ATOM        PEAK        modified coordinates   ',
     *   ' mate atom-  (input:)'/
     +   '  No. NAME    No. Integr.   x        y        z   ',
     *   '-DIST SHIFT  No. NAME' )
      IF (MSGS .GT. 0 .OR. NNOMAS + NNOMAW .GT. 0) WRITE (LIS1, 533)
  533 FORMAT (/'  No. ATOM    No. PEAK')
      IF (REN98 .OR. NOFREE) GOTO 592
      DO 570 I = 1, NATS
      ABSX = -1.
      II = IBOND(I)
      IF (II .EQ. 0) GOTO 553
      IF (ABS(ATXYZ(8,II)) .GT. 1.)
     *   ATXYZ(8,II) = ATXYZ(8,II) * ZSCAL / ZSCAL1
      IF (I .LE. 0) CALL KERROR ('kanniet', 545, 'DIRBF')
      J = LW(I)
      ABSX = ABS(ATXYZ(8,II))
      IF (II .EQ. I) THEN
         WRITE (LIS2, 545) II, ATNAMX(I), J, ABSX,
     *      (ATXYZ(K,II),K=1,3), DUM(I), XDUM(I)
  545    FORMAT (I5, 1X, A6, I4, F7.0, 1X, 3F9.5, 2F5.2, 1X, I4, 1X, A6)
      ELSE
         WRITE (LIS2, 545) II, ATNAME(II), J, ABSX,
     *      (ATXYZ(K,II),K=1,3), DUM(I), XDUM(I), I, ATNAMX(I)
         ENDIF
      IF (DUM(I) .GT. 0.23) THEN
         WRITE (LIS2, 547)
  547    FORMAT (53X, '----', 8X, 'WARNING')
         WRITE (LIS1, 549) II, ATNAME(I), J, ABSX, I, DUM(I)
  549    FORMAT (I5, 1X, A6, I5, F7.0, ' MATE FOR INPUT ATOM No.',
     *      I4, ' found at DIST = ', F4.2 )
         ENDIF
      IF (J .EQ. -4) THEN
         WRITE (CHOUT,
     *      FMT = '(12X, ''  This only input atom is not rejected'')')
         CALL SHOUT2
         GOTO 570
         ENDIF
      GOTO 570
  553 IF (J .EQ. -2) GOTO 564
      III = II
      IF (J .EQ. -1) III = 0
      WRITE (LIS2, 555) III, ATNAMX(I), I, ATNAMX(I)
      WRITE (LIS1, 555) III, ATNAMX(I), I, ATNAMX(I)
  555 FORMAT (I5, 1X, A6, '  WARNING: no mate found for input atom', I4/
     *      ' ---> ', A6, '  will be removed from',
     *        ' the list of original input atoms')
      IF (J .EQ. 0) CALL KERROR ('KANNIET', 555, 'DIRBF')
      GOTO 570
  564 WRITE (LIS1, 566) ATNAMX(I), ABSX, I
      WRITE (LIS2, 566) ATNAMX(I), ABSX, I
  566 FORMAT (4X, '0', 1X, A6, F6.0,
     *  ' : peak height is too low for input atom', I3)
  570 CONTINUE
      IF (NRECYR .LE. 1) THEN
         WRITE (LIS1, 572)
         WRITE (LIS2, 572)
  572    FORMAT (/' The input atoms are appended by the remaining ',
     *    '(unidentified) peaks,' / ' the tentative chemical ',
     *    ' assignment is based on the CRYSDA file' /
     *    ' and the peak height, not on any chemical argument.' / )
      ELSEIF (NRECYR .EQ. 2) THEN
         WRITE (LIS1, 575)
         WRITE (LIS2, 575)
  575    FORMAT (/' Note: the secondary input atoms, tentatively',
     *      ' assigned and added to' / ' the original list in',
     *      ' the foregoing cycle(s), are now merged with the '/
     *      ' new (unidentified) peaks, sorted',
     *      ' to peakheight, and reinterpreted.'/)
         ENDIF
  592 CALL KERNZI (0, LCELTY, 10)
      LCELTY(NTYPEZ) = LCELTY(NTYPEZ) + 9999
      L = 0
      IF (REN98 .OR. NOFREE) THEN
         CALL KERNZI (0, NCELIN, 10)
         ENDIF
      DO 600 N = 1, NTYPEZ
      LCELTY(N) = LCELTY(N) + NCELTY(N) - NCELIN(N)
      IF (LCELTY(N) .LT. 0) THEN
         LCELTY(N+1) = LCELTY(N+1) + LCELTY(N)
         LCELTY(N) = 0
         ENDIF
      IF (L .EQ. 0 .AND. LCELTY(N) .GT. 0) L = N
  600 CONTINUE
      LH = 3
      IF (ACELTY(NTYPEZ) .NE. 'H ') LH = 0
      NPC = NPIC - NATX
      NATXX = 99
      IF (NPC .NE. 0) NATXX = MAX0(5, (NATX-NATS)/7 +3)
      NAT = NATSN
      DO 603 J = 1, NPIC
      IF (X(4,J) .LE. 0.1) GOTO 603
      PEAKI = X(4,J) **2 * ZSCAL
      IF (PEAKI . LT. 50.) GOTO 604
      NAT = NAT + 1
      CALL KERNAB (X(1,J), ATXYZ(1,NAT), 3)
      ATXYZ(4,NAT) = 1.
      CALL KERNZA (0., ATXYZ(5,NAT), 6)
      ATXYZ(5,NAT) = XBTEMP(J)
      ATXYZ(6,NAT) = X(4,J)
      ATXYZ(8,NAT) = PEAKI
      ATXYZ(7,NAT) = FLOAT(J)
  603 CONTINUE
  604 CONTINUE
      NAT = MIN0 (NAT, NATX)
      NATNOH = NATSN
      DO 620 I = NATSN + 1, NAT
      CALL ATOMOC (0, ATXYZ(1,I), IOCC, 1)
      LLL = L
      IF (L .EQ. NTYPEZ.AND. LH .GT. 0 .AND. ATXYZ(8,I) .GT. 200.) THEN
         LLL = L - 1
         LH = LH - 1
         ENDIF
      CALL ATN4CN (ACELTY(LLL), I, I-1, ATNAME, I, ATNAME(I))
      J = NINT (ATXYZ(7,I))
      ATXYZ(7,I) = 0.
      N = JCON(J)
      KEY = 2
      IF (N .LE. 0) KEY = 3
      IF (REN98) KEY = 4
      IZAT(I) = NCELLZ(L)
      ATXYZ(10,I) = L
      CALL EXPEAK(KEY, I, I-1, IZOLD)
      IF (NRECYR .LE. 4 .AND. IZOLD .NE. 0) IRECY(5) = IRECY(5) + 1
      IF (N .GT. 0) THEN
         LW(N) = -3
         WRITE (LIS2, 545) I, ATNAME(I), J, ATXYZ(8,I), (X(K,J), K=1,3),
     *      DUM(N), XDUM(N), N, ATNAMX(N)
      ELSE
         WRITE (LIS2, 545) I, ATNAME(I), J, ATXYZ(8,I)
         ENDIF
      LIOCC = IMULT/IOCC(1)
      ATXYZ(9,I) = LIOCC
      LCELTY(L) = LCELTY(L) - LIOCC
      IF (IZAT(I) .NE. 1) NATNOH = NATNOH + 1
  615 IF (LCELTY(L) .LE. 0) THEN
         L = L + 1
         LCELTY(L) = LCELTY(L) + LCELTY(L-1)
         IF (LCELTY(L) .LE. 0) GOTO 615
         ENDIF
      IF (L .LT. NTYPEZ) GOTO 620
      IF (NCELIN(L).GT.0 .AND. LCELTY(L) .GT. 9999) GOTO 620
      IF (IZAT(I) .GT. 1) GOTO 620
      NATXX = NATXX - 1
      IF (NATXX .LE. 0) THEN
         NAT = I
         GOTO 640
         ENDIF
  620 CONTINUE
  640 CONTINUE
      I = MSGT - MSGS
      IF (I .GT. 0) WRITE (LIS1, 641) I
  641 FORMAT (' Of the atoms tentatively added in the forgoing cycle'/
     *       I5, ' were found at distances larger than 0.23 Angstrom.')
      IF (NATIN .LE. NATS) GOTO 666
      II = 0
      DO 643 I = NATS+1, NATIN
      IF (LW(I) .LT. 0) GOTO 643
      II = II + 1
  643 CONTINUE
      IF (II .EQ. 0) GOTO 666
      WRITE (LIS1, 644) II
  644 FORMAT (' Of the atoms tentatively added in the forgoing cycle'/
     *    I5, ' is(are) left out: no mate found or too weak to use.')
      WRITE (LIS2, FMT='('' Old No.  Old name   PEAK   HEIGHT'')')
      DO 647 I = NATS+1, NATIN
      IF (LW(I) .EQ. 0) WRITE (LIS2, 646) I, ATNAMX(I)
  646 FORMAT (I8, 4X, A6, I7, 4X, F6.0, F9.5)
      IF (LW(I) .LE. 0) GOTO 647
      J = LW(I)
      BSUM = ATXYZ(8,I)
      WRITE (LIS2, 646) I, ATNAMX(I), J, BSUM, X(4,J)
  647 CONTINUE
  666 IF (NQQQ.GT.0 .AND. NRECYR.LE.4) IRECY(6) = IRECY(6) + 1
      I = (IRECYD * 30)/ NAT
      IF (NAT .GT. 80) I = (IRECYD * 20)/ NAT
      IF (NAT .GT. 120) I = (IRECYD * 10)/ NAT
      I = MIN0 (19, I)
      IF ((NRECYR .LE. 6 .AND. I .GT. 0) .OR.
     *   (NRECYR .LE. 10 .AND. I .GE. 2)) IRECY(7) = I
      IF ((NRECYR.LE.6 .AND. ATXYZ(8,NAT).GT. 150.) .OR.
     *                       ATXYZ(8,NAT).GT. 250.) IRECY(8)=IRECY(8)+1
      IRECYS = 0
      DO 677 I = 1, 8
      IRECYS = IRECYS + IRECY(I)
  677 CONTINUE
      IF (IRECYS.GT.0 .AND. NRECYR.LE.8.AND.NRECYS.EQ.6) NRECYS=NRECYS-1
      IF (NRECYR .GE. 13) NRECYS = 9
      IF (NRECYS .GE. 9) THEN
         SWRECY = .FALSE.
         IPRY = LIS1
         ENDIF
      NRECY = NRECYR * 10 + NRECYS
      IF (.NOT. REN98 .AND. .NOT.NOFREE) CALL EXPEAK(-1, 0, 0, IDUMM)
      IF (.NOT. SWRECY) SWNOCY = .TRUE.
      IF (SWNOCY) GOTO 690
      PSQN = 1.0
      IF (PSQ .LT. 0.30) PSQN = AMIN1 (0.70, 3. * PSQ)
      IF (PSQ .LT. 0.07) PSQN = 0.20
      NNNN = 0
      DO 682 N = 1, NTYPE
  682 NNNN = NNNN + NCELTY(N) * NCELLZ(N) ** 2
      NNNR = FLOAT (NNNN) * PSQN * 1.00001
      N = 0
  685 N = N + 1
      IF (N .GT. NAT) GOTO 687
      ITYPE = NINT(ATXYZ(10,N))
      IF (ITYPE .EQ. 0) GOTO 687
      IF (NCELLZ(ITYPE) .EQ. 1) GOTO 687
      NNNR = NNNR - NINT(ATXYZ(9,N)) * NCELLZ(ITYPE) ** 2
      IF (NNNR .LT. 0) GOTO 687
      GOTO 685
  687 NATPSQ = N - 1
      NDIFF = MIN0 (MAX0 (2, NATX-NATIN) * 30 / 40, 50)
      NATDIF = NATPSQ
      IF (NDIFF .GT. 30) NATDIF = NATIN + NDIFF
      NATREC = MIN0 (NAT, NATPSQ, NATNOH, NATDIF)
      IF (NATNOH .LE. NATREC * 20 /18) NATREC = NATNOH
      NATREC = MAX0 (NATIN, NATREC)
      IF (NRECYR .LE. 2) GOTO 688
      IF (NRECYR .EQ. 3 .AND. NRECYS .LT. 5) NRECYS = 5
      IF (NATREC .EQ. NATNOH .AND. NRECYS .LT. 4) NRECYS = 5
      IF (NATREC .LE. NATIN .AND. NRECYS .LT. 5)  NRECYS = 5
      IF (R2X .LT. 0.20 .AND. NRECYS .LT. 4) NRECYS = 5
  688 NRECY = NRECYR * 10 + NRECYS
  690 CONTINUE
      CALL FILINQ (IDDS, 'DDSYST', 'FORMATTED', 'OUTPUT', KINQ)
      CALL FILINQ (ICOND, 'CONDA', 'FORMATTED', 'OUTPUT', KINQ)
      IF (SWNOCY) THEN
         WRITE (ICOND,
     *      FMT = '(''CONDA  '', A6, '' generated by FOUR'' )') CCODE
         IF (TITLE .NE. ' ')
     *      WRITE (ICOND, FMT = '(''TITLE '', A64 )') TITLE
         WRITE (ICOND, FMT = '(''PROGRAM NUTS AT2X'' / ''FINISH'')')
         WRITE (IDDS, FMT = '(''NUTS'' /  ''STOP'')')
         GOTO 785
         ENDIF
      WRITE (ICOND, FMT = '(''CONDA  '', A6, '' generated by FOUR'' )')
     *   CCODE
      IF (TITLE .NE. ' ')
     *   WRITE (ICOND, FMT = '(''TITLE '', A64 )') TITLE
      WRITE (ICOND, FMT = '(''PROGRAM DDMAIN          '')')
      WRITE (IDDS, FMT = '(''DDMAIN'')')
      IF (KPROG .GT. 0 .AND. KPROG .NE. 9 .AND. NRECYR.LE.3 ) THEN
         WRITE (ICOND, FMT = '(''OPTION 1 PHASEX 0 NRECY'', I4,
     *      '' NATL='', I5)') NRECY, NATIN
         WRITE (ICOND, FMT = '(''PROGRAM  PHASEX  '')')
         WRITE (ICOND, FMT = '(''PROGRAM  DDMAIN  '')')
         WRITE (ICOND, FMT = '(''OPTION 2 FOUR.PH 0 NRECY'',I4)') NRECY
         WRITE (IDDS,  FMT = '(''PHASEX'' / ''DDMAIN'')')
      ELSE
         WRITE (ICOND, FMT = '(''OPTION 3 FOUR 0 NRECY'', I4,
     *      '' NATL='', I5)') NRECY, NATIN
         ENDIF
      WRITE (ICOND, FMT = '(''PROGRAM  FOUR  NRECY '', I3)') NRECY
      WRITE (ICOND, FMT = '(''PROGRAM  NUTS AT2X'')')
      WRITE (ICOND, FMT = '(''FINISH'')')
      WRITE (IDDS, FMT = '(''FOUR'' /  ''NUTS'' / ''STOP'')')
  785 CALL FILCLO  (ICOND, 'KEEP')
      CALL FILCLO (IDDS, 'KEEP')
      RTYMAX = 0.0
      DO 800 I = 1, NAT
      CALL KERNAB (ATXYZ(1,I), X(1,I), 3)
      X(4,I) = ATXYZ(6,I)
      IF (ATXYZ(8,I) .GT. 0.0) GOTO 798
      X(4,J) = -ATXYZ(6,I)
      IF (ATXYZ(8,I) .GT. -1.1) X(4,J) = ATXYZ(8,I)
  798 N = NINT(ATXYZ(10,I))
      IF (N .EQ. 0) GOTO 800
      IF (RCELTY(N) .GT. RTYMAX) RTYMAX = RCELTY(N)
      ATXYZ(7,I) = RCELTY(N)
  800 CONTINUE
      NAT = MIN0 (NAT, NATNOH + 2 + NATNOH / 30)
      NPIC = NAT
      RTYMAX = 2.0 * RTYMAX
      IF (RTYMAX .GT. DMAXB .AND. ABS(1.95 - DMAXB) .LT. 0.001) THEN
         DMAXB = RTYMAX
         IF (DMAXB .GT. DMOUT) DMOUT = DMAXB
         ENDIF
      NNNN = 0
      DO 812 N = 1, NTYPE
      IF (NCELLZ(N) .NE. 1) NNNN = NNNN + NCELTY(N) * NCELLZ(N) ** 2
  812 CONTINUE
      ILAST = 1
      NNNR = 0
      DO 822 I =1, NAT
      IF (IZAT(I) .NE. 1) NNNR = NNNR + NINT(ATXYZ(9,I)) * IZAT(I) ** 2
      IF (IZAT(I) .NE. 1 .AND. NNNR .LE. NNNN) GOTO 820
      IF (NNNR .LE. NNNN*125/100) GOTO 820
      ATNAME(I) (1:1)  = 'Q'
      ATNAM = ATNAME(I) (3:6)
      CALL KERC2I (ATNAME(I)(2:2), LEND)
      IF (LEND.LT.0 .OR. LEND.GT.9) ATNAME(I)(2:6) = ATNAM
      IZAT(I) = 1
  820 CONTINUE
      IF (IFIX(ATXYZ(8,I)) .GE. 100) ILAST = I
      IF (NRECYR .LE. 2 .AND. IFIX(ATXYZ(8,I)) .GE. 60) ILAST = I
      IF (NRECYR .LE. 1 .AND. IFIX(ATXYZ(8,I)) .GE. 40) ILAST = I
  822 CONTINUE
      NAT = ILAST
      WRITE (LIS2, 833) NRECYO, IRECY, IRECYD, NRECY
  833 FORMAT (/ ' $TEMP NRECYO IRECY -D NRECY  ', I4, 1X, 8I3, I5, I4)
      IF (NRECYS .LT. 7) RETURN
      SUMB2Z = 0.
      ISUMZ = 0
      DO 901 I=1, NAT
      CALL EXPEAK(5, I, 0, IB100)
      FB100 = FLOAT (IB100) / 100.
      B2Z = FLOAT( IZAT(I) ) * FB100 * FB100
      IF (IB100 .LT. 0) B2Z = -B2Z
      SUMB2Z = SUMB2Z + B2Z
      ISUMZ = ISUMZ + IZAT(I)
  901 CONTINUE
      SUMB2Z = SUMB2Z / FLOAT (ISUMZ)
      SUMB2S = SQRT ( ABS (SUMB2Z) )
      IF (SUMB2Z .LT. 0.0) SUMB2S = - SUMB2S
      WRITE (LIS2, FMT='('' $TE AVERAGE B SHIFT IS:'', F7.3,
     *   '' reset !''  )') SUMB2S
      DO 905 I=1, NAT
      ATXYZ(5,I) = ATXYZ(5,I) - SUMB2S
  905 CONTINUE
      IF (SWNOCY .OR. NRECYS .LT. 8) RETURN
      WRITE (LIS2, FMT='('' $TEMPPP last line DIRBF '')')
      RETURN
      END
      SUBROUTINE ATN4CN (CHEM, NUMB, KEY, ATNAME, NAT, ATNEW)
      CHARACTER*2 CHEM
      CHARACTER*6 ATNAME(NAT), ATNEW, C, CNUMB
      ATNEW = ' '
      N = NUMB
  607 CALL KERI2C (N, CNUMB, 6)
      C(1:2) = CHEM
      C(3:6) = CNUMB
      IF (C(2:2) .EQ. ' ') C(2:6) = CNUMB
      IF (KEY .LE. 0) GOTO 999
      CALL KEREQ6 (C, ATNAME, KEY, KEND)
      IF (KEND .LE. 0) GOTO 999
      N = N + 100
      IF (N .LT. 300) N = N + 100
      IF (N .LT. 300) N = N + 100
      GOTO 607
  999 ATNEW = C
      RETURN
      END
      SUBROUTINE ATN2CN (ATNAM, CHEM, NUMB)
      CHARACTER*2 CHEM
      CHARACTER*6 ATNAM
      CHEM = ATNAM
      NUMB = 0
      IF (CHEM(2:2) .EQ. ' ') RETURN
      I = 1
      CALL KERC2I (ATNAM(2:2), LEND)
      IF (LEND.LT.0 .OR. LEND.GT.9) I=2
      IF (I .EQ. 1) CHEM(2:2) = ' '
      DO 101 N = I+1, 6
      CALL KERC2I (ATNAM(N:N), LEND)
      IF (LEND .EQ. 10) RETURN
      IF (LEND.GE.0 .AND. LEND.LE.9) THEN
         NUMB = 10*NUMB + LEND
      ELSE
         IF (NUMB .EQ. 0) NUMB = 999
         NUMB = - NUMB
         RETURN
         ENDIF
  101 CONTINUE
      RETURN
      END
      SUBROUTINE ATN24X (ATNAM, ATNAME, NAT, ATNEW)
      CHARACTER*6 ATNAME(NAT), ATNAM, ATNEW, C
      CHARACTER*2 CHEM
      C = ATNAM
      ATNEW = ' '
      CALL KEREQ6 (C, ATNAME, NAT, KEND)
      IF (KEND .LE. 0) THEN
         ATNEW = C
         RETURN
         ENDIF
      CALL ATN2CN (C, CHEM, NUMB)
      N = IABS(NUMB)
      CALL ATN4CN (CHEM, N, NAT, ATNAME, NAT, ATNEW)
      RETURN
      END
      SUBROUTINE EXPEAK(KEY, IAT, MIAT, IZOLD)
      COMMON /SYSTA/ IFILE(20), IFSTAT(20), ITIME(4),   KEYS(28),
     *               NFNUM,     NLIT,       NCOLN(32),  NCOLL(32),
     *               NFDOT(32), NFDOL(32),  NLUSER(32), FNUM(32),
     *               SWITCH(28)
      LOGICAL        SWITCH
      COMMON /SYSTB/ PROGNM,    PROSNM,     CCODE,      TITLE,
     *               CHIN,      LIT(32),    CHOUT
      CHARACTER      PROGNM *8, PROSNM *6,  CCODE *6,   TITLE *64,
     *               CHIN *80,  LIT *6,     CHOUT *72
      EQUIVALENCE (IDDL, IFILE(1))
      EQUIVALENCE (LIS1, IFILE(7)), (LIS2, IFILE(8))
      EQUIVALENCE (KEYS(27), IMAP)
      LOGICAL SWRECY
      EQUIVALENCE (SWITCH(17), SWRECY)
      COMMON /CRYSA/ CELL(6),  CELLSD(6),   RCELL(6),   VOLUM,
     *               WAVE,     CELALL(10),  AMOLW,      ZET,
     *               NELEC,    F000,        ABSMU,      ICENT,
     *               ILATT,    ISYST,       ILAUE,      IMULT,
     *               IUNIQ,    IPOLA,       NTYPE,      NSYMM,
     *         IRSYMM(3,3,24), TSYMM(3,24), NLATT,      TLATT(3,4),
     *         FRAC2C(3,3),    CART2F(3,3), RRMAT(3,3), SSMAT(3,3)
      COMMON /CRYSB/ SPGR,     WAVEAT,      CELATY(10)
      CHARACTER      SPGR *16, WAVEAT *2,   CELATY *2
      CHARACTER      CHR16 *16, CHR4 *4, CHR1 *1, CHR2 *2
      PARAMETER (MAXAT=993)
      COMMON /XATXYZ/ X(4,MAXAT), ATXYZ(10,MAXAT), IZAT(MAXAT)
      COMMON /ATNAMA/ ATNAME(MAXAT)
      CHARACTER ATNAME *6
      COMMON /SEARDA/ D2R, DMPIC, DMAXB, DMOUT, DMINB, ANGM(2), MCON,
     *        SEARDX, NPIC, NATIN, NAT, NATS, NATX, NATSN, BOV, IPRY,
     *        NRECY, NRECYR, NRECYS, PSQ, NATREC, KPROG, SCALEX, R2X
      COMMON /DIRBFA/ NCELTY(10), NCELLZ(10), NCELIN(10), NCELIX(10)
      COMMON /DIRBFB/ ACELTY(10)
      CHARACTER ACELTY *2
      DIMENSION PH(10,5), BFAC(5)
      DIMENSION SUM1(10,2), SUM2(10,2), SUM3(10,2), EXPP(10), EXPP2(10)
      DIMENSION RSUM(10,2), TSUM1(2), TSUM2(2), TSUM3(2)
      DIMENSION NCELIC(10), DBDP(10)
      CHARACTER *6 ATOLD
      CHARACTER *2 AC
      DATA ISKIP /0/, BDELTN /0./, TPH /8./, ATOLD /' '/, BR /3./
      DATA NUMB  /0/, XPPXX  /0.80/
      NTYPEZ = NTYPE
  101 IF (ACELTY(NTYPEZ) .EQ. ' ') THEN
         NTYPEZ = NTYPEZ - 1
         GOTO 101
         ENDIF
      IZOLD = 0
      MTYPE = NTYPEZ
      IF (ACELTY(NTYPE) .EQ. 'H ') MTYPE = MTYPE - 1
      IF (KEY .GT. 0) GOTO 203
      IF (KEY .LT. 0) GOTO 603
      CALL KERNZI (0, NCELIC, 10)
      CALL KERNZA (0., PH, 50)
      CALL KERNZA (0., BFAC, 5)
      IF (NRECY .LE. 11)  WRITE(LIS2,FMT='( / '' EXPH:'' ,
     *   '' Expected peak heights of atoms for five B values:'')')
      CALL LOGRD(IDDL, 'WRPEAK', KLOG)
      IF (KLOG .LT. 0) THEN
         CHOUT=' Peakheights not found on DDLOG file: rerun MERBIN !'
         CALL SHOUT
         GOTO 801
         ENDIF
      DO 150 IB = 1,5
      BFAC(IB) = FNUM(IB+1)
  150 CONTINUE
      IF (NRECY .LE. 11) WRITE (LIS2,FMT='(''    for B =   '',5F7.3 /
     *    '' Atom         '', 35(''-'')) ') BFAC
      DO 166 J=1,NTYPE
      READ (IDDL,FMT='(A16,A4,I3,A1,A2,F8.2,4F7.2)')CHR16,CHR4,J2,
     *      CHR1,CHR2,(PH(J,IH), IH=1,5)
      IF (CHR4 .NE. 'TYPE') THEN
         CHOUT = ' EXPEAK: DDLOG file incomplete! Rerun MERBIN '
         CALL SHOUT
         GOTO 801
         ENDIF
      IF (J2 .NE. J .OR. CELATY(J) .NE. CHR2) THEN
         CHOUT = ' EXPEAK: Atom types are inconsistent! Rerun MERBIN ?'
         CALL SHOUT
         GOTO 801
         ENDIF
      IF (NRECY .LE. 11)
     *   WRITE(LIS2,FMT='('' Type'',I2,3X,A2,1X,5F7.2)') J, CHR2,
     *     (PH(J,IH), IH=1,5)
  166 CONTINUE
      CALL LOGRD (IDDL, 'BP', KLOG)
      IF (KLOG .GT. 0) BR  = FNUM(4)
      CALL FILCLO (IDDL, 'KEEP')
      BDELTN = BFAC(2)-BFAC(1)
      BNMIN = AMAX1 (0.9, 0.9 * BFAC(1))
      BNMAX = AMIN1 (9.9, 1.1 * BFAC(5))
      ISKIP = 1
      CALL KERNZA (0.0, SUM1, 20)
      CALL KERNZA (0.0, SUM2, 20)
      CALL KERNZA (0.0, SUM3, 20)
      CALL KERNZA (0.0, RSUM, 20)
      CALL KERNZA (1.0, EXPP, 10)
      CALL KERNZA (1.0, EXPP2, 10)
      CALL KERNZA (0.1, DBDP, 10)
      DO 171 JTYPE = 1, NTYPE
      AC = ACELTY(JTYPE)
      IF   (AC.EQ.'P ' .OR. AC.EQ.'S ' .OR. AC.EQ.'CL' .OR. AC.EQ.'SE'
     * .OR. AC.EQ.'F ' .OR. AC.EQ.'BR' .OR. AC.EQ.'I ') EXPP(JTYPE)=.90
      IF   (AC.EQ.'O ' .OR. AC.EQ.'N ' .OR. AC.EQ.'C ') EXPP(JTYPE)=.80
      DBDP(JTYPE)=BDELTN/ AMIN1(-0.005, (PH(JTYPE,4)-PH(JTYPE,3))) /200.
  171 CONTINUE
      EXPP(NTYPEZ)= 0.70
      EXPP(MTYPE) = 0.70
      IF (NRECY .LE. 11) WRITE(LIS2,FMT=
     *   '(/'' EXPP atoms :'', 10(4X,A2))') (ACELTY(J), J=1,NTYPE)
      IF (NRECY .LE. 11) WRITE(LIS2,FMT=
     *   '( '' TEMP relax :'', 10F6.2 )') (EXPP(J), J=1,NTYPE)
      QQN = 0.
      QQS = 0.
      WQQN = 0.
      WQQS = 0.
      DO 172 I = 1, NTYPEZ
      IF (PH(I,3) .LT. 3.0) GOTO 172
      QQS = QQS + PH(I,3) * CELALL(I)
      QQN = QQN + CELALL(I)
      IF (CELATY(I) .EQ. 'C' .OR. CELATY(I) .EQ. 'N' .OR.
     *   (CELATY(I) .EQ. 'O' )) GOTO 172
      WQQS = WQQS + PH(I,3) * CELALL(I)
      WQQN = WQQN + CELALL(I)
  172 CONTINUE
      QQA = QQS/QQN
      WQQA = - 1.0
      IF (WQQN .GT. 9.) WQQA = WQQS / WQQN
      IQQ = NPIC * 95 / 100
      QQN = 0.
      QQS = 0.
      WQQN = 0.
      WQQS = 0.
      DO 174 J = 1, IQQ
      QQQ = ABS ( X(4,J) )
      KK = QQQ * 10000.
      KK = KK - IFIX ( QQQ ) * 10000
      KK = ( FLOAT(KK) )**2 / 4000.
      FKK = KK
      IF (FKK .LT. 250.) GOTO 174
      QQS = QQS + FKK / 100.
      QQN = QQN + 1.
      WQQS = WQQS + FKK / 100.
      WQQN = WQQN + 1.
  174 CONTINUE
      QQB = QQS/QQN
      QPPQQ = QQB / QQA
      WRITE (LIS2, FMT ='(
     *   '' $TEMP Aver.[PH/EXPH=XPPF] = QPPQQ ='', F 6.3)') QPPQQ
      IF (WQQA .LT. 0.0) GOTO 184
      IF (WQQN .LT. 17.) GOTO 184
      WQQB = WQQS/WQQN
      WPPQQ = WQQB / WQQA
      IZOLD = NINT (1000. * WPPQQ)
  184 XPPXX = QPPQQ * 0.90
      IF (NRECY .LE. 11) THEN
         XPPXX = QPPQQ * 0.60
         IF (NATS .LT. NPIC/5 .OR. PSQ .LT. 0.30) XPPXX = QPPQQ * 0.40
         IF (PSQ .LT. 0.20)  XPPXX = QPPQQ * 0.25
         ENDIF
      IF (PSQ .LT. 0.50)  XPPXX = QPPQQ * 0.50
      IF (NRECYR.EQ.2 .OR. NRECYR.EQ.4) XPPXX = QPPQQ * 0.60
      IF (NRECYR.EQ.3 .OR. NRECYR.EQ.5) XPPXX = QPPQQ * 0.70
      IF (NRECYR.GE.8 .OR. NRECYS.GE.8) XPPXX = QPPQQ * 0.75
      IF (IMAP .EQ. 1) XPPXX = 0.25
      DO 194 JTYPE = 1, NTYPE
      EXPP2(JTYPE) = EXPP(JTYPE) * XPPXX
  194 CONTINUE
      WRITE (LIS2, FMT='('' $TEMP EXPP2'', 10F6.2)')
     *   (EXPP2(IJ), IJ = 1, NTYPE)
      RETURN
  203 CONTINUE
      IF (ISKIP .EQ. 0) RETURN
      TT44 = 1.0
      IF (ABS(ATXYZ(4,IAT)) .LT. 0.9) THEN
         TT44 = ABS(ATXYZ(4,IAT))
         IF (IPRY .GT. 0) WRITE (LIS1, 205) IAT, ATNAME(IAT), TT44
         WRITE (LIS1, 205) IAT, ATNAME(IAT), TT44
  205    FORMAT (' Warning: reduced occupancy of atom', I4,
     *     ' = ', A6,'  is',F7.4)
         ENDIF
      KK = ATXYZ(6,IAT) * 10000.
      KK = KK - IFIX ( ATXYZ(6,IAT)) * 10000
      KK = ( FLOAT(KK) )**2 / 4000.
      FKK = KK
      JNAT = 1
      IF (KEY .EQ. 3) JNAT = 2
      MIATX = MIAT
      IF (KEY .GE. 4) MIATX = 0
      XPPXX2 = 0.90
      IF (JNAT .EQ. 2) XPPXX2 = 0.85
      IF (NRECYR .LE. 1 .AND. JNAT .EQ. 2) XPPXX2 = 0.80
      IZOLD = 0
      NAGAIN = 0
  209 N = NINT(ATXYZ(10,IAT))
      DO 210 JTYPE = 1,NTYPEZ
      IF (CELATY(JTYPE) .EQ. ACELTY(N)) GOTO 215
  210 CONTINUE
  215 CONTINUE
      BTEMP = ATXYZ(5,IAT)
      IF (BTEMP .LE. 0.0001) BTEMP = BR
      BTEMPP = 3.0 * BTEMP
      IF (BTEMP .LE. BFAC(1)) THEN
         TPH = PH(JTYPE,1)
      ELSEIF (BTEMP .GE. BFAC(5)) THEN
         TPH = PH(JTYPE,5)
      ELSE
         DO 300 IB = 2,5
         IF (BTEMP .LE. BFAC(IB)) THEN
            BRATIO = (BTEMP - BFAC(IB-1)) / BDELTN
            TPH = PH(JTYPE,IB-1)+((PH(JTYPE,IB)-PH(JTYPE,IB-1))*BRATIO)
            GOTO 325
            ENDIF
  300    CONTINUE
         ENDIF
  325 CONTINUE
      IF (TPH .LT. 0.1) TPH = 0.7
      EXPH = TPH * 100.
      EXPH = TT44 * EXPH
      XPPF = FKK/EXPH
      IF (KEY .EQ. 5) GOTO 551
      IF (NAGAIN .NE. 0 .OR. IPRY .EQ. 0) GOTO 1234
      IF (IZAT(IAT) .LE. 6) GOTO 1234
      IF (KEY .EQ. 1 .AND. IAT .GT. MIAT) GOTO 1234
      IF ((IZAT(IAT).GT.8  .AND. XPPF .LT. 0.5) .OR.
     *    (IZAT(IAT).GT.20 .AND. XPPF .LT. 0.6))
     *    WRITE (LIS2,350) IAT, ATNAME(IAT), FKK, EXPH
  350 FORMAT (' $TEMP WARNING PkHeight of atom',I4,
     *     ' = ', A6,' is',F6.0, ' Expected:', F6.0)
 1234 CONTINUE
      XPP = EXPP2(N) * XPPXX2
      IF (NAGAIN .LT. 0) GOTO 524
      IF (N .GE. MTYPE) GOTO 404
      IF (XPPF .GE. XPP) GOTO 404
      IF (NAGAIN .EQ. 0) THEN
         WRITE (LIS2,351) IAT, ATNAME(IAT), FKK, EXPH, IAT, XPPF
  351    FORMAT (' $TEMP For atom',I4, ' = ', A6,' PK height',F7.0,
     *      ' Expected Pk Height', F7.0 /
     *      ' $TEMP For atom',I4, '   PH/EXPH factor =', F5.2)
         RSUM(N,JNAT) = RSUM(N,JNAT) + 1.
         IZOLD = IZAT(IAT)
         ATOLD = ATNAME(IAT)
         NAGAIN = 1
         ENDIF
      N = N + 1
      IZAT(IAT) = NCELLZ(N)
      ATXYZ(10,IAT) = N
      IF (KEY .EQ. 1 .AND. MIAT.NE.0) THEN
         CALL ATN2CN (ATOLD, AC, NUMB)
         NUMB = IABS (NUMB)
         IF (NUMB .EQ. 0 .OR. NUMB .GE. 999) NUMB = IAT
      ELSE
         NUMB = IAT
         ENDIF
      CALL ATN4CN (ACELTY(N), NUMB, MIATX, ATNAME, MIATX, ATNAME(IAT))
      GOTO 209
  404 IF (NAGAIN .GT. 0) THEN
         WRITE (LIS2,407) IAT, ATNAME(IAT)
  407    FORMAT (' New name',I4, ' = ', A6)
         IF (IAT .LE. NATSN) WRITE (LIS1, 408) IAT, ATOLD, ATNAME(IAT)
  408    FORMAT (/' Original input atom ',I4, ' = ', A6,
     *      ' renamed to ', A6, ' LOW PEAK HEIGHT'/)
         ENDIF
      IF (KEY .EQ. 2) GOTO 422
      SUM1(N,JNAT) = SUM1(N,JNAT) + FKK
      SUM2(N,JNAT) = SUM2(N,JNAT) + EXPH
      SUM3(N,JNAT) = SUM3(N,JNAT) + 1.
  422 CONTINUE
      IF (IMAP .EQ. 1) GOTO 444
      IF (NRECYS .GE. 7) GOTO 444
      IF (NAGAIN .GE. 1 .OR. N .LE. 1 .OR. NRECYR .LE. 2
     *   .OR. PSQ .LE. 0.85 .OR. XPPF .LT. 0.90*EXPP(N)) GOTO 444
      IF (XPPF .LT. 1.20 * QPPQQ) GOTO 444
      XPPNEW = QPPQQ * EXPP2(N-1) / XPPXX
      DO 425 JT = 1,NTYPEZ
      IF (CELATY(JT) .EQ. ACELTY(N-1)) GOTO 426
  425 CONTINUE
  426 CONTINUE
      EXPHNW = EXPH * PH(JT,3) / PH(JTYPE,3)
      XPPFNW = FKK / EXPHNW
      IF (XPPFNW .GT. 1.10 * XPPNEW ) WRITE (LIS2,432)
     *   NRECYR, IAT, ATNAME(IAT), FKK, EXPHNW, XPPFNW, XPPNEW
  432 FORMAT (' $TEMP2 Cy', I3, ' At',I3,' ', A6,
     *  'FKK, N-1: EXPH XPPF lim', 2F6.0, 2F5.2)
      IF (XPPFNW .LT. 1.20 * XPPNEW) GOTO 444
      WRITE (LIS2,433) ATNAME(IAT), FKK, CELATY(JT), EXPHNW
  433 FORMAT (' $TEMP2 Cy  - Atom ', A6,' PkHeight', F6.0,
     *   '  EXpected PkHeight for ', A2, ' is', F6.0)
      III = MAX0 (2, NCELIX(N-1)/20)
      IF (NCELIC(N-1) .GE. III) GOTO 444
      IF ( NRECYR.EQ.3 .AND. XPPF.GT.2.5) GOTO 500
      IF ( NRECYR.EQ.3 .AND. XPPF.GT.1.8 .AND. R2X.LT.0.30) GOTO 500
      IF ( NRECYR.EQ.4 .AND. XPPF.GT.1.5 .AND. R2X.LT.0.20) GOTO 500
      IF ( NRECYR.EQ.5 .AND. XPPF.GT.1.4 .AND. R2X.LT.0.20) GOTO 500
      IF ( NRECYR.EQ.6 .AND. XPPF.GT.1.3 .AND. R2X.LT.0.30) GOTO 500
      IF ( NRECYR.GE.5 .AND. XPPF.GT.1.2 .AND. R2X.LT.0.15) GOTO 500
  444 RETURN
  500 IZOLD = IZAT(IAT)
      ATOLD = ATNAME(IAT)
      NAGAIN = -1
      WRITE (LIS2, FMT='('' $TE UPGRADE? Cy='',I3,'' at '',A6, 1X,
     *   ''PH XPPF R2'',F6.0, 2F5.2 )') NRECYR, ATOLD, FKK, XPPF, R2X
      N = N - 1
      IZAT(IAT) = NCELLZ(N)
      ATXYZ(10,IAT) = N
         IF (KEY .EQ. 1 .AND. MIAT.NE.0) THEN
         CALL ATN2CN (ATOLD, AC, NUMB)
         NUMB = IABS (NUMB)
         IF (NUMB .EQ. 0 .OR. NUMB .GE. 999) NUMB = IAT
      ELSE
         NUMB = IAT
         ENDIF
      CALL ATN4CN (ACELTY(N), NUMB, MIATX, ATNAME, MIATX, ATNAME(IAT))
      GOTO 209
  524 CONTINUE
      IF (XPPF .GE. XPP + 0.20) GOTO 544
      IF (NRECYR .GE. 6 .AND. XPPF .GE. XPP + 0.15) GOTO 544
      WRITE (LIS2, FMT='('' $TE UPGRADE? No: atom name not changed'')')
      IZAT(IAT) = IZOLD
      IZOLD = 0
      ATNAME(IAT) = ATOLD
      N = N + 1
      ATXYZ(10,IAT) = N
      GOTO 444
  544 WRITE (LIS2,547) IAT, ATNAME(IAT)
  547 FORMAT (' $TE UPGRADE! New name for big peak ',I4, ' = ', A6)
      IF (IAT .LE. NATSN) WRITE (LIS1, 548) IAT, ATOLD, ATNAME(IAT)
  548 FORMAT (/' Original input atom ',I4, ' = ', A6,
     *   ' renamed to ', A6, ' HIGH PEAK HEIGHT'/)
      RETURN
  551 CONTINUE
      DELB = DBDP(JTYPE) * ( FKK - QPPQQ * EXPH )
      BNEW = BTEMP + DELB
      IF (BNEW .LT. BNMIN) BNEW = BNMIN
      IF (BNEW .GT. BNMAX) BNEW = BNMAX
      BNEW = 0.25 * (BNEW + BTEMPP)
      ATXYZ(5,IAT) = BNEW
      DELB = BNEW - BTEMP
      IZOLD = NINT ( DELB * 100. )
      RETURN
  603 IF (ISKIP .EQ. 0) RETURN
      DO 614 JNAT = 1, 2
      TSUM1(JNAT) = 0.0
      TSUM2(JNAT) = 0.0
      TSUM3(JNAT) = 0.0
      DO 613 J = 1, NTYPEZ
      TSUM1(JNAT) = TSUM1(JNAT) + SUM1(J,JNAT)
      TSUM2(JNAT) = TSUM2(JNAT) + SUM2(J,JNAT)
      TSUM3(JNAT) = TSUM3(JNAT) + SUM3(J,JNAT)
      IF (SUM3(J,JNAT) .LT. 0.5) THEN
         SUM1(J,JNAT) = 0.0
      ELSE
         SUM1(J,JNAT) = SUM1(J,JNAT) / SUM3(J,JNAT)
         ENDIF
      IF (SUM3(J,JNAT) .LT. 0.5) THEN
         SUM2(J,JNAT) = 0.0
      ELSE
         SUM2(J,JNAT) = SUM2(J,JNAT) / SUM3(J,JNAT)
         ENDIF
  613 CONTINUE
      IF (TSUM3(JNAT) .LT. 0.5) THEN
         TSUM1(JNAT) = 0.0
         TSUM2(JNAT) = 0.0
      ELSE
         TSUM1(JNAT) = TSUM1(JNAT) / TSUM3(JNAT)
         TSUM2(JNAT) = TSUM2(JNAT) / TSUM3(JNAT)
         ENDIF
  614 CONTINUE
      ATOLD = 'all'
      WRITE (LIS2, 643)
  643 FORMAT (/' Statistics for PH (Peak height), EXPH (Expected PH),'
     *  /' (averaged) for input atoms (INPUT) and new peak (PEAKS) with'
     *  /' expected ratios PH/EXPH, and applied factor for Z-reduction')
      WRITE (LIS2, 644) (ACELTY(I), I=1,NTYPEZ), ATOLD
  644 FORMAT (/' For atom types: ', 11(3X, A3))
      WRITE (LIS2, 645) (SUM3(I,1), I=1,NTYPEZ), TSUM3(1)
  645 FORMAT ( ' nr.atoms: INPUT ', 11F6.0)
      WRITE (LIS2, 646) (SUM3(I,2), I=1,NTYPEZ), TSUM3(2)
  646 FORMAT ( '     + new PEAKS ', 11F6.0)
      WRITE (LIS2, 647) (SUM1(I,1), I=1,NTYPEZ), TSUM1(1)
  647 FORMAT (/' aver. PH: INPUT ', 11F6.0)
      WRITE (LIS2, 648) (SUM1(I,2), I=1,NTYPEZ), TSUM1(2)
  648 FORMAT ( '           PEAKS ', 11F6.0)
      WRITE (LIS2, 649) (SUM2(I,1), I=1,NTYPEZ), TSUM2(1)
  649 FORMAT (/' av. EXPH: INPUT ', 11F6.0)
      WRITE (LIS2, 650) (SUM2(I,2), I=1,NTYPEZ), TSUM2(2)
  650 FORMAT ( '           PEAKS ', 11F6.0)
      DO 664 JNAT = 1, 2
      DO 663 J = 1, NTYPEZ
      IF (SUM2(J,JNAT) .GT. 0.1) THEN
         SUM1(J,JNAT) = SUM1(J,JNAT) / SUM2(J,JNAT)
      ELSE
         SUM1(J,JNAT) = 0.0
         ENDIF
  663 CONTINUE
      IF (TSUM2(JNAT) .GT. 0.1) THEN
         TSUM1(JNAT) = TSUM1(JNAT) / TSUM2(JNAT)
      ELSE
         TSUM1(JNAT) = 0.0
         ENDIF
  664 CONTINUE
      WRITE (LIS2, 675) (SUM1(I,1), I=1,NTYPEZ), TSUM1(1)
  675 FORMAT (/'  PH/EXPH::INPUT ', 11F6.2)
      WRITE (LIS2, 676) (SUM1(I,2), I=1,NTYPEZ), TSUM1(2)
  676 FORMAT ( '           PEAKS ', 11F6.2)
      WRITE (LIS2, 691) (EXPP2(I), I=1,NTYPEZ)
  691 FORMAT (/' minimum PH/EXPH'  /' required: INPUT ', 11F6.2)
      RSUM1 = 0.
      RSUM2 = 0.
      DO 693 N = 1, NTYPEZ
      RSUM1 = RSUM1 + RSUM(N,1)
      RSUM2 = RSUM2 + RSUM(N,2)
      EXPP2(N) = XPPXX2 * EXPP2(N)
  693 CONTINUE
      WRITE (LIS2, 695) (EXPP2(I), I=1,NTYPEZ)
  695 FORMAT ( '   for new PEAKS ', 11F6.2)
      WRITE (LIS2, 697) (RSUM(I,1), I=1,NTYPEZ), RSUM1
  697 FORMAT (/' Nr of atoms with'
     *        /' reduced Z INPUT ', 11F6.0)
      WRITE (LIS2, 698) (RSUM(I,2), I=1,NTYPEZ), RSUM2
  698 FORMAT ( ' reduced Z PEAKS ', 11F6.0)
      WRITE (LIS2, FMT='(1X)')
      RETURN
  801 CHOUT=' Expected peakheights will not be looked at!'
      CALL SHOUT
      RETURN
      END
      SUBROUTINE DIRBD
      COMMON /SYSTA/ IFILE(20), IFSTAT(20), ITIME(4),   KEYS(28),
     *               NFNUM,     NLIT,       NCOLN(32),  NCOLL(32),
     *               NFDOT(32), NFDOL(32),  NLUSER(32), FNUM(32),
     *               SWITCH(28)
      LOGICAL        SWITCH
      EQUIVALENCE (LIS1,   IFILE(7))
      EQUIVALENCE (LIS2,   IFILE(8))
      WRITE (LIS1, 1)
      WRITE (LIS2, 1)
   1  FORMAT (/' WARNING !!, INPUT ATOMS ARE NOT AVERAGED WITH'/
     + ' PEAK POSITIONS:  LOOK AT YOUR HEAVY ATOM SHIFTS'/
     + ' ESPECIALLY WITH VERY STRONG "RESIDUAL" PEAKS...'/)
      RETURN
      END
      SUBROUTINE DIRBP
      COMMON /SYSTA/ IFILE(20), IFSTAT(20), ITIME(4),   KEYS(28),
     *               NFNUM,     NLIT,       NCOLN(32),  NCOLL(32),
     *               NFDOT(32), NFDOL(32),  NLUSER(32), FNUM(32),
     *               SWITCH(28)
      LOGICAL        SWITCH
      COMMON /SYSTB/ PROGNM,    PROSNM,     CCODE,      TITLE,
     *               CHIN,      LIT(32),    CHOUT
      CHARACTER      PROGNM *8, PROSNM *6,  CCODE *6,   TITLE *64,
     *               CHIN *80,  LIT *6,     CHOUT *72
      EQUIVALENCE (IATOMS, IFILE(1))
      EQUIVALENCE (LIS1,   IFILE(7)), (LIS2,  IFILE(8))
      EQUIVALENCE (KEYS(27), IMAP)
      COMMON /CRYSA/ CELL(6),  CELLSD(6),   RCELL(6),   VOLUM,
     *               WAVE,     CELALL(10),  AMOLW,      ZET,
     *               NELEC,    F000,        ABSMU,      ICENT,
     *               ILATT,    ISYST,       ILAUE,      IMULT,
     *               IUNIQ,    IPOLA,       NTYPE,      NSYMM,
     *         IRSYMM(3,3,24), TSYMM(3,24), NLATT,      TLATT(3,4),
     *         FRAC2C(3,3),    CART2F(3,3), RRMAT(3,3), SSMAT(3,3)
      PARAMETER (MAXAT=993)
      COMMON /XATXYZ/ X(4,MAXAT), ATXYZ(10,MAXAT), IZAT(MAXAT)
      COMMON /SEARDA/ D2R, DMPIC, DMAXB, DMOUT, DMINB, ANGM(2), MCON,
     *        SEARDX, NPIC, NATIN, NAT, NATS, NATX, NATSN, BOV, IPRY,
     *        NRECY, NRECYR, NRECYS, PSQ, NATREC, KPROG, SCALEX, R2X
      DIMENSION XS(3)
      IF (IMAP .EQ. 6 )  WRITE (LIS1, 1)
      IF (IMAP .EQ. 6 )  WRITE (LIS2, 1)
   1  FORMAT (' Patterson peak coordinates written to file ATOMS')
      WRITE (LIS1, 2)
   2  FORMAT (' List of ten highest Patterson peaks and their vector',
     + ' length',/,
     +        '    Peak  height ',4X,'x',7X,'y',7X,'z    length'/)
      WRITE (LIS2, 3)
   3  FORMAT (' List of Patterson peaks and their vector length'/
     +        '    Peak  height ',4X,'x',7X,'y',7X,'z    length'/)
      NPEAK = 1
      IF (IMAP .EQ. 6) THEN
         CHOUT = 'Patterson peaks (not atoms !)'
         CALL ATOMWA (IATOMS)
         ENDIF
      NAT = MIN0 (NPIC, NATX)
      DO 50 I = 1,NAT
      IF (X(4,I) .LE. 0.5) GOTO 333
      DDIST=999.
      DO 3040 II = NLATT,1,-1
      DO 30 L = 1,3
      XS(L) = X(L,I)-TLATT(L,II)
      IF (ABS(XS(L)) .GT. 0.5) XS(L) = XS(L) - SIGN(1.0,XS(L))
   30 CONTINUE
      DIST = SQRT ( QUAD2 (XS, XS) )
 3040 DDIST = AMIN1(DIST,DDIST)
      IF (I .LE. 10) WRITE (LIS1,  35) NPEAK, X(4,I), XS, DDIST
      WRITE (LIS2,  35) NPEAK, X(4,I), XS, DDIST
   35 FORMAT (I8, F8.0, 3F8.4, F8.2)
      IF (IMAP .EQ. 6) WRITE (IATOMS, 36) NPEAK, XS, X(4,I)
   36 FORMAT ('ATOM   Q', I3, 2X, 3F8.5, ' Peakheight ', F8.0)
      NPEAK = NPEAK + 1
   50 CONTINUE
  333 IF (IMAP .EQ. 6) WRITE (IATOMS, FMT='(''END'')')
      CALL FILCLO( IATOMS, 'KEEP')
      CALL KEPROX
      RETURN
      END
      SUBROUTINE DIRBB
      COMMON /SYSTA/ IFILE(20), IFSTAT(20), ITIME(4),   KEYS(28),
     *               NFNUM,     NLIT,       NCOLN(32),  NCOLL(32),
     *               NFDOT(32), NFDOL(32),  NLUSER(32), FNUM(32),
     *               SWITCH(28)
      LOGICAL        SWITCH
      LOGICAL SWRECY, REN98
      EQUIVALENCE (SWITCH(16), REN98)
      EQUIVALENCE (SWITCH(17), SWRECY)
      EQUIVALENCE (LIS1, IFILE(7)), (LIS2, IFILE(8))
      COMMON /CRYSA/ CELL(6),  CELLSD(6),   RCELL(6),   VOLUM,
     +               WAVE,     CELALL(10),  AMOLW,      ZET,
     +               NELEC,    F000,        ABSMU,      ICENT,
     +               ILATT,    ISYST,       ILAUE,      IMULT,
     +               IUNIQ,    IPOLA,       NTYPE,      NSYMM,
     +         IRSYMM(3,3,24), TSYMM(3,24), NLATT,      TLATT(3,4),
     +         FRAC2C(3,3),    CART2F(3,3), RRMAT(3,3), SSMAT(3,3)
      PARAMETER (MAXAT=993)
      COMMON /XATXYZ/ X(4,MAXAT), ATXYZ(10,MAXAT), IZAT(MAXAT)
      COMMON /ATNAMA/ ATNAME(MAXAT)
      CHARACTER ATNAME *6
      COMMON / / IFRAG(MAXAT), ISYM(MAXAT), IDUM(MAXAT),
     +           DUM(MAXAT),  IBOND(MAXAT*10), JBOND(MAXAT*10),
     +           KANG(64), LANG(64), KLANG(64)
      DIMENSION BLACOM(42000)
      EQUIVALENCE (BLACOM(1), IFRAG(1))
      DIMENSION JCON(MAXAT)
      EQUIVALENCE (JCON(1), ISYM(1))
      COMMON /SEARDA/ D2R, DMPIC, DMAXB, DMOUT, DMINB, ANGM(2), MCON,
     *        SEARDX, NPIC, NATIN, NAT, NATS, NATX, NATSN, BOV, IPRY,
     *        NRECY, NRECYR, NRECYS, PSQ, NATREC, KPROG, SCALEX, R2X
      DIMENSION XS(3), XST(3), XSHIFT(3)
      DIMENSION VEC(3,9)
      LOGICAL ID
      DATA ID / .FALSE. /
      IF (IPRY.NE.0) WRITE (IPRY, 100)
  100 FORMAT (/ ' Interatomic bonding distances  (Angstrom)',
     *   /'       Atom    peak',
     *   /'   No. name  integr.  N1 dist  N2 dist  N3 dist  ....'
     *   /'               x100'/)
      NMAX = 0
      IANG = 1
      DO 500 I = 1, NAT
      IF (IPRY.NE.0 .AND. I.EQ.NATSN+1 .AND. NRECYR .LE. 1)
     *   WRITE (IPRY, FMT = '('' New atoms:'')')
      IF (IPRY.NE.0 .AND. I.EQ.NATSN+1 .AND. NRECYR .GT. 1)
     *   WRITE (IPRY, FMT = '('' New atoms, including those added '',
     *      '' earlier (resorted/renamed):'')')
      RTYI = ATXYZ(7,I)
      N = 0
      II= 0
      DO 430 J = 1, NAT
      DAV2 = (ATXYZ(7,J) + RTYI) **2
      IF (I .EQ. J) ID = .TRUE.
      DO 400 IS = 1, NSYMM
      CALL SYMOP1 (IS, X(1,J), XS)
      DO 400 IC = 1, ICENT
      DO 400 IL = 1, NLATT
      IF (ID) THEN
         ID = .FALSE.
         GOTO 400
         ENDIF
      CALL SYMOP2 (IC, IL, XS, XST)
      CALL DISTSQ (ATXYZ(1,I), XST, DMAXB, XSHIFT, DIST2)
      IF (DIST2 .GT. DAV2) GOTO 400
      IF (I .NE. J) GOTO 350
      II = 1
      IF (DIST2 .LT. 0.04) GOTO 400
      II = 2
  350 N = N + 1
      IDUM(N) = J
      DUM(N) = SQRT(DIST2)
      IF (DUM(N) .LT. 0.85  .AND. IPRY.NE.0) THEN
         WRITE (IPRY, 355) I, J
         IF (I .LE. NATSN .AND. J .LE. NATSN) WRITE (IPRY, 356)
  355    FORMAT (I5, ' short contact with No.' , I4)
  356    FORMAT ('+', 30X, '... possible input error ... ')
         ENDIF
      IF (N .LE. 9) THEN
         CALL KERNAB (XSHIFT, VEC(1,N), 3)
      ELSE
         NMAX = 10
         ENDIF
  400 CONTINUE
  430 CONTINUE
      IFRAG(I) = N
      JCON(I) = 0
      IF (IPRY.NE.0) THEN
      IF (II .EQ. 1) WRITE (IPRY, 441) I
  441 FORMAT (/ I5,' lies on a symmetry element')
      IF (II .EQ. 2) WRITE (IPRY, 442) I
  442 FORMAT (/ I5,' is close to a symm. element')
      IF (N .EQ. 0) WRITE (IPRY, 443) I, ATNAME(I), ATXYZ(8,I)
  443 FORMAT (I5, 2X, A6, F6.0, ' *')
      ENDIF
      IF (N .LE. 0) GOTO 500
      IF (IPRY.NE.0)
     * WRITE (IPRY,445) I, ATNAME(I), ATXYZ(8,I), (IDUM(K),DUM(K),K=1,N)
  445 FORMAT (I5, 2X, A6, F6.0, 1X, 12(I4,F5.2)/ (20X, 12(I4,F5.2)) )
      IF (N .GT. 0) CALL GEOFOB (0, I, 0.0)
      IF (N .LE. 1) GOTO 500
      IF (IANG .GE. 10*MAXAT - 64) GOTO 500
      N = MIN0 (N, 9)
      IFRAG(I) = N
      IBAD= 10*I - 9
      CALL KERNAI (IDUM, IBOND(IBAD), N)
      JCON(I) = IANG
      DO 458 K = 1, N-1
      QK = DUM(K)
      QK2 = QK * QK
      DO 457 L = K+1, N
      QL = DUM(L)
      IF (QL .LT. 0.3 .OR. QK .LT. 0.3) THEN
         JBOND(IANG) = 0
         GOTO 453
         ENDIF
      DO 450 J = 1, 3
  450 XSHIFT(J) = VEC(J,L) - VEC(J,K)
      CALL VMATV1 (XSHIFT, RRMAT, XSHIFT, QJ2)
      JBOND(IANG) = NINT (
     *   ACOS ((QK2 + QL*QL - QJ2) / (2.00001 * QK*QL)) / D2R  )
  453 IANG = IANG + 1
  457 CONTINUE
  458 CONTINUE
  500 CONTINUE
      IF (IPRY.NE.0 .AND. NMAX .GT. 9) WRITE (IPRY, 502)
  502 FORMAT (/
     *  ' Note: angles are calculated only for the first 9 contacts')
      IANG = 1
      IF (.NOT.SWRECY) WRITE (LIS2, FMT='(  '' For distances and angles,
     * plots and COORDINATES see LIS1'')')
      IF (IPRY.NE.0) WRITE (IPRY, 510)
  510 FORMAT (/ ' Interatomic bonding angles'
     *   /'     Atom    peak'
     *   /'  No name  integr. N1-No-N2 ang  N1-No-N3 ang  ......'
     *   /'             x100'/)
      DO 800 I = 1, NAT
      IF (IPRY.NE.0 .AND. I .EQ. NATSN+1 .AND. NRECYR .LE. 1)
     *   WRITE (IPRY, FMT = '('' New atoms:'')')
      IF (IPRY.NE.0 .AND. I .EQ. NATSN+1 .AND. NRECYR .GT. 1)
     *   WRITE (IPRY, FMT = '('' New atoms, including those added '',
     *      '' earlier (resorted/renamed):'')')
      M = 0
      IF (JCON(I) .LE. 0) GOTO 790
      N = IFRAG(I)
      IBAD= 10*I - 9
      CALL KERNAI (IBOND(IBAD), IDUM, N)
      DO 758 K = 1, N-1
      DO 757 L = K+1, N
      IF (JBOND(IANG) .EQ. 0) GOTO 753
      M = M + 1
      KLANG(M) = JBOND(IANG)
      KANG(M) = IDUM(K)
      LANG(M) = IDUM(L)
  753 IANG = IANG + 1
  757 CONTINUE
  758 CONTINUE
  790 IF (M .LE. 0) THEN
         IF (IPRY.NE.0) THEN
           IF (NAT. LE. 99) WRITE (IPRY, 791) I, ATNAME(I), ATXYZ(8,I)
  791      FORMAT (I4, 1X, A6, F6.0, 1X, '*')
           IF (NAT. GT. 99) WRITE (IPRY, 792) I, ATNAME(I), ATXYZ(8,I)
  792      FORMAT (I5, 2X, A6, F6.0, 1X, '*')
           ENDIF
         GOTO 800
         ENDIF
      IF (NAT .LE. 99) THEN
         IF(IPRY.NE.0) WRITE (IPRY, 794) I, ATNAME(I), ATXYZ(8,I),
     *     (KANG(J), I, LANG(J), KLANG(J), J=1,M)
  794    FORMAT (I4, 1X, A6, F6.0, 4(I4,2I3,I4) / (17X, 4(I4,2I3,I4) ))
      ELSE
         IF(IPRY.NE.0) WRITE (IPRY, 795) I, ATNAME(I), ATXYZ(8,I),
     *     (KANG(J), LANG(J), KLANG(J), J=1,M)
  795    FORMAT (I5, 2X, A6, F6.0, 4(3I4, '.') / (19X, 4(3I4, '.') ))
         ENDIF
         CALL GEOFOB (M, I, 0.0)
  800 CONTINUE
      IF (REN98) CALL FCALR2
      RETURN
      END
      SUBROUTINE GEOFOB (KEYGEO, IGEO, ZSCAL)
      COMMON /SYSTA/ IFILE(20), IFSTAT(20), ITIME(4),   KEYS(28),
     *               NFNUM,     NLIT,       NCOLN(32),  NCOLL(32),
     *               NFDOT(32), NFDOL(32),  NLUSER(32), FNUM(32),
     *               SWITCH(28)
      LOGICAL        SWITCH
      LOGICAL REN98, SWRECY
      EQUIVALENCE (SWITCH(16), REN98), (SWITCH(17), SWRECY)
      EQUIVALENCE (LIS1, IFILE(7)), (LIS2, IFILE(8))
      COMMON /SYSTB/ PROGNM,    PROSNM,     CCODE,      TITLE,
     *               CHIN,      LIT(32),    CHOUT
      CHARACTER      PROGNM *8, PROSNM *6,  CCODE *6,   TITLE *64,
     *               CHIN *80,  LIT *6,     CHOUT *72
      COMMON /CRYSA/ CELL(6),  CELLSD(6),   RCELL(6),   VOLUM,
     +               WAVE,     CELALL(10),  AMOLW,      ZET,
     +               NELEC,    F000,        ABSMU,      ICENT,
     +               ILATT,    ISYST,       ILAUE,      IMULT,
     +               IUNIQ,    IPOLA,       NTYPE,      NSYMM,
     +         IRSYMM(3,3,24), TSYMM(3,24), NLATT,      TLATT(3,4),
     +         FRAC2C(3,3),    CART2F(3,3), RRMAT(3,3), SSMAT(3,3)
      COMMON /SEARDA/ D2R, DMPIC, DMAXB, DMOUT, DMINB, ANGM(2), MCON,
     *        SEARDX, NPIC, NATIN, NAT, NATS, NATX, NATSN, BOV, IPRY,
     *        NRECY, NRECYR, NRECYS, PSQ, NATREC, KPROG, SCALEX, R2X
      PARAMETER (MAXAT=993)
      COMMON /XATXYZ/ X(4,MAXAT), ATXYZ(10,MAXAT), IZAT(MAXAT)
      COMMON /ATNAMA/ ATNAME(MAXAT)
      CHARACTER ATNAME *6
      COMMON / / IFRAG(MAXAT), ISYM(MAXAT), IDUM(MAXAT),
     +           DUM(MAXAT),  IBOND(MAXAT*10), JBOND(MAXAT*10),
     +           KANG(64), LANG(64), KLANG(64),
     +           IFOB(3,1000), AFOB(1000), IFOBC(10,100),
     +           WFOB(MAXAT), WFOBA(MAXAT), WFOBC(MAXAT)
      DIMENSION BLACOM(42000)
      EQUIVALENCE (BLACOM(1), IFRAG(1))
      LOGICAL CONT
      DATA CONT / .FALSE. /
      DATA DMINX, DMIN, DMAX, DMAXX / 1.10, 1.20, 1.57, 1.70 /
      DATA EMINX, EMIN, EMAX, EMAXX /  90., 103., 126., 139. /
      DATA WANG / 0.5 /
      DATA DELMIN, DELMAX, ELMIN, ELMAX / .15, .15, 10., 10. /
      DATA LFOB, LFOBC, M / 0, 0, 0 /
      IF (REN98) THEN
         IF (KEYGEO .LT. 0) KEYS(24) = 0
         RETURN
         ENDIF
      IF (CONT) GOTO 111
      CONT = .TRUE.
      CALL KERNZI (0, IFOB, 3000)
      CALL KERNZA (0., AFOB, 1000)
      CALL KERNZI (0, IFOBC, 1000)
      LFOB = 0
      LFOBC = 0
      DELMIN = DMIN - DMINX
      DELMAX = DMAXX - DMAX
      ELMIN = EMIN - EMINX
      ELMAX = EMAXX - EMAX
  111 I = IGEO
      IF (KEYGEO) 517, 117, 317
  117 CONTINUE
      IF (IZAT(I) .GT. 9) RETURN
      M = IFRAG(I)
      IF (LFOB .EQ. 1000) GOTO 138
      DO 137 N = 1, M
      J = IDUM(N)
      IF (IZAT(J) .GT. 9) GOTO 137
      IF (J .LT. I) GOTO 137
      IF (DUM(N).GT.DMIN .AND. DUM(N).LT.DMAX) GOTO 137
      IF (DUM(N) .GT. 1.35) THEN
         FOB = AMIN1 (1.4, (DUM(N)-DMAX) / DELMAX)
      ELSE
         FOB = AMIN1 (1.4, (DMIN-DUM(N)) / DELMIN)
         ENDIF
      IF (FOB .LT. 0.1) GOTO 137
      LFOB = LFOB + 1
      IFOB(1,LFOB) = 0
      IFOB(2,LFOB) = I
      IFOB(3,LFOB) = J
      AFOB(LFOB) = FOB
  137 CONTINUE
  138 IF (M .LE. 4) RETURN
      IF (LFOBC .EQ. 100) RETURN
      M = MIN0 (M, 9)
      LFOBC = LFOBC + 1
      IFOBC(1, LFOBC) = I
      CALL KERNAI (IDUM(1), IFOBC(2, LFOBC), M)
      RETURN
  317 CONTINUE
      IF (IZAT(I) .GT. 9) RETURN
      IF (LFOB .EQ. 1000) RETURN
      M = KEYGEO
      DO 337 N = 1, M
      ANG = KLANG(N)
      IF (ANG.GT.EMIN .AND. ANG.LT.EMAX) GOTO 337
      K = KANG(N)
      IF (IZAT(K) .GT. 9) GOTO 337
      J = LANG(N)
      IF (IZAT(J) .GT. 9) GOTO 337
      IF (ANG .GT. 115) THEN
         FOB = AMIN1 (1.4, (ANG-EMAX) / ELMAX)
      ELSE
         FOB = AMIN1 (1.4, (EMIN-ANG) / ELMIN)
         ENDIF
      IF (FOB .LT. 0.1) GOTO 337
      LFOB = LFOB + 1
      IFOB(1,LFOB) = K
      IFOB(2,LFOB) = I
      IFOB(3,LFOB) = J
      AFOB(LFOB) = FOB * WANG
  337 CONTINUE
      RETURN
  517 CONTINUE
      IF (NAT .LE. 5) RETURN
      IF (NRECY .LE. 11) THEN
         MSKIP = (NAT - NATS + 3) / 2
      ELSE
         MSKIP = NAT / 2
         ENDIF
      SCFOB = SQRT(600./ZSCAL)
      NSKIP = 0
      NSKIPT = 0
  555 CALL KERNZA (0., WFOB, NAT)
      CALL KERNZA (0., WFOBA, NAT)
      CALL KERNZA (0., WFOBC, NAT)
      IOLD = 0
      IF (LFOB+LFOBC .EQ. 0) GOTO 901
      IF (LFOB .EQ. 0) GOTO 578
      DO 577 N = 1, LFOB
      K = IFOB(1,N)
      I = IFOB(2,N)
      J = IFOB(3,N)
      IF (I .EQ. 0 .OR. J .EQ. 0) GOTO 577
      IF (IZAT(I) .GT. 9) THEN
          IFOB(2,N) = 0
          GOTO 577
          ENDIF
      IF (K .EQ. 0) THEN
         WFOB(I) = SQRT(WFOB(I)**2 + AFOB(N)**2)
         WFOB(J) = SQRT(WFOB(J)**2 + AFOB(N)**2)
      ELSE
         WFOBA(I) = SQRT(WFOBA(I)**2 + AFOB(N)**2)
         WFOBA(J) = SQRT(WFOBA(J)**2 + AFOB(N)**2)
         WFOBA(K) = SQRT(WFOBA(K)**2 + AFOB(N)**2)
         IOLD = I
         ENDIF
  577 CONTINUE
  578 CONTINUE
      IF (LFOBC .EQ. 0) GOTO 598
      DO 587 N = 1, LFOBC
      I = IFOBC(1, N)
      IF (I .EQ. 0) GOTO 587
      IF (IZAT(I) .GT. 9) THEN
          IFOBC(1,N) = 0
          GOTO 587
          ENDIF
      WFOBC(I) = WFOBC(I) + 1.
      DO 581 NN = 2, 10
      NNN = IFOBC(NN, N)
      IF (NNN .EQ. 0) GOTO 581
      WFOBC(NNN) = WFOBC(NNN) + 1.
  581 CONTINUE
  587 CONTINUE
  598 CONTINUE
      XMAX = 0.
      IXMAX = 0
      DO 625 I = 1, NAT
      IF (IZAT(I).GT.9) GOTO 625
      XX = (WFOB(I) + WFOBA(I) + WFOBC(I))*SCFOB/X(4,I)
      IF (XX .GT. XMAX) THEN
         XMAX = XX
         IXMAX = I
         ENDIF
  625 CONTINUE
      I = IXMAX
      IF (I .EQ. 0) GOTO 901
      IF (XMAX .LT. 0.9) GOTO 901
      IF (SWRECY .AND. (NRECY .GT. 11 .OR. I .GT. NATSN)) THEN
         IF (ATNAME(I)(1:1) .NE. 'Q')
     *      WRITE (LIS2, FMT='('' Bad geometry: atom'',I4, 1X, A6,
     *      '' GEOFOB ='',F6.2,'' atom rejected :Q'')') I,ATNAME(I),XMAX
         IF (I .LE. NATSN) THEN
            NSKIP = NSKIP + 1
            WRITE (LIS1, FMT='('' Bad geometry: atom'',I4, 1X, A6,
     *          '' is rejected !'' )') I,ATNAME(I)
            ENDIF
         ATNAME(I)(1:1) = 'Q'
         CALL KERC2I (ATNAME(I)(2:2),LEND)
         IF ((LEND.LT.0).OR.(LEND.GT.9)) ATNAME(I)(2:6) = ATNAME(I)(3:6)
         IZAT(I) = 1
         IF (I .LE. NATREC + NSKIPT) NSKIPT = NSKIPT + 1
      ELSE
         IF (I .LE. NATSN .OR. I .LE. (NATREC*10)/11 )
     *      WRITE (LIS1, FMT='('' Bad geometry for atom'',I4, 1X, A6,
     *          ''  but atom is retained !'' )') I,ATNAME(I)
         WRITE (LIS2, FMT='('' Bad geometry: atom'',I4, 1X, A6,
     *   '' GEOFOB ='', F6.2, '' atom retained !!'' )') I,ATNAME(I),XMAX
         ENDIF
      IF (NSKIPT .GE. MSKIP) GOTO 901
      IF (LFOB .EQ. 0) GOTO 638
      NC = 0
      DO 637 N = 1, LFOB
      IF (IFOB(1,N).EQ.I .OR. IFOB(2,N).EQ.I .OR. IFOB(3,N).EQ.I)
     *   IFOB(2,N) = 0
      IF (IFOB(2,N) .EQ. 0) NC = NC + 1
  637 CONTINUE
      IF (NC .EQ. LFOB) LFOB = 0
  638 CONTINUE
      IF (LFOBC .EQ. 0) GOTO 698
      NC = 0
      DO 687 N = 1, LFOBC
      IF (IFOBC(1, N    ) .EQ. I) IFOBC(1, N    ) = 0
      IF (IFOBC(1, N    ) .EQ. 0) THEN
         NC = NC+ 1
         GOTO 687
         ENDIF
      NNC = 0
      DO 681 NN = 2, 10
      IF (IFOBC(NN, N    ) .EQ. I) IFOBC(NN, N    ) = 0
      IF (IFOBC(NN, N    ) .GT. 0) NNC = NNC + 1
  681 CONTINUE
      IF (NNC .GT. 4) GOTO 687
      IFOBC(1,LFOBC) = 0
      NC = NC + 1
  687 CONTINUE
      IF (NC .EQ. LFOBC) LFOBC = 0
  698 CONTINUE
      GOTO 555
  901 CONTINUE
      IF (NSKIPT .GT. NSKIP) THEN
         WRITE (CHOUT, 903) NSKIPT
  903    FORMAT (' Nr of peaks not output because of bad geometry:', I4)
         CALL SHOUT
         ENDIF
      KEYS(24) = NSKIP
      RETURN
      END
      SUBROUTINE FCALR2
      COMMON /SYSTA/ IFILE(20), IFSTAT(20), ITIME(4),   KEYS(28),
     *               NFNUM,     NLIT,       NCOLN(32),  NCOLL(32),
     *               NFDOT(32), NFDOL(32),  NLUSER(32), FNUM(32),
     *               SWITCH(28)
      LOGICAL        SWITCH
      COMMON /SYSTB/ PROGNM,    PROSNM,     CCODE,      TITLE,
     *               CHIN,      LIT(32),    CHOUT
      CHARACTER      PROGNM *8, PROSNM *6,  CCODE *6,   TITLE *64,
     *               CHIN *80,  LIT *6,     CHOUT *72
      EQUIVALENCE (ICRYS, IFILE(3))
      EQUIVALENCE (LIS1, IFILE(7)), (LIS2, IFILE(8))
      EQUIVALENCE (IBINFO, IFILE(11))
      COMMON /CRYSA/ CELL(6),  CELLSD(6),   RCELL(6),   VOLUM,
     *               WAVE,     CELALL(10),  AMOLW,      ZET,
     *               NELEC,    F000,        ABSMU,      ICENT,
     *               ILATT,    ISYST,       ILAUE,      IMULT,
     *               IUNIQ,    IPOLA,       NTYPE,      NSYMM,
     *         IRSYMM(3,3,24), TSYMM(3,24), NLATT,      TLATT(3,4),
     *         FRAC2C(3,3),    CART2F(3,3), RRMAT(3,3), SSMAT(3,3)
      COMMON /CRYSB/ SPGR,     WAVEAT,      CELATY(10)
      CHARACTER      SPGR *16, WAVEAT *2,   CELATY *2
      COMMON /FCALCA/ BP,       BR,       SCALE,    HKLMAX(3), STLMAX,
     *                IZTYPE(10), CELPAR(10), PSQQ, P1SQ,     ITRS(24),
     *        AMULT,  ASYMM,    ALATT,    ASYMCL,   NSYMC,    ASYMC,
     *                HKLX(3,24), IDHKL(24), HCODE, FOBS,     SIG,
     *                STL,      STL2,     ISS,      ENORM,
     *                FP,       PHIP,     FAP,      FBP,      EPSIL,
     *                EPSIL2,   SF2,      SF2P,     FPEXP(2,24)
      DIMENSION FITFO(3), FITFC(2)
      EQUIVALENCE (HCODE, FITFO(1)), (FP, FITFC(1))
      PARAMETER (MAXAT=993)
      COMMON /  / SICO(12500), FF(500,10),  EXPBP(500),   EXPBR(500),
     *            SUMF2(500),  SUMF2P(500), SFAC(13,10), ITAT(MAXAT),
     *            ITAT1(MAXAT), ITAT3(MAXAT),
     *            ACI(3,3,MAXAT), BCI(3,3,MAXAT)
      DIMENSION ECI4(3,3,MAXAT)
      DIMENSION ECI22(3,3,MAXAT)
      DOUBLEPRECISION ECI22, ECI4
      DIMENSION BNEWZ(3, MAXAT)
      EQUIVALENCE (BNEWZ(1,1), ECI22(1,1,1))
      DIMENSION BLACOM(42000)
      EQUIVALENCE (BLACOM(1), SICO(1))
      DIMENSION ACALC(MAXAT), BCALC(MAXAT), BNEWR2(MAXAT)
      COMMON /XATXYZ/ X(4,MAXAT), ATXYZ(10,MAXAT), IZAT(MAXAT)
      COMMON /ATNAMA/ ATNAME(MAXAT)
      CHARACTER ATNAME *6
      COMMON /SEARDA/ D2R, DMPIC, DMAXB, DMOUT, DMINB, ANGM(2), MCON,
     *        SEARDX, NPIC, NATIN, NAT, NATS, NATX, NATSN, BOV, IPRY,
     *        NRECY, NRECYR, NRECYS, PSQ, NATREC, KPROG, SCALEX, R2X
      COMMON /DIRBFA/ NCELTY(10), NCELLZ(10), NCELIN(10), NCELIX(10)
      COMMON /DIRBFB/ ACELTY(10)
      CHARACTER ACELTY *2
      DIMENSION NCELL1(10), NCELL3(10), IZTYP1(10), IZTYP3(10)
      DIMENSION BMINX(3), R2MINX(3)
      DOUBLEPRECISION SUMEO4, SUMTR1, SUMTR2, SUMTR3, SUMEC4, SUMEE4
      PARAMETER (MAXBUF = 198)
      DIMENSION BUFFO(MAXBUF)
      CHARACTER ATTEMP *6
      DIMENSION FFF(10), ADTRIG(24)
      DATA ADTRIG / 24*0.0 /
      PSQQ = PSQ
      SCALE = SCALEX
      WRITE (LIS2, FMT='('' $TEMP SCALE = '', F10.5)') SCALE
      BP = BOV
      BR = BOV
      DELB = 0.50 * BOV
      IF (DELB .LT. 1.) DELB = 1.
      IF (DELB .GT. 3.) DELB = 3.
      WRITE (LIS2, 101)
  101 FORMAT (/' Structure factor calculation for R2 refinement')
      WRITE (LIS2, 102) SCALEX, BOV, DELB
  102 FORMAT (' $TE Data used: Scale =', F9.5, ' Bov:', F6.3,
     *        ' delta-B:', F6.3 /)
      CALL BINIFF (1, IBINFO, 'BINFO', FITFO, NITFO, BUFFO, KENDFO)
      STLMAX = BUFFO(6)
      CALL KERNAB (BUFFO(7), HKLMAX, 3)
      KEYT = 2
      CALL FILINQ (ICRYS, 'CRYSDA', 'FORMATTED', 'INPUT', KINQ)
      CALL KERNZI (0, IZTYPE, 10)
      CALL FCALCI (KEYT, ATXYZ, IZAT, ITAT, NAT)
      CALL FILCLO (ICRYS, 'KEEP')
      BSYMCL = ASYMCL**2
      NREFL  = 0
      IREFL  = 0
      SUMTR1 = 0.
      SUMTR2 = 0.
      SUMTR3 = 0.
      SUMEO4 = 0.
      SUMEC4 = 0.
      SUMEE4 = 0.
      NAT9 = 9 * NAT
      DO 104 I3 = 1,NAT
      DO 104 I2 = 1,3
      DO 104 I1 = 1,3
      ECI22(I1,I2,I3) = 0.
      ECI4(I1,I2,I3) = -0.0001
  104 CONTINUE
      EFBP = 0.
      NCELL1(10) = 0
      NCELL3(1) = 0
      DO 105 M = 1, 9
      NCELL1(M) = NCELLZ(M+1)
      NCELL3(M+1) = NCELLZ(M)
  105 CONTINUE
      DO 106 M = 1, 10
      IZTYP1(M) = 0
      IZTYP3(M) = 0
      IF (NCELLZ(M) .LE. 1) THEN
         NCELL1(M) = 0
         NCELL3(M) = 0
         GOTO 106
         ENDIF
      IF (NCELL1(M) .EQ. 1) NCELL1(M) = 0
      IF (NCELL3(M) .EQ. 1) NCELL3(M) = 0
      IF (NCELLZ(M) .GT. 2 + NCELL1(M) .AND.
     *   10*NCELLZ(M) .GT. 12*NCELL1(M) ) NCELL1(M) = 0
      IF (NCELL3(M) .GT. 2 + NCELLZ(M) .AND.
     *   10*NCELL3(M) .GT. 12*NCELLZ(M) ) NCELL3(M) = 0
  106 CONTINUE
      IF (NRECYS .EQ. 7) WRITE (LIS2, 107) NCELLZ, NCELL1, NCELL3
  107 FORMAT (/' Atomic Z and next LOWER and HIGHER Z :'
     *        /' NCELLZ', 10I3/' NCELL1', 10I3/' NCELL3', 10I3)
      DO 114 J = 1, 10
      IZ = IZTYPE(J)
      IF (IZ .LE. 1) GOTO 114
      DO 110 M = 1, 10
      IF (NCELLZ(M) .EQ. IZ) GOTO 111
  110 CONTINUE
  111 CONTINUE
      IZTYP1(J) = NCELL1(M)
      IZTYP3(J) = NCELL3(M)
  114 CONTINUE
      IF (NRECYS .EQ. 7) WRITE (LIS2, 115) IZTYPE, IZTYP1, IZTYP3
  115 FORMAT (/' Corresponding atom type reference number : '
     *        /' IZTYPE', 10I3/' IZTYP1', 10I3/' IZTYP3', 10I3)
      IZ2 = 0
      NZ2 = 0
      DO 119 I = 1, NAT
      ITAT3(I) = 0
      IF (IZAT(I) .NE. IZ2) THEN
         IZ2 = IZAT(I)
         NZ2 = 1
      ELSE
         NZ2 = NZ2 + 1
         IF (NZ2 .GE. 10) GOTO 119
         ENDIF
      J = ITAT(I)
      IZ3 = IZTYP3(J)
      DO 117 M = 1, NTYPE
      IF (IZTYPE(M) .EQ. IZ3) THEN
         ITAT3(I) = M
         GOTO 119
         ENDIF
  117 CONTINUE
  119 CONTINUE
      IZ2 = 0
      NZ2 = 0
      DO 123 I = NAT, 1, -1
      ITAT1(I) = 0
      IF (IZAT(I) .NE. IZ2) THEN
         IZ2 = IZAT(I)
         NZ2 = 1
      ELSE
         NZ2 = NZ2 + 1
         IF (NZ2 .GE. 10) GOTO 123
         ENDIF
      J = ITAT(I)
      IZ1 = IZTYP1(J)
      DO 121 M = 1, NTYPE
      IF (IZTYPE(M) .EQ. IZ1) THEN
         ITAT1(I) = M
         GOTO 123
         ENDIF
  121 CONTINUE
  123 CONTINUE
      IF (NRECYS .GT. 7) GOTO 128
      NATLIM = MIN0 (NAT, 5)
      WRITE (LIS2, FMT='(/'' Atomic Z, JTYPE, J-, J+ ; top only :'')')
      DO 127 I = 1, NATLIM
      WRITE (LIS2, 125) IZAT(I), ITAT(I), ITAT1(I), ITAT3(I)
  125 FORMAT (' IZAT', I3,' ITAT', I3,' ITAT1', I3,' ITAT3', I3)
  127 CONTINUE
  128 CONTINUE
      STLX = 0.9 * STLMAX
  130 CALL BINIFF (0, IBINFO, 'BINFO', FITFO, NITFO, BUFFO, KENDFO)
      IF (KENDFO.LT.0) GOTO 600
      NREFL = NREFL + 1
      IF (FOBS .LT. 3.0 * SIG) GOTO 130
      CALL HKLC1U (HCODE, HKLX)
      CALL HKLSTL (HKLX, STL, STL2)
      IF (STL.GT.STLX .OR. HKLX(1,1).GT.HKLMAX(1) .OR.
     *      HKLX(2,1).GT.HKLMAX(2) .OR. HKLX(3,1).GT.HKLMAX(3))
     *      GOTO 130
      S = STL * 400. + 1.
      IS = IFIX(S)
      STLDEL = S - FLOAT(IS)
      ISS = NINT(S)
      SF2 = SUMF2 (ISS)
      GNORM = SQRT (SF2)
      EOBS2 = ( SCALE * FOBS / GNORM )**2
      EOBS4 = EOBS2**2
      IF (EOBS2 .LT. 0.01) GOTO 130
      IREFL = IREFL + 1
      DO 135 J=1,NTYPE
      FFF(J) = ( FF(IS,J) + (FF(IS+1,J)-FF(IS,J)) * STLDEL ) / GNORM
  135 CONTINUE
      Q5 = EXP (- DELB * STL2)
      Q4 = 1./ Q5
      CALL HKLEX1 (HKLX, HKLX)
      IF (NSYMM.EQ.1) GOTO 150
      DO 140 J=2,NSYMM
      IF (ITRS(J).EQ.0) GOTO 140
      ADTRIG(J) = HKLX(1,1)*TSYMM(1,J) + HKLX(2,1)*TSYMM(2,J) +
     *            HKLX(3,1)*TSYMM(3,J)
  140 CONTINUE
  150 FAP = 0.0
      FBP = 0.0
      CALL KERNZA (0.0, ACI, NAT9)
      IF (ICENT.EQ.1) CALL KERNZA (0.0, BCI, NAT9)
      DO 250 I=1,NAT
      A1 = 0.
      B2 = 0.
      DO 200 J=1,NSYMM
      TRIG = HKLX(1,J)*ATXYZ(1,I) + HKLX(2,J)*ATXYZ(2,I) +
     *       HKLX(3,J)*ATXYZ(3,I) + ADTRIG(J)
      IF (TRIG.LT.0.0) TRIG = TRIG - 0.00010
      ITRIG = MOD ( IFIX(TRIG * 10000. + 0.5), 10000)
      IF (ITRIG.LE.0) ITRIG = ITRIG + 10000
      A1 = A1 + SICO(ITRIG + 2500)
      IF (ICENT.EQ.1) B2 = B2 + SICO(ITRIG)
  200 CONTINUE
      IJ = ITAT(I)
      TF = ATXYZ(4,I) * EXP (- ATXYZ(5,I) * STL2)
      ACALC(I) = A1 * TF * FFF(IJ)
      FAP = FAP + ACALC(I)
      A1 = A1 * TF
      IF (ICENT .EQ. 1) THEN
         BCALC(I) = B2 * TF * FFF(IJ)
         FBP = FBP + BCALC(I)
         B2 = B2 * TF
         ENDIF
      DO 220 M = 1, 3
      IF (M .EQ. 1) IJ = ITAT1(I)
      IF (M .EQ. 2) IJ = ITAT(I)
      IF (M .EQ. 3) IJ = ITAT3(I)
      IF (IJ .EQ. 0) GOTO 220
      A1FFF = A1 * FFF(IJ)
      ACI(1,M,I) = A1FFF * Q4
      ACI(2,M,I) = A1FFF
      ACI(3,M,I) = A1FFF * Q5
      IF (ICENT .EQ. 2) GOTO 220
      B2FFF = B2 * FFF(IJ)
      BCI(1,M,I) = B2FFF * Q4
      BCI(2,M,I) = B2FFF
      BCI(3,M,I) = B2FFF * Q5
  220 CONTINUE
  250 CONTINUE
      EC2 = BSYMCL * (FAP*FAP + FBP*FBP)
      SUMTR1 = SUMTR1 + (0.99 * EOBS2 - EC2)**2
      SUMTR2 = SUMTR2 + (EOBS2 - EC2)**2
      SUMTR3 = SUMTR3 + (1.01 * EOBS2 - EC2)**2
      SUMEO4 = SUMEO4 + EOBS4
      SUMEC4 = SUMEC4 + EC2*EC2
      SUMEE4 = SUMEE4 + EC2*EOBS2
      DO 307 I = 1, NAT
      DO 305 M = 1, 3
      IF (M .EQ. 1 .AND. ITAT1(I) .EQ. 0) GOTO 305
      IF (M .EQ. 3 .AND. ITAT3(I) .EQ. 0) GOTO 305
      DO 303 L = 1, 3
      EFAP = FAP - ACALC(I) + ACI(L,M,I)
      IF (ICENT .EQ. 1) EFBP = FBP - BCALC(I) + BCI(L,M,I)
      ECI2 = BSYMCL * (EFAP*EFAP + EFBP*EFBP)
      ECI22(L,M,I) = ECI22(L,M,I) + EOBS2 * ECI2
      ECI4(L,M,I) = ECI4(L,M,I) + ECI2 * ECI2
  303 CONTINUE
  305 CONTINUE
  307 CONTINUE
      GOTO 130
  600 CONTINUE
      CALL FILCLO (IBINFO, 'KEEP')
      R1 = SUMTR1 / (SUMEO4 * 0.98)
      R2 = SUMTR2 / SUMEO4
      R3 = SUMTR3 / (SUMEO4 * 1.02)
      WRITE (LIS2, 603) R2, IREFL, R1, R2, R3
  603 FORMAT (//' $TEMP R2 =', F8.3, ' for ', I5, ' reflections'//
     *   ' For SCALE2 * [ 0.99 1.00 1.01 , R2 =:', 3F7.3 //)
      SCALE2 = SUMEC4 / SUMEE4
      R2NEW = 1. - SUMEE4**2 / SUMEO4 / SUMEC4
      WRITE (LIS2, FMT='('' $TEMP R2 : New Re-SCALE2 ='', F6.3,
     *   '' New R2 ='', F6.3)') SCALE2, R2NEW
      SC2K2 = 2. / SCALE2 / SUMEO4
      SCK4 = 1. / SCALE2**2 / SUMEO4
      RDIFF1 = 0.0
      RDIFF3 = 0.0
      SDIFF1 = 0.0
      SDIFF3 = 0.0
      BNEWRD = 0.0
      BNEWRA = 0.0
      DO 627 I = 1, NAT
      DO 615 M = 1, 3
      BMINX(M) = -1.
      R2MINX(M) = 999.
      IF (M .EQ. 1 .AND. ITAT1(I) .EQ. 0) GOTO 615
      IF (M .EQ. 3 .AND. ITAT3(I) .EQ. 0) GOTO 615
      DO 613 L = 1, 3
      ECI4(L,M,I) = 1. + SCK4 * ECI4(L,M,I) - SC2K2 * ECI22(L,M,I)
  613 CONTINUE
      RR1 = ECI4(1,M,I)
      RR2 = ECI4(2,M,I)
      RR3 = ECI4(3,M,I)
      CALL BOOTH (RR1, RR2, RR3, ATXYZ(5,I),
     *   DELB, BBOOTH, R2MIN)
      IF (M .EQ. 2) THEN
         BNEWR2(I) = BBOOTH
         BNEWRD = BNEWRD + BBOOTH - ATXYZ(5,I)
         BNEWRA = BNEWRA + ABS (BBOOTH - ATXYZ(5,I))
      ELSE
         IF (M .EQ. 1) THEN
            IF (RR1 .GT. R2MIN+0.0001) GOTO 6613
            RR12 = RR2 - RR1
            IF (RR12 .LT. 0.0002) GOTO 6613
            RR23 = RR3 - RR2
            IF (RR23 .LT. RR12) RR23 = RR12
            RR01A = AMAX1 (RR12 - RR23 + RR12, 0.0)
            RR01B = RR12 * RR12 / RR23
            RRDEL = - AMIN1 ( RR01A + RR01B, RR12 ) * 0.5
            WRITE (LIS2, FMT='('' $TE extrapol: Rold, A,B,-'', 4F8.4)')
     *         R2MIN, RR01A, RR01B, RRDEL
            R2MIN = R2MIN + RRDEL
         ELSEIF(M .EQ. 3) THEN
            IF (RR3 .GT. R2MIN+0.0001) GOTO 6613
            RR32 = RR2 - RR3
            IF (RR32 .LT. 0.0002) GOTO 6613
            RR21 = RR1 - RR2
            IF (RR21 .LT. RR32) RR21 = RR32
            RR01A = AMAX1 (RR32 - RR21 + RR32, 0.0)
            RR01B = RR32 * RR32 / RR21
            RRDEL = - AMIN1 ( RR01A + RR01B, RR32 ) * 0.5
            WRITE (LIS2, FMT='('' $TE extrapol: Rold, A,B,+'', 4F8.4)')
     *         R2MIN, RR01A, RR01B, RRDEL
            R2MIN = R2MIN + RRDEL
            ENDIF
         ENDIF
 6613 BMINX(M) = BBOOTH
      R2MINX(M) = R2MIN
      IZ = IZAT(I)
      IF (M .EQ. 1) II = ITAT1(I)
      IF (M .EQ. 3) II = ITAT3(I)
      IF (M .NE. 2) IZ = IZTYPE(II)
      IF (M .EQ. 2 .OR. R2MIN .LT. R2NEW - 0.001)
     *   WRITE (LIS2, 614) I, ATNAME(I), IZ, ( ECI4(L,M,I), L=1,3 ),
     *   ATXYZ(5,I), BBOOTH, R2MIN
  614 FORMAT (I4, 1X, A6, I2,' R2 for B123', 3F6.3,
     *   ' B,Bnew, R2', 2F5.2, F6.3)
  615 CONTINUE
      IF (ECI4(1,2,I) .GT. 0.0) THEN
         RDIFF1 = RDIFF1 - ECI4(1,2,I) + ECI4(2,2,I)
         SDIFF1 = SDIFF1 + 1.0
         ENDIF
      IF (ECI4(3,2,I) .GT. 0.0) THEN
         RDIFF3 = RDIFF3 + ECI4(2,2,I) - ECI4(3,2,I)
         SDIFF3 = SDIFF3 + 1.0
         ENDIF
      BNEWZ(1,I) = -999.
      IF (BMINX(1) .LT. 0. .AND. BMINX(3) .LT. 0.) GOTO 624
      M = 0
      IF (BMINX(1) .GT. 0. .AND. BMINX(3) .GT. 0.) THEN
         IF (R2MINX(1) .LT. R2MINX(3)) THEN
            IF (R2MINX(1) .LT. R2MINX(2)) M = 1
         ELSE
            IF (R2MINX(3) .LT. R2MINX(2)) M = 3
            ENDIF
      ELSE
         IF (R2MINX(1) .LT. R2MINX(2)) M = 1
         IF (R2MINX(3) .LT. R2MINX(2)) M = 3
         ENDIF
      IF (M .EQ. 0) GOTO 624
      IF (M .EQ. 1) II = ITAT1(I)
      IF (M .EQ. 3) II = ITAT3(I)
      R2DEL = ECI4(2,2,I) - R2MINX(M)
      IF (R2DEL .LT. 0.001) GOTO 624
      BNEWZ(1,I) = FLOAT(II)
      BNEWZ(2,I) = R2DEL
      BNEWZ(3,I) = BMINX(M)
      IZ = IZTYPE(II)
      WRITE (LIS2, FMT='('' $TE new Z? IZ, R2DEL, B'', I3, 2F7.3)')
     *   IZ, R2DEL, BMINX(M)
  624 CONTINUE
  627 CONTINUE
      IF (SDIFF1 .GT. 0.5) RDIFF1 = RDIFF1 / SDIFF1
      IF (SDIFF3 .GT. 0.5) RDIFF3 = RDIFF3 / SDIFF3
      WRITE (LIS2, FMT='('' $TE reduction in R2 by DELZ :'',2F8.4//)')
     *   RDIFF1, RDIFF3
      BNEWRD = BNEWRD / FLOAT(NAT)
      BNEWRA = BNEWRA / FLOAT(NAT)
      WRITE (LIS2, FMT='('' $TE Aver.delB., Av abs delB'', 2F8.4//)')
     *   BNEWRD, BNEWRA
      XXX = AMAX1 (0.1, 0.05*BOV)
      IF (XXX .GT. BNEWRA) XXX = BNEWRA
      RDIFF = 0.5 * (RDIFF1 + RDIFF3)
      WRITE (LIS1, FMT='('' '')')
      TDIFF1 = 0.00001
      TDIFF3 = 0.00001
      DO 628 I = 1, NAT
      IF (ECI4(1,2,I) .GT. 0.0)
     *   TDIFF1 = TDIFF1 + (- ECI4(1,2,I) + ECI4(2,2,I) - RDIFF1 )
      IF (ECI4(3,2,I) .GT. 0.0)
     *   TDIFF3 = TDIFF3 + (- ECI4(3,2,I) + ECI4(2,2,I) - RDIFF3 )
  628 CONTINUE
      IF (SDIFF1 .GT. 1.5) TDIFF1 =  (TDIFF1 / SDIFF1)
      IF (SDIFF3 .GT. 1.5) TDIFF3 =  (TDIFF3 / SDIFF3)
      WRITE (LIS2, FMT='('' $TE reduction: esd  by DELZ :'',2F8.4//)')
     *   TDIFF1, TDIFF3
      DO 649 I = 1, NAT
      IF (BNEWZ(1,I) .LT. 0.0) GOTO 630
      IF (BNEWZ(2,I) .LT. RDIFF) GOTO 630
      II = NINT(BNEWZ(1,I))
      IZ = IZTYPE(II)
      IF (IZ .LT. IZAT(I) .AND. BNEWZ(2,I) .GT. 2.9 * TDIFF1) GOTO 640
      IF (IZ .GT. IZAT(I) .AND. BNEWZ(2,I) .GT. 2.9 * TDIFF3) GOTO 640
  630 BDELR2 = (BNEWR2(I) - ATXYZ(5,I) - BNEWRD) * XXX / BNEWRA
      TEMPAB = ATXYZ(5,I)
      ATXYZ(5,I) = 0.25 * (3.*ATXYZ(5,I) + BDELR2 + BOV)
      WRITE (LIS2, 632) I, ATNAME(I), TEMPAB, BNEWR2(I), ATXYZ(5,I)
  632 FORMAT (I4, 1X, A6, ' B values: old, Booth, new =', 3F6.3)
      GOTO 649
  640 CONTINUE
      IZAT(I) = IZ
      ATXYZ(5,I) = 0.25 * (BNEWZ(3,I) + 3.* BOV )
      ATTEMP = ATNAME(I)
      CALL ATN4CN (CELATY(II), I, 0, ATNAME, 1, ATNAME(I))
      WRITE (LIS2, 643) I, ATTEMP, ATNAME(I), IZ, BNEWZ(2,I), BNEWZ(3,I)
  643 FORMAT (I4, 1X, 2A6, I2,' New name with R2-delta= ', F6.3,
     *   ' Bnew=', F5.2 )
      WRITE (LIS1, 644) I, ATTEMP, ATNAME(I), IZ
  644 FORMAT (' Note: Atom nr',I4,' ', A6, ' is renamed to ',
     *   A6, '(Z =',I2, ') based upon R2')
  649 CONTINUE
      CALL ATOMOC (2, ATXYZ, ITAT, NAT)
      RETURN
      END
      SUBROUTINE BOOTH (R1, R2, R3, X2, DELX, XM, RM)
      X1 = X2 - DELX
      X3 = X2 + DELX
      C = ( R1 +R3 -R2 -R2 ) / (2. * DELX * DELX)
      B = ( R3 - R1 ) / (2. * DELX) - C * 2. * X2
      IF (C .LT. 0.0001) THEN
         XM = X2
         RM = R2
         IF (B .GT. 0.0001 .AND. R1 .LT. R2) THEN
            XM = X1
            RM = R1
         ELSEIF (B .LT. -0.0001 .AND. R3 .LT. R2) THEN
            XM = X3
            RM = R3
            ENDIF
         GOTO 111
         ENDIF
      XM = - 0.5 * B / C
      IF (XM .GT. X3) THEN
         XM = X3
         RM = R3
      ELSEIF (XM .LT. X1) THEN
         XM = X1
         RM = R1
      ELSE
         RM = R2 + B * (XM-X2) + C * (XM*XM-X2*X2)
         ENDIF
  111 CONTINUE
      IF (XM .LT. 0.5) XM = 0.5
      RETURN
      END
      SUBROUTINE CLSTRS (LIS2)
      COMMON /CRYSA/ CELL(6),  CELLSD(6),   RCELL(6),   VOLUM,
     +               WAVE,     CELALL(10),  AMOLW,      ZET,
     +               NELEC,    F000,        ABSMU,      ICENT,
     +               ILATT,    ISYST,       ILAUE,      IMULT,
     +               IUNIQ,    IPOLA,       NTYPE,      NSYMM,
     +         IRSYMM(3,3,24), TSYMM(3,24), NLATT,      TLATT(3,4),
     +         FRAC2C(3,3),    CART2F(3,3), RRMAT(3,3), SSMAT(3,3)
      PARAMETER (MAXAT=993)
      COMMON /XATXYZ/ X(4,MAXAT), ATXYZ(10,MAXAT), IZAT(MAXAT)
      COMMON / / IFRAG(MAXAT), ISYM(MAXAT), IDUM(MAXAT),
     +           DUM(MAXAT),  IBOND(MAXAT*10), JBOND(MAXAT*10),
     +           XGEO
      DIMENSION BLACOM(42000)
      EQUIVALENCE (BLACOM(1), IFRAG(1))
      DIMENSION LW(MAXAT), LR(MAXAT)
      EQUIVALENCE (LW(1), IDUM(1)), (LR(1), DUM(1))
      COMMON /SIZEX/ KFRAG(20), NFRAG, LFRAG(20), NOFRAG, NNA
      COMMON /SEARDA/ D2R, DMPIC, DMAXB, DMOUT, DMINB, ANGM(2), MCON,
     *        SEARDX, NPIC, NATIN, NAT, NATS, NATX, NATSN, BOV, IPRY,
     *        NRECY, NRECYR, NRECYS, PSQ, NATREC, KPROG, SCALEX, R2X
      DIMENSION   IB(3), XS(3), X1(3), XSTOR(3)
      LOGICAL ISWFRM
      DATA    ISWFRM  /.FALSE./
      DFRG = 2.80
      IF (DMAXB .GT. DFRG) DFRG = DMAXB
      DFRG2 = DFRG * DFRG
      DMAXB2 = DMAXB * DMAXB
      ICHECK=0
      GOTO 1010
 1000 NAT = NPIC
 1010 NNA = 0
      MCON=0
      DO 1020 I=1,NAT
      IFRAG(I)=0
      LW(I) = 0
 1020 LR(I) = 0
      NFRAG=0
      NOFRAG=0
      NAT1=NAT-1
      DO 1180 II=1,NAT1
      IF (IFRAG(II).EQ.(-1000)) GOTO 1030
      IF (ICHECK.EQ.1.OR.IFRAG(II).NE.0) GOTO 1180
      IF (NOFRAG .LT. 20) THEN
         NOFRAG = NOFRAG + 1
         NFRAG = NFRAG + 1
         ENDIF
      IF (NOFRAG .EQ. 20 .AND. .NOT. ISWFRM) THEN
         WRITE (LIS2, FMT='('' Warning: 20 or more fragments '')')
         ISWFRM = .TRUE.
         ENDIF
 1030 ICHECK=0
      IFRAG(II)=NOFRAG
      KOUNT=1
      I=II
      IBEGIN=II
 1040 DO 1140 J=IBEGIN,NAT
      IF (IFRAG(J).LT.0.AND.IFRAG(J).NE.(-1000)) GOTO 1140
      JMOVE = 0
      KSYM=0
      DO 1120 K =1, IMULT
      IF (I.EQ.J.AND.K.EQ.1) GOTO 1120
      IF (JSYMM(I,J,K,IB,XS,X1).NE.0) GOTO 1120
      DIST2 = QUAD2 (X1, X1)
      IF (DIST2 .GT. DFRG2) GOTO 1120
      LR(J)=1
      IF (DIST2 .GT. DMAXB2) GOTO 1120
      IF (DIST2 .GT. (ATXYZ(7,I) + ATXYZ(7,J)) **2 ) GOTO 1120
      IF (I .EQ. J .AND. DIST2 .LT. 0.04) GOTO 1120
      IF (IABS(IFRAG(J)).EQ.NOFRAG) GOTO 1080
      IF (IFRAG(J).EQ.(-1000)) GOTO 1060
      JMOVE=1
      KSYM=K
      DO 1050 L=1,3
 1050 XSTOR(L)=XS(L)
 1060 IFRAG(J)=NOFRAG
      KOUNT=KOUNT+1
      IF (KSYM.EQ.K) GOTO 1100
 1080 JBND=100000*IB(1)+10000*IB(2)+1000*IB(3)+K
      IF (JBND.EQ.555001) GOTO 1100
      IF (NNA .EQ. MAXAT-1) WRITE (LIS2, 1085) NOFRAG, MAXAT-1
 1085 FORMAT (/' CLUSTER ',I3,' BONDS ',I4,' TIMES TO ITSELF')
      IF (NNA .LE. MAXAT-1) NNA = NNA + 1
      IF (I.LT.J) ISYM(NNA)=250000*NOFRAG+500*I+J
      IF (I .GE. J) ISYM(NNA) = 250000*NOFRAG + 500*J + I
      GOTO 1120
 1100 IF (I.EQ.J) GOTO 1120
      IDIST = 1000.0 * SQRT(DIST2) + 0.5
      MCON=MCON+1
      IBOND(MCON)=(512*I+J)*8192+IDIST
      MCON=MCON+1
      IBOND(MCON)=(512*J+I)*8192+IDIST
 1120 CONTINUE
      IF (JMOVE .EQ. 0) GOTO 1140
      DO 1125 L=1,3
      X(L,J) = XSTOR(L)
 1125 CONTINUE
 1140 CONTINUE
      IFRAG(I)=-IFRAG(I)
      DO 1160 I=1,NAT
      IF (IFRAG(I).EQ.NOFRAG)GOTO 1040
 1160 CONTINUE
      IF (NAT .EQ. NPIC) GOTO 1178
      DO 1175 I=1,NAT
      IF (KOUNT.LT.4) GOTO 1170
      IF (IABS(IFRAG(I)) .NE. NOFRAG) GOTO 1165
      IF (LW(I) .EQ. 0) GOTO 1170
      GOTO 1000
 1165 IF (LW(I) .EQ. 0) LW(I) = LR(I)
 1170 LR(I) = 0
 1175 CONTINUE
 1178 IF (KOUNT.EQ.1)GOTO 1179
      KFRAG(NOFRAG)=KOUNT
      GOTO 1180
 1179 IFRAG(II)=0
      NFRAG=NFRAG-1
      NOFRAG=NOFRAG-1
 1180 CONTINUE
      DO 1240 I=1,NAT
      IFRAG(I)=IABS(IFRAG(I))
 1240 CONTINUE
      IF (MCON.GT.0) CALL ISORT(IBOND,MCON)
      RETURN
      END
      SUBROUTINE PICTUR (NPROJ)
      COMMON /CRYSA/ CELL(6),  CELLSD(6),   RCELL(6),   VOLUM,
     *               WAVE,     CELALL(10),  AMOLW,      ZET,
     *               NELEC,    F000,        ABSMU,      ICENT,
     *               ILATT,    ISYST,       ILAUE,      IMULT,
     *               IUNIQ,    IPOLA,       NTYPE,      NSYMM,
     *         IRSYMM(3,3,24), TSYMM(3,24), NLATT,      TLATT(3,4),
     *         FRAC2C(3,3),    CART2F(3,3), RRMAT(3,3), SSMAT(3,3)
      COMMON /SYSTA/ IFILE(20), IFSTAT(20), ITIME(4),   KEYS(28),
     *               NFNUM,     NLIT,       NCOLN(32),  NCOLL(32),
     *               NFDOT(32), NFDOL(32),  NLUSER(32), FNUM(32),
     *               SWITCH(28)
      LOGICAL        SWITCH
      PARAMETER (MAXAT=993)
      COMMON /XATXYZ/ X(4,MAXAT), ATXYZ(10,MAXAT), IZAT(MAXAT)
      COMMON / / IFRAG(MAXAT), ISYM(MAXAT), IDUM(MAXAT),
     +           DUM(MAXAT),  IBOND(MAXAT*10), JBOND(MAXAT*10),
     +           XGEO
      DIMENSION BLACOM(42000)
      EQUIVALENCE (BLACOM(1), IFRAG(1))
      DIMENSION XA(3,MAXAT), XB(4,MAXAT)
      PARAMETER (MX7=3*MAXAT+1)
      EQUIVALENCE (XA(1,1), IBOND(1)), (XB(1,1), IBOND(MX7))
      COMMON /SEARDA/ D2R, DMPIC, DMAXB, DMOUT, DMINB, ANGM(2), MCON,
     *        SEARDX, NPIC, NATIN, NAT, NATS, NATX, NATSN, BOV, IPRY,
     *        NRECY, NRECYR, NRECYS, PSQ, NATREC, KPROG, SCALEX, R2X
      COMMON /SIZEX/ KFRAG(20), NFRAG, LFRAG(20), NOFRAG, NNA
      CHARACTER *124 CLINE1
      CHARACTER *6 CL
      DIMENSION B(3,3), V(3,3), IND(3), SUM(3), XMAX(3), XMIN(3)
      IF (IPRY.EQ.0)  RETURN
      NUM = KFRAG(NOFRAG)
      IF (NUM .LT. 4) RETURN
      ANUM = NUM
      COSW = (COS(D2R*CELL(5)) - COS(D2R*CELL(6))*COS(D2R*CELL(4))) /
     *                          (SIN(D2R*CELL(6))*SIN(D2R*CELL(4)))
      SINW=SQRT(1.0-COSW**2)
      A11=CELL(1)*SINW*SIN(D2R*CELL(6))
      A21=CELL(1)*COS(D2R*CELL(6))
      A22=CELL(2)
      A23=CELL(3)*COS(D2R*CELL(4))
      A31=CELL(1)*COSW*SIN(D2R*CELL(6))
      A33 = CELL(3) * SIN(D2R*CELL(4))
      CALL KERNZA (0.0, SUM, 3)
      CALL KERNZA (-10000.0, XMAX, 3)
      CALL KERNZA ( 10000.0, XMIN, 3)
      CALL KERNZA (0.0, B, 9)
      DO 108 I=1,NAT
      IF (IFRAG(I).NE.NOFRAG)GOTO 108
      XA(1,I)=X(1,I)*A11
      XA(2,I)=X(1,I)*A21+X(2,I)*A22+X(3,I)*A23
      XA(3,I)=X(1,I)*A31+X(3,I)*A33
      DO 106 J=1,3
      SUM(J)=SUM(J)+XA(J,I)
      DO 106 K=1,3
  106 B(J,K)=B(J,K)+XA(J,I)*XA(K,I)
  108 CONTINUE
      DO 120 J=1,3
      DO 120 K=1,3
  120 B(J,K)=B(J,K)-SUM(J)*SUM(K)/ANUM
      CALL EIGEN(B,V,IND)
      K=0
      DO 180 I=1,NAT
      IF (IFRAG(I).NE.NOFRAG) GOTO 180
      K=K+1
      XB(4,K)=I
      DO 140 J=1,3
      L=IND(J)
      XB(J,K)=XA(1,I)*V(1,L)+XA(2,I)*V(2,L)+XA(3,I)*V(3,L)
      XMAX(J)=AMAX1(XMAX(J),XB(J,K))
  140 XMIN(J)=AMIN1(XMIN(J),XB(J,K))
  180 CONTINUE
      NNN = 0
      N1 = 2
      N2 = 1
      IF (115.0/(XMAX(1)-XMIN(1)) .GE. 2.0/0.254) GOTO 190
      N1 = 1
      N2 = 2
  190 AMAX = AMAX1(XMAX(N2)-XMIN(N2), XMAX(N2+1)-XMIN(N2+1))
      SCALE = AMIN1(115.0/AMAX, 2.5/0.254)
      SCL = 0.254 * SCALE
  200 CALL SORT (XB, MAXAT, NUM, N1)
  210 NNN = NNN + 1
      WRITE (IPRY, 220) NOFRAG
  220 FORMAT (/// ' ', 72('+')// ' Cluster ',I3)
      GOTO (240, 280, 320), NNN
  240 WRITE (IPRY, 260) SCL
  260 FORMAT ('+', 17X, ' Plot on least squares plane,   scale  ='
     +,F6.2,'  cm /A' /)
      GOTO 360
  280 WRITE (IPRY, 300)
  300 FORMAT ('+', 17X, ' Plot on plane orthogonal to l.s. plane' /)
      GOTO 360
  320 WRITE (IPRY, 340)
  340 FORMAT ('+', 17X, ' Plot on most squares plane' )
  360 IX=0
      OFFSET = 2.5
      ALN = 6.
      CLINE1 = ' '
      DO 460 I =1, NUM
      IXREL = 0.1 * ALN * SCALE * (XMAX(N1)-XB(N1,I)) - FLOAT(IX) + 0.5
      IF (IXREL .LE. 0) GOTO 420
      WRITE (IPRY, 402) CLINE1
  402 FORMAT (A124)
      CLINE1 = ' '
      IX = IX + IXREL
      IF (IXREL .EQ. 1) GOTO 420
      DO 410 J = 1, IXREL-1
  410 WRITE (IPRY, 402) CLINE1
  420 IY=SCALE*(XB(N2,I)-XMIN(N2))+OFFSET
      K = XB(4,I) + 0.5
      CALL KERI2C (K, CL, 3)
      IYE = IY + 1
      IF (K .GT. 99) IYE = IY + 2
      IF (K .LE.  9) IYE = IY
      ICL = 1 + IYE - IY
      IF (CLINE1(IY:IYE) .NE. ' ') GOTO 430
      CLINE1(IY:IYE) = CL(1:ICL)
      GOTO 460
  430 CLINE1(IY:IY) = '*'
  460 CONTINUE
      WRITE (IPRY, 402) CLINE1
      IF (NPROJ .GT. NNN) GOTO 500
      IF (NPROJ .GT. 0 .OR. NNN .GE. 2) GOTO 600
      NTEMP = 4-NNN
      IND1 = IND(NTEMP)
      IND2 = IND(NTEMP-1)
      IF (B(IND2,IND2) .GT. 2.0*B(IND1,IND1)) GOTO 600
  500 IF (NNN .EQ. 2) GOTO 540
      IF (N1 .EQ. 2) GOTO 560
      N2 = 3
      GOTO 210
  540 N2 = 2
      IF (N1 .EQ. 3) GOTO 210
  560 N1 = 3
      GOTO 200
  600 CONTINUE
      RETURN
      END
      SUBROUTINE SCHOUT (KEYT, SCALAT, ZSCAL)
      COMMON /SYSTA/ IFILE(20), IFSTAT(20), ITIME(4),   KEYS(28),
     *               NFNUM,     NLIT,       NCOLN(32),  NCOLL(32),
     *               NFDOT(32), NFDOL(32),  NLUSER(32), FNUM(32),
     *               SWITCH(28)
      LOGICAL        SWITCH
      COMMON /SYSTB/ PROGNM,    PROSNM,     CCODE,      TITLE,
     *               CHIN,      LIT(32),    CHOUT
      CHARACTER      PROGNM *8, PROSNM *6,  CCODE *6,   TITLE *64,
     *               CHIN *80,  LIT *6,     CHOUT *72
      EQUIVALENCE (IATOMS, IFILE(1)), (IDDL, IFILE(2))
      EQUIVALENCE (ISPEK, IFILE(4))
      EQUIVALENCE (LIS1, IFILE(7)), (LIS2, IFILE(8))
      EQUIVALENCE (IATOLD, IFILE(10))
      EQUIVALENCE (IATX, KEYS(11))
      EQUIVALENCE (KEYS(27), IMAP)
      LOGICAL MOLEN, NOFREE
      EQUIVALENCE (SWITCH(2), MOLEN)
      EQUIVALENCE (SWITCH(14), NOFREE)
      LOGICAL DMAXCH, SWRECY
      EQUIVALENCE (SWITCH(28), DMAXCH), (SWITCH(17), SWRECY)
      COMMON /CRYSA/ CELL(6),  CELLSD(6),   RCELL(6),   VOLUM,
     *               WAVE,     CELALL(10),  AMOLW,      ZET,
     *               NELEC,    F000,        ABSMU,      ICENT,
     *               ILATT,    ISYST,       ILAUE,      IMULT,
     *               IUNIQ,    IPOLA,       NTYPE,      NSYMM,
     *         IRSYMM(3,3,24), TSYMM(3,24), NLATT,      TLATT(3,4),
     *         FRAC2C(3,3),    CART2F(3,3), RRMAT(3,3), SSMAT(3,3)
      COMMON /CRYSB/ SPGR,     WAVEAT,      CELATY(10)
      CHARACTER      SPGR *16, WAVEAT *2,   CELATY *2
      PARAMETER (MAXAT=993)
      COMMON /XATXYZ/ X(4,MAXAT), ATXYZ(10,MAXAT), IZAT(MAXAT)
      COMMON /ATNAMA/ ATNAME(MAXAT)
      CHARACTER ATNAME *6
      COMMON / / IFRAG(MAXAT), ISYM(MAXAT), IDUM(MAXAT),
     +           DUM(MAXAT),  IBOND(MAXAT*10), JBOND(MAXAT*10),
     +           XGEO
      DIMENSION BLACOM(42000)
      EQUIVALENCE (BLACOM(1), IFRAG(1))
      DIMENSION LW(MAXAT)
      EQUIVALENCE (LW(1), IDUM(1))
      COMMON /SIZEX/ KFRAG(20), NFRAG, LFRAG(20), NOFRAG, NNA
      COMMON /SEARDA/ D2R, DMPIC, DMAXB, DMOUT, DMINB, ANGM(2), MCON,
     *        SEARDX, NPIC, NATIN, NAT, NATS, NATX, NATSN, BOV, IPRY,
     *        NRECY, NRECYR, NRECYS, PSQ, NATREC, KPROG, SCALEX, R2X
      COMMON /DIRBFA/ NCELTY(10), NCELLZ(10), NCELIN(10), NCELIX(10)
      COMMON /DIRBFB/ ACELTY(10)
      CHARACTER ACELTY *2
      CHARACTER * 2   AQQQ
      DIMENSION DXYZM(3)
      CHARACTER *6 CFRAG, RFAC6
      CHARACTER *12 CFRAGX
      CHARACTER *4 RFAC
      LOGICAL SPEK
      SPEK = .FALSE.
      IF (IATX .EQ. 3 .OR. IATX .EQ. 5) SPEK = .TRUE.
      NRFAC = 0
      DMOUT2 = DMOUT * DMOUT
      DO 100 I = 1, 3
  100 DXYZM(I) = DMOUT * RCELL(I)
      IF (IPRY.NE.0) WRITE (IPRY, 101)
  101 FORMAT (/ 1X,71('-')/)
      IF (SWRECY) CALL DDRECY
      CALL FILINQ (ISPEK,  'SPF',   'FORMATTED', 'OUTPUT', KINQ)
      WRITE (ISPEK, 217) CCODE
  217 FORMAT ('TITL  : DIRDIF output for : ',A6)
      WRITE (ISPEK, 219) CELL
  219 FORMAT ('CELL  ',6F10.5)
      WRITE (ISPEK, 221) SPGR
  221 FORMAT ('SPGR  ',A16)
      IF (SPEK) WRITE (LIS2, FMT='(/
     *   '' Output atomic parameter file CCODE.SPF for PLUTON''/)')
      IF (SPEK) WRITE (LIS1, FMT='(/'' TEMP: '',
     *   '' Output atomic parameter file CCODE.SPF for PLUTON''/)')
      CALL FILINQ (IATOMS, 'ATOMS', 'FORMATTED', 'INPUT', KINQ)
      IT = MAX0 (0, 10000 * (ITIME(1)-1900) + 100 * ITIME(2) + ITIME(3))
      R2 = 9.99
      IF (R2X .GT. 0.01) R2 = R2X
      WRITE (CHOUT, 7102) CCODE, IT, KEYS(13), NRECYR, R2, SCALAT
 7102 FORMAT ('ATOMS ', A6, ' < FOUR date', I7,
     *   ' RUN', I4, ' CY=', I3, ' R2X=', F6.3, '  SC=', F10.6 )
      IF (SCALAT .LE. 0.01) CHOUT(60:72) = ' '
      WRITE (IATOMS, FMT = '(A72)') CHOUT
      WRITE (IATOMS, 799) NRECYR
  799 FORMAT ('REMARK OUTPUT FOURIER CYCLE', I3, ' [ R2X = R2(input) ]')
      BP = 0.
      BR = 0.
      CALL LOGRD (IDDL, 'BP', KLOG)
      IF (KLOG .GT. 0) THEN
         BP  = FNUM(3)
         BR  = FNUM(4)
         WRITE(LIS2, 1007) BP,BR
 1007    FORMAT (' Temp.factors from DDLOG: ',
     *        ' Bp =', F6.3, ' Br =', F6.3)
         ENDIF
      CALL FILCLO (IDDL, 'KEEP')
      WRITE(LIS2, 1008) BOV
 1008 FORMAT (' Overall temperat. factor ', ' Bov:', F6.3)
      WRITE(IATOMS, 1009) BOV
 1009 FORMAT ('REMARK Overall (Wilson) temp. factor, BOV=', F6.3)
      IF (BP .LT. .001) BP = BOV
      IF (BR .LT. .001) BR = BOV
      IF (KEYT .LE. 1) GOTO 105
      IF (NRECYR .LE. 1) GOTO 105
      NATT = 1
      IF (NATSN .GT. 0 .AND. NAT .GT. NATSN) NATT = NATSN
  105 IF (NATREC .EQ. 0) NATREC = NAT
      NTYPEZ = NTYPE
  106 IF (ACELTY(NTYPEZ) .EQ. ' ') THEN
         NTYPEZ = NTYPEZ - 1
         GOTO 106
         ENDIF
      NQQQ = 0
      AQQQ = ACELTY(NTYPEZ)
      IF (AQQQ .EQ. 'H') AQQQ = ACELTY(NTYPEZ-1)
      CALL GEOFOB (-1, 0, ZSCAL)
      IF (IPRY.NE.0) WRITE (IPRY, 102)
  102 FORMAT (/' COORDINATES of atoms and (interpreted) peaks', 15X,
     *  'cluster'/
     +  ' N# name PKintegr PKheight    x',8X,'y',8X,'z',6X,'B', 4X,
     +     'number'/ 13X, 'x100  x100'/)
      IILWM = 1.0 + VOLUM / 17.
      CALL ATOMOC (0, ATXYZ, LW, NAT)
      IILW = 0
      DO 127 I = 1, NAT
      LW(I) = (NSYMM * NLATT) / LW(I)
      IILW = IILW + ATXYZ(4,I) * LW(I)
      IF (IILW .GT. IILWM .AND. I .LT. NAT) THEN
         WRITE (LIS2, FMT='('' Limit NAT because at.vol. = 17 Ang'')')
         MAT = I
         GOTO 129
         ENDIF
  127 CONTINUE
      WIILW = VOLUM / FLOAT(IILW)
      WRITE(LIS2,FMT='('' $TE NAT, Volume/atom:'',I4,F5.1)') NAT,WIILW
  129 II = 0
      DO 180 I=1, NAT
      IFRAG(I)=IABS(IFRAG(I))
      CALL KERI2C (IFRAG(I), CFRAG, 2)
      IF (IFRAG(I) .EQ. 0) CFRAG = '0 '
      N=MIN0(IFRAG(I)+1,11)
      CFRAGX = ' '
      CFRAGX(N:N+1) = CFRAG
      RFAC = ' '
      IF (ATXYZ(9,I) .LT. 0.0001) THEN
         WRITE(LIS2, FMT = '('' occ. factor = 0: '', A6)') ATNAME(I)
      ELSE
         IRFAC = IMULT / NINT (ATXYZ(9,I))
         IF (IRFAC .GT. 1) THEN
            CALL KERI2C (IRFAC, RFAC6, 2)
            RFAC(2:4) = RFAC6
            RFAC(4:4) = 'R'
            IF (RFAC(3:3) .EQ. ' ') RFAC(3:3) = ':'
            NRFAC = NRFAC + 1
            ENDIF
         ENDIF
      K = X(4,I) **2 * ZSCAL
      KK = X(4,I) * 10000.
      KK = KK - IFIX ( X(4,I)) * 10000
      KK = ( FLOAT(KK) )**2 / 4000.
      IF (SWRECY .AND. I.LE.NATSN .AND. NRECYR.GE.3 .AND.
     *      IZAT(I).GE.5 .AND. IZAT(I).LE.9 .AND. K.LT.200) THEN
         IF (IPRY .GT. 0) WRITE (IPRY, FMT='('' Following weak atom '',
     *      A6, '' denoted Q !!!'')') ATNAME(I)
         IF (IPRY .EQ. 0) WRITE (LIS1, FMT='('' Initial weak atom '',
     *      A6, '' rejected, i.e. denoted Q !!!'')') ATNAME(I)
         ATNAME(I)(1:1) = 'Q'
         CALL KERC2I (ATNAME(I)(2:2),LEND)
         IF ((LEND.LT.0).OR.(LEND.GT.9)) ATNAME(I)(2:2) = '0'
         IZAT(I) = 1
         NQQQ = NQQQ + 1
         ENDIF
      IF (I.GT.NATSN .AND. NRECYR.GE.3 .AND. I*100/95.GE.NAT .AND.
     *   ATNAME(I)(1:1) .EQ. 'Q' .AND. K.GT.111 .AND. NQQQ.GE.-2) THEN
         IF (IPRY .GT. 0) WRITE (IPRY, FMT=
     *      '('' Following tail (Q) atom accepted  !!!'')')
         ATNAME(I)(1:1) = AQQQ(1:1)
         IF (AQQQ(2:2) .NE. ' ') THEN
            RFAC6 = ATNAME(I)
            ATNAME(I)(1:2) = AQQQ
            ATNAME(I)(3:6) = RFAC6(2:6)
            ENDIF
         IZAT(I) = NCELLZ(NTYPEZ)
         IF (IZAT(I) .EQ. 1) IZAT(I) = NCELLZ(NTYPEZ-1)
         NQQQ = NQQQ - 1
         ENDIF
      IF (KEYT .LE. 1) THEN
         IF (I .LE. NATSN) ATXYZ(5,I) = BP
         IF (I .GT. NATSN) ATXYZ(5,I) = BR
      ELSE
         IF (I.LE.NATSN .AND. ATXYZ(5,I).LT.0.0001) ATXYZ(5,I) = BP
         IF (I.GT.NATSN .AND. ATXYZ(5,I).LT.0.0001) ATXYZ(5,I) = BR
         ENDIF
      IF (I.GT.NATSN .AND. ATNAME(I)(1:1).EQ.'Q' .AND. K.LT.99) GOTO 136
      IF (K .LT. 1) GOTO 136
      IF (II .GE. NATREC) THEN
         ATNAME(I)(1:1) = 'Q'
         CALL KERC2I (ATNAME(I)(2:2),LEND)
         IF ((LEND.LT.0).OR.(LEND.GT.9)) ATNAME(I)(2:6) = ATNAME(I)(3:6)
         IZAT(I) = 1
         ENDIF
  136 IF (IPRY.NE.0) THEN
         WRITE (IPRY, 137) I,ATNAME(I),K,KK, RFAC,(X(J,I),J=1,3),
     *      ATXYZ(5,I), CFRAGX
  137    FORMAT (I4, 1X, A6, 2I6, A4, 3F9.5, F5.2, 1X, A12)
         IF (K .LT. 100 .AND. I .GT. NATSN) WRITE (IPRY, FMT=
     *       '(8X, A6, '' <---- not output ----'')') ATNAME(I)
         ENDIF
      IF (K .LT. 1) GOTO 180
      IF (K .GT. 100 .OR. I .LE. NATSN) THEN
         IF (IMAP .EQ. 1) ATXYZ(5,I) = 0.0
         IF (NOFREE) ATXYZ(5,I) = 0.0
         WRITE (IATOMS, 173)
     *      ATNAME(I), (X(J,I),J=1,3), (ATXYZ(J,I),J=4,5), K, CFRAG
  173    FORMAT ('ATOM   ', A6, 5F9.5, 3X, I6, '$', A6)
         IF (SPEK) WRITE (ISPEK, 174) ATNAME(I), (ATXYZ(J,I),J=1,3)
  174    FORMAT (A6,2X,3F10.5)
         IF (ATNAME(I)(1:1).NE.'Q') THEN
            II = II + 1
            IF (IPRY.NE.0 .AND. SWRECY .AND. II.EQ.NATREC) WRITE (IPRY,
     *         FMT='('' Following peaks not output for recycling'')')
            ENDIF
         ENDIF
  180 CONTINUE
      WRITE (IATOMS, FMT = '(''END'')')
      CALL COPY80 (IATOMS, 'ATOMS', IATOLD, 'ATOLD')
      CALL FILCLO (IATOMS, 'KEEP')
      IF (SPEK) CALL FILCLO (ISPEK, 'KEEP')
      IF (NRFAC .GT. 0 .AND. IPRY.NE.0) WRITE (IPRY, 236) NRFAC
  236 FORMAT (/ ' :R =', ' symmetry reduction factor for',
     *         I3, ' atoms at special positions ')
      NATSN = NATSN - KEYS(24)
      IF (NATSN.NE.NATS .OR. KEYS(23).EQ.0) THEN
         IRUN = KEYS(13)
         WRITE (CHOUT,FMT='(''RUN '',I3,'' NEW   NAT= '',I4,
     *       '' KPROG '', I3)') IRUN, NATSN, KPROG
         CALL LOGWR (IDDL)
         CALL FILCLO (IDDL, 'KEEP')
         ENDIF
      IF (SWRECY) THEN
         NATREC = II
         WRITE (CHOUT, 242) NATREC
  242    FORMAT (' Number of atoms output for Fourier recycling:', I4)
         CALL SHOUT
         ENDIF
      IF (NNA .LE. 0) GOTO 1900
      IF (NNA .GT. 1) CALL ISORT (ISYM, NNA)
      ISYM(NNA + 1) = 0
      JJ=0
      K=0
      DO 1890 I = 1, NNA
      II=ISYM(I)/250000
      IF (IPRY.NE.0 .AND. II.NE.JJ) WRITE (IPRY, 1860) II
 1860 FORMAT (/' Cluster', I3, '  joins to itself through the peak',
     *          'pairs)')
      JJ=II
      K=K+2
      LW(K-1)=MOD(ISYM(I),250000)/500
      LW(K)=MOD(ISYM(I),500)
      IF (K .GT. 2) THEN
         IF (LW(K) .EQ. LW(K-2)) K = K - 2
         ENDIF
      IF (K.LT.12.AND.II.EQ.ISYM(I+1)/250000) GOTO 1890
      IF (IPRY.NE.0) WRITE (IPRY, 1870) (LW(J),J=1,K)
 1870 FORMAT (1X, 6(I7,',',I3))
      K=0
 1890 CONTINUE
 1900 RETURN
      END
      SUBROUTINE DDRECY
      COMMON /SYSTA/ IFILE(20), IFSTAT(20), ITIME(4),   KEYS(28),
     *               NFNUM,     NLIT,       NCOLN(32),  NCOLL(32),
     *               NFDOT(32), NFDOL(32),  NLUSER(32), FNUM(32),
     *               SWITCH(28)
      LOGICAL        SWITCH,    SWRECY
      EQUIVALENCE   (SWITCH(17),SWRECY)
      COMMON /SYSTB/ PROGNM,    PROSNM,     CCODE,      TITLE,
     *               CHIN,      LIT(32),    CHOUT
      CHARACTER      PROGNM *8, PROSNM *6,  CCODE *6,   TITLE *64,
     *               CHIN *80,  LIT *6,     CHOUT *72
      EQUIVALENCE (LIS1, IFILE(7)), (LIS2, IFILE(8))
      EQUIVALENCE (IHELP, IFILE(10))
      COMMON /CRYSA/ CELL(6),  CELLSD(6),   RCELL(6),   VOLUM,
     *               WAVE,     CELALL(10),  AMOLW,      ZET,
     *               NELEC,    F000,        ABSMU,      ICENT,
     *               ILATT,    ISYST,       ILAUE,      IMULT,
     *               IUNIQ,    IPOLA,       NTYPE,      NSYMM,
     *         IRSYMM(3,3,24), TSYMM(3,24), NLATT,      TLATT(3,4),
     *         FRAC2C(3,3),    CART2F(3,3), RRMAT(3,3), SSMAT(3,3)
      PARAMETER (MAXAT=993)
      COMMON /XATXYZ/ X(4,MAXAT), ATXYZ(10,MAXAT), IZAT(MAXAT)
      COMMON /ATNAMA/ ATNAME(MAXAT)
      CHARACTER ATNAME *6
      COMMON / / IFRAG(MAXAT), ISYM(MAXAT), IDUM(MAXAT),
     +           DUM(MAXAT),  IBOND(MAXAT*10), JBOND(MAXAT*10),
     +           XGEO
      DIMENSION BLACOM(42000)
      EQUIVALENCE (BLACOM(1), IFRAG(1))
      COMMON /SEARDA/ D2R, DMPIC, DMAXB, DMOUT, DMINB, ANGM(2), MCON,
     *        SEARDX, NPIC, NATIN, NAT, NATS, NATX, NATSN, BOV, IPRY,
     *        NRECY, NRECYR, NRECYS, PSQ, NATREC, KPROG, SCALEX, R2X
      COMMON /DIRBFA/ NCELTY(10), NCELLZ(10), NCELIN(10), NCELIX(10)
      COMMON /DIRBFB/ ACELTY(10)
      CHARACTER ACELTY *2
      CHARACTER CAT *2
      CHARACTER*4 ATNAM
      LOGICAL FIRST
      DATA JZAT / 1 /
      FIRST = .FALSE.
      IF (NRECYS .GE. 9) RETURN
      IF (NRECYR.LE.1 .OR. (NRECYR.EQ.2 .AND. NATREC-NATSN .GT. 30))
     *   FIRST = .TRUE.
          IF (NRECYR .LE. 1) THEN
         CALL FILCLO (IHELP, 'KEEP')
         CALL XHELP (IHELP, LIS1, 760.)
         CALL FILCLO (IHELP, 'KEEP')
         ENDIF
      NTYPEZ = NTYPE
  106 IF (ACELTY(NTYPEZ) .EQ. ' ') THEN
         NTYPEZ = NTYPEZ - 1
         GOTO 106
         ENDIF
      NCAT = NTYPEZ
      IF (ACELTY(NCAT) .EQ. 'H' .AND. NCAT .GT. 1) NCAT = NCAT - 1
      IF (NCELLZ(NCAT) .GT. 6) RETURN
      CAT = ACELTY(NCAT)
      IRESET = 0
      IF (.NOT. FIRST) GOTO 500
      DO 300 I = NATSN + 1, NATREC
      IF (IRESET .GT. 0) GOTO 250
      IF (IZAT(I) .GT. 18) GOTO 300
      IF (IZAT(I) .LE. 14) GOTO 300
      IF (ATXYZ(8,I) .GT. 0.5 * FLOAT(IZAT(I)) ) GOTO 300
      IRESET = I
      WRITE (LIS2, FMT='('' $TEMP .... P .. S .. CL  low peak >?'')')
  250 IZAT(I) = NCELLZ(NCAT)
      CALL KERC2I (ATNAME(I) (2:2), KEND)
      ATNAM = ATNAME(I) (2:5)
      IF (KEND .GE. 10) ATNAM = ATNAME(I) (3:6)
      ATNAME(I) (1:2) = CAT
      ATNAME(I) (3:6) = ATNAM
      IF (ATNAME(I) (2:2) .EQ. ' ') ATNAME(I) (2:6) = ATNAM
      CALL ATN24X (ATNAME(I), ATNAME, NATREC, ATNAME(I))
  300 CONTINUE
  400 IF (IPRY.NE.0 .AND. IRESET .GT. 0)
     *  WRITE (IPRY, 402) IRESET, CAT
  402 FORMAT (/' Note about the recycling strategy:'/
     *        /' Starting from peak number' , I3,
     *  ' all peaks have been named ', A2/)
      RETURN
  500 NPK = 0
      PK = 0.
      NATREX = MIN0 (NATREC * 12/10 , NAT)
      DO 510 I = NATSN + 1, NATREC
      IF (IZAT(I) .GT. 10 .OR. IZAT(I) .LE. 1) GOTO 510
      NPK = NPK + 1
      PK = PK + ATXYZ(8,I)
      IF (NPK .GT. 1) GOTO 510
      JZAT = IZAT(I)
      CAT = ATNAME(I) (1:2)
      CALL KERC2I (ATNAME(I) (2:2), KEND)
      IF (KEND .LE. 10 .OR. KEND .GE. 37) CAT(2:2) = ' '
  510 CONTINUE
      IF (NPK .LE. 4) RETURN
      PK = PK / FLOAT(NPK)
      DO 600 I = NATSN + 1, NATREC
      IF (IZAT(I) .GT. 18 .OR. IZAT(I) .LT. 10) GOTO 600
      IF (ATXYZ(8,I) .GT. 1.5 * PK) GOTO 600
      IF (IRESET .EQ. 0) IRESET = I
      IZAT(I) = JZAT
      CALL KERC2I (ATNAME(I) (2:2), KEND)
      ATNAM = ATNAME(I) (2:5)
      IF (KEND .GE. 10) ATNAM = ATNAME(I) (3:6)
      ATNAME(I) (1:2) = CAT
      ATNAME(I) (3:6) = ATNAM
      IF (ATNAME(I) (2:2) .EQ. ' ') ATNAME(I) (2:6) = ATNAM
  600 CONTINUE
      GOTO 400
      END
      SUBROUTINE ISORT (N, M)
      DIMENSION N(M)
      INT=2
 1000 INT=INT+INT
      IF (INT.LT.M) GOTO 1000
      INT=MIN0(M,(3*INT)/4-1)
 1020 INT=INT/2
      IFIN=M-INT
      DO 1200 II=1,IFIN
      I=II
      J=I+INT
      IF (N(I).LE.N(J)) GOTO 1200
      IT=N(J)
 1080 N(J)=N(I)
      J=I
      I=I-INT
      IF (I.LE.0) GOTO 1090
      IF (N(I).GT.IT) GOTO 1080
 1090 CONTINUE
      N(J)=IT
 1200 CONTINUE
      IF (INT.NE.1) GOTO 1020
      RETURN
      END
      FUNCTION JSYMM(I,J,K,IB,XS,X1)
      DIMENSION IB(3),XS(3),X1(3)
      COMMON /CRYSA/ CELL(6),  CELLSD(6),   RCELL(6),   VOLUM,
     *               WAVE,     CELALL(10),  AMOLW,      ZET,
     *               NELEC,    F000,        ABSMU,      ICENT,
     *               ILATT,    ISYST,       ILAUE,      IMULT,
     *               IUNIQ,    IPOLA,       NTYPE,      NSYMM,
     *         IRSYMM(3,3,24), TSYMM(3,24), NLATT,      TLATT(3,4),
     *         FRAC2C(3,3),    CART2F(3,3), RRMAT(3,3), SSMAT(3,3)
      PARAMETER (MAXAT=993)
      COMMON /XATXYZ/ X(4,MAXAT), ATXYZ(10,MAXAT), IZAT(MAXAT)
      COMMON /SEARDA/ D2R, DMPIC, DMAXB, DMOUT, DMINB, ANGM(2), MCON,
     *        SEARDX, NPIC, NATIN, NAT, NATS, NATX, NATSN, BOV, IPRY,
     *        NRECY, NRECYR, NRECYS, PSQ, NATREC, KPROG, SCALEX, R2X
      DIMENSION DXYZM(3)
      LOGICAL CONT
      DATA CONT / .FALSE. /
      IF (CONT) GOTO 111
      CONT = .TRUE.
      DM = AMAX1 (2.8, DMAXB)
      DO 101 L = 1, 3
  101 DXYZM(L) = DM * RCELL(L)
  111 JSYMM = 0
      CALL OPER1 (K, XS, X(1,J))
      DO 1080 L=1,3
      IB(L)=5
 1060 X1(L)=X(L,I)-XS(L)
      IF (ABS(X1(L)).LE.0.5)GOTO 1070
      XS(L)=XS(L)+SIGN(1.0,X1(L))
      IB(L)=IB(L)+ISIGN(1,IFIX(2.0*X1(L)))
      GOTO 1060
 1070 IF (ABS(X1(L)) .GT. DXYZM(L)) GOTO 1100
 1080 CONTINUE
      RETURN
 1100 JSYMM = 1
      RETURN
      END
      SUBROUTINE EIGEN (B, V, IND)
      DIMENSION B(3,3),V(3,3),IND(3)
      DO 1020 I=1,3
      DO 1000 J=1,3
      V(I,J)=0.0
 1000 CONTINUE
      V(I,I)=1.0
 1020 CONTINUE
      NNN = 17
 1040 KNT=0
      NNN = NNN - 1
      IND1=1
      IND3=1
      DO 1500 I=1,2
      IP1=I+1
      DO 1460 J=IP1,3
      IF (ABS(B(I,J)) .LT. 0.000001 * B(IND1,IND1)) KNT = KNT + 1
      BIJ=B(I,I)-B(J,J)
      IF (ABS(B(I,J)).LT.ABS(BIJ))GOTO 1100
      IF (ABS(BIJ).EQ.0.0) THEN
         T = 1.
      ELSE
         T=SIGN(1.0,B(I,J)*BIJ)
         ENDIF
      GOTO 1180
 1100 CONTINUE
      T=B(I,J)/BIJ
 1180 CONTINUE
      IF (ABS(T) .GT. 1.E-6) THEN
         G=T/(2.0+2.0*T*T)
         SN=2.0*G/(1.0+G*G)
         CS=1.0-G*SN
      ELSE
         SN = 0.
         CS = 1.
         ENDIF
      DO 1200 K=1,3
      BIK=B(I,K)
      B(I,K)=CS*B(I,K)+SN*B(J,K)
      B(J,K)=CS*B(J,K)-SN*BIK
 1200 CONTINUE
      DO 1220 K=1,3
      BKI=B(K,I)
      B(K,I)=CS*B(K,I)+SN*B(K,J)
      B(K,J)=CS*B(K,J)-SN*BKI
      VKI=V(K,I)
      V(K,I)=CS*V(K,I)+SN*V(K,J)
      V(K,J)=CS*V(K,J)-SN*VKI
 1220 CONTINUE
 1460 CONTINUE
      IF (B(IP1,IP1).GT.B(IND1,IND1))IND1=IP1
      IF (B(IND3,IND3).GT.B(IP1,IP1))IND3=IP1
 1500 CONTINUE
      IF (KNT.LT.3 .AND. NNN.GT.0) GOTO 1040
      IND(1)=IND1
      IND(2)=1
      IND(3)=IND3
      IF (IND1.EQ.IND(2).OR.IND3.EQ.IND(2))IND(2)=2
      IF (IND1.EQ.IND(2).OR.IND3.EQ.IND(2))IND(2)=3
      DET = V(1,1)*(V(2,2)*V(3,3)-V(2,3)*V(3,2))+V(1,2)*(V(2,3)*V(3,1)
     1  -V(2,1)*V(3,3))+V(1,3)*(V(2,1)*V(3,2)-V(2,2)*V(3,1))
      IF (DET .GT. 0.0) RETURN
      DO 1520 I=1,3
      V(I,3) = -V(I,3)
 1520 CONTINUE
      RETURN
      END
      SUBROUTINE PP1
      COMMON /SYSTA/ IFILE(20), IFSTAT(20), ITIME(4),   KEYS(28),
     *               NFNUM,     NLIT,       NCOLN(32),  NCOLL(32),
     *               NFDOT(32), NFDOL(32),  NLUSER(32), FNUM(32),
     *               SWITCH(28)
      LOGICAL        SWITCH
      EQUIVALENCE (IPR1, IFILE(6)), (LIS2, IFILE(8))
      EQUIVALENCE (IFMAP, IFILE(17)), (ISCRA, IFILE(18))
      EQUIVALENCE (KEYS(27), IMAP),   (KEYS(28), IHALF)
      COMMON /FFTDA/ SCALE, MH(3), NPP(3), XLMIN(3), XLMAX(3)
      DIMENSION ITLE(20)
      EQUIVALENCE (SCALOR,ITLE(18))
      PARAMETER (KUSER2=30000, LUSER2=KUSER2/2)
      COMMON / / X(KUSER2)
      COMPLEX YCOM(LUSER2)
      EQUIVALENCE (X(1), YCOM(1))
      DIMENSION BLACOM(42000)
      EQUIVALENCE (BLACOM(1), X(1))
      INTEGER  P1, P2, R, SKIP, RECS,  D(5)
      NSIZE = KUSER2 - 1
      NX = NPP(1)
      NY = NPP(2)
      NZ = NPP(3)
      RECMAX = 32000.
      MAXSIZ = IFIX ( SQRT (RECMAX * 0.25 * FLOAT (NX*NY*NZ)) )
      NSIZE  = MIN0 (NSIZE, MAXSIZ)
      P1 = NSIZE / (2*NY*(MH(3) + 1))
      P2 = NSIZE / (NX*NZ)
      MPASS = MAX0 (1, NX / MAX0 (1, P1) )
      NBYTES = P1 * P2 * (MH(3)+1) * 8 + 8
      WRITE (LIS2, 101) NBYTES
  101 FORMAT (' Scratch file has max.', I6, ' bytes per record')
      IF (KEYS(26) .EQ. 0 .AND. MPASS .GT. 1) WRITE (IPR1, 102) MPASS
  102 FORMAT (' Approximate nr of intermediate ',
     *         'Fourier-transform passes needed:', I3)
      MPASS = MPASS / MAX0 (1, MPASS/4)
      NPASS = 0
      REWIND ISCRA
      R = -MH(1)
  111 IF (R+P1 .GT. MH(1)) P1 = MH(1) + 1 - R
      CALL RDHKL (YCOM, NY, MH(3)+1, P1, R)
      D(1) = 2*(NY*P1*(MH(3)+1))
      D(2) = 2
      D(3) = D(1)
      D(4) = D(1)
      D(5) = 2*NY
      CALL CMPLFT (X(1), X(2), NSIZE, NY, D)
      CALL WRITEY (YCOM, NY, MH(3)+1, P1, P2, ISCRA)
      NPASS = NPASS + 1
      IF (MPASS .GT. 1 .AND. KEYS(26) .LE. 11 .AND.
     *   MOD (NPASS, MPASS) .EQ. 0) WRITE (IPR1, FMT=
     *   '('' Intermediate transform,  pass'', I2)') NPASS
      R = R + P1
      IF (R .LE. MH(1)) GOTO 111
      IF (IMAP.EQ.5) THEN
         SCALE = 500.0 / SCALE
      ELSE
         SCALE = 3000.0 / SCALE
         ENDIF
      SCALOR = SCALE
      REWIND IFMAP
      NYNEW = NY
      WRITE (IFMAP) ITLE, IMAP, IHALF
      IF (IHALF.NE.0) NYNEW = MIN0 (NY-NY/2+3, NY)
      WRITE (IFMAP) NX,NZ,NYNEW,NY
      REWIND ISCRA
      SKIP = 0
      R = 0
      P1 = NSIZE / (2*NY*(MH(3) + 1))
      RECS = (NY - 1)/P2
      NPASS = 0
  200 IF (R+P2 .GT. NY) P2 = NY - R
      CALL READHL (YCOM, NX,NZ/2, P2,MH(1),MH(3),P1,SKIP,RECS,ISCRA)
      IF (R + P2 .LT. NY) REWIND ISCRA
      SKIP = SKIP + 1
      D(1) = NX * NZ * P2
      D(2) = NZ
      D(3) = NZ * NX
      D(4) = 2 * (MH(3) + 1)
      D(5) = 2
      CALL CMPLFT (X(1), X(2), NSIZE, NX, D)
      D(2) = 2
      D(3) = D(1)
      D(4) = D(1)
      D(5) = NZ
      CALL HERMFT (X(1), X(2), NSIZE, NZ/2, D)
      CALL OUTPUT (X, NZ, NX, P2, R, NY)
      NPASS = NPASS + 1
      IF (MPASS .GT. 1 .AND. KEYS(26) .LE. 11 .AND.
     *   MOD (NPASS, MPASS) .EQ. 0) WRITE (IPR1, FMT=
     *   '('' Final transform, pass'', I2)') NPASS
      R = R + P2
      IF (R .LT. NY) GOTO 200
      CALL FILCLO (ISCRA, 'DELETE')
      RETURN
      END
      SUBROUTINE RDHKL (X, NY, NZ, NX, HS)
      INTEGER HS
      COMPLEX X(NY,NZ,NX)
      COMMON /SYSTA/ IFILE(20), IFSTAT(20), ITIME(4),   KEYS(28),
     *               NFNUM,     NLIT,       NCOLN(32),  NCOLL(32),
     *               NFDOT(32), NFDOL(32),  NLUSER(32), FNUM(32),
     *               SWITCH(28)
      LOGICAL        SWITCH
      EQUIVALENCE (LIS2, IFILE(8)), (IBINFF, IFILE(16))
      EQUIVALENCE (KEYS(27), IMAP)
      COMMON /CRYSA/ CELL(6),  CELLSD(6),   RCELL(6),   VOLUM,
     *               WAVE,     CELALL(10),  AMOLW,      ZET,
     *               NELEC,    F000,        ABSMU,      ICENT,
     *               ILATT,    ISYST,       ILAUE,      IMULT,
     *               IUNIQ,    IPOLA,       NTYPE,      NSYMM,
     *         IRSYMM(3,3,24), TSYMM(3,24), NLATT,      TLATT(3,4),
     *         FRAC2C(3,3),    CART2F(3,3), RRMAT(3,3), SSMAT(3,3)
      COMMON /FFTDA/ SCALE, MH(3), NPP(3), XLMIN(3), XLMAX(3)
      PARAMETER (MAXBUF = 198)
      DIMENSION FITFFT(5), BUFFFT(MAXBUF)
      EQUIVALENCE (FITFFT(4),EI), (FITFFT(5),PHI)
      INTEGER H, HM, HL
      DIMENSION IHKLX(3,24), PHX(24), NSYMEX(24), GRI(2,24), TAB(15)
      DIMENSION IHKLI(3), ITSYMM(3,24), NRCO(5)
      DATA NEX / -1 /
      DATA NPASS / 0 /
      D2R = ATAN(1.0) / 45.0
      IF (NEX .GE. 0) GOTO 210
      NEX = 0
      SCALE = 0.0
      CALL KERNZI (0, NRCO, 5)
      DO 110 I=1,15
  110 TAB(I) = SIN (FLOAT(30*I) * D2R)
      DO 113 J=1,NSYMM
      DO 113 I=1,3
  113 ITSYMM(I,J) = NINT (TSYMM(I,J) * 12.0)
  210 CALL BINIFF (1, IBINFF, 'BINFFT', FITFFT, NITFFT, BUFFFT, NEND)
      NPASS = NPASS + 1
      HM = HS + NX - 1
      HL = HS
      IDISC = HL * HM
      IHMAX = MAX0 (IABS(HL), IABS(HM))
      IHMIN = MIN0 (IABS(HL), IABS(HM))
      DO 220 H=1,NX
      DO 220 L=1,NZ
      DO 220 K=1,NY
  220 X(K,L,H) = CMPLX(0.0,0.0)
  320 CALL BINIFF (0, IBINFF, 'BINFFT', FITFFT, NITFFT, BUFFFT, NEND)
      IF (NEND.LT.0) GOTO 500
      IF (NPASS .GT. 1) GOTO 330
      NRCO(1) = NRCO(1) + 1
      IF (EI .LT. 0.0) THEN
         NRCO(2) = NRCO(2) + 1
         GOTO 320
         ENDIF
  330 CALL KERF2I(FITFFT, IHKLI, 3)
      IF (ISYST.GT.4 .AND. ISYST.LT.8) GOTO 360
      INMAX = MAX0 (IABS(IHKLI(1)), IABS(IHKLI(2)), IABS(IHKLI(3)))
      INMIN = MIN0 (IABS(IHKLI(1)), IABS(IHKLI(2)), IABS(IHKLI(3)))
      IF (IDISC.LT.0) GOTO 360
      IF (INMIN.GT.IHMAX .OR. INMAX.LT.IHMIN) GOTO 320
  360 CALL FEXPAN (IHKLI, IHKLX, PHX, NSYMEX, NEXP)
      PHI = PHI * D2R
      EC  = EI * COS(PHI)
      ES  = EI * SIN(PHI)
      DO 380 J=1,NEXP
      H = IHKLX(1,J)
      IF (H.GT.HM .OR. H.LT.HL) GOTO 380
      IF (IABS(H) .GE. NPP(1)/2 .OR. IABS(H) .GT. MH(1)) THEN
         NRCO(3) = NRCO(3) + 1
         GOTO 380
         ENDIF
      K = IHKLX(2,J)
      IF (IABS(K) .GE. NPP(2)/2 .OR. IABS(K) .GT. MH(2)) THEN
         NRCO(4) = NRCO(4) + 1
         GOTO 380
         ENDIF
      L = IHKLX(3,J)
      IF (IABS(L) .GE. NPP(3)/2 .OR. IABS(L) .GT. MH(3)) THEN
         NRCO(5) = NRCO(5) + 1
         GOTO 380
         ENDIF
      NU = 0
      IF (IMAP.LE.0 .OR. IMAP.EQ.2 .OR. IMAP.GE.5) GOTO 370
      NSYMX = NSYMEX(J)
      NU = - IHKLI(1)*ITSYMM(1,NSYMX) - IHKLI(2)*ITSYMM(2,NSYMX)
     *     - IHKLI(3)*ITSYMM(3,NSYMX)
      NU = MOD(NU,12)
  370 IF (NU.LE.0) NU = NU + 12
      XS = TAB(NU)
      XC = TAB(NU+3)
      GRI(1,J) =  XC*EC - XS*ES
      GRI(2,J) = (XS*EC + XC*ES) * PHX(J)
      NEX = NEX + 1
      SCALE = SCALE + SQRT(GRI(1,J)*GRI(1,J)+GRI(2,J)*GRI(2,J))
      NOKO = 0
      IF (H.EQ.0 .AND. L.EQ.0 .AND. K.NE.0) NOKO = NY - K + 1
      H = H - HS + 1
      IF (K.LT.0) K = NY + K
      K = K + 1
      L = L + 1
      X(K,L,H) = CMPLX(GRI(1,J),GRI(2,J))
      IF (NOKO.NE.0) X(NOKO,L,H) = CONJG(X(K,L,H))
  380 CONTINUE
      GOTO 320
  500 CONTINUE
      IF (NX + HS .LE. MH(1)) RETURN
      CALL FILCLO (IBINFF, 'DELETE')
      WRITE (LIS2, 631) NPASS
  631 FORMAT (' Intermediate transforms required ', I3, ' passes')
      WRITE (LIS2, 690) NRCO(1)
  690 FORMAT (' Number of reflections from input file   =',I7)
      IF (NRCO(2) .GT. 0) WRITE (LIS2, 691) NRCO(2)
  691 FORMAT (' of which',I7,' were rejected'/)
      WRITE (LIS2, 692) NEX
  692 FORMAT (' Number of reflections in one hemisphere =',I7)
      IF (NEX .EQ. 0) CALL KERROR ('No reflections found', 0,'RDHKL')
      IF (NRCO(3).GT.0 .OR. NRCO(4).GT.0 .OR. NRCO(5).GT.0)
     *    WRITE (LIS2, 693)
  693 FORMAT (' not included in calculations, because: '/)
      IF (NRCO(3).GT.0) WRITE (LIS2, 694) MH(1), NRCO(3)
  694 FORMAT (8X,' having H greater than ',I3,'  were ',I6/)
      IF (NRCO(4).GT.0) WRITE (LIS2, 695) MH(2), NRCO(4)
  695 FORMAT (8X,' having K greater than ',I3,'  were ',I6/)
      IF (NRCO(5).GT.0) WRITE (LIS2, 696) MH(3), NRCO(5)
  696 FORMAT (8X,' having L greater than ',I3,'  were ',I6/)
      RETURN
      END
      SUBROUTINE FEXPAN (IHKLI, IHKLX, PHX, NSYMEX, NEXP)
      DIMENSION IHKLI(3), IHKLX(3,24), PHX(24), NSYMEX(24)
      COMMON /CRYSA/ CELL(6),  CELLSD(6),   RCELL(6),   VOLUM,
     *               WAVE,     CELALL(10),  AMOLW,      ZET,
     *               NELEC,    F000,        ABSMU,      ICENT,
     *               ILATT,    ISYST,       ILAUE,      IMULT,
     *               IUNIQ,    IPOLA,       NTYPE,      NSYMM,
     *         IRSYMM(3,3,24), TSYMM(3,24), NLATT,      TLATT(3,4),
     *         FRAC2C(3,3),    CART2F(3,3), RRMAT(3,3), SSMAT(3,3)
      DIMENSION IHKLS(3)
      NEXP = 0
      DO 3900 J=1,NSYMM
      CALL VXMATI (IHKLI, IRSYMM(1,1,J), IHKLS)
      B1 = 1.0
      IF (IHKLS(3)) 1590, 1570, 1600
 1570 IF (IHKLS(1)) 1590, 1580, 1600
 1580 IF (IHKLS(2)) 1590, 1600, 1600
 1590 B1 = -1.0
      DO 3610 I=1,3
 3610 IHKLS(I) = -IHKLS(I)
 1600 CONTINUE
      IF (J.EQ.1) GOTO 1630
      DO 3620 JJ=1, NEXP
      IF (IHKLX(1,JJ) .EQ. IHKLS(1) .AND. IHKLX(2,JJ) .EQ. IHKLS(2)
     *   .AND. IHKLX(3,JJ) .EQ. IHKLS(3)) GOTO 3900
 3620 CONTINUE
 1630 NEXP = NEXP + 1
      CALL KERNAI (IHKLS, IHKLX(1, NEXP), 3)
      PHX(NEXP) = B1
      NSYMEX(NEXP) = J
 3900 CONTINUE
      RETURN
      END
      SUBROUTINE WRITEY (X, NY, NZ, NX, SIZE, ISCRA)
      COMPLEX X(NY,NZ,NX)
      INTEGER SIZE, H, P, Q, R
      P = SIZE
      Q = 0
  100 R = Q + 1
      IF (Q+P .GT. NY) P = NY - Q
      Q = Q + P
      WRITE (ISCRA) (((X(K,L,H), H=1,NX), K=R,Q), L=1,NZ)
      IF (Q .LT. NY) GOTO 100
      RETURN
      END
      SUBROUTINE READHL (X,NX,NZ,NY,HMAX,LMAX,SIZE,SKIP,RECS,ISCRA)
      INTEGER HMAX, SIZE, SKIP ,RECS, H, HL, HU, P, Q
      COMPLEX X(NZ,NX,NY)
      LM = LMAX + 1
      P = SIZE
      HU = NX - HMAX
      IF (SKIP .LE. 0) GOTO 200
      DO 100 Q=1,SKIP
      READ (ISCRA)
  100 CONTINUE
  200 IF (HU + P .GT. NX) GOTO 400
      HL = HU + 1
      HU = HU + P
      READ (ISCRA) (((X(L,H,K), H=HL,HU), K=1,NY), L=1,LM)
      IF (RECS .LE. 0) GOTO 200
      DO 300 Q=1,RECS
      READ (ISCRA)
  300 CONTINUE
      GOTO 200
  400 IF (HU .NE. NX) GOTO 700
      HU = 0
  500 IF (HU+P .GT. HMAX+1) P = HMAX + 1 - HU
      HL = HU + 1
      HU = HU + P
      READ (ISCRA) (((X(L,H,K), H=HL,HU), K=1,NY), L=1,LM)
  550 IF (HU .EQ. HMAX + 1) GOTO 800
      IF (RECS .LE. 0) GOTO 500
      DO 600 Q=1,RECS
      READ (ISCRA)
  600 CONTINUE
      GOTO 500
  700 HL = HU + 1
      HU = HU + P - NX
      IF (HU .GT. HMAX+1) HU = HMAX + 1
      READ (ISCRA) (((X(L,H,K), H=HL,NX), (X(L,H,K), H=1,HU),
     *   K=1,NY), L=1,LM)
      GOTO 550
  800 DO 900 H=2,HU
      HL = NX + 2 - H
      DO 900 K=1,NY
      X(1,HL,K) = CONJG(X(1,H,K))
  900 CONTINUE
      HL = HMAX + 2
      HU = NX - HMAX
      IF (HU .LT. HL) GOTO 920
      DO 910 L=1,LM
      DO 910 K=1,NY
      DO 910 H=HL,HU
      X(L,H,K) = CMPLX(0.0,0.0)
  910 CONTINUE
  920 IF (LM .GE. NZ) GOTO 940
      P = LM + 1
      DO 930 K=1,NY
      DO 930 L=P,NZ
      DO 930 H=1,NX
      X(L,H,K) = CMPLX(0.0,0.0)
  930 CONTINUE
  940 RETURN
      END
      SUBROUTINE OUTPUT (X, NZ, NX, NY, Y, NYT)
      REAL X(NZ,NX,NY)
      INTEGER Y
      COMMON /SYSTA/ IFILE(20), IFSTAT(20), ITIME(4),   KEYS(28),
     *               NFNUM,     NLIT,       NCOLN(32),  NCOLL(32),
     *               NFDOT(32), NFDOL(32),  NLUSER(32), FNUM(32),
     *               SWITCH(28)
      LOGICAL        SWITCH
      COMMON /SYSTB/ PROGNM,    PROSNM,     CCODE,      TITLE,
     *               CHIN,      LIT(32),    CHOUT
      CHARACTER      PROGNM *8, PROSNM *6,  CCODE *6,   TITLE *64,
     *               CHIN *80,  LIT *6,     CHOUT *72
      EQUIVALENCE (LIS2, IFILE(8)), (IFMAP, IFILE(17))
      EQUIVALENCE (KEYS(28), IHALF)
      LOGICAL PRIMAP
      EQUIVALENCE (PRIMAP, SWITCH(11))
      COMMON /FFTDA/ SCALE, MH(3), NPP(3), XLMIN(3), XLMAX(3)
      INTEGER SEC, XL, XU
      INTEGER*2 INUT(250)
      DIMENSION LINE(50)
      DATA NCOL /36/
      IF (IHALF.NE.0)THEN
         NSECO = MIN0 (NYT-NYT/2+2, NYT)
         XLMAX(2) = AMIN1 (XLMAX(2), 0.5)
      ELSE
         NSECO = NYT
         ENDIF
      NXM = MIN0 (NX,  INT(XLMAX(1)*FLOAT(NX)) +1)
      NYM = MIN0 (NYT, INT(XLMAX(2)*FLOAT(NYT))+1)
      NZM = MIN0 (NZ,  INT(XLMAX(3)*FLOAT(NZ)) +1)
      ASC = 1000.0 / FLOAT(NX)
      DO 900 K=1,NY
      SEC = Y + K - 1
      NRSEC = SEC + 1
      DO 350 J=1,NZ
      DO 330 I=1,NX
      X(J,I,K) = X(J,I,K) * SCALE
      INUT4 = MIN0 (NINT (X(J,I,K)), 32767)
  330 INUT(I) = INUT4
      IF (NRSEC.LE.NSECO) WRITE (IFMAP) SEC, J, NX, (INUT(I),I=1,NX)
      IF (NRSEC.EQ.NYT .AND. IHALF.NE.0)
     +WRITE (IFMAP) SEC, J, NX, (INUT(I),I=1,NX)
  350 CONTINUE
      IF (.NOT. PRIMAP) GOTO 900
      XU = 0
  450 XL = XU + 1
      XU = XU + NCOL
      IF (XU.GT.NXM) XU = NXM
      DO 460 L=XL,XU
      INUT(L) = FLOAT((L-1))*ASC + 0.5
  460 CONTINUE
      NEC = 1000.*(FLOAT(SEC))/(FLOAT(NYT)) + 0.5
      WRITE (LIS2, 480) NEC, TITLE
  480 FORMAT ('1SECTION  Y = ', I4, 4X, A64 /)
      WRITE (LIS2, 500) (INUT(NN), NN=XL,XU,2)
  500 FORMAT (8X,'X =',I4,17I6)
      ILIM = XL + 1
      WRITE (LIS2, 520) (INUT(NN), NN=ILIM,XU,2)
  520 FORMAT (12X,18I6)
      WRITE (LIS2, 530)
  530 FORMAT ('0')
      DO 700 I=1,NZM
      LINE(1) = FLOAT((I-1))*1000./FLOAT(NZ) + 0.5
      L = 1
      DO 600 J=XL,XU
      L = L + 1
      LINE(L) = NINT (0.1 * X(I,J,K))
      IF (LINE(L) .LT. -99) LINE(L) = -99
  600 CONTINUE
      WRITE (LIS2, 630) (LINE(J), J=1,L)
  630 FORMAT ('0Z =',I4,' *  ',36I3)
  700 CONTINUE
      IF (XU.LT.NXM)GOTO 450
      IF (SEC.GE.NYM) PRIMAP = .FALSE.
  900 CONTINUE
      RETURN
      END
      SUBROUTINE CMPLFT (X, Y, NSIZE, N, D)
      REAL X(NSIZE), Y(NSIZE)
      INTEGER D(5),PMAX,PSYM,TWOGRP,FACTOR(15),SYM(15),UNSYM(15)
      PMAX   = 5
      TWOGRP = 4
      CALL SRFP (N, PMAX, TWOGRP, FACTOR, SYM, PSYM, UNSYM)
      CALL MDFTKD (N, FACTOR, D, X, Y, NSIZE)
      CALL DIPRP (N, SYM, PSYM, UNSYM, D, X, Y, NSIZE)
      RETURN
      END
      SUBROUTINE SRFP (PTS,PMAX,TWOGRP,FACTOR,SYM,PSYM,UNSYM)
      INTEGER PTS,PMAX,TWOGRP,PSYM, FACTOR(15), SYM(15), UNSYM(15)
      INTEGER PP(14), QQ(7), F, P, PTWO, Q, R
      N = PTS
      PSYM = 1
      F = 2
      P = 0
      Q = 0
  100 IF (N.LE.1) GOTO 500
      DO 200 J=F,PMAX
      IF (N.EQ.(N/J)*J) GOTO 300
  200 CONTINUE
      CALL KERNER (200, 'SRFP')
  300 F = J
      N = N / F
      IF (N.EQ.(N/F)*F) GOTO 400
      Q = Q + 1
      QQ(Q) = F
      GOTO 100
  400 N = N / F
      P = P + 1
      PP(P) = F
      PSYM = PSYM * F
      GOTO 100
  500 R = 1
      IF (Q.EQ.0) R = 0
      IF (P.LT.1) GOTO 700
      DO 600 J=1,P
      JJ = P + 1 - J
      SYM(J) = PP(JJ)
      FACTOR(J) = PP(JJ)
      JJ = P + Q + J
      FACTOR(JJ) = PP(J)
      JJ = P + R + J
      SYM(JJ) = PP(J)
  600 CONTINUE
  700 IF (Q.LT.1) GOTO 900
      DO 800 J=1,Q
      JJ = P + J
      UNSYM(J) = QQ(J)
      FACTOR(JJ) = QQ(J)
  800 CONTINUE
      SYM(P+1) = PTS / PSYM**2
  900 JJ = 2*P + Q
      FACTOR(JJ+1) = 0
      PTWO = 1
      J = 0
 1000 J = J + 1
      IF (FACTOR(J).EQ.0) GOTO 1200
      IF (FACTOR(J).NE.2) GOTO 1000
      PTWO = PTWO * 2
      FACTOR(J) = 1
      IF (PTWO .GE. TWOGRP) GOTO 1100
      IF (FACTOR(J+1).EQ.2) GOTO 1000
 1100 FACTOR(J) = PTWO
      PTWO = 1
      GOTO 1000
 1200 IF (P.EQ.0) R = 0
      JJ = 2*P + R
      SYM(JJ+1) = 0
      IF (Q.LE.1) Q = 0
      UNSYM(Q+1) = 0
      RETURN
      END
      SUBROUTINE DIPRP (PTS, SYM, PSYM, UNSYM, DIM, X, Y, NSIZE)
      REAL X(NSIZE), Y(NSIZE)
      INTEGER SYM(15), UNSYM(15), DIM(5), PTS, PSYM, DK, PUNSYM, TEST
      LOGICAL ONEMOD
      INTEGER SEP, DELTA, P, P0, P1, P2, P3, P4, P5, SIZE
      INTEGER V(14), MODULO(14), S(14), U(14)
      DATA MODS / 0 /
      NEST = 14
      NT = DIM(1)
      SEP = DIM(2)
      P2 = DIM(3)
      SIZE = DIM(4) - 1
      P4 = DIM(5)
      IF (SYM(1).EQ.0) GOTO 500
      DO 100 J=1,NEST
      U(J) = 1
      S(J) = 1
  100 CONTINUE
      N = PTS
      DO 200 J=1,NEST
      IF (SYM(J).EQ.0) GOTO 300
      JJ = NEST + 1 - J
      U(JJ) = N
      N = N / SYM(J)
      S(JJ) = N
  200 CONTINUE
  300 JJ = 0
      L = 1
      V(1) = 1
  310 L = L + 1
      V(L) = V(L-1)
  320 IF (L.LT.NEST) GOTO 310
      N = V(NEST)
      JJ = JJ + 1
      IF (JJ.GE.N) GOTO 400
      DELTA = (N-JJ) * SEP
      P1 = (JJ-1)*SEP + 1
      DO 350 P0=P1,NT,P2
      P3 = P0 + SIZE
      DO 350 P=P0,P3,P4
      P5 = P + DELTA
      T = X(P)
      X(P) = X(P5)
      X(P5) = T
      T = Y(P)
      Y(P) = Y(P5)
      Y(P5) = T
  350 CONTINUE
  400 V(L) = V(L) + S(L)
      IF (V(L).LE.U(L)) GOTO 320
      L = L - 1
      IF (L.NE.0) GOTO 400
  500 IF (UNSYM(1).EQ.0) GOTO 1900
      PUNSYM = PTS / PSYM**2
      MULT = PUNSYM / UNSYM(1)
      TEST = (UNSYM(1)*UNSYM(2)-1) * MULT * PSYM
      LK = MULT
      DK = MULT
      DO 600 K=2,NEST
      IF (UNSYM(K).EQ.0) GOTO 700
      LK = LK * UNSYM(K-1)
      DK = DK / UNSYM(K)
      U(K) = (LK-DK) * PSYM
      MODS = K
  600 CONTINUE
  700 ONEMOD = MODS.LT.3
      IF (ONEMOD) GOTO 900
      DO 800 J=3,MODS
      JJ = MODS + 3 - J
      MODULO(JJ) = U(J)
  800 CONTINUE
  900 MODULO(2) = U(2)
      JL = (PUNSYM-3) * PSYM
      MS = PUNSYM * PSYM
      DO 1800 J=PSYM,JL,PSYM
      K = J
 1000 K = K * MULT
      IF (ONEMOD) GOTO 1200
      DO 1100 I=3,MODS
      K = K - (K/MODULO(I))*MODULO(I)
 1100 CONTINUE
 1200 IF (K.GE.TEST) GOTO 1300
      K = K - (K/MODULO(2))*MODULO(2)
      GOTO 1400
 1300 K = K - (K/MODULO(2))*MODULO(2)+MODULO(2)
 1400 IF (K.LT.J) GOTO 1000
      IF (K.EQ.J) GOTO 1800
      DELTA = (K-J) * SEP
      DO 1600 L=1,PSYM
      DO 1500 M=L,PTS,MS
      P1 = (M+J-1)*SEP + 1
      DO 1500 P0=P1,NT,P2
      P3 = P0 + SIZE
      DO 1500 JJ=P0,P3,P4
      KK = JJ + DELTA
      T = X(JJ)
      X(JJ) = X(KK)
      X(KK) = T
      T = Y(JJ)
      Y(JJ) = Y(KK)
      Y(KK) = T
 1500 CONTINUE
 1600 CONTINUE
 1800 CONTINUE
 1900 RETURN
      END
      SUBROUTINE MDFTKD (N, FACTOR, DIM, X, Y, NSIZE)
      INTEGER FACTOR(15), DIM(5), F, P, R, S
      REAL X(NSIZE), Y(NSIZE)
      S = DIM(2)
      F = 0
      M = N
  100 F = F + 1
      P = FACTOR(F)
      IF (P.EQ.0) RETURN
      M = M / P
      R = M * S
      NDIR1 = NSIZE - R
      NDIR2 = NDIR1 - R
      NDIR3 = NDIR2 - R
      NDIR4 = NDIR3 - R
      GOTO (100, 200, 300, 400, 500), P
  200 CALL R2CFTK (N, M, X(1), Y(1), X(R+1), Y(R+1), DIM, NSIZE, NDIR1)
      GOTO 100
  300 CONTINUE
      CALL R3CFTK(N, M, X(1), Y(1), X(R+1), Y(R+1), X(2*R+1), Y(2*R+1)
     ., DIM, NSIZE, NDIR1, NDIR2)
      GOTO 100
  400 CALL R4CFTK (N, M, X(1), Y(1), X(R+1), Y(R+1), X(2*R+1), Y(2*R+1)
     ., X(3*R+1), Y(3*R+1), DIM, NSIZE, NDIR1, NDIR2, NDIR3)
      GOTO 100
  500 CALL R5CFTK (N, M, X(1), Y(1), X(R+1), Y(R+1), X(2*R+1), Y(2*R+1)
     ., X(3*R+1), Y(3*R+1), X(4*R+1), Y(4*R+1), DIM, NSIZE, NDIR1,
     +  NDIR2, NDIR3, NDIR4)
      GOTO 100
      END
      SUBROUTINE R2CFTK (N, M, X0, Y0, X1, Y1, DIM, NDIR0, NDIR1)
      INTEGER DIM(5), SIZE, SEP
      REAL X0(NDIR0), Y0(NDIR0), X1(NDIR1), Y1(NDIR1), IS, IU
      LOGICAL FOLD,ZERO
      DATA TWOPI / 6.2831853 /
      DATA C, S / 0.0, 0.0 /
      NT = DIM(1)
      SEP = DIM(2)
      L1 = DIM(3)
      SIZE = DIM(4) - 1
      K2 = DIM(5)
      NS = N * SEP
      M2 = M * 2
      FM2 = FLOAT(M2)
      MOVER2 = M/2 + 1
      MM2 = SEP * M2
      FJM1 = -1.0
      DO 600 J=1, MOVER2
      FOLD = J.GT.1 .AND. 2*J.LT.M+2
      K0 = (J-1)*SEP + 1
      FJM1 = FJM1 + 1.0
      ANGLE = TWOPI * FJM1 / FM2
      ZERO = ANGLE.EQ.0.0
      IF (ZERO) GOTO 200
      C = COS(ANGLE)
      S = SIN(ANGLE)
      GOTO 200
  100 FOLD = .FALSE.
      K0 = (M+1-J)*SEP + 1
      C = -C
  200 DO 500 KK=K0,NS,MM2
      DO 440 L=KK,NT,L1
      K1 = L + SIZE
      DO 420 K=L,K1,K2
      RS = X0(K) + X1(K)
      IS = Y0(K) + Y1(K)
      RU = X0(K) - X1(K)
      IU = Y0(K) - Y1(K)
      X0(K) = RS
      Y0(K) = IS
      IF (ZERO) GOTO 300
      X1(K) = RU*C + IU*S
      Y1(K) = IU*C - RU*S
      GOTO 420
  300 X1(K) = RU
      Y1(K) = IU
  420 CONTINUE
  440 CONTINUE
  500 CONTINUE
      IF (FOLD) GOTO 100
  600 CONTINUE
      RETURN
      END
      SUBROUTINE R3CFTK (N, M, X0, Y0, X1, Y1, X2, Y2, DIM,
     +     NDIR0, NDIR1, NDIR2)
      REAL X0(NDIR0),Y0(NDIR0),X1(NDIR1),Y1(NDIR1),X2(NDIR2),Y2(NDIR2),
     +      I0,I1,I2,IA,IB,IS
      LOGICAL FOLD,ZERO
      INTEGER DIM(5), SIZE, SEP
      DATA TWOPI / 6.2831853 /
      DATA A/-0.5/, B/0.866/
      DATA C1, S1, C2, S2 / 0.0, 0.0, 0.0, 0.0 /
      NT = DIM(1)
      SEP = DIM(2)
      L1 = DIM(3)
      SIZE = DIM(4) - 1
      K2 = DIM(5)
      NS = N * SEP
      M3 = M * 3
      FM3 = FLOAT(M3)
      MM3 = SEP * M3
      MOVER2 = M/2 + 1
      FJM1 = -1.0
      DO 600 J=1, MOVER2
      FOLD = J.GT.1 .AND. 2*J.LT.M+2
      K0 = (J-1)*SEP + 1
      FJM1 = FJM1 + 1.0
      ANGLE = TWOPI * FJM1 / FM3
      ZERO = ANGLE.EQ.0.0
      IF (ZERO) GOTO 200
      C1 = COS(ANGLE)
      S1 = SIN(ANGLE)
      C2 = C1*C1 - S1*S1
      S2 = S1*C1 + C1*S1
      GOTO 200
  100 FOLD = .FALSE.
      K0 = (M+1-J)*SEP + 1
      T = C1*A + S1*B
      S1 = C1*B - S1*A
      C1 = T
      T = C2*A - S2*B
      S2 = -C2*B - S2*A
      C2 = T
  200 DO 500 KK=K0,NS,MM3
      DO 440 L=KK,NT,L1
      K1 = L + SIZE
      DO 420 K=L,K1,K2
      R0 = X0(K)
      I0 = Y0(K)
      RS = X1(K) + X2(K)
      IS = Y1(K) + Y2(K)
      X0(K) = R0 + RS
      Y0(K) = I0 + IS
      RA = R0 + RS*A
      IA = I0 + IS*A
      RB = (X1(K)-X2(K)) * B
      IB = (Y1(K)-Y2(K)) * B
      IF (ZERO) GOTO 300
      R1 = RA + IB
      I1 = IA - RB
      R2 = RA - IB
      I2 = IA + RB
      X1(K) = R1*C1 + I1*S1
      Y1(K) = I1*C1 - R1*S1
      X2(K) = R2*C2 + I2*S2
      Y2(K) = I2*C2 - R2*S2
      GOTO 420
  300 X1(K) = RA + IB
      Y1(K) = IA - RB
      X2(K) = RA - IB
      Y2(K) = IA + RB
  420 CONTINUE
  440 CONTINUE
  500 CONTINUE
      IF (FOLD) GOTO 100
  600 CONTINUE
      RETURN
      END
      SUBROUTINE R4CFTK (N, M, X0, Y0, X1, Y1, X2, Y2, X3, Y3, DIM,
     +   NDIR0, NDIR1, NDIR2, NDIR3)
      REAL X0(NDIR0),Y0(NDIR0),X1(NDIR1),Y1(NDIR1),
     +     X2(NDIR2),Y2(NDIR2),X3(NDIR3),Y3(NDIR3)
      INTEGER DIM(5), SIZE, SEP
      LOGICAL FOLD,ZERO
      REAL I1,I2,I3,IS0,IS1,IU0,IU1
      DATA TWOPI / 6.2831853 /
      DATA C1, S1, C2, S2, C3, S3 / 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 /
      NT = DIM(1)
      SEP = DIM(2)
      L1 = DIM(3)
      SIZE = DIM(4) - 1
      K2 = DIM(5)
      NS = N * SEP
      M4 = M * 4
      FM4 = FLOAT(M4)
      MM4 = SEP * M4
      MOVER2 = M/2 + 1
      FJM1 = -1.0
      DO 600 J=1, MOVER2
      FOLD = J.GT.1 .AND. 2*J.LT.M+2
      K0 = (J-1)*SEP + 1
      FJM1 = FJM1 + 1.0
      ANGLE = TWOPI * FJM1 / FM4
      ZERO = ANGLE.EQ.0.0
      IF (ZERO) GOTO 200
      C1 = COS(ANGLE)
      S1 = SIN(ANGLE)
      C2 = C1*C1 - S1*S1
      S2 = S1*C1 + C1*S1
      C3 = C2*C1 - S2*S1
      S3 = S2*C1 + C2*S1
      GOTO 200
  100 FOLD = .FALSE.
      K0 = (M+1-J)*SEP + 1
      T = C1
      C1 = S1
      S1 = T
      C2 = -C2
      T = C3
      C3 = -S3
      S3 = -T
  200 DO 500 KK=K0,NS,MM4
      DO 440 L=KK,NT,L1
      K1 = L + SIZE
      DO 420 K=L,K1,K2
      RS0 = X0(K) + X2(K)
      IS0 = Y0(K) + Y2(K)
      RU0 = X0(K) - X2(K)
      IU0 = Y0(K) - Y2(K)
      RS1 = X1(K) + X3(K)
      IS1 = Y1(K) + Y3(K)
      RU1 = X1(K) - X3(K)
      IU1 = Y1(K) - Y3(K)
      X0(K) = RS0 + RS1
      Y0(K) = IS0 + IS1
      IF (ZERO) GOTO 300
      R1 = RU0 + IU1
      I1 = IU0 - RU1
      R2 = RS0 - RS1
      I2 = IS0 - IS1
      R3 = RU0 - IU1
      I3 = IU0 + RU1
      X2(K) = R1*C1 + I1*S1
      Y2(K) = I1*C1 - R1*S1
      X1(K) = R2*C2 + I2*S2
      Y1(K) = I2*C2 - R2*S2
      X3(K) = R3*C3 + I3*S3
      Y3(K) = I3*C3 - R3*S3
      GOTO 420
  300 X2(K) = RU0 + IU1
      Y2(K) = IU0 - RU1
      X1(K) = RS0 - RS1
      Y1(K) = IS0 - IS1
      X3(K) = RU0 - IU1
      Y3(K) = IU0 + RU1
  420 CONTINUE
  440 CONTINUE
  500 CONTINUE
      IF (FOLD) GOTO 100
  600 CONTINUE
      RETURN
      END
      SUBROUTINE R5CFTK (N, M, X0, Y0, X1, Y1, X2, Y2, X3, Y3, X4, Y4,
     *                   DIM, NDIR0, NDIR1, NDIR2, NDIR3, NDIR4)
      REAL X0(NDIR0),Y0(NDIR0),X1(NDIR1),Y1(NDIR1),X2(NDIR2),
     +     Y2(NDIR2),X3(NDIR3),Y3(NDIR3),X4(NDIR4),Y4(NDIR4),
     +     I0,I1,I2,I3,I4,IA1,IA2,IB1,IB2,IS1,IS2,IU1,IU2
      INTEGER DIM(5), SIZE, SEP
      LOGICAL FOLD,ZERO
      DATA TWOPI / 6.2831853 /
      DATA A1/0.30902/   ,B1/0.95106/   ,A2/-0.80902/   ,B2/0.58779/
      DATA C1, S1, C2, S2, C3, S3, C4, S4 / 0.,0.,0.,0.,0.,0.,0.,0. /
      NT = DIM(1)
      SEP = DIM(2)
      L1 = DIM(3)
      SIZE = DIM(4) - 1
      K2 = DIM(5)
      NS = N * SEP
      M5 = M * 5
      FM5 = FLOAT(M5)
      MM5 = SEP*M5
      MOVER2 = M/2 + 1
      FJM1 = -1.0
      DO 600 J=1, MOVER2
      FOLD = J.GT.1 .AND. 2*J.LT.M+2
      K0 = (J-1)*SEP + 1
      FJM1 = FJM1 + 1.0
      ANGLE = TWOPI * FJM1 / FM5
      ZERO = ANGLE.EQ.0.0
      IF (ZERO) GOTO 200
      C1 = COS(ANGLE)
      S1 = SIN(ANGLE)
      C2 = C1*C1 - S1*S1
      S2 = S1*C1 + C1*S1
      C3 = C2*C1 - S2*S1
      S3 = S2*C1 + C2*S1
      C4 = C2*C2 - S2*S2
      S4 = S2*C2 + C2*S2
      GOTO 200
  100 FOLD = .FALSE.
      K0 = (M+1-J)*SEP + 1
      T = C1*A1 + S1*B1
      S1 = C1*B1 - S1*A1
      C1 = T
      T = C2*A2 + S2*B2
      S2 = C2*B2 - S2*A2
      C2 = T
      T = C3*A2 - S3*B2
      S3 = -C3*B2 - S3*A2
      C3 = T
      T = C4*A1 - S4*B1
      S4 = -C4*B1 - S4*A1
      C4 = T
  200 DO 500 KK=K0,NS,MM5
      DO 440 L=KK,NT,L1
      K1 = L + SIZE
      DO 420 K=L,K1,K2
      R0 = X0(K)
      I0 = Y0(K)
      RS1 = X1(K) + X4(K)
      IS1 = Y1(K) + Y4(K)
      RU1 = X1(K) - X4(K)
      IU1 = Y1(K) - Y4(K)
      RS2 = X2(K) + X3(K)
      IS2 = Y2(K) + Y3(K)
      RU2 = X2(K) - X3(K)
      IU2 = Y2(K) - Y3(K)
      X0(K) = R0 + RS1+RS2
      Y0(K) = I0 + IS1+IS2
      RA1 = R0 + RS1*A1+RS2*A2
      IA1 = I0 + IS1*A1+IS2*A2
      RA2 = R0 + RS1*A2+RS2*A1
      IA2 = I0 + IS1*A2+IS2*A1
      RB1 = RU1*B1 + RU2*B2
      IB1 = IU1*B1 + IU2*B2
      RB2 = RU1*B2 - RU2*B1
      IB2 = IU1*B2 - IU2*B1
      IF (ZERO) GOTO 300
      R1 = RA1 + IB1
      I1 = IA1 - RB1
      R2 = RA2 + IB2
      I2 = IA2 - RB2
      R3 = RA2 - IB2
      I3 = IA2 + RB2
      R4 = RA1 - IB1
      I4 = IA1 + RB1
      X1(K) = R1*C1 + I1*S1
      Y1(K) = I1*C1 - R1*S1
      X2(K) = R2*C2 + I2*S2
      Y2(K) = I2*C2 - R2*S2
      X3(K) = R3*C3 + I3*S3
      Y3(K) = I3*C3 - R3*S3
      X4(K) = R4*C4 + I4*S4
      Y4(K) = I4*C4 - R4*S4
      GOTO 420
  300 X1(K) = RA1 + IB1
      Y1(K) = IA1 - RB1
      X2(K) = RA2 + IB2
      Y2(K) = IA2 - RB2
      X3(K) = RA2 - IB2
      Y3(K) = IA2 + RB2
      X4(K) = RA1 - IB1
      Y4(K) = IA1 + RB1
  420 CONTINUE
  440 CONTINUE
  500 CONTINUE
      IF (FOLD) GOTO 100
  600 CONTINUE
      RETURN
      END
      SUBROUTINE HERMFT (X, Y, NSIZE, N, DIM)
      REAL X(NSIZE), Y(NSIZE)
      INTEGER DIM(5), D2, D3, D4, D5
      DATA TWOPI / 6.2831853 /
      TWON = FLOAT(2*N)
      NT = DIM(1)
      D2 = DIM(2)
      D3 = DIM(3)
      D4 = DIM(4) - 1
      D5 = DIM(5)
      DO 100 I0=1,NT,D3
      I1 = I0 + D4
      DO 100 I=I0,I1,D5
      A = X(I)
      B = Y(I)
      X(I) = A + B
      Y(I) = A - B
  100 CONTINUE
      NOVER2 = N/2 + 1
      IF (NOVER2 .LT. 2) GOTO 500
      DO 400 I0 = 2, NOVER2
      ANGLE = TWOPI * FLOAT(I0-1) / TWON
      CO = COS(ANGLE)
      SI = SIN(ANGLE)
      K = (N + 2 - 2*I0)*D2
      K1 = (I0 - 1)*D2 + 1
      DO 300 I1=K1,NT,D3
      I2 = I1 + D4
      DO 200 I=I1,I2,D5
      J = I + K
      A = X(I) + X(J)
      B = X(I) - X(J)
      C = Y(I) + Y(J)
      D = Y(I) - Y(J)
      E = B*CO + C*SI
      F = B*SI - C*CO
      X(I) = A + F
      X(J) = A - F
      Y(I) = E + D
      Y(J) = E - D
  200 CONTINUE
  300 CONTINUE
  400 CONTINUE
      CALL CMPLFT (X, Y, NSIZE, N, DIM)
  500 RETURN
      END
      SUBROUTINE RDSECT (MAX, NNXP2, NNZ, NXZ3, LIN)
      PARAMETER (KUSER2=30000)
      COMMON / / NR3D
      INTEGER*2 NR3D(KUSER2)
      INTEGER*2 BLACOM(84000)
      EQUIVALENCE (BLACOM(1), NR3D(1))
      IF (MAX .GE. NXZ3) MAX = 0
      MX = MAX
      MAX = MAX - 2
      DO 1320 IZ=1,NNZ
      MIN = MAX + 3
      MAX = MAX + NNXP2
      READ (LIN) IYSEC, IZLIN, IXTOT, (NR3D(IX),IX=MIN,MAX)
      NR3D(MAX+1) = NR3D(MIN)
      NR3D(MAX+2) = NR3D(MIN+1)
 1320 CONTINUE
      MIN = MAX + 3
      MAX = MAX + NNXP2 + NNXP2 + 2
      DO 1340 IX=MIN,MAX
      MX = MX + 1
      NR3D(IX) = NR3D(MX)
 1340 CONTINUE
      RETURN
      END
      SUBROUTINE SORT (X, MAXAT, NAT, N)
      DIMENSION X(4,MAXAT), T(4)
      INT=2
 1000 INT=INT+INT
      IF(INT.LT.NAT)GO TO 1000
      INT = MIN0 (NAT, (3*INT)/4-1)
 1020 INT=INT/2
      IFIN=NAT-INT
      DO 1200 II=1,IFIN
      I=II
      J=I+INT
      IF (X(N,I) .GE. X(N,J)) GOTO 1200
      DO 1060 K=1,4
      T(K) = X(K,J)
 1060 CONTINUE
 1080 DO 1100 K=1,4
      X(K,J) = X(K,I)
 1100 CONTINUE
      J=I
      I=I-INT
      IF (I) 1140, 1140, 1120
 1120 IF (X(N,I) .LT. T(N)) GOTO 1080
 1140 DO 1160 K=1,4
      X(K,J) = T(K)
 1160 CONTINUE
 1200 CONTINUE
      IF(INT.NE.1)GO TO 1020
      RETURN
      END
      FUNCTION QUAD2 (X1, X2)
      DIMENSION X1(3), X2(3)
      COMMON /CRYSA/ CELL(6),  CELLSD(6),   RCELL(6),   VOLUM,
     +               WAVE,     CELALL(10),  AMOLW,      ZET,
     +               NELEC,    F000,        ABSMU,      ICENT,
     +               ILATT,    ISYST,       ILAUE,      IMULT,
     +               IUNIQ,    IPOLA,       NTYPE,      NSYMM,
     +         IRSYMM(3,3,24), TSYMM(3,24), NLATT,      TLATT(3,4),
     +         FRAC2C(3,3),    CART2F(3,3), RRMAT(3,3), SSMAT(3,3)
      CALL VMATV1 (X1, RRMAT, X2, DIST2)
      QUAD2 = DIST2
      RETURN
      END
      SUBROUTINE OPER1 (J, XN, XYZ)
      DIMENSION XN(3), XYZ(3)
      DIMENSION FS(3,3,24)
      COMMON /CRYSA/ CELL(6),  CELLSD(6),   RCELL(6),   VOLUM,
     *               WAVE,     CELALL(10),  AMOLW,      ZET,
     *               NELEC,    F000,        ABSMU,      ICENT,
     *               ILATT,    ISYST,       ILAUE,      IMULT,
     *               IUNIQ,    IPOLA,       NTYPE,      NSYMM,
     *         IRSYMM(3,3,24), TSYMM(3,24), NLATT,      TLATT(3,4),
     *         FRAC2C(3,3),    CART2F(3,3), RRMAT(3,3), SSMAT(3,3)
      LOGICAL CONT
      DATA CONT / .FALSE. /
      IF (CONT) GOTO 111
      CONT = .TRUE.
      CALL KERI2F (IRSYMM, FS, 9*NSYMM)
  111 ISYM = MOD(J,NSYMM)
      IF (ISYM .EQ. 0) ISYM = NSYMM
      IP = (J-1) / NSYMM
      ILAT = MOD(IP,NLATT) + 1
      IF (IP .LT. NLATT) THEN
         DO 120 I = 1,3
  120    XN(I) = XYZ(1) * FS(I,1,ISYM)
     *         + XYZ(2) * FS(I,2,ISYM)
     *         + XYZ(3) * FS(I,3,ISYM) + TSYMM(I,ISYM) + TLATT(I,ILAT)
      ELSE
         DO 130 I = 1,3
  130    XN(I) =-XYZ(1) * FS(I,1,ISYM)
     *         - XYZ(2) * FS(I,2,ISYM)
     *         - XYZ(3) * FS(I,3,ISYM) - TSYMM(I,ISYM) + TLATT(I,ILAT)
         ENDIF
      RETURN
      END
