U.S. DEPARTMENT OF THE INTERIOR 

U.S. GEOLOGICAL SURVEY

An interactive computer program to extrapolate the clay 

fraction distributions of truncated grain-size data

by

Poppe, L.J. 1 and Eliason, A.H/

Open-File Report 99-27

This report is preliminary and has not been reviewed for conformity with U.S. Geological Survey 

editorial standards (or with the North American Stratigraphic Code). Any use of trade, product, or 

firm names is for descriptive purposes only and does not imply endorsement by the U.S. 

Government.

1 Coastal and Marine Geology Program, USGS, Quissett Laboratories, Woods Hole, MA 02543

2 Eliason Data Services, 230 Meetinghouse Road, Mashpee, MA 02649
ABSTRACT

Since the mid-1960's, grain size analyses have been performed on marine sediment samples 

at the facilities of the Coastal and Marine Geology Program, U.S. Geological Survey in Woods Hole, 

Massachusetts. However, owing to analytical limitations, the procedures utilized commonly truncate 

the tails of the fine-fraction portion of the grain-size distributions. This report presents an interactive 

computer program that can be used to extrapolate incomplete grain-size distributions down to the 

colloidal clay fraction.

INTRODUCTION

Grain size is the most fundamental property of sediments. Geologists and sedimentologists 

use it to study trends in surface processes and sediment transport; engineers use it to study stability 

under load; geochemists use it to study kinetic reactions and the affinities of fine-grained particles 

and contaminants; and hydrologists use it when studying the movement of subsurface fluids (Blatt 

and others, 1972; McCave and Syvitski, 1991). Inasmuch as many geological observations consist 

of measurements made on a large number of specimens, the techniques and equipment used for 

particle-size analysis must be fast, accurate, and yield highly reproducible results.

The fine-grained fraction of a sediment is defined as silt (particles with diameters less than 

62 microns down to 4 microns) plus clay (particles with diameters less than 4 microns, with colloidal 

clay being less than 0.1 microns). Classical techniques to analyze the fine-grained portion of grainsize distributions are dominated by sedimentation methods, such as the hydrometer and pipette 

methods (Krumbein and Pettijohn, 1938; Folk, 1974). These methods involve preparing a dispersed, 

homogeneous suspension of the fine fraction, dilution with a dispersant solution (usually sodium 

hexametaphosphate) to 1000 ml, and allowing the particles to settle in a graduated cylinder. As the 

settling of sediment particles continues according to Stoke's Law, samples are either withdrawn 

(pipette) or measurements made (hydrometer) of the suspension at preset time intervals. However, 

many sedimentation laboratories have discontinued or limited the use of pipette and hydrometer 

techniques because of inherent problems with settling (e.g. Brownian motion), thermal convection, 

irregular particle shape, mass settling in density currents, rounding errors due to large multiplication 

factors, and, especially, the time necessary to extend an analysis down into the clay range. For 

example, over 24 hours are required to extend a pipette analysis down to 11(|) at 20°C, room 

temperature (Krumbein and Pettijohn, 1938).

Although originally designed to count and size blood cells (Coulter, 1957; Berg, 1958), 

electro-resistance multichannel particle-size analyzers have found other applications because of the 

short time and small amount of material required for analysis. These devices, which are currently 

manufactured by Coulter/Beckman (Coulter Counter) and Particle Data (Elzone), are also used in 

industry, in the biological sciences (Sheldon and Parsons, 1967), and in geology, where they are used 

to measure the grain size of sediments by measuring differences in electrical resistivity produced by 

sediment particles (McCave and Jarvis, 1973; Schiedler, 1976; Kranck and Milligan, 1979; 

Muerdter and others, 1981; Milligan and Kranck, 1991).

Because the aperture tubes used by electro-resistance multichannel particle-size analyzers 

can size particles of only 2-40% of the aperture diameter, at least two aperture tubes with 

overlapping ranges are required to determine the size distribution within the fine fraction of a marine

1
sediment. For example, during a standard multi-aperture analysis using a Coulter Counter, a 200 

micron aperture tube would be used to resolve the size distribution between 64 to 8 microns, and a 

30 micron aperture tube would be needed to resolve the size distribution between 12 microns down 

to approximately 0.6 microns. These individual analyses or passes through the respective aperture 

tubes are then mathematically combined to produce the fine-fraction (silt plus clay) distribution. The 

fine fraction is, in turn, combined with the coarse fraction data to generate a complete grain size 

distribution.

An obvious limitation of this method is that electro-resistance multichannel particle-size 

analyzers, even if calibrated to resolve down to about 0.6 microns (10.75(|>), can only resolve a 

portion of the clay size distribution. The fine clay in the 0.6 to 0.1 micron range (some of the 11(|> 

and all of the 12(|> and 13(|> fractions that extend down to the colloidal clay boundary) is not detected. 

Although an analysis performed down to 0.6 microns is adequate for most freshwater, estuarine, and 

shelf environments, the sediments from many deeper water marine environments (e.g. rise and 

abyssal plain) may contain significant material in the fine clay fraction. By combining a doublepoint pipette analysis with the data from a electro-resistance multichannel particle-size analysis 

details of the fine-fraction can be extended through the fine clay range, but this solution is less 

desirable because of the problems with pipette analyses cited above and because the resultant data 

are no longer at whole-(j) intervals.

Earlier efforts to extrapolate grain size data beyond the limits of the analysis used graphical 

methods (Schlee and Webster, 1967; Schlee, 1973). For example, if Schlee determined that more 

than 5% of the fine-grained sediment from a sample was less than 1 micron, additional data points 

were estimated for the finer sizes at whole-(|> intervals by projecting the grain-size curve as plotted 

on probability paper and following the slope of the line. While this is a practical solution, it is much 

more labor intensive than the computer program presented herein.

SOFTWARE

Data Input

The required format of the input file and the data dictionary for the fields of the input file are 

supplied in Appendices A and B, respectively. The fields within the input file are comma delimited; 

records are delimited by dollar signs. No imbedded spaces or commas are allowed in any of the 

fields of the input file.

This particular format was chosen for the input file because it is the same format as that of 

the output file (a pre-database file) created by GSTAT, the main data processing and statistics 

program used in the Sedimentation Laboratory of the Woods Hole Field Center, Coastal and Marine 

Geology Program, U.S. Geological Survey (Poppe and others, 1985). Although many of the fields 

that describe the sample history are unnecessary and may be left blank, we recommend their 

inclusion because they provide a more useful record. The necessary fields include the lab number 

and the phi class and associated cumulative frequency percentages. The lab number is important for 

two reasons. First, it uniquely identifies the sample; any number or character string can be used. 

Second, it is the first field; without it the program will crash because the record delimiter will appear 

out of order. The absence of any identifiers (e.g. navigation, sampling device), attributes (e.g. verbal
equivalent), or data (e.g. percent clay, statistics) in the remaining fields will not prevent the program 

from functioning..

Software Usage

The program CLAYEST, written in Microsoft QuickBasic, can be used to estimate the 

amount of clay present in a range between a user-defined upper size limit and 0.1 microns. It was 

designed to calculate that portion of the clay fraction not detected during a particle-size analysis by 

a Coulter Counter down to 13(j> (Coulter Electronics, 1988). Documentation is supplied in Appendix 

C; the software code is provided in Appendix D. The program, when specified, will create a 

hardcopy (Fig. 1) and/or a digital data file (Fig. 2). The optional hardcopy will contain: the sample 

identifiers; the original (data supplied with the input file) and revised (data generated by the 

program) frequency percentages; the original and revised percentages of gravel, sand, silt, and clay; 

a revised sediment classification (verbal equivalent); and the revised Method of Moments statistics. 

All sample identifiers of with a length greater than 11 characters will be truncated. The optional data 

file will be comma delimited and contain: the sample identifiers; the revised percentages of sand, 

gravel, silt, and clay; the revised verbal equivalent; the revised statistics; and the revised frequency 

percentages (-5 to 13(j>). The first line of this data file will contain a list of the field attributes.

The program CLAYEST may be run under DOS 3.0 or later. The program will also run in 

DOS under Windows 3.1/95/98, and directly under Windows by accessing it from the file manager 

or the Run command. The program is interactive and will prompt the operator for:

1.) A printer output on LPT1 (the default is no).

2.) Which to use: linear or exponential interpolation or the mean of both (the default 

is the mean). Linear interpolation may slightly over-estimate the amount of 

clay present; exponential interpolation may slightly under-estimate the 

amount of clay present.

3.) The smallest particle size in microns actually measured.

4.) The drive where the data file is located (the default is C:).

5.) The path name to the data file.

6.) What filename to read. If none is entered, the program will display a directory of 

the file names in the specified drive and path.

7.) The initial lab number. If no lab number is entered the program will search for 

the first lab number in the specified file.

8.) Whether to write the revised data to disk (the default is yes). The program will 

add an "E_", "L_", or "M_" depending on which interpolation is chosen to 

the beginning of the filename assigned to the output file. If the resultant 

filename is too long (over eight characters), the user will be warned and 

asked to enter a new file path and name.

9.) Whether the operator wishes to pause the scrolling screen after each record in

order to examine the data.

After processing is completed, the operator will be prompted to quit or continue.. If the 

