Line data Source code
1 : module phys_prop
2 :
3 : ! Properties of aerosols that are used by radiation and other parameterizations.
4 :
5 : ! *****N.B.*****
6 : ! This module is a utility used by the rad_constituents module. The properties stored
7 : ! here are meant to be accessed via that module. This module knows nothing about how
8 : ! this data is associated with the constituents that are radiatively active or those that
9 : ! are being used for diagnostic calculations. That is the responsibility of the
10 : ! rad_constituents module.
11 :
12 : use shr_kind_mod, only: r8 => shr_kind_r8
13 : use spmd_utils, only: masterproc
14 : use radconstants, only: nlwbands, nswbands, idx_sw_diag
15 : use ioFileMod, only: getfil
16 : use cam_pio_utils, only: cam_pio_openfile
17 : use pio, only: file_desc_t, var_desc_t, pio_get_var, pio_inq_varid, &
18 : pio_inq_dimlen, pio_inq_dimid , pio_nowrite, pio_closefile, &
19 : pio_seterrorhandling, PIO_BCAST_ERROR, PIO_INTERNAL_ERROR, PIO_NOERR
20 :
21 : use cam_logfile, only: iulog
22 : use cam_abortutils, only: endrun
23 :
24 : implicit none
25 : private
26 : save
27 :
28 : integer, parameter, public :: ot_length = 32
29 :
30 : public :: &
31 : physprop_accum_unique_files, &! Make a list of the unique set of files that contain properties
32 : ! This is an initialization step that must be done before calling physprop_init
33 : physprop_init, &! Initialization -- read the input datasets
34 : physprop_get_id, &! Return ID used to access the property data from the input files
35 : physprop_get ! Return data for specified ID
36 :
37 : ! Data from one input dataset is stored in a structure of type(physprop_type).
38 : type :: physprop_type
39 : character(len=256) :: sourcefile ! Absolute pathname of data file.
40 : character(len=ot_length) :: opticsmethod ! one of {hygro,nonhygro}
41 :
42 : ! for hygroscopic species of externally mixed aerosols
43 : real(r8), pointer :: sw_hygro_ext(:,:)
44 : real(r8), pointer :: sw_hygro_ssa(:,:)
45 : real(r8), pointer :: sw_hygro_asm(:,:)
46 : real(r8), pointer :: lw_hygro_abs(:,:)
47 :
48 : ! for nonhygroscopic species of externally mixed aerosols
49 : real(r8), pointer :: sw_nonhygro_ext(:)
50 : real(r8), pointer :: sw_nonhygro_ssa(:)
51 : real(r8), pointer :: sw_nonhygro_asm(:)
52 : real(r8), pointer :: sw_nonhygro_scat(:)
53 : real(r8), pointer :: sw_nonhygro_ascat(:)
54 : real(r8), pointer :: lw_abs(:)
55 :
56 : ! complex refractive index
57 : complex(r8), pointer :: refindex_aer_sw(:)
58 : complex(r8), pointer :: refindex_aer_lw(:)
59 :
60 : ! for radius-dependent mass-specific quantities
61 : real(r8), pointer :: r_sw_ext(:,:)
62 : real(r8), pointer :: r_sw_scat(:,:)
63 : real(r8), pointer :: r_sw_ascat(:,:)
64 : real(r8), pointer :: r_lw_abs(:,:)
65 : real(r8), pointer :: mu(:)
66 :
67 : ! for modal optics
68 : real(r8), pointer :: extpsw(:,:,:,:) ! specific extinction
69 : real(r8), pointer :: abspsw(:,:,:,:) ! specific absorption
70 : real(r8), pointer :: asmpsw(:,:,:,:) ! asymmetry factor
71 : real(r8), pointer :: absplw(:,:,:,:) ! specific absorption
72 : real(r8), pointer :: refrtabsw(:,:) ! table of real refractive indices for aerosols visible
73 : real(r8), pointer :: refitabsw(:,:) ! table of imag refractive indices for aerosols visible
74 : real(r8), pointer :: refrtablw(:,:) ! table of real refractive indices for aerosols infrared
75 : real(r8), pointer :: refitablw(:,:) ! table of imag refractive indices for aerosols infrared
76 :
77 : ! for core/shell optics
78 : real(r8), pointer :: extpsw2(:,:) ! specific extinction
79 : real(r8), pointer :: abspsw2(:,:) ! specific absorption
80 : real(r8), pointer :: asmpsw2(:,:) ! asymmetry factor
81 : real(r8), pointer :: absplw2(:,:) ! specific absorption
82 : real(r8), pointer :: corefrac(:) ! table of real refractive indices for aerosols visible
83 : integer :: nfraC ! number of Chebyshev coefficients
84 :
85 : ! for hygroscopic species of pure sulfate
86 : real(r8), pointer :: sw_hygro_ext_wtp(:,:)
87 : real(r8), pointer :: sw_hygro_ssa_wtp(:,:)
88 : real(r8), pointer :: sw_hygro_asm_wtp(:,:)
89 : real(r8), pointer :: lw_hygro_abs_wtp(:,:)
90 : real(r8), pointer :: wgtpct (:) ! weight percent!
91 : integer :: nwtp ! number weight percent
92 : ! for hygroscopic species of externally mixed aerosols
93 : real(r8), pointer :: sw_hygro_coreshell_ext(:,:,:,:,:)
94 : real(r8), pointer :: sw_hygro_coreshell_ssa(:,:,:,:,:)
95 : real(r8), pointer :: sw_hygro_coreshell_asm(:,:,:,:,:)
96 : real(r8), pointer :: lw_hygro_coreshell_abs(:,:,:,:,:)
97 : real(r8), pointer :: bcdust(:) ! table of bc-dust mass ratio
98 : real(r8), pointer :: kap(:) ! table of kappa
99 : real(r8), pointer :: relh(:) ! table of relative humidity
100 : integer :: nbcdust
101 : integer :: nkap
102 : integer :: nrelh
103 :
104 : ! microphysics parameters.
105 : character(len=32) :: aername ! for output of number concentration
106 : real(r8) :: density_aer ! density of aerosol (kg/m3)
107 : real(r8) :: hygro_aer ! hygroscopicity of aerosol
108 : real(r8) :: dryrad_aer ! number mode radius (m) of aerosol size distribution
109 : real(r8) :: dispersion_aer ! geometric standard deviation of aerosol size distribution
110 : real(r8) :: num_to_mass_aer ! ratio of number concentration to mass concentration (#/kg)
111 : ! *** Is this actually (kg/#) ???
112 : ! mode parameters
113 : integer :: ncoef ! number of Chebyshev coefficients
114 : integer :: prefr ! dimension in table of real refractive indices
115 : integer :: prefi ! dimension in table of imag refractive indices
116 : real(r8) :: sigmag ! geometric standard deviation of the number distribution for aerosol mode
117 : real(r8) :: dgnum ! geometric dry mean diameter of the number distribution for aerosol mode
118 : real(r8) :: dgnumlo ! lower limit of dgnum
119 : real(r8) :: dgnumhi ! upper limit of dgnum
120 : real(r8) :: rhcrystal ! crystalization relative humidity for mode
121 : real(r8) :: rhdeliques ! deliquescence relative humidity for mode
122 :
123 : endtype physprop_type
124 :
125 : ! This module stores data in an array of physprop_type structures. The way this data
126 : ! is accessed outside the module is via a physprop ID, which is an index into the array.
127 : integer :: numphysprops = 0 ! an incremental total across ALL clim and diag constituents
128 : type (physprop_type), pointer :: physprop(:)
129 :
130 : ! Temporary storage location for filenames in namelist, and construction of dynamic index
131 : ! to properties. The unique filenames specified in the namelist are the identifiers of
132 : ! the properties. Searching the uniquefilenames array provides the index into the physprop
133 : ! array.
134 : character(len=256), allocatable :: uniquefilenames(:)
135 :
136 : ! Number of evenly spaced intervals in rh used in this module and in the aer_rad_props module
137 : ! for calculations of aerosol hygroscopic growth.
138 : integer, parameter, public :: nrh = 1000
139 :
140 : !================================================================================================
141 : contains
142 : !================================================================================================
143 :
144 7680 : subroutine physprop_accum_unique_files(radname, type)
145 :
146 : ! Count number of aerosols in input radname array. Aerosols are identified
147 : ! as strings with a ".nc" suffix.
148 : ! Construct a cumulative list of unique filenames containing physical property data.
149 :
150 : character(len=*), intent(in) :: radname(:)
151 : character(len=1), intent(in) :: type(:)
152 :
153 : integer :: ncnst, i
154 : character(len=*), parameter :: subname = 'physprop_accum_unique_files'
155 : !------------------------------------------------------------------------------------
156 :
157 : ! Initial guess for number of files we need.
158 7680 : if (.not. allocated(uniquefilenames)) allocate(uniquefilenames(50))
159 :
160 7680 : ncnst = ubound(radname, 1)
161 :
162 53760 : do i = 1, ncnst
163 :
164 : ! check if radname is either a bulk aerosol or a mode
165 53760 : if (type(i) == 'A' .or. type(i) == 'M' .or. type(i) == 'B') then
166 :
167 : ! check if this filename has been used by another aerosol. If not
168 : ! then add it to the list of unique names.
169 33792 : if (physprop_get_id(radname(i)) < 0) then
170 19968 : numphysprops = numphysprops + 1
171 19968 : if (numphysprops > size(uniquefilenames)) then
172 0 : call double_capacity(uniquefilenames)
173 : end if
174 19968 : uniquefilenames(numphysprops) = trim(radname(i))
175 : endif
176 :
177 : endif
178 : enddo
179 :
180 : contains
181 :
182 : ! Simple routine to re-allocate an array with twice the size, but with
183 : ! the inital values being preserved.
184 0 : subroutine double_capacity(array)
185 : character(len=256), intent(inout), allocatable :: array(:)
186 :
187 0 : character(len=256), allocatable :: tmp(:)
188 : integer :: ierr
189 :
190 0 : allocate(tmp(size(array)*2), stat=ierr)
191 0 : if ( ierr /= 0 ) then
192 0 : call endrun('physprop_accum_unique_files: Allocation error.')
193 : end if
194 :
195 0 : tmp(:size(array)) = array
196 :
197 0 : deallocate(array, stat=ierr)
198 : if ( ierr /= 0 ) then
199 0 : call endrun('physprop_accum_unique_files: Deallocation error.')
200 : end if
201 :
202 0 : call move_alloc(tmp, array)
203 :
204 0 : end subroutine double_capacity
205 :
206 : end subroutine physprop_accum_unique_files
207 :
208 : !================================================================================================
209 :
210 1536 : subroutine physprop_init()
211 :
212 : ! Read properties from the aerosol data files.
213 :
214 : ! ***N.B.*** The calls to physprop_accum_unique_files must be made before calling
215 : ! this init routine. physprop_accum_unique_files is responsible for building
216 : ! the list of files to be read here.
217 :
218 : ! Local variables
219 : integer :: fileindex
220 : type(file_desc_t) :: nc_id ! index to netcdf file
221 : character(len=256) :: locfn ! path to actual file used
222 : character(len=32) :: aername_str ! string read from netCDF file -- may contain trailing
223 : ! nulls which aren't dealt with by trim()
224 :
225 : integer :: ierr ! error codes from mpi
226 :
227 : !------------------------------------------------------------------------------------
228 :
229 4608 : allocate(physprop(numphysprops))
230 :
231 21504 : do fileindex = 1, numphysprops
232 19968 : nullify(physprop(fileindex)%sw_hygro_ext)
233 19968 : nullify(physprop(fileindex)%sw_hygro_ssa)
234 19968 : nullify(physprop(fileindex)%sw_hygro_asm)
235 19968 : nullify(physprop(fileindex)%lw_hygro_abs)
236 :
237 19968 : nullify(physprop(fileindex)%sw_hygro_ext_wtp)
238 19968 : nullify(physprop(fileindex)%sw_hygro_ssa_wtp)
239 19968 : nullify(physprop(fileindex)%sw_hygro_asm_wtp)
240 19968 : nullify(physprop(fileindex)%lw_hygro_abs_wtp)
241 19968 : nullify(physprop(fileindex)%wgtpct)
242 :
243 19968 : nullify(physprop(fileindex)%sw_hygro_coreshell_ext)
244 19968 : nullify(physprop(fileindex)%sw_hygro_coreshell_ssa)
245 19968 : nullify(physprop(fileindex)%sw_hygro_coreshell_asm)
246 19968 : nullify(physprop(fileindex)%lw_hygro_coreshell_abs)
247 19968 : nullify(physprop(fileindex)%bcdust)
248 19968 : nullify(physprop(fileindex)%kap)
249 19968 : nullify(physprop(fileindex)%relh)
250 :
251 19968 : nullify(physprop(fileindex)%sw_nonhygro_ext)
252 19968 : nullify(physprop(fileindex)%sw_nonhygro_ssa)
253 19968 : nullify(physprop(fileindex)%sw_nonhygro_asm)
254 19968 : nullify(physprop(fileindex)%sw_nonhygro_scat)
255 19968 : nullify(physprop(fileindex)%sw_nonhygro_ascat)
256 19968 : nullify(physprop(fileindex)%lw_abs)
257 :
258 19968 : nullify(physprop(fileindex)%refindex_aer_sw)
259 19968 : nullify(physprop(fileindex)%refindex_aer_lw)
260 :
261 19968 : nullify(physprop(fileindex)%r_sw_ext)
262 19968 : nullify(physprop(fileindex)%r_sw_scat)
263 19968 : nullify(physprop(fileindex)%r_sw_ascat)
264 19968 : nullify(physprop(fileindex)%r_lw_abs)
265 19968 : nullify(physprop(fileindex)%mu)
266 :
267 19968 : nullify(physprop(fileindex)%extpsw)
268 19968 : nullify(physprop(fileindex)%abspsw)
269 19968 : nullify(physprop(fileindex)%asmpsw)
270 19968 : nullify(physprop(fileindex)%absplw)
271 19968 : nullify(physprop(fileindex)%refrtabsw)
272 19968 : nullify(physprop(fileindex)%refitabsw)
273 19968 : nullify(physprop(fileindex)%refrtablw)
274 19968 : nullify(physprop(fileindex)%refitablw)
275 :
276 19968 : nullify(physprop(fileindex)%extpsw2)
277 19968 : nullify(physprop(fileindex)%abspsw2)
278 19968 : nullify(physprop(fileindex)%asmpsw2)
279 19968 : nullify(physprop(fileindex)%absplw2)
280 19968 : nullify(physprop(fileindex)%corefrac)
281 :
282 19968 : call getfil(uniquefilenames(fileindex), locfn, 0)
283 19968 : physprop(fileindex)%sourcefile = locfn
284 :
285 : ! Open the physprop file
286 19968 : call cam_pio_openfile(nc_id, locfn, PIO_NOWRITE)
287 :
288 19968 : call aerosol_optics_init(physprop(fileindex), nc_id)
289 :
290 : ! Close the physprop file
291 21504 : call pio_closefile(nc_id)
292 :
293 : end do
294 1536 : end subroutine physprop_init
295 :
296 : !================================================================================================
297 :
298 67584 : integer function physprop_get_id(filename)
299 :
300 : ! Look for filename in the global list of unique filenames (module data uniquefilenames).
301 : ! If found, return it's index in the list. Otherwise return -1.
302 :
303 : character(len=*), intent(in) :: filename
304 : integer iphysprop
305 :
306 67584 : physprop_get_id = -1
307 574464 : do iphysprop = 1, numphysprops
308 574464 : if(trim(uniquefilenames(iphysprop)) == trim(filename) ) then
309 47616 : physprop_get_id = iphysprop
310 47616 : return
311 : endif
312 : enddo
313 :
314 : end function physprop_get_id
315 :
316 : !================================================================================================
317 :
318 0 : subroutine physprop_get(id, sourcefile, opticstype, &
319 : sw_hygro_ext, sw_hygro_ssa, sw_hygro_asm, lw_hygro_abs, &
320 : sw_hygro_ext_wtp, sw_hygro_ssa_wtp, sw_hygro_asm_wtp, lw_hygro_abs_wtp, &
321 : sw_nonhygro_ext, sw_nonhygro_ssa, sw_nonhygro_asm, &
322 : sw_nonhygro_scat, sw_nonhygro_ascat, lw_abs, &
323 : refindex_aer_sw, refindex_aer_lw, &
324 : r_sw_ext, r_sw_scat, r_sw_ascat, r_lw_abs, mu, &
325 : extpsw, abspsw, asmpsw, absplw, refrtabsw, &
326 : refitabsw, refrtablw, refitablw, &
327 0 : aername, density_aer, hygro_aer, dryrad_aer, dispersion_aer, &
328 : num_to_mass_aer, ncoef, prefr, prefi, sigmag, &
329 : dgnum, dgnumlo, dgnumhi, rhcrystal, rhdeliques, &
330 : extpsw2, abspsw2, asmpsw2, absplw2, corefrac, nfrac, &
331 : wgtpct, bcdust, kap, relh, &
332 : nkap, nwtp, nbcdust, nrelh, &
333 : sw_hygro_coreshell_ext, sw_hygro_coreshell_ssa, &
334 : sw_hygro_coreshell_asm, lw_hygro_coreshell_abs )
335 :
336 : ! Return requested properties for specified ID.
337 :
338 : ! Arguments
339 : integer, intent(in) :: id
340 : character(len=256), optional, intent(out) :: sourcefile ! Absolute pathname of data file.
341 : character(len=ot_length), optional, intent(out) :: opticstype
342 : real(r8), optional, pointer :: sw_hygro_ext(:,:)
343 : real(r8), optional, pointer :: sw_hygro_ssa(:,:)
344 : real(r8), optional, pointer :: sw_hygro_asm(:,:)
345 : real(r8), optional, pointer :: lw_hygro_abs(:,:)
346 : real(r8), optional, pointer :: sw_hygro_ext_wtp(:,:)
347 : real(r8), optional, pointer :: sw_hygro_ssa_wtp(:,:)
348 : real(r8), optional, pointer :: sw_hygro_asm_wtp(:,:)
349 : real(r8), optional, pointer :: lw_hygro_abs_wtp(:,:)
350 : real(r8), optional, pointer :: wgtpct(:)
351 : integer, optional, intent(out) :: nwtp
352 :
353 : real(r8), optional, pointer :: sw_hygro_coreshell_ext(:,:,:,:,:)
354 : real(r8), optional, pointer :: sw_hygro_coreshell_ssa(:,:,:,:,:)
355 : real(r8), optional, pointer :: sw_hygro_coreshell_asm(:,:,:,:,:)
356 : real(r8), optional, pointer :: lw_hygro_coreshell_abs(:,:,:,:,:)
357 : real(r8), optional, pointer :: kap(:)
358 : integer, optional, intent(out) :: nkap
359 : real(r8), optional, pointer :: bcdust(:)
360 : integer, optional, intent(out) :: nbcdust
361 : real(r8), optional, pointer :: relh(:)
362 : integer, optional, intent(out) :: nrelh
363 :
364 : real(r8), optional, pointer :: sw_nonhygro_ext(:)
365 : real(r8), optional, pointer :: sw_nonhygro_ssa(:)
366 : real(r8), optional, pointer :: sw_nonhygro_asm(:)
367 : real(r8), optional, pointer :: sw_nonhygro_scat(:)
368 : real(r8), optional, pointer :: sw_nonhygro_ascat(:)
369 : real(r8), optional, pointer :: lw_abs(:)
370 : complex(r8), optional, pointer :: refindex_aer_sw(:)
371 : complex(r8), optional, pointer :: refindex_aer_lw(:)
372 : real(r8), optional, pointer :: r_sw_ext(:,:)
373 : real(r8), optional, pointer :: r_sw_scat(:,:)
374 : real(r8), optional, pointer :: r_sw_ascat(:,:)
375 : real(r8), optional, pointer :: r_lw_abs(:,:)
376 : real(r8), optional, pointer :: mu(:)
377 : real(r8), optional, pointer :: extpsw(:,:,:,:)
378 : real(r8), optional, pointer :: abspsw(:,:,:,:)
379 : real(r8), optional, pointer :: asmpsw(:,:,:,:)
380 : real(r8), optional, pointer :: absplw(:,:,:,:)
381 : real(r8), optional, pointer :: refrtabsw(:,:)
382 : real(r8), optional, pointer :: refitabsw(:,:)
383 : real(r8), optional, pointer :: refrtablw(:,:)
384 : real(r8), optional, pointer :: refitablw(:,:)
385 : character(len=20), optional, intent(out) :: aername
386 : real(r8), optional, intent(out) :: density_aer
387 : real(r8), optional, intent(out) :: hygro_aer
388 : real(r8), optional, intent(out) :: dryrad_aer
389 : real(r8), optional, intent(out) :: dispersion_aer
390 : real(r8), optional, intent(out) :: num_to_mass_aer
391 : integer, optional, intent(out) :: ncoef
392 : integer, optional, intent(out) :: prefr
393 : integer, optional, intent(out) :: prefi
394 : real(r8), optional, intent(out) :: sigmag
395 : real(r8), optional, intent(out) :: dgnum
396 : real(r8), optional, intent(out) :: dgnumlo
397 : real(r8), optional, intent(out) :: dgnumhi
398 : real(r8), optional, intent(out) :: rhcrystal
399 : real(r8), optional, intent(out) :: rhdeliques
400 : ! for core/shell
401 : real(r8), optional, pointer :: extpsw2(:,:)
402 : real(r8), optional, pointer :: abspsw2(:,:)
403 : real(r8), optional, pointer :: asmpsw2(:,:)
404 : real(r8), optional, pointer :: absplw2(:,:)
405 : real(r8), optional, pointer :: corefrac(:)
406 : integer, optional, intent(out) :: nfrac
407 :
408 : ! Local variables
409 : character(len=*), parameter :: subname = 'physprop_get'
410 : !------------------------------------------------------------------------------------
411 :
412 5696288967 : if (id <= 0 .or. id > numphysprops) then
413 0 : write(iulog,*) subname//': illegal ID value: ', id
414 0 : call endrun('physprop_get: ID out of range')
415 : end if
416 :
417 5696288967 : if (present(sourcefile)) sourcefile = physprop(id)%sourcefile
418 5696288967 : if (present(opticstype)) opticstype = physprop(id)%opticsmethod
419 5696288967 : if (present(sw_hygro_ext)) sw_hygro_ext => physprop(id)%sw_hygro_ext
420 5696288967 : if (present(sw_hygro_ssa)) sw_hygro_ssa => physprop(id)%sw_hygro_ssa
421 5696288967 : if (present(sw_hygro_asm)) sw_hygro_asm => physprop(id)%sw_hygro_asm
422 5696288967 : if (present(lw_hygro_abs)) lw_hygro_abs => physprop(id)%lw_hygro_abs
423 5696288967 : if (present(sw_hygro_ext_wtp)) sw_hygro_ext_wtp => physprop(id)%sw_hygro_ext_wtp
424 5696288967 : if (present(sw_hygro_ssa_wtp)) sw_hygro_ssa_wtp => physprop(id)%sw_hygro_ssa_wtp
425 5696288967 : if (present(sw_hygro_asm_wtp)) sw_hygro_asm_wtp => physprop(id)%sw_hygro_asm_wtp
426 5696288967 : if (present(lw_hygro_abs_wtp)) lw_hygro_abs_wtp => physprop(id)%lw_hygro_abs_wtp
427 5696288967 : if (present(wgtpct)) wgtpct => physprop(id)%wgtpct
428 5696288967 : if (present(nwtp)) nwtp = physprop(id)%nwtp
429 5696288967 : if (present(sw_hygro_coreshell_ext)) sw_hygro_coreshell_ext => physprop(id)%sw_hygro_coreshell_ext
430 5696288967 : if (present(sw_hygro_coreshell_ssa)) sw_hygro_coreshell_ssa => physprop(id)%sw_hygro_coreshell_ssa
431 5696288967 : if (present(sw_hygro_coreshell_asm)) sw_hygro_coreshell_asm => physprop(id)%sw_hygro_coreshell_asm
432 5696288967 : if (present(lw_hygro_coreshell_abs)) lw_hygro_coreshell_abs => physprop(id)%lw_hygro_coreshell_abs
433 5696288967 : if (present(kap)) kap => physprop(id)%kap
434 5696288967 : if (present(nkap)) nkap = physprop(id)%nkap
435 5696288967 : if (present(bcdust)) bcdust => physprop(id)%bcdust
436 5696288967 : if (present(nbcdust)) nbcdust = physprop(id)%nbcdust
437 5696288967 : if (present(relh)) relh => physprop(id)%relh
438 5696288967 : if (present(nrelh)) nrelh = physprop(id)%nrelh
439 5696288967 : if (present(sw_nonhygro_ext)) sw_nonhygro_ext => physprop(id)%sw_nonhygro_ext
440 5696288967 : if (present(sw_nonhygro_ssa)) sw_nonhygro_ssa => physprop(id)%sw_nonhygro_ssa
441 5696288967 : if (present(sw_nonhygro_asm)) sw_nonhygro_asm => physprop(id)%sw_nonhygro_asm
442 5696288967 : if (present(sw_nonhygro_scat)) sw_nonhygro_scat => physprop(id)%sw_nonhygro_scat
443 5696288967 : if (present(sw_nonhygro_ascat)) sw_nonhygro_ascat => physprop(id)%sw_nonhygro_ascat
444 5696288967 : if (present(lw_abs)) lw_abs => physprop(id)%lw_abs
445 :
446 5696288967 : if (present(refindex_aer_sw)) refindex_aer_sw => physprop(id)%refindex_aer_sw
447 5696288967 : if (present(refindex_aer_lw)) refindex_aer_lw => physprop(id)%refindex_aer_lw
448 :
449 5696288967 : if (present(r_sw_ext)) r_sw_ext => physprop(id)%r_sw_ext
450 5696288967 : if (present(r_sw_scat)) r_sw_scat => physprop(id)%r_sw_scat
451 5696288967 : if (present(r_sw_ascat)) r_sw_ascat => physprop(id)%r_sw_ascat
452 5696288967 : if (present(r_lw_abs)) r_lw_abs => physprop(id)%r_lw_abs
453 5696288967 : if (present(mu)) mu => physprop(id)%mu
454 :
455 5696288967 : if (present(extpsw)) extpsw => physprop(id)%extpsw
456 5696288967 : if (present(abspsw)) abspsw => physprop(id)%abspsw
457 5696288967 : if (present(asmpsw)) asmpsw => physprop(id)%asmpsw
458 5696288967 : if (present(absplw)) absplw => physprop(id)%absplw
459 5696288967 : if (present(refrtabsw)) refrtabsw => physprop(id)%refrtabsw
460 5696288967 : if (present(refitabsw)) refitabsw => physprop(id)%refitabsw
461 5696288967 : if (present(refrtablw)) refrtablw => physprop(id)%refrtablw
462 5696288967 : if (present(refitablw)) refitablw => physprop(id)%refitablw
463 :
464 5696288967 : if (present(aername)) aername = physprop(id)%aername
465 5696288967 : if (present(density_aer)) density_aer = physprop(id)%density_aer
466 5696288967 : if (present(hygro_aer)) hygro_aer = physprop(id)%hygro_aer
467 5696288967 : if (present(dryrad_aer)) dryrad_aer = physprop(id)%dryrad_aer
468 5696288967 : if (present(dispersion_aer)) dispersion_aer = physprop(id)%dispersion_aer
469 5696288967 : if (present(num_to_mass_aer)) num_to_mass_aer = physprop(id)%num_to_mass_aer
470 :
471 5696288967 : if (present(ncoef)) ncoef = physprop(id)%ncoef
472 5696288967 : if (present(prefr)) prefr = physprop(id)%prefr
473 5696288967 : if (present(prefi)) prefi = physprop(id)%prefi
474 5696288967 : if (present(sigmag)) sigmag = physprop(id)%sigmag
475 5696288967 : if (present(dgnum)) dgnum = physprop(id)%dgnum
476 5696288967 : if (present(dgnumlo)) dgnumlo = physprop(id)%dgnumlo
477 5696288967 : if (present(dgnumhi)) dgnumhi = physprop(id)%dgnumhi
478 5696288967 : if (present(rhcrystal)) rhcrystal = physprop(id)%rhcrystal
479 5696288967 : if (present(rhdeliques)) rhdeliques = physprop(id)%rhdeliques
480 :
481 : ! For core/shell bins
482 5696288967 : if (present(extpsw2)) extpsw2 => physprop(id)%extpsw2
483 5696288967 : if (present(abspsw2)) abspsw2 => physprop(id)%abspsw2
484 5696288967 : if (present(asmpsw2)) asmpsw2 => physprop(id)%asmpsw2
485 5696288967 : if (present(absplw2)) absplw2 => physprop(id)%absplw2
486 5696288967 : if (present(corefrac)) corefrac => physprop(id)%corefrac
487 5696288967 : if (present(nfrac)) nfrac = physprop(id)%nfrac
488 :
489 5696288967 : end subroutine physprop_get
490 :
491 : !================================================================================================
492 : ! Private methods
493 : !================================================================================================
494 :
495 19968 : subroutine aerosol_optics_init(phys_prop, nc_id)
496 :
497 : ! Determine the opticstype, then call the
498 : ! appropriate routine to read the data.
499 :
500 : type(physprop_type), intent(inout) :: phys_prop ! data after interp onto cam rh mesh
501 : type(file_desc_t), intent(inout) :: nc_id ! indentifier for netcdf file
502 :
503 : integer :: opticslength_id, opticslength
504 : type(var_desc_t) :: op_type_id
505 : integer :: ierr ! mpi error codes
506 : character(len=ot_length) :: opticstype_str ! string read from netCDF file -- may contain trailing
507 : ! nulls which aren't dealt with by trim()
508 : !------------------------------------------------------------------------------------
509 :
510 19968 : ierr = pio_inq_dimid(nc_id, 'opticsmethod_len', opticslength_id)
511 19968 : ierr = pio_inq_dimlen(nc_id, opticslength_id, opticslength)
512 19968 : if ( opticslength .gt. ot_length ) then
513 0 : call endrun(" optics type length in "//phys_prop%sourcefile//" excedes maximum length of 32")
514 : endif
515 19968 : ierr = pio_inq_varid(nc_id, 'opticsmethod', op_type_id)
516 19968 : ierr = pio_get_var(nc_id, op_type_id,phys_prop%opticsmethod )
517 :
518 0 : select case (phys_prop%opticsmethod)
519 : case ('zero')
520 0 : call zero_optics_init(phys_prop, nc_id)
521 :
522 : case ('hygro')
523 0 : call hygro_optics_init(phys_prop, nc_id)
524 :
525 : case ('hygroscopic')
526 4608 : call hygroscopic_optics_init(phys_prop, nc_id)
527 :
528 : case ('hygroscopic_wtp')
529 0 : call hygroscopic_wtp_optics_init(phys_prop, nc_id)
530 :
531 : case ('hygroscopic_coreshell')
532 0 : call hygroscopic_coreshell_optics_init(phys_prop, nc_id)
533 :
534 : case ('nonhygro')
535 0 : call nonhygro_optics_init(phys_prop, nc_id)
536 :
537 : case ('insoluble')
538 4608 : call insoluble_optics_init(phys_prop, nc_id)
539 :
540 : case ('volcanic_radius','volcanic_radius1','volcanic_radius2','volcanic_radius3')
541 4608 : call volcanic_radius_optics_init(phys_prop, nc_id)
542 :
543 : case ('volcanic')
544 0 : call volcanic_optics_init(phys_prop, nc_id)
545 :
546 : case ('modal')
547 6144 : call modal_optics_init(phys_prop, nc_id)
548 :
549 : case ('sectional')
550 0 : call bin_optics_init(phys_prop, nc_id)
551 :
552 : case ('sectional_props')
553 0 : call bindef_optics_init(phys_prop, nc_id)
554 :
555 : ! other types of optics can be added here
556 :
557 : case default
558 : call endrun('aerosol_optics_init: unsupported optics type '//&
559 19968 : trim(phys_prop%opticsmethod)//' in file '//phys_prop%sourcefile)
560 : end select
561 :
562 5696288967 : end subroutine aerosol_optics_init
563 :
564 : !================================================================================================
565 :
566 0 : subroutine hygro_optics_init(phys_prop, nc_id)
567 :
568 : ! Read optics data of type 'hygro' and interpolate it to CAM's rh mesh.
569 :
570 : type (physprop_type), intent(inout) :: phys_prop ! data after interp onto cam rh mesh
571 : type (file_desc_t), intent(inout) :: nc_id ! indentifier for netcdf file
572 :
573 : ! Local variables
574 : integer :: ierr ! error flag
575 :
576 : integer :: rh_idx_id, lw_band_id, sw_band_id
577 : integer :: kbnd, krh
578 : integer :: rh_id, sw_ext_id, sw_ssa_id, sw_asm_id, lw_ext_id
579 : integer :: nbnd, swbands
580 :
581 : ! temp data from hygroscopic file before interpolation onto cam-rh-mesh
582 : integer :: nfilerh ! number of rh values in file
583 0 : real(r8), allocatable, dimension(:) :: frh
584 0 : real(r8), allocatable, dimension(:,:) :: fsw_ext
585 0 : real(r8), allocatable, dimension(:,:) :: fsw_ssa
586 0 : real(r8), allocatable, dimension(:,:) :: fsw_asm
587 :
588 : real(r8) :: rh ! real rh value on cam rh mesh (indexvalue)
589 : !------------------------------------------------------------------------------------
590 :
591 0 : allocate(phys_prop%sw_hygro_ext(nrh,nswbands))
592 0 : allocate(phys_prop%sw_hygro_ssa(nrh,nswbands))
593 0 : allocate(phys_prop%sw_hygro_asm(nrh,nswbands))
594 0 : allocate(phys_prop%lw_abs(nlwbands))
595 :
596 0 : ierr = pio_inq_dimid(nc_id, 'rh_idx', rh_idx_id)
597 :
598 0 : ierr = pio_inq_dimlen(nc_id, rh_idx_id, nfilerh)
599 :
600 0 : ierr = pio_inq_dimid(nc_id, 'lw_band', lw_band_id)
601 :
602 0 : ierr = pio_inq_dimid(nc_id, 'sw_band', sw_band_id)
603 :
604 0 : ierr = pio_inq_dimlen(nc_id, lw_band_id, nbnd)
605 :
606 0 : if (nbnd .ne. nlwbands) call endrun(phys_prop%sourcefile// &
607 0 : ' has the wrong number of lwbands')
608 :
609 0 : ierr = pio_inq_dimlen(nc_id, sw_band_id, swbands)
610 :
611 0 : if(swbands .ne. nswbands) call endrun(phys_prop%sourcefile// &
612 0 : ' has the wrong number of sw bands')
613 :
614 0 : ierr = pio_inq_varid(nc_id, 'rh', rh_id)
615 :
616 0 : ierr = pio_inq_varid(nc_id, 'ext_sw', sw_ext_id)
617 :
618 0 : ierr = pio_inq_varid(nc_id, 'ssa_sw', sw_ssa_id)
619 :
620 0 : ierr = pio_inq_varid(nc_id, 'asm_sw', sw_asm_id)
621 :
622 0 : ierr = pio_inq_varid(nc_id, 'abs_lw', lw_ext_id)
623 :
624 : ! specific optical properties on file's rh mesh
625 0 : allocate(fsw_ext(nfilerh,nswbands))
626 0 : allocate(fsw_asm(nfilerh,nswbands))
627 0 : allocate(fsw_ssa(nfilerh,nswbands))
628 0 : allocate(frh(nfilerh))
629 :
630 0 : ierr = pio_get_var(nc_id, rh_id, frh)
631 :
632 0 : ierr = pio_get_var(nc_id, sw_ext_id, fsw_ext)
633 :
634 0 : ierr = pio_get_var(nc_id, sw_ssa_id, fsw_ssa)
635 :
636 0 : ierr = pio_get_var(nc_id, sw_asm_id, fsw_asm)
637 :
638 0 : ierr = pio_get_var(nc_id, lw_ext_id, phys_prop%lw_abs)
639 :
640 : ! interpolate onto cam's rh mesh
641 0 : do kbnd = 1,nswbands
642 0 : do krh = 1, nrh
643 0 : rh = 1.0_r8 / nrh * (krh - 1)
644 0 : phys_prop%sw_hygro_ext(krh,kbnd) = &
645 0 : exp_interpol( frh, fsw_ext(:,kbnd) / fsw_ext(1,kbnd), rh ) &
646 0 : * fsw_ext(1, kbnd)
647 0 : phys_prop%sw_hygro_ssa(krh,kbnd) = &
648 0 : lin_interpol( frh, fsw_ssa(:,kbnd) / fsw_ssa(1,kbnd), rh ) &
649 0 : * fsw_ssa(1, kbnd)
650 0 : phys_prop%sw_hygro_asm(krh,kbnd) = &
651 0 : lin_interpol( frh, fsw_asm(:,kbnd) / fsw_asm(1,kbnd), rh ) &
652 0 : * fsw_asm(1, kbnd)
653 : enddo
654 : enddo
655 :
656 0 : deallocate (fsw_ext, fsw_asm, fsw_ssa, frh)
657 :
658 : ! read refractive index data if available
659 0 : call refindex_aer_init(phys_prop, nc_id)
660 :
661 : ! read bulk aero props
662 0 : call bulk_props_init(phys_prop, nc_id)
663 :
664 0 : end subroutine hygro_optics_init
665 :
666 : !================================================================================================
667 :
668 0 : subroutine zero_optics_init(phys_prop, nc_id)
669 :
670 : ! Read optics data of type 'nonhygro'
671 :
672 : type (physprop_type), intent(inout) :: phys_prop ! storage for file data
673 : type (file_desc_t), intent(inout) :: nc_id ! indentifier for netcdf file
674 :
675 : ! Local variables
676 : integer :: lw_band_id, sw_band_id
677 : integer :: sw_ext_id, sw_ssa_id, sw_asm_id, lw_ext_id
678 : integer :: swbands, nbnd
679 : integer :: ierr ! error flag
680 : !------------------------------------------------------------------------------------
681 :
682 : ! perhaps this doesn't even need allocated.
683 0 : allocate (phys_prop%sw_nonhygro_ext(nswbands))
684 0 : allocate (phys_prop%sw_nonhygro_ssa(nswbands))
685 0 : allocate (phys_prop%sw_nonhygro_asm(nswbands))
686 0 : allocate (phys_prop%lw_abs(nlwbands))
687 :
688 0 : phys_prop%sw_nonhygro_ext = 0._r8
689 0 : phys_prop%sw_nonhygro_ssa = 0._r8
690 0 : phys_prop%sw_nonhygro_asm = 0._r8
691 0 : phys_prop%lw_abs = 0._r8
692 :
693 0 : end subroutine zero_optics_init
694 :
695 : !================================================================================================
696 :
697 4608 : subroutine insoluble_optics_init(phys_prop, nc_id)
698 :
699 : ! Read optics data of type 'nonhygro'
700 :
701 : type (physprop_type), intent(inout) :: phys_prop ! storage for file data
702 : type (file_desc_t), intent(inout) :: nc_id ! indentifier for netcdf file
703 :
704 : ! Local variables
705 : integer :: lw_band_id, sw_band_id
706 : integer :: sw_ext_id, sw_ssa_id, sw_asm_id, lw_ext_id
707 : integer :: swbands, nbnd
708 : integer :: ierr ! error flag
709 : integer :: start(2), count(2)
710 : !------------------------------------------------------------------------------------
711 :
712 4608 : allocate (phys_prop%sw_nonhygro_ext(nswbands))
713 4608 : allocate (phys_prop%sw_nonhygro_ssa(nswbands))
714 4608 : allocate (phys_prop%sw_nonhygro_asm(nswbands))
715 4608 : allocate (phys_prop%lw_abs(nlwbands))
716 :
717 4608 : ierr = pio_inq_dimid(nc_id, 'lw_band', lw_band_id)
718 :
719 4608 : ierr = pio_inq_dimid(nc_id, 'sw_band', sw_band_id)
720 :
721 4608 : ierr = pio_inq_dimlen(nc_id, lw_band_id, nbnd)
722 :
723 4608 : if (nbnd .ne. nlwbands) call endrun(phys_prop%sourcefile// &
724 0 : ' has the wrong number of lwbands')
725 :
726 4608 : ierr = pio_inq_dimlen(nc_id, sw_band_id, swbands)
727 :
728 4608 : if (swbands .ne. nswbands) call endrun(phys_prop%sourcefile// &
729 0 : ' has the wrong number of sw bands')
730 :
731 : ! read file data
732 4608 : ierr = pio_inq_varid(nc_id, 'ext_sw', sw_ext_id)
733 4608 : ierr = pio_inq_varid(nc_id, 'ssa_sw', sw_ssa_id)
734 4608 : ierr = pio_inq_varid(nc_id, 'asm_sw', sw_asm_id)
735 4608 : ierr = pio_inq_varid(nc_id, 'abs_lw', lw_ext_id)
736 :
737 13824 : start = 1
738 13824 : count=(/1,swbands/)
739 :
740 4608 : ierr = pio_get_var(nc_id, sw_ext_id, start, count, phys_prop%sw_nonhygro_ext)
741 4608 : ierr = pio_get_var(nc_id, sw_ssa_id, start, count, phys_prop%sw_nonhygro_ssa)
742 4608 : ierr = pio_get_var(nc_id, sw_asm_id, start, count, phys_prop%sw_nonhygro_asm)
743 13824 : count = (/1,nbnd/)
744 4608 : ierr = pio_get_var(nc_id, lw_ext_id, start, count, phys_prop%lw_abs)
745 :
746 : ! read refractive index data if available
747 4608 : call refindex_aer_init(phys_prop, nc_id)
748 :
749 : ! read bulk aero props
750 4608 : call bulk_props_init(phys_prop, nc_id)
751 :
752 4608 : end subroutine insoluble_optics_init
753 :
754 : !================================================================================================
755 :
756 4608 : subroutine volcanic_radius_optics_init(phys_prop, nc_id)
757 :
758 : ! Read optics data of type 'volcanic_radius'
759 :
760 : type (physprop_type), intent(inout) :: phys_prop ! storage for file data
761 : type (file_desc_t), intent(inout) :: nc_id ! indentifier for netcdf file
762 :
763 : ! Local variables
764 : integer :: lw_band_id, sw_band_id, mu_id, mu_did
765 : integer :: sw_ext_id, sw_scat_id, sw_ascat_id, lw_abs_id
766 : integer :: swbands, nbnd, n_mu_samples
767 : integer :: ierr ! error flag
768 : !------------------------------------------------------------------------------------
769 :
770 4608 : ierr = pio_inq_dimid(nc_id, 'mu_samples', mu_did)
771 4608 : ierr = pio_inq_dimlen(nc_id, mu_did, n_mu_samples)
772 :
773 13824 : allocate (phys_prop%r_sw_ext(nswbands,n_mu_samples))
774 9216 : allocate (phys_prop%r_sw_scat(nswbands,n_mu_samples))
775 9216 : allocate (phys_prop%r_sw_ascat(nswbands,n_mu_samples))
776 13824 : allocate (phys_prop%r_lw_abs(nlwbands,n_mu_samples))
777 13824 : allocate (phys_prop%mu(n_mu_samples))
778 :
779 4608 : ierr = pio_inq_dimid(nc_id, 'lw_band', lw_band_id)
780 :
781 4608 : ierr = pio_inq_dimid(nc_id, 'sw_band', sw_band_id)
782 :
783 4608 : ierr = pio_inq_dimlen(nc_id, lw_band_id, nbnd)
784 :
785 4608 : if (nbnd .ne. nlwbands) call endrun(phys_prop%sourcefile// &
786 0 : ' has the wrong number of lwbands')
787 :
788 4608 : ierr = pio_inq_dimlen(nc_id, sw_band_id, swbands)
789 :
790 4608 : if (swbands .ne. nswbands) call endrun(phys_prop%sourcefile// &
791 0 : ' has the wrong number of sw bands')
792 :
793 : ! read file data
794 4608 : ierr = pio_inq_varid(nc_id, 'bext_sw', sw_ext_id)
795 4608 : ierr = pio_inq_varid(nc_id, 'bsca_sw', sw_scat_id)
796 4608 : ierr = pio_inq_varid(nc_id, 'basc_sw', sw_ascat_id)
797 4608 : ierr = pio_inq_varid(nc_id, 'babs_lw', lw_abs_id)
798 4608 : ierr = pio_inq_varid(nc_id, 'mu_samples', mu_id)
799 :
800 4608 : ierr = pio_get_var(nc_id, sw_ext_id, phys_prop%r_sw_ext)
801 4608 : ierr = pio_get_var(nc_id, sw_scat_id, phys_prop%r_sw_scat)
802 4608 : ierr = pio_get_var(nc_id, sw_ascat_id, phys_prop%r_sw_ascat)
803 4608 : ierr = pio_get_var(nc_id, lw_abs_id, phys_prop%r_lw_abs)
804 4608 : ierr = pio_get_var(nc_id, mu_id, phys_prop%mu)
805 :
806 : ! read bulk aero props
807 4608 : call bulk_props_init(phys_prop, nc_id)
808 :
809 4608 : end subroutine volcanic_radius_optics_init
810 :
811 : !================================================================================================
812 :
813 0 : subroutine volcanic_optics_init(phys_prop, nc_id)
814 :
815 : ! Read optics data of type 'volcanic'
816 :
817 : type (physprop_type), intent(inout) :: phys_prop ! storage for file data
818 : type (file_desc_t) , intent(inout) :: nc_id ! indentifier for netcdf file
819 :
820 : ! Local variables
821 : integer :: lw_band_id, sw_band_id
822 : integer :: sw_ext_id, sw_scat_id, sw_ascat_id, lw_abs_id
823 : integer :: swbands, nbnd
824 : integer :: ierr ! error flag
825 : !------------------------------------------------------------------------------------
826 :
827 0 : allocate (phys_prop%sw_nonhygro_ext(nswbands))
828 0 : allocate (phys_prop%sw_nonhygro_scat(nswbands))
829 0 : allocate (phys_prop%sw_nonhygro_ascat(nswbands))
830 0 : allocate (phys_prop%lw_abs(nlwbands))
831 :
832 0 : ierr = pio_inq_dimid(nc_id, 'lw_band', lw_band_id)
833 0 : ierr = pio_inq_dimid(nc_id, 'sw_band', sw_band_id)
834 :
835 0 : ierr = pio_inq_dimlen(nc_id, lw_band_id, nbnd)
836 :
837 0 : if (nbnd .ne. nlwbands) call endrun(phys_prop%sourcefile// &
838 0 : ' has the wrong number of lwbands')
839 :
840 0 : ierr = pio_inq_dimlen(nc_id, sw_band_id, swbands)
841 0 : if(masterproc) write(iulog,*) 'swbands',swbands
842 :
843 0 : if (swbands .ne. nswbands) call endrun(phys_prop%sourcefile// &
844 0 : ' has the wrong number of sw bands')
845 :
846 : ! read file data
847 0 : ierr = pio_inq_varid(nc_id, 'bext_sw', sw_ext_id)
848 0 : ierr = pio_inq_varid(nc_id, 'bsca_sw', sw_scat_id)
849 0 : ierr = pio_inq_varid(nc_id, 'basc_sw', sw_ascat_id)
850 0 : ierr = pio_inq_varid(nc_id, 'babs_lw', lw_abs_id)
851 :
852 0 : ierr = pio_get_var(nc_id, sw_ext_id, phys_prop%sw_nonhygro_ext)
853 0 : ierr = pio_get_var(nc_id, sw_scat_id, phys_prop%sw_nonhygro_scat)
854 0 : ierr = pio_get_var(nc_id, sw_ascat_id, phys_prop%sw_nonhygro_ascat)
855 0 : ierr = pio_get_var(nc_id, lw_abs_id, phys_prop%lw_abs)
856 :
857 : ! read bulk aero props
858 0 : call bulk_props_init(phys_prop, nc_id)
859 :
860 0 : end subroutine volcanic_optics_init
861 :
862 : !================================================================================================
863 :
864 4608 : subroutine hygroscopic_optics_init(phys_prop, nc_id)
865 :
866 : ! Read optics data of type 'hygroscopic' and interpolate it to CAM's rh mesh.
867 :
868 : type (physprop_type), intent(inout) :: phys_prop ! data after interp onto cam rh mesh
869 : type (file_desc_T), intent(inout) :: nc_id ! indentifier for netcdf file
870 :
871 : ! Local variables
872 : integer :: ierr ! error flag
873 :
874 : integer :: rh_idx_id, lw_band_id, sw_band_id
875 : integer :: kbnd, krh
876 : integer :: rh_id, sw_ext_id, sw_ssa_id, sw_asm_id, lw_ext_id
877 : integer :: nbnd, swbands
878 :
879 : ! temp data from hygroscopic file before interpolation onto cam-rh-mesh
880 : integer :: nfilerh ! number of rh values in file
881 4608 : real(r8), allocatable, dimension(:) :: frh
882 4608 : real(r8), allocatable, dimension(:,:) :: fsw_ext
883 4608 : real(r8), allocatable, dimension(:,:) :: fsw_ssa
884 4608 : real(r8), allocatable, dimension(:,:) :: fsw_asm
885 4608 : real(r8), allocatable, dimension(:,:) :: flw_abs
886 :
887 : real(r8) :: rh ! real rh value on cam rh mesh (indexvalue)
888 : character(len=*), parameter :: sub = 'hygroscopic_optics_init'
889 : !------------------------------------------------------------------------------------
890 :
891 4608 : allocate(phys_prop%sw_hygro_ext(nrh,nswbands))
892 4608 : allocate(phys_prop%sw_hygro_ssa(nrh,nswbands))
893 4608 : allocate(phys_prop%sw_hygro_asm(nrh,nswbands))
894 4608 : allocate(phys_prop%lw_hygro_abs(nrh,nlwbands))
895 :
896 4608 : ierr = pio_inq_dimid(nc_id, 'rh_idx', rh_idx_id)
897 4608 : ierr = pio_inq_dimlen(nc_id, rh_idx_id, nfilerh)
898 :
899 4608 : ierr = pio_inq_dimid(nc_id, 'lw_band', lw_band_id)
900 4608 : ierr = pio_inq_dimlen(nc_id, lw_band_id, nbnd)
901 4608 : if (nbnd .ne. nlwbands) call endrun(phys_prop%sourcefile// &
902 0 : ' has the wrong number of lwbands')
903 :
904 4608 : ierr = pio_inq_dimid(nc_id, 'sw_band', sw_band_id)
905 4608 : ierr = pio_inq_dimlen(nc_id, sw_band_id, swbands)
906 4608 : if(swbands .ne. nswbands) call endrun(phys_prop%sourcefile// &
907 0 : ' has the wrong number of sw bands')
908 :
909 4608 : ierr = pio_inq_varid(nc_id, 'rh', rh_id)
910 4608 : ierr = pio_inq_varid(nc_id, 'ext_sw', sw_ext_id)
911 4608 : ierr = pio_inq_varid(nc_id, 'ssa_sw', sw_ssa_id)
912 4608 : ierr = pio_inq_varid(nc_id, 'asm_sw', sw_asm_id)
913 4608 : ierr = pio_inq_varid(nc_id, 'abs_lw', lw_ext_id)
914 :
915 : ! specific optical properties on file's rh mesh
916 13824 : allocate(fsw_ext(nfilerh,nswbands))
917 9216 : allocate(fsw_asm(nfilerh,nswbands))
918 9216 : allocate(fsw_ssa(nfilerh,nswbands))
919 13824 : allocate(flw_abs(nfilerh,nlwbands))
920 13824 : allocate(frh(nfilerh))
921 :
922 4608 : ierr = pio_get_var(nc_id, rh_id, frh)
923 4608 : ierr = pio_get_var(nc_id, sw_ext_id, fsw_ext)
924 4608 : ierr = pio_get_var(nc_id, sw_ssa_id, fsw_ssa)
925 4608 : ierr = pio_get_var(nc_id, sw_asm_id, fsw_asm)
926 4608 : ierr = pio_get_var(nc_id, lw_ext_id, flw_abs)
927 :
928 : ! interpolate onto cam's rh mesh
929 69120 : do kbnd = 1,nswbands
930 64581120 : do krh = 1, nrh
931 64512000 : rh = 1.0_r8 / nrh * (krh - 1)
932 0 : phys_prop%sw_hygro_ext(krh,kbnd) = &
933 64512000 : exp_interpol( frh, fsw_ext(:,kbnd) / fsw_ext(1,kbnd), rh ) &
934 645120000 : * fsw_ext(1, kbnd)
935 0 : phys_prop%sw_hygro_ssa(krh,kbnd) = &
936 64512000 : lin_interpol( frh, fsw_ssa(:,kbnd) / fsw_ssa(1,kbnd), rh ) &
937 645120000 : * fsw_ssa(1, kbnd)
938 0 : phys_prop%sw_hygro_asm(krh,kbnd) = &
939 64512000 : lin_interpol( frh, fsw_asm(:,kbnd) / fsw_asm(1,kbnd), rh ) &
940 645184512 : * fsw_asm(1, kbnd)
941 : enddo
942 : enddo
943 78336 : do kbnd = 1,nlwbands
944 73806336 : do krh = 1, nrh
945 73728000 : rh = 1.0_r8 / nrh * (krh - 1)
946 0 : phys_prop%lw_hygro_abs(krh,kbnd) = &
947 73728000 : exp_interpol( frh, flw_abs(:,kbnd) / flw_abs(1,kbnd), rh ) &
948 737353728 : * flw_abs(1, kbnd)
949 : enddo
950 : enddo
951 :
952 4608 : deallocate (fsw_ext, fsw_asm, fsw_ssa, flw_abs, frh)
953 :
954 : ! read refractive index data if available
955 4608 : call refindex_aer_init(phys_prop, nc_id)
956 :
957 : ! read bulk aero props
958 4608 : call bulk_props_init(phys_prop, nc_id)
959 :
960 4608 : end subroutine hygroscopic_optics_init
961 :
962 : !================================================================================================
963 :
964 0 : subroutine nonhygro_optics_init(phys_prop, nc_id)
965 :
966 : ! Read optics data of type 'nonhygro'
967 :
968 : type (physprop_type), intent(inout) :: phys_prop ! storage for file data
969 : type (file_desc_t) , intent(inout) :: nc_id ! indentifier for netcdf file
970 :
971 : ! Local variables
972 : integer :: lw_band_id, sw_band_id
973 : integer :: sw_ext_id, sw_ssa_id, sw_asm_id, lw_ext_id
974 : integer :: swbands, nbnd
975 : integer :: ierr ! error flag
976 : !------------------------------------------------------------------------------------
977 :
978 0 : allocate (phys_prop%sw_nonhygro_ext(nswbands))
979 0 : allocate (phys_prop%sw_nonhygro_ssa(nswbands))
980 0 : allocate (phys_prop%sw_nonhygro_asm(nswbands))
981 0 : allocate (phys_prop%lw_abs(nlwbands))
982 :
983 0 : ierr = pio_inq_dimid(nc_id, 'lw_band', lw_band_id)
984 0 : ierr = pio_inq_dimid(nc_id, 'sw_band', sw_band_id)
985 :
986 0 : ierr = pio_inq_dimlen(nc_id, lw_band_id, nbnd)
987 :
988 0 : if (nbnd .ne. nlwbands) call endrun(phys_prop%sourcefile// &
989 0 : ' has the wrong number of lwbands')
990 :
991 0 : ierr = pio_inq_dimlen(nc_id, sw_band_id, swbands)
992 :
993 0 : if (swbands .ne. nswbands) call endrun(phys_prop%sourcefile// &
994 0 : ' has the wrong number of sw bands')
995 :
996 : ! read file data
997 0 : ierr = pio_inq_varid(nc_id, 'ext_sw', sw_ext_id)
998 0 : ierr = pio_inq_varid(nc_id, 'ssa_sw', sw_ssa_id)
999 0 : ierr = pio_inq_varid(nc_id, 'asm_sw', sw_asm_id)
1000 0 : ierr = pio_inq_varid(nc_id, 'abs_lw', lw_ext_id)
1001 :
1002 0 : ierr = pio_get_var(nc_id, sw_ext_id, phys_prop%sw_nonhygro_ext)
1003 0 : ierr = pio_get_var(nc_id, sw_ssa_id, phys_prop%sw_nonhygro_ssa)
1004 0 : ierr = pio_get_var(nc_id, sw_asm_id, phys_prop%sw_nonhygro_asm)
1005 0 : ierr = pio_get_var(nc_id, lw_ext_id, phys_prop%lw_abs)
1006 :
1007 : ! read refractive index data if available
1008 0 : call refindex_aer_init(phys_prop, nc_id)
1009 :
1010 : ! read bulk aero props
1011 0 : call bulk_props_init(phys_prop, nc_id)
1012 :
1013 0 : end subroutine nonhygro_optics_init
1014 :
1015 : !================================================================================================
1016 :
1017 9216 : subroutine refindex_aer_init(phys_prop, nc_id)
1018 :
1019 : ! Read refractive indices of aerosol
1020 :
1021 : type (physprop_type), intent(inout) :: phys_prop ! storage for file data
1022 : type (file_desc_T), intent(inout) :: nc_id ! indentifier for netcdf file
1023 :
1024 : ! Local variables
1025 : integer :: i
1026 : integer :: istat1, istat2, istat3 ! status flags
1027 : integer :: vid_real, vid_im ! variable ids
1028 9216 : real(r8), pointer :: ref_real(:), ref_im(:) ! tmp storage for components of complex index
1029 : character(len=*), parameter :: subname = 'refindex_aer_init'
1030 : !------------------------------------------------------------------------------------
1031 :
1032 : ! assume that the dimensions lw_band and sw_band have already been checked
1033 : ! by the calling subroutine
1034 :
1035 : ! Check that the variables are present before allocating storage and reading.
1036 : ! Since we're setting complex data values, both the real and imaginary parts must
1037 : ! be present or neither will be read.
1038 :
1039 : ! set PIO to return control to the caller when variable not found
1040 9216 : call pio_seterrorhandling(nc_id, PIO_BCAST_ERROR)
1041 :
1042 9216 : istat1 = pio_inq_varid(nc_id, 'refindex_real_aer_sw', vid_real)
1043 9216 : istat2 = pio_inq_varid(nc_id, 'refindex_im_aer_sw', vid_im)
1044 :
1045 9216 : if (istat1 == PIO_NOERR .and. istat2 == PIO_NOERR) then
1046 :
1047 9216 : allocate(ref_real(nswbands), ref_im(nswbands))
1048 :
1049 9216 : istat3 = pio_get_var(nc_id, vid_real, ref_real)
1050 9216 : if (istat3 /= PIO_NOERR) then
1051 0 : call endrun(subname//': ERROR reading refindex_real_aer_sw')
1052 : end if
1053 :
1054 9216 : istat3 = pio_get_var(nc_id, vid_im, ref_im)
1055 9216 : if (istat3 /= PIO_NOERR) then
1056 0 : call endrun(subname//': ERROR reading refindex_im_aer_sw')
1057 : end if
1058 :
1059 : ! successfully read refindex data -- set complex values in physprop object
1060 9216 : allocate(phys_prop%refindex_aer_sw(nswbands))
1061 138240 : do i = 1, nswbands
1062 0 : phys_prop%refindex_aer_sw(i) = cmplx(ref_real(i), abs(ref_im(i)),&
1063 138240 : kind=r8)
1064 : end do
1065 :
1066 9216 : deallocate(ref_real, ref_im)
1067 :
1068 : end if
1069 :
1070 9216 : istat1 = pio_inq_varid(nc_id, 'refindex_real_aer_lw', vid_real)
1071 9216 : istat2 = pio_inq_varid(nc_id, 'refindex_im_aer_lw', vid_im)
1072 :
1073 9216 : if (istat1 == PIO_NOERR .and. istat2 == PIO_NOERR) then
1074 :
1075 9216 : allocate(ref_real(nlwbands), ref_im(nlwbands))
1076 :
1077 9216 : istat3 = pio_get_var(nc_id, vid_real, ref_real)
1078 9216 : if (istat3 /= PIO_NOERR) then
1079 0 : call endrun(subname//': ERROR reading refindex_real_aer_lw')
1080 : end if
1081 :
1082 9216 : istat3 = pio_get_var(nc_id, vid_im, ref_im)
1083 9216 : if (istat3 /= PIO_NOERR) then
1084 0 : call endrun(subname//': ERROR reading refindex_im_aer_lw')
1085 : end if
1086 :
1087 : ! successfully read refindex data -- set complex value in physprop object
1088 9216 : allocate(phys_prop%refindex_aer_lw(nlwbands))
1089 156672 : do i = 1, nlwbands
1090 0 : phys_prop%refindex_aer_lw(i) = cmplx(ref_real(i), abs(ref_im(i)),&
1091 156672 : kind=r8)
1092 : end do
1093 :
1094 9216 : deallocate(ref_real, ref_im)
1095 :
1096 : end if
1097 :
1098 : ! reset PIO to handle errors internally
1099 9216 : call pio_seterrorhandling(nc_id, PIO_INTERNAL_ERROR)
1100 :
1101 9216 : end subroutine refindex_aer_init
1102 :
1103 : !================================================================================================
1104 :
1105 6144 : subroutine modal_optics_init(props, ncid)
1106 :
1107 : ! Read optics data for modal aerosols
1108 :
1109 : type (physprop_type), intent(inout) :: props ! storage for file data
1110 : type (file_desc_T), intent(inout) :: ncid ! indentifier for netcdf file
1111 :
1112 : ! Local variables
1113 : integer :: ierr
1114 : integer :: did
1115 : integer :: ival
1116 : type(var_desc_t) :: vid
1117 6144 : real(r8), pointer :: rval(:,:,:,:,:) ! temp array used to eliminate a singleton dimension
1118 :
1119 : character(len=*), parameter :: subname = 'modal_optics_init'
1120 : !------------------------------------------------------------------------------------
1121 :
1122 : ! Check dimensions for number of lw and sw bands
1123 :
1124 6144 : ierr = pio_inq_dimid(ncid, 'lw_band', did)
1125 6144 : ierr = pio_inq_dimlen(ncid, did, ival)
1126 6144 : if (ival .ne. nlwbands) call endrun(subname//':'//props%sourcefile// &
1127 0 : ' has the wrong number of lw bands')
1128 :
1129 6144 : ierr = pio_inq_dimid(ncid, 'sw_band', did)
1130 6144 : ierr = pio_inq_dimlen(ncid, did, ival)
1131 6144 : if (ival .ne. nswbands) call endrun(subname//':'//props%sourcefile// &
1132 0 : ' has the wrong number of sw bands')
1133 :
1134 : ! Get other dimensions
1135 6144 : ierr = pio_inq_dimid(ncid, 'coef_number', did)
1136 6144 : ierr = pio_inq_dimlen(ncid, did, props%ncoef)
1137 :
1138 6144 : ierr = pio_inq_dimid(ncid, 'refindex_real', did)
1139 6144 : ierr = pio_inq_dimlen(ncid, did, props%prefr)
1140 :
1141 6144 : ierr = pio_inq_dimid(ncid, 'refindex_im', did)
1142 6144 : ierr = pio_inq_dimlen(ncid, did, props%prefi)
1143 :
1144 : ! Allocate arrays
1145 : allocate( &
1146 : props%extpsw(props%ncoef,props%prefr,props%prefi,nswbands), &
1147 : props%abspsw(props%ncoef,props%prefr,props%prefi,nswbands), &
1148 : props%asmpsw(props%ncoef,props%prefr,props%prefi,nswbands), &
1149 : props%absplw(props%ncoef,props%prefr,props%prefi,nlwbands), &
1150 : props%refrtabsw(props%prefr,nswbands), &
1151 : props%refitabsw(props%prefi,nswbands), &
1152 : props%refrtablw(props%prefr,nlwbands), &
1153 153600 : props%refitablw(props%prefi,nlwbands) )
1154 :
1155 :
1156 : ! allocate temp to remove the mode dimension from the sw variables
1157 24576 : allocate(rval(props%ncoef,props%prefr,props%prefi,1,nswbands))
1158 :
1159 6144 : ierr = pio_inq_varid(ncid, 'extpsw', vid)
1160 6144 : ierr = pio_get_var(ncid, vid, rval)
1161 37079040 : props%extpsw = rval(:,:,:,1,:)
1162 :
1163 6144 : ierr = pio_inq_varid(ncid, 'abspsw', vid)
1164 6144 : ierr = pio_get_var(ncid, vid, rval)
1165 37079040 : props%abspsw = rval(:,:,:,1,:)
1166 :
1167 6144 : ierr = pio_inq_varid(ncid, 'asmpsw', vid)
1168 6144 : ierr = pio_get_var(ncid, vid, rval)
1169 37079040 : props%asmpsw = rval(:,:,:,1,:)
1170 :
1171 6144 : deallocate(rval)
1172 :
1173 : ! allocate temp to remove the mode dimension from the lw variables
1174 36864 : allocate(rval(props%ncoef,props%prefr,props%prefi,1,nlwbands))
1175 :
1176 6144 : ierr = pio_inq_varid(ncid, 'absplw', vid)
1177 6144 : ierr = pio_get_var(ncid, vid, rval)
1178 42375168 : props%absplw = rval(:,:,:,1,:)
1179 :
1180 6144 : deallocate(rval)
1181 :
1182 6144 : ierr = pio_inq_varid(ncid, 'refindex_real_sw', vid)
1183 6144 : ierr = pio_get_var(ncid, vid, props%refrtabsw)
1184 :
1185 6144 : ierr = pio_inq_varid(ncid, 'refindex_im_sw', vid)
1186 6144 : ierr = pio_get_var(ncid, vid, props%refitabsw)
1187 :
1188 6144 : ierr = pio_inq_varid(ncid, 'refindex_real_lw', vid)
1189 6144 : ierr = pio_get_var(ncid, vid, props%refrtablw)
1190 :
1191 6144 : ierr = pio_inq_varid(ncid, 'refindex_im_lw', vid)
1192 6144 : ierr = pio_get_var(ncid, vid, props%refitablw)
1193 :
1194 6144 : ierr = pio_inq_varid(ncid, 'sigmag', vid)
1195 6144 : ierr = pio_get_var(ncid, vid, props%sigmag)
1196 :
1197 6144 : ierr = pio_inq_varid(ncid, 'dgnum', vid)
1198 6144 : ierr = pio_get_var(ncid, vid, props%dgnum)
1199 :
1200 6144 : ierr = pio_inq_varid(ncid, 'dgnumlo', vid)
1201 6144 : ierr = pio_get_var(ncid, vid, props%dgnumlo)
1202 :
1203 6144 : ierr = pio_inq_varid(ncid, 'dgnumhi', vid)
1204 6144 : ierr = pio_get_var(ncid, vid, props%dgnumhi)
1205 :
1206 6144 : ierr = pio_inq_varid(ncid, 'rhcrystal', vid)
1207 6144 : ierr = pio_get_var(ncid, vid, props%rhcrystal)
1208 :
1209 6144 : ierr = pio_inq_varid(ncid, 'rhdeliques', vid)
1210 6144 : ierr = pio_get_var(ncid, vid, props%rhdeliques)
1211 :
1212 6144 : end subroutine modal_optics_init
1213 :
1214 : !================================================================================================
1215 :
1216 0 : subroutine bin_optics_init(props, ncid)
1217 :
1218 : ! Read optics data for modal aerosols
1219 :
1220 : type (physprop_type), intent(inout) :: props ! storage for file data
1221 : type (file_desc_T), intent(inout) :: ncid ! indentifier for netcdf file
1222 :
1223 : ! Local variables
1224 : integer :: ierr
1225 : integer :: did
1226 : integer :: ival
1227 : type(var_desc_t) :: vid
1228 :
1229 : character(len=*), parameter :: subname = 'bin_optics_init'
1230 : !------------------------------------------------------------------------------------
1231 :
1232 : ! Check dimensions for number of lw and sw bands
1233 :
1234 0 : ierr = pio_inq_dimid(ncid, 'lw_band', did)
1235 0 : ierr = pio_inq_dimlen(ncid, did, ival)
1236 0 : if (ival .ne. nlwbands) call endrun(subname//':'//props%sourcefile// &
1237 0 : ' has the wrong number of lw bands')
1238 :
1239 0 : ierr = pio_inq_dimid(ncid, 'sw_band', did)
1240 0 : ierr = pio_inq_dimlen(ncid, did, ival)
1241 0 : if (ival .ne. nswbands) call endrun(subname//':'//props%sourcefile// &
1242 0 : ' has the wrong number of sw bands')
1243 :
1244 : ! Get other dimensions
1245 0 : ierr = pio_inq_dimid(ncid, 'corefrac', did)
1246 0 : ierr = pio_inq_dimlen(ncid, did, props%nfrac)
1247 :
1248 :
1249 : ! Allocate arrays
1250 : allocate( &
1251 : props%extpsw2(props%nfrac,nswbands), &
1252 : props%abspsw2(props%nfrac,nswbands), &
1253 : props%asmpsw2(props%nfrac,nswbands), &
1254 : props%absplw2(props%nfrac,nlwbands), &
1255 0 : props%corefrac(props%nfrac) )
1256 :
1257 0 : ierr = pio_inq_varid(ncid, 'extpsw2', vid)
1258 0 : ierr = pio_get_var(ncid, vid, props%extpsw2)
1259 :
1260 0 : ierr = pio_inq_varid(ncid, 'abspsw2', vid)
1261 0 : ierr = pio_get_var(ncid, vid, props%abspsw2)
1262 :
1263 0 : ierr = pio_inq_varid(ncid, 'asmpsw2', vid)
1264 0 : ierr = pio_get_var(ncid, vid, props%asmpsw2)
1265 :
1266 0 : ierr = pio_inq_varid(ncid, 'absplw2', vid)
1267 0 : ierr = pio_get_var(ncid, vid, props%absplw2)
1268 :
1269 0 : ierr = pio_inq_varid(ncid, 'corefrac', vid)
1270 0 : ierr = pio_get_var(ncid, vid, props%corefrac)
1271 :
1272 0 : end subroutine bin_optics_init
1273 :
1274 :
1275 : !================================================================================================
1276 :
1277 0 : subroutine bindef_optics_init(props, ncid)
1278 :
1279 : ! Read optics data for modal aerosols
1280 :
1281 : type (physprop_type), intent(inout) :: props ! storage for file data
1282 : type (file_desc_T), intent(inout) :: ncid ! indentifier for netcdf file
1283 :
1284 : ! Local variables
1285 : integer :: ierr
1286 : integer :: did
1287 : integer :: ival
1288 : type(var_desc_t) :: vid
1289 :
1290 : character(len=*), parameter :: subname = 'bin_optics_init'
1291 : !------------------------------------------------------------------------------------
1292 :
1293 : ! Check dimensions for number of lw and sw bands
1294 :
1295 0 : ierr = pio_inq_dimid(ncid, 'lw_band', did)
1296 0 : ierr = pio_inq_dimlen(ncid, did, ival)
1297 0 : if (ival .ne. nlwbands) call endrun(subname//':'//props%sourcefile// &
1298 0 : ' has the wrong number of lw bands')
1299 :
1300 0 : ierr = pio_inq_dimid(ncid, 'sw_band', did)
1301 0 : ierr = pio_inq_dimlen(ncid, did, ival)
1302 0 : if (ival .ne. nswbands) call endrun(subname//':'//props%sourcefile// &
1303 0 : ' has the wrong number of sw bands')
1304 :
1305 0 : ierr = pio_inq_varid(ncid, 'density', vid)
1306 0 : ierr = pio_get_var(ncid, vid, props%density_aer)
1307 :
1308 0 : ierr = pio_inq_varid(ncid, 'hygroscopicity', vid)
1309 0 : ierr = pio_get_var(ncid, vid, props%hygro_aer)
1310 :
1311 : ! read refractive index data if available
1312 0 : call refindex_aer_init(props, ncid)
1313 :
1314 0 : end subroutine bindef_optics_init
1315 :
1316 : !================================================================================================
1317 :
1318 13824 : subroutine bulk_props_init(physprop, nc_id)
1319 :
1320 : ! Read props for bulk aerosols
1321 :
1322 : type (physprop_type), intent(inout) :: physprop ! storage for file data
1323 : type (file_desc_T), intent(inout) :: nc_id ! indentifier for netcdf file
1324 :
1325 : ! Local variables
1326 : integer :: ierr
1327 :
1328 : type(var_desc_T) :: vid
1329 :
1330 : logical :: debug = .false.
1331 :
1332 : character(len=*), parameter :: subname = 'bulk_props_init'
1333 : !------------------------------------------------------------------------------------
1334 :
1335 : ! read microphys
1336 13824 : ierr = pio_inq_varid(nc_id, 'name', vid)
1337 13824 : ierr = pio_get_var(nc_id, vid, physprop%aername)
1338 :
1339 : ! use GLC function to remove trailing nulls and blanks.
1340 : ! physprop%aername = aername_str(:GLC(aername_str))
1341 :
1342 13824 : ierr = pio_inq_varid(nc_id, 'density', vid)
1343 13824 : ierr = pio_get_var(nc_id, vid, physprop%density_aer)
1344 :
1345 13824 : ierr = pio_inq_varid(nc_id, 'sigma_logr', vid)
1346 13824 : ierr = pio_get_var(nc_id, vid, physprop%dispersion_aer)
1347 :
1348 13824 : ierr = pio_inq_varid(nc_id, 'dryrad', vid)
1349 13824 : ierr = pio_get_var(nc_id, vid, physprop%dryrad_aer)
1350 :
1351 13824 : ierr = pio_inq_varid(nc_id, 'hygroscopicity', vid)
1352 13824 : ierr = pio_get_var(nc_id, vid, physprop%hygro_aer)
1353 :
1354 13824 : ierr = pio_inq_varid(nc_id, 'num_to_mass_ratio', vid)
1355 13824 : ierr = pio_get_var(nc_id, vid, physprop%num_to_mass_aer)
1356 :
1357 : ! Output select data to log file
1358 13824 : if (debug .and. masterproc .and. idx_sw_diag > 0) then
1359 0 : if (trim(physprop%aername) == 'SULFATE') then
1360 0 : write(iulog, '(2x, a)') '_______ hygroscopic growth in visible band _______'
1361 0 : call aer_optics_log_rh('SO4', physprop%sw_hygro_ext(:,idx_sw_diag), &
1362 0 : physprop%sw_hygro_ssa(:,idx_sw_diag), physprop%sw_hygro_asm(:,idx_sw_diag))
1363 : end if
1364 0 : write(iulog, *) subname//': finished for ', trim(physprop%aername)
1365 : end if
1366 :
1367 13824 : end subroutine bulk_props_init
1368 :
1369 : !================================================================================================
1370 :
1371 138240000 : function exp_interpol(x, f, y) result(g)
1372 : ! Purpose:
1373 : ! interpolates f(x) to point y
1374 : ! assuming f(x) = f(x0) exp a(x - x0)
1375 : ! where a = ( ln f(x1) - ln f(x0) ) / (x1 - x0)
1376 : ! x0 <= x <= x1
1377 : ! assumes x is monotonically increasing
1378 : ! Author: D. Fillmore
1379 :
1380 : implicit none
1381 :
1382 : real(r8), intent(in), dimension(:) :: x ! grid points
1383 : real(r8), intent(in), dimension(:) :: f ! grid function values
1384 : real(r8), intent(in) :: y ! interpolation point
1385 : real(r8) :: g ! interpolated function value
1386 :
1387 : integer :: k ! interpolation point index
1388 : integer :: n ! length of x
1389 : real(r8) :: a
1390 :
1391 138240000 : n = size(x)
1392 :
1393 : ! find k such that x(k) < y =< x(k+1)
1394 : ! set k = 1 if y <= x(1) and k = n-1 if y > x(n)
1395 :
1396 138240000 : if (y <= x(1)) then
1397 : k = 1
1398 138101760 : else if (y >= x(n)) then
1399 1382400 : k = n - 1
1400 : else
1401 : k = 1
1402 289612800 : do while (y > x(k+1) .and. k < n)
1403 136719360 : k = k + 1
1404 : end do
1405 : end if
1406 :
1407 : ! interpolate
1408 138240000 : a = ( log( f(k+1) / f(k) ) ) / ( x(k+1) - x(k) )
1409 138240000 : g = f(k) * exp( a * (y - x(k)) )
1410 : return
1411 : end function exp_interpol
1412 :
1413 : !================================================================================================
1414 :
1415 129024000 : function lin_interpol(x, f, y) result(g)
1416 : ! Purpose:
1417 : ! interpolates f(x) to point y
1418 : ! assuming f(x) = f(x0) + a * (x - x0)
1419 : ! where a = ( f(x1) - f(x0) ) / (x1 - x0)
1420 : ! x0 <= x <= x1
1421 : ! assumes x is monotonically increasing
1422 : ! Author: D. Fillmore
1423 :
1424 : implicit none
1425 :
1426 : real(r8), intent(in), dimension(:) :: x ! grid points
1427 : real(r8), intent(in), dimension(:) :: f ! grid function values
1428 : real(r8), intent(in) :: y ! interpolation point
1429 : real(r8) :: g ! interpolated function value
1430 :
1431 : integer :: k ! interpolation point index
1432 : integer :: n ! length of x
1433 : real(r8) :: a
1434 :
1435 129024000 : n = size(x)
1436 :
1437 : ! find k such that x(k) < y =< x(k+1)
1438 : ! set k = 1 if y <= x(1) and k = n-1 if y > x(n)
1439 :
1440 129024000 : if (y <= x(1)) then
1441 : k = 1
1442 128894976 : else if (y >= x(n)) then
1443 1290240 : k = n - 1
1444 : else
1445 : k = 1
1446 270305280 : do while (y > x(k+1) .and. k < n)
1447 127604736 : k = k + 1
1448 : end do
1449 : end if
1450 :
1451 : ! interpolate
1452 129024000 : a = ( f(k+1) - f(k) ) / ( x(k+1) - x(k) )
1453 129024000 : g = f(k) + a * (y - x(k))
1454 : return
1455 : end function lin_interpol
1456 :
1457 : !================================================================================================
1458 :
1459 : subroutine aer_optics_log(name, ext, ssa, asm)
1460 :
1461 : ! Purpose:
1462 : ! write aerosol optical constants to log file
1463 :
1464 : ! Author: D. Fillmore
1465 :
1466 : character(len=*), intent(in) :: name
1467 : real(r8), intent(in) :: ext(:)
1468 : real(r8), intent(in) :: ssa(:)
1469 : real(r8), intent(in) :: asm(:)
1470 :
1471 : integer :: kbnd, nbnd
1472 : !------------------------------------------------------------------------------------
1473 :
1474 : nbnd = ubound(ext, 1)
1475 :
1476 : write(iulog, '(2x, a)') name
1477 : write(iulog, '(2x, a, 4x, a, 4x, a, 4x, a)') 'SW band', 'ext (m^2 kg^-1)', ' ssa', ' asm'
1478 : do kbnd = 1, nbnd
1479 : write(iulog, '(2x, i7, 4x, f13.2, 4x, f4.2, 4x, f4.2)') kbnd, ext(kbnd), ssa(kbnd), asm(kbnd)
1480 : end do
1481 :
1482 : end subroutine aer_optics_log
1483 :
1484 : !================================================================================================
1485 :
1486 :
1487 0 : subroutine aer_optics_log_rh(name, ext, ssa, asm)
1488 :
1489 : ! Purpose:
1490 : ! write out aerosol optical properties
1491 : ! for a set of test rh values
1492 : ! to test hygroscopic growth interpolation
1493 :
1494 : ! Author: D. Fillmore
1495 :
1496 : character(len=*), intent(in) :: name
1497 : real(r8), intent(in) :: ext(nrh)
1498 : real(r8), intent(in) :: ssa(nrh)
1499 : real(r8), intent(in) :: asm(nrh)
1500 :
1501 : integer :: krh_test
1502 : integer, parameter :: nrh_test = 36
1503 : integer :: krh
1504 : real(r8) :: rh
1505 : real(r8) :: rh_test(nrh_test)
1506 : real(r8) :: exti
1507 : real(r8) :: ssai
1508 : real(r8) :: asmi
1509 : real(r8) :: wrh
1510 : !------------------------------------------------------------------------------------
1511 :
1512 0 : do krh_test = 1, nrh_test
1513 0 : rh_test(krh_test) = sqrt(sqrt(sqrt(sqrt(((krh_test - 1.0_r8) / (nrh_test - 1))))))
1514 : enddo
1515 0 : write(iulog, '(2x, a)') name
1516 0 : write(iulog, '(2x, a, 4x, a, 4x, a, 4x, a)') ' rh', 'ext (m^2 kg^-1)', ' ssa', ' asm'
1517 :
1518 : ! loop through test rh values
1519 0 : do krh_test = 1, nrh_test
1520 : ! find corresponding rh index
1521 0 : rh = rh_test(krh_test)
1522 0 : krh = min(floor( (rh) * nrh ) + 1, nrh - 1)
1523 0 : wrh = (rh) *nrh - krh
1524 0 : exti = ext(krh + 1) * (wrh + 1) - ext(krh) * wrh
1525 0 : ssai = ssa(krh + 1) * (wrh + 1) - ssa(krh) * wrh
1526 0 : asmi = asm(krh + 1) * (wrh + 1) - asm(krh) * wrh
1527 0 : write(iulog, '(2x, f5.3, 4x, f13.3, 4x, f5.3, 4x, f5.3)') rh_test(krh_test), exti, ssai, asmi
1528 : end do
1529 :
1530 0 : end subroutine aer_optics_log_rh
1531 :
1532 :
1533 : !================================================================================================
1534 :
1535 0 : subroutine hygroscopic_coreshell_optics_init(phys_prop, nc_id)
1536 :
1537 : ! Read optics data of type 'hygroscopic_coreshell' and interpolate it to CAM's rh mesh.
1538 :
1539 : type (physprop_type), intent(inout) :: phys_prop ! data after interp onto cam rh mesh
1540 : type (file_desc_T), intent(inout) :: nc_id ! indentifier for netcdf file
1541 :
1542 : ! Local variables
1543 : integer :: ierr ! error flag
1544 :
1545 : integer :: rh_id, lw_band_id, sw_band_id, coreshell_id, dstbc_id, kap_id
1546 : integer :: sw_ext_id, sw_ssa_id, sw_asm_id, lw_abs_id
1547 : integer :: nbnd, swbands, did
1548 :
1549 : ! temp data from hygroscopic file before interpolation onto cam-rh-mesh
1550 : integer :: nrh ! number of rh values in file
1551 : integer :: nfrac ! number of core/shell ratio values in file
1552 : integer :: nbcdust,nkap
1553 :
1554 : real(r8) :: rh ! real rh value on cam rh mesh (indexvalue)
1555 : character(len=*), parameter :: sub = 'hygroscopic_coreshell_optics_init'
1556 : !------------------------------------------------------------------------------------
1557 :
1558 0 : if (masterproc) then
1559 0 : write(iulog,*) 'hygroscopic_coreshell_optics_init: Read file '//trim(phys_prop%sourcefile)
1560 : endif
1561 :
1562 0 : ierr = pio_inq_dimid(nc_id, 'lw_band', lw_band_id)
1563 0 : ierr = pio_inq_dimlen(nc_id, lw_band_id, nbnd)
1564 0 : if (nbnd .ne. nlwbands) call endrun(trim(phys_prop%sourcefile)// &
1565 0 : ' has the wrong number of lwbands')
1566 :
1567 0 : ierr = pio_inq_dimid(nc_id, 'sw_band', sw_band_id)
1568 0 : ierr = pio_inq_dimlen(nc_id, sw_band_id, swbands)
1569 0 : if(swbands .ne. nswbands) call endrun(trim(phys_prop%sourcefile)// &
1570 0 : ' has the wrong number of sw bands')
1571 :
1572 :
1573 0 : ierr = pio_inq_dimid(nc_id, 'coreshellratio', did)
1574 0 : ierr = pio_inq_dimlen(nc_id, did, phys_prop%nfrac)
1575 :
1576 0 : ierr = pio_inq_dimid(nc_id, 'dstbcratio', did)
1577 0 : ierr = pio_inq_dimlen(nc_id, did, phys_prop%nbcdust)
1578 :
1579 0 : ierr = pio_inq_dimid(nc_id, 'kap', did)
1580 0 : ierr = pio_inq_dimlen(nc_id, did, phys_prop%nkap)
1581 :
1582 0 : ierr = pio_inq_dimid(nc_id, 'rh_idx', rh_id)
1583 0 : ierr = pio_inq_dimlen(nc_id, rh_id, phys_prop%nrelh)
1584 :
1585 : allocate(phys_prop%sw_hygro_coreshell_ext(phys_prop%nrelh,nswbands, &
1586 0 : phys_prop%nfrac,phys_prop%nbcdust,phys_prop%nkap))
1587 : allocate(phys_prop%sw_hygro_coreshell_ssa(phys_prop%nrelh,nswbands, &
1588 0 : phys_prop%nfrac,phys_prop%nbcdust,phys_prop%nkap))
1589 : allocate(phys_prop%sw_hygro_coreshell_asm(phys_prop%nrelh,nswbands, &
1590 0 : phys_prop%nfrac,phys_prop%nbcdust,phys_prop%nkap))
1591 : allocate(phys_prop%lw_hygro_coreshell_abs(phys_prop%nrelh,nlwbands, &
1592 0 : phys_prop%nfrac,phys_prop%nbcdust,phys_prop%nkap))
1593 0 : allocate(phys_prop%corefrac(phys_prop%nfrac))
1594 0 : allocate(phys_prop%bcdust(phys_prop%nbcdust))
1595 0 : allocate(phys_prop%kap(phys_prop%nkap))
1596 0 : allocate(phys_prop%relh(phys_prop%nrelh))
1597 :
1598 0 : ierr = pio_inq_varid(nc_id, 'rh', rh_id)
1599 0 : ierr = pio_inq_varid(nc_id, 'coreshellratio', coreshell_id) ! modified by Pengfei for coreshell
1600 0 : ierr = pio_inq_varid(nc_id, 'dstbcratio', dstbc_id) ! modified by Pengfei for coreshell
1601 0 : ierr = pio_inq_varid(nc_id, 'kap', kap_id)
1602 :
1603 0 : ierr = pio_inq_varid(nc_id, 'ext_sw_coreshell', sw_ext_id)
1604 0 : ierr = pio_inq_varid(nc_id, 'ssa_sw_coreshell', sw_ssa_id)
1605 0 : ierr = pio_inq_varid(nc_id, 'asm_sw_coreshell', sw_asm_id)
1606 0 : ierr = pio_inq_varid(nc_id, 'abs_lw_coreshell', lw_abs_id)
1607 :
1608 0 : ierr = pio_get_var(nc_id, sw_ext_id, phys_prop%sw_hygro_coreshell_ext)
1609 0 : ierr = pio_get_var(nc_id, sw_ssa_id, phys_prop%sw_hygro_coreshell_ssa)
1610 0 : ierr = pio_get_var(nc_id, sw_asm_id, phys_prop%sw_hygro_coreshell_asm)
1611 0 : ierr = pio_get_var(nc_id, lw_abs_id, phys_prop%lw_hygro_coreshell_abs)
1612 0 : ierr = pio_get_var(nc_id, kap_id, phys_prop%kap)
1613 0 : ierr = pio_get_var(nc_id, rh_id, phys_prop%relh)
1614 0 : ierr = pio_get_var(nc_id, dstbc_id, phys_prop%bcdust)
1615 0 : ierr = pio_get_var(nc_id, coreshell_id, phys_prop%corefrac)
1616 :
1617 : ! read refractive index data if available
1618 0 : call refindex_aer_init(phys_prop, nc_id)
1619 :
1620 0 : end subroutine hygroscopic_coreshell_optics_init
1621 :
1622 : !================================================================================================
1623 :
1624 0 : subroutine hygroscopic_wtp_optics_init(phys_prop, nc_id)
1625 :
1626 : ! Read optics data of type 'hygroscopic' and interpolate it to CAM's rh mesh.
1627 :
1628 : type (physprop_type), intent(inout) :: phys_prop ! data after interp onto cam rh mesh
1629 : type (file_desc_T), intent(inout) :: nc_id ! indentifier for netcdf file
1630 :
1631 : ! Local variables
1632 : integer :: ierr ! error flag
1633 :
1634 : integer :: lw_band_id, sw_band_id, did
1635 : integer :: sw_ext_wtp_id, sw_ssa_wtp_id, sw_asm_wtp_id, lw_ext_wtp_id, wtp_id
1636 : integer :: nbnd, swbands
1637 :
1638 : real(r8) :: rh ! real rh value on cam rh mesh (indexvalue)
1639 : character(len=*), parameter :: sub = 'hygroscopic_wtp_optics_init'
1640 : !------------------------------------------------------------------------------------
1641 :
1642 : !st
1643 : ! Get other dimensions
1644 0 : ierr = pio_inq_dimid(nc_id, 'wgtpct', did)
1645 0 : ierr = pio_inq_dimlen(nc_id, did, phys_prop%nwtp)
1646 :
1647 :
1648 0 : allocate(phys_prop%sw_hygro_ext_wtp(phys_prop%nwtp,nswbands))
1649 0 : allocate(phys_prop%sw_hygro_ssa_wtp(phys_prop%nwtp,nswbands))
1650 0 : allocate(phys_prop%sw_hygro_asm_wtp(phys_prop%nwtp,nswbands))
1651 0 : allocate(phys_prop%lw_hygro_abs_wtp(phys_prop%nwtp,nlwbands))
1652 0 : allocate(phys_prop%wgtpct(phys_prop%nwtp))
1653 :
1654 0 : ierr = pio_inq_dimid(nc_id, 'lw_band', lw_band_id)
1655 0 : ierr = pio_inq_dimlen(nc_id, lw_band_id, nbnd)
1656 0 : if (nbnd .ne. nlwbands) call endrun(phys_prop%sourcefile// &
1657 0 : ' has the wrong number of lwbands')
1658 :
1659 0 : ierr = pio_inq_dimid(nc_id, 'sw_band', sw_band_id)
1660 0 : ierr = pio_inq_dimlen(nc_id, sw_band_id, swbands)
1661 0 : if(swbands .ne. nswbands) call endrun(phys_prop%sourcefile// &
1662 0 : ' has the wrong number of sw bands')
1663 :
1664 0 : ierr = pio_inq_varid(nc_id, 'ext_sw_wtp', sw_ext_wtp_id)
1665 0 : ierr = pio_inq_varid(nc_id, 'ssa_sw_wtp', sw_ssa_wtp_id)
1666 0 : ierr = pio_inq_varid(nc_id, 'asm_sw_wtp', sw_asm_wtp_id)
1667 0 : ierr = pio_inq_varid(nc_id, 'abs_lw_wtp', lw_ext_wtp_id)
1668 0 : ierr = pio_inq_varid(nc_id, 'wgtpct', wtp_id)
1669 :
1670 0 : ierr = pio_get_var(nc_id, sw_ext_wtp_id, phys_prop%sw_hygro_ext_wtp)
1671 0 : ierr = pio_get_var(nc_id, sw_ssa_wtp_id, phys_prop%sw_hygro_ssa_wtp)
1672 0 : ierr = pio_get_var(nc_id, sw_asm_wtp_id, phys_prop%sw_hygro_asm_wtp)
1673 0 : ierr = pio_get_var(nc_id, lw_ext_wtp_id, phys_prop%lw_hygro_abs_wtp)
1674 0 : ierr = pio_get_var(nc_id, wtp_id, phys_prop%wgtpct)
1675 :
1676 : ! read refractive index data if available
1677 0 : call refindex_aer_init(phys_prop, nc_id)
1678 :
1679 : ! read bulk aero props
1680 0 : call bulk_props_init(phys_prop, nc_id)
1681 :
1682 0 : end subroutine hygroscopic_wtp_optics_init
1683 :
1684 :
1685 0 : end module phys_prop
|
