CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
cires_ugwp_triggers.F90
1
3
5 contains
6!
7
9 subroutine slat_geos5_tamp_v0(im, tau_amp, xlatdeg, tau_gw)
10!=================
11! GEOS-5 & MERRA-2 lat-dependent GW-source function tau(z=Zlaunch) =rho*<u'w'>
12!=================
13 implicit none
14 integer :: im
15 real :: tau_amp, xlatdeg(im), tau_gw(im)
16 real :: latdeg, flat_gw, tem
17 integer :: i
18
19!
20! if-lat
21!
22 do i=1, im
23 latdeg = abs(xlatdeg(i))
24 if (latdeg < 15.3) then
25 tem = (latdeg-3.0) / 8.0
26 flat_gw = 0.75 * exp(-tem * tem)
27 if (flat_gw < 1.2 .and. latdeg <= 3.0) flat_gw = 0.75
28 elseif (latdeg < 31.0 .and. latdeg >= 15.3) then
29 flat_gw = 0.10
30 elseif (latdeg < 60.0 .and. latdeg >= 31.0) then
31 tem = (latdeg-60.0) / 23.0
32 flat_gw = 0.50 * exp(- tem * tem)
33 elseif (latdeg >= 60.0) then
34 tem = (latdeg-60.0) / 70.0
35 flat_gw = 0.50 * exp(- tem * tem)
36 endif
37 tau_gw(i) = tau_amp*flat_gw
38 enddo
39!
40 end subroutine slat_geos5_tamp_v0
41
43 subroutine slat_geos5_v0(im, xlatdeg, tau_gw)
44!=================
45! GEOS-5 & MERRA-2 lat-dependent GW-source function tau(z=Zlaunch) =rho*<u'w'>
46!=================
47 implicit none
48 integer :: im
49 real :: xlatdeg(im)
50 real :: tau_gw(im)
51 real :: latdeg
52 real, parameter :: tau_amp = 100.e-3
53 real :: trop_gw, flat_gw
54 integer :: i
55!
56! if-lat
57!
58 trop_gw = 0.75
59 do i=1, im
60 latdeg = xlatdeg(i)
61 if (-15.3 < latdeg .and. latdeg < 15.3) then
62 flat_gw = trop_gw*exp(-( (abs(latdeg)-3.)/8.0)**2)
63 if (flat_gw < 1.2 .and. abs(latdeg) <= 3.) flat_gw = trop_gw
64 else if (latdeg > -31. .and. latdeg <= -15.3) then
65 flat_gw = 0.10
66 else if (latdeg < 31. .and. latdeg >= 15.3) then
67 flat_gw = 0.10
68 else if (latdeg > -60. .and. latdeg <= -31.) then
69 flat_gw = 0.50*exp(-((abs(latdeg)-60.)/23.)**2)
70 else if (latdeg < 60. .and. latdeg >= 31.) then
71 flat_gw = 0.50*exp(-((abs(latdeg)-60.)/23.)**2)
72 else if (latdeg <= -60.) then
73 flat_gw = 0.50*exp(-((abs(latdeg)-60.)/70.)**2)
74 else if (latdeg >= 60.) then
75 flat_gw = 0.50*exp(-((abs(latdeg)-60.)/70.)**2)
76 end if
77 tau_gw(i) = tau_amp*flat_gw
78 enddo
79!
80 end subroutine slat_geos5_v0
81
83 subroutine init_nazdir_v0(naz, xaz, yaz)
84 use ugwp_common_v0 , only : pi2
85 implicit none
86 integer :: naz
87 real, dimension(naz) :: xaz, yaz
88 integer :: idir
89 real :: phic, drad
90 drad = pi2/float(naz)
91 if (naz.ne.4) then
92 do idir =1, naz
93 phic = drad*(float(idir)-1.0)
94 xaz(idir) = cos(phic)
95 yaz(idir) = sin(phic)
96 enddo
97 else
98! if (naz.eq.4) then
99 xaz(1) = 1.0 !E
100 yaz(1) = 0.0
101 xaz(2) = 0.0
102 yaz(2) = 1.0 !N
103 xaz(3) =-1.0 !W
104 yaz(3) = 0.0
105 xaz(4) = 0.0
106 yaz(4) =-1.0 !S
107 endif
108 end subroutine init_nazdir_v0
109 end module cires_ugwp_triggers
Define constants.