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