[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: [ccp4bb]: SFALL and anomalous dispersion
Somewhat laborious - we can calculat Fano from the anomalous scatterers
then add together the real and inmaginary SFS with appropriate
multiplication factors
eg F(+) = Ftot -f'*Fano +iF" Fano
F(-) = Ftot* -f'Fano* +if"Fano* where F* is the complex conjugate
It uses a jiffy called mtzadd.f
I will attach scripts - all pre GUI days!
Elranor
And mtzadd.f - not sure if it is in distribution..
Bernard Santarsiero wrote:
> *** For details on how to be removed from this list visit the ***
> *** CCP4 home page http://www.ccp4.ac.uk ***
>
> Is there a way to calculate F(+) and F(-) or DANO in SFALL or CPP4?
>
> Bernie Santarsiero
>
>
#!/bin/csh -f
#
# Calculate structure factors for all Atase atoms include Ses;
# standard formfactors
goto 2
sfall \
HKLIN /x/people/yao/atase/ecse3_sort_trun_l123sc.mtz \
HKLOUT /y/ccp4/work/1ecf_x0.5_4.mtz \
XYZIN ./1ecf_x0.5_4.pdb \
<< 'END-sfrkall'
TITL Structure factors calculed for all SD atoms constant formactor
MODE SFCALC XYZIN HKLIN
RESO 37 2.25
BINS 60
LABI FP=FSe1 SIGFP=SIGFSe1
LABO ALLIN FC=FC_x0.5_4 PHIC=AC_x0.5_4
END
'END-sfrkall'
#
#
# Get structure factors for the Se alone, using a Unitary structure factor
# Extract from Se pdb file
# CRYST1 116.900 157.500 106.300 90.00 90.00 90.00 C 2 2 21 16
# SCALE1 0.008554 0.000000 0.000000 0.00000
# SCALE2 0.000000 0.006349 0.000000 0.00000
# SCALE3 0.000000 0.000000 0.009407 0.00000
# ATOM 67 An MET A 10 1 27.688 55.343 19.523 1.00 27.19
# ATOM 455 An MET A 60 1 14.303 49.818 24.775 1.00 17.68
# Extract from atomsf.lib for the unitary formfactor
#
# An 2
# 102 102 0.000000
# 1.000000 0.000000 0.000000 0.000000
# 0.000000 0.000000 0.000000 0.000000
# 0.000000 1.000000 0.000000 1.000000
sfall \
HKLIN /y/ccp4/work/1ecf_x0.5_4.mtz \
HKLOUT /y/ccp4/work/1ecf_sd_x0.5.mtz \
XYZIN ./1ecf_sd_x0.5_4.pdb \
ATOMSF /y/ccp4/atase/atomsf.lib \
<< 'END-sfrkall'
TITL Structure factors calculed for all SD atoms constant formactor
MODE SFCALC XYZIN HKLIN
RESO 37 2.25
NOSCALE
BINS 60
LABI FP=FSe1 SIGFP=SIGFSe1
LABO ALLIN FC=FC_seu PHIC=AC_seu
END
'END-sfrkall'
#
#
#!/bin/csh -f
#
2:
/y/ccp4/cadd/junk/mtzadd \
hklin /y/ccp4/work/1ecf_sd_x0.5.mtz \
hklout /y/ccp4/work/junk1.mtz \
<<eof
LABI FC1=FC_x0.5_4 PHIC1=AC_x0.5_4 -
FC2=FC_seu PHIC2=AC_seu -
FC3=FC_seu PHIC3=AC_seu
LABO FSUM = FCL1_x0.5_4(+) PHISUM=PHICL1_x0.5_4(+)
SCALE FC1 1 0 # Scale Bfactor
SCALE FC2 -5.5 0 # Scale Bfactor
SCALE COMP FC3 1.3 0 # Complex ( ie advanvce phase by 90) Scale Bfactor
END
eof
#
/y/ccp4/cadd/junk/mtzadd \
hklin /y/ccp4/work/junk1.mtz \
hklout /y/ccp4/work/junk2.mtz \
<<eof
#LABI FC1=FC_x0.5_4 PHIC1=AC_x0.5_4 FC2=FC_seu PHIC2=AC_seu
LABI FC1=FC_x0.5_4 PHIC1=AC_x0.5_4 -
FC2=FC_seu PHIC2=AC_seu -
FC3=FC_seu PHIC3=AC_seu
LABO FSUM = FCL1_x0.5_4(-) PHISUM=PHICL1_x0.5_4(-)
SCALE FC1 1 0
SCALE FC2 -5.5 0
SCALE COMP FC3 -1.3 0
END
eof
#
/y/ccp4/cadd/junk/mtzadd \
hklin /y/ccp4/work/junk2.mtz \
hklout /y/ccp4/work/junk3.mtz \
<<eof
LABI FC1=FC_x0.5_4 PHIC1=AC_x0.5_4 -
FC2=FC_seu PHIC2=AC_seu -
FC3=FC_seu PHIC3=AC_seu
LABO FSUM = FCL2_x0.5_4(+) PHISUM=PHICL2_x0.5_4(+)
SCALE FC1 1 0
SCALE FC2 -4.5 0
SCALE COMP FC3 2.4 0
END
eof
#
/y/ccp4/cadd/junk/mtzadd \
hklin /y/ccp4/work/junk3.mtz \
hklout /y/ccp4/work/junk4.mtz \
<<eof
LABI FC1=FC_x0.5_4 PHIC1=AC_x0.5_4 -
FC2=FC_seu PHIC2=AC_seu -
FC3=FC_seu PHIC3=AC_seu
LABO FSUM = FCL2_x0.5_4(-) PHISUM=PHICL2_x0.5_4(-)
SCALE FC1 1 0
SCALE FC2 -4.5 0
SCALE COMP FC3 -2.4 0
END
eof
#
/y/ccp4/cadd/junk/mtzadd \
hklin /y/ccp4/work/junk4.mtz \
hklout /y/ccp4/work/junk5.mtz \
<<eof
LABI FC1=FC_x0.5_4 PHIC1=AC_x0.5_4 -
FC2=FC_seu PHIC2=AC_seu -
FC3=FC_seu PHIC3=AC_seu
LABO FSUM = FCL3_x0.5_4(+) PHISUM=PHICL3_x0.5_4(+)
SCALE FC1 1 0
SCALE FC2 -2.0 0
SCALE COMP FC3 1.4 0
END
eof
#
/y/ccp4/cadd/junk/mtzadd \
hklin /y/ccp4/work/junk5.mtz \
hklout /y/ccp4/work/junk6.mtz \
<<eof
LABI FC1=FC_x0.5_4 PHIC1=AC_x0.5_4 -
FC2=FC_seu PHIC2=AC_seu -
FC3=FC_seu PHIC3=AC_seu
LABO FSUM = FCL3_x0.5_4(-) PHISUM=PHICL3_x0.5_4(-)
SCALE FC1 1 0
SCALE FC2 -2.0 0
SCALE COMP FC3 -1.4 0
END
eof
#
/y/ccp4/cadd/junk/mtzMADmod hklin $SCRATCH/junk6.mtz \
hklout ./1ecf_x0.5_4_L1L2L3.mtz \
<<eof
LABI -
F1+=FCL1_x0.5_4(+) SIGF1+=SIGFSe(+)1 F1-=FCL1_x0.5_4(-) SIGF1-=SIGFSe(-)1 -
F2+=FCL2_x0.5_4(+) SIGF2+=SIGFSe(+)2 F2-=FCL2_x0.5_4(-) SIGF2-=SIGFSe(-)2 -
F3+=FCL3_x0.5_4(+) SIGF3+=SIGFSe(+)3 F3-=FCL3_x0.5_4(-) SIGF3-=SIGFSe(-)3
LABO -
F1=FCL1_x0.5_4 SIGF1=SIGFCL1_x0.5_4 D1=DCL1_x0.5_4 SIGD1=SIGDCL1_x0.5_4 -
F2=FCL2_x0.5_4 SIGF2=SIGFCL2_x0.5_4 D2=DCL2_x0.5_4 SIGD2=SIGDCL2_x0.5_4 -
F3=FCL3_x0.5_4 SIGF3=SIGFCL3_x0.5_4 D3=DCL3_x0.5_4 SIGD3=SIGDCL3_x0.5_4
END
eof
exit
C
C This code is distributed under the terms and conditions of the
C CCP4 licence agreement as `Part ii)' software. See the conditions
C in the CCP4 manual for a copyright statement.
C
C Authors: EJD, MDW
C
C mtzmnf.f
C STANDARD ITEM NAMES taken from $CDOC/mtzlib.doc
C Standard names are normally used for the following items:
C
C Name Item
C
C H, K, L Miller indices.
C .....
C I Intensity.
C SIGI sigma(I).
C .....
C FP Native 'F' value.
C FC Calculated 'F'.
C FPHn 'F' value for derivative 'n'.
C
C DP Anomalous difference for native data.
C DPHn Anomalous difference for derivative 'n'.
C
C SIGFP sigma(FP).
C SIGDP sigma(DP).
C SIGFPHn sigma(Fn).
C SIGDPHn sigma(DELn).
C
C PHIC Calculated Phase.
C PHIM Most probable phase.
C PHIB Phase.
C
C FOM figure of merit.
C WT weight
C
CCOLUMN TYPES
C The column types are as follows:
C
C H index h,k,l
C J intensity
C F structure amplitude, F
C D anomalous difference
C Q standard deviation of anything: J,F,D or other
C P phase angle in degrees
C W weight (of some sort )
C A phase probability coefficients (Hendrickson/Lattmann)
C .....
C I any other integer
C R any other real
C*************************************************************************************
C--------- description of Key_Worded Input --------------
C ================================
C
C
C (i) TITLE (optional)
C
C Replaces title in MTZ file
C
C (ii) LABIN (compulsory)
C
C Examples:
C LABI F1=FTOXD3 SIGF1=SIGFTOXD3 -
C F2=FAU20 SIGF2=SIGFAU20 -
C D2=ANAU20 SIGD2=SIGANAU20 -
C
C (iii) END (optional)
C
C Ends keyword input
C
PROGRAM MAIN
C ..
C .. Parameters ..
INTEGER MCOLS,NKEYWD,LSYM
PARAMETER (MCOLS=200,NKEYWD=200,LSYM=192)
C ..
C .. Local Arrays ..
REAL ADATA(MCOLS),BDATA(MCOLS),RSYML(4,4,LSYM),
+ SCL(MCOLS),B(MCOLS),COMPL(MCOLS)
INTEGER LOOKUP(MCOLS),IN(3),NCHANGES(MCOLS)
CHARACTER*1 CTPRGI(MCOLS),CTUSRI(MCOLS)*1,CTPRGO(2)
CHARACTER*30 LSPRGI(MCOLS),CLABS(MCOLS),LSUSRI(MCOLS)*30,
+ LSPRGO(2)
LOGICAL LOGMSS(MCOLS),LMSSOLD(MCOLS)
C ..
C .. Local Scalars ..
REAL S
INTEGER I,IC,IFAIL,IFIL1,IFLAG,IPRINT,JDO30,JDO35,JDO40,JDO60
INTEGER JDO65,JDO70,JDO80,JDO90,LUNOUT,NCOL,NF1,NLPRGI,NSPGRP
INTEGER NSYM,NSYMP,NUMFO,NUMFC,NUMI,NUMPH,NTOK,STARTFO,STARTFC
INTEGER STARTI,STARTPH
CHARACTER NAMSPG*10,PGNAM*10,TITLE*80
CHARACTER LTYPEX*1
LOGICAL ENDFIL,LTITLE
C ..
C .. Variables for Parser ..
REAL FVALUE(NKEYWD)
INTEGER IBEG(NKEYWD),IDEC(NKEYWD),IEND(NKEYWD),ITYPE(NKEYWD)
CHARACTER KEY*4,CVALUE(NKEYWD)*4,LINE*600,WORD*4
LOGICAL LEND
C ..
C .. External Routines ..
EXTERNAL CCPERR,CCPFYP,CCPRCS,CENTR,CENTRIC,EQUAL_MAGIC,LKYIN,
+ LROPEN,LRASSN,LRCLAB,LRCLOS,LRREFL,LRREFM,LRSYMI,LRSYMM,
+ LWOPEN,LWREFL,MTZINI,PARSER
C ..
C .. Data Statements ..
DATA NLPRGI/43/
DATA NLPRGO/2/
DATA LTITLE/.FALSE./
C Allow NUMFO sets of Fs and phases
C Up to NUMFC FCs
C Up to NUMI Is
DATA NUMFO,NUMFC,NUMI/5,5,5/
DATA LSPRGI/'H','K','L',
+ 'F1' ,'SIGF1' ,'PHIB1' ,'SIGD1' ,
+ 'F2' ,'SIGF2' ,'PHIB2' ,'SIGD2' ,
+ 'F3' ,'SIGF3' ,'PHIB3' ,'SIGD3' ,
+ 'F4' ,'SIGF4' ,'PHIB4' ,'SIGD4' ,
+ 'F5' ,'SIGF5' ,'PHIB5' ,'SIGD5' ,
+ 'F6' ,'SIGF6' ,'PHIB6' ,'SIGD6' ,
+ 'FC1','PHIC1' ,
+ 'FC2','PHIC2' ,
+ 'FC3','PHIC3' ,
+ 'FC4','PHIC4' ,
+ 'FC5','PHIC5' ,
+ 'I1' ,'SIGI1' ,
+ 'I2' ,'SIGI2' ,
+ 'I3' ,'SIGI3' ,
+ 'I4' ,'SIGI4' ,
+ 'I5' ,'SIGI5' ,
+ 153*' '/
DATA CTPRGI/'H','H','H',
+ 'F','Q','P','W',
+ 'F','Q','P','W',
+ 'F','Q','P','W',
+ 'F','Q','P','W',
+ 'F','Q','P','W',
+ 'F','Q','P','W',
+ 'F','P',
+ 'F','P',
+ 'F','P',
+ 'F','P',
+ 'F','P',
+ 'J','Q',
+ 'J','Q',
+ 'J','Q',
+ 'J','Q',
+ 'J','Q',
+ 153*' '/
DATA LSPRGO/'FSUM','PHISUM'/
DATA CTPRGO/'F','P'/
DATA LOOKUP/-1,-1,-1,0,0,0,0,0,0,0,0,189*0/
DATA NCHANGES/200*0/
C
C Starting indices of program label blocks
STARTFO=4
STARTFC=NUMFO*4 + STARTFO
STARTI=NUMFC*2 + STARTFC
LUNOUT=6
IFIL1 = 1
NREF = 0
PI=3.14159265
DTOR=PI/180.
CALL CCPFYP
CALL MTZINI
CALL CCPRCS(6,'MTZADD','$Date: 1996/04/15 15:36:46 $')
C ...Read keyword input
20 CONTINUE
KEY = ' '
LINE = ' '
C Reset number of words on line from 20 to NKEYWD
C Need lots for LABI lines
NTOK = NKEYWD
C
C **********************************************************
CALL PARSER(KEY,LINE,IBEG,IEND,ITYPE,FVALUE,CVALUE,IDEC,NTOK,LEND,
+ .TRUE.)
C **********************************************************
C
IF (LEND) GO TO 200
IF(KEY.EQ.'TITL') THEN
IF (NTOK.GE.2) THEN
LTITLE = .TRUE.
TITLE = LINE(IBEG(2) :IEND(NTOK))
WRITE(LUNOUT,'(10X,a,/,15X,a80,/)')
+ 'TITLE OF MTZ FILE REPLACED BY:',TITLE
ENDIF
GO TO 20
ELSE IF(KEY.EQ.'LABI') THEN
C---- LABIN i.e. set input labels
C *****************************************
CALL LKYIN(IFIL1,LSPRGI,NLPRGI,NTOK,LINE,IBEG,IEND)
C *****************************************
GO TO 20
ELSE IF(KEY.EQ.'LABO') THEN
C---- LABOUT i.e. set output labels
C *****************************************
CALL LKYOUT(IFIL1,LSPRGO,NLPRGO,NTOK,LINE,IBEG,IEND)
C *****************************************
GO TO 20
C---- SCALE [Complex] Label Scale B-factor
C ===========================
C
C F1 or F2 Scale and temperature factor
C
C
ELSE IF(KEY.EQ.'SCAL') THEN
ITOK = 2
COMPLEX = 0.0
WORD = LINE(IBEG(2):IBEG(2))
CALL CCPUPC(WORD)
IF(WORD(1:1) .EQ.'C' ) THEN
COMPLEX = 1.0
ITOK = 3
END IF
5 IF (ITOK .LE. NTOK .AND. ITOK .GT. 0) THEN
C
C---- Subroutine RDSCAL extracts from line, beginning at field ITOK,
C sets of 3 arguments:- Label (= F1 or F2); Scale; Bfactor
C
C Input: LSPRGI list of NLPRGI possible labels
C Output: ILPRGI label number found
C SCAL scale
C BB temperature factor
C ITOK next item to read, = 0 if finished, .lt. 0 if error
C
I = ITOK
C ***************************************************
CALL RDSCAL(ITOK,LINE,IBEG,IEND,ITYPE,FVALUE,NTOK,NLPRGI,
+ LSPRGI,ILPRGI,SCAL,BB)
C ***************************************************
IF (ITOK .EQ. I) THEN
C Allow for old version of RDSCAL which doesn't reset ITOK
ITOK = ITOK + 3
ELSEIF (ITOK .EQ. -1) THEN
WRITE (LUNOUT, '(A)')
$ ' **** No scale factor given ****'
KE = KE + 1
ELSEIF (ITOK .EQ. -2) THEN
WRITE (LUNOUT, '(A)')
$ ' **** Unrecognized program label ****'
KE = KE + 1
ENDIF
C
IF (ITOK .GE. 0) THEN
SCL(ILPRGI) = SCAL
B(ILPRGI) = BB
COMPL(ILPRGI) = COMPLEX
WRITE (LUNOUT,FMT=6004),ILPRGI,
$ LSPRGI(ILPRGI)(1:LENSTR(LSPRGI(ILPRGI))),
+ SCL(ILPRGI),B(ILPRGI),COMPL(ILPRGI)
6004 FORMAT (//
$ ' Scale factor and Temperature factor for',7X,
$ I4,A,4X,3F9.4)
GO TO 5
ENDIF
ENDIF
GO TO 20
C
ELSE IF(KEY.EQ.'END') THEN
WRITE(LUNOUT,'(/)')
GO TO 200
ELSE
Call ccperr(1,' Key word unrecognised')
END IF
C
C---- End of keyword input
200 CONTINUE
C
IPRINT = 2
C **********************************
CALL LROPEN(IFIL1,'HKLIN',IPRINT,IFAIL)
C **********************************
IF(IFAIL.NE.0) THEN
CALL CCPERR(1, ' No input file???????')
END IF
C
C ****************************************
CALL LRASSN(IFIL1,LSPRGI,NLPRGI,LOOKUP,CTPRGI)
C *****************************************
C
C =========================================
CALL LRCLAB(IFIL1,CLABS,CTUSRI,NCOL)
C =========================================
C
C Test correct column types.
NF1 = 0
NF2 = 0
DO 30 JDO30 = STARTFO,STARTFO+(NUMFO-1)*4,4
NF1 = NF1 + 1
NF2 = NF2 + LOOKUP(JDO30)
IF (LOOKUP(JDO30).NE.0 .NEQV. LOOKUP(JDO30+1).NE.0) THEN
WRITE (LUNOUT,FMT='(/A,I2,A,I2,A/)')
+ 'Unmatched F',NF1,' or SIGF',NF1,' assigned.'
CALL CCPERR(1,' Unmatched F/SIGF pair')
ENDIF
IF (LOOKUP(JDO30+2).NE.0 .NEQV. LOOKUP(JDO30+3).NE.0) THEN
WRITE (LUNOUT,FMT='(/A,I2,A,I2,A/)')
+ 'Unmatched PHIB',NF1,' or FOM',NF1,' assigned.'
CALL CCPERR(1,' Unmatched PHI/FOM pair')
ENDIF
IF (LOOKUP(JDO30+2).NE.0 .AND. LOOKUP(JDO30).EQ.0) THEN
WRITE (LUNOUT,FMT='(/A,I2,A,I2/)')
+ 'PHI',NF1,' assigned but no corresponding F',NF1
CALL CCPERR(1,' Unmatched F/PHI pair')
ENDIF
30 CONTINUE
NF1 = 0
DO 40 JDO40 = STARTFC,STARTFC+(NUMFC-1)*2,2
NF1 = NF1 + 1
NF2 = NF2 + LOOKUP(JDO40)
IF (LOOKUP(JDO40).NE.0 .NEQV. LOOKUP(JDO40+1).NE.0) THEN
WRITE (LUNOUT,FMT='(/A,I2,A,I2,A/)')
+ 'Unmatched FC',NF1,' or PHIC',NF1,' assigned.'
CALL CCPERR(1,' Unmatched FC/PHIC pair')
ENDIF
40 CONTINUE
NI1 = 0
NI2 = 0
DO 35 JDO35 = STARTI,STARTI+(NUMI-1)*2,2
NI1 = NI1 + 1
NI2 = NI2 + LOOKUP(JDO35)
IF (LOOKUP(JDO35).NE.0 .NEQV. LOOKUP(JDO35+1).NE.0) THEN
WRITE (LUNOUT,FMT='(/A,I2,A,I2,A/)')
+ 'Unmatched I',NI1,' or SIGI',NI1,' assigned.'
CALL CCPERR(1,' Unmatched I/SIGI pair')
ENDIF
35 CONTINUE
C
C
IF ( NF2 .GT.0 .AND. NI2 .GT.0)
+ CALL CCPERR(1,' Cannot mix Fs and Is')
C
C **********************
CALL LWOPEN(IFIL1,'HKLOUT')
C **********************
C
C Append output label for FSUM and PHISUM
C **********************
CALL LWASSN(IFIL1,LSPRGO,2,CTPRGO,1)
C **********************
C
C---- Iflag = 0 means replace old title.
C
ICHK = 1
IFLAG = 0
IF (LTITLE) CALL LWTITL(IFIL1,TITLE,IFLAG)
C
C---- Loop for reflections----------------------------------------------------
C
10 CONTINUE
C
C **************************
CALL LRREFL(IFIL1,S,ADATA,ENDFIL)
IF (ENDFIL) GO TO 100
CALL LRREFM(IFIL1,LOGMSS)
C **************************
C
CALL EQUAL_MAGIC(IFIL1,ADATA(NCOL+1),2)
IN(1) = NINT(ADATA(1))
IN(2) = NINT(ADATA(2))
IN(3) = NINT(ADATA(3))
C ICHK = MOD(NREF,1000)
FSUM = 0
PHISUM = 0.0
ASUM = 0
BSUM = 0.0
C Test whether all terms there - if not will have to omit reflection
C
C
C case (a): F, SIGF
DO 60 JDO60 = STARTFO,STARTFO+NUMFO*4-3,4
IF(LOOKUP(JDO60).EQ.0) GO TO 60
C ...if SIGF already flagged as "Missing"
IF(LOGMSS(LOOKUP(JDO60)) .OR. LOGMSS(LOOKUP(JDO60+1)) .OR.
+ LOGMSS(LOOKUP(JDO60+2)) .OR. LOGMSS(LOOKUP(JDO60+3))) GO TO 10
FOM = ADATA(LOOKUP(JDO60+3))
F1=FOM * ADATA(LOOKUP(JDO60))*SCL(JDO60)*EXP(B(JDO60)*S/4.0)
PHI1=ADATA(LOOKUP(JDO60+2))
IF(COMPL(JDO60).NE.0.0) PHI1 = PHI1 + COMPL(JDO60)*90
PHI1=DTOR*PHI1
ASUM = ASUM + F1*COS(PHI1)
BSUM = BSUM + F1*SIN(PHI1)
IF(ICHK.EQ.0)write(6,'(3i4,7f9.3)')in,ADATA(LOOKUP(JDO60)),
+ ADATA(LOOKUP(JDO60+2)),ADATA(LOOKUP(JDO60+3)),F1,PHI1,ASUM,BSUM
60 CONTINUE
C case (c): FC, SIGFC
DO 70 JDO70 = STARTFC,STARTFC+NUMFC*2-2,2
IF(LOOKUP(JDO70).EQ.0) GO TO 70
IF(LOGMSS(LOOKUP(JDO70)) .OR. LOGMSS(LOOKUP(JDO70+1)))
+ GO TO 10
F1= ADATA(LOOKUP(JDO70))*SCL(JDO70)*EXP(B(JDO70)*S/4.0)
PHI1=ADATA(LOOKUP(JDO70+1))
IF(COMPL(JDO70).NE.0.0) PHI1 = PHI1 + COMPL(JDO70)*90
PHI1=DTOR*PHI1
ASUM = ASUM + F1*COS(PHI1)
BSUM = BSUM + F1*SIN(PHI1)
IF(ICHK.EQ.0)write(6,'(3i4,6f9.3)')in,ADATA(LOOKUP(JDO70)),
+ ADATA(LOOKUP(JDO70+1)),F1,PHI1,ASUM,BSUM
70 CONTINUE
C case (d): I, SIGI
DO 80 JDO80 = STARTI+1,STARTI+NUMI*2-1,2
IF(LOOKUP(JDO80).EQ.0) GO TO 80
IF(LOGMSS(LOOKUP(JDO80)) .OR. LOGMSS(LOOKUP(JDO80+1)))
+ GO TO 10
F1= ADATA(LOOKUP(JDO80))*SCL(JDO80)*EXP(B(JDO80)*S/4.0)
PHI1=0.0
ASUM = ASUM + F1*COS(PHI1)
BSUM = BSUM + F1*SIN(PHI1)
80 CONTINUE
FSUM = SQRT(ASUM*ASUM + BSUM*BSUM)
PHISUM=0.0
IF (FSUM.gt.0.0001) PHISUM=ATAN2(BSUM,ASUM)/DTOR
ADATA(NCOL+1)=FSUM
ADATA(NCOL+2)=PHISUM
IF(ICHK.EQ.0)write(6,'(3i4,6f9.3)')in,ASUM,BSUM,FSUM,PHISUM
NREF = NREF + 1
CALL LWREFL(IFIL1,ADATA)
C ...Return for next reflection
GOTO 10
100 CONTINUE
CALL LRCLOS(IFIL1)
CALL LWCLOS(IFIL1,1)
WRITE(LUNOUT,'(//,A,/,A,//)')
+ ' CHANGES MADE TO MTZ FILE',
+ ' ========================'
WRITE(LUNOUT,'(/)')
CALL CCPERR(0,'Normal termination')
END