JUNE  1985
USER'S GUIDE FOR THE ADVANCED STATISTICAL TRAJECTORY
       REGIONAL AIR POLLUTION (ASTRAP)  MODEL
      ATMOSPHERIC SCIENCES RESEARCH LABORATORY
         OFFICE OF RESEARCH AND DEVELOPMENT
          ENVIRONMENTAL PROTECTION AGENCY
   RESEARCH TRIANGLE PARK, NORTH CAROLINA  27711

-------
USER'S GUIDE FOR THE ADVANCED STATISTICAL TRAJECTORY
       REGIONAL AIR POLLUTION (ASTRAP)  MODEL
                         by
                  Jack D.  Shannon
          Environmental  Research Division
            Argonne National  Laboratory
                 Argonne,  IL   60439
            Interagency Agreement Number
                    DW930060-01
                  Project Officer
                   Terry L.  Clark
        Meteorology and Assessment Division
      Atmospheric Sciences Research Laboratory
         Research Triangle Park,  NC  27711
      ATMOSPHERIC SCIENCES RESEARCH LABORATORY
         OFFICE OF RESEARCH AND DEVELOPMENT
          ENVIRONMENTAL PROTECTION AGENCY
   RESEARCH TRIANGLE PARK, NORTH CAROLINA  27711

-------
                                  NOTICE
     The information in this document has been funded by the United States
Environmental Protection Agency under Interagency Agreement Number DW930060-
01 to Argonne National Laboratory.  It has been subject to the Agency's
peer and administrative review, and it has been approved for publication
as an EPA document.  Mention of trade names or commercial products does
not constitute endorsement or recommendation for use.
                                    ii

-------
                                 ABSTRACT
     The Advanced Statistical Trajectory  Regional  Air Pollution (ASTRAP)
model simulates  long-range,  long-term  transport  and deposition  of  air
pollutants, primarily oxides of sulfur and nitrogen.  The ASTRAP model is
designed to combine ease of exercise with an appropriate detail of physical
processes for assessment applications related to acid deposition.

     The theoretical basis and the  computational structure of the ASTRAP
model are described.  The model  is  evaluated with  observed  data,  and an
illustrative exercise is provided.   Computer codes  of the model subpro-
grams are provided in appendices.

-------
                                 CONTENTS







Abstract	   i i i




Figures	    vi




Tables	   vii




Symbols and Abbreviations	  viii




Acronyms	    ix




Acknowledgements	     x




Executive Summary	     1




1.  Introduction	     3




2.  Data Requirements	     5




3.  Features and Limitations	     7




4.  Technical Description	     8




       Horizontal Dispersion	     8




       Vertical Dispersion.	    12




       Concentration and Deposition	    21




       Preprocessors and Postprocessors	    23




       Evaluation	    32




5.  Computer Aspects	    32




6.  Test Application	    39




References	    45




Appendix A:  Program HORZ listing	    48




Appendix B:  Program VERT listing	    59




Appendix C:  Program CONCDEP listing	    70
                                    v

-------
                                    FIGURES
Number                                                                     Page
   1   Grid systems for ASTRAP input and output	    5
   2   Time series of horizontal dry puff densities	    10
   3   Time series of horizontal wet puff densities.	    11
   4   Vertical eddy diffusivity for summer conditions in ASTRAP	    113
   5   Diurnal patterns of SCU dry deposition velocities in ASTRAP	    14
   6   Diurnal patterns of SO^ dry deposition velocities in ASTRAP	    15
   7   Diurnal patterns of N0/N02 dry deposition velocities in ASTRAP...    16
   8   Diurnal patterns of NO-j dry deposition velocities in ASTRAP	    17
   9   Diurnal patterns of the rate of transformation of SCU to
         S04 in ASTRAP	    18
  10   Diurnal patterns of the rate of transformation of NO/NC^ to
         N03 in ASTRAP	    19
  11   Comparison of ASTRAP simulations with observed annual wet
         sulfate deposition during 1980	    26
  12   Comparison of ASTRAP simulations with observed annual wet
         nitrate deposition during 1980	    27
  13   Comparison of ASTRAP simulations with observed average SOo
         concentrations for January and July 1978.	    28
  14   Comparison of ASTRAP simulations with observed average SO^
         concentrations for January and July 1978	    29
  15   Comparison of ASTRAP simulations with observed wet sulfate
         deposition for January and July 1978	    30
  16   ASTRAP program structure	    33

-------
                                    TABLES
Number                                                                    Page
  1   Early dispersion parameterization limits	    21
  2   ASTRAP wet deposition verification statistics	    31
  3   ASTRAP verification statistics for 1978 simulations	    31
  4   ASTRAP input	    34
  5   ASTRAP input/output files	    38
  6   ASTRAP test case S02 calculations	    40
  7   ASTRAP test case SO^ calculations	    41
  8   ASTRAP test case dry deposition calculations	    42
  9   ASTRAP test case wet deposition calculations	    43
 10   ASTRAP test case integrated deposition budget	    44
                                     VII

-------
                           SYMBOLS AND ABBREVIATIONS
V
H (x,y)
Hw (x,y)
R
P
N
M
b
L
E
c
C
n
m
d
D
W
S.D.
Ob
Mod
Res
Pom
 z
mean position of a puff on the horizontal grid
one standard deviation spread of puff in x and y
correlation between x and y in the puff
two-dimensional density function of horizontal dry puff
two-dimensional density function of wet-deposited puff
fractional bulk removal by precipitation
precipitation rate in cm/6 hr
number of source cells
fraction of unit emission still aloft
M in the case of no previous wet deposition
number of emission layers
emission rate in tonnes per six hours
one-dimensional surface concentration (yg m~ )
three-dimensional surface concentration (ug m  )
number of trajectories contributing to a dry deposition statistic
number of trajectories contributing to a wet deposition statistic
dry deposition increment per time step
                              _o
cumulative dry deposition (g m  )
                              _o
cumulative wet deposition (g m  )
standard deviation
observed value
modeled value
residual (Ob-Mod)
correlation between observations and modeled values
vertical eddy diffusivity in meters squared per second

-------
                                   ACRONYMS
ASTRAP      Advanced Statistical Trajectory Regional Air Pollution  model
VERT        VERTical dispersion subprogram
HORZ        HORZontal dispersion subprogram
CONCDEP     CONCentration and DEPosition subprogram
NMC         National Meteorological Center
AES         Atmospheric Environment Service of Canada
SURE        SUlfur Regional Experiment
MOI         Memorandum of Intent between the U. S. and Canada on  trans boundary
            air pollution
UTM    :    Universal Transverse Mercator map projection
LRT    :    Long-Range Transport and deposition
NAPAP  :    National Acid Precipitation Assessment Program

-------
                                ACKNOWLEDGMENTS

     This research has been  funded as part of the National  Acid Precipitation
Assessment  Program  by  the  U.  S.  Environmental  Protection  Agency,  under
Interagency  Agreement   DW  89930060-01   between   the  U.  S.   Environmental
Protection Agency and the U. S. Department of Energy.

     Dr. Barry Lesht has provided many helpful suggestions.  The assistance of
Mr. Steve Fine  in reorganizing the  coding was  extremely valuable.   Typing of
the user's  guide  by Mrs.  Rebecca  Spencer  and  Ms.  Anh Tran  is  gratefully
acknowledged.   Critical reviews  by Mr.  Terry  Clark  and  Ms.  Joan  Novak have
been quite useful.

-------
                               EXECUTIVE SUMMARY

     The Advanced Statistical Trajectory Regional Air Pollution  (ASTRAP)  model
calculates  atmospheric  concentrations  and  deposition of  oxides of sulfur  or
nitrogen  over spatial  scales  of  100  km to  continental, and  over  temporal
scales of one to three months.

     The  model  assumes linearity  between emissions  and deposition,  although
the linear rate itself varies diurnally and seasonally and can be  changed as a
user's option.  ASTRAP  is  not an appropriate model for  short-term or  episodic
conditions, nor for  near-source air quality.   Yearly simulations  can  be made
by  averaging  atmospheric  concentrations  and adding  deposition from  seasonal
simulations.

     ASTRAP is structured here  to treat as many as 730 source cells with  up  to
six emission  layers  and to  make calculations for a  51 x  45  grid  across the
United  States  and   Canada,   with  grid  spacing  100-130  km,   depending   upon
latitude.   Preprocessed emission fields  and time series  of  continental  wind
and precipitation analyses are  required  in order to exercise  the model.   If
the  user   desires   to  change  the  size  of   the  emission   grid   or the
concentration/deposition grid,  simple coding changes would be necessary.

     The deposition  fields are  integrated  across each  state and  province  in
turn,  but  in  this  version of  the code  all  sources  are  in a single inventory.
If  the  concentration and deposition from individual sources or source  areas
are to  be produced  separately  as  in a  source-receptor  transfer matrix,  some
coding changes would be required.

     The  subprograms in ASTRAP do   not  contain simple  options  to  select the
modeling alternatives  discussed above,  but  guidance  as  to particular  changes
that might be necessary is given in  the section on applications.

     Among the features of ASTRAP are the following:

        •  efficient long-term  simulations of continental  scale.

        •  diurnal  and  seasonal  variations  of  dry deposition   velocities,
           chemical  transformation rate, and vertical stability  profile.

        •  layers   in  the  vertical  to  allow  for  the   effect   of different
           effective emission heights.

        •  use  of   analyzed   synoptic  meteorology,  to  allow  examination  of
           meteorologically induced  variability.

        •  integrated deposition fields.

        •  leakage to the free  troposphere.

     The ASTRAP model consists of   three major  subprograms,  VERT, HORZ, and
CONCDEP, which  calculate  one-dimensional  vertical profiles,  two-dimensional
horizontal dispersion, and concentrations and deposition,  respectively.   On  an
IBM 3033,  VERT  requires  210K (bytes) storage and  5-10  minutes   of  computation
time,  plus  one  output  file,  for a  seasonal simulation.   HORZ  requires  500K
                                      -1-

-------
storage and  10 minutes  computation  time,  two input tape drives,  and  an  output
file for a seasonal  simulation.   CONCDEP  requires 500K storage and 10 minutes
computation  time,  an  emission  file,  the  output  files  from the  other  two
programs, a file used to  identify each receptor cell as to state or  province,
and an output file for  a  seasonal simulation.  If only the emission  inventory
or the receptor  grid is changed, only  CONCDEP  need be  rerun.   If the year  of
meteorology is changed,  both CONCDEP and HORZ must be rerun.   If  the  season  of
meteorology is changed,  all subprograms must be rerun.
                                      -2-

-------
                                1.  INTRODUCTION

     Acid  rain,  as  the  atmospheric deposition  of acidifying  substances  is
commonly  termed,  is only  one  of  many potential  environmental  and  health
problems   involving   long-range  transport   and   deposition   of   airborne
pollutants.   Investigations  of  the transport  and deposition of many  different
substances over  a  wide  range of spatial and temporal scales eventually may be
necessary, in order to develop and  improve assessments  of potential  actions
and  regulations.   A  common  approach  in  such  investigations  is  numerical
simulation of the  atmospheric  pathways.   The  Advanced Statistical Trajectory
Regional  Air  Pollution (ASTRAP)  model  described  here  has  been  developed  to
simulate  the  long-term (monthly  to  yearly),  regional-scale (resolution about
100 km) deposition of oxides of sulfur and nitrogen, the major contributors to
acid deposition.  Many  of the  modeling techniques in ASTRAP  are based  upon
previous work by Taylor (1921), Durst et al. (1959), Bolin  and Persson (1975),
and Sheih  (1977).   The  ASTRAP techniques can be extended to other  pollutants,
provided  that linear parameterizations of chemical transformation and  removal
processes  are suitable,  and  can be  extended  to other   spatial  scales,  if
appropriate meteorological data  can  be  obtained.   Unmodified extension  of  the
modeling techniques used in ASTRAP to shorter, episodic temporal  scales is  not
recommended,   because  of  certain statistical features of  the model, which  are
discussed in  following sections.

     The  ASTRAP  model consists  of three submodels and various  preprocessors
and  postprocessors.    Use  of   particular  preprocessors   and  postprocessors
depends  upon  application  and   data  availability.    The  submodels   are  the
vertical  diffusion program  (VERT),  the horizontal dispersion  program  (HORZ)
and  the concentration  and deposition  program (CONCDEP).   The  processes  of
vertical diffusion,  chemical  transformation,  and dry deposition are  simulated
in VERT through  parameterizations that are independent of  horizontal  location
or particular meteorological conditions.  In HORZ, the processes  of horizontal
advection, horizontal diffusion, and wet deposition are simulated with the  use
of  time series  of  wind  and  precipitation  analyses.   CONCDEP  combines  the
statistics produced  by the  first two  programs  with emission  inventories  to
produce fields of average atmospheric concentration and cumulative deposition.

     Data  requirements  for  ASTRAP  simulations  are  presented  in  Section  2,
while the  features  and  limitations of ASTRAP  are  sketched  in Section  3.   The
theoretical basis  of  the  ASTRAP model is  outlined in  Section 4.   The  primary
subprograms are described, preprocessors and postprocessors are discussed,  and
model verification statistics  are presented.   Specific  computer aspects  are
discussed in  Section  5  and a test simulation  is  presented  in Section 6.   The
computer  codes of  the subprograms, included in  appendices, are written  to  be
exercised with meteorological data from  1980,  as described  in Section  4.
                                     -3-

-------
                             2.   DATA REQUIREMENTS

     The data requirements for ASTRAP simulations are monthly to seasonal  time
series of  transport  wind and precipitation analyses,  and a source  inventory.
Ideally, wind and precipitation  analyses  should be at intervals of  six hours,
should have  a spatial resolution  of about 300  and  100  km, respectively,  and
should extend well  beyond the region of interest in  order to allow accurate
calculations  of  recurving  trajectories.    Unfortunately,  these  desirable
features are not always available.  Rawinsonde upper-air wind measurements  are
made routinely at intervals of 12 hours rather than six hours, and the spatial
coverage of  rawinsonde stations,  while adequate  over the  United  States  and
southern Canada,  is  sparse over  northern Canada and  essentially missing  over
oceanic areas.   This can cause  serious  difficulties  in  obtaining wind fields
suitable  for deposition  simulations  for  Florida or the  Canadian Maritime
Provinces,   for  instance.   A  similar  difficulty can arise  from   inadequate
spatial  coverage  of  precipitation  fields.    The   temporal  resolution   of
precipitation data  is usually  better  than for  upper air  winds,  since many
stations report hourly precipitation.

     Analysis of  a  one-layer transport wind field  requires definition of  the
depth of the  field.   There is no  universally  accepted transport layer, since
over  a  single  location  fresh pollutants might  be  transported in  a growing
planetary  boundary   layer  while   aged pollutants  might   be transported  in  a
deeper  layer  defined  by  the  previous  day's  maximum mixing  depth upstream.
Definition of a  mixed layer depth is adversely affected also by the practice
of making  routine soundings only  twice a day,  as there  is a typical diurnal
cycle of planetary boundary growth and decay.

     The objective  analysis  of  observed  winds to produce  a gridded analysis
does  filter  small-scale  wind  variations  from  the field,  but that does  not
appear  to  be a major  cause of  uncertainty  for  long-term simulations.    In
gridding  precipitation  observations,  however,  frequency of  precipitation
occurrence is increased as  the grid  increment  increases, since a cell average
is non-zero if any one of perhaps several observation  stations within the  cell
is  non-zero.    Without  adjustment by a  threshold  value,  this could  lead  to
overprediction   of   wet   deposition  by   ASTRAP,   since   the   wet  removal
parameterization  was developed  from observations  at a  point,  not across  a
grid, and is a function of both  precipitation  frequency and amount.

     Initial  preparation  of   wind  and  precipitation  fields  for  ASTRAP
simulations  has  been  performed  at  the University  of Michigan  by Professor
Perry Samson.   His wind  fields,  produced every  12 hours, are for 500 m layers
to 3000 m MSL for a  17 x 19 horizontal grid of  National  Meteorological Center
(NMC) spacing.   (The NMC polar  stereographic  projection, aligned along 80°  W
and true at  60°  N with 381  km spacing  there,  reduces to about 306  km spacing
at 30°  N.)   The wind grid is  the  coarser grid shown  in Figure 1, with values
of  west-east   and   south-north   components  at  the  cell  centroids   (i.e.,
components not  aligned along grid system).  The corners of the grid, too  far
from rawinsonde stations for the analysis method, are  coded as missing values,
as are  layers  that  were in  effect below the elevation of  the terrain.   Speed
is  given  in  meters  per  second.   Further processing has  been  necessary  at
Argonne before  the  analyses could be used  in  ASTRAP  simulations.   For  spring
and  summer,   the bottom 3  non-missing  values  have  been averaged  at   each
gridpoint;  for   fall and  winter the  bottom  2  non-missing values  have  been
                                     -4.

-------
Figure  1.  Grid  system  for  ASTRAP  input  and  output:  wind field  analyses
           (centers of  large cells);  precipitation analyses  (centered about
           intersections of  small cells with  complete  coverage  within dashed
           rectangle);  trajectory virtual  sources  (•);  emission field  and
           concentration/deposition  fields  (centers of small cells;  not  all
           cells contain  emissions).   North  Pole  (large  dot)  is  the (0,0)
           point  of  the  grid   system,  with  X  increasing  to   right  and  Y
           increasing downward.
                                     -5-

-------
averaged.  For grid  points  still with missing value indicators,  extrapolation
has been used to  produce  wind  values.  Finally, linear temporal  interpolation
has created fields at 6-hr intervals, and the wind analyses have  been  combined
into files  three months  long  (Dec  - Feb,  Mar  - May,  Jun -  Aug,  and Sep  -
Nov).   The wind  components are  converted  to  components along  the NMC  axes
internally in the HORZ subprogram.

     The precipitation analyses  produced  at the University of Michigan  are  on
a 50 x  45  grid  with 1/3  NMC spacing, as outlined in  Figure  1.  The  original
analyses are  hourly, with a special code  for  missing values.   Precipitation
amount  is  in  millimeters.   At Argonne the hourly fields have  been added  to
produce  6-hrly  fields, with  missing  values  filled  in  by  interpolation and
extrapolation.   The  precipitation data  are organized in monthly files.  Any
missing fields of either  wind  or precipitation have been replaced by  temporal
interpolations  in the  processing   at  Argonne,  so  there are  no  date-time
indicators written with the fields.

     For the version of ASTRAP in this user's guide, an SOo inventory  has  been
created in which  the seasonal  emissions in kilotonnes are given by effective
emission layer (0-100 m,  100-200  m,  200-300 m, 300-400 m, 400-600 m,  and  600-
800 m)  and  NMC  position  of  the lower  left corner of the  grid cell  (the NMC
coordinate  system is  0,0  at   the   north   pole).   No separate  SO/  emission
inventory  is  used,  so   primary  sulfate  emissions  factors   are applied  in
CONCDEP.  The primary  sulfate  emission factor for sources in  the bottom layer
is assumed to be 0.05 (i.e., one unit of SO™ equivalent emission  is  treated  as
0.95 units  of S02 and  1.5(0.05)  = 0.075 units  of  primary sulfate,  where the
1.5 factor arises because the ratio  of the  molecular weight of  sulfate to  that
of SOn  is 96/64); the  primary  sulfate emission factor for elevated  sources  is
assumed to be 0.03 in  Florida  and in the Northeast, and  0.015  elsewhere.  The
emission grid is  the same spacing as the precipitation grid,  but is irregular
because  the  inventory  is  organized  by  state  and province.    The  emission
information has  been derived from a  preliminary version of the  National  Acid
Precipitation Assessment  Program (NAPAP) inventory for  1980  produced  for the
1985 Assessment.   Wind,  precipitation,  and emission  data  are all written  in
binary format for efficiency.

     A  different  emission inventory for ASTRAP  would  need to be preprocessed
to the  form described above, or some of the model code would  require  changes.
If an NO   inventory  were to be  used,  this  version of ASTRAP would  expect the
inventory  to  be  expressed  as  equivalent  kilotonnes  of  N0~,   and  it would
require  a  different output  from  program  VERT (described  in  the  Vertical
Dispersion subsection of  the Technical Description section).
                                      -6-

-------
                         3.  FEATURES AND LIMITATIONS

     The ASTRAP model can be applied in emission policy assessments,  for  which
the  model  can be  exercised  to predict how  deposition fields would  change  in
response to  changes in the  pollutant  emission field  (Streets et al.,  1983).
This  assumes linearity between  emissions and  deposition.    In assessment  of
future actions, the specific meteorological  conditions and short-term emission
variations  that  will occur  are  not known,  and thus  an expected or  long-term
average  deposition  field,   such   as   provided   by   ASTRAP,  is   generally
appropriate.   The model can be  used  to estimate concentration and deposition
at  specific  receptor locations,  and can  estimate  the individual contribution
from different sources or source regions.

     The   following   major   simplifications   and   assumptions   have   been
incorporated into  the ASTRAP model:

    1.  Long-term  horizontal dispersion and long-term vertical dispersion can
        be simulated independently.

    2.  Long-term  horizontal diffusion can  be approximated  by  the  spread  of
        plume centerlines,  with  small-scale diffusion about  individual plumes
        ignored.

    3.  Chemical  transformation  can be parameterized  as  a linear  first-order
        process.

    4.  Wet  removal is  a  function of  the half-power  of  the  precipitation
        amount over 6 or 12 hours.

    5.  Transport  is two-dimensional.

    6.  The  dry  deposition  parameterization  is   horizontally  uniform   (not
        intrinsic  to  the ASTRAP approach  but a  part of the ASTRAP version
        described here).

     The separation of calculations of horizontal  and vertical dispersion  in
ASTRAP into  two-dimensional  horizontal  processes  and one-dimensional vertical
processes  is  much lower in  computational  requirements than  three-dimensional
Eulerian   numerical  modeling,  while  the  assumption   of   independence   of
horizontal and vertical dispersion processes over periods of  a month  or longer
appears  to  give  acceptable  results.    Because  of  the separation  of the
calculations   of   horizontal   and   vertical   dispersion,   however,   any
parameterization occurring solely in VERT has no horizontal variation, and any
parameterization  occurring  entirely in HORT is  independent  of  the  height  of
release.

     The  "statistical   trajectory"   in  ASTRAP  refers  to   the   treatment  of
horizontal  dispersion.   Rather  than  diffusion  of  individual plumes, the
distribution of the plume centerlines resulting from synoptic  scale variations
in  transport  wind  fields   is  calculated.    For  monthly  and  longer-term
dispersion  over   the  regional scale,  the  latter  process  is dominant.   For
episodic and shorter time  scales,  however,  ASTRAP  should  not be used without
modification to account for plume spreading.
                                     -7-

-------
                           4.   TECHNICAL DESCRIPTION

HORIZONTAL DISPERSION

     In  order  that  HORZ  need  be  exercised  only   once   for  a  particular
meteorological  period,   regardless   of  the  number   of  emission  scenarios
examined,  horizontal  dispersion  statistics  are  calculated  in  HORZ  for  a
virtual  source  grid covering  all of the  contiguous   states  and the  Canadian
provinces,  and  subsequently  interpolated  when  concentrations  and deposition
are calculated.  The grid is  of National Meteorological Center  (NMC)  spacing,
i.e., 381 km at 60° N, with 123 virtual sources as indicated in  Fig. 1.

     The meteorological data  are  time series of horizontal transport wind  and
precipitation fields, as described in the data requirements section.   The wind
fields  are  organized  in  quarterly  files,  while  the precipitation  data  are
organized in monthly files.

     Since  similar  calculations  are  made  for  each  virtual  source,  this
discussion  focuses  on a  single  virtual  source.   For a one-  to  three-month
sequence  of  meteorological  analyses,   simulated  tracers  of  unit mass  are
released from the source  at intervals of six hours.   Calculations in  HORZ  are
independent of  the height of  release.   The  tracers  are tracked  for  28 time
steps  (seven  days), unless  they  leave the  wind  grid.   At  each successive
tracer  position,  the  precipitation  field  is  checked  to  see  whether  there
should be wet removal.

     Positions are determined  for the midpoints of the six-hour segments,  and
thus the  first  trajectory segment  is  three hours  long,  while all subsequent
segments are  six hours long.    The  statistics  calculated in HORZ are  ensemble
statistics,   where  each   ensemble  represents  all   trajectory  positions
(endpoints) of a particular plume age (time since release) for each source  for
the length of the meteorological analysis.

     The  statistics  generated  for  a  puff  describing  the  density   of  the
ensemble of equal  age trajectory endpoints are the  coordinates y  and y   of
the mean position, the  standard deviations  o  and  
-------
concentration  and  deposition   are   calculated.     For  efficiency,   however,
extremely   low  contributions   outside   approximately  2.0  -   2.5   standard
deviations  of  plume spread are skipped.   The  contributions  at  the other points
are  scaled upward  by  1.046 to compensate  for the  mass  lying outside  the
threshold puff density.

     The wet  deposition parameterization is  based  upon the work  of  Hicks  and
Shannon (1979), but the coefficients  have been modified to  take  into  account a
different analysis  method for the precipitation fields.  The  parameterization
is
                      min