operator chooses to continue, the program will restart with the first prompt.
Figure 1. An example of a typical hardcopy of revised grain size data generated by the program 

CLAYEST. Each record is printed on a separate page.

USGS Sediment Size Analysis 

Clay Fraction Estimation

Lab # 

AM071

Device 

VV

Gravel

Sand

Silt

Clay

Field ID Proj. ID 

28-02 LIS

Location Depth 

CT 36

Phi Size 

-5

-4

-3

-2

-1

0

1

2

3

4

5

6

7

8

9

10

11

12

13

14

Original Revised

% 1.38

% 71.98

% 19.61

% 7.04 

Total:

Cruise ID 

SEAX96017

Top

0

Original 

Freq. % 

0.00

0.00

0.00

0.00

1.38

2.74

6.76

32.32

23.10

7.06

2.88

5.57

5.69

5.46

4.16

2.27

0.61

0.00

0.00

0.00

1.35

70.56

19.21

8.88 

100.00

Requestor Latitude Longitude 

HKNEBEL 41.02000 -73.05609

Bottom Sample Weight 

2 28.9046

Revised Linearly 

Freq. % 

0.00

0.00

0.00

0.00

1.38

2.69

6.63

31.68

22.64

6.92

2.82

5.46

5.58

5.35

4.08

2.23

1.29

0.86

0.43

0.00

Revised Classification: SILTY SAND

Revised Statistics

Median =

Mean =

Standard Deviation =

Skewness =

Kurtosis =

2.34

3.32

2.77

0.55

0.47
Figure 2. An example of a typical output file generated by the program CLAYEST. Each record 

contains the original sample identifiers and revised: fraction percentages, verbal equivalent, 

and statistics; extrapolation method; and revised frequency percents. The first record is a 

list of the field names; the remaining records are examples of typical data. The fields are 

delimited by commas and the records are delimited by dollar signs.

Lab Number,Field ID,Project ID,CruiseID,Requestor,Month,Day, Year, Latitude, Longitude, Sample 

Device, Location, Depth, Top, Bottom, Sample Weight, % Sand, % Gravel, % Silt,% Clay,Verbal 

Equivalent, Median, Mean, StandardDeviation, Skewness, Kurtosis, Extrapolation Method, 13 Phi 

Fq. %, 12 Phi Fq. %, 11 Phi Fq. %, 10 Phi Fq. %, 9 Phi Fq. %, 8 Phi Fq. %, 7 Phi Fq. %, 6 Phi 

Fq. %, 5 Phi Fq. %, 4 Phi Fq. %, 3 Phi Fq. %, 2 Phi Fq. %, 1 Phi Fq. %, 0 Phi Fq. %,-! Phi Fq. 

%,-2 Phi Fq. %,-3 Phi Fq. %,-4 Phi Fq. %,-5 Phi Fq. %,

AM307,ELIS-38,LIS,JD96-1,LPOPPE,8,27,96,41.28087,-72.22818,VV,CT,9.0,0,2,39.0394,93.52, 

0.00, 4.68, 1.80, SAND, 2.75, 2.97, 1.27, 1.70, 18.20, LINEAR, 0.13, 0.25, 0.38, 0.49, 0.56, 

0.54, 0.75, 1.29, 2.11,28.34,59.93, 2.90, 1.40, 0.93, 0.00, 0.00, 0.00, 0.00, 0.00,$
2

w 

U 

c* 

tu 

cu

^

U 

2

tu

D

a

tu 

c*

U2

g

H 

U 

<

£± 

U.

_J 

<

Q

O 

_j 

-J

O 

U

EXTRAPOLATED PORTION

OF THE GRAIN SIZE

DISTRIBUTION

AVERAGED SOLUTION

LINEAR 

SOLUTION ANALYZED 

OF THE GRA

DISTRIBU

EXPONENTIAL SOLU

0 0.1 0.6

INCREASING GRAIN SIZE (MICRONS)

Figure 3. Stylized plot showing that portion of the clay fraction typically truncated by electro-resistance particle 

are examples of linear, exponential, and average of the linear and exponential extrapolations to 0.1 

boundary by the program CLAYEST. Errors inherent in the linear and exponential solutions have bee

effect.
Comments

Figure 3 presents a stylized plot showing that portion of the grain-size distribution typically 

truncated by electro-resistance particle size analyzers. Also shown in this figure are examples of 

linear, exponential, and averaged extrapolations of the clay fraction to 0.1 microns, the clay-colloidal 

clay boundary. Errors inherent in the linear and exponential extrapolations have been exaggerated 

within this figure to more clearly demonstrate their effect on the grain-size distribution. The linear 

extrapolation tends to slightly overestimate the amount of clay present in a typical distribution and 

may be used by operators to account for material within the colloidal fraction. The exponential 

extrapolation tends to slightly underestimate the amount of clay present in a typical distribution and 

may be used by operators who want data that represents the minimal amount of clay present. 

Although an operator may select linear extrapolation, exponential extrapolation, or the mean of both 

extrapolations, the mean usually provides the most accurate estimate and is, therefore, the 

recommended solution for most applications.

ACKNOWLEDGMENTS

Interested parties can obtain copies of the documentation and compiled and uncompiled 

versions of the software on 3.5" diskettes from the Woods Hole Field Center of the Coastal and 

Marine Geology Program, U.S. Geological Survey, Woods Hole, MA. We thank J. Commeau and 

R. Rendigs for reviewing this report.

REFERENCES

Berg, R.H., 1958, Electronic size analysis of sub-sieve particles by flowing through a small liquid 

resistor: ASTM Special Technical Publication 234, p. 245-258.

Blatt, Harvey, Middleton, Gerard, and Murray, Raymond, 1972, Origin of Sedimentary Rocks: 

Englewood Cliffs, New Jersey, Prentice-Hall, 634 p.

Coulter, W.H., 1957, High speed automatic blood cell counter and cell size analyzer: proceedings 

of the National Electronics Conference, p. 1034-1042.

Coulter Electronics, Inc., 1988, Fine particle application notes: Instruction manual, 72 p. 

Folk, R.L., 1974, The petrology of sedimentary rocks: Hemphill Publishing Co., Austin, 182 p.

Kranck, K. and Milligan, T.G., 1979, The use of the Coulter Counter in studies of particle-size 

distributions in aquatic environments: Bedford Institute Report Series Bl-R, v. 79, p. 7.

Krumbein, W.C., and Pettijohn, F.J., 1938, Manual of Sedimentary petrology, New York, Appleton, 

Century, and Crofts, 549 p.

McCave, I.N., and Jarvis, J., 1973, Use of the Model T Coulter Counter in size analysis of fine to

7
coarse sand: Sedimentology, v. 20, p. 305-315.

McCave, IN., and Syvitski, J.P.M., 1991, Principles and methods of particle size analysis, In J.P.M. 

Syvitski (ed.), Principles, Methods, and Applications of Particle Size Analysis, New York, 

Cambridge University Press, p. 3-21.

Milligan, T.G., and Kranck, Kate, 1991, Electro-resistance particle size analyzers, In J.P.M. Syvitski 

(ed.), Principles, Methods, and Applications of Particle Size Analysis, New York, Cambridge 

University Press, p. 109-118.

Muerdter, D.R., Dauphin, J.P., and steele, G., 1981, An interactive computerized system for grain 

size analysis of silt using electro-resistance: Journal of Sedimentary Petrology, v. 51, p. 647- 

650.

Poppe, L.J., Eliason, A.H., and Fredericks, J.J., 1985, APSAS - An automated particle size analysis 

system: U.S. Geological Survey Circular 963, 77 p.

Schlee, J., and Webster, J., 1967, A computer program for grain-size data: Sedimentology, v. 8, p.

45-54.

Schlee, J., 1973, Atlantic continental shelf and slope of the United States sediment texture of the 

northeastern part: U.S. Geological Survey Professional Paper 529-L, 64 p.

Sheldon, R.W., and Parsons, T.R., 1967, A practical manual on the use of the Coulter Counter in 

marine research: Toronto, Coulter Electronics, 66 p.

Shideler, G.L., 1976, A comparison of electronic particle counting and pipet techniques in routine 

mud analysis: Journal of Sedimentary Petrology, v. 42, p. 122-134.
APPENDIX A - INPUT FILE FORMAT FOR AN INDIVIDUAL RECORD

LABNUMBER,FIELDID,PROJECTID,CRUISEID,REQUESTOR,MONTH,DAY,YEAR,LATI

TUDE,LONGITUDE,DEVICE,LOCATION,DEPTH,TOP,BOTTOM,SAMPLE,WEIGHT,%SA

ND,%GRAVEL,%SILT,%CLAY,SEDCLASS,MEDIAN,MEAN,STDDEV,SKEWNESS,KURT

OSIS,MODE1,M1F%,MODE2,,M2F%,MODE3,M3F%,#MODES,11(1),CFP11,10,CFP10,9(1),CF

