libsim Versione 7.1.11

◆ qcsummaryflagb()

elemental logical function, public qcsummaryflagb ( integer(kind=int_b), intent(in), optional  flag0,
integer(kind=int_b), intent(in), optional  flag1,
integer(kind=int_b), intent(in), optional  flag2,
integer(kind=int_b), intent(in), optional  flag3 
)

Check data validity based on multiple confidences.

Compute final decision boolean flag Quality control is complete if one of 3 conditions is verified: a) invalidated data b) gross error check failed c) tot variable less -1 Controllo di validita' del dato basato su test multipli. Per il calcolo della validita' del dato (flag booleano B33007), si prendono in considerazione 3 test; il dato risulta invalidato (flag booleano posto a false) se e solo se uno dei test risulta soddisfatto: a) il dato e' stato invalidato a mano (flag0=B33196=0) b) il dato non ha passato il gross erro check (flag1=B33192=0) c) la variabile tot risulta minore a -1 La variabile tot e' il risultato del confronto tra controllo climatologico (flag1, B33192), controllo temporale (flag2, B33193) e controllo spaziale (flag3, B33194). Ad ognuno di tali controlli e' stato attribuito un punteggio a seconda che ciascuno dei valori relativi ai flag di qualita' risulti inferiore od uguale-maggiore di 10. Nel dettaglio: se B33192 < 10 tot=-1; se B33192>=10 tot=0 se B33193 < 10 tot=-1; se B33193>=10 tot=1 se B33194 < 10 tot=-1; se B33194>=10 tot=1 Ogni dato e' controllato nei 3 flag di qualita' presenti, e viene valutata la somma risultante di tot. Se tot risulta inferiore a -1, qcsummaryflag e' posto a false ed il dato e' invalitato (B33007=0). Se tot risulta maggiore od uguale a -1 qcsummaryflag e' true ed il dato e' valido.

Definizione alla linea 1122 del file modqc.F90.