o
    0.5 P..1/2 M  .

    0.5 M        '
                                                      ±.
                                                           o.i
where R. . is the fractional pollutant deposition at  time  step  (plume age)  j  of
the  tracer  mass from  virtual  source i,  P^.  is the six-hour  precipitation  in
centimeters  encountered by  the tracer  during time step  j ,  and  M. .  is  the
initial  normalized  mass less the accumulated wet  removal before  time  step  j.
A  mass   equal  to the  wet  deposition  (or 1/2  of  the removal  from the  mixed
layer)  is shifted  to the free  troposphere,  where  it is  subject  to subsequent
wet  deposition processes but no  longer contributes  to  dry deposition nor  to
atmospheric   surface  concentrations.     Because   of   the   "min"   function,
precipitation  amounts greater than  1 cm have  the same effect as a 1-cm amount,
while because  of the  threshold,  amounts less  than  or equal  to  1 mm  have  no
removal  effect.

     As  previously  mentioned,   there are  separate  sets  of  statistics  for  wet
and  dry  tracer  ensembles.    The   same  trajectory  endpoints   are  used  in
calculation  for each  corresponding pair  of  wet  and dry  ensembles,  but  the
weights  are  different  (most of the wet ensemble  have  zero weight) , and  thus
the wet  and  dry ensemble puffs  differ.   The  ensemble puffs for dry deposition
and  surface  concentration  tend  to be  more  regular  in  the  trend of  their
overlapping  positions  than  do   the  wet  deposition  puffs.   This is  because  wet
removal  is   a  highly  irregular process  and  thus exhibits more  statistical
variation for  a single season.   The mathematical description representing  each
puff  j   from  source i  is  a   bivariant normal   density  function  (Mood  and
Graybill, 1963) :
 „,   ^
 H(x'y)
                I
              a   ,,  2
             x y /l-pxy
where x  and  y are any receptor  coordinates  and the terms y , y ,  a  ,  a  ,  and
P   are  functions of  source  and  plume age.  Examples of  sequences  of 28  puffs
are shown  in Figs.  2  and 3 for  dry and  wet  processes,  respectively, from  a
hypothetical  Oklahoma  source during  summer of  1980.   Weighting  of the mass
associated  with each  puff  is  not  shown  here;  older  puffs have less mass
because  of  wet and dry  removal  enroute.   Sometimes an  apparent  curl in  the
sequence of mean dry  puff centroids  during the last time steps is  produced as
a result of tracers beginning  to exit the boundaries of  the wind fields.   The
endpoint ensembles thus  become "slow"  subsets  of preceeding ensembles.   It is
clear that the dry puffs describe the mean dry  plume better than the  wet  puffs
                                     -9-

-------
           DRY  DEPOSITION  PUFFS    SUMMER  1980
             '-|N
              I   \
                     \

                     /•»
         ft

         x'l
        > J>
        /

       X

      x >

      /  -s
      /N

     / V
         ^/^
      I
                                       \A-
   x  '
   / / *• i   '
  * i  "^   '

N' '    ^<.
                               /
     >   /
\     >   /
\     l
»    »'x

 '
          /   >
         /   /
  /
  /

'>..
           /      7--*

         Y--J  '
/


/- ^
                                \      - ' "
                                 ^  .--* ;<•
                                 "    -' !'•

                                 /*^0.-  i  .
                                / *
                                 ' ~*
                                !'*.•*--    '-,
                               >, -.  V r -. , >

                              X* '  "  V--'
                              #  '  \ ,.  ,
                                                       ^S ' ,
                 \

                 \    t

                 \   t
                  N  |
    Figure  2. Time  series  of horizontal dry  puff densities from a hypothetical

            source in Oklahoma,  contoured  at  one standard deviation.  Puff

            centroids are indicated by (•).
                                 -10-

-------
          WET DEPOSITION PUFFS     SUMMER  1980
Figure 3. Time series of horizontal wet  puff densities from  a  hypothetical
         source  in  Oklahoma,  contoured  at one standard deviation.   Puff
         centroids are indicated by (o).
                               -11-

-------
describe the  wet  deposition "plume".   This is  a statistical feature arising
because dry  conditions  are much  more frequent.   The average of  a number  of
summers would produce a smoother wet plume.

     As in many atmospheric dispersion problems,  plume resolution  (here formed
by  overlapping  puffs)  is  most  critical near  the source.   The  basic model
calculations extend seven  days  in six-hour time  increments,  which means  that
one  of the  major computational  loops  in the  concentration  and deposition
program would normally  be exercised  28  times.    Postprocessing  statistical
trajectory  and vertical  integration  results  within  CONCDEP  before making
concentration and deposition calculations  combines the 28 six-hour steps  into
4 six-hour steps, 4 12-hour steps,  and 4 24-hour steps, extending  through the
same seven  days  of dispersion  but  reducing  computation time by  more than  a
factor of two.

VERTICAL DISPERSION

     The processes of turbulent vertical diffusion in the mixed  layer, leakage
to  the free  troposphere,  chemical  transformation,  and  dry  deposition   are
calculated in program VERT of  ASTRAP with a one-dimensional  version of  the
Gaussian Moment Conservation (QIC) technique (Shannon, 1979).  The  CMC version
here has nine layers plus a storage  layer representing the free  troposphere.

     The diurnal  variation or  "breathing" of  the planetary  boundary layer,
i.e.,  the   periodic   decoupling  and coupling  of  the  surface  and   elevated
portions of  the  mixing layer by  the cycle of  nocturnal inversion formation,
intensification,  lifting, and destruction by convective mixing from below,  has
been  characterized  in   field  investigations  (Hicks  et  al.,   1981)  and  is
parameterized  in  VERT  through  a  seasonally  and  diurnally  varying stability
profile (vertical eddy diffusivity K~ specified for each layer).   The  vertical
resolution is that of the layer thicknesses (0-100, 100-200, 200-300,  300-400,
400-600, 600-800,  800-1000, 1000-1400, and  1400-1800  m).  An example  of input
KZ  profiles  (here two-hourly  rather than  hourly as  in VERT)  is provided  in
Figure  4.     The  K   values   are   roughly  similar  to  the  eddy  viscosity
calculations of Yamada (1977), given the discrete  layered structure of ASTRAP,
which  acts  as a  filter  for some  turbulent eddies, and  given that the ASTRAP
parameterizations  represent  the   mix  of  sunny,  well-mixed  conditions   with
cloudy, more stable conditions.

     The variations   in  dry deposition  associated  with typical  diurnal  and
seasonal patterns  of  both  atmospheric stability  and  surface  resistances  are
parameterized through average diurnal  values  of dry deposition  velocities for
each season (Wesely and  Hicks,  1977; Wesely and Shannon, 1984;  Wesely et  al.,
1984).  The primary uncertainty  arises in extending the measurements  taken at
a few  high  quality,  well  characterized  research  sites  to  the  continent  as  a
whole.  The seasonal patterns currently used in ASTRAP for SOX and NOX species
are plotted in Figs. 5 through 8.

     A diurnal variation in the  linear first-order transformation  rate of the
primary pollutant  (S02  or NO/NO.,)  to a secondary pollutant (SO,   or NO.,"),
directly  or  indirectly  due   to  photochemical  activity  patterns,   is   also
included (Gillani et  al.,  1981).   The seasonal patterns, illustrated  in Figs.
9 and  10 for the SO  and NO  pollutant systems as  parameterized  in ASTRAP,  are
intended  to  include   all chemical  transformation  pathways  except  those
                                     -12-

-------
         SUMMER  VERTICRL  EDDY  DIFFUSIVITY
  8.
£8.
-C CD

CD


OJ
  CD
                                   10
                                                       n q
          300     600
900
 1200

TIME
1500    1800    2100    2400
                                                         2 — 1
  Figure 4.  Time-height cross section of vertical eddy diffusivity in m s   for


          summer conditions in ASTRAP.
                               -13-

-------
 SQ2  DEPOSITION  VELOCITIES
  oo
CJ
Q
   CD
    •
   o
   CM
    •
   O


   O
    •
   o
    0
8    12    16
 LOCfll TIME
20
24
Figure 5. Diurnal patterns of SCU dry deposition velocities of ASTRAP for
      winter (A), spring (+), summer (n), and autumn (o).
                    -14-

-------
 SP4  DEPOSITION  VELOCITIES
   oo
 CJ
   CD
    •
   O
   CM

   CD


   O
    •
   O
     0
8    12    16
 LOCRL TIME
20
Figure 6. Diurnal patterns of SO^ dry deposition velocities in ASTRAP for
      winter (A), spring (+),  summer (D), and autumn (o).
                    -15-

-------
NO/^02  DEPOSITION  VELOCITIES
    CO
  z:
  CJ
    CD

    O



    •^r
     •
    O
    o
     •
    o
      0
i
4
r
8
 r
12
 i
16
                 LOCRL TIME
20
 I
24
  Figure 7.  Diurnal patterns of NO/NCU dry deposition velocities in ASTRAP
        for winter (A), spring (+), summer (D), and autumn (o).
                      -16-

-------
 ND3
   o
   GO
 CJ
 Q
   CO
    •
   o
   o

   OJ
   O

   O
   o
DEPOSITION
           •e-e-
VELOCITIES
     0
        8     12    16
        LOCRL TIME
     20
24
Figure 8. Diurnal patterns of NO-, dry deposition velocities in ASTRAP for
      winter (A) spring and autumn (©) and summer (D).
                     -17-

-------
502-504  TRRNSFORMflTION  RRTE
5
    00-
    CD-I
    in-
  CL.
                8    12    16
                 LOCflL TIhE
 Figure 9. Diurnal patterns of the rate of transformation of S02 to SO,
       in ASTRAP for winter (A), spring (+), summer (Q), and autumn (o)
                     -18-

-------
MOX-
NITRRTE
         TRRNSFORMRTION  RRTES
       en -
       oo-
     z: C£>H
     r. i in
     O_
       o
         0
         8     12    16
         LOCRL TIME
                            20
24
    Figure 10.
Diurnal patterns of the rate of transformation of N0/N02 to NO., in
ASTRAP for winter (A), spring (+), summer (a), and autumn (o).
                        -19-

-------
associated with  precipitating clouds.   The rates of  SOj transformation,  for
example, are greater than normally found in clear-air experiments  (Meagher and
Bailey, 1983) because  the  rates  include the effects of chemical processing in
non-precipitating clouds (Hegg et al.,  1984.)

     Since most  emissions  from  minor  low-level sources  are  concentrated  in
urban areas, where polluted air  masses  typically are more reactive  chemically
(Alkezweeny and Powell,  1977)  than  in the less urbanized  surroundings  of  many
tall-stack  major  sources,  the  chemical  transformation rate   in VERT  is
increased during initial time steps for near-surface sources,

             ATxr = 0.5 B(3 - A)  (4 -  t)/(2 - 0.5A);  K3,  t<4       (3)
                        0, otherwise

where ATxr is  the  increase in transformation  rate  in  % hr  ,  B is  a constant
defined by  the  pollutant  system (1.0 for S02,  0.5  for  N02),  I is the model
layer,  t  is  the  plume  age in  hours,  and  A is  a seasonal  constant  (2  for
summer, 1 for spring and fall, 0 for winter).

     The dry deposition velocities and  low-level eddy diffusivities  during the
first  2.5  hours of  dispersion  are  not allowed to  decline below  arbitrarily
specified  seasonal  thresholds  (Table  1),  as  urban  areas are  less prone  to
strong nocturnal inversions than rural  areas  because of the local  heat island
and  aerodynamic roughness  effects  in  the  urban areas.   Likewise,  the  dry
deposition velocity  for SO^ during  the first 2.5 hours  of dispersion is  not
allowed to rise  above a specified seasonal maximum,  since urban  surfaces are
generally drier than  rural surfaces.   Elevated  releases  are less affected by
these  changes  than  are  low-level  releases,  because  the  diurnal stability
cycles  in layers  aloft  are  not  modified,  and the  pollutant  must   diffuse
downward before dry deposition can  take place.  The dry  deposition of  primary
particles is  set  to  1.0  cm s   during  the first 2.5  hours  to simulate  some
gravitational settling  of  larger particles from small sources, presumed to be
inefficiently controlled.

     Most  modelers  and  field  experimentalists  agree  that   there  must  be
"leakage"  from the  mixed  layer  to the  free  troposphere  through exchange
mechanisms such  as fair-weather cumuli,  but  no  definitive field  experiments
quantifying  such  vertical  transport  mechanisms  for pollutants  on  a regional
basis have been completed and analyzed.  A parameterization for  such a  loss is
now included  in ASTRAP, but the rates  have been set at relatively  low values
until such  time  as  more is learned  about the long-term  regional  significance
of  the  mechanism.    (Loss to  the  free troposphere  from  the  mixed layer  via
precipitation-related   processes  is  accomplished  in  subprogram HORZ  by
"depositing"  in  the  free  troposphere  1/2  of  the  bulk wet removal from  the
mixed  layer, where  it  is  subject  to subsequent  wet  removal  but  not  dry
removal.)

     Since time  evolution of  incremental  dry  deposition  and  of  the vertical
profiles of concentrations varies with  the  time  of release, simulated releases
are made  and diffused  for seven days  from each hour  of  the  day in turn and
then  averaged  separately for  emissions in each of  the first  six  layers,  and
stored in that  form.   Calculations  could be made also  for the  top three model
layers, but  normally  no  emissions  start  that high  (800 m).    In  an  earlier
version of ASTRAP (Shannon, 1981), a  diurnal emission variation  function was
                                     -20-

-------
              Table  1.  EARLY DISPERSION  PARAMETERIZATION LIMITS
Season
Summer
Autumn
Winter
Spring
Vd2 (max) '
0.50
0.40
0.40
0.40
7d2 (min)
0.30
0.15
0.15
0.20
Vd4 (min)
0.20
0.20
0.20
0.20
Kz (min)
3.0
3.0
3.0
3.0
7^2 and Vd4 are tne dry deposition velocities of  S02 and  SO^,  respectively,  in
cm s  .  KZ is the vertical eddy diffusity within the bottom  200  m in m2 s"1.


used  to weight  the  average,  but  experience  indicated  that  such  detailed
emission  information  was  almost  never  available  with  desired  accuracy  in
emission inventories or  scenarios,  and so that variation was  removed from the
model.    The   coding  here  includes  the parameterization,   but  there  is  no
variation in the rate.

     The  calculated  and  stored  normalized  statistics  from   the  vertical
dispersion  program  for  a  seasonal  simulation  include  the one-dimensional
surface concentration,  the  pollutant mass remaining  aloft,  and  the pollutant
mass deposited  by  dry processes  during the time increment  (six hours  except
the initial three-hour step).   Separate sets of  statistics are maintained for
S02,  primary  SO?"  and   secondary  SO?"  or  for  N0/N02,   primary  NOZ,   and
secondary NOo-

     The   eddy  diffusivity,   chemical  conversion,   and   dry   deposition
parameterizations  for  a  particular  season are  applied everywhere,  regardless
of  latitude  or  surface  vegetation.   This  is  an  obvious limitation of  the
ASTRAP algorithms, which include  some oversimplifications in order  to achieve
computational  simplicity.   A possible  improvement  would  be  to  calculate  the
vertical dispersion  program  with different  input parameterization  values  for
distinct source regions,  such as the  Deep South or  New England, with different
typical surface and meteorological conditions.  This still would  be  imprecise,
however, since the trajectories eventually pass out of the source  region.

CONCENTRATION AND DEPOSITION

     The program CONCDEP combines the statistics  produced by  programs HORZ and
VERT with an  emission  inventory to  calculate  primary  and secondary pollutant
surface atmospheric  concentrations  and  dry  and  wet deposition.   The version
described here produces those calculations for a  regularly spaced grid.

     Each emission source or  cell  is  nested within  the  coarser  trajectory
virtual source grid.   The emission  cells overlie the smaller cells  in Fig.  1,
although there are  no  man-made  emissions in  the  oceanic  cells nor  in some  of
the  continental  cells.   The bias  between the   emission location   and  the
                                     -21-

-------
appropriate  trajectory starting  point  is  added  to  all  ensemble  puff  mean
positions as  calculated from  that trajectory  virtual source.   The  standard
deviations of  each ensemble puff  are increased by  a term accounting  for the
emission and  deposition grid  resolution and,  for wet deposition  puffs,  the
precipitation  analysis grid resolution.   In  the  version  of  ASTRAP  provided
here,  all  of  the  grid increments  are 1/3  NMC,  and  the  assumed rectangular
distribution across a cell has a standard deviation  of  1/3 NMC unit//TI.

     A portion  of  the coding  in  CONCDEP is  designed to decrease computation
requirements by  ignoring  the effect of  puffs  outside their 2.5 a range.   The
mass associated with the puff is adjusted such that  there is no  loss  of mass.

     The  concentration and  deposition  fields are   calculated  by overlaying
puffs  (N source  cells  x 12  telescoped  time  steps) and adding their  densities
over  the  receptor  grid (illustrated  by the  smaller cells in Fig.  1).   Puff
superposition  avoids   computational  difficulties   often  encountered   with
continuous plume models in complex flows (Johnson,  1983) .  For  this  the  puffs
must be weighted by both  the emission rate per 6 hours for the  source and the
number of normalized tracer masses contributing to the  puff:

                       N   L       12      n
              C(x,y) = E   I   E   E       J-J     H(x.y)   c    ,   (4a)
                       i=l A=l  ^ j=l   o   3   i        J   J

                       N   L       12
              D(x,y) =» I   E   E   E   n   H(x,y) d    , and         (4b)
                       1=1 ji=i  1X> j=i  XJ         *J

                       N   L       12
              W(x,y) = E   I   E   E   m   H  (x,y) b   .            (4c)
                       1=1 £=i  ix- j=i  ^  w        x-3
where  C  is   the   long-term  average  concentration  of  either   gaseous   or
particulate  pollutant,  c.,.  is  the  corresponding  one-dimensional   surface
concentration of the  pollutrant  as  a function of emission  layer  and time step,
L is the number of  layers,  E^ is  the emission rate to  layer  H,  from source i,
n.f .   is  the  number  of  normalized  emission units  of  age j   from source  i
contributing to the statistic, n  is the number of simulated tracers released,
D  is  the  cumulative dry  deposition  of   pollutant,  d^,  is  the  deposition
increment of the corresponding pollutant during time step  j from layer  Si, W is
the  cumulative  wet  deposition  of  pollutant, rn^j is  the number of  normalized
emission units  from  source i  deposited during  time step  j,  HW  is   the  wet
analogue to  H, and b^,  is the  fraction of the original  pollutant mass  from
layer  i  still  airborne  at  time step  j if  there  had  been  no previous  wet
deposition.

     The NMC  polar stereographic  map projection defining the  various ASTRAP
input  and  output   fields  is  not   an  equal  area  projection.   In order  to
integrate the deposition field,  therefore,  the area represented  by  a grid cell
must be  computed.    This requires  a transformation of  the  NMC cell centroid
coordinates to  latitude  and then calculation  of  the  map  scale  factor  at  that
latitude.   The map scale  factor is then multiplied  by  the grid increment at
60°  N,  where the map  projection  is  true.
                                     -22-

-------
     Provided  that  deposition gradients are not  too  steep over the  region of
interest  relative  to  the  deposition grid  resolution,  the  simplest  way  to
define  an  irregular  area  of  deposition  integration  is  with  an  integer
indicator array,  which classifies  each  cell.   The  integrated deposition  can
then be  scaled by the ratio  of  the actual area  of  the region of  interest to
the sum  of  the areas  of  the cells  designated as being  within the  region of
interest.

PREPROCESSORS AND POSTPROCESSORS

     For ASTRAP use,  horizontal  locations  of emission sources must be  defined
in  NMC  coordinates,  or  in  units  convertible  to NMC coordinates,  such  as
latitude and  longitude  or UTM coordinates.   Effective plume height, or  stack
parameters   (height,   diameter,   exhaust   speed,   and  exhaust   temperature)
sufficient  to  calculate  the effective plume height  according to  a plume rise
formula such as that of Briggs (1971) must be  provided; ambient wind  speed  and
temperature can be estimated from climatology.  Since  information  on  fuel type
is often lost by the time an emission inventory is compiled,  the primary  SO^ ~
emission factor is  estimated  here  through  a parameterization that  varies with
location  (fuel  oil  combustion, most  common in the  Northeast,  emits a higher
fraction of sulfur  oxides  as  primary sulfate than  does  coal combustion)  and
layer  (surface sources, as  in residential  heating, are more  likely to  be oil-
fired).   If many  sources are  included  in  the  scenario,  it is   usually more
convenient  to  apportion  emissions  to a three-dimensional  grid,  expressed  in
polar  stereographic   coordinates  horizontally   and  expanding   coordinates
vertically.

     When  output  is  in  the  form of  gridded  concentration  and deposition
fields,  the normal mode  of graphic  production   is  to draw  isopleths  over  a
background  map of  the  region  of  interest.   Since  the  available graphics
packages vary considerably  from  one computer system to another, these  program
codes  are  not included here.   The output grid is defined in NMC  coordinates
which are (0,0) at the pole, X increasing to the  right, Y increasing  downward,
and oriented along 80° W longitude; conversions of NMC coordinates to latitude
and longitude, and latitude/longitude to NMC coordinates are  given below:

                 A - 6397/381
                 B = 3.1416/2
                 D = 57.29
                 C - A(l + sin(60,/D)l
               LAT - D(B - 2 tan~if(X/ -I- YZ)i/
-------
air  quality  networks  are  not  routinely operated,  although the  past  Sulfur
Regional  Experiment   (SURE)   (Environmental  Research  &  Technology,   1982)
provided a  valuable network  for much of  1977  and 1978.   Unfortunately, wet
deposition networks were  very  limited  then.   Current  observations  of wet
deposition measure  the  combined effect of  all emissions,  including  natural
emissions,   and  do  not  directly  evaluate  the  individual  source-receptor
relationships  actually  calculated  by  ASTRAP.    Observations  also   contain
"noise", in  terms  of  both  regional  unrepresentativeness  produced  by  local
source  effects  and   errors   associated  with  precipitation  undercatch  or
windblown dust contamination.   A model of  transport  and deposition represents
a relationship between emissions and deposition.  Model performance statistics
can  be  adversely affected  by errors  in emission  inventories,  as  well as by
errors or unrepresentativeness  of  observations.   Thus, even a "perfect"  model
of source-receptor relationships could never explain all observed variance.

     There are at  least  two approaches  to  "correcting"  observation data sets
that  contain apparent  local  effects  before  comparison  with model   output.
Detailed analysis of local sources near each observation site may allow one to
eliminate sites from the set subjectively, or the raw observations may  be used
to  produce an  objectively analyzed  field  on the  same  scale as  model  output,
with local effects and error presumably filtered to some extent.

     Verification   statistics   recommended   by  the   American  Meteorological
Society (Fox,  1981, 1984) apply  best to plume models,  i.e.,  near the  source
and  with  negligible travel  time for  removal  or  transformation  processes to
become  significant.   The  techniques offered  are  most useful  for individual
locations where  a  sensor is affected  by only  one  or a few  sources.    Also,  a
measurement is considered  both accurate and representative  of  an ensemble of
possible measurements for the particular conditions.   For long-range transport
(LRT) modeling the  verification techniques available  are not  so  direct.  Dry
deposition is clearly  important in the  sulfur  budget  and is likely important
in  pollutant  loading  to  ecological  systems,  but  it  is  not  monitored.
Horizontal  net  mass  fluxes  could  provide  an   important  check  for  model
consistency,  but are not measured.  Interpretation of monitored wet deposition
is  hampered  by  the  fact  that  various  monitoring  networks  have different
sampling and  analysis protocols.   Hence,  proper LRT model validation is  a
contentious and complicated matter.

     A  qualitative  evaluation of three models,  including ASTRAP,  developed as
part  of  the  Multistate  Atmospheric Power Production  Pollution  Study  (MAP3S)
program (Shannon  et al.,  1982)  gave a  subjective comparison  of  observed and
modeled  isopleth   patterns   of  atmospheric  concentrations,   and  presented
concurrent   simulations   of   wet    deposition   without   comparison  with
observations.  ASTRAP and seven other models were evaluated in studies  related
to  the  Memorandum  of Intent  (MOI)  on Transboundary  Mr Pollution  between the
United  States and  Canada  (Schiermeier and  Misra,  1982;  Ferguson and  Machta,
1982) .   However,   verification  statistics  were  inconclusive  because  large
discrepancies in  different emission inventories for  1978  were  never  resolved
and  because  the period  of comparison  was  rather limited—two months  in the
case of ASTRAP.   A comparison of ASTRAP simulations with selected high-quality
Canadian wet sulfate deposition observations was  given by Voldner  and  Shannon
(1983).
                                     -24-

-------
     More complete ASTKAP wet  deposition verification is given here.   Several
precipitation chemistry networks  in  the  United States and Canada provide  data
on  annual  wet deposition of  sulfate and  nitrate  at about  100  sites in  1980
(Niemann,  1982).   In this  comparison,  the contribution  from natural  sources
(not included in the emission inventory) is estimated in  terms of the  observed
rainfall at  the sites  by  subtracting  assumed concentrations  resulting  from
natural  sources from  the  observations.   Thus,  the  comparison  is  between
simulations  of  wet  deposition  of   anthropogenic  sulfate   and  nitrate   and
adjusted observations of the  same species.  Comparisons of ASTRAP  simulations
with observations are shown in  Fig.  11  for wet sulfate deposition  and in  Fig.
12  for wet  nitrate  deposition.   Statistics associated  with these  comparisons
are listed in  Table 2.   About  69% to 73%  of  the  observed variance of  sulfate
or  nitrate wet deposition  over the  continent is explained  by  the  model.    An
emission inventory for 1980, a  preliminary version of the NAPAP inventory,  has
been used.

     In  a  further  evaluation of  ASTRAP,  the  1980  SOX  emission inventory  has
also been used with meteorological data  from  1978 in recalculation of  S02  and
SOA ~  average  atmospheric  concentrations and  wet sulfate deposition   for
January  and  July   1978  and  for  the  entire  year   as  examined  in  the   MOI
studies.   The SURE January air quality observations were  for 10-31  January.
The appropriate 22 days of meteorological data were used, but that  is  probably
shorter  than  desirable for ASTRAP simulations.   Comparisons are  plotted  in
Figs. 13 through  15 and are summarized  in Table  3.   The arbitrary adjustment
for  the  effects of  natural sources at  wet  deposition receptors  (all in  the
east) is larger than that for the 1980 data; the 1978 observations  were mostly
from monthly  collection, whereas almost all of the  1980 were collected daily
or  weekly in  higher quality networks.   Model  performance is generally better
than  that  presented  for ASTRAP in  the MOI  model   comparisons,  even though
emissions are  for  a different  year.   The ASTRAP simulations  in  the MOI work
for annual  concentrations  and  deposition  were actually  an  average or sum  of
January and July results,  since those were the only months of meteorological
data  analysed  in model-usuable form for  ASTRAP  at the  time.   Here the  12
months of meteorological data  are December 1977  through  November of  1978,  so
there still is not a complete temporal overlap.
                                     -25-

-------
    WET  SULFRTE  RNNUflL  DEPOSITION   1980
  o-
  03
  o

,—, o
CJ
^ ^
in °.
\ o
rr •*•
O
en

CD
O
LJ
LJ
CO
 LJ
 E—
 CO
 CD °.
 tr o-
•  •

  %

'/  •'
   o

   o.
0.0
      10.0
                       20.0      30.0       40.0

                    hODELED (KG  S04/HECTRRE)
                                             50.0
                                                     60.0
  Figure 11. Comparison of ASTRAP simulations with observed annual wet sulfate

          deposition during 1980.  (Dashed lines are  for a perfect fit and

          for factor of two range about a perfect fit.)
                             -26-

-------
   WET  NITRRTE  flNNURL  DEPOSITION   1980
                  10.0     15.0     20.0     25.0     30.0     35.0
                  MODELED  (KG  N03/HECTRRE)
Figure 12.  Comparison of ASTRAP simulations with observed annual wet nitrate
         deposition during 1980.  (Dashed lines are for a perfect fit and
         for factor of two range about a perfect fit.)
                           -27-

-------
502  CONCENTRflTION   JflNURRY  RND  JULY   1978
  o-
  en
<-* o
E—
LJ
z:
'"""^ * ~r "
—' CD
O
LJ
CL.

cn °.
21 to-
ol "*•

CD
O

CJ
i—« o
2Z o.
•—• ro

Q
LJ
LJ o
CD  •
CD L0'
O ""
   a
   o
                                /
                                /
                               /
                              /
                              /
                /
                /
               /
                         /
                         /
              /
             /
             /
                      /
                      I
                      I
                     I
       •  %
         i
         i
        i
       i
     i
                 •
+  •     •+

•      /      4-
/•+ •
     0.0      15.0       30.0      45.0      60.0      75.0      90.0
             MODELED  (MICROGRRMS  PER CUBIC METER)
    Figure 13.  Comparison of ASTRAP simulations with observed average S02
             concentrations for  January (•) and July (+) 1978.  (Dashed lines
             are for a perfect fit and for factor of two range about a  perfect
             fit.)
                               -28-

-------
504  CONCENTRRTION  JflNURRY  RND  JULY  1978
  o
  rvj
25 «.
E—• cjD'
QJ —
O

CD
ID
LJ ~
Q_

CO
21
(H
01
CD
O o
a
LJ
>•
ce o
LJ ,,:•
en
CD
a
  o
  o
                                /
                                /
                               /
                             /
                             /
                         .+
     /
    /
    I
   I
   I
   I

   +
      I  f
      I '
                                    +
 /
/
/               +
           /  +
     0.0        4.0         3.0         12.0        16.0        20.0
             MODELED  (MICROGRflMS  PER CUBIC METER)
    Figure 14. Comparison of ASTRAP simulations with observed average SO,
            concentrations for January (•) and July (+) 1978.  (Dashed lines
            are for a perfect fit and for factor of two range about a perfect
            fit.)
                               -29-

-------
  SULFRTE  DEPOSITION  JflNURRY  RND  JULY  197E
  o
   •
  CO'
LJ
C£
cn o

«-" co-
CJ
LJ
x:
O
cn

CD
O
LJ
LJ
cn
CD
O

a
LJ
E— o
a
az
  a
   •
  a
£
     0.0
             2.0            4.0            6.0

               MODELED  (KG  504/HECTflRE)
8.0
  Figure 15.  Comparison of ASTRAP simulations with observed wet sulfate

           deposition for January (o) and July (+)  1978.  (Dashed lines are

           for a perfect fit  and  for factor of two range about a perfect
           fit.)
                             -30-

-------
       TABLE 2.  ASTRAP ANNUAL WET DEPOSITION VERIFICATION  STATISTICS.4

1980

1980
1980
1980
Caseb
sulfate

sulfate
nitrate
nitrate
#0bs .
95

95
95
95
OF
20.2

17.7
11.8
11.0
Mod
17.6

17.6
10.5
10.5
Res
2.5

0.1
1.2
0.4
S.D.
12.

12.
7.
7.
o
4

9
5
7
S.D
15

15
8
8
'm
.3

.3
.4
.4
S.D.r
8.3

8.0
4.3
4.2
O
om
0.83

0.84
0.85
0.86
aWet deposition units are kilograms/hectare.

 An asterisk  (*)  indicates that  observations have been  reduced by  estimated
natural source  contribution (SO^: 2 micromoles/liter  east of the Mississippi
or Manitoba and 7 micromoles/liter west; NOg:  1  micromole/liter east and  3.5
micromoles/liter west).
              TABLE 3.   ASTRAP VERIFICATION STATISTICS FOR 1978.a
Caseb
Jan S02
Jul S02
Annual
Jan SO,
July SO
Annual
Jan wet

Jul wet


so2

,
\
sulf ate^
•ft
sulfate
#0bs
17
20
3
29
47
9
7

13
Ob
40.7
27.3
30.3
6.8
11.5
9.1
1.5

3.1

Mod
35.5
21.8
32.8
6.7
10.0
9.2
1.1

2.1


Res
5
5
-2
0
1
-0
0

0
.2
.5
.4
.1
.5
.1
.4

.9
S.D
17
13
15
1
2
1
2

1
•o
.7
.9
.9
.7
.8
.9
.8

.8
«•»•.
16.5
13.5
13.2
1.6
3.4
1.5
0.9

1.7
S.D.r
17.6
13.5
5.0
2.3
2.5
1.9
2.2

0.8
Pom
0.44
0.49
0.64
0.06
0.69
0.35
0.67

0.82
Annual wet sulfate      8   28.0   14.7   13.2      15.6   13.1    6.8    0.79


aUnits are micrograms/cubic meter for air concentrations and kilograms/hectare
for wet deposition.

 An asterisk  (*)  indicates that  observations  have been  reduced by estimated
natural source contribution (5 micromoles/liter).
                                     -31-

-------
                              5.  COMPUTER ASPECTS

     As currently  structured,  ASTRAP consists of  three  programs, HORZ,  VERT,
and CONCDEP, the first  two  of  which must be run to create data files  required
by CONCDEP.   There is no prescribed order  in  running HORZ and VERT,  and  they
could be run simultaneously.   The  job  control language is appropriate for the
IBM 3033  mainframe  at  Argonne  National Laboratory,  and would  require  some
modifications  for  use  on other systems.    While  programming is in  standard
FORTRAN, input and output algorithms may require coding modification for other
computer systems.   A simplied  flow chart  for  the overall model is  shown  in
Fig.  16.    Obviously  there is   preprocessing of  raw wind  and   precipitation
observations in  order to produce  the  gridded fields, but  as the fields  were
initially processed elsewhere as discussed in the Data Requirements section no
discussion of their preprocessing is provided here.

     On an  IBM 3033, VERT  requires 210K (bytes) storage and 5-10 minutes  of
CPU time, plus one output file, for a seasonal simulation.  HORZ  requires  500K
storage and 10 minutes CPU time, two input tape drives, and an output  file for
a seasonal simulation.  CONCDEP requires 500K storage, 10-15  minutes CPU  time,
an  emission input  file,  the  output files  from  the  other  two subprograms,  a
file  used  to  identify  each  receptor  cell  as  to  state  or  province, and  an
output  file for a seasonal  simulation.   The  meteorological fields  on  tape
input  to HORZ  are written  in  binary,  as are  the output  files  from  HORZ and
VERT,  the emission file, and the output  file from CONCDEP.

     Quarterly  sets  of  wind  analyses   and  monthly sets   of   precipitation
analyses are written  on two separate tapes.   The  internal file number of the
nonlabeled  wind  file  (file  10  in  HORZ),   and  the  names and  internal  file
numbers of  three corresponding  precipitation data  sets (files 20-22 in HORZ),
as  well as the  tape number  for  different  years  must be  changed whenever  a
different  season is  simulated.  The names  of data sets  created  by the  three
programs should be changed every time a  different period  is simulated.

     In normal application ASTRAP is structured so as to  read an  emission  file
from  disk  storage.   Since  the emission inventory  is likely  to  be the  first
data set that a user would change,  however,  the appropriate code  is  "commented
out"  for  the  test example  and an  algorithm defining a single  hypothetical
source  is  added.   The write file is unit 9 for  each program; the input  files
in  CONCDEP must be  renamed when  different meteorological  periods  have  been
used in HORZ or VERT.

     The   output   file   from  CONCDEP   is   normally  stored  on  disk,   and
postprocessors  are later  used  to  display  the  results  graphically.    Since
available  graphics libraries will  vary  greatly  from one  computer  system  to
another and  are  often  partly  proprietary,  the postprocessor is  not  included
here.    Input  parameters for the ASTRAP  programs  are listed  in Table  4,  while
input/output  files  are  summarized  in  Table  5.    Some   changes   which  a
prospective user of  ASTRAP might wish to make  are discussed in  the following
section on the test application of ASTRAP.
                                      -32-

-------
              at
           Z CC
           O 01
           — x
           H 0.
           < CO
    Z     £ O
    O     as a.
    —     o o
    40     U. CC I—
    a     co H- iti
    u.     z    o
i—o.     < LU a
as — _i  oe uj 3
UJ Q <  H- CC 00
>•    >     U.
    _l O  _1    V)
    < £  < O co
    O UJ  O I- <
    -, a;  ~    £
    > Q
          uj coca
              o i
                                 to
                                 z
                                 o
                                UJ
                                CO
                                o
                                OS
                                O
                                 01
                                 _l
                                 M
                                 o.

                                 H
                                 3
                                 a.


                                 o
    CC
    Ol
    0.
    CO
OS H O
O Z £
1C O OJ
    N CC

    5 H
    O 01
    Z X
               /S
   V)
Q 01

Ol CC

U. CO


Q 01
Z £
 ~1 CO
 Ol Ol

 O. CC
    01
 O. CO

 (J Ol
 01 £
 CC —
a. H
co
z
o
CO

01
V)

u

u.

u
01
CL
CO

CC
o
u.

CO
Ol
                                0.

                                3
                                o
                                    v)  o
                                    z  —
                                    O  H-
                                    —  <
                                    H  as
                                    <  H-
                                    ce.  z
                                    H  UJ
                                    Z  U
                                    UJ  Z
                                    U  O
                                    z  u
                                    o
                                    u  u

                                    U  CC
                                    —  UJ
                                    OS  Z
                                CX. UJ  Q.
                                LU Z  (/)
                                a Q-  o
                                <-> v>  C
                                Z O  H-
                                O £  <
                                o t-
                                    <  I-

                                    I-  <
                                    Z  H
                                    <  3
                                    H-  _l Z
                                    3  -I O
                                    -I  O —
                                    _l  Q. h-
                                    O     —
                                    a.  > co
                                        ce o
                                    >-  < a.
                                    cc  a uj
                                    <  z a
                                    £  o
                                    «  o >
                                    CC  UJ CC
                                    a.  co Q
                                                              K V)
                                                              UJ >
                                                              a <
                                                           z Q -i
                                                           o 3 a.
                                                           — B5 CO
                                                           H-    «-
                                                           — Z Q
                                                           V) O
                                                           O — Q£
                                                           OL I- <
                                                           LLJ — _J
                                                           Q CO 3
                                                              O OQ
                                                           t- a. <
                                                           UJ UJ |-
                                                                              UJ
                                                                              _i

                                                                              u.

                                                                              H-
                                                                              3
                                                                              Q.
                                                                              H-
                                                                              3
                                                                              O
                                                                                           CO
                                                                                           UJ
                                                                                           (A
                                                                                           O
                                                                                          QL.
                                                                                                                        UJ
                                                                                                                        o:
                                                                                                                        a
                                                                                                                        3
                                                                                                                        CC
                                                                                                                        I—
                                                                                                                        CO
                                                                                                                        CC
                                                                                                                        CD
                                                                                                                        O
                                                                                                                        CC
                                                                                                                        a.
                                                                                                                       (-0
                                                                                                                        LU
                                                                                                                        CC
                                                 o
                                                 _1
                                                 UJ
                                                 l«4
                                                 u.
                                                 CO
                                                 CO
                                                 UJ
                                                         -33-

-------























Oi
Cd
H
i
2
INPUT
3
rn
^2

^^

Cd
rJ
33
H























Cd
SI
rH
>•
<(J
ftj
Sj


H
2
M
O
Cu



H
Z
r-3

O
H



5
o
oS
PH

CO
H
M
Z











Z
0

H
rH
M
«!
O
C/5
S





Cd

*5
Z



O
vO
»_^






O
•
vO
px
o




en

5
o


a-
Cd
S
O
CJ

•iH
B
cr
CO

o
^H








cu
4J
CO
4-1
CO
JS
o
CO
CU

O

CO
cu
^4


cd
3
<6



















































1







CO
M
>
>^ 4.J
CO tH
rH CO
<4H d
J2 0 CU
o -a
co d
cu o u
•H CU
>4H 4J >,
O O Cfl
Cfl rH
4-1 I4H 4J
iH CO
CO Cfl CU
d 3
cu co o
Q CO rH
Oi
O
CJ
Q

^^
vO
m
oo
CM

T3
CU
iJ
J-l
S
rJ
O
M-^
d
3

O
CM

4J
d
3


Ot
Cd
Q
CJ
Z
o
CJ
•a
cu
N
•H
rH
s
(^
0
z



d
>% 0

4-1 « Cfl CO
d Q« T3 co
cu cu d IH
en s 4J OB
cu co « o cu
d u co co
O CJ CU cfl CO T3
•H C S bO d
4J TH iH ^^ T3 CO
•H u a
CO CO 4J Cfl /—^
o co A d cu
a. co o cfl >M-H
CU B Cfl 4J P O
"O CU 3 Cfl -H M
Z rH B 4J CU
>, U rH iH iJ X
M M O O !-i Cfl CO
Q o <*H (X CX CU rH
H

W
Q








CU
4J
4J
§
fH
O
<4H
d
3

O
{**)

AJ
*r4
d
3


a,
Cd
a
CJ
Z
o
CJ



CJ
^





>4H
O

4J
d
cu

o
d
•H

rH
Cfl
4J
d
o

tH
l-l
O
3S


^
Q








































X

00
d
o
rH
cfl

CO
rH
rH
CU
o

d
o
tH
CO
co
•H
B
cu







CO
^^^ /



CU
4J
4J
a o
M •
o m
(4-4 pr |
£3 



CJ
"





<4H P
O >4 CU

4.) t>0 cfl
d d ^H
0) O
S rH -C
CU CO U
M CO
u co cu
d rH
•H rH UH
cu o
rH O
CO CO
4J d CO
do cu
o -H d
N CO -*
•H CO O
i-l iH -H
OS £
EC CU H


>< N
O Q


^
O
CM
r^


cu
4J
4J
a
M
O

C3
3

0
m

4-1
•iH
d



Cu
Cd
Q
CJ
z
0
CJ


f_l
CU
a.
H









cu
o
3
O
CO

^

CO
d
o
•H
CO
CO

B
Cd


2
Cd































d
o
CO
cfl
CU
CO











IJ
cu
CO
rH

•o
d
CO

rH
f^
cu
O







CM
t-H
^^






o

CM
Ct«
CN
rH



cn

5
o


a,
Q
CJ
Z
0
CJ



1








cu
.d a,
4J CU
4J
d co
•H
CU
iO S
d -H
tH 4J
4J
J3 T>
bo cu
•H a,
cu o
S U
CO
4-1 CU
4-1 rH
3 CU
ft,) iJ


3
Cd
                                           •o

                                           I
                                            o
                                            u
-34-

-------

















T)
cu
3
C
uj
T"I
4-1
c
o
o
-3-
Pd
J
05



























Pd
h-l
CO

>J
S2
H
5
3
o
PK


H
§
O

5
o
o





CO
H
Z









Z
o
M
H
P-i
M
O
CO
Pd
O





Pd

.,
T3

rH
CO
1-(
U
•H
c
M


Q
a
M

m
^~
^
in




CN
M
CO

c/>
o
0
P*
a
o
cj





1




J3
CO 4J
4J CU
C r-<
cu d«
S 0
cu co
r-l -H
CU
M
"0 0

a M
O 60 >%
T^ tH
4-1 CO CO 60
to cu ca a
r-4 S CO 1-1
3 -H cu S
S 4-1 O CO
•H 0) CU (-1
CO Q C T3


Q
z
M

In
'^
_^
"-'




CO
t-l
^

co
§
"
a,
Pd
Q
O
§
O





1








CU
o
a
•H
^
o
^4
a.
M
"O CU
C -H
CO 4-1
•H
CU 4-1
•u C
CO CU
4-1 -O
CO -H

Q
Q
Z
M











CO CO
M M

W CO
§ §
0 CJ
Ot
Pd
g
PS 21
Pd O
> 0





1 1



*
CO

II

i-H f-^ ^
II II S
fl
X $-1 3
o cu co
CO 4J
N_X C3 M
•H CN
4-1 3 ^
d /-v — ' n -d-
CO CN
4-> e so n
3 II 0 C
>H CO -H r-(
^H X CO U f-H
o o cu a, co
Q-i Z CO CO Pn

J
o
pj CO


^ 2
S
^D O^
v-x **
CN


r^
•
-* m
< to
00 CN
— H r-H

en co
§ §
u o
PM
Q
U H
z orf
O Pd
O >





CO
1 CN^
a





^ i
O 3 T3
cO cu «-i C
CU O 4-1 CQ
C -H

o > cu
O >, 5s

O tt T> iH
•H CU
•U -O >,
tO C i-l ^3
•H to CO l-i
> O >», 3
CU CU -H 4J O
!-l 4-1 4J -H JS
jQ CO l-i >
43 4J CU -H >.
<3 to > co .a

CU
H

M s2

S
^^
ON"
^




CO
M
2

CO
9
o

0
33





1

CU
O
M
3
O
CO

(— 4
CO
^j
S-l
•H
>

^
O
o
cu
1-1
CO
H
Pd
Q

<£
t-}
































G
0
T-l
4-1
CO
O

4-1
•H
C
cu
•a
•H

-o
t-l
M
bO



































J[
O
M

XI

5p
O
M

SO
C!
•H
cu
9)
CJ
O
U
0.








































4-1
4-1
CU


U
cu
3
0
rH

3
0
M






^
CN
"*
>-^


CD

CM
PE^
CN
t— 4

CO
9
o
Q-,
Pd
Q
CJ
O
CJ





I

1
CJ
1)
i—) O
CO 4-1
M a.
4-1 CIO CU
C3 4-*
rH -H CO
CO 4-1
•H 3 CU

a T-I -H
a> 1-1 4J
4J 4J
o c -a
a. o cu
o a
B o
3 CO 0
B 

z

^
S
                                           •a
                                           cu

                                           s
                                           o
                                           o
-35-

-------


























x"s
-o
0)
C
•H
C
o
o
• 33 O





1 1 1



CO
I-l
cu
^•t
« cd
CO i-H
A
-H u e
0
II M M -H
CU 3 CO
U U O CO
0) fl -C iH
S -H 1 8
a if ^o cu
3 ^r co
CO ~^ CN II O 0) O
J-)
C II bO ^ co M
o c cu cu
CO rH iH 43 01 43
co rH \-< a a a
CU CO CU 3 -H 3
CO fe CO Z 4J Z
Z
o
CO

Cd CO
Cft H t>3
Z Z Z

o\

M
r^
"-'
T3
CU
U
u
(rt
g
^4
O
4H
C
3

CM
£•4
M

4-1
•H
C
3

Oi
§





i








e
o

4-1
cfl
4-1
•H
Q.

O
CU
M
Cu

^4
3
O

1

Oi
M
CJ

2
0,



















x-s
CN
CN
1
0
"






























CO
rH
m

o
4-1







^
vO

cn
oo
Ci
TJ
cu
4-1
4-1


M
O
4H
e
3


0
CN

4-1
1-1
C
3
PU
CONCDE



T3
CU
N
•H
i
u
0
C







cu
e
^4
O M

M 4H
•H
CO CO
CO
csO cO
C S
iH
C Z
•H
CO M
e o
cu
a! co



Q
3
35
CO














































M
cu -d !-i
0) C 0)
4-1 Cfl >,
ca cd
" rH
0) 4J
ace
iH CO O
4J 4-1 iH
3 co
42 rH CO
CJ rH -H
co o a
01 Cu 01







^
vO

{*•>
co
CN
•O
cu
4J
4J
cd
a
S-l
o
4H
C
3


O
CN

4\J
•H
°

CONCDE



t3
cu
N
iH
rH
i-l
O
C








1
C
0
CJ

O
1
r-H

rH
cd
o
iH
4J
u
01




X!
o
CO










































CO «
cu cu
rH O)
•rH 4-1
4H CO
O T3 M
M 01 C CU
cu a « >>
iH CO
C 4J - rH
O 4-1
•H 43 C 0
4J CJ Cd O
CO CO JJ iH
t-t CU 3 Cfl
4-1 rH Cfl
C H rH -rj
CU O O S
CJ 4H CU CU







                                             -d
                                             cu
                                              O
                                              CJ
-36-

-------























•o
3
C
4_l
C
o
u
Ed
fO
<
























Ed
CSJ
M
CO
>*
^3
cm


H

o




H
Z


o
M

S
CJ
I




CO
H
M
Z
3





Z
0
1— 1
H
Pi
M
Crf
CJ
CO

Q





Ed

^
Z


CTi
^? rH
CN .J"
^ £
\-^
T3
0)

4-1
M
 33




CO
** ~gf
cu bo *
XI C W
4J O CO
1 rH t8
cO <4-) cO CU ^-•v
S CN 0 v-' 1 X!
MO 4J 4-1
O CO CO T3 CO 3
Un4 CO W C QJ O
en '4-1 O C -H 3 tn
C 0 Z 0) 5 1
co e co x!
MO) » O 4-1 'H 4-1
4-1 4-1 X O. M M
CO O S O X O
>, M z o a c
rH O CO -
M C - C 3 CO
3 O -* >< CO o -H
O i-l O -MO
03 4-1 CO X 4-1 CO >-i


^
erf
H 3











O
•
ul
ta



CO
Q
2§
CJ


0
03




1
CO
^
O 01
a M
CU 01
M J2
a.
4J CO
0) O C
s o. o
0 -H
rH M 4J
CO 4-1 U
a oi
O 0) >
•H cu e
4-1 M O
O M-l CJ
CO
M 0 X
en
o
erf
H
£3

^•^
CN
CN
A
CN




"*
-^*
ft,
CN
"



CO
a
3



Ed





CO

• M CD
^t cu S
4-1 S 1-f
•H 3 4-1
O rH
O Cu «
rH 4J
0) 4-1 C
> O CO
4-1
C3 rH 3
O CO iH
•H e rH
4J O O
•H 1-1 a >,
CO 4-1 CO
O O •> T3
CU fi 0)
CU 3 bflu-i
Q 14-1 CO O



Q
t>
,-y
CO
CN
r-H
CO
CN
5
T3
CU
4-1
4J
M
O
4-1
a
3

O


4-1
1-1
C
13
cu
Ed
a
o
o
CJ




1


M-! CO
>4-l U
3 f-l
Q. 4-1
CO
C f-l
O 4J
•H CO
4J 4-1
T-t CO
CO
O 0)
CX rH
0) XI
T3 S
0)
4-1 CO
1 g



Q
3


3
CN
y^

T3
CU
4-1
4-1
g
O
*4-(
c
a

o
c*O

4-1
•H
C3
3
Oi
Ed
a
CJ
z
0
CJ




CJ





"4-1
o
CO
CO rH
CU rH
4-1 CU
cfl O
C
i-l C
-0 O
rl fl
O CO
O CO
O i-l
X cu



CO
X


s
CN
^^


0)
4J
4-1
g
O
(4^
C3


O
CO

4-1
•H
a

Oi
Ed
P
CJ
Z
o
CJ




CJ





<4_l
o
CO
CO rH
0) rH
4-1 a)
CO U
c
•H C
T3 0
(-1 1-1
O CO
O CO
O iH
S
>< CU



CO
>•*
-37-

-------
                     TABLE 5.  ASTRAP INPUT/OUTPUT FILES
            INPUT
 UNIT #     OUTPUT  PROGRAM
                    PARAMETER NAMES
                                     COMMENTS
            OUTPUT   HORZ
            OUTPUT   CONCDEP
10
10
INPUT    HORZ
INPUT    CONCDEP
11
OUTPUT   VERT
20
INPUT
ITP         INPUT
(20,21,22)
CONCDEP
         HORZ
30
50
INPUT
                   CONC(6,28,123),
                   WD(6,28,123),
                   LANDE(19,16),  NF,
                   NXEM,  NYEM,  XEMLL,
                   YEMLL, DXEM, DYEM,
                   NTS

                   S02(51,45),
                   S04(51,45),
                   DDEP(51,45),
                   WDEP(51,45),
                   SUMM,  ADRY(60),
                   AWET(60)
                  ,  V(17,19)
          C0(6,28,123),
          WD(6,28,123),
          LANDE(19,16),
          NF

          AVSFC(28,3)
          AVDD(28,3)
          AVBGT(28,3)
SOX(28,3,6),
DEPT(28,3,6),
SBUD(28,3,6)

PRECIP(50,45)
CONCDEP   EM(720,6),
          XS(720) YS(720),
          NS, NZ, DX, DY
INPUT    CONCDEP   INDD(51,45)
                               Output from HORZ read into
                               CONCDEP containing season-
                               al simulation of horizon-
                               tal statistics
Final output file from
CONCDEP containing all
concentration and
deposition data, usually
stored on disk to be used
later by postprocessors

Nonlabeled seasonal 6-hour
wind file, requires input
tape drive.

Results of HORZ, contain-
ing seasonal simulation of
horizontal trajectory
statistics

Output file from VERT
read into CONCDEP
containing 1-D surface
concentration data, dry
deposition increment and
total airborne mass for
each emission layer

Results of VERT for all
emission layers
                               Three monthly 6-hour
                               precipitation files,
                               require input tape drive
                     Emissions in kT S0x as SO™
                     by source cell and layer
                     for each season
                               Input file which contains
                               state and province
                               geographical information
                                     -38-

-------
                             6.  TEST APPLICATION

     The  computer  codes attached  as appendices  to this  report simulate  the
ASTRAP exercise  described  here.   In order to make the test example  simple,  a
section  of  code defining  a  single source  is  inserted.    In  more  typical
simulations  a preprocessed  emission inventory, gridded  on a 1/3 NMC  scale,
provides  seasonal  SO-  emissions  for  layered  virtual sources.    This  more
relevant  code  remains   in   the program,  but  is  skipped  by  a  GO  TO  4444
statement.   In normal ASTRAP  simulations certain  randomizations  are done  to
avoid biases  associated with discrete  source and deposition grids.   These are
temporarily  commented out  for this test example  (see statements containing
RANF  in  CONCDEP).    An  emission  of   1000  kT  SOX as SO,  per year  from  a
hypothetical  source in Oklahoma (same  location as source of ensemble  puffs in
Figs. 2 and 3) is combined with summer  vertical  dispersion  statistics  and with
trajectory statistics for  summer (June through  August) of 1980 in CONCDEP to
produce  the  concentration and  deposition fields printed  in  Tables 6  through
9.   The  printed fields  begin in  the  llth small cell  column  of Fig. 1.   The
integrated deposition amounts  for each state  and  province are  shown  in  Table
10.   (Deposition to  the Great  Lakes  is  totaled  without  regard  to  state  or
province.)

     Obviously,  a  user of  ASTRAP would  wish to make  changes from this  test
case of model coding.  The most likely  changes would be in  either  the  emission
field, the concentration  and deposition  field,  or  both,  in CONCDEP.  The key
to defining  successfully  a  different emission  field  is to describe the  field
horizontally  in  terms  of  NMC  coordinates and  NMC spacing,  if gridded,  and
vertically   according   to   the   layers   in   VERT.      A   change   in   the
concentration/deposition  grid  can  be  easily   accomplished  by  changing  the
variables describing the dimension, spacing, and lower  left corner.  To  change
to an  irregular  set  of  n  receptors,   however, the  two-dimensional  receptor
grids would need to be changed to  (n x  1)  fields, and the X and  Y arrays  would
need  to  be   described   in   NMC coordinates  and  dimensioned  to  (n).     The
integration of the deposition field should then  be skipped.

     For   the  more  sophisticated  user   who  might  wish   to   change   a
parameterization, the default values included  in the  code as input data  and
illustrated in Figs.  4  through 10  could  be replaced,  as well as some default
values (such as minimum deposition velocities during initial dispersion)  which
are part of  the  program  code.   A change  in the  input meteorological  fields to
different  spatial  scales   can  be  fairly  straightforward,  as  long  as  the
dimension,  spacing,  and  corner  of the various  arrays  are  changed  in  a
consistent  fashion.   A  change  in the   temporal  scale  would  involve  more
complicated changes.

     This user's  guide  is incomplete without  the  accompanying  program  codes,
which have been  augmented  by numerous comment statements.   It is suggested
that any potential user read the code in  detail.
                                     -39-

-------
    a« aa a> a> a< « a a>« a» e» a a a> a«(ea to s « o> s  a> «  a a a  a> tas « «a « e* « <*
o
~d
o:«
UJ —
                                                                40

-------
   ««<• «>« ««»«» a> s IB «« tg « 

                                                                      41

-------
  Bjaaaaaaasjaa«aaa
  in	
  o>aa"a«sj»aaaia«iaa«siaa«aasiaa

  *
  ujaaaaaasiaaaaaaaa««aasiaaaaa«iaaaaaaasiaaa«iaa«aaaa

                                        asiaa

        	siaaaa
  •»Bia<»si«aasiaa5siaaaaaaaaaaaaaaaaaaa a a a a s a

  CM
                              sj a a a a a a a a a a — — — ^j <^j — — a a a a a a oj a
                                                  aaa

                                                  — — • (M
                                esaaasiaaaaa — Nf>-r-»t-j — aaaaoasja
                                                 J — CVJ en Ul r- at — «» «J
omaaaaaaaaaaaaaaaaaaaaasjaaaaaaaasaaasioaaaaaaaaa
 .

                                  aaaaaaaaaa — (M-»oi^
z — aaaaa»aaasiaaaaaaaaaaa«)aa=>aaaaa —
OM
—                                 a a a aa a a a a a — — cn CT< •» a a a« a a a

MM ...........
u>                                 a a a a a aa ani a a a — a a «j a s a a as
oajQaaa«aaaa«jsjaaaaaaaaaas)aaaaaaaaaaaaaaaaaiaaaaj

                                                      59 a a «
 M                                       aj ttj eg si ta a ca a sj a *a a a sj GJ a «> ta
 h-inaaaaaaaaaainaaaaaaaaaaiaaaaaciaaaaaaasaaaaaaaoaa
 =>                                           »a a a a a a a a a a a a a a
 uma«>aaaaaaaaaaaaaaata

-------
o

o
o

z
o
o
cu
Ul
o
                                     J 01 eo ui — 5» 
 z

 >
 o
                          — 
-------
REFERENCES:

Alkezweeny, A.  S.  and D. C.  Powell,  1977:  Estimation  of  transformation rate
    of S0~  and  SO,  from atmospheric  concentration data.   Atmos.  Environ.  11,
    179-182.

Bolin,  B.  and   C.  Persson,  1975:    Regional dispersion  and  deposition  of
    atmospheric  pollutants  with particular  application to  sulfur  pollution
    over western Europe.  Tellus  24,  281-310.

Briggs, G.  A.,  1971:   Some  recent analyses of  plume  rise observations.  Proc.
    2nd International Clean Air  Congress (Edited by  H. M. England  and  W.  T.
    Baery).  Academic Press, New  York.

Durst, C. S., A.  F.  Crossley  and N. E. Davies,  1959:   Horizontal  diffusion in
    the atmosphere as determined by geostrophic trajectories.   J.  Fluid  Mech.
    _6, 401-422

Ferguson,  H.  L.  and L.  Machta,  1982:   Atmospheric Sciences and Analysis Work
    Group  2  Final  Report.    Submitted  to  Coordinating   Committee  of  the
    Memorandum of Intent on Transboundary Air  Pollution between Canada and the
    United  States.

Environmental Research  & Technology,  1982:   The Sulfate Regional  Experiment:
    Data  Base  Inventory and  Summary  of  Major  Index  File  Programs.   Electric
    Power Research Institute Report EA-1904.

Fox, D. G., 1981:  Judging air quality  model performance.   Bull. Amer. Meteor.
    Soc. 62, 599-609.

Fox, D. G.,  1984:   Uncertainty  in  air quality modeling.   Bull. Amer. Meteor.
    Soc. 65, 27-36.

Gillani, N. V.,  S.  Kohli and W.  E.  Wilson, 1981:  Gas-to-particle conversion
    of sulfur in  power  plant  plumes  -  I.  Parameterization of  the conversion
    rate for dry, moderately polluted ambient  conditions.   Atmos.  Environ.  15,
    2293-2313.

Hegg,   D.    A.,  L.  F.   Radke  and  P.   V.  Hobbs,   1984:    Measurements  of
    transformations   in  the  physical  and  chemical   properties   of clouds
    associated  with  onshore  flow  in  Washington  State.   J.  Climate  Appl.
    Meteor. 23,  979-984.

Hicks, B.  B., G.  D. Hess, M.  L.  Wesely,  T. Yamada, P.  Frenzen,  R.  L. Hart,  D.
    L. Sisterson, P. E.  Hess,  F.  C. Kulhanek,  R. C. Lipschutz and  G.  A. Zerbe,
    1981:    The  Sangamon  Field  Experiments:    Observations  of  the  Diurnal
    Evolution of the Planetary Boundary Layer  Over Land.  Report ANL/RER/81-1,
    Argonne National Laboratory.

Hicks, B.  B. and J. D. Shannon,   1979:   A  method for modeling the deposition of
    sulfur  by precipitation  over  regional  scales.     J.  Appl.  Meteor.  18,
    1415-1420.
                                     -45-

-------
Johnson, W. B.,  1983:   Interregional exchanges of air pollution:   model types
    and applications.  J. Air Poll. Assoc. 33, 563-574.

Meagher,  J.   F.  and  E.  M.  Bailey,  1983:    The seasonal  variation  of  the
    atmospheric SC^ to S0^2~ conversion rate.  J. Geophys. Res. 88,  1525-1527.

Mood,  A.  M.   and  F.  A.  Graybill,  1963:    Introduction to  the  Theory  of
    Statistics (2nd ed.).  McGraw-Hill, p. 198.

Niemann, B.  L.,  1982:   The  1980 data set for further evaluation of  regional
    air  quality/acid  deposition simulation  models.    Informal  report to  Acid
    Rain Staff, U.S. Environmental Protection  Agency.

Schiermeier,   F.  A.  and P. K. Misra,  1982:   Regional Modeling  Subgroup  Report
    2F-M.  Submitted to Work Group 2:  Atmospheric Sciences and Analysis.

Shannon, J.  D.,  1979:   A Gaussian  moment-conservation diffusion  model.   J.
    Appl. Meteor. 18, 1406-1414.

Shannon,  J.   D.,  1981:    A  model  of  regional  long-term   average  sulfur
    atmospheric  pollution, surface removal,  and  net  horizontal flux.   Atmos.
    Environ.   15,  689-701.

Shannon,  J.   D.,   L.   Kleinman,  C.  Benkovitz,   and  C.  Berkowitz,   1982:
    Intercomparison  of  MAP3S models  of  long-range  transport  and  deposition.
    Preprints, 3rd  Joint Conf.  on Applications  of  Air Pollution Meteorology,
    American Meteorological Society, Boston, pp.  29-32.

Sheih,  C.  M., 1977:   Application  of a  statistical trajectory  model  to  the
    simulation of  sulfur pollution  over  northeastern United  States.   Atmos.
    Environ.   11,  173-178.

Streets, D. G., D. A. Knudson and J. D. Shannon,  1983:  Selected  strategies to
    reduce acidic deposition in the U.S. Environ. Sci. Technol. 17,  474A-485A.

Taylor,  G.  I., 1921:   Diffusion by continuous  movements.  Proc. Math.  SL20,
    196-212.

Voldner,  E.   C.   and  J. D.  Shannon,  1983:    Evaluation  of  predicted  wet
    deposition  fields  of  sulfur  in  eastern   Canada.    Transactions,   Air
    Pollution  Control  Association Specialty Conference  on The Meteorology of
    Acid Deposition, 387-399.

Wesely, M. L.  and  B.  B.  Hicks,  1977:  Some  factors  that affect the deposition
    rates  of  sulfur dioxide  and similar  gases  on vegetation.   J.  Air Poll.
    Control Assoc. 27, 1110-1116.

Wesely,  M.  L. and  J.  D. Shannon,  1984:    Improved  estimates  of  sulfate  dry
    deposition in eastern North America.  Environ. Progress _3_,  78-81.

Wesely, M. L., D. R. Cook, R. L. Hart and R. E.  Spear, 1984:   Measurements and
    parameterization  of  particulate  sulfur  dry deposition  over  grass.    J.
    Geophys.   Res. 90, 2131-2143.
                                     -46-

-------
Yaraada, T.,  1977:  Analyses of the Sangamon data.  Argonne National Laboratory
    Radiological and Environmental  Research Division Annual Report ANL-77-65,
    Part IV, 49-52.
                                     -47-

-------
Appendix A:  Program HORZ Listing
                48

-------
 1.     //HORZ 003
 2.     //*MAIN ORC
 3.     // EXEC FGICLG
 i..     //SYSIN DD *
 5.     C ... USING WIND AND PRECIP ANALYSES FROM U MICH, THIS' PROGRAM WILL
 5.     C     PRODUCE THE TRAJECTORY STATISTICS REQUIRED BY ASTRAP.
 7.     C
 8.     C                  VARIABLES USED IN MORE THAN ONE SUBROUTINE
 9.     C
l.'J.     C      CONC :  STATISTICS DESCRIBING THE AIRBORNE PORTION OF THE
11.     C          TRACERS.
12.     C       CONCd,*,*)  AND  CONC(2,*,*> :  USED TO CALCULATE THE MEAN
13.     C           LOCATION IN THE X AND Y DIRECTIONS, RESPECTIVELY,
14.     C           OF THE AIRBORNE TRACER FOR EACH TIME STEP/SOURCE
15.     C           COMBINATION.
IS.     C       CONC(3,*,*)  AND  CONCU,*,*) :  USED TO CALCULATE THE STANDARD
17.     C           DEVIATION IN THE X AND Y DIRECTIONS OF THE AIRBORNE
13.     C           TRACER FOR EACH TIME STEP/SOURCE COMBINATION.
19.     C       CONC<5,*,*> :  THE SUM OF THE AMOUNTS OF AIRBORNE TRACER FOR
2.'j .     C           THE TIME STEP/SOURCE COMBINATIONS.
21.     C       CONC(S,*,*> : USED TO DETERMINE THE CORRELATION BETWEEN THE X
22.     C        AND Y POSITIONS IN THE ENSEMBLE PUFFS.
23.     C      DC :   THE NUMBER OF MEASUREMENTS CONTRIBUTING TO THE
24.     C          RESULTS FOR EACH TIME STEP/SOURCE COMBINATION.
25.     C      DELT :  THE NUMBER OF SECONDS IN A TIME STEP.
26.     C      DELX  AND  DELY :  LENGTH OF AN NMC GRID INCREMENT IN METERS
27.     C          AT 6£f DEGREES NORTH LATITUDE.
25.     C      DXP  AND  DYP :  THE LENGTHS IN NMC UNITS OF A PRECIPITATION GRID
29.     C          CELL IN THE X AND Y DIRECTIONS.
30.     C      DXWI   AND  DYWI :  THE LENGTH IN NMC UNITS OF A WIND GRID CELL
31.     C          IN THE X AND Y DIRECTIONS, RESPECTIVELY.
32.     C      FAC3 :  A SEASONAL FACTOR USED TO DETERMINE THE SPEED OF THE
33.     C          NOCTURNAL JET.
34.     C      INCPR :  THE TIME INCREMENT FOR V/IND/PRECIP DATA.
35.     C      IPER :  A NUMBER FROM 1 TO 4 REPRESENTING THE SIX-HOUR PERIOD
33.     C          OF THE DAY (GMT) TO WHICH THE WIND DATA APPLY.
37.     C      LANDE:  REPRESENTS THE TRAJECTORY CALCULATION GRID.  CELLS WHICH
33.     C          HAVE NONZERO VALUES WILL BE USED AS VIRTUAL SOURCES.
39.     C          THESE CELLS WILL BE ASSIGNED SEQUENTIAL IDENTIFICATION
4J.     C          NUMBERS PROCEEDING ROW-BY-ROW FROM THE LOWER LEFT CORNER.
'-:.     C      MEAS :  THE NUMBER OF MEASUREMENTS READ.
^2.     C      NDAY :  THE NUMBER OF DAYS THE PROGRAM SHOULD MODEL.
43.     C      NF :   THE NUMBER OF DATA SETS FOR THE DESIRED PERIOD.
44.     C      NS :   THE NUMBER OF SOURCE (LAND-CONTAINING) CELLS.
45.     C      NTS :  THE NUMBER OF TIME STEPS TO BE CALCULATED.
iJ.     C      NX?  AND  MYP :  THE NUMBER OF CELLS IN THE X AND Y DIRECTIONS IN
47.     C          THE PRECIPITATION GRID.
43.     C      NXWI   AND  NYWI :  THE NUMBER OF CELLS IN THE X AND Y DIRECTIONS
49.     C          IN THE WIND GRID.
5.-J.     C      PRECIP :  24-HOUR PRECIPITATION TOTALS.
51.     C      ?Z :   SET TO PI DIVIDED BY 2.
52.     C      RO :   A CONSTANT USED TO CONVERT LATITUDE AND LONGITUDE TO NMC
£3.     C          UNITS.
E-,.     C      TR :   THE FRACTION OF THE TRACER LEFT FOR EACH TIME STEP/SOURCE
55.     C           COMBINATION.
5.:.     C      U  AND  V :  X AND Y COMPONENTS Or THE TRANSPORT WIND.
57.     C        ALONG 30 DEGREES WEST LONGITUDE X IS WEST-EAST AND V IS
53.     C          NORTH-SOUTH.
59.     C      UTROP :  THE FRACTION OF A TRACER LOST TO THE UPPER TROPOSPHERE
c.I .     C          QY CONVECTION. '
        C      WC :   THE NUMBER OF WET DEPOSITION EVENTS FOR EACH TIME STEP/
        C          SOURCE COMB I NATION.
        C      WD :   WET DEPOSITION DATA
                                       49

-------
 54.

 So.
 67.
 63.
 59.
 72.
 11 .
 72 .
 73 .
 74.
 75 .
 73 .
 77.
 73.
 79.
 a .y .
 21 .
 o o .
 SS.
 57.
 S3 .
 S3.
 50.

 92.
 93 .
 94 ,
 95.
 96 ,
 97 .
 C £
 99.
1
l.C'3 .
l/i'9 .
1 1 £',
Hi .
112.
113.
11'.
115,
116.
1 17.
11C.
C
r
C
C
C
C
C
c
C
c
c
c

c
c

c
c
C
C
c
c
c
c
c
c
c
c
c
c
       WD(1,*,*)   AND WD(2,*,*) :   USED TO COMPUTE THE MEAN
           LOCATION IN THE X AND Y DIRECTIONS, RESPECTIVELY,
           OF  WET DEPOSITION FOR EACH TIME STEP/SOURCE COMBINATION.
       WD(3,*,*)   AND  WD<4,*,*> :  USED TO COMPUTE THE STANDARD
           DEVIATION IN THE X AND  Y DIRECTIONS OF WET DEPOSITION FOR
           EACH TIME STEP/SOURCE COMBINATION.
       WD<5,*,*>  :   THE AMOUNT OF  MASS REMOVED BY WET DEPOSITION FOR
           EACH TIME STEP/SOURCE COMBINATION.
       WD(6,*,*)  :  USED TO DETERMINE THE CORRELATION BETWEEN THE X AND
           Y POSITIONS IN THE ENSEMBLE PUFFS.

   .  READ PARAMETERS AND INITIALIZE VARIABLES
     CALL I NIT

   ,  READ DATA AND  COMPUTE TRAJECTORIES
     CALL MODEL

   ,  COMPUTE STATISTICS AND PRINT  RESULTS.
     CALL RESULT
     STOP
     END

     SUBROUTINE INIT

   .  READ USER-SPECIFIED PARAMETERS AND INITIALIZE VARIABLES.

                         VARIABLES USED ONLY  IN INIT

      DXEM  AND  DYEM :  NMC LENGTH OF AM TRAJECTORY GRID CELL  IN  NMC
          UNITS IN  THE X AND Y DIRECTIONS, RESPECTIVELY.
      NXEM  AND  NYEM :  NUMBER OF CELLS IN THE X AND Y DIRECTIONS,
          RESPECTIVELY, ON THE TRAJECTORY GRID.

     COMMON/WIND/IK 17,19),V<17,19),DXWI,DYWI,XWILL,YWILL,DELT,X( 17),
    X  Y( 19),DELX,DELY,NXWI,NYWI,MEAS -
     COMMON/PREC/PRECIP(5#,45),WD<6,23,123),DXP,DYP,XPLL,
    %  YPLL,UTROP,WC<23,123>,NXP,NYP
     COMMON/STAT/ XC(29,123 ) ,YC(29,123 ) ,TR(29,123>,CONC(6,28,123),
    Y,  F ACS, DC ( 28, 1 23 ) , P2 , RO , XEMLL , TRL < 29 , 123 >,
    %  YEMLL,DXEM,DYEM,NS,NTS,IPER,LANDE(19,16),NXEM,NYEM
     COMMON /CONTRL/ NDAY,NF,ITP,  IOUTST,IOUTEN,IDD(4)

   .  DEFINE MISCELLANEOUS CONSTANTS
     P2=3.14159/2.
               R0=6397./381
                              CALCULATE NMC COORDINATES  FOR  SOURCE
125 .
125.
127.
 ...  READ DATA FOR  LANDE
     CELLS.
     NXEM=19
     NYEM=16
     XEMLL = -11 .5
     YEMLL=20.5
     DXEM=1 .
     DYEM=1.
     READ 2023,( ,!=!,NXEM),0=1,NYEM)
2023 FORMAT(1913)
     NS=0
     DO 1030 0 = 1 , NYEM
       DO 1030 1 = 1, NXEM
         IF  (LANDE(I,0>.EQ.0> GO TO 1030
         NS = MS +• 1
         IF  (NS.GT.123) GO TO 40
         LANDE = MS
         XC< 1 , NS ) = XEi1LL •> DXEX"< 1-1 )
                                       50

-------
128.                YC< 1,NS)=YEMLL - D YEM* ( 0- 1 >- 1 6 .3
129.       1030     CONTINUE
_;4 .            GO TO 60
131 .         42 PRINT 2050
132.       2050 FORMAT* /'^'WARNING:  TOO MANY SOURCE CELLS.  ONLY 123 USED #*'/)
103.            NS=123
134.         60 PRINT 2070, NS
135.       2070 FORMAT* 15,'  SOURCE CEILS WERE READ1)
135.            PRINT 2080, ( ')
1~3.            PRINT 2100,XPLL,YPLL
174.            PRINT 2110,DXP,DYP
175.            XPLL=XPLL-DXP/2.0
173.            YPLL=YPLL+DYP/2.0
177.            YPLL=YPLL-15.0
173.      C
179.      C ... INPUT MISCELLANEOUS PARAMETERS.
135.            READ 2140,UTROP
1C: .       2140 FORMAT(F5.0)
132.            PRINT 2150,UTROP
133.       2150 FORMAT*/1  FRACTION OF  SULFUR LOST TO THE" UPPER TROPOSPHERE BY CONV
              XECTIVE PRECIPITATION = ',F6.3>
         C
         C ...  READ DAY/TIME INFORMATION.
               READ 2220, I DO , MDAY , NTS
          2220  FORMAT* 4* IX, 12), 14. IX, 12 )
               INCPR=6
           240  NF=4rtNDAY
                                       51

