-
Notifications
You must be signed in to change notification settings - Fork 12
/
Copy pathhydromet_pdf_parameter_module.F90
112 lines (83 loc) · 4.19 KB
/
hydromet_pdf_parameter_module.F90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
!---------------------------------------------------------------------------
! $Id$
!===============================================================================
module hydromet_pdf_parameter_module
! Description:
! This module defines the derived type hydromet_pdf_parameter.
! References:
! None
!-------------------------------------------------------------------------
use clubb_precision, only: &
core_rknd ! Variable(s)
implicit none
private ! Default scope
public :: hydromet_pdf_parameter, & ! Variable type
init_hydromet_pdf_params ! Procedure
integer, parameter, private :: &
max_hydromet_dim = 8
type hydromet_pdf_parameter
real( kind = core_rknd ), dimension(max_hydromet_dim) :: &
hm_1, & ! Mean of hydrometeor, hm (1st PDF component) [un vary]
hm_2, & ! Mean of hydrometeor, hm (2nd PDF component) [un vary]
mu_hm_1, & ! Mean of hm (1st PDF component) in-precip (ip) [un vary]
mu_hm_2, & ! Mean of hm (2nd PDF component) ip [un vary]
sigma_hm_1, & ! Standard deviation of hm (1st PDF comp.) ip [un vary]
sigma_hm_2, & ! Standard deviation of hm (2nd PDF comp.) ip [un vary]
corr_w_hm_1, & ! Correlation of w and hm (1st PDF component) ip [-]
corr_w_hm_2, & ! Correlation of w and hm (2nd PDF component) ip [-]
corr_chi_hm_1, & ! Correlation of chi and hm (1st PDF component) ip [-]
corr_chi_hm_2, & ! Correlation of chi and hm (2nd PDF component) ip [-]
corr_eta_hm_1, & ! Correlation of eta and hm (1st PDF component) ip [-]
corr_eta_hm_2 ! Correlation of eta and hm (2nd PDF component) ip [-]
real( kind = core_rknd ), dimension(max_hydromet_dim,max_hydromet_dim) :: &
corr_hmx_hmy_1, & ! Correlation of hmx and hmy (1st PDF component) ip [-]
corr_hmx_hmy_2 ! Correlation of hmx and hmy (2nd PDF component) ip [-]
real( kind = core_rknd ) :: &
mu_Ncn_1, & ! Mean of Ncn (1st PDF component) [num/kg]
mu_Ncn_2, & ! Mean of Ncn (2nd PDF component) [num/kg]
sigma_Ncn_1, & ! Standard deviation of Ncn (1st PDF component) [num/kg]
sigma_Ncn_2 ! Standard deviation of Ncn (2nd PDF component) [num/kg]
real( kind = core_rknd ) :: &
precip_frac, & ! Precipitation fraction (overall) [-]
precip_frac_1, & ! Precipitation fraction (1st PDF component) [-]
precip_frac_2 ! Precipitation fraction (2nd PDF component) [-]
end type hydromet_pdf_parameter
contains
!=============================================================================
subroutine init_hydromet_pdf_params( hydromet_pdf_params )
! Description:
! Initialize the elements of hydromet_pdf_params.
! References:
!-----------------------------------------------------------------------
use constants_clubb, only: &
zero ! Constant(s)
implicit none
! Output Variable
type(hydromet_pdf_parameter), intent(out) :: &
hydromet_pdf_params ! Hydrometeor PDF parameters [units vary]
! Initialize hydromet_pdf_params.
hydromet_pdf_params%hm_1 = zero
hydromet_pdf_params%hm_2 = zero
hydromet_pdf_params%mu_hm_1 = zero
hydromet_pdf_params%mu_hm_2 = zero
hydromet_pdf_params%sigma_hm_1 = zero
hydromet_pdf_params%sigma_hm_2 = zero
hydromet_pdf_params%corr_w_hm_1 = zero
hydromet_pdf_params%corr_w_hm_2 = zero
hydromet_pdf_params%corr_chi_hm_1 = zero
hydromet_pdf_params%corr_chi_hm_2 = zero
hydromet_pdf_params%corr_eta_hm_1 = zero
hydromet_pdf_params%corr_eta_hm_2 = zero
hydromet_pdf_params%corr_hmx_hmy_1 = zero
hydromet_pdf_params%corr_hmx_hmy_2 = zero
hydromet_pdf_params%mu_Ncn_1 = zero
hydromet_pdf_params%mu_Ncn_2 = zero
hydromet_pdf_params%sigma_Ncn_1 = zero
hydromet_pdf_params%sigma_Ncn_2 = zero
hydromet_pdf_params%precip_frac = zero
hydromet_pdf_params%precip_frac_1 = zero
hydromet_pdf_params%precip_frac_2 = zero
return
end subroutine init_hydromet_pdf_params
!===============================================================================
end module hydromet_pdf_parameter_module