NBSIR 76-1102

Computer Programs for

Structural Chemistry: Status.

A Fortran Program for

Statistical Analysis of

Crystallographic Quantities

LeRoy W. Schroeder* and Brian Dickens

*Research Associate of the American Dental Assoc.

Health Foundation Research Unit at the National

Bureau of Standards

January 1976

Institute for Materials Research

National Bureau of Standards

Washington, D. C. 20234

Issued July, 1 976

U. S. DEPARTMENT OF COMMERCE

NATIONAL BUREAU OF STANDARDS
-
NBSIR 76-1102

COMPUTER PROGRAMS FOR

STRUCTURAL CHEMISTRY: STATUS.

A FORTRAN PROGRAM FOR

STATISTICAL ANALYSIS OF

CRYSTALLOGRAPHIC QUANTITIES

LeRoy W. Schroeder* and Brian Dickens

’Research Associate of the American Dental Assoc.

Health Foundation Research Unit at the National

Bureau of Standards

Institute for Materials Research

National Bureau of Standards

Washington, D. C. 20234

January 1 976

Issued July, 1 976

U S. DEPARTMENT OF COMMERCE, Elliot L. Richardson, Secretary

Edward O. Vetter, Under Secretary

Dr. Betsy Ancker-Johnson, Assistant Secretary for Science and Technology

NATIONAL BUREAU OF STANDARDS, Ernest Ambler, Acting Director

This investigation was supported in part byGrant DE 00572 to the American Dental Assoc.

Health Foundation from the National Institute

of Dental Research and is part of the dental

research program conducted by the National

Bureau of Standards in cooperation with the

American Dental Association Health Foundation.

Copies of this report and NBS Magnetic Tape 11

containing card images of the program described

herein can both be obtained from:

National Technical Information Service

Department of Commerce

5285 Park Royal Road

Springfield, Virginia 22151

Abstract

This report describes a FORTRAN computer program for

evaluation of (i) the results of crystallographic leastsquares refinements by examination of the residuals

6 Ri = (Fo-Fc) /a (Fo) , (ii) the differences in sets of data

collected by different methods from the same crystal and

the appropriateness of the assigned standard errors by

examination of the statistics

6m i = [FOi (set 1) - kFo^set 2)1 /[Cj 3 (set 1) + k

2 a } 2 (set 2 ) ] s

, (iii) differences in data sets collected by the same method

from different crystals of the same material again using

the 6m statistics, and (iv) the significance of differences

in parameters in different models representing the crystal

structure of the same material by examination of the

statistics