-------
192.            ITP=20
1-Z.       2250 FORMAT(4(IX, 12) )
1?-' .            PRINT 2250,17?
Ii5.       2260 FORMAT = 9999.
210.                WC<1,0)=0.0
211.                DC(I,0)=0.0
212.                TRLCI+i,0)=0.
213.                TR(1 + 1,0 )=0.
21-;.                DO 1330 K=l ,o
215.                  CONC,NXP,NYP
2-!C.            COMMON/STAT/  XC ( 29 , 1 23 ) , YC< 29 , 1 23 ) , TR( 29 , 1 23 ) , CONC( 6 , 28 , 123),
2-; 4.           %  FAC3,DC(28, 1 23 ) , P 2 , RO , XEMLL , TRL ( 29 , 123),
245.           %  YEMLL,DXEM,DYEM,NS,NTS,IPER,LAilDE(19,16),NXEM,NYEM
2 5.            COMMON /CONTRL/ NDAY,NF,ITP,IOUTGT,IOUTEN,IDD(4 )
Z-\l.            DIMENSION IDPU)
         C
               NEWMON=1
               MEAS=0
               NF=0
         C
253.      C ... COMPUTE  TRAOECTORIES
254.            DO 1222  ITI = I ,,IDAY
255.              DO 1222 IPER=1,4
                                        52

