Line data Source code
1 : module aerosol_optics_cam
2 : use shr_kind_mod, only: r8 => shr_kind_r8
3 : use shr_kind_mod, only: cl => shr_kind_cl
4 : use cam_logfile, only: iulog
5 : use radconstants, only: nswbands, nlwbands, idx_sw_diag, idx_uv_diag, idx_nir_diag
6 : use radconstants, only: get_lw_spectral_boundaries
7 : use phys_prop, only: ot_length
8 : use physics_types,only: physics_state
9 : use physics_buffer,only: physics_buffer_desc
10 : use ppgrid, only: pcols, pver
11 : use physconst, only: rga, rair
12 : use cam_abortutils, only: endrun
13 : use spmd_utils, only: masterproc
14 : use rad_constituents, only: n_diag, rad_cnst_get_call_list
15 : use cam_history, only: addfld, add_default, outfld, horiz_only, fieldname_len
16 : use cam_history_support, only: fillvalue
17 :
18 : use tropopause, only : tropopause_findChemTrop
19 : use wv_saturation, only: qsat
20 :
21 : use aerosol_properties_mod, only: aerosol_properties
22 : use modal_aerosol_properties_mod, only: modal_aerosol_properties
23 : use carma_aerosol_properties_mod, only: carma_aerosol_properties
24 :
25 : use aerosol_state_mod, only: aerosol_state
26 : use modal_aerosol_state_mod,only: modal_aerosol_state
27 : use carma_aerosol_state_mod,only: carma_aerosol_state
28 :
29 : use aerosol_optics_mod, only: aerosol_optics
30 : use refractive_aerosol_optics_mod, only: refractive_aerosol_optics
31 : use hygrocoreshell_aerosol_optics_mod, only: hygrocoreshell_aerosol_optics
32 : use hygrowghtpct_aerosol_optics_mod, only: hygrowghtpct_aerosol_optics
33 :
34 : implicit none
35 :
36 : private
37 :
38 : public :: aerosol_optics_cam_readnl
39 : public :: aerosol_optics_cam_init
40 : public :: aerosol_optics_cam_final
41 : public :: aerosol_optics_cam_sw
42 : public :: aerosol_optics_cam_lw
43 :
44 : type aero_props_t
45 : class(aerosol_properties), pointer :: obj => null()
46 : end type aero_props_t
47 : type aero_state_t
48 : class(aerosol_state), pointer :: obj => null()
49 : end type aero_state_t
50 :
51 : type(aero_props_t), allocatable :: aero_props(:) ! array of aerosol properties objects to allow for
52 : ! multiple aerosol representations in the same sim
53 : ! such as MAM and CARMA
54 :
55 : ! refractive index for water read in read_water_refindex
56 : complex(r8) :: crefwsw(nswbands) = -huge(1._r8) ! complex refractive index for water visible
57 : complex(r8) :: crefwlw(nlwbands) = -huge(1._r8) ! complex refractive index for water infrared
58 : character(len=cl) :: water_refindex_file = 'NONE' ! full pathname for water refractive index dataset
59 :
60 : logical :: carma_active = .false.
61 : logical :: modal_active = .false.
62 : integer :: num_aero_models = 0
63 : integer :: lw10um_indx = -1 ! wavelength index corresponding to 10 microns
64 : real(r8), parameter :: lw10um = 10._r8 ! microns
65 :
66 : character(len=4) :: diag(0:n_diag) = (/' ','_d1 ','_d2 ','_d3 ','_d4 ','_d5 ', '_d6 ','_d7 ','_d8 ','_d9 ','_d10'/)
67 :
68 : type out_name
69 : character(len=fieldname_len), allocatable :: name(:) ! nbins
70 : end type out_name
71 :
72 : type(out_name), allocatable :: burden_fields(:) ! num_aero_models
73 : type(out_name), allocatable :: aodbin_fields(:)
74 : type(out_name), allocatable :: aoddust_fields(:)
75 : type(out_name), allocatable :: burdendn_fields(:) ! num_aero_models
76 : type(out_name), allocatable :: aodbindn_fields(:)
77 : type(out_name), allocatable :: aoddustdn_fields(:)
78 :
79 : contains
80 :
81 : !===============================================================================
82 1024 : subroutine aerosol_optics_cam_readnl(nlfile)
83 : use namelist_utils, only : find_group_name
84 : use spmd_utils, only : mpicom, masterprocid, mpi_character, mpi_success
85 :
86 : character(len=*), intent(in) :: nlfile ! filepath for file containing namelist input
87 :
88 : integer :: unitn, ierr
89 : character(len=cl) :: errmsg
90 : character(len=*), parameter :: subname = 'aerosol_optics_cam_readnl'
91 :
92 : ! ===================
93 : ! Namelist definition
94 : ! ===================
95 : namelist /aerosol_optics_nl/ water_refindex_file
96 :
97 : ! =============
98 : ! Read namelist
99 : ! =============
100 1026 : if (masterproc) then
101 2 : open( newunit=unitn, file=trim(nlfile), status='old' )
102 2 : call find_group_name(unitn, 'aerosol_optics_nl', status=ierr)
103 2 : if (ierr == 0) then
104 0 : read(unitn, aerosol_optics_nl, iostat=ierr)
105 0 : if (ierr /= 0) then
106 0 : write(errmsg,'(2a,i10)') subname,':: ERROR reading namelist, error code: ',ierr
107 0 : call endrun(errmsg)
108 : end if
109 : end if
110 2 : close(unitn)
111 : end if
112 :
113 : ! ============================
114 : ! Broadcast namelist variables
115 : ! ============================
116 1024 : call mpi_bcast(water_refindex_file, len(water_refindex_file), mpi_character, masterprocid, mpicom, ierr)
117 1024 : if (ierr/=mpi_success) then
118 0 : call endrun(subname // ':: ERROR mpi_bcast '//trim(water_refindex_file))
119 : end if
120 :
121 1024 : if (masterproc) then
122 2 : write(iulog,*) subname,': water_refindex_file = ',trim(water_refindex_file)
123 : end if
124 :
125 1024 : end subroutine aerosol_optics_cam_readnl
126 :
127 : !===============================================================================
128 0 : subroutine aerosol_optics_cam_init
129 : use rad_constituents, only: rad_cnst_get_info
130 : use phys_control, only: phys_getopts
131 : use ioFileMod, only: getfil
132 :
133 : character(len=*), parameter :: prefix = 'aerosol_optics_cam_init: '
134 : integer :: nmodes=0, nbins=0, iaermod, istat, ilist, i
135 :
136 : logical :: call_list(0:n_diag)
137 : real(r8) :: lwavlen_lo(nlwbands), lwavlen_hi(nlwbands)
138 : integer :: m, n
139 :
140 : character(len=fieldname_len) :: fldname
141 : character(len=128) :: lngname
142 : logical :: history_aero_optics ! output aerosol optics diagnostics
143 : logical :: history_amwg ! output the variables used by the AMWG diag package
144 : logical :: history_dust ! output dust diagnostics
145 :
146 : character(len=cl) :: locfile
147 :
148 : call phys_getopts(history_amwg_out = history_amwg, &
149 : history_aero_optics_out = history_aero_optics, &
150 0 : history_dust_out = history_dust )
151 :
152 0 : num_aero_models = 0
153 :
154 0 : call rad_cnst_get_info(0, nmodes=nmodes, nbins=nbins)
155 0 : modal_active = nmodes>0
156 0 : carma_active = nbins>0
157 :
158 : ! count aerosol models
159 0 : if (modal_active) then
160 0 : num_aero_models = num_aero_models+1
161 : end if
162 0 : if (carma_active) then
163 0 : num_aero_models = num_aero_models+1
164 : end if
165 :
166 0 : if (num_aero_models>0) then
167 0 : allocate(aero_props(num_aero_models), stat=istat)
168 0 : if (istat/=0) then
169 0 : call endrun(prefix//'array allocation error: aero_props')
170 : end if
171 : end if
172 :
173 0 : iaermod = 0
174 :
175 0 : if (modal_active) then
176 0 : iaermod = iaermod+1
177 0 : aero_props(iaermod)%obj => modal_aerosol_properties()
178 : end if
179 0 : if (carma_active) then
180 0 : iaermod = iaermod+1
181 0 : aero_props(iaermod)%obj => carma_aerosol_properties()
182 : end if
183 :
184 0 : if (water_refindex_file=='NONE') then
185 0 : call endrun(prefix//'water_refindex_file must be specified')
186 : else
187 0 : call getfil(water_refindex_file, locfile)
188 0 : call read_water_refindex(locfile)
189 : end if
190 :
191 0 : call get_lw_spectral_boundaries(lwavlen_lo, lwavlen_hi, units='um')
192 0 : do i = 1,nlwbands
193 0 : if ((lwavlen_lo(i)<=lw10um) .and. (lwavlen_hi(i)>=lw10um)) then
194 0 : lw10um_indx = i ! index corresponding to 10 microns
195 : end if
196 : end do
197 0 : call rad_cnst_get_call_list(call_list)
198 :
199 0 : do ilist = 0, n_diag
200 0 : if (call_list(ilist)) then
201 0 : call addfld ('EXTINCT'//diag(ilist), (/ 'lev' /), 'A','/m',&
202 0 : 'Aerosol extinction 550 nm, day only', flag_xyfill=.true.)
203 0 : call addfld ('EXTINCTUV'//diag(ilist), (/ 'lev' /), 'A','/m',&
204 0 : 'Aerosol extinction 350 nm, day only', flag_xyfill=.true.)
205 0 : call addfld ('EXTINCTNIR'//diag(ilist), (/ 'lev' /), 'A','/m',&
206 0 : 'Aerosol extinction 1020 nm, day only', flag_xyfill=.true.)
207 0 : call addfld ('ABSORB'//diag(ilist), (/ 'lev' /), 'A','/m',&
208 0 : 'Aerosol absorption, day only', flag_xyfill=.true.)
209 0 : call addfld ('AODVIS'//diag(ilist), horiz_only, 'A',' ', &
210 0 : 'Aerosol optical depth 550 nm', flag_xyfill=.true.)
211 0 : call addfld ('AODVISst'//diag(ilist), horiz_only, 'A',' ', &
212 0 : 'Stratospheric aerosol optical depth 550 nm, day only', flag_xyfill=.true.)
213 0 : call addfld ('AODNIRst'//diag(ilist), horiz_only, 'A',' ', &
214 0 : 'Stratospheric aerosol optical depth 1020 nm, day only',flag_xyfill=.true.)
215 0 : call addfld ('AODUVst'//diag(ilist), horiz_only, 'A',' ', &
216 0 : 'Stratospheric aerosol optical depth 350 nm, day only', flag_xyfill=.true.)
217 0 : call addfld ('AODUV'//diag(ilist), horiz_only, 'A',' ', &
218 0 : 'Aerosol optical depth 350 nm, day only', flag_xyfill=.true.)
219 0 : call addfld ('AODNIR'//diag(ilist), horiz_only, 'A',' ', &
220 0 : 'Aerosol optical depth 1020 nm, day only',flag_xyfill=.true.)
221 0 : call addfld ('AODABS'//diag(ilist), horiz_only, 'A',' ', &
222 0 : 'Aerosol absorption optical depth 550 nm, day only', flag_xyfill=.true.)
223 0 : call addfld ('AODxASYM'//diag(ilist), horiz_only, 'A',' ', &
224 0 : 'Aerosol optical depth 550 * asymmetry factor, day only', flag_xyfill=.true.)
225 0 : call addfld ('EXTxASYM'//diag(ilist), (/ 'lev' /), 'A',' ', &
226 0 : 'extinction 550 nm * asymmetry factor, day only', flag_xyfill=.true.)
227 0 : call addfld ('AODTOT'//diag(ilist), horiz_only, 'A','1',&
228 0 : 'Aerosol optical depth summed over all sw wavelengths', flag_xyfill=.true.)
229 :
230 0 : call addfld ('EXTINCTdn'//diag(ilist), (/ 'lev' /), 'A','/m',&
231 0 : 'Aerosol extinction 550 nm, day only', flag_xyfill=.true.)
232 0 : call addfld ('EXTINCTUVdn'//diag(ilist), (/ 'lev' /), 'A','/m',&
233 0 : 'Aerosol extinction 350 nm, day only', flag_xyfill=.true.)
234 0 : call addfld ('EXTINCTNIRdn'//diag(ilist), (/ 'lev' /), 'A','/m',&
235 0 : 'Aerosol extinction 1020 nm, day only', flag_xyfill=.true.)
236 0 : call addfld ('ABSORBdn'//diag(ilist), (/ 'lev' /), 'A','/m',&
237 0 : 'Aerosol absorption, day only', flag_xyfill=.true.)
238 0 : call addfld ('AODVISdn'//diag(ilist), horiz_only, 'A',' ', &
239 0 : 'Aerosol optical depth 550 nm', flag_xyfill=.true.)
240 0 : call addfld ('AODVISstdn'//diag(ilist), horiz_only, 'A',' ', &
241 0 : 'Stratospheric aerosol optical depth 550 nm, day only', flag_xyfill=.true.)
242 0 : call addfld ('AODNIRstdn'//diag(ilist), horiz_only, 'A',' ', &
243 0 : 'Stratospheric aerosol optical depth 1020 nm, day only', flag_xyfill=.true.)
244 0 : call addfld ('AODUVstdn'//diag(ilist), horiz_only, 'A',' ', &
245 0 : 'Stratospheric aerosol optical depth 350 nm, day only', flag_xyfill=.true.)
246 0 : call addfld ('AODUVdn'//diag(ilist), horiz_only, 'A',' ', &
247 0 : 'Aerosol optical depth 350 nm, day only', flag_xyfill=.true.)
248 0 : call addfld ('AODNIRdn'//diag(ilist), horiz_only, 'A',' ', &
249 0 : 'Aerosol optical depth 1020 nm, day only', flag_xyfill=.true.)
250 0 : call addfld ('AODABSdn'//diag(ilist), horiz_only, 'A',' ', &
251 0 : 'Aerosol absorption optical depth 550 nm, day only', flag_xyfill=.true.)
252 0 : call addfld ('AODxASYMdn'//diag(ilist), horiz_only, 'A',' ', &
253 0 : 'Aerosol optical depth 550 * asymmetry factor, day only', flag_xyfill=.true.)
254 0 : call addfld ('EXTxASYMdn'//diag(ilist), (/ 'lev' /), 'A',' ', &
255 0 : 'extinction 550 nm * asymmetry factor, day only', flag_xyfill=.true.)
256 0 : call addfld ('AODTOTdn'//diag(ilist), horiz_only, 'A','1',&
257 0 : 'Aerosol optical depth summed over all sw wavelengths, day only')
258 :
259 0 : if (lw10um_indx>0) then
260 0 : call addfld('AODABSLW'//diag(ilist), (/ 'lev' /), 'A','/m',&
261 0 : 'Aerosol long-wave absorption optical depth at 10 microns')
262 : end if
263 0 : call addfld ('TOTABSLW'//diag(ilist), (/ 'lev' /), 'A',' ', &
264 0 : 'LW Aero total abs')
265 :
266 0 : if (ilist>0 .and. history_aero_optics) then
267 0 : call add_default ('EXTINCT'//diag(ilist), 1, ' ')
268 0 : call add_default ('ABSORB'//diag(ilist), 1, ' ')
269 0 : call add_default ('AODVIS'//diag(ilist), 1, ' ')
270 0 : call add_default ('AODVISst'//diag(ilist), 1, ' ')
271 0 : call add_default ('AODABS'//diag(ilist), 1, ' ')
272 : end if
273 :
274 : end if
275 : end do
276 :
277 0 : if (num_aero_models>0) then
278 :
279 0 : allocate(burden_fields(num_aero_models), stat=istat)
280 0 : if (istat/=0) then
281 0 : call endrun(prefix//'array allocation error: burden_fields')
282 : end if
283 0 : allocate(aodbin_fields(num_aero_models), stat=istat)
284 0 : if (istat/=0) then
285 0 : call endrun(prefix//'array allocation error: aodbin_fields')
286 : end if
287 0 : allocate(aoddust_fields(num_aero_models), stat=istat)
288 0 : if (istat/=0) then
289 0 : call endrun(prefix//'array allocation error: aoddust_fields')
290 : end if
291 :
292 0 : allocate(burdendn_fields(num_aero_models), stat=istat)
293 0 : if (istat/=0) then
294 0 : call endrun(prefix//'array allocation error: burdendn_fields')
295 : end if
296 0 : allocate(aodbindn_fields(num_aero_models), stat=istat)
297 0 : if (istat/=0) then
298 0 : call endrun(prefix//'array allocation error: aodbindn_fields')
299 : end if
300 0 : allocate(aoddustdn_fields(num_aero_models), stat=istat)
301 0 : if (istat/=0) then
302 0 : call endrun(prefix//'array allocation error: aoddustdn_fields')
303 : end if
304 :
305 0 : do n = 1,num_aero_models
306 :
307 0 : allocate(burden_fields(n)%name(aero_props(n)%obj%nbins()), stat=istat)
308 0 : if (istat/=0) then
309 0 : call endrun(prefix//'array allocation error: burden_fields(n)%name')
310 : end if
311 0 : allocate(aodbin_fields(n)%name(aero_props(n)%obj%nbins()), stat=istat)
312 0 : if (istat/=0) then
313 0 : call endrun(prefix//'array allocation error: aodbin_fields(n)%name')
314 : end if
315 0 : allocate(aoddust_fields(n)%name(aero_props(n)%obj%nbins()), stat=istat)
316 0 : if (istat/=0) then
317 0 : call endrun(prefix//'array allocation error: aoddust_fields(n)%name')
318 : end if
319 :
320 0 : allocate(burdendn_fields(n)%name(aero_props(n)%obj%nbins()), stat=istat)
321 0 : if (istat/=0) then
322 0 : call endrun(prefix//'array allocation error: burdendn_fields(n)%name')
323 : end if
324 0 : allocate(aodbindn_fields(n)%name(aero_props(n)%obj%nbins()), stat=istat)
325 0 : if (istat/=0) then
326 0 : call endrun(prefix//'array allocation error: aodbindn_fields(n)%name')
327 : end if
328 0 : allocate(aoddustdn_fields(n)%name(aero_props(n)%obj%nbins()), stat=istat)
329 0 : if (istat/=0) then
330 0 : call endrun(prefix//'array allocation error: aoddustdn_fields(n)%name')
331 : end if
332 :
333 0 : do m = 1, aero_props(n)%obj%nbins()
334 :
335 0 : write(fldname,'(a,i2.2)') 'BURDEN', m
336 0 : burden_fields(n)%name(m) = fldname
337 0 : write(lngname,'(a,i2.2)') 'Aerosol burden bin ', m
338 0 : call addfld (fldname, horiz_only, 'A', 'kg/m2', lngname, flag_xyfill=.true.)
339 0 : if (history_aero_optics) then
340 0 : call add_default (fldname, 1, ' ')
341 : end if
342 :
343 0 : fldname = 'AOD_'//trim(aero_props(n)%obj%bin_name(0,m))
344 0 : aodbin_fields(n)%name(m) = fldname
345 0 : lngname = 'Aerosol optical depth, day only, 550 nm, '//trim(aero_props(n)%obj%bin_name(0,m))
346 0 : call addfld (aodbin_fields(n)%name(m), horiz_only, 'A', ' ', lngname, flag_xyfill=.true.)
347 0 : if (history_aero_optics) then
348 0 : call add_default (fldname, 1, ' ')
349 : end if
350 :
351 0 : write(fldname,'(a,i2.2)') 'AODDUST', m
352 0 : aoddust_fields(n)%name(m) = fldname
353 0 : write(lngname,'(a,i2,a)') 'Aerosol optical depth, day only, 550 nm mode ',m,' from dust'
354 0 : call addfld (aoddust_fields(n)%name(m), horiz_only, 'A', ' ', lngname, flag_xyfill=.true.)
355 0 : if (history_aero_optics) then
356 0 : call add_default (fldname, 1, ' ')
357 : end if
358 :
359 0 : write(fldname,'(a,i2.2)') 'BURDENdn', m
360 0 : burdendn_fields(n)%name(m) = fldname
361 0 : write(lngname,'(a,i2)') 'Aerosol burden, day night, bin ', m
362 0 : call addfld (burdendn_fields(n)%name(m), horiz_only, 'A', 'kg/m2', lngname, flag_xyfill=.true.)
363 0 : if (history_aero_optics) then
364 0 : call add_default (fldname, 1, ' ')
365 : end if
366 :
367 0 : fldname = 'AODdn_'//trim(aero_props(n)%obj%bin_name(0,m))
368 0 : aodbindn_fields(n)%name(m) = fldname
369 0 : lngname = 'Aerosol optical depth 550 nm, day night, '//trim(aero_props(n)%obj%bin_name(0,m))
370 0 : call addfld (aodbindn_fields(n)%name(m), horiz_only, 'A', ' ', lngname, flag_xyfill=.true.)
371 0 : if (history_aero_optics) then
372 0 : call add_default (fldname, 1, ' ')
373 : end if
374 :
375 0 : write(fldname,'(a,i2.2)') 'AODdnDUST', m
376 0 : aoddustdn_fields(n)%name(m) = fldname
377 0 : write(lngname,'(a,i2,a)') 'Aerosol optical depth 550 nm, day night, bin ',m,' from dust'
378 0 : call addfld (aoddustdn_fields(n)%name(m), horiz_only, 'A', ' ', lngname, flag_xyfill=.true.)
379 0 : if (history_aero_optics) then
380 0 : call add_default (fldname, 1, ' ')
381 : end if
382 :
383 : end do
384 :
385 : end do
386 :
387 : end if
388 :
389 : call addfld ('AODDUST', horiz_only, 'A',' ', 'Aerosol optical depth 550 nm from dust, day only', &
390 0 : flag_xyfill=.true.)
391 : call addfld ('AODSO4', horiz_only, 'A',' ', 'Aerosol optical depth 550 nm from SO4, day only', &
392 0 : flag_xyfill=.true.)
393 : call addfld ('AODPOM', horiz_only, 'A',' ', 'Aerosol optical depth 550 nm from POM, day only', &
394 0 : flag_xyfill=.true.)
395 : call addfld ('AODSOA', horiz_only, 'A',' ', 'Aerosol optical depth 550 nm from SOA, day only', &
396 0 : flag_xyfill=.true.)
397 : call addfld ('AODBC', horiz_only, 'A',' ', 'Aerosol optical depth 550 nm from BC, day only', &
398 0 : flag_xyfill=.true.)
399 : call addfld ('AODSS', horiz_only, 'A',' ', 'Aerosol optical depth 550 nm from seasalt, day only', &
400 0 : flag_xyfill=.true.)
401 : call addfld ('AODABSBC', horiz_only, 'A',' ', 'Aerosol absorption optical depth 550 nm from BC, day only',&
402 0 : flag_xyfill=.true.)
403 : call addfld ('BURDENDUST', horiz_only, 'A','kg/m2', 'Dust aerosol burden, day only' , &
404 0 : flag_xyfill=.true.)
405 : call addfld ('BURDENSO4', horiz_only, 'A','kg/m2', 'Sulfate aerosol burden, day only' , &
406 0 : flag_xyfill=.true.)
407 : call addfld ('BURDENPOM', horiz_only, 'A','kg/m2', 'POM aerosol burden, day only' , &
408 0 : flag_xyfill=.true.)
409 : call addfld ('BURDENSOA', horiz_only, 'A','kg/m2', 'SOA aerosol burden, day only' , &
410 0 : flag_xyfill=.true.)
411 : call addfld ('BURDENBC', horiz_only, 'A','kg/m2', 'Black carbon aerosol burden, day only', &
412 0 : flag_xyfill=.true.)
413 : call addfld ('BURDENSEASALT', horiz_only, 'A','kg/m2', 'Seasalt aerosol burden, day only' , &
414 0 : flag_xyfill=.true.)
415 : call addfld ('SSAVIS', horiz_only, 'A',' ', 'Aerosol single-scatter albedo, day only', &
416 0 : flag_xyfill=.true.)
417 :
418 : call addfld ('AODDUSTdn', horiz_only, 'A',' ', 'Aerosol optical depth 550 nm from dust, day night', &
419 0 : flag_xyfill=.true.)
420 : call addfld ('AODSO4dn', horiz_only, 'A',' ', 'Aerosol optical depth 550 nm from SO4, day night', &
421 0 : flag_xyfill=.true.)
422 : call addfld ('AODPOMdn', horiz_only, 'A',' ', 'Aerosol optical depth 550 nm from POM, day night', &
423 0 : flag_xyfill=.true.)
424 : call addfld ('AODSOAdn', horiz_only, 'A',' ', 'Aerosol optical depth 550 nm from SOA, day night', &
425 0 : flag_xyfill=.true.)
426 : call addfld ('AODBCdn', horiz_only, 'A',' ', 'Aerosol optical depth 550 nm from BC, day night', &
427 0 : flag_xyfill=.true.)
428 : call addfld ('AODSSdn', horiz_only, 'A',' ', 'Aerosol optical depth 550 nm from seasalt, day night', &
429 0 : flag_xyfill=.true.)
430 : call addfld ('AODABSBCdn', horiz_only, 'A',' ', 'Aerosol absorption optical depth 550 nm from BC, day night',&
431 0 : flag_xyfill=.true.)
432 : call addfld ('BURDENDUSTdn', horiz_only, 'A','kg/m2', 'Dust aerosol burden, day night' , &
433 0 : flag_xyfill=.true.)
434 : call addfld ('BURDENSO4dn', horiz_only, 'A','kg/m2', 'Sulfate aerosol burden, day night' , &
435 0 : flag_xyfill=.true.)
436 : call addfld ('BURDENPOMdn', horiz_only, 'A','kg/m2', 'POM aerosol burden, day night' , &
437 0 : flag_xyfill=.true.)
438 : call addfld ('BURDENSOAdn', horiz_only, 'A','kg/m2', 'SOA aerosol burden, day night' , &
439 0 : flag_xyfill=.true.)
440 : call addfld ('BURDENBCdn', horiz_only, 'A','kg/m2', 'Black carbon aerosol burden, day night', &
441 0 : flag_xyfill=.true.)
442 : call addfld ('BURDENSEASALTdn', horiz_only, 'A','kg/m2', 'Seasalt aerosol burden, day night' , &
443 0 : flag_xyfill=.true.)
444 : call addfld ('SSAVISdn', horiz_only, 'A',' ', 'Aerosol single-scatter albedo, day night', &
445 0 : flag_xyfill=.true.)
446 :
447 0 : if (history_amwg) then
448 0 : call add_default ('AODDUST01' , 1, ' ')
449 0 : call add_default ('AODDUST03' , 1, ' ')
450 0 : call add_default ('AODDUST' , 1, ' ')
451 0 : call add_default ('AODVIS' , 1, ' ')
452 : end if
453 :
454 0 : if (history_dust) then
455 0 : call add_default ('AODDUST01' , 1, ' ')
456 0 : call add_default ('AODDUST02' , 1, ' ')
457 0 : call add_default ('AODDUST03' , 1, ' ')
458 : end if
459 :
460 0 : if (history_aero_optics) then
461 0 : call add_default ('AODDUST01' , 1, ' ')
462 0 : call add_default ('AODDUST03' , 1, ' ')
463 0 : call add_default ('ABSORB' , 1, ' ')
464 0 : call add_default ('AODVIS' , 1, ' ')
465 0 : call add_default ('AODUV' , 1, ' ')
466 0 : call add_default ('AODNIR' , 1, ' ')
467 0 : call add_default ('AODABS' , 1, ' ')
468 0 : call add_default ('AODABSBC' , 1, ' ')
469 0 : call add_default ('AODDUST' , 1, ' ')
470 0 : call add_default ('AODSO4' , 1, ' ')
471 0 : call add_default ('AODPOM' , 1, ' ')
472 0 : call add_default ('AODSOA' , 1, ' ')
473 0 : call add_default ('AODBC' , 1, ' ')
474 0 : call add_default ('AODSS' , 1, ' ')
475 0 : call add_default ('BURDEN01' , 1, ' ')
476 0 : call add_default ('BURDEN02' , 1, ' ')
477 0 : call add_default ('BURDEN03' , 1, ' ')
478 0 : call add_default ('BURDENDUST' , 1, ' ')
479 0 : call add_default ('BURDENSO4' , 1, ' ')
480 0 : call add_default ('BURDENPOM' , 1, ' ')
481 0 : call add_default ('BURDENSOA' , 1, ' ')
482 0 : call add_default ('BURDENBC' , 1, ' ')
483 0 : call add_default ('BURDENSEASALT', 1, ' ')
484 0 : call add_default ('SSAVIS' , 1, ' ')
485 0 : call add_default ('EXTINCT' , 1, ' ')
486 0 : call add_default ('AODxASYM' , 1, ' ')
487 0 : call add_default ('EXTxASYM' , 1, ' ')
488 :
489 0 : call add_default ('AODdnDUST01' , 1, ' ')
490 0 : call add_default ('AODdnDUST03' , 1, ' ')
491 0 : call add_default ('ABSORBdn' , 1, ' ')
492 0 : call add_default ('AODVISdn' , 1, ' ')
493 0 : call add_default ('AODUVdn' , 1, ' ')
494 0 : call add_default ('AODNIRdn' , 1, ' ')
495 0 : call add_default ('AODABSdn' , 1, ' ')
496 0 : call add_default ('AODABSBCdn' , 1, ' ')
497 0 : call add_default ('AODDUSTdn' , 1, ' ')
498 0 : call add_default ('AODSO4dn' , 1, ' ')
499 0 : call add_default ('AODPOMdn' , 1, ' ')
500 0 : call add_default ('AODSOAdn' , 1, ' ')
501 0 : call add_default ('AODBCdn' , 1, ' ')
502 0 : call add_default ('AODSSdn' , 1, ' ')
503 0 : call add_default ('BURDENdn01' , 1, ' ')
504 0 : call add_default ('BURDENdn02' , 1, ' ')
505 0 : call add_default ('BURDENdn03' , 1, ' ')
506 0 : call add_default ('BURDENDUSTdn' , 1, ' ')
507 0 : call add_default ('BURDENSO4dn' , 1, ' ')
508 0 : call add_default ('BURDENPOMdn' , 1, ' ')
509 0 : call add_default ('BURDENSOAdn' , 1, ' ')
510 0 : call add_default ('BURDENBCdn' , 1, ' ')
511 0 : call add_default ('BURDENSEASALTdn', 1, ' ')
512 0 : call add_default ('SSAVISdn' , 1, ' ')
513 0 : call add_default ('EXTINCTdn' , 1, ' ')
514 0 : call add_default ('AODxASYMdn' , 1, ' ')
515 0 : call add_default ('EXTxASYMdn' , 1, ' ')
516 : end if
517 :
518 0 : call addfld( 'SULFWTPCT', (/ 'lev' /), 'I', '1', 'Sulfate Weight Percent' )
519 :
520 0 : end subroutine aerosol_optics_cam_init
521 :
522 : !===============================================================================
523 0 : subroutine aerosol_optics_cam_final
524 :
525 : integer :: iaermod
526 :
527 0 : do iaermod = 1,num_aero_models
528 0 : if (associated(aero_props(iaermod)%obj)) then
529 0 : deallocate(aero_props(iaermod)%obj)
530 0 : nullify(aero_props(iaermod)%obj)
531 : end if
532 : end do
533 :
534 0 : if (allocated(aero_props)) then
535 0 : deallocate(aero_props)
536 : endif
537 :
538 0 : end subroutine aerosol_optics_cam_final
539 :
540 : !===============================================================================
541 0 : subroutine aerosol_optics_cam_sw(list_idx, state, pbuf, nnite, idxnite, tauxar, wa, ga, fa)
542 :
543 : ! calculates aerosol sw radiative properties
544 :
545 : integer, intent(in) :: list_idx ! index of the climate or a diagnostic list
546 : type(physics_state), intent(in), target :: state ! state variables
547 :
548 : type(physics_buffer_desc), pointer :: pbuf(:)
549 : integer, intent(in) :: nnite ! number of night columns
550 : integer, intent(in) :: idxnite(nnite) ! local column indices of night columns
551 :
552 : real(r8), intent(inout) :: tauxar(pcols,0:pver,nswbands) ! layer extinction optical depth
553 : real(r8), intent(inout) :: wa(pcols,0:pver,nswbands) ! layer single-scatter albedo
554 : real(r8), intent(inout) :: ga(pcols,0:pver,nswbands) ! asymmetry factor
555 : real(r8), intent(inout) :: fa(pcols,0:pver,nswbands) ! forward scattered fraction
556 :
557 : character(len=*), parameter :: prefix = 'aerosol_optics_cam_sw: '
558 :
559 : integer :: ibin, nbins
560 : integer :: iwav, ilev
561 : integer :: icol, istat
562 : integer :: lchnk, ncol
563 :
564 0 : type(aero_state_t), allocatable :: aero_state(:) ! array of aerosol state objects to allow for
565 : ! multiple aerosol representations in the same sim
566 : ! such as MAM and CARMA
567 :
568 : class(aerosol_optics), pointer :: aero_optics
569 :
570 : real(r8) :: dopaer(pcols)
571 : real(r8) :: mass(pcols,pver)
572 : real(r8) :: air_density(pcols,pver)
573 :
574 0 : real(r8), allocatable :: pext(:) ! parameterized specific extinction (m2/kg)
575 0 : real(r8), allocatable :: pabs(:) ! parameterized specific absorption (m2/kg)
576 0 : real(r8), allocatable :: palb(:) ! parameterized single scattering albedo
577 0 : real(r8), allocatable :: pasm(:) ! parameterized asymmetry factor
578 :
579 : real(r8) :: relh(pcols,pver)
580 : real(r8) :: sate(pcols,pver) ! saturation vapor pressure
581 : real(r8) :: satq(pcols,pver) ! saturation specific humidity
582 : real(r8) :: sulfwtpct(pcols,pver) ! sulf weight percent
583 :
584 : character(len=ot_length) :: opticstype
585 : integer :: iaermod
586 :
587 : real(r8) :: aodvis(pcols) ! extinction optical depth in vis
588 : real(r8) :: aoduv(pcols) ! extinction optical depth in uv
589 : real(r8) :: aodnir(pcols) ! extinction optical depth in nir
590 : real(r8) :: absorb(pcols,pver)
591 : real(r8) :: aodabs(pcols) ! absorption optical depth
592 :
593 : real(r8) :: aodabsbc(pcols) ! absorption optical depth of BC
594 :
595 : real(r8) :: aodtot(pcols)
596 :
597 : real(r8) :: extinct(pcols,pver)
598 : real(r8) :: extinctnir(pcols,pver)
599 : real(r8) :: extinctuv(pcols,pver)
600 :
601 : real(r8) :: asymvis(pcols) ! asymmetry factor * optical depth
602 : real(r8) :: asymext(pcols,pver) ! asymmetry factor * extinction
603 :
604 : real(r8) :: wetvol(pcols,pver)
605 : real(r8) :: watervol(pcols,pver)
606 :
607 : real(r8) :: vol(pcols)
608 : real(r8) :: dustvol(pcols)
609 :
610 : real(r8) :: scatdust(pcols)
611 : real(r8) :: absdust(pcols)
612 : real(r8) :: dustaodbin(pcols)
613 :
614 : real(r8) :: scatbc(pcols)
615 : real(r8) :: absbc(pcols)
616 :
617 : real(r8) :: scatpom(pcols)
618 : real(r8) :: abspom(pcols)
619 :
620 : real(r8) :: scatsslt(pcols)
621 : real(r8) :: abssslt(pcols)
622 :
623 : real(r8) :: scatsoa(pcols)
624 : real(r8) :: abssoa(pcols)
625 :
626 : real(r8) :: scatsulf(pcols)
627 : real(r8) :: abssulf(pcols)
628 :
629 : real(r8) :: burden(pcols)
630 : real(r8) :: burdendust(pcols), burdenso4(pcols), burdenbc(pcols), &
631 : burdenpom(pcols), burdensoa(pcols), burdenseasalt(pcols)
632 :
633 : real(r8) :: hygrodust(pcols), hygrosulf(pcols), hygrobc(pcols), &
634 : hygropom(pcols), hygrosoa(pcols), hygrosslt(pcols)
635 :
636 : real(r8) :: aodbin(pcols)
637 :
638 0 : complex(r8), pointer :: specrefindex(:) ! species refractive index
639 :
640 : class(aerosol_state), pointer :: aerostate
641 : class(aerosol_properties), pointer :: aeroprops
642 :
643 : real(r8) :: specdens
644 : character(len=32) :: spectype ! species type
645 0 : real(r8), pointer :: specmmr(:,:)
646 : real(r8) :: hygro_aer !
647 :
648 : real(r8) :: scath2o, absh2o, sumscat, sumabs, sumhygro
649 :
650 : real(r8) :: aodc ! aod of component
651 :
652 : ! total species AOD
653 : real(r8) :: dustaod(pcols), sulfaod(pcols), bcaod(pcols), &
654 : pomaod(pcols), soaaod(pcols), ssltaod(pcols)
655 :
656 : real(r8) :: aodvisst(pcols) ! stratospheric extinction optical depth
657 : real(r8) :: aoduvst(pcols) ! stratospheric extinction optical depth in uv
658 : real(r8) :: aodnirst(pcols) ! stratospheric extinction optical depth in nir
659 : real(r8) :: ssavis(pcols)
660 : integer :: troplev(pcols)
661 :
662 : integer :: i, k
663 :
664 0 : nullify(aero_optics)
665 :
666 0 : lchnk = state%lchnk
667 0 : ncol = state%ncol
668 :
669 : !REMOVECAM - no longer need this when CAM is retired and pcols no longer exists
670 0 : troplev(:) = 0
671 : !REMOVECAM_END
672 0 : call tropopause_findChemTrop(state, troplev)
673 :
674 0 : mass(:ncol,:) = state%pdeldry(:ncol,:)*rga
675 0 : air_density(:ncol,:) = state%pmid(:ncol,:)/(rair*state%t(:ncol,:))
676 :
677 0 : aodvis = 0._r8
678 0 : aodnir = 0._r8
679 0 : aoduv = 0._r8
680 0 : aodabs = 0._r8
681 0 : absorb = 0._r8
682 0 : aodtot = 0._r8
683 0 : tauxar = 0._r8
684 0 : extinct = 0._r8
685 0 : extinctnir = 0._r8
686 0 : extinctuv = 0._r8
687 0 : asymvis = 0.0_r8
688 0 : asymext = 0.0_r8
689 0 : ssavis = 0.0_r8
690 0 : aodvisst = 0.0_r8
691 0 : aoduvst = 0.0_r8
692 0 : aodnirst = 0.0_r8
693 :
694 0 : burdendust = 0.0_r8
695 0 : burdenso4 = 0.0_r8
696 0 : burdenbc = 0.0_r8
697 0 : burdenpom = 0.0_r8
698 0 : burdensoa = 0.0_r8
699 0 : burdenseasalt = 0.0_r8
700 :
701 0 : aodabsbc = 0.0_r8
702 0 : dustaod = 0.0_r8
703 0 : sulfaod = 0.0_r8
704 0 : pomaod = 0.0_r8
705 0 : soaaod = 0.0_r8
706 0 : bcaod = 0.0_r8
707 0 : ssltaod = 0.0_r8
708 :
709 0 : if (num_aero_models<1) return
710 :
711 0 : allocate(aero_state(num_aero_models), stat=istat)
712 0 : if (istat/=0) then
713 0 : call endrun(prefix//'array allocation error: aero_state')
714 : end if
715 :
716 0 : iaermod = 0
717 0 : if (modal_active) then
718 0 : iaermod = iaermod+1
719 0 : aero_state(iaermod)%obj => modal_aerosol_state( state, pbuf )
720 : end if
721 0 : if (carma_active) then
722 0 : iaermod = iaermod+1
723 0 : aero_state(iaermod)%obj => carma_aerosol_state( state, pbuf )
724 : end if
725 :
726 0 : allocate(pext(ncol), stat=istat)
727 0 : if (istat/=0) then
728 0 : call endrun(prefix//'array allocation error: pext')
729 : end if
730 0 : allocate(pabs(ncol), stat=istat)
731 0 : if (istat/=0) then
732 0 : call endrun(prefix//'array allocation error: pabs')
733 : end if
734 0 : allocate(palb(ncol), stat=istat)
735 0 : if (istat/=0) then
736 0 : call endrun(prefix//'array allocation error: palb')
737 : end if
738 0 : allocate(pasm(ncol), stat=istat)
739 0 : if (istat/=0) then
740 0 : call endrun(prefix//'array allocation error: pasm')
741 : end if
742 :
743 0 : aeromodel: do iaermod = 1,num_aero_models
744 :
745 0 : aeroprops => aero_props(iaermod)%obj
746 0 : aerostate => aero_state(iaermod)%obj
747 :
748 0 : nbins=aeroprops%nbins(list_idx)
749 :
750 0 : sulfwtpct(:ncol,:pver) = aerostate%wgtpct(ncol,pver)
751 0 : call outfld('SULFWTPCT', sulfwtpct(1:ncol,:), ncol, lchnk)
752 :
753 0 : binloop: do ibin = 1, nbins
754 :
755 0 : dustaodbin(:) = 0._r8
756 0 : burden(:) = 0._r8
757 0 : aodbin(:) = 0.0_r8
758 :
759 0 : call aeroprops%optics_params(list_idx, ibin, opticstype=opticstype)
760 :
761 0 : select case (trim(opticstype))
762 : case('modal') ! refractive method
763 : aero_optics=>refractive_aerosol_optics(aeroprops, aerostate, list_idx, ibin, &
764 0 : ncol, pver, nswbands, nlwbands, crefwsw, crefwlw)
765 : case('hygroscopic_coreshell')
766 : ! calculate relative humidity for table lookup into rh grid
767 0 : call qsat(state%t(:ncol,:), state%pmid(:ncol,:), sate(:ncol,:), satq(:ncol,:), ncol, pver)
768 0 : relh(:ncol,:) = state%q(1:ncol,:,1) / satq(:ncol,:)
769 0 : relh(:ncol,:) = max(1.e-20_r8,relh(:ncol,:))
770 : aero_optics=>hygrocoreshell_aerosol_optics(aeroprops, aerostate, list_idx, &
771 0 : ibin, ncol, pver, relh(:ncol,:))
772 : case('hygroscopic_wtp')
773 : aero_optics=>hygrowghtpct_aerosol_optics(aeroprops, aerostate, list_idx, &
774 0 : ibin, ncol, pver, sulfwtpct(:ncol,:))
775 : case default
776 0 : call endrun(prefix//'optics method not recognized')
777 : end select
778 :
779 0 : if (associated(aero_optics)) then
780 :
781 0 : wetvol(:ncol,:pver) = aerostate%wet_volume(aeroprops, list_idx, ibin, ncol, pver)
782 0 : watervol(:ncol,:pver) = aerostate%water_volume(aeroprops, list_idx, ibin, ncol, pver)
783 :
784 0 : wavelength: do iwav = 1, nswbands
785 :
786 0 : vertical: do ilev = 1, pver
787 :
788 0 : call aero_optics%sw_props(ncol, ilev, iwav, pext, pabs, palb, pasm )
789 :
790 0 : call init_diags
791 :
792 0 : column: do icol = 1,ncol
793 0 : dopaer(icol) = pext(icol)*mass(icol,ilev)
794 0 : tauxar(icol,ilev,iwav) = tauxar(icol,ilev,iwav) + dopaer(icol)
795 0 : wa(icol,ilev,iwav) = wa(icol,ilev,iwav) + dopaer(icol)*palb(icol)
796 0 : ga(icol,ilev,iwav) = ga(icol,ilev,iwav) + dopaer(icol)*palb(icol)*pasm(icol)
797 0 : fa(icol,ilev,iwav) = fa(icol,ilev,iwav) + dopaer(icol)*palb(icol)*pasm(icol)*pasm(icol)
798 :
799 0 : call update_diags
800 :
801 : end do column
802 :
803 : end do vertical
804 : end do wavelength
805 :
806 : else
807 0 : call endrun(prefix//'aero_optics object pointer not associated')
808 : end if
809 :
810 0 : deallocate(aero_optics)
811 0 : nullify(aero_optics)
812 :
813 0 : call output_bin_diags
814 :
815 : end do binloop
816 : end do aeromodel
817 :
818 0 : call output_tot_diags
819 :
820 0 : deallocate(pext)
821 0 : deallocate(pabs)
822 0 : deallocate(palb)
823 0 : deallocate(pasm)
824 :
825 0 : do iaermod = 1,num_aero_models
826 0 : deallocate(aero_state(iaermod)%obj)
827 0 : nullify(aero_state(iaermod)%obj)
828 : end do
829 :
830 0 : deallocate(aero_state)
831 :
832 : contains
833 :
834 : !===============================================================================
835 0 : subroutine init_diags
836 0 : dustvol(:ncol) = 0._r8
837 0 : scatdust(:ncol) = 0._r8
838 0 : absdust(:ncol) = 0._r8
839 0 : hygrodust(:ncol) = 0._r8
840 0 : scatsulf(:ncol) = 0._r8
841 0 : abssulf(:ncol) = 0._r8
842 0 : hygrosulf(:ncol) = 0._r8
843 0 : scatbc(:ncol) = 0._r8
844 0 : absbc(:ncol) = 0._r8
845 0 : hygrobc(:ncol) = 0._r8
846 0 : scatpom(:ncol) = 0._r8
847 0 : abspom(:ncol) = 0._r8
848 0 : hygropom(:ncol) = 0._r8
849 0 : scatsoa(:ncol) = 0._r8
850 0 : abssoa(:ncol) = 0._r8
851 0 : hygrosoa(:ncol) = 0._r8
852 0 : scatsslt(:ncol) = 0._r8
853 0 : abssslt(:ncol) = 0._r8
854 0 : hygrosslt(:ncol) = 0._r8
855 0 : end subroutine init_diags
856 :
857 : !===============================================================================
858 0 : subroutine update_diags
859 :
860 : integer :: ispec
861 :
862 0 : if (iwav==idx_uv_diag) then
863 0 : aoduv(icol) = aoduv(icol) + dopaer(icol)
864 0 : extinctuv(icol,ilev) = extinctuv(icol,ilev) + dopaer(icol)*air_density(icol,ilev)/mass(icol,ilev)
865 0 : if (ilev<=troplev(icol)) then
866 0 : aoduvst(icol) = aoduvst(icol) + dopaer(icol)
867 : end if
868 :
869 0 : else if (iwav==idx_sw_diag) then ! vis
870 0 : aodvis(icol) = aodvis(icol) + dopaer(icol)
871 0 : aodabs(icol) = aodabs(icol) + pabs(icol)*mass(icol,ilev)
872 0 : extinct(icol,ilev) = extinct(icol,ilev) + dopaer(icol)*air_density(icol,ilev)/mass(icol,ilev)
873 0 : absorb(icol,ilev) = absorb(icol,ilev) + pabs(icol)*air_density(icol,ilev)
874 0 : ssavis(icol) = ssavis(icol) + dopaer(icol)*palb(icol)
875 0 : asymvis(icol) = asymvis(icol) + dopaer(icol)*pasm(icol)
876 0 : asymext(icol,ilev) = asymext(icol,ilev) + dopaer(icol)*pasm(icol)*air_density(icol,ilev)/mass(icol,ilev)
877 :
878 0 : aodbin(icol) = aodbin(icol) + dopaer(icol)
879 :
880 0 : if (ilev<=troplev(icol)) then
881 0 : aodvisst(icol) = aodvisst(icol) + dopaer(icol)
882 : end if
883 :
884 : ! loop over species ...
885 :
886 0 : do ispec = 1, aeroprops%nspecies(list_idx,ibin)
887 : call aeroprops%get(ibin, ispec, list_ndx=list_idx, density=specdens, &
888 0 : spectype=spectype, refindex_sw=specrefindex, hygro=hygro_aer)
889 0 : call aerostate%get_ambient_mmr(list_idx, ispec, ibin, specmmr)
890 :
891 0 : burden(icol) = burden(icol) + specmmr(icol,ilev)*mass(icol,ilev)
892 :
893 0 : vol(icol) = specmmr(icol,ilev)/specdens
894 :
895 0 : select case ( trim(spectype) )
896 : case('dust')
897 0 : dustvol(icol) = vol(icol)
898 0 : burdendust(icol) = burdendust(icol) + specmmr(icol,ilev)*mass(icol,ilev)
899 0 : scatdust(icol) = vol(icol) * specrefindex(iwav)%re
900 0 : absdust(icol) =-vol(icol) * specrefindex(iwav)%im
901 0 : hygrodust(icol)= vol(icol)*hygro_aer
902 : case('black-c')
903 0 : burdenbc(icol) = burdenbc(icol) + specmmr(icol,ilev)*mass(icol,ilev)
904 0 : scatbc(icol) = vol(icol) * specrefindex(iwav)%re
905 0 : absbc(icol) =-vol(icol) * specrefindex(iwav)%im
906 0 : hygrobc(icol)= vol(icol)*hygro_aer
907 : case('sulfate')
908 0 : burdenso4(icol) = burdenso4(icol) + specmmr(icol,ilev)*mass(icol,ilev)
909 0 : scatsulf(icol) = vol(icol) * specrefindex(iwav)%re
910 0 : abssulf(icol) =-vol(icol) * specrefindex(iwav)%im
911 0 : hygrosulf(icol)= vol(icol)*hygro_aer
912 : case('p-organic')
913 0 : burdenpom(icol) = burdenpom(icol) + specmmr(icol,ilev)*mass(icol,ilev)
914 0 : scatpom(icol) = vol(icol) * specrefindex(iwav)%re
915 0 : abspom(icol) =-vol(icol) * specrefindex(iwav)%im
916 0 : hygropom(icol)= vol(icol)*hygro_aer
917 : case('s-organic')
918 0 : burdensoa(icol) = burdensoa(icol) + specmmr(icol,ilev)*mass(icol,ilev)
919 0 : scatsoa(icol) = vol(icol) * specrefindex(iwav)%re
920 0 : abssoa(icol) = -vol(icol) * specrefindex(iwav)%im
921 0 : hygrosoa(icol)= vol(icol)*hygro_aer
922 : case('seasalt')
923 0 : burdenseasalt(icol) = burdenseasalt(icol) + specmmr(icol,ilev)*mass(icol,ilev)
924 0 : scatsslt(icol) = vol(icol) * specrefindex(iwav)%re
925 0 : abssslt(icol) = -vol(icol) * specrefindex(iwav)%im
926 0 : hygrosslt(icol)= vol(icol)*hygro_aer
927 : end select
928 : end do
929 :
930 0 : if (wetvol(icol,ilev)>1.e-40_r8 .and. vol(icol)>0._r8) then
931 :
932 0 : dustaodbin(icol) = dustaodbin(icol) + dopaer(icol)*dustvol(icol)/wetvol(icol,ilev)
933 :
934 : ! partition optical depth into contributions from each constituent
935 : ! assume contribution is proportional to refractive index X volume
936 :
937 0 : scath2o = watervol(icol,ilev)*crefwsw(iwav)%re
938 0 : absh2o = -watervol(icol,ilev)*crefwsw(iwav)%im
939 0 : sumscat = scatsulf(icol) + scatpom(icol) + scatsoa(icol) + scatbc(icol) + &
940 0 : scatdust(icol) + scatsslt(icol) + scath2o
941 0 : sumabs = abssulf(icol) + abspom(icol) + abssoa(icol) + absbc(icol) + &
942 0 : absdust(icol) + abssslt(icol) + absh2o
943 0 : sumhygro = hygrosulf(icol) + hygropom(icol) + hygrosoa(icol) + hygrobc(icol) + &
944 0 : hygrodust(icol) + hygrosslt(icol)
945 :
946 0 : scatdust(icol) = (scatdust(icol) + scath2o*hygrodust(icol)/sumhygro)/sumscat
947 0 : absdust(icol) = (absdust(icol) + absh2o*hygrodust(icol)/sumhygro)/sumabs
948 :
949 0 : scatsulf(icol) = (scatsulf(icol) + scath2o*hygrosulf(icol)/sumhygro)/sumscat
950 0 : abssulf(icol) = (abssulf(icol) + absh2o*hygrosulf(icol)/sumhygro)/sumabs
951 :
952 0 : scatpom(icol) = (scatpom(icol) + scath2o*hygropom(icol)/sumhygro)/sumscat
953 0 : abspom(icol) = (abspom(icol) + absh2o*hygropom(icol)/sumhygro)/sumabs
954 :
955 0 : scatsoa(icol) = (scatsoa(icol) + scath2o*hygrosoa(icol)/sumhygro)/sumscat
956 0 : abssoa(icol) = (abssoa(icol) + absh2o*hygrosoa(icol)/sumhygro)/sumabs
957 :
958 0 : scatbc(icol)= (scatbc(icol) + scath2o*hygrobc(icol)/sumhygro)/sumscat
959 0 : absbc(icol) = (absbc(icol) + absh2o*hygrobc(icol)/sumhygro)/sumabs
960 :
961 0 : scatsslt(icol) = (scatsslt(icol) + scath2o*hygrosslt(icol)/sumhygro)/sumscat
962 0 : abssslt(icol) = (abssslt(icol) + absh2o*hygrosslt(icol)/sumhygro)/sumabs
963 :
964 :
965 0 : aodabsbc(icol) = aodabsbc(icol) + absbc(icol)*dopaer(icol)*(1.0_r8-palb(icol))
966 :
967 :
968 :
969 0 : aodc = (absdust(icol)*(1.0_r8 - palb(icol)) + palb(icol)*scatdust(icol))*dopaer(icol)
970 0 : dustaod(icol) = dustaod(icol) + aodc
971 :
972 0 : aodc = (abssulf(icol)*(1.0_r8 - palb(icol)) + palb(icol)*scatsulf(icol))*dopaer(icol)
973 0 : sulfaod(icol) = sulfaod(icol) + aodc
974 :
975 0 : aodc = (abspom(icol)*(1.0_r8 - palb(icol)) + palb(icol)*scatpom(icol))*dopaer(icol)
976 0 : pomaod(icol) = pomaod(icol) + aodc
977 :
978 0 : aodc = (abssoa(icol)*(1.0_r8 - palb(icol)) + palb(icol)*scatsoa(icol))*dopaer(icol)
979 0 : soaaod(icol) = soaaod(icol) + aodc
980 :
981 0 : aodc = (absbc(icol)*(1.0_r8 - palb(icol)) + palb(icol)*scatbc(icol))*dopaer(icol)
982 0 : bcaod(icol) = bcaod(icol) + aodc
983 :
984 0 : aodc = (abssslt(icol)*(1.0_r8 - palb(icol)) + palb(icol)*scatsslt(icol))*dopaer(icol)
985 0 : ssltaod(icol) = ssltaod(icol) + aodc
986 :
987 : end if
988 0 : else if (iwav==idx_nir_diag) then
989 0 : aodnir(icol) = aodnir(icol) + dopaer(icol)
990 0 : extinctnir(icol,ilev) = extinctnir(icol,ilev) + dopaer(icol)*air_density(icol,ilev)/mass(icol,ilev)
991 :
992 0 : if (ilev<=troplev(icol)) then
993 0 : aodnirst(icol) = aodnirst(icol) + dopaer(icol)
994 : end if
995 :
996 : end if
997 :
998 0 : aodtot(icol) = aodtot(icol) + dopaer(icol)
999 :
1000 0 : end subroutine update_diags
1001 :
1002 : !===============================================================================
1003 0 : subroutine output_bin_diags
1004 :
1005 : integer :: icol
1006 :
1007 0 : if (list_idx == 0) then
1008 :
1009 0 : call outfld(burdendn_fields(iaermod)%name(ibin), burden, pcols, lchnk)
1010 0 : call outfld(aoddustdn_fields(iaermod)%name(ibin), dustaodbin, pcols, lchnk)
1011 0 : call outfld(aodbindn_fields(iaermod)%name(ibin), aodbin, pcols, lchnk)
1012 :
1013 0 : do icol = 1, nnite
1014 0 : burden(idxnite(icol)) = fillvalue
1015 0 : aodbin(idxnite(icol)) = fillvalue
1016 0 : dustaodbin(idxnite(icol)) = fillvalue
1017 : end do
1018 :
1019 0 : call outfld(burden_fields(iaermod)%name(ibin), burden, pcols, lchnk)
1020 0 : call outfld(aoddust_fields(iaermod)%name(ibin), dustaodbin, pcols, lchnk)
1021 0 : call outfld(aodbin_fields(iaermod)%name(ibin), aodbin, pcols, lchnk)
1022 :
1023 : endif
1024 :
1025 0 : end subroutine output_bin_diags
1026 :
1027 : !===============================================================================
1028 0 : subroutine output_tot_diags
1029 :
1030 : integer :: icol
1031 :
1032 0 : call outfld('AODUVdn'//diag(list_idx), aoduv, pcols, lchnk)
1033 0 : call outfld('AODVISdn'//diag(list_idx), aodvis, pcols, lchnk)
1034 0 : call outfld('AODABSdn'//diag(list_idx), aodabs, pcols, lchnk)
1035 0 : call outfld('AODNIRdn'//diag(list_idx), aodnir, pcols, lchnk)
1036 0 : call outfld('AODTOTdn'//diag(list_idx), aodtot, pcols, lchnk)
1037 0 : call outfld('EXTINCTUVdn'//diag(list_idx), extinctuv, pcols, lchnk)
1038 0 : call outfld('EXTINCTNIRdn'//diag(list_idx), extinctnir, pcols, lchnk)
1039 0 : call outfld('EXTINCTdn'//diag(list_idx), extinct, pcols, lchnk)
1040 0 : call outfld('ABSORBdn'//diag(list_idx), absorb, pcols, lchnk)
1041 0 : call outfld('EXTxASYMdn'//diag(list_idx), asymext, pcols, lchnk)
1042 0 : call outfld('AODxASYMdn'//diag(list_idx), asymvis, pcols, lchnk)
1043 :
1044 0 : call outfld('AODVISstdn'//diag(list_idx), aodvisst,pcols, lchnk)
1045 0 : call outfld('AODUVstdn'//diag(list_idx), aoduvst, pcols, lchnk)
1046 0 : call outfld('AODNIRstdn'//diag(list_idx), aodnirst,pcols, lchnk)
1047 :
1048 0 : do icol = 1, nnite
1049 0 : aodvis(idxnite(icol)) = fillvalue
1050 0 : aodnir(idxnite(icol)) = fillvalue
1051 0 : aoduv(idxnite(icol)) = fillvalue
1052 0 : aodabs(idxnite(icol)) = fillvalue
1053 0 : aodtot(idxnite(icol)) = fillvalue
1054 0 : extinct(idxnite(icol),:) = fillvalue
1055 0 : extinctnir(idxnite(icol),:) = fillvalue
1056 0 : extinctuv(idxnite(icol),:) = fillvalue
1057 0 : absorb(idxnite(icol),:) = fillvalue
1058 0 : asymext(idxnite(icol),:) = fillvalue
1059 0 : asymvis(idxnite(icol)) = fillvalue
1060 0 : aodabs(idxnite(icol)) = fillvalue
1061 0 : aodvisst(idxnite(icol)) = fillvalue
1062 0 : aoduvst(idxnite(icol)) = fillvalue
1063 0 : aodnirst(idxnite(icol)) = fillvalue
1064 : end do
1065 :
1066 0 : call outfld('AODUV'//diag(list_idx), aoduv, pcols, lchnk)
1067 0 : call outfld('AODVIS'//diag(list_idx), aodvis, pcols, lchnk)
1068 0 : call outfld('AODABS'//diag(list_idx), aodabs, pcols, lchnk)
1069 0 : call outfld('AODNIR'//diag(list_idx), aodnir, pcols, lchnk)
1070 0 : call outfld('AODTOT'//diag(list_idx), aodtot, pcols, lchnk)
1071 0 : call outfld('EXTINCTUV'//diag(list_idx), extinctuv, pcols, lchnk)
1072 0 : call outfld('EXTINCTNIR'//diag(list_idx), extinctnir, pcols, lchnk)
1073 0 : call outfld('EXTINCT'//diag(list_idx), extinct, pcols, lchnk)
1074 0 : call outfld('ABSORB'//diag(list_idx), absorb, pcols, lchnk)
1075 0 : call outfld('EXTxASYM'//diag(list_idx), asymext, pcols, lchnk)
1076 0 : call outfld('AODxASYM'//diag(list_idx), asymvis, pcols, lchnk)
1077 0 : call outfld('AODVISst'//diag(list_idx), aodvisst,pcols, lchnk)
1078 0 : call outfld('AODUVst'//diag(list_idx), aoduvst, pcols, lchnk)
1079 0 : call outfld('AODNIRst'//diag(list_idx), aodnirst,pcols, lchnk)
1080 :
1081 : ! These diagnostics are output only for climate list
1082 0 : if (list_idx == 0) then
1083 0 : do icol = 1, ncol
1084 0 : if (aodvis(icol) > 1.e-10_r8) then
1085 0 : ssavis(icol) = ssavis(icol)/aodvis(icol)
1086 : else
1087 0 : ssavis(icol) = 0.925_r8
1088 : endif
1089 : end do
1090 0 : call outfld('SSAVISdn', ssavis, pcols, lchnk)
1091 :
1092 0 : call outfld('BURDENDUSTdn', burdendust, pcols, lchnk)
1093 0 : call outfld('BURDENSO4dn' , burdenso4, pcols, lchnk)
1094 0 : call outfld('BURDENPOMdn' , burdenpom, pcols, lchnk)
1095 0 : call outfld('BURDENSOAdn' , burdensoa, pcols, lchnk)
1096 0 : call outfld('BURDENBCdn' , burdenbc, pcols, lchnk)
1097 0 : call outfld('BURDENSEASALTdn', burdenseasalt, pcols, lchnk)
1098 :
1099 0 : call outfld('AODABSBCdn', aodabsbc, pcols, lchnk)
1100 :
1101 0 : call outfld('AODDUSTdn', dustaod, pcols, lchnk)
1102 0 : call outfld('AODSO4dn', sulfaod, pcols, lchnk)
1103 0 : call outfld('AODPOMdn', pomaod, pcols, lchnk)
1104 0 : call outfld('AODSOAdn', soaaod, pcols, lchnk)
1105 0 : call outfld('AODBCdn', bcaod, pcols, lchnk)
1106 0 : call outfld('AODSSdn', ssltaod, pcols, lchnk)
1107 :
1108 :
1109 0 : do icol = 1, nnite
1110 :
1111 0 : ssavis(idxnite(icol)) = fillvalue
1112 0 : asymvis(idxnite(icol)) = fillvalue
1113 :
1114 0 : burdendust(idxnite(icol)) = fillvalue
1115 0 : burdenso4(idxnite(icol)) = fillvalue
1116 0 : burdenpom(idxnite(icol)) = fillvalue
1117 0 : burdensoa(idxnite(icol)) = fillvalue
1118 0 : burdenbc(idxnite(icol)) = fillvalue
1119 0 : burdenseasalt(idxnite(icol)) = fillvalue
1120 0 : aodabsbc(idxnite(icol)) = fillvalue
1121 :
1122 0 : dustaod(idxnite(icol)) = fillvalue
1123 0 : sulfaod(idxnite(icol)) = fillvalue
1124 0 : pomaod(idxnite(icol)) = fillvalue
1125 0 : soaaod(idxnite(icol)) = fillvalue
1126 0 : bcaod(idxnite(icol)) = fillvalue
1127 0 : ssltaod(idxnite(icol)) = fillvalue
1128 :
1129 : end do
1130 :
1131 0 : call outfld('SSAVIS', ssavis, pcols, lchnk)
1132 0 : call outfld('AODxASYM', asymvis, pcols, lchnk)
1133 0 : call outfld('BURDENDUST', burdendust, pcols, lchnk)
1134 0 : call outfld('BURDENSO4' , burdenso4, pcols, lchnk)
1135 0 : call outfld('BURDENPOM' , burdenpom, pcols, lchnk)
1136 0 : call outfld('BURDENSOA' , burdensoa, pcols, lchnk)
1137 0 : call outfld('BURDENBC' , burdenbc, pcols, lchnk)
1138 0 : call outfld('BURDENSEASALT', burdenseasalt, pcols, lchnk)
1139 0 : call outfld('AODABSBC', aodabsbc, pcols, lchnk)
1140 0 : call outfld('AODDUST', dustaod, pcols, lchnk)
1141 0 : call outfld('AODSO4', sulfaod, pcols, lchnk)
1142 0 : call outfld('AODPOM', pomaod, pcols, lchnk)
1143 0 : call outfld('AODSOA', soaaod, pcols, lchnk)
1144 0 : call outfld('AODBC', bcaod, pcols, lchnk)
1145 0 : call outfld('AODSS', ssltaod, pcols, lchnk)
1146 :
1147 : end if
1148 :
1149 0 : end subroutine output_tot_diags
1150 :
1151 : end subroutine aerosol_optics_cam_sw
1152 :
1153 : !===============================================================================
1154 0 : subroutine aerosol_optics_cam_lw(list_idx, state, pbuf, tauxar)
1155 :
1156 : ! calculates aerosol lw radiative properties
1157 :
1158 : integer, intent(in) :: list_idx ! index of the climate or a diagnostic list
1159 : type(physics_state), intent(in), target :: state ! state variables
1160 :
1161 : type(physics_buffer_desc), pointer :: pbuf(:)
1162 :
1163 : real(r8), intent(inout) :: tauxar(pcols,pver,nlwbands) ! layer absorption optical depth
1164 :
1165 :
1166 : real(r8) :: dopaer(pcols)
1167 : real(r8) :: mass(pcols,pver)
1168 :
1169 : character(len=*), parameter :: prefix = 'aerosol_optics_cam_lw: '
1170 :
1171 : integer :: ibin, nbins
1172 : integer :: iwav, ilev
1173 : integer :: ncol, icol, istat
1174 :
1175 0 : type(aero_state_t), allocatable :: aero_state(:) ! array of aerosol state objects to allow for
1176 : ! multiple aerosol representations in the same sim
1177 : ! such as MAM and CARMA
1178 :
1179 : class(aerosol_optics), pointer :: aero_optics
1180 : class(aerosol_state), pointer :: aerostate
1181 : class(aerosol_properties), pointer :: aeroprops
1182 :
1183 0 : real(r8), allocatable :: pabs(:)
1184 :
1185 : real(r8) :: relh(pcols,pver)
1186 : real(r8) :: sate(pcols,pver) ! saturation vapor pressure
1187 : real(r8) :: satq(pcols,pver) ! saturation specific humidity
1188 : real(r8) :: sulfwtpct(pcols,pver) ! sulf weight percent
1189 :
1190 : character(len=32) :: opticstype
1191 : integer :: iaermod
1192 :
1193 : real(r8) :: lwabs(pcols,pver)
1194 0 : lwabs = 0._r8
1195 0 : tauxar = 0._r8
1196 :
1197 0 : nullify(aero_optics)
1198 :
1199 0 : allocate(aero_state(num_aero_models), stat=istat)
1200 0 : if (istat/=0) then
1201 0 : call endrun(prefix//'array allocation error: aero_state')
1202 : end if
1203 :
1204 0 : iaermod = 0
1205 0 : if (modal_active) then
1206 0 : iaermod = iaermod+1
1207 0 : aero_state(iaermod)%obj => modal_aerosol_state( state, pbuf )
1208 : end if
1209 0 : if (carma_active) then
1210 0 : iaermod = iaermod+1
1211 0 : aero_state(iaermod)%obj => carma_aerosol_state( state, pbuf )
1212 : end if
1213 :
1214 0 : ncol = state%ncol
1215 :
1216 0 : mass(:ncol,:) = state%pdeldry(:ncol,:)*rga
1217 :
1218 0 : allocate(pabs(ncol), stat=istat)
1219 0 : if (istat/=0) then
1220 0 : call endrun(prefix//'array allocation error: pabs')
1221 : end if
1222 :
1223 0 : aeromodel: do iaermod = 1,num_aero_models
1224 :
1225 0 : aeroprops => aero_props(iaermod)%obj
1226 0 : aerostate => aero_state(iaermod)%obj
1227 :
1228 0 : nbins=aero_props(iaermod)%obj%nbins(list_idx)
1229 :
1230 0 : sulfwtpct(:ncol,:pver) = aerostate%wgtpct(ncol,pver)
1231 :
1232 0 : binloop: do ibin = 1, nbins
1233 :
1234 0 : call aeroprops%optics_params(list_idx, ibin, opticstype=opticstype)
1235 :
1236 0 : select case (trim(opticstype))
1237 : case('modal') ! refractive method
1238 : aero_optics=>refractive_aerosol_optics(aeroprops, aerostate, list_idx, ibin, &
1239 0 : ncol, pver, nswbands, nlwbands, crefwsw, crefwlw)
1240 : case('hygroscopic_coreshell')
1241 : ! calculate relative humidity for table lookup into rh grid
1242 0 : call qsat(state%t(:ncol,:), state%pmid(:ncol,:), sate(:ncol,:), satq(:ncol,:), ncol, pver)
1243 0 : relh(:ncol,:) = state%q(1:ncol,:,1) / satq(:ncol,:)
1244 0 : relh(:ncol,:) = max(1.e-20_r8,relh(:ncol,:))
1245 : aero_optics=>hygrocoreshell_aerosol_optics(aeroprops, aerostate, list_idx, &
1246 0 : ibin, ncol, pver, relh(:ncol,:))
1247 : case('hygroscopic_wtp')
1248 : aero_optics=>hygrowghtpct_aerosol_optics(aeroprops, aerostate, list_idx, &
1249 0 : ibin, ncol, pver, sulfwtpct(:ncol,:))
1250 : case default
1251 0 : call endrun(prefix//'optics method not recognized')
1252 : end select
1253 :
1254 0 : if (associated(aero_optics)) then
1255 :
1256 0 : wavelength: do iwav = 1, nlwbands
1257 :
1258 0 : vertical: do ilev = 1, pver
1259 0 : call aero_optics%lw_props(ncol, ilev, iwav, pabs )
1260 :
1261 0 : column: do icol = 1, ncol
1262 0 : dopaer(icol) = pabs(icol)*mass(icol,ilev)
1263 0 : tauxar(icol,ilev,iwav) = tauxar(icol,ilev,iwav) + dopaer(icol)
1264 0 : lwabs(icol,ilev) = lwabs(icol,ilev) + pabs(icol)
1265 : end do column
1266 :
1267 : end do vertical
1268 :
1269 : end do wavelength
1270 :
1271 : else
1272 0 : call endrun(prefix//'aero_optics object pointer not associated')
1273 : end if
1274 :
1275 0 : deallocate(aero_optics)
1276 0 : nullify(aero_optics)
1277 :
1278 : end do binloop
1279 : end do aeromodel
1280 :
1281 0 : call outfld('TOTABSLW'//diag(list_idx), lwabs(:,:), pcols, state%lchnk)
1282 :
1283 0 : if (lw10um_indx>0) then
1284 0 : call outfld('AODABSLW'//diag(list_idx), tauxar(:,:,lw10um_indx), pcols, state%lchnk)
1285 : end if
1286 :
1287 0 : deallocate(pabs)
1288 :
1289 0 : do iaermod = 1,num_aero_models
1290 0 : deallocate(aero_state(iaermod)%obj)
1291 0 : nullify(aero_state(iaermod)%obj)
1292 : end do
1293 :
1294 0 : deallocate(aero_state)
1295 :
1296 0 : end subroutine aerosol_optics_cam_lw
1297 :
1298 : !===============================================================================
1299 : ! Private routines
1300 : !===============================================================================
1301 :
1302 0 : subroutine read_water_refindex(infilename)
1303 : use cam_pio_utils, only: cam_pio_openfile
1304 : use pio, only: file_desc_t, var_desc_t, pio_inq_dimlen, pio_inq_dimid, pio_inq_varid, &
1305 : pio_get_var, PIO_NOWRITE, pio_closefile, pio_noerr
1306 :
1307 :
1308 : ! read water refractive index file and set module data
1309 :
1310 : character*(*), intent(in) :: infilename ! modal optics filename
1311 :
1312 : ! Local variables
1313 :
1314 : integer :: i, ierr
1315 : type(file_desc_t) :: ncid ! pio file handle
1316 : integer :: did ! dimension ids
1317 : integer :: dimlen ! dimension lengths
1318 : type(var_desc_t) :: vid ! variable ids
1319 : real(r8) :: refrwsw(nswbands), refiwsw(nswbands) ! real, imaginary ref index for water visible
1320 : real(r8) :: refrwlw(nlwbands), refiwlw(nlwbands) ! real, imaginary ref index for water infrared
1321 :
1322 : character(len=*), parameter :: prefix = 'read_water_refindex: '
1323 : !----------------------------------------------------------------------------
1324 :
1325 : ! open file
1326 0 : call cam_pio_openfile(ncid, infilename, PIO_NOWRITE)
1327 :
1328 : ! inquire dimensions. Check that file values match parameter values.
1329 :
1330 0 : ierr = pio_inq_dimid(ncid, 'lw_band', did)
1331 0 : if (ierr /= pio_noerr ) then
1332 0 : call endrun(prefix//'pio_inq_dimid lw_band')
1333 : end if
1334 0 : ierr = pio_inq_dimlen(ncid, did, dimlen)
1335 0 : if (ierr /= pio_noerr ) then
1336 0 : call endrun(prefix//'pio_inq_dimlen lw_band')
1337 : end if
1338 0 : if (dimlen /= nlwbands) then
1339 0 : write(iulog,*) 'lw_band len=', dimlen, ' from ', infilename, ' ne nlwbands=', nlwbands
1340 0 : call endrun(prefix//'bad lw_band value')
1341 : endif
1342 :
1343 0 : ierr = pio_inq_dimid(ncid, 'sw_band', did)
1344 0 : if (ierr /= pio_noerr ) then
1345 0 : call endrun(prefix//'pio_inq_dimid sw_band')
1346 : end if
1347 0 : ierr = pio_inq_dimlen(ncid, did, dimlen)
1348 0 : if (ierr /= pio_noerr ) then
1349 0 : call endrun(prefix//'pio_inq_dimlen sw_band')
1350 : end if
1351 0 : if (dimlen /= nswbands) then
1352 0 : write(iulog,*) 'sw_band len=', dimlen, ' from ', infilename, ' ne nswbands=', nswbands
1353 0 : call endrun(prefix//'bad sw_band value')
1354 : endif
1355 :
1356 : ! read variables
1357 0 : ierr = pio_inq_varid(ncid, 'refindex_real_water_sw', vid)
1358 0 : if (ierr /= pio_noerr ) then
1359 0 : call endrun(prefix//'pio_inq_varid refindex_real_water_sw')
1360 : end if
1361 0 : ierr = pio_get_var(ncid, vid, refrwsw)
1362 0 : if (ierr /= pio_noerr ) then
1363 0 : call endrun(prefix//'pio_get_var refrwsw')
1364 : end if
1365 :
1366 0 : ierr = pio_inq_varid(ncid, 'refindex_im_water_sw', vid)
1367 0 : if (ierr /= pio_noerr ) then
1368 0 : call endrun(prefix//'pio_inq_varid refindex_im_water_sw')
1369 : end if
1370 0 : ierr = pio_get_var(ncid, vid, refiwsw)
1371 0 : if (ierr /= pio_noerr ) then
1372 0 : call endrun(prefix//'pio_get_var refiwsw')
1373 : end if
1374 :
1375 0 : ierr = pio_inq_varid(ncid, 'refindex_real_water_lw', vid)
1376 0 : if (ierr /= pio_noerr ) then
1377 0 : call endrun(prefix//'pio_inq_varid refindex_real_water_lw')
1378 : end if
1379 0 : ierr = pio_get_var(ncid, vid, refrwlw)
1380 0 : if (ierr /= pio_noerr ) then
1381 0 : call endrun(prefix//'pio_get_var refrwlw')
1382 : end if
1383 :
1384 0 : ierr = pio_inq_varid(ncid, 'refindex_im_water_lw', vid)
1385 0 : if (ierr /= pio_noerr ) then
1386 0 : call endrun(prefix//'pio_inq_varid refindex_im_water_lw')
1387 : end if
1388 0 : ierr = pio_get_var(ncid, vid, refiwlw)
1389 0 : if (ierr /= pio_noerr ) then
1390 0 : call endrun(prefix//'pio_get_var refiwlw')
1391 : end if
1392 :
1393 : ! set complex representation of refractive indices as module data
1394 0 : do i = 1, nswbands
1395 0 : crefwsw(i) = cmplx(refrwsw(i), abs(refiwsw(i)), kind=r8)
1396 : end do
1397 0 : do i = 1, nlwbands
1398 0 : crefwlw(i) = cmplx(refrwlw(i), abs(refiwlw(i)), kind=r8)
1399 : end do
1400 :
1401 0 : call pio_closefile(ncid)
1402 :
1403 0 : end subroutine read_water_refindex
1404 :
1405 0 : end module aerosol_optics_cam
|
