libsim Versione 7.2.4

◆ vdgeb()

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

Data gross error check.

Parametri
[in]flagconfidenza

Definizione alla linea 1071 del file modqc.F90.

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