EPA-650/2-74-113
JULY 1974
Environmental Protection Technology Series
iiiiiiiiiii
3)
o
V
UJ
O
-------
RESEARCH REPORTING SERIES
Research reports of the Office of Research and Development, U.S. Environ-
mental Protection Agency, have been grouped into series. These broad
categories were established to facilitate further development and applica-
tion of environmental technology. Elimination of traditional grouping was
consciously planned to foster technology transfer and maximum interface
in related fields. These series are:
1. ENVIRONMENTAL HEALTH EFFECTS RESEARCH
2. ENVIRONMENTAL PROTECTION TECHNOLOGY
3. ECOLOGICAL RESEARCH
4. ENVIRONMENTAL MONITORING
5. SOCIOECONOMIC ENVIRONMENTAL STUDIES
6. SCIENTIFIC AND TECHNICAL ASSESSMENT REPORTS
9. MISCELLANEOUS
This report has been assigned to the ENVIRONMENTAL PROTECTION
TECHNOLOGY series. This series describes research performed to
develop and demonstrate instrumentation, equipment and methodology
to repair or prevent environmental degradation from point and non-
point sources of pollution. This work provides the new or improved
technology required for the control and treatment of pollution sources
to meet environmental quality standards.
Copies of this report are available free of charge to Federal employees,
current contractors and grantees, and nonprofit organizations - as
supplies permit - from the Air Pollution Technical Information Center,
Environmental Protection Agency, Research Triangle Park, North
Carolina 27711; or, for a fee, from the National Technical Information
Service, Springfield, Virginia 22161.
-------
EPA-650/2-74-113
REMOTE SENSING
OF POLLUTANTS
COMPUTERIZED REDUCTION OF LONG-PATH ABSORPTION DATA
by
V. E. Derr, M. H. Ackley,
M. J . Post, and R. F. Calfee
Wave Propagation Laboratory
Environmental Research Laboratories
NOAA
Boulder, CO 80302
IAG No, EPA-IAG-077(D)
ROAP TTo. 26AAP
Program Element No. 1AA010
EPA Project Officer: Dr. H. M. Barnes, Jr.
Chemistry and Physics Laboratory
National Environmental Research Center
Research Triangle Park, North Carolina 27711
Prepared for
OFFICE OF RESEARCH AND DEVELOPMENT
U.S. ENVIRONMENTAL PROTECTION AGENCY
WASHINGTON, D.C. 20460
July 1974
-------
This report has been reviewed by the Office of Research and Development, Environ-
mental Protection Agency, and approved for publication. Approval does not
signify that the contents necessarily reflect the views and policies of the Environ-
mental Protection Agency, nor does mention of trade names or commercial products
constitute endorsement or recommendation for use.
Publication No. EPA-650/2-74-113
-------
PREFACE
This report is arranged in a flexible digital format. The
sections and subsections are designated by adding a decimal point and
a number. For example, 5.1.2 is the second sub-subsection of the first
subsection of section 5. References to sections are given without the
use of the word section; for example, (7.3) refers the reader to the third
subsection of section 7. Equations are invariably and uniquely
referred to by a colon notation; that is, (3:7) refers the reader to
equation seven of section 3. Figures and tables have a unique reference
also; for example, F3.1 and T4.7 are respectively the first figure of
section 3 and the seventh table of section 4. Appendixes are indicated by
a capital A; that is, A2.1 is the first section of the second appendix.
Figures, tables, and equations of appendixes are denoted respectively by
FA3.3, TA2.4, and A3:2.
The International System of Units is used throughout except where
the traditional unit still holds so strongly that confusion would result
from absolute purity. Pollutants for example are often measured in parts
per million. Conversion factors are given in A4.
This report is the product of a joint effort by EPA and NOAA. It
conforms to EPA format rather than that of NOAA.
iii
-------
TABLE OF CONTENTS
Abstract x
1. Introduction 1
1.1 The Physical Problem 4
1.2 The Mathematical Problem 9
2. Present Status of the Library of Absorption Spectra 19
3. Best-Fit Calculations of Concentrations From a 22
Simulated Complex System
3.1 Accuracy of Concentration Calculations 22
3.1.1 Calculation of Concentrations for 32
Spectra Without Noise
3.1.2 Effect of Finite Signal-to-Noise Ratio (S/N) 33
3.1.3 Effects of "Maverick" Gases in the Sample 34
4. Calculation of the Selected Gas Sample Concentrations 39
From Spectrophotometer Data
4.1 Effect of Mismatched Resolution of Experimental 42
Data and Library Spectra
4.2 Calculation of the Concentration of CO 50
4,3 Calculation of the Concentration of SO- in the 56
Presence of CH.
4.4 Calculation of the Concentration of CH. 57
5. Conclusions and Recommendations 63
6. Acknowledgements 65
7. List of Symbols 66
8. References 67
Al. Sources and Contents of the Library of Spectra 71
Al.l ' Motivation for the Library 71
A1.2 Sources and Content 71
A1.3 Programs Used in Processing Spectra 91
-------
Page
A2. Program EPAGAS Documentation and Analysis 111
A2.1 Program Summary 111
A2.1.1 List of Fortran Symbols and Descriptions 112
A2.1.2 List of Program and Subprogram Descriptions 121
A2.2 Simulation of Gaussian Noise 170
A2.3 Constituent and Transmission Error Formulation 171
A2.4 Interpolation Formulas Used in Transmission Library 173
A2.5 Program Execution Techniques 177
A2.5.1 Theoretical Gas Library Spectra 177
A2.5.2 Experimental Transmission Data 178
A2.5.3 Minimization of "Goodness of Fit". Controls 179
A2.5.4 Error Parameters 180
A2.5.5 Miscellaneous Execution Techniques 182
A3. Characteristics of the ROSE Long-Path Spectrophotometer 185
System
A4. Convenient Conversion Factors 189
vi
-------
FIGURES
Page
Fl.l(a-c) Infrared Spectra of Selected Atmospheric Constituents n
F1.2 Infrared Spectrum of the Normal Atmosphere 14
F1.3 A Portion of the Transmission Spectra of CH., HNO_, CH., and
H20 and Their Composite Spectra 15
F3.1 Principal Regions of Absorption of Selected Normal and
Pollutant Atmospheric Constituents OQ
F4.1 Spectrum of S02 and CH4 Obtained by FTS (1300.8-1310.5 cm'1) 40
F4.2 Spectrum of S02 and CH4 Obtained by FTS (1310.6-1320.4 cm"1) 41
F4.3 Spectrum of S02 and CH4 Obtained by FTS (1320.5-1330.5 cm"1) 41
F4.4 Comparison of a Theoretical Spectrum (High Resolution) of Water
Vapor and the Spectrum Observed Through a Spectrometer With
1.0 cm"1 Slit Width 43
F4.5 Experimental Spectrum, CO, Nominal Resolution 2.5 cm" Obtained
on ROSE System 44
F4.6 Calculated Spectrum, CO, Resolution 3.0 cm" 44
F4.7 Calculated Spectrum, CO, Resolution 1.0 cm" 45
F4.8 Calculated Spectrum, CO, Resolution 2.3 cm" 45
F4.9(a-c) Experimental Spectrum of Combination of SO and CH., Nominal
-1
Slit Width 1.0 cm (ROSE System) 47
F4.10(a-c) Calculated Spectrum, S02 and CH., Resolution 1.0 cm" 4g
F4.ll Spectrum of CO Calculated by Program DEGRADE (Resolution
2.3 cm ) 52
F4.12 Spectrum of CO, Experimental Data From ROSE Spectrophotometer ^
F4.13 Computed Best-Fit Transmission of CO (EPAGAS Program),
Resolution 3.0 cm"1 c/
F4.14 Computed Best-Fit Transmission of CO (EPAGAS Program),
Resolution 2.3 cm"1 55
F4.15 Change in Spectrum of CO for Small Change of Concentration
(Calculated by Program DEGRADE) 56
F4.16 Computed Best-Fit Transmission of CH Compared with the Given
Experimental Transmission 60
vii
-------
FA1.1 Infrared Absorbing Regions of Atmospheric Gases 79
FA3.1 Resolution vs. Wavelength for Various Slit Width Settings
for ROSE Spectrophotometer 187
viii
-------
TABLES
Tl.l Annual Average Concentrations of Important Pollutants in
Cities (1968) 4
T2.1 Gases Available in the Spectral Library 21
T3.1 Calculation of Gas Concentrations From Simulated Composite
Spectra (Set 1) (Infinite Signal-to-Noise Ratio) 23
T3.2 Calculation of Gas Concentrations From Simulated Composite
Spectra (Set 2) (Infinite Signal-to-Noise Ratio) 24
T3.3(a-c) Calculation of Gas Concentrations From Simulated Composite
Spectra (Set 1) (S/N=10), (S/N=30), (S/N=100) 25
T3.4 Calculation of Gas Concentrations From Simulated Composite
Spectra (Set 2) (S/N = 10, 30, 100) 28
T3.5(a-b) Gas Concentrations From Simulated Composite Spectra (Set 1)
(Maverick Gases) 36
T3.6 Gas Concentrations From Simulated Composite Spectra (Set 2)
(Maverick Gases) 38
T4.1 Calculation of the Concentration of a Sample of CO by
Application of EPAGAS to an Experimental Spectrum at Several
Spectral Regions and Two Instrumental Resolutions 51
T4.2 Calculation of the Concentration of a Sample of S02 in the
Presence of CH. by Application of EPAGAS to an Experimental
Spectrum 58
T4.3 Calculation of the Concentration of a Sample of CH4 and H20
by Application of EPAGAS to an Experimental Spectrum 62
TA1.1 Spectra From the WPL Data Set 72
TA1.2 Spectra Processed From the Scientific Literature 75
ix
-------
ABSTRACT
Atmospheric gaseous pollutants are very numerous in industrial
regions. It is estimated that 25 or more pollutant molecules may be
found in the atmosphere in significant quantities. The measurement
of the concentration of each gas from the complex spectrum obtained by
a long-path infrared spectrophotometer requires the fitting of trial
spectra composed from a library of spectra. The fitting procedure
adjusts the concentrations of the trial spectra until a "best fit" in a
least-squares sense is produced. This report is a description of the
physical, mathematical, and calculational principles and procedures for
the use of a digital computer program to determine concentrations of
atmospheric gases in a path of a few kilometers. Detailed instructions
for the computer program and a library of spectra are provided.
-------
1. INTRODUCTION
The Chemistry and Physics Laboratory of the Environmental Protection
Agency (EPA) in Raleigh, NC, has been studying the feasibility of various
methods of measuring the gaseous pollutants of the earth's atmosphere. One
of the most important and promising of these methods is the use of long-path
infrared absorption by gases of the atmosphere. This method, useful only
for gas concentration measurement and not for measuring the concentration
of aerosols, employs a beam of radiation of wide bandwidth from a stable
source strong in infrared radiation. The collimated beam passes through
a part of the atmosphere to be studied, typically about 1 to 4 km. It
may be received at the end of that path or returned by reflectors to a
receiver near the transmitter. The receiver contains an infrared spectrometer
Many methods other than long-path infrared absorption have been
considered for the analysis of atmospheric gases. Much effort has
been expended to measure accurately the concentrations of pollutants near
cities and in the ambient atmosphere far from cities to obtain background
concentrations. Other methods include Raman spectroscopy, the use of
tunable laser sources, and observation of the ultraviolet electronic spectra
of molecules. These methods may either provide the entire analysis or
contribute partially. However, in this report, we confine our attention
to the use of long-path infrared absorption spectra and attempt to exploit
this method as fully as time and resources allow. It is the authors'
belief that other methods will prove important in the future, but the long-
established science of infrared absorption spectroscopy allows important
progress to be made very quickly so that at least preliminary data may be
obtained on the pollutants in the atmospheres of cities. It is probable
that combinations of several methods will be utilized.
The analysis of the composite infrared spectra prduced by many atmos-
pheric constituents is laborious and difficult. It requires a large computer
and a library of spectra. The Wave Propagation Laboratory of the National
Oceanic and Atmospheric Administration (NOAA) in Boulder, CO, has worked
-------
for many years in a joint program with EPA to aid in the analysis of such
composite spectra. The EPA financed portion of the work has been done under
contract number EPA-IAG-077(D) entitled, "Remote Sensing of Pollutants."
This report is specialized in several ways. Firstly, it makes
the assumption that the instrument used in obtaining the data is that
described in 1.1. This assumption is necessary in order to be specific
concerning the parameters used in our computational method. Without
such an assumption, it would be necessary to consider wide ranges of
parameters of a long-path spectrophotometer. The complexity caused by
this would not clarify the method and would complicate the analysis. It
is anticipated that this specialization will not hamper the'extension of
the method to other specific devices. This has been accomplished by
letting the input parameters to the computer program which analyzes the
complex composite spectra to be at the discretion of the user. Thus, the
choice of a different instrument would merely be reflected in different
input parameters to the calculation. A second way in which the work has
been limited is by the choice of spectra placed in the library. This
choice has been dictated by the impossibility of obtaining a complete
library of spectra with the resources available.
Fortunately, for the immediate application of the method, a relatively
small library of spectra is probably quite adequate. The spectra of the
gases of importance in the Los Angeles, CA, atmosphere, the Raleigh, NC,
atmosphere, and the atmosphere near refineries have btfen given preference
whenever possible. However, the library of spectra has been limited by the
unavailability of suitable spectra of some gases. Spectra were obtained
from the scientific literature whenever possible. In certain cases, spectra
were determined by laboratory experiment when they were not available in
the literature. The entries and sources of the spectral library are given
in Al.
The major normal atmospheric gases, fLO, CO-, N20, and CH., are part
of the population. (Oxygen and nitrogen have no significant spectra in the
frequency range of interest.) Normal constituents and pollutants are
discussed in 2. The listed gases were chosen because of their presence
-------
in the atmosphere and because detailed spectroscopic parameters are available
for them, including either line positions, strengths, and line-broadening
constants, or high-resolution experimental spectra.
The spectra of molecules have long been used as a definitive
identifier, and when the spectrophotometer is properly calibrated, molecular
concentrations can also be measured. In many cases, however, the method
has been applied to relatively large concentrations in the laboratory,
to pure gases, or at most to a mixture of a few constituents. It is
true that laboratory spectrometers have been used to obtain low-level
concentrations of a few species in the midst of large concentrations of
other gases, but this method can only be successfully performed by simple
observation when the regions of significant spectral absorption are well
separated. Because most organic pollutant molecules have rich infrared
spectra, this spectral region may be used for detection and measurement, but
the potential for interference between overlapping spectral lines is high.
Thus, the measurement of atmospheric pollutants is complicated by low
concentrations and by the simultaneous presence of many gases in temporally
and geographically varying amounts. Because low concentrations may be
significant for biological well-being, the method of analysis must be
sensitive. Table 1.1 (Burriss et al., 1972)* shows typical concentrations
of some pollutants in several typical cities. The sensitivity must be
sufficient to determine such concentrations and, perhaps, a factor of 10
less. The interference due to aerosols is not considered in this report.
Because hydrocarbons (HC) are a large percentage of the total atmos-
pheric pollutant content of many cities, they need special attention. In
Chicago, IL, in 1968, the average HC concentration is 5380 yg/m , compared
3
with CO at 7750 yg/m . The molecules making up the total HC content of
the atmosphere for any location are not clearly known because many current
analyzes separate only methane and non-methane HC. However, if there is
no outstanding single HC and the concentration of HC is distributed among
many different molecules, the weak absorption of each will provide
* References are to authors and year of publication.
-------
Table 1.1
ANNUAL AVERAGE CONCENTRATIONS
OF IMPORTANT POLLUTANTS IN CITIES (1968)
City
Chicago, IL
Cincinnatti, OH
Philadelphia, PA
Denver, CO
St. Louis, MO
Washington, DC
CO
7750
7000
-
6750
5750
4250
HC*
5380
4220
2690
5380
7300
2690
NO,
98
59
78
78
39
98
NO
94
-
67
54
40
54
so2
340
57
230
29
86
110
*HC = hydrocarbons
3
Units: pg/nf
(From Burriss et al., 1972)
effectively a small absorption background which can be ignored. This
problem can be solved definitively only after a detailed analysis of
constituents is obtained for a particular region, and the spectra of these
constituents are known.
Thus, the considerations discussed above require the analysis of a
complex spectrum of many gases of relatively low amounts. This report is
a description of the physical, mathematical, and computational procedures
necessary to obtain the concentrations of many gases contained in an
atmospheric sample by long-path infrared spectroscopy. The types of
gases present in the sample are assumed known; the concentrations are
unknown.
1.1 THE PHYSICAL PROBLEM
We consider an atmosphere composed, perhaps inhomogeneously, of
normal and pollutant constituents. The pressure is assumed to be one
f 7
atmosphere (1.02 x 10 dynes/cm ). We assume that the constituents do not
-------
significantly change during the course of an observation. The observation
is through a portion of the atmosphere approximately 61 cm in diameter and
1 to 4 km long. (These dimensions arise from the characteristics of a
long-path spectrophotometer (see A3).) To focus on a specific problem,
we choose the Research Triangle Area of North Carolina and a representative
concentration of pollutants for that area will be used when available.
For those pollutants whose average concentrations in Raleigh are not known,
the known values in other comparable areas will be used. If no reliable
values of trace materials can be obtained, they will be assumed to lie in
the few parts per billion range. Further information may be found in Al,
where TA1.1 and TA1.2 show the lowest concentrations assumed.
Although the method to be described is quite general, for clarity
we will use the parameters of a specific system for remote sensing of
pollutants (Streiff and Claysmith, 1972). This system, when used in an
absorption mode, is a field-type infrared (IR) scanning spectrophotometer
with a blackbody light source collimated by a telescope for long-path
operation. Some system parameters of interest are given in A3.
The equipment can be used as a receiving spectrophotometer to study
emissions from any sufficiently hot source, since it was designed to measure
the emission spectra of gases in smokestacks. The methods developed herein
are applicable to the study of emissions by substituting a library of
emission spectra in place of the library of absorption spectra.
Thus, by means of the long-path spectrophotometer, we obtain the
superimposed spectra of the atmospheric constituents in the spectral ranges
3 to 5 urn (3330-1820 cm"1) and 7 to 13.5 urn (1430-740 cm"1). The 5 to
7 urn spectral range is not used because of the very strong absorption of
water vapor. It is the purpose of the method described in this report to
obtain from the superimposed spectra the concentrations of the molecular
species present in the path.
Each of the gases present in the path acts independently of the
others according to the Bouguer-Beers law. That is, a small thickness of
-------
gas absorbs a fractional amount, = = -K'dx, which implies
I - Ioe~K'X (1:1)
where I is the initial intensity, K1 is the absorption coefficient, and
x is the absorber thickness. (Physically, Kf is the reciprocal of the
distance at which the intensity is 1/e of the initial intensity.) If the
first absorber is followed by a second, equal to the first, it absorbs an
equal fraction of the radiation incident upon it. If we consider a plane
wave impinging on "slabs" of absorbing gases, the effect is the same
whether the slabs follow one another or the gases are all placed in one
slab. If the gases have absorption coefficients KJ, K' etc, then the
intensity after passage through the absorbers is I = I e~^ 1 2+*"'x.
The quantities KJ, K', ... are complicated functions of the frequency
-1
(herein designated as v and measured in cm ) and of the temperature,
pressure, and concentrations of the gases. Hence I is also a complicated
function of the same variables. The long-path spectrophotometer measures
I/I as a function of frequency and, from that data, it is necessary to
determine the concentration of the gases.
To see how the concentration enters the problem, we write a K'x as
Kpx = KW = aN x, where N is the number of absorbing molecules per unit
volume, p is the density of absorbing gas per unit area of the radiation
path, a is the molecular absorption coefficient, W is the area density of
gas, and K is the absorption coefficient. A very wide range of units has
been used for these quantities in the literature. The reader should refer to
Deutschman and Calfee (1967) and Calfee (1971) for transformations between
commonly used sets of units. Table A1.4 lists a summary of these trans-
formations. In this report, W is in units of molecules per square centimeter
and I is in units of centimeters squared per molecule.
If the spectrophotometer measures I/I as a function of wavelength
after suitable adjustment and calibration, the data available (by taking
the natural logarithm of each side of 1:1) is the product N
KW = £n-=- , where KW is a function of frequency. But KW = .Z^K.W.,
I ' M 7 j^j ! x»
o
-------
where N is the total number of absorbing gases, K. is the absorption
coefficient, and V\L is the (area) density of the i gas. If we have the
necessary spectral information about the gases, the K. is known, and we
need to determine the concentrations (determinable from W,) of the gases
present. The ROSE* system allows for an I/I output by having an internal
blackbody that can be set at the same temperature as the remote blackbody.
Their signals fall on the same detector with different chopping frequencies,
They are separated and amplified by two phase sensitive locking amplifiers,
and then ratioed. The problem that occurs (spurious variation of I/I ) is
presumably due to the differences in the optical paths of the two beams.
Thus a true I/I is not provided.
If the instrument only provides relative transmission, I, rather
than I/I , it is necessary to scale the data before the calculation of
concentrations. The scaling must be done so that the input data to the
computer program lies between zero and one, just as does I/I . This can
be accomplished by simple scaling only if the response of the spectrophoto-
meter is the same over all frequency ranges. If not, the data must be
scaled in sections. In either case, only relative concentrations can be
obtained. If absolute concentrations are required, and I/I is not
provided by the spectrophotometer, they may be obtained by adding the
spectrum of a known amount of calibrating gas. There are several ways of
doing this in practice. If, for example, the humidity is known along the
path, this could provide an absolute scale. The accuracy would suffer,
however, for obtaining low-level concentrations in comparison with the
relatively larger concentration of water vapor, and it would be more
accurate to choose an independent measurement of a gas of lower
concentration. Another method would be to insert a test cell containing
a gas of known concentration, adjusted so that its absorption falls in a
spectral range unlikely to interfere with the unknown spectra. The concen-
tration should be adjusted so as to fall within the range of most of the
unknown gases, by trial and error, if necessary. The material of the cell
* Remote Optical Sensing of Emissions
-------
must be nonabsorbent or accounted for by observing the effect of an empty
cell. I can also be obtained by recording a spectrum upwind of an extended
source or in a clean environment. The quantity can then be obtained down-
wind at the same path length and reference blackbody temperature. The
resulting I/I will then be an absolute spectrum.
The problem presented by the infrared continuum can under some
circumstances be obviated by this method. The infrared continuum, a general
increase of all absorption curves over a large part of the spectral range
being considered here, occurs when long paths though the atmosphere are used.
At present, only empirical methods of dealing with it are known. Over the
relatively short path (2-4 km) planned for the spectrophotometer described
in this report, the problem will often arise. (See Shapiro and Gush (1966),
Burch (1970), and McClatchey et al., (1971).) However, in contamination such
as occurs in Los Angeles, the continuum is often a significant contribution
to the spectrum.
Determining concentrations is complicated by the inescapably finite
signal-to-noise ratio (S/N) of the spectrophotometer and the noise due to
atmospheric scintillation. Thus, the signal from the spectrophotometer
has a noise component. It is probably not additive, but has been non-
linearly mixed by the system detector. However, a complete noise analysis
of the detector is not available, and the character of the noise is
dependent on the particular system. Therefore, we will make an assumption,
for convenience, that the noise on the signal is additive Gaussian noise,
that is, the output of the system is composed of a signal and noise, (I/I +N),
where the mean value of N is zero. The bandwidth of the Gaussian noise is
determined by the time constant of the circuitry; it has essentially the
same spectral bounds as the signal I/I . The precision of the measurement
is dependent on the S/N ratio and it, together with the resolution of the
ratio digital voltmeter, determine the error. The ratio digital voltmeter
4
of our system has a resolution of 1 part in 10 . However, the
S/N ratio depends on atmospheric turbulence and is probably never as large
4
as 10 (in the absorption mode). Hence the signal-to-noise ratio is usually
the limiting factor in the precision of concentration measurement.
8
-------
Streiff and Claysmith (1972) discuss the signal and noise
characteristics of the spectrophotometer system. Their data and a
private communication from M. Streiff indicate that practical working
ranges of S/N run from 10 to 100. Without further tests of the system
under working conditions, a complete noise analysis cannot be performed.
To offer guidance on the degradation suffered in precision under finite
signal-to-noise ratios, the computer analysis has been performed with an
added simulated noise (see section 3.1.2).
In a different sense than the noise of atmospheric fluctuations
and the spectrophotometer system discussed above, the presence of water
vapor and to some extent that of natural carbon dioxide in the atmosphere
constitute important sources of noise or interference. Water vapor
concentration varies considerably in time and space, and its very strong
absorption spectrum is found through most of the frequency region of
interest. It interferes least with other spectra in the region from 2400
to 2600 cm" . Unfortunately, only a few gases of interest have observable
absorptions there, and thus water vapor must always be one of the
unknown gases whose concentration is determined by the procedure described
herein. If it is measurable otherwise, it may be desirable to consider
removing its effect from the data by preprocessing. When the
absolute humidity and temperature of the sample are known, the water vapor
spectrum may be computed (1.2) and divided from the absorption spectrum
to eliminate its effect. Carbon dioxide is also a naturally variable
constituent, but its variation and absorption are less than those of water;
and it usually is not an important interferant.
1.2 THE MATHEMATICAL PROBLEM
The transmission spectrum of an atmospheric constituent I/I is a
unique function of frequency v, temperature T, pressure p, and concentration
W, that is, (I/I ). = f-(v, W., P, T), where the subscript indexes the
-------
specific gas. Figure 1.1 (a-c) shows typical spectra of normal constituents
and pollutants, F1.2 shows the normal atmosphere. If the transmission
path contains two different gases, the resulting transmission spectrum is
obtained by application of the Bouguer-Beers law (1:1). According to the
discussion of 1.1, the spectrum due to two gases is obtained by
multiplying the relative transmissions, that is,
-(K W +K W 1
(I/IQ) = (i/i D! (I/I0)2 = e 2 ' ^i:2)
and in general, for N gases, N
(I/IQ) = n (I/I^ = e ^ x A , (1:3)
where (I/I ) is the transmission of the composite atmosphere.
It should be noted that in whatever part of the frequency range one
of the gases has zero (or small) transmission, the composite spectrum
(I/I ) is also zero. In any part of the frequency range where the
transmission is one (no absorption) for all gases but one, the spectrum
will be identical to that one. Otherwise, the spectrum will be the
product of the constituent spectra and will generally present a complex
problem to the experimenter attempting to identify the constituents and
measure their concentrations. Figure 1.3 shows spectra of two gases,
CH. (methane) and HNO_(nitric acid), and their combined spectrum.
4 o
The spectra of the gases (I/IQ)i are functions of v, T, P, and W.
We assume that only the concentrations W. are unknown, and the compound
spectrum recorded as a function of frequency v is the output of the
spectrophotometer system.
In practice, (I/I ). as found in the library is accurate for only
those values of P and T for which the spectrum was determined. Departures
of (I/I ). from true (I/I ), however, are negligible for the range of
temperatures and pressures found on the earth's surface.
10
-------
16
A
O
It
14
4
14
U
HNOj
HjO
C02
HjO
HN03
SOj
SOj
03
NH,
CH4
1400
NH3
C02
CO
CO
NjO
1800 2200
FREQUENCY
HjO
HCL
HCL
CH4
2600
3000
Figure 1.la
INFRARED SPECTRA OF SELECTED ATMOSPHERIC CONSTITUENTS
(Bands are identified on the overlay)
-Ifl
U)
V>
0
N
-------
u '« -"
J
u
** '"
. ill
i "~ " - ' .*; -- "' ...: - ;:
u ;.,' .'».,. u.,. ,5; is .
j . ''*' " > *
u ' *''".;;..;" ' - " t
M - ^ u
'J " *
u ^F ' ^^
1 . ; " - ;
Ji*
'». ---
u
Figure l.lb
INFRARED SPECTRA OF SELECTED ATMOSPHERIC CONSTITUENTS
-------
r 1*
i"
u ., . -'mini" | j
,u
j
Iu
i
COKEITMTIOKWOLECI/UICVl ^^ ,
: W NITMIC OXIDE " ~ ' II"I-IIBI" «"l"«""' - - U
Ui ^ ";' H»mi>e« niLIIDE 114til" IHJMin
J ti«« iT»m« MI inn' lumin. imnii
IM^II m' i.H x n't' IH i tin
U
J
.
mum »' IK i nil. t« ,,»
i arinnari
i dP »tm *mtM< j
U
J
U
I- -. -- ~._ «i niiniiiiniwit IHIII" iwiii'i
_.._ ' MI iiTiacinnmiDE i< n". IH > n'l
u «I,C»,IIC« METOK IN nil CHI IM« ||1I. IHI ||H
I | ll»liH CM I I H » ll«. I 11|H
1 mmm "* " H to FIMU "»m.oi-Iui«ii».i«iii» iU
J . ' "' ' ' J
. ' ~ ' ~ ^ "
u
*->- ' ' *
^ y-'*>^'' - L
>«"r~ Si S" i«i ' nii ai ' m - ~tm
mmm «-'>
Figure 1.1c
INFRARED SPECTRA OF SELECTED ATMOSPHERIC CONSTITUENTS
-------
I*
'llf-l
H-V-T-, - Y-
VI , V 'I i I -Y^
*
asmu,
P : I
o.t |M-i)
Figure 1 .2
INFRARED SPECTRUM OF THE NORMAL ATMOSPHERE
14
-------
COMPOSITE
SPECTRUM
Figure 1.3
A PORTION OF THE TRANSMISSION SPECTRA OF
CH4, HN03> CH4, AND H20 AND THEIR COMPOSITE SPECTRA
15
-------
Therefore, in this application p and T are assumed constants, with
p = 1 atm(1.02 x 10 dynes/cm ) and T = 296°K. For special circumstances,
accurate (I/I ). may be calculated for any temperature and pressure only
for molecules whose line parameters are fully known. This is accomplished
using program DEGRADE (Deutschman and Calfee, 1967) . For spectra taken
experimentally, there is no known transformation to a new set of P and T
without a detailed study of the molecule.
We call the output of the spectrophotometer S(v). We must choose
the concentrations W. for all the gases which are in the sample path in
such a way that the spectrum (I/I ) of (1:3) is the same as, or
as close as possible to, the experimental data S(v). Simply put, this is
performed by choosing a trial set of W., calculating (I/I ), and then
comparing (I/I ) with S(v). The comparison is made by calculating an
error function:
n
(I/I j W^) - S(V)
1:4
Here the symbol (I/I ; W^) signifies the composite transmission spectrum,
O 1 *.-!
calculated from the spectral library under the assumption that the i
»
gaseous constituent has the trial concentration W? . Of course (I/Io;Wn)
is a function of frequency. The number of terms in the sum is n. The
summation over v is somewhat arbitrary; maximum accuracy will be achieved
by using all statistically independent measurements. A good compromise
for the number of points was found by trial and error to be 5 per
resolution element. (In the calculation described in this report, the
data S(v) and the spectra of the individual gases are stored at increments
of Av = 0.2 or 0.6 cm" , depending on the frequency range (see section 2).)
The trial spectra of individual gases are obtained by interpolation
and extrapolation from a set of spectra given in the library for specific
concentrations which span the range of concentrations expected in the
atmospheric sample.
16
-------
The first estimate of the concentrations, W^!, is generally not correct,
of course. After the calculation of the first error, E^, the concentrations
are varied in an orderly stepwise manner, calculating successive
corresponding errors E., until the set of concentrations is obtained
that minimize E.. These concentrations are then the "best" estimate of
the concentrations in the atmospheric path in which the sample data
S(v) were obtained.
The calculation is performed by the computer subroutine MINMYZD,
developed by Slutz and Winkelman (1964) and applied by Lawrence and
Hallenbeck (1965). MINMYZD is described in detail in A2.
This computer routine, although in principle as simple as described
above, requires careful design to prevent costly, long computer runs and
inaccuracies. The problem is to find the minimum value of a
multidimensional function, E, whose variables are the W.. Care must be
taken to avoid local minima in this function and to find the absolute
minimum. The subroutine MINMYZD has features which practically assure that
result in applications such as herein discussed.
It is not possible and not desirable to discuss here the problems
of uniqueness involved in minimization problems such as this one.
Instead, the long experience with this program and its successful
application to many physical problems gives assurance of its usefulness
and accuracy.
Arbitrarily chosen initial concentrations prolong the computer
running time and may cause inaccuracy in the final results. Hence initial
concentrations are selected for each constituent by computing a "best fit"
in a very limited spectral range where the constituent has a prominent
spectrum and the interference from other gases is small. These selected
concentrations are then employed as initial, approximate choices in the
computation over the whole spectral range. When the number of gases is
small, these initial values may be sufficient for the experimenter's needs.
It is true that the optimization or best-fit computation here
employed is not the only way the concentrations could be obtained from the
experimental data. Many other methods exist, under the names of best-fit,
17
-------
minimization, or inversion methods. No attempt was made to evaluate
which of these would be the best method. The MINIMYZD routine was selected
because of its long history of successful efficient operation. Methods such
as that of Backus and Gilbert (Westwater and Cohen, 1973) offer some
advantages of error estimation. Matthes (1973) has devised a method to
circumvent problems arising in least-squares fitting which results in
negative concentration estimates. This problem may be prevented in the
present program by placing a lower bound of zero on any concentration
estimate, but it is probably better practice to eliminate such absolute
constraints on the MINMYZD subroutine to prevent nonlinear effects on
the search process.
18
-------
2. PRESENT STATUS OF THE LIBRARY OF ABSORPTION SPECTRA
It is essential in the identification of pollutant gases by
long-path absorption spectroscopy to possess accurate spectra of the
gases expected to occur in the atmosphere under examination. These
spectra must be complete in the spectral range covered by the long-path
spectrophotometer because "maverick" gases (i.e., gases present in the
sample but not in the spectral library) cause inaccuracy in the
calculation of the concentrations. Some examples of such inaccuracies
are presented in 3.
The spectra must be known in such a way that they may be determined
as a function of temperature, pressure, and concentrations to suit the
conditions in the test atmosphere. If, for example, the line position,
the line strength, and the broadening parameters of a gas as functions of
temperature and pressure are known, the spectra may be computed for any
concentration by computer routines (Deutschman and Calfee, 1967, and
McClatchey et al., 1973). If these parameters are unknown, but accurately
measured experimental spectra are available for concentrations near those
of interest, we may extrapolate from them to other concentrations.
To make use of such experimental data, it is necessary that
the spectra be obtained at a temperature and pressure near those normally
encountered over polluted cities, that resolution is comparable to or
better than that of the spectrophotometer used, and that the ordinate of
the spectra is calibrated as an absolute transmittance (or, that the line
strengths in some way are known).
Whether such spectra are available depends on whether researchers
have been interested in these specific molecules. No concerted effort
has been made to the authors' knowledge to obtain detailed spectra of
all atmospheric constituents including pollutants. However, an effort
has been made by a loosely knit group of researchers* to obtain detailed
* The Group on Atmospheric Transmission Studies (GOATS)
19
-------
spectra on the normal atmospheric constituents. These data are summarized
by McClatchey et al., (1973).
In general, however, the search of literature on the spectra of
normal and pollutant atmospheric constituents reveals a wasteland
punctuated by a few oases of excellent work. Spectroscopists, academic
and industrial, have not been supported in obtaining the high-resolution
(1 cm" ) infrared spectra with a well calibrated intensity scale at
atmospheric pressure needed for concentration measurements in an
atmosphere containing many interfering gases.
The chief limitation of the application of this (or any)
spectroscopic technique to the detection of atmospheric constituents lies
in the ignorance of the spectra of many atmospheric constituents,
particularly pollutants. Just as the criminal's fingerprints cannot be
used to identify him unless they are in police files, so observation of
a gas spectrum only confuses as if it is not in our library. Notably,
spectra are unavailable for many hydrocarbons. Also spectra of excited
(high temperature) molecules such as found in smoke stacks are not in the
literature.
The sources of the spectral library are summarized in Al and
detailed descriptions are given there of the parameters of the
spectra. Also described are manipulations performed in converting spectra
found in the literature into computer compatible form. We summarize
in T2.1 the gases available in the spectra library. Spectra for
selected concentrations for all these gases are displayed in Fl.l(a-c).
20
-------
Table 2.1
GASES AVAILABLE IN THE SPECTRAL LIBRARY
Gas
Water
Carbon dioxide
Ozone
Nitrous oxide
Methane
Carbon monoxide
Hydrochloric acid
Sulfur dioxide
Nitric acid
Ammonia
Nitric oxide
Hydrogen sulfide
Ethylene
Nitrogen dioxide
Acetone
Formaldehyde
Ethane
Normal -butane
I so -butane
Propane
Pentane
Formula
H20
co2
°3
N20
CH4
CO
HCL
so2
HN03
NH3
NO
H2S
C2H4
N02
CH3CH2HCO
H2CO
C2H6
n-C.H.,.,
4 10
i-C H
4 10
C3H8
CCH10
5 12
Accuracy of Spectrum *
1
1
1
1
1
1
1
2
2
2
1
2
2
2
3
2
2
3
3
3
3
*Accuracy of spectrum: (1) Line information highly trustworthy;
(2) Good experimental spectrum available, trust-
worthy; and
(3) Low resolution only, some uncertainty.
21
-------
3. BEST-FIT CALCULATIONS OF CONCENTRATIONS
FROM A SIMULATED COMPLEX SYSTEM
The application of the MINMYZD fitting routine to specific problems
requires a controlling program to direct progress from the input of data
through the minimization process to the final gas constituent concentrations
and error estimates. This program, called EPAGAS, is described in detail
in A2.
Using this program, the entire computation was tested by preparing
trial spectra, simulating the data obtained from a long atmospheric path
spectrophotometer experiment, and then obtaining the concentrations. The
tests were made "blind", that is, the programmer was presented with the
simulated sample data without knowledge of the concentrations until the
computer results were compared with the known composition. In this
section, a number of such computations are reviewed. The objectives were
to "debug" the program, to test the success of determining the concentrations
of gases from their composite spectrum, to determine the accuracy achievable,
and to estimate the degradation of accuracy due to finite signal-to-noise
ratio in the data. Further, the effects of a lack of a one-to-one
correspondence of the library of spectra with the atmospheric constituents
are examined.
In these preliminary tests, library spectra of H20, N20, CH4, S02,
CO, C0_, and 0_ were employed in various concentrations of combinations to
produce a complex spectrum.
3.1 ACCURACY OF CONCENTRATION CALCULATIONS
Tables 3.1 through 3.4 list the trial runs made on the
combinations of gases indicated at the top of the table. Each fitting
was made over a narrow frequency range to determine the suitability
of that range for detection of the gas under investigation and to observe
inaccuracies arising when the gas has no absorption in the region.
It should be noted here, and it is emphasized and clarified in
A2, that the analysis of complex spectra by the method which is
the subject of this report should not be done blindly or automatically.
22
-------
t»
Table 3.1
CALCULATION OF GAS CONCENTRATIONS
FROM SIMULATED COMPOSITE SPECTRA (SET 1)
(INFINITE SIGNAL-TO-NOISE RATIO)
Code
8888 a
b
c
d
e
f
jr
h
i
Kaver.unher
Range
(cm-1)
1250.
1250,
12--C.
1290.
1310.
1330.
1350.
1370.
1380.
0-1249.8
0-1269.8
0-128P.S
0-1309. g
0-1329.8
0-1349.8
0-1369.8
0-1389. S
2-1400.0
3
3
_;
3
3
3
3
3
3
H20 (1
Kg = ?.73
We
-3x10"
.73xlf>22
22
.73xl022
-»T
.73x10"
.-3xlC22
.73xl022
.C8xl022
.6'xlO22
x 10"
AK
.00
.00
.00
.00
.00
.CO
.00
.05
.06
Kg = l.PS x 101 8
We tv:
1.09X1P18
l.OPxlP18
l.OSxlP18
1.10xl01S
1.23X1018
1.21X101*
1.21xl018
1.2U1018
1.21X1P18
-.01
-,P1
.00
-.02
-.15
-.13
-.13
-.13
-.13
a: (S
Kg = l.PS x
Vc
1 Q
1.08x10
1 Q
1.08x10
TO
1.08x10"
l.OSxlO19
1.06X1019
1.08X1019
1.06xlP19
1.4SxlP19
10
1.4SX1T1-
AioJ»
tx
.OP
.OP
.00
.00
.02
.00
.02
-.40
-.40
SO 2 (20.2)
h'g = 3.?4 x IP17
We iW
ggxlo17
17.88X101'
17.S6XK1'
P.75X1P17
4.29X101''
3.51X1017
3.32xl017
4.SSX1P1'
S.lOxIP1'
-14.34
-14.34
-14.34
2.79
- 0-75
0.03
0.22
- 1.34
- 1.56
HN03 (30.2)
hg = 3.54 x 1017
he Ah
0.0'xlO1'
O.lOxlP1'
S.SSxlO1
3.45X1017
3.42X1017
3.56X101'
4.32X101'
lfi.3Sx!0]
18.23.xl01'
3.47
3.44
- 0.01
0.09
0.12
- 0.02
- 0.78
-14.84
-14.69
Final Transmission
Error
RMS
K td. Unwtd.
.03
.06
.03
.05
.03
.02
.03
.03
.03
6.60
6.60
6.41
1.25
0.35
0.06
0.37
6.68
6.61
Wg » given concentration in simulated sample, molecules cm"
_ 7
Kc « calculated concentration, molecules cm
AW » Wg-Wc (the factor of a power of ten
is omitted), molecules cm"2
See 3.1 for definition of the weighted and
unweighted final transmission errors (last
two columns). Numbers following gas
formulae are computer codes.
A composite sample spectrum was obtained by
adding the absorption coefficients as a
function of frequency and calculating the
transmission (See 1.2). The sample spectrum
was then analyzed by MINMYZD to obtain the
calculated concentrations. The error is the
difference between the calculated and given
transmission. (See 3.1),
-------
to
Table 3.2
CALCULATION OF GAS CONCENTRATIONS
FROM SIMULATED COMPOSITE SPECTRA (SET 2)
(INFINITE SIGNAL-TO-NOISE RATIO)
Code
9999 a
b
c
d
Kavenumber
Range
Cera-1)
2100.0-2159.4
2160.0-2219.4
2220.0-2279.4
2240.4-2299.8
H20 Cl
Wg = 3.73
Kc
3.74xl021
3,72xl021
3.68xl021
3.69xl021
x 1021
AW
-.01
.01
.05
.04
COj (2B1)
Wg » 1.08 x 1021
We
l.OSxlO21
l.OSxlO21
l.OSxlO21
l.OSxlO21
AW
.00
.03
.00
.00
0, (3B1)
Wg « 2.69
We
2.84xl017
5.94xl017
S.94xl017
5.94xl01T
x 1017
AK
-0.15
-3.25
-3.25
-3.25
N20 (4B1)
Kg * 1.08 x 101 e
Kc
1.22X1018
l.OSxlO18
l.OSxlO18
l.OSxlO18
AK
-.14
.00
.00
.00
CO (6B1)
Wg = 1.88
Kc
l.SSxlO19
l.SSxlO19
l.SSxlO19
2.98xl019
x 10"
AK
0.00
0.00
0.00
-1.10
Final Transmission
Error
RMS
Ktd.
.03
.03
.03
.03
Unwtd .
0.09
1.45
1.45
1.53
Wg given concentration in simulated sample, molecules cm"
We calculated concentration, molecules cm"
AW m Kg-Wc (the factor of a power of ten
is omitted), molecules cm"2
See 3.1 for definition of the weighted and
unweighted final transmission errors (last
two columns). Numbers following gas
formulae are computer codes.
A composite sample spectrum was obtained by
adding the absorption coefficients as a
function of frequency and calculating the
transmission (See 1.2). The sample spectrum
was then analyzed by M1NMYZD to obtain the
calculated concentrations. The error is the
difference between the calculated and given
transmission. (See 3.1),
-------
cn
Table 3.3a
CALCULATION OF GAS CONCENTRATIONS
FROM SIMULATED COMPOSITE SPECTRA (SET 1)
Code
S/N = 10
8888 a
b
c
d
e
f
g
h
i
Kavenumber
Range
(cm-1)
1230.0-1249.8
1250.0-1269.8
1270.0-1289.8
1290.0-1309.8
1310.0-1329.8
1330.0-1349.8
1350.0-1369.8
1370.0-1389.8
1380.2-1400.0
H2n (i
Kg = 3.73
We
3.80xl022
3.74xl022
3.S3xl022
3.98X1022
3.98xl022
3.T2xl022
3.7Sxl022
3.71xl022
3.65X1022
x 1022
tx
-.07
-.01
.20
-.25
-.25
.01
-.02
.02
.08
K20 (4A1)
Kg = 1.08 x 101 *
We AW
O.SOxlO18
0.82xl018
1.04xl018
1.03X1018
0.75X1018
i.aixio18
1.21xl018
1.21X1018
1.21X1018
.28
.26
.04
.05
.33
-.13
-.13
-.13
-.13
Q! * (5A1)
Wg = 1.08 X IP19
Kc AN
1.04xl019
1.13xl019
1.06xl019
l.OSxlO18
1.02xl019
l.llxlO19
0.99xl019
1.48xl019
1.48xl019
.04
-.05
.02
.00
.06
-.03
.09
-.40
-.40
S02 CO. 2)
Wg * 3.54 x 10ir
Kc AK
17.88X1017
17.88xl017
17.88xl017
0.28xl017
2.50X101'
2.75X1017
4.15xl017
4.24X101'
5.60X101'
-14,34
-14.34
-14.34
3.26
1.05
0.79
- 0.61
- 0.70
- 2.06
HN03 (30.2)
Wg = 3.54 x 1017
We AW
- 1.89xl017
0.68X1017
4.81xl017
3.54xl017
3.08xl017
3.9SX1017
4.69X1017
18.17xl01;
18.22X101'
5.43
2.86
-1.27
0.00
0.46
- 0.41
- 1.15
-14.63
-14.68
Final Transmission
Error
RMS
Ktd . Unwtd .
2.08
2.28
1.86
1.50
0.76
0.80
0.28
0.14
0.14
6.86
6.54
6.44
1.46
0.54
0.40
0.59
6.55
6.63
Wg » given concentration in simulated sample, molecules
We « calculated concentration, molecules cm"
AW Wg-Wc (the factor of a power of ten
is omitted), molecules cm"2
See 3.1 for definition of the weighted and
unweighted final transmission errors (last
two columns). Numbers following gas
formulae are computer codes.
-2
A composite sample spectrum was obtained by
adding the absorption coefficients as a
function of frequency and calculating the
transmission (See 1.2). The sample spectrum
was then analyzed by MINMYZD to obtain the
calculated concentrations. The error is the
difference between the calculated and given
transmission, (See 3,1,2).
-------
Code
S/N =30
8888 a
b
c
d
e
f
g
h
i
Kavenumber
Range
(en-1)
1230.0-1249.8
12SO. 0-1269. 8
1270.0-1289.8
1290.0-1309.8
1310.0-1329.8
1330.0-1349.8
' 1350.0-1369.8
1370.0-1389.8
1380.2-1400.0
HiO (1A1)
Kg = 3.73 x 10
We
3.76X1022
3.66X1022
3.68X1022
3. 73x1 022
3.7SX1022
3.73X1022
3. 76x1 O22
3.68xl022
3.69x!0:2
22
Oi
03
07
OS
00
02
00
03
05
04
N20 (4A1)
Kg - l.OR x 10
Kc
O.SOxlO18
1.23xl01R
l.OSxlO18
l.llxlO18
l.SSxin18
1.21X1018
1.21xl018
1.21xl018
1.21X1018
i t
fiK
58
15
03
03
30
13
13
13
13
Table
(S/JJ
CIU
Kg - 1.08
Kc
l.OSxlO19
10
1. 07x10*-
l.osxin19
1C
i.07xin
0.97xl01?
1Q
i.nsxio
1.07xl019
IP
1.48x10
1.48X1019
3,3t>
30)
(5A1)
x 10"
AK
.00
.01
.03
.01
.11
.03
.01
-.40
-.40
SO; (20.
Kg ' 3.54 x
Kc
17.88X101'
17.88xl01?
17.88xl017
l.OOxlO1'
4.35xl017
3. 64x1 O17
3. 06x1 O17
4.83xl017
4,54xl017
2)
10"
Ui
-14.34
-14.34
-14.34
2.P4
- 0.81
- 0.10
0.48
- 1.29
- 1.00
HNO,
Kg - 3.54
Kc
O.lSxlO17
-0.69X1017
4.11xl017
3.50X1017
3.40xl017
3.54X1017
4. 58x1 O17
18.20xl017
18.24X1017
(30.2)
x 10lr
AW
3.39
4.23
- 0-57
0.04
0.14
0.00
- 1.04
-14.66
-14.70
Final
RMS
Ktd.
.65
.75
.65
.46
.25
.26
.10
.05
.06
Transmission
Error
Unwtd .
6.59
6.69
6.42
1.14
0.39
0.07
0.52
6.58
6.59
Kg given concentration in simulated sample, molecules
Kc « calculated concentration , molecules cm"2
flK Kg-Kc (the factor of a power of ten
is omitted), molecules cm"^
See 3.1 for definition of the weighted and
unweighted final transmission errors (last
two columns). Numbers following gas
formulae are computer codes.
cnr
A composite sample spectrum was obtained by
adding the absorption coefficients as a
function of frequency and calculating the
transmission (See 1.2). The sample spectrum
was then analyzed by MINMYZD to obtain the
calculated concentrations. The error is the
difference between the calculated and given
transmission. (See 3.1.2).
-------
Table 3.3c
(S/N-100)
Code
S/N = 100
8888 a
b
c
d
e
f
g
h
i
h'avenumber
Range
(cm-1)
1230.0-1249.8
1250.0-1269.8
1270.0-1289.8
1290.0-1309.8
1310.0-1329.8
1330.0-1349.8
1350.0-1369.8
'1370.0-1389.8
1380.2-1400.0
HjO (1A1)
Wg = 3.73
Kc
3.74X1022
3.71X1022
3.73xl022
3.72xl022
3.72xl022
3.74X1022
3.74xl022
3. 68x1 O22
3.67xl022
x 1022
AW
-.01
.02
.00
.01
.01
-.01
-.01
.05
.06
NjO (4A.1)
Wg « 1.08
h'c
l.llxlO18
1.13xl018
1.09xl018
i.oaxio18
1.21X1018
1.21X1018
1.21X1018
l.ZlxlO18
1.21X1018
x 10"
ex
-.03
-.05
-.01
-.01
-.13
-.13
-.13
-.13
-.13
CHW (5A1)
Wg = 1.08
Kc
1.07X1019
1.07X1019
1.09X1019
10
1.08x10
l.OSxlO19
1.07X1019
1.06X1019
1.48X1019
1.48X1019
x 10"
AW
.01
.01
-.01
.00
.00
.01
.02
-.40
-.40
SO 2 (20.2)
Wg 3.54 x 1017
We
17.88xl017
17.88xl017
17.88xl017
0.74xl017
4.11xl017
3.57xl017
3.25X1017
4.84X1017
S.19xl017
AW
-14.34
-14.34
-14.34
2.80
- 0.57
- 0.03
0.29
- 1.30
- 1.65
HNOj (30.2)
Wg = 3.54 x 1017
Kc
- O.OSxlO1'
O.OSxlO1'
3. 49x1 O17
3.SOX101'
3.SOX1017
3-SlxlO17
4.37X1017
18.37X1017
18.23xl017
AW
3.57
3.46
0.05
0.04
0.04
0.03
- 0.83
-14.83
-14.69
Final Transmission
Error
RMS
Ktd.
.25
.22
.20
.17
.07
.09
.04
.03
.03
Unwtd.
6.61
6.60
6.41
1.25
0.26
0-06
0.40
6.66
6.61
Wg » given concentration in simulated sample, molecules cm
We » calculated concentration , molecules cm"
AW » Wg-Kc (the factor of a power of ten
is omitted) , molecules cm-2
See 3.1 for definition of the weighted and
unweighted final transmission errors (last
two columns). Numbers following gas
formulae are computer codes.
A composite sample spectrum was obtained by
adding the absorption coefficients as a
function of frequency and calculating the
transmission (See 1.2). The sample spectrum
was then analyzed by MINMYZD to obtain the
calculated concentrations. The error is the
difference between the calculated and given
transmission. (See 3.1.2).
-------
Table 3.4
CALCULATION OF GAS CONCENTRATIONS
FROM SIMULATED COMPOSITE SPECTRA (SET 2)
(S/N 10, 30, 100)
00
Code
S/N 10
9999 a
b
c
A
S/N » 30
9999 a
b
c
d
S/N « 100
9999 a
b
c
d
Kavenunber
Range
Ccn.-1)
2100.0-2159.4
2160.0-2219.4
2220.0-2279.4
2240.4-2299.8
2100.0-2159.4
2160.0-2219.4
2220.0-2279.4
2240.4-2299.8
2100.0-2159.4
2160.0-2219.4
2220.0-2279.4
2240.4-2299.8
HjO (1B1)
Wg » 3.73 x in21
he
3. 66x1 O21
4.n3xl021
S.lOxlO21
-0.93xl021
3.73X1021
3.89xl021
4.72xl021
2. 67x1 O21
3.72xl021
3.68xl021
3. 68x1 O21
3. 78x1 O21
AK
0.07
-0.35
-1.37
4.66
0.00
-n.ie
-0.99
1.06
0.01
0.05
o.ps
-0.05
C02 (2B1)
Kg 1.08 x 10J1
Kc
1.07xl021
0.35X1021
l.llxlO21
l.lOxlP21
1.14xl021
1.27X1021
1.09xl021
1.08xl021
1.07xl021
l.OSxlO21
l.OSxlO21
l.OSxlO21
AK
.01
.73
-.03
-.02
-.06
-.19
-.01
.00
.01
.00
.00
.00
0, C3B1J
Wg - 2.69 x 1017
Kc
1.79X1017
5.94xl017
5.94xl017
5. 94x1 O17
2.S3xl017
S.94X1017
5.94xl017
5.94X101'
S.OlxlO17
S.94xl017
5.94xl017
5.94X1017
AW
0.90
-3.25
-3.25
-3.25
0.16
-3.25
-3.25
-3.25
- .32
-3.25
-3.25
-3.25
N20 (4B1)
Kg » 1.08 x 10l 8
Kc
1.34X1018
1.07xl018
l.llxlO18
1.06X1018
1.19xl018
l.OSxlO18
l.OSxlO18
1.07xl018
1.19xl018
l.OSxlO18
l.OSxlO18
l.OSxlO18
AK
-.26
.01
-.03
.02
-.11
.00
-.01
.01
-.11
.00
.00
.00
CO (6B1)
Kg = 1.88 x 10"
We
1.91X1019
l.SSxlO19
1.87xl019
1.33xl019
l.SSxlO19
l.BSxlO19
1.89X1019
1.49xl019
l.SSxlO19
l.SSxlO19
1.87xl019
2.59xl019
AK
-.03
.00
.01
.55
.00
.00
-.01
.39
.00
.00
.01
-.71
Final Transmission
Error
RMS
Ktd.
1.09
1.28
2.64
3.54
0.36
0.44
0.86
1.19
0.10
0.14
0.26
0.31
Unwtd.
0.42
1.50
1.58
2.55
0.09
1.46
1.52
1.54
0.15
1.4S
1.45
1.49
Wg - given concentration in simulated sample, molecules
We calculated concentration, molecules cm"2
AW Kg-Kc (the factor of a power of ten
is omitted), molecules cm'2
See 3.1 for definition of the weighted and
unweighted final transmission errors (last
two columns). Numbers following gas
formulae are computer codes.
-2
A composite sample spectrum was obtained by
adding the absorption coefficients as a
function of frequency and calculating the
transmission (See 1.2). The sample spectrum
was then analyzed by MINMYZD to obtain the
calculated concentrations. The error is the
difference between the calculated and given
transmission. (See 3.1.2).
-------
To apply directly the method to a sample of unknown composition, utilizing
the full library of test spectra, would be an example of how much must be
paid (in computer costs) for total ignorance of the sample. Maximum
information available on the composition of the unknown sample should be
utilized in selecting the library of test spectra. In the examples
discussed below (3.1.3), we will see the penalty on the accuracy of the
calculation for omitting a spectrum from the library when the gas is in
the sample.
Further, we note that suitable choices of frequency regions to
investigate the presence of a given gas must be made. Obviously, if at
concentration levels of interest, a gas shows no (or weak) absorption,
its concentration may not be measured there. In F3.1 (and Fl.la)
consider the two spectral regions of absorption by HNO_, namely
-1
those at 880 and 1300 cm . The higher frequency band has strong potential
interference from H20, S02, N20, CH4, and all HC. Thus, if HN03 is
suspected in a sample, the concentration should be calculated in the
lower frequency band. Therefore, the analyst should consult the
positioning of spectra in F3.1 before beginning his analysis.
Suppose an analyst had unlimited free computer time. Should he
not then simplify the exercise of judgment and insert his entire library
into the computer and analyze the whole frequency range? But here he will
find a requirement for careful judgment also. For, as we will see below,
in regions where absorption is effectively zero, the results of the
calculation are small but essentially random. Therefore, the analyst must
choose from among the calculated concentrations obtained in all regions
only those where the absorption is substantial for the gas of interest.
Such a weighting of values could b~ added to the computer program. It
is uncertain whether the advantages outweigh the added complexity in
the program.
In T3.1 through T3.3, a uniform notation is used. Two separate
sets of gases are used. Set 1 consists of H.O, N_0, CH., SCL, and
HN03, and has the Code 8888. Set 2 consists of H20, CO 0 , N_0, and CO,
29
-------
OJ
o
VI/I/IfI C;Hi ////////I
n r
\iiiiiiiJiiiiiiHiiiiitiifiiiiiiiiiiiiiiiiiiiiiiiiiiiiiii mjtiiiiinninnniiT7Tiiniiiiinii[fHiHiiiii\
\.,,,.}.,.,,. i., ,,,;.--,,,.,.}..., ^
v nil ililfll /////// in in EH- ///////////////////////*l
I///////////////////////// Cl«
i///////////////////////////// He i ////iiiinniiiiii/iiiii/k
(
r//////////////7 « ll/llllllliillih
\ *!
f////////////////
/////////////////I I//////////1
//////////I
777//////////////
i/r/in/itifiinirrnn
i//////////////////////////////////////////////////////'
^///////////////////////////CjC1
FIEWilCT SCII 1)
Figure 3.1
PRINCIPAL REGIONS OF ABSORPTION OF
SELECTED NORMAL AND POLLUTANT ATMOSPHERIC CONSTITUENTS
(The region from 1460 to 1800 cnr1 has
strong water vapor absorption)
! T
tl ' i ' '
til loU 1200 14
_^_
It UN
J
llll
I
2lL
22
1
N 24
I U
1 i ' 1 1 '
N 2600 2800 30W 3^00 J«
-------
and has the Code 9999. H20 is included in both sets because it is always
present in the atmosphere; C0_, though always present, is not as variable.
The other gases were chosen to be representative spectra in the upper and
lower frequency ranges (see Fl (a-c) and F3.1). The concentration (molecules
-2
cm ) is denoted for each gas by Wg. The concentrations are chosen to
match representative values in the Los Angeles area. The frequency range
for which the calculation was made is designated by a small letter. The
total range was divided into subranges to observe the effect of spectral
strength and structure on the accuracy of the computation. The calculated
concentration is denoted by We. The difference (Wg-Wc)xlO~^ is a measure
of the accuracy (y is the exponent in Wg). Two measures of rms error are
given, the weighted error and the unweighted error. The weighted error is
of use to the programmer. Here we will refer to the unweighted error.
The penultimate columns of the tables contain those final transmission
errors defined as
1 N (f T I2
N £ ) I(True Transmission),-(Calculated Transmission)k CFj>
where CF. is a "confidence factor" introduced to weight in favor of high
concentrations because they produce, generally speaking, less error.
(See A2.5.) The transmissions involved are those due to the total
sample of all gases. The summation is over each experimental data value
of the frequency range given in the tables. It thus is a figure of merit
for "goodness of fit" over that range. The unweighted error is obtained by
letting all CF =1.
Jv
The unweighted error presented in the final column is a different
sort of measure. It is simply the rms of the fractional errors (Wg-Wc)/Wg,
on a line of the table. In the error column, factors of 10 have been omitted.
-------
Tables 3.1 and 3.2 are based on an assumption of infinite
signal-to-noise ratio. Tables 3.3 (a-c) and 3.4 illustrate the effect
of finite S/N. Tables 3.5 (a, b) and 3.6 illustrate the effect of
"maverick" gases. All these cases are discussed in the following sections.
3.1.1 Calculation of Concentrations^or Spectra Without Noise
Consider the calculations coded as 8888 (a-i) in T3.1 (Set 1).
A simulated sample was prepared, containing concentrations of the gases
listed, unknown to the analyst. Using the simulated sample as data,
concentrations were computed by the MINMYZD method. In the frequency
ranges given in the second column, a concentration of water of
22
Wg = 3.73 x 10 is relatively strong and the small errors CAW) under
HO are to be expected when, as here, the S/N is assumed infinite.
The inaccuracy lies between 0 and 2 percent, the larger error occurring
toward the strong "black-out" absorption band of H90 lying between 1430
-1
and 1800 cm . Interference from the other constituents of the sample
has little effect on the calculated concentrations of H-0.
The calculated concentration of N_0 from the data is done with an
error not exceeding 14 percent. The larger error occurs where the absorption
is very small, and such regions should be lightly weighted in a final
averaging. By such weighting, we may estimate a maximum error of 5 percent
in this case.
Methane (CH.) has a strong absorption spectrum and very distinctive
spectral structure in part of the frequency range covered in T3.1, and
the error in the calculated concentration is small except where the
absorption is very weak in the range 1370 to 1400 cm~ , being 37 percent.
17 -2
The given concentration of S02, Wg = 3.54 x 10 molecules cm ,
is relatively small. Better results than shown can be obtained by using
the stronger absorption spectrum between 1050 and 1200 cm~ (see Fl.la).
Here, the errors, excluding the regions from 1230 to 1310 cm , range from
1 to 44 percent. By using the results only from the strongest absorbing
region, the errors are quite small, averaging 9 percent.
The last gas considered in this set is HNOj. Its spectrum is strong,
32
-------
but the concentration is relatively low. Again we see small errors where
the absorption is strong, and large errors amounting to 100 percent where
the absorption is weak.
The calculations listed under Codes 9999 (a-d) in T3.2 give
a similar set of calculations for a different set of gases, namely H.O,
CCL, 0-, N20, and CO. The spectra of these gases are shown in Fl.la,
except for 0_. Its absorption, in the frequency range 2100 to
-1
2300 cm , is too weak to show on the scale of the graph; it is listed
however in the library of spectra.
In this set of gases, the errors, except for 0 , are quite small.
nd ti.jy become considerable in some cases at
S/N = 10. Generally, the calculated concentrations are within an order of
magnitude of the true sample concentrations (except in the case of SO-
and HNO, in those regions of the spectrum where absorption is weak).
3 '
Whether these errors are tolerable or not depends on the needs of the
experimenter. In many cases, order of magnitude results are all that are
needed for preliminary assessments of a polluted atmosphere. Greater
33
-------
accuracy could be achieved by longer paths and longer integration times,
provided that the instrument is stable against long period drift and the
sample is stationary.
Even at S/N = 10, the errors may be kept low by choosing regions where
the spectrum is strong. It is certainly true, however, that as
S/N becomes much less than 10, the accuracy will decrease significantly.
Table 3.4 may be compared with T3.2 for the effect of noise on
concentrations calculated in the spectral range from 2100-2300 cm" . The gases
CSet 2) H20, C02, 03, N20, and CO are present in the simulated sample
(Codes 9999 (a-d)). Here conclusions parallel those of the above
discussion; a signal-to-noise ratio more than 10 is required to keep
errors below 20 percent, and S/N = 100 produces negligible errors.
The reader should not be perturbed because a few low S/N cases
appear more accurate than those with higher or infinite S/N. This occurs
for two reasons: (1) the noise added may accidently provide an occasional
better fit; and (2) the repetition of the algorithm seeking the best fit is
at the discretion of the operator. In the case of low S/N, the search was
allowed to proceed for more iterations than in the less noisy cases in order
to provide the best opportunity of overcoming the noise.
3.1.3 Effects of "Maverick" Gases in the Sample
In analysis of a spectrum from a specific geographical location, it
is important to have in the library of spectra, utilized in the computer
analysis, a spectrum of each of the gases expected in the area. However, the
library should be kept as small as possible to minimize computer storage
problems. Thus, the library data entered into the computer memory are
chosen to match as closely as possible the set of gases that are in the
atmosphere under measurement. For a larger than necessary library (i.e.,
some gases in the library are not in the test atmosphere), the only
important penalties are increased computer memory required and the
additional costs due to the necessity of breaking the frequency spectrum
into separately treated sections, each run separately on the computer, if
the limited computer memory requires it.
34
-------
To test the effect of having gases present in the sample, but
not present in the computer library (a "maverick" gas), runs were
made on both previously used sets of gases. From each set of calculations,
spanning for Set 1, from 1230 to 1400 cnf in nine steps, each gas,
except H90, was successively deleted from the library (T3.5a,b).
^ _i
Similar procedures were followed in the 2100 to 2300 cm range
(T3.6).
The gas deleted from the library was included in the simulated
data, but was not among those which the computer attempted to fit to the
sample data by best choice of concentration. It thus acted as a
localized noise or perturbation of the spectrum.
Comparing T3.1 with T3.5 (a,b) and T3.2 with T3.6, we see that
appreciable errors may arise by having an insufficient library. The
greatest of these appear where the spectrum of the maverick gas overlaps
a weak spectrum of a gas whose spectrum is in the library. This result
indicates that as a practical measure it is better to add a gas to the
library if its presence in the sample is uncertain.
An examination of the spectra of Sets 1 and 2 in Fl.1 (a-c) and
F3.1 shows that the spectral structure must act as an appreciable noise
background where the gas is a maverick. The MINMYZD method, by utilizing
goodness of fit over a wide spectral region, partially eliminates this
noise effect.
When a gas spectrum is present in the library, but near or at zero
concentration in the sample, the calculation occasionally produces negative
concentrations, small in absolute value. They could be prevented by a
simple change in the program, but it was judged preferable not to affect
the convergence process by puttiug a barrier at zero concentrations. Such
negative concentrations should be considered equal to zero. This situation
is discussed in A2.5.
35
-------
Table 3.5a
GAS CONCENTRATIONS FROM SIMULATED COMPOSITE SPECTRA (SET 1)
Final Transmission
Code
8888 a,B
b,B
c,B
d,B
e,B
f,B
gjB
h,B
i,B
1811 a,Bl
b.ll
c.ll
d.Bl
«,B1
f.Bl
I.B1
h.Bl
i.Bl
h'avenumber
Range
(cm-1)
1230.0-1249.8
1250.0-1269.8
1270.0-1289.8
1290.0-1309.8
1310.0-1329.8
1330.0-1349.8
1350.0-1369.8
1370.0-1389.8
1380.2-1400.0
1230.1-1249.8
1250.0-1269.8
1270.0-1289.1
1290.0-1309.8
1310.0-1329.8
1330.0-1349.8
1350.0-1369.1
1370.0-1389.8
1380.2-1400.0
HzO (1
Wg - 3.73
We
3.78xl022
4-OOxlO22
3.92xl022
4.04X1022
3.69xl022
3.73xl022
3. 73x1 O22
3.68xl022
3.67X1022
4.««tcll22
3.3t*1022
2.90xl022
2.0«xl022
4.S6K1022
4.09K1022
4.33X1022
3. 6»xl O22
3.67xl022
x 10"
AW
-0.05
-0.27
-0.19
-0.31
0.04
0.00
0.00
0.05
0.06
-0-»3
O.J4
0.83
1.67
-0.83
-0.36
-0.60
0.05
0.06
NjO
Wg - 1.
We
(4A1)
08 x 10l '
AW
N.O deleted from
library
8888a,B
8888i,B
S.JlxlO18
4.I5X1011
-O.SSxlO11
5.52xl018
3.71xl011
1.21xl018
1.21X1018
1.21xl018
1.21X1018
for runs
through
-7.43
-3.77
1.67
-4.44
-2.63
-0.13
-0.13
-0.13
-0.13
CH fc (5A1)
Wg 1.08 x
We
l.ioxio19
1.18X1019
i.r>2xio19
1.34xlP19
l.SlxlO19
1°
1.08x10 '
1.06X1019
1.48xl019
1.48X1019
CM. deleted
4
library for
10"
AW
-.02
-.10 .
..06
-.26
-.23
,00
.02
-.40
-.40
from
rum
888la,Bl through
88881, Bl
SUj (20,
Wg * 3.54 x
We
17.88xl017
17.88xl017
17.88xl017
- O.lSxlO17
0.42xl017
3.51X1017
3.32xl017
4.88xl017
S.lOxlO17
17.88xl«17
17.8«xl017
17.Mxl017
- 1.12xlP17
Il.34xl017
0.32xl017
I.96X1017
4.8Ixl017
S.lOxlO17
.2)
in"
AW
-14.34
-14.34
-14.34
3.72
3.12
0.03
0.22
-1.34
- 1.56
-14.34
-14.34
-14.34
4.46
- 6.80
3.22
2.58
- 1.34
- 1.56
HN03 (30.2)
Wg « 3.54 x 10'7
We AK
- 0.07xl017 3.61
24.58xl017 -21.04
S.lOxlO17 - 4.56
S.SlxlO17 - 2.27
4.85xl017 - 1.31
3.56xl017 - 0.02
4.32xl017 - 0.78
18.38X1017 -14.84
18.23xl017 -14.«
M.WxlO17 -21.77
-17.MX1017 21.44
25. OSxlO17 -21.51
- 2.27xl017 5.11
- 0.2TX101* 3.11
9-llxlO17 - 5.64
24.01xl017 -20.47
18.3«xl017 -14.84
18.23X1017 -14.69
Error
RMS
Wtd.
0.35
0.74
1.76
0.60
0.11
0.02
0.03
O.03
0.03
T.lf
f.42
T.79
1«.3%
0.72
t.71
.74
8.03
a 03
Unwtd .
7.39
12.73
7.52
2.19
1.71
0.02
0.41
7.41
7.3»
!».*
11.
12. m
4.41
4. IS
3.21
10.32
7.41
7.39
-2
Wg given concentration in simulated sample, molecules cm"2
He calculated concentration, molecules cm
AW Wg-Wc (the factor of a power of ten
is omitted), molecules cnr2
See 3.1 for definition of the weighted and
unweighted final transmission errors (last
two columns). Numbers following gas
formulae are computer codes.
A composite sample spectrum IMS obtained by
adding the absorption coefficients as a
function of frequency and calculating the
transmission (See 1.2). The sample spectrum
was then analyzed by MINMYZD to obtain the
calculated concentrations. The error is the
difference between the calculated and given
transmission. (See 3.1.3).
Table 3.5b
-------
Table 3.5b
Final Transmission
Wavenumber H20 (1A1)
Range Wg = 3.73
Code (cm-1) Kc
8888 a,B2 1230.0-1249.8 3.73xl022
2?
b,B2 1250.0-1269.8 3.73x10
c,B2 1270.0-1289.8 3.73xl022
d,B2 1290.0-1309.8 3,73xl022
e,B2 1310.0-1329.8 3.64xl022
77
f,B2 1330.0-1349.8 4.05x10^
g,B2 1350.0-1369.8 4.08xl022
h,B2 1370.0-1389.8 4.04xl022
i,B2 1380,2-1400.0 4.04xl022
8888 a,B3 1230.0-1249.8 3.73xl022
b,B3 1250.0-1269.8 3.73xl022
c,B3 1270.0-1289.8 4.0nxl022
d,B3 1290.0-1309.8 Z.rjxlO22
e,B3 1310.0-1329.8 4.34xl022
f,B3 1330.0-1349.8 4.00xl022
g,B3 1350.0-1369.8 3.64xl022
??
h,B3 1370.0-1389.8 3.68x10
i,B3 1380.2-1400.0 3.67xl022
x 1022
AK
.00
.00
.00
.00
.09
-.32
-.35
-.31
-.31
.00
.00
-.32
.20
-.61
-.27
.09
.05
.06
N20 (4A1)
Wg = 1 . 08
We
1.09xl018
1 0
1.09xlOiB
l.OSxlO18
l.OSxlO18
0.31X1018
1 R
1.21xl018
1.21X1018
l.ZlxlO18
1.21X1018
1.09xl018
1.09xl018
1.32xl018
3.02xl018
2.93X1018
1.21xl018
1.21xl018
1R
1.21X1018
1.21X1018
x in1 8
Aw
-0.01
-0.01
0.00
o.oo
0.77
-0.13
-0.13
-0.13
-0.13
-0.01
-0.01
-0.24
-1.94
-1.85
-0.13
-0.13
- .13
- .13
CH i. (5A1)
Wg = 1,08
Kc
l.OSxlO19
1 Q
1.08x10 M
l.OSxlO19
l.OSxlO19
1.16X1019
1 Q
1.06x10
1.09X1019
1.48X1019
1.48X1019
i.nsxin19
1.09xl019
1.22X1019
0.86X1019
0.82X1019
1.42xl019
1.07xl019
1Q
1.48X101'
1.48X1019
x 1019
AK
.00
.00
.00
.00
-.08
.02
-.01
-.40
-.40
.00
.00
-.14
.22
.26
-.34
.01
-.40
-.40
S02 (20.2)
Wg = 3.54 x 101'
We
SO. deleted
library for
AW
from
runs
8888a,B2, through
8888i,B2
17.88xl017
17.88xl017
17.88xl017
14.24xl017
ll.SOxlO17
lO.lOxlO17
6.14xl017
1 7
. 4.86x10
S.lOxlO17
-14.34
-14.34
-14.34
-10.70
- 8.26
- 6.56
- 2.60
- 1.32
- 1.56
HNOj (30.2)
Wg = 3.54 x
We
0.07xl017
i 7
0.10x10
3.55xl017
3.54X1017
S.19xl017
1 7
4-SOxlO1'
7.82X1017
18.25X1017
18.23X1017
HMO. deleted
1017
AW
3.47
3.44
- 0.01
0.00
- 1.65
- 0.96
- 4.28
-14,71
-14.69
from
j
library for runs
8888a,B3 through
8888i,B3
Error
RMS
Wtd.
0.03
0.06
0.03
0.03
0.11
0.36
0.06
0.02
0.02
0.03
0.06
0.72
1.76
0.19
0.81
0.04
0.03
0.03
Unwtd.
1.74
1.72
0.01
0.00
0.91
0.51
2.15
7.36
7.35
7.17
7.17
7.17
5.44
4.25
3.29
1.30
0.69
0.81
Wg » given concentration in simulated sample, molecules rnT
We calculated concentration, molecules cm"2
AW « Wg-Wc (the factor of a power of ten
is omitted), molecules cm-2
See 3.1 for definition of the weighted and
unweighted final transmission errors (last
two columns). Numbers following gas
formulae are computer codes.
A composite sample spectrum was obtained by
adding the absorption coefficients as a
function of frequency and calculating the
transmission (See 1.2). The sample spectrum
was then analyzed by MINMYZD to obtain the
calculated concentrations. The error is the
difference between the calculated and given
transmission. (See 3.1.3).
-------
Table 3.6
GAS CONCENTRATIONS FROM SIMULATED COMPOSITE SPECTRA (SET 2)
00
Final Transmission
Code
9999 a, A
b,A
c,A
d,A
9999 a,Al
b.Al
c,Al
d.Al
9999 a,A2
b.Al
C.A2
d,A2
9999 a, A3
b,A3
C.A3
d.A3
Wavenuraber
Range
(cm-1)
2100.0-21S9.4
2160.0-2219.4
2220.0-2279.4
2240.4-2299.8
2100.0-2159.4
2160.0-2219.4
2220.0-2279.4
2240.4-2299.8
2100.0-2159.4
2160.0-2219.4
2220.0-2279.4
2240.4-2299.8
2100.0-2159.4
2160.0-2219.4
2220.0-2279,4
2240.4-2299.8
Wg given concentration
H20 (IB]
Wg 3.73 x
We
3. 65x1 021
3.73X1021
49.64xl021
11. 74x1 O21
3.68xl02i
3.72xl021
3.68xl021
3.69xl021
3.72xl021
l.SSxlO21
8! -xlO21
11 74x1 O21
12.55xl021
36.19xl021
3.59xl021
3.61xl021
in simulated
1) C02 (2B1)
1021
AW
0.08
0.00
-45.91
- 8.01
0.05
0.01
0.05
0.04
0.01
2.15
-77.51
- 8.01
- 8.82
-32.46
0.14
0.12
sample
Kg - 1.08 x 1021
We
CO- deleted
library for
AW
from
runs
9999a,A through
9999d,A
1.29xl021
l.OSxlO21
l.OBxlO21
l.OSxlO21
l.OSxlO21
5.76xl012
O.SOxlO21
0.92xl021
ll.llxlO21
-5.34xl021
l.OSxlO21
l.OSxlO21
- 0.21
0.03
0.00
0.00
0.03
- 4.68
0.28
0.16
-10.03
6.42
0.00
0.00
03 (3B1)
Wg = 2.69 x 1017
We AW
5.07xl017 - 2.38
5.94X1017 - 3.25
5.94xl017 - 3.25
4.40xl017 - 1.71
0, deleted from
library for runs
9999a,Al through
9999d,Al
2.64xl017 O.OS
5.94xl017 - 3.25
5.94xl017 - 3.25
4.40X1017 - 1.71
45.48xl017 -42.79
S.94xl017 - 3.25
5.94xl017 - 3.25
5.94X1017 - 3.25
NjO (4B1)
Wg - 1.08 x 101 *
We AK
0.64xl018 o-44
l.OSxlO18 o-OO
0.93xl018 0-15
O.SlxlO18 0.57
O.SlxlO18 0.27
l.OSxlO18 0.00
l.OSxlO18 0.00
l.OSxlO18 0.00
N_0 deleted from
£.
library runs
9999a,A2 through
99994, A2
44.14x10 -43.06
1.79xl018 - 0.71
1.14xl018 -0.06
l.OSxlO18 0.00
CO (6B1)
Wg « 1.88
We
1.93xl019
1.88X1019
1.78xl019
2.98xl019
1.91xl019
1.88xl019
1Q
1.88x10
2.98xl019
l.SSxlO1^
3.88X1019
14.11xl019
2.98X1019
CO deleted
library for
x 1011
AW
- 0.05
0.00
0.10
- 1.10
- 0.03
0.00
0.00
- 1.10
- 0,01
- 2.00
-12.23
- 1.10
from
runs
9999a,A3 through
9999d,A3
Error
PMS
Ktd .
0.43
0.03
24.22
63.50
0.20
0.03
0.03
0.03
0.04
24.13
33.67
29.44
16. SS
18.60
2.06
0,03
Unwtd .
1.21
1.63
23.01
4.14
0.17
0.02
0.03
0,55
0.03
3.20
39.27
4.13
31.08
16.63
1.63
1.63
' niolecules cm A composite sample spectrum was obtained by
We » calculated concentration, molecules cm"*
adding the absorption coefficients as a
AW « Wg-Wc (the factor of a power of ten
is omitted), molecules cm-2
See 3.1 for definition of the weighted and
unweighted final transmission errors (last
two columns). Numbers following gas
formulae are computer codes.
function of frequency and calculating the
transmission (See 1.2). The sample spectrum
was then analyzed by MINMYZD to obtain the
calculated concentrations. The error is the
difference between the calculated and given
transmission. (See 3.1.3).
-------
4. CALCULATION OF THE SELECTED GAS SAMPLE CONCENTRATIONS
FROM SPECTROPHOTOMETER DATA
The results presented in 3 were obtained from simulated data
samples produced by calculating composite spectra from selected gas
spectra for specified concentrations. While these are an adequate test
of the computer program, they do not test the overall system of
spectrophotometer and computer algorithm. Accordingly, composite spectra
were prepared by Dr. Wm. Herget, EPA, Research Triangle Park, NC, on
the ROSE* system described in 1.1 and A3. These data were recorded
on digital magnetic tape, and this information became the input to the
EPAGAS program described in A2. In addition, for comparison, some data
were obtained on a Fourier transform spectrometer (FTS).
The data were not obtained from a long atmospheric path because of
the difficulty of obtaining, in that case, an independent measure of the
gaseous concentrations. Instead, a cell containing approximately measured
concentrations of selected gases was placed in the path. The data were
normalized, approximately, by removing the cell and running over the same
spectral range to determine the spectra of ambient H_0 and C0~.
Several problems arose in connection with normalization and
wavelength calibration. These are discussed in 1.1 and A2.5, respectively.
The result of normalization and wavelength uncertainties in the spectro-
photometer output and the added effect of mismatched resolutions of the
spectrophotometer and the spectra library result in some small inaccuracies
in the calculated concentrations. These problems can all be solved by
better calibration methods, and effort is underway to do so. The ROSE
spectrophotometer resolution is a function of wavelength, which is
characteristic of grating spectrometers. It can be prevented by a non-
linear slit drive mechanism or by using a Fourier Transform Spectrometer
(FTS). The FTS can be expected to have a linear frequency scale, but the
*Remote Optical Sensing of Emissions.
39
-------
data may need a shift of frequency to correct for an uncertain zero calibra-
tion. The ROSE spectrophotometer data and the data from the FTS correspond
well in general. However, in a number of regions, the correspondence is
sufficiently poor to cause errors in calculated concentrations. Data with
such discrepancies may be seen by comparing F4.1, F4.2, and F4.3 with
F4.9 (a-c).
However, in spite of these difficulties, th« computational method
performs well, and as shown below, the results of the calculation are of good
accuracy when precautions are taken to calibrate the data.
Several calculations from experimental data are described in the
following sections, utilizing different gases in different frequency ranges.
1300.S 1302.S 1304.7 130C.C 13084
FREQUENCY (cm*1)
Figure 4.1
SPECTRUM OF SOo AND CH4 OBTAINED
BY FTS (1300.8-1310.5 cm-1)
13104
40
-------
zo
us
us
IjO
0.8
0.4
0.2
0.0
1310.6 1312.6 1314.5 13164 1318.4 1320.4
FREQUENCY (cm'1)
Figure 4.2
SPECTRUM OF S02 AND CH4 OBTAINED
BY FTS (1310.6-1320.4 cnr1)
2.0
US
US
o
»«
I 0-8
0.6
0.4
0.2
0.0
1320.5 1322.5 1324.5 132&5 1328.5
FREQUENCY (cm'1)
Figure 4.3
SPECTRUM OF SOe AND CH4 OBTAINED
BY FTS (1320.5-1330.5 cm-1)
1330-5
41
-------
4.1 EFFECT OF MISMATCHED RESOLUTION OF EXPERIMENTAL DATA
AND LIBRARY SPECTRA
The resolution of spectroscopic instruments has a large effect on the
appearance of molecular spectra. It is well known that instrumental
resolution can change the shape of a theoretical spectrum so thoroughly
that recognition is very difficult. Appearance of the spectrum is of
great importance in any best-fit program because it compares two curves,
the one obtained from the spectral library with the one from the data.
Figure 4.4 shows an example of a theoretical spectrum degraded by an
instrumental bandwidth (Calfee, 1966).
The gross changes of recognizable features seen in F4.4 are not
usually a problem in spectral analysis of pollutants because the mismatch in
resolution is not as large. However, considerable error in analysis can
be introduced by such mismatch. Figure 4.5 is a position of the reconstructed
spectrum of carbon monoxide reproduced from data obtained on the ROSE
19 -2
system. The concentration of the sample was 1.5 x 10 molecules cm and
the spectral slit width was 3 cnf according to the accompanying data.
Difficulty was encountered in analyzing the spectrum, hence a series of
computed spectra were made with various amounts of gas and with different
slit widths. Figure 4.6 was computed using a spectral slit width
of 3 cm"1 and with 1.5 x 1018, 7.5 x 1018, and 1.5 x 1019 molecules cm"2 to
see if a different amount of gas in the cell could be the cause of the
discrepancy. This was not the case because none of the computed spectra
match the data.
Next, a calculation was made with the same quantities of gas, but
with a narrower slit width (2.75 cm" ). Figure 4.7 shows the result of
this calculation. This result approaches ;loser to the experimental
spectrum. Figure 4.8 was made by using 2.3 cm" as the spectral slit
width. Here the agreement between the computed spectrum and the experimental
19 -2
curve is very good for a concentration of 1.5 x 10 molecules cm . This
analysis points up the significance of the spectral slit width in the
structure of the spectrum and the best-fit procedure. (The slit width as a
function of frequency is shown in A3.)
42
-------
THEORETICAL SPECTRUM
P= 0-1 Qtm
00
90
801
z.
o
<7i JO
01
3j 60
j
£ 50
50 i
20
1C
3900
/I
w = OOlO P' cm H.,0
80 70
FREQUENCY cm"1
DEGRADED SPECTRUM
= Ola1m
UPPER CURVE w = OOOI pr cm H20 LOWER CURVE w = OQIO pr cm H20 SPECTRAL SLIT (20)= 1-0 cm'1
r\
3850
FREQUENCY cm"
Figure 4.4
COMPARISON OF A THEORETICAL SPECTRUM (HIGH RESOLUTION)
OF WATER VAPOR AND THE SPECTRUM OBSERVED THROUGH A
SPECTROMETER WITH 1.0 cnr1 SLIT WIDTH (IN MANY REGIONS
THE COMPARISON IS VERY DIFFICULT)
43
-------
2.0
LB
US
kj 1-2
d
§ °-8
H
0.6
OA
C.2
0.0
ZO
L3
o
0.8
0.6
a4
0.0
t. 1.
a. T.
3. 1.5 « w
1*
2068.92072.9 2076.9 2080.9 208$^) 208941
FREQUENCY (cm'1)
Figure 4.5
EXPERIMENTAL SPECTRUM, CO,
NOMINAL RESOLUTION 2.5 cm-1
OBTAINED ON ROSE SYSTEM
(Compare with F4.6)
2068.9 2072.9 2076.9 208O.9 208&0 2089*
FREQUENCY (cm*1)
Figure 4.6
CALCULATED SPECTRUM, CO,
RESOLUTION 3.0 cm-1
(Compare with F4.5)
44
-------
zo
US
US
IA
"
0.6
0.4
0.2
0.0
1. 1.1 a U
J. T.»«10W
>. Lliiol*
lM*l«cnl«> cm'2)
206&8 2072.8 2076.9 2080.9 2085.0 2089-0
FREQUENCY (cm'1)
Figure 4.7
CALCULATED SPECTRUM, CO,
RESOLUTION 1;0 Ciir1
(Compare with F4.5)
2.0
US
us
< 1.2
t J
0-8
0.6
0.4
C.2
0.0
20684 2072.8 2076.9 2080.9 2085.0 2O89.O
FREQUENCY (cm*1)
Figure 4.8
CALCULATED SPECTRUM, CO,
RESOLUTION 2.3 cm-1
(Compare with F4.5)
45
-------
The spectra shown in F4.9 (a-c) were reconstructed from EPA
ROSE data of a combination of SO (7.56 x 1018 molecules cm"2) and CH.
19 -2 -1
(1.23 x 10 molecules cm ) with a spectral slit width of 1 cm . Figure
4.10 (a-c) shows the computed spectra for the same regions with the same
amounts of gas and the same slit width used for the calculations. The
most obvious disagreement appears between the curves of F4.9b (experimental)
and F4.10b. It would appear that some gas other than SO and CH. (or noise)
(see 4.4) was in the sample to produce the additional structure. An
examination of the entire series shows a general agreement, but again the
experimental curves show an indication of a higher resolution than that
used for the computed spectra. It would also appear that the use of a
lower value of concentration would produce a computed spectra which would
come closer to matching the experimental curves. Probably both a lower
concentration and a higher resolution would produce a spectrum which would
agree with the experimental spectrum.
This discussion illustrates the necessity for an accurate
normalization procedure and precision in determining instrumental
resolution to avoid ambiguity in spectral analysis. The resolution
of spectrometers generally changes over their spectral range. This makes
it necessary that either the experimenter or the analyst must be able to
compensate for the effect.
46
-------
2.0
US
U6
J- 1.0
d
< 0.8
0.6
OA
0.2
0.0
1300.8 1302.8 1304.7 1306.6 1303.6 1310*
FREQUENCY (cm'1)
Figure 4.9a
EXPERIMENTAL SPECTRUM OF
COMBINATION OF S02 AND CH,,
NOMINAL SLIT WIDTH 1.0 cm-1
(ROSE SYSTEM) (1300.8-1310.5 cm-1)
(Compare with F4.10a,b,c)
2.0
18
US
1/4
1.0
0.6
0.4
0.2
0.0
1310.6 1312.6 1314.5 13164 1313.4 1320.4
FREQUENCY (cm'1)
Figure 4.9b
EXPERIMENTAL SPECTRUM OF
COMBINATION OF S02 AND CH4,
NOMINAL SLIT WIDTH 1.0 cm-1
(ROSE SYSTEM) (1310.6-1320.4 cnT1)
(Compare with F4.10a,b,c)
-------
00
2.0
18
U6
o
H 1>2
',.* *j<>
" 'f
c °'8
H
0.6
0.2
0.0
1320.5 1322.5 1324.5 13263 1328.5
FREQUENCY (cm'1)
Figure 4.9c
EXPERIMENTAL SPECTRUM OF
COMBINATION OF SO? AND CH4,
NOMINAL SLIT WIDTH 1.0 cm-1
(ROSE SYSTEM) (1320.5-1330.5 cm-1)
(Compare with F4.10a,b,c)
1330.5
2.0
L8
1J6
1.0
as
0.6
0.4
0.2
0.0
1300.3 1302.8 1304.7 130S.6 1300.6 13103
FREQUENCY (cm'1)
Figure 4.10a
CALCULATED SPECTRUM, S02 AND CH4
(RESOLUTION 1.0 cm-1)
(1300.8-1310.5 cm-1)
(Compare with F4.9a,b,c)
-------
2.0
U8
U6
-u
«£>
0.6
O4
0.0
2.0
B -
Z
SJ 0.8
0.6
0.4
0.2
0.0
1310.6 1312.6 1314.5 1316.5 1318.4 132O.4
FREQUENCY (cm-1)
Figure 4.1Ob
CALCULATED SPECTRUM, S02 AND CH4
(RESOLUTION 1.0 cm-1)
(1310.6-1320.4 cm-1)
(Compare with F4.9a,b,c)
1320.5 1322.5 1324.5 13263 1328.5 1330.5
FREQUENCY (cm'1)
Figure 4.10c
CALCULATED SPECTRUM, SO? AND CH/
(RESOLUTION 1.0 cm-1)
(1320.5-1330.5 cm-1)
(Compare with F4.9a,b,c)
-------
4.2 CALCULATION OF THE CONCENTRATION OF CO
A composite spectrum of H20, CO CO, CH4, and SO was prepared on
the ROSE spectrophotometer over the spectral region from 1250 to
1330 cm" and in the range 2050 to 2150 cm'1. Because only a relative
intensity was obtained, normalization proceeded by obtaining a background
run containing only spectra of E^O and CO-. The ratio of the two spectra
gave the normalized spectrum I/I . The EPAGAS program was then used to
calculate the concentration of CO in three spectra regions, 2073-2094,
2094-2114, and 2114-2132 cm"1. T4.1 shows the results of the calculation.
The spectral regions are listed in the second column. Several initial
trial concentrations were used to determine the mathematical stability of
the process. As can be seen for each region, different starting estimates
made little change in final calculated concentrations. However, as
discussed in 4.1, the resolution used in the spectral library strongly
affects the accuracy. Part of the calculations summarized in T4.1 used
the nominal spectral library resolution (3.0 cm) for this region. However,
as pointed out in 4.1, the resolution used in preparing the data was
actually closer to 2.3 cm" . Nine values of concentration were calculated
assuming a resolution of 2.3 cm" . These values average 1.26 x 1019
_2
molecules cm ; the sample concentration measured chemically was
19 2
1.5 x 10 molecules cm" . It is believed that this agreement is excellent
in view of the calibration and normalization problems.
Note that the spectral resolution was carefully determined in the
region 2073-2094 cm" and that the errors are the least there, increasing
as the calculation is performed for more distant regions.
By using program DEGRADE (Deutschman and Calfee, 1967) before the
fitting routine in EPAGAS, the actual r;-olution of the sample was deter-
mined to be 2.3 cm . This choice of resolution was verified by the
results given in T4.1.
50
-------
Table 4.1
CALCULATION OF THE CONCENTRATION OF A SAMPLE OF CO BY
EPAGAS TO AN EXPERIMENTAL SPECTRUM AT SEVERAL SPECTRAL
Frequency
Region Correction
Code
la
b
c
d
c
2a
b
c
d
e
3a
b
c
d
e
Average
Average
Spectral
2073.0-2094.0
2073.0-2094.0
2073.0-2094.0
2073.0-2094.0
2073.0-2094.0
2094.0-2114.0
2094.0-2114.0
2094.0-2114.0
2094.0-2114.0
2094.0-2114.0
2114.0-2134.0
2114.0-2134.0
2114.0-2134.0
2114.0-2134.0
2114. 0-2. 34.'0
CO Concentration
CO Concentration
Av
0.0
0.0
c.o
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
and Error
and Error
Resolution
3.0
3.0
2.3
2.3
2.3
3.0
3.0
2.3
2.3
2.3
3.0
3.0
2.3
2.3
2.3
cn
Trial (Est.)
Concentration
(molecules cm~
l.SPxlO19
I Q
1.88x10
1 O
5.38x10
l.SOxlO19
1 R
1.88x10
l.SOxlO19
1 P
1.88x10
1 8
5.38x10
1Q
1.50x10
1.88xl019
l.SOxlO19
1.88xl01P
5.38xl018
l.SOxlO19
1.88X1019
for 3.0 resolution: 8.44xl01S
for 2.3 resolution: 1.26xl019
CO Calculated
Concentration
2) (molecules cm-2)
7.33xl018
1 Q
7.32x10
1 Q
1.18x10
1.18xl019
10
1.18x10
9.27xl018
1 D
9.29x10
1.26X1019
1.26xl019
1°
1.26x10
8.71xl018
s.rixio18
1.35xl01P
10
l.SSxlO1-
1.35xl01P
; Error 2.24.
; Error 1.87.
Final
Transmi ssions
2.51
2.51
1.25
1 .-25
1.25
2.24
2.24
2.18
2.18
2.18
1.96
1.96
2.19
2.19
2.19
Average" for a
Given Resolution
and Region
Cal.
Cone. Error
1 a
7.33x10 2.51
10
1.18x10 1.25
1 C
9.28x10 2.24
1.26xl019 2.18
I Q
8.71xl0lc 1.96
1.35xl019 2.19
Sample Concentration: l.SOxlO15 mol cm"2 (Obtained by Chemical Analysis).
The calculated concentration with the least error (1.26xl019) compares well with the sample concentration
NOTE: Powers of ten suppressed in errors quoted above. Concentrations are in molecules cm--.
-------
Figures 4.11 and 4.12 show, respectively, the data produced by the
ROSE spectrophotometer and the best fitting calculated spectrum. The
slight differences in this case are probably due to noise. Although the
shapes of these spectra match quite well, the absolute transmission
differences account for the slight discrepancy in calculated: and experimental
concentration. These absolute transmission differences are unresolvable
without better normalization.
c;
1^
<
2.0
IB
1A
IjO
Q&
0.6
0.4
0.2
0.0
2.0
_L
_L
J.
J_
J I I L
Figure 4.11
SPECTRUM OF CO CALCULATED BY
PROGRAM DEGRADE (RESOLUTION
2.3 cnr1)'
(Compare with £4.12)
2068.92072.9 2076.9 2080.9 20854) 20894)
FREQUENCY (cm*1)
n 14
5
H "
0.4
012
0.0
206&9 2072.9 2076.9 2080.9 20854)
FREQUENCY (cm'1)
Figure 4.12
SPECTRUM OF CO, EXPERIMENTAL
DATA FROM-ROSE
SPECTROPHOTOMETER
(Compare with F4.ll)
208*0
52
-------
The effect of resolution on the goodness of fit in EPAGAS may be
seen in F4.13. In the first of these, the dashed curve represents the
given experimental data, the full line is the calculated value. The
results of the calculation are given in (resolution 3.0 cm" ) T4.1. These
curves should be contrasted with the first calculation listed in T4.1 using
2.3 cm" (line 7), shown in F4.14. The much better fit of given and
calculated data results in a more accurate value of concentration and once
again shows the importance of resolution matching.
The efficiency of the program in performing a "best fit" may be seen
in F4.14. Notice that the two curves are not sufficiently alike to be
fitted exactly. However, the optimization of the adjustable parameters
of the calculated spectrum has clearly resulted in an excellent fit in
a least-squares sense. An understanding of the change of a typical spectrum
with a small change of concentration may be seen in F4.15. Clearly good
accuracy and low noise of the instrument data are required for accurate
analysis.
53
-------
Tilf-?1 ^
»*-uii ^
OX)
0.2
0.4
1*1
-2073.4
-2075.6
2077.6
2079.8
2081.8
2083.9
tj -2085.9
2087.9
-2084.9
2O92.0
2094.0
0.6
I
0.8
1.0
Experimental
Computed
Gas: 6
Amount: 7.317 + 01C
Percent: 1OO
Figure 4.13
COMPUTED BEST-FIT TRANSMISSION OF CO (EPAGAS PROGRAM)
RESOLUTION 3.0 cm-1
54
-------
0.0
0.2
0.4
0.6
i
0.8
1.0
2073.4
2075.6
- Experimental
^ Computed
l&l
2077.6
- 2079.8
2081.8
2083.9
-2085.9
IL.
2087.9
2089.9
-2092.0
Gas: 6
Amount: 1.183 +019
Percent: 100
-2093.9
Figure 4.14
COMPUTED BEST-FIT TRANSMISSION OF CO (EPAGAS PROGRAM)
RESOLUTION 2.3 cm-1
55
-------
CO CMC 2.SKM+01*
CO CMC 00+01*
20004 2040.4 2OSO3 2120.3 216O2
FREQUENCY (cm'1)
Figure 4.15
CHANGE IN SPECTRUM OF CO
FOR SMALL CHANGE OF CONCENTRATION
(Calculated by Program DEGRADE)
2200.2
4.3 CALCULATION OF THE CONCENTRATION OF S02 IN THE PRESENCE OF CH4
The gas sample cell used in obtaining data for C02, as discussed in
4.2, contained also SO- and CH.. These gases have significant spectra in
-1
the region from 1150 to 1420 cm (see Fl.l (a-c)). However, the data have
good signal-to-noise ratio only in the ranges 1300-1330 cm" and 1250-1300
cm" . We first consider the 1300-1330 cm" range. In this range,
the spectrophotometer data shows evidence of incomplete elimination of
water spectra or noise, hence the data from the Fourier transform
spectrometer were utilized. In this ri..nee, the spectral structure of CH4
does not lend itself to fitting uniquely, and we here seek only to find
the concentration of S02 in the presence of interference due to CH4.
56
-------
Table 4.2 shows the results of concentration calculations by MINMYZD for
several different conditions of initial trial concentration and frequency
corrections. The stability of the calculated concentration of S0_,
regardless of starting point, is apparent. The frequency corrections listed
for the various regions are linear frequency shifts due to an uncertainty
of the zero point in the Fourier transform.
As expected, the concentrations calculated for CH. here vary widely
and are untrustworthy. However, the presence of CH. does not radically
affect the value of S0_ concentration computed in any of the spectral
-2
regions. The average value of S02 concentration molecules cm is very
close to the chemically determined concentration in spite of the
many uncertainties in the data parameters. It is particularly encouraging
to note that the concentration of the sample lay outside the span of library
concentrations, but the program EPAGAS, by use of the subroutine FINDTRN,
was able to extrapolate properly. The library does not include such
large concentrations because it is not expected to be so large in field
experiments, even in highly polluted geographical regions.
4.4 CALCULATION OF THE CONCENTRATION OF CH,
4
An attempt was made to analyze the spectra for the concentration of
CH.. The results of this calculation were not as successful as the
4
determination of the concentrations of CO and SO . The cause of the
discrepancy is either noise or unknown interfering spectra (or both) in
the data. The problem is illustrated in F4.16. There we see that the
optimum fit of CH. to the data still allows large errors because structure
is found in the data which does not come from the CH. molecule.
The source of the structure very likely lies among the following:
a. The library spectra of CH4 is incorrect.
b. There may be small SO,, lines from an adjacent band not
found in the library.
c. The sample contained an unknown gas.
d. The spectrum contains instrumental noise.
e. Water lines are incompletely cancelled.
57
-------
Table 4.2
CALCULATION OF THE CONCENTRATION OF THE SAMPLE OF S02 IN THE PRESENCE OF CH
BY APPLICATION OF EPAGAS TO AN EXPERIMENTAL SPECTRUM
in
oo
Code
la*
b
c*
d
e
f
2a*
b*
c
d
e
f
3a
b
c
d
e
f
Spectral
Region
1300.
1300.
1300.
1300.
1300.
1300.
1310.
1310.
1310.
1310.
1310.
1310.
1320.
1320.
1320.
1320.
1320.
1320.
0-1310.
0-1310.
0-1310.
0-1310.
0-1310.
0-1310.
0-1320.
0-1320.
0-1320.
0-1320.
0-1320.
0-1320.
0-1330.
0-1330.
0-1330.
0-1330.
0-1330.
0-1330.
Frequency
Correction
Av
0
0
0
0
0
0
-------
Table 4.2 (cont.)
in
t£>
Code
4a
b
c
d
e
f
5a*
b*
c
d
e
f
Average
for all
Samp 1 e
Spectral
Region
1300.0-1330.0
1300.0-1330.0
1300.0-1330.0
1300.0-1330.0
1300.0-1330.0
1300.0-1330.0
1310.0-1330.0
1310.0-1330.0
1310.0-1330.0
1310.0-1330.0
1310.0-1330.0
1310.0-1330.0
Concentrations
Regions
Concentrations:
Trial (Est.) Calculated
Frequency Concentrations Concentrations
Correction (molecules cm- ) (molecules cm-2)
Av
0.0
0.0
-0.9
-0.9
-1.2
-1.2
0.0
0.0
-,0.9
-0.9
-1.2
-1.2
CH i, S02 CH i, S02
7.56X1018 1.23xl019 S.llxlO18 6.54xl018
1.23xl019 7.56xl018 S.lSxlO18 6.56xl018
7.56xl018 1.23xl019 S.OOxlO18 6.84xl018
1.23xl019 7.56X1018 S.OOxlO18 6.84xl018
7.56xl018 1.23xl019 S.OSxlO18 7.04xl018
1.23X1019 7.56xl018 S.OSxlO18 7.04xl018
7.56X1018 1.23xl019 1.71xl018 6.20xl018
1.23xl019 7.56xl018 1.91xl018 6.20xl018
7.56xl018 1.23xl019 9.62xl018 6.38xl018
1.23xl019 7.56xl018 9.62xl018 6.38xl018
7.S6xl018 1.23X1019 1.06X1019 6.52xl018
1.23xl019 7.56X1018 l.OSxlO19 6.52xl018
Averages for a
Final Given Frequency
Transmission and Range
Error
20.05
20.05
11.49
11.49
9.22
9.22
9.77
9.77
6.34
6.34
6.05
6.05
CHi, S02 Error
S.lSxlO18 6.55xl018 20.05
1 p I Q
8.00x10 6.84x10 11.49
S.OSxlO18 7.04X1018 9.22
-
9.62xl018 6.38X1018' 6.34
1.07xl019 6.52xl018 6.05
Resolution CHA SO, Error
0.
-0.
-1.
CH4 = 1.
SO, = 7.
0 4.88xl018 7.00xl018 20.60
9 6.88xl018 7.16xl018 9.55
2 5.39xl01S 6.75xl018 8.07
, __ . _
23x10 molecules cm (Obtained by Chemical Analysis)
18 -2
56x10 molecules cm
* Mot used in Averages
-------
0.0
0.2
0.4
0.6
1.0
it.
-1279.8
-1280.9
- -1281.9
-1283.0
. .1284.0
2£ 1285.0
Experimental
. Computed
1286.0
-1287.1
--1288.1
12894
1290.1
G«s:S
Amount: 8.244+018
Percent: 100
Figure 4.16
COMPUTED BEST-FIT TRANSMISSION OF CH4
COMPARED WITH THE GIVEN EXPERIMENTAL TRANSMISSION
60
-------
We comment on these possibilities:
a. The library spectrum of CH has been checked by many
experimenters. It is possible, but considered unlikely,
that it is incorrect.
b. A chemical analysis by EPA of the sample did not
».
indicate any other gases. The analysis method should
be reviewed to see if it would reveal other gases
than those known to be in the sample.
c. The spectrum furnished by EPA from the ROSE
spectrophotometer clearly has noise in some regions
of the spectrum. The section analyzed here appeared
less noisy than any other by visual inspection, but if
the noise were low frequency drifts, it would not be seen
as obvious noise.
d. An initial suspicion that water absorptions were
responsible was discounted by the investigation
described below.
In summary, the reason for the discrepancy between the value of
1.23 x 10 9 molecules per cm2 of CH4 and the value of 8.22 x 1018 calculated
by EPAGAS is not known at this time. Both instrumental uncertainties and
possible library inadequacies must be investigated. It should be noted
that even with these difficulties, the values differ by a factor of only 1.5.
Table 4.3 shows a summary of calculations done to determine the
concentrations of CH. and H^O in the spectrum. The water values are low
(for water) and are scattered, when compared in different frequency
regions, indicating that the spectra contains very little water absorption.
The CH. values are, on the other Land, quite consistent. Figure 4,16
shows the impossibility of fitting the given spectrum with a CH.-only
absorption spectrum. Structure in the given spectrum cannot be matched by
any concentration of the CH. spectrum.
61
-------
Table 4.3
CALCULATION OF THE CONCENTRATION OF A SAMPLE
OF CH4 AND H20 BY APPLICATION OF EPAGAS TO
AN EXPERIMENTAL SPECTRUM
Frequency
Region Correction
Code Spectral Av
la
b
c
d
2a
b
c
d
1250.
1260.
1270.
1280.
1250.
1260.
1270.
1280.
Averages for
and a given
Correction
0-1260.0
0-1270.0
0-1280.0
0-1290.0
0-1260.0
0-1270.0
0-1280.0
0-1290.0
all regions
frequency
-0.6
-0.6
-0.6
-0.6
-0.9
-0.9
-0.9
-0.9
Av
-0.6
-0.9
3
3
3
3
3
3
3
3
Trial (fist.)
Concentrations
(molecules cm-2)
H20 CHi,
.73xl022
.73xl022
2?
.73x10
.73xl022
22
.73x10
22
.73x10^
22
.73x10
22
.73x10
1.
1.
1.
1.
1.
1.
1.
1.
23xl019
23xl019
23xl019
23xl019
23xl019
23xl019
23xl019
23xl019
-2
6
-1
5
-9
7
4
1
-2
1
Calculated Final
Concentrations Transmission
H20 OU Error
,17xl021
.82xl020
.44xl020
.04xl02°
.14X1020
.09xl02°
.42X1020
.SOxlO20
H2° 20
.82x10^
.04xl02°
7.79xl018
8.46xl018
1 R
8.67x10
I O
7.96x10
1 R
8.03x10
1 R
8.17x10
I Q
6.85x10
1 O
7.90x10
OU
8.22x10
7.74X1018
6.07
6.88
6.73
7.14
5.41
6.14
8.68
10.20
Error
6.71
7.61
-------
5. CONCLUSIONS AND RECOMMENDATIONS
Many geophysical studies require the measurement of constituent
concentrations in situations where the identifying signals from many
components are mixed in the spectrometer. In remote measurement of
aerosols by electromagnetic scatter, in measurement of ocean constituents,
in the measurement of pollutants in the water and in the air, the same
complex combination of identifying spectra occur. Means of separating
and measuring a large number of components of a mixture usually require
analysis by a large computer because the data set is quite large.
Not only must the field data be large enough to give sufficient
statistically independent information, but the identifiers (spectra)
require a large computer memory to avoid breaking the problem into too many
pieces.
The calculational algorithms utilized in the program EPAGAS have
been shown in this report to determine efficiently and accurately the
concentrations of gas from a complex long-path infrared spectrum. Two
important requirements are:
1. The spectral library must be accurate in wavelength,
pressure and temperature broadening parameters,
and spectral line or band intensity.
2. The spectrophotometer used to obtain the field data
must accurately measure the spectral intensity (at
least relatively), and the wavelength and the spectral
resolution must match the resolution of the
spectra in the library.
The first requirement is met (with qualifications detailed in A2)
by the spectral library assembled for this report (Fl.l(a-c)), but
adequate spectral data analysis of complex spectra is rare in the
scientific literature, as discussed in Al. Other gases than those
considered in this report are important for pollution monitoring. The
technological changes in the United States resulting from the fuel
shortage and the expected mineral shortage will require new industrial
63
-------
developments such as oil shale processing. These can be expected to
result in new effluents whose spectra will most likely not be well
known. Indeed, only a small portion of atmospheric pollutants presently
in industrial areas are well enough known spectroscopically.
Thus, it is emphatically recommended that continuing support be
given to competent spectroscopists to obtain the basic spectroscopic
data needed for air and water pollution. Absorption spectra, Raman
spectra, and emission spectra are needed for the various techniques
presently in development.
It is further recommended that vigorous efforts be made to obtain
field data under controlled conditions to develop experience with
the whole system. Initially, the experiments must include, along the
path, point samplers to establish the levels of gaseous pollutants for
comparison with the spectral analysis. Remote sensing systems must not
be considered complete until proven in comparison with more conventional
methods.
64
-------
6. ACKNOWLEDGMENTS
The authors would like to thank Ms. K. Yoshida and Mr. D. Gregory
for their meticulous assistance in processing the data for this report.
Mr. R. Lawrence has provided wise counsel in the use of the MINMYZD
program. Dr. R. Schwiesow has made essential contributions to the library
of spectra and intelligent aid in all phases of the task. Dr. W. Herget
of EPA (Research Triangle Park) has clarified many important problems and
supplied experimental data. The authors express their gratitude to
Ms. B. Stahl for her skill and monumental patience in producing the
manuscript.
65
-------
7. LIST OF SYMBOLS
I Radiation Intensity
K1 Absorption coefficient (cm" )
K Absorption coefficient
(cm molecules" )
W Area density (or concentration)
(cm'2)
P Density absorbing gas per
unit area of radiation
(cm'3)
N Noise power
v Spatial frequency (wave-
number) (cm" )
T Temperature ( K)
_2
P Pressure (dynes-cm or atm)
S(v) Spectrophotometer output
E Error function
Wg Given concentration
(molecules cm )
We Computed concentration
(molecules cm )
CF Confidence factor
X Wavelength (cm)
C Speed of light in vacuum
(cm sec" )
f Frequency (sec" )
S/N Signal-to-noise ratio
66
-------
8. REFERENCES
Abels, L. L., and J. H. Shaw (1966), Widths and strengths of vibration
- rotation lines in the fundamental band of nitric oxide, J_. Mol.
Spectry., 20, 11-28.
Box, G. E. P., and Mervin E. Muller (1958), A note on the generation of
random normal deviates, Ann. Math. Statist., 29, 610-611.
Burch, D. E. (1970), Investigation of the absorption of infrared radiation
by atmospheric gases, Semiannual Technical Report, Jul. 51, 1969 -
Jan. 31, 1970, Philco Ford Corp., Aeronutronic Div., Newport
Beach, CA., Contract F19628-69-C-0263, U4784.
Burriss, R. C., E. Griggs, T. Connolly, and W. J. Veigele (1972),
Land Use Planning for Air Quality in the Pikes Peak Area, Kaman
Sciences Corp., Colorado Springs, CO, Contract USE-CO-08-00-0002,
K-72-946-FR.
Calfee, R. F. (1966), Infrared absorption properties of atmospheric gases,
J_. Quant. Spectrosc. Radiat. Transfer, jS, 221-228.
Calfee, R. F. (1971), A note on terminologies used in gaseous absorption
processes, NOAA Tech. Rept. ERL 211-WPL15, Boulder, CO.
Cossee, P., and J. H. Schachtschneider (1966), Vibrational analysis of
acetone, acetaldehyde, and formaldehyde, J_. Chem. Phys., 44, (1),
97-111.
Deutschman, E. M., and R. F. Calfee (1967), Two computer programs to
produce theoretical absorption spectra of water vapor and carbon
dioxide, ESSA Tech. Rept., IER 31-ITSA 51, Boulder, CO.
Lawrence, R. S., and M. Hallenbeck (1965), A method for obtaining the
parameters of electron-density profiles from topside ionograms,
NBS Tech. Note No. 515, Boulder, CO.
Levin, A., and C. F. Meyer, The infrared absorption spectra of acetylene,
ethylene and ethane, J_. Opt. Soc. Am., 1£, (3), 157-164.
Matthes, W. (1973), "Unfolding of Composite Spectra by Linear Regression,"
Private communication of a paper to be published.
67
-------
McClatchey, R. A., R. W. Fenn, F. E. Volz, J. S. Garing, and J. E. A.
Selby (1971), Optical properties of the atmosphere, AFCRL
Environmental Research Paper #354.
McClatchey, R. A., W. S. Benedict, S. A. Clough, D. E. Burch, R. F. Calfee,
and K. Fox (1973), AFCRL atmospheric absorption line parameters
compilation, AFCRL Environmental Research Paper #434, AFCRL-TR-73-
0096, A.F. Sys. Comm., Hanscom Field, Bedford, MA.
Milne, J. B., and A. Ruoff (1967), Coriolis coupling constants in sulfur
trioxide, J_. Mol. Spectry., 23, 408-415.
Nakagawa T., and V. Morino (1971), Coriolis Interactions in the v. and vfi
band of formaldehyde, J_. Mol. Spectros., 58, 84-106.
Nielsen, H. H., and E. F. Barker (1931), Infrared absorption bands in
hydrogen sulphide, Phys. Rev. 57, 727-732.
Pierson, R. H., A. N. Fletcher, and E. St. C. Gantz (1956), Catalog of
infrared spectra for qualitative analysis of gases, Anal. Chem.,
2£, 1218-1239.
Shapiro, M. M., and H. P. Gush (1966), The collision-induced fundamental
and first overtone bands of oxygen and nitrogen, Can. J_. Phys.,
4£, 949-963.
Shaw, J. H. (1956), Nitric oxide fundamentals, £. Chem. Phys., 24, (2),
399-402.
Slutz, R. J., and J. R. Winkelman (1964), Shape of the magnetospheric
boundary under solar wind pressure, J_. Geophys. Res., 69, (23),
4933-4948.
Smith, L. G. (1940), Fine structure of the 3.35 u band of ethylene,
£. Chem. Phys., £, 798-799.
Smith, W. L., and I. M. Mills (1964), Coriolis perturbation in the infrared
spectrum of ethylene, J_. Chem. Phys.. 40, (8), 2095-2109.
Streigg, M. L., and C. R. Claysmith (1972), Design and construction of a
system for remote optical sensing of emissions, Final Rept., Contract
EPA-CPA 22-69-142, General Dynamics, Convair Aerospace Div.,
San Diego, CA.
68
-------
U. S. Naval Res. Lab., American Petroleum Institute Research Project 44,
(in-house project).
Westwater, E. R., and A. Cohen (1973), Application of Backus-Gilbert
inversion technique to determination of aerosol size distributions
from optical scattering measurement, Appl. Opt., 12, (6), 1340-1348,
69
-------
APPENDIX 1
TABLE OF CONTENTS
Page
Al. Sources and Contents of the Library of Spectra 71
Al.l Motivation for the Library 71
A1.2 Sources and Content 71
A1.3 Programs Used in Processing Spectra 91
70
-------
APPENDIX 1
Al. SOURCES AND CONTENTS OF THE LIBRARY OF SPECTRA
Al.l MOTIVATION FOR THE LIBRARY
As outlined in the body of the report and explained in detail in
A2, EPAGAS is a program which employs a least-squares fitting method
to match an experimentally taken spectrum, say of a polluted atmosphere,
with an internally generated reference spectrum. Through constructing
reference spectra by adjusting the concentrations of gases in its
library, EPAGAS fits the data spectrum to a reference spectrum and
thus determines the concentration of the gases in the data spectrum.
The closer that the library comes to containing spectra of all the
molecules found in the experimental spectra, the more accurate are the
results of EPAGAS. (See 3.1.3.) This is the motivation for assembling
a large and accurate library of spectra.
Physically, the library consists of a set of condensed computer cards
which contain values of the absolute transmittance of each gas at each
concentration at each frequency used by EPAGAS-a total of about 256,000
points.
A1.2 SOURCES AND CONTENT
Spectra for the 10 gases shown in Fl.l were assembed by this
laboratory* over the past 10 years and are summarized in TAl.l. The
spectra of the first seven gases listed in the table (set numbers 1A1
through 7B1) were calculated theoretically using line parameters and
self-broadening constants. SO-, HNO,, and NH_ spectra were taken
experimentally. All assumed room temperature (296 K).
*Wave Propagation Laboratory, ERL/NOAA. These data were assembled from
theoretical computations and experimental data from many sources. See
McClatchey et al., 1973.
71
-------
Table Al.l
SPECTRA FROM TUB WPL DATA SET
N)
SET NUMBER
1A1
1B1
2A1
2B1
3A1
3B1
4A1
4B1
5A1
5B1
6B1
7B1
20.1
20.2
30.1
30.2
40.1
GAS
H20-Water
H20-tfater
C02-Carbon Dioxide
C02 -Carbon Dioxide
Oj-Ozone
0, -Ozone
NjO-Nitrous Oxide
N20-Nitrous Oxide
CH4 -Methane
CH^-Methane
CO-Carbon Monoxide
HCL-Hydrochloric Acid
S02-Sulfur Dioxide
S02-SulfuT Dioxide
HNOj-Nitric Acid
HNOj-Nitric Acid
NH, -Ammonia
WAVENUMBER
REGION
740.
1800.
710.
1810.
970.
2040.
1120.
2130.
1180.
2488.
2000.
2400.
1074.
1299.
824.
1234.
1112,
0-1450.0
0-3200.4
0-1100.0
2-2499.0
0-1100.0
0-2139.0
0-1360.0
6-2599.8
0-1400.0
2-3199.8
4-22S9.6
6-3169.8
8-1231.2
2-1409.2
8- 951.2
8-1383.8
0-1211.4
INCREMENT
(cm'1)
0.2
0.6
0.2
0.6
0.2
0.6
0.2
0.6
0.2
0.6
0.6
0.6
0.2
0.2
0.2
0.2
0.2
RESOLUTION
1.0
3.0
1.0
3.0
1,0
3.0
1.0
3.0
1.0
3.0
3.0
3.0
1.0
1.0
1.0
1.0
1.0
NUMBER OF
POINTS PER
CONCENTRATION
3551
2335
1951
1149
651
166
1201
783
1101
1187
433
1283
783
551
633
746
498
7.
6.
7.
6.
8.
8.
1.
1.
1.
5.
5.
2.
2.
2.
2.
5.
5.
2.
2.
1.
3.
3.
3.
1.
3.
46 xlO21,
77 xlO22
46 xlO20,
77 XlO21
61 xlO20,
61 XlO20,
08 xlO ,
08 xlO18
08 xlO17,
38 xlO17,
38 xlO17,
69 xlO18,
69 xlO19
69 xlO18,
69 xlO19
38 xlO 8,
38 XlO19
69 xlO16,
69 xlO18
06 xlO18,
54 xlO16,
54 xlO18
54 xlO17,
06 xlO17,
54 xlO17,
CONCENTRATIONS
(molecules cm-2)
2.24 xlO22,
2.24 xlO21,
1.08 xlO21,
1.08 xlO21,
1.88 xlO17,
1.88 xlO17,
1.08 xlO18,
1.08 xlO18,
5.38 xlO18,
5.38 xlO18,
1.08 xlO19,
1.08 xlO17,
3.54 xlO18,
1.06 xlO17,
1.06 xlO18,
3.54 xlO17,
1.06 xlO18,
3.73 xlO22,
3.73 xlO21,
1.29 xlO21
1.29 xlO21
2.69 xlO17,
2.69 xlO17,
1.88 xlO18
1.88 xlO18
1.08 xlO19,
1.08 xlO19,
1.88 xlO19,
2.69 xlO17,
1.06 xlO19,
3.54 xlO17,
3.54 xlO18
1.06 xlO18,
3.54 xlO18,
5.22 xlO22,
5.22 xlO21,
5.38 xlO17,
5.38 xlO17
1.08 xlO
1.88 xlO19,
1.88 xlO19,
2.69 xlO19,
1,08 xlO18,
3.54 xlO19
1.06 xlO18,
3.54 xlO18
1.06 xlO19
-------
The theoretical gas spectra were calculated using program DEGRADE
(Deutschman and Calfee, 1967) for specific resolutions, temperatures,
concentrations, and data increments. In other applications, it is
possible to calculate new theoretical spectra for differing sets of these
parameters. For the remaining experimental spectra in TA1.1 and all of
those in TA1.2, except NO, it is only possible to generate lower
resolution spectra of different concentrations using programs LOWRES and
MORTRAN. See A1.4 for brief explanations and listings of these and
other programs used in processing the library of spectra.
In TA1.1, "increment" refers to the spacing on the wavenumber
axis between transmittance values. Wavenumber v is related to frequency
f 1
f and wavelength X by v= = y , where c is the speed of light. Note
° -1 -1
that in the lower region, 700 to 1450 cm , the increments are 0.2 cm .
In the upper region, 1800 to 3200 cm , the increments are 0.6 cm and
they are spaced evenly from 1800.0 cm" . This standardization of spacing
is necessitated by the manner in which EPAGAS compares spectra. In
the lower region, most spectra in the library were of 1.0 cm resolution,
while in the upper region, they were typically of 5.0 cm~ resolution.
Thus the 0.2 and 0.6 cm" spacings give five data points per resolution
element, providing essentially complete information about the spectra.*
*Trying to compare higher or lower resolution experimental spectra with
library spectra will introduce errors, but effects on EPAGAS outputs are
presently unknown. It should be noted that most commercial spectrometers
change resolution while making a large frequency scan. The adjustments
to the program can be made quite simply by changing the input data to
correspond to the library requirements.
73
-------
To increase the library, the literature was searched for usable
spectra of 23 molecules. These molecules have been found in polluted
atmospheres, but no high resolution spectra had been collected for use
in this application. Usable data for only 11 gases, those listed in
TA1.2, could be found with the time and resources available. For
the remaining 12 gases, either no spectra were available or published
information was unusable because no absolute scale of transmittance or
absorption was given or because the concentrations could neither be calcu-
lated nor estimated with available information. Except in a few special
cases (SO ), it was not possible to obtain spectra experimentally due to
the expense and time involved with the procurement and handling of samples
and the proper ratioing of data to obtain absolute values of transmittance.
A summary sheet for each of the 11 molecules for which useful
information was found follows this discussion of general techniques
(A1.3). Note that gas 4 is missing. This gas (S03) was deleted, after
the coding system had been finalized, when it was determined that the
spectra were in error. The important parameters required for the conversion
of the published data to a form usable by program EPAGAS are listed in seven
columns. First, a photocopy of the region of interest (column i) was made
from the literature and photo-enlarged to a size suitable for scaling on a
digitizer. Absorption values for even increments of wavenumber (column 4)
were punched onto cards, which in turn were used to produce one card for
each frequency together with its corresponding value of transmittance, plus
peripheral identifiers such as the resolution (column 2) and concentration
(column 3) found in the literature.
It was necessary to produce these cards to manipulate the spectra
with existing in-house programs. SPTR, a simple read-punch
program, was used for this purpose. As needed, the spectra were then
interpolated linearly to finer spacings (column 5) using program INTERP
and/or degraded to a lower resolution (column 6) using program LOWRES.
The latter program employs a sliding triangular slit function as a low
pass filter to smooth high-resolution data into low-resolution data.
74
-------
Table A1.2
SPECTRA PROCESSED FROM THE SCIENTIFIC LITERATURE
01
SET NUMBER
1
2
3a
3b
'3c
3d
3e
3f
5
6a
6b
6c
6d
7b
7c
7e
8a
8b
8c
GAS
NO-Nitric Oxide
H,S -Hydrogen Sulfide
C-H -Ethylene
C2H4-Ethylene
C2H4-Ethylene
CH -Ethylene
CH -Ethylene
C-H -Ethylene
N02-Nitrogen Dioxide
CH3CH2HCO-Acetone
CH3CH2HCO-Acetone
CH3CH2HCO-Acetone
CH3CH2HCO-Acetone
H2CO- Formaldehyde
H2CO- Formaldehyde
H?CO-Formaldehyde
C0H,-Ethane
i 6
C-H, -Ethane
2. t>
C~ti -Ethane
REGION
1805.4
2623.8
800.0
1380.0
1800.0
1970.4
2936.4
3050.4
2853.0
700.0
1150.0
1325.0
2900.4
851.0
1077.0
1239.0
792.6
1360.0
2863.2
- 1935.0
- 2755.2
- 1098.0
- 1450.0
- 1955.4
- 2101.8
- 3049.8
- 3200.4
- 2967.0
- 1149.8
- 1315.8
- 1450.0
- 3099.0
-- 1071.6
- 1237.0
- 1439.0
- 870.0
- 1450.0
- 3050.4
INCREHEiii"
r\ ,-
U . Q
0.6
0.2
0.2
0.6
0.6
0.6
0.6
0.6
0.2
0.2
0.2
0.6
"0.2
0.2
'0.2
0.2
0.2
0.6
NUMBER OF
POINTS PER
CONCENTRATION
217
220
1491
351
260
220
190
251
191
2250
830
626
332
1104
801
1001
388
451
313
LOWEST CONCENTRATION
(molecules cm-2}
1.06 x 1018
20
1.06 x 10
18
1.06 x 10
18
3.54 x 10
18
3.34 x 10
19
3.54 x 10
3.54 * 1018
18
1.06 x 10
17
3.54 x 10
18
1.06 x 10
17
3.54 x 10
17
3.54 x 10
18
1.06 x 10
19
1.06 x 10
3.54 x 1018
1 0
3.54 x 10
3.54 x 1018
3.54 x 1018
3.54 x 1017
-------
Table A1.2 (cont.)
ON
SET NUMBER
9a
9b
lOa
lOb
lla
lib
12a
12b
12c
GAS
n -C.H.Q-Normal Butane
n -C.H.n-Normal Butane
i -C,H]0-Iso-Butane
i -C.H.--Iso-Butane
4 10 .
CTH .-Propane
-> o
C_HQ-Propane
J o
C-H --Pentane
C H -Pentane
C5H 2 -Pentane '
WAVENUMBER
REGION
(cm-1)
700.
1800.
700.
1800.
700.
1800.
700.
1800.
2544.
0 -
0 -
0 -
o -
0 -
o -
0 -
o -
o -
1450.0
3200.4
1450.0
3200.4
1450.0
3200.4
1450.0
1999.8
3200.4
INCREMENT
(cm-1)
0
0
0
0
0
0
0
0
0
.2
.6
.2
.6
.2
.6
.2
.6
.6
NUMBER OF
POINTS PER
CONCENTRATION
3751
2335
3751
2335
3751
2335
3751
334
1095
LOWEST CONCENTRATION
(molecules cm-2)
3.
3.
1.
1.
1.
1.
3.
1.
3.
54
54
06
06
06
06
54
06
54
X
X
X
X
X
X
X
X
X
io18
18
10IH
i n
io19
19
10 *
10
io19
19
ioiy
1 8
io18
19
10 *
18
10iB
Note: Insufficient data was found for 0
2,
S03, HCN, ^2°2' C3H4°' C2H4°' C4H8°' C3H8» C7H8' and
N_-N2 collisional spectra. In most cases, the spectrum was of too low resolution or unknown concentration.
-------
Finally program MORTRAN, using Bouguer-Beerfs law, calculated
transmittance values for five concentrations from each data point for
storage in the library. Any particular problems encountered or special
techniques employed for a molecule are summarized in the final
paragraph(s) on the individual summary sheets.
In several cases, molecular bands were found that could not be
processed, either because absolute transmittance scales were not used by the
author, resolution was extremely low, or no value of concentration could be
deduced. Such bands were listed as missing if other bands of the same
molecule could be processed. It is possible, though unlikely, that bands
exist which the literature search did not find. Such bands were neither
processed nor listed as missing.
The accuracy of the spectra in the library was controlled primarily
by the accuracy of the blowup process. Typically, a 2 x 6 in. journal
spectrum was enlarged to 18 x 54 in. and data points were taken every
0.04 in. on the frequency axis of the enlargement. The uncertainty in
reproducing the position of a spectrum on the frequency axis was about
0.12 in., or -3 of the data-taking intervals (column 4). The uncertainty
in controlling the length of the transmittance axis is -5 percent of full
scale. Both of these uncertainties assume the reported spectra were
ratioed and calibrated.
The spectra for several gases at concentrations expected in the
atmosphere were nearly straight lines. That is, more than 90 percent of the
band had transmittance values above 0.97. Such spectra are nearly useless
for use in EPAGAS. Thus, the first concentration which caused at least
10 percent of the band to fall under 97 percent transmission, and the next
four higher concentrations in the standard sequence were calculated and used
in the library. The standard sequence is . . . , 1.06 x 10 , 3.54 x 10 ,
1.06 x 10n*1, 3.54 x 10n+1, 1.06 x 10n+2, . . . (molecules cm'2), where n is
an integer. The "lowest concentration" in the last column on table A1.2
refers to the first number selected from this sequence and implies the
next four which were used in the library.
77
-------
Concentrations on TA1.1 and TA1.2 and on the summary sheets are
18
given in molecules per square centimeter for a 1-km path. 1.06 x 10
molecules cm" corresponds to 3 torr partial pressure at 296 K. For
conversion to other units, consult A4 (Calfee, 1971) which lists
many of the standard units used in studying gaseous absorption processes
and summarizes their interrelationships. See the reference for detailed
explanations.
Figures l.lb and l.lc display selected spectra for each of the gases
processed. The overlays to Fl.la and Fl.lb depict band designations.
Where disagreement occurs between literature sources, the band was not
identified, but the disagreement was noted on the summary sheet. Due to
the complicated vibration and rotation modes of the larger molecules in
Fl.lc, individual bands are not easily identifiable and therefore no
overlay is provided.
Figure Al.l summarizes the state of knowledge in the infrared of
absorbing gases found in the atmosphere.
78
-------
woe »soe »ooc «soc oooc IMOO 1000 i'»c 1000 -MOO .woe MOO
Figure Al.l
INFRARED ABSORBING REGIONS OF ATMOSPHERIC GASES
-------
GAS 1: NO-Nitric Oxide
00
o
1
Frequency
Range
(cm-1)
1805.4-1935.
1. J. Chem.
2. J. Molec
2
Original
Resolution
(cm-1)
0
Phys. , Vol.
. Spec. , Vol.
3
Original
Concentration
(mol. cm~2)
24, #2, p. 399,
20, p. 11, 1966
4
Increment
Data Taken
(cm-1)
1956, Shaw.
, Abels § Shaw.
5
Converted
Increment
(cm-1)
0.6
6
Converted
Resolution
(cm-1)
3.0
7
Interpolated
Concentrations
(mol. cr.r2) Source
1 Q 0(1
1.06x10 -1.06x10 u 1,2
The NO spectra were calculated from theoretical line positions
(Shaw, 1956), theoretical line strengths (Abels 5 Shaw, 1966), and
experimental half-widths (Abels § Shaw, 1966) using program DEGRADE.
Where not reported, positions were computed using formulas presented
in Shaw, 1956, and half widths were estimated from adjacent lines.
Line parameter cards containing this unreported information are
identifiable by a missing last significant figure in the field
describing the position and/or half-width.
-------
GAS 2: H9S-Hydrogen Sulfide
1
Frequency
Range
(cm-1)
2
Original
Resolution
(cm-1)
3
Original
Concentration
(mol. cm"2)
4
Increment
Data Taken
(cm-1)
5
Converted
Increment
(cm-1)
6
Converted
Resolution
(cm-1)
7
Interpolated
Concentrations
(mol. cm- 2)
Source
2623.8-2755.2
1.2
2.7 x 10
20
0.6
1.06xl02°-1.06xl022
oo
1. Phys. Review, Vol. 37, p. 728, 1931, Nielsen § Barker
-1
A strong band, 1100-1400cm~", is missing. The resolution
claimed by Nielsen and Barker seems high. The apparent resolution
of their spectrum is close to 3.0cm-1, the resolution required for
this region. Therefore, LOWRES was not used in processing this
spectrum for the library.
-------
GAS
C,ll -Ethylcnc
oo
to
i : 5
Frequency Original Original
Range Resolution Concentration
(cm-1) (era-1) (nol. cm"2)
a
b
c
d
e
t"
800
1380
1800
1970
2936
3050
.0-1098.0
.0-1450.0
.0-1955.4
.4-2101.8
.4-3049.8
.4-3200.4
0.8
1.5
1.6
1.5
0.6
1.5
1.1 x
3.0 x
1.0 x
2.5 x
3.4 x
2.4 x
1C)
lO1'
io19
io20
io20
io19
IO19
4
Increment
Data Taken
(en-1)
0.
0.
0.
0.
0.
0.
4
2 ,
6
6
6
6
5 (-
Converted Converted
Increment Resolution
(ca-1) (cs-1)
u.: - i.
3.
5 .
5 .
5.0 5.
1.
Interpolated
Concentrations
(n:ol. c:r.--) Source-
IS
ObxlO -1
54xl01S-3
54xl018-3
54xl0ly-5
54xl018-3
06xl018-l
""0
. 54.X1020
. 54xl020
. 54xl021
.5-txlO20
.06xl020
1
2
2
2
3
-
1. J. Chem. Phys., Vol. 40, «S, p. 2096, 1964 Smith (, Mills.
2. J. Opt. Soc. Amor., Vol. 16, "3, p. 155, 1928, Levin f, Meyer.
3. J. Chem. Phys., Vol. 8, p. 799, 1940, Smith.
No accurate concentrations were reported by Levin and Meyer.
The original concentrations gix'en in column 3 were first estimated
from Sadler low resolution band strengths and then refined by
comparison with neighboring spectra where possible.
For instance, the original spectra in the 2936-3049 cm"1 and
3050-3200cm"l regions were assumed to have been taken at concen-
trations of 2.7 x IO1 mol-cm and 3.0 x IO19 mol-cm2, respectively
When converted to a concentration of 5.54 x IO18 mol-cm'2, a
noticeable mismatch occurred between the two spectra at 3050cm"1.
To cause the two spectra to coincide at 3050cm"1, the assumed
concentration for the first spectrum had to be increased to
3.4 x IO19 mol-cm~2 (25.3% increase) ai'.d the second lowered to
2.4 x IO19 mol-cm"2 (19.1% decrease). These refined figures were
the ones reported on the summary sheet, and they reflect a 50uo
uncertainty in the reported concentration.
Once again, the reported resolution seems too high i if all
the Levin and Meyer spectra. Apparent resolution is about 3.0cm"1
Note that the apparent resolution is satisfactory for the
1800-3200cm"1 region but too low for the TOO-lSOOcnT1 recion.
-------
GAS 5: NO--Nitrogen Dioxide
1
Frequency
Range
(cm-I)
2853.0-2967.0
2
Original
Resolution
(cm-1)
0.6
3
Original
Concentration
(mol . cm~^)
1.1 x 1019
4
Increment .
Data Taken
(cm-1)
0.6
5
Converted
Increment
(cm-1)
6
Converted
Resolution
(cm-1)
1.0
7
Interpolated
Concentrations
(mol. cm- 2)
3.54xl017-3.54xl019
Source
1
oo
1. Private communication from Goldman, University of Denver, 1969.
Only about one half of the band is available from the band
center to one wing. Another band, centered at 1617cm~* is
available, but was omitted from the library due to strong f^O
absorption in this region.
-------
GAS 6:
oo
CH3CH2HCO-Acetone
a
b
c
d
1
Frequency
Range I
(cm-1)
700.0-1149.8
1150.0-1315.8
1325.0-1450.0
2900.4-3099.0
2
Original
Resolution
(cm-1)
(IR-7)
(IR-7)
_7
(IR 7)
(IR-7)
3
Original
Concentration
(raol. cm"2)
8.2 x 10
3.2 x 10
18
9.9 x 10
1 0
9.9 x 10
4
Increment
Data Taken
(cm-1)
1.0
1.0
1.0
1.0
5
Converted
Increment
(cm-1)
0.2
0.2
0.2
0.6
6 7
Converted Interpolated
Resolution Concentrations
(cm-1) (mol. cm-2)
1.06xl018-1.06xl020
3.54xl017-3.54xl019
17 1Q
3.54x10 -3.54x10 y
1.06xl018-1.06xl020
Source
1
1
1
1
1. J. Chem. Phys., Vol. 44, #1, p. 97, 1966, Cossee 6 Schachtschneider.
A band exists from 1680-1825cm as reported by Cassee and
Schachtschneider, but it was omitted due to strong H20 absorption
in this region.
The resolution is lower than desired (approximately 7cm )
due to the low resolution of the Beckman IR-7 spectrometer. Since
a range of resolutions may be preset on this instrument, one needs
to consult the workers to determine the exact resolution of the
various runs.
-------
oo
tn
GAS 7: H2CO-Formaldehyde
1 2 3
Frequency Original Original
Range Resolution Concentration
(cm-l) (cm-1) (mol. cm-2)
a
b
c
d
e
851.
997.
1077.
1239.
1321.
0- 997.2
4-1071.6
0-1237.0
0-1325.0
0-1439.0
0.3
0.3
0.3
0.3
0.3
7.1 x
3.5 x
4.7 x
2.4 x
5.3 x
1 Q
ioiy
1 Q
1019
1Q
ioiy
10
ioiy
1Q
ioiy
4 5
Increment Converted
Data Taken Increment
(cm-1) (cm"1)
0.
0.
0.
0.
0.
2
2
2
2
2
6 7
Converted Interpolated
Resolution Concentrations
(cm-1) (mol. cm- 2) Source
1.
1.
1.
1.
1.
0
0
0
0
0
1Q 71
1.06x10-1.06x10
iq ? I
1.06x10 -1. 06x10
18 20
3.54x10-3.54x10
ic 20
3. 54X101 -3.54x10
18 20
3.54x10-3.54x10
1
1
1
1
1
1. J. Molec. Spec., Vol. 38, p.. 90, 1971, Nakagawa § Morino.
Temperatures of the original spectra range from 150°C
to 25°C, accounting for small discrepancies which may occur
between regions. See Nakagawa § Morino for details.
Note an error on page 90, Figure le. The reported cell
length is 100cm but should read 10cm.
Nakagawa and Morino report the band head at 1167cm~
to be v., while Herzberg claims it is v
6*
-------
00
GAS 8: C,H, -Ethane
/ o
1234
Frequency Original Original Increment
Range Resolution Concentration Data Taken
(cm-1) (cm-1) (mol. cm'2) (cm-1)
1Q
a 792.6- 870.0 0.6 6.0 x 10 0.2
1Q
b 1360.0-1450.0 1.3 4.0 x 10iy 0.2
1Q
c 2863.2-3050.4 1.4 6.0 x 10 0.6
56 7
Converted Converted Interpolated
Increment Resolution Concentrations
(cm"1) (cnr1) (mol. cm~2) Source
in 20
3.54x10-3.54x10 1
18 20
3.54x10-3.54x10 1
17 1Q
3.54x10-3.54x10 1
1. J. Opt. Soc. Amer., Vol. 16, #3, p. 155, 1928, Levin 6 Meyer.
No accurate concentrations were reported. The original concen-
trations given in column 2 were estimated from Sadler low resolution
band strengths.
The reported resolution seems too high. The apparent resolution
is about 1.0cm"1 in the 793-870cm"1 region and about 3.0cm"1 in the
remaining two regions.
-------
00
a
b
1
Frequency
Range
(crn-I)
700.0-1450.0
1800.0-3200.4
2
Original
Resolution
(cm-1)
3
Original
Concentration
(mol. cm~2)
2.5 x 102°
2.5 x 102°
GAS 9:
4
Increment
Data Taken
(cm-1)
5.0
5.0, 10.0
C4H1Q-N-Butane
5
Converted
Increment
(cm-1)
0.2
0.6
6
Converted
Resolution
Ccm-1) -
7
Interpolated
Concentrations
(mol. err 2)
3.54xl018-3.54xl02°
3.54xl018-3.54xl020
Source
1
1-
1. Anal. Chem., Vol. 28, p. 1218, 1956, Pierson, Fletcher, Grantz.
The resolution for all regions is unknown, but is
approximately 10-15cm~l. The instrument used was a Perkin-Elmer
Model 21.
-------
GAS 10:
00
00
C4H1ft-I-Butane
a
b
1
Frequency
Range
(cm-1)
700.0-1450.0
1800.0-3200.4
2
Original
Resolution
(cm-1)
.
3
Original
Concentration
(mol. cm" 2)
2.5 x 10
2.5 x 1020
4
Increment
Data Taken
(cm-1)
5.0
5.0, 10.0
5
Converted
Increment
(cm-1)
0.2.
0.6
6
Converted
Resolution
(cm-1)
7
Interpolated
Concentrations
(mol. cm" 2)
1.06xl019-1.06xl021
1.06xl019-1.06xl021
Source
1
1
1. Anal. Chem., Vol. 28, p. 1218, 1956, Pierson, Fletcher, Grantz.
The resolution for all regions is unknown, but is
approximately lO-lScnr1. The instrument used was a Perkin-
Elmer Model 21.
-------
00
GAS 11: CjHg-Propane
a
b
1 2
Frequency Original
Range Resolution
(cm-1) (cm-1)
700.0-1450.0
1800.0-3200.4
3
Original
Concentration
(mol.
2.5
2.5
cm-2)
xlO20
xlO20
4
Increment
Data Taken
(cm-1)
5.0
5.0, 10.0
5
Converted
Increment
(cm-1)
0.2
0.6
6 7
Converted Interpolated
Resolution Concentrations
(cm-1) (mol.
1.06xl019
1.06xl019
cm-2) Source
-1.06xl021 1
-1.06xl021 L
1. Anal. Chem., Vol. 28, p. 1218, 1956, Pierson, Fletcher, Grantz.
The resolution for all regions is unknown, but is
approximately 10-lScm"1. The instrument used was a Perkin-
Elmer Model 21.
-------
GAS 12;
to
o
C5H12-Pentane
a
b
c
1
Frequency
Range
(cm-1)
700.0-U50.0
1800.0-1999.8
2544.0-3200.4
2
Original
Resolution
(cm-1)
-
3
Original
Concentration
(mol. era"2)
9.1 x 10 9
19
9.1 X 10
9.1 x 1019
4
Increment
Data Taken
(on-1)
5.0
5.0
10.0
5
Converted
Increment
(cm-1)
0.2
0.2
0.6
6
Converted
Resolution
(cm-1)
-
7
Interpolated
Concentrations
(raol. en" 2)
3.54xl018-3.54xl020
iq 21
1.06x10-1.06x10
3.54xl018-3.54xl02°
Source
1
1
1
1. U.S. Naval Research Lab., American Petroleum Institute Research Project 44.
The original spectrum for this molecule was for a liquid
sample 0.174mm thick, and therefore it may differ significantly
from the gaseous spectrum. Equivalent gaseous concentrations
were figured from the density of liquid pentane and the thickness
of the sample.
The resolution for all regions is unknown, but is approximately
10-15cm'1.
-------
A1.3 PROGRAMS USED IN PROCESSING SPECTRA
FORTRAN PROGRAM
CDC 3800
SPTR
M. J. Post
Aug. 17, 1973
Given absorption values A expressed as an integer number between
999 and 0, corresponding to 99.9 and 0.0 percent absorption, SPTR punches one
card for each absorption value. On this card are punched wavenumber and
the value of transmittance.
Card A. - A
--- i
Column
1-80
Description
Absorption values A(I)
Output
Card A.
1-10
11-20
54-58
59-70
73-80
F10.2
F10.3
F5.2
E12.3
A8
Wavenumber V.
Transmittance T.
Resolution RES
Concentration CONG
Molecule name MOL
Note; 1. An EOF card is required at the end of the data cards.
2. RES, CONC, and MOL are entered on source cards, no^, data cards
3. V. are generated internally by changing source cards.
91
-------
PROGRAM SPTR
THIS PROGRAM TAKES COMPACTED SPECTRAL ABSORPTION VALUFS (FORMAT 2014)
AND OUTPUTS CARDS. ONE FREQUENCY AND ONE TRANSMISSION VALUE PER.
CARD, PLUS PERIPHERAL IDENTIFIERS.
IT IS USEFUL IN CONVERTING DIGITIZED ANALOG SPECTRAL INFORMATION
(FROM STRIP CHARTS, FOR EXAMPLE) TO CARDS IN A FORMAT USED IN
ASSOCIATED PROGRAMS.
PPOG3AM SPTR
DIMENSION A(20>
i FORMAT <20F4.4>
2 FORMAT
-------
FORTRAN PROGRAM
CDC 3800
LOWRES
M. J. Post
Aug. 17, 1973
Given a transmission spectrum T(N) of high resolution, LOWRES
degrades it to a spectrum TF of lower resolution A.
Used in the computations are the given increment DV and the output
increment DELV. A triangular slit function of half-width A is used
to compute each new TF.
Input
Card A
Column
1-10
11-20
21-25
26-30
31-35
Card B
11-20
F10.3
F5.2
F5.2
F5.2
F10.3
Description
Beginning wavenumber
of given spectrum, VI
End wavenumber of
given spectrum, V2
Increment given
spectrum, DV
Lower resolution
desired, A
Desired output
increment, DELV
Given spectrum
93
-------
Output
LOWRES - (cont.)
Card A.
Column
1-10
11-20
54-58
50-70
73-80
F10.5
F5.2
E12.3
A8
Description
Wavenumber of output
transmittance, VI.
Transmittance TF. of
output spectrum
Resolution of output
spectrum RES (=A)
Concentration of output
spectrum CONG
Molecular name MOL
NOTE: An EOF card is required at the end of the data cards. DELV must
be an exact multiple of DV and A must be an odd multiple of DV.
CONG and MOL are entered on source cards, not data cards.
94
-------
PROGRAM LOWRES
THIS PROGRAM TAKES A SPECTRUM OF HIGHER RESOLUTION GIVEN IN INCREMENTS
OF DV AND DEGRADES IT TO A SPECTRUM OF RESOLUTION A, GIVEN IN INCREMENTS
OF DELV, POSITIONED AT THE ORIGINAL FREQUENCIES.
DELV MUST BE AN EXACT MULTIPLE OF DV.
A MUST BE AN ODD MULTIPLE OF DV.
DIMENSION T(2000)
1 FORMAT (2F10.3.3F5.0)
2 FORMAT (10X,F10.3)
3 FORMAT (F10.3,F10.5,33X,F5.2,E12.3,2X,A8)
RES=3.0
CONC=1.05E+19
MOL=8H N02
READ 1,V1,V2,DV,A,DELV
MULT=DELV/DV + .00001
SLIT=A*2.
NMX=SLIT/DV-.99999
SLTFTR=DV/A**2
VEND=V2-A+DV
N-l
20 READ 2,T(N)
IF (EOF,60) 30,25
25 N=N+1
GO TO 20
30 NJ=1
VI=V1+A-DV
V=VI+DV
50 SUMT=0.
DO 150 NN=NJ,NMX
TT=(A-ABS(V-VI))* T(NN)
SUMT=SUMT+TT
150 V=V+DV
TF=SUMT*SLTFTR
PRINT 3,VI,TF,RES,CONC,MOL
PUNCH 3,VI,TF,RES,CONC,MOL
NJ=NJ+MULT
NMX=NMX+MULT
VI=VI+DV*MULT
IF (VI. GT. VEND) CALL' EXIT
V=VI-A+DV '
GO TO 50
END
95
-------
FORTRAN PROGRAM INTERP
CDC 3800 M. J. Post
Aug. 17, 1973
Given a spectrum spaced in increments of DELV, INTERP linearly
interpolates it to finer increments DIV and outputs it on cards.
Input
Card A
Column Field Description
1-5 F5.2 Input spectral
increment DELV
6-10 F5.2 Output spectral
increment DIV
Card B
1-10 F10.3 Beginning wavenumber
of spectrum VI
11-20 F10.5 First transmittance
value of spectrum Tl
54-58 F5.3 Resolution of spectrum
RES
59-70 E12.3 Concentration of
spectrum CONG
73-80 A8 Molecular name MOL
Card C
11-20 F10.5 Transmittance values
T2.
96
-------
Output
Card A.
INTERP - (cont.)
Column
1-10
11-20
54-58
59-70
73-80
F10.5
F5.3
E12.3
A8
Description
Wavenumber of
transmittance value,
VI.
Transmittance value T.
Resolution of spectrum
RES
Concentration of
spectrum CONC
Molecular name MOL
97
-------
C T-JIS PROGRAM TAKES A SPECTRUM SPACED I ^ INCREMENTS OF DELV AND LINEARLY
c INTERPOLATES IT TO FINFR INCREMENTS DIY.
c
C DFLV MUST BE AS INTEGRAL MULTIPLE OF DIV.
C
1 FORMAT (2F5.2)
? FORMAT UOX.FIO.S)
3 FORMAT (F10.3.P10.t;»33X»fS.3,pi2.3,2x,Afi)
RFfiD l.DFLV.OIV
RF6D 3«V1»T1»RES«CC)NC»MOL
N=HELV/OIV * .001
VI=V1-DIV
10 READ 2,12
IF (rOF«f.O) ^0»12
12 PART =(T2-T1)/>J
00 15 K=1»N
T=T1«PART*
98
-------
For a detailed discussion of program DEGRADE, see Deutschman and
Calfee, 1967.
99
-------
OEOPADE
«Mi)«F(5) , AR<5)«TT<5)«TF<5),SUMT(5),TRAP$UM<5> tFMTH9)t
«UJ<9) »FHT4<9)
COMMON/ 1/uNU (3000 >.SA( 3000 ).AA ( 3000) < CA Y ( 1 > /?/T
t» FnPMAT4H K1CHDELV =F7.^»3X7HBOUNO «F7
9 FORMAT (///33HODATA FOR PLOTTING ARE ON TAPF c-
10 FOHMAT(I«»)
11 Fn»MAT (60H1TOO MANY LINE CARDS«EXCEEDS DIMENSION. LAST CARD RF.AD
IAS =F9.?)
12 Fn«MAT(19HH? TOO LARGE. 17 «I5»SX9HGNJ ( 17) =F9.2.SX12HV1-A*BOUND
13 RFAD i«(FMTKi>*i-i*9)
14 3FAD l.(FMT2«I=1.K1)
C
C SFT UP TEMPERATURE CORRECTION CONSTANTS.
C
1810 C51 = (TEMPO-TEMP)/(TEMPO»TFMP».6950)
1820 C2= (TEMPO/TEM») *»RX
1830 CA=< (TEMPO/TEMP)»«CX)»P
CftLL 09 EXJN
C RfAD IN Ll^E DATA AND COJNT CARDS
C
VLO=VI-A-BOUND-I .$ \/Hi=v?»A*8ouso
20 1=0
^i 1=1*1
C
C
25 RFAD FMTI *GV tS»ALF»F3P«MOL
IF(GM .I.T.VLO) GO TO 25
lr?7.27?
27? 1F(MOL.F3.1)GO 70 ?74
273 A| F = .l
274 AA(I)=(ALF*CA)*«2
G\.lUI)=GN
2f>01 S«(I)=S»ALF»CA»CS2»EXP(-?PP»CSU $ GO TO 21
27 11=1-1
28 FORMAT (2F10.J. 15)
2«« GNJ(I)tGMtMIl)»Il
100
-------
C CHECK THAT G^U
3? XPOUND = VI -A BOUND
33 lF(GNU(I7i> .<3T. XBOUND)34»36
34 P9INT 12.17,&MJH7) .XROUND
35 GO TO 1001
3ft PPINT 7.PtTEMP,A,DV,DEl.V,ROUND
37 P°TNT 3R.(ri(K),K=l,Kl)
38 Fr>RMAT<26X,*rfl =»E10.3«5X»»W2 *»E10.3«5X»«3 =»F10.3«5X *W4 =»r!0.3
l.SX«rt5 =»E10.3/)
39 PPINT 5
C
C SFT UP CALCULATION COUNTERS.
C
40 IS=1
41 SLIT=2.«A
46 lniV=DELV/DV».01
47 NJn001 = iniV»l
48 Nonoa = (SLIT/DV) * 1.00001
49 MOOD3 = NJODD2 - IDIV
53 *J = l
bs SLTFTP, = ov/JD APPLICATION or VARIABLE PORTION OF SLIT
c FACTOR INTEGRAL
c
130 DO 133 K=l,Kl
134 T(MN.K) = EXPF(-rf(ff)»CAY)
136 T7 509»164tl66
164 TPAPSUM«) = T=JAPSUM(K) »0.5«E(K)
165 GO TO 174
101
-------
TOAPSUMdO = T3APSUM(K) « EITE(2. 6 )V»(TF(K),K=1»K1>
1B1 V=V»DELV
190 On 191 K=1,K1
191 5UMT(K)aO.O
200 IF(V-V2>201,201.235
f*
C SAVING OF AL^EADIT CALCULATED TRANSMISSION VALUES Fn» NEXT STEP
C IM THE CALCULATION
C
201 DO 205 K=l,Kl
202 DO 204 NN = NOOD1»»400D2
203 KK- = NN-If)IV
20* T(KK.K)=T(NN«K)
205 CONTINUE
210 vr=v
211 V=V-A
21S 00 219 NN = 1.SODD3
216 DO 218 K=l.Kl
217 TT(K) = (A-AiU<;F(^-\/m»T(N^.K)
21B S(|MT(K)=SUMT(K)»TT(K)
219 V=V»DV
225 NJ = NOD03 « 1
229 GO TO 79
C
c TPAPAZOIDAL «ULE CALCULATION OF TOTAL ABSORPTANCE
c
235 BPINT 3»P»A.TFMP
236 Do 238 K=l,Kl
237 AR(K) = DELV»(TRAPSUM(K)-0.5»E(K»)
238 POINT R*AB(K)*K(K)
239 *FAD 240.V1.V2
240 Fr,»MAT(2Fl0.2)
241 1F 250.^.P
250 GO TO(251»10tl).\P
251 Ei>-n FILE 2
252 3F*IND ?
2,55 GO TO 1001
C
C -UCKOUT STATEMENTS
C
500 N* = 31
501 GO TO 1000
509 \v = Iftl
1000 POINT 4,MA,I1,VUU1)
1001 CALL EXIT
102
-------
SURROUTIML A8SCOEF (V.BOU>JD« li
COMMON/1 /GNU ( 3000 ) «SA (3000 ) .AA (3000) .CAY ( 1 )/2/T 13000.5) . 15
C
C DETERMINATION OF INDEXING VALUES FQR ABSORPTION COEFFICIENT
C
10 00 14 I =15*11
11 lF(V-BOUNl>-GNU12»12»14
1? TS = I
13 GO TO 15
14 CONTINUE
15 00 19 K = IStll
1ft IF
5? RFTURN
END
SCOPE
ILDAD
2,10000
(I2.I4,2F6.3,2FR.«;.F5.2.E10.3.I2.F7.2)
103
-------
SUBROUTINE ABSCOEF(V,BOUND,II,P)
COMMON/1/GNU(3000),SA<3000),AA (3000),CAY(l)/2/T(3000,5),15
c
c DETERMINATION OF INDEXING VALUES FOR ABSORPTION COEFFICIENT CALC.
c
10 DO 14 I =15.11
11 IF(V-BOUND-GNU(I))12,12,14
12 15 = I
13 GO TO 15
14 CONTINUE
15 DO 19 K = 15,11
16 IF(V+BOUND-GNU(K))17,17,19
17 16 = K
18 GO TO 25
19 CONTINUE
c
c CALCULATION OF ABSORPTION COEFFICIENT WITH LORENTZ LINE SHAPE
c
25 CAY1 = CAY2 =0.0
30 DO 46 I = 15,16
32 Y = ABSF(V-GNU(I))
34 IF(Y-2.)36,36,42
36 SUM1=SA(I)/(Y**2+AA(I))
38 CAY1 = CAY1 + SUM!
40 GO TO 46
42 SUM2=SA(I)/Y**2
44 CAY2 = CAY2 + SUM2
46 CONTINUE
50 CAY = 0.3183*(CAY1+CAY2)
52 RETURN
END
SINGLE-BANK COMPILATION.
104
-------
FORTRAN PROGRAM
CDC 3800
MORTRAN
Dana Gregory
June 29, 1973
Given the original gas concentration, W , the corresponding
theoretical transmission, T(W , v.). as a function of the wavenumber v.,
and the original gas concentration, W , MORTRAN will compute theoretical
transmissions, Tl(W.,v.), for the desired concentrations W. .
Used in the computation is Beer's Formula,
.,.
where TlfW.v.) is a function of the set of transmission values for each
i J
desired concentration.
Card A
Card B
Column
1-10
11-20
21-30
31-40
46-50
1-10
11-20
15
E10.3
E10.3
Description
Beginning wavenumber
Ending wavenumber
Wavenumber increment
Original concentration,
Wo
NW, number of desired
concentrations to be
used in signal run
1st desired
concentration, W.
2 desired
concentration, W_
105
-------
Column
Field
MORTRAN - (cont.)
Description
71-80
Card C.-C
i n
1-10
11-20
E10.3
F10.0
F10.0
i desired
concentration, W.
Wavenumber, v.
Transmission,
T(Wo,v.)
Output
Card A.
1-10
11-20
21-30
31-40
11-15
F10.2
F10.2
F10.2
E10.3
Card Bi-Bn
1-5
6-10
15
3PF5.0
3PF5.0
Beginning wavenumber
Ending wavenumber
Wavenumber increment
Concentration, W.
i = 0; i = 1, NW
Sequence number
Transmission,
T(Wi,vj)
Transmission,
76-80
3PF5.0
Transmission,
106
-------
MORTRAN - (cont.)
SUBROUTINES USED: None.
MORTRAN: Computes additional transmissions and outputs in packed format,
107
-------
01
9.?
j i
3«+
05
06
065
07
OS
09
21
22
c-l ON <!?' (1 i ) , T (51 1, 01 ,
Get (U/, F ',"V,W,n;j
IF (ECf ,6C>
( yr
") 0 L f I = 1 , f'T
=?eflr: Of , v(T» ,T (u ,V.nv)'';:»,OR
X=V(I-1)+CV
Y=V(I)
I P ( A PS ( Y- X » . I f..CJ01)C9tM1
CO^TTNUF
IF (NCtcOS*lc.xllr.NT)?l ,22
NCAftCS
13=1
PUNCH
JO
PRINT 30J2,K,
CONTINUE
00 WO J=liNWB
EXPT=WN(J)/W
00 ID 1=1, NT
Tl (I)=T(I)**Py;'T
PUNCH 20, ev^v/.
j>
PRINT 3031, P7,EV,OV,WN
10
20
3001
00 i»0 Ksl
PUNCH 30» Kf 2tK,(Ti(u ,L=IG,IE:)
3002 FORHAT«2X,If ,(15 (3PF5.0J ) )
I3=IE*15
IE=IE+15
J.O
850
851
) »Tl (3LOO)
CONTINUE
GO TO Cl
PRINT 851
FQPMAT(12) ,''1X,*SOC'PY, THt«?E ARF NO MORF CARDS')
108
-------
I MOPTRflN
CALL EXIT
900 PRINT 901, NT,'lV,fV,riV
901 FOPMAT(20X»Fc-ny IN C Af?CS*//7X« NT = * ,Tt , *rV=*,F13 . 2 , ?X*EV=* , F13 . 2
5 3X,*DV=*,F1P.'>
CALL EXIT
910 *
-------
APPENDIX 2
TABLE OF CONTENTS
A2. Program EPAGAS Documentation and Analysis 111
A2.1 Program Summary 111
A2.1.1 List of Fortran Symbols and Descriptions 112
A2.1.2 List of Program and Subprogram Descriptions 121
EPAGAS CALC
CLOSER CONST
CXMIZ EROR1
EROR2 FINDIT
FINDTRN GAUSSL
GRAPHS IBINSER
MINMYZD PARTIN
PLOTD PLOTT
PRINT1 PRINT2
PRINTS READALL
READ1 READ2
READ3 READ4
READ6 READ?
READS RMSERR
THEOTRN
A2.2 Simulation of Gaussian Noise 170
A2.3 Constituent and Transmission Error Formulation 171
A2.4 Interpolation Formulas Used in Transmission Library 173
A2.5 Program Execution Techniques 177
A2.5.1 Tbeoretical Gas Library Spectra 177
A2.5.2 Experimental Transmission Data 178
A2.5.3 Minimization of "Goodness of Fit" Controls 179
A2.5.4 Error Parameters 180
A2.5.5 Miscellaneous Execut. Techniques 182
110
-------
APPENDIX 2
A2. PROGRAM EPAGAS DOCUMENTATION ANALYSIS
A2.1 PROGRAM SUMMARY
In this appendix, the main minimization program EPAGAS and all of the
subprograms needed by EPAGAS are described. All coding was done in the
FORTRAN applicable to the CDC 3800 computer. Each description lists
what that particular subprogram does, what input parameters are needed,
and what output parameters are produced.
The computer program must be given: (1) a library of theoretical
transmissions for each gas constituent which might possibly be present
in the measured total transmissions. For each constituent, the theoretical
transmissions are known functions of both wavenumber and gas
concentration. The library must be large enough to bound the regions of
interest in both wavenumber and concentrations to be expected. See
THEOTRN description for exact preparation of the library table. Due
to possible computer memory limitation, the user may wish to "tailor-
make" separate libraries for different experimental situations. However,
care must be taken in selecting proper wavenumber spacings and individual
concentrations so that double linear interpolation within the transmission
table will result in correct transmissions values for wavenumbers and
concentrations lying between grid points of the library table.
(2) Next, a set of measured experimental total transmissions TG
as a function only of known wavenumber v must be provided. This spectrum,
obtained from the spectrophotometer, will be compared to a trial spectrum
TC which is a function of gas constituent, gas concentration W, and
wavenumber v as determined from the library of theoretical transmissions.
By adjusting the amounts of W for each gas constituent as a function of wave-
number, the program tries to get TC to agree as best as possible with TG,
thus determining the gas composition and the amounts of each gas
constituent that were present in the sample producing the total
transmission, TG. READ7 and/or READS descriptions provide specific
information on the preparation details for inputting TG.
Ill
-------
(3) The program must also be provided with a starting set of gas
concentrations W for each gas constituent that is suspected to be present
in the sample producing the given experimental total transmission TG.
Also given is a corresponding set of concentration increments, WI, which
are used when the minimization process tries to vary the individual W.
by an amount WI. to find the best agreement between TG and TC. See the
READ6 description for preparation of W and WI.
(4) Along with the initial values of W and WI, the program also
needs error bounds for each of the individual gas constituent
concentrations. See EROR1 description for detail error formulas. The
errors defined on the gas concentrations are used primarily to help
guarantee an acceptable solution by keeping the gas concentrations
within physically reasonable bounds. READ2 description provides details
for inputting the error parameters.
(5) Finally, various program controls must be provided. Some of
these control parameters are used in directing the minimization process
and are listed in detail in the READS description. The remaining control
parameters are used in directing printing and plotting of intermediate
results and in providing certain convergence criterion. Details for these
parameters are found in the READ4 description.
A2.1.1 List of Fortran Symbols and Descriptions
Single-Valued Variables
Variable Description
CFM A constant confidence factor to be applied to the
experimental transmission data if individual ones
are not specified.
GLB The minimum root-mean-square (rms) error permitted
before using the minimization routine MINMYZD again.
112
-------
Variable
IPLOT
IPRINT
MK
N
NG
NT
Q
RMS
RMSTP
RMST
RMSW
Description
A plotting control which, when set to m, permits a
graph of the given and currently computed
transmissions, as functions of wavenumber, to be
plotted after every m use of MINMYZD.
A printing control which, when set to n, permits
printing of intermediate results after every n
use of MINMYZD.
The maximum number of uses of MINMYZD permitted for
a particular data set.
N = NG + NT. The total number of errors in the
error array ERR.
The total number of gas constituents to be found
in the library of theoretical transmissions.
The number of experimental transmission values used.
The minimum percent by which each try of MINMYZD
must improve the mean square deviation.
Current root-mean-square (rms) error computed
from the error array ERR.
Lowest root-mean-square (rms) error computed
from the error array ERR.
Lowest root-mean-square (rms) error computed
from the NG+1 through NG+NT elements of the error
array ERR. This rms value corresponds to
the errors based on the difference between the
given experimental transmission values and the
theoretical transmission values currently computed.
Lowest root-mean-square (rms) error computed
from the 1 through NG elements of the error
array ERR. This rms value corresponds to the
errors defined on the individual gas concentrations.
113
-------
Variably
UB
UE
UI
UMAX
UMIN
ZEP
Description
Beginning wavenumber corresponding to the given
transmission values if they are equally spaced in
wavenumber.
End wavenumber corresponding to the given
transmission values if they are equally spaced
in wavenumber.
Wavenumber increment corresponding to the given
transmission values if they are equally spaced in
wavenumber.
Maximum wavenumber value for which the given
transmission values will range over.
Minimum wavenumber value for which the given
transmission values will range over.
A number used by MINMYZD which tells how much of the
predicted final step to take.
Array Variables
Variable
CF
ERR
Maximum
dimension
150
160
Description
Confidence factor, O-CF(K), K=l, NT,
corresponding to the NT given
transmission values.
Errors based on the individual gas
concentrations and the weighted
differences between the given and
computed transmission values.
Elements 1 through NG are the errors
corresponding to the gas concentrations,
and elements NG+1 through NG+NT are
the errors corresponding to the
transmission errors.
114
-------
Variable
IDGAS
Maximum
dimension
10
IPW
10
KX
10
Description
Integer numerical identifiers (12)
for each of the gases to be found in
the gas library of theoretical
transmissions. IDGAS(K), K=l, NG
correspond to the K gas.
One of the error parameters used in
computing the errors based on the
individual gas concentrations. This
parameter is an odd positive integer
used as a power in the error
function. IPW controls the rate of
increase of the error as the current
gas concentration ventures outside
the acceptable range for that
particular concentration. IPW(K),
K=l, NG correspond to the K* gas.
The first six elements of the array KX
are input control parameters used by
the main minimization routine MINMYZD.
The last four elements of the array are
output parameters of MINMYZD and
provide various convergence information.
KX(1) is the maximum number of new sets
of errors (minimization cycles) to use.
If KX(1)=0, then one cycle is used. If
KX(1)<0, then each gas concentration
increment WI is reduced by 10 percent after
each complete use of MINMYZD, and
|KX(1)| is used for the maximum number
of cycles to be used in one complete
use of MINMYZD.
115
-------
Variable
Maximum
dimension
Description
KX(2) is the number of cycles to count
for a bad cycle, that is, one which did
not result in enough improvement,
KX(2)-1.
KX(3) is the percent improvement in
tenths in the rms error for the
cycle to be counted as a good cycle.
Any cycle resulting in less than this
causes KX(2) cycles to be counted
rather than one cycle, KX(3)-0.
KX(4) is the limit as to the multiple
of a WI increment which may be taken
in one full step (ZEP=1). If KX(4)=0,
then 100 is used. If KX<0, then no
step will be permitted in the -WI
direction for the first complete use
of MINMYZD. Future uses of MINMYZD
uses [lCX(4)| .
KX(5) is the minimum ZEP in 0.001's
permitted. KX>0.
KX(6) is the minimum total percent
improvement in tenths for ZEP not to
be changed. KX(6)>0.
KX(7): If KX(7)=1, then at least one
single step (even if it were a step
of WI.) resulted in an improvement of
at least KX(3); otherwise1, KX(7)=0.
KX(8): If KX(8)=1, then there was some
improvement sometime; otherwise,
KX(8)=0.
116
-------
Maximum
Variable dimension Description
KX(9): If KX(9)=-1, then all WI's
were too small. Convergence has
probably been reached. If KX(9)=0,
then the solution matrix was
singular. If KX(9)=1, then the matrix
was normal.
KX(10) is the actual total percent
improvement in tenths.
NU 10 Number of wavenumber values used for
each gas in the library of theoretical
transmissions. NU(K) is the number
of wavenumber values used for the
Kth gas, K=l, NG.
NW 10 Number of concentration values used
for each gas in the library of
theoretical transmissions. NW(K) is
the number of concentrations used for
the Kthgas, K»l, NG.
S 5 An array used by READALL and its
associated reading routines. If S(K)=0,
then the input information was correct.
If S(K)j*0, then a reading error occurred.
TC 150 Computed theoretical transmission
values, TC(K), K=l, NT based on the
NG current gas concentrations being
used. TC(K) corresponds to the K
value. TC(K) = j{G TCI(I,K).
TCI 10 Individual compu- ~ ted theoretical
transmission values TCI(I,K) for a
given wavenumber U(K), K=l, NT and
given gas concentrations W(I), 1=1, NG.
117
-------
Variable
TG
TRAN
Maximum
dimension
150
(5000,5)
U
150
UR
(10,3)
WHR
10
Description
Given experimental transmission values
TG(K), K=l, NT as a function of
known wavenumber U(K), K=l, NT.
A two-dimensional array which contains
all of the theoretical transmission
values that make up the gas library.
The transmission values are a
function of known wavenumbers,
known gas constituents, and known
concentrations for each of the
specified gas constituents. See
THEOTRN description for exact
formulation of the array.
Known wavenumbers, U(K), K=l, NT
corresponding to the array TG of the
NT given experimental transmission
values.
A two-dimensional array UR(I,J)
containing the wavenumber range and
spacing used for each of the gas
constituents in the gas library.
I corresponds to the I gas. J=l
corresponds to the beginning wave-
number, J-2 corresponds to the end
wavenumber, and J=3 corresponds to
the wavenumber increment.
One of the error parameters used in
computing the errors based on the
individual gas concentrations.
This parameter is the half-range of
permitted variation in concentration
for the K gas constituent and is
118
-------
Maximum
Variable dimension Description
simply the difference between the
mean concentration WM(K) and the
lowest concentration of the K gas
constituent.
WIO 10 Original increments to be applied to
the starting gas concentration WO(K),
K-=l, NG. These are needed by MINMYZD
when it tries to adjust the original
concentrations WO to get the best
agreement between the given and
computed transmission values. These
values may be modified in the
minimization process, and hence
become WI(K), K=l, NG.
WI 10 Current incremental values of gas
concentrations, WI(K), K-l, NG.
WM 10 One of the error parameters used in
computing the errors based on the
individual gas concentrations. This
parameter is the mean concentration
value of the desired upper and lower
bounds of the concentrations of the
Kth gas constituent, K=l, NG.
WO 10 Original gas concentrations to be used
in the minimization process. WO(K)
correspoi
K=l, NG.
corresponds to the K gas constituent,
119
-------
Maximum
Variable dimension Description
WR (10,5) A two-dimensional array WR(I,J)
containing the actual gas
concentrations used in the gas
library. I corresponds to the I
gas constituent, and J corresponds to
the J gas concentration of the I
gas constituent.
W 10 Current gas concentrations being used
in the minimization process. W(K),
K=l, NG corresponds to the K gas
constituent.
WW 10 One of the error parameters used in
computing the errors based on the
individual gas concentrations. This
parameter is a "normalization"
factor and is used to scale the gas
concentration errors. WW(K), K=l,
NG are the factor for the K gas
constituent.
120
-------
A.2.1.2 List of Program and Subprogram Descriptions
FORTRAN PROGRAM EPAGAS
CDC 3800 Margot Ackley
July 30, 1973
Program EPAGAS is a computer program which is designed to determine
atmospheric composition from measured total transmissions at known wave-
numbers. The program utilizes a set of theoretical transmission spectra
of selected gas constituents at known values of concentrations and at known
wavenumbers, and a set of observed total transmissions at known wave-
numbers only. By employing a particular minimization technique, EPAGAS
tries to fit the observed transmissions to the theoretical transmissions
by varying the concentrations of the individual constituents.
Assuming a certain concentration of each gas constituent in an
atmospheric path, a corresponding transmission can be found from a large
table of known transmissions as a function of gas constituent, wavenumber,
and gas concentration. The total transmission at a specified wavenumber
is then computed and compared to the observed transmission value obtained
from a spectrophotometer. The difference between the computed and
observed transmissions is then minimized by varying the concentrations of
the individual gas constituents. The minimization process tries to reduce
the rms error of an error function defined on the errors of the given
and computed transmissions and on the concentrations. (See A2.3.) The
concentration errors are used mainly to help keep the concentrations within
such a range as to provide physically realistic solutions.
The program uses various controls to direct its progress in the
minimization process. The initial and final gas constituent concentrations,
computed transmissions, and errors are always printed together with a plot
of the given and computed transmissions as a function of wavenumber.
Intermediate results may also be printed and plotted.
For preparation of the theoretical gas library table, see the THEOTRN
description. For preparation of the given transmission data and all
121
-------
EPAGAS - (cont.)
other input needed by the program, see the READALL description. It is
imperative that all descriptions pertaining to the program be read and
understood.
USES:
CALC
CLOSER
CONST
CXMIZ
EROR1
EROR2
FINDIT
FINDTRN
GAUSS L
GRAPHS
IBINSER
MINMYZD
PARTIN
PLOTD
PLOTT
PRINT1
PRINT2
PRINTS
READALL
READ1
READ2
READS
READ4
READ6
READ?
READ8
RMSERR
THEOTRN
122
-------
FORTRAN SUBROUTINE CALC
CDC 3800 Margot Ackley
(Internal to EPAGAS) July 30, 1973
CALC (W, ERR, RMS)
Given the array W of NG current gas concentrations and the selected
arguments from the two COMMON blocks, COMMON/100/ and COMMON/101/, CALC
computes an error array ERR of length NG + NT, where NT = the number
of computed and given transmissions and NG ~ the number of gas
concentrations. CALC also computes the rms error RMS needed by
MINMYZD, and the individual rms errors RMSW for the error based only
on the concentrations, and RMST for the error based only on the
transmissions. (See A2.3.)
The two COMMON blocks utilized are as follows: COMMON/100/TRAN
(5000,5), NU(10), NW(10), UR(10,3), WR(10,5), TG(150), U(150), CF(150),
TC(150), WM(10), WHR(IO), WW(10), IPW(IO); COMMON/101/NG, NT, RMSTP,
RMST, RMSW.
The first NG elements of the array ERR are the errors based on
errors from the NG gas concentrations. (See EROR1 description.) The
NG + 1 elements through the NG + NT elements of the array ERR are the
errors based on errors from the NT given and currently computed
transmissions. (See EROR2 description.)
Before computing the individual transmission errors, CALC must
first calculate the NT computed transmissions TC based on the current NG
gas concentrations W given to CALC. CALC uses subroutine FINDTRN to find
the individual computed transmission TCI(i,j) for a given wavenumber
U(j) and a given gas concentration W(i). The total transmission for the
particular wavenumber U(j) is then the product of all the transmissions
for all the gas concentrations W(i) for the particular U(j) where:
NG
TC(j) - n
123
-------
CALC - (cont.)
If the current use of CALC produces an RMS error which is smaller
than the previous rms error, then RMSTP = RMS, and the array TC will
then contain the computed transmissions producing the lower RMS. If
RMS is greater, then the values RMSTP and TC will contain the last set
of values resulting in an improved rms error.
For definitions of arguments and variables used, see below
descriptions:
EPAGAS: W
READ7(8): NT, U, TG, CF
THEOTRN: NG
USES; EROR1, EROR2, FINDTRN(IBINSER), and RMSERR.
CALC: Computes concentration and transmission errors for MINMYZD.
124
-------
FORTRAN FUNCTION CLOSER
CDC 3800 J. R. Winkelman
(Internal to EPAGAS) June 15, 1965
CLOSER (A, X, B)
CLOSER is the value of X in the qlosed region A-X-B, taking on
the values of the end points for X outside the range. More exactly,
if X-A, the value of CLOSER will be A. If X-B, the value of CLOSER will
be B. Otherwise, CLOSER=X.
USES: No subroutines.
CLOSER: Restrict X to a closed region.
125
-------
FORTRAN SUBROUTINE CONST .
CDC 3800 Margot Ackley
(Internal to EPAGAS) July 30, 1973
CONST (S, UMIN, UMAX, ZEP, KX, Q, MK, GLB, CFM, IPRINT, IPLOT)
CONST sets certain constants needed by program EPAGAS and its
associated subroutines. These values are as follows:
S(l) = S(2) = -1. S(3) = S(4) = S(5) = 0.
WM(1) = 3.237E22 WHR(l) = 3.437E22
WM(2) = 1.235E20 WHR(2) = 1.235E20
WM(3) = 1.235E20 WHR(3) = 1.235E20
WM(4) = 1.235E20 WHR(4) = 1.235E20
WM(5) = 1.245E20 WHR(5) = 1.245E20
WM(6) = 0.715E20 WHR(6) = 0.715E20
WW(1) = WW(2) = WW(3) = WW(4) = WW(5) = WW(6) = 1
IPWC1) = IPW(2) = IPW(3) = IPW(4) = IPW(5) = IPW(6) = 1
UMIN = 1000. UMAX = 1250.
KX(1) = 10 ZEP = 0.6
KX(2) =3 Q = 10.0
KX(3j =50 GLB = 3.0
KX(4) =30 CFM = 10.0
KX(5) = 100 MK = 10
KX(6) =100 IPRINT = 1
KX(7) = 0 IPLOT = 1
KX(8) = KX(9) = KX(10) = 0
126
-------
CONST - (cont.)
The definitions of the above arguments are given in the below
subroutine descriptions.
S READALL
WM, WHR, WW, IPW READ2
UMIN, UMAX READ1
ZEP, KX READS
Q, GLB, CFM, MK READ4
IPRINT, IPLOT
CONST utilizes COMMON/100/ where: COMMON/100/TRAN(5000,5) ,
NU(10), NW(10), UR(10,3), WR(10,5), TGC150), UC150), CF(150),
TCC150), WM(10), WHR(IO), WW(10), IPW(IO).
USES: No subroutines.
CONST: Set values for EPAGAS.
127
-------
FORTRAN FUNCTION CXMIZ
CDC 3800 Margot Ackley
(Internal to EPAGAS) July 30, 1973
CXMIZ (KX, G, GLB, MK, RMS)
Given the array KX (see READ3 description) , the arguments Q and
GLB (see READ4 description) , MK the remaining uses of MINMYZD to be
permitted (see READ4 description), and RMS, the current rms error,
CXMIZ checks to see if MINMYZD is to be used again.
If KX(9) = -1, then CXMIZ =2.0 and returns. If KX(9) = 0, then
CXMIZ = 3.0 and returns. If KX(9) = 1, then CXMIZ checks the following:
if MK-0, if QY >KX(10), or if RMS-GLB, then CXMIZ = 2.0.
If none of the above occur, then CXMIZ = 1.0 and returns.
Where:
CXMIZ =1, use MINMYZD again,
CXMIZ = 2, do not use MINMYZD again,
CXMIZ = 3, singular matrix, halt computation.
USES: No subroutines.
CXMIZ: Checks for future uses of MINMYZD in EPAGAS.
128
-------
FORTRAN SUBROUTINE EROR1
CDC 3800 Margot Ackley
(Internal to EPAGAS) July 30, 1973
EROR1 (W, ERR)
Given NG values of W, where W(I) is the concentration of the I gas
constituent and corresponding error parameters WM, WHR, WW, and IPW, EROR1
computes the NG errors ERR as defined below:
IPW(I)
I wrri -WMfT** 1
ERR(I) = WW(I)'
(W(I)-WM(I)\
WHR(I) I
I = 1, NG
WM, WHR, WW, and IPW are provided through COMMON/100/: COMMON/100/
TRAN(5000,5), NU(10), NW(10), UR(10,3), WR(10,5), TG(150), U(150), CF(150),
TC(150), WM(10), WHR(IO), WW(10), IPW(IO).
USES: No subroutines.
EROR1_: Computes gas constituent errors for EPAGAS.
129
-------
FORTRAN FUNCTION EROR2
CDC 3800 Margot Ackley
(Internal to EPAGAS) July 30, 1973
EROR2 (TG, TC, CF)
EROR2 is the error computed from the discrepancy between the given
transmission TG and the computed transmission TC, weighted by a
confidence factor CF, where
EROR2 = CF*(TG-TC).
USES: No subroutines.
EROR2: Computes transmission errors for EPAGAS.
130
-------
FORTRAN FUNCTION FINDIT
CDC 3800 J. R. Winkelman
(Internal to EPAGAS) December 18, 1962
FINDIT (A, L, M, N)
The value of FINDIT will be that A(L) which would be obtained if
the array A was sorted between A(M) and A(N). A(L) has its correct
new value, that is, A(K)-ACL) for M-K-L, A(K)-A(L) for L-K-N. Thus,
the array is in better sort on exit from FINDIT. The array A outside of
the range M, N is unchanged.
Example: To find the 20 percentile and median of the array from A(10)
to A(20),
PMED = FINDIT (A, 15, 10, 20)
P20 = FINDIT (A, 12, 10, 14) or FINDIT (A, 12, 10, 20).
USES: PARTIN.
FINDIT: Find the Lth sorted value.
131
-------
FORTRAN SUBROUTINE FINDTRN
CDC 3800 Margot Ackley
(Internal to EPAGAS) July 30, 1973
FINDTRN (W, UP, TT)
Given the arguments W and UP and selected variables from the two
COMMON blocks COMMON/100/TRAN(5000,5), NU(10), NW(10), UR(10,3), WR(10,5),
TG(150), U(150), CF(150), TC(ISO) , WM(10) , WHR(IO) , WW(10) , IPW(IO) ; and
COMMON/101/NG, NT, RMSTP, RMST, RMSWFINDTRN finds the corresponding
theoretical transmissions TT.
The array W contains the NG current values of gas concentrations, and
UP is some given wavenumber within the range of interest. The computed
theoretical transmissions are placed in the array TT, which has a
maximum length of 10, where TT(I) corresponds to the theoretical
transmissions from the gas library as a function of wavenumber UP and
concentration W(I) of the I1" gas.
Because the current concentrations W and wavenumber UP may not
necessarily correspond to those wavenumbers and concentrations in the
gas library table TRAN (see THEOTRN description) , FINDTRN uses double
linear interpolation within the two-dimensional table TRAN to determine
the transmission TT(I) corresponding to W(I) and UP.
Wavenumber UP will never lie outside the range of the table.
However, the concentrations W could possibly lie outside the range due
to the minimization process which varies only the array W. In the event
that the individual W(I) lies above the range of the concentrations of
the I gas, the highest two concentrations of that gas are used to
determine a transmission value TT(I) for W(I) and UP. In the case where
W(I) lies below the range of the table, the smallest concentration in
the table and a concentration of zero are used in determining the
corresponding transmission. The interpolation formula for this case
utilizes the fact that, by definition, a transmission value for zero
concentration is 1. (See A2.4.) The interpolation process can possibly
132
-------
FINDTRN - (cont.)
result in transmissions outside the range 0 to 1, but the structure of
the error functions tends to discourage selections of those concentrations
W which make this occur.
USES: IBINSER.
FINDTRN: Finds theoretical transmissions as a function of wavenumber
and concentrations for EPAGAS.
133
-------
FORTRAN SUBROUTINE GAUSSL
CDC 3800 L. David Lewis
(Internal to EPAGAS) January 20, 1970
GAUSSL (A, NRD, NR, NC, D)
Given an augmented matrix A, its row dimension NRD (from the
DIMENSION statement defining A), the number of rows NR (equals the order of
the matrix), and the number of columns NC (equals NR plus the number of
vectors);
Find D, the determinant of the matrix and the solution vectors
of the indicated linear system. The given matrix is destroyed, and the
given vectors are replaced by the solution vectors.
Restrictions: A must be stored columnwise;
NC-NR, NR-NRD. If NC = NR, only D is computed.
Method: Gauss elimination with pivotal condensation.
Error conditions: D = 0, the matrix was singular; the augmented
matrix has been destroyed.
USES: No subroutines.
GAUSSL: Determinant of a matrix by Gauss elimination and solution.
134
-------
FORTRAN SUBROUTINE GRAPHS
CDC 3800 Margot Ackley
(Internal to EPAGAS) July 30, 1973
GRAPHS CA, B, Y, NGC, NXC, NUC, NOC, NBC, NQ)
Edits an array NQ of 101 alphanumeric characters (Rl format) to
allow computer plotting of functions.
For A
-------
GRAPHS - (cont.)
In use, to plot n functions, GRAPHS would be called n times, with
the n function and scale values, then a line would be printed. A
suitable print statement might be: WRITE OUTPUT TAPE 61, 99, X, (NQ(I),
1=1, 101), 99 FORMAT (!HbE16.9, 2X, 101R1), where X is the abscissa value
or perhaps one of the function values.
Restrictions: A must not equal B; zero (0) cannot appear as a plotted
character.
USES: No subroutines.
GRAPHS: Plot graphs on the line printer.
136
-------
FORTRAN FUNCTION IBINSER
CDC 3800 Margot Ackley
(Internal to EPAGAS) July 30, 1973
IBINSER (X, N, P)
Given an array X of length N containing monotonically increasing
numbers and given a number P, IBINSER performs a binary search in array
X to find the integer value I for which X(I)-P>X(I+1). IBINSER is then
set equal to I. that is, 1-IBINSER-N. If N=2, IBINSER*!.
USES; No subroutines.
IBINSER: Binary search of a monotonically increasing array.
137
-------
FORTRAN SUBROUTINE MINMYZD
CDC 3800 R. Gregg Merrill
(J. R. Winkelman)
(Internal to EPAGAS) March 20, 1970
MINMYZD (W, WI, NG, CALC, ERR, N, ZEP, KX, PRNT, RMS, XTM)
This subroutine is used to minimize a function of several variables
when the true value of the function is known at several points. In other
words, it fits an arbitrary function of NG variables to N known points.
For example, many positions of a satellite over a few days may be known
and the elements (period, major axis, etc.) are desired which will give
an orbit for which the rms distance from the given position is a
minimum.
Because the convergence procedures and flexibility are much greater
than would be required on "nice functions," a suggested set of control
factors KX will be given so that the average programmer does not need to
become familiar with all the options. (See summary.)
In detail, the subroutine starts with a given set of NG-10 parameters
W and an increment WI for each W. Within the WI array, a zero implies
that the corresponding W is not to be changed. Also given are the N-160
errors ERR (an error at each of the given points where these are of no
concern to MINMYZD) and an rms error RMS where RMS may be obtained by
RMS = RMSERR(ERR,N). NG may be 1, but N must be -NG.
MINMYZD makes a step in W. of WI. and looks at the resultant errors
obtained by calling CALC(W,ERR,RMS). If the RMS for this new set is smaller,
the slopes of the individual errors with respect to VI ^ are calculated and
stored. If larger, a step of WI. is made in the other direction and slopes
are calculated. If neither RMS is smaller, WI.. is halved for future use, as
the step size is probably too large. If no errors changed at all, WIi is
set to zero (normally because WI was halved out of existence and is no
longer important to the minimization process).
138
-------
MINMYZD - (cont.)
After all W's have been varied, the matrix
(BERR^SW )*(AW.) = -(ERR/) (1)
is solved for the changes AW. in W. The new parameters are then
Wi = Wi+ZEP*AWi If the important 3ERR./8W. is nearly constant, ZEP
should be 1. As a minimum is approached, ZEP must be made smaller. This
is done internally under control of the control values KX.
After the new ERR's and RMS have been calculated, PRNT is checked.
If PRNT is zero, no internal printing is done. If PRNT ? 0, a line is
printed specifying the ZEP used and the RMS corresponding to that ZEP.
It also prints a ZEP and RMS which would have theoretically minimized
P P
the curve (slopes are known for ZEP = 0, and the theoretical value RMS
for any ZEP can be compared against the actual RMS evaluated at the given
ZEP to predict the second derivatives and the corresponding minimum).
If at any time a new RMS is less than the old one, the W and ERR
arrays are replaced by the better values and the necessary adjustments
made to KX (7 and 8).
The KX array contains six input control parameters and four output
parameters. These control the action of MINMYZD. None of the controls
KX(l-6) are altered by the program.
KX(1) The total number of cycles to take, even if conditions are
still good.
KX(2) The number of cycles to count for a cycle which does not cause
enough improvement.
KX(3) The acceptable improvement in 0.1 percent (105 is 10.5 percent)
Any cycle resulting in less than this causes KX(2) cycles to
be counted rather than one cycle.
139
-------
MINMYZD - (cont.)
KD(4) The limit as to the multiple of a WI increment which may be
taken in one full step (ZEP=1). A zero will be used as 100.
A small slope might predict a step of 10 which is clearly
too large. If KX(4) is negative, no step is made in the
-WI direction. It is normally faster to have KX(4)<0 for
the first few uses of MINMYZD.
KX(5) Minimum ZEP to be used in 0.001 (10 is 0.010). If ZEP is allowed
to go too small (alas, the predictions are not always good)
the step does not go anywhere. If ZEP minimum is too large,
ZEP cannot be made small enough to account for the
nonlinearity in ERR for final convergence.
KX(6) If good improvements are made, it is natural for ZEP to become
smaller as MINMYZD moves away from the region in which the
derivatives were calculated. If the total improvement
KX(10) is greater than KX(6) in 0.1 percent, then the ZEP at
exit will be the same as the ZEP on entry, that is, things are
going well so stick to the same ZEP. Otherwise, exit with
the last predicted ZEP .
P
KX(7) through KX(10) are used to convey certain information to the
user.
KX(7) If 1, at least one single step (even if it were a setup of
WI.) resulted in an improvement of at least KX(3); otherwise,
it is zero.
KX(8) If 2, there was some improvement sometime; otherwise, zero.
KX(9) If -1, the WI's are all too small; the system has probably
converged. If 0, the matrix was singular; this is not expected
to occur. +1 is normal.
KX(10) Total improvement in 0.1 percent relative to the RMS input.
750 would mean the new RMS was % of the old RMS.
After PRNT is checked, the improvement is checked against KX(3).
If sufficient, one is subtracted from the cycle count set up from KX(1).
140
-------
MINMYZD - (cent.)
A ZliP gives an error RMS: (X =ai(KX(3)))
No
Any improvement? i.e.,
RMS^RMS0 where RMS0
is lowest rms error?
Does predicted ZDP_ give an
error RMSpSRMSo-XVRMSn?
Yes
Count as a bad cycle
but continue using ZEPp
Yes
Is RMS*RMS0-X%*RMS0? N»_
Nn
Count as a good
cycle. Use ZEPp
regardless of
Yes Is
Count as a bad cycle
but continue using ZEP*
where ZI-P*= Zl-P -ZOP
P
V
141
-------
MINMYZD - (cont.)
If the count is -1, the new ERR is put into equation (1) and the system
continues from there using the predicted ZEP . Any time the cycle count
reaches zero, an exit is made.
If the improvement is not good enough, KX(2) is subtracted from the
cycle count. If the count -0 or the predicted RMS is also not good
enough, an exit is made. Otherwise, a calculation is made using a ZEP
based on the old ZEP and the predicted ZEP , hoping for an improvement.
If no improvement is made but the predicted RMS is good enough, then
calculations are made using the predicted ZEP ; otherwise, an exit is made.
(See flow chart for termination of MINMYZD.)
Helpful hints: If KX(2) is one less than KX(1), then, after a "good"
cycle, one bad cycle will end it. If a "bad" cycle occurs immediately,
then one more try with the adjusted ZEP will usually give good results.
On complicated functions (almost anything will work on simple ones), it
may be necessary to set KX(5) higher from 100 to 300 for a few steps if ZEP is
being predicted too small. The XTM matrix, currently (160,11), is the
partial derivatives matrix XTM(I,K) = 3ERR./3W,, 1=1, N, K=l, NG, where K
1 K.
is the Kth non-zero WI that is, there are no zero columns.
If one or more rows of XTM completely dominate, it may result in
a singular matrix. This condition is checked, and the matrix is weighted
if necessary. The weights are stored in the column (KK) following the
normal information. This number XTM(I,KK) will normally be 1. If it is
less than 1, then the I row has been multiplied by XTM(I,KK) and appears
smaller than when the row was actually calculated.
SUMMARY:
Given: Arrays W and WI, NG words long where W contains the parameters,
WI contains the suggested increments, and a zero WI means that the corres-
ponding W is not to be varied; a subroutine CALC, which given W, calculates
an error array ERR (positive and/or negative), N long, and the rms error
RMS; a variable ZEP; and values for PRNT, KX(1), KX(2), KX(3), KX(4),
KX(5), and KX(6).
142
-------
MINMYZD - (cont.)
MINMYZD exits with improved values of W; WI may be slightly altered;
the array ERR of new errors associated with W; RMS correspond to ERR; ZEP
may be a better value (must be smaller toward the end); KX(7)=1, at least
one single step exceeded 0.1 percent improvement (from KX(3)) or else is
zero; KX(8)=1 some improvement occurred sometime or else is zero; KX(9)=1),
all is well, =0 unexpected, =-1 all WI's reduced to zero; you are done,
congratulations! KX(10) the relative improvement for this one step of
MINMYZD is 0.1 percent.
Normally KX(9) should be tested, and then test KX(7,8,or 10) or
RMS to determine whether to call MINMYZD again or to quit as good enough.
Note that all values needed for input are proper as returned from MINMYZD,
so only a call is necessary, no further work.
USES: CALC, CLOSER, internal MINMYZD (FINDIT (PARTIN)), and GAUSSL.
MINMYZD: To minimize an arbitrary function of several variables.
143
-------
FORTRAN SUBROUTINE PARTIN
CDC 3800 J. R. Winkelman
(Internal to EPAGAS) December 18, 1962
PARTIN (K, L, M, N, I, J)
Given the array K, a location in it K(L), and the limits of the
area of concern K(M) to K(N), PARTIN shuffles the array until the following
inequalities hold:
K(U) - K(L) for M - IJ - J
K(IJ) - K(L) for I - IJ - N
where M-J-I-N. I and J are output.
PARTIN uses fixed-point arithmetic internally so that it handles
fixed and floating arrays without prejudice.
USES: No subroutines.
PARTIN: Partition (partial sort) an array.
144
-------
FORTRAN SUBROUTINE PLOTD
CDC 3800 Margot Ackley
(Internal to EPAGAS) July 30, 1973
PLOTD (W, IDGAS)
Given the above arguments and selected arguments from the two COMMON
blocks, COMMON/100/ and COMMON/101/, PLOTD plots an on-line graph of the
difference between the given transmissions TG and final computed
transmissions TC versus the wavenumber U as a function of the final set of
gas concentrations W. Scale for the transmission differences is set
from -0.10 to +0.10.
The two COMMON blocks utilized are as follows: COMMON/100/TRAN
(5000,5), NU(10), NW(10), UR(10,3), WR(10,5), TG(150), U(150), CF(150),
TC(150), WM(10), WHR(IO), WW(10), IPW(IO); and COMMON/101/NG, NT, RMSTP,
RMST, RMSW.
For definitions of arguments and variables used, see below
descriptions:
EPAGAS: W
THEOTRN: IDGAS, NG
READ?(8): NT, U, TG
CALC: TC
PLOTD starts a new page.
USES: GRAPHS
PLOTD: Plots difference of given and final computed transmissions versus
wavenumber.
145
-------
FORTRAN SUBROUTINE PLOTT
CDC 3800 Margot Ackley
(Internal to EPAGAS) July 30, 1973
PLOTT (W, IDGAS, MM)
Given the above arguments and selected arguments from the two COMMON
blocks, COMMON/100/ and COMMON/101/, PLOTT plots an on-line graph of the
given transmissions TG and corresponding computed transmissions TC versus the
wavenumber U as a function of the current values of the gas concentrations
W. Scale for transmissions is set from 0 to 1.
The two COMMON blocks utilized are as follows: COMMON/100/TRAN
(5000,5), NU(10), NW(10), UR(10,3) , WR(10,5), TG(150), U(150), CF(150),
TC(150), WM(10), WHR(IO), WW(10), IPW(IO); and COMMON/101/NG, NT, RMSTP,
RMST, RMSW.
For definitions of arguments and variables used, see below
descriptions:
EPAGAS: W
THEOTRN: IDGAS, NG
READ?(8): NT, U, TG
CALC: TC
If NM=0, a title line stating that the plot is before use of MINMYZD
is printed; if NM>0, a title line stating that the plot is after NX
uses of MINMYZD is printed; and if NNKO, a title line stating that the plot
is after final use of MINMYZD is printed.
PLOTT starts a new page.
USES: GRAPHS
PLOTT: Plot curves of given and currently computed transmission versus
wavenumber.
146
-------
FORTRAN SUBROUTINE PRINT1
CDC 3800 Margot Ackley
(Internal to EPAGAS) July 30, 1973
PRINT1 (Q, GLB, CFM, MK, IPRINT, IPLOT, UMIN, UMAX, ZEP, KX, IDGAS, W,
WI)
Given all of the above arguments and selected arguments from the
two COMMON blocks, COMMON/100/ and COMMON/101/, PRINT1 prints all the
input data for program EPAGAS.
The two COMMON blocks are as follows: COMMON/100/TRAN(5000,5),
NU(10), NW(10), UR(10,3), WR(10,5), TG(150), U(150), CF(150) , TC(150),
WM(10), WHR(IO), WW(10), IPW(IO); and COMMON/101/NG, NT, RMSTP, RMST,
RMSW.
For definitions of arguments and variables printed, see below
descriptions:
READ4: Q, GLB, CFM, MK, IPRINT, IPLOT
READ1: UMIN, UMAX
READS: ZEP, KX
THEOTRN: IDGAS, NG
READALL: W, WI
READ2: WM, WHR, WW, IPW
READ?(8): NT, U, CF, TG
PRINT1 starts a new page.
USES: No subroutines.
PRINT1: Prints all input for EPAGAS.
147
-------
FORTRAN SUBROUTINE PRINT2
CDC 3800 Margot Ackley
(Internal to EPAGAS) July 30, 1973
PRINT2 (NM, Q, GLB, MK, ZEP, KX, RMS, IDGAS, WO, W, WI, ERR)
Given the above arguments and selected arguments from the two
COMMON blocks, COMMON/100/ and COMMON/101/, PRINT2 prints the intermediate
data and errors before the first use of MINMYZD and after any successive
uses of MINMYZD.
The two COMMON blocks are as follows: COMMON/100/TRAN(5000,5) ,
NU(10), NW(10), UR(10,3), WR(10,5), TG(150), U(150), CF(150), TC(150),
WM(10), WHR(IO), WW(10), IPW(IO); and COMMON/101/NG, NT, RMSTP, RMST,
RMSW.
For definitions of arguments and variables printed, see below
descriptions:
READ4: Q, GLB, MK
READS: ZEP, KX
CALC: RMS
THEOTRN: IDGAS, NG
READ6: WO
READALL: W, WI
CALC: ERR, RMST, RMSW, TC
READ2: WM, WHR, WW, IPW
READ?(8): NT, U, CF, TG
If NM=0, a title line stating that the data is before use of MINMYZD
is printed. If NM>0, a title line stating that the data is after NM uses
of MINMYZD is printed.
PRINT2 starts a new page.
USES: No subroutines.
PRINT2: Prints intermediate data and errors for EPAGAS.
148
-------
FORTRAN SUBROUTINE PRINTS
CDC 3800 Margot Ackley
(Internal to EPAGAS) July 30, 1973
PRINTS (IDGAS, WO, W, ERR, RMS)
Given the above arguments and selected arguments from the two
COMMON blocks, COMMON/100/ and COMMON/101/, PRINTS prints the final data
and errors.
The two COMMON blocks are as follows: COMMON/100/TRAN(5000,5),
NU(10), NW(10), UR(10,3), WR(10,5), TG(150), U(150) , CF(150), TC(150),
WM(10), WHR(IO), WW(10), IPW(IO); and COMMON/100/NG, NT, RMSTP, RMST, RMSW.
For definitions of the arguments and variables printed, see below
descriptions:
THEOTRN: IDGAS, NG
READ6: WO
READALL: W
CALC: ERR, RMS, RMST, RMSW, TC
READ2: WM, WHR, WW, IPW
READ7(8): NT, U, CF, TG
PRINTS starts a new page.
USES: No subroutines.
PRINTS: Prints final data and errors for EPAGAS.
149
-------
FORTRAN SUBROUTINE READALL
CDC 3800 Margot Ackley
(Internal to EPAGAS) July 30, 1973
READALL CS, UMIN, UMAX, ZEP, KX, Q, GLB, MK, CFM, IPRINT, IPLOT, W, WI,
WO, WIO)
Given the array S of length 5 (see CONST description) , READALL reads
in and checks input data to be used by EPAGAS. Each of the seven different
types of data must be preceded by the appropriate header card containing
the card type number in column 1.
Error array
Card type element FORTRAN reading subprogram
1 S(3) READ1 (UMIN, UMAX)
2 S(4) READ2 (FAKEARG)
3 none READS (ZEP, KX)
4 none READ4 (Q, GLB, CFM, MK, IPRINT,
IPLOT)
6 S(l) READ6 (WO, WIO)
7 S(2) READ7 (UMIN, UMAX, CFM)
8 S(2) READS (UMIN, UMAX, CFM)
S(K) = 0 indicates no error in reading and S(K) f 0 indicates that
an error has occurred. If a card type 5 is read in, S(5) = -1 and
indicates an error. See above subprogram descriptions for preparation of
data. If 1, 2, 3, 4 type data is not read in, then those values given by
CONST will be used in computation. Type 6 data must initially be read in,
and type 7 or_ 8 (but not both) must also be initially read in. As many
sets of data as desired may be processed. However, a type 9 card must
follow every set of data (9 in column 1). This card will cause checking
150
-------
READALL - (cont.)
of data and begins computation if no reading errors occur. READALL
provides printed error analysis for any errors.
Arrays W and WI are set equal to the arrays WO and WIO. These
are the gas concentrations W with corresponding increments WI, and the
original concentrations WO with corresponding original increments WIO.
If no concentrations are read in after the initial set, then each new
set of data will be processed, starting with the original set of
concentrations WO and increments of WIO. If one wishes to use the final
computed set of concentrations of the last data set as the starting
concentrations of the next data set, this may be accomplished by using
READ6 with the appropriate header card, followed by a card with a 2 in
column 80. This last set of concentrations and increments will then
become the original concentrations and increments for all future sets
of data.
After the initial set of data, sets of data need only to contain
those values which one wants to vary. The unchanged values which will
be used in processing the new data set will be those which had last been
read in.
READALL utilizes two blocks of COMMON: COMMON/100/TRAN(5000,5),
NU(IO), NW(10), UR(10,3), WR(10,5), TG(150), UC150), CFC150), TC(150),
WM(10), WHR(IO), WW(10), IPW(IO); and COMMON/101/NG, NT, RMSTP, RMST,
RMSW.
USES: READ1, READ2, READ3, READ4, READ6, READ7, and READS.
READALL: Reads in and checks data for EPAGAS.
151
-------
FORTRAN FUNCTION RE ADI
CDC 3800 Margot Ackley
(Internal to EPAGAS) July 30, 1973
READ1 (UMIN, UMAX)
READ1 reads in the minimum wavenumber value UMIN and the maximum
wavenumber value UMAX for which the given transmissions will range over.
If UMIN-0 and/or UMIN-UMAX, then an error has orrcurred and READ1
returns equal to -1. If no error occurs, then READ1 = 0.
CARD FORMATS:
CONTROL CARD-read by READALL
Col. Field Description
1 Integer Card type - use 1
DATA CARD-read by READ1
1-10 Floating UMIN>0
11-20 Floating UMZX>UMIN
USES; No subroutines.
READ1: Reads in wavenumber range for EPAGAS.
152
-------
FORTRAN FUNCTION
CDC 3800
(Internal to EPAGAS)
READ2 (FAKEARG)
READ2
Margot Ackley
July 30, 1973
Given NG, the number of gas constituents, READ2 reads in the
necessary gas concentration error parameters to be used by EROR1. NG
values each of WM, WHR, WW, and IPW are read in.
If WW(K)<0, WM(K)<0, WHRQQ-0, and/or IPW(K)-0 for K-l, NG, then a
reading error has occurred and READ2 returns equal to -1. If no error
occurs, then READ2=0. See EROR1 description for detailed use of WW, WM,
WHR, and IPW. FAKEARG is a fake argument and has no real use in the
function. It is set to zero.
READ2 utilizes two blocks of COMMON: COMMON/100/TRAN(5000,5),
NU(10), NW(10), UR(10,3), WR(10,5), TG(ISO), U(150), CF(150), TC(150),
WM(10), WHR(IO), WW(10), IPW(10); and COMMON/101/NG, NT, RMSTP, RMST, RMSW.
CARD FORMATS:
CONTROL CARD-read by READALL
Col.
Field
Description
Integer
Card type - use 2
1-13
21-33
41-53
69-70
DATA CARD-read by READ 2
E13.6
E13.6
E13.6
Integer
NMQQ-0
WHR(K)>0
WW(K)-0
IPW(K)>0
. K=l, NG
USES: No subroutines.
READ2: Reads in gas concentration error parameters,
153
-------
FORTRAN SUBROUTINE
CDC 3800
(Internal to EPAGAS)
READ3 (ZEP, KX)
READ5
Margot Ackley
June 30, 1973
READS reads in the necessary controls needed for subroutine MINMYZD
which are ZEP and the first six elements of the array KX.
CARD FORMATS:
CONTROL CARD-read by READALL
Col.
Field
Description
Integer
Card type - use 3
DATA CARD-read by READ3
1-10
11-20
Floating
Integer
21-30
Integer
ZEP (0
-------
READ5 - (cont.)
Col.
31-40
41-50
Integer
51-60
61-70
Integer
Integer
Description
KX(3), the percent improvement in
tenths in the rms error to be
accepted as a good cycle. Any cycle
resulting in less than this causes
KX(2) cycles to be counted rather
than one cycle, KX(3)-0.
KX(4), the limit as to the multiple
of a WI increment which may be taken
in one full step CZEP=1). If set to
zero, then 100 is used. If KX(4) is
negative, then no step will be
permitted in the -WI direction for the
first complete try of MINMYZD.
Future tries of MINMYZD, use KX(4) .
KX(5), minimum ZEP in 0.001's,
KX(5)>0.
KX(6), total percent improvement
in tenths, KX(6)>0.
USES: No subroutines.
READS: Reads controls for MINMYZD.
155
-------
FORTRAN SUBROUTINE READ4
CDC 3800 Margot Ackley
(Internal to EPAGAS} July 30, 1973
READ4 (Q, GLB, CFM, MK, IPRINT, IPLOT)
READ4 reads in the necessary program controls for EPAGAS: Q is
the number of percent by which each complete try of MINMYZD must improve
the mean square deviation when using CALC, GLB is the minimum rms
error permitted before using MINMYZD again, MK is the maximum number of
uses of MINMYZD permitted, and CFM is a constant confidence factor to
be applied to the given experimental transmission data if individual
ones are not provided.
IPRINT=n permits printing of intermediate results using PRINT2
after every nth use of MINMYZD, and IPLOT=m permits a graph of the given
and computed transmissions as functions of wavenumber to be plotted
after every mth use of MINMYZD. If either IPRINT and/or IPLOT is zero,
then no intermediate printing and/or plotting will occur.
CARD FORMATS:
CONTROL CARD-read by READALL
Col. Field Description
1 Integer Card type - use 4
DATA CARD-read by READ4
1-10
11-20
21-30
31-40
41-50
51-60
Floating
Floating
Floating
Integer
Integer
Integer
Q>0
GLB-0
CFM>0
MK>0
IPRINT-0
IPLOT^O
156
-------
READ4 - (cont.)
USES: No subroutines.
READ4: Read program controls for EPAGAS.
157
-------
FORTRAN FUNCTION READ6
CDC 3800 Margot Ackley
(Internal to EPAGAS) July 30, 1973
READ6 (WO, WIO)
Given NG, the number of gas constituents, READ6 reads in NG
original gas concentrations WO(K) and their corresponding increments
WIO(K), K=l, NG. If the first data card read in contains a 2 in column
80 instead of W0(l) and WIO(l), then the original start concentrations
and increments will be set to the final concentrations and increments
of the previous data run and READ6 is set to -1. If not, then READ6
reads in the rest of WO and WIO values and returns equal to zero.
READ6 utilizes one COMMON block where: COMMON/101/NG, NT, RMSTP,
RMST, RMSW.
CARD FORMATS:
CONTROL CARD-read by READALL
Col. Field Description
1 Integer Card type - use 6
DATA CARD-read by READ6
1-13 E13.6 WO(K)
21-33 E13.6 WIO(K) K=l, NG
USEjS: No subroutines.
READ6: Reads in starting gas concentrations and corresponding increments.
158
-------
FORTRAN FUNCTION READ7
CDC 3800 Margot Ackley
(Internal to EPAGAS) July 30, 1973
READ? (UMIN, UMAX, CFM)
Given UMIN, UMAX, and CFM, where UMIN and UMAX are, respectively, the
minimum and maximum wavenumbers for the corresponding given transmissions
and CFM is a constant confidence factor to be used for each individual
transmission value if separate factors are not provided; READ? reads
in individual transmissions TG(K) and corresponding confidence factors
CF(K). The transmissions are assumed to be equally spaced in wavenumber
starting at UB, stepping in UI increments, and ending at UE. NT number
of pairs of TG(K) and CF(K) are read in where NT=1 + (UE-UB)/UT. The first
data card (a.) read in gives the values UB, UE, and UI. The following
data cards (b.) contain TG and CF.
If any of the following conditions occur: UBUMAX, UI-0,
NT>150, or TG(K)<0, TG(K)>1, CF(K)<0, K=l, NT, then an error has occurred
and READ? returns equal to -1. Otherwise, if no error occurs, then
READ7=0 upon return. If CF(K)=0, then CF(K) is set to CFM for all
appropriate K=l, NT.
Before return, the array U(K), K=l, NT is filled where U(K)=UB +
READ? utilizes two blocks of COMMON: COMMON/100/TRAN(5000,5) ,
NU(10), NW(10), UR(10,3), WR(10,5), TG(150), U(150) , CF(150), TC(150),
WMC10), WHR(IO), WWC10), IPW(IO); and COMMON/ 1 01 /NG, NT, RMSTP, RMST,
RMSW.
159
-------
CARD FORMATS:
READ7 - (cont.)
CONTROL CARD-read by READALL
Col.
Field
Description
Integer
Card type - use 7
DATA CARD (a.) - read by READ7
1-10
11-20
21-30
Floating
Floating
Floating
UB
UE
UI
UB-UMIN
<
UE-UMAX
UI>0
DATA CARD (b.) - read by READ7
1-5
6-10
11-15
16-20
21-25
26-30
31-35
36-40
F5.4
F5.5
F5.4
F5.5
F5.4
F5.5
F5.4
F5.5
TG(j)
CF(j)
TG(j+l)
CF(j-H)
TG(j+2)
CF(j+2)
TGCJ+3)
CFCJ+3)
j-1. 9, 17, ...
-------
READ? - (cont.)
USES: No subroutines.
READT^ Reads in given transmissions and corresponding confidence factors
for EPAGAS.
161
-------
FORTRAN FUNCTION READS
CDC 3800 Margot Ackley
(Internal to EPAGAS) July 30, 1973
READS (UMIN, UMAX, CFM)
Given UMIN, UMAX, and CFM, where UNIM and UMAX are, respectively, the
minimum and maximum wavenumbers for the corresponding given transmissions
and CFM is a constant confidence factor to be used for each individual
transmission value if separate factors are not provided, READS reads in
the individual wavenumber values U(K), the corresponding given
transmission TG(K), and its confidence factor CF(K), K=l, NT. The first
data card (a.) read in gives the value of NT, the number of given
transmissions to be read. The following cards (b.) contain TG, U,
and CF.
If any of the following conditions U(K)UMAX, TG(K)<-0.2,
or TG(K)>1.01, or CF(K)<0, K=l, NT, or NT>150 occur, then an error has
been encountered and READS returns equal to -1. Otherwise, no error
occurred and READ returns equal to 0. If CF(K)=0, then CF(K) is set to
CFM for all appropriate K=l, NT.
READS utilizes two blocks of COMMON: COMMON/101/TRANC5000,S) ,
NU(10), NW(10), UR(10,3), WR(10,5), TG(150), U(150), CF(150), TC(150) ,
WM(10), WHR(IO), WW(10), IPW(IO); and COMMON/101/NG, NT, RMSTP, RMST,
RMSW.
CARD FORMATS:
CONTROL CARD-read by READALL
Col. Field Description
1 Integer Card type - use 8
162
-------
READS - (cont.)
DATA CARD
(a.) - read by READS
Col.
Field
Description
1-3
Integer
NT
(KNT-150
DATA CARD
(b.) - read by READS
1-10
11-20
21-30
Floating
Floating
Floating
TG(K)
U(K)
CF(K)
K=l, NT
USES: No subroutines.
READS; Reads in given transmission and corresponding wavenumbers and
confidence factors for EPAGAS.
163
-------
FORTRAN FUNCTION RMSERR
CDC 3800 J- R- Winkelman
(Internal to EPAGAS) June 15, 1965
RMSERR (A, N)
Given an array A of length N, RMSERR is the root mean square of
the array A where
h
N j |
RMSERR =I £
USES: No subroutines.
RMSERR: Root mean square of an array.
164
-------
FORTRAN FUNCTION THEOTRN
CDC 3800 Margot Ackley
(Internal to EPAGAS) July 30, 1973
THEOTRN (IDGAS)
THEOTRN reads in the gas library of theoretical transmissions
for a given number of specified gases as a function of both gas
concentration and wavenumber. If any errors occur in reading or certain
variables are outside permitted ranges, THEOTRN returns immediately
equal to one. If all reading is correct, then THEOTRN returns equal to
zero.
Before the actual gas library is read in, THEOTRN first reads in
a set of data cards that specify which gases are in the library, which
concentrations for each gas appear, and which wavenumber region is used for
each gas.
Variables appearing on the data cards with their corresponding
dimensions and restrictions, if applicable, are given below:
NG: Total number of different gases to be found in
the library, (KNG-10.
IDGAS: Array of maximum length 10 which contains a
numerical (12) identifier for each of the gases,
where IDGAS(I) « identifier for the Ith gas, 1=1,
NG.
NU: Array of maximum length 10 which contains the
number of wavenumber values used for each gas,
where NU(I) = number of wavenumber values for the
Ith gas, 1=1, NG.
NW: Array of maximum length 10 which contains the number
of concentrations used for each gas, where NW(I) =
number of concentrations for the I gas, 1=1, NG;
0
-------
THEOTRN - (cont.)
UR: A two-dimensional array UR(I,J) dimensioned (10,3).
UR contains the wavenumber range and spacing used
for each of the gases. I corresponds to the Ith
gas. J=l corresponds to the beginning wavenumber,
J=2 corresponds to the end wavenumber, and J=3
corresponds to the wavenumber increment. Note
that 1. + (UR(I,2) - UR(I,1))/UR(I,3) should be
equal to NU(I).
WR: A two-dimensional array WR(I,J) dimensioned (10,5).
WR contains the actual concentrations used for each
gas. I corresponds to the I gas. J corresponds
to the J concentration for the particular I
gas and must be given in strictly increasing order
where l.E-20
-------
THEOTRN - (cont.)
Col.
Field
Description
IDGAS(NG)
Card 3.
1-4
5-8
Integer
Integer
NU(1)
NU(2)
NU(NG)
Card 4.
1-4
5-8
Integer
Integer
NW(1)
NW(2)
NW(NG)
Card 5.
1-8
9-16
17-24
25-32
33-40
41-48
49-56
57-64
65-72
73-80
Floating
Floating
Floating
Floating
Floating
Floating
Floating
Floating
Floating
UR(1,2)
UR(1,3)
UR(2,1)
URC2.2)
UR(2,3)
UR(3,1)
URC3,2)
UR(3,3)
Blank
Use successive cards up to a maximum of four.
167
-------
THEOTRN - (cont.)
Card 5. +J, l-J-13
Col.
1-10
11-20
Field
E10.3
E10.3
Description
WR(I,1D
WR(I,2) 1=1, NG
WR(I,NW(I))
Use successive NG cards, one for each gas.
THEOTRN next reads in the actual gas library.
LIBRARY CARD FORMAT:
Description
NAME
ICON
NN
T(J)
T(J-H)
T(J+2)
T(J+3)
Col.
1-2
3-4
5-8
9-12
13-16
17-20
21-24
Field
Integer
Integer
Integer
F4.3
F4.3
F4.3
F4.3
77-80
F4.3
T(J+17)
The gas with identification IDGAS(l) must come first, then the gas
with IDGAS(2) is next, etc. Gas with identification IDGAS(NG) comes last.
NAME on each gas library card corresponds to IDGAS of the particular
gas being read in. Within a particular gas, the concentrations W must be
increasing. The first concentration corresponds to ICON=1, second
concentration corresponds to ICON=2, and last concentration corresponds to
168
-------
THEOTRN - (cont.)
. . - - -~ - - -
ICON=NW(I) for the I gas. Within a given gas and a given concentration,
the cards are serialized by NN where NN=1 for the first card, NN=2 for
the second card, etc. The maximum value for NN within a set is NU(I)/18
for the I gas. T(J) corresponds to the theoretical transmission for
the particular Ith gas and its particular concentration. T(l) is the
transmission corresponding to first wavenumber UR(I,1), T(2) is the
transmission corresponding to the wavenumber UR(K,1) + UR(I,3), and
the last transmission, T(NU), corresponds to the last wavenumber
UR(I,2) for the I gas. Within a given gas and particular concentration,
all cards have 18 values of transmission except for the last card which
may have less.
THEOTRN does extensive checking of the card order when reading in
the library. If the cards are out of order, an error has occurred and
THEOTRN terminates equal to one.
As the individual transmission T is read in, it is stored in
the large transmission array TRAN (I,J) which has been dimensioned
(5000,5). J corresponds to the concentration number ICON, where l-J-5
for all gases. 1=1 corresponds to the transmission of the first gas,
first wavenumber. I=NU(1) corresponds to the transmission of the
first gas, last wavenumber. I=NU(1) + 1 corresponds to the transmission
of the second gas, first wavenumber. I = NU(1) + NU(2) corresponds
NG
to the transmission of the second gas, last wavenumber. T _
< * ~
-5000 corresponds to the transmission of the last gas,
last wavenumber. Unless all gases have the same number of concentrations,
the actual amount of storage utilized will not be rectangular.
THEOTRN utilizes two COMMON blocks COMMON/ 100/ and COMMON/ 101/
where COMMON/100/TRAN(5000,5) , NU(10), NW(10), UR(10,3), WR(10,5),
TG(150), U(150), CF(150), TC(150), WM(10) , WHR(IO), WW(10) , IPW'(10) ; and
COMMON/101/NG, NT, RMSTP, RMST, RMSW.
USES: No subroutines.
THEOTRN : Reads in gas library of theoretical transmissions for EPAGAS.
169
-------
A2.2 SIMULATION OF GAUSSIAN NOISE
To study the effects of noise in the transmission data, a method
of obtaining random Gaussian noise was developed. Many large-scale
computers (CDC 3800 in particular) provide a library routine which
furnishes the user with a sequence of uniformly distributed random
numbers ranging from 0 to 1. Given two such numbers A. and B., a
Gaussian random variable X. can be obtained from a transformation
provided by Box and Muller (1958) , where
X., = (-2 logeA.)^ cos (2JIB.) .
Three different signal-to-noise ratios 10, 30, and 100 were
incorporated in this report. Letting SNR be the signal-to-noise ratio
desired, then SNR is defined as
SNR =
n
where o, is the standard deviation of the data and a is the standard
d n
deviation of the Gaussian noise for the set of data to be analyzed. The
noise N. used to modify a particular transmission value T. is
defined as
N. = CT X. ,
i n i '
where the noise over the set has mean zero and standard deviation a .
n
The new transmission values TN. incorporating the noise are then
defined as
TN. = T. + N. .
111
170
-------
A2.3 CONSTITUENT AND TRANSMISSION ERROR FORMULATION
To direct its progress, the minimization process uses a
total rms error E as a measure of the "goodness of fit." By varying
the amounts of each of the gas constituents, the process tries to
minimize the total rms error E which is both a function of the rms
error ET from the individual transmission errors and the rms error
EW from the errors of the individual constituent amounts. (See EROR2
and EROR1 descriptions, respectively.)
For the benefit of the user, all three rms errors E, ET, and E
i w
are standard printed variables. In this program, the total rms
error E and the constituent rms error E are calculated. From these
w
two, the transmission rms error E_, can then be calculated as follows.
Let NG = the number of gas constituents,
NT = the number of transmission values,
A. = the individual errors from the gas constituents,
i=l, NG,
A. = the individual errors from the transmissions,
i=NG+l, NG+NT,
N = NG + NT.
But
Then,
E =
N
NG
E2 =
NG
NG
I*
1 NG 7
* V A2
NG A* i
1=1
NG+NT=N
*
N
.2 NT
i N~
^F * £
NT
i=NG+l
A2
171
-------
and
Therefore,
N * E2 - NG *_E?
NT
172
-------
A2.4 INTERPOLATION FORMULAS USED IN TRANSMISSION LIBRARY
Given a wavenumber v and a gas concentration W for a specified
constituent, a corresponding transmission T(v,W) can be found from the
transmission library. Because the library contains transmissions at
discrete values of wavenumber and concentration, a form of double
linear interpolation must be utilized.
Assuming the v and W lie within the region of the transmission
table, that is,
W. - W, - W - HL - W
mm 1 2 max
v . < < < <
mm - V - v - V -
then by interpolating first in the W direction, we have,
and
W -
T(V2,W) -
W -
Therefore,
W-W
and
W-W
TCv2,w) = _i *
2 W2 W2
Interpolating next in the v direction, we have
T(V,W) -
= T(v2,W) -
V - V,
(A2:2)
(A2:3)
173
-------
Therefore,
T(v.W) -
" T(V1»W)
(A2:4)
By substituting the expressions for TCv^W) and T(v ,W) from (A2:2)
into equation (A2:4), we have the following interpolation formula
v-v
W ^
!'V "
-T(vrW2)
(A2:5)
Because of the structure of the program, the wavenumber v will
always lie within the region of the transmission table. However, the
concentration W for a particular gas constituent could possibly lie
outside the region due to the minimization process which varies only
the W's. In the event W is outside the region, we have two possible
cases:
W < W
max
W < W .
nun
For the first case, let
W , = W.
-------
T(vrW) -
W-W,
w-w.
(A2:6)
and
T(v,W) - T(v
T(v2,W) - TCv^
v-
(A2:7)
By referring to (A2:l) and (A2:3), we can then use
(A2:5) for the final formula with W, = W , and W0 = W
l max-l 2 max
For the second case, where
, = 0 - W < W . = W_
1 - rain 2
or
Wl ° < Wmin '
we may utilize the fact that
T(v,0) = T(v,W ) = 1 for all v.
(A2:8)
Applying (A2:8) to (A2:5), we have
v-v.
r2
W-Q *
^ *[ T(v2,W2)-l+l-T(v1,W9)
1.
1-1
(A2:9)
175
-------
Therefore,
T(V,W) = 1 + ^ *
T(V1,W2) -
where
v-v
*
W2 = Wmin'
(A2:10)
In summary, given v and W to find T(v,W), where v . - v., - v - v_ - v :
< <, < < mm 1 2 max
For W . - W, - W - W_ - W
mm 1 2 max
For
For
V-V
T - ^ *
-T(vrW2)J + Tj
W-W
,,W2) -
W < W,
max '
T(V,W) =
v-v.
Vvi
w-w
max-l
W -W ,
max max-l
* T(Vn,W ) - TfV0-W ,) +
1 v 2* max' v 2 max-l'
,,W J - TfV^.W ) T(V..W .) - T(V.,W .
1' max-l' v 1' max-* | ^2' max-l' * 1' max-l'
W-W
max-l
W -W .
max max-1
W < W . ,
min'
T(v,W) =
T(V.,W ) - T(V.,W .
v 1 maxj * 1* max-l
1f ,
lf max-l
W
min
1-1,
176
-------
A2.5 PROGRAM EXECUTION TECHNIQUES
Although program EPAGAS is a totally self-contained program, its
proper utilization requires art as well as science if successful results
are desired. Before attempting a run, the user should thoroughly acquaint
himself with all the program and subprogram descriptions. The following
categories require careful checking before the program is run and are a
trouble-shooting checklist in case of improper results.
A2.5.1 Theoretical Gas Library Spectra
The user cannot hope to have good results if the library of spectra
does not contain all of the gas constituents which are (or are suspected)
to occur in the experimental data. Also the wavenumber range for the
gas constituents may not be sufficient for the range of the experimental
data. If this occurs, an error message will be printed. The immediate
solution to this problem is either to increase the wavenumber range of the
library spectra or to analyze a'smaller region of the experimental data.
The range of the concentrations of a particular gas constituent is
another item that should be checked. The library should have a range of
concentrations that covers the expected concentration of a particular
gas constituent. If, during the process of trying to find a solution,
EPAGAS tries to use a concentration outside the region of given
concentrations, a method of linear extrapolation is used to derive a
particular gas concentration. For concentrations less than 10% of the
total range of the region lying beyond the outer limits of the library,
there generally will not be problems. However, concentrations greater
than 10% could cause the program to give erroneous results.
The spacing of the library spectra in both wavenumber and
concentration is also very important. The program uses linear
interpolation in both wavenumber and concentration to determine
transmissions for a given wavenumber and a trial concentration. (See
FINDTRN decription.) If the grid spacing is not fine enough to give
reasonably correct results with linear interpolation, the program may not
be able to give satisfactory results. The fineness of the library is
177
-------
often limited by available computer memory and by the total number of gas
constituents present. In the case where the user has experimental data
containing many gas constituents, it might be necessary to limit the wave-
number region to a fairly narrow band or to pick regions where only a few
gas constituents are prominent.
The program has extensive built-in checks to make sure that the
library spectra and associated wavenumber and concentration information
are properly read in. If any of the library cards are out of order or
if inconsistencies occur, an error message is printed and further execution
is terminated. In the event this happens, carefully check the THEOTRN
description for proper input of the library spectra.
A2.5.2 Experimental Transmission Data
Unsuccessful results can often be caused by experimental data which
have an intolerable noise level. Signal-to-noise ratio, S/N, should be on
the order of 100 or better. For S/N<100, the user should select wavenumber
regions that have only a single gas constituent or regions which have
little interference from other gases.
Before utilizing EPAGAS, the user may find that it is necessary to
do some "preprocessing" of the data. Due to a variety of reasons,
instrumental drift in particular, the data may have to be corrected
either linearly or nonlinearly for both wavenumber and/or transmission.
Absolute transmission values of I/I without noise should range from 0 to
1. Past experience in the use of program EPAGAS has shown that wavenumber
drift is particularly detrimental to obtaining successful results. Wave-
number and intensity calibrations can be obtained from careful comparison
with known spectra.
Another very important item to check is the instrumental resolution.
Resolution should be the same or very nearly the same as the resolution
used in producing the library spectra. It is recommended that instrumental
resolution should deviate from that of the library by no more than 10 percent.
Also the resolution should not fluctuate as a function of wavenumber which
can be a major problem for some types of spectrometers unless it is
corrected.
178
-------
The temperature and pressure of the experimental data should also
correspond as close as possible to those used in producing the library
spectra.
In order that the experimental data is well defined, spacing in wave-
number should be sufficiently fine. If the spacing is too coarse,
important spectral line features may be lost. If this is the case, a
suitable solution would almost be impossible to find. Wavenumber spacings
of at least 1/5 of the resolution are strongly recommended.
Two other important aspects of the experimental data should also
be taken into consideration. The slit function of the instrument should
closely approximate that of the library. The digital resolution in the
conversion of the analog I/IQ information to digital form has a definite
effect on the final accuracy of determinations of gas concentrations.
The higher the digital resolution, the more accurately EPAGAS can
determine the concentrations.
A2.5.3 Minimization or "Goodness of Fit" Controls
Often less than desirable results are obtained because EPAGAS has
not been permitted to "work hard enough" at finding a satisfactory
solution. A thorough understanding of the descriptions of MINMYZD, READ3,
READ4, and CXMIZ would be extremely helpful (along with some experience
in using EPAGAS) in remedying this problem. MINMYZD description gives
the main minimization procedures and presents a discussion of the con-
vergence controls that direct its progress toward a solution. READS des-
cription provides information on the controls, used by MINMYZD, which are
labelled ZEP and array KX. Particular attention should be paid to the first
six elements of the array KX. The description of READ4 gives a discussion
of the variables Q, GLB, and MK which are used in conjunction with the 9th
and 10th elements of the array KX by CXMIZ to determine if another try by
MINMYZD to find a solution should be made. It is recommended that when
first using the program on a new set of experimental data, the
controls should be set as loosely as computer costs will permit.
179
-------
A2.5.4 Error Parameters
Proper formulation of the errors defined on the individual gas
concentrations is most helpful for successful utilization of EPAGAS. The
user though must keep in mind that the term "error" is used rather loosely
in connection with the gas concentrations. These errors are mainly
a way in which to guarantee that the minimization process does not try
to select individual gas concentrations which are known to be physically
unsuitable. The errors (or better, restraints) on the individual gas
concentrations are of the general form,
i* A
rIPW(I)
W(I)-WM(I) "
WHR(I)
where I = i gas constituent. W is the current trial concentration being
used by the minimization process, and all the other variables are input
error parameters. (See EROR1 description.) WM is considered to be the
mean acceptable concentration. For most cases (unless additional
knowledge exists of the amount of the gas constituent present), WM is simply
the mean concentration value based on the highest and lowest library
concentrations of the particular gas constituent. WHR is the half range
of permitted variation of the particular concentration and is simply the
difference between the highest library concentration and the mean
concentration. As long as W remains within the concentration range of
the library, the absolute value of
in a fairly small error. If the
W-WM
WHR
will be between 0 and 1, resulting
minimization process tries to
select W's outside of the library region, the error will become rapidly
larger as W gets farther and farther away from the acceptable region. The
variable IPW serves the purpose of determining just how fast the errors
will grow as W moves farther and farther away. IPW should always be an
odd, positive number to preserve the "direction" of the error. The larger
IPW is, the faster the error will increase once W has entered the unacceptable
regions. From past experience, IPW=3 usually is a good starting value.
WW is a weighting or normalization factor and is used to keep individual
180
-------
gas constituent errors in line with each other and also in proportion
with the individual transmission errors. Sometimes several runs have to
be made before a proper determination of the WW's can be made; however,
unity is always good for an initial try.
It is generally best to keep the gas constituent error parameters
fairly loose; however, if prior knowledge of the concentration of a
particular gas constituent in the experimental data exists, then one may
wish to narrow the acceptable region around the suspected concentration.
This would be done by setting WM equal to the suspected concentration,
by considerably decreasing WHR, by increasing IPW up to 9 or 11, and by
setting WW to 10 or 20 times the values of WW used for the other gas
constituents.
The errors, depending on the computed and given transmission values,
are of the general form,
ERR(I) = CF(I)*(TG(I)-TC(I)),
where I = i value corresponding to the i wavenumber. TG is the given
experimental transmission, TC is the currently computed transmission based
on the current gas concentrations as determined from the minimization
process, and CF is a confidence or weighting factor. The role of CF is
twofold. If some regions of the experimental data are more reliable
than others, then the corresponding CF's should be larger than those
applied to transmissions in the more uncertain regions. This has an effect
on the minimization process so as to let it somewhat "ignore" the less
reliable data and not to spend its efforts seeking some erroneous solution
to very noisy or bad data. Also, the CF's play a similar role for the
transmission errors that the WW's do for the gas concentration errors. Care
should be taken to select the CF's so that the transmission errors are
proportional to the concentration error. The user must always keep in
mind that the minimization process does not involve the individual gas
concentration or transmission errors, but rather tries to minimize the
composite rms error of all the errors. Hence, it is vitally important
181
-------
that the effects of the gas concentration errors and the transmission
errors on the rms error are properly balanced.
A2.5.5 Miscellaneous Execution Techniques
To commence its minimization process, program EPAGAS needs to be
furnished with a set of beginning gas constituent concentrations. If
no prior knowledge of the concentrations exists, either the mean or
the center concentration is always a good choice. As mentioned in the
previous section, if the concentration of a particular gas constituent is
approximately known, then this value should be the one selected for the
beginning gas concentration. Also, in this case, the gas error parameters
should be tightened for this particular gas constituent. From prior
experience, improper selection of beginning trial gas concentrations has
generally had no disastrous effects on the final solution except to increase
computation time. In most cases, a unique solution is found regardless
of the starting point. However, in some rare cases, a change in the
starting positions or a change in the order of the gas constituents can
sometimes be helpful in seeking a suitable solution. The size of the
initial increments by which the individual concentrations are varied by
the minimization process must also be carefully selected. If the
increments are too large, the first change in a particular gas concentration
may completely "jump" over the actual true concentration. On the other
hand, if the increments are too small, a proper solution may not be reached
because applying the increment to the original gas concentration results in
an insignificant change. When in doubt, use an increment on the larger
size because the minimization process will decrease the size of the
increment, but will never increase it.
Program EPAGAS is equipped to furnish intermediate printing of
current gas concentrations, of currently computed transmissions, and of
corresponding errors, along with a graph showing the given and currently
computed transmissions as a function of wavenumber. Careful studying of
these intermediate results can often furnish the user with a better feeling
of how certain changes in the individual gas concentrations affect the
overall rms error.
182
-------
The intermediate printing also provides the user with the current
values of selected elements of the KX array, KX(7) through KX(10).
These values also can give the user an insight into the success of each
cycle or try of the minimization process. (See MINMYZD description.)
Last, but not least, there are some sets of experimental data which
do not lend themselves properly to the minimization technique of
determining atmospheric gas constituents and their corresponding
concentrations from their superimposed infrared spectra. This may occur
whenever the data is too noisy or too imprecise, or whenever the scales
(I/IQ or v) must be corrected nonlinearly. Without good data, good results
cannot be achieved, but if the nonlinear corrections are known, they may
be applied in the computer program.
183
-------
APPENDIX 3
TABLE OF CONTENTS
Page
A3. Characteristics of the ROSE Long-Path Spectrophotometer 185
System
184
-------
APPENDIX 3
A3. CHARACTERISTICS OF THE ROSE LONG-PATH SPECTROPHOTOMETER SYSTEM
Design range (transmission) 0.4 to 4.0 km
(0.25 to 2.5 mi.)
Spectral regions 3 to 5,5 y (1820 to 3330 cm" and
7 to 13.5 y (740 to 1430 cm"1)
Resolution to 0.01 y
Scan time 2-min. minimum (either spectral region),
130-min. maximum
(either spectral region)
Display wavenumber, I, I/I (see section 1.1)
Source and Reference Blackbodies
Operating temperature 1800°K
Temperature stability -4°C
Emissivity 0.99 - 0.01
Aperture 12.7 x 6.3 mm (0.50 x 0.25 in.)
Telescopes
Focal length 304.8 cm (120 in.)
Clear aperture 61.0 cm (24 in.)
Field of view 0.25 deg. (-0.125 deg.)
Blur circle diameter 150 y
Primary f/no 2.0
Secondary magnification 2.5
Source Modulator
Operating frequency 570 Hz
Reference Chopper^
Frequency 330 Hz
Stability -0.3 Hz
185
-------
Monochromator
Mfgr. and model
Type
Gratings
Slits
Wavenumber drive
Min. scan time
Detector
Mfgr. and model
Type
Operating temp.
Size
Field of view
D* (Am, 630)
Detector Cooler
Mfgr. and model
Type
Operating temp, range
Operating temperature
Cooldown time
Perkin Elmer (Norwalk, CT)
Model 21OB
Littrow grating with linear wavenumber
drive
3 to 5.5 y; 240 I/mm 20° 3400 cm"1
(2.94 y)
7 to 13.5 y; 101 I/mm 22°2' 1333 cm'1
(7.5 y)
2
masked to 15.0 cm net area (nominal)
width, 0 to 2.0 mm micrometer-
controlled; height, 0 to 12.0 mm
micrometer-controlled
3/32 to 6 rpm in steps of 2x
approximately 2 minutes
Santa Barbara Research Center
(Goleta, CA) 8679'1
Mercury-doped Germanium
26°K (58 psia HJ
0.2 x 2.0 mm (0.008 x 0.08 in.)
90°
0.84 x 1010m Hz^/W
Cryogenic Technology Inc.
(Waltham, MA) Model 20 with
3833-005 temperature controller
closed-cycle He refrigerator,
consisting of compressor and cold
head connected by flexible hoses
approximately 10 to 28 K
26°K (58 psia ty
approximately 15 min.
186
-------
15
10
240 L/mm
Silt Width
6
101 L/mm
Silt Width ->i
i i i
1000 2000 3000
Wavenumber (cm )
4000
600
1000 1400
Wavenumber (cm"1)
1800
Figure A3.1
RESOLUTION VS. WAVENUMBER FOR VARIOUS SLIT WIDTH SETTINGS
FOR ROSE SPECTROPHOTOMETER
-------
APPENDIX 4
TABLE OF CONTENTS
Page
A4. Convenient Conversion Factors 189
188
-------
APPENDIX 4
A4. CONVENIENT CONVERSION FACTORS
U (atm cm) STP X
_2
U (gm cm ) X
U (atm cm) STP X
U (gm cm~2)H20 X
K(v)(atm cm)"1 STP X
KCv)(gm cm"2)"1 X
K(v)(atm cm)'1 STP X
K(v)(gm cm"2)"1 X
K(v)Catm cm)"1 STP X
K(v)(atm cm)"1 STP X
K(v)(gm cm"2)"1 X
1.219 x 10
273
82.06 x ±*
2.689 x 10
3.34 x 1022
19
y
82.06 x ^
M
1.219 x 10
-2
3.72 x 10
A
M 0
356.3 x £
M
4.343
4.343
-20
M
273
gm cm
(atm cm) STP
molecules cm"
molecules cm
(gm cm" )~
(atm cm)" STP
-2 -1
(molecules cm )
-2 -1
(molecules cm )
_2
db/(gm cm )
db/(atm cm) STP
_2
db/(gm cm )
atm cm
1 \
m I
STP
82.06 x ^
M
cm
-1
gm cm
-2
gm cm
1.219 x 10~2 X
M
273
cm
-1
atm cm
STP
atm cm
STP
3.72 x 10
-20
cm
-1
(molecules cm
cm
ppm X
0 = temperature deg. K,
M = molecular weight
M
I
44.96 M
yg m
189
-------
oo
A = Avogadro's number = 6.0225 x 10"
S = Line intensity
K = Absorption coefficient
U = Optical path
190
-------
TECHNICAL REPORT DATA
(Please read Instructions on the reverse before completing)
1. REPORT NO.
EPA-650/2-74-113
3. RECIPIENT'S ACCESSION-NO.
4. TITLE AND SUBTITLE
Remote Sensing of Pollutants: Computerized Reduction
of Long-path Absorption Data
5. REPORT DATE
July 1974
6. PERFORMING ORGANIZATION CODE
7. AUTHOR(S)
V.E. Derr, M.H. Ackley, M.J. Post, and R.F. Calfee
8. PERFORMING ORGANIZATION REPORT NO.
9. PERFORMING ORGANIZATION NAME AND ADDRESS
Environmental Research Laboratories
National Oceanic and Atmospheric Administration
Boulder, Colorado 80302
10. PROGRAM ELEMENT NO.
1AA010; ROAP 26AAP
11. CONTRACT/GRANT NO.
EPA-IAG-077(D)
12. SPONSORING AGENCY NAME AND ADDRESS
National Environmental Research Center
U.S. Environmental Protection Agency
Research Triangle Park, North Carolina 27711
13. TYPE OF REPORT AND PERIOD COVERED
Final
14. SPONSORING AGENCY CODE
15. SUPPLEMENTARY NOTES
16. ABSTRACT
Atmospheric gaseous pollutants are very numerous in industrial regions. It is
estimated that 25 or more pollutant molecules may be found in the atmosphere in
significant quantities. The measurement of the concentration of each gas from the
complex spectrum obtained by a long-path infrared spectrophotometer requires the
fitting of trial spectra composed from a library of spectra. The fitting procedure
adjusts the concentrations of the trial spectra until a "best fit" in a least-
squares sense is produced. This report is a description of the physical, mathemati-
cal, and calculational principles and procedures for the use of a digital computer
program to determine concentrations of atmospheric gases in a path of a few kilo-
meters. Detailed instructions for the computer program and a library of spectra
are provided.
17. KEY WORDS AND DOCUMENT ANALYSIS
a. DESCRIPTORS
Gaseous pollutants
Infrared spectrophotometry
Computer program
13. DISTRIBUTION STATEMENT
Unlimited
b.lDENTIFIERS/OPEN ENDED TERMS .
19. SECURITY CLASS (This Report)
Unclassified
20. SECURITY CLASS (This page)
Unclassified
c. COSATI l-'icld/Group
21. NO. OF PAGES
200
22. PRICE
iPA Form 2220-1 (9-73)
191
------- |