-------
256.     C
257.     C ...     SET SEASONAL VARIABLES IF A NEW MONTH WILL BE STARTED.
258.               IF (NEWMON.NE.i> GO TO S3
259.               NEWMON=0
263-.               ITEMP=IDD(1>/3+l
251.               GO T0(40,10,20,30,40),ITEMP
252.     C ...     SPRING (MARCH-MAY)
253.        10     FAC3=1.I5
254.               GO TO 50
255.     C ...     SUMMER (JUNE-AUGUST )
235.        20     FAC3=1.25
257.               GO TO 50
258.     C ...     FALL (SEPTEMBER-NOVEMBER)
259.        30     FAC3=1.15
270.               GO TO 50
271.     C ...     WINTER (DECEMBER-FEBRUARY )
272.        40     FAC3=1.05
273.        50     MEAS=MEAS+1
274.     C
275.     C . .  .     READ WIND DATA.
275.           READ <10,END=996) U,V
277.     C
27G.     C ...      READ PRECIPITATION DATA.
279.       998 READ (ITP,END=999> PRECIP
23.0.     C
231.           GO TO 1220
282.       999 REWIND ITP
283.           ITP=ITP+1
2C4.           IF (ITP.GT.22) GO TO 996
255.           GO TO 998
215.      1220     CALL TRAJ
237.      1222 NF=NF+1
288.       995 CONTINUE
289.           RETURN
290.           END
291.
292.
293.           SUBROUTINE TRAJ
294.     C
295.     C ... TRAJ TRANSPORTS TRACERS FOR A TIME STEP AND COLLECTS
295.     C     STATISTICS.
297.     C
298.     C                         VARIABLES USED ONLY IN TRAJ
299.     C
300.     C      DA ,   DB ,   DD ,  AND  DE :   RECIPROCALS OF THE. DISTANCES FROM
3^1.     C          THE TRACER'S LOCATION TO THE FOUR CLOSEST WIND MEASUREMENTS.
3.S'2.     C      Dl ,   D2 ,   D3 ,  AND  D4 :   PERPENDICULAR DISTANCES FROM THE
303.     C          TRACER'S LOCATION TO THE FOUR CLOSEST WIND CELLS.
304.     C      DT :   NUMBER OF SECONDS IN A CENTER-DIFFERENCED TIME STEP.
305.     C      FRMOV :  THE FRACTION OF A TRACER REMOVED--BY WET DEPOSITION
306.     C          AND CLOUD CONVECTIVE EXCHANGE.
307.     C      FRMAI1 :  THE FRACTION OF THE TRACER REMAINING AFTER WET
300.     C          DEPOSITION.
309.     C      FRMOV1 :  THE FRACTION OF THE TRACER DEPOSITED 3Y WET DEPOSITION.
310.     C      IM :   THE TIME STEP BEING PROCESSED.
311.     C      I :   THE NEXT TIME STEP.
312.     C      IQ  AND  JQ  :  THE COLUMN AND ROW, RESPECTIVELY, OF THE TRACER'S
313.     C          LOCATION ON THE TRAJECTORY CALCULATION GRID.
314.     C      IX :   COLUMN Of THE WIND GRID TO THE LEFT OF THE POINT OF
315.     C          INTEREST.
315.     C      IY :   ROW OF THE WIND GRID BELOW THE POINT OF INTEREST.
317.     C      IXX :  COLUMN OF THE WIND GRID TO THE RIGHT OF THE POINT OF
318.     C          INTEREST.
319.     C      IYY :  ROW OF THE WIND GRID  ABOVE THE POINT OF INTEREST.
                                       53

-------
320.      C      JX  AND  JY :   COLUMN AND ROW ON THE PRECIPITATION GRID FOR THE
321.      C          TRACER'S NSW LOCATION.
322.      C      MM :   INITIAL  VALUE OF THE INNER LOOP INDEX SO ONLY TIME
323.      C          STEPS THAT HAVE VALUES WILL BE PROCESSED.
324.      C      PP :   PRECIPITATION IN THE TRACER'S CELL.
325.      C      UU  AMD  VV :   THE ESTIMATIONS OF THE U AND V COMPONENTS OF
326.      C          THE WIND IN THE TRANSPORT LAYER.
327.      C
328.            COMMON/WIND/U<17,19),V(17,1 9 ) , DXWI,DYW!,XWILL,YWILL,DELT,X{17>,
329.           %  Y<19),DELX,DELY,NXWI,NYWI,MEAS
330.            COMMON/PREC/PRECIP(50,45),WD<5,23,123>,DXP,DYP,XPLL ,
331.           %  YPLL.UTROP,WC(23,123),NXP,NYP
332.            COMMON/STAT/ XC<29,123),YC(29,1Z3),TR<29,123),COMC(5,28,123),
333.           %  FAC3,DC(28,123),P2,RO,XEMLL,TRL(29,123),
334.           %  YEMLL,DXEM,DYEM,NS,NTS,IPER,LANDE<19,16 ) ,NXEM,NYEM
335.      C
336.      C ... COMPUTE TRAJECTORIES FOR THE TIME STEP/SOURCE COMBINATIONS.
337.      C
333.      C ... UNTIL  THERE ARE DATA FOR EACH TIME STEP, ONLY LOOK AT
339.      C     THE TRAJECTORIES FOR WHICH THERE ARE DATA.
340.            MM=MAX0(1,(NTS+1-MEAS) )
341 .            DO 1040 ITS = MM,NTS
342.              DO 1040 J = l ,MS
343.                IM=NTS + 1  -ITS
344.                I=IM+1
3-15.      C
346.                DT=DELT
347.      C         IF THIS IS  THE FIRST TIME STEP, FORWARD DIFFERENCING IS DONE
340.      C         SO  DT  IS  DIVIDED BY  2.
349.                IF (IM.EQ.l ) DT=DELT*.5
350.      C ...     CHECK TO SEE IF THE TRACER IS OUTSIDE THE WIND GRID ( XC  AND
351.      C          YC  ARE SET TO 9999 WHEN A TRACER IS OUTSIDE THE WIND GRID).
352.                IF (XC(IM.J>.GE.1000.> GO TO 30
353.      C
354.      C ...     USE BILINEAR INTERPOLATION TO ESTIMATE THE WIND AT  THE CURRENT
355.      C         TRACER LOCATION.
356.      C
357.      C         FIRST OBTAIN THE FOUR  CLOSEST GRIDDED WIND VALUES.
358.                IX = -XWILL>/DXWI + 1.
359.                IY = (YWILL - YC(IM,J ) )/DYWI + 1.
360.            IXX=IX+1
3S1.            IYY=!Y+1
362.            IF (IXX.LT.l) IXX=1
363.            IF UYY.LT.l) IYY=1
364.            IF -V( IY )
374.                D3 = XC( IM,0  )-)(< IXX >
375.                D4 = YC(IM,0>-Y< IYY)
375.      C
377.      C ...     OBTAIN WEIGHTINGS.
378.                DA=1 ./< SCRT(D1*DH-D2*D2) + 1 .E-10)
379.                DB=1 ./< SQRT*1 . E-10)
381 .                DD=1 ./( SQRT(D3*D3+D4*D4)+1 . E-1.9)
382.                DENOM=1./
-------
334.     C ...     COMPUTE SPEEDS.
385.               UU=(DA*U< IX,1Y)-H)B*U( IX , I YY )+DE*U( I XX , IY >
386.          %     +DD*U< IXX,IYY) )*DENOM
387.               VV=(DA*V< IX, IY>+DB*V< IX, IYY)+DE*V< IXX,IY>
388.          %     +DD*V
405.               UMAP=UMAP*SIG/DELX
405.               VMA?=VMAP-SIG/DELY
A 07 .     C
403.     C ...     TRANSPORT THE TRACER  FOR  A TIME  STEP.
409.               XC(I,J)=XC(IM,J>+UMAP"DT
410.               YC( I ,J ) = YC( IM,J >+VMAP*DT
411.               TR(I,J )=TR(IM,J )
.112.               TRL(I.J)=TRL(IM,J)
413.     C ...     HAVE TRAJECTORIES DEPARTED THE WIND  GRID?
4U.               IF   GO TO  30
416.               IF (YC( I,J).GT.{Y(1 ) + 2. »  GO TO  30
417.               IF (YC< I,J ).LT.(Y(NYWI )-2. >> GO  TO  30
418.     C
419.               FRMAIl=1.0
420.     C
421.     C ... IF THE TRACER  IS NOT WITHIN  THE  PRECIPITATION  ANALYSIS,
422.     C     EXTEND THE LAST CELL ON THE  ROW  OUTWARD.
"23.               JX = (XC( I ,0 5-XPLD/DXP  +  1.
424.               JY=(YPLL - YC(I,J))/DYP  +  1.
425.               IF  JX = NXP
428.           IF (JX.LT.1 ) JX=1
A29.     C
430.     C ...     CALCULATE WET DEPOSITION  IF  THERE  IS PRECIPITATION.
431.               PP=PRECIP(JX,JY>
432.               IF < PP.LE.1.0) GO TO  20
433.               FRMOV = SQRT( . 1»PP )
434.     C ...     DURING THE WINTER { FAC3  =1.05 IN  THE  WINTER)  IN NORTHERN
435.     C         AREAS, WET DEPOSITION  IS  REDUCED BECAUSE THERE  IS SNOW
436.     C         INSTEAD OF RAIN.
437.               IF (FAC3.EQ.1.05  .AND.  YC(I,J).LT.0) FRMOV=FRMOV*.5
430.               IF !FRMOV.GT.1.0) FRHOV=1.0
'139 .               FRMOV1 = < 1 .0-UTROP )*FRMOV
^'40.               FRMAI1 = 1 .0-FRMOV1
                   TRL( I,J )=TRL< !M,J >*( 1.-FRHOV)
                   TR( I,J >=TR( IM,J )-FRt1AIl
         C
"-M.     C ...     GATHER STATISTICS IF  ENOUGH  STEPS  HAVE  ELAPSED.
'- J 5 .               WTF RAC = TR ! IM , J ) '-• F RMOV 1
446.               V/D < 5 , If j, J ) =\/D { 5 , I M , J ) +WTFRAC
447.               WD< 1 , IM.J )=>/D{  1 , IF1,J ) + XC( I ,0 --WTFRAC
                                       55

-------
448.
449.
450.
451.
452.
453.
454.
455.
456.
457.
458.
459.
460.
461 .
462.
453.
464 .
465.
456.
467.
463 .
469.
''-70.
471.
472.
473.
474.
475.
475.
477.
478.
479.
480.
481 .
482.
433.
484.
•1- i i. ,
AS3.
4 £ 4 .
=WD(2,IM,
         WD(3,IM,0)=WD(3,IM,
         WD<4,IM,0 )=WD(4,IM,
         WD(5,IM,0 )=WD<5,IM,
         WC(IM,0)=WC
                             0)+YC{I,J}*WTFRAC
                             0)+XC(I,0)*XC + YC + XC< I,J >*YC=CONC(2,
                     >=CONC(3,
                     >=CONC(4,
          CONC<6,IM,0)=CONC<6
          DC( IM,0 ) = DC( IM.O> + l
          GO TO 1040
          IM,
          IM,
          IM,
J )=CONC(5, IM,
         ,IM,
         .0
0 ) •»• XC( I,0)*DRFRAC
0) + YC(1,0 >*DRFRAC
0)*XC(I,0 )*XC( I ,
0)+YC(1,0 >*YC( I,
0)+DRFRAC
0) + YC*DRFRAC
0)*DRFRAC

0>*DRFRAC
c
c
c
c
c
c
c
:  ...      TAG  TRACERS  OUTSIDE GRID
   33      XC(I,0 ) = 9999.
          YC< I,J) = 9999.
          TR(1,0)=0.0
          TRL( I,0 )=0.0

 1040      CONTINUE
      RETURN
      END

      SUBROUTINE  RESULT

    .  FINISH  CALCULATION. OF STATISTICS AND OUTPUT RESULTS.

                         VARIABLES USED ONLY IN RESULT

       TITLE  ,  WDLBL  ,   AND  REMLBL :  LABELS FOR OUTPUT.

      COMMON/PREC/PRECIP<50,45),WD(6,28,123),DXP,DYP,XPLL,
    %  YPLL.UTROP,WC(28,123),NXP,NYP
      COMMON/STAT/  XC ( 29 , 1 23 >, YC ( 29 , 12.3 ), TR{ 29 , 1 23 ), CONC( 6 , 28 ,123 ),
    %  FAC3,DC<20,123),P2,RO,XEMLL,TRL(29,123),
    %  YEMLL,DXEM,DYEM,NS,NTS,IPER,LANDE(19,15),NXEM,NYEM
      COMMON  /CONTRL/  NDAY , TJF , ITP , IOUTST, IOUTEN , I DD( 4 )
      DOUBLE  PRECISION TITLEt7,6 ) ,WDL3L(5 ) ,REMLBL(5>,NUL3L(7)
      DATA
                          , ' MEAN LO1
                          1 MEAN LO1,
                             S.D. I ' ,
                             S.D. I'.
                              SUM 0',
                          , ' X-Y COR'
                          OSITION ','
    '/.  REMLBL/1  RZMAINI ' , 'NG TRACE1,
    %  NULBL/1   NUMBER1,'  OF MEAS','
                           FOR'/
        TITLE/3*'
             3* .
             4*'
             4*'
             4*.
              4*'
        WDLBL/' WET DEP:
        REMLBL/1 REMAIN I ' , '
        NULBL/'  NUMBER','
           1 THE RES','ULTS
, 'CATION I
'CATION I '
1 N X DIRE '
1 N Y DIRE '
'F THE MA'
, 'RELATION
EXTRACT I ' ,
'R FRACTI '
UREMENTS' ,
1 , 'N X DIRE
, ' N Y DIRE1
, 'CTION OF '
, 'CTION OF '
, 'SSES FOR'
'ON
, 'ONS EXTR1
1 CONTRIB1 ,
1 , 'CTION OF1 ,
, 'CTION OF ' ,
' / ,
'/,
, 'ACTION '/,
'UTING TO1 ,
      COMPLETE THE CALCULATIONS  FOR   CONC  ,
      DO 1030 1=1,NTS
        DO 1030 0 = 1 ,NS
          IF (CONC(5,I,0J.LE.0.)  GO  TO 10
          COMC( 1 , I .0 )=CONC< 1,1,0 )/COt!C(5, I
          CONC(2, I , J >=COMC(2, I ,0 )/COi\!C<5, I
          IF (DC( I ,J ).LT.3. )  GO  TO  21
          DD = DCiI.0)/(DC(I,0 >-l.0)
                                                     0 )
                                                     0 )
          DQ = DC' I ,0 ;/'( DC( I ,0 )-2.
          CONC(3, I ,0 > = DD«((CONC(
            CONC' 1,1,0))
                                   I,0 )/CONC<5,1,0) )-CONC(1,I ,0 >*
                                        56

-------
512,
513.
514.
515.
515.
517.
518.
519.
520.
521 .
522.
523.
524.
525.
525.
527.
•523 .
32S.
 -.
5 -' 4 .
545 .
   .
525 .
5"£G.
£57 .
5£1 .
 5- o
 j t..
T- ~ O
O SO .
3 O 4 .
5 5 o .
572 ,
573 .
574 ,
575 .
CONG (4, !,J) = DD*« CONG (4
  CONC( 2, ! ,0 ) )
CONCC6, I,0)=DQ*((CONC<6
 CONC<2, 1,0 )>
IF =CCNC< 2,1,0)+ 15.0
CONC( 3,1,0 >=SQRT(CONC<3, 1,0))
CONC( 4,1,0 )=SQRT(CONC<4, 1,0))
GO TO  1030
                                   1 ,0 }/CONC( 5 , I ,0 ) )-CONC( 2 , I ,0)*

                                   I,J)/CONC<5, I,0))-CONC< 1 ,1 .0)*