Spi = [pi (set 1) - pi (set 2)] / [rr 2 pi(set 1) + a

2

p 4 (set 2)1^.

Procedure (i) provides diagnostic tests of the overall fit

between observed and calculated crystallographic quantities,

procedure (ii) focuses attention on experimental methods,

procedure (iii) focuses attention on sample homogeneity,

and procedure (iv) examines the sensitivity of refinement

models to sample characteristics and experimental techniques.

Procedures (iii) and (iv) have been discussed by Abrahams

and Keve (1971) , who suggested that the examination be

accomplished by means of plots of residuals or statistics

against the expected normal distribution quantiles. Here

we have programmed their procedure and have extended the

treatment to two additional kinds of plots. One kind of

plot compares residuals with the independent variable d*

and the other kind compares residuals with the calculated

variable F&/sin20. In addition, we have used the Miller

indices to divide the data into classes and octants so that

the possibility of anisotropic effects can be examined.

The program will handle up to 7200 data points in each of

the two sets compared and is oriented specifically toward

examining diffraction data.
The FORTRAN coding has been designed to minimize

changes necessary to adapt this program to other computers.

All plots are made on the line printer for convenience and

speed. Only four input cards for program direction are

required in most cases. This report includes a description

of the general procedure, data specifications, program

logic, a listing of the FORTRAN code, and samples of input

and output.

Key words: Crystallographic data; errors; probability

plots; residuals; statistics; uncertainties

IV
Table of Contents

, Page

I. General

^ I. A. Introduction, Disclaimers and

Acknowledgments 1

I.B. Raison d

1 etre for STATUS 2

I.C. Application of Statistical Analysis to

Crystallography 5

1. Use of the 5mi Statistic 5

2 . Use of the Sr* Statistic 9

3. Use of the 5p< Statistic qq

I.D. General Description of STATUS 11

1. Overall Procedure 11

2 . Crystallographic Data Required as Input . . 13

3. Format of the Output Plot 14

4. Program Configuration and Core Requirements 15

5. Subroutine Hierarchy and Call Sequence . . 17

6 . Description of Subroutines 21

I. E. References 2 5

II. Program Manual 2 6

II. A. Card Input 26

II. B. Data Files 30

1. Input Files 30

2. Output Files 31

3. Temporary Files 31

4. File Assignments 31

II. C. Examples of Output Plots 32

1. Normal Probability Plot 32

2 , d* Plot 3 5

3. Log (Font's in2 8 ) Plot 42

4. Half-Normal Probability Plot 48

III. Detailed Description of Main Algorithms 51

III. A. Main Routine STATUS 51

1. Common Blocks and Arrays 52

2 . Flow Diagram of STATUS 55

v
Page

III. B. Subroutine ANISO 57

III. C. Subroutine DRUMRD and DRUMRT 59

III. D. Subroutine FITCHK 60

III. E. Subroutine FREAD 63

III. F. Subroutine INFOFC . 6 5

III. G. Subroutine LSFIT 67

III. H. Subroutine PLOTEM 69

III. I. Subroutine REREAD.. 72

III. J. Subroutine SCALE 7 4

IV. FORTRAN Code. ........... 77

IV. A. Introduction 77

IV. B. Features Specific to UNIVAC 1108 FORTRAN V • 77

IV. C. Listings 79

V. DESCRIPTION OF NBS MAGNETIC TAPE 11 144

vi
Computer Programs for Structural Chemistry: STATUS.

A FORTRAN Program for Statistical Analysis

of Crystallographic Quantities

LeRoy W. Schroeder and Brian Dickens

I. General

I. A. Introduction, Disclaimers and Acknowledgments

The text of this report describes the main features

of the FORTRAN program STATUS identified in Line 50 of the

main routine as the January 1975 version. The program elements

of STATUS as given in this report may be consulted for other

details. STATUS was written by LeRoy W. Schroeder and

Brian Dickens at the National Bureau of Standards, Washington,

D. C. 20234. This report is also intended as a program

manual to aid in the use and implementation of STATUS.

Questions concerning STATUS not covered in this report may

be addressed to the authors at the above address.

The authors have extensively tested the program on

the UNIVAC 1108 Exec 8 system but, of course, cannot

guarantee its performance.

LWS would like to acknowledge several helpful discussions

with Dr. James J. Filliben of the Statistical Engineering

Laboratory at NBS . Acknowledgment is due Dr. Filliben as author of subroutines SORT and UNIMED, Dr. W. V. Loebenstein

as author of SIMLEQ , R. J. Arms for subroutine PINV and

S. Peavy as author of PLOTS. The authors would also like to

acknowledge the able assistance of Pamela Kingsbury in

preparing this report.
I.B. Raison d'etre for STATUS

A recent editorial (Science, 1972) by J. R. MacDonald,

Chairman of the Numerical Data Advisory Board of the National

Research Council, was concerned with the question of the

trustworthiness of experimental data. MacDonald points out

that measures of uncertainty either are generally not given or are

unaccompanied by any estimate of their reliability.

Crystallographic structural studies provide details

of the average atomic arrangements, bond distances and angles,

impurity distributions and site populations, and can even be

used as a means, albeit expensive, of chemical analysis. From

this, one obtains basic structural information about molecular

and ionic geometries, ion packings, structural relationships

and so on. In addition, the fine details in the parameters

resulting from a structure determination provide information

about the extent of positional disorder in the structure,

ranges of impurity substitution in a given atomic arrangement,

and the means by which such substitution can affect the

physical and chemical properties of the host structure. Along

more interpretive lines, structural parameters have been used

as "data" in such areas of chemistry as crystal energy calculations, investigations of epitactical relationships and

elucidation of hydrogen bonding by combining diffraction and

spectroscopic results.

The structural parameters comprising the "model" of

the crystal structure are quantities estimated by non-linear

least-squares procedures. These procedures also provide

estimates of the uncertainties associated with the structural

parameters. However, the actual values of the derived parameters and their associated uncertainties depend on the magnitudes and uncertainties in the primary crystallographic data,

i.e., the structural amplitudes derived from measurements of

diffraction intensities from the specimen, assumed here to be

a single crystal. Two considerations therefore arise:

(1) The observed amplitudes should be made as free

of systematic error as possible and the uncertainties should

be derived from a variance that combines the mean-square

random error and any remaining mean- square systematic error,

and

(2) The derived model should be as free of

systematic effects as possible, and uncertainties in the

parameters of the model should include both random

error and any bias resulting from systematic effects.

An example of such an effect is thermal diffuse scattering,

for which corrections are difficult to apply.

2
Consideration (1) has generally been dealt with by

applying corrections and checking agreement between equivalent reflections. Crystallographically equivalent reflections

provide "repeat points" if anisotropic effects are neglected.

Abrahams (1964) has proposed a method for evaluating the

independent variances entering the measurement. Some of these

are apparatus dependent and not easily evaluated.

Consideration (2) has been partially dealt with by

referring to the standard and weighted factors,

defined as

R = S | | Fo | - | Fc | | /£ | Fo | and R = (£w ( | Fo| - |Fc | ) 2 /Sw Fo2

) 77

respectively, as criteria by which to judge the fit to the

diffraction data of quantities calculated from the model. The

R factors are insufficient for assessing (i) estimates of

errors in the model parameters and (ii) to what degree the

model has compensated for systematic errors in the measured

structural amplitudes. An attempt at (i) could be made by . comparing supposedly identical chemical entities such as a C-C bond length as calculated in different structures after

refinement. Such an attempt must assume that crystal packing

effects are negligible.

The statistical procedures outlined in the Abstract

and described further in Section I.C provide a means of

dealing with these considerations. The statistical

analysis follows standard procedures and might be handled by

a general program such as OMNITAB II (Hogben, Peavy and

Varner, 1971), designed for such analyses. However, the

generality of OMNITAB makes for awkward handling of large

quantities of data such as are found in crystallography.

Analysis of experimental crystallographic data requires at

least four numbers, h, k, i and a statistic, for each

reflection. A typical data set of 3000 reflections would

therefore require handling of about 12,000 numbers.

OMNITAB II was developed around the "worksheet" concept

in that columns of data ( i.e . , large arrays for this

application) are entered and various operations performed

on these columns. Difficulties arise in cataloging

statistics by Miller indices using OMNITAB because a separate array is required for each Miller index, i.e .

,

210 belongs to the h=2 layer, the k=l layer, the 1=0 layer

and the +++octant.

A program such as STATUS written specifically for

crystallographic applications removes these difficulties.

It also can make use of crystallographic data files produced by other crystallographic programs such as the X-RAY

3
SYSTEM (Stewart et al., 1972) . It minimizes the amount of

manipulation of the data required and facilitates statistical analysis during such stages of structure refinement

as in comparing two data sets where the reflection lists

must be searched for matches. It calculates relevant crystallographic quantities and breaks the data into appropriate

subsets with minimum user intervention.

4
I.C. Application of Statistical Analysis to Crystallography

• t. Abrahams and Keve (1971) point out that normal probability plots overcome many of the shortcomings of the

conventional R factors because they make use of the individual

residuals, (Fo-Fc) /a (F0 ) and are useful in comparing sets of

observed data (Fo) from different crystals of the same

material. Structural parameters obtained from different

samples of the same material can also be compared. We have

extended the procedures of Abrahams and Keve to cover (i)

data collected on the same sample by different experimental

procedures to test experimental methods, (ii) data collected

by the same experimental method oh different samples of the

same material to test sample characteristics, and (iii)

sensitivity of models to experimental methods by examining

parameters obtained using data from case (i)

The statistics or residuals can be examined without

assumptions about their distribution. The residuals,

(Fo-Fc) /a (F ) , contain all available information about the

manner in which the fitted model fails to properly explain

the observed variation in Fo. Generally, residuals are

plotted against time sequence, the fitted variable Fc, and

independent variables such as the Miller indices and the

magnitude of the reciprocal lattice vector. Examination of

trends among the residuals provides a detailed basis for

judging how well the calculated quantities derived from the

assumed model fit the observed data. Plots of various subsets of residuals focus on specific parts of the assumed

model.

Although the indices are the only independent variables

free from error (by virtue of their being integers) here we

take the magnitude of the reciprocal lattice vector as being

free from systematic error because random and systematic

errors in the cell parameters are usually insignificant

when compared with errors in Fo. Crystallographic leastsquares refinements generally make the same assumption.

Three types of statistics are of use in crystallography.

They are

x (i) 6mt = [Fo (A) i -kFo(B)t] / [<? 2 (A) t +k2 a

2 (B)i] 2 which involves two sets of observed quantities, Fo (A)

and Fo (B) i

,

together with their associated standard

deviations, a (A) i and a (B) i . Fo (A) i and Fo(B)i are

related by the scale factor k, chosen so that £6mi

is a minimum. This statistic is naturally 1 independent of any model of the crystal structure

and serves to focus attention on experimental

procedures and sample characteristics such as

homogeneity.

5
(ii) 6Rt = [

| Fo | i - | Fc | i ] / 9(Fo)i

where the observed and calculated quantities, Fo and

Fc, are related by a scale factor which is usually

determined by the crystallographic least-squares

procedure by minimizing

E6Ri 2 = E [ (| Fo | i - | Fc | i ) / a(Fo)i] 2.

1 i

crfFo^i is the standard deviation of Fo. The 6R t statistic,

which is actually a residual, measures the discrepancy

between observed structural amplitudes and those calculated from the model parameters.

( iii) 6pi = |p(A)i - p (B) i | / [a p 2 (A)i + © p 2 (B) 1 ] 2

where the model parameters p (A) if p (B) t , and associated standa

deviations <? p (A) i and cr p (B) i are derived from leastsquares refinements using two data sets obtained from

(a) different crystals of the same material or slightly

different materials (one sample being more impure than

the other, for example) or (b) the same crystal but by

different experimental procedures. This statistic

measures the differences between models and hence

examines their sensitivity to differences in samples

and experimental procedures.

The information and trends contained in sets of these

statistics is revealed best by plots. Principal types of

plots and their diagnostic capabilities will be discussed

for each of the three statistics in turn.

1. Use of the Sm* Statistic

The statistic 6mi can be used in three major types of

plots, as discussed below. In addition, a subset 6.^^ of

the 6mi , chosen on the basis of reflection

hOO, hkO, etc., or layer, or octant can be

type of plot. This enables both isotropic

effects to be examined.

class, i .

e

.

,

used in each

and anisotropic

6
(a) Order statistics can be used to relate the 6m t (ordered in increasing magnitude) to the quantiles

Xi expected for a normal distribution. Meaningful

conclusions that can be drawn from plots of this

type have been discussed by Abrahams and Keve (1971)

and a summary based on their discussion follows. The

distribution of 6mi , calculated from two Fo values

and their associated er values, is normal (Gaussian)

if the 6m i values contain only random error. In that

case a linear 6mi vs X< plot arises and the deviation

of its slope from unity indicates the scale factor to

be applied to the pooled a value used in obtaining the

6rrii . A markedly non-linear plot or a plot with nonzero intercept indicates that the two data sets differ

systematically or that one or both sets of ‘(Fo)

systematically mis-estimated . Thus, the two data

sets cannot be related by a single scale factor and

cannot come from the same normal population. A

search for the source of discrepancy must then be made.

If the distribution of 6mi is identical to the

normal distribution the slope will be 1.0. A slope

greater than 1.0 indicates that the average value of

O’ (Fo) is too small. This situation may arise in two

ways: (1) The estimate of the random error as

derived from counting statistics in Fo is too small,

for example, if significant coincidence losses in

counting occur, or (2) a more likely effect is that

the Fo values contain systematic differences which

are comparable or larger than the random error.

Hence, Fo(A)-Fo(B) is greater than AFC based on random

error only. Situations (1) and (2) may be distinguished by such plots as those involving 6pi if

models have been refined for the two sets, or by plots

of 6mi against d* and log (Fo~/sin2 6 ) . Trends in the

latter plots indicate systematic differences in the

two data sets.

A slope less than 1.0 when the 6mi are plotted

against the normal distribution indicates that the

a (Fq ) are on average too large.

An appropriate way to adjust the scale of cr(Fo) is

to apply the observed slope in inverse proportion to

the average sigma for each data set. Systematic

differences may exist in the two data sets in addition

to the random error. Although 6nu plots used alone

cannot give an idea of the relative magnitudes of these

systematic differences (unless c(Fo) is known to be

correct) , a combination of 6mi and 6pi plots can.

Diagnosis using this combination is amply discussed

in the Appendix of the paper by Abrahams and Keve.

7
(b) The 6mt statistics can be plotted against d* for

all the data or for various reflection classes. Such

plots should show no trends if there are no systematic differences between the two data sets. A hornshaped plot symmetrical about the abscissa would

indicate, mis-estimation of cr(Fo) varying in a systematic way with d* and sin ©. Plots involving

different reflection classes or layers can be used to

check on the anisotropy of error in an effect such as

correction for absorption, which is a function of

d* , i . e , , sin I .

(c) Plots of 6m j. against log (F§/sin2 1 ) for all the

data and various reflection classes permit some

estimate of the importance of extinction effects.

The variable Fo/sin2© was chosen because it is

proportional to crystallographic Q, ( i . e .

,

X

3

F

2/V2 sin2 © ) , scaled by V /X°

.

This scaling is of no consequence

here because it enters into the scale factor when the

two data sets are scaled together. Both primary and

secondary extinction effects are proportional to Q

and so these plots provide tests for differences in

extinction between the two data sets. Differences in

slope for reflection classes would be indicative of

differences in anisotropic extinction. Horn-shaped

plots are indicative of mis-estimation of a(Fo) in a systematic manner with Fo . One might suspect simultaneous diffraction if large Situ values tended to be

associated with small values of Fo

.

8
Systematic differences in the observed data sets result

from such effects as:

(i) errors in the measurement of Fot (set 1) and Foi (set 2),

(ii) physical differences in the samples used, e ,q . , differences in the degree of primary and secondary extinction

or in absorption corrections because differences in crystal

shapes were not fully taken into account, and

(iii) chemical differences between the samples used for

collecting data. These chemical differences may also lead

to physical differences, such as in extinction or diffuse

scattering, depending on the degree of disorder.

2 . Use of the 6Ri Statistic

The 6Ri statistic, which is the residual obtained in

least-squares analysis, can be plotted in the same manner

as the 6rru statistic.

(a) The 6Ri can be formed into order statistics and

plotted against the expected normal distribution. A

linear plot with zero intercept (within limits corresponding to the error in the least-squares scale factor)

indicates that the Fo contain no systematic error or

trends that cannot be accounted for in some manner by

the refinement model. Conversely, a non-linear plot

indicates the presence of effects not taken account of

by the model. A linear plot with a slope less than

1.0 would indicate overestimation of a(Fo) , or if the

a (Fo) values are known to be correct, which is rarely

the case, that the model has been able to partially

account for some systematic effect in the Fo set. A

slope greater than one would indicate an inadequate

model or underestimated ct(Fo) values. These various

situations may be clarified by further plots as

described below.

9
(b) The residuals 6 R 1 may be plotted against d* for

all the data and for various reflection classes. If

the model is adequate, a linear plot of uniform

scatter and with a slope approximately zero should be

obtained. A curved or sloped plot indicates an effect

in the observed data not adequately accounted for by

the model. In addition, if the scatter is not uniform, i .

e

. , large | 6 Ri | tend to be at small d* while

small | 6 Ri

|

tend to be at large d* , giving a hornshaped plot about the abscissa, a mis-estimation of

a (Fo) in a systematic manner with d* is indicated.

Differences between plots involving reflection

classes such as layers may suggest experimental or

sample orientation effects. Plots showing uniform

scatter with a non-zero slope would indicate a model

with inadequate occupancy factors or thermal

parameters

.

A linear plot indicates that the residuals belong

to a normal distribution and justifies the use of

variance tests based on the F-distribution, such as

Hamilton's R-factor ratio test (Hamilton, 1964).

(c) The 6Ri statistics may be plotted against

log (Fc e /sin2 fi) , an appropriate function of the

" fitted variable" . As indicated previously in

section 1 where this type plot was discussed in

connection with 6m j, these plots can be used to test

for extinction effects, which are proportional to

Fc ci /sin2 9. Plots involving different reflection

classes and octants permit tests for anisotropy when

using unmerged data. Again, plots should show uniform scatter and no trends. A regular horn-shaped

linear plot symmetrically placed about the abscissa

would indicate that the model adequately fits the

data, but that the values of a(Fo) have been misestimated in some systematic way. If the largest

positive 6R t values tend to be associated with small

values of Fc, one might check for simultaneous

diffraction, which would serve to increase Fo

considerably for small values of Fo.

10
3. Use of the 6pi Statistic

The third type of statistic, 6pi , is plotted against

the corresponding half-normal distribution quantiles because the order of the parameters p (A) t and p(B)i is not

significant and the sign of 6pi is meaningless. No scaling

is involved because the parameters are on the same scale.

A linear plot with a slope of 1.0 and intercept approximately 0 should result. A non-linear plot indicates the

presence in data sets of effects which cause the parameter

values derived from the two sets to systematically differ

from one another. A linear plot with a slope less than 1.0

indicates the presence in the data sets of some systematic

errors or effects comparable in size to the random errors

that the refinement model cannot account for, a situation

which results in large values for a(p). { . Such effects may

also be indicated by trends in the residuals, e .

g

. , in

various 6 Ri plots.

A slope greater than 1.0 indicates that the average

a (p) i is too small. This may happen if the set of Fo contains

systematic trends which the model parameters have been able

to absorb, and should be confirmed by the absence of noticeable trends in 6 Ri plots. Any systematic bias should

ideally be removed from the model parameters. Alternatively,

a(p)j_ for each parameter set may be multiplied by the slooe

so that the resulting a (p) t values will provide an estimate

of the variance due to both random errors and the meansquare errors associated with the systematic bias.

Abrahams and Keve provide in the Appendix and figure 6 of

their paper a detailed explanation of how a combination of

6mi and 6pi plots may be used to estimate the relative

magnitude of random and systematic errors. 6Ri plots help

to formulate corrections and improve the model.

I.D. General Description of Program STATUS

1. Overall Procedure

The uses of the STATUS program have been described in

Section I.C. Because statistical analyses have

been applied specifically to crystallographic

quantities in STATUS, we have taken care to write the program in a manner which minimizes data manipulation and

sorting while allowing maximum flexibility. This is

accomplished by keeping the four arrays EM, X, LOGFC and

LogHKL in core at all times. Array DM holds the statistics

under analysis. Array X contains the normal quantiles,

d* values, or log (Fc“/sin2 0) values. Cataloging of the

I-th statistic into its reflection classes and octant is

accomplished by non-zero bits in the associated computer

word LOGHKL (I) . For a 36 bit word, as available on the

11
UNIVAC 1108 , this allows up to 33 classes (3 bits are required to specify the octant) while requiring one word per

statistic. The LOGFC(I) word holds the address on mass

storage of the reflection information associated with the

I-th statistic. Other routines such as DMDQ, DMCALC, and

PRYNT use this mass storage address to locate the reflection

information (stored on disk or drum)

.

Analysis of residuals from least-squares refinement

is straightforward. Reflection information consisting of

the Miller indices, Fo, Fc, sigma and, as an option, the

mean path length or the extinction coefficient are read in.

The statistic is calculated and stored. The reflection

information is written on word-addressable mass storage and

its address is stored in the LOGFC array. Fc must be available if plots involving the fitted variable, Fc 2 /sin26 are

to be made. Note that most files produced by crystallographic programs [e.g., files used as input for a Fourier

calculation and the X-RAY SYSTEM (Stewart, et al. 1972;

hereafter referred to as XRAY) binary file] contain the

Miller indices, Fo and Fc, rather than the statistics

(Fo-Fc) /ef (Fo) even though least-squares programs calculate

this quantity or an equivalent residual.

Treatment of two sets of Fo data is a little more

complex in that the two values of Fo and their associated

s (Fo) must be matched up unless the input files are the

same length and in the same order. Also the scale factor

relating the two data sets must be found before the

statistics can be calculated. Subroutine SCALE matches

up the two sets of Fo and optionally writes the Miller

indices and the two Fo and or (Fo) values on unit NR. SCALE

also finds the value of the scale factor. Optionally,

subroutine SPEC can be called to refine the value of the

scale factor to ensure that the sum of the squaredstatistics is a minimum. Once generated on unit NR by

STATUS, the file of matched values can be read by subroutine REREAD allowing one to make subsequent runs using

different values of the scale factor without requiring

the initial time-consuming sort. This feature would not be

available if only the indices and statistic were saved on

an output file for subsequent runs.

The order of data treatment was developed so that at

most three sorts are required. In order to produce normal

probability plots the incoming statistics must be sorted

in order of increasing magnitude to form order statistics.

The corresponding X, LOGFC, and LOGHKL arrays are also

rearranged to preserve relative indexing. Generation of

d* and Fc 2 /sin2S values is accomplished by subroutine

12
DMDQ using the addresses stored in array LOGFC to locate on

mass storage the indices and Fc values that correspond to a given statistic. Subroutine CATLOG is also called at that

time to classify the statistic because the indices are

currently available. The d* values stored in the X

array are then sorted in order of increasing magnitude and

the corresponding statistics in DM and the catalog information in LOGFC and LOGHKL are similarly rearranged. Various

subsets of the statistics, such as those belonging to the

hkO layer, can then be selected and plotted. In this manner the variation of the residuals with g*, for example,

may be examined as a function of reflection class, layer,

or octant with only one sorting of the data.

STATUS can also be used for statistical analysis of

structural parameters. They are read in by subroutine

INFING if the Fourier file written by least-squares program RFINE4 is available or in other cases by subroutine

USER which must be written by the user. Subroutine MTCHEK

finds all atoms with two sets of parameters and STATUS

calls subroutine DPCALC to calculate the statistic. Onlv

a half-normal probability plot is required in this case

since the order of equivalent parameters is irrelevant.

Subroutine UNIMED obtains the medians for the I-th order

statistics and these are used together with the percentage

points of the normal distribution to obtain half-normal

quantiles. This avoids the errors associated with small

samples (Hamilton and Abrahams, 1972). All structural parameters, statistics and half-normal quantiles are printed out.

2 . Crystallographic Data Required as Input

In keeping with the design objective, STATUS makes use

of crystallographic data files generated by most programs

used for crystallographic calculations. A comparison of

two sets of observed data requires the Miller indices; h,

k, i, Fo (set 1) , Fo(set 2) , Q (Fo) (set 1) , and a (Fo) (set 2) , quantiles available on most crystallographic data files.

Analysis of least-squares residuals requires Fc instead of

the second set of Fo. Fc is available on files used as

input to Fourier analysis.

Input of two sets of observed data is accomplished by

subroutine FREAD which reads a BCD file. Optionally, STATUS

will produce a binary file of matched Fo and er values which

can be used as input for subsequent runs. Quantities for

analysis of least-squares residuals may be read from the

XRAY binary file by subroutine INFOFC or from the Fourier

file written by RFINE4 (Finger and Prince, 1975) by subroutine FINGFO. The specifications of these input files

will be discussed in the section on Input Files.

13
Input of atomic parameters is treated somewhat differently because of the wide variety of sources possible.

Presently, they may be read from the Fourier file written by

RFINE4 or from cards or another file by subroutine USER

which must be written by the user. Subroutine USER allows

for structural parameter input from many different sources

including published works.

Only a small amount of information needs to be input

on cards. This includes a title to identify the job, the

unit cell parameters, the wavelength at which the data were

collected, and the form of the input data (Fo or Fo 2

) . If

the optional spherical absorption correction is requested

the values of pr for the two data sets are required.

Formats for card input will be discussed in the section on

Card Input.

3. Format of the Output Plots

All plots are output on the line printer so that no

special plotting devices are required. Thus the user has

a visual representation of the statistical analysis as

soon as it is completed. A detailed description including

examples of actual output will be discussed in the section

on Output Plots.

The title of each output plot indicates the type of

plot, amount of data in the set or subset, number of points

in the plot and the octant or the reflection class the

statistics belong to if a subset is being plotted. The

actual quantities plotted are listed under each plot so

information corresponding to the residuals plotted can be

referred to.

In producing normal probability plots, the program

ignores those statistics which lie on the extremes of the

plot ( i . e

.

, with | Xi | > 2.0, where Xi is the normal probability plot quantile) so that they will not influence the

calculated slope and intercept out of proportion to their

importance. The slope, intercept and goodness of fit are

calculated from the functions (i) y = mx + c, (ii) y

= mx

2 + nx + c, and (iii) y = mx

3 + nx

2 + px + c . If the

plot is non-linear (somewhat kinked) the fit of the cubiccontaining form is appreciably better as revealed by the

goodness of fit. However, in this case the slope, given

as the coefficient of x in the equation, is that of the

central portion, and is probably smaller than it will be

when the scaling error has been corrected. Thus, the nonlinear functions serve mainly as diagnostic services to

aid in judgment of linearity rather than to obtain estimates

14
of slopes. The 6mi quantities, the least-squares line and the

cubic-containing function are plotted on the same plot

using the symbols • + and *, respectively. Multiple points

are indicated by a digit signifying the number of superpositions. If there are more than 50 points in the data

set, every (N/50 + 1) th point is plotted starting with the

point halfway through the first range. When residuals are

being analyzed the variation of scatter of the residuals

with the independent variable is important. If (N/50 + 1)

is six or more ( i .

e

. , each plotted point represents an

interval of six or more data points) the average positive

and negative deviation from the values calculated by the

linear function for the interval is plotted. These values

are indicated by the symbols U and L, respectively, and aid

in judging the variation of the scatter with respect to

the independent variable.

4. Program Configuration and Requirements

The program STATUS consists of a control routine

(also called STATUS) and 35 subroutines which perform the

actual functions. This modular form allows use of standard

documented subroutines when possible and makes it more

convenient to modify for specific situations if necessary.

It also allows for overlays and mapping for conservation of

core if needed. The present configuration requires about

10.000 (decimal) words (9,000 if mapped) for the code and

40.000 for data. This allotment enables treatment of 7200

statistics or 1000 structural parameters (100 atoms in the

asymmetric unit) which is more than adequate for most

crystallographic situations where the use of STATUS is

warranted. The amount of storage reserved for data can

easily be changed as described in the section on the FORTRAN

code

.

STATUS requires two special functions in addition to

the standard library functions. They are NTRAN, a UNIVAC

routine handling transfers of information between central

memory and mass storage and FLD for manipulation of specified bits in a computer word. The section on Features

Specific to UNIVAC 1108 FORTRAN V describes these functions

and where they are called in more detail. Users can substitute appropriate routines for their computing systems.

There are a few FORTRAN statements which are not American

Standard FORTRAN and they are discussed in the section on

Features Specific to UNIVAC 1108 FORTRAN V.

15
A maximum of six logical units including the card

reader (or remote terminal) and line printer is required

for I/O. One unit (NDRUM) refers to word-addressable random access mass storage and is used by the program for

temporary storage. Two units (NTAPEA, NTAPEB) can refer

to any convenient devices holding the input data files.

An optional unit (NR) may also refer to any convenient mass

storage device (tape, drum or disk) but not a card unit

since information in that file is buffered to save I/O time.

The next section presents the subroutine hierarchy

and call sequence to aid those who wish to map the program.

The present mapping scheme is as follows:

Main segment (always in core)

Routines: STATUS, DRUMRD, DRUMRT, POSITN, MAXHKL, HKLGEN

First sub-segment-starts at the end of the main segment

Routines: PRIME, CELL, FREAD, SCALE, ABSORB

Second sub-segment-overlays the first.

Routines: REREAD, SPEC, DMCALC , NFACT

Third sub- segment-over lays the first and second

Routines: INFOFC

Fourth sub-segment-overlays the previous segments

Routines: SORT, PROB, PINV, UNIMED, PLOTEM, LSFIT, FITCHK,

PLOTS, SIMLEQ, FINGFO

Fifth sub-segment- star ts at end of the fourth

Routines: DMDQ, ANISO, PRYNT, CATLOG

Sixth sub-segment'Over lays the fifth.

Routines: DPCALC, INFING, MTCHEK, PRYNTT, USER

16
5. Subroutine Hierarchy and Call Sequence

As an aid to program modification by users who

may need to replace machine specific routines or to

reduce the size of the program the order in which

subroutines are called is given below in diagramatic

form.

i STATUS

t [prime

| CELL JPOSITN

NTRAN

I ABSORB | |NTRAlT~j

17
FINGFO INFOFC

Imaxhkl DRUMRT

|

NTRAN]

© ©

18
® -

___ [user;

1 ^—' ®

'© [

1—

®

If ingfo

j

i infingI

Idrumrt

1 NTRAN

1

19
20
6. Description of Subroutines

ABSORB - Calculates the transmission factors for a spherical crystal given the Miller indices and reciprocal cell constants.

ANiSO - Calls PLOTEM to set up plots for various reflection

classes according to variable NCL . Sets flag II so that

LSFIT, FITCHK and PLOTEM take only those 6mi values and

reflections whose catalog word LOGHKL has the II-th bit

non-zero

.

CATLOG - Catalogs a reflection according to class (hOO,

OkO, etc.), layer, or octant, by setting equal to 1 the

corresponding bit in word in the LOGHKL array which is

associated with that reflection. Requires the indices and

reflection number.

CELL - Reads in real cell, calculates reciprocal cell and

stores results in common block G.

CHANGE - non-existent at present but will contain cards

which will allow user to keep track of any changes made

to present program.

CHANGE - Contains comment cards which allow user to keep

track of changes made to program. Also contains general

update information and version date.

DMCALC - Reads off mass storage using the addresses stored

in LOGFC the Fo and a values for two corresponding reflections. Calculates 6mi (unit deviate) which is stored in

array DM, (6m), e.s.d. of 6mi /and a weighted R factor over

the two data sets. Prints out these quantities to give a characterization of the 6m distribution.

DMDQ - Reads from mass storage using the addresses stored

in array LOGFC the indices and other information pertaining

to the reflection associated with a given 6mi value.

Depending on the variables ND and NQ in its call argument

list, DMDQ calculates d* or crystallographic Q = (F0^sin2@)

and stores the result in the X array so that 6m can be

plotted against whatever is in X. DMDQ also calls subroutine CATLOG to catalog the reflections by classes if

desired.

DPCALC - Calculates 6p t for two corresponding atomic

occupancy, positional or thermal parameters.

DRUMRD - Reads a record containing the indices, f and g for

one reflection from word-accessible mass storage (drum or

disk)

.

DRUMRT - Writes one reflection record on word-accessible

mass storage.

21
FINGFO - Specialized subroutine for reading indices, Fo,

o

1 (Fo) and Fc for each reflection from end of a Fourier file

written by program RFINE4 (Finger and Prince, 1975). Calculates the corresponding 6Ri value and stores it in array DM.

Can also be used to obtain two sets of Fo and er (Fo) from

different Fourier files containing different data sets obtained from the same compound.

FITCHK - Calculates points from the least-squares line and

cubic-containing functions obtained by LSFIT and stores

them, together with the original 6M^ values , for combined

plotting against abscissa values stored in array X. Also

calculates the goodness-of-fit for. the various least-squares

functions

.

FREAD - Reads reflection records containing indices, Fo and

<7 (Fo) from input units NTAPEA and NTAPEB. Stores the indices, Fo and er on disk or drum and keeps the value of the

slowest varying index in array ID for use in the matching

procedure in subroutine SCALE. In this way SCALE is able

to read in only those reflections with the same value of

that index, conserve memory space and match efficiently

using small batches. The position of each reflection in

the mass storage file is known from the position of its

slowest varying index in the ID array, since both the mass

storage file and the ID array were written incrementally as

the reflections were read in.

HKLGEN - Unpacks the Miller indices from word JKL according

to the order specified by JJ, KK and LL.

INFING - Specialized routine to read parameters from the

end of the Fourier file written by program RFINE4

.

INFOFC - Reads h,k,4, Fo, Fc, and a (Fo) from XRAY binary file.

Calculates 6m, writes h,k, l, Fo, and Fc on mass storage,

stores address in array LOGFC.

LSFIT - Calculates least-squares line and quadratic- and

cubic-containing functions relating variables X and Y.

Y is the complete set of statistics, SMj_, or a subset thereof.

X may be the normal quantiles, d* values, Q values or subsets

if reflection classes are considered.

MAXHKL - Finds maximum and minimum values of Miller indices

and stores values in common block D for further use.

22
MTCHEK - Matches up parameters for two parameter sets,

checks for unequal numbers of atoms ( e .

g

. , one set may

include hydrogens, etc.). Atom matching is based only on

the alphanumerical name given for each atom, so that corresponding atoms must have the same name, e .

g

. Ca 1, with an

identical number of embedded blanks.

NFACT - Computes N factoral which is used in routine SPEC.

PINV - Inverse probability function (percentage-point) used

to obtain the normal quantile corresponding to i-th value

of P (X) for the j-ordered statistics.

PLOTEM - Sets up plot titles, calls LSFIT for least-squares

fits, writes out coefficients of fitted functions (lines

and curves) , calls FITCHK to calculate goodness-of-fits,

prints out titles for plots, and calls PLOTS to do the

actual plotting.

PLOTS - Routine to plot up to five curves stored in X and

Y arrays, using the line printer. Keeps track of multiple

points and finds its own limits for the axes. Does not

call a new page and does not label the axes.

POSITN - Dummy routine at present. User may add statements

to position logical unit NTAPE so that NFILE is ready for

I/O.

PRIME - Initializes the logical unit variables: IN - card

reader, NOUT - printer, NDRUM - mass storage. Separate

routine allows user flexibility and easy modification.

PROB - Calculates i-th value of P (X) for the j-ordered

DM values so that PINV can be used to obtain the corresponding normal quantile.

PRYNT - Prints out all SMj values and corresponding reflection information or only those greater than 2.0 if desired.

Reads values from mass storage and requires LOGFC word

which tells where reflection record corresponding to SMj

value is located. Calls HKLGEN to obtain the indices.

PRYNTT - Prints out 6p t values for parameters and corresponding atom names, etc. Requires LOG word which tells

what atom and type of parameter (occupancy, positional,

etc.) go with a given 6pi value.

REREAD - Reads the quantities Fo and a (Fo) for data sets

one and two from a file written on unit NR by routine SCALE.

Thus, the scale factor may be adjusted without re-sorting

the original lists of data to find the two Fo and Q values

having common indices. New 6mi values are calculated.

23
SCALE - Reads the indices, Fo and ex (Fo) for the two data

sets from mass storage and finds reflection information

having common indices. Calls MAXHKL to get range of indices,

applies absorption corrections if needed. Also buffers up

the matched reflection list and outputs this list on unit NR

if desired. Calculates the scale factor required to put the

two data sets on a common basis by forming (set 1) /

£f01 ( set 2 ) .

1

iSIMLEQ - solves simultaneous equations in matrix form

AX=B

.

SORT - Sorts the 6M^ values in ascending order so that thev

form order statistics and keeps the LOGFC words in the same

order as the sorted 6M^ values. The LOGFC word tells where

reflection information corresponding to the i-th 6Mi value

is located on mass storage and enables reflection information associated with the i-th 6M^ value to be accessed easily

by routine PRYNT. For atomic parameters, the 6P^ values are

stored in array DM and the ordinal numbers of corresponding

atoms are stored in LOG. The SORT subroutine orders the

6Pj[ values and keeps the LOG values in correspondence with

the 6Pi values. Parameter information is printed via

PRYNTT

.

SPEC - Obtains the scale factor needed to make SSrrg

3 a minimum. This may be different from that obtained in a routine SCALE.

UNIMED - Computes an approximation to the median of the

i-th order statistic. Used here to generate quantities

needed to obtain quantiles for the half-normal plots.

USER - Routine to allow user to input atomic parameters from

various devices and formats. Must be user written.

24
I.E. References

Abrahams, S. C. (1964) Acta Cryst. 17, 1327.

Abrahams, S. C. and Keve, E. T. (1971) Acta Cryst. A27, 157.

Hamilton, W. C. (1964) Statistics in Physical Science,

Estimation, Hypotheses Testing, and Least Squares . New York: Ronald Press.

Hamilton, W. C. and Abrahams, S. C. (1972) . Acta Cryst.

A28, 215.

Hogben, D. , Peavy, S. T. and Varner, R. N. (1971) OMNI TAB II

User's Reference Manual , NBS Technical Note 552

.

MacDonald, J. R. (1972) Science 176.

Finger, L. W. and Prince, E. (1975) A System of FORTRAN IV

Computer Programs for Crystal Structure Computations , NBS Technical Note 854.

Stewart, J. M., Kruger, G. J. , Ammon, H. L., Dickinson, C.,

and Hall, S. R. (1972) The X-RAY System of Crystallographic

Programs . Technical Report TR-192, University of Maryland,

College Park, Maryland.
II. Program Manual

II. A. Card Input

Although this section is called Card Input, STATUS

could be run from a remote terminal if the printer output

can be directed elsewhere. Such details depend on the

particular operating system and are not discussed here. The

following description refers to cards (or card images) that

are read by the program STATUS itself. Card input is referenced by logical unit IN and is presently set to 5 by

subroutine PRIME.

Card 1 Title Card - read by main routine STATUS

FORMAT (12A6)

Col 1-72 any title

Card 2 Cell Card - read by subroutine CELL

FORMAT (13X, 3F8. 3, 3F9. 5) (Compatible with XRAY FORMAT)

Col 14-21 A cell dimension

22-29 B cell dimension

30-37 C cell dimension

38-46 Alpha or Cos Alpha

47-55 Beta or Cos Beta

56-64 Gamma or Cos Gamma

Card 3 Control Parameters - read by main routine STATUS

FORMAT (A3, 1215, 13, 1612)

Col 1-3 Punch CTL - program checks third card read for

this label

5 NTAPEA - first source of input data-reflections

and/or atomic parameters.

10 JJA - sort order of indices of reflections

15 KKA

2 0 LLA

The sort order governs the order in which H,K and

L are packed together in the " JKL" word, ( i .

e

.

,

- JJA=2 , KKA=1, LLA=3 would pack

the indices in the order K, H, L . ) The reflections

are broken into batches with the first index

constant in the packed word (i.e.-in the example

above, one batch will have all K=0, another will

have all K=l, etc.)

The sort will run fastest

for a given number of reflections when there are

the most batches, so the sort order should give

1 to the index associated with the largest (real

space) dimension of the unit cell.

25 NTAPEB - second source of input data.

26
Col 30 NFTEST 0

1

no plotting of reflection data,

read input from RFINE4 Fourier tape

and compare Fo and Fc

2 read in two different reflection

sets, each from a different RFINE4

3

Fourier tape.

read reflections from two different

tapes in XRAY BCD input FORMAT

e.g.. Col 15 less than indicator

19 H

23 K

27 L

. 28-37 Fo

4

38-47 a (Fo)

read tape previously written by

subroutine SCALE. The tape is

designated NTAPEA (Col 5) , the

input file on this tape is designated by NFILEA (Col 60)

5 read Fo and Fc from XRAY binary

file assigned to unit NTAPEA.

35 NPTEST 0

1

no atomic or group parameters,

atomic parameters from RFINE4

Fourier tape.

2 parameters read in by special subroutine USER supplied by the user.

40 NFILEA No. of input file on NTAPEA, 1 for first

file, 2 for second file, etc.

45 NFILEB Same information for unit NTAPEB

Following section selects output options.

48 NPR I NT 0 no printing of statistics (DM) and

associated values.

1

2

print all DM and associated values,

print only those with absolute value

greater than 2 . This is about

5 percent of the data. About 55 DM

will be printed per page of output.

50 NO 1 for plots of statistics vs normal

0

quantiles

no such plots

52 NOC Governs reflection type for normal quantile

0

1

plot.

all data regardless of class or octant

standard classes H00, 0K0, 00L, HKO

,

HOL, OKL, HHH , HHL, HKK, HKH, and HKL

(non-zero indices)

27
2 H layers, up to 10

3 K layers, up to 10

4 L layers, up to 10

5 does NOC = 1 and NOC = 2

6 does NOC = 1 and NOC = 3

7 dees NOC = 1 and NOC = 4

54 NLAYO Number

NOC.GT.

of layers

1 . for normal plots

56 NOCTO Option specifying octants to plot.

0 all octants plotted together for

each class or layer. Use NOCTO = 0

for merged data.

1 HKL and -H,-K,-L, i .

e

. , +++ and

(check triclinic)

2 +++, and +-+ (useful for monoclinic Y axis unique)

.

3 +++, , H

b, and +H— (useful for

monoclinic Z axis unique)

.

4 H—l

f- ,

, H—l— , H

b, —I

b (useful

orthorhombic equivalents)

.

7 all octants plotted separately for

each class or layer.

Note the following concerning pages of output.

Number of plots generated - 1 if NOC = 0

plus 11 * (NOCTO + 1) if NOC = 1

plus (NLAYO) (NOCTO +1) if NOC = 2,3,4

plus (11+NLAYO) (NOCTO+1) if NOC = 5,6,7

with two pages output generated per plot.

58 ND 1 for plots of statistics vs d*

(2sin (9) /wavelength)

0 no such plots

60 NDC Governs reflection type for d* plots

values as above for NOC.

62 NLAYD Number of layers for d* plots

if NDC.GT. 1.

64 NOCTD As for NOCTO above.

Estimate number of plots by replacing variables in above

formulae with NDC, NLAYD , NOCTD.

66 NQ 1 for plots of statistics vs Q (Fo 2 /sin (2 © )

)

0 no such plots

28
68 NQC Values as above for NOC

70 NLAYQ Number of layers for Q plots if layers

being done.

72 NOCTQ As for NOCTQ above

Calculate number of plots as above using NQC,

NLAYQ, NOCTQ.

It is possible to produce a maximum of 588

plots or approximately 1200 pages. Therefore

do not abuse the flexibility of this program.

Program checks to see how many plots will be

generated and stops if more than 120 would be

produced in one run.

74 NR Output of matched reflections (pair of

Fo, etc.) on unit NR allows for SCALE

adjustments without resorting.

76 NFILER File to be used on unit NR.

78 NABS Greater than zero for spherical absorption corrections.

80 NINDM Greater than zero for adjustment of

scale factor to ensure that £6m2 is a minimum.

4 Read by main routine STATUS only if NQ or NABS ^ 0.

FORMAT (2F10.7, Il,4F5.2)

Col 1-10 Wavelength—needed for Q plots and absorption

correction

11-20 Scale factor to be applied to second set of

structure factors (needed only if matched

reflection being read by subroutine REREAD

from file written by SCALE)

21 0 for Fo data, 1 for Fo 2 data input.

Following 2 quantities needed only if NABS.GT.l

22-25 PR for first data set.

26-30 PR for second data set.

31-35 Delta pR (max PR - min pR) for first data set

(may be left blank)

36-40 Delta pR for second data set (may be left blank)

29
II. B. Data Files

1. Input Files

A new run will require at least one and more likely two

input files which are referenced by logical units NTAPEA and

NTAPEB . These input files may have different formats depending on their source. Control parameter NFTEST allows

for selection of input type (setting NFTEST = 0 means no

input is expected) . Different types will be discussed in

turn

.

a. Comparison of two sets of observed structure

amplitudes

.

If NFTEST = 2 subroutine FINGFO will read the Fourier

file written by program RFINE4. FINGFO was specifically written for use with RFINE4.

If NFTEST = 3, subroutine FREAD will be called to read

a file having values of h, k, i, Fo and a for each

reflection. The present form of the read statements

allows for a BCD file of card images compatible with

the format (13X, 12, 314, F10.2, F10.4), which is

compatible with the XRAY reflection data format.

However, a different format or a binary file may be

read after the read statement in FREAD has been changed.

A description of the procedure and flow chart for FREAD

has been included in Section III.E to aid in modification .

If NFTEST = 4, subroutine REREAD is called to read a file of two sets of corresponding Fo 1 s and their a values. This file must have been previously written by

subroutine SCALE and is a binary data file consisting

of blocks of 490 words each (this is a logical record

for a FORTRAN write operation) . Each block consists of

70 reflection records. A reflection record consists of

the seven quantities h, k, i, Fo(set 1), Fo(set 2),

a (set 1) and a (set 2) . Logical unit NR refers to this

file. No modification should be necessary unless the

user wishes to optimize storage (such optimization requires a detailed knowledge of the local computing

system)

.

30
b. Comparison of a set of observed and calculated

structure amplitudes . If NFTEST = 1, subroutine FINGFO

will read Fo and Fc values from the Fourier file

written by program RFINE4.

If NFTEST = 5, subroutine INFOFC will read in indices,

Fo, Fc, and g from the XRAY system binary data file.

A complete description of this file is given in the

report describing the XRAY system. Subroutine INFOFC

can easily be modified by users wishing to read from

a different type of file.

A description of the procedure and flow chart for

INFOFC is given in Section III.F. This description taken together with the listing should aid in

modification

2

Output Files

All output except the optional file of matched reflections already discussed is on the line printer. Logical

unit NOUT references the printer and is set to 6 in subroutine PRIME.

3.

Temporary File

During execution, STATUS refers to a random-access

word addressable file which is written by subroutine DRUMRT

and read by DRUMRD. Presently, this is a binary file and

is written and read in blocks of four words: the packed

index word, JKL; the reflection flag, NT; Fo, and g. It is

addressed relative to the beginning of the file. If Fo and

Fc are being considered, the file structure has JKL, NT, Fo,

g followed immediately by JKL, NT, Fc, absorption path

length or extinction coefficient. If two sets of observed

data are being treated the corresponding values may not

follow one another but both their addresses on the file are

stored in one LOGFC word.

4.

File Assignments

A description of actual file assignments to the physical

devices is not given here since this depends on the particular

computing system in use. STATUS allows the user to specify

the actual values of the logical units either as control

parameters on input card 3 or by simple modification of subroutine PRIME. This should provide sufficient flexibility

regardless of the particular form of file assignment.

31
II. C. Examples of Output Plots

Several plots from an actual production run have been

selected to provide examples of the types of plots possible

and their uses. The purpose of the production run was to

compare agreement among intensities equivalent by symmetry.

The unit cell of a tin phosphate appeared to be monoclinic,

but the intensities of some of the equivalent reflections

were not consistent with monoclinic symmetry even after an

absorption correction was applied (with some difficulty

because the crystal was an irregularly shaped fragment)

.

The statistics, 1 6mhke = (Fhke - hk£ + ° hk i

S

and ern^-j = (Fhkl - FHk£ ) / (o|kJj + 9^-5 )® were formed by

reading in hk£ and hk£ reflections as two separate files for the

first case and hkl and hk 9, reflections for the second case.

The results were analyzed in various ways in hopes of detecting systematic effects. A re-establishment of the scale

factor by minimizing Z6m^ ^ gave a value of 1.019. This is

satisfactorily close to the expected value of 1.000

because our scaling of ‘the data through remeasurement of

standard reflections has a precision of 1% on intensities.

1. Normal Probability Plot

The first line of the plot title (Fig. 1) indicates that

all statistics regardless of reflection class or octant were

plotted against the expected normal quantiles. The designation DM (VERT) reminds us that the ordinate axis refers to

the statistics. The line printer does not allow resolution

of every point if there are more than 100. Plotting 50 or

less eliminates this problem and allows for easier reading.

Thus, the second line states that every 106th point of the

total 5275 in the data set was plotted. The next line gives

the equation of the least-squares line and quadratic that

best relate the statistics to the normal quantiles. The

third line gives the cubic equation that relates the statistics to the normal quantiles. Finally, the fourth line gives

the goodness-of-fit, i

S = (£ (6 i -6 calc ) S / N-P) ?

, for the

various curves. These quantities are intended to aid in

judging the variation. The definitions of the plot symbols

are given below the plot. An accompanying table shows the

values corresponding to the 50 points that were plotted.

The cubic function curves downward more rapidly than

the line for negative statistics and upward more rapidly

for positive statistics. This is also revealed directly by

the lower value of the goodness-of-fit parameter for the

cubic function when compared with those for the quadratic

function line. Thus, the distribution of statistics has

"thicker" tails than the normal, ( i . e

.

, there are more extreme values than expected) and the distribution tends to

be symmetrical about its mean.

32
Figure 1 » a iI 't i I l l l l l + i i I i i oo<

2:

ccQ

CL I O I 2! I

c_5 k-

<T

O

CD | O I O I

Cl

Q CL

— mo K OO CO * LL - CC *

LL) I X >0 > h- U * >3* IZI-O- • o *- < nj

LL • LL l/' Q O K < z z>

h-

— o CO o J Q. I/) Z

Q_ LL < O a. in < <M O 3 * — o * k- c0 X <

1/3 CL Z k* # C c co CO <

1— < ao

O

k- LUO

O

—

LUa

x

co *

U.

O + LL

OO O' LU -i

X O' OH- • * O O uj eg m

_J CO st 1

i/>

m

02

m

—

— O'

-I z

LL

—

LU

O

cla

11 j ini

-1 r- o < cmin co z

4.

—

-j

O'm 1/3 CM -Jo 1 U- II - . ffl I co *-« 11 in | 2; m Q

1/3 :

« LU U Z ~ CO k nh a: « z _j o — CJ K O co — a lu co u.

cc LU

X -a

a ll k- 3 < O o 0 O 1/1 O u"

co l/; —< k- UU > CO l- Z

c. < CO

Q

UJ LU < O > -J LU C

LU _JO # CM# CM # • + & CM

eg * + cm

<\j # + cmCM tt + CM + • #

C\J # + • *

CM * + •

CM #

fr CM # 1 1 1 1 1 1 1 + 1 1 1 1 i 1 1 1 1 1 1 1 r 1 1 1 1 1 + 1 1 1 1 1 1 1 1 1 + 1 1 1 1 1 1 o+

CC

C\J O'm

do 1 o I o I 2o i * IIIIIIIII +I ?

xO

cr> ! 2
Figure 1

contd.

—I

23

I

li

CJ

r0 L

J

<

-• :d

-H r-t r-t MX r-M o o O O O o O o o O o o O O o o r-t ^-t r-t M-i o O O O O O O O O o O o O O O o H r-l r-H r-l •—

!

H

—J L 05 o O O o o o O o o o o o o o o o o o o c o o o o o o o O O o O O O o O o o o o O o o o o o O o o o

3 4 4 4 4 4 4 4- 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 1 1 1 1 4 4 4 4 4 4 4 4 4 + 4 4 4 4 4 4 4 4 4 4 4 4 4

• LA a K 4 a 4- A CO 4 4- tO A •0 A CM K 4^ (\J CM CM CM A A A A 4- 't A 4- (A ro m 4" A O 4* r-t O' o CO O' O A C- A A 0- O h- —

4

lA c 00 >0 4" A c\J *“•4 o >0 O' f\l A O' a K (M A *—

i

A r-t A i—l A •H A CX5 O' h- CM h» CM K cm r- csl CO ro O' 4- >—

1

K A r—i A o r—CM r. A A A

t -4 -* r-t «-4 O' CO CO K >0 45 A A <t 4 a a cm CM A H CM f- -< i—

t

CM CM ro CO 4- 4' A A A K h- 00 A A •—

A

r—i rH -* r-t

-J

1

t—uX

1 1 i 1 1 1 1 1 1 1 i 1 i 1 1 i 1 i 1 1 i 1 1 1 1

#

«\j S

Cj 2*—

a

O' O' A h- r-t O' K A CO A r-t re. A K A O O' 4- •Nl A C- A CO CO c- o 4" A CM CM 4* c- CM ro* 4- r—

1

CO CM O CO A A A CM A CO r- r\l O

11 A o o oo A cm CO o A O O' *—

t

A o A A r—i O O O' O —

t

CM CO 4- r\J o CO A 4- CM O CO A 4- A r\l *—

<

r—t r—

1

CM 4- A O O 4* A A O O' A A A Os) A O >0 a O CO A 4* cm -* O' CO r- A A a A CM r-4 O O r—

1

CM CM m 4- A A A CO O' o *-* CM CO 4- A A CO O' A A A CM 00 K

-It- •A

1 41 A

1 a

1

c\J

1 cm1 CM

1 1 1 1 i 1

i i l i i i i 1 1 i -1 9 1 l - - -1 - - CM CM CM CM re. CO 4

1

• LL a o

- o

_J 4- A A A 4- (M 0 A o- O <A -4 A o 4 o A A r~4 4- A CM A A ro A A CM r—

1

ro Is- A A 4- K A CM A O A A CM r*- CO

-j on At C- CM 4- O A O CM A r-t A r0 A c~ 4- A r-t o iO K h- A A A A A o A C- A A O r—

1

CM 4* A A -t A A A A A CM 4" A A A A

O A 4 CO A v'O tM O A h- A 4* CM —

t

O A CO r- A 4- A CM O O CM A j- A A A A O CM A 0 K A O CM 4* A A O A K <\J < A A cr» CM CM CM f\l CM r-t r-t r—

t

r-H

i 1 i i t i | i 1 rH r-t •—t ^t —^ r—

t

r-4 r—

t

CM CM CM CM CM A A A 4" > k l i l i i l i l i I i i l i i A ian

LU l/) CO H (T1 H

t— LU <T CM f—t CM < oc ^ (nj r^ OHlTiHOlTih

Or-tOOG-'MCCvt

<! <h co 'D iTi rn (v A cm A a A K A

r^- h o h ^

r-* O A CO I" vO in Q >0 nj O h o -i

cm cm <f- a oo w* o

•t ro (\J H O G

^ co o in co m o

r-i o o 0s k k cCM A 4 4* A A K

m <i A4-ACMrMCM.-tr-i—t.—

t

o > I I I I I I I I t I QC. C <r ^ to

oo •-•OK — OO 2T Z CO m

Q < O

a

ii II < C:

^ ^ I I I I I I I I I I

I I

nj fCi «t in o h A A O r-l <—I r\| ;\J njin<f-inxnr*-coa'0^

no cm cm c\j cm c\i cm cm A A

cmA

n^inin0'Hc\jGina'fCi0'inHa'O't>j- AAAAAOAAOCMA—•cmAOK'-^A

cc o o -i(v^if\ >o t:0 (J'Onj>t'0 (nMNin

cn st in >o f*- co ^

cn m r^i m in rn m

o h cm a vj* lo

4* 4* 4- 4- 4* 4- vOMD

•4* >1* *3" O O4” A
The non-linear character of the plot indicates the presence

of non-random effects. This inference is confirmed by plots

of the deviations against various independent variables such

as d*. The coefficients of the linear equation describe

the mean and standard deviation of a normal distribution that

is closest (in the least-squares sense) to the actual distribution . The slope indicates that the average AF/a is 2.2 rather

than 1. This may result from Af being greater than AR, the

deviation due to random errors, or from a being underestimated, or from a combination of these effects. Further

progress can be made if one has independent or prior knowledge

about the correctness of the assigned a. Our experience with

our normal data collection procedure suggests that the

assigned a based on counting statistics tends to

represent certain errors satisfactorily. These errors

we associate with short-term (time required to record

a data point) instrument and environmental fluctuations

and are taken as random. Since AF/a = AR/ct + AS/a = 2.2 and Ar/ct

1, AS/a ~1.2 (AS is the deviation resulting

from systematic effects in the two F values) . AS = (s(l) - S(2) ), where S(l) and S(2) are vectors since they

have both a magnitude and direction. Thus o< AS/a < | S (1) | /

n + S(2)/a, and |S(1) \ / a + |s(2) | / cr^ 1.2. The systematic effects, |

S

( 1 ) | and |s(2) | are at least 0.6a for the two

equivalent sets of reflections and may be much larger if

significant cancellation of systematic effects has occurred.

However, we note that the slope for the cubic near a normal

quantile of zero is 1.65 which shows that the actual

distribution envelopes the normal distribution and suggests

that systematic effects have not completely cancelled each

other out. A slope of less than one would have indicated

significant cancellation of systematic effects.

2. d* Plot

The second plot (Fig. 2) shows the variation of statistics from reflections belonging to the (lk£) layer with the

independent variable d*. The " +++OCTANT (S) " indicates that

all indices are positive which in this special case (see

section II. C., Introduction) means statistics Ahk£ formed

from an<3 F^ ?

. The plot heading has already been

explained. The plot shows that the deviations vary with d* . and in particular, that negative deviates occur below 0.50A-1

(Bragg angles <15 degrees) indicating that F^^ is larger

than the equivalent F^^ . At higher angles the trend of

35
the deviates is toward zero as it should be. The plot

symbols U and L represent the positive and negative

deviations of the statistics from the values given by the

line as averaged over the 14 statistics that each interval

represents. Thus U and L show the scatter of the statistics about their linear trend. We see that the scatter is

not uniform as it should be, but that the largest Af/ct

values occur for d* values below 0.50A-1 (the lower third of

the range) . This implies that the a values for the

observed F values in this region are underestimated and that

systematic effects are significant because a t is proportional to

Ft + S where S is the systematic effect and F T , c T are the

"true" values of F and n. There will be cases for F T + S<FT and thus a<a T especially when S>R.

36
«.+

nCT

ANT

(S»1KL

REFLECTIONS

PLOT

OF

DMIVERT)

AGAINST

I)

ST

AK

VFRY

1A

TH

POINT

IN

PLOTS

695

POINTS

IN

THIS

SET

OF

REFLECTIONS

50

POINTS

IN

THIS

PLOT

LEAST

SQUARES

LINE

IS

DM*

-1.585

l.331*X —

LEAST

SQUARES

QUADRATIC

IS

DM=

-3.553 +

5.893

*X+

-2.237

*X** 2

EAST

SQUARES

CUWIC

IS

DM=

-7.A69 *

21.305

*X+

-19.36A*

X**2 +

5.

663*

X** 3

OODNESS

OF

FIT

IS

3.031

FOR

LS

LINF,

3.011

FOR

LS

QUADRATIC

AND

2.990

FOR

CUBIC

FIT

IF

MURE

THAN A

PTS

IN

PLOT

Figure 2 + =+ I I o+

ccO'

ro

>c I I I l + i I l l I * * i i r i CM i i i i i-i ii-i + ii i-iriii

CM fSJ * + 3 # 3 #

fr 4- * + * 4* « 4- * 4-

CM CM * + * 4- # # 4 4- * + #

CM *# • * CM + CM

CM fO CM + -S- • 4 #

4- * 4 CM 4 * 4 *

CM# 4 * 4 I I I I I I I I I I I | | | I | I I I 4 I I I I I I I I I +

37 4 44 Oo

roOo

COoo-

~v

rCi CM I

12775+02+
=

DMC

>I

(ORDINATE*,

SYMBOL

IS

.

B

=

DM(I>-

L.

S.

LINE,

SYMBOL

IS

*

C

=

DM(

I>

-

L.

S.

CUBIC,

SYMBOL

IS

=

ABSCISSA

VALUES

(DSTAR,

LOG

F**2/SIN(2*THETA)

,

OR

QUANTILE

SYMBOLS

L,

U

ARE

FOR

AVERAGE

-

AND

+

DEVIATIONS

(A-B

VALUES

FOR

14

POINTS

IN

THE

INTERVAL.*

POINT

ABC

D

1

-7.405

-1.272

-3.457

.235+00

Figure 2

contd.

+

o O o o o o o o o o o o o o o o o o o O o o H r—4 o O o o o o o o o o o o o o o o o O o o o o o o o o

4 4 •f 4 4 4 + + 4 4 4 4 4 4 4 4 4 -f 4 4 4 4 4 4 4 4

sO CM ro in in 4" in o CO r- CM O' CO O' in o o CO 00 O' 4- >t O <—

(

O' CM in c—

1

4" 00 CO r- o in in o co f~( CO 4- 00 O' o O pH

CM cn 4* <r in in in vO vO h- 00 CO 00 O' O' O' O' O' r—A r~* i H

oooooooooooooooooooooo + + + + + + + + + + + + + + + + ++ + + + + (Mjfi'OcDiNj'Oco-irorfir-0(NJ>00'<NJ-+a'<M'+R-f'J

)^OvOvt--''4-o^f\j f-*vO«iooa'covONf(,n^tr--o^a'm'OOf\jOLr>vOvOco>oc7'NtLnof\j--0'f^r^ r-(NO'0'Oa'ONa''OcoHrvjLA^'^'J, Lnir\Nj'ff»n-<a'^'Oin^fNjHp^fnLf\'ONO'0'OLn't'-'»-*ina'M^Ln^o r^a'LnoocovOLnm(\j^ooooOf*---(^^^r-<^r^oooooooooooooooooooocNi(Nm>tr-

(VI ,-4 *-4 ^4 I I I I I I I

I I I I

i t r i i i i i i i

HvOvtff10C'CONCOOJ

OCTCvjCOo^OsO^iCD

r-tooo'o^cocor-r^vo ^ -1 H I

I I I I I I

I l l

<tvtmocomo'>t-o

^cor-inOr-.-Hr^co

vOini^mir\'t’tf^rn

I I I I I I I I I

‘^a'<NJtr'Os '4'(NJrO'tsO c\jsOsOOin(Nj^(^r->t

I I I I I \ I I I ^

COCOCDsffOvOO'tO'H Hnco(\nncoCM^NH OOOn^pjHMfMn

I

NsOf^io^ro^f^oco

ino'c0soc0'00''froo n co 't in in in vor^oo

>0 ao cm in ^ <n 0'Hffi'Oir\ffth-

•4* cm O' o m ~4 O o in c\j O' in c\j r- O co in n rg

'O N ro m h

O' H f\J IT* 00 >0

*-* >fr —* *-* •-» LH

O' fC| O 0s M O o ro cm it* O' m cm

O' CO fT) CO sO •—» o

(Fl M h- O sf H n

»0 n 'O in *h co co

h- (m O' m- r- a*4 *-*

f\j in in so eg —

h

in (\j in co ^ o -i

c\j o >o in fh <4

- *-< O' r—H nO O '“l »0

rg r-i vO in o

m —* r- *-h rn in

in —« m o ^ cm m

in co so 'J' r\i —i nj 1)11 I >t | >t cm so

I

r\j m r-* o

i r i cm 4* c\j | in

I I I

cm

I cm cmI

cm rf\ *t in vO O —I CM ro vt in

CM CM <M CNJ CM CM

vO N CD O' o

CM CM CM CM m h fM m >t

ro ro m co in so n- co O' o m m m m m <}- —• cm cn s^- 4* 4" 4- 4* in so

4- 4-

Is- CO O'

4* 4* 4" o

in

< o

33

151

.180+01
Figure 3 is a plot of statistics from the h = 2 layer

against d* . The "++- OCTANT ( S ) V and the reflection class

indicates the statistics are those of the form

Again the statistics show a variation with d*, and the

linear slope of 3.9 indicates a more severe trend than that

for the lk 4 statistics. Absorption could well be the cause

since it is known to have a larger effect on those reflections occurring at larger d* values (higher 2© angles)

.

Non-uniform scatter of the statistics is also indicated.

39
Figure 3 *x#

(NJ

O

r*- -J 4- CL

(NJ z I ~

00 4- | K 4 4 CL * X 1 * vt 1 in Z

ro < 1 1 —! X 1 « X 1 O 1 r-U UJ 4

CL CL 1 < a 1 K 2: 1 0O 1 o 4 LL 11 (NJ 1 CO X 1 x 1 X • LL 4

00 CO 1 Z 1 LJ 11 < CO 1 o X 3 1 < O II O 1 U 2* 1 CL O ac 1 —

ro

Q j X 00 00 # LL 4

CL k—t X # r LU X X ro 1 > X O # Ch 1 w —«

CvJ f0

21 z X oj « i O M < s0 in 1 CL • 1 X (/) O H 1 0 H < —1 1 Z 3 I 4 X — Of 1 0 O 3 1 -J CL 00 z 1 CL UJ 4 < 1 O' CL 1 < (NJ LU 1 3 * •—« 1 a # K 1 00 X < 1 00 CL 4 z X tt O 1 0 lO in < 1 < s0 3 1 K LU CM <3 1 (3 -J e 1 LU s0 00 1 -J 1 CO .J 1 LL 1 1 LU 1 a. X CL 4 # O 1 LL <5 4 LL 1 O0 O O 1 2T X ro 1 O K • # r~l 1 >—1

LU <0 >4 1 X 00 <0 • 1 O <0 in 1 LU 00 in 1 U to* • 4

LLX 0 1 LU X 4 m i a: l 1 z in 1 NX 0 LU 1 ro z 1 00 0 4 1 J X in -J 1 Z 1 x 1 nn co 00 4 * 0 -J 1 CL # 1

11 CL 1 3

(NJ CO 2:

O 1 CO 0 LL I in 1 00 0—1 1 00 tax

II CM 1 X 2: >4 1 O Q • 4 u in 1 CL 00 l w 1 Z 1 t—i

LU 0 l Z x 00 l X X CO NX 1 Z -J 3 1 n-1 O h- 1 (/) O 00 BH 4

ox CL LU 00 LL 4 X CL UJ _t z X < CL IX O

<.

X 3 < O 4 X O 3 CM O f\l 00 O 00 (O O 00 00 h- X X 00 1 > () X z CD 4 CL <1 t/> Q 9 4 LU LU < O > -J LU O

LU -J

O i i i * i i i \ + i r i i i i i i I + i r i i i * 3 • 4 _l • f\| 4- _| * (NJ -J •

(O + —

J

(NJ cm .j

f0 3 -J m • -j

<NJ <\J _J e N'D J

CM (NJ

ro * _j 3 4# • _l

(NJ 4 * _J

(NJ # CM

(NJ # (NJ t (NJ # —I 3 • m (nj m

D(0 -J D + N J 3 4 # • -J 3 (NJ # -J

ID 4- (NJ _J • 3 4* fNJ ZD (NJ (NJ 3 • -J (NJ

(nj ro 3 • 4 # -J 3 4 (NJ —I 3 • m 3 3

ID

(NJ

J

(NJ (NJ (NJ _j 3

(Nl 4 # 4 • 3 3 (NJ • -J # 4 # 4 # 4

(NJ # 4

(NJ 4

(NJ I I I I I I I • _J 4 * i i + i i i i i o

4* O'm

COoo’ 43 I I i + I I I i * 4

(NJ O4 I

(NJ 4DOI I I I f 4 I ^ I I I \ t I I 4 I I I I + I I » CMO4r—* 't0

CNJ CMt1 I I I t I \ +

(NJ o4O

CMo0

CO1
=

DM

(II

(ORDINATE),

SYMBOL

IS

.

B

=

DM<!)-

L.

S.

LINE,

SYMBOL

IS

*

C

=

DM

(

)I

-

L.

S.

CUBIC,

SYMBOL

IS

=

ABSCISSA

VALUES

(DSTAR,

LOG

F**2/SIN(2*THETA)

,

OR

QUANTILE

SYMBOLS

L,

U

ARE

FOR

AVERAGE

-

AND

+

DEVIATIONS

(A-B

VALUES

FOR

12

POINTS

IN

THE

INTERVAL.)

POINT

ABC

D

1

2.315

-3.484

-3.093

.466+00

Figure 3

contd.

o O o o o o O o o o o o o o o o _ p_4 pH ,—H ph pH pH pH pH pH pH H pH pH pH pH pH pH pH pH pH pH o o o o o o o o o o o o o o O o o o o o o o o o o o o o o o a o o o o o o o o o o o o o o o o o

4> + + f + 4* + •f + 4- + + + + + 4- + + + 4- + + 4* 4- + 4- + 4- + + 4- 4- 4- 4» + 4- 4» 4- + 4- 4- •f 4- CO X co m co 4- X 4* pH o CM 4- O' ao m in o CM X o pH 4* 0 X o cnj m o m in >0 X O' <\j •+ X O' m X O 4- X O m

CO CNJ >0 4- 00 o m x 4* in X o (NJ in X o o o o —

h

r—4 f—

1

+—

i

CNJ CM Cvj CsJ m m m m ro m 4* sh 4- 4- n in m X -0 X X r- h-

<4* in in X X X h- r*- co X CO O' O' O' O' H pH —

pH — pH ph -h p=< pH *“H pH pH -H p=i pH pH pH pH pH pH rH pH pH pH PH pH pH pH pH

cn vO o (\i 0s st S CO fO N © ^ pH CM m m (Nj p< -H

cm co cm 't oo >fr ro —• ph 4* r*- in O CT' IT) (\J O O' ^ O' h O CO ^

x in cm cm m O' 0s

Is- in >t ro o r*^ o

o o co r- oo X O'. X h- 4* m

4 NhOOh

f\j CD o m o <t m

O in o 4 N «t m m >0 h- x X

^ h ^ iTi in cn ^ N O O' 4 H

O' O' O O ^ O' O'

CD O' 4 o m O O' O MNJ N (tl O

CD N >0 (O M O 4

fri m m ro m m ^ (\j cm nj r\j nj h h h ^ | | | | | | -h pH I |

I I I I I I I I I I I I I I I I I I I

0'ri>O 0'n^aDC0 f\jincDnj

O' 4cooh^4rDmc\jsin mcMOo*r->oin4, moo'X

mcnmcMCMcMcMCMCMnjpHpH

r i r i i

i i i i i \ i

morncoovJcrfO'fvj't'fcD't

incccoir\cor\J'£>4-cOvOmcr'

f>->oiri>j-roro^-(oc7s cor»NO

pH H H pH pHI H pH | | | ^

i i r i i i i i

vor--a'inor^vocor^^4-r-i nj'tcoin(\iO''tOco^Mm Din4^(\OOOOHfMfCl

I I I I I I I I

O'O'inxr^oinpHxcNJinx inoin^njinco^^mino' 4ins0fwcO0'O'Hnjfn4, in

coi^inr^rvjo'CDOoO' com {MXin'OO'O'incMmXpHO

in^op^mtT'O'incDLnr^in^

l irincvjinom^^M I

minr^vOcNjO'O'ininmpHcn ino'NOrnirifrKNjhO'coin

cm

I I I

m<\jommaor^-^>t-*o0't ino^'tO'coOino'^cOH oooinfoinoiAOo^f'Jco

(NJ | pH | |H H (\| >-4

| <\J

I I I

•4- mQ'i^'fr'00> 0'00^''—

i

4-m^4-coc\jO'h-a'r-om 04-c\jooocNJomr^sOro

*—

4

F—I »—( CNJ | | -H pH pH pH <\j

I I I

r\j m m o r*- co cm m 4- in X r- x O^O—iCNJcn^^vor^co

phimcmcmcmcmcmcmcmcm

O' o h (\i co 4 in

cvj m m m m m cn

>0 r- m m

GO O' m m

o ^ cnj m 4* 4- 4^ 4- 4^ 4" in

*4" X

4*

00 O'

4“ 4-

41 < Q
3. Log (Fo^sin2 9 ) Plot

Figure 4 shows a plot of all statistics against

log (Fo/sin29) . in this particular case we are testing two

equivalent sets of Fo values for anisotropy in extinction.

This type of plot could also show absolute effects of extinction ( i . e

.

, between Fo and Fc) if least-squares residuals

were used, because in this case the abscissa values would be

log (Fc 3 /sin2 9 ) . The plot shows that negative statistics tend

to occur at large abscissa values. The plotted U and L

symbols show that the largest |af|/cj(Fo) values tend to occur

at large values of Fo/sin28. However, sin2 0 = cos§ sin9

and thus l/sin28 = 2/(Xd*cos8) . Thus the variations of the

statistics with d* and Fo^sin28 are not independent.

We need the 51og[2/ ( Xd*cos8) ] / ^d* in order to relate

1 = 36 mhk/9d *

• Sd*/^logf2/ ( Xd*cos9 ) ] . Then ^log (2/Xd*cos 9) /Sd*

= -1/2.3 d* and its inverse is -2.3 d* , showing that the

inverse relation between the two variations varies with d*

.

Thus the slopes of plots against d* should be of different

sign from slopes of plots against log (Fo/sin2 8) . Here the

d* plot shows a positive slope (negative statistics tend

to occur at small d* values) while the log (Fo^sin2 8) plots

shows a negative trend as expected. However, if the

statistics vary in a systematic manner with Fo 3 we would

expect a more pronounced trend in the log (Fo^sin2 8) plot

because the range of Fo 3 values is 50-100 times greater than

the range of l/sin29. The ordering in Fo^sin28 plots is

primarily due to Fo 2

„ A comparison of the d* plot (Fig. 2)

with this one shows that the cubic curve reaches an extreme

of -4 for low d* values while the corresponding value here

is -2.6. This indicates that the disagreement between

intense equivalent reflections occurring at low d* values

is not worse than the disagreement between equivalent

reflections occurring at low d* regardless of their intensity.

the variations as ^Sm^^Slogr (2/\d*cos9)

42
ALL

REFLECTIONS

INCLUDED

PLOT

OF

DM(VERT)

AGAINST

L

OGF

**

2/S

I

N

(

2T)

FVERYI

1

1

TH

POINT

IN

PLOTS5527

POINTS

IN

THIS

SET

OF

REFLECTIONS

50

POINTS

IN

THIS

PLOT

LEAST

SQUARES

LINE

IS

0M=

2.350

+

-.83l*X

—

LEAST

SQUARES

QUADRATIC

IS

DM=

-2.559

«-

2.387

*X+

-.507

*X**2

L

FAS T

SQUARES

CUBIC

IS

DM=

15.222

+

-15.048

*X

+

5.004*

X**2

*

-,561*X**3

GOODNESS

OF

FIT

IS

3.557

FOR

LS

LINE,

3.549

FOP

LS

QUADRATIC

AND

3.542

FOR

CUBIC

FIT

IF

MORE

THAN

4

PTS

IN

PLOT

Figure 4

+ I l l I I I I I I + l I I I I l » I l + 1111 CM

! i I ii + i»itiiiKi + iiiirriii +

+ *

3 3 3

* 3

33

3

3 33

3

3

CM 3

3 3

3

3

3

3

+ *

+ #

+ *

+ *

CM

CMCM # +

• CM CM * +

<N

(Minm

• CM

CM * + # + 3 # +

O'

*4

*

I I I I I I I I + » I I I I I I I I O

+

I | I I I I \ V I + I I I I I » I I I

54297+01+
=

DM

(II

IORDINATE),

SYMBOL

IS

.

B

=

DMC

I»

—

L.

S.

LINE,

SYMBOL

IS

*

C

=

DM(I)

-

L.

S.

CUBIC,

SYMBOL

IS

=

ABSCISSA

VALUES

CDSTAR,

LOG

F**2/SIN(2*THETA)

,

OR

QUANTILESYMBOLS

L,

U

APE

FOR

AVERAGE

-

AND

+

DEVIATIONS

(A-B

VALUES

FOR

III

POINTS

IN

THE

INTERVAL.)

POINT

ABC

D

1

.175

.755

.807

.192+01

Figure 4

contd.

+

OOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOO + + + + + + + + + + + + ++ + + + + + ++ + + + + + + + + + + + + + + + + + + + + + + ++ + + + + som.-4>0'-*4-co(Mincoominr*-ocM4-r-onj4-r*-a'(\j4‘'00'<\J4'r^-cMino>(MsO*-iinoina'm-*»-4*-<—

<

4

- 000 ^

o»-<(\j(\jcOfOfo>tNt'tinLnLninsOsO'0'Or>-r^h-r^-r-cocooocoa'as ON ooo^^(\jc\irofnrO'4-iAvor*-coO''-, c\jco (NOsj<Nf^(NfNrvj(Njrgrvjrvj(\jf^(\j(NjcNjrvif\jojf^c\jc^OsJfsjrsjrvj(NfNJOsj(Nrnrn(^rOf^fnf^rnfnrorrirnrornrorn't>t't

m ffi co (\j ro

in m in o in

in *4

- ro m c\j n 'O m

nj o «o n h

CM »—

4

*—4 r—

4

^-4

00 CO O 'O O

0s co o o o o o o

rsj r\j in co cm vO in 4* rO <\J cm »-h O O O O o o o in ^ ^ ^

o o o o

o o o o o ^ co rn o 4} h h n tn m

o o o o o

cO O' h- O (\j 4“ 4- in r- co o o o o o

in in ^ o

cn h rno o O H H fH (\|

O' ^ t\j N O' O O O' CO O' nj m >£) O'

C\l 4* co co

4- hI I I I I I I I I I I I I I I I I I

I I I I -« TM

I I

n-cocoinrovomrocoo m^rH^rnon-'t^o

•Oinin4-4-4roromc\j

OinonJNnjNl-HO

h- 4-m—<aof^inco*H

C\| C\| (\| (\J <—I »—I r—I i—

I

r4 l-if-iHH 43 cofnin't OM^cnnonjinhO oooooooooo

I I I I I

^inm'Om'OO^ON

cm in co —• 4" r*- —i in co m

r t i I t i i I i t

^HsO'tm-iNinOvt

vtmLnsOh-aoO'OCM'O

h I I I I I I -i -«

I I I COOOOOOM^MO invOnjHO'O'Onj'Oco

r-'OsOc\ioo 4-m<\ifor^

I ^ ^ tvj

I * ,

I

'OhNO'>oinco*HCO

in 4 inHcor-r-) 0'CD

Nj-osjO'rHtNjcNjvocO'cn

I f-l M H I

I

<MinQ' 4*^c0 tf'C7

N O vOO inn^oM^'OcninO'-i

CM C\J (\J CM CM r-4 r-4

l III

c\j(\jc\Jincoc7'fs-'t '0 O^nKiN'tnjniri njcOHoOinsOOm

h m h h c\i r-4

I I

oofoomNtnomcD

covD 4 coo'n- 0' 4 ino HHinocD'OhOOnj

cm cm m i Nin i so III l I

cm m 4- in \0 n- co o o C\] CO 4- m s0 h- co O' O h nj m 4

r-1 CM (M CsJ (M C\J m

CM CMCMfMfMfOrOfOfO 4 in 'O h co O' oro co m co ro co 4

* >—

1

cm ro 4" in

4- 4- 4^ 4^ 4-

vO r4* 4- co O' O

4- 4- in

< o

44
Figure 5 shows the variation of statistics from the

h2 l layer with log (Fo/^si^ 0) . The statistics, ^>m

]r

-L 2

' are

from the F, _ „ and F - . layers so we are

h2 1 h2 l

comparing the two octants of the Fh, quadrant. The

statistics tend towards negative values for large Fo2/sin2 0

values with the cubic curve showing that the trend is most

evident for the 63 largest statistics (7 plot points x 9

reflections per point = 63)

.

This is about 15% of the data.

The overall trend is very similar to that observed in the

overall plot except for the very large and very small

Fo// sin2S values where it is more severe.

45
+++

OCTANTIS)

H

2

L

REFLECTIONS

PLOT

OF

DMIVEPT}

AGAINST

LOGF**2/S

IN1

2T)

EVERY

9

TH

POINT

IN

PLOTS

438

POINTS

IN

THIS

SET

OF

REFLECTIONS

49

POINTS

IN

THIS

PLOT

LEAST

SQUARES

LINE

IS

DM=

2.309

+

-.829*X

—

LEAST

SQUARES

QUADRATIC

IS

DM=

-9.459

+

6.754

*X+

-1.174

*X**2

LEAST

SQUARES

CUBIC

IS

DM=

49.253

+

-50.463

*X

+

16.615*

X**2

+

-1.825*X**3

GOODNESS

OF

FIT

IS

3.006

FOP

LS

LINE,

2.955

FOP

LS

QUADRATIC

AND

2.881

FOR

CUBIC

FIT

IF

MORE

THAN

4

PTS

IN

PLOT

Figure 5

+ I I i i i I I I i + tliliii 3

I I + I I I I I I I I I I I I I I I I I I + I I I I I • +

I I I +

1

1 3 • + * 3 4»

CXJ * 3 + m

1 0 3 3 (XI C\| (XJ 4- 3 4-

1 • 3 4" * -J

1 3 4> *

1

1 3 3 (XI (XJ

1

I 3 4- (XI mJ

1

1 3 4* # -J • +

1

1 3 4- * -J

1 e 3 c 3 >fr «=J -J

1 3 m

1 3 (XI

1 3 « (NJ

1 3 3 • (Xj CXI -J O

1 4- 3 fr f • —1

1

j 3 + (XI

1 3 3 exj (XI (XJ (XJ

1 (XI 4* 0

1 (XJ 4 • -J

1 3 0 * + «=J

1 eg 4

1 • 3 * 4= i

rn 4* _J

1 (XJ * 4- -J

• 3 # 4-

I 3 •fr 4s (NJ +

4* CN

+ 1111 r—4 o

4- o

>*

CXI

till i i i » i \ i i r +

+ i *—) o

4*

Is- 4-

3 * + 9 _J 3 m

•4* -j

o

•f O'

<*

Is- sO

46

19206+01

.24839+01

.30472+01

.36105+01

.41738+01

.47371+01
=

DM(

l

)

(ORDINATE),

SYMBOL

IS

.

B

=

DM(I»-

L.

S.

LINE,

SYMBOL

IS

*

C

=

DM<

I

)

-

L.

S.

CUBIC,

SYMBOL

IS

=

ABSCISSA

VALUES

(DSTAR,

LOG

F**2/S

I

N(

2*THET

IA

,

HR

QUANTILE

SYMBOLS

L,

U

ARE

FOR

AVERAGE

-

AND

+

DEVIATIONS

(A-B

VALUES

FOR

9

POINTS

IN

THE

INTERVAL.)

POINT

ABC

0

1

.175

.716

1.432

.192+01

Figure 5

contd.

P-4 r-R —

,

P-4 rR ^4 i—4 —4 r-R *-» i—4 _4 —4 P—4 rR _4 _4 1—4 rR p—4 —4 rR p-i r-R P—4 rR r^R r-R r—4 —4 ^4 rR r-R rR rR ,_4 rR r—

1

P—4 rR r—4 o o o O O O c o o o o o O O o o o o O o O o o o o o Q o o o o o o o o o o o o o o o o o o o o o

4 4- 4 4* 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 + 4 4 4 4 4 4 4 4 4 4 4

CM >0 in — O >r >0 —

<

m 00 O CM >a* n- 00 'S’ O' m >0 O >0 CM rR •"4 m r- O' m «n >0 co r^ m r-R O' m 00 O' »—

4

m O n0 >* -R —4 <M m m 'J* in in in O >0 >0 >0 v0 r- 00 CO co O' O' O' O' O H rR CM CM CM m m m m vJ* in in -0 r- co o r-4 CM in

CNJ CNJ c\l <M CNJ CM CM CM CM CM CM CM CM CM CM CM CM CM CM CM CM CM CM CM CM m m m m m m m m m m m m m m m m m m >a- «•*

fO O O fVJ

r*- ^ in go ^ m o o

co <js cm in m cm r- o' CM CM CM o r\j m

cm m m in m m m m

m rn o

cn ir\ >}- sf m m m m mooCD n0 m cm cm

oo

>*<M

*C >*co

<M rH

r\j cr 'J* co cn ^ h o 4- in 4 ^ o o

N -< (\l h ff) >0 n n o m in co

r4 r- ^ C\J (\J (\J

>0 lA >fr in 'O o h n >t in m m m m mO O' N O' cohinoo

m m m cm o

O' in m o ^

vO h vO m

m O' nj cm >t

i i i i i i r i i i i i i i i i i i i i ii | -4 >0

I I I

n?- vO vj- m co co in h ^ o h m m ^ n ^ (T

oo o' Vs

Is- m

co o (nj o O n

CM cm CNJ f\J -4 O O' CM in r"- >t O' co m o

-R -4 o o o O

I

inr^'Om'Ovt

CMrOsOCT'CNi^ o o o o —< -r

f I I l l I

in in co O' nj co in O' >o in v0 h h M di m fo

I I I I I I

O O' O *4 r-4 »4 o h in n co 0s *4-

I I I I I I

o co O' O' h in

nj m O' vO o (\i n in in -o h- od

t r i i i i

^ —i (m r— r— (\j

| r-4 rR "R r-4

I I I I I

Ovomr^onm^in>j'm'0^ocor«-rw otn'4'>ocommmminav '4*nj^'4'^'-*r^aocom+-*c\)NOr^C''4-mosjaN (Nj rg>4'm-<a'«o-‘r^cocO'aDr^>oh-inr^Oc\icooao'i’oooor^mcr'Onjinc\jNO<Ni+-HcomOsCCNjincNj^Hmcnmr^in '4-sOm«omcocoinocoocMf^f^-4->o-oomcostr^om^ocop-i(Mmin^G^^nm^inrr. voroma'OdiniNj^ascNjvt

I \ | —i | | | | | r-4 r-H 4 | | | ^(Ohminj | m | nmcom

II I I I I I I I I I I I I I I

c\jm>tmsor-cocro — c\jm>tcnvor^coa'0 sOh-co(j'0^c\jm'Xinvor-cocro^TMm'Tin>or^cou'0'-<c >jm'rin HHHH(\j(\j(\jnJN(\j(\jnj(NjN(n(nri(nfn(nfnfnffifn>j-'l4>t>t'tN co O'

>*• ^ ^

47 < c
4. Half-Normal Probability Plot

Figure 6 is a plot of statistics, J>Pt = p t (set 1) - Pi (set 2) / [s^Pi (set 1) - o

2 pi (set 2)] 3 against halfnormal quantiles. Half-normal quantiles are used because

the order of subtraction of sets 1 and 2 is not meaningful.

The Pi are fractional coordinates from structural refinements of data sets from two separate crystals of

3 -Ca 3 (P04 ) 2 . Each crystal contained about 5 wt% of Mg.

Thus, this plot tests the sensitivity of structural

parameters to sample differences since the same experimental procedure was used in both cases. The plot shows

that the actual 6p t (represented by •) wiggle about the

least-squares line (represented by *) . This is mostly

due to the paucity of data (43 points) . Indeed the plot

shows that the cubic curve fits the variation only slightly

better. The 6pi represents mainly random differences and

shows the samples are nearly identical structurally. The

slope of 1.275 for the line suggests that errors in the

parameters, n(p)iare too small by 20%. These errors

are estimates from the least-squares and are likely to

be too small if the parameters in the model account for

errors in the measured Fo.

48
DM

QUANTITIES

BASED

ON

X.Y.2*S

PLOT

OF

DM(VERT)

AGAINST 1/

2N0RMALQU

ANT

Figure 6 + i i i I i

ID + § i + u. 1 O H | CM•1 lL i i U 1 0 1 H Z> 1 3 II V 1 -1 2 1 a O a: 1 n O 1 co c

* LL * 1 z X 0 1 K V * O' 1 M O 1 Z K n — < N 1 a • 1 <0 O i 1 K < 1 z 3 +

N*

O 1 o O 1 a co Z | UJ + < 1 m CL 1 < CM <J 1 D * M O * hCO X <a H * G

CO

X < < ** D

UJ

N

C

J • CO 11 J X cr • * o 0 + u. N

CM

X 0 • * ro

vO • CM

P0• + •

CO UJ o z • + 1 J o CO -1 o J 1 • 1

II 1 cr l X O 1 G U. 1

CO 1 o | C/) II 0 | K X 1 O O • + Ja CO 1| z 1 CM UJ u 1 Z M CO I K M 0 ~ 1 Z J D 1 ** U H 1 0 CO « a UJ (0 lL +

CL UJ **

X < a LL o K D < O

<3 ->

O

CO o CO ro

CO (0 CM

k UJ

n > CO H z CM a < (0 Q •

LLi UJ < O > -1 UJ O

LU -i I I I + 111 + I I I + i i i i i r i r + i i i i i i i i i + i i i r i i i • * • + * • +• + *+ *+ * • + * • + *+ <M+ C\J

<M + • # +

CM +* # +* +* +* + * +* + * + *

C\l + • CM • I I I I I I I I + i i i r i i i i + i » i i i » • i • + i i i i i i i i i i i i i i i

CM

CM *

CM *

CM I I <

56614 - 01

+4
DM(

I)

(ORDINATE).

SYMBOL

IS

.

B

=

OM(

I

)-

L.

S.

LINE,

SYMBOL

IS

*

C

=

DM(

I)

-

L.

S.

CUBIC,

SYMBOL

IS

:

ABSCISSA

VALUES

(DSTAR,

LOG

F

**

2/S I

N<

2*T

HET

)A

,

OR

QUANTILE

POINT

ABC

D

1

.011

-.057

.006

.193-01

Figure 6

contd.

+

o o o o o o O O o o o o o o O o o o O o o o O O O o o o **•1 -* rH -« f-t •”4 _1 o o o o o o o o o o o o o o o o o o o o O o o o o o o o o o o o o o o o o o o o o o

1 1 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + 4- + + +

CO 4 OJ o CO yQ 4 (M rt O' 00 CO CO co O' ro y£) o 4 O' in (M o o o ro ro CD ro O' 4 N in ro ro 4 O'

vD 4 o ro in co —> 4 N O'' (M in CO *-« 4 N o 4 r^ o 4 r^ o 4 GO CM M) o 4 co o o T-l (M ro ro 4 in vD N GO 4 N •* (M (M (M CVJ ro ro ro 4 4 4 in in in vD v£) yQ N N r- GO CO O' O' O' ** -* -* **

oo ro o o» -« id ^ ro in vo O' 4 O O Q O O H H

v00' r0v0o4C0CM ^^(McMrororo4

oo in ro ro »£ o \0 4 vo^vo^'OcvJNm

^lfilf)'O\ONNC0

0'in^coincMO' voroo cDoioonW^n^ in rovOr^44(\jO''0(M oomcMO'inocM^cM

in\£)NNG0(7'C7'O'00

»-« 4 O' 4 OvO^CO

(M'-*4ao(MinO'CM oooo«—*—

I

4

, H(Di/)nHo(7'0' 'OoroN ^ in O' cm vo cm ro ro ro 4 4 4 in in

0'o<Mino'roo''0 oino'DNW'DvOvOvONh-CCOOO

vO~*vOCMf^rOO' vO O'oo^^-.cMCMro <tmin'ONooO'-*ro

O' 4 4 O' o n <c (M 4 (M CM 4 vO GO o o — ^in O' m *-< oo (\j O' O' o ro 4 »n vo ^ -• C\l CVJ (M CNJ CMfnO'oO'-fr-a-N 'Drno'CM'OciN

cm ro ro 4 m 'O N

oo cm in cm 4 4 o

GO 00 o ^ (\J ro

O' in h- n cm o» O'

ro ^ ro co CD GO o

ro 4 4 4 in in \o ^ N O ^ W h

VO CO CM CM o in © © S N © O'

(MIO^©©N©0'0 CM ro 4 in VO Is- CD 0'OH(Mm<tm©s©0'o

•^CMCM(M(M<Mc\]CMCMC\J(MfO ^ CM

ro ro

ro 4

ro ro in vo

ro ro ro

GO O' o ^ (M

ro ro 4 4 4

ro4

ii li < Q

50

COO
III. Detailed Description of Main Algorithms

A general description of the program STATUS and its

subroutines has already been given. Because STATUS consists of

a main control routine and 35 subroutines, it is not practical

to describe the logic of each routine in great detail. Instead

we have provided macro-flow charts for the control routines

and those routines likely to be modified by a user, typically

those dealing with input and output. A description of each

routine includes its purpose, its procedure, a definition of

the important variables and a macro-flow chart detailing the

logical procedure. These aids together with the FORTRAN

listing should enable a user to make any required modification.

III. A. Main Routine STATUS

Purpose

Controls the execution of a job by reading from cards

or their images the control parameters discussed in the

section on Card Inputs.

Procedure

Reads job title card (stops on EOF) , a card with cell

constants, then control card and finally a card with experimental parameters. Checks control parameters for consistency

and stops if some are out of bounds and cannot be corrected.

Calls subroutines to read in data, calculate statistics, etc.

When job is finished attempts to read another job card.

Important Variables:

See section on Common Blocks and Arrays.

See section on Card Inputs for description of control parameters.

JC - total number of statistics, JA in set one, JB in set two,

JC = JA + JB.

JZ - pointer for drum address of Fc or absorption path length.

M - packed word indicating class, layer and octant.

M - 1000 X class type + lOx number of layers + octant option.

II - number of bit in word LOGHKL to be examined in cataloging

of reflections into classes.

LA - number of atomic parameters in first set.

LB - number of atomic parameters in second set.

IP - number of parameters common to the two sets.

KIND - type of parameter; 1 - occupancy factor, 2 - coordinate,

3 - thermal.

51
1. Common Blocks and Arrays

a. Parameter Variables — these set the

array limits.

NN - Maximum number of 6m or residuals

(7200 at present)

.

MM - Maximum number of atoms, 100 at present

(equivalent to ~1000 parameters in

least-squares)

«

N - Maximum number of points in plots - restricted to 50 by PLOTS.

NNA - Maximum number of reflection records in

storage

.

NNB=NNA + 1 - Boundary for overlay by equivalence

statement

.

NNC=NN/2 + 1 - Boundary for overlay by equivalence

statement

.

NND=NN/2 + NNA + 1 - Boundary for overlay by

equivalence statement.

b. Common Blocks

Block/D MAXH, MINH, MAXK, MINK, MAXL , MINL

Maximum and minimum values of the Miller

indices

.

Block/H/DUMMY (NN) - Used in equivalencing arrays.

Block/G/ASTAR, BSTAR, CSTAR, COSAST, COSBST, COSGST

Reciprocal cell values.

Block/lO/lN, NOUT, NDRUM, NTAPEA, NTAPEB

.

IN - Logical unit for card reader.

NOUT - Logical unit for line printer.

NDRUM - Logical unit for mass storage.

NTAPEA - Logical unit for input data set 1.

NTAPEB - Logical unit for input data set 2.

52
Block/ORDER/JJA, KKA, LLA, NR, NFILER, NABS

JJA, KKA, and LLA give order in which

Miller indices are to be packed in JKL

word, i.e., 2, 1, 3 for JJA, KKA and LLA

gives JKL = 1,000,000 (k+200) + 1000

(h+200) + ( i,+200) . The indices are packed into the JKL word to save

memory storage. Variable index order is an

option which is helpful when reflections are being

matched efficiently.

NR - Output unit from SCALE

NFILER - File number on unit NR

NABS - Flag for optional absorption

correction

.

Block/TI TLE/ATI TLE (14)

ATITLE - 14 Hollerith words set by STATUS

and PLOTEM to provide headings for

the various plots.

c. Arrays listed in groups requiring similar dimensions.

1. Arrays holding reflection information:

Dimension NN; IDA, IDB, NNA ; FA, FB,

SGA, SGB, JKLA, JKLB, KA, KB and LOGFC

.

2. Arrays holding atom parameter information:

Dimension MM; NAMEA, NAMES, OCA, OCB, XYZA,

XYZB, SIGXA, SIGXB, BA, BB, SIGBA, SIGBB,

SIGOA, SIGOB, MATCH, LOG.

3. Arrays holding statistic information:

Dimension NN; X, LOGHKL.

4. Arrays holding plotting information:

N; DMSMAL , XSMALL

.

53
d. Equivalences to conserve storage:

1. LOGFC, LOG: LOGFC is written over when

parameters treated.

2. DM, IDA, X, IDB: DM and X are written over

IDA and IDB.

3. LOGHKL and JKLA, JKLB, KA, KB.

Note: These equivalences work only if MM ^ 100

i o e . , there is a maximum of 100 atoms

per asymmetric unit at present.

54
iii.a.2. FLOW DIAGRAM OF STATUS

(j5 TART J— Initialize counters and plot titles for X

CALL PRIME to assign logical unit numbers (SMread TITLE card for job

-( STOP

read cell parameters in cell and

calculate metric tensor.

read program control parameters, check for errors, and print out what program will attempt to do.

rp or INFOFC to I

, Fo, Fc, g-frcT

CALL SPEC to minimize 2. Hi.1

-UiU- I ,

| CALL SORT to order statistics

|4

CALL PROB to calculate normal quantiles corresponding to the <f # CALL PRYNTT

. - |.

to print —out desired statistics and associated reflection information.

CALL PLOTEM to produce normal probability plot

CALL ANISO to produce normal probability plots for specified reflection classes and octants.

CALL DMDQ to obtain d* values and

catalog reflections if required.

CALL SORT to sort d* in order of increasing value.

I CALL PLOTEM to produce plot of statistics against d* values .

CALL ANISO to produce plots of statistics against d* value for

reflection classes and octants.

CALL SORT to sort O values in order of increasing magnitude.

CALL PLOTEM to produce a plot of statistics against Q values.

'reflection' classes

? Y

CALL ANISO to produce plots of statistics against Q values for

reflection classes and octants.

55
56
III. B. Subroutine ANISO

Purpose

Controls the selection of statistics for plotting

subsets of the complete data set.

Procedure

Each statistic has associated with it the LOGHKL word

that specifies all classes, layers and octant that the statistic belongs to. These may be plotted separately so comparisons

can be made. Subroutine ANISO "tells" subroutine PLOTEM what

to plot in the following way. If the variable MCLT in the call

to PLOTEM is positive PLOTEM will include the statistic in a plot if the Ilth bit of the LOGHKL word is non-zero. If an

octant is specified via MCLT, PLOTEM checks bits 33, 34 and 35

to see if the statistic belongs to the specified octant.

Variables

M - a packed word which specifies class, number of

layers and octant option.

MT=M/1000; = 1 standard classes, = 2 hOO layers, = 3 OkO

layers, = 4 004 layers, = 5 standard classes

plus hOO layers, etc.

MLAY - number of layers to be plotted, up to 10.

MOCT - octant options.

= 0 plot all octants together.

= 1 plot the +++ (h positive, k positive, 4 positive)

and octants separately.

= 2 plot the +++, , and +-+ octants.

= 3 plot the +++, , ++-, and +-+ octants.

= 4 plot the +++, ++-, +-+, -++, and octants.

= 7 plot all eight octants.

NOCT - is the octant number, one through eight, where octant

1 is the +++ octant, octant 2 is the ++- octant, etc.

Octants are represented in the LOGHKL word by a 3 bit

configuration with zero for positive indices and 1 for

negative indices, i .

e

. , -H— is 001.

MCLT - "hundreds" digit is the octant number, NOCT, "units"

digit is the class or layer type, MT = 1, 2, 3, or 4.

II - number of the bit in the LOGHKL word to be tested,

0th through 20th.

IBEG - first bit to be examined (the 11th if layers only)

.

ISTOP- the last bit to be examined.

57
Subroutine ANISO

58
III. C. Subroutine DRUMRD

Purpose

:

* Reads a reflection record from word-addressable drum

and unpacks the four words comprising the record.

Procedure

:

Sets a pointer, JPOINT, then calls UNIVAC routine

NTRAN to read in the record. The pointer indexes from the

beginning of the drum files and is equal to (1-1) *4, where

the I-th record is wanted. The listing is simple enough

for no flow chart to be needed.

Important variables:

JHKL - packed word containing the indices.

F - value of Fo

SG - value of corresponding sigma.

NT - flag indicating unobserved or some other qualityindicator .

I - the reflection record address.

Subroutine DRUMRT

Purpose

:

Writes a reflection record on word-addressable mass

storage

.

Procedure

:

Packs the indices into one word and places this packed

word together with NT, F and SG into the buffer FDRUM and

calls UNIVAC routine . NTRAN to write out buffer.

Important Variables:

N - flag indicating whether reflection is "observed" or not.

IJ - array holding the indices h, k, and l.

F, S - Fo and corresponding sigma.

Modification Aids for DRUMRD and DRUMRT:

Any procedure that operates on random-access, word-orrecord-addressable mass storage is desirable for these two

routines. Specific procedures depend on local computing

facilities

.

59
III. D. Subroutine FI TCHK

Purpose

Calculates values corresponding to the least-squares

line, quadratic, and cubic fits for statistics or selected

subsets of the statistics. If the plots of the statistics

against d* or Q are being made, obtain the scatter of

the statistics about the least-squares line.

Procedure

Examine input flags, if abscissa is normal quantiles

leave out the tails. Select statistics according to octant

as specified by NOCT and class or layer as specified by II.

The plotting subroutine handles 50 points conveniently,

accordingly there are K/50 + 1 points per interval where

K is the number of statistics in the subset. If there are

at least 6 points in the interval the positive and negative

deviations of the value calculated from the least-squares

line are computed. The average positive and negative deviations can be calculated and plotted with the statistic,

least-squares line, quadratic, and cubic values for the

midpoint of the interval.

Important Variables:

LHALF - one-half of L, the number of points in the interval.

DM - array holding JC statistics for the data set.

X - array holding JC abscissa values.

DMSMAL - array holding subset of statistics to be plotted.

XSMALL - array holding corresponding subset of abscissa values

K - number of statistics in the subset.

M - "hundreds" digit specifies octant number, 1 through 8

"units" digit positive if classes or layers wanted

(subsets of statistics)

.

II - number of the bit in the LOGHKL word to be examined.

ALINE, BLINE - coefficients for linear equation.

ACURVE, BCURVE, CCURVE - coefficients for quadratic

equation

60
ACUBIC, BCUBIC , CCUBIC, DCUBIC - coefficients for cubic

, equation.

SUMLN, SUMCV, SUMCUB - goodness-of-fit for line, quadratic

and cubic equation.

NOCT - octant number, 1 through 8.

61

III. E. Subroutine FREAD

Purpose:

Read in an observed data set, check for input errors,

call MAXHKL to establish largest indices, store (for

efficient sorting and reflection matching) slowest varying

index for later identification of members of batches based

on this index, establish the minimum and maximum values of

the slowest varying index, and write out reflection information on drum.

Procedure

:

The read-in is accomplished by looping over the maximum

number of statistics allowed so that array overflow cannot

occur. An end of file on the input file stops the read-in

before the loop is completed. Reflections with read errors,

zero values of Fo, with the unobserved indicator set to 1,

and with indices all zero are skipped. A count is kept of

the accepted reflections. Subroutine MAXHKL is called to

establish the maximum value of the indices and the slowest

varying index is stored in the ID array for each reflection.

Modification Aids:

This routine is the most likely one to be modified by

users to fit their particular situation. The only part of

read statement which should not be changed is the index

read-in of the form IJ (JJ) , IJ (KK) , IJ (LL) , which allows

the order to be specified. IJ (1) is the slowest varying

index which is stored in the ID array for each reflection.

Important Variables:

ID - array holding the value of the slowest-varying

index for each reflection;

J - total number of usable reflections.

JJ, KK, LL - order of storing the indices, i .

e

.

,

JJ=2, KK=1,

LL=3 stores them in the order k, h, l.

IJMIN, UMAX - minimum and maximum values of the slowest

varying index.

63
Subroutine FREAD

START

loop over storage

NN maximum locations

FIN

read reflection Rag,

indices. Fo 3 O from

input file on unit NTA.PE

t

w print out total

! No. of reflections,

1 input unit No.

(^RETURN)

Ite, increment count ^ call MAXHKL

! store slowest

varying index

replace

Index

replace UMAX

with index

write, flag,

indices, Fo 5 O on drum

64
III. F. Subroutine INFOFC

Purpose

Obtain reflection information from the binary data

file generated by the X-RAY SYSTEM. Calculate the

statistics and store corresponding reflection information

on drum for further reference during the run.

Procedure

Reads physical records from the binary data file until

it finds the first one (logical record 15) containing the

reflection information. Users are advised to consult the

manual of the X-RAY SYSTEM for a complete description of

the binary file. Each record is read into a buffer and

reflection records are unblocked from this buffer. Array

LOGFC contains the drum address for reflection information

corresponding to each statistic- Statistics are stored in

array DM. The FORTRAN listing gives further explanation of

the procedure.

Important Variables:

JY - total number of usable statistics.

JZ - drum address for Fc and T-bar.

DM - array holding the JY statistics.

LOGFC - array holding drum addresses for information

corresponding to the statistics for each reflection.

NN - maximum number of statistics allowed by size of

DM array.

NFTEST- extra flag for miscellaneous purposes.

IPT - unblocking pointer-points to the last word unblocked.

NWORD - number of words in the block.

IJ - array holding the indices.

65
Subroutine INFOFC

66
III. G. Subroutine LSFIT

Purpose

Select the required statistics and express their

variation with whatever is selected as the abscissa (i.e.,

normal quantiles, d* values) in terms of linear,

quadratic and cubic equations.

Procedure

Examine input flags, if abscissa is normal quantiles

leave out the tails, i.e ., those values of statistics for

which the absolute values of the corresponding normal

quantiles are greater than 2.0. Select statistics according to specifications indicated by variables M and II.

Form the required sums for a least-squares fit of linear,

quadratic and cubic equations.

Important Variables:

X - array holding a total of J abscissa values.

Y - array holding a total of J statistics.

ALINE - intercept for the linear equation.

BLINE - slope for the linear equation.

ACURVE, BCURVE , CCURVE - coefficients for quadratic equation.

ACUBIC, BCUBIC , CCUBIC, DCUBIC - coefficients for cubic equation

M - "hundreds" digit specified octant number, 1 through 8.

"units" digit positive if classes or layers wanted.

II - number of the bit in the LOGHKL word to be tested

to divide the reflections into classes.

K - number of statistics in the subset.

NOCT - octant number, 1 through 8.

67
Subroutine LSFIT

63
III. H. Subroutine PLOTEM

Purpose

:

Set up the proper title for the plot, acquire all the

information required for the plot and print out the information actually plotted.

Procedure

:

Examines control variables to determine what type of

plot is wanted. Gets the corresponding title words and

stores them in array ATITLE. Call subroutine LSFIT to

obtain least-squares lines and curves for the statistics or

subset. Calls subroutine FITCHK to calculate the values

predicted by the least-squares lines and curves and establishes

the scatter about the line if d* or log (Fc 3 /sin2 9) plots

are being made. Subroutine PLOTS is called to produce the

actual plot on the line printer after a new page has been

titled properly by PLOTEM. After the plot has been made,

PLOTEM prints out the plot symbols for identification of

the various curves and lists the values that have been

plotted. Subroutine PLOTS sets its own axial limits from

the input data.

Important Variables:

JC - total number of statistics stored in DM array.

X - array holding abscissa values , normal or half-normal

quantiles, d* values, or log CFc

3 /sin2 9) values.

DMSMAL- array holding 50 values from all the statistics or

subset to be plotted.

XSMALL- array holding the corresponding abscissa values.

M - control variable? "hundreds" digit is the octant

number, "tens" digit indicates half-normal plot

(structural parameters, etc.), "units" digit =1

indicates standard classes; =2, h layers; =3, k layers;

=4, H layers.

LAYT - layer type as specified by "units" digit of M.

II - number 0 through 10 for bit standard class; bit

number 11 through 20 for layers 1-10.

J - number of points actually plotted.

69
K - number of points in reflection class or layer.

L - number of points in interval, i.e., L = K/50 + 1;

if K/50 < 1 every point will be plotted.

ALINE, BLINE - coefficients for linear equation.

ACURVE, BCURVE , CCURVE - coefficients for quadratic.

ACUBIC, BCUBIC, CCUBIC, DCUBIC - coefficients for cubic

equation

.

SUMLN, SUMCV, SUMCUB - goodness-of-fit for line, quadratic,

and cubic equations.

NPT - number of curves to plot, 3 or 5

.

NUMBER - number of points in each curve.

ATITLE - array holding the plot title.

70
Subroutine PLOTEM

71
III. I. Subroutine REREAD

Purpose

:

Reads the file of matched reflections (two sets of

observed data) previously generated by subroutine SCALE,

calculates and stores the corresponding statistic using the

value for the scale factor SKALE (one of REREAD input arguments) . Writes the reflection pair on drum file for further

reference and stores the location in catalog word LOGFC

.

Also characterizes the distribution of statistics.

Procedure

The input file is buffered into a buffer of 490 words

which contain information on 7 pairs of matched reflections.

Each pair consists of the indices and two values of Fo and

two values of a (Fo) . These are unblocked until the buffer

is empty, when a new read-in takes place. The end of the

data is marked by a value of 99 for the first index. After

unblocking, counters JY and JZ are incremented, the statistic is calculated, and the indices and associated values

for the pair are written on the drum file. The location of

this information is stored in the LOGFC (I) word corresponding

to the I-th statistic. Upon completion of input the mean

and e.s.d. of the statistics are calculated and printed out

to provide a characterization of their distribution.

Important Variables:

BUFF - array holding 7 pairs of matched reflections.

NBUFF -size of array BUFF.

JY - total number of statistics calculated.

JZ - number of Fo, «r values, used for pointing to

corresponding drum locations.

DM - array holding JY statistics.

LOGFC- array holding locations of information on drum

corresponding to the statistics.

JJ, KK, LL - order in which to unpack the indices for

matched pair.

SKALE -scale factor relating the two observed data sets.

SUMM- mean of the distribution of statistics.

SSUM- standard deviation of distribution of statistics.

72
Subroutine REREAD

print number

of statistics,

input unit file

number

calc. <j>

loop over

statistics

add (S-&?)

i to sum

FIN

i—

calc. ESD for

statistics, print

out mean, ESD

c

RETURN

73
III. J. Subroutine SCALE

Purpose

Match reflections from two sets of observed data,

store drum addresses of the pair on LOGFC word and calculate

the scale factor relating the two observed sets of data.

Optional tasks include applying spherical absorption corrections to the data sets and writing out matched pairs on

unit NR.

Procedure

Reads the packed index word, Fo and a (Fo) values for

each data set from the file written on drum by subroutine

FREAD. This information is read in by batches, a batch

being all those reflections with the same value of the

slowest varying Miller index (batching reflections saves

time in the sort procedure) . Corresponding reflections

are matched by comparing the index words JKLA and JKLB and

matched indices, Fo and a (Fo) values are optionally copied

into a buffer which is then output on unit NR. The buffer

holds information corresponding to 70 pairs of reflections.

The corresponding drum addresses of the pair are stored in

LOGFC (I) for the I-th statistic. If requested, a spherical

absorption correction is applied before pairs are written

out on unit NR. Subroutine ABSORB is called to calculate

the transmission factor and its error for each reflection.

The sum of Fo for each data set is accumulated. Finally,

the scale factor is calculated and printed out for reference

together with the total number of statistics, the number of

reflections in each data set and the input unit number.

Important Variables:

IDA, IDB - arrays holding values of slowest varying index

for the data sets.

JA, JB - number of reflections in data set one and two.

LOGFC - array holding drum address for reflections making

up the statistic.

SKALE - scale factor relating the two data sets.

IJMINA, IJMAXA - minimum and maximum values of slowest

varying index for set one.

JKLA - packed indices for data set one.

74
JKLB - packed indices for data set two.

KA, KB- arrays holding drum addresses of the matched pair.

FA, FB -arrays holding Fo values for the two data sets.

SGA, SGB - arrays holding sigma values for the two data sets.

NBUFF - size of buffer BUFF.

75
Subroutine SCALE

start)

initialize counters

position input

to correct file

loop over batches^

loop over refIections\FIN

in first data set

'

increment count read Fo , sigma, index word from drum, store

in working arrays.

store drum address

loop over reflections^

in 2nd. set

increment count read Fo , sigma, index word from drum, store

in working arrays.

store drum address

of information

<

>p ov

In batch for first set r

loop over reflections

in corresponding batch .second set >

.EEL

FIN

increment count of statistics store drum address in LOGFC

add Fo(l) and Fo(2) to sums for scale factor

70, 80 end of loop over reflections In sets 1 and 2 for present batch <

get another batch>-®

Y

calculate seal sets. K= 21 Fo

e factor relating to two

(set 1)/X Fo (set 2)

print out input unit numbers, scale

factor, no. of statistics, no. of reflections in each set.

\

'

76
IV. FORTRAN Code

IV. A. Introduction

STATUS was written using UNIVAC 1108 FORTRAN V code

and is generally compatible with FORTRAN IV. The section

on Common Blocks and Arrays describes common block storage

and its purpose. Arrays requiring similar dimensions and

the parameter variables which set the dimensions are also

given. A number of equivalences used to conserve storage

are indicated along with the names of the arrays that are

later written over. Storage allotment can be changed by

changing the value of the parameter variables. Users

whose compilers will not allow the parameter statement can

remove it and dimension the arrays explicitly in the main

routine, STATUS and, if required, in the subroutines also.

The section on Features Specific to UNIVAC 1108 FORTRAN V

indicates where these statements and all other features of

the code that are not American Standard FORTRAN are located

in STATUS. The main control parameters are described in the

section on Input Parameters.

IV. B. Features Specific to UNIVAC 1108 FORTRAN V

NTRAN - a UNIVAC routine which handles transfers of information between central memory and mass storage or tape devices.

It performs such functions as reading blocks of data into

arrays, skipping files, closing files, etc.

NTRAN is called by the following subroutines:

DRUMRD Line number 17

DRUMRT " " 23

FLD - a UNIVAC FORTRAN V function for manipulation of

specified bits in a computer word.

Called by the following subroutines:

CATLOG Line number 38 - 47, 61 -63

FITCHK II ll 44

INFING II II 33 - 34

LSFIT II ll 52

PARAMETER statement - Assigns numerical values to variables

specifying array dimensions and enables one to change program

storage by redimensioning arrays conveniently.

Used in the following routines:

77
STATUS Line number

FINGFO

FITCHK

FREAD

INFING

LSFIT

REREAD

SCALE

USER

29

19

10

11

18

12

13

17

4

I/O statements with END = statement number optional in the

following routines:

STATUS Line number 98

FREAD " " 17

Dimension statements including PARAMETER variables:

STATUS Line number 55

Statement involving array ATITLE (variable = 6H text)

STATUS Line number 84

" " 2 54

" " 276

" " 305

" " 37 8

" " 453

" " 47 6

PLOTEM Check entire subroutine

- 88

- 256

- 278

- 307

- 386

- 454

- 477

78
12

3

4

5

6

7

8

9

10

1 1 12

1 3

1 4

15

1 & 17

18

19

2 C 21

22

23

24

25

26

27

28

29

30

3 1 32

33

34

35

36

37

38

39

40

4 1 42

43

44

45

46

47

48

49

5C

51

52

53

54

55

56

57

IV. C. LISTINGS

C STATUS MAIN PROGRAM

C THE PROCEDURE IN THIS PROGRAM IS THAT DESCRIBED IN THE PAPER • NORMAL

C PROBABILITY PLOT ANALYSIS OF ERROR IN MEASURED AND DERIVED QUANTITIES

C STANDARD DEVIATIONS 1 BY S.C. ABRAHAMS AND E.T.KEVE ACTA CRYST.C1971)

C A 27, P, 157-165.

C THIS PROCEDURE IS TO COMPARE (DELTA F)/( SIGMA F) WITH A NORMAL

C DISTRIBUTION. THE PLOT SHOULD IDEALLY BE A STRAIGHT LINE WITH A C SLOPE =1.0.

C

C RESIDUAL ANALYSIS DEVELOPED by L.W .SChROEDER BASED ON N.R. DRAPER

C AND H. SMITH, ‘APPLIED REGRESSION ANALYSIS' (1966)WI|_EY. CHAPTER 3.

C NOTE THAT DELTA F/SIGMA F IS THE RESIDUAL IN 'UNIT DEVIATE* FORM.

C PLOTS OF RESIDUALS VS INDEPENDENT VARIABLES ( D -ST A R » ETC . ) SHOULD NOT

C SHOW ANY TREND, ALSO, THE SCATTER SHOULD BE UNIFORM WHEN CORRECT

C WEIGHTS HAVE BEEN APPLIED.

C

C PROGRAM WRITTEN BY B. DICKENS AND L.W. SCHROEDER (311.05) AT THE

C NATIONAL BUREAU OF STANDARDS, WASHINGTON, D.C. 20234.

C PHONE (301) 926-2455.

C

C USE AT YOUR OWN RISK PROGRAM ONLY TESTED ON UNIVAC 11C8-EXEC-8.

C

c

C ******** *********************************** #*****#*44**#** ********

c

C PARAMETER STATEMENT ALLOWS EASY ADJUSTMENT OF ARRAY SIZES.

C QUANTITIES TO DO WITH AVAILABLE CORE AND ARRAY DIMENSIONS

C CAUTION DON , T MAKE MM MORE THAN 100 WITHOUT ALTERING EQUIVALENCES

PARAMETER N=50 , NNA=1R00, NN=7200, NNB=180 1, NNC=360 1» NND=5401»

2 MM= 60 * MM M= 4 0 0* NBUFF=70C, MM MM=3

2

C I N I T I ALL I ZE MAXIMUM AND MINIMUM VALUES OF MILLER INDICES TO

C RIDICULOUS VALUES

DATA MAXH,MAXK,MAXL,MINH,MINK,MINL /-1CCC, -1000, -1000,10 00 . 1 COO . 2 1000/

COMMON /D/ MAXH , MI NH ,MAXK, MI NK, MAXL . MI NL

COMMON /TITLE/ AT I TL E ( 1 4

)

C I NITIALLl ZE HIGHEST AND LOWEST VALUES OF MILLER INDEX KEPT IN CORE IN

C SAME WAY

DATA I JM 4XA , I JMAXB , I JMI NA , I JMI NB /- 1 0 00 , - 1 0 00 . 1 0 00 t 1 0 C 0/

COMMON /IO/ IN, NOUT ,NDRUM, NT APEA, NTAPEB

C DUMMY IS USED IN EQUIVALENCE STATEMENTS SOME SUBROUTINES TO SAVE

C SPACE

COMMON /h/ DUMMY ( N N

)

C JJA.KKA AND LLA ARE THE ORDER IN WHICH THE MILLER INDICES ARE PACKED

C INTO THE JKL WORD. THEY ARE IN COMMON WITH NR AND NABS BECAUSE

C THESE QUANTITIES ARE ALL USED AT THE END OF THE SCALE SUBROUTINE.

COMMON /ORDER/ J J A , KK A , LLA , NR , NF I LE R , N ABS

DI MENS ION TITLE(12)»X(NN), DMSMAL (N.5).XSMALL(N,5),DM(NN),

2 LOGHKL(NN) »LOGFC( NN

)

DIMENSION DOC ( 2 ) / 6 HJ AN • 6H 1975 /

C THESE ARE THE RECIPROCAL CELL VALUES

CCMMON /G/ ASTAR.BSTAR ,CSTAR , C 0 SA ST , COSB ST , COSGST

COMMON / EX PT/ W A VE , S K AL E, UR ( 2 ) , DU R ( 2 ) , NF

C THESE ARE USED IN SUBROUTINE SCALE

DIMENSION IDA(NN), IDB(NN),FA(NNA),FB(NNA)»SGA(NNA) ,SGB(NNA) , 2 J KLA ( NNA) , J KLE ( NNA ) , KA(NNA) , KB ( NNA

)

C THESE ARE USED FOR THE ATOMIC PARAMETER TESTING.

79
58

59

60

61

62

63

64

65

66

67

68

69

70

71

72

73

74

75

76

77

78

79

80

81

82

83

84

85

86

87

88

89

90

91

92

93

94

95

96

97

98

99

ICO

10 1 10 2

1 03

1 C 4

105

1 C 6

107

1 08

109

1 10

1 1 1

1 1 2

1 13

1 1 4

1 15

D I ME NS I ON NAMEA< MM) , OC A< MM) , XYZA< 3, MM) , S IGX A ( 3, MM ) , BA ( 6, MM ) , 2 S I GBA ( 6 » MM ) * NAMEB ( MM ) , OCB( MM) ,XYZB(3,MM),SIGXB(3,MM),BB(6,MM),

3 SI GBB (6 »M M ) , SI GOA (MM) , S I GOB (MM) « LOG ( NN ) , MATCH ( MM

)

C LOTS OF EQUIVALENCES TO SAVE SPACE. SO THAT LOTS OF REFLECTIONS CAN

C BE TREATED.

EQUIVALENCE ( DUMMY .LOGHKL

)

EQUIVALENCE (DM, IDA), (X.IDB), ( LOGHKL ( 1 ) » J KL A ( 1 ) ) , ( LOGHKL ( NNB ) , 2 JKLE? ( 1 ) ) » ( LOGHKL (NNC ) ,KA( 1 )) , ( LO GHKL (NND ),KB( 1) ) EQUIVALENCE (LOGHKLC2CO ) ,OCA(l )

)

EQUIVALENCE

EQUI VALENCE

EQUIVALENCE

EQUIVALENCE

EQUI VALENCE

EQUIVALENCE

EQUI VALENCE

EQUI VALENCE

EQUIVALENCE

EQUIVALENCE

EQUI VALENCE

EQUIVALENCE

EQUI VALENCE

EQUIVALENCE

C I

1 0

20

EQUIVALENCE

N I T I ALL I SE SOME

DO 20 1=1, 14

AT I TLE ( I ) = 6H

AT I TLE( 8 )=6H PLOT

ATITLE(9 ) = 6HOF DM(

AT I TLE ( 1 0 ) =6HVERT

)

AT I TLE ( 1 1

) = 6H AGAIN

ASSIGN UNIT NUMBERS

CALL PRIME

WRITE ( NOUT , 370 ) DOC

( LOGHKL( 300 ) ,OCB(l )

)

(LOGHKL(4CO ),SIGOA( 1 )

)

(LOGHKL (500) ,SIGOB(l )

)

( LOGHKL ( 600 ) , XYZA( 1,1))

( LOGHKL (9C0 ),XYZB(1,1))

( LOGHKL ( 1200 ) , SIGXAC 1,1))

( LOGHKL ( 1500 ), SIGXB( 1, 1 )

)

( LOGHKL ( 1 80 0 ) , B A ( l , 1 ) ) (LOGHKL( 2400) , 8B( 1 , 1 )

)

(LOGHKL( 3000 ), S IGBA( 1 , 1 )

)

( LOGHKL ( 36C 0 ) ,SIGBB( 1 , 1 )

)

( LOGHKL ( 4200 ) , NAMEA( 1 )

)

( LOGHKL ( 4300 ) , NAMEB ( 1 )

)

( LOGHKL ( 4400) , MATCH( 1 )

)

(LOGFCC 1 ) ,LOG( 1 ) ) OF title of plots

NDRUM

c

c

c

cc30

c

cc

***$$$$* $$4:* 4c4c *4 * * 4c 4c 4* 4c4c 4c 4c 4c 4c 4c 4: 4c 4c 4c 4c 4c 4c 4c 4c 4c 4c 4c ********* ** 4c4c4c4c4c4c4c4c

READ AND WRITE TITLE

CONTINUE

READ ( IN.460 ,END=340 ) TITLE

WRITE (NOUT, 470) TITLE

READ IN UNIT CELL

CALL CELL

READ IN REFLECTION TAPE NUMBER AnD ORDER IN WHICH INDICES VARY FOR

TAPE A AND TAPE B.ALSO READ FLAGS FOR READ REFLECTIONS AND READ

PARAMETERS, ETC.

READ ( IN, 480 ) CARD , NTAPE A, J J A , KKA , LLA , NT A PE B , NF TEST , NP TEST , NF ILE A

,

2NFI LEB, NPRI NT, NO , NOC , NL AYO , NOCTO , ND , NDC, NL A YD , NOCTD , NQ , NQC , NL A YQ ,

N

30CTQ.NR , NABS, MI NDM , NF

WRITE (NOUT, 380) NT APE A , NF ILEA.NTAPEB, NF ILEB,JJA,KKA, LLA , NFTE ST , NP

2TE ST , NPR I NT

IF (CARD.NE.3HCTL) GO TO 330

READ 4TH CARD ONLY IF NECESSARY, I.E. FOR SCALING, Q-PLOTS, OR

IF ABSORPTION CORRECTION WANTED.

IF(NQ.LE.O . AND .NABS.LE.O.AND.NFTEST .NE.4) GO TO 35

READ (IN, 390 ) W AVE , SK ALE , UR , DUR

35 CONTINUE

80
116

117

118

1 19

120

121

122

123

124

125

1 26

127

1 28

1 29

130

131

132

133

1 34

135

136

137

138

1 39

140

141

142

143

144

145

146

147

148

149

1 50

151

152

153

154

155

156

157

158

159

160

161

162

163

164

165

166

167

168

169

170

171

172

173

C

C END OF CARD READ-IN, PRINT HEADINGS

IF (NR.GT.O) WRITE (N0ut,400) NR

IF (NABS.GT.O) WRITE (NOUT.410)

IF (MINDM.GT.O) WRITE (NOUT.420)

C

C CHECK INPUT PARAMETERS FOR ERRORS - STOP IF NECESSARY.

IF( JJA.LT.O .OR. KKA.LT. 0 .OR.LLA.LT. 0 ) GO TO 330

IF ( J J A . GT • 3 • OR . KK A .GT . 3.0R.LLA.GT • 3) GO TO 330

IF ( NFTEST.LT. 0. OR. NFTEST. GT. 5) GO TO 330

IF ( NPTEST .LT .0 .OR .NPTEST. GT .2 ) GO TO 330

IF ( N0C.GT.7.0R.NDC.GT.7.0R.NQC.GT.7) GO TO 330

IF (NABS .LT.O. OR.NR.LT. 0) GO TO 330

IF (WAVE. LT.-.000001 .OR. WAVE. GT. 3.0 ) GO TO 33C

IF ( SKALE.LT.-.O 1 ) GO TO 330

IF (NLAYO.LT.O .OR.NLAYO.GT.1C) GO TO 330

IF ( NOCTO.LT.O.OR.NOCTO.GT.7 ) GO TO 330

IF ( NOCTO.EQ • 5 .OR • NOCTO .EQ .6 ) N0CT0=7

IF (NLAYD.LT.O .OR.NLAYD.GT.10> GO TO 33^

IF ( NOCTD.LT.O .OR. NOCTD. GT .7 ) GO To 330

IF (NOCTD. EQ. 5. OR. NOCTD .EQ. 6) N0CTD=7

IF ( NLAYQ.LT. 0. OR. NLAYQ.GT. 10 ) GO TO 330

IF (NOCTQ.LT.O.OR.NOCTQ.GT.T ) GO TO 330

IF ( NOCT Q. EQ • 5 • OR . NOCT Q • EQ .6 ) N0CTQ=7

IF(NFLAG.LT.-Ol) GO TO 330

C

C PRINT OUT WHAT THE PROGRAM WILL ATTEMPT TO DO

IF (NO.EQ. 1) WRITE (N0UT.430) NOC , NLAYO, NOC TO

IF (ND.EQ.l) WRITE (NOUT.440) NDC , NLAYD , NOC TD

IF (NQ.EQ.l) WRITE (NOUT,450> NQC, NLAYQ, NOC TQ

C COMPUTE THE NUMBER OF PLOTS PROGRAM WILL ATTEMPT TO GENERATE

NOPLT S=0

ELE VO= 1

1

ELEVD= 1 1 ELEV Q= 1

1

IF ( N0C.EQ.2.0R.N0C.EQ.3.0R.N0C.EQ.4) ELEVQ = 0 IF (NOC.LE.l) NLAYO=0

IF (N0C.EQ.2.OR. NDC • EQ .3.0R.NDC.EQ.4) ELEVD=0

IF (NDC.LE.l) NL AY D=0

IF ( NQC. EQ .2 .OR .NQC.EQ .3. OR .NQC. EQ .4) ELEVQ=0

IF ( NOC • LE. 0 ) ELEVO=0

IF (NDC .LE.O )ELEVD=0

IF ( NQC »L E • 0 ) ELEVQ= 0

IF (NQC.LE.l) NL AY Q=0

IF (NO.EQ.l) NOPLT S=NOPLTS+ l + ( ELEVO+NLAYO ) * ( NOC TO+ 1

)

IF (ND.EQ.l) NOPLTS=NOPLTS+l+(ELEVD+NLAYD)*(NOCTD+l

)

IF (NQ.EQ.l) NOPLTS=NOPLTS+ l+( ELEVQ+NLAYQ) * ( NOCTQ+1

)

IF (NPTEST. GT • 0 ) NOPLTS=NOPLTS + 3

NPAGES=2*NOPLTS + 5

WRITE ( NOUT , 350 ) N OPLT S , NP AGE

S

IF ( NPAGES.GT. 250 ) WRITE (NOUT.360)

IF ( NPAGES. GT. 250) STOP

C ABSORPTION CORRECTION VALID FOR MU R=0.0 TO 1.0

DO 40 JC = 1 ,

2

IF ( UR( JC) .LT.-.0001 . OR.UR( JC) .GT. 1 .0 ) GO TO 330

40 IF (DUR( JC) .GT.UR( JC) .OR.DUR(JC) .LT.-.0001 ) GO TO 330

C GET READY FOR ASSIGNED GC TO STATEMENTS LATER

NFTEST=NFTEST+

1

81
174

175

176

177

178

179

180

181

182

183

184

1 85

1 86

187

188

1 89

190

191

192

193

194

195

196

197

198

199

200

20 1 202

203

204

205

206

20 7

208

209

210

21 1 212

213

214

215

216

217

218

2 S 9

220

221

222

223

224

225

226

2 27

228

229

230

231

NPTE ST=NPTEST+

1

C POSITION INPUT TAPES NT APE A AND NTAPEB TO NFILEA AND NFILEB

IF ( NT APEA.GT .6) CALL POSITN ( NT APE A, NF ILEA ) jc=o

JZ=-1

JA = 0

JB= 0

C IF FO AND FC BEING READ, COPY OVER SORT ORDER,

JJB= JJA

KKB=KKA

LLB = LL A C SET FLAG FOR CATALOGING-lATER TURNED OFF IF NO CLASSES NEEDED.

NCV=1

C IN PROGRAM, _ NFTEST— 1 NO F,S

C NFTEST = 2 READ FOBS AND FC FROM FINGER FOURIER TAPE

C NFTEST = 3 READ ONE SET OF FOBS FROM FIRST FINGER FOURIER TAPE AND

C SET FROM SECOND FINGER FOURIER TAPE

C NFTEST = 4 F, S FROM XRAY TYPE INPUT, 2 DIFFERENT TAPES,

2

C SETS OF F , S » C NFTEST=5 READ FO 9 S FROM TAPE PREVIOUSLY WRITTEN BY SUBROUT I NESC

A

C

C NFTEST=6 READ Fo AND FC FROM X-RAY 70 BINARY DATA FILE OP USER

C FILE IF INFOFC SUITABLY MODIFIED.

GO TO C 1 40 ,50 , 50,60, 80, 90), NFTEST

C READ F,S FROM FINGER FOURIER TAPE

50 CALL FINGFO ( J A, LO GFC , DM , J J A , KKA , LLA , NFTE ST , I D A , I JMI NA , I JM A X A , NT AP

2 EA » J Z

)

C ADD THESE REFLECTIONS TO TOTAL

JC=JC+JA

C IF WORKING WITH FO AND FC THEY NEED NOT BE SORTED AND SCALED. JUMP

C TO SORTING OF DM VALUES, WHICH WERE CALCULATED IN FINGFO FOR THIS

C CASE.

IF C NFTEST. EQ®2 J GO TO 100

C IF ANOTHER FO SET IS NEEDED, READ IT NOW

IF ( NT AP EB , GT • 6 ) CALL POSITN C NT APEB, NF SLEB

I

CALL FfNGFO I J B , LOGF C * D M , J JB » KKB , LL B, NFT ES T ,IDB,IJMINB, I JMAXB, NT AP

2EB, JZ)

C ADD THESE REFLECTIONS TO TOTAL

JC=JC+J8

GO TO 70

C READ FIRST REFLECTION SET FROM XRAY67 TYPE INPUT

60 CALL FREAD I NT AREA , IDA, JA, JJA, KKA, LLA, IJMINA, IJMAXA)

C READ IN SECOND REFLECTION SET

IF ( NTAPEB. GT.6 ) CALL POSITN { NT APEB, NF ILEB

)

CALL FREAD ( NT APEB , I DB , JB, J J B, K KB , LLB , I JM I NB , I J MAXB

)

C FIND COMMON REFLECTIONS AND SCALE FACTOR

70 CONTINUE

CALL SCALE C I DA « J A , I DB , JB , JC ,LOGFC, SKALE ,TJMINA,IJMAXA,FA,FB»SGA,S

2 GB, JKLA , JKLB « K A » KB

)

C CALL SPEC FOR MINIMUM OF SUM«DM**2> IF REQUESTED

IF CMINDM.GT.O) CALL SPEC C JC , SK AL E , DM , LOGF C , NN

)

C CALCULATE DM QUANTITIES FROM F, S AND SIGMAS

IF IMINDM.LE.O) CALL DMCALC I JC , SK ALE , DM , LO GFC , NN

)

GO TO ICO

C

80 CALL REREAD C JC , LOGFC , D M , J J A , KK A , LL A , NTAPE A , J Z , SK ALE

)

C CALL SPEC FOR MINIMUM OF SUM { DM442 ) IF REQUESTED

IF (MINDM.GT.O) CALL SPEC ( JC , SK ALE , DM , LOGF C , NN

)

82
23 2 233

234

235

236

237

238

239

240

241

242

243

244

245

246

247

248

249

25C

251

252

25 3

254

255

256

257

258

259

26 0 261

262

263

264

26 5

266

267

268

269

270

271

27 2 273

27 4

275

276

277

278

279

280

281

282

28 3

284

285

286

287

28 8

289

GO TO ICO

C

90 CALL INFOFC ( NT Ap E A , JC , JZ', NFTEST , I J M I N A , I JM AXA , LOGFC , D M , NN ) C

C

c ******** ********** ************************* ***********************

c

C SORT DM VALUES In ORDER OF INCREASING MAGNITUDE

100 CALL SORT ( DM , X . LO GFC , LOGHKL • JC , NN

)

C CALCULATE NORMAL PROBABLITY QUANTILES IF NEEDED.

IF ( NPR I NT • GT • O • OR .NO.GE.l ) CALL PPOB (DM,X,JC,NN)

C PRINT DM VALUES AS REQUESTED VIA NPRINT

IF (NPRINT. GT.0) CALL PRYNT ( D M . LOGFC , LOGHKL , JC , J J A , KK A . LL A , J JB , KK

2B.LLB, NPRI NT, X, NFTEST, NN

)

C

C

c ******************************************************************

c

C THIS SECTION FOR PLOTS OF RESIDUALS. VS NORMAL QUANTITLES

IF (NO.NE.l) GO TO 110

M=C

C SET UP REMAINDER OF TITLE.

ATITLE( 12)=6HST NOR

ATITLE( 1 3) =6HMAL QU

AT ITLE( 14) =6 HA NT IL E WRITE (NOUT.550)

C CATALOG REFLECTIONS IF REQUIRED.

IF (NOC.GT.O) CALL DMDQ ( D M , LOGF C , J C . J J A , KK A , LL A » X , LOGHKL , 0 , 0 , NOC

,

2NN ) C PLOT DM VS THE NORMAL QUANTILES.

CALL PLOTEM ( J C» X » DM. DMSMAL. XSMALL »M» I I • LOGHKL , NN » N

)

IF (NOC.LE.C) GO TO 110

C PLOT REFLECTION CLASSES

M=NOC* 1C C0+NLAYO*1 0 + NO CTO

CALL ANISO ( JC, LOGHKL, DMSMAL, XSMALL, DM, X,M, NN,N)

lie continue

c

c

c ** **************** ********** **************************************

c

C SKIP OVER PLOTTING OF DM AGAINST D^StAR IF THESE PLOTS NOT REQUESTED.

IF (ND.NE.l) GO TO 120

M=0

C FIX UP TITLE

ATITLEC 1 2)=6HST

ATI TLE( 13) =6HDST AR

ATITLE( 1 4 ) =6H

WRITE (NOUT.550)

C OBTAIN DSTAR VALUES AND CATALOG IF REQUIRED.

CALL DMDQ ( D M , LOGF C , JC , J J A , K K A , LL A , X , L OGHKL . 1 . 0 , NDC , NN

)

C DO OVERALL PLOT

C SORT X ARRAY HOLDING DSTAR**2 VALUES AND REARRANGE CORRESPONDING ARRAY

CALL SORT (X .DM, LOGFC, LOGHKL , JC, NN)

CALL PLOTEM ( JC , X , DM , DM SM AL , X SMALL , M , I I , LOG HKL , NN , N

)

C SKIP PLOTTING OF REFLECTION CLASSES IF SO REQUESTED

IF (NDC.LE.C) GO TO 120

M=NDC* 1 C 00 + NLAYD*l O + NOCTD

C PLOT DM AS FUNCTION OF H , K , L OF REFLECTION. THE NECESSARY INFORMATION

83
290

291

29 2 29 3

29 4

295

296

297

29 8

299

3C 0 30 1 3C2

3C 3

30 4

30 5

306

307

30 8

309

310

31 1 312

31 3

31 4

31 5

316

317

318

319

320

321

322

323

324

325

326

327

328

329

33 r 331

332

333

334

335

336

337

338

339

340

341

34 2

343

344

345

346

347

C HAS ALREADY BEEN STORED IN THE LOGHKL ARRAY IN THE DMD SUBROUTINE.

CALL ANISO ( JC. LOGHKL. DMSMAL . XSMALL , DM, X, M, NN,N

5

120 CONTINUE

C

C

C

C ****   ** ********** ** ******** ******** **** * **** *******************

c

IF C NQ.NE. 1) GO TO 130

C THIS SECTION FOR PLOTS OF RESIDUALS VS . CRYSTAL Q ( F**2/SI N ( 2THET A )

5

C GET Q VALUES NOW

M=0

C CALCULATE Q VALUES AND STORE IN X ARRAY

CALL DMDQ ( DM ,L0GF C • JC , JJ A , K KA , LL A , X , L OGHKL » 0 t 1 , NQC , NN

)

C FIX UP GENERAL TITLE

ATI TLE ( 1 25 =6HST LOG

ATITLEI 135= 6HF * 4 2/

S

AT I TLE « 1 4)=6HIN( 2T

)

WRITE ( NOUT . 550 5 C DO OVERALL PLOT FIRST.

C SORT X ARRAY HOLDING Q VALUES.

CALL SORT (X ,DM,LOGFC,LOGHKL, JC, NN)

CALL PLOTEM C JC , X , DM , DMS MAL , XSMALL , M , I I , LO G HKL , NN , N

)

C PLOTS FOR CLASSES IF DESIRED.

IF <NQC.LE.05 GO TO 130

C PLOT RESIDUAL VS Q FOR VARIOUS CLASSES OF REFLECTIONS.

M=NQC* If 00-8-NLAYQ*! O + NOCTQ

CALL ANISO < JC.LOGHKL, DMSMAL, XSMALL, DM, X.M, NN,N)

13C CONTINUE

C

cc ****************************** ************* ***********************

c

C TREAT ATOMIC PARAMETERS IN HALF NORMAL PLOTS

C

C

C NPTEST = 1 no parameters

C NPTEST = 2 PARAMETERS FROM FINGER FOURIER TAPE

C NPTEST = 3 PARAMETERS READ IN.

C

140 GO TO <30,160,1505, NPTEST

C SPECIAL USER INPUT TO BE SPECIFIED BY USER.

150 CALL USER < N A ME A , OC A , S I GOA , X YZ A , SI GX A, B A , S I GB A , L A , NTA PE A 5 CALL USER < NAMES, OCB, SI GOB, XYZB, SIGXB, BB.SI GBB . LB » NT APEB

)

GO TO 180

C READ PARAMETERS FROM END OF FINGER FOURIER TAPE

C SKIP OVER REFLECTION SET 1 IF NECESSARY

16C IF < NFTEST.EQ .35 GO TO 170

CALL FINGFO (JA.LOGFC, DM, JJA»KKA,LLA, NFTEST , IDA, I JMINA, IJMAXA ,NTAP

2EA, JZ)

C SKIP OVER REFLECTION SET 2 IF NECESSARY (HAVE TO DO SAME AS FOR

C SET 15

CALL FINGFO < J B, LO GFC , DM , J JB , KKB , LLB , NFTE ST , I DB , I JMI NB , I JMAXB , NT AP

2EB, JZ)

170 CONTINUE

C READ ATOM SET 1.

CALL INFING < N AME A , OCA • S IGO A , X YZ A , S I GX A , B A , SI GB A , L A , NTA PE A

)

C READ ATOM SET 2 84
348

34 9

350

35 1 352

353

354

355

356

357

358

359

360

361

362

363

364

365

366

367

368

369

370

371

37 2 373

374

37 5

376

377

378

379

380

381

382

383

384

385

386

387

388

389

390

391

39 2 393

394

395

396

397

398

399

40 0 401

40 2 403

404

405

CALL INFING ( N A ME3 , OC B * S IGOB , X Y ZB . S I GX B * BB * S I GBB . LB , NT APEB

)

C

I cc ******************************************************************

c

C PRELIMINARY CHECK ON ATOM SETS-EQUAL NUMBERS OF ATOMS EXPECTED AT

C FIRST SIGHT.

180 IF (LA.NE.LB) WRITE (NOUT.490) LA, LB

C CHECK THERE ARE AT LEAST TWO SIMILAR ATOMS IN THE 2 SETS

C FIND CORRESPONDING ATOMS IN THE 2 ATOM SETS

CALL MTCHEK ( N A ME A , L A, N AMEB , LB , M AT CH, I P , MM

)

C TEST FOR MORE DATA IF LESS THAN 3 ATOMS IN COMMON BETWEEN ATOMIC SETS

IF ( IP ,LE. 2) GO TO 10

C THERE MUST BE AT LEAST TWO EQUVALENT ATOMS IF WE GET THIS FAR

C CALCULATE DP VALUES AND ASSOCIATED QUANTITIES FOR OCCUPANCY , PO S IT IONAL

C PARAMETERS AND THERMAL PARAMETERS SEPARATELY.

C CHECK WHETHER STANDARD DEVIATIONS ARE GREATER THAN ZERO TO SEE WHICH

C OF THESE QUANTITIES HAVE BEEN VARIED.

Cc

c ******** ********** ** ******** *************** ************* **********

c

C HALF-NORMAL PROBABILITY PLOTS DONE IN FOLLOWING SECTION.

C DO OCCUPANCIES HERE.

NWH AT=0

C TELLS SUBROUTINE PLOTEM HOW TO MAKE REST OF TITLE FOR THIS PART

C NO SPECIAL REFLECTION CLASSES

M=0

11 = 0

C FIX UP TITLE

ATI TLE( 1 )=6H DM QU

AT I TLE ( 2 ) = 6H AnT I TI

ATITLE(3)=6HES BAS

ATITLEI 4 )=6HED ON

AT I TLE ( 5 ) =6H0CCUPA

ATI TLE ( 6) =6HNCI ES

ATITLE( 1 2)=6HST 1/2

AT I TL E ( 13) = 6HNOR MA L

AT I TLE ( 14) =6hQU A NT

.

C TEST FOR ZERO IN FLOATING POINT.

TEST=.0C0C01

WRITE ( NOUT , 5C 0 ) C CHECK STANDARD DEVIATIONS OF OCCUPANCY TO SEE IF THEY ARE NON ZERO

JC=0

DO 190 1=1, MMM

C LOG WILL TELL PRYNTT SUBROUTINE WHICH PARAMETERS ON WHICH ATOMS HAVE

C WHICH DM VALUES

19C LOG ( I

) =0

C TELL PRYNTT SUBROUTINE THESE ARE OCCUPANCIES

KI ND=

1

DO 200 1=1, IP

C GET THE TWO ATOMIC CATALOGUE NUMBERS.

J=MATC H( I ) /ICO

K=M AT CH ( I)-J*1C0

C TEST FOR NON-ZERO SIGMAS

IF( S IGOA ( J ) .LT . TEST. OR . SIGOB ( K) . LT . TEST ) GO TO 200

C CALCULATE DM VALUE IF APPROPRIATE

CALL DPCALC ( D , OCA ( J ) , OC 8( K) , S I GOA ( J ) , S I GO B ( K ) )

85
406

407

4C8

409

410

41 1 412

413

414

415

416

417

4 18

419

420

421

42 2

423

424

425

426

427

428

429

430

431

432

433

434

435

436

437

438

439

440

441

442

443

444

445

446

447

448

449

450

451

452

45 3

454

455

456

457

458

459

460

46 1 462

463

JC=JC+

1

C STORE IT AWAY

DM(JC) =D

C STORE WHERE IT CAME FROM

LOGC JC )=MATCH( I ) *100+1

2 00 CONTINUE

C CONTINUE PROCESS IF THERE ARE ENOUGH PARAMETERS FOUND

210 IF (JC.GE.2) GO TO 220

WRITE ( NOUT ,510) JC

C GO TO NEXT STAGE

GO TO 240

C ARRAYS X AND ID B ARE DUMMYS AT PRESENT.

220 CALL SORT ( DM » X ,LO G, IDB , JC . NN

)

C CALCULATE HALF NORMAL PROBABILITY DISTRIBUTION FOR THESE PARAMETERS.

C GET APPROXIMATION TO MEDIAN OF I-TH ORDER STATISTIC

CALL UNIMED (JC«X)

C

C GET EXPECTED VALUES VIA PERCENTAGE POINTS AND MEDIAN.

C REFERENCE DANIEL , TECH NUME TRI CS , 19 59 , P A GES 311-341

DO 230 LS=1 , JC

Q=X( LS

)

Q=C1 . 0—Q ) /2

.

X(LS)=PINV(Q)

230 CONTINUE

C PRINT OUT ALL DM VALUES FOR ATOMIC PARAMETERS

NPR I NT= 1 IF ( NPR I NT . GT • 0 ) CALL PRYNTT ( JC , D M , X , LOG, NPR I NT , K I N D, OC A, OCB , S I GO

2A» SI GOB, XYZ A »XYZB» SI GX A»SIGXB»BA»BB,SI GB A , S IGBB, NAME A, NAMES, NN, MM

)

C PLOT DM VALUES AGAINST X VALUES, ETC.

WRITE C NOUT « 550

)

CALL PLOTEMI JC,X .DM, DMS MAL , X SMALL , M , I I , LOGHKL. NN, N

)

C GO TO NEXT STAGE

240 NWHAT=NWHAT+1

C

C RE- IN IT I ALL I ZE

JC=A DO 250 1=1,400

250 LOG ( I ) =0

C OFF WE GO

C

GO TO (260,290,3201, NWHAT

C DO POSITIONAL PARAMETERS HERE

C

260 WRITE (NOUT, 520)

C TELLS PRYNTT SUBROUTINE WE ARE NOW DOING POSITIONAL PARAMETERS

KI ND=2

C FIX UP TITLE

AT I TLE (5 ) = 6HX,Y,Z'

ATI TLE (6 ) =6 HS

C GET DM® S IF APPROPRIATE (NON-ZERO SIGMAS OF PARAMETERS) AND PLOT THEM

DO 280 1=1 , IP

J=MATCH ( I ) /I 00

K=M AT CH ( I S- J*10C

DO 270 L=1 .

3

IF ( SIGXA( L. J) .LT. TEST. OR. SI GXB( L, K) .LT .TEST ) GO TO 270

CALL DPCALC (D.XYZACL, J ) ,XYZB(L,K ) , SIGXA(L , J) , SIGXB( L,K)

)

JC=JC+

1

DM( JC ) =D

86
464

465

466

467

468

469

470

471

472

473

474

475

476

477

478

479

480

481

482

483

484

485

486

487

48 8

489

490

491

492

49 3

494

495

496

497

498

499

500

50 1 502

50 3

5C4

505

5C6

5 07

50 8

509

510

511

512

513

514

515

516

517

518

519

520

521

LOG ( JC) = MATCH! I ) *1 00+L+1

270 CONTINUE

280 CONTINUE

GO TO 210

C

C DO THERMAL PARAMETERS HERE! ASSUMED ANISOTROPIC)

c290 WRITE ! NOUT * 5 30

)

C THERMAL PARAMETERS NOW, 0 SUBROUTINE PR YNTT ( PR YNTTE , 3RD DECLENSION

C (CAN'T WASTE 6 YEARS OF LATIN, GOT TO WORK IT IN SOMEWHERE))

KIND=3

C FIX UP TITLE

ATITLE! 5 ) = 6HBETAS

AT ITLE! 6 ) = 6H

C AGAIN, GET DM VALUES, PLOT THEM.

DO 310 1=1, IP

J=MATCH! I ) / 1 C 0

K=MATCH( I)-J*10C

DO 300 L=1.6

IF ( SIGBA(L,J) .LT. TEST .OR. SIGBBIL, K) .LT. TEST) GO TO 30-' CALL DPCALC ( D , B A ( L , J ) , BB ( L , K ) , S I GB A ( L » J ) , S IG BB ( L , K ) ) JC=JC+1

DM( JC) =D

LOG ( JC ) = MATCH( I ) 71 00+L+4

300 CONTINUE

310 CONTINUE

GO TO 210

C FINISHED THIS PART

320 WRITE (NOUT, 540)

C CHECK FOR MORE DATA, BEGINNING WITH TITLE.

GO TO 30

Cc

c t****************************************** ***********************

c330 CONTINUE

C COMES HERE WHEN INPUT PARAMETERS FAULTY.

WRITE (NOUT, 560)

340 STOP

C

c

c

c

c350 FORMAT !1H0,I5,15H PLOTS POSS I BLE , I 5 , 1 5H PAGES OF PLOTS)

360 FORMAT (1H0.49H TOO MUCH OUTPUT POSSIBLE-EXECUTION STOPP ED-TH I NK

)

370 FORMAT (1H1,98HLIST SUBROUTINE CHANGE TO KEEP UP TO DATE ON THE Ev

20LUTI0N OF THIS PROGRAM. THIS IS THE VERSION OF , 2 A 6 , // 9H ** ***UN I T

,

314, 52HMUST BE ASSIGNED FOR INTERMEDIATE SCRATCH F I LE . 77 777 //

)

38C FORMAT ( 1HC.28HFIRST INPUT DATA SET ON UNIT, 13, 6H FILE , I5.30H SEC

20ND INPUT DATA SET ON UNIT, 13, 6H FILE , I 5 / 1 HO , 1 OH SORT ORDER, 312, 8H

3 NFTEST=,I2,8H NPT EST= , I 2 , 1 3 H PRINT OPTION, 12)

390 FORMAT! 2F1C. 7 , 1 X ,4F5. 2)

400 FORMAT ! 1 HO , 34 H MAT CHED REFLECTIONS OUTPUT ON UNIT, 14)

410 FORMAT ! 1H0, 51 HABSORPT I ON CORRECTIONS WILL BE APPLIED TO DATA SETS

2 ) 4 20 FORMAT ( 1H0 , 66HSUBP0UTI NE SPEC WILL BE CALLED TO VARY SCALE TO MIN

2IMIZE SUM DM772)

87
522 4 30 FORMAT ( 1H0 » 47H STATISTICS VS NORMAL QUANTILES WILL BE PLO TTED/5 H

523 2N0C= * I

2

,7H NLA YO=, 13 , 1 5H OCTANT OPTION=,I3)

524 4 40 FORMAT ( 1H0,62H STATISTICS VS D- ST AR ( 2 S I N ( THET A ) / W AV ELENGTH ) WILL

525 2 be plotted /sh ndc=,I2,7h nlayd=,

i

3 , ish octant option=,I3)

526 4 50 FORMAT < 1 H 0 » 57H STATISTICS VS CRYST . Q< F** 2/ S I N ( 2*THET A ) ) WILL BE P

527 2LOTTED/5H NQC=*I2.7H NL A YQ= , I 3 , 1 5H OCTANT OPTTON=,I3)

528 4 60 FORMAT ( 1 2 A 6 ) 529 4 70 FORMAT l 1H1,5X,12A6)

530 480 FORMAT C A3, 12,815, 13,1612)

531 490 FORMAT ( 1 4H *** -WAPNI NG-, 16 , 1 5H ATOMS IN SET 1,I6,9H IN SET 2)

532 500 FORMAT ( 60 H 1 TEST OCCUPANCIES IF THEY HAVE NON-ZERO STANDARD DEV I AT

53 3 2IONNS)

534 5 10 FORMAT ( 64 H LESS THAN 2 COMMON ATOMIC PARAMETERS VARIED, GO ON TO N

535 2 EXT STAGED

536 5 20 FORMAT (27H1TEST POSITIONAL PARAMETERS)

537 5 30 FORMAT (24H1TEST THERMAL PARAMETERS)

538 540 FORM AT ( 23H *** JOB COMPLETED *** ) 539 550 FORMAT ( 1 HI D 540 560 FORMAT 1 1H0 » 57HCQNTR OL PARAMETER I N CORR ECT- CHECK INPUT-EXECUTION S

541 2 TOPPED

)

54 2 C

54 3 END

88
123

4

5

6

7

8

9

1C

1 1

1 2

1 3

1 4

15

1 6

17

IB

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

48

49

50

5 1 52

53

54

55

56

57

SUBROUTINE ABSORB ( J , K , L * A A * AB . D A A , O AB

)

C

c *#**:£#***###***#:<£ 4 4*4 ** 4*** *********** ******* ****** **********

cC SUBROUTINE ABSORB CALCULATES THE ABSORPTION CORRECTIONS TO BE APPLIED

C TO THE F VALUES IF REQURED. IT IS SET UP FOR SPHERICAL CRYSTALS AT

C PRESENT.

C RECIPROCAL CELL CONSTANTS

COMMON /G/ ASTAR,BSTAR,CSTAR,COSAST,COSBST, COSGST

C 10 UNITS

COMMON /IO/ IN, NOUT.NDRUM. NTAPE A , NTAPEB

COMMON /EXPT/ W A VE , SK AL E , U R ( 2 ) , DUR ( 2 ) , NF

DIMENSION A( 2) « DAD UR ( 2) «ANORM( 2)

C NF= FLAG FOR F OR F**2, 0 FOR F, 1 FOR F**2

C SPHERICAL ABSORPTION CORRECTIONS FOR 1) CRYSTAL A 2) CRYSTAL B C DEGREE TO RADIAN CONVERSION PARAMETER

RAD= 3 . 1416/180.

IF (NTIME.EQ.l) GO TO 40

DO 10 L= 1 i

2

IF ( DUR(L) .LT. ( .02*UR(L) ) ) DUR ( L ) = . 0 2*UR ( L

)

C NORMALIZE THE CORRECTION TO FIRST DATA SET SO SCALE DOES NOT CHANGE

C VERY MUCH. COMPUTE FACTOR FOR 2-THETA =0 CASE.

C ABSORPTION CORRECTION EQUATION FROM R OUS E , C 00 PER, ET AL . ACTA CRYST.

C ( 1970) A 26 , P. 682-691 . A( L) =EXP (-1 . 5108 *UR (L )+ . 0951 *UR (L ) **2)

10 CONTINUE

ANORM ( 1 ) = 1 .0

ANORM ( 2 ) =A( 1 ) / A( 2)

WRITE ( NOUT , 30 ) NF , WA VE , UR , D UR , A NOR

M

C SET FLAG TO GO TO STATEMENT 40 ON THE NEXT CALL.

NTI ME =

1

40 CONTINUE

FJ = J FK=K

FL = L

C CALCULATE D SPACING AND THEN 2THETA FOR J,K,L REFLECTION

TEMPT H= ( FJ4ASTAR ) * 4 2 + ( FK *B ST AR ) * * 2+ ( FL*CS T A R ) 4 * 2+ 2 . 7FK7FL* BST AR*CS

2TAR*COSAST+2 .*FL*FJ*AST AR *CST AR *CO S BST + 2 . *F J*FK * AST AR *B STA R*C 0 SGST

D=1 . /SQRT( TEMPTH)

TEMPTH= ( ASIN( WAvE/( 2 .*D ) ) )/RAD*2.

TH2=TEMPTH

C CALCULATE ABSORPTION CORRECTION VIA ROUSE EQUATION.

FACT 1 = 1 . 510 8- .031 5*SIN< TH2 )**2

FACT 2=-. 0951 -.2898 *SIN( TH2 )**2

C CALCULATE TRANSMISSION FACTOR FOR TWO DATA SETS.

DO 20 L= 1,2

A ( L ) =E XP ( — 1 . 0*FACT1*UR(L)-1 . 0*F ACT2*UR ( L ) ** 2 ) *A NOR M< L

)

DADUR(L)=A(L)*(-1. 04FACT1-2. 0*F ACT2*UR ( L )

)

20 CONTINUE

A A = A ( 1 ) A B= A ( 2)

IF (NF.NE.l) AA=SQRT ( A A

)

IF (NF.NE.l) AB=SQRT(AA)

DAA=DADUR( 1 )*DUR< 1 ) DAB=DADUR< 2 )*DUR ( 2

)

C

RETURN

89
58

59

6 o 61

62

63

64

65

C

C

c30

c

FORMAT ( 2X.3HFXX, I 1, 5HWAVE= . F 1 0 . 7, 7HMUR( 1 )= ,F5. 2,7HMUR( 2)=,F5.2»10

2HDELMUR { 1 ) = ,F5. 2 , 1 OHDELMUR ( 2 )= »F5. 2 , 22 HNORM AL I Z AT I ON FACTO PS= , 2F5.

33)

END

90
12

3

4

5

6

T

8

9

1C

1 1 12

13

14

1 5

16

17

1 8

19

2C

2 1 22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

4 1 42

43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

SUBROUTINE ANISO ( JC , L OGHKL , DM S M AL * X SM ALL * D M, X , M , NN , N

)

C

C **************************** ** ************* *************** ********

c

C SUBROUTINE ANISO MAKES PLOTEM PLOT EACH CLASS OF REFLECTIONS IN TURN.

C DM IS WHOLE SAMPLE POPULATION FOR ORDINATE IN PLOTS, X IS SAME FOR

C ABSCISSA, DMSMAL IS THE SUBSET OF DM WHICH IS ACTUALLY PLOTTED,

C X SMALL IS SAME FOR X, JC IS TOTAL NUMBER IN DM POPULATION, LOGHKL

C IS ONE WORD FOR EACH DM WHICH TELLS WHICH CLASSES THE REFLECTION

C FALLS IN . DIMENSION DM(NN) ,X(NN) , D MS MALI N , 5 ) ,XSMALL(N,5) , LOGHKL ( NN

)

C M IS FLAG INDICATING CLASS TYPE

C SET UP THE BIT NUMBER ,11, WHICH GOVERNS CATALOGUE LOOK-UP FOR

C REFLECTION CLASS.

C MAKE SURE M . GT . C

.

IF (M.LT.O) RETURN

C NON-ZERO VALUE OF M TELLS REFLECTION PLOTTING AND CURVE FITTING

C ROUTINES TO DIVIDE REFLECTIONS UP INTO CLASSES.

IBEG=0

I STOP= 1

0

II = 0

C SET OPTION SWITCH =0 FOR PLOTTING ALL OCTANTS TOGETHER IN ONE PLOT

C PER CLASS.

NOC T=0

C UNPACK M TO ESTABLISH WHAT WILL BE PLOTTED.

C GET CLASS OPTION

MT=M/1CC0

C GET NUMBER OF LAYERS TO BE PLOTTED.

MLAY=(M-MT* 10C0 )/l

0

C GET OCTANT OPTION WHICH IS 0,1, 2, 3, 4, OR 7

moct=(m-mt*ioco-mlay*io

)

C DO STANDARD CLASSES FIRST, IF WANTED.

IF ( MT . EQ. 1 .OR .MT. GE. 5) GO TO 20

10 CONTINUE

C SET UP MT VARIABLE FOR LAYER TYPE,1=STD CL ASSES , 2= H, 3=K , 4=L LAYERS.

IF (MT.GE.5) MT= MT —

3

C SET UP RANGE OF BITS TO BE CHECKED.

IBEG= 1

1

I STOP= IBEG+MLAY

C LOOP THROUGH CLASSES AND LAYERS SETTING BIT AND OCTANT NUMBERS.

C

20 CONTINUE

C INITIALIZE II

I I = I BEG

C INITIALIZE OCTANT NUMBER IF SEPARATE OCTANTS ARE TO BE CONSIDERED.

IF (MOCT.GT.O) NOCT=l

30 CONTINUE

C II IS THE REFLECTION CLASS BEING CURVE-FITTED AND PLOTTED

C CALL PLOTTING AND CURVE FITTING ROUTINES

C VARIABLE MCLT TELLS LAYER TYPE AND OCTANT WANTED.

MCLT=NOCT* 1C0+MT

CALL PLOTEM ( JC, X, DM, DMSMAL , XSMALL « MCLT , I I , LO GHKL , NN , N

)

IF (MOCT.LE.O) GO TO 50

C 0 TH OPTION SPECIFIES ALL OCTANTS.

IF (N0CT.EQ.8) GO TO 50

C INCREMENT OCTANT NUMBER

40 NOCT=NOC T+

1

91
58

59

60

6 1 62

63

64

65

66

67

68

69

70

71

72

73

74

75

76

77

78

7<3 80

C CHECK IF OCTANT OPTION SPECIFIES THIS OCTANT.

IF ( MOCT .EQ . 1 . AND. NOCT .LT. 8) GO TO 40

C FIRST OPTION SPECIFIES +H.+K.+L (OCTANT 1) AND -H,-K,-L (OCTANT 8)

IF ( MOCT.EQ. 2. AND. (NOCT. GT.4 .AND. NOCT.LT .8 ) ) GO TO 40

IF (MOCT. EQ. 2. AND. NOCT. EQ. 2) GO TO 40

C SECOND OPTION SPECIFIES OCTANTS 1.3 AND 8.

IF (MOCT .EQ .3 .AND . ( NOCT .GT. 3 .AND. NOCT.LE .7) ) GO TO 40

C THIRD OPTION SPECIFIES OCTANTS 1,2,3 AND 8.

IF (MOCT. EQ. 4. AND. (NOCT. EQ. 4. OR. NOCT. EQ. 6. OR. NOCT. EQ. 7) ) GO TO 40

C FOURTH OPTION SPECIFIES ALL OCTANTS EXCEPT 4,6, AND 7

C LAST OPTION SPECIFIES ALL OCTANTS.

C NOW THAT WE HAVE THE CORRECT OCTANT NUMBER-PRODUCE PLOT.

GO TO 30

C

C INCREMENT BIT NUMBER - BEGIN A NEW CLASS OR LAYER.

50 11=11+1

C RESET OCTANT NUMBER TO FIRST OCTANT.

IF (N0CT.EQ.8) NOCT=l

IF ( I I .LT. ISTOP1 GO TO 30

IF (MT.GE.51 GO TO 10

PETUPN

C

END

92
12

3

4

5

6

7

8

9

1 0

1 1 1 2

1 3

14

1 5

1 6

17

1 8

19

20

21

22

23

24

25

26

27

28

29

30

3 1 32

33

34

35

36

37

38

39

40

4 1 42

43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

<JU

SUBROUTINE CATLOG ( J.K.L'.LOGHKL. IiLAY.NN)

C

c

c

cc

c

c

c

cc

c

c

c

c

c

ccc

c

c

cc

c

c

c

cc

c

******************************************* ***********************

SUBROUTINE CATLOG PACKS UP THE LOGHKL WORD ACCORDING TO WHICH OF THE

CLASSES BELOW THE REFLECTION J,K,L FALLS IN.

MAXIMUM AND MINIMUM VALUES OF MILLER INDICES IN THIS DATA SET

COMMON /D/ MAXH, MI NH .MAXK, MI NK , MAXL, MI NL

DIMENSION LOGHKL ( N N

)

SET BITS IN LOGHKL WORD

REFLECTION WITH INDICES

DEPENDING ON WHICH

J * K » L FALLS IN.

REFLECTION classes the

BIT 0

123

4

5

6

7

8

9

10

BITS FOR

1 1

1 2

1 3

ETC.

20

H . 0 * 0

0 * K * 0

0,0,

L

H, K, 0

H, 0 ,L

0, K,L

H , K , L ( H , K OR L NOT EQUAL TO ZERO)

H, H, H

H , H , L

H, K ,K

H , K , H

LAYER DATA

1 KL , OR HI L OR HK 1 LAYER

2 KL , OR H2L OR HK2 LAYER

3KL , OR H3L OR HK 3 LAYER

1 OKL , OR H10L OR HK 1 0 LAYER

BITS 21 THROUGH 33 FOR FUTURE USE

BIT 33 0 IF

BIT 34 0 IF

BIT 35 0 IF

THUS CCO OR

LAY IS CLASS

DO STANDARD

H IS PLUS, 1 IF

K IS PLUS, 1 IF

L IS PLUS, 1 IF

OCTAL 0 WILL BE

OR LAYER TYPE

CLASSES

IS

IS

IS

MINUS

MINUS

MINUS

THE +++ OCTANT

IF ( K.EQ.O . AND.L.EQ.C ) FLD (0,1. LOGHKL ( I ) ) = 1 IF ( J .EQ.O .AND.L.EQ.O ) FLD ( 1 , 1 , LOGHKL! I ) )=

1

IF ( J. EQ.O. AND. K.EQ.O) FLD ( 2 , l , L OGHKL ( I )

) =

1

IF (L.EQ.O) FLD( 3, 1 , LOGHKL (I ) ) =

1

IF (K.EQ.O) FLD( 4 , 1 .LOGHKL ( I )

) =

1

IF (J.EQ.O) FLD ( 5 , 1 .LOGHKL { I ) )=

1

IF ( J .NE .0. AND.K. NE. 0. AND. L. NE.O ) FLD( 6 , 1 , L OGHKL ( I ) ) = 1 IF ( J .EQ .K .AND.K. EQ .L) FLD ( 7 , 1 , LOGHKL ( I )

) =

1

IF (J.EQ.K) FLD(8,1, LOGHKL (I ) ) =

1

IF (K.EQ.L) FLD<9, l,LOGHKL(I ) )=1

CLASSIFY THE APPROPRIATE LAYER

IF (LAY.LE.l ) (LAY. EQ.2

)

( LAY. EQ. 3)

(LAY . EQ .4

)

10 N= 1,10

NB= N+ 1 0

IF

IF

IF

DO

GO TO

NH= J NH=K

NH=L

20

1 0

(NH.EQ.N.OR.NH.EQ.-N) FLD(NB, 1 , L 0 GHKL ( I) )=1

CONTINUE

IF

93
58 C

59 20 CONTINUE

60 C ASSIGN REFLECTION TO ITS PROPER 0

6 1 IF (J.LT.O) FLD( 33 * 1 * LO GHKL ( I ) ) = 1 62 IF (K.LT.O) FL D ( 34 r 1 » LOGHKL ( I ) ) = 1 63 IF IL.LT.O) FLDC35 * 1 « LO GHKL ( I ) )=

1

64 RETURN

65 c 66 END

94
12

3

4

5

6

7

8

9

1 0

1 1

1 2

1 3

14

15

16

17

1 8

19

2 ? 21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

3 T

38

39

40

41

42

43

44

45

U

U

SUBROUTINE CELL

*****:(£ *:£:£* ******** ******** ***************** *** * 4c ** * * * * * * ** * * * * * * * * READ THE CELL CARD. TRANSFORM THE ANGLES IF NECESSARY,

REAL CELL, AND CALCULATE THE RECIPROCAL CELL

SUBROUTINE TO

WRITE OUT THE

I/O UNITS

COMMON /IO/ I N, NOUT, NDRUM, NT APEA ,NT APEB

RECIPROCAL CELL PARAMETERS

COMMON /G/ ASTAR, BSTAR , CSTAR»COS AST, COSBST, COSGST

PARAMETER REQUIRED TO CHANGE DEGREES TO RADIANS

RAD=2.*3.1416/36C.

READ CELL, TRANSFORM ANGLES IF NECE SS ARY , WR I T E OUT CELL

READ (IN, 10) A,B,C ,AL,BE,GA

ANGLES MAY BE DEGREES OR COSINES

IF (AL.GT.l.) COSA=COS ( ALTRAD)

CO S A= AL

AL=A COS ( AL ) / R AD

COSB=COS( 3E*RAD>

C0SB=8E

SE=ACOS( BE)/RAD

COSG=COS( GA*RAD)

COSG=GA

G A= ACOS ( GA) /RAD

A , B , C , AL , B E , G A

IF

IF

IFIF

IF

IFIF

IF

( AL .LE.l . )

( AL.LE. 1 .

)

( BE .GT. 1 .

)

( BE.LE. 1 .

)

( BE.LE. 1 .

)

( GA .GT. 1 . )

( GA .LE.l .

)

(GA.LE. 1.

)

WRITE (6,20)

AL= AL* RAD

BE=BE*RAD

GA=G A*R A D

S=( AL+BE+GA ) /2

.

VOL= 2 • *A *B*C* SORT ( SIN( S ) * S I N ( S- AL ) * SI N ( S-BE ) *S I N ( S- GA ) ) CALCULATE RECIPROCAL CELL

ASTAR = B*C*SI N( AD/VOL

BSTAR=C*A*SIN( BE )/ VOL

CST AR= A* B* SIN(GA)/VOL

COSAST=( CO SB*COS G— COS A ) /(SIN(BE)*SIN(GA))

COSBST= ( COSG*COSA— COSB )/(SIN(GA)*SIN(AL))

COSGST=( COSA*COSB— COSG ) / (S IN ( AL ) *S I N( BE)

)

RETURN

1C FORMAT ( 13X ,3F8.3, 3F9. 5)

20 FORMAT ( 32H THE CELL DIMENSIONS READ IN ARE , 3F 1 0 . 4 , 3F 1 0 . 2

)

C

END

95
123

4

5

6

7

8

9

10

1 1 12

13

14

15

16

17

18

1 9

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

48

49

50

51

52

53

n

n

SUBROUTINE DMCALC { JC , S K ALE , DM , L OGF C » NN

}

C

C X^**’*'**'**'*'********##*'*** ** ## * # #** # ## 5jc *:$: ********$*************$

C

C ROUTINE CALCULATES DM*T HE 'UNIT DEVIATE*.

C ALSO THE MEAN WHICH SHOULD BE CLOSE TO ZERO.

C ALSO THE E.S.D. WHICH SHOULD BE CLOSE TO ONE.

C

C WE NOW HAVE JC COMMON REFLECT I ONS . THE T R CATALOGUE NUMBERS ARE PACKED

C INTO THE LOGFC ARRAY.

COMMON / 10/ IN, NOUT

C NOW CALCULATE DM(1),S

DIMENSION DMINN) ,LOGFC( NN)

RNUM=0 .0

RDEN= 0 .0

SUMM=0 . SUM=0.

SS = SKALE **

2

DO 10 1=1, JC

C FIND THE SEQUENCE NUMBERS OF THE TWO CORRESPONDING DEFLECTIONS, ONE

C IN EACH DATA SET

J=LOGFC I I)/l ceooo

K=LOGFC( I>-10C000*

J

C READ THEM FROM THE DRUM

CALL DRUMRD < JHKL , F A , SGA , NT , J

)

CALL DRUMRD ( JHKL , FB , SGB , N T , K

)

C CALCULATE THE DM QUANTITY

DMC I ) = CFA-SKALE*FB J/SQRTI S GA ** 2 + SS *SGB*4 2

)

SUM=SUM+DM( I )4*2

SUMM=SUMM+DMC I

)

RNUM=RNUM+ABS( DMC I ) ) RDEN=RDEN+F A/SGA

C CALCULATE MEAN AND ESD OF THE DM DISTRIBUTION

10 CONTINUE

SUMMA=SUMM/JC

ESD=0.

0

DO 20 1=1, JC

ESD=ESD+CDM( I }-SUMMA)**2

20 CONTINUE

ESD=SQRT (ESD/( JC-1 ) ) R=RNUM/RDEN

WRITE ( NOUT , 30 ) JC , SUM , SUM M , SUMM A , ES

D

WRITE C NOUT , 40 ) R

RETURN

30 FORMAT C1H0.16H SUMMARY OF THE ,I5,14H DM S T AT I ST I CS// 1 7H SUM OF D

2M**2 IS ,E10.4»14H SUM OF DM IS ,E10.4//28H MEAN OF DM DISTRIBUTIO

3N = , F8 • 3 , 27H ESD. OF DM DISTRIBUTION - ,F6.3)

40 FORMAT (1H0,6H WR = ,F8.4,51H FOR COMPARISON WITH WR OBTAINED IN

2 REFINEMENTS)

C

END

96
123

4

5

6

7

8

9

10

1 1 12

13

1 4

15

16

17

1 8

19

20

21

22

23

24

25

26

27

28

29

30

3 1 32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

58

SUBROUTINE DMDQ ( 0 M , LOGFC , JC , J J , KK * LL , X , LOG HKL , ND . NQ , NCL , N N

)

I C

C ******** I*********************************** *************** ********

cC SUBROUTINE TO READ EACH REFLECTION IN TURN FROM THE DRUM

C GENERATE THE MILLER INDICES CALCULATE THE D VALUE AND PREPARE FOR

C D-STAR PLOT OR Q PLOTS.

C DM IS ORDINATE VALUES. X IS ABSCISSA. SEQUENCE NUMBERS OF

C CORRESPONDING REFLECTIONS ARE PACKED INTO LOGFC, INFORMATION ON

C REFLECTIONS CLASSES EACH REFLECTION BELONGS TO I S PACKED INTO LOGHKL

C JJ.KK AND LL ARE PACKING ORDER OF MILLER INDICES IN JKL WORD, ND IS

C A PARAMETER WHICH CONTROLS WHAT IS TO BE PLOTTED AGAINST WHAT.

C RECIPROCAL CELL PARAMETERS

COMMON /G/ ASTAR,BSTAR,CSTAR,COSAST,COSBST, COSGST

COMMON /EXPT/ WA VE , SK AL E , UR ( 2 ) , D UR ( 2 ) , NF

DIMENSION DM(NN) , LOGFC (NN) , X ( NN ), LOGhKL < NN

)

C FIRST I N I T I ALL I ZE LOGHKL AWAY SO THAT WE HAVE ALL BITS ZERO.

DO 10 1=1, JC

LOGHKL ( I

) =0

10 CONTINUE

DO 30 1=1, JC

c obtain indices of reflections associated with dm value.

C READ REFLECTION INDICES AND FCALC FROM THE DRUM

LB=LOGFC ( I )/10CC0 r' L A=L OGF C ( I) -LB* 100 00 0 CALL DRUMRD ( J KL , F , S G, NT , L A

)

C GENERATE the MILLER INDICES

CALL HKLGEN ( J KL , J , K , L , J J , KK , LL

)

C FLOAT MILLER INDICES.

H J= J HK=K

HL=L

C CALCULATE D-STAR VALUE IF WANTED.

IF (ND.NE. 1 .AND.NQ.NE.l ) GO TO 20

TEMPT H= (HJ*AST AR ) * *2 + ( HK*BST AR)**2+(HL*CSTAR) ** 2+ 2 . *HK*HL*B STA R *C

S

2TAR*CO SAST+2 .*HL*H J*AST AR*CSTAR*COSBST+2 .* H J* HK* AST A R *BST AR*COS GST

C STORE D VALUE IN X ARRAY WHICH IS LATER PLOTTED AGAINST DM ARRAY

X( I

)

SQR T { TEMPTH)

C CHECK IF Q (F**2/SIN( 2THETA) ) VALUES ARE WANTED.

IF ( NQ .NE. 1 ) GO TO 20

SIN2T=. 5*WAVE*X( I

)

C STORE LOG Q IN ARRAY X C CHECK IF DATA F OR F**2

IF (NF.LE.O) F=F**2

X( I )=ALOG10(F/SIN2T)

20 CONTINUE

C CHECK IF REFLECTION CLASSES ARE TO BE PLOTTED SEPARATELY. IF NOT,

C return.

C OBTAIN CLASS OR LAYER TYPE

NL AY = NCL

IF ( NLAY.LE.O ) GO TO 30

IFINLAY. GE .5 ) NL A Y=NLA Y—

3

C CATALOG REFLECTIONS ACCORDING TO CLASS OR LAYER AND OCTANT.

CALL CATLOG(J,K,L, LOGHKL, I , NLAY , NN

)

30 CONTINUE

RETURN

C

END 97
cm

ro

ifl

vo

I

1

8

9

SUBROUTINE DPCALC COP.PA.PB^SA.SB)

cc *#**** ** $$$$$$$***$$*$:***$$$$$**$$*$$* ###*#**:*

C

C STANDARD CALCULATION OF DP QUANTITY FOR HALF NORMAL PLOT

DP= ABS(ABS(PA)— ABS(PB) )/SQRT( SA**2+SB**2)

RETURN

C

END

98
123

4

5

6

7

8

9

1 9

1 1

1 2 13

14

15

15

17

1 8

19

20

21

22

23

24

SUBROUTINE DRUMRD ( JHKL ,

F

? SG , NT , I

)

C

C ******** *$$***$$***4: ******** **************************************

C

C SUBROUTINE TO READ RECORD (REFLECTION) FROM RANDOM ACCESS MASS

C STORAGE (DISK, DRUM).

C READS ONE RECORD PER CALL

C MAY BE EASILY MODIFIED FOR DIFFERENT I/O.

COMMON /IO/ IN.NOUT.NDRUM.NTAPEA.NTAPEB

C NEXT TWO STATEMENTS ARE A MECHANISM TO READ BOTH FLOATING AND FIXED

C POINT QUANTITIES FROM THE DRUM IN ONE OPERATION.

DIMENSION LDRUM( 4) ,FDRUM( 4)

EQUIVALENCE ( LDRUM.FDRUM)

C J POINT INDICATES LOCATION OF REFLECTION ON DRUM

JPO INT= ( I- 1 ) *4

C READ FROM DRUM, BUT FIRST COUNT FROM BEGINNING OF DRUM.

CALL NTRAN ( NDRUM , 6 , - 1 0 0C 0 0 , 6 , J POI NT , 2 , 4 , LDRU M, LOPT , 22

)

JHKL=LDRUM( 1

)

NT=LDRUM (2

)

F=FDRUM( 3)

SG=FDRUM (4

)

RETURN

C

END

99
123

4

5

6

7

a9

1 0

1

1

1 2

1 3

14

15

1 6

17

18

19

20

21

22

23

24

25

26

SUBROUTINE DRUMRT (N,IJ,F,S)

C

C **** afc******^** ************ **************** ********** **************

C

C SUBROUTINE TO WRITE REFLECTION AND RANDOM ACCESS MASS STORAGE (DISK,

C DRUM)®

C MAY BE EASILY MODIFIED FOR OTHER TY PR I/O.

C N IS FLAG FOR UNOBSERVED REFLECTIONS.

C IJ ARRAY IS MILLER INDICES, F IS STRUCTURE FACTOR, S IS S I G( FO

)

COMMON / 10/ IN, NOUT ,NDRUM» NTAPEA,NT APEB

DIMENSION I J( 3) ,LDRUM( 4 ) , F DR UM ( 4

)

C LDRUM AND F DRUM ARE EQUIVALENT SO THAT FIXED AND FLOATING POINT

C QUANTITIES CAN BE WRITTEN ON THE DRUM IN ONE OPERATION, WHICH SAVES

C TIME.

EQUIVALENCE (LDRUM ,FDRUM)

C PACK MILLER INDICES INTO ONE WORD TO SAVE COMPUTER SPACE AND WRITING

C TIME.

LDRUM ( 1)=( IJ( 2) + 100)*l 0 0 0+ IJ(3)+1C0 + (IJ( 1 ) + 100)*1000C00

LDRUMC 2 ) =N

FDR UM( 3) =F

FDRUM ( 4 )=S

C WRITE THE DRUM.

CALL NTRAN (NDRUM, 1.4, LDRUM, LOP, 22)

RETURN

C

END
123

4

5

6

7

8

9

10

1 1 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

39

40

41

42

43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

SUBROUTINE FINGFO ( JY « L OGF C.DM.JJ.KK.LL, NFT EST , ID* I J M I N » UMAX

2 NT APE * J Z

)

C ************************ ** ******** ************* ******* **   ** ****

c

C SUBROUTINE TO READ REFLECTIONS FROM FIN GL S FOURIER TAPE

C FINGLS FOURIER TAPE IS WRITTEN BY PROGRAM RFINE 2.

C RFINE 2 WRITTEN BY L. W. F I NGER , GEOPHYS I C AL LABORATORY , CARNEG IE C INSTITUTE OF WASHI NGTON , D . C

C NTAPE A AND NTAPEB NOT IN COMMON HERE BECAUSE WE WANT TO USE SAME

C SUBROUTINE FOR BOTH TAPES. THEY ARE TRANSFERRED IN VIA NT APE ARGUMENT

C IN CALL STATEMENT.

INTEGER A COMMON /TO/ I N , NOUT , NDR UM

C ID ARRAY IS ONE MILLER INDEX WHICH IS USFD TO PULL REFLECTIONS OFF

C DRUM IN BATCHES OF (SAY) CONSTANT K. THE SORTING ON H«L GOES MUCH

C FASTER WHEN ONLY REFLECTIONS WITH ONE INDEX THE SAME ARE SORTED AT

C THE SAME TIME.

PARAMETER N=5C, NN A= 18 0 9* NN=7200» NNB=1801, NNC=3601, NND = 54^1.

2 MM=60, MMM=4CC, NBUFF=70C, MMMM=32

DIMENSION ID(NN)

C DUMMY IS USED TO OVERLAP IN CORE SOME QUANTITIES WHICH ARE NOT

C NEEDED AT SAME TIME.

C FBUF AND 1FBUF ARE BUFFERS USED IN READING FINGLS FOURIER TAPE F,

COMMON /H/ DUMMY ( N N

)

C I J ( 3 ) IS MILLER INDICES. DM IS ORDINATE IN PLOTS. LOGFC IS WORD

C CONTAINING PACKED SEQUENCE NUMBERS ON DRUM OF CORRESPONDING

C REFLECTIONS. RMAT AND TRANS ARE SYMMETRY ELEMENT COMPONENTS WHICH ARE

C NOT NEEDED HERE.

DIME NS ION FBUF ( 504 ) , IFBUF( 5<'4),IJ(3),DM(NN) .LOGFC ( NN ) , RM AT ( 3 , 3, 24 ) 2. TRANS (3.24)

C EQUIVALENCES TO CONSERVE STORAGE

EQUIVALENCE ( FBUF ( 1 ) , IFBUF( 1 ) . D U MM Y ( 1 ) ,VOL. RMAT ( 1.1,1),

2 TRANS (1,1))

C READ FIRST RECORD, CONTAINING SYMMETRY ELEMENT INFORMATION, ON FINGLS

C FOURIER TAPE.

READ ( NT APE ) NSYM.VOL, ( ((RMAT(I,J,K).I=1,3) .J=1.3),K=1,NSYM),((TRA

2NS ( J , I ) , J=1 ,3), 1=1 , N SYM

)

JX = C

C BEGIN BY READING REFLECTIONS INTO BUFFER

GO TO 20

C CHECK WHETHER ALL BUFFER PROCESSED.

10 IF (NFBR-500) 30.20.2C

C READ INTO BUFFER.

20 READ (NT APE) FBUF

C ZERO COUNTER

NFBR=0

C GET MILLER INDICES

30 IJ( JJ )=IFBUF(NFBR+1

)

C CHECK FOR END OF LIST, INDICATED BY J=99.

IF ( I J ( J J ) —9 9 ) 40,60,40

40 IJ(KK)=lFBUF(NFBR+2)

I J(LL)=IFBUF(NFBR+3)

C FIND MAXIMUM H K AND L

CALL MAXHKL ( I J . J J . KK , LL

)

FOBS=FBUF( NFBR+4)

FC=FBUF( NFBR+5)

101
58

59

60

6 1 62

63

64

65

66

67

68

69

70

71

72

73

74

75

76

77

78

79

80

81

82

83

84

85

86

87

88

89

90

91

92

93

94

95

96

97

98

99

ICO

10 1

1 02

10 3

104

105

1 C 6

1 C7

108

109

110

1 1 1

1 1 2

1 1 3

1 14

1 1 5

C LESS THAN INDICATOR.

A= IFBUFC NFBR+8

)

C SIGCFOS

SGA=FBUF(NFBR+9)

C EXTINCTION CORRECTION FOR THIS REFLECTION

EXY=FBUF CNFBR+1 0

)

NFBR= INF BR + 10

)

DIMENSION NSAMPC2Q 1 , FEXY C 20

1

NEXY=( EX Y+ .02) * 20

IFCNEXY.LE.O } NEXY=1

NS AM PC NEXY)=NSAMP< NEXY) +1

NTOT AL=NTOT AL + 1 C REJECT REFLECTIONS REJECTED FOR SEVERE EXTINCTION IN RFINE.

IFC EXY .LE. .7} GO TO 10

C REJECT LESS-THANS

IF ( A . EQ • 1 ) GO TO 10

C INCREMENT VARIOUS COUNTERS. JX IS NUMBER OF REFLECTIONS READ IN OFF

C THIS TAPE. JY IS TOTAL NUMBER OF REFLECTIONS READ IN, JZ IS TO SPACE

C REFLECTIONS IN PAIRS ON DRUM BECAUSE FO AND FC ARE TREATED AS

C SEPARATE REFLECTIONS, BUT THEY ARE NOT SCALED BECAUSE IT IS ASSUMED

C THAT FINGLS SCALING IS MORE APPROPRIATE.

JX= JX+

1

JY= JY+

JZ = J Z+2

C DO NOT OVERWRITE REFLECTION ARRAYS.

IF CJX.GE.NN) GO TO 60

C IF TWO SEPARATE SETS OF FO* S ARE BEING READ IN, WE HAVE TO SAVE THE

C VALUES USED LATER IN THE SORT ROUTI NE I . E . , THE MAXIMUM AND MINIMUM

C VALUES OF THE INDEX KEPT IN CORE.

IF C NFTEST.NE.3l GO TO 5C

C GET MAXIMUM AND MINIMUM INDICES FOR LATER PLOTTING OF DM VALUES

C AGAINST RANGES OF INDICES.

CALL MAX HK L ( IJ, JJ.KK.LL)

C THE ID ARRAY LATER GIVES THE POSITION OF THE REFLECTION ON THE DRUM,

C THE CONTENTS OF THE ID ARRAY IS THE MOST SLOWLY VARYING INDEX OF THE

C REFLECTION.

IDC JY) = I JC 1

5

C IJMIN AMD UMAX GIVE THE RANGE IN THE VALUES OF THE MOST SLOWLY VaRYC ING INDEX.

IF C IJ( 1 ).LT. IJMIN) I J M I N= I J ( 1 ) IF C XJ( 1 >.GT.I JMAX) I JMAX=IJ ( 1

)

C WRITE REFLECTION ON DRUM.

CALL DRUMRT ( A, I J, FOBS, S GA

)

GO TO 10

C IF FO AND FC ARE BEING COMPARED, GET DM VALUE HERE.

5 0 DM( J Y) = ( ABS C FOBS ) -* ABS IFC) )/SGA

C STORE LOCATIONS OF FO AND FC ON DRUM.

LOGFCI JY )=JZ*1CC000+JZ+1

C WRITE REFLECTION ON DRUM WITH FOBS FOR FIRST WRITE.

CALL DRUMRT C A , I J , FOBS , SGA

)

C WRITE SECOND REFLECTION ON DRUM OR IN THIS CASE SECOND HALF OF FIRST

c REFLECTION

C EXY IS EXTINCTION FACTOR. STICK IT IN PLACE OTHERWISE USED FOR SECOND

C S I G MA CFO).

CALL DRUMRT C A , I J, FC , EXY

)

GO TO 10

C WRITE HOW MANY REFLECTIONS READ IN AND WHICH TAPE THEY ARE READ FROM

60 WRITE ( NOUT ,90 ) NTOTAL , NTAPE , JX

102
116 DO 70 1=1.20

117 70 FE X Y ( I } = I/20

.

118 WR ITE (6.100)

1 1 9 WRITE (6,110) (FEXY( I ) , NS AMP ( I ) , 1=1 , 20

)

1 20 DO 80 1=1,20

121 FEX Y ( I )=0.

122 80 NSAMPC I ) =0

123 RETURN

124 C

125 c 126 c 127 90 FORMAT (I10.28H REFLECTIONS READ FROM UNIT .I6.9H OF THESE, 16, 63H

128 2 HAD EXTINCTION COEFFICIENTS GREATER THAN 0.7 AND WERE 0BSERVED/37H

1 29 3 ONLY THESE WERE SAVED FOR LATER USE.)

1 30 1 00 FORMAT ( 49H DISTRIBUTION OF EXTINCTION FACTORS IN INPUT LIST/40H

131 2 TOP OF RANGE NUMBER OF REFLECTIONS)

132 1 1C FORMAT (F12.2, 8X, 110)

133 c 134 END

103
I

123

4

5

6

7

8

9

1 0

1 1

1 2

1 3

1 4

1 5

1 6

1 7

1 8

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

4 2 43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

SUBROUTINE FIT CHK C DM, X * JC , D MSM AL , X SMALL , L , K, AL INE, BL I NE , A CUR VE

,

2 BCURVE, CCUR VE, SUMLN , SU MC V , M , I I , LOGHKL « J , ACU BI C , BCU BT C, CCUBIC,

3 DCUBI C , SUM CUB > C

c *******4 ****************** *: ******* ********************** **********

c

C SUBROUTINE TO CALCULATE POINTS FOR LEAST SQUARES LINE AND CURVE AN4

C TO ALSO CALCULATE GOODNESS OF FITS FOR L.S. LINE, QUADRATIC AND

C CUBIC CURVES.

PARAMETER N=5C , NNA=1800, NN=72CC, NNB=1801, NNC=36C1, NND=54Fl,

2 MM= 60 , MM M=40 0 , NBUFF=7C0, MMMM=32

DI MENS ION DM(NN) , XCNN) , DMSMALI N, 5) , XSMALLI N ,5) , LOGHKL ( NN)

CC J=1 INDICATES THAT X ARRAY CONTAINS NORMAL QUANTILES, SET FLAG NQ

NQ=0

IF CJ.EQ.l) NQ=1

C RESET J SO IT CAN KEEP COUNT OE THE NUMBER OF POINTS IN PLOT

j=e

JK=0

C LHALF IS CONCERNED WITH PICKING THE PLOTTED POINT FROM THE MIDDLE OF

C THE BATCH WHEN EVERY LT H POINT IS PLOTTED.

C PLOTTING ROUTINES ARE SET UP FOR 50 POINTS MAXIMUM.

LHALF=L/2

SUMLN=0

SUMC V=0

SUMCUB=C

SNEG=0 .

0

SP0S=0 .0

NPQS=0

NEG=0

N I N=0

C LOOK OVER ALL QUANT IT IES REFLECTIONS, PARAMETERS, ETC.

DO 80 1=1, JC

C LEAVE OUT TAILS OF DISTRIBUTION IF X ARRAY CONTAINS NORMAL QUANTILES

IF (NQ.NE.l) GO TO 10

IF C ABS( X( I ) ) oGT.2.0) GO TO 80

10 CONTINUE

C PULL OUT OCTANT NUMBER, 0 THROUGH 8, 0 = AL L OCT ANTS , 1 = + + +, ET C

.

NOCT=M/l 00

IF ( ( M—NOCT* 1 CO I.LE.O) GO TO 20

C IF ONLY PARTIAL DATA IS TO BE PLOTTED, SELECT DATA HERE

C THIS PULLS OUT THE APPROPRIATE REFLECTIONS BASED ON THE BITS PLACED

C IN THE LOGHKL WORD BY SUBROUTINE CATLOG

NTEST= FLD( I I , 1 , LOGHKL c I > ) NTEST=ABSCNTEST)

IF ( NTEST • EQ » 0 ) GO TO 80

20 CONTINUE

C CHECK IF REFLECTION BELONGS TO OCTANT WANTED

IFC NOCTeLE.O ) GO TO 25

NTEST= FLD(33,3, LOGHKL C I I ) NTE ST=ABS! NTEST)

IF( ( NQCT-1 ). NE. NTEST) GO TO 80

25 CONTINUE

C CALCULATE POINTS ON LEAST SQUARES LINES

D= AL INE+BLINE*X( I

)

C CALCULATE POINTS ON LEAST SQUARES CURVES

E=ACUR VE + BCUR VE*X C I ) + CCURVE*X« I )**2

> 58

» 59

5 60

i 61

i 62

t 63

tj 6

1

6 65

( 6 !

j 61

|

61

I 6'

V

7

7 7

" 7

• 7

• 7

"

7

i

5 !

5 f

17 I

!» I

I

I

I

[

t

104
58

59

6 C 6 1 62

63

64

65

66

67

68

69

7*

71

72

73

74

75

76

77

78

79

80

81

82

83

84es86

87

88

89

90

91

92

93

94

95

96

97

98

99

ICO

101

10 2

1 03

104

105

106

107

1 C 8

109

no

11

1

112

113

114

115

IF (K.LT.4) GO TO 3" C CALCULATE POINTS ON L.S. CUBIC CURVE IF ENOUGH POINTS IN PLOT.

F=ACUBIC+BCUBIC*X< I )+CCUBIC*X( I )**2+DCUB IC*X< I ) **3

C SUM UP FOR GOODNESS OF FIT PARAMETER.

SUMCUB= SUM CUB+ ( F— DM( I ) ) **2

30 CONTINUE

SUMLN=SUMLN+ (D—DM( I ) )** 2 SUMCV=SUMCV+( E— DM( I ) )**2

C COPY ONLY REQUIRED POINTS INTO PLOTTING ARRAYS

C CALC QUANTITIES FOR INDICATION OF DM SCATTER.

C FIRST CHECK IF NECESARY CONDITIONS ARE MET.

C OMIT IF NORMAL QUANTILES OR LESS THAN 6 POINTS IN INTERVAL

IF (L .LT .6 .OR .NQ.EQ. 1 ) GO TO 70

C CALC DIFFERENCE FOR POINT ON LS LINE

DEL=DM( I )-D

IF (DEL) 40,50,50

C SUMM UP NEGATIVE DEVIATIONS

40 SNEG=SNEG+DEL

NEG=NEG+1

GO TO 60

C SUMM UP POSITIVE DEVATIONS

50 SPOS=SPOS+DEL

NPOS=NPOS+

1

60 NIN=NIN+1

IF (NIN.LT.L) GO TO 70

I F ( N EG .L E • 0 ) NEG=

1

IF(NPOS.LE.O )NPOS=l

C HAVE ALL POINTS IN INTERVAL NOW

C STORE QUANTITIES IN DS MALL ARRAY

DMSM AL <J»4)=SNE G/NEG

DMSMAL ( J ,5)=SPOS/NPOS

C RESET COUNTERS AND SUMS

SPOS=0 .0

SNEG= 0 .

0

NPOS=0

NEG = 0 N I N= 0

70 CONTINUE

C JK COUNTS HOW MANY POINTS WE HAVE SO FAR

JK= JK+

1

C NTEST IS THE POINT IN THE MIDDLE OF THIS RANGE.

NTE ST= ( JK/L) *L+LHALF

IF ( JK ,NE. NTEST ) GO TO 80

C IF THIS IS THE POINT WE WANT, COME HERE.

J=J+ 1 C STORE ORIGINAL DM VALUE IN PLOTTING ARRAYS. THESE ARE DMSMAL AND

C X SMALL

•

DM SM AL ( J , 1 )=DM( I

)

C STORE VALUE FOR LINE

DMSMAL! J , 2) =D

C STORE VALUE FOR CUBIC CURVE

DMSMAL ( J , 3 ) =F

C STORE DISTRIBUTION QUANTILE FOR PLOTTING ARRAYS. THEY HAVE TO BE TOLD

C 3 TIMES.

X SMALL ( J ,1 )=X( I ) XSMALLI J ,2 ) =X( I

)

X SM ALL ( J ,3 ) = X( I ) XSMALL ( J ,4 >=X ( I

)

105
1 16

1 1 7

1 10

119

120

121

122

123

124

X SM ALL (

J

,5 >=X( I > 80 CONTINUE

c calculate goodness of fit of least squares LINE AND CURVES to data.

SUMLN=SQRTCSUMLN/( K-2)

)

SUMC V= SQRT ( SUMC V/( K-3) ) IF (K.GT.4) SUMCUB = SQRT ( SUMCUB/< K-4 ) ) RETURN

C

END

106
12

3

4

5

6

7

8

9

10

1 1

1 2

1 3

1 4

15

16

17

1 8

19

20

21

22

23

24- 25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

4 2 43

44

45

46

47

48

49

50

SUBROUTINE FREAD ( NT AP E , ID , J . J J , KK , LL , I JMI N , I JM AX

)

C i C I**:************************** *************** ***********************

C

C SUBROUTINE TO READ THE INPUT REFLECTIONS FROM N TAPE TO WRITE THOSE R

C NS ON THE DRUM AND TO KEEP IN STORE IN THE ID ARRAY THE SLOWLY VARYIN

C AS A CATALOGUE NUMBER

C MAY EASILY BE MODIFIED FOR UNFORMATTED I/O

INTEGER A COMMON /IO/ IN.NOUT.NDRUM.NTAPEA.NTAPEB

PARAMETER N=50 * NNA=1800, NN=7200, NNB=1801 , NNC=360 1 , NND=540l,

2 MM=60. MMM=4C0t NBUFF=700 , MMMM=32

DIMENSION I J ( 3 ) » ID( NN

)

J=0

DO 10 1= 1 » NN

C READ IN REFLECTION. STOP READING ON END OF FILE.

READ (NTAPE,3C,END=20,ERR=10) A , I J ( J J ) . I J ( K K ) , I J ( L L ) , F ,

S

C DISREGARD ANY ZERO VALUES OF F

IF (E.LT..001) GO TO 10

C THROW OUT UNOBSERVED REFLECTIONS

IF (A.E0.2) GO TO 10

C MAKE SURE WE DON'T GET 0.0,0. FROM BLANK CARD OR SOME OTHER GOOF.

IF ( I J < 1 ) . EQ .0 . AND . IJ( 2 ) .EQ. 0. AND. I J( 3) . EQ . 0 ) GO TO 10

C WE HAVE GOT A GOOD REFLECTION.

J=J+ 1 C GET MAXIMUM AND MINIMUM INDICES FOR LATER PLOTTING OF DM VALUES

C AGAINST RANGES OF INDICES.

CALL MAXHKL ( I J , J J , KK , LL

)

C THE ID ARRAY LATER GIVES THE POSITION OF THE REFLECTION ON THE DRUM,

C THE CONTENTS OF THE ID ARRAY IS THE MOST SLOWLY VARYING INDEX OF THE

C REFLECTION.

I D ( J)=IJ( 1

)

C IJMIN AMD UMAX GIVE ThE RANGE IN THE VALUES OF THE MOST SLOWLY VARYC ING INDEX.

IF ( I J ( 1 ) .LT. I JM IN ) I J M I N= I J ( 1 ) IF ( IJ( 1 ) • GT. IJMAX) I JMAX= I J ( 1

)

C WRITE REFLECTION ON DRUM

CALL DRUMRT (A.IJ.F.S)

10 CONTINUE

C READING OF REFLECTIONS COMPLETED

20 WRITE (NOUT.40) J.NTAPE

RETURN

C

C X-RAY BCD FORMAT

C

C

30 FORMAT ( 13X, I2.3I4.F1C. 2.F10.4)

40 FORMAT (I10.27H REFLECTIONS READ FROM TAPE, 15)

C

END

107
12

3

4

5

6

7

8

9

1C

1 1 12

1 3

1 4

15

SUBROUTINE HKLGEN (JKL»J»K»L»JJ»KK»LL)

C

c ****************** ************************************************

c

C SUBROUTINE TO UNPACK THE MILLER INDICES

DIMENSION I J ( 3

)

I JC 1 )= JKL/1C00000- 100

I J(2)=C JKL-I IJ(l)+100)*lC000O0)/lO00-100

IJ< 3) = JKL-{ IJ(1)+100)*1C0C 00 0-(IJ(2)+100)*100 0-100

J= I J ( J J

)

K= I J ( KK

)

L=I JILL)

RETURN

C

END

108
12

3

4

5

6

7 ' 8

9

1C

1 1

1 2

13

14

15

1 6

17

18

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

48

49

50

SUBROUTINE INFING C NAME , OCC A , S I GO* X YZ, S I GX Y Z , BBET A , S I GB * NA TOM

,

2 NT APE)

C

c ******************************************************************

c

C SUBROUTINE TO READ PARAMETERS INFORMATION FROM END OF FINGER FOURIER

C NAME IS THE ATOM NAME, TO BE WRITTEN IN A6

C OCCA IS THE ATOM OCCUPANCY

C SIGO IS THE STANDARD DEVIATION OF THE OCCUPANCY

C X YZ ARE THE POSITIONAL PARAMETERS

C SIGXYZ ARE THE STANDARD DEVIATIONS ON THE POSITIONAL PARAMETERS.

c bbeta are the anisotropic thermal parameters.

c sigb are the standard deviations on bbeta.

C TAG IS THE FINGLS REPRESENTATION OF THE ATOM NAME.

C ISOT, ISCAT, OCCUP AND SITE ARE CONCERNED WITH TOTAL SITE

C OCCUPANCY, NUMBERS OF SCATTERING FACTORS, ETC. AND ARE NOT USED HERE.

COMMON /IO/ IN, NOUT

PARAMETER N=50, NNA=18C0, NN=72CQ, NNB=1801, NNC=3601, NND=5401,

2 MM=60, MMM=4C0» NBUFF=70C, MMMM=32

COMMON /H/ DUMMY ( N N

)

DIMENSION NAME (MM) , OCC A ( MM ) , S I GO ( MM ) , X YZ ( 3 , MM ) , S I GX Y Z ( 3 , MM )

,

2 BBETA (6 ,MM ) , SIGB< 6, MM ) , TAG( 2» MM ) . I SOT ( MM) ,ISCAT(2,MM), OCC UP ( MM) , 3 SITE (MM)

C WRITE UNWANTED ATOM PARAMETERS ON TOP OF ONE ANOTHER

EQUI VALENCE (ISOT(l),ISCAT(l,l), OCCUP ( 1 ) ,S ITE( 1 ) ,DUMMY( 1 )

)

C READ END OF TAPE F AS WRITTEN BY FINGLS. REFLECTIONS HAVE ALREADY

C BEEN READ.

READ ( NT APE) NATOM , ( C T A G( J , I ) , J= 1 , 2 ) , I S OT ( I ) , ( I SC AT ( J , I ) , J = 1 , 2 ) , OC

2C A ( I ) ,OCCUP( I),SITE( I), ( XYZ( J,I),J=1,3),(BBETA(J,I),J=1,6),I = 1,NAT

3 OM ),(SIGO(I),(SI GX Y Z ( J , I),J=1,3),(SIGB(J,I) ,J=1 ,6) , 1 = 1 , NATOM)

C PACK ATOM NAME INTO ONE LOCATION

DO ID 1=1, NATOM

FLD ( 0,18, NAME ( I ) )=FLD(0 , 1 8 , T AG ( 1 , I )

)

FLD( 18, 1 8, NAME ( I ) ) =FLD (0,18,TAG(2,I))

1C CONTINUE

WRITE (NOUT, 20) NT APE

WRITE (NOUT , 30 ) (NAME(I), OCCA ( I ) . S I GO ( I ) , ( X YZ ( J , I ) , SIGXYZ( J , I ) , J = 1 2,3), ( BBETA ( J, I ) , SI GB ( J . I ) » J= 1 . 6 ) , I = 1 .NATOM)

RETURN

C

C

C

20 FORMAT ( 36 H 1 ATOMIC PARAMETERS AS READ FROM UNIT.I4/99H QUANTITIES

2 ARE E OCCUPANCIES, SIGMA OF OCCUPANCIES, X, SIG(X), Y, SIG(Y), Z,

3S I G ( Z ) ON FIRST LINE./55H AND BETAS, EACH FOLLOWED BY ITS SIGMA,

4 ON SECOND L INE/55H ORDER OF BETAS IS (1,1) (2,2) (3,3) (1.2) (1,3

5) (2,3)/)

30 FORMAT ( 2X A6 , 2 X 2F9 . 6 , 6F l 0 . 6/8 X 1 2F 1 0 . 7

)

C

END

109
SUBROUTINE INFOFC ( NT AP E * JY

,

JZ . NFTEST , I JM

I

1 N , I J M AX , LOGFC , DM , NN ) 2

3

4

5

6

7

8

9

10

1 1 1 2

1 3

1 4

1 5

1 6

17

18

1 9

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

C

Cc

c

cc

c

cc

c

c

cc

cc

cc

c

c

c

c

cc

cc

cc

c

c

cc

cc

cc

c

10 c

c

c

c

c20

**** ******** ** * ***   ** ** ** **************************************

ROUTINE TO READ IN FOBS AND FCALC V ALUES, CALC U LATE THE CORRESPONDING

DM, WRITE OUT FOBS AND F C ALC IN PAIRS ON DRUM.

PRESENTLY SET UP FOR READING THE X-RAY 70 SYSTEM BINARY DATA FILE.

x-ray 70 system of computer programs, edited by j ,m .stewapt, un iversity

OF MARYLAND, COLLEGE PARK, MARYLAND . 20742.

NT APE IS LOGICAL UNIT FOR INPUT FILE.

JY IS COUNT OF REFLECTIONS READ.

JZ IS POINTER TO SPACE FOBS AND FCALC IN PAIRS.

LOGFC ARRAY HOLDS LOCATIONS OF FOBS AND FCALC ON DRUM.

REPLACE STATEMENTS BETWEEN MO CONTINUE TO 20 CONTINUE* AND *30

CONTINUE TO 50 CONTINUE* FOR YOUR OWN READ IN

COMMON /IO/ I N, NOUT , NDRUM

COMMON /H/ DUMMY (501)

DIMENSION DM ( NN ) »LOGFC( NN) , I J ( 3)

DIMENSION IOBUF (50 0 ) ,FOBUFC 500

)

EQUIVALENCE < IOBUF( 1 > ,FOBUF( 1 » , DUMMY! 1 )

)

INITIALIZE COUNTERS

JY=Q

JZ= -1

NEXT SECTION CONTAINS READ IN STATEMENTS FOR X-RAY 70 BINARY FILE.

THIS FILE CONTAINS 25 LOGICAL RECORDS, THE FIRST FOUR WORDS OF EACH

RECORD CHARACTERIZE T hE RECORD TYPE.

WORD 1 IS COUNT OF WORDS IN THE 'PHYSICAL* RECORD

FORTRAN WRITE).

(PRODUCED BY THE

WORD 2 NUMBER OF THE PHYSICAL RECORD BELONGING TO

=0 IF LOGICAL RECORD = PHYSICAL RECORD.

LOGICAL RECORD,

WORD 3 IDENTIFICATION OF LOGICAL RECORD TYPE, I.E.

CONSTANTS, R EFLECT I ONS , ETC

.

HI STORY ,CELL

WORD 4 NUMBER OF WORDS PER BLOCK IN A BLOCKED RECORD, E.G.

OF WORDS PER REFLECTION IN LOGICAL RECORD 15.

NUMBER

BEGIN READING THE BINARY FILE,

CONTINUE

READ ( NT APE ) IOBUF 1 ,( IOBUF ( J ), J= 2, I OBUF 1

)

SKIP OVER FIRST 14 LOGICAL RECORDS AS REFLECTION INFORMATION IS IN

THE 15TH.

IF ( IOBUF! 31 «LT. 15 ) GO TO 10

IF ( I0BUF(3) .GT. 15) GO TO 30

INIALIZE UNBLOCKING POINTER, AND GET NUMBER OF WORDS IN BLOCK.

I PT= 4

NWORD= I OBUF ( 4

)

UNBLOCK THE REFLECTION RECORDS FROM THE BUFFER.

CONTINUE

IJ! 1 )= I OBUF ( IPT+1

)

110
58

59

60

61

62

63

64

65

66

67

68

69

70

71

72

73

74

75

76

77

78

79

80

81

82

83

84

85

86

87

88

89

90

91

92

93

94

95

IJ< 2 )= I O BUF ( IPT+2)

IJ ( 3 )= I 0 BUF( IPT + 3)

JCODE=IOBUF( IPT+9)

FOB S=FOB UF ( IPT+ 1 2)

S I G= FOBU E( I PT + 13 1 FCALC=FOBUF< IPT+14)

TBAR=FOBUF( I PT+21)

C TBAR IS CARRIED ALONG FOR CHECKING POSSIBLE EXTINCTION EFFECTS.

C

C STATEMENTS FROM *20 CONTINUE TO 30 CONTINUE' PROVIDE COMPAT AB IL I TY

C WITH OTHER ROUTINES IN STATUS.

C

C INCREMENT REFLECTION COUNT

JY=JY+1

JZ=jZ+2

C DO NOT OVER WRITE ARRAYS

IF (JY.GT.NN) GO TO 30

C STORE LOCATIONS OF FO AND FC ON THE DRUM

LOGFC ( JY )=JZ*1COOOO+JZ+1

C CALCULATE DM VALUE

DM( JY )=( FOBS-FCALC )/SI

G

C WRITE REFLECTION ON DRUM

CALL DRUMRT ( JCODE, I J, FOBS, S IG)

CALL DRUMRT ( JCODE , I J , FCAL C , TBAR

)

IPT=I PT+NWORD

C CHECK IF MORE REFLECTIONS IN RECORD

IF ( IPT. LT . I OBUF1 ) GO TO 20

C CHECK IF MORE PHYSICAL RECORDS IN LOGICAL RECORD 15.

IF ( IOBUFC 2) .NE.0) GO TO 10

30 CONTINUE

C OUTPUT SUMMARY OF REFLECTION READ-IN

WRITE (NOUT.40) JY.NTAPE

RETURN

C

c40 FORMAT ( I 1 0 , 27 H REFLECTIONS READ FROM UNIT, 15)

C

END

\ 111
12

3

4

5

6

7

8

9

SUBROUTINE LSFIT (Y,X,J,ALINE,BlINE,M,II * LOGHKL » K » ACURVE » B CURVE » 2 CCURVE, ACUBIC,BCUBIC,CCUBIC,DCUB!C) ******************************************* ***********************

SUBROUTINE TO CALCULATE EQUATIONS OF LEAST SQUARES LINE, AND

QUADRATIC AND CUBIC CURVES WHICH RELATE VARIABLES X AND Y. LOGHKL

TELLS WHICH REFLECTIONS ARE TO BE INCLUDED IF REFLECTIONS ARE BEING

DONE BY CLASS,

1 c C VARIABLE M WILL BE ZERO IF NO CLASSES.

1

1

COMMON /IO/ IN. NOUT , NDRUM , NT APE A * NT APEB

12 PARAMETER N=50, NNA=1800, NN=7200» NNB=1801, NNC= 3601, NND

13 2 MM=60, MMM=4C0, NBUFF=700, MMMM=32

1 4 DIMENSION Y( NN).» X( NN) , LOGHKL CNN)

1 5 C K =1 WHEN X ARRAY CONTAINS NORMAL QUANTILES, SET FLAG NQ

16 NQ= 0

17 IF (K.EQ.1J NQ=1

1 8 C RESET K TO ZERO SO IT CAN COUNT THE NUMBER OF DM IN THE FIT.

1 9 K=Q

20 c INIT I ALL I ZE

21 SI GY=0

22 SIGYX=n

23 SIGYXX=0

24 S1GYX3=0

25 STGX=0

26 S I GX X=C

27 S I G X3= 0

28 SIGX4=0

29 SI GX5=0

3 n SIGX6=0

31 ACUB I C=0

32 BCUB IC=0

33 CCUB IC = 0

34 DCU B I C=0

35 DO 30 1=1,

J

36 c LEAVE TAILS OUT OF DISTRIBUTION IF X ARRAY CONTAINS NORMAL QU

37 IF C NQ oN E ® 1 I GO TO 10

38 IF (ABS(XCI)), GT .2.0! GO TO 30

30 10 CONT INUE

40 c SKIP TO STATEMENT 20 IF ALL REFLECTIONS ARE TO BE INCLUDED

41 c PULL OUT OCTANT NUMBER, 0 THROUGH 8, 0 =ALL OCTANTS, 1=+++, ETC

42 N0CT=M/1 00

43 c CHECK IF CLASSES OR LAYERS

44 IF( (M-N0CT*1C 05 .LE.O) GO TO 20

45 c 46 c IF ONLY SELECTED REFLECTIONS ARE TO BE INCLUDED, FILTER THEM i 47 NTEST=FLD( I I , I , LOGHKLC I 1 ) 48 IF C NTESTo EQ.O ) GO TO 30

49 20 CONTINUE

50 c CHECK IF REFLECTION BELONGS TO OCTANT IF SPECIFIED

51 IFI NOCT» LE.O ) GO TO 25

52 NTEST =FLD (33,3, LOGHKL I I )

I

53 NTEST = ABSCNTEST1

54 IF( ( NOCT-ll . NE. NTEST) GO TO 30

55 25 CONTINUE

56 C SUM UP QUANTITIES FOR LINE, QUADRATIC AND CUBIC.

57 XX=X(

I

) **2

OUT HERE

112
58

59

60

61

62

63

64

65

66

67

68

69

70

71

72

73

74

75

76

77

78

79

80

81

82

83

84

85

86

87

88

89

90

91

92

93

94

95

96

97

98

99

ICO

101

102

103

104

105

1 06

107

108

109

110

1 1 1

1 1 2

113

114

SIGY=SIGY+Y( I

)

SIGYX=S IGYX+YC I ) *X< I

)

SI GY XX=S I G YX X+Y ( I )*XX

S I GYX3=S IGYX3+Y ( I )*XX*X< I

)

SI G X=S I G X+ X ( I ) SIGXX=SIGXX+XX

SIGX3=SIGX3+XX*X< I

)

SIGX4=SI GX4+XX*XX

SIGX5=S IGX5+XX*XX*X( I

)

SI G X6= SIGX6+I X X* X ( I ) )**2

C COUNT HOW MANY IN SUMMATION

K=K+ 1 30 CONTINUE

C MAKE SURE WE HAVE ENOUGH

IF (K.LT.3) RETURN

DIMENSION A( 6, 6) ,C( 6) , B( 6)

C LS LIN5

C( 1 )=SIGY

C( 2 )=S IGYX

A(1 ,1 )=K

A ( 1 , 2) =SIGX

A( 2, 1 ) =S IGX

A ( 2*2) =S IGXX

CALL SIMLEQ ( 2) A*C*B)

ALINE=B( 1 ) BLI NE=B ( 2)

C SOLVE FOR LS QUADRATIC, REMEMBER WE HAVE SOME QUANTITIES ALREADY

C STORED FOR LINE CALCULATION,

C(3)=SIGYXX

A ( 1 ,3)=S IGXX

A( 2,3) =S I GX3

A ( 3, 1 )=S IGXX

A ( 3 , 2 ) =S I GX3

A ( 3, 3) =SIGX4

CALL SIMLEQ <3,A,C,B>

ACURVE=B(1 ) PCURVE=B( 2)

CCUPVE=B< 3)

IF (K.GT.4) GO TO 40

RETURN

C SOLVE FOR LS CUBIC

40 C( 4 ) =S IGYX3

A(1,4)=SIGX3

A ( 2,4) =S I GX4

A(3,4)=SIGX5

A( 4, 1 ) = S IGX3

A(4,2)=SIGX4

A ( 4, 3)=SIGX5

A(4,4)=SIGX6

CALL SIMLEQ <4,A,C,B)

ACUB IC=B ( 1 ) BCU BI C=B ( 2

)

CCUB I C=8 ( 3

)

DCUB I C= B ( 4

)

RETURN

C

END

113

\ ' •
£•

(jj

12

5

6

7

8

9

10

1 1 1 2

1 3

1 4

1 5

16

SUBROUTINE MAXHKL ( I J * J J . KK , LL

)

C

c ******** ******* ************* ** ************************** ** *** ***^

c

C SUBROUTINE TO OBTAIN THE MAXIMUM AND MINIMUM VALUES AT THE MILLER I! DIMENSION I J ( 3

)

COMMON /D/ MAXH. MI NH »MAXK, MI NK»MAXL,MINL

IF C I J ( J J) -LT .MINH ) MIN H= I J ( J J ) IF ( I J( JJ) .GT.MAXH! maxh=i J I J J

1

IF ( I J (KK) .LT .M INK ) MINK=IJ(KK)

IF ( IJ(KK) ® GT »MAXK ) MAX K= S J { KK ) IF < I J( LL) .LT.MINL ) MINL=SJ(LL)

IF ( I J (LL ) .GT .MAXL ) MAXL= IJ(LL)

RETURN

C

END

114
12

3

4

5

6

7

8

9

10

1 1 12

1 3

1 4

1 5

1 6

17

18

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

n

n SUBROUTINE MTCHEK { NAM EA , L A . N AME B, L B , M ATCH, KOUNT , MM

)

***********:***   ** ***:***** ****** ********* *************** ********

i C

C SUBROUTINE TO CHECK NUMBER OF COMMON ATOMS BETWEEN TWO ATOMIC SETS

C ONE ATOMIC SET MAY HAVE MOREiE.G* HYDROGEN ATOMS, ETC.

COMMON /1 0/ INtNOUT ,NDRUM,NT APEA, NT APEB

C NAMEA IS FIRST SET OF ATOM NAMES, NAmEB IS SECOND SET, MATCH IS A C WORD INTO WHICH THE SEQUENCE NUMBERS OF CORRESPONDING ATOMS ARE

C PACKED.

DIMENSION NAME A ( MM ), NAMEB ( MM ), MATCH (MM

)

KOUNT= 0

C DO OVER BOTH ATOM SETS, LA ATOMS IN SET 1, LB IN SET 2.

DO 20 1=1, LA

DO 10 J= 1 , LB

IF ( NAME A( I ) .NE. NAMEB( J ) ) GO TO 10

KOUNT= KO UNT+

C STORE CATALPGUE NUMBER OF MATCHED ATOMS

C I ISPOSITION OF FIRST ATOM IN FIRST SET, J IS POSITION OF SEC ND ATO

C SECOND SET

MATCH( KOUNT )=I *100+J

10 CONTINUE

20 CONTINUE

WRITE ( NOUT , 30 ) KOUNT

RETURN

Cc

c30 FORMAT (10H THERE ARE,I4,46H ATOMS IN COMMON BETWEEN ATOMIC SETS 1 2 AND 2/61H IF THIS IS LESS THAN 3 ATOMS , A BANDON PARAMETER CHECKI

3NG PART)

C

END

115
123

4

5

6

7

8

9

1 0

1 1

1 2

1 3

FUNCTION NF ACT (N)

C

C ******************** *********************** ***********************

c

C INTEGER FUNCTION TO CALCULATE N FACTORIAL.

NT= N —

1

NF AC T=N

DO 10 J= 1 NT

NF AC T=NF ACT* ( N-J)

10 CONTINUE

RETURN

C

END

116
123

4

5

6

7

8

9

10

1 1

1 2

1 3

14

1 5

16

IT

18

19

20

21

22

23

24

25

26

27

28

29

30

31

FUNCTION PINV <Q)

Cc ****************** ****** *^ **************** * ***********************

cC INVERSE PROBABILITY FUNCTION BY R.J.ARMS 2/27/67

C Q BETWEEN 0. AND 1.

C ERROR LESS THAN 4.5E-4 . SEE NBS HANDBOOK 26.2.23 . DATA C0.C1.C2. D1,D2,D3.E1,E2,K /2.S15517 , .802853. .010328,

2 1*432788, .1 89 269 , .0 0 1 30 8 , 1 . E-06 , 1 . E- 3 5, 4/

IF (Q.LT. 1 .+E1 .AND.Q.GT.-El ) GO TO 10

P IN V=0 . IF (K.EQ.O) RETURN

K=K- 1 WRITE (6,30) Q

RETURN

C

10 F=l.

R=AMA X 1 ( 0. , AM I N 1 ( 1 • , Q)

)

IF (Q.LT. .5) GO TO 20

F=-l . R=1 . -R

20 a=SQRT (-2. *ALOG( AMAXK E 2. R ) ) ) PINV=F*( A- ( (C2*A+C 1 ) *A+C0 )/ ( ( ( D3TA + D2) *A +D 1 )*A+1. ) ) RETURN

C

cc30 FORMAT ( 1 HO , 48H I NV ER SE PRO 3 I L I TY FUNCTION HAS INVALID A R G U M E N T= E 1

6

2.8, /IX, 1 8H VALUE SET TO ZERO.)

C

END

117
12

3

4

5

6

7

8

9

10

1 1

1 2 13

1 4

15

16

1 7

18

19

20

2 1 22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

SUBROUTINE PLOTEM ( J C, X , DM , D MS M AL , X SMALL , M , I I , L OGHKL , NN , N

)

C

C **************************** 44 4* *4*444**44*** *********** **********

C

C SUBROUTINE TO CALCULTE LEAST SQUARES LINE AND LEAST SQUARES QUADRATIX

C FOR DM AGAINST X AND THEN PLOT DM AGAINST X. LEAST-SQUARES

c line and cubic curve included in 3rd. plot.

COMMON /IO/ IN* NOUT , NDRUM , NT APE A * NTAPEB

C ROUTINE SETS UP PLOT TITLES, CALLS ROUTINES TO CALCULATE LEAST—SQUARE?

C LINE* QUADRATIC AND CUBIC.

COMMON /TITLE/ AT I TLEC 14 5 C DM AND X ARE ALL QUANTITIES. DMSMAL AND XS M ALL ARE SUBSETS OF DM AND

C X SELECTED FOR PLOTS. NUMBER IS NUMBER OF POINTS IN PLOT.

DI MENS ION DM( NN) ,X (NN) . DMSMALC N, 5 ) , X SMALL ( N , 5) . NUMBER < 5)

DIMENSION ALAY(lO)/6H 1 . 6H 2 ,6H 3 ,6H 4 , 6H 5 , 16H 6 , 6H 7 , 6H 8 ,6H 9 * 6H 10 /

DIMENSION OCT (9)/6H +++ »6H -&• + - , 6H +-+ . 6H + , 6H -++ . 6H

1-+- . 6H —+ . 6H —-- ,6H ALL /

C L OGHKL TELLS REFLECTION CLASSES. IT IS NEEDED FOR TRANSFER TO FITCHf

C A ND L SF I T DIMENSION LOGHKL(NN)

K=0

C GET EQUATIONS FOR LEAST-SQUARES LINE, QUADRATIC, CUBIC.

C THESE ARE DM = ALINE + BLINE4X

C DM = ACURVE + BCURVE4X + CCURV E4X**2

C AND DM = ACUBIC + BCUB IC* X + CCUBIC*X**2 + DCUBIC*X**3

C

C IN THE ARGUMENT LIST, M IS NON-ZERO IF SOME SUBSET OF POINTS IS TO

C BE SOUGHT AND II IS THE NUMBER OF THE BIT IN THE L OGHKL WORD TO

C EXAMINE ( 0= IGNORE , 1= INCLUDE) IF REFLECTION CLASSES ARE BEING PICKED

C SIGNAL =LSFIT= AND =FITCHK= IF X ARRAY HOLDS NORMAL QUANTILES.

IF ( ATI TL F( 12) . EQ . 6 HST NOR . OR . C AT I TL E< 1 2 ) . EQ . 6HST 1/2))K = 1 J=K

C

C K WILL CONTAIN NUMBER OF POINTS IN THE SET AFTER WEEDING OUT IS

C DONE IN =LSFIT=.

C JC I S TOTAL NUMBER OF POINTS CONSIDERED.

CALL LSF IT ( DM, X, JC * AL I NE , BL I NE , M, I I , LOGHKL • K , ACURVE , BCURV E, C CUR V! 2 , A CU BI C, BCUBIC, C CU B I C , D CU B I C ) C MAKE SURE WE HAVE ENOUGH POINTS TO PROCEED.

IF (K.GT.4) GO TO 10

WRITE C NOUT, 240 I K

RETURN

C

10 L=K/N+1

C CALCULATE GOODNESS OF FITS AND FILL OUT PLOTTING ARRAYS DMSMAL AND

C X SMALL • N IS DIMENSION OF PLOTTING ARRAYS, EVERY L'TH POINT WILL BE

C PLOTTED, J IS NUMBER OF POINTS IN FINAL PLOT, SUMLN, SUMCV AND SUMCUf

C ARE GOODNESS OF FITS FOR LS LINE, QUADRATIC AND CUBIC,

C RESPECTIVELY.

CALL FIT CHK ( DM, X, JC, DMSMAL, XSM ALL , L , K , ALI NE , BL I NE , ACURVE , BCURVE,<

2CURVE, SUMLN , SUMCV, M, II, LOGHKL, J, ACUBIC, BCU B IC * C CUB I C , DC UB I C , SUMCUf

3) c tell plots how many points to plot for line, curve and actual data.

DO 20 1=1,3

number c I )=

j

20 CONTINUE

118
58

59

60

6 1 62

63

64

65

66

67

68

3 69

70

71

' 72

73

74

75

76

77

78

- 79

80

8 1 82

83

84

85

86

87

88

89

E 90

91

92

93

94

95

96

97

98

99

100

10 1

1 02

103

104

105

106

1 C 7

. 108

109

1 1 0

1 1 1

1 1 2

113

1 1 4

115

NPT = 3

IF( ATITLEC 12) .EQ.6HST 1/2) GO TO 200

IF ( ATITLEC 12) .EQ.6HST NOR . 0 R . L . LT . 6 ) GO TO 30

C SET UP FOR PLOTS INVOLVING DM SCATTER

NPT=5

NUMBER! 4)=J

NUMBER ( 5 ) = J 30 CONTINUE

C THE NEXT PART GETS THE APPROPRIATE TITLE.

IF (M.LE.0) GO TO 180

C GET OCTANT NUMBER AND CLASS OR LAYER TYPE,(1=STD CLASSES , 2=H, 3=K , 4=L

>

NOCT=M/l 00

layt=m-noct*ioo

IF ( NOCT • EQ • 0 ) NOCT=9

C SET UP TITLE

ATI TLE<2)=6H

AT ITLE( 3)=6H REFL

ATI TLE(4 )=6HECT ION

ATITLEC 5)=6HS

ATITLEC 6 ) = 6H

ATI TLE (7 )=6H

IF ( I I . GT. 10 ) GO TO 150

NTEST= I I+l

GO TO (40,50.60,70.80, 90,100,110,120,130, 140), NTEST

C FOLLOWING SECTION FOR REFLECTION CLASSES

40 ATITLE( 1 )=6HH,0 ,0

GO TO 190

C

50 ATITLEC 1 )=6HC,K,C

