libsim Versione 7.2.1

◆ vdgec()

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

Data gross error check.

Parametri
[in]flagconfidenza

Definizione alla linea 1250 del file modqc.F90.

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