P9,8(|>,CI¥8Jc|),CI¥7,6(|>,CFP6,5(|>,C^

CFP-lr2(j>,CFP-2r3(|),CFP-3,-4(|>,CFP-4,-5(j),CFP-5,$
APPENDIX B-DATA DICTIONARY

The data dictionary contains an index of, and definitions for, the abbreviations and the sediment fiel

(Appendix A).

Short Field Name Full-Length Field Name 

LABNUMBER Unique sample identifier

FIELDID 

PROJECTID

CRUISEID

REQUESTOR

MONTH

DAY

YEAR

LATITUDE

LONGITUDE

DEVICE

LOCATION

DEPTH

TOP

BOTTOM

WEIGHT

%GRAVEL

%SAND

%S1LT

%CLAY

SEDCLASS 

MEDIAN

MEAN

Field identifier 

Project name

Cruise identifier

Requestor

Month

Day

Year

Latitude

Longitude

Sampling device

Area

Depth or sounding in meters

Top depth

Bottom depth

Sample weight

Percent gravel

Percent sand 

Percent silt 

Percent clay

Lithology

Median

Mean

Field Description

Unique identification number (consisting of 2 characters a

sedimentation laboratory of the USGS's Woods Hole F

Texture Database

Sample name or number as assigned in the field 

Project under which samples were taken, or data generated (i.e. 

indicates more specific area 

Cruise on which samples were colloected

Name of Principle Investigator or designee requesting analysis 

Month sample collected/analyzed 

Day sample collected/analyzed 

Calender year sample collected/analyzed

Latitude in decimal-degrees (south latitudes are depicted by ne

Longitude in decimal-degrees (west longitudes are depicted by 

Device used to collect the sample

General geographic area (general enough to easily locate on a s

Measured depth of water overlying sediment at sample time, no

Top depth of the sample below the sediment-water interface, in 

Bottom depth of the sample below the sediment-water interface

Weight of sample in grams 

Gravel content in percent dry weight of the sample (particles wi

mm; -l(j) and larger) 

Sand content in percent dry weight of the sample (particles with 

but greater than or equal to 0.0625 mm ; 0$ through 4

Silt content in percent dry weight of the sample (particles with 

mm, but greater than or equal to 0.004 mm; 5$ throug

Clay content in percent dry weight of the sample (particles with 

mm; 9(j> and smaller) 

General lithologic description 

Middle point in the grain size distribution in phi units 

Average value in the grain size distribution in phi units

10
STDDEV

SKEWNESS

KURTOSIS

MODE1

M1F%

MODE2

M2F%

MODE3

M3F%

#MODE

114)

CFP11

10(j>

CFP10

94)

CFP09

84)

CFP08

74)

CFP07

64)

CFP06

54)

CFP05

44)

CFP04

34)

CFP03

24)

CFP02

14>

CFP01

04)

CFPOO

-14>

CFP-1

-24)

CFP-2

Standard deviation (sorting)

Skewness

Kurtosis

Mode 1 class

Mode 1 strength

Mode 2 class

Mode 2 strength

Mode 3 class

Mode 3 strength

Number of modes

11 phi

Cfpll

10 phi

CfplO

9 phi

Cfp09

8 phi

CfpOS

7 phi

Cfp07

6 phi

Cfp06

5 phi

Cfp05

4 phi

Cfp()4

3 phi

Cfp03

2 phi

Cfp02

1 phi

CfpOl

Ophi

CfpOO

-Iphi

Cfp-1

-2 phi

Cfp-2

Standard deviation (root mean square of the deviations) of the 

Skewness (deviation from symmetrical form) of the grain size d

Kurtosis (degree of curvature near the mode) of the grain size d

First mode (particle size that occurs the most number of times) 

Modal strength of the first mode in percent

Second mode in phi units

Modal strength of the second mode in percent

Third mode in phi units

Modal strength of the third mode in percent

Number of modes

11 phi fraction

Cumulative frequency percent of the 11 (J> fraction and coarser 

10 phi fraction

Cumulative frequency percent of the 10(J) fraction and coarser 

9 phi fraction

Cumulative frequency percent of the 9(J> fraction and coarser (l

8 phi fraction

Cumulative frequency percent of the 8(J> fraction and coarser (l

7 phi fraction

Cumulative frequency percent of the 7(J> fraction and coarser (l

6 phi fraction

Cumulative frequency percent of the 6(J> fraction and coarser (l

5 phi fraction

Cumulative frequency percent of the 5(J> fraction and coarser (l

4 phi fraction

Cumulative frequency percent of the 4(J> fraction and coarser (

3 phi fraction

Cumulative frequency percent of the 3cJ> fraction and coarser (

2 phi fraction

Cumulative frequency percent of the 2cJ> fraction and coarser (

1 phi fraction

Cumulative frequency percent of the IcJ) fraction and coarser (

0 phi fraction

Cumulative frequency percent of the 0(J> fraction and coarser (

-1 phi fraction

Cumulative frequency percent of the-lcj) fraction and coarser (

-2 phi fraction

Cumulative frequency percent of the-2(J> fraction and coarser (

11
-3({) -3 phi -3 phi fraction

CFP-3 Cfp-3 Cumulative frequency percent -3([> fraction and coarser (less th

-4(|> -4 phi -4 phi fraction

CFP-4 Clp-4 Cumulative frequency percent -4(|> fraction and coarser (less th

-5((> -5 phi -5 phi fraction

CFP-5 Cfp-5 Cumulative frequency percent -5([> fraction and coarser

12
APPENDIX C PROGRAM DOCUMENTATION

Program Name: CLAYEST.EXE 

Type; Main program

Purpose; This program will calculate that portion of the clay fraction not detected during a particlesize analysis by a Coulter Counter (down to 13.0 phi). The operator may select linear or 

exponential extrapolation, or the mean of both.

Operating System; DOS version 3.0 or later, or DOS under Windows 3.1 or Windows 95

Source Language: Microsoft QuickBASIC version 4.0 or later 

Library-standard

Source Code; CLAYEST.BAS 

Program Category; Data processing

Input; The pre-database file created by GSTAT, the main data processing and statistics program 

used in the Sedimentation Laboratory of the Coastal and Marine Geology Program, Woods 

Hole Section. The format for this pre-database file and an example of a typical file is shown 

in Appendix A.

Output; When specified, this program will generate a hardcopy (Fig. 1) and/or a data file (Fig. 2).

The optional hardcopy will contain: the sample identifiers; the original and revised frequency 

percents; the original and revised percentages of gravel, sand, silt, and clay; a revised 

sediment classification (verbal equivalent); and the revised statistics. All sample 

identifiers of length greater than 11 characters will be truncated.

The optional data file will be comma delimited and contain: the sample identifiers; the 

revised percentages of sand, gravel, silt, and clay; the revised verbal equivalent; the 

revised statistics; and the revised frequency percentages (-5 to 13 phi). The first line 

of this data file will contain a list of the field attributes.

Usage; To use the program from DOS, change to the working directory and type: CLAYEST

To use the program from Windows 95/98, access it from the file manager or the Run 

command.

The program CLAYEST is interactive and will prompt the operator for: 

1.) A printer output on LPT1 (the default is no).

13
2.) Which to use: linear or exponential interpolation or the mean of both (the default 

is the mean). Linear interpolation may slightly over-estimate the amount of 

clay present; exponential interpolation may slightly under-estimate the 

amount of clay present.

3.) The smallest particle size in microns actually measured.

4.) The drive where the data file is located (the default is C:).

5.) The path name to the data file.

6.) What filename to read. If none is entered, the program will display a directory of 

the file names in the specified drive and path.

7.) The initial Lab Number. If no lab number is entered the program will search for 

the first lab number in the specified file.

8.) Whether to write the revised data to disk (the default is yes). The program will 

add an "E_", "L_", or "M_" depending on which interpolation is chosen to 

the beginning of the filename assigned to the output file. If the resultant 

filename is too long (over eight characters), the user will be warned and 

asked to enter a new file path and name.

9.) Whether the operator wishes to pause the scrolling screen after each record in 

order to examine the data.

Author; Eliason Data Services, Mashpee, MA 02649 (508)477-1400 

Maintenance: L. Poppe, USGS, Woods Hole, MA 02543

14
APPENDIX D - PROGRAM SOFTWARE

1 ** Name: CLAYEST97.BAS **

'** Type: MAIN [ X ] FUNCTION [ ] SUB [ ]  

' **  Language: Microsoft QuickBASIC 4.00 or above **

' ** Library: Standard **

' ** Arguments: None -MainModule **

1 ** Programmer: A. Eliason- Eliason Data- (508)477-1400 **

' ** Date: Version # in main module includes date **

1 ** as: MMDDYYYY **

' ** Description: Please see comments within this module **

I  

' **  Usage: If no info here please see the line **

' ** comments. **

1 Clay Extrapolation Program - CLAYES97.BAS

1 Last rev 97-06-24 AHE

DEFINT I-N 'use Fortran conventions

DEFDBL A-H, O-Z

DIM PHIP(32), FP(32), FPcor(32), FPcorlOO(32), Iden(4), F(2002)

DIM XFit(4), YFit(4)

DIM s(32), index(32)

DIMItext$(1000)

DIM Lable$(70)

DECLARE SUB CPI (PHIPO, FP(), N, DP, F(), IP, DPHI)

DECLARE SUB BPLOT (TMIN, TMAX, X(), N) 'EDS sub for FOPLOT

DECLARE SUB PressAnyKey ()

DECLARE SUB CSCOEF (N, X(), Y(), s(), index())

DECLARE SUB CSFITl (N, X(), Y(), s(), index(), XI, sx)

DECLARE SUB CSFIT2 (N, Xl(), Yl(), YPI, YPN, X, Y)

DECLARE SUB FFIT (N, X(), Y(), IType, A, B)

DECLARE FUNCTION U2Phi (US)

Version$ = "062497.01"

Diag% = 0

FileWrite% = 1

Cont% = 0

InterpFlag = 1

ISWFlag = 0 'no Schlee-Webster

15
Lprnt = 0

PhiNudge = .045 ' this factor is added to conversions from micron dia to Phi 

1 in order to adjust for the offset between whole micron 

1 diamaters and phi boundries eg: 1 micron actually= 9.965..Phi 

1 so we add .045 to Phi in order to cross the 10 Phi boundry.

NextFile:

CLS

PRINT

PRINT " USGS Sediment Size Analysis"

PRINT " Clay Fraction Estimation Version "; Version$

PRINT

PRINT "OPTIONS:"

INPUT " Do you want printed output (Y/N)? [N] ", A$

IF UCASE$(LEFT$(A$, 1)) = "Y" THEN Lprnt = 1

PRINT

ExType = 4 ' default to MEAN

InterpFlag = 1

1 lets not kill this yet rev 6-18-97

'INPUT " Do you want to use Extrapolation (Y/N)? [Y] ", A$

'IF UCASE$(LEFT$(A$, 1)) = "N" THEN InterpFlag = 0

' If NO is selected then only the original data are printed/output 

IF InterpFlag THEN

PRINT "Please READ the following carefully:"

PRINT

PRINT " You may elect to estimate clay content from the smallest measured size to"

PRINT " 13 phi. This may be done by LINEAR or EXPONENTIAL extrapolation, or the "

PRINT " MEAN of both sets of results. "

PRINT " Linear extrapolation may tend to over-estimate the clay present, while"

PRINT " exponential may under-estimate."

PRINT " It is up to the operator to determine the best method based"

PRINT " upon the original data."

PRINT

INPUT " Do you wish to use LINEAR, EXPONENTIAL, or MEAN (L/E/M)? [M] ", A$ 

IF UCASE$(LEFT$(A$, 1)) = "L" THEN ExType = 1 

IF UCASE$(LEFT$(A$, 1)) = "E" THEN ExType = 3 

1 else ExType = the default, MEAN, = 4 

PRINT

PRINT " The program needs to know the smallest measured size based upon" 

PRINT " the Coulter calibration. Normally this will be a value between" 

PRINT " 0.5 and 0.8 microns. You MUST enter a value." 

Size:

16
PRINT 

DO

INPUT " Smallest measured micron diameter"; UD 

LOOP WHILE UD <= 0 

IF UD < .5 OR UD > .8 THEN

PRINT USING " You have entered a value which corresponds to ###.## Phi"; U2Phi(UD)

INPUT " Is this correct (Y/N)? ", A$

IF UCASE$(LEFT$(A$, 1)) <> "Y" THEN GOTO Size: 

END IF

PhiMin = U2Phi(UD) ' actulal phi size 

PhiMin = PhiMin + PhiNudge ' correct to USGS boundry 

IntPhi = INT(PhiMin) ' Integer portion of smallest measured 

IntPhiFrac = IntPhi + 1 ' Fraction containing the smallest meas. 

NdxPhiCor = IntPhiFrac + 6 ' index to Fraction Array 

PhiCorFac = PhiMin - IntPhi ' divisor for correction 

PctPhiCor = 100 - INT(PhiCorFac * 10000) / 100 ' just for information

PRINT USING " The value you entered corresponds to ##.## Phi"; PhiMin

PRINT USING " The ## Phi Fraction will be increased in order to account"; IntPhiFrac

PRINT USING " for the ##.## percent perceived to be undetected in this fraction."; PctPhiCor

IF Daig% THEN PRINT " Value will be divided by "; PhiCorFac

PRINT

END IF

PressAnyKey

1 Read the input field lables

FOR I = 1 TO 68

READ Lable$(I)

IF Diag% = 1 THEN PRINT Lable$(I); " ~ "; 

NEXT I 

PRINT

1 new file read routine 9-16-96 ############ AHE

FilePath:

INPUT "Enter Drive Letter for data file - [C:] ", Drv$

IF Drv$ = "" THEN Drv$ = "C:"

IF LEN(Drv$) = 1 THEN Drv$ = Drv$ + ":"

IF LEN(Drv$) > 2 THEN GOTO FilePath

PRINT

PRINT "You will now be asked for a PATH then a FILENAME. If you don't know, leave"

PRINT "one or both blank and you will be shown a Directory screen."

INPUT "File PATH to read? - eg '\DATAV [\] ", FPN$

IF FPN$ = "" THEN FPN$ = "\"

17
IF LEFT$(FPN$, 1) <> "\" THEN FPN$ = "\" + FPN$ 

IF RIGHT$(FPN$, 1) <> T THEN FPN$ = FPN$ + T 

PRINT "Data Path = "; Drv$ + FPN$ 

INPUT "FILENAME to read? "; FLN$

IF FLN$ = "" THEN ' Guess s/he doesn't know 

FILES Drv$ + FPN$

PRINT "Please start again. When asked, enter a FILENAME from the list above." 

PRINT "You MUST include the extension (eg AAAAAA.PG)" '6-97 AHE 

PRINT "This program WILL NOT accept FILENAMES of more than 6 characters nor ext-" 

PRINT "tions greater than 3 characters" 

GOTO FilePath

END IF

LabNo:

INPUT "Initial LAB NUMBER eg AHE001' ? ", LabNo$

LabNo$ = UCASE$(LabNo$)

LnLen = LEN(LabNo$)

FoundFlag% = 0

Drv$ = UCASE$(Drv$)

FPN$ = UCASE$(FPN$)

FLN$ = UCASE$(FLN$)

PRINT "Attempting to open "; Drv$ + FPN$ + FLN$ 

OPEN Drv$ + FPN$ + FLN$ FOR INPUT AS #1

DO WHILE NOT EOF(l)

INPUT #1,A$

IF (UCASE$(LEFT$(A$, LnLen)) = LabNo$) OR (UCASE$(MID$(A$, 2, LnLen)) = LabNo$) 

THEN

IF LEFT$(A$, 1) = "$" THEN A$ = MID$(A$, 2, LnLen)

FoundFlag% = 1

EXIT DO 

END IF

IF Diag% THEN PRINT A$; " "; 

LOOP

IF FoundFlag% THEN 

PRINT 

PRINT " Found "; A$

LN$ = A$

18
ELSE

BEEP

PRINT

PRINT " Did NOT find "; LabNo$

PRINT " Please try again with another Lab Number." 

END IF 

IF FoundFlag% = 0 THEN

CLOSE #1

PRINT

GOTO LabNo 

END IF

INPUT "Do you wish to write extrapolated data to Disk? [Y]"; A$ 

IF UCASE$(LEFT$(A$, 1)) = "N" THEN FileWrite% = 0 

IF FileWrite% = 1 THEN

RevType$ = "M_" 'default- mean of both

IF ExType = 1 THEN RevType$ = "L_" 'linear

IF ExType = 3 THEN RevType$ = "E_" 'exponential

OutFile$ = Drv$ + FPN$ + RevType$ + FLN$ 

IF LEN(FLN$) > 10 THEN 

PRINT 

PRINT A 1V.LL l -L "nmnnmiMmmmnnHMmmn" ........................................

PRINT "CAN NOT Write output data to "; OutFile$ 

PRINT "!!!!!!!!!! File Name too Long !!!!!!!!!!" 

PRINT "You must specify a new output file Path/Name!" 

BEEP 

ELSE

PRINT "Output data will be written to the following Path/File," 

PRINT "unless you enter another complete pathname." 

END IF

PRINT "Output Filename [ "; OutFileS;" ] ";

INPUT "", A$

IF UCASE$(A$) = "Y" THEN A$ = ""

IF A$ <> "" THEN OutFileS = AS

PRINT "Attempting to open "; OutFileS; " for output."

OPEN OutFileS FOR OUTPUT AS #2

PRINT "Open sucessful."

' Write the attribute names to disk

FOR I = 1 TO 26

PRINT #2, Lable$(I); ",";

NEXT I

' attribute name # 27 is type of extrapolation

PRINT #2, "Extrapolation Method,";

19
' kill this LJP 10-24-96 PRINT #2, "Mode 1 Fq.%, Mode 2 Fq.%, Mode 3 Fq.%,Number of 

Modes,";

FOR I = 13 TO -5 STEP -1 ' -4 to -5 10/96 ARE

PRINT #2,1; "Phi Fq.%,"; 

NEXT I

PRINT #2,""

END IF 

PRINT

Cont% = 0 ' Correct error ARE 6-24-97

INPUT "Do you wish to pause and examine the data for each record? [Y]", A$

IF UCASE$(LEFT$(A$, 1)) = "N" THEN Cont% = 1 True

'PressAnyKey 

NextSamp:

Itext$(l) = LN$

1 = 2

DO

A$ = "": B$ = "" 

DO

DO ' ignore CR & LF

B$ = INPUT$(1,#1)

LOOP WHILE B$ = CHR$(10) OR B$ = CHR$(13) 

IF B$ = "," THEN EXIT DO ' end of field - exit 

A$ = A$ + B$ ' concatenate

IF B$ = "$" THEN EXIT DO ' because $ is not comma delimited 

LOOP

Itext$(I) = A$ 

1 = 1+1

IF Diag% THEN PRINT A$; "~"; 

LOOP UNTIL A$ = "$"

IF Diag% THEN 

PRINT : PRINT 

FOR I = 1 TO 68

PRINT USING "\ \&"; Lable$(I); Itext$(I) 

IF (I MOD (22) = 0) THEN PressAnyKey 

NEXT I 

PressAnyKey 

END IF

20
' check for proper field placement 

IF Itext$(68) <> "$" THEN

PRINT : PRINT

PRINT "End of Record Character ($) NOT IN PROPER LOCATION"

PRINT "Aborting Program"

BEEP: BEEP: BEEP

PRINT : PRINT

END 

END IF

' Now convert Cumulative Fq% to Fq%

FORj = 1TO32

PHIP(j) = j - 6: FP(j) = 0: FPcor(j) = 0: FPcorlOO(j) = 0 'init the arrays 

NEXTj

Cumul = 0: FqPct = 0: j = 1

FOR I = 67 TO 35 STEP -2 ' step from largest to smallest size

FqPct = VAL(Itext$(I)) - Cumul

FqPct = INT(FqPct * 100 + .499999999#) / 100' round to two places

Cumul = Cumul + FqPct

PHIP(j) = VAL(Itext$(I - 1)) ' Phi Size from input data

FP(j) = FqPct ' Actual phi percentage

IF Diag% THEN PRINT USING "## \ \ ##.## ###.##"; PHIP(j); Itext$(I); FP(j); Cumul

j=j + l

'*** PHIP(j) now contains -5phi at 1, -4 at 2,... 11 phi at 17 *** 

NEXT I 

PRINT

N = 20 ' -5 to 14 phi FRACTION inclusive 

DPHI = .2

XSand = VAL(Itext$(17)) 

Grav = VAL(Itext$(18)) 

Silt = VAL(Itext$(19)) 

clay = VAL(Itext$(20))

' test for total percent = 100 

CFP = 0 

FORj= 1 TON

CFP = CFP + FP(j) 

NEXTj

21
IF CFP > 101 OR CFP < 99 THEN

PRINT "Cumulative Frequency Percent is not 100 +/- 1"

PRINT

IF Lprnt THEN

LPRINT "Cumulative Frequency Percent is not 100 +/- 1"

LPRINT 

END IF 

Press Any Key 

GOTO NextSamp 

END IF

1 test if Phi intervals are even

DP = PHDP(2) - PHIP(l) '=.5

FORj = 3TON

Test2 = PHIPG) - PHIP(j - 1) '=.5

Test3 = ABS(DP - Test2) '=0

IF Diag% THEN PRINT DP, Test2, Test3: PRINT

IF (Test3 >.001) THEN

PRINT Iden$; " Phi intervals are not even"

PRINT

IF Lprnt THEN

PRINT Iden$; " Phi intervals are not even" 

PRINT 

END IF 

PressAnyKey 

GOTO NextSamp 

END IF 

DP = Test2 

NEXT]

'PressAnyKey

RealStuff:

IF InterpFlag = 0 THEN GOTO PrintOut

PRINT

PRINT "Doing the initial Extrapolation"

FOR I = 1 TO N ' copy the array

FPcor(I) = FP(I) 

NEXT I

FOR I = N TO 1 STEP -1 ' find the smallest non-zero fraction 

IFFP(I)>OTHEN

NZ = I ' NZ= index to smallest FRACTION with data

22
EXIT FOR 

END IF

NEXT I

' Make the percentage correction to the fraction containg the smallest

' measurable data

imin = NdxPhiCor

IF FP(imin) > 0 THEN

FPcor(imin) = FP(imin) / PhiCorFac

PRINT USING "The smallest measurable fraction is ## Phi with a value of ##.##%"; IntPhiFrac; 

FP(imin)

PRINT USING "##.## has been changed to ##.##"; FP(imin); FPcor(imin)

PRINT 

ELSE

PRINT USING "The ## Phi fraction was zero so it needed no change."; IntPhiFrac

PRINT 

END IF

' Now do Extrapolation from the smallest non zero fraction to 14 FRACTION 

1 which is zero. Don't forget the 13 FRACTION contains the 12.001 to 13 data

PRINT USING "The smallest non-zero fraction is ## Phi with a value of ##.##%"; PHIP(NZ);

FPcor(NZ)

IF ExType = 1 THEN PRINT "Linear ";

IF ExType = 3 THEN PRINT "Exponential ";

IF ExType = 4 THEN PRINT "Mean of L./E. ";

IF ExType <> 4 AND ExType <> 3 AND ExType <> 1 THEN PRINT "error at ExType. Program

halted.": STOP

PRINT "Extrapolation will now be performed to 13 Phi.(14 Phi Fraction=0)"

PRINT

' Set up for calls to FFTT (Linear and Exponential)

index = NZ

XFit(l) = index

YFit(l) = FPcor(index)

XFit(2) = index 

YFit(2) = FPcor(index)

XFit(3) = 20 

YFit(3) = .0000001

CALL FFTT(3, XFitQ, YFitQ, 1, Alin, Blin) 

CALL FFTT(3, XFitQ, YFitQ, 3, Aexp, Bexp)

23
index = NZ

'Npoints = 20 - NZ

'CLinFac = FPcor(index) / Npoints

FOR I = index+1 TO 19

IF ExType = 1 THEN FPcor(I) = Alin + Blin * I ' linear

IF ExType = 3 THEN FPcor(I) = Aexp * EXP(Bexp * I) ' logarithmic

IF ExType = 4 THEN FPcor(I) = ((Alin + Blin * I) + (Aexp * EXP(Bexp * I))) / 2

NEXT I

1 Re-normalize to 100 percent

TotalPct = 0

FOR I = 1 TO 20 ' N?

TotalPct = TotalPct + FPcor(I) 

NEXT I 

FOR I = 1 TO 20

FPcor 100(1) = 100 * FPcor(I) / TotalPct 

NEXT I

Printout:

PRINT

IF Cont% = 0 THEN PressAnyKey

1 Print the header data

CLS

PRINT " USGS Sediment Size Analysis"

PRINT " Clay Fraction Estimation Version "; Version$

PRINT

PRINT "Lab # Field ID Proj. ID Cruise ID Requestor Latitude Longitude"

PRINT

PRINT USING "\ \"; Itext$(l); Itext$(2); Itext$(3); Itext$(4); Itext$(5); Itext$(9); Itext$(10)

PRINT

PRINT

PRINT "Device Location Depth Top Bottom Samp. Weight"

PRINT

PRINT USING "\ \"; Itext$(l 1); Itext$(12); Itext$(13); Itext$(14); Itext$(15); Itext$(16)

PRINT

IF Cont% = 0 THEN PressAnyKey

PRINT

PRINT" Original Revised";

IF ExType = 1 THEN PRINT " Linearly"

24
IF ExType = 3 THEN PRINT " Exponentially"

IF ExType = 4 THEN PRINT " by Mean of Linear and Exponential"

PRINT " Phi Size Freq. % Freq. %" 

FORj = 1 TON

PRINT USING " ## ##.## ##.##"; PHIP(j); FP(j); FPcorlOO(j) 

NEXTj 

PRINT 

IF Lprnt THEN

' Print the original data to LPT1

LPRINT " USGS Sediment Size Analysis"

LPRINT " Clay Fraction Estimation Version "; Version$; ""

LPRINT

LPRINT" Lab# Field ID Proj. ID Cruise ID Requestor Latitude Longitude"

LPRINT" ";

LPRINT USING "\ \"; Itext$(l); Itext$(2); Itext$(3); Itext$(4); Itext$(5); Itext$(9); Itext$(10)

LPRINT

LPRINT" Device Location Depth Top Bottom Samp. Weight"

LPRINT" ";

LPRINT USING "\ \"; Itext$(l 1); Itext$(12); Itext$(13); Itext$(14); Itext$(15); Itext$(16)

LPRINT

LPRINT" Original Revised";

IF ExType = 1 THEN LPRINT " Linearly"

IF ExType = 3 THEN LPRINT " Exponentially"

IF ExType = 4 THEN LPRINT " by Mean of Linear and Exponential"

LPRINT " Phi Size Freq. % Freq. %"

LPRINT

FORj = 1TON

LPRINT USING" ## ##.## ##.##"; PfflP(j); FP(j); FPcorlOO(j) 

NEXTj 

LPRINT 

END IF

' test for total percent = 100 

CFP = 0 

FORj = lTON

CFP = CFP + FP(j) 

NEXTj 

IF CFP > 101 OR CFP < 99 THEN

25
PRINT "Cumulative Frequency Percent is not 100 +/- 1"

PRINT

IF Lprnt THEN

LPRINT "Cumulative Frequency Percent is not 100 +/- 1"

LPRINT 

END IF 

PressAnyKey 

GOTO NextSamp 

END IF

1 test if Phi intervals are even

DP = PHIP(2) - PHBP( 1) '=.5

FORj = 3TON

Test2 = PHBPG) - PHBP(j - 1) '=.5

Test3 = ABS(DP - Test2) =0

IF Diag% THEN PRINT DP, Test2, Test3: PRINT

IF (Tests >.001) THEN

PRINT Iden$; " Phi intervals are not even"

PRINT

IF Lprnt THEN

PRINT Iden$; " Phi intervals are not even" 

PRINT 

END IF 

PressAnyKey 

GOTO NextSamp 

END IF 

DP = Test2 

NEXTj

IF Cont% = 0 THEN PressAnyKey

Empl = PHBP(l) - .5 * DP 'Find center of 1st class interval

' to be used at end for statistics. 

IP = (PHBP(N) - PHBP( 1)) / DPHI +1.5 'Calculate # of Extrapolated points

'(Span/Interp_Interval)+1.5 =# points 

IF ISWFlag = 0 THEN

'PRINT "No S-W Extrapolation Requested"

'PRINT

'IF Lprnt THEN

' LPRINT "No S-W Extrapolation Requested"

' LPRINT

'END IF

IP = N ' Extrapolated Points = Original

DPHI = DP ' Delta Phi (interp) = Original

FOR I = 1 TO N

26
F(I) = FPcor 100(1) ' Extrapolated Data = Original 

'changed FP(i) to FPcor(i) 10-24-96 AHE

FP(I) = FPcorlOO(I) 'reload FP() 10-25-96 

NEXT I 

ELSE

PRINT "S-W Extrapolated Values will be used"

PRINT

IF Lprnt THEN

LPRINT "S-W Extrapolated Values were used"

LPRINT 

END IF 

FOR I = 1 TO N

FP(I) = FPcorlOO(I) 

NEXT I

GOTO SchleeWebster

CSCOEF ***************************************************

PRINT "Running CSCOEF...";

CALL CSCOEF(N, PHIPO, FPcor(), s(), indexQ)

PRINT " DONE"

FOR I = 1 TO N

PRINT s(I)

NEXT I 

Xl=PfflP(l) 

FOR I = 1 TO IP

'CALL CSFTT1(N, PHIP(), FPcor(), s(), indexQ, XI, SX)

'CALL CSFIT2(N, PHIPO, FPcor(), 0, 0, XI, SX)

F(I) = sx

XI =X1 +DPHI 

NEXT I

' Re-normalize to Phi Percentages 

DPI = DPHI/DP 

FOR I = 1 TO IP

IF F(I) < 0 THEN F(I) = 0

F(I) = F(I) * DPI

IF Diag% THEN 

PRINT I, F(I) 

IF (I MOD (15) = 0) THEN PressAnyKey

END IF 

NEXT I

'Scale the data for plotting 

Ymax = 0 

FOR I = 1 TO IP

IF F(I) > Ymax THEN Ymax = F(I)

27
NEXT I

IF Lprnt THEN LPRINT "Plot Using CSFIT1 - Natural Cubic Spline" 

'IF Lprnt THEN LPRINT "Plot Using CSFIT2 - Clamped Cubic Spline" 

CALLBPLOT(0!, Ymax + 1.5, F(), IP) 

IF Lprnt THEN LPRINT CHR$( 12) 

1 Ends test 

SchleeWebster:

1 move all the original values up by 1 in the arrays but leave first 

1 element as it was. 

FOR I = 1 TO N

'IF Diag% THEN PRINT FP(J + 1), PHIP(J + 1)

j=N-I

FP(j + 2) = FP(j + 1)

PHIP(j + 2) = PHIP(j + 1) 

NEXT I

PRINT

PRINT USING "### Points at ##.## Phi will be fit to ### points at #.### Phi"; N; DP; IP; DPHI

IF Lprnt THEN

LPRINT USING "### Points at ##.## Phi will be fit to ### points at #.### Phi"; N; DP; IP; 

DPHI

LPRINT 

END IF

CALL CPI(PHIP(), FPQ, N, DP, F(), IP, DPHI) 

' Re-normalize to Phi Percentages 

DPI = DPHI/DP 

FOR I = 1 TO IP

IF F(I) < 0 THEN F(I) = 0

F(I) = F(I) * DPI

IF Diag% THEN 

PRINT I, F(I) 

IF (I MOD (15) = 0) THEN PressAnyKey

END IF

NEXT I 

PRINT 

END IF

'Here is a cut at the Lloyd Breslaus Sediment Class Program

28
'IF Cont% = 0 THEN PressAnyKey 'REM'd out to elininate double PAC 6-97

PRINT " Calculating Descriptive Classification Percentages"

PRINT

' Get the new percentages

Gravel = 0

Sand = 0

Silt = 0

clay = 0

FOR I = 1 TO 5 'Gravel

Gravel = Gravel + FP(I) 

NEXT I 

FOR 1 = 6 TO 10 'Sand

Sand = Sand + FP(I) 

NEXT I 

FOR I =11 TO 14 'Silt

Silt = Silt + FP(I) 

NEXT I 

FOR I =15 TO 20 'Clay

clay = clay + FP(I) 

NEXT I

IF Sand = 0 THEN Sand - .001 

IF Gravel = 0 THEN Gravel = .001 

IF Silt = 0 THEN Silt = .001 

IF clay = 0 THEN clay = .001 

DO

IF Gravel >= 50 THEN C% = 1: EXIT DO

IF Gravel >= 10 THEN C% = 2: EXIT DO

XSand = Sand

Sand = Sand + Gravel

IF Sand >= 75 THEN C% = 3: EXIT DO

IF Silt >= 75 THEN C% = 4: EXIT DO

IF clay >= 75 THEN C% = 5: EXIT DO

' Sand, Silt, & Clay are all less than 75% 

SanSil = Sand / Silt 

ClySnd = clay/Sand 

SilCly = Silt / clay

IF Sand <= 20 THEN

IF SanSil > 1 THEN C% = 6: EXIT DO 

IF SilCly < 1 THEN C% = 7: EXIT DO 

IF ClySnd > 1 THEN

29
C% = 8: EXIT DO 

ELSE C% = 9: EXIT DO 

END IF 

ELSE ' Sand> 20%

IF clay <= 20 THEN

IF SanSil < 1 THEN C% = 9: EXIT DO 

IF SilCly > 1 THEN C% = 10: EXIT DO 

IF ClySnd > 1 THEN 

C% = 6: EXIT DO 

ELSE

C%= 11: EXIT DO 

END IF 

ELSE ' Clay > 20%

IF Silt > 20 THEN C% = 12: EXIT DO 

IF ClySnd > 1 THEN 

C% = 6: EXIT DO 

ELSE

C%= 11: EXIT DO 

END IF 

END IF 

END IF

PRINT "ERROR IN BRESLAU ROUTINE": BEEP: END ' ya can't get heah from theya 

LOOP

SELECT CASE C% 

CASE 1

C$ = "GRAVEL" 

CASE 2

C$ = "GRAVELLY SEDIMENT" 

CASES C$ = "SAND" 

CASE 4

C$ = "SILT" 

CASES C$ = "CLAY" 

CASE 6

C$ = "SANDY CLAY" 

CASE?

C$ = "SELTY CLAY" 

CASE 8

C$ = "CLAYEY SILT" 

CASE 9

C$ = "SANDY SILT"

30
CASE 10

C$ = "SILTY SAND" 

CASE 11

C$ = "CLAYEY SAND" 

CASE 12

C$ = "SAND SILT CLAY" 

CASE ELSE

C$ = "Error in classification routine" 

END SELECT

PRINT" Original Revised"

PRINT USING "Gravel % ##.## ##.##"; VAL(Itext$( 18)); Gravel 

PRINT USING "Sand % ##.## ##.##"; VAL(Itext$(17)); XSand 

PRINT USING "Silt % ##.## ##.##"; VAL(Itext$(19)); Silt 

PRINT USING "Clay % ##.## ##.##"; VAL(Itext$(20)); clay 

Total = Gravel + XSand + Silt + clay 

PRINT USING " Total = ###.##"; Total 

PRINT

PRINT "Revised Classification: "; C$ 

Verbal$ = C$ 

PRINT

IF Lprnt THEN

LPRINT" Original Revised"

LPRINT USING " Gravel % ##.## ##.##"; VAL(Itext$(18)); Gravel

LPRINT USING" Sand %##.## ##.##"; VAL(Itext$( 17)); XSand

LPRINT USING " Silt % ##.## ##.##"; VAL(Itext$(19)); Silt

LPRINT USING " Clay % ##.## ##.##"; VAL(Itext$(20)); clay

LPRINT USING " Total = ###.##"; Total

LPRINT

LPRINT " Revised Classification: "; C$

LPRINT

END IF

IF Cont% = 0 THEN PressAnyKey

'Scale the data for plotting 

Ymax = 0 

FOR I = 1 TO IP

IF F(I) > Ymax THEN Ymax = F(I) 

NEXT I

31
CALL BPLOT(0!, Ymax + 1.5, F(), IP)

STATS:

' now do the statistics via Schlee/Webster 1966

' Start with MODE

'IP = N ' need if not using schlee-webster

Dell=-.l

PRINT "Revised Statistics"

IF Lprnt THEN

LPRINT

LPRINT " Revised Statistics" 

END IF

FOR I = 2 TO IP

Del2 = F(I) - F(I - 1)

IF ((Del2 * Dell) < 0) AND (Del2 < 0) AND ((F(I - 1) - 5! * DPHI) > 0) THEN 

Xi = I - 2

Xem = Empl+Xi*DPHI 

Xf=F(I- 1) 

1 = 1+ 1

PRINT" Mode "; Xem, Xf 

IF Lprnt THEN 

'LPRINT

'LPRINT " Mode "; Xem, Xf 

END IF 

END IF 

Dell =Del2 

NEXT I

' MEDIAN

Sum = 0

FOR I = 1 TO IP

Sum = F(I) + Sum

IF (Sum - 50!) >= 0 THEN EXIT FOR 

NEXT I 

Xi = I - 1

Ct = Empl+Xi*DPHI

CtMed = Ct - (Sum - 50!) * DPHI / F(I) + .5 * DPHI 

Fmt$ \ \ =#####.##" 

PRINT USING FmtS; "Median"; CtMed 

IF Lprnt THEN LPRINT USING Fmt$; "Median"; CtMed

32
' Moments about the Mean

Sum = 0: Suml = 0: Sum2 = 0: Sum3 = 0: Sum4 = 0

FOR I = 1 TO IP

Xi = I - 1

Ct = Empl+Xi*DPHI

Sum = Sum + F(I)

Suml = Suml+F(I)*Ct

Sum2 = Sum2 + F(I) * Ct A 2

Sum3 = Sum3 + F(I) * Ct A 3

Sum4 = Sum4 + F(I) * Ct A 4 

NEXT I

Enl = Suml / Sum 

En2 = Sum2 / Sum 

EnS = SumS / Sum 

En4 = Sum4 / Sum 

Zm2 = En2-Enl A 2

ZmS = EnS - 3! * En2 * Enl + 2! * Enl A 3 

Zm4 = En4 + Enl * (-4 * En3 -f 6! * Enl * En2 - 3! * Enl A 3) 

DPHI2 = DPMI * DPMI

Zm4 = Zm4 - .5 * DPHI2 * Zm2 -f .02916667# * DPHI2 * DPHI2 

Zm2 = Zm2-DPHI2/12! 

'PRINT Enl; En2; En3; En4; Zm2; ZmS; Zm4

' Standard Deviation

Sigma = SQR(Zm2)

' Skewness

Skew = .5 * ZmS / (Sigma * Zm2)

' Kurtosis

ZKurt = Zm4/Zm2 A 2-3!

' Print those puppies

PRINT USING Fmt$; "Mean"; Enl

PRINT USING Fmt$; "Standard Deviation"; Sigma

PRINT USING Fmt$; "Skewness"; Skew

PRINT USING Fmt$; "Kurtosis"; ZKurt

PRINT

IF Lprnt THEN

LPRINT USING Fmt$; "Mean"; Enl 

LPRINT USING Fmt$; "Standard Deviation"; Sigma 

LPRINT USING Fmt$; "Skewness"; Skew 

LPRINT USING Fmt$; "Kurtosis"; ZKurt 

LPRINT CHR$( 12) ' Form Feed

END IF

33
' write the corrected data to disk 

IF FileWrite% THEN

FORI= 1TO 16

PRINT #2, Itext$(I); ","; ' Lab# - Samp Wt.

NEXT I

PRINT #2, USING "###.##,"; Sand; Gravel; Silt; clay;

PRINT #2, Verbal$;",";

PRINT #2, USING "###.##,"; CtMed; Enl; Sigma; Skew; ZKurt;

'Print #2, "9999,9999,9999,0," ' NEEDS WORK AHE - Killed LJP 10/96

IF ExType = 1 THEN PRINT #2, "LINEAR,";

IF ExType = 3 THEN PRINT #2, "EXPONENTIAL,";

IF ExType = 4 THEN PRINT #2, "MEAN_L_E,";

FORI= 19 TO 1 STEP-1

PRINT #2, USING "###.##,"; FPcorlOO(I); 'round to 2 dp'

NEXT I 

PRINT #2, "$" 

END IF

IF Cont% = 0 THEN PRINT "For next Sample -"; : PressAnyKey 

' See if the file is at end of data

DO WHILE NOT EOF(l)

INPUT #1,LN$

IF LN$ <> "" THEN GOTO NextSamp ' start next record 

LOOP 

CLS

CLOSE #1 

CLOSE #2 

PRINT

PRINT "End of File"

INPUT "Enter 'C to Continue or 'Q' to quit program. [Q] ", A$ 

IF UCASE$(LEFT$(A$, 1)) = "C" THEN RESTORE: GOTO NextFile

END

' Label Data

DATA Lab Number,Field ID,Project ID,Cruise ID,Requestor,Month,Day,Year

DATA Latitude,Longitude,Sample Device,Location,Depth,Top,Bottom,Sample Weight

DATA % Sand,% Gravel,% Silt,% Clay,Verbal Equivalent,Median,Mean

DATA Standard Deviation,Skewness,Kurtosis,Mode 1 ,Mode 1 Frequency %

DATA Mode 2,Mode 2 Frequency %,Mode 3,Mode 3 Frequency %,Number of Modes

34
DATA "11 Phi"," 11 Phi Cumulative Fq %","10 Phi","10 Phi Cummulative Fq %" 

DATA " 9 Phi"," 9 Phi Cumulative Fq %"," 8 Phi"," 8 Phi Cummulative Fq %" 

DATA " 7 Phi"," 7 Phi Cumulative Fq %"," 6 Phi"," 6 Phi Cummulative Fq %" 

DATA " 5 Phi"," 5 Phi Cumulative Fq %"," 4 Phi"," 4 Phi Cummulative Fq %" 

DATA " 3 Phi"," 3 Phi Cumulative Fq %"," 2 Phi"," 2 Phi Cummulative Fq %" 

DATA " 1 Phi"," 1 Phi Cumulative Fq %"," 0 Phi"," 0 Phi Cummulative Fq %" 

DATA "-1 Phi","-l Phi Cumulative Fq %","-2 Phi","-2 Phi Cummulative Fq %" 

DATA "-3 Phi","-3 Phi Cumulative Fq %","-4 Phi","-4 Phi Cummulative Fq %" 

DATA "-5 Phi","-5 Phi Cumulative Fq %",End Flag 

DATA ""

1 ** Name: BPLOT **

1 ** Type: MAIN [ ] FUNCTION [ ] SUB [ X ]  

1 **  Language: Microsoft QuickBASIC 4.00 or above **

' ** Library : Standard **

' ** Arguments: **

'** Programmer: A. Eliason- Eliason Data- (508)477-1400 **

' ** Date: Version # in main module includes date **

1 ** as: MMDDYYVV **

1 ** Description: Please see comments within this module **

 ## **

' ** Usage: If no info here please see the line ** 

' ** comments. **

1 ** Name: BPLOT **

'** Type: MAIN [ ] FUNCTION [ ] SUB [ X ]  

1 **  Language: Microsoft QuickBASIC 4.00 or above **

'** Programmer: A. Eliason- Eliason Data- (508)477-1400 **

1 ** Date: Version # in main module includes date **

1 ** as: MMDDYYVV **

1 ** Description: Please see comments within this routine **

 ## **

1 ** Usage: If no info here please see the line ** 

' ** comments. **

SUB BPLOT (TMIN, TMAX, X(), N) 

SHARED Lprnt, Cont%

35
M = 0

SPAN = TMAX - TMIN 

UNIT = SPAN / 70 

FOR I = 1 TO N

M = M+1

IF M = 20 THEN M = 0: IF Cont% = 0 THEN PressAnyKey

ISPC = (X(I) - TMIN) / UNIT

Fract = (I- l)/5-5

IFract = Fract - .4999

'PRINT USING "## ##.##### "; IFract; X(i);

PRINT USING "##.##### "; X(I);

PRINT SPC(ISPC); "*" '; ' SPC(70 - ISPC);

'IF Lprnt THEN

' LPRINT USING "## ##.##### " ; IFract; X(i);

' LPRINT SPC(ISPC); "*" '; ' SPC(70 - ISPC);

'END IF 

NEXT I 

PRINT 

END SUB

SUB CPl (H(), Q(), N, DELH, F(), M, DELHC) 

SHARED Diag%

= 0 

Q(N + 2) = 0

FORj=2TON

DELQ = (Q(j + 1) - Q(j)) / DELH

Y = 3!*Q<j)-3!*Q<j + l)

IMAX = DELH / DELHC '= number of Extrapolated points per interval 

FOR I = 1 TO MAX

JP = (Z - H(2)) / DELHC + 1 .4999 ' Prevent rounding error

'PRINT "*********** jp= " jp- " " INPUT A$

XK = (Z-H(j))/DELH

X = XK+1!

A = Q(j-l) + XK*(Y-Q(j-l))

FC = A + X * (Q(j) - A + .5 * (Q(j + 1) - 2! * Q(j) + A) * XK) 

IF (j - 2) = 0 THEN GOTO Contl2 

IF (j - N) o 0 THEN GOTO ContS 

Contl2:

' Test for Extrapolated value less than zero

'PRINT FC

IFFC<OTHEN

Z = H(j) + DELHC

36
IF G - 2) = 0 THEN

PRINT "Exponential Extrapolation at Beginning" 

GOTO ContV 

END IF

IF 0 - N) = 0 THEN

PRINT "Exponential Extrapolation at End" 

ELSE

GOTO ContS 

END IF 

ContT:

IF Q(j + 1) = 0 THEN Q(j + 1) = .000001 'prevent div by zero 

IF Q(j) = 0 THEN Q(j) = .000001 

SLPE = LOG(Q(j) / Q(j + 1)) / DELH 

FOR H = 2 TO MAX

JP = (Z - H(2)) / DELHC + 1 .4999

ANTI = SLPE * (Z - H(j))

F(JP) = Q(j) / EXP(ANTI)

Z = Z + DELHC 

NEXTH 

GOTOContl: 

END IF 

ContS: F(JP) = FC

Cont2: Z = Z + DELHC 

NEXT I

Contl:

NEXTj

F(JP+1) = 0 

IF Diag% THEN PRINT "CPI WAS RUN. EXITING CPI."

END SUB

SUB CSCOEF (N, X(), Y(), s(), index()) STATIC 

1 Copyright (c) 1993 Crescent Software 

' Cublic spline coefficients subroutine

37
'Input

' N = number of X and Y data points

' X() = array of X data points (N rows by 1 column)

' Y() = array of Y data points (N rows by 1 column)

1 Output

1 S() = array of cubic spline coefficients

' (N rows by 1 column)

' INDEXQ = array of indices (N rows by 1 column)

REDIM RHO(N), TAU(N)

NM1=N- 1

FOR I = 1 TO N

index(I) = I 

NEXT I

' ascending order data sort 

FORI= 1TONM1

FORj = IPlTON 

H = index(I) 

IJ = index(j)

IF (X(II) > X(IJ)) THEN 

ITEMP = index(I) 

index(I) = index(j) 

index(j) = ITEMP 

END IF 

NEXTj 

NEXT I

NM2 = N - 2

RHO(2) = 0#

TAU(2) = 0#

FORI = 2TONM1 

EM1 =index(I- 1) 

H = index(I) 

IIP1 =index(I+ 1)
Hi = x(npi)-x(n)

Temp = (HIMl / HI) * (RHO(I) + 2#) + 2# 

RHO(I+ !) = -!#/Temp

D = 6# * ((Y(IIP1) - Y(H)) / HI - (Y(H) - Y(IIM1)) / HIMl) / HI 

TAU(I + 1) = (D - HIMl * TAU(I) / HI) / Temp 

NEXT I

s(l) = 0#

s(N) = 0#

' compute cubic spline coefficients

FOR I = 1 TO NM2

EB = N-I

s(EB) = RHO(EB + 1) * s(EB + 1) + TAU(EB + 1) 

NEXT I

ERASE RHO, TAU

END SUB

SUB CSFITl (N, X(), Y(), s(), index(), XI, sx) STATIC 

' Copyright (c) 1993 Crescent Software

' Cublic spline Extrapolation subroutine 

' Input

' N = number of X and Y data points 

' X() = array of X data points (N rows by 1 column) 

' Y() = array of Y data points (N rows by 1 column) 

1 SQ = array of cubic spline coefficients

(N rows by 1 column)

1 INDEX() = array of indices (N rows by 1 column) 

' XI = X data value to fit

1 Output

' SX = cubic spline Extrapolated value for X

FOR I = 2 TO N

H = index(I)

39
IF (XI <= X(II)) THEN EXIT FOR 

NEXT I

L = I- 1

IL = index(L) 

HJP1 =index(L+ 1)

A = X(ILP1)-X1 

B = XI - X(IL)

HL = X(ILP1)-X(IL)

sx = A * s(L) * (A * A / HL - HL) / 6# + B * s(L + 1) * (B * B / HL - HL) / 6# + (A * Y(IL) + B 

* Y(ILP1))/HL

END SUB

SUB CSFTT2 (N, Xl(), Yl(), YPI, YPN, X, Y) STATIC 

1 Copyright (c) 1989, 1990 Crescent Software

1 Clamped cubic spline Extrapolation subroutine 

'Input

' N = number of X and Y data points 

' Xl() = vector of X data points ( N rows ) 

' Yl() = vector of Y data points ( N rows ) 

' YPI = derivative at data point 1

YPN = derivative at data point N 

' X = X data point to fit

' Output

' Y = Extrapolated Y data point

1 NOTE: Xl(l) < Xl(2) < Xl(3) < ... < X1(N)

REDEVI U(N), YPP(N)

YPP(l) = -.5

= (3#/ (XI (2) - Xl(l))) * ((Yl(2) - Yl(l)) / (Xl(2) - Xl(l)) - YPI)

FOR I = 2 TO N - 1

SIG = (X1(I)-X1(I- 1))/(X1(I+ 1)-X1(I- 1))

40
P = SIG* YPP(I- l) + 2# 

YPP(I) = (SIG- 1#)/P

U(I) = (6# * ((Y1(I + 1) - Y1(I)) / (XI(I + 1) - X1(I)) - (Y1(I) - Y1(I - 1)) / (XI(I) - X1(I - 1))) 

/ (XI(I + 1) - X1(I - 1)) - SIG * U(I - 1)) / P 

NEXT I

QN = .5

UN = (3# / (X1(N) - X1(N - 1))) * (YPN - (Y1(N) - Y1(N - 1)) / (X1(N) - X1(N - 1)))

YPP(N) = (UN - QN * U(N - 1)) / (QN * YPP(N - 1) + 1#)

FORK = N-1TO1 STEP-1

YPP(K) = YPP(K) * YPP(K + 1) + U(K) 

NEXTK

' Extrapolate

KLO= 1 

KHI = N

WHILE (KHI - KLO > 1) 

K = (KHI + KLO) / 2# 

IF(X1(K)>X)THEN

KHI = K 

ELSE

KLO = K 

END IF 

WEND

H = X1(KHI)-X1(KLO)

A = (X1(KHI)-X)/H 

B = (X-X1(KLO))/H

Y = A * Yl(KLO) + B * Yl(KHI) + ((A A 3 - A) * YPP(KLO) + (B A 3 - B) * YPP(KHI)) * (H

A 2) / 6#

ERASE U, YPP 

END SUB 

SUB FFIT (N, X(), Y(), IType. A, B) STATIC

1 Function curve fit subroutine

41
1 1nput

1 N = number of data points (N >= 3)

' X() = vector of X data points ( N rows )

1 Y() = vector of Y data points ( N rows )

1 FTYPE = type of curve fit

1 1 = linear Y = A + B * X

1 2 = logarithmic Y = A + B * LOG(X)

1 3 = exponential Y = A * EXP(B * X)

' Output

1 A, B = coefficients of curve fit

Xl=0# 

Yl=0# 

Z = 0# 

X2 = 0# 

Y2 = 0#

SELECT CASE IType 

CASE1

' Linear

FOR I = 1 TO N

X1=X1+X(I) 

Yl = Yl + Y(I) 

X2 = X2 + X(I) * X(I) 

Y2 = Y2 + Y(I) * Y(I) 

Z = Z + X(I) * Y(I) 

NEXT I

X1=X1/N 

Y1=Y1/N

B = (Z - N * XI * Yl) / (X2 - N * XI * XI) 

A = Y1-B*X1 

CASE 2

' Logarithmic 

FOR I = 1 TO N

XLOGX = LOG(X(I))

XI =X1 + XLOGX

Y1=Y1+Y(I)

X2 = X2 + XLOGX * XLOGX

Y2 = Y2 + Y(I) * Y(I)

Z = Z + XLOGX * Y(I)

42
NEXT I

X1 = X1/N 

Y1=Y1/N

B = (Z - N * XI * Yl) / (X2 - N * XI * XI) 

A = Y1 -B *X1 

CASES

' Exponential 

FOR I = 1 TO N

XLOGY = LOG(Y(I))

X1=X1+X(I)

X2 = X2 + X(I) * X(I)

Y1=Y1+ XLOGY

Y2 = Y2 + XLOGY * XLOGY

Z = Z + X(I) * XLOGY 

NEXT I

X1=X1/N 

Y1=Y1/N

B = (Z - N * XI * Yl) / (X2 - N * XI * XI) 

A = EXP(Y1-B*X1) 

END SELECT

END SUB

DEFSNG A-H, O-Z

SUB PressAnyKey

X% = POS(O)

Y% = CSRLIN

PRESfT "Press any key to continue.";

DO

LOOP WHILE INKEYS = ""

LOCATE Y%, X%

PRINT" ";

LOCATE Y%, X%

END SUB

DEFDBL A-H, O-Z 

FUNCTION U2Phi (Usize) 

' convert micron size to phi

phi = -2.12720466241656D-16 + (-1.44269504088896#) * LOG(Usize / 1000)

43
'PRINT USING "###.######### mm = "; Usize; 

'PRINT USING "###.### phi"; Phi 

U2Phi = phi 

END FUNCTION

44