libsim Versione 7.1.11
|
◆ map_distinct_ttr_mapper()
map distinct Definizione alla linea 1120 del file stat_proc_engine.F90. 1121! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
1122! authors:
1123! Davide Cesari <dcesari@arpa.emr.it>
1124! Paolo Patruno <ppatruno@arpa.emr.it>
1125
1126! This program is free software; you can redistribute it and/or
1127! modify it under the terms of the GNU General Public License as
1128! published by the Free Software Foundation; either version 2 of
1129! the License, or (at your option) any later version.
1130
1131! This program is distributed in the hope that it will be useful,
1132! but WITHOUT ANY WARRANTY; without even the implied warranty of
1133! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
1134! GNU General Public License for more details.
1135
1136! You should have received a copy of the GNU General Public License
1137! along with this program. If not, see <http://www.gnu.org/licenses/>.
1138#include "config.h"
1139
1146IMPLICIT NONE
1147
1148! per ogni elemento i,j dell'output, definire n elementi di input ad
1149! esso contribuenti (arrayof_ttr_mapper) con le seguenti informazioni
1150TYPE ttr_mapper
1151 INTEGER :: it=imiss ! k
1152 INTEGER :: itr=imiss ! l
1153 INTEGER :: extra_info=imiss ! dtratio for intervals, extreme for point in time
1154 TYPE(datetime) :: time=datetime_miss ! time for point in time
1155END TYPE ttr_mapper
1156
1157INTERFACE OPERATOR (==)
1158 MODULE PROCEDURE ttr_mapper_eq
1159END INTERFACE
1160
1161INTERFACE OPERATOR (/=)
1162 MODULE PROCEDURE ttr_mapper_ne
1163END INTERFACE
1164
1165INTERFACE OPERATOR (>)
1166 MODULE PROCEDURE ttr_mapper_gt
1167END INTERFACE
1168
1169INTERFACE OPERATOR (<)
1170 MODULE PROCEDURE ttr_mapper_lt
1171END INTERFACE
1172
1173INTERFACE OPERATOR (>=)
1174 MODULE PROCEDURE ttr_mapper_ge
1175END INTERFACE
1176
1177INTERFACE OPERATOR (<=)
1178 MODULE PROCEDURE ttr_mapper_le
1179END INTERFACE
1180
1181#undef VOL7D_POLY_TYPE
1182#undef VOL7D_POLY_TYPES
1183#undef ENABLE_SORT
1184#define VOL7D_POLY_TYPE TYPE(ttr_mapper)
1185#define VOL7D_POLY_TYPES _ttr_mapper
1186#define ENABLE_SORT
1187#include "array_utilities_pre.F90"
1188
1189#define ARRAYOF_ORIGTYPE TYPE(ttr_mapper)
1190#define ARRAYOF_TYPE arrayof_ttr_mapper
1191#define ARRAYOF_ORIGEQ 1
1192#define ARRAYOF_ORIGGT 1
1193#include "arrayof_pre.F90"
1194! from arrayof
1195
1196
1197CONTAINS
1198
1199! simplified comparisons only on time
1200ELEMENTAL FUNCTION ttr_mapper_eq(this, that) RESULT(res)
1201TYPE(ttr_mapper),INTENT(IN) :: this, that
1202LOGICAL :: res
1203
1204res = this%time == that%time
1205
1206END FUNCTION ttr_mapper_eq
1207
1208ELEMENTAL FUNCTION ttr_mapper_ne(this, that) RESULT(res)
1209TYPE(ttr_mapper),INTENT(IN) :: this, that
1210LOGICAL :: res
1211
1212res = this%time /= that%time
1213
1214END FUNCTION ttr_mapper_ne
1215
1216ELEMENTAL FUNCTION ttr_mapper_gt(this, that) RESULT(res)
1217TYPE(ttr_mapper),INTENT(IN) :: this, that
1218LOGICAL :: res
1219
1220res = this%time > that%time
1221
1222END FUNCTION ttr_mapper_gt
1223
1224ELEMENTAL FUNCTION ttr_mapper_lt(this, that) RESULT(res)
1225TYPE(ttr_mapper),INTENT(IN) :: this, that
1226LOGICAL :: res
1227
1228res = this%time < that%time
1229
1230END FUNCTION ttr_mapper_lt
1231
1232ELEMENTAL FUNCTION ttr_mapper_ge(this, that) RESULT(res)
1233TYPE(ttr_mapper),INTENT(IN) :: this, that
1234LOGICAL :: res
1235
1236res = this%time >= that%time
1237
1238END FUNCTION ttr_mapper_ge
1239
1240ELEMENTAL FUNCTION ttr_mapper_le(this, that) RESULT(res)
1241TYPE(ttr_mapper),INTENT(IN) :: this, that
1242LOGICAL :: res
1243
1244res = this%time <= that%time
1245
1246END FUNCTION ttr_mapper_le
1247
1248#include "arrayof_post.F90"
1249#include "array_utilities_inc.F90"
1250
1251
1252! common operations for statistical processing by differences
1253SUBROUTINE recompute_stat_proc_diff_common(itime, itimerange, stat_proc, step, &
1254 otime, otimerange, map_tr, f, keep_tr, time_definition, full_steps, &
1255 start)
1256TYPE(datetime),INTENT(in) :: itime(:)
1257TYPE(vol7d_timerange),INTENT(in) :: itimerange(:)
1258INTEGER,INTENT(in) :: stat_proc
1259TYPE(timedelta),INTENT(in) :: step
1260TYPE(datetime),POINTER :: otime(:)
1261TYPE(vol7d_timerange),POINTER :: otimerange(:)
1262INTEGER,ALLOCATABLE,INTENT(out) :: map_tr(:,:,:,:,:), f(:), keep_tr(:,:,:)
1263INTEGER :: nitr
1264LOGICAL,ALLOCATABLE :: mask_timerange(:)
1265INTEGER,INTENT(in) :: time_definition
1266LOGICAL,INTENT(in),OPTIONAL :: full_steps
1267TYPE(datetime),INTENT(in),OPTIONAL :: start
1268
1269INTEGER :: i, j, k, l, dirtyrep
1270INTEGER :: steps
1271LOGICAL :: lfull_steps, useful
1272TYPE(datetime) :: lstart, pstart1, pstart2, pend1, pend2, reftime1, reftime2, tmptime
1273TYPE(vol7d_timerange) :: tmptimerange
1274TYPE(arrayof_datetime) :: a_otime
1275TYPE(arrayof_vol7d_timerange) :: a_otimerange
1276
1277! compute length of cumulation step in seconds
1279
1280lstart = datetime_miss
1281IF (PRESENT(start)) lstart = start
1282lfull_steps = optio_log(full_steps)
1283
1284! create a mask of suitable timeranges
1285ALLOCATE(mask_timerange(SIZE(itimerange)))
1286mask_timerange(:) = itimerange(:)%timerange == stat_proc &
1287 .AND. itimerange(:)%p1 /= imiss .AND. itimerange(:)%p2 /= imiss &
1288 .AND. itimerange(:)%p1 >= 0 &
1289 .AND. itimerange(:)%p2 > 0
1290
1291IF (lfull_steps .AND. steps /= 0) THEN ! keep only timeranges defining intervals ending at integer forecast steps or analysis timeranges
1292 mask_timerange(:) = mask_timerange(:) .AND. &
1293 (itimerange(:)%p1 == 0 .OR. & ! all analyses
1296ENDIF
1297! mask_timerange includes all candidate timeranges
1298
1299nitr = count(mask_timerange)
1300ALLOCATE(f(nitr))
1301j = 1
1302DO i = 1, nitr
1303 DO WHILE(.NOT.mask_timerange(j))
1304 j = j + 1
1305 ENDDO
1306 f(i) = j
1307 j = j + 1
1308ENDDO
1309
1310! now we have to evaluate time/timerage pairs which do not need processing
1311ALLOCATE(keep_tr(nitr, SIZE(itime), 2))
1312CALL compute_keep_tr()
1313
1314ALLOCATE(map_tr(nitr, SIZE(itime), nitr, SIZE(itime), 2))
1315map_tr(:,:,:,:,:) = imiss
1316
1317! scan through all possible combinations of time and timerange
1318DO dirtyrep = 1, 2
1319 IF (dirtyrep == 2) THEN ! dirty and expensive trick for sorting descriptors
1324 ENDIF
1325 DO l = 1, SIZE(itime)
1326 DO k = 1, nitr
1327 CALL time_timerange_get_period(itime(l), itimerange(f(k)), &
1328 time_definition, pstart2, pend2, reftime2)
1329
1330 DO j = 1, SIZE(itime)
1331 DO i = 1, nitr
1332 useful = .false.
1333 CALL time_timerange_get_period(itime(j), itimerange(f(i)), &
1334 time_definition, pstart1, pend1, reftime1)
1335 tmptimerange = vol7d_timerange_new(timerange=stat_proc)
1336
1337 IF (reftime2 == pend2 .AND. reftime1 == pend1) THEN ! analysis
1338 IF (pstart2 == pstart1 .AND. pend2 > pend1) THEN ! =-|
1339 CALL time_timerange_set_period(tmptime, tmptimerange, &
1340 time_definition, pend1, pend2, reftime2)
1341 IF (lfull_steps) THEN
1343 useful = .true.
1344 ENDIF
1345 ELSE
1346 useful = .true.
1347 ENDIF
1348
1349 ELSE IF (pstart2 < pstart1 .AND. pend2 == pend1) THEN ! -=|
1350 CALL time_timerange_set_period(tmptime, tmptimerange, &
1351 time_definition, pstart2, pstart1, pstart1)
1352 IF (lfull_steps) THEN
1354 useful = .true.
1355 ENDIF
1356 ELSE
1357 useful = .true.
1358 ENDIF
1359 ENDIF
1360
1361 ELSE IF (reftime2 == reftime1) THEN ! forecast, same reftime
1362 IF (pstart2 == pstart1 .AND. pend2 > pend1) THEN ! |=-
1363 CALL time_timerange_set_period(tmptime, tmptimerange, &
1364 time_definition, pend1, pend2, reftime2)
1365 IF (lfull_steps) THEN
1367 useful = .true.
1368 ENDIF
1369 ELSE
1370 useful = .true.
1371 ENDIF
1372
1373 ELSE IF (pstart2 < pstart1 .AND. pend2 == pend1) THEN ! |-=
1374 CALL time_timerange_set_period(tmptime, tmptimerange, &
1375 time_definition, pstart2, pstart1, reftime2)
1376 IF (lfull_steps) THEN
1378 useful = .true.
1379 ENDIF
1380 ELSE
1381 useful = .true.
1382 ENDIF
1383
1384 ENDIF
1385 ENDIF
1386 useful = useful .AND. tmptime /= datetime_miss .AND. &
1387 tmptimerange /= vol7d_timerange_miss .AND. tmptimerange%p2 == steps
1388
1389 IF (useful) THEN ! add a_otime, a_otimerange
1390 map_tr(i,j,k,l,1) = append_unique(a_otime, tmptime)
1391 map_tr(i,j,k,l,2) = append_unique(a_otimerange, tmptimerange)
1392 ENDIF
1393 ENDDO
1394 ENDDO
1395 ENDDO
1396 ENDDO
1397ENDDO
1398
1399! we have to repeat the computation with sorted arrays
1400CALL compute_keep_tr()
1401
1402otime => a_otime%array
1403otimerange => a_otimerange%array
1404! delete local objects keeping the contents
1407
1408#ifdef DEBUG
1409CALL l4f_log(l4f_debug, &
1414CALL l4f_log(l4f_debug, &
1417CALL l4f_log(l4f_debug, &
1419CALL l4f_log(l4f_debug, &
1421CALL l4f_log(l4f_debug, &
1423CALL l4f_log(l4f_debug, &
1425#endif
1426
1427CONTAINS
1428
1429SUBROUTINE compute_keep_tr()
1430
1431keep_tr(:,:,:) = imiss
1432DO l = 1, SIZE(itime)
1433 DO k = 1, nitr
1434 IF (itimerange(f(k))%p2 == steps) THEN
1435 CALL time_timerange_get_period(itime(l), itimerange(f(k)), &
1436 time_definition, pstart2, pend2, reftime2)
1437 useful = .false.
1438 IF (reftime2 == pend2) THEN ! analysis
1441 useful = .true.
1442 ENDIF
1443 ELSE IF (lfull_steps) THEN
1445 useful = .true.
1446 ENDIF
1447 ELSE
1448 useful = .true.
1449 ENDIF
1450 ELSE ! forecast
1451 IF (lfull_steps) THEN
1453 useful = .true.
1454 ENDIF
1455 ELSE
1456 useful = .true.
1457 ENDIF
1458 ENDIF
1459 IF (useful) THEN
1460! CALL time_timerange_set_period(tmptime, tmptimerange, &
1461! time_definition, pstart2, pend2, reftime2)
1462 keep_tr(k,l,1) = append_unique(a_otime, itime(l))
1463 keep_tr(k,l,2) = append_unique(a_otimerange, itimerange(f(k)))
1464 ENDIF
1465 ENDIF
1466 ENDDO
1467ENDDO
1468
1469END SUBROUTINE compute_keep_tr
1470
1471END SUBROUTINE recompute_stat_proc_diff_common
1472
1473
1474! common operations for statistical processing by metamorphosis
1475SUBROUTINE compute_stat_proc_metamorph_common(istat_proc, itimerange, ostat_proc, &
1476 otimerange, map_tr)
1477INTEGER,INTENT(in) :: istat_proc
1478TYPE(vol7d_timerange),INTENT(in) :: itimerange(:)
1479INTEGER,INTENT(in) :: ostat_proc
1480TYPE(vol7d_timerange),POINTER :: otimerange(:)
1481INTEGER,POINTER :: map_tr(:)
1482
1483INTEGER :: i
1484LOGICAL :: tr_mask(SIZE(itimerange))
1485
1486IF (SIZE(itimerange) == 0) THEN ! avoid segmentation fault in case of empty volume
1487 ALLOCATE(otimerange(0), map_tr(0))
1488 RETURN
1489ENDIF
1490
1491! useful timeranges
1492tr_mask(:) = itimerange(:)%timerange == istat_proc .AND. itimerange(:)%p2 /= imiss &
1493 .AND. itimerange(:)%p2 /= 0 ! .AND. itimerange(:)%p1 == 0
1494ALLOCATE(otimerange(count(tr_mask)), map_tr(count(tr_mask)))
1495
1496otimerange = pack(itimerange, mask=tr_mask)
1497otimerange(:)%timerange = ostat_proc
1498map_tr = pack((/(i,i=1,SIZE(itimerange))/), mask=tr_mask)
1499
1500END SUBROUTINE compute_stat_proc_metamorph_common
1501
1502
1503! common operations for statistical processing by aggregation
1504SUBROUTINE recompute_stat_proc_agg_common(itime, itimerange, stat_proc, tri, &
1505 step, time_definition, otime, otimerange, map_ttr, dtratio, start, full_steps)
1506TYPE(datetime),INTENT(in) :: itime(:)
1507TYPE(vol7d_timerange),INTENT(in) :: itimerange(:)
1508INTEGER,INTENT(in) :: stat_proc
1509INTEGER,INTENT(in) :: tri
1510TYPE(timedelta),INTENT(in) :: step
1511INTEGER,INTENT(in) :: time_definition
1512TYPE(datetime),POINTER :: otime(:)
1513TYPE(vol7d_timerange),POINTER :: otimerange(:)
1514TYPE(arrayof_ttr_mapper),POINTER :: map_ttr(:,:)
1515INTEGER,POINTER,OPTIONAL :: dtratio(:)
1516TYPE(datetime),INTENT(in),OPTIONAL :: start
1517LOGICAL,INTENT(in),OPTIONAL :: full_steps
1518
1519INTEGER :: i, j, k, l, na, nf, n
1520INTEGER :: steps, p1, maxp1, maxp2, minp1mp2, dstart, msteps
1521INTEGER(kind=int_ll) :: stepms, mstepms
1522LOGICAL :: lforecast
1523TYPE(datetime) :: lstart, lend, pstart1, pstart2, pend1, pend2, reftime1, reftime2, tmptime
1524TYPE(arrayof_datetime) :: a_otime
1525TYPE(arrayof_vol7d_timerange) :: a_otimerange
1526TYPE(arrayof_integer) :: a_dtratio
1527LOGICAL,ALLOCATABLE :: mask_timerange(:) ! improve !!!!
1528TYPE(ttr_mapper) :: lmapper
1529CHARACTER(len=8) :: env_var
1530LOGICAL :: climat_behavior
1531
1532
1533! enable bad behavior for climat database (only for instantaneous data)
1534env_var = ''
1535CALL getenv('LIBSIM_CLIMAT_BEHAVIOR', env_var)
1536climat_behavior = len_trim(env_var) > 0 .AND. .NOT.PRESENT(dtratio)
1537
1538! compute length of cumulation step in seconds
1540
1541! create a mask of suitable timeranges
1542ALLOCATE(mask_timerange(SIZE(itimerange)))
1543mask_timerange(:) = itimerange(:)%timerange == tri &
1544 .AND. itimerange(:)%p1 /= imiss .AND. itimerange(:)%p1 >= 0
1545
1546IF (PRESENT(dtratio)) THEN
1547 WHERE(itimerange(:)%p2 > 0 .AND. itimerange(:)%p2 /= imiss) ! for avoiding FPE MOD(steps, 0)
1548 mask_timerange(:) = mask_timerange(:) .AND. mod(steps, itimerange(:)%p2) == 0
1549 ELSEWHERE
1550 mask_timerange(:) = .false.
1551 END WHERE
1552ELSE ! instantaneous
1553 mask_timerange(:) = mask_timerange(:) .AND. itimerange(:)%p2 == 0
1554ENDIF
1555
1556#ifdef DEBUG
1557CALL l4f_log(l4f_debug, &
1558 '(re)compute_stat_proc_agg, number of useful timeranges before choosing analysis/forecast: '// &
1559 t2c(count(mask_timerange)))
1560#endif
1561
1562! euristically determine whether we are dealing with an
1563! analysis/observation or a forecast dataset
1564na = count(mask_timerange(:) .AND. itimerange(:)%p1 == 0)
1565nf = count(mask_timerange(:) .AND. itimerange(:)%p1 > 0)
1566lforecast = nf >= na
1567#ifdef DEBUG
1568CALL l4f_log(l4f_debug, &
1570#endif
1571
1572! keep only timeranges of one type (really necessary?)
1573IF (lforecast) THEN
1574 CALL l4f_log(l4f_info, &
1575 'recompute_stat_proc_agg, processing in forecast mode')
1576ELSE
1577 mask_timerange(:) = mask_timerange(:) .AND. itimerange(:)%p1 == 0
1578 CALL l4f_log(l4f_info, &
1579 'recompute_stat_proc_agg, processing in analysis mode')
1580ENDIF
1581
1582#ifdef DEBUG
1583CALL l4f_log(l4f_debug, &
1584 '(re)compute_stat_proc_agg, number of useful timeranges: '// &
1585 t2c(count(mask_timerange)))
1586#endif
1587
1588IF (SIZE(itime) == 0 .OR. count(mask_timerange) == 0) THEN ! avoid segmentation fault in case of empty volume
1589 ALLOCATE(otime(0), otimerange(0), map_ttr(0,0))
1590 IF (PRESENT(dtratio)) ALLOCATE(dtratio(0))
1591 RETURN
1592ENDIF
1593
1594! determine start and end of processing period, should work with p2==0
1595lstart = datetime_miss
1596IF (PRESENT(start)) lstart = start
1597lend = itime(SIZE(itime))
1598! compute some quantities
1599maxp1 = maxval(itimerange(:)%p1, mask=mask_timerange)
1600maxp2 = maxval(itimerange(:)%p2, mask=mask_timerange)
1601minp1mp2 = minval(itimerange(:)%p1 - itimerange(:)%p2, mask=mask_timerange)
1602IF (time_definition == 0) THEN ! reference time
1603 lend = lend + timedelta_new(sec=maxp1)
1604ENDIF
1605! extend interval at the end in order to include all the data in case
1606! frac_valid<1; must use < and not <= in the DO WHILE loops some lines
1607! below in order to exclude the last full interval which would be empty
1608lend = lend + step
1609IF (lstart == datetime_miss) THEN ! autodetect
1610 lstart = itime(1)
1611! if autodetected, adjust to obtain real absolute start of data
1612 IF (time_definition == 0) THEN ! reference time
1613 lstart = lstart + timedelta_new(sec=minp1mp2)
1614 ELSE ! verification time
1615! go back to start of longest processing interval
1616 lstart = lstart - timedelta_new(sec=maxp2)
1617 ENDIF
1618! full_steps is effective only in analysis mode and when start is not
1619! specified (start by itself in analysis mode implies full_steps with
1620! respect to start instead of absolute full steps)
1621 IF (optio_log(full_steps) .AND. .NOT.lforecast) THEN
1623 ENDIF
1624ENDIF
1625
1626#ifdef DEBUG
1627CALL l4f_log(l4f_debug, &
1629#endif
1630
1631! create output time and timerange lists
1632
1633IF (lforecast) THEN ! forecast mode
1634 IF (time_definition == 0) THEN ! reference time
1636
1637! apply start shift to timerange, not to reftime
1638! why did we use itime(SIZE(itime)) (last ref time)?
1639! CALL getval(lstart-itime(SIZE(itime)), asec=dstart)
1641! allow starting before first reftime but restrict dtstart to -steps
1642! dstart = MAX(0, dstart)
1644 DO p1 = steps + dstart, maxp1, steps
1645 CALL insert_unique(a_otimerange, vol7d_timerange_new(stat_proc, p1, steps))
1646 ENDDO
1647
1648 ELSE ! verification time
1649
1650! verification time in forecast mode is the ugliest case, because we
1651! don't know a priori how many different (thus incompatible) reference
1652! times we have, so some assumption of regularity has to be made. For
1653! this reason msteps, the minimum step between two times, is
1654! computed. We choose to compute it as a difference between itime
1655! elements, it could be also computed as difference of itimerange%p1
1656! elements. But what if it is not modulo steps?
1657 mstepms = steps*1000_int_ll
1658 DO i = 2, SIZE(itime)
1660 IF (stepms > 0_int_ll .AND. stepms < mstepms) THEN
1661 msteps = stepms/1000_int_ll
1663 ENDIF
1664 ENDDO
1665 msteps = mstepms/1000_int_ll
1666
1667 tmptime = lstart + step
1668 DO WHILE(tmptime < lend) ! < since lend has been extended
1669 CALL insert_unique(a_otime, tmptime)
1670 tmptime = tmptime + step
1671 ENDDO
1672
1673! in next loop, we used initially steps instead of msteps (see comment
1674! above), this gave much less combinations of time/timeranges with
1675! possible empty output; we start from msteps instead of steps in
1676! order to include partial initial intervals in case frac_valid<1;
1677! however, as a gemeral rule, for aggregation of forecasts it's better
1678! to use reference time
1679 DO p1 = msteps, maxp1, msteps
1680 CALL insert_unique(a_otimerange, vol7d_timerange_new(stat_proc, p1, steps))
1681 ENDDO
1682
1683 ENDIF
1684
1685ELSE ! analysis mode
1686 tmptime = lstart + step
1687 DO WHILE(tmptime < lend) ! < since lend has been extended
1688 CALL insert_unique(a_otime, tmptime)
1689 tmptime = tmptime + step
1690 ENDDO
1691 CALL insert_unique(a_otimerange, vol7d_timerange_new(stat_proc, 0, steps))
1692
1693ENDIF
1694
1697otime => a_otime%array
1698otimerange => a_otimerange%array
1701! delete local objects keeping the contents
1704
1705#ifdef DEBUG
1706CALL l4f_log(l4f_debug, &
1707 'recompute_stat_proc_agg, output time and timerange: '//&
1709#endif
1710
1711IF (PRESENT(dtratio)) THEN
1712! count the possible i/o interval ratios
1713 DO k = 1, SIZE(itimerange)
1714 IF (itimerange(k)%p2 /= 0) &
1715 CALL insert_unique(a_dtratio, steps/itimerange(k)%p2) ! guaranteed to be integer
1716 ENDDO
1718 dtratio => a_dtratio%array
1720! delete local object keeping the contents
1722
1723#ifdef DEBUG
1724 CALL l4f_log(l4f_debug, &
1726 ' possible aggregation ratios, from '// &
1728#endif
1729
1730 ALLOCATE(map_ttr(SIZE(otime),SIZE(otimerange)))
1731 do_itimerange1: DO l = 1, SIZE(itimerange)
1732 IF (.NOT.mask_timerange(l)) cycle do_itimerange1
1733 do_itime1: DO k = 1, SIZE(itime)
1734 CALL time_timerange_get_period(itime(k), itimerange(l), &
1735 time_definition, pstart1, pend1, reftime1)
1736 do_otimerange1: DO j = 1, SIZE(otimerange)
1737 do_otime1: DO i = 1, SIZE(otime)
1738 CALL time_timerange_get_period_pop(otime(i), otimerange(j), step, &
1739 time_definition, pstart2, pend2, reftime2)
1740 IF (lforecast) THEN
1741 IF (reftime1 /= reftime2) cycle do_otime1
1742 ENDIF
1743
1744 IF (pstart1 >= pstart2 .AND. pend1 <= pend2 .AND. &
1746 lmapper%it = k
1747 lmapper%itr = l
1748 lmapper%extra_info = steps/itimerange(l)%p2 ! dtratio, guaranteed to be integer
1749 n = append(map_ttr(i,j), lmapper)
1750 cycle do_itime1 ! can contribute only to a single interval
1751 ENDIF
1752 ENDDO do_otime1
1753 ENDDO do_otimerange1
1754 ENDDO do_itime1
1755 ENDDO do_itimerange1
1756
1757ELSE
1758
1759 ALLOCATE(map_ttr(SIZE(otime),SIZE(otimerange)))
1760 do_itimerange2: DO l = 1, SIZE(itimerange)
1761 IF (.NOT.mask_timerange(l)) cycle do_itimerange2
1762 do_itime2: DO k = 1, SIZE(itime)
1763 CALL time_timerange_get_period(itime(k), itimerange(l), &
1764 time_definition, pstart1, pend1, reftime1)
1765 do_otimerange2: DO j = 1, SIZE(otimerange)
1766 do_otime2: DO i = 1, SIZE(otime)
1767 CALL time_timerange_get_period_pop(otime(i), otimerange(j), step, &
1768 time_definition, pstart2, pend2, reftime2)
1769 IF (lforecast) THEN
1770 IF (reftime1 /= reftime2) cycle do_otime2
1771 ENDIF
1772
1773 IF (climat_behavior .AND. pstart1 == pstart2) cycle do_otime2
1774 IF (pstart1 >= pstart2 .AND. pend1 <= pend2) THEN ! useful
1775 lmapper%it = k
1776 lmapper%itr = l
1777 IF (pstart1 == pstart2) THEN
1778 lmapper%extra_info = 1 ! start of interval
1779 ELSE IF (pend1 == pend2) THEN
1780 lmapper%extra_info = 2 ! end of interval
1781 ELSE
1782 lmapper%extra_info = imiss
1783 ENDIF
1784 lmapper%time = pstart1 ! = pend1, order by time?
1785 n = insert_sorted(map_ttr(i,j), lmapper, .true., .true.)
1786! no CYCLE, a single input can contribute to multiple output intervals
1787 ENDIF
1788 ENDDO do_otime2
1789 ENDDO do_otimerange2
1790 ENDDO do_itime2
1791 ENDDO do_itimerange2
1792
1793ENDIF
1794
1795END SUBROUTINE recompute_stat_proc_agg_common
1796
1797
1798SUBROUTINE compute_stat_proc_agg_sw(vertime, pstart, pend, time_mask, &
1799 max_step, weights)
1800TYPE(datetime),INTENT(in) :: vertime(:)
1801TYPE(datetime),INTENT(in) :: pstart
1802TYPE(datetime),INTENT(in) :: pend
1803LOGICAL,INTENT(in) :: time_mask(:)
1804TYPE(timedelta),OPTIONAL,INTENT(out) :: max_step
1805DOUBLE PRECISION,OPTIONAL,INTENT(out) :: weights(:)
1806
1807INTEGER :: i, nt
1808TYPE(datetime),ALLOCATABLE :: lvertime(:)
1809TYPE(datetime) :: half, nexthalf
1810INTEGER(kind=int_ll) :: dt, tdt
1811
1812nt = count(time_mask)
1813ALLOCATE(lvertime(nt))
1814lvertime = pack(vertime, mask=time_mask)
1815
1816IF (PRESENT(max_step)) THEN
1817! new way
1818! max_step = timedelta_0
1819! DO i = 1, nt - 1
1820! IF (lvertime(i+1) - lvertime(i) > max_step) &
1821! max_step = lvertime(i+1) - lvertime(i)
1822! ENDDO
1823! IF (lvertime(1) - pstart > max_step) max_step = lvertime(1) - pstart
1824! IF (pend - lvertime(nt) > max_step) max_step = pend - lvertime(nt)
1825! old way
1826 IF (nt == 1) THEN
1827 max_step = pend - pstart
1828 ELSE
1829 half = lvertime(1) + (lvertime(2) - lvertime(1))/2
1830 max_step = half - pstart
1831 DO i = 2, nt - 1
1832 nexthalf = lvertime(i) + (lvertime(i+1) - lvertime(i))/2
1833 IF (nexthalf - half > max_step) max_step = nexthalf - half
1834 half = nexthalf
1835 ENDDO
1836 IF (pend - half > max_step) max_step = pend - half
1837 ENDIF
1838
1839ENDIF
1840
1841IF (PRESENT(weights)) THEN
1842 IF (nt == 1) THEN
1843 weights(1) = 1.0
1844 ELSE
1846 half = lvertime(1) + (lvertime(2) - lvertime(1))/2
1848 weights(1) = dble(dt)/dble(tdt)
1849 DO i = 2, nt - 1
1850 nexthalf = lvertime(i) + (lvertime(i+1) - lvertime(i))/2
1852 weights(i) = dble(dt)/dble(tdt)
1853 half = nexthalf
1854 ENDDO
1856 weights(nt) = dble(dt)/dble(tdt)
1857 ENDIF
1858ENDIF
1859
1860END SUBROUTINE compute_stat_proc_agg_sw
1861
1862! get start of period, end of period and reference time from time,
1863! timerange, according to time_definition.
1864SUBROUTINE time_timerange_get_period(time, timerange, time_definition, &
1865 pstart, pend, reftime)
1866TYPE(datetime),INTENT(in) :: time
1867TYPE(vol7d_timerange),INTENT(in) :: timerange
1868INTEGER,INTENT(in) :: time_definition
1869TYPE(datetime),INTENT(out) :: reftime
1870TYPE(datetime),INTENT(out) :: pstart
1871TYPE(datetime),INTENT(out) :: pend
1872
1873TYPE(timedelta) :: p1, p2
1874
1875
1876p1 = timedelta_new(sec=timerange%p1) ! end of period
1877p2 = timedelta_new(sec=timerange%p2) ! length of period
1878
1880! (timerange%p1 > 0 .AND. timerange%p1 < timerange%p2) .OR. &
1881 timerange%p1 < 0 .OR. timerange%p2 < 0) THEN ! is this too pedantic and slow?
1882 pstart = datetime_miss
1883 pend = datetime_miss
1884 reftime = datetime_miss
1885 RETURN
1886ENDIF
1887
1888IF (time_definition == 0) THEN ! time == reference time
1889 reftime = time
1890 pend = time + p1
1891 pstart = pend - p2
1892ELSE IF (time_definition == 1) THEN ! time == verification time
1893 pend = time
1894 pstart = time - p2
1895 reftime = time - p1
1896ELSE
1897 pstart = datetime_miss
1898 pend = datetime_miss
1899 reftime = datetime_miss
1900ENDIF
1901
1902END SUBROUTINE time_timerange_get_period
1903
1904
1905! get start of period, end of period and reference time from time,
1906! timerange, according to time_definition. step is taken separately
1907! from timerange and can be "popular"
1908SUBROUTINE time_timerange_get_period_pop(time, timerange, step, time_definition, &
1909 pstart, pend, reftime)
1910TYPE(datetime),INTENT(in) :: time
1911TYPE(vol7d_timerange),INTENT(in) :: timerange
1912TYPE(timedelta),INTENT(in) :: step
1913INTEGER,INTENT(in) :: time_definition
1914TYPE(datetime),INTENT(out) :: reftime
1915TYPE(datetime),INTENT(out) :: pstart
1916TYPE(datetime),INTENT(out) :: pend
1917
1918TYPE(timedelta) :: p1
1919
1920
1921p1 = timedelta_new(sec=timerange%p1) ! end of period
1922
1924! (timerange%p1 > 0 .AND. timerange%p1 < timerange%p2) .OR. &
1925 timerange%p1 < 0 .OR. timerange%p2 < 0) THEN ! is this too pedantic and slow?
1926 pstart = datetime_miss
1927 pend = datetime_miss
1928 reftime = datetime_miss
1929 RETURN
1930ENDIF
1931
1932IF (time_definition == 0) THEN ! time == reference time
1933 reftime = time
1934 pend = time + p1
1935 pstart = pend - step
1936ELSE IF (time_definition == 1) THEN ! time == verification time
1937 pend = time
1938 pstart = time - step
1939 reftime = time - p1
1940ELSE
1941 pstart = datetime_miss
1942 pend = datetime_miss
1943 reftime = datetime_miss
1944ENDIF
1945
1946END SUBROUTINE time_timerange_get_period_pop
1947
1948
1949! set time, timerange%p1, timerange%p2 according to pstart, pend,
1950! reftime and time_definition.
1951SUBROUTINE time_timerange_set_period(time, timerange, time_definition, &
1952 pstart, pend, reftime)
1953TYPE(datetime),INTENT(out) :: time
1954TYPE(vol7d_timerange),INTENT(inout) :: timerange
1955INTEGER,INTENT(in) :: time_definition
1956TYPE(datetime),INTENT(in) :: reftime
1957TYPE(datetime),INTENT(in) :: pstart
1958TYPE(datetime),INTENT(in) :: pend
1959
1960TYPE(timedelta) :: p1, p2
1961INTEGER(kind=int_ll) :: dmsec
1962
1963
1964IF (time_definition == 0) THEN ! time == reference time
1965 time = reftime
1966 p1 = pend - reftime
1967 p2 = pend - pstart
1968ELSE IF (time_definition == 1) THEN ! time == verification time
1969 time = pend
1970 p1 = pend - reftime
1971 p2 = pend - pstart
1972ELSE
1973 time = datetime_miss
1974ENDIF
1975
1976IF (time /= datetime_miss) THEN
1978 timerange%p1 = int(dmsec/1000_int_ll)
1980 timerange%p2 = int(dmsec/1000_int_ll)
1981ELSE
1982 timerange%p1 = imiss
1983 timerange%p2 = imiss
1984ENDIF
1985
1986END SUBROUTINE time_timerange_set_period
1987
1988
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 |