C FO'S ABSORPTION CORRECTION FOLLOWING THE METHOD OF C WALKER -STUART. C Original source code obtained from Babu Varghese, IIT/M, C (author ?) 1991 C----------------------- modified mvr/sept.,'97------------------ C C Input Files: C CR : Orientation matrix(reciprocal axis matrix) C Diffraction geometry (Philips option 1) C etc. Sample file at the end of this code C LIST3.DAT : SHELXL93 LIST 3 output (option 1) C or C LIST4.DAT : SHELXL93 LIST 4 (w) output scaled by w C if necessary - see SHELXL93 manual C (option 2) C Output Files: C ABS.LST : correction details C NEWF.HKL : for option 1 : This should be given as C HKLF 3 input for SHELXL93 C or C NEWF2.HKL : for option 2 : This should be given as C HKLF 4 input for SHELXL93 C-------------------------------------------------------------------- C C MAIN PROGRAM INTEGER LP,SHX,CR,SCR,SC1,SC2,NEW,SC3,SCZ CHARACTER*72 TITLE REAL*8 UB(3,3) DIMENSION MNR1(34), MNR2(34), MN(2,34) CHARACTER*1 FLAG REAL*8 ROT(3,3),ROTINV(3,3),UBP(3,3) REAL*8 PQ(34,34),TN(34),RS(34,34),TN1(34) REAL*8 SCALA,OME,FOQ,SCMN(34),SCM(12),SMN,CMN REAL*8 SUMFO,SUMFM,SUMFC,DELFC,DELFM,RINIZ,RFINAL DATA LP,SHX,CR,SCR,SC1,SC2,SC3,NEW,SCZ/11,12,13,14,18,19,21,22,15/ DATA IH0,IK0,IL0,LCO,THRESL/0,0,0,500,1.0E-18/ C C DATA (MN(I,J),J=1,34),I=1,2)/ C 9 0,0,0,0,0,0,1,2,3,4,1,1,1,1,2,3,4, C 9 0,0,0,0,0,0,1,2,3,4,1,1,1,1,2,3,4, C 9 1,2,3,4,5,6,0,0,0,0,1,2,3,4,2,3,4, C 9 1,2,3,4,5,6,0,0,0,0,1,2,3,4,2,3,4/ C DATA MNR1/ 0,0,0,0,0,0,1,2,3,4,1,1,1,1,2,3,4, 9 0,0,0,0,0,0,1,2,3,4,1,1,1,1,2,3,4/ DATA MNR2/ 1,2,3,4,5,6,0,0,0,0,1,2,3,4,2,3,4, 9 1,2,3,4,5,6,0,0,0,0,1,2,3,4,2,3,4/ DO 501 J = 1, 34 MN(1,J) = MNR1(J) MN(2,J) = MNR2(J) 501 CONTINUE c MLSNE=0 c WRITE(*,502) MLSNE 502 FORMAT(1H+,30X,I3) C WRITE(*,701) 701 FORMAT(1X,' Input mode ? 1 for LIST3, 2 for LIST4') READ(*,*) INMOD IF (INMOD.EQ.1) THEN CALL IN3 ELSEIF (INMOD.EQ.2) THEN CALL IN4 ELSE STOP ENDIF C OPEN(CR,FILE='CR',STATUS='OLD') OPEN(SHX,FILE='SHX',STATUS='OLD') OPEN(SCR,FILE='SCR',STATUS='NEW') OPEN(LP, FILE='ABS.LST',STATUS='NEW') WRITE(LP,1000) READ(CR,1010)TITLE WRITE(LP,1020)TITLE READ(CR,*)IGEOM GO TO (10,12),IGEOM 10 WRITE(LP,1030) GO TO 13 12 WRITE(LP,1040) 13 CONTINUE READ(CR,*)WAVE WRITE(LP,1050)WAVE READ(CR,*)((UB(I,J),J=1,3),I=1,3) WRITE(LP,1052)((UB(I,J),J=1,3),I=1,3) READ(CR,*)NCICLI ICICLO=1 NUM=0 SUMFO=0. SUMFC=0. C READ SHELX LIST 5 OUTPUT AND CREATE THE FILE SCRATCH SCR NCO=0 WRITE(*,505) 505 FORMAT(1X) 35 READ(SHX,101,END=50)IH,IK,IL,FO,FC,SIGFO NCO=NCO+1 IF(MOD(NCO,LCO).EQ.0)WRITE(*,506)NCO 506 FORMAT(1H+,10X,I5) IF(FO.LT.3.0*SIGFO) GO TO 35 WRITE(SCR,1070)IH,IK,IL,FO,SIGFO,FC NUM=NUM+1 SUMFO=SUMFO+FO SUMFC=SUMFC+ABS(FC) GO TO 35 50 CONTINUE WRITE(*,506)NCO C MLSNE = MLSNE+1 C WRITE(*,502) MLSNE C EVALUATE SCALE FACTOR 2000 WRITE(LP,1200)ICICLO WRITE(*,1200) ICICLO SCALA=SUMFC/SUMFO WRITE(LP,1080)NUM,SCALA APSMIN=1.0E10 APSMAX=-1.0E10 ATMIN=1.0E10 ATMAX=-1.0E10 C BUILD THE 34 X 34 LINEAR SYSTEM IN POLAR ANGLES PHI AND MU DO 100 I=1,34 TN(I)=0. DO 100 J=1,34 100 PQ(I,J)=0. NRIF=0 OPEN(SCR,FILE='SCR') NCO=0 WRITE(*,505) 110 READ(SCR,1070,END=350)IH,IK,IL,FO,SIGFO,FC NCO=NCO+1 IF(MOD(NCO,LCO).EQ.0)WRITE(*,507)NCO 507 FORMAT(1H+,15X,I5) FCC=ABS(FC) IF(FCC.LT.3.0*SCALA*SIGFO.OR.FCC.GT.2.0*SCALA*FO) GO TO 110 NRIF=NRIF+1 CALL ANG(IGEOM,IH,IK,IL,UB,WAVE,TET,CHI,PHI) CALL PHIMU(WAVE,TET,CHI,PHI,PHIP,PHIS,U) OME=1.0/(SCALA*SIGFO) FOQ=FO*FO*OME DO 150 I=1,17 CALL SENOMN(MN(1,I),MN(2,I),PHIP,PHIS,U,SCMN(I)) CALL COSMN(MN(1,I+17),MN(2,I+17),PHIP,PHIS,U,SCMN(I+17)) 150 CONTINUE DO 200 I=1,17 DO 200 J=1,34 PQ(I,J)=PQ(I,J)+FOQ*SCMN(J)*SCMN(I) PQ(I+17,J)=PQ(I+17,J)+FOQ*SCMN(J)*SCMN(I+17) 200 CONTINUE DO 210 N=1,34 TN(N)=TN(N)+FO*OME*SCMN(N)*(ABS(FC)/SCALA-FO) 210 CONTINUE GO TO 110 350 CONTINUE WRITE(*,507)NCO WRITE(LP,1072)NRIF C SOLVE THE 34 X 34 LINEAR EQUATIONS CALL GAUSJO(PQ,TN,34) WRITE(LP,1090)(MN(1,I),MN(2,I),TN(I),I=1,34) C BUILD THE 12 X 12 LINEAR SYSTEM IN THETA DO 360 I=1,12 TN1(I)=0. DO 360 J=1,12 360 RS(I,J)=0. C MLSNE = MLSNE+1 C WRITE(*,502) MLSNE OPEN(SCR,FILE='SCR') IRIF=0 NCO=0 WRITE(*,505) 370 READ(SCR,1070,END=600)IH,IK,IL,FO,SIGFO,FC NCO=NCO+1 IF(MOD(NCO,LCO).EQ.0)WRITE(*,508)NCO 508 FORMAT(1H+,20X,I5) FCC=ABS(FC) IF(FCC.LT.3.0*SCALA*SIGFO.OR.FCC.GT.2.0*SCALA*FO) GO TO 370 IRIF=IRIF+1 CALL ANG(IGEOM,IH,IK,IL,UB,WAVE,TET,CHI,PHI) CALL PHIMU(WAVE,TET,CHI,PHI,PHIP,PHIS,U) APS=1. DO 400 I=1,17 CALL SENOMN(MN(1,I),MN(2,I),PHIP,PHIS,U,SMN) CALL COSMN(MN(1,I+17),MN(2,I+17),PHIP,PHIS,U,CMN) APS=APS+SMN*TN(I)+CMN*TN(I+17) 400 CONTINUE IF(APS.LT.0.)CALL CORR(1,IH,IK,IL,TET,CHI,PHI,FO,APS) FM=SCALA*FO*APS OME=1.0/(SCALA*SIGFO) FMM=FM*FM*OME DO 450 M=1,6 SCM(M)=SIN(M*TET) SCM(M+6)=COS(M*TET) 450 CONTINUE DO 500 I=1,6 DO 500 J=1,12 RS(I,J)=RS(I,J)+FMM*SCM(J)*SCM(I) RS(I+6,J)=RS(I+6,J)+FMM*SCM(J)*SCM(I+6) 500 CONTINUE DO 510 I=1,12 TN1(I)=TN1(I)+FM*OME*SCM(I)*(ABS(FC)-FM) 510 CONTINUE GO TO 370 600 CONTINUE WRITE(*,508)NCO WRITE(LP,1073)IRIF C SOLVE THE 12 X 12 LINEAR SYSTEM CALL GAUSJO(RS,TN1,12) WRITE(LP,1100) DO 615 I=1,12 JJJ=I IF(JJJ.GT.6)JJJ=JJJ-6 615 WRITE(LP,1110)JJJ,TN1(I) C FO ABSORPTION CORRECTION FROM FILE SHX CREATE THE FILE NEW C AND PREPARE THE FILE SC1 FOR NEW CYCLE OF ABSORPTION CORRECTIONS C MLSNE = MLSNE+1 C WRITE(*,502) MLSNE OPEN(SHX,FILE='SHX') OPEN(NEW,FILE='NEW',STATUS='NEW') OPEN(SC1,FILE='SC1',STATUS='NEW') OPEN(SC2,FILE='SC2',STATUS='NEW') NCO=0 WRITE(*,505) 610 PAR=0. ISET=0 611 READ(SHX,101,END=700)IH,IK,IL,FO,FC,SIGFO NCO=NCO+1 IF(MOD(NCO,LCO).EQ.0)WRITE(*,509)NCO 509 FORMAT(1H+,25X,I5) IF(IH.EQ.IK.AND.IK.EQ.IL.AND.IL.EQ.0) GO TO 700 SIGMA=SIGFO CALL ANG(IGEOM,IH,IK,IL,UB,WAVE,TET,CHI,PHI) CALL PHIMU(WAVE,TET,CHI,PHI,PHIP,PHIS,U) APS=1. DO 620 I=1,17 CALL SENOMN(MN(1,I),MN(2,I),PHIP,PHIS,U,SMN) CALL COSMN(MN(1,I+17),MN(2,I+17),PHIP,PHIS,U,CMN) APS=APS+SMN*TN(I)+CMN*TN(I+17) 620 CONTINUE IF(APS.LT.0.)CALL CORR(1,IH,IK,IL,TET,CHI,PHI,FO,APS) FM=SCALA*FO*APS ATET=1. DO 650 M=1,6 SCM(M)=SIN(M*TET) SCM(M+6)=COS(M*TET) ATET=ATET+SCM(M)*TN1(M)+SCM(M+6)*TN1(M+6) 650 CONTINUE IF(ATET.LT.0.)CALL CORR(2,IH,IK,IL,TET,CHI,PHI,FO,ATET) FM=FM*ATET SIGFO=SIGFO*SCALA IF(ICICLO.NE.NCICLI) GO TO 690 WRITE(NEW,1060)IH,IK,IL,FM,SIGFO,ISET,PAR 690 WRITE(SC2,101)IH,IK,IL,FM,FC,SIGFO IF(FO.LT.3.0*SIGMA)GO TO 611 WRITE(SC1,1070)IH,IK,IL,FM,SIGFO,FC IF(APS.GT.APSMAX)APSMAX=APS IF(APS.LT.APSMIN)APSMIN=APS IF(ATET.GT.ATMAX)ATMAX=ATET IF(ATET.LT.ATMIN)ATMIN=ATET GO TO 610 700 CONTINUE WRITE(*,509)NCO IF(ICICLO.NE.NCICLI) GO TO 705 WRITE(NEW,1060) IH0,IK0,IL0,FM,SIGFO,ISET,PAR 705 CONTINUE C MLSNE = MLSNE+1 C WRITE(*,502) MLSNE C EVALUATE THE R FACTOR AS R=SUM(FO-FC)/SUM(FO) AFTER ABSORPTION C CORRECTION AND BEFORE ABSORPTION CORRECTION. C EVALUATE NEW SUM(FM) AND SUM(FC) TO CALCULATE NEW SCALE FACTOR C FOR NEXT CYCLE NNN=0 SUMFC=0. SUMFM=0. DELFC=0. DELFM=0. OPEN(SC1,FILE='SC1') OPEN(SCR,FILE='SCR') OPEN(SCZ,FILE='SCZ',STATUS='NEW') NCO=0 WRITE(*,505) 850 READ(SCR,1070,END=800)IH,IK,IL,FO,SIGFO,FC READ(SC1,1070,END=970)IHN,IKN,ILN,FM,SIGFO,FC NCO=NCO+1 IF(MOD(NCO,LCO).EQ.0)WRITE(*,710)NCO 710 FORMAT(1H+,30X,I5) IF(IH.NE.IHN.OR.IK.NE.IKN.OR.IL.NE.ILN)GO TO 970 NNN=NNN+1 AAA=FO*SCALA IF (ICICLO.GT.1) GO TO 333 WRITE (SCZ,1070) IHN,IKN,ILN,AAA,SIGFO 333 DELFC=DELFC+ABS(AAA-ABS(FC)) SUMFO=SUMFO+FO*SCALA DELFM=DELFM+ABS(FM-ABS(FC)) SUMFM=SUMFM+FM SUMFC=SUMFC+ABS(FC) GO TO 850 800 RINIZ=DELFC/SUMFO RFINAL=DELFM/SUMFM WRITE(LP,1130)NNN,RINIZ,RFINAL WRITE(*,710)NCO C PRINT MAX AND MIN VALUES OF THE ABSORPTION COEFFICIENTS. WRITE(LP,1140)APSMAX,APSMIN,ATMAX,ATMIN C SET IN NEW CYCLE OF LEAST SQUARES ABSORPTION CORRECTION IF(ICICLO.EQ.NCICLI)GOTO 2001 C MLSNE = MLSNE+1 C WRITE(*,502) MLSNE OPEN(SCR,FILE='SCR') OPEN(SC1,FILE='SC1') NCO=0 WRITE(*,505) 940 READ(SC1,1070,END=950)IH,IK,IL,FM,SIGFO,FC NCO=NCO+1 IF(MOD(NCO,LCO).EQ.0)WRITE(*,711)NCO 711 FORMAT(1H+,35X,I5) WRITE(SCR,1070)IH,IK,IL,FM,SIGFO,FC GO TO 940 950 CONTINUE WRITE(*,711)NCO OPEN(SC2,FILE='SC2') OPEN(SC3,FILE='SC3',STATUS='NEW') NCO=0 WRITE(*,505) 990 READ(SC2,101,END=980)IH,IK,IL,FM,FC,SIGFO NCO=NCO+1 IF(MOD(NCO,LCO).EQ.0)WRITE(*,512)NCO 512 FORMAT(1H+,40X,I5) WRITE(SC3,101)IH,IK,IL,FM,FC,SIGFO GO TO 990 980 CONTINUE WRITE(*,512)NCO SUMFC=SUMFC SUMFO=SUMFM ICICLO=ICICLO+1 GO TO 2000 970 WRITE(LP,1210) C MLSNE = MLSNE+1 C WRITE(*,502) MLSNE STOP 2001 GOTO(2002,2003),INMOD 2002 CALL OUTF STOP 2003 CALL OUTF2 STOP 101 FORMAT(3I4,3F8.2) C105 FORMAT(6A1,1X,6A1) 1000 FORMAT(1H1,20(1H.),'PROGRAM ABSORB ',20(1H.)) 1010 FORMAT(A72) 1020 FORMAT(/1X,A72) 1030 FORMAT(/1X,'PHILIPS PW 1100 DIFFRACTION GEOMETRY ') 1040 FORMAT(1X,'SIEMENS A.E.D. DIFFRACTION GEOMETRY') 1050 FORMAT(1X,'WAVELENGTH:',F8.6) 1052 FORMAT(//1X,'UB(3,3) ORIENTATION MATRIX ',(/1X,3F15.7)) 1060 FORMAT(3I4,2F8.2,I4,6F8.5) 1070 FORMAT(3I4,3F8.2) 1072 FORMAT(///,'NUMBER OF REFLECTIONS USERD TO CREATE ',1X, 1 /'THE LINEAR SYSTEM OF ORDER IN THE POLAR ANGLES ',1X, 2/'PHI AND MU :',I5) 1073 FORMAT(1H1,'NUMBER OF REFLECTIONS USED TO CREATE',1X, 1/'THE LINEAR SYSTEM OF ORDER 12 IN THE SCATTERING ANGLE',1X, 2/'THETA :',I5) 1080 FORMAT(//1X,'SCALE FACTOR EVALUATED ON',I5,'REFLECTIONS',F15.7) 1090 FORMAT(//1X,'COEFFICIENTS FOR FOURIER SERIES EXPANSION', 1/1X,'IN PHI AND MU OF THE ABSORPTION CORRECTIONS:', 2/5X,'M',5X,'N',13X,'PQ',/,(1X,I5,1X,I5,F15.6)) 1100 FORMAT(1X,'COEFFICIENTS FOR FOURIER SERIES EXPANSION', 1/1X,'IN THETA OF THE ABSORPTION CORRECTIONS :', 2/5X,'N',13X,'RS') 1110 FORMAT(1X,I5,F15.6) 1130 FORMAT(1H1,'TEST OF FO-FC AGREEMENT ON',I5,'REFLECTIONS:', 2//1X,'BEFORE FO''S ABSORPTION CORRECTIONS', 3/5X,'R=SUM(ABS(K*FO-ABS(FC)))/SUM(FO)=',F15.6, 4//1X,'AFTER FO',1H','S ABSORPTION CORRECTIONS', 5/5X,'R=SUM(ABS(FM-ABS(FC)))/SUM(FM)=',F15.6 6/1X,'WHERE K IS THE SCALE FACTOR GIVEN ABOVE , AND', 7 / 'FM ARE THE FO',1H','S CORRECTED FOR ABSORPTION.') 1140 FORMAT(////1X,'MAX AND MIN VALUES OF THE ABSORPTION', 1/1X,'CORRECTIONS', 2/1X,'EVALUATED ON THE REFLECTIONS WHICH SATISFY:', 2'FO>3.0*SIGMA(FO)',/1X, 3'CORRECTION IN PHI AND MU',10X,'MAX=',F12.6,10X,'MIN=',F12.6, 4//1X,'CORRECTION IN THETA ',13X,'MAX=',F12.6,10X,'MIN=',F12.6) 1200 FORMAT(1H1,'CYCLE',I3/) 1210 FORMAT(1X,'SEQUENCE ERROR PLEASE RERUN WITH NCYCLE=1') END SUBROUTINE ANG(IGEOM,IH,IK,IL,UB,WAVE,TET,CHI,PHI) C IGEOM=1 PHILIPS PW1100 GEOMETRY C IGEOM=2 SIEMENSA.E.D.GEOMETRY REAL*8 UB(3,3) REAL NORM C* LP=2 LP=11 PIG=3.14159265 X=UB(1,1)*IH+UB(1,2)*IK+UB(1,3)*IL Y=UB(2,1)*IH+UB(2,2)*IK+UB(2,3)*IL Z=UB(3,1)*IH+UB(3,2)*IK+UB(3,3)*IL RAD=SQRT(X*X+Y*Y+Z*Z) NORM=SQRT(X*X+Y*Y) GO TO (100,200),IGEOM 100 SNCHI=-Z/RAD SENT=RAD/2.0 IF(SENT.LT.1.0) GO TO 110 SENT=1.0 WRITE(LP,50)IH,IK,IL 110 CHI=ATAN2(SNCHI,SQRT(1.0-SNCHI*SNCHI)) IF(CHI.GT.PIG)CHI=CHI-2.0*PIG TET=ATAN2(SENT,SQRT(1.0-SENT*SENT)) SNPHI=X/NORM COSPHI=Y/NORM PHI=ATAN2(SNPHI,COSPHI) RETURN 200 SENT=RAD*WAVE/2.0 IF(SENT.LT.1.0) GO TO 210 SENT=1.0 WRITE(LP,50)IH,IK,IL 210 TET=ATAN2(SENT,SQRT(1.0-SENT*SENT)) SNCHI=Z/RAD COSCHI=NORM/RAD CHI=ATAN2(SNCHI,COSCHI) SNPHI=Y/NORM COSPHI=X/NORM PHI=ATAN2(SNPHI,COSPHI) IF((CHI*57.29577).GT.-5.0)RETURN CHI=-CHI PHI=PHI+PIG IF(PHI.GE.2.0*PIG)PHI=PHI-(2.0*PIG) RETURN 50 FORMAT(1X,'CONFLICTING BETWEEN MILLER INDICES',3I6, 1/1X,'AND UB ORIENTATION MATRIX') END SUBROUTINE PHIMU(WAVE,TET,CHI,PHI,PHIP,PHIS,U) PIG=3.14159265 ARG=SIN(TET)*SIN(CHI) UU=ATAN2(ARG,(1.0-ARG*ARG)) U=SQRT(ABS(UU)) IF(UU.LT.0)U=-U U=U*7.0*WAVE TG=TAN(TET) CS=COS(CHI) ARG=1.0/SQRT(1.0+TG*TG*CS*CS) DPHI=ATAN2(SQRT(1.0-ARG*ARG),ARG) PHIP=PHI+DPHI PHIS=PHI+PIG-DPHI RETURN END SUBROUTINE SENOMN(M,N,PHIP,PHIS,U,SMN) REAL*8 SMN SMN=SIN(M*U+N*PHIP)+SIN(M*U+N*PHIS) RETURN END SUBROUTINE COSMN(M,N,PHIP,PHIS,U,CMN) REAL*8 CMN CMN=COS(M*U+N*PHIP)+COS(M*U+N*PHIS) RETURN END SUBROUTINE GAUSJO(A,B,N) C THIS ROUTINE SOLVES A SYSTEM OF N LINEAR EQUATIONS C WITH N LESS OR EQUAL TO 34 C USING THEGAUSS JORDAN METHOD C DOUBLE PRECISION A,B,RP,PV,T,RIS C DIMENSION A(34,34),B(34),RP(35) REAL*8 A(34,34),B(34),RP(35),PV,T,RIS T=0. NP1=N+1 NM1=N-1 DO 10 L=1,N PV=0. LM=N-L+1 DO 1 K=1,LM IF(DABS(A(K,1)).LT.PV) GO TO 1 PV=DABS(A(K,1)) KK=K 1 CONTINUE IF(PV.GT.T) GO TO 2 RETURN 2 DO 3 J=1,N RIS=A(KK,J) A(KK,J)=A(1,J) A(1,J)=RIS 3 CONTINUE RIS=B(KK) B(KK)=B(1) B(1)=RIS DO 4 J=1,NM1 4 RP(J)=A(1,J+1)/A(1,1) RP(N)=B(1)/A(1,1) RP(NP1)=1./A(1,1) DO 6 I=2,N RIS=A(I,2)-A(I,1)*RP(1) DO 5 J=2,NM1 5 A(I,J)=A(I,J+1)-A(I,1)*RP(J) A(I,N)=B(I)-A(I,1)*RP(N) B(I)=-A(I,1)*RP(NP1) A(I,1)=RIS 6 CONTINUE DO 8 J=1,N DO 7 I=1,NM1 7 A(I,J)=A(I+1,J) A(N,J)=RP(J) 8 CONTINUE DO 9 I=1,NM1 9 B(I)=B(I+1) B(N)=RP(NP1) 10 CONTINUE DO 12 I=1,N RIS=A(I,1) DO 11 J=1,NM1 11 A(I,J)=A(I,J+1) A(I,N)=B(I) B(I)=RIS 12 CONTINUE RETURN END SUBROUTINE PROD(A,B,C) REAL*8 A(3,3),B(3,3),C(3,3) DO 1 I=1,3 DO 1 J=1,3 1 C(I,J)=0. DO 2 I=1,3 DO 2 J=1,3 DO 2 K=1,3 C(I,J)=A(I,K)*B(K,J)+C(I,J) 2 CONTINUE RETURN END C SUBROUTINE CORR(N,IH,IK,IL,TET,CHI,PHI,FO,A) C* LP=2 LP=11 RAD=57.29577 T=TET*RAD C=CHI*RAD P=PHI*RAD IF(N.EQ.2)GO TO 77 WRITE(LP,100)IH,IK,IL,T,C,P,FO,A 100 FORMAT(1X,'H',I3,3X,'K',I3,3X,'L',I3,3X,'THETA',F7.3,3X, 1'CHI',F7.3,3X,'PHI',F7.3,3X,'FO',F7.2,3X,'APS ',F8.4) GO TO 444 77 WRITE(LP,109)IH,IK,IL,T,C,P,FO,A 109 FORMAT(1X,'H',I3,3X,'K',I3,3X,'L',I3,3X,'THETA',F7.3,3X, 1'CHI',F7.3,3X,'PHI',F7.3,3X,'FO',F7.2,3X,'ATET',F8.4) 444 A=-A RETURN END SUBROUTINE INV(A,B) REAL*8 A(3,3),B(3,3),DET C* LP=2 LP=11 B(1,1)=A(2,2)*A(3,3)-A(2,3)*A(3,2) B(1,2)=-(A(1,2)*A(3,3)-A(1,3)*A(3,2)) B(1,3)=A(1,2)*A(2,3)-A(1,3)*A(2,2) B(2,1)=-(A(2,1)*A(3,3)-A(2,3)*A(3,1)) B(2,2)=A(1,1)*A(3,3)-A(1,3)*A(3,1) B(2,3)=-(A(1,1)*A(2,3)-A(1,3)*A(2,1)) B(3,1)=A(2,1)*A(3,2)-A(2,2)*A(3,1) B(3,2)=-(A(1,1)*A(3,2)-A(1,2)*A(3,1)) B(3,3)=A(1,1)*A(2,2)-A(1,2)*A(2,1) DET=0. DO 90 L=1,3 90 DET=DET+A(1,L)*B(L,1) IF(DABS(DET).GT.1.0D-12) GO TO 100 WRITE(LP,120)DET 120 FORMAT(1X,'SUBR. INV:DET=0. ',E15.8) STOP 100 DO 95 M=1,3 DO 95 L=1,3 95 B(M,L)=B(M,L)/DET RETURN END C SUBROUTINE IN4 CHARACTER*1 FLAG DATA IH0,IK0,IL0,LCO,THRESL/0,0,0,500,1.0E-18/ WRITE(*,503) 503 FORMAT(3X,'Reading LIST4...'/) OPEN(16,FILE='LIST4.DAT',STATUS='OLD') OPEN(12,FILE='SHX',STATUS='NEW') NCO=0 45 READ(16,102,END=40)IH,IK,IL,FCS,FOS,SIGFOS,FLAG NCO=NCO+1 IF(MOD(NCO,LCO).EQ.0) WRITE(*,504)NCO 504 FORMAT(1H+,5X,I5) 102 FORMAT(3I4,2F12.2,F10.2,1X,A1) FO=SQRT(ABS(FOS)) FC=SQRT(FCS) IF(FO.LT.THRESL) THEN SIGFO=0.01 ELSE SIGFO=SIGFOS/(2.*FO) ENDIF 511 IF(FOS.LT.0.) FO=-FO WRITE(12,101)IH,IK,IL,FO,FC,SIGFO GO TO 45 40 CONTINUE WRITE(*,504)NCO C MLSNE = MLSNE+1 C WRITE(*,502) MLSNE 101 FORMAT(3I4,3F8.2) RETURN END C SUBROUTINE IN3 DATA IH0,IK0,IL0,LCO,THRESL/0,0,0,500,1.0E-18/ WRITE(*,503) 503 FORMAT(3X,'Reading LIST3...'/) OPEN(16,FILE='LIST3.DAT',STATUS='OLD') OPEN(12,FILE='SHX',STATUS='NEW') NCO=0 45 READ(16,*,END=40)IH,IK,IL,FO,SIGFO,A,B NCO=NCO+1 IF(MOD(NCO,LCO).EQ.0) WRITE(*,504)NCO 504 FORMAT(1H+,5X,I5) 102 FORMAT(3I4,4F8.2) FC=SQRT(A*A + B*B) IF(SIGFO.LT.THRESL) SIGFO = 0.01 WRITE(12,101)IH,IK,IL,FO,FC,SIGFO GO TO 45 40 CONTINUE WRITE(*,504)NCO C MLSNE = MLSNE+1 C WRITE(*,502) MLSNE 101 FORMAT(3I4,3F8.2) RETURN END C SUBROUTINE OUTF WRITE(*,503) 503 FORMAT(3X,'Writing F-data in NEW.HKL...'/) OPEN(22,FILE='NEW') OPEN(3,FILE='NEWF.HKL',STATUS='NEW') NCO=0 20 READ(22,1060,END=100) IH,IK,IL,FM,SIGFO,ISET,PAR NCO=NCO+1 IF(MOD(NCO,LCO).EQ.0) WRITE(*,504)NCO 504 FORMAT(1H+,35X,I5) 1060 FORMAT(3I4,2F8.2,I4,6F8.5) WRITE(3,10)IH,IK,IL,FM,SIGFO GOTO 20 10 FORMAT(3I4,2F8.2) 100 WRITE(*,504)NCO RETURN END C SUBROUTINE OUTF2 WRITE(*,503) 503 FORMAT(3X,'Writing F2-data in NEW.HKL...'/) OPEN(22,FILE='NEW') OPEN(3,FILE='NEWF2.HKL',STATUS='NEW') NCO=0 20 READ(22,1060,END=100) IH,IK,IL,FM,SIGFO,ISET,PAR NCO=NCO+1 IF(MOD(NCO,LCO).EQ.0) WRITE(*,504)NCO 504 FORMAT(1H+,35X,I5) 1060 FORMAT(3I4,2F8.2,I4,6F8.5) F2=FM*FM SIGF2=SIGFO*2.0*ABS(FM) IF (FM.LT.0.)F2=-F2 WRITE(3,10)IH,IK,IL,F2,SIGF2 GOTO 20 10 FORMAT(3I4,2F8.2) 100 WRITE(*,504)NCO RETURN END C C------------------------------------------------------------------- C Sample input file CR C TITLE CI2 IN C2/c :title C 1 :IGEOM (Philip's) : format-free C 0.71073 :wave-length : ,, C 0.041418 0.032050 -0.002130 :orientation matrix: ,, C -0.032674 0.038381 0.010599 : : C 0.017842 -0.004142 0.081845 : : C 1 :number of cycles : ,, C--------------------------------------------------------------------