source: LMDZ6/branches/LMDZ-ECRAD/libf/phylmd/ecrad/radiation_ice_optics_baran.F90 @ 3880

Last change on this file since 3880 was 3880, checked in by idelkadi, 3 years ago

Online implementation of the radiative transfer code ECRAD in LMDZ.

  • Inclusion of the ecrad directory containing the sources of the ECRAD code
  • Adaptation of compilation scripts (CPP_ECRAD keys)
  • Call of ecrad in radlwsw_m.F90 under the logical key iflag_rrtm = 2
File size: 3.8 KB
Line 
1! radiation_ice_optics_fu.F90 - Scheme for ice optical properties adapted from Baran's data
2!
3! (C) Copyright 2014- ECMWF.
4!
5! This software is licensed under the terms of the Apache Licence Version 2.0
6! which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
7!
8! In applying this licence, ECMWF does not waive the privileges and immunities
9! granted to it by virtue of its status as an intergovernmental organisation
10! nor does it submit to any jurisdiction.
11!
12! Author:  Robin Hogan
13! Email:   r.j.hogan@ecmwf.int
14!
15
16module radiation_ice_optics_baran
17
18  implicit none
19  public
20
21  ! The number of ice coefficients depends on the parameterization
22  integer, parameter :: NIceOpticsCoeffsBaran = 9
23  integer, parameter :: NIceOpticsCoeffsBaran2016 = 5
24
25contains
26
27 
28  !---------------------------------------------------------------------
29  ! Compute ice-particle scattering properties using a
30  ! parameterization as a function of ice water mixing ratio only
31  subroutine calc_ice_optics_baran(nb, coeff, ice_wp, &
32       &  qi, od, scat_od, g)
33
34    use parkind1, only : jprb
35    !use yomhook,  only : lhook, dr_hook
36
37    ! Number of bands
38    integer, intent(in)  :: nb
39    ! Coefficients read from a data file
40    real(jprb), intent(in) :: coeff(:,:)
41    ! Ice water path (kg m-2) and mixing ratio (kg kg-1)
42    real(jprb), intent(in) :: ice_wp, qi
43    ! Total optical depth, scattering optical depth and asymmetry factor
44    real(jprb), intent(out) :: od(nb), scat_od(nb), g(nb)
45
46    !real(jprb) :: hook_handle
47
48    !if (lhook) call dr_hook('radiation_ice_optics:calc_ice_optics_baran',0,hook_handle)
49
50    od  = ice_wp * (coeff(1:nb,1) + coeff(1:nb,2) &
51         &      / (1.0_jprb + qi*coeff(1:nb,3)))
52    scat_od = od * (coeff(1:nb,4) + coeff(1:nb,5) &
53         &      / (1.0_jprb + qi*coeff(1:nb,6)))
54    ! To apply the simple parameterization in Baran et al. (2014), use the
55    ! following instead, but note that it overestimates shortwave absorption:
56    !    od = ice_wp * coeff(1:nb,1)
57    !    scat_od = od * coeff(1:nb,4)
58    g = coeff(1:nb,7) + coeff(1:nb,8) / (1.0_jprb + qi*coeff(1:nb,9))
59
60    !if (lhook) call dr_hook('radiation_ice_optics:calc_ice_optics_baran',1,hook_handle)
61
62  end subroutine calc_ice_optics_baran
63
64
65  !---------------------------------------------------------------------
66  ! Compute ice-particle scattering properties using a
67  ! parameterization as a function of ice water mixing ratio and
68  ! temperature
69  subroutine calc_ice_optics_baran2016(nb, coeff, ice_wp, &
70       &  qi, temperature, od, scat_od, g)
71
72    use parkind1, only : jprb
73    !use yomhook,  only : lhook, dr_hook
74
75    ! Number of bands
76    integer, intent(in)  :: nb
77    ! Coefficients read from a data file
78    real(jprb), intent(in) :: coeff(:,:)
79    ! Ice water path (kg m-2) and mixing ratio (kg kg-1)
80    real(jprb), intent(in) :: ice_wp, qi
81    ! Temperature (K)
82    real(jprb), intent(in) :: temperature
83    ! Total optical depth, scattering optical depth and asymmetry factor
84    real(jprb), intent(out) :: od(nb), scat_od(nb), g(nb)
85   
86    ! Powers of temperature, some multiplied by qi
87    real(jprb) :: qi_T, T2, qi_over_T4
88   
89    !real(jprb) :: hook_handle
90
91    !if (lhook) call dr_hook('radiation_ice_optics:calc_ice_optics_baran2016',0,hook_handle)
92
93    T2 = temperature * temperature
94
95    if (qi < 1.0e-3_jprb) then
96      qi_T = qi * temperature
97      qi_over_T4 = 1.0_jprb / (T2 * T2)
98    else
99      qi_T = 1.0e-3_jprb * temperature
100      qi_over_T4 = 1.0_jprb / (T2 * T2)
101    end if
102
103    od      = ice_wp * coeff(1:nb,1) * qi_over_T4
104    scat_od = od * (coeff(1:nb,2) + coeff(1:nb,3) * qi_T)
105    g       = coeff(1:nb,4) + coeff(1:nb,5) * qi_T
106
107    !if (lhook) call dr_hook('radiation_ice_optics:calc_ice_optics_baran2016',1,hook_handle)
108
109  end subroutine calc_ice_optics_baran2016
110
111end module radiation_ice_optics_baran
Note: See TracBrowser for help on using the repository browser.