Line data Source code
1 : module cam_history_buffers
2 : use shr_kind_mod, only: r8 => shr_kind_r8
3 : use cam_history_support, only: max_chars, field_info, hentry, dim_index_2d
4 : use cam_abortutils, only: endrun
5 : use pio, only: var_desc_t
6 :
7 : implicit none
8 :
9 :
10 : contains
11 1495368 : subroutine hbuf_accum_inst (buf8, field, nacs, dimind, idim, flag_xyfill, fillvalue)
12 : !
13 : !-----------------------------------------------------------------------
14 : !
15 : ! Purpose: Accumulate instantaneous values of field in 2-D hbuf.
16 : ! Set accumulation counter to 1.
17 : !
18 : !-----------------------------------------------------------------------
19 : !
20 : real(r8), pointer :: buf8(:,:) ! 2-D history buffer
21 : integer, pointer :: nacs(:) ! accumulation counter
22 : integer, intent(in) :: idim ! Longitude dimension of field array
23 : logical, intent(in) :: flag_xyfill ! non-applicable xy points flagged with fillvalue
24 : real(r8), intent(in ) :: field(idim,*) ! real*8 array
25 : type (dim_index_2d), intent(in ) :: dimind ! 2-D dimension index
26 : real(r8), intent(in) :: fillvalue
27 : !
28 : ! Local indices
29 : !
30 : integer :: ieu, jeu ! number of elements in each dimension
31 : integer :: i, k ! loop indices
32 :
33 : logical :: bad ! flag indicates input field fillvalues not applied consistently
34 : ! with vertical level
35 :
36 1495368 : call dimind%dim_sizes(ieu, jeu)
37 :
38 2990736 : do k=1,jeu
39 26464536 : do i=1,ieu
40 24969168 : buf8(i,k) = field(i,k)
41 : end do
42 : end do
43 :
44 1495368 : if (flag_xyfill) then
45 0 : do i=1,ieu
46 0 : if (field(i,1) == fillvalue) then
47 0 : nacs(i) = 0
48 : else
49 0 : nacs(i) = 1
50 : end if
51 : end do
52 : else
53 1495368 : nacs(1) = 1
54 : end if
55 :
56 1495368 : return
57 : end subroutine hbuf_accum_inst
58 : !#######################################################################
59 :
60 97880040 : subroutine hbuf_accum_add (buf8, field, nacs, dimind, idim, flag_xyfill, fillvalue)
61 : !
62 : !-----------------------------------------------------------------------
63 : !
64 : ! Purpose: Add the values of field to 2-D hbuf.
65 : ! Increment accumulation counter by 1.
66 : !
67 : !-----------------------------------------------------------------------
68 : !
69 : real(r8), pointer :: buf8(:,:) ! 2-D history buffer
70 : integer, pointer :: nacs(:) ! accumulation counter
71 : integer, intent(in) :: idim ! Longitude dimension of field array
72 : logical, intent(in) :: flag_xyfill ! non-applicable xy points flagged with fillvalue
73 : real(r8), intent(in ) :: field(idim,*) ! real*8 array
74 : type (dim_index_2d), intent(in ) :: dimind ! 2-D dimension index
75 : real(r8), intent(in) :: fillvalue
76 : !
77 : ! Local indices
78 : !
79 : integer :: ieu, jeu ! number of elements in each dimension
80 : integer :: i,k ! indices
81 :
82 97880040 : call dimind%dim_sizes(ieu, jeu)
83 :
84 97880040 : if (flag_xyfill) then
85 1498464 : do k=1,jeu
86 13259664 : do i=1,ieu
87 12510432 : if (field(i,k) /= fillvalue) then
88 5880600 : buf8(i,k) = buf8(i,k) + field(i,k)
89 : end if
90 : end do
91 : end do
92 : !
93 : ! Ensure input field has fillvalue defined invariant in the z-direction, then increment nacs
94 : !
95 749232 : call check_accum (field, idim, ieu, jeu, fillvalue)
96 12510432 : do i=1,ieu
97 12510432 : if (field(i,1) /= fillvalue) then
98 5880600 : nacs(i) = nacs(i) + 1
99 : end if
100 : end do
101 : else
102 1166096016 : do k=1,jeu
103 17946363816 : do i=1,ieu
104 17849233008 : buf8(i,k) = buf8(i,k) + field(i,k)
105 : end do
106 : end do
107 97130808 : nacs(1) = nacs(1) + 1
108 : end if
109 :
110 97880040 : return
111 : end subroutine hbuf_accum_add
112 :
113 : !#######################################################################
114 :
115 0 : subroutine hbuf_accum_variance (hbuf, sbuf, field, nacs, dimind, idim, flag_xyfill, fillvalue)
116 : !
117 : !-----------------------------------------------------------------------
118 : !
119 : ! Purpose: accumulate a running variance for standard deviations
120 : !
121 : ! The method of computing variance is from http://www.johndcook.com/blog/standard_deviation/
122 : ! and
123 : ! 1962 paper by B. P. Welford and is presented in Donald Knuth's Art of
124 : ! Computer Programming, Vol 2, page 232, 3rd edition.
125 : !-----------------------------------------------------------------------
126 : !
127 : real(r8), pointer :: hbuf(:,:) ! 2-D history buffer
128 : real(r8), pointer :: sbuf(:,:) ! 2-D history buffer
129 : integer, pointer :: nacs(:) ! accumulation counter
130 : integer, intent(in) :: idim ! Longitude dimension of field array
131 : logical, intent(in) :: flag_xyfill ! non-applicable xy points flagged with fillvalue
132 : real(r8), intent(in) :: field(idim,*) ! real*8 array
133 : type (dim_index_2d), intent(in ) :: dimind ! 2-D dimension index
134 : real(r8), intent(in) :: fillvalue
135 : !
136 : ! Local indices
137 : !
138 : integer :: ieu, jeu ! number of elements in each dimension
139 : integer :: i,k ! indices
140 : real(r8) :: tmp
141 :
142 0 : call dimind%dim_sizes(ieu, jeu)
143 :
144 0 : if (flag_xyfill) then
145 :
146 : !
147 : ! Ensure input field has fillvalue defined invariant in the z-direction, then increment nacs
148 : !
149 0 : call check_accum (field, idim, ieu, jeu, fillvalue)
150 0 : do i=1,ieu
151 0 : if (field(i,1) /= fillvalue) then
152 0 : nacs(i) = nacs(i) + 1
153 : end if
154 : end do
155 :
156 0 : do k=1,jeu
157 0 : do i=1,ieu
158 0 : if (field(i,k) /= fillvalue) then
159 0 : if (nacs(i)==1) then
160 0 : hbuf(i,k) = field(i,k)
161 0 : sbuf(i,k) = 0._r8
162 : else
163 0 : tmp = hbuf(i,k)
164 0 : hbuf(i,k) = hbuf(i,k) + (field(i,k)-hbuf(i,k))/nacs(i)
165 0 : sbuf(i,k) = sbuf(i,k) + (field(i,k)-hbuf(i,k))*(field(i,k)-tmp)
166 : endif
167 : end if
168 : end do
169 : end do
170 :
171 : else
172 :
173 : ! increment counter before variance calculation
174 0 : nacs(1) = nacs(1) + 1
175 :
176 0 : do k=1,jeu
177 0 : do i=1,ieu
178 0 : if (nacs(1)==1) then
179 0 : hbuf(i,k) = field(i,k)
180 0 : sbuf(i,k) = 0._r8
181 : else
182 0 : tmp = hbuf(i,k)
183 0 : hbuf(i,k) = hbuf(i,k) + (field(i,k)-hbuf(i,k))/nacs(1)
184 0 : sbuf(i,k) = sbuf(i,k) + (field(i,k)-hbuf(i,k))*(field(i,k)-tmp)
185 : endif
186 : end do
187 : end do
188 :
189 : end if
190 :
191 0 : end subroutine hbuf_accum_variance
192 :
193 : !#######################################################################
194 :
195 0 : subroutine hbuf_accum_add00z (buf8, field, nacs, dimind, idim, flag_xyfill, fillvalue)
196 : !
197 : !-----------------------------------------------------------------------
198 : !
199 : ! Purpose: Add the values of field to 2-D hbuf, only if time is 00z.
200 : ! Increment accumulation counter by 1.
201 : !
202 : !-----------------------------------------------------------------------
203 : !
204 : use time_manager, only: get_curr_date
205 :
206 : real(r8), pointer :: buf8(:,:) ! 2-D history buffer
207 : integer, pointer :: nacs(:) ! accumulation counter
208 : integer, intent(in) :: idim ! Longitude dimension of field array
209 : logical, intent(in) :: flag_xyfill ! non-applicable xy points flagged with fillvalue
210 : real(r8), intent(in ) :: field(idim,*) ! real*8 array
211 : type (dim_index_2d), intent(in ) :: dimind ! 2-D dimension index
212 : real(r8), intent(in) :: fillvalue
213 : !
214 : ! Local indices
215 : !
216 : integer :: ieu, jeu ! number of elements in each dimension
217 : integer :: i,k ! indices
218 : integer :: yr, mon, day, tod
219 :
220 : ! get the time of day, return if not 00z
221 0 : call get_curr_date (yr,mon,day,tod)
222 0 : if (tod /= 0) return
223 :
224 0 : call dimind%dim_sizes(ieu, jeu)
225 :
226 0 : if (flag_xyfill) then
227 0 : do k=1,jeu
228 0 : do i=1,ieu
229 0 : if (field(i,k) /= fillvalue) then
230 0 : buf8(i,k) = buf8(i,k) + field(i,k)
231 : end if
232 : end do
233 : end do
234 : !
235 : ! Ensure input field has fillvalue defined invariant in the z-direction, then increment nacs
236 : !
237 0 : call check_accum (field, idim, ieu, jeu, fillvalue)
238 0 : do i=1,ieu
239 0 : if (field(i,1) /= fillvalue) then
240 0 : nacs(i) = nacs(i) + 1
241 : end if
242 : end do
243 : else
244 0 : do k=1,jeu
245 0 : do i=1,ieu
246 0 : buf8(i,k) = buf8(i,k) + field(i,k)
247 : end do
248 : end do
249 0 : nacs(1) = nacs(1) + 1
250 : end if
251 :
252 : return
253 0 : end subroutine hbuf_accum_add00z
254 :
255 : !#######################################################################
256 :
257 1489176 : subroutine hbuf_accum_max (buf8, field, nacs, dimind, idim, flag_xyfill, fillvalue)
258 : !
259 : !-----------------------------------------------------------------------
260 : !
261 : ! Purpose: Accumulate the maximum values of field in 2-D hbuf
262 : ! Set accumulation counter to 1.
263 : !
264 : !-----------------------------------------------------------------------
265 : !
266 : real(r8), pointer :: buf8(:,:) ! 2-D history buffer
267 : integer, pointer :: nacs(:) ! accumulation counter
268 : integer, intent(in) :: idim ! Longitude dimension of field array
269 : logical, intent(in) :: flag_xyfill ! non-applicable xy points flagged with fillvalue
270 : real(r8), intent(in ) :: field(idim,*) ! real*8 array
271 : type (dim_index_2d), intent(in ) :: dimind ! 2-D dimension index
272 : real(r8), intent(in) :: fillvalue
273 : !
274 : ! Local indices
275 : !
276 : integer :: ieu, jeu ! number of elements in each dimension
277 : integer :: i, k
278 :
279 1489176 : call dimind%dim_sizes(ieu, jeu)
280 :
281 1489176 : if (flag_xyfill) then
282 0 : do k=1,jeu
283 0 : do i=1,ieu
284 0 : if (nacs(i) == 0) then
285 0 : buf8(i,k) = -huge (buf8)
286 : end if
287 0 : if (field(i,k) > buf8(i,k) .and. field(i,k) /= fillvalue) then
288 0 : buf8(i,k) = field(i,k)
289 : end if
290 : end do
291 : end do
292 : else
293 2978352 : do k=1,jeu
294 26354952 : do i=1,ieu
295 23376600 : if (nacs(1) == 0) then
296 7824600 : buf8(i,k) = -huge (buf8)
297 : end if
298 24865776 : if (field(i,k) > buf8(i,k)) then
299 7824600 : buf8(i,k) = field(i,k)
300 : end if
301 : end do
302 : end do
303 : end if
304 :
305 1489176 : if (flag_xyfill) then
306 0 : call check_accum (field, idim, ieu, jeu,fillvalue)
307 0 : do i=1,ieu
308 0 : if (field(i,1) /= fillvalue) then
309 0 : nacs(i) = 1
310 : end if
311 : end do
312 : else
313 1489176 : nacs(1) = 1
314 : end if
315 :
316 1489176 : return
317 0 : end subroutine hbuf_accum_max
318 :
319 : !#######################################################################
320 :
321 1489176 : subroutine hbuf_accum_min (buf8, field, nacs, dimind, idim, flag_xyfill, fillvalue)
322 : !
323 : !-----------------------------------------------------------------------
324 : !
325 : ! Purpose: Accumulate the minimum values of field in 2-D hbuf
326 : ! Set accumulation counter to 1.
327 : !
328 : !-----------------------------------------------------------------------
329 : !
330 : real(r8), pointer :: buf8(:,:) ! 2-D history buffer
331 : integer, pointer :: nacs(:) ! accumulation counter
332 : integer, intent(in) :: idim ! Longitude dimension of field array
333 : logical, intent(in) :: flag_xyfill ! non-applicable xy points flagged with fillvalue
334 : real(r8), intent(in ) :: field(idim,*) ! real*8 array
335 : type (dim_index_2d), intent(in ) :: dimind ! 2-D dimension index
336 : real(r8), intent(in) :: fillvalue
337 : !
338 : ! Local indices
339 : !
340 : integer :: ieu, jeu ! number of elements in each dimension
341 : integer :: i, k
342 :
343 1489176 : call dimind%dim_sizes(ieu, jeu)
344 :
345 1489176 : if (flag_xyfill) then
346 0 : do k=1,jeu
347 0 : do i=1,ieu
348 0 : if (nacs(i) == 0) then
349 0 : buf8(i,k) = +huge (buf8)
350 : end if
351 0 : if (field(i,k) < buf8(i,k) .and. field(i,k) /= fillvalue) then
352 0 : buf8(i,k) = field(i,k)
353 : end if
354 : end do
355 : end do
356 : else
357 2978352 : do k=1,jeu
358 26354952 : do i=1,ieu
359 23376600 : if (nacs(1) == 0) then
360 7824600 : buf8(i,k) = +huge (buf8)
361 : end if
362 24865776 : if (field(i,k) < buf8(i,k)) then
363 7824600 : buf8(i,k) = field(i,k)
364 : end if
365 : end do
366 : end do
367 : end if
368 :
369 1489176 : if (flag_xyfill) then
370 0 : call check_accum (field, idim, ieu, jeu, fillvalue)
371 0 : do i=1,ieu
372 0 : if (field(i,1) /= fillvalue) then
373 0 : nacs(i) = 1
374 : end if
375 : end do
376 : else
377 1489176 : nacs(1) = 1
378 : end if
379 :
380 1489176 : return
381 : end subroutine hbuf_accum_min
382 :
383 0 : subroutine hbuf_accum_addlcltime (buf8, field, nacs, dimind, idim, flag_xyfill, fillvalue, c , decomp_type,&
384 : lcltod_start, lcltod_stop)
385 : !
386 : !-----------------------------------------------------------------------
387 : !
388 : ! Purpose: Add the values of field to 2-D hbuf, only if the local time
389 : ! is in the range specified.
390 : ! Increment accumulation counter by 1.
391 : !
392 : !-----------------------------------------------------------------------
393 : !
394 : use time_manager, only: get_curr_date
395 : use phys_grid, only: get_rlon_all_p, phys_decomp
396 : use physconst, only: pi
397 : use phys_grid, only: get_ncols_p
398 : use dyn_grid, only: dyn_grid_get_elem_coords
399 :
400 : type (dim_index_2d), intent(in ) :: dimind ! 2-D dimension index
401 : real(r8), pointer :: buf8(:,:) ! 2-D history buffer
402 : integer, pointer :: nacs(:) ! accumulation counter
403 : integer, intent(in) :: idim ! Longitude dimension of field array
404 : logical, intent(in) :: flag_xyfill ! non-applicable xy points flagged with fillvalue
405 : real(r8), intent(in) :: field(idim,*) ! real*8 array
406 : integer, intent(in) :: c ! chunk (physics) or latitude (dynamics) index
407 :
408 : integer, intent(in) :: decomp_type, lcltod_start, lcltod_stop
409 : real(r8), intent(in) :: fillvalue
410 :
411 : !
412 : ! Local indices
413 : !
414 : integer :: ieu, jeu ! number of elements in each dimension
415 : integer :: i,k ! indices
416 : integer :: yr, mon, day, tod
417 : integer :: ncols
418 :
419 0 : integer, allocatable :: lcltod(:) ! local time of day (secs)
420 0 : logical, allocatable :: inavg(:) ! is the column in the desired local time range?
421 0 : real(r8), allocatable :: rlon(:) ! column longitude (radians)
422 0 : integer, allocatable :: cdex(:) ! global column index
423 :
424 0 : call dimind%dim_sizes(ieu, jeu)
425 :
426 0 : allocate( inavg(1:ieu) , lcltod(1:ieu) )
427 0 : lcltod(:) = 0
428 :
429 : !
430 : ! Get the time of day and longitude and compute the local time.
431 : !
432 0 : call get_curr_date (yr,mon,day,tod)
433 :
434 0 : if ( decomp_type == phys_decomp ) then
435 :
436 0 : ncols = get_ncols_p(c)
437 0 : ieu = ncols
438 0 : allocate( rlon(ncols) )
439 0 : call get_rlon_all_p(c, ncols, rlon)
440 0 : lcltod(1:ncols) = mod((tod) + (nint(86400._r8 * rlon(1:ncols) / 2._r8 / pi)), 86400)
441 :
442 : else
443 :
444 0 : ncols = ieu
445 0 : allocate(rlon(ncols),cdex(ncols))
446 0 : call dyn_grid_get_elem_coords( c, rlon=rlon, cdex=cdex )
447 0 : lcltod(:) = -999999
448 0 : where( cdex(:)>0 ) lcltod(1:ieu) = mod((tod) + (nint(86400._r8 * rlon(1:ncols) / 2._r8 / pi)), 86400)
449 :
450 : endif
451 :
452 :
453 :
454 : !
455 : ! Set a flag to indicate that the column is in the requested local time range.
456 : ! If lcltod_stop is less than lcltod_start, then the time is wrapping around 24 hours.
457 : !
458 0 : inavg(:) = .false.
459 :
460 0 : if (lcltod_stop < lcltod_start) then
461 : ! the ".and.(lcltod(:)>0" condition was added to exclude the undefined (-999999) columns
462 0 : where((lcltod(:) >= lcltod_start) .or. ((lcltod(:) < lcltod_stop).and.(lcltod(:)>0))) inavg(:) = .true.
463 0 : else if (lcltod_stop == lcltod_start) then
464 0 : where(lcltod(1:ieu) == lcltod_start) inavg(1:ieu) = .true.
465 : else
466 0 : where((lcltod(:) >= lcltod_start) .and. (lcltod(:) < lcltod_stop)) inavg(:) = .true.
467 : end if
468 :
469 0 : if (flag_xyfill) then
470 0 : do k=1,jeu
471 0 : do i=1,ieu
472 0 : if (inavg(i) .and. (field(i,k) /= fillvalue)) then
473 0 : buf8(i,k) = buf8(i,k) + field(i,k)
474 : end if
475 : end do
476 : end do
477 : !
478 : ! Ensure input field has fillvalue defined invariant in the z-direction, then increment nacs
479 : !
480 0 : call check_accum (field, idim, ieu, jeu, fillvalue)
481 0 : do i=1,ieu
482 0 : if (inavg(i) .and. (field(i,1) /= fillvalue)) then
483 0 : nacs(i) = nacs(i) + 1
484 : end if
485 : end do
486 : else
487 0 : do k=1,jeu
488 0 : do i=1,ieu
489 0 : if (inavg(i)) then
490 0 : buf8(i,k) = buf8(i,k) + field(i,k)
491 : end if
492 : end do
493 : end do
494 : !
495 : ! NOTE: Because of the local time, some columns in the chunk may not be used in the
496 : ! average, so nacs need to be dimension with the number of columns unlike the other
497 : ! averaging techniques.
498 : !
499 0 : do i=1,ieu
500 0 : if (inavg(i)) nacs(i) = nacs(i) + 1
501 : end do
502 : end if
503 :
504 0 : deallocate( inavg , rlon, lcltod )
505 0 : if (allocated(cdex)) deallocate(cdex)
506 :
507 0 : return
508 0 : end subroutine hbuf_accum_addlcltime
509 :
510 : !#######################################################################
511 :
512 :
513 : !#######################################################################
514 :
515 749232 : subroutine check_accum (field, idim, ieu, jeu, fillvalue)
516 : !
517 : integer, intent(in) :: idim
518 : real(r8), intent(in) :: field(idim,*) ! real*8 array
519 : integer, intent(in) :: ieu,jeu ! loop ranges
520 :
521 : logical :: bad
522 : integer :: i,k
523 : real(r8), intent(in) :: fillvalue
524 : !
525 : ! For multilevel fields ensure that all levels have fillvalue applied consistently
526 : !
527 749232 : bad = .false.
528 749232 : do k=2,jeu
529 749232 : do i=1,ieu
530 0 : if (field(i,1) == fillvalue .and. field(i,k) /= fillvalue .or. &
531 0 : field(i,1) /= fillvalue .and. field(i,k) == fillvalue) then
532 0 : bad = .true.
533 : end if
534 : end do
535 : end do
536 :
537 749232 : if (bad) then
538 0 : call endrun ('CHECK_ACCUM: inconsistent level values')
539 : end if
540 :
541 749232 : return
542 0 : end subroutine check_accum
543 :
544 : end module cam_history_buffers
|
