libsim Versione 7.2.6
|
◆ volgrid6d_compute_vert_coord_var()
Method for building a volume containing the vertical coordinate as a variable. This method produces a volgrid6d volume, derived from this, containing a single variable, horizontally constant, on the same input levels, which describes the vertical coordinate in the form of a physical variable. The grid, time and timerange metadata are the same as for the original volume. Only a single vertical level type, the one matching the level argument, is converted to a variable. The level argument can also indicate the layer between two surfaces of the same type, in that case the variable representing the vertical coordinate will be set to the value of the midpoint between the two layers. If something goes wrong, e.g. no level matches level argument or the level canot be converted to a physical value, an empty volume is returned.
Definizione alla linea 1038 del file volgrid6d_class_compute.F90. 1039! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
1040! authors:
1041! Davide Cesari <dcesari@arpa.emr.it>
1042! Paolo Patruno <ppatruno@arpa.emr.it>
1043
1044! This program is free software; you can redistribute it and/or
1045! modify it under the terms of the GNU General Public License as
1046! published by the Free Software Foundation; either version 2 of
1047! the License, or (at your option) any later version.
1048
1049! This program is distributed in the hope that it will be useful,
1050! but WITHOUT ANY WARRANTY; without even the implied warranty of
1051! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
1052! GNU General Public License for more details.
1053
1054! You should have received a copy of the GNU General Public License
1055! along with this program. If not, see <http://www.gnu.org/licenses/>.
1056#include "config.h"
1057!> Extension of volgrid6d_class with methods for performing simple
1058!! statistical operations on entire volumes of data. This module
1059!! includes methods for performing operations such as (re)averaging
1060!! and (re)cumulation in time of entire volgrid6d_class::volgrid6d
1061!! objects.
1062!!
1063!! \ingroup volgrid6d
1071IMPLICIT NONE
1072
1073CONTAINS
1074
1075!> General-purpose method for computing a statistical processing on
1076!! data in a volgrid6d object already processed with the same statistical
1077!! processing, on a different time interval specified by \a step and
1078!! \a start. This method tries to apply all the suitable specialized
1079!! statistical processing methods according to the input and output
1080!! statistical processing requested. The argument \a stat_proc_input
1081!! determines which data will be processed, while the \a stat_proc
1082!! argument determines the type of statistical process to be applied
1083!! and which will be owned by output data.
1084!!
1085!! The possible combinations are:
1086!!
1087!! - \a stat_proc_input = 254
1088!! - \a stat_proc = 0 average instantaneous observations
1089!! - \a stat_proc = 2 compute maximum of instantaneous observations
1090!! - \a stat_proc = 3 compute minimum of instantaneous observations
1091!! - \a stat_proc = 4 compute difference of instantaneous observations
1092!!
1093!! processing is computed on longer time intervals by aggregation,
1094!! see the description of volgrid6d_compute_stat_proc_agg()
1095!!
1096!! - \a stat_proc_input = *
1097!! - \a stat_proc = 254 consider statistically processed values as
1098!! instantaneous without any extra processing
1099!!
1100!! see the description of volgrid6d_decompute_stat_proc()
1101!!
1102!! - \a stat_proc_input = 0, 1, 2, 3, 4, 200
1103!! - \a stat_proc = \a stat_proc_input recompute input data on
1104!! different intervals
1105!!
1106!! the same statistical processing is applied to obtain data
1107!! processed on a different interval, either longer, by
1108!! aggregation, or shorter, by differences, see the description of
1109!! volgrid6d_recompute_stat_proc_agg() and
1110!! volgrid6d_recompute_stat_proc_diff() respectively; it is also
1111!! possible to provide \a stat_proc_input \a /= \a stat_proc, but
1112!! it has to be used with care.
1113!!
1114!! - \a stat_proc_input = 0
1115!! - \a stat_proc = 1
1116!!
1117!! a time-averaged rate or flux is transformed into a
1118!! time-integrated value (sometimes called accumulated) on the same
1119!! interval by multiplying the values by the length of the time
1120!! interval in seconds, keeping constant all the rest, including
1121!! the variable; the unit of the variable implicitly changes
1122!! accordingly, this is supported officially in grib2 standard, in
1123!! the other cases it is a forcing of the standards.
1124!!
1125!! - \a stat_proc_input = 1
1126!! - \a stat_proc = 0
1127!!
1128!! a time-integrated value (sometimes called accumulated) is
1129!! transformed into a time-averaged rate or flux on the same
1130!! interval by dividing the values by the length of the time
1131!! interval in seconds, see also the previous description of the
1132!! opposite computation.
1133!!
1134!! If a particular statistical processing cannot be performed on the
1135!! input data, the program continues with a warning and, if requested,
1136!! the input data is passed over to the volume specified by the \a
1137!! other argument, in order to allow continuation of processing. All
1138!! the other parameters are passed over to the specifical statistical
1139!! processing methods and are documented there.
1140SUBROUTINE volgrid6d_compute_stat_proc(this, that, stat_proc_input, stat_proc, &
1141 step, start, full_steps, frac_valid, max_step, weighted, clone)
1142TYPE(volgrid6d),INTENT(inout) :: this !< volume providing data to be recomputed, it is not modified by the method, apart from performing a \a volgrid6d_alloc_vol on it
1143TYPE(volgrid6d),INTENT(out) :: that !< output volume which will contain the recomputed data
1144INTEGER,INTENT(in) :: stat_proc_input !< type of statistical processing of data that has to be processed (from grib2 table), only data having timerange of this type will be processed, the actual statistical processing performed and which will appear in the output volume, is however determined by \a stat_proc argument
1145INTEGER,INTENT(in) :: stat_proc !< type of statistical processing to be recomputed (from grib2 table), data in output volume \a that will have a timerange of this type
1146TYPE(timedelta),INTENT(in) :: step !< length of the step over which the statistical processing is performed
1147TYPE(datetime),INTENT(in),OPTIONAL :: start !< start of statistical processing interval
1148LOGICAL,INTENT(in),OPTIONAL :: full_steps !< if \a .TRUE. apply processing only on intervals starting at a forecast time or a reference time modulo \a step
1149REAL,INTENT(in),OPTIONAL :: frac_valid !< minimum fraction of valid data required for considering acceptable a recomputed value, default=1.
1150TYPE(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
1151LOGICAL,INTENT(in),OPTIONAL :: weighted !< if provided and \c .TRUE., the statistical process is computed, if possible, by weighting every value with a weight proportional to its validity interval
1152LOGICAL , INTENT(in),OPTIONAL :: clone !< if provided and \c .TRUE. , clone the gaid's from \a this to \a that
1153
1154INTEGER :: dtmax, dtstep
1155
1156
1157IF (stat_proc_input == 254) THEN
1158 CALL l4f_category_log(this%category, l4f_info, &
1159 'computing statistical processing by aggregation '//&
1161
1162 CALL volgrid6d_compute_stat_proc_agg(this, that, stat_proc, &
1163 step, start, full_steps, max_step, clone)
1164
1165ELSE IF (stat_proc == 254) THEN
1166 CALL l4f_category_log(this%category, l4f_error, &
1167 'statistical processing to instantaneous data not implemented for gridded fields')
1168 CALL raise_error()
1169
1170ELSE IF (stat_proc_input == stat_proc .OR. &
1171 (stat_proc == 0 .OR. stat_proc == 2 .OR. stat_proc == 3)) THEN
1172! avg, min and max can be computed from any input, with care
1173
1174 IF (count(this%timerange(:)%timerange == stat_proc_input) == 0) THEN
1175 CALL l4f_category_log(this%category, l4f_warn, &
1177! return an empty volume, without signaling error
1179 CALL volgrid6d_alloc_vol(that)
1180
1181 ELSE
1182! euristically determine whether aggregation or difference is more suitable
1183 dtmax = maxval(this%timerange(:)%p2, &
1184 mask=(this%timerange(:)%timerange == stat_proc))
1186
1187#ifdef DEBUG
1188 CALL l4f_category_log(this%category, l4f_debug, &
1190#endif
1191
1192 IF (dtstep <= dtmax) THEN
1193 CALL l4f_category_log(this%category, l4f_info, &
1194 'recomputing statistically processed data by difference '// &
1196 CALL volgrid6d_recompute_stat_proc_diff(this, that, stat_proc, step, &
1197 full_steps, start, clone)
1198 ELSE
1199 CALL l4f_category_log(this%category, l4f_info, &
1200 'recomputing statistically processed data by aggregation '// &
1202 CALL volgrid6d_recompute_stat_proc_agg(this, that, stat_proc, step, start, &
1203 full_steps, frac_valid, clone, stat_proc_input)
1204 ENDIF
1205 ENDIF
1206
1207ELSE ! IF (stat_proc_input /= stat_proc) THEN
1208 IF ((stat_proc_input == 0 .AND. stat_proc == 1) .OR. &
1209 (stat_proc_input == 1 .AND. stat_proc == 0)) THEN
1210 CALL l4f_category_log(this%category, l4f_info, &
1211 'computing statistically processed data by integration/differentiation '// &
1213 CALL volgrid6d_compute_stat_proc_metamorph(this, that, stat_proc_input, &
1214 stat_proc, clone)
1215 ELSE
1216 CALL l4f_category_log(this%category, l4f_error, &
1218 ' not implemented or does not make sense')
1219 CALL raise_error()
1220 ENDIF
1221
1222ENDIF
1223
1224END SUBROUTINE volgrid6d_compute_stat_proc
1225
1226
1227!> Specialized method for statistically processing a set of data
1228!! already processed with the same statistical processing, on a
1229!! different time interval. This method performs statistical
1230!! processing by aggregation of shorter intervals.
1231!!
1232!! The output \a that volgrid6d object contains elements from the original volume
1233!! \a this satisfying the conditions
1234!! - timerange (vol7d_timerange_class::vol7d_timerange::timerange)
1235!! of type \a stat_proc (or \a stat_proc_input if provided)
1236!! - any p1 (analysis/observation or forecast)
1237!! - p2 > 0 (processing interval non null, non instantaneous data)
1238!! and equal to a multiplier of \a step
1239!!
1240!! Output data will have timerange of type \a stat_proc and p2 = \a
1241!! step. The supported statistical processing methods (parameter \a
1242!! stat_proc) are:
1243!!
1244!! - 0 average
1245!! - 1 accumulation
1246!! - 2 maximum
1247!! - 3 minimum
1248!! - 4 difference
1249!! - 200 vectorial mean
1250!!
1251!! The start of processing period can be computed automatically from
1252!! the input intervals as the first possible interval modulo \a step,
1253!! or, for a better control, it can be specified explicitely by the
1254!! optional argument \a start. Be warned that, in the final volume,
1255!! the first reference time will actually be \a start \a + \a step,
1256!! since \a start indicates the beginning of first processing
1257!! interval, while reference time (for analysis/oservation) is the end
1258!! of the interval.
1259!!
1260!! The purpose of the optional argument \a stat_proc_input is to allow
1261!! processing with a certain statistical processing operator a dataset
1262!! already processed with a different operator, by specifying the
1263!! latter as stat_proc_input; this is useful, for example, if one
1264!! wants to compute the monthly average of daily maximum temperatures;
1265!! however this has to be used with care since the resulting data
1266!! volume will not carry all the information about the processing
1267!! which has been done, in the previous case, for example, the
1268!! temperatures will simply look like monthly average temperatures.
1269SUBROUTINE volgrid6d_recompute_stat_proc_agg(this, that, stat_proc, &
1270 step, start, full_steps, frac_valid, clone, stat_proc_input)
1271TYPE(volgrid6d),INTENT(inout) :: this !< volume providing data to be recomputed, it is not modified by the method, apart from performing a \a volgrid6d_alloc_vol on it
1272TYPE(volgrid6d),INTENT(out) :: that !< output volume which will contain the recomputed data
1273INTEGER,INTENT(in) :: stat_proc !< type of statistical processing to be recomputed (from grib2 table), only data having timerange of this type will be recomputed and will appear in the output volume
1274TYPE(timedelta),INTENT(in) :: step !< length of the step over which the statistical processing is performed
1275TYPE(datetime),INTENT(in),OPTIONAL :: start !< start of statistical processing interval
1276LOGICAL,INTENT(in),OPTIONAL :: full_steps !< if \a .TRUE. and \a start is not provided, apply processing only on intervals starting at a forecast time or a reference time modulo \a step
1277REAL,INTENT(in),OPTIONAL :: frac_valid !< minimum fraction of valid data required for considering acceptable a recomputed value, default=1.
1278LOGICAL, INTENT(in),OPTIONAL :: clone !< if provided and \c .TRUE. , clone the gaid's from \a this to \a that
1279INTEGER,INTENT(in),OPTIONAL :: stat_proc_input !< to be used with care, type of statistical processing of data that has to be processed (from grib2 table), only data having timerange of this type will be recomputed, the actual statistical processing performed and which will appear in the output volume, is however determined by \a stat_proc argument
1280
1281INTEGER :: tri
1282INTEGER i, j, n, n1, ndtr, i3, i6
1283TYPE(arrayof_ttr_mapper),POINTER :: map_ttr(:,:)
1284INTEGER,POINTER :: dtratio(:)
1285REAL :: lfrac_valid
1286LOGICAL :: lclone
1287REAL,POINTER :: voldatiin(:,:), voldatiout(:,:)
1288
1289
1290NULLIFY(voldatiin, voldatiout)
1291IF (PRESENT(stat_proc_input)) THEN
1292 tri = stat_proc_input
1293ELSE
1294 tri = stat_proc
1295ENDIF
1296IF (PRESENT(frac_valid)) THEN
1297 lfrac_valid = frac_valid
1298ELSE
1299 lfrac_valid = 1.0
1300ENDIF
1301
1303! be safe
1304CALL volgrid6d_alloc_vol(this)
1305
1306! when volume is not decoded it is better to clone anyway to avoid
1307! overwriting fields
1308lclone = optio_log(clone) .OR. .NOT.ASSOCIATED(this%voldati)
1309! initialise the output volume
1311CALL volgrid6d_alloc(that, dim=this%griddim%dim, ntimerange=1, &
1312 nlevel=SIZE(this%level), nvar=SIZE(this%var), ini=.false.)
1313that%level = this%level
1314that%var = this%var
1315
1316CALL recompute_stat_proc_agg_common(this%time, this%timerange, stat_proc, tri, &
1317 step, this%time_definition, that%time, that%timerange, map_ttr, &
1318 dtratio=dtratio, start=start, full_steps=full_steps)
1319
1320CALL volgrid6d_alloc_vol(that, decode=ASSOCIATED(this%voldati))
1321
1322do_otimerange: DO j = 1, SIZE(that%timerange)
1323 do_otime: DO i = 1, SIZE(that%time)
1324
1325 DO n1 = 1, SIZE(dtratio)
1326 IF (dtratio(n1) <= 0) cycle ! safety check
1327
1328 DO i6 = 1, SIZE(this%var)
1329 DO i3 = 1, SIZE(this%level)
1330 CALL volgrid_get_vol_2d(that, i3, i, j, i6, voldatiout)
1331 ndtr = 0
1332 DO n = 1, map_ttr(i,j)%arraysize
1333 IF (map_ttr(i,j)%array(n)%extra_info == dtratio(n1)) THEN
1334 ndtr = ndtr + 1
1335 CALL volgrid_get_vol_2d(this, i3, map_ttr(i,j)%array(n)%it, &
1336 map_ttr(i,j)%array(n)%itr, i6, voldatiin)
1337
1338 IF (ndtr == 1) THEN
1339 voldatiout = voldatiin
1340 IF (lclone) THEN
1342 map_ttr(i,j)%array(n)%itr,i6), that%gaid(i3,i,j,i6))
1343 ELSE
1344 that%gaid(i3,i,j,i6) = this%gaid(i3, map_ttr(i,j)%array(n)%it, &
1345 map_ttr(i,j)%array(n)%itr,i6)
1346 ENDIF
1347
1348 ELSE ! second or more time
1349 SELECT CASE(stat_proc)
1350 CASE (0, 200, 1, 4) ! average, vectorial mean, accumulation, difference
1352 voldatiout(:,:) = voldatiout(:,:) + voldatiin(:,:)
1353 ELSEWHERE
1354 voldatiout(:,:) = rmiss
1355 END WHERE
1356 CASE(2) ! maximum
1358 voldatiout(:,:) = max(voldatiout(:,:), voldatiin(:,:))
1359 ELSEWHERE
1360 voldatiout(:,:) = rmiss
1361 END WHERE
1362 CASE(3) ! minimum
1364 voldatiout(:,:) = min(voldatiout(:,:), voldatiin(:,:))
1365 ELSEWHERE
1366 voldatiout(:,:) = rmiss
1367 END WHERE
1368 END SELECT
1369
1370 ENDIF ! first time
1371 ENDIF ! dtratio(n1)
1372 ENDDO ! ttr
1373
1374#ifdef DEBUG
1375 CALL l4f_log(l4f_debug, &
1376 'compute_stat_proc_agg, ndtr/dtratio/frac_valid: '// &
1378#endif
1379 IF (ndtr > 0) THEN ! why this condition was not here before?
1380 IF (real(ndtr)/real(dtratio(n1)) >= lfrac_valid) THEN ! success
1381 IF (stat_proc == 0) THEN ! average
1383 voldatiout(:,:) = voldatiout(:,:)/ndtr
1384 END WHERE
1385 ENDIF
1386 CALL volgrid_set_vol_2d(that, i3, i, j, i6, voldatiout)
1387#ifdef DEBUG
1388 CALL l4f_log(l4f_debug, &
1389 'compute_stat_proc_agg, coding lev/t/tr/var: '// &
1391#endif
1392 ELSE
1393! must nullify the output gaid here, otherwise an incomplete field will be output
1394 IF (lclone) THEN
1396 ELSE
1398 ENDIF
1399#ifdef DEBUG
1400 CALL l4f_log(l4f_debug, &
1401 'compute_stat_proc_agg, skipping lev/t/tr/var: '// &
1403#endif
1404 ENDIF
1405 ENDIF ! ndtr > 0
1406
1407 ENDDO ! level
1408 ENDDO ! var
1409 ENDDO ! dtratio
1411 ENDDO do_otime
1412ENDDO do_otimerange
1413
1414DEALLOCATE(dtratio, map_ttr)
1415
1416END SUBROUTINE volgrid6d_recompute_stat_proc_agg
1417
1418
1419!> Method for statistically processing a set of instantaneous data.
1420!! This method performs statistical processing by aggregation of
1421!! instantaneous data.
1422!!
1423!! The output \a that volgrid6d object contains elements from the original volume
1424!! \a this satisfying the conditions
1425!! - timerange (vol7d_timerange_class::vol7d_timerange::timerange)
1426!! of type 254 (instantaeous)
1427!! - any p1 (analysis/observation or forecast)
1428!! - p2 = 0 (processing interval null, instantaneous data)
1429!!
1430!! Output data will have timerange of type \a stat_proc, and p2 = \a
1431!! step. The supported statistical processing methods (parameter \a
1432!! stat_proc) are:
1433!!
1434!! - 0 average
1435!! - 2 maximum
1436!! - 3 minimum
1437!! - 4 difference
1438!!
1439!! A maximum distance in time for input valid data can be assigned
1440!! with the optional argument \a max_step, in order to filter datasets
1441!! with too long "holes".
1442SUBROUTINE volgrid6d_compute_stat_proc_agg(this, that, stat_proc, &
1443 step, start, full_steps, max_step, clone)
1444TYPE(volgrid6d),INTENT(inout) :: this !< volume providing data to be recomputed, it is not modified by the method, apart from performing a \a volgrid6d_alloc_vol on it
1445TYPE(volgrid6d),INTENT(out) :: that !< output volume which will contain the recomputed data
1446INTEGER,INTENT(in) :: stat_proc !< type of statistical processing to be recomputed (from grib2 table), only data having timerange of this type will be recomputed and will appear in the output volume
1447TYPE(timedelta),INTENT(in) :: step !< length of the step over which the statistical processing is performed
1448TYPE(datetime),INTENT(in),OPTIONAL :: start !< start of statistical processing interval
1449LOGICAL,INTENT(in),OPTIONAL :: full_steps !< if \a .TRUE. and \a start is not provided, apply processing only on intervals starting at a forecast time or a reference time modulo \a step
1450TYPE(timedelta),INTENT(in),OPTIONAL :: max_step !< maximum allowed distance in time between two contiguougs valid data within an interval, for the interval to be eligible for statistical processing
1451LOGICAL , INTENT(in),OPTIONAL :: clone !< if provided and \c .TRUE. , clone the gaid's from \a this to \a that
1452
1453INTEGER :: tri
1454INTEGER i, j, n, ninp, i3, i6
1455TYPE(arrayof_ttr_mapper),POINTER :: map_ttr(:,:)
1456TYPE(timedelta) :: lmax_step
1457LOGICAL :: lclone
1458REAL,POINTER :: voldatiin(:,:), voldatiout(:,:)
1459
1460
1461NULLIFY(voldatiin, voldatiout)
1462tri = 254
1463IF (PRESENT(max_step)) THEN
1464 lmax_step = max_step
1465ELSE
1466 lmax_step = timedelta_max
1467ENDIF
1468
1470! be safe
1471CALL volgrid6d_alloc_vol(this)
1472
1473! when volume is not decoded it is better to clone anyway to avoid
1474! overwriting fields
1475lclone = optio_log(clone) .OR. .NOT.ASSOCIATED(this%voldati)
1476! initialise the output volume
1478CALL volgrid6d_alloc(that, dim=this%griddim%dim, ntimerange=1, &
1479 nlevel=SIZE(this%level), nvar=SIZE(this%var), ini=.false.)
1480that%level = this%level
1481that%var = this%var
1482
1483CALL recompute_stat_proc_agg_common(this%time, this%timerange, stat_proc, tri, &
1484 step, this%time_definition, that%time, that%timerange, map_ttr, &
1485 start=start, full_steps=full_steps)
1486
1487CALL volgrid6d_alloc_vol(that, decode=ASSOCIATED(this%voldati))
1488
1489do_otimerange: DO j = 1, SIZE(that%timerange)
1490 do_otime: DO i = 1, SIZE(that%time)
1491 ninp = map_ttr(i,j)%arraysize
1492 IF (ninp <= 0) cycle do_otime
1493
1494 IF (stat_proc == 4) THEN ! check validity for difference
1495 IF (map_ttr(i,j)%array(1)%extra_info /= 1 .OR. &
1496 map_ttr(i,j)%array(ninp)%extra_info /= 2) THEN
1498 cycle do_otime
1499 ENDIF
1500 ELSE
1501! check validity condition (missing values in volume are not accounted for)
1502 DO n = 2, ninp
1503 IF (map_ttr(i,j)%array(n)%time - map_ttr(i,j)%array(n-1)%time > &
1504 lmax_step) THEN
1506 cycle do_otime
1507 ENDIF
1508 ENDDO
1509 ENDIF
1510
1511 DO i6 = 1, SIZE(this%var)
1512 DO i3 = 1, SIZE(this%level)
1513 CALL volgrid_get_vol_2d(that, i3, i, j, i6, voldatiout)
1514
1515 IF (stat_proc == 4) THEN ! special treatment for difference
1516 IF (lclone) THEN
1518 map_ttr(i,j)%array(1)%itr,i6), that%gaid(i3,i,j,i6))
1519 ELSE
1520 that%gaid(i3,i,j,i6) = this%gaid(i3, map_ttr(i,j)%array(1)%it, &
1521 map_ttr(i,j)%array(1)%itr,i6)
1522 ENDIF
1523! improve the next workflow?
1524 CALL volgrid_get_vol_2d(this, i3, map_ttr(i,j)%array(ninp)%it, &
1525 map_ttr(i,j)%array(ninp)%itr, i6, voldatiin)
1526 voldatiout = voldatiin
1527 CALL volgrid_get_vol_2d(this, i3, map_ttr(i,j)%array(1)%it, &
1528 map_ttr(i,j)%array(1)%itr, i6, voldatiin)
1529
1531 voldatiout(:,:) = voldatiout(:,:) - voldatiin(:,:)
1532 ELSEWHERE
1533 voldatiout(:,:) = rmiss
1534 END WHERE
1535
1536 ELSE ! other stat_proc
1537 DO n = 1, ninp
1538 CALL volgrid_get_vol_2d(this, i3, map_ttr(i,j)%array(n)%it, &
1539 map_ttr(i,j)%array(n)%itr, i6, voldatiin)
1540
1541 IF (n == 1) THEN
1542 voldatiout = voldatiin
1543 IF (lclone) THEN
1545 map_ttr(i,j)%array(n)%itr,i6), that%gaid(i3,i,j,i6))
1546 ELSE
1547 that%gaid(i3,i,j,i6) = this%gaid(i3, map_ttr(i,j)%array(n)%it, &
1548 map_ttr(i,j)%array(n)%itr,i6)
1549 ENDIF
1550
1551 ELSE ! second or more time
1552 SELECT CASE(stat_proc)
1553 CASE (0, 1) ! average, accumulation
1555 voldatiout(:,:) = voldatiout(:,:) + voldatiin(:,:)
1556 ELSEWHERE
1557 voldatiout(:,:) = rmiss
1558 END WHERE
1559 CASE(2) ! maximum
1561 voldatiout(:,:) = max(voldatiout(:,:), voldatiin(:,:))
1562 ELSEWHERE
1563 voldatiout(:,:) = rmiss
1564 END WHERE
1565 CASE(3) ! minimum
1567 voldatiout(:,:) = min(voldatiout(:,:), voldatiin(:,:))
1568 ELSEWHERE
1569 voldatiout(:,:) = rmiss
1570 END WHERE
1571 END SELECT
1572
1573 ENDIF ! first time
1574 ENDDO
1575 IF (stat_proc == 0) THEN ! average
1577 voldatiout(:,:) = voldatiout(:,:)/ninp
1578 END WHERE
1579 ENDIF
1580 ENDIF
1581 CALL volgrid_set_vol_2d(that, i3, i, j, i6, voldatiout)
1582 ENDDO ! level
1583 ENDDO ! var
1585 ENDDO do_otime
1586ENDDO do_otimerange
1587
1588DEALLOCATE(map_ttr)
1589
1590
1591END SUBROUTINE volgrid6d_compute_stat_proc_agg
1592
1593
1594!> Specialized method for statistically processing a set of data
1595!! already processed with the same statistical processing, on a
1596!! different time interval. This method performs statistical
1597!! processing by difference of different intervals. Data with both
1598!! analysis/observation or forecast timerange are processed.
1599!!
1600!! The output \a that volgrid6d object contains elements from the original volume
1601!! \a this satisfying the conditions
1602!! - timerange (vol7d_timerange_class::vol7d_timerange::timerange)
1603!! of type \a stat_proc
1604!! - any p1 (analysis/observation or forecast)
1605!! - p2 > 0 (processing interval non null, non instantaneous data)
1606!! and equal to a multiplier of \a step if \a full_steps is \c .TRUE.
1607!!
1608!! Output data will have timerange of type \a stat_proc and p2 = \a
1609!! step. The supported statistical processing methods (parameter \a
1610!! stat_proc) are:
1611!!
1612!! - 0 average
1613!! - 1 accumulation
1614!! - 4 difference
1615!!
1616!! Input volume may have any value of \a this%time_definition, and
1617!! that value will be conserved in the output volume.
1618SUBROUTINE volgrid6d_recompute_stat_proc_diff(this, that, stat_proc, step, full_steps, start, clone)
1619TYPE(volgrid6d),INTENT(inout) :: this !< volume providing data to be recomputed, it is not modified by the method, apart from performing a \a volgrid6d_alloc_vol on it
1620TYPE(volgrid6d),INTENT(out) :: that !< output volume which will contain the recomputed data
1621INTEGER,INTENT(in) :: stat_proc !< type of statistical processing to be recomputed (from grib2 table), only data having timerange of this type will be recomputed and will appear in the output volume
1622TYPE(timedelta),INTENT(in) :: step !< length of the step over which the statistical processing is performed
1623LOGICAL,INTENT(in),OPTIONAL :: full_steps !< if provided and \a .TRUE., process only data having processing interval (p2) equal to a multiple of \a step
1624TYPE(datetime),INTENT(in),OPTIONAL :: start !< if provided, together with \a full_steps, processes data on intervals starting at \a start +- an integer amount of \a step intervals
1625LOGICAL,INTENT(in),OPTIONAL :: clone !< if provided and \c .TRUE. , clone the gaid's from \a this to \a that
1626INTEGER :: i3, i4, i6, i, j, k, l, nitr, steps
1627INTEGER,ALLOCATABLE :: map_tr(:,:,:,:,:), f(:), keep_tr(:,:,:)
1628REAL,POINTER :: voldatiin1(:,:), voldatiin2(:,:), voldatiout(:,:)
1629!LOGICAL,POINTER :: mask_timerange(:)
1630LOGICAL :: lclone
1631TYPE(vol7d_var),ALLOCATABLE :: varbufr(:)
1632
1633
1634! be safe
1635CALL volgrid6d_alloc_vol(this)
1636! when volume is not decoded it is better to clone anyway to avoid
1637! overwriting fields
1638lclone = optio_log(clone) .OR. .NOT.ASSOCIATED(this%voldati)
1639! initialise the output volume
1641CALL volgrid6d_alloc(that, dim=this%griddim%dim, &
1642 nlevel=SIZE(this%level), nvar=SIZE(this%var), ini=.false.)
1643that%level = this%level
1644that%var = this%var
1645
1646! compute length of cumulation step in seconds
1648
1649! compute the statistical processing relations, output time and
1650! timerange are defined here
1651CALL recompute_stat_proc_diff_common(this%time, this%timerange, stat_proc, step, &
1652 that%time, that%timerange, map_tr, f, keep_tr, &
1653 this%time_definition, full_steps, start)
1654nitr = SIZE(f)
1655
1656! complete the definition of the output volume
1657CALL volgrid6d_alloc_vol(that, decode=ASSOCIATED(this%voldati))
1658! allocate workspace once
1659IF (.NOT.ASSOCIATED(that%voldati)) THEN
1660 ALLOCATE(voldatiin1(this%griddim%dim%nx, this%griddim%dim%ny), &
1661 voldatiin2(this%griddim%dim%nx, this%griddim%dim%ny), &
1662 voldatiout(this%griddim%dim%nx, this%griddim%dim%ny))
1663ENDIF
1664
1665! copy the timeranges already satisfying the requested step, if any
1666DO i4 = 1, SIZE(this%time)
1667 DO i = 1, nitr
1669 l = keep_tr(i, i4, 1)
1670 k = keep_tr(i, i4, 2)
1671#ifdef DEBUG
1672 CALL l4f_category_log(this%category, l4f_debug, &
1675#endif
1676 DO i6 = 1, SIZE(this%var)
1677 DO i3 = 1, SIZE(this%level)
1679 IF (lclone) THEN
1681 ELSE
1682 that%gaid(i3,l,k,i6) = this%gaid(i3,i4,f(i),i6)
1683 ENDIF
1684 IF (ASSOCIATED(that%voldati)) THEN
1685 that%voldati(:,:,i3,l,k,i6) = this%voldati(:,:,i3,i4,f(i),i6)
1686 ELSE
1687 CALL volgrid_get_vol_2d(this, i3, i4, f(i), i6, voldatiout)
1688 CALL volgrid_set_vol_2d(that, i3, l, k, i6, voldatiout)
1689 ENDIF
1690 ENDIF
1691 ENDDO
1692 ENDDO
1693 ENDIF
1694 ENDDO
1695ENDDO
1696
1697! varbufr required for setting posdef, optimize with an array
1698ALLOCATE(varbufr(SIZE(this%var)))
1699DO i6 = 1, SIZE(this%var)
1700 varbufr(i6) = convert(this%var(i6))
1701ENDDO
1702! compute statistical processing
1703DO l = 1, SIZE(this%time)
1704 DO k = 1, nitr
1705 DO j = 1, SIZE(this%time)
1706 DO i = 1, nitr
1708 DO i6 = 1, SIZE(this%var)
1709 DO i3 = 1, SIZE(this%level)
1710
1713! take the gaid from the second time/timerange contributing to the
1714! result (l,f(k))
1715 IF (lclone) THEN
1717 that%gaid(i3,map_tr(i,j,k,l,1),map_tr(i,j,k,l,2),i6))
1718 ELSE
1719 that%gaid(i3,map_tr(i,j,k,l,1),map_tr(i,j,k,l,2),i6) = &
1720 this%gaid(i3,l,f(k),i6)
1721 ENDIF
1722
1723! get/set 2d sections API is used
1724 CALL volgrid_get_vol_2d(this, i3, l, f(k), i6, voldatiin1)
1725 CALL volgrid_get_vol_2d(this, i3, j, f(i), i6, voldatiin2)
1726 IF (ASSOCIATED(that%voldati)) &
1727 CALL volgrid_get_vol_2d(that, i3, &
1728 map_tr(i,j,k,l,1), map_tr(i,j,k,l,2), i6, voldatiout)
1729
1730 IF (stat_proc == 0) THEN ! average
1732 voldatiout(:,:) = &
1733 (voldatiin1(:,:)*this%timerange(f(k))%p2 - &
1734 voldatiin2(:,:)*this%timerange(f(i))%p2)/ &
1735 steps
1736 ELSEWHERE
1737 voldatiout(:,:) = rmiss
1738 END WHERE
1739 ELSE IF (stat_proc == 1 .OR. stat_proc == 4) THEN ! acc, diff
1741 voldatiout(:,:) = voldatiin1(:,:) - voldatiin2(:,:)
1742 ELSEWHERE
1743 voldatiout(:,:) = rmiss
1744 END WHERE
1745 IF (stat_proc == 1) THEN
1746 CALL vol7d_var_features_posdef_apply(varbufr(i6), voldatiout)
1747 ENDIF
1748 ENDIF
1749
1750 CALL volgrid_set_vol_2d(that, i3, &
1751 map_tr(i,j,k,l,1), map_tr(i,j,k,l,2), i6, voldatiout)
1752
1753 ENDIF
1754 ENDDO
1755 ENDDO
1756 ENDIF
1757 ENDDO
1758 ENDDO
1759 ENDDO
1760ENDDO
1761
1762IF (.NOT.ASSOCIATED(that%voldati)) THEN
1763 DEALLOCATE(voldatiin1, voldatiin2, voldatiout)
1764ENDIF
1765
1766END SUBROUTINE volgrid6d_recompute_stat_proc_diff
1767
1768
1769!> Specialized method for statistically processing a set of data
1770!! by integration/differentiation.
1771!! This method performs statistical processing by integrating
1772!! (accumulating) in time values representing time-average rates or
1773!! fluxes, (\a stat_proc_input=0 \a stat_proc=1) or by transforming a
1774!! time-integrated (accumulated) value in a time-average rate or flux
1775!! (\a stat_proc_input=1 \a stat_proc=0). Analysis/observation or
1776!! forecast timeranges are processed. The only operation performed is
1777!! respectively multiplying or dividing the values by the length of
1778!! the time interval in seconds.
1779!!
1780!! The output \a that volgrid6d object contains elements from the
1781!! original volume \a this satisfying the conditions
1782!! - timerange (vol7d_timerange_class::vol7d_timerange::timerange)
1783!! of type \a stat_proc_input (0 or 1)
1784!! - any p1 (analysis/observation or forecast)
1785!! - p2 > 0 (processing interval non null, non instantaneous data)
1786!!
1787!! Output data will have timerange of type \a stat_proc (1 or 0) and
1788!! p1 and p2 equal to the corresponding input values. The supported
1789!! statistical processing methods (parameter \a stat_proc) are:
1790!!
1791!! - 0 average
1792!! - 1 accumulation
1793!!
1794!! Input volume may have any value of \a this%time_definition, and
1795!! that value will be conserved in the output volume.
1796SUBROUTINE volgrid6d_compute_stat_proc_metamorph(this, that, stat_proc_input, stat_proc, clone)
1797TYPE(volgrid6d),INTENT(inout) :: this !< volume providing data to be recomputed, it is not modified by the method, apart from performing a \a volgrid6d_alloc_vol on it
1798TYPE(volgrid6d),INTENT(out) :: that !< output volume which will contain the recomputed data
1799INTEGER,INTENT(in) :: stat_proc_input !< type of statistical processing of data that has to be processed (from grib2 table), only data having timerange of this type will be processed, the actual statistical processing performed and which will appear in the output volume, is however determined by \a stat_proc argument
1800INTEGER,INTENT(in) :: stat_proc !< type of statistical processing to be recomputed (from grib2 table), data in output volume \a that will have a timerange of this type
1801LOGICAL , INTENT(in),OPTIONAL :: clone !< if provided and \c .TRUE. , clone the gaid's from \a this to \a that
1802
1803INTEGER :: j, i3, i4, i6
1804INTEGER,POINTER :: map_tr(:)
1805REAL,POINTER :: voldatiin(:,:), voldatiout(:,:)
1806REAL,ALLOCATABLE :: int_ratio(:)
1807LOGICAL :: lclone
1808
1809NULLIFY(voldatiin, voldatiout)
1810
1811! be safe
1812CALL volgrid6d_alloc_vol(this)
1813! when volume is not decoded it is better to clone anyway to avoid
1814! overwriting fields
1815lclone = optio_log(clone) .OR. .NOT.ASSOCIATED(this%voldati)
1816
1817IF (.NOT.((stat_proc_input == 0 .AND. stat_proc == 1) .OR. &
1818 (stat_proc_input == 1 .AND. stat_proc == 0))) THEN
1819
1820 CALL l4f_category_log(this%category, l4f_warn, &
1821 'compute_stat_proc_metamorph, can only be applied to average->accumulated timerange and viceversa')
1822! return an empty volume, without signaling error
1824 CALL volgrid6d_alloc_vol(that)
1825 RETURN
1826ENDIF
1827
1828! initialise the output volume
1830CALL volgrid6d_alloc(that, dim=this%griddim%dim, ntime=SIZE(this%time), &
1831 nlevel=SIZE(this%level), nvar=SIZE(this%var), ini=.false.)
1832that%time = this%time
1833that%level = this%level
1834that%var = this%var
1835
1836CALL compute_stat_proc_metamorph_common(stat_proc_input, this%timerange, stat_proc, &
1837 that%timerange, map_tr)
1838
1839! complete the definition of the output volume
1840CALL volgrid6d_alloc_vol(that, decode=ASSOCIATED(this%voldati))
1841
1842IF (stat_proc == 0) THEN ! average -> integral
1843 int_ratio = 1./real(that%timerange(:)%p2)
1844ELSE ! cumulation
1845 int_ratio = real(that%timerange(:)%p2)
1846ENDIF
1847
1848DO i6 = 1, SIZE(this%var)
1849 DO j = 1, SIZE(map_tr)
1850 DO i4 = 1, SIZE(that%time)
1851 DO i3 = 1, SIZE(this%level)
1852
1853 IF (lclone) THEN
1855 ELSE
1856 that%gaid(i3,i4,map_tr(j),i6) = this%gaid(i3,i4,j,i6)
1857 ENDIF
1858 CALL volgrid_get_vol_2d(this, i3, i4, map_tr(j), i6, voldatiin)
1859 CALL volgrid_get_vol_2d(that, i3, i4, j, i6, voldatiout)
1861 voldatiout = voldatiin*int_ratio(j)
1862 ELSEWHERE
1863 voldatiout = rmiss
1864 END WHERE
1865 CALL volgrid_set_vol_2d(that, i3, i4, j, i6, voldatiout)
1866 ENDDO
1867 ENDDO
1868 ENDDO
1869ENDDO
1870
1871
1872END SUBROUTINE volgrid6d_compute_stat_proc_metamorph
1873
1874!> Method for building a volume containing the vertical coordinate as
1875!! a variable.
1876!! This method produces a volgrid6d volume, derived from \a this,
1877!! containing a single variable, horizontally constant, on the same
1878!! input levels, which describes the vertical coordinate in the form
1879!! of a physical variable. The grid, time and timerange metadata are
1880!! the same as for the original volume. Only a single vertical level
1881!! type, the one matching the \a level argument, is converted to a
1882!! variable. The \a level argument can also indicate the layer between
1883!! two surfaces of the same type, in that case the variable
1884!! representing the vertical coordinate will be set to the value of
1885!! the midpoint between the two layers. If something goes wrong,
1886!! e.g. no level matches \a level argument or the level canot be
1887!! converted to a physical value, an empty volume is returned.
1888SUBROUTINE volgrid6d_compute_vert_coord_var(this, level, volgrid_lev)
1889TYPE(volgrid6d),INTENT(in) :: this !< volume with the vertical levels
1890TYPE(vol7d_level),INTENT(in) :: level !< vertical level to be converted to variable, only the type(s) of level are used not the value(s)
1891TYPE(volgrid6d),INTENT(out) :: volgrid_lev !< output volume with the variable describing the vertical coordinate
1892
1893INTEGER :: nlev, i, ii, iii, iiii
1894TYPE(grid_id) :: out_gaid
1895LOGICAL,ALLOCATABLE :: levmask(:)
1896TYPE(volgrid6d_var) :: lev_var
1897
1899IF (.NOT.ASSOCIATED(this%gaid)) THEN
1900 CALL l4f_log(l4f_error, 'volgrid6d_compute_vert_coord_var: input volume not allocated')
1901 RETURN
1902ENDIF
1903! if layer, both surfaces must be of the same type
1905 CALL l4f_log(l4f_error, 'volgrid6d_compute_vert_coord_var: requested (mixed) layer type not valid')
1906 RETURN
1907ENDIF
1908
1909! look for valid levels to be converted to vars
1910ALLOCATE(levmask(SIZE(this%level)))
1911levmask = this%level%level1 == level%level1 .AND. &
1912 this%level%level2 == level%level2 .AND. c_e(this%level%l1)
1914nlev = count(levmask)
1915IF (nlev == 0) THEN
1916 CALL l4f_log(l4f_error, 'volgrid6d_compute_vert_coord_var: requested level type not available')
1917 RETURN
1918ENDIF
1919
1920out_gaid = grid_id_new()
1921gaidloop: DO i=1 ,SIZE(this%gaid,1)
1922 DO ii=1 ,SIZE(this%gaid,2)
1923 DO iii=1 ,SIZE(this%gaid,3)
1924 DO iiii=1 ,SIZE(this%gaid,4)
1927 EXIT gaidloop
1928 ENDIF
1929 ENDDO
1930 ENDDO
1931 ENDDO
1932ENDDO gaidloop
1933
1934! look for variable corresponding to level
1935lev_var = convert(vol7d_var_new(btable=vol7d_level_to_var(level)), &
1936 grid_id_template=out_gaid)
1938 CALL l4f_log(l4f_error, 'volgrid6d_compute_vert_coord_var: no variable corresponds to requested level type')
1939 RETURN
1940ENDIF
1941
1942! prepare output volume
1944 time_definition=this%time_definition) !, categoryappend=categoryappend)
1945CALL volgrid6d_alloc(volgrid_lev, ntime=SIZE(this%time), nlevel=nlev, &
1946 ntimerange=SIZE(this%timerange), nvar=1)
1947! fill metadata
1948volgrid_lev%time = this%time
1949volgrid_lev%level = pack(this%level, mask=levmask)
1950volgrid_lev%timerange = this%timerange
1951volgrid_lev%var(1) = lev_var
1952
1953CALL volgrid6d_alloc_vol(volgrid_lev, decode=.true.)
1954! fill data
1955DO i = 1, nlev
1957 volgrid_lev%voldati(:,:,i,:,:,:) = real(volgrid_lev%level(i)%l1 + &
1958 volgrid_lev%level(i)%l2)* &
1959 vol7d_level_to_var_factor(volgrid_lev%level(i))/2.
1960 ELSE
1961 volgrid_lev%voldati(:,:,i,:,:,:) = real(volgrid_lev%level(i)%l1)* &
1962 vol7d_level_to_var_factor(volgrid_lev%level(i))
1963 ENDIF
1964ENDDO
1965! fill gaid for subsequent export
1967 DO i=1 ,SIZE(volgrid_lev%gaid,1)
1968 DO ii=1 ,SIZE(volgrid_lev%gaid,2)
1969 DO iii=1 ,SIZE(volgrid_lev%gaid,3)
1970 DO iiii=1 ,SIZE(volgrid_lev%gaid,4)
1972 ENDDO
1973 ENDDO
1974 ENDDO
1975 ENDDO
1977ENDIF
1978
1979END SUBROUTINE volgrid6d_compute_vert_coord_var
1980
Restituiscono il valore dell'oggetto nella forma desiderata. Definition datetime_class.F90:322 Costruttori per le classi datetime e timedelta. Definition datetime_class.F90:311 Functions that return a trimmed CHARACTER representation of the input variable. Definition datetime_class.F90:349 Restituiscono il valore dell'oggetto in forma di stringa stampabile. Definition datetime_class.F90:327 Make a deep copy, if possible, of the grid identifier. Definition grid_id_class.F90:336 Apply the conversion function this to values. Definition volgrid6d_var_class.F90:396 This module defines an abstract interface to different drivers for access to files containing gridded... Definition grid_id_class.F90:249 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. Definition stat_proc_engine.F90:211 Extension of volgrid6d_class with methods for performing simple statistical operations on entire volu... Definition volgrid6d_class_compute.F90:214 This module defines objects and methods for managing data volumes on rectangular georeferenced grids. Definition volgrid6d_class.F90:240 Class for managing physical variables in a grib 1/2 fashion. Definition volgrid6d_var_class.F90:218 |