GO TO 190

C

6C ATITLE( 1 )=6H0 ,0 ,L

GO TO 190

C

70 ATITLEC 1 )=6HH,K,0

GO TO 190

C

80 ATI TLE C 1 )=6HH,0 ,

L

GO TO 190

C

90 ATITLE ( 1 )=6H0,K,L

GO TO 190

C

100 ATITLEC 1 )=6HH,K,L

ATITLEC 5 ) = 6HH, K, AN

ATI TLE (6 )=6HD L NO

AT ITLEC 7)=6HT = 0.

GO TO 190

C

110 ATITLEC 1 )=6HH, H,H

GO TO 190

C

120 ATITLEC 1 ) = 6HH. H,L

GO TO 190

C

130 ATITLEC 1 )=6HH,K,K

GO TO 190

C

119
1 16

1 1 7

1 18

1 19

12 ° 121

122

123

1 24

125

126

127

128

1 29

130

131

132

133

134

1 35

136

137

138

1 39

140

1 4 1 142

143

144

145

146

147

148

149

1 50

151

152

153

154

1 55

1 56

157

158

1 59

160

161

162

163

164

165

166

167

168

1 69

170

17 1 172

173

140 ATITLEC 1 )=6HH,K,H

GO TO 190

C TITLES FOR LAYERS FOLLOW.

