BUGSrad
 All Classes Files Functions Variables
comscp1.f
Go to the documentation of this file.
1 
2 
3 ! CVS: $Id: comscp1.F,v 1.3 2001/04/30 08:48:56 norm Exp $
4 ! CVS: $Name: $
5 
6 !-----------------------------------------------------------------------
7 
8  subroutine comscp1
9  + ( ncol , nlm , taer , tcldi
10  +, tcldw , tgm , tray , waer
11  +, wcldi , wcldw , wray , asyaer
12  +, asycldi , asycldw , tccld1 , tcclr1
13  +, asycld , asyclr , fwcld , fwclr
14  + )
15 
16  use kinds
17 
18  implicit none
19 
20 !-----------------------------------------------------------------------
21 ! REFERENCES:
22 ! comscp1 combines single scattering properties due to Rayleigh absorp
23 ! tion, aerosols, water continuum, gray gaseous absorption, ice crystals
24 ! and water droplets.
25 ! Laura D. Fowler/slikrock (08-12-97).
26 
27 ! send comments to laura@slikrock.atmos.colostate.edu and
28 ! partain@atmos.colostate.edu
29 
30 ! MODIFICATIONS:
31 ! * changed declarations to adapt the code from BUGS4 to BUGS5.
32 ! Laura D. Fowler/slikrock (02-01-00).
33 
34 ! SUBROUTINES CALLED:
35 ! none.
36 
37 ! FUNCTIONS CALLED:
38 ! none.
39 
40 ! INCLUDED COMMONS:
41 ! none.
42 
43 ! ARGUMENT LIST VARIABLES:
44 ! All arrays indexed as nlm correspond to variables defined in the
45 ! middle of layers. In this subroutine, all the arrays are defined as
46 ! local arrays in BUGSswr.
47 
48 ! INPUT ARGUMENTS:
49 ! ----------------
50  integer (kind=int_kind), intent(in)::
51  & ncol !Length of sub-domain..
52  &, nlm !Number of layers.
53 
54  real (kind=dbl_kind), intent(in), dimension(ncol,nlm)::
55  & asyaer !Asymmetry factor of aerosols (-).
56  &, asycldi !Asymmetry factor of ice clouds (-).
57  &, asycldw !Asymmetry factor of water clouds (-).
58  &, taer !Optical depth of aerosols (-).
59  &, tcldi !Optical depth of ice clouds (-).
60  &, tcldw !Optical depth of water clouds (-).
61  &, tgm !Optical depth of water vapor continuum (-).
62  &, tray !Optical depth due to Rayleigh absorption (-).
63  &, waer !Single scattering albedo of aerosols (-).
64  &, wcldi !Single scattering albedo of ice clouds (-).
65  &, wcldw !Single scattering albedo of water clouds (-).
66  &, wray !Single scattering albedo due to Rayleigh
67 ! absorption (-).
68 
69 ! OUTPUT ARGUMENTS:
70 ! -----------------
71  real (kind=dbl_kind), intent(out), dimension(ncol,nlm)::
72  & asyclr !Clear-sky asymmetry factor (-).
73  &, asycld !All-sky asymmetry factor (-).
74  &, tcclr1 !Clear-sky optical depth (-).
75  &, tccld1 !All-sky optical depth (-).
76  &, fwclr !Total clear-sky single-scattering albedo (-).
77  &, fwcld !Total cloudy single-scattering albedo (-).
78 
79 ! LOCAL LIST VARIABLES:
80 
81  integer (kind=int_kind)::
82  & i !Horizontal index.
83  &, l !Vertical index.
84 
85  real (kind=dbl_kind)::
86  & wwray,wwaer,wwcldi,wwcldw
87 
88 !-----------------------------------------------------------------------
89 
90  do l = 1, nlm
91  do i = 1, ncol
92 
93  tcclr1(i,l) = tgm(i,l) + tray(i,l) + taer(i,l)
94  tccld1(i,l) = tcclr1(i,l)+ tcldi(i,l) + tcldw(i,l)
95 
96  wwray = wray(i,l)*tray(i,l)
97  wwaer = waer(i,l)*taer(i,l)
98  wwcldi = wcldi(i,l)*tcldi(i,l)
99  wwcldw = wcldw(i,l)*tcldw(i,l)
100 
101  fwclr(i,l) = wwray+wwaer
102  fwcld(i,l) = fwclr(i,l)+wwcldi+wwcldw
103 
104  if(fwclr(i,l).gt.1.e-10) then
105  asyclr(i,l) = (asyaer(i,l)*wwaer)/fwclr(i,l)
106  else
107  asyclr(i,l) = 1.
108  endif
109 
110  if(fwcld(i,l).gt.1.e-10) then
111  asycld(i,l) = (asyaer(i,l)*wwaer+asycldi(i,l)*wwcldi
112  + + asycldw(i,l)*wwcldw)
113  + / fwcld(i,l)
114  else
115  asycld(i,l) = 1.
116  endif
117 
118  enddo
119  enddo
120 
121  return
122  end subroutine comscp1
123 
124 !-----------------------------------------------------------------------