»
) )
                                            GO
                                            GO
                                          TO
                                          TO
                                                   20
                                                   20
   10

   21
   20
     CONCt1,
     CONC(2,
     CONC<6,
     CONC(3.
     CONC(4,
     CONC(2,
                  0 )=0.
                  0)=0.
                  0)=0.0
                  0)=DKP/SQRT(12,
                  0)=DYP/SQRT(12.
                   0 )=CONC(2,I,0 ) + 16.0
 1030
     CONTINUE
C
C
               COMPLETE THE CALCULATIONS
               DO 1050  1=1,NTS
                 DO 1060 0=1,NS
                           FOR  WD
                   IF (WD<5,I ,0).LE.£T. )  GO  TO 40
                   WD(1,I,0)=WD(1,1,0)/WD(5,1,0)
                   WD(2,I,0)=WD(2,I,0)/WD(5,I,0 )
                   IF (WC( I,0 ).LT.3. >  GO TO 51
                   WW=WC(I,0)/(WC(I,0)-1.0)
                   WQ=WC(I,0 )/(WC(I,0)-2. )
                   WD < 3,I,0)=WW* « WD(3,I,0)/WD(5,
                            )=WW*((WD(4,I,0 >/WD(5,
                            >=WQ*((WD(6,I,0 )/WD(5,
                            ,I,0 ).LE.(DXP*DXP/12. !
                            ,I,0 ).LE.(DYP*DYP/12. )
          WD(4,1 ,0
          WD(5 , I , 0
          IF (WD(3
          I r ( '.;D ( 4
          WD(3,I,0
          WD(4,I
          WD(2, I
          GO TO
                                   I,0»-WD( 1
                                   I ,0 »-WD(2
                                   1,0)>-WD(1
                                   > GO TO 50
                                   ) GO TO 50
                                                   I,0)*WD(1
                                                   I,0>*WD(2
                                                   I,0)*WD(2
                                                   ,I,J
                                                   !i!o
                       ))
                       »
                       »
            40

            51
            50
WD( 1 .
WD(2,
WD(6,
WD(3,
WD(4,
WD(2,
I ,0
I ,0
I ,0
1,0
I ,0
I ,0
             >=SQRT(WD(3,I ,0) )
           ,0 )=SQRT(WD(4, 1,0))
           ,0 )=WD(2, I ,0 ) + 16.0
           1060
             )=0.
             }=3.
             )=0.0
             ) = DXP/SQRT( 12. )
             ) = DYP/SQRT( 12
                                          >
                            )=WD(2, I ,0 ) + 16.0
          1060
          CONTINUE
 1080

 2093
                                      :>,L=1,7),R£MLBL,(J,
-------
577.
573 .
579.
530.
531 .
532.
583.
584.
535 .
585.
587.
588.
539.
590.
591.
592.
593.
534.
595 .
59sl
597.
59C.
599.
6."' rj .
5 C' 1 .
6.: 2.
      END
//GO.FT09F001 DD DSN = B254S8 .TSUM80HC , UNI7=ALLPFRM,
// DISP=(NEW,CATLG),DCB=(LRECL=X,BLKSIZE=6144,RECFM=VBS),
// SPACE=(TRK,(10,5),RISE)
//GO.FT10F001 DD LABEL = ( 19,NL, ,IN ),UNIT=READ5250,OISP = , UN IT=READ625£T,
// DISP=,DCB=,
// VOL=SER=201436
//GO.?T21FJr01 DD DSN=OUL80. PF8HR , LABEL = { 8 , SL , , I N> , UNIT = READ6250,
// DISP = (OLD,KEEP > , DCB = (LRECL = X,BLKSIZE = 6136,RECFM=VBS),
// VOL=SER=201435
//GO.FT22F001 DD DSN=AUG80. PF6HR , LABEL = { 9 , SL., , IN >, UN IT = READ6250,
// DISP=(OLD,KEEP),DCB=,
// VOL=
'SER=201436









//GO.SYSIN DD *
0
0
0
0
0
0
0
0
0
0
0
vj
0
0
0
0
0
0
0
0
0000000
0 0123 0000
0 0120121122 0 0
0114115116117118119
01041051061071081091101
87
71
55
42
0
0
0
0
2T
0
0
88 89 90 91 92 93 94
72 73 74 75 76 77 78
57 58 59 50 61 62 63
43 44 45 45 47 48 49
30 31 32 33 34 35 35
0 19 20 21 22 23 24
0 0 3 11 12 13 14
0004567
0000200
0000000
0000000
0.500
6
MO
/*
1
DA

80 0 92 28
YR TI NDY NT

0 SS
0 0
0 0
0 0
11 0
95 95
79 80
64 65
50 51
37 38
25 25
15 16
8 9
0 0
0 0
0 0
UTROP
MO , DA


0
0
0
0
0
97
81
56
52
39
27
17
10
3
1
0

,YR


0
0
0
0
0
98
32
57
53
40
28
18
3
0
0
0

,HR


0
0
0
0
0
99
83
68
54
41
29
0
0
0
0
0

0 0
0 0
0 0
0 0
01.121
0
0
0
0
13
0
0
0
0
0
100101102103.
84 85
69 70
55 0
0 0
0 5
0 0
0 0
0 0
0 0
0 0

86
0
0
0
0
0
0
0
0
0

OF START, #






0
0
0
0
0
0
0
0
0
0

DAYS


0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0

*


INDICATOR
ARRAY FOR
VIRTUAL
SOURCES OF
TRAJECTORY
STATISTICS











#• TIMESTETS


                                        58

-------
Appendix B:  Program VERT Listing
                59

-------
  1.      //VERT JOB  REGION»ZU3rK,TIME=ljer.,CLASS-X
  2.      //"MAIN ORG=RM080
  3.      /'/  EXEC FGICLG
  4.      //SYSIN DD  *
  5.      C       THIS PROGRAM  CALCULATES THE  HORIZONTALLY  INTEGRATED VERTICAL
  5.      C      PROFILES  OF  GAS,  PRIMARY AND  SECONDARY  PARTICLES
  7.      C      AS A  FUNCTION  OF  TIME  AFTE3.  RELEASE  (PLUME  AGE).
  3.      C
  9.      C       DIURNAL  AND VERTICAL  VARIATIONS  OF  PARAMETERS ARE SIMULATED,
 10.      C      BUT NOT HORIZONTAL  VARIATIONS.  TYPICAL PATTERNS OF VERTICAL EDDY
 11.      C      DIFFUSIVITY  ARE  SPECIFIED BY  LAYER AND  TIME OF DAY.  THE RATE OF
 12.      C      TRANSFORMATION FROM GAS  TO PARTICLE  IS. A FUNCTION OF TIME OF DAY.
 13.      C      DEPOSITION VELOCITIES  FOR S02 OR  NO/N02 AND S04 OR N03
 14.      C      ARE SPECIFIED  BY  TIME  OF DAY.
 15.      C
 15.      C       THE  TYPICAL STABILITY CYCLE  IMPLIED 3Y THE EDDY DIFFUSIVITIES IS
 17.      C      A NOCTURNAL  INVERSION, FORMING  AROUND SUNSET,  DEEPENING DURING
 13.      C      THE NIGHT, AND RISING  AND ERODING FROM  BELOW DURING THE MORNING.
 19.      C      DURING THE AFTERNOON THE ATMOSPHERE  IS  UNSTABLE UP TO THE
 23.      C      SYNOPTIC-SCALE MIXING  HEIGHT.
 21.      C
 22.      C       THE  MODEL CAM HAVE A  MAXIMUM OF  9 LAYERS  WITHOUT REDIMENSIONING,
 23.      C      WITH  THE  NUM3ER  OF  LAYERS AND THE LAYER THICKNESSES BEING INPUT.
 24.      C      CONCENTRATIONS ARE. ALSO  CALCULATED AT THE  SURFACE.
 25.      C
 2o.      C      THROUGHOUT THE PROGRAM,  S02/N02 IS POLLUTANT 1, PRIMARY S04/N03
 27,      C      IS POLLUTANT 2,  AND SECONDARY S04/N03 IS POLLUTANT 3.
 23.      C
 29.      C
 3~.      C                          VARIABLES IN  LABELLED  COMMON1
 31.      C
 32.      C       ABVML  :   AVERAGE MASS OF S  OR N  ABOVE  THE  MIXING LAYER
 33.      C           FOR  EACH  TIME STEP  AND  POLLUTANT.
 34.      C       A123 ,.  A?  ,  AND  AS  :   PARTIAL  CALCULATIONS  FOR EACH LAYER.
 35.      C       CONV :   FRACTION OF A MOLECULE COMPOSED OF S  OR N (BY MASS) FOR
~3G.      C           EACH POLLUTANT.
 37.      C       DCOR :   AIR DENSITY AT  EACH  LEVEL AS A FRACTION. OF SURFACE
         C           DENSITY.
         C       DEPT :   DRY DEPOSITION  MASS  INCREMENT FOR EACH TIME STEP AND
         C           POLLUTANT,  FOR A  UNIT EMISSION.
 4i.      C       02 : THICKNESS  OF EACH LAYER IN METERS.
 .12.      C       DZD  ,  DZU  ,  AND DZZ  :   PARTIAL  CALCULATIONS  FOR EACH LAYER.
 43.      C       EMR  :  RELATIVE  DIURNAL EMISSION RATE.
 4-1.      C       FAC  :  SEASONAL  FACTOR  AFFECTING THE TRANSFORMATION RATE.
 45.      C       KZ : VERTICAL  EDDY DIFFUSIVITY  BY  LAYER AMD  HOUR.
 45.      C       NTS  :  NUMBER OF TIME STEPS.
 47,      C       NZ : NUMBER  OF  LAYERS.
 4S.      C       NZ1  :  ONE  LESS  THAN  THE NUMBER  OF  LAYERS.
 49.      C       SBUD :   REMAINING AIRBORNE  MASS  FROM A UNIT EMISSION OF
 53.      C           ATMOSPHERIC  S OR  N, BY  TIME  STEP AND  POLLUTANT.
 51.      C       SOXS :   AVERAGE  ONE-DIMENSIONAL  SURFACE CONCENTRATION OF
 52.      C           S OR N  FOR A UNIT EMISSION.
 53.      C       SRCLVL  : SOURCE LEVEL  POINTER.
 54.      C       TR : HOURLY  S02-S04  OR N02-N03  TRANSFORMATION RATE.
 E5.      C       USTAR  :   FRICTION VELOCITY  IM M/SEC.
 55.      C       VD : DEPOSITION VELOCITIES  3Y AG£  CF  THE  PLUME, HOUR OF
 57.      C           THE  DAY,  AND POLLUTANT.
 Sa.      C       VK : VON KARMAN CONSTANT.
 59.      C       Zl : HEIGHT  FOR SURFACE CALCULATIONS.
 £j.j.      C       22 : HEIGHT  OF  THE CONSTANT FLU!! LAYER.
 6:.      C       IPOL : POLLUTANT SYSTEM INDICATOR <1=SOX,  2=t!OX)
 5i!.      C
 53.      C
                                       60

-------
 54.      C                             VARIABLES ONLY  IN  MAIM
 55.      C
 56.      C      AVADML :   AVERAGE MASS CF SOX OR NOX AB'OVI  THE  MIXING LAYER
 67.      C          CORRESPONDING TO TRAJECTORY TIME STEPS.
 GG.      C      A'-'BGT :   AVERAGE MASS OF AIR30RME  SOX  OR  NOX COERTSPOND! NG TO
 53.      C          TRAJECTORY TIME STEPS.
 73.      C      AVDD :  AVERAGE DRY DEPOSITION CORRESPONDING TO TRAJECTORY
 71.      C          TIME STEPS.
 72.      C      AVSFC :   AVERAGE SURFACE CONCENTRATIONS CORRESPONDING TO
 73.      C          TRAJECTORY TIME STEPS.
 74.      C      NSEASN :   CODE FOR THE SEASON <1=SUMMER,  2-AUTUMN,  3=WINTER,
 73.      C          4=SPRING>.
 7,5.      C      NVI :  NUMBER OF VICTUAL SOURCE LEVELS.
 77.      C
 73.      C
 79.            COMMON /DP/ DZ<10),A123(9),A7(9>,A3(9 ) ,DZD<9 ), DZU(9 ), DZZC 9 )
 ZiJ.            COMMON /R/ DEPT( 163,3), 30XS( 158,3), SBUD( 168,3), ABVML( 168,3), Z1..Z2,
 31. .           y.  VK,USTAR,KZ(2,9,24),VD{2,24,2>,TR(24),DCOR{ 10 ) , COMV( 4 ) , EMR< 9 , 24 )
 32.            COMMON /INT/ NTS , fiZ, NZ1 , SRCLVL, FAC , I POL
 T3.            DOU3LE PRECISION DZ, Al 23 , A7 , A8 , DZD , DZU., DZZ
 •2-\.            DOU3LE PRECISION A1,A2,A3
 Co.            DOUBLE PRECISION SEASON(4 ) ,IP0(2 )
 e3.            REAL AVSFC (28, 3), AVDD (23, 3) , AV3GT< 28 , 3.) , AVABML ( 28 ,3 ) , KZ
 37.            INTEGER SRCLVL
 S3.            DATA SEASON/' SUMMER ',' AUTUMN ','  WINTER ',' SPRING  '/
 89.            DATA IPO/' SULFUR ' ,'NITROGEN ' /
 90.      C
 91 .      C ... INPUT DATA
 32.            READ 2010,NSEASN,IPOL
 93.       2010 FORWAT(2I3)
 34.            NZ=9
 95.            NTS-158
 i'S.            NVI =6
 37.            Z1 = 2.0
 30.            22=35.0
 39.            PRINT 2021,IPO(IPOL)
!£?.       2J021 FORMATS  CALCULATIONS FOR QXIDTS OF ' ,A8 )
131.            PRINT 2020,NZ,NVI,SEASCfU NSEASN>, NTS,Zl,Z2
132.       2020 FORMAT(I2,' LAYERS WITH THE  BOTTOM1,12,' BEING SOURCE  LEVELS.1/
103.           %  '  MODEL',AS,'WITH',14,' TIME STEPS.'/' SURFACE CONCENTRATIONS',
104.           X  '  CALCULATED AT',F6.1,'M AND CONSTANT FLUX  LAYER AT',
1J5,           %  76 . 1 , 'M. ' >
10S.      C
I.i7.      C
153.      C ... SET CONSTANTS.
               USTAR=0.4
i.l,            CONV<1)=32./54.
i;2.            COMV(2>=32./96.
•13.            CONV<3>=32./96.
114.            CONV<4)=96./64.
115.            IF (IPOL.EQ. 1 ) GO TO  10
iio.            COMV<1 ) = 14./46.
117.            CCMV(2)=14./G2.

ll^J.            CO;IV(4> = 62./46.
125.         10 CONTINUE
'.2;..      C ... CHECK PARAMETERS.
122.            IF riZ.LE.9 .AND. NVL.LE.9  .AND.  NSEASJH ^L£. 4)  GO TO 50
12?.            PRINT 2.750
:24.       2JET50 FORMATi'1 A USEP.-SUPPL.I^D PARAffETUR  IS  TOO  LARG.E  AND MAKES'/
125.           %  '  EXECUTION  IMPOSSIBLE.   EXECUTION HALTED.')
::ro.            STOP
127,         50 IF (iiZ.ST.l .AMD. NVI.GT.-Cf  .AND.  MSIASti. GT .3 .AND.  Z1.GT.0)
                                        61

-------
128.
129.
135.
131.
132.
133.
134.
135.
135.
137.
133.
139.
1 40 .
141.
142.
143.
144.
145.
145.
147.
143.
1-19.
i30.
I 51 .
152.
153.
154.
155.
156.
157.
1 58 .
1 = 3.
160.
161 .
162.
153.
134.
1G5.
166.
167.
153.
i i 3 .
1-0.
1 - [
1-2 .
1 7C .
1-4.
175.
176.
177.
178.
179.
130.
1 G 1 ,
132.
133.
1S4.
105.
135.
i'37.
183.
139.
Isfl.
191.
?

2070
2

C
80
C
C . . .

2090

2100


2110
C
C . . .


2121
2120


2122



2123
2130

2124


2125

1140
C
C


9934
C
r*
c
c
c . . .
c
c
c

c . . .
150



C . . .
160



C . . .
170

  GO TO 80
PRINT 2070
FORMAT* '  A USER-SUPP1_X£D PARAMETER  IS TOO SMALL ACT  MAKES1/
  1  EXECUTION IMPOSSIBLE.  EXECUTION HALTED.')
STOP
NZL=NZ-1
                 AMD DENSITY
                ,I=1,MZ)
READ THICKNESSES
READ 2090,
FORMAT*//1 LAYER THICKNESSES
READ 2090,,1=1,NZ)
FORMAT, 0 = 1 , 24), I =-1, NZ )
              1 )
              IPOL
              1> GO TO 2125
  READ
  READ
  READ
  READ
  FORMAT* 12F5
 DO 2124 KQ-1
 READ 2120, TR
  IF (IPOL.EQ
  GO TO  1140
  READ 2 120, WASTE
  READ 2120, WASTE
CONTINUE

INPUT DUMMY VALUES FOR DIURNAL EMISSION  VARIATION.
DO 9934  1=1 ,9
 DO 9934 0 = 1 ,24
EMR( 1,0 )=1.0

SET DEPOSITION VELOCITIES  FOR THE  FIRST  2.5  HOURS EQUAL TO
THOSE FOR AFTER 2.5 HOURS,  UNLESS  THEY ARE LE.SS  THAN A
SEASONAL MINIMUM.
IN SUBROUTINE  INTGRT THE PRIMARY  SULFATE. DEPOSITION VELOCITY
IS SET TO 1.0  CM PER SECOND  DURING THE FIRST 2.5 HOURS TO
ALLOW FOR GRAVITATIONAL SETTLING  OF LARGE PARTICLES.

GO TO < 150, 160,170,1S0),NSEASM
SUMMER
S02MIN=.0030
S02MAX= .0050
GO TO  190
FALL
S02MIN=.0015
S02 MAX =.0040
S04MIM=..a020
GO TO  190
WINTER
S02MIM-.
                         62

-------
192.           S04M!N=
193.           GO TO  190
194.     C ... SPRING
195.       130 S02MIN=.002S
196.           S02MAX=.0040
197.           S04MIN=.0020
133.       190 DO 1200 1=1,24
199.             VD!2,I,1 )=VD(1,1,1)
200.             IF ;VD(2,I,1 >.LT.S:02MIN> VDC2-,-! , 1 )=S02M1N
251.             IF ;VD(2,I ,1 ).GT.S02MAX> VD( 2,-L, 1 ) = S02MAX
202.             VD;2,I,2)=VD;1,1,2)
253.      1200   IF ; VD(2., I , 2 ) . LT . S04MIN > VD< 2,1, 2 ) = S04MI N
254.     C
205.     C ... PRINT DEPOSITION VELOCITIES.
255.           PRINT 2210
207.      2210 FORMAT;- FIRST 2.5 HOURS OF DISPERSION')
208.           PRINT 2220,>
2:5.           PRINT 2220,,I=1,24)
211.           PRINT 2210
2:2.           PRINT 2230,(VD;2,I,2),1=1,24>
213.      2230 FORMAT;' PARTICLE  DEPOSITION VELOCITY  IN  M/S',/,z(5X,i2F7.4,/>>
21.;.           PRINT 2230,;VD( 1,1,2),1 = 1,24)
215.     C
216.     C ... SET EDDY DIFFUSIVITIES FOR THE FIRST 2.5  HOURS  EQUAL  TO THOSE
217.     C     FOR AFTER THE FIRST  2.5  HOURS, EXCEPT  FOU THE FIRST
210.     C     AND SECOND LAYERS.
219.           DO 1240 1=1,NZ
225.             DO 1240 J=l,24
221.      1240     KZ<2,1,0)=KZ(1,1,0)
222.           DO 1250 J = l ,24
223.      1250   IF ;:
235.           PRINT 2260,<
231 .           PRINT 2270,TR
232.      2270 FORMAT;' HOURLY  TRANSFORMATION RATE',/,2;5X,i2F7.4,/)>
233.     C
234.     C ... CALCULATE, FOR LATER  USE, FACTORS  RELATED TO THE THICKNESSES
235.     C     OF THE LAYERS.
235.           DO 2280 K=l,NZ
237.             !*
245.             A2 = *DZ< K )*DZ< K >/( 4.*DZD(K) )
245.             A3=*DZU/4. )
Z-".7.             A123;X)=A1-A2->A3
243.             A7<;<) = ( DZD
2-'9.      2280   AS-vK) = DZ( K)*DZ{ K )/( 4 . *DZD( K ) )
205.           FAC=2-.5rt:ABS;NSEASN-3)
251 .     C
2B2.     C ... SET DUMMY VALUES FOR  USE  IN COMPUTING  DIFFUSION UP  ACROSS THE
252.     C     TOP OF THE MIXING  LAYER.
25-1.           oz(NZ+i ) = DZ;NZ>
               DCCR
