PROGRAM TRACOR ************************** File: TRACOR.NEW: 21 July '97 TRACOR.FOR * Last update: 11 Nov. 1998 * Source: ROG 1986 , SGG 1987 * After 'CSUBFILE FOUR2' follows the second part of TRACOR = Fourier pgm ************************************************************************ C wensen: - FOUR groter in Tracor - FOM voor BOTS ! *TRACOR LOG of recent modifications (last on top 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 (ICRYS, IFILE(3)) EQUIVALENCE (IBINFF, IFILE(16)) 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) DIMENSION SEDEL2(50,2) CALL KEPROG('TRACOR') WRITE (LIS2, FMT = '('' Last TRACOR update: 11 Nov. 1999'')') CALL RDCRYS (ICRYS) IF (IPOLA .EQ. 7) CALL KERROR * ('TRACOR not applicable to space group P1' , 0, 'TRACIN') CALL FCALIN CALL TRACIN CALL COTRA1 IF (IPOLA .LE. 0) THEN CALL D3FOUR CALL FILCLO (IBINFF, 'DELETE') ELSE CALL D12FOU ENDIF CALL COTRA3 CALL R2TRAC (SEDEL2) CALL SHIFTR (SEDEL2) CALL KEPROX WRITE (LIS2, FMT='(/'' Test: DDOKA exit MAIN SUBPROGRAM ''/)') STOP 99 END SUBROUTINE FCALIN 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)), (IATOMS, IFILE(2)) EQUIVALENCE (IDDS, IFILE(1)), (ICRYS,IFILE(3)) EQUIVALENCE (ICON, IFILE(4)), (LIS1, IFILE(7)), (LIS2, IFILE(8)) EQUIVALENCE (IBINFO, IFILE(11)), (IBINFC, IFILE(12)) 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), PSQ, 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), FITFC2(51) EQUIVALENCE (HCODE, FITFO(1)), (EPSIL2, FITFC2(1)) COMMON /COTRA/ MULS(48,48), INDX(3,3), IDIMF, TT(3,48), * EFST(48), FHST(48), PHST(48), CELLN(6), MAXHKL(3), * P111, P77, EMIN, EO2AV, SCX, BOVX, DAMPX, SMM, IOMAPS, SMAX, * IRR(3,3,48), IORG, BBB, D2R, R2D, * ATEM1, ATEM2, ATEM3, ATEM4, NTEM5, ATEM6 PARAMETER (MAXAT = 993) PARAMETER (MAXBUF = 198) COMMON / / SICO(12500), FF(500,10), EXPBP(500), EXPBR(500), * SUMF2(500), SUMF2P(500), SFAC(13,10), * ATXYZ(10,MAXAT), IZAT(MAXAT), ITAT(MAXAT), NAT, * BUFFO(MAXBUF), BUFFC(MAXBUF), BUFBUF(MAXBUF) DIMENSION BLACOM(42000) EQUIVALENCE (BLACOM(1), SICO(1)) COMMON /ATNAMA/ ATNAME(MAXAT) CHARACTER * 6 ATNAME PARAMETER (LCMAX = 15) CHARACTER * 6 LCONDA(LCMAX) DATA LCONDA / 'TRACOR', 'EMIN' , 'SCSG' , 'BHSG' , 'DAMP' , * 'SCALE' , 'BBB' , 'STLMAX', 'SMM' , 'MAXHKL', * 'PRIMAP', 'PSQMAX', 'XXXXXX', 'WILSON', 'PRINT' / WRITE (LIS2, 101) 101 FORMAT (' Structure factor calculation for TRACOR:' * / ' EXPAND data to P1 symmetry (or centered equivalent)') CALL RDCRYS (ICRYS) DO 110 I=1,NTYPE 110 BUFFO(I) = CELALL(I) / ZET I = NTYPE J = NINT(ZET) WRITE (LIS2, 114) J, (CELATY(K), BUFFO(K), K=1,I) 114 FORMAT (' Z:', I3 / ' FORMUL:', 6(2X,A2,F6.1) / * ( 8X, 6(2X,A2,F6.1))) CALL LOGRD (IDDL, 'MERBSC', KLOG) IF (KLOG .LT. 0 .OR. LIT(2) .NE. 'SCALE') THEN CHOUT=' MERBIN SCALE not found on DDLOG file: what happened?' CALL SHOUT CALL KERROR ('You cleared too much, probably!', 114, 'FCALIN') ENDIF SCALE = FNUM(2) BOV = FNUM(3) BP = BOV BR = BOV WRITE (LIS2, 116) SCALE, BOV 116 FORMAT (' Data from DDLOG: Scale =', F9.5, ' Bov:', F6.3) CALL FILCLO (IDDL, 'KEEP') CALL BINIFF (1, IBINFO, 'BINFO', FITFO, NITFO, BUFFO, KENDFO) STLMAX = BUFFO(6) CALL KERNAB (BUFFO(7), HKLMAX, 3) CALL KERF2I (HKLMAX, MAXHKL, 3) EMIN = 0.0 SCX = SCALE BOVX = BOV DAMPX = BP SMM = 0.00001 IOMAPS = 0 P77 = 1. 120 CALL RDCOND (ICON, LCONDA, LCMAX, KEND) GOTO (120, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 120, 120,120), KEND IF (KEND.LE.0) GOTO 140 2 EMIN = FNUM(1) WRITE (LIS1, FMT = '('' Input value for Emin: '', F10.2)') EMIN GOTO 120 3 SCSG = FNUM(1) IF(SCSG .LT. 0.001) SCSG = 1.0000 SCX = SCALE * SCSG WRITE (LIS1, FMT = '('' SCALE '', F10.5,'' = SCALE '', F10.5, * '' multiplied by extra scale: '', F10.5)') SCX, SCALE, SCSG GOTO 120 4 BHSG=FNUM(1) BOVX = BOV + BHSG WRITE (LIS1, FMT = '('' Input: additional Bp'', F10.5 )') BHSG WRITE (LIS1, FMT = '('' Now to be used: Bov'', F10.5 )') BOVX GOTO 120 5 DAMP = 0.5 * FNUM(1) DAMPX = DAMPX - DAMP WRITE (LIS1, FMT = '('' Supplied: damping parameter '', F10.5)') * FNUM(1) GOTO 120 6 IF (FNUM(1) .LT. 0.0001) GOTO 120 SCALE = FNUM(1) WRITE (LIS1, 125) SCALE 125 FORMAT (' Scale from CONDA file: Scale =', F9.5) GOTO 120 7 IF (NFNUM.LT.1) CALL KERNER (7, 'FCALIN') IF (FNUM(1).GT.0.0001) THEN BOV = FNUM(1) BP = BOV BR = BOV ENDIF IF (FNUM(2).GT.0.0001) BP = FNUM(2) IF (FNUM(3).GT.0.0001) BR = FNUM(3) WRITE (LIS1, 126) FNUM(1), FNUM(2), FNUM(3) 126 FORMAT (' Temp. factors from CONDA: Bov=', * F6.3, ' Bp =', F6.3, ' Br =', F6.3) IF (FNUM(2).LE.0.0001 .OR. FNUM(3).LE.0.0001) * WRITE (LIS1, 127) BOV, BP, BR 127 FORMAT (' Temp. factors to be used: Bov=', * F6.3, ' Bp =', F6.3, ' Br =', F6.3) GOTO 120 8 IF (NFNUM.NE.1) CALL KERNER (8, 'FCALIN') IF (FNUM(1) .LT. 0.0001 .OR. FNUM(1) .GE. STLMAX) GOTO 120 STLMAX = FNUM(1) WRITE (LIS1, 128) STLMAX 128 FORMAT (' Skip reflections if sin(th)/lambda >', F8.4) GOTO 120 9 SMM =FNUM(1) IF(SMM.LT.0.0001) SMM=0.000001 WRITE (LIS1, 129) SMM 129 FORMAT (' Input factor (for testruns) SMM = ', F8.4) GOTO 120 10 J = 0 DO 137 I = 1, 3 IF (FNUM(I).LT.0.0001) GOTO 137 K = NINT(FNUM(I)) IF (K .GE. MAXHKL(I)) GOTO 137 J = 1 MAXHKL(I) = K 137 CONTINUE IF (J .EQ. 0) GOTO 120 WRITE (LIS1, 138) MAXHKL 138 FORMAT (' Skip reflections if indices exceed MAXHKL =', 3I3) CALL KERI2F (MAXHKL, HKLMAX, 3) GOTO 120 11 IOMAPS = 1 GOTO 120 12 IF (FNUM(1) .GE. 0.01) P77 = FNUM(1) GOTO 120 140 CALL FILCLO (ICON, 'KEEP') SMAX = STLMAX NSET = 0 CALL FILINQ (IATOMS, 'ATOMS', 'FORMATTED', 'INPUT', KINQ) IF (KINQ.EQ.-1) CALL KERROR (' No ATOMS file found', * 140, 'FCALIN') CALL KERINA (IATOMS, LIT, 1, LEND) IF (LIT(1) .NE. 'ATOMS') CALL KERROR * (' Incorrect header on ATOMS file', 143, 'FCALIN') BACKSPACE IATOMS CALL ATOMIN (IATOMS, ATXYZ, ATNAME, IZAT, MAXAT, NAT, KEYT) NSET = NSET + 1 IF (KEYT .GT. 1) WRITE (LIS1, FMT='( * '' Temperature factors present in ATOMS file ignored'')') KEYT = 1 NATQ = NAT NATH = 0 N = 1 143 CONTINUE IF (ATNAME(N)(1:1).EQ.'H' .AND. IZAT(N).EQ.1) NATH = NATH + 1 IF (ATNAME(N)(1:1).EQ.'H' .OR. IZAT(N).EQ.1 .OR. * ATNAME(N)(1:1) .EQ. 'Q') THEN IF (N .EQ. NAT) GOTO 148 DO 146 N1 = N, NAT - 1 CALL KERNAB (ATXYZ(1,N1+1), ATXYZ(1,N1), 4) ATNAME(N1) = ATNAME(N1+1) 146 IZAT(N1) = IZAT(N1+1) 148 NAT = NAT - 1 N = N - 1 ENDIF CALL KERNZA (0.0, ATXYZ(5,N), 6) N = N + 1 IF (N .LE. NAT) GOTO 143 IF (NATH .NE. 0) WRITE (LIS1, FMT= * '('' Number of H atoms rejected:'', I3)') NATH N = NATQ - NAT - NATH IF (N .GT. 0) WRITE (LIS1, FMT= * '('' Nr of Q-atoms (= peaks) rejected:'', I3)') N IF (NAT .LE. 0) CALL KERROR ('.... No atoms left!', 0, 'FCALIN') CALL ATOMPR (LIS1, 7, ATXYZ, ATNAME, IZAT, NAT) CALL KERNAB (BUFFO(5), BUFFC(5), MAXBUF - 4) BUFFC(16) = SCALE BUFFC(17) = BOV BUFFC(18) = SCALE BUFFC(19) = BP BUFFC(20) = BR BUFFC(21) = STLMAX CALL KERNAB (HKLMAX, BUFFC(22), 3) CALL FCALCI (KEYT, ATXYZ, IZAT, ITAT, NAT) BUFFC(25) = NAT BUFFC(26) = P1SQ BUFFC(27) = PSQ CALL FCALCX RETURN END SUBROUTINE FCALCX 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 (IPR1, IFILE(6)), (LIS1, IFILE(7)), (LIS2, IFILE(8)) EQUIVALENCE (IBINFO, IFILE(11)), (IBINFC, IFILE(12)) 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), PSQ, 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), FITFC2(51) EQUIVALENCE (HCODE, FITFO(1)), (EPSIL2, FITFC2(1)) PARAMETER (MAXAT = 993) PARAMETER (MAXBUF = 198) COMMON / / SICO(12500), FF(500,10), EXPBP(500), EXPBR(500), * SUMF2(500), SUMF2P(500), SFAC(13,10), * ATXYZ(10,MAXAT), IZAT(MAXAT), ITAT(MAXAT), NAT, * BUFFO(MAXBUF), BUFFC(MAXBUF), BUFBUF(MAXBUF) DIMENSION BLACOM(42000) EQUIVALENCE (BLACOM(1), SICO(1)) COMMON /ATNAMA/ ATNAME(MAXAT) CHARACTER * 6 ATNAME WRITE (LIS1, FMT='('' Calculate structure factors'')') FICENT = FLOAT (ICENT) WRITE (LIS1, 183) P1SQ 183 FORMAT (' Scattering fraction of known part:'/ * ' excluding symmetry related molecules: P1**2 =', F6.3) NITFC = 3 + 2 * NSYMM CALL BINOFF (27, IBINFC,'BINFC2', FITFC2, NITFC, BUFFC, KENDFC) NREFL = 0 IREFL = 0 SUMNR2 = 0. R2SUMA = 0. R2SUMB = 0. ALATT2 = NLATT * NLATT CALL BINIFF (1, IBINFO, 'BINFO', FITFO, NITFO, BUFFO, KENDFO) NITFO = 3 200 CALL BINIFF (0, IBINFO, 'BINFO', FITFO, NITFO, BUFFO, KENDFO) IF (KENDFO.LT.0) GOTO 220 IREFL = IREFL + 1 CALL HKLC1U (HCODE, HKLX) CALL HKLSTL (HKLX, STL, STL2) IF (STL.GT.STLMAX .OR. HKLX(1,1).GT.HKLMAX(1) .OR. * HKLX(2,1).GT.HKLMAX(2) .OR. HKLX(3,1).GT.HKLMAX(3)) THEN FITFC2(1) = -999. GOTO 215 ENDIF CALL FCALX1 (ATXYZ, ITAT, NAT) SUMNR2 = SUMNR2 + ASYMM / EPSIL2 SF2 = SUMF2 (ISS) SF2P= SUMF2P(ISS) NREFL = NREFL + 1 G2 = EXPBR(ISS) * EXPBR(ISS) * SUMF2(ISS) * ALATT2 EO2 = ( SCALE * FOBS ) **2 / G2 EO4 = EO2 * EO2 DO 207 I = 1, NSYMM IF (FPEXP(1,I) .LT. 0.) GOTO 207 R2SUMA = R2SUMA + (EO2 - FPEXP(1,I)*FPEXP(1,I)/G2)**2 R2SUMB = R2SUMB + EO4 207 CONTINUE 215 CONTINUE CALL BINOFF (0, IBINFC, 'BINFC2',FITFC2, NITFC, BUFFC, KENDFC) GOTO 200 220 CONTINUE CALL BINOFF (-1, IBINFC, 'BINFC2',FITFC2, NITFC, BUFFC, KENDFC) IF (IREFL.NE.NREFL) WRITE (LIS1, 242) NREFL 242 FORMAT (' Note:', I6,' input reflections accepted ') WRITE (LIS1, 254) NINT(SUMNR2) 254 FORMAT (' Number of reflections after expansion:', I12) WRITE (IPR1, 258) 258 FORMAT (' Structure factor calculation finished') R2 = R2SUMA / R2SUMB R2E = 1. - P1SQ WRITE (LIS1, FMT='(/'' R2 for the model in SpGr P1:'', F5.2, * '' ( expected:'', F5.2, '' )''/ 34(''-'')/ )') R2, R2E RETURN END SUBROUTINE FCALX1 (ATXYZ, ITAT, NAT) DIMENSION ATXYZ(10,NAT), ITAT(NAT) 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 /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 /FCALCA/ BP, BR, SCALE, HKLMAX(3), STLMAX, * IZTYPE(10), CELPAR(10), PSQ, 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) COMMON / / SICO(12500), FF(500,10), EXPBP(500), EXPBR(500), * SUMF2(500), SUMF2P(500), SFAC(13,10) DIMENSION BLACOM(42000) EQUIVALENCE (BLACOM(1), SICO(1)) DIMENSION FFF(10) S = STL * 400. + 1. IS = IFIX(S) STLDEL = S - FLOAT(IS) ISS = NINT(S) DO 110 J=1,NTYPE IF (CELPAR(J).LE.0.0) GOTO 110 FFF(J) = FF(IS,J) + (FF(IS+1,J)-FF(IS,J)) * STLDEL 110 CONTINUE CALL HKLEX1 (HKLX, HKLX) CALL HKLEX2 (HKLX, IDHKL, IEPS, IEPS2) EPSIL = IEPS EPSIL2 = IEPS2 DO 600 J=1,NSYMM IF (IDHKL(J).EQ.0) GOTO 200 K = IABS(IDHKL(J)) FPEXP(1,J) = -K FPEXP(2,J) = FPEXP(2,K) IF (FPEXP(2,J).LT.0.0001) FPEXP(2,J)=0.0001 IF (IDHKL(J).LT.0) FPEXP(2,J)=-FPEXP(2,J) GOTO 600 200 FAP = 0.0 FBP = 0.0 DO 400 I=1,NAT TRIG = HKLX(1,J)*ATXYZ(1,I) + HKLX(2,J)*ATXYZ(2,I) + * HKLX(3,J)*ATXYZ(3,I) IF (TRIG.LT.0.0) TRIG = TRIG - 0.00010 ITRIG = MOD ( IFIX(TRIG * 10000. + 0.5), 10000) IF (ITRIG.LE.0) ITRIG = ITRIG + 10000 IJ = ITAT(I) FAP = FAP + SICO(ITRIG + 2500) * FFF(IJ) FBP = FBP + SICO(ITRIG) * FFF(IJ) 400 CONTINUE FP = ALATT * SQRT (FAP*FAP + FBP*FBP) * EXPBP(ISS) PHIP = 0.0 IF (FP.GT.0.001) PHIP = ATAN2(FBP,FAP) / 0.0174532925 IF (PHIP.LT.0.0) PHIP = PHIP + 360. FPEXP(1,J) = FP FPEXP(2,J) = PHIP 600 CONTINUE RETURN END SUBROUTINE TRACIN 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 (ICON, IFILE(4)) EQUIVALENCE (LIS2, IFILE(8)) EQUIVALENCE (IBINFO, IFILE(11)) EQUIVALENCE (IBINFC, IFILE(12)) LOGICAL SWPRI EQUIVALENCE (SWPRI, SWITCH(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) COMMON /CRYSB/ SPGR, WAVEAT, CELATY(10) CHARACTER SPGR *16, WAVEAT *2, CELATY *2 PARAMETER (KFA = 20671) PARAMETER (MAXBUF = 198) COMMON / / FA(KFA), FB(KFA), + BUFFC2(MAXBUF), BUFFO(MAXBUF), BUFFFT(MAXBUF), + FITFC2(51), FITFO(3), FITFFT(5), + NITFC2, KENDFC, NITFO, KENDFO DIMENSION BLACOM(42000) EQUIVALENCE (BLACOM(1), FA(1)) PARAMETER (MXP=50, MXP1=MXP + 1) COMMON /PIEK/ RHO(MXP1), SH(3,MXP1), IBOTS(MXP1), BOTS(MXP1), NPK COMMON /COTRA/ MULS(48,48), INDX(3,3), IDIMF, TT(3,48), * EFST(48), FHST(48), PHST(48), CELLN(6), MAXHKL(3), * P1SQ, PSQ, EMIN, EO2AV, SCX, BOVX, DAMPX, SMM, IOMAPS, SMAX, * IRR(3,3,48), IORG, BBB, D2R, R2D, * AMULT, ASYMM, ALATT, ASYMCL, NSYMC, ASYMC DIMENSION IRP(3,3) DIMENSION EO2A(4,2) DATA EO2A /1.000,1.250,2.000,3.250,1.000,1.5705,2.5251,3.9080/ D2R = ATAN(1.0) / 45.0 R2D = 1.0 / D2R CALL KERNZA(0.0, RHO, MXP1) CALL KERNZA(0.0, SH, 3*MXP1) NPK = 0 NR = 0 AMULT = FLOAT (IMULT) ASYMM = FLOAT (NSYMM) ALATT = FLOAT (NLATT) ASYMCL= FLOAT (ICENT*NLATT) NSYMC = NSYMM * ICENT ASYMC = FLOAT (NSYMC) CALL BINIFF (1,IBINFC,'BINFC2',FITFC2,NITFC2,BUFFC2,KENDFC) SC = BUFFC2(18) BH = BUFFC2(19) BOV = BH PSQMAX = PSQ P1SQ = BUFFC2(26) PSQ = BUFFC2(27) IF (PSQMAX .GE. 0.01) PSQ = AMIN1(PSQ, PSQMAX) SMAX = BUFFC2(21) CALL KERNAB (BUFFC2, BUFFFT, 27) BUFFFT(28) = 5. CALL BINIFF (1,IBINFO,'BINFO',FITFO,NITFO,BUFFO,KENDFO) IEMIN = IFIX (EMIN * 2. + 1.5) IF (IEMIN.GT.4) IEMIN=4 EMINX = FLOAT(IEMIN)/2. -0.5 IF (ABS(EMINX - EMIN) .GT. 0.01) THEN EMIN = EMINX WRITE (CHOUT, FMT = '('' Reset value for Emin: '', F10.2)') EMIN CALL SHOUT ENDIF EO2AV = EO2A(IEMIN,ICENT) BBB = EO2AV * (1. - PSQ * P1SQ) - (PSQ - PSQ * P1SQ) IF (SWPRI) WRITE (LIS2, 334) P1SQ, PSQ, EMIN, EO2AV, BBB 334 FORMAT (4X, 'P1SQ PSQ EMIN ', + ' EO2AV BBB' / 5F10.5 ) IF (SWPRI) WRITE (LIS2, 336) SCX, BOVX, DAMPX, SMM 336 FORMAT (3X,'Scale',6X,'BH',6X,' Damp-F',6X, * 'SMM sharpening params' / 6F10.5 ) DO 402 I=1,NSYMM L=I+NSYMM DO 402 J=1,3 TT(J,I)=TSYMM(J,I) IF(ICENT.EQ.2) THEN TTJL=AMOD(1.0-TT(J,I),1.0) TT(J,L)=TTJL ENDIF DO 402 K=1,3 IF(ICENT.EQ.2) IRR(K,J,L) = -IRSYMM(K,J,I) 402 IRR(K,J,I) = IRSYMM(K,J,I) WRITE (LIS2, FMT='('' Symmetry operations for sp.gr. '', A16)') * SPGR DO 419 I=1,NSYMC WRITE (LIS2, 407) I, ((IRR(L,M,I),M=1,3), TT(L,I), L=1,3) 407 FORMAT (I4, ' :', 3(' / ', 3I3, F8.3)) DO 419 J=1,NSYMC CALL MATAXI (IRR(1,1,I), IRR(1,1,J), IRP) DO 418 M=1,NSYMC IFIT=0 DO 415 K=1,3 DO 415 L=1,3 IF (IRR(K,L,M) .EQ. IRP(K,L)) IFIT=IFIT+1 415 CONTINUE IF (IFIT .NE. 9) GOTO 418 MULS(I,J)=M GOTO 419 418 CONTINUE CALL KERROR ('Symmetry error' , 0, 'SYTRA') 419 CONTINUE IF (SWPRI) THEN WRITE (LIS2,722) (I,I=1,NSYMC) 722 FORMAT('0Multiplication table for symmetry operations'/ + (5X,24I3/)) DO 730 I=1,NSYMC 730 WRITE (LIS2, 732) I, (MULS(I,J),J=1,NSYMC) 732 FORMAT(I3, 2X, 24I3) ENDIF RETURN END SUBROUTINE COTRA1 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 (IPR1, IFILE(6)), (LIS1, IFILE(7)), (LIS2, IFILE(8)) EQUIVALENCE (IBINFO, IFILE(11)), (IBINFC, IFILE(12)) EQUIVALENCE (ITS14, IFILE(14)), (ITS15, IFILE(15)) EQUIVALENCE (IBINFF, IFILE(16)) LOGICAL SWPRI EQUIVALENCE (SWPRI, SWITCH(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 (KFA = 20671) PARAMETER (MAXBUF = 198) COMMON / / FA(KFA), FB(KFA), + BUFFC2(MAXBUF), BUFFO(MAXBUF), BUFFFT(MAXBUF), + FITFC2(51), FITFO(3), FITFFT(5), + NITFC2, KENDFC, NITFO, KENDFO DIMENSION BLACOM(42000) EQUIVALENCE (BLACOM(1), FA(1)) DIMENSION BUF(2,24) EQUIVALENCE (EPSIL2, FITFC2(1)), (SFTOT2, FITFC2(2)), + (SFPAR2, FITFC2(3)), (BUF(1,1),FITFC2(4)) PARAMETER (MXP=50, MXP1=MXP + 1) COMMON /PIEK/ RHO(MXP1), X1(3,MXP1), IBOTS(MXP1), BOTS(MXP1), NPK COMMON /COTRA/ MULS(48,48), INDX(3,3), IDIMF, TT(3,48), * EFST(48), FHST(48), PHST(48), CELLN(6), MAXHKL(3), * P1SQ, PSQ, EMIN, EO2AV, SCX, BOVX, DAMPX, SMM, IOMAPS, SMAX, * IRR(3,3,48), IORG, BBB, D2R, R2D, * AMULT, ASYMM, ALATT, ASYMCL, NSYMC, ASYMC DIMENSION SHKL(3,24), ISHKL(3,48), LHK(3), M3D(3), IORGX(14) DIMENSION MHKL(3), M2HKL(3), KM3D(3) LOGICAL SWHKL, SWHKL2 DIMENSION MSS(4) LOGICAL LMSS(4), LMSS4, LMSS5, LMSS6, LMSS7 EQUIVALENCE (LMSS(1), LMSS4), (LMSS(2), LMSS5), * (LMSS(3), LMSS6), (LMSS(4), LMSS7) EQUIVALENCE (MSS(1), MSS4), (MSS(2), MSS5), * (MSS(3), MSS6), (MSS(4), MSS7) DATA ZERO /0.0 / , IZRO, NNN / 0, 999/ DATA IORGX /1,1,1,2, 2, 5, 5, 3, 3, 3, 3, 3, 4, 4/ DATA I1, I2 / 0, 0/ IORG = IORGX(ILAUE) CALL KERNZI (0, INDX, 9) CALL KERNZI (0, KM3D, 3) CALL KERNZI (0, M2HKL, 3) CALL KERNAI (MAXHKL, MHKL, 3) CALL KERNZI (0, MAXHKL, 3) IF (IPOLA .NE. 8) GOTO 102 DO 101 I = 1, 3 CELLN(I) = 1.4 * CELL(I) CELLN(I+3) = 90. M2HKL(I) = 2 * MHKL(I) 101 M3D(I) = 6. * SMAX * CELL(I) CELLN(6) = 120. GOTO 110 102 DO 103 I=1,3 M3D(I) = SMAX * 2.31 * CELL(I) M2HKL(I) = MHKL(I) IF (ILAUE .EQ. 3) M3D(I) = SMAX * 2.0 * CELL(I) IF (IUNIQ .EQ. I) M3D(I) = SMAX * 2.0 * CELL(I) CELLN(I+3)= CELL(I+3) 103 CELLN(I) = CELL(I)/2.0 GOTO (110,104,108,104,104), IORG 104 CELLN(1) = CELL(1) CELLN(2) = CELLN(1) M3D(1) = SMAX * 4.0 * CELL(1) M3D(2) = M3D(1) M3DMMM = SQRT ( 0.5 * FLOAT( KFA - 300) ) KM3D(1) = MIN0 (M3D(1), M3DMMM) KM3D(2) = KM3D(1) KM3D(3) = MHKL(3) M2HKL(1) = 2 * MHKL(1) M2HKL(2) = 2 * MHKL(1) IF (IORG .EQ. 2) GOTO 110 CELLN(3) = CELL(1) M3D(3) = M3D(1) M2HKL(3) = MHKL(1) IF (IORG .EQ. 4) GOTO 110 M3D(1) = SMAX * 4.62 * CELL(1) M3D(2) = M3D(1) M3D(3) = M3D(1) M2HKL(1) = 2 * MHKL(1) M2HKL(2) = 2 * MHKL(1) M2HKL(3) = 2 * MHKL(1) GOTO 110 108 CELLN(1) = CELL(1) CELLN(2) = CELLN(1) M3D(1) = SMAX * 4.62 * CELL(1) M3DMMM = SQRT ( 0.5 * FLOAT( KFA - 300) ) M3D(1) = MIN0 (M3D(1), M3DMMM) M3D(2) = M3D(1) M2HKL(1) = 2 * MHKL(1) M2HKL(2) = 2 * MHKL(1) M2HKL(3) = 2 * MHKL(1) 110 WRITE (LIS2, 112) CELL, VOLUM, CELLN, IORG, M3D, M2HKL, KM3D 112 FORMAT (/' CELL: ', 3F8.3, 3F7.2, ' VOLUME:' ,F9.2/ * ' CELLN:', 3F8.3, 3F7.2, ' transformed cell'/ * ' (Origin keys: IORG =',I2, ')' / * ' (Expected maximum transformed indices: ', 3I3,' )'/ * ' ( [ TEMP: M2HKL=',3I3,' TEMP: KM3D=', 3I3,' )' /) IF (IPOLA .EQ. 0) THEN IDIMF=3 CALL BINOFF (28, IBINFF, 'BINFFT', FITFFT, 5, BUFFFT, KENDFF) GOTO 130 ENDIF WRITE (LIS2, 119) 119 FORMAT(' POLAR SPACE GROUP - NO 3-D SEARCH') IDIMF=2 IF (IPOLA.EQ.3 .OR. IPOLA.EQ.5 .OR. IPOLA.EQ.6) IDIMF=1 CALL FILINQ (ITS15, 'SCRA15', 'UNFORMATTED', 'SCRATCH', KINQ) GOTO (122, 120, 123, 121, 122, 121, 120, 121), IPOLA 120 INDX(1,1) = 1 INDX(2,3) = 1 INDX(3,2) =-1 M3D(2)=0 GOTO 130 121 I=1 GOTO 125 122 I=2 GOTO 125 123 I=3 125 J=MOD(I,3)+1 K=6-I-J INDX(1,I) = 1 INDX(2,J) = 1 INDX(3,K) = 1 M3D(K)=0 IF(IDIMF.EQ.1) M3D(J)=0 130 CONTINUE CALL KERNZA (0.0, FA, KFA) CALL KERNZA (0.0, FB, KFA) KLAAR = 0 131 M2 = 2 * M3D(3) + 1 M1 = M2 * (2 * M3D(2) + 1) M3 = M3D(2) * M2 NEED = M3D(1) * M1 + M3D(2) * M2 + M3D(3) WRITE (LIS2, 132) M3D,NEED 132 FORMAT(' Max indices and storage for general scan: ',3I3,I8) M5 = 2 * M3D(2) + 1 MS1 = (M3D(2) + 1) * M2 MS2 = (M3D(1) + 1) * M2 MS3 = (M3D(1) + 1) * M5 MS12= MS1 + MS2 MS123= MS12 + MS3 KFA1 = KFA - 100 IF (MS123 .GT. KFA1) * WRITE (LIS2, 133) MS1, MS2, MS3, MS123, KFA1 133 FORMAT(/' Storage 2D scans: BC AC AB Sum LIMIT:'/ * ' -----------------', 3I6, I7, I8 /) KLAAR = KLAAR + 1 IF (MS123 .GT. KFA1 .AND. KLAAR .LE. 4) THEN KLA = NINT ( 100. * SQRT ( FLOAT( MS123) / FLOAT( KFA1 ) ) ) DO 134 I = 1, 3 M3D(I) = 100 * M3D(I) / KLA 134 CONTINUE GOTO 131 ENDIF IF (MS123 .GT. KFA1) THEN WRITE (LIS1, 133) MS1, MS2, MS3, MS123, KFA1 WRITE (LIS1, FMT='(/'' Rerun TRACOR with lower sinth/lam'')') CALL KERROR ('Storage problems: KFA1.', 133, 'COTRA1') ENDIF MSTOT = MS123 NKFA = MIN0 (NEED, KFA) IDIMFX = IDIMF IF (ILAUE .EQ. 1) IDIMFX = 4 IF (NEED.LE.KFA .AND. IDIMFX.EQ.3) IDIMFX = 5 MSS4 = MS123 + MS2 MSS5 = MSS4 + MS2 MSS6 = MSS5 + MS2 MSS7 = MSS6 + MS2 DO 135 M = 1,4 LMSS(M) = .TRUE. IF (MSS(M) .LE. KFA - 100) THEN LMSS(M) = .FALSE. MSTOT = MSS(M) ENDIF 135 CONTINUE NFOUR = 0 KBOT = 1 KTOP = NEED NPASS = 1 IPASS = 0 IF (NEED .LE. KFA .OR. ILAUE .EQ. 1) GOTO 140 CALL FILINQ (ITS14, 'SCRA14', 'UNFORMATTED', 'SCRATCH', KINQ) KTOP = KFA - MSTOT NPASS = (NEED - 1 + MSTOT) / KFA + 1 140 NACC = 0 NRIN = 0 NRINP = 0 SEOBS = 0.0 SES2 = 0.0 TOTREF= 0.0 NUNOBS = 0 IADDM = -9999 NREF14 = 0 SWHKL = .FALSE. SWHKL2 = .FALSE. IF (.NOT. SWPRI) SWHKL2 = .TRUE. ITREF = 0 ITREFE = 0 ITREFS = 0 IF (SWPRI) WRITE (LIS2, 150) 150 FORMAT('0 h k l Eobs eps SF2 SF2P', * ' Fp phase Fp phase .....') 160 CALL BINIFF (0,IBINFC,'BINFC2',FITFC2,NITFC2,BUFFC2,KENDFC) IF (KENDFC .LT. 0) GOTO 500 NSREF = NSYMM / NINT (EPSIL2) NRIN = NRIN + NSREF NRINP = NRINP + 1 CALL BINIFF (0,IBINFO,'BINFO',FITFO,NITFO,BUFFO,KENDFO) IF (KENDFO .LT. 0) CALL KERNER(-2, 'COTRA1') IF (FITFC2(1) .LT. 0.) THEN ITREFS = ITREFS + 1 GOTO 160 ENDIF CALL HKLC1U (FITFO(1), SHKL) CALL HKLSTL (SHKL, STL, STL2) IF (STL .GT. SMAX) GOTO 160 IF (FITFO(2) .LE. 6.01*FITFO(3)) THEN NUNOBS = NUNOBS + 1 FITFO(2) = FITFO(2) / 4. ENDIF EF = SCX * FITFO(2) * EXP(BOVX * STL2) / SQRT(SFTOT2 * ALATT) IF (EF.GT.3.0) EF=3.0 IF (EF.LT.EMIN) THEN ITREFE = ITREFE + 1 GOTO 160 ENDIF NACC = NACC + NSREF SXX=STL/SMM IF(SXX.GT.SMM) SXX=1.0 SXX = SQRT(SXX) * EXP(DAMPX * STL2) / SQRT(SFPAR2 * ALATT) CALL HKLEX1 (SHKL, SHKL) CALL HKLEX2 (SHKL, ISHKL, IEPS, IEPS2) EPSIL = IEPS CALL KERF2I (SHKL, ISHKL, 3*NSYMM) IF (SWHKL2) GOTO 172 IF (SWHKL) THEN SWHKL2 = .TRUE. SWHKL = .FALSE. GOTO 172 ENDIF IF (ISHKL(1,1).EQ.0 .OR. ISHKL(2,1).EQ.0 .OR. * ISHKL(3,1).EQ.0) GOTO 172 IF (ISHKL(1,1).EQ.ISHKL(2,1) .OR. ISHKL(2,1).EQ.ISHKL(3,1) .OR. * ISHKL(3,1).EQ.ISHKL(1,1)) GOTO 172 IF (ISHKL(1,1).EQ.-ISHKL(2,1) .OR. ISHKL(2,1).EQ.-ISHKL(3,1) .OR. * ISHKL(3,1).EQ.-ISHKL(1,1)) GOTO 172 SWHKL = .TRUE. 172 IF (SWPRI .AND. (SWHKL .OR. NRINP.LT.22)) WRITE * (LIS2,180) (ISHKL(IJJ,1),IJJ=1,3), EF, (FITFC2(I),I=1,NITFC2) 180 FORMAT (3I4, F5.2, F3.0, 2F7.0, 1X, 7(F7.2,F6.1) * / 35X, 7(F7.2,F6.1)) SUMES=0.0 DO 190 IS=1,NSYMM PHST(IS) = BUF(2,IS) IF (PHST(IS) .LT. 0.0) PHST(IS) = 360. + PHST(IS) IF (BUF(1,IS) .GE. 0.0) THEN EFST(IS) = BUF(1,IS) * SXX SUMES = SUMES + EFST(IS)**2 ELSE ISK = NINT (-BUF(1,IS)) EFST(IS) = EFST(ISK) ENDIF 190 CONTINUE IF (ICENT.EQ.1) GOTO 196 DO 195 IS = NSYMM+1 , NSYMC DO 194 I=1,3 194 ISHKL(I,IS)= - ISHKL(I,IS-NSYMM) EFST(IS) = EFST(IS-NSYMM) 195 PHST(IS) = 360. - PHST(IS-NSYMM) 196 CONTINUE TOTREF = TOTREF + 1. SUMES = SUMES / FLOAT(NSREF) SES2 = SES2 + SUMES EF2 = EF**2 SEOBS = SEOBS + EF2 EF2 = EF2 - PSQ * SUMES - BBB DO 400 I=2,NSYMC DO 390 II=1,NSYMM PTRA = SHKL(1,II)*TT(1,I)+ SHKL(2,II)*TT(2,I)+ SHKL(3,II)*TT(3,I) IJ = MULS(II,I) PTRA = PHST(II) - PHST(IJ) - 360.*PTRA EFTRA = EF2 * EFST(II) * EFST(IJ) / EPSIL IF (EFTRA .LT. 0.) THEN EFTRA = ABS(EFTRA) PTRA = PTRA + 180. ENDIF PTRA = AMOD(-PTRA, 360.) DO 301 J=1,3 301 LHK(J) = ISHKL(J,II) - ISHKL(J,IJ) GOTO (302, 303, 303, 305, 304), IORG 302 LHK(1)=LHK(1)/2 LHK(2)=LHK(2)/2 303 LHK(3)=LHK(3)/2 GOTO 305 304 IF (IPOLA .NE. 8) GOTO 305 LHK(1) = LHK(1) - LHK(2) LHK(2) = LHK(2) - LHK(3) LHK(3) = 0 305 DO 306 J=1,3 IF (IABS(LHK(J)) .LE. M3D(J)) GOTO 306 ITREF = ITREF + 1 GOTO 390 306 CONTINUE DO 308 J=1,3 308 MAXHKL(J) = MAX0 (MAXHKL(J), IABS(LHK(J))) IF(LHK(1)) 320, 310, 330 310 IF(LHK(2)) 320, 315, 330 315 IF(LHK(3)) 320, 390, 330 320 DO 322 J=1,3 322 LHK(J)=-LHK(J) PTRA = - PTRA 330 IADD = LHK(3) + LHK(2) * M2 + LHK(1) * M1 IF (PTRA .LT. 0.) PTRA = PTRA + 360. IF (SWHKL) WRITE (LIS2, 331) (ISHKL(L,I),L=1,3), * (ISHKL(L,II),L=1,3), (ISHKL(L,IJ),L=1,3), LHK, IADD 331 FORMAT (' HKL for: I,II,IJ,M; IADD: ', 4(3I3,2X), I6) GOTO (338, 338, 343, 341, 338), IDIMFX CALL KERROR ('Impossible ...', 338, 'COTRA1') 338 FA(IADD) = FA(IADD) + EFTRA * COS(PTRA * D2R) FB(IADD) = FB(IADD) + EFTRA * SIN(PTRA * D2R) GOTO 390 341 FITFFT(1) = LHK(1) FITFFT(2) = LHK(2) FITFFT(3) = LHK(3) FITFFT(4) = EFTRA FITFFT(5) = PTRA CALL BINOFF(0, IBINFF, 'BINFFT', FITFFT, 5, BUFFFT, KENDFF) GOTO 390 343 CONTINUE IF (LHK(1) .EQ. 0) THEN IADD = LHK(2) * M2 + LHK(3) GOTO 338 ENDIF IF (LHK(2) .EQ. 0) THEN IADD = MS1 + LHK(1) * M2 + LHK(3) GOTO 338 ENDIF IF (LHK(3) .EQ. 0) THEN IADD = MS12 + LHK(1) * M5 + LHK(2) GOTO 338 ENDIF IF (ILAUE .LE. 3) GOTO 341 IF (LMSS4) GOTO 345 IF (LHK(1) .EQ. LHK(2)) THEN IADD = MS123 + LHK(1) * M2 + LHK(3) GOTO 338 ENDIF IF (LMSS5) GOTO 345 IF (LHK(1) .EQ. -LHK(2)) THEN IADD = MSS4 + LHK(1) * M2 + LHK(3) GOTO 338 ENDIF IF (LMSS6) GOTO 345 IF (LHK(2) .EQ. -2 * LHK(1)) THEN IADD = MSS5 + LHK(1) * M2 + LHK(3) GOTO 338 ENDIF IF (LMSS7) GOTO 345 IF (2 * LHK(2) .EQ. -LHK(1)) THEN IADD = MSS6 - LHK(2) * M2 + LHK(3) GOTO 338 ENDIF 345 IF (IADD .GT. KTOP) GOTO 380 IADD = IADD + MSTOT GOTO 338 380 WRITE (ITS14) IADD, EFTRA, PTRA IADDM = MAX0(IADDM,IADD) NREF14 = NREF14 + 1 390 CONTINUE 400 CONTINUE GOTO 160 500 IF (ITREF .GT. 0) WRITE (LIS1, FMT='( * '' Nr of refl. bypassed (HKL > M3D)'', I6)') ITREF IF (ITREFE .GT. 0) WRITE (LIS1, FMT='( * '' Nr of refl. bypassed (E < EMIN) '', I6)') ITREFE IF (ITREFS .GT. 0) WRITE (LIS1, FMT='( * '' Nr of refl. bypassed (S > SMAX) '', I6)') ITREFS IF (NREF14 .GT. 0) THEN WRITE(LIS2, FMT = '('' max. indices-address : '', I10, * '' (calc:'',I10,'')'',/, * '' Nr. of refl. to file ITS14: '', I10,/)') * IADDM, NEED, NREF14 NEED = IADDM WRITE (LIS1, FMT='('' Max. transformed hkl: '', 3I3)') MAXHKL ENDIF WRITE(LIS2, FMT = '('' Nr. of weak reflections : '', I10)') * NUNOBS WRITE (LIS2, FMT='('' Max. transformed hkl: '', 3I3)') MAXHKL IF (NPASS .GT. 1) WRITE (IPR1, FMT= '('' Approximately'', I3, * '' passes needed for output/input scratch file: '')') NPASS ITREF = 0 GOTO 630 600 IDIMFX = 6 KBOT=KTOP+1 KTOP = MIN0 (KTOP + KFA, NEED) NKFA = KTOP - KBOT + 1 MSTOT = 0 CALL KERNZA (0.0, FA, NKFA) CALL KERNZA (0.0, FB, NKFA) REWIND ITS14 NREFFI = 0 NREFFO = 0 620 READ (ITS14) IADD, F, P NREFFI = NREFFI + 1 IF (IADD.LT.1) GOTO 625 IF (IADD.LT.KBOT .OR. IADD.GT.KTOP) GOTO 620 NREFFO = NREFFO + 1 J = IADD - KBOT + 1 FA(J) = FA(J) + F * COS(P * D2R) FB(J) = FB(J) + F * SIN(P * D2R) GOTO 620 625 CONTINUE IF (IPASS .EQ. 0) WRITE(LIS2, FMT = * '('' Nr. of refl. from file ITS14:'',I9 )') NREFFI WRITE(LIS2, FMT = * '('' Nr. of reflections stored :'',I9 )') NREFFO 630 IF (IDIMFX .EQ. 4) GOTO 650 DO 640 K = 1, NKFA J = K + KBOT - 1 F = SQRT (FA(K)**2 + FB(K)**2) IF (F .LT. 0.000001) GOTO 640 IF (F .LT. 0.01) THEN ITREF = ITREF + 1 GOTO 640 ENDIF NFOUR = NFOUR + 1 GOTO (632, 632, 635, 631, 632, 632), IDIMFX 631 CALL KERROR ('Impossible ...', 631, 'COTRA1') 632 IADD = J 633 I1 = (IADD + M3 + M3D(3)) / M1 IADD = IADD - I1 * M1 I2 = (IADD + M3 + M3D(3)) / M2 - M3D(2) I3 = IADD - I2 * M2 IF (IDIMF.EQ.3) GOTO 639 I = I1 * INDX(1,1) + I2 * INDX(1,2) + I3 * INDX(1,3) I2= I1 * INDX(2,1) + I2 * INDX(2,2) + I3 * INDX(2,3) I1=I IF (IDIMF.EQ.1) WRITE(ITS15) I1,FA(K),FB(K) IF (IDIMF.EQ.2) WRITE(ITS15) I1,I2,FA(K),FB(K) GOTO 640 635 IF (J. GT. MSTOT) GOTO 638 IF (J .LE. MS1) GOTO 632 IF (J .LE. MS12) THEN I1 = (J - MS1 + M3D(3)) / M2 I2 = 0 I3 = J - MS1 - I1 * M2 GOTO 639 ENDIF IF (J .LE. MS123) THEN I1 = (J - MS12 + M3D(2)) / M5 I2 = J - MS12 - I1 * M5 I3 = 0 GOTO 639 ENDIF IF (ILAUE .LE. 3) CALL KERROR (' Kan niet ', 635, 'COTRA') IF (J .LE. MSS4) THEN I1 = (J - MS123 + M3D(3)) / M2 I2 = I1 I3 = J - MS123 - I1 * M2 GOTO 639 ENDIF IF (J .LE. MSS5) THEN I1 = (J - MSS4 + M3D(3)) / M2 I2 = -I1 I3 = J - MSS4 - I1 * M2 GOTO 639 ENDIF IF (J .LE. MSS6) THEN I1 = (J - MSS5 + M3D(3)) / M2 I2 = -2 * I1 I3 = J - MSS5 - I1 * M2 GOTO 639 ENDIF IF (J .LE. MSS7) THEN I2 = (J - MSS6 + M3D(3)) / M2 I1 = 2 * I2 I3 = J - MSS6 - I2 * M2 I2 = - I2 GOTO 639 ENDIF CALL KERROR (' Kan niet ', 635, 'COTRA') 638 IADD = J - MSTOT GOTO 633 639 FITFFT(1) = I1 FITFFT(2) = I2 FITFFT(3) = I3 FITFFT(4) = F P = ATAN2(FB(K),FA(K)) * R2D IF (P .LT. 0.0) P = P + 360. FITFFT(5) = P CALL BINOFF(0, IBINFF, 'BINFFT', FITFFT, 5, BUFFFT, KENDFF) 640 CONTINUE IF (NREF14 .EQ. 0) GOTO 650 IPASS = IPASS + 1 IF (NPASS .GT. 1) THEN IF (IPASS .EQ. 1) WRITE (ITS14) IZRO,ZERO,ZERO WRITE (IPR1, FMT='('' Pass'', I3, '' completed'')') IPASS WRITE (LIS2, FMT='('' Pass'', I3, '' completed'')') IPASS IF (KTOP .LT. NEED) GOTO 600 ENDIF 650 IF (IDIMF.EQ.1) WRITE(ITS15) NNN,ZERO,ZERO IF (IDIMF.EQ.2) WRITE(ITS15) NNN,IZRO,ZERO,ZERO IF (IDIMF.NE.3) THEN REWIND ITS15 ELSE CALL BINOFF(-1, IBINFF, 'BINFFT', FITFFT, 5, BUFFFT, KENDFF) ENDIF IF (ITREF .GT. 0) WRITE (LIS1, FMT='( * '' Nr of refl. bypassed (F < 0.01) '', I6)') ITREF IF (SWPRI) THEN WRITE (LIS2, 902) NRIN,NACC, NINT(TOTREF), NFOUR 902 FORMAT (' Nr of reflections input:', I9 / * ' Nr of reflections accepted:', I6,' with', I6, ' independent.' * /' Nr of terms that contributed to the Fourier map:', I8) WRITE (LIS2, 903) MAXHKL 903 FORMAT (' Maximum indices for Fourier synthesis: ', 3I3) ENDIF CALL FILCLO (ITS14, 'DELETE') SEOBS = SEOBS / TOTREF SES2 = SES2 / TOTREF PPES = P1SQ * EO2AV + (1.-P1SQ) EOOR = EO2AV - PSQ * PPES - BBB ENUL = SEOBS - PSQ * SES2 - BBB WRITE (LIS2, 910) EO2AV,SEOBS,PPES,SES2,EOOR,ENUL 910 FORMAT (/' Average of theory observed' / +' EOBS**2 ',F10.4,5X,F10.4 / +' ES**2 ',F10.4,5X,F10.4 / +' EO2-IMV-OR ',F10.4,5X,F10.4 //) RETURN END SUBROUTINE D12FOU 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 (ITS15, IFILE(15)) 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 (NPKM2=100,MNGR2=163, MNGR1=200,LC24=24) COMMON / / SICO(1250), RHOSUM(MNGR2,MNGR2), RHOMAX(NPKM2), + IC1(MNGR1), IC2(MNGR1), IMAX(NPKM2), JMAX(NPKM2), + XMAX(NPKM2), YMAX(NPKM2), RHOP(NPKM2), IRHO(LC24) DIMENSION BLACOM(42000) EQUIVALENCE (BLACOM(1), SICO(1)) PARAMETER (MXP=50, MXP1=MXP + 1) COMMON /PIEK/ RHO(MXP1), SH(3,MXP1), IBOTS(MXP1), BOTS(MXP1), NPK COMMON /COTRA/ MULS(48,48), INDX(3,3), IDIMF, TT(3,48), * EFST(48), FHST(48), PHST(48), CELLN(6), MAXHKL(3), * P1SQ, PSQ, EMIN, EO2AV, SCX, BOVX, DAMPX, SMM, IOMAPS, SMAX, * IRR(3,3,48), IORG, BBB, D2R, R2D, * AMULT, ASYMM, ALATT, ASYMCL, NSYMC, ASYMC DIMENSION AL(3) CALL SICOT (SICO, 1250) IF (IPOLA .EQ. 8) THEN CALL KERNAB (CELLN, AL, 3) ANGLE = 120. GOTO 135 ENDIF DO 112 L=1,2 AL(L)=0. DO 110 I=1,3 X=0. DO 108 J=1,3 108 X = X + INDX(L,J) * RRMAT(I,J) 110 AL(L) = AL(L) + X * INDX(L,I) 112 AL(L)=SQRT(AL(L)) AL(3)=0. DO 122 I=1,3 X=0. DO 121 J=1,3 121 X = X + (INDX(1,J) - INDX(2,J)) * RRMAT(I,J) 122 AL(3) = AL(3) + X * (INDX(1,I) - INDX(2,I)) AL(3)=(AL(1)*AL(1)+AL(2)*AL(2)-AL(3))/AL(1)/AL(2) ANGLE = ACOS(0.5 * AL(3)) * R2D IF (IDIMF .NE. 1) GOTO 135 NGR1=AL(1)/0.21 IF (NGR1.GT.200) NGR1=200 WRITE (LIS1, 127) AL(1),NGR1 127 FORMAT('0 FOURIER IS CALCULATED ALONG AXIS OF',F8.3,' A' / +' NUMBER OF GRID POINTS IN THIS DIRECTION (QX):',I5 ) CALL D1FOUR (IOMAPS,NGR1,SCALE) GOTO 137 135 NGR1=AL(1)/0.25 IF (NGR1.GT.163) NGR1=163 NGR2=AL(2)/0.25 IF (NGR2.GT.163) NGR2=163 WRITE (LIS1, 136) AL(1), AL(2), ANGLE, NGR1, NGR2 136 FORMAT('0FOURIER IS CALCULATED IN PLANE OF AXES OF', /, 6X, +F8.3,' AND',F8.3,' A, ANGLE = ',F8.3 / +' NUMBER OF GRID POINTS IN FIRST DIRECTION (QX):',I5/ +' NUMBER OF GRID POINTS IN SECOND DIRECTION (QY):',I5) CALL D2FOUR (IOMAPS,NGR1,NGR2,SCALE) 137 WRITE (LIS1, 138) SCALE 138 FORMAT(' FOURIER SCALEFACTOR = ',F10.4) CALL FILCLO (ITS15, 'DELETE') CALL KERNZA (0.0, SH, 3 * MXP1) RHOMIN = 0.5 * RHOP(1) DO 206 I=1,NPK RHO(I) = RHOP(I) IF (RHO(I) .LT. RHOMIN) THEN NPK=I-1 RETURN ENDIF DO 204 K=1,3 IF (INDX(1,K) .NE. 0) SH(K,I) = XMAX(I) 204 CONTINUE IF (IDIMF .EQ. 1) GOTO 206 DO 205 K=1,3 IF (INDX(2,K) .NE. 0) SH(K,I) = YMAX(I) 205 CONTINUE 206 CONTINUE RETURN END SUBROUTINE COTRA3 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 (LIS1, IFILE(7)), (LIS2, IFILE(8)) EQUIVALENCE (ICRYS, IFILE(3)), (ITS15, IFILE(15)) 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 (MXP=50, MXP1=MXP + 1) COMMON /PIEK/ RHO(MXP1), SH(3,MXP1), IBOTS(MXP1), BOTS(MXP1), NPK DIMENSION ISAME(MXP) COMMON /COTRA/ MULS(48,48), INDX(3,3), IDIMF, TT(3,48), * EFST(48), FHST(48), PHST(48), CELLN(6), MAXHKL(3), * P1SQ, PSQ, EMIN, EO2AV, SCX, BOVX, DAMPX, SMM, IOMAPS, SMAX, * IRR(3,3,48), IORG, BBB, D2R, R2D, * AMULT, ASYMM, ALATT, ASYMCL, NSYMC, ASYMC DIMENSION ORGVEC(3,8), SHD(3), SHDT(3), SHDTX(3), PLATT(3,4) LOGICAL SWORG DATA NORIG, SWORG /1 , .FALSE. / DO 238 I=1,NPK IF (IPOLA .EQ. 8) THEN SH(3,I) = - SH(2,I) SH(2,I) = SH(2,I) - SH(1,I) SH(3,I) = SH(3,I) - SH(1,I) SH(2,I) = SH(2,I) - SH(1,I) SH(1,I) = 0.0000001 CALL KERNZA (0., PLATT, 3) CALL KERNZA (0., ORGVEC, 3) ENDIF DO 211 J=1,3 SHJI = AMOD(2.0+SH(J,I),1.0) 211 SH(J,I) = SHJI GOTO (212,220,215,230,230), IORG 212 SH(1,I)=SH(1,I)/2.0 SH(2,I)=SH(2,I)/2.0 215 SH(3,I)=SH(3,I)/2.0 GOTO 238 220 IF (SH(1,I).LT.0.5) GOTO 215 SH(1,I)=SH(1,I)-0.5 SH2I = AMOD(SH(2,I)+0.5, 1.0) SH(2,I) = SH2I GOTO 215 230 IF (SH(1,I).LT.0.5) GOTO 238 DO 235 J=1,3 SHJI = AMOD(SH(J,I)+0.5, 1.0) SH(J,I) = SHJI 235 CONTINUE 238 CONTINUE IF (IPOLA .EQ. 8) GOTO 242 CALL RDCRYB (ICRYS, 'NORIG', KEND) IF (KEND .EQ. -1) * CALL KERROR (' No NORIG record on CRYSDA file', 238, 'COTRA3') READ(CHIN, FMT = '(10X, I10)') NORIG DO 239 N = 1,NORIG 239 READ (ICRYS, FMT = '(10X, 3F10.7)') (ORGVEC(J,N), J= 1,3) CALL FILCLO (ICRYS, 'KEEP') IF (NORIG * NLATT .EQ. 1) GOTO 248 DO 241 N = 1, NLATT PLATT(1,N) = TLATT(1,N) PLATT(2,N) = TLATT(2,N) PLATT(3,N) = TLATT(3,N) IF (IPOLA.EQ.1 .OR. IPOLA.EQ.3 .OR. IPOLA.EQ.5) PLATT(1,N) = 0. IF (IPOLA.EQ.2 .OR. IPOLA.EQ.3 .OR. IPOLA.EQ.6) PLATT(2,N) = 0. IF (IPOLA.EQ.4 .OR. IPOLA.EQ.5 .OR. IPOLA.EQ.6) PLATT(3,N) = 0. 241 CONTINUE 242 DO 246 N1 = 1,NPK-1 IF (RHO(N1) .LT. 0.0) GOTO 246 DO 245 N2 = N1+1, NPK IF (RHO(N2) .LT. 0.0) GOTO 245 DO 243 J = 1,3 243 SHD(J) = SH(J,N1) - SH(J,N2) DO 244 IL = 1, NLATT SHDT(1) = SHD(1) + PLATT(1,IL) SHDT(2) = SHD(2) + PLATT(2,IL) SHDT(3) = SHD(3) + PLATT(3,IL) DO 244 N = 1,NORIG CALL DISTSQ (SHDT, ORGVEC(1,N), 0.4, SHDTX, DIST2) IF (DIST2 .LT. 0.16) THEN IF (.NOT. SWORG) * WRITE (LIS2, FMT = '('' TEMP: Removed (close contact):'')') WRITE (LIS2,FMT='('' shiftvector'',I3,'' ='',3F8.4,'' ='',I3)') * N2, (SH(J,N2),J=1,3), N1 RHO(N2) = -9999. SWORG = .TRUE. GOTO 245 ENDIF 244 CONTINUE 245 CONTINUE 246 CONTINUE IF (.NOT. SWORG) GOTO 248 NPK1 = 0 DO 247 I = 1,NPK IF (RHO(I) .LT. 0.0) GOTO 247 NPK1 = NPK1 + 1 CALL KERNAB (SH(1,I), SH(1,NPK1), 3) RHO(NPK1) = RHO(I) 247 CONTINUE NPK = NPK1 248 CONTINUE CALL DISCHK (NA, ISAME, MXP) NBOTS = 0 NBOTSP = 0 DO 270 I=1,NPK IF (IBOTS(I) .NE. 0) NBOTS = NBOTS + 1 IF (BOTS(I) .LE. 3.00) NBOTSP = NBOTSP + 1 270 CONTINUE IF (NBOTS .EQ. 0) WRITE (LIS1, FMT='(/'' No close contacts''/)') IF (NBOTS .GT. 0) WRITE (LIS1, 275) NBOTS 275 FORMAT (I4, ' of the shift vectors lead to bad contacts.'/ * ' Number of contacts (nr) and shortests distance:'/ * 42X, 'nr dist') IF (NBOTS .EQ. 0 .AND. NBOTSP .GT. 0) WRITE (LIS1, FMT='( * 33X, ''Shortest distance:'')') RHOX = RHO(1) * .5 WRITE(LIS1, FMT='('' Possible shift vectors:''/ * '' TR= FOM Tx Ty Tz'' )') NPKXX = NPK DMAX=2.999 MSAME = NA / 3 DO 280 I = 1, NPKXX IF (I .LE. 2 .OR. RHO(I) .GT. RHOX) NPK = I IF (IBOTS(I) .GT. 0) THEN IF (I .LE. 2 .OR. RHO(I) .GT. RHOX) * WRITE (LIS1,276) I,RHO(I),(SH(J,I),J=1,3), IBOTS(I), BOTS(I) 276 FORMAT (I4, F7.0, 3F9.5, I6, F7.2) ELSEIF (BOTS(I) .GT. DMAX) THEN IF (I .LE. 2 .OR. RHO(I) .GT. RHOX) * WRITE (LIS1, 276) I, RHO(I), (SH(J,I),J=1,3) ELSE IF (I .LE. 2 .OR. RHO(I) .GT. RHOX) * WRITE (LIS1, 277) I, RHO(I), (SH(J,I),J=1,3), BOTS(I) 277 FORMAT (I4, F7.0, 3F9.5, 6X, F7.2) ENDIF IF (ISAME(I).GT.MSAME) * WRITE (LIS1,FMT='('' Symmetry ?? shift nr. '', I3)') I IBOTS(I) = MIN0 ( IBOTS(I), 79) 280 CONTINUE RHO(NPK+1)=-9999.0 RETURN END SUBROUTINE R2TRAC (SEDEL2) DIMENSION SEDEL2(50,2) 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 (IPR1, IFILE(6)), (LIS1, IFILE(7)), (LIS2, IFILE(8)) EQUIVALENCE (IBINFO, IFILE(11)), (IBINFC, IFILE(12)) 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 (MAXBUF = 198) DIMENSION BUFFC2(MAXBUF), BUFFO(MAXBUF), + FITFC2(51), FITFO(3) INTEGER NITFC2, KENDFC, NITFO, KENDFO DIMENSION BUF(2,24) EQUIVALENCE (SFTOT2, FITFC2(2)), + (SFPAR2, FITFC2(3)), (BUF(1,1),FITFC2(4)) PARAMETER (MXP=50, MXP1=MXP + 1) COMMON /PIEK/ RHO(MXP1), X1(3,MXP1), IBOTS(MXP1), BOTS(MXP1), NPK COMMON /COTRA/ MULS(48,48), INDX(3,3), IDIMF, TT(3,48), * EFST(48), FHST(48), PHST(48), CELLN(6), MAXHKL(3), * P1SQ, PSQ, EMIN, EO2AV, SCX, BOVX, DAMPX, SMM, IOMAPS, SMAX, * IRR(3,3,48), IORG, BBB, D2R, R2D, * AMULT, ASYMM, ALATT, ASYMCL, NSYMC, ASYMC DIMENSION SHKL(3,24), SORTR2(50) WRITE (LIS1, 120) 120 FORMAT (/' Calculate R2 values for all shift vectors'/ * ' (ignore possible overlap between molecules)') CALL BINIFF (1, IBINFO, 'BINFO', FITFO, NITFO, BUFFO, KENDFO) CALL BINIFF (1, IBINFC, 'BINFC2',FITFC2,NITFC2,BUFFC2,KENDFC) TOTREF = 0. SEOBS2 = 0. SEFS2 = 0. SEOBS4 = 0. CALL KERNZA (0.0, SEDEL2, 100) MPK = MIN0 (50, NPK) 160 CALL BINIFF (0,IBINFC,'BINFC2',FITFC2,NITFC2,BUFFC2,KENDFC) IF (KENDFC .LT. 0) GOTO 500 CALL BINIFF (0,IBINFO,'BINFO',FITFO,NITFO,BUFFO,KENDFO) IF (KENDFO .LT. 0) CALL KERNER(-2, 'R2TRAC') IF (FITFC2(1) .LT. 0.) GOTO 160 CALL HKLC1U (FITFO(1), SHKL) CALL HKLSTL (SHKL, STL, STL2) IF (STL .GT. SMAX) GOTO 160 FNORM = EXP(BOVX * STL2) / SQRT(SFTOT2 * ALATT) EF = SCX * FITFO(2) * FNORM IF (EF.GT.3.0) EF=3.0 IF (EF.LT.EMIN) GOTO 160 EF2 = EF**2 SEOBS2 = SEOBS2 + EF2 SEOBS4 = SEOBS4 + EF2**2 CALL HKLEX1 (SHKL, SHKL) DO 190 IS=1,NSYMM PHST(IS) = BUF(2,IS) IF (BUF(1,IS) .GE. -0.01) THEN EFST(IS) = BUF(1,IS) ELSE ISK = NINT (-BUF(1,IS)) EFST(IS) = EFST(ISK) ENDIF PTRA = SHKL(1,1)*TT(1,IS)+ SHKL(2,1)*TT(2,IS)+ SHKL(3,1)*TT(3,IS) PHSTX = PHST(IS) + 360. * PTRA PHST(IS) = AMOD (PHSTX, 360.) 190 CONTINUE IF (ICENT.EQ.1) GOTO 196 DO 195 IS = NSYMM+1 , NSYMC EFST(IS) = EFST(IS-NSYMM) 195 PHST(IS) = 360. - PHST(IS-NSYMM) 196 CONTINUE TOTREF = TOTREF + 1. DO 215 I = 1, MPK FA1 = 0. FB1 = 0. DO 210 IS=1,NSYMM PTRA =SHKL(1,IS)*X1(1,I)+ SHKL(2,IS)*X1(2,I)+ SHKL(3,IS)*X1(3,I) PTRA = 360.*PTRA + PHST(IS) PTRA = AMOD(PTRA, 360.) IF (PTRA .LT. 0.) PTRA = PTRA + 360. FA1 = FA1 + EFST(IS) * COS(PTRA * D2R) IF (ICENT.EQ.1) FB1 = FB1 + EFST(IS) * SIN(PTRA * D2R) 210 CONTINUE IF (ICENT.EQ.2) FA1 = FA1 + FA1 EFS2 =(FA1**2 + FB1**2) * FNORM * FNORM IF (EFS2.GT.9.0) EFS2=9.0 SEDEL2(I,1) = SEDEL2(I,1) + ( EF2 - EFS2 ) **2 IF (I .EQ. 1) SEFS2 = SEFS2 + EFS2 215 CONTINUE GOTO 160 500 CONTINUE CALL FILCLO (IBINFC,'DELETE') CALL FILCLO (IBINFO,'KEEP') SEOBS2 = SEOBS2 / TOTREF SEFS2 = SEFS2 / TOTREF WRITE (LIS2, 910) SEOBS2,SEFS2 910 FORMAT (/' Average of observed' / +' EOBS**2 ',F10.4 / +' EFS**2 ',F10.4 ,' for first shift only') DO 919 I = 1, MPK SEDEL2(I,1) = SEDEL2(I,1) / SEOBS4 919 CONTINUE WRITE (LIS1, 923) (I, SEDEL2(I,1), I=1, MPK) WRITE (LIS2, 923) (I, SEDEL2(I,1), I=1, MPK) 923 FORMAT (/' R2 for all shift vectors '/ * ' nr R2 nr R2 ......'/ * 5(I4, F6.3) ) WRITE (LIS1, FMT='( '' R2 used for output selection.'' )') CALL KERNZA (0.0, SORTR2, 50) SORTR2(1) = SEDEL2(1,1) DO 933 I = 2, MPK DO 930 L = I, 2, -1 IF (SEDEL2(I,1) .LT. SORTR2(L-1)) THEN SORTR2(L) = SORTR2(L-1) IF (L .EQ. 2) SORTR2(L-1) = SEDEL2(I,1) ELSE SORTR2(L) = SEDEL2(I,1) GOTO 933 ENDIF 930 CONTINUE 933 CONTINUE WRITE (LIS2, 943) ( SORTR2(I), I=1, MPK) 943 FORMAT (/' R2 sorted decreasing ... '/ * ' R2 R2 R2 R2 R2'/ * 5(4X, F6.3) ) R2MIN = SORTR2(1) XXX = ( 1. + R2MIN ) /2. IF (XXX .GT. 0.91) XXX = R2MIN * 1.10 XXXR2 = XXX - R2MIN DO 950 I = 1, MPK SEDEL2(I,2) = ( XXX - SEDEL2(I,1) ) / XXXR2 IF (SEDEL2(I,2) .LT. -9.99) SEDEL2(I,2) = -9.99 IF (SEDEL2(I,2) .LE. 0.) GOTO 950 IF (SEDEL2(I,1) .GT. 2. * R2MIN) SEDEL2(I,2) = -9.99 950 CONTINUE IF (MPK .LE. 9) RETURN R2MAX = SORTR2(9) DO 953 I = 1, MPK IF (SEDEL2(I,1) .GT. R2MAX) SEDEL2(I,1) = 999. + SEDEL2(I,1) 953 CONTINUE RETURN END SUBROUTINE SHIFTR (SEDEL2) DIMENSION SEDEL2(50,2) 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 LOGICAL SWPRI EQUIVALENCE (SWPRI, SWITCH(10)) EQUIVALENCE (IATOMS, IFILE(2)), (IATOLD, IFILE(10)) EQUIVALENCE (IPR1, IFILE(6)), (LIS1, IFILE(7)), (LIS2, IFILE(8)) PARAMETER (KUSER2=30000) PARAMETER (MAXAT=993) COMMON / / XX(KUSER2), ATXYZ(10,MAXAT), NATX DIMENSION BLACOM(42000) EQUIVALENCE (BLACOM(1), XX(1)) COMMON /ATNAMA/ ATNAME(MAXAT) CHARACTER * 6 ATNAME PARAMETER (MXP=50, MXP1=MXP + 1) COMMON /PIEK/ RHO(MXP1), SH(3,MXP1), IBOTS(MXP1), BOTS(MXP1), NPK DIMENSION XLOCK(3) DATA DMAXTR / 0.431 / IF (SWPRI) WRITE (LIS2, FMT = ' ('' The new coordinates are: '' / * '' Atom name x y z'' /)') DO 123 N = 1,NATX ATXYZ(8,N) = ATXYZ(1,N) ATXYZ(9,N) = ATXYZ(2,N) 123 ATXYZ(10,N) = ATXYZ(3,N) REWIND IATOMS NOTR2 = 0 DO 540 I = 1, NPK INEW = I - NOTR2 DO 330 N = 1,NATX ATXYZ(1,N) = ATXYZ(8,N) + SH(1,I) ATXYZ(2,N) = ATXYZ(9,N) + SH(2,I) ATXYZ(3,N) = ATXYZ(10,N) + SH(3,I) CALL LOCKIN (ATXYZ(1,N), DMAXTR, XLOCK, DISTTR, NPOSTR) IF (NPOSTR .GT. 1) THEN IF (I .EQ. 1) WRITE (LIS1, 226) N, ATNAME(N) 226 FORMAT (' Shift nr 1, Atom', I4,' = ', A6, * ' lies on a special position') CALL KERNAB (XLOCK, ATXYZ(1,N), 3) IF (I .EQ. 1 .AND. DISTTR .GT. 0.2) WRITE (LIS1, FMT =' * ('' Atom locked-in: dist ='', F5.2, '' ang.'')') DISTTR IF (DISTTR .GT. 0.2) WRITE (LIS2, FMT =' * ('' Atom locked-in: dist ='', F5.2, '' ang.'')') DISTTR ENDIF IF (SWPRI) WRITE (LIS2, FMT = ' (4X, A6, 3F10.5)' ) * ATNAME(N), (ATXYZ(J,N), J=1,3) 330 CONTINUE R299 = SEDEL2(I,1) IF (R299 .GT. 999.) R299 = R299 - 999. IF (SEDEL2(I,2) .LE. 0.) THEN NOTR2 = NOTR2 + 1 WRITE (LIS1, 332) I, R299, SEDEL2(I,2) GOTO 540 ELSEIF (SEDEL2(I,1) .LT. 99.) THEN WRITE (LIS1, 331) I, R299, SEDEL2(I,2), INEW 331 FORMAT ( ' TRACOR shift #',I2, ' R2 =', F6.3, * ' Q2=',F5.2, ' accepted: TR=', I2) ELSEIF (I .LT. 9) THEN WRITE (LIS1, 1332) I, R299, SEDEL2(I,2), INEW 1332 FORMAT ( ' TRACOR shift #',I2, ' R2 =', F6.3, * ' Q2=',F5.2, ' accepted: TR=',I2, ' R2=BIG !') ELSE NOTR2 = NOTR2 + 1 WRITE (LIS1, 332) I, R299, SEDEL2(I,2) 332 FORMAT ( ' TRACOR shift #',I2, ' R2 =', F6.3, * ' Q2=',F5.2, ' rejected R2 is too BIG') GOTO 540 ENDIF WRITE (CHOUT, 333) INEW, (SH(J,I),J=1,3), * SEDEL2(I,1), RHO(I), IBOTS(I) 333 FORMAT ('TR=', I3, ' Shift:', 3F8.4, ' R2=', F6.3, * ' FOM=', F6.0, ' X= ', I2) IFOM = NINT(RHO(I)) IF (SEDEL2(I,1) .GT. 99.) SEDEL2(I,1) = SEDEL2(I,1) - 999. CALL ATOMWT (IATOMS, ATXYZ, ATNAME, NATX, INEW, * SEDEL2(I,1), IFOM, IBOTS(I)) 540 CONTINUE INEW = NPK - NOTR2 CALL FILCLO (IATOMS, 'KEEP') WRITE (IPR1, 544) NPK 544 FORMAT (' Number of shift vectors found by TRACOR :', I3) IF (INEW .LT. NPK) WRITE (IPR1, 545) INEW 545 FORMAT (' Number of vectors accepted after R2 test:', I3) IF (IBOTS(1) .GT. 0) THEN WRITE (IPR1, 550) WRITE (LIS1, 550) 550 FORMAT(/' NOTE: the first shift vector brings the fragment in' * /' collision with a symmetry related molecule!'/) MSAME = NATX / 3 IF (IBOTS(1) .GT. MSAME) THEN WRITE (IPR1, 554) WRITE (LIS1, 554) 554 FORMAT(' The fragment is probably put on a symmetry element:' * /' check if its internal symmetry is correct.'/) ENDIF ENDIF IF (INEW .GT. 1) WRITE (LIS1, 120) 120 FORMAT (/' All shift vectors have been applied and tested. '/ * ' All accepted parameter sets are written to the ATOMS file,' / * ' and transferred to TRAVEC'/ * ' but TRAVEC is going to reorder the atom sets.'/) RETURN END SUBROUTINE ATOMWT (IATOMS, ATXYZ,ATNAME,NAT, NNNN, R2,IFOM,IBOTS) DIMENSION ATXYZ(10,NAT), ATNAME(NAT) CHARACTER *6 ATNAME 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)) CHARACTER *6 CBOTS, TBOTS IF (NNNN .EQ. 1) WRITE (LIS2, FMT='(/ * '' Output ATOMS header lines''/)') IF (IBOTS .GT. 0) THEN CALL KERI2C (IBOTS, TBOTS, 2) CBOTS(5:6) = TBOTS(1:2) CBOTS(1:4) = ' X= ' ELSE CBOTS = ' ' ENDIF IT = MAX0 (0, 10000 * (ITIME(1)-1900) + 100 * ITIME(2) + ITIME(3)) WRITE (IATOMS, 102) CCODE, IT, KEYS(13), NNNN, R2, IFOM, CBOTS 102 FORMAT ('ATOMS ', A6, ' < TRACOR :' * , I7, ' RUN', I4, ' TR=', I3, ' R2= ', F5.3, ' FOM=', I5, A6) WRITE (LIS2, 1102) CCODE, IT, KEYS(13), NNNN, R2, IFOM, CBOTS 1102 FORMAT (/' ', 'ATOMS ', A6, ' < TRACOR :' * , I7, ' RUN', I4, ' TR=', I3, ' R2= ', F5.3, ' FOM=', I5, A6) IF (CHOUT .NE. ' ') THEN WRITE (IATOMS, FMT = '(''REMARK '', A65)') CHOUT(1:65) WRITE (LIS2, FMT = '('' REMARK '', A65)') CHOUT(1:65) CHOUT = ' ' ENDIF DO 109 I = 1, NAT WRITE (IATOMS, 104) ATNAME(I), (ATXYZ(J,I), J=1,3) 104 FORMAT ('ATOM ', A6, 2X, 3F9.5) 109 CONTINUE WRITE (IATOMS, FMT = '(''END'')') RETURN END SUBROUTINE D1FOUR (IO,NGR1,SCALE) 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)) EQUIVALENCE (ITS15, IFILE(15)) 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 (MNGR1=200, MNGR2=163, NPKM1=50, NPKM2=100, LC24=24) COMMON / / SI(1250), RHOSUM(MNGR1), JUNK(MNGR2*MNGR2-MNGR1), + RHOMAX(NPKM2), IC1(MNGR1), IC2(MNGR1), IMAX(NPKM2), + JMAX(NPKM2), XMAX(NPKM2), YMAX(NPKM2), RHOP(NPKM2), + IRHO(LC24) DIMENSION BLACOM(42000) EQUIVALENCE (BLACOM(1), SI(1)) PARAMETER (MXP=50, MXP1=MXP + 1) COMMON /PIEK/ RHO(MXP1), TX(3,MXP1), IBOTS(MXP1), BOTS(MXP1), NPK DIMENSION CO(1000) EQUIVALENCE (SI(251),CO(1)) CALL KERNZA (0.,RHOSUM,MNGR1) PA = 1000. / FLOAT(NGR1) DO 40 I1=1,NGR1 40 IC1(I1) = FLOAT(I1) * PA + 0.5 60 READ(ITS15) IH,FA,FB IF (IH.EQ.999) GOTO 70 DO 50 I1=1,NGR1 II = MOD (IH * IC1(I1), 1000) IF (II .LE. 0) II = II + 1000 RHOSUM(I1)=RHOSUM(I1)+FA*CO(II)+FB*SI(II) 50 CONTINUE GOTO 60 70 CALL KERNZA (0.,RHOMAX,NPKM1) CALL KERNZI (0,IMAX,NPKM1) SCALE=0.0 DO 100 I1=1,NGR1 SCALE=SCALE+RHOSUM(I1)**2 DO 110 I3=1,NPKM1 IF (RHOSUM(I1).GT.RHOMAX(I3)) GOTO 130 110 CONTINUE GOTO 100 130 DO 120 I4=1,NPKM1 IF ((NPKM1+1-I4).EQ.I3) GOTO 210 IMAX(NPKM1+1-I4)=IMAX(NPKM1-I4) RHOMAX(NPKM1+1-I4)=RHOMAX(NPKM1-I4) 120 CONTINUE 210 IMAX(I3)=I1 RHOMAX(I3)=RHOSUM(I1) 100 CONTINUE DO 170 I1=1,NPKM1-1 IF(RHOMAX(I1).LE. 0.01) GOTO 170 I4=I1+1 DO 180 I2=I4,NPKM1 IF (RHOMAX(I2).LE. .01) GOTO 180 IDIF1=IABS(IMAX(I1)-IMAX(I2)) IDIF1=MIN0(IDIF1,NGR1-IDIF1) IF (IDIF1.LT.3 ) GOTO 200 GOTO 180 200 RHOMAX(I2)=0. 180 CONTINUE 170 CONTINUE SCALE = 100. * SQRT (FLOAT(NGR1) / SCALE) IF (IO.LE.0) GOTO 330 WRITE (LIS2,1003) 1003 FORMAT (' ***1D.-FOURIER MAP*** '/) ITW2=1 ITW =1 WRITE (LIS2,1008) 1008 FORMAT(' U= 1 , 2 , 3....ETC'/) GOTO 310 340 ITW2=ITW2+1 ITW=LC24*(ITW2-1)+1 310 NPAG=LC24*ITW2 IF (NGR1.LT.NPAG) NPAG=NGR1 I3=0 DO 1010 I1=ITW,NPAG I3=I3+1 IRHO(I3)=IFIX(RHOSUM(I1)*SCALE+.5) 1010 CONTINUE WRITE (LIS2,1000) (IRHO(I1),I1=1,I3) 1000 FORMAT (1X,24I5) IF (NGR1.GT.NPAG) GOTO 340 330 N=0 DO 250 I1=1,NPKM1 IF (RHOMAX(I1).LE. .01) GOTO 250 IX1=IMAX(I1)-1 IX2=IMAX(I1) IX3=IMAX(I1)+1 IXX3=IX3 IXX1=IX1 IF (IX2.EQ.NGR1) IXX3=IX3-NGR1 IF (IX2.EQ.1) IXX1=IX1+NGR1 RHO1=RHOSUM(IXX1) RHO2=RHOMAX(I1) RHO3=RHOSUM(IXX3) IF (RHO1.GT.RHO2 .OR. RHO3.GT.RHO2) GOTO 250 CALL INTPOL(RHO1,RHO2,RHO3,IX1,IX3,X,NGR1) XMAX(I1)=AMOD(X,1.0) RHOP(I1)=RHO2*SCALE N=N+1 RHOP(N)=RHOP(I1) XMAX(N)=XMAX(I1) IF (N .EQ. MXP) GOTO 260 250 CONTINUE 260 RHOP(N+1)=-9999. XMAX(N+1) = 0. NPK = N RETURN END SUBROUTINE D2FOUR (IO, NGR1, NGR2, SCALE) 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)) EQUIVALENCE (ITS15, IFILE(15)) 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 (MNGR1=200, MNGR2=163, NPKM2=100, LC24=24) COMMON / / SI(1250),RHOSUM(MNGR2,MNGR2),RHOMAX(NPKM2),IC1(MNGR1), + IC2(MNGR1), IMAX(NPKM2), JMAX(NPKM2), XMAX(NPKM2), + YMAX(NPKM2),RHOP(NPKM2),IRHO(LC24) DIMENSION BLACOM(42000) EQUIVALENCE (BLACOM(1), SI(1)) PARAMETER (MXP=50, MXP1=MXP + 1) COMMON /PIEK/ RHO(MXP1), TX(3,MXP1), IBOTS(MXP1), BOTS(MXP1), NPK DIMENSION CO(1000),AH(MNGR1),BH(MNGR1) EQUIVALENCE (SI(251),CO(1)),(AH(1),IMAX(1)),(BH(1),XMAX(1)) CALL KERNZA (0.,RHOSUM,MNGR2*MNGR2) PA1 = 1000. / FLOAT(NGR1) PA2 = 1000. / FLOAT(NGR2) DO 40 I1=1,NGR1 40 IC1(I1) = FLOAT(I1) * PA1 + 0.5 DO 42 I2=1,NGR2 42 IC2(I2) = FLOAT(I2) * PA2 + 0.5 IHO=-1 48 READ(ITS15) IH,IK,FA,FB IF (IHO.EQ.-1) GOTO 53 IF (IH.EQ.IHO) GOTO 54 DO 52 I1=1,NGR1 IHX = MOD (IHO * IC1(I1), 1000) IF (IHX .LE. 0) IHX = IHX + 1000 DO 52 I2=1,NGR2 52 RHOSUM(I1,I2)=RHOSUM(I1,I2)+AH(I2)*CO(IHX)+BH(I2)*SI(IHX) IF (IH.EQ.999) GOTO 70 53 CALL KERNZA(0.,AH,MNGR1) CALL KERNZA(0.,BH,MNGR1) IHO=IH 54 DO 55 I2=1,NGR2 KY = MOD (IK * IC2(I2), 1000) IF (KY .LE. 0) KY = KY + 1000 AH(I2)=AH(I2) + FA*CO(KY) + FB*SI(KY) 55 BH(I2)=BH(I2) + FB*CO(KY) - FA*SI(KY) GOTO 48 70 CONTINUE CALL KERNZA (0.,RHOMAX,NPKM2) CALL KERNZI (0,IMAX,NPKM2) CALL KERNZI (0,JMAX,NPKM2) SCALE=0.0 DO 100 I1=1,NGR1 DO 100 I2=1,NGR2 SCALE=SCALE+RHOSUM(I1,I2)**2 IF (RHOSUM(I1,I2).LE.RHOMAX(NPKM2)) GOTO 100 DO 110 I3=1,NPKM2 IF (RHOSUM(I1,I2).GT.RHOMAX(I3)) GOTO 130 110 CONTINUE GOTO 100 130 DO 120 I4=1,NPKM2 IF ((NPKM2+1-I4).EQ.I3) GOTO 210 IMAX(NPKM2+1-I4)=IMAX(NPKM2-I4) JMAX(NPKM2+1-I4)=JMAX(NPKM2-I4) RHOMAX(NPKM2+1-I4)=RHOMAX(NPKM2-I4) 120 CONTINUE 210 IMAX(I3)=I1 JMAX(I3)=I2 RHOMAX(I3)=RHOSUM(I1,I2) 100 CONTINUE DO 170 I1=1,NPKM2-1 IF (RHOMAX(I1) .LE. .01) GOTO 170 I4=I1+1 DO 180 I2=I4,NPKM2 IF (RHOMAX(I2) .LE. .01) GOTO 180 IDIF1=IABS(IMAX(I1)-IMAX(I2)) IDIF1=MIN0(IDIF1,NGR1-IDIF1) IDIF2=IABS(JMAX(I1)-JMAX(I2)) IDIF2=MIN0(IDIF2,NGR2-IDIF2) IF (IDIF1.GE.3 .OR. IDIF2.GE.3) GOTO 180 RHOMAX(I2)=0. 180 CONTINUE 170 CONTINUE SCALE = 100. * SQRT (FLOAT(NGR1 * NGR2) / SCALE) IF (IO.LE.0) GOTO 330 ITW2=1 ITW=1 GOTO 310 340 ITW2=ITW2+1 ITW=LC24*(ITW2-1)+1 310 NPAG=LC24*ITW2 IF (NGR1.LT.NPAG) NPAG=NGR1 WRITE (LIS2,1012) (I,I=ITW,NPAG) 1012 FORMAT ('1',3X,'U=',24(2X,I3)) WRITE (LIS2,1013) 1013 FORMAT(' ','V=') DO 300 I1=1,NGR2 I3=0 DO 302 I2=ITW,NPAG I3=I3+1 IRHO(I3)=IFIX(RHOSUM(I2,I1)*SCALE+.5) 302 CONTINUE WRITE (LIS2,1000) I1,(IRHO(I2),I2=1,I3) 1000 FORMAT ('0',I3,2X,24I5) 300 CONTINUE IF (NGR1.GT.NPAG) GOTO 340 330 N=0 DO 250 I1=1,NPKM2 IF (RHOMAX(I1).LE. .01) GOTO 250 IX1=IMAX(I1)-1 IX2=IMAX(I1) IX3=IMAX(I1)+1 IY1=JMAX(I1)-1 IY2=JMAX(I1) IY3=JMAX(I1)+1 IXX3=IX3 IXX1=IX1 IF (IX2.EQ.NGR1) IXX3=IX3-NGR1 IF (IX2.EQ.1) IXX1=IX1+NGR1 RHO1=RHOSUM(IXX1,IY2) RHO2=RHOMAX(I1) RHO3=RHOSUM(IXX3,IY2) IF (RHO1.GT.RHO2 .OR. RHO3.GT.RHO2) GOTO 250 CALL INTPOL(RHO1,RHO2,RHO3,IX1,IX3,X,NGR1) IM1=IMAX(I1) XMAX(I1)=AMOD(X,1.0) RHOP(I1)=RHO2 IYY1=IY1 IYY3=IY3 IF (IY2.EQ.NGR2) IYY3=IY3-NGR2 IF (IY2.EQ.1) IYY1=IY1+NGR2 RHO1=RHOSUM(IX2,IYY1) RHO2=RHOMAX(I1) RHO3=RHOSUM(IX2,IYY3) IF (RHO1.GT.RHO2 .OR. RHO3.GT.RHO2) GOTO 250 CALL INTPOL (RHO1,RHO2,RHO3,IY1,IY3,X,NGR2) YMAX(I1)=AMOD(X,1.0) RHOP(I1)=(RHOP(I1)+RHO2)/2.*SCALE N=N+1 RHOP(N)=RHOP(I1) XMAX(N)=XMAX(I1) YMAX(N)=YMAX(I1) IF (N .EQ. MXP) GOTO 260 250 CONTINUE 260 RHOP(N+1)=-999. XMAX(N+1) = 0. YMAX(N+1) = 0. NPK = N IF (NPK .LE. 1) RETURN IF (ILAUE .NE. 8) RETURN RHM = 0.45 * RHOP(1) NPK = 0 DO 550 I1 = 1, N IF (RHOP(I1) .LT. RHM) GOTO 550 NPK = NPK + 1 RHOP(NPK) = RHOP(I1) XMAX(NPK) = XMAX(I1) YMAX(NPK) = YMAX(I1) IF (I1 .EQ. N) GOTO 550 RHM1 = 0.83 * RHOP(I1) DO 540 I2 = I1+1, N IF (RHOP(I2) .LT. RHM1) GOTO 540 XM1 = XMAX(I1) - XMAX(I2) IF (XM1 .LT. 0.0) XM1 = XM1 + 1.0 YM1 = YMAX(I1) - YMAX(I2) IF (YM1 .LT. 0.0) YM1 = YM1 + 1.0 IF (ABS(XM1-0.333).LE.0.02 .AND. ABS(YM1-0.667).LE.0.02) GOTO 530 IF (ABS(XM1-0.667).GT.0.02 .OR. ABS(YM1-0.333).GT.0.02) GOTO 540 530 RHOP(I2) = 0.0 540 CONTINUE 550 CONTINUE RETURN END SUBROUTINE INTPOL(RHO1,RHO2,RHO3,IX1,IX3,X,IGRID) S=1. IF (RHO1.LT.RHO3) GOTO 10 RHO33=RHO1-RHO3 RHO2=RHO2-RHO3 RHO11=0. IX1=IX3 S=-1. GOTO 20 10 RHO33=RHO3-RHO1 RHO2=RHO2-RHO1 RHO11=0. 20 TEL=(RHO33/RHO2)-4. A=(.5*RHO33)-RHO2 B=RHO2-A DENOM=(2.*RHO33/RHO2)-4. X=TEL/DENOM RHOM = (A*X+B)*X IF (NINT(S) .EQ. -1) GOTO 40 RHOM=RHOM+RHO1 GOTO 50 40 RHOM=RHOM+RHO3 50 X=(FLOAT(IX1)+S*X)/FLOAT(IGRID) RHO2=RHOM RETURN END SUBROUTINE D3FOUR 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 (ICRYS, IFILE(3)), (ICON, IFILE(4)) EQUIVALENCE (IBINBI, IFILE(18)) 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 (KUSER2=30000, KUSER1=KUSER2/3) PARAMETER (MXP=50, MXP1=MXP + 1) COMMON /PIEK/ RHO(MXP1), X(3,MXP1), IBOTS(MXP1), BOTS(MXP1), NPK COMMON /COTRA/ MULS(48,48), INDX(3,3), IDIMF, TT(3,48), * EFST(48), FHST(48), PHST(48), CELLN(6), MAXHKL(3), * P1SQ, PSQ, EMIN, EO2AV, SCX, BOVX, DAMPX, SMM, IOMAPS, SMAX, * IRR(3,3,48), IORG, BBB, D2R, R2D, * AMULT, ASYMM, ALATT, ASYMCL, NSYMC, ASYMC COMMON /FFTDA/ SCALE, MH(3), NPP(3), XLMIN(3), XLMAX(3) CALL KERNAI (MAXHKL, MH, 3) CALL KERNAB (CELLN, CELL, 6) NSYMM = 1 IF (ILAUE .EQ. 5) THEN NLATT = 2 IMULT = 2 TLATT(1,2) = 0.5 TLATT(2,2) = 0.5 TLATT(3,2) = 0. ELSE NLATT = 1 IMULT = 1 ENDIF ICENT = 1 ILAUE = 1 ISYST = 1 ILATT = 1 CALL RCELLR (CELL, VOLUM, RCELL) CALL CELLRR (CELL, RRMAT) CALL MATF2C (CELL, FRAC2C) KEYS(27) = 5 CALL FFTRIN (KUSER1) CALL PP1 CALL SEART CALL RDCRYS (ICRYS) RETURN END SUBROUTINE FFTRIN (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)), (ICON, IFILE(4)) EQUIVALENCE (LIS1, IFILE(7)), (LIS2, IFILE(8)) EQUIVALENCE (IBINFF, IFILE(16)), (IFMAP, IFILE(17)) EQUIVALENCE (ISCRA, IFILE(18)) EQUIVALENCE (KEYS(27), IMAP), (KEYS(28), IHALF) LOGICAL SWPRI, PRIMAP 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 /SEARTA/ D2R, DMPIC, DMAXB, DMOUT, DMINB, ANGM(2), MCON, * SEARDX, NPIC, NATIN, NAT, NATS, NRECY, PSQ, NATREC, KPROG PARAMETER (MAXBUF = 198) DIMENSION FITFFT(5), BUFFFT(MAXBUF), IGM3(3) DIMENSION MAXHKL(3) FACTOR = 0.3 CALL KERNZA (0.0, XLMIN, 3) CALL KERNZA (1.0, XLMAX, 3) IGM3(1) = 1 IGM3(2) = 1 IGM3(3) = 2 CALL BINIFF (1, IBINFF, 'BINFFT', FITFFT, NITFFT, BUFFFT, NEND) PSQ = BUFFFT(27) CALL FILINQ (ISCRA, 'BINBIG', 'UNFORMATTED', 'OUTPUT', KINQ) CALL FILINQ (IFMAP, 'FMAPT', 'UNFORMATTED', 'OUTPUT', KINQ) CALL KERF2I (BUFFFT(22), MAXHKL, 3) IF (MAXHKL(1) .EQ. 0) CALL KERF2I (BUFFFT(7), MAXHKL, 3) 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 IHALF = 0 IHALF = 1 IF (ILATT .EQ. 1 .OR. ILATT .EQ. 3) IHALF = 0 FACTOR = 999. DO 219 I=1,3 NPP(I) = 2 * MH(I) + 2 219 FACTOR = AMIN1 (FACTOR, 0.90 * CELL(I)/NPP(I)) IF (IHALF.NE.0 .AND. IGM3(3).EQ.1) IGM3(2) = 2 GOTO 242 240 DO 241 I=1,3 241 NPP(I) = CELL(I) / FACTOR + 0.5 242 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, 282) FACTOR 282 FORMAT (' The GRID spacing is approximately', F6.3, ' Angstrom') 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 SEART 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 (ICON, IFILE( 4)) EQUIVALENCE (LIS1, IFILE( 7)), (LIS2, IFILE(8)) EQUIVALENCE (IFMAP, IFILE(17)) LOGICAL SWPRI EQUIVALENCE (SWPRI, SWITCH(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 (KUSER2=30000) PARAMETER (MAXAT=993) COMMON / / NSTORE(KUSER2), X(4,MAXAT) DIMENSION NBLACO(42000) EQUIVALENCE (NBLACO(1), NSTORE(1)) PARAMETER (MXP=50, MXP1=MXP + 1) COMMON /PIEK/ RHO(MXP1), XT(3,MXP1), IBOTS(MXP1), BOTS(MXP1), NPK COMMON /SEARTA/ D2R, DMPIC, DMAXB, DMOUT, DMINB, ANGM(2), MCON, * SEARDX, NPIC, NATIN, NAT, NATS, NRECY, PSQ, NATREC, KPROG DIMENSION ITLE(20) D2R = ATAN(1.0) / 45.0 DMPIC = 0.5 REWIND IFMAP READ (IFMAP) ITLE,IMAP,IHALF WRITE (LIS1, 192) 192 FORMAT(/' FFT PEAK SEARCH FOR TRACOR === VERSION JUN 1987') NPIC = MXP NAT = MXP CALL PKSRT (X, MAXAT, IHALF, IFMAP) CALL FILCLO (IFMAP, 'DELETE') IF (SWPRI) WRITE (LIS2, 445) 445 FORMAT ('0Untransformed peaks from 3D Fourier map'/ * ' PEAK',7X,'HEIGHT',13X,'X',9X,'Y',9X,'Z') RHOMIN = 0.50 * X(4,1) DO 520 I=1,NPIC IF (X(4,I) .LT. RHOMIN) GOTO 530 DO 515 J=1,3 515 XT(J,I)=X(J,I) RHO(I)=X(4,I) IF (SWPRI) WRITE (LIS2, 517) I,RHO(I),(XT(J,I),J=1,3) 517 FORMAT (' ' ,3X,I2,7X,F6.0,10X,F7.4,3X,F7.4,3X,F7.4) IF (I .EQ. 50) GOTO 550 520 CONTINUE I = NPIC + 1 530 I = I - 1 550 DO 555 J=1,3 555 XT(J,I+1)=0. RHO(I+1)=-9999. NPK = I RETURN END SUBROUTINE DISCHK (NA, ISAME, MXP2) DIMENSION ISAME(MXP2) 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(2)) 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 (KUSER2=30000, KUSER1=KUSER2/3) PARAMETER (MAXAT=993) COMMON / / XX(KUSER2), ATXYZ(10,MAXAT), NATX DIMENSION IZAT(MAXAT), WORK(3,KUSER1) EQUIVALENCE (WORK(1,1), IZAT(1), XX(1)) DIMENSION BLACOM(42000) EQUIVALENCE (BLACOM(1), XX(1)) COMMON /ATNAMA/ ATNAME(MAXAT) CHARACTER * 6 ATNAME PARAMETER (MXP=50, MXP1=MXP + 1) COMMON /PIEK/ RHO(MXP1), SH(3,MXP1), IBOTS(MXP1), BOTS(MXP1), NPK COMMON /COTRA/ MULS(48,48), INDX(3,3), IDIMF, TT(3,48), * EFST(48), FHST(48), PHST(48), CELLN(6), MAXHKL(3), * P1SQ, PSQ, EMIN, EO2AV, SCX, BOVX, DAMPX, SMM, IOMAPS, SMAX, * IRR(3,3,48), IORG, BBB, D2R, R2D, * AMULT, ASYMM, ALATT, ASYMCL, NSYMC, ASYMC DIMENSION TL(3,4) DIMENSION COSA(3), SINA(3), Q(3,3), SSH(3), XY(3), XGRAV(3) DIMENSION RATIO(MXP) CALL KERNZA(0.0, TL, 12) GOTO (37,31,31,31,31,34,33),ILATT 31 DO 32 I=1,3 32 TL(I,2)=0.5 IF(ILATT.LT.5) TL(ILATT-1,2)=0.0 GOTO 37 33 TL(1,2)=1.0/3.0 TL(2,2)=2*TL(1,2) TL(3,2)=TL(2,2) TL(1,3)=TL(2,2) TL(2,3)=TL(1,2) TL(3,3)=TL(1,2) GOTO 37 34 DO 36 I=2,4 DO 35 J=1,3 35 TL(J,I)=0.5 36 TL(I-1,I)=0.0 37 CONTINUE DMAX=3.0 DMIN = 2.4 DO 110 I=1,3 COSA(I) = COS(CELL(I+3) * D2R) 110 SINA(I) = SQRT(1.0- COSA(I)**2 ) COSGS=(COSA(1)*COSA(2)-COSA(3)) / SINA(1)/SINA(2) Q(1,1)=CELL(1)*SINA(2)*SQRT(1.0-COSGS**2) Q(1,2)=0.0 Q(1,3)=0.0 Q(2,1)=-CELL(1)*COSGS*SINA(2) Q(2,2)=CELL(2)*SINA(1) Q(2,3)=0.0 Q(3,1)=CELL(1)*COSA(2) Q(3,2)=CELL(2)*COSA(1) Q(3,3)=CELL(3) CALL FILINQ (IATOMS, 'ATOMS', 'FORMATTED', 'INPUT', KINQ) IF (KINQ.EQ.-1) CALL KERROR (' No ATOMS file found', 0, 'DISCHK') CALL ATOMIN (IATOMS, ATXYZ, ATNAME, IZAT, MAXAT, NA, KEYT) NATX = NA CALL KERNZA(0.0, XGRAV, 3) DO 130 N=1,NA DO 130 I=1,3 130 XGRAV(I) = XGRAV(I) + ATXYZ(I,N) DO 131 I=1,3 131 XGRAV(I) = XGRAV(I) / FLOAT(NA) DO 140 IS=1,NPK DO 140 I=1,3 IF (XGRAV(I) + SH(I,IS) .GT. 1.0) SH(I,IS) = SH(I,IS) - 1.0 140 CONTINUE WRITE (LIS2, 149) 149 FORMAT (/' Shifts, with distance check'/ ' Peak height ', 3X, * 'Tx',6X,'Ty',6X,'Tz', 8X, 'contacts by symmetry: Angstrom') WRITE (LIS2, 150) DMAX 150 FORMAT (63X, 'max: ', F4.2) K=0 DO 168 ISS=2,NSYMC DO 167 IL=1,NLATT DO 166 N=1,NA K=K+1 DO 165 I=1,3 WORK(I,K)=TL(I,IL)+TT(I,ISS) DO 164 J=1,3 164 WORK(I,K)=WORK(I,K)+ATXYZ(J,N)*IRR(I,J,ISS) 165 WORK(I,K)=AMOD(4.0+WORK(I,K),1.0) 166 CONTINUE 167 CONTINUE 168 CONTINUE RHOMIN = 0.5 * RHO(1) NPRINT = 0 ISAVE = 0 NBADX = 0 DO 200 IS=1,NPK IBOTS(IS)=0 BOTS(IS)=999. ISAME(IS) = 0 IF (RHO(IS) .LT. RHOMIN) THEN IF (ISAVE .EQ. 0) THEN ISAVE = IS - 1 WRITE (LIS2, FMT='(/ * '' The following shift(s) are rejected:'' /)') ENDIF WRITE (LIS2,2) IS,RHO(IS),(SH(J,IS),J=1,3) RHO(IS) = - 1.0 RHOMIN = RHO(1) + 1. GOTO 200 ENDIF RATIO(IS) = 1. IF (IS .GT. 1) THEN IF (RHO(IS) .GT. RHO(IS-1)) RHO(IS) = RHO(IS-1) RATIO(IS) = RHO(IS-1) / RHO(IS) IF (RATIO(IS) .GT. 1.25) RHOMIN = RHO(IS) * 0.99 ENDIF RATI2 = 0.0 IF (IS .GT. 2) RATI2 = RHO(IS-2) / RHO(IS) IF (RATI2 .GT. 1.30) RHOMIN = RHO(IS) * 0.99 RATIO(1) = AMAX1 (RATIO(1), RATIO(IS)) IF (IS .GE. 5 ) THEN TQ = RATIO(IS) + RATIO(IS-1) + RATIO(IS-2) - 3.0 IF (RATIO(1)-1.0 .GT. 4.0*TQ .AND. RHO(IS) .LT. 0.7*RHO(1) * .AND. IS .LE. 10) RHOMIN = RHO(IS) * 0.99 ENDIF NPRIN = 0 DMAX=3.0 DMAX2 = DMAX * DMAX NBAD = 0 WRITE (LIS2,2) IS,RHO(IS),(SH(J,IS),J=1,3) 2 FORMAT(I4, F7.0, 3F8.4) DO 197 N=1,NA DO 170 I=1,3 170 XY(I)=ATXYZ(I,N)+SH(I,IS) IFI=0 DO 188 ISS=2,NSYMC DO 175 I=1,3 SSH(I)=0.0 DO 172 J=1,3 172 SSH(I)=SSH(I)+SH(J,IS)*IRR(I,J,ISS) 175 SSH(I)=AMOD(4.5+SSH(I),1.0)-0.5 DO 186 IL=1,NLATT IST=IFI+N IFI=IFI+NA DO 184 J=IST,IFI D1=WORK(1,J)+SSH(1)-XY(1) D1=(AMOD(4.5+D1,1.0)-0.5) D4=ABS(D1)*Q(1,1) IF (D4.GT.DMAX) GOTO 184 D2=WORK(2,J)+SSH(2)-XY(2) D2=(AMOD(4.5+D2,1.0)-0.5) D5=ABS(D2*Q(2,2)+D1*Q(2,1)) IF (D5.GT.DMAX) GOTO 184 D3=WORK(3,J)+SSH(3)-XY(3) D3=(AMOD(4.5+D3,1.0)-0.5) D6=ABS(D3*Q(3,3)+D2*Q(3,2)+D1*Q(3,1)) IF (D6.GT.DMAX) GOTO 184 D = D4*D4+D5*D5+D6*D6 IF (D.GT.DMAX2) GOTO 184 D = SQRT(D) BOTS(IS) = AMIN1(BOTS(IS), D) M=MOD(J-1,NA)+1 IF ((NPRINT .GT. 50 .AND. NPRIN .EQ. 8) .OR. * (NPRINT .GT. 70 .AND. NPRIN .EQ. 12)) THEN DMAX = AMAX1 (BOTS(IS), DMIN) DMAX2 = DMAX * DMAX WRITE (LIS2, 150) DMAX IF (D .GT. DMAX) GOTO 184 ENDIF IF (NPRINT .GT. 99 .AND. NPRIN .EQ. 3) THEN DMAX = DMIN DMAX2 = DMAX * DMAX WRITE (LIS2, 150) DMAX IF (D .GT. DMAX) GOTO 184 ENDIF NPRINT = NPRINT + 1 NPRIN = NPRIN + 1 IF (NLATT.GT.1) GOTO 180 WRITE (LIS2,3) ATNAME(N),ATNAME(M),ISS,D 3 FORMAT(35X, A6, ' ---- ', A6, ' (Sym.', I2, ') =', F6.2) GOTO 182 180 WRITE (LIS2,4) ATNAME(N),ATNAME(M),ISS,IL,D 4 FORMAT(35X, A6, ' --- ', A6, ' (Sym.', 2I2, ') =', F6.2) 182 CONTINUE IF (D .LT. 0.63) ISAME(IS) = ISAME(IS) + 1 IF (D .GE. 0.43 .AND. D .LE. 0.83) GOTO 193 IF (N .NE. M) THEN IF (MIN0 (IZAT(N), IZAT(M) ) .GE. 7 .AND. D .LE. 2.2) GOTO 193 IF (D .LE. DMIN) GOTO 193 ELSE WRITE (LIS2, FMT='(63X, ''symmetry?'')') ENDIF 184 CONTINUE 186 CONTINUE 188 CONTINUE GOTO 197 193 IBOTS(IS) = IBOTS(IS) + 1 197 CONTINUE IF (IBOTS(IS) .EQ. 0) GOTO 200 WRITE (LIS2, FMT='(66X, ''? skip'')') NBAD = NBAD + 1 IF (IS .LE. 5) THEN NBADX = NBADX + 1 WRITE (LIS1, 2252) IS, RHO(IS) 2252 FORMAT (' Warning: shift nr.',I3, ' (pk =', F5.0, * ') leads to bad contacts') IF (ISAME(IS) .GT. 2) WRITE (LIS1, 2253) IS, ISAME(IS) 2253 FORMAT (' Note for shift nr.',I3, ':', I4, * ' atom pairs coincide by symmetry !') ELSE IF (NBADX.GT.0) WRITE (LIS1,FMT='(9X,''Also for others ..'')') NBADX = -99 ENDIF 200 CONTINUE IF (ISAVE .GT. 0) NPK = ISAVE RETURN END SUBROUTINE PKSRT (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 /SEARTA/ D2R, DMPIC, DMAXB, DMOUT, DMINB, ANGM(2), MCON, * SEARDX, NPIC, NATIN, NAT, NATS, NRECY, PSQ, NATREC, KPROG PARAMETER (KUSER2=30000) COMMON / / NR3D INTEGER*2 NR3D(KUSER2) INTEGER*2 ILACOM(84000) EQUIVALENCE (ILACOM(1), NR3D(1)) DIMENSION XS(3), X1(3), IDIFF(19), B(19) DIMENSION DXYZM(3) 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, 'PKSRT') DX = 1.0 / FLOAT(NNX) DY = 1.0 / FLOAT(NNYT) DZ = 1.0 / FLOAT(NNZ) LEVEL = 0 LIMIT = MIN0(MAXAT, 2*NAT) 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) READ(LIN) 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) B(10) = AMAX1(B(10), BINT) 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 = NPIC LEVEL = X(4,NPIC) + 0.5 1980 CONTINUE 2000 CONTINUE IF (IY .GE. NNY) GO TO 2100 IF (IY - NNYOLD + 2) 1400, 1200, 1400 2100 CALL SORT (X, MAXAT, NO, 4) NNN = MIN0 (NO, NPIC) IF (NNN .EQ. NPIC) RETURN LEVEL = LEVEL - 100 IF(LEVEL.GE.(-200)) GO TO 1100 NPIC=NO 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 (' Appr. number of 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.LT.3.OR.IMAP.GT.4) 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 ILACOM(84000) EQUIVALENCE (ILACOM(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