150 CONTINUE

NTEST= I I -1 0

IFC LAYT. NE. 4) GO TO 160

ATITlEI 1 ) = 6H H K ATI TLE ( 2 ) = AL AY { NTEST

)

GO TO 190

C

160 CONTINUE

C TEST IF K LAYERS WANTED

IF(LAYT .NE. 3) GO TO 170

ATI TLE ( 1 1=6 H H

AT I TLE ( 2)=ALAY( NTEST!

ATITLEC 3 )=6HL REFL

GO TO 190

C

170 CONTINUE

C DO H LAYERS HERE

IFILAYT.NE. 2) GO TO 190

ATITLEI 1 ) = AL AY { NTEST

)

ATITLEC 2)=6HK L

GO TO 190

C

180 ATITLEI 1)=6H ALL R

AT I TLE ( 2 )=6HEFLECT

ATI TLE ( 3 >=6HI ONS I AT I TLE ( 4 ) = 6HNCLUDE

ATITLE(5)=6HD

AT I TLE ( 6 ) =6H

AT ITLE( 7 )=6H

GO TO 20 0

C WRITE OUT HEADING

190 WRITE(NOUT,260)OCTCNQCT)

C WRITE OUT HEADING

200 WRITE (NOUT. 270) ATITLE

C WRITE OUT EQUATIONS OF LINE ETC. AND GOODNESS OF FITS

