PROGRAM VIEW C C ROBERT H. BLESSING C MEDICAL FOUNDATION OF BUFFALO C 73 HIGH STREET C BUFFALO, NEW YORK 14203, USA C TELEPHONE: (716) 856-9600 C ELECTRONIC MAIL: Blessing@MFB.Buffalo.Edu C C JULY 1994 C CHARACTER LABEL*2 COMMON /VTKB/ IVT,IKB,IC COMMON /SIZE/ NXVT,NYVT,NXLP,NYLP COMMON /PLOT/ ITYPE,LU,NX,NY,ILOG,FSCALE,ISMOOTH COMMON /LABL/ LABEL(4) COMMON /RFLN/ IEND,II,IH,IK,IL,ANGLES(4),WIDTH,SPEED,FBG, & BG1,SCAN,BG2,YNET,SIGY,XTIME COMMON /QLIM/ QMIN,QMAX COMMON /BG/ A0,A1,E0,E1,B0,B1,J1,J2 C C COMMON BLOCKS LABELED /BLOCKI/ ARE FROM PROGRAM BGLP FOR VIEWING C PROCESSED DATA. C PARAMETER (NMAX=100) CHARACTER ATIME*8,ADATE*9,TITLE*80,FILE1*80,FILE2*80 COMMON /BLOCK0/ LU1,LU2,LU3,ILP COMMON /BLOCKC/ NCOND,ICOND(38) COMMON /BLOCK1/ ATIME,ADATE,TITLE,FILE1,FILE2 COMMON /BLOCK2/ GINV(3,3),WLE,WL1,WL0,WL2 COMMON /BLOCK3/ WINDOW,CUTOFF,NSMOOTH,II1,II2,XT1,XT2,SC1,SC2, & XX0,XW1,XW2 COMMON /BLOCK4/ C0KIN,C1KIN,C2KIN,C0DYN,C1DYN,C2DYN,FRACTD,VARFD, & INEUTRON COMMON /BLOCK5/ IDIFF,U(3,3),MODEL1,MODEL2,Q1(3,3),Q2(3,3),T1,T2, & F,SIGU(3,3),SIGQ1(3,3),SIGQ2(3,3),SIGT1,SIGT2,TANTHM COMMON /BLOCK6/ VARTHE,VARTIM,TAU,VARTAU,ATT,VARATT COMMON /BLOCK7/ NSTEP,XX(NMAX),YY(NMAX),VX(NMAX),VY(NMAX) COMMON /BLOCK8/ TH0,X0,SIGX0,XE,X1,X2,W1,W2,I1,I2,L1,L2,R,SIGR,EPS CALL STARTUP CALL DRIVER END C----------------------------------------------------------------------- BLOCK DATA CHARACTER LABEL*2 COMMON /VTKB/ IVT,IKB,IC COMMON /SIZE/ NXVT,NYVT,NXLP,NYLP COMMON /PLOT/ ITYPE,LU,NX,NY,ILOG,FSCALE,ISMOOTH COMMON /LABL/ LABEL(4) COMMON /RFLN/ IEND,II,IH,IK,IL,ANGLES(4),WIDTH,SPEED,FBG, & BG1,SCAN,BG2,YNET,SIGY,XTIME COMMON /QLIM/ QMIN,QMAX COMMON /BG/ A0,A1,E0,E1,B0,B1,J1,J2 C C COMMON BLOCKS LABELED /BLOCKI/ ARE FROM PROGRAM BGLP FOR VIEWING C PROCESSED DATA. C PARAMETER (NMAX=100) CHARACTER ATIME*8,ADATE*9,TITLE*80,FILE1*80,FILE2*80 COMMON /BLOCK0/ LU1,LU2,LU3,ILP COMMON /BLOCKC/ NCOND,ICOND(38) COMMON /BLOCK1/ ATIME,ADATE,TITLE,FILE1,FILE2 COMMON /BLOCK2/ GINV(3,3),WLE,WL1,WL0,WL2 COMMON /BLOCK3/ WINDOW,CUTOFF,NSMOOTH,II1,II2,XT1,XT2,SC1,SC2, & XX0,XW1,XW2 COMMON /BLOCK4/ C0KIN,C1KIN,C2KIN,C0DYN,C1DYN,C2DYN,FRACTD,VARFD, & INEUTRON COMMON /BLOCK5/ IDIFF,U(3,3),MODEL1,MODEL2,Q1(3,3),Q2(3,3),T1,T2, & F,SIGU(3,3),SIGQ1(3,3),SIGQ2(3,3),SIGT1,SIGT2,TANTHM COMMON /BLOCK6/ VARTHE,VARTIM,TAU,VARTAU,ATT,VARATT COMMON /BLOCK7/ NSTEP,XX(NMAX),YY(NMAX),VX(NMAX),VY(NMAX) COMMON /BLOCK8/ TH0,X0,SIGX0,XE,X1,X2,W1,W2,I1,I2,L1,L2,R,SIGR,EPS DATA IKB,IVT /5, 6/ DATA LU1,LU2,LU3,ILP /1, 2, 3, 4/ DATA NXVT,NYVT,NXLP,NYLP /48, 20, 96, 50/ DATA WINDOW,CUTOFF,NSMOOTH,II1,II2,XT1,XT2,SC1,SC2 /0.08333, 3, 2, & 0, 999999, 0, 999999, 0, 1/ DATA XX0,XW1,XW2 /0.5, 0, 0/ DATA U,MODEL1,MODEL2,Q1,Q2,T1,T2,F,SIGU,SIGQ1,SIGQ2,SIGT1,SIGT2, & TANTHM /9*0, 0, 0, 9*0, 9*0, 0, 0, 0, 9*0, 9*0, 9*0, 0, 0, 0/ DATA VARTHE,VARTIM,TAU,VARTAU,ATT,VARATT /0.001, 0, 0, 0, 1, 0/ END C----------------------------------------------------------------------- SUBROUTINE STARTUP CHARACTER ATIME*8,ADATE*9,TITLE*80,FILE1*80,FILE2*80, & BTIME*8,BDATE*9,HCOND*18,W*3,LABEL*2,ANS*1,DIFF*4 COMMON /VTKB/ IVT,IKB,IC COMMON /SIZE/ NXVT,NYVT,NXLP,NYLP COMMON /PLOT/ ITYPE,LU,NX,NY,ILOG,FSCALE,ISMOOTH COMMON /LABL/ LABEL(4) COMMON /BLOCK0/ LU1,LU2,LU3,ILP COMMON /BLOCKC/ NCOND,ICOND(38) COMMON /BLOCK1/ ATIME,ADATE,TITLE,FILE1,FILE2 COMMON /BLOCK3/ WINDOW,CUTOFF,NSMOOTH,II1,II2,XT1,XT2,SC1,SC2, & XX0,XW1,XW2 COMMON /BLOCK5/ IDIFF,U(3,3),MODEL1,MODEL2,Q1(3,3),Q2(3,3),T1,T2, & F,SIGU(3,3),SIGQ1(3,3),SIGQ2(3,3),SIGT1,SIGT2,TANTHM COMMON /BLOCK6/ VARTHE,VARTIM,TAU,VARTAU,ATT,VARATT CALL TIME (ATIME) CALL DATE (ADATE) c ATIME='hh:mm:ss' c ADATE='dd-mmm-yy' C C OPEN (UNIT=IVT,FILE='SYS$OUTPUT',STATUS='NEW') C OPEN (UNIT=IKB,FILE='SYS$INPUT',STATUS='NEW') C WRITE (IVT,1001) ATIME,ADATE 1001 FORMAT ('0PROGRAM VIEW/BGLP. ',A,1X,A) C C 10 WRITE (IVT,1002) C 1002 FORMAT (' SCREEN WIDTH? 80 OR 132 (132). ',$) C READ (IKB,1000) W C 1000 FORMAT (A) C IF (W.EQ.' ') W='132' C IF (W.NE.'80'.AND.W.NE.'132') GO TO 10 C IF (W.EQ.'80') NXVT=48 C IF (W.EQ.'132') NXVT=96 C 20 WRITE (IVT,1003) C 1003 FORMAT (' VIEW RAW OR PROCESSED DATA? R OR P (P). ',$) C READ (IKB,1000) ANS C IF (ANS.EQ.' ') ANS='P' C IF (ANS.NE.'R'.AND.ANS.NE.'P') GO TO 20 C IF (ANS.EQ.'R') ITYPE=0 C IF (ANS.EQ.'P') ITYPE=1 C IF (ANS.EQ.'P'.AND.W.NE.'132') THEN C WRITE (IVT,1009) C STOP C END IF C 1009 FORMAT ('0SET SCREEN WIDTH TO 132 TO VIEW PROCESSED DATA.') C 10 WRITE (IVT,1002) 1002 FORMAT (' LOCK KEYBOARD FOR UPPER CASE LETTERS.'/' WINDOW AREA MAX &IMIZED (AT LEAST 125 X 55 CHARACTERS)? Y OR N (Y). ',$) READ (IKB,1000) ANS 1000 FORMAT (A) IF (ANS.EQ.' ') ANS='Y' IF (ANS.NE.'Y'.AND.ANS.NE.'N') GO TO 10 IF (ANS.EQ.'Y') THEN NXVT=96 NYVT=40 ELSE STOP 'LOCK UPPER CASE, MAXIMIZE WINDOW, AND RESTART PROGRAM.' END IF 20 WRITE (IVT,1003) 1003 FORMAT (' VIEW RAW OR PROCESSED DATA? R OR P (P). ',$) READ (IKB,1000) ANS IF (ANS.EQ.' ') ANS='P' IF (ANS.NE.'R'.AND.ANS.NE.'P') GO TO 20 IF (ANS.EQ.'R') ITYPE=0 IF (ANS.EQ.'P') ITYPE=1 C C IF ITYPE = 0, VIEW RAW DATA PROFILES. C IF (ITYPE.EQ.0) THEN 11 WRITE (IVT,1004) READ (IKB,1000) FILE1 OPEN (UNIT=LU1,FILE=FILE1,STATUS='OLD',FORM='UNFORMATTED', & ERR=11) WRITE (IVT,1005) READ (IKB,1000) DIFF OPEN (UNIT=ILP,FILE='view.lp',STATUS='NEW', & CARRIAGECONTROL='FORTRAN') END IF 1004 FORMAT ('0ENTER RAW REFLECTION DATA FILE NAME.') 1005 FORMAT ('0ENTER DIFFRACTOMETER TYPE. CAD4, P3, OR OTHER.') C C IF ITYPE = 1, VIEW PROCESSED DATA PROFILES USING SUBROUTINES FROM C PROGRAM BGLP AND THE BGLP.DAT PARAMETERS FROM PROGRAM REFPK. C IF (ITYPE.EQ.1) THEN CALL INPUT (BTIME,BDATE) WRITE (IVT,1006) FILE1 WRITE (IVT,1007) BTIME,BDATE,TITLE IF (IDIFF.EQ.4) DIFF='CAD4' IF (IDIFF.EQ.3) DIFF='P3' WRITE (ILP,1001) ATIME,ADATE WRITE (ILP,1006) FILE1 WRITE (ILP,1007) BTIME,BDATE,TITLE IF (NCOND.NE.0) THEN WRITE (IVT,1010) DO 15 I=1,NCOND CALL HKLCOND (HCOND,ICOND(I)) WRITE (IVT,1011) HCOND 15 CONTINUE END IF END IF 1006 FORMAT ('0RAW DATA FILE: '/' ',A) 1007 FORMAT ('0BGLP.DAT PARAMETERS FROM PROGRAM REFPK JOB:'/ &' TIME DATE TITLE'/' ',A,1X,A,1X,A) 1010 FORMAT ('0CONDITIONS LIMITING POSSIBLE REFLECTIONS:'/) 1011 FORMAT (' ',A) C C STORE LABELS FOR DIFFRACTOMETER ANGLES FOR PROGRAM VIEW OUTPUT. C IF (DIFF.EQ.'CAD4') THEN LABEL(1)='TH' LABEL(2)='PH' LABEL(3)='OM' LABEL(4)='KA' ELSE IF (DIFF.EQ.'P3') THEN LABEL(1)='TT' LABEL(2)='OM' LABEL(3)='PH' LABEL(4)='CH' ELSE LABEL(1)='A1' LABEL(2)='A2' LABEL(3)='A3' LABEL(4)='A4' END IF RETURN END C----------------------------------------------------------------------- SUBROUTINE DRIVER PARAMETER (NMAX=100) CHARACTER ATIME*8,ADATE*9,TITLE*80,FILE1*80,FILE2*80 COMMON /VTKB/ IVT,IKB,IC COMMON /RFLN/ IEND,II,IH,IK,IL,ANGLES(4),WIDTH,SPEED,FBG, & BG1,SCAN,BG2,YNET,SIGY,XTIME COMMON /SIZE/ NXVT,NYVT,NXLP,NYLP COMMON /PLOT/ ITYPE,LU,NX,NY,ILOG,FSCALE,ISMOOTH COMMON /BLOCK0/ LU1,LU2,LU3,ILP COMMON /BLOCK1/ ATIME,ADATE,TITLE,FILE1,FILE2 COMMON /BLOCK7/ NSTEP,XX(NMAX),YY(NMAX),VX(NMAX),VY(NMAX) COMMON /BLOCK8/ TH0,X0,SIGX0,XE,X1,X2,W1,W2,I1,I2,L1,L2,R,SIGR,EPS CHARACTER A*2,SEARCH*2 SAVE LU=IVT NX=NXVT NY=NYVT 2 WRITE (IVT,999) 999 FORMAT ('0TYPE "?" OR "O" TO LIST OPTIONS.') 1 L=ILOG F=FSCALE S=ISMOOTH ILOG=0 FSCALE=1 ISMOOTH=0 C C ENTER OPTION. C WRITE (IVT,1001) 1001 FORMAT (' OPTION? ',$) READ (IKB,1000) A 1000 FORMAT (A) IF (A.EQ.'I'.OR.A.EQ.'H'.OR.A.EQ.'HA'.OR.A.EQ.'T'.OR.A.EQ.'S'.OR. & A.EQ.'Q'.OR.A.EQ.'PK'.OR.A.EQ.'BG'.OR.A.EQ.'ER') THEN SEARCH=A IC=0 ELSE IF (A.EQ.'C') THEN A=SEARCH IC=1 END IF IF (A.EQ.'?') GO TO 11 IF (A.EQ.'O') GO TO 11 IF (A.EQ.'R') GO TO 12 IF (A.EQ.'A') GO TO 13 IF (A.EQ.'B') GO TO 14 IF (A.EQ.'L') GO TO 15 IF (A.EQ.'I') GO TO 21 IF (A.EQ.'H') GO TO 22 IF (A.EQ.'HA') GO TO 23 IF (A.EQ.'T') GO TO 24 IF (A.EQ.'S') GO TO 25 IF (A.EQ.'Q') GO TO 26 IF (A.EQ.'PK') GO TO 31 IF (A.EQ.'BG') GO TO 32 IF (A.EQ.'ER') GO TO 33 IF (A.EQ.'LN') GO TO 50 IF (A.EQ.'SM') GO TO 51 IF (A.EQ.'1') GO TO 52 IF (A.EQ.'2') GO TO 52 IF (A.EQ.'5') GO TO 52 IF (A.EQ.'M') GO TO 53 IF (A.EQ.'W') GO TO 54 IF (A.EQ.'P') GO TO 55 IF (A.EQ.'DE') GO TO 56 IF (A.EQ.'MD') GO TO 57 IF (A.EQ.'E') GO TO 99 IF (A.EQ.'X') GO TO 99 IF (A.NE.' ') GO TO 2 C C DEFAULT: DISPLAY THE NEXT PROFILE. C 10 CALL READR CALL PLOT0 GO TO 1 C C LIST USER'S OPTIONS. C 11 WRITE (IVT,1011) 1011 FORMAT ('0OPTIONS ARE:'/ &'0LIST USER''S OPTIONS O or ? REWIND FILE R'/ &' ADVANCE N REFLNS. A BACKSPACE N REFLNS. B'/ &' LIST NEXT 20 REFLNS. L SEARCH SERIAL NO. I'/ &' SEARCH HKL H SEARCH SYST. ABSNT. HKL HA'/ &' SEARCH TWO-THETA T SEARCH SIN(TH)/L S'/ &' SEARCH Q = I/SIGMA(I) Q SEARCH DISPLCD. PKS. PK'/ &' SEARCH IRRGLR. BKGDS. BG SEARCH SCAN ERRORS ER'/ &' CONTINUE SEARCH C CHANGE TO LOG SCALE LN'/ &' SMOOTH PROFILE PLOT SM DISPLAY ON SCALE 1X 1'/ &' MULTIPLY SCALE 2X 2 MULTIPLY SCALE 5X 5'/ &' MOVE PK/BG POINT(S) M WRITE REVISED REFLN. W'/ &' PRINT PROFILE PLOT P DELETE MEASUREMENT DE'/ &' MODIFY SCAN PROFILE MD EXIT E or X &'/ &' DEFAULT: DISPLAY NEXT PROFILE') GO TO 1 C C REWIND FILE 1. C 12 REWIND LU1 GO TO 1 C C ADVANCE N REFLECTIONS. C 13 WRITE (IVT,1013) 1013 FORMAT (' ENTER NUMBER OF REFLECTIONS TO ADVANCE. ',$) READ (IKB,2013,ERR=13) N 2013 FORMAT (I10) IF (N.LE.0) GO TO 1 DO 113 I=1,N READ (LU1,END=12) 113 CONTINUE GO TO 1 C C BACKSPACE N REFLECTIONS. C 14 WRITE (IVT,1014) 1014 FORMAT (' ENTER NUMBER OF REFLECTIONS TO BACKSPACE. ',$) READ (IKB,2014,ERR=14) N 2014 FORMAT (I10) IF (N.LT.0) GO TO 1 DO 114 I=1,N+1 BACKSPACE (LU1,ERR=1) 114 CONTINUE GO TO 1 C C LIST THE NEXT 20 REFLECTIONS. C 15 WRITE (IVT,1015) ANGLE 1015 FORMAT (' I H K L ',4(5X,A2),' Y SIGMA(Y)', &' XTIME'/' - - - - -- -- -- -- ', &' - -------- -----') DO 115 I=1,20 CALL READR IF (IEND.NE.0) GO TO 12 IF (ITYPE.EQ.0) THEN CALL BGPKBG (0.0,NSTEP,XX,YY,VX,VY,BG1,SCAN,BG2,R,SIGR) ELSE CALL PROCESS END IF WRITE (IVT,2015) II,IH,IK,IL,ANGLES,R,SIGR,XTIME 115 CONTINUE 2015 FORMAT (' ',I6,2X,3I3,2X,4F7.2,3F10.2) GO TO 1 C C SEARCH FOR SERIAL NUMBER(S). C 21 CALL ISEARCH GO TO 1 C C SEARCH FOR HKL(S). C 22 CALL HSEARCH GO TO 1 C C SEARCH FOR SYSTEMATICALLY ABSENT HKL(S). C 23 CALL ASEARCH GO TO 1 C C SEARCH FOR TWO-THETA VALUE(S). C 24 CALL TSEARCH GO TO 1 C C SEARCH FOR SIN(THETA)/LAMBDA VALUE(S). C 25 CALL SSEARCH GO TO 1 C C SEARCH FOR NET INTENSITY Q = I/SIGMA(I) VALUE(S). C 26 CALL QSEARCH GO TO 1 C C SEARCH FOR PROFILES WITH DISPLACED PEAK. C 31 CALL PSEARCH GO TO 1 C C SEARCH FOR PROFILES WITH IRREGUALR BACKGROUND. C 32 CALL BSEARCH GO TO 1 C C SEARCH FOR SCANS WITH DISCONTINUITY ERRORS. C 33 CALL ESEARCH GO TO 1 C C CHANGE TO LOGARITHMIC COUNT-RATE SCALE. C 50 ILOG=1 IF (ITYPE.EQ.0) THEN CALL PLOT0 ELSE CALL PLOT1 END IF GO TO 1 C C SMOOTH PROFILE PLOT. C 51 ILOG=L FSCALE=F WRITE (IVT,1051) 1051 FORMAT (' ENTER NUMBER OF SMOOTHING PASSES (DEFAULT: N = 2). ',$ &) READ (IKB,2051,ERR=51) ISMOOTH 2051 FORMAT (I2) IF (ISMOOTH.EQ.0) ISMOOTH=2 IF (ISMOOTH.LT.0) ISMOOTH=0 IF (ITYPE.EQ.0) THEN CALL PLOT0 ELSE CALL PLOT1 END IF GO TO 1 C C MULTIPLY (LINEAR) COUNT-RATE SCALE. C 52 IF (A.EQ.'1') FSCALE=1 IF (A.EQ.'2') FSCALE=2*F IF (A.EQ.'5') FSCALE=5*F ISMOOTH=S IF (ITYPE.EQ.0) THEN CALL PLOT0 ELSE CALL PLOT1 END IF GO TO 1 C C MOVE PEAK/BACKGROUND POINTS. C 53 ILOG=L FSCALE=F ISMOOTH=S CALL MOVE GO TO 1 C C WRITE REVISED REFLECTION RECORD TO REVISED.BLP FILE. C 54 WRITE (LU3) II,IH,IK,IL,ANGLES,R,SIGR,XTIME GO TO 1 C C PRINT PROFILE PLOT. C 55 ILOG=L FSCALE=F ISMOOTH=S LU=ILP NX=NXLP NY=NYLP WRITE (LU,1055) ATIME,ADATE,FILE1 1055 FORMAT ('1'/'0PROGRAM VIEW/BGLP. ',A,1X,A,'. ',A) IF (ITYPE.EQ.0) THEN CALL PLOT0 ELSE CALL PLOT1 END IF LU=IVT NX=NXVT NY=NYVT GO TO 1 C C DELETE THE DISPLAYED PROFILE MEASUREMENT. C 56 WRITE (IVT,1056) 1056 FORMAT (' DELETE THE DISPLAYED PROFILE MEASUREMENT? Y OR N (N).', &' ',$) READ (IKB,1000) A IF (A.EQ.'Y') THEN SIGR=-SIGR GO TO 54 END IF GO TO 1 C C MODIFY THE SCAN PROFILE AND RE-EVALUATE THE MEASUREMENT. C 57 IF (ITYPE.EQ.0) THEN WRITE (IVT,1057) ELSE CALL MODIFY END IF 1057 FORMAT (' MODIFY COMMAND IS EXECUTED ONLY FOR PROCESSED PROFILES.' &) GO TO 1 C C TERMINATE SESSION. C 99 IF (ITYPE.EQ.1) THEN REWIND (LU3) READ (LU3,END=95) WRITE (IVT,1099) READ (IKB,1000,ERR=99) A IF (A.NE.'N') CALL REVISE END IF 1099 FORMAT (' KEEP THE REVISED DATA.BLP FILE? Y OR N (Y). ',$) 95 REWIND ILP DO 81 I=1,8 READ (ILP,1000,END=89) A 81 CONTINUE GO TO 88 89 CLOSE (UNIT=ILP,STATUS='DELETE') 88 STOP 'PROGRAM VIEW FINIS' END C----------------------------------------------------------------------- SUBROUTINE BGPKBG (F,N,X,Y,U,V,B1,PK,B2,YNET,SIGY) DIMENSION X(1),Y(1),U(1),V(1) IF (F.EQ.0) THEN Q=1/6.0 ELSE Q=F END IF Q=Q*(X(N)-X(1)) X1=X(1)+Q X2=X(N)-Q B1=0 PK=0 B2=0 VB=0 VP=0 DO 1 I=2,N-1 DX=0.50*(X(I+1)-X(I-1)) VD=0.25*(U(I+1)+U(I-1)) YI=Y(I)*DX VI=V(I)*DX**2+Y(I)**2*VD IF (X(I).LT.X1) THEN B1=B1+YI VB=VB+VI ELSE IF (X1.LE.X(I).AND.X(I).LE.X2) THEN PK=PK+YI VP=VP+VI ELSE B2=B2+YI VB=VB+VI END IF 1 CONTINUE Q=(X2-X1)/(X(N)-X(1)) Q=Q/(1-Q) YNET=PK-Q*(B1+B2) SIGY=SQRT(VP+Q**2*VB) RETURN END C----------------------------------------------------------------------- SUBROUTINE ISEARCH COMMON /VTKB/ IVT,IKB,IC COMMON /RFLN/ IEND,II,IH,IK,IL,ANGLES(4),WIDTH,SPEED,FBG, & BG1,SCAN,BG2,YNET,SIGY,XTIME COMMON /BLOCK0/ LU1,LU2,LU3,ILP SAVE IF (IC.EQ.1) GO TO 1 IMIN=0 IMAX=0 5 WRITE (IVT,1001) 1001 FORMAT (' ENTER SERIAL NUMBER(S). I OR IMIN, IMAX. ',$) READ (IKB,1002,ERR=5) IMIN,IMAX 1002 FORMAT (2I10) IF (IMIN.EQ.0.AND.IMAX.EQ.0) RETURN IMIN=ABS(IMIN) IMAX=ABS(IMAX) IF (IMAX.LT.IMIN) IMAX=IMIN 1 CALL READR IF (IEND.NE.0) RETURN I=ABS(II) IF (I.LT.IMIN.OR.I.GT.IMAX) GO TO 1 CALL PLOT0 RETURN END C----------------------------------------------------------------------- SUBROUTINE HSEARCH COMMON /VTKB/ IVT,IKB,IC COMMON /RFLN/ IEND,II,IH,IK,IL,ANGLES(4),WIDTH,SPEED,FBG, & BG1,SCAN,BG2,YNET,SIGY,XTIME COMMON /BLOCK0/ LU1,LU2,LU3,ILP CHARACTER H1*1,H2*1,H3*1,A*1 SAVE IF (IC.EQ.1) GO TO 1 I1=0 I2=0 I3=0 H1=' ' H2=' ' H3=' ' 5 WRITE (IVT,5001) 5001 FORMAT (' ENTER H, K, L. ',$) READ (IKB,5002,ERR=5) I1,I2,I3 5002 FORMAT (3I5) IF (I1.EQ.0.AND.I2.EQ.0.AND.I3.EQ.0) RETURN C C ZONES HK0, H0L, OR 0KL? OR LEVELS HKn, HnL, OR nKL? C C AXES H00, 0K0, OR 00L? OR ROWS Hkl, hKl, OR hkL? C IF (I1.EQ.99) H1='H' IF (I2.EQ.99) H2='K' IF (I3.EQ.99) H3='L' C C 100? C IF (I1.EQ.1.AND.I2.EQ.0.AND.I3.EQ.0) THEN WRITE (IVT,1100) READ (IKB,1000) A IF (A.EQ.'Y') THEN H1='H' H2='0' H3='0' ELSE H1='1' H2='0' H3='0' END IF END IF 1000 FORMAT (A) 1100 FORMAT (' SEARCH FOR ALL AXIAL H00? Y OR N (N). ',$) C C 010? C IF (I1.EQ.0.AND.I2.EQ.1.AND.I3.EQ.0) THEN WRITE (IVT,1010) READ (IKB,1000) A IF (A.EQ.'Y') THEN H1='0' H2='K' H3='0' ELSE H1='0' H2='1' H3='0' END IF END IF 1010 FORMAT (' SEARCH FOR ALL AXIAL 0K0? Y OR N (N). ',$) C C 001? C IF (I1.EQ.0.AND.I2.EQ.0.AND.I3.EQ.1) THEN WRITE (IVT,1001) READ (IKB,1000) A IF (A.EQ.'Y') THEN H1='0' H2='0' H3='L' ELSE H1='0' H2='0' H3='1' END IF END IF 1001 FORMAT (' SEARCH FOR ALL AXIAL 00L? Y OR N (N). ',$) C C 110, HH0, OR HHL? C IF (I1.EQ.1.AND.I2.EQ.1.AND.I3.EQ.0) THEN WRITE (IVT,1110) READ (IKB,1000) A IF (A.EQ.'Y') THEN H1='H' H2='H' H3='L' ELSE WRITE (IVT,2110) READ (IKB,1000) A IF (A.EQ.'Y') THEN H1='H' H2='H' H3='0' END IF END IF END IF 1110 FORMAT (' SEARCH FOR ALL ZONAL HHL? Y OR N (N). ',$) 2110 FORMAT (' SEARCH FOR ALL AXIAL HH0? Y OR N (N). ',$) C C 101, H0H, OR HKH? C IF (I1.EQ.1.AND.I2.EQ.0.AND.I3.EQ.1) THEN WRITE (IVT,1101) READ (IKB,1000) A IF (A.EQ.'Y') THEN H1='H' H2='K' H3='H' ELSE WRITE (IVT,2101) READ (IKB,1000) A IF (A.EQ.'Y') THEN H1='H' H2='0' H3='H' END IF END IF END IF 1101 FORMAT (' SEARCH FOR ALL ZONAL HKH? Y OR N (N). ',$) 2101 FORMAT (' SEARCH FOR ALL AXIAL H0H? Y OR N (N). ',$) C C 011, 0KK, OR HKK? C IF (I1.EQ.0.AND.I2.EQ.1.AND.I3.EQ.1) THEN WRITE (IVT,1011) READ (IKB,1000) A IF (A.EQ.'Y') THEN H1='H' H2='K' H3='K' ELSE WRITE (IVT,2011) READ (IKB,1000) A IF (A.EQ.'Y') THEN H1='0' H2='K' H3='K' END IF END IF END IF 1011 FORMAT (' SEARCH FOR ALL ZONAL HKK? Y OR N (N). ',$) 2011 FORMAT (' SEARCH FOR ALL AXIAL 0KK? Y OR N (N). ',$) C C 111 OR HHH? C IF (I1.EQ.1.AND.I2.EQ.1.AND.I3.EQ.1) THEN WRITE (IVT,1111) READ (IKB,1000) A IF (A.EQ.'Y') THEN H1='H' H2='H' H3='H' ELSE H1='1' H2='1' H3='1' END IF END IF 1111 FORMAT (' SEARCH FOR ALL AXIAL HHH? Y OR N (N). ',$) C C SIGNED INDICES? STANDARD REFERENCE REFLECTIONS? C ISIGN=0 IOMIT=1 IF ((H1.NE.'0'.AND.H1.NE.'H').OR. & (H2.NE.'0'.AND.H2.NE.'H'.AND.H2.NE.'K').OR. & (H3.NE.'0'.AND.H3.NE.'H'.AND.H3.NE.'K'.AND.H3.NE.'L')) THEN WRITE (IVT,1003) READ (IKB,1000) A IF (A.EQ.'N') ISIGN=1 WRITE (IVT,1004) READ (IKB,1000) A IF (A.EQ.'N') IOMIT=0 END IF 1003 FORMAT (' IGNORE INDEX SIGNS AND SEARCH ABSOLUTE VALUES? Y OR N ( &Y). ',$) 1004 FORMAT (' OMIT OCCURRENCES AS STANDARD REFERENCE REFLECTIONS? Y O &R N (Y). ',$) IF (ISIGN.EQ.0) THEN IF (H1.NE.'H') I1=ABS(I1) IF (H2.NE.'H'.AND.H2.NE.'K') I2=ABS(I2) IF (H3.NE.'H'.AND.H3.NE.'K'.AND.H3.NE.'L') I3=ABS(I3) END IF I=0 C C ZONES HK0, H0L, OR 0KL? OR LEVELS HKn, HnL, nKL? C IF (H1.EQ.'H'.AND.H2.EQ.'K'.AND.H3.NE.'L') I=1 IF (H1.EQ.'H'.AND.H2.NE.'K'.AND.H3.EQ.'L') I=2 IF (H1.NE.'H'.AND.H2.EQ.'K'.AND.H3.EQ.'L') I=3 C C ZONES HHL, HKH, OR HKK? C IF (H1.EQ.'H'.AND.H2.EQ.'H'.AND.H3.EQ.'L') I=4 IF (H1.EQ.'H'.AND.H2.EQ.'K'.AND.H3.EQ.'H') I=5 IF (H1.EQ.'H'.AND.H2.EQ.'K'.AND.H3.EQ.'K') I=6 C C AXES H00, 0K0, OR 00L? OR ROWS Hkl, hKl, OR hkL? C IF (H1.EQ.'H'.AND.H2.NE.'K'.AND.H3.NE.'L') I=7 IF (H1.NE.'H'.AND.H2.EQ.'K'.AND.H3.NE.'L') I=8 IF (H1.NE.'H'.AND.H2.NE.'K'.AND.H3.EQ.'L') I=9 C C AXES HH0, H0H, 0KK, OR HHH? OR ROWS HHn, HnH, OR nKK? C IF (H1.EQ.'H'.AND.H2.EQ.'H'.AND.H3.NE.'L') I=10 IF (H1.EQ.'H'.AND.H2.NE.'K'.AND.H3.EQ.'H') I=11 IF (H1.NE.'H'.AND.H2.EQ.'K'.AND.H3.EQ.'K') I=12 IF (H1.EQ.'H'.AND.H2.EQ.'H'.AND.H3.EQ.'H') I=13 C C Q = I/SIGMA(I) SEARCH LIMITS? C CALL QLIMITS C C SEARCH. C 1 CALL READR IF (IEND.NE.0) RETURN IF (II.LT.0.AND.IOMIT.EQ.1) GO TO 1 J=IH K=IK L=IL IF (ISIGN.EQ.0) THEN J=ABS(J) K=ABS(K) L=ABS(L) END IF IF (I.EQ.0.AND.(J.NE.I1.OR.K.NE.I2.OR.L.NE.I3)) GO TO 1 IF (I.EQ.1.AND.L.NE.I3) GO TO 1 IF (I.EQ.2.AND.K.NE.I2) GO TO 1 IF (I.EQ.3.AND.J.NE.I1) GO TO 1 IF (I.EQ.4.AND.ABS(J).NE.ABS(K)) GO TO 1 IF (I.EQ.5.AND.ABS(J).NE.ABS(L)) GO TO 1 IF (I.EQ.6.AND.ABS(K).NE.ABS(L)) GO TO 1 IF (I.EQ.7.AND.(K.NE.I2.OR.L.NE.I3)) GO TO 1 IF (I.EQ.8.AND.(J.NE.I1.OR.L.NE.I3)) GO TO 1 IF (I.EQ.9.AND.(J.NE.I1.OR.K.NE.I2)) GO TO 1 IF (I.EQ.10.AND.(ABS(J).NE.ABS(K).OR.L.NE.I3)) GO TO 1 IF (I.EQ.11.AND.(ABS(J).NE.ABS(L).OR.K.NE.I2)) GO TO 1 IF (I.EQ.12.AND.(ABS(K).NE.ABS(L).OR.J.NE.I1)) GO TO 1 IF (I.EQ.13.AND.(ABS(J).NE.ABS(K).OR.ABS(K).NE.ABS(L))) GO TO 1 CALL QTEST (IQ) IF (IQ.EQ.0) GO TO 1 RETURN END C----------------------------------------------------------------------- SUBROUTINE ASEARCH COMMON /VTKB/ IVT,IKB,IC COMMON /RFLN/ IEND,II,IH,IK,IL,ANGLES(4),WIDTH,SPEED,FBG, & BG1,SCAN,BG2,YNET,SIGY,XTIME COMMON /BLOCK0/ LU1,LU2,LU3,ILP CHARACTER*18 H SAVE IF (IC.EQ.1) GO TO 1 5 WRITE (IVT,1001) 1001 FORMAT (' ENTER CONDITION FOR POSSIBLE REFLECTION (E.G., H0L,L+H=2 &N).') READ (IKB,1002,ERR=5) H 1002 FORMAT (A) IF (H.EQ.' ') RETURN I=0 CALL HKLCOND (H,I) IF (I.EQ.0) THEN WRITE (IVT,1005) GO TO 5 END IF 1005 FORMAT (' CONDITION MUST BE ONE OF THE FOLLOWING. NO LEADING OR E &MBEDDED BLANKS.'// &' HKL,H+K=2N HKL,K+L=2N HKL,L+H=2N'/ &' HKL,H+K,K+L,L+H=2N HKL,H+K+L=2N HKL,-H+K+L=3N'/ &' HKL,H-K+L=3N HKL,H-K=3N HK0,H=2N'/ &' HK0,K=2N HK0,H+K=2N HK0,H+K=4N'/ &' 0KL,K=2N 0KL,L=2N 0KL,K+L=2N'/ &' 0KL,K+L=4N H0L,L=2N H0L,H=2N'/ &' H0L,L+H=2N H0L,L+H=4N HH(-2H)L,L=2N'/ &' H(-H)0L,L=2N HHL,L=2N(R-AXES) HHL,L=2N'/ &' HKH,K=2N HKK,H=2N HHL,2H+L=4N'/ &' HKH,2H+K=4N HKK,2K+H=4N H00,H=2N'/ &' H00,H=4N 0K0,K=2N 0K0,K=4N'/ &' 00L,L=2N 00L,L=4N 000L,L=2N'/ &' 000L,L=3N 000L,L=6N') CALL QLIMITS 1 CALL READR IF (IEND.NE.0) RETURN IF (II.LT.0) GO TO 1 CALL HKLTEST (IH,IK,IL,1,I,IABSENT) IF (IABSENT.EQ.0) GO TO 1 CALL QTEST (IQ) IF (IQ.EQ.0) GO TO 1 RETURN END C----------------------------------------------------------------------- SUBROUTINE TSEARCH CHARACTER LABEL*2 COMMON /LABL/ LABEL(4) COMMON /VTKB/ IVT,IKB,IC COMMON /RFLN/ IEND,II,IH,IK,IL,ANGLES(4),WIDTH,SPEED,FBG, & BG1,SCAN,BG2,YNET,SIGY,XTIME COMMON /PLOT/ ITYPE,LU,NX,NY,ILOG,FSCALE,ISMOOTH COMMON /BLOCK0/ LU1,LU2,LU3,ILP SAVE DATA TOL /0.01/ IF (IC.EQ.1) GO TO 1 TMIN=0 TMAX=0 5 WRITE (IVT,1001) 1001 FORMAT (' ENTER TWO-THETA VALUE(S). TT OR TTMIN, TTMAX. ',$) READ (IKB,1002,ERR=5) TMIN,TMAX 1002 FORMAT (2F10.0) IF (TMIN.EQ.0.AND.TMAX.EQ.0) RETURN IF (TMAX.EQ.0) TMAX=TMIN IF (TMAX.LT.TMIN) THEN T=TMIN TMIN=TMAX TMAX=T END IF IF (LABEL(1).EQ.'TH') THEN TMIN=0.5*TMIN TMAX=0.5*TMAX END IF TMIN=TMIN-TOL TMAX=TMAX+TOL CALL QLIMITS 1 CALL READR IF (IEND.NE.0) RETURN IF (II.LT.0) GO TO 1 T=ANGLES(1) IF (T.LT.TMIN.OR.T.GT.TMAX) GO TO 1 CALL QTEST (IQ) IF (IQ.EQ.0) GO TO 1 RETURN END C----------------------------------------------------------------------- SUBROUTINE SSEARCH COMMON /VTKB/ IVT,IKB,IC COMMON /PLOT/ ITYPE,LU,NX,NY,ILOG,FSCALE,ISMOOTH COMMON /RFLN/ IEND,II,IH,IK,IL,ANGLES(4),WIDTH,SPEED,FBG, & BG1,SCAN,BG2,YNET,SIGY,XTIME COMMON /BLOCK0/ LU1,LU2,LU3,ILP COMMON /BLOCK2/ GINV(3,3),WLE,WL1,WL0,WL2 DIMENSION H(3) SAVE DATA TOL /0.001/ IF (ITYPE.EQ.0) THEN WRITE (IVT,1010) RETURN END IF 1010 FORMAT ('0THE SIN(THETA)/LAMBDA SEARCH REQUIRES THE METRIC MATRI', &'X FROM THE ''BGLP.DAT'''/' FILE FROM PROGRAM REFPK. WITH RAW PRO &FILE DATA YOU CAN DO A SEARCH ON THE'/' DIFFRACTOMETER TWO-THETA V &ALUES.') IF (IC.EQ.1) GO TO 1 SMIN=0 SMAX=0 5 WRITE (IVT,1001) 1001 FORMAT (' ENTER SIN(THETA)/LAMBDA VALUE(S). S OR SMIN, SMAX. ',$ &) READ (IKB,1002,ERR=5) SMIN,SMAX 1002 FORMAT (2F10.0) IF (SMIN.EQ.0.AND.SMAX.EQ.0) RETURN IF (SMAX.EQ.0) SMAX=SMIN SMIN=ABS(SMIN) SMAX=ABS(SMAX) IF (SMAX.LT.SMIN) THEN T=SMIN SMIN=SMAX SMAX=T END IF SMIN=SMIN-TOL SMAX=SMAX+TOL CALL QLIMITS 1 CALL READR IF (IEND.NE.0) RETURN IF (II.LT.0) GO TO 1 H(1)=IH H(2)=IK H(3)=IL S=SINTHL(H,GINV) IF (S.LT.SMIN.OR.S.GT.SMAX) GO TO 1 CALL QTEST (IQ) IF (IQ.EQ.0) GO TO 1 RETURN END C----------------------------------------------------------------------- SUBROUTINE QSEARCH COMMON /VTKB/ IVT,IKB,IC COMMON /RFLN/ IEND,II,IH,IK,IL,ANGLES(4),WIDTH,SPEED,FBG, & BG1,SCAN,BG2,YNET,SIGY,XTIME IF (IC.EQ.0) CALL QLIMITS 1 CALL READR IF (IEND.NE.0) RETURN IF (II.LT.0) GO TO 1 CALL QTEST (IQ) IF (IQ.EQ.0) GO TO 1 RETURN END C----------------------------------------------------------------------- SUBROUTINE QLIMITS COMMON /VTKB/ IVT,IKB,IC COMMON /RFLN/ IEND,II,IH,IK,IL,ANGLES(4),WIDTH,SPEED,FBG, & BG1,SCAN,BG2,YNET,SIGY,XTIME COMMON /PLOT/ ITYPE,LU,NX,NY,ILOG,FSCALE,ISMOOTH COMMON /QLIM/ QMIN,QMAX 5 WRITE (IVT,1001) 1001 FORMAT ('0ENTER MINIMUM AND MAXIMUM Q = I/SIGMA(I) LIMITS.'/' ENTE &R QMIN, QMAX (DEFAULT: NO LIMITS). ',$) READ (IKB,1002,ERR=5) QMIN,QMAX 1002 FORMAT (2F10.0) IF (QMIN.EQ.0.AND.QMAX.EQ.0) RETURN IF (QMAX.LT.QMIN) QMAX=9E9 IF (ITYPE.EQ.0) THEN 6 WRITE (IVT,1003) READ (IKB,1002,ERR=6) FBG IF (FBG.EQ.0) FBG=1/6.0 END IF 1003 FORMAT ('0ENTER THE DECIMAL FRACTION OF THE SCAN WIDTH TO BE TAKEN & AS BACKGROUND'/' AT EACH END OF THE SCAN.'/' ENTER F (DEFAULT: F & = 1/6 = 0.1667). ',$) RETURN END C----------------------------------------------------------------------- SUBROUTINE QTEST (IQ) PARAMETER (NMAX=100) COMMON /RFLN/ IEND,II,IH,IK,IL,ANGLES(4),WIDTH,SPEED,FBG, & BG1,SCAN,BG2,YNET,SIGY,XTIME COMMON /PLOT/ ITYPE,LU,NX,NY,ILOG,FSCALE,ISMOOTH COMMON /QLIM/ QMIN,QMAX COMMON /BLOCK7/ NSTEP,XX(NMAX),YY(NMAX),VX(NMAX),VY(NMAX) COMMON /BLOCK8/ TH0,X0,SIGX0,XE,X1,X2,W1,W2,I1,I2,L1,L2,R,SIGR,EPS Q=0 IQ=0 IF (QMIN.EQ.0.AND.QMAX.EQ.0) THEN CALL PLOT0 IQ=1 ELSE IF (ITYPE.EQ.0) THEN CALL BGPKBG (FBG,NSTEP,XX,YY,VX,VY,BG1,SCAN,BG2,YNET,SIGY) IF (SIGY.NE.0) Q=YNET/SIGY IF (QMIN.LE.Q.AND.Q.LE.QMAX) THEN CALL PLOT0 IQ=1 END IF ELSE CALL PROCESS IF (SIGR.NE.0) Q=R/SIGR IF (QMIN.LE.Q.AND.Q.LE.QMAX) THEN CALL PLOT1 IQ=1 END IF END IF RETURN END C----------------------------------------------------------------------- SUBROUTINE PSEARCH PARAMETER (NMAX=100) COMMON /VTKB/ IVT,IKB,IC COMMON /RFLN/ IEND,II,IH,IK,IL,ANGLES(4),WIDTH,SPEED,FBG, & BG1,SCAN,BG2,YNET,SIGY,XTIME COMMON /PLOT/ ITYPE,LU,NX,NY,ILOG,FSCALE,ISMOOTH COMMON /BLOCK0/ LU1,LU2,LU3,ILP COMMON /BLOCK3/ WINDOW,CUTOFF,NSMOOTH,II1,II2,XT1,XT2,SC1,SC2, & XX0,XW1,XW2 COMMON /BLOCK7/ NSTEP,XX(NMAX),YY(NMAX),VX(NMAX),VY(NMAX) COMMON /BLOCK8/ TH0,X0,SIGX0,XE,X1,X2,W1,W2,I1,I2,L1,L2,R,SIGR,EPS SAVE IF (IC.EQ.1) GO TO 1 DMIN=0 DMAX=0 5 WRITE (IVT,1001) 1001 FORMAT ('0ENTER A SEARCH RANGE FOR PEAK DISPLACEMENTS FROM THE SCA &N MID-POINT AS'/' DECIMAL FRACTIONS OF THE SCAN WIDTH.'/' ENTER DM &IN, DMAX (DEFAULTS: DMIN = 0.1, DMAX = 0.5). ',$) READ (IKB,1002,ERR=5) DMIN,DMAX 1002 FORMAT (2F10.0) IF (DMIN.EQ.0) DMIN=0.1 IF (DMAX.EQ.0) DMAX=0.5 1 CALL READR IF (IEND.NE.0) RETURN IF (ITYPE.EQ.0) THEN CALL CENTROID (WINDOW,CUTOFF,XX,YY,VX,VY,1,NSTEP,X0,SIGX0) ELSE CALL PROCESS END IF D=ABS(X0-0.5*(XX(1)+XX(NSTEP)))/(XX(NSTEP)-XX(1)) IF (D.LT.DMIN.OR.D.GT.DMAX) GO TO 1 IF (ITYPE.EQ.0) THEN CALL PLOT0 ELSE CALL PLOT1 END IF RETURN END C----------------------------------------------------------------------- SUBROUTINE BSEARCH PARAMETER (NMAX=100) COMMON /VTKB/ IVT,IKB,IC COMMON /RFLN/ IEND,II,IH,IK,IL,ANGLES(4),WIDTH,SPEED,FBG, & BG1,SCAN,BG2,YNET,SIGY,XTIME COMMON /PLOT/ ITYPE,LU,NX,NY,ILOG,FSCALE,ISMOOTH COMMON /BLOCK0/ LU1,LU2,LU3,ILP COMMON /BLOCK3/ WINDOW,CUTOFF,NSMOOTH,II1,II2,XT1,XT2,SC1,SC2, & XX0,XW1,XW2 COMMON /BLOCK7/ NSTEP,XX(NMAX),YY(NMAX),VX(NMAX),VY(NMAX) COMMON /BLOCK8/ TH0,X0,SIGX0,XE,X1,X2,W1,W2,I1,I2,L1,L2,R,SIGR,EPS DIMENSION X(NMAX),Y(NMAX) SAVE IF (IC.EQ.1) GO TO 1 5 WRITE (IVT,1001) 1001 FORMAT ('0ENTER A DECIMAL FRACTION OF THE SCAN WIDTH FOR THE LOCAL & BACKGROUND'/' SAMPLES.'/' ENTER F (DEFAULT: F = 0.0833 = 1/12)', &'. ',$) READ (IKB,1002,ERR=5) F 1002 FORMAT (F10.0) IF (F.EQ.0) F=1/12.0 6 WRITE (IVT,1003) 1003 FORMAT ('0ENTER A SIGNIFICANCE THRESHOLD (Q = DELTA/SIGMA(DELTA)', &') FOR LOCAL'/' COUNT RATE DIFFERENCES.'/' ENTER QMIN (DEFAULT: Q &MIN = 3). ',$) READ (IKB,1002,ERR=6) Q IF (Q.LE.0) Q=3 7 WRITE (IVT,1004) 1004 FORMAT ('0ENTER A MAXIMUM PERMITTED RATIO [R = Y(1)/Y(2) OR Y(N)/Y &(N-1)]'/' FOR ADJACENT BACKGROUND SAMPLES AT THE TWO ENDS OF THE S &CAN.'/' ENTER RMAX (DEFAULT: RMAX = 2). ',$) READ (IKB,1002,ERR=7) R IF (R.EQ.0) R=2 1 CALL READR IF (IEND.NE.0) RETURN N=NINT(F*NSTEP) Y1=0 Y2=0 V1=0 V2=0 DO I=1,N Y1=Y1+YY(I) Y2=Y2+YY(I+N) V1=V1+VY(I) V2=V2+VY(I+N) END DO IF (Y1-Y2.GT.Q*SQRT(V1+V2).AND.Y1.GT.R*Y2) GO TO 2 Y1=0 Y2=0 V1=0 V2=0 DO I=NSTEP-N+1,NSTEP Y1=Y1+YY(I) Y2=Y2+YY(I+N) V1=V1+VY(I) V2=V2+VY(I+N) END DO IF (Y2-Y1.GT.Q*SQRT(V1+V2).AND.Y2.GT.R*Y1) GO TO 2 GO TO 1 2 CALL PLOT0 RETURN END C----------------------------------------------------------------------- SUBROUTINE ESEARCH PARAMETER (NMAX=100) COMMON /VTKB/ IVT,IKB,IC COMMON /RFLN/ IEND,II,IH,IK,IL,ANGLES(4),WIDTH,SPEED,FBG, & BG1,SCAN,BG2,YNET,SIGY,XTIME COMMON /BLOCK0/ LU1,LU2,LU3,ILP COMMON /BLOCK7/ NSTEP,XX(NMAX),YY(NMAX),VX(NMAX),VY(NMAX) SAVE IF (IC.EQ.1) GO TO 1 5 WRITE (IVT,1001) 1001 FORMAT ('0ENTER A SIGNIFICANCE THRESHOLD (Q = ABS(DELTA)/SIGMA(DEL &TA)) FOR LOCAL'/' COUNT RATE DIFFERENCES.'/' ENTER QMIN (DEFAULT', &': QMIN = 3). ',$) READ (IKB,1002,ERR=5) Q 1002 FORMAT (F10.0) IF (Q.EQ.0) Q=3 6 WRITE (IVT,1003) 1003 FORMAT ('0ENTER A THRESHOLD FACTOR [F = ABS(Y(I) - Y(I-1))/(YMAX - & YMIN)] FOR'/' TESTING FOR SCAN DISCONTINUITY ERRORS.'/' ENTER FMI &N (DEFAULT: FMIN = 0.75). ',$) READ (IKB,1002,ERR=6) F IF (F.EQ.0) F=0.75 1 CALL READR IF (IEND.NE.0) RETURN C C TEST FOR FLAT OR BLANK SCAN. C YMIN=YY(1) YMAX=YY(1) DO I=2,NSTEP IF (YY(I).LT.YMIN) THEN YMIN=YY(I) VMIN=VY(I) END IF IF (YY(I).GT.YMAX) THEN YMAX=YY(I) VMAX=VY(I) END IF END DO IF (YMAX-YMIN.LE.SQRT(VMIN+VMAX)) THEN CALL PLOT0 RETURN END IF C C TEST FOR STATISTICALLY FLAT SCAN. C SUMY=0 SUMV=0 DO I=1,NSTEP SUMY=SUMY+YY(I) SUMV=SUMV+VY(I) END DO YMEAN=SUMY/NSTEP ESD=SQRT(SUMV/NSTEP) SUMV=0 DO I=1,NSTEP SUMV=SUMV+(YY(I)-YMEAN)**2 END DO RMSD=SQRT(SUMV/(NSTEP-1)) IF (RMSD.LT.Q*ESD) GO TO 1 VMEAN=RMSD**2/NSTEP C C A SCAN WITH A SIGNIFICANT PEAK SHOULD HAVE A SIGNIFICANTLY SKEWED C DISTRIBUTION OF COUNT RATES. C SUMY=0 SUMV=0 DO I=1,NSTEP D=YY(I)-YMEAN SUMY=SUMY+D**3 SUMV=SUMV+D**4*(VY(I)+VMEAN) END DO SKEW=SUMY/NSTEP ESD=3*SQRT(SUMV)/NSTEP IF (ABS(SKEW).LT.Q*ESD) GO TO 1 C C TEST FOR DISCONTINUITY. C T=F*(YMAX-YMIN) DO I=2,NSTEP D=ABS(YY(I)-YY(I-1)) S=SQRT(VY(I)+VY(I-1)) IF (D.GT.Q*S.AND.D.GT.T) THEN CALL PLOT0 RETURN END IF END DO GO TO 1 END C----------------------------------------------------------------------- SUBROUTINE PLOT0 C C CALLED FROM DRIVER AND SEARCH SUBROUTINES. C PARAMETER (NMAX=100) DIMENSION XTEMP(NMAX),YTEMP(NMAX),A(100,50) COMMON /BLOCK7/ NSTEP,XX(NMAX),YY(NMAX),VX(NMAX),VY(NMAX) COMMON /PLOT/ ITYPE,LU,NX,NY,ILOG,FSCALE,ISMOOTH COMMON /RFLN/ IEND,II,IH,IK,IL,ANGLES(4),WIDTH,SPEED,FBG, & BG1,SCAN,BG2,YNET,SIGY,XTIME IF (ITYPE.EQ.0) THEN C C RAW PROFILE DATA C CALL BGPKBG (FBG,NSTEP,XX,YY,VX,VY,BG1,SCAN,BG2,YNET,SIGY) CALL ZDATA (NMAX,NSTEP,XX,YY,NX,XTEMP,YTEMP,ISMOOTH) CALL ADATA (ILOG,FSCALE,NX,NY,XTEMP,YTEMP,A,0,0,0,0,0) CALL RPLOT (A,NSTEP) ELSE C C PROCESSED PROFILE DATA C CALL PROCESS CALL PLOT1 END IF RETURN END C----------------------------------------------------------------------- SUBROUTINE PLOT1 C C CALLED FROM SUBROUTINE DRIVER, PLOT0, QTEST, MOVE, OR MODIFY. C PARAMETER (NMAX=100) DIMENSION XTEMP(NMAX),YTEMP(NMAX),A(100,50) COMMON /BLOCK2/ GINV(3,3),WLE,WL1,WL0,WL2 COMMON /BLOCK7/ NSTEP,XX(NMAX),YY(NMAX),VX(NMAX),VY(NMAX) COMMON /BLOCK8/ TH0,X0,SIGX0,XE,X1,X2,W1,W2,I1,I2,L1,L2,R,SIGR,EPS COMMON /PLOT/ ITYPE,LU,NX,NY,ILOG,FSCALE,ISMOOTH CALL ZDATA (NMAX,NSTEP,XX,YY,NX,XTEMP,YTEMP,ISMOOTH) CALL ADATA (ILOG,FSCALE,NX,NY,XTEMP,YTEMP,A,I1,I2,L1,L2,NSTEP) IF (I1.EQ.0.AND.I2.EQ.0.AND.L1.EQ.0.AND.L2.EQ.0) THEN I0=0 IE=0 ELSE I0=NINDEX(X0) IF (WLE.EQ.0) THEN IE=0 ELSE IE=I0+NINT(NSTEP*(XE-X0)/(XX(NSTEP)-XX(1))) IE=MAX(IE,-99) END IF END IF CALL PPLOT (A,NSTEP,I1,IE,L1,I0,L2,I2,R,SIGR) RETURN END C----------------------------------------------------------------------- SUBROUTINE ZDATA (M,N,X,Y,NN,XX,YY,NSMOOTH) DIMENSION X(1),Y(1),XX(1),YY(1) IF (NN.LE.NINT(N/3.0)) THEN DO 3 I=1,NN J=3*I XX(I)=(X(J-2)+X(J-1)+X(J))/3 YY(I)=(Y(J-2)+Y(J-1)+Y(J))/3 3 CONTINUE ELSE IF (NN.LE.NINT(N/2.0)) THEN DO 2 I=1,NN J=2*I XX(I)=(X(J-1)+X(J))/2 YY(I)=(Y(J-1)+Y(J))/2 2 CONTINUE ELSE DO 1 I=1,MIN(N,NN) XX(I)=X(I) YY(I)=Y(I) 1 CONTINUE END IF J=MIN(N,NN) IF (NSMOOTH.GT.0) CALL SMOOTH (NSMOOTH,J,YY) IF (J.LT.M) THEN DO 9 I=J+1,M XX(I)=0 YY(I)=0 9 CONTINUE END IF RETURN END C----------------------------------------------------------------------- SUBROUTINE ADATA (ILOG,FSCALE,M,N,X,Y,A,I1,I2,L1,L2,NSTEP) C C BUILD THE DISPLAY ARRAY FOR A REFLECTION PROFILE. C CHARACTER A*1 DIMENSION A(100,50),X(1),Y(1) COMMON /BG/ A0,A1,E0,E1,B0,B1,J1,J2 YMAX=1 DO 10 I=1,M YMAX=MAX(YMAX,Y(I)) 10 CONTINUE SCALE=FSCALE*N/YMAX IF (ILOG.NE.0) SCALE=N/LOG(1+YMAX) DO 1 I=1,M YI=Y(I) IF (ILOG.NE.0) YI=LOG(1+YI) YI=NINT(SCALE*YI) DO 1 J=1,N A(I,J)=' ' IF (YI.GE.N-J+1) A(I,J)='X' 1 CONTINUE C C DUMMY ARGUMENTS I1 = I2 = L1 = L2 = 0 FLAG A RAW PROFILE. C IF (I1.EQ.0.AND.I2.EQ.0.AND.L1.EQ.0.AND.L2.EQ.0) RETURN DO 2 I=1,M IF (I1.LE.I.AND.I.LE.I2) THEN XI=X(I) IF (E0.EQ.0) THEN BI=A0+A1*XI ELSE IF (I.LT.J1) BI=A0+A1*XI IF (J1.LE.I.AND.I.LT.L1) BI=0 IF (L1.LE.I.AND.I.LE.J2) BI=E0+E1*XI IF (I.GT.J2) BI=B0+B1*XI END IF IF (BI.LT.0) BI=0 IF (ILOG.NE.0) BI=LOG(1+BI) J=NINT(SCALE*BI) IF (1.LE.J.AND.J.LE.N) THEN J=N-J+1 IF (A(I,J).NE.' ') A(I,J)='*' END IF END IF 2 CONTINUE I=MIN(M,NSTEP) DO 3 J=1,N IF (I1.GT.1.AND.A(I1,J).NE.' ') A(I1,J)='*' IF (I2.LT.I.AND.A(I2,J).NE.' ') A(I2,J)='*' IF (L1.GE.1.AND.A(L1,J).NE.' ') A(L1,J)='*' IF (L2.LE.I.AND.A(L2,J).NE.' ') A(L2,J)='*' 3 CONTINUE RETURN END C----------------------------------------------------------------------- SUBROUTINE RPLOT (A,NSTEP) C C PLOT A RAW REFLECTION PROFILE AND THE DIFFRACTOMETER DATA. C CHARACTER A*1 DIMENSION A(100,50) CHARACTER LABEL*2 COMMON /LABL/ LABEL(4) COMMON /PLOT/ ITYPE,LU,NX,NY,ILOG,FSCALE,ISMOOTH COMMON /RFLN/ IEND,II,IH,IK,IL,ANGLES(4),WIDTH,SPEED,FBG, & BG1,SCAN,BG2,YNET,SIGY,XTIME WRITE (LU,9000) 9000 FORMAT (////) IF (ILOG.NE.0) THEN WRITE (LU,3001) II, (A(I,1),I=1,NX) ELSE IF (FSCALE.NE.1) THEN WRITE (LU,2001) II,NINT(FSCALE),(A(I,1),I=1,NX) ELSE WRITE (LU,1001) II, (A(I,1),I=1,NX) END IF IF (ISMOOTH.NE.0) THEN WRITE (LU,2000) ISMOOTH,(A(I,2),I=1,NX) ELSE WRITE (LU,1000) (A(I,2),I=1,NX) END IF WRITE (LU,1002) IH,IK,IL,(A(I,3),I=1,NX) WRITE (LU,1000) (A(I,4),I=1,NX) WRITE (LU,1003) LABEL(1),ANGLES(1),(A(I,5),I=1,NX) WRITE (LU,1003) LABEL(2),ANGLES(2),(A(I,6),I=1,NX) WRITE (LU,1003) LABEL(3),ANGLES(3),(A(I,7),I=1,NX) WRITE (LU,1003) LABEL(4),ANGLES(4),(A(I,8),I=1,NX) WRITE (LU,1000) (A(I, 9),I=1,NX) WRITE (LU,1004) WIDTH,(A(I,10),I=1,NX) WRITE (LU,1005) SPEED,(A(I,11),I=1,NX) WRITE (LU,1006) BG1, (A(I,12),I=1,NX) WRITE (LU,1007) SCAN, (A(I,13),I=1,NX) WRITE (LU,1008) BG2, (A(I,14),I=1,NX) WRITE (LU,1000) (A(I,15),I=1,NX) WRITE (LU,1009) YNET, (A(I,16),I=1,NX) WRITE (LU,1010) SIGY, (A(I,17),I=1,NX) WRITE (LU,1011) XTIME,(A(I,18),I=1,NX) WRITE (LU,1000) (A(I,19),I=1,NX) WRITE (LU,1000) (A(I,20),I=1,NX) IF (NY.GT.20) THEN DO 1 J=21,NY WRITE (LU,1000) (A(I,J),I=1,NX) 1 CONTINUE END IF IF (NX.LE.50) THEN IF (NSTEP.LE.50) THEN WRITE (LU,1051) ELSE WRITE (LU,1052) END IF ELSE WRITE (LU,1100) END IF RETURN 1000 FORMAT (' ',20X,100A1) 1001 FORMAT (' I ',I10,5X,100A1) 2001 FORMAT (' I ',I10,I4,'x',100A1) 3001 FORMAT (' I ',I10,' LN Y',100A1) 2000 FORMAT (' ',10X,I5,'x SM ',100A1) 1002 FORMAT (' HKL ',3I3,5X,100A1) 1003 FORMAT (' ',A2,' ',F10.2,5X,100A1) 1004 FORMAT (' WIDTH',F10.2,5X,100A1) 1005 FORMAT (' SPEED',F10.2,5X,100A1) 1006 FORMAT (' BG1 ',F10.0,5X,100A1) 1007 FORMAT (' SCAN ',F10.0,5X,100A1) 1008 FORMAT (' BG2 ',F10.0,5X,100A1) 1009 FORMAT (' YNET ',F10.2,5X,100A1) 1010 FORMAT (' SIGY ',F10.2,5X,100A1) 1011 FORMAT (' XTIME',F10.2,5X,100A1) 1051 FORMAT (' ',19X,'+----+----+----+----+----+----+----+----+----+--- &-+'/' ',19X,'0 10 20 30 40 50') 1052 FORMAT (' ',19X,'+----+----+----+----+----+----+----+----+----+--- &-+'/' ',19X,'0 10 20 30 40 50 60 70 80 90 100') 1100 FORMAT (' ',19X,'+----+----+----+----+----+----+----+----+----+--- &-+----+----+----+----+----+----+----+----+----+----+'/' ',19X, &'0 10 20 30 40 50 60', &' 70 80 90 100') END C----------------------------------------------------------------------- SUBROUTINE PPLOT (A,NSTEP,I1,IE,L1,I0,L2,I2,R,SIGR) C C PLOT A PROCESSED REFLECTION PROFILE. C CHARACTER A*1 DIMENSION A(100,50) CHARACTER LABEL*2 COMMON /LABL/ LABEL(4) COMMON /PLOT/ ITYPE,LU,NX,NY,ILOG,FSCALE,ISMOOTH COMMON /RFLN/ IEND,II,IH,IK,IL,ANGLES(4),WIDTH,SPEED,FBG, & BG1,SCAN,BG2,YNET,SIGY,XTIME WRITE (LU,9000) 9000 FORMAT (////) IF (I1.EQ.0.AND.L1.EQ.0.AND.L2.EQ.0.AND.I2.EQ.0) WRITE (LU,9001) 9001 FORMAT ('0REFLN. OUTSIDE ''BGLP.DAT'' SER. NO. OR XTIME RANGE, O', &'R SYMMETRY FORBIDDEN.') IF (ILOG.NE.0) THEN WRITE (LU,3001) II, (A(I,1),I=1,NX) ELSE IF (FSCALE.NE.1) THEN WRITE (LU,2001) II,NINT(FSCALE),(A(I,1),I=1,NX) ELSE WRITE (LU,1001) II, (A(I,1),I=1,NX) END IF IF (ISMOOTH.NE.0) THEN WRITE (LU,2000) ISMOOTH,(A(I,2),I=1,NX) ELSE WRITE (LU,1000) (A(I,2),I=1,NX) END IF WRITE (LU,1002) IH,IK,IL,(A(I,3),I=1,NX) WRITE (LU,1000) (A(I,4),I=1,NX) WRITE (LU,1003) LABEL(1),ANGLES(1),(A(I,5),I=1,NX) WRITE (LU,1003) LABEL(2),ANGLES(2),(A(I,6),I=1,NX) WRITE (LU,1003) LABEL(3),ANGLES(3),(A(I,7),I=1,NX) WRITE (LU,1003) LABEL(4),ANGLES(4),(A(I,8),I=1,NX) WRITE (LU,1000) (A(I, 9),I=1,NX) WRITE (LU,1004) WIDTH,I1,(A(I,10),I=1,NX) WRITE (LU,1005) SPEED,IE,(A(I,11),I=1,NX) WRITE (LU,1006) L1, (A(I,12),I=1,NX) WRITE (LU,1007) I0, (A(I,13),I=1,NX) WRITE (LU,1008) L2, (A(I,14),I=1,NX) WRITE (LU,1009) I2, (A(I,15),I=1,NX) WRITE (LU,1000) (A(I,16),I=1,NX) WRITE (LU,1010) R, (A(I,17),I=1,NX) WRITE (LU,1011) SIGR, (A(I,18),I=1,NX) WRITE (LU,1012) XTIME,(A(I,19),I=1,NX) WRITE (LU,1000) (A(I,20),I=1,NX) IF (NY.GT.20) THEN DO 1 J=21,NY WRITE (LU,1000) (A(I,J),I=1,NX) 1 CONTINUE END IF IF (NX.LE.50) THEN IF (NSTEP.LE.50) THEN WRITE (LU,1051) ELSE WRITE (LU,1052) END IF ELSE WRITE (LU,1100) END IF DO 2 I=1,NX A(I,1)=' ' 2 CONTINUE IF (I1.GT. 1) A(I1,1)='1' IF (I2.LT.NX) A(I2,1)='2' IF (L1.GE. 1) A(L1,1)='1' IF (L2.LE.NX) A(L2,1)='2' IF (IE.GT. 0) A(IE,1)='E' A(I0,1)='0' WRITE (LU,1000) (A(I,1),I=1,NX) RETURN 1000 FORMAT (' ',20X,100A1) 1001 FORMAT (' I ',I8,7X,100A1) 2001 FORMAT (' I ',I8,I5,'x ',100A1) 3001 FORMAT (' I ',I8,' LN Y ',100A1) 2000 FORMAT (' ',13X,I3,'x SM',100A1) 1002 FORMAT (' HKL ',3I3,7X,100A1) 1003 FORMAT (' ',A2,' ',F8.2,7X,100A1) 1004 FORMAT (' WIDTH',F8.2,' I1',I3,1X,100A1) 1005 FORMAT (' SPEED',F8.2,' IE',I3,1X,100A1) 1006 FORMAT (' ',13X,' L1',I3,1X,100A1) 1007 FORMAT (' ',13X,' I0',I3,1X,100A1) 1008 FORMAT (' ',13X,' L2',I3,1X,100A1) 1009 FORMAT (' ',13X,' I2',I3,1X,100A1) 1010 FORMAT (' YNET ',3X,E10.3,2X,100A1) 1011 FORMAT (' SIGY ',3X,E10.3,2X,100A1) 1012 FORMAT (' XTIME',F8.2,7X,100A1) 1051 FORMAT (' ',19X,'+----+----+----+----+----+----+----+----+----+--- &-+'/' ',19X,'0 10 20 30 40 50') 1052 FORMAT (' ',19X,'+----+----+----+----+----+----+----+----+----+--- &-+'/' ',19X,'0 10 20 30 40 50 60 70 80 90 100') 1100 FORMAT (' ',19X,'+----+----+----+----+----+----+----+----+----+--- &-+----+----+----+----+----+----+----+----+----+----+'/' ',19X, &'0 10 20 30 40 50 60', &' 70 80 90 100') END C----------------------------------------------------------------------- SUBROUTINE MOVE PARAMETER (NMAX=100) COMMON /VTKB/ IVT,IKB,IC COMMON /BLOCK7/ NSTEP,XX(NMAX),YY(NMAX),VX(NMAX),VY(NMAX) COMMON /BLOCK8/ TH0,X0,SIGX0,XE,X1,X2,W1,W2,I1,I2,L1,L2,R,SIGR,EPS CHARACTER POINT*2 DIMENSION POINT(6),DELTA(6) I1OLD=I1 I2OLD=I2 X0OLD=X0 XEOLD=XE L1OLD=L1 L2OLD=L2 DO 1 I=1,6 POINT(I)=' ' DELTA(I)=0 1 CONTINUE 5 WRITE (IVT,1000) 1000 FORMAT (' ENTER POINT LABELS AND SHIFTS, E.G., IE-3,L1-5,L2+2. ', &$) READ (IKB,1001,ERR=5) (POINT(I),DELTA(I),I=1,6) 1001 FORMAT (6(A2,F4.0)) IF (DELTA(1).EQ.0) RETURN DX=(XX(I2)-XX(I1))/(I2-I1+1) DO 2 I=1,6 IF (POINT(I).EQ.'I1') I1=I1+DELTA(I) IF (POINT(I).EQ.'I2') I2=I2+DELTA(I) IF (POINT(I).EQ.'I0') THEN XE=XE+DELTA(I)*DX L1=L1+DELTA(I) L2=L2+DELTA(I) X0=X0+DELTA(I)*DX GO TO 3 END IF IF (POINT(I).EQ.'IE') XE=XE+DELTA(I)*DX IF (POINT(I).EQ.'L1') L1=L1+DELTA(I) IF (POINT(I).EQ.'L2') L2=L2+DELTA(I) 2 CONTINUE 3 CONTINUE IF (I1.LT.1.OR.I2.GT.NSTEP.OR. & L1.LT.MIN(L1OLD,I1).OR.L2.GT.MAX(L2OLD,I2).OR. & X0.LT.XX(I1).OR.X0.GT.XX(I2).OR.XE.GT.X0) THEN WRITE (IVT,1009) I1=I1OLD I2=I2OLD X0=X0OLD XE=XEOLD L1=L1OLD L2=L2OLD CALL PLOT1 RETURN END IF 1009 FORMAT ('0SHIFT OUT OF SCAN RANGE. RE-ENTER MOVE COMMAND AFTER '' &OPTION?'' PROMPT.') IF (WLE.EQ.0) THEN CALL BGLP1 ELSE EPS=0.5*(W1+W2) J1=NINDEX(XE-EPS) J2=NINDEX(XE+EPS) IF (J2.LT.L1-1) THEN I1=MAX(I1,J2) CALL BGLP1 ELSE CALL BGLP2 END IF END IF CALL PLOT1 RETURN END C----------------------------------------------------------------------- SUBROUTINE MODIFY PARAMETER (NMAX=100) COMMON /VTKB/ IVT,IKB,IC COMMON /BLOCK2/ GINV(3,3),WLE,WL1,WL0,WL2 COMMON /BLOCK7/ NSTEP,XX(NMAX),YY(NMAX),VX(NMAX),VY(NMAX) COMMON /BLOCK8/ TH0,X0,SIGX0,XE,X1,X2,W1,W2,I1,I2,L1,L2,R,SIGR,EPS DIMENSION K(5),L(5) DO 10 I=1,5 K(I)=0 L(I)=0 10 CONTINUE 5 WRITE (IVT,1000) 1000 FORMAT (' ENTER STEP NUMBERS FOR FIVE OR FEWER STEPS TO BE REPLACE &D BY A LOCAL'/' AVERAGE COUNT RATE.'/' ENTER I1, I2, I3, I4, I5.', &' ',$) READ (IKB,1001,ERR=5) K 1001 FORMAT (5I5) IF (K(1).LE.0) RETURN 6 WRITE (IVT,1002) 1002 FORMAT (' ENTER STEP NUMBERS FOR FIVE OR FEWER STEPS FOR CALCULATI &NG A LOCAL'/' AVERAGE REPLACEMENT COUNT RATE.'/' ENTER J1, J2, J3, & J4, J5. ',$) READ (IKB,1001,ERR=6) L IF (L(1).LE.0) RETURN DO 15 I=1,5 IF (K(I).LT.0.OR.K(I).GT.NSTEP) GO TO 5 IF (L(I).LT.0.OR.L(I).GT.NSTEP) GO TO 6 15 CONTINUE N=0 Y=0 V=0 DO 1 I=1,5 IF (L(I).NE.0) THEN N=N+1 Y=Y+YY(L(I)) V=V+VY(L(I)) END IF 1 CONTINUE Y=Y/N V=V/N DO 2 I=1,5 IF (K(I).NE.0) THEN YY(K(I))=Y VY(K(I))=V END IF 2 CONTINUE IF (WLE.EQ.0) THEN CALL BGLP1 ELSE EPS=0.5*(W1+W2) J1=NINDEX(XE-EPS) J2=NINDEX(XE+EPS) IF (J2.LT.L1-1) THEN I1=MAX(I1,J2) CALL BGLP1 ELSE CALL BGLP2 END IF END IF CALL PLOT1 RETURN END C----------------------------------------------------------------------- SUBROUTINE REVISE CHARACTER ATIME*8,ADATE*9,TITLE*80,FILE1*80,FILE2*80 COMMON /VTKB/ IVT,IKB,IC COMMON /BLOCK0/ LU1,LU2,LU3,ILP COMMON /BLOCK1/ ATIME,ADATE,TITLE,FILE1,FILE2 PARAMETER (M=100) DIMENSION JJ(M),JH(3,M),AJ(4,M),YJ(M),SJ(M),TJ(M),IH(3),AI(4) C C TEST FOR EMPTY FILE OF REVISED REFLECTIONS. C REWIND LU3 READ (LU3,END=21) GO TO 22 21 WRITE (IVT,6000) 6000 FORMAT ('0NO REVISED REFLECTIONS WERE RECORDED.') CLOSE (UNIT=LU3,STATUS='DELETE') RETURN 22 CONTINUE C C TEST FOR PRESENCE OF OLD 'DATA.BLP' FILE. C CLOSE (UNIT=LU1,STATUS='KEEP') OPEN (UNIT=LU1,FILE=FILE2,STATUS='OLD',FORM='UNFORMATTED',ERR=23) GO TO 24 23 CONTINUE C C IF THERE IS NO OLD 'DATA.BLP' FILE, WRITE A NEW 'DATA.BLP' FILE, C OMITTING ANY DELETED REFLECTIONS. C WRITE (IVT,6001) FILE2 6001 FORMAT ('0THERE IS NO OLD ''DATA.BLP'' FILE: ',A/' A NEW '' &DATA.BLP'' FILE WILL BE WRITTEN.') WRITE (ILP,6010) ATIME,ADATE,FILE1,FILE2 6010 FORMAT ('1'/'0PROGRAM VIEW/BGLP. ',A,1X,A/'0DATA.RAW FILE: ',A/' &0NEW DATA.BLP FILE: ',A/' --------'/'0 I H K L & A1 A2 A3 A4 Y SIGY XTIME'/ &' - - - - -- -- -- -- -', &' ---- -----') OPEN (UNIT=LU1,FILE=FILE2,STATUS='NEW',FORM='UNFORMATTED') REWIND LU3 51 READ (LU3,END=59) II,IH,AI,YI,SI,TI C C REFLECTIONS TO BE DELETED ARE FLAGGED BY A NEGATIVE SIGMA(Y). C IF (SI.GT.0) THEN WRITE (LU1) II,IH,AI,YI,SI,TI WRITE (ILP,6002) II,IH,AI,YI,SI,TI ELSE WRITE (ILP,6003) II,IH,AI END IF 6002 FORMAT ('0NEW ',I6,2X,3I3,2X,4F7.2,3F10.2) 6003 FORMAT ('0MEASUREMENT DELETED IN ''NEW.BLP'' FILE:'/' OLD ',I6, &2X,3I3,2X,4F7.2) GO TO 51 59 CONTINUE CLOSE (UNIT=LU1,STATUS='KEEP') CLOSE (UNIT=LU3,STATUS='DELETE') RETURN 24 CONTINUE C C OPEN NEW 'DATA.BLP' FILE. C OPEN (UNIT=LU2,FILE='revised.blp',STATUS='NEW',FORM='UNFORMATTED') WRITE (ILP,1000) ATIME,ADATE,FILE1,FILE2,FILE2 1000 FORMAT ('1'/'0PROGRAM VIEW/BGLP. ',A,1X,A/'0RAW DATA FILE: ',A/ &'0OLD.BLP FILE: ',A/'0NEW.BLP FILE: ',A/' -------'/'0 ', &' I H K L A1 A2 A3 A4 Y SIGY & XTIME'/' - - - - -- -- -- -- & - ---- -----') C C STORE REVISED REFLECTION DATA IN MEMORY. C REWIND LU3 11 N=0 DO 10 J=1,M READ (LU3,END=15) JJ(J),(JH(I,J),I=1,3),(AJ(I,J),I=1,4), & YJ(J),SJ(J),TJ(J) N=N+1 10 CONTINUE 15 IF (N.EQ.0) GO TO 9 C C READ OLD AND WRITE NEW 'DATA.BLP' FILES. C REWIND LU1 1 READ (LU1,END=8) II,IH,AI,YI,SI,TI C C WORK BACKWARDS THROUGH THE LIST OF MODIFIED REFLECTIONS IN ORDER TO C SELECT THE LAST MODIFICATION IN ANY CASE(S) MODIFIED MORE THAN ONCE. C DO 5 J=N,1,-1 IF (JJ(J).NE.II) GO TO 5 IF (TJ(J).NE.TI) GO TO 5 DO 3 I=1,3 IF (JH(I,J).NE.IH(I)) GO TO 5 3 CONTINUE DO 4 I=1,4 IF (AJ(I,J).NE.AI(I)) GO TO 5 4 CONTINUE IF (SJ(J).GT.0) THEN WRITE (ILP,1001) II,IH,AI,YI,SI,TI,YJ(J),SJ(J) WRITE (LU2) II,IH,AI,YJ(J),SJ(J),TI ELSE WRITE (ILP,1002) II,IH,AI,YI,SI,TI END IF GO TO 1 5 CONTINUE 1001 FORMAT ('0OLD ',I6,2X,3I3,2X,4F7.2,3F10.2/' NEW ',47X,2F10.2) 1002 FORMAT ('0MEASUREMENT DELETED IN ''NEW.BLP'' FILE:'/' OLD ',I6, &2X,3I3,2X,4F7.2) WRITE (LU2) II,IH,AI,YI,SI,TI GO TO 1 8 IF (N.EQ.M) GO TO 11 9 CLOSE (UNIT=LU3,STATUS='DELETE') RETURN END C----------------------------------------------------------------------- C C THE FOLLOWING SUBROUTINES - SOME SLIGHTLY MODIFIED FOR USE IN C PROGRAM VIEW - ARE FROM PROGRAM BGLP. C C----------------------------------------------------------------------- SUBROUTINE BGLP1 C C LINEAR BACKGROUND C PARAMETER (NMAX=100) DIMENSION X(NMAX),Y(NMAX) COMMON /BLOCK7/ NSTEP,XX(NMAX),YY(NMAX),VX(NMAX),VY(NMAX) COMMON /BLOCK8/ TH0,X0,SIGX0,XE,X1,X2,W1,W2,I1,I2,L1,L2,R,SIGR,EPS COMMON /BG/ A0,A1,E0,E1,B0,B1,J1,J2 E0=0 E1=0 B0=0 B1=0 DATA RAD /0.017453293/ C C FIT A STRAIGHT LINE TO THE PROFILE IN THE TWO BACKGROUND REGIONS, C REQUIRING THAT, FOR THE PURPOSE OF FITTING THE LINE, EACH C BACKGROUND SAMPLE INCLUDES AT LEAST TWO POINTS. C J1=MAX(L1-1,I1+1) J2=MIN(L2+1,I2-1) N=0 DO 1 I=I1,J1 N=N+1 X(N)=XX(I) Y(N)=YY(I) 1 CONTINUE DO 2 I=J2,I2 N=N+1 X(N)=XX(I) Y(N)=YY(I) 2 CONTINUE CALL FITLINE (N,X,Y,A0,A1,V0,V1,COV) C C SUBTRACT LINEAR BACKGROUND POINT-BY-POINT, DIVIDE OUT LORENTZ AND C POLARIZATION FACTORS POINT-BY-POINT, AND INTEGRATE THE PEAK C PROFILE BETWEEN PEAK LIMITS. C R=0 VR=0 DO 9 I=MAX(I1+1,L1),MIN(L2,I2-1) C C BACKGROUND SUBTRACTION C XI=XX(I) YI=YY(I)-A0-A1*XI VI=VY(I)+V0+V1*XI**2+2*XI*COV C C LORENTZ AND POLARIZATION CORRECTIONS C TH=(XX(I)+TH0)*RAD VT=VX(I)*RAD**2 CALL LZ (TH,VT,FL,VL) CALL PZ (TH,VT,FP,VP) CV=SQRT(VL*VP) YJ=YI IF (YI.EQ.0) YI=0.5*SQRT(VI) VI=(VI/YI**2+VL/FL**2+VP/FP**2+CV/(FL*FP))*(YI/(FL*FP))**2 YI=YJ/(FL*FP) C C INTEGRATED NET INTENSITY BY TRAPEZOIDAL INTEGRATION C DX=0.50*(XX(I+1)-XX(I-1)) VD=0.25*(VX(I+1)+VX(I-1)) R=R+YI*DX VR=VR+YI**2*VD+DX**2*VI 9 CONTINUE SIGR=SQRT(VR) RETURN END C----------------------------------------------------------------------- SUBROUTINE BGLP2 C C STRUCTURED BACKGROUND UNDER LOW ANGLE REFLECTIONS DUE TO THE C ABSORPTION EDGE OF THE BETA FILTER AND THE PEAK IN THE WHITE C RADIATION BACKGROUND C C R. J. NELMES (1975). ACTA CRYST. A31, 273-279. C PARAMETER (NMAX=100) DIMENSION X(NMAX),Y(NMAX) COMMON /BLOCK7/ NSTEP,XX(NMAX),YY(NMAX),VX(NMAX),VY(NMAX) COMMON /BLOCK8/ TH0,X0,SIGX0,XE,X1,X2,W1,W2,I1,I2,L1,L2,R,SIGR,EPS COMMON /BG/ A0,A1,E0,E1,B0,B1,J1,J2 DATA RAD /0.017453293/ C C FIT STRAIGHT LINE SEGMENTS TO THE STRUCTURED BACKGROUND. C A0=0 A1=0 VA0=0 VA1=0 COVA=0 E0=0 E1=0 VE0=0 VE1=0 COVE=0 B0=0 B1=0 VB0=0 VB1=0 COVB=0 J1=NINDEX(XE-EPS) J2=NINDEX(XE+EPS) IF (J2.LE.I1) GO TO 99 C C FIT ONE STRAIGHT LINE TO THE LOW ANGLE BACKGROUND BELOW THE TAIL C OF THE BETA FILTER ABSORPTION EDGE. C J1=MIN(J1,L1-1) J1=MAX(J1,I1) N=0 DO 1 I=I1,J1 N=N+1 X(N)=XX(I) Y(N)=YY(I) 1 CONTINUE IF (N.GE.3) THEN CALL FITLINE (N,X,Y,A0,A1,VA0,VA1,COVA) ELSE IF (N.EQ.2) THEN A0=0.5*(Y(1)+Y(2)) VA0=0.5*A0 ELSE A0=Y(1) VA0=A0 END IF C C FIT A SECOND STRAIGHT LINE TO THE ABSORPTION EDGE. C J2=MIN(J2,L1-1) J2=MAX(J2,J1) N=0 DO 2 I=J1,J2 N=N+1 X(N)=XX(I) Y(N)=YY(I) 2 CONTINUE IF (N.GE.3) THEN CALL FITLINE (N,X,Y,E0,E1,VE0,VE1,COVE) ELSE IF (N.EQ.2) THEN E0=0.5*(Y(1)+Y(2)) VE0=0.5*E0 ELSE E0=Y(1) VE0=E0 END IF C C THE ABSORPTION EDGE SHOULD HAVE POSITIVE INTERCEPT AND SLOPE. C IF (E0.LT.0.OR.E1.LT.0) GO TO 99 C C FIT A THIRD STRAIGHT LINE TO THE HIGH ANGLE, WHITE RADIATION C BACKGROUND ABOVE THE TAIL OF THE BRAGG PEAK. C J2=MAX(J2,L2+1) J2=MIN(J2,I2) N=0 DO 3 I=J2,I2 N=N+1 X(N)=XX(I) Y(N)=YY(I) 3 CONTINUE IF (N.GE.3) THEN CALL FITLINE (N,X,Y,B0,B1,VB0,VB1,COVB) ELSE IF (N.EQ.2) THEN B0=0.5*(Y(1)+Y(2)) VB0=0.5*B0 ELSE B0=Y(1) VB0=B0 END IF C C THE HIGH ANGLE BACKGROUND SHOULD HAVE POSITIVE INTERCEPT AND C NEGATIVE SLOPE. C IF (B0.LT.0.OR.B1.GT.0) GO TO 99 C C EVALUATE THE LOW ANGLE BACKGROUND AT XE - EPSILON AND THE HIGH C ANGLE BACKGROUND AT XE + EPSILON. C J1=NINDEX(XE-EPS) J2=NINDEX(XE+EPS) J1=MAX(J1,I1) J2=MIN(J2,I2) X1=XX(J1) X2=XX(J2) Y1=A0+A1*X1 Y2=B0+B1*X2 V1=VA0+VA1*X1**2+X1*COVA V2=VB0+VB1*X2**2+X2*COVB C C THE BACKGROUND AT XE - EPSILON SHOULD BE LOWER THAN THE C BACKGROUND AT XE + EPSILON. C IF (Y1.GE.Y2) GO TO 99 C C CALCULATE THE COEFFICIENTS OF THE LINE JOINING THE LOW ANGLE C BACKGROUND AT X1 = XE - EPSILON TO THE HIGH ANGLE BACKGROUND AT C X2 = XE + EPSILON. C IF (J1.LT.L1.AND.L1.LT.J2) THEN C C THE LOW ANGLE PEAK LIMIT FALLS ON THE ABSORPTION EDGE. C X1=XX(L1-1) Y1=E0+E1*X1 V1=VE0+VE1*X1**2+X1*COVE END IF E1=(Y2-Y1)/(X2-X1) E0=Y1-E1*X1 VE1=(V1+V2)/(X2-X1)**2 VE0=(V1*X2**2+V2*X1**2)/(X2-X1)**2 COVE=-(V1*X2+V2*X1)/(X2-X1)**2 C C SUBTRACT BACKGROUND POINT-BY-POINT, DIVIDE OUT LORENTZ AND C POLARIZATION FACTORS POINT-BY-POINT, AND INTEGRATE THE PEAK C PROFILE BETWEEN PEAK LIMITS. C R=0 VR=0 DO 9 I=MAX(I1+1,L1),MIN(L2,I2-1) C C BACKGROUND SUBTRACTION C XI=XX(I) IF (XI.LT.X1) THEN C0=A0 C1=A1 V0=VA0 V1=VA1 COV=COVA ELSE IF (X1.LE.XI.AND.XI.LE.X2) THEN C0=E0 C1=E1 V0=VE0 V1=VE1 COV=COVE ELSE C0=B0 C1=B1 V0=VB0 V1=VB1 COV=COVB END IF YI=YY(I)-C0-C1*XI VI=VY(I)+V0+V1*XI**2+2*XI*COV C C PREVENT TRUNCATION OF THE LOW ANGLE TAIL OF THE PEAK BY THE C ESTIMATED ABSORPTION EDGE. C IF (XI.LT.0.AND.YI.LT.0) YI=0 C C LORENTZ AND POLARIZATION CORRECTIONS C TH=(XX(I)+TH0)*RAD VT=VX(I)*RAD**2 CALL LZ (TH,VT,FL,VL) CALL PZ (TH,VT,FP,VP) CV=SQRT(VL*VP) YJ=YI IF (YI.EQ.0) YI=0.5*SQRT(VI) VI=(VI/YI**2+VL/FL**2+VP/FP**2+CV/(FL*FP))*(YI/(FL*FP))**2 YI=YJ/(FL*FP) C C INTEGRATED NET INTENSITY BY TRAPEZOIDAL INTEGRATION C DX=0.50*(XX(I+1)-XX(I-1)) VD=0.25*(VX(I+1)+VX(I-1)) R=R+YI*DX VR=VR+YI**2*VD+DX**2*VI 9 CONTINUE SIGR=SQRT(VR) IF (R.LT.CUTOFF*SIGR) GO TO 99 RETURN 99 CONTINUE C C IF STRUCTURED BACKGROUND MODEL IS NOT APPROPRIATE, RESORT TO C LINEAR BACKGROUND. C J2=MIN(NINDEX(XE+EPS),L1-1) I1=MAX(I1,J2) CALL BGLP1 RETURN END C----------------------------------------------------------------------- SUBROUTINE CENTROID (WINDOW,CUTOFF,XX,YY,VX,VY,I1,I2,X0,SIGX0) C C INTENSITY-WEIGHTED CENTROID OF PEAK ABOVE LINEAR BACKGROUND C DIMENSION XX(1),YY(1),VX(1),VY(1) PARAMETER (NMAX=100) DIMENSION X(NMAX),Y(NMAX),V(NMAX) EQUIVALENCE (V(1),X(1)) C C FIND THE SMALLEST X AND LARGEST X FOR WHICH THE LOCAL AVERAGE Y(X) C INCREASES SIGNIFICANTLY ABOVE THE LOCAL AVERAGE BACKGROUND. C N=0 DO 11 I=I1,I2,+1 N=N+1 Y(N)=YY(I) V(N)=VY(I) 11 CONTINUE CALL FINDL (N,Y,V,WINDOW,CUTOFF,L) IF (L.EQ.0) GO TO 9 L1=I1+L N=0 DO 12 I=I2,I1,-1 N=N+1 Y(N)=YY(I) V(N)=VY(I) 12 CONTINUE CALL FINDL (N,Y,V,WINDOW,CUTOFF,L) IF (L.EQ.0) GO TO 9 L2=I2-L IF (L1.GE.L2) GO TO 9 C C FIT A STRAIGHT LINE TO THE BACKGROUND NEAR THE PEAK LIMITS. C L=NINT(WINDOW*N) N=0 DO 1 I=L1-L,L2+L IF (I.LT.L1.OR.I.GT.L2) THEN N=N+1 X(N)=XX(I) Y(N)=YY(I) END IF 1 CONTINUE CALL FITLINE (N,X,Y,C0,C1,V0,V1,COV) C C FIND THE INTENSITY-WEIGHTED CENTROID BETWEEN THE PEAK LIMITS. C SUMY=0 SUMYX=0 SUMX2V=0 SUMXV=0 SUMV=0 DO 2 I=L1,L2 XI=XX(I) YI=YY(I)-C0-C1*XI VI=VY(I)+V0+V1*XI**2+2*COV*XI DX=0.50*(XX(I+1)-XX(I-1)) VD=0.25*(VX(I+1)+VX(I-1)) VI=YI**2*VD+DX**2*VI YI=YI*DX SUMY=SUMY+YI SUMYX=SUMYX+YI*XI SUMX2V=SUMX2V+XI**2*VI SUMXV=SUMXV+XI*VI SUMV=SUMV+VI 2 CONTINUE C C TEST FOR A SIGNIFICANTLY POSITIVE PEAK ABOVE BACKGROUND BETWEEN C THE TRUNCATED PEAK LIMITS. C IF (SUMY.LT.CUTOFF*SQRT(SUMV)) GO TO 9 X0=SUMYX/SUMY SIGX0=SQRT(ABS(SUMX2V-2*X0*SUMXV+X0**2*SUMV)/SUMY**2) RETURN 9 CONTINUE C C RETURN X0 = 90 DEGREES TO FLAG AN ERROR CONDITION. C X0=90 RETURN END C----------------------------------------------------------------------- SUBROUTINE FINDL (N,Y,V,WINDOW,CUTOFF,L) C C FIND THE X AT WHICH THE LOCAL AVERAGE Y(X) FIRST INCREASES C SIGNIFICANTLY ABOVE THE LOCAL AVERAGE BACKGROUND. C DIMENSION Y(1),V(1) L=NINT(WINDOW*N) DO 1 I=L+1,N-L YB=0 YI=0 VB=0 VI=0 DO 2 J=0,L YB=YB+Y(I-J) VB=VB+V(I-J) YI=YI+Y(I+J) VI=VI+V(I+J) 2 CONTINUE IF (YI-YB.GT.CUTOFF*SQRT(VI+VB)) THEN L=I RETURN END IF 1 CONTINUE L=0 RETURN END C----------------------------------------------------------------------- SUBROUTINE FITLINE (N,X,Y,C0,C1,V0,V1,COV) C C LEAST-SQUARES STRAIGHT LINE BACKGROUND C C GIVE EQUAL WEIGHT TO ALL BACKGROUND POINTS BECAUSE: (1) IT IS C EXPECTED THAT THE BACKGROUND SHOULD BE APPROXIMATELY CONSTANT, AND C THUS CONSTANT WEIGHTS SHOULD BE APPROPRIATE; AND (2) WEIGHTS C PROPORTIONAL TO THE RECIPROCALS OF THE VARIANCES OF THE C EXPERIMENTAL MEASUREMENTS, ALTHOUGH APPROPRIATE FOR NORMALLY C DISTRIBUTED DATA, ARE NOT APPROPRIATE FOR POISSON DISTRIBUTED C DATA. WEIGHT = 1/VARIANCE WOULD BIAS THE FIT TOWARD THE C BACKGROUND POINTS WITH THE LOWEST COUNT RATES. (SEE P. F. PRICE C (1979). ACTA CRYST. A35, 57-60.) C DIMENSION X(1),Y(1) SUMX=0 SUMY=0 SUMX2=0 SUMXY=0 DO 1 I=1,N XI=X(I) YI=Y(I) SUMX=SUMX+XI SUMY=SUMY+YI SUMX2=SUMX2+XI*XI SUMXY=SUMXY+XI*YI 1 CONTINUE DET=N*SUMX2-SUMX*SUMX C0=(SUMX2*SUMY-SUMX*SUMXY)/DET C1=(-SUMX*SUMY+N*SUMXY)/DET V0=SUMX2/DET V1=N/DET COV=-SUMX/DET CHISQ=0 DO 2 I=1,N CHISQ=CHISQ+(Y(I)-C0-C1*X(I))**2 2 CONTINUE CHISQ=CHISQ/(N-2) V0=CHISQ*V0 V1=CHISQ*V1 COV=CHISQ*COV RETURN END C----------------------------------------------------------------------- SUBROUTINE GEOM (IDIFF,ANGLES,TWOTH,OMEGA,CHI,PHI) C C DIFFRACTOMETER GEOMETRIES C DIMENSION ANGLES(4) REAL KAPPA C C KAPPA AXIS INCLINATION ANGLE C C ALPHA = PI/3.6 = 50 DEGREES C DATA SINA, COSA /0.766044, 0.642788/ DATA DEG, RAD /57.2957795, 0.017453293/ IF (IDIFF.EQ.1) THEN C C HAMILTON'S DIFFRACTOMETER AXES C C WALTER C. HAMILTON (1974). INTERNATIONAL TABLES FOR X-RAY C CRYSTALLOGRAPHY, VOL. IV, PP. 273-284. C TWOTH = ANGLES(1) OMEGA = ANGLES(2) CHI = ANGLES(3) PHI = ANGLES(4) ELSE IF (IDIFF.EQ.2) THEN C C BUSING'S AND LEVY'S DIFFRACTOMETER AXES C C W. R. BUSING AND H. A. LEVY (1967). ACTA CRYST. 22, 457-464. C TWOTH = ANGLES(1) OMEGA = ANGLES(2) CHI = ANGLES(3) PHI = ANGLES(4) C C TRANSFORM TO SETTING ANGLES AS DEFINED BY HAMILTON. C OMEGA = -OMEGA CHI = -CHI PHI = -PHI ELSE IF (IDIFF.EQ.3) THEN C C SIEMENS (NEE NICOLET, NEE SYNTEX) P3 DIFFRACTOMETER C TWOTH = ANGLES(1) OMEGA = ANGLES(2) PHI = ANGLES(3) CHI = ANGLES(4) C C TRANSFORM TO SETTING ANGLES AS DEFINED BY HAMILTON. C THETA = TWOTH/2 OMEGA = OMEGA - THETA OMEGA = -OMEGA ELSE IF (IDIFF.EQ.4) THEN C C ENRAF-NONIUS CAD4 DIFFRACTOMETER KAPPA ANGLES C THETA = ANGLES(1) PHI = ANGLES(2) OMEGA = ANGLES(3) KAPPA = ANGLES(4) C C TRANSFORM FROM KAPPA TO EULERIAN SETTING ANGLES. C KAPPA = KAPPA*RAD SINX = SINA*SIN(KAPPA/2) COSX = SQRT(COSA**2 + SINA**2*COS(KAPPA/2)**2) CHI = 2*ARCTAN(SINX,COSX) C C KAPPA IS MECHANICALLY LIMITED TO THE RANGE 0 TO 180 DEGREES. C WITH ALPHA = 50 DEGREES, THIS LIMITS CHI TO THE RANGE -100 TO C +100 DEGREES. C SINX = COSA*SIN(KAPPA/2)/COS(CHI/2) COSX = COS(KAPPA/2)/COS(CHI/2) DELTA = ARCTAN(SINX,COSX)*DEG OMEGA = OMEGA + DELTA PHI = PHI + DELTA CHI = CHI*DEG C C TRANSFORM TO SETTING ANGLES AS DEFINED BY HAMILTON. C TWOTH = 2*THETA OMEGA = OMEGA - THETA OMEGA = -OMEGA CHI = -CHI PHI = -PHI END IF C C RETURN ANGLES -180 .LE. A .LE. +180. C TWOTH = A180(TWOTH) OMEGA = A180(OMEGA) PHI = A180(PHI) CHI = A180(CHI) RETURN END C----------------------------------------------------------------------- FUNCTION ARCTAN(SINX,COSX) IF (SINX.LT.-1) SINX=-1 IF (SINX.GT.+1) SINX=+1 IF (COSX.LT.-1) COSX=-1 IF (COSX.GT.+1) COSX=+1 ARCTAN=ATAN2(SINX,COSX) RETURN END C----------------------------------------------------------------------- FUNCTION A180(A) A=MOD(A,360.0) IF (A.LT.-180) A=A+360 IF (A.GT.+180) A=A-360 A180=A RETURN END C----------------------------------------------------------------------- FUNCTION A360(A) A=MOD(A,360.0) IF (A.LT.0) A=A+360 A360=A RETURN END C----------------------------------------------------------------------- SUBROUTINE HKLCOND (H,I) C C CONDITIONS LIMITING POSSIBLE REFLECTIONS C CHARACTER*18 H,A DIMENSION A(38) DATA A / & 'HKL,H+K=2N ','HKL,K+L=2N ','HKL,L+H=2N ', & 'HKL,H+K,K+L,L+H=2N','HKL,H+K+L=2N ','HKL,-H+K+L=3N ', & 'HKL,H-K+L=3N ','HKL,H-K=3N ','HK0,H=2N ', & 'HK0,K=2N ','HK0,H+K=2N ','HK0,H+K=4N ', & '0KL,K=2N ','0KL,L=2N ','0KL,K+L=2N ', & '0KL,K+L=4N ','H0L,L=2N ','H0L,H=2N ', & 'H0L,L+H=2N ','H0L,L+H=4N ','HH(-2H)L,L=2N ', & 'H(-H)0L,L=2N ','HHL,L=2N(R-AXES) ','HHL,L=2N ', & 'HKH,K=2N ','HKK,H=2N ','HHL,2H+L=4N ', & 'HKH,2H+K=4N ','HKK,2K+H=4N ','H00,H=2N ', & 'H00,H=4N ','0K0,K=2N ','0K0,K=4N ', & '00L,L=2N ','00L,L=4N ','000L,L=2N ', & '000L,L=3N ','000L,L=6N '/ IF (I.EQ.0) THEN DO 1 I=1,38 IF (H.EQ.A(I)) RETURN 1 CONTINUE I=0 ELSE H=A(I) END IF RETURN END C----------------------------------------------------------------------- SUBROUTINE HKLTEST (H,K,L,NCOND,ICOND,IABSENT) C C CONDITIONS LIMITING POSSIBLE REFLECTIONS C C ( 1) HKL,H+K=2N C ( 2) HKL,K+L=2N C ( 3) HKL,L+H=2N C ( 4) HKL,H+K,K+L,L+H=2N C ( 5) HKL,H+K+L=2N C ( 6) HKL,-H+K+L=3N C ( 7) HKL,H-K+L=3N C ( 8) HKL,H-K=3N C ( 9) HK0,H=2N C (10) HK0,K=2N C (11) HK0,H+K=2N C (12) HK0,H+K=4N C (13) 0KL,K=2N C (14) 0KL,L=2N C (15) 0KL,K+L=2N C (16) 0KL,K+L=4N C (17) H0L,L=2N C (18) H0L,H=2N C (19) H0L,L+H=2N C (20) H0L,L+H=4N C (21) HH(-2H)L,L=2N C (22) H(-H)0L,L=2N C (23) HHL,L=2N (RHOMBOHEDRAL AXES) C (24) HHL,L=2N C (25) HKH,K=2N C (26) HKK,H=2N C (27) HHL,2H+L=4N C (28) HKH,2H+K=4N C (29) HKK,2K+H=4N C (30) H00,H=2N C (31) H00,H=4N C (32) 0K0,K=2N C (33) 0K0,K=4N C (34) 00L,L=2N C (35) 00L,L=4N C (36) 000L,L=2N C (37) 000L,L=3N C (38) 000L,L=6N C C SET IABSENT = 0 FOR AN ALLOWED REFLECTION, C IABSENT = 1 FOR A FORBIDDEN REFLECTION. C INTEGER H DIMENSION ICOND(38) IABSENT=0 IF (NCOND.EQ.0) RETURN DO 90 I=1,NCOND GO TO (1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22, & 23,24,25,26,27,28,29,30,31,32,33,34,35,36,37,38) ICOND(I) 1 IF (MOD(H+K,2).NE.0) GO TO 91 GO TO 90 2 IF (MOD(K+L,2).NE.0) GO TO 91 GO TO 90 3 IF (MOD(L+H,2).NE.0) GO TO 91 GO TO 90 4 IF (MOD(H+K,2).NE.0.OR.MOD(K+L,2).NE.0.OR.MOD(L+H,2).NE.0) &GO TO 91 GO TO 90 5 IF (MOD(H+K+L,2).NE.0) GO TO 91 GO TO 90 6 IF (MOD(-H+K+L,3).NE.0) GO TO 91 GO TO 90 7 IF (MOD(H-K+L,3).NE.0) GO TO 91 GO TO 90 8 IF (MOD(H-K,3).NE.0) GO TO 91 GO TO 90 9 IF (L.EQ.0.AND.MOD(H,2).NE.0) GO TO 91 GO TO 90 10 IF (L.EQ.0.AND.MOD(K,2).NE.0) GO TO 91 GO TO 90 11 IF (L.EQ.0.AND.MOD(H+K,2).NE.0) GO TO 91 GO TO 90 12 IF (L.EQ.0.AND.MOD(H+K,4).NE.0) GO TO 91 GO TO 90 13 IF (H.EQ.0.AND.MOD(K,2).NE.0) GO TO 91 GO TO 90 14 IF (H.EQ.0.AND.MOD(L,2).NE.0) GO TO 91 GO TO 90 15 IF (H.EQ.0.AND.MOD(K+L,2).NE.0) GO TO 91 GO TO 90 16 IF (H.EQ.0.AND.MOD(K+L,4).NE.0) GO TO 91 GO TO 90 17 IF (K.EQ.0.AND.MOD(L,2).NE.0) GO TO 91 GO TO 90 18 IF (K.EQ.0.AND.MOD(H,2).NE.0) GO TO 91 GO TO 90 19 IF (K.EQ.0.AND.MOD(L+H,2).NE.0) GO TO 91 GO TO 90 20 IF (K.EQ.0.AND.MOD(L+H,4).NE.0) GO TO 91 GO TO 90 21 IF (K.EQ.H.AND.MOD(L,2).NE.0) GO TO 91 GO TO 90 22 IF (K.EQ.-H.AND.MOD(L,2).NE.0) GO TO 91 GO TO 90 23 IF (K.EQ.H.AND.MOD(L,2).NE.0) GO TO 91 GO TO 90 24 IF (ABS(K).EQ.ABS(H).AND.MOD(L,2).NE.0) GO TO 91 GO TO 90 25 IF (ABS(L).EQ.ABS(H).AND.MOD(K,2).NE.0) GO TO 91 GO TO 90 26 IF (ABS(L).EQ.ABS(K).AND.MOD(H,2).NE.0) GO TO 91 GO TO 90 27 IF (ABS(K).EQ.ABS(H).AND.MOD(2*H+L,4).NE.0) GO TO 91 GO TO 90 28 IF (ABS(L).EQ.ABS(H).AND.MOD(2*H+K,4).NE.0) GO TO 91 GO TO 90 29 IF (ABS(L).EQ.ABS(K).AND.MOD(2*K+H,4).NE.0) GO TO 91 GO TO 90 30 IF (K.EQ.0.AND.L.EQ.0.AND.MOD(H,2).NE.0) GO TO 91 GO TO 90 31 IF (K.EQ.0.AND.L.EQ.0.AND.MOD(H,4).NE.0) GO TO 91 GO TO 90 32 IF (H.EQ.0.AND.L.EQ.0.AND.MOD(K,2).NE.0) GO TO 91 GO TO 90 33 IF (H.EQ.0.AND.L.EQ.0.AND.MOD(K,4).NE.0) GO TO 91 GO TO 90 34 IF (H.EQ.0.AND.K.EQ.0.AND.MOD(L,2).NE.0) GO TO 91 GO TO 90 35 IF (H.EQ.0.AND.K.EQ.0.AND.MOD(L,4).NE.0) GO TO 91 GO TO 90 36 IF (H.EQ.0.AND.K.EQ.0.AND.MOD(L,2).NE.0) GO TO 91 GO TO 90 37 IF (H.EQ.0.AND.K.EQ.0.AND.MOD(L,3).NE.0) GO TO 91 GO TO 90 38 IF (H.EQ.0.AND.K.EQ.0.AND.MOD(L,6).NE.0) GO TO 91 90 CONTINUE RETURN 91 IABSENT=1 RETURN END C----------------------------------------------------------------------- SUBROUTINE INPUT (BTIME,BDATE) C C READ CONTROL DATA AND EVALUATE VARIABLES. C CHARACTER ATIME*8,ADATE*9,TITLE*80,FILE1*80,FILE2*80,HCOND*18, & BTIME*8,BDATE*9 COMMON /VTKB/ IVT,IKB,IC COMMON /BLOCK0/ IO1,IO2,IO3,ILP COMMON /BLOCKC/ NCOND,ICOND(38) COMMON /BLOCK1/ ATIME,ADATE,TITLE,FILE1,FILE2 COMMON /BLOCK2/ GINV(3,3),WLE,WL1,WL0,WL2 COMMON /BLOCK3/ WINDOW,CUTOFF,NSMOOTH,II1,II2,XT1,XT2,SC1,SC2, & XX0,XW1,XW2 COMMON /BLOCK4/ C0KIN,C1KIN,C2KIN,C0DYN,C1DYN,C2DYN,FRACTD,VARFD, & INEUTRON COMMON /BLOCK5/ IDIFF,U(3,3),MODEL1,MODEL2,Q1(3,3),Q2(3,3),T1,T2, & F,SIGU(3,3),SIGQ1(3,3),SIGQ2(3,3),SIGT1,SIGT2,TANTHM COMMON /BLOCK6/ VARTHE,VARTIM,TAU,VARTAU,ATT,VARATT DIMENSION CELL(6) CHARACTER TYPE1*4,TYPE2*4 PI=ACOS(-1.0) C C READ CONTROL DATA. C OPEN (UNIT=IO1,FILE='bglp.dat',STATUS='OLD') READ (IO1,500) BTIME READ (IO1,500) BDATE READ (IO1,500) TITLE READ (IO1,500) FILE1 READ (IO1,500) FILE2 READ (IO1,501) CELL READ (IO1,502) IDIFF,INEUTRON READ (IO1,501) THM,RHOM,FRACTD,SIGFD READ (IO1,501) WLE,WL1,WL0,WL2 READ (IO1,501) SIGTHE,SIGTIM READ (IO1,501) TAU,SIGTAU,ATT,SIGATT READ (IO1,501,END=8) WINDOW READ (IO1,501,END=8) CUTOFF READ (IO1,502,END=8) NSMOOTH READ (IO1,505,END=8) SC1,SC2 READ (IO1,503,END=8) II1,II2,XT1,XT2 READ (IO1,504,END=8) ((U(I,J),J=1,3),I=1,3) READ (IO1,504) ((SIGU(I,J),J=1,3),I=1,3) READ (IO1,500) TYPE1 READ (IO1,504) ((Q1(I,J),J=1,3),I=1,3) READ (IO1,504) ((SIGQ1(I,J),J=1,3),I=1,3) READ (IO1,504) T1 READ (IO1,504) SIGT1 READ (IO1,500) TYPE2 READ (IO1,504) ((Q2(I,J),J=1,3),I=1,3) READ (IO1,504) ((SIGQ2(I,J),J=1,3),I=1,3) READ (IO1,504) T2 READ (IO1,504) SIGT2 READ (IO1,504,END=9) F READ (IO1,505,END=9) XX0,XW1,XW2 500 FORMAT (1X,A) 501 FORMAT (1X,6F10.5) 502 FORMAT (1X,2I10) 503 FORMAT (1X,2I10,2F10.3) 504 FORMAT (3(1X,3E15.8/)) 505 FORMAT (1X,4F10.2) GO TO 9 8 XX0=0.5 XW1=0.25 XW2=0.75 9 CLOSE (UNIT=IO1,STATUS='KEEP') IF (TYPE1.EQ.'LRNZ') MODEL1=1 IF (TYPE1.EQ.'GAUS') MODEL1=2 IF (TYPE1.EQ.'LRNZ') MODEL1=1 IF (TYPE1.EQ.'PLMC') MODEL1=3 IF (TYPE2.EQ.'LRNZ') MODEL2=1 IF (TYPE2.EQ.'GAUS') MODEL2=2 IF (TYPE2.EQ.'PLMC') MODEL2=3 C C CALCULATE RECIPROCAL SPACE METRIC TENSOR. C CALL METRIC (CELL,GINV) C C EVALUATE MONOCHROMATOR VARIABLES. C THM=THM*PI/180 RHOM=RHOM*PI/180 COSSQ=COS(RHOM)**2 SINSQ=SIN(RHOM)**2 C0KIN=COS(2*THM)**2 C1KIN=COSSQ+SINSQ*C0KIN C2KIN=SINSQ+COSSQ*C0KIN C0DYN=ABS(COS(2*THM)) C1DYN=COSSQ+SINSQ*C0DYN C2DYN=SINSQ+COSSQ*C0DYN IF (FRACTD.NE.0.AND.SIGFD.EQ.0) SIGFD=0.05*FRACTD TANTHM=TAN(THM) C C CONVERT FROM MICROSECONDS TO SECONDS. C TAU=TAU*1E-6 SIGTAU=SIGTAU*1E-6 SIGTIM=SIGTIM*1E-6 C C CALCULATE VARIANCE VALUES. C VARFD=SIGFD**2 VARTAU=SIGTAU**2 VARATT=SIGATT**2 VARTHE=SIGTHE**2 VARTIM=SIGTIM**2 C C READ CONDITIONS LIMITING POSSIBLE REFLECTIONS. C OPEN (UNIT=IO1,FILE='hklcond.dat',STATUS='UNKNOWN') NCOND=0 5 READ (IO1,506,END=6) HCOND 506 FORMAT (A) I=0 CALL HKLCOND (HCOND,I) NCOND=NCOND+1 ICOND(NCOND)=I GO TO 5 6 IF (NCOND.EQ.0) THEN CLOSE (UNIT=IO1,STATUS='DELETE') ELSE CLOSE (UNIT=IO1,STATUS='KEEP') END IF C C OPEN INPUT AND OUTPUT REFLECTION DATA FILES. C OPEN (UNIT=IO1,FILE=FILE1,STATUS='OLD',FORM='UNFORMATTED') OPEN (UNIT=IO3,STATUS='SCRATCH',FORM='UNFORMATTED') OPEN (UNIT=ILP,FILE='view.lp',STATUS='NEW') RETURN END C----------------------------------------------------------------------- SUBROUTINE LZ (TH,VT,FL,VL) C C LORENTZ FACTOR. EQUATORIAL DIFFRACTION GEOMETRY C U=COS(2*TH) V=SIN(2*TH) FL=1/V VL=VT*(-2*U/V**2)**2 RETURN END C----------------------------------------------------------------------- SUBROUTINE METRIC (A,GB) C C CALCULATE RECIPROCAL LATTICE PARAMETERS AND DIRECT AND C RECIPROCAL LATTICE METRIC TENSORS FROM DIRECT LATTICE C PARAMETERS. C DIMENSION A(6),B(6),GA(3,3),GB(3,3) PI=ACOS(-1.0) DEG=180/PI RAD=PI/180 CA=COS(A(4)*RAD) CB=COS(A(5)*RAD) CG=COS(A(6)*RAD) SA=SIN(A(4)*RAD) SB=SIN(A(5)*RAD) SG=SIN(A(6)*RAD) V=A(1)*A(2)*A(3)*SQRT(1-CA**2-CB**2-CG**2+2*CA*CB*CG) B(1)=A(2)*A(3)*SA/V B(2)=A(1)*A(3)*SB/V B(3)=A(1)*A(2)*SG/V B(4)=(CB*CG-CA)/(SB*SG) B(5)=(CA*CG-CB)/(SA*SG) B(6)=(CA*CB-CG)/(SA*SB) GA(1,1)=A(1)*A(1) GA(1,2)=A(1)*A(2)*CG GA(1,3)=A(1)*A(3)*CB GA(2,1)=GA(1,2) GA(2,2)=A(2)*A(2) GA(2,3)=A(2)*A(3)*CA GA(3,1)=GA(1,3) GA(3,2)=GA(2,3) GA(3,3)=A(3)*A(3) GB(1,1)=B(1)*B(1) GB(1,2)=B(1)*B(2)*B(6) GB(1,3)=B(1)*B(3)*B(5) GB(2,1)=GB(1,2) GB(2,2)=B(2)*B(2) GB(2,3)=B(2)*B(3)*B(4) GB(3,1)=GB(1,3) GB(3,2)=GB(2,3) GB(3,3)=B(3)*B(3) B(4)=ACOS(B(4))*DEG B(5)=ACOS(B(5))*DEG B(6)=ACOS(B(6))*DEG RETURN END C----------------------------------------------------------------------- FUNCTION NINDEX(X) C C FIND THE INDEX I OF THE ABSCISSA XX(I) NEAREST THE VALUE X. C C THE ABSCISSAE IN THE ARRAY XX(NSTEP) MUST BE ORDERED SUCH THAT XX(I) C INCREASES WITH I. C PARAMETER (NMAX=100) COMMON /BLOCK7/ NSTEP,XX(NMAX),YY(NMAX),VX(NMAX),VY(NMAX) IF (X.LT.XX(1)) THEN NINDEX=0 RETURN END IF IF (X.GT.XX(NSTEP)) THEN NINDEX=NSTEP+1 RETURN END IF DO 1 I=2,NSTEP IF (XX(I).GE.X) THEN IF (X-XX(I-1).LT.XX(I)-X) THEN NINDEX=I-1 ELSE NINDEX=I END IF RETURN END IF 1 CONTINUE END C----------------------------------------------------------------------- SUBROUTINE PROCESS PARAMETER (NMAX=100) CHARACTER ATIME*8,ADATE*9,TITLE*80,FILE1*80,FILE2*80 COMMON /VTKB/ IVT,IKB,IC COMMON /RFLN/ IEND,II,IH,IK,IL,ANGLES(4),WIDTH,SPEED,FBG, & BG1,SCAN,BG2,YNET,SIGY,XTIME COMMON /BG/ A0,A1,E0,E1,B0,B1,J1,J2 COMMON /BLOCK0/ IO1,IO2,IO3,ILP COMMON /BLOCKC/ NCOND,ICOND(38) COMMON /BLOCK1/ ATIME,ADATE,TITLE,FILE1,FILE2 COMMON /BLOCK2/ GINV(3,3),WLE,WL1,WL0,WL2 COMMON /BLOCK3/ WINDOW,CUTOFF,NSMOOTH,II1,II2,XT1,XT2,SC1,SC2, & XX0,XW1,XW2 COMMON /BLOCK5/ IDIFF,U(3,3),MODEL1,MODEL2,Q1(3,3),Q2(3,3),T1,T2, & F,SIGU(3,3),SIGQ1(3,3),SIGQ2(3,3),SIGT1,SIGT2,TANTHM COMMON /BLOCK6/ VARTHE,VARTIM,TAU,VARTAU,ATT,VARATT COMMON /BLOCK7/ NSTEP,XX(NMAX),YY(NMAX),VX(NMAX),VY(NMAX) COMMON /BLOCK8/ TH0,X0,SIGX0,XE,X1,X2,W1,W2,I1,I2,L1,L2,R,SIGR,EPS DIMENSION H(3),YYTEMP(NMAX),VYTEMP(NMAX) PI=ACOS(-1.0) DEG=180/PI RAD=PI/180 C C SKIP REFLECTIONS THAT ARE NOT WITHIN THE SERIAL NUMBER OR EXPOSURE C TIME LIMITS. C IF (ABS(II).LT.II1.OR.ABS(II).GT.II2) GO TO 9 IF (XTIME.LT.XT1.OR.XTIME.GT.XT2) GO TO 9 C C SKIP SYMMETRY FORBIDDEN REFLECTIONS. C CALL HKLTEST (IH,IK,IL,NCOND,ICOND,IABSENT) IF (IABSENT.NE.0) GO TO 9 C C STORE INDICES IN REAL ARRAY. C H(1)=IH H(2)=IK H(3)=IL C C CALCULATE SIN(THETA)/LAMBDA. C S=SINTHL(H,GINV) C C CALCULATE SPECTRAL THETA VALUES. C THE=ASIN(S*WLE)*DEG TH1=ASIN(S*WL1)*DEG TH0=ASIN(S*WL0)*DEG TH2=ASIN(S*WL2)*DEG C C TRANSFORM TO SETTING ANGLES AS DEFINED BY WALTER HAMILTON, C INTERNATIONAL TABLES FOR X-RAY CRYSTALLOGRAPHY, VOL. IV, 1974, PP. C 273-284. C CALL GEOM (IDIFF,ANGLES,TWOTH,OMEGA,CHI,PHI) C C SET SCAN LIMITS (DEFAULT VALUES ARE SC1 = 0, SC2 = 1). C DX=XX(NSTEP)-XX(1) I1=NINDEX(XX(1)+SC1*DX) I2=NINDEX(XX(1)+SC2*DX) C C SMOOTH PROFILE FOR PEAK CENTROID LOCATION. C DO 11 I=1,NSTEP YYTEMP(I)=YY(I) VYTEMP(I)=VY(I) 11 CONTINUE IF (NSMOOTH.NE.0) THEN CALL SMOOTH (NSMOOTH,NSTEP,YYTEMP) CALL SMOOTH (NSMOOTH,NSTEP,VYTEMP) END IF C C FIND INTENSITY-WEIGHTED CENTROID OF PEAK ABOVE BACKGROUND. C CALL CENTROID (WINDOW,CUTOFF,XX,YYTEMP,VX,VYTEMP,I1,I2,X0,SIGX0) C C TREAT SCANS THAT HAVE NO SIGNIFICANT PEAK ABOVE BACKGROUND. C IF (X0.EQ.90) THEN IF (XX0.NE.0) THEN X0=XX(I1)+(XX(I2)-XX(I1))*XX0 SIGX0=0 ELSE CALL XCALC (XX,I1,I2,H,U,SIGU,GINV,TWOTH,OMEGA,CHI,PHI, & X0,SIGX0) END IF WEAK=1 ELSE WEAK=0 END IF C C CALCULATE POSITIONS OF SPECTRAL FEATURES. C XE=X0+(THE-TH0) X1=X0+(TH1-TH0) X2=X0+(TH2-TH0) IF (WEAK.EQ.0) CALL XCHECK (XE,X1,X0,X2,XX,YYTEMP,I1,I2) C C CALCULATE PEAK WIDTHS. C IF (XW1.NE.0.OR.XW2.NE.0) THEN W1=X1-(XX(I1)+(XX(I2)-XX(I1))*XW1) W2=(XX(I1)+(XX(I2)-XX(I1))*XW2)-X2 SIGW1=0 SIGW2=0 ELSE TANTH=TAN(TH0*RAD) CALL WCALC (XX,I1,I2,CHI,PHI,MODEL1,Q1,T1,F,TANTH,TANTHM, & SIGQ1,SIGT1,W1,SIGW1) CALL WCALC (XX,I1,I2,CHI,PHI,MODEL2,Q2,T2,F,TANTH,TANTHM, & SIGQ2,SIGT2,W2,SIGW2) END IF C C LIMIT THE ESTIMATED UNCERTAINTIES IN PEAK POSITION AND PEAK WIDTHS TO C ONE STEP. C DX=(XX(I2)-XX(I1))/(I2-I1+1) SIGX0=MIN(SIGX0,DX) SIGW1=MIN(SIGW1,DX) SIGW2=MIN(SIGW2,DX) C C CALCULATE PEAK LIMITS. C L1=NINDEX(X1-SIGX0-W1-SIGW1) L2=NINDEX(X2+SIGX0+W2+SIGW2) C C PERFORM BACKGROUND SUBTRACTION AND LORENTZ AND POLARIZATION C CORRECTIONS. C IF (WLE.EQ.0.OR.WEAK.EQ.1) THEN CALL BGLP1 ELSE C C EPSILON IS AN ESTIMATE OF THE HALF-WIDTH OF THE BETA-FILTER ABSORPTION C EDGE BROADENED DUE TO FINITE INSTRUMENTAL RESOLUTION. C EPS=0.5*(W1+W2) J1=NINDEX(XE-EPS) J2=NINDEX(XE+EPS) IF (J2.LT.L1-1) THEN I1=MAX(I1,J2) CALL BGLP1 ELSE CALL BGLP2 END IF END IF RETURN 9 I1=0 I2=0 L1=0 L2=0 X0=0 XE=0 A0=0 A1=0 E0=0 R=0 SIGR=0 RETURN END C----------------------------------------------------------------------- SUBROUTINE PZ (TH,VT,FP,VP) C C POLARIZATION FACTOR. WITH OR WITHOUT INCIDENT BEAM CRYSTAL C MONOCHROMATOR C COMMON /BLOCK4/ C0K,C1K,C2K,C0D,C1D,C2D,FD,VF,INEUTRON C C C0K = COS(2*THM)**2 KINEMATIC (MOSAIC) MONOCHROMATOR CRYSTAL C C0D = ABS(COS(2*THM)) DYNAMIC (PERFECT) MONOCHROMATOR CRYSTAL C C C1 = COS(RHO)**2 + SIN(RHO)**2*C0 C C2 = SIN(RHO)**2 + COS(RHO)**2*C0 C C RHO = 0 PARALLEL MONOCHROMATOR-DIFFRACTOMETER GEOMETRY C RHO = 90 DEG PERPENDICULAR GEOMETRY C C FD = 0 ALL KINEMATIC, NO DYNAMIC MONOCHROMATOR DIFFRACTION C FD = 1 ALL DYNAMIC, NO KINEMATIC MONOCHROMATOR DIFFRACTION C C VF = VAR(FD), THE VARIANCE OF THE DYNAMIC DIFFRACTION FRACTION C C FOR AN UNPOLARIZED INCIDENT BEAM (I.E., NO MONOCHROMATOR), C THM = 0, C0 = C1 = C2 = 1, AND FD = VF = 0. C C FOR NEUTRON DIFFRACTION DATA, INEUTRON = 1. C IF (INEUTRON.EQ.1) THEN FP=1 VP=0 ELSE U=COS(2*TH) V=SIN(2*TH) PK=(C1K+C2K*U**2)/(1+C0K) VK=VT*(-4*C2K*U*V/(1+C0K))**2 IF (FD.EQ.0) THEN FP=PK VP=VK ELSE PD=(C1D+C2D*U**2)/(1+C0D) VD=VT*(-4*C2D*U*V/(1+C0D))**2 FP=FD*PD+(1-FD)*PK VP=FD**2*VD+(1-FD)**2*VK+(PD**2+PK**2)*VF END IF END IF RETURN END C----------------------------------------------------------------------- SUBROUTINE READ1 (IFILE,IEND,II,IH,IK,IL,ANGLES,WIDTH,SPEED,XTIME, & NSTEP,YY) DIMENSION ANGLES(4),YY(1) INTEGER*2 JH,JK,JL,JY(96) READ (IFILE,END=9) II,JH,JK,JL,ANGLES,WIDTH,SPEED,(X,I=1,5),XTIME, & JY IH=JH IK=JK IL=JL NSTEP=96 DO 1 I=1,NSTEP YY(I)=JY(I) 1 CONTINUE IEND=0 RETURN 9 IEND=1 RETURN END C----------------------------------------------------------------------- SUBROUTINE READR C C SUBROUTINE READR MUST RETURN: C C INTEGER*4 VARIABLES: C C II - MEASUREMENT SERIAL NUMBER C IH,IK,IL - REFLECTION INDICES C NSTEP - NUMBER OF SCAN STEPS C C REAL*4 ARRAYS: C C ANGLES(4) - DIFFRACTOMETER SETTING ANGLES C XX(NSTEP) - ROCKING-ANGLE STEPS C YY(NSTEP) - COUNT-RATE STEPS C VX(NSTEP) - ROCKING-ANGLE VARIANCES C VY(NSTEP) - COUNT-RATE VARIANCES C C REAL*4 VARIABLE: C C XTIME - RADIATION EXPOSURE TIME C C THE DIFFRACTOMETER ANGLES MUST BE IN DEGREES, AND THEIR ORDER IS C IMPORTANT: C C DIFF. ANGLES(1) (2) (3) (4) C C INT.TAB. TWO-THETA OMEGA CHI PHI C BUS.LEV. TWO-THETA OMEGA CHI PHI C P3 TWO-THETA OMEGA PHI CHI C CAD4 THETA PHI OMEGA KAPPA C C THE ROCKING-ANGLE STEPS MUST INCREASE (FROM NEGATIVE TO POSITIVE) C IN ORDER OF INCREASING ABSOLUTE VALUE OF THE BRAGG ANGLE, THETA. C C RECOMMENDED UNITS: C C ROCKING ANGLE - DEGREES (THETA) C COUNT RATES - COUNTS PER SECOND C EXPOSURE TIME - HOURS C PARAMETER (NMAX=100) COMMON /VTKB/ IVT,IKB,IC COMMON /BLOCK0/ LU1,LU2,LU3,ILP COMMON /BLOCK7/ NSTEP,XX(NMAX),YY(NMAX),VX(NMAX),VY(NMAX) COMMON /RFLN/ IEND,II,IH,IK,IL,ANGLES(4),WIDTH,SPEED,FBG, & BG1,SCAN,BG2,YNET,SIGY,XTIME IEND=0 CALL READ1 (LU1,IEND,II,IH,IK,IL,ANGLES,WIDTH,SPEED,XTIME,NSTEP, & YY) IF (IEND.NE.0) THEN REWIND LU1 WRITE (IVT,1000) RETURN END IF 1000 FORMAT (' END OF FILE. FILE REWOUND.') C C CONVERT FROM STEP-SCAN COUNT PROFILE TO COUNT-RATE VERSUS ROCKING- C ANGLE PROFILE. C CALL XYDATA (WIDTH,SPEED,NSTEP,XX,YY,VX,VY) RETURN END C----------------------------------------------------------------------- FUNCTION SINTHL(H,GINV) DIMENSION H(3),GINV(3,3) Q=0 DO 1 I=1,3 DO 1 J=1,3 Q=Q+H(J)*GINV(I,J)*H(I) 1 CONTINUE SINTHL=0.5*SQRT(Q) RETURN END C----------------------------------------------------------------------- SUBROUTINE SMOOTH (NSMOOTH,N,Y) C C SMOOTHING OF EQUALLY SPACED DATA BY FOURTH DIFFERENCES C C EQUIVALENT TO A FIVE-POINT LEAST-SQUARES FIT OF A PARABOLA, C C Y = A0 + A1*X + A2*X**2, C C TO THE (I - 2)ND, (I - 1)ST, ITH, (I + 1)ST, AND (I + 2)ND C POINTS, I.E., TO X = -2, -1, 0, 1, AND 2. C C C. LANCZOS (1956). APPLIED ANALYSIS, PP. 316-320. PRENTICE- C HALL, ENGLEWOOD CLIFFS, NEW JERSEY. C DIMENSION Y(N),T(-2:+2) IF (NSMOOTH.LE.0.OR.N.LT.6) RETURN DO 9 I=1,NSMOOTH DO 1 J=-2,+2 T(J)=Y(3+J) 1 CONTINUE D3=T(1)-3*T(0)+3*T(-1)-T(-2) D4=T(2)-4*T(1)+6*T(0)-4*T(-1)+T(-2) Y(1)=Y(1)+D3/5+3*D4/35 Y(2)=Y(2)-2*D3/5-D4/7 Y(3)=Y(3)-3*D4/35 DO 2 J=4,N-2 DO 3 K=-2,+1 T(K)=T(K+1) 3 CONTINUE T(2)=Y(J+2) D4=T(2)-4*T(1)+6*T(0)-4*T(-1)+T(-2) Y(J)=Y(J)-3*D4/35 2 CONTINUE D3=T(2)-3*T(1)+3*T(0)-T(-1) Y(N-1)=Y(N-1)+2*D3/5-D4/7 Y(N)=Y(N)-D3/5+3*D4/35 DO 9 J=1,N IF (Y(J).LT.0) Y(J)=0 9 CONTINUE RETURN END C----------------------------------------------------------------------- SUBROUTINE WCALC (X,I1,I2,CHI,PHI,MODEL,Q,T,F,TANTH,TANTHM, & SIGQ,SIGT,W,SIGW) C C CALCULATED PEAK WIDTHS C DIMENSION X(1),Q(3,3),SIGQ(3,3),Z(3) DATA RAD /0.017453293/ C C GET THE COMPONENTS ALONG CARTESIAN AXES FIXED IN THE CRYSTAL OF A UNIT C VECTOR, Z, NORMAL TO THE EQUATORIAL PLANE OF THE INCIDENT AND C DIFFRACTED BEAM VECTORS. C COSCH=COS(CHI*RAD) SINCH=SIN(CHI*RAD) COSPH=COS(PHI*RAD) SINPH=SIN(PHI*RAD) Z(1)=SINPH*SINCH Z(2)=COSPH*SINCH Z(3)=COSCH C C CALCULATE PEAK WIDTHS. C U=0 V=0 DO 10 I=1,3 DO 10 J=1,3 U=U+Z(I)*Z(J)*Q(I,J) V=V+(Z(I)*Z(J)*SIGQ(I,J))**2 10 CONTINUE IF (MODEL.EQ.1) THEN C C LORENTZIAN CONVOLUTION MODEL C W=SQRT(U)+T*TANTH SIGW=SQRT(V/(2*U)**2+(SIGT*TANTH)**2) ELSE IF (MODEL.EQ.2) THEN C C GAUSSIAN CONVOLUTION MODEL C W=SQRT(U+(T*TANTH)**2) SIGW=SQRT(V+(2*T*TANTH**2*SIGT)**2)/(2*W) ELSE IF (MODEL.EQ.3) THEN C C GAUSSIAN CONVOLUTION MODEL FOR PARALLEL MONOCHROMATOR GEOMETRY C W=SQRT(ABS(U-2*F*T*TANTH/TANTHM+T*(TANTH/TANTHM)**2)) SIGW=SQRT(V/(2*U)**2+(2*F*SIGT*TANTH/TANTHM)**2 & +(SIGT*(TANTH/TANTHM)**2)**2) END IF RETURN END C----------------------------------------------------------------------- SUBROUTINE XCALC (X,I1,I2,H,U,SIGU,GINV,TWOTH,OMEGA,CHI,PHI, & X0,SIGX0) C C EXPECTED PEAK POSITION C DIMENSION X(1),H(3),U(3,3),SIGU(3,3),GINV(3,3),Y(3),V(3) DATA RAD,DEG /0.017453293, 57.2957795/ C C GET THE COMPONENTS ALONG CARTESIAN AXES FIXED IN THE CRYSTAL OF A C UNIT VECTOR, Y, PARALLEL TO THE RECIPROCAL LATTICE VECTOR, DSTAR, C IN THE DIFFRACTING CONDITION. C DO 1 I=1,3 Y(I)=0 V(I)=0 DO 1 J=1,3 Y(I)=Y(I)+H(J)*U(J,I) V(I)=V(I)+(H(J)*SIGU(J,I))**2 1 CONTINUE DSTAR=2*SINTHL(H,GINV) DO 2 I=1,3 Y(I)=Y(I)/DSTAR V(I)=V(I)/DSTAR**2 2 CONTINUE IF (TWOTH.LT.0) THEN DO 3 I=1,3 Y(I)=-Y(I) 3 CONTINUE END IF C C CALCULATE EXPECTED PEAK DISPLACEMENT, DELTA(OMEGA). C COSPH=COS(PHI*RAD) SINPH=SIN(PHI*RAD) COSCH=COS(CHI*RAD) SINCH=SIN(CHI*RAD) SINOM=Y(1)*COSPH-Y(2)*SINPH COSOM=(Y(1)*SINPH+Y(2)*COSPH)*COSCH-Y(3)*SINCH DELTA=ARCTAN(SINOM,COSOM)*DEG-OMEGA DELTA=A180(DELTA) IF (TWOTH.LT.0) DELTA=-DELTA C C THE POSITIVE SENSE OF THE OMEGA ROTATION HAS BEEN DEFINED SUCH C THAT DELTA(X0) = -DELTA(OMEGA). C X0=0.5*(X(I1)+X(I2))-DELTA VARSIN=V(1)*COSPH**2+V(2)*SINPH**2 VARCOS=(V(1)*SINPH**2+V(2)*COSPH**2)*COSCH**2+V(3)*SINCH**2 OMEGA=ARCTAN(SINOM,COSOM) VAROM=(COS(OMEGA)**2/COSOM)**2*(TAN(OMEGA)**2*VARCOS+VARSIN) SIGX0=SQRT(VAROM)*DEG RETURN END C----------------------------------------------------------------------- SUBROUTINE XCHECK (XE,X1,X0,X2,X,Y,I1,I2) C C FOR WEAK, HIGH-ANGLE REFLECTIONS, SUBROUTINE CENTROID SOMETIMES C FINDS X0 TOO CLOSE TO THE ALPHA-ONE PEAK. ENSURE THAT X1 IS AT C THE ALPHA-ONE PEAK. C DIMENSION X(1),Y(1) D=(X(I2)-X(I1))/(I2-I1+1) IF (X0-X1.LE.D) RETURN J1=NINDEX(X1) J0=NINDEX(X0) J1=MAX(J1,I1+1) J0=MIN(J0,I2-1) IF (J0-J1.LE.1) RETURN YMAX=0 DO 1 I=J1,J0 YI=Y(I-1)+Y(I)+Y(I+1) IF (YI.GT.YMAX) THEN YMAX=YI IMAX=I END IF 1 CONTINUE D=X(IMAX)-X1 IF (D.GT.0) THEN XE=XE+D X1=X1+D X0=X0+D X2=X2+D END IF RETURN END C----------------------------------------------------------------------- SUBROUTINE XYDATA (WIDTH,SPEED,NSTEP,XX,YY,VX,VY) C C CONVERT FROM STEP-SCAN COUNT PROFILE TO COUNT-RATE VERSUS ROCKING- C ANGLE PROFILE, AND APPLY COINCIDENCE CORRECTIONS FOR COUNTING LOSSES. C C RETURN: C C ROCKING-ANGLE STEPS XX(NSTEP) IN DEGREES (THETA) INCREASING C (FROM NEGATIVE TO POSITIVE) IN ORDER OF INCREASING ABSOLUTE C VALUE OF THE BRAGG ANGLE, THETA; C C COUNT-RATE STEPS YY(NSTEP) IN COUNTS PER SECOND; C C VARIANCE STEPS VX(NSTEP) AND VY(NSTEP). C C THIS VERSION OF SUBROUTINE XYDATA TREATS OMEGA, OMEGA/THETA, OR C OMEGA/TWO-THETA STEP SCANS WITH: C C NSTEP EQUALLY SPACED STEPS, C ----- C TOTAL SCAN WIDTH IN DEGREES OF OMEGA ROTATION, C ----- C SCAN SPEED IN DEGREES (OMEGA) PER MINUTE. C ----- C C NEGATIVE SCAN WIDTH IS A FLAG THAT THE ATTENUATOR WAS USED. C DIMENSION XX(1),YY(1),VX(1),VY(1) COMMON /BLOCK6/ VARTHE,VARTIM,TAU,VARTAU,ATT,VARATT C C ROCKING-ANGLE ABSCISSAE C X0=0.5*(1+NSTEP) DX=ABS(WIDTH)/NSTEP VARD=2*VARTHE/NSTEP**2 DO 1 I=1,NSTEP XX(I)=(I-X0)*DX VX(I)=(I-X0)**2*VARD 1 CONTINUE C C COUNT-RATE ORDINATES C SPEED=SPEED/60 T=(ABS(WIDTH)/NSTEP)/SPEED VART=2*(2*VARTHE/SPEED**2+VARTIM)/NSTEP**2 DO 2 I=1,NSTEP Y=YY(I) C C ALLOW FOR STEPS WITH VERY FEW OR NO COUNTS. C Y=Y+0.5 VARY=Y X=Y/T VARX=(VARY+X**2*VART)/T**2 C C DEAD TIME CORRECTION C Y=X/(1-TAU*X) VARY=Y**4*(VARX/X**4+VARTAU) C C ATTENUATOR CORRECTION C IF (WIDTH.LT.0) THEN X=Y VARX=VARY Y=ATT*X VARY=X**2*VARATT+ATT**2*VARX END IF YY(I)=Y VY(I)=VARY 2 CONTINUE RETURN END C-----------------------------------------------------------------------