libsim Versione 7.1.11

◆ invalidatedc()

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

Data invalidated check.

Parametri
[in]flagattributo di invalidazione del dato

Definizione alla linea 1270 del file modqc.F90.

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