CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
cires_tauamf_data.F90
1
3
5
6 use machine, only: kind_phys
7!...........................................................................................
8! tabulated GW-sources: GRACILE/Ern et al., 2018 and/or Resolved GWs from C384-Annual run
9!...........................................................................................
10implicit none
11
12 integer :: ntau_d1y, ntau_d2t
13 real(kind=kind_phys), allocatable :: ugwp_taulat(:)
14 real(kind=kind_phys), allocatable :: tau_limb(:,:), days_limb(:)
15 logical :: flag_alloctau = .false.
16 character(len=255):: ugwp_taufile = 'ugwp_limb_tau.nc'
17
18 public :: read_tau_amf, cires_indx_ugwp, tau_amf_interp
19
20contains
21
23 subroutine read_tau_amf(me, master, errmsg, errflg)
24
25 use netcdf
26 integer, intent(in) :: me, master
27 integer :: ncid, iernc, vid, dimid, status
28 integer :: k
29
30 character(len=*), intent(out) :: errmsg
31 integer, intent(out) :: errflg
32!
33
34 iernc=nf90_open(trim(ugwp_taufile), nf90_nowrite, ncid)
35
36 if(iernc.ne.0) then
37 write(errmsg,'(*(a))') "read_tau_amf: cannot open file_limb_tab data-file ", &
38 trim(ugwp_taufile)
39 print *, 'cannot open ugwp-v1 tau-file=',trim(ugwp_taufile)
40 errflg = 1
41 return
42 else
43
44
45 status = nf90_inq_dimid(ncid, "lat", dimid)
46! if (status /= nf90_noerr) call handle_err(status)
47!
48 status = nf90_inquire_dimension(ncid, dimid, len =ntau_d1y )
49
50 status = nf90_inq_dimid(ncid, "days", dimid)
51 status = nf90_inquire_dimension(ncid, dimid, len =ntau_d2t )
52
53 if (me == master) print *, ntau_d1y, ntau_d2t, ' dimd of tau_ngw ugwp-v1 '
54 if (ntau_d2t .le. 0 .or. ntau_d1y .le. 0) then
55 print *, 'ugwp-v1 tau-file=', trim(ugwp_taufile)
56 print *, ' ugwp-v1: ', 'ntau_d2t=',ntau_d2t, 'ntau_d2t=',ntau_d1y
57 stop
58 endif
59
60 if (.not.allocated(ugwp_taulat)) allocate (ugwp_taulat(ntau_d1y ))
61 if (.not.allocated(days_limb)) allocate (days_limb(ntau_d2t))
62 if (.not.allocated(tau_limb)) allocate (tau_limb(ntau_d1y, ntau_d2t ))
63
64 iernc=nf90_inq_varid( ncid, 'DAYS', vid )
65 iernc= nf90_get_var( ncid, vid, days_limb)
66 iernc=nf90_inq_varid( ncid, 'LATS', vid )
67 iernc= nf90_get_var( ncid, vid, ugwp_taulat)
68 iernc=nf90_inq_varid( ncid, 'ABSMF', vid )
69 iernc= nf90_get_var( ncid, vid, tau_limb)
70
71 iernc=nf90_close(ncid)
72
73 endif
74
75 end subroutine read_tau_amf
76
78 subroutine cires_indx_ugwp (npts, me, master, dlat,j1_tau,j2_tau, w1_j1tau, w2_j2tau)
79
80 use machine, only: kind_phys
81
82 implicit none
83
84 integer, intent(in) :: npts, me, master
85 real(kind=kind_phys) , dimension(npts), intent(in) :: dlat
86
87 integer, dimension(npts), intent(inout) :: j1_tau, j2_tau
88 real(kind=kind_phys) , dimension(npts), intent(inout) :: w1_j1tau, w2_j2tau
89
90!locals
91
92 integer :: i,j, j1, j2
93!
94 do j=1,npts
95 j2_tau(j) = ntau_d1y
96 do i=1,ntau_d1y
97 if (dlat(j) < ugwp_taulat(i)) then
98 j2_tau(j) = i
99 exit
100 endif
101 enddo
102
103
104 j2_tau(j) = min(j2_tau(j),ntau_d1y)
105 j1_tau(j) = max(j2_tau(j)-1,1)
106
107 if (j1_tau(j) /= j2_tau(j) ) then
108 w2_j2tau(j) = (dlat(j) - ugwp_taulat(j1_tau(j))) &
109 / (ugwp_taulat(j2_tau(j))-ugwp_taulat(j1_tau(j)))
110 else
111 w2_j2tau(j) = 1.0
112 endif
113 w1_j1tau(j) = 1.0 - w2_j2tau(j)
114 enddo
115 return
116 end subroutine cires_indx_ugwp
117
119 subroutine tau_amf_interp(me, master, im, idate, fhour, j1_tau,j2_tau, ddy_j1, ddy_j2, tau_ddd)
120 use machine, only: kind_phys
121 implicit none
122
123!input
124 integer, intent(in) :: me, master
125 integer, intent(in) :: im, idate(4)
126 real(kind=kind_phys), intent(in) :: fhour
127
128 real(kind=kind_phys), intent(in), dimension(im) :: ddy_j1, ddy_j2
129 integer , intent(in), dimension(im) :: j1_tau,j2_tau
130!ouput
131 real(kind=kind_phys), dimension(im) :: tau_ddd
132!locals
133
134 integer :: i, j1, j2, it1, it2 , iday
135 integer :: ddd
136 real(kind=kind_phys) :: tx1, tx2, w1, w2, fddd
137!
138! define day of year ddd ..... from the old-fashioned "GFS-style"
139!
140 call gfs_idate_calendar(idate, fhour, ddd, fddd)
141
142 it1 = 2
143 do iday=1, ntau_d2t
144 if (fddd .lt. days_limb(iday) ) then
145 it2 = iday
146 exit
147 endif
148 enddo
149
150 it2 = min(it2,ntau_d2t)
151 it1 = max(it2-1,1)
152 if (it2 > ntau_d2t ) then
153 print *, ' Error in time-interpolation for tau_amf_interp '
154 print *, ' it1, it2, ntau_d2t ', it1, it2, ntau_d2t
155 print *, ' Error in time-interpolation see cires_tauamf_data.F90 '
156 stop
157 endif
158
159 w2 = (fddd-days_limb(it1))/(days_limb(it2)-days_limb(it1))
160 w1 = 1.0-w2
161
162 do i=1, im
163 j1 = j1_tau(i)
164 j2 = j2_tau(i)
165 tx1 = tau_limb(j1, it1)*ddy_j1(i)+tau_limb(j2, it1)*ddy_j2(i)
166 tx2 = tau_limb(j1, it2)*ddy_j1(i)+tau_limb(j2, it2)*ddy_j2(i)
167 tau_ddd(i) = tx1*w1 + w2*tx2
168 enddo
169
170 end subroutine tau_amf_interp
171
173 subroutine gfs_idate_calendar(idate, fhour, ddd, fddd)
174
175 use machine, only: kind_phys
176 implicit none
177! input
178 integer, intent(in) :: idate(4)
179 real(kind=kind_phys), intent(in) :: fhour
180!out
181 integer, intent(out) :: ddd
182 real(kind=kind_phys), intent(out) :: fddd
183!
184!locals
185!
186 real(kind=kind_phys) :: rjday
187 integer :: jdow, jdoy, jday
188 real(8) :: rinc(5)
189 real(4) :: rinc4(5)
190
191 integer :: iw3jdn
192 integer :: jd1, jddd
193
194 integer idat(8),jdat(8)
195
196
197 idat(1:8) = 0
198 idat(1) = idate(4)
199 idat(2) = idate(2)
200 idat(3) = idate(3)
201 idat(5) = idate(1)
202 rinc(1:5) = 0.
203 rinc(2) = fhour
204!
205 call w3movdat(rinc, idat,jdat)
206! jdate(8)- date and time (yr, mo, day, [tz], hr, min, sec)
207 jdow = 0
208 jdoy = 0
209 jday = 0
210 call w3doxdat(jdat,jdow, ddd, jday)
211 fddd = float(ddd) + jdat(5) / 24.
212 end subroutine gfs_idate_calendar
213
214end module cires_tauamf_data