vvEPA
United States
Environmental Protection
Agency
Environmental Sciences Research
Laboratory
Research Triangle Park NC 27711
EPA-600 4-79-068
November
Research and Development
Long-Range
Transport and
Transformation of
SO2 and Sulfate
-------
RESEARCH REPORTING SERIES
Research reports of the Office of Research and Development, U S Environmental
Protection Agency, have been grouped into nine series These nine broad cate-
gories were established to facilitate further development and application of en-
vironmental technology. Elimination of traditional grouping was consciously
planned to foster technology transfer and a maximum interface in related fields
The nine 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 (STAR)
7 Interagency Energy-Environment Research and Development
8. "Special" Reports
9 Miscellaneous Reports
This report has been assigned to the ENVIRONMENTAL MONITORING series
This series describes research conducted to develop new or improved methods
and instrumentation for the identification and quantification of environmental
pollutants at the lowest conceivably significant concentrations It also includes
studies to determine the ambient concentrations of pollutants in the environment
and/or the variance of pollutants as a function of time or meteorological factors.
This document is available to the public through the National Technical Informa-
tion Service, Springfield, Virginia 22161
-------
EPA-600/4-79-068
November 1979
LONG-RANGE TRANSPORT AND TRANSFORMATION OF S02 AND SULFATE
By
Teizi Henmi and Elmar R. Reiter
Department of Atmospheric Science
Colorado State University
Fort Collins, Colorado 80523
Grant No. R 805271
Project Officer
George C. Holzworth
Meteorology and Assessment Division
Environmental Sciences Research Laboratory
Research Triangle Park, NC 27711
ENVIRONMENTAL SCIENCES RESEARCH LABORATORY
OFFICE OF RESEARCH AND DEVELOPMENT
U.S. ENVIRONMENTAL PROTECTION AGENCY
RESEARCH TRIANGLE PARK, NC 27711
-------
DISCLAIMER
This report has been reviewed by the Environmental Science Research
Laboratory, U.S. Environmental Protection Agency, and approved for publica-
tion. Approval does not signify that the contents necessarily reflect the
views and policies of the U.S. Environmental Protection Agency, nor does
mention of trade names or commercial products constitute endorsement or
recommendation for use.
-------
ABSTRACT
This report provides a technical description and program listing of
computational models used in the calculation of long-range transport and
transformation of SOp and sulfate. Two computer schemes have been developed
by which the concentrations of SOp and sulfate along trajectories, as well as
the average concentration distribution over gridded regions, can be calcu-
lated.
In Model A, the effect of the diurnal change of the mixing layer height
is not taken into consideration. Trajectories of mean winds in a single
layer between the surface and the top of the daytime mixing layer are cal-
culated. The thickness of the layer is kept constant. On the other hand, in
Model B, trajectories of the mean winds of the nighttime ground-based stable
layer and of the daytime mixing layer, and trajectories of the mean wind of
the layer between daytime mixing height and nighttime stable layer height are
taken into consideration.
Model B is applied to study the effects of major point sources of SOp on
the air quality over the area between 35°N and 45°N and between 75°W and
95°W, encompassing the Ohio River basin. The calculated concentrations of
SOp and of sulfate are compared with data of observed concentrations. It is
found that the distribution patterns of calculated concentrations of both,
S02 and sulfate are statistically correlated with those of observed concen-
trations.
The observational data on transformation rates of SOp to sulfate are
scrutinized and an empirical formula for the transformation rate as a func-
tion of relative humidity is deduced. Using this formula, calculations are
made of the regional residence times of SOp and sulfate over the region of
the United States east of 105°W longitude, where major industrial sources are
located. The regional residence times of S0? are in the range between 15 and
30 hours for the cold season, and between 15 and 40 hours for the warm sea-
son.
The regional residence time of sulfate is about one order of magnitude
greater than that of SOp, ranging between 150 and 450 hours for the cold
season and between 200 and 500 hours for the warm season.
m
-------
Using a cumulus model, sulfate content in rainwater is calculated. It
is found that aerosol capture by cloud water through microphysical processes
is efficient enough to produce observed levels of sulfate content in rain-
water.
Future research needs in further development of this study are also
described.
This report was submitted in fulfillment of Grant No. R 805271 by Colo-
rado State University under the sponsorship of the U.S. Environmental Pro-
tection Agency. This report covers the period May 1, 1977, to April 30,
1979. Work was completed as of April 30, 1979.
IV
-------
CONTENTS
ABSTRACT iii
FIGURES vi i
TABLES xi
LIST OF SYMBOLS xii
ACKNOWLEDGMENT xv
1. Introduction 1
2. Conclusion 3
3. Recommendations 5
3.1. Refinement of the Model 5
3.2. Verification of the Model 6
3.3 Application of the Model to Study Future Air Quality .... 6
4. Explanation of Models 7
4.1 Calculations of Trajectories, Dispersion Parameter and Other
Parameters 7
4.2 Calculations of the Concentrations of S02 and Sulfate ... 12
4.3 Output of Models 30
4.4 Preliminary Comparison of the Output of Model A and Model B. 31
5. Detailed Structure of the Models 43
5.1 Trajectory Program 43
5.2 Concentration Calculation Program 49
-------
6. Applications of the Model 56
6.1 Introduction 56
6.2 Application (1) 56
6.3 Application (2) 67
6.4 Summary 82
7. Regional Residence Times of S09 and Sulfate Over the Eastern
United States L 83
7.1 Introduction 83
7.2 Data Used 83
7.3 Approach 84
7.4 Results 92
7.5 Sensitivity Analysis 97
7.6 Summary 97
8. Numerical Study on Sulfate Content in Rainwater 99
8.1 Introduction 99
8.2 The Cloud Model 100
8.3 Sulfate Aerosols in the Atmosphere 100
8.4 Results 102
8.5 Conclusion . 104
References 106
Appendices
A. Fortran Listing of Model A 109
B. Fortran Listing of Model B 125
VI
-------
FIGURES
Number Page
1 Schematic diagram showing how trajectories starting at
different times are calculated in Model B. For each K the
upper diagram shows how the horizontal projection of tra-
jectories may split and the lower diagram shows how tra-
jectories appear in vertical layers 8
2 Examples of trajectories originating at 38.4°N and 90.1°W,
starting on May 1, 1974. K=l, 2, 3, and 4 correspond to the
starting times of 00, 06, 12, and 18 GMT, respectively. M=l
correspond to day number. Number(s) at the end of trajectories
are for identification. Where trajectories have split, the
smaller number designates the layer between surface and night-
time stable layer height 10
3 Similar to Figure 2, except trajectories starting on May 2,
1974 11
4 Example of I. along trajectories for M=l, K=l. S.L. indicates
"Single-Layer Model" A; 1, 2, 3 and 4 are trajectory numbers
in Model B as shown in Fig. 2 13
5 Example of Ih along trajectories for M=l, K=2 14
6 Example of I. along trajectories for M=l, K=3 15
7 Example of I, along trajectories for M=l, K=4 16
8 Definition of a, (3, d and t. (For explanation see text.) . . 17
9 Schematic diagram of the models 19
10 Schematic diagram of determination of removal rate and trans-
formation rate along a trajectory 23
11 Example of the concentrations of SOp and sulfate along tra-
jectory 1 for M=l, K=4 24
12 Example of the concentration of SOp and sulfate along tra-
jectory 2 for M=2, K=3 25
vn
-------
Number
13
14
15
16
17
18
Example of the concentration of S0« and sulfate along tra-
jectory 2 for M=2, K=4 .
Definition of angles a, p, a, &, 6-^ and 62; and distances
d and t
0
Distribution of average concentrations of SCL (ug/m).
Model A
0
Distribution of average concentrations of sulfate (ug/m ).
Model A
2
Distribution of deposition amounts of SCL (Kg/Km ). Model A.
2
Distribution of deposition amounts of sulfate (Kg/Km ).
Model A
Page
26
28
32
33
34
35
19 Distribution of average concentrations of SCL at the surface
(ug/m3). Model B 36
20 Distribution of average concentrations of SCL above the
3
nighttime mixing height (ug/m ). Model B 37
21 Distribution of average concentrations of sulfate at the
surface (ug/m3). Model B 38
22 Distribution of average concentrations of sulfate above
the nighttime boundary layer (ug/m ). Model B 39
2
23 Distribution of deposition amounts of S(L (Kg/Km ). Model B. 40
2
24 Distribution of deposition amounts of sulfate (Kg/Km).
Model B 41
25 Flow diagram of trajectory program of Model A 44
26 Flow diagram of trajectory program of Model B 45
27 Flow chart of concentration calculations by Model A 52
28 Flow chart of concentration calculations by Model B 53
29 Distribution of removal velocity (cm/sec) due to precipitation
for warm season (May-October), 1974 57
30 Distribution of average ground level relative humidity
(percent) for warm season (May-October), 1974 58
vm
-------
Number Page
0
31 Locations and emission intensities (x 10 ton/year) of major
S(L point source ...................... 59
o
32 Distribution of average SCL concentrations (ng/m ) at the
surface for the period of May 1 through 10, 1974 ...... 61
33 Distribution of average S02 concentrations (yg/m ) above
the nighttime mixing height for the period of May
1 through 10, 1974 ..................... 62
o
34 Distribution of average sulfate concentration (pg/m ) at the
surface for the period of May 1 through 10, 1974 ...... 63
35 Distribution of average sulfate concentrations
the above nighttime mixing height for the period May 1
through 10, 1974 ...................... 64
2
36 Total deposition amounts (Kg/Km ) of S02 due to dry
deposition and precipitation scavenging for the period
of May 1 through 10, 1974 .................. 65
2
37 Total deposition amounts of sulfate (Kg/Km ) due to dry
deposition and precipitation scavenging for the period of
May 1 through 10, 1974 ................... 66
3
38 Locations and intensities of S02 emission sources (x 10
ton/year) .......................... 68
39 Scheme of calculating 24-hour average concentrations .... 70
3
40 Distribution of S0? concentrations (ng/m at the surface level
on May 11, 1974 71
3
41 Distribution of S0? concentrations (ng/m ) above the night-
time stable layer 72
42 Distribution of sulfate concentrations ((jg/m ) at the surface
level on May 11, 1974 73
3
43 Distribution of sulfate concentrations (ng/m ) above the
nighttime stable layer, on May 11, 1974 74
2
44 Distribution of S02 deposition amounts (Kg/Km ) due to dry
deposition plus precipitation of May 11, 1974 75
2
45 Distribution of sulfate deposition amounts (Kg/Km ) due to
dry deposition plus precipitation on May 11, 1974 76
-------
Number
46
47
48
49
50
51
52
53
F
2
Distribution of S0? deposition amounts (Kg/Km ) due to pre-
cipitation only on May 11, 1974
2
Distribution of sulfate deposition amounts (Kg/Km ) due to
precipitation only on May 11, 1974
Relationship between calculated and observed S02 concentration
Relationship between calculated and observed sulfate concen-
tration
Seasonally averaged relative humidity over the eastern United
States for the cold season
Seasonally averaged relative humidity over the eastern United
States for the warm season
Schematic plot of the concentrations of S0? and sulfate as a
function of time
Transformation rate of S00 to sulfate as a function of rela-
?age
77
78
80
Rl
85
8fi
88
tive humidity 90
54 The regional residence time of SO^ in hours, TSQ , for the
cold season 93
55 The regional residence time of S02 in hours, T<-Q , for the
warm season 94
56 The regional residence time of sulfate in hours, TSM,, for
the cold season 95
57 The regional residence time of sulfate in hours, TSM,, for
the warm season 96
58 Size distributions of typical urban and background aerosols
(Whitby, 1978) 101
59 Sulfate content of rainwater as a function of L.W.C. of the
cloud 103
60 (-)v as a function of L.W.C. of the cloud 105
-------
TABLES
Number
1
2
3
4
5
6
7
8
9
Input parameters for the trajectory program of Model A ...
Output parameters of Program A
Input parameters for the trajectory program of Program B . .
Output parameters for the trajectory program of Model B . . -.
Input data for the concentration calculations of Model A . .
Output parameters of the concentration calculation program
of Model A ...
Input data for the concentration calculations of Model B . .
Output parameters of the concentration calculation program
of Model B
Sensitivity analyses for oarameters
Page
47
48
50
51
54
54
54
55
98
-------
LIST OF SYMBOLS
C-, concentration of S00
1 L
C-.Q initial concentration of S02
C-j, concentration of S02 after time interval t-,
Cy concentration of sulfate
C2Q initial concentration of sulfate
Cp, concentration of sulfate after time interval t-.
€2 maximum concentration of sulfate
d length of trajectory segment
f probability of precipitation
H-, nighttime stable layer height
H2 daytime mixing layer height
H mean mixing layer height
h upper value of integral with respect to height
h bottom value of integral with respect to height
K transformation rate of S02 to sulfate
k , depletion rate of S02 due to precipitation
k 2 depletion rate of sulfate due to precipitation
L mean width of pollutant plume
L width of the original source area
L-, width of pollutant plume at time step I
N number of days multiplied by the number of trajectories per
day
xii
-------
P precipitation rate
P mean precipitation rate
RH relative humidity
R-, quantity representing S02 removal terms
Rp quantity representing sulfate removal terms
T^Q residence time of S02
SDL residence time of sulfate
t distance defined as a function of d and I.
t time
t time interval n
u mean wind speed
Vd vertical velocity through the top of the mixing layer
Vg deposition velocity
Vg-, deposition velocity of SO,,
Vg2 deposition velocity of sulfate
V precipitation scavenging velocity
V x-component of mean wind
/\
V ycomponent of mean wind
X-. distance of trajectory segment at time interval t-^
z distance between grid point and trajectory endpoint
a angle between a latitude circle and a trajectory segment
P angle between distance d and distance t
a angle extending between the line Z and a latitude circle
6-, angle defined as a + n/2
6? angle defined as a + 3n/2
e angle defined as a function of d and I.
xm
-------
6 -, longitude of the endpoint of trajectory segment at time step
n"i n-1
6 longitude of the endpoint of trajectory segment at time step n
K concentration of a pollutant in rainwater
X. depletion rate of S09 due to dry deposition
°1 i
\, depletion rate os sulfate due to dry deposition
d2
\ scavenging rate due to precipitation
I. mean width of pollutant plume
a dispersion parameter
h
a standard deviation of x-component wind
x
a standard deviation of y-component wind
y
i. duration of dry period
T duration of wet period
X concentration of pollutant in air
iL -, latitude of the endpoint of trajectory segment at time step
n"i n-1
4* latitude of the endpoint of trajectory segment at time step n
xiv
-------
ACKNOWLEDGMENT
The trajectory calculation portion of our models in this report is an
expansion and modification of an original scheme developed by Heffter and
Taylor (1975).
Programming help by Ms. Ann M. Starr and Mr. Dan R. Westhoff are greatly
appreciated.
The distributions of 10-day average concentrations ot S0« and sulfate,
were calculated on the EPA computer. Mr. Adrian Busse was very helpful in
running our program on the computer. Thanks are also due to the project
officer, Mr. G.C. Holzworth, for giving us constructive suggestions and
providing necessary data.
xv
-------
CHAPTER 1
INTRODUCTION
Growing reliance upon coal as a fuel for power generation has increased
the international and national concern about S02 and sulfate pollution.
Recent observational evidence in the U.S. and Europe indicates that exposure
to sulfate aerosols and ozone is possible over large regions extending for
hundreds of kilometers downwind from major sources.
In our previous studies under EPA Grant No. R803685, we have developed a
long-range transport model suitable for keeping track of pollutants down-
stream of large industrial complexes (Henmi et al., 1977, 1978). Transport
modelling plays an important role in understanding the relationship between
anthropogenic S02 emission and atmospheric concentrations of SO^ and sulfate,
and serves as a basis for the simulation of present and future changes in
pollutant concentrations and in atmospheric parameters such as visibility,
turbidity and acidity of rainwater.
In order to estimate the fate of pollutants in the atmosphere, the
transport in each of the following layers must be studied:
(1) Nighttime stable layer and daytime mixing layer containing a near-
uniform vertical distribution of pollutants;
(2) The layer between daytime mixing height and nighttime stable-layer
height, in which pollutants are trapped during the nighttime;
(3) The free atmosphere containing pollution which has been carried
up, beyond the range of boundary-layer mixing, by large-scale
ascent or by penetrative convection.
Taking into consideration the transports in the layers described in (1)
and (2), we have developed long-range transport and transformation models of
S0? and sulfate. The details of our models are described in Chapters 4 and
5. Comparisons between the observed concentrations of SOp and sulfate and
those calculated by the model are described in Chapter 6.
In previous studies (Henmi et al., 1977, 1978), we considered the de-
termination of the residence time of S02 in the mixing layer. Assuming that
-------
precipitation, deposition and the transformation into sulfate are the
mechanisms of removal, the regional residence time of SOp was calculated for
the region of the United States east of 105°W latitude. However, the short-
coming of this study was that the fate of sulfate had not been considered.
During the present contract period, we have extended this study to include
the residence time of sulfate. The results are described in Chapter 7.
The concentration of sulfate in precipitation is a consequence of
several microphysical processes occurring within and beneath the clouds. The
microphysical processes include Brownian motion, inertial impaction and
nucleation. In Chapter 8, taking into consideration these processes, a
numerical study is performed on the sulfate concentration in rainwater.
-------
CHAPTER 2
CONCLUSION
One of our principal objectives under Research Grant No. R 805271 has
been to develop computer models for the study of long-range transport and
transformation of S02 and sulfate. These models are intended for use with
input parameters obtainable routinely from the national meteorological ser-
vices. We have developed two computational schemes by which the concen-
trations of S02 and sulfate along trajectories, as well as the average con-
centration distributions over gridded regions, can be calculated. In Model
A, trajectories of the mean wind in a single layer between the surface and
the top of the daytime mixing layer are calculated. In Model B, we take into
consideration the trajectories of the mean wind in (a) the layer between the
surface and the top of the daytime mixing layer, (b) the layer between the
surface and the top of the nighttime stable layer (i.e. the height of the
ground-based inversion layer in which pollutants are essentially "trapped"),
and (c) the layer between the nighttime stable layer and the daytime mixing
height. Because of the more realistic assumptions on the diurnal variation
of the mixing height, Model B is superior to Model A, however it requires
more computer time to run.
The characteristics of these models can be summarized as follows:
1. As input data, wind speed and direction which are observed by ra-
diosonde, and relative humidity and precipitation data observed at
the surface are utilized.
2. The transformation rate of S02 to sulfate is expressed as a func-
tion of relative humidity, and the wet removal rate is parameter-
ized as a function of precipitation rate and precipitation prob-
ability.
3. The program can be applied to any region over the northern hemi-
sphere where input data are available.
Model B was applied to calculate the geographical distributions of
average concentrations of SOp and sulfate over the region between 35°N and
45°N and between 75°W and 95°W which encompasses the Ohio River basin. The
distribution patterns of the calculated concentrations of S0? and of sulfate
were statistically correlated with the patterns of observed concentrations.
-------
The regional residence times of S02 and sulfate were calculated for the
region of the United States east of 105°W longitude. The regional residence
times of S02 were in the range between 15 and 30 hours for the cold season
and in the range between 15 and 40 hours for the warm season. The regional
residence time of sulfate was about one order of magnitude greater than that
of S02, ranging between 150 and 450 hours for the cold season and between 200
and 500 hours for the warm season.
Using a cumulus model, sulfate content in rainwater was calculated.
Although the study has shown that microphysical processes such as condensa-
tion, Brownian diffusion, attachments due to thermophoresis and diffusio-
phoresis can be efficient enough to attain the sulfate content of rainwater
which has been observed, this conclusion should be regarded as tentative.
The combined effects of the microphysical processes and of S02 oxidation to
sulfate in water must be studied in the future.
-------
CHAPTER 3
RECOMMENDATIONS
3.1. REFINEMENT OF THE MODEL
Although our efforts of developing the model have reached the stage
where the distribution of average concentrations, as well as deposition
amounts of S0? and sulfate due to a multitude of sources can be calculated,
there are still several aspects to be improved.
3.1.1. Saving of Computing Time
The model (Model B) in the present version requires, on the EPA computer
(UNIVAC 1100A), about five hours of computing time to calculate 10-day aver-
age concentrations of S02 and sulfate emanating from 60 sources, and the
deposition amounts of these two pollutants over an area of 10 degrees lati-
tude by 20 degrees longitude. The results of the calculations are described
in Chapter 6. In the model, more than 90 percent of the computing time is
spent in interpolation schemes of concentrations along the trajectories into
grid point values. It seems that improvement of this interpolation scheme is
the key to utilizing computer time in a more economic way.
3.1.2. Improvement of the Expression of Transport Process Near the Source
In the present model, the assumption is made that pollutants are in-
stantaneously and uniformly mixed into the mixing layer after being emitted
from the sources. In order to improve the performance of the model, more
realistic formulas to express concentrations, particularly for the first 6 to
12 hours, must be considered.
3.1.3. Use of Daily Radiosonde Data for the Calculation of the Mixing
Height
Climatological averages of the heights of the daytime mixing layer and
of the nighttime stable layer determined by Holzworth (1972) have been used
in the model. For the calculation of the concentration distributions for a
short time period, it is desirable to determine the mixing heights based on
twice daily radiosonde data.
3.1.4. Modification of the Model for Climatological Data
It is possible to calculate long-term (seasonal or annual) average
concentrations, using the present model repeatedly. However, due to the
-------
large amount of computing time required for this process, it is not practical
to do so. Statistical studies of trajectories, dispersion parameters, mixing
height, precipitation and relative humidity over the period of interest must
be carried out. The model will have to be modified to handle statistical
distributions of these parameters.
3.2. VERIFICATION OF THE MODEL
In Chapter 6, we describe our preliminary efforts to verify the model.
Further tests should be conducted using actual concentration data of SOp
and sulfate observed at stations within the area of interest. There exists a
strong correlation between aerosol mass concentrations in the size range of
0.1 to 1 micron in diameter and the extinction of light by scattering
(Charlson, 1969; Hidy et a!., 1976). Furthermore, measurements in and near
St. Louis, MO, have shown that sulfates often account for 50 percent or more
of the aerosol mass in this size range (Charlson et al., 1974). Therefore,
comparisons of the geographical distribution of visibility with the sulfate
distribution derived from our calculations seems a useful means to verify the
model.
3.3. APPLICATION OF THE MODEL TO STUDY FUTURE AIR QUALITY
Once a reliable model performance is assured, the model shall be applied
to study the effects of future emission sources on regional air quality, dry
deposition and acid rain climatology. The model application should focus on
such questions as the impact of dramatically increasing coal use in the
Midwest, the impact of the Ohio River basin emissions on the New England
states, etc. The results of such studies will serve as useful tools in
formulating a national pollution control in the face of the urgently needed
development of new energy sources.
-------
CHAPTER 4
EXPLANATION OF MODELS
Models A and B each consist of two separate programs:
(1) A program that calculates the trajectories of plumes leaving from
any chosen source area(s), dispersion parameters, and the para-
meters necessary for concentration calculations along trajectory
segments using observed wind data.
(2) A program that calculates the concentrations of sulfur dioxide
and sulfate along trajectory segments and interpolates the plume
concentrations to grid intersections.
A program that calculates trajectories was originally developed by
NOAA's Air Resources Laboratories (Heffter and Taylor, 1975) and has been
substantially modified for our purpose. In the following, the details of our
program will be described using examples.
4.1 CALCULATIONS OF TRAJECTORIES, DISPERSION PARAMETER AND OTHER PARAMETERS
4.1.1 Trajectories
In Model A, as mentioned previously, the effect of the diurnal change of
mixing layer height is not taken into consideration. Trajectories of mean
winds in a single layer between the surface and the top of the daytime mixing
layer are calculated. The thickness of the layer is kept constant. On the
other hand, in Model B, trajectories of mean winds of the nighttime stable
layer and daytime mixing layer, and trajectories of the mean winds of the
layer between daytime mixing height and nighttime stable layer height are
taken into consideration. In Figure 1 we show how trajectories of the mean
wind starting at different times are calculated by Model B. The program is
capable of calculating trajectories either four times a day originating at
00, 06, 12, and 18 GMT, or two times a day originating at 00, and 12 GMT. It
is designated in the program that the time between 00 GMT (18 CST) and 12 GMT
(6 CST) is nighttime and the time between 12 GMT (6 CST) and 00 GMT (18 CST)
is daytime. In Fig. 1, H, is the nighttime mixing height and H~ is the day-
time mixing height. Solid lines and broken lines indicate daytime and night-
time trajectory segments, respectively. K = 1 describes the case in which
the air trajectory leaves the origin at 00 GMT (18 CST). In this case we
assume that the pollutants emanating from the surface are contained and
transported in the nighttime mixing layer. At 12 GMT (6 CST), the pollutants
-------
t
\
1
t
t
to
"
v'
.3-
II
*
/
/
l
t
1
1
\
\
\
V
\
\
\
\
\
\
1
t
\
\
\
\
\
t
/
/
/
/
/
/
/
/
/
/
/
/
/
*
/
/
•
II
l
1
1
T
1
? 3
t
1
1
1
if
if
l
1
I
=• c
I1
I
i
r <
N
N
N
O
N
IN
N
0
N
CJ
>~
M
o
N
M
N
O
N
(M
m
o
D
1
\
\
\
\
\ 1
t
u
X.
\
\
(N
it
^
i
i
i
i
\
\
\
\
\
\
V
I
/
/
/
/
r
\
\
\
\
\
\
\
/
i
i
\
\
\
\
\
\
\
\
/
/
i
»
i
\
\
\
\
t
/
t
/
/
/
i
i
i
i
i
i
i
/
/
/
^
/
TT Tl
IttI
It
1
I 1
t
e i c
l * + 1
t t III
1 1
1 1
l
1 i
N - '
r x
-0 ">
QJ O)
+J •!-
ra s-
r- 0
3 -P
c u
_— OJ
r^ •*'
fO '""»
£ U «
2 j_
S- (O r-
0) -P.^
M- C
"- 2 *-
5 '-STS
£ 13^1
g.|a
t*fe
IE^
in -l~>
•» ro
l" 0 S-
0> ^ -P
'£ " >
° E ?
o 2
•^ss
oj "3 5
?3^^=
•P in
5 g|§
^|^>
° o,"3.^
c o> -a
•i— x:
3E-P i.
O CD
M -«=i<£ S
M in o
x: •—
m ro »
^"5
5 «fe-g
^^ m
U -P
•r** * *^~
M -P00--
<3 to CL.
- E i— in
ai m
H J= T3 >,
* U O 03
O co s: E
(U
S-
-------
contained in this nighttime stable layer become mixed uniformly into the
daytime mixing layer due to convective mixing spawned by radiational heating.
At 00 GMT (18 CST) of the second day, the pollutants distributed uniformly
below the daytime mixing height are separated into two layers, the layer be-
tween the surface and the nighttime stable layer height, and the layer be-
tween the nighttime stable layer height and daytime mixing height. The air
masses in the two layers take different tracks. At 12 GMT (6 CST) of the
second day, the pollutants which followed different tracks during the pre-
vious night, are again mixed into the daytime mixing layer. K = 2, 3 and 4
correspond to the cases where trajectories leave at 06 GMT (24 CST), 12 GMT
(6 CST) and 18 GMT (12 CST), respectively. The figure shows that, for two-
day periods, a single trajectory at the beginning separates into four dif-
ferent trajectories.
Figures 2 and 3 are examples of calculated trajectories starting on May
1 and 2, 1974, respectively. The origin of these trajectories is located at
38.4°N and 90.1°W. In these figures, trajectories calculated using Model A
are drawn with thick broken lines. Trajectories drawn with solid lines rep-
resent daytime segments and those with broken lines represent nighttime seg-
ments of trajectories.
4.1.2 Dispersion Parameter
The dispersion parameter av, defined as
av
(1)
is calculated along each trajectory segment. V and V are the x- and y-
x y
components of the mean wind, u is the mean wind speed, and av and av are
defined as
" h
h
0
.
CVx[z] -
h -
V2 <*
h
0
(2)
h
/ (
h
0
Vy[Z] - Vy)2 dZ
[ h'ho
(3)
Here, h and h are the top and the bottom of the layer. Details of deriva-
tion of av, have been described by Henmi et al. (1977, 1978). After the
-------
I O) O)
S- J- £
O CU -T-
U JC -U
IS -P
*
/
'
/
i ''
i
\
, i
X*
•
/
f
'
/
\
-\
*
j
/\
\
\
\\ i
\
/
jr
\
i
m
\:
1
t
tr
UJ
i
i
i
•"'*
i
i F
t- t-
i 1
S;
,.'
i
>
f
!
-^- • ro to
r- >> u
Oi i — -i— (1)
r-H O) M- O
-•r- j_) 4—
rH -P C S-
>> O) T3 ">
tO CX.r-
>- 2: to c
c ^ ^ S
O M— 3t
., f« S-
-P ra S-
S- 00 Ol
tO i— 1 "* >>
-P OJ ro
i/) "O '<~ i —
« ro 9
01 - S- -P
uo -P ro
-a o c
c M- en
ro » o T-
o tn
^ o -a eu
0 c "0
•Sf <•- OJ
. 0 i-
00 n, 0»
CO */) ^- jD
\
ro
»_
^
r
'>
?
^
\
i\
•^ •
;
f
\
^
^
r
^
*
\
i
^
i-x
^»
14
t
T
"*«&
!
?\
•^
S/
<
\
en "J t-
c en
a) • -P
•.-05--,-
j_ +j a> i—
0 JO Q. •
-> Q. > ••-
ro m v_ ro a)
s- eu fa jc jc
-P s- -n
s_ in S-
to «d- o i —
a) -a -P
i— TD C U OJ
a. c o eu i—
ro in ro n3
X - O) S- -P
LU CO S- -P Ul
CNJ
ai
s~
Z!
cn
10
-------
7
f
ro
cvj S-
£
CM
«
CD
C
to
01
s
0
u
01
re
S-
Q.
OJ
U
X
OJ
C\J
OJ
o
-l->
CO
CO
O)
S-
CJ>
11
-------
trajectory analysis has been repeated for n time steps of equal length, the
mean width of the pollutant plume L or Z, along the mean wind is given by
L = Z, = 2(ov. i • t + 0v. 0 • t + + av. „ -t) + L (4)
n h,l h,2 h,n o
where L is the width of the original source area perpendicular to the mean
wind.
Figures 4, 5, 6 and 7 show examples of Zh with L = 0. These Z^ values
were calculated along trajectories exhibited in Figure 2. M in Figure 4 in-
dicates the number of the day (e.g. 1, 2, 3, etc.) of the trajectory calcula-
tion and K indicates the number of the six-hour interval within that day,
during which the trajectory was started. The thick solid lines in these
figures represent Z.-values calculated using Model A. The width of the
plume, Z. , is given in degrees of latitude (1 degree latitude = 111.137 Km).
4.1.3 Other Parameters
Other parameters which are necessary for concentration calculations in
the second program are the angles a and p, the distances d and t, and the
direction indicator of the trajectory, UK. These parameters are defined in
Figure 8.
The distance d is the length of the trajectory segment at time step t .
The distance t is defined as
• + zh(tn)2
where Z. (t ) is the dispersion parameter given by Eq. (4).
a is the angle between a latitude circle and the trajectory segment
during time step t ; p is the angle between d and t.
The direction indicator UK gives the quadrant into which the trajectory
points. East has been chosen as the reference direction. Thus, UK is
defined as
UK
UK
UK
UK
=
=
=
=
1
2
3
4
for
0
90
180
270
o
o
o
o
<
<
<
<
a
a
a
a
<
<
<
<
90
180
270
360
o
o
o
o
4.2 CALCULATIONS OF THE CONCENTRATIONS OF S02 AND SULFATE
4.2.1 Basic Equation
Realizing that our knowledge about the transformation rate of sulfur
dioxide to sulfate at the present does not justify incorporation of anything
more complex than the parameterization of linearized kinetics, the following
two equations for the calculation of concentrations of sulfur dioxide and
sulfate are used:
12
-------
S.L.
NIGHT
TIME(hrs)
Figure 4. Example of I. along trajectories for M = 1, K = 1. S.L. indicates
"Single-Layer Model" A; 1, 2, 3, and 4 are trajectory numbers in
Model B as shown in Fig. 2.
13
-------
NIGHT
TIME(hrs)
NIGHT
Figure 5. Example of I. along trajectories for M = 1, K = 2.
14
-------
DAY
TIMEChra)
Figure 6. Example of I. along trajectories for M = 1, K = 3.
15
-------
Figure 7. Example of 1 along trajectories for M = 1, K = 4.
16
-------
!JK=2
UK =3
IJK =
Figure 8. Definition of a, p, d and t. (For explanation see text.)
17
-------
Let us assume that the pollutants are transported with the mean wind
speed u, as shown in Fig. 9, from an initial vertical plane at location "0"
in which the concentrations C-,/, and C?n are uniformly distributed. C-.Q and
C™ are, respectively, the initial concentrations of sulfur dioxide and sul-
fate. The width of the source area perpendicular to the mean wind is taken
as L , its height is h, equivalent to the height of the mixed layer. After
time t = x/u the pollutants cross a vertical plane P at a distance x, from
the source area. The width of the plume at P is defined as
L(x) = LQ + 2avh • t (5)
where av, is the dispersion parameter.
If C-, and Cp are the concentrations of sulfur dioxide and of sulfate
at distance x, the rate at which pollutants_cross the vertical plane P are
U'h-L(x)'C-| and u*lvL(x)'C?. Assuming that u and h are constant during the
time interval required to travel the distance dx, the net rates of loss
between planes at x and x + dx are given by u-h'-r- (LC-,)-dx and U'h-T-
(LC2)«dx. These amounts are balanced by removal processes such as dry depo-
sition, precipitation scavenging and loss from the top of the mixed layer,
and by transformation processes.
The mass balance equations for sulfur dioxide C-, and sulfate C2 are
given by
(6)
(7)
where V is the deposition velocity, V . is the vertical velocity through the
top of the mixing layer, h is the height of the mixing layer, f is the
probability of precipitation during the time x/u, V is the precipitation
scavenging velocity and K is the transformation rate of sulfur dioxide to
sulfate. Subscripts I and 2 stand for sulfur dioxide and sulfate, respec-
o =
tively. The ratio •*- is the ratio of molecular weights of sulfate (SO,) to
sulfur dioxide (S0?).
The solutions to Eqs. (6) and (7) are given by
Rl . xl (8)
= -[v '
- -tv <
h fVsl + Vdll •
• LCj - KCj • Lh
• LC2 + | KCj • Lh
L + 2avh
u
18
-------
TOP OF MIXING LAYER
Figure 9. Schematic diagram of the models.
19
-------
C2 = -*X-» • exp
L LQ + 2avh *
il
exp -V* 9-1 (9)
where R^ = V - + fV$1 + Vdl + K-h, and R2 = V 2 + fVg2 + Vd2. Equations (8)
and (9) can be applied along a trajectory as follows: After the first inter-
val, t-., of trajectory analysis, the concentrations of the pollutants at the
endpoint of a trajectory segment, x-., are calculated as C-.-. and C^-i and are
regarded as the concentrations of the new hypothetical line source with the
width of !_-, = L +2 v.-, t-,. By repeating this calculation process for
successive trajectory elements, the concentrations along the trajectory at
any distance from the source area can be calculated.
4.2.2 Removal Mechanisms and Transformation of S02 to Sulfate
4.2.2a Removal Due to Precipitation
The removal of pollutants due to precipitation can be expressed in terms
of scavenging velocity, V (Henmi et al., 1977, 1978). The parameter is
defined as
V = (-) • P (10)
x v
where K is the concentration of a pollutant in rainwater, x the concentration
of a pollutant in air, and P the precipitation rate during a precipitation
i/
event. The subscript v means that the ratio (-) is formed on a volume basis.
Since precipitation is episodic, the probability of a precipitation event f
must be multiplied to V . The probability of a precipitation event is de-
fined as
f =
where, t. and t are the duration of dry and wet periods, respectively.
In application (1) described in Chapter 6, and in the calculation of
regional residence times of S02 and sulfate discussed in Chapter 7, the pre-
cipitation rate P and the probability of a precipitation event f based on
climatological data are used. On the other hand, in application (2) in
20
-------
Chapter 6 where 24-hour average concentrations of S02 and sulfate are
calculated, P and f are computed for every 12 hours at 81 weather stations
and are interpolated for individual trajectory segments.
The actual values of the ratio (K/X) for S0~ and sulfate are listed in
Chapters 6 and 7.
4.2.2b Removal Due to Dry Deposition
Removal of SCL due to dry deposition has been reviewed in detail in our
previous reports (Henmi et al., 1977 and 1978). Dry deposition of sulfate
has not been studied extensively. The only known measurements of sulfate
deposition velocity in the atmosphere yielded approximately 0.1 cm/sec
(Engelman and Sehmel, 1976).
Realizing the difference of atmospheric conditions between daytime and
nighttime, different values of dry deposition velocities are assigned for
daytime and nighttime. Actual values of dry deposition velocities used for
the concentration calculations of S00 and sulfate are described in Chapter
6. l
4.2.2c Removal Due to Vertical Transport of Pollutants into the Free Atmo-
sphere
Pollutants released into the mixing layer are dispersed into the free
atmosphere by several processes, such as turbulent diffusion, convective
activity, large-scale ascent from horizontal convergence and upgliding of
warm air over cold air. The relative effectiveness of these processes de-
pends on meteorological conditions. With improved knowledge of these pro-
cesses we hope to be able to parameterize their impact on vertical mixing by
chosing appropriate values for V,, the vertical velocity through the top of
the mixing layer. However, in the present models, we assume that V, = 0.
4.2.2d Transformation of S00 to Sulfate
The following empirical relationship between the transformation rate and
relative humidity is obtained:
K = 3.304xlO~4-exp (0.063-RH) hr"1 (12)
K is the transformation rate of S02 to sulfate and RH is the relative humid-
ity in percent. The derivation of this equation is described in Chapter 7.
The removal rate due to precipitation scavenging and the transformation
rate of SO- to sulfate along trajectory segments are interpolated using the
21
-------
data of meteorological stations which are surrounding the midpoint of a
trajectory segment. An interpolation routine from NCAR's Software Support
Library was chosen and adapted in the model. The interpolation routine along
with its attendant subroutines uses a given set of points and their coordi-
nates (latitude and longitude) to set up a grid of triangles using those
points as vertices (see Fig. 10). A linear interpolation from the triangle
sides surrounding the midpoint of a trajectory segment is performed. If
points outside the grid, (x , y ), must be interpolated, this routine will
extrapolate from the edge of the grid to those points.
4.2.4 Interpretation of Plume Concentrations to Grid Intersections
Models A and B are capable of calculating the average concentrations of
sulfur dioxide and sulfate downstream from source areas over gridded regions.
The procedure can be summarized as follows:
(i) The concentrations along trajectories are calculated using Eqs.
(8) and (9). Latitudes and longitudes of endpoints of trajectory
segments, and dispersion parameters are used as input data. Ex-
amples of the concentrations of SCL and sulfate along trajectories
are shown in Figs. 11 to 13. In these figures, the source con-
q
centrations of S00 and of sulfate are assumed to be 30 ug/m and 0
3
ug/m , respectively. Concentrations were calculated by Model B.
They are different along different trajectories due to different
values of dispersion parameters and of removal terms.
(ii) Each grid point over the region is examined as to whether or not
it is located within the plume. If a grid point is located outside
the plume, the concentration at that point is assigned the value 0.
If a grid point is located within the plume, the concentrations
based on Eqs. (8) and (9) are assigned to that point.
The criteria of determining whether or not grid points are located
within the plume are discussed in the next section.
(iii) For each trajectory segment, step (ii) is repeated, and the values
of concentrations at affected grid points are accumulated.
(iv) The calculations of concentrations for trajectories from different
source areas are carried out repeating steps (i) and (ii). In
Model B, the calculations of concentrations for trajectories of the
winds in different layers are also conducted repeating steps (i)
and (ii). The concentrations due to different sources and due to
different trajectories are accumulated for each grid point.
(v) After step (iv) is completed, the average concentrations at grid
points are calculated dividing the accumulated concentrations by
N. Here N is defined as the number of days multiplied by the
number of trajectories per day.
22
-------
Figure 10. Schematic diagram of determination of removal rate and transfor-
mation rate along a trajectory.
23
-------
(cW/0Tl)ONOO
s-
o
o
4J
u
0)
•f— 3
to
s_
O)
E
O
'ro
O)
-a
c.
re
CM
o
to
in
c.
o
(O
s_
OJ
o
c
o
o
0)
-C
-l->
M-
O
Ol
'o.
tt3
X
ai
s-
3
O5
-------
UW/9rl)'ONOO
CM
S_
o
CM
o
-M
(J
O)
•I-J
ro
-p
O5
O)
-P
3
IO
T3
C.
to
CM
O
o
c
o
(0
S-
HI
o
C
o
(J
Ol
CD
'o.
to
X
CM
OJ
S-
01
25
-------
CM
S-
o
CM
u
o>
'T~>
03
Ol
-p
ro
T3
C
03
o
to
-t-J
ro
Ol
u
o
u
O
O)
Q.
(0
X
ro
01
26
-------
4.2.5 Criteria of Determining Locations of Grid Points Relative to the
Plume
Referring to Fig. 14, the following parameters are defined:
8n-i> ^n-r longitude and latitude of the endpoint of trajectory
segment at time step n-1.
6 , ib : longitude and latitude of the endpoint of trajectory
segment at time step n.
d: the distance between (0 -,, »b -,) and (0 , ib ) and de-
fined as n~l "~l n n
" + y in nautical miles
where
= (6n-l - 6n) ' cos (Vl ' n/180) ' 60 (13)
y = (tb - ib -,) • 60.
J VYn Tn-l
t: the distance defined as
71* (15)
z: the distance between the grid point and the point (6n_i, ^n-r
a: the angle defined as
a = tan"1 £ (16)
P: the angle defined as
P = tan"1 ^ (17)
y: the angle extending between the line z and a latitude circle.
s: the angle defined as
e = tan"1 f ™
*h
5-,: the angle defined as
dl = a + n/2 (19)
52: the angle defined as
62 = a + 27t (2Q)
27
-------
3aniiivi
c.
to
in
O
00.
01
-------
With these parameter definitions, the criteria for determining whether or not
grid points are located within the trapezoidal plume segment are given such
that z and y must satisfy the following conditions:
^n - Vl and 6n - Vl ^IJK = ^
A. p > a
1. a < y < a + p 2. In + (a-p) < y < 2n
d d
z - cos (y-a) z - cos (ya)
3. 0 < a < a H. or + p < -^ «, 6-,
. , d z < . cos
_^
- cos (y-a) - cos (crye)
5. 62 < y < 2n + (a-p)
a < t cos j£+P) ,
cos (y-cre)
B. p < or
1. a-p i);n_1 and 6n > 6n_1 (UK = 2) and
[III] »|»n < ^^ and 6n < 6n.1 (UK = 3)
1. a-p 2rt
1. a-p
-------
3. 6^^ > y a + p - 2n
< . cos (e+p)
- cos (ors-y)
B. a + p < 2n
4. 69 < y < a - p
2
t cos (s+p)
"
<
- cos
2. 69 < y < a - p
z
cos (s+p)
3.
- cos (a-y)
a + p < y < 2n
cos (s+p)
z < t
- cos
4.3 OUTPUT OF MODELS
J2-
z - cos (y-cr£)
4. 0 < y < 6-^
7 < t cos JL+P)
- cos (cry-£)
Output parameters calculated for each grid point by the models are as
follows:
Program A
1. Average concentration of
S02 (|jg/m3)
2. Average concentration of
sulfate (|jg/m)
Program B
1. Average concentration of SQ0 at the
3
surface (fjg/m )
2. Average concentration of S02 above
the nighttime stable layer height
3. Deposition amount of SO,, 3.
(Kg/Km2)
4. Deposition amount of sulfate 4.
(Kg/Km2)
7.
8.
Average concentration of sulfate at
o
the surface (pg/m )
Average concentration of sulfate
above the nighttime stable layer
2
height (ng/m )
Total deposition amount of S09
2
(Kg/Km) due to dry deposition and
precipitation scavenging
Total deposition amount of sulfate
2
(Kg/Km ) due to dry deposition and
precipitation scavenging
Deposition amount of S09 due to pre-
2
cipitation (Kg/Km )
Deposition amount of sulfate due to
2
precipitation (Kg/Km)
30
-------
Here, the deposition amounts of SO,, and sulfate on the surface are calculated
from (V i+fV i)xC-, and (V 2+fV 2)xC2- *n Model B, the deposition amounts due
to precipitation only are also calculated.
4.4 PRELIMINARY COMPARISON OF THE OUTPUT OF MODEL A AND MODEL B
In order to compare the performances of Models A and B, the follow-
ing preliminary calculations were conducted. The source was assumed to be
located at 38.4°N and 90.1°W. The region between 30°N and 45°N and between
90°W and 105°W was considered. The source concentrations of S09 and sulfate
33
were assumed as 30 ug/m and 0 ug/m , respectively.
It was also assumed that no precipitation fell over the area. Realizing
the difference of atmospheric conditions between daytime and nighttime,
different values of removal terms were assigned for these times. The follow-
ing values were used for this calculation.
Daytime Nighttime
V -, 2 cm/sec 1 cm/sec
V - 0-2 cm/sec 0.1 cm/sec
K 6xlO"6 sec"1 3xlO"6 sec"1
Vs-^ 0 cm/sec 0 cm/sec
V 2 0 cm/sec 0 cm/sec
V.-, 0 cm/sec 0 cm/sec
V >2 0 cm/sec 0 cm/sec
The calculation of trajectories and dispersion parameters were made
using observed wind data supplied by NOAA's Air Resources Laboratory at
Silver Spring, Maryland. The calculations were conducted for two days, from
May 1, 1974 to May 2, 1974. Trajectories were initiated every 6 hours. The
daytime mixing height H? was assumed to be 2000 m, and the nighttime mixing
height H^ was assumed to be 500 m. For this test-run the background concen-
trations of both S0? and sulfate in the region under consideration were
assumed to be zero.
Trajectories and dispersion parameters for this example have been shown
in Figs. 2 and 3 and Figs. 4 to 7, respectively.
The output of Model A is shown in Figs. 15 to 18. Figures 19 to 24
give the output of Model B. It can be seen from a comparison of these
figures that S02 concentrations obtained from both Model A and Model B are
31
-------
5 ; • • ; -I • •
-. \ \ * =. Y ' '
?•?••/!•=
O>
-o
o
CM
O
CO
m
s_
-P
c
O)
u
c
o
(J
O)
O5
(C
s-
OJ
>
to
o
c
o
-t->
3
si
o
Lf)
HI
S-
3
D5
32
-------
• i
OJ
T3
O
O)
IL
OJ
+J
(O
O
1/5
O
03
S-
-t->
c
cu
u
c
O
u
o>
D5
03
s_
QJ
ro
O
-p
01
s_
D5
M } M M M M ; i 3 i i ; ;• : ; z 2
33
-------
01
TD
O
CNJ
cn
CM
O
to
o
in
-P
c:
3
o
re
c.
o
o
Q.
O)
C.
o
in
o
-------
O)
-o
o
CXJ
CD
CD
-P
(0
o
ra
c
a
in
O
Q.
c
o
S-
-p
to
00
o>
O5
35
-------
CO
i i
»'
.
o>
TD
O
O)
U
s-
3
in
(0
O
CO
in
c
O
-p
m
s_
0)
U
c:
o
O
0)
O5
(O
S-
d)
>
(D
o
-p
S-
in
•i—
a
en
O)
s-
3
O5
-------
£
on
QJ
_c
01
c
x
E
Ol
O5
C
O)
^:
-t->
OJ
>
o
XI
to
CM
O
I/)
c
o
to
S-
01
u
c
o
(J
QJ
CD
tO
s_
CD
>
to
c:
o
s_
Q O
O
CSJ
O)
S-
02
37
-------
01
T3
O
tS»tttttlt
"•••*.? ? ' Xl
tr>
01
(J
(0
0)
+J
03
in
c
o
ro
S-
O)
(J
o
u
o>
O5
(O
S-
o>
>
(XJ
_Q
si
-l->
to
o
H H M M
CSJ
01
38
-------
O)
S-
QJ
i i i ;
03
-D
C
o
_a
01
en
c.
0)
£:
-(->
O)
o
JD
a
to
C
o
(T3
O)
O
C
o
u
ai
D5
03
ai
o
c
o
M H H
CQ
S_ T3
-*-> O
O
CM
CSJ
a>
O)
39
-------
t:;;tt«::«:
» « » ! ! i * • • S !
CQ
0)
-a
o
05
CO
o
3
O
1
o
+3
a
o.
01
T3
O
c
o
•r—
•4->
3
s
+J
l/l
CO
CM
01
S-
40
-------
O)
-a
o
05
O)
-fJ
fO
o
in
c
o
rO
O
in
O
a.
01
TJ
c
o
3
_a
s-
CsJ
O)
CD
41
-------
very similar. However, the sulfate distribution calculated by Model A shows
different features from that calculated by Model B. Because of the prelimi-
nary nature of these test-runs, conclusions from the above results will not
be drawn.
The computing time required by Model B is about three times longer
than that required by Model A. Interpolation schemes of concentrations
along trajectory into grid points use more than 90 percent of the computing
time. The computing time is proportional to the number of grid points,
number of sources, number of trajectories and number of days for which
calculations are carried out.
42
-------
CHAPTER 5
DETAILED STRUCTURE OF THE MODELS
In order to calculate the distributions of SOp and sulfate concentra-
tions due to multiple sources, Model A as well as Model B are divided into
two major routines, a trajectory program and a concentration program. The
present versions of either Model A or Model B can calculate the concentra-
tions due to five sources at a time. The concentrations due to more than
five sources can be calculated by applying the model repeatedly.
Both models can be used to calculate the average concentrations of S02
and sulfate for time periods longer than 24 hours.
5.1 TRAJECTORY PROGRAM
5.1.1 Flow Charts
Flow charts of the trajectory calculation routines for Model A and Model
B are shown in Figs. 25 and 26.
5.1.2 Description
The function of the trajectory program is to use upper-air multi-level
wind speed and direction observations from a randomly distributed network of
reporting stations to compute trajectories, dispersion parameters and other
parameters along the trajectories within layers of specified thickness.
The major input parameters to the program are latitude and longitude
values of stations, wind direction and speed data as a function of height
from rawinsonde observations. These data are on magnetic tapes.
The major output consists of latitude and longitude values of the end-
points of trajectory segments, mean wind speeds, dispersion parameters and
other parameters which have been described in the previous chapter.
Model A
The program of trajectory calculations consist of a main program and
nine subroutines, plus seven subroutines which are used in the transmittal of
information from the data.
43
-------
O5
s_
o
X
O)
0.
O)
4J
X
QJ
Read governing
parameters from
card
Print governing
parameter values
Calculate
u", v", ov , av
A y
for each radiosonde
station
observed
wind data
Calculate
• latitude and longitude of endpoint of trajectory
segment
• mean wind and dispersion parameter along trajectory
• other parameters required in the second program
Print above parameters
in file
Stop
Figure 25. Flow diagram of trajectory program of Model A.
44
-------
i.
o
X
-------
TRAJET is the main program in which the loops of number of days, number
of origins, number of trajectories per day, and number of trajectory segments
are defined. The output parameters are written in this program.
INPUT is the subroutine to read control parameters from input cards and
to define other parameters. The input parameters for the program are given
in Table 1.
DTABKO is the subroutine to make data blocks of averaged wind and dis-
persion parameters, calling RDAVQ in which wind data are read in, and x- and
y-components of the mean wind are calculated for each station.
AWIND is called from RDAVO to average winds and to calculate dispersion
parameters for each station.
DISQ is the subroutine in which latitude and longitude of trajectory
segments, dispersion parameters, mean wind speeds and other parameters are
calculated using the data block calculated in DTABKO. Subroutine ITSIW is
called from DISO and calculates the weighting factors and the contributions
of individual data to a trajectory segment calculation. Subroutine XSZ is
to calculate angles a and p, and distance parameters t and d which are used
in the second program of concentration calculations.
KKMM is the subroutine to determine proper trajectory segments, and time
and day in a data block.
ALTDTA provides alternate data for the output of the program when
outputs are not calculated because of unavailability of input data.
CLSTM is a subroutine to provide an alternate set of input data within
the closest synoptic time interval for the calculation of output parameters.
Subroutines BOUT, POSTP, RDDATE, RDSTHD RDWIND, RDDATA and ATAPE are
used to read the data from the input tapes.Output parameters of the trajec-
tory program of Model A are given in Table 2. These parameters for different
sources (origin) are stored in different files, since they are used as input
data for the concentration computation program.
Model B
The structure of this program is similar to that of Model A. The major
differences of Model B from Model A are:
(i) A temporary file to store wind data for repeated use is added.
(ii) Subroutines DIS02 and DIS03 are added. The roles of these sub-
routines are similar to those of DISO, but DIS02 is for the wind of
the layer between the surface and the nighttime stable layer
height, and DIS03 is for the wind of the layer between the night-
time stable layer height and the daytime mixing height.
46
-------
TABLE 1. INPUT PARAMETERS FOR THE TRAJECTORY PROGRAM OF MODEL A
Card No.
1
2*6
Parameter
Names
NO
NDY
OID
OLAT
OLON
Definition
Number of origins
Number of days for which computations are
desired
Name of origins (sources)
Latitude of origin (sources)
Longitude of origin (sources)
7*11
-12
13
14
15
16
17
18
WIOROT
IBDY
IMO
IBYR
IHOUR
IBMO
IDIR
NDYDUR
NDYDTA
WTYPE
Width of sources
Day \
Month ( of the beginning of
Year I
Hour ) computation
Month in logical card
Direction of trajectory, forward (FORW)
Duration in days
Number of days of wind input data
Type of wind input data, observed (OBS)
of boundaries for
observed wind input
data
of transport layer in meters
boundaries for map
DTABT
DTABB
DTABL
DTABR
LBAAT
LTAAT
ALATT
ALATB
ALONL
Top latitude
Bottom latitude
Left longitude
Right longitude
Base | ,
Top I of
Top )
Bottom > b
Left )
47
-------
Parameter
Names
TSLAT
TSLON
UBAR
SPREAD
ALPHA
BETA
T
R2
UK
ICNT
SICNT
TABLE 2. OUTPUT PARAMETERS OF PROGRAM A
Definition
Latitude of the endpoint of trajectory segment
Longitude of the endpoint of trajectory segment
Mean wind speed along trajectory segment
Dispersion parameter
Angle a defined in Figure 8
Angle p defined in Figure 8
Distance t defined in Figure 8
Distance d defined in Figure 8
Direction indicator
Parameters necessary for mapping
48
-------
(iii) Subroutine CONVER is added in which output parameters are re-
placed by differently-named parameters.
(iv) The main program TRAJET becomes longer to calculate simultaneously
four trajectories and parameters along these trajectories.
Governing input parameters for Model B are read in a subroutine
INPUT and are shown in Table 3.
Output parameters of Model B are given in Table 4. Since trajectories
of mean winds in three different layers are computed, there are four trajec-
tories for two-day periods. Output parameters for each trajectory are stored
separately in a different file.
5.2 CONCENTRATION CALCULATION PROGRAM
5.2.1 Flow Charts
Flow charts of the concentration calculation programs of Models A and
B are shown in Figs. 27 and 28.
5.2.2 Description
The function of the concentration calculation program is to use input
data generated in the trajectory program for the computation of average
concentrations and deposition amounts of S0? and sulfate over gridded areas.
Model A
The program consists of a main program and two subroutines.
In the main program TRAJET, governing parameters which have been given
originally in the trajectory program are transferred into the concentration
calculation program. Input data generated in the trajectory program are also
read from files in the main program.
Subroutine CONCAL is the major part of the concentration calculation
program, in which"tfie concentration calculations along trajectories and
interpolations of the plume concentrations to grid points are made.
terms
polations of the plume concentrations to grid points are made.
Subroutine REMOVE is called from CONCAL to give the values of removal
according to the number of loops which represent the time of day.
Input data which must be read in from cards are source concentrations of
S0? and of sulfate, and the width of the source. These data and are shown
in Table 5.
Model B
In addition to the main program and two subroutines used in Model A,
Model B calls upon the subroutines DETCON, PRECIP, TRIANG, CONHUL, BADGET,
TMESH, TMESHI, SORTI, PVEC, and STEST, and the functions TRINTR and LTEST.
49
-------
TABLE 3. INPUT PARAMETERS FOR THE TRAJECTORY PROGRAM OF PROGRAM B
Card No.
1
2^6
4
5
6
7
8
9
10
Parameter
Names
NO
NDY
OID
OLAT
OLON
WIORT
ERS00
L
ERSUL
IBDY
IMO
IBYR
IHOUR
IBMO
TDIR
NDYDUR
NDYDTA
WTYPE
DTABT
DTABB
DTABL
DTABR
LBAAT
LTAAT
LMAAT
ALATT
ALATB
ALONL
Definition
Number of sources
Number of days for which computations are
desired
Name \
[»
Width )
Emission rate
of S02
Emission rate from sources
of sulfate
Day \
Month ( of the beginning of
Year ( computation
Hour J
Month in logical cord
Direction of trajectory, forward (FORW)
Duration in days
Number of days of wind input data
Type of wind input data, observed (OBS)
Top latitude \
Bottom latitude I uounaaries tor
Left longitude ( observed wind i
Right longitude )
Bottom of mixing height
Daytime mixing height
Nighttime mixing height
Top )
Bottom \ boundaries for map
Left )
nput
50
-------
TABLE 4. OUTPUT PARAMETERS FOR THE TRAJECTORY PROGRAM OF MODEL B
Parameter
Names
TTLAT1
TTLON1
UBAR1
SPREA1
ALPHA!
BETA1
Tl
R21
IJK1
TTLAT2
TTLON2
UBAR2
SPREA2
ALPHA2
BETA2
T2
R22
IJK2
TTLAT3
TTLON3
UBAR3
SPREA3
ALPHAS
BETA3
T3
R23
IJK3
TTLAT4
TTLON4
UBAR4
SPREA4
ALPHA4
BETA4
T4
R24
IJK4
Definition
Latitude of endpoint of trajectory 1
Longitude of endpoint of trajectory 1
Mean wind speed along trajectory 1
Dispersion parameter along trajectory 1
Angle a along trajectory 1, defined in Fig. 8
Angle p along trajectory 1, defined in Fig. 8
Distance t along trajectory 1, defined in Fig. 8
Distance d along trajectory 1, defined in Fig. 8
Direction indicator along trajectory 1
Latitude of endpoint of trajectory 2
Longitude of endpoint of trajectory 2
Mean wind speed along trajectory 2
Dispersion parameter along trajectory 2
Angle a along trajectory 2
Angle p along trajectory 2
Distance t along trajectory 2
Distance d along trajectory 2
Direction indicator along trajectory 2
Latitude of endpoint of trajectory 3
Longitude of endpoint of trajectory 3
Mean wind speed along trajectory 3
Dispersion parameter along trajectory 3
Angle a along trajectory 3
Angle p along trajectory 3
Distance t along trajectory 3
Distance d along trajectory 3
Direction indicator along trajectory 3
Latitude of endpoint of trajectory 4
Longitude of endpoint of trajectory 4
Mean wind speed along trajectory 4
Dispersion parameter along trajectory 4
Angle a along trajectory 4
Angle p along trajectory 4
Distance t along trajectory 4
Distance d along trajectory 4
Direction indicator along trajectory 4
51
-------
O)
u
3
O
GO
X
O)
Read governing parameters
Read input data generated
in trajectory program
Calculate concentration
of S02 and sulfate
along trajectory
i
Writo input data
Calculate the concentrations
and deposition amount
on grid point
Accumulate the concentrations
and deposition amount
Calculate the average
concentrations due to present
number of sources
Source concentrations
of S02 and sulfate
and source width
Write the average concentration
and deposition amounts
due to present number of sources
Stop
Figure 27. Flow chart of concentration calculations by Model A.
52
-------
O
4J
O
Read governing parameters
File of governing
parameters
Read input data generated
in trajecotry program
File of input data
Calculate concentration of
SOp and sulfate along tra-
jectory
* .
Precipitation rate
and relative humidity
data
I
Write concentrations along trajectory
Calculate the concentrations and
deposition amounts on grid points
I
Accumulate the concentrations and
deposition amounts in file
Stop
Read file
Calculate the average concentrations
Write the average concentrations and
deposition amounts
Stop
Figure 28. Flow chart of concentration calculations by Model B.
53
-------
TABLE 5. INPUT DATA FOR THE CONCENTRATION CALCULATIONS OF MODEL A
Parameter
Name Definition
ORCON Source concentration of SOp
SORCON Source concentration of sulfate
WIORO Width of source area
Output parameters of this program are shown in Table 6. These out-
put parameters are written for every grid point.
TABLE 6. OUTPUT PARAMETERS OF THE CONCENTRATION CALCULATIONS OF MODEL A
Parameter
Name Definition
o
XXJ Average concentration of S00 (ug/m )
3
SXXJ Average concentration of sulfate (|jg/m )
DXXJ Deposition amount of SO, (Kg/Km2)
2
DSXXJ Deposition amount of sulfate (Kg/Km )
The role of the subroutine DETCON is to redefine the concentration value
according to the number of trajectories and the number of loops. The
subroutine PRECIP calls upon TRIANG and its attendant routines to interpolate
precipitation and relative humidity data for the midpoints of trajectory
segments. CONHUL, TMESH, TMESHI, SORTRI. PVEC, and STEST, and the functions
TRINTR and LTEST are the attendant routines to TRIANG. The subroutine
BADGET calculates the sulfur budget in the area of interest.
Input data read from cards are given in Table 7.
TABLE 7. INPUT DATA FOR THE CONCENTRATION CALCULATIONS OF MODEL B
Parameter
Name Definition
1ST Station number
LAT Latitude of station
LON Longitude of station
RH Relative humidity at stations
PS Precipitation rate at stations
Output parameters of Model B are shown in Table 8. In addition to the
various parameters at grid points, the concentrations along trajectories are
given as output.
54
-------
TABLE 8. OUTPUT PARAMETERS OF THE CONCENTRATION CALCULATION PROGRAM OF
MODEL B
Parameter
Name Definition
I Number of trajectory
TSLAT Latitude of trajectory endpoint
TSLON Longitude of trajectory endpoint
CONC Concentration of S02 along trajectory
SCONC Concentration of sulfate along trajectory
M Subscript for the loop of day
K Subscript for the loop of trajectory per day
J Subscript for the loop of trajectory segment
XXJ Average surface concentration of SO,, (ug/m )
XXJ2 Average concentration of S09 above the nighttime
3
mixing height (ug/m ) 3
SXXJ Average surface concentration of sulfate (ug/m )
SXXJ2 Average concentration of sulfate above the nighttime
o
mixing height (ug/m ) 9
DXXJ Deposition amount of S09 (Kg/Km)
2
DSXXJ Deposition amount of sulfate (Kg/Km ) ,,
PXXJ Precipitation deposition amount of S09 (Kg/Km)
2
PSXXJ Precipitation deposition amount of sulfate (Kg/Km )
FORTRAN programs of Model A and Model B are given in Appendix A
and B, respectively.
55
-------
CHAPTER 6
APPLICATIONS OF THE MODEL
6.1 INTRODUCTION
In this chapter we describe the results of model applications. From
Model B we calculated the geographical distributions of average concentra-
tions of S02 and of sulfate over the region between 35°N and 45°N and between
75°W and 95°W, which encompasses the Ohio River basin. Two applications were
made. In application (1), 10-day averaged concentrations and deposition
amounts of SOp and sulfate were calculated. In application (2), 24-hour
averaged concentrations and deposition amounts over the region were calcu-
lated. The computed concentrations at the surface level were compared with
the concentrations observed at stations located in the region.
6.2 APPLICATION (1)
6.2.1 Input Data
In this application we calculated 10-day averaged concentrations
for the period from May 1 to May 10, 1974. Averaged precipitation scavenging
velocities and average relative humidities based on climatological data were
used to calculate removal rates of SOp and sulfate and transformation rates
of S02 to sulfate. Figures 29 and 30 were obtained by interpolating the
data which were observed at sixty weather stations in the eastern United
States. The data were used for the period between May and October, 1974.
c
Sixty point sources of S02 whose emission rate is more than 10 ton/year
were taken into consideration. The geographical locations and the emission
rates are shown in Fig. 31 (Clark, 1978). Other input data used for the
calculations, in addition to precipitation scavenging velocity, relative
humidity and emission rate, are as follows:
Daytime mixing height 1600 m
Nighttime stable layer height 250 m
Dry deposition velocity of SOp
daytime 2 cm/sec
nighttime 1 cm/sec
56
-------
£-
O)
-Q
O
+J
(J
O
I
rO
O)
in
S-
rd
S-
O
c
O
-p
fO
-p
u
HI
J-
Q.
O)
u
o>
to
(J
o
o
o
TO
>
O
E
O)
S-
W) I—
••- cn
Q i-l
CSJ
ID
S-
cn
57
-------
c:
o
to
ro
u
S-
O)
O-
01
>
OJ
O)
>
Ol
o
S-
D)
O)
CD
CD
O)
>
rd
3
-Q
to o
StJ
o
o
OO
O)
S-
O5
58
-------
1-
CO
GO
CO
00
vS
-rO
X
CM
X
V \ CM ~
JSV> i e
l
-------
Dry deposition velocity of sulfate
daytime 0.2 cm/sec
nighttime 0.1 cm/sec
The ratio (K/X) for S09 5xl04
4
The ratio (K/X) for sulfate 5x10
(K/X) are the ratios of the concentration of S0? or sulfate in rainwater
to the concentrations of those pollutants in air. The heights of the day-
time mixing layer and of the nighttime stable layer were based on data by
Holzworth (1972).
Calculations were conducted using wind data provided by the Air Re-
sources Laboratory, NOAA at Maryland. Trajecotries were initiated from each
source every 12 hours.
6.2.2 Results and Discussion
The results of our computations are shown in Figs. 32 through 37. The
average concentrations of $62 at the surface level and above the nighttime
stable layer height are shown in Figs. 32 and 33, respectively. As expected,
high concentrations of S02 were found in the area surrounding strong sources.
The eastern part of Ohio, the western part of Pennsylvania and West Virginia
were under a high concentration of S02, where several point sources with high
emission rates were located.
Distributions of sulfate were more uniform than those of S02, as can be
seen in Figs. 34 and 35. This is due to the fact that sulfate is formed with
time, and that a maximum concentration of sulfate is attained several hours
after S0~ emission from the source. The computed amounts of sulfate concen-
trations were several factors greater than the values observed normally.
The concentrations vary sensitively with changes in the removal terms,
as discussed in Chapter 7.
The distribution patterns at the surface level (Figs. 32 and 34) and
at the level above the nighttime stable layer (Figs. 33 and 35) were only
slightly different. This is probably due to the fact that during this 10-day
period the wind shear with height was not large enough to produce a signifi-
cant difference between the concentrations at the surface level and those at
the level above the nighttime stable layer.
The deposition amounts of SOp and sulfate, for the 10-day period, are
shown in Figs. 36 and 37. It can be seen that the patterns of contour lines
of SOp as well as of sulfate were similar to the patterns of concentrations
of the respective pollutants at the surface level. Since the deposition
amount is calculated as the product of surface level concentration times
60
-------
CM
-a
o
s-
O)
Q.
O>
_c
+J
s_
o
1/5
C
O
HI
u
c
o
-------
s-
o
O)
>
o
X)
re
/—>
en
E
cn
^.
*~-s
in
c
o
-p
re
01 ,
g!
u
CM <
O
CO
O)
en
c:
o
-l-> CL
in
• r- CU
ro
co
O)
S-
D)
62
-------
£_
-t->
03
in
c
o
s_
4->
c
1-
3
O)
63
-------
J-
o
OJ
CJ
01
+J
+J
_c
O5
c
Ol
-C
-(->
O)
>
o
.0
in
c
o
-13
nj
Ol
u •
c «*•
o r^
O CT)
1-H
Ol
-P »
re
s_
W> O>
55
LD
CO
01
3
O)
64
-------
05
c
en
c
O)
U
U)
to
u
QJ
S-
Q.
(0
in
o
CL
O)
"O
-a
o
-P
QJ
3
T3
CSJ
O
OO I
c\j
05
O5
to
•*->
c
3
O
O "
5 o
w 'si
li
T3
CLI
ro
-(j
o
OO
Ol
i.
13
CT)
65
-------
re
-p
u
CD
Q.
-a
s=.
03
c.
o
in
o
o.
OJ
T3
cr>
Q) r-l
3
-a -
o
/—.T-l
CM
=1
O
O)
•a
> o
3 OJ
o a.
to a>
o
+j o
O O5
O- C
•!-
T3 O5
C
i— OJ
m >
-i-> nj
o o
t— wi
ro
a>
S-
Z3
05
66
-------
deposition velocity, including dry and wet deposition, the similarity
mentioned above could be expected.
6.3 APPLICATION (2)
In this application we calculated the distributions of 24-hour average
concentrations and deposition amounts of SC^ and sulfate over the region
between 35°N and 45°N and between 75°W and 95°W. The concentrations of S02
and of sulfate were measured most extensively over the region on May 11,
1974. In order to verify our model we calculated the 24-hour average con-
centrations for the same day.
6.3.1 Input Data
From application (1), we learned that the concentration values of sul-
fate calculated over the region were several factors greater than the values
observed normally. In order to obtain more realistic and dependable results,
we used the following assumptions for the removal terms:
Dry deposition velocity of SOp
Daytime 2 cm/sec
Nighttime 2 cm/sec
Dry deposition velocity of sulfate
Daytime 0.4 cm/sec
Nighttime 0.4 cm/sec.
Because of the difficulty in obtaining hourly data of relative humidity over
the region for the past, the following values were used for the transforma-
tion rate of S02 to sulfate:
Daytime K = 0.1 hr'1
Nighttime K = 0.01 hr"1
Such a difference in the magnitudes of transformation rates between daytime
and nighttime is reported by Husar et al. (1977).
Another change made in the present calculations is in the value of the
dimensionless ratio (K/X)V for SOp and sulfate. Here K is the pollutant
content of rainwater and x is the concentration in air. In the previous
A
application, it was assumed that (K/X) for both SO,, and sulfate is 5 x 10 .
f- 5
In the present application, we assume that (K/X) is 5 x 10 . Eventually the
assumption of equal values of this ratio for both, S02 and sulfate, will have
to yield as more observational evidence becomes available. In Fig. 38 the
locations and the emission rates of major S0? sources used for the calcula-
tions are shown. As was the case in application (1), the SO^ sources whose
67
-------
S-
«j
O)
o
-p
CO
o
in
0)
u
3
O
in
in
in
O
to
in
cu
in
c
Ol
-a
c
(T3
C
O
u
o
ao
ro
a>
s-
^
CD
68
-------
n
emission rates are more than 10 ton/year are taken into consideration. The
emission of S0? from these sources contains about 90 percent of the total
emission over the area.
Clark's (1978) original emission rate estimate of S02 was revised after
application (1) was carried out. For the present purpose, gridded emission
maps of S0? for the United States and Canada east of the Rocky Mountains
supplied by EPA (Clark, 1979) were used.
The data of precipitation rate and precipitation probability for 81
stations located in the region were used to calculate the deposition amounts.
The following heights were used:
Daytime mixing height 1600 m
Nighttime stable layer . 250 m
6.3.2 Method of Calculation
Trajectories from each source area were calculated every 12 hours,
starting at 6 and 18 CST, and each trajectory was pursued for 48 hours.
Therefore, in order to calculate 24-hour average concentrations, we must take
into consideration the contributions of trajectories which started up to 2
days earlier. In Fig. 39, the scheme for calculating 24-hour concentrations
is shown. Trajectory segments drawn in solid lines were used to calculate
the 24-hour average concentrations. For this calculation, 1-degree grid
spacing was used to save computing time.
6.3.3 Results and Discussion
The distributions of 24-hour average concentrations of SC^ at the sur-
face and above the nighttime stable layer are shown in Figs. 40 and 41. It
can be seen that high concentrations of SC^ were calculated over the area
surrounding the sources.
Figures 42 and 43 contain the distributions of sulfate concentration at
the surface level and at the level above the nighttime stable layer. Both
figures show the very widespread area covering several states under high con-
centration of sulfate on this day. Clear air covered the western part of
this region (Iowa and Missouri).
The distributions of deposition amounts of S0? and sulfate due to both,
dry deposition and precipitation, are shown in Figs. 44 and 45. The dis-
tributions of deposition amounts of SOp and sulfate due to precipitation
only are shown in Figs. 46 and 47. On May 11, 1974, some of the eastern
states did not receive precipitation, thus no deposition due to precipitation
was recounted there.
69
-------
CM._
o" T
- T L L
OsJ.. 1
- T "
I
.. I - L
CM.. i -
".. i I
— CM —
II II II
T
T
1 T
T
T T T
. .
. 1
CM — CVJ
it it ii
I 1
T 1
T i
T 1
T 1
T I
T
1
ro
Q
CO
1
>-
§
— CM
n ii
^
-------
cr>
c
o
CD
>
Ol
HI
u
S-
3
I/)
Ol
1/5
C
O
c
Ol
o
u
CM
o
to
Ol
3
O5
71
-------
CNJ
oi
S
CO
-P
in
O)
-t->
0)
>
o
.0
to
c.
o
(O
S-
(J
c
o
CSJ
o
oo
o
•13
^
-Q
si
+J
en
O)
s_
72
-------
U>
C
o
O)
>
QJ
OJ
U
O)
TO
I/)
c
o
o>
u
c
o
U
OJ
+J
03
o
c
o
CM
«d-
Ol
s_
3
05
73
-------
S-
Ol
01
!n
(0
.u
in
Oi
-P
O)
-l-i
o>
>
o
_Q
CO
/—-s
co
c
o
(0
-(->
c:
o>
u
c:
o
u
Ol
-l->
(O
c.
o
(fl
CO
-------
c
o
c
o
rO
-t->
O-
U
-------
CSJ
03
-P
U
Ol
i.
O-
tO
IS)
o
a.
o>
TD
T3
O
4->
Ol
3
•o
O)
3
O
C
O
in
o
O.
a>
-a
o>
CVI
C. T
O
-!-> rH
to
•r- C
O O
LD
O)
S-
3
CT
76
-------
c
o
>>
c
o
c
o
•r—
+J
re
+j
•r~-
Q.
0)
Q.
-------
ra
-P
u
CD
CL
O)
3
-o
esj
1
O)
o
(O
c
o
-J3
•r—
o
Q.
3
O5
78
-------
The sulfur budget over the region was calculated as follows:
Sulfur emission 24,546 tons
Removal by wet deposition
as S02 4,293 tons
as sulfate 3,005 tons
7,298 tons
Removal by dry deposition
as S02 10,176 tons
as sulfate 2,771 tons
12,947 tons
Total deposition of sulfur 20,246 tons
Amount exported from the boundaries
of the region 4,300 tons
According to these calculations, about 30 percent of sulfur are removed by
precipitation, about 53 percent are removed by dry deposition and only 17
percent are exported from the boundary of the region. Of course, these
estimates are based on the removal parameters used.
In order to verify the performance of the model, the concentrations of
S02 and sulfate calculated for the surface level were compared with concen-
trations observed at stations located in the region. Observed data supplied
by EPA were used.
Numerous stations observed S02 concentrations on May 11, 1974. In order
to avoid data points containing local effects, we used only those data taken
at stations classified as "rural", "remote", and "suburban-residential" (EPA,
1974). Data taken at 73 stations defined as above, were compared with the
calculated concentrations of S02. On the other hand, the number of stations
which observed the sulfate concentrations is limited. Sulfate data taken at
41 stations, regardless of their classifications, were compared.
The observed concentrations of S02 are plotted in Fig. 48 against the
calculated ones. Similar plots for sulfate are shown in Fig. 49. Notice
that the scales are logarithmic in both figures. A few points in both dia-
grams are overlapping.
The correlation coefficients between the mean deviation of the observed
concentrations and mean deviation of the calculated concentrations were cal-
culated as follows:
for S02 concentrations, r = 0.635, r = 0.302
and for sulfate concentrations, r = 0.487, r = 0.393.
79
-------
fee
i- e
2 \
UJ
CM
LJ
a:
LU
CO
GO
o
10
1.0
O.I
0.1
*»«
H*'hk K'M' kik
LJIlJ
1.0
10
CALCULATED CONCENTRATION OF
Figure 48. Relationship between calculated and observed S02 concentration.
80
-------
If
£l
tt: 10
H e
LU CP
O
-
o:
LU
en
CD
o
10 r
CALCULATED CONCENTRATION OF
SULFATE (/xg/m3) —
Figure 49. Relationship between calculated and observed sulfate concentra-
tion.
81
-------
Here, r is the correlation coefficient, r is the critical value for cor-
relation coefficients with 99.9% confidence levels. Since r is larger than
rc in both SCL and sulfate, we can conclude that the distribution patterns of
calculated concentrations of SOp and sulfate are statistically related to
those of observed concentrations with the correlation coefficients r = 0.635
and 0.487, respectively.
6.4 SUMMARY
The model was applied to calculate the distributions of SCL and sulfate
concentrations over the region between 35°N and 45°N, and between 75°W and
95°W. In the first example, using climatological data of relative humidity
and precipitation rate, the 10-day average concentrations of SCL and sulfate
were calculated. Although the distribution patterns seem qualitatively
reasonable, it was noticed that the magnitudes of sulfate concentration were
several factors greater than those observed normally. The reason for this
discrepancy is that the magnitudes of removal terms used were inappropriate.
For this reason, the magnitudes of these removal terms were changed in the
second application, where the distribution of 24-hour average concentrations
were calculated. The concentrations of both SCL and sulfate thus obtained
were compared with those observed. The distribution patterns of calculated
concentrations of SC^ and sulfate were statistically related to those of
observed concentrations with the correlation coefficients r = 0.635 and
0.487, respectively.
With further refinements, we feel confident that correlation coef-
ficients with even larger values can be obtained. In order to save computing
time, we calculated the concentration values for 1 degree grid spacing in
the second application. If 0.5 degree grid spacing had been used, as in the
first application, we might have obtained better results.
82
-------
CHAPTER 7
REGIONAL RESIDENCE TIMES OF S02 AND SULFATE OVER
THE EASTERN UNITED STATES
7.1 INTRODUCTION
A parameter which adequately characterizes the fate of pollutants over
long time- and space-scales is the "residence time" or "turnover time" in the
atmosphere. The residence time of a pollutant can be defined as the time
required to decrease the burden of that pollutant in the atmosphere over a
certain region to a value which is 1/e of the original value, assuming that
no further input of the pollutant into the atmosphere from below or across
the lateral boundaries of the region occurs. The turnover time is the time
equal to the total mass of a pollutant in the atmosphere over a certain
region divided by the removal flux of the pollutant.
We calculated the residence time of SOp over the eastern U.S. (Henmi et
al., 1977 and 1978). However, there were shortcomings in those earlier
studies: (i) the transformation rate of SOp to sulfate was assumed to be
constant and (ii) the residence time of sulfate was not considered.
In the present study, the observational data on transformation rates
were scrutinized and an empirical formula for the transformation rate as a
function of relative humidity was deduced. The residence time of sulfate was
defined mathematically.
Calculations were made of the regional residence times of S02 and sul-
fate over the region of the United States east of 105°W longitude where major
industrial activities are located. The year was divided into the cold season
(January - April, and November and December) and the warm season (May -
October), and the regional residence times for each of these two seasons were
calculated.
7.2 DATA USED
7.2.1 Mean Mixing-Layer Height, H, and Precipitation Data
Mean mixing-layer height and precipitation data were analyzed in a
previous study (Henmi et al. , 1977 and 1978). The same data were used here.
83
-------
7.2.2 Relative Humidity Data
Since the transformation rate of S0? to sulfate is derived as a function
of relative humidity, data on relative humidity are necessary as input. The
monthly average relative humidities of the year of 1974 were obtained from
the Climatological Data, National Summary (U.S. Department of Commerce,
1974), and the seasonal average values for 61 stations located in the study
area were used for the computation of transformation data.
Figures 50 and 51 show the seasonal average relative humidities over
the eastern U.S. for the cold season and the warm season, respectively.
7.3 APPROACH
7.3.1 Residence Times of S02 and Sulfate
In order to calculate the residence times of S02 and sulfate, we assume
that they are distributed uniformly vertically through the mixing layer, that
no significant net transports of the pollutants occur across the boundaries
of the region, that the pollutants are removed from the layer by dry deposi-
tion and precipitation scavenging, and that SOp is transformed into sulfate
by a first-order reaction. Under these assumptions, we can write
dC-,
- + k
Pl
dC2 3
-af = -(xd2 + VC2 + I Kci
C, and C0 are the concentrations of S00 and of sulfate, A, and A, are the
12 2 ' d-. dp
depletion rates due to dry deposition, k and k are the depletion rates
pl P2
due to precipitation. Subscripts 1 and 2 stand for S00 and sulfate. The
3
transformation rate from S02 to sulfate is represented by K. The ratio j is
the ratio of molecular weights of sulfate (SOT) to sulfur dioxide (SOp). The
solutions to these equations are given by
C1 = C1Q exp[-(Ad + kp + K)t] (23)
84
-------
Figure 50. Seasonally averaged relative humidity over the eastern United
States for the cold season.
85
-------
Figure 51. Seasonally averaged relative humidity over the eastern United
States for the warm season.
86
-------
~\f\
kp2)t]
exp[-(\. + k + K)t]j (24)
dl pl '
where C-.Q and C^g are initial concentrations of SCL and sulfate, respective-
ly. Since we assumed that sulfate is produced only by the transformation of
S02, CpQ = 0. Therefore,
K)
exp[-(\ri + k + K)t]| (25)
°1 pl '
The concentrations C-, and C2 computed from Eqs. (23) and (25) are shown
in Fig. 52. According to this figure, the e-folding residence time of S0?
T<.Q , can be defined as
S09 K . + k + K
2 d P
(26)
The time T at which the concentration of sulfate, Cp, reaches a maximum,
2 ,
C , can be calculated from
dC
-i«
Therefore,
= 0
-
K) -
+ k + K)-exp[-(\. + k + K)t]J
1 "l 1 "l '
ln
(27)
The maximum concentration of sulfate, Cp-, is given by
87
-------
0
Figure 52.
Tm Ts
S02
TIME
Tc
SUL
Schematic plot of the concentrations of
function of time.
and sulfate as a
-------
3 K
C2m = 1 (Ad + k + K) - (Ad + k ) C]
k_ + K
•exp
In
(A , + A + K) - (A. + A J A, + A + K
1 Pi 9 Po 1 Pi
exp
-i
(28)
We define the residence time of sulfate 1 -. as the time elapsed until the
concentration of sulfate, Cp, becomes e 'C*-. Therefore,
2 _
exp -(A, + k + K) T
sulj '
-1
exp
~ P
~\ + V + K
-1
•
exp
(A. + kn + K\
dT P! \
n -1 1
\*\ )
- 1
(2<
0
Equation (29) is solved numerically with iteration methods.
7.3.2 Chemical Transformation
There are several observational studies available on the transformation
rate of SOp to sulfate. Unfortunately, the studies conducted in the field
reveal a variability of the transformation rate over three orders of magni-
tude, depending on the concentrations of the plume, relative humidity, etc.
In this study we are interested in the fates of SCL and of sulfate after both
pollutants have been uniformly mixed in the mixing layer. For this reason,
the observational data of the concentrations of S0? and sulfate taken at
least one hour after emissions were scrutinized, and the transformation rates
were plotted in Fig. 53 as a function of relative humidity. From this
figure, the transformation rate K is derived as
K = 3.304 x 10"4 • exp(0.063 • RH) hr"1
where RH is the relative humidity in percent.
(30)
In previous studies (Henmi et al., 1977 and 1978), a constant value of
1 x 10 sec was used for the transformation rate. In the present study,
the relationship of Eq. 30 was used.
89
-------
1.0
Smith and Jeffrey (1974)
University of Utah (1974)
Lusis and Phillips (1977)
Gartrell (1963)
.1
LU
OC
z
o
cr
O
P. 01
.001
40
60
Figure 53.
80
Transformation rate of S02 to sulfate as a function of relative
humidity.
90
-------
7.3.3 Dry Deposition and Precipitation Scavenging
Deposition velocities for SOp have been observed by a number of in-
vestigators. In the previous study, we reviewed several reports and assumed
that, for our purpose, the dry deposition velocity of SCL is 1 cm/sec re-
gardless of the region and the season. Little is known about the deposition
velocity of sulfate. The only known measurements of sulfate deposition
velocity in the atmosphere yielded approximately 0.1 cm/sec (Engelman and
Sehmel, 1976).
Accordingly, the values
V = 1 cm/sec
gl
Vg2 =0.1 cm/sec
were used in the present study. Here, V , is the deposition velocity of
S00 and Vrt is that of sulfate. 9
2 g2
The depletion rate A . due to dry deposition can be written as
H
where H is the mean mixing layer height.
The scavenging rate due to precipitation, \ , can be expressed as
A = (*) . ! = !§
p x v H H
where K is the concentration of S02 or SOT in rainwater, x the concentrations
of those pollutants in air, and P the mean precipitation rate during a
precipitation event. The subscript v means that the rat-'o (K/X) is formed
on a volume basis. The depletion rate due to precipitation, k , is given by
where f is the probability of precipitation events. It can be expressed by
t
f =
i is the length of the average wet period and i. is the length of the aver-
age dry period. Data of t and T, were shown by Henmi et al. (1977 and
1978). p d
91
-------
Values of (-) for SCL have been scrutinized in our previous report, and
4
the value of 5 x 10 has been chosen as typical. We have little knowledge
about the value of (-) for sulfate. Therefore, in this study we assume,
tentatively, the value of 5 x 10 for that pollutant to be identical with
that for S02.
7.4 RESULTS
The regional residence times of S02 in hours, T<-Q , for the cold and
warm seasons are shown in Figs. 54 and 55. The regional residence times
of sulfate in hours, T^y. , for the cold and warm season are given in Figs.
56 and 57.
The following can be seen from these diagrams:
1. Both T^Q and T<.y, are longer in the warm season than in the cold
season over the whole region studied. This is due to shorter dry periods
and shallower mixing layers in the cold season than in the warm season.
2. Short residence times characterize the region surrounding the
Great Lakes and the coastal regions.
3. Values for T<-Q and Tq,, are found to be relatively long in the
western parts of the area under investigation, where the mixing layer height
is large, the precipitation frequency is small and the relative humidity is
small.
4. The regional residence times of S0«, T™ , are in the range between
15 and 30 hours for the cold season and in the range between 15 and 40 hours
for the warm season. T<.Q in the present study is shorter than the values
derived in the previous study. The reason for this difference is that in
the previous study a constant value of 1 x 10 sec for the transforma-
tion rate was used.
5. The regional residence time of sulfate, T, , defined in the pre-
sent study is about one order of magnitude greater than T™ . Values for
Tgy, are in the range between 150 and 450 hours for the cold season and in
the range between 200 and 500 hours for the warm season.
92
-------
Figure 54. The regional residence time of SCL in hours,
season.
, for the cold
93
-------
40
Figure 55. The regional residence time of S00 in hours, Tcn , for the warm
season. L bU2
94
-------
Figure 56. The regional residence time of sulfate in hours,
cold season.
, for the
95
-------
Figure 57. The regional residence time of sulfate in hours,
warm season.
, for the
96
-------
7.5 SENSITIVITY ANALYSIS
As can be seen from Eqns. (26) and (29), the residence times of SOp and
sulfate are dependent on several physical parameters whose values are deter-
mined empirically or experimentally. These parameters include the dry de-
position velocity of S0?, Vg-,, the dry deposition velocity of sulfate, Vg^,
the transformation rate of SQ~ to sulfate, K, and the dimensionless ratios
of SO, and sulfate contents in rainwater and air, (-)cn and (-)Cn • The
^ X ->U2 X ->u 4
actual values of these parameters may be different from those used in our
calculations of residence time. Thus, it is necessary to perform sensitivity
analyses by varying systematically the values of these parameters.
For each parameter, calculations were made with values varying about the
reference value, whereas the other parameters remain at their reference val-
ues. Changes in residence times of S0? and sulfate are given in Table 9.
The changes of residence times given in this table are based on climatologi-
cal data for the cold season obtained at Charleston, West Virgina. Calcula-
tions for several stations show a similar range of changes in residence
times.
This table also indicates that the residence times change substantially
with the changes of these parameters. Therefore, we must regard our results
as approximate until better knowledge about these parameters is obtained.
7.6 SUMMARY
In this study, we have extended a previous investigation of the regional
residence times of SCL, using an empirical relationship between the trans-
formation rate K and the relative humidity. The residence time of sulfate,
T<-y, , was defined mathematically and T^,,, was also calculated over the region
of the United States east of 105°W longitude where major industrial activi-
ties are located.
The residence times of S02 and sulfate were calculated over the region
of the United States east of 105°W longitude where major industrial activi-
ties are located. The regional residence times of S02 are in the range be-
tween 15 and 30 hours for the cold season and in the range between 15 and 40
hours for the warm season.
The regional residence time of sulfate is about one order of magnitude
greater than that of SO^, ranging between 150 and 450 hours for cold sea-
son and between 200 and 500 hours for the warm season.
97
-------
CO
oc
LU
I—
LU
s:
—
1—
1 — 1
l>»
1 — 1
\—
CO
z
LU
CO
,
cn
LU
_J
OQ
£
O>
XT
_4_3
o
in
<+- 1-
O CD
-P
tn o>
OJ E
rs re
t— t-
re re
> Q-
Ol
r-
O)
**• E
O -r-
CO I—
C 01
•i- U
C
~>
CO CO MIX MIX
> > *•—' v— •*
S3 &?
CM
CM ro
+ i
a>
CM E
0 T-
CO r-
c ai
•r- U
t~~
ai a>
cn -o
C T-
re cn
-C CD
C_5 Q£
^^ ^^
f~^i f~^
co r^-
+ i
£-
x:
a> TD
=> 0)
i — in
> "^
\ s-
r~ x:
LO ^^
CM r^»
CD LO
CD CM
O)
U
C
at
Qi
j-
x:
^»
r-
LO
CM
CD
CD
E 0
O -P
• r—
-P CM
ro CD
£-
O)
-P
0)
E
re
re
CL.
E CO
S-
0 t-
^» o
tn
C. O)
re -p ^d~
s- re o
1— S- CO ^£.
u
m
o in
0) \ t.
in E -c
"^v U ^^
E r^^
U i— 1 LO
CM
rH 0 O
II II o'
rH CM ||
cn cn
> > ^.
"&Q ^Q
csj i_n
. ,
00 U5
0 00
CO
|
-J-
^^
to ^^
• CD
CD
VO OT
p*^
+ i
^~ ID
0 O
i-H rH
X X
LO LO
LO
CD
rH -3-
X 0
LO CO
T3 C
C «r-
re re i^*
in S- o
tn CM co
0) O C. , — .
i— CO T- MIX
C "^t
O M- in
•i- 0 -P S-
in c -r— CM
c o a> re o
cu T- -P cn
'E -p c -a /->
•r- re o c: MIX
a s- o re ^
"^" LO
U CD O S-
a> T-H i— i x:
in x X \
^v LO LO r-~
E LO
O II II CM
i-H 0
CM *3" •
II 000
CO CO
T-\ ^ <~^ \\
cn M| x M| x
>
cn r^. to CM
• • * •
CM LO CO i—
1 II r—
1
^^ ^^ ^? ^^
CD CD CD CD
• • * •
CD O CD O
U
a>
tn
'*j
•P | *
• r— *f*.
in o CM
o o cn
a. i— >
Ol O)
o >
u
O) ^d" LD
tn O CD
-\^ f«— fl T~H
E X X
U LO LO
rH
II II
0
CM ^J
II 0 0
CO CO
CM /—s *~^
Cn Ml x Mix
> ^ ^
LO r*.
. .
00 LO
1 rH
^? fe^
0 0
* ,
co r*^«
CM CO
1 i
U U
a) 01
tn in
E E
U U
CM CO
(J
•1 — '1 —
CO O rH
o o cn
0. r~ >
0) CU
a >
i-
x:
**\
r^
LO
CM
o
,
o
II
-------
CHAPTER 8
NUMERICAL STUDY ON SULFATE CONTENT IN RAINWATER
8.1 INTRODUCTION
The extent, severity, and causes of acid precipitation are becoming
topics of major concern in the United States. In Europe it was reported many
years ago that the acidification of lakes and rivers in Scandinavia was
caused mainly by long-range transport of pollutants from Biitain and from the
European continent. Acid rain in the northeastern United States is associ-
ated with air parcels that have travelled through major industrial areas in
the midwestern and east-central states.
Acid rains, associated with sulfate pollution, cause structural, agri-
cultural and ecological damage. Because of this potential damage, it is
important to predict accurately air and water quality under future use of
fuels with high sulfur content.
Reliable continuous precipitation chemistry data for the United States
are sparse. However, several attempts have been made to determine the chemi-
cal composition of rain. Junge and Werby (1958) measured the concentrations
of chloride, sodium, potassium, calcium, and sulfate in rainwater over the
United States. Average sulfate concentrations over the United States were in
the range of 1 mg/liter to 5 mg/liter. Recent studies in a region surround-
ing the St. Louis area (Hales and Dana, 1979) show the extremely wide vari-
ability in chemical content from one event to another. However, the average
sulfate content of rainwater is in the range between 4 mg/liter and 17
mg/liter. Sulfate concentrations in rainwater observed in Ithaca, N.Y., and
Aurora, N.Y., are within the range of 1 to 8 mg/liter (Dittenhoefer and
Dethier, 1976). Thus, it may be concluded that, although the sulfate content
of rainwater varies with storm type, trajectory of air mass, geopraphical lo-
cation, etc., it is within the range of 1 to 20 mg/liter.
Recently, Scott (1978) has proposed that sufficient sulfate aerosol
exists to account for observed levels of sulfate in rainwater. According to
Scott, sulfate formation from sulfur dioxide in cloud water and rainwater is
not fast enough to produce appreciable quantities of sulfate during the
lifetime of the cloud particles contributing to precipitation.
The concentration of sulfate in precipitation is a consequence of sever-
al microphysical processes occurring within and beneath the clouds. Micro-
physical processes such as Brownian motion, attachments due to thermophoresis
99
-------
and diffusiophoresis, inertia! impaction and nucleation serve to remove the
sulfate aerosol from the air.
The purpose of this chapter is to examine the validity of Scott's idea.
Taking into consideration the microphysical processes mentioned above, model
calculations of rainout processes of sulfate aerosols were carried out.
8.2 THE CLOUD MODEL
In the case of a precipitating cumulus cloud, the discrete particle
_0
sizes of importance range from aerosol particles of the order of 10 micro-
meters in diameter to precipitation drops as large as 5x10 micrometers in
diameter. Therefore, in order to study the interaction of pollutant par-
ticles and cloudwater droplets, a model must be chosen in which the micro-
physical processes are incorporated in as much detail as possible. One such
model is the one-dimensional cumulus model by Cotton (1972a, b). This model
involves upward integration with height following the rise of the convective
bubble or plume.
A numerical model which is composed of cloud dynamics and microphysics
of cloud particles and aerosol has been developed (Henmi et al., 1977, 1978).
In this model, the aerosols in the planetary boundary layer are transported
up into the cumulus cloud. The cloud is envisaged as an assembly of cloud-
water droplets and rainwater droplets intermingled with aerosol pollutants,
some of which are free-floating in the cloud air and some are collected by
the droplets. This model is used for the present calculations to estimate
the amount of aerosol pollutants captured by cloudwater and rainwater drop-
lets.
8.3 SULFATE AEROSOLS IN THE ATMOSPHERE
Atmospheric aerosol is composed of "fine" particles having a diameter
smaller than 2 micrometers and coarse" particles having a diameter greater
than 2 micrometers. Those fine and coarse modes are usually quite different
in chemical properties. The physical separation of the fine and coarse modes
is justified because the condensation process produces fine particles whereas
mechanical processes produce mostly coarse particles. Practically all of the
sulfate found in atmospheric aerosol is in the fine mode (Whitby, 1978).
According to Whitby (1978), on the average, the sulfate content of con-
tinental aerosol is 15 to 30 percent of the total aerosol mass, whereas
marine aerosols contain sulfate in the range of 30 to 60 percent.
In Fig. 58 the size distributions of typical urban aerosol and back-
ground aerosol are shown. The mass concentrations of urban aerosol and
3 3
background aerosol given in this figure are 99.6 |jg/m and 48.7 ug/m , re-
3 3
spectively. The accumulation modes contain 61.3 ug/m and 18.7 ugm , of
urban aerosol and background aerosol. For the present study, we assume that
50 percent of the mass of aerosol in the accumulation mode is composed of
sulfate.
100
-------
I05E
ro
i
o
a.
Q
£
<
URBAN
BACKGROUND
0.001
100
Figure 58. Size distributions of typical urban and background aerosols
(Whitby, 1978).
101
-------
8.4 RESULTS
Following the procedures described by Henmi et al. (1977, 1978), numeri-
cal calculations of cumulus clouds are performed. In order to generate
different clouds, different input data of sounding and entrainment constants
are used. In the present calculation, the concentration of cloudwater drop-
_0
lets is assumed to be 300 cm , which is typical for continental clouds. The
theoretical calculations by Howell (1949), Mordy (1960), and Neiburger and
Chi en (1960) show that the concentration of cloudwater droplets is determined
by the number of effective condensation nuclei in the cloud base region.
Furthermore, the particles effective as condensation nuclei are, in general,
hygroscopic. Sulfate particles are known to be very hygroscopic (Fletcher,
1966). Therefore, two cases are studied. In the first case, the simple
assumption is made that the 300 largest aerosol particles per cm in the
accumulation mode are consumed as condensation nuclei within the first layer
of integration above the cloud base. In the second case, condensation nuclei
are assumed to be supplied by other sources. Sulfate particles are captured
by cloudwater droplets and rainwater droplets through microphysical mecha-
nisms such as Brownian diffusion, inertial impaction and attachment due to
thermophoresis and diffusiophoresis.
Using the aerosol spectra shown in Fig. 58, the sulfate content of
rainwater is calculated. In Fig. 59, the sulfate content of rainwater is
shown as a function of the average liquid water content (LWC) of the cloud.
This average liquid water content is calculated as a value averaged over the
cloud height. In general, LWC is proportional to the height of the cloud and
to the precipitation amount. In Fig. 59 the broken lines represent the
3
result of the first case where the 300 largest aerosol particles per cm in
the accumulation mode are consumed as condensation nuclei. The solid lines
represent the results of the second case where sulfate particles are not used
as condensation nuclei.
From this figure, the following can be seen:
i) The sulfate content in rainwater is inversely proportional to LWC.
ii) The condensation process is a dominant mechanism for determining
the sulfate content in rainwater. However, about 60 percent of the
sulfate mass in rainwater is captured by other microphysical mech-
anisms.
iii) The concentration of aerosol is a determining factor of the sulfate
content in rainwater. Rainwater which has collected urban aerosol
contains more sulfate than rainwater which has collected background
aerosol.
iv) The sulfate content of rainwater from these simulated clouds is
within the range of 1 and 8 mg/1. These values are well within the
range suggested by observations.
102
-------
-------
Figure 60 shows the relationship between (K/X) and LWC, where (K/X) is
the dimension!ess ratio of pollutant concentrations in rainwater and air.
(K/X) is inversely proportional to LWC. Hales and Dana (1979) derived a
similar relationship as
(*) = 1.9xl06A-°'88
A.
where, A is the rain amount (cm).
tf
Similarly, an inverse relationship between (-) and precipitation amount,
A
precipitation rate or LWC has been found by Engelmann (1971). The observa-
tions by Scott and Lanlainen (1979) show that (K/X) for sulfate was 7.9x10
and 2.03x10 in two different rain events. It seems that (K/X) values from
our simulated clouds are close to the observed values.
8.5 CONCLUSION
Although the present study has shown that microphysical processes such
as condensation, Brownian diffusion, attachment due to thermophoresis and
diffusiophoresis and inertial impaction can be efficient enough to attain the
sulfate content of rainwater, this conclusion should be regarded as tenta-
tive. The sulfate content of aerosol in the present study was assumed to be
3 3
about 30 ug/m and 9 ug/m for urban and background aerosols. A recent study
by Hales and Dana (1979) is contradictory to the study by Scott (1978) and to
our present study. According to Hales and Dana, a rapid oxidation of S02 to
sulfate occurs in cloud systems in warm, polluted environments and leads to a
possible explanation for the observed sulfate content of rainwater.
The combined effects of the microphysical processes and of S02 oxidation
to sulfate in water must be studied in the future.
104
-------
URBAN
AEROSOL
BACKGROUND
AEROSOL
3.0 r
t
2.0
in
O
1.0
WITH CONDENSATION
WITHOUT CONDENSATION
WITH CONDENSATION
WITHOUT CONDENSATION
0.5 1.0 1.5 2.0
L.W.C. (g/m3)
Figure 60. (-) as a function of L.W.C. of the cloud.
105
-------
REFERENCES
Charlson, R.J. 1969: Atmospheric visibility related to aerosol mass
concentration. Environ. Sci. Tech., 3, 915-918.
, R.J. Vanderpol, D.S. Covert, A.P. Waggoner and N.C. Ahlquist, 1976:
Sulfuric acid and ammonium sulfate aerosols: Optical detection in the
St. Louis region. Science, 184, p. 156.
Clark, T.L., 1978: Private communication from Mr. Holzworth, EPA.
, T.L., 1979: Gridded annual pollutant emissions east of the Rockies.
Private communication.
Cotton, W.R., 1972: Numerical simulation of precipitation development in
supercooled cumuli - Part I. Mon. Wea. Rev., 100, 757-763.
, 1972: Numerical simulation of precipitation development in super-
cooled cumuli - Part II. Mon. Wea. Rev., 100, 764-784.
Dittenhoeffer, A.C., and B.E. Dethier, 1976: The precipitation chemistry of
western New York State: A meteorological interpretation. Research
Project Technical Completion Report OWRT, Project No. A-044-NY. Agree-
ment No. 14010001-3532, Cornell University Water Resources and Marine
Science Center, Ithaca, New York.
Engelmann, R.G. and Schmel, G.A., 1976: Atmosphere-sulfate exchange of par-
ticulate and gaseous pollutants. ERDA Document CONF-740921, U.S. Energy
Research and Development Administration, Oak Ridge, TN.
, 1971: Scavenging prediction using ratios of concentrations in air
and precipitation. J. Appl. Met., 10, 493-497.
Gartrell, F.E., F.W. Thomas, and S.B. Carpenter, 1963: Atmospheric oxidation
of S09 in coal-burning power plant plumes. American Industrial Hyg.
Assoc/J., 24, 113-120.
Fletcher, N.H., 1966: The Physics of Rainclouds. Cambridge at University
Press, pp. 390.
Hales, J.M. and M.T. Dana, 1979: Precipitation scavenging of urban pollut-
ants by convective storm systems. Journal of App1ied Meteorology, 18,
294-316.
106
-------
Heffter, J.L., and A.D. Taylor, 1975: A regional-continental scale trans-
port, diffusion and deposition model. NOAA Technical Memorandum, ERL
ARL-50. Air Resources Laboratories, Silver Springs, Maryland.
Henmi, T. , E.R. Renter, and R. Edson, 1977 Residence time of atmospheric
pollutants and long-range transport. Environmental Research Paper, No.
]2_ Colorado State University, Fort Collins, Colorado.
, and , 1978: Residence time of atmospheric pollutants and
long-range transport. EPA-600/4-78-003, Environmental Monitoring Series.
Environmental Sciences Research Laboratory, Office of Research and De-
velopment, U.S. Environmental Protection Agency, Research Triangle Park,
North Carolina 27711.
Hidy, G.M. et al., 1974: Characterization of aerosols in California, Final
Report IV, California Air Resources Board.
Holzworth, G.C., 1972: Mixing heights, wind speeds, and potential for urban
air pollution throughout the contiguous United States. Ap-101, Office
of Air Programs, Environmental Protection Agency.
Howe11, W.E., 1949: The growth of cloud drops in uniformly cooled air. J.
Met., 6, 134-149.
Junge, C.E., and R.T. Werby, 1958: The concentration of chloride, sodium,
potassium, calcium, and sulfate in rainwater over the United States, J.
Meteorol., 15, 417-425.
Lusis, M.A., and C.R. Phillips, 1977: The oxidation of S02 to sulfates
in dispersing plumes. Atmospheric Environment, jj, 239-241.
Mordy, W.A., 1959: Computations of the growth by condensation of a popula-
tion of cloud droplets. Tell us, 11, 16-44.
Neiburger, M. , and C.W. Chein, 1960: Computations of the growth of cloud
drops by condensation using an electronic digital computer. Physics of
Precipitation, Geophysical Monograph Series, 5, 191-208, AGU, Washing-
ton, D.C.
Scott, B.C., 1978: Parameterization of sulfate removal by precipitation.
Journal of Applied Meteorology, 17, 1375-1396.
and N.S. Lanlainen, 1979: On the concentration of sulfate in precipi-
tation. Journal £f Applied Meteorology, 18, 138-147.
Smith, F.B., and G.H. Jeffrey, 1975: Airborne transport of sulphur dioxide
from the U.K. Atmospheric Environment, 9, 643-659.
University of Utah, 1974: Cited in "Power plant stack plumes in complex
terrain, an appraisal of current research" by R.C. Koch et al.
EPA-600/7-77-020. Environmental Sciences Research Laboratory Office
of Research and Development. U.S. Environmental Protection Agency.
107
-------
United States Department of Commerce, 1974: Climatological Data National
Survey.
United States Environmental Protection Agency, 1976: Directory of air
quality monitoring sites active in 1974. EPA 450/2-76-008. U.S. EPA.
Monitoring and Data Analysis Division, Research Triangle Park, NC 27711.
Whitby, K.T., 1978: The physical characteristics of sulfur aerosols, At-
mospheric Environment, 12, 135-159.
108
-------
APPENDIX - A
FORTRAN LISTING OF MODEL A
Trajectory Calculation Program (Model A)
PROGRAM THAJET
COMMON /STAT/ SLAT(100.4.8), SLON(100.4.8), XAw(100.4.8).
Z YAWllOO.4,8) t XSIGHA(100.4,8), YSIGMAUOO, 4,8)
••FORWARD AND BACKWARD TRAJECTORY MODEL
••USING OBSERVED AND ANALYZED INPUT WINDS
LOGICAL SICNT
••FOR OBSERVED INPUT WINDS ONLY»*****»«*»*«*********»*»**********»***
DIMENSION SLATT(IOO), SLONT(IOO), XAWTUOO). YAWT(lOO)
DIMENSION NSTAI4. 8)
DIMENSION XXSIGUOO). YYSIG(IOO)
INTEGER SlU(100t 4)
INTEGER TSLAT.TSLON
••SEE DATA MxNSTA
••*••*«*•«*•*••••*•«••••••**«•«««*«*•«••«••••••****«•«•••••••••••••*•
»•••••*•»•»»•••••«»«»••••<
LOGICAL SJTAG15).SICNTT
INTEGER xiPod), YIP061 .XTSIPOOI »YTSIPOO) ,wip<30)
DIMENSION 010(5)iOLAT(5)iOLON(5)
DIMENSION wioHOT(5)
INTEGER DXSIGP130),DYSIGP(30)
_ ._. ._.. .IBYR.NDY.TDIR,
COMMON/INPyTl/IBDY.lBMO.IBYR.f
A POLLAT.JHOUR.IMO.
A NDYDUR,NOYOTA,tfTYPE.(
...... .......... _,DTAST»DTA8B,OTA8LfDTABRt
A LBAAT,LTAAT,ALATT,ALATB,ALONL
••MAXIMUM NUMBER OF ORIGINS
DATA MXNO/S/
••MAXIMUM NUMBER OF STATIONS
DATA MXNSTA/100/
DATA OBStANA/3HOBS,3HANA/
DATA BACK/4HBACK/
DATA SJTAG/1H , 1H. , IH*, 1H«, 1H-X
DATA ITS.ITPDY,IMM/20, 4, 8/
0
Z NTS, NMM)
••SET ITAB FORWARD/BACKWARD TO i FOR FORWARDI 2 FOR BACKWARD
ITABF8«1
IF (TDIR .EQ. BACK) ITABFB«2
••LOOP THROUGH NUMBER OF DAYS
DO 400 M=1.NDY
CALL DTABKO ( M, MXNSTA, NTPDY» NMM» NSTA, ITABF8, NDYDTA,
z iTPDYt IMM, DTABT, DTABB, DTABL, DTABR, LBAAT, LTAAT,
Z OLON(l), SID)
••LOOP THROUGH NUMBER OF ORIGINS
DO 399 IO»1,NO
••LOOP THROUGH NUMBER OF TRAJECTORIES PER DAY
DO 300 K-l.NTPDY
IF (ITABFB ,EQ. 2) GO TO 32
MM«1
KK«1
GO TO 33
32 MM«NMM
KK = 4
••LOOP THROUGH NUMBER OF TRAJECTORY SEGMENTS
33 00 299
DO 299 J = 1, NTS
IF(J.EQ.l) 80 TO 40
IF (TSLAT .LE. 900) GO TO 42
••TRAJECTORY HAS BEEN TERMINATED
••DUE TO INSUFFICIENT AMOUNT OF INPUT DATA-
••TERMINATION CODE INFORMATION
••REPLACES ALL REMAINING TRAJECTORY SEGMENT LAT,LON VALUES
109
-------
GO TO 200
40 TSLAT1»OLAT(IO)
TSLON1«OLON(IO)
SPRED»0.
SPREDl'O.
42 JTAG*0
CALL KKMM i j. «• M. KK, MM, NMM, NDYDTA, NTSPDY»ITABFB,IRTNS>
••TRAJECTORY TERMINATED DUE TO NO DATA READ INTO DATA-BLOCK(DTA8K)
••OR CLOSEST-TIME(CLSTM) CALCULATION OUTSIDE DATA-BLOCK
TSLAT
TSL
ON»9996
UBAR»0.
SPREAD=0.
ALPHA=99.
B£TA*99.
T=99.
R2=99.
IJK«0
ICNT»0
SICNT=SJTAG<1)
GO TO 200
90 CALL OISO { J, K, M, KK, MM, MXNSTA, NSTA, ITABFB.
A ITS. ITPOY, IMM, RADIUS. TSDUR* MBIP. MEIP, SLATT, SLONT, TSLATi
A, TSLON1, XAWT, YAWT. XXSIG, YYSIG, KKT, MMT, SJTAG, JTAG. ICNTMN,
A ICNTMX, XIP. YIP, XTSIP, YTSIP, DXSIGP. OYSIGP, WIP, TSLATE,
A TSLON2, ISNAP, wtOROT(lO), UBART. SIGMAH, ICNTT, SICNTT, SPRED1,
A SPREDZ, ALPHAA, BETAA, TA, RSA, IJKA. IRTNSI
IF (IRTNS .EQ. 1TO) GO TO 170
CALL ALTDTA ( SJTAG, JTAG, TSLAT2, TSLON2, UBART, SIGMAH,
2 ICNTT, SICNTT, ALPHAA, BETAA, TA, R2A, IJKA, IRTNS)
IF (IRTNS .EQ. 170) GO TO 170
UO CALL CLSTM ( J, KK, MM, JTAG, NTPDY, ITABFB, M, NDYDTA, NMM,
Z KKT, MMT, IRTNS)
IF (IRTNS .EQ. 90) GO TO 90
IF(IRTNS .EQ. 66) GO TO 86
170 TSLAT»TSLAT2»10,«,5
TSLON=INT(TSLON2« 10.* 1800,5) -1800
UBAR=UBART«.5
ICNT»ICNTT
SICNT«SICNTT
TSLAT1»TSLAT2
TSLON1»TSLON2
SPRED1»SP»ED2
SPRED=SPRED« SIGMAH
SPREAD»SPRED/60.
AlPHAaALPHAA
BETA*BETAA
T»TA
R2>R2A
IJK=IJKA
200 KUNIT*20*IO
WRITE(KUNlT) TSLAT, TSLON, USAR, SPREAD, ALPHA, BETA, T, R2, UK,
Z ICNTi SICNT
IF(TSLON
299 CONTINUE
IFJTSLON .EQ. 9996) GO TO 300
300 CONTINUE
399 CONTINUE
400 CONTINUE
••END OF TIME AND ORIGIN LOOPS
WRITE TO A PERMANENT FILE THOSE VARIABLES WHICH
WILL BE NEEDED BY THE SECOND PROGRAM.
DO 475 I » 1, 5
WRITElU, 455) OlDtl), OLATII), OLONJI)
475 CONTINUE
WRIT£(14. 480) NO. NDY. NTS, NTPDY, LTAAT
WRITEI14, 490) TSDUR, ALATT, ALONL
455 FORMAT
480 FORMATI5I10)
490 FORMAT13F15.5)
SUBROUTINE ALTOTA(SJTAG, JTAG, TSLAT2, TSLON2, UBART, SIGMAH,
Z ICNTT, SICNTT, ALPHAA, BETAA, TA, R2A, IJKA, IRTNS)
••ALTERNATE DATA SELECTION
LOGICAL SJTAG(S),SICNTT
DATA OSS,ANA/3HOBS,3HAN
IRTNS»0
JTAG«JTA6*1
IF (JTAG .LE. 4) RETURN
NA/
110
-------
F (JTAG .LE. 4) RETURN
w,-
TSLAT2S99.9
TSLON2=999.9
UBART«0.
SIGMAH*0,
ALPHAA*99.
BETAA-99.
TA«99,
R2A*99.
IJKA»0
ICNTT«--..
SICNTT'SJ
IRTNS-170
NTT«ICNTT
SICNTTǤJTAO(JTAGT)
REJ
ENC
RETURN
"NO
SUBROUTINE AHlNOlSHTt SATH, 01, USAAT, LTAAT, WHT« XW, YMi NLVL,
Z XAWt YAWi XSIGMA, YSIGMA)
••AVERAGE WINDS
INTEGER SlD.SHT,SATHtBTH,TTH,TH
INTEGER SHT,SATHiBTH,TTH,TH
INTEGER WHT(50)
DIMENSION XW(50) ,YW(50)
DIMENSION MHT(IOO)
B«SHT*LBAAT
F (LBAAT ,NE. 0 .OR. DI .GT. 60. .OR. SHT ,GT. SATH) |_B = 3ATH»LBAAT
UT*SATH*LTAAT
STHXW=0.
STHYW«0.
L
I
.
••CALCULATE M HEIGHTS FOR WIND LEVELS
IF(NLVL.EO.l) GO TO 21
MLVL=NLVL«2-1
DO 12 M=1,MLVL,2
MHT(M)*WHT( (M*l)/2)
12 CONTINUE
••CALCULATE M HEIGHTS FOR MID LEVELS
DO 14 M«2iMLVL»2
MHT(M)=(MHT(M-I)*MHT(M*1) ) /2
H CONTINUE
••DETERMINE M LAYER FOR LB
IF(LB.GE.MHTtl)) GO TO IB
GO TO 30
18 DO 20 M»2.MLVL
IF(LB.LT.MHT(M» GO TO 30
20 CONTINUE
21 XAW399.
YAW*99.
XSIGMAX99.
YSIGMA=99.
RETURN
••CALCULATE WEIGHT SUMS AT LB
30 MLBoM
BTH=MHT(M)-LB
Mlz(M»l)/2
STHXW»STHXW*BTH»XW (Ml )
STHYW»STHYW*BTH»YW (Ml )
STH=STH+BTH
••DETERMINE M LAyER FOR LT
IF(LT.LE.MHT(MLVD) GO TO 48
MM-MLVL
60 TO 70
MLVLMI*ML .. .
DO so M«I,MLVLMI
48 MLVLMl«MLVL.l
MM*MLVL«M
IF(LT.GT.MHT(MM)) GO TO 70
50 CONTINUE
XAW'99.
YAW*99.
XSIGMA»99.
YSIGMA=99.
RETURN
••CALCULATE WEIGHT SUMS AT LT
70 MMLTM1=MM-1
TTH=LT-MHT(MM)
MM2=(MM»^)/2
STHXW»STHXW*TTH»XW(MM2)
STHYW=STHYW»TTH«YW(MM2)
STH=STH*TTH
IFtMLB.GT.MMLTMl) GO TO 300
••CALCULATE INTERMEDIATE WEIGHT SUMS
DO 100 M»MLB,MMLTM1
111
-------
TH»MHT(M*i)-MHT(M)
M2«(M*2)/2
STHXW=STHXW»TH«XW-MHT(M)
XSlGMA=XSlGMA*TH»(XW«M*2)/2)-XAW)««2
YSlGMA»YSl6MA*TH»(YW((M*2)/2)-YAW)«»2
101 CONTINUE
IFIXSIGMA.LT.O.O.OR.YSIGMA.LT.O.O) 30 TO 500
IFtSTH.LEtO.O) GO TO 500
XSIGMA»SgRTUS10MA/STH)
YSIGMA*SQ«t(YSIGMA/STH)
RETURN
500 XAW-99.
YAW'99.
XSIGMA»99.
YS!9MA«99.
ACTURN
ENO
POUTI**
)
IFjIJTAG.EQ.*) GO TC
:OND CLOSEST TIME, w
IF(ITABFB.EQ,2) GO
1* TRAJECTORIES
IF(MOD(J,2),EQ.O) GO TO 30
GO TO 10
SUBROUTINt CLSTMU, KK, MM, JTAG, NTPDY, ITABFB, M, NDYDT*. NMM,
2 KKT, «MT. IRTNS)
••CLOSEST TIME FOR ALTERNATE DATA
IRTNS*86
IFjIJTAG.EQ.4) GO TO 8
••SECOND CLOSEST TIME, JTAG«2
>F.ORIF(ITA|FB.EQ42, GO TO 3
2
•• BACK TRAJECTORIES
3 IF(MOD(Jt£>«EQ.O> GO TO 10
GO TO 30
. IRD CLOSEST TIME, JTAO*4(OBS ONLY)
8 IF(ITABFB.EQ.2) GO TO 2
GO TO 3
10 KKT*KK*1
IF (KKT.LE.NTPDY) GO TO 20
KKT»1
MMT»MM*1
GO TO 60
20 MMTsMM
GO TO 60
IF(KKT4SE.1> GO TO 40
30 KKT-KK-1
IFIKKT.GL
KKT»NTPDY
MMT*MM-1
GO TO 60
40 MMT-MM
60 IF (MMT ,LT, 1. .OR. MMT .GT. NMM) RETURN
IF(ITABFB.£5.2> §0 TO 62
IF (M * MMT -1 .GT. NDYDTA) RETURN
IHTNS»90
RETURN
62 IF (M«MMT-l .LT. NMM) RETURN
IRTNS»90
RETURN
END
SUBROUTINE DISOfJ, K, M, KK, MM, MXNSTA, NSTA, ITABFB,
A ITS, ITPOY, IMM. RADIUS, TSDuR. MBIP, MEIP, SLATT. SLONT, TSLATi
A, TSLON1, XAWT, YAWT, XXSIG, YYSIG, KKT, MMT. SJTAG. JTAGf ICNTMNI
A ICNTMX, XIP, YIP, XTSIP, YTSIP, DXSIGP, DYS1GP, HIP, TSLAT2,
A TSLON2, ISNAP, WIOROT , UBART, SIGMAH, ICNTT, SICNTT, SPRED1,
A SPRE02. ALPHAA, BETAA, TA, R2A, IJKA, IRTNSI
112
-------
••DISPLACEMENT CALCULATION USING OBSERVED WINDS
LOGICAL SJTAG(5).SICNTT
DIMENSION NSTA(ITPDY.IMM)
DIMENSION XXSIGIMXNSTA),YYSIG,YAWT(MXNSTA>
DIMENSION SLATT(MXNSTA),SLONT(MXNSTA)
INTEGER XlPI,NSTADO
SLATT=77\
YYSIG(I)«77.
96 CALL ITSIM ( JTAG, RADIUS. SLAT < I ,KK,MM) . SLON ( I .KK.MM) ,
, SLONT(I), XAW ( I ,KK,MM) , YAW ( I .KK.MM) , XSIG
. KK. MM), XAMT(I). YAWT < I ) , XXSIG ( I ) , YY
WI« IORT
IF (IJRTNS ,EQ. 150) GO TO 150
SLATT(I), SLONT(I), XAW ( I ,KK,MM) , YAW ( I .KK.MM) , XSIGMA ( I ,KK .MM) ,
Z YSIGMA(I. KK. MM), XAMT(I). YAWT < I ) , XXSIG ( I ) , YYSIG(I), TSDUR,
Z TSLATli TSLONl, ITABFB, NEAR, XI, YI, XTSl, YTSI, DXSIGI, DYSlGI.
Z WI« IORTNS)
IF (IJRTNS ,
ICNTT"ICNTT»1
XIP(ICNTTI«XI*,5
YIP(ICNTT)"YI».5
XTSIP(ICNTT)«xfsi/3.*.5
YTSIP (ICNTT)»YTSI/3.».5
DXSIGP(ICNTT)»DXSIG!/3.«0.5
OYSIGP(ICNTT)=DYSIGI/3.*0,5
W IP (ICNTT)=WI» 1000000.
SWXTSI=SWXTSI+WI»XTSI
SDXSIG»SOXSIG*Wl»QXSjGI
SDYSIG*SDYSIG*Wl»DYSIGI
SMl'SWUHl
150 CONTINUE
IFdCNTT .LT. ICNTMN .AND, NEAR ,NE. 1) RETURN
••NOT SUFFICIENT fiflSERVED WINDS
••TRAJECTORY SEGMENT DISPLACEMENT**************************************
IF (SWl .EQ. 0.) SWI=.0001
XTS=SWXTSI/SWl
YTS=SWYTSI/S«I
XSDX=SOXSIG/SWI
YSOY«SOYSIG/SWI
XTSP=XTS/3.*,5
YTSP=YTS/3.*,5
TSLAT2»TSLAT1«YTS/60,
TSLAT2»TSLAT1«YTS/60,
TSLON2«TSLONl-XTS/(60.*COS(TSLAT1*P1/180.))
UBART«,515»SQRT(XTS»XTS*YTS«YTS)/TSDUR
UBARTT«UBART/0,515
SIGMAH«(AdS(XTS)*YSDY*ABS(YTS)»XSOX)/(UBARTT»TSDUR)
SPRED2SSPRED1»SIGMAH
CALL XSZ(TSLAT1.TSLON1.TSLAT2.TSLON2»XTS.YTS,SPRE02.ALPHAA,BETAA,
IJKA.R2A,TA,WIOROT)
lCNTT«SjtAG(JTAG*l)
A
S
IF (M ,GE. MBIP .AND. M ,LE. MEIP) GO TO 172
IRTNS»170
RETURN
172 IFUSNAP ,EO, 1) GO TO 175
••PRINT TRAJECTORY SEGMENT DISPLACEMENT ANO»»»»»»»««»»*»»*»»»**«»«»««««
••INDIVIDUAL STATION RELATIVE L-OCATIONS AND DISPLACEMENTS**************
113
-------
IF(ISNAP.EO.l) GO TO 175
ISNAP-1
173 FORMATUH »i j « M TAG CNT XTS vrst/
A 1H ,i XI YI XTSIYTSI Ml',
A i XI YI XTSIYTSI MI XI YI XTSIYTSI wl't
A « XI YI XTSIYTSI HIM
175 WRITE(6,1HO) J,K,M, JTAG, ICNTT»XTSP»YTSP,
A (XlP(IP),YlP,XTSIP(IP),YTSIP, XSI6AA(io0.4i8), YSlGMAtiOOt 4
COMMON /OATERO/IOMO.IOYR.tDOY.IOHR.NRSTiiNftEc
DIMENSION SI (100), SLA JlOOtSLO (100), XA(l60),Y
XAM (
4,8)
lM
ilM
XA(l60),YA( 100) tXS(lOO),
IF(M.GT.l) 00 TO 5
YS(IOO)
GT.
..
MRIT£(6!4)
FORMAT (1H ,27X,» NUMBER OF NUMBER NUMBER OF»/
B IM ,27X,»REPOHTIN6 OF STATIONS IN*/
STAT10NS "ECOROS gOUN^ARIES t ) '
S DO 200 MM" 1, NMM
GO TO (l,2)!lTABFB
1 IF(M.Ea.I) 80 TO 8
1FXAW(I,K,MM»1)
YAW(I.K.MM)«YAW(I,K,MM*1)
XSlGMA(i,K,HM)«XSIGMA(I,K(MM»l)
YSIGMAU,K,MM)*YSIGMA(I,K,MM»1)
YSI8MA(I
9 CONTINUE
10 CONTINUE
GO TO 200
8 DO 199 K»1»NTPDY
A
LL POSJP(M,MM,K)
D0~l6i T«lil66
ALL RDAVO* ( MXNSTA, DTABT, DTABB. DTABL* DTABR, LBAAT, LTAAT.
OLAT,OLON,SI,SLA,SLO,XA,YA,XS,YS,NSTA(K,MM))
o 101 v • ---
SID(I.>
5LAT(Ii
5LON(I,
SLAT(I»K,MM)aSLA(|)
SLON(I,K,MM)»SLO(I)
XAW(I,K,MM)*XA(I)
YAWil.K,MM)«YA(I)
XSlGMA(I,KfMM)*XS(I)
YSlGMA(i,K,MM)»YS(I)
101 CONTINUE
199 CONTINUE
200 CONTINUE
RETURN
END
SUBROUTINE I-NPUT(NTPDY. TSDUR, RADIUS, ICNTMN, ICNTMX, MBIP, MEIP.
I !"!XNO,..Np, OID, OLAT, OLON, WlOROT, IHRINP, NLVLA, NTSPDY,
NTS, NMM)
•INPUT INFORMATION
114
-------
DIMENSION OID(MXNO).OLAT(MXNO).OLON(MXNO)
DIMENSION wioROTisi
COMMON/INPUT1/IBOY»IBMO,18YR.NDYtTOIR,
A POLLAT,IHOUR,1MO»
A NDYDUR,NDYDTA,WTYPE,DTABT,DT4BB,DTABL,DTABR,
A LBAAT.LTAAT.ALATT,ALATB,ALONL
»•»««••••»»»«•••»»»»•«««•••••••••••»«•»•«•••••»••••••••••••••*••••••••«
••INPUT DATA CARDS*************************** •****•**•*********•»••••*«
••CARD 0*»
••NUMBER OF ORIGINS
READ(StlOO) NO.NDY
100 FORMAT(I2ilX,I2)
••CARD 1»»
••ORIGIN IDENTIFIED BY LETTERS TDIRtNDYDUR
3 FORMAT(A4tlXiI2)
"NUMBER*OF DAYS OF WIND INPUT DATAJNDYDTAJ
READ(5.4) NDYDTA
4 FORMAT(I2)
••CARD 5»»
AND ANALYZED (O.A,
5
••CARD 6»»
••BOUNDARIES FOR OBSERVED WIND INPUT DATA
••TOP LATITUDE (DTABT)t BOTTOM LATITUDE (DT
••LEFT LONOITUOEJDTABL). RIGHT LONGlTyDE (DTAR>
. BEJp(5,6) OTABT.pTiBBjDTABL.DTABfi
OM LATITUDEtDTABBJlj,
GHT LONGITUDE
-------
MBIP=0
MEIP=0
••PRINTING iNTtRVAL IN HQUHS
,
••NUMBER OF ANALYSIS LEVELS
• •END INPUT PANAMETERS*«»**»«»««»»»*»»*«»**«*»«*»»»»»»»»«*»»«**»**»»»»*
«»«««««»<>««»•«•««*»««»•«»«»««»««»»»««»«<»«««•**««««*«»»•«««•«««*«»«»»•»«
••NUMBER OF TRAJECTORY SEGMENTS PER DAY
NTSPDY=24/TSDUR
••NUMBER OF TRAJECTORY SEGMENTS
••NUMBER*OF DAYS IN DATA BLOCK (DTABK)
NMM»NDYDUR*1
• •PRINT^INPUT^iJATA CARD INFORMATION
[NT U
WRITE
(Oini 10)«OLAT(IO),OLON(IO),IOal,NO)
WRITE(6,2D IBOY,IBMO.IHYH.NDY,
3 TOIR.NnYDUK,
4 NOYDTA,
I DTABTlt)TABB,DTABL«DTA8H,
208FORMAT(1H1,16XA?INPUTAFOR TRAJECTORY COMPUTATIONS*/
A1H ,«INPUT«/1H ,» DATA*/1H ,i CARD*/
11H 4X.41H1 ORIGIN »••«»••••»•••*»••«•*«•»«••*•««« «
A A3XM'.F4;itlX,Fft.l,.)./10(46X,A3,»<«,F4.1,lX.F6.1,«)'/»
C~TE THAT COMPUTATIONS BEGIN ••»•••••» ,12,IX,A3,IX,I
B1H ,4X,41H NUMBER OF DAYS COMPUTATIONS DESIRED »• »I2/
31H J4X,41H3 DIRECTION IN TIME ••«««•«•••••• A4/
AIM ,4X»41rt DURATION IN DAYS «••• ** «12/
41H !4xt41H4 NUMBER OF DAYS OF WIND INPUT DATA •••• ,I2/
51H I4X.41H5 TYPE OF WIND INPUT DATA •«•••••••••«•• ,A3/
61H ,4X,41H6 BOUNDARIES FOR OBSERVED WIND INPUT DATA/
A1H j4x!41H TOP AND BOTTOM LATITUDES ••••••«•••••• ,2X,F4.1,
SIM ^X^I^'LEF^AND RIGHT LONGITUDES ••» F6.i.u,F6.i/
2 s ;tx:^tt mir^l m ffil&M&ii i4,ix,i4,
~-15 - 4-- ALATT,ALATB,ALONL
BOUNDARIES FOR MAPS IN SUBROUTINES/
AIM i4X,41« TOP LATITUDE •••»•••»••••••»•••••»•*•• .2X.F4.1/
81H iJxUlM BOTTOM LATITUDE «•«••••••••••••••»••«* ,2X,F4.1/
C1H j4x!41H LEFT LONGITUDE «••«•»•••••••• •• ,F6.l>
WRITE(6,3U)
30 FORMAT(lHl)
RETURN
END
SUBROUTINE ITSIWUTAG, RADIUS, SL^'vlkrr' 5femT,D* T^T! *Tcl'om-*
z XSIGMA, YSIGMA, XAWT, YAWT, xxsis, YYSIG, TSOUR, TSLATI, TSLONI,
z ITABFH! NEAR, x, Y, XTS, YTS, DXSIG. DYSIG, w, IJRTNS)
LAfSH-JCATf 6° T° 3
ILONSH-SLONT
3 SLATSH-SLAT
L*T
LON
H)«60.»COS(SLATSH«PI/180,)
***0-
%
11 In5fAl:§S?5^§;To^R''1
XAMSH»XAWT
YAWSH«YAW!.
XXSIGS»XXSIG
YYSIGS»YYSIG
60 TO 18
20 XAWSH«XAM
116
-------
YAWSH-YAW
XXSIGS-XSIGMA
YYSIGS-YSIGMA
18 GO TO <21»22),ITABFB
21 XTS»1.9*«XAWSH»TSDUR
YTS«I.94«YAWSH«TSDUR
DXSIG«1.9*«XXSIGS»TSOUR
OYSIG"1,9*»YYSIGS«TSDUR
GO TO 23
22 XTS«-1.94«XAWSH»TSOUR
YTS*-1.94*YAWSH«TSDUR
OXSIG«-1,94«XXSIGS»TSDUR
DYSIG«-i.*4»YYSlGS«TSDUR
23 CONTINUE
DTS«SQRT(XTS»XTS*YTS«YTS>
XW=-X*XTS/2,
YW« Y-YTS/2,
OWSQ)>XW»XM»YW«YW
IF(OWSQ.EU.O.) OWSQ=.OOOl
DISTW-l./DWSQ
AI_ I Ny • 1
IFID ,N£. 0. .AND. DTS .NE. 0.) AUINW«1,-,5«ABS((YTS«X-XTS«Y)/
(DTS»0* .0001))
W«OISTW«ALINW
IJRTNSS0
REJURN
w
U.LVL)
XW(LVL)=-«SPD»S:N(WC .
YW(LVL>5-*SPO»COStw[>lR»Pl/180.)
50 CONTINUE
IF(LBAAT.NE.o) GO TO 60
XI«(OLON-iLON(I))«fiO.»COS(SLAT(I)«PI/lflO.)
YI»(SLAT(I)-OLAT)»60.
OI»SQRT(Xl*XI«YI*YI)
60 CALL AwlNQ ( SHT, SATH, DI, LbAAT, LTAAT, *HT, XW, YW. NLVL.
Z XAW(I), YAW(I), XSIGMA(I), YSIGMA(D)
loo CONTINUE
110 «RITF. (6,120) MSTA
'20 FORMAT(1H«,54X,13)
KE'IJRN
END
SUBHOUTINE ROl)ATA(I)
COMMON /iNuAT/onAr ISIOJ.IEND.NXTSTA.LSTSTA
isw=o
100 CALL TAPLX (iwINDUATA7,l,l,IOST,blO.OOAT,LEN)
iFnOST.EQ.O.OR.IOST.EQ.*) GOTO 1J5
IF(IOST.EQ.l) GOTO 116
GOTO 117
116 IF(ISW.NE.O) GO TO 98
GO TO 100
117 *RITE<6,2002)
3002 FOHMATO FR RDDATA. I-O,02»)
115 ISW«0
IF (LtN.EQ.510) GO TO 150
IF (LEN.EQ.510) GO TO 150
WRITE (6,2003) LF.N
2003 FORMAT!' FH RODATA, BAD LENiIlO)
150 IaI-170
NXTSTA*NXTSTA-170
GO TO 99
98 IEND=10
99 RETURN
EM)
SUBROUTINE RDDATE
FOR CSU TAPES
COMMON /DATERO/IDMO,IDYR»IDDY»IDHR,NRSTA,NREC,IM
COMMON /INUAT/ODAT<5loi,IENO,NXTSTA,LSTSTA
DATA 1ZEK/OX
ISW»0
CALL TAPtX( • «i I MUD AT A i ,1»1,IOST,510»ODAT,LEN)
IF(IOST.EQ.O.OR.IOST,E0.4) GO TO lib
IF(IOST.E'l.l) GO TO 116
GO TO 11r
117
-------
116 IFUSrt.NEtO) GO TO 98
ISW*10
GO TO ino
117
FORMATC FR RDOATE, I-0«,02)
lib ISW=0
EN
IF(LEN.<5T.20> GO TO 100
IDMO=OOAT(1)
IDDY=OOAT<3)
IDHRsQOAT (4)
NRSTAsOOAT (5)
NRECaODAT (6)
IM=ODAT(7)
NXTSTA=17l
LSTSTA=171
GO TO 99
98 IENO=10
99 RETURN
END
SUBROUTINE RDSTHD
COMMON /STHORD/ISTNOiALAT,AUON,ISTH,ISATH,NLVL
COMMON XI'MUAT/ODAT(510) , lENDtNXTSTA ,I_STSTA
DIMENSION OAT(3iI70)
EQUIVALENCE (ODATtn .DAT
-------
INTEGER SIIMMXNSTAI ,SIOT»SHT,SATH
INTEGER W*T(50)
DIMENSION SLAT(MXNSTA) , SLON (MXNST A )
DIMENSION XAw(MXNSTA),YAW(MXNSTA)
DIMENSION XSIGMAIMXNSTAI.YSIGMAIMXNSTA)
DIMENSION xw (bco ,vwtbo>
DIMENSION MTH(ja)
COMMON /DATERU/IUMO, IDYW, IDOY? IDHH,'\IWSTA,NREC, IM
COMMON /STHURO/ISTNO»ALAT»ALON, I STh , I _ATH, NLVL
COMMON XINOAT/ODAT (510) i I END ,NXTST A ,LSTST A
COMMON /WINDRD/HDS (3«5o>
DATA PI/3.14159/
UATA MTH /•JAN«f«FEb«,iMARi,iAPR',»MAY«. • JUN« , « JUL' , » AUG ' , ISEP« ,
A ,*
EQUIVALENCE (ODAT ( 1 ) ,DAT ( 1 ) )
IST=LSTSTA*2
IF(IST,GT.170) CALL RDDATA(IST)
IENsIST*
HDS(2,K)«OAT(a,I)
HDS(3,KisUAT(3,I)
10 CONTINUE
IF (JEN.EQ.IEN) GO TO 99
CALL ROOATA(IEN)
DO 20 1=1, IEN
K=K*1
HDSd ,K) sOAT( 1,1)
HDS(?,K)sUAT(2,I)
HDS(3,K)»UAT(3,I)
20 CONTINUE
99 RETURN
END
119
-------
SUBROUTINE XSZ(TSLAT1tTSLONl,TSLAT2,TSLON2.XTS,YTS.SPRED2,ALPHAA
A ,BETAA,IJKA,R2A,TA,WIOROT)
DETERMINES THE DISTANCE,R2,AND THE ANGLES,ALPHAA.BETWEEN TRAJECTORY
END POIMTS» AND THE DISTANCE TA, AND THE ANGLE,BETAA, BETWEEN TRAJECTORY
ENDPOINT AND SPREAD END POINT— ._..—....
PIB3. 14151*
WORO»WIOROT/1852.0
R2-XTS«XTS*YTS»YTS
R2A=SQRT(R2)
S2«(WORO/2,»SPRED2)«(WORO/2.*SPRED2»
T2«R2»S2
TA»SQRT(T2)
iF(R2,EQ,0.) GO TO 101
IFtXTS.EQ.O.) 60 TO 1
F(R2.EQ.O.) GO TO 1
F0,
?F(t§lJf2tGE.TSLATl.AND,TSLON2,LE.TSLONl) GO TO 10
F(TSLAT2.GE;TSLATl.ANDlTSLON2jGT.TSLONi) GO TO 20
IF TSLAT2.LT.TSLAT1.AND.TSLON2IGT.TSLON1) GO TO 30
IF TSLAT2.LTITSLAT1IANDITSLON2ILE.TSLONI) GO TO 40
10 |JKA»1
ALPHAA«ALPHAA
RETURN
20 IJKA»2
RETURN
30 IJKA=3
ALPHAA-PI*ALPHAA
RETURN
40 IJKA*4
ALPHAA»2,»PI-ALPHAA
RETURN
101 ALPHAA-0.
BETAA-0.
•A«0
URN
Program of Concentration Calculation (Model A)
PROGRAM TRAJET
. ..__ .-j_, -IN DATA COMPUTED FROM THE PREVEOUS
TO COMMUTE TM£ RESULTANT CONCENTRATIONS ETC.
COMMON /OUTP1/ TSLAT(20,4,8), TSLON (20,4,8) , UBAH(20,4,8)
, SPREAO(20,4»8)
COMMON /OUTP2/ ALPHA(20,4,8)( BETA (20,4,8) , IJK(20(4,8)
, T(20,4»8). R2(20,4,6)
COMMON/6RIO/BLAT ( 31 ) ,8LON(31) , CONG (31,31 ) ,SCONG (31,31 ) ,XXJ (31,31 ),
^1?^
i^ «£s?:
NSION OID(
0 I * 1, 5
REAO(14, 1
IMENSION OID(5),OLAT(5),OLON(5»
DO 2
-- -,~ -. 10) OID(I), OLAT(I), OLON(I)
10 FORMAT (A3, F4,l, F6.1)
20 CONTINUE
READtU, 30) NO, NOY, NTS, NTPDY, LTAAT
30 FORMAT(SIlO)
READ(14, 40 ) TSOUR* ALATT, ALONL
40 FORMATOF15.5)
00 500_IO"1,NO
- •••- - 'SO
IFINO.EQ.
WRITE(6,4
416 FORMAT(IHJ,A3,«
) GO TO 450
6) 010(10),OLAT(IO)»OLON(IO)
KUNIT-20*
DO 425 M-
0
,NDY
DO"420 >' ^PP?
0 J • 1, NTS
READ (KUNIT1 TSLAT(J.K.M), TSLON(J,K,M), UBARIJ,K,M),
SPREAO(JtK,M). ALPHA(JfK,M). BETA(J,K,M ), T(J,K,M), R2(0,K,M),
, l{NT(J,K,M). SICNT(J.K.M)
120
-------
420 CONTINUE
422 CONTINUE
425 CONTINUE
450 CONTINUE
CALL CONCAL (IO,NO,OLAT,OLONt NTS
A .NTPDYiNDY,TSOUR,ALATT,ALONL.LTAAT)
SUBROUTINE CONCALCIOINO.OLAT.OLON,
ANTSiNTPDy.NDY.TSDURtALATTiALONL.LTAAT)
-CALCULATE THE CONCENTRATIONS ON GRID POINTS--—- ---- -——.--
COMMON /OUTPl/ TSLAT(20f4,8) , TSLON (20 .4 ,8) , U8AR(20»4,8)
Z , SPREAD(20,4,8)
COMMON /OUTP2/ ALPHA(20t4,8) , BETA (20 |4 1 8) , IJK(20i4,8)
Z , T<20,4}8), R2(20,4,8)
COMMON/GRlD/6LAT(31) iBLON(31) »CONC> < 31 • 31 > ,SCONG ( 31 . 31 ) »XX-J(31f 31)
Z SXXJ(3li31)tDXXJ(3l,31).DSXXJ(3i,31)
DIMENSION VG1 (20.4) ,VG2(20,4) ,CONVR(20t4) ,H(20,4) ,RAMOA1 (2of4) .
Z RAMDA2(20,«)
INTEGER TSLAT.TSLON
DIMENSION WOR(20,4,8), CONC (20 ,4, 8) , SCONC ( 20 t 4, 8)
DIMENSION GAMMAISLS!) ,c(3i,3i)
DIMENSION o«coN(5> ,SORCON
DIMENSION OLAT< SI.OLONI s>
DIMENSION^INDEX(31,31)«INDEX2(31»31)
00 48oO Kal,NTPnY
DO 3000 Jrl,NTS
CALL REMOVE (K, J.VGA, VGB,TRAR,HGTtRAMAAtRAMBB,LTAAT)
VG1 (0,K)=VGA
CONVR( Jt K)=TRAR
H(J,K)«HGT
RAMOA1 ( J,K)=RAMAA
RAMOA2(J,K)iRAMR8
3000 CONTINUE
4000 CONTINUE
READI5.7T) ORCON(IO) .SORCON ( 10) , WIORO ( 10)
77 FORMATI3F10.1)
*ORO=WlORO(IO)
DO 87 M = 1, 8
DO B6 K = 1 , 4
00 85 J * 1, 20
*OR(Jf K, M)aO.
CONClJt K, M)ao.
SCONCIJt Ki M)aO.
85 CONTINUE
86 CONTINUE
87 CONTINUE
IF (IO.NE.U GO TO 101
DO 100 TI»1,31
DO 99 JJ » 1, 31
III=II-1
JJK*JJ-1
BLAT(II)=ALATT-0.5»III
BLON(JJ) »ALONL-0.5«JJK
CONG(II|JJ)»0,0
XXvJ(II»JJ>»0.0
SCONG(II.JJ) =0.0
SXXJ(II,J>J)=0.0
DXXJ(II,JJ)=0,0
DSXXJ(IItJJ)«0.0
99 CONTINUE
100 CONTINUE
101 CONTINUE ,
00 800 MsliNDY
DO 700 K«1«NTPDY
DO 600 J=1»NTS
REMVl^VGl ( J»K) *RAMDA1 ( J,K)«H(JfK) «CONVR ( J . K ) «H ( J ,K )
REMV2=VG2(J.K) *RAMDA2 ( J,K) «H < J,K)
DREMV=REMV1-REMV2
VG1H=VG1 (JtK)/H(J,K)
VG2H»VG2(J»K)/H(J,K)
DA»EXP(-VG1H«TSDUR*360Q.)
OB=EXP ( -V62H*TSDUR«3600. )
DC=EXPJDREMVH»TSDUR«3600.0>
DPlaEXP(-RAMDAl (JtK)»TSDUR«3600.)
DP2«EXP(-RAMDA2(J«K)«TSDUR»3600.)
DCON=EXP(-CONVR(J,K)»TSDUR»3600.)
SWIDTH»WO«0*2.»1852.«SPREAD(J,K.M)«60.
IF(J.NE.l) GO TO 601
WOR(J,K,M)»WORO
CONC(J.KtM)a( (ORCON(IO)«WOR( JtKiM) ) /SWIDTH) »DA»DCON«DP1
SCONC(J,K»M)a (SORCON(IO)»WOR(J,KiM)/SWlDTH)«DB»DP2
A * ( 1,5«CONVR(J»K)«H< J,K ) /OREMV) »CONC (J»KtM)«DC
GO TO 600
601 WOR(JtKfM) =«OHO*2,«60.«l852.«SPREAD(J.lfK,M)
121
-------
CONC(J«KiM)=( (CONC(J-1.K.M)»WOR(J,K,M) ) /SWIDTH) «DA»DCON«DPl
SCONC(J»KtM)«( (SCONC(J-1»K,M)»WOR(J,K,M) ) /SWIOTH) »DB»DP2
A * (1.5»CONVR1) ALPHA tBET A ( JiK |M) . IJK ( J,K ,M) » T (si i K ,M)
AiR2(J£K.M) t SPREAD! j,K,M) .UBAR ( J,K,M) .TSLAT ( J,K,M) , TSLON ( J tK , M)
651 FORMAT<2F10.3,I5,4F10.3,2I10)
IF ( J ,EQ. 1) 60 TO 599
IF (TSUATU, K, Ml/10. ,6T. BLAT(l) .AND.
Z TSLAT/10. .GE. BLAT(l) .OR.
Z TSLAT(0. K, MJ/10. .LT. BLATOll .AND,
Z TSLATlO-1. K, M)/IO. ,LE. BLAT (31 ) ) 60 TO 700
IF (TSLONUi K, M)/10. ,GT. BLON(l) .AND.
Z TSLON(J-1» K, Ml/10. ,GE, BLON(l) .OR.
Z TSLON(J» K, M>/10. .LT. 8LONO1) .AND,
Z TSLON1J-1, K, M)/10. .LE. BLON(31» GO TO 700
ON
599 CONTINUE ,
DO 10 II-U31
DO 9 JJ « 1. 31
) 60
ic
IF(J.EQ.l)"6o"fo 50
TTUAT»TSLAT(J-1,K,M)/10.
mON*TSLON(J-l,K,M)/10.
60 TO 51
50 TTLAT»OLAT(|0)
TTLON-OLON(IO)
51 CONTINUE
Gx«(TTLON-BLON(JJM»60,«COS(TTLAT»PI/lBO,)
GY»(BLAT(II)-TTLAT)«60.
IF(SLON(JJ).EQ.TTLON) GO TO 5
IFJBLATJID.EQ.TTLAT) GO TO 6
RATIO»ABS(GY/GX)
6AMMA1II.JJ)sATAN(RATIO)
IF.GE.TTLAT.AND.8LON(JJ),UE.TTLON) GO TO
IF(BLAT(II).GE.TTLAT.AND,8LON(JJ).GT.TTLON) GO TO
IF(BLAT(II),UT.TTLAT.AND.BLON(JJ).GT.TTUON) GO TO
iF(8LAT(Il>.LT.TTLAT.AND.8LON(JJ).LE.TTLON) GO TO
21 GAMMA(il«JJ)*GAMMA(II,JJ)
GO TO 30
22 GAMMA(II.JJ)ePI-GAMMA(Il,jJ)
GO TO 30
23 GAMMA(Il.JJ)*PI«GAMMA(II«JJ)
GO TO 30
24 GAMMA(II,JJ)»2.*PI-GAMMA(II,JJ)
GO TO 30
5 IF(BLAT(I1).GE,TTUAT) GO TO 4
6AMMA(It,JJ)*l.S*PI
GO TO 30
4 GAMMA(ll.JJ)iO,5«PI
GO TO 30
6 IF(BLON(JJ).LE.TTLON) 60 TO 7
GAMMA(IIiJJ)«Pl
GO TO 30
7 GAMMA(11,JJ)»0.
30 C(II.JJ>'SQRT(GX»GX«6Y*GY)
9 CONTINUE
10 CONTINUE
OELTA1«ALPHA»0.5»PI
DELTA2«ALHHA(J,K,M)»1.5»Pf
PH*Il«AUPHA(j,K,M)*B|TA(J»K,M)
INDEX2I,
CONTINUE
CONTINUt
DO 250 11*1.31
110
DO 246 JJ « I, 31
AIF(J.EQjI^ANO^BLAT(II).EQ.OLAT(IO).AND.BUON(JJ).EO,OLON(IO))
IF(J,EQ.1.AND,C(II.JJ).LE,WORO/(2.«1852.) . AND. GAMMA (II ,JJ| ,EQ,
A DELTA1.0R. GAMMA ( II, JJ) .EQ.DELTA2) GO T^ 355
EE » T(J. K, M)» T(J, K, M) - R2(J, K, M)«R2(J, K, M)
E«SQRT (EE)
F*R2(Jt K. M)»E
(II. JJ) ,6T. F
J!UJ,K,M):EQ:I)
JK(J,K,M).EQ.2)
JK(J.K,M).EQ.3)
JK(J,K,M),EQ.4) JO
FdJK(J.KtM) ,|Q 0) GO TO 140
F(BETA(J.KfM),LT,ALPHAjJ,K,M)) GO TO 111
F (GAMMA < I I, JJ) . GE. ALPHA (J.K.M) , AND. GAMMA (II. JJ) .LE.PHAI 1 ,OR.
i,jJK6E:<2.»PI»PHAI2).AN5.GAMMi(lTiJJ),L^,2.ip!,oA.
l.JJj.GTtd. . AND. GAMMA ( II ,JJ) .LT.AtPHA (S.KjM? ) * *
R2(Jt K. M)»E
F(C(II. JJ) ,6T. F) GO TO
F
-------
A «(1.5*CONVrt »CONC(J»K,M)»D6
INDEX(II,JJ)«2
GO TO 249
350 CONG(II,JJ)*ORCON(IO)
SCONG(II,JJ)*SORCON(IO)
INDEXlII»JJ)=2
GO TO 249
162 CONG(II,JJ)=CONC(J,K,M)
SCONG(II,JJ)=SCONC(J,K,M)
INDEX«2
GO TO 249
140 CONG(II,JJ)=0.
SCONG(II,JJ)aO.O
INDEX(II,JJ)=1
249 1F(INDEX2(II,JJ),£0,2.AND.INDEX(II.JJ).EQ.2) GO TO 251
GO TO 252
SCONGdl, JJ) *0.0
252 XXJtII»JJ>»XXJ(II,JJ)»CONG(lI«JJ>
SxxJ(ll,JJ)=«SXXJ(II,vJJ)»SCONG(II,JJ)
DCONG=CON6(II,JJ)«(V61(J,K)+RAMDA1(J,K))»H(J.K) «3, 6*6.0
DSCONG*SCONG>0,0
SCON6(II»JJ)«6.0
GAMMA(II,JJ)«0.0
C(IIiJJl*0.0
660 CONTINUE
661 CONTINUE
650 CONTINUE
IF(IO.NE.l) GO TO 700
CCCCCCC JJJ « 4 TIMES THE NUMBER OF DAYS CONSIDERED
JJJ=28
700 CONTINUE
800 CONTINUE
DO 802 11*1,31
DO 801 JJ • 1, 31
XXJ(II,JJ>»XXJ(II.JJ)/FLOAT(JJJ)
SXXJ(n,JO)«SXXJ(II,JJ»/FLOAT(JJJ)
801 CONTINUE
802 CONTINUE ,
WRITE(6,914)
914 FORMATUHU8X,' CONCENTRATION OF S02 «////)
WRITE(6,915) (BLON(JJ),JJ»1,16)
915 FOftMAT(8X,l6Ffl.l/////)
DO 903 II»1,31
WRITE(6,904) BLAT(II),(XXJ(II,JJ)tJJ«l.!6)
90* FORMAT(F8,l,l6E8.a////)
903 CONTINUE
«RITE(6.925) (BLON(JJ),JJ»17,31)
925 FORMAT(lHl,ax»l*F8.1/V///)
DO 913 11*1,31
WRITE(6,941) 8LAT(II),(XXJ(IIiJJ)tJJ=17,31)
941 FORMAT(F8.1,15E8.2////)
913 CONTINUE
WRITE(6,9i>0)
950 FORMATUHliSX, ' CONCENTRATION OF S04 •////)
WRITE(6,915) (8LON(JJ),JJ»1,16)
00 951 II*U3J
WRITE(6,904) BLAT(11),(SXXJ(11,JJ),JJ = 1,16)
951 CONTINUE
WRITE(6,925) (8LON(JJ),JJ»17,31)
DO 952 11*1,31
WRITE(6,941) 8LAT(II),(SXXJ(II,JJ) ,JJ«17,31)
952 CONTINUE
WRITE(6,953)
FORMAT(1H1,8x,'utposiiION or aoi
WRITE(6,915) (BLON(JJ),JJ«1,16)
953 FORMATUHl,8X,'DEPOSITION OF S02 (KG/KM2) • ////i
WRITE(6,915) (F - . . - - -
00 954 11*1,31
WRITE (6, 904) BLATl II) , (DXXJ(II.JJ) ,JJ=1,16)
954 CONTINUE
WRlTEf6,925) (BLON t JJ) , JJ»17 , 31 )
00 955 11*1,31
WRITE (6, 94 if SLAT (I I) , (OXXJ(II,JJ) ,JJ«17,31)
955 CONTINUE
956 FQRMAT(1H1,8X, 'DEPOSITION OF S04 (KG/KM2) •////)
WRITE(6,915) (BLON(JJ) ,JJ«1,16)
DO 957 11*1,31
WRITE (6, 904) BLAT(II> , (DSXXJ ( 1 1 » JJ> , JJ»1 , 16)
957 CONTINUE
123
-------
WPITE<6,925) ( BLON(JJ),JJ«17»31)
Do 958 II*1»31
WRITE(6»9*1) BLAT(II),(OSXXJ(11,JJ)»JJ«17,31)
958 CONTINUE
DO 920 11*1,31
DO 9J9 JJ = 1, 31
XXJ(IIiJJ>»XXJ(II,JJ)»FLOAT(JJJ)
sxxj(iiiJJ>
919 CONTINUE
920 CONTINUE
RETURN
END
110
120
130
SUBROUTINE REMOVE(K,J,VGAfVGBtTRARtHGTtRAMAA,RAMBB,LTAAT>
IF(K.E«.1> GO TO " "
lF.2> GO TO ...
lF(K.Ea,3 GO TO 130
!F(K.EO.*> GO TO 140
IF(J.LE,4) GO TO
IFulLE^e) GO TO
IFU.LE.l2) GO TO
IFIJiLEllb) GO TO
GO TO 111
IF(J.LE.2)
111
ttt
m
GO TO
111
112
GO TO 111
GO TO 112
GO TO 111
JFIJ.LE.6) GO TO
IF(J.LE.IO)
IF{J.LE,14)
IF(J.Lt.lB)
GO TO 112
lF(J.Lt.4) GO TO
IFIJ.LE.8) GO TO
IF(jItEll2) GO TO
IF(J.LE.16) GO TO
GO TO 112
112
111
112
111
140 IF(J.LE.2) GO TO 112
IFIJ.LE.6) GO TO 111
.LE.10)
• LE.14
IFIJ.LE.
GO TO
GO TO
GO TO
IF(J.LE.IB)
GO TO 111
VGAsQ.Ol
VGBsO.OOl
TRAR»3.0E-6
HGT=FLOAT
-------
APPENDIX - B
FORTRAN LISTING OF MODEL B
Trajectory Calculation Program (Model B)
PROGRAM TRAJET
1 (OUTPUT, DAT A 1,TAPFS=DATA1,TAPF6=OUT^UT,TAPE1=10"2H,TAPE2'
2 TAPE21 = l002b,TA°F22=10G/e.B,TAPF23 = 1002B,TAWE2<»=l 10?B,
3 TAPE25=lUO?b,TApE2t>=lreyH,TAPF2r=1002B,TAPE2fl=l 002H,
4 TAPE29= 100 ?B,T APT 30 = 100 2B, TAPE 31 = 100 2B, TAPE 32= tOO 2B,TAPt33=1002R,
5 TAPE3<*=lOO?B,TAPr35 = no 2*, TAPE 36 = 1 002B,TAPE37=tno2H,TAPt 38 = 10024,
6 T APE3<5= 1 00?H »TAPF*0 = innaB, TAPE 1^ = 1 00 RH)
DIMENSION F.hiSr>2(S) ft HSUL (S)
DIMENSION OIDCi) ,OLAT (b) ,OLON(5) , ImlOHOT (»
Z YAW (flO,2,5) .XSlGMAteO.^.S) , YSIGKft I 80 ,2 ,5)
COHMON/STAT2/SLA1?(80,?,5) , SLON2 ( 80 ,2 ,5) , XAW2 ( «0,2»5>,
Z YAW2<80.2,5> , *SIKM2<80,2,5) • YSI GM2 ( 40 , 2.5)
COMMON/STATS/SLAT'HtHO.^.S) , SL"N3 ( 80 , 2 ,5 ) , XA W3 (
Z YAw3(8Q»2,'S),XSIPM3(ao,2.S),Y?IGM3(BO,2»S)
C"»FORWARD AND HACR^'ARU TRAJECTORY MODEL
C*»USING OBSERVED ANO ANALYZED INPUT '.VINOS
LOGICAL SICNT1.SICMT?, SICNT3, SICNT4.SICNT
C»»FOR 08SERVEU INPUT IsIrJDS OfJLY«»»»»»«»»i>-»(**<»»»«»**»
DIMENSION SLATT (HP) ,SLONT(HO) ,XAWT(dO) ,YAWT{80)
DIMENSION NSTA(?i[;>iNSTA^(?,5).NSTA3(2i5)
DIMENSION XXSI<3(aO) lYYSlG(flO)
INTEGER SID(80.2)
INTEGER SID2(aO. 2) ,5103(80.2)
INTEGER TSLAT.TSL""J
INTEGER TTLAT1.TTI.ON1,TTLAT?.,TTLON3,TTLAT3,TTLON3,TTLAT4»TTLON4
C»«SEE DATA
LOGICAL SJTAGCS) .SICNTT
INTEGER XIP(3())iYTP(.3Q)tXTSIP(10),YTSIP(3n),WlP(30)
INTEGER OXSIGD(30) ,OYSIGP(30)
COMMON/ 1 NPun/I3DY,IBMO,lHYRtNOYtTDIRt
A POLLST,IHOUR,IhO,
A NDYOt/^,NnY01A,wl YPE,PTAB1 , DT ABB ,DTABL .OTABK,
A LBAAT,LMAAT,LTAAT,ALATT, ALATH,ALONL
COMMON/SOS T/TSI ATS, TSLONS , SpREDl ,
A SPREO,S^«EAU,ALPHA,riETA,T,R?,IJK
C»»MAXIMUM NUMBER OF STATIONS
DATA MXNSTA/ RO/
DATA OBS, ANA/3H03S,3HANA/
DATA BACK/4H8ACK/
DATA SJT*G/1H , 1H. , 1H+ , 1H» , 1H-X
DATA ITS«lTPOY«lMM/rt,2,5/
ISNAP=0
CALL INPUT ( NTPDY, TSDIJR, RADIUS, ICNTMN, ICNTMX, MrilP, MElP,
z MXNO, NO, nin. OLAT, OLON, WIOROT, IHRIMP, NLVI.A, NTSPOY,
Z NTS, NMM,ERSO?,EKSUL)
C»»SET ITAB FO«WAHD/HACKWAHD TO 1 FOR FORWARD, 2 F0« BACKWARD
1TABFB=1
IF (TOIR .EQ, PACK) ITABFB = ?.
C«»LOOP THROUGH NUHPF.:? OF DAyS
00 400 "=1,NOY
CALL DTA8KO ( M, MXNSTA, MTPDY, NMM, NST A,NST A2 ,NST A3 ,
A IT58FR,NOYnTA,
Z ITPOY, IMM, DTA8T, DTAB^, OTABL, DTABR, LBAAT, LMAAT.LTAAT,
Z OLAT lOLON, SI 0 , SID2, SID3)
C«»LOOP THROUG^ NUMBER OF DAYS
DO 399 10=1, NO
C
C
C»«LOOP THROUGH NUMBER OF TRAJECTORIES PER DAY
DO 300 N=l »4,NTPOY
125
-------
.l) K=l
IF(N.EQ.J) K = 2
K K T — K
IF (ITABK4 .EQ. GO TO 42
C««TRAJECTORV HAS BEEN TF.RHIMATEO
C««OUE TO INSUFFICIENT AMOUNT OF INPUT DATA-
C»»TERMINATION CODE INFORMATION
C»«REPL ACES ALL REMAINING TRAJECTORY SEGMENT LAT«LON VALUES
GO TO 20
C ,»«*««»««»»»«» INI T I AL 1 7 AT I ON«n*»»*«»»**»»*»*«»*«»»»*
=L.Af
GO TO 200
««»««»»»«»
40 TLATll=OL.AT(If»
TLAT12=OLAT(IO)
TLAT13=OLAT(IO)
TLAT14=OLAT
TLON12=OLON( 10)
TLON13=OUON(IO)
TLON14=OLON(IO)
SPPD11=0.
SPR02l=Ot
SPRD3UO.
SPRD41=0.
SPRE1=0.
SPRF2=0.
SP»E3=0.
SPRE4=0.
42 JTAG=0
JTG2=0
JTG3=0
JTG4=0
CALL KKMM { J, N, H, KK. MMt NMM, NDYDTA, NTSPDY » I T A«FB ,
IF (IRTNS .FQ. -JO) GO TO 9Q
C««TRAJECTORY TErtHINATEO OUE TO MO PATA READ INTO DATA-HLOCK (0 r ABK )
C»»OR CLOSEST-riME(CLSTM) CALCULATION OUTblDE DATA-BLOCK
^6 TTLAT1=9V6
TTLON3=9996
TTLON4s9996
UHAP1=0.0
U8AR3=0.0
U«AR4=0.0
SPREA1=0.0
SPREA2=0.0
SPREA3=0.0
$PREA4=0.0
ALPHA2=99.
ALPHA3=99.
ALPHA4=99.
8ETA1=99.
BETA2=
-------
SICNT1=SJTAG(1>
SICNT2=SJT4G l)
SICNT3=SJTAG 1)
SICNT4=SJTAG(1)
GO TO 200
90 IF(K.EQ.l) GO !'•> 1100
IF(K.EO.Z) GO TO 13QO
IF(J.GE.3.ANO.J.LF. 4) GO TO 1120
IFfJ.GE.b.ANO.J.LF.fa) GO TO 1130
IF(J.GF..7.ANO.J.LF.B> GO TO 1140
1110 CALL DIS°2(J,K,y,KK,MM,i*!XNSTA;NSTA2,ITABF8, ITS,
ZITPDY»TMM»RADIL'S,TSDUR,MBIP,MFIP,SLATT.SLONT, TLAT 1 1 . TL ONI 1, XAWT,
Z YAWT,XX3lG,YYSIG.KK2.MM2«SJTAG,JTG2,ICNTMN,ICNTMX,XlP,YlP,
Z XTSIP»YTSIP,DXSIfiPfDYSlGPfWIP.TSLAT?,TSLON2, I SNAP, W IOP.OT,
Z UBAHT,SlGMAH,ICNTT,SICNTT,SPRnil,SPRE02,ALPHAA,BETAA,TA,R2A,IJKA,
ZIRTNS)
IF2,ALPHAA,BETAA,TA,R2A,IJKA,
Z IRTNS)
IF (IRTMS«t.Q,170) GO TO 170
CALL ALT
-------
TLATl3=TSLAT5
TLATl^TSLATS
TLON11=TSLON5
TLON12=TSLON5
TLON14=TSLON5
SPRD11=SPHED1
SPRD31=SPKE01
SPRD41=SPKED1
SPRE2=SPHtD
SPREA1=SP«EAD
SPPEA2=SPHEAD
SPREA3=SPREAD
60 TO 200
87 TTLAT1=9?6
TTLAT3=996
TTLON1=9996
TTLON3*9996
UBAR1=0.
UBAH3=0.
SPREA1=0.
SPRtA3=0.
ALPHA1=S°]TAG°JTG?, TSLAT2, TSLON2 ..JBART .SIGMAH, ICNTT,
Z SICNTT.ALPHAA,HETAA,TA,q2A,IJKA,IHTNS)
iFdRTNSjEQ.]- -- -- '-•
141 C ------
171 CALL CONVEX(TSLAT2|TSLON?,UBART,ICNTT,SICNTT,SPRE1«SIGMAH,
Z SHRE02,ALPHAA,BETAA,TA,R2A,IJKA,TSLAT,TSLON,
.Z UHAR.ICNT.SICNT)
TTLAT1=TSLAT
TTLAT3»TSLAT
TTLONl»TbLON
TTLON3=TSLON
UBAR3allBAR
ALPHA1«ALPHA
ALPHA3«ALPHA
BETAl=flETA
BETA3»RETA
T1»T
R2l=R2
IJK1=UK
IJK3=IJK
ICNT1=TCNT
ICNT3*ICNT
SICNTUSICNT
SICNT3=SICNT
TLAT11=TSLAT5
.
TLONIUTSLONS
SPRE3=SP«tn
SPREAUSPREAO
GO TO 1321
88
128
-------
TTLON2=9996
TTLON4=9996
UBAR2=0.
U0AR4=0.
ALPHA2=99.
ALPriA4=99.
BETA2=99.
BETA4=99.
T2=99.
T4=99.
R22=99.
R24=99.
IJK2=0
IJK4=0
ICNT?=0
ICNT4=P
GO TO POO
93 CONTINUE
1321 CALL DIS03tJ,K,M,KK,MM,MXNSTA,NSTA3»ITABFB, ITS,
ZITP[>YiIMM»RADJIJS,TSDURt.18IP,MEIP,SLATTtSLONT,TLATl2tTLONl2.XAWTt
Z YAWTtXXiIG,YYSlG,KK3»M>n»SJTAG, JT03, ICNTMN , 1CN TMX , XIP t Y t P,
Z XTSIP,YTsiP,OXSIGP,QY5IGP,'-iIP,TSL4Ta,TSLON2,ISNAP,WlOROV,
Z UBART.SlGMAHi JCMTTiSIChlTT,SPRD2l»SP«EOafALHHAA,BETAA,TA,R2A,IJKAt
IF(IRTNS.EQ.170) GO TO 172
CALL *LTOTA( SJTAG,jrG3. TSL AT2, TSLON2.UBART, SIGMAH, ICNTTt
Z SlCNTT«*LPHAA,3c;TArt,TAtR2A,IJKA,IRTNS)
IF(IRTNS.EQ.17fl) GO TO 172
CALL ClSTM< J,KK,MM,JTG3tNTPDYfITABFHtvi,NDYDTAtNMM,KK3»MM3,IRTNS)
IF (IRTNS.tQ.90) GO TO 93
IF( IRTNS.tQ.86) GO TO 88
172 CALL COKlVEK(TSLAT?fTSLON?tUBARTiICNTT«SlCNTT,SPRF.2»SIGHAH,
Z SPRED2,ALPHAA,8ETAAiTA,R2A, IJKA.TSLAT. TSLON,
ZTTLAT2=TSUAT U.AR , ICNT, SICNT,
TTLON2=TSLON
ALPHA2=ALPHA
ALPHA4=ALPHA
BETA2=PETA
BETA4=8ETA
R22=P2
R24=R2
IJK2=UK
IJK4 = UK
ICNT2=ICNT
SlCNT2=SfCNT
TLAT12= T5LAT5
TLAT14= TSLAT5
SPR021=SHHEDl
SPRD41=SPHED1
.
GO TO 200
89 TTLAT1=996
TTLON1=9996
TTLON3=9996
UBARl^f).
UBAR3=0.
SPREA1=0.
SPREA3=0.
ALPHA3=99.
BtTAl=99.
BETA3=99.
Tl=99.
T3=99.
R21=99.
R23=99.
IJK1=0
IJK3=0
ICNT3=n
SICNT1=SJTAG(1 )
SICNT3=SJTAG(1)
fiO TO 13J1
94 CONTINUE
1330 CALL OISO ( J , K ,M ,nK , MM ,MxNST A^-iSTl ,ITArtFB, ITS.
Z I TPP Y.I MM. RAO I US, rSi)UR,i1UIP,Mf TP,SLATT,SLONT,TLAT11 ,TLONll ,XAwT,
129
-------
7 YAWT,XXi>IG,YYSir,,KKT.MMT,SJTAr;,JTAC>, ICNTMM,ICMT«X,X IP»Y1 P.
7 xTSIP,YTSTP,nxsir.P,OYSlGP,^I",TSLAT?,TSLON?,t=;MAP,^IOflOr,
I USAHT,SlGl5AH,l2riTT!stC'VrT;SPRRll,SPRtD2,ALPMAft,HETAA,TA,«?A,IJKA,
7 Th?TNSl
, TSL ATP> , TSLON2 .U8ART , Si G*AH , ICNTT ,
2 SICNTT,flLPHAA.8F.TAA,TA,M2A,IJKA,lWTNS)
,
IF (IRTNS. EG. 90) GO TO 9ft
173 CALLRTNS'tQCONVE^TSLA?I,TSLON?,UHART, ICNTT, bICNTT.SPREl .SIGMAH,
lfJL«LL v.u >5pwE0^;ALpHA«,dt
7 'JHAR.I'CNT.SICNT)
TTLAT1=TSLAT
TTLAT3=T3LAT
TTLOM3=T3LON
ALPHA1=ALPHA
ALPHA3=ALPHA
BETAl=RErA
BETA3=RETA
T1 = T
T3 = T
R2l=M2
R23=R2
ICNT3=ICNT
SICNT1=S£CNT
SICNT3=SlCNT
TLATll=TSLATb
TLAT13=TSLATS
TLON11=TSLON5
TLON13=TSLON5
SPRD31=S^HE01
GO TO 1331
61 TTUAT2=996
TTLON2=9996
TTLON4=9
UBAR2=0.
SPREA2=0.
SPREA4=0.
BETA2=99,
BETAft=99.
T2=99.
T*=99.
R22=99.
R24=99,
IJK2=0
IJK4=0
ICNT2=n
ICNT4=0
SlCNT2=SJT4fi (1)
SICNT4=SJTAG( I)
GO TO 200
95 CONTINUE
1331 CALL OIS" (JtK,M,KKiMM,MxNSTA,NSTA ,ITARFB, ITS.
ZITPOYtTMM,RftOIl|S,TSOUH,MBIP,MFIP,SLATT»SLONT,TLATl?,TLONl2«XAWT,
Z YAWT.XXSIG,YY5IG,KK2,MM2,SJTAG,JTG2,ICNTMN,ICNTMX,XlP,YlP,
7. XTSIP.YTSIP.OXSIOP.DYSlGP.WIP.TSLATZtTSLONZilSNAP.HlOROr,
Z UHART,SiliMAH, tCNTT . SI CMTT , SPHD21 , SPREU2, ALPHAA , BETA A , TA , H^A , IJK A ,
Z IRTNS)
lF(IRTNS.tQ.170) GO TO 174
CALL ALTiJTAI SJTAG.JTG2. TSLAT2, TSLON2,uBART .SIGMAH, I CNTT .
Z S ICNTT t ALPHA A.BETA A, TAtR2A, IJK A, IRTNS)
iFtlRTNS.tQ,] 70) GO TO 174
CALL CLSTM(J,KK,MM,JTG?,NTPDYiITABFb,M,NOYDTA,NMM,KK2tMM2, IRTNS)
IF ( IRTNS. E.Q. 90! GO TO 95
iFdRTMS.tQ.Hb) GO TO hi
174 CALL CONVER(TSLflT2.TSLON?,UBART.ICNTT,SICNTT,SPRE2.SIGMAH,
Z SPHE02,ALPHAA.dETAAtTA,R2A, UKA,TSLAT» TSLON,
Z UBAh.ICNT.SICNT)
TTLAT2=TSLAT
TTLAT4=TSLAT
TTLON2=TSLCN
TTLON4=TSLON
ALPHA2=ALPHA
ALPHA4sALPHA
130
-------
UK2=IJK
IJKAsIJK
ICNT2=ICNT
ICHT<»=1C'NT
SICNT2=SICNT
SICM* = SlC\T
TLAT12=TSLATS
TLAT1*=TSLAT5
TLON12=TSLON5
60 TO 200
62
SPRFA1=0.
IJK1=0
SICNT1=S-JTAG(1)
GO TO 13*1
96 CONTINUE
1340 CALL DIS02(J,K.M»KK.MM,f«XNSTA,MSTA2tITABFi3, ITS,
ZITPOY«IMM,RAUI!JS,TSOUH,MKIP,MFIP,SL&TT,SLONT,TL
ZITPOY«IMM,RAUI!JS,TSOUH,MKIP,MFIP,SL&TT,SLONT,TLATH,TLONH,XA«T,
Z YAWT.XXSlQ,YYSIG,KKTiMMT«SJTAG»JTAG,ICMrMN,ICNrMX,XlPfYtP,
2 XTSIH.YjSIP.OxsIGP.DYSISP.WIPf TSLAT2. TSLON2i I SNAP, W IOROT ,
Z UBART,SIGMAH,1CNTT»SICNTT,SPRD11»SPHED2,ALPHAA,HETAA,TA.H2A,IJKA«
Z IRTNS)
IF(IRTNS.EQ.170I GO TO 1 7S
CALL ALT>JTA( SJTAG.JTAG, TSL AT2, TSLON2 , U8ART t S IGMAH, ICNTT »
z SICNTT.ALPHAA,I^:TAA,TA»«?A,IJKA,IRTNS)
CALL CLSTM(J,KK,MM,JTAG«NTPOY, ITA6FB ,M,NOYOTA,NMM,KKT»MMT , IRTNSl
IF(IRTMS.EQ.90) GO TO 96
IF(IRTNS.EQ.86) 60 TO 62
175 CALL CONVER(TSLAT2,TSLON2iUBAKT.ICNTT,SICNTT.SPREl»SIGMAH,
Z SPRE02,ALPHAA,HETAA,TA,R2A,IJKA,TSLAT, TSLON,
Z UBAK.KNT.SICNT)
510 TTLAT1=T^LAT ^
TTLON1=T5LON
UBAR1=UHAR
ALPHA1=ALPHA
BETA1=8ETA
T1*T
R21=R2
ICNT1=ICNT
SICNT1=SICNT
TLAT11=TSLAT5
SPRD11=SP«ED1
SPRE1=SPHED
SP«EA1=SHREAD
GO TO 13*1
63 TTLAT2=9*6
TLON2=9996
AR2=0.
TT
U8
ALPHA2=99.
BETA2=99.
P22=99.
ICNT2=0
SICNT2=SOTAG(1)
GO TO 13*2
97 CONTINUE
1341 CALL DIS02(J,K,M,KK»MM,MXNSTA,NSTA2,ITABFB, ITS.
ZITPDY»tMM,RAai!JS,TSUUR,M,*IP,M£lP,SLATT,SLONT,TLAT12,TLONl2,XAwT,
Z
Z WT»XX,YYS'3,KK
-------
,86) GO TO 63
CGNVtR(TSLAT2,TSLOM:>,UBART,ICNTT,SICNTT,SPRE2»SIGMAH,
SPREO«d,ALPHAA,BETAA,TA,R2A,IJKA,TSLAT. fSLON,
IF ( IRTNS,tQ.B6)
176 CALL
2
Z (JBAR, ICNTt SICNT)
512 TTLAT2=TSLAT
ALPHA3=ALPHA
BtTA2=BETA
T2 = T
R22=R2
IJK?=IJK
ICNT2=ICNT
S1CNT2»SICNT
TSLAT-5
SPRD2l=SPHEni
SPREA2=SP«EAO
GO TO 13*2
TTLON3=9996
SPREA3aO.
ALPHA3s99.
IJK3=0
ICNT3=n
SICNT3=SJTAG(1)
<30 TO 13*3
1342 CALL DTS03< J,K,M,KK»MM,MXNSTA.NSTA3»1TABFR, ITS.
ZITPDY.TMM,RAUH1S,TSDUP,'1HIP,MF:iP,SLATT,SLONT,TLAT13,TLONl3,XAWT,
2 YAWTtXXSIG,YYSlG,KK3,HM3.SJTAr,1JTG3,ICNJMN,ICNTMX,XlP,YlP,
" "
Z XlblK.YISlH.ilXbUjH^UYIiltjK.WIK.IbLAI^.IbLUN^.IcNAP.WIOKOI,
Z UBART.SlGMAH, ICiMTT, SIC'xITT ,SPR031 , SPKED2 , ALPHA A , SET A A , TA , R2A , UK A,
Z IRTNS)
IFflRTNS.EQ. 170) GO TO 177
CALL ALTUTA( SJTAG,jrG3. TSLAT2»TSION2,UBART,SIGMAH,ICNTT,
Z SICNTT,ALPHAA,BETAA,TA,R2A,IJKA,IHTNb>
IFdRTNS.t-Q.170) GO TO 177
CALL CLSTM(J,KK,MM,JTG3»NTPDY»ITABFB,M,NDYIHA,NMM,KK3,MM3,jRTNS)
IF(IRTNS.EQ.90) GO TO 9H
IF(IRTNS.EQ.H^) 50 TO 64
177 CALL CONVE.R(TSLAT2,TSLON2,UBAHT,ICNTT,SICNTT,SPRE3»SIGrtAH,
Z SPREU?,ALPHAAtdETAA,TA,R2A,IJKA,TSLAT,TSLON,
Z llBAH.ICNT, SICNT)
514 TTLAT3aTSLAT
TTLON3=TSLON
ALPHA3aALPHA
BETA3=BETA
T3 = T
R23=R2
IJK3=IJK
ICNT3=ICNT
SICNT3=SICNT
TLATl3=TSLAT5
TLON13=TSLON5
SPPD3USPHED1
GO TO 13*3
65 TTLAT«=9*6
TTLON4=9996
UBAR4=0.
SPRPA4=0
ALPHA4
=0.
=99.
T4=99.
R24=99.
ICNT4=0
GO TO 200
1343 CALL Ols03o ro 65
17S CALL CONVt:H(TSLAT2,TSLON?,UBAr(T,ICNTT,SICNTT,SPRE4.SIGMAH,
132
-------
Z SPRED/',»LPHAA,8ETAA,TA,R2A, IJKA , TSLAT . I"SLON,
2 I|HAR,ICNT,SIC GO TO 300
IF(TTL*Tl.GT.^5n.ftND.'JTl.ftT2.6T.*50.AND.TTLATi.GT.450.AND.
Z TTLAT4.GT.450) no TO ?50
IF(TTLATl.UT.3Sn.AND.TTLAT2.LT.350.AND.TTLAT3.LT.350.ANO.
Z TTLAT4.LT.350) GO TO 250
IFM.LT.750.AND.TTLON2.LT.750.AND.TTLON3.LT.750.AND.
Z TTLON4.LT.350) ^0 TO 250
IF(TTI_ONl.GT.'J5n.AND.TTLON2.nT.950.AND.TTUON3.GT.<350.AND.
Z TTLON4.GT.950) GO 10 250
GO TO
250 IAL
255 WRITEIKU^IT) TTL AT 1 , TTl ONl , UH API , SPKEA 1 ,ALPMA1 ,HETA1,T1.W21.IJK1
Z , IALL
WH1TF-(6.<^500) TTL ATI , TTL ON] ,uR»«l t SPRtAl , ALPHA ] ,HETAl,Tl«R21t
Z IJKI.KU'NIT
2500 FORMAT (lX,?HO,fiKl'S.3, TTL AT4 , TTLON4 , UHAH4 , SPRE A4 , ALPhA4 . BET A4, T4.R24,
Z IJK4,KUNIT
299 CONTINUE
300 CONTINUE
399 CONTINUE
400 CONTINUt
DO 475 1=1,5
WRITE (14) 010 ( r ) ,OLAT ( 1) ,OLON( I) ,WIUROT ( I ) ,EHS02 (I ) ,ERSUt-( I)
475 CONTINUE
WRITE (14) NO,NL)Y,MTS,"JTPDY,LTAATiLMAAT,LRAAT
WRITE(14) TSOUR,ALATT,ALONL
STOP
END
SUBROUTINE ALFPTA (SJTAG, JTAG, TSLAT^, TSLOM^, UBART, SI&MAH,
Z ICNTT, SICNTT. ALPHA A, t^ETAA, TA, H2A, IJKA,
C
C»«ALTERNATE DATA SELECTION
C
LOGICAL SJTAG(^) ,SICNTT
DATA ORS, ANA/3HOH?,3HAi\JA/
IRTNS=n
JTAG=JTAG»1
IF (JTAG .LE. 4) RETURN
JTAGT=5
TSLAT2=94.q
TSLON2=999.9
UBART=0.
SIGMAH=0,
133
-------
JA=99.
R2A=99.
I JKA = 0
ICNTT=ICNTT
SICNTT = SJTAG(JTA'3T)
IRTNS=170
RETURN
END
SUBROUTINE AWINf) (SHT, 5ATH, DI , LBAATt LTAAT. WHT, XW» Y*, NLVL,
Z XAW, YAW, xSlfiMA, YSK.MA)
C
C»«AVERAGE WIN^S
c
INTEGER SHT,SATN,BTH,TTH,TH
INTEGER rfHT(50)
DIMENSION XW<5<1) ,YW(50)
DIMENSION MHTI ! ') *> )
|_8 = SHT*LHAAT
IF (LBAAT .NE. 0 .OR. DI .GT. f<0. .OR. SHT ,GT. SATH) LB=SATH+L3AAT
LT=SATH+LTAAT
STHXW=O.
STHYW=0.
STH=0.
C«»CALCULATE M HEIGHTS FOR WIN') LEVELS
IF(NLVL.EU.I) GO TO 21
MLVL=NLVL»2-1
DO 13 M=1,MLVL,2
MDM»(M»1) /^
MHT (M) =WHT (MOM)
12 CONTINUE
C««CALCULATE M HEIGHTS FOR MI.) LEVELS
DO l
-------
00 iol M=ML*,MMLTMI
TH=MHT (M+l )-MHT(M)
MDM=(M +
YSIGMA=SUKT(YSIGMA/STH)
RETURN
500 XAW=99.
YAW=99.
XSIGMA39?.
RETURN
END
SUBROUTINE BIN(LENtS,$,$)
COMMON/INDATX ODAT1510).NXTSTA.LSTSTA
R|AD^2>2LEN,
RETURN 2
END
SUBROUTINE CLST^1(J, KKt "IM. JTAG. NTPDYt ITA8FB, M, NOYDTA, NMM,
7. KKT, MMT, IRTNS)
C
C»»CLOSEST TIME FOR 4LT?:RNATE OATA
IRTNS = (?h
IF (JTAG. ta. 4) SO TO 1
C»»S£COND CLOSEST TIMt, JTAG=?
IF(IT*RFH.E0.2) GO TO 3
C#»FORw TRAJECTORIES
2 IF(MOD(J»2) .EQ.O) GO TO 30
GO TO 10
C«* BACK TRAJECTORIES
3 IF(MOD(J«2) .EQ.3) GO TO 10
GO TO TO
C««THIRD CLOSEST TIME, JTAG=4(0»S ONLY)
8 lF(ITARFrf.EQ.2) GO TO Z
GO TO 3
10 KKT=KK+1
IF (KKT.L£..NTPrJV) 00 TO 20
KKT = 1
MMT=MM*1
GO TO 60
20 MMT=MM
GO TO f>0
30 KKT=KK-1
IF(KKT.GE.I) GO TO 40
KKT=NTPOY
GO TO fiO
40 MMT=MM
60 IF (MMT .LT. 1. .OR. MMT ,GT. NMM) HETUKN
lF(ITABFc*.EO.H) ijO TO 62
IF (M * MMT -1 ,«T. iMOYOfA) RFTURN
IRTNS=90
RETURN
62 IF (M+MMT-1 .LT. NMM) RtTURN
IRTNS=ALPHAA"JETAAfTAtR2AfIJKA'TSLATtTsLnNt
LOGICAL SICNT,SICNTT
INTEGER TSLAT,TSLON
COMMON/SOST/TSl-ATS.TSLO»'i,SPRFDl,
A SPREO, SPREAD, ALPHA, HET«,T,R2, UK
TSLAT = TSLAT2«l'i.*,5
TSLON=INT (TSLON2»10.*lfl00.5) -1300
ICNT=ICNTT
SICNT=SICNTT
135
-------
TSLAT5=T5LAT?
TSLON5=TSLON2
SPREDUSPUE02
ALPHA=ALPMAA
BETA=BETAA
T = TA
R2=R2A
IJK=IJKA
RETURN
END
FHC0025 31 PAGES PRINTED.
SUBROUTINE niSO(J, K, M» ^K, MM, MKNSTA, NSTA, ITABF8,
A ITS. ITPOY, IMM, RADIUS, TSDUR, MBIP, MEIP, SLATT, SLONT, TSLAfl
A, TSLONlt XAWT» YAWT, XXSIG, YYSIG, KKT, MMT, SJTAG, JTAG, ICNTMN,
A ICNTMX, XIP, YIP, XTSIP, YTSIP, DXSIGP, OYSIGP, WIP, TSLAT?,
A TSLON2, ISNAP, VIOHOT , IJRART, SIGMAH, ICNTT, SICNTT. SPRE01,
A SPRED2, ALPHAA, RETAA, TA, R?A, IJKA, IRTNSI
c
C»«DISPLACEMENT CALCULATION USING OBSERVED WINDS
C
LOGICAL SJTAG<5) tSICNTT
DIMENSION NSTAUTPriY.IMM)
DIMENSION XXSlG(f'XNSTA) , YYSIG (MXNSTA)
DIMENSION XAWT(^XNSTA) , YAt»T {MXNSTA )
DIMENSION SLATT(MXNSTA) ,Sl ONT(MXNSTA)
INTEGER XIP(ICNIMX) , YIP (ICNTMX) , XTSIP ( I CNTMX) , vTSIP ( ICNTMX) ,
A WIP(ICNTMX) ,XTSP,VTSP
COMMON/STAT/SLAT ( 60,2,?i), SLON (80,2,5),XAw ( 8n.2»5),
Z YAW (80f2,S) ,XSK'MA(80,?.B) , YSIGMA (80 ,2,5)
INTEGER UXSIGPdCNTMX) ,DYSIGP( ICNTMX)
DATA PI/J.14159/
IRTNS=0
SWXTSI=0.
SWYTSIaQ.
SDXSlGsO.
SDYSIGaO.
SWI=0.
ICNTT=0
NEAR=0
IF(JTAG.EQ.O) SO TO <»0
IF (NSTA(KKT,MMT) .FO.O) f-FTURN
NSTAOO=Ni>TA(KKT,MMT)
DO 80 I=l,NSTAnO
SLATT(I)=SLAT( I,KKT»MMT)
SLONT in =SLON( T,KKT,MMTI
XAWT(I)=XAW(I»KKT,MMT)
YAWT(I)»YAw(I,KKT.MMT)
XXSIG.. VAWT(I),XXSIG(1,,YYSI(5(I), TSDUR, '
I ISLAI11SrM2LONlt IT4HfrR« NtAR, XI, YI, XTSI, YTSI, DXSIGI, OYSIGI,
L W I , I JH I "NS )
IF (IJRTNS .EQ. 1^0) 00 TO 1?0
ICNTT=ICMTT*J
XlP(lCNTT)sxl+.5
YlP(ICNTn=Yl».S
.
VTSIP(ICNTT)=yTSI/3.*.5
DXSIGP(ICNTD=nxSIGI/-i. + n.S
DYSlGPf ICNTT) =OvsiGI/3.* 0.5
WIP( I CNTT)=wl»] 000000.
SWXTSI=S"XTSI*WI»XTSI
SDXSIG=SOXSIG+w!«PxSIGT
SWI=SWI*WI
150 CONTINUE
IFdCNTT .LT. ICiNlTMN .AND. NEftR ,NE. 1) RETURN
C»«NOT SllFFICItNT 0«St"RVEO WINDS
C»» TRAJECTORY SEGMENT 01 SPLICE ME NT »»•»»» «•»»»»•»»«»««».»« »««•»«»»«»««»«»«««««
IF (swi .EQ. n.) swi=.noni
XTS=SWXTSI/SvJl
136
-------
xsox=snxsiG/sv,i
XTSP=XTS/3.+.S
YTSP=YTS/J.+.5
TSLAT2=TSLATl+vTS/60.
TSLON2=T3LON1-XTSX(60,*COS(TSLATI*P 1/180.))
U8A«T=.51b«SORT(XTS«XTS+rTS«YTS>/TSDUR
SIGMAH=(Ae(S(XT5) *YSDY+AHS / (UB/\RTT*TSOUR V
CAUL *SZ < TSLAT1 « I SLON1 , TSLAT2,TSLON2,XTS,YTS,SPKED?,ALPHAA,RETAA,
A IJKAf RJ?A,TA,WIOROT)
SICNTT=SJTAG(JTAQ*1)
IF (M .Gt. MBIP .AND. M .LE. MFIP) GO TO 172
IRTNS = 17')
RETURN
172 IF(ISNAP ,EQ. 1) '"O TO 17?
C»*PRINT TRAJECTORY Sfc'PMENT n I SPL ACFMENT AND««»«#«*#»«««*«»««»«*»»»«»«««
C««INOIVIDUAL STATION KFLATIVf LOCATIONS AND D I SPLACEMENTS*«*»*»*»»*«<»»»
IF ( ISNAP.tQ.l) GO TO 175
ISNAP=1
WRITER, 173)
173 FORMATUH ,« j K M TA.5 CNT XTS YTS»/
A 1H ,« XI YI XTSIYTSI WI»,
A « XI YI XTSIYTSI WI XI YI XTSIYTSI Wl«t
A » XI VI XTSIYTSI WI»)
175 WRITE(6tldO) J«^ t'1, JTAfi, TCNTT, XTSP, YTSP,
A (XIP(IP),YIP(IP).XTSIP(IP)»YTSIP(IP),WIP(IP),IP=1,ICNTT)
180 FORMAT)!* ,3I3»*I4/
A 1H ,25x.*(lX»2It5»3-H)/1H t25X,4(3X»2l5»3I<»)/
A 1H ,?5X,4(3Xt2I5i il4)/!H , 25X.4 (3X , 21 5, 3 I f ) /
1H ,
A 1H ,25X,4(1Xi2IS, 3UI/JH , 2bX, 4 ( 3X , £15. 3 I 4)
A 1H ,2bX.*(3Xt/
IRTNS=0
SWXTSI=0.
SWYTSI=0.
SDXSIG=0.
SDYSIG=0.
SWI=0.
ICNTT=0
NEAR=0
IF( JTAG.tQ.O) GO TO "»0
IF(NSTA2(KKT,MMT) .EQ.01 RETURN
NSTADO=NSTA2(KKT,MMT)
DO 80 I=1.NSTADO
SLATT( I) =SLAT2( I, KKT. MMT)
SLONT(I) =SLON2(I , KKT, MMT )
XAWT(I) =XAW2(I,K|
-------
IF (IJ9TNS .EO. ISO) GO TO
1CNTT=ICNTT»1
XlP(ICNTT)=xI+.5
YlP(ICNT1)=Yl + .S
XTSIPUCNTT)=XTSI/3. + .S
YTSIP(ICNTT)=YT'JI/1, + .5
. .
OYSIGPf ICNTT) =r>YSir, 1/3. +0.5
wTP( iCNTT) =wi«iononoo.
SWYTSI=SWYTSI+WI*YTSI
150 CONTINUE
IF(ICNTT .LT. IC.NTMN ,ANO. NEAR .NE. 1) RETUHN
C««NOT SUFFICItNT O'lSc^ ~"
C<»«TRAJECTORY StGMEMT
IF (SWI «tQ. 0.)
XTS=SWXTSI/SWI
YTS=SWYTSI/SWl
YSDY=SDYSIG/S
175 WRITE(6tl80) J,K ,M, JTAGt ICNTT t XTSP t YTSPt
A (XIO(H') ,YIP(IP),XTSIP(IP) .YTSIP(IP) ,WIP(IP) ,1^=1, ICNTT)
180 FORMAT (1^ »3I3»"»I4/
A 1-t »i?5Xt4 (3X»2I5. jI/»)/lH ,25X,4(3X,2I5,3I4)/
A IN ,25X,<»(TX,2I5»3U)/1H t 25X, 4 ( 3X , 2 15,31 *) X
A 1* t25X.4<3x»215.3U)/lH » 25X, 4 (3X, 2 15, 314) /
A I* ,25X,4(.TX,2I5,3I4)/1H ,25X,4 (3X, 2 15 , 31*1 )
IRTNS=17U
RETURN
END
SUBROUTINE OI503(J,K,M,KK,MM,MXNSTA,NSTA3, ITAHFR,
A ITS, ITPDY, IMM, RADIUS, TSDUR, MBJP, MEIP, SLATT, SLONT, TSLATl
A, TSLON1* XArtT, Y4WT, XXSlG, YYSIG, KKT, MMT, SJTAG, JTAli, ICNTMNt
A ICNTMX, XIP, VIP, XTStf, YTSIP, DXSIGP, DYSIGP, WIP, TSLAT2,
A TSLON?, ISNAP, *IOHOT , URAUT, SIGMAH, ICNTT, SICNTT, SPRE01,
A SPREO?, ALPHAA, BETAA, TA, R?A, IJKA, IRTNSI
DIMENSION NSTA3(ITPUY,IMM)
DIMENSION XXSHHMXNSTA) , YYSIG (MXNSTA)
DIMENSION XAWT fMXNSTA) , Y AWT (MXNSTA )
DIMENSION SLATT(MXNSTA) , SLONT (MXNSTA)
INTEGER *IP( ICNTMX) ,YIP< 1CNTMX) , XTS IP ( ICNTMX) , YTSIP (
A «ap=-<'»W3( I, KKT, MMT)
138
-------
XXSIG(I)=XSIGM'»(I,KKT,MMT)
YYSrG(I)=YSIGM3(l,KKT,MMT)
80 CONTINUE
GO TO 95
90 IF (NSTA3
YlP(ICNTDsYl*.S
XTSIP(TC'NTT)=XTS!/3.*.5
YTSlP(IC^TT)=YTSI/3.*.5
.
DYSIGP(ICNTT) =DYSIGl/3.+0.5
WlP(lCNTT)=WI*1000000.
SWXTSI=S"XTSI+WI*XTSI
SWYTSI=S*YTSI+WI«YTSI
SDYSIG=SL>YSIG+WI«DYSIGI
150 CONTINUE
IFflCNTT .LT, TCNTMN ,A" J,K,M,JTAG, ICNTT , X TSP , YTSP,
A Iri ,a5X,4(3Xi2I5,3l4)/lH ,25X,4 (3X , 215,
A 1H ,25X,*(TXi2I5.3U)/lH , ^5X , 4( 3X , 21 5 , 3 U )
A 1H ,25X,4(3X«2I5,3U)/1H , 25X,<* ( 3X , ? 15, 31 4) /
A 1H ,?5X,'»(3X»2I5,3U)/1H ,25X,4(3X,2I5,3I4) )
IRTNS=170
RETURN
END
SUBROUTINE OTARKO(H, MXNSTA, NTPDY, NMM, NST A,NST A3 ,NST A3,
1>
2 OLON^'siu^siO?*SID3?TAB9t DTAflL' OTABR« L8AAT, LMAAT.LTAAT.OLAT,
C
C««OATA BLOCK »-OR OBSERVED INFORMATION
INTFGtR JlC? («>.f'<;T.\. IT"Li\ ) , Sl^t ( ."XNS T A , I TPU Y )
INTtGER 5 I (hO) iSr? (80) , r. I 3 (HO)
DIMENSION SLA (80»,SLO (00), xA (80), YA (HO),XS (Hn),YS («0)
DIMENSION SLA? (HO) , SLO? (rtO) , XA?(80) »YA?(SO) ,XS?(HO) ,YS?(8o
DIMENSION SLA3(fiO),SL.(J3{fO),XA3(HO)fYAj(80) tXS3 (HO) ,YS3(«0)
COMMON/STAT/SLAT ( 60,?, 5), SI ON < HO ,
-------
2 YAW <30«<2»5) «XSIG*,b) . XAW2 ( 80,2,5),
Z YAw2(flO«2,5) ,XSIRM2(flO,i;»S) ,YSI 6^2(80,2, 5)
COMMON/STAT3/SLAT3(80.2»'J> , SI ON3 < HO . 2 , 5) , XAW3 < 80,2,5),
Z YAW3(flO»^tS),XSIC-*'3(BOt2iSl,Y
NSTA2(K,'v'f ) =NSTAT?
DO 11 I=1,NSTAT2
SLAT2(I,K,MM)=SLAT2(I,K,MM»l)
SLON2(I,K»MM) =SLON2(I »KiMM«l)
XAW2(I»K,MM)=X'iW2( I»K,MM+l)
YAW2(I,K,MM)=YA*2(I,K,MM+1)
XSIGM2( I ,K,MM) =XSIGM2(t ,K,MM+1 )
YSIGM2(I«K,MM)=YSIGM2(I,K,MM»1)
11 CONTINUE
NSTAT3=NSTA3(K,MM+1)
NSTA3(K,MM)=NSTAn
DO 12 I=1,NSTAT3
SLAT3(ItK»MM)=SLAT3(I,K,MM+l)
.
XAW3(I,K,MM)=XAH3(I »K,MM*1)
YAW3(I,K,MM)=YAW3(I,K,MM*1)
XSIGM3(I»K,MM) =XSIPM3(I,K,MM*1)
YSIGM3( I»K,MM) =YSK,W3( I tK,MM*l)
12 CONTINUE
10 CONTINUE
60 TO 200
8 DO 199 N=1,4,MTPI)Y
IF(N.EQ.l) K=\
IF(N.E0.3) K=?
IFLAG=1
CALL PnsfPtM.MM.K;)
CALL ^R^AVO ( rtXNSTAt DTABT, DTAB8, DTABL, DTARR, LBAAj, LTAAT,
Z OLAT.OLON.SI ,SLA ,SLO tXA ,YA ,XS ,YS ,NSTA (K,MM))
DO 101 I^*»80
SID (I,K)=SI(I)
SLAT (IiK,MM)'SLA(I)
SLON (I,K,MM)=SLO(I)
XAW (I,K,MM)=<4(I)
YAW (I,K,MM)=YA(I)
XSIGMA(I,K,MM) =X5(I)
YSIGMA( ^K.MM) =YS(I)
101 CONTINUE
IFLAG=2
REWIND 2
CALL POSTP (M,MM,K)
CALL RUAVO(MX'V9TA,OTARr.DTABR,DTAHL,OTABP,LBAAT,LMAAT,
2 OLAT,OLON,SI?,SLA2.SL02,XA2,YA2,XS
-------
YAW3(I,K«MM)=YA3(I )
XSIGM3(I«K,MM)=XS-< (I)
YSIGM3(I.K,MfM = YS3 < I)
103 CONTINUE
REWIND 2
1^9 CONTINUE
200 CONTINUE
RETljHN
END
SUBROUTINE INPUT (NTPDY, TSOUR, RADIUS* ICNTMN. ICNTMX, WHIP, MElP,
Z MXNO. NO. 010, 01 AT, OLON, WIOKOT. IHRINP, NLVLA, NTSPOY,
I NTS. NMM.ERSOP.EPSUL)
C«*INPUT INFORMATION
DIMENSION ERSO?(5) .t'RSULC?)
DIMENSION OID ,EHSUL(!)
1 FORMAT (A-Mx.K4.1,lX,F>'i.l,lX.F10.1,lX»F10.0.1X,F10.0)
101 CONTINUE
C««CARO 2«»
C«»DATE COMPUTATIONS HEGIN: OAYdBDY). MONTH(IRHO), YEAR(IBYR)
C»«NUMBER OF DAYS COMPUTATIONS DESIRFDINDY)
RE AD (5, 3) I80Y.IMO, IBYC? . I HOUR , I HMO
2 FORMAT(4(I2.1X) .£X,A3)
C»»CARD 3«*
C«»DIRECTION IN TiiEiTniq) , FORwARD(FOrtv«) . BACKWARO(RACK)
C»«DUHATION IN OAYSJNIIYDIJR)
PEAO{5,3) TDP?,NOYOUR
3 FORMAT(A^»1X,I?)
C««CARD 4«»
C««NUM8ER OF DAYS OF ml NO INPUT DATA (NOYOTA )
4
C»*CARD _
C*«TYPF OF WIM> INPUT DATA(V»TYPE)
C«*08SERVED(OaS), ANALYZED (ANA) ,OBSF«VED AND ANALYZED(0,A)
READ<5,5> WTYfJF
5 FORMAT(*3)
C**CARD 6«»
C»»BOUND4RIES ''OR OBSERVED WIND INPUT DATA
C**TOP LATlTUDt(DTA8T!, BOTTOM LATITUDE(OT48B),
C»»LEFT LONGITUDE(OTAHL), HIGnr LONGITUDE(OTABR)
READ(5,6) QTABT«OTABB,QTABL»OTftHR
f> FORMAT(F4.1,lX,M.l.lX,fL^.l.lX.F6,l)
C»»CARD fl«»
C«*TRANSPaRT L*YER BftSE(lBAAT), AND TOP(LTAAT) IN METERS
C«»ABOVE AVERAGE TERRAIN
READJ5.8) LBAAT.LTAAT.LMAAT
C»*CARO 0»«
•C»*BOUNDAPIES K)R MAPS IN SIJBMOUTI NFS:
C«» LATITUDE 0*- MAP TOP BOUNDARY ( ALATT)
C»« LATITUDE OF MAP BOTTOM BOUNDARY(ALATB)
C»«LONGlTUDE OF MAP LFFT BOUNDARY(ALONL)
C«»(MAP WIDTH IN LONGITUDINAL DEGREES IS OBTAINED
C SEPARATE TABULATION)
READ(5,9) ALATT..M ATB,ALONL
C*«END INPUT DATA
C
C
Co» INPUT PARAMETERS4****"**4*4***'**"**0 ****<**<>*** ************************
C««NUMBER OF TWAJECT^KIFS PEP DAY
NTPDY=2
C»«TRAJECTORY SEGMENT DURATION IN HOURS
C»»RAOIUS IN NAUTICAL MILES FOK INCIUOING OBSERVED WINDS
C*«IN DISPLACEMENT CALCULATIONS
C««MINIMUM NUMWER OF STATIONS *ITHIH HADIUS
C»«FOH A DISPLACEMENT CALCULATION
iCNTMNr?
141
-------
C*»MAXIMUM NUMBER OF STATIONS '*'ITHIM RADIUS
C*»FOR A DISPLACEMENT CALCULATION
ICNTMX=JO
C*«LATITUD£ ABOVE WHICH DISTANCES ARE CALCULATED
c IN A POLAR COORDINATE SYSTEM
POLl AT=60.
c«»HEGiNNiNG AND ENDING -1 IN M-I,NOY LOOP FOR PRINTING
C«»TRAJECTORY SEGMENT DISPLACEMENTS AND INDIVIDUAL STATION
C»»RELATIVt LOCATIONS AMD DISPLACEMENTS
MBIPsO
MEIP=0
C«»PRINTING INTERVAL IN HOURS
C»«NUMBER OF ANALYSIS LEVELS
NLVLA=4
C««F.ND INPUT PARAMETERS*****
C********************************* **************************************
C
C»«NUMBER OF TRAJECTORY SEGMENTS PER DAY
NTSPDY=2*/TSDUW
C»»NUMBEH OF TRAJECTORY SEGMENTS
NTS=8
C«*NUMBER OF DAYS IN DATA SLOCK (DTABKJ
NMM=NDYDUR»1
C«»PRINT INPUT OATA CARD iwo^MATiorf
WRITE (6.20) (OTU(IO) .OLAT(IO) .OLON(IO) ,10=1, NO)
WRITE (f. ,21) IBOY
3 TDIH.NDYDUR,
4 NDYOTA,
5 WTYPE,
8 L«AAT,LTAAT
30 FORMAT(lrtl,16X,*INPUT FOR TRAJECTORY COMPUTATIONS*/
AIM ,*INPUT*/1H ,« DATA»/1H ,* CARD*/
A A3, «{«,F4. 1, IX, Ff>. 1, ») */10(46X, A3, *(»,F4.1,lx,F6. It «)«/))
31 FORMAT (
21H ,4X,4lHa DATE THAT COMPUTATIONS BEGIN «««*»»»»* , 1 2, IX, A3, 1 X, I
B?H ,4X,4lH NUMBER OF DAYS COMPUTATIONS DESIRED »« ,I2/
31H ,4X,41M3 DIRECTION IN TIME »»»««»*«»»*»»«»«*»»» ,A4/
A1H ,4X,4lH DURATION IN OAyS *«**«**»*»«*«»«*«»»** ,\Z/
41H ,4x,4iH4 NUMBER OF DAYS OF WIND INPUT DATA «««» ,ia/
51H ,4X,4lH5 T^PE OF i^Ii-lD INPUT DATA »»«»»»*»»«»»»» ,A3/
6JH ,4X,4IH6 BOUNDARIES FOR OBSERVED WIND INPUT DATA/
AlH ,4X,4lH TOP AND BOTTOM LATITUDES ««««»««»«»««» ,2X,F*.l,
A 3x,f-4.1/
81H ,4X,4lH LEFT AND RIGHT LONGITUDES ««»»»««««*«» ,F6. 1 » 1 X ,F6. 1 /
B1H ,4X,30H8 TRANSPORT LAYER BASE AND TOP/
AlH ,4X,4lH IN METERS ABOVE AVERAGE TERRAIN »*»»«» ,14, IX, 14)
WRITE<6,<25) ALATT,ALAT^,ALONL
25 FORMAT (
91H ,4X,3&Hq HOUNOARIES FOR MAPS IN SUBROUTINES/
AlH ,4X,41H TOP LATITUDE «««*«»***•»*»*»•»»»««««*«* ,?X,F4.1/
B1H ,4X,4lri BOTTOM LATITUDE »*»*»««»*»««»»»*«*»*»» ,?X,F*.l/
C1H ,4X,4lrt LEFT LONGITUDE »*»*»»««««««««*»««««««« ,F6.1>
WRITE(6,JO)
30 FORMAT(lHl)
RETURN
END
SUBROUTINE ITSI«< ( JTAO, RADIUS, SLAT, SLON, SLATT, SLONT, XAW, YAw,
Z XSIGMA, YSIGMA, XAWT, YAttT, 9.) RETURN
GO TO 12
10 IF (XAW .E<3. 9<).) RETURN
12 IFfO.LE. RADIUS/?.) NEAR=1
IF ( JTAP.tQ.O) GO TO 20
XAWSH=XA«T
YAWSH=YA»(T
XXSIGS=XXSir,
YYSIGS=YYSIG
GO TO 18
20
142
-------
XXSIGSsXSIGMA
YYSIGS=Y5IGMA
18 GO TO (21 .22) , ITAHFB
21 XTS = 1 ,94»XAwSH«TSr>U»«
YTS*l.UK
DXSIG*! . V
GO TO ?3
22 XTS = -1 ,^'
YTS=-1 . 9<*
DXSIG=-1.<»4»XXS1GS»TSIHIR
DYSIR=-1 • 'M«YYSlGS»Tsnijh
23 CONTINUE
OTS=SQPT (XTS«YTb*YTS"YTS)
YW= Y-YTS/2.
IF([>WSO.E«.n.)
DISTW=1./UWSQ
ALINW=1 .
IF (0 .NE. 0. .AND. GTS ,NE. 0.) AL I NW= 1 . - . <5*A8s ( ( YTS«X-XTs«Y) /
I (DTS*D+ ,000\))
,
W=DISTW»ALINW
IJHTNS=0
RETURN
END
SUBROUTINE KKMM(J, K, *l, KK, MM, NMM, MOYDTA, NTSPOY , »T ABFB , IRTNs)
C»»OETERMINE PHQPER TRAJECTORY SEGMFMT
C««TIME OF DAY(KK) ANO OAY(MM) IN DATA BLOCK DTABK
IRTNS=0
KKS=KK
IF (ITA6FB .EG. 2) GO TO 52
KKT = 1. + .S">FI.OAT{J) *K-1
KK=MOD(KKT-1«4) +1
It (KK.I.T.KKS) MM=MM+1
IF(M*MM-l.LE.Nr)YQTA) IRTNS = QO
RETURN
52 KKT=5.5-.5#FLOAT(MOO( J-l .NTSPDY) +1 ) +K-1
IF(KK.GT.KKS) MM=MM-1
IF(M*MM-1.GE.NMM)
RETURN
END
SUBROUTINE POSTP(M,MM,K)
C««POSITION TAPES FOR STARTING RUN AT DESIGNATED BEGINNING
COMMON/INPUT 1/IPDY, I BMOtI8YR,NDYtTOIRf
A POLLAT,IHOUR»lMO,
A NO YOUR t NOYDTA,WTYPF.|OTAHT»OTABB,DTABL,OTABR»
A LBAAT»LMAAT,LTAAT,ALATT,ALATH,ALONL
COMMON /UATERO/IOMO,IDYM,IODY«IOHR,NRSTA,NREC,IM
COMMON /I NO AT/00 AT (510) • IENO.NXTST A « LSTSTA
IEND=0
5 CALL ROOATE
IF (IENO.Nt.0) STOP 333
IF (M.NF.l.OR.MM.NF. i .OR.K.NE. i ) GO TO 100
IF (IM.NE.IMO) (JO TO 5
IF ( IUDY.ME. IHDY) GO TO 5
IF(IDHR.NE.IHOIIH) GO TO 5
100 RETURN
END
SUBROUTINE RUAVODAT(SlO)»!£NO,NXTSTA,LSTSTA
COMMON X*INDRO/H|)S (3,bn)
DATA Pt/J.141?;^/
WRITE (6. 15) IDYHiIDMOt IDDY, IOHR,NRSTA,NHEC
IS FORMATflH ,9X,I*»lXiA3,lX,l2,«-»,I2,«Z»t6x,U.f«X,I5)
NSTA=0
1SOVER=MXNSTA
143
-------
IF(NRSTA.EQ.O) 00 TO
DO 100 I« = l ,l>" ~
56 CALL RDSTHD
SIDTslSTNO
SUATT=AUAT
IF < SUATT.UE.DTA6T.AND.SU ATT. GF.DTAHB. AND.
A SLONT.UE.OT ABU.AND.SLONT.GE.OTAUR) GO TO 40
CALU RDWlNO
GO TO 100
40 NSTA=NSTA+1
IFINSTA.LE.MXNSTA) GO TO 44
1SOVER*IS>OVER*1
WRITE (6, 42) ISOVER,MXNSTA»SIDT»SUATT,SUONT
42 FORMATflH »70X.«lSTAOVER=»,l3,^X.»MXNSTA=*tI3,
A 3X, 15, M«,F4.1, IX, Fft. !,»)»)
NSTA=MXNSTA
44 CONTINUE
I=NSTA
SID(I)=SIDT
SUON< i) =SUONT
CALL RDWIND
DO 50 UVU=1,NLVL
WHT(LVU)=HOS<1,UVL)
WSPD=HDS(3,LVU)
XW(UVU)=-WSPD»SIM(WOIR»PI/1BO.)
YW(LVU) =-WSPD»COS('«/DIR«PI/iHO.)
50 CONTINUE
Xl=(OUON-SLON(T) > *fiO.«COS < SLAT f I ) «P 1/180 .)
Yl=,200a> IOST
2002 FORMATC FR RDDATE, i-oi»02»
IFUFUAG.EQ.2) GO TO 120
GO TO 130
H4 CALL BOUT(LEN,$115,$10l)
5 ISM»0
IF(LEN.GT.20) GO TO 100
IDMO»ODAT<1)
IDYR»OOAT(2)
IDDY»Ot)AT< J)
IDHR«ODAT<4)
NRSTA»OOAT<5)
NREC»ODAT<6)
IM»ODAT<7)
NXTSTA«171
LSTSTA«171
GO TO 99
98 IEND»10
GO TO 99
101 STOP 404
99 RETURN
END
144
-------
SUBROUTINE RDOATAM)
COMMON XINOAT/00AT(510>,IEND,NXTSTA,LSTSTA
COMMON/FLAG/IFLAG
100 iniFLAG.EQ.2) GO TO 120
130 CALL TAPEX (IWINOOATA
,!_!_ ,Kr^n ».w4.,u«^,- •,1,1,IOST,S10,ODAT,LEN>
IFUOST.EQ.O.OR.IOST.EQ.4) GO TO 114
IFdOST.EQ.l) GO TO 116
GO TO 117
120 CALL BIN (LEN,$115,$116,$117)
116 IF(ISW.NE.O) GO TO §8
ISW*10
GO TO 100
117 WRITE(6,2002) 10ST
2002 FORMATt' FR HDOATA
115
ORMATM fH HDOATA, I-0«,02)
F(IFLAG.ta.2) GO TO 120
0 TO 130
CALL BOUT(LEN,$115,5101)
ir (LEN.EQ.5lO) GO TO 150
WRITE(6,2003) LEN
2003 FORMAT!* ^R RODATA, BAD LEN*I10)
150 I»I-170
NXTSTA=NXTSTA-170
GO TO 99
98 IENDMO
GO TO 99
101 STOP M
99 RETURN
END
SUBROUTINE HDSTHO
C FOR CSU TAPE5
COMMON /5THDHO/I3TMO, ALM »ALOM,1STH,ISATH,NLVL
COMMON /INOAT/O'JAT (blO).lEND,NXTSTA,LSTSTA
COMMON/FLAQ/IFLflG
DIMENSION DATO»171)
EQUIVALENCE (OPAT(1),DA1(1))
DATA IOU«/0/
IF(NXTSTA.GT.l'O) CALL HOOATA(IDUM)
ISTNO=OAT(1,N»TSTA)
ALAT=DAT(^,MXTC.T«)
ALON=OAT(3,NXT^TAJ
I=NXTSTA»1
IF(I.GT.l'O) CALL SODATA(I)
ICNT=ICNT+l
ISTH=OAT(1,1)
ISATH=OAT(2,1)
NLVL=IMT(3,I)
LSTSTA=NXT«;TA
NXTSTA=N*TSTA* (Ni.wL*2j
RETURN
END
SUBROUTINE ROrtlNO
C FOR CSU TAPES
COMMON /STHDHO/ISTNO,ALAT,ALOM,IST^.ISATH.NLVL
COMMON /INOAT'OOAT(510) , lEND,NXTSTA,LSrSTA
COMMON /HINDRD/HOS(3,50)
DIMENSION l)ATn»170)
EQUIVALENCE(OQAT(1),DAT(1))
IST=LSTSTA*?
IF( IST.GT.iro> CALL RODATA(IST)
IEN=IST*(NLVL-1)
150
10
20
99
JEN=I£N
K = 0
IF(IEN.LT.171) GO TO ISO
JEN=1 fn
DO 10 I = 1ST. JEN
K=K*1
HOS ( 1 ,K) =OAT (I.I)
HDS(2,K) =OAT(^, t )
HDS ( .3»K) =UAT (1,1)
CONTINUE
IF{ JEN.EQ.IEN) GO TO 49
CALL RPOATAlUN)
DO 20 I = U IhN
K = K*1
HOS(1,K) =OAT( 1,1)
HDS(2»K) =U4T (S, I)
HOS(3,K)=UftT (J.I)
CONTINUE
END
145
-------
§UBROUTINE RDWIND
SU TAPES
COMMON /STHORD/ISTNO»ALAT,ALON,IST8,ISATH,NLVL
COMMON /INOAT/OnAT(510),IEND,NXTSTA,LSTSTA
COMMON /WlNDRO/HDS(3.50)
DIMENSION UATI3.170)
EQUIVALENCE(ODATd) ,DAT(1) )
IST=LSTSTA*2
IF(IST.GT.170) CALL RDDATA(IST)
IEN»IST*(NLVL-1)
JEN=IEN
IFIIEN.LT.171) GO TO 150
JEN=170
150 DO 10 I = 1ST. JEN
K»K*1
HDSd.K)=OAT(l,I)
HDS(2.K)=UAT(2,I)
HDS(3.K)=UAT(3,I)
10 CONTINUE
IF(JEN.EQ.IEN) GO TO 99
CALL RODATA(IEN)
DO 20 I=1.IEN
K«K+1
HOS(ltK)«UAT(l,I)
HDS(2.K)»OAT(2,I)
HDS(3.K)»L)AT(3,I)
20 CONTINUE
99 RETURN
END
SUgROyTINE XSZ(TSLAT1.TSLONUTSLAT2,TSLON2,XTS,YTS,SPRED2,ALPHAA
A .BETAA.IJKA.f?A,TA»WIOHUT)
c DETERMINES THE DISTANCE.^,AND THE ANGLES,ALPHAA.BETWEEN TRAJECTORY
PI=3. 14159
WORO=XTOKOT/1HS2.')
R2=XTS»XTS+YTS»YTS
R2A=SQRT(R2)
T2=R2*S2
TA=SQHT(T2)
IF(R2.EQ.O.) RO TO 101
IF(XTS.EU.O.) GO TO 1
Zl=ABS(YTS/XTS)
ALPHAA=ATAN(71 )
GO TO 5
1 ALPHAA=.5«PI
5 IF (HaA.E'J.O. ) GO TO <5>
72=SORT(S2/R2)
BETAA»AT*N(72)
GO TO 7
6 BETAA=0.
7 CONTINUE
IF(TSLAT^.GE.TSLAT1.ANO.TSLON?.LE.TSLON1) GO TO 10
lF(TSLAT«;.GE.T?;LATt.flNO.TSI.Oh?.GT.TSLONl) GO TO 20
IKITSLAT^.LT.TSLATI . ANO.TSLONP.GT.TSLON1 ) GO TO 30
IF(TSLAT«J.LT.TSLAT1.AND.TSLON2.LE.TSLON1) GO TO *0
10 IJKA=1
ALPHAA=ALPHAA
RETURN
20 IJKA=2
ALPHAAsPI-ALPHArt
RETURN
30 IJKA=3
ALPHAAaPl+ALPHAA
RETURN
40 IJKA-54
ALPHAA=2«»PI-ALPHAA
RKTURN
101 ALPHAAsO.
BETAA=0.
IJKAsO
RETURN
END
146
-------
Program o4- Concentration Calculation (Model B)
PROGRAM CAl CON
1 (OATA2,l)ATA1,(U,TH'llT,TAPtc> = PATA?, TAPE *> = OUTPUT ,
A TAPE7=DATA3,
2 TAPE?l = lnr)2fi.TAPF.22=1002B,TAPE23=1002B,TAPfc.24 = in02R,
3TAPf:25=1002fl,TAPE?ft=1002H,TAPF?7 = lnO/;B,TAPE2H=i 00?fl , T APT. 29 = 1 002B .
5 TAPE35=1002B, TAPE 16 = 1002 rt.TAPF 17=100 2B, TAPE 3 H = inO?H,TAPe.J9=10Q 28,
6T*PE40= 10028, T'iPEl 4 = 1 002^, T APESO , T APE51 = 1 002B , TAPF.52 )
C THIS PROGRAM wn_l_ READ IN DATA COMPUTEU FROM THE PPEVEOUS PHQGRAM
C TO COMPUTE THE RESULTANT CONCENTRATIONS ETC.
COMMON /UUTPUT1/ TSLAT(8,2,5) ,TSLON(8,2,5> ,UBAR<8,2,5) ,
Z SPREAD (a»2»e>>
COMMON /OUTPUT?/ ALPHA<8,2,5) ,BETA < 8 , 2, S> ,IJK(R,2, ) /2.
X(I)=tALONL-LnN(I) )»COS(ABtON»PI/180.)»60.
Y(I) = (ALATT-LAT(I))»f,0.
N=f4* 1
495 CONTINUE
DO 500 IP=1,NO
WHITE(f>,M6) OID(TO), OLAT(IO) .OLON(IO)
416 F ORM A T ( 1 H 1 , A3 . » C " , F* . 1 , ] X » F 6 . 1 • " > " )
IF(IO.EQ.I) GO TO 1
IF(IO,FQ.2) GO TO 7
IF(IO.FQ,3) GO TO 3
IF(IO,EQ.*.) GO TO 4
IF(IO.E0.5) GO TO 5
1 111=21
KKK=24
GO TO 6
2 III=2S
KKK=28
GO TO 6
3 111=29
KKK=32
GO TO 6
4 111=33
KKK=36
GO TO 6
5 111=37
KKK=40
6 DO 100 I=III.KKK
KUNIT=I
DO 450 M=1,NOY
DO 451 K = 1,NTP'1Y
ICOUNT(M,K)=0
451 CONTINUE
450 CONTINUE
00 425 M=1,NDY
DO 422 K = 1, fJTPDY
00 420 vl = 1. NTS
HEAD (KUNIT) TSLAT(J,K,M) , TSLON ( J,K ,M) , UBAR(J.K.M),
Z SPREAOf J,K,M) , ALPHA( J,K,M) , BETA(0,K,M ), T(J.K,M), R2(J,K,M),
Z IJK(J,K.M) ,IALL
IFdALt. £0.999) GO TO 422
ICOUNT(M,K)=ICOUNT(M,K) *1
420 CONTINUE
422 CONTINUE
425 CONTINUE
CALL CONCAL t i, OL AT, OLON, NTS, NTPDY.NOY.TSOUR.AL ATT, ALONL. LTAAT.
Z LMAAT,LBAAT,EPS02.fc.rtSUL,WIOROT,IO,ICOUNT,N)
147
-------
100 CONTINUE
500 CONTINUE
STOP
END
SUBROUTINE y*OGET
Z (REMV1»REMV2,DPEMV,AVGl,'
2 A-i4MOAl,4HAMOA?,ACON\/RiAH,AUBAR,
Z AR£,AWOR,I,J,K,M,CS02»CS04,OS02,
Z • 1SULF»PSOa,PSULF,TSULF>
C»tnn>«»»THIS SUBPROGRAM CALCULATES THE DEPOSITION AMOUNTS OF
C 502 AND S04 DlJF TO DRY AfjD WET DEPOSITIONS ALONG
C TRAJECTORY AND ACCUMULATES THE VALUES FOR f»f. AREA
C OVER THE E.'ITIKE PERIOD
COFF = 4.3?*1.0£.-5
C»»«*»««»»««
100
goo
300
400
110
111
TF{M.
IF(M
IF(M
IF(M
IF(M
IF(M
IF(M
IF(M
IF(H
•^ — »
Eg.
.tu
.EQ
. tu
.tu
.to
.EQ
IEQ
.EQ
l.ANO.
.l.ANO
.2.ANU
.2. AND
.2.ANLI
.3. AND
.3.ANO
.3. AND
.4. AND
?K
.K
.K
.K
.K
.K
.K
.K
EQ.l
.Eu.
.E'l.
.EQ.
.EQ.
.EO.
42-
.EQ.
,EQ.
^
1
?
?
i
i
2
GO
.AND
.AND
.AND
.AND
.AND
.AND
.AND
.AND
>
TO
.J
.J
.J
.J
.J
:d
.j
1000
.LE.7)
.LE.5)
.LE.3)
.F.Q.8)
.EU.l)
.GE.6)
.LE.4)
.GE.2)
GO
GO
GO
GO
GO
GO
GO
GO
GO
TO
TO
TO
TO
TO
TO
TO
TO
TO
1000
1000
1000
1000
1000
1000
1000
1000
1000
SULn=CS02
SUFT=CS04
IJsMOU; (I-?]) .4)*1
GO TO (IPO,?00,300,400),IJ
GO TO (110,130),K
GO TO (210,230),K
GO TO (310,330)»K
GO TO (410,430)»K
IF(J.LE,2) GO TO 111
IF(J.LE.4) GO TO 112
IF(J.LE.6) GO TO 113
IF(J.LE.S) GO TO
CSUL'J=SULD
CSUFT = SLIFT
GO TO 2000
GO TO HI
GO TO 112
IFU.UE.2) c,o TO 112
IF(J.(-E.4) (50 TO 111
1F(J.LE.6) GO TO 112
IFU.LE.8) GO TO 111
~ '- -- 21?
it- (J.Lt.O) (jU ID 211
IF(J.LE.B) GO TO 112
CSULO=SULO
CSUFT=SUFT
GO TO 2000
CSULD»O.O
CSUFT=0.0
GO TO 2000
1FU.LE.2) GO TO 212
1F(J.LE.4) QO TO 211
IF(J.LE.6> GO TO 112
GO TO 111
113
11*
130
210
211
212
230
310
330
Ho
1000
2000
C»»«»»»«» CALCULATION JF DRY DEPOSITION AMOUNTS»«««««***«
GO TO
IF(J.LE.6)
GO TO 212
GO TO 211
GO TO 212
IF(J.LE.b) CO TO 212
GO TO 211
CSULDsO.O
CSUFTsO.O
IF(REMV
-------
PSULFSUSULF1»ABAMDA2
TSULF=(DSO?*PSO?)/2.0+ ( DSUI.F*PSllLF )/3. 0
KKK=1
GO TO J500
3000 HS02=0.0
(>SULF=O.O
PS02=0.0
PSULF=0.0
TSULF=0.0
3500 CONTINUE
RETURN
END
FHCQ025 2b PAGCS
SUBROUTINE CONCAL (I.OUAT.OLON.NTS. NTPOY.NDY.TSDUR.ALATT.ALONL,
Z LTAAT,
Z LMAAT»L.l*AAT.KHS03.tMSUL»WIORO,IOi ICOUNT «N>
o— — CALCULATE. THE CONCENTRATIONS ON GRID POINTS-- ------- - —
COMMON /OUTPUT!/ TSLAT(H,2,5) »TSLON(H,2.!3) ,U6AR(8,2»5) ,
Z SP«t:AO(8,2,5)
COMMON /OUTPUT?/ ALPHAIH.S.SI ,RETA(8,2»5) ,UK<8,2,5) ,
Z T(8.2»^ ) .R2(8,2.5>
DIMENSION PXXJ(11.21),PSXXJ (11.21)
DIMtNSlON PSI(13,100)
DIMENSION P«;A ( 1.1, ion)
COMMON/INPUT!/ PS(IOO)
COMMON/t>RID/*LAT
DIMENSION ORCONC5) .SORCO^(S) ,WIORO(5)
DIMENSION OLAT f'5) .OLON(S)
DIMENSION ERSO?(^>) ,EKSUL<5)
DIMENSION ICOLIr-JT (S,2)
DO ,
PSI (i.KK) =O.U
445 CONTINUE
HEWING 7
DO 555 11=2,7
DO 55b KK=1.8l
READ( '.5571 PSI ( II, KK)
557 FORMAT(F4.3)
556 CONTINUE
555 CONTINUE
DO 558 II=8,U
DO 559 KK=l.ril
PSI ( H.KK) =0.0
559 CONTINUE
558 CONTINUE
REWIND 7
DO 5000 *=1, NDV
DO 4000 K = l,NTPl)Y
ICC=ICOUNT(M,K)
DO 4444 J=l,lcr,2
DO 4S55 KK=1.81
IF(M.tO.l.AND.K.EQ.l) r,0 TO llll
If {M. tQ. LAND. K, KQ, 2) GO TO 1112
IF(M.tQ.2.AND.K.E0.1) GO TO 1113
lF(M.tQ.2.AND.K.EO.?) GO TO 1114
IF (M.tQ.3.«ND.K.EQ.l) GO TO 111S
IF (M.tQ.3.AND.K.EQ.2) GO TO 111&
lF(M.tQ.4.ANU.K.EQ.]) GO TO 1117
IFtM.EQ.4.AND.K.etl.2) GO TO 1118
mi jj=(j+n/2
GO TO 1121
1112 JJ=(J*D/2+l
GO TO 1121
1113 JJ=(J*l>/2+2
GO TO 1121
1114 JJ=(J+l)/2*3
GO TO 1121
1115 JJ=(J+l>/2*4
GO TO 1121
1116 JJ=(J*l>/2»5
GO TO 1121
1117 JJ=(J*l>/2*6
GO TO 1121
1118 JJ=(J*l)/2+7
1121 PSA(J,KK)=PS1 ( JJ,KK)
PSA(J*1«KK)=PSI ( JJ.KK)
4555 CONTlNUt
4444 CONTINUE.
149
-------
DO 3000 J=1,ICC
IF(TSLAT(J,K,n.FU.990,
AY=(ALATT-AY) «60.
60 TO 299H
2999 AX30.5«(ULON( 10) * TSLON ( J , K , M) /in.)
AY = 0.5«(OLAT(I1)*TSLAT »COS( A6LON»PI/iaO.) »60.
AY=( ALATT-AY) *f.O.
2998 CONTINUE
DO 399I3 KK=1.R1
PS(KK) =PSA ( J.KK)
-JQQR rnNTTNUE
J ^ "CALL REMOVE (I,M,K,J, vGA,vGfl,T«AR,HGT,RAMAA,RftMaa,
Z LTAAT,LMAAT,L*AAT,AX,AY,AtATT»ALONL,M
VG1 ( J»K»M) =V/GA
VG2( JtK.M) =VGH
COHVRI J.K,M)=tRAR
H(J,K,M)=MGT
RAMDAl (J«K,M) soAMAA
RAMDA2(J«K,M) =HAMBR
3000 CONTINUE
00 399-J J=1,ICC
DO 399H KK=1 ,H1
PSA(J|KK)=0.0
PS(KK)=0.0
3998 CONTINUE
3999 CONTINUE
4000 CONTINUE
5000 CONTINUE
WORO=WIO*0(IO)
DO 87 M = 1, NOY
DO 86 K = 1, NTPDY
ICC=ICOUNT(MtK)
DO 8S J = 1, ICC
WORY
DO 700 K=ltNTPOY
ICC=ICOUNT(M,K)
DO 600 J=l tICC
IF(TSLAT(J,K,M) . F.U. 999. OR. TSLON ( J,K .M) .EQ.999Q)
Z GO TO 6f6
REMV1=VG1 (J, K.M) *BAMDA1 (JiK.M) +CONVH (J,K,M)«H(J»K,M)
( J,K,M) *RAMDA2( J»K,M)
"
VG1H=VG1
VG2H = Vfi2
-------
SCONC(J,K,M)s (SORCON(IO)«WOH(J,K?M)/SWIDTH)«()R«OP2
A *( 1 ,5«CONVR<>I»K,M)»H(J,K,M)/HrtEMV )»CONC(J,K,M)*OC
GO TO 577
601 WOR(J,K,M) =WORO + ,f.»60.*l()b2.«r,PRt:An(J-l,K,M)
CONC (J,K»M) = ( (CONC (J-UK,M)«WOP( J,K«M) ) /S WIDTH) *[)A«OCON«Opl
SCONC(J,K,M)r((SCONC(J-1,K,M)»WOH(J,K,M))/SWIDTH)»Drt»DP2
A * < l.b»CONVR(,.UK,M>*H( J,K?M)/OREMV )
A *CnNC(J,K,M)"DC
577 WRITE(ft,1001) I»TSLAT(J,K,M),TSLON(J, ,M),CONC(J,K,M),SCONC(J ,K ,M)
A f M , K i J
IF (J.E'}.1) GO TO 3333
CS02sCONC(J-l,K.M)
CS04=SCONC(J-1,K,M)
GO TO 5555
3333 CS02=OHCON(IO)
CS04*SOHCON(IO)
5555 CALL «AfJ6ET
Z (REMVl,RFMV2inRfMv,VGl (J.K.M) ,VG?(J,K,M) ,
Z RAMDAl(J,KtM),RAMOA2(J,KtM),CONVR(J»K.M)»
Z H ( J»K,M) ,UBAR(.i,K,M) tR2 (J»K,M) ,WOR(J,K,M)
Z I
.
APSULF=APSULF*PSULF
ATSUL=ATSUL+TSULF
1001 FORMAT(10X,3no«?F20.5t3l7)
600 CONTINUE
666 CONTINUE
oo ana 11=1.11
=1,?J
DO 221 JJ
INDEx
-------
_^ = AL^HA(J,-<,M).HETA(J,K,M)
EE=T(J,K,M)«T!J.K,EQ.l. \NO.J.GE.2) GO TO 150
IF(M.tQ.4.AMD.K.EQ.2) GO TO 150
IF(J.EQ.1.AN0.9LAT(II).tfJ.OLAT { 10).AMU.BLON(JJ).EQ.OLON(IO))
A GO TO 350
1F(J.E,LE,WORO/(2.«lb52.) .AND.GAMMA(11,Jj).EQ.
A DELTA1.OR.GAMMMII,JJ) .EQ.DELTAS) GO TO 350
IF(C(II, JJ) .GT. F) GO TO 150
IF (IJK(J,K,M) .EQ.l) GO TO HO
IF (IJK(J,K,M).FU.2) GO TO 120
IF(IJK(J«K,M).K«.3> GO TO 120
lF(IJK(JtK,M» .F.Q.4) 60 TO 130
IF(IJK(J,K,M) .F.lJ.O) GO TO 140
110 IF(BETA(J,K,M),LT.ALPHA(J,K,M)) GO TO HI
IF(GAMMA!II.JJ).G£.«LPHMJ.K»H).ANO.GAMMAjIT.JJ).LE.PHAI1.OR.
A GAMMA(I I,JJ) .OF.(2.«f'I*PHAI?).AND.GAMMA(11,JJ).LE,2.»PI.OR.
A GAMMAdl.JJl .OT.O. . AND.GAMMA (I I , JJ) ,LT. ALPHA (J, K ,M) )
A GO TO IIP
GO TO 113
112 D1=R2(J,K,M)/COS(Al.PHA(J,K,M)-GAMMA(II,JJ))
D1=ABS(01)
IF(C(II,JJ).LE.U1) GO TO 165
GO TO 149
111 IF(GAMMA(II, JJ| ,GE.PHA]2.AND,r,AMMA(II,JJ) .LE.PHAI1) GO TO 112
GO TO 114
113 1F(GAMMA(II,JJ).GT.PHAIl.AND.GAMMA(IIiJJ).LE.DELTA1) GO TO 170
IF(GAMMA(I I,JJ) .GE.DELTA2.AND.RAMMA(II,JJ),LT.(a.*PI*PHAl£))
A GO TO 1/0
GO TO 150
11V IF(GAMMAtII,JJ) ,GT.PHAI1.AND.GAMMA(II,JJ).LE.DELTA1) GO TO 170
IF (GAMMA (I I, JJ) ,GF..DELTA2.AND.r,AMMA(II»JJ) ,LT.2.*PI .OR.
A GAMMAdl.JJ) .GE.O. . AND.GAMMA (II , JJ) .LT.PHA I ?) GO TO 170
GO TO 150
120 IF(GAMMA(II,JJ>,GE.PHAI2.AND.GAMMA(II,JJ).LE.PHAI1) GO TO H2
IF (GAMMA ( II, JJ) .GE.DELTA2.ANO.r,AMMA(II,JJ) ,l,T . PHAI 2. OR,
A GAMMAdl.JJI ,r,T.PHAIl.ANO.GAMMA0.
IFlDOX.GT. WIDTH) GO TO 150
GO TO 160
) *COX»(SPPFAO(J,K,M)«60.-WORO/(2.«lflS2.) )
c«*j
IF(DDX.GT. WIDTH) fit) TO 150
O
..
GO TO 160
]49 AA = C(JIfJJ)»C(TI.JJ>*^(JjK,'"M<>H<'MJtK,M)-?.*C(IT , J J ) »R2 ( J ,K , M )
A «COS(ALHHA(J,K,M)-GAMMM(II,JJ) )
A=SQRT(AA)
lF(A.LF..ti GO TO liS2
150 CONG(II,JJ)=0.0
SCONGI II , JJ) =0.0
CONG2(II« JJ)=0.0
SCONG2 ( I t , JJ) =0.0
INDEX ( II' JJ) =1
GO TO ?4V
165 CDX = C(II»JJ)«COS(GAMMA( 1 1 , J J) -ALPHA ( J ,K ,M) )
152
CONTINUE
-------
JF (COx.tiE.R? ( J.K.M) ) GO TO 150
,K,M) *CONV«(J,K,M)»H(J»K,M)
DHEMVH=DHEMV/H(J,K,M)
VG1H=VG1 ( J.K.M) /M ( J,K,M)
VG?H = \/G2( J.K.M! /H( J,K,M)
D4sEXP<-VGlH»Cnx»1852./llriAR(J,K.M) )
D5»EXP(-VG2H»CrX'M852./UHAp/H(,J,K,M) )
IF(U8AR( J,K,M) .hQ.0.0) GC TO 140
IF(R2(J,K,M) .Eo.o.O) GO TO 140
IF(J.EO.l) GO TO 161
SlGH = 2.»60.»lfir>2.» (SPREAO(J-1,K,M) * (COX/H2 IJ»K ,M) )«( SPREAD ( J.K.M)
A -SPREAOfJ-l.h.M) ) )
CONG! I! t JJ) = ( (COMC( J-l,K,M)»wOR( J-l »K.M) ) /(WQRO*SIGH) ) *04«07«DP4
SCQNGdl «JJ) = ( (SCONCl J-l ,K,M) *WO« ( J-l t K,M) ) /( -JOPO + SIGH) ) »D5»DPS
A * (1 ,5«(.ONVR(J,K,M) »H(J,K,M) /DREMV)
A »CONC(J,K.M) <
INDEX(II»JJ)=2
GO TO 34V
161 SlGH=(2.*60.«lf«52.»SPREAD( J,K,M) »COX/R2 ( J ,K , M) )
CONGf II, JJ) =( (OHCON(IO) "WOHO) /(WORO+SIGH) ) «04»07»f)P4
SCONG(II»JJ) =( (SOBCON(IO)*WOHO) /(WOrtO»SlGHn*05»nP5
*{I.5»CONVR(J,KiM)»H(J,K,M)/DHEMV)»CONC(JfK.M)«D6
GO TO
350 CONG(IItJJ)=ORCON(lO)
SCONG(IItJJ) =SOHCON( 10)
GO TO 24 V
162 CONG(II,JJ)=CONC(.J,«,M)
SCONG(IItJJ) =?CCNC(J,K,M)
INDEX(II»JJ)=2
GO TO ?<*^
140 CONG(II,JJ)=0.
SCONGf tltjj)=0.0
INDEX(II»JJ)=1
249 IF(INDEX
2*8 CONTINUE
250 CONTINUE
DO 661 11=1.11
00 660 JJ = 1. ?l
CONG( I It JJ) =0.0
SCONGf II »JJ) =0.0
CONG2(IItJJ)=0.0
SCONG2(IItJJ)=0.0
6AMMA ( I I» JJ) =0.0
C(II.JJI=U.O
660 CONTINUE
661 CONTINUE
650 CONTINUE
700 CONTINUE
800 CONTINUE
REWIND 51
DO 795 Il3l.lt
WRITE (51) BLAT(U)
795 CONTINUt
DO 796 JJ=1,21
WRITE(^l) HLON(JJ)
796 CONTINUE
REwINU
SO
00 802 11=1,11
WRITE (50) (XXJ(IT.JJ) ,SXXJ(IT,JJ) ,XXJ2(II,JJ) ,SXX J2 ( 1 1 « JJ) .
1 OxxJI 1 I »JJ) »OSXXJ( I I , JJ) ,PXXJ(II, JJ) tPSXXJdl. JJ) « JJ = 1 ,21 )
802 CONTINUt
153
-------
RE«rlNf) 5?
WRITE C>2) ADSO?,AOSULF,ftPS02»APSULF*ATSUL
SUBROUTINE CONHUL ( X , > ,NPTS, IP. IP? t IER) „ „ „
COMMON/T«IAN1/I1.NP,ASCALE,YSCALE,XMIN,YMIN,XMAX.YMAX,NL«NT
LOGICAL EQUAL ,START,l.AH
DIMENSION X0 TO 111
DO 107 I r«IL = l«N»TS
IF(IPdTAIL) .f.Q.O) GO TO 108
107 CONTINUE
lOfl CONTINUE
109 I2 = IP(ID
JF( I2.EiJ.O) RO TO 111
IF(Y(I1).EQ. Y(I2) ,AND.<(I1) .EQ.Xdaj I GO TO HO
11 = 13
GO TO 10V
110 NP=NP-1
IER=1
1500
IP(I1)=I^(I2>
IP(I2)=0
IP(ITAIL)=I2
ITAIL=I2
GO TO 10*
111 CALL PVEC(IP,IH)
IF(NP,GE.4) GO TO 113
113 CONTINUE
1501
RETURN
113 1=IP(1)
XMIN=X(I)
I=IP(NP)
XMAX=X(I)
YMIN=Y(D
00 114 I=
YMAX*AMAX1 (YMAX.Y ( I ) )
YMINsAMlNl (YMIM,Y ( I) )
114 CONTINUE
DX=XMAX-XMIN
DY=YMAX-YMIN
IF(DX.NE.O. .AND. OY.NE.O.) GO TO 115
IEP=33
WRITERS, 1502)
1502 FORMAT(1X,»IER=33»)
RETURN
115 XSCALE=. '07/OX
YSCALE=. /07/OX
00 116 1=2, NP
J=IP(I)
IF (X ( J) .NE.XMIN) GO TO 117
1 16 CONTINUE
IER=34
1S03 FORMAT ( IX, »IER=")'»
-------
120 IF (ABS(DX2) + ABS (DY?)-AHS (DX3) -AHS(DYJ) ) 12<>»122,121
121 IP GO TO 124
11*11*1
IF { .NOT.STAPT.OR. II .LT.3) GO TO 123
START=. FALSE.
IP(NPTS+1>=IP(NP*1 )
IP(NP*1)=IP(1)
NP=NP*1
123 IF(NP-I1.GE.2) fin TO 118
1ER=35
WRITE<6,1504)
1504 FORMAT <1<<.»IER=35*)
RETURN
124 CONTINUE
RETURN
END
SUBROUTINE DETCON ( I , J, K , SOLO, SUF T. SUL02. SUFT2, CSULD,CSUFT,CSULO;>,
Z CSUFT?)
IJ=MOP< (1-21) .4) *1
GO TO (100,200.300.400) ,IJ
100 IF(K.E0.1 ) GO TO 110
tF(K.EQ.2 ) GO TO 130
200 IF(K.E(J.l ) GO TO 210
IF(K.EO. 2) GO TO 230
300 IF ) GO TO 212
GO TO 211
410 GO TO 212
430 IF(J.LE. 6 > GO TO 212
GO TO 2H
END
COMMON /INpuri/ PSflOO)
CALL TRIANG(X, Y.N.I TR 1,1 SUP, IUSED,IE«)
RAMAA=FI
RAMB8=F1
RETURN
END
155
-------
INTEGER riEAO
N=0
I=HEAD
101 IF(I.tQ.O) GO TO 102
N = N» 1
J=LIST ( I)
LIST (I) =<*
I*J
GO TO 101
102 IF(N.EO.J) RETURN
J = N
K=LIST(J>
IF(K.E<}. N) GO 10 JOS
103 I=LIST(J)
IF(I.GT.O) GO TO \04
LIST(J)=-I
GO TO 105
104 LlST(v))=-K
K = J
J=I
IF(J.NF..N) GO TO 103
105
GO TO 102
END
SUBROUTINE REMOVE (I tM,K » U . VGA, VGH , TH'AH, HGT, KAMA 4 , RAMBR.
z LTAAT,LMAAT,LHAAT,AX,/\Y, A|_AT "•
IJ = MOD«I-21)«M*1
GO TO COO. ?UO, 300, 400) ,IJ
100 IF(K.EQ.l) GO TO 110
IF(K.tQ. GO TO 111
IF(J.UE.4) GO TO U2
IF(J.LE.6 ) GO TO 113
GO TO 11*
130 !F(J.LE.«i) GO TO 131
IF(J.I_F..4) GO TO 132
IF(J.LE.6 ) GO TO 133
GO TO 134
210 IF GO TO
IF(J.LF.4) GO TO
IF(J.LE.6 ) GO TO 213
GO TO 214
230 IF(J.LE.«J) GO TO 231
IF(J.LE.4) GO TO 232
IF(J.LE.6 ) GO TO 233
GO TO 234
310 IF(J.LF..2) GO TO 311
312
GO TO
IF(J.|_E.6 ) GO TO 313
GO TO 31*
330 iFfj.LF.O GO TO 331
IF(J.LE.4) GO TO 332
IF(J.LF..O ) GO TO 333
GO TO 334
410 IFU.LE.O GO TO 411
IF(J.LE.4) GO TO 4U
IF(J.LE.6 ) GO TO 413
GO TO 414
430 IF(J.LE.
-------
132 GO TO 111
133 GO TO 112
134 GO TO 111
211 GO TO 111
212 GO TO 112
213 VGA = f).Q
VGH=0.0
-
CALL PRtClP(AX,AY,A|_ATT.ALOW ,RAMAA ,RAMBB, N)
RAMAA = RAMAA«0.11'>1 , OE5/j6flO.
RAMBB=HAMHB«0.01»1 . OL5/3frnO.
TRAW=0.01»(1
RETURN
314 GO TO 112
231 GO TO 11
INTEGER DIRECT
LOGICAL LAB
ASSIGN 101 TO NEXT
LA8=.FALSE.
RETURN
ENTRY SO*T2
GO TO NEXT,(101,103.121)
101 CONTINUE
M=l
DIRECT=Q
1 = 1
JHEAO=1
LI5T(1)=0
IHEAD(1)=0
102 IF(I.GE.N) GO TO 114
11=1+1
L1 = I
L2=I1
ASSIGN 103 TO NEXT
RETURN
103 IF(ITE«5T) 105,104,111
104 IF(OIRF.CD 113,106,107
105 IFfOIKECn 110.10tMl07
106 DIRECT=1
JHEAD=I
107 LIST(D=I1
IHEAD(I)=0
108 1=11
GO TO 102
109 LIST(D=0
IH€AO(I)=0
110 IHEAD(M)=JHEAO
M = M»1
DIRECT=0
GO TO 10H
111 IF (DIRECT) 113,H2,10
-------
117 M=MP
IF(M.NE.I) GO TO 115
L1=IH£AU (1)
IHEAO(l) =0
RETURN
118
IF(l.EQ.l) GO TO 117
IHEAD(I)=0
GO TO 11'
119 J1=IHEAU(I)
J2sIHEAD(I*l)
IHEAD(I)=0
IHEAD(I*U=0
I2 = J2
NCHAIN=0
120 L)»J1
ASSIGN 121 TO NEXT
RETURN
121 IF(ITEST.LE.O) PO TO 1?2
IF(NCHAIN.EQ. 0) IHF.AO>
122 IF(NCHAIN.EQ.O) IHEAD(HP)sJl
IF(NCHAIN.EQ.?) Llsr(I?)=Jl
NCHAINal
JlsLIST(Il)
IF(JI.NE.O) GO TO 120
LlST(Il)=J2
GO TO 116
END
SUBROUTINE STEST ( X . Y,LlNE,ITRI,LL«FLIP)
COMMON/T^IANIX II«NP,XsCftLE,YSCALE,XMIN,YMIN,XMAX,YMAx.NU,NT
COMMON/TNIAN2/ L<^,L3t L*«L5
DIMENSION X(l)fY(l)»LINE(5,l),ITBI(3,l)«E(5)iINEXT(3)
LOGICAL f'LIP
DATA lNEXT/2,3,1/
DATA FAC/1.E13/
DATA FACI/l.E-13/
OX(I.J)aXSCALE»(X(J)-X(I))
DY( I . J)=YSCALE»(Y(J)-Y(I))
EL(I)=FACI«FLOAT(IA6S
J2=LINE(2»L1)
DO 101 I=H,3
IF(IABS(ITRI (I-l,KTl) ) .EC.L1) GO TO 102
101 CONTINUE
1 = 1
102 L2=ITRI (I.KT1)
lalNEXT(I)
L3=ITRI ( I.KT1)
IF(L2.LT.O) GO TO 103
GO TO 10*
103 L"-L2
J3=LINE(1.L)
DO 10S 1=2,3
IF(IAHs]lTRI(I-1,KT2)) .EQ.L1) GO TO 106
104 DO 10S 1=2,3
IF(IAHS(
105 CONTINUE
151
106 L4=ITRI(I»KT2)
L5=ITRI (I»KT2)
GO TO 107
GO TO 10H
107 L=-L4
108
.
IFfAl.LE.O.) RETURN
A2 = DX(J4«J2)
IF(A2.LE.Oi.)
158
-------
LL2=IABS(L2)
E<2)=ELsDX(Jl,J2)
DY1?=OY(J1,J2)
S3«(DX12»UY(J1.J3)-OY12»OX(J1,J3) ) / Hf-TUHN
LINE(1«H) =J4
LINE(2,L1)=J3
LINE(3,LI)=KTl
LINE(b,H)=E43»FAC
1=3
IF(L2.LT. 0) I=t
LINE(I.LL2)=KT2
IF(L4.LT.O) 1=4
LlNE(I.H.4)sKTl
ITRI(1»KTI)=L1
ITRI(2,KT1)=L3
ITRI (3,KT1)=L4
.ITRI (3.KT2)sL2
FLI
RETURN
END
.
T=»SX«DY-UX«SY
IF (T.fST. TOU! bO TO 10S
IF (T.LT.-TOL ) GO TO 1C3
IF(A8S(SY) .GT. APStSX)) GO TO 101
R=OX/OY
60 TO 102
101 R=DY/SY
10^ IF(n..(_E.H .AND. R.LE.l.) GO TO 104
103 LTEST=-1
RETURN
104 LTEST=0
RETURN
105 LTEST=1
RETURN
END
SUBROUTINE TMESH ( X . Y.NPTS, IP, LINE . 1 TRI , IER)
COMMON/T-
-------
i
f (I.EQ.I1)
TRI(1,I)*I
ITRI <3,I)=-(I*I1-1>
IF 1TRT(3,I)=-*|L
103 CONTINUE
DO 10^ 1=1, ML
01=LINE(1»I)
J2*LINE(2,n
LINE(5, I)=FAC»
1 SORT ( (XSCALE»(X(JD-X(J2) ) )«»2* (YSCALE»(Y ( J1)-Y(J2
2 M>*»2>
104 CONTINUE
DO 105 L=I1P1,NL
CALL STEST < X, Y ,L INE, ITRI »L .r L I PI)
105 CONTINUE
Kl=Kl*l
IFlKltGT. NP) r-,0 TO 14?
DO 141 K=K1,NP
DO 11J KT=1,NT
,
L13IT«I (l.KT)
IF (Ll.t-T. 0) r,0 TO 106
Jl»LINE(ltLl)
GO TO 10f
106 L=-L1
107 DXl = XSCALt«(X(.J)-X(Jl) )
OY1=YSCAUE«(Y(J)-Y(J1) )
SXl=XSCALt« (X U12) -X(Jl) )
SYl = YSCALt« ( Y( J<») -V( Jl) )
LT=LT£ST(SXl,SYl»Dxl»UYl )
IF(LT) 113, 10*. Id1?
108 L2=ITRI <-?«KT>
IF(L2.LT.O) GO TO
J3=LINF.(£,L2)
GO TO 110
10^ L=-L2
J3 = LINE(UL)
110 IF(LT) 113,114,111
111
0X2=0X1-5X1
DY2=DY1-SY1
IF(LTEST(SX2,SY2,nx2,DY2) ) 113.120,112
112 SX3=-(SX1*SX2)
SY3=-(SYl*SY2)
OX3=SX3*UX1
DY3=SY3»UY1
IER=36
WRITEC6. 15101
1S10 FOHMAT(1X,*IER=
RETURN
114 LINE(1,NL+1)=J1
LINE(2,NL*1)=J
LINEO,NL*1)=NT
LINE(4,NL+l) =KT
LINE(5,NL*1)=FAC«SORT (0X3 »
LINF. (2.NC + 2) =J
LlNE(3,NL+2) =KT
LINE(4,NL*2)=MT»1
LlNE(5,NL»2)=FAC*SORT(DX2«»2*nY2«»2)
LINE (l,NL+3) =J3
LINE(3iNL*J)=NT*l
LINE(4,Nl*3)=NT*2
IF(L2.LT. 0) fcO TO
GO TO 116
115 L=-L2
LINE(4,L)=NT*1
116 L3=ITRI(3,KT)
IF(LS.LT.O) GO TO 117
LlHE(3,LJ)=NT+3
GO TO 11*
117 L=-L3
LlNE(4,L)=NT+2
118 CONTINUE
ITRI(2,KT)=NL*?
ITRI (3.KT|s.(fJL + l )
ITRI (lfNT*n=L?
ITRI(2,NT*l)=Nl. + 3
ITRI(3,NT+l)s.(NL+2)
ITRI(1,NT»2)=L3
ITRI (3,NT*2) =- fNl.t
160
-------
CALL STEST(X,Y.LIMF»ITPI,L l,FI.TPl)
CALL STEST(X,Y,LIME,ITRI,L?,FLIP1>
CALL STEST(XtY.LINEiITHl,L3,FLIHl>
GO TO 140
119 L3*ITRI(J»KT)
GO TO 122
120 L3=L1
L1=L2
L2=ITRI<3»KT>
JT=J3
J1=J2
J2=JT
GO TO 122
121 L3=L2
L2=L1
L1=ITRI(3.KT)
JT=J2
J2 = J1
J1=J3
122 CONTINUE
IF(Ll.LT.O) GO TO 123
KT2=LINE(4,L11
GO TO 124
123 L=-L1
KT2=LINE(3,L)
124 IF(KT2.NE.O) GO TO 125
IER = 2
WRITE(6tl5ll)
1511 FORMAT(1X,»IEW = ?«')
GO TO 141
125 CONTINUE
IF(L1.EU.-ITRI(1,KT2)I JO TO t?6
IF (Ll.EQt-ITRT(2.KT2)) PO TO 1?7
IF (Ll.EQ.-ITRI(3.KT2)) GO TO l?fi
1512
RETURN
126 V>=ITRI
L5=ITRI (3,KT?)
GO TO 12^
127 L4 = ITRI(-3»Kr2)
L5=ITRI (1 tKT2)
GO TO 12^
12S L-+=I
12^ IF(L*.LT.O) GO TO 130
J4 = LINF. () ))»»?)
133 LINE(1,NL»2) =J3
LINE(2tNL*2) =J
LINE(3,NL+2)=KT
LINE (*,NL*2) =NT»1
1 SORT ( (^SCALE«{>( { J3) -X ( J) ) ) «»2+ ( Y SCALE* ( Y ( J3) -Y(J) ) )««2)
LINE (1,NL+1) =J1
LINE(2,NL*1)=J
LINE(3,NL*1)=NT«1
IF (KT2.EU.O) GO TO 134
GO TO
134 LINE (4,NL+1) =0
135 LINE(5,NI.+ 1)=FAC«
1 SQRTt (XSCALE* (X(J1)-X{J)))»«2*(YSCALE»(Y(J1)-Y(J))
IF(L3.LT«0) GO TO 136
LINE (3.L3) =NT+1
GO TO 137
136 L=-L3
LINE(4,L)=NT»1
137 ITRI(1,KT)=L1
ITRI (2,KT) =L2
ITRI (3,KT)=NL + ?
ITRI (l.NT*l)=NU + l
ITRI (2,NT + 1)=-(NL*?)
ITRI <3,NT+i)=L3
LINE(1,NL»3)=J4
161
LINE (4,NL
-------
IF(L4.LT« 0) GO TO 138
LINE <3.L<*1 = NT*?
GO TO 13^
138 L = -L*
LINF. <4,L>=NT + ?
139 ITRI (l,NT+2)=-(NL*l)
ITRI (2,NT+2)=L4
ITRI(3,NT*2)=NL»3
ITRI <1,K12)=-LI
ITRI <2,KT2)=-(ML+3)
CALL'STEST (xt?, LIMP. ITRI. L3,FLi Pi)
L
CALL STEST ( X , Y.L T NE , I T*I ,L4 .
'
" Y, L
CALL STEST(X,Y,LINF.,ITPI,U1,FLIP1) ,
CALL STE^T (X,Y.LTNEtIT*I.ML+l.FLIPl>
5 Y, LINE, ITRI, NL*3,FLIP1>
140 CONTINUE
NL=NL+3
NT = NT*2
141 CONTINUE
142 CONTINUE
RETURN
END
SUHROUTI^E TMESHI(XtY,LlN(i,lTRIfIEK)
COMMON/T^IANl/NO.MP.XSCAl.E.YSCALE.XMIN.YMIN.XMAX.YMAX.NL.NT
COMMON/TWJAN3/ LL(4)
DIMENSION X(l) ,Y (1 ) tLINEfb,!) , IT" I (3,1 )
LOGICAL I-LIP
LM=N8*1
KOUNT=3«NL
101 CONTINUE
LMIN=LM
LM = NL
00 103 L=LMIN.NL
It ( LINE(S,L) .LT.O) GO TO 10T
CALL STEST
KOUNT=KOUNT-1
iFfKOUNT.tQ. 0) PO TO 104
LINE(5,L)=ISIGN(LINE(5,I.).-I)
IF(. NOT. CLIP) GO TO 103
DO 102 1=1,4
LI = IABS(LL( I) )
IF(LINEC*,LI) .EQ.n.OH.LlNF,(3.1 I) ,EO.O) GO TO 10?
LINE<5,LI)=IA3S(UNE(S,LI) )
IF(LI.LT.U) UMrMIMO(LM,Lt)
102 CONTINUE
103 CONTINUE
IF(LM.LT.NL) GO TO 101
GO TO 10b
104 CONTINUE
WRITEI6.1S1S)
1515 FORMAT! 1X»MER = 3»)
105 CONTINUE
LMIN=NH* I
DO 106 L=LMIN,NL
LINE(!s,L)=IAHS(LINF. (b.D)
106 CONTINUE
RETURN
END
SUBROUTINE TRIANG ( X. Y , N, I TRI . I SUP , IUSER, IER )
DIMENSION X(M) .Y(N), ITRI(ISHP)
COMMON/T-iIANl/Nie,MPTS,XSCALE,YSCALE«XMIN,YMIN,xMAX,YMAX,ML,NT
CALL CO^HUL (X.Y.'M. ITRI (?) , ITRI (N*3) ,ItR)
IF (IF.R .O.E. 32) RETURN
ITRI (1) = NB
NL = 2«Nd*3«(NPTS-Nd-2)
NT =NtU
-------
DO 102 1=1, NL
J1=ITHI(IU)
J2=ITRI (IL+l)
DX=X(J2)-X(J1)
DY=Y(J2) -Y(J1)
DXX=XX-X (Jl )
DYY=YY-Y
IFtD.LT.U. .OR. O.GT.OP) <*0 TO 101
D=DXX»DXX+DYY»OYY-G<»D/OH
IF(D.LT.U. ) 0=0.
IF (D.GE.UMIN.ANi'.pMlN.GE.O. ) GO TO 101
DMIN=0
ILSAV=IL
101 IL=IL*5
103 CONTINUE
K=ITRI <2>
0*)«
TRINTR«Z(J1) »K1+Z( J2) *F?+Z( J3)»F3
RETURN
103 T= (DX«DXX*DY»OYY)/ (DX«DX+DY»DY)
RETURN
104 CONTINUE
TRINTR=Z»XXJ2UI,JJ),SXXJ2(II,JJ),
163
-------
1 UXXJ ( T I »JJ) ,nsXKJ ( II »JJ! »»XXJdIiJJ) »PSXXJdI»JJ) « JJ=1 ,
200 CONTINUE
DO «02 11=1,11
DO 801 JJ = I. 21
XXJ( II ,J>J) =XXJ ( II.JJ) /FLOAT (JJ.I)
SXXJ(II«JJ)=SXX.J (Hi JJ I/FLO AT(JJJ)
XXJ2(IIfJJ)=XXJ2(IIiJJ)/FLOAT/2,0 + nsXXJdI, JJl/3.0
03 = l. (BLON(JJ) , JJ=17,?1)
925 FORMAT ( 1^1 ,8X,5f"8. 1/////1
DO 913 1 1 = 1 ,11
WRITE (ft, *41) ALATdl) , (XXJ(II,JJ) ,JJ = 17,21)
941 FORMAT (8LON(JJ) ,JJ=17,21)
00 952 11=1,11
WHITE! ft.^*l) BL AT(II) , (SXX.JtII ,JJ) ,JJ=17,21)
952 CONTINUE
953 FORMAT ( 1^1 ,HX, »OF.POS I T I ON OF «;02 ( KG/KM2 )»////)
WRITE(6»l)l5) (PLON(JJ) ,JJ=l,lf,)
00 P54 11=1,1]
WRITEX6iy04) HLAT (II) , (UXXJI II , JJ) ,JJ=1, 16)
ITE(6,^25) (t'LON(JJ) ,JvJ=17.21)
955 11=1,11
ITE<6.*41) HLAT(II),(DXXJ(II,JJ) ,JJS17,21)
955 CONTINUE
956 FORMAT (im, fiX, ^DEPOSITION OF S04 (KG/KM2) *////)
WRITE (^,^15) (BLON( JJ) ,J,J = l,lft)
00 957 11=1,11
WHITE (6,904) RLATdI>»(OsXxJ(II»JJ)»JJ = l»lf>)
957 CONTINUE
WRITE(6,^25) ( HLON(JJ) ,JJ = 17,?1)
DO 958 11=1,11
WRITE (ft, 'Ml) HI. AT (I I) , (OSXXJdl »JJ) ,JJ=17,21)
95H CONTINUE
WRITE (ft,^7Q)
970 FORMAT(lHltHX,«CONCENTH«TION OFS02 ABOVE B. LA YER (M IGRAM/M3) »
Z////)
WRITE (ft, ^15) (4LON( JJ) , JJ=1 , 16)
DO 971 11=1,11
WRITE (ft, ^04) f»I..U(II),(XXJ?(ITfJJ)»JJ=l,\6)
971 CONTINUE
WRITE (ft, '
-------
'/.MITE (
980 hORMAT ( 1^1 , 8X , ^COMCENTnM f I0'j OF SULFATE ABOVE B . L AYEP ( Ml GRAM/M3 ) <
981 ' T ' ' ' ' SX *J? ' " ' JJ ' ' JJ= l '
goI3i^iii?!irLON'JJ=17';'1>
982 CONTINue"411 "l "T ' l l' ' ' SXX J? ( 11 ' JJ > ' JJ= 1 7 ''2l '
WRITF (-5,700)
700 FORMAT dm ,ax, !>nF.posisiON OF 502 DUE TO HAI,
00 710 I 1 = 1,11
riL*T(!I),;PXXJ(II,JJ),JJ=l,16)
710
gg'Tinj
1 ^LATtII),(PxXJ(IItJJ),JJ=17,21)
Ht
CONTTNUc
WRITF (6,730)
730 FOHMAT(lHl,»Jx,*nEPOSISION OF SULFATE DUE TO PA IN (KR/KM^J #////)
740 CONTINUE9041 HLAT«"'"pSXXJ(II,JJ),JJ = l
WRITE (6.925) (RLON(JJ) . JJ=17,21)
750
.
850 1'1D'<>A'0llNT °r SULFUR(KG) IN THt
851
WRITF.(6tH«;0)
852 FORMAT(10X.»RAIN OEP03ISTOK OF SULFUR ( KG) =» ,F 1 5.
REMIND 5?
AOS02,AOSULr , APS02,APSULF,ATSUL
WRITE(6.1111) ADSO?
WRITE(6,l]l2) AOSULF
WRITE(6,m3) APS02
TE (b, l ] IS) ATSUL
1111 FORMAT MH] , IPX, «UMY DEPOSITION OF SO? ( KG ) =» . E 1 •=, . 5 )
}J}| £OW^I OX,«nRY nEPCSlTION OF SCH
-------
APSO?=0.0
APSULF*0.0
_
REWIND 52
WRITE(i>2) ADSOa.ADSULf , APSO? , APSULF .ATSUL
STOP
END
166
-------
TECHNICAL REPORT DATA
(Please read Instructions on the reverse before completing)
1 REPORT NO
EPA-600/4-79-068
3 RECIPIENT'S ACCESSI ON-NO.
4. TITLE AND SUBTITLE
LONG-RANGE TRANSPORT AND TRANSFORMATION OF SO,
AND SULFATE i
5 REPORT DATE
November 1979
6. PERFORMING ORGANIZATION CODE
7. AUTHOR(S)
Teizi Henmi and Elmar R. Reiter
8. PERFORMING ORGANIZATION REPORT NO.
9. PERFORMING ORGANIZATION NAME AND ADDRESS
Colorado State University
Fort Collins, Colorado 80523
10. PROGRAM ELEMENT NO.
1AA603A AE-009 (FY-78)
11. CONTRACT/GRANT NO.
R-805271
12. SPONSORING AGENCY NAME AND ADDRESS
Environmental Sciences Research Laboratory - RTP, NC
Office of Research and Development
U.S. Environmental Protection Agency
Research Triangle Park, North Carolina 27711
13. TYPE OF REPORT AND PERIOD COVERED
Final 5/77 - 4/79
14. SPONSORING AGENCY CODE
EPA/600/09
15. SUPPLEMENTARY NOTES
16. ABSTRACT
Technical descriptions and computer programs are presented for two models that
calculate long-range transport, diffusion, transformation of S02 to sulfate,
and dry and precipitation deposition of initially emitted S02. One model treats the
mixing layer height as constant; the other (at the expense of computer time) varies
the mixing layer height diurnally and tracks pollutants in three layers—the daytime
mixing layer, the nocturnal ground-based stable layer, and the daytime mixed layer
that remains above the nocturnal stable layer. Application of the multi-layer model
over a region encompassing the Ohio River Basin produced patterns of S02 and sulfate
concentrations that are statistically correlated with observed concentrations.
An empirical formula for the transformation rate of S02 to sulfate is derived
and used in calculations of regional residence times of S02 (TS02) and sulfate (TSUL)
for the U.S. east of 105°W longitude. TS02 ranges 15-30 and 15-40 hours for the cold
and warm seasons, respectively; TSUL ranges 150-450 and 200-500 hours for the cold
and warm seasons, respectively.
Using a cumulus cloud model, results showed that sulfate aerosol capture by
cloud water through microphysical processes is sufficient to produce observed levels
of sulfate in rain water.
17.
KEY WORDS AND DOCUMENT ANALYSIS
DESCRIPTORS
b. IDENTIFIERS/OPEN ENDED TERMS
c. COS AT I Field/Group
Air pollution
Sulfur dioxide
* Atmospheric diffusion
Conversion
Distance
Sulfates
Computer programs
13B
07B
04A
09B
13. DISTRIBUTION STATEMENT
RELEASE TO PUBLIC
19. SECURITY CLASS (This Report)
UNCLASSIFIED
21. NO. OF PAGES
183
20. SECURITY CLASS (Thispage)
UNCLASSIFIED
22. PRICE
EPA Form 2220-1 (9-73)
167
------- |