libsim Versione 7.1.11

◆ vdgec()

elemental logical function vdgec ( character(len=vol7d_cdatalen), intent(in)  flag)
private

Data gross error check.

Parametri
[in]flagconfidenza

Definizione alla linea 1256 del file modqc.F90.

1257! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
1258! authors:
1259! Davide Cesari <dcesari@arpa.emr.it>
1260! Paolo Patruno <ppatruno@arpa.emr.it>
1261
1262! This program is free software; you can redistribute it and/or
1263! modify it under the terms of the GNU General Public License as
1264! published by the Free Software Foundation; either version 2 of
1265! the License, or (at your option) any later version.
1266
1267! This program is distributed in the hope that it will be useful,
1268! but WITHOUT ANY WARRANTY; without even the implied warranty of
1269! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
1270! GNU General Public License for more details.
1271
1272! You should have received a copy of the GNU General Public License
1273! along with this program. If not, see <http://www.gnu.org/licenses/>.
1274#include "config.h"
1275
1278
1425module modqc
1426use kinds
1429use vol7d_class
1430
1431
1432implicit none
1433
1434
1436type :: qcpartype
1437 integer (kind=int_b):: att
1438 integer (kind=int_b):: gross_error ! special valuer for "*B33192" when gross error check failed
1439 integer (kind=int_b):: invalidated ! special valuer for "*B33196" when manual invalidation happen
1440end type qcpartype
1441
1443type(qcpartype) :: qcpar=qcpartype(10_int_b,0_int_b,1_int_b)
1444
1445integer, parameter :: nqcattrvars=4
1446CHARACTER(len=10),parameter :: qcattrvarsbtables(nqcattrvars)=(/"*B33196","*B33192","*B33193","*B33194"/)
1447
1448type :: qcattrvars
1449 TYPE(vol7d_var) :: vars(nqcattrvars)
1450 CHARACTER(len=10) :: btables(nqcattrvars)
1451end type qcattrvars
1452
1454interface init
1455 module procedure init_qcattrvars
1456end interface
1457
1459interface peeled
1460 module procedure peeledrb, peeleddb, peeledbb, peeledib, peeledcb &
1461 ,peeledri, peeleddi, peeledbi, peeledii, peeledci &
1462 ,peeledrr, peeleddr, peeledbr, peeledir, peeledcr &
1463 ,peeledrd, peeleddd, peeledbd, peeledid, peeledcd &
1464 ,peeledrc, peeleddc, peeledbc, peeledic, peeledcc
1465end interface
1466
1467
1469interface vd
1470 module procedure vdi,vdb,vdr,vdd,vdc
1471end interface
1472
1474interface vdge
1475 module procedure vdgei,vdgeb,vdger,vdged,vdgec
1476end interface
1477
1479interface invalidated
1480 module procedure invalidatedi,invalidatedb,invalidatedr,invalidatedd,invalidatedc
1481end interface
1482
1483private
1484
1485public vd, vdge, init, qcattrvars_new, invalidated, peeled, vol7d_peeling
1486public qcattrvars, nqcattrvars, qcattrvarsbtables
1487public qcpar, qcpartype, qcsummaryflagb ! ,qcsummaryflagi
1488
1489contains
1490
1491
1492! peeled routines
1493#undef VOL7D_POLY_SUBTYPE
1494#undef VOL7D_POLY_SUBTYPES
1495#undef VOL7D_POLY_ISC
1496#define VOL7D_POLY_SUBTYPE REAL
1497#define VOL7D_POLY_SUBTYPES r
1498
1499#undef VOL7D_POLY_TYPE
1500#undef VOL7D_POLY_TYPES
1501#undef VOL7D_POLY_ISC
1502#undef VOL7D_POLY_TYPES_SUBTYPES
1503#define VOL7D_POLY_TYPE REAL
1504#define VOL7D_POLY_TYPES r
1505#define VOL7D_POLY_TYPES_SUBTYPES rr
1506#include "modqc_peeled_include.F90"
1507#include "modqc_peel_util_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 dr
1514#include "modqc_peeled_include.F90"
1515#include "modqc_peel_util_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 ir
1522#include "modqc_peeled_include.F90"
1523#include "modqc_peel_util_include.F90"
1524#undef VOL7D_POLY_TYPE
1525#undef VOL7D_POLY_TYPES
1526#undef VOL7D_POLY_TYPES_SUBTYPES
1527#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1528#define VOL7D_POLY_TYPES b
1529#define VOL7D_POLY_TYPES_SUBTYPES br
1530#include "modqc_peeled_include.F90"
1531#include "modqc_peel_util_include.F90"
1532#undef VOL7D_POLY_TYPE
1533#undef VOL7D_POLY_TYPES
1534#undef VOL7D_POLY_TYPES_SUBTYPES
1535#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1536#define VOL7D_POLY_TYPES c
1537#define VOL7D_POLY_ISC = 1
1538#define VOL7D_POLY_TYPES_SUBTYPES cr
1539#include "modqc_peeled_include.F90"
1540#include "modqc_peel_util_include.F90"
1541
1542
1543#undef VOL7D_POLY_SUBTYPE
1544#undef VOL7D_POLY_SUBTYPES
1545#undef VOL7D_POLY_ISC
1546#define VOL7D_POLY_SUBTYPE DOUBLE PRECISION
1547#define VOL7D_POLY_SUBTYPES d
1548
1549#undef VOL7D_POLY_TYPE
1550#undef VOL7D_POLY_TYPES
1551#undef VOL7D_POLY_TYPES_SUBTYPES
1552#define VOL7D_POLY_TYPE REAL
1553#define VOL7D_POLY_TYPES r
1554#define VOL7D_POLY_TYPES_SUBTYPES rd
1555#include "modqc_peeled_include.F90"
1556#undef VOL7D_POLY_TYPE
1557#undef VOL7D_POLY_TYPES
1558#undef VOL7D_POLY_TYPES_SUBTYPES
1559#define VOL7D_POLY_TYPE DOUBLE PRECISION
1560#define VOL7D_POLY_TYPES d
1561#define VOL7D_POLY_TYPES_SUBTYPES dd
1562#include "modqc_peeled_include.F90"
1563#undef VOL7D_POLY_TYPE
1564#undef VOL7D_POLY_TYPES
1565#undef VOL7D_POLY_TYPES_SUBTYPES
1566#define VOL7D_POLY_TYPE INTEGER
1567#define VOL7D_POLY_TYPES i
1568#define VOL7D_POLY_TYPES_SUBTYPES id
1569#include "modqc_peeled_include.F90"
1570#undef VOL7D_POLY_TYPE
1571#undef VOL7D_POLY_TYPES
1572#undef VOL7D_POLY_TYPES_SUBTYPES
1573#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1574#define VOL7D_POLY_TYPES b
1575#define VOL7D_POLY_TYPES_SUBTYPES bd
1576#include "modqc_peeled_include.F90"
1577#undef VOL7D_POLY_TYPE
1578#undef VOL7D_POLY_TYPES
1579#undef VOL7D_POLY_TYPES_SUBTYPES
1580#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1581#define VOL7D_POLY_TYPES c
1582#define VOL7D_POLY_TYPES_SUBTYPES cd
1583#include "modqc_peeled_include.F90"
1584
1585
1586#undef VOL7D_POLY_SUBTYPE
1587#undef VOL7D_POLY_SUBTYPES
1588#undef VOL7D_POLY_ISC
1589#define VOL7D_POLY_SUBTYPE INTEGER
1590#define VOL7D_POLY_SUBTYPES i
1591
1592#undef VOL7D_POLY_TYPE
1593#undef VOL7D_POLY_TYPES
1594#undef VOL7D_POLY_TYPES_SUBTYPES
1595#define VOL7D_POLY_TYPE REAL
1596#define VOL7D_POLY_TYPES r
1597#define VOL7D_POLY_TYPES_SUBTYPES ri
1598#include "modqc_peeled_include.F90"
1599#undef VOL7D_POLY_TYPE
1600#undef VOL7D_POLY_TYPES
1601#undef VOL7D_POLY_TYPES_SUBTYPES
1602#define VOL7D_POLY_TYPE DOUBLE PRECISION
1603#define VOL7D_POLY_TYPES d
1604#define VOL7D_POLY_TYPES_SUBTYPES di
1605#include "modqc_peeled_include.F90"
1606#undef VOL7D_POLY_TYPE
1607#undef VOL7D_POLY_TYPES
1608#undef VOL7D_POLY_TYPES_SUBTYPES
1609#define VOL7D_POLY_TYPE INTEGER
1610#define VOL7D_POLY_TYPES i
1611#define VOL7D_POLY_TYPES_SUBTYPES ii
1612#include "modqc_peeled_include.F90"
1613#undef VOL7D_POLY_TYPE
1614#undef VOL7D_POLY_TYPES
1615#undef VOL7D_POLY_TYPES_SUBTYPES
1616#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1617#define VOL7D_POLY_TYPES b
1618#define VOL7D_POLY_TYPES_SUBTYPES bi
1619#include "modqc_peeled_include.F90"
1620#undef VOL7D_POLY_TYPE
1621#undef VOL7D_POLY_TYPES
1622#undef VOL7D_POLY_TYPES_SUBTYPES
1623#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1624#define VOL7D_POLY_TYPES c
1625#define VOL7D_POLY_ISC = 1
1626#define VOL7D_POLY_TYPES_SUBTYPES ci
1627#include "modqc_peeled_include.F90"
1628
1629
1630#undef VOL7D_POLY_SUBTYPE
1631#undef VOL7D_POLY_SUBTYPES
1632#undef VOL7D_POLY_ISC
1633#define VOL7D_POLY_SUBTYPE INTEGER(kind=int_b)
1634#define VOL7D_POLY_SUBTYPES b
1635
1636#undef VOL7D_POLY_TYPE
1637#undef VOL7D_POLY_TYPES
1638#undef VOL7D_POLY_TYPES_SUBTYPES
1639#define VOL7D_POLY_TYPE REAL
1640#define VOL7D_POLY_TYPES r
1641#define VOL7D_POLY_TYPES_SUBTYPES rb
1642#include "modqc_peeled_include.F90"
1643#undef VOL7D_POLY_TYPE
1644#undef VOL7D_POLY_TYPES
1645#undef VOL7D_POLY_TYPES_SUBTYPES
1646#define VOL7D_POLY_TYPE DOUBLE PRECISION
1647#define VOL7D_POLY_TYPES d
1648#define VOL7D_POLY_TYPES_SUBTYPES db
1649#include "modqc_peeled_include.F90"
1650#undef VOL7D_POLY_TYPE
1651#undef VOL7D_POLY_TYPES
1652#undef VOL7D_POLY_TYPES_SUBTYPES
1653#define VOL7D_POLY_TYPE INTEGER
1654#define VOL7D_POLY_TYPES i
1655#define VOL7D_POLY_TYPES_SUBTYPES ib
1656#include "modqc_peeled_include.F90"
1657#undef VOL7D_POLY_TYPE
1658#undef VOL7D_POLY_TYPES
1659#undef VOL7D_POLY_TYPES_SUBTYPES
1660#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1661#define VOL7D_POLY_TYPES b
1662#define VOL7D_POLY_TYPES_SUBTYPES bb
1663#include "modqc_peeled_include.F90"
1664#undef VOL7D_POLY_TYPE
1665#undef VOL7D_POLY_TYPES
1666#undef VOL7D_POLY_TYPES_SUBTYPES
1667#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1668#define VOL7D_POLY_TYPES c
1669#define VOL7D_POLY_ISC = 1
1670#define VOL7D_POLY_TYPES_SUBTYPES cb
1671#include "modqc_peeled_include.F90"
1672
1673
1674#undef VOL7D_POLY_SUBTYPE
1675#undef VOL7D_POLY_SUBTYPES
1676#undef VOL7D_POLY_ISC
1677#define VOL7D_POLY_SUBTYPE CHARACTER(len=vol7d_cdatalen)
1678#define VOL7D_POLY_SUBTYPES c
1679
1680#undef VOL7D_POLY_TYPE
1681#undef VOL7D_POLY_TYPES
1682#undef VOL7D_POLY_TYPES_SUBTYPES
1683#define VOL7D_POLY_TYPE REAL
1684#define VOL7D_POLY_TYPES r
1685#define VOL7D_POLY_TYPES_SUBTYPES rc
1686#include "modqc_peeled_include.F90"
1687#undef VOL7D_POLY_TYPE
1688#undef VOL7D_POLY_TYPES
1689#undef VOL7D_POLY_TYPES_SUBTYPES
1690#define VOL7D_POLY_TYPE DOUBLE PRECISION
1691#define VOL7D_POLY_TYPES d
1692#define VOL7D_POLY_TYPES_SUBTYPES dc
1693#include "modqc_peeled_include.F90"
1694#undef VOL7D_POLY_TYPE
1695#undef VOL7D_POLY_TYPES
1696#undef VOL7D_POLY_TYPES_SUBTYPES
1697#define VOL7D_POLY_TYPE INTEGER
1698#define VOL7D_POLY_TYPES i
1699#define VOL7D_POLY_TYPES_SUBTYPES ic
1700#include "modqc_peeled_include.F90"
1701#undef VOL7D_POLY_TYPE
1702#undef VOL7D_POLY_TYPES
1703#undef VOL7D_POLY_TYPES_SUBTYPES
1704#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1705#define VOL7D_POLY_TYPES b
1706#define VOL7D_POLY_TYPES_SUBTYPES bc
1707#include "modqc_peeled_include.F90"
1708#undef VOL7D_POLY_TYPE
1709#undef VOL7D_POLY_TYPES
1710#undef VOL7D_POLY_TYPES_SUBTYPES
1711#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1712#define VOL7D_POLY_TYPES c
1713#define VOL7D_POLY_ISC = 1
1714#define VOL7D_POLY_TYPES_SUBTYPES cc
1715#include "modqc_peeled_include.F90"
1716
1717
1718subroutine init_qcattrvars(this)
1719
1720type(qcattrvars),intent(inout) :: this
1721integer :: i
1722
1723this%btables(:) =qcattrvarsbtables
1724do i =1, nqcattrvars
1725 call init(this%vars(i),this%btables(i))
1726end do
1727
1728end subroutine init_qcattrvars
1729
1730
1731type(qcattrvars) function qcattrvars_new()
1732
1733call init(qcattrvars_new)
1734
1735end function qcattrvars_new
1736
1737
1745SUBROUTINE vol7d_peeling(this, data_id, keep_attr, delete_attr, preserve, purgeana)
1746TYPE(vol7d),INTENT(INOUT) :: this
1747integer,INTENT(inout),pointer,OPTIONAL :: data_id(:,:,:,:,:)
1748CHARACTER(len=*),INTENT(in),OPTIONAL :: keep_attr(:)
1749CHARACTER(len=*),INTENT(in),OPTIONAL :: delete_attr(:)
1750logical,intent(in),optional :: preserve
1751logical,intent(in),optional :: purgeana
1752
1753integer :: inddativar,inddatiattrinv,inddatiattrcli,inddatiattrtem,inddatiattrspa,inddativarattr
1754type(qcattrvars) :: attrvars
1755
1756INTEGER(kind=int_b),pointer :: invbb(:,:,:,:,:),clibb(:,:,:,:,:),tembb(:,:,:,:,:),spabb(:,:,:,:,:)
1757INTEGER,pointer :: invbi(:,:,:,:,:),clibi(:,:,:,:,:),tembi(:,:,:,:,:),spabi(:,:,:,:,:)
1758REAL,pointer :: invbr(:,:,:,:,:),clibr(:,:,:,:,:),tembr(:,:,:,:,:),spabr(:,:,:,:,:)
1759DOUBLE PRECISION,pointer :: invbd(:,:,:,:,:),clibd(:,:,:,:,:),tembd(:,:,:,:,:),spabd(:,:,:,:,:)
1760CHARACTER(len=vol7d_cdatalen),pointer :: invbc(:,:,:,:,:),clibc(:,:,:,:,:),tembc(:,:,:,:,:),spabc(:,:,:,:,:)
1761
1762call l4f_log(l4f_info,'starting peeling')
1763
1764call init(attrvars)
1765
1766! generate code per i vari tipi di dati di v7d
1767! tramite un template e il preprocessore
1768
1769
1770#undef VOL7D_POLY_SUBTYPE
1771#undef VOL7D_POLY_SUBTYPES
1772#define VOL7D_POLY_SUBTYPE REAL
1773#define VOL7D_POLY_SUBTYPES r
1774
1775#undef VOL7D_POLY_TYPE
1776#undef VOL7D_POLY_TYPES
1777#define VOL7D_POLY_TYPE REAL
1778#define VOL7D_POLY_TYPES r
1779#include "modqc_peeling_include.F90"
1780#undef VOL7D_POLY_TYPE
1781#undef VOL7D_POLY_TYPES
1782#define VOL7D_POLY_TYPE DOUBLE PRECISION
1783#define VOL7D_POLY_TYPES d
1784#include "modqc_peeling_include.F90"
1785#undef VOL7D_POLY_TYPE
1786#undef VOL7D_POLY_TYPES
1787#define VOL7D_POLY_TYPE INTEGER
1788#define VOL7D_POLY_TYPES i
1789#include "modqc_peeling_include.F90"
1790#undef VOL7D_POLY_TYPE
1791#undef VOL7D_POLY_TYPES
1792#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1793#define VOL7D_POLY_TYPES b
1794#include "modqc_peeling_include.F90"
1795#undef VOL7D_POLY_TYPE
1796#undef VOL7D_POLY_TYPES
1797#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1798#define VOL7D_POLY_TYPES c
1799#include "modqc_peeling_include.F90"
1800
1801
1802#undef VOL7D_POLY_SUBTYPE
1803#undef VOL7D_POLY_SUBTYPES
1804#define VOL7D_POLY_SUBTYPE DOUBLE PRECISION
1805#define VOL7D_POLY_SUBTYPES d
1806
1807#undef VOL7D_POLY_TYPE
1808#undef VOL7D_POLY_TYPES
1809#define VOL7D_POLY_TYPE REAL
1810#define VOL7D_POLY_TYPES r
1811#include "modqc_peeling_include.F90"
1812#undef VOL7D_POLY_TYPE
1813#undef VOL7D_POLY_TYPES
1814#define VOL7D_POLY_TYPE DOUBLE PRECISION
1815#define VOL7D_POLY_TYPES d
1816#include "modqc_peeling_include.F90"
1817#undef VOL7D_POLY_TYPE
1818#undef VOL7D_POLY_TYPES
1819#define VOL7D_POLY_TYPE INTEGER
1820#define VOL7D_POLY_TYPES i
1821#include "modqc_peeling_include.F90"
1822#undef VOL7D_POLY_TYPE
1823#undef VOL7D_POLY_TYPES
1824#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1825#define VOL7D_POLY_TYPES b
1826#include "modqc_peeling_include.F90"
1827#undef VOL7D_POLY_TYPE
1828#undef VOL7D_POLY_TYPES
1829#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1830#define VOL7D_POLY_TYPES c
1831#include "modqc_peeling_include.F90"
1832
1833
1834#undef VOL7D_POLY_SUBTYPE
1835#undef VOL7D_POLY_SUBTYPES
1836#define VOL7D_POLY_SUBTYPE INTEGER
1837#define VOL7D_POLY_SUBTYPES i
1838
1839#undef VOL7D_POLY_TYPE
1840#undef VOL7D_POLY_TYPES
1841#define VOL7D_POLY_TYPE REAL
1842#define VOL7D_POLY_TYPES r
1843#include "modqc_peeling_include.F90"
1844#undef VOL7D_POLY_TYPE
1845#undef VOL7D_POLY_TYPES
1846#define VOL7D_POLY_TYPE DOUBLE PRECISION
1847#define VOL7D_POLY_TYPES d
1848#include "modqc_peeling_include.F90"
1849#undef VOL7D_POLY_TYPE
1850#undef VOL7D_POLY_TYPES
1851#define VOL7D_POLY_TYPE INTEGER
1852#define VOL7D_POLY_TYPES i
1853#include "modqc_peeling_include.F90"
1854#undef VOL7D_POLY_TYPE
1855#undef VOL7D_POLY_TYPES
1856#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1857#define VOL7D_POLY_TYPES b
1858#include "modqc_peeling_include.F90"
1859#undef VOL7D_POLY_TYPE
1860#undef VOL7D_POLY_TYPES
1861#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1862#define VOL7D_POLY_TYPES c
1863#include "modqc_peeling_include.F90"
1864
1865
1866#undef VOL7D_POLY_SUBTYPE
1867#undef VOL7D_POLY_SUBTYPES
1868#define VOL7D_POLY_SUBTYPE INTEGER(kind=int_b)
1869#define VOL7D_POLY_SUBTYPES b
1870
1871#undef VOL7D_POLY_TYPE
1872#undef VOL7D_POLY_TYPES
1873#define VOL7D_POLY_TYPE REAL
1874#define VOL7D_POLY_TYPES r
1875#include "modqc_peeling_include.F90"
1876#undef VOL7D_POLY_TYPE
1877#undef VOL7D_POLY_TYPES
1878#define VOL7D_POLY_TYPE DOUBLE PRECISION
1879#define VOL7D_POLY_TYPES d
1880#include "modqc_peeling_include.F90"
1881#undef VOL7D_POLY_TYPE
1882#undef VOL7D_POLY_TYPES
1883#define VOL7D_POLY_TYPE INTEGER
1884#define VOL7D_POLY_TYPES i
1885#include "modqc_peeling_include.F90"
1886#undef VOL7D_POLY_TYPE
1887#undef VOL7D_POLY_TYPES
1888#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1889#define VOL7D_POLY_TYPES b
1890#include "modqc_peeling_include.F90"
1891#undef VOL7D_POLY_TYPE
1892#undef VOL7D_POLY_TYPES
1893#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1894#define VOL7D_POLY_TYPES c
1895#include "modqc_peeling_include.F90"
1896
1897
1898
1899#undef VOL7D_POLY_SUBTYPE
1900#undef VOL7D_POLY_SUBTYPES
1901#define VOL7D_POLY_SUBTYPE CHARACTER(len=vol7d_cdatalen)
1902#define VOL7D_POLY_SUBTYPES c
1903
1904#undef VOL7D_POLY_TYPE
1905#undef VOL7D_POLY_TYPES
1906#define VOL7D_POLY_TYPE REAL
1907#define VOL7D_POLY_TYPES r
1908#include "modqc_peeling_include.F90"
1909#undef VOL7D_POLY_TYPE
1910#undef VOL7D_POLY_TYPES
1911#define VOL7D_POLY_TYPE DOUBLE PRECISION
1912#define VOL7D_POLY_TYPES d
1913#include "modqc_peeling_include.F90"
1914#undef VOL7D_POLY_TYPE
1915#undef VOL7D_POLY_TYPES
1916#define VOL7D_POLY_TYPE INTEGER
1917#define VOL7D_POLY_TYPES i
1918#include "modqc_peeling_include.F90"
1919#undef VOL7D_POLY_TYPE
1920#undef VOL7D_POLY_TYPES
1921#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1922#define VOL7D_POLY_TYPES b
1923#include "modqc_peeling_include.F90"
1924#undef VOL7D_POLY_TYPE
1925#undef VOL7D_POLY_TYPES
1926#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1927#define VOL7D_POLY_TYPES c
1928#include "modqc_peeling_include.F90"
1929
1930
1931
1932IF (.NOT.PRESENT(keep_attr) .AND. .NOT.PRESENT(delete_attr) .and. .not. optio_log(preserve)) THEN ! destroy all attributes
1933 IF (ASSOCIATED(this%voldatiattrr)) DEALLOCATE(this%voldatiattrr)
1934 IF (ASSOCIATED(this%voldatiattrd)) DEALLOCATE(this%voldatiattrd)
1935 IF (ASSOCIATED(this%voldatiattri)) DEALLOCATE(this%voldatiattri)
1936 IF (ASSOCIATED(this%voldatiattrb)) DEALLOCATE(this%voldatiattrb)
1937 IF (ASSOCIATED(this%voldatiattrc)) DEALLOCATE(this%voldatiattrc)
1938
1939 CALL delete(this%datiattr)
1940 CALL delete(this%dativarattr)
1941END IF
1942
1943IF (PRESENT(keep_attr)) THEN ! set to missing non requested attributes and reform
1944
1945 if (optio_log(preserve)) call l4f_log(l4f_warn,"preserve parameter ignored: keep_attr passed")
1946 CALL keep_var(this%datiattr%r)
1947 CALL keep_var(this%datiattr%d)
1948 CALL keep_var(this%datiattr%i)
1949 CALL keep_var(this%datiattr%b)
1950 CALL keep_var(this%datiattr%c)
1951 CALL qc_reform(this,data_id, miss=.true., purgeana=purgeana)
1952
1953ELSE IF (PRESENT(delete_attr)) THEN ! set to missing requested attributes and reform
1954
1955 if (optio_log(preserve)) call l4f_log(l4f_warn,"preserve parameter ignored: delete_attr passed")
1956 CALL delete_var(this%datiattr%r)
1957 CALL delete_var(this%datiattr%d)
1958 CALL delete_var(this%datiattr%i)
1959 CALL delete_var(this%datiattr%b)
1960 CALL delete_var(this%datiattr%c)
1961 CALL qc_reform(this,data_id, miss=.true., purgeana=purgeana)
1962
1963ELSE IF (PRESENT(purgeana)) THEN
1964
1965 CALL qc_reform(this,data_id, purgeana=purgeana)
1966
1967ENDIF
1968
1969
1970CONTAINS
1971
1972
1974subroutine qc_reform(this,data_id,miss, purgeana)
1975TYPE(vol7d),INTENT(INOUT) :: this
1976integer,INTENT(inout),pointer,OPTIONAL :: data_id(:,:,:,:,:)
1977logical,intent(in),optional :: miss
1978logical,intent(in),optional :: purgeana
1979
1980integer,pointer :: data_idtmp(:,:,:,:,:)
1981logical,allocatable :: llana(:)
1982integer,allocatable :: anaind(:)
1983integer :: i,j,nana
1984
1985if (optio_log(purgeana)) then
1986 allocate(llana(size(this%ana)))
1987 llana =.false.
1988 do i =1,size(this%ana)
1989 if (associated(this%voldatii)) llana(i)= llana(i) .or. any(c_e(this%voldatii(i,:,:,:,:,:)))
1990 if (associated(this%voldatir)) llana(i)= llana(i) .or. any(c_e(this%voldatir(i,:,:,:,:,:)))
1991 if (associated(this%voldatid)) llana(i)= llana(i) .or. any(c_e(this%voldatid(i,:,:,:,:,:)))
1992 if (associated(this%voldatib)) llana(i)= llana(i) .or. any(c_e(this%voldatib(i,:,:,:,:,:)))
1993 if (associated(this%voldatic)) llana(i)= llana(i) .or. any(c_e(this%voldatic(i,:,:,:,:,:)))
1994
1995#ifdef DEBUG
1996 if (.not. llana(i)) call l4f_log(l4f_debug,"remove station"//t2c(i))
1997#endif
1998
1999 end do
2000
2001 nana=count(llana)
2002
2003
2004 allocate(anaind(nana))
2005
2006 j=0
2007 do i=1,size(this%ana)
2008 if (llana(i)) then
2009 j=j+1
2010 anaind(j)=i
2011 end if
2012 end do
2013
2014
2015 if(present(data_id)) then
2016 allocate(data_idtmp(nana,size(data_id,2),size(data_id,3),size(data_id,4),size(data_id,5)))
2017 data_idtmp=data_id(anaind,:,:,:,:)
2018 if (associated(data_id))deallocate(data_id)
2019 data_id=>data_idtmp
2020 end if
2021
2022 call vol7d_reform(this,miss=miss,lana=llana)
2023
2024 deallocate(llana,anaind)
2025
2026else
2027
2028 call vol7d_reform(this,miss=miss)
2029
2030end if
2031
2032end subroutine qc_reform
2033
2034
2035SUBROUTINE keep_var(var)
2036TYPE(vol7d_var),intent(inout),POINTER :: var(:)
2037
2038INTEGER :: i
2039
2040IF (ASSOCIATED(var)) THEN
2041 if (size(var) == 0) then
2042 var%btable = vol7d_var_miss%btable
2043 else
2044 DO i = 1, SIZE(var)
2045 IF (all(var(i)%btable /= keep_attr(:))) THEN ! n.b. ALL((//)) = .TRUE.
2046 var(i)%btable = vol7d_var_miss%btable
2047 ENDIF
2048 ENDDO
2049 end if
2050ENDIF
2051
2052END SUBROUTINE keep_var
2053
2054SUBROUTINE delete_var(var)
2055TYPE(vol7d_var),intent(inout),POINTER :: var(:)
2056
2057INTEGER :: i
2058
2059IF (ASSOCIATED(var)) THEN
2060 if (size(var) == 0) then
2061 var%btable = vol7d_var_miss%btable
2062 else
2063 DO i = 1, SIZE(var)
2064 IF (any(var(i)%btable == delete_attr(:))) THEN ! n.b. ANY((//)) = .FALSE.
2065 var(i) = vol7d_var_miss
2066 ENDIF
2067 ENDDO
2068 end if
2069ENDIF
2070
2071END SUBROUTINE delete_var
2072
2073END SUBROUTINE vol7d_peeling
2074
2075
2076end 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.