libsim  Versione 7.2.6
volgrid6d_class_compute.F90
1 ! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
2 ! authors:
3 ! Davide Cesari <dcesari@arpa.emr.it>
4 ! Paolo Patruno <ppatruno@arpa.emr.it>
5 
6 ! This program is free software; you can redistribute it and/or
7 ! modify it under the terms of the GNU General Public License as
8 ! published by the Free Software Foundation; either version 2 of
9 ! the License, or (at your option) any later version.
10 
11 ! This program is distributed in the hope that it will be useful,
12 ! but WITHOUT ANY WARRANTY; without even the implied warranty of
13 ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 ! GNU General Public License for more details.
15 
16 ! You should have received a copy of the GNU General Public License
17 ! along with this program. If not, see <http://www.gnu.org/licenses/>.
18 #include "config.h"
19 
30 USE grid_id_class
32 USE simple_stat
33 IMPLICIT NONE
34 
35 CONTAINS
36 
102 SUBROUTINE volgrid6d_compute_stat_proc(this, that, stat_proc_input, stat_proc, &
103  step, start, full_steps, frac_valid, max_step, weighted, clone)
104 TYPE(volgrid6d),INTENT(inout) :: this
105 TYPE(volgrid6d),INTENT(out) :: that
106 INTEGER,INTENT(in) :: stat_proc_input
107 INTEGER,INTENT(in) :: stat_proc
108 TYPE(timedelta),INTENT(in) :: step
109 TYPE(datetime),INTENT(in),OPTIONAL :: start
110 LOGICAL,INTENT(in),OPTIONAL :: full_steps
111 REAL,INTENT(in),OPTIONAL :: frac_valid
112 TYPE(timedelta),INTENT(in),OPTIONAL :: max_step ! maximum allowed distance in time between two single valid data within a dataset, for the dataset to be eligible for statistical processing
113 LOGICAL,INTENT(in),OPTIONAL :: weighted
114 LOGICAL , INTENT(in),OPTIONAL :: clone
115 
116 INTEGER :: dtmax, dtstep
117 
118 
119 IF (stat_proc_input == 254) THEN
120  CALL l4f_category_log(this%category, l4f_info, &
121  'computing statistical processing by aggregation '//&
122  trim(to_char(stat_proc_input))//':'//trim(to_char(stat_proc)))
123 
124  CALL volgrid6d_compute_stat_proc_agg(this, that, stat_proc, &
125  step, start, full_steps, max_step, clone)
126 
127 ELSE IF (stat_proc == 254) THEN
128  CALL l4f_category_log(this%category, l4f_error, &
129  'statistical processing to instantaneous data not implemented for gridded fields')
130  CALL raise_error()
131 
132 ELSE IF (stat_proc_input == stat_proc .OR. &
133  (stat_proc == 0 .OR. stat_proc == 2 .OR. stat_proc == 3)) THEN
134 ! avg, min and max can be computed from any input, with care
135 
136  IF (count(this%timerange(:)%timerange == stat_proc_input) == 0) THEN
137  CALL l4f_category_log(this%category, l4f_warn, &
138  'no timeranges of the desired statistical processing type '//t2c(stat_proc)//' available')
139 ! return an empty volume, without signaling error
140  CALL init(that)
141  CALL volgrid6d_alloc_vol(that)
142 
143  ELSE
144 ! euristically determine whether aggregation or difference is more suitable
145  dtmax = maxval(this%timerange(:)%p2, &
146  mask=(this%timerange(:)%timerange == stat_proc))
147  CALL getval(step, asec=dtstep)
148 
149 #ifdef DEBUG
150  CALL l4f_category_log(this%category, l4f_debug, &
151  'stat_proc='//t2c(stat_proc)//' dtmax='//t2c(dtmax)//' dtstep='//t2c(dtstep))
152 #endif
153 
154  IF (dtstep <= dtmax) THEN
155  CALL l4f_category_log(this%category, l4f_info, &
156  'recomputing statistically processed data by difference '// &
157  t2c(stat_proc_input)//':'//t2c(stat_proc))
158  CALL volgrid6d_recompute_stat_proc_diff(this, that, stat_proc, step, &
159  full_steps, start, clone)
160  ELSE
161  CALL l4f_category_log(this%category, l4f_info, &
162  'recomputing statistically processed data by aggregation '// &
163  t2c(stat_proc_input)//':'//t2c(stat_proc))
164  CALL volgrid6d_recompute_stat_proc_agg(this, that, stat_proc, step, start, &
165  full_steps, frac_valid, clone, stat_proc_input)
166  ENDIF
167  ENDIF
168 
169 ELSE ! IF (stat_proc_input /= stat_proc) THEN
170  IF ((stat_proc_input == 0 .AND. stat_proc == 1) .OR. &
171  (stat_proc_input == 1 .AND. stat_proc == 0)) THEN
172  CALL l4f_category_log(this%category, l4f_info, &
173  'computing statistically processed data by integration/differentiation '// &
174  t2c(stat_proc_input)//':'//t2c(stat_proc))
175  CALL volgrid6d_compute_stat_proc_metamorph(this, that, stat_proc_input, &
176  stat_proc, clone)
177  ELSE
178  CALL l4f_category_log(this%category, l4f_error, &
179  'statistical processing '//t2c(stat_proc_input)//':'//t2c(stat_proc)// &
180  ' not implemented or does not make sense')
181  CALL raise_error()
182  ENDIF
183 
184 ENDIF
185 
186 END SUBROUTINE volgrid6d_compute_stat_proc
187 
188 
231 SUBROUTINE volgrid6d_recompute_stat_proc_agg(this, that, stat_proc, &
232  step, start, full_steps, frac_valid, clone, stat_proc_input)
233 TYPE(volgrid6d),INTENT(inout) :: this
234 TYPE(volgrid6d),INTENT(out) :: that
235 INTEGER,INTENT(in) :: stat_proc
236 TYPE(timedelta),INTENT(in) :: step
237 TYPE(datetime),INTENT(in),OPTIONAL :: start
238 LOGICAL,INTENT(in),OPTIONAL :: full_steps
239 REAL,INTENT(in),OPTIONAL :: frac_valid
240 LOGICAL, INTENT(in),OPTIONAL :: clone
241 INTEGER,INTENT(in),OPTIONAL :: stat_proc_input
242 
243 INTEGER :: tri
244 INTEGER i, j, n, n1, ndtr, i3, i6
245 TYPE(arrayof_ttr_mapper),POINTER :: map_ttr(:,:)
246 INTEGER,POINTER :: dtratio(:)
247 REAL :: lfrac_valid
248 LOGICAL :: lclone
249 REAL,POINTER :: voldatiin(:,:), voldatiout(:,:)
250 
251 
252 NULLIFY(voldatiin, voldatiout)
253 IF (PRESENT(stat_proc_input)) THEN
254  tri = stat_proc_input
255 ELSE
256  tri = stat_proc
257 ENDIF
258 IF (PRESENT(frac_valid)) THEN
259  lfrac_valid = frac_valid
260 ELSE
261  lfrac_valid = 1.0
262 ENDIF
263 
264 CALL init(that)
265 ! be safe
266 CALL volgrid6d_alloc_vol(this)
267 
268 ! when volume is not decoded it is better to clone anyway to avoid
269 ! overwriting fields
270 lclone = optio_log(clone) .OR. .NOT.ASSOCIATED(this%voldati)
271 ! initialise the output volume
272 CALL init(that, griddim=this%griddim, time_definition=this%time_definition)
273 CALL volgrid6d_alloc(that, dim=this%griddim%dim, ntimerange=1, &
274  nlevel=SIZE(this%level), nvar=SIZE(this%var), ini=.false.)
275 that%level = this%level
276 that%var = this%var
277 
278 CALL recompute_stat_proc_agg_common(this%time, this%timerange, stat_proc, tri, &
279  step, this%time_definition, that%time, that%timerange, map_ttr, &
280  dtratio=dtratio, start=start, full_steps=full_steps)
281 
282 CALL volgrid6d_alloc_vol(that, decode=ASSOCIATED(this%voldati))
283 
284 do_otimerange: DO j = 1, SIZE(that%timerange)
285  do_otime: DO i = 1, SIZE(that%time)
286 
287  DO n1 = 1, SIZE(dtratio)
288  IF (dtratio(n1) <= 0) cycle ! safety check
289 
290  DO i6 = 1, SIZE(this%var)
291  DO i3 = 1, SIZE(this%level)
292  CALL volgrid_get_vol_2d(that, i3, i, j, i6, voldatiout)
293  ndtr = 0
294  DO n = 1, map_ttr(i,j)%arraysize
295  IF (map_ttr(i,j)%array(n)%extra_info == dtratio(n1)) THEN
296  ndtr = ndtr + 1
297  CALL volgrid_get_vol_2d(this, i3, map_ttr(i,j)%array(n)%it, &
298  map_ttr(i,j)%array(n)%itr, i6, voldatiin)
299 
300  IF (ndtr == 1) THEN
301  voldatiout = voldatiin
302  IF (lclone) THEN
303  CALL copy(this%gaid(i3, map_ttr(i,j)%array(n)%it,&
304  map_ttr(i,j)%array(n)%itr,i6), that%gaid(i3,i,j,i6))
305  ELSE
306  that%gaid(i3,i,j,i6) = this%gaid(i3, map_ttr(i,j)%array(n)%it, &
307  map_ttr(i,j)%array(n)%itr,i6)
308  ENDIF
309 
310  ELSE ! second or more time
311  SELECT CASE(stat_proc)
312  CASE (0, 200, 1, 4) ! average, vectorial mean, accumulation, difference
313  WHERE(c_e(voldatiin(:,:)) .AND. c_e(voldatiout(:,:)))
314  voldatiout(:,:) = voldatiout(:,:) + voldatiin(:,:)
315  ELSEWHERE
316  voldatiout(:,:) = rmiss
317  END WHERE
318  CASE(2) ! maximum
319  WHERE(c_e(voldatiin(:,:)) .AND. c_e(voldatiout(:,:)))
320  voldatiout(:,:) = max(voldatiout(:,:), voldatiin(:,:))
321  ELSEWHERE
322  voldatiout(:,:) = rmiss
323  END WHERE
324  CASE(3) ! minimum
325  WHERE(c_e(voldatiin(:,:)) .AND. c_e(voldatiout(:,:)))
326  voldatiout(:,:) = min(voldatiout(:,:), voldatiin(:,:))
327  ELSEWHERE
328  voldatiout(:,:) = rmiss
329  END WHERE
330  END SELECT
331 
332  ENDIF ! first time
333  ENDIF ! dtratio(n1)
334  ENDDO ! ttr
335 
336 #ifdef DEBUG
337  CALL l4f_log(l4f_debug, &
338  'compute_stat_proc_agg, ndtr/dtratio/frac_valid: '// &
339  t2c(ndtr)//'/'//t2c(dtratio(n1))//'/'//t2c(lfrac_valid))
340 #endif
341  IF (ndtr > 0) THEN ! why this condition was not here before?
342  IF (real(ndtr)/real(dtratio(n1)) >= lfrac_valid) THEN ! success
343  IF (stat_proc == 0) THEN ! average
344  WHERE(c_e(voldatiout(:,:)))
345  voldatiout(:,:) = voldatiout(:,:)/ndtr
346  END WHERE
347  ENDIF
348  CALL volgrid_set_vol_2d(that, i3, i, j, i6, voldatiout)
349 #ifdef DEBUG
350  CALL l4f_log(l4f_debug, &
351  'compute_stat_proc_agg, coding lev/t/tr/var: '// &
352  t2c(i3)//'/'//t2c(i)//'/'//t2c(j)//'/'//t2c(i6))
353 #endif
354  ELSE
355 ! must nullify the output gaid here, otherwise an incomplete field will be output
356  IF (lclone) THEN
357  CALL delete(that%gaid(i3,i,j,i6))
358  ELSE
359  CALL init(that%gaid(i3,i,j,i6)) ! grid_id lacks a nullify method
360  ENDIF
361 #ifdef DEBUG
362  CALL l4f_log(l4f_debug, &
363  'compute_stat_proc_agg, skipping lev/t/tr/var: '// &
364  t2c(i3)//'/'//t2c(i)//'/'//t2c(j)//'/'//t2c(i6))
365 #endif
366  ENDIF
367  ENDIF ! ndtr > 0
368 
369  ENDDO ! level
370  ENDDO ! var
371  ENDDO ! dtratio
372  CALL delete(map_ttr(i,j))
373  ENDDO do_otime
374 ENDDO do_otimerange
375 
376 DEALLOCATE(dtratio, map_ttr)
377 
378 END SUBROUTINE volgrid6d_recompute_stat_proc_agg
379 
380 
404 SUBROUTINE volgrid6d_compute_stat_proc_agg(this, that, stat_proc, &
405  step, start, full_steps, max_step, clone)
406 TYPE(volgrid6d),INTENT(inout) :: this
407 TYPE(volgrid6d),INTENT(out) :: that
408 INTEGER,INTENT(in) :: stat_proc
409 TYPE(timedelta),INTENT(in) :: step
410 TYPE(datetime),INTENT(in),OPTIONAL :: start
411 LOGICAL,INTENT(in),OPTIONAL :: full_steps
412 TYPE(timedelta),INTENT(in),OPTIONAL :: max_step
413 LOGICAL , INTENT(in),OPTIONAL :: clone
414 
415 INTEGER :: tri
416 INTEGER i, j, n, ninp, i3, i6
417 TYPE(arrayof_ttr_mapper),POINTER :: map_ttr(:,:)
418 TYPE(timedelta) :: lmax_step
419 LOGICAL :: lclone
420 REAL,POINTER :: voldatiin(:,:), voldatiout(:,:)
421 
422 
423 NULLIFY(voldatiin, voldatiout)
424 tri = 254
425 IF (PRESENT(max_step)) THEN
426  lmax_step = max_step
427 ELSE
428  lmax_step = timedelta_max
429 ENDIF
430 
431 CALL init(that)
432 ! be safe
433 CALL volgrid6d_alloc_vol(this)
434 
435 ! when volume is not decoded it is better to clone anyway to avoid
436 ! overwriting fields
437 lclone = optio_log(clone) .OR. .NOT.ASSOCIATED(this%voldati)
438 ! initialise the output volume
439 CALL init(that, griddim=this%griddim, time_definition=this%time_definition)
440 CALL volgrid6d_alloc(that, dim=this%griddim%dim, ntimerange=1, &
441  nlevel=SIZE(this%level), nvar=SIZE(this%var), ini=.false.)
442 that%level = this%level
443 that%var = this%var
444 
445 CALL recompute_stat_proc_agg_common(this%time, this%timerange, stat_proc, tri, &
446  step, this%time_definition, that%time, that%timerange, map_ttr, &
447  start=start, full_steps=full_steps)
448 
449 CALL volgrid6d_alloc_vol(that, decode=ASSOCIATED(this%voldati))
450 
451 do_otimerange: DO j = 1, SIZE(that%timerange)
452  do_otime: DO i = 1, SIZE(that%time)
453  ninp = map_ttr(i,j)%arraysize
454  IF (ninp <= 0) cycle do_otime
455 
456  IF (stat_proc == 4) THEN ! check validity for difference
457  IF (map_ttr(i,j)%array(1)%extra_info /= 1 .OR. &
458  map_ttr(i,j)%array(ninp)%extra_info /= 2) THEN
459  CALL delete(map_ttr(i,j))
460  cycle do_otime
461  ENDIF
462  ELSE
463 ! check validity condition (missing values in volume are not accounted for)
464  DO n = 2, ninp
465  IF (map_ttr(i,j)%array(n)%time - map_ttr(i,j)%array(n-1)%time > &
466  lmax_step) THEN
467  CALL delete(map_ttr(i,j))
468  cycle do_otime
469  ENDIF
470  ENDDO
471  ENDIF
472 
473  DO i6 = 1, SIZE(this%var)
474  DO i3 = 1, SIZE(this%level)
475  CALL volgrid_get_vol_2d(that, i3, i, j, i6, voldatiout)
476 
477  IF (stat_proc == 4) THEN ! special treatment for difference
478  IF (lclone) THEN
479  CALL copy(this%gaid(i3, map_ttr(i,j)%array(1)%it,&
480  map_ttr(i,j)%array(1)%itr,i6), that%gaid(i3,i,j,i6))
481  ELSE
482  that%gaid(i3,i,j,i6) = this%gaid(i3, map_ttr(i,j)%array(1)%it, &
483  map_ttr(i,j)%array(1)%itr,i6)
484  ENDIF
485 ! improve the next workflow?
486  CALL volgrid_get_vol_2d(this, i3, map_ttr(i,j)%array(ninp)%it, &
487  map_ttr(i,j)%array(ninp)%itr, i6, voldatiin)
488  voldatiout = voldatiin
489  CALL volgrid_get_vol_2d(this, i3, map_ttr(i,j)%array(1)%it, &
490  map_ttr(i,j)%array(1)%itr, i6, voldatiin)
491 
492  WHERE(c_e(voldatiin(:,:)) .AND. c_e(voldatiout(:,:)))
493  voldatiout(:,:) = voldatiout(:,:) - voldatiin(:,:)
494  ELSEWHERE
495  voldatiout(:,:) = rmiss
496  END WHERE
497 
498  ELSE ! other stat_proc
499  DO n = 1, ninp
500  CALL volgrid_get_vol_2d(this, i3, map_ttr(i,j)%array(n)%it, &
501  map_ttr(i,j)%array(n)%itr, i6, voldatiin)
502 
503  IF (n == 1) THEN
504  voldatiout = voldatiin
505  IF (lclone) THEN
506  CALL copy(this%gaid(i3, map_ttr(i,j)%array(n)%it,&
507  map_ttr(i,j)%array(n)%itr,i6), that%gaid(i3,i,j,i6))
508  ELSE
509  that%gaid(i3,i,j,i6) = this%gaid(i3, map_ttr(i,j)%array(n)%it, &
510  map_ttr(i,j)%array(n)%itr,i6)
511  ENDIF
512 
513  ELSE ! second or more time
514  SELECT CASE(stat_proc)
515  CASE (0, 1) ! average, accumulation
516  WHERE(c_e(voldatiin(:,:)) .AND. c_e(voldatiout(:,:)))
517  voldatiout(:,:) = voldatiout(:,:) + voldatiin(:,:)
518  ELSEWHERE
519  voldatiout(:,:) = rmiss
520  END WHERE
521  CASE(2) ! maximum
522  WHERE(c_e(voldatiin(:,:)) .AND. c_e(voldatiout(:,:)))
523  voldatiout(:,:) = max(voldatiout(:,:), voldatiin(:,:))
524  ELSEWHERE
525  voldatiout(:,:) = rmiss
526  END WHERE
527  CASE(3) ! minimum
528  WHERE(c_e(voldatiin(:,:)) .AND. c_e(voldatiout(:,:)))
529  voldatiout(:,:) = min(voldatiout(:,:), voldatiin(:,:))
530  ELSEWHERE
531  voldatiout(:,:) = rmiss
532  END WHERE
533  END SELECT
534 
535  ENDIF ! first time
536  ENDDO
537  IF (stat_proc == 0) THEN ! average
538  WHERE(c_e(voldatiout(:,:)))
539  voldatiout(:,:) = voldatiout(:,:)/ninp
540  END WHERE
541  ENDIF
542  ENDIF
543  CALL volgrid_set_vol_2d(that, i3, i, j, i6, voldatiout)
544  ENDDO ! level
545  ENDDO ! var
546  CALL delete(map_ttr(i,j))
547  ENDDO do_otime
548 ENDDO do_otimerange
549 
550 DEALLOCATE(map_ttr)
551 
552 
553 END SUBROUTINE volgrid6d_compute_stat_proc_agg
554 
555 
580 SUBROUTINE volgrid6d_recompute_stat_proc_diff(this, that, stat_proc, step, full_steps, start, clone)
581 TYPE(volgrid6d),INTENT(inout) :: this
582 TYPE(volgrid6d),INTENT(out) :: that
583 INTEGER,INTENT(in) :: stat_proc
584 TYPE(timedelta),INTENT(in) :: step
585 LOGICAL,INTENT(in),OPTIONAL :: full_steps
586 TYPE(datetime),INTENT(in),OPTIONAL :: start
587 LOGICAL,INTENT(in),OPTIONAL :: clone
588 INTEGER :: i3, i4, i6, i, j, k, l, nitr, steps
589 INTEGER,ALLOCATABLE :: map_tr(:,:,:,:,:), f(:), keep_tr(:,:,:)
590 REAL,POINTER :: voldatiin1(:,:), voldatiin2(:,:), voldatiout(:,:)
591 !LOGICAL,POINTER :: mask_timerange(:)
592 LOGICAL :: lclone
593 TYPE(vol7d_var),ALLOCATABLE :: varbufr(:)
594 
595 
596 ! be safe
597 CALL volgrid6d_alloc_vol(this)
598 ! when volume is not decoded it is better to clone anyway to avoid
599 ! overwriting fields
600 lclone = optio_log(clone) .OR. .NOT.ASSOCIATED(this%voldati)
601 ! initialise the output volume
602 CALL init(that, griddim=this%griddim, time_definition=this%time_definition)
603 CALL volgrid6d_alloc(that, dim=this%griddim%dim, &
604  nlevel=SIZE(this%level), nvar=SIZE(this%var), ini=.false.)
605 that%level = this%level
606 that%var = this%var
607 
608 ! compute length of cumulation step in seconds
609 CALL getval(step, asec=steps)
610 
611 ! compute the statistical processing relations, output time and
612 ! timerange are defined here
613 CALL recompute_stat_proc_diff_common(this%time, this%timerange, stat_proc, step, &
614  that%time, that%timerange, map_tr, f, keep_tr, &
615  this%time_definition, full_steps, start)
616 nitr = SIZE(f)
617 
618 ! complete the definition of the output volume
619 CALL volgrid6d_alloc_vol(that, decode=ASSOCIATED(this%voldati))
620 ! allocate workspace once
621 IF (.NOT.ASSOCIATED(that%voldati)) THEN
622  ALLOCATE(voldatiin1(this%griddim%dim%nx, this%griddim%dim%ny), &
623  voldatiin2(this%griddim%dim%nx, this%griddim%dim%ny), &
624  voldatiout(this%griddim%dim%nx, this%griddim%dim%ny))
625 ENDIF
626 
627 ! copy the timeranges already satisfying the requested step, if any
628 DO i4 = 1, SIZE(this%time)
629  DO i = 1, nitr
630  IF (c_e(keep_tr(i, i4, 2))) THEN
631  l = keep_tr(i, i4, 1)
632  k = keep_tr(i, i4, 2)
633 #ifdef DEBUG
634  CALL l4f_category_log(this%category, l4f_debug, &
635  'volgrid6d_recompute_stat_proc_diff, good timerange: '//t2c(f(i))// &
636  '->'//t2c(k))
637 #endif
638  DO i6 = 1, SIZE(this%var)
639  DO i3 = 1, SIZE(this%level)
640  IF (c_e(this%gaid(i3,i4,f(i),i6))) THEN
641  IF (lclone) THEN
642  CALL copy(this%gaid(i3,i4,f(i),i6), that%gaid(i3,l,k,i6))
643  ELSE
644  that%gaid(i3,l,k,i6) = this%gaid(i3,i4,f(i),i6)
645  ENDIF
646  IF (ASSOCIATED(that%voldati)) THEN
647  that%voldati(:,:,i3,l,k,i6) = this%voldati(:,:,i3,i4,f(i),i6)
648  ELSE
649  CALL volgrid_get_vol_2d(this, i3, i4, f(i), i6, voldatiout)
650  CALL volgrid_set_vol_2d(that, i3, l, k, i6, voldatiout)
651  ENDIF
652  ENDIF
653  ENDDO
654  ENDDO
655  ENDIF
656  ENDDO
657 ENDDO
658 
659 ! varbufr required for setting posdef, optimize with an array
660 ALLOCATE(varbufr(SIZE(this%var)))
661 DO i6 = 1, SIZE(this%var)
662  varbufr(i6) = convert(this%var(i6))
663 ENDDO
664 ! compute statistical processing
665 DO l = 1, SIZE(this%time)
666  DO k = 1, nitr
667  DO j = 1, SIZE(this%time)
668  DO i = 1, nitr
669  IF (c_e(map_tr(i,j,k,l,1))) THEN
670  DO i6 = 1, SIZE(this%var)
671  DO i3 = 1, SIZE(this%level)
672 
673  IF (c_e(this%gaid(i3,j,f(i),i6)) .AND. &
674  c_e(this%gaid(i3,l,f(k),i6))) THEN
675 ! take the gaid from the second time/timerange contributing to the
676 ! result (l,f(k))
677  IF (lclone) THEN
678  CALL copy(this%gaid(i3,l,f(k),i6), &
679  that%gaid(i3,map_tr(i,j,k,l,1),map_tr(i,j,k,l,2),i6))
680  ELSE
681  that%gaid(i3,map_tr(i,j,k,l,1),map_tr(i,j,k,l,2),i6) = &
682  this%gaid(i3,l,f(k),i6)
683  ENDIF
684 
685 ! get/set 2d sections API is used
686  CALL volgrid_get_vol_2d(this, i3, l, f(k), i6, voldatiin1)
687  CALL volgrid_get_vol_2d(this, i3, j, f(i), i6, voldatiin2)
688  IF (ASSOCIATED(that%voldati)) &
689  CALL volgrid_get_vol_2d(that, i3, &
690  map_tr(i,j,k,l,1), map_tr(i,j,k,l,2), i6, voldatiout)
691 
692  IF (stat_proc == 0) THEN ! average
693  WHERE(c_e(voldatiin1(:,:)) .AND. c_e(voldatiin2(:,:)))
694  voldatiout(:,:) = &
695  (voldatiin1(:,:)*this%timerange(f(k))%p2 - &
696  voldatiin2(:,:)*this%timerange(f(i))%p2)/ &
697  steps
698  ELSEWHERE
699  voldatiout(:,:) = rmiss
700  END WHERE
701  ELSE IF (stat_proc == 1 .OR. stat_proc == 4) THEN ! acc, diff
702  WHERE(c_e(voldatiin1(:,:)) .AND. c_e(voldatiin2(:,:)))
703  voldatiout(:,:) = voldatiin1(:,:) - voldatiin2(:,:)
704  ELSEWHERE
705  voldatiout(:,:) = rmiss
706  END WHERE
707  IF (stat_proc == 1) THEN
708  CALL vol7d_var_features_posdef_apply(varbufr(i6), voldatiout)
709  ENDIF
710  ENDIF
711 
712  CALL volgrid_set_vol_2d(that, i3, &
713  map_tr(i,j,k,l,1), map_tr(i,j,k,l,2), i6, voldatiout)
714 
715  ENDIF
716  ENDDO
717  ENDDO
718  ENDIF
719  ENDDO
720  ENDDO
721  ENDDO
722 ENDDO
723 
724 IF (.NOT.ASSOCIATED(that%voldati)) THEN
725  DEALLOCATE(voldatiin1, voldatiin2, voldatiout)
726 ENDIF
727 
728 END SUBROUTINE volgrid6d_recompute_stat_proc_diff
729 
730 
758 SUBROUTINE volgrid6d_compute_stat_proc_metamorph(this, that, stat_proc_input, stat_proc, clone)
759 TYPE(volgrid6d),INTENT(inout) :: this
760 TYPE(volgrid6d),INTENT(out) :: that
761 INTEGER,INTENT(in) :: stat_proc_input
762 INTEGER,INTENT(in) :: stat_proc
763 LOGICAL , INTENT(in),OPTIONAL :: clone
764 
765 INTEGER :: j, i3, i4, i6
766 INTEGER,POINTER :: map_tr(:)
767 REAL,POINTER :: voldatiin(:,:), voldatiout(:,:)
768 REAL,ALLOCATABLE :: int_ratio(:)
769 LOGICAL :: lclone
770 
771 NULLIFY(voldatiin, voldatiout)
772 
773 ! be safe
774 CALL volgrid6d_alloc_vol(this)
775 ! when volume is not decoded it is better to clone anyway to avoid
776 ! overwriting fields
777 lclone = optio_log(clone) .OR. .NOT.ASSOCIATED(this%voldati)
778 
779 IF (.NOT.((stat_proc_input == 0 .AND. stat_proc == 1) .OR. &
780  (stat_proc_input == 1 .AND. stat_proc == 0))) THEN
781 
782  CALL l4f_category_log(this%category, l4f_warn, &
783  'compute_stat_proc_metamorph, can only be applied to average->accumulated timerange and viceversa')
784 ! return an empty volume, without signaling error
785  CALL init(that)
786  CALL volgrid6d_alloc_vol(that)
787  RETURN
788 ENDIF
789 
790 ! initialise the output volume
791 CALL init(that, griddim=this%griddim, time_definition=this%time_definition)
792 CALL volgrid6d_alloc(that, dim=this%griddim%dim, ntime=SIZE(this%time), &
793  nlevel=SIZE(this%level), nvar=SIZE(this%var), ini=.false.)
794 that%time = this%time
795 that%level = this%level
796 that%var = this%var
797 
798 CALL compute_stat_proc_metamorph_common(stat_proc_input, this%timerange, stat_proc, &
799  that%timerange, map_tr)
800 
801 ! complete the definition of the output volume
802 CALL volgrid6d_alloc_vol(that, decode=ASSOCIATED(this%voldati))
803 
804 IF (stat_proc == 0) THEN ! average -> integral
805  int_ratio = 1./real(that%timerange(:)%p2)
806 ELSE ! cumulation
807  int_ratio = real(that%timerange(:)%p2)
808 ENDIF
809 
810 DO i6 = 1, SIZE(this%var)
811  DO j = 1, SIZE(map_tr)
812  DO i4 = 1, SIZE(that%time)
813  DO i3 = 1, SIZE(this%level)
814 
815  IF (lclone) THEN
816  CALL copy(this%gaid(i3,i4,map_tr(j),i6), that%gaid(i3,i4,j,i6))
817  ELSE
818  that%gaid(i3,i4,map_tr(j),i6) = this%gaid(i3,i4,j,i6)
819  ENDIF
820  CALL volgrid_get_vol_2d(this, i3, i4, map_tr(j), i6, voldatiin)
821  CALL volgrid_get_vol_2d(that, i3, i4, j, i6, voldatiout)
822  WHERE (c_e(voldatiin))
823  voldatiout = voldatiin*int_ratio(j)
824  ELSEWHERE
825  voldatiout = rmiss
826  END WHERE
827  CALL volgrid_set_vol_2d(that, i3, i4, j, i6, voldatiout)
828  ENDDO
829  ENDDO
830  ENDDO
831 ENDDO
832 
833 
834 END SUBROUTINE volgrid6d_compute_stat_proc_metamorph
835 
850 SUBROUTINE volgrid6d_compute_vert_coord_var(this, level, volgrid_lev)
851 TYPE(volgrid6d),INTENT(in) :: this
852 TYPE(vol7d_level),INTENT(in) :: level
853 TYPE(volgrid6d),INTENT(out) :: volgrid_lev
854 
855 INTEGER :: nlev, i, ii, iii, iiii
856 TYPE(grid_id) :: out_gaid
857 LOGICAL,ALLOCATABLE :: levmask(:)
858 TYPE(volgrid6d_var) :: lev_var
859 
860 CALL init(volgrid_lev) ! initialise to null
861 IF (.NOT.ASSOCIATED(this%gaid)) THEN
862  CALL l4f_log(l4f_error, 'volgrid6d_compute_vert_coord_var: input volume not allocated')
863  RETURN
864 ENDIF
865 ! if layer, both surfaces must be of the same type
866 IF (c_e(level%level2) .AND. level%level1 /= level%level2) THEN
867  CALL l4f_log(l4f_error, 'volgrid6d_compute_vert_coord_var: requested (mixed) layer type not valid')
868  RETURN
869 ENDIF
870 
871 ! look for valid levels to be converted to vars
872 ALLOCATE(levmask(SIZE(this%level)))
873 levmask = this%level%level1 == level%level1 .AND. &
874  this%level%level2 == level%level2 .AND. c_e(this%level%l1)
875 IF (c_e(level%level2)) levmask = levmask .AND. c_e(this%level%l2)
876 nlev = count(levmask)
877 IF (nlev == 0) THEN
878  CALL l4f_log(l4f_error, 'volgrid6d_compute_vert_coord_var: requested level type not available')
879  RETURN
880 ENDIF
881 
882 out_gaid = grid_id_new()
883 gaidloop: DO i=1 ,SIZE(this%gaid,1)
884  DO ii=1 ,SIZE(this%gaid,2)
885  DO iii=1 ,SIZE(this%gaid,3)
886  DO iiii=1 ,SIZE(this%gaid,4)
887  IF (c_e(this%gaid(i,ii,iii,iiii))) THEN ! conserve first valid gaid
888  CALL copy(this%gaid(i,ii,iii,iiii), out_gaid)
889  EXIT gaidloop
890  ENDIF
891  ENDDO
892  ENDDO
893  ENDDO
894 ENDDO gaidloop
895 
896 ! look for variable corresponding to level
897 lev_var = convert(vol7d_var_new(btable=vol7d_level_to_var(level)), &
898  grid_id_template=out_gaid)
899 IF (.NOT.c_e(lev_var)) THEN
900  CALL l4f_log(l4f_error, 'volgrid6d_compute_vert_coord_var: no variable corresponds to requested level type')
901  RETURN
902 ENDIF
903 
904 ! prepare output volume
905 CALL init(volgrid_lev, griddim=this%griddim, &
906  time_definition=this%time_definition) !, categoryappend=categoryappend)
907 CALL volgrid6d_alloc(volgrid_lev, ntime=SIZE(this%time), nlevel=nlev, &
908  ntimerange=SIZE(this%timerange), nvar=1)
909 ! fill metadata
910 volgrid_lev%time = this%time
911 volgrid_lev%level = pack(this%level, mask=levmask)
912 volgrid_lev%timerange = this%timerange
913 volgrid_lev%var(1) = lev_var
914 
915 CALL volgrid6d_alloc_vol(volgrid_lev, decode=.true.)
916 ! fill data
917 DO i = 1, nlev
918  IF (c_e(level%level2)) THEN
919  volgrid_lev%voldati(:,:,i,:,:,:) = real(volgrid_lev%level(i)%l1 + &
920  volgrid_lev%level(i)%l2)* &
921  vol7d_level_to_var_factor(volgrid_lev%level(i))/2.
922  ELSE
923  volgrid_lev%voldati(:,:,i,:,:,:) = real(volgrid_lev%level(i)%l1)* &
924  vol7d_level_to_var_factor(volgrid_lev%level(i))
925  ENDIF
926 ENDDO
927 ! fill gaid for subsequent export
928 IF (c_e(out_gaid)) THEN
929  DO i=1 ,SIZE(volgrid_lev%gaid,1)
930  DO ii=1 ,SIZE(volgrid_lev%gaid,2)
931  DO iii=1 ,SIZE(volgrid_lev%gaid,3)
932  DO iiii=1 ,SIZE(volgrid_lev%gaid,4)
933  CALL copy(out_gaid, volgrid_lev%gaid(i,ii,iii,iiii))
934  ENDDO
935  ENDDO
936  ENDDO
937  ENDDO
938  CALL delete(out_gaid)
939 ENDIF
940 
941 END SUBROUTINE volgrid6d_compute_vert_coord_var
942 
943 END MODULE volgrid6d_class_compute
Distruttori per le 2 classi.
Restituiscono il valore dell'oggetto nella forma desiderata.
Costruttori per le classi datetime e timedelta.
Functions that return a trimmed CHARACTER representation of the input variable.
Restituiscono il valore dell'oggetto in forma di stringa stampabile.
Make a deep copy, if possible, of the grid identifier.
Apply the conversion function this to values.
Classi per la gestione delle coordinate temporali.
This module defines an abstract interface to different drivers for access to files containing gridded...
Module for basic statistical computations taking into account missing data.
Definition: simple_stat.f90:25
This module contains functions that are only for internal use of the library.
Extension of volgrid6d_class with methods for performing simple statistical operations on entire volu...
This module defines objects and methods for managing data volumes on rectangular georeferenced grids.
Class for managing physical variables in a grib 1/2 fashion.
Class for expressing an absolute time value.
Class for expressing a relative time interval.
Derived type associated to a block/message/record/band of gridded data coming from a file-like object...
Derived type defining a dynamically extensible array of TYPE(ttr_mapper) elements.
Object describing a rectangular, homogeneous gridded dataset.
Definition of a physical variable in grib coding style.

Generated with Doxygen.