WR ITE ( NOUT, 280) L .K.J.ALI NE » BL l NE . ACURVE, B CURV E, C CURVE, ACUBI C. BCU

2BIC.CCUBIC .DCUBIC. SUMLN, SUMCV, SUMCUB

C PLOT X SMALL ARRAYS AGAINST DM SMALL ARRAYS

CALL PLOTS ( NPT, XSMALL , DMSMAL, NUMBER, N)

C TELL USER WHAT SYMBOLS IN PLOTS MEAN.

WRITE C NOUT . 290 ) IF (ATITLEC 12) .EG. 6HST NOR . OR. L . LT . 6 ) GO TO 210

WRITE (NOUT, 220) L

2 1 C CONTINUE

WRITE (NOUT, 230)

C WRITE OUT NUMBERS USED IN PLOT.

WRITE ( NOUT, 300) { ( I , DMSMAL ( I , 1 ) , DMSMAL ( I , 2 ) , DMSMAL ( 1,3) , XSMALL ( I

,

21 ) ) , 1= 1, J)

RETURN

C

C

c220 FORMAT ( 1H , 64H SYMBOLS L, U ARE FOR AVERAGE - AND + DEVIATIONS (A

2 —B VALUES FOR, I 5 » 2 5H POINTS IN THE INTERVAL.))

