LCOV - code coverage report
Current view: | top level - physics/cam - aerosol_optics_cam.F90 (source / functions) |  |  | Hit | Total | Coverage
Test: | coverage.info |  | Lines: | 554 | 672 | 82.4 %
Date: | 2024-12-17 22:39:59 |  | Functions: | 9 | 14 | 64.3 %

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

Generated by: LCOV version 1.14