-------
256.
257.
253.
259.
230.
231 .
252.
253.
254.
2S5.
265.
267.
253.
259.
270.
271.
272.
273.
274.
275.
275.
277.
278.
279.
280.
281 .
232.
233.
234.
2S5.
233.
287.
283.
289.
290.
231.
292.
293.
294.
293.
295.
297.
298.
299.
300 .
301 .
302.
303.
304 .
305.
305.
307.
3.V8.
309.
310.
311 .
312.
313.
314.
315.
316 .
317.
313 .
313.
C
C

C






C
C
C

C
C
C










C








C

























C
C
    .  DO FOR EACH SOURCE. LEVEL.
      DO 1390 SRCLVL=1,NVI

        DO 1290 1=1,168
          DO 1290 0=1,3
            ASVMLtI,J)=0.
            S3UD<  I,0)=0.
            SOXS(I,0)=0.
 12.90       DEPT=AVSFC( 11,0 ) + St>XS( 1,0 5/6.JET
    AVDD( II,0)=AVDD< II,0)+DEPT< 12.2 ,0 >/Z.0
    IF (II.ME.l) AVDD< II ,0 >=AVDD( II ,0)-DEPT< 11,0 5/2.JET
    AVBGT( II,0)=, '  PRIMARY PARTICLE', 2< T21_ 14-( 1
     XX,F6.5>/>,' SECONDARY  PARTICLE',2(T21,14(1X,F5.5)/))
        PRINT 2350
 2360   FCRMATi ' MASS OF DRY  CrEPOSITTON: ' >
        PRINT 2350,AVDD
        PRINT 2370
 2370   FORMAT< ' MASS OF S OR N  ABOVE  THE  MIXING  LAYER'.1)
        PRINT 2350,AVABML
        PRINT 2380
 2380   FORMAT(' MASS OF AIRBORNE  S OR  N:1)
        PRINT 2350,AVBGT
        WRITE (11) AVSFC,AVDD,AVBGT
 1390 CONTINUE
      STOP
      END
      SUBROUTINE  INTGRT

C ... DIFFUSE UNIT ZMISSIONS  OT  SOX  OR XOX ORIG1NATTN-G •' IN A GIVEN
                               64

-------
223.     C     SOURCE LAYER.
321 .     C
322.     C
323.     C                           VARIABLE'S  ONLY  IN  INTCRT
324.     C
325.     C      AML  :  AMOUNT ABOVE THE MI X ING...LAYER FOtt  A  1-HOUR
326.     C          TIME STEP.
327.     C      CP  :  THE DISPLACEMENT  IN  METETtS  OF  THE CENTER  OF MASS  OF  A
328.     C          PUFF AFTER DIFFUSION.
329.     C      DEP  :  AMOUNT OF DEPOSITION  IN  A  1-HOUR TIME STEP.
330.     C      DTT  :  LENGTH OF A SU3STEP  LN SECONDS.
331.     C      FL2  AND  FL4 :  FLUXES OF GAS  AND PARTICLE.
332.     C      ID  :  EQUALS 1 AFTER THE FIRST  2.5 HOURS.   EQUALS 2
333.     C          DURING THE FIRST 2.5 HOURS.
334.     C      L  :  HOUR OF DAY SEING  MODELLED.
335.     C      Ml  ,  M2 , AND  M3 :  THE  FRACTION OF THE MASS  ORIGINALLY1 IN
336.     C          THE CURRENT LAYER WHICH  HAS DIFFUSED  TO THE  LOWER
337.     C          ADJACENT, THE CURRENT, AND  THE UPPER  ADJACENT LAYERS,
333.     C          RESPECTIVELY.
339.     C      NTI  :  NUMBER OF SUBSTEPS  IN AN HOUR.
340.     C      ORIGHR :  HOUR OF THE DAY  THE PLUME  ORIGINATED.
341.     C      P  :  FRACTIONAL ALLOTMENT  OF POLLUTANT  PUFFS FROM EACH  LAYER
342.     C          TO THE LOWER, CURRENT, AND  UPPER LAYERS FOR  EACH POLLUTANT.
343.     C      RATIO :  RATIO BETWEEN  THE CONCENTRATIONS AT THE  TOP OF THE
344.     c          CONSTANT FLUX LAYER AND  THE SURFACE.
345.     C      SOX  :  VERTICAL PROFILES OF  SOX OR NOX.
346.     C      TRR  :  FRACTION OF GAS  NOT TRANSFORMED  TO PARTICLE.
347.     C      SX  :  HOLDS NEW PROFILE OF SOX  OR NOX.
343.     C      TX  :  LENGTH OF A TIME  STEP  IN  SECONDS.
349.     C      VAR  :  THE SPATIAL VARIANCE  OF  THE POLLUTANTS'  CONCENTRATIONS.
350.     C
351 .     C
352.           COMMON /DP/ DZ< 113), A 1.23 ( 9 ), A7< 9 > ,A8( 9 ), DZD< 9 ), DZU( 9 >, DZZ< 9 )
353.           COMMON /R/ DEPT< 153 , 3 ), SOXS< 168, 3 ), SBUD< 1 S3 , 3 ), ABVML ( 1 63 , 3 ) ,21,22.,
354.          X  VK,USTAR,KZ{2,9,24>,VD(2,24,2),TR<24>,DCOR<10),CONV(4>,EMR<9,24>
355.           COMMON /INT/ NTS,NZ,NZi,SRCLVL,FAC,I POL
355.           DOUBLE PRECISION DZ, A123 ,A7 ,AS , DZD , DZU., DZZ
337.           DOUBLE PRECISION SOX ( \S3, 3 ), F 1 < 9 ), F3( 9 ) , AA , AAA , TRR
358.           DOUBLE PRECISION A12,A21,P(9,3,3 ) ,A1,A2,A3,AB,D,DD
3£9.           DOUBLE PRECISION RATIO*6),SX(9,3>,VDD,VDU,CP,VAR
3 0.0.           DOUBLE PRECISION Ml,M2,M3,DEP<3 ) ,AML<3>
361.           REAL KZ
362.           INTEGER SRCLVL,ORIGHR
3S3.     C
354.     C ... SET DUMMY VALUES FOR USE IN CALCULATING  DIFFUSION  ACROSS
365.     C     THE TOP OF THE MIXING LAYER ASSOCIATED. WITH  CLOUD  PENETRATION.
366.           SOX(NZ+1,1>=0.
357.           SOX(NZ+1,2)=0.
358.           SOX(NZ+1,3)=0.
369.     C
370.           DO  1190 ORIGHR=1,24
371.             AML(1)=0.
372.             AML(2>=0.
373.             AML(3)=0.
374.     C
375.     C ...   KIPUT UNIT EMISSIONS OF GAS AND  PRIMARY  PARTICLE INTO
375.     C       THE SOURCE LAYERS.  PUT A VERY SMALL AMOUNT OF  SOX OR  NOX INTO
377.     C       THE OTHER LAYERS TO ELIMINATE  DIVISION BY  ZERO.
373.             SOX< SRCLVL,1 ) = !.*EMR(SRCLVL,ORIGHR)
379.             SOX{SRCLVL,25 = 1.*EMR(SRCLVL,ORIGHR )
20,'j.             DO 1U10 0 = 1 ,3
331.               DO 1010 1 = 1, NZ
332.                 17 ( I.EQ.SRCLVL .AND. J.LT.3) GO TO  101.3
333.                 SCX! I , J ) = OZ( I )*1 . E-ljCT
                                       65

-------
384.
385.
385.
387.
388.
339.
390.
391 .
392.
393.
394.
395.
395.
397.
398.
399.
400.
401 .
402 .
403.
 1013
407 ,
4,19
412.
413.
414.
415.
415.
417.
413.
419.
420.
421 .
422.
423.
424.
42S.
428.

423.
429.
430.
431 .
'132.
433.
434.
435.
435.
437.
433.
439.
44£f,
441 .
442.
443 .
444.
445.
445.
4 4 7 ,
C
C
C
C
C
C
 1020
C
C
C
          1021
C
C
C
 1030
 1040
C
C  . . .
    50
    CONTINUE

PRODUCE INTEGRATION TO CORRESPOND TO THE TRAJECTORY
CALCULATIONS.
DO 1190 ITT=1,NTS
  ID=1
  IF (ITT.LT.4) ID=2

  DEP(1>=0.
  DEP<2)=0.
  DEP(3)=0.
  L=MOD((ITT+ORIGftR-2),24)+l

  COMPUTE RATIOS BETWEEN CONCENTRATIONS AT THE'TOP OF  THE
  CONSTANT FLUX LAYER AND AT THE SURFACE HEIGHT.
  DO 1020 0=1,2
    RR=VK*USTAR/VD=RR/(RR+ALOG = RATIO(5)

  COMPUTE THE LARGEST ALLOWABLE TIME INCREMENT OF  INTEGRATION
  FOR A SUBSTEP.
  TX=3600.
  DTT=TX
  RS=10000.
  DO 1030 1 = 1,NZ.
    RR=DZ/KZ(ID,I,L)
    IF (RR.LT.RS)  RS=RR

  AAA=0.
  DO 1040 0=1,3
    DO 1040  I = l,NZ.l
      SS = SOX+A2*DZ(I >
      AA=A12/A21
      IF (AA.GT.AAA)  AAA=AA
      CONTINUE

  DTT=0. i250*RS/AAA
  IF (DTT.GT.TX) DTT=TX
  NTI=TX/DTT + 0.9999
  DTT=TX/NTI

  CHECK FLUXES.  IF EITHER IS  TOO  LARGE,  REDUCE THE  TIME  STEP,
  FL2 = VD(!D,L,1}*DTT/DZ<1 )
  FL4=VD
-------
448.     C ...     IF FIRST TOUR, REDUCE  THE TIME STEP  SINCE DOING FORWARD
449.     C         INSTEAD OF CENTERED  FINITE  DIFFERENCING.
450.        60     IF UTT.EQ.l)  DTT=DTT/2.
451.     C
452.     C ...     MODIFY TRANSFORMATION  RATE  FOR SOURCES IN EOTTOM TWO LAYERS
453.     C         DURING THE FIRST  2.5 HOURS.
454.               AFAC=0.
455.               IF (ITT.GT.3  .OR.  SRCLVL.GT.2) GO TO 70
456.               AFAC=0..25*(3.-SRCLVL)'-<4.-FLOAT( ITT) )/( FAC*FLOAT( IPOD)
457.        70     TRR=1 .0-(TR< L )+AFAC* .01 )
458.     C
459.     C ...     MAKE A PORTION OF  THE  DIFFUSIVITV CALCULATIONS.
450.               DO 1080 K=2,NZ.
451.                 F1(IO=2D0*KZ< ID.K-1 ,L)*DTT
462.      1080       F3(K)=2D0*KZ< ID,K,L)*DTT
453.               Fl(l>=0.0
464.               F3(1)=F1(2>
455.               F3(NZ)=F3*DCOR< KU ) )
434.                     A3=SOX(KD,I >/< DZ( KD >*DCOR< KD ) )
405.     C
435.                     VDU=F3+A3*DZ )
433.                     CP=(VDU-VDD>*.5
',09.                     VAR=(DZZC<)  -i-  5D-r-'/A123 M3=0.
436.     C ...           FOR BOTTOM  LAYER
497.                      IF (K.EQ. 1 )  Ml=0.
493.                     M2=1.-M1-M3
499.                     P=0.
          1110           CONTINUE
         C
5J9.     C ...       REALLOCATE  DISPERSED  PUFFS.
51J.                 DO 1120  0=1 ,3
511.                   SX
-------
512.
513.
514.
515.
515.
517.
513.
519.
520.
521.
522.
523.
524.
525.
526.
527.
528.
529.
530.
531.
532.
533.
534.
535.
536.
537.
533.
539.
540.
541.
542.
543.
544.
545.
5 -IS.
547.
540.
549.
550.
551 .
532 .
553!
534.
555.
555.
557 .
553.
559.
550.
561 .
562.
563.
554.
C '^ C
566.
567.
568.
559.
570.
571 .
572.
573.
574.
575.





1120
t
/
C
C . . .






1
c
c



1130
C
C . . .


1140
1150
C
C . . .


1150





C
C . . .





1170
1130
1190
C
C . . .
c





1210


//GO.:
// UN
// SP-
//GO.:
1
        SX< 1,0)=P( U2.,0)*SOX(UO)+4>(2,1,J)*SOX{2,0)
        AML<0)=-AML.<0)+P*SOX(NZ,0 )
        DO 1120 1=2,NZ1
          11=1-1
          12=1+1
          SX 00=2
        IF 0*DTT*RATIO{ 0 + 3)*SXC1,0)
         /DZ(1 )
     ALLOW NO MORE THAN  3%  OF THE S  OR  N TO BE  DEPOSITED IN
        ONE SUBSTEP.
        IF (D.GT.(SX( 1 ,0 )*.03) )  D=SX( 1 ,0 )*0.,03.00
        DEP<0)=DE?<0)+D
        SX( 1,0) = SX<1 ,0 )-D
        CONTINUE

      TRANSFER  NEW,  VALUES TO ARRAY   SOX .
      DO 1140  1 = 1,NZ
        DO 1140 0=1,3
          SOX<1,0 > = SX(1,0)
      CONTINUE

    CALCULATE CHEMICAL TRANSFORMATION  FROM  GAS  TO PARTICLE.
    DO 1160  I=1,NZ
      SOX=TRR*SOX<1,1)
    AML(3)=AML< 3)+CONV(4)*(1.0-TRR >*AML(1)
    ABVML(ITT,3)=A3VML(ITT,3>+AML<3>
    ABVML
SUM FIELDS.
DO 1180 0=1,3
  SOXS
-------
575.
577.
570.
573.
530.
531 .
532.
533.
5C4.
505 .
535.
507.
533.
589.
590.
591 .
592.
593.
534.
595 .
596.
597.
590.
599.
60,0' .
601 .
602.
603 .
60- .
5^5.
60'6.
657 .
60S ,
509.
610.
611.
612.
613.
614.
615.
613.
617.
613.
619.
620.
62 i .
622,
623.
624.
625 .
625.
627.
623.
629.
630.
631 .
632.
633.
634.
635.
636.
6 C 7 .
630.
639.
100
1 .
10
90
5
45
3
50
16
120
1 .
10.
1 .
10.
1.
10.
0.5"
10.
1 .
10.
1 .
10.
1 .
10.
1.
10.
1.
10.
40
400
100
1000
10
65
5
35
3
40
8
50
1 .
10.
1 .
10.
1.
10.
.75
10.
1.
10.
1 .
10.
1.
10.
1 .
1 .
1 ,
1 .
24
200
50
485 .
10
70
100 100
.988 .975 .
10 10 10
85 30 75
555
45 45 45
333
55 50 45
16 16 16
110 100 90
1 . 1.
10. 10.
1. 1.
10. 10.
1. 1.
10. 10.
0.5 0.5
10. 10.
1. 1 .
10. 10.
1 . 1 .
10. 10.
1 . 1 .
10. 10.
1 . 1 .
10. 10.
1 . 1 .
10. 10.
40 40 40
380 360 300
100 100 100
950 900 750
10 10 10
55 50 45
555
36 36 30
388
35 30 25
888
55 50 45
1. 1.
10. 10.
1 . 1 .
10. 10.
1 . 1 .
10. 10.
.75 .75
10. 10.
1 . 1 .
10. 10.
1 . 1 .
10. 10.
1 . 1 .
10. 10.
1 . 1.
1 . 1 .
1 . 1 .
1 . 1 .
24 24 24
180 L60 140
50 50 60
''•50 400 350
10 10 10
60 55 50
100
.965
15
70
5
40
10
40
20
80
1.
10.
1 .
10.
1.
10.
0.5
10.
1.
10.
1 .
10.
1 .
10.
1 .
10.
1 .
10.
40
220
100
550
10
40
5
25
8
20
3
40
1 .
10.
1 .
10.
1 .
10.
.75
10.
1 .
10.
1 .
10.
1 .
10.
1 .
1 .
1 .
1 .
24
112
60
280
13
40
200 200
.950 .930 .
25
55
10
30
15
35
30
70
2.
10.
1 .
10.
1.
10.
45 55
30 15
15 25
20 10
30 45
25 10
60 90
50 20
4.
10.
1 .
10.
1 .
10.
0.S 0.5
10.
1 .
10.
i
10!
i .
10.
i .
10.
i.
10.
60
140
150
350
20
35
8
20
15
15
15
35
1 .
10.
1 .
10.
1 .
10.
.75
10.
1 .
5.
1 .
5.
1 .
3.
1 .
1 .
1 .
1 .
24
30
50
200
10
35
10.
1 .
5.
1 .
5 .
1 .
5.
1 .
5.
1 .
5.
80 140
80 60
200 3B0
200 150
30 40
20 10
12 16
10 5
25 30
10 a
30 45
20 B
2.
5.
1 .
5.
1 .
5 .
.75
5.
1 .
1 .
1 .
1 .
1 .
1 .
1.
1 .
1 .
1 .
5-6 80
55 24
140 200
1 .10 30
10 30
10 10
200
910 .
80
10
35
5
55
8
110
16
8.
6.
3.
5.
1 .
5.
0.5
5.
1 .
1.
1 .
1.
1 .
1.
1 .
1 .
1 .
1.
220
40
550
100
50
10
20
5
35
8
55
8
5.
1.
2.
1.
1 .
1 .
.75"
.75
1 .
1 .
1 .
1 .
1 .
1 .
1.
1 .
1 .
1 .
112
24
280
50
50
10
400
.830
90
10
40
5
60
8
120
15
10.
2.
7.
i .
i!
i .
3.
0.5
1.
1 .
1 .
1.
1 .
1 .
1 .
1 .
1 .
1.
300
40
T50
100
60
10
25
5
40
8
60
8
3.
1 .
5.
1.
3.
1 .
.75
.75
1 .
1 .
i .
1 .
1 .
1 .
i .
1 .
1 .
1 .
U0
24
350
60
60
10
400 LAYER TWICKNFSS 
.340 DENSITY CORRECTION
90 90
10 10
45 45
5 5
60 60
8 8
120 120
16 16
10. 10
1. 1
10. 10
1 . 1
13. 10
1 . 1
10. 10
0.5 0.
5. 10
1. 1
1. 5
1. 1
1. 5
1 . 1
1 . 1
1 . 1
1. 1
1. 1
360 3SJ0T
40 40
900 9-50
100 100
65 65
10 10
30 36
5 5
40 40
8 3
60 60
8 8
10. 10
1. 1
8. 10
1. 1
7. 10
1. 1
3. 10
S02 DRY
*
S04 DRY

NO/N02

DEP
100
DEP

S.UMMFR



DRY DEP

NITRATE/HNO3

10.
1.
. 10.
1 .
. 10.
1 .
13.
5 0.5
. 10.
1 .
. 10.
1 .
. 10.
1.
5.
1 .
1 .
1 .
S02 OXID
TIMES

10.
1.
10.
1 .
10.
1 .
10.
0.5
10.
1.
10.
1 .
10.
1 .
10.
1 .
5.
1 .

DRY DEP

DIURNAL
PATTERN OF
VERTICAL
EDDY
DIFFUSIV.
IN METERS
SQUARED
PER S
BY HR AND
LAYER








. RATE IN X/ttR
100
NOX OXLDATION









. 10.
1.
10.
1 .
. 10.
1.
10.
.75 .75 .75
1 . 3
1 . 1
1 . 1
1 . 1
1 . 1
1 . 1
1. 1
1 . 1
1 . 1
1 . 1
160 150
24 24
400 450
50 60
73 70
10 10
10.
1 .
1 .
1 .
1 .
1 .
1 .
1 .
1 .
1 .















10.
1 .
10.
1 .
10.
1 .
10.
.75
10.
1 .
10.
1 .
10.
1 .
1 .
1 .
1 .
1 .







RATE

FALL































69

-------
Appendix C:  Program CONCDEP Listing
                 70

-------
640.        ;3   5   5   5    5    5    5   12   15  18  21   24                  WINTER
641.
642.
643.
644.
5-5.
645 .
647.
643.
649.
650.
651 .
652.
653.
654.
655.
655 .
657.
658.
659.
663.
651 .
632.
663.
654 .
6G5 .
66S.
567.
5S3.       10  10  10   10   15   25   35   55   70  80  80   80                  STRING
G59 .
67.2?.
671.
672.
673.
6 7 -I.
675.
675.
6 7 7 .
673.
679.
63.1?.
o» 1 .
632.
633.
684.
6S5.
6S6.
637.
683.
639.
695.
691 .
6S2.
693.
694.
6'J5.
653.
697.
693.
;3
25
9
45
15
3£?
1.
10.
1.
10.
1 .
ISS.
.75
10.
i .
12.
I.
10.
1.
1.
1 .
1 .
1 .
1 .
16
112
40
280
10
80
5
40
a
50
a
60
1 .
10.
1.
10.
1 .
10.
.75
10.
1 .
IS.
1 .
12.
1 .
10.
i .
10.
1 .
10.
32
325
35
800
/*
5
25
9
35
15
80
1
i *
10.
1 .
10.
1 .
10.
.75
10.
1 .
10.
1 .
10.
1.
1 .
1 .
1 .
1 .
1 .
16
96
40
5 5
25 25
9 9
30 25
15 15
75 70
1.
10.
1 .
10.
1 .
10.
.75
10.
1 .
10.
1.
10.
1.
1 .
1 .
1 .
1 .
1.
15 16
30 64
40 40
240 Z00 1S0
10
75
5
40
3
45
3
55
1 .
10.
1 .
10.
1 .
10.
.75
IS.
1 .
10.
1 .
10.
1 .
10.
1 .
10.
1 .
10.
32
230 2
80
10 10
70 65
5 5
40 40
8 8
40 35
8 B
50 45
1.
10.
1 .
li?.
1 .
10.
.75
10.
1 .
10.
1 .
10.
1.
10.
1 .
15.
1 .
10.
32 32
40 200
80 80
750 500 500


5
16
9
20
15
65
1 .
10.
1 .
10.
1 .
10.
.75
10.
1 .
10.
1 .
10.
1.
1 .
1 .
1 .
1 .
1 .
16
48
40
L20
15
60
5
35
3
30
10
40
1 .
10.
1 .
10.
1 .
10.
.75
10.
1 .
10.
1 .
10.
1 .
10.
1 .
10.
1 .
10.
32
160
80
400

5
12
9
15
15
40
1.
5.
1.
5.
1 .
5.
.75
5.
1 .
1 .
1 .
1 .
1.
1 .
1 .
1 .
1 .
1.
16
32
40
80
25
45
10
25
15
25
15
35
1.5
10.
1 .
10.
1 .
10.
.75'
10.
1 .
10.
1 .
10.
1 .
10.
1 .
10.
1 .
5 .
56
1.25
140
300

5 12
10 5
9 30
9 9
15 45
20 15
1 .
1 .
1.
1.
1 .
1.
.75 .
.75 .
1 .
1 .
1 .
1 .
1.
1 .
1 .
1 .
1 .
1 .
16 32
16 16
40 80
40 40
35 55
30 15
15 20
15 10
25 35
20 10
30 45
25 10
3.
10.
1 .5
10.
1 .
10.
.75 .
10.
1 .
10.
1 .
10.
1.
10.
1.
5.
1 .
i
80 12.B
80 56
200 300
200 140

15
5
35
9
70
15
1 .
1.
1.
1.
1.
1 .
75
75
1 .
1 .
1 .
1 .
1 .
1.
1 .
1 .
1.
1 .
48
16
120
40
70
10
30
5
45
3
55
3
7.
5.
4.
5.
3.
5.
75
5.
1 .
5.
1 .
5.
1.
5.
1.
1.
1 .
1 .
18
5
40
9
80
15
5 .
I!
i .
i.
i .
i .
.75
.75
1 .
1 .
1 .
1.
1.
1 .
1 .
1.
1 .
1 .
64
16
L60
40
80
10
35
5
50
8
60
8
1.0.
1 .
3.
1 .
6.
1 .
5.
.75
1 .
1 .
1 .
1 .
1 .
1 .
1 .
1 .
1 .
1 .
1 60'. 200
32
400
80

32
•500
80

21
5
45
9
80
15
10.
1.
5.
1 .
3.
1 .
.75
.75
1 .
1.
1.
1.
1.
1 .
1 .
1 .
1 .
1.
80
16
2JST0
40
80
10
40
5
50
3
60
3
10.
1 .
10.
1 .
10.
1 .
10.
.75
5.
1 .
1 .
1 .
1.
1 .
1 .
1.
1 .
1 .
240
32
500
80

24
5
45
9
80
15
10.
1.
10.
1.
10.
1 .
10.
.75
1 .
1 .
1 .
1 .
1.
1 .
1 .
1 .
1.
1 .
9-6
16
240
40
80
10
40
5
50
3
60
8
10.
1 .
10.
1 .
10.
1 .
10.
.75
10.
1 .
5.
1 .
5.
1 .
1 .
1 .
1 .
1 .
28J2T
32
7'00
80







10.
1.
L0r.
1 .
1J3.
1.
10.
.75
10.
1 .
1 .
1 .
1.
1.
1 .
1.
1.
1.












LSr.
l.
10.
1 .
10.
l .
10.
.75
10.
1 .
10.
1 .
10.
1.
5.
1.
i[
I .











10.
1.
10.
1 .
10.
1 .
10.
.75
10.
1 .
10.
1 .
1 .
1.
1 .
1 .
1.
1 .












10.
1 .
10.
1.
10.
i *
10.
.75
10.
1 .
10.
1 .
10.
1 .
10.
1.
5.
1 .





                                        71

-------
 3.
 9.
10.
11 .
12.
13.
14.
15.
15.
17.
33,
O *v ,
C 5
3 5 ,

33.
42 .
44
//CONCDEP JOB
//*VAIN ORG=RM080
// EXEC FGICLG
//SYS IN DD *
      DIMENSION
      DInZMSIOM
      DIMENSION
      DIMENSION
      DIME MS I ON
      DIMENSION
      DIMENSION
      REAL MAXM
      R0 = 5397./381 .
      R018=l.866*RO
      P2=3.14159/2.
                          = 1,REGION=5J2TJ3TK,CLASS=W
C
C

C
C
C
        C
        C
        C
        C
        C
                                                                  ,45)
                CO(6,28,123),WD(6,28,123),LANDE< 19,16)
                S02<51 ,45), 204(51 ,45>,DDEP(51,45),WDEP(51,
                X(51 ),Y(45), INDO(51 ,45)
                SOX<28,3,6),DEPT{28,3,6),S8UD(28,3,6),MAXN( 12)
                EM(720,6),XS{720),YS(720),WF(6>,EW( 12)
                ISTPR<60>, IND(51 ,45)
                AREA(60),AREAA(60),ADRY(60),AWET(60)
      NSP
      MSP =
         333:
                  IS
                  12
THE NUMBER OF TELESCOPED TIME STEPS,
              SBXX AMD SBYY ARE THE VARIANCES OF A UNIFORM DISTRIBUTION ACROSS
               A 1/3 NMC
              SBXX=(1,/3.
              SBYY = ( 1 ./3
                 CELL
                 )*( 1 ,
                 )*( 1
                               AND ARE USED IN RESOLUTION LIMITATIONS,
      C!<1 AND OUTSID ARE USED IN REDISTRIBUTING  PUFF MASS  BEYOND  .
       2.5 STANDARD DEVIATIONS.
      OUTSID=1.045
      CK1=0.0175*0.3939*2.0*3.14159
      CK1=ALOG(CK1 )

      FAC1 AMD FAC2 ARE USED  IN CALCULATING CONCENTRATION  AND  DEPOSITION
       FIELDS.
      FAC1=1.E06*OUTSID
      FAC2=365."OUTSID
      UTROP IS THE FRACTION OF SULFUR MASS  REMOVED
       THE PBL 3Y PRECIPITATION TRANSFERRED TO  THE
        TROPOSPHERE.
      UTROP=l./2.
      UT = UTROP/(1.-UTROP )
                                                    FROM
                                                    FREE
      INPUT PARAMETERS
      NXC=51
      NYC = .-t5
      X<1)=-10.83
      DXC=1 ./3.
                               DESCRIBING CONCENTRATION/DEPOSITION  GRIDS.
      V ( 1 ) = 1 9 . S 3 3
      DYC=1./3.
      YS=Y< 1 ) + DYC/2.0
      DO  2333 1=2,NXC
      X( I > = X! 1-1 J + DXC
      DO  3332 0=2,NYC
      YCO ) = Y(0-1 )-DYC

      FACX=38l.*381.*2.*3.1416/(1.866*1.366)
      FACXX = 381 .-'381 .*DXC*DYC/( 1 .366*1 .865)

      DEFINE PARAMETERS DESCRIBING  TRAJECTORY  PROGRAM FIELDS.
      NG=19
      MG=>5
      X L = - 1 1 . 5
      YL=20.5
                                        72

-------
 64.            GX=1.0
 S 5 .            G Y = 1 .0
 IS.            NTS=2S

 63.            DO 901 1=1,NXC
 69.            DO 901 0=1,NYC
 70.            S02*I,0)=0.0
 71.            S04(I,0)=0.0
 72.            DDEP
 84.       7003 FORMAT*1812)
 S5.       7002 FORMAT*33I2)
 CJ.      C
 S7.      C ... INDD  IDENTIFIES EACH STATE AND PROVINCE.
 £3.            DO 34.J 00 = 1 ,45
 S3.            0=46-00
 20.        940 READ  (50,7005) *1.51*1.61*10.0
         C
iZ'7.      C ... MAXN  IS  USED IN CALCULATING THE MAXIMUM POTENTIAL TRAJECTORIES
103.      C      CONTRIBUTING TO A STATISTIC IN THE TELESCOPED TIME STEPS.
l.-:9.            READ  8020, MAXN
lit/.       3020 FORMAT* 12F6. 1 )
111.      C
112.      C ... EW IS  USED TO WEIGHT PUFFS IN THE TELESCOPED TIME STEPS.
113.            READ  8021,EW
114.       8021 FORMAT*12F2.0)
115.      C
ue.      c
117.      C ... READ  RESULTS OF TRAOECTORY PROGRAM.
113.      C     IN THE TRAOECTORY PROGRAM ARRAY CO IS CALLED CONC.
113.            READ  <10)  CO,WD,LANDE,NF
12.0'.            XNF = NF
12..            DO 700 0 = 1, 123
122.            DO 701 1=1,23
i:3.            CO*6, I ,0)=CO<5,I,0)/*CO*3,1,0)*CO(d,!,J) )
124.        701 WD(5,I,0)=WO(6,I,0)/*WD(3,I,0 )*WD*4,I,0) )
125.      C
125.      C
127.      C
                                       73

-------
12S.      C
129.      C ... POST-PROCE'SSING OF TRAJECTORY STATISTICS TO COMPRESS TIMESTEPS,
1?7.            DO 500 K=3,4
131.            DO 502 1=1,4
132.            IA=2*I+4
133.            IB=IA-1
134.            IC=I+4
135.            W1=CO<5,IA,0 )
135.            W2 = CO<5,IB,J )
137.            WA=W1+W2
133.            IF (WA.LE.0. ) WA=1.
139.            X1=UD(5,IA,J)
140.            X2=WD<5,IB,J)
141.            XA=X1+X2
142.            IF (XA.LE.0.) XA=1.
143.            Y1=CO(K-2,IA,0)
1A4.            Y2 = CO(K-2,IB,0 )
145.            Z1=WD(K-2,IA,0)
145.            Z2=WD(K-2,IB,J)
147.            YY=\W1*Y1+W2*Y2)/WA
140.            ZZ=(X1*ZUX2«Z2)/XA
143.            CC=(W1*CO(K,IA.O>*CO(K,IA,0)+W2*CO(K,IB,0)*CO* + X2*(Z2-ZZ>*(Z2-ZZ»/XA
153.            C0(K,IC,0 ) = SQRT(CC>
154.            WD(K,IC,0>=SQRT(DD)
155.        5J?2 CONTINUE
156.            DO 503 1=1,4
137.            IA=4*I+12
lie.            IB=IA-1
               IC=IA-2
               ID=IA-3
               IE=I+8
               W1=CO(5,IA,0 )
153.            W2 = CO<5,18,0 )
164.            W3=CO(5,IC,J)
165.            W4 = CO(5,ID,J )
166.            WA=W1+W2+W3+W4
157.            IF (WA.LE.0. ) WA=1.
153.            X1=WD<5,IA,0>
:.-;.            X2=WD(5, IB.O )
i7,.            X3=WD(5,IC,J)
1"'..            X4=WD<5, ID,0 )
1/2.            XA=:;i-x2-;:3-M4
173.            IF < XA.LZ.Cf. )  XA=1 .
174.            Yl=CO(K-2,IA,0)
175.            Y2 = CO(!<-2, IB,0 )
175.            Y3 = CO(K-2,IC,0 )
177.            Y4=CO(K-2,ID,0>
178.            Zl=WD(K-2,IA,0>
179.            Z2=WD(K-2,IB,J>
180.            Z3=WD
181.            Z4=WD( !<-2, ID,J >
102.            YY = (Wl*Yl+W2*Y2-s-W3*Y3-t-W4*Y4)/WA
133.            ZZ=< X1*Z1+X2*Z2 + X3*Z3 + X4*Z4>/XA
1S4.            CC = (Wl-*CO(K,IA,0 ) + ( Yl-YY )*( Y 1-YY » +
1G5.           1 W2-«CO + ( Y2-YY >*( Y2-YY ) ) +
1SS.           2 W3*(CO(I<:, IC,0 )?;CO(K, IC,0 > + ( Y3-YY)*( Y3-YY) > +
187.           3 W4*/WA
133.            DD = < Xla(WD*(Zl-ZZ ) )+
1S9.           1 X2~(V/D(K, IB,0 )*WD(K, IB,0 ) + ( 22-ZZ )'-; { Z2-ZZ ) ) +
190.           2 X3*(WD(K,IC,0)-WD(K,IC,0)-K Z3-ZZ)*(Z3-ZZ ) ) +
191 .           3 X4-'/XA
                                        74

-------
192.
193.
194.
135.
1S6.
197.
193.
199.
200.
201 .
202.
203.
204.
205.
206.
207.
208.
259.
210.
21 1 .
212.
213.
214.
215.
216.
217.
213.
21S.
220.
221 .
222.
223.
224.
225.
226.
227.
223.
229.
230.
231 .
232.
233 .
234.
235.
235.
237.
233.
239 .
240.
241.
242.
243.
2 -' 4 .
245 .
246.
247.
243.
249.
250.
251 .
252.
253.
254.
255.
CO(K, 11, J )=SQRT
503 CONTINUE
500 CONTINUE
DO 710 !<=! ,5
IF (K.EQ.3) GO TO 7
IF (K.EQ.4) GO TO 7
DO 702 1=1,4
IA=2*I+4
IB=IA-1
IC=I+4
W1=CO(5, IA,0 )
W2=CO<5, IB,0 )
WA=W1+W2
X1=WD(5, IA,0 )
X2=WD(5, IB,0 )
XA=X1+X2
IF (WA.LE.0. ) WA=1 .
IF ( XA.LE .0. ) XA=1 .
C0( K, IC,0 )=(W1*CO
W4 = CO(5, ID,0 )
WA=W1+W2+W3+W4
X1=WD(5, IA,0 )
X2=WD(5, IB,0 )
X3=WD<5, IC,0)
X4=WD<5, ID,0 >
XA=X1+X2+X3+X4
IF ( XA.LE. 0. ) XA=1 .
IF (WA.LE.0. ) WA=1 .
CO(K, IE,0 >= C0<5, I ,0 ) = -0.99
, I ,0 >*CO(3, I ,0 )+2.*SBXX)
, I ,0 )*CO< 4,1,0 )+2.*SBYY)
, I ,0 )*WD< 3,1,0 )+3.*SBXX )
, I ,0 )*WD< 4,1,0 ) + 3.*SBYY>