230 FORMAT ( 1H ,52H POINT ABC D

120
DATA IN DATA SET NOT ENOUGH

174

175

176

177

178

179

18C

18 1 182

183

184

185

186

187

188

189

190

191

192

2 ) 240 FORMAT (I4.59H

2TO PLOT)

260 FORMAT! 1H1 , A6 ,9H0CTANT< S),)

270 FORMAT ( 1H+ , 26X, 1 4A6

)

280 FORMAT (6H EVERY, 13, 18H TH POINT IN PLOTS, 14, 35H POINTS IN THIS SE

2T , 1 4 , 20 H POINTS IN THIS PL0T/32H LEAST SQUARES LI

3NE IS DM=,F9.3,3H + , F9 . 3 , 2H* X , 34 H LEAST SQUARES QUADRATIC

4 IS DM=,F9.3,3H + ,F9.3,4H *X+,F9.3,6H *X**2/27H LEAST SQUARES CUB

5 I C IS DM=,F9.3,3H + .F9.3.6H *X + ,F9.3,9H* X**2 + , F9 . 3 , 5 H* X ** 3/

69H GOODNESS OF FIT IS,F9.3,14H FOR LS LINE, ,F9.3,18H FOR LS QUAD

7RATIC,5H AND ,F9.3,41H FOR CUBIC FIT IF MORE THAN 4 PTS IN PLOT)

