BUGSrad
 All Classes Files Functions Variables
cloudg.f
Go to the documentation of this file.
1 
2 
3 ! CVS: $Id: cloudg.F,v 1.8 2006/11/16 19:54:12 norm Exp $
4 ! CVS: $Name: $
5 
6 !-----------------------------------------------------------------------
7 
8  subroutine cloudg
9  + ( ncol , nlm , mb , ib
10  +, pp , tt , wcont , re
11  +, pdist , cnrw , cniw , cnri
12  +, cnii , xlam , tcld , wcld
13  +, asycld , flag
14  + )
15 
16  use kinds
17 
18  implicit none
19 
20 !-----------------------------------------------------------------------
21 ! REFERENCES:
22 ! Cleaned up version of cloud.f from G. Stephens. Computes the cloud
23 ! optical properties.
24 
25 ! This routine has been modified to include the the corrections to
26 ! ADT for spherical particles based upon the work of David Mitchell
27 ! DRI. All the derivations have been carried out for the modified
28 ! gamma distribution assuming that m=0.5 (em in program), a
29 ! parameter in eqn (5) of Mitchell (1994).
30 
31 ! tcld, wcld, asycld are the optical depth, single scattering albedo,
32 ! and asymmetry parameter of cloud particles based on the use of
33 ! ADT theory as used by Stephens et al (1990). Effective radius re
34 ! is input (in microns) and the water content is in g/m3. The logical
35 ! variable flag is .false. for water and .true. for ice.
36 
37 ! send comments to laura@slikrock.atmos.colostate.edu and
38 ! partain@atmos.colostate.edu.
39 
40 ! MODIFICATIONS:
41 ! * changed declarations to adapt the code from BUGS4 to BUGS5.
42 ! Laura D. Fowler/slikrock (02-01-00).
43 
44 ! SUBROUTINES CALLED:
45 ! none.
46 
47 ! FUNCTIONS CALLED:
48 ! none
49 
50 ! INCLUDED COMMONS:
51 ! none.
52 
53 ! ARGUMENT LIST VARIABLES:
54 ! INPUT ARGUMENTS:
55 ! ----------------
56  logical (kind=log_kind), intent(in)::
57  & flag !If true, computes optical properties of ice clouds, of
58 ! of water clouds otherwise.
59 
60  integer (kind=int_kind), intent(in)::
61  & ncol !Length of sub-domain.
62  &, nlm !Number of layers.
63  &, mb !Total number of spectral intervals.
64  &, ib !Index of spectral interval.
65 
66  real (kind=dbl_kind), intent(in), dimension(mb)::
67  & cnrw !Real part of refractive index (Water clouds).
68  &, cniw !Imaginary part of refractive index (Water clouds).
69  &, cnri !Real part of refractive index (Ice clouds).
70  &, cnii !Imaginary part of refractive index (Ice clouds).
71  &, xlam !Center of spectral band.
72 
73  real (kind=dbl_kind), intent(in)::
74  & pdist !
75 
76  real (kind=dbl_kind), intent(in), dimension(ncol,nlm)::
77  & tt !Temperature (K).
78  &, wcont !Cloud water/ice content (g/m^-3).
79  &, re !Cloud effective radius (mu).
80 
81  real (kind=dbl_kind), intent(in), dimension(ncol,nlm+1)::
82  & pp !Level pressure (hPa).
83 
84 ! OUTPUT ARGUMENTS:
85 ! -----------------
86  real (kind=dbl_kind), intent(out), dimension(ncol,nlm)::
87  & tcld !Cloud optical depth (-).
88  &, wcld !Cloud single scattering albedo (-).
89  &, asycld !Cloud asymmetry factor (-).
90 
91 ! LOCAL VARIABLES:
92  complex (kind=dbl_kind)::
93  & cm,um
94 
95  integer (kind=int_kind)::
96  & i, l
97 
98  real (kind=dbl_kind)::
99  & abs , area , c0 , c1
100  &, cnr , cni , dz , eps
101  &, ext , f2 , f3 , no
102  &, p0 , p1 , p2 , pi
103  &, rm , xm , vm , rho_water
104 
105 !-----------------------------------------------------------------------
106 
107 !---- initialize local and output arrays:
108 
109  tcld(:,:) = 0.
110  wcld(:,:) = 0.
111  asycld(:,:) = 1.
112 
113 !-- initialize miscellaneous constants and indices of refraction:
114 
115  pi = acos(-1.)
116  eps = 1.e-5
117 
118  if(flag) then
119  !Ice case
120  cnr = cnri(ib)
121  cni = cnii(ib)
122  rho_water = 1.0e6 !gm^-3
123  else
124  cnr = cnrw(ib)
125  cni = cniw(ib)
126  rho_water = 1.0e6 !gm^-3
127  endif
128 
129 !-- constants depending upon the characteristic width of the distribu
130 ! tion.(these may be made to vary with hydrometeor species and thus
131 ! pdist could be made to depend upon level and column numbers).
132 
133 ! p0 = 0.
134  p0 = pdist
135  p1 = p0 + 1.
136  p2 = p0 + 2.
137  f2 = p1 * p0
138  f3 = p2 * f2
139 
140 !---- calculate cloud optical properties:
141 
142  do l = 1, nlm
143  do i = 1, ncol
144  if(wcont(i,l) .gt. eps) then
145  dz=29.286*log(pp(i,l+1)/pp(i,l)) * tt(i,l)
146  rm = re(i,l)/p2
147  no = wcont(i,l) / ( (4.*pi/3.)*f3*rho_water*rm**3 ) !Particles per cubic micrometer
148  area = pi*f2*no*rm**2*1.0e6 !The factor converts inverse micrometers,
149  !i.e., micrometer^2/micrometer^3, to inverse meters
150  c0 = 2.*area
151  c1 = c0/f2
152  xm = 2.*pi*rm/xlam(ib)
153  cm = cmplx(cnr,-cni)
154 
155  if (ib .eq. 1) then
156  !For band 1 (0.5 um) only compute the extinction.
157  um = 2.*xm*(cnr-1.)*cmplx(0.d0,1.d0)
158  ext = c0 + 2.*c1*real(p0/(um*(um+1.)**p1)
159  + + 1./(um**2*(um+1.)**p0)-1./um**2)
160  tcld(i,l) = ext*dz
161  wcld(i,l) = 0.999999
162  asycld(i,l) = 0.85
163  else
164  !Compute both extinction and absorption coefficients for all other bands.
165  um = 2.*xm*(cm-1.)*cmplx(0.d0,1.d0)
166  ext = c0 + 2.*c1*real( p0/(um*(um+1.)**p1)
167  + + 1./(um**2*(um+1.)**p0)-1./um**2)
168  vm = 4.*xm*cni
169  abs = area + c1*sngl( p0/(vm*(vm+1.)**p1)
170  + + 1./(vm**2*(vm+1.)**p0) - 1./vm**2 )
171 ! abs = area + c1*sngl(
172 ! + p0/(dble(vm)*(dble(vm)+1.)**dble(p1))
173 ! + + 1./(dble(vm)**2*(dble(vm)+1.)**dble(p0))
174 ! + - 1./dble(vm)**2)
175  tcld(i,l) = ext*dz
176 
177  if (ext.lt.abs) ext = abs
178  wcld(i,l) = (ext-abs)/ext
179 
180  if(wcld(i,l) .lt. 0.) then
181  print *,wcld(i,l),ext,abs,wcont(i,l)
182  print *,pp(i,l),pp(i,l+1)
183  print *,tt(i,l)
184  print *,re(i,l)
185  stop
186  endif
187  endif
188  asycld(i,l)=0.85 !default, overridden in bugs_swr(), bugs_lwr()
189  endif
190  enddo
191  enddo
192 
193  return
194  end subroutine cloudg
195 !-----------------------------------------------------------------------