C ... READ RESULTS OF VERTICAL DISPERSION PROGRAM.
C NVI IS THE NUMBE
C VERTICAL INTEGRA
DO 10 0=1,5
R OF EFFECTIVE EMISSION LEVELS FOR WHICH
TIONS ARE MADE.

75

-------
256.
257.
258.
259.
260.
261.
262.
263.
254.
265.
266.
267.
268.
269.
270.
271 .
272.
273.
274.
275.
276.
277 .
278.
279.
230.
281.
282.
283.
284.
285.
286.
287.
288.
289.
290.
291 .
292.
293.
294.
295.
296.
297.
298.
299.
300.
301 .
302.
303.
304.
305.
306.
307.
308.
309.
310.
311 .
312.
313.
314.
315.
316.
317.
318.
319.




C
C
C
C
C
C





















C
C
C
C

C







C
C
C
C
C
C









C
C
C
C
C
READ (20) «SOX(I,K,0),I = 1,NTS>,K=1,3),«DEPT(I,K,0),I = 1,NTS>,
# K=1,3),«SBUD,K=1 ,3)
10 CONTINUE
REWIND 22




... POST-PROCESSING OF VERTICAL INTEGRATION STATISTICS TO COMPRESS
TIMESTEPS.
DO 706 K=1,NVI
DO 706 J=l,3
DO 707 1=1 ,4
IA=2*I+4
IB=IA-1
IC=I+4
SOX( IC,J,K)=(SQX< IA,J,K)+SOX< IB,J,K)>/2.
DEPT +DEPT( IB,J,K) >/2.
707 SBUD< IC,0,K) = (SGUD< IA ,0 , K)+SBUD< IB,J,K) 5/2.
DO 708 1=1 ,4
IA=4*I+12
IB=IA-1
IC=IA-2
ID=IA-3
IE=I+8
SOX( IE,J,K>=(SOX( IA,0,IO + SOX< I2,J,K)+SOX( 1C , J , K >+SOX( ID ,0 , K ) )/4.
DEPT+DEPT< 1C , 0 , K >+DEPT< ID,J,K»
S /4 .
708 SBUD + SBUD(ID,0,!<))
# /4 .
706 CONTINUE




GO TO 4444
... IS DEFINES THE APPROPRIATE SEASONAL EMISSION FIELD.
DO 4460 L = l , IS
READ (30) NS,NZ,( (EM< I ,0 >, 1 = 1 ,NS),0 = 1 ,NZ),(XS( I ),
1 I=1,NS),(YS(I),I=1,MS),DX,DY
4460 CONTINUE
DO 4488 1=1, MS
XS( I ) = XS< I ) + DX
4488 YS< I ) = YS( I )-DY