290 FORMAT ( 1H ,107HA = DM( I) (ORDINATE), SYMBOL IS . B = DM( I

) — L. S

2. LINE. SYMBOL IS 7 C = DM(I) - L. S. CUBIC, SYMBOL IS +/74H D = 3ABSCISSA VALUES (DSTAR, LOG F**2/S I N< 2*THET A ) . OR QUANTILE

4 ) 300 FORMAT ( I 1 0 . 3F 1 0 • 3 , El 5 . 3

C

END

121
1 SUBROUTINE PLOTS ( NARGS, X,Y, NRMX.NROW)

23

45

6

7

8

P

1 0

1

1

1 2

13

14

1 5

16

1 7

1 8

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

4 1 42

43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

C

C

c

ccc

c

c

c

c

cc

c

c

c

c

c

c

c

cc

c

c

cc

c

c

cc

c

c

c

c

c

#*#*** ** ********** **** * ********** 4  

4 ** **** ****** *** *** ********

PLOTS UP TO 5 PAIRS OF COLUMNS. THE NOTATION IS AS FOLLOWS.

X - A MATRIX HAVING UP TO 5 COLUMNS AND NROW ROWS

Y - A MATR I

X

HAVING UP TO 5 COLUMNS AND NROW ROWS

THE FIRST COLUMN OF Y IS PLOTTED VERSUS THE FIRST COLUMN OF X SECOND Y SECOND X ETC

NROW - THE NUMBER OF ROWS IN X AND Y IN THEIR DIMENSION STATEMENTS

IN THE MAIN PROGRAM. E.G. IF X AND Y APE DIMENSIONED

AS X( 130,4) , Y( 130,4) IN THE MAIN PROGRAM THEN SET NROW=

NARGS - THE NUMBER OF PAIRS OF COLUMNS TO BE PLOTTED.

NRMX - A VECTOR CONTAINING THE LENGTHS OF THE DATA TO BE

PLOTTED IN THE SEVERAL COLUMN PAIRS. THAT IS,

NRMX ( 1 ) IS NUMBER OF POINTS TO BE PLOTTED IN FIRST COLUMN

NRMX! 2)

NRMX I 3)

ETC.

THE SYMBOLS USED ARE FIRST SECOND

© 4

IF MORE THAN ONE POINT FALLS ON THE

THE NUMBER IS PRINTED.

LIMITS FOR THE PLOT ARE FIGURED OUT

AND LARGEST DATA POINTS.

SECOND

THIRD

THIRD FOURTH FIFTH COLUMN PA

+ L U

SAME POSITION A TALLEY IS KEPT AN

AUTOMATICALLY BASED ON THE SMALLE

C4444

C  DOES NOT CALL A NEW PAGE.

C

C** 4* WRITTEN BY S. PEAVY 5/1/67

C COMMENTS BY B.L. JOINER 6/4/69

C

c

DI MENSION NRMX (5) , X C NR OW , 5 ) ,Y(NROW,5)

DIMENSION PRINT! 10 1 ) , XP( 6) , BOOLI 5) , IDGT( 9)

INTEGER PRINTsBLANK

EQUIVALENCE (XO.XMIN), (Xl.XMAX), ( YO , YM IN ) , ( Y 1 , YMAX

)

INTEGER BOOL

data BOOLI 1 ) , BOOLI 2) , BOOL ( 3 ) , BOOL (4 ) , BOOL ( 5 ) / 1 H. , 1 H* , 1 H+ , 1 HL

,

2 1 HU/ , COL X /6HCOLUMN/ .BLANK /1H /

DATA IDGTI 1 ), IDGTI 2 ) ,IDGTI 3),IDGTI4),!DGTI5),IDGTI6), IDGTI 7 )

,

2 IDGTI 8) . IDGTI9 ) / 1 H2 . 1 H3 , 1 H4 , 1 H5, 1 H6, 1 H7 , 1 H8. 1 H9 , 1 HX/

IPR I NT =6

X1 = XI 1 , 1 ) X0 = X 1 DO 10 ix=l, NARGS

K2=NRMXI IX

)

DO 10 1=1, K2

122
58

59

60

61

62

63

64

65

66

6 7

68

69

70

71

72

73

74

75

76

77

78

79

80

8 1 82

83

84

85

86

8^

88

89

90

9 1 92

93

94

95

96

97

98

99

ICO

10 1

1 C 2 103

104

105

106

107

1 08

1 C 9

110

1 1 1

1 1 2

113

1 14

115

IF ( XI .LT.XJ I ,1 X) ) X1=X(I.IX)

IF ( XO .GT.X( I, IX) ) X0=X(I,IX)

10 CONTINUE

Y1=Y< 1,1)

YO = Yl

DO 20 J = 1 , N ARGS

K4=NRMX( J)

DO 20 1=1, K4

IF < Y1 ,LT. Y< I , J) ) Y1=Y(I,J)

20 IF (YO.GT.YC I, J) ) YO=Y(I,J)

WRITE ( I PR I NT , 250

)

C**** DETERMINE X AND Y INCREMENTS FOR PLOT

YDELTA=( YMAX-YMIN ) /50

XDELTA= ( XMAX — XMIN) /100

YL=YM AX-YDELTA/2

YT=Y MAX

C  THE I LOOP CONTROLS THE 5 DIVISIONS OF THE Y ORDINATE

DO 190 1=1,6

L= 1 C ** ** THE J LOOP IS FOR EACH LINE OF PRINT WITHIN THE DIVISIONS

DO 190 J=l,10

C  BLANK OUT PRINT BUFFER LINE.

DO 30 K= 1,101

30 PRINT(K) = BLANK

C **  * * THE KK INDEX IS FOR EACH CURVE. KK LESS THAN 6.

DO 150 KK=1,NARGS

K4 = NR MX ( KK

)

K5= 1 c**** this determines if y(k) value is on the present print line

DO 1 40 K=1 , K4

IF (Y(K.KK)-YT) 40,40.140

40 IF (Y(K,KK) —YL ) 140,140,50

C *** * YES. Y ( K ) BELONGS ON THIS PRINT LINE

C**** THEREFOPE DETERMIND WHERE ALL THE X(K5) fall ON THE X-AXIS

50 XL= XMIN

XT= X MI N+ XDELT A/2 • DO 130 KA= 1,101

IF ( X( K5 ,KK) -XL) 120,60,60

60 IF ( X( K5,KK )— XT) 70,120,120

70 IF ( PRINT(KA)— BLANK) 90,80,90

80 PRl NT( KA ) =BOOL( KK)

GO TO 140

C**4* IF MORE THEN ONE POINT FALLS ON THE PRINT POSITION, TALLY THE NO.

C**** OF POINTS.

90 DO 100 K KK= 1,9

IF ( PR I NT ( K A

) — I DGT ( KKK ) ) 100,110,100

100 CONTINUE

PRINT(KA)=IDGT( 1

)

GO TO 140

1 10 IF (PRINT(KA) .NE.IDGT(9)) PR I NT ( K A ) = I DGT ( K K K+ 1

)

GO TO 140

120 XL= XT

130 XT= XT+ XD ELT A 140 K5=K5+1

150 CONTINUE

YT=YL

YL=YL—YDELT A GO TO (160,180), L

123
1 16

1 17

1 1 8

1 19

120

121

122

123

1 24

125

1 26

127

128

1 29

130

13 1 132

133

1 34

135

1 36

1 37

138

139

140

141

1 42

143

144

160 IF ( 1 — 5) 170*170*2 00

17C L=2

C**** THIS PATH IS EXECUTED ONCE IN EVERY DIVISION OF THE Y-AXIS.

C**** TENTH LINE, STARTING WITH ZERO LINE

YP=YT -9-YDELT A/ 2

.

WRITE C I PR I NT « 220 ) YP, PRINT

GO TO 190

180 WRITE C I PR 1 NT * 2 30 ) PRINT

C**** PRINTS LINE

190 CONTINUE

200 WRITE ( I PR I NT * 2 20 1 Y MI N* PR I NT

C **** LAST LINE OF PRINT OUT PLUS X VALUES ALONG X-AXIS.

WRITE ( I PR I NT « 2 50

)

XPC 1 1=XMIN

XPC 6 )= XMAX

XR=20 . *XDELTA

DO 210 1=2,5

210 XPC I ) = XPC 1-1 )+XR

WRITE C I PR I NT , 240 1 XP

RETURN

C

c

c220 FORMAT { 1 X , E 1 2 . 5 » 1 H+, 1 0 1 A 1 , 1 H+

)

230 FORMAT ( 1 3 X , 1 H— , 1 0 1 A 1 *1

H

—

)

240 FORMAT C 6C 7X » E 1 3. 5 ?

1

250 FORMAT C 14 X , 1 H+ , 1 0 « 1 0 H

C

EVERY

END
UUUUUU

U

V

uu

1 SUBROUTINE POSITN { NT APE * NF I LE

**************************** **************************** ** ******** SUBROUTINE TO POSITION ANY TAPE =NTAPE= WITH FILE =NF I LE = READY

SEPARATE ROUTINE ALLOWS EASY MODIFICATION FOP LOCAL COMPUTER CONVENT IONS

COMMON / IO/ IN.NOUT

MAYBE THE USER DOESN'T WANT TO POSITION HIS TAPE. IF SO, RETURN.

IF (NFILE.EQ.O) RETURN

PUT IN THE NECESSARY STATEMENTS HERE.

WRITE (NOUT.IO) NTAPE.NFILE

RETURN

16

IT

1 8

10

C

FORMAT ( 5H TAPE,l4»19H POSITIONED TO FILE, 14)

END

125

l
uu

1 SUBROUTINE PRIME

4’5

6

7

8

9

1C

1 1

1 2

1 3

1 4

1 5

16

cCC ALLOWS USER TO SET UNITS NUMBERS FOR ADAPTION TO LOCAL COMPUTER

C CONVENTIONS.

COMMON / 10/ IN* NOUT * NDRUM. NT APE A . NT APEB

C CAPO READER

IN=5

C PRINTED OUTPUT

NOUT =6

C NORM IS DRUM UNIT NUMBER

NDRUM=27

RETURN

C

END

126

_
123

4

5

6

7

8

9

1C

1 1

1 2

1 3

1 4

1 5

1 6

17

1 8

19

20

21

22

23

24

25

26

2?

28

29

30

SUBROUTINE PROB (DM.X.JC.NN)

Cc ************************:£*** *********************************** ***

c DIMENSION DM(NN) ,X(NN)

C WE NOW HAVE JC VALUES OF DM ARRANGED IN INCREASING ORDER . COMPARE

C DISTRIBUTION OF DM,S WITH NORMAL PROBABILITY FUNCTION.

C PROBABILITY OF INDIVIDUAL DM

C DM VALUES GO ON BOTH SIDES OF ZERO

C FIND HALFWAY

NTEST=JC/2

FJC= JC

DO 10 1=1, JC

C MAKE PROBABILITY SUITABLE FOR HALF NORMAL PLOT PROGRAM P INV ( Q ) , SUPPL IC ED BY R.J.ARMS OF NBS.

FI=I

C CALCULATE PROBABILITY IN FLOATING POINT ARITHEMATIC

PRO B= A BS ( ( F J C— 2.TFI + 1 )/(FJC)

)

C MAKE PROB SUITABLE FOR USE IN ARM'S ROUTINE

IF ( PROB ,GT . 1 .0 ) PR08=1.

Q=< 1 .-PROB >/2

.

C OBTAIN CORRESPONDING NORMAL QUANTILE VIA PINV(Q).

X ( I

) =P I NV { Q

)

C CONVERT TO FULL NORMAL PLOT VALUE.

C END TEST GOVERNS SIGN OF X VALUE

IF (I.LE.NTEST) X(I)=-X(I>

1C CONTINUE

RETURN

C

END

12 7
12

3

4

5

6

7

8

9

in

1

1

1 2

1 3

14

15

16

1 7

1 8

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

4 1

SUBROUTINE PRY NT ( DM ,LOGFC .LOGhKL, JC, J J A, KK A,LL A, J JB , KKB ,LLB

,

2 NPRINT, X, NFTEST.NN)

C

C ******************************************* ***********************

cC SUBROUTINE TO PRINT ALL DM VALUES WITH ATTACHED H , K , L AND X, OR

C PRINT THOSE VALUES WITH/DM/ GREATER THAN 2.

COMMON /IO/ IN, NOUT * NDRUM, NT APE A » NT APEB

DIMENSION DM ( NN ) , LOGFC (NN)»X(NN) .LOGHKL(NN)

NOUT=2 0

WRITE C NOUT. 30)

IF ( NFTEST .EQ. 2 .OR .NFTEST.EQ .6) WRITE (NOUT, 20)

DO 10 1=1, JC

C SELECT DM VALUES AS REQUESTED IF NPRI NT=2

IF ( NPRINT. EQ.2. AND. ABS ( DM< I ) ) ,LT. 2. ) GO TO 10

C GET REFLECTION LOCATIONS ON DRUM

J=LOGFC( I )/l 0C000

K=LOGFC( I) -1 00 COO*

J

C READ REFLECTIONS OFF DRUM.

CALL DRUMRD ( JHKLA »FA» SGA.NA.J

)

CALL DRUMRD ( JHKLB , F B, S GB , NB , K

)

C FIND HKL FOR EACH REFLECTION.

CALL HKL GEN ( J HKL A , J A, K A , L A, J J A , KK A , LL A

)

CALL HKLGEN ( JHKLB , JB , KB ,L B , JJB , KKB , LLB

)

C WRITE OUT QUANTITIES.

WRITE ( NOUT, 40) J A , K A , L A , J B . KB , L B, N A , NB , F A , FB , S GA , SGB , DM ( I ) , X

(

I ) ,

J

2 , K 10 CONTINUE

NOUT=6

END FILE 20

RETURN

C

Cc20 FORMAT C1H0.63HF(2) IS FCALC AND SQ(2) GIVES EXTINCTION FACTOR, OR

2 • TBAR.ETC.)

30 FORMAT (90H1 HI K1 Ll H2 K2 L2 LT F(l) F(2) SG(1)

2SG( 2 ) DM(I) X(I) LIST PLACES//)

40 FORMAT ( 6 1 4 , 2 I 2, 2F 8 . 3, 2F8 . 4 , 2F 8 . 3 , 2 I 6)

C

END

128
123

4

5

6

7

8

9

1 0

1 1

1 2

13

14

1 5

16

17

18

19

20

2 1 22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

4 1 42

43

44

45

46

47

48

49

50

51

SUBROUTINE PRYNTT ( J C , D M, X , L OG , NPR I NT , K I ND , OC A , OCB , S I GO A , S I GOB

,

2 XYZA,XYZ8,SIGXA,SIGXB,BA,BB,SIGBA,SIGBB,NAMEA,NAMEB,NN,MM)

c

c ******************************************* ***********************

c

C SUBROUTINE PRINTS DM VALUES TOGETHER WITH

C CORRESPONDING PARAMETERS AND THEIR NAMES.

COMMON /IO/ IN.NOUT

DIME NS ION DM( NN) , X( NN) , LOG ( NN ) , OCA C MM) , 0CB( MM),SIGOA(MM),SIGOB(MM)

2 ,XYZA( 3,MM) , XYZB(3, MM ) , SIGXA (3, MM) , SIGXB( 3, MM ) ,BA( 6, MM) . BB (6,MM)

,

3 SI GBA ( 6 iMM ) « SI GBB (6* MM ) .NAME A ( MM) . NAMEB ( MM

)

DIMENSION PAR AM ( 10 ) DATA PAR AMI 1 ) , PARAMI 2) . PAR AMI 3) * PAR AM ( 4) ,PARAM( 5) ,PARAM< 6) . 2 PARAMI7) ,PARAM( 8) .PARAMC9) ,PARAM(10) /6 HOCCUP ,6H X ,6H Y , 3 6H Z »6HB( 1. 1) . 6HB (2»2)»6HB(3»3)» 6HB ( 1. 2) »6HB( 1 .3) ,6HB( 2.3)/

WRITE (NOUT.60)

DO 50 1=1. JC

IF ( NPRINT.EQ.2 . AND. ABS (DM( I ) ) .LT. 2. ) GO TO 50

C FIND NUMBERS OF ATOMS ASSOCIATED WITH EACH DM VALUE

MATCH=LOG( I ) / 1 0 C C FIND WHICH PARAMETER OF THESE ATOMS

NPARA=LOG( I

) -MATCH* 100

NA=MATCH/1 00

NB=MATCH-NA*1 CO

C WRITE OUT PART OF LINE

WRITE (NOUT.70) DM ( I ) , X (

) , N AME A ( N A ) , N A MEB ( NB

)

C WRITE OUT ACTUAL PARAMETERS. BASED ON VALUE OF KIND =1 FOR OCCUPANCY

C 2 FOR POSITIONAL. 3 FOR THERMAL

GO TO (10,20,30), KIND

10 WRITE ( N OUT , 80 ) OC A ( NA ) , OC B ( NB ) , S I GOA ( NA ) , S I GOB ( N B

)

GO TO 40

20 NP=NPARA-1

WR ITE ( NOUT, 80 ) XYZ A ( NP , N A ) , X Y ZB ( NP , NB ) , S I G X A ( NP , N A ) , S I GXB ( NP , NB

)

GO TO 40

30 NP=NPARA—

WRITE (NOUT, 80) BA ( NP , N A ) , BB ( NP , NB ) , Si G B A ( N P , N A ) , S I G BB ( NP, NB

)

C WRITE OUT HOLLERITH LITERAL TO TELL EXPLICITLY TO THE OUTSIDE WORLD

C WHICH PARAMETER THAT DM IS FOR.

40 WRITE (NOUT, 90 PARAM(NPARA)

50 CONTINUE

RETURN

C

C

c60 FORMAT

2 AND ARD

80 FORMAT

70 FORMAT

90 FORMAT

C

END

( 7X , 79 H DM X ATOM 1 ATOM 2 PARAMETERS ST

DEVIATIONS VARIABLE)

( 1H+.36X, 2F 8.5, 5 X , 2F 8.5)

( 7X , 2F7. 3.2 (2XA6 )

)

( 1H+,80X,A6)

129
123

4

5

6

7

8

9

1 0

1 1

1 2

1 3

1 4

15

16

17

1 8

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

SUBROUTINE REREAD ( JY i L OGF C»DM»JJ,KK,LL»NT APE, J Z » SKALE)

C

c ******************************************************************

C IN ORDER ON UNIT NR BY SCALE ROUTINE IN STATUS PROGRAM.

C

C SUBROUTINE TO REREAD REFLECTION INFORMATION PREVIOUSLY WRITTEN

C WRITING out at that TIME AND READING IN FROM UNIT NT APE HERE SAVES

C THE TIME REQUIRED FOR THE SEARCH FOR MATCHES IN 2 LISTS OF REFLECTION

C NTAPEA AND NTAPEB NOT IN COMMON HERE BECAUSE WE WANT TO USE SAME

C SUBROUTINE FOR BOTH TAPES. THEY ARE TRANSFERRED IN VIA NT APE ARGUMENT

C IN CALL STATEMENT.

COMMON /TO/ IN.NOUT.NDRUM

PARAMETER N = 50, NNA=1800, NN=7200, NNB=18C1, NNC=36C 1 , NND=54C1,

2 MM=60, MMM=4C0, NBUFF=490

C DUMMY IS USED TO OVERLAP IN CORE SOME QUANTITIES WHICH ARE NOT

C NEEDED AT SAME TIME.

COMMON /H/ DUMMY I N N ) C I J ( 3 ) IS MILLER INDICES. DM IS ORDINATE IN PLOTS, LOGFC IS WORD

C CONTAINING PACKED SEQUENCE NUMBERS ON DRUM OF CORRESPONDING

C REFLECTIONS. BUFF AND I BUFF ARE INPUT BUFFERS FOR THE TAPE READING

DIMENSION BUFF(NBUFF) , I B UF F(NBUFF) , IJI3) .DM(NN) , LOGFC INN)

C EQUIVALENCES TO CONSERVE STORAGE

EQUI VALENCE IBUFFC 1 ) , I BUFF C 1 ) , DUMMY ( 1 ) ) WRITE (NOUT.60) SKALE

SS=SKALE**2

jx=o

NT= 0

SUM=0

.

SUM M = Q

•

C BEGIN BY READING REFLECTIONS INTO BUFFER

1 0 NBR=0

READ C NT APE ) BUFF

20 CONTINUE

I J ( J J ) =1 BUFF I NBR+1 ) IF C I J C J J ) .EQ.99) GO TO 30

IF CJX.GT.NN) GO TO 30

I J( KK) =1 BUFF ( NBR+2

)

I JILL ) = I BUFFC NBR+3

)

FA=BUFF ( NBR+4

)

FB = BUFF« NBR+ 5

)

SGA=BUFF ( NBR+6

)

SGB =BUFF (NBR+7)

NBR=N^R+7

C INCREMENT VARIOUS COUNTERS • JX IS NUMBER OF REFLECTIONS READ IN OFF

C THIS TAPE. JY IS TOTAL NUMBER OF REFLECTIONS READ IN.JZ IS TO SPACE

C REFLECTIONS IN PAIRS ON DRUM BECAUSE FO AND FC ARE TREATED AS

C SEPARATE REFLECTION.

JX= JX+

1

JY= JY+

JZ= JZ+2

C GET MAXIMUM AND MINIMUM INDICES FOR LATER PLOTTING OF DM VALUES

C AGAINST RANGES OF INDICES.

CALL MAXHKL I I J , J J , KK , LL

)

C CALCULATE DM VALUE AND RECORD REFLECTION LOCATIONS ON DRUM IN LOGFC

C WORD.

DM( JY) =( FA-SKALE+FB) /( SORT ( S GA** 2+SS*S GB** 2 )

)

LOGFC I JY )=JZ*1CCC00+JZ+1

130
58

59

60

61

62

63

64

65

66

67

68

69

70

71

72

73

74

75

76

77

78

79

8 C 81

82

83

84

85

86

87

O'n

SUM= SUM+ DM ( JY)**2

SUMM=SUMM+DM( JY

)

C WRITE REFLECTION ON DRUM.