1123! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
1124! authors:
1125! Davide Cesari <dcesari@arpa.emr.it>
1126! Paolo Patruno <ppatruno@arpa.emr.it>
1127
1128! This program is free software; you can redistribute it and/or
1129! modify it under the terms of the GNU General Public License as
1130! published by the Free Software Foundation; either version 2 of
1131! the License, or (at your option) any later version.
1132
1133! This program is distributed in the hope that it will be useful,
1134! but WITHOUT ANY WARRANTY; without even the implied warranty of
1135! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
1136! GNU General Public License for more details.
1137
1138! You should have received a copy of the GNU General Public License
1139! along with this program. If not, see <http://www.gnu.org/licenses/>.
1140#include "config.h"
1141
1144
1291module modqc
1292use kinds
1295use vol7d_class
1296
1297
1298implicit none
1299
1300
1302type :: qcpartype
1303 integer (kind=int_b):: att
1304 integer (kind=int_b):: gross_error ! special valuer for "*B33192" when gross error check failed
1305 integer (kind=int_b):: invalidated ! special valuer for "*B33196" when manual invalidation happen
1306end type qcpartype
1307
1309type(qcpartype) :: qcpar=qcpartype(10_int_b,0_int_b,1_int_b)
1310
1311integer, parameter :: nqcattrvars=4
1312CHARACTER(len=10),parameter :: qcattrvarsbtables(nqcattrvars)=(/"*B33196","*B33192","*B33193","*B33194"/)
1313
1314type :: qcattrvars
1315 TYPE(vol7d_var) :: vars(nqcattrvars)
1316 CHARACTER(len=10) :: btables(nqcattrvars)
1317end type qcattrvars
1318
1320interface init
1321 module procedure init_qcattrvars
1322end interface
1323
1325interface peeled
1326 module procedure peeledrb, peeleddb, peeledbb, peeledib, peeledcb &
1327 ,peeledri, peeleddi, peeledbi, peeledii, peeledci &
1328 ,peeledrr, peeleddr, peeledbr, peeledir, peeledcr &
1329 ,peeledrd, peeleddd, peeledbd, peeledid, peeledcd &
1330 ,peeledrc, peeleddc, peeledbc, peeledic, peeledcc
1331end interface
1332
1333
1335interface vd
1336 module procedure vdi,vdb,vdr,vdd,vdc
1337end interface
1338
1340interface vdge
1341 module procedure vdgei,vdgeb,vdger,vdged,vdgec
1342end interface
1343
1345interface invalidated
1346 module procedure invalidatedi,invalidatedb,invalidatedr,invalidatedd,invalidatedc
1347end interface
1348
1349private
1350
1351public vd, vdge, init, qcattrvars_new, invalidated, peeled, vol7d_peeling
1352public qcattrvars, nqcattrvars, qcattrvarsbtables
1353public qcpar, qcpartype, qcsummaryflagb ! ,qcsummaryflagi
1354
1355contains
1356
1357
1358! peeled routines
1359#undef VOL7D_POLY_SUBTYPE
1360#undef VOL7D_POLY_SUBTYPES
1361#undef VOL7D_POLY_ISC
1362#define VOL7D_POLY_SUBTYPE REAL
1363#define VOL7D_POLY_SUBTYPES r
1364
1365#undef VOL7D_POLY_TYPE
1366#undef VOL7D_POLY_TYPES
1367#undef VOL7D_POLY_ISC
1368#undef VOL7D_POLY_TYPES_SUBTYPES
1369#define VOL7D_POLY_TYPE REAL
1370#define VOL7D_POLY_TYPES r
1371#define VOL7D_POLY_TYPES_SUBTYPES rr
1372#include "modqc_peeled_include.F90"
1373#include "modqc_peel_util_include.F90"
1374#undef VOL7D_POLY_TYPE
1375#undef VOL7D_POLY_TYPES
1376#undef VOL7D_POLY_TYPES_SUBTYPES
1377#define VOL7D_POLY_TYPE DOUBLE PRECISION
1378#define VOL7D_POLY_TYPES d
1379#define VOL7D_POLY_TYPES_SUBTYPES dr
1380#include "modqc_peeled_include.F90"
1381#include "modqc_peel_util_include.F90"
1382#undef VOL7D_POLY_TYPE
1383#undef VOL7D_POLY_TYPES
1384#undef VOL7D_POLY_TYPES_SUBTYPES
1385#define VOL7D_POLY_TYPE INTEGER
1386#define VOL7D_POLY_TYPES i
1387#define VOL7D_POLY_TYPES_SUBTYPES ir
1388#include "modqc_peeled_include.F90"
1389#include "modqc_peel_util_include.F90"
1390#undef VOL7D_POLY_TYPE
1391#undef VOL7D_POLY_TYPES
1392#undef VOL7D_POLY_TYPES_SUBTYPES
1393#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1394#define VOL7D_POLY_TYPES b
1395#define VOL7D_POLY_TYPES_SUBTYPES br
1396#include "modqc_peeled_include.F90"
1397#include "modqc_peel_util_include.F90"
1398#undef VOL7D_POLY_TYPE
1399#undef VOL7D_POLY_TYPES
1400#undef VOL7D_POLY_TYPES_SUBTYPES
1401#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1402#define VOL7D_POLY_TYPES c
1403#define VOL7D_POLY_ISC = 1
1404#define VOL7D_POLY_TYPES_SUBTYPES cr
1405#include "modqc_peeled_include.F90"
1406#include "modqc_peel_util_include.F90"
1407
1408
1409#undef VOL7D_POLY_SUBTYPE
1410#undef VOL7D_POLY_SUBTYPES
1411#undef VOL7D_POLY_ISC
1412#define VOL7D_POLY_SUBTYPE DOUBLE PRECISION
1413#define VOL7D_POLY_SUBTYPES d
1414
1415#undef VOL7D_POLY_TYPE
1416#undef VOL7D_POLY_TYPES
1417#undef VOL7D_POLY_TYPES_SUBTYPES
1418#define VOL7D_POLY_TYPE REAL
1419#define VOL7D_POLY_TYPES r
1420#define VOL7D_POLY_TYPES_SUBTYPES rd
1421#include "modqc_peeled_include.F90"
1422#undef VOL7D_POLY_TYPE
1423#undef VOL7D_POLY_TYPES
1424#undef VOL7D_POLY_TYPES_SUBTYPES
1425#define VOL7D_POLY_TYPE DOUBLE PRECISION
1426#define VOL7D_POLY_TYPES d
1427#define VOL7D_POLY_TYPES_SUBTYPES dd
1428#include "modqc_peeled_include.F90"
1429#undef VOL7D_POLY_TYPE
1430#undef VOL7D_POLY_TYPES
1431#undef VOL7D_POLY_TYPES_SUBTYPES
1432#define VOL7D_POLY_TYPE INTEGER
1433#define VOL7D_POLY_TYPES i
1434#define VOL7D_POLY_TYPES_SUBTYPES id
1435#include "modqc_peeled_include.F90"
1436#undef VOL7D_POLY_TYPE
1437#undef VOL7D_POLY_TYPES
1438#undef VOL7D_POLY_TYPES_SUBTYPES
1439#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1440#define VOL7D_POLY_TYPES b
1441#define VOL7D_POLY_TYPES_SUBTYPES bd
1442#include "modqc_peeled_include.F90"
1443#undef VOL7D_POLY_TYPE
1444#undef VOL7D_POLY_TYPES
1445#undef VOL7D_POLY_TYPES_SUBTYPES
1446#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1447#define VOL7D_POLY_TYPES c
1448#define VOL7D_POLY_TYPES_SUBTYPES cd
1449#include "modqc_peeled_include.F90"
1450
1451
1452#undef VOL7D_POLY_SUBTYPE
1453#undef VOL7D_POLY_SUBTYPES
1454#undef VOL7D_POLY_ISC
1455#define VOL7D_POLY_SUBTYPE INTEGER
1456#define VOL7D_POLY_SUBTYPES i
1457
1458#undef VOL7D_POLY_TYPE
1459#undef VOL7D_POLY_TYPES
1460#undef VOL7D_POLY_TYPES_SUBTYPES
1461#define VOL7D_POLY_TYPE REAL
1462#define VOL7D_POLY_TYPES r
1463#define VOL7D_POLY_TYPES_SUBTYPES ri
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 DOUBLE PRECISION
1469#define VOL7D_POLY_TYPES d
1470#define VOL7D_POLY_TYPES_SUBTYPES di
1471#include "modqc_peeled_include.F90"
1472#undef VOL7D_POLY_TYPE
1473#undef VOL7D_POLY_TYPES
1474#undef VOL7D_POLY_TYPES_SUBTYPES
1475#define VOL7D_POLY_TYPE INTEGER
1476#define VOL7D_POLY_TYPES i
1477#define VOL7D_POLY_TYPES_SUBTYPES ii
1478#include "modqc_peeled_include.F90"
1479#undef VOL7D_POLY_TYPE
1480#undef VOL7D_POLY_TYPES
1481#undef VOL7D_POLY_TYPES_SUBTYPES
1482#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1483#define VOL7D_POLY_TYPES b
1484#define VOL7D_POLY_TYPES_SUBTYPES bi
1485#include "modqc_peeled_include.F90"
1486#undef VOL7D_POLY_TYPE
1487#undef VOL7D_POLY_TYPES
1488#undef VOL7D_POLY_TYPES_SUBTYPES
1489#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1490#define VOL7D_POLY_TYPES c
1491#define VOL7D_POLY_ISC = 1
1492#define VOL7D_POLY_TYPES_SUBTYPES ci
1493#include "modqc_peeled_include.F90"
1494
1495
1496#undef VOL7D_POLY_SUBTYPE
1497#undef VOL7D_POLY_SUBTYPES
1498#undef VOL7D_POLY_ISC
1499#define VOL7D_POLY_SUBTYPE INTEGER(kind=int_b)
1500#define VOL7D_POLY_SUBTYPES b
1501
1502#undef VOL7D_POLY_TYPE
1503#undef VOL7D_POLY_TYPES
1504#undef VOL7D_POLY_TYPES_SUBTYPES
1505#define VOL7D_POLY_TYPE REAL
1506#define VOL7D_POLY_TYPES r
1507#define VOL7D_POLY_TYPES_SUBTYPES rb
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 DOUBLE PRECISION
1513#define VOL7D_POLY_TYPES d
1514#define VOL7D_POLY_TYPES_SUBTYPES db
1515#include "modqc_peeled_include.F90"
1516#undef VOL7D_POLY_TYPE
1517#undef VOL7D_POLY_TYPES
1518#undef VOL7D_POLY_TYPES_SUBTYPES
1519#define VOL7D_POLY_TYPE INTEGER
1520#define VOL7D_POLY_TYPES i
1521#define VOL7D_POLY_TYPES_SUBTYPES ib
1522#include "modqc_peeled_include.F90"
1523#undef VOL7D_POLY_TYPE
1524#undef VOL7D_POLY_TYPES
1525#undef VOL7D_POLY_TYPES_SUBTYPES
1526#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1527#define VOL7D_POLY_TYPES b
1528#define VOL7D_POLY_TYPES_SUBTYPES bb
1529#include "modqc_peeled_include.F90"
1530#undef VOL7D_POLY_TYPE
1531#undef VOL7D_POLY_TYPES
1532#undef VOL7D_POLY_TYPES_SUBTYPES
1533#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1534#define VOL7D_POLY_TYPES c
1535#define VOL7D_POLY_ISC = 1
1536#define VOL7D_POLY_TYPES_SUBTYPES cb
1537#include "modqc_peeled_include.F90"
1538
1539
1540#undef VOL7D_POLY_SUBTYPE
1541#undef VOL7D_POLY_SUBTYPES
1542#undef VOL7D_POLY_ISC
1543#define VOL7D_POLY_SUBTYPE CHARACTER(len=vol7d_cdatalen)
1544#define VOL7D_POLY_SUBTYPES c
1545
1546#undef VOL7D_POLY_TYPE
1547#undef VOL7D_POLY_TYPES
1548#undef VOL7D_POLY_TYPES_SUBTYPES
1549#define VOL7D_POLY_TYPE REAL
1550#define VOL7D_POLY_TYPES r
1551#define VOL7D_POLY_TYPES_SUBTYPES rc
1552#include "modqc_peeled_include.F90"
1553#undef VOL7D_POLY_TYPE
1554#undef VOL7D_POLY_TYPES
1555#undef VOL7D_POLY_TYPES_SUBTYPES
1556#define VOL7D_POLY_TYPE DOUBLE PRECISION
1557#define VOL7D_POLY_TYPES d
1558#define VOL7D_POLY_TYPES_SUBTYPES dc
1559#include "modqc_peeled_include.F90"
1560#undef VOL7D_POLY_TYPE
1561#undef VOL7D_POLY_TYPES
1562#undef VOL7D_POLY_TYPES_SUBTYPES
1563#define VOL7D_POLY_TYPE INTEGER
1564#define VOL7D_POLY_TYPES i
1565#define VOL7D_POLY_TYPES_SUBTYPES ic
1566#include "modqc_peeled_include.F90"
1567#undef VOL7D_POLY_TYPE
1568#undef VOL7D_POLY_TYPES
1569#undef VOL7D_POLY_TYPES_SUBTYPES
1570#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1571#define VOL7D_POLY_TYPES b
1572#define VOL7D_POLY_TYPES_SUBTYPES bc
1573#include "modqc_peeled_include.F90"
1574#undef VOL7D_POLY_TYPE
1575#undef VOL7D_POLY_TYPES
1576#undef VOL7D_POLY_TYPES_SUBTYPES
1577#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1578#define VOL7D_POLY_TYPES c
1579#define VOL7D_POLY_ISC = 1
1580#define VOL7D_POLY_TYPES_SUBTYPES cc
1581#include "modqc_peeled_include.F90"
1582
1583
1584subroutine init_qcattrvars(this)
1585
1586type(qcattrvars),intent(inout) :: this
1587integer :: i
1588
1589this%btables(:) =qcattrvarsbtables
1590do i =1, nqcattrvars
1591 call init(this%vars(i),this%btables(i))
1592end do
1593
1594end subroutine init_qcattrvars
1595
1596
1597type(qcattrvars) function qcattrvars_new()
1598
1599call init(qcattrvars_new)
1600
1601end function qcattrvars_new
1602
1603
1611SUBROUTINE vol7d_peeling(this, data_id, keep_attr, delete_attr, preserve, purgeana)
1612TYPE(vol7d),INTENT(INOUT) :: this
1613integer,INTENT(inout),pointer,OPTIONAL :: data_id(:,:,:,:,:)
1614CHARACTER(len=*),INTENT(in),OPTIONAL :: keep_attr(:)
1615CHARACTER(len=*),INTENT(in),OPTIONAL :: delete_attr(:)
1616logical,intent(in),optional :: preserve
1617logical,intent(in),optional :: purgeana
1618
1619integer :: inddativar,inddatiattrinv,inddatiattrcli,inddatiattrtem,inddatiattrspa,inddativarattr
1620type(qcattrvars) :: attrvars
1621
1622INTEGER(kind=int_b),pointer :: invbb(:,:,:,:,:),clibb(:,:,:,:,:),tembb(:,:,:,:,:),spabb(:,:,:,:,:)
1623INTEGER,pointer :: invbi(:,:,:,:,:),clibi(:,:,:,:,:),tembi(:,:,:,:,:),spabi(:,:,:,:,:)
1624REAL,pointer :: invbr(:,:,:,:,:),clibr(:,:,:,:,:),tembr(:,:,:,:,:),spabr(:,:,:,:,:)
1625DOUBLE PRECISION,pointer :: invbd(:,:,:,:,:),clibd(:,:,:,:,:),tembd(:,:,:,:,:),spabd(:,:,:,:,:)
1626CHARACTER(len=vol7d_cdatalen),pointer :: invbc(:,:,:,:,:),clibc(:,:,:,:,:),tembc(:,:,:,:,:),spabc(:,:,:,:,:)
1627
1628call l4f_log(l4f_info,'starting peeling')
1629
1630call init(attrvars)
1631
1632! generate code per i vari tipi di dati di v7d
1633! tramite un template e il preprocessore
1634
1635
1636#undef VOL7D_POLY_SUBTYPE
1637#undef VOL7D_POLY_SUBTYPES
1638#define VOL7D_POLY_SUBTYPE REAL
1639#define VOL7D_POLY_SUBTYPES r
1640
1641#undef VOL7D_POLY_TYPE
1642#undef VOL7D_POLY_TYPES
1643#define VOL7D_POLY_TYPE REAL
1644#define VOL7D_POLY_TYPES r
1645#include "modqc_peeling_include.F90"
1646#undef VOL7D_POLY_TYPE
1647#undef VOL7D_POLY_TYPES
1648#define VOL7D_POLY_TYPE DOUBLE PRECISION
1649#define VOL7D_POLY_TYPES d
1650#include "modqc_peeling_include.F90"
1651#undef VOL7D_POLY_TYPE
1652#undef VOL7D_POLY_TYPES
1653#define VOL7D_POLY_TYPE INTEGER
1654#define VOL7D_POLY_TYPES i
1655#include "modqc_peeling_include.F90"
1656#undef VOL7D_POLY_TYPE
1657#undef VOL7D_POLY_TYPES
1658#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1659#define VOL7D_POLY_TYPES b
1660#include "modqc_peeling_include.F90"
1661#undef VOL7D_POLY_TYPE
1662#undef VOL7D_POLY_TYPES
1663#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1664#define VOL7D_POLY_TYPES c
1665#include "modqc_peeling_include.F90"
1666
1667
1668#undef VOL7D_POLY_SUBTYPE
1669#undef VOL7D_POLY_SUBTYPES
1670#define VOL7D_POLY_SUBTYPE DOUBLE PRECISION
1671#define VOL7D_POLY_SUBTYPES d
1672
1673#undef VOL7D_POLY_TYPE
1674#undef VOL7D_POLY_TYPES
1675#define VOL7D_POLY_TYPE REAL
1676#define VOL7D_POLY_TYPES r
1677#include "modqc_peeling_include.F90"
1678#undef VOL7D_POLY_TYPE
1679#undef VOL7D_POLY_TYPES
1680#define VOL7D_POLY_TYPE DOUBLE PRECISION
1681#define VOL7D_POLY_TYPES d
1682#include "modqc_peeling_include.F90"
1683#undef VOL7D_POLY_TYPE
1684#undef VOL7D_POLY_TYPES
1685#define VOL7D_POLY_TYPE INTEGER
1686#define VOL7D_POLY_TYPES i
1687#include "modqc_peeling_include.F90"
1688#undef VOL7D_POLY_TYPE
1689#undef VOL7D_POLY_TYPES
1690#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1691#define VOL7D_POLY_TYPES b
1692#include "modqc_peeling_include.F90"
1693#undef VOL7D_POLY_TYPE
1694#undef VOL7D_POLY_TYPES
1695#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1696#define VOL7D_POLY_TYPES c
1697#include "modqc_peeling_include.F90"
1698
1699
1700#undef VOL7D_POLY_SUBTYPE
1701#undef VOL7D_POLY_SUBTYPES
1702#define VOL7D_POLY_SUBTYPE INTEGER
1703#define VOL7D_POLY_SUBTYPES i
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#undef VOL7D_POLY_SUBTYPE
1733#undef VOL7D_POLY_SUBTYPES
1734#define VOL7D_POLY_SUBTYPE INTEGER(kind=int_b)
1735#define VOL7D_POLY_SUBTYPES b
1736
1737#undef VOL7D_POLY_TYPE
1738#undef VOL7D_POLY_TYPES
1739#define VOL7D_POLY_TYPE REAL
1740#define VOL7D_POLY_TYPES r
1741#include "modqc_peeling_include.F90"
1742#undef VOL7D_POLY_TYPE
1743#undef VOL7D_POLY_TYPES
1744#define VOL7D_POLY_TYPE DOUBLE PRECISION
1745#define VOL7D_POLY_TYPES d
1746#include "modqc_peeling_include.F90"
1747#undef VOL7D_POLY_TYPE
1748#undef VOL7D_POLY_TYPES
1749#define VOL7D_POLY_TYPE INTEGER
1750#define VOL7D_POLY_TYPES i
1751#include "modqc_peeling_include.F90"
1752#undef VOL7D_POLY_TYPE
1753#undef VOL7D_POLY_TYPES
1754#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1755#define VOL7D_POLY_TYPES b
1756#include "modqc_peeling_include.F90"
1757#undef VOL7D_POLY_TYPE
1758#undef VOL7D_POLY_TYPES
1759#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1760#define VOL7D_POLY_TYPES c
1761#include "modqc_peeling_include.F90"
1762
1763
1764
1765#undef VOL7D_POLY_SUBTYPE
1766#undef VOL7D_POLY_SUBTYPES
1767#define VOL7D_POLY_SUBTYPE CHARACTER(len=vol7d_cdatalen)
1768#define VOL7D_POLY_SUBTYPES c
1769
1770#undef VOL7D_POLY_TYPE
1771#undef VOL7D_POLY_TYPES
1772#define VOL7D_POLY_TYPE REAL
1773#define VOL7D_POLY_TYPES r
1774#include "modqc_peeling_include.F90"
1775#undef VOL7D_POLY_TYPE
1776#undef VOL7D_POLY_TYPES
1777#define VOL7D_POLY_TYPE DOUBLE PRECISION
1778#define VOL7D_POLY_TYPES d
1779#include "modqc_peeling_include.F90"
1780#undef VOL7D_POLY_TYPE
1781#undef VOL7D_POLY_TYPES
1782#define VOL7D_POLY_TYPE INTEGER
1783#define VOL7D_POLY_TYPES i
1784#include "modqc_peeling_include.F90"
1785#undef VOL7D_POLY_TYPE
1786#undef VOL7D_POLY_TYPES
1787#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1788#define VOL7D_POLY_TYPES b
1789#include "modqc_peeling_include.F90"
1790#undef VOL7D_POLY_TYPE
1791#undef VOL7D_POLY_TYPES
1792#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1793#define VOL7D_POLY_TYPES c
1794#include "modqc_peeling_include.F90"
1795
1796
1797
1798IF (.NOT.PRESENT(keep_attr) .AND. .NOT.PRESENT(delete_attr) .and. .not. optio_log(preserve)) THEN ! destroy all attributes
1799 IF (ASSOCIATED(this%voldatiattrr)) DEALLOCATE(this%voldatiattrr)
1800 IF (ASSOCIATED(this%voldatiattrd)) DEALLOCATE(this%voldatiattrd)
1801 IF (ASSOCIATED(this%voldatiattri)) DEALLOCATE(this%voldatiattri)
1802 IF (ASSOCIATED(this%voldatiattrb)) DEALLOCATE(this%voldatiattrb)
1803 IF (ASSOCIATED(this%voldatiattrc)) DEALLOCATE(this%voldatiattrc)
1804
1805 CALL delete(this%datiattr)
1806 CALL delete(this%dativarattr)
1807END IF
1808
1809IF (PRESENT(keep_attr)) THEN ! set to missing non requested attributes and reform
1810
1811 if (optio_log(preserve)) call l4f_log(l4f_warn,"preserve parameter ignored: keep_attr passed")
1812 CALL keep_var(this%datiattr%r)
1813 CALL keep_var(this%datiattr%d)
1814 CALL keep_var(this%datiattr%i)
1815 CALL keep_var(this%datiattr%b)
1816 CALL keep_var(this%datiattr%c)
1817 CALL qc_reform(this,data_id, miss=.true., purgeana=purgeana)
1818
1819ELSE IF (PRESENT(delete_attr)) THEN ! set to missing requested attributes and reform
1820
1821 if (optio_log(preserve)) call l4f_log(l4f_warn,"preserve parameter ignored: delete_attr passed")
1822 CALL delete_var(this%datiattr%r)
1823 CALL delete_var(this%datiattr%d)
1824 CALL delete_var(this%datiattr%i)
1825 CALL delete_var(this%datiattr%b)
1826 CALL delete_var(this%datiattr%c)
1827 CALL qc_reform(this,data_id, miss=.true., purgeana=purgeana)
1828
1829ELSE IF (PRESENT(purgeana)) THEN
1830
1831 CALL qc_reform(this,data_id, purgeana=purgeana)
1832
1833ENDIF
1834
1835
1836CONTAINS
1837
1838
1840subroutine qc_reform(this,data_id,miss, purgeana)
1841TYPE(vol7d),INTENT(INOUT) :: this
1842integer,INTENT(inout),pointer,OPTIONAL :: data_id(:,:,:,:,:)
1843logical,intent(in),optional :: miss
1844logical,intent(in),optional :: purgeana
1845
1846integer,pointer :: data_idtmp(:,:,:,:,:)
1847logical,allocatable :: llana(:)
1848integer,allocatable :: anaind(:)
1849integer :: i,j,nana
1850
1851if (optio_log(purgeana)) then
1852 allocate(llana(size(this%ana)))
1853 llana =.false.
1854 do i =1,size(this%ana)
1855 if (associated(this%voldatii)) llana(i)= llana(i) .or. any(c_e(this%voldatii(i,:,:,:,:,:)))
1856 if (associated(this%voldatir)) llana(i)= llana(i) .or. any(c_e(this%voldatir(i,:,:,:,:,:)))
1857 if (associated(this%voldatid)) llana(i)= llana(i) .or. any(c_e(this%voldatid(i,:,:,:,:,:)))
1858 if (associated(this%voldatib)) llana(i)= llana(i) .or. any(c_e(this%voldatib(i,:,:,:,:,:)))
1859 if (associated(this%voldatic)) llana(i)= llana(i) .or. any(c_e(this%voldatic(i,:,:,:,:,:)))
1860
1861#ifdef DEBUG
1862 if (.not. llana(i)) call l4f_log(l4f_debug,"remove station"//t2c(i))
1863#endif
1864
1865 end do
1866
1867 nana=count(llana)
1868
1869
1870 allocate(anaind(nana))
1871
1872 j=0
1873 do i=1,size(this%ana)
1874 if (llana(i)) then
1875 j=j+1
1876 anaind(j)=i
1877 end if
1878 end do
1879
1880
1881 if(present(data_id)) then
1882 allocate(data_idtmp(nana,size(data_id,2),size(data_id,3),size(data_id,4),size(data_id,5)))
1883 data_idtmp=data_id(anaind,:,:,:,:)
1884 if (associated(data_id))deallocate(data_id)
1885 data_id=>data_idtmp
1886 end if
1887
1888 call vol7d_reform(this,miss=miss,lana=llana)
1889
1890 deallocate(llana,anaind)
1891
1892else
1893
1894 call vol7d_reform(this,miss=miss)
1895
1896end if
1897
1898end subroutine qc_reform
1899
1900
1901SUBROUTINE keep_var(var)
1902TYPE(vol7d_var),intent(inout),POINTER :: var(:)
1903
1904INTEGER :: i
1905
1906IF (ASSOCIATED(var)) THEN
1907 if (size(var) == 0) then
1908 var%btable = vol7d_var_miss%btable
1909 else
1910 DO i = 1, SIZE(var)
1911 IF (all(var(i)%btable /= keep_attr(:))) THEN ! n.b. ALL((//)) = .TRUE.
1912 var(i)%btable = vol7d_var_miss%btable
1913 ENDIF
1914 ENDDO
1915 end if
1916ENDIF
1917
1918END SUBROUTINE keep_var
1919
1920SUBROUTINE delete_var(var)
1921TYPE(vol7d_var),intent(inout),POINTER :: var(:)
1922
1923INTEGER :: i
1924
1925IF (ASSOCIATED(var)) THEN
1926 if (size(var) == 0) then
1927 var%btable = vol7d_var_miss%btable
1928 else
1929 DO i = 1, SIZE(var)
1930 IF (any(var(i)%btable == delete_attr(:))) THEN ! n.b. ANY((//)) = .FALSE.
1931 var(i) = vol7d_var_miss
1932 ENDIF
1933 ENDDO
1934 end if
1935ENDIF
1936
1937END SUBROUTINE delete_var
1938
1939END SUBROUTINE vol7d_peeling
1940
1941
1942end module modqc
Variables user in Quality Control.
Definition: modqc.F90:392
Test di dato invalidato.
Definition: modqc.F90:417
Remove data under a defined grade of confidence.
Definition: modqc.F90:397
Check data validity based on single confidence.
Definition: modqc.F90:407
Check data validity based on gross error check.
Definition: modqc.F90:412
Definition of constants to be used for declaring variables of a desired type.
Definition: kinds.F90:251
Definitions of constants and functions for working with missing values.
Utilities and defines for quality control.
Definition: modqc.F90:363
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:374

Generated with Doxygen.