libsim Versione 7.2.0

◆ vdb()

elemental logical function vdb ( integer(kind=int_b), intent(in)  flag)

Data validity check for confidence.

Parametri
[in]flagconfidenza

Definizione alla linea 1057 del file modqc.F90.

1058! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
1059! authors:
1060! Davide Cesari <dcesari@arpa.emr.it>
1061! Paolo Patruno <ppatruno@arpa.emr.it>
1062
1063! This program is free software; you can redistribute it and/or
1064! modify it under the terms of the GNU General Public License as
1065! published by the Free Software Foundation; either version 2 of
1066! the License, or (at your option) any later version.
1067
1068! This program is distributed in the hope that it will be useful,
1069! but WITHOUT ANY WARRANTY; without even the implied warranty of
1070! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
1071! GNU General Public License for more details.
1072
1073! You should have received a copy of the GNU General Public License
1074! along with this program. If not, see <http://www.gnu.org/licenses/>.
1075#include "config.h"
1076
1079
1226module modqc
1227use kinds
1230use vol7d_class
1231
1232
1233implicit none
1234
1235
1237type :: qcpartype
1238 integer (kind=int_b):: att
1239 integer (kind=int_b):: gross_error ! special valuer for "*B33192" when gross error check failed
1240 integer (kind=int_b):: invalidated ! special valuer for "*B33196" when manual invalidation happen
1241end type qcpartype
1242
1244type(qcpartype) :: qcpar=qcpartype(10_int_b,0_int_b,1_int_b)
1245
1246integer, parameter :: nqcattrvars=4
1247CHARACTER(len=10),parameter :: qcattrvarsbtables(nqcattrvars)=(/"*B33196","*B33192","*B33193","*B33194"/)
1248
1249type :: qcattrvars
1250 TYPE(vol7d_var) :: vars(nqcattrvars)
1251 CHARACTER(len=10) :: btables(nqcattrvars)
1252end type qcattrvars
1253
1255interface init
1256 module procedure init_qcattrvars
1257end interface
1258
1260interface peeled
1261 module procedure peeledrb, peeleddb, peeledbb, peeledib, peeledcb &
1262 ,peeledri, peeleddi, peeledbi, peeledii, peeledci &
1263 ,peeledrr, peeleddr, peeledbr, peeledir, peeledcr &
1264 ,peeledrd, peeleddd, peeledbd, peeledid, peeledcd &
1265 ,peeledrc, peeleddc, peeledbc, peeledic, peeledcc
1266end interface
1267
1268
1270interface vd
1271 module procedure vdi,vdb,vdr,vdd,vdc
1272end interface
1273
1275interface vdge
1276 module procedure vdgei,vdgeb,vdger,vdged,vdgec
1277end interface
1278
1280interface invalidated
1281 module procedure invalidatedi,invalidatedb,invalidatedr,invalidatedd,invalidatedc
1282end interface
1283
1284private
1285
1286public vd, vdge, init, qcattrvars_new, invalidated, peeled, vol7d_peeling
1287public qcattrvars, nqcattrvars, qcattrvarsbtables
1288public qcpar, qcpartype, qcsummaryflagb ! ,qcsummaryflagi
1289
1290contains
1291
1292
1293! peeled routines
1294#undef VOL7D_POLY_SUBTYPE
1295#undef VOL7D_POLY_SUBTYPES
1296#undef VOL7D_POLY_ISC
1297#define VOL7D_POLY_SUBTYPE REAL
1298#define VOL7D_POLY_SUBTYPES r
1299
1300#undef VOL7D_POLY_TYPE
1301#undef VOL7D_POLY_TYPES
1302#undef VOL7D_POLY_ISC
1303#undef VOL7D_POLY_TYPES_SUBTYPES
1304#define VOL7D_POLY_TYPE REAL
1305#define VOL7D_POLY_TYPES r
1306#define VOL7D_POLY_TYPES_SUBTYPES rr
1307#include "modqc_peeled_include.F90"
1308#include "modqc_peel_util_include.F90"
1309#undef VOL7D_POLY_TYPE
1310#undef VOL7D_POLY_TYPES
1311#undef VOL7D_POLY_TYPES_SUBTYPES
1312#define VOL7D_POLY_TYPE DOUBLE PRECISION
1313#define VOL7D_POLY_TYPES d
1314#define VOL7D_POLY_TYPES_SUBTYPES dr
1315#include "modqc_peeled_include.F90"
1316#include "modqc_peel_util_include.F90"
1317#undef VOL7D_POLY_TYPE
1318#undef VOL7D_POLY_TYPES
1319#undef VOL7D_POLY_TYPES_SUBTYPES
1320#define VOL7D_POLY_TYPE INTEGER
1321#define VOL7D_POLY_TYPES i
1322#define VOL7D_POLY_TYPES_SUBTYPES ir
1323#include "modqc_peeled_include.F90"
1324#include "modqc_peel_util_include.F90"
1325#undef VOL7D_POLY_TYPE
1326#undef VOL7D_POLY_TYPES
1327#undef VOL7D_POLY_TYPES_SUBTYPES
1328#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1329#define VOL7D_POLY_TYPES b
1330#define VOL7D_POLY_TYPES_SUBTYPES br
1331#include "modqc_peeled_include.F90"
1332#include "modqc_peel_util_include.F90"
1333#undef VOL7D_POLY_TYPE
1334#undef VOL7D_POLY_TYPES
1335#undef VOL7D_POLY_TYPES_SUBTYPES
1336#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1337#define VOL7D_POLY_TYPES c
1338#define VOL7D_POLY_ISC = 1
1339#define VOL7D_POLY_TYPES_SUBTYPES cr
1340#include "modqc_peeled_include.F90"
1341#include "modqc_peel_util_include.F90"
1342
1343
1344#undef VOL7D_POLY_SUBTYPE
1345#undef VOL7D_POLY_SUBTYPES
1346#undef VOL7D_POLY_ISC
1347#define VOL7D_POLY_SUBTYPE DOUBLE PRECISION
1348#define VOL7D_POLY_SUBTYPES d
1349
1350#undef VOL7D_POLY_TYPE
1351#undef VOL7D_POLY_TYPES
1352#undef VOL7D_POLY_TYPES_SUBTYPES
1353#define VOL7D_POLY_TYPE REAL
1354#define VOL7D_POLY_TYPES r
1355#define VOL7D_POLY_TYPES_SUBTYPES rd
1356#include "modqc_peeled_include.F90"
1357#undef VOL7D_POLY_TYPE
1358#undef VOL7D_POLY_TYPES
1359#undef VOL7D_POLY_TYPES_SUBTYPES
1360#define VOL7D_POLY_TYPE DOUBLE PRECISION
1361#define VOL7D_POLY_TYPES d
1362#define VOL7D_POLY_TYPES_SUBTYPES dd
1363#include "modqc_peeled_include.F90"
1364#undef VOL7D_POLY_TYPE
1365#undef VOL7D_POLY_TYPES
1366#undef VOL7D_POLY_TYPES_SUBTYPES
1367#define VOL7D_POLY_TYPE INTEGER
1368#define VOL7D_POLY_TYPES i
1369#define VOL7D_POLY_TYPES_SUBTYPES id
1370#include "modqc_peeled_include.F90"
1371#undef VOL7D_POLY_TYPE
1372#undef VOL7D_POLY_TYPES
1373#undef VOL7D_POLY_TYPES_SUBTYPES
1374#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1375#define VOL7D_POLY_TYPES b
1376#define VOL7D_POLY_TYPES_SUBTYPES bd
1377#include "modqc_peeled_include.F90"
1378#undef VOL7D_POLY_TYPE
1379#undef VOL7D_POLY_TYPES
1380#undef VOL7D_POLY_TYPES_SUBTYPES
1381#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1382#define VOL7D_POLY_TYPES c
1383#define VOL7D_POLY_TYPES_SUBTYPES cd
1384#include "modqc_peeled_include.F90"
1385
1386
1387#undef VOL7D_POLY_SUBTYPE
1388#undef VOL7D_POLY_SUBTYPES
1389#undef VOL7D_POLY_ISC
1390#define VOL7D_POLY_SUBTYPE INTEGER
1391#define VOL7D_POLY_SUBTYPES i
1392
1393#undef VOL7D_POLY_TYPE
1394#undef VOL7D_POLY_TYPES
1395#undef VOL7D_POLY_TYPES_SUBTYPES
1396#define VOL7D_POLY_TYPE REAL
1397#define VOL7D_POLY_TYPES r
1398#define VOL7D_POLY_TYPES_SUBTYPES ri
1399#include "modqc_peeled_include.F90"
1400#undef VOL7D_POLY_TYPE
1401#undef VOL7D_POLY_TYPES
1402#undef VOL7D_POLY_TYPES_SUBTYPES
1403#define VOL7D_POLY_TYPE DOUBLE PRECISION
1404#define VOL7D_POLY_TYPES d
1405#define VOL7D_POLY_TYPES_SUBTYPES di
1406#include "modqc_peeled_include.F90"
1407#undef VOL7D_POLY_TYPE
1408#undef VOL7D_POLY_TYPES
1409#undef VOL7D_POLY_TYPES_SUBTYPES
1410#define VOL7D_POLY_TYPE INTEGER
1411#define VOL7D_POLY_TYPES i
1412#define VOL7D_POLY_TYPES_SUBTYPES ii
1413#include "modqc_peeled_include.F90"
1414#undef VOL7D_POLY_TYPE
1415#undef VOL7D_POLY_TYPES
1416#undef VOL7D_POLY_TYPES_SUBTYPES
1417#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1418#define VOL7D_POLY_TYPES b
1419#define VOL7D_POLY_TYPES_SUBTYPES bi
1420#include "modqc_peeled_include.F90"
1421#undef VOL7D_POLY_TYPE
1422#undef VOL7D_POLY_TYPES
1423#undef VOL7D_POLY_TYPES_SUBTYPES
1424#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1425#define VOL7D_POLY_TYPES c
1426#define VOL7D_POLY_ISC = 1
1427#define VOL7D_POLY_TYPES_SUBTYPES ci
1428#include "modqc_peeled_include.F90"
1429
1430
1431#undef VOL7D_POLY_SUBTYPE
1432#undef VOL7D_POLY_SUBTYPES
1433#undef VOL7D_POLY_ISC
1434#define VOL7D_POLY_SUBTYPE INTEGER(kind=int_b)
1435#define VOL7D_POLY_SUBTYPES b
1436
1437#undef VOL7D_POLY_TYPE
1438#undef VOL7D_POLY_TYPES
1439#undef VOL7D_POLY_TYPES_SUBTYPES
1440#define VOL7D_POLY_TYPE REAL
1441#define VOL7D_POLY_TYPES r
1442#define VOL7D_POLY_TYPES_SUBTYPES rb
1443#include "modqc_peeled_include.F90"
1444#undef VOL7D_POLY_TYPE
1445#undef VOL7D_POLY_TYPES
1446#undef VOL7D_POLY_TYPES_SUBTYPES
1447#define VOL7D_POLY_TYPE DOUBLE PRECISION
1448#define VOL7D_POLY_TYPES d
1449#define VOL7D_POLY_TYPES_SUBTYPES db
1450#include "modqc_peeled_include.F90"
1451#undef VOL7D_POLY_TYPE
1452#undef VOL7D_POLY_TYPES
1453#undef VOL7D_POLY_TYPES_SUBTYPES
1454#define VOL7D_POLY_TYPE INTEGER
1455#define VOL7D_POLY_TYPES i
1456#define VOL7D_POLY_TYPES_SUBTYPES ib
1457#include "modqc_peeled_include.F90"
1458#undef VOL7D_POLY_TYPE
1459#undef VOL7D_POLY_TYPES
1460#undef VOL7D_POLY_TYPES_SUBTYPES
1461#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1462#define VOL7D_POLY_TYPES b
1463#define VOL7D_POLY_TYPES_SUBTYPES bb
1464#include "modqc_peeled_include.F90"
1465#undef VOL7D_POLY_TYPE
1466#undef VOL7D_POLY_TYPES
1467#undef VOL7D_POLY_TYPES_SUBTYPES
1468#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1469#define VOL7D_POLY_TYPES c
1470#define VOL7D_POLY_ISC = 1
1471#define VOL7D_POLY_TYPES_SUBTYPES cb
1472#include "modqc_peeled_include.F90"
1473
1474
1475#undef VOL7D_POLY_SUBTYPE
1476#undef VOL7D_POLY_SUBTYPES
1477#undef VOL7D_POLY_ISC
1478#define VOL7D_POLY_SUBTYPE CHARACTER(len=vol7d_cdatalen)
1479#define VOL7D_POLY_SUBTYPES c
1480
1481#undef VOL7D_POLY_TYPE
1482#undef VOL7D_POLY_TYPES
1483#undef VOL7D_POLY_TYPES_SUBTYPES
1484#define VOL7D_POLY_TYPE REAL
1485#define VOL7D_POLY_TYPES r
1486#define VOL7D_POLY_TYPES_SUBTYPES rc
1487#include "modqc_peeled_include.F90"
1488#undef VOL7D_POLY_TYPE
1489#undef VOL7D_POLY_TYPES
1490#undef VOL7D_POLY_TYPES_SUBTYPES
1491#define VOL7D_POLY_TYPE DOUBLE PRECISION
1492#define VOL7D_POLY_TYPES d
1493#define VOL7D_POLY_TYPES_SUBTYPES dc
1494#include "modqc_peeled_include.F90"
1495#undef VOL7D_POLY_TYPE
1496#undef VOL7D_POLY_TYPES
1497#undef VOL7D_POLY_TYPES_SUBTYPES
1498#define VOL7D_POLY_TYPE INTEGER
1499#define VOL7D_POLY_TYPES i
1500#define VOL7D_POLY_TYPES_SUBTYPES ic
1501#include "modqc_peeled_include.F90"
1502#undef VOL7D_POLY_TYPE
1503#undef VOL7D_POLY_TYPES
1504#undef VOL7D_POLY_TYPES_SUBTYPES
1505#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1506#define VOL7D_POLY_TYPES b
1507#define VOL7D_POLY_TYPES_SUBTYPES bc
1508#include "modqc_peeled_include.F90"
1509#undef VOL7D_POLY_TYPE
1510#undef VOL7D_POLY_TYPES
1511#undef VOL7D_POLY_TYPES_SUBTYPES
1512#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1513#define VOL7D_POLY_TYPES c
1514#define VOL7D_POLY_ISC = 1
1515#define VOL7D_POLY_TYPES_SUBTYPES cc
1516#include "modqc_peeled_include.F90"
1517
1518
1519subroutine init_qcattrvars(this)
1520
1521type(qcattrvars),intent(inout) :: this
1522integer :: i
1523
1524this%btables(:) =qcattrvarsbtables
1525do i =1, nqcattrvars
1526 call init(this%vars(i),this%btables(i))
1527end do
1528
1529end subroutine init_qcattrvars
1530
1531
1532type(qcattrvars) function qcattrvars_new()
1533
1534call init(qcattrvars_new)
1535
1536end function qcattrvars_new
1537
1538
1546SUBROUTINE vol7d_peeling(this, data_id, keep_attr, delete_attr, preserve, purgeana)
1547TYPE(vol7d),INTENT(INOUT) :: this
1548integer,INTENT(inout),pointer,OPTIONAL :: data_id(:,:,:,:,:)
1549CHARACTER(len=*),INTENT(in),OPTIONAL :: keep_attr(:)
1550CHARACTER(len=*),INTENT(in),OPTIONAL :: delete_attr(:)
1551logical,intent(in),optional :: preserve
1552logical,intent(in),optional :: purgeana
1553
1554integer :: inddativar,inddatiattrinv,inddatiattrcli,inddatiattrtem,inddatiattrspa,inddativarattr
1555type(qcattrvars) :: attrvars
1556
1557INTEGER(kind=int_b),pointer :: invbb(:,:,:,:,:),clibb(:,:,:,:,:),tembb(:,:,:,:,:),spabb(:,:,:,:,:)
1558INTEGER,pointer :: invbi(:,:,:,:,:),clibi(:,:,:,:,:),tembi(:,:,:,:,:),spabi(:,:,:,:,:)
1559REAL,pointer :: invbr(:,:,:,:,:),clibr(:,:,:,:,:),tembr(:,:,:,:,:),spabr(:,:,:,:,:)
1560DOUBLE PRECISION,pointer :: invbd(:,:,:,:,:),clibd(:,:,:,:,:),tembd(:,:,:,:,:),spabd(:,:,:,:,:)
1561CHARACTER(len=vol7d_cdatalen),pointer :: invbc(:,:,:,:,:),clibc(:,:,:,:,:),tembc(:,:,:,:,:),spabc(:,:,:,:,:)
1562
1563call l4f_log(l4f_info,'starting peeling')
1564
1565call init(attrvars)
1566
1567! generate code per i vari tipi di dati di v7d
1568! tramite un template e il preprocessore
1569
1570
1571#undef VOL7D_POLY_SUBTYPE
1572#undef VOL7D_POLY_SUBTYPES
1573#define VOL7D_POLY_SUBTYPE REAL
1574#define VOL7D_POLY_SUBTYPES r
1575
1576#undef VOL7D_POLY_TYPE
1577#undef VOL7D_POLY_TYPES
1578#define VOL7D_POLY_TYPE REAL
1579#define VOL7D_POLY_TYPES r
1580#include "modqc_peeling_include.F90"
1581#undef VOL7D_POLY_TYPE
1582#undef VOL7D_POLY_TYPES
1583#define VOL7D_POLY_TYPE DOUBLE PRECISION
1584#define VOL7D_POLY_TYPES d
1585#include "modqc_peeling_include.F90"
1586#undef VOL7D_POLY_TYPE
1587#undef VOL7D_POLY_TYPES
1588#define VOL7D_POLY_TYPE INTEGER
1589#define VOL7D_POLY_TYPES i
1590#include "modqc_peeling_include.F90"
1591#undef VOL7D_POLY_TYPE
1592#undef VOL7D_POLY_TYPES
1593#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1594#define VOL7D_POLY_TYPES b
1595#include "modqc_peeling_include.F90"
1596#undef VOL7D_POLY_TYPE
1597#undef VOL7D_POLY_TYPES
1598#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1599#define VOL7D_POLY_TYPES c
1600#include "modqc_peeling_include.F90"
1601
1602
1603#undef VOL7D_POLY_SUBTYPE
1604#undef VOL7D_POLY_SUBTYPES
1605#define VOL7D_POLY_SUBTYPE DOUBLE PRECISION
1606#define VOL7D_POLY_SUBTYPES d
1607
1608#undef VOL7D_POLY_TYPE
1609#undef VOL7D_POLY_TYPES
1610#define VOL7D_POLY_TYPE REAL
1611#define VOL7D_POLY_TYPES r
1612#include "modqc_peeling_include.F90"
1613#undef VOL7D_POLY_TYPE
1614#undef VOL7D_POLY_TYPES
1615#define VOL7D_POLY_TYPE DOUBLE PRECISION
1616#define VOL7D_POLY_TYPES d
1617#include "modqc_peeling_include.F90"
1618#undef VOL7D_POLY_TYPE
1619#undef VOL7D_POLY_TYPES
1620#define VOL7D_POLY_TYPE INTEGER
1621#define VOL7D_POLY_TYPES i
1622#include "modqc_peeling_include.F90"
1623#undef VOL7D_POLY_TYPE
1624#undef VOL7D_POLY_TYPES
1625#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1626#define VOL7D_POLY_TYPES b
1627#include "modqc_peeling_include.F90"
1628#undef VOL7D_POLY_TYPE
1629#undef VOL7D_POLY_TYPES
1630#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1631#define VOL7D_POLY_TYPES c
1632#include "modqc_peeling_include.F90"
1633
1634
1635#undef VOL7D_POLY_SUBTYPE
1636#undef VOL7D_POLY_SUBTYPES
1637#define VOL7D_POLY_SUBTYPE INTEGER
1638#define VOL7D_POLY_SUBTYPES i
1639
1640#undef VOL7D_POLY_TYPE
1641#undef VOL7D_POLY_TYPES
1642#define VOL7D_POLY_TYPE REAL
1643#define VOL7D_POLY_TYPES r
1644#include "modqc_peeling_include.F90"
1645#undef VOL7D_POLY_TYPE
1646#undef VOL7D_POLY_TYPES
1647#define VOL7D_POLY_TYPE DOUBLE PRECISION
1648#define VOL7D_POLY_TYPES d
1649#include "modqc_peeling_include.F90"
1650#undef VOL7D_POLY_TYPE
1651#undef VOL7D_POLY_TYPES
1652#define VOL7D_POLY_TYPE INTEGER
1653#define VOL7D_POLY_TYPES i
1654#include "modqc_peeling_include.F90"
1655#undef VOL7D_POLY_TYPE
1656#undef VOL7D_POLY_TYPES
1657#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1658#define VOL7D_POLY_TYPES b
1659#include "modqc_peeling_include.F90"
1660#undef VOL7D_POLY_TYPE
1661#undef VOL7D_POLY_TYPES
1662#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1663#define VOL7D_POLY_TYPES c
1664#include "modqc_peeling_include.F90"
1665
1666
1667#undef VOL7D_POLY_SUBTYPE
1668#undef VOL7D_POLY_SUBTYPES
1669#define VOL7D_POLY_SUBTYPE INTEGER(kind=int_b)
1670#define VOL7D_POLY_SUBTYPES b
1671
1672#undef VOL7D_POLY_TYPE
1673#undef VOL7D_POLY_TYPES
1674#define VOL7D_POLY_TYPE REAL
1675#define VOL7D_POLY_TYPES r
1676#include "modqc_peeling_include.F90"
1677#undef VOL7D_POLY_TYPE
1678#undef VOL7D_POLY_TYPES
1679#define VOL7D_POLY_TYPE DOUBLE PRECISION
1680#define VOL7D_POLY_TYPES d
1681#include "modqc_peeling_include.F90"
1682#undef VOL7D_POLY_TYPE
1683#undef VOL7D_POLY_TYPES
1684#define VOL7D_POLY_TYPE INTEGER
1685#define VOL7D_POLY_TYPES i
1686#include "modqc_peeling_include.F90"
1687#undef VOL7D_POLY_TYPE
1688#undef VOL7D_POLY_TYPES
1689#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1690#define VOL7D_POLY_TYPES b
1691#include "modqc_peeling_include.F90"
1692#undef VOL7D_POLY_TYPE
1693#undef VOL7D_POLY_TYPES
1694#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1695#define VOL7D_POLY_TYPES c
1696#include "modqc_peeling_include.F90"
1697
1698
1699
1700#undef VOL7D_POLY_SUBTYPE
1701#undef VOL7D_POLY_SUBTYPES
1702#define VOL7D_POLY_SUBTYPE CHARACTER(len=vol7d_cdatalen)
1703#define VOL7D_POLY_SUBTYPES c
1704
1705#undef VOL7D_POLY_TYPE
1706#undef VOL7D_POLY_TYPES
1707#define VOL7D_POLY_TYPE REAL
1708#define VOL7D_POLY_TYPES r
1709#include "modqc_peeling_include.F90"
1710#undef VOL7D_POLY_TYPE
1711#undef VOL7D_POLY_TYPES
1712#define VOL7D_POLY_TYPE DOUBLE PRECISION
1713#define VOL7D_POLY_TYPES d
1714#include "modqc_peeling_include.F90"
1715#undef VOL7D_POLY_TYPE
1716#undef VOL7D_POLY_TYPES
1717#define VOL7D_POLY_TYPE INTEGER
1718#define VOL7D_POLY_TYPES i
1719#include "modqc_peeling_include.F90"
1720#undef VOL7D_POLY_TYPE
1721#undef VOL7D_POLY_TYPES
1722#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1723#define VOL7D_POLY_TYPES b
1724#include "modqc_peeling_include.F90"
1725#undef VOL7D_POLY_TYPE
1726#undef VOL7D_POLY_TYPES
1727#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1728#define VOL7D_POLY_TYPES c
1729#include "modqc_peeling_include.F90"
1730
1731
1732
1733IF (.NOT.PRESENT(keep_attr) .AND. .NOT.PRESENT(delete_attr) .and. .not. optio_log(preserve)) THEN ! destroy all attributes
1734 IF (ASSOCIATED(this%voldatiattrr)) DEALLOCATE(this%voldatiattrr)
1735 IF (ASSOCIATED(this%voldatiattrd)) DEALLOCATE(this%voldatiattrd)
1736 IF (ASSOCIATED(this%voldatiattri)) DEALLOCATE(this%voldatiattri)
1737 IF (ASSOCIATED(this%voldatiattrb)) DEALLOCATE(this%voldatiattrb)
1738 IF (ASSOCIATED(this%voldatiattrc)) DEALLOCATE(this%voldatiattrc)
1739
1740 CALL delete(this%datiattr)
1741 CALL delete(this%dativarattr)
1742END IF
1743
1744IF (PRESENT(keep_attr)) THEN ! set to missing non requested attributes and reform
1745
1746 if (optio_log(preserve)) call l4f_log(l4f_warn,"preserve parameter ignored: keep_attr passed")
1747 CALL keep_var(this%datiattr%r)
1748 CALL keep_var(this%datiattr%d)
1749 CALL keep_var(this%datiattr%i)
1750 CALL keep_var(this%datiattr%b)
1751 CALL keep_var(this%datiattr%c)
1752 CALL qc_reform(this,data_id, miss=.true., purgeana=purgeana)
1753
1754ELSE IF (PRESENT(delete_attr)) THEN ! set to missing requested attributes and reform
1755
1756 if (optio_log(preserve)) call l4f_log(l4f_warn,"preserve parameter ignored: delete_attr passed")
1757 CALL delete_var(this%datiattr%r)
1758 CALL delete_var(this%datiattr%d)
1759 CALL delete_var(this%datiattr%i)
1760 CALL delete_var(this%datiattr%b)
1761 CALL delete_var(this%datiattr%c)
1762 CALL qc_reform(this,data_id, miss=.true., purgeana=purgeana)
1763
1764ELSE IF (PRESENT(purgeana)) THEN
1765
1766 CALL qc_reform(this,data_id, purgeana=purgeana)
1767
1768ENDIF
1769
1770
1771CONTAINS
1772
1773
1775subroutine qc_reform(this,data_id,miss, purgeana)
1776TYPE(vol7d),INTENT(INOUT) :: this
1777integer,INTENT(inout),pointer,OPTIONAL :: data_id(:,:,:,:,:)
1778logical,intent(in),optional :: miss
1779logical,intent(in),optional :: purgeana
1780
1781integer,pointer :: data_idtmp(:,:,:,:,:)
1782logical,allocatable :: llana(:)
1783integer,allocatable :: anaind(:)
1784integer :: i,j,nana
1785
1786if (optio_log(purgeana)) then
1787 allocate(llana(size(this%ana)))
1788 llana =.false.
1789 do i =1,size(this%ana)
1790 if (associated(this%voldatii)) llana(i)= llana(i) .or. any(c_e(this%voldatii(i,:,:,:,:,:)))
1791 if (associated(this%voldatir)) llana(i)= llana(i) .or. any(c_e(this%voldatir(i,:,:,:,:,:)))
1792 if (associated(this%voldatid)) llana(i)= llana(i) .or. any(c_e(this%voldatid(i,:,:,:,:,:)))
1793 if (associated(this%voldatib)) llana(i)= llana(i) .or. any(c_e(this%voldatib(i,:,:,:,:,:)))
1794 if (associated(this%voldatic)) llana(i)= llana(i) .or. any(c_e(this%voldatic(i,:,:,:,:,:)))
1795
1796#ifdef DEBUG
1797 if (.not. llana(i)) call l4f_log(l4f_debug,"remove station"//t2c(i))
1798#endif
1799
1800 end do
1801
1802 nana=count(llana)
1803
1804
1805 allocate(anaind(nana))
1806
1807 j=0
1808 do i=1,size(this%ana)
1809 if (llana(i)) then
1810 j=j+1
1811 anaind(j)=i
1812 end if
1813 end do
1814
1815
1816 if(present(data_id)) then
1817 allocate(data_idtmp(nana,size(data_id,2),size(data_id,3),size(data_id,4),size(data_id,5)))
1818 data_idtmp=data_id(anaind,:,:,:,:)
1819 if (associated(data_id))deallocate(data_id)
1820 data_id=>data_idtmp
1821 end if
1822
1823 call vol7d_reform(this,miss=miss,lana=llana)
1824
1825 deallocate(llana,anaind)
1826
1827else
1828
1829 call vol7d_reform(this,miss=miss)
1830
1831end if
1832
1833end subroutine qc_reform
1834
1835
1836SUBROUTINE keep_var(var)
1837TYPE(vol7d_var),intent(inout),POINTER :: var(:)
1838
1839INTEGER :: i
1840
1841IF (ASSOCIATED(var)) THEN
1842 if (size(var) == 0) then
1843 var%btable = vol7d_var_miss%btable
1844 else
1845 DO i = 1, SIZE(var)
1846 IF (all(var(i)%btable /= keep_attr(:))) THEN ! n.b. ALL((//)) = .TRUE.
1847 var(i)%btable = vol7d_var_miss%btable
1848 ENDIF
1849 ENDDO
1850 end if
1851ENDIF
1852
1853END SUBROUTINE keep_var
1854
1855SUBROUTINE delete_var(var)
1856TYPE(vol7d_var),intent(inout),POINTER :: var(:)
1857
1858INTEGER :: i
1859
1860IF (ASSOCIATED(var)) THEN
1861 if (size(var) == 0) then
1862 var%btable = vol7d_var_miss%btable
1863 else
1864 DO i = 1, SIZE(var)
1865 IF (any(var(i)%btable == delete_attr(:))) THEN ! n.b. ANY((//)) = .FALSE.
1866 var(i) = vol7d_var_miss
1867 ENDIF
1868 ENDDO
1869 end if
1870ENDIF
1871
1872END SUBROUTINE delete_var
1873
1874END SUBROUTINE vol7d_peeling
1875
1876
1877end module modqc
Variables user in Quality Control.
Definition: modqc.F90:386
Test di dato invalidato.
Definition: modqc.F90:411
Remove data under a defined grade of confidence.
Definition: modqc.F90:391
Check data validity based on single confidence.
Definition: modqc.F90:401
Check data validity based on gross error check.
Definition: modqc.F90:406
Definition of constants to be used for declaring variables of a desired type.
Definition: kinds.F90:245
Definitions of constants and functions for working with missing values.
Utilities and defines for quality control.
Definition: modqc.F90:357
Module for quickly interpreting the OPTIONAL parameters passed to a subprogram.
Classe per la gestione di un volume completo di dati osservati.
Definisce il livello di attendibilità per i dati validi.
Definition: modqc.F90:368

Generated with Doxygen.