libsim Versione 7.1.11
|
◆ map_inv_distinct_ttr_mapper()
map inv distinct Definizione alla linea 1216 del file stat_proc_engine.F90. 1218! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
1219! authors:
1220! Davide Cesari <dcesari@arpa.emr.it>
1221! Paolo Patruno <ppatruno@arpa.emr.it>
1222
1223! This program is free software; you can redistribute it and/or
1224! modify it under the terms of the GNU General Public License as
1225! published by the Free Software Foundation; either version 2 of
1226! the License, or (at your option) any later version.
1227
1228! This program is distributed in the hope that it will be useful,
1229! but WITHOUT ANY WARRANTY; without even the implied warranty of
1230! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
1231! GNU General Public License for more details.
1232
1233! You should have received a copy of the GNU General Public License
1234! along with this program. If not, see <http://www.gnu.org/licenses/>.
1235#include "config.h"
1236
1243IMPLICIT NONE
1244
1245! per ogni elemento i,j dell'output, definire n elementi di input ad
1246! esso contribuenti (arrayof_ttr_mapper) con le seguenti informazioni
1247TYPE ttr_mapper
1248 INTEGER :: it=imiss ! k
1249 INTEGER :: itr=imiss ! l
1250 INTEGER :: extra_info=imiss ! dtratio for intervals, extreme for point in time
1251 TYPE(datetime) :: time=datetime_miss ! time for point in time
1252END TYPE ttr_mapper
1253
1254INTERFACE OPERATOR (==)
1255 MODULE PROCEDURE ttr_mapper_eq
1256END INTERFACE
1257
1258INTERFACE OPERATOR (/=)
1259 MODULE PROCEDURE ttr_mapper_ne
1260END INTERFACE
1261
1262INTERFACE OPERATOR (>)
1263 MODULE PROCEDURE ttr_mapper_gt
1264END INTERFACE
1265
1266INTERFACE OPERATOR (<)
1267 MODULE PROCEDURE ttr_mapper_lt
1268END INTERFACE
1269
1270INTERFACE OPERATOR (>=)
1271 MODULE PROCEDURE ttr_mapper_ge
1272END INTERFACE
1273
1274INTERFACE OPERATOR (<=)
1275 MODULE PROCEDURE ttr_mapper_le
1276END INTERFACE
1277
1278#undef VOL7D_POLY_TYPE
1279#undef VOL7D_POLY_TYPES
1280#undef ENABLE_SORT
1281#define VOL7D_POLY_TYPE TYPE(ttr_mapper)
1282#define VOL7D_POLY_TYPES _ttr_mapper
1283#define ENABLE_SORT
1284#include "array_utilities_pre.F90"
1285
1286#define ARRAYOF_ORIGTYPE TYPE(ttr_mapper)
1287#define ARRAYOF_TYPE arrayof_ttr_mapper
1288#define ARRAYOF_ORIGEQ 1
1289#define ARRAYOF_ORIGGT 1
1290#include "arrayof_pre.F90"
1291! from arrayof
1292
1293
1294CONTAINS
1295
1296! simplified comparisons only on time
1297ELEMENTAL FUNCTION ttr_mapper_eq(this, that) RESULT(res)
1298TYPE(ttr_mapper),INTENT(IN) :: this, that
1299LOGICAL :: res
1300
1301res = this%time == that%time
1302
1303END FUNCTION ttr_mapper_eq
1304
1305ELEMENTAL FUNCTION ttr_mapper_ne(this, that) RESULT(res)
1306TYPE(ttr_mapper),INTENT(IN) :: this, that
1307LOGICAL :: res
1308
1309res = this%time /= that%time
1310
1311END FUNCTION ttr_mapper_ne
1312
1313ELEMENTAL FUNCTION ttr_mapper_gt(this, that) RESULT(res)
1314TYPE(ttr_mapper),INTENT(IN) :: this, that
1315LOGICAL :: res
1316
1317res = this%time > that%time
1318
1319END FUNCTION ttr_mapper_gt
1320
1321ELEMENTAL FUNCTION ttr_mapper_lt(this, that) RESULT(res)
1322TYPE(ttr_mapper),INTENT(IN) :: this, that
1323LOGICAL :: res
1324
1325res = this%time < that%time
1326
1327END FUNCTION ttr_mapper_lt
1328
1329ELEMENTAL FUNCTION ttr_mapper_ge(this, that) RESULT(res)
1330TYPE(ttr_mapper),INTENT(IN) :: this, that
1331LOGICAL :: res
1332
1333res = this%time >= that%time
1334
1335END FUNCTION ttr_mapper_ge
1336
1337ELEMENTAL FUNCTION ttr_mapper_le(this, that) RESULT(res)
1338TYPE(ttr_mapper),INTENT(IN) :: this, that
1339LOGICAL :: res
1340
1341res = this%time <= that%time
1342
1343END FUNCTION ttr_mapper_le
1344
1345#include "arrayof_post.F90"
1346#include "array_utilities_inc.F90"
1347
1348
1349! common operations for statistical processing by differences
1350SUBROUTINE recompute_stat_proc_diff_common(itime, itimerange, stat_proc, step, &
1351 otime, otimerange, map_tr, f, keep_tr, time_definition, full_steps, &
1352 start)
1353TYPE(datetime),INTENT(in) :: itime(:)
1354TYPE(vol7d_timerange),INTENT(in) :: itimerange(:)
1355INTEGER,INTENT(in) :: stat_proc
1356TYPE(timedelta),INTENT(in) :: step
1357TYPE(datetime),POINTER :: otime(:)
1358TYPE(vol7d_timerange),POINTER :: otimerange(:)
1359INTEGER,ALLOCATABLE,INTENT(out) :: map_tr(:,:,:,:,:), f(:), keep_tr(:,:,:)
1360INTEGER :: nitr
1361LOGICAL,ALLOCATABLE :: mask_timerange(:)
1362INTEGER,INTENT(in) :: time_definition
1363LOGICAL,INTENT(in),OPTIONAL :: full_steps
1364TYPE(datetime),INTENT(in),OPTIONAL :: start
1365
1366INTEGER :: i, j, k, l, dirtyrep
1367INTEGER :: steps
1368LOGICAL :: lfull_steps, useful
1369TYPE(datetime) :: lstart, pstart1, pstart2, pend1, pend2, reftime1, reftime2, tmptime
1370TYPE(vol7d_timerange) :: tmptimerange
1371TYPE(arrayof_datetime) :: a_otime
1372TYPE(arrayof_vol7d_timerange) :: a_otimerange
1373
1374! compute length of cumulation step in seconds
1376
1377lstart = datetime_miss
1378IF (PRESENT(start)) lstart = start
1379lfull_steps = optio_log(full_steps)
1380
1381! create a mask of suitable timeranges
1382ALLOCATE(mask_timerange(SIZE(itimerange)))
1383mask_timerange(:) = itimerange(:)%timerange == stat_proc &
1384 .AND. itimerange(:)%p1 /= imiss .AND. itimerange(:)%p2 /= imiss &
1385 .AND. itimerange(:)%p1 >= 0 &
1386 .AND. itimerange(:)%p2 > 0
1387
1388IF (lfull_steps .AND. steps /= 0) THEN ! keep only timeranges defining intervals ending at integer forecast steps or analysis timeranges
1389 mask_timerange(:) = mask_timerange(:) .AND. &
1390 (itimerange(:)%p1 == 0 .OR. & ! all analyses
1393ENDIF
1394! mask_timerange includes all candidate timeranges
1395
1396nitr = count(mask_timerange)
1397ALLOCATE(f(nitr))
1398j = 1
1399DO i = 1, nitr
1400 DO WHILE(.NOT.mask_timerange(j))
1401 j = j + 1
1402 ENDDO
1403 f(i) = j
1404 j = j + 1
1405ENDDO
1406
1407! now we have to evaluate time/timerage pairs which do not need processing
1408ALLOCATE(keep_tr(nitr, SIZE(itime), 2))
1409CALL compute_keep_tr()
1410
1411ALLOCATE(map_tr(nitr, SIZE(itime), nitr, SIZE(itime), 2))
1412map_tr(:,:,:,:,:) = imiss
1413
1414! scan through all possible combinations of time and timerange
1415DO dirtyrep = 1, 2
1416 IF (dirtyrep == 2) THEN ! dirty and expensive trick for sorting descriptors
1421 ENDIF
1422 DO l = 1, SIZE(itime)
1423 DO k = 1, nitr
1424 CALL time_timerange_get_period(itime(l), itimerange(f(k)), &
1425 time_definition, pstart2, pend2, reftime2)
1426
1427 DO j = 1, SIZE(itime)
1428 DO i = 1, nitr
1429 useful = .false.
1430 CALL time_timerange_get_period(itime(j), itimerange(f(i)), &
1431 time_definition, pstart1, pend1, reftime1)
1432 tmptimerange = vol7d_timerange_new(timerange=stat_proc)
1433
1434 IF (reftime2 == pend2 .AND. reftime1 == pend1) THEN ! analysis
1435 IF (pstart2 == pstart1 .AND. pend2 > pend1) THEN ! =-|
1436 CALL time_timerange_set_period(tmptime, tmptimerange, &
1437 time_definition, pend1, pend2, reftime2)
1438 IF (lfull_steps) THEN
1440 useful = .true.
1441 ENDIF
1442 ELSE
1443 useful = .true.
1444 ENDIF
1445
1446 ELSE IF (pstart2 < pstart1 .AND. pend2 == pend1) THEN ! -=|
1447 CALL time_timerange_set_period(tmptime, tmptimerange, &
1448 time_definition, pstart2, pstart1, pstart1)
1449 IF (lfull_steps) THEN
1451 useful = .true.
1452 ENDIF
1453 ELSE
1454 useful = .true.
1455 ENDIF
1456 ENDIF
1457
1458 ELSE IF (reftime2 == reftime1) THEN ! forecast, same reftime
1459 IF (pstart2 == pstart1 .AND. pend2 > pend1) THEN ! |=-
1460 CALL time_timerange_set_period(tmptime, tmptimerange, &
1461 time_definition, pend1, pend2, reftime2)
1462 IF (lfull_steps) THEN
1464 useful = .true.
1465 ENDIF
1466 ELSE
1467 useful = .true.
1468 ENDIF
1469
1470 ELSE IF (pstart2 < pstart1 .AND. pend2 == pend1) THEN ! |-=
1471 CALL time_timerange_set_period(tmptime, tmptimerange, &
1472 time_definition, pstart2, pstart1, reftime2)
1473 IF (lfull_steps) THEN
1475 useful = .true.
1476 ENDIF
1477 ELSE
1478 useful = .true.
1479 ENDIF
1480
1481 ENDIF
1482 ENDIF
1483 useful = useful .AND. tmptime /= datetime_miss .AND. &
1484 tmptimerange /= vol7d_timerange_miss .AND. tmptimerange%p2 == steps
1485
1486 IF (useful) THEN ! add a_otime, a_otimerange
1487 map_tr(i,j,k,l,1) = append_unique(a_otime, tmptime)
1488 map_tr(i,j,k,l,2) = append_unique(a_otimerange, tmptimerange)
1489 ENDIF
1490 ENDDO
1491 ENDDO
1492 ENDDO
1493 ENDDO
1494ENDDO
1495
1496! we have to repeat the computation with sorted arrays
1497CALL compute_keep_tr()
1498
1499otime => a_otime%array
1500otimerange => a_otimerange%array
1501! delete local objects keeping the contents
1504
1505#ifdef DEBUG
1506CALL l4f_log(l4f_debug, &
1511CALL l4f_log(l4f_debug, &
1514CALL l4f_log(l4f_debug, &
1516CALL l4f_log(l4f_debug, &
1518CALL l4f_log(l4f_debug, &
1520CALL l4f_log(l4f_debug, &
1522#endif
1523
1524CONTAINS
1525
1526SUBROUTINE compute_keep_tr()
1527
1528keep_tr(:,:,:) = imiss
1529DO l = 1, SIZE(itime)
1530 DO k = 1, nitr
1531 IF (itimerange(f(k))%p2 == steps) THEN
1532 CALL time_timerange_get_period(itime(l), itimerange(f(k)), &
1533 time_definition, pstart2, pend2, reftime2)
1534 useful = .false.
1535 IF (reftime2 == pend2) THEN ! analysis
1538 useful = .true.
1539 ENDIF
1540 ELSE IF (lfull_steps) THEN
1542 useful = .true.
1543 ENDIF
1544 ELSE
1545 useful = .true.
1546 ENDIF
1547 ELSE ! forecast
1548 IF (lfull_steps) THEN
1550 useful = .true.
1551 ENDIF
1552 ELSE
1553 useful = .true.
1554 ENDIF
1555 ENDIF
1556 IF (useful) THEN
1557! CALL time_timerange_set_period(tmptime, tmptimerange, &
1558! time_definition, pstart2, pend2, reftime2)
1559 keep_tr(k,l,1) = append_unique(a_otime, itime(l))
1560 keep_tr(k,l,2) = append_unique(a_otimerange, itimerange(f(k)))
1561 ENDIF
1562 ENDIF
1563 ENDDO
1564ENDDO
1565
1566END SUBROUTINE compute_keep_tr
1567
1568END SUBROUTINE recompute_stat_proc_diff_common
1569
1570
1571! common operations for statistical processing by metamorphosis
1572SUBROUTINE compute_stat_proc_metamorph_common(istat_proc, itimerange, ostat_proc, &
1573 otimerange, map_tr)
1574INTEGER,INTENT(in) :: istat_proc
1575TYPE(vol7d_timerange),INTENT(in) :: itimerange(:)
1576INTEGER,INTENT(in) :: ostat_proc
1577TYPE(vol7d_timerange),POINTER :: otimerange(:)
1578INTEGER,POINTER :: map_tr(:)
1579
1580INTEGER :: i
1581LOGICAL :: tr_mask(SIZE(itimerange))
1582
1583IF (SIZE(itimerange) == 0) THEN ! avoid segmentation fault in case of empty volume
1584 ALLOCATE(otimerange(0), map_tr(0))
1585 RETURN
1586ENDIF
1587
1588! useful timeranges
1589tr_mask(:) = itimerange(:)%timerange == istat_proc .AND. itimerange(:)%p2 /= imiss &
1590 .AND. itimerange(:)%p2 /= 0 ! .AND. itimerange(:)%p1 == 0
1591ALLOCATE(otimerange(count(tr_mask)), map_tr(count(tr_mask)))
1592
1593otimerange = pack(itimerange, mask=tr_mask)
1594otimerange(:)%timerange = ostat_proc
1595map_tr = pack((/(i,i=1,SIZE(itimerange))/), mask=tr_mask)
1596
1597END SUBROUTINE compute_stat_proc_metamorph_common
1598
1599
1600! common operations for statistical processing by aggregation
1601SUBROUTINE recompute_stat_proc_agg_common(itime, itimerange, stat_proc, tri, &
1602 step, time_definition, otime, otimerange, map_ttr, dtratio, start, full_steps)
1603TYPE(datetime),INTENT(in) :: itime(:)
1604TYPE(vol7d_timerange),INTENT(in) :: itimerange(:)
1605INTEGER,INTENT(in) :: stat_proc
1606INTEGER,INTENT(in) :: tri
1607TYPE(timedelta),INTENT(in) :: step
1608INTEGER,INTENT(in) :: time_definition
1609TYPE(datetime),POINTER :: otime(:)
1610TYPE(vol7d_timerange),POINTER :: otimerange(:)
1611TYPE(arrayof_ttr_mapper),POINTER :: map_ttr(:,:)
1612INTEGER,POINTER,OPTIONAL :: dtratio(:)
1613TYPE(datetime),INTENT(in),OPTIONAL :: start
1614LOGICAL,INTENT(in),OPTIONAL :: full_steps
1615
1616INTEGER :: i, j, k, l, na, nf, n
1617INTEGER :: steps, p1, maxp1, maxp2, minp1mp2, dstart, msteps
1618INTEGER(kind=int_ll) :: stepms, mstepms
1619LOGICAL :: lforecast
1620TYPE(datetime) :: lstart, lend, pstart1, pstart2, pend1, pend2, reftime1, reftime2, tmptime
1621TYPE(arrayof_datetime) :: a_otime
1622TYPE(arrayof_vol7d_timerange) :: a_otimerange
1623TYPE(arrayof_integer) :: a_dtratio
1624LOGICAL,ALLOCATABLE :: mask_timerange(:) ! improve !!!!
1625TYPE(ttr_mapper) :: lmapper
1626CHARACTER(len=8) :: env_var
1627LOGICAL :: climat_behavior
1628
1629
1630! enable bad behavior for climat database (only for instantaneous data)
1631env_var = ''
1632CALL getenv('LIBSIM_CLIMAT_BEHAVIOR', env_var)
1633climat_behavior = len_trim(env_var) > 0 .AND. .NOT.PRESENT(dtratio)
1634
1635! compute length of cumulation step in seconds
1637
1638! create a mask of suitable timeranges
1639ALLOCATE(mask_timerange(SIZE(itimerange)))
1640mask_timerange(:) = itimerange(:)%timerange == tri &
1641 .AND. itimerange(:)%p1 /= imiss .AND. itimerange(:)%p1 >= 0
1642
1643IF (PRESENT(dtratio)) THEN
1644 WHERE(itimerange(:)%p2 > 0 .AND. itimerange(:)%p2 /= imiss) ! for avoiding FPE MOD(steps, 0)
1645 mask_timerange(:) = mask_timerange(:) .AND. mod(steps, itimerange(:)%p2) == 0
1646 ELSEWHERE
1647 mask_timerange(:) = .false.
1648 END WHERE
1649ELSE ! instantaneous
1650 mask_timerange(:) = mask_timerange(:) .AND. itimerange(:)%p2 == 0
1651ENDIF
1652
1653#ifdef DEBUG
1654CALL l4f_log(l4f_debug, &
1655 '(re)compute_stat_proc_agg, number of useful timeranges before choosing analysis/forecast: '// &
1656 t2c(count(mask_timerange)))
1657#endif
1658
1659! euristically determine whether we are dealing with an
1660! analysis/observation or a forecast dataset
1661na = count(mask_timerange(:) .AND. itimerange(:)%p1 == 0)
1662nf = count(mask_timerange(:) .AND. itimerange(:)%p1 > 0)
1663lforecast = nf >= na
1664#ifdef DEBUG
1665CALL l4f_log(l4f_debug, &
1667#endif
1668
1669! keep only timeranges of one type (really necessary?)
1670IF (lforecast) THEN
1671 CALL l4f_log(l4f_info, &
1672 'recompute_stat_proc_agg, processing in forecast mode')
1673ELSE
1674 mask_timerange(:) = mask_timerange(:) .AND. itimerange(:)%p1 == 0
1675 CALL l4f_log(l4f_info, &
1676 'recompute_stat_proc_agg, processing in analysis mode')
1677ENDIF
1678
1679#ifdef DEBUG
1680CALL l4f_log(l4f_debug, &
1681 '(re)compute_stat_proc_agg, number of useful timeranges: '// &
1682 t2c(count(mask_timerange)))
1683#endif
1684
1685IF (SIZE(itime) == 0 .OR. count(mask_timerange) == 0) THEN ! avoid segmentation fault in case of empty volume
1686 ALLOCATE(otime(0), otimerange(0), map_ttr(0,0))
1687 IF (PRESENT(dtratio)) ALLOCATE(dtratio(0))
1688 RETURN
1689ENDIF
1690
1691! determine start and end of processing period, should work with p2==0
1692lstart = datetime_miss
1693IF (PRESENT(start)) lstart = start
1694lend = itime(SIZE(itime))
1695! compute some quantities
1696maxp1 = maxval(itimerange(:)%p1, mask=mask_timerange)
1697maxp2 = maxval(itimerange(:)%p2, mask=mask_timerange)
1698minp1mp2 = minval(itimerange(:)%p1 - itimerange(:)%p2, mask=mask_timerange)
1699IF (time_definition == 0) THEN ! reference time
1700 lend = lend + timedelta_new(sec=maxp1)
1701ENDIF
1702! extend interval at the end in order to include all the data in case
1703! frac_valid<1; must use < and not <= in the DO WHILE loops some lines
1704! below in order to exclude the last full interval which would be empty
1705lend = lend + step
1706IF (lstart == datetime_miss) THEN ! autodetect
1707 lstart = itime(1)
1708! if autodetected, adjust to obtain real absolute start of data
1709 IF (time_definition == 0) THEN ! reference time
1710 lstart = lstart + timedelta_new(sec=minp1mp2)
1711 ELSE ! verification time
1712! go back to start of longest processing interval
1713 lstart = lstart - timedelta_new(sec=maxp2)
1714 ENDIF
1715! full_steps is effective only in analysis mode and when start is not
1716! specified (start by itself in analysis mode implies full_steps with
1717! respect to start instead of absolute full steps)
1718 IF (optio_log(full_steps) .AND. .NOT.lforecast) THEN
1720 ENDIF
1721ENDIF
1722
1723#ifdef DEBUG
1724CALL l4f_log(l4f_debug, &
1726#endif
1727
1728! create output time and timerange lists
1729
1730IF (lforecast) THEN ! forecast mode
1731 IF (time_definition == 0) THEN ! reference time
1733
1734! apply start shift to timerange, not to reftime
1735! why did we use itime(SIZE(itime)) (last ref time)?
1736! CALL getval(lstart-itime(SIZE(itime)), asec=dstart)
1738! allow starting before first reftime but restrict dtstart to -steps
1739! dstart = MAX(0, dstart)
1741 DO p1 = steps + dstart, maxp1, steps
1742 CALL insert_unique(a_otimerange, vol7d_timerange_new(stat_proc, p1, steps))
1743 ENDDO
1744
1745 ELSE ! verification time
1746
1747! verification time in forecast mode is the ugliest case, because we
1748! don't know a priori how many different (thus incompatible) reference
1749! times we have, so some assumption of regularity has to be made. For
1750! this reason msteps, the minimum step between two times, is
1751! computed. We choose to compute it as a difference between itime
1752! elements, it could be also computed as difference of itimerange%p1
1753! elements. But what if it is not modulo steps?
1754 mstepms = steps*1000_int_ll
1755 DO i = 2, SIZE(itime)
1757 IF (stepms > 0_int_ll .AND. stepms < mstepms) THEN
1758 msteps = stepms/1000_int_ll
1760 ENDIF
1761 ENDDO
1762 msteps = mstepms/1000_int_ll
1763
1764 tmptime = lstart + step
1765 DO WHILE(tmptime < lend) ! < since lend has been extended
1766 CALL insert_unique(a_otime, tmptime)
1767 tmptime = tmptime + step
1768 ENDDO
1769
1770! in next loop, we used initially steps instead of msteps (see comment
1771! above), this gave much less combinations of time/timeranges with
1772! possible empty output; we start from msteps instead of steps in
1773! order to include partial initial intervals in case frac_valid<1;
1774! however, as a gemeral rule, for aggregation of forecasts it's better
1775! to use reference time
1776 DO p1 = msteps, maxp1, msteps
1777 CALL insert_unique(a_otimerange, vol7d_timerange_new(stat_proc, p1, steps))
1778 ENDDO
1779
1780 ENDIF
1781
1782ELSE ! analysis mode
1783 tmptime = lstart + step
1784 DO WHILE(tmptime < lend) ! < since lend has been extended
1785 CALL insert_unique(a_otime, tmptime)
1786 tmptime = tmptime + step
1787 ENDDO
1788 CALL insert_unique(a_otimerange, vol7d_timerange_new(stat_proc, 0, steps))
1789
1790ENDIF
1791
1794otime => a_otime%array
1795otimerange => a_otimerange%array
1798! delete local objects keeping the contents
1801
1802#ifdef DEBUG
1803CALL l4f_log(l4f_debug, &
1804 'recompute_stat_proc_agg, output time and timerange: '//&
1806#endif
1807
1808IF (PRESENT(dtratio)) THEN
1809! count the possible i/o interval ratios
1810 DO k = 1, SIZE(itimerange)
1811 IF (itimerange(k)%p2 /= 0) &
1812 CALL insert_unique(a_dtratio, steps/itimerange(k)%p2) ! guaranteed to be integer
1813 ENDDO
1815 dtratio => a_dtratio%array
1817! delete local object keeping the contents
1819
1820#ifdef DEBUG
1821 CALL l4f_log(l4f_debug, &
1823 ' possible aggregation ratios, from '// &
1825#endif
1826
1827 ALLOCATE(map_ttr(SIZE(otime),SIZE(otimerange)))
1828 do_itimerange1: DO l = 1, SIZE(itimerange)
1829 IF (.NOT.mask_timerange(l)) cycle do_itimerange1
1830 do_itime1: DO k = 1, SIZE(itime)
1831 CALL time_timerange_get_period(itime(k), itimerange(l), &
1832 time_definition, pstart1, pend1, reftime1)
1833 do_otimerange1: DO j = 1, SIZE(otimerange)
1834 do_otime1: DO i = 1, SIZE(otime)
1835 CALL time_timerange_get_period_pop(otime(i), otimerange(j), step, &
1836 time_definition, pstart2, pend2, reftime2)
1837 IF (lforecast) THEN
1838 IF (reftime1 /= reftime2) cycle do_otime1
1839 ENDIF
1840
1841 IF (pstart1 >= pstart2 .AND. pend1 <= pend2 .AND. &
1843 lmapper%it = k
1844 lmapper%itr = l
1845 lmapper%extra_info = steps/itimerange(l)%p2 ! dtratio, guaranteed to be integer
1846 n = append(map_ttr(i,j), lmapper)
1847 cycle do_itime1 ! can contribute only to a single interval
1848 ENDIF
1849 ENDDO do_otime1
1850 ENDDO do_otimerange1
1851 ENDDO do_itime1
1852 ENDDO do_itimerange1
1853
1854ELSE
1855
1856 ALLOCATE(map_ttr(SIZE(otime),SIZE(otimerange)))
1857 do_itimerange2: DO l = 1, SIZE(itimerange)
1858 IF (.NOT.mask_timerange(l)) cycle do_itimerange2
1859 do_itime2: DO k = 1, SIZE(itime)
1860 CALL time_timerange_get_period(itime(k), itimerange(l), &
1861 time_definition, pstart1, pend1, reftime1)
1862 do_otimerange2: DO j = 1, SIZE(otimerange)
1863 do_otime2: DO i = 1, SIZE(otime)
1864 CALL time_timerange_get_period_pop(otime(i), otimerange(j), step, &
1865 time_definition, pstart2, pend2, reftime2)
1866 IF (lforecast) THEN
1867 IF (reftime1 /= reftime2) cycle do_otime2
1868 ENDIF
1869
1870 IF (climat_behavior .AND. pstart1 == pstart2) cycle do_otime2
1871 IF (pstart1 >= pstart2 .AND. pend1 <= pend2) THEN ! useful
1872 lmapper%it = k
1873 lmapper%itr = l
1874 IF (pstart1 == pstart2) THEN
1875 lmapper%extra_info = 1 ! start of interval
1876 ELSE IF (pend1 == pend2) THEN
1877 lmapper%extra_info = 2 ! end of interval
1878 ELSE
1879 lmapper%extra_info = imiss
1880 ENDIF
1881 lmapper%time = pstart1 ! = pend1, order by time?
1882 n = insert_sorted(map_ttr(i,j), lmapper, .true., .true.)
1883! no CYCLE, a single input can contribute to multiple output intervals
1884 ENDIF
1885 ENDDO do_otime2
1886 ENDDO do_otimerange2
1887 ENDDO do_itime2
1888 ENDDO do_itimerange2
1889
1890ENDIF
1891
1892END SUBROUTINE recompute_stat_proc_agg_common
1893
1894
1895SUBROUTINE compute_stat_proc_agg_sw(vertime, pstart, pend, time_mask, &
1896 max_step, weights)
1897TYPE(datetime),INTENT(in) :: vertime(:)
1898TYPE(datetime),INTENT(in) :: pstart
1899TYPE(datetime),INTENT(in) :: pend
1900LOGICAL,INTENT(in) :: time_mask(:)
1901TYPE(timedelta),OPTIONAL,INTENT(out) :: max_step
1902DOUBLE PRECISION,OPTIONAL,INTENT(out) :: weights(:)
1903
1904INTEGER :: i, nt
1905TYPE(datetime),ALLOCATABLE :: lvertime(:)
1906TYPE(datetime) :: half, nexthalf
1907INTEGER(kind=int_ll) :: dt, tdt
1908
1909nt = count(time_mask)
1910ALLOCATE(lvertime(nt))
1911lvertime = pack(vertime, mask=time_mask)
1912
1913IF (PRESENT(max_step)) THEN
1914! new way
1915! max_step = timedelta_0
1916! DO i = 1, nt - 1
1917! IF (lvertime(i+1) - lvertime(i) > max_step) &
1918! max_step = lvertime(i+1) - lvertime(i)
1919! ENDDO
1920! IF (lvertime(1) - pstart > max_step) max_step = lvertime(1) - pstart
1921! IF (pend - lvertime(nt) > max_step) max_step = pend - lvertime(nt)
1922! old way
1923 IF (nt == 1) THEN
1924 max_step = pend - pstart
1925 ELSE
1926 half = lvertime(1) + (lvertime(2) - lvertime(1))/2
1927 max_step = half - pstart
1928 DO i = 2, nt - 1
1929 nexthalf = lvertime(i) + (lvertime(i+1) - lvertime(i))/2
1930 IF (nexthalf - half > max_step) max_step = nexthalf - half
1931 half = nexthalf
1932 ENDDO
1933 IF (pend - half > max_step) max_step = pend - half
1934 ENDIF
1935
1936ENDIF
1937
1938IF (PRESENT(weights)) THEN
1939 IF (nt == 1) THEN
1940 weights(1) = 1.0
1941 ELSE
1943 half = lvertime(1) + (lvertime(2) - lvertime(1))/2
1945 weights(1) = dble(dt)/dble(tdt)
1946 DO i = 2, nt - 1
1947 nexthalf = lvertime(i) + (lvertime(i+1) - lvertime(i))/2
1949 weights(i) = dble(dt)/dble(tdt)
1950 half = nexthalf
1951 ENDDO
1953 weights(nt) = dble(dt)/dble(tdt)
1954 ENDIF
1955ENDIF
1956
1957END SUBROUTINE compute_stat_proc_agg_sw
1958
1959! get start of period, end of period and reference time from time,
1960! timerange, according to time_definition.
1961SUBROUTINE time_timerange_get_period(time, timerange, time_definition, &
1962 pstart, pend, reftime)
1963TYPE(datetime),INTENT(in) :: time
1964TYPE(vol7d_timerange),INTENT(in) :: timerange
1965INTEGER,INTENT(in) :: time_definition
1966TYPE(datetime),INTENT(out) :: reftime
1967TYPE(datetime),INTENT(out) :: pstart
1968TYPE(datetime),INTENT(out) :: pend
1969
1970TYPE(timedelta) :: p1, p2
1971
1972
1973p1 = timedelta_new(sec=timerange%p1) ! end of period
1974p2 = timedelta_new(sec=timerange%p2) ! length of period
1975
1977! (timerange%p1 > 0 .AND. timerange%p1 < timerange%p2) .OR. &
1978 timerange%p1 < 0 .OR. timerange%p2 < 0) THEN ! is this too pedantic and slow?
1979 pstart = datetime_miss
1980 pend = datetime_miss
1981 reftime = datetime_miss
1982 RETURN
1983ENDIF
1984
1985IF (time_definition == 0) THEN ! time == reference time
1986 reftime = time
1987 pend = time + p1
1988 pstart = pend - p2
1989ELSE IF (time_definition == 1) THEN ! time == verification time
1990 pend = time
1991 pstart = time - p2
1992 reftime = time - p1
1993ELSE
1994 pstart = datetime_miss
1995 pend = datetime_miss
1996 reftime = datetime_miss
1997ENDIF
1998
1999END SUBROUTINE time_timerange_get_period
2000
2001
2002! get start of period, end of period and reference time from time,
2003! timerange, according to time_definition. step is taken separately
2004! from timerange and can be "popular"
2005SUBROUTINE time_timerange_get_period_pop(time, timerange, step, time_definition, &
2006 pstart, pend, reftime)
2007TYPE(datetime),INTENT(in) :: time
2008TYPE(vol7d_timerange),INTENT(in) :: timerange
2009TYPE(timedelta),INTENT(in) :: step
2010INTEGER,INTENT(in) :: time_definition
2011TYPE(datetime),INTENT(out) :: reftime
2012TYPE(datetime),INTENT(out) :: pstart
2013TYPE(datetime),INTENT(out) :: pend
2014
2015TYPE(timedelta) :: p1
2016
2017
2018p1 = timedelta_new(sec=timerange%p1) ! end of period
2019
2021! (timerange%p1 > 0 .AND. timerange%p1 < timerange%p2) .OR. &
2022 timerange%p1 < 0 .OR. timerange%p2 < 0) THEN ! is this too pedantic and slow?
2023 pstart = datetime_miss
2024 pend = datetime_miss
2025 reftime = datetime_miss
2026 RETURN
2027ENDIF
2028
2029IF (time_definition == 0) THEN ! time == reference time
2030 reftime = time
2031 pend = time + p1
2032 pstart = pend - step
2033ELSE IF (time_definition == 1) THEN ! time == verification time
2034 pend = time
2035 pstart = time - step
2036 reftime = time - p1
2037ELSE
2038 pstart = datetime_miss
2039 pend = datetime_miss
2040 reftime = datetime_miss
2041ENDIF
2042
2043END SUBROUTINE time_timerange_get_period_pop
2044
2045
2046! set time, timerange%p1, timerange%p2 according to pstart, pend,
2047! reftime and time_definition.
2048SUBROUTINE time_timerange_set_period(time, timerange, time_definition, &
2049 pstart, pend, reftime)
2050TYPE(datetime),INTENT(out) :: time
2051TYPE(vol7d_timerange),INTENT(inout) :: timerange
2052INTEGER,INTENT(in) :: time_definition
2053TYPE(datetime),INTENT(in) :: reftime
2054TYPE(datetime),INTENT(in) :: pstart
2055TYPE(datetime),INTENT(in) :: pend
2056
2057TYPE(timedelta) :: p1, p2
2058INTEGER(kind=int_ll) :: dmsec
2059
2060
2061IF (time_definition == 0) THEN ! time == reference time
2062 time = reftime
2063 p1 = pend - reftime
2064 p2 = pend - pstart
2065ELSE IF (time_definition == 1) THEN ! time == verification time
2066 time = pend
2067 p1 = pend - reftime
2068 p2 = pend - pstart
2069ELSE
2070 time = datetime_miss
2071ENDIF
2072
2073IF (time /= datetime_miss) THEN
2075 timerange%p1 = int(dmsec/1000_int_ll)
2077 timerange%p2 = int(dmsec/1000_int_ll)
2078ELSE
2079 timerange%p1 = imiss
2080 timerange%p2 = imiss
2081ENDIF
2082
2083END SUBROUTINE time_timerange_set_period
2084
2085
Restituiscono il valore dell'oggetto nella forma desiderata. Definition: datetime_class.F90:328 Functions that return a trimmed CHARACTER representation of the input variable. Definition: datetime_class.F90:355 Quick method to append an element to the array. Definition: stat_proc_engine.F90:371 Destructor for finalizing an array object. Definition: stat_proc_engine.F90:384 Method for inserting elements of the array at a desired position. Definition: stat_proc_engine.F90:362 Method for packing the array object reducing at a minimum the memory occupation, without destroying i... Definition: stat_proc_engine.F90:394 This module contains functions that are only for internal use of the library. Definition: stat_proc_engine.F90:217 Classe per la gestione di un volume completo di dati osservati. Definition: vol7d_class.F90:279 |