... THE GO TO 4444 STATEMENT AND THE FOLLOWING 11 STATEMENTS
... AND COMMENT CARDS SHOULD BE DELETED IN REAL SIMULATIONS.
4444 NS=1
XS( 1 ) = -4.5
YS< 1 ) = 15.5
EM< 1,5)=250.0
EM< 1 , 1 )=0.0
EM< 1 ,2)=0.ff
EM( 1 ,3)=0.0
EM( 1 ,4 >=0.0
EM( 1 ,6)=0.0
... EMISSION SOURCE IS LOCATED IN EASTERN OKLAHOMA.
... EMISSION RATE IS 250 KT S02 PER SEASON.



76

-------
320.     C
321 .     C
322.           SUMM=0.0
323.           DO 771 1=1,NS
324.           DO 771 0=1,NVI
325.           SUMM=SUMM+EM<1,0 )/365.
326.       771 EM(I,0>=1.E03*EM- GO TO 9999
343.           IF (Il.GT.MG) GO TO 9999
344.           IF (Ol.LT.l ) GO TO 9999
345.           IF (Ol.GT.MG) GO TO 9999
346.           IF (LANDE(II,015.EQ.0) GO TO 779
347.           GO TO 773
348.       779 13=11-1
349.           14=11+1
350.           03=01-1
351.           04=01+1
352.           IF ( I3.LT.1 ) 13 = 1
353.           IF (I4.GT.NG) I4=NG
354.           IF (03.LT.1 ) 03=1
355.           IF (04.GT.MG) 04=MG
356.           DO 777 i<=13, 14
357.           DO 777 L=03,04
358.           IF .GT.0) GO TO 7708
359.       777 CONTINUE
360.           GO TO 9999
361.      7708 I1=K
362.           01=L
363.       778 CONTINUE
364.           XBI=XS- )
355.           YBI=YS(IO)-
356.           I=LANDE( II ,J1 )
337.      8000 FORMAT('  TOTAL EMISSIONS = ',F8.1,' KT S PER SIMULATION1)
o < p      p
369.     C ... PMF IS THE  FRACTION OF TOTAL SULFUR EMITTED AS
370.     C      PRIMARY SULFATE.  IT IS HIGHER IN THE NORTHEAST, FLORIDA, AND
371.     C      IN THE  FIRST EMISSION LAYER, UNDER THE ASSUMPTION
372.     C      THAT MOST  EMISSIONS THERE ARE FROM OIL COMBUSTION.
373.           PMF=0.015
374.           IF (XS( 10 KGT.0.5) PMF = 1.33
375.           IF (YS< 10 ).GT. 17.5 .AnO. XS( I 0>.GT.-3.0 ) PMF=0.03
376 .     C
377.           DO 86  M=l,NVI
373.        86 WF(M)=0.£f
379.     C
3S0.     C ... DO-LOOP  FOR THE  TELESCOPED TIME STEPS.
381.           DO 100 0=1,NSP
382.           00=0-1
383.           IF (O.EQ.1 ) 00 = 1
                                        77

-------
334.     C
385.     C ... CALCULATE MASS WEIGHTING  OF  DRY  AND WET PUFFS.
325.           TW=CO(5,0,I>/  GO  TO  95
396.           PMQ=PMF
397.           IF (M.EQ.l)  PMQ=0.05
398.           EMQ=1.0-PMQ
393.           PMQ=1.5*PMQ
4.;/0.           V2=V2+EW(0)*EM(IJ,M)*SOX(0,1,M)*EMQ
401 .           V4=V4+EW< J )*EM< 10 , M ) * < SOX < 0 , 2 , M ) *PMCH-SOX ( 0 , 3 , M >*EMQ )
4.J2.           IF (O.EQ. 1 )  GO TO  995
403.           DO 85 L=l,00
4.74.        85 WF(M>=WF(M>+(UT*WD(5,L,I) *(EMQ*(DEPT(0J,1,M>
^75.          1 + DEPT(00,3,M)) + PMQ-/( UT*WD < 5 , L , I )+CO(5,L, I ) >>
<;75.           WETF=WETF + EW< J )•'£;•!( 10 , M )*< WF < M )-M SBUD< J J , 1 , M > + SBUD< 0 0 . 3 , M > >*EMQ
407.          # + PMQ*SBUD(00,2,M>)
458.           WQ = UTi'WD(5,JO, I >+CO(5,00, I >
409.           IF (WQ.LE.0.> GO  TO  996
410.           WF < M)=WF < M)*(1.0-WD(5,0,1)/WQ >
411.           GO TO 996
412.       995 WETF=WETF+0.5*EW<0)*EM*EM
416.        95 CONTINUE
417.     C
418.     C ... ADD BIAS  BETWEEN  TRAJECTORY  VIRTUAL SOURCE AND EMISSION  CELL TO
419.     C      PUFF CENTROID, AND  DETERMINE RECEPTOR GRID  RECTANGLE CONTAINING
420.     C      2.5 STANDARD  DEVIATION PUFF FOR CONCENTRATIONS AMD DRY  DEP.
421 .           IF (TV/.LE.0. > GO  TO  601
422.           XB = XBI+CO< 1,0,1 )
423.           YA = YBI-i-CO( 2,0,1)
''•2-\.           SCALE=CO< 3,0,1 )-vCO(4,0, I )*FACX
•!i5.           R01=SG!rf{ 1 .-C0< 5,0,1 )-CO( 5,0,1))
425.           R02 = -2.'-'R01"ROi
427.           ROW=-2.*CO<6,0,I>
423.           R01=R01*SCALE
429.           TWR=TW/R01
43,0'.           I1 = (XB-2.5"CO< 3,0, I )-X0)/DXC + 1.
431.           !2 = (XB + 2.5*CO(3,0,I )-X0)/DXC + 1.
432.           IF { II .LT.1 )  11 = 1
433.           IF ( I2.GT.NXO I2 = NXC
434.           01 = ( Y0-YA-2.5*CO<4,0, I ) >/DYC +• 1.
435.           02=
-------
443.           PH!=P2-2.*A7AN-t-Y(L)*Y )/R013)
449.           SIG=1.0+SIN-(RANFt-1 5-0.5 )*DYC
434.           DLX=ABS(X5-XB )/CO(3,J,I )
455.           DLY=ABS-X0)/DXC + 1.0
491.           IF (II.LT.1 )  11 = 1
<192.           IF (I2.GT.NXC)  I2 = NXC
493.           01 = (Y0-B-2.5*WD(4,0,I) )/DYC + 1.0
434.           02=(Y0-B+2.5*WD(4,0,I))/DYC + 1.0
495.           IF (01.LT.1>  01=1
495.           IF (02.GT.NYC)  02=NYC
497.           DO 92 K=I1 , 12
49S.           DO 91 L=01 .02
4C9.           IF (IND(K,L).EQ.0) GO  TO 91
5.-0.           MI = 1
5 " 1 .           F I = 1 .0
5,:2.           IF (O.GT.l ) GO  TO 54
53.           N I = 13
5.-'.           FI=0.1
5.5.        64 CONTINUE
: .:.           ?HI=?2-2.*ATAN(SQRT(X
-------
512.      C
513.      C
5.4.      C
515.      C ... FOLLOWING 2 STATEMENTS SHOULD Be DELETED AND THE 2 STATEMENTS
516.      C ... FOLLOWING THEM SHOULD BE ACTIVATED  IN REAL  SIMULATIONS.
517.            X5=X(K>
518 .            Y5 = Y( L )
519.      C     X5 = X< K)-(RANF(-1 )-0.5)*DXC
520.      C     Y5 = Y(L )-(RANF(-l )-0.5)*DYC
521.            WLX=ABS(X5-A)/WD(3,0,I )
522.            WLY=ABS< V5-B >/WD( 4 , J , I )
523.            B B = < WL X * WL X +WRCW*WL X *WL Y + WL Y * WL Y > /WR2
524.            IF (BB.LT.CK1) GO TO 93
525.            WCON=FI*WTWW*EXP=WDEP(K,L >+WCON*WETF
527.         93 CONTINUE
5:3.         91 CONTINUE
529.         92 CONTINUE
S3£T.        600 CONTINUE
531 .        100 CONTINUE
532.       9999 CONTINUE
533.            DO 9230 1=1, NXC
534.            DO 9233 0 = 1 ,NYC
535.            S02(  I ,0 ) = S02( I ,0 )*FAC1
536.            S04(  I ,0 )=S04< I ,0 )*FAC1
537.            DDEP< 1,0 )=DDEP( 1,0 >*FAC2
53S.       9200 WDEPf I ,0 )=WDEP( I ,0 )*FAC2
329.            PRINT 20S0
2 ::.            PRKIT 2373
5,;.       2053 FORMATi i H 1 ,//. ' TABLE 6:  AVERAGE CONCENTRATION OF'
              1  , '  S02 IN MICROGRAMS PER CUBIC  METER1)
               DO 173 00 = 1 ,MYC
           170 PRINT 20S1 ,0,(S02( 1,0 ), 1=1 1,50)
          2061 FORMAT( !3, IX, 40F3.0)
G<7.      2070 FORMAT;'      11 12 13  14  15  16  17  13  19  23  21  22  23  24  25  26  27 23
543.          S29 30 31  32  33 34 35  36  37  38  39  40  41  42  43  44  45  46  47  48  49 50
549.          % '  >
550.           PRINT 2062
551.           PRINT 2070
552.      2062 FORMAT; i HI ,//,' TASLE T.  AVERAGE CONCENTRATION  OF •
553.          2  ,' S04 IN MICROGRAMS PER CUBIC METER')
5=J.           DO 171 00=!,NVC
dc5,           0 = M\'C + 1-J3
5 = 0.       1/1 PRI;;T 2051 ,o , i S04( i ,o ) , 1 = 11 ,50)
5S7.           PRI-IT 2063
5?8.           PRINT 2070
5E9.      2063 FORMAT( 1H1 ,/, ' TABLE  8:  CUMULATIVE  DRY  DEPOSITION  OF'
              3  ,' TOTAL SULFUR  IN GRAMS PER  SQUARE  METER SIMULATION1)
               DO 172 00 = 1 ,NYC
               0=NYC+1-00
553.       172 PRINT 2064 , 0 , < DDEP < I , 0 > , I = 1 1 , 50 )
5G4.      2064 FORMAT< 13, 1X.40F3.2)
555.           PRINT 2065
555.           PRINT 2070
557.      2065 FORMAT; ii-ii ,/,' TABLE  9:  CUMULATIVE  WET DEPOSITION  OF'
5S3.          4  ,' TOTAL SULFUR  IN GRAMS PER  SQUARE  METER SIMULATION1)
Ec9.           DO 173 00=1 ,NYC
570.           0=NYC+1-00
571.       173 PRINT 2064 . 0 , ( WDEP ( I , 0 ) , I = 1 1 , 50 )
572.           COMDRY=0.0
573.           CONWET=0.0
574.           DO 200 1 = 1 , NXC
573.           DO 200 0 = 1 , MYC
                                       80

-------
S7S.           IF ( INDD< !,OKIE.JEM GO TO  283
577.           PHI=P2-2.*ATAN< < XU >*X< I > + Y< 0 >*Y< 0 »**0. 5/R018 >
578.           SIG=1.+SIN(PHI>
579.           SF=SIG*SIG*FACXX
5C0.           IQ=INDD(I,J)
531.           IF ( IQ.EQ.0) GO TO  2.882
582.           AREAA.LE.0.) GO TO  921
590.           ADRY*1.E-3*AREA(I)/AREAA(I>
591.           AWET*1.E-3*AREA,AWET(I>,! = !,60)
613.      2044 FORMAT(10X,A4,2F12.3 )
314.     C
315.           STOP
615.           END
j.T'.     //GO.FT10F7J1  DD DSN = B2545S . TSUMS0HC , D I SP=OLD
£i3.     //GO.FT50F001  DD DSN=B25458.ITC.DATA,DISP=OLD
3:9.     //GO.FT20F001  DD DSN=B2546S.VSUMMER,DISP=OLD
S:.J.     //GC.FT30F001  DD DSN = 325468 . EMI 980SU . AU84 , D I SP =OLD
62i.     //GO.FT£r9F;j.ai  DD D3n = B25468 . CONDTE ST , D I SP = < NEW, CATLG ),
622.     // UNIT=LONGSHRD,DC3=,
523.     // SPACE = (TRK,{10,5),RLSE )
624.     //GO.SYSIN  DD *
325.       3   IS IDENTIFIES  EMISSION SEASON  ( 1=W,2 = SP,3 = SU,4 = F )
626.
627.
628.
629 .
o.. o ,
c"o
537
623.
639 ,
3 IS IDENTIF1
0
33333
33331133
22111111
22111111
22111111
2 2 2 1 1 1 i i
22211111
22211111
22211111
EES EMISSION
3 3
3 3 1
3 3 1
3311
3311
33331 1
311111
i 1 1 1 i 1
111111
111111
1 i i i 1 1
1 i 1 1 1 1
111111
111111
SEASON ( 1=W
3 3
1 3 3
1 3 3
1 3 3
1 3 3
0 0
3 3
3 3
1 3 3
1133
11133
111133
111113
111133
                                       81

-------
540.
641 .
542 .
643.
644.
6^5.
646.
647.
5 i>. 8 .
649.
222111111
222111111
222211111
222211111
222111111
222111111
222111111
279711111
£££11111
277711111
£££lllii
279711111
£££11111
277711111
£££11111

C £ C £""£"1 — 1-1

2 — 7 — 7 — 7 — 7 — 7—1 — 1 — 1
£££££i 1 1
2-2-2-2-2-2-2-1-1
2—7—7—7—7—7-7—7— 1
£££u£££ I
.7—7-7—7—7—7—7-7—7
£ £ £ £ £ £ £ £ £
•2-2-2-2-2-2-2-2-2
,7-7 — 7— 9 — 9 — 7 — 9 — 9 — o
£ £ £ £ £ £ £ £ t.
2 — 9 — ""^ •— 7 — 7 — ""* — 7 — 7 — 7
£££££_ £ £ £
•2-2-2-2-2-2-2-2-2
•2-2-2-2-2-2-2-2-2
•2-2-2-2 33333
•2-2 333
•233
3 3
0
0
S3
0
1
1
1
1
1
1
1
1
1
1
i



" 1 -

— 1 —
1
-1-
— 1 —
1
— 1 —
-2


3
2
3







1
1
1
1
1
1
1
1
1
1
1






1
1
1
1
1
3

o
O
3
3








1 1 1
1 1 1
1 1 1
1 1 1
1 1-1
1 1 - 1 -
1 1-1-

1—1—1






— 1 — 1 — 1
ill
-1-1-1
— 1 — 1 — 1
ill
— 1 77
i O %3
333












i
1
1
1
•i-
- 1 -

i







3-

o
o
3
0

3-
0








1
1
1
i
1
1
1

1




i


1


3

i
i
i
3








1
1
1
1
1
1
-1

1







"I

1
-1


1
-1
3








1 I
i i
i i
i i
i i
-i i

— 1—1







-1-1
- 1 — 1
11
-1-1

— 1 — 1
-1-1
3 3








3
3
1
1
1
i
1

_ 1







-1

1
-1

1
3
3








r>
o
3
3 3
1 3
1 1
1 1
I 1
1 1
i 1
1. 1




-1-1


-1-1

11
-1-1

1 — 1
— 1 — 1
11
-1 3
3 3








3
3
3
3
1
i
1
1
1







-1

1
-1
1

3









3
3
3
. 1
1
_ 1







-1

1
-1












_ 1
i
— 1
i
_ 1







-1

i
-1













. 1
1
— 1







-1

i
3
M
J











— 1
i
— 1





0
w
-1

i
3














1




3 —


-1

i
3













_
^
o






3

















^
  INDICATES EASTERN U.S.A.
0



3
3
3
3
3
3
3
3
3





3
3
3
3
2
2
2
3
o
J
3
3
3


3
3
2
2
2
2
2
2
2
2
9
3
3
3
3
3
2
2
2
2
2
2
2
2
2
2
2
2
3
3
2
2
2
2
2
2
2
2
2
2
2
2
2
3
2
2
2
2
O
£
2
2
2
2
2
2
2
2
2
2
3
3
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
3
2
2
2
2
2
2
2
9
2
2
2
2
2
2
2
2
3
3
3
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
3
3
3
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
3
3
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
3
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
3
3
3
3
3
2
2
2
2
2
2
2
2
2
O
£.
2
O
2
2
2
9
2
2
--)
2
3
3
3
3
3
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
7
2
9
2
3
3
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
3
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
9
2-
3
3
3
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
3
3
3
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
2
-2
                                       82

-------
704,
725.
705.
737 .
7Z3.
739.
710.
711 .
712.
713.
714.
715.
716.
717.
718.
719.
720.
721 .
722.
723.
724.
725.
72G.
727.
728.
  5160 11390
  8360  5640
  5820  8410
  4960  5270
  4220 26730
 25530 35630
  AL  AZ  AR
  MD  MA  MI
  PA  RI  SC
NFLD  NS ONT
    00    10
 111122
/*
322 2-2-2
3-2-2-2-2-2'
3-2-2-2-2-2
3-2-2-2-2-2-
3-2-2-2-2-2'
3 3-2-2-2-2'
  3-2-2-2-2-
  3 3-2-2-2
    3 3 2-2-
    3 3 2-2-
    3 3 2-2-
    3322-
 5310 15870
 3630  5630
       5970
       4120
        961
       2S00
     CO  CT
     MS  MO
 4770
 7070
 8490
24650
 CA
 MN
 SD  TN  TX
PEI QUESASK
   20    30
224444
•2-2-2-2-2-2
•2-2-2-2-2-2
•2-2-2-2-2-2
•2-2-2-2-2-2
•2-2-2-2-2-2
•2-2-2-2-2-2
•2-2-2-2-2-2
-2-2-2-2-2  3
•2-2-2-2  3  3
•2-2-2  3  3
•2-2  3  3
•2-2  3  3
 10420    501
  3230  4040
 14710  7720
  6990  9700
  4080  6820
 15540  2110
  DE  DC  FL
  MT  MB  NV
  UT  VT  VA
                                      (2) INDICATES WESTERN U.S.A.
                                      (3> INDICATES NEAR-SHORE  REGION
                                          NECESSARY FOR PLOTTING
  206    07
 4850  3320
11050   930
 4530   121
 2420  5610
41260   218
 GA  ID  IL
 NH  NO  NM
 WA  WV  WI
     5860  7300
     1060   825
      784 12170
     3110  7700
     9790  9470
    59490 21570
     IN  IA  KS
     NY  NC  ND
     WYG LKALTA
           STATE/PROV.
           AREA IN TENS
           OF SQ.  MILES
               45
          65
   85
105
135
175
KY  LA
OH  OK
BC MAM

 215
 ME
 OR
 NB

255
                                       83

-------
                                   TECHNICAL REPORT DATA
                           (Please read Instructions on the reverse before completing)
 REPORT NO.
                             2.
                                                           3. RECIPIENT'S ACCESSION NO.
 TITLE AND SUBTITLE
                                                           5. REPORT DATE
USER'S GUIDE FOR THE ADVANCED STATISTICAL  TRAJECTORY
AIR POLLUTION (ASTRAP)  MODEL
             6. PERFORMING ORGANIZATION CODE
 AUTHOR(S)

 Jack D. Shannon
                                                           8. PERFORMING ORGANIZATION REPORT NO.
 PERFORMING ORGANIZATION NAME AND ADDRESS
 Environmental Research Division
 Argonne National Laboratory
 Argonne, IL  60439
             10. PROGRAM ELEMENT NO.
                CCVN1A/05-3102 (FY-85)
             11. CONTRACT/GRANT NO.

               DW930060-01
2. SPONSORING AGENCY NAME ANO ADDRESS
Atmospheric Sciences Research Laboratory - RTF ,NC
Office of Research  and Development
U.S.  Environmental  Protection Agency
Research Triangle Park, NC  27711
             13. TYPE OF REPORT AND PERIOD COVERED
                Final   (10/82 - 6/85)	
             14. SPONSORING AGENCY CODE
                 EPA/600/09
5. SUPPLEMENTARY NOTES
6. ABSTRACT
      The Advanced  Statistical Trajectory Regional Air Pollution (ASTRAP) model
 simulates long-range,  long-term transport and deposition of air pollutants,  primarily
 oxides of sulfur and nitrogen.  The  ASTRAP model is designed to combine ease of
 exercise with an appropriate detail  of physical processes for assessment applications
 related to acid deposition.

      The theoretical basis and the computational structure of the ASTRAP model are
 described.  The model is evaluated with observed data,  and an illustrative exercise
 is provided.  Computer codes of the  model subprograms are provided in  appendices.
7.
                                KEY WORDS ANO DOCUMENT ANALYSIS
                  DESCRIPTORS
                                              b.lDENTIFIERS/OPEN ENDED TERMS
                              COSATI Field/Group
18. DISTRIBUTION STATEMENT

     RELEASE TO  PUBLIC
19. SECURITY CLASS (This Report I
    UNCLASSIFIED
                                                                          21. NO. OF PAGES
20. SECURITY CLASS (This page)
    UNCLASSIFIED
22. PRICE
EPA Form 2220-1 (R«v. 4-77)   PREVIOUS EDITION is OBSOLETE

-------
-------
                                  K,a a a
                                 >aaaaasiaaaaaaaaQsi3iaa«iaaaaaaa3(aaaasi3iaaaaaaaaa«i«>

                            s'aaaaaaaaaaa
  uiaaaaaaaaasiaaaaaaaaaaaaaaaasiaaaaeiaaasiaaaaaaaaisi

                         ' a a a « a a a a a a a a a


                         'a'aa'aaaaasiasiaa
                         iaaaaaaaa«aa«a<

                        a'siaaaaaaaiaaaaaaaa


                      aiaaaaaaaaaaiaaaaaaasi
o
M CO Q
^ rv«
3rt'««	«	
                                                               •aaaaaa
                    a a a a a a a a a a a a— — — — — — — — —— — a a a a a a       si


U4M>-*-	•	
                    —   	   ... —                        .isa!a«aaa>saaa
    Bi»a«8)«««
                  ai« a> ai ai ai is a> a> — -
x                sissi«BaiB«>i»ia«Biaiaia--iMtn-t-»uiu)(«)(y — s a s a s « a                   i-


  (SJ .............................................       _1
c£                  st * st «» « «j  «» «* — — N n ••» « at in u» t* M — ttt t4 si « 9 « «t                   »



                    «,ai ai a ai a> a; aj ai ai Q a 3> — — (M n -»ui uj uj •»  ai            '
                                              t8<9si3i a ai a> » a> — — C
Otv< .............................................       u.
                            «i a aj ai ai ai aia a a a> a> — — —  ai a a ai                       o
—                            s aaia a a si ai aa aa a — — — — .» .» a « a a                         o
i_aaaaaaaaaaaaaasi«aaaaaa«aaia«aaaaaaa«iaaaaa3iaaaaa       —
— N .............................................       >-
u>                               si a aa si si aa a « aa s> a a a a a a aa
OaiStStS)9«IBSSlSIS)««BIS

u "                               a a a a a a a a a a a o a a a « a a a
o —
                                           a a a a a a a Q a a s
                                                   «   s
 ^^aasaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaai
 S—  .............................................
 o
 unaaas>aaaaas>aaaasiaaaasiaaaaaaaaaaS)aaaaaaaaaaaaaa

                                                        4-3