CALL DRUMRT { NT , I J , FA . S QA

>

CALL DRUMRT ( NT , I J . F B, S GB

)

IF ( NBR. EQ.NBUFF ) Go TO 10

GO TO 20

C WRITE HOW MANY REFLECTIONS READ IN AND WHICH TAPE THEY WERE READ FROM

30 WRITE (NOUT.70) JY.NTAPE

C CALCULATE MEAN OF DM DISTRIBUTION

SUMM=SUMM/JY

C CALCULATE E.S.D. OF DM DISTRIBUTION.

SSUM=0 .0

DO 40 JX=1 . JY

SSUM=SSUM+ ( DM ( JX l-SUMM ) **2

40 CONTINUE

SSUM=SQRT( SSUM/( JY-1 )

)

WRITE (NOUT.50) SU M, SUM M , SSU

M

RETURN

C

Cc50 FORMAT (1H0.13H SUM DM**2 = .E10.4.11H MEAN DM = .F6.3.31H E.S.D.

20F DM DISTRIBUTION = ,F6.3)

0 FORMAT ( 8 1 H SCALE FACTOR APPLIED TO SECOND SET OF STRUCTURE FACTO

2RS IN SUBROUTINE REREAD IS.F10.3)

70 FORMAT (I10.27H REFLECTIONS READ FROM TAPE. 16)

C

END

131
12

3

4

5

6

7

8

9

1 0

1 1

1 2 13

14

15

1 6

17

18

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

4 2 43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

SUBROUTINE SCALE C ID A , J A , I DB , JB , JC , LOGFC * SK AL E , I J MI N A , I J MA X A, F A , FB

2 . SG A , SGB » JKL A , J KLB ,KA,KB)

C

C ******************************************************************

C

C SUBROUTINE TO CALCULATE THE LEAST SQUARES SCALE BETWEEN THE TWO SETS

C REFLECTIONS

C ALSO APPLIES CORRECTION FOR ABSORPTION IF NEEDED.

COMMON /TO/ IN. NOUT , NDR UM * NT APE A » NT APEB

C JJ.KK, AND LL TELL PACKING ORDER OF H K AND L IN JKL WORD.

C NR, IF NON-ZERO, IS NUMBER OF UNIT SORTED VALUES ARE TO BE WRITTEN ON

C (SAVES SORTING NEXT TIME). NABS IS NON-ZERO IF ABSORPTION

C CORRECTIONS ARE TO BE MADE.

COMMON /ORDER/ J JA » KKA , LLA , NR, NF ILER , NABS

C FA ETC ARE THE ARRAYS WHICH ARE ACTUALLY SORTED. THE SORTING IS DONE

C IN BATCHES WITH ONE CONSTANT INDEX TO SAVE TIME.

PARAMETER N = 50 , NNA=18C0, NN=7200, NNB=18C1, NNC=3601, NND=5401,

2 MM= 60 , MMM=4CC» NBUFF=490

DIMENSION FA(NNA> ,FB(NNA) .SGA(NNA) ,SGB(NNA) , JKLA(NNA), JKLB(NNA)

,

2 KA(NNA) , KB ( NNA

)

C IDA AND IDB ARE THE VALUES OF ONE OF THE MILLER INDICES. THEY ARE

C USED TO SELECT REFLECTIONS FROM THE DRUM. THE LOGFC ARRAY TELLS

C WHEPE THE SELECTED REFLECTIONS ARE ON THE DRUM.

D I ME NS ION LOGFC ( NN ) , I D A ( NN ) , ID B ( NN

)

C THIS IS A BUFFER FOR WRITING OUT ON UNIT NR IF REQUESTED.

DIMENS ION BUFF ( NBUFF ) , I BUFF ( NBUFF

)

EQUIVALENCE (BUFF, IBUFF)

SUM A = 0

SUM B=Q

NBR=0

JC = 0 A A= 1 . AB= 1 . IF (NR.GT.6) CALL POSXTN (NR.NFILER)

C BEGIN LOOKING FOR WHICH REFLECTIONS ARE READ AND SORTED

NJKL— I JMINA—

1

C LOOP OVER REFLECTION BATCHES BEGINS HERE.

10 NJKL=NJKL+1

C STOP, LOOK IF LAST HKL RANGE HAS BEEN PROCESSED ON THE COMPUTER

IF ( NJKL.GT. IJMAXA ) Go TO 90

KOUNTA=0

KOUNTB=0

C READ REFLECTIONS OFF DRUM IN BATCHES WITH SAME SLOWLY VARYING INDEX,

C DO FIRST REFLECTION SET HERE

DO 20 1=1, JA

IF ( IDA< I) .NE.NJKL) GO TO 2"

KOUNT A=K OUNT A+

1

C PREVENT POSSIBLE ARRAY OVERFLOW

IF( KOUNTA.GT.NNA) GO TO 25

CALL DRUMRD ( J, FT, SGT, NT, I

)

FA ( KOUNTA) =FT

SG A ( KOUNTA ) = SGT

JKL A ( K OUNT A ) = J KA ( KOUNTA)=I

20 CONTINUE

C DO SECOND REFLECTION SET HERE

25 DO 30 1=1, JB

132
I

58

59

60

61

62

63

64

65

66

67

68

69

70

71

72

73

74

75

76

77

78

79

80

8 1 82

83

84

85

86

87

88

89

90

9 1 92

93

94

95

96

97

98

99

ICO

1C1

102

1 03

10 4

105

106

107

1 C 8

109

no

i i i 1 12

1 13

1 14

1 1 5

IF < IDB( I ) .NE .NJKL ) GO TO 30

KQUNTB=K0UNT8+1

IF( KOUNTB. GT .NNA ) GO TO 35

I J A= I + J

CALL DRUMRD ( J , FT, SGT, NT , I JA

FB{ KOUNTB)=FT

SGB ( KOUNTB) =SGT

JKLB ( KOUNTB ) =

KB ( KOUNT B ) = I J A 30 CONTINUE

C CHECK THAT BOTH SETS INCLUDE SOME REFLECTIONS

35 IF( KOUNTA. EQ.O . OR. KOUNT B. EQ.O ) GO TO in

C WE NOW HAVE KOUNTA REFLECTIONS IN THE F A » SGA ,ETC. ARRAYS, AND

C REFLECTIONS IN THE CORRESPONDING FB, SGB, ARRAYS,

c take each jkla in turn and search jklb list for match.

C ASSUME EACH REFLECTION IS UNIQUE IN EACH SET

DO 80 1=1, KOUNTA

DO 70 M=l, KOUNTB

IF (JKLA(I)— JKLB(M)) 70,40,70

C CALCULATE QUANTITIES NEEDED FOR SCALE FACTOR ESTIMATION IF

C REFLECTIONS MATCH.

40 INDEX=JKLA( I

CALL HKL GEN( INDEX , J , K , L , J J A, KK A , LL A

IF (NABS.LE.O) GO TO 50

C CALCULATE ABSORPTION CORRECTIONS IF NABS IS NON-ZERO.

C APPLY ABSORPTION CORRECTIONS

CALL ABSORB ( J,K,L,AA,AB,DAA,DAB)

C PROPAGATE ERROR IN ABSORPTION CORRECTION INTO SIGMAS

VAR A= ( SGA< I ) /F A ( I)4 42+ ( DAA/AA) 442

VARB=(SGB( M)/FB( M) ) 4*2+ < DAB/ AB ) 4*2

F A ( I ) =F A ( I ) / A A FBI M)=FB(M)/AB

SG A ( I

) =F A( I)4SQRT( VARA)

SGB ( M) =F B( M ) 4 SORT ( VARB)

50 IF (NR. EQ.O) GO TO 60

C FILL UP BUFFER

I BUFF ( NBR+1 )=J

I BUFF ( NB R+2 ) = K

I B UFF ( NB R+3 ) =L

BUFF(NBR+4)=FA( I)

BU FF ( NBR +5 ) = FB(M)

BUFF ( NBR+6) =SGA{ I)

BUFF (NBR +7 )=SGB( M)

NBR=NBR+7

IF ( NBR.LT. NBUFF) GO TO 60

C WRITE OUT BUFFER CONTAINING 70 REFLECTIONS

WRITE (NR) BUFF

C REZERO BUFFER COUNTER.

NBR=0

C STORE RE-ORDERED INFORMATION FOR DM CALCULATIONS LATER.

60 JC=JC+1

IF(JC.GT.NN) GO TO 90

C LOGFC CONTAINS THE POSITION OF THE REFLECTIONS ON DRUM

LOGFC( JC )=K A ( I )4100C00+KB( M)

SUMA=SUMA+FA( I

SUMB=SUMB+FB( M)

70 CONTINUE

80 CONTINUE

KOUNTB

133
1 1 6

1 17

HR

119

1 20

121

122

123

1 24

1 25

1 26

127

128

1 29

13C

131

132

133

134

135

136

137

138

139

140

141

142

1 43

144

145

146

C NOW THAT ALL REFLECTIONS FOR GIVEN SLOWLY VARYING INDEX HAVE BEEN READ

C OFF THE DRUM AND SORTED* I NCREASE THE SLOWLY VARYING INDEX BY 1 AND

C SEE IF THERE ARE ANY MORE REFLECTIONS To READ AND SORT. CONTINUE UNTIL

CREFLECTIONS FOR HIGHEST VALUE OF SLOWLY VARYING INDEX HAVE BEEN RE AO

C AND SORTED.

GO TO 10

C CHECK IF WE DO HAVE SOME COMMON REFLECTIONS.

90 IF (JC.GT.O) GO TO 100

WRITE (NOUT,120)

STOP

C CALCULATE THE SCALE FACTOR BETWEEN TWO SETS OF REFLECTIONS

IOC CONTINUE

SKALE=SUMA/SUMB

WRITE (NOUT.130) NTA PE A » NT APE B , SK ALE

,

J C

,

J A

,

JB

IF ( NR.EQ.O) GO TO 1 10

C SET END OF DATA FLAG ( H=99 ) AND WRITE OUT LAST PART OF DATA ON UNIT

C NR .

I BUFF ( NBR+1 ) =99

WP f TE (NR) BUFF

END FILE NR

110 RE TURN

C

C

c120 FORMAT (49H STOPPED IN SCALE ROUTINE,- NO COMMON REFLECTIONS)

130 FORMAT (46H THE SCALE FACTOR BETWEEN DATASETS 1 FROM TAPE, 16, 1 7H A 2ND 2 FROM TAPE, 16, 3H IS,F8.3/11H THERE WErE,I6,38H REFLECTIONS CO

3MM0N TO BOTH SETS , W IT H , I 6/ 2 5H REFLECTIONS IN SET 1 AND,I6,2CH REF

4LECTIONSIN SET 2)

C

END

134
123

4

5

6

7

8

9

10

1 1 12

1 3

14

15

16

17

18

19

20

21

22

23

24

25

26

27

2e

29

30

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

SUBROUTINE SIMLEQ (N,A,C,X)

Cc ******** ********** ************************************** ** ********

c

C GENERAL SOLUTION OF A SET OF LINEAR SIMULTANEOUS EQUATIONS

C SIGMA( A( I * J) *X< J) ) =C( I > BY WILLIAM V. LOEBENSTEIN

DIMENSION A( 6 » 6 ) «B

(

616 ) ,F( 6 , 6 ),X< 6 ) ,XR( 6 ) . BR( 6 . 6 ) ,D{ 6 ) . DR( 6 ) ,C( 6 ) D( 1 )=C( 1

)

DO 10 J= 1 » N

10 B( 1. J)=A( 1 , J)

N1 = N —

DO 20 I=1,N1

20 F( 1, I)=-A( 1 + 1, 1 )/B( 1, 1 ) DO 50 L= 2 » N

D(L)=C(L)

L1=L-1

DO 30 K=1,L1

30 D(L)=0(L)+F(K,L1)*D(K)

DO 40 J=L,N

B( L » J) =A ( L * J ) DO 40 K= 1 » L 1 40 B(L, J)=B(L, J)+F(K,Ll )*B<K, J)

IF (L.GT.N1) GO TO 60

DO 50 I=L.N1

FCL, I ) =- A( 1+1 »L)/B(L,L)

DO 50 K=1,L1

50 F(L,I)=F(L.I)-F(K, I)*B(K,L)/B(L.L)

60 DO 70 1=1,

K=N+ 1 — I 70 DR(K)=D(I)

DO 80 J= 1 ,

DO 80 1=1,

K=N+ 1-

L=N+ 1 — J 80 BR( K,L) =B( I , J

)

XR ( 1 ) = DR( 1 )/BR< 1,1)

DO 90 K= 2 ,

XR(K)=DR(K)/BR(K,K)

K 1 = K—

DO 90 1=1, K1

90 XR(K) = XR(K)-BR(K,I)*XR( I)/BR(K»K)

DO 100 K= 1 » N

I=N+ 1 —

100 X(I)=XR(K)

RETURN

C

END

135
123

4

5

6

7

8

9

1 0

1 1

1 2

l 3

1 4

15

16

1 7

16

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

46

47

48

49

50

51

52

53

54

55

56

57

C

cc

c

ccc

c

c

c

cc

cc

c

c

c

cc

cc

cc

cc

ccc

cc

c

ccc

c

SUBROUTINE SORT ( Y , Y D 1 ,YD2»YD3*N» NL

)

******** ************** ************ ********* ************* **********

THIS ROUTINE SORTS THE ELEMENTS OF THE INPUT VECTOR Y.

Y IS A SINGLE PRECISION VECTOR OF (UNSORTED) OBSERVATIONS,

THE INTEGER VALUE N = NUMBER OF ELEMNETS . NL IS THE

MAXIMUM NUMBER OF ELEMENTS DIMENSIONED FOR IN THE CALLING PROGRAM,

THE OUTPUT FROM THIS ROUTINE IS THE SINGLE PRECISION VECTOR Y INTO

THE SORTED OBSERVATIONS HAVE BEEN PLACED.

RESTRICTIONS ON THE MAXIMUM ALLOWABLE VALUE OF N—THE DIMENSIONS

OF VECTORS IU AND IL (DEFINED AND USED INTERNALLY WITHIN THIS ROUT

DETERMINE THE MAXIMUM ALLOWABLE VALUE OF N FOR THIS

ROUTINE. IF IU AND IL EACH HAVE DIMENSION K, THEN N MAY NOT EXCEE

2**(K+1) - 1. FOR THIS ROUTINE AS WRITTEN, THE DIMENSIONS OF IU A HAVE BEEN SET TO 36, THUS THE MAXIMUM ALLOWABLE VALUE OF N IS

APPROXIMATELY 137 BILLION, SINCE THIS EXCEEDS THE MAXIMUM ALLOWAB

VALUE FOR AN INTEGER VARIABLE IN MANY COMPUTERS, AND SINCE A SORT

BILLION ELEMENTS IS PRESENTLY IMPRACTICAL AND UNLIKELY, THEREFORE

TEST FOR WHETHER THE INPUT SAMPLE SIZE N EXCEEDS 137 BILLION HAS B INCORPORATED INTO THIS ROUTINE. IT IS THUS ASSUMED THAT THERE IS

(PRACTICAL) RESTRICTION ON THE MAXIMUM VALUE OF N FOR THIS ROUTINE

PRINTING NONE UNLESS AN ERROR CONDITION EXISTS

THIS ROUTINE IS SINGLE PRECISION IN INTERNAL OPERATION.

SUBROUTINES NFEDED —NONE

SORTING METHOD BINARY SORT

REFERENCE CACM MARCH 1969, PAGE 186 (BINARY SORT ALGORITHM BY RIC

C. SINGLETON.

PAGE 54.

PAGE 624.

PAGE 41.

STATISTICAL ENGINEERING LABORATORY ( WASHINGTON, D.C. 20234 JUNE 19

WRITTEN

NATI ONAL

-'-C ACM JANUARY 1970,

CACM OCTOBER 1970,

J AC M JANUARY 1961,

BY JAMES J. FILL IBEN

,

BUREAU OF STANDARDS,

TIME FOR SORTING IS PROPORTIONAL TO N LN N IN CONTRAST TO N**2

FOR THE BUBBLE SORT.

DIMENSION Y ( 1 ) , YD1 ( 1 ) , YD2( !),Yd3C1)

DIMENSION IUI 36) . I L( 36)

C

C

C

10

c20

c30

I PR = 6

CHECK THE INPUT ARGUMENTS FOR ERRORS

IF (N.LT.l) GO TO 20

IF (N .EQ .1 ) GO TO 30

HOLD=Y t 1 ) DO 10 1=2,

N

IF ( Y ( I ) . NE • HOLD ) GO TO 40

CONTINUE

WRITE ( I PR ,180) HOLD

RETURN

WRITE ( I PR , 1 90

)

WRITE ( I PR ,210) N

RETURN

WRITE ( I PR , 200

)

136
58

59

60

61

62

63

64

65

66

6T

68

69

7 C 71

72

73

74

75

76

77

78

79

80

81

82

83

84

85

86

87

88

89

90

91

92

93

94

95

96

9 7

98

99

ICC

101

102

1 C 3

1 C 4

105

106

1 C 7

108

109

110

1 1 1 1 1 2

1 1 3

1 1 4

1 1 5

c*n

cn

nnnnnnr>4»n

RETURN

0 CONTINUE

I SLIGHT MODIFICATION BY L • W. SCHROEDER, ( NB S » WASH I NGTON , D . C. ) SO THAT

DEPENDENT ARRAYS ARE REARRANGED AS IS ARRAY Y TO PRESERVE RELATIVE

INDEXING.

CHECK TO SEE IF THE INPUT VECTOR IS ALREADY SORTED

NM 1 = N— 1 DO 50 1=1, NM1

IP1=I+1

IF ( Y( I ) .LE. Y< IP 1 ) ) GO TO 50

GO TO 60

0 CONTINUE

RETURN

C M=1

1= 1 J=N

70 IF CI.GE.J) GO TO 140

8 0 K=I

MID=(I+J)/2

AM ED= Y ( M ID ) C SAVE CORRESPONDING VALUES FOR DEPENDENT ARRAYS

C0RYD1=YD1 (MID)

CO RY D2=YD2(MID)

COR YD3=Y D3 (MID)

IF ( Y( I ) .LE.AMED) GO TO 90

Y( MI D) =Y ( I ) Y ( I ) =A MED

AMED=Y (MID)

C EXTEND EXCHANGE TO DEPENDENT ARRAYS

YD 1 ( MID) =YD 1(1)

Y02 ( MID )=YD2 ( I ) YD 3 ( MID) =YD3< I

)

YD 1(1 ) = COR YD 1 YD2 ( I ) =CORYD2

YD3 ( I ) =CORYD3

C0RYD1=YD1 (MID)

CORY D2=Y D2 (MID)

CORYD3=YD3 ( MID)

90 L=J

IF (Y( J) .GE.AMED) GO TO 110

Y ( MID) = Y ( J

)

Y ( J )=AMED

AMED=Y( M ID ) C EXTEND EXCHANGE TO DEPENDENT ARRAYS

YD 1 ( MID )=YD1 ( J ) YD2( MID) =YD2 ( J

)

YD3(MI D)=YD3( J

)

YD1 ( J ) = C ORYD

1

YD 2 ( J ) =C ORYD 2 YD3 ( J )=C0RYD3

CORYDl =YD1 ( MID)

C0RYD2=YD2( MID)

CORY D3=Y D3 ( MID)

137
1 1 6 IF ( Y( I ) .LE. AMED) GO TO 110

1 1 7 Y( MID) =Y ( I )

1 18 Y ( I ) =A ME D

119 AMED=Y (MID)

1 20 C EXTEND EXCHANGE TO DEPENDENT ARRAYS

121 YOl (MID ) = YD 1 ( I ) 122 YD2 ( MI D) =YD2 ( I ) 127 YD3 ( MID) =YD3( I

)

124 YD1 ( I

) =CORYD

1

125 YD 2(1) =CORYD2

1 26 YD3 ( I )=CORYD3

1 27 CORYD1 =YD1 ( M ID)

1 28 CORY D2=YD2 (MID)

1 29 CORYD3=YD3( MID)

130 GO TO 110

131 c 132 1 00 Y( L )=Y( K ) 133 Y ( K ) =T T 134 c SET VALUES FOR CORRESPONDING ELEMENTS IN DEPENDENT ARRAYS

1 35 YD 1 ( L ) = YD1 ( K

)

136 YD2 (L)=YD2f K)

137 YD3( L)=YD3(K)

138 YDl ( K ) =T YD 1 139 YD 2 ( K ) =T YD2

140 YD3 ( K ) =T YD 3

141 1 10 L=L-1

142 IF ( Y( L) .GT. AMED) GO TO 110

143 TT = Y(L )

1 44 c SAVE VALUES FOR DEPENDENT ARRAYS

145 TYD 1=YD 1 ( L

)

1 46 TYD2=YD2 (L

)

147 TYD 3=YD3 ( L

)

148 1 20 K=K + 1

1 49 IF (Y(K) .LT. AMED) GO TO 120

1 50 IF ( K oLE.L) GO TO 104

151 L M I =L—

I

1 52 JMK=J-K

153 IF (LMI.LE.JMK) GO TO 130

1 54 XL( M)~I

155 IU( M)=L

156 I=K

157 M=M+ 1 158 GO TO 150

159 c

1 60 1 30 IL ( M ) = K

161 IU( M ) = J 162 J=L

163 M=M + 1 164 GO TO 150

165 C

l 66 1 40 M= M— 1 167 IF ( M.EQ.O ) RETURN

168 1= I L ( M )

1 69 J=IU(M)

170 150 JMI=J-I

171 IF (JMI.GE.ll) GO TO 80

172 IF ( I • EQ • 1 ) GO TO 70

173 1=1-1

138
174

175

176

177

178

179

1 80

181

182

183

184

185

1 86

187

188

1 89

190

191

192

193

194

195

196

197

198

199

200

201

20 2

203

204

205

206

207

U

U

160 1 = 1+1

IF ( I .EQ .J ) GO TO 140

AMED=Y C 1 + 1 ) C SAVE CORRESPONDING VALUES FOR DEPENDENT ARRAYS

CO RY D 1 =Y D 1 ( 1+1 ) COR YD2 =YD2 ( I + 1

)

C O R Y D 3= Y D 3 ( 1+ 1 ) IF (Y( I ) .LE. AMED) GO TO 160

K= I 170 Y(K+1)=Y(K)

C REARRANGE DEPENDENT ARRAYS IN LIKE MANNER

YD 1( K+ 1 ) =YD 1 ( K)

YD2IK+1 )=YD2(K)

YD 3 ( K+l )=YD3(K)

K = K— 1 IF ( AMED .LT .Y(K ) ) GO TO 170

Y( K+l ) =A MED

C SET VALUES FOR CORRESPONDING ELEMENTS IN DEPENDENT ARRAYS

YD 1 ( K+ 1 ) =CORYD

YD2 ( K+ 1 ) =COR YD 2 YD3( K+ 1

) =C0RYD3

GO TO 160

1 80

1 9C

200

210

C

FORMAT ( 1H » 1 0 8 H*4 *** NON-FATAL DIAGNOSTIC THE FIRST INPUT ARGUM

2ENT (A VECTOR) TO THE SORT SUBROUTINE HAS ALL ELEMENTS =»E15.8,6

3 H 4***4)

FORMAT (1H « 58H FATAL ERROR - NUMBER OF ELEMENTS To BE SORTED IS N

2EGATIVE)

FORMAT (1H « 64 H DI AGONISTIC - ONLY ONE ELEMENT TO BE SORTED BY SOR

2 T SUBROUTINE)

FORMAT ( 1H »35H***4* THE VALUE OF THE ARGUMENT IS »I8,6H *#4*4)

END

139
123

4

5

6

7

8

9

1 0

1 1

1 2

13

1 4

15

1 6

17

18

19

20

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

37

38

39

40

41

42

43

44

45

4 6

47

48

49

50

5 1 52

53

54

55

56

57

SUBROUTINE SPEC C JC , SK A l_E , DM »LOGFC » NN ) C

C 44 4*4* 44 4444*444444444444444 4*44 444444444444444444*44444444**4444*

C

DI MENS ION DMC NN) ,LOGFC( NN ) » SKU 1 ), DY( 1 1 , 1 1 ) , SDMSC 11). SDMI 11)

DIMENSION NSIGN(ll)

C FIND THE SCALE FACTOR WHICH MAKES THE SUM OVER DM442 A MINIMUM.

N=5

DX=0 .006

NTRY=1

10 IF (NTRY.EQ.2) N= 1

1

NF=N/2+l

NC H=0

C

C CALCULATE SCALE FOR POINTS ON EITHER SIDE OF TRIAL VALUE TO MAKE A C DIFFERENCE TABLE.

C INITIALIZE SUMS

DO 20 1 = J P 1 1 SDMS C J ) =0 .0

20 SDMI j) =C . A DO AO I = 1 » N

FI=I-NF

SKI I ) =SKALE*( 1 .+DX4F I ) PRINT 180 * SKI I ) CALL DMC ALC ( JC , SK ( I ) » D M . LOGFC , NN

5

DO 30 J=1,JC

SDMI I ) =SDM C I )+DM( J ) 30 SDM S( I ) — SDMS I I )+DMC J)442

40 CONTINUE

C FORM FIRST DIFFERENCES <K=1) FOR DIFFERENCE TABLE.

NN=N-1

DO 50 J= 1 » NN

DYC 1 » J ) = SDMS! J+ 1 ) — SDMS I J)

C CHECK TO SEE THAT SDMS GOES THRU A MINIMUM.

NSI GNI J ) =S I G N < 1 * DY ( 1 « J ) ) IF « NSIGNC J) .NE.NSIGNl 1 ) ) NCH=1

50 CONTINUE

C CHECK IF SIGN CHANGED. IF NOT START OVER

IF INCH.EQ.l) GO TO 60

IF (NTRY.EQ.2) PRINT 190

IF (NTRY.EQ.2) RETURN

NTRY=NTR Y+

1

GO TO 10

C

60 CONTINUE

C PRINT OUT TABLE HEADING

PRINT 20

0

C FORM K TH ORDER ADJUSTED DIFFERENCES

DO 70 K=2»NN

L=N-K

DO 70 J=1,L

7 0 DY(K,J) = (DYCK-l.J+l)-DY(K-l , J) ) /NFACTIK

)

C WRITE OUT DIFFERENCES.

DO 80 K= 1 » N

80 PRINT 210, SKI I ) .SDMI I ) .SDMS Cl ) . CDY ( J. I) , J= 1 » N)

C FIND MINIMUM VIA 0.001 STEPS IN SC AL E , ST AR T I NG FROM MINIMUM CALCULATE

SMI N=1 04410

140
58

59

60

61

62

63

64

65

66

67

6 R 69

70

71

72

73

74

75

76

77

78

79

80

81

82

83

84

85

86

87

88

89

9C

91

92

93

94

95

96

97

98

99

ICO

101

102

103

1 C 4

1 C 5

106

107

108

U

U

90

1 00

c

c

1 10

1 20

1 30

140

C

1 50

C

1 60

C

1 70

180

1 9C

200

210

2 20

C

DO 100 J = 1,N

IF < SDMS( J)-SMIN) 90,100,100

SMIN=SDMS( J

)

JMIN=J , CONTINUE

INTERPOLATE USING ADJUSTED DIFFERENCES via NEWTONS FOR MULA ( 3RD. ORDER

)

SC = SK( JMIN)

NTRY=0

TRY MAKING SCALE LARGER

NEX=2

SC=SK( JMIN

)

SSAVE=SK( JMIN)

CONTINUE

SC=SC+.00 14 ( - 1 . 0**NEX)

NT RY=NTRY+

MO=( SC-SK( 1 ) +D X ) /D X IF (MO) 170, 170, 120

IF (MO.GT.N) GO TO 170

IF (N-3-MO) 130,140,140

MO= N—

3

F MO= MO—

U=( SC— SK(JMIN) —FMO 4D X) /DX

SI N=( ( ( U — 2 • ) *DY (3, MO) + DY< 2, MO ) )*<U- 1 . ) +DY( 1 , MO) )*U+SDMS( MO)

IF ( NTRY.EQ. 1 . AND. SI N. GT.SMI N) GO TO 150

IF ( SIN.GT.SMIN) GO TO 160

SSAVE=SC

SSIN=SIN

GO TO 110

SC=SK( JMIN

)

NEX= 1 GO TO 1 1

0

PRINT 220 , SKALE ,SSAVE,SSIN

SKALE=SSAVE

CALL DMCA^C ( JC , SK AL E , DM , LOGFC , NN

)

RETURN

PRINT 190

RETURN

FORMAT ( 1H0, 15H TRIAL SCALE = ,F10.4,25H SET BY SUBROUTINE =SPEC=)

FORMAT ( 1 HO , 43 H SCALE MINIMIZATION FAILED, SCALE NOT CHANGED)

FORMAT (1H0.52H SCALE SUM DM SUM DM**2 DIFFERENCE

2 S ) FORMAT (2X.14E9.4)

FORMAT ( 1 HO , 1 5HPRE V I OUS SC AL E= , F 1 0 . 4, 1 0 HNE W SC AL E= , F 1 0 . 4 , 9HSUM D M S

2Q=,E9.4)

END

141
123

4

5

6

7

a9

1 0

1

1

1 2

1 3

14

1 5

16

17

18

1 9

29

21

22

23

24

25

26

27

28

29

30

31

32

33

34

35

36

C

C

c

c

c

cc

ccc

cc

c

c

c

ccc

cc

cc

cc

c

cc

SUBROUTINE UNI MED <N,X1 ******************************************* *********************** STATISTIC (FOR I = 1,2 > N ) FROM A UNIFORM DISTRIBUTION (ON THE

INTERVAL (C,l)»«

THIS IS IDENTICAL TO THE MEDIAN OF THE BETA DISTRIBUTION WITH PARA

I AMD N-I + l FOR 1=1,2, ...,N.

THE INPUT TO THIS ROUTINE IS THE DES ipED INTEGER SAMPLE SIZE N

AND AN EMPTY SINGLE PRECISION VECTOR X (OF DIMENSION AT LEAST N) I WHICH THE N GENERATED UNIFORM ORDER STATISTIC MEDIANS WILL BE PLAC

THE OUTPUT FROM THIS ROUTINE IS THE SINGLE PRECISION VECTOR X INTO WHICH THE N GENERATED UNIFORM ORDER STATISTIC MEDIANS

HAVE BEEN PLACED.

ALL OF THE PROBABILITY PLOT ROUTINES MAKE USE OF THIS ROUTINE.

JUSTIFICATION AND ACCURACY OF THE ALGORITHM USED IS FOUND IN AN

UNPUBLISHED JJF MANUSCRIPT.

THERE IS NO RESTRICTION ON THE MAXIMUM VALUE OF N FOR THIS ROUTINE

PRINTING NONE UNLESS AN ERROR CONDITION EXISTS

THIS ROUTINE IS SINGLE PRECISION IN INTERNAL OPERATION

SUBROUTINES NEEDED NONE

REFERENCE UNPUBLISHED JJF MANUSCRIPT

WRITTEN BY JAMES J. FILLIBEN, STATISTICAL ENGINEERING LABORATORY ( NATIONAL BUREAU OF STANDARDS, WASHINGTON, D.C. 20234 JUNE 19

DIMENSION X( 1}

A N= N

IPR=6

CHECK THE INPUT ARGUMENTS FOR ERRORS

IF < N.LT. 1 > GO TO 1C

IF (N.EQ.l) GO TO 29

GO TO 30

37 10 WRITE (IPR

38 WRITE (IPR

39 RETURN

40 c

41 20 WRITE (IPR

42 30 CONTI NUE

43 C

44 G AM= •

3

45 DO 40 1=1,

46 A I = I 47 X( I ) = ( AI-G

48 40 continue

49 RETURN

50 C

51 c 52 50 FORMAT ( 1H

53 2 UNIMED SUI 54 60 FORMAT ( 1H

55 2 ENT TO the

56 70 FORMAT ( 1H

57 C

58 END

,9lH***** FATAL ERROR THE FIRST

5ROUTINE IS NON-POSITIVE  )

,1C0H** *** NON-FATAL DI AGNOSTIC — — THE FIRST

UNIMED SUBROUTINE HAS THE VALUE 1 *****)

,35H***** THE VALUE OF THE ARGUMENT IS ,18

142

INPUT ARGUMENT TO THE

INPUT ARGUM

6H  )
nn SUBROUTINE USER ( N A » O, S O . X , S X i B » SB , L , NT APE ) ** ****** ************ **************************************** ****** PARAMETER N = 5<*. NN A= 1800* NN=72C0» NNB=1801, NNC=36'M. NN0=54C1»

2 MM=60» MMM=4CC'» NBUFF=70C, MMMM=32

COMMON /H/ DUMMY (NN)

DT MENS ION NA (MM) , 0( MM) , SO( MM) , X( 3,MM) , SX ( 3, MM) , B(6, MM ) , SB( 6,MM

C PUT YOUR READ IN STATEMENTS HERE FOR SPECIAL USER INPUT OF PARAMETERS

C AND STANDARD DEVIATIONS.

C NA = NAME OF ATOM

C 0= OCCUPANCY

C SO = ERROR ON OCCUPANCY

C X= POSITIONAL PARAMETERS

C SX = ERROR ON POSITIONAL PARAMETERS

C 8= THERMAL PARAMETERS

C SB= ERROR ON THERMAL PARAMETERS

C L= NUMBER OF ATOMS IN LIST

C NT APE= UN IT WITH PARAMETERS ON IT.

C USE EQUIVALENCE STATEMENT WITH DUMMY (1) AS IN SUBROUTINE I NF I NG TO

C SAVE STORAGE SPACE.

PETUPN

C

END

143
V. Description of NBS Magnetic Tape 11

This tape contains one file of card images, one 80

character card per record, while comprise the FORTRAN program

STATUS described in this report.

The file contains the STATUS main routine and 35

subroutines in the following order: ABSORB, ANISO, CATLOG,

CELL, CHANGE, DMCALC , DMDQ , DPCALC, DRUMRD , DRUMRT , FINGFO,

FITCHK , FREAD, HKLGEN , INFING, INFOFC, LSFIT , MAXHKL,

MITCHEK, NFACT , PINV, PLOTEM, PLOTS, POSITN, PRIME, PROB,

PRYNT, PRYNTT, REREAD, SCALE, SIMLEQ , SORT, SPEC, UNIMED,

and USER.

The file of card images was recorded on a certified 2400

foot reel in BCD characters, 556 bpi density, even parity and

in 7-track mode. The National Technical Information Service

can generate tapes in other modes if so requested.

144
SBS-1 14A (REV. 7-73) P UTs7 DEPT. OF COMM.

BIBLIOGRAPHIC DATA

SHEET

1. PUBLICATION OR REPORT NO.

76- 1102

2. Gov’t Accession

No.

3. Recipient’s Accession No.

4. TITLE AND SUBTITLE

Computer Programs for Structural Chemistry! STATUS.

A FORTRAN Program for Statistical Analysis of

Crystallographic Quantities

5. Publication Date July, 1976

6. Performing Organization Code

7. AUTHOR(S)

LeRoy W. Schroeder and Brian Dickens

8. Performing Organ. Report No.

7 6-

9. PERFORMING ORGANIZATION NAME AND ADDRESS

NATIONAL BUREAU OF STANDARDS

DEPARTMENT OF COMMERCE

WASHINGTON, D.C. 20234

10. Project/Task/Work Unit No. 3110585

11. Contract/Grant No. DE 00572

12. Sponsoring Organization Name and Complete Address (Street, City, State, ZIP)

American Dental Association

211 East Chicago Avenue

Chicago, Illinois 60611

13. Type of Report & Period Covered

14. Sponsoring Agency Code

15. SUPPLEMENTARY NOTES

The program is stored on NBS Magnetic Tape 11, available from NTIS.

16. ABSTRACT (A 200-word or less factual summary of most significant information. If document includes a significant

bibliography or literature survey, mention it here.)

This Report describes a FORTRAN computer program for evaluation of (i) the

results of crystallographic least-squares refinements by examination of the residuals,

(ii) differences in sets of data collected by different methods from the same crystal,

(iii) differences in data sets collected by the same method from different crystals

of the same material, and (iv) the differences in parameters in different models

representing the crystal structure of the same material.

Part of the evaluation is accomplished by plots of residuals against the expected

normal distribution quantiles. Additional plots compare residuals with the independent

variable, d*, and with the calculated variable Fc2/sin28. The Miller indices can be

used to divide a data set into various classes and octants so that the possibility

of anisotropic effects can be examined.

The program will treat up to 7200 data points in each of two experimental

data sets or 1000 structural parameters, and is oriented specifically towards

examining crystallographic data.

17. KEY WORDS (six to twelve entries; alphabetical order; capitalize only the first letter of the first key word unless a proper name; separated by semicolons

)

Crystallographic data; errors; probability plots; residuals; statistics; uncertainties

18. AVAILABILITY [£ Unlimited

^ For Official Distribution. Do Not Release to NTIS

_ Order From Sup. of Doc., U.S. Government Printing Office Washington, D.C. 20402, SD Cat. No. 03

19. SECURITY CLASS

(THIS REPORT)

UNCL ASSIFIED

X Order From National Technical Information Service (NTIS)

Springfield, Virginia 22151 Manuscript & tape available.

20. SECURITY CLASS

(THIS PAGE)

UNCLASSIFIED

21. NO. OF PAGES

151

22. Price

$ 6.75

USCOMM-DC 29042-P74
.