libsim Versione 7.1.11

◆ vdgeb()

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

Data gross error check.

Parametri
[in]flagconfidenza

Definizione alla linea 1077 del file modqc.F90.

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